Setting up ODE in Matlab for ODE 45 solver

Hi, I am having difficulty setting this problem up for the ODE45 solver. It is a kinetics equation, relating R (rate constants as a function of concentrations C), T (Temperature) and h (heat of formation) all of which vary with distance x. There are 3 species that take part in the reaction.
How do I formulate the function that calculates all these parameters to call in the ODE45 solver ( something like dTdx = .... ) How to plot x,T and x,C (what would be the correct syntax to plot these). Any help with a pseudo code will be greatly appreciated.

7 commentaires

Torsten
Torsten le 6 Mar 2017
The molar balances run fine now ?
https://de.mathworks.com/matlabcentral/answers/327213-trouble-with-ode45-for-an-array-of-values
Best wishes
Torsten.
DIP
DIP le 6 Mar 2017
Yes Torsten. I am a noob. There were discrepancies in the units which I somehow missed. I will update the answers in a bit.
Torsten
Torsten le 6 Mar 2017
Just add the equation for T to the three molar balances and solve the four ODEs simultaneously.
Best wishes
Torsten.
function dTdx=hw2_5(x,C,T)
% Constants
A = 2.2e10;
Ea = 130000;
Ru = 8.3145;
rho=1.3;
cp=1200;
u=1.5;
% Rate Constants and heat of formation
h(1)=-110530;
h(2)=0;
h(3)=-393520;
R(1) = -A * (exp(-Ea /(Ru*T))) * C(1) * sqrt(C(2));
R(2) = 0.5 * R(1);
R(3) = -R(1);
dTdx=-(sum(R(1)*h(1)+R(2)*h(2)-R(3)*h(3)))/(cp*rho*u);
% Solution 4: Chemical Species as a Function of Axial Distance
clc
clear
% Constants
L = 0.12;
N = 39;
delx = L/(N-1);
xspan = 0:delx:0.12;
C0 = [0.52 0.99 0.0 750];
% ODE 45 Code
[x,C]=ode45('hw2_5',xspan,C0);
% Plots
% subplot(2,1,1);
plot(x,C(:,1),'b--o',x,C(:,2),'g--+',x,C(:,3),'r--s')
legend('C_{CO}','C_{O2}','C_{CO2}')
xlabel('Axial (x) Direction [m]')
ylabel('Concentrations [mol/m^3]')
title('Molar Concentration vs Distance')
DIP
DIP le 7 Mar 2017
Hi Torsten, I am getting matrix dimension errors.
DIP
DIP le 7 Mar 2017
Modifié(e) : DIP le 7 Mar 2017
hi Torsten, there are only 2 ODEs
dC/dx=R
dT/dx=RH
R and H are functions of T and C. R and C are arrays.

Connectez-vous pour commenter.

 Réponse acceptée

Torsten
Torsten le 7 Mar 2017
function dydx = hw2_5(x,y)
C(1) = y(1);
C(2) = y(2);
C(3) = y(3);
T = y(4);
% Constants
A = 2.2e10;
Ea = 130000;
Ru = 8.3145;
rho=1.3;
cp=1200;
u=1.5;
% Rate Constants and heat of formation
h(1)=-110530;
h(2)=0;
h(3)=-393520;
R(1) = -A * (exp(-Ea /(Ru*T))) * C(1) * sqrt(C(2));
R(2) = 0.5 * R(1);
R(3) = -R(1);
dTdx=-(R(1)*h(1)+R(2)*h(2)+R(3)*h(3))/(cp*rho*u);
dydx = [R(1);R(2);R(3);dTdx]
Best wishes
Torsten.

1 commentaire

DIP
DIP le 8 Mar 2017
Modifié(e) : DIP le 8 Mar 2017
That worked awesome. Thank You Torsten. Could you recommend any ODE book that has examples solved in matlab ???? Thank you again

Connectez-vous pour commenter.

Plus de réponses (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by