numerical solution of a system of ODE which is not in standard form
2 次查看(过去 30 天)
显示 更早的评论
Hi
I' m wondering can MATLAB solve a system of ODE of these forms numerically?
H11*dydt+H12*dxdt=H13
H21*dydt+H22*dxdt=H23
I can change it to standard form but I have very large numbers which MATLAB has problem with it and I receive this massage:
Output truncated. Text exceeds maximum line length for Command Window display.
. I will upload my output. I want to solve these ODE numerically so I need the coefficient of all of terms in each eqaution. I utilized the collect cammand but it did not work. How cam I control and arrenge these large coeffisient. how can i determine each coeffisient correctly?
I would be gratefull if some one can help me.
thank you in advance
raha
this is some terms of output differential equation
(930441627500546658162730981425505183132667665707339718275447859267760748743794621
3848815439656366790278984387099311194969091496593187779728810864687259359515581251
8191929791546041229963635314657368394923219128924360408416905997843608302300038558
79197879071442125658901144......................................................0
0000000000000000000*Iass^8*IBss^5*Ipss^4 - 191746860148563438549216527115219808387369783
2286756806355856960741923824185970189......
采纳的回答
Star Strider
2020-7-4
4 个评论
raha ahmadi
2020-7-5
编辑:raha ahmadi
2020-7-5
Dear Star Strider
Thank you for your answer.
Its better to explain my code .
I m going to solve numerically this two coupled ODE :
H11*dIa/dIb+H12* dIp/dIB=H13
H21*dIa/dIB+H22*dIp/dIB=H23
My parameters (system states) are Iass , IBss and Ipss. Also in this system dIadIB is equal to dIass/dIBss and dIp/dIB is dIpss/dIBss.
I think MATLAB need explicit form of ODE( I mean : y’=f(x,y,t) ) to solve . so I change it to explicit form
dIadIb=(H12*H23-H22*H13)/(H12*H21-H22*H11);
dIpdIb=(H11*H23-H21*H13)/(H11*H22-H21*H12);
beacause I had to calculate some mathematical eqautins to reach this ODE system, I prefered to start symbolically because I must plot Iass and Ipss in terms of Ibss after finding coefficient define its in vector form and write its function and solve with ode15s:
but I couldn’t extract coefficient.after sending my quaestin I found that I must solve it with ode15i. Is it right? I will search about what you suggested
I don’t know how write this code. (sorry it is very long)
I really appreciate your help
Best regards
Raha
syms Iass IBss Ipss
%%
%constants
e=1.6*1e-19;
nA=3.2;
nB=3.18;
nP=3.19;
T0=300;
mili=1e-3;
c=3*1e8;
sigma=0.75;
R1=24;
R2=7;
R3=9.4;
micro=1e-6;
lA=300*micro;
wA=1.2*micro;
tA=0.15*micro;
vA=lA*wA*tA;
aA=lA*wA;
lP=250*micro;
wP=1.5*micro;
tP=0.25*micro;%tickness
vP=lP*wP*tP;
aP=lP*wP;
%Bragg section
lB=300*micro;
leffB=135*micro;
wB=1.5*micro;
tB=0.25*micro;%tickness
vB=lB*wB*tB;
aB=lB*wB;
RHOv=4.825e3;
V=[vA,vP,vB,vA,vP,vB];
Cv=3.124e2;
Cc=2.46e-7;
Cs=5.32e-5;
Cj=Cv*RHOv.*V;
C=mean(Cj);
RF=0.35;
RLambda=0.35;
Gamma=0.8;
Ba=-5.97*1e-27;
Bb=-5.97*1e-27;
Bp=-5.97*1e-27;
Ca=1.1*1e-4;
CB=1.6*1e-4;
Cp=1.6*1e-4;
rhoN=[9.24e-5,9.24e-5,2.19e-3,2.229e-3];
dna=[0.15,0.15,2.5,0.5];
dna=dna.*micro;
Ea=0.83;
a1a=sigma*Ea;
a2a=sigma^2*sum(rhoN.*dna)/aA;
Rna=sum(rhoN.*dna)/aA;
dnB=[0.25,0.15,2.5,0.5];
dnB=micro.*dnB;
EB=0.95;
a1B=sigma*EB;
a2B=sigma^2*sum(rhoN.*dnB)/aB;
RnB=sum(rhoN.*dnB)/aB;
dnp=[0.25,0.15,2.5,0.5];
dnp=micro.*dnp;
Ep=0.95;
a1p=sigma*Ep;
a2p=sigma^2*sum(rhoN.*dnp)/aP;
Rnp=sum(rhoN.*dnp)/aP;
l0=lA+lB+lP;
rho0=2.229e-3;
d0=100*micro;
w0=400*micro;
rc=(rho0*d0)/(w0*l0);
a1=0.2*1.1*1e8;
a2=1.1*1e-16-(0.2*2.7)*1e-17;
a3=(4*1e-41-2.7*1e-41);
Tw=302;
T0=300;
b1=1.24*0.28*1e8;
b2=(0.28*1e-16-0.48*1e-17);
b3=(1+1.6*1e-2)*(Tw-T0);
p1=b1;
p2=b2;
p3=b3;
ngA=nA+(0.1*nA);
ngA=3.20+0.32;
dgdNa=4.83*1e-21;
eta=1e-4;
Ntr=1.37*1e24;
mili=1e-3;
nm=1e-9;
Ia0=75*mili;
lambda0=1553*nm;
GammaB=2.4342*1e-7;
zetaA=0.325;
zetaB=0.488;
zetaP=0.488;
VgA=c/ngA;
g=zetaA*dgdNa*VgA;
alphaA=2*1e3; %m^-1
alphaP=1.5*1e3;
alphaB=1.5*1e3;
alphaL=5.46*1e-6 ;
dalphadN=1.8*1e-21;
AprimA=dalphadN;
AprimP=dalphadN;
AprimB=dalphadN;
tauBbar=2.6037e-23; %2.603667095500614e-23
tauPbar =2.3057e-23;
v1=sigma/(e*vA);
v2=sigma/(e*vP);
v3=sigma/(e*vB);
tauA=1/(a1+2*a2*Ntr+3*a3*Ntr^2);
zeta0=zetaA*alphaA*lA+zetaB*alphaB*lB+zetaP*alphaP*lP;
lc0=nA*lA+nB*leffB+nP*lP;
gammaP0=c/lc0*(zeta0-0.5*log(Gamma^2*RF*RLambda));
Ja=(sigma*Ia0)/(e*vA);
S0=zeros(1,4);
S0(1)=-g^2*gammaP0; S0(2)=Ja*g^2-(1/tauA)*Ntr*g^2-2*gammaP0*g/tauA+...
eta*a2*Ntr^2*g^2;
S0(3)=Ja*g/tauA-(1/tauA)^2*Ntr*g-gammaP0*(1/tauA)^2+2*a2*g*Ja*eta*Ntr;
S0(4)=+eta*a2*Ja^2;
rut=roots(S0);
II=find(rut>0);
Sss0=rut(II);
Nss0=(Ja+zetaA*g*VgA*Ntr*Sss0)/(1/tauA+zetaA*g*VgA*Sss0);
tauAnew=1/(tauA^-1+g*Sss0);
nAss0=nA+Ba*(Nss0-Ntr);
lc0=nAss0*lA+nB*leffB+nP*lP;
gammaP0=c/lc0*(zeta0-0.5*log(Gamma^2*RF*RLambda));
mNum=2*lc0/(1553*1e-9);
LambdaB0=1553*1e-9/(2*nB);
alphaAss=AprimA*Iass;
alphaBss=AprimB*IBss;
alphaPss=AprimP*Ipss;
m11=-1/(R2*Cc);
m12=-m11;
m13=1/(Cc*R1);
m21=1/(R2*Cs);
m22=-(1/R2+1/R3)*(1/Cs);
m23=0;
m31=6/(R2*Cc);
m32=-m31;
m33=-1/(R1*C)-6/(Cc*R1);
M=[m11,m12,m13;m21,m22,m23;m31,m32,m33];
detM=det(M);
P=a1a*Ia0+a2a*Ia0^2;
Pc=rc*(1-sigma)^2*Ia0^2;
U1=Pc/Cc;
U2=T0/(R3*Cs);
U3=P/C-(6*Pc)/Cc;
m=[-U1,m12,m13;-U2,m22,m23;-U3,m32,m33];
detm=det(m);
A1=(1/detM)*(-(1/Cc)*(1/R2+1/R3)*...
(6/(R1*Cc+1/(R1*C)))*(1/Cs)+(6/(C*Cc*Cs))*(1/R1)*(1/R2+1/R3));
p=a1a*Ia0+a2a*Ia0^2;
iaas=[a2a a1a -p/3];
iaa=roots(iaas);
iaa1=iaa(imag(iaa)==0);
iaa2=find(iaa1>0);
Iaass=iaa(iaa2);
%I_BB*
ibbs=[a2B a1B -p/3];
ibb=roots(ibbs);
ibb1=iaa(imag(ibb)==0);
ibb2=find(ibb1>0);
Ibbss=ibb(ibb2);
%I_PP*
ipps=[a2p a1p -p/3];
ipp=roots(ipps);
ipp1=ipp(imag(ipp)==0);
ipp2=find(ipp1>0);
Ippss=ipp(ipp2);
Itot=Iaass+Ibbss+Ippss;
Iast=Iaass-Iass;
IBst=Ibbss-IBss;
Ipst=Ippss-Ipss;
LAss=lA*(1+alphaL)*4*Rna*A1*a2a*(Iast-Iass)/(a1a+2*a2a*Iast)*Itot;
LBss=lB*(1+alphaL)*4*RnB*A1*a2B*(IBst-IBss)/(a1B+2*a2B*IBst)*Itot;
LPss=lP*(1+alphaL)*4*Rnp*A1*a2p*(Ipst-Ipss)/(a1p+2*a2p*Ipst)*Itot;
zetass=zetaA*alphaAss*LAss+zetaB*alphaBss*LBss+zetaP*alphaPss*LPss;
Tc=detm/detM
omega=A1*rc*(1-sigma)^2*(Itot);
Aa=omega*(4*a2a*(Iast-Iass))/(a1a+2*a2a*Iast);
AB =omega*(4*a2B*(IBst-IBss))/(a1B+2*a2B*IBst);
Ap= omega*(4*a2p*(Ipst-Ipss))/(a1p+2*a2p*Ipst);
nssA=Ba*Iass*tauAnew+Ca*(Aa*Iass+AB*IBss+Ipss);
nssB=Bb*IBss*tauBbar+CB*(Aa*Iass+AB*IBss+Ipss);
nssP=Ba*Ipss*tauPbar+Cp*(Aa*Iass+AB*IBss+Ipss);
LambdaBss=LambdaB0*(1+alphaL)*(Aa*Iass+AB*IBss+Ap*Ipss);
Lcss=nssA*LAss+nssB*LBss+nssP*LPss;
gammaPss=c/(Lcss)*(zetass-0.5*log(Gamma^2*RF*RLambda));
%%
%initial value
Iass0=75*mili;
IBss0=0;
Ipss0=0;
%%
D1=Ba*tauAnew*v1*LAss+Ca*lA*Aa+nA*LAss*alphaL*Aa+CB*LBss*Aa+nB*LBss*alphaL*Aa+...
Cp*LPss*Aa;
D2=LAss*AB*(Ca+nA*alphaL)+Bb*tauBbar*v2*LBss+CB*LBss*AB+...
AB*nB*LBss*alphaL+Cp*LPss*AB;
D3=Ca*LAss*Ap+nA*LAss*alphaL*Ap+CB*LBss*Ap+nB*LBss*alphaL*Ap+...
Bp*LPss*tauPbar*v3+Cp*LPss*Ap;
X=-2/mNum*(LAss*Ca+nA*LAss*alphaL)+...
2*LambdaBss*(CB+alphaL*Bb)-...
2/mNum*(nB*LBss*alphaL+CB)-2/mNum*(Cp*LPss+nP*LPss*alphaL);
Y=zetaA*alphaA*LAss*alphaL+zetaB*alphaB*LBss*alphaL+zetaP*alphaP*LPss*alphaL;
%%
H11= tauAnew*v1*(-2/mNum*LAss*Ba)+Aa*X;
H12= tauPbar*v3*(-2/mNum*LPss*Bp)+Ap*X;
H13= -(tauBbar*v2*(2*LambdaBss*Bb-2/mNum*LBss*Bb)+AB*X);
H21=2*gammaPss*tauAnew*v1-g*Iass*tauAnew*tauAnew*v1+2*g*Ntr*tauAnew*v1+Iass*tauAnew*c/(Lcss^2)...
*D1*(zetass-0.5*log(Gamma^2*RF*RLambda))+c/Lcss*(zetaA*AprimA*tauAnew*v2+Aa*Y);
H22= c*Iass*tauAnew/(Lcss^2)*D3*(zetass-0.5*log(Gamma^2*RF*RLambda))+(c/Lcss)*(zetaP*AprimP...
*tauPbar*v3+Ap*Y);%Lcss= ss value of length
H23= -((Iass*tauAnew*c)/(Lcss^2)*D2*(zetass-0.5*log(Gamma^2*RF*RLambda))...
+c/(Lcss)*(zetaB*AprimB*tauBbar*v2+AB*Y));
%% A system of two ODE
dIadIb=(H12*H23-H22*H13)/(H12*H21-H22*H11);
dIpdIb=(H11*H23-H21*H13)/(H11*H22-H21*H12);
ddd=collect (dIadIb,[Iass IBss Ipss])
Star Strider
2020-7-5
I find your code difficult to follow, so in part I will defer to John’s analysis.
I still do not understand what you are doing, however the magnitudes of your constants are significantly greater than MATLAB’ floating-point limits, defined as realmax (1.797693134862316e+308), so I doubt if you could solve this numerically and get meaningful results, regardless.
I added these lines at the end of your posted code:
ddd = simplify(ddd, 'Steps',1000);
ddd = vpa(ddd, 5)
ddd_fcn = matlabFunction(ddd)
producing a long but tractable vpa result:
ddd =
-(2.0313e+127*(1.0485e+34*Iass - 6.8648e+32)*(1.4e+32*IBss - 1.0907e+31)*(1.8017e+565*IBss^6*Iass^5*Ipss^4 - 3.713e+564*IBss^6*Iass^5*Ipss^3 + 2.7034e+563*IBss^6*Iass^5*Ipss^2 - 7.9594e+561*IBss^6*Iass^5*Ipss + 7.3457e+559*IBss^6*Iass^5 + 3.6843e+597*IBss^6*Iass^4*Ipss^4 - 2.1412e+596*IBss^6*Iass^4*Ipss^3 - 4.8022e+595*IBss^6*Iass^4*Ipss^2 + 4.8364e+594*IBss^6*Iass^4*Ipss - 1.1825e+593*IBss^6*Iass^4 - 9.1201e+596*IBss^6*Iass^3*Ipss^4 + 5.2034e+595*IBss^6*Iass^3*Ipss^3 + 1.1894e+595*IBss^6*Iass^3*Ipss^2 - 1.1841e+594*IBss^6*Iass^3*Ipss + 2.8658e+592*IBss^6*Iass^3 + 8.4375e+595*IBss^6*Iass^2*Ipss^4 - 4.7133e+594*IBss^6*Iass^2*Ipss^3 - 1.101e+594*IBss^6*Iass^2*Ipss^2 + 1.0819e+593*IBss^6*Iass^2*Ipss - 2.5875e+591*IBss^6*Iass^2 - 3.4562e+594*IBss^6*Iass*Ipss^4 + 1.884e+593*IBss^6*Iass*Ipss^3 + 4.5132e+592*IBss^6*Iass*Ipss^2 - 4.369e+591*IBss^6*Iass*Ipss + 1.0304e+590*IBss^6*Iass + 5.2862e+592*IBss^6*Ipss^4 - 2.8e+591*IBss^6*Ipss^3 - 6.9085e+590*IBss^6*Ipss^2 + 6.5727e+589*IBss^6*Ipss - 1.5244e+588*IBss^6 + 3.2995e+565*IBss^5*Iass^6*Ipss^4 - 6.7997e+564*IBss^5*Iass^6*Ipss^3 + 4.9508e+563*IBss^5*Iass^6*Ipss^2 - 1.4576e+562*IBss^5*Iass^6*Ipss + 1.3452e+560*IBss^5*Iass^6 + 3.6034e+565*IBss^5*Iass^5*Ipss^5 + 7.3686e+597*IBss^5*Iass^5*Ipss^4 - 4.2824e+596*IBss^5*Iass^5*Ipss^3 - 9.6045e+595*IBss^5*Iass^5*Ipss^2 + 9.6728e+594*IBss^5*Iass^5*Ipss - 2.365e+593*IBss^5*Iass^5 + 1.4815e+597*IBss^5*Iass^4*Ipss^5 + 1.9238e+600*IBss^5*Iass^4*Ipss^4 - 1.4547e+599*IBss^5*Iass^4*Ipss^3 - 2.5453e+598*IBss^5*Iass^4*Ipss^2 + 2.9382e+597*IBss^5*Iass^4*Ipss - 7.7342e+595*IBss^5*Iass^4 - 3.6673e+596*IBss^5*Iass^3*Ipss^5 - 5.0391e+599*IBss^5*Iass^3*Ipss^4 + 3.7113e+598*IBss^5*Iass^3*Ipss^3 + 6.6741e+597*IBss^5*Iass^3*Ipss^2 - 7.5633e+596*IBss^5*Iass^3*Ipss + 1.9633e+595*IBss^5*Iass^3 + 3.3928e+595*IBss^5*Iass^2*Ipss^5 + 4.9085e+598*IBss^5*Iass^2*Ipss^4 - 3.5168e+597*IBss^5*Iass^2*Ipss^3 - 6.5091e+596*IBss^5*Iass^2*Ipss^2 + 7.2344e+595*IBss^5*Iass^2*Ipss - 1.8492e+594*IBss^5*Iass^2 - 1.3898e+594*IBss^5*Iass*Ipss^5 - 2.1074e+597*IBss^5*Iass*Ipss^4 + 1.4659e+596*IBss^5*Iass*Ipss^3 + 2.7989e+595*IBss^5*Iass*Ipss^2 - 3.0462e+594*IBss^5*Iass*Ipss + 7.6508e+592*IBss^5*Iass + 2.1256e+592*IBss^5*Ipss^5 + 3.3638e+595*IBss^5*Ipss^4 - 2.2645e+594*IBss^5*Ipss^3 - 4.4757e+593*IBss^5*Ipss^2 + 4.7595e+592*IBss^5*Ipss - 1.1712e+591*IBss^5 + 1.4978e+565*IBss^4*Iass^7*Ipss^4 - 3.0867e+564*IBss^4*Iass^7*Ipss^3 + 2.2474e+563*IBss^4*Iass^7*Ipss^2 - 6.6168e+561*IBss^4*Iass^7*Ipss + 6.1066e+559*IBss^4*Iass^7 + 3.2995e+565*IBss^4*Iass^6*Ipss^5 + 3.6843e+597*IBss^4*Iass^6*Ipss^4 - 2.1412e+596*IBss^4*Iass^6*Ipss^3 - 4.8022e+595*IBss^4*Iass^6*Ipss^2 + 4.8364e+594*IBss^4*Iass^6*Ipss - 1.1825e+593*IBss^4*Iass^6 + 1.8017e+565*IBss^4*Iass^5*Ipss^6 + 1.4815e+597*IBss^4*Iass^5*Ipss^5 + 1.923e+600*IBss^4*Iass^5*Ipss^4 - 1.453e+599*IBss^4*Iass^5*Ipss^3 - 2.5466e+598*IBss^4*Iass^5*Ipss^2 + 2.9387e+597*IBss^4*Iass^5*Ipss - 7.7348e+595*IBss^4*Iass^5 - 2.2028e+597*IBss^4*Iass^4*Ipss^6 - 1.1526e+600*IBss^4*Iass^4*Ipss^5 - 1.0992e+600*IBss^4*Iass^4*Ipss^4 + 1.6892e+599*IBss^4*Iass^4*Ipss^3 + 3.4654e+596*IBss^4*Iass^4*Ipss^2 - 8.2247e+596*IBss^4*Iass^4*Ipss + 2.6311e+595*IBss^4*Iass^4 + 5.4528e+596*IBss^4*Iass^3*Ipss^6 + 3.0168e+599*IBss^4*Iass^3*Ipss^5 + 2.1096e+599*IBss^4*Iass^3*Ipss^4 - 3.8071e+598*IBss^4*Iass^3*Ipss^3 + 8.9002e+596*IBss^4*Iass^3*Ipss^2 + 9.797e+595*IBss^4*Iass^3*Ipss - 3.7412e+594*IBss^4*Iass^3 - 5.0447e+595*IBss^4*Iass^2*Ipss^6 - 2.9378e+598*IBss^4*Iass^2*Ipss^5 - 1.8457e+598*IBss^4*Iass^2*Ipss^4 + 3.5122e+597*IBss^4*Iass^2*Ipss^3 - 1.1063e+596*IBss^4*Iass^2*Ipss^2 - 6.3108e+594*IBss^4*Iass^2*Ipss + 2.7242e+593*IBss^4*Iass^2 + 2.0664e+594*IBss^4*Iass*Ipss^6 + 1.2612e+597*IBss^4*Iass*Ipss^5 + 7.6012e+596*IBss^4*Iass*Ipss^4 - 1.465e+596*IBss^4*Iass*Ipss^3 + 5.0078e+594*IBss^4*Iass*Ipss^2 + 2.1785e+593*IBss^4*Iass*Ipss - 9.9377e+591*IBss^4*Iass - 3.1606e+592*IBss^4*Ipss^6 - 2.013e+595*IBss^4*Ipss^5 - 1.1968e+595*IBss^4*Ipss^4 + 2.292e+594*IBss^4*Ipss^3 - 7.9158e+592*IBss^4*Ipss^2 - 3.1297e+591*IBss^4*Ipss + 1.4245e+590*IBss^4 - 4.3976e+564*IBss^3*Iass^7*Ipss^4 + 9.0127e+563*IBss^3*Iass^7*Ipss^3 - 6.521e+562*IBss^3*Iass^7*Ipss^2 + 1.9072e+561*IBss^3*Iass^7*Ipss - 1.7529e+559*IBss^3*Iass^7 - 9.6875e+564*IBss^3*Iass^6*Ipss^5 - 1.6977e+597*IBss^3*Iass^6*Ipss^4 + 2.2314e+596*IBss^3*Iass^6*Ipss^3 - 1.7171e+594*IBss^3*Iass^6*Ipss^2 - 7.1691e+593*IBss^3*Iass^6*Ipss + 2.282e+592*IBss^3*Iass^6 - 5.2899e+564*IBss^3*Iass^5*Ipss^6 - 6.8266e+596*IBss^3*Iass^5*Ipss^5 - 8.9071e+599*IBss^3*Iass^5*Ipss^4 + 1.3231e+599*IBss^3*Iass^5*Ipss^3 - 1.812e+597*IBss^3*Iass^5*Ipss^2 - 4.3063e+596*IBss^3*Iass^5*Ipss + 1.5054e+595*IBss^3*Iass^5 + 1.015e+597*IBss^3*Iass^4*Ipss^6 + 5.3302e+599*IBss^3*Iass^4*Ipss^5 + 1.9783e+599*IBss^3*Iass^4*Ipss^4 - 4.3147e+598*IBss^3*Iass^4*Ipss^3 + 1.5303e+597*IBss^3*Iass^4*Ipss^2 + 7.2037e+595*IBss^3*Iass^4*Ipss - 3.3966e+594*IBss^3*Iass^4 - 2.5059e+596*IBss^3*Iass^3*Ipss^6 - 1.3913e+599*IBss^3*Iass^3*Ipss^5 - 1.6394e+598*IBss^3*Iass^3*Ipss^4 + 5.9998e+597*IBss^3*Iass^3*Ipss^3 - 3.2204e+596*IBss^3*Iass^3*Ipss^2 - 2.1576e+594*IBss^3*Iass^3*Ipss + 2.9989e+593*IBss^3*Iass^3 + 2.3114e+595*IBss^3*Iass^2*Ipss^6 + 1.3507e+598*IBss^3*Iass^2*Ipss^5 + 6.4777e+596*IBss^3*Iass^2*Ipss^4 - 4.3794e+596*IBss^3*Iass^2*Ipss^3 + 2.8648e+595*IBss^3*Iass^2*Ipss^2 - 2.1318e+593*IBss^3*Iass^2*Ipss - 1.3588e+592*IBss^3*Iass^2 - 9.4357e+593*IBss^3*Iass*Ipss^6 - 5.7777e+596*IBss^3*Iass*Ipss^5 - 1.3643e+595*IBss^3*Iass*Ipss^4 + 1.6382e+595*IBss^3*Iass*Ipss^3 - 1.1602e+594*IBss^3*Iass*Ipss^2 + 1.468e+592*IBss^3*Iass*Ipss + 3.4077e+590*IBss^3*Iass + 1.4375e+592*IBss^3*Ipss^6 + 9.1832e+594*IBss^3*Ipss^5 + 1.5417e+593*IBss^3*Ipss^4 - 2.455e+593*IBss^3*Ipss^3 + 1.7602e+592*IBss^3*Ipss^2 - 2.4732e+590*IBss^3*Ipss - 3.9273e+588*IBss^3 + 4.8236e+563*IBss^2*Iass^7*Ipss^4 - 9.8239e+562*IBss^2*Iass^7*Ipss^3 + 7.0565e+561*IBss^2*Iass^7*Ipss^2 - 2.0479e+560*IBss^2*Iass^7*Ipss + 1.873e+558*IBss^2*Iass^7 + 1.0626e+564*IBss^2*Iass^6*Ipss^5 + 2.4934e+596*IBss^2*Iass^6*Ipss^4 - 4.0675e+595*IBss^2*Iass^6*Ipss^3 + 1.7677e+594*IBss^2*Iass^6*Ipss^2 + 9.0941e+591*IBss^2*Iass^6*Ipss - 1.3333e+591*IBss^2*Iass^6 + 5.8024e+563*IBss^2*Iass^5*Ipss^6 + 1.0026e+596*IBss^2*Iass^5*Ipss^5 + 1.3174e+599*IBss^2*Iass^5*Ipss^4 - 2.3746e+598*IBss^2*Iass^5*Ipss^3 + 1.1432e+597*IBss^2*Iass^5*Ipss^2 + 3.7586e+594*IBss^2*Iass^5*Ipss - 8.8645e+593*IBss^2*Iass^5 - 1.4908e+596*IBss^2*Iass^4*Ipss^6 - 7.8801e+598*IBss^2*Iass^4*Ipss^5 - 1.7197e+598*IBss^2*Iass^4*Ipss^4 + 4.9782e+597*IBss^2*Iass^4*Ipss^3 - 2.8263e+596*IBss^2*Iass^4*Ipss^2 + 6.8425e+593*IBss^2*Iass^4*Ipss + 1.7569e+593*IBss^2*Iass^4 + 3.6632e+595*IBss^2*Iass^3*Ipss^6 + 2.0471e+598*IBss^2*Iass^3*Ipss^5 - 6.7652e+596*IBss^2*Iass^3*Ipss^4 - 3.6825e+596*IBss^2*Iass^3*Ipss^3 + 2.8874e+595*IBss^2*Iass^3*Ipss^2 - 3.0544e+593*IBss^2*Iass^3*Ipss - 1.2063e+592*IBss^2*Iass^3 - 3.3609e+594*IBss^2*Iass^2*Ipss^6 - 1.9765e+597*IBss^2*Iass^2*Ipss^5 + 2.0049e+596*IBss^2*Iass^2*Ipss^4 + 1.1162e+595*IBss^2*Iass^2*Ipss^3 - 1.597e+594*IBss^2*Iass^2*Ipss^2 + 3.123e+592*IBss^2*Iass^2*Ipss + 3.2599e+590*IBss^2*Iass^2 + 1.3637e+593*IBss^2*Iass*Ipss^6 + 8.4009e+595*IBss^2*Iass*Ipss^5 - 1.0446e+595*IBss^2*Iass*Ipss^4 - 1.1557e+593*IBss^2*Iass*Ipss^3 + 4.92e+592*IBss^2*Iass*Ipss^2 - 1.2796e+591*IBss^2*Iass*Ipss - 2.5064e+588*IBss^2*Iass - 2.063e+591*IBss^2*Ipss^6 - 1.3252e+594*IBss^2*Ipss^5 + 1.7188e+593*IBss^2*Ipss^4 + 1.6792e+590*IBss^2*Ipss^3 - 6.6153e+590*IBss^2*Ipss^2 + 1.8539e+589*IBss^2*Ipss - 9.4607e+585*IBss^2 - 2.3415e+562*IBss*Iass^7*Ipss^4 + 4.7344e+561*IBss*Iass^7*Ipss^3 - 3.3722e+560*IBss*Iass^7*Ipss^2 + 9.6969e+558*IBss*Iass^7*Ipss - 8.8166e+556*IBss*Iass^7 - 5.158e+562*IBss*Iass^6*Ipss^5 - 1.4907e+595*IBss*Iass^6*Ipss^4 + 2.6797e+594*IBss*Iass^6*Ipss^3 - 1.5334e+593*IBss*Iass^6*Ipss^2 + 2.4892e+591*IBss*Iass^6*Ipss + 1.5882e+589*IBss*Iass^6 - 2.8166e+562*IBss*Iass^5*Ipss^6 - 5.9941e+594*IBss*Iass^5*Ipss^5 - 7.9326e+597*IBss*Iass^5*Ipss^4 + 1.5618e+597*IBss*Iass^5*Ipss^3 - 9.6562e+595*IBss*Iass^5*Ipss^2 + 1.6773e+594*IBss*Iass^5*Ipss + 1.0672e+592*IBss*Iass^5 + 8.9125e+594*IBss*Iass^4*Ipss^6 + 4.7445e+597*IBss*Iass^4*Ipss^5 + 7.7321e+596*IBss*Iass^4*Ipss^4 - 2.7382e+596*IBss*Iass^4*Ipss^3 + 1.9209e+595*IBss*Iass^4*Ipss^2 - 3.5488e+593*IBss*Iass^4*Ipss - 1.9443e+591*IBss*Iass^4 - 2.1762e+594*IBss*Iass^3*Ipss^6 - 1.2248e+597*IBss*Iass^3*Ipss^5 + 1.0527e+596*IBss*Iass^3*Ipss^4 + 1.0974e+595*IBss*Iass^3*Ipss^3 - 1.2873e+594*IBss*Iass^3*Ipss^2 + 2.8092e+592*IBss*Iass^3*Ipss + 1.134e+590*IBss*Iass^3 + 1.9822e+593*IBss*Iass^2*Ipss^6 + 1.1738e+596*IBss*Iass^2*Ipss^5 - 1.7939e+595*IBss*Iass^2*Ipss^4 + 4.7853e+593*IBss*Iass^2*Ipss^3 + 2.9742e+592*IBss*Iass^2*Ipss^2 - 1.0697e+591*IBss*Iass^2*Ipss - 1.871e+588*IBss*Iass^2 - 7.9756e+591*IBss*Iass*Ipss^6 - 4.9455e+594*IBss*Iass*Ipss^5 + 8.6232e+593*IBss*Iass*Ipss^4 - 4.1144e+592*IBss*Iass*Ipss^3 + 4.9156e+589*IBss*Iass*Ipss^2 + 2.1635e+589*IBss*Iass*Ipss - 2.4916e+586*IBss*Iass + 1.1947e+590*IBss*Ipss^6 + 7.7198e+592*IBss*Ipss^5 - 1.3768e+592*IBss*Ipss^4 + 7.1481e+590*IBss*Ipss^3 - 6.1926e+588*IBss*Ipss^2 - 2.1249e+587*IBss*Ipss + 7.1254e+584*IBss + 4.2413e+560*Iass^7*Ipss^4 - 8.5046e+559*Iass^7*Ipss^3 + 5.9975e+558*Iass^7*Ipss^2 - 1.7057e+557*Iass^7*Ipss + 1.5397e+555*Iass^7 + 9.3433e+560*Iass^6*Ipss^5 + 3.1501e+593*Iass^6*Ipss^4 - 5.9515e+592*Iass^6*Ipss^3 + 3.7981e+591*Iass^6*Ipss^2 - 8.8293e+589*Iass^6*Ipss + 4.2059e+587*Iass^6 + 5.1019e+560*Iass^5*Ipss^6 + 1.2667e+593*Iass^5*Ipss^5 + 1.6876e+596*Iass^5*Ipss^4 - 3.4765e+595*Iass^5*Ipss^3 + 2.3793e+594*Iass^5*Ipss^2 - 5.815e+592*Iass^5*Ipss + 2.8083e+590*Iass^5 - 1.8834e+593*Iass^4*Ipss^6 - 1.0093e+596*Iass^4*Ipss^5 - 1.4371e+595*Iass^4*Ipss^4 + 5.6961e+594*Iass^4*Ipss^3 - 4.4255e+593*Iass^4*Ipss^2 + 1.1211e+592*Iass^4*Ipss - 5.4565e+589*Iass^4 + 4.5613e+592*Iass^3*Ipss^6 + 2.5843e+595*Iass^3*Ipss^5 - 2.6844e+594*Iass^3*Ipss^4 - 1.5628e+593*Iass^3*Ipss^3 + 2.4895e+592*Iass^3*Ipss^2 - 7.3032e+590*Iass^3*Ipss + 3.6627e+588*Iass^3 - 4.1152e+591*Iass^2*Ipss^6 - 2.4528e+594*Iass^2*Ipss^5 + 4.1508e+593*Iass^2*Ipss^4 - 1.7687e+592*Iass^2*Ipss^3 - 1.7954e+590*Iass^2*Ipss^2 + 1.7233e+589*Iass^2*Ipss - 9.7854e+586*Iass^2 + 1.6372e+590*Iass*Ipss^6 + 1.0214e+593*Iass*Ipss^5 - 1.935e+592*Iass*Ipss^4 + 1.1566e+591*Iass*Ipss^3 - 2.0727e+589*Iass*Ipss^2 - 4.6589e+586*Iass*Ipss + 8.5934e+584*Iass - 2.4194e+588*Ipss^6 - 1.5716e+591*Ipss^5 + 3.02e+590*Ipss^4 - 1.8883e+589*Ipss^3 + 4.0881e+587*Ipss^2 - 1.9022e+585*Ipss))/(- 2.0096e+791*IBss^8*Iass^6*Ipss^4 + 3.7741e+790*IBss^8*Iass^6*Ipss^3 - 2.5129e+789*IBss^8*Iass^6*Ipss^2 + 6.8477e+787*IBss^8*Iass^6*Ipss - 6.0318e+785*IBss^8*Iass^6 + 1.7188e+825*IBss^8*Iass^5*Ipss^4 - 3.228e+824*IBss^8*Iass^5*Ipss^3 + 2.1493e+823*IBss^8*Iass^5*Ipss^2 - 5.8568e+821*IBss^8*Iass^5*Ipss + 5.159e+819*IBss^8*Iass^5 - 5.1334e+824*IBss^8*Iass^4*Ipss^4 + 9.5499e+823*IBss^8*Iass^4*Ipss^3 - 6.2947e+822*IBss^8*Iass^4*Ipss^2 + 1.699e+821*IBss^8*Iass^4*Ipss - 1.4873e+819*IBss^8*Iass^4 + 6.1112e+823*IBss^8*Iass^3*Ipss^4 - 1.1252e+823*IBss^8*Iass^3*Ipss^3 + 7.3354e+821*IBss^8*Iass^3*Ipss^2 - 1.9591e+820*IBss^8*Iass^3*Ipss + 1.7032e+818*IBss^8*Iass^3 - 3.6245e+822*IBss^8*Iass^2*Ipss^4 + 6.5992e+821*IBss^8*Iass^2*Ipss^3 - 4.2502e+820*IBss^8*Iass^2*Ipss^2 + 1.1221e+819*IBss^8*Iass^2*Ipss - 9.6811e+816*IBss^8*Iass^2 + 1.0708e+821*IBss^8*Iass*Ipss^4 - 1.926e+820*IBss^8*Iass*Ipss^3 + 1.224e+819*IBss^8*Iass*Ipss^2 - 3.191e+817*IBss^8*Iass*Ipss + 2.7301e+815*IBss^8*Iass - 1.2606e+819*IBss^8*Ipss^4 + 2.2373e+818*IBss^8*Ipss^3 - 1.4015e+817*IBss^8*Ipss^2 + 3.604e+815*IBss^8*Ipss - 3.0552e+813*IBss^8 - 6.0289e+791*IBss^7*Iass^7*Ipss^4 + 1.1322e+791*IBss^7*Iass^7*Ipss^3 - 7.5386e+789*IBss^7*Iass^7*Ipss^2 + 2.0543e+788*IBss^7*Iass^7*Ipss - 1.8095e+786*IBss^7*Iass^7 + 3.9344e+790*IBss^7*Iass^6*Ipss^5 + 5.1565e+825*IBss^7*Iass^6*Ipss^4 - 9.684e+824*IBss^7*Iass^6*Ipss^3 + 6.4478e+823*IBss^7*Iass^6*Ipss^2 - 1.7571e+822*IBss^7*Iass^6*Ipss + 1.5477e+820*IBss^7*Iass^6 - 3.3651e+824*IBss^7*Iass^5*Ipss^5 + 9.2289e+827*IBss^7*Iass^5*Ipss^4 - 1.895e+827*IBss^7*Iass^5*Ipss^3 + 1.3466e+826*IBss^7*Iass^5*Ipss^2 - 3.8452e+824*IBss^7*Iass^5*Ipss + 3.4784e+822*IBss^7*Iass^5 + 1.005e+824*IBss^7*Iass^4*Ipss^5 - 2.8907e+827*IBss^7*Iass^4*Ipss^4 + 5.8644e+826*IBss^7*Iass^4*Ipss^3 - 4.1165e+825*IBss^7*Iass^4*Ipss^2 + 1.1622e+824*IBss^7*Iass^4*Ipss - 1.044e+822*IBss^7*Iass^4 - 1.1965e+823*IBss^7*Iass^3*Ipss^5 + 3.5932e+826*IBss^7*Iass^3*Ipss^4 - 7.1989e+825*IBss^7*Iass^3*Ipss^3 + 4.9875e+824*IBss^7*Iass^3*Ipss^2 - 1.3911e+823*IBss^7*Iass^3*Ipss + 1.2401e+821*IBss^7*Iass^3 + 7.0961e+821*IBss^7*Iass^2*Ipss^5 - 2.2171e+825*IBss^7*Iass^2*Ipss^4 + 4.3837e+824*IBss^7*Iass^2*Ipss^3 - 2.9947e+823*IBss^7*Iass^2*Ipss^2 + 8.2437e+821*IBss^7*Iass^2*Ipss - 7.2878e+819*IBss^7*Iass^2 - 2.0965e+820*IBss^7*Iass*Ipss^5 + 6.7921e+823*IBss^7*Iass*Ipss^4 - 1.3243e+823*IBss^7*Iass*Ipss^3 + 8.9108e+821*IBss^7*Iass*Ipss^2 - 2.4183e+820*IBss^7*Iass*Ipss + 2.1184e+818*IBss^7*Iass + 2.468e+818*IBss^7*Ipss^5 - 8.2664e+821*IBss^7*Ipss^4 + 1.588e+821*IBss^7*Ipss^3 - 1.0512e+820*IBss^7*Ipss^2 + 2.8095e+818*IBss^7*Ipss - 2.4366e+816*IBss^7 - 6.0289e+791*IBss^6*Iass^8*Ipss^4 + 1.1322e+791*IBss^6*Iass^8*Ipss^3 - 7.5386e+789*IBss^6*Iass^8*Ipss^2 + 2.0543e+788*IBss^6*Iass^8*Ipss - 1.8095e+786*IBss^6*Iass^8 + 7.8688e+790*IBss^6*Iass^7*Ipss^5 + 5.1565e+825*IBss^6*Iass^7*Ipss^4 - 9.684e+824*IBss^6*Iass^7*Ipss^3 + 6.4478e+823*IBss^6*Iass^7*Ipss^2 - 1.7571e+822*IBss^6*Iass^7*Ipss + 1.5477e+820*IBss^6*Iass^7 + 1.6847e+791*IBss^6*Iass^6*Ipss^6 - 6.7302e+824*IBss^6*Iass^6*Ipss^5 + 1.8467e+828*IBss^6*Iass^6*Ipss^4 - 3.7917e+827*IBss^6*Iass^6*Ipss^3 + 2.6944e+826*IBss^6*Iass^6*Ipss^2 - 7.6935e+824*IBss^6*Iass^6*Ipss + 6.9595e+822*IBss^6*Iass^6 - 1.4409e+825*IBss^6*Iass^5*Ipss^6 - 1.0634e+828*IBss^6*Iass^5*Ipss^5 + 1.3089e+829*IBss^6*Iass^5*Ipss^4 - 2.9093e+828*IBss^6*Iass^5*Ipss^3 + 2.2688e+827*IBss^6*Iass^5*Ipss^2 - 7.0448e+825*IBss^6*Iass^5*Ipss + 6.6945e+823*IBss^6*Iass^5 + 4.3033e+824*IBss^6*Iass^4*Ipss^6 + 3.3278e+827*IBss^6*Iass^4*Ipss^5 - 4.3935e+828*IBss^6*Iass^4*Ipss^4 + 9.6818e+827*IBss^6*Iass^4*Ipss^3 - 7.4816e+826*IBss^6*Iass^4*Ipss^2 + 2.3032e+825*IBss^6*Iass^4*Ipss - 2.1777e+823*IBss^6*Iass^4 - 5.1231e+823*IBss^6*Iass^3*Ipss^6 - 4.1347e+826*IBss^6*Iass^3*Ipss^5 + 5.7642e+827*IBss^6*Iass^3*Ipss^4 - 1.2611e+827*IBss^6*Iass^3*Ipss^3 + 9.6689e+825*IBss^6*Iass^3*Ipss^2 - 2.9533e+824*IBss^6*Iass^3*Ipss + 2.7796e+822*IBss^6*Iass^3 + 3.0385e+822*IBss^6*Iass^2*Ipss^6 + 2.55e+825*IBss^6*Iass^2*Ipss^5 - 3.7465e+826*IBss^6*Iass^2*Ipss^4 + 8.1387e+825*IBss^6*Iass^2*Ipss^3 - 6.1906e+824*IBss^6*Iass^2*Ipss^2 + 1.8756e+823*IBss^6*Iass^2*Ipss - 1.7566e+821*IBss^6*Iass^2 - 8.9768e+820*IBss^6*Iass*Ipss^6 - 7.8071e+823*IBss^6*Iass*Ipss^5 + 1.2095e+825*IBss^6*Iass*Ipss^4 - 2.6082e+824*IBss^6*Iass*Ipss^3 + 1.9673e+823*IBss^6*Iass*Ipss^2 - 5.9076e+821*IBss^6*Iass*Ipss + 5.5029e+819*IBss^6*Iass + 1.0568e+819*IBss^6*Ipss^6 + 9.4945e+821*IBss^6*Ipss^5 - 1.5526e+823*IBss^6*Ipss^4 + 3.3221e+822*IBss^6*Ipss^3 - 2.4829e+821*IBss^6*Ipss^2 + 7.3823e+819*IBss^6*Ipss - 6.834e+817*IBss^6 - 2.0096e+791*IBss^5*Iass^9*Ipss^4 + 3.7741e+790*IBss^5*Iass^9*Ipss^3 - 2.5129e+789*IBss^5*Iass^9*Ipss^2 + 6.8477e+787*IBss^5*Iass^9*Ipss - 6.0318e+785*IBss^5*Iass^9 + 3.9344e+790*IBss^5*Iass^8*Ipss^5 + 1.7188e+825*IBss^5*Iass^8*Ipss^4 - 3.228e+824*IBss^5*Iass^8*Ipss^3 + 2.1493e+823*IBss^5*Iass^8*Ipss^2 - 5.8568e+821*IBss^5*Iass^8*Ipss + 5.159e+819*IBss^5*Iass^8 + 1.6847e+791*IBss^5*Iass^7*Ipss^6 - 3.3651e+824*IBss^5*Iass^7*Ipss^5 + 9.2272e+827*IBss^5*Iass^7*Ipss^4 - 1.8947e+827*IBss^5*Iass^7*Ipss^3 + 1.3464e+826*IBss^5*Iass^7*Ipss^2 - 3.8446e+824*IBss^5*Iass^7*Ipss + 3.4779e+822*IBss^5*Iass^7 - 7.1838e+790*IBss^5*Iass^6*Ipss^7 - 1.4409e+825*IBss^5*Iass^6*Ipss^6 - 1.0634e+828*IBss^5*Iass^6*Ipss^5 + 1.3055e+829*IBss^5*Iass^6*Ipss^4 - 2.9024e+828*IBss^5*Iass^6*Ipss^3 + 2.2639e+827*IBss^5*Iass^6*Ipss^2 - 7.0313e+825*IBss^5*Iass^6*Ipss + 6.6825e+823*IBss^5*Iass^6 + 6.1443e+824*IBss^5*Iass^5*Ipss^7 + 3.063e+827*IBss^5*Iass^5*Ipss^6 - 7.5885e+828*IBss^5*Iass^5*Ipss^5 - 6.1238e+828*IBss^5*Iass^5*Ipss^4 + 1.6085e+828*IBss^5*Iass^5*Ipss^3 - 1.307e+827*IBss^5*Iass^5*Ipss^2 + 4.1106e+825*IBss^5*Iass^5*Ipss - 3.9223e+823*IBss^5*Iass^5 - 1.835e+824*IBss^5*Iass^4*Ipss^7 - 9.5761e+826*IBss^5*Iass^4*Ipss^6 + 2.5508e+828*IBss^5*Iass^4*Ipss^5 + 1.194e+828*IBss^5*Iass^4*Ipss^4 - 3.4471e+827*IBss^5*Iass^4*Ipss^3 + 2.8476e+826*IBss^5*Iass^4*Ipss^2 - 8.9731e+824*IBss^5*Iass^4*Ipss + 8.5576e+822*IBss^5*Iass^4 + 2.1846e+823*IBss^5*Iass^3*Ipss^7 + 1.1891e+826*IBss^5*Iass^3*Ipss^6 - 3.3519e+827*IBss^5*Iass^3*Ipss^5 - 1.2382e+827*IBss^5*Iass^3*Ipss^4 + 3.7613e+826*IBss^5*Iass^3*Ipss^3 - 3.1212e+825*IBss^5*Iass^3*Ipss^2 + 9.802e+823*IBss^5*Iass^3*Ipss - 9.321e+821*IBss^5*Iass^3 - 1.2957e+822*IBss^5*Iass^2*Ipss^7 - 7.3294e+824*IBss^5*Iass^2*Ipss^6 + 2.1821e+826*IBss^5*Iass^2*Ipss^5 + 7.1678e+825*IBss^5*Iass^2*Ipss^4 - 2.2296e+825*IBss^5*Iass^2*Ipss^3 + 1.8459e+824*IBss^5*Iass^2*Ipss^2 - 5.7582e+822*IBss^5*Iass^2*Ipss + 5.4509e+820*IBss^5*Iass^2 + 3.8279e+820*IBss^5*Iass*Ipss^7 + 2.2424e+823*IBss^5*Iass*Ipss^6 - 7.0556e+824*IBss^5*Iass*Ipss^5 - 2.1881e+824*IBss^5*Iass*Ipss^4 + 6.8509e+823*IBss^5*Iass*Ipss^3 - 5.6343e+822*IBss^5*Iass*Ipss^2 + 1.7412e+821*IBss^5*Iass*Ipss - 1.6385e+819*IBss^5*Iass - 4.5064e+818*IBss^5*Ipss^7 - 2.7247e+821*IBss^5*Ipss^6 + 9.0711e+822*IBss^5*Ipss^5 + 2.745e+822*IBss^5*Ipss^4 - 8.5571e+821*IBss^5*Ipss^3 + 6.9678e+820*IBss^5*Ipss^2 - 2.1281e+819*IBss^5*Ipss + 1.9877e+817*IBss^5 + 6.9778e+790*IBss^4*Iass^9*Ipss^4 - 1.2948e+790*IBss^4*Iass^9*Ipss^3 + 8.5111e+788*IBss^4*Iass^9*Ipss^2 - 2.2911e+787*IBss^4*Iass^9*Ipss + 2.0021e+785*IBss^4*Iass^9 - 1.3661e+790*IBss^4*Iass^8*Ipss^5 - 5.9681e+824*IBss^4*Iass^8*Ipss^4 + 1.1074e+824*IBss^4*Iass^8*Ipss^3 - 7.2796e+822*IBss^4*Iass^8*Ipss^2 + 1.9596e+821*IBss^4*Iass^8*Ipss - 1.7124e+819*IBss^4*Iass^8 - 5.8495e+790*IBss^4*Iass^7*Ipss^6 + 1.1684e+824*IBss^4*Iass^7*Ipss^5 - 3.2305e+827*IBss^4*Iass^7*Ipss^4 + 6.5588e+826*IBss^4*Iass^7*Ipss^3 - 4.6002e+825*IBss^4*Iass^7*Ipss^2 + 1.2969e+824*IBss^4*Iass^7*Ipss - 1.1638e+822*IBss^4*Iass^7 + 2.4944e+790*IBss^4*Iass^6*Ipss^7 + 5.0031e+824*IBss^4*Iass^6*Ipss^6 + 3.7116e+827*IBss^4*Iass^6*Ipss^5 - 4.9724e+828*IBss^4*Iass^6*Ipss^4 + 1.098e+828*IBss^4*Iass^6*Ipss^3 - 8.4933e+826*IBss^4*Iass^6*Ipss^2 + 2.615e+825*IBss^4*Iass^6*Ipss - 2.4725e+823*IBss^4*Iass^6 - 2.1334e+824*IBss^4*Iass^5*Ipss^7 - 1.0656e+827*IBss^4*Iass^5*Ipss^6 + 2.8949e+828*IBss^4*Iass^5*Ipss^5 + 1.2366e+828*IBss^4*Iass^5*Ipss^4 - 3.6546e+827*IBss^4*Iass^5*Ipss^3 + 3.0356e+826*IBss^4*Iass^5*Ipss^2 - 9.585e+824*IBss^4*Iass^5*Ipss + 9.1484e+822*IBss^4*Iass^5 + 6.3343e+823*IBss^4*Iass^4*Ipss^7 + 3.3117e+826*IBss^4*Iass^4*Ipss^6 - 9.6699e+827*IBss^4*Iass^4*Ipss^5 - 9.8531e+826*IBss^4*Iass^4*Ipss^4 + 5.1858e+826*IBss^4*Iass^4*Ipss^3 - 4.6882e+825*IBss^4*Iass^4*Ipss^2 + 1.52e+824*IBss^4*Iass^4*Ipss - 1.4643e+822*IBss^4*Iass^4 - 7.493e+822*IBss^4*Iass^3*Ipss^7 - 4.0853e+825*IBss^4*Iass^3*Ipss^6 + 1.2643e+827*IBss^4*Iass^3*Ipss^5 + 1.0446e+825*IBss^4*Iass^3*Ipss^4 - 4.1203e+825*IBss^4*Iass^3*Ipss^3 + 4.0507e+824*IBss^4*Iass^3*Ipss^2 - 1.3417e+823*IBss^4*Iass^3*Ipss + 1.3011e+821*IBss^4*Iass^3 + 4.4134e+821*IBss^4*Iass^2*Ipss^7 + 2.4997e+824*IBss^4*Iass^2*Ipss^6 - 8.19e+825*IBss^4*Iass^2*Ipss^5 + 2.4064e+824*IBss^4*Iass^2*Ipss^4 + 1.9616e+824*IBss^4*Iass^2*Ipss^3 - 2.0593e+823*IBss^4*Iass^2*Ipss^2 + 6.9011e+821*IBss^4*Iass^2*Ipss - 6.7073e+819*IBss^4*Iass^2 - 1.2941e+820*IBss^4*Iass*Ipss^7 - 7.5865e+822*IBss^4*Iass*Ipss^6 + 2.6345e+824*IBss^4*Iass*Ipss^5 - 1.1867e+823*IBss^4*Iass*Ipss^4 - 5.3023e+822*IBss^4*Iass*Ipss^3 + 5.771e+821*IBss^4*Iass*Ipss^2 - 1.9347e+820*IBss^4*Iass*Ipss + 1.8754e+818*IBss^4*Iass + 1.5112e+818*IBss^4*Ipss^7 + 9.1372e+820*IBss^4*Ipss^6 - 3.3682e+822*IBss^4*Ipss^5 + 1.685e+821*IBss^4*Ipss^4 + 6.2595e+820*IBss^4*Ipss^3 - 6.8555e+819*IBss^4*Ipss^2 + 2.2749e+818*IBss^4*Ipss - 2.1886e+816*IBss^4 - 9.6373e+789*IBss^3*Iass^9*Ipss^4 + 1.7645e+789*IBss^3*Iass^9*Ipss^3 - 1.1432e+788*IBss^3*Iass^9*Ipss^2 + 3.0351e+786*IBss^3*Iass^9*Ipss - 2.6283e+784*IBss^3*Iass^9 + 1.8868e+789*IBss^3*Iass^8*Ipss^5 + 8.2428e+823*IBss^3*Iass^8*Ipss^4 - 1.5092e+823*IBss^3*Iass^8*Ipss^3 + 9.7778e+821*IBss^3*Iass^8*Ipss^2 - 2.5959e+820*IBss^3*Iass^8*Ipss + 2.248e+818*IBss^3*Iass^8 + 8.079e+789*IBss^3*Iass^7*Ipss^6 - 1.6138e+823*IBss^3*Iass^7*Ipss^5 + 4.4948e+826*IBss^3*Iass^7*Ipss^4 - 9.0124e+825*IBss^3*Iass^7*Ipss^3 + 6.2292e+824*IBss^3*Iass^7*Ipss^2 - 1.7311e+823*IBss^3*Iass^7*Ipss + 1.5392e+821*IBss^3*Iass^7 - 3.4451e+789*IBss^3*Iass^6*Ipss^7 - 6.91e+823*IBss^3*Iass^6*Ipss^6 - 5.1508e+826*IBss^3*Iass^6*Ipss^5 + 7.3958e+827*IBss^3*Iass^6*Ipss^4 - 1.6238e+827*IBss^3*Iass^6*Ipss^3 + 1.2468e+826*IBss^3*Iass^6*Ipss^2 - 3.8075e+824*IBss^3*Iass^6*Ipss + 3.5821e+822*IBss^3*Iass^6 + 2.9466e+823*IBss^3*Iass^5*Ipss^7 + 1.4746e+826*IBss^3*Iass^5*Ipss^6 - 4.3133e+827*IBss^3*Iass^5*Ipss^5 - 1.3742e+827*IBss^3*Iass^5*Ipss^4 + 4.3715e+826*IBss^3*Iass^5*Ipss^3 - 3.6632e+825*IBss^3*Iass^5*Ipss^2 + 1.1541e+824*IBss^3*Iass^5*Ipss - 1.0987e+822*IBss^3*Iass^5 - 8.6917e+822*IBss^3*Iass^4*Ipss^7 - 4.5526e+825*IBss^3*Iass^4*Ipss^6 + 1.4334e+827*IBss^3*Iass^4*Ipss^5 - 2.6055e+824*IBss^3*Iass^4*Ipss^4 - 4.3666e+825*IBss^3*Iass^4*Ipss^3 + 4.3657e+824*IBss^3*Iass^4*Ipss^2 - 1.4531e+823*IBss^3*Iass^4*Ipss + 1.4117e+821*IBss^3*Iass^4 + 1.0209e+822*IBss^3*Iass^3*Ipss^7 + 5.5748e+824*IBss^3*Iass^3*Ipss^6 - 1.8662e+826*IBss^3*Iass^3*Ipss^5 + 1.7496e+825*IBss^3*Iass^3*Ipss^4 + 1.8888e+824*IBss^3*Iass^3*Ipss^3 - 2.7664e+823*IBss^3*Iass^3*Ipss^2 + 1.0015e+822*IBss^3*Iass^3*Ipss - 1.0008e+820*IBss^3*Iass^3 - 5.9665e+820*IBss^3*Iass^2*Ipss^7 - 3.3831e+823*IBss^3*Iass^2*Ipss^6 + 1.2036e+825*IBss^3*Iass^2*Ipss^5 - 1.5681e+824*IBss^3*Iass^2*Ipss^4 - 2.3767e+822*IBss^3*Iass^2*Ipss^3 + 1.0244e+822*IBss^3*Iass^2*Ipss^2 - 4.1214e+820*IBss^3*Iass^2*Ipss + 4.2455e+818*IBss^3*Iass^2 + 1.7348e+819*IBss^3*Iass*Ipss^7 + 1.0175e+822*IBss^3*Iass*Ipss^6 - 3.8532e+823*IBss^3*Iass*Ipss^5 + 5.5876e+822*IBss^3*Iass*Ipss^4 - 5.4909e+820*IBss^3*Iass*Ipss^3 - 2.2278e+820*IBss^3*Iass*Ipss^2 + 9.8445e+818*IBss^3*Iass*Ipss - 1.0362e+817*IBss^3*Iass - 2.0072e+817*IBss^3*Ipss^7 - 1.2132e+820*IBss^3*Ipss^6 + 4.9004e+821*IBss^3*Ipss^5 - 7.3062e+820*IBss^3*Ipss^4 + 1.2552e+819*IBss^3*Ipss^3 + 2.3106e+818*IBss^3*Ipss^2 - 1.0627e+817*IBss^3*Ipss + 1.1229e+815*IBss^3 + 6.616e+788*IBss^2*Iass^9*Ipss^4 - 1.1933e+788*IBss^2*Iass^9*Ipss^3 + 7.6063e+786*IBss^2*Iass^9*Ipss^2 - 1.9881e+785*IBss^2*Iass^9*Ipss + 1.7038e+783*IBss^2*Iass^9 - 1.2953e+788*IBss^2*Iass^8*Ipss^5 - 5.6587e+822*IBss^2*Iass^8*Ipss^4 + 1.0206e+822*IBss^2*Iass^8*Ipss^3 - 6.5057e+820*IBss^2*Iass^8*Ipss^2 + 1.7004e+819*IBss^2*Iass^8*Ipss - 1.4573e+817*IBss^2*Iass^8 - 5.5463e+788*IBss^2*Iass^7*Ipss^6 + 1.1079e+822*IBss^2*Iass^7*Ipss^5 - 3.1074e+825*IBss^2*Iass^7*Ipss^4 + 6.1449e+824*IBss^2*Iass^7*Ipss^3 - 4.1779e+823*IBss^2*Iass^7*Ipss^2 + 1.1425e+822*IBss^2*Iass^7*Ipss - 1.0053e+820*IBss^2*Iass^7 + 2.365e+788*IBss^2*Iass^6*Ipss^7 + 4.7437e+822*IBss^2*Iass^6*Ipss^6 + 3.5515e+825*IBss^2*Iass^6*Ipss^5 - 5.4525e+826*IBss^2*Iass^6*Ipss^4 + 1.1899e+826*IBss^2*Iass^6*Ipss^3 - 9.064e+824*IBss^2*Iass^6*Ipss^2 + 2.7429e+823*IBss^2*Iass^6*Ipss - 2.5659e+821*IBss^2*Iass^6 - 2.0228e+822*IBss^2*Iass^5*Ipss^7 - 1.0138e+825*IBss^2*Iass^5*Ipss^6 + 3.1855e+826*IBss^2*Iass^5*Ipss^5 + 8.7515e+825*IBss^2*Iass^5*Ipss^4 - 2.8942e+825*IBss^2*Iass^5*Ipss^3 + 2.4251e+824*IBss^2*Iass^5*Ipss^2 - 7.5886e+822*IBss^2*Iass^5*Ipss + 7.1894e+820*IBss^2*Iass^5 + 5.9239e+821*IBss^2*Iass^4*Ipss^7 + 3.1069e+824*IBss^2*Iass^4*Ipss^6 - 1.0532e+826*IBss^2*Iass^4*Ipss^5 + 4.4882e+824*IBss^2*Iass^4*Ipss^4 + 2.234e+824*IBss^2*Iass^4*Ipss^3 - 2.4357e+823*IBss^2*Iass^4*Ipss^2 + 8.2399e+821*IBss^2*Iass^4*Ipss - 8.0339e+819*IBss^2*Iass^4 - 6.9028e+820*IBss^2*Iass^3*Ipss^7 - 3.7732e+823*IBss^2*Iass^3*Ipss^6 + 1.3651e+825*IBss^2*Iass^3*Ipss^5 - 1.8221e+824*IBss^2*Iass^3*Ipss^4 - 1.8207e+822*IBss^2*Iass^3*Ipss^3 + 1.0991e+822*IBss^2*Iass^3*Ipss^2 - 4.4905e+820*IBss^2*Iass^3*Ipss + 4.6456e+818*IBss^2*Iass^3 + 3.9994e+819*IBss^2*Iass^2*Ipss^7 + 2.2687e+822*IBss^2*Iass^2*Ipss^6 - 8.7632e+823*IBss^2*Iass^2*Ipss^5 + 1.4832e+823*IBss^2*Iass^2*Ipss^4 - 5.6996e+821*IBss^2*Iass^2*Ipss^3 - 1.8406e+820*IBss^2*Iass^2*Ipss^2 + 1.3129e+819*IBss^2*Iass^2*Ipss - 1.5187e+817*IBss^2*Iass^2 - 1.1518e+818*IBss^2*Iass*Ipss^7 - 6.7531e+820*IBss^2*Iass*Ipss^6 + 2.791e+822*IBss^2*Iass*Ipss^5 - 5.1164e+821*IBss^2*Iass*Ipss^4 + 2.6925e+820*IBss^2*Iass*Ipss^3 - 9.415e+817*IBss^2*Iass*Ipss^2 - 2.0887e+817*IBss^2*Iass*Ipss + 2.8607e+815*IBss^2*Iass + 1.3189e+816*IBss^2*Ipss^7 + 7.9611e+818*IBss^2*Ipss^6 - 3.5288e+820*IBss^2*Ipss^5 + 6.5901e+819*IBss^2*Ipss^4 - 3.7219e+818*IBss^2*Ipss^3 + 4.0358e+816*IBss^2*Ipss^2 + 1.6473e+815*IBss^2*Ipss - 2.6121e+813*IBss^2 - 2.257e+787*IBss*Iass^9*Ipss^4 + 4.0032e+786*IBss*Iass^9*Ipss^3 - 2.5053e+785*IBss*Iass^9*Ipss^2 + 6.4345e+783*IBss*Iass^9*Ipss - 5.4497e+781*IBss*Iass^9 + 4.4187e+786*IBss*Iass^8*Ipss^5 + 1.9304e+821*IBss*Iass^8*Ipss^4 - 3.4239e+820*IBss*Iass^8*Ipss^3 + 2.1428e+819*IBss*Iass^8*Ipss^2 - 5.5034e+817*IBss*Iass^8*Ipss + 4.6611e+815*IBss*Iass^8 + 1.892e+787*IBss*Iass^7*Ipss^6 - 3.7793e+820*IBss*Iass^7*Ipss^5 + 1.0673e+824*IBss*Iass^7*Ipss^4 - 2.0783e+823*IBss*Iass^7*Ipss^3 + 1.3872e+822*IBss*Iass^7*Ipss^2 - 3.7256e+820*IBss*Iass^7*Ipss + 3.2395e+818*IBss*Iass^7 - 8.0681e+786*IBss*Iass^6*Ipss^7 - 1.6183e+821*IBss*Iass^6*Ipss^6 - 1.2163e+824*IBss*Iass^6*Ipss^5 + 1.9976e+825*IBss*Iass^6*Ipss^4 - 4.3298e+824*IBss*Iass^6*Ipss^3 + 3.2684e+823*IBss*Iass^6*Ipss^2 - 9.7864e+821*IBss*Iass^6*Ipss + 9.0941e+819*IBss*Iass^6 + 6.9006e+820*IBss*Iass^5*Ipss^7 + 3.4612e+823*IBss*Iass^5*Ipss^6 - 1.169e+825*IBss*Iass^5*Ipss^5 - 2.992e+824*IBss*Iass^5*Ipss^4 + 1.0034e+824*IBss*Iass^5*Ipss^3 - 8.3541e+822*IBss*Iass^5*Ipss^2 + 2.5861e+821*IBss*Iass^5*Ipss - 2.4327e+819*IBss*Iass^5 - 2.0047e+820*IBss*Iass^4*Ipss^7 - 1.0521e+823*IBss*Iass^4*Ipss^6 + 3.8445e+824*IBss*Iass^4*Ipss^5 - 2.2644e+823*IBss*Iass^4*Ipss^4 - 6.654e+822*IBss*Iass^4*Ipss^3 + 7.6288e+821*IBss*Iass^4*Ipss^2 - 2.5849e+820*IBss*Iass^4*Ipss + 2.5132e+818*IBss*Iass^4 + 2.3155e+819*IBss*Iass^3*Ipss^7 + 1.266e+822*IBss*Iass^3*Ipss^6 - 4.9587e+823*IBss*Iass^3*Ipss^5 + 7.3814e+822*IBss*Iass^3*Ipss^4 - 1.0768e+821*IBss*Iass^3*Ipss^3 - 2.6055e+820*IBss*Iass^3*Ipss^2 + 1.1897e+819*IBss*Iass^3*Ipss - 1.2617e+817*IBss*Iass^3 - 1.3285e+818*IBss*Iass^2*Ipss^7 - 7.5334e+820*IBss*Iass^2*Ipss^6 + 3.1665e+822*IBss*Iass^2*Ipss^5 - 5.8259e+821*IBss*Iass^2*Ipss^4 + 3.0901e+820*IBss*Iass^2*Ipss^3 - 1.295e+818*IBss*Iass^2*Ipss^2 - 2.3031e+817*IBss*Iass^2*Ipss + 3.1828e+815*IBss*Iass^2 + 3.7854e+816*IBss*Iass*Ipss^7 + 2.2165e+819*IBss*Iass*Ipss^6 - 1.0025e+821*IBss*Iass*Ipss^5 + 1.9795e+820*IBss*Iass*Ipss^4 - 1.2736e+819*IBss*Iass*Ipss^3 + 2.6427e+817*IBss*Iass*Ipss^2 + 6.1427e+814*IBss*Iass*Ipss - 3.8786e+813*IBss*Iass - 4.2841e+814*IBss*Ipss^7 - 2.5796e+817*IBss*Ipss^6 + 1.259e+819*IBss*Ipss^5 - 2.5216e+818*IBss*Ipss^4 + 1.6908e+817*IBss*Ipss^3 - 4.1198e+815*IBss*Ipss^2 + 1.9691e+813*IBss*Ipss + 2.1283e+811*IBss + 3.0603e+785*Iass^9*Ipss^4 - 5.3278e+784*Iass^9*Ipss^3 + 3.2667e+783*Iass^9*Ipss^2 - 8.229e+781*Iass^9*Ipss + 6.878e+779*Iass^9 - 5.9915e+784*Iass^8*Ipss^5 - 2.6175e+819*Iass^8*Ipss^4 + 4.5569e+818*Iass^8*Ipss^3 - 2.794e+817*Iass^8*Ipss^2 + 7.0382e+815*Iass^8*Ipss - 5.8827e+813*Iass^8 - 2.5655e+785*Iass^7*Ipss^6 + 5.1245e+818*Iass^7*Ipss^5 - 1.4567e+822*Iass^7*Ipss^4 + 2.7888e+821*Iass^7*Ipss^3 - 1.8233e+820*Iass^7*Ipss^2 + 4.8007e+818*Iass^7*Ipss - 4.1194e+816*Iass^7 + 1.094e+785*Iass^6*Ipss^7 + 2.1943e+819*Iass^6*Ipss^6 + 1.6551e+822*Iass^6*Ipss^5 - 2.9104e+823*Iass^6*Ipss^4 + 6.259e+822*Iass^6*Ipss^3 - 4.675e+821*Iass^6*Ipss^2 + 1.3822e+820*Iass^6*Ipss - 1.274e+818*Iass^6 - 9.3568e+818*Iass^5*Ipss^7 - 4.6935e+821*Iass^5*Ipss^6 + 1.7061e+823*Iass^5*Ipss^5 + 4.2498e+822*Iass^5*Ipss^4 - 1.4212e+822*Iass^5*Ipss^3 + 1.1698e+821*Iass^5*Ipss^2 - 3.5674e+819*Iass^5*Ipss + 3.3235e+817*Iass^5 + 2.6944e+818*Iass^4*Ipss^7 + 1.4139e+821*Iass^4*Ipss^6 - 5.5779e+822*Iass^4*Ipss^5 + 3.5364e+821*Iass^4*Ipss^4 + 8.8699e+820*Iass^4*Ipss^3 - 1.0251e+820*Iass^4*Ipss^2 + 3.427e+818*Iass^4*Ipss - 3.2994e+816*Iass^4 - 3.0817e+817*Iass^3*Ipss^7 - 1.684e+820*Iass^3*Ipss^6 + 7.1549e+821*Iass^3*Ipss^5 - 1.0911e+821*Iass^3*Ipss^4 + 2.2845e+819*Iass^3*Ipss^3 + 3.0509e+818*Iass^3*Ipss^2 - 1.4497e+817*Iass^3*Ipss + 1.5403e+815*Iass^3 + 1.7492e+816*Iass^2*Ipss^7 + 9.9059e+818*Iass^2*Ipss^6 - 4.5414e+820*Iass^2*Ipss^5 + 8.4915e+819*Iass^2*Ipss^4 - 4.7913e+818*Iass^2*Ipss^3 + 5.1819e+816*Iass^2*Ipss^2 + 2.0939e+815*Iass^2*Ipss - 3.3148e+813*Iass^2 - 4.925e+814*Iass*Ipss^7 - 2.8771e+817*Iass*Ipss^6 + 1.4279e+819*Iass*Ipss^5 - 2.8537e+818*Iass*Ipss^4 + 1.9024e+817*Iass*Ipss^3 - 4.567e+815*Iass*Ipss^2 + 2.0052e+813*Iass*Ipss + 2.5668e+811*Iass + 5.5019e+812*Ipss^7 + 3.3009e+815*Ipss^6 - 1.7787e+817*Ipss^5 + 3.5944e+816*Ipss^4 - 2.4737e+815*Ipss^3 + 6.6123e+813*Ipss^2 - 5.6818e+811*Ipss)
and:
ddd_fcn =
function_handle with value:
@(IBss,Iass,Ipss)((Iass.*1.048540363105328e+34-6.864833896770217e+32).* ...
that will not fit so I truncated it. That will allow you to substitute the variables and get the array you need.
raha ahmadi
2020-7-6
Dear Star Strider
I’m so grateful for your help. I couldnt handle my ODE coeffisients but now I can. it seems its implicit form and very different numbers ( numbers as large as 1e24 and also as small as 1e-9 make it difficult to solve. They are also nolinear. Odes15s seems does not work anymore ( In fact I m new in soling ODE.) I'll use your and John' s points. I' ll try one more time .
I hope best wishes for you and john
Raha
Star Strider
2020-7-6
As always, our pleasure!
The ‘dde_fcn’ with the appropriate arguments should allow you to calculate ‘ddd’. I did not try that myself because I have no idea what the arguments should be.
更多回答(1 个)
John D'Errico
2020-7-5
Your question seems to be one of solving a system of ODEs with a mass matrix. For example in the help for ODE45, we see this option:
ode45 can solve problems M(t,y)*y' = f(t,y) with mass matrix M that is
nonsingular. Use ODESET to set the 'Mass' property to a function handle
MASS if MASS(T,Y) returns the value of the mass matrix. If the mass matrix
is constant, the matrix can be used as the value of the 'Mass' option. If
the mass matrix does not depend on the state variable Y and the function
MASS is to be called with one input argument T, set 'MStateDependence' to
'none'. ODE15S and ODE23T can solve problems with singular mass matrices.
That seems to be the question you have posed.
H11*dydt+H12*dxdt=H13
H21*dydt+H22*dxdt=H23
So all you need to do is pass in the mass matrix as
M = [H11, H12;H21, H22];
If you compute these elements in symbolic form but still just numbers, then use double to convert them to double precsion numbers.
1 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!发生错误
由于页面发生更改,无法完成操作。请重新加载页面以查看其更新后的状态。
您也可以从以下列表中选择网站:
如何获得最佳网站性能
选择中国网站(中文或英文)以获得最佳网站性能。其他 MathWorks 国家/地区网站并未针对您所在位置的访问进行优化。
美洲
- América Latina (Español)
- Canada (English)
- United States (English)
欧洲
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom(English)
亚太
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)