Introduction
Convective diffusion equation is one kind of basic motion equation that linearizes the nonlinear equation of viscous fluid. It can be used to describe river pollution, air pollution,distribution of pollutants in nuclear waste pollution, fluid flow and heat conduction in fluid, mass, heat transport process [1-26]. Considering the classification of equations, it belongs to parabolic or elliptic equations, but it also presents the basic properties of hyperbolic equations caused by convection-dominant. Therefore, it is of great theoretical and practical significance to construct a numerical method for three dimensional convection-diffusion equation which can reflect its characteristic properties.
[13,22,26-32], used traditional finite difference and [8,12] adopted finite element method, however numerical solutions produces serious numerical oscillations and numerical dissipation phenomena. For the purpose of eliminating the numerical oscillation caused by convection dominance [12,31] proposed the characteristic finite element method and characteristic finite difference method, which are effective in numerical calculation and applications. There are a lot of papers have discussed the numerical solution of convection-diffusion equation with high accuracy, good stability and for small diffusion coefficient in recent years as [16,17]. The above method improves the traditional method, but it also has many insurmountable effects. The streamline diffusion method reduces the numerical diffusion, but artificially imposes the streamline direction. The modified finite element method can be flexible with large time step without reducing the precision of approximation and has high stability at the front of flow front, which eliminates the numerical diffusion phenomenon, but too dissipated. The mixed finite element method for the convection dominant diffusion equation for the convection part of the equation is discretion by the characteristic difference method, the diffusion part of the equation is discretion by the mixed element method. For the larger time step, they all have a non-physical oscillation. [12,31] pointed out that the classical upwind format does not produce numerical vibration or oscillations, however discrete convection terms with only first order accuracy. [17] gives Modified Upwind Format is second-order accuracy, but for the strong convective dominance problem, the format is still only one order accuracy, and occurs serious numerical dissipation phenomenon. These methods above usually used to linear interpolation on convection terms. It does not oscillate, but the accuracy is low.
It is difficult to solve the numerical solutions for strong convection-dominant problems, [16,18,19,24-27,32] proposed corresponding improvement measures. The above method improves the traditional method, but it also has many insurmountable defects: the streamline dissipation method reduces the numerical diffusion, the disadvantage is: artificially impose the direction of the streamline. The modified finite element method can adopt a large time step without reducing the rescission of the approximation, has high stability, eliminates the numerical diffusion phenomenon, but also reflects the dissipation phenomenon. Because of the characteristics of the equation itself, it brings some difficulties to establish an accurate and effective numerical solution method.
We give an implicit difference scheme, which is proper to diffusion dominant cases. The difference scheme overcomes the disadvantages of the characteristic line direction method, which is a low-order difference in the time layer. Difference scheme is third ordered accuracy in the time layer and avoids non-physical oscillations. It overcomes the disadvantages of least square fitting method, in which diffusion phenomenon occurs because of higher polynomial fitting on diffusion terms. We theoretically and numerically verify that the implicit finite difference scheme is unconditionally stable and second order convergent in both space and time layers.The algorithms we propose in this paper can improve the efficiency without reducing the approximation accuracy and avoid the numerical diffusion phenomenon.
For the convection dominant diffusion problem and diffusion-dominant problem, the low-order scheme has serious numerical dissipation, and the high-order scheme is prone to numerical dispersion and non-physical oscillation. Therefore, we construct a numerical method with high precision, stable and suitable for small diffusion coefficients. It reflects the characteristic properties of hyperbolic equation.
1 Problem
We consider a three dimensional convection-diffusion equation
in which a, b, c, v1, v2, v3 constants and vi > 0 for i = 1, 2, 3. If vi < 0 i = 1, 2, 3, the initial-value problem is called not well-posed. If the a,b,c is smaller than v1, v2, v3, that is, the convection effect is relatively weak. In such problems, diffusion dominates, and the equations are elliptic or parabolic. If Pe is large, that is, the diffusion of solute molecules is slow relative to the fluid velocity. In such problems, the convection is dominant and the equations have the characteristics of hyperbolic equations.
1.1 If
then is called convection-dominant, if
the problem (1.1) is strong convection dominance.
We discussing the convection-dominant problem in section 2 and section 3, and the considering diffusion-dominant problem in section 4. The initial-boundary value of (1.1) is
which describes the diffusion of a substance in a medium that is moving with speed vi, i = 1, 2, 3. The solution of (1.1) is the concentration of the diffusing substance. Analytic solution of the initial-value problem (1.1) and (1.2) is given in [4]. Generally it is difficult to write down a formula for a classical solution. References [9] has considered four techniques of solving partial differential equations, It is too complicated and asks too much initial-boundary conditions [33-38]. In classical theory of heat and the energy equation of incompressible fluid flow satisfies three dimensional convection-diffusion equation. Solving the heat equation, energy equation is a hot topic in physics area. Since we give an implicit difference scheme for numerical solution of (1.1). We will give its compatibility, stability and convergence consequently.
1.1 Implicit difference scheme
Here we establish the implicit difference scheme
with initial condition
The difference scheme (1.3) is a two-layers implicit difference scheme which involves the fourteen different points on the time layers n and n+1
See Figure 1.
Figure 1: The points that the Implicit difference scheme (1.3) used.
Rewrite the initial-value problem as
and rewrite the difference scheme (1.1) as
Now we trying to calculate it's Truncated error.
1.2 Truncated error
Assume now u(x, y, z, t) is sufficient smooth and for t third order differential, for space variable x,y and z are fourth order differential. For calculation is easy, we choose the same step length h of x,y and z axis as
, the step length of the time space is τ. u(x, y, z, t) is the solution of initial-valued problem (1.4) and
is a solution of the difference-initial problem (1.5). Here we give Labels for later calculation
Gives Taylor Expanding of each blocks of (1.5) at point
we have
(1.7) - (1.8), we have
and likewise (1.13), (1.8)
Let
.
then we have
The truncated error of implicit difference scheme (1.5) is
delete the high level item, it is easy to see that the truncated error of implicit difference scheme (1.5) is
.
1.3 Compatibility
1.2 The differential scheme (1.5) is compatible with the differential equation (1.4).
Proof: Let us now check compatibility of (1.5)
According to the definition of compatibility, when the time step and spatial step to zero, the truncated error tends to zero, thus the difference equation (1.5) tends to initial-valued problem (1.6), then the compatibility conditions hold.
1.3.1 Precision: Rewrite the difference scheme (1.5), we have
Thus, difference scheme is second ordered. Truncated error of (1.5)
The difference scheme (1.5) is second-order precision for h, for τ respectively.
1.4 Stability
Generally use the Fourier Transformation Method for stability of the constant coefficient problem.
1.3 The difference scheme (1.5) is stable.
Proof: For finding growth factor of the difference scheme, we usually use Fourier transforming Method. Rearrange the items (1.26) we have
in which
The last item is third ordered corresponding τ, eliminating
here the difference scheme
Let
in which
Let F is stands for Fourier transformation,then each items in (1.31) after Fourier transformation is below
Likely we can write
Likely we can write
and
likely we can write
likely we can write
Then
and pay attention to minus items of numerator, then we have
in which
the implicit difference scheme (1.5) is unconditionally stable [39-45].
The stability of three dimensional problem is difficult than one and two dimensional problem.Even though (1.5) is unconditional stable theoretically, in practice restricted by mesh grid ratio and convection-diffusion coefficients.
It is worth to show that the explicit difference scheme for three dimensional problem double times constrictions about time step, that is why we propose the second ordered implicit difference scheme in this section. To solve the three dimensional convection-diffusion equation by using implicit difference scheme directly is difficult. So the numerical methods becomes hot topic in recent decades [46-50].
1.5 Convergance
Lax and Richmyer (1956) gives the Lax Equivalence Theorem which helps to determine the convergence when we dont know the exact solution. It is always used to constant coefficient problem. For variable coefficient problem, we have used to Energy Inequality Method.
1.4 According to the Lax Equivalence Theorem, The implicit difference scheme (1.5) is converges to three dimensional convection-diffusion equation.
The Lax equivalence theorem can be used when the initial value problem is linear, well-posed and has the periodic initial and boundary condition. The problem (1.1) is a well-posed the First-Dirichlet Initial - Boundary problem with periodic initial condition. It is worth to show that the explicit difference scheme for three dimensional problem wants triple times constrictions about time τ. So we propose the second ordered implicit difference scheme in this paper [51-60]. Here we give an algorithm of numerical solution of the three dimensional convection-diffusion equation which using (1.5).
1.6 Stepwise alternating direction implicit difference method
The three dimensional problem
We rewrite implicit method (1.5), we need initial and boundary condition
We want special step
for calculation. The The first,second and third terms are the diffusion terms at the right end and the second, third terms at the left end are the convection term. Because of the characteristics of the equation itself, it is difficult to establish an accurate, effective method. Convection velocity corresponding x,y and z are constants, v1,v2,v3 are diffusion constant coefficients. If the convection coefficients a,b and c are small, the convection effect is relatively weak, and diffusion dominates, equations are elliptic or parabolic [61-65].
If the number of a,b and c are large, the diffusion of solute molecules is slow relative to the fluid velocity. In such problems, the convection is dominant, the equation has the characteristics of hyperbolic equations. The problem (1.1) and initial-boundary conditions (1.2), the diffusion of a substance in a medium that is moving with speeds to x,y and z directions. The unknown function is the concentration of the diffusing substance. The conventional Galerkin Finite element method is used to solve the convection dominant problem.
Peaceman, Rachford and Douglas proposed Alternating Direction Implicit Method.
we have use the implicit method (1.69) for Lu, as
we do implicit difference to x direction to L1u of
and the
is a high ordered residual deleting it,we have
This is a Fractional Step method, second ordered unconditionally stable. The numerical solution is need initial -boundary conditions. Now we give an Alternating Direction Implicit method (ADIM) algorithm for three dimensional problem (??) which using the implicit difference scheme (??) inderectly. It is necessary to show that, now we propose the details which change the implicit to explicit, and improve the accuracy as well, reducing the computational work [66-70].
Algorithm of ADIM
step one: step 1: input initial condition
and boundary condition
and function u(x, y, t), and constant coefficient a, b, v, step length l λ,η.
step 2: Pretreatment: calculating
beyond x direction and
beyond y direction, and
are for right boundary condition.
step 3: calculate the
by using
calculate the
using
step 4: iterating step 2 and step 3, stop when n=N-1.
Because the difference scheme proposed step 2 is unconditionally stable that will not affect the convergence of (??) and (??). This is the lively usage of Fractional Step Methods(FSM) and Alternative
1.7 Greedy algorithm
1.5 The problem (1.1) is uniquely solved by (1.5) directly.
Proof: According to the time layers, splitting (1.1) on three direction, we have Below we use the each blocks of implicit method (1.5) on each equation of (1.37).
1.7.1 Implicit split on x axis
Rewrite (1.5) on x direction, we have the iterative linear system below
in which
Let,
then (1.51) represents a linear system as
The coefficient matrices and corresponding vectors in linear system (1.52)
1.6 The sufficient conditions to find a unique solution of the linear system (1.52) using the catch-up method is these inequalities hold on
are the tri-diagonal elements of Ax, or the coefficients of
in (1.51) like
If
is a diagonally dominant, inverse of A exist.
1.7 When
the coefficient matrix Ax in is diagonally dominant.
Proof: Because if (1.59) hold on, then we have
the coefficient matrix is diagonally dominant. The sufficient condition that (1.52) has a unique solution is the inequalities (1.59) hold on.
1.8 If (1.59) holds on, then (1.58), holds on too, then the linear system (1.52) has a unique solution.
1.9 If
and
then A- exist, and
Proof: The linear system (1.52) is an iteration form on time layers
If
then the difference iteration method (1.52) is convergent.Because
in which
is a diagonal matrix, and
is a positive real number, suppose kx is a real number, and
,let
, then (1.61) becomes
1.10 When
then
The solution of linear system (1.52) is convergent on the x direction.
1.11 The difference scheme (1.51) has unique solution on x direction, requires these inequities established
1) Ax is a diagonally dominant matrix. 2)
1.7.2 Implicit Split On y Axis
We can rewrite (1.5) beyond y directions
Let l=1,2,...,L-1 ,then (1.62) represents a linear system,as
1.12 If the inequality system
holds on, then the coefficient matrix in (1.63) is diagonally dominant.
1.13 When inequality system (1.64) is hold on,then the linear system (1.63) has a unique solution on y direction.
Proof is same as Theorem 2.10, here we omit.
In linear system (1.63), the coefficient matrices and corresponding vectors are
1.14 When
is a diagonally dominant,
exist, and (1.63) has a unique solution on y direction.
1.7.3 Implicit split on z Axis
The same way, we rewrite (1.5) on z directions, we have
Let m=1,2,...,M-1, then (1.70) represents a linear system as
In which, the coefficient matrices and corresponding vectors are
1.15 If the inequality system
holds on, then the coefficient matrix in (1.72) is diagonally dominant.
1.16 When inequality system (1.77) hold on, the linear system (1.71) has a unique solution on z direction.
Implicit difference scheme (1.5) for three dimensional convection-diffusion equation (1.1) constructed by three one dimensional linear systems [71-80]. Here the implicit difference method is second ordered on time layers and on x,y,z directions, second ordered on convection terms,second ordered on diffusion terms.The convergence better than Finite element method, characteristic line method,and mesh-less method which is only first ordered on time layers and first ordered on convection terms, see [18,19,24].
2. Implicit Split Alternating Direction Method (ISADM)
The numerical solution of (1.37) on x,y,z direction separately have no physical meaning. Solving the three dimensional convection-diffusion equation, we give an Implicit Split Alternating Direction Difference Method. This algorithm based on (1.5), which solves the problem (1.1) globally, six ordered the convergence on time layers, at-least second ordered on convection terms.Equation (1.1) is a tri-linear system, coefficients are three dimensional tensor. The main idea of this method is splits (1.1) on the x,y and z directions,each blocks of it use implicit format blocks of (1.5), then iterating (1.52), (1.63), (1.71) alternately. Split equation (1.1)
into the sum of three one-dimensional equations such as:
Here we use implicit difference schemes in each partial differential equation in (2.1) on x,y and z directions separately. At the time layers, we use
step once, then we have
It is clear that adding (2.2), (2.2) (2.3) and (2.4) is equal to the implicit difference scheme (1.5). Time t is a curvature on fourth order space. Three dimensional problem (1.1) solved by iterating (2.2), (2.3), (2.4), in each step, work only on one direction, then iterate alternately. By this way, we can use implicit difference scheme (1.5) indirectly, solve the three dimensional convection diffusion equation (1.1) globally, and greatly improves accuracy, numerical solution quickly converges to analytic solution, and with no unnecessary vibration [81-88].
This Implicit Split Alternating Direction Method(ISADM) can also usedto the convection dominant,diffusion-dominant problem.
2.1 Convergences of ISADM
Rewrite (2.2), (2.3) and (2.4) we have
When
, three difference equation (2.5), (2.5), (2.6), (2.5) (2.5) (2.4) becomes three-dimensional linear systems
The coefficient matrices and corresponding vectors are same as (1.52), (1.52) (1.63), (1.52) (1.63) (1.71) Calculating the values on level n+1 by using (2.8).
2.1 If inequalities system
hold on in (2.8), then
exist, and we have
in which
same as
in Theorem 2.10.
2.2 (1.1) can be uniquely solved by linear systems (2.8) in each x,y and z direction alternately.
2.3 Solution of three dimensional linear system (2.8) is convergent to the analytic solution of (1.1)
Proof: We check the convergence of (2.8). Input the first linear system into the second one in (2.8), we have
We put (2.11). into the third linear system in (2.8), we have
Let
then (2.12) becomes
because we have
then we know (2.8) at least six ordered convergent.Problem (2.12) are inverse problem, three-dimensional tri-linear system (2.8) has a unique solution when the inequality systems (2.10) hold on. So (2.8) convergent. The speed of convergence rate far exceeds the Super-Relaxation Method (SOR), quicker than general Alternating Direction Implicit Method. Here we give the Algorithm of Alternating Direction Implicit Split Difference Method below [89,90].
2.2 Pretreatment for initial-boundary condition
Numerical solution of (1.1) by using (2.8) always depends on initial boundary conditions. If the initial-boundary condition is continuous, use implicit split scheme (2.8), numerical solution converges to analytic solution so quickly Figure 2.
Figure 2: The pretreatment of left and right initial-boundary condition on x axis.
If the initial-boundary conditions is discontinuous, then there appears oscillation near the discontinuous initial-boundary points. Hence we show the pretreatment for discontinuous initial-boundary points. Before we calculate the numerical solution, do this pretreatment first. See Figure 3, for pretreatment on x axis. The pretreatment on y and z axis is the same. We omit their Figure here. When j = 1, pretreatment for left boundary condition on x direction is
Figure 3: The pretreatment on y and z axis is the same.
see the left side of Figure 3. When l = 1, pretreatment for left boundary condition on y direction is
When m = 1, pretreatment for left boundary condition on z direction is
While j = J - 1, the pretreatment for right boundary condition on x axis direction is
see the right of Figure 3. While l = L-1, the pretreatment for right boundary condition on y axis direction is
While m = M-1, the pretreatment for right boundary condition on z axis direction is
Inside of definition domain, as
, we use (2.5), (2.5) (2.6), (2.5) (2.6) (2.7) alternately. We calculate all value of
from
Now we will give the algorithm of ISADM, see Algorithm 1. Now we give the algorithm of Stepwise Alternating Direction Implicit Method(SADIM)for three dimensional problem (1.1), which the implicit difference scheme
.
When
,
becomes
2.2.1 Numerical example
Three dimensional convection diffusion equation
and initial- boundary condition given as
The initial-boundary problem is a steady state blocks on four dimensional space, w hen t=0, then it is a super symmetric 3d objects which formed by the intersection of six paraboloids. Here we show its one block on axis y. See Figure 4. When t not equal to zero, then t is curves on four dimensional space. The super symmetric 3D objects moves along the t curve, and her whole body and endogenous changes when time t changes [91-96].
Figure 4: When t not equal to zero,then t is curves on four dimensional space.
2.2.2 Other effects of stability
The main idea of this section is, the difference scheme is established to ensure that its solution satisfies the properties of mass, momentum and total energy conservation on the whole solution region and even on each grid. In a word, the study on the conservation of difference scheme should include two contents: How to construct the appropriate conservation scheme; that the results of numerical solution should be tested when the conservation scheme is difficult to design [97-100]. It is obvious the conservation of this practical calculation result is closely related to the grid scale effect, differential remainder effect and numerical boundary effect and conservation of format mentioned above. Therefore, conservation should be used as a test of the calculation process,according to which the grid scale should be checked and dissipating dispersion and numerical boundaries to ensure the stability of the calculated results.
The difference scheme (1.5) proper to strong diffusion dominant problem, and the problem of time dominance. For the diffusion dominant case, when the v is too big or too small, or the time step chosen too large or too small, the numerical solution dispersive [101-132]. Though the difference scheme is unconditionally stable in theoretically, in practice, the stability of the numerical solution affected grid scale effects, residual effect of difference scheme, numerical boundary effect, and conservation effect of scheme, etc.
Difference Equation (1.5) can easily be extended to more than three spatial dimension, and the analytical solutions for several dimensions are also available. Thus, they can be used to test numerical methods in one, two, three and more dimensional convection-diffusion equations and still unconditionally stable and higher accuracy on time step and convection terms.This is the lightning of our paper.
Conclusion
This paper mainly discussing a numerical solution for tree dimensional convective-diffusion equation with initial-boundary condition. An implicit difference scheme constructed, its truncated error, stability, convergence are analyzed. For solving convergent numerical solution, we give an Implicit Split Alternating Direction Method and it's Algorithms. If
, the problem (1.1) is well-posed, and can solve (1.1) uniquely by the implicit difference scheme (1.5). The implicit method is unconditional stable, second-ordered convergent on x,y, and z direction. Theoretical analysis and experimental results show that when the grid ratio is properly selected, the difference scheme is quite stable.
However, the convection coefficient is too large, there appears non-physical vibration. Use Algorithm we proposed reduces the non-physical oscillation and eliminates the dispersion effectively.
The problem with discontinuous initial-boundary condition, much more numerical solution appears perturbation near the discontinuous points. We have use the Saus scheme before use algorithm, retreating discontinuous initial-boundary conditions. Since the implicit difference scheme performs better than the standard Galerkin finite difference scheme, and quicker than SOR iteration method.
The implicit difference scheme (1.5) used in Algorithm alternately, eliminating discontinuous initial-boundary conditions one by hand, solving three dimensional problem, the algorithm almost six ordered accuracy. It can be extended to three dimensional variable coefficient convective-diffusion equation directly. The characteristic line method and Finite Element methods uses all of the initial-boundary conditions and the three degrees differential values, so the computational work is much more. The implicit difference method we proposed are over coming this deficiency. For three dimensional problem, Implicit method works three directions at a time, thus tripling the calculation at a point. The parallel computational work easily processing this algorithm in this paper, six ordered convergent, unconditionally Von-Neumann stable.For convection dominance and dig divergence dominance, the number of meshes should be increased and the length reduced.
However, no more in-depth study was undertaken in this area. We comparing with the characteristic difference method and the characteristic finite element method, the method we proposed is effective, simple in contraction and have a good convergence order. Therefore, exploring the numerical solutions of these problems using implicit difference format provides a feasible numerical model and a convenient, quick operation method for both theorists and practitioners, which has certain theoretical value and application value. Authors are grateful for the stimulating discussions and enlightenment's from from my teacher.