Application Center - Maplesoft

App Preview:

Lesson 23 -- Application: Chemical Reactions

You can switch back to the summary page by clicking here.

Learn about Maple
Download Application


 

Lesson23.mw

ORDINARY DIFFERENTIAL EQUATIONS POWERTOOL

Lesson 23 -- Application: Chemical Reactions

Prof. Douglas B. Meade

Industrial Mathematics Institute

Department of Mathematics

University of South Carolina

Columbia, SC 29208

URL:   http://www.math.sc.edu/~meade/

E-mail: meade@math.sc.edu

Copyright  2001  by Douglas B. Meade

All rights reserved

-------------------------------------------------------------------

>

Outline of Lesson 23

23.A First-Order Reactions

23.B Second-Order Reactions

23.C Example:  2NO + O``[2]  -`>`  2NO``[2]

>

Initialization

> restart;

> with( DEtools ):

> with( plots ):

Warning, the name changecoords has been redefined

>

23.A First-Order Reactions

An example of a first-order (chemical) reaction is the conversion of t -butyl chloride into t -butyl alcohol, a reaction expressed chemically with the notation

(CH``[3] )``[3] CCl + NaOH  -`>`  (CH``[3] )``[3] COH + NaCl

A chemical reaction is first-order if the molecules of a substance A decompose into smaller molecules at a rate proportional to the amount of substance A remaining at any time. Thus, if A(t) denotes the amount of substance A at time t , then the amount of substance A can be modeled with the first-order linear (and separable) ODE

> ode1 := diff( A(t), t ) = -k * A(t);

ode1 := diff(A(t), t) = -k*A(t)

>

where the reaction rate k is positive and has units of 1/time.

The solution to this ODE, with initial condition

> ic1 := A(0) = A[0];

ic1 := A(0) = A[0]

>

can be written down on sight, but it is simpler to use Maple to enter the result in a Maple worksheet

> sol1 := dsolve( { ode1, ic1 }, A(t) );

sol1 := A(t) = A[0]*exp(-k*t)

>

While first-order chemical reactions are easy to model and analyze, they are not very common.

>

23.B Second-Order Reactions

An example of a second-order chemical reaction is

CH``[3] Cl + NaOH  -`>`  CH``[3] OH + NaCl

in which one molecule of methyl alcohol and one molecule of sodium hydroxide combine to form one molecule of methyl hydroxide and one molecule of sodium chloride. This reaction proceeds at a rate proportional to the product of the remaining concentrations of methyl chloride and sodium hydroxide.

The general form for a second-order reaction is

A + B  -`>`  C

To model a second-order reaction, let alpha and beta denote the initial amounts of chemicals A and B and denote the amount of chemical C at time t by C(t) . The amount of chemical A remaining at time t is alpha-C(t) ; likewise, beta-C(t) is the remaining amount of chemical B . Hence, the amount of C is modeled by by the ODE

> ode2 := diff( C(t), t ) = k * ( alpha - C(t) ) * ( beta - C(t) );

ode2 := diff(C(t), t) = k*(alpha-C(t))*(beta-C(t))

>

where k is a reaction constant. This first-order ODE is nonlinear (but separable). With initial condition

> ic2 := C(0) = C[0];

ic2 := C(0) = C[0]

>

the amount of chemical C is found to be

> sol2 := dsolve( { ode2, ic2 }, C(t) );

sol2 := C(t) = (-exp(k*t*alpha-k*t*beta)*beta*C[0]+exp(k*t*alpha-k*t*beta)*alpha*beta-alpha*beta+C[0]*alpha)/(-exp(k*t*alpha-k*t*beta)*C[0]+exp(k*t*alpha-k*t*beta)*alpha-beta+C[0])

>

provided alpha <> beta . If alpha = beta , the solution becomes

> simplify( dsolve( { eval( ode2, beta=alpha ), ic2 }, C(t) ) );

C(t) = (t*k*alpha^2-t*k*C[0]*alpha+C[0])/(k*t*alpha-t*k*C[0]+1)

>

Note that it is not possible to obtain this solution simply by substituting beta = alpha into the general solution, as the following failure attests,

> eval( sol2, beta=alpha );

Error, numeric exception: division by zero

>

but is obtained in the limit as beta approaches alpha , as the following limit shows.

> map( limit, sol2, beta=alpha );

C(t) = (t*k*alpha^2-t*k*C[0]*alpha+C[0])/(k*t*alpha-t*k*C[0]+1)

>

23.C Example: 2NO + O``[2] -> 2NO``[2]

Two molecules of nitrous oxide combine with one molecule of oxygen to form two molecules of nitrogen dioxide, according to the chemical equation

2NO + O``[2]  -`>`  2NO``[2]

At room temperature this reaction can be modeled with the IVP

> ode3 := diff( NO2(t), t ) = k * ( alpha - NO2(t) )^2 * ( beta - NO2(t)/2 );

> ic3 := NO2(0) = A;

ode3 := diff(NO2(t), t) = k*(alpha-NO2(t))^2*(beta-1/2*NO2(t))

ic3 := NO2(0) = A

>

where NO2(t) is the concentration of nitrogen dioxide at time t , alpha is the initial concentration of nitrogen oxide, beta is the initial concentration of oxygen, and A is the initial concentration of nitrogen dioxide.

This is not a second-order reaction!  This is, in fact, a third-order reaction.  The general chemical equation

m*A+n*B  -`>`  k*C

where m, n , and k , are positive integers, gives rise to the ODE

dC/dt = K*(alpha-``(m/k)*C)^m*(beta-``(n/z)*C)^y

where alpha and beta are the initial amounts of A and B , respectively.  The reaction is of order m+n because m molecules of substance A must simultaneously meet up with n molecules of substance B to form k molecules of substance C .

The solution to the given IVP

> ode3;
ic3;

diff(NO2(t), t) = k*(alpha-NO2(t))^2*(beta-1/2*NO2(t))

NO2(0) = A

>

with a separable ODE is

> sol3 := dsolve( { ode3, ic3 }, NO2(t) );

sol3 := NO2(t) = exp(RootOf(2*alpha^2*ln(-alpha+A)-2*alpha^2*ln(-2*beta+A)-2*alpha^2*ln(exp(_Z)+2*beta-alpha)-8*beta^2+4*t*k*beta^2*A*exp(_Z)+4*alpha*beta-12*t*k*beta^2*alpha*A+t*k*alpha^2*A*exp(_Z)+6...sol3 := NO2(t) = exp(RootOf(2*alpha^2*ln(-alpha+A)-2*alpha^2*ln(-2*beta+A)-2*alpha^2*ln(exp(_Z)+2*beta-alpha)-8*beta^2+4*t*k*beta^2*A*exp(_Z)+4*alpha*beta-12*t*k*beta^2*alpha*A+t*k*alpha^2*A*exp(_Z)+6...sol3 := NO2(t) = exp(RootOf(2*alpha^2*ln(-alpha+A)-2*alpha^2*ln(-2*beta+A)-2*alpha^2*ln(exp(_Z)+2*beta-alpha)-8*beta^2+4*t*k*beta^2*A*exp(_Z)+4*alpha*beta-12*t*k*beta^2*alpha*A+t*k*alpha^2*A*exp(_Z)+6...sol3 := NO2(t) = exp(RootOf(2*alpha^2*ln(-alpha+A)-2*alpha^2*ln(-2*beta+A)-2*alpha^2*ln(exp(_Z)+2*beta-alpha)-8*beta^2+4*t*k*beta^2*A*exp(_Z)+4*alpha*beta-12*t*k*beta^2*alpha*A+t*k*alpha^2*A*exp(_Z)+6...

>

This "explicit" solution is of almost no use.

Note, however, that the qualitative behavior of solutions can be determined directly from the ODE. There are two equilibria. The equilibrium at NO2 = 2*beta is stable and, because of the quadratic term, NO2 = alpha is semi-stable.

Assuming numeric values for the parameters are available, a numerical solution is easy to obtain.  For example, at a temperature of 25C, the rate constant is 7130 liter``^2 /(mole``^2 second). Assume that the initial concentration of NO is alpha = 0.003 mole/liter, the initial concentration of O``[2] is beta = 0.0021 mole/liter and the initial concentration of NO``[2] is A = 0 mole/liter.  These parameter values are then

> param3 := { k=7.13*10^3, alpha=0.003, beta=0.002, A=0 };

param3 := {k = 7130.00, alpha = 0.3e-2, A = 0, beta = 0.2e-2}

>

Figure 23.1, a plot of the numerical solution to this IVP for the first hour (3600 seconds) of the reaction, is obtained with

> sol3n := dsolve( eval( { ode3, ic3 }, param3 ), NO2(t), numeric ):

> odeplot( sol3n, [t,NO2(t)], 0..3600, view=[0..3600,0..0.003], title="Figure 23.1" );

[Plot]

>

Note that the concentration is increasing to the equilibrium solution NO2 = 0.003.

Slightly more information can be obtained from Figure 23.2 where the semi-stable solution curve (in green) is superimposed on the direction field for this equation.  In addition, the stable equilibrium solution NO2(t) = 0.004 is drawn in blue.

> display(DEplot( eval(ode3,param3), NO2(t), t=0..3600, [[0,0],[0,.004]], NO2=0..0.01, stepsize=25, arrows=medium, dirgrid=[15,15], linecolor=[green,blue]), tickmarks=[4,[.003,.004,.01]], title="Figure 23.2" );

[Plot]

>

It also helps to plot the right-hand side of the (autonomous) ODE as a function of the dependent variable, as done in Figure 23.3.

> plot( eval(rhs(ode3),param3 union {NO2(t)=NO2}), NO2=0.001..0.005, labels=[NO2,"NO2 '     "], title="Figure 23.3" );

[Plot]

>

From Figure 23.3, it is seen that the rate of change of the concentration of NO``[2] is positive when the concentration of NO``[2] is in the interval [0, 0.003) or (0.003, 0.004) and the concentration of NO``[2] decreases only when the concentration exceeds 0.004.  The equilibrium solutions correspond to the intercepts on the horizontal axis.  

For the stable equilibrium at NO2(t) = 0.004, the derivative is positive when the concentration of NO``[2] is less than 0.004, so the concentration will increase towards 0.004.  The derivative is correspondingly negative when the concentration of NO``[2] is greater than 0.004, so the concentration will decrease towards 0.004.

For the semi-stable equilibrium at NO2(t) = 0.003, the derivative is positive when the concentration of NO``[2] is less than 0.003, so the concentration will increase towards 0.003, but the derivative is also positive when the concentration of NO``[2] is greater than 0.003, so the concentration will increase away from 0.003.

>

[Back to ODE Powertool Table of Contents]

>