Invalid indexing or function definition.
Afficher commentaires plus anciens
clc;
clear all;
close all;
syms y1(t) y2(t) y3(t) y4(t) y5(t) y6(t) y7(t) y8(t)
C = 3*10^10;
KT = 0.086173324; % MeV
M0 = 1.0078; % u
M1 = 1.0078; % u
M2 = 2.014102; % u
M3 = 1.0078; % u
M4 = 2.014102; % u
M5 = 2.014102; % u
M6 = 3.016030; % u
M7 = 2.014102; % u
M8 = 3.016030; % u
M9 = 3.016030; % u
M10 = 3.016030; % u
M11 = 4.002603; % u
Z0 = 1;
Z1 = 1;
Z2 = 1;
Z3 = 1;
Z4 = 1;
Z5 = 1;
Z6 = 2;
Z7 = 1;
Z8 = 2;
Z9 = 2;
Z10 = 2;
Z11 = 2;
Eat = 0.98951013;
m01 = (M0*M1)/(M0+M1);
m23 = (M2*M3)/(M2+M3);
m45 = (M4*M5)/(M4+M5);
m67 = (M6*M7)/(M6+M7);
m89 = (M8*M9)/(M8+M9);
m1011 = (M10*M11)/(M10+M11);
cont = (931.494)./(3*10^10)^2;
cont1 = 10^-24;
sy0 = 0.1;
Q2 = (4.08*10^-15)*((y8./KT)^(-2/3))*exp(-3.381*((y8./KT)^(-1/3)))*(1+(3.82*y8./KT)+1.51*((y8./KT)^2)+0.144*(y8./KT)^3-1.14*(10^-2)*(y8./KT)^4)./(6.023*10^23);
Q4 =(2.58*10^3)*((y8./KT)^(-2/3))*exp(-3.721*((y8./KT)^(-1/3)))*(1+(3.96*y8./KT)+0.116*(y8./KT)^2)./(6.023*10^23);
Q5 = ((5.59*10^10)*(((y8./KT))^(-2/3))*exp(-12.277*((y8./KT)^(-1/3)))*(1-0.135*(y8./KT)+2.54*0.01*((y8./KT)^2)-1.29*0.001*((y8./KT)^3)))./(6.023*10^23);
D = (2.07*10^14)*((y8./KT)^(3/2))*(1-0.860*((y8./KT)^(1/2))+0.429*(y8./KT))*exp(-25.82/(y8./KT));
nB = (0.6*10^-9)*(((2*1.202)./(pi.^2))*(((10^6)*sy0.*8065.73).^3));
nP = 0.75*nB;
nN = 0.25*nB;
nEn = ((((3*1.202)./(2*pi.^2))*(((10^6)*sy0.*8065.73).^3)))+ nP;
nEp = ((((3*1.202)./(2*pi.^2))*(((10^6)*sy0.*8065.73).^3)))*(exp((-0.511)./(sy0)));
Wtn2 = 1000;
Wconst2 = ((1.239*10^-10)^3);
Wb2 = ((7.56*10^-16)*((1.239*10^-6)^3))./((1.6*10^-19)*((8.617*10^-5)^4));
WV2 = (3*10^10)*sqrt(1-(1./((y8/0.511)+1))^2);
WF2 = (((1.166*10^-23)./(8065.73*10^6))^2)*((y8)^2)*(exp(-(0.783)./y8));
WF1 = (((1.166*10^-23)./(8065.73*10^6))^2)*((y8)^2);
L1 = (y3./Wtn2)+(y5*D)-(y1*y4*WF2*WV2)-(y1*y2*WF1*WV2);
L2 = (((y4)^2)*Q2*0.5)-(y2*y3*WF1*WV2)-(y1*y2*WF1*WV2);
L3 = (y1*y4*WF2*WV2)+(y5*D)-(y2*y3*WF1*WV2)-(y3./Wtn2);
L4 = (y3./Wtn2)+(y2*y3*WF1*WV2)+(y5*D)+(((y6)^2)*Q5)-(y1*y4*WF2*WV2)-(((y4)^2)*Q2*0.5)-(y4*y5*Q4);
L5 = (((y4)^2)*Q2*0.5)-(y5*D)-(y4*y5*Q4);
L6 = (y4*y5*Q4)-(((y6)^2)*Q5*0.5);
L7 = (((y6)^2)*Q5*0.5);
S1 = (((((y4)^2)*Q2*0.5)*1.44)-((y5*D)*2.22)+((y4*y5*Q4)*5.49)+((((y6)^2)*Q5*0.5)*12.86));
S2 = (((y2*y3*WF1*WV2))+((y1*y4*WF2*WV2))+(2*(y1*y2*WF1*WV2))+((y3./Wtn2))+((((y4)^2)*Q2*0.5)))*y8;
L8 = ((Wconst2).*(S1-S2))./(4*Wb2.*((y8)^3));
ode1 = diff(y1)==L1;
ode2 = diff(y2)==L2;
ode3 = diff(y3)==L3;
ode4 = diff(y4)==L4;
ode5 = diff(y5)==L5;
ode6 = diff(y6)==L6;
ode7 = diff(y7)==L7;
ode8 = diff(y8)==L8;
odes = [ode1;ode2;ode3;ode4;ode5;ode6;ode7;ode8];
% S = dsolve(odes)
cond1 = y1(0)==nEn;
cond2 = y2(0)==nEp;
cond3 = y3(0)==nN;
cond4 = y4(0)==nP;
cond5 = y5(0)==0;
cond6 = y6(0)==0;
cond7 = y7(0)==0;
cond8 = y8(0)==0.1;
conds = [cond1;cond2;cond3;cond4;cond5;cond6;cond7;cond8];
[y1Sol(t),y2Sol(t),y3Sol(t),y4Sol(t),y5Sol(t),y6Sol(t),y7Sol(t),y8Sol(t)] = dsolve(odes,conds);
t = linspace(10^-25,10^-22,10);
loglog(t,y1Sol(t),'k-','linewidth',2)
Réponses (1)
Cris LaPierre
le 4 Nov 2021
Remove the (t) from the output variables for dsolve. It is treating t as an index.
[y1Sol,y2Sol,y3Sol,y4Sol,y5Sol,y6Sol,y7Sol,y8Sol] = dsolve(odes,conds);
2 commentaires
mohamed elkasemy
le 4 Nov 2021
Cris LaPierre
le 4 Nov 2021
The reason is because of the warning message that is still there (Unable to find symbolic solution). t is a 1x10 vector, but y1Sol is an empty vector.
Catégories
En savoir plus sur Mathematics dans Centre d'aide et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!