Collection and analysis of rate data - 2013

From Introduction to Reactor Design: 3K4
Jump to navigation Jump to search
The printable version is no longer supported and may have rendering errors. Please update your browser bookmarks and please use the default browser print function instead.
Class date(s): 28 February to 06 March
Download video: Link (plays in Google Chrome) [198 M]

Download video: Link (plays in Google Chrome) [370 M]

Download video: Link (plays in Google Chrome) [230 M]

Textbook references

  • F2011: Chapter 7
  • F2006: Chapter 5

Suggested problems

F2011 F2006
Problem 7-7 (a) Problem 5-6 (a)
Problem 7-8 (a) Problem 5-7 (a)
Problem 7-15 Not in this edition

Class materials

28 February 2013 (07C)

04 March 2013 (08A)

06 March 2013 (08B-1)

Source code: R

R code for the example covered in class:

t <- c(0, 50, 100, 150, 200, 250)    # [minutes]
# or you can write t <- seq(0, 250, 50)

CA <- c(154, 114, 87, 76, 70, 58)  # [mol/m^3]

# Plot the data first. See the expected trend. It it decreasing nonlinearly, so it is either first or second order (but not zeroth order)
plot(t, CA)
grid()

# Assume it is a first order, rA = -k CA
# So plot log(CA0/CA) against t
plot(t, log(CA[1]/CA))
grid()

# Does not quite look like a straight line, but let us fit a least-squares model anyway

mod.first <- lm( log(CA[1] / CA) ~ t + 0)  # force intercept to be zero
summary(mod.first)

#Call:
#lm(formula = log(CA[1]/CA) ~ t + 0)
#
#Residuals:
#         1          2          3          4          5          6 
# 1.874e-16  8.988e-02  1.493e-01  7.361e-02 -5.502e-02 -7.784e-02 
#
#Coefficients:
#   Estimate Std. Error t value Pr(>|t|)    
#t 0.0042174  0.0002555   16.51 1.49e-05 ***
#---
#Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
#
#Residual standard error: 0.09474 on 5 degrees of freedom
#Multiple R-squared: 0.982,	Adjusted R-squared: 0.9784 
#F-statistic: 272.5 on 1 and 5 DF,  p-value: 1.489e-05 

# Now assume a second order system, rA = -k CA^2
# So plot 1/CA against t
plot(t, 1/CA)
grid()
abline(mod.first)

# Definitely more linear. Let us fit the least squares model to check:
mod.second <- lm( I(1/CA) ~ t)
summary(mod.second)
#Call:
#lm(formula = I(1/CA) ~ t)
#
#Residuals:
#         1          2          3          4          5          6 
#-2.751e-04 -5.219e-05  6.146e-04  2.227e-04 -7.051e-04  1.951e-04 
#
#Coefficients:
#             Estimate Std. Error t value Pr(>|t|)    
#(Intercept) 6.769e-03  3.692e-04   18.33 5.21e-05 ***
#t           4.111e-05  2.439e-06   16.86 7.26e-05 ***
#---
#Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
#
#Residual standard error: 0.0005101 on 4 degrees of freedom
#Multiple R-squared: 0.9861,	Adjusted R-squared: 0.9826 
#F-statistic: 284.2 on 1 and 4 DF,  p-value: 7.259e-05

<rst> <rst-options: 'toc' = False/> <rst-options: 'reset-figures' = False/>

  1. . :math:`\displaystyle \int{ \frac{1}{(1-X)^2} \,dX} =`

</rst>

So we got a higher \(R^2\) on the second order model. The reaction rate can be considered second order. What is the reaction rate constant, \(k_A\) for this system?

Source code: Using MATLAB or Octave

t = [0, 50, 100, 150, 200, 250];    % [minutes]
CA = [154, 114, 87, 76, 70, 58];    % [mol/m^3]

% Plot the data first. See the expected trend. It it decreasing nonlinearly, so 
% it is either first or second order (but not zeroth order)
figure
plot(t, CA, 'o')
grid on

% Some template code for the linear regression
x = t;
y = log(CA(1)./CA);
n = numel(y);

% Assume it is a first order, rA = -k CA; so plot log(CA0/CA) against t
figure
plot(x, y, 'o')
grid on

% Does not quite look like a straighTSt line, but let us fit a least-squares model anyway
% Fit the model y = b x + 0 (i.e. force the intercept to zero) 
X = [x(:)]; 
y = y(:);              % make sure it is a column vector       
b = inv(X'*X)*X'*y;    % recall this from the 3E4 course 
%b = regress(y, X);     % only if you have the Statistics Toolbox in MATLAB
 
% Plot the data along with the fitted line:
figure
plot(x, y, 'o')
hold('on')
grid('on')
plot(x, X*b, 'r')
xlabel('x values')
ylabel('y values')
legend({'Original data', 'Fitted line'}, 'Location', 'Best')
 
% Additional calculations
predicted = X*b;
resids = y - predicted;            % resids = e = y - Xb
RSS = resids' * resids;            % residual sum of squares
TSS = sum((y - mean(y)).^2);       % total sum of squares
R2_first_order = 1 - RSS/TSS


% Now try a second order model: -rA = k CA^2; so plot 1/CA against t
figure
plot(t, 1./CA, 'o')
grid('on')

% Definitely more linear. Let us fit the least squares model,
% 1/CA = 1/CA0 + kt, or y = b0 + b1 x, where y = 1/CA and x = t, and there is an intercept
% Some template code for the linear regression
x = t;
y = 1./CA;
n = numel(y);
X = [ones(n,1) x(:)]; 
y = y(:);              % make sure it is a column vector       
b = inv(X'*X)*X'*y;    % recall this from the 3E4 course 
% b = regress(y, X);     % only if you have the Statistics Toolbox in MATLAB

% Plot the data along with the fitted line:
figure
plot(x, y, 'o')
hold('on')
grid('on')
plot(x, X*b, 'r')
xlabel('x values')
ylabel('y values')
legend({'Original data', 'Fitted line'}, 'Location', 'Best')
 
% Additional calculations
predicted = X*b;
resids = y - predicted;            % resids = e = y - Xb
RSS = resids' * resids;            % residual sum of squares
TSS = sum((y - mean(y)).^2);       % total sum of squares
R2_second_order = 1 - RSS/TSS