Résolution 2ème ordre ODE, Matlab- l'accélération de l'équation a besoin de sa propre valeur afin d'inclure un autre terme différent

j'ai cette 2ème ordre ODE à résoudre dans Matlab:Résolution 2ème ordre ODE, Matlab- l'accélération de l'équation a besoin de sa propre valeur afin d'inclure un autre terme différent

(a + f(t))·(dx/dt)·(d²x/dt²) + g(t) + ((h(t) + i(t)·(d²x/dt² > b·(c-x)))·(dx/dt) + j(t))·(dx/dt)² + k(t)·(t > d) = 0

où

a, b, c, d sont connues constantes
f(t), g(t), h(t), i(t), j(t), k(t) sont connus des fonctions dépendantes de la t
x est la position
dx/dt est la vitesse
d²x/dt² est l'accélération

et noter les deux conditions que

i(t) est introduit dans l'équation si (d²x/dt² > b·(c-x))
k(t) est introduit dans l'équation si (t > d)

Ainsi, le problème pourrait être résolu avec une structure similaire à Matlab comme cet exemple:

[T,Y] = ode45(@(t,y) [y(2); 'the expression of the acceleration'], tspan, [x0 v0]);

où

T est le vecteur temporel, Y est le vecteur de position (colonne 1 as y(1)) et la vitesse (colonne 2 comme y(2)).
ode45 est le solveur ODE, mais un autre pourrait être utilisé.
tspan, x0, v0 sont connus.
the expression of the acceleration signifie une expression pour d²x/dt², mais vient ici le problème, car il est dans la condition pour i(t) et « en dehors » en même temps la multiplication (a + f(t))·(dx/dt). Ainsi, l'accélération ne peut pas être écrit en Matlab comme d²x/dt² = something

Quelques questions qui pourraient aider:

une fois l'état (d²x/dt² > b·(c-x)) et/ou (t > d) est satisfaite, le terme respectif i(t) et/ou k(t) sera être introduit jusqu'à la fin de l'heure déterminée au tspan.
pour la condition (d²x/dt² > b·(c-x)), le terme d²x/dt² pourrait être écrit comme la différence des vitesses, comme y(2) - y(2)', si y(2)' est la vitesse de l'instant précédent, divisé par le temps de l'étape définie dans tspan. Mais je ne sais pas comment accéder à la valeur précédente de la vitesse pendant la résolution du ODE

Merci à l'avance!

Source

2017-09-24 Sergio

Vous pouvez essayer un solveur implicite, comme [ 'ode15i'] (https://uk.mathworks.com/help/matlab/ref/ode15i.html) – Steve

@ Steve: Et qu'est-ce que cela va résoudre? – Wrzlprmft

@Wrzlprmft Cela résout le fait qu'il s'agit d'un ODE implicite, où 'x_tt: = d^2x/dt^2' est * effectivement * un argument implicite comme' H (x_tt - b (c * x)) 'où' H (t) 'est la fonction d'escalier Heaviside' H (t) = 1 si t> = 0, 0 si t <0'. Si la solution existe ou a un sens est une autre question à regarder. – Steve

Tout d'abord, vous devriez reduce your problem to a first-order differential equation, en remplaçant dx/dt par une variable dynamique pour la vélocité. C'est quelque chose que vous devez faire de toute façon pour résoudre l'ODE et ainsi vous n'avez pas besoin d'accéder aux valeurs précédentes de la vélocité. Pour ce qui est de la réalisation de vos conditions, il suffit de modifier la fonction que vous transmettez à ode45 pour tenir compte de cela. A cet effet, vous pouvez utiliser d²x/dt² dans la partie droite de votre ODE. Gardez à l'esprit que les solveurs ODE n'aiment pas les discontinuités, vous pouvez donc vouloir lisser l'étape ou simplement redémarrer le solveur avec une fonction différente, une fois que la condition est remplie (crédit Steve).

Source

2017-09-25 08:58:10 Wrzlprmft

Je suis entièrement d'accord avec le premier paragraphe, mais en ce qui concerne la seconde, 'ode45' relis sur que vous pouvez exprimer l'équation en termes de la dérivée d'ordre le plus élevé, ce qui n'est pas possible en raison du' i (t) * (x_tt> b (cx)) * x_t * x_tt' terme – Steve

De toute façon, vous ne pouvez pas évaluer explicitement cette condition puisqu'elle reste activée une fois que vous l'avez activée et vous devez de toute façon compter sur l'absence de paradoxe. (Voir aussi mon édition.) – Wrzlprmft

+1 D'accord, cela me semble plus logique maintenant. J'ai commencé à écrire la partie "redémarrer le solveur avec une fonction différente", donc je vais finir et publier, mais cette réponse semble l'avoir couverte. – Steve

Le deuxième terme conditionnel k(t)*(t>d) devrait être assez simple à implémenter, donc je vais passer à côté de ça.

Je diviser votre équation en deux parties:

F1(t,x,x',x'') := (a+f(t))x'x'' + g(t) + (h(t)x'+j(t))x'' + k(t)(t>d), 
F2(t,x,x',x'') := F1(t,x,x',x'') + i(t)x'x'',

où les premiers dérivés de temps désignent. Comme suggéré dans this other answer

[...] ou tout simplement redémarrer le solveur avec une autre fonction

vous pouvez résoudre le ODE F1 pour t \in [t0, t1] =: tspan. Ensuite, vous trouverez la première fois tstar où x''> b(c-x) et les valeurs x(tstar) et x'(tstar), et de résoudre F2 pour t \in [tstar,t1] avec x(tstar), x'(tstar) comme conditions de départ. Après avoir dit tout cela, la mise en œuvre correcte de cela devrait utiliser events, comme suggéré dans LutzL's comment.

Source

2017-09-25 09:37:18 Steve

Alors, je devrais utiliser quelque chose qui ressemble à ceci:

function [value,isterminal,direction] = ODE_events(t,y,b,c) 

value = d²x/dt² - b*(c-y(1)); % detect (d²x/dt² > b·(c-x)), HOW DO I WRITE d²x/dt² HERE? 
isterminal = 0;     % continue integration 
direction = 0;     % zero can be approached in either direction

puis inclure dans le fichier (.m), où mon ode est ceci:

refine = 4; % I do not get exactly how this number would affect the results 
options = odeset('Events',@ODE_events,'OutputFcn',@odeplot,'OutputSel',1, 'Refine',refine); 

[T,Y] = ode45(@(t,y) [y(2); (1/(a + f(t))*(y(2)))*(- g(t) - ((h(t) + i(t))·(y(2)) - j(t)·(y(2))² - k(t)*(t > d)) ], tspan, [x0 v0], options);

Comment puis-je gérer i(t)? parce que i(t)*(d²x/dt² > b*(c-y(1))))*(y(2)) doit être inclus d'une manière ou d'une autre.

Merci encore

Source

2017-09-30 09:06:23 Sergio

Résolution 2ème ordre ODE, Matlab- l'accélération de l'équation a besoin de sa propre valeur afin d'inclure un autre terme différent

Répondre

Questions connexes