Hairer & Wanner’s Brusselator problem

A good test of the Rosenbrock methods can be applied to the PDE from Hairer & Wanner. The PDE’s for the diffusion equation (also known as the Brusselator in one spatial coordinate) are:

 

The initial & boundary conditions are:

The first step in the solution of this PDE is to convert it into equivalent ODE’s, with the use of the finite difference technique & the MOL. Here is the procedure to solve the PDE’s:

• Choose the number of grid points N to be used for the finite difference method.

• Calculate the discrete location of the points on the grid:

•  Calculate the spacing between the grid points 

•  Use the finite difference technique to represent the 2^nd order PDE’s as ODE’s

 

•  Use the MOL to convert the PDE’s to ODE’s

•  Solve the ODE’s in the above step. The Rosenbrock technique with the Kaps-Rentrop GRK4A coefficients was employed in this case, with N= 40.

Plotted below are the 3D solutions of vectors u & v. 

 

MOL vs Finite Difference

During a recent hallway conversation with my esteemed colleague Dr. Grigoriev, we were discussing the Method of Lines and he had posed the question of why it was referred to as the MOL and not the Finite Difference method.

Numerical solution of any differential equation (ordinary or partial) requires marching one step at a time, given the initial conditions at time t=0. The basic principle behind this marching goes back to the fundamental principle of the limit of a function. Given a function f(x), the derivative of the function can be represented as

 

The smaller the value of h, the more accurate is the derivative function. This principle of the limit, represented typically as the Taylor’s series expansion, is the “discretized” version of the ODE. Of course, “discretization” is a misnomer for an ODE since there is only one independent variable.

In the case of PDE’s, there exists the spatial dependence in addition to a time variable. The discretization thus has to be performed in both the spatial and the time domains. Once discretized separately as time and spatial domains, the independent variable becomes just one. The PDE has thus been converted into an ODE. This process of “replacement of the spatial derivatives with ODE’s” is referred to as the Method of Lines (ref: Hamdi et al).

So, how does one achieve this conversion from PDE to ODE? That is where the finite difference method comes in handy. The finite difference is thus the tool that discretizes the spatial domain, while the process of performing this change is the MOL.

Stated yet another way, the finite difference method is used to convert the 2^nd order partial equation by using the following equality,  

 

followed by the Method of Lines, which uses this knowledge to replace the spatial derivatives with the ODE’s, i.e.,

 

PDE solution with the Shampine coefficients.

Presented below is a solution to the heat equation using the Shampine coefficients for the Rosenbrock algorithm, with the following partial differential equation:

   

The initial & boundary conditions are as follows:

                               

The method of lines was used to convert the PDE’s to a familiar form of the ODE’s so as to allow easy implementation with the generalized Rosenbrock formulation. The slab of length 1 unit was divided into 20 segments. Plotted below is the temperature distribution over the slab, as a function of time. 

Rosenbrock – Shampine’s parameters

Shampine (Implementation of Rosenbrock Methods) derived a set of coefficients for the Rosenbrock methods, with a gamma value of 0.5. Presented below are the coefficients for Shampine’s method, with the scaling accounted for direct implementation in the generalized method (ref: here). The generalized Rosenbrock algorithm can be easily implemented with these coefficients.

         A21=2.
         A31=48./25.
         A32=6./25.
         C21=-8.
         C31=372./25.

         C32=12./5.
         C41=-112./125.
         C42=-54./125.
         C43=-2./5.
         B1=19./9.
         B2=1./2.
         B3=25./108.
         B4=125./108.
         E1=17./54.
         E2=7./36.
         E3=0.
         E4=125./108.
         GAMMA=.5
         C2= 0.1000000000000000E01
         C3= 0.6000000000000000
         D1= 0.5000000000000000
         D2=-0.1500000000000000E01
         D3= 0.2420000000000000E01
         D4= 0.1160000000000000