Numerical Integration                                                                 (last updated  9/9/99)

The information in this tutorial is located in the Section 16.4 of the MATLAB manual,  The Student Edition of MATLAB, Version 5, The MATH WORKS Inc. Prentice-Hall, Englewood Cliffs, NJ 07632, 1997.

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

This tutorial contains the following sections;

Introduction

Integration Using a Set of Discrete Data

Example Problem - Using a Set of Discrete Data

Integration When a Function is Available

Example Problem - When a Function is Available
 
=============================================================================================================

Introduction

Integration arises in a number of situations in chemical engineering. MATLAB has the capability of performing integrations both symbolically and numerically. Symbolically either definite or indefinite integration can be performed. Symbolic integration is not addressed in this tutorial. See Section 6.4 of the MATLAB manual for information regarding symbolic integration. Numerical integration is addressed in Section 5.18.4.

We are interested in numerically determining the area under the curve over a specific given interval. That is we want to determine an estimate for

(1) Integration using a set of discrete data.
(2) Integration when a function is available.

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

Integration Using a Set of Discrete Data    Section 16.4

When a set of discrete data is available, MATLAB will perform the numerical integration using the MATAB command trapz, which uses the trapezoidal rule. The data must be loaded into arrays. Given arrays x and y, the syntax of the command is ( >> I = trapz(x,y) ). The integration will b performed over the array x.

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

Example Problem - Using a Set of Discrete Data

Estimation of Total Mass From Flow and Concentration Data
 

Here is the m-file which solves the example problem.
______________________________________________________________________

%    Darin Ridgway
%   ChE XXX
%   July 11, 1998

%   Enter the data
time = [0:120:1440];                %   time in minutes
Q = [4.2 3.5 2.8 5.6 8.4 12.6 12.4 9.7 14.6 10.8 7.0 5.0 3.4];        %   Q is in m^3/min
C = [12.43 11.80 11.27 15.82 18.22 11.97 8.48 6.79 18.65 22.58 16.42 6.15 11.86];      %    C is in mg/m^3

%    Calculate the mass flow rate in mg/min
qc = Q.*C;

%   Call the trapz routine to integrate the data
mass = trapz(time,qc);       % This is the total mass in miiligrams

%   Create output
disp('    Darin Ridgway')
disp('    ChE XXX')
disp('    July 11, 1998')

fprintf('\n\n   The total mass is %6.2e mg \n', mass)
____________________________________________________________________________

Here is the response in the Command Window
____________________________________________________________________________

    Darin Ridgway
    ChE XXX
    July 11, 1998

The total mass is 1.6133e+005  mg
_____________________________________________________________

Notes:

(a) The more data you have to use the better your estimate will be. This is especially true where the integrand is changing quickly.

(b) If your limits of integration are not directly in the data set, you must create them from the data set by interpolating or extrapolating. For the data set above this would happen if I had desired the total mass between 3 P.M. and 9 P.M.

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

Integration When a Function is Available    Section 16.4

When a function is available instead of a set of discrete data a different MATLAB command is utilized. There are two choices;

(a) Use the function to create a distinct array of points, then use the trapz command.
(b) Use the MATLAB command quad or quad8. The only difference between them is quad8 is more rigorous in using higher order approximations.

The syntax is quad('function', lower limit, upper limit), where 'function' is a function m-file name containing the function to be integrated, and lower and upper limt are the limits of inegration.

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

Example Problem - When a Function is Available

Change in Enthalpy for Hydrogen Cyanide.

Calculate the change in enthalpy for 300 mol/h of hydrogen cyanide being heated from 100 to 300°C.

The heat capacity in J/mol is given as (35.3 + 2.908e-2*T + 1.092e-5*T^2) where T is in degrees C.

(Heat capacity data from Felder and Rousseau, Elementary Principles of Chemical Processes, 2nd ed., Wiley and Sons Inc., 1986.)
_________________________________________________________________________________

Here is the m-file which solves the example problem.
______________________________________________________________________

%    Darin Ridgway
%   ChE XXX
%   July 11, 1998

%    Call the quad routine to solve the system.  The function is in file HCN.m
delH = quad('HCN',100,300);

%   Create output
disp('    Darin Ridgway')
disp('    ChE XXX')
disp('    July 11, 1998')

fprintf('\n\n   The change in enthalpy is %6.2e J/h \n', delH) ________________________________________________________________________________________

Here is the function file
________________________________________________________________________________________

%    Darin Ridgway
%   ChE XXX
%   July 11, 1998

function delH = HCN(T)

specH = 35.3 + 0.02908*T + 1.092e-5*T.^2;    %   The units on specH are J/mol
delH = 300 * specH;    % The units on delH are J/h
______________________________________________________________________________________

Here is the response in the Command Window
______________________________________________________________________________________

%    Darin Ridgway
%   ChE XXX
%   July 11, 1998

The change in enthalpy is 2.4954e+006 J/h