An efficient finite element model for dynamic analysis of gravity dam-reservoir-foundation interaction problems

In this paper, a finite element (FE)-based model that efficiently evaluates the dynamic behavior of dam-reservoir-foundation interaction (DRFI) problem was proposed including the radiation of waves to the unbounded rock and reservoir domains. Lagrangian fluid elements were used to discretize the near-field reservoir domain, while the presented infinite fluid elements were used to discretize the far-field reservoir domain. The fully coupled equation of motion for DRFI problem was solved by direct method. A two-dimensional (2D) plane-strain FE formulation of the problem is written in FORTRAN 90 programming language. Investigations were conducted on the effect of near-field domain size (length and depth) on the dynamic behavior of DRFI, dam-foundation interaction (DFI), and dam-reservoir interaction (DRI) problems. The results of this study demonstrate that the proposed model outperforms many other models that have been evolved in the literature in terms of accuracy and speed. The reflected hydrodynamic pressures at the far-field reservoir domain were efficiently absorbed by the suggested infinite fluid elements. The near-field domains size has a noticeable impact on the dynamic behavior of the dam. Making an exact choice about the size is more challenging. However, it was observed that the size of 1.5H is the physically appropriate response.

In the modeling of the DRI, DRFI problems using Lagrangian aproach, the coupled equation of motion is generally solved by the methods such as Mode superposition, Newmark, Wilson-θ and HHT-α method.The far-field reservoir domain effect as well as foundation effect is not included in the FE modeling of the problem as radiation damping.Furthermore, the material damping mechanism of the system is generally selected as viscous damping (Bayraktar, Akkose et al. 2005a, Calayir and Karaton 2005a, Calayir and Karaton 2005b, Akköse, Adanur et al. 2007, Hacıefendioğlu, Başağa et al. 2007, Akköse, Adanur et al. 2008, Akkose, Bayraktar et al. 2008, Bayraktar, Altunişik et al. 2009, Wang, Zhang et al. 2014).In these solution methods, the dynamic response of the system is generally obtained by contribution of modes.The number of modes in FE modeling of DRFI and DRI problems is seen to be large in amount.Among this large amount of modes, the contribution of important modes on the dynamic responses of the problems is difficult to be estimated.Many tests and experimental analysis is needed to estimate the number of prominent modes.However, if the contribution of all modes is taken into account, it will take a long time to analyze the problem, which is inefficient and time-consuming.Another drawback with these solution technics is their stability, which needs further work.All of these drawbacks could be removed by solving the coupled equation of motion in the frequency domain using the direct method.In this method, the coupled equation of motion is transformed to the frequency domain.There is no need for further tests to find the number of important modes anymore.The radiation of reflected waves could easily be included by utilizing infinite elements in the model which is the most suitable approach in the FE method.As material damping, either hysteretic or Rayleigh damping models can be included in the coupled equation of motion.The whole system can be solved at the same time.This method is also faster than many other methods.But, the transformation of frequency domain responses to the time domain is the only further work in this method.
The dynamic analysis of roller compacted concrete gravity dam has been investigated by (Ghaedi, Hejazi et al. 2018) considering DRFI effects and flexibility of foundation.It is reported that the flexibility of foundation domain has a great effect on the responses.The effect of base-rock characteristics on the dynamic response of DFI problem was investigated by (Bayraktar, Hançer et al. 2005b) considering rigid base foundation input, massless foundation input and deconvolved bedrock foundation models.In this paper, the responses obtained are compared with each other.By comparing the mean responses, it is reported that the massless-foundation input is inadequate model, but it is a practical technique.Similarly, the dynamic analysis of DRFI problem is investigated by (Ghaemian, Noorzad et al. 2019, Salamon, Wood et al. 2021) incorporating rigid, massless, and massed foundation models and it is reported that the massless foundation model overestimates the dynamic response of dam.The seismic responses of DRFI problem are evaluated by (Gorai and Maity 2019) under near-fault and far-fault ground motion data.In the FE modeling of the problem, the width and depth of the foundation domain is assumed to be 3b and 1.5b, respectively, where b is the width of dam structure.The seismic performance of DRFI problem is investigated by (Chen, Yang et al. 2019).In this FE model, an artificial viscous boundary condition is used at the truncated foundation boundaries with 3 times dam height in each direction.The seismic response evaluation of gravity dams including the effects of foundation mass is implemented by (Asghari, Taghipour et al. 2018).In this FE analysis, viscous dampers are utilized at the truncated boundaries of the foundation.The dimension of foundation in the model is varied from 2-10H, where H is the height of dam.In this paper the length of reservoir domain is assumed to be constant which is 2H in all conditions.And it is reported that the distance 5H is relatively appropriate dimension for truncation of boundaries for foundation domain.The effect of reservoir length on the seismic performance of concrete gravity dams under near and far-fault ground motion is investigated by (Bayraktar, Türker et al. 2010).The depth of foundation (H) in the model is assumed to be constant.It is concluded that the length of 3 and 4H is the sufficient length for seismic performance evaluation of dams, where H is the height of dam structure.(Kartal, Cavusli et al. 2017) demonstrated that the 3H is an appropriate length for seismic response evaluation of concrete gravity dams considering concrete nonlinearity.The same result is also reported by (Hariri-Ardebili, Seyed-Kolbadi et al. 2016, Wang, Zhang et al. 2021), while the length of 5H is reported by (Pelecanos, Kontoe et al. 2013).While considering all these research papers, the researchers have investigated the most prominent topics regarding the dynamic behavior of gravity dams using FE method.Many different numerical modelings for dynamic analysis of this problem are proposed either in time or frequency domain.So far, very little attention has been paid to the impact of near-field domains size on the dynamic responses.There is no clear and consistent result reported on the precise dimension of near-field domains either for reservoir or foundation domain and both.
This study presents a novel and effective FE model for the dynamic analysis of DRI, DFI, and DRFI.The main goal of this study is determination of a consistent size regarding the near-field domains.

DESCRIPTION OF THE FE MODEL
The finite-infinite element model of DRI, DFI, and DRFI problems developed in this study is explained here.In this numerical model, the dam structure and near-field rock foundation domains were discretized by 8-noded iso-parametric quadrilateral finite elements.The reservoir domain was divided into near-field and far-field domains.The near-field reservoir domain includes the sloshing behavior of water with its compressibility.The near-field reservoir domain was discretized by 9-noded Lagrangian fluid elements.The free-surface region of water simulating the sloshing behavior was discretized by 1D 3-noded linear elements.These finite element types are well-known in the literature, and therefore, the basic mathematical formulas that lead to the interpolation functions of the elements are not given.The interpolation functions of the finite elements were taken from (Klaus-Jurgen 1982).
The far-field reservoir domain which absorbs the reflected hydrodynamic waves was discretized by 7-noded fluid infinite elements.This infinite element model was taken from (Yerli, Temel et al. 1998) and then, it was modified to use in the unbounded reservoir domain.The foundation of concrete gravity dam and reservoir domain was assumed as bedrock.The far-field foundation domain which radiates the reflected waves was discretized by 7-noded solid infinite elements developed by (Yerli, Temel et al. 1998).The equations of motion for each domain (dam, reservoir and foundation) were obtained by the energy method.Then, the coupled equation of motion for the system was derived using frictionless contact elements at the interfaces between mediums (Akkas, Akay et al. 1979, Hacıefendioğlu 2009, Hacıefendioğlu, Soyluk et al. 2012).The consistent mass matrices of the elements in the dam, near-field reservoir domain, and near-field foundation domain was calculated by 3x3 Gauss numerical integration method.The stiffness matrices of the elements in the dam and foundation domain were calculated by 3x3 Gauss numerical integration method, while the stiffness matrix of the elements in the near-field reservoir domain and the free surface of the reservoir was calculated by 2x2 reduced Gauss numerical integration method as recommended by (Wilson and Khalvati 1983).The mass and stiffness matrices of the infinite elements in the far-end boundaries of the reservoir and foundation domain were calculated using the Newton-Cotes numerical integration method (Bettess and Zienkiewicz 1977).
The coupled equation of motion was transformed into the Laplace domain and then solved by the direct method.In this method, the damping mechanism of the system can be assumed as either frequency-independent (hysteretic damping) or frequency dependent (Rayleigh damping) parameter.In this paper, the damping parameter was assumed as hysteretic.The coupled system was solved at the same time.The complex responses computed in the Laplace domain were finally transformed into time-domain responses by using Durbin's method (Durbin 1974).

Finite Element Formulation of the Problem Using Lagrangian Approach
The coupled equation of motion for the DRFI problem was derived using the Lagrangian approach.In this approach, the fluid elements were assumed as linearly elastic, nonrotating, and non-viscous.The basic equations for the fluid domain in the Lagrangian approach were given here according to (Akköse and Şimşek 2010).In the dynamic motion of the dam, foundation and the reservoir domain, the displacements are defined as unknown parameters.Therefore, the Latin American Journal of Solids and Structures, 2022, 19(6), e459 4/21 compatibility condition at the dam-reservoir and reservoir foundation interfaces was automatically satisfied.The stressstrain relationships of a two-dimensional fluid element in matrix form can be given as follows: where P, 11 C and v ε represents the pressure, bulk modulus and volumetric strain ( ) of the fluid element, respectively.Since the fluid element was assumed to be irrotational, a constraint parameter can be included in the stressstrain Eq. (1) using penalty methods (Klaus-Jurgen 1982, Zienkiewicz, Taylor et al. 2000).In this equation, z P , 22 C and z W denotes the rotational stress, constraint parameter and the rotation of the fluid element around z axis ( z normal to the (x-y) plane, respectively.The equation of motion for the fluid domain was derived using the energy principle; the finite element formulation of the total strain energy for the fluid domain can be written as follows: where f U and f K are nodal displacement vector and stiffness matrix of fluid domain, respectively.An important behavior of water is the capability of displacements without altering its volume.This behavior can be observed at the surface of reservoirs and storage tanks in the form of sloshing waves that move up and down.The effects of sloshing waves at the surface, which are an important aspect, can be considered in Lagrangian approach.It is feasible to consider the effect of free surface waves in the terms of fluid potential energy.The potential energy of fluid system can be written as: where sf U and f S are the up-down nodal displacement vector and stiffness matrix of surface fluid elements, respectively.Moreover, the kinetic energy of the fluid domain can be written as: where f U  and f M are the nodal velocity vector and mass matrix of fluid domain, respectively.The following set of equations that describe the motion of fluid domain can be obtained substituting Eqs.(2-4), into the Lagrange's (Clough and Penzien 1975) equation as: where f U  and  f K are the nodal acceleration vector and the stiffness matrix of fluid domain including the free surface stiffness matrix, respectively.The stiffness matrix of the entire system can be calculated by summing the stiffness matrices of all 2D elements ( with the stiffness matrices of all 1D surface elements ( Furthermore, the mass matrix of the whole system calculates by the summation of the mass matrices of all 2D elements ( ) exist in fluid domain.The stiffness and mass matrices of one element can be calculated as: Solids and Structures, 2022, 19(6), e459 5/21 where, e f B and f C in Eq. ( 6) denotes the strain-displacement matrix and elasticity matrix (Eq. ( 1) of a fluid element, respectively.Furthermore, f ρ , g , e f h and e f H expresses the mass density of water, gravitational constant and a vector containing interpolation functions of 1D surface element and 2D finite fluid element, respectively.The same procedure can be adopted deriving the FE equation of motion for dam structure.The coupled equation of motion for dam-reservoir interaction problem may be obtained by imposing the boundary conditions at the interface between reservoir and dam structure.For this purpose, the nodal displacements perpendicular to the interface between two mediums must be equal: This boundary condition can be implemented using truss elements (Akkas, Akay et al. 1979, Hacıefendioğlu 2009, Hacıefendioğlu, Soyluk et al. 2012) between interface nodes or penalty methods (Klaus-Jurgen 1982).In this paper, short and axially almost rigid truss elements were used as interface link elements.Therefore, the coupled equation of motion for dam-structure interaction problem subjected to ground acceleration can be written as follows: where, M, K, and g a represent the mass matrix, stiffness matrix of the coupled system and time varying ground acceleration vector, respectively.U  and U express the nodal acceleration vector and nodal displacement vector of the coupled system, respectively.By including the hysteretic damping parameter into the coupled Laplace transformed Eq. ( 8), a set of linear algebraic equations can be obtained as: where, s,  U , and (s) represent complex Laplace variable, Laplace transformed displacement vector and Laplace transformed ground acceleration vector, respectively.Also, η and i denote the hysteretic damping factor and imaginary number ( i = -1 ).The hysteretic damping factor in this equation can be written as η=2ζ (Gutierrez and Chopra 1978).
The solution of Eq. ( 9) at each time step can easily be implemented using Cholesky method.

Mathematical Formulations of the Infinite Element Model
Only the necessary interpolation functions for infinite elements are given here.The mapping functions of 7-noded infinite element are given as: Latin American Journal of Solids and Structures, 2022, 19(6), e459 6/21 ( ) The horizontal and vertical displacement fields of 7-noded infinite element in the Laplace domain can be calculated as: where  i N denotes the displacement shape functions of the infinite element, and can be written as: In Eq. ( 12),  ( ) k P ξ, s expresses the wave-propagation function in the infinite element.This type of function is originally developed by (Chuhan and Chongbin 1987) in the Laplace domain such as, The propagation of compression (P), shear (S) and Rayleigh (R) waves in the mediums are included in the infinite element model developed by (Yerli, Temel et al. 1998).Thus, the propagation function is written as, where α and c k are the displacement decay parameter and corresponding wave velocity, respectively.The ξ and L parameters in the function show the infinite direction and characteristic length of the element.Furthermore, the displacement behavior in many unbounded domains could be the form of 1/x function.By selecting an adequate decay coefficient α, it is possible to develop an exponential propagation function as (Bettess and Zienkiewicz 1977), Latin American Journal of Solids and Structures, 2022, 19(6), e459 7/21 where X 1 and X 4 are the x coordinate of the element at the point 1 and 4 (Figure 1).The constants a, b, and c given in the function can be calculated by using the Eq. ( 11).For this purpose, the following relationships were obtained by considering the nodes 1, 4, and 6 (or 3, 5, and 7) in the infinite element: Thus, the constants a, b and c in this equation can easily be obtained by assuming    and their satisfaction with Eq. ( 11) in three different cases (Yerli, Temel et al. 1998).The infinite element stiffness and mass matrices were calculated numerically.The typical Gauss-Legendre integration technique was applied for the finite direction.However, Newton-Cotes numerical integration scheme was used to calculate infinite integrals.The typical expression of Newton-Cotes numerical integration scheme is shown in Eq. ( 17): ( ) where F(ζ) , i β , W k and ξ k represent the displacement interpolation function, the type of propagated waves at the infinite direction weight values and abscissa values, respectively (for details: Bettess and Zienkiewicz, 1977).

VALIDATIONS OF THE DEVELOPED PROGRAMS
Many examples from previous studies were selected, and comparative studies were conducted to verify the results of the developed programs in FORTRAN 90.Generally, three types of analyzes were carried out in these investigations.Modal analysis, static, and dynamic analysis were implemented.In the modal analysis section, the natural frequencies of rigid fluid tank (Akköse and Şimşek 2010), Koyna dam (Chopra andChakrabarti 1971, Nayak andMaity 2013), Sariyar dam-reservoir interaction (Köseoğlu 2007), Sariyar dam-reservoir-foundation interaction (Köseoğlu 2007) problems were analyzed.The static analysis of Sariyar dam-foundation interaction (Bayraktar 1991) problem was carried out in the static analysis part.Finally, the dynamic analysis of Koyna dam, Koyna dam-reservoir interaction (Chopra andChakrabarti 1971, Chopra andChakrabarti 1973) and soil-structure interaction (Düzgün 2007) problems were carried out in the dynamic analysis part.In this paper, only the comparative results for natural frequencies of a rigid fluid tank (Akköse and Şimşek 2010), the dynamic analysis of Koyna dam (Chopra andChakrabarti 1971, Chopra andChakrabarti 1973) and soilsturucture interaction (Düzgün 2007)  solved by mode superposition method, while the equation of motion for soil-structure interaction problem was solved by direct method in the Laplace domain.The finite element model of rigid fluid tank is shown in Figure 2. Several studies (Lamb 1924, Olson and Bathe 1983, Wilson and Khalvati 1983, Akköse and Şimşek 2010) have examined the modal analysis of this tank.The first 10 natural frequencies of the tank were calculated and compared in Table 1.According to Table 1, the first 5 modes reflect the sloshing modal frequencies of the fluid tank, which are in the range of 0.3567-0.7653Hz.The first sloshing frequency of the tank given in Table 1 was analytically calculated by (Lamb 1924, Wilson andKhalvati 1983).Also, the sixth modal frequency given in Table 1 was the first compressive modal frequency of the rigid tank analytically calculated by (Olson and Bathe 1983).It is clearly seen that the results obtained from the developed program have an excellent agreement with analytical and numerıcal results given in (Lamb 1924, Olson and Bathe 1983, Wilson and Khalvati 1983, Akköse and Şimşek 2010).3).Detailed information about this problem can be obtained from (Chopra andChakrabarti 1971, Chopra andChakrabarti 1973).As a result of the analysis, the horizontal displacement time-histories of the crest point at the dam are comparatively given in Figure 4.The discrepancies between the maximum displacement responses given in (Chopra andChakrabarti 1971, Chopra andChakrabarti 1973) and calculated results are seen to be 1.32% for dam and 0.07% for dam-reservoir problem.The dynamic results calculated by the developed program appear to be in good agreement with the results given in (Chopra andChakrabarti 1971, Chopra andChakrabarti 1973).
Finally, the dynamic analysis of soil-structure interaction problem implemented in (Düzgün 2007) is reanalyzed here in order to control the results of direct method and the correct work of the 7-noded infinite element model in unbounded soil medium.The detailed information about this problem can be obtained from (Düzgün 2007).The comparative timehistory displacement response of the point selected at the top of the structure is given in Figure 5.The discrepancy between the maximum displacement responses is seen to be about 5%.By considering the comparative results given in Table 1, Figure 4 and Figure 5 it can be concluded that the developed programs in FORTRAN 90 give quite accurate results.It has also been seen that the coupling of finite and infinite elements works accurately to absorb the reflected dynamic waves in the unbounded mediums.

PARAMETRIC STUDIES
The dynamic responses of three different problems were examined here.Furthermore, the effect of near-field domain size on the dynamic response of DRI, DFI, and DRFI problems was evaluated by including and not including the unbounded domains in the model.The dynamic responses of these problems were analyzed under two different conditions.The boundaries of the near-field domain in rock foundation were constrained in x and y directions, only the x direction was constrained in the reservoir domain with the name WithOut Infinite Domain (WOID); in the second condition, the infinite elements were used in the far-field boundaries of the model.This condition was named With Infinite Domain (WID).These conditions are further explained in Figures 6-7.
The length and depth of the near-field foundation domain as well as the length of the near-field reservoir domain were assumed to be from 1-15 times dam height.The dimensions of domains in the problems were given in Figures 6-7, and Table 2.The suitable finite element mesh of the coupled problems was obtained by calculating the natural frequencies of the system preparing many different mesh schemes.The material properties of dam, water and rock foundation are given in Table 2.The length and stiffness of interaction elements used at the interfaces between waterdam and water-bedrock were assumed to be 0.001 m and 2x10 14 , respectively.This assumption correctly satisfies the dynamic interaction between these domains (Akkas, Akay et al. 1979).The hysteretic damping factor for the coupled system was assumed as η=10% .This damping factor is equivalent to 5% viscous damping, which was experimentally measured in several dams (Chopra 2020).The displacement decay parameter used in the infinite elements was selected as 4.08 (Bettess and Bettess 1984).The horizontal acceleration time-histories of the Koyna earthquake were imposed to the models as dynamic loads (Figure 3).In the Laplace domain analysis, the number of Laplace parameters (m) is assumed to be 10.The α parameter in the Laplace transformation is assumed to be αT = 6 (Durbin 1974) where T is the time duration and equal 10.24 sec.The time step in the analysis was selected as 0.01 sec.In each dynamic analysis, there were a total of 1,024 steps.The constraint parameter in the water was assumed to be 100 times its bulk modulus (Wilson and Khalvati 1983).The velocity of compressive (P), shear (S) and Rayleigh (R) waves in the rock foundation was assumed to be 2,300, 950 and 450 m/sec, respectively.The velocity of compressive waves in the water was assumed to be 2,100 m/sec (Schumann, Stipp et al. 2014).The flow charts for expressing the research methodology and demonstration of the steps in the developed programs are given in Figures 8-9, respectively.

RESULTS AND DISCUSSIONS
In this section, the applicability of the numerical model is assessed with dynamic responses of DRI, DFI, and DRFI problems under horizontal component of 1967 Koyna earthquake.The dynamic responses of DRI, DFI, and DRFI problems are evaluated in individual subsections.In the last part, the results of these problems are compared and discussed.

Dynamic Responses of DRI Problem
In the dynamic analysis of the DRI problem, the length of the near-field reservoir domain is varied from 1-7H, where H is the height of the dam.As mentioned before, the boundary condition at the far-field domain was modeled in two conditions; (i) constrained in the x-direction (WOID) and (ii) infinite fluid elements (WID).The foundation was assumed Latin American Journal of Solids and Structures, 2022, 19(6), e459 12/21 to be rigid (Figures 6a, 7a).The dynamic response obtained from condition (i) is expected to be inaccurate in contrast with the response of condition (ii) due to the reflection of hydrodynamic waves in the bounded reservoir domain.Condition (ii) is more realistic because the transmission of reflected hydrodynamic waves to infinity is accurately provided by using the proposed infinite fluid elements.The distribution of maximum hydrodynamic pressures at the upstream face of the dam model for different reservoir lengths is compared in Figure 10.According to Figure 10, it can be seen that the hydrodynamic pressure at the upstream face of the dam is varied by variation of the reservoir length.In general, it is seen to be augmented by increasing the reservoir length.The distribution of hydrodynamic pressure at the dam face given in Figure 10a is significantly higher than that response given in Figure 10b.The behavior of the response is seen to be not reasonable due to reflection of hydrodynamic waves in the bounded reservoir domain.In Figure 10b, the response is seen to be more reasonable and the difference between each length response is smaller than those responses given in Figure 10a.The results also reveal that the absorption capacity of the far-field reservoir domain significantly reduces the hydrodynamic pressures at the heel of the dam.Furthermore, the corresponding hydrodynamic pressure approaches infinity at the resonant frequency when the far-field reservoir boundary is constrained, absorbing energy in the far-field reservoir boundary always results in the smallest hydrodynamic response of the dam.For instance, the peak value of hydrodynamic pressure at the heel of the dam for (WOID, α = 1) and (WID, α = 1) conditions are seen to be 0.675 and 0.246 MPa, respectively.The decrease in the responses was founded at about 63%.It can be concluded that the hydrodynamic pressure will be overestimated by constraining the far-field reservoir boundaries (without infinite elements) and the hydrodynamic pressure will be underestimated by including the absorbing far-field reservoir boundary (with infinite fluid elements) in the model.
The horizontal crest displacement time histories of dam structure for reservoir lengths 1H, 4H, and 7H are given in Figure 11.According to Figure 11, it is seen that the displacement time histories are different in each reservoir length.The displacement response for (WOID, α = 1) condition is seen to be the maximum, while in (WID, α = 1) condition it is the minimum response over 1H, 4H, and 7H lengths.Similarly, the absolute maximum displacement responses at the upstream face of the dam are given in Figure 12.According to this figure, it is seen that by extending the length of the near-field reservoir, the crest displacements are also increased.For example, an increase in the near-field reservoir length for (WID, α = 1) and (WID, α = 3) conditions leads to an increase in the absolute maximum displacements from 0.0162 to 0.0178 m with a difference about 9%.With increasing the length of the near-field reservoir domain, the hydrodynamic pressure increases, and thus, the crest displacements are amplified.The maximum compressive and tensile stress responses at the base of the dam model are compared and given in Figure 13.According to this figure, the same behavior is also observed while investigating the stress responses.The stress responses at the base of the dam structure for (WID) condition are seen to be nearly the same for all reservoir lengths.Moreover, the stress responses are smaller than those responses given in (WOID) condition.Thus, it should be noted that the infinite fluid elements used in the far-field reservoir domain effectively absorb the reflected waves and they can represent reality.

Dynamic Responses of DFI Problem
The dynamic responses of the DFI problem was investigated in this section.The length and depth of rock foundation domain were varied from 1-15H in both upstream and downstream sides of the dam structure (Figures 6b, 7b).As mentioned in the previous section, two conditions (WOID) and (WID) were also considered in this problem.The displacement time histories of the crest point of the dam structure for both conditions with α = 1, 4, 7 are given in Figure 14.According to this figure, it is seen that the displacement responses are different for different α values.Moreover, by increasing the size of near-field domains in the problem, the crest displacement is also increased.While including the radiation of waves in the foundation domain, the displacement responses are significantly decreased (Figure 14b).By imposing constrained boundary conditions to the far-field foundation domain, the responses are increased due to the Latin American Journal of Solids and Structures, 2022, 19(6), e459 13/21 reflection of waves in the bounded rock foundation domain.The absolute maximum displacement responses at the upstream face of the dam structure are given in Figure 15.Only the responses for α =1-7 are given.According to Figure 15, it is seen that by increasing the near-field domain size, the displacement responses significantly increase in both conditions.But, the responses of (WID) condition for all α values are the smallest due to absorption of the reflected wave in the unbounded foundation domain.The maximum compressive and tensile stress responses at the base of the dam structure is given in Figure 16; in this figure, only the responses for α =1-7 in both conditions are given.It is also observed that by increasing the near-field domain size, the stress responses at the base of the dam significantly increase.This trend in (WID) condition is much slower and smaller than (WOID) condition.

Dynamic Responses of DRFI Problem
The dynamic responses of the DRFI problem are investigated in this section.The length of the near-field reservoir as well as the length and depth of the near-field rock foundation at upstream face and downstream face sides are varied from 1-15H (Figures 6c, 7c).Only the responses for α = 1-7 are given here.The hydrodynamic pressure at the dam face in (WOID) and (WID) conditions are given in Figure 17.According to this figure, it is seen that the hydrodynamic pressure at the dam face increases by increasing the dimension of near-field domains.The behavior of hydrodynamic response in (WOID) and (WID) conditions can be seen to be different due to different boundary conditions imposed on the far-field boundaries.The response in (WID) condition is seen to be significantly bigger than (WOID) condition due to the reflection of waves in the bounded domains.The hydrodynamic pressure response in (WOID) condition is far from reality because the transmission of reflected waves to unbounded domains (water and rock foundation) is not provided.It can be seen that in the DRI and DRFI problems, the hydrodynamic pressure increases due to including the effect of rock foundation in the analysis.For example, the peak of hydrodynamic pressure at the heel of the dam for DRI and DRFI problems in (WID, α = 3) condition is 0.311 and 0.972 MPa, respectively.In other words, the peak of hydrodynamic pressure is increased by about 68% for the dam with the rock foundation domain effect.The displacement time-history response for the crest point of the dam structure is given in Figure 18.Only, the responses for α = 1, 4, and 7 are compared.It is observed that the responses are generally increased by increasing the size of near-field domains.Furthermore, the responses for (WID) condition are seen to be smaller than the responses for (WOID) condition.The same behavior is also observed in the absolute maximum displacement responses at the upstream face of the dam structure which is given in Figure 19.The significant stress responses in the dam body generally occur at the heel, toe, and some regions at the below of the neck point.The stress concentration generally occurs in these regions of the dam body.Therefore, the variation of compressive and tensile stress responses at the heel of dam structure from α = 1-15 in (WOID) and (WID) conditions are comparatively investigated here.The stress responses of DRI, DFI, and DRFI problems are individually evaluated in Figures 21a,  21b, 21c, and finally the responses of these problems are compared with each other in Figure 21d.
According to Figure 21a, it is seen that the stress responses vary with the variation of near-field domain size in (WOID) and (WID) conditions.Moreover, the responses after α = 4.5 for both conditions are seen to be nearly the same.The maximum tensile stress for condition (WOID) is calculated as 3.159 MPa and for condition (WID), it is calculated as 3.138 MPa.The difference between responses is seen to be 0.66%.The same evaluation is also conducted on DFI (Figure 21b) and DRFI (Figure 21c) problems.And it is seen that the stress responses after α = 12 are nearly the same and the variations are negligible.The maximum compressive stress responses at the heel of the dam in the DFI problem for (WOID, α = 12) and (WID, α = 12) conditions are seen to be -5.7 and -5.409 MPa, respectively.The difference between responses is calculated as 5.1%.The same evaluation for the DRFI problem is conducted and the difference is found to be 0.24%.The stress responses at the heel of dam structure for DRI, DFI, and DRFI problems are comparatively evaluated in Figure 21d.According to this figure, the stress response of the DRI problem is significantly smaller than the DRFI problem.Thus, the assumption of a rigid foundation in the DRI problem gives inaccurate response.To further clarify this condition, the compressive stress responses at the heel of the dam in the DRI and DRFI problems for (WID, α = 4.5) condition are seen to be -2.953 and -38.563MPa, respectively.The difference between responses of the DRFI and DRI problems is calculated as 92.34%.The effect of reservoir on the stress responses at the heel of the dam is also important.To demonstrate the importance, the stress response in the DRFI and DFI problems is evaluated considering (WID, α = 4.5) condition.The tensile stress responses at the heel of the dam in the DRFI problem can be seen to be 44.452MPa and it is 21.966 MPa in the DFI problem.The rise of stress at that point is seen to be 50.5%.According to Figure 21d, it can be seen that the radiation effect of waves in the bounded domains for the DFI and DRFI problems can effectively be provided if the near-field domains size would be assumed as 12H.By taking the size of near-field domains as 12 times dam height, there is no need for including the far-field domains in the model.Because the near-field domain size is sufficiently big to absorb the reflected dynamic waves.And thus, the truncated boundaries can be constrained in the x and y direction for foundation and only the x-direction for water.At the truncated boundaries, there is no need for infinite elements or transmitting boundaries anymore.The stress responses in (WOID, α = 12) and (WID, α = 12) conditions are seen to be approximately the same.The truncation of boundaries at this size is not economic because the analysis process takes a long time and it is time-consuming.Furthermore, if the stress level in (WOID, α = 12) or (WID, α = 12) condition is the correct response, the same stress level can be obtained by assuming the near-field domain size of 1.5H.This condition is schematically explained in Figure 22 for tensile stress variation at the heel of the dam structure.The same result can also be obtained while observing the compressive stress response.
The near-filed domain size in the DRFI problem is concluded to be three times (Chen, Yang et al. 2019, Gorai andMaity 2019) and five times (Asghari, Taghipour et al. 2018) dam height.The size which gives the maximum averaging dam responses is concluded as an appropriate near-field domain size in these papers.This evaluation for determining the size of near-filed domains is inaccurate and has no physical interpretation.The size in which the responses in both (WOID) and (WID) conditions are nearly the same or constant can be the appropriate near-field domain size.This evaluation has a physical meaning and it is naturally correct.Moreover, significant errors can occur in the responses when the near-field domain size is assumed to be 3 or 5 times dam height.For example, the tensile stress responses at the heel of the dam in the DRFI problem for (WID, α = 1.5) and (WID, α = 3) conditions are seen to be 16.116 and 43.307 MPa, respectively.The error is calculated as 62.78%.The tensile stress at the heel of the dam for (WID, α = 1.5) and (WID, α = 5) conditions are seen to be 16.116 and 25.995 MPa, respectively.The error is calculated as 38%.
In this paper, the dynamic analysis of a DRI problem selected from (Küçükarslan 2004b) took 150 seconds using a Latitude E6430 x64-based computer, while the same problem took 367 seconds to solve in (Küçükarslan 2004b).However, the time duration for dynamic analysis of this problem depends on the number of nodes, elements, time steps in the analysis, and method of analysis.By comparing the time durations regarding the solution of the DRI problem, it was seen that the developed numerical model is faster than many other approaches.

CONCLUSIONS
In the Lagrangian modeling of DRI, DFI, and DRFI problems an efficient numerical model was developed in FORTRAN 90 programming language, which is more accurate and faster than many approaches.The coupled equation of motion which gives the complex displacement responses of these problems is solved by the direct method in the Laplace domain.Subsequently, the complex displacement responses were transformed into time domain using Durbin's method.Other dynamic responses were derived utilizing the reel displacement responses.Furthermore, the effect of the near-field domains sizes on the dynamic responses of these problems were investigated including and not including the far-field domains, and an appropriate size for the near-field domains was recommended.The further results can be listed as follow: • The proposed numerical model is seen to be more accurate and faster than many other methods developed in the literature.

•
The size of the near-field domains affects the stress responses.The stress responses for (WOID) and (WID) circumstances, however, were seen to be the same or constant at a certain size.This size was determined to be 4.5 times dam height for the DRI problem and 12 times dam height for the DFI and DRFI problems.This indicates that if the size of the near-field domains is large enough, the reflected waves in the bounded domains can be totally absorbed.

•
If the length of the reservoir domain is assumed to be 4 times the height of the dam, there is no need for transmitting boundaries at the far-field reservoir domain in the DRI problem.Similarly, if the size of the near-field domains is assumed to be 12 times dam height, there is no need for transmitting boundaries at the far-field rock domain in the DFI and DRFI problems.

•
The stress level obtained at the heel of dam structure for the DFI and DRFI problems in (WOID, α = 12) or (WID, α = 12) conditions is seen to be equivalent to the stress level for (WID, α = 1.5) condition.Furthermore, the responses obtained in these conditions are the appropriate response.
• Significant errors have been detected in the stress responses of the heel point in the DFI and DRFI problems when the near-field domain size is assumed to be 3 or 5 times dam height.The errors in these problems for sizes 3 and 5 were seen to be 62.78% and 38%, respectively.

•
The dynamic responses of the handled problems in this study are amplified by increasing the size of the near-field domains.This increasing trend is bigger in (WOID) condition than in (WID) condition.
• By restricting the far-field reservoir boundaries in the DRI problem, the hydrodynamic pressure will be overestimated, while including the proposed infinite fluid elements at far-field reservoir boundaries the hydrodynamic pressure will be underestimated.The difference between responses at the heel of the dam for these conditions was 63%.

•
In the DRI and DRFI problems, the hydrodynamic pressure increases due to inclusion of rock foundation in the model.The assumption of rigid base foundation gives inaccurate hydrodynamic pressure response.The difference between hydrodynamic responses at the heel of the dam for these problems was 68%.
In this paper, the dynamic loads generated by horizontal component of Koyna earthquake data are applied to all FE nodes of the model.The static and the hydrostatic loads were not included in the analysis.Only the severe loads are applied to model.As a future work, the findings of this paper can be combined by considering all types of loads in the problem, as well as spatially varying ground motion data.

Figure 1 :
Figure 1: 7-noded infinite element model used far-field reservoir and rock foundation boundaries, a) dimensional coordinates, b) non-dimensional coordinates

Figure 2 :
Figure 2: Finite element mesh of rigid fluid tank

Figure 4 :Figure 5 :
Figure 4: Validation of the developed program with the dynamic results given in(Chopra andChakrabarti 1971, Chopra andChakrabarti 1973)

Figure 6 :
Figure 6: A concrete gravity dam model interacts with, a) bounded reservoir, b) bounded rock foundation, c) bounded reservoir and rock foundation at the same time

Figure 7 :
Figure 7: A concrete gravity dam model interacts with, a) unbounded reservoir, b) unbounded rock foundation, c) unbounded reservoir and rock foundation at the same time

Figure 8 :Figure 9 :
Figure 8: Flow chart for research methodology

Figure 10 :
Figure 10: The effect of near-field reservoir length on the distribution of hydrodynamic pressure at upstream face of the dam in DRI problem, a) without infinite domain, b) with infinite domain

Figure 11 :
Figure 11: The effect of near-field reservoir length on the horizontal crest displacement of the dam in DRI problem, a) without infinite domain, b) with infinite domain

Figure 12 :
Figure 12: The effect of near-field reservoir length on the horizontal displacement at upstream face of the dam in DRI problem, a) without infinite domain, b) with infinite domain problems Ahmad Yamin Rasa et al.Latin American Journal of Solids and Structures, 2022, 19(6), e459 14/21

Figure 13 :Figure 14 :Figure 15 :
Figure 13: The effect of near-field length on the maximum compressive and tensile stresses at the base of dam in DRI problem, a) without infinite domain, b) with infinite domain

Figure 16 :
Figure 16: The effect of near-field rock foundation size on the maximum tensile and compressive stresses at the base of dam in DFI problem, a) without infinite domain, b) with infinite domain

Figure 17 : 21 Figure 18 :Figure 19 :Figure 20 :Figure 21
Figure 17: The effect of near-field domains size on the distribution of hydrodynamic pressure at upstream face of the dam in DRFI problem, a) without infinite domain, b) with infinite domain

Figure 22 :
Figure 22: The variation of stress level at the heel of the dam in DRFI, DFI and DRI problems according to near-field domains size problems were given.The coupled equation of motion for Koyna dam problem was Latin American Journal of Solidsand Structures, 2022, 19(6), e4598/21

Table 1
First 10 natural frequencies of rigid fluid tank (rad/sec)

Table 2
Material properties and dimension of the problem