## Reaction Engineering

## Reaction Engineering

(OP)

Hi

I am new with Matlab and trying to solve this problem.

The problem is a mutli-reaction system of the form

dc1/dt=kc1+kc2

dc2/dt=kc2-kc3

... etc.

where k=reaction constant

Here is a code I wrote in Matlab:

function dc=bio(t,c)

T=298

E=13145 %constant

R=1.98 % constant

A=1.4e5 %constant

k=A*e^(-E/RT) %(note K is a function of T)

dc=[dc1;dc2]

Then I made another m.file

t=0

tf=20

initial=[1 0]

[t,c]=ode45(@bio,[t0 tf],initial);

So far everything is fine and the model works.

Now I like to change the T (over a wide range)and

run the model again. I tried a few things but none worked. Can someone please tell me how I can do that?

I am new with Matlab and trying to solve this problem.

The problem is a mutli-reaction system of the form

dc1/dt=kc1+kc2

dc2/dt=kc2-kc3

... etc.

where k=reaction constant

Here is a code I wrote in Matlab:

function dc=bio(t,c)

T=298

E=13145 %constant

R=1.98 % constant

A=1.4e5 %constant

k=A*e^(-E/RT) %(note K is a function of T)

dc=[dc1;dc2]

Then I made another m.file

t=0

tf=20

initial=[1 0]

[t,c]=ode45(@bio,[t0 tf],initial);

So far everything is fine and the model works.

Now I like to change the T (over a wide range)and

run the model again. I tried a few things but none worked. Can someone please tell me how I can do that?

## RE: Reaction Engineering

=====================================

Eng-tips forums: The best place on the web for engineering discussions.

## RE: Reaction Engineering

The problem was when I put the loop in the the "solver" file (I mean the one which I call the ode45), it takes only the first value of T and repeats the calculation. As I said I am new with Matlab so I may be making a stupid mistake!

## RE: Reaction Engineering

T=298

in your bio file (derivative file).

You have to do it in your main file. Something like:

global T

t=0

tf=20

initial=[1 0]

T=298

for index = 1:100

T=T+0.1*index

[t,c]=ode45(@bio,[t0 tf],initial);

end

Then your bio file doesn't change T but instead pulls it in as a global variable. Bio.m/ Something lik

dc=bio(t,c)

global T

E=13145 %constant

R=1.98 % constant

A=1.4e5 %constant

k=A*e^(-E/RT) %(note K is a function of T)

dc=[dc1;dc2]

Note it looks like you may have left out a few steps in defining dc1 and dc2 but I assume they can be calculated in the dc=bio() derivative procedure as a funciton of T.

=====================================

Eng-tips forums: The best place on the web for engineering discussions.

## RE: Reaction Engineering

% file Main.m

global T

t0=0

tf=20

initial=[1 1]

T=298

for index = 1:3

T=T+10*index

[t,c]=ode45('bio',[t0 tf],initial);

plot(t,c)

pause

end

%file bio.m

function dc=bio(t,c)

global T

E=13145 %constant

R=1.98 % constant

A=1.e5 %constant

k=A*exp(-E/R/T) %(note K is a function of T)

d1=k*c(1)+k*c(2)

d2= k*c(1)-k*c(2)

dc=[d1;d2]

=====================================

Eng-tips forums: The best place on the web for engineering discussions.

## RE: Reaction Engineering

Thanks again. It works!

## RE: Reaction Engineering

plot(t,c(:,1))

title(['C1, T=' num2str(T)])

pause

plot(t,c(:,2))

title(['C2, T=' num2str(T)])

Eng-tips forums: The best place on the web for engineering discussions.

## RE: Reaction Engineering

Eng-tips forums: The best place on the web for engineering discussions.