errors in calling function
이전 댓글 표시
Hi,I'm trying to call the myequation in another script ,here is my first script:
clear all,clc
%constantes :::
%temperature
T=17000;
%constante k: recombinaison radiative
k1=2.36 *(1e-12)*(T/300)^(-0.29)*exp(17.60/T);
k2 =3.24 *(1e-12) *(T/300)^(-0.66);
k3=3.50 *(1e-12)*(T/300)^(-0.53)* exp(3.20/T);
%constante k::Reconmbinaison
k4=1*10^(-10);
k5=4.54*10^(-10);
k6=9.10*10^(-10);
k7=8.98*10^(-9)*(T/300)^(-0.5);
k8=5.56*10^(-11)*(T/300)^(0.41)*exp(26.90/T);
k9=4.98*10^(-10)*exp(-18116/T);
k10=6*10^(-11)*(T/300)^(-0.16);
k11=8.69*10^(-11)*exp(-22600/T);
k12=1*10^(-10);
k13=2.94*10^(-11)*(T/300)^(0.5)*exp(-58025/T);
k14=1*10^(-9);
k15=4.80*10^(-10)
k16=1.00*10^(-9)*(T/300)^(-0.5);
k17=8.30*10^(-10);
k18=2.42*10^(-12)*(T/300)^(-0.21)*exp(44/T);
k19=1.18*10^(-11)*exp(-20413/T);
k20=9.82*10^(-12)*(T/300)^(-0.21)*exp(-5.20/T);
k21=1.66*10^(-10)*exp(-14100/T);
k22=1.15*10^(-10)*exp(-13400/T);
k23=2.51*10^(-10)*exp(-38602/T);
k24=1.00*10^(-10);
k25=5.37*10^(-11)*exp(-13800/T);
k26=5.00*10^(-11)*exp(-200/T);
k27=5.00*10^(-12)*exp(-900/T);
k28=6.00*10^(-12);
k29=2.00*10^(-10)*(T/300)^(-0.12);
k30=1.30*10^(-10);
k31=5.00*10^(-10);
k32=2.90*10^(-10);
k33=3.10*10^(-10);
k34=3.66*10^(-11);
k35=2.26*10^(-12)*(T/300)^(0.86)*exp(-3134.0/T);
k36=1.80*10^(-10);
k37=3.38*10^(-11)*(T/300)^(-0.17)*exp(2.80/T);
k38=3*10^(-12);
k39=1.00*10^(-10)*(T/300)^(0.4);
k40=6.10*10^(-10);
k41=1*10^(-13);
k42=5*10^(-11);
k43=5*10^(-11);
k44=1.15*10^(-10);
k45=4*10^(-11);
k46=1*10^(-10);
k47=1.50*10^(-11)*exp(-4300/T);
k48=5.3137*10^(-10);
k49=4.10*10^(-10);
k50=8*10^(-10);
k51=1.32*10^(-12);
k52=2.8*10^(-12)*exp(23400/T);
k53=1*10^(-10)*exp(-55200/T);
k54=8.60*10^(-11);
k55=5.99*10^(-12)*exp(-24075/T);
k56=2.02*10^(-11)*(T/300)^(-0.19)*exp(31.90/T);
k57=2.51*10^(-11)*exp(-30653/T);
k58=3.47*10^(-11)*(T/300)^(-1.33)*exp(-242/T);
%constante k::Recombinaison dissociative
k59=3*10^(-7)*(T/300)^(-0.5);
k60=1.70*10^(-7)*(T/300)^(-0.3);
%constante k::Association radiative
k61=4.01*10^(-18)*(T/300)^(0.17)*exp(-101.50/T);
k62=1.08*10^(-18)*(T/300)^(0.07)*exp(-57.50/T);
k63=3.14*10^(-18)*(T/300)^(-0.15)*exp(-68/T);
k64=4.00*10^(-14)*(T/300)^(-1);
k65=4.69*10^(-19)*(T/300)^(1.52)*exp(50.50/T);
k66=3*10^(-16)*(T/300)^(-1);
k67=5*10^(-10)*(T/300)^(-3.7)*exp(-800/T);
k68=4.36*10^(-18)*(T/300)^(0.35)*exp(-161.30/T);
k69=5.72*10^(-19)*(T/300)^(0.37)*exp(-51/T);
k70=4.90*10^(-20)*(T/300)^(1.58);
k71=3.71*10^(-18)*(T/300)^(0.24)*exp(-26.10/T);
%constante k::Transfert de charge
k72=7.05*10^(-10)*(T/300)^(-0.03)*exp(16.70/T);
k73=1*10^(-10);
k74=1.10*10^(-10);
k75=1.10*10^(-10);
k76=1.10*10^(-10);
k77=5.20*10^(-11);
k78=4.80*10^(-10);
k79=4.90*10^(-12)*(T/300)^(0.5)*exp(-4580/T);
k80=6.30*10^(-10);
k81=1.90*10^(-11);
k82=1.00*10^(-11);
k83=1.40*10^(-10);
k84=7.30*10^(-10)*exp(-890/T);
k85=1.00*10^(-9);
k86=3.11*10^(-10);
k87=1.00*10^(-9);
k88=8.25*10^(-10);
k89=1.10*10^(-9)*(T/300)^(-0.50);
k90=4.51*10^(-10);
k91=1* 10^(-11);
k92=8.40*10^(-10);
k93=8.50^(-10);
k94=3.40*10^(-10);
k95=5*10^(-11);
k96=1.20*10^(-10);
k97=6.60*10^(-10);
k98=4.60*10^(-10);
k99=1.00*10^(-10)*(T/300)^(-0.5);
k100=7.40*10^(-11);
k101=6.30*10^(-10);
k102=3.30*10^(-10);
k103=5.70*10^(-10);
%constante k::Neutralisation
k104=7.51*10^(-8)*(T/300)^(-0.5);
k105=7.51*10^(-8)*(T/300)^(-0.5);
k106=7.51*10^(-8)*(T/300)^(-0.5);
k107=7.51*10^(-8)*(T/300)^(-0.5);
k108=7.51*10^(-8)*(T/300)^(-0.5);
k109=7.51*10^(-8)*(T/300)^(-0.5);
k110=7.51*10^(-8)*(T/300)^(-0.50);
k111=7.51*10^(-8)*(T/300)^(-0.50);
k112=7.51*10^(-8)*(T/300)^(-0.50);
k113=7.51*10^(-8)*(T/300)^(-0.50);
k114=7.51*10^(-8)*(T/300)^(-0.5);
k115=7.51*10^(-8)*(T/300)^(-0.5);
k116=7.51*10^(-8)*(T/300)^(-0.5);
k117=7.51*10^(-8)*(T/300)^(-0.5);
%constante k::Photoionisation
k118=3.10*10^(-10)*exp(-3.3/T);
k119=4.90*10^(-8)*exp(-0.50/T);
k120=1.09*10^(-8)*exp(-0.5/T);
k121=4.10*10^(-10)*exp(-3.80/T);
k122=1*10^(-11)*exp(-1.70/T);
k123=2.96*10^(-9)*exp(-2/T);
k124=5*10^(-9)*exp(-2.1/T);
k125=6.88*10^(-9)*exp(-1.5/T);
k126=7.90*10^(-10)*exp(-2.1/T);
k127=6.10*10^(-9)*exp(-0.5/T);
k128=3.50*10^(-11)*exp(-2/T);
k129=2.30*10^(-10)*exp(-3.90/T);
k130=2.10*10^(-10)*exp(-3.50/T);
k131=4.70*10^(-10)*exp(-2.10/T);
%constante k::Attachement radiatif
k132=2.25*10^(-15);
k133=1.50*10^(-15);
k134=2.00*10^(-15)*(T/300)^(-0.5);
k135=1.70*10^(-14)*(T/300)^(-0.5);
%constante k::Détachement associatif
k136=5*10^(-10);
k137=1*10^(-9);
k138=1*10^(-9);
k139=5*10^(-10);
k140=5*10^(-10);
k141=5*10^(-11);
k142=5*10^(-10);
k143=1.90*10^(-10);
k144=2.90*10^(-10);
k145=3.10*10^(-10)*(T/300)^(-0.83);
k146=6.50*10^(-10);
k147=2.20*10^(-10);
k148=1.00*10^(-9);
k149=1*10^(-9);
%constante k::Excitation/ Désexcitation
k150=2.0425*10^(-7)*T^(-0.3)*exp(-29489.55/T);
k151=2.6*10^(-10)*T^(0.89)*exp(-11372.9/T);
k152=1.2*10^(-10)*T^(0.59)*exp(-18916.15/T);
k153=6.0278*10^(-7)*T^(-0.4)*exp(-19880.36/T);
k154=1.6306*10^(-17);
k155=3.6532*10^(-11);
%constante k::Détachement électronique
k156=7*10^(-10);
k157=2*10^(-10);
k158=3.6*10^(-10);
%constante k::Dissociation par impact électronique
k159= 6.2*10^(-10) * T^(1.75) * exp(-113032.7/T );
%constante k::Ionisation par impact électronique
k160= 7.1*10^(-9) * T^(0.76) * exp(-181038/T );
%constante k::Attachement électronique cn(1)6/s
k161=1*10^(-31);
k162=1*10^(-31);
%constante k::Processus de recombinaison cn(1)6/s
k163=5.70*10^(-34)*T^(-2.6);
k164=2.5919*10^(-80);
k165=1.31*10^(-31)*T^(-1.5);
function dndt=myequations(t,n,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,k18,k19,k20,k21,k22,k23,k24,k25,k26,k27,k28,k29,k30,k31,k32,k33,k34,k35,k36,k37,k38,k39,k40,k41,k42,k43,k44,k45,k46,k47,k48,k49,k50,k51,k52,k53,k54,k55,k56,k57,k58,k59,k60,k61,k62,k63,k64,k65,k66,k67,k68,k69,k70,k71,k72,k73,k74,k75,k76,k77,k78,k79,k80,k81,k82,k83,k84,k85,k86,k87,k88,k89,k90,k91,k92,k93,k94,k95,k96,k97,k98,k99,k100,k101,k102,k103,k104,k105,k106,k107,k108,k109,k110,k111,k112,k113,k114,k115,k116,k117,k118,k119,k120,k121,k122,k123,k124,k125,k126,k127,k128,k129,k130,k131,k132,k133,k134,k135,k136,k137,k138,k139,k140,k141,k142,k143,k144,k145,k146,k147,k148,k149,k150,k151,k152,k153,k154,k155,k156,k157,k158,k159,k160,k161,k162,k163,k164,k165)
dndt=zeros(43,1)
%Chemin C
dndt(1)= k1*n(2)*n(42)+k15*n(5)*n(9)+k16*n(5)*n(26)+k25*n(4)*n(26)...
+k29*n(4)*n(9)+k33*n(4)*n(11)+k39*n(7)*n(26)+k40*n(7)*n(27)...
+k43*n(7)*n(9)+k44*n(7)*n(12)+k48*n(9)*n(9)+k59*n(11)*n(42)...
+k72*n(2)*n(24)+k104*n(3)*n(5)+k105*n(3)*n(8)+k106*n(3)*n(25)...
+k117*n(18)*n(2)+k119*n(3)*n(43)+k122*n(11)*n(43)+k124*n(13)*n(43)...
+k130*n(22)*n(43)-k4*n(1)*n(37)-k8*n(1)*n(15)-k9*n(1)*n(26)...
-k10*n(1)*n(24)-k11*n(1)*n(20)-k12*n(1)*n(35)-k13*n(1)*n(22)...
-k61*n(2)*n(1)-k64*n(1)*n(13)-k65*n(1)*n(4)-k66*n(1)*n(9)...
-k67*n(1)*n(5)-k68*n(1)*n(1)-k69*n(1)*n(7)-k73*n(1)*n(27)...
-k74*n(1)*n(11)-k75*n(1)*n(23)-k76*n(1)*n(21)-k77*n(1)*n(19)...
-k118*n(1)*n(43)-k132*n(1)*n(42)-k136*n(1)*n(3)-k137*n(1)*n(12)...
-k138*n(1)*n(14)-k139*n(1)*n(6)-2*k164*n(7)*n(1)*n(1)...
+k164*n(7)*n(1)^(2)-k163*n(1)*n(15)*n(4)+k163*n(1)*n(15)*n(4)...
-k165*n(1)*n(32)*n(4)+k165*n(1)*n(32)*n(4);
%Chemin C+:
dndt(2)= -k1*n(2)*n(42)-k5*n(2)*n(15)-k6*n(2)*n(30)-k7*n(2)*n(36)...
-k61*n(2)*n(1)-k62*n(2)*n(7)-k63*n(2)*n(4)-k72*n(2)*n(24)...
-k117*n(18)*n(2)+k45*n(7)*n(11)+k73*n(1)*n(27)...
+k74*n(1)*n(11)+k75*n(1)*n(23)+k76*n(1)*n(21)...
+k77*n(1)*n(19)+k118*n(1)*n(43)+k122*n(11)*n(43);
%Chemin C-:
dndt(3)= k32*n(4)*n(12)+k132*n(1)*n(42)-k14*n(3)*n(24)-k104*n(3)*n(5)...
-k105*n(3)*n(8)-k106*n(3)*n(25)-k119*n(3)*n(43)...
-k136*n(1)*n(3)-k140*n(3)*n(4)-k141*n(3)*n(15)-k142*n(3)*n(7);
%Chemin O:
dndt(4)= k2*n(5)*n(42)+k8*n(1)*n(15)+k10*n(1)*n(24)+k13*n(1)*n(22)...
+k14*n(3)*n(24)+k35*n(7)*n(15)+k36*n(7)*n(19)...
+k37*n(7)*n(24)+k38*n(7)*n(32)+k52*n(15)*n(24)...
+k53*n(15)*n(20)+k55*n(15)*n(22)+k56*n(15)*n(26)...
+k78*n(5)*n(9)+k79*n(5)*n(22)+k80*n(5)*n(30)...
+k81*n(5)*n(15)+k84*n(6)*n(15)+k85*n(6)*n(26)...
+k104*n(3)*n(5)+2*k107*n(5)*n(6)+k108*n(5)*n(18)...
+k109*n(6)*n(8)+k120*n(6)*n(43)+2*k126*n(15)*n(43)...
+k128*n(19)*n(43)+k130*n(22)*n(43)+k131*n(24)*n(43)...
+k162*n(42)*n(4)*n(15)-k19*n(4)*n(24)-k20*n(4)*n(32)...
-k21*n(4)*n(30)-k22*n(4)*n(30)-k23*n(4)*n(20)...
-k24*n(4)*n(36)-k25*n(4)*n(26)-k26*n(4)*n(26)...
-k27*n(4)*n(13)-k28*n(4)*n(35)-k29*n(4)*n(9)...
-k30*n(4)*n(21)-k31*n(4)*n(14)-k32*n(4)*n(12)...
-k33*n(4)*n(11)-k63*n(2)*n(4)-k65*n(1)*n(4)...
-2*k70*n(4)*n(4)-k82*n(4)*n(21)-k83*n(4)*n(23)...
-k133*n(4)*n(42)-k140*n(3)*n(4)-k143*n(4)*n(6)...
-k144*n(4)*n(12)-k161*n(42)*n(4)*n(15)...
-k162*n(42)*n(4)*n(15)-k163*n(1)*n(15)*n(4)-k165*n(1)*n(32)*n(4);
%Chemin O+:
dndt(5)= k5*n(2)*n(15)+k34*n(8)*n(15)+k82*n(4)*n(21)...
+k83*n(4)*n(23)+k128*n(19)*n(43)-k19*n(4)*n(24)...
-k20*n(4)*n(32)-k21*n(4)*n(30)-k22*n(4)*n(30)...
-k23*n(4)*n(20)-k24*n(4)*n(36)-k25*n(4)*n(26)...
-k26*n(4)*n(26)-k27*n(4)*n(13)-k28*n(4)*n(35)...
-k29*n(4)*n(9)-k30*n(4)*n(21)-k31*n(4)*n(14)...
-k32*n(4)*n(12)-k33*n(4)*n(11)-k63*n(2)*n(4)...
-k65*n(1)*n(4)-2*k70*n(4)*n(4)-k82*n(4)*n(21)...
-k83*n(4)*n(23)-k133*n(4)*n(42)-k140*n(3)*n(4)...
-k143*n(4)*n(6)-k144*n(4)*n(12)-k161*n(42)*n(4)*n(15)...
-k162*n(42)*n(4)*n(15)-k163*n(1)*n(15)*n(4)-k165*n(1)*n(32)*n(4);
%Chemin O-:
dndt(6)= k133*n(4)*n(42)+k161*n(42)*n(4)*n(15)-k84*n(6)*n(15)...
-k85*n(6)*n(26)-k107*n(5)*n(6)-k109*n(6)*n(8)-k120*n(6)*n(43)...
-k139*n(1)*n(6)-k143*n(4)*n(6)-k145*n(6)*n(24)...
-k146*n(6)*n(22)-k147*n(6)*n(7);
%Chemin N:
dndt(7)= k3*n(8)*n(42)+k9*n(1)*n(26)+k11*n(1)*n(20)...
+k18*n(5)*n(20)+k19*n(4)*n(24)+k23*n(4)*n(20)...
+k26*n(4)*n(26)+k30*n(4)*n(21)+2*k60*n(21)*n(42)...
+k86*n(8)*n(15)+k87*n(8)*n(9)+k88*n(8)*n(22)...
+k89*n(8)*n(26)+k90*n(8)*n(24)+k105*n(3)*n(8)...
+k109*n(6)*n(8)+k110*n(8)*n(18)+k111*n(8)*n(12)...
+k112*n(8)*n(28)+k113*n(8)*n(14)+2*k129*n(20)*n(43)...
+k131*n(24)*n(43)+2*k159*n(20)*n(42)-k35*n(7)*n(15)...
-k36*n(7)*n(19)-k37*n(7)*n(24)-k38*n(7)*n(32)...
-k39*n(7)*n(26)-k40*n(7)*n(27)-k41*n(7)*n(13)...
-k42*n(7)*n(14)-k43*n(7)*n(9)-k44*n(7)*n(12)...
-k45*n(7)*n(11)-k46*n(7)*n(35)-k62*n(2)*n(7)...
-k69*n(1)*n(7)-k71*n(8)*n(7)-k91*n(7)*n(21)...
-k142*n(3)*n(7)-k147*n(6)*n(7)-k164*n(7)*n(1)*n(1);
%Chemin N+:
dndt(8)=k91*n(7)*n(21)-k3*n(8)*n(42)-k34*n(8)*n(15)...
-k71*n(8)*n(7)-k86*n(8)*n(15)-k87*n(8)*n(9)...
-k88*n(8)*n(22)-k89*n(8)*n(26)-k90*n(8)*n(24)...
-k105*n(3)*n(8)-k109*n(6)*n(8)-k110*n(8)*n(18)...
-k111*n(8)*n(12)-k112*n(8)*n(28)-k113*n(8)*n(14);
%Chemin C2:
dndt(9)= k9*n(1)*n(26)+k12*n(1)*n(35)+k13*n(1)*n(22)...
+k27*n(4)*n(13)+k41*n(7)*n(13)+k68*n(1)*n(1)...
+k74*n(1)*n(11)+k94*n(11)*n(24)+k111*n(8)*n(12)...
+k115*n(12)*n(25)+k121*n(9)*n(43)+k123*n(12)*n(43)...
+k124*n(13)*n(43)+k136*n(1)*n(3)+k153*n(10)*n(42)...
+k154*n(10)^(2)+k155*n(10)*n(1)-k15*n(5)*n(9)...
-k29*n(4)*n(9)-k43*n(7)*n(9)-k47*n(9)*n(15)...
-2*k48*n(9)*n(9)-k49*n(9)*n(19)-k66*n(1)*n(9)...
-k78*n(5)*n(9)-k87*n(8)*n(9)-k92*n(9)*n(23)...
-k93*n(9)*n(27)-k121*n(9)*n(43)-k134*n(9)*n(42)...
-k148*n(9)*n(12)-k150*n(9)*n(42);
%Chemin C2(d)
dndt(10)= k150*n(9)*n(42)-k153*n(10)*n(42)-k154*n(10)^(2)-k155*n(10)*n(1);
%Chemin C2+:
dndt(11)= k61*n(2)*n(1)+k78*n(5)*n(9)+k87*n(8)*n(9)...
+k92*n(9)*n(23)+k93*n(9)*n(27)-k33*n(4)*n(11)...
-k45*n(7)*n(11)-k50*n(11)*n(15)-k59*n(11)*n(42)...
-k74*n(1)*n(11)-k94*n(11)*n(24)-k122*n(11)*n(43);
%Chemin C2-:
dndt(12)= -k32*n(4)*n(12)-k44*n(7)*n(12)-k111*n(8)*n(12)...
-k115*n(12)*n(25)-k123*n(12)*n(43)-k137*n(1)*n(12)...
-k144*n(4)*n(12)-k148*n(9)*n(12)-k149*n(12)*n(13)...
+k31*n(4)*n(14)+k42*n(7)*n(14)+k134*n(9)*n(42);
%Chemin C3:
dndt(13)= k48*n(9)*n(9)+k66*n(1)*n(9)+k113*n(8)*n(14)...
+k114*n(14)*n(25)+k125*n(14)*n(43)...
+k137*n(1)*n(12)-k27*n(4)*n(13)-k41*n(7)*n(13)...
-k64*n(1)*n(13)-k124*n(13)*n(43)...
-k135*n(13)*n(42)-k149*n(12)*n(13);
%Chemin C3-:
dndt(14)= k135*n(13)*n(42)-k31*n(4)*n(14)...
-k42*n(7)*n(14)-k113*n(8)*n(14)...
-k114*n(14)*n(25)-k125*n(14)*n(43)-k138*n(1)*n(14);
%Chemin O2:
dndt(15)= k17*n(5)*n(32)+k19*n(4)*n(24)+k20*n(4)*n(32)...
+k21*n(4)*n(30)+k57*n(24)*n(24)+k70*n(4)*n(4)...
+k77*n(1)*n(19)+k97*n(19)*n(32)+k98*n(19)*n(24)...
+k108*n(5)*n(18)+k110*n(8)*n(18)+k116*n(18)*n(25)...
+k117*n(18)*n(2)+k127*n(18)*n(43)+k143*n(4)*n(6)...
+2*k156*n(18)*n(15)+2*k157*n(18)*n(16)...
+2*k158*n(18)*n(17)+k161*n(42)*n(4)*n(15)...
-k5*n(2)*n(15)-k8*n(1)*n(15)-k34*n(8)*n(15)...
-k35*n(7)*n(15)-k47*n(9)*n(15)-k50*n(11)*n(15)...
-k51*n(15)*n(37)-k52*n(15)*n(24)-k53*n(15)*n(20)...
-k54*n(15)*n(27)-k55*n(15)*n(22)-k56*n(15)*n(26)...
-k81*n(5)*n(15)-k84*n(6)*n(15)-k86*n(8)*n(15)...
-k95*n(15)*n(21)-k96*n(15)*n(23)-k126*n(15)*n(43)...
-k141*n(3)*n(15)-k151*n(15)*n(42)-k152*n(15)*n(42)...
-k156*n(18)*n(15)-k161*n(42)*n(4)*n(15)...
-k162*n(42)*n(4)*n(15)-k163*n(1)*n(15)*n(4);
%Chemin O2(a):
dndt(16)= k151*n(42)*n(15)-k157*n(16)*n(18);
%Chemin O2(b):
dndt(17)= k152*n(42)*n(15)-k158*n(17)*n(18);
%Chemin O2-:
dndt(18)= k84*n(6)*n(15)+k162*n(42)*n(4)*n(15)...
-k108*n(5)*n(18)-k110*n(8)*n(18)...
-k116*n(18)*n(25)-k117*n(18)*n(2)-k127*n(18)*n(43)...
-k156*n(18)*n(15)-k157*n(18)*n(16)-k158*n(18)*n(17);
%Chemin O2+:
dndt(19)= k81*n(5)*n(15)+k86*n(8)*n(15)+k95*n(15)*n(21)...
+k96*n(15)*n(23)-k36*n(7)*n(19)-k49*n(9)*n(19)-k77*n(1)*n(19)...
-k97*n(19)*n(32)-k98*n(19)*n(24)-k128*n(19)*n(43);
%Chemin N2:
dndt(20)= k21*n(4)*n(30)+k37*n(7)*n(24)+k39*n(7)*n(26)...
+k57*n(24)*n(24)+k76*n(1)*n(21)+k82*n(4)*n(21)...
+k91*n(7)*n(21)+k95*n(15)*n(21)+k99*n(21)*n(26)...
+k100*n(21)*n(22)-k11*n(1)*n(20)-k18*n(5)*n(20)...
-k23*n(4)*n(20)-k53*n(15)*n(20)...
-k129*n(20)*n(43)-k159*n(20)*n(42)-k160*n(20)*n(42);
%Chemin N2+:
dndt(21)= k40*n(7)*n(27)+k71*n(8)*n(7)+k160*n(20)*n(42)...
-k30*n(4)*n(21)-k60*n(21)*n(42)-k76*n(1)*n(21)...
-k82*n(4)*n(21)-k91*n(7)*n(21)...
-k95*n(15)*n(21)-k99*n(21)*n(26)-k100*n(21)*n(22);
%Chemin CO:
dndt(22)= k4*n(1)*n(37)+k5*n(2)*n(15)+k8*n(1)*n(15)...
+k24*n(4)*n(36)+k26*n(4)*n(26)+k27*n(4)*n(13)...
+k28*n(4)*n(35)+k29*n(4)*n(9)+k31*n(4)*n(14)...
+k32*n(4)*n(12)+2*k47*n(9)*n(15)+k49*n(9)*n(19)...
+k50*n(11)*n(15)+k54*n(15)*n(27)+k58*n(24)*n(37)...
+k65*n(1)*n(4)+k75*n(1)*n(23)+k83*n(4)*n(23)...
+k92*n(9)*n(23)+k96*n(15)*n(23)+k102*n(24)*n(23)...
+k139*n(1)*n(6)+k140*n(3)*n(4)-k13*n(1)*n(22)...
-k55*n(15)*n(22)-k79*n(5)*n(22)-k88*n(8)*n(22)...
-k100*n(21)*n(22)-k101*n(22)*n(27)...
-k130*n(22)*n(43)-k146*n(6)*n(22);
%Chemin CO+:
dndt(23)= k7*n(2)*n(36)+k15*n(5)*n(9)+k33*n(4)*n(11)...
+k49*n(9)*n(19)+k50*n(11)*n(15)+k63*n(2)*n(4)...
+k67*n(1)*n(5)+k79*n(5)*n(22)+k88*n(8)*n(22)...
+k100*n(21)*n(22)+k101*n(22)*n(27)-k75*n(1)*n(23)...
-k83*n(4)*n(23)-k92*n(9)*n(23)...
-k96*n(15)*n(23)-k102*n(24)*n(23);
%Chemin NO:
dndt(24)= k20*n(4)*n(32)+k22*2*n(4)*n(30)+k23*n(4)*n(20)...
+k24*n(4)*n(36)+k25*n(4)*n(26)+k34*n(8)*n(15)...
+k35*n(7)*n(15)+k51*n(15)*n(37)+k106*n(3)*n(25)...
+k114*n(14)*n(25)+k115*n(12)*n(25)+k116*n(18)*n(25)...
+k147*n(6)*n(7)-k10*n(1)*n(24)-k14*n(3)*n(24)-k19*n(4)*n(24)...
-k37*n(7)*n(24)-k52*n(15)*n(24)-k57*2*n(24)*n(24)...
-k58*n(24)*n(37)-k72*n(2)*n(24)-k90*n(8)*n(24)...
-k94*n(11)*n(24)-k98*n(19)*n(24)-k102*n(24)*n(23)...
-k103*n(24)*n(27)-k131*n(24)*n(43)-k145*n(6)*n(24);
%Chemin NO+:
dndt(25)= k6*n(2)*n(30)+k16*n(5)*n(26)+k17*n(5)*n(32)...
+k18*n(5)*n(20)+k30*n(4)*n(21)+k36*n(7)*n(19)...
+k54*n(15)*n(27)+k72*n(2)*n(24)+k90*n(8)*n(24)...
+k94*n(11)*n(24)+k98*n(19)*n(24)+k102*n(24)*n(23)...
+k103*n(24)*n(27)-k106*n(3)*n(25)-k114*n(14)*n(25)...
-k115*n(12)*n(25)-k116*n(18)*n(25);
%Chemin CN:
dndt(26)= k4*n(1)*n(37)+k6*n(2)*n(30)+k7*n(2)*n(36)...
+k10*n(1)*n(24)+k11*n(1)*n(20)+k12*n(1)*n(35)...
+k28*n(4)*n(35)+k41*n(7)*n(13)+k42*n(7)*n(14)...
+k43*n(7)*n(9)+k45*n(7)*n(11)+k46*2*n(7)*n(35)...
+k69*n(1)*n(7)+k73*n(1)*n(27)+k93*n(9)*n(27)...
+k101*n(22)*n(27)+k103*n(24)*n(27)...
+k112*n(8)*n(28)+k142*n(3)*n(7)...
+k164*n(7)*n(1)*n(1)-k9*n(1)*n(26)-k16*n(5)*n(26)...
-k25*n(4)*n(26)-k26*n(4)*n(26)-k39*n(7)*n(26)...
-k56*n(15)*n(26)-k85*n(6)*n(26)...
-k89*n(8)*n(26)-k99*n(21)*n(26);
%Chemin CN+:
dndt(27)= k62*n(2)*n(7)+k89*n(8)*n(26)+k99*n(21)*n(26)...
-k40*n(7)*n(27)-k54*n(15)*n(27)...
-k73*n(1)*n(27)-k93*n(9)*n(27)...
-k101*n(22)*n(27)-k103*n(24)*n(27);
%Chemin CN-:
dndt(28)= k14*n(3)*n(24)+k44*n(7)*n(12)+k85*n(6)*n(26)-k112*n(8)*n(28);
%Chemin O3:
dndt(29)= k163*n(1)*n(15)*n(4);
%Chemin N2O:
dndt(30)= k38*n(7)*n(32)+k53*n(15)*n(20)+k58*n(24)*n(37)...
-k6*n(2)*n(30)-k21*n(4)*n(30)...
-k22*n(4)*n(30)-k80*n(5)*n(30);
%Chemin N2O+:
dndt(31)= k80*n(5)*n(30);
%Chemin NO2:
dndt(32)= k52*n(15)*n(24)+k145*n(6)*n(24)-k17*n(5)*n(32)...
-k20*n(4)*n(32)-k38*n(7)*n(32)...
-k97*n(19)*n(32)-k165*n(1)*n(32)*n(4);
%Chemin NO2+:
dndt(33)= k97*n(19)*n(32);
%Chemin CO2:
dndt(34)= k51*n(15)*n(37)+k55*n(15)*n(22)...
+k141*n(3)*n(15)+k146*n(6)*n(22);
%Chemin C2N:
dndt(35)= -k12*n(1)*n(35)-k28*n(4)*n(35)-k46*n(7)*n(35);
%Chemin CNO:
dndt(36)= -k7*n(2)*n(36)-k24*n(4)*n(36);
%Chemin OCN
dndt(37)= -k4*n(1)*n(37)-k51*n(15)*n(37)...
-k58*n(24)*n(37)+k56*n(15)*n(26);
%Chemin NO3:
dndt(38)= k165*n(1)*n(32)*n(4);
%Chemin C2O:
dndt(39)= k144*n(4)*n(12);
%Chemin C4:
dndt(40)= k64*n(1)*n(13)+k138*n(1)*n(14)+k148*n(9)*n(12);
%Chemin C5:
dndt(41)= k149*n(12)*n(13);
%Chemin e:
dndt(42)= k118*n(1)*n(43)+k119*n(3)*n(43)+k120*n(6)*n(43)...
+k121*n(12)*n(43)+k123*n(12)*n(43)+k125*n(14)*n(43)...
+k127*n(18)*n(43)+k136*n(1)*n(3)+k137*n(1)*n(12)...
+k138*n(1)*n(14)+k139*n(1)*n(6)+k140*n(3)*n(4)...
+k141*n(3)*n(15)+k142*n(3)*n(7)+k143*n(4)*n(6)...
+k144*n(4)*n(12)+k145*n(6)*n(24)+k146*n(6)*n(22)...
+k147*n(6)*n(7)+k148*n(9)*n(12)+k149*n(12)*n(13)...
+k150*n(9)*n(42)+k151*n(15)*n(42)+k152*n(15)*n(42)...
+k153*n(10)*n(42)+k156*n(18)*n(15)+k157*n(18)*n(16)...
+k158*n(18)*n(17)+k159*n(20)*n(42)...
+2*k160*n(20)*n(42)-k1*n(2)*n(42)-k2*n(5)*n(42)...
-k3*n(8)*n(42)-k59*n(11)*n(42)-k60*n(21)*n(42)...
-k132*n(1)*n(42)-k133*n(4)*n(42)-k134*n(9)*n(42)...
-k135*n(13)*n(42)-k150*n(9)*n(42)-k151*n(15)*n(42)...
-k152*n(15)*n(42)-k153*n(10)*n(42)-k159*n(20)*n(42)...
-k160*n(20)*n(42)-k161*n(42)*n(4)*n(15)...
-k162*n(42)*n(4)*n(15);
%Chemin hv:
dndt(43)= -k118*n(1)*n(43)-k119*n(3)*n(43)-k120*n(6)*n(43)...
-k121*n(9)*n(43)-k122*n(11)*n(43)-k123*n(12)*n(43)...
-k124*n(13)*n(43)-k125*n(14)*n(43)-k126*n(15)*n(43)...
-k127*n(18)*n(43)-k128*n(19)*n(43)-k129*n(20)*n(43)...
-k130*n(22)*n(43)-k131*n(24)*n(43)+k1*n(2)*n(42)...
+k2*n(5)*n(42)+k3*n(8)*n(42)+k61*n(2)*n(1)...
+k62*n(2)*n(7)+k63*n(2)*n(4)+k64*n(1)*n(13)...
+k65*n(1)*n(4)+k66*n(1)*n(9)+k67*n(1)*n(5)...
+k68*n(1)*n(1)+k69*n(1)*n(7)+k70*n(4)*n(4)...
+k71*n(8)*n(7)+k132*n(1)*n(42)+k133*n(4)*n(42)...
+k134*n(9)*n(42)+k135*n(13)*n(42);
end
Here is my second script:
function [t,n]= call_myequations()
%condition intial
%IC=IC1=IC2=....=IC43=10^6cm^3
IC =10^19*ones(1,43);
tic,[t,n] = ode45(@(t,n)myequations(t,n,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,k18,k19,k20,k21,k22,k23,k24,k25,k26,k27,k28,k29,k30,k31,k32,k33,k34,k35,k36,k37,k38,k39,k40,k41,k42,k43,k44,k45,k46,k47,k48,k49,k50,k51,k52,k53,k54,k55,k56,k57,k58,k59,k60,k61,k62,k63,k64,k65,k66,k67,k68,k69,k70,k71,k72,k73,k74,k75,k76,k77,k78,k79,k80,k81,k82,k83,k84,k85,k86,k87,k88,k89,k90,k91,k92,k93,k94,k95,k96,k97,k98,k99,k100,k101,k102,k103,k104,k105,k106,k107,k108,k109,k110,k111,k112,k113,k114,k115,k116,k117,k118,k119,k120,k121,k122,k123,k124,k125,k126,k127,k128,k129,k130,k131,k132,k133,k134,k135,k136,k137,k138,k139,k140,k141,k142,k143,k144,k145,k146,k147,k148,k149,k150,k151,k152,k153,k154,k155,k156,k157,k158,k159,k160,k161,k162,k163,k164,k165),[0 10e-3],IC)
toc
plot(t,n)
end
call_myequations
Unrecognized function or variable 'k1'.
Error in
call_myequations>@(t,n)myequations(t,n,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,k18,k19,k20,k21,k22,k23,k24,k25,k26,k27,k28,k29,k30,k31,k32,k33,k34,k35,k36,k37,k38,k39,k40,k41,k42,k43,k44,k45,k46,k47,k48,k49,k50,k51,k52,k53,k54,k55,k56,k57,k58,k59,k60,k61,k62,k63,k64,k65,k66,k67,k68,k69,k70,k71,k72,k73,k74,k75,k76,k77,k78,k79,k80,k81,k82,k83,k84,k85,k86,k87,k88,k89,k90,k91,k92,k93,k94,k95,k96,k97,k98,k99,k100,k101,k102,k103,k104,k105,k106,k107,k108,k109,k110,k111,k112,k113,k114,k115,k116,k117,k118,k119,k120,k121,k122,k123,k124,k125,k126,k127,k128,k129,k130,k131,k132,k133,k134,k135,k136,k137,k138,k139,k140,k141,k142,k143,k144,k145,k146,k147,k148,k149,k150,k151,k152,k153,k154,k155,k156,k157,k158,k159,k160,k161,k162,k163,k164,k165)
(line 6)
tic,[t,n] =
ode45(@(t,n)myequations(t,n,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,k18,k19,k20,k21,k22,k23,k24,k25,k26,k27,k28,k29,k30,k31,k32,k33,k34,k35,k36,k37,k38,k39,k40,k41,k42,k43,k44,k45,k46,k47,k48,k49,k50,k51,k52,k53,k54,k55,k56,k57,k58,k59,k60,k61,k62,k63,k64,k65,k66,k67,k68,k69,k70,k71,k72,k73,k74,k75,k76,k77,k78,k79,k80,k81,k82,k83,k84,k85,k86,k87,k88,k89,k90,k91,k92,k93,k94,k95,k96,k97,k98,k99,k100,k101,k102,k103,k104,k105,k106,k107,k108,k109,k110,k111,k112,k113,k114,k115,k116,k117,k118,k119,k120,k121,k122,k123,k124,k125,k126,k127,k128,k129,k130,k131,k132,k133,k134,k135,k136,k137,k138,k139,k140,k141,k142,k143,k144,k145,k146,k147,k148,k149,k150,k151,k152,k153,k154,k155,k156,k157,k158,k159,k160,k161,k162,k163,k164,k165),[0
10e-3],IC)
Error in odearguments (line 90)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 115)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
Error in call_myequations (line 6)
tic,[t,n] =
ode45(@(t,n)myequations(t,n,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,k18,k19,k20,k21,k22,k23,k24,k25,k26,k27,k28,k29,k30,k31,k32,k33,k34,k35,k36,k37,k38,k39,k40,k41,k42,k43,k44,k45,k46,k47,k48,k49,k50,k51,k52,k53,k54,k55,k56,k57,k58,k59,k60,k61,k62,k63,k64,k65,k66,k67,k68,k69,k70,k71,k72,k73,k74,k75,k76,k77,k78,k79,k80,k81,k82,k83,k84,k85,k86,k87,k88,k89,k90,k91,k92,k93,k94,k95,k96,k97,k98,k99,k100,k101,k102,k103,k104,k105,k106,k107,k108,k109,k110,k111,k112,k113,k114,k115,k116,k117,k118,k119,k120,k121,k122,k123,k124,k125,k126,k127,k128,k129,k130,k131,k132,k133,k134,k135,k136,k137,k138,k139,k140,k141,k142,k143,k144,k145,k146,k147,k148,k149,k150,k151,k152,k153,k154,k155,k156,k157,k158,k159,k160,k161,k162,k163,k164,k165),[0
10e-3],IC)
>> call_myequations
Unrecognized function or variable 'k1'.
Error in
call_myequations>@(t,n)myequations(t,n,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,k18,k19,k20,k21,k22,k23,k24,k25,k26,k27,k28,k29,k30,k31,k32,k33,k34,k35,k36,k37,k38,k39,k40,k41,k42,k43,k44,k45,k46,k47,k48,k49,k50,k51,k52,k53,k54,k55,k56,k57,k58,k59,k60,k61,k62,k63,k64,k65,k66,k67,k68,k69,k70,k71,k72,k73,k74,k75,k76,k77,k78,k79,k80,k81,k82,k83,k84,k85,k86,k87,k88,k89,k90,k91,k92,k93,k94,k95,k96,k97,k98,k99,k100,k101,k102,k103,k104,k105,k106,k107,k108,k109,k110,k111,k112,k113,k114,k115,k116,k117,k118,k119,k120,k121,k122,k123,k124,k125,k126,k127,k128,k129,k130,k131,k132,k133,k134,k135,k136,k137,k138,k139,k140,k141,k142,k143,k144,k145,k146,k147,k148,k149,k150,k151,k152,k153,k154,k155,k156,k157,k158,k159,k160,k161,k162,k163,k164,k165)
(line 6)
tic,[t,n] =
ode45(@(t,n)myequations(t,n,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,k18,k19,k20,k21,k22,k23,k24,k25,k26,k27,k28,k29,k30,k31,k32,k33,k34,k35,k36,k37,k38,k39,k40,k41,k42,k43,k44,k45,k46,k47,k48,k49,k50,k51,k52,k53,k54,k55,k56,k57,k58,k59,k60,k61,k62,k63,k64,k65,k66,k67,k68,k69,k70,k71,k72,k73,k74,k75,k76,k77,k78,k79,k80,k81,k82,k83,k84,k85,k86,k87,k88,k89,k90,k91,k92,k93,k94,k95,k96,k97,k98,k99,k100,k101,k102,k103,k104,k105,k106,k107,k108,k109,k110,k111,k112,k113,k114,k115,k116,k117,k118,k119,k120,k121,k122,k123,k124,k125,k126,k127,k128,k129,k130,k131,k132,k133,k134,k135,k136,k137,k138,k139,k140,k141,k142,k143,k144,k145,k146,k147,k148,k149,k150,k151,k152,k153,k154,k155,k156,k157,k158,k159,k160,k161,k162,k163,k164,k165),[0
10e-3],IC)
Error in odearguments (line 90)
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ode45 (line 115)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
Error in call_myequations (line 6)
tic,[t,n] =
ode45(@(t,n)myequations(t,n,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,k18,k19,k20,k21,k22,k23,k24,k25,k26,k27,k28,k29,k30,k31,k32,k33,k34,k35,k36,k37,k38,k39,k40,k41,k42,k43,k44,k45,k46,k47,k48,k49,k50,k51,k52,k53,k54,k55,k56,k57,k58,k59,k60,k61,k62,k63,k64,k65,k66,k67,k68,k69,k70,k71,k72,k73,k74,k75,k76,k77,k78,k79,k80,k81,k82,k83,k84,k85,k86,k87,k88,k89,k90,k91,k92,k93,k94,k95,k96,k97,k98,k99,k100,k101,k102,k103,k104,k105,k106,k107,k108,k109,k110,k111,k112,k113,k114,k115,k116,k117,k118,k119,k120,k121,k122,k123,k124,k125,k126,k127,k128,k129,k130,k131,k132,k133,k134,k135,k136,k137,k138,k139,k140,k141,k142,k143,k144,k145,k146,k147,k148,k149,k150,k151,k152,k153,k154,k155,k156,k157,k158,k159,k160,k161,k162,k163,k164,k165),[0
10e-3],IC)
>> call_myequations
Error: File: call_myequations.m Line: 10 Column: 1
All functions in a script must be closed with an 'end'.
>> call_myequations
Error: File: call_myequations.m Line: 2 Column: 17
Local function name must be different from the script name.
답변 (2개)
[t,n] = call_myequations();
function [t,n]= call_myequations()
%condition intial
%IC=IC1=IC2=....=IC43=10^6cm^3
IC =10^19*ones(1,43);
obj = init_fun;
tic,[t,n] = ode45( obj, [0 10e-3],IC);
toc
plot(t,n)
end
function fun = init_fun
%constantes :::
%temperature
T=17000;
%constante k: recombinaison radiative
k1=2.36 *(1e-12)*(T/300)^(-0.29)*exp(17.60/T);
k2 =3.24 *(1e-12) *(T/300)^(-0.66);
k3=3.50 *(1e-12)*(T/300)^(-0.53)* exp(3.20/T);
%constante k::Reconmbinaison
k4=1*10^(-10);
k5=4.54*10^(-10);
k6=9.10*10^(-10);
k7=8.98*10^(-9)*(T/300)^(-0.5);
k8=5.56*10^(-11)*(T/300)^(0.41)*exp(26.90/T);
k9=4.98*10^(-10)*exp(-18116/T);
k10=6*10^(-11)*(T/300)^(-0.16);
k11=8.69*10^(-11)*exp(-22600/T);
k12=1*10^(-10);
k13=2.94*10^(-11)*(T/300)^(0.5)*exp(-58025/T);
k14=1*10^(-9);
k15=4.80*10^(-10);
k16=1.00*10^(-9)*(T/300)^(-0.5);
k17=8.30*10^(-10);
k18=2.42*10^(-12)*(T/300)^(-0.21)*exp(44/T);
k19=1.18*10^(-11)*exp(-20413/T);
k20=9.82*10^(-12)*(T/300)^(-0.21)*exp(-5.20/T);
k21=1.66*10^(-10)*exp(-14100/T);
k22=1.15*10^(-10)*exp(-13400/T);
k23=2.51*10^(-10)*exp(-38602/T);
k24=1.00*10^(-10);
k25=5.37*10^(-11)*exp(-13800/T);
k26=5.00*10^(-11)*exp(-200/T);
k27=5.00*10^(-12)*exp(-900/T);
k28=6.00*10^(-12);
k29=2.00*10^(-10)*(T/300)^(-0.12);
k30=1.30*10^(-10);
k31=5.00*10^(-10);
k32=2.90*10^(-10);
k33=3.10*10^(-10);
k34=3.66*10^(-11);
k35=2.26*10^(-12)*(T/300)^(0.86)*exp(-3134.0/T);
k36=1.80*10^(-10);
k37=3.38*10^(-11)*(T/300)^(-0.17)*exp(2.80/T);
k38=3*10^(-12);
k39=1.00*10^(-10)*(T/300)^(0.4);
k40=6.10*10^(-10);
k41=1*10^(-13);
k42=5*10^(-11);
k43=5*10^(-11);
k44=1.15*10^(-10);
k45=4*10^(-11);
k46=1*10^(-10);
k47=1.50*10^(-11)*exp(-4300/T);
k48=5.3137*10^(-10);
k49=4.10*10^(-10);
k50=8*10^(-10);
k51=1.32*10^(-12);
k52=2.8*10^(-12)*exp(23400/T);
k53=1*10^(-10)*exp(-55200/T);
k54=8.60*10^(-11);
k55=5.99*10^(-12)*exp(-24075/T);
k56=2.02*10^(-11)*(T/300)^(-0.19)*exp(31.90/T);
k57=2.51*10^(-11)*exp(-30653/T);
k58=3.47*10^(-11)*(T/300)^(-1.33)*exp(-242/T);
%constante k::Recombinaison dissociative
k59=3*10^(-7)*(T/300)^(-0.5);
k60=1.70*10^(-7)*(T/300)^(-0.3);
%constante k::Association radiative
k61=4.01*10^(-18)*(T/300)^(0.17)*exp(-101.50/T);
k62=1.08*10^(-18)*(T/300)^(0.07)*exp(-57.50/T);
k63=3.14*10^(-18)*(T/300)^(-0.15)*exp(-68/T);
k64=4.00*10^(-14)*(T/300)^(-1);
k65=4.69*10^(-19)*(T/300)^(1.52)*exp(50.50/T);
k66=3*10^(-16)*(T/300)^(-1);
k67=5*10^(-10)*(T/300)^(-3.7)*exp(-800/T);
k68=4.36*10^(-18)*(T/300)^(0.35)*exp(-161.30/T);
k69=5.72*10^(-19)*(T/300)^(0.37)*exp(-51/T);
k70=4.90*10^(-20)*(T/300)^(1.58);
k71=3.71*10^(-18)*(T/300)^(0.24)*exp(-26.10/T);
%constante k::Transfert de charge
k72=7.05*10^(-10)*(T/300)^(-0.03)*exp(16.70/T);
k73=1*10^(-10);
k74=1.10*10^(-10);
k75=1.10*10^(-10);
k76=1.10*10^(-10);
k77=5.20*10^(-11);
k78=4.80*10^(-10);
k79=4.90*10^(-12)*(T/300)^(0.5)*exp(-4580/T);
k80=6.30*10^(-10);
k81=1.90*10^(-11);
k82=1.00*10^(-11);
k83=1.40*10^(-10);
k84=7.30*10^(-10)*exp(-890/T);
k85=1.00*10^(-9);
k86=3.11*10^(-10);
k87=1.00*10^(-9);
k88=8.25*10^(-10);
k89=1.10*10^(-9)*(T/300)^(-0.50);
k90=4.51*10^(-10);
k91=1* 10^(-11);
k92=8.40*10^(-10);
k93=8.50^(-10);
k94=3.40*10^(-10);
k95=5*10^(-11);
k96=1.20*10^(-10);
k97=6.60*10^(-10);
k98=4.60*10^(-10);
k99=1.00*10^(-10)*(T/300)^(-0.5);
k100=7.40*10^(-11);
k101=6.30*10^(-10);
k102=3.30*10^(-10);
k103=5.70*10^(-10);
%constante k::Neutralisation
k104=7.51*10^(-8)*(T/300)^(-0.5);
k105=7.51*10^(-8)*(T/300)^(-0.5);
k106=7.51*10^(-8)*(T/300)^(-0.5);
k107=7.51*10^(-8)*(T/300)^(-0.5);
k108=7.51*10^(-8)*(T/300)^(-0.5);
k109=7.51*10^(-8)*(T/300)^(-0.5);
k110=7.51*10^(-8)*(T/300)^(-0.50);
k111=7.51*10^(-8)*(T/300)^(-0.50);
k112=7.51*10^(-8)*(T/300)^(-0.50);
k113=7.51*10^(-8)*(T/300)^(-0.50);
k114=7.51*10^(-8)*(T/300)^(-0.5);
k115=7.51*10^(-8)*(T/300)^(-0.5);
k116=7.51*10^(-8)*(T/300)^(-0.5);
k117=7.51*10^(-8)*(T/300)^(-0.5);
%constante k::Photoionisation
k118=3.10*10^(-10)*exp(-3.3/T);
k119=4.90*10^(-8)*exp(-0.50/T);
k120=1.09*10^(-8)*exp(-0.5/T);
k121=4.10*10^(-10)*exp(-3.80/T);
k122=1*10^(-11)*exp(-1.70/T);
k123=2.96*10^(-9)*exp(-2/T);
k124=5*10^(-9)*exp(-2.1/T);
k125=6.88*10^(-9)*exp(-1.5/T);
k126=7.90*10^(-10)*exp(-2.1/T);
k127=6.10*10^(-9)*exp(-0.5/T);
k128=3.50*10^(-11)*exp(-2/T);
k129=2.30*10^(-10)*exp(-3.90/T);
k130=2.10*10^(-10)*exp(-3.50/T);
k131=4.70*10^(-10)*exp(-2.10/T);
%constante k::Attachement radiatif
k132=2.25*10^(-15);
k133=1.50*10^(-15);
k134=2.00*10^(-15)*(T/300)^(-0.5);
k135=1.70*10^(-14)*(T/300)^(-0.5);
%constante k::Détachement associatif
k136=5*10^(-10);
k137=1*10^(-9);
k138=1*10^(-9);
k139=5*10^(-10);
k140=5*10^(-10);
k141=5*10^(-11);
k142=5*10^(-10);
k143=1.90*10^(-10);
k144=2.90*10^(-10);
k145=3.10*10^(-10)*(T/300)^(-0.83);
k146=6.50*10^(-10);
k147=2.20*10^(-10);
k148=1.00*10^(-9);
k149=1*10^(-9);
%constante k::Excitation/ Désexcitation
k150=2.0425*10^(-7)*T^(-0.3)*exp(-29489.55/T);
k151=2.6*10^(-10)*T^(0.89)*exp(-11372.9/T);
k152=1.2*10^(-10)*T^(0.59)*exp(-18916.15/T);
k153=6.0278*10^(-7)*T^(-0.4)*exp(-19880.36/T);
k154=1.6306*10^(-17);
k155=3.6532*10^(-11);
%constante k::Détachement électronique
k156=7*10^(-10);
k157=2*10^(-10);
k158=3.6*10^(-10);
%constante k::Dissociation par impact électronique
k159= 6.2*10^(-10) * T^(1.75) * exp(-113032.7/T );
%constante k::Ionisation par impact électronique
k160= 7.1*10^(-9) * T^(0.76) * exp(-181038/T );
%constante k::Attachement électronique cn(1)6/s
k161=1*10^(-31);
k162=1*10^(-31);
%constante k::Processus de recombinaison cn(1)6/s
k163=5.70*10^(-34)*T^(-2.6);
k164=2.5919*10^(-80);
k165=1.31*10^(-31)*T^(-1.5);
fun = @(t,n)myequations(t,n,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,k18,k19,k20,k21,k22,k23,k24,k25,k26,k27,k28,k29,k30,k31,k32,k33,k34,k35,k36,k37,k38,k39,k40,k41,k42,k43,k44,k45,k46,k47,k48,k49,k50,k51,k52,k53,k54,k55,k56,k57,k58,k59,k60,k61,k62,k63,k64,k65,k66,k67,k68,k69,k70,k71,k72,k73,k74,k75,k76,k77,k78,k79,k80,k81,k82,k83,k84,k85,k86,k87,k88,k89,k90,k91,k92,k93,k94,k95,k96,k97,k98,k99,k100,k101,k102,k103,k104,k105,k106,k107,k108,k109,k110,k111,k112,k113,k114,k115,k116,k117,k118,k119,k120,k121,k122,k123,k124,k125,k126,k127,k128,k129,k130,k131,k132,k133,k134,k135,k136,k137,k138,k139,k140,k141,k142,k143,k144,k145,k146,k147,k148,k149,k150,k151,k152,k153,k154,k155,k156,k157,k158,k159,k160,k161,k162,k163,k164,k165);
end
function dndt = myequations(t,n,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,k14,k15,k16,k17,k18,k19,k20,k21,k22,k23,k24,k25,k26,k27,k28,k29,k30,k31,k32,k33,k34,k35,k36,k37,k38,k39,k40,k41,k42,k43,k44,k45,k46,k47,k48,k49,k50,k51,k52,k53,k54,k55,k56,k57,k58,k59,k60,k61,k62,k63,k64,k65,k66,k67,k68,k69,k70,k71,k72,k73,k74,k75,k76,k77,k78,k79,k80,k81,k82,k83,k84,k85,k86,k87,k88,k89,k90,k91,k92,k93,k94,k95,k96,k97,k98,k99,k100,k101,k102,k103,k104,k105,k106,k107,k108,k109,k110,k111,k112,k113,k114,k115,k116,k117,k118,k119,k120,k121,k122,k123,k124,k125,k126,k127,k128,k129,k130,k131,k132,k133,k134,k135,k136,k137,k138,k139,k140,k141,k142,k143,k144,k145,k146,k147,k148,k149,k150,k151,k152,k153,k154,k155,k156,k157,k158,k159,k160,k161,k162,k163,k164,k165)
dndt = zeros(43,1);
%Chemin C
dndt(1)= k1*n(2)*n(42)+k15*n(5)*n(9)+k16*n(5)*n(26)+k25*n(4)*n(26)...
+k29*n(4)*n(9)+k33*n(4)*n(11)+k39*n(7)*n(26)+k40*n(7)*n(27)...
+k43*n(7)*n(9)+k44*n(7)*n(12)+k48*n(9)*n(9)+k59*n(11)*n(42)...
+k72*n(2)*n(24)+k104*n(3)*n(5)+k105*n(3)*n(8)+k106*n(3)*n(25)...
+k117*n(18)*n(2)+k119*n(3)*n(43)+k122*n(11)*n(43)+k124*n(13)*n(43)...
+k130*n(22)*n(43)-k4*n(1)*n(37)-k8*n(1)*n(15)-k9*n(1)*n(26)...
-k10*n(1)*n(24)-k11*n(1)*n(20)-k12*n(1)*n(35)-k13*n(1)*n(22)...
-k61*n(2)*n(1)-k64*n(1)*n(13)-k65*n(1)*n(4)-k66*n(1)*n(9)...
-k67*n(1)*n(5)-k68*n(1)*n(1)-k69*n(1)*n(7)-k73*n(1)*n(27)...
-k74*n(1)*n(11)-k75*n(1)*n(23)-k76*n(1)*n(21)-k77*n(1)*n(19)...
-k118*n(1)*n(43)-k132*n(1)*n(42)-k136*n(1)*n(3)-k137*n(1)*n(12)...
-k138*n(1)*n(14)-k139*n(1)*n(6)-2*k164*n(7)*n(1)*n(1)...
+k164*n(7)*n(1)^(2)-k163*n(1)*n(15)*n(4)+k163*n(1)*n(15)*n(4)...
-k165*n(1)*n(32)*n(4)+k165*n(1)*n(32)*n(4);
%Chemin C+:
dndt(2)= -k1*n(2)*n(42)-k5*n(2)*n(15)-k6*n(2)*n(30)-k7*n(2)*n(36)...
-k61*n(2)*n(1)-k62*n(2)*n(7)-k63*n(2)*n(4)-k72*n(2)*n(24)...
-k117*n(18)*n(2)+k45*n(7)*n(11)+k73*n(1)*n(27)...
+k74*n(1)*n(11)+k75*n(1)*n(23)+k76*n(1)*n(21)...
+k77*n(1)*n(19)+k118*n(1)*n(43)+k122*n(11)*n(43);
%Chemin C-:
dndt(3)= k32*n(4)*n(12)+k132*n(1)*n(42)-k14*n(3)*n(24)-k104*n(3)*n(5)...
-k105*n(3)*n(8)-k106*n(3)*n(25)-k119*n(3)*n(43)...
-k136*n(1)*n(3)-k140*n(3)*n(4)-k141*n(3)*n(15)-k142*n(3)*n(7);
%Chemin O:
dndt(4)= k2*n(5)*n(42)+k8*n(1)*n(15)+k10*n(1)*n(24)+k13*n(1)*n(22)...
+k14*n(3)*n(24)+k35*n(7)*n(15)+k36*n(7)*n(19)...
+k37*n(7)*n(24)+k38*n(7)*n(32)+k52*n(15)*n(24)...
+k53*n(15)*n(20)+k55*n(15)*n(22)+k56*n(15)*n(26)...
+k78*n(5)*n(9)+k79*n(5)*n(22)+k80*n(5)*n(30)...
+k81*n(5)*n(15)+k84*n(6)*n(15)+k85*n(6)*n(26)...
+k104*n(3)*n(5)+2*k107*n(5)*n(6)+k108*n(5)*n(18)...
+k109*n(6)*n(8)+k120*n(6)*n(43)+2*k126*n(15)*n(43)...
+k128*n(19)*n(43)+k130*n(22)*n(43)+k131*n(24)*n(43)...
+k162*n(42)*n(4)*n(15)-k19*n(4)*n(24)-k20*n(4)*n(32)...
-k21*n(4)*n(30)-k22*n(4)*n(30)-k23*n(4)*n(20)...
-k24*n(4)*n(36)-k25*n(4)*n(26)-k26*n(4)*n(26)...
-k27*n(4)*n(13)-k28*n(4)*n(35)-k29*n(4)*n(9)...
-k30*n(4)*n(21)-k31*n(4)*n(14)-k32*n(4)*n(12)...
-k33*n(4)*n(11)-k63*n(2)*n(4)-k65*n(1)*n(4)...
-2*k70*n(4)*n(4)-k82*n(4)*n(21)-k83*n(4)*n(23)...
-k133*n(4)*n(42)-k140*n(3)*n(4)-k143*n(4)*n(6)...
-k144*n(4)*n(12)-k161*n(42)*n(4)*n(15)...
-k162*n(42)*n(4)*n(15)-k163*n(1)*n(15)*n(4)-k165*n(1)*n(32)*n(4);
%Chemin O+:
dndt(5)= k5*n(2)*n(15)+k34*n(8)*n(15)+k82*n(4)*n(21)...
+k83*n(4)*n(23)+k128*n(19)*n(43)-k19*n(4)*n(24)...
-k20*n(4)*n(32)-k21*n(4)*n(30)-k22*n(4)*n(30)...
-k23*n(4)*n(20)-k24*n(4)*n(36)-k25*n(4)*n(26)...
-k26*n(4)*n(26)-k27*n(4)*n(13)-k28*n(4)*n(35)...
-k29*n(4)*n(9)-k30*n(4)*n(21)-k31*n(4)*n(14)...
-k32*n(4)*n(12)-k33*n(4)*n(11)-k63*n(2)*n(4)...
-k65*n(1)*n(4)-2*k70*n(4)*n(4)-k82*n(4)*n(21)...
-k83*n(4)*n(23)-k133*n(4)*n(42)-k140*n(3)*n(4)...
-k143*n(4)*n(6)-k144*n(4)*n(12)-k161*n(42)*n(4)*n(15)...
-k162*n(42)*n(4)*n(15)-k163*n(1)*n(15)*n(4)-k165*n(1)*n(32)*n(4);
%Chemin O-:
dndt(6)= k133*n(4)*n(42)+k161*n(42)*n(4)*n(15)-k84*n(6)*n(15)...
-k85*n(6)*n(26)-k107*n(5)*n(6)-k109*n(6)*n(8)-k120*n(6)*n(43)...
-k139*n(1)*n(6)-k143*n(4)*n(6)-k145*n(6)*n(24)...
-k146*n(6)*n(22)-k147*n(6)*n(7);
%Chemin N:
dndt(7)= k3*n(8)*n(42)+k9*n(1)*n(26)+k11*n(1)*n(20)...
+k18*n(5)*n(20)+k19*n(4)*n(24)+k23*n(4)*n(20)...
+k26*n(4)*n(26)+k30*n(4)*n(21)+2*k60*n(21)*n(42)...
+k86*n(8)*n(15)+k87*n(8)*n(9)+k88*n(8)*n(22)...
+k89*n(8)*n(26)+k90*n(8)*n(24)+k105*n(3)*n(8)...
+k109*n(6)*n(8)+k110*n(8)*n(18)+k111*n(8)*n(12)...
+k112*n(8)*n(28)+k113*n(8)*n(14)+2*k129*n(20)*n(43)...
+k131*n(24)*n(43)+2*k159*n(20)*n(42)-k35*n(7)*n(15)...
-k36*n(7)*n(19)-k37*n(7)*n(24)-k38*n(7)*n(32)...
-k39*n(7)*n(26)-k40*n(7)*n(27)-k41*n(7)*n(13)...
-k42*n(7)*n(14)-k43*n(7)*n(9)-k44*n(7)*n(12)...
-k45*n(7)*n(11)-k46*n(7)*n(35)-k62*n(2)*n(7)...
-k69*n(1)*n(7)-k71*n(8)*n(7)-k91*n(7)*n(21)...
-k142*n(3)*n(7)-k147*n(6)*n(7)-k164*n(7)*n(1)*n(1);
%Chemin N+:
dndt(8)=k91*n(7)*n(21)-k3*n(8)*n(42)-k34*n(8)*n(15)...
-k71*n(8)*n(7)-k86*n(8)*n(15)-k87*n(8)*n(9)...
-k88*n(8)*n(22)-k89*n(8)*n(26)-k90*n(8)*n(24)...
-k105*n(3)*n(8)-k109*n(6)*n(8)-k110*n(8)*n(18)...
-k111*n(8)*n(12)-k112*n(8)*n(28)-k113*n(8)*n(14);
%Chemin C2:
dndt(9)= k9*n(1)*n(26)+k12*n(1)*n(35)+k13*n(1)*n(22)...
+k27*n(4)*n(13)+k41*n(7)*n(13)+k68*n(1)*n(1)...
+k74*n(1)*n(11)+k94*n(11)*n(24)+k111*n(8)*n(12)...
+k115*n(12)*n(25)+k121*n(9)*n(43)+k123*n(12)*n(43)...
+k124*n(13)*n(43)+k136*n(1)*n(3)+k153*n(10)*n(42)...
+k154*n(10)^(2)+k155*n(10)*n(1)-k15*n(5)*n(9)...
-k29*n(4)*n(9)-k43*n(7)*n(9)-k47*n(9)*n(15)...
-2*k48*n(9)*n(9)-k49*n(9)*n(19)-k66*n(1)*n(9)...
-k78*n(5)*n(9)-k87*n(8)*n(9)-k92*n(9)*n(23)...
-k93*n(9)*n(27)-k121*n(9)*n(43)-k134*n(9)*n(42)...
-k148*n(9)*n(12)-k150*n(9)*n(42);
%Chemin C2(d)
dndt(10)= k150*n(9)*n(42)-k153*n(10)*n(42)-k154*n(10)^(2)-k155*n(10)*n(1);
%Chemin C2+:
dndt(11)= k61*n(2)*n(1)+k78*n(5)*n(9)+k87*n(8)*n(9)...
+k92*n(9)*n(23)+k93*n(9)*n(27)-k33*n(4)*n(11)...
-k45*n(7)*n(11)-k50*n(11)*n(15)-k59*n(11)*n(42)...
-k74*n(1)*n(11)-k94*n(11)*n(24)-k122*n(11)*n(43);
%Chemin C2-:
dndt(12)= -k32*n(4)*n(12)-k44*n(7)*n(12)-k111*n(8)*n(12)...
-k115*n(12)*n(25)-k123*n(12)*n(43)-k137*n(1)*n(12)...
-k144*n(4)*n(12)-k148*n(9)*n(12)-k149*n(12)*n(13)...
+k31*n(4)*n(14)+k42*n(7)*n(14)+k134*n(9)*n(42);
%Chemin C3:
dndt(13)= k48*n(9)*n(9)+k66*n(1)*n(9)+k113*n(8)*n(14)...
+k114*n(14)*n(25)+k125*n(14)*n(43)...
+k137*n(1)*n(12)-k27*n(4)*n(13)-k41*n(7)*n(13)...
-k64*n(1)*n(13)-k124*n(13)*n(43)...
-k135*n(13)*n(42)-k149*n(12)*n(13);
%Chemin C3-:
dndt(14)= k135*n(13)*n(42)-k31*n(4)*n(14)...
-k42*n(7)*n(14)-k113*n(8)*n(14)...
-k114*n(14)*n(25)-k125*n(14)*n(43)-k138*n(1)*n(14);
%Chemin O2:
dndt(15)= k17*n(5)*n(32)+k19*n(4)*n(24)+k20*n(4)*n(32)...
+k21*n(4)*n(30)+k57*n(24)*n(24)+k70*n(4)*n(4)...
+k77*n(1)*n(19)+k97*n(19)*n(32)+k98*n(19)*n(24)...
+k108*n(5)*n(18)+k110*n(8)*n(18)+k116*n(18)*n(25)...
+k117*n(18)*n(2)+k127*n(18)*n(43)+k143*n(4)*n(6)...
+2*k156*n(18)*n(15)+2*k157*n(18)*n(16)...
+2*k158*n(18)*n(17)+k161*n(42)*n(4)*n(15)...
-k5*n(2)*n(15)-k8*n(1)*n(15)-k34*n(8)*n(15)...
-k35*n(7)*n(15)-k47*n(9)*n(15)-k50*n(11)*n(15)...
-k51*n(15)*n(37)-k52*n(15)*n(24)-k53*n(15)*n(20)...
-k54*n(15)*n(27)-k55*n(15)*n(22)-k56*n(15)*n(26)...
-k81*n(5)*n(15)-k84*n(6)*n(15)-k86*n(8)*n(15)...
-k95*n(15)*n(21)-k96*n(15)*n(23)-k126*n(15)*n(43)...
-k141*n(3)*n(15)-k151*n(15)*n(42)-k152*n(15)*n(42)...
-k156*n(18)*n(15)-k161*n(42)*n(4)*n(15)...
-k162*n(42)*n(4)*n(15)-k163*n(1)*n(15)*n(4);
%Chemin O2(a):
dndt(16)= k151*n(42)*n(15)-k157*n(16)*n(18);
%Chemin O2(b):
dndt(17)= k152*n(42)*n(15)-k158*n(17)*n(18);
%Chemin O2-:
dndt(18)= k84*n(6)*n(15)+k162*n(42)*n(4)*n(15)...
-k108*n(5)*n(18)-k110*n(8)*n(18)...
-k116*n(18)*n(25)-k117*n(18)*n(2)-k127*n(18)*n(43)...
-k156*n(18)*n(15)-k157*n(18)*n(16)-k158*n(18)*n(17);
%Chemin O2+:
dndt(19)= k81*n(5)*n(15)+k86*n(8)*n(15)+k95*n(15)*n(21)...
+k96*n(15)*n(23)-k36*n(7)*n(19)-k49*n(9)*n(19)-k77*n(1)*n(19)...
-k97*n(19)*n(32)-k98*n(19)*n(24)-k128*n(19)*n(43);
%Chemin N2:
dndt(20)= k21*n(4)*n(30)+k37*n(7)*n(24)+k39*n(7)*n(26)...
+k57*n(24)*n(24)+k76*n(1)*n(21)+k82*n(4)*n(21)...
+k91*n(7)*n(21)+k95*n(15)*n(21)+k99*n(21)*n(26)...
+k100*n(21)*n(22)-k11*n(1)*n(20)-k18*n(5)*n(20)...
-k23*n(4)*n(20)-k53*n(15)*n(20)...
-k129*n(20)*n(43)-k159*n(20)*n(42)-k160*n(20)*n(42);
%Chemin N2+:
dndt(21)= k40*n(7)*n(27)+k71*n(8)*n(7)+k160*n(20)*n(42)...
-k30*n(4)*n(21)-k60*n(21)*n(42)-k76*n(1)*n(21)...
-k82*n(4)*n(21)-k91*n(7)*n(21)...
-k95*n(15)*n(21)-k99*n(21)*n(26)-k100*n(21)*n(22);
%Chemin CO:
dndt(22)= k4*n(1)*n(37)+k5*n(2)*n(15)+k8*n(1)*n(15)...
+k24*n(4)*n(36)+k26*n(4)*n(26)+k27*n(4)*n(13)...
+k28*n(4)*n(35)+k29*n(4)*n(9)+k31*n(4)*n(14)...
+k32*n(4)*n(12)+2*k47*n(9)*n(15)+k49*n(9)*n(19)...
+k50*n(11)*n(15)+k54*n(15)*n(27)+k58*n(24)*n(37)...
+k65*n(1)*n(4)+k75*n(1)*n(23)+k83*n(4)*n(23)...
+k92*n(9)*n(23)+k96*n(15)*n(23)+k102*n(24)*n(23)...
+k139*n(1)*n(6)+k140*n(3)*n(4)-k13*n(1)*n(22)...
-k55*n(15)*n(22)-k79*n(5)*n(22)-k88*n(8)*n(22)...
-k100*n(21)*n(22)-k101*n(22)*n(27)...
-k130*n(22)*n(43)-k146*n(6)*n(22);
%Chemin CO+:
dndt(23)= k7*n(2)*n(36)+k15*n(5)*n(9)+k33*n(4)*n(11)...
+k49*n(9)*n(19)+k50*n(11)*n(15)+k63*n(2)*n(4)...
+k67*n(1)*n(5)+k79*n(5)*n(22)+k88*n(8)*n(22)...
+k100*n(21)*n(22)+k101*n(22)*n(27)-k75*n(1)*n(23)...
-k83*n(4)*n(23)-k92*n(9)*n(23)...
-k96*n(15)*n(23)-k102*n(24)*n(23);
%Chemin NO:
dndt(24)= k20*n(4)*n(32)+k22*2*n(4)*n(30)+k23*n(4)*n(20)...
+k24*n(4)*n(36)+k25*n(4)*n(26)+k34*n(8)*n(15)...
+k35*n(7)*n(15)+k51*n(15)*n(37)+k106*n(3)*n(25)...
+k114*n(14)*n(25)+k115*n(12)*n(25)+k116*n(18)*n(25)...
+k147*n(6)*n(7)-k10*n(1)*n(24)-k14*n(3)*n(24)-k19*n(4)*n(24)...
-k37*n(7)*n(24)-k52*n(15)*n(24)-k57*2*n(24)*n(24)...
-k58*n(24)*n(37)-k72*n(2)*n(24)-k90*n(8)*n(24)...
-k94*n(11)*n(24)-k98*n(19)*n(24)-k102*n(24)*n(23)...
-k103*n(24)*n(27)-k131*n(24)*n(43)-k145*n(6)*n(24);
%Chemin NO+:
dndt(25)= k6*n(2)*n(30)+k16*n(5)*n(26)+k17*n(5)*n(32)...
+k18*n(5)*n(20)+k30*n(4)*n(21)+k36*n(7)*n(19)...
+k54*n(15)*n(27)+k72*n(2)*n(24)+k90*n(8)*n(24)...
+k94*n(11)*n(24)+k98*n(19)*n(24)+k102*n(24)*n(23)...
+k103*n(24)*n(27)-k106*n(3)*n(25)-k114*n(14)*n(25)...
-k115*n(12)*n(25)-k116*n(18)*n(25);
%Chemin CN:
dndt(26)= k4*n(1)*n(37)+k6*n(2)*n(30)+k7*n(2)*n(36)...
+k10*n(1)*n(24)+k11*n(1)*n(20)+k12*n(1)*n(35)...
+k28*n(4)*n(35)+k41*n(7)*n(13)+k42*n(7)*n(14)...
+k43*n(7)*n(9)+k45*n(7)*n(11)+k46*2*n(7)*n(35)...
+k69*n(1)*n(7)+k73*n(1)*n(27)+k93*n(9)*n(27)...
+k101*n(22)*n(27)+k103*n(24)*n(27)...
+k112*n(8)*n(28)+k142*n(3)*n(7)...
+k164*n(7)*n(1)*n(1)-k9*n(1)*n(26)-k16*n(5)*n(26)...
-k25*n(4)*n(26)-k26*n(4)*n(26)-k39*n(7)*n(26)...
-k56*n(15)*n(26)-k85*n(6)*n(26)...
-k89*n(8)*n(26)-k99*n(21)*n(26);
%Chemin CN+:
dndt(27)= k62*n(2)*n(7)+k89*n(8)*n(26)+k99*n(21)*n(26)...
-k40*n(7)*n(27)-k54*n(15)*n(27)...
-k73*n(1)*n(27)-k93*n(9)*n(27)...
-k101*n(22)*n(27)-k103*n(24)*n(27);
%Chemin CN-:
dndt(28)= k14*n(3)*n(24)+k44*n(7)*n(12)+k85*n(6)*n(26)-k112*n(8)*n(28);
%Chemin O3:
dndt(29)= k163*n(1)*n(15)*n(4);
%Chemin N2O:
dndt(30)= k38*n(7)*n(32)+k53*n(15)*n(20)+k58*n(24)*n(37)...
-k6*n(2)*n(30)-k21*n(4)*n(30)...
-k22*n(4)*n(30)-k80*n(5)*n(30);
%Chemin N2O+:
dndt(31)= k80*n(5)*n(30);
%Chemin NO2:
dndt(32)= k52*n(15)*n(24)+k145*n(6)*n(24)-k17*n(5)*n(32)...
-k20*n(4)*n(32)-k38*n(7)*n(32)...
-k97*n(19)*n(32)-k165*n(1)*n(32)*n(4);
%Chemin NO2+:
dndt(33)= k97*n(19)*n(32);
%Chemin CO2:
dndt(34)= k51*n(15)*n(37)+k55*n(15)*n(22)...
+k141*n(3)*n(15)+k146*n(6)*n(22);
%Chemin C2N:
dndt(35)= -k12*n(1)*n(35)-k28*n(4)*n(35)-k46*n(7)*n(35);
%Chemin CNO:
dndt(36)= -k7*n(2)*n(36)-k24*n(4)*n(36);
%Chemin OCN
dndt(37)= -k4*n(1)*n(37)-k51*n(15)*n(37)...
-k58*n(24)*n(37)+k56*n(15)*n(26);
%Chemin NO3:
dndt(38)= k165*n(1)*n(32)*n(4);
%Chemin C2O:
dndt(39)= k144*n(4)*n(12);
%Chemin C4:
dndt(40)= k64*n(1)*n(13)+k138*n(1)*n(14)+k148*n(9)*n(12);
%Chemin C5:
dndt(41)= k149*n(12)*n(13);
%Chemin e:
dndt(42)= k118*n(1)*n(43)+k119*n(3)*n(43)+k120*n(6)*n(43)...
+k121*n(12)*n(43)+k123*n(12)*n(43)+k125*n(14)*n(43)...
+k127*n(18)*n(43)+k136*n(1)*n(3)+k137*n(1)*n(12)...
+k138*n(1)*n(14)+k139*n(1)*n(6)+k140*n(3)*n(4)...
+k141*n(3)*n(15)+k142*n(3)*n(7)+k143*n(4)*n(6)...
+k144*n(4)*n(12)+k145*n(6)*n(24)+k146*n(6)*n(22)...
+k147*n(6)*n(7)+k148*n(9)*n(12)+k149*n(12)*n(13)...
+k150*n(9)*n(42)+k151*n(15)*n(42)+k152*n(15)*n(42)...
+k153*n(10)*n(42)+k156*n(18)*n(15)+k157*n(18)*n(16)...
+k158*n(18)*n(17)+k159*n(20)*n(42)...
+2*k160*n(20)*n(42)-k1*n(2)*n(42)-k2*n(5)*n(42)...
-k3*n(8)*n(42)-k59*n(11)*n(42)-k60*n(21)*n(42)...
-k132*n(1)*n(42)-k133*n(4)*n(42)-k134*n(9)*n(42)...
-k135*n(13)*n(42)-k150*n(9)*n(42)-k151*n(15)*n(42)...
-k152*n(15)*n(42)-k153*n(10)*n(42)-k159*n(20)*n(42)...
-k160*n(20)*n(42)-k161*n(42)*n(4)*n(15)...
-k162*n(42)*n(4)*n(15);
%Chemin hv:
dndt(43)= -k118*n(1)*n(43)-k119*n(3)*n(43)-k120*n(6)*n(43)...
-k121*n(9)*n(43)-k122*n(11)*n(43)-k123*n(12)*n(43)...
-k124*n(13)*n(43)-k125*n(14)*n(43)-k126*n(15)*n(43)...
-k127*n(18)*n(43)-k128*n(19)*n(43)-k129*n(20)*n(43)...
-k130*n(22)*n(43)-k131*n(24)*n(43)+k1*n(2)*n(42)...
+k2*n(5)*n(42)+k3*n(8)*n(42)+k61*n(2)*n(1)...
+k62*n(2)*n(7)+k63*n(2)*n(4)+k64*n(1)*n(13)...
+k65*n(1)*n(4)+k66*n(1)*n(9)+k67*n(1)*n(5)...
+k68*n(1)*n(1)+k69*n(1)*n(7)+k70*n(4)*n(4)...
+k71*n(8)*n(7)+k132*n(1)*n(42)+k133*n(4)*n(42)...
+k134*n(9)*n(42)+k135*n(13)*n(42);
end
댓글 수: 3
Imene Yed
2021년 5월 31일
Walter Roberson
2021년 5월 31일
The code is right there in my post,
function fun = init_fun
immediately after the end of the code for call_myequations
The k-values are not known in call_myequations. So you must pass them to the function:
function [t,n] = call_myequations(k1,k2,...)
( I think it's time to make k a vector: k(1),k(2),... instead of k1,k2,...)
Maybe you want to create another function parameters where you define the k vector:
function k = parameters()
T = ...;
k(1) = ...;
k(2) = ...;
...
end
And a function main that calls all these other functions in order:
function main
k = parameters();
[t,n] = call_myequations(k);
postprocess(t,n);
end
댓글 수: 20
Imene Yed
2021년 5월 31일
Imene Yed
2021년 5월 31일
Imene Yed
2021년 5월 31일
Torsten
2021년 5월 31일
oh, come on.
5 scripts:
main.m, parameters.m, call_myequations.m, myequations.m and postprocess.m.
These are the typical routines in problem solving:
Driver, model parameters, solver calling, problem formulation and solution postprocessing.
Imene Yed
2021년 5월 31일
Walter Roberson
2021년 5월 31일
Sigh. 5 minutes with a decent editor. Most of the time was spent reviewing to be sure that everything was okay.
See attached.
Imene Yed
2021년 5월 31일
Walter Roberson
2021년 5월 31일
So take the code I attached, find the function lines, and write each function to a different file.
Walter Roberson
2021년 6월 1일
Attached is the code broken into individual files.
Stephen23
2021년 6월 1일
"I have to replace all the k1 k2... with k(1) k(2) it will take a huge time ."
Not at all.
With notepad++ it takes just six seconds to type the following regular expressions into the "Replace" dialog box:
Find what = k(\d+)
Replace with = K\($1\)
and then click the big "Replace All" button.
Walter Roberson
2021년 6월 1일
It was a bit more work in vi; I had to use
1,$s/k\(\d\+\)/k(\1)/g
After that it was a matter of deleting out the k variables from the parameter lists.
Imene Yed
2021년 6월 1일
Walter Roberson
2021년 6월 1일
You already had function call_myequations and function myequations. Those were moved into their own files and the k calling sequence was changed to a vector.
You already had a script that initialized all the constants, but did not run the other code.
Your code to call myequations assumed that the constants existed even though it did not run the code that created the constants.
What changed is that I took the script that creates the constants and made it into a function that returns an anonymous function. call_myequations now executes the function that initializes the constants, and calls fminunc using the function handle. And there is a trivial script that starts it all going so that I could run it here for you.
So in terms of your design the change is to move initializing the constants into a function that is expected to return a function handle. Otherwise it is the same code structure you already had.
Imene Yed
2021년 6월 1일
Walter Roberson
2021년 6월 1일
Your function is
function [t,n]= call_myequations()
%condition intial
%IC=IC1=IC2=....=IC43=10^6cm^3
IC =10^19*ones(1,43);
tic,[t,n] = ode45( t,n,k,[0 10e-3],IC);
toc
plot(t,n)
end
You pass k into call_myequations but call_myequations does not expect to receive any inputs, but it expects to use k.
You call ode45() passing in t as the first parameter to it. There is no t defined inside the function, and t is not a parameter to the function, and is not global, and this is not a nested function with shared variables: so in order for that to be valid, t would have to be a function with no parameters that executed to give the function handle that ode45 needs as its first parameter. But t is not a function that you have defined.
The second parameter you pass in to ode45 is n. However, as above, n is not defined in context. It would have to be a function with no parameters that returned a time interval vector... but n is not a function that you have defined.
The third parameter you pass in to ode45 is k. However, as above, k is not defined in context. It would have to be a function with no parameters that returned a vector or array of initial conditons... but k is not a function that you have defined.
The fourth parameter you pass in to ode45 is [0 10e-3] . The fourth parameter to ode45() must be [] or an options structure, never numeric.
The fifth parameter you pass in to ode45 is IC. ode45 has no documented fifth parameter. For historical reasons, parameters after the fourth parameter have their value taken and passed as additional parameters to the ode function. It is not at all recommended to rely on this, as there are cases (especially with some options set) where such parameters would be used for other purposes. MATLAB makes no attempt to be fully backwards compatible for those cases, since the ability to pass extra parameters to function has not been documented for close to 20 years.
Your function myequations expects parameters t, n, and k. If you were to define t as a function that returned a handle to myequations, then IC would be passed in to myequations in the k position.
It would be a lot easier for you to recode call_myequations as
function [t,n]= call_myequations(k)
%condition intial
%IC=IC1=IC2=....=IC43=10^6cm^3
IC =10^19*ones(1,43);
tic,[t,n] = ode45( @(t,n)myequations(t,n,k),[0 10e-3],IC);
toc
plot(t,n)
end
Imene Yed
2021년 6월 2일
Imene Yed
2021년 6월 2일
Walter Roberson
2021년 6월 2일
You would need to return them from function main and assign them to variables.
카테고리
도움말 센터 및 File Exchange에서 Ordinary Differential Equations에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
