Solving system of ODEs and then get the optimum value for certain parameters.

3 views (last 30 days)
I want to optimize fro R1 and R2 with the objective function (t_f-1.2e-6).^2+(t_t-50e-6).^2 == 0. Hereis my code that takes R1 and R2 to solve two differential eqns and give me t_f and t_t.
function [t_f,t_t] = fun1(R1,R2)
C1 = 25e-9;
C2 = 1200e-12;
V0 = 10e3;
F = @(t,V) [((1/(R1*C1))*(V(2) - V(1)));((V(1)/(R1*C2))-(V(2)/(R1*C2))-(V(2)/(R2*C2)))];
[t1, V1]=ode45(F,[0 300e-6], [V0 0]);
Vp = max(V1(:,2));
for n=1:size((V1(:,2)),1)
if V1(n,2) == Vp
t_f = t1(n);
elseif ((V1(n,2)-(Vp/2))*1e-3)^2 < 1e-5
t_t =t1(n);
end
end
end
Someone please help me in using optimization tools of matlab like fminsearch or ga to get the optimal values of R1 and R2.

Accepted Answer

Eduard Reitmann
Eduard Reitmann on 6 Aug 2018
I am not sure what type of input arguments the "fun1" needs, but it gives the following error when I try fun1(1,1):
Output argument "t_t" (and maybe others) not assigned during call to
"test>fun1".
Both output arguments (t_f & t_t) need to be assigned in the "fun1" function in order to solve the objective function. I made the following changes to "fun1" to assign both outputs (assume zero if not assigned):
if V1(n,2) == Vp
t_f = t1(n);
t_t = 0;
elseif ((V1(n,2)-(Vp/2))*1e-3)^2 < 1e-5
t_t = t1(n);
t_f = 0;
end
I you are happy with the modification above, this code should work theoretically.
fun = @(x) objfun(x(1),x(2));
x0 = [1;1];
options.Display = 'iter';
x = fminsearch(fun,x0,options);
R1 = x(1)
R2 = x(2)
Just remember to define this objective function (after script or separate file, depending on MATLAB version):
function fval = objfun(R1,R2)
[t_f,t_t] = fun1(R1,R2);
fval = (t_f-1.2e-6).^2 + (t_t-50e-6).^2;
end
I got the following answer (you can play with the "options" parameter to increase the accuracy of the answer):
R1 =
1.2505e+03
R2 =
1.3598e+03
  3 Comments
Eduard Reitmann
Eduard Reitmann on 6 Aug 2018
I would rewrite your fun1 code to:
function [t_f,t_t,t,V] = fun1(R1,R2)
C1 = 25e-9;
C2 = 1200e-12;
V0 = 10e3;
F = @(t,V) [((1/(R1*C1))*(V(2) - V(1)));((V(1)/(R1*C2))-(V(2)/(R1*C2))-(V(2)/(R2*C2)))];
[t, V]=ode45(F,[0 300e-6], [V0 0]);
V = V(:,2);
[Vmax,imax] = max(V);
t_f = t(imax);
ihalf = find(V>=Vmax/2,1,'last');
t_t = t(ihalf);
end
The execution code should then be:
fun = @(x) objfun(x(1),x(2));
x0 = [100;100];
options.Display = 'iter';
options.TolX = 1e-20;
x = fminsearch(fun,x0,options);
R1 = x(1)
R2 = x(2)
[t_f,t_t,t,V] = fun1(R1,R2);
figure
plot(t,V,t_f,max(V),'*',t_t,max(V)/2,'*')
t_f
t_t
With the objective function:
function [fval,t,V] = objfun(R1,R2)
[t_f,t_t,t,V] = fun1(R1,R2);
fval = (t_f-1.2e-6).^2 + (t_t-50e-6).^2;
end
Using initial conditions R1=100 and R2=100, I get:
R1 =
2.6172e+03
R2 =
191.7561
This gives a function value of 3.84394e-36 instead of 1.6731e-14 for R1=199 and R2=2.5e3, i.e. more optimized answer.

Sign in to comment.

More Answers (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!