Nonlinear Algebraic Systems 1 Iterative solution methods 2

Nonlinear Algebraic Systems 1. Iterative solution methods 2. Fixed-point iteration 3. Newton-Raphson method 4. Secant method 5. Matlab tutorial 6. Matlab exercise

Single Nonlinear Equation l Nonlinear algebraic equation: f(x) = 0 » Analytical solution rarely possible » Need numerical techniques » Multiples solutions may exist l Iterative solution » » » l Start with an initial guess x 0 Algorithm generates x 1 from x 0 Repeat to generate sequence x 0, x 1, x 2, … Assume sequence convergences to solution Terminate algorithm at iteration N when: Many iterative algorithms available » Fixed-point iteration, Newton-Raphson method, secant method

Fixed-Point Iteration l Formulation of iterative equation » A solution of x = g(x) is called a fixed point l Convergence » The iterative process is convergent if the sequence x 0, x 1, x 2, … converges: » Let x = g(x) have a solution x = s and assume that g(x) has a continuous first-ord » A function satisfying theorem is called a contraction mapping: » K determines the rate of convergence

Newton-Raphson Method l Iterative equation derived from first-order Taylor series ex l Algorithm » Input data: f(x), df(x)/dx, x 0, tolerance (d), maximum number of iter » Given xn, compute xn+1 as: » Continue until | xn+1 -xn| < d|xn| or n = N

Convergence of the Newton-Raphson Method l Order » Provides a measure of convergence rate » Newton-Raphson method is second-order l l Assume f(x) is three times differentiable, its first- and second-order derivatives are non-zero at the solution x = s & x 0 is sufficiently close to s, then the Newton method is second-order & exhibit quadratic converge to s Caveats » The method can converge slowly or even diverge for poorly chosen x 0 » The solution obtained can depend on x 0 » The method fails if the first-order derivative becomes zero (singularity)

Secant Method l Motivation » Evaluation of df/dx may be computationally expensive » Want efficient, derivative-free method l Derivative approximation l Secant algorithm l Convergence » Superlinear: » Similar to Newton-Raphson (m = 2)

Matlab Tutorial l l Solution of nonlinear algebraic equations with Matlab FZERO – scalar nonlinear zero finding » Matlab function for solving a single nonlinear algebraic equation » Finds the root of a continuous function of one variable » Syntax: x = fzero(‘fun’, xo) – ‘fun’ is the name of the user provided Matlab m-file function (fun. m) that evaluates & returns the LHS of f(x) = 0. – xo is an initial guess for the solution of f(x) = 0. » Algorithm uses a combination of bisection, secant, and inverse quadratic interpolation methods.

Matlab Tutorial cont. l l l Solution of a single nonlinear algebraic equation: Write Matlab m-file function, fun. m: Call fzero from the Matlab command line to find the solution: >> xo = 0; >> fzero('fun', xo) ans = 0. 5376 l Different initial guesses, xo, can give different solutions: >> fzero('fun', 1) ans = 1. 2694 >> fzero('fun', 4) ans = 3. 4015

Nonisothermal Chemical Reactor l Reaction: A B l Assumptions » Pure A in feed » Perfect mixing » Negligible heat losses » Constant properties (r, Cp, DH, U) » Constant cooling jacket temperature (Tj) l Constitutive relations » Reaction rate/volume: r = kc. A = k 0 exp(-E/RT)c. A » Heat transfer rate: Q = UA(Tj-T)

Model Formulation l Mass balance l Component balance l Energy balance

Matlab Exercise l Steady-state model l Parameter values » » » l k 0 = 3. 493 x 107 h-1, E = 11843 kcal/kmol (-DH) = 5960 kcal/kmol, r. Cp = 500 kcal/m 3/K UA = 150 kcal/h/K, R = 1. 987 kcal/kmol/K V = 1 m 3, q =1 m 3/h, CAf = 10 kmol/m 3, Tf = 298 K, Tj = 298 K. Problem » Find the three steady-state points:

Matlab Tutorial cont. l FSOLVE – multivariable nonlinear zero finding » Matlab function for solving a system of nonlinear algebraic equations » Syntax: x = fsolve(‘fun’, xo) – Same syntax as fzero, but x is a vector of variables and the function, ‘fun’, returns a vector of equation values, f(x). » Part of the Matlab Optimization toolbox » Multiple algorithms available in options settings (e. g. trustregion dogleg, Gauss-Newton, Levenberg-Marquardt)

Matlab Exercise: Solution with fsolve l Syntax for fsolve » » l x = fsolve('cstr', xo, options) 'cstr' – name of the Matlab m-file function (cstr. m) for the CSTR model xo – initial guess for the steady state, xo = [CA T] '; options – Matlab structure of optimization parameter values created with the optimset function Solution for first steady state, Matlab command line input and output: >> xo = [10 300]'; >> x = fsolve('cstr', xo, optimset('Display', 'iter')) Norm of step First-order Trust-region Iteration Func-count f(x) optimality radius 0 3 1. 29531 e+007 1. 76 e+006 1 1 6 8. 99169 e+006 1 1. 52 e+006 1 2 9 1. 91379 e+006 2. 5 7. 71 e+005 2. 5 3 12 574729 6. 25 6. 2 e+005 6. 25 4 15 5605. 19 2. 90576 7. 34 e+004 6. 25 5 18 0. 602702 0. 317716 776 7. 26 6 21 7. 59906 e-009 0. 00336439 0. 0871 7. 26 7 24 2. 98612 e-022 3. 77868 e-007 1. 73 e-008 7. 26 Optimization terminated: first-order optimality is less than options. Tol. Fun. x = 8. 5637 311. 1702

Matlab Exercise: cstr. m function f = cstr(x) ko = 3. 493 e 7; E = 11843; H = -5960; rho. Cp = 500; UA = 150; R = 1. 987; V = 1; q = 1; Caf = 10; Tf = 298; Tj = 298; Ca = x(1); T = x(2); f(1) = q*(Caf - Ca) - V*ko*exp(-E/R/T)*Ca; f(2) = rho. Cp*q*(Tf - T) + -H*V*ko*exp(-E/R/T)*Ca + UA*(Tj-T); f=f';