Search

While the literature on numerical methods (e.g. finite elements and, to a certain extent, virtual elements) concentrates on convex elements, there is a need to probe beyond this limiting constraint so that the field can be further explored and developed. Thus, in this paper, we employ the virtual element method for non-convex discretizations of elastodynamic problems in 2D and 3D using an explicit time integration scheme. In the formulation, a diagonal matrix-based stabilization scheme is proposed to improve performance and accuracy. To address efficiency, a critical time step is approximated and verified using the maximum eigenvalue of the local (rather than global) system. The computational results demonstrate that the virtual element method is able to consistently handle general nonconvex elements and even non-simply connected elements, which can lead to convenient modeling of arbitrarily-shaped inclusions in composites.

We present a stable and accurate moving particle semi-implicit (MPS) for free-surface flows. For simulating free-surface flows, previous formulations using renormalized procedure is only applicable to interior particles. Poor accuracy of the free-surface particle in turn contaminate the overall solution accuracy. In this study, dummy particle is introduced to stabilize the renormalized Laplacian model on free surface and to enforce the free-surface condition. The proposed renormalized Laplacian model enables accurate pressure calculation for both internal and free-surface particles. Modification also includes a modified gradient model in which linear consistency is guaranteed and stability is enhanced by introducing dummy particle. The proposed MPS scheme provides a more accurate and stable simulation of incompressible free surface flows. Several two-dimensional numerical simulations are presented to demonstrate the enhanced performance of the proposed MPS method.

Starting from a highly continuous isogeometric analysis discretization, we introduce hyperplanes that partition the domain into subdomains and reduce the continuity of the discretization spaces at these hyperplanes. As the continuity is reduced, the number of degrees of freedom in the system grows. The resulting discretization spaces are finer than standard maximal continuity IGA spaces. Despite the increase in the number of degrees of freedom, these finer spaces deliver simulation results faster with direct solvers than both traditional finite element and isogeometric analysis for meshes with a fixed number of elements. In this work, we analyze the impact of continuity reduction on the number of Floating Point Operations (FLOPs) and computational times required to solve fluid flow and electromagnetic problems with structured meshes and uniform polynomial orders. Theoretical estimates show that for sufficiently large grids, an optimal continuity reduction decreases the computational cost by a factor of O ( p 2 ) . Numerical results confirm these theoretical estimates. In a 2D mesh with one million elements and polynomial order equal to five, the discretization including an optimal continuity pattern allows to solve the vector electric field, the scalar magnetic field, and the fluid flow problems an order of magnitude faster than when using a highly continuous IGA discretization. 3D numerical results exhibit more moderate savings due to the limited mesh sizes considered in this work.

We present a method of CutFEM type for the Poisson problem with either Dirichlet or Neumann boundary conditions. The computational mesh is obtained from a background (typically uniform Cartesian) mesh by retaining only the elements intersecting the domain where the problem is posed. The resulting mesh does not thus fit the boundary of the problem domain. Several finite element methods (XFEM, CutFEM) adapted to such meshes have been recently proposed. The originality of the present article consists in avoiding integration over the elements cut by the boundary of the problem domain, while preserving the optimal convergence rates, as confirmed by both the theoretical estimates and the numerical results.

In this paper, we propose a diffuse interface model and finite element approximation for two-phase magnetohydrodynamic (MHD) flows with different viscosities and electric conductivities. An energy stable scheme, which is based on the finite element method for the spatial discretization and first order semi-implicit scheme combined with convex splitting method for the temporal discretization, is proposed to solve this new model. The numerical scheme is proved to be mass-conservative and energy law preserving. By Leray–Schauder fixed point theorem, the existence of solutions to the numerical scheme is shown. The uniqueness of the numerical solutions is obtained. Utilizing the stability of the numerical scheme and the compactness method, the existence of the weak solutions to the two-phase MHD model is established as well. Furthermore, given more regularity on the weak solution, the convergence of the numerical scheme is derived. Finally, numerical experiments are provided to verify the theoretical results and validate the proposed model.

Suspension of an oxytactic bacteria (e.g. the species Bacillus subtilis) placed in a container with its upper surface open to the atmosphere results in the formation of complex bioconvection patterns. The bacteria consume the oxygen diluted in the water, thereby causing the decrease of oxygen concentration everywhere except on the free surface. Through the free surface, which is in direct contact with the air, oxygen diffuses into the water. Slightly denser than water, the oxytactic bacteria are able to swim towards the higher concentration of oxygen (i.e. upwards) and they concentrate in a thin layer below the free surface. This causes the change of the suspension density and Rayleigh–Taylor type instabilities to occur. The chemotaxis phenomenon has been successfully modeled within continuum mechanics approach under certain simplifications. The set of (non-linearly) coupled equations describing the process involves the Boussinesq approximation of the Navier–Stokes equations governing the fluid motion and two convection–diffusion type equations governing the bacteria and oxygen concentrations. One of the simplifications that might significantly influence numerical simulations is the boundary condition for fluid equation on the free surface. This condition ensures that the vertical component of the velocity is zero, thus keeping the position of free surface fixed. This assumption significantly simplifies numerical procedure since the non-linearly coupled system can then be solved on stationary grid. However, allowing the motion of the free surface and completing the system with appropriate boundary conditions on contact line (liquid–solid–gas interface), a more realistic model is derived and new insights on nonlinear dynamics of the chemotaxis phenomenon are obtained. Our aims in this paper are to upgrade the currently available model into a more realistic one in both two and three dimensions, to propose a numerical procedure to deal with the new system (now posed on time-dependent domain) and, finally, to show the difference between this new model and the previous simplified one.

Fractures play an important role in the geothermal energy recovery from enhanced geothermal systems (EGS) involving interacted multi-physical fields. A strongly coupled thermo-hydro-mechanical (THM) model is proposed to simulate the process of the long-term geothermal production in three-dimensional geothermal reservoirs containing arbitrary discrete fracture networks (DFNs). By introducing a strong discontinuity concept into the DFN model, the aperture variation of each fracture induced by the fluid pressure, external stresses and thermal expansion in the period of production can be captured. Non-isothermal fluid flow in DFN and the local thermal non-equilibrium (LTNE) between fluid and rock matrix are formulated and coupled with the fracture deformation model. Verification is carried out against an analytical solution, followed by a sensitivity and convergence analysis concerning time step and mesh size. This approach is then applied to Habanero EGS project in Australia to evaluate the geothermal productivity and efficiency for a period of 20 years with different injection and production pressures. The results demonstrate that the proposed model is robust and effective to simulate the THM coupled process in 3D fractured reservoirs.

Two triangular G 1 -conforming elements, based on the triangular Gregory patch, suitable for the analysis of the Kirchhoff plate model are presented. Both have cubic normal derivative along the sides, so that they can be effectively used in combination with generalized Hermitian elements. The Gregory patch consists in a rational enhancement of the base-polynomial spaces useful to design G 1 -conforming elements on general C 0 -conforming unstructured meshes. Because of the presence of the rational functions, the second derivatives at the corners of the element present finite discontinuities, that prevent the elements from passing the bending patch test. The discontinuities are removed using a constrained version of the Gregory patch, with Lagrange multipliers. In this way, the rational conforming space collapses into a conforming rearrangement of the original polynomial interpolant spaces. The proposed formulation design elements that pass the bending patch test and present optimal rate of convergence on general unstructured meshes.

This work proposes a novel metric based algorithm for quadrilateral mesh generating. Each quad-mesh induces a Riemannian metric satisfying special conditions: the metric is a flat metric with cone singularities conformal to the original metric, the total curvature satisfies the Gauss–Bonnet condition, the holonomy group is a subgroup of the rotation group { e i k π ∕ 2 } , there is cross field obtained by parallel translation which is aligned with the boundaries, and its streamlines are finite geodesics. Inversely, such kind of metric induces a quad-mesh. Based on discrete Ricci flow and conformal structure deformation, one can obtain a metric satisfying all the conditions and obtain the desired quad-mesh. This method is rigorous, simple and automatic. Our experimental results demonstrate the efficiency and efficacy of the algorithm.

This paper presents a new high-order mixed scheme of polygonal finite element method (PFEM) for incompressible flow computation. The system of polygonal basis functions with bubble-enriched functions is applied to represent the fluid velocity field, while the pressure field is approximated by the set of polygonal basis functions. The new element is named as MINIPe. From a practical point of view, the present element naturally satisfies the well-known inf–sup condition. The numerical tests of the incompressible fluid flow problems governed by Stokes equation systems are presented to verify the performance quality of proposed method with respect to the stability and convergence rate.

This paper proposes a new multiscale topology optimization method for the concurrent design of multiphase composite structures under a certain range of excitation frequencies. Distinguished from the existed studies, a general concurrent design formulation for the dynamic composite structures with more than two material phases is developed. The macrostructureand its microstructures with multiple material phases are optimized simultaneously. The integral of the dynamic compliances over an interval of frequencies is formulated as the optimization objective, so as to minimize the frequency response within the concerned excitation range. The effective properties of the multiphase microstructures are evaluated by using the numerical homogenization method, which actually serves as a link to bridge the macro and micro finite element analyses. Furthermore, to describe the boundaries of multiple material phases for the microstructure, a parametric color level set method (PCLSM) is developed by using an efficient interpolation scheme. In this way, L level set functions can represent at most 2 L material phases without any overlaps. Moreover, these “color” level sets are updated by directly using the well-established gradient-based algorithm, which can greatly facilitate the proposed method to solve the multi-material optimizations with multiple design constraints. Several 2D and 3D numerical examples are used to demonstrate the effectiveness of the proposed method in the concurrent design of the dynamic composite structures under the excitation frequency ranges.

The computational modelling of the heart motion within a cardiac cycle is an extremely challenging problem due to (a) the complex multi-scale interaction that takes place between the electrophysiology and electrochemistry at cellular level and the macro-scale response of the heart muscle, and (b) the large deformations and the strongly anisotropic and quasi-incompressible behaviour of the myocardium. These pose an extreme challenge to the scalability of electro-mechanical solvers due to the size and conditioning of the system of equations required to obtain accurate solutions, both in terms of wall deformation and transmembrane potential propagation. In the search towards an efficient modelling of electro-activation, this paper presents a coupled electromechanical computational framework whereby, first, we explore the use of an efficient stabilised low order tetrahedral Finite Element methodology and compare it against a very accurate super enhanced mixed formulation previously introduced by the authors in Garcia-Blanco et al. (2019) and, second, we exploit the use of tailor-made staggered and staggered linearised solvers in order to assess their feasibility against a fully monolithic approach. Through a comprehensive set of examples, culminating in a realistic ventricular geometry, we aim to put forward some suggestions regarding the level of discretisation and coupling required to ensure sufficiently reliable results yet with an affordable computational time.

Markov Chain Monte Carlo (MCMC) simulation has a considerable computational burden when the target probability density function (PDF) evaluation involves a black-box, potentially computationally-expensive, numerical model. A novel framework to accelerate MCMC is developed here for such applications. It leverages a Kriging surrogate model approximation to the target PDF to improve computational efficiency, while preserves convergence properties to the exact target PDF, avoiding potential accuracy problems introduced through the surrogate model error. The approach relies on the delayed-rejection (DR) scheme and combines the rapid exploration characteristics of global (independent) proposals with the local search robustness of random walk proposals under the MCMC setting. The global proposal is chosen as the Kriging-based target PDF approximation. This proposal may resemble, depending on the surrogate model error, the actual target PDF very well, and therefore can provide significant computational benefits, including fast convergence of the Markov chain to its stationary distribution and, more importantly, low correlation between samples. At each MCMC step a candidate sample is generated from the independent proposal. If rejected, DR allows an extra random walk around the current sample. The DR step guarantees convergence to the actual target PDF, circumventing robustness problems that may arise when the Kriging-based independent proposal offers a poor approximation to (i.e., underestimates) the actual target PDF. The overall computational efficiency is further improved through an adaptive updating of the Kriging surrogate model during the MCMC sampling phase by extracting information from candidate samples whose inclusion in the training database can substantially enhance Kriging’s accuracy. The computational efficiency and robustness of the established algorithm, termed Adaptive Kriging Delayed Rejection Adaptive Metropolis algorithm(AK-DRAM), are demonstrated in two analytical benchmark problems and two engineering problems focusing on conditional failure sample simulation and Bayesian inference.

This paper presents an approach to the automatic enrichment of finite elements in the vicinity of a stress singularity. The enrichment consists of semi-analytical singular modes constructed using the Scaled Boundary Finite Element Method (SBFEM). In contrast to analytical methods, the SBFEM provides modes for inhomogeneous and anisotropic materials without additional effort. The finite element basis can be of arbitrary order and remains unaltered by the enrichment. The approach requires enrichment in only one layer of elements around a node. Due to the compatibility of SBFEM with FEM, there is no need for transitional elements, and there are no parasitic terms. The approach is tested for several benchmark problems. The stress intensity factors are computed based on techniques inspired by the SBFEM. The proposed procedure is compared to a standard finite element implementation and shows a significant improvement in the error of the displacement field for problems involving singular stresses.

Traditional sensitivity analysis methods for the model with correlated inputs and univariate output fail to provide satisfactory results for multivariate output. In this work, we first establish a reasonable contribution classification for the univariate output with the correlated input. Then the covariance decomposition method is extended to the case of correlated inputs as a reference, and the vector projection sensitivity index is extended to aggregate the correlated and uncorrelated contributions of the input to multiple outputs. The definition of the new sensitivity index is based on the vector projection, which can take into account both uncertainties and correlations among multiple outputs by projecting the conditional variance vector (built by the full marginal variance contributions) on the unconditional variance vector (built by unconditional variance magnitudes and correlation of the multiple outputs). The mathematical properties of the extended vector projection sensitivity index are discussed and its relations with other existing sensitivity indices are highlighted. Two numerical examples and two engineering examples about an aircraft structure are employed to illustrate the validity and potential benefits of the extended vector projection sensitivity index.

The Immersed Boundary Method has been used to simulate a range of fundamental flows and turbulence. Such studies have demonstrated the method’s promising applicability for engineering analysis. However, to the authors’ knowledge, flows in coupled components or in scenarios with coupled physics, such as rotor–stator interaction, fan–intake interaction, aeroelastics, aeroacoustics, etc., are still rarely investigated using high-fidelity methods. Due to its high computational costs, the complexity of geometry meshing process and the requirement for moving boundaries limit the investigation of flows in such environments. Previous research suggests that high-fidelity simulations with an acceptable geometry modelling strategy may tackle these issues and provide useful insights. There exists a hierarchy of geometry modelling methods which includes the conventional Directly Mesh Resolving (DMR) method, the Immersed Boundary Method (IBM), and the IBM with Smeared Geometry (IBMsg, or eIBMg). The present research proposes an alternative to these approaches in the form of the Euler IBM with local force (eIBMl) by imposing a distribution function generated from blade configuration. Compared to the eIBMg, this method can include more realistic flow physics within each blade passage without smearing its geometry. This method is applied to the study of fan–intake interaction focusing on the transport of inlet distortion through blade passages and pressure wave propagation.

As an efficient uncertainty quantification (UQ) methodology for moment propagation and probability analysis of quantities of interest, polynomial chaos (PC) expansions have received broad and sustained attentions. However, the exponentially increasing cost of building PC representations with increasing dimension of uncertainty, i.e., the curse of dimensionality, seriously restricts the practical application of PC at the industrial level. Some efficient strategies applying adaptive basis selection algorithm for sparse optimization (or l 1 -minimization) of PC show great potential compared to the classical full PC. However, these strategies mainly focus on forward selection algorithms, which are incapable of correcting any error made by these algorithms. Hence, this paper develops a novel adaptive forward–backward selection (AFBS) algorithm for reconstructing sparse PC. The proposed algorithm by a reasonable combination of forward selection and adaptively backward elimination technique can efficiently correct mistakes made by earlier forward selection steps, which retains the most significant PC terms and discards redundant or insignificant ones. The accuracy of built PC metamodel is checked by a cross-validation procedure. As a consequence, the most significant PC terms are detected sequentially and corresponding PC coefficients are accurately recovered. It largely enhances the sparsity of PC and improves the prediction accuracy compared to the popular forward selection algorithms, e.g., least angle regressions (LARs). To validate the efficiency of the proposed algorithm, a complex analytical function with Gaussian distribution inputs and two challenging aerodynamic applications including a sonic boom propagation analysis considering atmospheric uncertainty and a natural-laminar-flow (NLF) airfoil computation under both geometrical and operational uncertainties are elaborated. With an in-depth comparison with some popular PC reconstruction methodologies, the performance of the devised AFBS method is assessed comprehensively. The results demonstrate that the proposed AFBS method can efficiently identify the significant PC contributions describing the problems, and reproduce sparser PC metamodel and more accurate UQ results than those classical full PC and LARs methods.

Compressive-shear fracture is commonly observed in rock-like materials. However, this fracture type cannot be captured by current phase field models (PFMs), which have been proven an effective tool for modeling fracture initiation , propagation, coalescence , and branching in solids. The existing PFMs also cannot describe the influence of cohesion and internal friction angle on load–displacement curve during compression tests. Therefore, to develop a new phase field model that can simulate well compressive-shear fractures in rock-like materials, we construct a new driving force in the evolution equation of phase field. Strain spectral decomposition is applied and only the compressive part of the strain is used in the new driving force with consideration of the influence of cohesion and internal friction angle. For ease of implementation, a hybrid formulation is established for the phase field modeling. Then, we test the brittle compressive-shear fractures in uniaxial compression tests on intact rock-like specimens as well as those with a single or two parallel inclined flaws. All numerical results are in good agreement with the experimental observation, validating the feasibility and practicability of the proposed PFM for simulating brittle compressive-shear fractures.

The uncertainty analysis of structures generally involves uncertain parameters of different types. In order to derive predictions regarding uncertain structural responses, it is crucial to represent the uncertainty appropriately according to the underlying information available. This paper presents a unified framework for uncertainty analysis under probabilistic, evidence, fuzzy and interval uncertainties, by which the quantities with sufficient data, sparse data, and subjective information can be simultaneously considered. A Taylor expansion-based unified uncertainty analysis (T-UUA) method is first proposed for small uncertainty problems. By temporarily neglecting the evidence, fuzzy and interval variables, the probability-evidence-interval-fuzzy model is degraded into a random problem, in which the expectations and variances of responses can be obtained as functions in terms of evidence, interval and fuzzy uncertainties. Then, through dealing with the evidence variables, the previous expectations and variances are further expressed as a summation of functions in terms of fuzzy and interval variables with basic probability assignments (BPAs). The fuzziness is then discretized by using α -cut technique and thus the expectations and variances are further expressed as functions of only intervals. Afterwards, by reconsidering the interval uncertainties, the bounds of the expectations and variances are computed via combining Taylor expansion with interval arithmetic. In addition, a dimensional reduction (DR)/efficient global optimization (EGO)-based unified uncertainty analysis (DR/EGO-UUA) method is also presented to solve the large uncertainty problems. The framework of DR/EGO-UUA is similar as T-UUA. However, in DR/EGO-UUA, the second moments of responses are computed by DR integrations, and their upper and lower bounds are calculated by the EGO. Finally, three numerical examples are investigated to demonstrate the effectiveness of the proposed methods.

The aim of this paper is to accurately solve short wave scattering problems governed by the Helmholtz equation using the Bernstein–Bezier Finite Element method (BBFEM), combined with a conformal perfectly matched layer (PML). Enhanced PMLs, where curved geometries are represented by means of the blending map method of Gordon and Hall, are numerically investigated. In particular, the performance of radial and elliptical shaped PMLs, with a parabolic absorption function, are assessed and compared in terms of accuracy against second order Bayliss–Gunzberger–Turkel (BGT 2 ) based local absorbing boundary conditions. Numerical results dealing with problems of Hankel source radiation and wave scattering by a rigid cylinder show that the radial shaped PML, with the h and p versions of BBFEM, enables the recovery of the predicted algebraic and exponential convergence rates of a high order finite element method (FEM). Furthermore, radial shaped BGT 2 and PML have a similar performance, as long as the wave is not sufficiently well resolved. But, BGT 2 performs poorly as the wave resolution increases. Additionally, the effect of harmonics of higher modes on accuracy is examined. The study reveals that the PML outperforms BGT 2 for almost all propagating modes. However, a similar performance is achieved with both methods either with higher modes or a low wave resolution. Results from a multiple scattering benchmark problem provide evidence of the good performance of the proposed PMLs and the benefit of elliptical shaped PMLs in reducing significantly the size of the computational domain, without altering accuracy. The choice of the PML parameters ensuring optimal performance is also discussed.

Topology optimization of shell–infill structures has been a hot topic in the optimization community. The crucial issue of such design problem is how to accurately describe the material interfaces, which is often regarded as a relatively difficult problem in density-based topology optimization. This paper presents an erosion-based interface identification method to handle this difficulty by defining the different parts of the original and eroded structures as interfaces. The theoretical relation for determining the erosion parameters is derived to control the interface thickness accurately. It turns out that the erosion-based method offers excellent performance even for the extreme cases where the ratios of member sizes to a prescribed thickness are small, and at the same time, it is easy to implement since we only need to consider the frequently-used filtering and projection processes in the applied erosion operation. Then, we provide an improved SIMP-based topology optimization method for shell–infill structures based on the new idea of defining interfaces. In this paper, a tweak to the existing two-step filtering/projection process is made to separate the shell and infill, and a corresponding interpolation function is developed to model the whole composite objects. Enhanced by the worst-case based robust formulation, the minimum length scale of optimized structures can be controlled. Several 2D and 3D compliance optimization examples are provided to illustrate the effectiveness of the proposed method.

When a multi-material structure is optimized by using the level set method, besides the topological changes that are familiar to people (i.e., vanishing of a hole, nucleation of a hole, merging of holes, and breaking apart of a region), another topological change, i.e., the conversion between a pair of materials, is also possible to happen but receives little attention. In this case, a hole is no longer a region of void but becomes occupied by a different solid material. In our present study, such a topological change is called the generalized hole nucleation, and the material conversion scheme of the BESO is integrated into the level set based framework to deal with it. A strategy determining where and when to conduct generalized hole nucleation is presented. In addition, the shape derivatives of the structure’s traction free boundary and the interfaces between material regions are carefully analyzed. The results of several numerical examples prove that the proposed method is effective and efficient.

The combination of tetrahedral and hexahedral elements in a single conformal mesh requires pyramids or prisms to make the transition between triangular and quadrilateral faces. This paper presents high order exact sequences of finite element approximations in H 1 ( Ω ) , H ( curl , Ω ) , H ( div , Ω ) , and L 2 ( Ω ) based on such kind of three dimensional mesh configurations. The approach is to consider composite polynomial approximations based on local partitions of the pyramids into two or four tetrahedra. The traces associated with triangular faces of these tetrahedral elements are constrained to match the quadrilateral shape functions on the quadrilateral face of the pyramid, in order to maintain conformity with shared neighboring hexahedron, or prism. Two classes of composite exact sequences are constructed, one using classic Nedelec spaces of first kind, and a second one formed by enriching these spaces with properly chosen higher order functions with vanishing traces. Projection-based interpolants satisfying the commuting diagram property are presented in a general form for each type of element. The interpolants are expressed as the sum of linearly independent contributions associated with vertices, edges, faces, and volume, according to the kind of traces appropriate to the space under consideration. Furthermore, we study applications to the mixed formulation of Darcy’s problems based on compatible pairs of approximations in { H ( div , Ω ) , L 2 ( Ω ) } for such tetrahedral–hexahedral–prismatic–pyramidal meshes. An error analysis is outlined, showing same (optimal) orders of approximation in terms of the mesh size as one would obtain using purely hexahedral or purely tetrahedral partitions. Enhanced accuracy for potential and flux divergence variables is obtained when enriched space configurations are applied. The predicted convergence orders are verified for some test problems.

Coupling the finite volume method (FVM) and the moving-particle semi-implicit (MPS) method, a conservative hybrid method is proposed for simulation of incompressible interfacial flow. In this hybrid method, bulk flow is resolved using the multi-moment constrained interpolation profile (CIP) FVM while the interface region is rendered using moving particles calculated by several enhanced MPS schemes. A new boundary condition is developed for particles located on the boundary of the interface. Numerical particles in the MPS method are used to determine the density of each mesh cell and the conservation of mass is guaranteed by a volume-fraction remedy procedure. The approach allows accurate interface tracking without diminishing numerical efficiency. Several benchmark simulations are conducted to demonstrate the accuracy and efficiency of the proposed hybrid method.

This paper proposes a stress-constrained topology optimization method for structures integrating truss-like microstructures and a solid material. By utilizing the representative volume element (RVE) method, the truss-like microstructures are equivalent to homogeneous materials. Then, a power-law material interpolation scheme for bi-material, including the solid material and the equivalent material, is constructed. Based on the mechanical analysis, a yield index g is introduced to indicate the yielding situation of the truss-like microstructures. In addition, the yielding of truss-like microstructures is controlled implicitly by constraining the von Mises stress of the equivalent materials. By means of two p-norm clustered stress measures, the stress of the solid material and the equivalent material are constrained. Finally, based on the adjoint method, the sensitivity of stress and displacement constraints is obtained and the optimization can be solved by a gradient-based optimization algorithm. Two numerical examples are presented to demonstrate the validity and effectiveness of the proposed methodology.

In this paper, we propose accurate and efficient finite difference methods to discretize the two- and three-dimensional fractional Laplacian $(-\Delta)^{\frac{\alpha}{2}}$ ($0 < \alpha < 2$) in hypersingular integral form. The proposed finite difference methods provide a fractional analogue of the central difference schemes to the fractional Laplacian, and as $\alpha \to 2^-$, they collapse to the central difference schemes of the classical Laplace operator $-\Delta$. We prove that our methods are consistent if $u \in C^{\lfloor\alpha\rfloor, \alpha-\lfloor\alpha\rfloor+\epsilon}({\mathbb R}^d)$, and the local truncation error is ${\mathcal O}(h^\epsilon)$, with $\epsilon > 0$ a small constant and $\lfloor \cdot \rfloor$ denoting the floor function. If $u \in C^{2+\lfloor\alpha\rfloor, \alpha-\lfloor\alpha\rfloor+\epsilon}({\mathbb R}^d)$, they can achieve the second order of accuracy for any $\alpha \in (0, 2)$. These results hold for any dimension $d \ge 1$ and thus improve the existing error estimates for the finite difference method of the one-dimensional fractional Laplacian. Extensive numerical experiments are provided and confirm our analytical results. We then apply our method to solve the fractional Poisson problems and the fractional Allen-Cahn equations. Numerical simulations suggest that to achieve the second order of accuracy, the solution of the fractional Poisson problem should {\it at most} satisfy $u \in C^{1,1}({\mathbb R}^d)$. One merit of our methods is that they yield a multilevel Toeplitz stiffness matrix, an appealing property for the development of fast algorithms via the fast Fourier transform (FFT). Our studies of the two- and three-dimensional fractional Allen-Cahn equations demonstrate the efficiency of our methods in solving the high-dimensional fractional problems., Comment: 24 pages, 6 figures, and 6 tables

We develop and implement a new nonlocal crack-bridging model for evaluating toughening and strengthening mechanisms of fiber reinforced cementitious composites . The quasibrittle fracture of cement matrix is characterized by XFEM-cohesive crack framework. To explore fiber-bridging mechanism, the horizontal fiber phases corresponding to different embedding length are superimposed on the near-crack region of the matrix phase , while an effective homogeneous Cauchy medium is assumed outside the multiphase region. The sliding interface for each fiber phase is described by a slip field, leading to a new nonlocal slip model being proposed which can be used to treat arbitrary fiber distributions. This model is formulated in a variationally consistent framework. In this study, it is verified that the toughening effect is induced by the interfacial shear stress rather than the fiber tensile stress. Prior to evaluate new problems, the capability of the proposed model is validated to reproduce consistent fracture behaviors with the experimentations. It is concluded from our systematic studies that a higher fiber modulus leads to a higher loading capacity, but not a higher composite toughness, while unaffects its residual strength. It is found that the enhancing interfacial strength has improved the composites toughness dramatically, but has limited contribution to the strengthening effect .

In this work, we propose an adaptive phase field method (PFM) to simulate quasi-static brittle fracture problems. The phase field equations are solved using the scaled boundary finite element method (SBFEM). The adaptive refinement strategy is based on an error indicator evaluated directly from the solutions of the SBFEM without any need for stress recovery techniques. Quadtree meshes are adapted to perform mesh refinement. The polygons with hanging nodes in the quadtree decomposition are treated as n − sided polygons within the framework of the SBFEM and do not require any special treatment in contrast to the conventional finite element method. Several benchmark problems are used to demonstrate the robustness and the efficacy of the proposed technique. The adaptive refinement strategy reduces the mesh burden when adopting the PFM to model fracture. Numerical results show an improvement in the computational efficiency in terms of the number of elements required in the standard PFM without compromising the accuracy of the solution.

With the rapid development of reliability analysis theory, the safety assessment for engineering systems with hybrid epistemic uncertainties has received increasing attentions. To overcome the imperfection of traditional single-uncertainty modeling methods, this paper introduces evidence variables and fuzzy variables simultaneously to describe the epistemic uncertain parameters. A novel dual-stage reliability analysis framework is presented, where the first stage incorporates the evidence information by the belief and plausibility measures; the second stage incorporates the fuzzy information by a membership function-like formula. To improve the computational efficiency of response bound prediction associated with various joint focal elements and λ k -cut interval variables, a universal metamodel is constructed by radial basis functions in the global design space. The Latin Hypercube design method without overlapped data is adopted as the sampling strategy. Finally, two numerical examples verify the effectiveness of the proposed method in both mathematical theory and engineering application.

Peridynamics is a nonlocal reformulation of solid mechanics which can reproduce the solution of classical continuum mechanics as its horizon approaches zero. Taking this idea into consideration, we proposed a new coupling approach between meshfree peridynamics and finite element method. An intermediate variable horizon element called “morakkab” is introduced to transfer force and momentum between peridynamic and FEM regions. In the original peridynamic formulation, artificial ghost forces will arise if the horizon varies with position; The creation of these artifacts is avoided by utilizing a modified peridynamic equation of motion based on the concept of co-family. The accuracy of the method is illustrated in 2D benchmark tests. Using the modified equation of motion, compared to the original peridynamic formulation, dramatically decreases the relative errors with respect to the analytical solution.

A localizing gradient damage model with micro inertia effect is proposed for the dynamic fracture of quasi-brittle materials. The objective is to achieve mesh independent solutions, and to avoid spurious effects associated with the conventional nonlocal enhancement. The proposed localizing gradient damage model closely resembles the conventional gradient enhancement, albeit with an interaction domain that decreases with damage, complemented by a micro inertia effect. We first consider a classical crack branching problem, where the localizing gradient damage model is shown to resolve the mesh sensitivity issue, as well as to correctly reproduce the crack profile. Moreover, the micro inertia effect is observed to retard the crack velocity. Next, the tensile loading of a Polymethyl Methacrylate plate is considered. It is shown that the proposed model effectively captures the experimentally observed transition of crack profiles as the loading rate increases, i.e. from a straight crack propagation, to sub-branching, and finally to macro branching. Numerical results in terms of crack patterns, crack velocities, and fracture energies are in good agreement with the experimental data. To furthermore demonstrate the superior performance of the localizing gradient damage model, the macro branching problem is solved using the conventional gradient enhancement with micro inertia. It is shown that a spurious damage growth and an erroneous interaction between closely spaced cracks suppress the development of macro branching, even though reasonable values are obtained for the fracture energy and crack velocity. The localizing gradient damage model is able to fully resolve these issues.

The paper presents a new two-scale methodology for modelling of damage in ductile heterogeneous materials employing homogenization approach over the two microstructural models. The first microstructural model comprises elastoplastic behaviour without softening, while only microstructural softening is homogenized in the second model. The scale transition of the state variables is based on the first-order computational homogenization scheme. The microstructural damage response is governed by the nonlocal variable embedded into the gradient elastoplastic formulation using the implicit localization approach. The averaging computation over the microstructural softening zones is driven by the damage variable and the local equivalent strain. Therein discretization is performed by using the mixed finite element formulation including the nonlocal equivalent plastic strain interpolation, while the standard displacement based finite elements, employing the simple constitutive relation with damage variable, are used at the macrolevel. The capability of modelling material damage responses at both the micro and macroscale is demonstrated by the benchmark examples.

The boundary element method (BEM) is a powerful tool in computational acoustics, because the analysis is conducted only on structural surfaces, compared to the finite element method (FEM) which resorts to special techniques to truncate infinite domains. The isogeometric boundary element method (IGABEM) is a recent progress in the category of boundary element approaches, which is inspired by the concept of isogeometric analysis (IGA) and employs the spline functions of CAD as basis functions to discretize unknown physical fields. As a boundary representation approach, IGABEM is naturally compatible with CAD and thus can directly perform numerical analysis on CAD models, avoiding the cumbersome meshing procedure in conventional FEM/BEM and eliminating the difficulty of volume parameterization in isogeometric finite element methods. The advantage of tight integration of CAD and numerical analysis in IGABEM renders it particularly attractive in the application of structural shape optimization because (1) the geometry and the analysis can be interacted, (2) remeshing with shape morphing can be avoided, and (3) an optimized solution returns a CAD geometry directly without postprocessing steps. In the present paper, we apply the IGABEM to structural shape optimization of three dimensional exterior acoustic problems, fully exploiting the strength of IGABEM in addressing infinite domain problems and integrating CAD and numerical analysis. We employ the Burton–Miller formulation to overcome fictitious frequency problems, in which hyper-singular integrals are evaluated explicitly. The gradient-based optimizer is adopted and shape sensitivity analysis is conducted with implicit differentiation methods. The design variables are set to be the positions of control points which directly determine the shape of structures. Finally, numerical examples are provided to verify the algorithm.

This paper applies high-order discontinuous Galerkin finite element solutions of the shallow water equations to realistic coastal hydrodynamics problems. The goal is to obtain results on par with the accuracy of current low-order models, while significantly surpassing their performance. Accomplishing this requires developing mesh discretizations with boundary and parameter field representations that are consistent with the high-order accuracy of the numerical methods. These developments allow coarse-mesh, high-order, solutions to provide comparable accuracy to today’s state-of-the-art high-resolution, low-order coastal models at reduced computational expense.

Stiffness degradation due to matrix cracks is the main initial form of damage in composite laminates. This paper presents a framework to model geometrically nonlinear large deformation behaviors of matrix cracked hybrid composite double-curved deep shell containing carbon nanotube reinforced composite (CNTRC) layers. Two types of structures, namely Structure-I and Structure-II, are investigated. The CNTRC layers in Structure-I are considered with CNTs arranged in uniformly distributions, while Structure-II is arranged in functionally graded distributions. The degraded stiffness of cracked layers is modeled via the self-consistent model (SCM) micromechanical framework. To describe the geometrically nonlinear large deflection behaviors and account for deep and moderate thick shells, the von Karman geometric nonlinearity assumptions and the term 1/(1+ ς /R) are considered in the relationship between displacement and strain. The IMLS-Ritz method is employed to discretize the non-linear partial differential equations. The modified Newton–Raphson method in combination with the arc-length iteration technique is adopted to solve the discretized equations. Comparison studies indicate that the proposed predictive model can furnish very accurate results for linear and nonlinear behaviors of thin to moderately thick as well as shallow and deep laminated doubly-curved shells. Parametric studies on the effect of CNT distribution, matrix crack density, load type, length-to-thickness ratio, radius-to-length ratio, aspect ratio, boundary condition, and fiber ply-angle on the geometrically nonlinear large deformation behaviors of spherical hybrid composite shells are investigated.

A 3D fluid–structure coupling between Smoothed Particle Hydrodynamics (SPH) and Finite Element (FE) methods is proposed in this paper, with its application to complex tire hydroplaning simulations on rough ground. The purpose of this work is to analyze the SPH–FE coupling capabilities for modeling efficiently such a complex phenomenon. On the fluid side, the SPH method is able to handle the three complex interfaces of the hydroplaning phenomenon: free-surface, ground/fluid and fluid/tire interfaces. On the solid side, the FE method is used for its ability to treat tire–ground contact. A new algorithm dedicated to such SPH–FE coupling strategies is proposed to optimize the computational efficiency through the use of differed time steps between fluid and solid solvers. This way, the number of calls to the FE solver is minimized while maintaining the accuracy and stability of the coupling. The ratio between these respective time steps relies on a control procedure based on pressure loading. The present 3D SPH–FE model is first validated with different academic test cases and experimental data before considering the complex problem of the 3D hydroplaning simulations. Hydroplaning simulations are performed and analyzed on 3D configurations involving both smooth and rough grounds.

Composite shells with complex geometry are widely used in aerospace structures. Due to the complexity of geometry and curvilinear fiber path, the analysis and optimization based on finite element analysis (FEA) for complex variable-stiffness (VS) shells is extremely time-consuming. By comparing with FEA, isogeometric analysis (IGA) exhibits higher prediction efficiency of buckling load. In this work, the formula of geometric stiffness matrix for complex VS shells is derived for the first time based on degenerated shell method using IGA, which is the basis of performing linear buckling analysis. Then, a new variable curvature quasi-linear function (VCQLF) to describe curvilinear fiber path is proposed, which can further expand the design space of VS shells. After that, two frameworks for shape optimization of complex shells are put forward and then compared, and it is found that the one based on LOFT function can provide representative control variables of shape and effectively reduces the number of design variables for complex shells. Finally, a novel collaborative optimization framework of fiber path and shell shape using IGA is established. By comparison of traditional methods, it is demonstrated that the proposed framework can greatly improve the efficiency of optimization and fully explore the buckling load of complex VS shells.

We present a computational framework for the shape optimization of a Horizontal-Axis Wind Turbine (HAWT) rotor blade under uncertainty. Our framework integrates aerodynamic simulations based on the blade element method which utilizes reduced order models of the blade structure and wind load with design sensitivity analysis and nonlinear programming. The wind velocity is modeled as a stochastic process to account for variations in time and space. An additional stochastic process accounts for uncertainties in the structural material properties. The uncertainty propagation is based on a non-intrusive polynomial chaos expansion that allows accurate estimation of stochastic performance metrics such as generated power and structural compliance. Sensitivities of cost and constraint functions with respect to shape parameters namely twist angles are computed with an efficient scheme to enable gradient based optimization. To demonstrate the effect of uncertainty, designs obtained from optimization under uncertainty are compared to those obtained from deterministic optimization.

This paper introduces an algorithm for reconstructing three-dimensional (3D) heterogeneous media, in which the randomly distributed individual entities, such as aggregates, pores or other material phases with different sizes and shapes can be explicitly represented. The reconstruction procedure achieves high efficiency and reliability by integrating random process theory and stochastic reconstruction procedure. Specifically, the random process theory is applied to generate an initial model based on the statistics obtained from the digital image of heterogeneous material sample and then this coarse model is evolved by the simulated annealing algorithm to match the same statistical descriptors of sample structure. Two 3D models of different sizes are reconstructed. The autocorrelation function is applied to validate the accuracy of the reconstructed models. The comparison between the reconstructed models and the real model shows good agreement. The proposed method is suitable for reconstructing large size 3D random media for its high computational efficiency.

Lattice structures have been widely studied due to their advantage of low stiffness-to-weight ratio or sometimes auxetic properties. This paper presents a topology optimization method for structures with functionally-graded infill lattices with buckling constraints, which minimizes compliance while ensuring a prescribed level of structural stability against buckling failures. To realize topologically-optimized structures filled with functionally-graded lattices, Helmholtz PDE-filter with a variable radius is applied on the density field in Solid Isotropic Material with Penalization (SIMP) method. Buckling load factors based on the linear buckling analysis is employed as buckling constraints. Numerical examples show that proposed method can generate stiff structures comparable to the ones by the SIMP, with functionally-graded infill lattices that improve the structural stability by avoiding long, slender features under compression.

International audience; This work aims at proposing a new procedure for parametric problems whose separated representation has been considered difficult, or whose SVD compression impacted the results in terms of performance and accuracy. The proposed technique achieves a fully separated representation for layered domains with interfaces exhibiting waviness or – more generally – deviating from planar surfaces, parallel to the coordinate plane. This will make possible a simple separated representation, equivalent to others, already analyzed in some of our former works. To prove the potentialities of the proposed approach, two benchmarks will be addressed, one of them involving an efficient space–time separated representation achieved by considering the same rationale.

A block coordinate descent algorithm within a discontinuous Galerkin finite element formulation is proposed for both explicit and implicit energy minimization in solid bodies in which rigid-cohesive fractures initiate and propagate. The crack opening dependent Barenblatt type surface cohesive energy function is non-differentiable at the origin. In each iteration or each explicit update, the algorithm minimizes the potential with respect to each crack opening displacement unknown and with respect to the block of deformation unknowns, sequentially. This decomposes minimization of the full non-differentiable problem into a number of “small” non-differentiable sub-problems that must be solved locally at the inter-element boundaries of the finite element mesh and an “easy” differentiable sub-problem that characterizes global equilibrium. As a result, non-differentiability is effectively treated locally and the algorithm can be easily incorporated into standard finite element codes. On the basis of a convexity analysis of the proposed non-differentiable energy functional, we obtain a minimum cohesive process zone resolution criterion, known empirically in the previous literature as a requirement for capturing the amount of dissipated fracture energy correctly. The method is free of any regularization parameters and preserves the discrete nature of fracture without introducing a stress discontinuity at initiation of decohesion. Robustness of the method is shown through several numerical examples of fragmentation and branching and through comparisons with existing numerical and experimental results.

We derive a thermodynamically-consistent, three-dimensional, rate form-based finite-deformation constitutive theory and computational approach for damage & fracture in nonlinear viscoelastic materials. The key ingredient that allows us to develop a physical criterion for fracture is the use of a Gibbs potential-based multinetwork formulation of viscoelasticity. An approach that is based on the use of the criticality of the averaged Gibbs potential over a fracture process zone is used for the initiation & propagation of cracks in a body. The rate form-based constitutive theory and fracture criterion are implemented into the Abaqus (2018) finite-element program through a user-material subroutine interface. Crack propagation is modeled as the failure of elements leading to loss of mechanical resistance to deformation. The tearing & fracture response of a viscoelastic material using a notch-in-plate sample deformed in simple tension is also simulated. By comparing the simulation results for a local damage criterion versus a truly nonlocal damage criterion, we show how the spread of Gibbs free energy around the crack tip influences the nature of the crack growth. As expected, the material fracture response described by a local damage criterion exhibits pathological mesh dependence whereas the response obtained using the nonlocal damage criterion is mesh objective regardless of mesh density, element type and orientation. This work presents two major computational advances: (1) the conversion of the rate form-based finite-deformation constitutive theory into an objective time-integration procedure is straightforward because the material time derivative-based constitutive equations are frame-invariant, and (2) the physically-observed crack initiation & propagation process in solids can be accurately & robustly simulated using simple numerical techniques such as the element failure method if a truly nonlocal fracture criterion is utilized. Finally, we benchmark results obtained from our computational framework and numerical simulations for modeling crack initiation & propagation to physical experimental data of mixed-mode crack propagation in a viscoelastic material.

Stress intensity factors (SIFs) are used in continuum fracture mechanics to quantify the stress fields surrounding a crack in a homogeneous material in the linear elastic regime. Critical values of the SIFs define an intrinsic measure of the resistance of a material to propagate a crack. At atomic scales, however, fracture occurs as a series of atomic bonds breaking, differing from the continuum description. As a consequence, a formal analog of the continuum SIFs calculated from atomistic simulations can have spatially localized, microstructural contributions that originate from varying bond configurations. The ability to characterize fracture at the atomic scale in terms of the SIFs offers both an opportunity to probe the effects of chemistry, as well as how the addition of a microstructural component affects the accuracy. We present a novel numerical method to determine SIFs from molecular dynamics (MD) simulations. The accuracy of this approach is first examined for a simple model, and then applied to atomistic simulations of fracture in amorphous silica. MD simulations provide time and spatially dependent SIFs, with results that are shown to be in good agreement with experimental values for fracture toughness in silica glass.

Isogeometric Kirchhoff–Love elements have received an increasing attention in geometrically nonlinear analysis of elastic shells. Nevertheless, some difficulties still remain. Among the others, the highly nonlinear expression of the strain measure, which leads to a complicated and costly computation of the discrete operators, and the existence of locking, which prevents the use of coarse meshes for slender shells and low order NURBS, are key issues that need to be addressed. In this work, exploiting the hypothesis of small membrane strains, we propose a simplified strain measure with a third order polynomial dependence on the displacement variables which allows an efficient evaluation of the discrete quantities. Numerical results show practically no difference to the original model, even for very large displacements and composite structures. Patch-wise reduced integrations are then investigated to deal with membrane locking in large deformation problems. An optimal integration scheme for third order C 2 NURBS, in terms of accuracy and efficiency, is identified. Finally, the recently proposed Newton method with mixed integration points is used for the solution of the discrete nonlinear equations with a great reduction of the iterative burden with respect to the standard Newton scheme.

In order to exclude spurious failure in compression the consideration of tension compression asymmetry (TCA) in damage models is of great interest in research as well as in industry. This paper presents an efficient and unifying algorithmic procedure for the incorporation of TCA into damage models. With the presented approach the implementation effort is drastically reduced since the equations and subroutines of the model without TCA can be reused. The strategy is first introduced and illustrated for an isotropic damage model. Afterwards, the procedure is applied to a gradient-extended anisotropic damage model with a second order damage tensor recently developed by the authors. For this model studies at integration point level as well as structural examples are presented. Two different aspects are included in the presented TCA approach and are demonstrated with the results of the simulations: (i) different damage evolution in tension and compression and (ii) stiffness recovery (crack closure) under compression. Due to the utilized gradient-extended formulation, being very efficient since only one additional degree of freedom is introduced, the finite element computations are shown to deliver mesh-objective results.

We develop superconvergent postprocessing techniques and accurate a-posteriori error estimators for variational multiscale formulations (VMS) of partial differential equations . We show that including VMS fine-scale information in the Quantity of Interest (QoI) approximation can provide a superior rate of convergence , including exact recovery for a non-trivial class of QoIs. Further improved postprocessing is achieved by reinterpreting the standard VMS fine-scale approximation via a local bubble-function projection. Decomposing the error in this modified QoI evaluation into an element-local component which depends only on the stabilization and a non-local dual-weighted residual component allows us to obtain accurate QoI error estimates. This analysis also leads to a new adjoint-based expression for the classic stabilization parameter, the intrinsic time ‘ τ ’. Numerical experiments on challenging boundary layer problems demonstrate the superior convergence of the modified QoI evaluation and the reliability of the adjoint error estimator for various QoIs. We also provide computational confirmation of exact recovery of an average value QoI using our estimator with the canonical ‘ τ ’ [Hughes (1995)], as predicted by our analysis.

This study presents a novel multi-scale computational method to analyze the dynamic thermo-mechanical performance of heterogeneous shell structures with orthogonal periodic configurations. The heterogeneities of heterogeneous shell structures are taken into account by periodic layouts of unit cells on the microscale in orthogonal curvilinear coordinate system. The new second-order two-scale approximate solutions for these multi-scale problems are constructed based on the multi-scale asymptotic analysis. Furthermore, the error estimates for the second-order two-scale (SOTS) solutions are obtained under some hypotheses. And then, a novel SOTS numerical algorithm based on finite element method (FEM), finite difference method (FDM) and decoupling method is brought forward in detail. Finally, some numerical examples are presented to verify the feasibility and validity of our multi-scale computational method. They also demonstrate that our multi-scale computational method can accurately capture the micro-scale dynamic thermo-mechanical responses in heterogeneous block structure, plate, cylindrical and doubly-curved shallow shells. In this paper, a unified multi-scale computational framework is established for dynamic thermo-mechanical problems of heterogeneous materials and structures with orthogonal periodic configurations. The asymptotic homogenization theory in Cartesian coordinate system and cylindrical coordinate system can be directly obtained based on the results in this paper.

A new 3-D numerical approach for the time dependent wave and heat equations as well as for the time independent Laplace equation on irregular domains with the Dirichlet boundary conditions has been developed. Trivial Cartesian meshes and simple 27-point uniform and nonuniform stencil equations are used for 3-D irregular domains. The calculation of the coefficients of the stencil equations is based on the minimization of the local truncation error of the stencil equations and yields the optimal order of accuracy of the new technique. Very small distances ( 0 . 1 h − 1 0 − 9 h where h is the grid size) between the grid points of a Cartesian mesh and the boundary do not worsen the accuracy of the new technique. At similar 27-point stencils, the accuracy of the new approach is two orders higher than that for the linear finite elements. The numerical results for irregular domains show that at the same number of degrees of freedom, the new approach is even much more accurate than the high-order (up to the fifth order) tetrahedral finite elements with much wider stencils. The wave and heat equations can be uniformly treated with the new approach. The order of the time derivative in these equations does not affect the coefficients of the stencil equations of the semi-discrete systems. The new approach can be directly applied to other partial differential equations.

In this work, an entirely automatized adaptive semi-explicit/explicit numerical approach is proposed to analyse nonlinear dynamic models. In this methodology, the time integrators are spatially and temporally locally evaluated, allowing different values to be assigned to each element of the spatial discretization and to each time step of the analysis. Thus, a very flexible procedure is enabled, permitting the technique to properly locally adapt to the properties of the model, as well as to the computed results. The novel algorithm is very accurate and extremely efficient, since it deals with reduced systems of equations and it requires no iterative procedures for solution. The main features of the novel methodology can be summarized as follows: (i) it is simple; (ii) it is locally defined; (iii) it enables stabilized analyses; (iv) it is an efficient non-iterative single-solution procedure; (v) it provides enhanced accuracy; (vi) it enables advanced controllable algorithmic dissipation in the higher modes; (vii) it considers a link between the temporal and the spatial discretizations; (viii) it stands as a single-solve framework based on reduced systems of equations; (ix) it provides an updating criterion for the tangent nonlinear matrices; (x) under special conditions, it may reproduce the Central Difference Method; (xi) it is entirely automated, requiring no expertise from the user.