MATLAB: No difference in return values when using ode45

MATLABode45

Hi,
I've been trying to model the reaction of a Water Gas Shift reaction using some kinetics I found by using ode45. I have used it before but I must be doing something wrong because the values of x(1)-x(5) are not changing. My code looks like this:
% Script calling Moe Kinetics Function
%

% Ptot=8.85; % (bar) ptot
% xco =.0499; % Mole Fractions
% xco2=.0844;
% xh20=.2913;
% xh2 =.4532;
% xch4=.1212;
%
% Pco =.0501*8.85 = x(1) % Partial Pressures (bar)

% Pco2=.0830*8.85 = x(2)

% Ph2o=.2824*8.85 = x(3)

% Ph2 =.4485*8.85 = x(4)

% Pch4=.1360*8.85 = x(5)

% T = 628.35 K = x(6)
R = 8.3144598*(10^-5); % (m^3*bar)/(K*mol)

% Initial condition Partial Pressures (Pa) and Initial temperature (K)
ic=[.1212*8.85,.4532*8.85,.2913*8.85,.0844*8.85,.0499*8.85,628.35];
Vspan=[0:.001:10];
[V, x] = ode45('project_wgs',Vspan,ic);
% Concentrations of species
x1=x(:,1);
x2=x(:,2);
x3=x(:,3);
x4=x(:,4);
x5=x(:,5);
T=x(:,6); % Temperature (K)
% Converting concentration to partial pressure (bar)
x11 = x1.*R.*T;
x22 = x2.*R.*T;
x33 = x3.*R.*T;
x44 = x4.*R.*T;
x55 = x5.*R.*T;
P = x11+x22+x33+x44+x55; % Total Pressure
function [dCdV] = project_wgs (V,x)
R = 8.3144598*(10^-5); % (m^3*bar)/(K*mol)
q = 5.888*(10^5)/60; % (m^3/min)
% Pco =.0501*8.85 = x(1) % Partial Pressures (bar)
% Pco2=.0830*8.85 = x(2)
% Ph2o=.2824*8.85 = x(3)
% Ph2 =.4485*8.85 = x(4)
% Pch4=.1360*8.85 = x(5)
% T = 628.35 K = x(6)
% Rate equation in (mol/(m^3*min))
r = ((7.26*exp(-15429/R*x(6))*x(1)*x(3)) - (551.6*exp(-53482/R*x(6))*x(2)*x(4)))*1000*7633.65;
dCdV(1,1) = r/q;
dCdV(2,1) = r/q;
dCdV(3,1) = r/q;
dCdV(4,1) = r/q;
dCdV(5,1) = r/q;
dCdV(6,1) = x(6);

Best Answer

One part of your calculation for r is exp(-15429/R*x(6)). With R around 1e-5 and x(6) initially around 600, this term in your expression for r is roughly exp(-1e11). That underflows to 0. Similarly, the other part of your calculation for r, exp(-53482/R*x(6)), also underflows. So r is equal to 0 and so is r/q.