Journal of Water Resource and Protection
Vol. 1  No. 5 (2009) , Article ID: 867 , 6 pages DOI:10.4236/jwarp.2009.15043

Two Modified QUICK Schemes for Advection-Diffusion Equation of Pollutants on Unstructured Grids

Linghang XING

Hydraulic Research Laboratory, Changjiang River Scientific Research Institute, Wuhan,China


Received June 28, 2009; revised August 17, 2009; accepted August 26, 2009

Keywords: Unstructured Grids, Q-QUICK/UQ-QUICK, Numerical Computation, Advection-Diffusion Equation of Pollutants


In this paper, two modified QUICK schemes, namely Q-QUICK and UQ-QUICK, for improving the precision of convective flux approximation are verified in advection-diffusion equation of pollutants on unstructured grids. The constructed auxiliary nodes for Q-QUICK/UQ-QUICK are composed of two neighboring nodes plus the next upwind node, the later node is generated from intersection of the line of current neighboring nodes and their corresponding interfaces. 2D unsteady advection-diffusion equation of pollutants is conducted for their verifications on unstructured grids. The numerical results show that Q-QUICK and UQ-QUICK have similar computational accuracy to the central difference scheme and similar numerical stability to upwind difference scheme after applying the deferred correction method. In addition, their corresponding CPU times are approximately equivalent to those of traditional difference schemes and their abilities for adapting high grid deformation are robust.

1. Introduction

The physical processes of pollutant transportation in flowing water are mainly consisting of advection and diffusion, which are usually governed by advection-diffusion equation of pollutants. Generally, the action of advection process would dominate the transportation of pollutants. Thereby, it is so important to improve the precision of convective flux approximation. Till to now, many works have been done for numerical simulation of advection-diffusion equation of pollutants on structured grids [1−8]. The natural computational boundaries are irregular. If structured grids applied, the workload of CFD pretreatment would be increased and the numerical precision in the boundaries would be decreased. Unstructured grid can produce arbitrary geometry and can well fit to complex physical boundary. Presently, the computation based on unstructured grid becomes more and more popular. However, many high-precision schemes on uniform grid can not be applied to unstructured grids directly. It is significative to extend them to unstructured grid computation. In the past decade, a number of difference schemes to calculate convective flux were developed for incompressible flow simulation. They include upwind difference scheme

(UDS), central differencing scheme (CDS), hybrid differencing scheme (HDS), the quadratic upstream, quadratic upstream extended and quadratic upstream extended revised difference schemes (QUICK, QUDSE, QUDSER as modified by Pollard and Siu [9]), the locally exact scheme (LEDS) (1972), and the power difference schemes (PDS) of Pantankar [10]. The unconditionally convergent schemes UDS/HDS/LEDS/PDS can be significantly inaccurate under coarse grids, thus they require considerable grid refinement to produce acceptable results. This makes them expensive. Moreover, they implicitly introduce the numerical diffusion term and distort the solution. In terms of accuracy and computational efficiency, it appears that the QUICK/QUDSE/QUDSER may offer the best compromise [11]. In uniform grids, they can have over second-order precision for convective flux approximation. However, QUICK/ QUDSE/ QUDSER need two nodes upstream. It is not so easy to apply these high-order schemes to unstructured grid directly, especially in three dimensional problems. Moreover, to find the exact locations of the next upwind node would require a very complex pointer system and consume more memory and CPU time. Davidson L. 1996 [12] introduced one method where the next upwind node is constructed by intersection from the line of two adjacent central nodes and its corresponding interface. Presentlyit is named Q-QUICK. According to the method, another similar scheme named UQ-QUICK is also introduced for comparative investigation. As for time discretization, the first-order and second-order schemes are taken into account. To this end, many compound schemes are formed and their corresponding numerical performances including numerical precision, stability and CPU time are fully demonstrated in 2D unsteady advection-diffusion equation of pollutants. In order to accelerate convergence speed for linear equations, generalized Minimum Residual (GMRES) [13] method with the Incomplete LU (ILUT) precondition is used.

2. Q-Quick and UQ-Quick

QUICK [14] is the third order approximation of the convection term. However, this high-order scheme is not easy to apply to unstructured grid directly. To find the exact location of the next upwind node would increase the geometrical complexity and consume relative more memory and CPU time. For a compromise, DAVIDSON L. introduced a modified QUCIK scheme in the end of his paper [12]. However, little work has been done later for this recommended scheme. Thus, it is necessary for fully studying its numerical characteristics including numerical precision, convergent stability and CPU time. As for comparison, another modified QUICK scheme is also introduced, namely UQ-QUICK.

2.1. QUICK

In a uniform grid (Figure 1 for definitions of points WW, W, P, E, EE etc.), the quadratic upstream difference scheme of Leonard [14] at the east cell-face can be written as


It is a third order approximation of the convection term.

2.2. Q-QUICK

Considering Figure 2(a), the next upwind node U is constructed by intersection of line and its corresponding interface. A better way is probably to use reconstruction schemes, namely to compute the gradient in node P and use Taylor expansion to obtain the value at point U.


For present study, the second order approximation is considered. Thus the first two terms from the right of the equation are kept. And then, it is imperative to estimate

Figure 1. QUICK scheme in uniform grid.

Figure 2. The next upwind node reconstruction for Q-QUICK.

Figure 3. The next upwind node reconstruction for UQ-QUICK.

the face value at interface.In doing this, it is assumed that the flow direction is from right to left. To this end,the face value can be interpolated by QUICK method as the same way acted on the structured grid. So the normal face value is derived according as



If the flow direction is inverse (see Figure 2(b)), the same Formula (3) can be derived.


The reconstruction method of NQ-QUICK has no intrinsic difference to that of Q-QUICK except that the three auxiliary nodes, and (see Figure 3) in a line are perpendicular to interface. In addition, the node values of and should be obtained first before interpolation by QUICK method. However, a tough problem here is that the relative positions of nodes, and are not as simple as those of Q-QUICK. In present study, two cases are summarized as follows.

l Case 1 If (see Figure 3(a))



If (see Figure 3(b)), similar to that of case 1 ().

l Case 2 If (see Figure 3(c))


where,;, andhave the same formulas as written in (4).

On high twisted grid, the coefficients of and in (5) may be much larger, which can cause the face value abnormal and result in the whole computation failure. So some flexible method should be implemented for the whole computational continuity. A better way for mending the defect may use the linear interpolation method for assuring a second-order precision at least. Then, the interface value can be written as


where, is the linear interpolation factor.

To this end, the coefficients of andare as follows

; (7)

If (see Figure 3(d)), similar to that of case 2().

3. Numerical Verification

In present study, a two-dimensional advection-diffusion equation of pollutants is discretized by compounds of different time schemes and different convective flux schemes. The former includes UDS and Crank-Nicolson and the later consists of UDS, CDS, HDS, PDS, Q-QUICK and UQ-QUICK. To this end, comprehensive comparisons for their numerical performances are investigated, including relative errors, CPU time and numerical stability.

3.1. Governing Equation and Initial Condition

In general, a two-dimensional advection-diffusion equation of pollutants can be written as


where, C is concentration of pollutants, u and v are velocities along x and y directions; Dx and Dy are diffusive coefficients;SØ is the source item.

The initial computational condition is governed by a unit Gauss impulse within a rectangular plane. It is depicted as follow:

Table 1. Definition of.


Its corresponding resolution is [15]


where, Dx = Dy= 0.01m2/s; u=v=0.8m/s; x and y are defined as m, m (see Figure 4(a)); the time step is set as s.

3.2. Numerical Discretization

Generally, the advection-diffusion equation of pollutants can be first integrated over a control volume and then using the Gauss’s divergence theorem, finally the common form over a CV is


where, is listed for time discretization;

, standing for the convection and diffusion fluxes through one interface f respectively, andis the surface outwardly normal vector.

In the following section, the first-order and second-order schemes for time discretization are both considered and the corresponding common equations for convenience of program compiling are also illuminated. In addition, over - relaxed approach [16−18] is adopted for cross derivative term approximation in the numerical computation.

First-order time discretization (UDS)

A common form for first-order time discretization along with different schemes for convective flux discretization can be derived as follow:


where,;, andare the splitting vectors;; has different form for each scheme and its definition is tabulated in Table 1.

Second-order time discretization (Crank-Nicolson)

A common form for Crank-Nicolson can be written as


where, different coefficients are as a follow

The definition of is the same as above, superscript ‘0’ denotes the last time step.

3.3. Numerical Results

The verification of Q-QUICK/UQ-QUICK has been investigated on three values of unstructured mesh numbers: 454, 1770 and 3955 with quadrilateral grids only. As it can be seen from Figure 4, the computational Grid 1 and Grid 2 are intentionally twisted along the pollutant downstream. Generalized Minimum Residual (GMRES (30)) method with the Incomplete LU (ILUT) precondition is used to accelerate the convergence of the linear equation. The numerical precision is indicated by relative errors defined as follow:


If t=1.25s, the relative errors and CPU time of different schemes are listed in Table 2, concentrations of each scheme along line y=x are also illuminated in Figure 5.

• Numerical precision: At the same scheme of time discretization, the numerical precision of Q-QUICK/ UQQUICK/ CDS are much better than UDS/HDS/PDS, as for Q-QUICK/ UQ-QUICK/ CDS themselves, Q-QUICK/ UQ-QUICK exhibit a little lower relative errors that those of CDS especially on sparse grids. In addition, Q-QUICK and UQ-QUICK show a similar numerical precision. Furthermore, with the increasing mesh numbers, QQUICK/ UQ-QUICK/CDS can fit well to the benchmarks quickly and on Grid 3 they can perform a perfectly match.

• CPU time: Q-QUICK and UQ-QUICK consume a slight longer CPU times than those of other schemes. As for comparison by them, Q-QUICK seems a little economical. This is reasonable for its little auxiliary nodes in scheme reconstruction.

• Numerical stability: After applying the deferred correction method for Q-QUICK/UQ-QUICK and over-relaxed approach for cross derivative term approximation, most of the schemes can keep a good numerical stability except that CDS has a little vibration on relative sparse grids. In addition, these

Figure 4. Computational domain and grids.

Table 2. Relative errors and CPU time.

schemes are in sensitive to high grid deformation.

For further clearly representing the numerical precision by each scheme, the distributions of concentration at 1.25s within the whole computational domain are also illustrated in Figure 6 (Grid 1, Crank-Nicolson). As referenced to analytic solution, the suggested scheme Q-QUICK or UQ-QUICK can keep pollutant cloud shape well more similar to that of benchmark, UDS and HDS numerically dissipate a lot and the corresponding shape almost distort, UDS exhibits better than the former two schemes, however, its shape shows not so matching especially in the left part.

4. Conclusions

In present study, two modified QUICK namely Q-QUICK and UQ-QUICK are introduced and their verifications were conducted on 2D unsteady advection-diffusion equation of pollutants. The main conclusions are as follow:

1) For comparisons in relative error and computational time amongst various schemes on unstructured grids, it can be seen that Q-QUICK and UQ-QUICK have similar computational accuracy to the central difference scheme and similar numerical stability to upwind difference scheme after applying the deferred correction method.

2) Q-QUICK and UQ-QUICK’s corresponding CPU times are approximately equivalent to those of traditional difference schemes and their abilities for adapting high grid deformation are robust.

Figure 5. Illustration of comparisons of concentration along line y=x at t=1.25s.

Figure 6. Illustration of comparisons of concentration in the plane at t=1.25s. (Grid1, Crank-Nicolson).

3) It is so promising to apply the suggested schemes to simulate pollutant transportation in shallow water with the merit of higher numerical precision.

5. Acknowledgments

This research was sponsored by Research Initiation Funds for Doctor (No.YJJ0906/SL01) of Changjiang River Scientific Research Institute ; 2009 National Public Research Institutes for Basic R & D Operating Expenses Special Project (No.YWF0901 and No.YWF0905); National Basic Research Program of China (No.2007CB714100, No. 2007CB714106).


  1. R. S~mkiewiez, “Oscillation-free solution of shallow water equations for nonstaggered grid,” Journal of Hydraulic Engineering, ASCE, Vol. 119, No. 10, pp. 1118−1137, 1993.
  2. J. Fletcher, “Computational techniques for fluid dyrtaralcs,” 8pfinger-Vedag, Berlin, No. 1, 1991.
  3. J. Shi and E. F. Toro, “Fully discrete high-order shock-eaptunng numerical schemes,” Internationa1 Journal for Numerical Methods in Fluids, Vol. 23, pp. 241−269, 1996.
  4. S. Sankaranarayanan, N. J. Shankar, and H. F. Cheoag, “Three-dimensional finite difference model for transport conservative pollutants,” Ocean Engneering, Vol. 25, No. 6, pp. 425−442, 1998.
  5. R. J. Sobey, “Fraetional step algorithm for estruaine mass transport,” International Journal for Numerical Methods in Fluids, Vol. 3, pp. 567−581, 1983.
  6. B. P. Leonard, “Simple high accuracy resolution program for convective modelling of discontinuities” International Journal for Numerical Methods in Fluids, Vol. 8, pp. 1291−1318, 1988.
  7. B. J. Noye and H. H. Tan, “A third-order semi-implicit finite difference method for solving one-dimensional eonvectlon-dlffusion equation,” International Journal for Numerical Method in Engineering, Vol. 26, pp. 1615− 1629, 1988.
  8. B. J. Noye and H. H. Tan, “Finite difference method for the two-dimensional advetion difusion equation” International Journal for Numerical Methods in Fluids, Vol. 9, pp. 75−98, 1989.
  9. A. Pollard, A. L. Siu, and L. W. Siu, “The calculation of some laminar flows using various discretization schemes,” Comp. Meth. Appl. Mech. Eng., Vol. 35, pp. 293−313, 1982.
  10. S. V. Pantankar, Numerical Heat Transfer, McGraw−Hill, New York, 1980.
  11. M. K. Patel and N. C. Markatos, “An evaluation of eight discretization schemes for two-dimensional convection-diffusion equations,” International Journal for Numerical Methods in Fluids, Vol. 6, pp. 129−154, 1986.
  12. L. Davidson, “A pressure correction method for unstructured meshes with arbitrary control volumes,” International Journal for Numerical Methods in Fluids, Vol. 22, pp. 265−281, 1996.
  13. Y. Saad and M. H. Schultz, “Gmres: A generalized minimal residual algorithm for solving nonsymmetric linear systems,” SIAM. J. Sci. Stat. Comput., Vol. 7, 1986.
  14. B. P. Leonard, “Third-order finite-difference method for steady two-dimensional convection,” Numerical Methods in Laminar and Turbulent Flow, pp. 807−819, 1978.
  15. S. I. Karaa and J. Zhan, “High order ADI method for solving unsteady convection diffusion problems,” Journal
  16. of Computational Physics, Vol. 198, No. 1, pp. 1−9, 2004.
  17. H. Jasak, “Error analysis and estimation for the finite volume method with applications to fluid flows,” Ph.D. thesis, Imperial College, University of London, London, 1996.
  18. B. Basara, “Employment of the second-moment turbulence closure on arbitrary unstructured grids,” International Journal of Numerical Methods in Fluids, Vol. 44, pp. 377−407, 2001.
  19. Y. Y. Tsui and Y. F. Pan, “A pressure-correction method for incompressible flows using unstructured methes,” Numerical Heat Transfer, Part B, Vol. 49, pp. 43−65, 2006.