The paper extends the formulation of a 2D geometrically exact beam element proposed in our previous paper [1] to curved elastic beams. This formulation is based on equilibrium equations in their integrated form, combined with the kinematic relations and sectional equations that link the internal forces to sectional deformation variables. The resulting first-order differential equations are approximated by the finite difference scheme and the boundary value problem is converted to an initial value problem using the shooting method. The paper develops the theoretical framework based on the Navier-Bernoulli hypothesis but the approach could be extended to shear-flexible beams. The initial shape of the beam is captured with high accuracy, for certain shapes including the circular one even exactly. Numerical procedures for the evaluation of equivalent nodal forces and of the element tangent stiffness are presented in detail. Unlike standard finite element formulations, the present approach can increase accuracy by refining the integration scheme on the element level while the number of global degrees of freedom is kept constant. The efficiency and accuracy of the developed scheme are documented by five examples that cover circular and parabolic arches and a spiral-shaped beam. It is also shown that, for initially curved beams, a cross effect in the relations between internal forces and deformation variables arises, i.e., the bending moment affects axial stretching and the normal force affects the curvature.
In this paper we propose a novel approach to discretize linear port-Hamiltonian systems while preserving the underlying structure. We present a finite element exterior calculus formulation that is able to mimetically represent conservation laws and cope with mixed open boundary conditions using a single computational mesh. The possibility of including open boundary conditions allows for modular composition of complex multi-physical systems whereas the exterior calculus formulation provides a coordinate-free treatment. Our approach relies on a dual-field representation of the physical system that is redundant at the continuous level but eliminates the need of mimicking the Hodge star operator at the discrete level. By considering the Stokes-Dirac structure representing the system together with its adjoint, which embeds the metric information directly in the codifferential, the need for an explicit discrete Hodge star is avoided altogether. By imposing the boundary conditions in a strong manner, the power balance characterizing the Stokes-Dirac structure is then retrieved at the discrete level via symplectic Runge-Kutta integrators based on Gauss-Legendre collocation points. Numerical experiments validate the convergence of the method and the conservation properties in terms of energy balance both for the wave and Maxwell equations in a three dimensional domain. For the latter example, the magnetic and electric fields preserve their divergence free nature at the discrete level.
In our previous work [SIAM J. Sci. Comput. 43(3) (2021) B784-B810], an accurate hyper-singular boundary integral equation method for dynamic poroelasticity in two dimensions has been developed. This work is devoted to studying the more complex and difficult three-dimensional problems with Neumann boundary condition and both the direct and indirect methods are adopted to construct combined boundary integral equations. The strongly-singular and hyper-singular integral operators are reformulated into compositions of weakly-singular integral operators and tangential-derivative operators, which allow us to prove the jump relations associated with the poroelastic layer potentials and boundary integral operators in a simple manner. Relying on both the investigated spectral properties of the strongly-singular operators, which indicate that the corresponding eigenvalues accumulate at three points whose values are only dependent on two Lam\'e constants, and the spectral properties of the Calder\'on relations of the poroelasticity, we propose low-GMRES-iteration regularized integral equations. Numerical examples are presented to demonstrate the accuracy and efficiency of the proposed methodology by means of a Chebyshev-based rectangular-polar solver.
Targeting simulations on parallel hardware architectures, this paper presents computational kernels for efficient computations in mortar finite element methods. Mortar methods enable a variationally consistent imposition of coupling conditions at high accuracy, but come with considerable numerical effort and cost for the evaluation of the mortar integrals to compute the coupling operators. In this paper, we identify bottlenecks in parallel data layout and domain decomposition that hinder an efficient evaluation of the mortar integrals. We then propose a set of computational strategies to restore optimal parallel communication and scalability for the core kernels devoted to the evaluation of mortar terms. We exemplarily study the proposed algorithmic components in the context of three-dimensional large-deformation contact mechanics, both for cases with fixed and dynamically varying interface topology, yet these concepts can naturally and easily be transferred to other mortar applications, e.g. classical meshtying problems. To restore parallel scalability, we employ overlapping domain decompositions of the interface discretization independent from the underlying volumes and then tackle parallel communication for the mortar evaluation by a geometrically motivated reduction of ghosting data. Using three-dimensional contact examples, we demonstrate strong and weak scalability of the proposed algorithms up to 480 parallel processes as well as study and discuss improvements in parallel communication related to mortar finite element methods. For the first time, dynamic load balancing is applied to mortar contact problems with evolving contact zones, such that the computational work is well balanced among all parallel processors independent of the current state of the simulation.
A singularly perturbed parabolic problem of convection-diffusion type with a discontinuous initial condition is examined. An analytic function is identified which matches the discontinuity in the initial condition and also satisfies the homogenous parabolic differential equation associated with the problem. The difference between this analytical function and the solution of the parabolic problem is approximated numerically, using an upwind finite difference operator combined with an appropriate layer-adapted mesh. The numerical method is shown to be parameter-uniform. Numerical results are presented to illustrate the theoretical error bounds established in the paper.
In this work, we present a modification of explicit Runge-Kutta temporal integration schemes that guarantees the preservation of any locally-defined quasiconvex set of bounds for the solution. These schemes operate on the basis of a bijective mapping between an admissible set of solutions and the real domain to strictly enforce bounds. Within this framework, we show that it is possible to recover a wide range of methods independently of the spatial discretization, including positivity preserving, discrete maximum principle satisfying, entropy dissipative, and invariant domain preserving schemes. Furthermore, these schemes are proven to recover the order of accuracy of the underlying Runge-Kutta method upon which they are built. The additional computational cost is the evaluation of two nonlinear mappings which generally have closed-form solutions. We show the utility of this approach in numerical experiments using a pseudospectral spatial discretization without any explicit shock capturing schemes for nonlinear hyperbolic problems with discontinuities.
Computations of incompressible flows with velocity boundary conditions require solution of a Poisson equation for pressure with all Neumann boundary conditions. Discretization of such a Poisson equation results in a rank-deficient matrix of coefficients. When a non-conservative discretization method such as finite difference, finite element, or spectral scheme is used, such a matrix also generates an inconsistency which makes the residuals in the iterative solution to saturate at a threshold level that depends on the spatial resolution and order of the discretization scheme. In this paper, we examine inconsistency for a high-order meshless discretization scheme suitable for solving the equations on a complex domain. The high order meshless method uses polyharmonic spline radial basis functions (PHS-RBF) with appended polynomials to interpolate scattered data and constructs the discrete equations by collocation. The PHS-RBF provides the flexibility to vary the order of discretization by increasing the degree of the appended polynomial. In this study, we examine the convergence of the inconsistency for different spatial resolutions and for different degrees of the appended polynomials by solving the Poisson equation for a manufactured solution as well as the Navier-Stokes equations for several fluid flows. We observe that the inconsistency decreases faster than the error in the final solution, and eventually becomes vanishing small at sufficient spatial resolution. The rate of convergence of the inconsistency is observed to be similar or better than the rate of convergence of the discretization errors. This beneficial observation makes it unnecessary to regularize the Poisson equation by fixing either the mean pressure or pressure at an arbitrary point. A simple point solver such as the SOR is seen to be well-convergent, although it can be further accelerated using multilevel methods.
In this paper we are concerned with Trefftz discretizations of the time-dependent linear wave equation in anisotropic media in arbitrary space dimensional domains $\Omega \subset \mathbb{R}^d~ (d\in \mathbb{N})$. We propose two variants of the Trefftz DG method, define novel plane wave basis functions based on rigorous choices of scaling transformations and coordinate transformations, and prove that the corresponding approximate solutions possess optimal-order error estimates with respect to the meshwidth $h$ and the condition number of the coefficient matrices, respectively. Besides, we propose the global Trefftz DG method combined with local DG methods to solve the time-dependent linear nonhomogeneous wave equation in anisotropic media. In particular, the error analysis holds for the (nonhomogeneous) Dirichlet, Neumann, and mixed boundary conditions from the original PDEs. Furthermore, a strategy to discretize the model in heterogeneous media is proposed. The numerical results verify the validity of the theoretical results, and show that the resulting approximate solutions possess high accuracy.
The aim of this work is to provide the first strong convergence result of numerical approximation of a general time-fractional second order stochastic partial differential equation involving a Caputo derivative in time of order $\alpha\in(\frac 12; 1)$ and driven simultaneously by a multiplicative standard Brownian motion and additive fBm with Hurst parameter $H\in(\frac 12, 1)$, more realistic to model the random effects on transport of particles in medium with thermal memory. We prove the existence and uniqueness results and perform the spatial discretization using the finite element and the temporal discretization using a fractional exponential integrator scheme. We provide the temporal and spatial convergence proofs for our fully discrete scheme and the result shows that the convergence orders depend on the regularity of the initial data, the power of the fractional derivative, and the Hurst parameter $H$.
We present a theoretical and computational framework based on fractional calculus for the analysis of the nonlocal static response of cylindrical shell panels. The differ-integral nature of fractional derivatives allows an efficient and accurate methodology to account for the effect of long-range (nonlocal) interactions in curved structures. More specifically, the use of frame-invariant fractional-order kinematic relations enables a physically, mathematically, and thermodynamically consistent formulation to model the nonlocal elastic interactions. In order to evaluate the response of these nonlocal shells under practical scenarios involving generalized loads and boundary conditions, the fractional-Finite Element Method (f-FEM) is extended to incorporate shell elements based on the first-order shear-deformable displacement theory. Finally, numerical studies are performed exploring both the linear and the geometrically nonlinear static response of nonlocal cylindrical shell panels. This study is intended to provide a general foundation to investigate the nonlocal behavior of curved structures by means of fractional order models.
In this paper, a theoretical scheme is proposed for shape-programming of thin hyperelastic plates through differential growth. First, starting from the 3D governing system of a hyperelastic (neo-Hookean) plate, a consistent finite-strain plate equation system is formulated through a series-expansion and truncation approach. Based on the plate equation system, the problem of shape-programming is studied under the stress-free assumption. By equating the stress components in the plate equations to be zero, the explicit relations between growth functions and geometrical quantities of the target shape of the plate are derived. Then, a theoretical scheme of shape-programming is proposed, which can be used to identify the growth fields corresponding to arbitrary 3D shapes of the plate. To demonstrate the efficiency of the scheme, some typical examples are studied. The predicted growth functions in these examples are adopted in the numerical simulations, from which the target shapes of the plate can be recovered completely. The scheme of shape-programming proposed in the current work is applicable for manufacture of intelligent soft devices.