Spatio-temporal plug flow reactor simulation with method of lines
33 views (last 30 days)
Show older comments
I'm trying out and tinkering a few Matlab codes and techniques in my spare time and encountered this particular problem. Note that I'm fairly new to Matlab.
The code and problem statement that I started from is from the below link, which is transient (time-dependent) simulation of a plug flow reactor with a single component.
What I want is that another component, say B, is specified as product of the reaction. Hence, it has two components involved with this reaction scheme:
A --> B
with additional rate of B formation as -vo*diff(C)./diff(V)-k*CA(2:end).^2
The problem is that the vectorization allows only representation of C(i) as concentration of A along the volume/length of the reactor, instead of type of reacting species. As result, the type of reacting species cannot be specified in the code as matrix of C. The code I wrote so far and the error message is given as follows.
function PFR_t_L_sim_v2_1
% function method_of_lines
clear all; clc; close all
function dCdt = reaction(t,C)
r_1 = k*C(1,2:end).^2;
% we use vectorized operations to define the odes at each node
% point.
dCdt_1 = -vo*diff(C)./diff(V)- r_1;
dCdt_2 = -vo*diff(C)./diff(V)+ r_1;
dCdt = [0; dCdt_1; dCdt_2];
dCdt = dCdt(:);
end
Ca0 = 2; % Entering concentration
Cb0 = 1; % Entering concentration
vo = 2; % volumetric flow rate
volume = 20; % total volume of reactor, spacetime = 10
k = 1; % reaction rate constant
N = 100; % number of points to discretize the reactor volume on
init = zeros(N,1); % Concentration in reactor at t = 0
init(1) = Ca0; % concentration component 1 at entrance
init(2) = Cb0; % concentration component 1 at entrance
V = linspace(0,volume,N)'; % discretized volume elements, in column form
tspan = [0 20];
[t,C] = ode15s(@reaction,tspan,init);
plot(t,C(:,end))
xlabel('time')
ylabel('Concentration at exit')
legend('C_A', 'C_B','Location','N')
end
Error message
Error using odearguments (line 95)
PFR_T_L_SIM_V2_1/REACTION returns a vector of length 1, but the length of initial conditions vector is 100. The vector returned by PFR_T_L_SIM_V2_1/REACTION and
the initial conditions vector must have the same number of elements.
Could anyone help me to sort this out?
Thanks in advance.
0 Comments
Accepted Answer
Torsten
on 13 May 2022
Edited: Torsten
on 13 May 2022
function PFR_t_L_sim_v2_1
% function method_of_lines
clear all; clc; close all
function dCdt = reaction(t,C)
C1 = C(1:N);
C2 = C(N+1:2*N);
dC1dt = zeros(N,1);
dC2dt = zeros(N,1);
dC1dt(2:N) = -vo*(C1(2:N)-C1(1:N-1))./(V(2:N)-V(1:N-1)) - k*C1(2:N).^2;
dC2dt(2:N) = -vo*(C2(2:N)-C2(1:N-1))./(V(2:N)-V(1:N-1)) + k*C1(2:N).^2;
dCdt = [dC1dt;dC2dt];
%r_1 = k*C(1,2:end).^2;
% we use vectorized operations to define the odes at each node
% point.
%dCdt_1 = -vo*diff(C)./diff(V)- r_1;
%dCdt_2 = -vo*diff(C)./diff(V)+ r_1;
%dCdt = [0; dCdt_1; dCdt_2];
%dCdt = dCdt(:);
end
Ca0 = 2; % Entering concentration
Cb0 = 1; % Entering concentration
vo = 2; % volumetric flow rate
volume = 20; % total volume of reactor, spacetime = 10
k = 1; % reaction rate constant
N = 100; % number of points to discretize the reactor volume on
init = zeros(2*N,1); % Concentration in reactor at t = 0
init(1) = Ca0; % concentration component 1 at entrance
init(N+1) = Cb0; % concentration component 1 at entrance
%init = zeros(N,1); % Concentration in reactor at t = 0
%init(1) = Ca0; % concentration component 1 at entrance
%init(2) = Cb0; % concentration component 1 at entrance
V = linspace(0,volume,N)'; % discretized volume elements, in column form
tspan = [0 20];
[t,C] = ode15s(@reaction,tspan,init);
plot(t,[C(:,N),C(:,2*N)])
%plot(t,C(:,end))
xlabel('time')
ylabel('Concentration at exit')
legend('C_A', 'C_B','Location','N')
end
More Answers (0)
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!