Received 12 May 2016; accepted 9 July 2016; published 12 July 2016
The governing equation for a convective-dispersive transport phenomenon is given by:
where c is concentration; v is the flow velocity; is the dispersive tensor; t is the time variable and is the rate of solute mass exchange. Equation (1) yields satisfactory results in controlled environments in the absence of media heterogeneities. However, it exhibits certain limitations in modeling real life phenomenon. Its deviation from Fickian theory, results in an infinite speed of propagation of the scalar field (Cattaneo  Curtin and Pipkin1  , Ferreira and Oliveira  , Hristov  ). In order to overcome this unphysical property, Fick’s law is replaced by a flux accompanied by a memory term. It was proposed (Cattaneo’s  ) that the flux at a certain point x and at time t should be a consequence of scalar variation at some point x but after an elapsed time. Based on this, Ficks equation is modified to read:
where is a delay parameter. Applying a first order Taylor series expansion to the flux we obtain the so-called Cattaneo’s flux term.
Without any loss in generality, Equation (1) can be recast to read
where are the reaction and source/sink terms. The source/sink term accounts for the injection or extraction of mass into or out of the system. are the diffusion coefficient and velocity variables. The term inside the integral term accounts for the memory of the system. Equation (4) is a one-dimensional integro-partial differential equation which describes the transport of quantities such as heat, vorticity, energy and mass and is very important in many physical systems especially those involving the environment. It comes with the following constraints.
An initial condition of the type
and boundary conditions
where. In addition, are given constants.
Fundamental solutions generated from “free-space” Greens functions are usually applied in boundary integral analysis of transport equations. Efforts to make them more applicable to boundary element method (BEM) formulation, especially for non-homogeneous material, can be found in (Shaw  ) and (Shaw and Manolis  ). These have not only met with mixed results but have also been found to be inadequate for a good number of problems of engineering interest.
For example, adopting the classical boundary element approach to viscoelastic materials or to cases which incorporate the memory term requires fundamental solutions that are yet not well known except for the simple Maxwell material whose analytical solution is available (Gaul and Schanz  ). Numerical attempts to deal with this include that of Schanz and Antes  . They studied a BEM viscoelastic formulation based on the so called convolution integral method proposed by Lubich  in 1988. In this work, a quadrature formula whose weights are determined by the Laplace transform of the fundamental solution is applied to evaluate a convolution inte- gral and the entire formulation numerically handled by a multi-step procedure.
The hybrid boundary-integral-finite element procedure adopted herein employs the Green’s function of the Laplace operator, and retains the time marching feature of classical BEM, but its spatial discretization enlists a local support feature that is similar to that of the finite element. This guarantees some attractive outcomes; it facilitates the hybridization process especially for schemes that are element or nodal based (Onyejekwe   ) and establishes a hybrid computational method which exploits the strengths of its composite parts. This is the key motivation for this study.
2. Numerical Formulation
A typical boundary integral discretization requires both the complementary differential equation and the free- space Green’s function.
where is the source node and x is the field node and by the same token both t and represent the source and field nodes in the temporal coordinates. The flexibility offered by a hybridization procedure allows the integral kernel in equation (4) to be evaluated at each node.
We consider defined in the region where
We adopt a composite weighted trapezoidal rule in the temporal coordinate for the memory term.
Equation (4) is finally converted to its integral analog by the Green’s second identity to yield:
where are given and I is the “memory” or the hereditary’ term; can be a constant or a function of
the dependent variable,
This hybrid formulation unlike a typical boundary integral formulation requires that the computational domain be discretized into elements over which the scalar variable will assume some functional distributions. This is represented for the dependent variable as
where represents the current time level, while denotes previous time level. and are interpolating functions in space and time respectively. The element equation can now be written as:
Equation (14) is put in a compact matrix form
3. Numerical Test Examples
Example 1 Consider the following special case of Equation (4) with reaction term of the nonlinear variety (Khuri and Sayfy  ),
where subject to the initial condition
and boundary conditions
Next we account for the source term where
The exact solution is given by
The numerical solutions are obtained with respectively. The closeness of Figure 1 to that of Khuri and Sayfy  and the magnitude of absolute errors as a result of comparison between analytic and numerical results in Table 1 demonstrate the reliability of the numerical formulation developed herein.
Example 2 We consider a special case of Equation (4) given by:
Figure 1. 3D plot of scalar profiles for Ex. 1.
Table 1. Comparison of numerical and analytical results at time = 0.025.
At time = 0.025000, L2-Norm = 0.000178, L-Infinity Norm = 0.000463, Root Mean Square Error = 0.000899.
where, are the velocity and diffusion coefficient s respectively.
The initial and boundary conditions are given as:
Equation (21) is a transport equation that involves a convective-diffusive process coupled with a memory term. It is one of the most important equations of mathematical physics and has been used to describe the transport of such quantities like mass, heat, vorticity, energy and momentum. In the absence of the memory term, the analytical solution can be determined by method of separation of variables. This is expressed as:
The numerical solution has been found to develop a steep front or boundary layers at x = 1 (Feng and Tian 2006). Profiles of the numerical solution of equation (21) and the corresponding analytical solution without the memory term (equation (24)) are displayed in Figure 2(a) and Figure 2(b) for the following problem parameters:. While some instabilities can be detected in Figure 2(b), none can be found in 2a. This confirms an earlier work by Feng and Tian  . In addition, this observation validates the
Figure 2. 3D numerical scalar profiles (with memory term) for. 3D analytic scalar profiles (without the memory term) for.
interesting physics brought to play by the introduction of the memory term in the C-D equation. The modification of the flux term which takes care of the infinite speed of species propagation introduces an additional diffusion term. Again a close look at Equation (21) shows that we have two elliptic operators whose second order derivatives are multiplied by positive parameters. These derivatives model the diffusion process while the only first order derivative is associated with the convective process. Diffusion has always been known to play a significant role in the region of rapid changes (the boundary layers) where the scalar profile exhibits a steep profile. When a standard numerical procedure is applied to a C-D problem, if the diffusion term is relatively smaller than the convective term, the computed solution profile is often oscillatory. Whereas if there exists an “excess” (superfluous) diffusion term, the computed profile is smeared as demonstrated in Figure 2(a). It will also be interesting to compare the effects on the scalar profiles introduced by different magnitudes and discretizations of the convection term. However, the details of this will be the theme of a forthcoming paper.
The main thrust of this paper is a seamless incorporation of non-Fickian flux and memory into a localized element-based boundary element method. The advantage of this technique lies in the possibility of admitting realistic experimental information into a governing differential equation without complicating BEM formulation. For example, it offers a feasible alternative for applying BEM to describe transport processes in viscoelastic media where the inclusion of memory term is vital.