Runge-Kutta Time Stepping Research Articles

Three-Dimensional Dynamic Model of TEHD Tilting-Pad Journal Bearing—Part II: Parametric Studies

This paper presents a new analysis method for a thermo-elasto-hydro-dynamic (TEHD) tilting pad journal bearing (TPJB) system to reach a static equilibrium condition adopting nonlinear transient dynamic solver, whereas earlier studies have used iteration schemes such as Newton–Raphson method. The theoretical TPJB model discussed in Part I of this research is combined into a newly developed algorithm to perform a bearing dynamic analysis and present dynamic coefficients. In the nonlinear transient dynamic solver, physical and modal coordinates coexist for computational efficiency, and transformation between modal and physical coordinate is performed at each numerical integration time step. Variable time step Runge–Kutta numerical integration scheme is adopted for a reliable and fast calculation. Nonlinear time transient dynamic analysis and steady thermal analysis are combined to find the static equilibrium condition of the TPJB system, where the singular matrix issue of flexible pad finite element (FE) model is resolved. The flexible pad TPJB model was verified by comparison with other numerical results. Simulation results corresponding with the theoretical model explained in Part I are presented and discussed. It explains how the TPJB dynamic behavior is influenced by a number of eigenvector of flexible pad FE model and pad thickness. Preload change under fluid and thermal load is examined.

Numerical simulation on hyperbolic diffusion equations using modified cubic B-spline differential quadrature methods

In this article, a modified cubic B-spline differential quadrature method (MCB-DQM) is proposed to solve a hyperbolic diffusion problem in which flow motion is affected by both convection and diffusion. One dimensional hyperbolic non-homogeneous heat, wave and telegraph equations are also considered along with two dimensional hyperbolic diffusion problem. The method reduces the hyperbolic problem into a system of nonlinear ordinary differential equations. The system is then solved by the optimal four stage three order strong stability-preserving time stepping Runge–Kutta (SSP-RK43) scheme. The reliability and efficiency of the method has been tested on seven examples. The stability of the method is also discussed and found to be unconditionally stable.

Stabilisation of discrete steady adjoint solvers

A new implicit time-stepping scheme which uses Runge–Kutta time-stepping and Krylov methods as a smoother inside FAS-cycle multigrid acceleration is proposed to stabilise the flow solver and its discrete adjoint counterpart. The algorithm can fully converge the discrete adjoint solver in a wide range of cases where conventional point-implicit methods fail due to either physical or numerical instability. This enables the discrete adjoint to be applied to a much wider range of flow regimes. In addition, the new algorithm offers improved efficiency when applied to stable cases for which the conventional Block–Jacobi solver can fully converge. Both stable and unstable cases are presented to demonstrate the improved robustness and performance of the new scheme. Eigen-analysis is presented to outline the mechanism of the adjoint stabilisation effect.

A step-size selection strategy for explicit Runge–Kutta methods based on Lyapunov exponent theory

Using a stability result for variable time-step Runge–Kutta methods applied to nonautonomous real linear scalar test problem that decays exponentially fast, a step-size selection algorithm is devised. The step-size selection algorithm is based partly on stability information obtained by estimating the discrete Lyapunov exponent of a Runge–Kutta method applied to a nonautonomous linear scalar problem. The utility of the approach is illustrated in numerical experiments that demonstrate how the new algorithm performs against a standard accuracy based step-size selection strategy.

Gas cavity–body interactions: Efficient numerical solution

The paper investigates the interactions occurring between a gas cavity, the surrounding liquid and the nearby structures. In more detail, focus is in the characterization of the various dynamical phases (e.g. “acoustic phase” and “gas bubble phase”) and in the design of a modelling approach aimed at minimizing the computational efforts needed to analyse the cases of gas cavities spatially close to the target structure. Hence, a domain decomposition (DD) strategy is proposed which enables efficient computations. Hyperbolic flow equations govern the flow evolution and, while the inner domain 1D solution is calculated by means of an HLL scheme for the fluxes and a 1st order time stepping, the outer domain 3D solution is achieved on the basis of a MUSCL scheme coupled with a 3rd order Runge–Kutta time stepping. Various comparative tests, based on use of the full-scale experimental data by Smith (1975), have been used to test the DD strategy. A simplified approach is, finally, proposed to be used for the complex case of multiple cavity explosions. Use of the approach reveals that the worst load scenario for the target structure occurs when all cavities explode simultaneously.

Experimental validation of a numerical model for predicting the trajectory of blood drops in typical crime scene conditions, including droplet deformation and breakup, with a study of the effect of indoor air currents and wind on typical spatter drop trajectories

Bloodstain Pattern Analysis (BPA) provides information about events during an assault, e.g. location of participants, weapon type and number of blows. To extract the maximum information from spatter stains, the size, velocity and direction of the drop that produces each stain, and forces acting during flight, must be known.A numerical scheme for accurate modeling of blood drop flight, in typical crime scene conditions, including droplet oscillation, deformation and in-flight disintegration, was developed and validated against analytical and experimental data including passive blood drop oscillations, deformation at terminal velocity, cast-off and impact drop deformation and breakup features. 4th order Runge–Kutta timestepping was used with the Taylor Analogy Breakup (TAB) model and Pilch and Erdman's (1987) expression for breakup time. Experimental data for terminal velocities, oscillations, and deformation was obtained via digital high-speed imaging. A single model was found to describe drop behavior accurately in passive, cast off and impact scenarios.Terminal velocities of typical passive drops falling up to 8m, distances and times required to reach them were predicted within 5%. Initial oscillations of passive blood drops with diameters of 1mm<d<6mm falling up to 1.5m were studied. Predictions of oscillating passive drop aspect ratio were within 1.6% of experiment. Under typical crime scene conditions, the velocity of the drop within the first 1.5m of fall is affected little by drag, oscillation or deformation.Blood drops with diameter 0.4–4mm and velocity 1–15m/s cast-off from a rotating disk showed low deformation levels (Weber number<3). Drops formed by blunt impact 0.1–2mm in diameter at velocities of 14–25m/s were highly deformed (aspect ratios down to 0.4) and the larger impact blood drops (∼1–1.5mm in diameter) broke up at critical Weber numbers of 12–14. Most break-ups occurred within 10–20cm of the impact point. The model predicted deformation levels of cast-off and impact blood drops within 5% of experiment. Under typical crime scene conditions, few cast-off drops will break up in flight. However some impact-generated drops were seen to break up, some by the vibration, others by bag breakup.The validated model can be used to gain deep understanding of the processes leading to spatter stains, and can be used to answer questions about proposed scenarios, e.g. how far blood drops may travel, or how stain patterns are affected by winds and draughts.

Polynomial differential quadrature method for numerical solutions of the generalized Fitzhugh–Nagumo equation with time-dependent coefficients

In this paper, polynomial differential quadrature method (PDQM) is applied to find the numerical solution of the generalized Fitzhugh–Nagumo equation with time-dependent coefficients in one dimensional space. The PDQM reduces the problem into a system of first order non-linear differential equations. Then, the obtained system is solved by optimal four-stage, order three strong stability-preserving time-stepping Runge–Kutta (SSP-RK43) scheme. The accuracy and efficiency of the proposed method are demonstrated by three test examples. The numerical results are shown in max absolute errors (L∞), root mean square errors (RMS) and relative errors (L2) forms. Numerical solutions obtained by this method when compared with the exact solutions reveal that the obtained solutions are very similar to the exact ones.

Numerical solution of Burgers’ equation with modified cubic B-spline differential quadrature method

In this paper, a new numerical method, “modified cubic-B-spline differential quadrature method (MCB-DQM)” is proposed to find the approximate solution of the Burgers’ equation. The modified cubic-B-spline basis functions are used in differential quadrature to determine the weighting coefficients. The MCB-DQM is used in space, and the optimal four-stage, order three strong stability-preserving time-stepping Runge–Kutta (SSP-RK43) scheme is used in time for solving the resulting system of ordinary differential equations. To check the efficiency and accuracy of the method, four examples of Burgers’ equation are included with their numerical solutions, L2 and L∞ errors and comparisons are done with the results given in the literature. The proposed method produces better results as compared to the results obtained by almost all the schemes available in the literature, and approaching to the exact solutions. The presented method is seen to be easy, powerful, efficient and economical to implement as compared to the existing techniques for finding the numerical solutions for various kinds of linear/nonlinear physical models.

A parallel local timestepping Runge–Kutta discontinuous Galerkin method with applications to coastal ocean modeling

Geophysical flows over complex domains often encompass both coarse and highly resolved regions. Approximating these flows using shock-capturing methods with explicit timestepping gives rise to a Courant–Friedrichs–Lewy (CFL) timestep constraint. This approach can result in small global timesteps often dictated by flows in small regions, vastly increasing computational effort over the whole domain. One approach for coping with this problem is to use locally varying timesteps. In previous work, we formulated a local timestepping (LTS) method within a Runge–Kutta discontinuous Galerkin framework and demonstrated the accuracy and efficiency of this method on serial machines for relatively small-scale shallow water applications. For more realistic models involving large domains and highly complex physics, the LTS method must be parallelized for multi-core parallel computers. Furthermore, additional physics such as strong wind forcing can effect the choice of local timesteps. In this paper, we describe a parallel LTS method, parallelized using domain decomposition and MPI. We demonstrate the method on tidal flows and hurricane storm surge applications in the coastal regions of the Western North Atlantic Ocean.

An unstructured finite volume numerical scheme for extended 2D Boussinesq-type equations

We present a high-order well-balanced unstructured finite volume (FV) scheme on triangular meshes for modeling weakly nonlinear and weakly dispersive water waves over slowly varying bathymetries, as described by the 2D depth-integrated extended Boussinesq equations of Nwogu, rewritten here in conservation law form. The FV scheme numerically solves the conservative form of the equations following the median dual node-centered approach, for both the advective and dispersive part of the equations. For the advective fluxes, the scheme utilizes an approximate Riemann solver along with a well-balanced topography source term upwinding. Higher order accuracy in space and time is achieved through a MUSCL-type reconstruction technique and through a strong stability preserving explicit Runge–Kutta time stepping. Special attention is given to the accurate numerical treatment of moving wet/dry fronts and boundary conditions. The model is applied to several examples of non-breaking wave propagation over variable topographies and the computed solutions are compared to experimental data. The presented results indicate that the presented FV model is robust and capable of simulating wave transformations from relatively deep to shallow water, providing accurate predictions of the wave's propagation, shoaling and runup.

A high-order adaptive time-stepping TVD solver for Boussinesq modeling of breaking waves and coastal inundation

We present a high-order adaptive time-stepping TVD solver for the fully nonlinear Boussinesq model of Chen (2006), extended to include moving reference level as in Kennedy et al. (2001). The equations are reorganized in order to facilitate high-order Runge–Kutta time-stepping and a TVD type scheme with a Riemann solver. Wave breaking is modeled by locally switching to the nonlinear shallow water equations when the Froude number exceeds a certain threshold. The moving shoreline boundary condition is implemented using the wetting–drying algorithm with the adjusted wave speed of the Riemann solver. The code is parallelized using the Message Passing Interface (MPI) with non-blocking communication. Model validations show good performance in modeling wave shoaling, breaking, wave runup and wave-averaged nearshore circulation.

A path-conservative method for a five-equation model of two-phase flow with an HLLC-type Riemann solver

Compressible multi-phase flows are found in a variety of scientific and engineering problems. The development of accurate and efficient numerical algorithms for multi-phase flow simulations remains one of the challenging issues in computational fluid dynamics. A main difficulty of numerical methods for multi-phase flows is that the model equations cannot always be written in conservative form, though they may be hyperbolic and derived from physical conservation principles. In this work, assuming a hyperbolic model, a path-conservative method is developed to deal with the non-conservative character of the equations. The method is applied to solve the five-equation model of Saurel and Abgrall for two-phase flow. As another contribution of the work, a simplified HLLC-type approximate Riemann solver is proposed to compute the Godunov state to be incorporated into the Godunov-type path-conservative method. A second order, semi-discrete version of the method is then constructed via a MUSCL reconstruction with Runge–Kutta time stepping. Moreover, the method is then extended to the two-dimensional case by directional splitting. The method is systematically assessed via a series of test problems with exact solutions, finding satisfactory results.

Nonuniform time-step Runge–Kutta discontinuous Galerkin method for Computational Aeroacoustics

With many superior features, Runge–Kutta discontinuous Galerkin method (RKDG), which adopts Discontinuous Galerkin method (DG) for space discretization and Runge–Kutta method (RK) for time integration, has been an attractive alternative to the finite difference based high-order Computational Aeroacoustics (CAA) approaches. However, when it comes to complex physical problems, especially the ones involving irregular geometries, the time step size of an explicit RK scheme is limited by the smallest grid size in the computational domain, demanding a high computational cost for obtaining time accurate numerical solutions in CAA. For computational efficiency, high-order RK method with nonuniform time step sizes on nonuniform meshes is developed in this paper. In order to ensure correct communication of solutions on the interfaces of grids with different time step sizes, the values at intermediate-stages of the Runge–Kutta time integration on the elements neighboring such interfaces are coupled with minimal dissipation and dispersion errors. Based upon the general form of an explicit p-stage RK scheme, a linear coupling procedure is proposed, with details on the coefficient matrices and execution steps at common time-levels and intermediate time-levels. Applications of the coupling procedures to Runge–Kutta schemes frequently used in simulation of fluid flow and acoustics are given, including the third-order TVD scheme, and low-storage low dissipation and low dispersion (LDDRK) schemes. In addition, an analysis on the stability of coupling procedures on a nonuniform grid is carried out. For validation, numerical experiments on one-dimensional and two-dimensional problems are presented to illustrate the stability and accuracy of proposed nonuniform time-step RKDG scheme, as well as the computational benefits it brings. Application to a one-dimensional nonlinear problem is also investigated.

High-order finite-volume methods for the shallow-water equations on the sphere

This paper presents a third-order and fourth-order finite-volume method for solving the shallow-water equations on a non-orthogonal equiangular cubed-sphere grid. Such a grid is built upon an inflated cube placed inside a sphere and provides an almost uniform grid point distribution. The numerical schemes are based on a high-order variant of the Monotone Upstream-centered Schemes for Conservation Laws (MUSCL) pioneered by van Leer. In each cell the reconstructed left and right states are either obtained via a dimension-split piecewise-parabolic method or a piecewise-cubic reconstruction. The reconstructed states then serve as input to an approximate Riemann solver that determines the numerical fluxes at two Gaussian quadrature points along the cell boundary. The use of multiple quadrature points renders the resulting flux high-order. Three types of approximate Riemann solvers are compared, including the widely used solver of Rusanov, the solver of Roe and the new AUSM +-up solver of Liou that has been designed for low-Mach number flows. Spatial discretizations are paired with either a third-order or fourth-order total-variation-diminishing Runge–Kutta timestepping scheme to match the order of the spatial discretization. The numerical schemes are evaluated with several standard shallow-water test cases that emphasize accuracy and conservation properties. These tests show that the AUSM +-up flux provides the best overall accuracy, followed closely by the Roe solver. The Rusanov flux, with its simplicity, provides significantly larger errors by comparison. A brief discussion on extending the method to arbitrary order-of-accuracy is included.

DNS of gas bubbles behaviour using an improved 3D front tracking model—Model development

In recent years CFD has proven to be a valuable and powerful tool to advance our understanding of complex multiphase flow systems arising in industrial applications. However, the predictive capabilities of this tool are determined by many factors of physical and numerical origin but in particular by the quality of the closures adopted for the description of the interface forces. The objective of this study is to improve the front tracking method in order to compute such forces with sufficient accuracy. This paper describes the further development of a 3D front tracking model to achieve improved volume conservation and circumvent problems related to the representation of surface tension. First, we have included a method to handle the pressure jump at the interface. This causes the spurious currents, observed in conventional front tracking, to decrease with two orders of magnitude. Also the advection scheme has been adapted, using higher order velocity interpolation (using cubic splines), and Runge–Kutta time-stepping, in order to prevent considerable volume changes of the dispersed phase. Test simulations involving a stationary bubble, a standard advection test and an oscillating droplet, demonstrate the effect of these improvements. The implementation of these procedures enlarged the computational window and in particular enabled the simulation of very small bubbles, where large surface forces dominate, without any significant spurious currents or volume loss.

An unstructured node-centered finite volume scheme for shallow water flows with wet/dry fronts over complex topography

An unstructured high-resolution, node-centered finite volume algorithm for triangular grids is developed in order to simulate unsteady, two-dimensional, shallow water flows over arbitrary topography with wetting and drying. The algorithm utilizes Roe’s approximate Riemann solver, to compute the numerical fluxes, while second-order spatial accuracy is achieved with a MUSCL reconstruction technique, using a slope limiter to control the total variation of the reconstructed field, and an explicit four stage Runge–Kutta time stepping. The novel aspects of the algorithm include the extension to second-order of the topography source term treatment and the wet/dry front treatment, using the MUSCL slope limiting method, within the node-centered finite volume formulation. The developed wet/dry treatment is found to ensure absolute mass conservation. At the beginning of the procedure a triangular grid is locally refined (performing h-enrichment) at regions with steep variations of the topography in order to obtain a better approximation of irregular terrain characteristics and improve the accuracy of the results in the flow for those regions in which a more complex flow is associated to more abrupt geometry without an excess in computational cost. The proposed numerical scheme is extensively validated, with respect to its effectiveness and robustness, against benchmark test cases and experimental data related to propagation and run-up of long waves over arbitrary topographies.

Dissipativity of Runge–Kutta methods for neutral delay differential equations with piecewise constant delay

The numerical approximation of neutral dissipative initial value problems with piecewise constant delay by fixed time stepping Runge–Kutta methods is considered and two sufficient conditions for the dissipativity of the system are given. The concept of (weak) E ( λ ) -dissipativity is introduced. It is shown that under some assumptions a DJ -irreducible, algebraically stable Runge–Kutta method is (weakly) E ( λ ) -dissipative.

Phase Space Stability Error Control with Variable Time‐stepping Runge‐Kutta Methods for Dynamical Systems

AbstractWe consider a phase space stability error control for numerical simulation of dynamical systems. We illustrate how variable time‐stepping algorithms perform poorly for long time computations which pass close to a fixed point. A new error control was introduced in [9], which is a generalization of the error control first proposed in [8]. In this error control, the local truncation error at each step is bounded by a fraction of the solution arc length over the corresponding time interval. We show how this error control can be thought of either a phase space or a stability error control. For linear systems with a stable hyperbolic fixed point, this error control gives a numerical solution which is forced to converge to the fixed point. In particular, we analyze the forward Euler method applied to the linear system whose coefficient matrix has real and negative eigenvalues. We also consider the dynamics in the neighborhood of saddle points. We introduce a step‐size selection scheme which allows this error control to be incorporated within the standard adaptive algorithm as an extra constraint at negligible extra computational cost. Theoretical and numerical results are presented to illustrate the behavior of this error control. (© 2004 WILEY‐VCH Verlag GmbH &amp; Co. KGaA, Weinheim)

Simulation of reactive sputtering kinetics in real in-line processing chambers

We present a dynamic, coupled macroscopic model (DOGMA) for the process of reactive sputtering suited for description of real in-line processing. In our model we divide the complex volume of a processing chamber into simple cells, e.g. parallelepipeds. For each cell gettering kinetics of reactive gas at metallized surfaces as well as sputtering kinetics of the target are calculated using a Runge–Kutta time step method. The inert and reactive gas flow is obtained via flow conductances defined for each connection between adjacent cells. For deposition of sputtered particles on chamber and substrate surfaces, a pressure dependent distribution matrix is used, which is obtained from single-particle Monte Carlo calculations. The flow conductances for different pressure conditions are calculated using the well known ‘Direct Simulation Monte Carlo’ method implemented on a Linux cluster, where for each volume cell an individual calculation task can be spawned. Tuned with flow conductances and particle distribution factors from Monte Carlo calculations our coupled model is capable to describe real in-line sputtering processes on a single-processor PC at almost real time speed. All numerical models are implemented within a simulation system RIG-VM developed at Fraunhofer IST, which enables easy communication between separate modules, e.g. during modeling of a controlling cycle and ensures interchangeability and reusability by using XML datasheets and web browser based user interfaces. In this paper we apply our model to reactive sputtering of SiO 2 with λ-probe based oxygen partial pressure stabilization. The model is adapted to experiment via material parameter variation. The resulting simulation system can be applied for optimization of process parameters with respect to process stability and film homogeneity.

Efficient Cell-Centered Multigrid Scheme for the Three-Dimensional Navier-Stokes Equations

A cell-centered scheme for three-dimensional Navier-Stokes equations, which is based on central-difference approximations and Runge-Kutta time stepping, is described. By using local time stepping, implicit residual smoothing, a multigrid method, and carefully controlled artificial dissipative terms, good convergence rates are obtained for two- and three-dimensional flows. The emphases are on the implicit smoothing and artificial dissipative terms with locally variable coefficients which depend on cell aspect ratios. The computational results for two-dimensional subsonic airfoil flows and three-dimensional transonic C-D nozzle flows are essentially coincident with other experimental and calculated results.