Optimization of Plug Flow Reactor by Varying Temperature and Pressure

95 views (last 30 days)
I'm trying to optimize a complex reaction system to have the largest output of desired product C by varying the temperature and pressure of the reaction in an adiabatic PFR. I've created the design equation of the reaction:
Function dYdV = Designequation(v,y,T,P)
Fa = y(1);
Fb = y(2);
Fc = y(3);
Fd = y(4);
% Explicit equations
T = T(1);
P = P(1);
Cto = P / (8.314 * 10^-5) / T;
E1 = 15000;
E2 = 17500;
Ft = Fa + Fb + Fc + Fd;
Ca = Cto * Fa / Ft;
Cb = Cto * Fb/Ft;
Cc = Cto * Fc/Ft;
Cd = Cto * Fd/Ft;
k1 = 0.075 * exp(E1 /1.987 * (1/300 - (1/T)));
k2 = 0.0015 * exp(E2 / 1.987 * (1 / 300 - (1 / T)));
Fao = 50;
ra = -k1*Ca*Cb;
rb = -k2*Cb*Cc;
% Differential equations
dFadV = ra;
dFbdV = (2*ra)+rb;
dFcdV = rb-ra;
dFddV = -rb;
dYdV = [dFadV; dFbdV; dFcdV; dFddV];
The design equation is solved for volume and the initial molar flow rates of A and B using ode113 with an initial temperature of 300K and 1 atm:
Vspan = [0 2]; % Range for the volume of the reactor
y0 = [50; 25; 0; 0]; % Initial values for the dependent variables i.e Fa,Fb,and Fc, Fd
T = 300;
P = 1;
[v, y] =ode113(@(v,y) Designequation(v,y,T,P),Vspan,y0);
dFadT = y(:,1);
dFbdT =y(:,2);
dFcdT =y(:,3);
dFddT = y(:,4);
plot(v,dFcdT); hold on
Now I want to find the largest value of the molar flow rate of C (dFcdT), by varying the temperature and pressure, but I haven't found a way that I can implement that based on searches of the MatLab forum.
I apologize if this has been answered in other forum posts, but I'm desperate.
  2 Comments
Torsten
Torsten on 16 Apr 2019
So you want to maximize dFcdT(end) ? Or integral_{v=0}^{v=2} dFcdT dv ? Or something else ?
Sean Powers
Sean Powers on 16 Apr 2019
I want to maximice dFcdT(end) while keeping the vspan 0:2 and y0 the same.

Sign in to comment.

Accepted Answer

Torsten
Torsten on 16 Apr 2019
Edited: Torsten on 16 Apr 2019
Although I guess that there are bounds on P and T, here is a code that will give you a start.
In case you want to constrain P and T, use "fmincon" instead of "fminsearch".
function main
format long
T = 300;
P = 1;
x0 = [T P];
options = optimset('MaxFunEvals',100000,'MaxIter',100000)
x = fminsearch(@fun,x0,options)
end
function obj = fun(x)
T = x(1);
P = x(2);
Vspan = [0 2]; % Range for the volume of the reactor
y0 = [50; 25; 0; 0]; % Initial values for the dependent variables i.e Fa,Fb,and Fc, Fd
[v, y] =ode15s(@(v,y) Designequation(v,y,T,P),Vspan,y0);
obj = -y(end,3)
end
function dYdV = Designequation(v,y,T,P)
Fa = y(1);
Fb = y(2);
Fc = y(3);
Fd = y(4);
% Explicit equations
Cto = P / (8.314 * 10^-5) / T;
E1 = 15000;
E2 = 17500;
Ft = Fa + Fb + Fc + Fd;
Ca = Cto * Fa / Ft;
Cb = Cto * Fb/Ft;
Cc = Cto * Fc/Ft;
Cd = Cto * Fd/Ft;
k1 = 0.075 * exp(E1 /1.987 * (1/300 - (1/T)));
k2 = 0.0015 * exp(E2 / 1.987 * (1 / 300 - (1 / T)));
Fao = 50;
ra = -k1*Ca*Cb;
rb = -k2*Cb*Cc;
% Differential equations
dFadV = ra;
dFbdV = (2*ra)+rb;
dFcdV = rb-ra;
dFddV = -rb;
dYdV = [dFadV; dFbdV; dFcdV; dFddV];
end
  9 Comments
Torsten
Torsten on 23 Apr 2019
@Sean Powers:
And did the exaustive search lead to a different result than the usual optimization with "fminsearch" ?

Sign in to comment.

More Answers (0)

Community Treasure Hunt

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

Start Hunting!