Temporal Discretization Errors Research Articles

Greedy trial subspace selection in meshfree time-stepping scheme with applications in coupled bulk-surface pattern formations

Combining kernel-based collocation methods with time-stepping methods to solve parabolic partial differential equations can potentially introduce challenges in balancing temporal and spatial discretization errors. Typically, using kernels with high orders of smoothness on some sufficiently dense set of trial centers provides high spatial approximation accuracy that can exceed the accuracy of finite difference methods in time. The paper proposes a greedy approach for selecting trial subspaces in the kernel-based collocation method applied to time-stepping to balance errors in both well-conditioned and ill-conditioned scenarios. The approach involves selecting trial centers using a fast block-greedy algorithm with new stopping criteria that aim to balance temporal and spatial errors. Numerical simulations of coupled bulk-surface pattern formations, a system involving two functions in the domain and two on the boundary, illustrate the effectiveness of the proposed method in reducing trial space dimensions while maintaining accuracy.

Verification and validation of a coupled CFD–FEM approach with overset structured grid through free decay and regular wave tests

This paper presents a comprehensive verification and validation study applied to the hydrodynamic analysis of a semi-submersible Floating Offshore Wind Turbine (FOWT) platform. The numerical simulation is conducted using a coupled Computational Fluid Dynamics-Finite Element Method (CFD–FEM) approach. The CFD program incorporates an overset grid interpolation method, facilitating a more systematic grid refinement in the estimation of discretization uncertainties. A mooring analysis program, based on the FEM scheme, is integrated with the CFD program to simulate the dynamic responses of mooring system. The necessity of the nonlinear mooring model is demonstrated by comparing solutions with a comparative test using a linear mooring model. The verification study, based on a least-squares-based extrapolation method, quantifies the spatial and temporal discretization errors of the numerical model. The numerical solutions based on the baseline model are further validated against the model test measurements. A pitch free decay model test and a regular wave model test conducted within the OC5 project are utilized for the verification and validation study. Results from both tests demonstrate the reliability and accuracy of the coupled method.

K-space dispersion error compensators for the fractional spatial derivatives based constant-Q viscoelastic wave equation modeling

Unlike the fractional temporal derivatives based constant-Q viscoelastic wave equation, the fractional spatial derivatives based one is preferred for the wavefield simulation due to its tractable memory consumption. This viscoelastic wave equation is often solved by the staggered-grid pseudospectral method, which has the spatial spectral accuracy but temporal low order accuracy. The k-space method is one of the potential remedies for the temporal discretization errors. However, it is difficult to fully remove the temporal dispersion errors by directly applying the elastic k-space error compensators to the viscoelastic wave equation, which contains additional temporal derivatives associated with the amplitude decay terms of the P- and S- waves. To tackle the issue, we have derived a new wavenumber-time domain viscoelastic k-space equation by the eigenvalue decomposition method and by the solution to the matrix differential equations. The error compensation terms inside our k-space equation can simultaneously correct the discretization errors of different types of temporal derivatives and are wave mode adaptive. The accuracy of the k-space wave equation is verified by the viscoelastic Green's function's solution. To further boost the simulation efficiency of the k-space method in heterogeneous media, we have proposed an improved k-space method, which decreases the number of the costly mixed-domain operators. The simulation results from the simple and complex models demonstrate that our two k-space methods are superior to the conventional staggered-grid pseudospectral method for the viscoelastic wave modeling.

Adaptive Timestepping, Linearization, and A Posteriori Error Control for Multiphase Flow of Immiscible Fluids in Porous Media with Wells

Summary This work focuses on the development of adaptive timesteps, stopping criteria, and error control strategies for reservoir simulations with fully implicit (FIM) solvers. Using a rigorous error control framework, an adaptive time selector combined with nonlinear stopping criteria is used to control nonlinear iterations as well as to balance accuracy and robustness for challenging nonlinear simulations. In reservoir simulation, efficiently solving a system of nonlinear equations arising from the FIM method can be computationally burdensome for complex recovery processes. Theoretically, an FIM reservoir simulator has no stability limit on the timestep size. In practice, standard Newton’s method often fails to converge for large timestep sizes and must therefore cut the timestep multiple times to achieve convergence, resulting in a large number of unnecessary iterations. Another cause of nonlinear convergence issues is the presence of wells, which are often presented as singular point/line sources strongly coupled to the reservoir model, posing additional restrictions on the timestep choice. Here, we use a posteriori error estimators to avoid unnecessary nonlinear iterations and timestep cuts when solving immiscible multiphase flow. First, we estimate error components (e.g., spatial, temporal, and nonlinear) and then apply these to balancing criteria, providing us with dynamic and adaptive strategies to control timestep and nonlinear iterations. The error estimators are fully and locally computable, inexpensive to use, and target the various error components, including well singularities. The method provides an adaptive criterion for stopping the nonlinear iteration process whenever the linearization error does not significantly affect the overall error. Simultaneously, timesteps are adapted to maintain a constant size of the temporal discretization error with respect to the total error. Altogether, this avoids using unnecessary linearization iterations, wasteful timestep cuts, and too small timesteps. To demonstrate the effectiveness of these adaptive features, we present results for a suite of cases, covering both standard benchmarks and conceptual problems incorporating highly heterogeneous media with multiple wells. The proposed timestep selector cooperates with the new stopping criteria to improve nonlinear solver performance and increases robustness for cases with high nonlinearity. Perhaps most important, the adaptive features ensure balanced temporal and spatial errors while maintaining sufficiently small nonlinear errors, which ensures solution accuracy by accurately reproducing saturation fronts, production plateau, and breakthrough times.

Optimal convergence analysis of a linearized second-order BDF-PPIFE method for semi-linear parabolic interface problems

The article proposes and analyzes the optimal error estimates of a second-order backward difference formula (BDF2) numerical scheme for the semi-linear parabolic interface problems. The partially penalized immersed finite element (PPIFE) methods are used for the spatial discretization to resolve discontinuity of the diffusion coefficient across the interface. The classical extrapolation method is adopted to treat the nonlinear term, which effectively avoids the complicated numerical calculation of the nonlinearity. Our error analysis is based on the corresponding time-discrete system, which neatly splits the error into two parts: the temporal discretization error and the spatial discretization error. Since the spatial discretization error is independent of time step size τ, we can unconditionally derive the optimal error estimates in both L2 norm and semi-H1 norm, while previous works always require the coupling condition of time step and space size. Numerical experiments are given to confirm the theoretical analysis.

NTFA-enabled goal-oriented adaptive space–time finite elements for micro-heterogeneous elastoplasticity problems

In this work, we establish a goal-oriented space–time finite element method for a class of dissipative heterogeneous materials. Those materials are modeled on both micro- and macroscale, with a scale transition of volume averaging type satisfying the Hill–Mandel condition. A nonuniform transformation field analysis is performed on the microscopic inelastic strain fields for a model reduction. Reduced variables are deduced from a space–time decomposition of those inelastic strain fields. Closed-form constitutive relations are derived from some dissipative considerations, thus resulting into a reduced order homogenization problem. The resulting model error is sufficiently small for the considered class of materials, thus leaving the discretization error of the finite element method as a main error source. For ease of error estimate, we rewrite the reduced order problem in a multifield formulation. Based on duality techniques, a backward-in-time dual problem is derived from a Lagrange method, rendering error representations of a user-defined quantity of interest. Combining a patch recovery technique, a computable error estimator is developed to quantify both spatial and temporal discretization errors. By means of a localization technique, local error estimators are used to drive a greedy adaptive refinement algorithm in space and time. The effectiveness of the resulting algorithm is confirmed by several numerical examples w.r.t. a prototype model.

On the properties of high-order least-squares finite-volume schemes

High-order finite-volume schemes based on polynomial least-squares methods are an active research topic for the discretization of hyperbolic equations as they allow to obtain high-order spatial discretization schemes in arbitrary meshes. However, few studies have analyzed their performance in good-quality/near-to-uniform meshes, which are commonly used as a meshing strategy in zones where turbulent effects are important. In this paper, the theoretical numerical properties of commonly used least-squares (LSQ) k-exact high-order finite volume schemes are studied in one-dimensional and in several two-dimensional meshes (with some remarks regarding their properties in three-dimensional meshes). These results are compared to those obtained using fully-constrained polynomial reconstructions only compatible with structured meshes. The numerical properties of the schemes are investigated through the von Neumann analysis methodology applied to the one-dimensional and two-dimensional finite volume formulation, including temporal discretization errors. This analysis is also extended to non-uniform and unstructured two-dimensional meshes. At last, the schemes are tested with several numerical experiments using the linear advection, the Euler equations and the Navier-Stokes equations. The analysis of both studies yields similar conclusions regarding the numerical errors and stability of the different studied schemes showing that the high-order least-squares finite volume schemes yield stable and robust results across different uniform and non-uniform unstructured meshes. However, their performance is heavily degraded in simulations with low mesh resolution compared to schemes specially catered to structured meshes. On the other hand, the latter schemes lack stability and robustness in general structured meshes and its formulation cannot be straightforwardly extended to unstructured meshes. Moreover, this work shows that the use of weighted-LSQ can drastically improve the results of LSQ schemes in under-resolved simulations.

A novel and efficient sequential-explicit technique for poroelasticity problems

This paper presents a novel sequential-explicit coupling technique to solve hydromechanical problems in poroelastic media. The proposed approach combines appropriate coupling terms defined in a rate form and an automatic time-stepping scheme that controls the temporal discretization error through a single parameter. Such technique has been implemented in a framework that manages sequential simulations, exchanges of information between models and adapts the time marching in a multirate approach. We present several numerical tests in which we evaluate the new scheme’s stability, accuracy, and performance through comparisons with the fully implicit solution. We observe that the proposed technique is stable, robust and accurate, and can be used in engineering applications for large scale and long time durations.

Two Dimensional Analysis of Hybrid Spectral/Finite Difference Schemes for Linearized Compressible Navier\u2013Stokes Equations

This study seeks to compare different combinations of spatial dicretization methods under a coupled spatial temporal framework in two dimensional wavenumber space. The aim is to understand the effect of dispersion and dissipation on both the convection and diffusion terms found in the two dimensional linearized compressible Navier–Stokes Equations (LCNSE) when a hybrid finite difference/Fourier spectral scheme is used in the x and y directions. In two dimensional wavespace, the spectral resolution becomes a function of both the wavenumber and the wave propagation angle, the orientation of the wave front with respect to the grid. At sufficiently low CFL number where temporal discretization effects can be neglected, we show that a hybrid finite difference/Fourier spectral schemes is more accurate than a full finite difference method for the two dimensional advection equation, but that this is not so in the case of the LCNSE. Group velocities, phase velocities as well as numerical amplification factor were used to quantify the numerical anisotropy of the dispersion and dissipation properties. Unlike the advection equation, the dispersion relation representing the acoustic modes of the LCNSE contains an acoustic terms in addition to its advection and viscous terms. This makes the group velocity in each spatial direction a function of the wavenumber in both spatial directions. This can lead to conditions for which a hybrid Fourier spectral/finite difference method can become less or more accurate than a full finite difference method. To better understand the comparison of the dispersion properties between a hybrid and full FD scheme, the integrated sum of the error between the numerical group velocity V^{*}_{grp,full} and the exact solution across all wavenumbers for a range of wave propagation angle is examined. In the comparison between a hybrid and full FD discretization schemes, the fourth order central (CDS4), fourth order dispersion relation preserving (DRP4) and sixth order central compact (CCOM6) schemes share the same characteristics. At low wave propagation angle, the integrated errors of the full FD and hybrid discretization schemes remain the same. At intermediate wave propagation angle, the integrated error of the full FD schemes become smaller than that of the hybrid scheme. At large wave propagation angle, the integrated error of the full FD schemes diverges while the integrated error of the hybrid discretization schemes converge to zero. At high reduced wavenumber and sufficiently low CFL number where temporal discretization error can be neglected, it was found that the numerical dissipation of the viscous term based on the CDS4, DRP4, CCOM6 and isotropy optimized CDS4 schemes (hbox {CDS4}_{{opt}}) schemes was lower than the actual physical dissipation, which is only a function of the cell Reynolds number. The wave propagation angle at which the numerical dissipation of the viscous term approaches its maximum occurs at pi /4 for the CDS4, DRP4, CCOM6 and hbox {CDS4}_{{opt}} schemes.

A Time-Domain Preconditioner for the Helmholtz Equation

Time-harmonic solutions to the wave equation can be computed in the frequency or in the time domain. In the frequency domain, one solves a discretized Helmholtz equation, while in the time domain, the periodic solutions to a discretized wave equation are sought, e.g. by simulating for a long time with a time-harmonic forcing term. Disadvantages of the time-domain method are that the solutions are affected by temporal discretization errors and that the spatial discretization cannot be freely chosen, since it is inherited from the time-domain scheme. In this work we address these issues. Given an indefinite linear system satisfying certain properties, a matrix recurrence relation is constructed, such that in the limit the exact discrete solution is obtained. By iterating a large, finite number of times, an approximate solution is obtained, similarly as in a time-domain method for the Helmholtz equation. To improve the convergence, the process is used as a preconditioner for GMRES, and the time-harmonic forcing term is multiplied by a smooth window function. The construction is applied to a compact-stencil finite-difference discretization of the Helmholtz equation, for which previously no time-domain solver was available. Advantages of the resulting solver are the relative simplicity, small memory requirement and reasonable computation times.

Flexible goal-oriented adaptivity for higher-order space–time discretizations of transport problems with coupled flow

In this work, a flexible higher-order space–time adaptive finite element approximation of convection-dominated transport with coupled fluid flow is developed and studied. Convection-dominated transport is a challenging subproblem in poromechanics in which coupled transport with flow, chemical reaction and mechanical response in porous media is considered. Key ingredients are the arbitrary degree discontinuous Galerkin time discretization of the primal and dual problems for the Dual Weighted Residual (DWR) approach, an a posteriori error estimation for the transport problem coupled with flow and its implementation in an advanced software architecture. The error estimate allows the separation of the temporal and spatial discretization error contributions which facilitates the simultaneous adjustment of the time and space mesh. The performance of the approach and its software implementation is studied by numerical convergence examples as well as an example of physical interest for convection-dominated cases.

Analysis of transient wave propagation in inhomogeneous media using edge-based gradient smoothing technique and bathe time integration method

As is known to all, the classical finite element method (FEM) always fails to effectively solve the wave problems in the relatively large wave numbers due to the troublesome numerical error issue. In this paper, the standard FEM with edge-based gradient smoothing technique and Bathe time integration method is developed to analyze the transient wave propagation problems in inhomogeneous media. We explicitly show that the numerical error of the calculated solutions for transient wave propagation problems consists of two different parts, namely the spatial discretization error and the temporal discretization error. Due to the edge-based gradient smoothing technique and the appropriate numerical dissipation effects from Bathe time integration scheme, it is found that the total numerical error can be significantly suppressed and more accurate numerical solutions can be obtained. Several typical numerical examples have been conducted to examine the capacity of the proposed method in solving transient wave propagation problems in inhomogeneous media.

Validation Studies of Linear Oscillating Compressor Cascade and Use of Influence Coefficient Method

Abstract Advanced predictions of blade flutter have been continually pursued. It is noted however that validation cases of unsteady CFD methods against experimental cases with detailed 3D unsteady pressures are still rather lacking. The main objectives of the present work are two-folds. First, validate and understand the characteristics of blade tip clearance, as well as a bubble-type flow separation for an unsteady CFD solver against a 3D oscillating cascade experiment. And second, examine the applicability of the influence coefficient method (ICM) as widely used in an oscillating linear cascade setup. In the first part, the capability of a widely used commercial solver (CFX) for unsteady flows induced by a 3D oscillating compressor cascade is examined. The present computations have shown consistently a destabilizing effect of increasing blade tip clearance, in agreement with the experiment. More remarkably, the computational analyses reveal a distinctive interplay between the inlet endwall boundary layer and the tip clearance in relation to the aerodynamic damping. Different inlet endwall boundary layer thicknesses are shown to lead to qualitatively different aeroelastic stability characteristics in relation to tip clearance. The aero-damping variation with the tip clearance under the influence of the inlet endwall boundary layer seems to correlate closely to a balancing act between the passage vortex and the tip leakage vortex. The tip clearance aeroelastic behavior seems also in line with a simple quasi-steady analysis. On the other hand, the mid-chord laminar bubble separation on suction surface, though with a clear signature in the local aero-damping, has negligible effects on the overall stability. The second part aims to examine computationally the applicability of the influence coefficient method in a linear cascade setup. The comparison between the cascade-based ICM data and a baseline “tuned cascade” shows that the differences in the sensitivity to the far-field treatment can be significant, depending on inter-blade phase angles. On the other hand, non-linearity effects closely relevant to the basic linear assumption of the ICM are shown to only have a small influence. The present results suggest that extra caution should be exercised when comparing a CFD-based tuned cascade model with a finite cascade-based ICM model, at conditions close to acoustic resonance. The resultant discrepancies may well arise from the inherently different far-field sensitivities between the two models, rather than those typical numerical and physical modeling aspects of interest (e.g., meshing, spatial and temporal discretization errors as well as turbulence modeling).

Sequential local mesh refinement solver with separate temporal and spatial adaptivity for non-linear two-phase flow problems

Convergence failure and slow convergence rates are among the biggest challenges with solving the system of non-linear equations numerically. Although mitigated, such issues still linger when using strictly small time steps and unconditionally stable fully implicit schemes. The price that comes with restricting time steps to small scales is the enormous computational load, especially in large-scale models. To address this problem, we introduce a sequential local mesh refinement framework of temporal and spatial adaptivity to optimize convergence rate and prevent convergence failure, while not restricting the whole system to small time steps, thus improving computational efficiency. Two types of error estimators are introduced to estimate the spatial discretization error, the temporal discretization error separately. These estimators provide a global upper bounds on the dual norm of the residual and the non-conformity of the numerical solution for non-linear two phase flow models. The mesh refinement algorithm starts from solving the problem on the coarsest space-time mesh, then the mesh is refined sequentially based on the spatial error estimator and the temporal error estimator. After each refinement, the solution from the previous mesh is used to estimate the initial guess of unknowns on the current mesh for faster convergence. Numerical results are presented to confirm accuracy of our algorithm as compared to the uniformly fine time step and fine spatial discretization solution. We observe around 25 times speedup in the solution time by using our algorithm.

Stable and contact‐free time stepping for dense rigid particle suspensions

SummaryWe consider suspensions of rigid bodies in a two‐dimensional viscous fluid. Even with high‐fidelity numerical methods, unphysical contact between particles occurs because of spatial and temporal discretization errors. We extend a time stepping method that avoids overlap by imposing a minimum separation distance between all pairs of bodies. In its original form, the method discretizes interactions between different particles explicitly. Therefore, to avoid stiffness, a large minimum separation distance is used. In this paper, we introduce a new implicit time stepping method that is able to simulate dense suspensions with large time step sizes and a small minimum separation distance. The method is tested on various unbounded and bounded flows, and rheological properties of the resulting suspensions are computed.

Convergence and error analysis of an automatically differentiated finite volume based heat conduction code

Purpose This paper aims to investigate the convergence and error properties of a finite volume-based heat conduction code that uses automatic differentiation to evaluate derivatives of solutions outputs with respect to arbitrary solution input(s). A problem involving conduction in a plane wall with convection at its surfaces is used as a test problem, as it has an analytical solution, and the error can be evaluated directly. Design/methodology/approach The finite volume method is used to discretize the transient heat diffusion equation with constant thermophysical properties. The discretized problem is then linearized, which results in two linear systems; one for the primary solution field and one for the secondary field, representing the derivative of the primary field with respect to the selected input(s). Derivatives required in the formation of the secondary linear system are obtained by automatic differentiation using an operator overloading and templating approach in C++. Findings The temporal and spatial discretization error for the derivative solution follows the same order of accuracy as the primary solution. Second-order accuracy of the spatial and temporal discretization schemes is confirmed for both primary and secondary problems using both orthogonal and non-orthogonal grids. However, it has been found that for non-orthogonal cases, there is a limit to the error reduction, which is concluded to be a result of errors in the Gauss-based gradient reconstruction method. Originality/value The convergence and error properties of derivative solutions obtained by forward mode automatic differentiation of finite volume-based codes have not been previously investigated.

Goal-oriented h-type adaptive finite elements for micromorphic elastoplasticity

Generalized continuum theories become necessary for many scenarios like strain localization phenomena or size effects. In this work, we consider a class of higher order continua, namely micromorphic continua, where the kinematics is enhanced by means of a microstructure undergoing an affine micro deformation. The increasing number of the degrees of freedom in such a theory clearly motivates an application of adaptive methods. For linear elastic micromorphic problems, we have shown the consistency of the resulted dual problem, ensuring an optimal convergence order, in Ju and Mahnken (2017). In this work, we extend the framework of goal-oriented adaptivity to time-dependent problems, i.e. plasticity problems, where a backwards-in-time dual problem is crucial to account for the error accumulation over time caused by both error generation and error transport. Our focus is limited to an adaptive control of spatial discretization errors of the finite element method (FEM), while the temporal discretization errors are not considered for simplicity. Based on duality techniques, we derive exact error representations. For practice, four computable error estimators are proposed, where two different ways to obtain enhanced solutions are considered, and where additionally an approximate forwards-in-time dual problem neglecting error transport is introduced. By means of certain localization techniques, these estimators are used to drive an adaptive mesh refinement algorithm. Their effectiveness is shown by several numerical examples based on a prototype model.

A fully analytic space-time depletion benchmark for fixed-source and eigenvalue calculations and its application to Monte Carlo simulation

Formal analysis of the convergence and stability of coupled neutron transport/depletion calculations is difficult due to the complicated nature of the nonlinear system of equations that must be solved. In this paper, we present a new benchmark problem for coupled neutron transport/depletion calculations along with its fully analytic solution. The benchmark problem – which we call the “Candlestick benchmark” because it models a progressive burnup from one end of a slab to another – is straightforward to implement and can be used to quantify the accuracy and stability of coupled transport/depletion calculations. Descriptions and solutions are provided for both a fixed-source and eigenvalue version of the benchmark. The Candlestick benchmark is a first-of-its-kind analytic solution to a nonlinear spatial depletion problem, and it is useful for the purpose of validating the solution methods used by reactor depletion solvers. As a demonstration of the usefulness of the Candlestick benchmark, we have used it to assess some of the convergence behaviors of Monte Carlo depletion calculations, including the effects of spatial discretization, temporal discretization, and stochastic error on the results of the calculation.

POD reduced-order modeling for evolution equations utilizing arbitrary finite element discretizations

The main focus of the present work is the inclusion of spatial adaptivity for the snapshot computation in the offline phase of model order reduction utilizing proper orthogonal decomposition (POD-MOR) for nonlinear parabolic evolution problems. We consider snapshots which live in different finite element spaces, which means in a fully discrete setting that the snapshots are vectors of different length. From a numerical point of view, this leads to the problem that the usual POD procedure which utilizes a singular value decomposition of the snapshot matrix, cannot be carried out. In order to overcome this problem, we here construct the POD model/basis using the eigensystem of the correlation matrix (snapshot Gramian), which is motivated from a continuous perspective and is set up explicitly, e.g., without the necessity of interpolating snapshots into a common finite element space. It is an advantage of this approach that the assembly of the matrix only requires the evaluation of inner products of snapshots in a common Hilbert space. This allows a great flexibility concerning the spatial discretization of the snapshots. The analysis for the error between the resulting POD solution and the true solution reveals that the accuracy of the reduced-order solution can be estimated by the spatial and temporal discretization error as well as the POD error. Finally, to illustrate the feasibility of our approach, we present a test case of the Cahn–Hilliard system utilizing h-adapted hierarchical meshes and two settings of a linear heat equation using nested and non-nested grids.

A Time-Spectral Algorithm for Fractional Wave Problems

This paper develops a high-accuracy algorithm for time fractional wave problems, which employs a spectral method in the temporal discretization and a finite element method in the spatial discretization. Moreover, stability and convergence of this algorithm are derived, and numerical experiments are performed, demonstrating the exponential decay in the temporal discretization error provided the solution is sufficiently smooth.