DAE Problem for implicit method. Discrete dynamic system.
1 Ansicht (letzte 30 Tage)
Ältere Kommentare anzeigen
Cesar García Echeverry
am 8 Mai 2021
Bearbeitet: Cesar García Echeverry
am 11 Mai 2021
%Parameters
Nest=100;
alfa = 1.2;
beta=0.75;
nu=alfa-1;
h = 0.9;
f = 0;
%Initial values
X=zeros(Nest,1);
Y=zeros(Nest,1);
X(1)=f;
Y(Nest)=h;
%Initial conditions to the DAE Method
init1=[X;Y];
tspan=[0 30];
rpar=[alfa,beta,nu,f,h,Nest];
M=[0 1 ;0 0];
options=odeset('Refine',5,'MassSingular','yes','Mass',M,'RelTol',1e-4,'AbsTol',1e-6);
[t,var]=ode15s(@diffX,tspan,init1,options,rpar);
function [vard]=diffX(t,var,rpar)
alfa=rpar(1);
beta=rpar(2);
nu=rpar(3);
f=rpar(4);
h=rpar(5);
Nest=rpar(6);
% var=reshape(var,[Nest,2]);
X=var(1:Nest);
Y=var(Nest+1:end);
for i=2:Nest-1
vard(i,1)=beta.*X(i-1)+Y(i+1)-beta.*X(i)-Y(i);
vard(i,2)=(alfa.*X(i))./(1+nu.*X(i-1));%Y(i+1)+nu.*X(i-1).*Y(i)-alfa.*X(i);
end
vard=reshape(vard,[2*(Nest-1),1]);
vard=[0;vard;Y(Nest)];
% vard=vard';
end
2 Kommentare
Akzeptierte Antwort
Walter Roberson
am 8 Mai 2021
X=zeros(Nest,1);
Y=zeros(Nest,1);
init1=[X;Y];
So your initial conditions are 2*Nest x 1 and your function must return 2*Nest x 1.
For a mass matrix M then M*y_prime = f(x, y) and f must return 2*Nest x 1 and y_prime must be that size. In order to have M*y_prime be the same size as f(x, y) then M must be a square matrix with rows (and columns) the same size as the number of rows in y_prime. Which is 2*Nest
But does M have that size? No. M is 2x2 and so can only be used if Nest is 1.
14 Kommentare
Walter Roberson
am 10 Mai 2021
PDE are often dealt with by dividing an area into discrete domains.
One of your recent questions mentions distillation columns. Physically those are continuous columns -- but subdividing into discrete parts to model behaviour through the column would be a common approach. It is also the most common approach to PDE.
Weitere Antworten (0)
Siehe auch
Kategorien
Mehr zu Boundary Conditions finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!