Of PDE solvers

Schnepf, A., Schrefl, T., Wenzel, W.W., 2002. The suitability of pde-solvers in rhizosphere modeling, exemplified by three mechanistic rhizosphere models. J. Plant Nutr. Soil Sci. 165, 713-718. 

Within esqjlicit schemes the computational effort to obtain the solution at the new time step is very small the main effort lies in a multiplication of the old solution vector with the coeflicient matrix. In contrast, implicit schemes require the solution of an algebraic system of equations to obtain the new solution vector. However, the major disadvantage of explicit schemes is their instability [84]. The term stability is defined via the behavior of the numerical solution for t —> . A numerical method is regarded as stable if the approximate solution remains bounded for t —> oo, given that the exact solution is also bounded. Explicit time-step schemes tend to become unstable when the time step size exceeds a certain value (an example of a stability limit for PDE solvers is the von-Neumann criterion [85]). In contrast, implicit methods are usually stable. 

Access to powerful computers and to commercial partial-differential-equation (PDE) solvers has facilitated modeling of the impedance response of electrodes exhibiting distributions of reactivity. Use of these tools, coupled with development of localized impedance measurements, has introduced a renewed emphasis on the study of heterogenous surfaces. This coupling provides a nice example for the integration of experiment, modeling, and error analysis described in Chapter 23. 

To overcome the inconsistence between the LBM and other PDE solvers on a same grid set, an alternative solution is to modify the boundary condition treatments of the PDE solver for the electric potential distribution to be consistent with the LBM bounded by the bounce-back rale. In this contribution, the halfway... 

The solvers available in COMSOL are extremely customizable and capable of handling a wide array of PDEs, geometries, and initial/boundary conditions. Here, we present an example steady-state solver setup, with multiple optional features enabled. 

The perturbed state variables are complex and the coefficients use the steady state solution as functions of x. Instead of a time solution, the PDE solver now has a series of steady state solutions with co as a varying parameter. 

The boundary condition implementations play a very critical role in the accuracy of the numerical simulations. The hydrodynamic boundary conditions for the LBM have been studied extensively. The conventional bounce-back rule is the most popular method used to treat the velocity boundary condition at the solid-fluid interface due to its easy implementation, where momentum from an incoming fluid particle is bounced back in the opposite direction as it hits the wall [20]. However, the conventional bounce-back rule has two main disadvantages. First, it requires the dimensionless relaxation time to be strictly within the range (0.5,2), otherwise the prediction will deviate from the correct result. Second, the nonslip boundary implemented by the conventional bounce-back rule is not located exactly on the boundary nodes, as mentioned before, which will lead to inconsistence when coupling with other partial differential equation (PDE) solvers on a same grid set [17]. 

Figure 10.30 shows a 3D plot for u (i.e., 0) as a function of x and f, using 200 mesh points for each dimension. Notice that the general behavior or topology of the surface is very similar to that of the anal5d ical solution shown earlier in Fig. 10.24. Keep in mind that the solution given by the MATLAB PDE solver is approximate.

The initial conditions and the spatial array of x points are defined on lines 19 and 20 of the listing. The call to the PDE solver on line 22 then returns a set of solutions from 0 to tmax (0.01 in this example) with 5 time solutions returned and with 5 sub time intervals used in solving the equations. Finally line 23 plots three of the solutions and line 24 saves all the ealeulated results. In order to monitor the progress of the solution set, the nprint = 1 parameters are set on tine 5 for both the odefd() and pdeivbv() functions. This produces printed output such as that shown in the selected output listing. It can be noted in the selected output that the solutions at each time increment are requiring 2 iterations of the basic Newton method used to solve the sets of equations. It is noted that the corrections on the second iteration are always within the machine aceuracy limits. This is because the equa-... 

A few other bookkeeping sections in the code convert the input solution array from an x-y labeled two index matrix into a single column array if needed on lines 454 through 459 and a final loop converts the single column solution array into a j,i two index matrix of solution values on lines 486 through 488. The type of solution desired may be specified as a table of values in order to pass additional parameters to the PDE solver. If a table is specified for the typsola parameter, two additional parameters, rx and ry, can be extracted from the input list as seen on line 444. For a direct sparse matrix solution with the SPM parameter these are not used. For the COE and SOR methods, this feature can be used to specify a spectral radius value with only the rx parameter used. For the ADI method these two values can be used to specify wx and wy parameters. In all cases, default values of these parameters are available, if values are not input to the various functions. 

The balance equations described in the previous sections include both space and time derivatives. Apart from a few simple cases, the resulting set of coupled partial differential equations (PDE) cannot be solved analytically. The solution (the concentration profiles) must be obtained numerically, either using self-developed programs or commercially available dynamic process simulation tools. The latter can be distinguished in general equation solvers, where the model has to be implemented by the user, or special software dedicated to chromatography. Some providers are given in Tab. 6.3. 

Using the boundary conditions (equations (5.54) and (5.55)) the boundary values uo and un+i can be eliminated. Hence, the method of lines technique reduces the nonlinear parabolic PDE (equation (5.48)) to a nonlinear system of N coupled first order ODEs (equation (5.52)). This nonlinear system of ODEs is integrated numerically in time using Maple s numerical ODE solver (Runge-Kutta, Gear, and Rosenbrock for stiff ODEs see chapter 2.2.5). The procedure for using Maple to solve nonlinear parabolic partial differential equations with linear boundary conditions can be summarized as follows ... 

In section 5.2.4, a stiff nonlinear PDE was solved using numerical method of lines. This stiff problem was handled by calling Maple s stiff solver. The temperature explodes after a certain time. The numerical method of lines (NMOL) technique was then extended to coupled nonlinear parabolic PDEs in section 5.2.5. By comparing with the analytical solution, we observed that NMOL predicts the behavior accurately. 

FEMLAB is a partial differential equation solver (PDE) available commercially from COMSOL. Inc, Included with this text is a special version of FEMLAB that has been prepared to solve problems intolving tubular reactors. Specifically, one can solve CRE problems with heat elTccts involving both axial and radial gradients in concentration and temperature simply by loading the FEMLAB CD on one s computer and running the program. One can also use it to solve isothermal CRE problem,s with reaction and diffusion. 

The mathematical model forms a system of coupled hyperbolic partial differential equations (PDEs) and ordinary differential equations (ODEs). The model could be converted to a system of ordinary differential equations by discretizing the spatial derivatives (dx/dz) with backward difference formulae. Third order differential formulae could be used in the spatial discretization. The system of ODEs is solved with the backward difference method suitable for stiff differential equations. The ODE-solver is then connected to the parameter estimation software used in the estimation of the kinetic parameters. More details are given in Chapter 10. The comparison between experimental data and model simulations for N20/Ar step responses over RI1/AI2O3 (Figure 8.8) demonstrates how adequate the mechanistic model is. 

A few examples will be demonstrated in the following section. Each feature was incorporated in the software because it has been found necessary or useful in some modelling project. All the modelling tasks that previously required tailor made solutions in each project can now be solved in a unified manner in the ModEst environment. ModEst is able to deal with explicit algebraic, implicit algebraic (systems of nonlinear equations) and ordinary differential equations. As the standard way to handle PDE systems, the Numerical Method of Lines, which transforms a PDE system to a number of ODE components is used. In addition, any model with a solver provided may be dealt with as an algebraic system. The basic numerical tools are contained in the well tested public domain software (Bias, Linpack, Eispack, LSODE). 

---