Use of ODE45 for concentration plot help
10 次查看(过去 30 天)
显示 更早的评论
Currently writing a code for the modelling of a chemical reactor for teh reaction of NH3+HCl -> NH4Cl by solving a first order ODE to gain a concentration plot.
I have tried using ODE45 but I cannot get the code to work to consider all three components in the system.
My boundary conditions are initial reactant concentrations are both 0.0109 mol/m3 and my constant k = 2.02*10^10 m3/s. Final reactant concentration = 0 assuming 100% converion.
The kinetic equation is just -ra = k*Cnh3*Chcl
Any help is appreciated, my code is below
%% Matlab code for Irreversible first order ammonium chloride synthesis reaction
%% Using ode45
function xdot = f(t,x)
% Rate constants
k1 = 2.02;
% Reactants
cA = x(1);
cB = x(2);
% Reaction rates
r1 = k1*cA*cB;
% Differential equations
xdot(1) = -r1;
xdot(2) = r1;
xdot = [xdot(1); xdot(2)];
%% Initial conditions
cA0 = 0.0109;
cB0 = 0.0109;
x0 = [cA0; cB0];
%% Time span
tspan = [0 6];
%% Solving ODE
[t,x] = ode45(@f,tspan,x0);
%% Plotting
plot(t,x(:,1),'r-',t,x(:,2),'b-');
title('Ammonium Chloride Synthesis (Irreversible)');
xlabel('time');
end
0 个评论
回答(1 个)
Bjorn Gustavsson
2023-2-14
To the best of my understanding you should include all three species in the ode-function. Perhaps something like this:
function xdot = f(t,x)
% Rate constants
k1 = 2.02;
% Reactants
cA = x(1);
cB = x(2);
cC = x(3);
% Reaction rates
r1 = k1*cA*cB;
% Differential equations
xdot = zeros(3,1);
xdot(1) = -r1;
xdot(2) = -r1;
xdot(3) = r1;
Then you specify the initial conditions for all 3 species:
%% Initial conditions
cA0 = 0.0109;
cB0 = 0.0109;
cC0 = 0;
x0 = [cA0; cB0; cC0];
%% Time span
tspan = [0 160]; % It seems necessary to increase time-span - slow reaction?
%% Solving ODE
[t,x] = ode45(@f,tspan,x0);
%% Plotting
ph = plot(t,x);
set(ph(1),'color','r');
set(ph(2),'color','b')
set(ph(3),'color','g','linewidth',2)
title('Ammonium Chloride Synthesis (Irreversible)');
xlabel('time');
Then I get curves that seem familiar to me - but I'm not an expert in this type of chemistry so I don't know if reaction-response should be on these time-scales or if they are completely off.
HTH
5 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!