Boundary element modelling of dynamic behavior of piecewise homogeneous anisotropic elastic solids

A traditional direct boundary integral equations method is applied to solve three-dimensional dynamic problems of piecewise homogeneous linear elastic solids. The materials of homogeneous parts are considered to be generally anisotropic. The technique used to solve the boundary integral equations is based on the boundary element method applied together with the Radau IIA convolution quadrature method. A numerical example of suddenly loaded 3D prismatic rod consisting of two subdomains with different anisotropic elastic properties is presented to verify the accuracy of the proposed formulation.


Introduction
The analytical solutions of the dynamic problems of three-dimensional anisotropic and linear elasticity are severely restricted to a limited number of specific geometrical configurations and boundary conditions. More sophisticated problems, e.g., when nonhomogeneous regions are considered, can only be solved by suitable numerical tools.
Nowadays, the most popular domain-type numerical method is the finite element method (FEM). However, the boundary element method (BEM) is an interesting and highly-accurate alternative to FEM, due to its semi-analytical nature and the reduction of the problem dimension by one. In the recent decades, application of BEM for elastic wave propagation in isotropic elastic solids and media have received considerable attention. Two-and three-dimensional, timedomain, frequency-and Laplace-domain BEM formulations for transient isotropic elastodynamic problems are reported in abundance in the literature (see, e.g., [1][2][3][4][5][6][7]). Certain difficulties of development of the boundary element techniques for anisotropic and linear elastic solids are usually associated with the mathematical complexity of the corresponding dynamic fundamental solutions. They are not available in the explicit closed forms.
In this paper, a multi-domain three-dimensional boundary element formulation is applied for transient dynamic anisotropic and linear elastic multi-material, i.e., piecewise homogeneous, problems. The employed boundary element technique is formulated in the Laplace domain. The transient time-domain response is calculated via an inverse transformation by Radau IIA Runge-Kutta convolution quadrature method (CQM). Three-dimensional anisotropic elastodynamic fundamental solutions are expressed as surface integrals over half a unit sphere. Different interpolation functions are used to approximate the boundary variables. To demonstrate the reliability and accuracy of the proposed boundary element formulation, a numerical example is presented.
In equations (1)-(6), σ ij denote the Cauchy stress tensor, u i are displacement components, ρ is the mass density, C ijkl is the elastic stiffness tensor, and ε ij is the linear strain tensor; u * i and t * i denote the known displacement and traction vectors on the external boundaries Γ u and Γ t , respectively. For the outward unit normal vector n i , the traction vector t i is expressed as For subsequent derivations, we employ the Laplace transform in order to eliminate the time variablef with s denoting the Laplace transform parameter. Applying Laplace transform (8) to equations (1)-(7) and combining the corresponding transforms of equations (1)-(3), we obtain a set of decoupled Laplace domain problems To solve the boundary-value problem defined in equations (9)-(10), the problem is reformulated as boundary integral equations (BIEs) for the displacements: where r = y − x, with y and x being the field and source points, respectively; c ij (x) are the jump terms, p.v. stands for the Cauchy principal-value integral,ḡ ij andh jk are the Laplace domain fundamental solutions for displacement and traction, respectively. Following the usual boundary element procedure, the external boundary Γ = ∂Ω is discretized into quadrangular boundary elements. Quadratic shape functions are adopted to approximate the geometry of the elements. For the spatial discretization of the displacements and tractions over each boundary element, the linear and constant interpolation functions are used, respectively. The nodal collocation scheme is employed for the spatial discretization of boundary integral equations (12). For a piecewise homogeneous solid, which consists of a finite number distinct homogeneous subdomains, the discretization process of the BIEs (12) is performed for each homogeneous subdomain. Under the assumption of continuity of the displacements and equilibrium of tractions on the each interface of adjacent subdomains, the resulting complex-valued system of linear algebraic equations is obtained for a fixed value of Laplace transform parameter sĀ (s)f (s) =b(s), whereĀ(s) is the coefficient matrix, the vectorf(s) combines unknown boundary field variables (displacements and tractions), andb(s) is a known vector.

Fundamental solutions
The elastodynamic fundamental solutions for the anisotropic and homogeneous solids can be split into the sum of a singular static and a regular dynamic part For the three-dimensional displacement fundamental solution, g S ij (r) andḡ D ij (r, s) have the following expressions [8,9]

Numerical inversion of the Laplace transform
The convolution quadrature method was developed by Lubich [10,11] for numerical approximation of the convolution integral, and it can also be used as a numerical inverse Laplace transform method. In this paper, we use CQM based on a Radau IIA Runge-Kutta method [12] given by Butcher' where ∆t is the time step size, N is total number of time steps, and 0 < R < 1 is a CQM parameter.

Numerical example
A three-dimensional piecewise homogeneous rod of total size 1 × 1 × 3 m, as shown in figure 1, is considered. The rod is subjected to a Heaviside-type loading t 3 = t * 3 H(t), t * 3 = −10 5 Pa at one end and fixed at another. The rest of the surface is free of tractions.

Conclusions
In this study, we applied a Radau IIA Runge-Kutta CQM based boundary element method for transient dynamic analysis of generally anisotropic and linear elastic piecewise homogeneous three-dimensional solids. For the spatial discretization, a collocation method with mixed boundary elements is employed. The singular and regular parts of the anisotropic dynamic fundamental solutions are expressed as integrals over a unit circle and a half of a unit sphere, respectively. The comparison of the obtained boundary element results with the corresponding finite element solutions of the provided numerical example shows that the present boundary element formulation for multi-domain anisotropic elastodynamic problems is accurate and stable.