ODE solving ERROR with 5 eq
1 visualización (últimos 30 días)
Mostrar comentarios más antiguos
Hi,
I can't understant the error message from a long but simple code....
clc
clear all
close all
format long
c=3e8;
r=30e-6;
Aeff=pi.*r.*r;
E=1e-6;
lambda0=1030e-9;
om0=2.*c./lambda0;
C0=0;
eta=1E15; %Value can be changed
b2=4e-28;
b3=1.5e-45;
T0=100./om0;
Tr=20./om0;
n2=1e-23;
gamma=2.*pi.*n2./(lambda0.*Aeff);
syms Tp(z) C(z) T(z) Om(z) phi(z)
ode1 = diff(Tp) == (b2-b3.*Om).*C./Tp
ode2 = diff(C) == (b2-b3.*Om).*( (1+C.*C)./Tp.^2 ) + gamma.*E./(sqrt(2.*pi).*Tp) .*(1- Om./om0 )
ode3 = diff(T) == -b2.*Om + b3./2 .*(Om.*Om + (1+C.*C)./2.*Tp.^2 ) + 3.*gamma.*E./(2.*sqrt(2.*pi).*om0.*Tp )
ode4 = diff(Om) == gamma.*E./(sqrt(2.*pi).*Tp.^3) .*(Tr-C./om0) - eta.*E./(sqrt(2.*pi).*Tp)
ode5 = diff(phi) == 0.5.*b2.*(1./(Tp).^2 - Om.^2) + b3.*Om./3 .*(Om.*Om + 3./4 .* (C.^2 -1)./ Tp.^2 ) + 3.*gamma.*E.*(1+Om./om0) ./(4.*sqrt(2.*pi).*Tp) - 0.5.*eta.*E;
odes = [ode1; ode2; ode3; ode4; ode5]
cond1 = Tp(0) == T0;
cond2 = C(0) == 0;
cond3 = T(0) == 0;
cond4 = Om(0) == 0;
cond5 = phi(0) == 0;
conds = [cond1; cond2; cond3; cond4; cond5];
[TpSol(z), CSol(z), TSol(z), OmSol(z), phiSol(z)] = dsolve(odes,conds)
error message =
Error using sym/subsindex (line
685)
Invalid indexing or function
definition. When defining a
function, ensure that the body of
the function is a SYM object. When
indexing, the input must be
numeric, logical or ':'.
Error in ODE (line 45)
[TpSol(z), CSol(z), TSol(z),
OmSol(z), phiSol(z)] =
dsolve(odes,conds)
Do you have an idea?
0 comentarios
Respuestas (1)
Star Strider
el 8 de Ag. de 2019
The error was with:
[TpSol(z), CSol(z), TSol(z), OmSol(z), phiSol(z)] = dsolve(odes,conds)
since MATLAB assumed that ‘TpSol(z)’ and the rest were either function calls or that ‘z’ is an index.
However your system does not have an analytical solution. You will have to integrate it numerically.
Try this:
[VF,Sbs] = odeToVectorField(odes)
odefcn = matlabFunction(VF, 'Vars',{T,Y})
tspan = linspace(0, 100);
ics = zeros(1,4)+eps;
[t,y] = ode15s(odefcn, tspan, ics);
figure
plot(t, y)
grid
lgndc = sprintfc('%s', Sbs);
legend(lgndc, 'Location','E')
8 comentarios
Star Strider
el 12 de Ag. de 2019
@Torsten — Thank you! It definitely could. I didn’t see that (early here).
Ver también
Categorías
Más información sobre Ordinary Differential Equations en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!