Department of Chemical Engineering
MATLAB Tutorial
 

Solving a System of Linear Equations                                         (last updated  9/9/99)

The information in this tutorial is located in the MATLAB manual. Any page or section numbers refer to the following:

The Student Edition of MATLAB, Version 5, The MATH WORKS Inc. Prentice-Hall, 1997.

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

This tutorial contains the following sections;

Introduction

Using Matrix Algebra Commands   -  (the preferred method)

Using the Symbolic Algebra solve Command

Using the Numerical Equation Solver, ie. the fsolve Command

Example Problem - Linear System of Equations

Example Problem Solution Using Matrix Algebra Commands

Example Problem Solution Using the Symbolic Algebra solve Command

Example Problem Solution Using the Numerical fsolve Command
 
=============================================================================================================

Introduction

There are a number of common situations in chemical engineering where systems of linear equations appear. Any system of 'n' linear equations in 'n' unknowns can be expressed in the general form:

where xj is the jth variable,
aij is the constant coefficient on the jth variable in the ith equation, and
bj is the constant term.
 

There are at least three ways in MATLAB of solving these systems;

(1) using matrix algebra commands,   (the preferred method)
(2) using the symbolic algebra solve command,
(3) using the numerical equation solver, ie. the fsolve command.

All three will be explained and demonstrated. Remember, as always you should work in an m-file.

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

Using Matrix Algebra Commands      This is the preferred method

The system of equations given above can be expressed in the following matrix form.

In short-hand notation this is listed as Ax = b.

Remember that the rules of matrix multiplication dictate that the 'A' matrix be consistent with the 'x' and 'b' vectors. The first column of 'A' contains the coefficients for x1. The first row of 'A' and element of 'b' contain the coefficients for the first equation.

To determine the variables contained in the column vector 'x', complete the following steps.

(a) Create the coefficient matrix 'A'. Remember the semi-colon notation for moving to the next row when entering data into a matrix. Also remember to include zeroes where an equation doesn't contain a variable.

(b) Create the right-hand-side column vector 'b' containing the constant terms from the equation. This must be a column vector, not a row.

(c) Calculate the values in the 'x' vector by left dividing 'b' by 'A', by typing x = A\b

Note that this is different from x = b/A when matrices are involved.

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

Using the Symbolic Algebra solve Command

You need to understand symbolic algebra if you are going to use this command.  If you do not understand what it means to do symbolic manipulation , it is not recommended that you use this method.

MATLAB is able to solve equations symbolically. For example the equation 3x + 4y = 6 can be solved for y. First enter the equation symbolically using single quotes, then execute the solve command. You must tell MATLAB which variable to solve for. The syntax is:

solve('equation(s) to be solved','variable(s) to solve for').
 

Multiple equations and variables are separated by commas.

e = '3*x + 4*y = 6'
e =
    3*x + 4*y = 6

y = solve(e,'y')
y =
     -3/4*x+3/2

Note that quotes aren't needed if the equation has been set equal to a symbolic variable.
 

Alternately, the equation can be entered into the solve command directly:

y = solve('3*x +4*y = 6','y')
y =
    -3/4*x+3/2
 

An example with multiple equations with a symbolic solution is given here.

f1 = 'x - 9*y = 180'
f1 =
    x - 9*y = 180

f2 = '32*x + 6*y = 300'
f2 =
    32*x + 6*y = 300

solve(f1,f2,'x,y')
ans =
     y = -130/7, x = 90/7
 
These are symbolic answers. To convert to numeric answers use the eval command.

eval(ans)
y =
    -18.5714
x =
    12.8571

If MATLAB cannot determine a unique symbolic solution it will try to solve it numerically.

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

Using the Numerical Equation Solver, ie. The fsolve Command

fsolve is the numerical, non-linear equation solver built into MATLAB. In almost every case it will successfully solve a linear system of equations.

Note: fsolve is part of the Optimization Toolbox and therefore can only be used on a machine with the professional version. It is not part of the Student Edition.
 

To use fsolve, the following steps must be executed.

(a) Create a function m-file containg the equations to be solved. Remember, the equations must be in the form of f(x) = 0.

(b) Make an initial guess as to the roots of the equation. The better your initial guess the more likely that the roots will be located.

(c) From the MATLAB prompt, execute the fsolve command. The syntax is;

           q = fsolve('function m-file name',initial guess).

For further information on using the fsolve command, see the tutorial on solving non-linear equations.

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

Example Problem - Linear System of Equations

The following series of reversible reactions occurs;

                    A <=====> B

                    A <=====> D

                    B <=====> C

                    B <=====> D

                    C <=====> D

Assume that the reactions are carried out in a batch reactor at constant temperature and pressure, and that the reactor initially contains only component A at a concentration of 1 mole/liter. Given the first-order rate constants for each reaction, determine the equilibrium concentration of each component.

Given an example set of rate constants, and applying a mass balance on three of the components and the constraint that the concentration of all four species must sum to one, results in the following system of linear equations.

                    CA + CB + CC + CD = 1

                    - 0.3CA + 0.02CB + 0.05CD = 0

                    0.1CA - 0.82CB + 0.1CC + 0.1CD = 0

                    0.5CB - 0.11CC + 0.1CD = 0

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

Solution Using Matrix Algebra Commands        The Preferred Method

The objective is to solve the following matrix equation for the vector C.

where C(1) through C(4) represent CA through CD.
________________________________________________________________________

Here is the m-file used to solve the problem
________________________________________________________________________

%   Darin Ridgway
%   Chemical Engineering
%   Last Updated  July 7, 1998

%   Create the coefficient matrix, 'A'.
A = [  1       1       1       1
        -0.3  0.02     0    0.05
         0.1 -0.82   0.1    0.1
          0     0.5   -0.11  0.1  ];

%   Create the right hand side of knowns.  Be sure it is a column
b = [ 1    0    0    0 ]' ;

%   Use matrix division to solve the system
C = A\b;

%   create the output
disp('    Darin Ridgway')
disp('   Chemical Engineering')
disp('   July 8, 1998')
disp( '  Example Problem')
fprintf('\n\n\n  The concentrations in mole/lit are\n\n')
fprintf(' Species     Conc\n')
fprintf('     A            %5.2f \n', C(1))
fprintf( '    B            %5.2f \n', C(2))
fprintf( '    C            %5.2f \n', C(3))
fprintf( '    D            %5.2f \n', C(4))

________________________________________________________________________

Here is the response in the Command Window
________________________________________________________________________

   Darin Ridgway
   Chemical Engineering
   July 8, 1998
   Example Problem

  The concentrations in mole/lit are

 Species     Conc
     A           0.04
     B           0.11
     C           0.66
     D           0.19
___________________________________________________________________________
 
(d) A quick check is to multiply the coefficient matrix 'A' by the calculated vector'C'. This should result in the 'b' vector being returned.

>> check = A*C
check =
   1.0000
   0.0000
   0.0000
   0.0000

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

Solution Using the Symbolic Algebra solve Command

(a) Enter the four equations as symbolic equations by giving the equation a name (in this case e1, e2, e3, and e4) and using the single quotes around the algebraic expression.
______________________________________________________________

Here is the m-file used to solve the example peroblem using the solve command
______________________________________________________________

%   Darin Ridgway
%   Chemical Engineering
%   Last updated July 8, 1998

%   Enter the equations.  Be sure to put the equations in single quotes
e1 = 'CA + CB + CC + CD = 1' ;
e2 = '-0.3*CA + 0.02*CB + 0.05*CD = 0' ;
e3 = '0.1*CA - 0.82*CB + 0.1*CC + 0.1*CD = 0' ;
e4 = '0.5*CB - 0.11*CC + 0.1*CD = 0' ;
 
 %   Use the solve command to find a solution
solve(e1,e2,e3,e4,'CA,CB,CC,CD') ;

fprintf('\n  All concentrations are in moles/liter \n\n')
______________________________________________________________

Here is the response in the command Window.  This is an example of solve finding a numerical solution
______________________________________________________________

CC = .6648682957791177, CD = .1878768644874642, CB = .1086956521739131,

CA = .3855918755950492e-1

All concentrations are in moles/liter

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

Solution Using the Numerical Equation Solver, ie. The fsolve Command

_______________________________________________________________

Here is the function file that contains the equations to be solved.  They must be in the form f(C) = 0.
_______________________________________________________________

%   Darin Ridgway
%   Chemical Engineering
%   Last Updated  -  July 8, 1998

% This is function m-file reaction.m

function z = reaction(C)

%     C(1) = the concentration of species A
%     C(2) = the concentration of species B
%     C(3) = the concentration of species C
%     C(4) = the concentration of species D

z(1) = C(1) + C(2) + C(3) + C(4) - 1;
z(2) = -0.3*C(1) + 0.02*C(2) +0.05*C(4);
z(3) = 0.1*C(1) - 0.82*C(2) + 0.1*C(3) + 0.1*C(4);
z(4) = 0.5*C(2) - 0.11*C(3) + 0.1*C(4);
_________________________________________________________________

Here is the m-file that sets the initial guess and calls fsolve
_________________________________________________________________

%   Darin Ridgway
%   Chemical Engineering
%   Last Updated  -  July 8, 1998

%   Make an initial guess of the roots
C0 = [ 0.25 0.25 0.25 0.25 ] ;

%   Execute the fsolve command, giving it the m-file name of equations and the initial guess.
C = fsolve('reaction',C0);

%   Generate the output

disp('    Darin Ridgway')
disp('   Chemical Engineering')
disp('   July 8, 1998')
disp( '  Example Problem')
fprintf('\n\n\n  The concentrations in mole/lit are\n\n')
fprintf(' Species     Conc\n')
fprintf('     A            %5.2f \n', C(1))
fprintf( '    B            %5.2f \n', C(2))
fprintf( '    C            %5.2f \n', C(3))
fprintf( '    D            %5.2f \n', C(4))
____________________________________________________________________

Here is the response in the Command Window
____________________________________________________________________

   Darin Ridgway
   Chemical Engineering
   July 8, 1998
   Example Problem

The concentrations in mole/lit are

 Species     Conc
     A           0.04
     B           0.11
     C           0.66
     D           0.19
_____________________________________________________________________