Introduction to pde1dMg¶
The section describes an algorithm and a corresponding integrator function pde1dMg for solving a system of partial differential equations
Equation 1
\(u_t \equiv \frac{\partial u}{\partial t} = f (u, x, t), x_L < x < x_R, t > t_0\)
This software is a one-dimensional differential equation solver. It requires the user to provide initial and boundary conditions in addition to a function for the evaluation of \(u_t\). The integration method is noteworthy due to the maintenance of grid lines in the space variable, \(x\). Details for choosing new grid lines are given in Blom and Zegeling, (1994). The class of problems solved with pde1dMg is expressed by Equation 1 and given in more detail by:
Equation 2
The vector \(u \equiv \left[ u^1, \ldots, u^{\mathit{NPDE}} \right]^T\) is the solution. The integer value \(\mathit{NPDE} \ge 1\) is the number of differential equations. The functions \(R_j\) and \(Q_j\) can be regarded, in special cases, as flux and source terms. The functions \(u, C_{j,k}, R_j, Q_j\) are expected to be continuous. Allowed values for the integer \(m\) are any of \(m = 0, 1, 2\). These are respectively for problems in Cartesian, cylindrical or polar, and spherical coordinates. In the two cases with \(m > 0\), the interval \(\left[ x_L, x_R \right]\) must not contain \(x = 0\) as an interior point.
The boundary conditions have the master equation form
Equation 3
In the boundary conditions the functions \(\beta_j\) and \(\gamma_j\) are continuous. In the two cases with \(m > 0\), with an endpoint of \(\left[ x_L, x_R \right]\) at 0, the finite value of the solution at \(x = 0\) must be ensured. This requires the specification of the solution at \(x = 0\), or it implies that \(R_j | _{x = x_L} = 0\) or \(R_j | _{x = x_R} = 0\). The initial values satisfy \(u \left( x, t_0 \right) = u_0 (x), \phantom{...} x \in \left[ x_L, x_R \right]\) , where \(u_0\) is a piece-wise continuous vector function of \(x\) with \(\mathit{NPDE}\) components.
The user must pose the problem so that mathematical definitions are known for the functions
These functions are provided to the function
pde1dMg in the form of two user-supplied
functions. This form of the usage interface is explained below and
illustrated with several examples. \(u_0\) can be supplied as the input
argument u
or by an optional user-supplied function. Users comfortable
with the description of this algorithm may skip directly to the
Examples section.
Description Summary¶
Equation 1 is approximated at \(N = \mathit{ngrids}\) time-dependent grid values \(x_L = x_0 < x_1 < \ldots < x_i(t) < \ldots < x_{N+1} = x_R\). Using the total differential \(\tfrac{du}{dt} = u_t + u_x \tfrac{dx}{dt}\) transforms the differential equation to the form
Using central divided differences for the factor \(u_x\) leads to the system of ordinary differential equations in implicit form
The terms \(U_i, F_i\) respectively represent the approximate solution to the partial differential equation and the value of \(f (u, x, t)\) at the point \((u, x, t) = \left( u_i, x_i(t), t \right)\). The truncation error from this approximation is second-order in the space variable \(x\). The above ordinary differential equations are underdetermined, so additional equations are added for determining the time-dependent grid points. These additional equations contain parameters that can be adjusted by the user. Often it will be necessary to modify these parameters to solve a difficult problem. For this purpose the following quantities are needed:
The values \(n_i\) are the so-called point concentration of the grid. The parameter \(\kappa \geq 0\) denotes a spatial smoothing value. Now the grid points are defined implicitly so that
The parameter \(\tau \geq 0\) denotes a time-smoothing value. If the value \(\tau\) is chosen to be large, this results in a fixed spatial grid. Increasing \(\tau\) from its default value avoids the error condition where grid lines cross. The divisors are defined by
The value \(\kappa\) determines the level of clustering or spatial smoothing of the grid points. Decreasing \(\kappa\) from its default values also decreases the amount of spatial smoothing. The parameters \(M_i\) approximate arc length and help determine the shape of the grid or \(x_i\) distribution. The parameter \(\tau\) prevents the grid movement from adjusting immediately to new values of the \(M_i\) , thereby avoiding oscillations in the grid that cause large relative errors in the solution. This is important when applied to solutions with steep gradients.
The discrete form of the differential equation and the smoothing equations are combined to yield the implicit system of differential equations
This is usually a stiff differential-algebraic system. It is solved using
the integrator differentialAlgebraicEqs. If
differentialAlgebraicEqs
is needed during the evaluations of the
differential equations or boundary conditions, it must be done in a separate
thread to avoid possible problems with pde1dMg’s
internal use of differentialAlgebraicEqs
. The only options for
differentialAlgebraicEqs
set by pde1dMg
are the Maximum BDF Order,
and the absolute and relative error values, documented as maxBdfOrder
,
and atolRtolScalars
.