DAE Problem for implicit method. Discrete dynamic system.
5 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
Cesar García Echeverry
le 8 Mai 2021
Modifié(e) : Cesar García Echeverry
le 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 commentaires
Réponse acceptée
Walter Roberson
le 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 commentaires
Walter Roberson
le 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.
Plus de réponses (0)
Voir également
Catégories
En savoir plus sur Boundary Conditions dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!