2nd order nonlinear differential equations as "State Space Representation" for ode45()

17 views (last 30 days)
Hi Community,
in recent days I read a lot about solving differential equations with the numerical solver ode45(..). But after all it is not clear to me how to handle the following problem, which is only schematically described, due to too high complexity of the actual problem (stated below).
Consider a system of two nonlinear differential equations with two unknowns to be solved for. The first one is of second order in both unknowns and the second equation is of first order in both unknows. The unknowns are y1 & y2 .
Eq1:
f(y1(t)'', y2(t)'', y1(t)', y2(t)', y1(t), y2(t))
e.g. a*y1(t)'' + b*y2(t)'' + c*y2(t)' + d*y1(t) = e*y1(t)'
Eq2:
f(y1(t)', y2(t)', y1(t), y2(t)) e.g.: g*y1(t)' - h*y2(t)' - k*y2(t) = 0
To solve these equations with ode45(...) it is necessary to write them as first order diff. eqns ( State Space Model ). Therefore the equations should be solved for y1(t)'' & y2(t)'', but I don't know how to do this WITHOUT differentiating the second equation and solving the system with backsubstitution. I want to avoid this because it makes things way more complicated due to all the nonlinearities.
Anyone an idea?
For those who are interested: This is my electrical circuit I want to solve for xk, xs & xv. Where
vk = xk', vs = xs' & vv = xv'.
Moreover the electro-mechanical transform has a nonlinearity of the form:
Bl(xk) = bl0 + bl1*xk + bl2*xk^2
The resistances Rs & Rv look like:
Rs = rs0 + rs1 * vs + rs2 * vs^2
Rv = rv0 + rv1 * vv + rv2 * vv^2
and the capacitances Cs & Cv:
Cs = cs0 + cs1 * xs + cs2 * xs^2
Cv = cv0 + cv1 * xv + rv2 * xv^2
and the generator Ug(t) is of the form:
Ug(t) = U*cos(2*pi*f*t)

Accepted Answer

Torsten
Torsten on 3 Mar 2015
Use ODE15I instead of ODE45.
Best wishes
Torsten.
  1 Comment
Helge
Helge on 3 Mar 2015
Edited: Helge on 4 Mar 2015
The documentation of ode15i( ) states that it solves differential equations of the form: f(t, y, y') = 0.
In my situation it is like f(t, y, y', y'', y''') = 0. Can I handle it the same way or do I have to do some substitutions first?

Sign in to comment.

More Answers (1)

Torsten
Torsten on 4 Mar 2015
Set y1=y, y2=y', y3=y''. Then solve the system
y1'=y2
y2'=y3
f(t,y1,y2,y3,y3')=0
Best wishes
Torsten.
  1 Comment
Helge
Helge on 5 Mar 2015
Hi Torsten,
I defined my function the following way, where xs & xv are my two unknows and xsd & xvd the first derivative, xsdd & xvdd the second and finally xsddd & xvddd the third. The rest are just fixed parameters.
function out = myode15func(t, y, yp, U_g, R_e, L_e, M_m, bl_0, bl_1,...
bl_2, r_s_0, r_s_1, r_s_2, s_s_0, s_s_1, s_s_2, r_v_0, r_v_1,...
r_v_2, s_v_0, s_v_1, s_v_2, f0)
xs = y(1);
xsd = y(2);
xsdd = y(3);
xv = y(4);
xvd = y(5);
xvdd = y(6);
xsddd = yp(3);
xvddd = yp(6);
out = zeros(6, 1);
out(1) = xsd - yp(1);
out(2) = xsdd - yp(2);
out(3) = bl_0.*(xsd + xvd) - U_g.*cos(2.*pi.*f0.*t) + L_e.*((r_s_0.*xsdd + s_s_0.*xsd + M_m.*(xsddd + xvddd) + 3.*r_s_2.*xsd.^2.*xsdd + 3.*s_s_2.*xs.^2.*xsd + 2.*r_s_1.*xsd.*xsdd + 2.*s_s_1.*xs.*xsd)./(bl_0 + bl_1.*(xs + xv) + bl_2.*(xs + xv).^2) - ((bl_1.*(xsd + xvd) + 2.*bl_2.*(xs + xv).*(xsd + xvd)).*(s_s_2.*xs.^3 + s_s_1.*xs.^2 + s_s_0.*xs + r_s_2.*xsd.^3 + r_s_1.*xsd.^2 + ...;
out(4) = xvd - yp(4);
out(5) = xvdd - yp(5);
out(6) = 1.0./(bl_0+bl_1.*xs+bl_1.*xv+bl_2.*xs.^2+bl_2.*xv.^2+bl_2.*xs.*xv.*2.0).^2.*(R_e.*bl_0.*s_s_0.*xs-R_e.*bl_0.*s_v_0.*xv-L_e.*bl_1.*r_s_0.*xsd.^2-L_e.*bl_1.*r_s_1.*xsd.^3-L_e.*bl_1.*r_s_2.*xsd.^4+L_e.*bl_1.*r_v_0.*xvd.^2+L_e.*bl_1.*r_v_1.*xvd.^3+L_e.*bl_1.*r_v_2.*xvd.^4+R_e.*bl_0.*r_s_1.*xsd.^2+R_e.*bl_0.*r_s_2.*xsd.^3-R_e.*bl_0.*r_v_1.*xvd.^2-...;
end
Solving this with ode15s () seems to work at the first glance.
Thanks a lot, I wasnt aware of the ode15i()-solver and its syntax!

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!