Equations Of Linear Elasticity Research Articles

Gradient preserving Operator Inference: Data-driven reduced-order models for equations with gradient structure

Hamiltonian Operator Inference has been introduced in Sharma et al. (2022) to learn structure-preserving reduced-order models (ROMs) for Hamiltonian systems. This approach constructs a low-dimensional model using only data and knowledge of the Hamiltonian function. Such ROMs can keep the intrinsic structure of the system, allowing them to capture the physics described by the governing equations. In this work, we extend this approach to more general systems that are either conservative or dissipative in energy, and which possess a gradient structure. We derive the optimization problems for inferring structure-preserving ROMs that preserve the gradient structure. We further derive an a priori error estimate for the reduced-order approximation. To test the algorithms, we consider semi-discretized partial differential equations with gradient structure, such as the parameterized wave and Korteweg–de-Vries equations, and equations of three-dimensional linear elasticity in the conservative case and the one- and two-dimensional Allen-Cahn equations in the dissipative case. The numerical results illustrate the accuracy, structure-preservation properties, and predictive capabilities of the gradient-preserving Operator Inference ROMs.

A global space–time Trefftz DG scheme for the time-dependent isotropic elastic wave equations

In this paper we are concerned with Trefftz discretizations of the time-dependent linear elastic wave equation in arbitrary space dimensional domains Ω⊂Rd(d∈N). We propose a Trefftz discontinuos Galerkin method, construct Trefftz basis functions and prove that the corresponding approximate solutions possess optimal-order error estimates with respect to the meshwidth. A space–time domain decomposition preconditioner is introduced for the system generated by the proposed variational formula. Besides, we propose the global Trefftz DG method combined with local DG methods to solve the nonhomogeneous model. The numerical results verify the validity of the theoretical results, and show that the resulting approximate solutions possess high accuracy, and the proposed preconditioner is effective.

A finite element framework for fluid–structure interaction of turbulent cavitating flows with flexible structures

We present a finite element framework for the numerical prediction of cavitating turbulent flows interacting with flexible structures. The vapor-fluid phases are captured through a homogeneous mixture model, with a scalar transport equation governing the spatio-temporal evolution of cavitation dynamics. High-density gradients in the two-phase cavitating flow motivate the use of a positivity-preserving Petrov–Galerkin stabilization method in the variational framework. A mass transfer source term introduces local compressibility effects arising as a consequence of phase change. The turbulent fluid flow is modeled through a dynamic subgrid-scale method for large eddy simulations. The flexible structure is represented by a set of eigenmodes, obtained through a modal decomposition of the linear elasticity equations. While a partitioned iterative approach is adopted to couple the structural dynamics and cavitating fluid flow, the deforming flow domain is described by an arbitrary Lagrangian-Eulerian frame of reference. We establish the fidelity of the proposed framework by comparing it against experimental and numerical studies for both rigid and flexible hydrofoils in cavitating flows. Under unstable partial cavitating conditions, we identify specific vortical structures leading to cloud cavity collapse. We further explore features of cavitating flow past a rigid body such as re-entrant jet and turbulence-cavity interactions during cloud cavity collapse. Based on the validation study conducted over a flexible NACA66 rectangular hydrofoil, we elucidate the role of cavity and vortex shedding in governing the structural dynamics. Subsequently, we identify a broad spectrum frequency band whose central peak does not correlate to the frequency content of the cavitation dynamics or the natural frequencies of the structure, indicating the induction of unsteady flow patterns around the hydrofoil. Finally, we discuss the coupled fluid–structure dynamics during a cavitation cycle and the underlying mechanism associated with the promotion and mitigation of cavitation.

A thin-plate approximation for ocean wave interactions with an ice shelf

A variational principle is proposed to derive the governing equations for the problem of ocean wave interactions with a floating ice shelf, where the ice shelf is modelled by the full linear equations of elasticity and has an Archimedean draught. The variational principle is used to form a thin-plate approximation for the ice shelf, which includes water–ice coupling at the shelf front and extensional waves in the shelf, in contrast to the benchmark thin-plate approximation for ocean wave interactions with an ice shelf. The thin-plate approximation is combined with a single-mode approximation in the water, where the vertical motion is constrained to the eigenfunction that supports propagating waves. The new terms in the approximation are shown to have a major impact on predictions of ice shelf strains for wave periods in the swell regime.

A Model for Simulating the Upward Flow of a Viscous Fluid in a Fracture Network

Fluid migration in a fracture network plays an important role in the oil accumulation mechanism and hence is key to oil exploration. In this study, we build a model by combining one-dimensional (1D) Navier–Stokes equations, linear elastic equations, and energy equations, and validate the model by reproducing the thickness profile of a fluid-driven crack measured in an experiment. We employ this model to simulate the upward flow of viscous fluid in a single fracture during hydrocarbon migration. The simulation suggests that the parameters of both the fluid and the surrounding rock matrix, as well as the boundary condition imposed on the fracture outlet, affect the upward flow in the fracture. We then extend our model from the single fracture to the bifurcated fracture and the fracture network by maintaining homogeneous pressure and mass conservation at the connection of the channels. We find that the increase in network complexity leads to an increase in the inlet pressure gradient and inlet speed, and a decrease in the outlet pressure gradient and outlet speed. The effective area where the fluid is driven upward from the inlet to the outlet is offset toward the inlet. More importantly, the main novelty of our model is that it allows us to evaluate the effect of inconsistencies in individual branch parameters, such as matrix stiffness, permeability, temperature, and boundary conditions, on the overall upward flow of viscous fluid. Our results suggest that the heterogeneity enforces the greater impact on the closer branches.

An ANS/ATFs‐based unsymmetric solid‐shell finite element algorithm for high‐quality finite deformation analysis of hyper‐elastic shell

AbstractAn effective updated Lagrangian (UL) algorithm is designed for extending the recent distortion‐tolerant unsymmetric 8‐node, 24‐DOF hexahedral solid‐shell element, US‐ATFHS8, to finite deformation analysis of hyper‐elastic shell structures. The distinguishing feature of this unsymmetric element is that two different interpolation schemes are employed for virtual displacement and real stress calculations, respectively. The assumed natural strain (ANS) method with shell assumptions, referring to the current configuration, is introduced to modify the strain tensors derived from the assumed virtual displacement fields in terms of isoparametric coordinates, thereby mitigating shear locking and trapezoidal locking. On the other hand, the analytical trial functions (ATFs) derived from the general solutions of homogenous governing equations for linear elasticity are updated in each increment step to obtain the incremental deformation gradient, which is then utilized for calculating the real stresses for curing the numerical difficulties in large deformation problems. Numerical examples show that the proposed algorithm enables the hyper‐elastic solid‐shell element US‐ATFHS8 to exhibit excellent performance in both regular and distorted meshes and yield considerable results even when other models cannot work.

Penalty‐free discontinuous Galerkin method

AbstractIn this article, we present a new high‐order discontinuous Galerkin (DG) method, in which neither a penalty parameter nor a stabilization parameter is needed. We refer to this method as penalty‐free DG. In this method, the trial and test functions belong to the broken Sobolev space, in which the functions are in general discontinuous on the mesh skeleton and do not meet the Dirichlet boundary conditions. However, a subset can be distinguished in this space, where the functions are continuous and satisfy the Dirichlet boundary conditions, and this subset is called admissible. The trial solution is chosen to lie in an augmented admissible subset, in which a small violation of the continuity condition is permitted. This subset is constructed by applying special augmented constraints to the linear combination of finite element basis functions. In this approach, all the advantages of the DG method are retained without the necessity of using stability parameters or numerical fluxes. Several benchmark problems in two dimensions (Poisson equation, linear elasticity, hyperelasticity, and biharmonic equation) on polygonal (triangles, quadrilateral, and weakly convex polygons) meshes as well as a three‐dimensional Poisson problem on hexahedral meshes are considered. Numerical results are presented that affirm the sound accuracy and optimal convergence of the method in the norm and the energy seminorm.

Adaptive [formula omitted]-matrix computations in linear elasticity

This article deals with the efficient numerical treatment of the Lamé equations. The equations of linear elasticity are considered as boundary integral equations and solved in the setting of the boundary element method (BEM). Using BEM, one is faced with the solution of a system of equations with a fully populated system matrix, which is in general very costly. In order to overcome this difficulty, adaptive and approximate algorithms based on hierarchical matrices and the adaptive cross approximation are proposed. These new methods rely on error estimators and refinement techniques known from adaptivity but are not used here to improve the mesh. We apply these new techniques to both, the efficient solution of Lamé equations and to the multiplication with given data.

Numerical investigation of liquid sloshing in 2D flexible tanks subjected to complex external loading

Due to the complexity of fluid–structure interactions (FSI), the majority of studies in the literature dealing with the sloshing problem are restricted to rigid tanks. This paper is devoted to a numerical investigation of the liquid sloshing behavior in a flexible tank subjected to external loading. A numerical methodology is proposed, taking into account the FSI problem by coupling two open-source codes: OpenFOAM for the fluid and FEniCS for the solid, using the preCICE library, a free library for fluid–structure interaction. The Arbitrary Lagrangian–Eulerian formulation is used for the two-phase flow system to solve the Navier–Stokes equations in the fluid domain using the finite volume method. Simultaneously, the linear-elastic equation of the structure is solved using the finite element method. An implicit coupling scheme is considered at the fluid–structure interface. The numerical methodology is validated by the results given in literature for harmonic excitation at different frequencies. Subsequently, an analysis of complex external loading, such as Gabor wavelets and earthquake ground motion, is conducted to highlight the significant impact of the wall flexibility on sloshing, as well as the influence of hydrodynamic forces on the structure’s deformation. The proposed coupling methodology is robust and effective, it can be applied to all types of liquids and materials. A dataset of one of the studied cases is given as a supplement to the paper (Kha et al., 2024).

Discovering interpretable physical models using symbolic regression and discrete exterior calculus

Computational modeling is a key resource to gather insight into physical systems in modern scientific research and engineering. While access to large amount of data has fueled the use of machine learning to recover physical models from experiments and increase the accuracy of physical simulations, purely data-driven models have limited generalization and interpretability. To overcome these limitations, we propose a framework that combines symbolic regression (SR) and discrete exterior calculus (DEC) for the automated discovery of physical models starting from experimental data. Since these models consist of mathematical expressions, they are interpretable and amenable to analysis, and the use of a natural, general-purpose discrete mathematical language for physics favors generalization with limited input data. Importantly, DEC provides building blocks for the discrete analog of field theories, which are beyond the state-of-the-art applications of SR to physical problems. Further, we show that DEC allows to implement a strongly-typed SR procedure that guarantees the mathematical consistency of the recovered models and reduces the search space of symbolic expressions. Finally, we prove the effectiveness of our methodology by re-discovering three models of continuum physics from synthetic experimental data: Poisson equation, the Euler’s elastica and the equations of linear elasticity. Thanks to their general-purpose nature, the methods developed in this paper may be applied to diverse contexts of physical modeling.

An analytic model for the flow induced in syringomyelia cavities.

A simple two-dimensional fluid-structure-interaction problem, involving viscous oscillatory flow in a channel separated by an elastic membrane from a fluid-filled slender cavity, is analyzed to shed light on the flow dynamics pertaining to syringomyelia, a neurological disorder characterized by the appearance of a large tubular cavity (syrinx) within the spinal cord. The focus is on configurations in which the velocity induced in the cavity, representing the syrinx, is comparable to that found in the channel, representing the subarachnoid space surrounding the spinal cord, both flows being coupled through a linear elastic equation describing the membrane deformation. An asymptotic analysis for small stroke lengths leads to closed-form expressions for the leading-order oscillatory flow, and also for the stationary flow associated with the first-order corrections, the latter involving a steady distribution of transmembrane pressure. The magnitude of the induced flow is found to depend strongly on the frequency, with the result that for channel flow rates of non-sinusoidal waveform, as those found in the spinal canal, higher harmonics can dominate the sloshing motion in the cavity, in agreement with previous in vivo observations. Under some conditions, the cycle-averaged transmembrane pressure, also showing a marked dependence on the frequency, changes sign on increasing the cavity transverse dimension (i.e. orthogonal to the cord axis), underscoring the importance of cavity size in connection with the underlying hydrodynamics. The analytic results presented here can be instrumental in guiding future numerical investigations, needed to clarify the pathogenesis of syringomyelia cavities.

Optimal surrogate boundary selection and scalability studies for the shifted boundary method on octree meshes

The accurate and efficient simulation of Partial Differential Equations (PDEs) in and around arbitrarily defined geometries is critical for many application domains. Immersed boundary methods (IBMs) alleviate the usually laborious and time-consuming process of creating body-fitted meshes around complex geometry models (described by CAD or other representations, e.g., STL, point clouds), especially when high levels of mesh adaptivity are required. In this work, we advance the field of IBM in the context of the recently developed Shifted Boundary Method (SBM). In the SBM, the location where boundary conditions are enforced is shifted from the actual boundary of the immersed object to a nearby surrogate boundary, and boundary conditions are corrected utilizing Taylor expansions. This approach allows choosing surrogate boundaries that conform to a Cartesian mesh without losing accuracy or stability. Our contributions in this work are as follows: (a) we show that the SBM numerical error can be greatly reduced by an optimal choice of the surrogate boundary, (b) we mathematically prove the optimal convergence of the SBM for this optimal choice of the surrogate boundary, (c) we deploy the SBM on massively parallel octree meshes, including algorithmic advances to handle incomplete octrees, and (d) we showcase the applicability of these approaches with a wide variety of simulations involving complex shapes, sharp corners, and different topologies. Specific emphasis is given to Poisson’s equation and the linear elasticity equations.

Eringen’s model via linearization of nonlocal hyperelasticity

We consider Riesz’ fractional gradient and a truncated version of it. The equations of nonlocal nonlinear elasticity based on those gradients are known. We perform a formal linearization and arrive at the equations of linear elasticity based on those nonlocal operators. We prove the existence of solutions of the linear equations, notably, by a nonlocal version of Korn’s inequality. Finally, we show that the linearizations obtained are particular cases of Eringen’s model with singular kernels.

Parameter-free preconditioning for nearly-incompressible linear elasticity

It is well known that via the augmented Lagrangian method, one can solve Stokes' system by solving the nearly incompressible linear elasticity equation. In this paper, we show that the converse holds, and approximate the inverse of the linear elasticity operator with a convex linear combination of parameter-free operators. In such a way, we construct a uniform preconditioner for linear elasticity for all values of the Lamé parameter λ∈[0,∞). Numerical results confirm that by using inf-sup stable finite-element spaces for the solution of Stokes' equations, the proposed preconditioner is robust in λ.

Quantitative osteoporosis diagnosis of porous cancellous bone using poroelastodynamic modal analysis.

Osteoporosis is a skeletal disease characterized by reduced bone mass and microarchitectural deterioration, leading to increased fragility. This study presents a novel three-dimensional poroelastodynamic model for analyzing cancellous bone free vibration responses. The model incorporates the Navier-Stokes equations of linear elasticity and the Biot theory of porous media, allowing the investigation of osteoporosis-related changes. The analysis considers parameters like porosity, density, elasticity, Poisson ratio, and viscosity of bone marrow within the porous medium. Our findings indicate that natural frequencies of cancellous bone play a crucial role in osteoporosis prediction. By incorporating experimental data from 12 mouse femurs, we unveil insights into osteoporosis prediction. Increased porosity reduces bone stiffness, lowering natural frequencies. However, it also increases bone mass loss relative to stiffness, leading to higher frequencies. Therefore, the natural frequencies of osteoporotic bone are always higher than the natural frequencies of normal bone. Additionally, an increase in bone marrow within the pores, while increasing damping effects, also increases natural frequencies, which is another indication of osteoporosis growth in bone. The presence of bone marrow within the pores further influences natural frequencies, providing additional insights into osteoporosis growth. Thinner and smaller bones are found to be more susceptible to osteoporosis compared to larger and bigger bones due to their higher natural frequencies at equivalent porosity levels.

Analysis‐aware defeaturing of complex geometries with Neumann features

AbstractLocal modifications of a computational domain are often performed in order to simplify the meshing process and to reduce computational costs and memory requirements. However, removing geometrical features of a domain often introduces a non‐negligible error in the solution of a differential problem in which it is defined. In this work, we extend the results from Buffa et al. (Math Models Methods Appl Sci. 2022; 32(02):359–402.) by studying the case of domains containing an arbitrary number of distinct Neumann features, and by performing an analysis on Poisson's, linear elasticity, and Stokes' equations. We introduce a simple, computationally cheap, reliable, and efficient a posteriori estimator of the geometrical defeaturing error. Moreover, we also introduce a geometric refinement strategy that accounts for the defeaturing error: starting from a fully defeatured geometry, the algorithm determines at each iteration step which features need to be added to the geometrical model to reduce the defeaturing error. These important features are then added to the (partially) defeatured geometrical model at the next iteration, until the solution attains a prescribed accuracy. A wide range of two‐ and three‐dimensional numerical experiments are finally reported to illustrate this work.

Isogeometric Tearing and Interconnecting solvers for linear elasticity in multi-patch Isogeometric Analysis with theory for two dimensional domains

We consider the linear elasticity equation, discretized with multi-patch Isogeometric Analysis. A standard discretization error analysis is based on Korn’s inequality, which degrades for certain geometries, such as long and thin cantilevers. This phenomenon is known as geometry locking. We observe that high-order methods, like Isogeometric Analysis is beneficial in such a setting. The main focus of this paper is the construction and analysis of a domain decomposition solver, namely an Isogeometric Tearing and Interconnecting (IETI) solver, where we prove that the convergence behavior does not depend on the constant of Korn’s inequality for the overall domain, but only on the corresponding constants for the individual patches. Moreover, our analysis is explicit in the choice of the spline degree. Numerical experiments are provided which demonstrates the efficiency of the proposed solver.

Deformation of an elastic second gradient spherical body under equatorial line density of dead forces

We consider deformations of an elastic body having initially a spherical shape. Assumed deformation energy depends on the first and second gradient of displacements. We apply an equatorial line density of dead loads, that are forces per unit line length directed in radial direction and applied along the equator of the sphere. We restrict ourselves our analysis to the case of linearized second strain gradient isotropic elasticity (for which the more general energy was determined by Mindlin) with only one characteristic length. Differently to what happens in first gradient continua, i.e. in classic linear elasticity, we show that for the particular class second gradient continua considered here these forces do not determine infinite displacements in the direction of applied dead line forces. Instead, using a series method for the solution of the considered boundary-value problem, we demonstrate that the displacements are finite. So in the deformed configuration there is not the formation of an edge at the material points where the forces are applied. Further investigations are therefore needed for establishing if this elastic-regime edge formation is made possible: (I) either in the case of more general linear elastic constitutive equations or (II) only when large deformations are considered or (III) if non-elastic phenomena are involved.

Deep Ritz method with adaptive quadrature for linear elasticity

In this paper, we study the deep Ritz method for solving the linear elasticity equation from a numerical analysis perspective. A modified Ritz formulation using the H1/2(ΓD) norm is introduced and analyzed for linear elasticity equation in order to deal with the (essential) Dirichlet boundary condition. We show that the resulting deep Ritz method provides the best approximation among the set of deep neural network (DNN) functions with respect to the “energy” norm. Furthermore, we demonstrate that the total error of the deep Ritz simulation is bounded by the sum of the network approximation error and the numerical integration error, disregarding the algebraic error. To effectively control the numerical integration error, we propose an adaptive quadrature-based numerical integration technique with a residual-based local error indicator. This approach enables efficient approximation of the modified energy functional. Through numerical experiments involving smooth and singular problems, as well as problems with stress concentration, we validate the effectiveness and efficiency of the proposed deep Ritz method with adaptive quadrature.

Physics-informed neural networks for mesh deformation with exact boundary enforcement

In this work, we have applied physics-informed neural networks (PINN) for solving mesh deformation problems. We used the collocation PINN method to capture the new positions of the vertex nodes while preserving the connectivity information. We use linear elasticity equations for mesh deformation. To prevent vertex collisions or edge overlap, the mesh movement in this work is conducted in steps with relatively small movements. For moving boundary problems, the exact position of the boundary is essential for having an accurate solution. However, PINNs are frequently unable to satisfy Dirichlet boundary conditions exactly. To overcome this issue, we have used hard boundary condition enforcement to automatically satisfy Dirichlet boundary conditions. Specifically, we first trained a PINN with soft boundary conditions to obtain a particular solution. Then, this solution was tuned with exact boundary positions and a proper distance function by using a new PINN considering only the equation residual. To assess the accuracy of our approach, we used the classical translation and rotation tests and compared them with a proper mesh quality metric considering the change in the element area and shape. The results show the accuracy of this approach is comparable with that of finite element solutions. We also solved different moving boundary problems, resembling commonly used fluid–structure interaction problems. This work provides insight into using PINN for mesh-deformation problems without needing a discretization scheme with reasonable accuracy.