A spectral collocation method for compressible,non-similar boundary layers

An efficient and highly accurate algorithm based on a spectral collocation method is developed for numerical solution of the compressible, two-dimensional and axisymmetric boundary layer equations. The numerical method incorporates a fifth-order, fully implicit marching scheme in the streamwise (timelike) dimension and a spectral collocation method based on Chebyshev polynomial expansions in the wall-normal (spacelike) dimension. The discrete governing equations are cast in residual form and the residuals are minimized at each marching step by a preconditioned Richardson iteration scheme which fully couples energy, momentum and continuity equations. Preconditioning on the basis of the finite difference analogues of the governing equations results in a computationally efficient iteration with acceptable convergence properties. A practical application of the algorithm arises in the area of compressible linear stability theory, in the investigation of the effects of transverse curvature on the stability of flows over axisymmetric bodies. The spectral collocation algorithm is used to derive the non-similar mean velocity and temperature profiles in the boundary layer of a ‘fuselage’ (cylinder) in a high-speed (Mach 5) flow parallel to its axis. The stability of the flow is shown to be sensitive to the gradual streamwise evolution of the mean flow and it is concluded that the effects of transverse curvature on stability should not be ignored routinely.     

On the efficiency of semi‐implicit and semi‐Lagrangian spectral methods for the calculation of incompressible flows

Classical semi‐implicit backward Euler/Adams–Bashforth time discretizations of the Navier–Stokes equations induce, for high‐Reynolds number flows, severe restrictions on the time step. Such restrictions can be relaxed by using semi‐Lagrangian schemes essentially based on splitting the full problem into an explicit transport step and an implicit diffusion step. In comparison with the standard characteristics method, the semi‐Lagrangian method has the advantage of being much less CPU time consuming where spectral methods are concerned. This paper is devoted to the comparison of the ‘semi‐implicit’ and ‘semi‐Lagrangian’ approaches, in terms of stability, accuracy and computational efficiency. Numerical results on the advection equation, Burger's equation and finally two‐ and three‐dimensional Navier–Stokes equations, using spectral elements or a collocation method, are provided. Copyright © 2001 John Wiley & Sons, Ltd.     

Single‐degree of freedom Hermite collocation for multiphase flow and transport in porous media

The classical collocation method using Hermite polynomials is computationally expensive as the dimensionality of the problem increases. Because of the use of a C1‐continuous basis, the method generates two, four and eight unknowns per node for one, two and three‐dimensional problems, respectively. In this paper we propose a numerical strategy to reduce the nodal unknowns to a single degree of freedom at each node. The reduction of the unknowns is due to the use of Lagrangian polynomials to approximate the first‐order derivatives over the minimal compact stencil surrounding each node. For the solvability of the problem the reduction of the number of collocation equations is done by a nodal weighting strategy. We have applied the proposed approach to enhance the efficiency of a collocation‐based multiphase flow and transport simulator. Benchmark cases illustrate the higher performance of the new methodology when compared to classical Hermite collocation. Copyright © 2004 John Wiley & Sons, Ltd.     

A Fourier–Chebyshev spectral collocation method for simulating flow past spheres and spheroids

An accurate Fourier–Chebyshev spectral collocation method has been developed for simulating flow past prolate spheroids. The incompressible Navier–Stokes equations are transformed to the prolate spheroidal co‐ordinate system and discretized on an orthogonal body fitted mesh. The infinite flow domain is truncated to a finite extent and a Chebyshev discretization is used in the wall‐normal direction. The azimuthal direction is periodic and a conventional Fourier expansion is used in this direction. The other wall‐tangential direction requires special treatment and a restricted Fourier expansion that satisfies the parity conditions across the poles is used. Issues including spatial and temporal discretization, efficient inversion of the pressure Poisson equation, outflow boundary condition and stability restriction at the pole are discussed. The solver has been validated primarily by simulating steady and unsteady flow past a sphere at various Reynolds numbers and comparing key quantities with corresponding data from experiments and other numerical simulations. Copyright © 1999 John Wiley & Sons, Ltd.     

配点型广义节点无网格法的基本原理及其应用

借鉴流形方法思想,引入广义节点的概念,对传统的无网格法进行了改进,建立了可具有任意高阶多项式插值函数的广义节点无网格方法。同时采用径向插值函数构造具有插值特性的逼近函数;采用配点法建立系统的离散方程。在阐述了这种方法基本原理的同时,针对线弹性力学问题给出了这种方法的数值计算列式。与传统无网格方法相比,这种方法更具有一般性;同时由于采用了配点法而不需要背景积分网格,所以可以认为这种方法是某种真正意义上的无网格法。当选取0阶广义节点位移插值函数时便可得到传统的无网格法;在不增加支持域内节点数目的条件下,通过选取高阶广义节点位移插值函数可以提高计算精度。最后通过算例分析,对0阶、1阶及2阶广义节点无网格法与现有的有关解答进行了对比,论证了其合理性。     

Three‐dimensional transient Navier–Stokes solvers in cylindrical coordinate system based on a spectral collocation method using explicit treatment of the pressure

A spectral collocation method is developed for solving the three‐dimensional transient Navier–Stokes equations in cylindrical coordinate system. The Chebyshev–Fourier spectral collocation method is used for spatial approximation. A second‐order semi‐implicit scheme with explicit treatment of the pressure and implicit treatment of the viscous term is used for the time discretization. The pressure Poisson equation enforces the incompressibility constraint for the velocity field, and the pressure is solved through the pressure Poisson equation with a Neumann boundary condition. We demonstrate by numerical results that this scheme is stable under the standard Courant–Friedrichs–Lewy (CFL) condition, and is second‐order accurate in time for the velocity, pressure, and divergence. Further, we develop three accurate, stable, and efficient solvers based on this algorithm by selecting different collocation points in r‐, ? ‐, and z‐directions. Additionally, we compare two sets of collocation points used to avoid the axis, and the numerical results indicate that using the Chebyshev Gauss–Radau points in radial direction to avoid the axis is more practical for solving our problem, and its main advantage is to save the CPU time compared with using the Chebyshev Gauss–Lobatto points in radial direction to avoid the axis. Copyright © 2010 John Wiley & Sons, Ltd.     

A Chebyshev collocation method for the Navier–Stokes equations with application to double-diffusive convection

A Chebyshev collocation method for solving the unsteady two-dimensional Navier–Stokes equations in vorticity–streamfunction variables is presented and discussed. The discretization in time is obtained through a class of semi-implicit finite difference schemes. Thus at each time cycle the problem reduces to a Stokes-type problem which is solved by means of the influence matrix technique leading to the solution of Helmholtz-type equations with Dirichlet boundary conditions. Theoretical results on the stability of the method are given. Then a matrix diagonalization procedure for solving the algebraic system resulting from the Chebyshev collocation approximation of the Helmholtz equation is developed and its accuracy is tested. Numerical results are given for the Stokes and the Navier–Stokes equations. Finally the method is applied to a double-diffusive convection problem concerning the stability of a fluid stratified by salinity and heated from below.     

Fifth‐order Hermitian schemes for computational linear aeroacoustics

We develop a class of fifth‐order methods to solve linear acoustics and/or aeroacoustics. Based on local Hermite polynomials, we investigate three competing strategies for solving hyperbolic linear problems with a fifth‐order accuracy. A one‐dimensional (1D) analysis in the Fourier series makes it possible to classify these possibilities. Then, numerical computations based on the 1D scalar advection equation support two possibilities in order to update the discrete variable and its first and second derivatives: the first one uses a procedure similar to that of Cauchy–Kovaleskaya (the ‘Δ‐P5 scheme’); the second one relies on a semi‐discrete form and evolves in time the discrete unknowns by using a five‐stage Runge–Kutta method (the ‘RGK‐P5 scheme’). Although the RGK‐P5 scheme shares the same local spatial interpolator with the Δ‐P5 scheme, it is algebraically simpler. However, it is shown numerically that its loss of compactness reduces its domain of stability. Both schemes are then extended to bi‐dimensional acoustics and aeroacoustics. Following the methodology validated in (J. Comput. Phys. 2005; 210 :133–170; J. Comput. Phys. 2006; 217 :530–562), we build an algorithm in three stages in order to optimize the procedure of discretization. In the ‘reconstruction stage’, we define a fifth‐order local spatial interpolator based on an upwind stencil. In the ‘decomposition stage’, we decompose the time derivatives into simple wave contributions. In the ‘evolution stage’, we use these fluctuations to update either by a Cauchy–Kovaleskaya procedure or by a five‐stage Runge–Kutta algorithm, the discrete variable and its derivatives. In this way, depending on the configuration of the ‘evolution stage’, two fifth‐order upwind Hermitian schemes are constructed. The effectiveness and the exactitude of both schemes are checked by their applications to several 2D problems in acoustics and aeroacoustics. In this aim, we compare the computational cost and the computation memory requirement for each solution. The RGK‐P5 appears as the best compromise between simplicity and accuracy, while the Δ‐P5 scheme is more accurate and less CPU time consuming, despite a greater algebraic complexity. Copyright © 2007 John Wiley & Sons, Ltd.     

On the accurate prediction of the wall-normal velocity in compressible boundary-layer flow

We consider a problem which arises in the numerical solution of the compressible two-dimensional or axisymmetric boundary-layer equations. Numerical methods for the compressible boundary-layer equations are facilitated by transformation from the physical (x, y) plane to a computational (ξ, η) plane in which the evolution of the flow is ‘slow’ in the time-like ξ direction. The commonly used Levy-Lees transformation results in a computationally well-behaved problem, but it complicates interpretation of the solution in physical space. Specifically, the transformation is inherently non-linear, and the physical wall-normal velocity is transformed out of the problem and is not readily recovered. Conventional methods extract the wall-normal velocity in physical space from the continuity equation, using finite-difference techniques and interpolation procedures. The present spectrally accurate method extracts the wall-normal velocity directly from the transformation itself, without interpolation, leaving the continuity equation free as a check on the quality of the solution. The present method for recovering wall-normal velocity, when used in conjunction with a highly accurate spectral collocation method for solving the compressible boundary-layer equations, results in a discrete solution which satisfies the continuity equation nearly to machine precision. As demonstration of the utility of the method, the boundary layers of three prototypical high-speed flows are investigated and compared: the flat plate, the hollow cylinder, and the cone. An important implication for classical linear stability theory is also briefly discussed.     

含肿瘤皮肤组织传热分析的广义有限差分法

生物传热分析在低温外科手术、肿瘤热疗、病热诊断等临床医学治疗和诊断中有着广泛的应用. 由于健康皮肤组织内肿瘤的存在使得肿瘤附近区域的温度会明显升高, 这一特性常被用于检测皮肤组织内的肿瘤生长, 因此有必要开展生物传热数值分析的研究. 本文以含肿瘤的皮肤组织为研究对象, 将一种新型区域型无网格配点法——广义有限差分法应用于能描述含肿瘤皮肤组织传热过程的Pennes方程求解. 广义有限差分法利用泰勒展开式与移动最小二乘法将计算区域内的每个离散点上的物理量导数表示成其与邻近点物理量及权重系数的线性组合, 进而构建得到仅含各离散点未知物理量的线性方程组. 该方法不仅具有无需划分网格、避免数值积分等无网格配点法的优点, 同时还克服了大多数无网格配点法中插值矩阵高度病态的问题, 为此类方法在大规模工程数值计算中的应用提供了可能性. 本文首先介绍了模拟含肿瘤皮肤组织传热过程的广义有限差分法离散模型, 随后通过不含肿瘤与含规则形状肿瘤的基准算例, 检验广义有限差分法的计算精度与收敛性; 在此基础上, 通过数值模拟研究不同肿瘤形状及肿瘤位置分布对皮肤组织内温度分布的影响.      

Scattering of plane monochromatic waves from a heterogeneous inclusion of arbitrary shape in a poroelastic medium: An efficient numerical solution

Scattering of plane longitudinal monochromatic waves from a heterogeneous inclusion of arbitrary shape in an infinite poroelastic medium is considered. Wave propagation in the medium is described by Biot’s equations of poroelasticity. The scattering problem is formulated in terms of the volume integral equations for displacements of the solid skeleton and fluid pressure in the pore space in the region occupied by the inclusion. An efficient numerical method is applied to solve these equations. In the method, Gaussian approximating functions are used for discretization of the problem. For regular node grids, the matrix of the discretized problem has Toeplitz’s properties, and the Fast Fourier Transform technique can be used for the calculation of matrix–vector products. The latter accelerates substantially the process of iterative solution of the discretized problem. For material parameters of typical sedimentary rocks, the system of differential equations of poroelasticity contains a differential operator with a small parameter. As the result, the wave field in the inclusion region is split up into a slowly changing part, and boundary layer functions concentrated near the inclusion interface. The method of matched asymptotic expansions is used for the numerical solution in this case. For a spherical inclusion, the results of the numerical and matched asymptotic expansion methods are compared with a semi-analytical series solution. For a non-spherical heterogeneous inclusion, an example of the numerical solution is presented.     

自由面势流问题的域外奇点边界元法及其数值误差分析

讨论了域外奇点边界元法在自由面势流问题计算中的作用，并以连续及离散Fourier分析对该方法(就m阶面元的一般情况)进行数值误差分析，导出了计及面元阶数、奇点至自由面垂向距离、配置点移动、差分格式等因素影响的数值误差一般表达式。从理论上证明了自由面势流问题计算中采用域外奇点法可改善离散产生的数值色散误差并能结合配置点前移(向上游)等方法以数值满足辐射条件。     

边界节点法的计算稳定性研究

边界节点法利用满足控制方程的非奇异通解作为基函数,半解析边界数值离散偏微分方程,具有精度高、收敛快、易编程等优点,是一种纯无网格配点方法.但是在求解具体问题时,随着节点数的增加,边界节点法经常得到严重病态的插值矩阵.本文利用有效条件数评价边界节点法求解Helmholtz问题线性方程组的计算稳定性;然后利用三种正则化方法处理其病态的线性方程组,并与高斯消元法比较计算精度和收敛性.通过数值实验,本文研究了有效条件数、误差和正则化方法之间的关系.     

Parallel ILU preconditioning,a priori pivoting and segregation of variables for iterative solution of the mixed finite element formulation of the Navier–Stokes equations

A parallel ILU preconditioning algorithm for the incompressible Navier–Stokes equations has been designed, implemented and tested. The computational mesh is divided into N subdomains which are processed in parallel in different processors. During ILU factorization, matrices and vectors associated with the nodes on the interface between the subdomains are communicated to the equation matrices to the adjacent subdomain. The bases for the parallel algorithm are an appropriate node ordering scheme and a segregation of velocity and pressure degrees of freedom. The inner nodes of the subdomain are numbered first and then the nodes on the interface between the subdomains. To avoid division by zero during the ILU factorization, the equations corresponding to the velocity degrees of freedom are assembled first in the global equation matrix, followed by the equations corresponding to the pressure degrees of freedom. Copyright © 2003 John Wiley & Sons, Ltd.     

Acoustic simulation using a novel approach for reducing dispersion error

It is well‐known that the traditional finite element method (FEM) fails to provide accurate results to the Helmholtz equation with the increase of wave number because of the ‘pollution error’ caused by numerical dispersion. In order to overcome this deficiency, a gradient‐weighted finite element method (GW‐FEM) that combines Shepard interpolation and linear shape functions is proposed in this work. Three‐node triangular and four‐node tetrahedral elements that can be generated automatically are first used to discretize the problem domain in 2D and 3D spaces, respectively. For each independent element, a compacted support domain is then formed based on the element itself and its adjacent elements sharing common edges (or faces). With the aid of Shepard interpolation, a weighted acoustic gradient field is then formulated, which will be further used to construct the discretized system equations through the generalized Galerkin weak form. Numerical examples demonstrate that the present algorithm can significantly reduces the dispersion error in computational acoustics. Copyright © 2016 John Wiley & Sons, Ltd.     

Collocated discrete least squares meshless (CDLSM) method for the solution of transient and steady‐state hyperbolic problems

A collocated discrete least squares meshless method for the solution of the transient and steady‐state hyperbolic problems is presented in this paper. The method is based on minimizing the sum of the squared residuals of the governing differential equation at some points chosen in the problem domain as collocation points. The collocation points are generally different from nodal points, which are used to discretize the problem domain. A moving least squares method is employed to construct the shape functions at nodal points. The coefficient matrix is symmetric and positive definite even for non‐symmetric hyperbolic differential equations and can be solved efficiently with iterative methods. The proposed method is a truly meshless method and does not require numerical integration. Advantages of the collocation points are shown to be threefold: First, the collocation points are shown to be responsible for stabilizing the method in particular when problems with shocked solution are attempted. Second, the collocation points are also shown to improve the accuracy of the solution even for problems with smooth solutions. Third, the collocation points are shown to contribute to the efficiency of the method when solving steady‐state problems via faster convergence of the resulting algorithm. The ability of the method and in particular the effect of collocation points are tested against a series of one‐dimensional transient and steady‐state benchmark examples from the literature and the results are presented. A sensitivity analysis is also carried out to investigate the effect of the base polynomials on the accuracy and convergence characteristics of the method in solving steady‐state problems. The results show the ability of the proposed method to accurately solve difficult hyperbolic problems considered. The method is also shown to be particularly stable for problems with shocked solution due to the inherent stabilizing mechanism of the method. Copyright © 2008 John Wiley & Sons, Ltd.     

A numerical scheme for solving a periodically forced Reynolds equation

A modified Reynolds equation is used to model the air‐film in a high‐speed squeeze‐film bearing. The axial position of the bearing stator is prescribed as a finite amplitude periodic oscillation. A numerical approach is considered for solving the uncoupled and coupled periodic problems associated with this model. The uncoupled problem requires the computation of the squeeze‐film dynamics when the rotor is held at a fixed axial position and the coupled problem incorporates the additional air–rotor interaction since the rotor position is unknown and modelled as a spring‐mass‐damper system. The details of a Fourier spectral collocation scheme are provided for the reduction of the modified Reynolds equation to a system of non‐linear, first‐order ordinary differential equations in space. Using the Matlab boundary value problem solver bvp4c this system of equations is solved to give the periodic pressure distributions and rotor heights. The high degree of accuracy in the spectral collocation scheme is demonstrated through comparison with an appropriate analytical solution. Further analysis indicates that the direct periodic solver is at least 10 times faster than the equivalent Crank–Nicholson finite‐difference scheme. For changing values of a selected physical parameter the method of arc‐length continuation is employed to track branches of solutions computed using the spectral collocation scheme. A selection of results is presented to demonstrate the range of accessible solutions and the robust nature of the numerical scheme. Copyright © 2010 John Wiley & Sons, Ltd.     

Runge–Kutta methods for a semi-analytical prediction of milling stability

On the basis of Runge–Kutta methods, this paper proposes two semi-analytical methods to predict the stability of milling processes taking a regenerative effect into account. The corresponding dynamics model is concluded as a coefficient-varying periodic differential equation with a single time delay. Floquet theory is adopted to predict the stability of machining operations by judging the eigenvalues of the state transition matrix. This paper firstly presents the classical fourth-order Runge–Kutta method (CRKM) to solve the differential equation. Through numerical tests and analysis, the convergence rate and the approximation order of the CRKM is not as high as expected, and only small discrete time step size could ensure high computation accuracy. In order to improve the performance of the CRKM, this paper then presents a generalized form of the Runge–Kutta method (GRKM) based on the Volterra integral equation of the second kind. The GRKM has higher convergence rate and computation accuracy, validated by comparisons with the semi-discretization method, etc. Stability lobes of a single degree of freedom (DOF) milling model and a two DOF milling model with the GRKM are provided in this paper.     

An Approach to Transport in Heterogeneous Porous Media Using the Truncated Temporal Moment Equations: Theory and Numerical Validation

In the last decade, the characterization of transport in porous media has benefited largely from numerical advances in applied mathematics and from the increasing power of computers. However, the resolution of a transport problem often remains cumbersome, mostly because of the time-dependence of the equations and the numerical stability constraints imposed by their discretization. To avoid these difficulties, another approach is proposed based on the calculation of the temporal moments of a curve of concentration versus time. The transformation into the Laplace domain of the transport equations makes it possible to develop partial derivative equations for the calculation of complete moments or truncated moments between two finite times, and for any point of a bounded domain. The temporal moment equations are stationary equations, independent of time, and with weaker constraints on their stability and diffusion errors compared to the classical advection–dispersion equation, even with simple discrete numerical schemes. Following the complete theoretical development of these equations, they are compared firstly with analytical solutions for simple cases of transport and secondly with a well-performing transport model for advective–dispersive transport in a heterogeneous medium with rate-limited mass transfer between the free water and an immobile phase. Temporal moment equations have a common parametrization with transport equations in terms of their parameters and their spatial distribution on a grid of discretization. Therefore, they can be used to replace the transport equations and thus accelerate the achievement of studies in which a large number of simulations must be carried out, such as the inverse problem conditioned with transport data or for forecasting pollution hazards.     

On the practical importance of the SSP property for Runge–Kutta time integrators for some common Godunov‐type schemes

We investigate through analysis and computational experiment explicit second and third‐order strong‐stability preserving (SSP) Runge–Kutta time discretization methods in order to gain perspective on the practical necessity of the SSP property. We consider general theoretical SSP limits for these schemes and present a new optimal third‐order low‐storage SSP method that is SSP at a CFL number of 0.838. We compare results of practical preservation of the TVD property using SSP and non‐SSP time integrators to integrate a class of semi‐discrete Godunov‐type spatial discretizations. Our examples involve numerical solutions to Burgers' equation and the Euler equations. We observe that ‘well‐designed’ non‐SSP and non‐optimal SSP schemes with SSP coefficients less than one provide comparable stability when used with time steps below the standard CFL limit. Results using a third‐order non‐TVD CWENO scheme are also presented. We verify that the documented SSP methods with the number of stages greater than the order provide a useful enhanced stability region. We show by analysis and by numerical experiment that the non‐oscillatory third‐order reconstructions used in (Liu and Tadmor Numer. Math. 1998; 79 :397–425, Kurganov and Petrova Numer. Math. 2001; 88 :683–729) are in general only second‐ and first‐order accurate, respectively. Copyright © 2005 John Wiley & Sons, Ltd.