INTRODUCTION
Modeling behavior of concrete materials in structural components is one of
the most interesting fields in structural engineering. Complex models are required
in order to capture the nonlinear behavior that arises as the structure is loaded.
Several researchers have introduced models in order to simulating the nonlinear
behavior and also failure of concrete under static and dynamic loads. The methods
based on Continuum Crack Model (CCM), Discrete Crack Model (DCM), Interface
Crack Approach (ICA) and the method of constrains are samples of most favorite
methods for numerical simulation of concrete behavior. Utilizing these methods
in conjunction with analytical solutions depends on the complexities of the
structures, material properties and the boundary conditions. Numerical models
such as finite element; finite difference; finite volume; extended finite element
and mesh free approaches are common methods for simulation of structural components
and their behavior (Yu et al., 2008).
The DCM requires monitoring the response and modifying topology of the mesh
corresponding to the current crack configurations at each state of loading.
However, this approach explicitly represents the crack as a separation of nodes
which is a more realistic representation of the opening crack. This model is
also useful when the location and direction of cracks is recognizable before
loading the structure. Hohberg (1990) studied the mechanism
of joint elements under water pressure using DCM. Ahmadi
and Razavi (1992) represented a finite element model of discrete cracks
for modeling joints. They considered perfectly elasticplastic behavior for
joint in tension and linear elastic behavior in compression and shear. Ahmadi
et al. (2001) introduced a nonlinear joint element with coupled tensionshear
behavior for analysis of arch dams. Lotfi and Espandar (2004)
used discrete crack method, nonorthogonal smeared crack and combination of
them for seismic analysis of dams. Du and Tu (2007)
combined explicit finite element method with transmitting boundary to study
the effects of contraction joint opening in concrete dams.
The methods based on CCM are divided to two major groups, i.e., damage mechanics approach and smeared crack approach. In the smeared crack approach cracks and joints are modeled in an average sense by appropriately modifying the material properties at the integration points of regular finite elements. Smeared cracks are convenient when the crack orientations are not known beforehand, because the formation of a crack involves no remeshing or new degrees of freedom.
Some researchers applied continuum damage mechanics to study the concrete behavior
such as Sumarac et al. (2003), Labadi
and Hannachi (2005), Contrafatto and Cuomo (2006),
Grassl and Jirasek (2006), He et
al. (2006), Cicekli et al. (2007) and
Khan et al. (2007). Also, Mirzabozorg
et al. (2004) utilized damage mechanics approach to conduct seismic
nonlinear analysis of concrete gravity dams in 2D space including damreservoir
interaction effects. Ardakanian et al. (2006)
considered nonlinear seismic behavior of mass concrete in 3D space which is
based on an anisotropic damage mechanics model.
Failure based on the smeared crack approach and classification of its branches
were studied by Malvar and Fourney (1990), Weihe
et al. (1998), Mosler and Meschke (2004)
and Phama et al. (2006). Mirzabozorg
and Ghaemian (2005) developed a model based on smeared crack approach in
3D space. In their study, they analyzed 3D models including damreservoir interaction
effects and considered nonlinear behavior of the structure. In addition, Mirzabozorg
et al. (2007) investigated nonuniform cracking in smeared crack
approach for 3D analysis of concrete dams and Mirzabozorg
et al. (2010) studied nonlinear behavior of concrete dams under nonuniform
earthquake ground motion records.
In the present study, a coaxial smeared crack model was introduced in 3D space for crack analysis of concrete structural components under static and dynamic conditions. For this purpose an advanced failure criterion is utilized. Finite element model of a set of beams are studied for constant and variable shear transfer coefficient conditions. In addition, capability of the proposed method is investigated for considering the behavior of coupled structurefluidsoil system under various load types.
CONSTITUTIVE LAW FOR CONCRETE
The proposed method for crack analysis of concrete structural components should be able to simulate the behavior of the element in various states as following; Presoftening behavior; fracture energy conservation; nonlinear behavior during the softening phase and finally crack closing/reopening behavior. The following subsections represent a brief review on general concepts of these stages.
Presoftening phase: Generally, the relationship of the stress and strain vectors at the presoftening phase is given as:
where, [D]_{elastic} is the elastic modulus matrix; {σ} and {ε} are the vector of stress and strain components. The modulus matrix in elastic condition maybe defined for isotropic, orthotropic and anisotropic materials.
Softening phase: During the softening phase, the elastic stressstrain relationship is substituted with an anisotropic modulus matrix which corresponds to the stiffness degradation level in the three principal directions. In the present study, the secant modulus stiffness approach, SMS, is unitized for the stiffness formulation in which the constitutive relation is defined in terms of total stresses and strains, shown in Fig. 1.
The stiffness modulus matrix based on the smeared crack propagation model is given in Eq. 2. It is worth noting that the extracted modulus matrix is coaxial with the principal strains in the considered location within the cracked element:
Where:
where, η_{1}, η_{2} and η_{3} are the
ratio of the softened Young’s modulus in the three principal directions and
the initial isotropic elastic modulus and β_{12}, β_{23}
and β_{13} are shear retension factors corresponding to the principal
directions given as:

Fig. 1: 
SMS formulation of stiffness modulus matrix 
The constitutive matrix given in Eq. 2 is transformed to the global coordinate system as following:
where, [T] is the strain transformation matrix in 3D space. Based on the maximum
strain reached in each principal direction, the secant modulus matrix is determined.
Increasing the normal strain in each principal direction leads to reduction
of the corresponding softened Young’s modulus. Finally, when the maximum
strain reaches the fracture strain, the considered Gaussian point within the
element in the corresponding direction is fully cracked and its contribution
in the stiffness matrix of the element is eliminated. Based on Eq.
2 to 5, any change in principal strains or their directions
in each Gaussian point leads to an update requirement of the global constitutive
matrix, .
Satisfying the energy conservation principle in each Gaussian point leads to
the fracture strain under static and dynamic loads:
where, h_{c} is the characteristic dimension of the considered Gaussian
point and is assumed equal to the cubic root of the Gaussian point’s volume contribution;
σ_{0} is the stress corresponding to the softening strain and G_{f}
is the specific fracture energy. The primed quantities show the dynamic constitutive
parameters.
The strainrate sensitivity of fracture energy is applied through a dynamic magnification factor DMF_{f}, so that:
Failure criterion: The strength of concrete under multiaxial stresses
is a function of the state of stress and cannot be predicted by limitation of
simple tensile, compressive and shearing stresses independently of each other.
In the elasticity based models, a suitable failure criterion is incorporated
for a complete description of the ultimate strength surface. Criteria such as
yielding, load carrying capacity and initiation of cracking have been used to
define failure (Babu et al., 2005). Many failure
criterions have been proposed for brittle material as well as concrete. Some
of more familiar criterions are MohrCoulomb criteria, DruckerPrager, Chen
and Chen (1975), Ottosen (1977), Hsieh
et al. (1982), Willam and Warnke (1974), Menetrey
and William (1995), Sankarasubramanian and Rajasekaran
(1996) and Fan and Wang (2002).
In the present study an advanced failure criterion is used for initiation and
propagation of cracks in concrete structural components (Menetrey
and William, 1995). In the general form this yield function can be expressed
as follow:
where, ξ is hydrostatic stress invariant, ρ is deviatory stress invariant, θ is deviatory polar angle, r (θ, e) is an elliptic function, e describes the shape of the deviatory trace, m represents the frictional resistance of material, c is cohesion of material and f'_{c} is uniaxial compressive strength of concrete. The main parameters in the above equation are defined as follow:
where, I_{1} is the first invariant of the Cauchy stress tensor; J_{2} and J_{3} are the second and third invariants of the deviatory part of the Cauchy stress tensor; σ_{ii} is principal stress; S_{ij}, S_{ji}, S_{jk }and S_{ki} are deviatory stresses. In addition, the elliptic function is in the form of Eq. 10, as following:
FORMULATION OF FLUIDSTRUCTUREINTERACTION
The general coupled equation of motion for Fluidstructureinteraction (FSI)
under dynamic loads is discussed in this section. Considering the coupled damreservoirfoundation
system, the governing equation in the reservoir medium is Helmholtz equation
from the Euler’s equation given as Ghaemian et
al. (2003):
where, p, C and t are the hydrodynamic pressure, pressure wave velocity in
liquid and time, respectively. Boundary conditions required to apply on the
reservoir medium to solve Eq. 11 are explained HaririArdebili
and Mirzabozorg (2010). The equations of the damfoundation (as the structure)
and the reservoir take the form:
where, [M], [C] and [K] are the mass, damping and stiffness matrices of the
structure including the dam body and its foundation media and [G], [C’]
and [K’] are matrices representing the mass, damping and stiffness equivalent
matrices of the reservoir, respectively. The matrix [Q] is the coupling matrix;
{f_{1}} is the vector including both the body and the hydrostatic force;
{P} and {U} are the vectors of hydrodynamic pressures and displacements, respectively
and {Ug} is the ground acceleration vector.
FINITE ELEMENT MODELING AND CASE STUDIES
In order to investigate the ability of proposed method in static and dynamic
analysis of concrete structural components two kinds of models were developed.
The first one is a simple notched beam under threepoint bending test and the
second one is coupled prototype gravity damreservoirfoundation system under
a combination of various load types. Figure 2 shows the model
used for the threepoint bending test. It’s a squaresection concrete beam
with an initial notch depth of 51 mm in its center while the depth of the beam
is 102 mm (Malvar and Warren, 1988). The material properties
are as follow; E = 21.7 GPa, v = 0.2, f’_{t} = 2.4 MPa, f’_{c}
= 29.0 MPa and G_{f }= 35 N m^{1}. Finite element model of
half of concrete beam with fine mesh is also depicted in Fig.
2.
Figure 3 shows the finite element model and also dimension
of prototype gravity damreservoirfoundation system. The model consists of
solid elements for simulation of body and massless foundation and also fluid
elements for modeling the reservoir.

Fig. 2: 
Geometry and finite element model of singleedgenotched beam
subjected to three point bending test 

Fig. 3: 
Finite element model, dimensions and boundary conditions in
coupled system 
Boundary conditions are also depicted in this Fig. 3. Material
properties of mass concrete in static and dynamic conditions are as follow;
Modulus of elasticity is 40 and 46 GPa, Poisson’s ratio is 0.2 and 0.14,
density is 2640 kg m^{3}, tensile strength is 2 MPa and dynamic tensile
stress is 3 MPa. Modulus of elasticity in foundation area is 30 GPa and Poisson’s
ratio is 0.2. Sound velocity in water is assumed to be 1460 m sec^{1}.
Three types of dynamic loads were used in present study, i.e., intensifying
sinusoidal loading, increasing step loading and real earthquake ground motion.
The two first load types were selected to investigation the pattern, route and
the type of cracking of huge mass concrete specimen under dynamic loads including
fluidstructure interaction. Acceleration timehistories of all three loads
are depicted in Fig. 4ac. Damping was assumed
to be 10% of critical damping in all cases.
RESULTS AND DISCUSSION
Verification of proposed model: In Fig. 5, the load
versus loadpoint deflection curves from the finite element model is compared
with the experimental result of Malvar and Warren (1988).

Fig. 4(ac): 
Acceleration timehistories applied to the coupled system
(a) Intensifying sinusoidal load, (b) Increasing step load and (c) Real
earthquake ground motion 

Fig. 5: 
Load versus loadpoint deflection compared with the experimental
result for constant and variable shear transfer coefficients 
Two type of smeared crack were assumed for each beams, i.e., in the first approach
shear transfer coefficient is assumed to have two constant parameters for open
and closed cracks and in the second case variable shear transfer coefficient
is take into account in smeared crack model. In the constant coefficient approach
the preassumed values in open and cracked conditions are 0.1 and 0.9, respectively.

Fig. 6: 
Final crack profile in concrete beam under the concentrated
load at the middle span 
There is good agreement between predictions of the two approaches through all
loading stages: elastic, hardening and softening. The numerical results also
agree very well with the experimental result which demonstrates the soundness
of the present algorithm. The final crack profiles of concrete beam are shown
in Fig. 6. Although, the general pattern of crack profiles
based on constant and variable shear transfer coefficients are close to each
other, some differences can be shown in final crack profile. The crack is started
from upper parts of the notch near the corner and extended almost linearly toward
the upper edge of the beam.
Concrete dam model and FSI problem: In this section application of proposed smeared crack model is investigated for the problems consisting fluidstructureinteraction. For this purpose the introduced coupled system in previous section is selected as case study. The modeled concrete structure is assumed to be a prototype of a concrete gravity dam, its support assumed to be modeling of foundation rock and finally the fluid is considered as reservoir water.
So the damreservoirfoundation system is subjected to various types of loading as mentioned before. In each case the nonlinear responses of the system are compared with reference linear responses. Three nonlinear models are considered for each loading in which the first two ones use constant shear transfer coefficients for open and closed cracks and the third one uses proposed rotating smeared crack model with ability of updating shear transfer coefficient in each load step of dynamic analysis.
Intensifying sinusoidal loading: In this section the results of coupled
system under an intensifying sinusoidal loading are investigated. The values
of shear transfer coefficients for open and close cracks were assumed to be
0.1 and 0.9 in the first case and 0.3 and 0.7 in the second one. Figure
7 shows the timehistory of crest displacement in stream direction for various
linear and nonlinear models. As it is clear using nonlinear model based on smeared
crack model lead to higher displacements than to linear one. Because of continuously
intensifying nature of this loading the behavior of system gradually goes from
linear elastic range to nonlinear phase and finally leads to global instability
of the coupled system.
Based on this Fig. 7 up to time t = 1.36 sec all models show
the same behavior which shows that models are in linear elastic range. After
t = 1.36 sec the nonlinear behavior of models are begin in the form of concrete
cracking.

Fig. 7: 
Crest displacement under intensifying sinusoidal dynamic load
for linear and nonlinear models 
All three models have different behavior showing the importance of shear transfer
coefficient in nonlinear analyses. The βconstant (0.1/0.9 values) model
fails earliest at the time t = 4.14 sec showing that in the current model, assuming
value near unit for shear transfer coefficient in closed cracks condition and
value near zero in open cracks condition can lead to unrealistic results. The
βconstant (0.3/0.7 values) model is the second one which fails under the
same loading at time t = 4.28 sec. Although this models endure more than previous
one under intensifying loading, it has limitations for fully capture the nonlinear
behavior of concrete under dynamic loading. The third model in which β
varies based on the proposed model, experiences some large deformations up to
28 mm in first 5 sec of dynamic load. However, there is no global instability
in this case.
Figure 8ad represent the nonconcurrent
envelope of the first principal stress within the dam body. As it is shown,
considering that there is no limitation in linear model its stresses are increased
to 12.7 MPa which is mainly concentrated near the toe and heel of base. Using
nonlinear models concentrate the area with high tensile stress in damfoundation
interface and also some areas around the neck. Consequently theses areas have
great capability for cracking. The similarity of βvariable model to βconstant
(0.1/0.9 values) is more than βconstant (0.3/0.7 values). Figure
9ac show the propagation of crack profile by increasing
the intensity of loading for nonlinear models.

Fig. 8(ad): 
Nonconcurrent envelope of first principal stresses in dam
body (a) Linear model, (b) βvariable model, (c) βconstant(0.1/0.9)
and (d) βconstant(0.3/0.7) 

Fig. 9(ac): 
Crack propagation in dam body and final crack profile in failure
time (a) βvariable model, (b) βconstant(0.1/0.9) and (c) βconstant(0.3/0.7) 
Also the final crack profiles are shown in this figure. As can be seen at
t = 2.00 sec almost all models have close crack profiles and at t = 3.00 sec
βvariable and βconstant (0.1/0.9 values) models have considerable
cracked element in neck area. For the final crack profiles all models show complete
cracking in base of the dam and in neck area.
Increasing step loading: This subsection is provides the responses
of damreservoirfoundation system under increasing steptype loading. This
load type is similar to previous one by this way that both is increasing and
is different from that one considering that it applies acceleration to system
which is not damps rapidly. Figure 10 shows timehistory
of crest displacement in stream direction for linear and nonlinear models. As
can be seen using linear model generate uniform displacement curve similar to
input acceleration while the damping of acceleration in the system is clearly
recognizable.

Fig. 10: 
Crest displacement under increasing steptype dynamic load
for linear and nonlinear models 
Based on these results the coupled system damps the input acceleration in
time interval about 0.33 sec. Comparing the nonlinear models shows that βvariable
model fails at t = 6.14 sec while using βconstant(0.1/0.9) leads to failure
of the system at t = 8.08 sec, respectively. In addition based on βconstant(0.1/0.9)
model, the damping of input acceleration in nonlinear phase is done very slowly.
Figure 11ac represent the nonconcurrent
envelope of first principal stress within the dam body under spettype loading.
As it is shown considering that there is no limitation in linear model its stresses
are increased to 11.92 MPa which is mainly concentrated near the toe and heel
of base and also in middle parts of upstream face. The stress envelopes of variable
and constant shear transfer coefficients are very close to each other, while
βconstant(0.1/0.9) model generates more high stress areas in base of
dam. Figure 12a and b show the propagation
of crack profile at different times for nonlinear models. As can be seen the
final crack profile for the βvariable model at time 8.94 sec and βconstant
(0.1/0.9 values) model at t = 8.54 sec are almost the same. In addition for
the lower times βconstant model generates more damage in base of dam.
There is no crack in neck area under this type of loading.
Real ground motion loading: Figure 13 shows timehistory
of crest displacement for linear and nonlinear models using real ground motion.
As can be seen because moderate intensity of selected ground motion none of
nonlinearbased models were failed.

Fig. 11(ac): 
:Nonconcurrent envelope of first principal stresses in dam
body (a) Linear model; (b) βvariable model and (c) βconstant(0.1/0.9) 

Fig. 12(ab): 
Crack propagation in dam body and final crack profile in failure
time (a) variable model and (b) constant(0.1/0.9) 

Fig. 13: 
Crest displacement under real ground motion for linear and
nonlinear models 

Fig. 14(ac): 
Nonconcurrent envelope of first principal stresses in dam
body (a) Linear model, (b) βvariable model and (c) βconstant(0.1/0.9) 

Fig. 15(ab): 
Crack propagation in dam body and final crack profile in failure
time (a) βvariable model and (b) βconstant(0.1/0.9) 
Although, both nonlinear models generate larger values for displacement that
to linear model, using βconstant(0.1/0.9) model leads to some larger
displacement that to βvariable model and two models have great consistency
with each other.
Figure 14ac represent the nonconcurrent
envelope of first principal stress within the dam body under real ground motion.
Linear model experiences up to 15.53 MPa tensile stress in toe of dam. The general
templates of two nonlinear models are very close each other and overstressed
area is concentrated in base of dam at vicinity of foundation. Figure
15a and b shows the propagation of crack profile for
nonlinear models under real ground motion. In spite of the previous load types,
there is no continuously increasing loading in this example. Like the Fig.
13 and 14, there is great similarity between crack profiles
not only for the final profiles but also the propagation of crack based on two
approaches have almost same trend.
CONCLUSION
In the present study a coaxial rotating smeared crack model was introduced for nonlinear behavior of mass concrete in 3D space. The advantages of proposed model are; using variable shear transfer coefficient which is updated in each load step; utilizing an advanced failure criterion for concrete in 3D space and ability of this method in simulation of cracking process in concrete with high accuracy. The proposed model was verified using the threepoint bending test of concrete beam. Results show very good consistency between numerical analysis and experimental test. On the other hand, finite element model of prototype gravity damreservoirfoundation system was provided in order to investigation the nonlinear dynamic behavior of large concrete specimens considering fluidstructureinteraction. Three types of dynamic loading were selected for this purpose; i.e., two increasing loads and one real ground motion record. The results were compared for constant and proposed variable shear transfer coefficients. It was found that generally shear transfer coefficient affects the results of dynamic analysis of concrete structures meaningfully. In high intensity levels of dynamic load, constant shear transfer coefficient leads to early failure of structure than to variable coefficient model. Also both final crack profile and the propagation of cracks within concrete specimens are highly affected by shear transfer coefficient. Finally it was concluded that the proposed model can be used for static and dynamic crack analysis of concrete structural components considering the effects of fluidstructureinteraction and hydrodynamic pressure on specimen.