Hello, I know this is a while after you posted the question, but...
I have looked through the code and tried to comment where possible, my interpretation of the code (see the end of this answer).
First and foremost, there are a few spelling errors in the actual code (e.g. 'lable' not 'label') - it is worth noting that MATLAB will not run if functions aren't spelt correctly. I have corrected these. I am also coming into this with no idea what any of the equations you have in your MATLAB file are - I have tried to guess what some mean, but you will have a better idea than I will.
Another important thing is to ensure that in your main script file (the first one you provided), the function (you called this 'Code for 10.6 B') you are calling in the ode45 line is spelt exactly how the file is saved. I have changed the name to make it work - but as a rule, the function name must be the same as the filename under which it is saved - so 'CSTR_Gas' must be the function name and so must the filename.
Anyway, below is the code - I have written comments and tidied it up a bit, keeping tidy code is a better way to make it easy to read and understandable.
Main script (call this whatever you want)
clear all % I recommend adding this to the start of all your main scripts
close all % to help clear any variables, close any figures and clear the
clc % command window - helps keep things tidy.
global Cv P_D MWA R T rowo Fo CAo k1 k2 V P F y CB
%%Defining Variables
% Variable % % Description %
Cv = 15.0; % Valve Coefficient / some kind of friction factor ??
P_D = 10.0; % Downstream pressure ??
MWA = 50.0; % Molecular mass of component A
MWB = 50.0; % Molecular mass of component A
R = 0.082057; % Molar gas constant
T = 350.0; % Temperature
Po = 30.0; % Initial pressure ??
rowo= (Po*MWA)/(R*T); % Initial density of gas ??
Fo = 10.0; % Inlet flow rate ??
CAo = rowo/MWA; % Initial inlet concentration ??
k1 = 0.5; % Rate constant for forward reaction ??
k2 = 0.05; % Rate constant for backward reaction ??
V = 1000.0; % Reactor volume
tf = 30.0; % Final time?
% Setting the Relative/Absolute Tolerances (Rel/AbsTol) to 1 * 10 ^{-6}
% Very basically, this is how accurate you want your solution:
options = odeset('RelTol',1e-6,'AbsTol',[1e-6 1e-6]);
% This is the ode45 solver (an iterative solver based on 4th and 5th order
% Runge-Kutta numerical methods). It is calling the 'CSTR_Gas' function file
% and integrating the solution between 0 and tf. [rowo CA0] are the initial conditions of the ODE
%'Options' is passing the above defined tolerances to the ode45 solver:
[t,var] = ode45(@CSTR_Gas, [0 tf], [rowo CAo], options);
row = var(:,1); % I assume row is density, it's spelt 'rho'.
CA = var(:,2); % Concentration of A.
% These two lines extract the density and concentration
% of A from the 'var' array created in the function
% file.
P=row*R*T/MWA; % Pressure at each time point
F=Cv*sqrt((P-P_D)./row); % Flow rate ??
y=CA*R*T./P; % Some kind of partial pressure ratio or mole fraction ??
XA=(CAo-CA)/CAo; % Convesion of A
CB=(P/(R*T)).*(1.0-y); % Concentration of B, uses the ratio y above to calculate CB
%%Plotting Solutions
subplot(411)
% This would make more sense to write: subplot(4,1,1) - creates a 4 x 1 grid,
% places the subsequent plot in the 1st position.
plot (t,CA, 'g') % Plots the concentration against time in the first plot position.
ylabel(['\bf C_A','mol/L']); % labels the y axis with the quoted text.
% Same as above - specifying where to place the plot in the subplot - 4 x 1
% grid, 2nd position etc...
subplot (412) %advance to plot #2
plot (t,XA, 'r') %plot #2
ylabel('\bf Conversion, X_A');
% Same as above - specifying where to place the plot in the subplot - 4 x 1
% grid, 3rd position etc...
subplot(413) % advance to plot #3
plot (t,y, 'b') %plot #3
ylabel('\bf y');
% Same as above - specifying where to place the plot in the subplot - 4 x 1
% grid, 4th position etc...
subplot(414) % advance to plot #4
plot (t,CB, 'c') %plot #4
xlabel('\bf time unit'), ylabel('\bf C_B, mol/L');
Function File (call this "CSTR_Gas" - letter-for-letter)
function [deriv] = CSTR_Gas(t,var)
global Cv P_D MWA R T rowo Fo CAo k1 k2 V P F y CB
% This defines the function CSTR_Gas which has two input arguments which
% are:
% t and var and the output is the derivative of var, deriv.
%Since we have two ode's, we initialize dy as 2 by 1 matrix.
row=var(1); % I am confused as to what 'row' is referring to, perhaps 'rho' - i.e. density?
CA=var(2); % But these two lines allocate spaces to put the solutions
deriv = zeros(2,1); % Creates an empty vector to store the ODE solutions in
P=row*R*T/MWA; %row=P*MW/(RT); % These lines just calculate
F=Cv*sqrt((P-P_D)/row); %flow through a valve % each variable (P, F, y, CB) at each
y=CA*R*T/P; %yA=PA/P; yA=CA*RT/P; % iteration of the ODE solver
CB=(P/(R*T))*(1-y); %CB+P*(1-yA)/(RT); %
deriv(1) = (rowo*Fo/V)-(row*F/V); %
deriv(2) =(Fo*CAo/V)-(F*CA/V)-k1*CA+k2*CB; %
Now, open up MATLAB, start a new script file, copy and paste each of these files below into two separate scripts. Save the main script as whatever you like and save the function file as "CSTR_Gas", make sure they are in the same folder. Once they are saved, click on the main script file and click on the "Run" button - if you click "Run" on the function file it won't work. It should produce the attached file.
Hope that helps,
Ash
