Optimal upwinding of N+2 Petrov-Galerkin method based upon frequency fitting

Abstract

The N+2 Petrov-Galerkin method employs weighting functions which are two degree order higher than basis functions for trial solution. Recent studies have shown that N+2 Petrov-Galerkin method yields promising results, however, determination of an appropriate upwinding level is still problematic. In this paper, optimal upwinding of N+2 Petrov-Galerkin method is investigated systematically by using fitting algorithm. Introduction It is quite well known that standard finite element method yields poor solutions when it is applied to convection-dominatedtransport problem. A great deal of effort has been made to overcome such numerical difficulty (Hughes, 1978 ; Hughes, 1979 ; Celia et al., 1990), however, problem is still under being challenged. Since numerical errors arise from symmetric treatment of convection terms in conventional finite element method, weighting functions are modified in order to give more weight in upwind direction as commonly done in finite difference Transactions on Ecology and Environment vol 17, © 1998 WIT Press, www.witpress.com, ISSN 1743-3541 666 Computer Methods in Water Resources XII This technique is referred to as the Petrov-Galerkin (PG) method. The N+& PG method uses weighting functions which are k polynomial degree higher than basis functions for trial solution. The N+l PG method, having quadratic perturbation term in weighting functions, has been successfully applied to steady-state convection-dominated transport problem (Christie et al., 1976 ; Heinrich et al., 1977). However, it failed to improve numerical solutions for time-dependent problem, which necessitated introduction of N+2 PG The weighting functions of N+2 PG method include both asymmetric quadratic and symmetric cubic perturbation terms. Although it is true that N+2 PG method alleviates dispersive errors without introducing too much dissipation, a weakness of this method lies in determining appropriate levels of upwinding. Specially, numerical experiments have shown that symmetric cubic term rather than asymmetric quadratic term in weighting functions plays a major role in reducing numerical errors due to time dependency (Dick, 1983 ; Weterink and Shea, 1989 ; Bouloutas and Celia, 1991). That is, using both quadratic and cubic terms in N+2 PG method deteriorates numerical solutions in time-dependent problem. Tezduyar and Ganjoo (1986) proposed perturbation functions dependent upon Courant number in streamline upwind PG Westerink and Shea (1989) presented a graphical method to determine upwinding parameters for cubic perturbation function for pure convection problem. Miller and Cornew (1992) showed that standard deviation of a Gaussian source affects optimal upwinding level in N+2 PG Carrano and Yeh (1994) introduced a spectral weighting technique to optimize phase error in N+2 PG In this paper, choosing appropriate upwinding parameters in N+2 PG method is investigated. By using Fourier expansion, enlargement factors of analytical and numerical solutions are compared to find optimal values for upwinding parameters. This procedure is termed as frequency fitting algorithm in Bouloutas and Celia (1991). This study will reveal role of each perturbation term in N+2 weighting function explicitly, and help in choosing appropriate upwinding levels in N+2 PG