Ordinary differential equations numerical solutions

Batch and bio process modeling often results in a set of ordinary differential equations. Numerical solution techniques for solving these differential equations normally involve approximating the differential equations by difference equations that are solved in a step-by-step marching fashion. 

In this text all numerical problems involve integration of simultaneous ordinary differential equations or solution of simultaneous algebraic equations. You should have no trouble finding ways to solve algebraic equations with a calculator, a spreadsheet, a personal computer, etc. 

The three chapters in this volume deal with various aspects of singular perturbations and their numerical solution. The first chapter is concerned with the analysis of some singular perturbation problems that arise in chemical kinetics. In it the matching method is applied to find asymptotic solutions of some dynamical systems of ordinary differential equations whose solutions have multiscale time dependence. The second chapter contains a comprehensive overview of the theory and application of asymptotic approximations for many different kinds of problems in chemical physics, with boundary and interior layers governed by either ordinary or partial differential equations. In the final chapter the numerical difficulties arising in the solution of the problems described in the previous chapters are discussed. In addition, rigorous criteria are proposed for... 

Shampine S 1994 Numerical Solutions of Ordinary Differential Equations (New York Chapman and Hall)... 

L.R. Petzold, L.O. Jay, and J. Yen. Numerical solution of highly oscillatory ordinary differential equations. Acta Numerica, pages 437-484, 1997. 

Numerical Solution of Ordinary Differential Equations as Initial... 

Lapidus, L., and J. Seinfeld. Numerical Solution of Ordinary Differential Equations, Academic, New York (1971). 

NUMERICAL SOLUTION OF ORDINARY DIFFERENTIAL EQUATIONS AS INITIAL VALUE PROBLEMS... 

Solving the gas dynamics expressions of Kuhl et al. (1973) requires numerical integration of ordinary differential equations. Hence, the Kuhl et al. paper was soon followed by various papers in which Kuhl s numerical exact solution was approximated by analytical expressions. 

In the general case of a piston flow reactor, one must solve a fairly small set of simultaneous, ordinary differential equations. The minimum set (of one) arises for a single, isothermal reaction. In principle, one extra equation must be added for each additional reaction. In practice, numerical solutions are somewhat easier to implement if a separate equation is written for each reactive component. This ensures that the stoichiometry is correct and keeps the physics and chemistry of the problem rather more transparent than when the reaction coordinate method is used to obtain the smallest possible set of differential... 

A marching-ahead solution to a parabolic partial differential equation is conceptually straightforward and directly analogous to the marching-ahead method we have used for solving ordinary differential equations. The difficulties associated with the numerical solution are the familiar ones of accuracy and stability. 

Unsteady behavior in an isothermal perfect mixer is governed by a maximum of -I- 1 ordinary differential equations. Except for highly complicated reactions such as polymerizations (where N is theoretically infinite), solutions are usually straightforward. Numerical methods for unsteady CSTRs are similar to those used for steady-state PFRs, and analytical solutions are usually possible when the reaction is first order. 

In this chapter we concentrate on dynamic, distributed systems described by partial differential equations. Under certain conditions, some of these systems, particularly those described by linear PDEs, have analytical solutions. If such a solution does exist and the unknown parameters appear in the solution expression, the estimation problem can often be reduced to that for systems described by algebraic equations. However, most of the time, an analytical solution cannot be found and the PDEs have to be solved numerically. This case is of interest here. Our general approach is to convert the partial differential equations (PDEs) to a set of ordinary differential equations (ODEs) and then employ the techniques presented in Chapter 6 taking into consideration the high dimensionality of the problem. 

Another kind of situation arises when it is necessary to take into account the long-range effects. Here, as a rule, attempts to obtain analytical results have not met with success. Unlike the case of the ideal model the equations for statistical moments of distribution of polymers for size and composition as well as for the fractions of the fragments of macromolecules turn out normally to be unclosed. Consequently, to determine the above statistical characteristics, the necessity arises for a numerical solution to the material balance equations for the concentration of molecules with a fixed number of monomeric units and reactive centers. The difficulties in solving the infinite set of ordinary differential equations emerging here can be obviated by switching from discrete variables, characterizing macromolecule size and composition, to continuous ones. In this case the mathematical problem may be reduced to the solution of one or several partial differential equations. 

This equation must be solved for yn +l. The Newton-Raphson method can be used, and if convergence is not achieved within a few iterations, the time step can be reduced and the step repeated. In actuality, the higher-order backward-difference Gear methods are used in DASSL [Ascher, U. M., and L. R. Petzold, Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations, SIAM, Philadelphia (1998) and Brenan, K. E., S. L. Campbell, and L. R. Petzold, Numerical Solution of Initial-Value Problems in Differential-Algebraic Equations, North Holland Elsevier (1989)]. 

Other methods can be used in space, such as the finite element method, the orthogonal collocation method, or the method of orthogonal collocation on finite elements. One simply combines the methods for ordinary differential equations (see Ordinary Differential Equations—Boundary Value Problems ) with the methods for initial-value problems (see Numerical Solution of Ordinary Differential Equations as Initial Value Problems ). Fast Fourier transforms can also be used on regular grids (see Fast Fourier Transform ). 

Absorption columns can be modeled in a plate-to-plate fashion (even if it is a packed bed) or as a packed bed. The former model is a set of nonlinear algebraic equations, and the latter model is an ordinary differential equation. Since streams enter at both ends, the differential equation is a two-point boundary value problem, and numerical methods are used (see Numerical Solution of Ordinary Differential Equations as Initial-Value Problems ). 

Dynamic simulations are also possible, and these require solving differential equations, sometimes with algebraic constraints. If some parts of the process change extremely quickly when there is a disturbance, that part of the process may be modeled in the steady state for the disturbance at any instant. Such situations are called stiff, and the methods for them are discussed in Numerical Solution of Ordinary Differential Equations as Initial-Value Problems. It must be realized, though, that a dynamic calculation can also be time-consuming and sometimes the allowable units are lumped-parameter models that are simplifications of the equations used for the steady-state analysis. Thus, as always, the assumptions need to be examined critically before accepting the computer results. 

Analytical solutions of the self-preserving distribution do exist for some coalescence kernels, and such behavior is sometimes seen in practice (see Fig. 40). For most practical applications, numerical solutions to the population balance are necessary. Several numerical solution techniques have been proposed. It is usual to break the size range into discrete intervals and then solve the series of ordinary differential equations that result. A geometric discretization reduces the number of size intervals (and equations) that are required. Litster, Smit and Hounslow (1995) give a general discretized population balance for nucleation, growth and coalescence. Figure 41 illustrates the evolution of the size distribution for coalescence alone, based on the kernel of Ennis Adetayo (1994). 

Kunii and Levenspiel(1991, pp. 294-298) extend the bubbling-bed model to networks of first-order reactions and generate rather complex algebraic relations for the net reaction rates along various pathways. As an alternative, we focus on the development of the basic design equations, which can also be adapted for nonlinear kinetics, and numerical solution of the resulting system of algebraic and ordinary differential equations (with the E-Z Solve software). This is illustrated in Example 23-4 below. 

An implementation of this algorithm, using the sequential procedure within the MATLAB environment, was proposed by Figueroa and Romagnoli (1994). To solve step 2, the constr function from the MATLAB Optimization Toolbox has been used. The numerical integration necessary in this step has been performed via the function ode45 for the solution of ordinary differential equations. 

Figure 3.6 Numerical solution of ordinary differential equations sketch of the four steps of the Runge-Kutta method to the order four giving the n+1 th estimate y(n+1) from the nth estimate y ". ...

Digital simulation is a powerful tool for solving the equations describing chemical engineering systems. The principal difficulties are two (1) solution of simultaneous nonlinear algebraic equations (usually done by some iterative method), and (2) numerical integration of ordinary differential equations (using discrete finite-difference equations to approximate continuous differential equations). 

As discussed in the introduction to this chapter, the solution of ordinary differential equations (ODEs) on a digital computer involves numerical integration. We will present several of the simplest and most popular numerical-integration algorithms. In Sec, 4.4.1 we will discuss explicit methods and in Sec. 4.4.2 we will briefly describe implicit algorithms. The differences between the two types and their advantages and disadvantages will be discussed. 

The process inputs are defined as the heat input, the product flow rate and the fines flow rate. The steady state operating point is Pj =120 kW, Q =.215 1/s and Q =.8 1/s. The process outputs are defined as the thlrd moment m (t), the (mass based) mean crystal size L Q(tK relative volume of crystals vr (t) in the size range (r.-lO m. In determining the responses of the nonlinear model the method of lines is chosen to transform the partial differential equation in a set of (nonlinear) ordinary differential equations. The time responses are then obtained by using a standard numerical integration technique for sets of coupled ordinary differential equations. It was found that discretization of the population balance with 1001 grid points in the size range 0. to 5 10 m results in very accurate solutions of the crystallizer model. 

---