Solving and plotting ode's in matlab

%reaction data
E1=42; % forward activation energy 
Cn2=1; % feed concentration
R =8.314;
A1=.04;
A2=.1;
Co2=1;
Fa0= Cn2 + Co2;
T0=500 ;  % feed temperature
Cpn2 =29.7198; Cpo2=29.051; Cpno =29.3018; %j/mol.K
deltaCp = 2*Cpno-Cpo2-Cpno; 
Hrxn = 180.5; % kj/mol at refrence temp 300K
Ua =0.01;
Ta0 =550; 
mc = 0.011
Cpc= 34.5

syms X(V) T(V) Ta(V) 

ra = -A1*exp(-E1/(8.314*T))*(Cn2^2)*((1-X)^2)*((T0/T)^2);


eqn1 = diff(X, V) == ra/Fa0; %X

eqn2 = diff(T, V) == (Ua*(Ta-T)+ra*Hrxn)/(Fa0*(Cn2*Cpn2+Co2*Cpo2)); %T

eqn3 = diff(Ta, V) == Ua*(Ta-T)/(mc*Cpc);


[odes, vars] = odeToVectorField(eqn1, eqn2, eqn3);
fun = matlabFunction(odes,'Vars',{'t','Y'});
x0 = [0 500 550 ];
tspan = [0 550];
[t, sol] = ode45(fun,tspan,x0);
X = sol(:,1);
T = sol(:,2);
Ta = sol(:,3);

plot(t,sol(:,1))
hold
grid
title('X vs Volume')
xlabel('Volume')
ylabel('X')
hold off

plot(t,sol(:,2),t,sol(:,3)) %plotting Ta on same graph
hold 
grid
title('T vs Volume')
xlabel('Volume')
ylabel('T')

hold off

for i = 1:501
    Ra(i) = Fa0*(sol(i,2)/t(i,1));
    plot(t,Ra)
    hold
    grid
    title('-ra vs V')
    xlabel('Volume')
    ylabel('-ra')
    hold off
    
end

• This is my new code. Although i think i properly indexed everything the output is not as expected.
• Except the initialized values all values of sol array is NaN for some reason.

[Walter Roberson]

We need all of your code in text form in order to test it ourselves.

[Vv,Yv]=ode45('varT',[0:500],[0;500]);

That is two initial conditions. You cannot have 3 outputs for two initial conditions.

[Sujay Bharadwaj]

I forgot to attach the edited version but i did add 3 initial conditions to my code. I edited my question and added all my code in text form. i can link the .m files if necessary

[Sujay Bharadwaj]

Moreover i wrote a similar code which executed successfully but it had only two ode's while this one has 3. Should i give you access to our matlab drive just for reference?

[Walter Roberson]

You must pass in one initial condition for each output -- you need

[Sujay Bharadwaj]

I just tweaked the code with 2 ordinary differential equation's (2 outputs) a bit to get this one. If you could tell me where to make the changes to iitialize the new variable Ta0.

Do i change something here ?

[Vv,Yv]=ode45('varT',[0:500],[0;500],[0,500]);

[Sujay Bharadwaj]

i changed the entire code and made it a lot more straightforward and concise now alothough i am not getting any errors the plots are not as expected. I am editing the question with the revised code.