Borrar filtros
Borrar filtros

How can I solve a second degree DAE ?

4 visualizaciones (últimos 30 días)
Zoé Cord'homme
Zoé Cord'homme el 8 de Jul. de 2022
Comentada: Zoé Cord'homme el 13 de Jul. de 2022
Hi !
For a project, I am currently needing to solve a second degree (or of index 2, I am not too familiar with those) DAE.
The equations are :
d²x/dt² = 1/m * (- f2(x-x1, dx/dt - dx1/dt))
0 = f2(x-x1, dx/dt - dx1/dt) - f1(x1,dx1/dt)
Where m is a scalar and f1 and f2 are functions looking like this :
function force=f1(x-x1, dx/dt - dx1/dt)
%Coefficients
p1 = 1.4347e+07;
p2 = 1.9757e+06;
p3 = 3.1841e+05;
if vrel<=0 %COMPRESSION
force = p1*(x-x1)^2 + p2*(x-x1) + p3;
else %DETENTE
force = p1*(x-x1)^2 + p2*(x-x1) + p3 - 260000;
end
So far, I have coded this :
m_materiel= ...;
M=[1 0 0 0
0 1 0 0
0 0 0 0
0 0 0 1]; % mass matrix
V0 =...; %initial condition
y0=[0 0 V0 V0];
dt=0.01;
tf=2;
tspan=0:dt:tf;
options = odeset('Mass',M,'Vectorized','on');
[t,Y]=ode15s(@(t,Y) f(t,Y,m_materiel),tspan,y0,options);
%-----------------------------------------
function out=f(t,Y,m_materiel)
out =[Y(3,:);
Y(4,:);
f_MI7984(Y(1,:),Y(3,:),data1,data2)-fQS_continue_MI20(Y(2,:)-Y(1,:),Y(4,:)-Y(3,:));
1/m*(fQS_continue_MI20(Y(2,:)-Y(1,:),Y(4,:)-Y(3,:)))];
I have used a similar technique than with classical second order differential equations,
Here "out" is dY where Y = [x1; x; dx1/dt; dx/dt]
When I run the code, I get multiple errors :
Error using vertcat
Dimensions of matrices being concatenated are not consistent.
Error in interp1q (line 31)
[~, j] = sort([x;xxi]);
Error in f (line 2)
out =[Y(3,:);
Error in @(t,Y)f(t,Y,m_materiel)
Error in odenumjac (line 143)
Fdel = feval(F,Fargs_expanded{:});
Error in daeic12 (line 37)
[DfDy,Joptions.fac,nF] = odenumjac(fun, {t0,y,args{:}}, f, Joptions);
Error in ode15s (line 310)
[y,yp,f0,dfdy,nFE,nPD,Jfac] = daeic12(odeFcn,odeArgs,t,ICtype,Mt,y,yp0,f0,...
Error in attelageeqVSmur (line 24)
[t,Y]=ode15s(@(t,Y) f(t,Y,m_materiel),tspan,y0,options);
Can you help me out ? Do you know how to solve this kind of equation ? Thank you :)

Respuesta aceptada

Torsten
Torsten el 8 de Jul. de 2022
Editada: Torsten el 8 de Jul. de 2022
Write your equations as
dx/dt = x2
dx2/dt = 1/m * (- f2(x-x1, x2 - dx1/dt))
0 = f2(x-x1, x2 - dx1/dt) - f1(x1,dx1/dt)
Are you able to solve the third equation for dx1/dt ?
If yes, insert the expression for dx1/dt in the call to f2 in equation 2 and use ode15s, if no, use ode15i.
  6 comentarios
Torsten
Torsten el 13 de Jul. de 2022
Please include the equations you try to solve and the code you use.
Zoé Cord'homme
Zoé Cord'homme el 13 de Jul. de 2022
found the issue thx ! my initial conditions were inconsistent !

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Ordinary Differential Equations en Help Center y File Exchange.

Productos


Versión

R2014a

Community Treasure Hunt

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

Start Hunting!

Translated by