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 件のコメント

Torsten
Torsten 2017 年 3 月 6 日
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 2017 年 3 月 6 日
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 2017 年 3 月 6 日
Just add the equation for T to the three molar balances and solve the four ODEs simultaneously.
Best wishes
Torsten.
DIP
DIP 2017 年 3 月 7 日
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);
DIP
DIP 2017 年 3 月 7 日
% 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 2017 年 3 月 7 日
Hi Torsten, I am getting matrix dimension errors.
DIP
DIP 2017 年 3 月 7 日
編集済み: DIP 2017 年 3 月 7 日
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.

サインインしてコメントする。

 採用された回答

Torsten
Torsten 2017 年 3 月 7 日

1 投票

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 件のコメント

DIP
DIP 2017 年 3 月 8 日
編集済み: DIP 2017 年 3 月 8 日
That worked awesome. Thank You Torsten. Could you recommend any ODE book that has examples solved in matlab ???? Thank you again

サインインしてコメントする。

その他の回答 (0 件)

カテゴリ

質問済み:

DIP
2017 年 3 月 6 日

編集済み:

DIP
2017 年 3 月 8 日

Community Treasure Hunt

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

Start Hunting!

Translated by