Effacer les filtres
Effacer les filtres

Why do I receive this error msg when solving a simple DAE with odefunction?

1 vue (au cours des 30 derniers jours)
Hainan Wang
Hainan Wang le 9 Juin 2021
Commenté : Hainan Wang le 20 Juin 2021
When I solve this simple DAE with three variables y1(t), y2(t), y3(t),
y1' = y2' + y1;
y2' = y1;
y3 = y1 + y2;
y1(0) = 1.0;
y2(0) = 0.5;
y3(0) = 1.5;
syms y1(t) y2(t) y3(t)
eqs = [diff(y1(t),t) == diff(y2(t),t) + y1(t),...
diff(y2(t),t) == y1(t),...
y3(t) == y1(t) + y2(t)];
vars = [y1(t) y2(t) y3(t)];
initConditions = [1 0.5 1.5];
[M,F] = massMatrixForm(eqs,vars);
M = odeFunction(M,vars);
F = odeFunction(F,vars);
opt = odeset('mass',M);
ode15s(F, [0 1], initConditions, opt)
I receive this error msg:
Need a better guess y0 for consistent initial conditions.
This DAE is simple, if subs y2' = y1 into first equation, then y1' = 2y1. With y1(0) = 1, y2(0) = 0.5, can solve y1 = e^2t and y2 = e^2t/2. And y3(0) = y1(0)+y2(0)=1.5. What's wrong with my guess on initial condition, most likely for y3?

Réponses (1)

Harshitha Kallam
Harshitha Kallam le 15 Juin 2021
Yes, the most likely reason for error is the initial value of y3 that is resulting in inconsistency.
A possible workaround can be to first solve for y1(t), y2(t) and then calculate y3(t) which is a summation of those two. Alternatively, the below modification can be done.
eqs = [diff(y1(t),t) == diff(y2(t),t) + y1(t),...
diff(y2(t),t) == y1(t),...
diff(y3(t),t) == diff(y1(t) + y2(t),t)]; %changing the last eqn
  1 commentaire
Hainan Wang
Hainan Wang le 20 Juin 2021
Thank you for your response!
However, I figure out that I can solve it with below code while keeping the 3rd equation as y3(t)=y1(t)+y2(t).
M = [1 -1 0; 0 1 0; 0 0 0];
optODE = odeset('Mass', M);
initConditions = [1 0.5 1.5]';
ode15s(@(t,y)dyneqn(t,y), [0 1], initConditions, optODE)
function dydt = dyneqn(t,y)
dydt = zeros(3,1);
dydt(1) = y(1);
dydt(2) = y(1);
dydt(3) = y(1) + y(2) - y(3);
end

Connectez-vous pour commenter.

Catégories

En savoir plus sur Programming dans Help Center et File Exchange

Tags

Produits


Version

R2021a

Community Treasure Hunt

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

Start Hunting!

Translated by