Higher-order accurate implicit time integration schemes for transport problems

The present paper is concerned with the numerical solution of transient transport problems by means of spatial and temporal discretization methods. The generalized initial boundary value problem of various nonlinear transport phenomena like heat transfer or mass transport is discretized in space by p-finite elements. After finite element discretization, the resulting first-order semidiscrete balance has to be solved with respect to time. Next to the classical generalized-α integration method predicated on the Newmark approach and the evaluation at a generalized midpoint also implicit Runge–Kutta time integration schemes, are presented. Both families of finite difference-based integration schemes are derived for general first-order problems. In contrast to the above-mentioned algorithms, temporal discontinuous and continuous Galerkin methods evaluate the balance equation not at a selected time instant within the timestep, but in an integral sense over the whole time step interval. Therefore, the underlying semidiscrete balance and the continuity of the primary variables are weakly formulated within time steps and between time steps, respectively. Continuous Galerkin methods are obtained by the strong enforcement of the continuity condition as special cases. The introduction of a natural time coordinate allows for the application of standard higher-order temporal shape functions of the p-Lagrange type and the well-known Gau?–Legendre quadrature of associated time integrals. It is shown that arbitrary order accurate integration schemes can be developed within the framework of the proposed temporal p-Galerkin methods. Selected benchmark analyses of calcium diffusion demonstrate the properties of all three methods with respect to non-smooth initial or boundary conditions. Furthermore, the robustness of the present time integration schemes is also demonstrated for the highly nonlinear reaction–diffusion problem of calcium leaching, including the pronounced changes of the reaction term and non-smooth changes of Dirichlet boundary conditions of calcium dissolution.     

Acoustoelastic effect simulation by time–space finite element formulation based on quadratic interpolation of the acceleration

Acoustoelastic effect describes the change of ultrasound velocity due to the initial stress. Its simulation involves a numerical analysis of nonlinear elastodynamics and requires high accuracy in the time domain. A time–space finite element formulation, derived from the quadratic interpolation of the acceleration within a time segment, is proposed for an accurate simulation of the acoustoelastic effect in the present study. Ten different integration schemes are generated based on this formulation and nine of them are found to be conditionally stable. Among the nine stable schemes, one is found to obtain a spectral radius of one when the normalized step ratio is less than 5.477, indicating no numerical dissipation or numerical divergence. Compared with integration schemes from previous studies, this integration scheme demonstrates better performance in calculation accuracy and energy conservation. A two-stage approach, namely the static stage and the dynamic stage, has been employed in the simulation of the acoustoelastic effect. The former stage is adopted to obtain the initial stress and the latter stage, where the proposed integration scheme is implemented, is adopted to simulate the ultrasound propagation in an initial stress state. The simulation results of the dynamic stage show that the ultrasound velocity increases in a compression stress state and decreases in a tension stress state for aluminum alloy, which is in good agreement with previous experimental studies. Together with the simulation result of the static stage, it is conjectured that the acoustoelastic effect results from the stress-dependent elastic modulus.     

A finite element method for high Reynolds number viscous fluid flow using two step explicit scheme

This paper presents the finite element method for the analysis of unsteady viscous flow of fluid at high Reynolds numbers. The method is based on the explicit numerical integration scheme in time and uses three node triangular finite elements. For the convenience of the formulation, slight compressibility is considered. For the explicit scheme, the selective lumping two step scheme has been successfully employed. Vortex shedding behind a cylinder has been computed and compared with the conventional experimental results. The results agree favourably when both schemes are compared.     

精细积分方法的评估与改进

详细分析了结构动力分析的精细积分方法的稳定性、计算精度，在此基础上提出了对现有精细积分方法的改进策略。算例证实了本文对精细积分方法改进的科学性与可行性。     

瞬态热传导问题的精细积分数值流形方法研究

数值流形方法（NMM）因其特有的双覆盖系统（数学覆盖和物理覆盖）在域离散方面具有独特的优势，而精细时间积分法则具有精度高、无条件稳定、无振荡以及计算结果不依赖于时间步长等特点。发展了用于研究二维瞬态热传导问题的精细积分NMM。结合待求问题的控制方程和边界条件，并基于修正变分原理导出了NMM的总体方程，给出了求解此类时间相依方程的精细时间积分及空间积分策略，选取了两个典型算例对方法的有效性进行了验证，结果表明本文方法可以高效高精度地求解瞬态热传导问题。     

THE COMPACTION OF TIME-DEPENDENT VISCOUS GRANULAR MATERIALS CONSIDERING INERTIAL FORCES

In the present paper, compactions of time-dependent viscous granular materials are simulated step by step using the automatic adaptive mesh generation schemes. Inertial forces of the viscous incompressible aggregates are taken into account. The corresponding conservation equations, the weighted-integral formulations, and penalty finite element model are investigated. The fully discrete finite element equations for the simulation are derived. Polygonal particles of aggregates are simplified as mixed three-node and four-node elements. The automatic adaptive mesh generation schemes include contact detection algorithms, and mesh upgrade schemes. Solutions of the numerical simulation are in good agreement with some results from literatures. With minor modification, the proposed numerical model can be applied in several industries, including the pharmaceutical, ceramic, food, and household product manufacturing.     

无条件稳定的更新精细积分方法

提出将Pade逼近与精细积分方法中的指数矩阵运算技巧结合起来,建立了精细积分法的更新形式及计算过程,对该更新精细积分方法的稳定性进行了论证与探讨.结果表明,该更新精细积分方法是无条件稳定的,整个积分方法的精度取决于所取Pade逼近的阶数与高斯积分点的数量.数值例题也显示了该方法的高效率及其可行性.     

Spectral/hp penalty least‐squares finite element formulation for unsteady incompressible flows

In this paper, we present spectral/hp penalty least‐squares finite element formulation for the numerical solution of unsteady incompressible Navier–Stokes equations. Pressure is eliminated from Navier–Stokes equations using penalty method, and finite element model is developed in terms of velocity, vorticity and dilatation. High‐order element expansions are used to construct discrete form. Unlike other penalty finite element formulations, equal‐order Gauss integration is used for both viscous and penalty terms of the coefficient matrix. For time integration, space–time decoupled schemes are implemented. Second‐order accuracy of the time integration scheme is established using the method of manufactured solution. Numerical results are presented for impulsively started lid‐driven cavity flow at Reynolds number of 5000 and transient flow over a backward‐facing step. The effect of penalty parameter on the accuracy is investigated thoroughly in this paper and results are presented for a range of penalty parameter. Present formulation produces very accurate results for even very low penalty parameters (10–50). Copyright © 2008 John Wiley & Sons, Ltd.     

精细积分法在迷宫密封转子系统非线性动力学特性计算中的应用

基于Muszynska密封力模型，建立了迷宫密封转子系统的非线性动力学模型，将精细积分法推广应用于非线性情况，计算了迷宫密封不平衡转子系统的动力学特性，依据Floquet理论讨论其分岔特性。研究表明：在2^N类算法计算指数矩阵基础上提出的精细积分法和传统的数值计算方法相比，其精度高，在分析中通过取不同步长计算对比，表明该方法在某些情况下可以采取较大时间步长，有效提高了计算速度。     

结构动力方程的更新精细积分方法

将高斯积分方法与精细积分方法中的指数矩阵运算技巧结合起来,建立了精细积分法的更新形式及计算过程,对该更新精细积分方法的稳定性进行了论证与探讨。在实施精细积分过程中不必进行矩阵求逆,整个积分方法的精度取决于所选高斯积分点的数量。这种方法理论上可实现任意高精度,计算效率较高,其稳定性条件极易满足。数值例题也显示了这种方法的有效性。     

二维定常不可压缩粘性流动N-S方程的数值流形方法

将流形方法应用于定常不可压缩粘性流动N-S方程的直接数值求解,建立基于Galerkin加权余量法的N-S方程数值流形格式,有限覆盖系统采用混合覆盖形式,即速度分量取1阶和压力取0阶多项式覆盖函数,非线性流形方程组采用直接线性化交替迭代方法和Nowton-Raphson迭代方法进行求解.将混合覆盖的四节点矩形流形单元用于阶梯流和方腔驱动流动的数值算例,以较少单元获得的数值解与经典数值解十分吻合.数值实验证明,流形方法是求解定常不可压缩粘性流动N-S方程有效的高精度数值方法.     

二维CE/SE算法改进及其高阶精度格式

通过提出一种新的守恒元和解元划分方式对二维时-空守恒元解元算法(CE/SE)进行了改进,推导了改进CE/SE算法的一、二阶精度计算格式,并给出了更高阶精度计算格式的构造方法。利用得到的改进CE/SE格式对激波反射问题、后台阶扰流问题及隔墙坑道传播问题进行了数值模拟。数值结果表明,对CE/SE算法的改进是成功的。改进CE/SE算法有诸多优点,值得在数值模拟中推广使用。     

精细时程积分法及其数值衍生格式应用评估

旋翼气动弹性耦合动力学方程本质上是一组刚性比较大的非线性偏微分方程。在有限元结构离散后，可改写为非齐次微分方程组，其中非齐次项是桨叶运动量（位移与速度）和气动载荷的函数。针对这类方程，本文尝试引入精细积分法及其衍生格式，借助数值方法计算Duhamel积分项。从积分精度与数值稳定性方面比较研究具有代表性的精细库塔法和高精度直接积分法。结合隐式积分算法，评估精细积分法应用于旋翼动力学方程的可行性。算例表明，精细积分法对矩形直桨叶动力学方程具有足够的求解精度。     

Hopping numerical approximations of the hyperbolic equation

Polynomial functions can be used to derive numerical schemes for an approximate solution of hyperbolic equations. A conventional derivation technique requires a polynomial to pass through every node values of a continuous computational stencil, leading to severe manifestation of the Gibbs phenomenon and strict time‐step limitation. To overcome the problem, this paper introduces polynomials that skip regularly (‘hop’ over) one or more nodes from the computational grid. Polynomials hopping over odd and even nodes yield a series of explicit numerical schemes of a required accuracy, with Lax–Friedrichs method being a particular simplest case. The schemes have two times wider stability interval compared to conventional continuous‐stencil explicit methods. Convex combinations of odd‐ and even‐node‐based updates improve further accuracy and stability of the method. Out of considered combinations (up to third‐order accuracy), derived odd‐order methods are stable for the Courant number ranging from 0 to 3, and even‐order ones from 0 to 5. A 2‐D extension of the hopping polynomial method exhibits similar properties. Copyright © 2007 John Wiley & Sons, Ltd.     

基于精细积分技术的非线性动力学方程的同伦摄动法

将精细积分技术（PIM）和同伦摄动方法（HPM）相结合，给出了一种求解非线性动力学方程的新的渐近数值方法。采用精细积分法求解非线性问题时，需要将非线性项对时间参数按Taylor级数展开，在展开项少时，计算精度对时间步长敏感；随着展开项的增加，计算格式会变得越来越复杂。采用同伦摄动法，则具有相对筒单的计算格式，但计算精度较差，应用范围也限于低维非线性微分方程。将这两种方法相结合得到的新的渐近数值方法则同时具备了两者的优点，既使同伦摄动方法的应用范围推广到高维非线性动力学方程的求解，又使精细积分方法在求解非线性问题时具有较简单的计算格式。数值算例表明，该方法具有较高的数值精度和计算效率。     

A consistent numerical scheme for the von Mises mixed-hardening constitutive equations

Analytical solutions are derived for the von Mises mixed-hardening elastoplastic model under rectilinear strain paths, and the concept of response subspace is introduced such that the original five-dimensional problem in deviatoric stress space is reduced to a more economic two-dimensional problem, of which two coordinates (x,y) suffice to determine normalized active stress. Furthermore, in this subspace a Minkowski spacetime can be endowed, on which the group action is found to be a proper orthochronous Lorentz group SO_o(2,1). The existence of a fixed point attractor in the normalized active stress space is demonstrated by the long-term behavior deduced from the analytical solutions, which together with the response stability is further verified by Lyapunov's direct method. Two numerical schemes based on a nonlinear Volterra integral equation and on a group symmetry are derived, the latter of which exactly preserves the consistency condition for every time step. The consistent scheme is stable, accurate and efficient, because it updates the stress point automatically on the yield surface at each time step without any iteration. For the purpose of comparison and contrast, numerical results calculated by the above two schemes as well as by the radial return method were displayed for several loading examples.     

Numerical study of the atmospheric flow over a coastal cliff

This work presents a numerical method suitable for the study of the development of internal boundary layers (IBL) and their characteristics for flows over various types of coastal cliffs. The IBL is an important meteorological occurrence for flows with surface roughness and topographical step changes. A two‐dimensional flow program was used for this study. The governing equations were written using the vorticity–velocity formulation. The spatial derivatives were discretized by high‐order compact finite differences schemes. The time integration was performed with a low storage fourth‐order Runge–Kutta scheme. The coastal cliff (step) was specified through an immersed boundary method. The validation of the code was done by comparison of the results with experimental and observational data. The numerical simulations were carried out for different coastal cliff heights and inclinations. The results show that the predominant factors for the height of the IBL and its characteristics are the upstream velocity, and the height and form (inclination) of the coastal cliff. Copyright © 2010 John Wiley & Sons, Ltd.     

Exponential time integration using Krylov subspaces

The application of exponential integrators based on Krylov techniques to large‐scale simulations of complex fluid flows with multiple time‐scales demonstrates the efficiency of these schemes in reducing the associated time‐step restrictions due to numerical stiffness. Savings of approximately 50% can be achieved for simulations of the three‐dimensional compressible Navier–Stokes equations while still maintaining a truncation error typical of explicit time‐stepping schemes. Exponential time integration techniques of this type are particularly advantageous for fluid flows with a wide range of temporal scales such as low‐Mach number, reactive or acoustically dominated flows. Copyright © 2008 John Wiley & Sons, Ltd.     

An automatic constraint violation stabilization method for differential/ algebraic equations of motion in multibody system dynamics

ＩｎｔｒｏｄｕｃｔｉｏｎＴｈｅｄｙｎａｍｉｃｅｑｕａｔｉｏｎｓｏｆｍｏｔｉｏｎｏｆｍｕｌｔｉｂｏｄｙｓｙｓｔｅｍｓｗｉｔｈｃｏｎｓｔｒａｉｎｔｓａｒｅｔｈｅｆｏｌｌｏｗｉｎｇｄｉｆｆｅｒｅｎｔｉａｌ／ａｌｇｅｂｒａｉｃｅｑｕａｔｉｏｎｓ,ｉ．ｅ．,Ｅ...