Vectorizing DAE with strong State Dependance Mass matrix
显示 更早的评论
Hello,
Please help me to vectorize my DAE, I want Vectorized='on'
Thank you
Andrii
a=0.4; b=0.1; c=0.9;
A = @(y) c/y(1);
Ms=@(t,y) [A(y) 0 -a/y(2); 0 b 0; 0 0 0];
Ds=@(t,y) [y(1)*y(3)-y(2)+1; y(1)-1; a*y(1)+y(2)+y(3)];
Y0=[1, 1, -2];
opts = odeset('MStateDependence','strong','MassSingular','yes','Mass',Ms,Vectorized='off');
%opts = odeset('MStateDependence','strong','MassSingular','yes','Mass',Ms,Vectorized='on');
[t,Y] = ode15s(Ds,[0 10],Y0,opts);
plot(t,Y,"-x")
采纳的回答
The Vectorized option in odeset means that your DAEs can be writen as 
, apparently your DAEs can not be expressed as separable form.
, apparently your DAEs can not be expressed as separable form.[a,b,c] = deal(.4,.1,.9);
Ms = @(t,y) [c./y(1),0,-a./y(2);0,b,0;zeros(1,3)];
Ds = @(t,y) [y(1).*y(3)-y(2)+1;y(1)-1;a.*y(1)+y(2)+y(3)];
Y0=[1;1;-2];
opts = odeset('MStateDependence','strong','MassSingular','yes','Mass',Ms,Vectorized='off');
[t,Y] = ode15s(Ds,[0 10],Y0,opts);
plot(t,Y)
7 个评论
Chuguang Pan
2025-7-9
编辑：Chuguang Pan
2025-7-9
@Andrii There is an example which use Vectorized option=="on" for solving stiff ODEs, it maybe helpful.
You should supply consistent initial values Y0 for the unknowns. Since your Y0 vector does not satisfy the last algebraic equation in t = 0, you will get almost arbitrary solutions for y(1) - y(3). The reason is the following:
If Y0 does not satisfy the algebraic equation, MATLAB must modify one (or more) of the given initial values in Y0=[1;1;-2]. But which one ?
If MATLAB assumes that y(1) and y(2) are differential variables and y(3) is an algebraic variable, it will start with y(1) = y(2) = 1 and y(3) = -a*1 - 1 = -1.4.
If MATLAB assumes that y(1) and y(3) are differential variables and y(2) is an algebraic variable, it will start with y(1) = 1, y(3) = -2 and y(2) = -a*1 + 2 = 1.6.
If MATLAB assumes that y(2) and y(3) are differential variables and y(1) is an algebraic variable, it will start with y(2) = 1, y(3) = -2 and y(1) = (-1 + 2)/a = 2.5.
In all three cases, you will get different solutions.
Looking at what MATLAB returns for y at time t = 0 shows that the chosen Y0 is even more complicated:
[a,b,c] = deal(.4,.1,.9);
Ms = @(t,y) [c./y(1),0,-a./y(2);0,b,0;zeros(1,3)];
Ds = @(t,y) [y(1).*y(3)-y(2)+1;y(1)-1;a.*y(1)+y(2)+y(3)];
Y0=[1;1;-2];
opts = odeset('MStateDependence','strong','MassSingular','yes','Mass',Ms,Vectorized='off');
[t,Y] = ode15s(Ds,[0 10],Y0,opts);
Y(1,:)
The vectorized option for a system of only three equations really doesn't matter.
Which is not obvios why not.
If MATLAB would accept 'Vectorize','on' in case of a mass matrix, you would have to return a matrix which has vectors as elements, thus a 3d- instead of a 2d - object.
If you want to try to vectorize the code, you could:
- Differentiate the algebraic equation
- Don't define a mass matrix, but just multiply f by M^(-1) in the function defining the ODEs
Then, without vectorization, your code would be like below. I'm not able to vectorize it - I had to use a for-loop in function "fun" to deal with matrix inputs for y.
Note that you get different results compared to the code with algebraic equation. The reason again is that your Y0 vector is not consistent.
[a,b,c] = deal(.4,.1,.9);
Y0=[1;1;-2];
options = odeset('Vectorized','on');
[t,Y] = ode15s(@(t,y)fun(t,y,a,b,c),[0 10],Y0,options);
plot(t,Y)
function dydt = fun(t,y,a,b,c)
dydt = zeros(3,size(y,2));
for i = 1:size(y,2)
M = [c/y(1,i),0,-a/y(2,i);0,b,0;a,1,1];
f = [y(1,i)*y(3,i)-y(2,i)+1;y(1,i)-1;0];
dydt(:,i) = M\f;
end
end
Andrii
2025-7-11
Torsten
2025-7-11
Lets see if vectorization will save me some calculation time.
But the code from above is not vectorized - you are still left with the need to vectorize the loop.
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Numerical Integration and Differential Equations 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
