# Coupled Numerical Methods to Analyze Interacting Acoustic-Dynamic Models by Multidomain Decomposition Techniques

Authors: Delfim Soares
Journal: Mathematical Problems in Engineering
 (2011)

---

## Abstract

In this work, coupled numerical analysis of interacting acoustic and dynamic models is focused. In this context, several numerical methods, such as the finite difference method, the finite element method, the boundary element method, meshless methods, and so forth, are considered to model each subdomain of the coupled model, and multidomain decomposition techniques are applied to deal with the coupling relations. Two basic coupling algorithms are discussed here, namely the explicit direct coupling approach and the implicit iterative coupling approach, which are formulated based on explicit/implicit time-marching techniques. Completely independent spatial and temporal discretizations among the interacting subdomains are permitted, allowing optimal discretization for each sub-domain of the model to be considered. At the end of the paper, numerical results are presented, illustrating the performance and potentialities of the discussed methodologies.


---

## Body

## 1. Introduction

Usually, an engineer is faced with the analysis of a problem where two or more different physical systems interact with each other, so that the independent solution of any one system is impossible without simultaneous solution of the others. Such systems are known as coupled, and the intensity of such coupling is dependent on the degree of interaction [1]. Numerical algorithms consider that coupled systems may interact by means of common interfaces and/or overlapped subdomains. The former, usually referred to as interface-coupling, considers that coupling occurs on domain interfaces via the boundary conditions imposed there. Generally, distinct domains describe different physical situations, but it is possible to consider coupling between domains that are physically similar in which different discretization processes have been used. In the second case (problems in which the various domains totally or partially overlap), coupling occurs through the differential governing equations, describing different physical phenomena. In this work, only interface coupling problems are considered, and the interactions between acoustic fluids and elastodynamic solids are focused. In this context, one can mention a number of different applications: interaction between fluids such as air, water, or lubricants and structural elements such as buildings, dams, offshore structures, mechanical components pressure vessels, and so forth, systems composed of the same medium, with subdomains discretized by different numerical methods (finite difference, finite elements, boundary elements, etc.) and/or different refinement levels, and so forth.In the present work, several numerical methods are considered to discretize the different subdomains of the global model, taking into account interface coupled analyses. Although nowadays there are several powerful numerical techniques available, none of them can be considered most appropriate for all kinds of analysis, and, usually, the coupling of different numerical methodologies is necessary to analyze complex problems more effectively. In this context, the coupling of different numerical methods is recommended, in order to profit from their respective advantages and to evade their disadvantages. Two basic coupling algorithms are discussed here, considering multidomain decomposition techniques. In the first algorithm, explicit time-marching procedures are employed for wave propagation analysis at some subdomains of the model. Since explicit algorithms allow the computation of the current time-step response as function of only previous time-steps information; those subdomains can be independently analyzed directly, at each time step, allowing the development of an explicit direct coupling approach (ExDCA). On the other hand, when implicit time-marching procedures are considered, the computation of the current time-step response depends on the current time-step information, and interacting subdomains modeled by these techniques cannot be independently analyzed directly, being an iterative procedure necessary to analyze these coupled subdomains, once multidomain decomposition techniques are regarded. For this case, a second coupling algorithm is discussed here, referred to as implicit iterative coupling approach (ImICA).Taking into account an explicit direct or an implicit iterative multidomain decomposition technique, the coupling of several numerical procedures is carried out here. In this work, the coupling of the finite difference method (FDM), finite element method (FEM), boundary element method (BEM), and meshless local Petrov-Galerkin method (MLPG) is focused. In the last decades, these methodologies have been intensively applied to model acoustic-dynamic coupled models, taking into account different coupling strategies and time- and frequency-domain analyses. Considering the FDM, Vireaux [2] employed staggered grids to analyze acoustic-dynamic models in the 80s; nowadays, several advanced techniques are available based on the FDM, including those based on coupled methods [3–5]. In fact, it did not take long to couple numerical methods to analyze interacting acoustic-dynamic models, and most of these procedures are based on FEM-BEM coupling techniques [6–17] although there are several other procedures based on different numerical methodologies [18–30].When time-domain acoustic-dynamic coupled analyses are focused, the coupling of media with different properties (high properties contrast) and/or the coupling of numerical procedures with different spatial/temporal behavior may lead to inaccurate results or, even worse, instabilities. Thus, it is important to develop robust discretization techniques that not only are able to provide accurate and stable analyses, but also are computationally efficient. In this work, a multilevel time-step procedure is presented, as well as nonmatching interface nodes techniques are referred, allowing each subdomain of the model to be independently and optimally discretized, efficiently improving the accuracy and the stability of the analyses.The paper is organized as follows: first, basic equations concerning acoustic and dynamic models are presented, as well as interface interacting relations; in the sequence, numerical modeling of the acoustic/dynamic subdomains is briefly addressed taking into account domain- and boundary-discretization techniques. In Section4, coupling algorithms are discussed, focusing on explicit direct and implicit iterative procedures. At the end of the paper, three numerical applications (taking into account several different configurations) are presented, illustrating the performance and potentialities of the discussed methodologies.


## 2. Governing Equations

In this section, acoustic and elastic wave equations are briefly presented. Each one of these wave propagation models is used to mathematically describe different subdomains of the global problem. At the end of the section, basic equations concerning the coupling of acoustic and dynamic subdomains are described.


### 2.1. Acoustic Subdomains

The scalar wave equation is given by(2.1)(κp,i),i-ρp̈-ξṗ+S=0,
where p(X,t) stands for hydrodynamic pressure distribution and S(X,t) for body source terms. Inferior commas (indicial notation is adopted) and over dots indicate partial space (p,i=∂p/∂xi) and time (ṗ=∂p/∂t) derivatives, respectively. ξ(X) stands for the viscous damping coeficient; ρ(X) is the mass density, and κ(X) is the bulk modulus of the medium. In homogeneous media, ρ and κ are constant, and the classical wave equation (disregarding damping) can be written as(2.2)p,ii-p̈c2+s=0,
where c=κ/ρ is the wave propagation velocity. The boundary and initial conditions of the problem are given by(i) boundary conditions (t>0, X∈Γ where Γ=Γ1∪Γ2):(2.3a)p(X,t)=p̅(X,t)forX∈Γ1,(2.3b)q(X,t)=p,j(X,t)nj(X)=q̅(X,t)forX∈Γ2,(ii) initial conditions (t=0,X∈Γ∪Ω):(2.4a)p(X,0)=p¯0(X),(2.4b)ṗ(X,0)=p¯̇0(X),where the prescribed values are indicated by over bars, and q represents the flux along the boundary whose unit outward normal vector components are represented by nj. The boundary of the model is denoted by Γ (Γ1∪Γ2=Γ and Γ1∩Γ2=0) and the domain by Ω.


### 2.2. Dynamic Subdomains

The elastic wave equation for homogenous media is given by(2.5)(cd2-cs2)uj,ji+cs2ui,jj-üi-ζu̇i+bi=0,
where ui and bi stand for the displacement and the body force distribution components, respectively. The notation for time and space derivatives employed in (2.1) is once again adopted. In (2.5), cd is the dilatational wave velocity and cs is the shear wave velocity; they are given by cd2=(λ+2μ)/ρ and cs2=μ/ρ, where ρ is the mass density, and λ and μ are the Lamé’s constants. ζ stands for viscous damping-related parameters. Equation (2.5) can be obtained from the combination of the following basic mechanical equations (proper to model heterogeneous media):(2.6a)σij,j-ρüi-ρζu̇i+ρbi=0,(2.6b)σij=λδijεkk+2μεij,(2.6c)εij=ui,j+uj,i2,where σij and εij are, respectively, stress and strain tensor components, and δij is the Kronecker delta (δij=1, for i=j, and δij=0, for i≠j). Equation (2.6a) is the momentum equilibrium equation; (2.6b) represents the constitutive law of the linear elastic model, and (2.6c) stands for kinematical relations. The boundary and initial conditions of the elastodynamic problem are given by(i) boundary conditions (t>0, X∈Γ where Γ=Γ1∪Γ2): (2.7a)ui(X,t)=u̅i(X,t)forX∈Γ1,(2.7b)τi(X,t)=σij(X,t)nj(X)=τ̅i(X,t)forX∈Γ2,(ii) initial conditions (t=0, X∈Γ∪Ω):(2.8a)ui(X,0)=u̅i0(X),(2.8b)u̇i(X,0)=u̅̇i0(X),where the prescribed values are indicated by over bars, and τi denotes the traction vector along the boundary (nj, as indicated previously, stands for the components of the unit outward normal vector).


### 2.3. Acoustic-Dynamic Interacting Interfaces

On the acoustic-dynamic interface boundaries, the dynamic subdomain normal (normal to the interface) accelerations (ün) are related to the acoustic subdomain fluxes (q), and the acoustic subdomain hydrodynamic pressures (p) are related to the dynamic subdomain normal tractions (τn). These relations are expressed by the following equations:(2.9a)ün-(1ρ)q=0,(2.9b)τn+p=0,where in (2.9a) and (2.9b) the sign of the different subdomain outward normal directions is taken into account (outward normal vectors on the same interface point are opposite for each subdomain). In (2.9a), ρ is the mass density of the interacting acoustic subdomain medium.


## 2.1. Acoustic Subdomains

The scalar wave equation is given by(2.1)(κp,i),i-ρp̈-ξṗ+S=0,
where p(X,t) stands for hydrodynamic pressure distribution and S(X,t) for body source terms. Inferior commas (indicial notation is adopted) and over dots indicate partial space (p,i=∂p/∂xi) and time (ṗ=∂p/∂t) derivatives, respectively. ξ(X) stands for the viscous damping coeficient; ρ(X) is the mass density, and κ(X) is the bulk modulus of the medium. In homogeneous media, ρ and κ are constant, and the classical wave equation (disregarding damping) can be written as(2.2)p,ii-p̈c2+s=0,
where c=κ/ρ is the wave propagation velocity. The boundary and initial conditions of the problem are given by(i) boundary conditions (t>0, X∈Γ where Γ=Γ1∪Γ2):(2.3a)p(X,t)=p̅(X,t)forX∈Γ1,(2.3b)q(X,t)=p,j(X,t)nj(X)=q̅(X,t)forX∈Γ2,(ii) initial conditions (t=0,X∈Γ∪Ω):(2.4a)p(X,0)=p¯0(X),(2.4b)ṗ(X,0)=p¯̇0(X),where the prescribed values are indicated by over bars, and q represents the flux along the boundary whose unit outward normal vector components are represented by nj. The boundary of the model is denoted by Γ (Γ1∪Γ2=Γ and Γ1∩Γ2=0) and the domain by Ω.


## 2.2. Dynamic Subdomains

The elastic wave equation for homogenous media is given by(2.5)(cd2-cs2)uj,ji+cs2ui,jj-üi-ζu̇i+bi=0,
where ui and bi stand for the displacement and the body force distribution components, respectively. The notation for time and space derivatives employed in (2.1) is once again adopted. In (2.5), cd is the dilatational wave velocity and cs is the shear wave velocity; they are given by cd2=(λ+2μ)/ρ and cs2=μ/ρ, where ρ is the mass density, and λ and μ are the Lamé’s constants. ζ stands for viscous damping-related parameters. Equation (2.5) can be obtained from the combination of the following basic mechanical equations (proper to model heterogeneous media):(2.6a)σij,j-ρüi-ρζu̇i+ρbi=0,(2.6b)σij=λδijεkk+2μεij,(2.6c)εij=ui,j+uj,i2,where σij and εij are, respectively, stress and strain tensor components, and δij is the Kronecker delta (δij=1, for i=j, and δij=0, for i≠j). Equation (2.6a) is the momentum equilibrium equation; (2.6b) represents the constitutive law of the linear elastic model, and (2.6c) stands for kinematical relations. The boundary and initial conditions of the elastodynamic problem are given by(i) boundary conditions (t>0, X∈Γ where Γ=Γ1∪Γ2): (2.7a)ui(X,t)=u̅i(X,t)forX∈Γ1,(2.7b)τi(X,t)=σij(X,t)nj(X)=τ̅i(X,t)forX∈Γ2,(ii) initial conditions (t=0, X∈Γ∪Ω):(2.8a)ui(X,0)=u̅i0(X),(2.8b)u̇i(X,0)=u̅̇i0(X),where the prescribed values are indicated by over bars, and τi denotes the traction vector along the boundary (nj, as indicated previously, stands for the components of the unit outward normal vector).


## 2.3. Acoustic-Dynamic Interacting Interfaces

On the acoustic-dynamic interface boundaries, the dynamic subdomain normal (normal to the interface) accelerations (ün) are related to the acoustic subdomain fluxes (q), and the acoustic subdomain hydrodynamic pressures (p) are related to the dynamic subdomain normal tractions (τn). These relations are expressed by the following equations:(2.9a)ün-(1ρ)q=0,(2.9b)τn+p=0,where in (2.9a) and (2.9b) the sign of the different subdomain outward normal directions is taken into account (outward normal vectors on the same interface point are opposite for each subdomain). In (2.9a), ρ is the mass density of the interacting acoustic subdomain medium.


## 3. Numerical Modelling

Several numerical methods can be applied to discretize each subdomain of the coupled acoustic-dynamic model, according to their properties and advantages/disadvantages. In the following sub-sections, some numerical methods are briefly discussed, addressing their basic characteristics.


### 3.1. Domain-Discretization Methods

In the numerical methods based on domain discretization, the whole domain of the model is discretized into basic structures (elements, cells, points, etc.), and the spatial treatment of the governing equations is carried out considering these basic structures. In this case, matrix system of equations, as indicated in (3.1), is usually obtained, where the mass (M), damping (C), and stiffness (K) matrices, as well as the load vector (F), are computed according to the spatial discretization techniques being employed (3.1)MẌ(t)+CẊ(t)+KX(t)=F(t).In (3.1), X(t) stands for the pressure/displacement vector (X≡P or X≡U for acoustic or dynamic formulations, respectively) at time t (spatial and temporal discretizations are considered separately). In the present work, the finite difference method (FDM), the finite element method (FEM), and the meshless local Petrov-Galerkin method (MLPG) are focused, taking into account domain-discretization techniques.The FDM was one of the first methods developed to analyze complex problems governed by differential equations [31, 32]. It is easy to implement and considerably efficient; however, it may become extremely restricted when complex geometries are considered, because it is usually based on a regular distribution of points. The FEM, on the other hand, is well suited to analyze complex geometries, requiring in counterpart a considerably amount of input data [1, 33–36]. It is also quite an efficient technique, being the most popular method available nowadays to analyze intricate engineering problems. It is easy to implement and can be generalized to analyze complex models quite easily. Its main disadvantages are related to modelling unbounded domains and high gradient variations, as well as difficulties related to mesh generation. In the past few years, meshless methods have emerged essentially stimulated by these difficulties related to mesh generation [37, 38]. Mesh generation is delicate in many situations, for instance, when the domain has complicated geometry; when the mesh changes with time, as in crack propagation, and remeshing is required at each time step; when a Lagrangian formulation is employed, especially with nonlinear PDEs, and so forth. In addition, the need for flexibility in the selection of approximating functions (e.g., the flexibility to use nonpolynomial approximating functions) has played a significant role in the development of meshless methods (many meshless approximations give continuous variation of the first- or higher-order derivatives of a primitive function in counterpart to classical polynomial approximation where secondary fields have a jump on the interface of elements. Therefore, meshless approximations are leading to more accurate results in many cases). The main disadvantages of meshless methods are still their high computational costs and, in some cases, their lack of stability.Once the spatial treatment of the governing equations is carried out by a domain-discretization technique and (3.1) is obtained, its time domain analysis must also be considered. In this case, finite difference techniques are usually applied, rendering an algebraic system of equations, as described in (3.2), which must be solved at each time step n(3.2)AXn=Bm.In (3.2), A and B stand for the effective matrix and vector of the model, respectively, and the entries of X stand for the unknown variables. One should observe that vector B accounts for boundary prescribed conditions and domain sources, as well as some other previous step contributions (previous to m). Taking into account explicit time-marching techniques, m=n-1, whereas, for implicit time-marching techniques, m=n. In this work, several explicit and implicit techniques are considered. The central difference method and the Green-Newmark method [39–41], for instance, are explicit techniques that are here considered associated with the FDM and the FEM. Similarly, the Houbolt method [42] and the Newmark method [43] are implicit techniques that are here considered associated with the MLPG and the FEM.


### 3.2. Boundary-Discretization Methods

In boundary-discretization methods, just the boundary of the model is discretized, taking into account once again some basic structure, such as elements and point distributions. In this case, transient fundamental solutions are employed, and mixed approaches are focused, rendering numerical procedures based on more than one field incognita. The matrix system of equations that arises considering this kind of discretization can be written as(3.3)AXn=BYn+Zn,
where the entries of X and Y stand for the unknown and known (i.e., prescribed conditions) variables at the boundary of the model, respectively. A and B are effective matrices related to these variables at the current time step, and Z accounts for eventual domain-discretized terms (body sources, initial conditions etc.) and time convolution contributions.In the present work, the boundary element method (BEM) is focused as a boundary-discretization technique [44–47]. As it is well known, the BEM is well suited to analyze unbounded domains and to model high gradient variations, once it is based on fundamental solutions that satisfy the Sommerfeld radiation condition and that can properly deal with singularities in the model. The BEM is also flexible and efficient, allowing the discretization of complex geometries, as long as homogeneous media are considered. For heterogeneous media (or other more complex models, such as those considering anisotropy and nonlinear behavior), the BEM may be considered an inappropriate numerical tool, since, in these cases, its formulation may become very complex and prohibitive.There are also some “hybrid” formulations that are difficult to classify as a domain-or a boundary-discretization technique. This is the case, for instance, for some meshless techniques that are based on local boundary discretization (see, e.g., the LBIE—local boundary integral equation method [37]). In these meshless techniques, only boundary discretization is considered; however, the boundaries in focus are those of fictitious domains inside the real domain and, as a consequence, the whole real domain is in fact being discretized. Another hybrid formulation that is focused here is the domain boundary element method (DBEM) [40, 48, 49]. In this approach, nontransient fundamental solutions are considered, and the matrix system of equations that arises is a mix of (3.1) and (3.3), with some matrices being computed based on boundary discretizations and others being computed based on domain discretizations. Analogously as described in the previous sub-section, the DBEM also requires time-marching techniques to treat the time domain ordinary differential matrix equation that arises. Here, the Houbolt method is considered as such numerical technique.


## 3.1. Domain-Discretization Methods

In the numerical methods based on domain discretization, the whole domain of the model is discretized into basic structures (elements, cells, points, etc.), and the spatial treatment of the governing equations is carried out considering these basic structures. In this case, matrix system of equations, as indicated in (3.1), is usually obtained, where the mass (M), damping (C), and stiffness (K) matrices, as well as the load vector (F), are computed according to the spatial discretization techniques being employed (3.1)MẌ(t)+CẊ(t)+KX(t)=F(t).In (3.1), X(t) stands for the pressure/displacement vector (X≡P or X≡U for acoustic or dynamic formulations, respectively) at time t (spatial and temporal discretizations are considered separately). In the present work, the finite difference method (FDM), the finite element method (FEM), and the meshless local Petrov-Galerkin method (MLPG) are focused, taking into account domain-discretization techniques.The FDM was one of the first methods developed to analyze complex problems governed by differential equations [31, 32]. It is easy to implement and considerably efficient; however, it may become extremely restricted when complex geometries are considered, because it is usually based on a regular distribution of points. The FEM, on the other hand, is well suited to analyze complex geometries, requiring in counterpart a considerably amount of input data [1, 33–36]. It is also quite an efficient technique, being the most popular method available nowadays to analyze intricate engineering problems. It is easy to implement and can be generalized to analyze complex models quite easily. Its main disadvantages are related to modelling unbounded domains and high gradient variations, as well as difficulties related to mesh generation. In the past few years, meshless methods have emerged essentially stimulated by these difficulties related to mesh generation [37, 38]. Mesh generation is delicate in many situations, for instance, when the domain has complicated geometry; when the mesh changes with time, as in crack propagation, and remeshing is required at each time step; when a Lagrangian formulation is employed, especially with nonlinear PDEs, and so forth. In addition, the need for flexibility in the selection of approximating functions (e.g., the flexibility to use nonpolynomial approximating functions) has played a significant role in the development of meshless methods (many meshless approximations give continuous variation of the first- or higher-order derivatives of a primitive function in counterpart to classical polynomial approximation where secondary fields have a jump on the interface of elements. Therefore, meshless approximations are leading to more accurate results in many cases). The main disadvantages of meshless methods are still their high computational costs and, in some cases, their lack of stability.Once the spatial treatment of the governing equations is carried out by a domain-discretization technique and (3.1) is obtained, its time domain analysis must also be considered. In this case, finite difference techniques are usually applied, rendering an algebraic system of equations, as described in (3.2), which must be solved at each time step n(3.2)AXn=Bm.In (3.2), A and B stand for the effective matrix and vector of the model, respectively, and the entries of X stand for the unknown variables. One should observe that vector B accounts for boundary prescribed conditions and domain sources, as well as some other previous step contributions (previous to m). Taking into account explicit time-marching techniques, m=n-1, whereas, for implicit time-marching techniques, m=n. In this work, several explicit and implicit techniques are considered. The central difference method and the Green-Newmark method [39–41], for instance, are explicit techniques that are here considered associated with the FDM and the FEM. Similarly, the Houbolt method [42] and the Newmark method [43] are implicit techniques that are here considered associated with the MLPG and the FEM.


## 3.2. Boundary-Discretization Methods

In boundary-discretization methods, just the boundary of the model is discretized, taking into account once again some basic structure, such as elements and point distributions. In this case, transient fundamental solutions are employed, and mixed approaches are focused, rendering numerical procedures based on more than one field incognita. The matrix system of equations that arises considering this kind of discretization can be written as(3.3)AXn=BYn+Zn,
where the entries of X and Y stand for the unknown and known (i.e., prescribed conditions) variables at the boundary of the model, respectively. A and B are effective matrices related to these variables at the current time step, and Z accounts for eventual domain-discretized terms (body sources, initial conditions etc.) and time convolution contributions.In the present work, the boundary element method (BEM) is focused as a boundary-discretization technique [44–47]. As it is well known, the BEM is well suited to analyze unbounded domains and to model high gradient variations, once it is based on fundamental solutions that satisfy the Sommerfeld radiation condition and that can properly deal with singularities in the model. The BEM is also flexible and efficient, allowing the discretization of complex geometries, as long as homogeneous media are considered. For heterogeneous media (or other more complex models, such as those considering anisotropy and nonlinear behavior), the BEM may be considered an inappropriate numerical tool, since, in these cases, its formulation may become very complex and prohibitive.There are also some “hybrid” formulations that are difficult to classify as a domain-or a boundary-discretization technique. This is the case, for instance, for some meshless techniques that are based on local boundary discretization (see, e.g., the LBIE—local boundary integral equation method [37]). In these meshless techniques, only boundary discretization is considered; however, the boundaries in focus are those of fictitious domains inside the real domain and, as a consequence, the whole real domain is in fact being discretized. Another hybrid formulation that is focused here is the domain boundary element method (DBEM) [40, 48, 49]. In this approach, nontransient fundamental solutions are considered, and the matrix system of equations that arises is a mix of (3.1) and (3.3), with some matrices being computed based on boundary discretizations and others being computed based on domain discretizations. Analogously as described in the previous sub-section, the DBEM also requires time-marching techniques to treat the time domain ordinary differential matrix equation that arises. Here, the Houbolt method is considered as such numerical technique.


## 4. Coupling Procedures

In this work, the global model is divided in different subdomains, and each subdomain is analysed independently (as an uncoupled model), taking into account the numerical discretization techniques discussed in Section3. The interactions between the different subdomains of the global model are considered taking into account the accelerations/tractions and fluxes/pressures at the common interfaces, as well as the continuity equations (2.9a) and (2.9b). Two coupling procedures are discussed here, namely, (i) an explicit direct coupling approach (ExDCA); (ii) an implicit iterative coupling approach (ImICA).In the first procedure (i.e., the ExDCA), explicit time-marching schemes (e.g., the central difference method, the Green-Newmark method, etc.) are employed in some of the subdomains that are analyzed by domain-discretization methods. In the second procedure (ImICA), implicit time-marching schemes are considered within the subdomains. Since the ImICA is based on implicit algorithms (m=n in (3.2)), successive renewals of variables at common interfaces are considered in the coupling analysis (iterative coupling process), until convergence is achieved. On the other hand, the ExDCA is based on explicit algorithms (m=n-1 in (3.2)), and, as consequence, a direct coupling procedure can be developed, as it is described in the subsections that follow.For both explicit direct and implicit iterative coupling procedures, it is appropriate to consider different temporal discretizations within each subdomain. This is the case since optimal time steps are usually quite different taking into account dynamic and acoustic models, as well as different discretization techniques (especially taking into account some time-marching schemes that are conditionally stable). For instance, as it has been extensively reported in the literature, for small time steps, the time-domain BEM may become unstable, whereas, for large time-steps, excessive numerical damping may occur [44, 45]. Thus, in order to ensure stability and/or accuracy, usually a much smaller FEM time-step is required when coupled BEM-FEM analyses are considered (especially if the central difference method is employed associated to the FEM, which requires a low critical time-step). This situation may be amplified if subdomains with considerably different wave propagation velocities are interacting. In the next subsection, the adoption of different temporal discretizations within each subdomain of the global model is briefly discussed. In the sequence, the ExDCA and the ImICA are described.


### 4.1. Multilevel Time-Step Discretization

In order to consider different time steps in each subdomain, interpolation/extrapolation procedures along time are performed. Here, several schemes are considered for this temporal data manipulation, according to the discretization techniques involved. For instance, when the BEM is considered discretizing an interacting subdomain, temporal interpolation and extrapolation procedures are carried out based on the BEM time interpolation functions. In this case, time extrapolation procedures can be applied with confidence since they are consistent with the time-domain BEM formulation. Once time interpolation and extrapolation techniques are being employed, coupled implicit subdomains can be easily independently analysed (ImICA) taking into account different time steps. If explicit subdomains are considered (ExDCA), a subdomain solution can be computed independently of the current time step. As a consequence, just time interpolation procedures, associated with subcycling techniques, may be necessary if different time steps are required. Using these temporal data manipulations, optimal modelling in each subdomain may be achieved, which is very important regarding flexibility, efficiency, accuracy, and stability.


### 4.2. Explicit Direct Coupling

In the explicit direct coupling (as well as in the implicit iterative coupling), natural boundary conditions are prescribed at the acoustic and at the dynamic subdomains common interfaces. Two explicit direct coupling approaches are discussed here, the first one considering acoustic explicit subdomains and the second one considering dynamic explicit subdomains. For both approaches, the acoustic subdomain time steps are considered larger than the dynamic subdomain time steps (when different time-steps are regarded), since the wave propagation velocities in solids are usually higher than in acoustic fluids.In the first explicit direct coupling algorithm discussed here, the pressures related to the acoustic subdomains are computed directly, since their evaluation only takes into account results corresponding to previous time steps (m=n-1 in (3.2)). Once the acoustic pressures are evaluated, they are employed to compute tractions which are employed as prescribed interface boundary conditions (natural boundary condition) for the dynamic subdomains, and the displacements/velocities/accelerations of the model are computed by analysing these subdomains. The accelerations are then employed to evaluate the acoustic fluxes which are applied as prescribed interface boundary conditions (natural boundary condition) for the acoustic subdomains. If necessary, the time derivatives of the acoustic pressures can be computed. The next time-step computations are then initiated, repeating the above-described procedures.The detailed algorithm for this first ExDCA is presented in Table1, taking into account different temporal discretizations for the acoustic and for the dynamic subdomains (tp and tu, resp.—βj and ζj stand for time interpolation/extrapolation coefficients in the tables that follow). Space interpolation procedures may also be adopted in order to consider independent subdomain spatial discretizations (i.e., disconnected interface nodes); this can be accomplished by considering proper interface interpolating functions Nu(·) and Np(·), which are based on relations (2.9a) and (2.9b).Table 1
ExDCA-1 algorithm.
Time-step loop (based ontp)(1) Acoustic subdomains analyses: evaluation ofPtp.(2) Subcycling (untiltu=tp):(2.1) pressure temporal interpolation:Ptu=∑j=0JuβjP(t-jΔt)p,(2.2) force-pressure compatibility (spatial interpolation):Ftu=Nu(Ptu),(2.3) dynamic subdomains analyses: evaluation ofUtu,(2.4) evaluation of time derivatives ofUtu: U̇tu (if necessary), Ütu.(3) Flux-acceleration compatibility (spatial interpolation):Qtp=Np(Ütp).(4) Evaluation of time derivatives ofPtp: Ṗtp, P̈tp (if necessary).In this work, this first algorithm is employed associated to FEM-FEM coupled procedures in which the acoustic subdomains are modelled considering the Green-Newmark method (explicit technique), and the dynamic subdomains are modelled considering the Newmark method (implicit technique), as well as to FEM-FEM, and FEM-FDM coupled procedures in which all subdomains are modelled considering the central difference method (explicit technique).In the second explicit direct coupling algorithm focused here, the displacements related to the dynamic subdomains are computed directly, since their evaluation only takes into account results corresponding to previous time steps (m=n-1 in (3.2)). Once the displacements are evaluated, they are employed to compute the accelerations and, as a consequence, the acoustic fluxes, which are employed as prescribed interface boundary conditions (natural boundary condition) for the acoustic subdomains. The acoustic subdomains are then analyzed, and the acoustic pressures are computed. The pressures are then employed to evaluate the normal tractions which are applied as prescribed interface boundary conditions (natural boundary condition) for the dynamic subdomains. If necessary, the velocities of the model are computed. The next time-step computations are then initiated, repeating the above-described procedures.The detailed algorithm for this second ExDCA is presented in Table2, taking into account different temporal and spatial discretizations for the acoustic and for the dynamic subdomains. In this work, this methodology is considered applied to FEM-BEM coupled procedures in which acoustic subdomains are modelled by the BEM, and dynamic subdomains are modelled by the FEM associated to the Green-Newmark method (explicit technique).Table 2
ExDCA-2 algorithm.
Time-step loop (based ontu)(1) Dynamic subdomains analyses: evaluation ofUtu.(2) Evaluation of time derivatives ofUtu: Ütu.(3) Acceleration temporal extrapolation:Ütp=∑j=0JpζjÜ(t-jΔt)u.(4) Flux-acceleration compatibility (spatial interpolation):Qtp=Np(Ütp).(5) Acoustic subdomains analyses: evaluation ofPtp.(6) Pressure temporal interpolation:Ptu=∑j=0JuβjP(t-jΔt)p.(7) Force-pressure compatibility (spatial interpolation):Ftu=Nu(Ptu).(8) Evaluation of time derivatives ofUtu: U̇tu (if necessary).(9) Evaluation of time derivatives ofPtp: Ṗtp, P̈tp (if necessary, when tu=tp).


### 4.3. Implicit Iterative Coupling

In the implicit iterative approach, each subdomain of the model is analysed independently (as in the ExDCA), and a successive renewal of the variables at the common interfaces is performed, until convergence is achieved. In order to maximize the efficiency and robustness of the iterative coupling algorithm, the evaluation of an optimised relaxation parameter is introduced, taking into account the minimisation of a square error functional.Initially, in the ImICA, the dynamic subdomains are analysed and the displacements at the common interfaces are evaluated, as well as its time derivatives. A relaxation parameterα is introduced in order to ensure and/or to speed up convergence, such that (superscript k stands for the iterative step)(4.1)Xt(k+1)=(α)Xt(k+α)+(1-α)Xt(k),
where the relaxation parameter can be introduced associated to the displacement variable (X≡U) or to the acceleration variable (X≡Ü). Once the accelerations are computed, they are employed to calculate the acoustic fluxes, which are prescribed as interface boundary conditions (natural boundary condition) for the acoustic subdomains. The acoustic subdomains are then analyzed, and the pressures of the model are computed, which are employed to evaluate dynamic forces at the common interfaces (natural boundary condition). The dynamic subdomains are then once again analyzed, repeating the whole process until convergence is achieved. Once convergence is achieved, the next time-step computations are initiated, repeating the above-described procedures.The algorithm representing the ImICA is presented in Table3, taking into account different temporal and spatial discretizations within each subdomain of the model. In this work, this algorithm is employed associated to FEM-BEM, BEM-BEM (as well as DBEM-BEM), and MLPG-MLPG coupling procedures (for MLPG-MLPG coupled analyses, different time-steps techniques are not considered here).Table 3
ImICA algorithm.
Time-step loop (based ontu)(1) Iterative analysis (until convergence):(1.1) dynamic subdomains analyses: evaluation ofUtu(k+α),(1.2 or 1.3) evaluation of time derivatives ofUtu(k+λ):Ütu(k+λ),(1.3 or 1.2) adoption of a relaxation parameter:Xtu(k+1)=αXtu(k+α)+(1-α)Xtu(k),(1.4) acceleration temporal extrapolation:Ütp(k+1)=ζ0Ütu(k+1)+∑j=1JpζjÜ(t-jΔt)u,(1.5) flux-acceleration compatibility (spatial interpolation):Qtp(k+1)=Np(Ütp(k+1)),(1.6) acoustic subdomains analyses: evaluation ofPtp(k+1),(1.7) pressure temporal interpolation:Ptu(k+1)=β0Ptp(k+1)+∑j=1JuβjP(t-jΔt)p,(1.8) force-pressure compatibility (spatial interpolation):Ftu(k+1)=Nu(Ptu(k+1)),(2) Evaluation of time derivatives ofUtu: U̇tu (if necessary).(3) Evaluation of time derivatives ofPtp: Ṗtp, P̈tp (if necessary, when tu=tp).The effectiveness of the iterative coupling methodology is intimately related to the relaxation parameter selection; an inappropriate selection forα can drastically increase the number of iterations in the analysis or, even worse, make convergence unfeasible. Once appropriate α values are considered, convergence is usually achieved in quite few iterative steps, providing an efficient and robust iterative coupling technique. In order to evaluate an optimal relaxation parameter, the following square error functional is here minimized:(4.2)f(α)=‖Xt(α)(k+1)-Xt(k)(α)‖2.Substituting (4.1) into (4.2) yields(4.3)f(α)=‖αWt(k+α)+(1-α)Wt(k)‖2=(α2‖Wt(k+α)‖2+2α(1-α)(Wt(k+α),Wt(k))+(1-α)2‖Wt(k)‖2),
where the inner product definition is employed (e.g., (W,W)=||W||2) and new variables, as defined in (4.4),  are considered(4.4)Wt(k+λ)=Xt(k+λ)-Xt(k+λ-1).To find the optimalα that minimizes the functional f(α), (4.3) is differentiated with respect to α, and the result is set to zero, yielding(4.5)α=(Wt(k),Wt(k)-Wt(k+α))‖Wt(k)-Wt(k+α)‖2,
which is an efficient and easy to implement expression that provides an optimal value for the relaxation parameter α, at each iterative step. It is important to note that the relation 0<α≤1 must hold. In the present work, the optimal relaxation parameter is evaluated according to (4.5) and if α∉(0.01;1.00), the previous iterative-step relaxation parameter is adopted. For the first iterative step, α=0.5 is selected.


## 4.1. Multilevel Time-Step Discretization

In order to consider different time steps in each subdomain, interpolation/extrapolation procedures along time are performed. Here, several schemes are considered for this temporal data manipulation, according to the discretization techniques involved. For instance, when the BEM is considered discretizing an interacting subdomain, temporal interpolation and extrapolation procedures are carried out based on the BEM time interpolation functions. In this case, time extrapolation procedures can be applied with confidence since they are consistent with the time-domain BEM formulation. Once time interpolation and extrapolation techniques are being employed, coupled implicit subdomains can be easily independently analysed (ImICA) taking into account different time steps. If explicit subdomains are considered (ExDCA), a subdomain solution can be computed independently of the current time step. As a consequence, just time interpolation procedures, associated with subcycling techniques, may be necessary if different time steps are required. Using these temporal data manipulations, optimal modelling in each subdomain may be achieved, which is very important regarding flexibility, efficiency, accuracy, and stability.


## 4.2. Explicit Direct Coupling

In the explicit direct coupling (as well as in the implicit iterative coupling), natural boundary conditions are prescribed at the acoustic and at the dynamic subdomains common interfaces. Two explicit direct coupling approaches are discussed here, the first one considering acoustic explicit subdomains and the second one considering dynamic explicit subdomains. For both approaches, the acoustic subdomain time steps are considered larger than the dynamic subdomain time steps (when different time-steps are regarded), since the wave propagation velocities in solids are usually higher than in acoustic fluids.In the first explicit direct coupling algorithm discussed here, the pressures related to the acoustic subdomains are computed directly, since their evaluation only takes into account results corresponding to previous time steps (m=n-1 in (3.2)). Once the acoustic pressures are evaluated, they are employed to compute tractions which are employed as prescribed interface boundary conditions (natural boundary condition) for the dynamic subdomains, and the displacements/velocities/accelerations of the model are computed by analysing these subdomains. The accelerations are then employed to evaluate the acoustic fluxes which are applied as prescribed interface boundary conditions (natural boundary condition) for the acoustic subdomains. If necessary, the time derivatives of the acoustic pressures can be computed. The next time-step computations are then initiated, repeating the above-described procedures.The detailed algorithm for this first ExDCA is presented in Table1, taking into account different temporal discretizations for the acoustic and for the dynamic subdomains (tp and tu, resp.—βj and ζj stand for time interpolation/extrapolation coefficients in the tables that follow). Space interpolation procedures may also be adopted in order to consider independent subdomain spatial discretizations (i.e., disconnected interface nodes); this can be accomplished by considering proper interface interpolating functions Nu(·) and Np(·), which are based on relations (2.9a) and (2.9b).Table 1
ExDCA-1 algorithm.
Time-step loop (based ontp)(1) Acoustic subdomains analyses: evaluation ofPtp.(2) Subcycling (untiltu=tp):(2.1) pressure temporal interpolation:Ptu=∑j=0JuβjP(t-jΔt)p,(2.2) force-pressure compatibility (spatial interpolation):Ftu=Nu(Ptu),(2.3) dynamic subdomains analyses: evaluation ofUtu,(2.4) evaluation of time derivatives ofUtu: U̇tu (if necessary), Ütu.(3) Flux-acceleration compatibility (spatial interpolation):Qtp=Np(Ütp).(4) Evaluation of time derivatives ofPtp: Ṗtp, P̈tp (if necessary).In this work, this first algorithm is employed associated to FEM-FEM coupled procedures in which the acoustic subdomains are modelled considering the Green-Newmark method (explicit technique), and the dynamic subdomains are modelled considering the Newmark method (implicit technique), as well as to FEM-FEM, and FEM-FDM coupled procedures in which all subdomains are modelled considering the central difference method (explicit technique).In the second explicit direct coupling algorithm focused here, the displacements related to the dynamic subdomains are computed directly, since their evaluation only takes into account results corresponding to previous time steps (m=n-1 in (3.2)). Once the displacements are evaluated, they are employed to compute the accelerations and, as a consequence, the acoustic fluxes, which are employed as prescribed interface boundary conditions (natural boundary condition) for the acoustic subdomains. The acoustic subdomains are then analyzed, and the acoustic pressures are computed. The pressures are then employed to evaluate the normal tractions which are applied as prescribed interface boundary conditions (natural boundary condition) for the dynamic subdomains. If necessary, the velocities of the model are computed. The next time-step computations are then initiated, repeating the above-described procedures.The detailed algorithm for this second ExDCA is presented in Table2, taking into account different temporal and spatial discretizations for the acoustic and for the dynamic subdomains. In this work, this methodology is considered applied to FEM-BEM coupled procedures in which acoustic subdomains are modelled by the BEM, and dynamic subdomains are modelled by the FEM associated to the Green-Newmark method (explicit technique).Table 2
ExDCA-2 algorithm.
Time-step loop (based ontu)(1) Dynamic subdomains analyses: evaluation ofUtu.(2) Evaluation of time derivatives ofUtu: Ütu.(3) Acceleration temporal extrapolation:Ütp=∑j=0JpζjÜ(t-jΔt)u.(4) Flux-acceleration compatibility (spatial interpolation):Qtp=Np(Ütp).(5) Acoustic subdomains analyses: evaluation ofPtp.(6) Pressure temporal interpolation:Ptu=∑j=0JuβjP(t-jΔt)p.(7) Force-pressure compatibility (spatial interpolation):Ftu=Nu(Ptu).(8) Evaluation of time derivatives ofUtu: U̇tu (if necessary).(9) Evaluation of time derivatives ofPtp: Ṗtp, P̈tp (if necessary, when tu=tp).


## 4.3. Implicit Iterative Coupling

In the implicit iterative approach, each subdomain of the model is analysed independently (as in the ExDCA), and a successive renewal of the variables at the common interfaces is performed, until convergence is achieved. In order to maximize the efficiency and robustness of the iterative coupling algorithm, the evaluation of an optimised relaxation parameter is introduced, taking into account the minimisation of a square error functional.Initially, in the ImICA, the dynamic subdomains are analysed and the displacements at the common interfaces are evaluated, as well as its time derivatives. A relaxation parameterα is introduced in order to ensure and/or to speed up convergence, such that (superscript k stands for the iterative step)(4.1)Xt(k+1)=(α)Xt(k+α)+(1-α)Xt(k),
where the relaxation parameter can be introduced associated to the displacement variable (X≡U) or to the acceleration variable (X≡Ü). Once the accelerations are computed, they are employed to calculate the acoustic fluxes, which are prescribed as interface boundary conditions (natural boundary condition) for the acoustic subdomains. The acoustic subdomains are then analyzed, and the pressures of the model are computed, which are employed to evaluate dynamic forces at the common interfaces (natural boundary condition). The dynamic subdomains are then once again analyzed, repeating the whole process until convergence is achieved. Once convergence is achieved, the next time-step computations are initiated, repeating the above-described procedures.The algorithm representing the ImICA is presented in Table3, taking into account different temporal and spatial discretizations within each subdomain of the model. In this work, this algorithm is employed associated to FEM-BEM, BEM-BEM (as well as DBEM-BEM), and MLPG-MLPG coupling procedures (for MLPG-MLPG coupled analyses, different time-steps techniques are not considered here).Table 3
ImICA algorithm.
Time-step loop (based ontu)(1) Iterative analysis (until convergence):(1.1) dynamic subdomains analyses: evaluation ofUtu(k+α),(1.2 or 1.3) evaluation of time derivatives ofUtu(k+λ):Ütu(k+λ),(1.3 or 1.2) adoption of a relaxation parameter:Xtu(k+1)=αXtu(k+α)+(1-α)Xtu(k),(1.4) acceleration temporal extrapolation:Ütp(k+1)=ζ0Ütu(k+1)+∑j=1JpζjÜ(t-jΔt)u,(1.5) flux-acceleration compatibility (spatial interpolation):Qtp(k+1)=Np(Ütp(k+1)),(1.6) acoustic subdomains analyses: evaluation ofPtp(k+1),(1.7) pressure temporal interpolation:Ptu(k+1)=β0Ptp(k+1)+∑j=1JuβjP(t-jΔt)p,(1.8) force-pressure compatibility (spatial interpolation):Ftu(k+1)=Nu(Ptu(k+1)),(2) Evaluation of time derivatives ofUtu: U̇tu (if necessary).(3) Evaluation of time derivatives ofPtp: Ṗtp, P̈tp (if necessary, when tu=tp).The effectiveness of the iterative coupling methodology is intimately related to the relaxation parameter selection; an inappropriate selection forα can drastically increase the number of iterations in the analysis or, even worse, make convergence unfeasible. Once appropriate α values are considered, convergence is usually achieved in quite few iterative steps, providing an efficient and robust iterative coupling technique. In order to evaluate an optimal relaxation parameter, the following square error functional is here minimized:(4.2)f(α)=‖Xt(α)(k+1)-Xt(k)(α)‖2.Substituting (4.1) into (4.2) yields(4.3)f(α)=‖αWt(k+α)+(1-α)Wt(k)‖2=(α2‖Wt(k+α)‖2+2α(1-α)(Wt(k+α),Wt(k))+(1-α)2‖Wt(k)‖2),
where the inner product definition is employed (e.g., (W,W)=||W||2) and new variables, as defined in (4.4),  are considered(4.4)Wt(k+λ)=Xt(k+λ)-Xt(k+λ-1).To find the optimalα that minimizes the functional f(α), (4.3) is differentiated with respect to α, and the result is set to zero, yielding(4.5)α=(Wt(k),Wt(k)-Wt(k+α))‖Wt(k)-Wt(k+α)‖2,
which is an efficient and easy to implement expression that provides an optimal value for the relaxation parameter α, at each iterative step. It is important to note that the relation 0<α≤1 must hold. In the present work, the optimal relaxation parameter is evaluated according to (4.5) and if α∉(0.01;1.00), the previous iterative-step relaxation parameter is adopted. For the first iterative step, α=0.5 is selected.


## 5. Numerical Aspects and Applications

In the following sub-sections, some numerical applications are presented, illustrating the performance and potentialities of the discussed coupling methodologies. In the first application, a multidomain column is analyzed, considering several geometrical and physical configurations, as well as coupling procedures. In this case, acoustic-acoustic, acoustic-dynamic, and dynamic-dynamic coupled models are discussed, taking into account axisymmetric, two-dimensional, and three-dimensional configurations. In the second application, a dam-reservoir system is analyzed, considering once more several coupling techniques. In this case, a two-dimensional model is focused, and some advanced analyses are carried out, such as the modeling of nonlinear behavior and infinite media. In the last application, a tube of steel, submerged in water, is analyzed. In this case, axisymmetric models are focused, and, once again, several geometric and numeric configurations are considered. Along the applications discussed here, a large scope of coupling procedures is presented, namely: (i) for the ExDCA—FEM-FEM, FEM-BEM and FEM-FDM coupling procedures; (ii) for the ImICA—FEM-BEM, DBEM-BEM (which is referred to here as BEM-BEM 1), BEM-BEM (which is referred to here as BEM-BEM 2) and MLPG-MLPG coupling procedures. In this way, the reader can compare and better visualize some benefits and drawbacks of each methodology, considering an ample range of configurations.


### 5.1. Multidomain Column

The first example is that of a prismatic body behaving like a one-dimensional column. Initially, the column is analysed as an acoustic model [50]. It is fixed at one end (p̅(t)=0) and subjected to a unitary Heaviside type forcing function acting at the opposite end (q̅(t)=H(t)). A sketch of the model is shown in Figure 1(a). The material properties of the column are c=1.0m/s and ρ=1.0kg/m3. The geometry of the model is defined by L=1.0m. As depicted in Figure 1(a), 28 boundary elements of equal length and 40 quadrilateral finite elements are employed in the coupled mesh. Regarding time discretization, three different cases of analysis are considered here, namely, (i) ΔtF=1.0ΔtB, (ii) ΔtF=0.2ΔtB, and (iii) ΔtF=0.1ΔtB; where ΔtB=0.06s, and the subscripts F and B are related to the FEM and to the BEM, respectively.Figure 1
Column model: (a) FEM-BEM acoustic-acoustic two-dimensional model; (b) MLPG-MLPG and FEM-FEM acoustic-dynamic two-dimensional model; (c) FEM-BEM acoustic-acoustic axisymmetric model; (d) FEM-FEM and FEM-FDM acoustic-acoustic/acoustic-dynamic/dynamic-dynamic three-dimensional model; (e) sketch of the three-dimensional model spatial discretization.In Figure2, time history results are depicted, at points A and B (see Figure 1(a)), taking into account the ExDCA and the ImICA. Potential (pressure) and flux results are presented considering the three different cases of analysis, and they are compared to the analytical solution [51], plotted as a solid line. As can be seen, a higher level of accuracy is observed when different time steps are considered within each subdomain, regarding their optimal temporal discretization. The robustness of the multilevel time-step algorithm must be highlighted: as illustrated in the present application, the algorithm deals properly with highly different subdomain temporal discretizations.Time-history results at points A and B taking into account FEM-BEM coupling procedures and different temporal discretizations for each subdomain: (a) explicit direct coupling analysis; (b) implicit iterative coupling analysis.
(a)
ExDCA(b)
ImICAIn a second approach for the column model, the acoustic-dynamic coupled problem is focused (fluid-solid column [27, 28]). A sketch of the problem is depicted in Figure 1(b). The geometry of the model is defined once again by L=1.0m, and the column is submitted to a time Heavisite force acting at one of its ends. The physical properties of the media are (i) fluid subdomain: κ=100N/m2 (bulk modulus) and ρ=1kg/m3 (density); (ii) solid subdomain: E=100N/m2 (Young modulus), ν=0 (Poisson rate), and ρ=1kg/m3 (density).Two spatial-temporal MLPG discretizations are considered to analyse the model, namely: (i) discretization 1—153 nodes are employed to spatially discretize each subdomain, and the time step adopted isΔt=0.0025s; (ii) discretization 2—561 nodes are employed to spatially discretize each subdomain, and the time-step adopted is Δt=0.00125s. In Figure 3, displacement time-history results at points A and B of the solid subdomain and hydrodynamic pressure time-history results at points B and C of the fluid subdomain are plotted, considering discretizations 1 and 2 and the MLPG-MLPG ImICA. Analytical time histories are also depicted in Figure 3, highlighting the good accuracy of the numerical results.Time-history results at points A, B, and C taking into account MLPG-MLPG coupling procedures (ImICA) and different refinement levels: (a) discretization 1; (b) discretization 2.
(a)(b)The same fluid-solid column is analysed considering FEM-FEM coupled procedures based on the Green-Newmark method and on the Newmark method (ExDCA). In this case, 200 square finite elements are employed to discretize each subdomain of the model, and the time discretization is specified byΔt=0.002s. The obtained results are depicted in Figure 4(a). In this analysis, the same wave propagation velocities are considered within the solid and fluid subdomains, that is, cd=c=10m/s. Taking into account a steel-water column, the physical properties of the media are (i) fluid subdomain: κ=2.3175·109N/m2, ρ=1030kg/m3; (ii) solid subdomain: E=2.1·1011N/m2, ν=0.3, and ρ=7700kg/m3; in this case, cd=6059m/s and c=1500m/s. In Figure 4(b), the displacement and hydrodynamic pressure time-history results at point B are depicted, considering Δtf=Δts and Δtf=4Δts, where Δts=3.33·10-6s, and the subscripts s and f are related to the solid and fluid subdomains, respectively. As one can observe, good agreement between the results is observed, in spite of the different time discretizations adopted within each subdomain.Time-history results for the solid-fluid column at points A, B, and C taking into account FEM-FEM coupling procedures (ExDCA) and different physical models: (a) homogeneous wave propagation velocities; (b) heterogeneous wave propagation velocities.
(a)(b)In a third approach for the column model, the propagation of acoustic waves through a prismatic circular column is analysed (axisymmetric model [16]). A sketch of the problem is depicted in Figure 1(c). The downside half of the column is modelled by the BEM, and the upside half is modelled by the FEM. The properties of the media are (i) Medium1 (BEM)—ρ=1.0kg/m3 and c=1.0m/s; (ii) Medium2 (FEM)—ρ=1.0kg/m3 and c=2.0m/s. The spatial discretization is considered as follows: 32 linear boundary elements discretize Medium1, and 64 linear-square finite elements discretize Medium2. The time steps selected are ΔtB=0.60s and ΔtF=0.15s. Time history acoustic pressures at points A and B (see Figure 1(c)) are depicted in Figure 5. The numerical results are compared with those of a finite difference solution (based on the work presented by Cohen and Joly [52]) and, as one may observe, the time responses are quite similar considering these two different methodologies.Time-history results for acoustic pressures taking into account FEM-BEM coupling procedures (ImICA) applied to a heterogeneous axisymmetric model: (a) results at point A; (b) results at point B.
(a)(b)In a fourth approach for the column model, three-dimensional analyses are considered, taking into account FEM-FEM and FEM-FDM explicit coupling approaches (ExDCA) based on the central difference time-marching method [24]. A sketch for the three-dimensional column is depicted in Figure 1(d). The geometrical dimensions of the column are: 10 m  ×  10 m  ×  50 m. Two media (of equal length and cross-section) compose the column; the physical properties of each medium are (null Poisson rate is adopted for elastodynamic subdomains) (i) Medium1—ρ=1.0kg/m3 and c(d)=10m/s; (ii) Medium2—ρ=1.0kg/m3 and c(d)=5m/s.Three different numerical models are considered to simulate this problem, taking into account different coupling procedure combinations. A sketch of the three models adopted is presented in Figure1(e). Details about each numerical model are given by (i) Model 1: elastodynamic FEM formulations are employed. Two independent FEM meshes are adopted, the first one with 2600 linear hexahedral elements and the other one with 2500 linear hexahedral elements (100 elements, i.e., one “element layer”, are used for mesh superposition, as described in [24]); (ii) Model 2: elastodynamic and acoustic FEM formulations are employed, as well as acoustic FDM formulation. 2500 linear hexahedral elements are adopted for the FEM elastodynamic mesh, and 1000 linear hexahedral elements are adopted for the FEM acoustic mesh. 2178 grid points are employed by the space fourth-order FDM discretization (grid points for mesh superposition included); (iii) Model 3: acoustic FEM and FDM formulations are employed. Two independent FEM meshes are adopted, each one with 2000 linear hexahedral elements. 1815 grid points are employed by the space fourth-order FDM discretization (grid points for mesh superposition included).Two numerical analyses are considered, namely: (i) homogeneous analysis, where the entire column is considered composed byMedium1; (ii) heterogeneous analysis, where half of the column is considered composed by Medium1, and the other half by Medium 2. The results achieved for the three different numerical models described above are depicted in Figure 6. The heterogeneous analysis considers two different time steps, namely, Δt1=0.05s (Medium 1) and Δt2=0.10s (Medium 2). For the homogeneous analysis, Δt=0.05s is adopted for the entire domain. In Figure 6, reference results are also depicted; these results correspond to a standard FEM simulation with 5000 linear hexahedral elements and Δt=0.05s (homogeneous and heterogeneous analyses). As one can see, results for all simulations are in good agreement.Time-history results (displacement/pressure x time) at the interface of the three-dimensional column taking into account FEM-FEM and FEM-FDM coupling procedures (ExDCA): (a) homogeneous analysis; (b) heterogeneous analysis.
(a)(b)Considering this first example, the advantages of the discussed multidomain decomposition procedures may be highlighted under several aspects: different time steps are easily adopted for each subdomain and, as a consequence, the algorithm becomes quite robust even when considering media with high properties contrast; moreover, less systems of equations need to be solved along the time-marching process; not all subdomains need to be considered at initial time steps, the activation/initialisation of different subdomains may be controlled based on the properties of the model (wave propagation velocities, etc.), saving most of the computational effort of the first time steps, and so forth.


### 5.2. Dam-Reservoir System

In this second example, a dam-reservoir system, as depicted in Figure7, is analyzed. The structure is subjected to a sinusoidal distributed vertical load, acting on its crest with an angular frequency of ω=18 rad/s.  The material properties of the dam are ν=0.25; E=3.437·109 N/m2; ρ=2000 kg/m3. The adjacent water is characterized by a mass density ρ=1000 kg/m3 and a wave velocity c=1436 m/s. The model is analyzed considering water levels defined by H=50 m and H=35 m.Figure 7
Sketch of the dam-reservoir system.Several ImICA and ExDCA are employed to analyze the dam-reservoir system. Taking into account the ImICA, the following discretizations are considered: (i) FEM-BEM: in this case, 93 quadrilateral finite elements are employed to discretize the dam, and the fluid is discretized by constant-length boundary elements (ℓ=5m). The time-steps adopted are Δtf=0.00350s and Δts=0.00175s [14]; (ii) MLPG-MLPG: 113 nodes are employed to discretize the dam, and the fluid is discretized by a regular equally spaced (horizontally sufficiently extended) distribution of nodes. The time-step adopted for the analyses is Δt=0.002s [27]; (iii) BEM-BEM 1: 34 linear boundary elements of equal length, and 102 linear triangular cells are employed to model the dam, and the fluid is discretized by constant-length boundary elements (ℓ=5m). The time-steps adopted are Δtf=0.003s and Δts=0.001s [23]; (iv) BEM-BEM 2: same as before, without the domain triangular cell mesh. The time-step adopted for the analyses is Δt=0.003s [25]. Results for these ImICA are depicted in Figure 8, taking into account H=50 m and H=35 m.Time-history results for the dam-reservoir system considering the ImICA: (a)H=35 m; (b) H=50 m.
(a)(b)Taking into account the ExDCA, the following discretizations are considered: (i) FEM-BEM—same as in the ImICA. The time-steps adopted areΔtf=0.00350s and Δts=0.000875s [15]; (ii) FEM-FEM—77 quadrilateral elements are employed to discretize the dam, and the fluid is discretized by a regular (horizontally sufficiently extended) distribution of square elements. The time-steps adopted for the analyses are Δtf=0.0014s and Δts=0.0007s [28]. Results for these ExDCA are depicted in Figure 9, taking into account H=50 m and H=35 m.Time-history results for the dam-reservoir system considering the ExDCA: (a)H=35 m; (b) H=50 m.
(a)(b)In this example, the advantages of employing different discretization procedures to analyze different subdomains of the global model can be explored. For instance, for the semi-infinite fluid domain, the BEM can be regarded as an appropriate discretization technique (infinite domain analysis), whereas domain-discretization methods can be applied to model the dam and consider some eventual more complicate behavior. In Figure10, results are depicted (H=50 m) considering linear and nonlinear behavior (elastoplastic analysis, von Mises yield criterion) for the dam and an FEM-BEM discretization [17] for the model (ImICA).Time-history results considering linear and nonlinear material behavior and FEM-BEM implicit iterative coupling analyses.
(a)(b)The results presented so far are obtained taking into account a closed-domain dam, (null displacements are prescribed at the base of the dam and null fluxes are prescribed at the base of the storage lake). As is well known, boundary element formulations are an extremely elegant tool to model infinite media. As a consequence, in the present BEM-BEM 2 coupling context, analyses considering an opened-domain dam (acoustic-dynamic coupling also being carried out at the base of the storage lake) can be provided very easily. For the opened-domain dam case, time-history results are depicted in Figure11 (H=50 m), considering the BEM-BEM 2 [25]. In Figure 12, some time snap shots are depicted, describing the displacement evolution of the closed/opened-domain dam.Time-history results for the opened-domain dam-reservoir system considering a BEM-BEM implicit iterative coupling analysis.
(a)(b)Scaled displacement results for the dam (H=50 m) along time: (a) closed-domain dam; (b) opened-domain dam.
(a)(b)


### 5.3. Tube of Steel Submerged in Water

In this application, two analyses of a tube of steel submerged in water (axisymmetric models) are carried out. A sketch of the first model is depicted in Figure13(a) [16]. A punctual source is located at the centre of the tube (axisymmetric axis), and it emits a signal characterized by three time-sinusoidal cycles (s(t)=sin(wt), where w=10kHz). The properties of the media are (i) water: ρ=1000kg/m3 and c=1500m/s; (ii) steel: ρ=7700kg/m3, E=2.1·1011N/m2, and ν=0.3. In an FEM-BEM ImICA, the boundary of the tube (water cavity) is discretized by acoustic linear boundary elements with length ℓ=0.02m. The tube itself is discretized by elastodynamic linear-square finite elements. The time steps selected are Δtf=8.0·10-6s and Δts=2.0·10-6s. In an FEM-FDM ExDCA, analogous discretization is adopted for the tube, and part of the fluid is discretized by a sufficiently extended FDM mesh. Time-history hydrodynamic pressures at points A and B and displacements at point C (see Figure 13(a)) are depicted in Figure 14, for the FEM-BEM ImICA and for the FEM-FDM ExDCA. As one may observe, the time responses of these two different methodologies are quite similar.Sketch of the tube submerged in water: (a) first case of analysis; (b) second case of analysis; (c) detail of the FEM mesh adopted to model the neighbourhood of the spherical source.
(a)(b)(c)Time-history results for the tube of steel submerged in water (first case of analysis): (a) pressures at point A; (b) pressures at point B; (c) horizontal displacements at point C.
(a)(b)(c)A sketch of the second model is depicted in Figure13(b) [24]. In this case, most of the domain is modelled by the FDM acoustic formulation (water). The metallic tube (marine riser) is modelled by the FEM elastodynamic formulation. A thin water layer surrounding the tube is also modelled by the FEM (acoustic formulation). Two different modelling procedures are adopted to simulate the source: (i) the source is considered punctual, and an excitation term is introduced in the correspondent grid point of the FDM mesh; (ii) the source is considered spherical (radius 0.03429 m), and an FEM mesh is introduced to properly model its neighbourhood (this mesh is depicted in Figure 13(c)). Results obtained from a laboratory experiment [53], as illustrated in Figure 16, are used to validate the numerical response. The marine riser is characterized by ϕ=410mm (external diameter), t=12mm (thickness), and h=4.7m (height). The source produces a sinusoidal excitation with frequency of 20kHz and duration of 3.0·10-4s. The time steps adopted for each subdomain of the numerical model are Δt1=2.0·10-7s (FDM mesh and FEM spherical source mesh); Δt2=1.0·10-7s (FEM mesh surrounding the tube); Δt3=0.5·10-7s (FEM tube mesh).The results achieved for the hydrodynamic pressure at the receiver (hydrophone) are depicted in Figures15(a) and 15(b), for the punctual and spherical source cases, respectively. The results obtained by the experimental analysis are depicted in Figure 15(c). As one may observe, good agreement between experimental and numerical (spherical source case) simulations is obtained (the scale of the graphics should be ignored, since the source intensities adopted in each analysis are different). Comparing the results depicted in Figures 15(a) and 15(b), one can clearly observe the energy dissipation in the source-receptor direction, due to the scattering induced by the spherical source. Figure 17 depicts three snap shots (FDM mesh) of the numerical analysis (punctual source case) and shows some interesting and important features related to the present wave propagation configuration as, for instance, wave fronts (head waves) arising from the faster propagation through the tube (Figure 17(b)) generate a reinforcement of amplitude at the wave front region close to the tube (Figure 17(c)). If one interprets the phenomenon thinking on rays (ray tracing theory), one may be led to erroneously interpret this reinforcement of the amplitude at oblique incidence.Numerical results (pressure x time) at the receiver (second case of analysis) considering FEM-FDM explicit direct coupling procedures and (a) punctual and (b) spherical sources. (c) Experimental results at the receiver.
(a)(b)(c)Figure 16
Photos of the experiment: (a) tube installation through the water tank input gate; (b) tank facilities; (c) tube inside the tank (view through the gate); (d) acoustic transductor ITC 1032 (source/receptor).Pressure distribution for the punctual source case (FDM mesh) at three different moments: (a) begin of propagation; (b) wave fronts due to the faster propagation through the tube wall (head waves); (c) reinforcement of amplitude.
(a)(b)(c)


## 5.1. Multidomain Column

The first example is that of a prismatic body behaving like a one-dimensional column. Initially, the column is analysed as an acoustic model [50]. It is fixed at one end (p̅(t)=0) and subjected to a unitary Heaviside type forcing function acting at the opposite end (q̅(t)=H(t)). A sketch of the model is shown in Figure 1(a). The material properties of the column are c=1.0m/s and ρ=1.0kg/m3. The geometry of the model is defined by L=1.0m. As depicted in Figure 1(a), 28 boundary elements of equal length and 40 quadrilateral finite elements are employed in the coupled mesh. Regarding time discretization, three different cases of analysis are considered here, namely, (i) ΔtF=1.0ΔtB, (ii) ΔtF=0.2ΔtB, and (iii) ΔtF=0.1ΔtB; where ΔtB=0.06s, and the subscripts F and B are related to the FEM and to the BEM, respectively.Figure 1
Column model: (a) FEM-BEM acoustic-acoustic two-dimensional model; (b) MLPG-MLPG and FEM-FEM acoustic-dynamic two-dimensional model; (c) FEM-BEM acoustic-acoustic axisymmetric model; (d) FEM-FEM and FEM-FDM acoustic-acoustic/acoustic-dynamic/dynamic-dynamic three-dimensional model; (e) sketch of the three-dimensional model spatial discretization.In Figure2, time history results are depicted, at points A and B (see Figure 1(a)), taking into account the ExDCA and the ImICA. Potential (pressure) and flux results are presented considering the three different cases of analysis, and they are compared to the analytical solution [51], plotted as a solid line. As can be seen, a higher level of accuracy is observed when different time steps are considered within each subdomain, regarding their optimal temporal discretization. The robustness of the multilevel time-step algorithm must be highlighted: as illustrated in the present application, the algorithm deals properly with highly different subdomain temporal discretizations.Time-history results at points A and B taking into account FEM-BEM coupling procedures and different temporal discretizations for each subdomain: (a) explicit direct coupling analysis; (b) implicit iterative coupling analysis.
(a)
ExDCA(b)
ImICAIn a second approach for the column model, the acoustic-dynamic coupled problem is focused (fluid-solid column [27, 28]). A sketch of the problem is depicted in Figure 1(b). The geometry of the model is defined once again by L=1.0m, and the column is submitted to a time Heavisite force acting at one of its ends. The physical properties of the media are (i) fluid subdomain: κ=100N/m2 (bulk modulus) and ρ=1kg/m3 (density); (ii) solid subdomain: E=100N/m2 (Young modulus), ν=0 (Poisson rate), and ρ=1kg/m3 (density).Two spatial-temporal MLPG discretizations are considered to analyse the model, namely: (i) discretization 1—153 nodes are employed to spatially discretize each subdomain, and the time step adopted isΔt=0.0025s; (ii) discretization 2—561 nodes are employed to spatially discretize each subdomain, and the time-step adopted is Δt=0.00125s. In Figure 3, displacement time-history results at points A and B of the solid subdomain and hydrodynamic pressure time-history results at points B and C of the fluid subdomain are plotted, considering discretizations 1 and 2 and the MLPG-MLPG ImICA. Analytical time histories are also depicted in Figure 3, highlighting the good accuracy of the numerical results.Time-history results at points A, B, and C taking into account MLPG-MLPG coupling procedures (ImICA) and different refinement levels: (a) discretization 1; (b) discretization 2.
(a)(b)The same fluid-solid column is analysed considering FEM-FEM coupled procedures based on the Green-Newmark method and on the Newmark method (ExDCA). In this case, 200 square finite elements are employed to discretize each subdomain of the model, and the time discretization is specified byΔt=0.002s. The obtained results are depicted in Figure 4(a). In this analysis, the same wave propagation velocities are considered within the solid and fluid subdomains, that is, cd=c=10m/s. Taking into account a steel-water column, the physical properties of the media are (i) fluid subdomain: κ=2.3175·109N/m2, ρ=1030kg/m3; (ii) solid subdomain: E=2.1·1011N/m2, ν=0.3, and ρ=7700kg/m3; in this case, cd=6059m/s and c=1500m/s. In Figure 4(b), the displacement and hydrodynamic pressure time-history results at point B are depicted, considering Δtf=Δts and Δtf=4Δts, where Δts=3.33·10-6s, and the subscripts s and f are related to the solid and fluid subdomains, respectively. As one can observe, good agreement between the results is observed, in spite of the different time discretizations adopted within each subdomain.Time-history results for the solid-fluid column at points A, B, and C taking into account FEM-FEM coupling procedures (ExDCA) and different physical models: (a) homogeneous wave propagation velocities; (b) heterogeneous wave propagation velocities.
(a)(b)In a third approach for the column model, the propagation of acoustic waves through a prismatic circular column is analysed (axisymmetric model [16]). A sketch of the problem is depicted in Figure 1(c). The downside half of the column is modelled by the BEM, and the upside half is modelled by the FEM. The properties of the media are (i) Medium1 (BEM)—ρ=1.0kg/m3 and c=1.0m/s; (ii) Medium2 (FEM)—ρ=1.0kg/m3 and c=2.0m/s. The spatial discretization is considered as follows: 32 linear boundary elements discretize Medium1, and 64 linear-square finite elements discretize Medium2. The time steps selected are ΔtB=0.60s and ΔtF=0.15s. Time history acoustic pressures at points A and B (see Figure 1(c)) are depicted in Figure 5. The numerical results are compared with those of a finite difference solution (based on the work presented by Cohen and Joly [52]) and, as one may observe, the time responses are quite similar considering these two different methodologies.Time-history results for acoustic pressures taking into account FEM-BEM coupling procedures (ImICA) applied to a heterogeneous axisymmetric model: (a) results at point A; (b) results at point B.
(a)(b)In a fourth approach for the column model, three-dimensional analyses are considered, taking into account FEM-FEM and FEM-FDM explicit coupling approaches (ExDCA) based on the central difference time-marching method [24]. A sketch for the three-dimensional column is depicted in Figure 1(d). The geometrical dimensions of the column are: 10 m  ×  10 m  ×  50 m. Two media (of equal length and cross-section) compose the column; the physical properties of each medium are (null Poisson rate is adopted for elastodynamic subdomains) (i) Medium1—ρ=1.0kg/m3 and c(d)=10m/s; (ii) Medium2—ρ=1.0kg/m3 and c(d)=5m/s.Three different numerical models are considered to simulate this problem, taking into account different coupling procedure combinations. A sketch of the three models adopted is presented in Figure1(e). Details about each numerical model are given by (i) Model 1: elastodynamic FEM formulations are employed. Two independent FEM meshes are adopted, the first one with 2600 linear hexahedral elements and the other one with 2500 linear hexahedral elements (100 elements, i.e., one “element layer”, are used for mesh superposition, as described in [24]); (ii) Model 2: elastodynamic and acoustic FEM formulations are employed, as well as acoustic FDM formulation. 2500 linear hexahedral elements are adopted for the FEM elastodynamic mesh, and 1000 linear hexahedral elements are adopted for the FEM acoustic mesh. 2178 grid points are employed by the space fourth-order FDM discretization (grid points for mesh superposition included); (iii) Model 3: acoustic FEM and FDM formulations are employed. Two independent FEM meshes are adopted, each one with 2000 linear hexahedral elements. 1815 grid points are employed by the space fourth-order FDM discretization (grid points for mesh superposition included).Two numerical analyses are considered, namely: (i) homogeneous analysis, where the entire column is considered composed byMedium1; (ii) heterogeneous analysis, where half of the column is considered composed by Medium1, and the other half by Medium 2. The results achieved for the three different numerical models described above are depicted in Figure 6. The heterogeneous analysis considers two different time steps, namely, Δt1=0.05s (Medium 1) and Δt2=0.10s (Medium 2). For the homogeneous analysis, Δt=0.05s is adopted for the entire domain. In Figure 6, reference results are also depicted; these results correspond to a standard FEM simulation with 5000 linear hexahedral elements and Δt=0.05s (homogeneous and heterogeneous analyses). As one can see, results for all simulations are in good agreement.Time-history results (displacement/pressure x time) at the interface of the three-dimensional column taking into account FEM-FEM and FEM-FDM coupling procedures (ExDCA): (a) homogeneous analysis; (b) heterogeneous analysis.
(a)(b)Considering this first example, the advantages of the discussed multidomain decomposition procedures may be highlighted under several aspects: different time steps are easily adopted for each subdomain and, as a consequence, the algorithm becomes quite robust even when considering media with high properties contrast; moreover, less systems of equations need to be solved along the time-marching process; not all subdomains need to be considered at initial time steps, the activation/initialisation of different subdomains may be controlled based on the properties of the model (wave propagation velocities, etc.), saving most of the computational effort of the first time steps, and so forth.


## 5.2. Dam-Reservoir System

In this second example, a dam-reservoir system, as depicted in Figure7, is analyzed. The structure is subjected to a sinusoidal distributed vertical load, acting on its crest with an angular frequency of ω=18 rad/s.  The material properties of the dam are ν=0.25; E=3.437·109 N/m2; ρ=2000 kg/m3. The adjacent water is characterized by a mass density ρ=1000 kg/m3 and a wave velocity c=1436 m/s. The model is analyzed considering water levels defined by H=50 m and H=35 m.Figure 7
Sketch of the dam-reservoir system.Several ImICA and ExDCA are employed to analyze the dam-reservoir system. Taking into account the ImICA, the following discretizations are considered: (i) FEM-BEM: in this case, 93 quadrilateral finite elements are employed to discretize the dam, and the fluid is discretized by constant-length boundary elements (ℓ=5m). The time-steps adopted are Δtf=0.00350s and Δts=0.00175s [14]; (ii) MLPG-MLPG: 113 nodes are employed to discretize the dam, and the fluid is discretized by a regular equally spaced (horizontally sufficiently extended) distribution of nodes. The time-step adopted for the analyses is Δt=0.002s [27]; (iii) BEM-BEM 1: 34 linear boundary elements of equal length, and 102 linear triangular cells are employed to model the dam, and the fluid is discretized by constant-length boundary elements (ℓ=5m). The time-steps adopted are Δtf=0.003s and Δts=0.001s [23]; (iv) BEM-BEM 2: same as before, without the domain triangular cell mesh. The time-step adopted for the analyses is Δt=0.003s [25]. Results for these ImICA are depicted in Figure 8, taking into account H=50 m and H=35 m.Time-history results for the dam-reservoir system considering the ImICA: (a)H=35 m; (b) H=50 m.
(a)(b)Taking into account the ExDCA, the following discretizations are considered: (i) FEM-BEM—same as in the ImICA. The time-steps adopted areΔtf=0.00350s and Δts=0.000875s [15]; (ii) FEM-FEM—77 quadrilateral elements are employed to discretize the dam, and the fluid is discretized by a regular (horizontally sufficiently extended) distribution of square elements. The time-steps adopted for the analyses are Δtf=0.0014s and Δts=0.0007s [28]. Results for these ExDCA are depicted in Figure 9, taking into account H=50 m and H=35 m.Time-history results for the dam-reservoir system considering the ExDCA: (a)H=35 m; (b) H=50 m.
(a)(b)In this example, the advantages of employing different discretization procedures to analyze different subdomains of the global model can be explored. For instance, for the semi-infinite fluid domain, the BEM can be regarded as an appropriate discretization technique (infinite domain analysis), whereas domain-discretization methods can be applied to model the dam and consider some eventual more complicate behavior. In Figure10, results are depicted (H=50 m) considering linear and nonlinear behavior (elastoplastic analysis, von Mises yield criterion) for the dam and an FEM-BEM discretization [17] for the model (ImICA).Time-history results considering linear and nonlinear material behavior and FEM-BEM implicit iterative coupling analyses.
(a)(b)The results presented so far are obtained taking into account a closed-domain dam, (null displacements are prescribed at the base of the dam and null fluxes are prescribed at the base of the storage lake). As is well known, boundary element formulations are an extremely elegant tool to model infinite media. As a consequence, in the present BEM-BEM 2 coupling context, analyses considering an opened-domain dam (acoustic-dynamic coupling also being carried out at the base of the storage lake) can be provided very easily. For the opened-domain dam case, time-history results are depicted in Figure11 (H=50 m), considering the BEM-BEM 2 [25]. In Figure 12, some time snap shots are depicted, describing the displacement evolution of the closed/opened-domain dam.Time-history results for the opened-domain dam-reservoir system considering a BEM-BEM implicit iterative coupling analysis.
(a)(b)Scaled displacement results for the dam (H=50 m) along time: (a) closed-domain dam; (b) opened-domain dam.
(a)(b)


## 5.3. Tube of Steel Submerged in Water

In this application, two analyses of a tube of steel submerged in water (axisymmetric models) are carried out. A sketch of the first model is depicted in Figure13(a) [16]. A punctual source is located at the centre of the tube (axisymmetric axis), and it emits a signal characterized by three time-sinusoidal cycles (s(t)=sin(wt), where w=10kHz). The properties of the media are (i) water: ρ=1000kg/m3 and c=1500m/s; (ii) steel: ρ=7700kg/m3, E=2.1·1011N/m2, and ν=0.3. In an FEM-BEM ImICA, the boundary of the tube (water cavity) is discretized by acoustic linear boundary elements with length ℓ=0.02m. The tube itself is discretized by elastodynamic linear-square finite elements. The time steps selected are Δtf=8.0·10-6s and Δts=2.0·10-6s. In an FEM-FDM ExDCA, analogous discretization is adopted for the tube, and part of the fluid is discretized by a sufficiently extended FDM mesh. Time-history hydrodynamic pressures at points A and B and displacements at point C (see Figure 13(a)) are depicted in Figure 14, for the FEM-BEM ImICA and for the FEM-FDM ExDCA. As one may observe, the time responses of these two different methodologies are quite similar.Sketch of the tube submerged in water: (a) first case of analysis; (b) second case of analysis; (c) detail of the FEM mesh adopted to model the neighbourhood of the spherical source.
(a)(b)(c)Time-history results for the tube of steel submerged in water (first case of analysis): (a) pressures at point A; (b) pressures at point B; (c) horizontal displacements at point C.
(a)(b)(c)A sketch of the second model is depicted in Figure13(b) [24]. In this case, most of the domain is modelled by the FDM acoustic formulation (water). The metallic tube (marine riser) is modelled by the FEM elastodynamic formulation. A thin water layer surrounding the tube is also modelled by the FEM (acoustic formulation). Two different modelling procedures are adopted to simulate the source: (i) the source is considered punctual, and an excitation term is introduced in the correspondent grid point of the FDM mesh; (ii) the source is considered spherical (radius 0.03429 m), and an FEM mesh is introduced to properly model its neighbourhood (this mesh is depicted in Figure 13(c)). Results obtained from a laboratory experiment [53], as illustrated in Figure 16, are used to validate the numerical response. The marine riser is characterized by ϕ=410mm (external diameter), t=12mm (thickness), and h=4.7m (height). The source produces a sinusoidal excitation with frequency of 20kHz and duration of 3.0·10-4s. The time steps adopted for each subdomain of the numerical model are Δt1=2.0·10-7s (FDM mesh and FEM spherical source mesh); Δt2=1.0·10-7s (FEM mesh surrounding the tube); Δt3=0.5·10-7s (FEM tube mesh).The results achieved for the hydrodynamic pressure at the receiver (hydrophone) are depicted in Figures15(a) and 15(b), for the punctual and spherical source cases, respectively. The results obtained by the experimental analysis are depicted in Figure 15(c). As one may observe, good agreement between experimental and numerical (spherical source case) simulations is obtained (the scale of the graphics should be ignored, since the source intensities adopted in each analysis are different). Comparing the results depicted in Figures 15(a) and 15(b), one can clearly observe the energy dissipation in the source-receptor direction, due to the scattering induced by the spherical source. Figure 17 depicts three snap shots (FDM mesh) of the numerical analysis (punctual source case) and shows some interesting and important features related to the present wave propagation configuration as, for instance, wave fronts (head waves) arising from the faster propagation through the tube (Figure 17(b)) generate a reinforcement of amplitude at the wave front region close to the tube (Figure 17(c)). If one interprets the phenomenon thinking on rays (ray tracing theory), one may be led to erroneously interpret this reinforcement of the amplitude at oblique incidence.Numerical results (pressure x time) at the receiver (second case of analysis) considering FEM-FDM explicit direct coupling procedures and (a) punctual and (b) spherical sources. (c) Experimental results at the receiver.
(a)(b)(c)Figure 16
Photos of the experiment: (a) tube installation through the water tank input gate; (b) tank facilities; (c) tube inside the tank (view through the gate); (d) acoustic transductor ITC 1032 (source/receptor).Pressure distribution for the punctual source case (FDM mesh) at three different moments: (a) begin of propagation; (b) wave fronts due to the faster propagation through the tube wall (head waves); (c) reinforcement of amplitude.
(a)(b)(c)


## 6. Conclusions

The present paper discusses multidomain decomposition techniques to model the propagation of interacting acoustic-elastic waves considering several coupling procedures. Two basic algorithms are presented here, namely, the ExDCA (explicit direct coupling approach) and the ImICA (implicit iterative coupling approach), which are based on explicit and implicit time-marching schemes, respectively, and multidomain decomposition coupling procedures. Within the context of these two basic algorithms, several coupled numerical methods are presented along the paper, such as FEM-FEM, FEM-FDM, FEM-BEM, BEM-BEM, DBEM-BEM, and MLPG-MLPG. Independent temporal and spatial (i.e., no matching nodes in common interfaces) discretizations within interacting subdomains are also discussed in the paper, being several applications of the discussed multilevel time-step algorithm presented along Section5, illustrating its good performance and potentialities. As a matter of fact, in Section 5, several numerical applications are considered (e.g., acoustic-acoustic/acoustic-dynamic/dynamic-dynamic wave propagation problems, two-dimensional/three-dimensional/axisymmetric models, different coupled numerical procedures, etc.), illustrating as a whole the good flexibility, accuracy, stability, and robustness of the discussed methodologies.



---

*Source: 245170-2011-09-18.xml*