Finite Element Analysis of Dam-Reservoir Interaction Using High-Order Doubly Asymptotic Open Boundary

The dam-reservoir system is divided into the near field modeled by the finite element method, and the far field modeled by the excellent high-order doubly asymptotic open boundary DAOB . Direct and partitioned coupled methods are developed for the analysis of dam-reservoir system. In the direct coupled method, a symmetric monolithic governing equation is formulated by incorporating the DAOB with the finite element equation and solved using the standard timeintegration methods. In contrast, the near-field finite element equation and the far-field DAOB condition are separately solved in the partitioned coupled methodm, and coupling is achieved by applying the interaction force on the truncated boundary. To improve its numerical stability and accuracy, an iteration strategy is employed to obtain the solution of each step. Both coupled methods are implemented on the open-source finite element code OpenSees. Numerical examples are employed to demonstrate the performance of these two proposed methods.


Introduction
The coupled analysis of dam-reservoir interaction has great significance for the design and safety evaluation of concrete dams under earthquakes.The finite element method and substructure method are often applied for the analysis of dam-reservoir system Figure 1 .The dam structure as well as the near-field reservoir with irregular geometry is discretized with finite elements.The remaining part of the reservoir, called the far field, is simplified as a semi-infinite layer with constant depth.On the truncated boundary, which separates the near and far field, the equations of motion and continuity should be satisfied simultaneously.In the early studies, frequency-domain analysis methods 1, 2 are often used for linear problems.However, for a nonlinear analysis of the dam-reservoir system, it is necessary to develop a direct time domain analysis.Zienkiewicz  radiation condition 5 .Tsai and his coworkers 6-8 established the time-domain models for the dynamic interaction analysis of dam-reservoir system by using the transmitting boundary.This approach is temporally global; that is, it requires the evaluation of convolution integrals.Boundary element method is undoubtedly a powerful numerical tool for analysis of problems involving unbounded domain.When the boundary element method 9-13 is applied to the direct time-domain analysis of dam-reservoir interaction, the formulation is temporally and spatially global.As a result, great numerical effort is required to solve the transient problems.This hinders its application to long-time analysis of large-scale engineering problems.The scaled boundary finite element method is a semianalytical method which is efficient for the analysis of dam-reservoir interaction.Details of this approach are given, for example, by Li et al. 14 and Lin et al. 15 .
The high-order open boundaries are promising alternatives for the simulation of semiinfinite reservoir in the analysis of dam-reservoir system.The high-order open boundaries 16, 17 are of increasing accuracy as the order of approximation increases.Moreover, the formulations are temporally local so that they are computationally efficient.As demonstrated by Prempramote et al. 18 , these high-order open boundaries are singly asymptotic at the high-frequency limit and are appropriate for the radiative fields, where virtually all of the field energy is propagating out to infinity 19 .However, in some classes of application, such as a semi-infinite reservoir with constant depth also known as a wave guide , a cutoff frequency exists.When the excitation frequency is close to or below the cutoff frequency, the wave field is largely nonradiative.In such cases, the high-order transmitting boundaries break down at late times in a time domain analysis 18 .To model an unbounded domain with the presence of nonradiative wave fields, one advance is the introduction of the doubly asymptotic boundaries 19-22 .Thus, the dynamic stiffness is exact at both the highfrequency and the low-frequency limit i.e., statics , with its formulation spatially global.However, the highest order denoting the accuracy reported in the literature is three 23 .
Recently, a novel high-order doubly asymptotic open boundary for one-dimensional scalar wave equation is proposed by Prempramote et al. 18 by extending the work in 24 .This high-order doubly asymptotic open boundary is capable of accurately mimicking the unbounded domain over the entire frequency i.e., from zero to infinity .This open boundary condition is constructed by using the continued fraction solution of dynamic stiffness matrix without explicitly evaluating its solution at discrete frequencies.When applied for a semiinfinite layer with a constant depth, the constants of the continued fraction solutions with any order are determined explicitly and recursively.Excellent accuracy and stability for longtime transient analysis are reported.In this paper, two coupled numerical methods for dam-reservoir interaction analysis based on the excellent high-order doubly asymptotic open boundary 25 are proposed.In the direct coupled method, the high-order doubly asymptotic open boundary is directly incorporated with the near-field finite element equation.A monolithic governing equation for the whole dam-reservoir system is formulated with sparse and symmetric coefficients matrices, which can be solved using the standard time-integration methods.In the partitioned coupled method, the near-field finite element equation and the far-field high-order doubly asymptotic boundary condition are separately solved.The high-order doubly asymptotic open boundary is split into two parts.The first part is proved to be the Sommerfeld radiation boundary, which can be included in the damping matrix of the near-field finite element equation.The second part includes all the high-order terms and is governed by a system of firstorder ordinary differential equations.These two sets of equations are solved by a sequential staggered implicit-implicit partition algorithm.To improve the stability and accuracy of this partitioned coupled method, an iteration strategy is employed to obtain the solution of each step.Both of these two coupled methods are numerically implemented on the open-source finite element code OpenSees to analyze the gravity dam-reservoir and arch dam-reservoir interaction.
This paper is organized as follows.In Section 2, the finite element formulation of dam-reservoir system is addressed.In Section 3, the scaled boundary finite element method is applied to 3-dimensional semi-infinite layer with constant depth.The governing equation on the truncated boundary is obtained.In Section 4, the scaled boundary finite element equation is decoupled by modal decomposition.Based on the continued fraction solutions of the dynamic stiffness, a high-order doubly asymptotic open boundary is constructed by introducing auxiliary variables.In Section 5, two coupled numerical methods are presented: the direct coupled method and the partitioned coupled method.Both numerical methods are implemented on the open source finite element code OpenSees.In Section 6, numerical examples of a gravity dam and an arch dam are presented.In the final section, conclusions are stated.

Modeling of Dam-Reservoir System
A typical dam-reservoir system is shown in Figure 1.The reservoir is divided into two parts: the near field with irregular geometry and the far field extending to the infinity with constant depth.The dam structure and near field reservoir are dicretized with finite elements, and the far field reservoir is modeled with high-order doubly asymptotic boundary.
The hydrodynamic pressure is denoted as p p x, y, z, t , and the acceleration of water particle is denoted as { ü} { ü x, y, z, t }.Assuming the water in the reservoir to be compressible, inviscid, and irrotational with small amplitude movements, the hydrodynamic pressure p in the reservoir satisfies the scalar wave equation where Δ is the Laplace operator and c is the compression wave velocity where K is the bulk modulus of water and ρ is the mass density.Without considering the material damping, the finite element formulation for damreservoir system can be partitioned as where M , K , and Q are the mass matrix, static stiffness matrix, the coupling matrix between structure and acoustic fluid, respectively, and {f} is the external force vector.
Subscript s denotes the degrees of freedom of the dam structure; subscript b denotes the degrees of freedom on the truncated boundary; subscript f denotes the degrees of freedom of the near-field water except for those on the truncated boundary.The interaction load applied to the semi-infinite reservoir far field by the near-field reservoir is denoted as {r}, so the external load applied to the near-field reservoir on the truncated boundary is equal to −{r}.
The mass and stiffness matrices of water treated as acoustic fluid are expressed as where N denotes the shape function of finite elements.To solve 2.7 for the dam-reservoir system, the relationship between the interaction load {r} and the hydrodynamic pressure {p} of the semi-infinite reservoir is formulated in the following section.

Summary of the Scaled Boundary Finite Element Method for Semi-Infinite Reservoir with Constant Depth
The scaled boundary finite element method is a semianalytical method developed to model unbounded domains with arbitrary shape 27 .The scaled boundary finite element formulation for the two-dimensional semi-infinite reservoir with constant depth was described in detailed 25 .For the sake of completeness, a brief summary of the equations necessary for the interpretation of high-order doubly asymptotic open boundary for hydrodynamic pressure is presented in this section.
To facilitate the coupling with the finite elements of the near-field reservoir, the truncated boundary is discretized by elements that have the same nodes and shape function as the finite elements.The derivation is summarized for three-dimensional semi-infinite reservoir with a vertical boundary Figure 2 .Streamlined expressions are presented as follows.
For an acoustic fluid, the relationship between acceleration and hydrodynamic pressure is equivalent to that between stress and displacement in stress analysis and is expressed as where {L} {∂/∂x ∂/∂y ∂ /∂z} is the differential operator.The equation of continuity without considering the volumetric stress-strain relationship of compressible water is written as

Mathematical Problems in Engineering
The Cartesian coordinates of a point x, y, z and inside the semi-infinite reservoir are expressed as where the coordinate ξ is equal to 0 on the vertical boundary.The Jacobian matrix for coordinate transformation from x, y, z to ξ, η, ζ is expressed as For a three-dimensional problem, where |J| is the determinant of the Jacobian matrix.The partial differential operator defined in 3.1 is expressed as Along horizontal lines passing through a node on the boundary, the nodal hydrodynamic pressure functions {p} {p ξ, t } are introduced.On the boundary, the nodal hydrodynamic pressure follows as {p b } {p ξ 0, t }.One isoparametric element S e is shown in Figure 2

3.11
Applying Galerkin's weighted residual technique in the circumferential directions η, ζ to 3.2 , the scaled boundary finite element equation for the three-dimensional semi-infinite reservoir with constant depth is obtained as 12 where E 0 , E 2 , and M 0 are coefficient matrices

3.13
The coefficient matrices of 3.12 are evaluated by standard numerical integration techniques in the finite element method.Similar to the static stiffness and mass matrices in the finite element method, both of the coefficients E 0 and E 2 are sparse and positive definite.This formulation is the same as that of two-dimensional case 25 .It is applicable for Mathematical Problems in Engineering 3-dimensional case with vertical truncated boundary of arbitrary geometry, such as the arch dam-reservoir system.
The acoustic nodal load vector r {r ξ, t } on a vertical boundary with a constant ξ is obtained based on virtual work principle and is expressed as

3.14
Assuming the time-harmonic behavior {p ξ, t } {P ξ, ω }e iωt and {r ξ, t } {R ξ, ω }e iωt ω is the excitation frequency with the amplitudes of the hydrodynamic pressure {P } {P ξ, ω }, 3.12 is transformed into the frequency domain as and the amplitudes of the acoustic nodal load {R} {R ξ, ω } are expressed as 3.16

High-Order Doubly Asymptotic Open Boundary for Hydrodynamic Pressure
The derivation of high-order doubly asymptotic open boundary for hydrodynamic pressure is implemented based on the modal expression of 3.15 .The streamlined expressions in 25 are summarized in this section.

Modal Decomposition of Scaled Boundary Finite Element Equation
Following the procedure in detailed 25 , the system of ordinary differential equations in 3.15 can be decoupled by a modal transformation.The modes are obtained from the following generalized eigenvalue problem • stands for a diagonal matrix : where λ 2 j is the diagonal matrix of positive eigenvalues, h is a characteristic length e.g., the depth of the semi-infinite layer to nondimensionlize the eigenvalues, and Φ are the matrix of eigenvectors representing the modes, which are normalized as As a result, the inverse of the eigenvector matrix can be obtained by the matrix multiplication The relationship between amplitudes of the hydrodynamic pressure and amplitudes of the modal hydrodynamic pressure { P } { P ξ, ω } is defined as Substituting 4.5 into 3.15 premultiplied with Φ T and using 4.2 and 4.3 lead to a system of decoupled equations with the dimensionless frequency where index j indicates the modal number.Substituting 4.5 into 3.16 , the acoustic nodal force vector is expressed as The amplitude of the modal nodal force vector or R j −h P j,ξ .4.9 Premultiplying 4.8 and using 4.2 and 4.9 yield R h Φ T {R}.

4.10
This equation transforms the amplitude of the acoustic nodal force vector to the amplitude of the modal force vector.The modal dynamic stiffness coefficient S j a 0 is defined as R j S j a 0 P j .

4.11
By eliminating R j and P j from 4.6 , 4.9 , and 4.11 , an equation for the modal dynamic stiffness coefficient is derived as

Doubly Asymptotic Continued Fraction Solution for Modal Dynamic Stiffness
Based on a doubly asymptotic solution of the modal dynamic stiffness coefficient S j a 0 , a temporally local open boundary 18 is constructed for a single mode of wave propagation.The solution of 4.12 is expressed as a doubly asymptotic continued fraction.An order M H high-frequency continued fraction is constructed first as

4.13
It is demonstrated by Prempramote et al. 18 that the high-frequency continued fraction solution does not converge when the excitation frequency is below the cutoff frequency.To determine a valid solution over the whole frequency range, an M L order low-frequency continued fraction solution is constructed for the residual term Y M H 1 j a 0 .Denoting the residual term for mode j as 4.14 The continued fraction solution for Y L,j a 0 at the low frequency limit is expressed as

4.15
The doubly asymptotic continued fraction solution is constructed by combining the high-frequency continued fraction solution in 4.13 with the low-frequency solution in 4.15 using 4.14 .

High-Order Doubly Asymptotic Open Boundary
Following the procedure developed for the modal space 18 , the acoustic force-pressure relationship in the time domain is formulated by using the continued fraction solution of the dynamic stiffness coefficient and introducing auxiliary variables.A system of firstorder ordinary differential equations with symmetric coefficient matrices is obtained, which represents a temporally local open boundary condition.
Substituting the continued fraction solution of the modal dynamic stiffness coefficient 4.13 -4.15 into 4.11 and applying the inverse Fourier transform, the time-domain highorder doubly asymptotic open boundary in the modal space is expressed as For an order M H M L doubly asymptotic continued fraction solution, the residual term p M L 1 L,j 0 applies.Equations 4.16 -4.20 constitute a system of first order ordinary differential equations relating the interaction load { r j }, hydrodynamic pressure { p j } and the modal auxiliary variables p in the modal space.This system of ordinary differential equations is a temporally local high-order doubly asymptotic open boundary condition for the semi-infinite reservoir with constant depth, which is directly established on the nodes of a vertical boundary.This boundary condition can be coupled seamlessly with finite element method.

Coupled Numerical Methods for Dam-Reservoir Interaction Analysis
Based on the high-order doubly asymptotic open boundary, two coupled numerical methods will be presented in this section: the direct coupled method and the partitioned coupled method.Both coupled numerical methods are implemented on the open-source objectoriented finite element code OpenSees, that is, the open system for earthquake engineering simulation.

The Direct Coupled Method
In the direct coupled method, by incorporating the high-order doubly asymptotic open boundary with the near-field finite element equation, a monolithic governing equation for the near-field structure and auxiliary variables representing the far-field reservoir is formulated.
To derive a symmetric monolithic formulation, modal transform is applied to the auxiliary variables defined in modal space { p i } and { p i L } as where {p i } and {p i L } are the auxiliary variables defined in real space.Using 4.5 , 4.10 , 5.1 , and left-multiplying 4.16 -4.20 by with the symmetric and positive definite matrix

5.7
The residual term E 0 { ṗ with the block matrices

14 Mathematical Problems in Engineering
The global load vector f g and unknowns u g are expressed as the semicolon ";" stands for the vertical concatenation of vectors

5.11
This system of linear equations describes the complete dam-reservoir system taking account of the influence of the semi-infinite reservoir.Similar to the near-field finite element formulation 2.7 , the coefficient matrices of this formulation are sparse and symmetric.This dynamic system can be solved using a standard time-integration method, such as the Newmark family schemes.Equation 5.8 is of order where N s is the number of degrees of freedom of dam structure, N f is the number of degrees of freedom of nearfield reservoir except for those of the truncated boundary, and N fb is the number of degrees freedom of truncated boundary.Compared with the near-field finite element equation 2.7 , additional auxiliary M H M L 1 N fb degrees of freedom are introduced in 5.8 .From the point view of numerical implementation, 5.8 can be implemented on existing finite element code without any special treatment and solved by existing finite element solver.The block coefficient matrices E 0 and A in 5.10 are evaluated only once in the analysis.However, the formulation of the direct coupled method involves additional auxiliary variables, which requires more computational effort and memory to solve this dynamic system.

The Partitioned Coupled Method
In the partitioned coupled method, the near-field finite element equation and the far-field high-order doubly asymptotic boundary condition are separately solved.They are coupled by the interaction force on the truncated boundary.The deviation of the partitioned coupled method is detailed in 25 .Streamlined expressions are presented as follows.
Using 4.3 and 4.10 , 4.16 left-multiplied by Φ −T is rewritten as

5.12
Substituting 5.12 into 2.7 , the finite element formulation of the dam-reservoir system considering the interaction between the near-field and the semi-infinite reservoir is expressed as

5.13
Mathematical Problems in Engineering 15 The high-order doubly asymptotic open boundary is split into two parts.The first part is the additional damping term on left-hand side of 5.13 , which is equivalent to the Sommerfeld radiation boundary as demonstrated in 25 .The second part is the coupling term Φ −T λ j { p 1 }/h on the right-hand side of 5.13 representing the contribution of the high-order terms of the doubly asymptotic boundary.It can be regard as an external load acted on the truncated boundary.For efficiency consideration in the numerical implementation, the coupling term Φ −T λ j { p 1 }/h is evaluated in the modal space.
Assembling 4.17 -4.20 multiplied by h/λ j , a system of ordinary differential equation for the modal auxiliary variables is formulated as where the coefficient matrices are expressed as

5.15
The unknown vector {z A,j t } consists of the modal auxiliary variables of mode j the semicolon ";" stands for the vertical concatenation of vectors z A,j t p

Mathematical Problems in Engineering
The only nonzero entry on the right-hand side is the modal hydrodynamic pressure p j obtained from 4.5 f A,j t − p j ; 0; • • • ; 0; 0; 0; • • • ; 0 .

5.17
It is demonstrated by Wang et al. 25 that 5.14 is stable up to the order M H M L 100 at least.Equations 5.13 for the near field and 5.14 for far field are coupled by the auxiliary variables { p 1 } defined in the modal space.Using the same time integration scheme, such as the trapezoidal rule, these two sets of equations are solved by a sequential staggered implicit-implicit partitioned procedure proposed in 28, 29 .The value of the modal auxiliary variables { p 1 } at time station t n 1 is determined by the last-solution extrapolation predictor 28, 29

5.18
The auxiliary variables { p 1 } are obtained by integrating 5.14 for prescribed hydrodynamic pressure {p}.The algorithm proposed in 25 to solve 5.13 and 5.14 proceeds as follows.
Step 2. Extract { p 1 } n from {z A,j } n 0 of each mode and assigned to { p 1 } p n 1 in 5.18 .
Step 3. Form the right-hand term of 5.13 , compute {u} n 1 and {p} n 1 by using an implicit method.
Step 4. Calculate modal hydrodynamic pressure { p} n 1 and form the right-hand term of 5.14 .
Step 5. Compute {z A,j } n 1 for each mode by using an implicit method.
Step 6. Increment n to n 1 and go to Step 2.
Since the predicted vector { p 1 } p n 1 has been used rather than { p 1 } n 1 in solving 5.13 and 5.14 , this algorithm may lead to a numerical instability or poor accuracy.To avoid these, the solution algorithm for the partitioned coupled method in this paper is modified.The solution process within one step given by Step 2 to Step 5 is repeated a number of times in an iterative manner 26 .In each additional cycle, { p 1 } p n 1 in Step 2 is extracted by the last computed {z A,j } n 1 in Step 5. Numerical experiments demonstrate that stability and accuracy are improved by performing a few additional iterations.
The order of 5.14 is only M H M L 1, and little computational effort is required to solve 5.14 .Consequently, additional memory required for the solution of partitioned method is less than that of the direct coupled method.
This partitioned coupled method and the sequential staggered implicit-implicit partition algorithm 25 are both implemented in the open-source finite element code OpenSees rather than ABAQUS; as a result, the computational efficiency is greatly improved without any time costing restart analysis in ABAQUS.

Numerical Examples
Two numerical examples are analyzed to evaluate the accuracy and efficiency of the two present coupled numerical methods.The first example is a gravity dam with an irregular near field reservoir.The open boundary is employed to represent the regular semi-infinite reservoir of a constant depth.The time integration for these two coupled numerical methods is performed by the trapezoidal integration.
It is demonstrated by Wang et al. 25 that an order M H M L 10 high-order doubly asymptotic open boundary condition is of excellent accuracy.So, the order M H M L 10 open boundary condition is chosen for these two numerical examples.

Gravity Dam
A typical gravity dam-reservoir system with an irregular near field is shown in Figure 3  The finite element mesh is shown in Figure 3 b .The system is divided into three parts: the dam body, the near-field reservoir, and the far-field semi-infinite reservoir with constant depth.The dam body is discretized with eight-node solid elements, the near-field reservoir with 156 eight-node acoustic fluid elements, and the dam-reservoir interface with 13 threenode interface elements.The far-field reservoir is modeled by 13  elements on the truncated boundary, which share the same nodes and are compatible with those of the near-field acoustic fluid elements.The total number of nodes of the whole model is 653.
To verify the results, an extended mesh covering a far-field reservoir region of 7200 m is analyzed.The size of extended mesh is sufficiently large to avoid the pollution of the dam response by the waves reflected on the truncated boundary for a time duration of 2 × 7200/1438.7 ≈ 10 s.
The El Centro earthquake ground motion see Figure 4 is imposed as the horizontal acceleration at the base of the dam.The time step is chosen as 0.002 s during which the pressure wave travels about one third of the distance between two adjacent nodes.The partitioned coupled method is performed without any additional iteration; that is, the solution process is the same as that in 25 .The responses of hydrodynamic pressure at heel of the gravity dam of the first 20 s are plotted in Figure 5. Excellent agreement between the solutions of both coupled numerical methods and the extended mesh solution is observed during the first 10 s.To demonstrate the improvement of the modified solution algorithm for partitioned coupled method, the time step is increased to 0.02 s.The responses of the first 20 s are plotted in Figure 6.As it is shown, the results of the solution algorithm proposed by Wang et al. 25 tend to be divergent and inaccurate.However, both the solutions of direct coupled method and partitioned coupled method agree with the solution of extended mesh very well.The number of additional iterations within each step for partitioned coupled method recorded is usually one and no more than four in this example.
The computer time for the gravity dam example list in Table 1 is recorded on a PC with a 2.93 GHz Intel Core i7 CPU.High computational efficiency of both coupled methods is observed.The computer time of the direct coupled method is about one tenth of that of the extended mesh.The computational effort of the partitioned coupled method is directly associated with the number of additional iterations.When there is no additional iteration, the computer time of the partitioned coupled method is nearly equal to that of the direct coupled method.

Arch Dam
An arch dam-reservoir system is shown in Figure 7.The physical properties of dam body and water are the same as that in the example of gravity dam.The arch dam is of a height of 22 m and the near-field reservoir covers a region of the dam height.The dam body is discretized with 272 twenty-node hexahedron solid elements, the near-field reservoir with 1088 twenty- node hexahedron acoustic fluid elements, and the dam-reservoir interface with 136 eightnode interface elements.The far-field reservoir is modeled by 136 eight-node quadratic elements on the truncated boundary, which share the same nodes and are compatible with those of the near-field acoustic fluid elements.The total number of nodes of the whole model is 6652.To verify the results, a similar extended mesh covering the region of 7200 m is analyzed.
Similar to the gravity dam example, the El Centro earthquake ground motion is imposed as the horizontal Y direction acceleration at the base of the arch dam.The time step is chosen as 0.02 s.The responses of the hydrodynamic pressure at the heel of the arch dam are plotted in Figure 8.The solution of the direct coupled method agrees with solution of the extended mesh and is long time stable.However, as for the partitioned coupled method, numerical divergence is observed at the early time of the analysis.Numerical instability of such coupling strategy is also reported in 12 .It can be expected as the partitioned coupled method is conditional stable; that is, the integration time step is limited by stability limits.When the time step is greater than the stability limit, numerical instability may occur, for example, the results of the arch dam example.As the dam-reservoir system is quite complicated, it is difficult to determine the stability limits of different application cases.When the partitioned coupled method is applied, smaller time step should be used.
The computer time of the direct coupled method recorded is 1463 s, which is about one ninth of that of the extended mesh, that is, 12547 s.It also demonstrates the high computational efficiency of the direct coupled method.

Conclusions
Two coupled numerical methods were developed for the dam-reservoir interaction analysis by incorporating the finite element method with the excellent high-order doubly asymptotic open boundary.The dam-reservoir system is divided into the near-field with irregular geometry and the far-field by the truncated boundary.In the direct coupled method, a global monolithic equation for the whole dam-reservoir system is formulated with sparse and symmetric coefficient matrices, which can be solved by the standard finite element solver.In the partitioned coupled method, the near-field finite element equation and the high-order open boundary condition are separately solved.They are coupled by the interaction force applied on the truncated boundary.The partitioned coupled method is achieved by using a sequential staggered implicit-implicit procedure.To improve the numerical stability and accuracy of the algorithm, an iteration strategy was employed to obtain the solution of each step.
Both of the two coupled numerical methods are implemented on the open-source finite element code OpenSees.Numerical experiments demonstrated the high efficiency and accuracy of both coupled numerical methods.The memory required for the solution of the partitioned method is less than that of the direct coupled method.Although the numerical stability and accuracy of the partitioned coupled method can be improved by additional iterations within each step, the partitioned coupled method is conditionally stable yet.Its stability is also related to the predictor.Further research should be carried out to improve the stability of the partitioned coupled method.In contrast, the direct coupled method is unconditionally stable if only an unconditionally stable time integration algorithm such as trapezoidal integration is chosen.Consequently, larger time steps can be used in the direct coupled method than that in the partitioned coupled method.

Figure 3 :
Figure 3: A gravity dam-reservoir system with irregular near field: a Geometry; b Mesh.
a , which is the same as the flexible dam example in 25 .The dam body has a modulus of elasticity E 35 GPa, Posisson's ration v 0.2, and mass density ρ 2400 kg/m 3 .The water in the reservoir has a pressure wave velocity c 1438.7 m/s and mass density ρ 1000 kg/m 3 .

Figure 6 :
Figure 6: Hydrodynamic pressure at heel of the gravity dam under El Centro ground motion with a time step of 0.02 s.
Wang et al. 25 extend the high-order doubly asymptotic open boundary condition by Prempramote at al. 18 to the analysis of the hydrodynamic pressure of a semiinfinite reservoir with constant depth.By applying the sequential staggered implicit-implicit partition algorithm, the high-order doubly asymptotic open boundary is coupled with the general purpose finite element software ABAQUS to analyze the gravity dam-reservoir interaction system.Numerical examples demonstrate the high accuracy and long-time stability of this proposed technique.Givoli et al. 26 proposed a new finite element scheme for the solution of timedependent semi-infinite wave-guide problems by incorporating with a high-order open boundary.Two versions of finite element formulation, namely, the augmented version and the split version, are proposed.Good performance of this scheme is demonstrated in numerical examples, but the global mass and stiffness matrices of the augmented formulation are nonsymmetric.
denotes global.M g , C g , and K g are the global mass matrix, global damping matrix, and global stiffness matrix, which are expressed as

Table 1 :
Computer time for the gravity dam example.