How to find the differential equation ? ode solver

So I am very new to matlab ; now wish to use ode solver. I am a little stuck with a thing ; If I have an ouput plot as my objective to plot ; and I have the solutions of the differential equations - is there a way to know he differential equation which can be put in the ode solver in order to reach the objective plot ?

Answer 1

[Torsten]

If the output plot function is f(t) and f(t0)=f0, then the corresponding ODE reads

dy/dt = f'(t) with y(t0) = f0.

This ODE should reproduce y = f.

Best wishes

Torsten.

[STP]

% Constants
beta=5;
alfa = 2.*beta/(beta+1);
tau1=4.4;
tau2=5;
Tc=1.24;
gamma=alfa.*(2-exp(-tau1));
% Time frames
t1=0:0.01:tau1;
t11 = tau1:0.01:8;
t2 = tau1:0.01:tau2;
t22 = tau2:0.01:8;
t3 = tau2:0.01:8;

% Part A
EeA = @(t) -alfa .*exp(-t./Tc) + alfa;
EeA1 = EeA(t1);
EeA2 = EeA(t11);
plot(t1,EeA1,'-b',t11,EeA2,'--c','lineWidth',2)
hold on
Ee1 = EeA(tau1);
scatter(tau1,Ee1,'ro','lineWidth',2,'MarkerFaceColor','r')

% Part B
EeB = @(t) gamma.*exp(-((t-tau1)./Tc)) -alfa;
EeB1 = EeB(t2);
EeB2 = EeB(t22);
plot(t2, EeB1,'-b',t22, EeB2,'--c','lineWidth',2)
Ee2 = EeB(tau2);
scatter(tau2,Ee2,'ro','lineWidth',2,'MarkerFaceColor','r')

% Part C
EeC = @(t) Ee2.*exp(-((t3-tau2)./Tc));
EeC1 = EeC(t3);
plot(t3, EeC1,'-b','lineWidth',2)
hold off
 xlabel('t');
 ylabel('E_e');
 xticklabels({'0', '4.4', '5'})
 xticks([ 0 4.4 5])
 xticklabels({'0', '4.4', '5'})
 xticklabels({'0', 't1', 't2'})

Above mentioned is my normal code with ode solver to plot the curve shown : an the equations are the solutions of the equation I need to put in ode solver in order to get the same curve ; How to go about it in such equations ? Many thanks!

[Torsten]

dEeA = @(t,y) alfa/Tc * exp(-t/Tc);
tspan = t1;
y0 = 0;
[T1 EeA1] = ode45(dEeA,tspan,y0);
tspan = t11;
y0 = EeA1(end);
[T11 EeA2] = ode45(dEeA,tspan,y0);
plot(t1,EeA1,'-b',t11,EeA2,'--c','lineWidth',2)

[STP]

Super, thanks! Although can you tellme how did u come to dEeA; A tool or function? or you solve it manually ? ? I have a few more equations ; and I wish to be able to do it myself :)

[Torsten]

syms alfa t Tc
f = -alfa*exp(-t/Tc) + alfa;
df = diff(f,t)

[STP]

Super; thanks a lot!

I have this solution :

EeB = @(t) gamma.*exp(-((t-tau1)./Tc)) -alfa;
EeB1 = EeB(t2);
EeB2 = EeB(t22);
plot(t2, EeB1,'-b',t22, EeB2,'--c','lineWidth',2)
Ee2 = EeB(tau2);
scatter(tau2,Ee2,'ro','lineWidth',2,'MarkerFaceColor','r')

And I found out the df the way you did; but when I use the ode solver; I get error -

dEeB = @(t,y) (-gamma *exp(-(t - tau1)/Tc)/Tc) ;
tspan2 = t2;
y1 = EeA1(end);
[T2 EeB1] = ode45(dEeB,tspan2,y1);

plot(t2,EeB1,'-b','lineWidth',2)

[Torsten]

Please show the complete code you are using and the error message you get.

[STP]

% Constants
beta = 5;
alfa = 2.*beta./(beta+1);
Tc = 0.20;
gamma = alfa.*(2-exp(-t./Tc));
t1=0:0.01:4.4;
t2=4.4:0.01:5;
t3=5:0.01:12



dEeA = @(t,y) alfa/Tc * exp(-t./Tc);
tspan1 = t1;
y0 = 0;
[T1 EeA1] = ode15s(dEeA,tspan1,y0);

plot(t1,EeA1,'-b','lineWidth',2)



dEeB = @(t,y) -gamma *exp(-(t - tau1)/Tc)/Tc ;
tspan2 = t2;
y1 = EeA1(end);
[T2 EeB1] = ode15s(dEeB,tspan2,y1);

plot(t2,EeB1,'-b','lineWidth',2)

What I wish to replicate is this plot done via this code(without ode solver) :

% Part A
EeA = @(t) -alfa .*exp(-t./Tc) + alfa;
EeA1 = EeA(t1);
EeA2 = EeA(t11);
plot(t1,EeA1,'-b',t11,EeA2,'--c','lineWidth',2)
hold on
Ee1 = EeA(tau1);
scatter(tau1,Ee1,'ro','lineWidth',2,'MarkerFaceColor','r')

% Part B
EeB = @(t) gamma.*exp(-((t-tau1)./Tc)) -alfa;
EeB1 = EeB(t2);
EeB2 = EeB(t22);
plot(t2, EeB1,'-b',t22, EeB2,'--c','lineWidth',2)
Ee2 = EeB(tau2);
scatter(tau2,Ee2,'ro','lineWidth',2,'MarkerFaceColor','r')

% Part C
EeC = @(t) Ee2.*exp(-((t3-tau2)./Tc));
EeC1 = EeC(t3);
plot(t3, EeC1,'-b','lineWidth',2)
hold off

Error message is ::

Error using odearguments (line 113)
Inputs must be floats, namely single or double.

Error in ode15s (line 150)
    odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0,
    options, varargin);
 

[Torsten]

% Constants

beta=5;
alfa = 2.*beta/(beta+1);
tau1=4.4;
tau2=5;
Tc=1.24;
gamma=alfa.*(2-exp(-tau1));
% Time frames
t1=0:0.01:tau1;
t11 = tau1:0.01:8;
t2 = tau1:0.01:tau2;
t22 = tau2:0.01:8;
t3 = tau2:0.01:8;

dEeA = @(t,y) alfa/Tc * exp(-t/Tc);
tspan = t1;
y0 = 0;
[T1 EeA1] = ode45(dEeA,tspan,y0);
tspan = t11;
y0 = EeA1(end);
[T11 EeA2] = ode45(dEeA,tspan,y0);
plot(t1,EeA1,'-b',t11,EeA2,'--c','lineWidth',2)

dEeB = @(t,y)-gamma *exp(-(t - tau1)/Tc)/Tc ;
tspan = t2;
y0 = gamma-alpha;
[T1 EeB1] = ode45(dEeB,tspan,y0);
tspan = t22;
y0 = EeB1(end);
[T11 EeB2] = ode45(dEeB,tspan,y0);
plot(t2, EeB1,'-b',t22, EeB2,'--c','lineWidth',2)
Ee2 = EeB2(end);
scatter(tau2,Ee2,'ro','lineWidth',2,'MarkerFaceColor','r')

dEeC = @(t,y) -Ee2*exp(-(t-tau2)/Tc)/Tc;
tspan = t3;
y0 = Ee2;
[T1 EeC1] = ode45(dEeC,tspan,y0);
plot(t3, EeC1,'-b','lineWidth',2)
hold off
 xlabel('t');
 ylabel('E_e');
 xticklabels({'0', '4.4', '5'})
 xticks([ 0 4.4 5])
 xticklabels({'0', '4.4', '5'})
 xticklabels({'0', 't1', 't2'})

[STP]

Hi, But If I use the above: I get -

Error using alpha Too many output arguments

.

[Torsten]

I did not conform to your misspelling of "alpha" :-)

Use

y0 = gamma-alfa;

instead of

y0 = gamma-alpha;

[STP]

Yes, I realised it aftre posting! Thanks :)

Although the third plot is off ; its not like it should be; doesnt match with the ouput plot. Does it match for you?

dEeC = @(t,y) -Ee2*exp(-(t-tau2)/Tc)/Tc;

the third equation solution is : ( gamma.*exp(-((t-tau1)./Tc)) -alfa).*exp(-((t3-tau2)./Tc));

time span be [4.4 5]

[Torsten]

But you prescribed

EeC = @(t) Ee2.*exp(-((t3-tau2)./Tc));

which does not depend on t.

And now you want to use

EeC = @(t) ( gamma.*exp(-((t-tau1)./Tc)) -alfa).*exp(-((t3-tau2)./Tc));

?

[STP]

Oh it does. Maybe this picture below will make it clear; Ee2 is nothing but the second equation we already plot with the exponential part multiplied to it.

EeA = @(t) -alfa .*exp(-t./Tc) + alfa;

EeB = @(t) gamma.*exp(-((tau2-tau1)./Tc)) -alfa;

EeC = @(t) ( gamma.*exp(-((tau2-tau1)./Tc)) -alfa).*exp(-((t3-tau2)./Tc));

t1= [0 4.4]

t2 = [4.4 5]

t3 = [5 ~]

alfa = 5

tau1 = 4.4

tau2= 5

etc

[Torsten]

According to the equations,

EeC = @(t) ( gamma.*exp(-(tau2-tau1)./Tc) -alfa).*exp(-(t-tau2)./Tc);

which means that you get

% Constants

beta=5;
alfa = 2.*beta/(beta+1);
tau1=4.4;
tau2=5;
Tc=1.24;
gamma=alfa.*(2-exp(-tau1));

% Time frames
t1=0:0.01:tau1;
t11 = tau1:0.01:8;
t2 = tau1:0.01:tau2;
t22 = tau2:0.01:8;
t3 = tau2:0.01:8;

dEeA = @(t,y) alfa/Tc * exp(-t/Tc);
tspan = t1;
y0 = 0;
[T1 EeA1] = ode45(dEeA,tspan,y0);
tspan = t11;
y0 = EeA1(end);
[T11 EeA2] = ode45(dEeA,tspan,y0);
plot(t1,EeA1,'-b',t11,EeA2,'--c','lineWidth',2)

dEeB = @(t,y)-gamma *exp(-(t - tau1)/Tc)/Tc ;
tspan = t2;
y0 = gamma-alfa;
[T1 EeB1] = ode45(dEeB,tspan,y0);
tspan = t22;
y0 = EeB1(end);
[T11 EeB2] = ode45(dEeB,tspan,y0);
plot(t2, EeB1,'-b',t22, EeB2,'--c','lineWidth',2)
Ee2 = EeB1(end);
scatter(tau2,Ee2,'ro','lineWidth',2,'MarkerFaceColor','r')

dEeC = @(t,y) -Ee2*exp(-(t-tau2)/Tc)/Tc;
tspan = t3;
y0 = Ee2;
[T1 EeC1] = ode45(dEeC,tspan,y0);
plot(t3, EeC1,'-b','lineWidth',2)
hold off
 xlabel('t');
 ylabel('E_e');
 xticklabels({'0', '4.4', '5'})
 xticks([ 0 4.4 5])
 xticklabels({'0', '4.4', '5'})
 xticklabels({'0', 't1', 't2'})

[STP]

Thanks a ton! :)

[STP]

Hi, sorry to bother again, I found out the initial conditions; and wanted to check whether they after ode solver and the above back drafting method you showed me leads to the same output; but somehow it doesnt; Can you make out whatsa wrong?

Tc = 2.0, alfa = beta=5; alfa = 2.beta/(beta+1); tau1= 2; tau2= 2.4; gamma=alfa.(2-exp(-tau1));*

Tc.*dEe/dt + Ee = -alfa.*Ek ... 1

initial condions be : for time span1= [0 t1] or eg [0 4.2] Ek = -1

for time span2 = [t1 t2] or [4.2 5] Ek = 1 How does solving eq.1 lead to below

EeA = @(t) -alfa .*exp(-t) + alfa; EeB = @(t) gamma.*exp(-(t-tau1)) - alfa;

My output comes nowhere like it should; its a straight line than the exponential curve it should be

beta = 5;
Ek = -1;
Tc = 2.0;
tau1 = 2;
tau2 = 2.4;
gamma=alfa.*(2-exp(-tau1));
alfa = 2.*beta/(beta+1);
dEedt  = @(t,y) ((-alfa.*Ek) - Ee)./Tc
t1 = [0 4.2]
y0 = 0;
[t1, EeA] = ode15s(dEedt,t1, y0);
plot(t1,EeA)