代码如下:前四行都是会用到的参数,odefun中有3个微分方程与1个代数方程
lamdap=975e-9;lamdas=1064e-9;t=1e-3;ap=2100e-27;ep=2200e-27;A=2.826e-7;
R1=0.6;as=5e-27;es=331e-27;alfap=8e-3;alfal=6e-3;tp=0.8e-3;N=2.1e28;R2=1;
h=6.626e-34;c=3e8;vp=c/lamdap;vs=c/lamdas;ts=0.8e-3;Pp0=0.2;
Pssat=h*vs*A/(t*tp*(es+as));Ppsat=h*vp*A/(t*tp*(ep+ap));
odefun=@(z,y)[-((N-y(4))*ap-y(4)*ep)*tp*y(1)-alfap*y(1);
(y(4)*ep-(N-y(4))*ap*y(2)-alfal*y(2));
(y(4)*ep-(N-y(4))*ap*y(3)-alfal*y(3));
((tp*y(1)*ap/(Ppsat*(ap+ep))+((y(2)+y(3))*as)/Pssat*(as+es))/(y(1)/Ppsat+1+(y(2)+y(3))/Pssat))-y(4)/N];%微分与代数方程组
M=[1 0 0 0;0 1 0 0;0 0 1 0;0 0 0 0];%质量矩阵
options=odeset('mass',M);%对于DAE问题,mass属性必须设置
y0=[1;0;0.1;N];%初值
[z,y]=ode15s(odefun,[0:0.01:0.2],y0,options);
figure('numbertitle','off','name','DAE demo—by Matlabsky')
plot(z,y,'*')
legend('y1(z)','y2(z)','y3(z)','y4(z)')
用matlab运行后程序没有错误,但是出现了如下的情况:
警告: 矩阵为奇异值、接近奇异值或缩放错误。结果可能不准确。RCOND = NaN。
> In ode15s at 589
In Untitled4 at 13
警告: 在 t=0.000000e+00 处失败。在时间 t 处,若不将步长降至允许的最小值(7.905050e-323)以下,积分公差要求无法满足。
> In ode15s at 668
In Untitled4 at 13
问题就出在ode15ss那里。方程是肯定有解的,不知道是不是因为很多参数的值很大的原因。