On damping created by heterogeneous yielding in the numerical analysis of nonlinear reinforced concrete frame elements
OON DAMPING CREATED BY HETEROGENEOUS YIELDING INTHE NUMERICAL ANALYSIS OF NONLINEAR REINFORCEDCONCRETE FRAME ELEMENTS
Pierre Jehel ∗ and R´egis Cottereau Laboratoire MSSMat / CNRS-UMR 8579, ´Ecole Centrale Paris, Grande voiedes Vignes, 92295 Chˆatenay-Malabry Cedex, France Department of Civil Engineering and Engineering Mechanics, ColumbiaUniversity, 630 SW Mudd, 500 West 120th Street, New York, NY, 10027, USA
March 2 nd , 2015: Paper accepted for publication in Computers and Structures
DOI: 10.1016/j.compstruc.2015.03.001
Abstract.
In the dynamic analysis of structural engineering systems, it iscommon practice to introduce damping models to reproduce experimentallyobserved features. These models, for instance Rayleigh damping, account forthe damping sources in the system altogether and often lack physical basis. Wereport on an alternative path for reproducing damping coming from materialnonlinear response through the consideration of the heterogeneous characterof material mechanical properties. The parameterization of that heterogeneityis performed through a stochastic model. It is shown that such a variabilitycreates the patterns in the concrete cyclic response that are classically regardedas source of damping.
Keywords: damping ; concrete ; nonlinear constitutive relation ; material hetero-geneity ; stochastic field 1.
Introduction
In the last few decades, a great deal of attention was paid to the comprehensionand modeling of damping mechanisms in inelastic time-history analyses (ITHA) ofconcrete and reinforced concrete (RC) structures [2, section 2.4]. Figure 1, adaptedfrom [33], shows the uniaxial cyclic compressive strain-stress ( E -Σ) response mea-sured on a concrete test specimen. Throughout this paper, the term “uniaxial”implies that there is only one loading direction and that the stress, respectivelystrain, of interest is the normal component of the stress, respectively strain, vectorin the loading direction. In other words, when it comes to constitutive relation be-tween stress and strain, the work presented thereafter is developed in a 1D setting.In figure 1, the so-called backbone curve, which is the envelope of the response(dashed line), shows the following phases: (i) an inelastic phase with positive slope( E ≤ . × − for that particular example, where E is the measured strain),and (ii) an inelastic phase with negative slope before the specimen collapses. Forconcrete, no elastic phase can really be identified, and hysteresis loops appear in ∗ Corresponding author: pierre.jehel[at]centralesupelec.fr a r X i v : . [ c s . C E ] M a r unloading-reloading cycles even for limited strain amplitudes. Other salient fea-tures include: (i) a residual deformation after unloading, and (ii) a progressivedegradation of the stiffness (slope of the unloading-loading segments). The hys-teresis loops are one of the sources of the damping that is observed in free vibrationrecordings of concrete beams. Other sources include friction at joints [27] or at theconcrete-steel interface in reinforced concrete [14]. These other sources of dampingwill not be discussed in this paper, where we will concentrate on material damping. unloading-reloading cycles even for limited strain amplitudes. Other salient fea-tures include: (i) a residual deformation after unloading, and (ii) a progressivedegradation of the stiffness (slope of the unloading-loading segments). The hys-teresis loops are one of the sources of the damping that is observed in free vibrationrecordings of concrete beams. Other sources include friction at joints [27] or at theconcrete-steel interface in reinforced concrete [14]. These other sources of dampingwill not be discussed in this paper, where we will concentrate on material damping. Figure 1.
Strain-stress concrete experimental response in pseudo-static cyclic uniaxial compressive loading (adapted from [33]). Σand E are the homogeneous compression stress and strain in theconcrete test specimen that are measured in the loading direction.Σ and E are spatial mean quantities in the sense that Σ is com-puted as the load in the hydraulic cylinder of the testing machinedivided by the area of the specimen cross section, and E is com-puted as the displacement of the cylinder divided by the length ofthe concrete specimen.Most classical uniaxial constitutive models of concrete for numerical simulationdo not dissipate any energy in unloading-reloading cycles (see e.g. figure 2 [top left]).It is then common practice to add a viscous damping model (Rayleigh damping)to the inelastic structural model, to reproduce phenomena that are experimentallyobserved at the structural level (decreasing amplitude of displacements in free vi-bration). However, Rayleigh damping is well known to lack physical justification,even when care is taken to avoid generating spurious damping forces [7, 18].Another class of approach aims at reproducing more precisely the features ofFigure 1 through elaborate inelastic constitutive relations. Figure 2 shows typi-cal examples of uniaxial relations found in the literature. The relation describedin [9] [top left] defines different response phases for different strain intensities, withan additional coefficient to control the loss of stiffness. The constitutive relationdescribed in [23] [bottom right] comes from a formulation developed in the frame-work of thermodynamics with internal variables. It reproduces damping featuresreasonably well, in particular for higher amplitudes, but requires the identifica-tion of a rather large number of parameters. These first two types of relations Figure 1.
Strain-stress concrete experimental response in pseudo-static cyclic uniaxial compressive loading (adapted from [33]). Σand E are the homogeneous compression stress and strain in theconcrete test specimen that are measured in the loading direction.Σ and E are spatial mean quantities in the sense that Σ is com-puted as the load in the hydraulic cylinder of the testing machinedivided by the area of the specimen cross section, and E is com-puted as the displacement of the cylinder divided by the length ofthe concrete specimen.Most classical uniaxial constitutive models of concrete for numerical simulationdo not dissipate any energy in unloading-reloading cycles (see e.g. figure 2 [top left]).It is then common practice to add a viscous damping model (Rayleigh damping)to the inelastic structural model, to reproduce phenomena that are experimentallyobserved at the structural level (decreasing amplitude of displacements in free vi-bration). However, Rayleigh damping is well known to lack physical justification,even when care is taken to avoid generating spurious damping forces [7, 18].Another class of approach aims at reproducing more precisely the features ofFigure 1 through elaborate inelastic constitutive relations. Figure 2 shows typi-cal examples of uniaxial relations found in the literature. The relation describedin [9] [top left] defines different response phases for different strain intensities, withan additional coefficient to control the loss of stiffness. The constitutive relationdescribed in [23] [bottom right] comes from a formulation developed in the frame-work of thermodynamics with internal variables. It reproduces damping featuresreasonably well, in particular for higher amplitudes, but requires the identifica-tion of a rather large number of parameters. These first two types of relationsare somehow defined by parts for different loading regimes. They hence require a wider set of parameters and seem to contradict the seemingly smooth transition be-tween regimes observed experimentally. The constitutive relation described in [48][top right] is heuristically defined from a database of experiments. It reproducesunloading-reloading hysteresis mechanisms, but lacks a theoretical basis. Finally,the relation described in [32] [bottom left] is based on a physical model of damageand friction. It manages to dissipate energy in unloading-reloading cycles, but thelack of obvious physical meaning for some parameters can render their identificationdifficult. are somehow defined by parts for different loading regimes. They hence require awider set of parameters and seem to contradict the seemingly smooth transition be-tween regimes observed experimentally. The constitutive relation described in [48][top right] is heuristically defined from a database of experiments. It reproducesunloading-reloading hysteresis mechanisms, but lacks a theoretical basis. Finally,the relation described in [32] [bottom left] is based on a physical model of damageand friction. It manages to dissipate energy in unloading-reloading cycles, but thelack of obvious physical meaning for some parameters can render their identificationdifficult. Figure 2.
Typical strain-stress relations in pseudo-static cycliccompressive loading for different models: [9] [top left], [48] [topright], [32] [bottom left], and [23] [bottom right].The main purpose of this paper is to present a multi-scale stochastic nonlinearconcrete model that can be accommodated in an efficient structural frame element(fiber element), and that participates to the overall structural damping in dynamicloading. In particular, this implies the developed concrete model be capable of rep-resenting hysteresis loops in unloading-reloading cycles at macro-scale (the scalewhere such behavior as in Figure 1 can be observed). In this work, this is achievedby the introduction, at an underlying meso-scale, of spatial variability in the pa-rameters. At meso-scale, an elasto-plastic response with linear kinematic hardeningand heterogeneous yield stress is considered. This choice is mainly driven by itssimplicity and its relevance is illustrated in the numerical applications.The main issue with modeling the heterogeneity of the yield stress lies in the pa-rameterization. On the one hand, local information on the heterogeneity of concreteis available at a scale that we wish to avoid (because of the associated computa-tional costs). On the other hand, identification becomes extremely difficult whenvery fine models are considered. We therefore choose to model the heterogeneity ofthe yield stress through a stochastic model. Hence, only three parameters controlthat heterogeneity: a mean value, a variance, and a correlation length. The choice ofparameterizing the fluctuating field of constitutive parameters by statistical quan-tities means that there might be fluctuations in the quantities of interest measuredfor different realizations of the random model. However, as will become apparentin the examples, some sort of homogenization comes in and these fluctuations canrightfully be ignored.
Figure 2.
Typical strain-stress relations in pseudo-static cycliccompressive loading for different models: [9] [top left], [48] [topright], [32] [bottom left], and [23] [bottom right].The main purpose of this paper is to present a multi-scale stochastic nonlinearconcrete model that can be accommodated in an efficient structural frame element(fiber element), and that participates to the overall structural damping in dynamicloading. In particular, this implies the developed concrete model be capable of rep-resenting hysteresis loops in unloading-reloading cycles at macro-scale (the scalewhere such behavior as in Figure 1 can be observed). In this work, this is achievedby the introduction, at an underlying meso-scale, of spatial variability in the pa-rameters. At meso-scale, an elasto-plastic response with linear kinematic hardeningand heterogeneous yield stress is considered. This choice is mainly driven by itssimplicity and its relevance is illustrated in the numerical applications.The main issue with modeling the heterogeneity of the yield stress lies in the pa-rameterization. On the one hand, local information on the heterogeneity of concreteis available at a scale that we wish to avoid (because of the associated computa-tional costs). On the other hand, identification becomes extremely difficult whenvery fine models are considered. We therefore choose to model the heterogeneity ofthe yield stress through a stochastic model. Hence, only three parameters controlthat heterogeneity: a mean value, a variance, and a correlation length. The choice of parameterizing the fluctuating field of constitutive parameters by statistical quan-tities means that there might be fluctuations in the quantities of interest measuredfor different realizations of the random model. However, as will become apparentin the examples, some sort of homogenization comes in and these fluctuations canrightfully be ignored.Several authors in the literature have considered random models of fluctuatingnonlinear materials [16, 20, 1, 35], in particular for concrete [21, 34, 26, 49, 45,50, 30] or in the context of dynamic analysis [36, 28, 44]. We consider here amodeling framework that is a combination of ingredients found in several previouspapers [29, 5, 21, 30], with a fluctuating yield stress modeled as a random field withnon-zero correlation length. However, the objective in these papers was to assessthe influence of parameter uncertainty on some quantity of interest. An objectivewith the current paper is to observe the effect of randomness at a meso-scale on thenonlinear stress-strain relation at macro-scale. The work herein presented shouldtherefore be seen as an innovative proposal for parameterization of a nonlinearstress-strain relation.In Section 2, we recall the theoretical formulation of the inelastic beam modelthat will be used throughout this paper. The stochastic multi-scale constitutive re-lation developed to represent concrete cyclic behavior in reinforced concrete frameelements is introduced in section 3. Concrete behavior at macro-scale is retrievedfrom the description of a meso-scale where elasto-plastic response with linear kine-matic hardening and spatially variable yield stress is assumed. In particular, weemphasize in sections 3.3 and 3.4 the heterogeneity of the yield stress and the pa-rameterization of that heterogeneity through a random model. In section 4, wereport on the limiting case of vanishing correlation length and monotonic loading,for which several results can be derived analytically. Section 5 presents numericalapplications of the model in the context of dynamic structural analysis of reinforcedconcrete frame elements.2.
2D continuum Euler-Bernoulli inelastic beam
Classical displacement-based formulation has been retained here although othermixed formulations can in certain cases show better performances [46]. For thesake of conciseness, we present the beam element in the 2D case, extension to 3Dis straightforward.2.1.
Euler-Bernoulli kinematics.
We define the continuum beam B = { x ∈ R | x ∈ [0 , L ]; x ∈ [ − h/ , h/ x ∈ [ − w/ , w/ } . Such a beam has length L and uniform rectangular cross-section S of size w × h . We consider an orthonormalbasis ( i , i , i ) of R , so that any material point in space x = (cid:80) i =1 x i i i . In the 2Dsetting adopted here, Euler-Bernoulli kinematics can be written at any point x ∈ B and at any time t ∈ [0 , T ] as(1) u ( x , t ) = (cid:18) u ( x , t ) = u S ( x , t ) − x θ S ( x , t ) u ( x , t ) = u S ( x , t ) (cid:19) with θ S ( x , t ) = ∂u S ( x , t ) ∂x .u and u are the longitudinal and transversal components of the displacementvector u at any point in the beam. u S and u S are rigid body translations and θ S is rigid body rotation of section S at position x along the beam axis. Thus, u S , u S and θ S only depend on x .For small transformations, strain tensor reads E = (cid:0) D ( u ) + D T ( u ) (cid:1) , where D ( · ) = (cid:80) i =1 ∂ · ∂x i ⊗ i i , with ⊗ the tensor product and · T the transpose operation.Then, defining the axial strain (cid:15) S = ∂u S /∂x and the curvature χ S = ∂ u S /∂x ,it comes:(2) E ( x , t ) = E ( x , t ) i ⊗ i with(3) E ( x , t ) = (cid:15) S ( x , t ) − x χ S ( x , t ) . Variational formulation.
Suppose B is loaded with forces per unit lengthof beam b and concentrated forces applied at beam ends. A variational form of theproblem of finding u such that beam equilibrium is satisfied is: find u such that(4) 0 = (cid:90) L (cid:26)(cid:90) S (cid:0) δ(cid:15) S − x δχ S (cid:1) Σ d S (cid:27) dx − (cid:90) L δ u · b dx − δ Π bc , where δ u is any kinematically admissible displacement field Σ = Σ · i ⊗ i with Σ the stress tensor, · a matrix product here, and δ Π bc the potential for the forces atbeam ends.Introducing the normal and bending forces in the beam cross-sections as(5) N = (cid:90) S Σ d S and M = − (cid:90) S x Σ d S , equation (4) can then be rewritten as(6) 0 = (cid:90) L δ e S · q dx − (cid:90) L δ u · b dx − δ Π bc , where q = ( N, M ) T and e S = (cid:0) (cid:15) S , χ S (cid:1) T .2.3. Inelastic constitutive behavior.
Cross-section inelastic constitutive response q ( e S ; t ) is thereafter represented using uniaxial material constitutive response Σ( E ; x , t )integrated over the cross-section, rather than using a direct relation between sectiondisplacements and forces. This approach leads to what is often referred to as fiberbeam element.With ∆ denoting an increment of some quantity, we introduce the tangent mod-ulus D as ∆Σ = D × ∆ E , that is, from equation (3),(7) ∆Σ( x , t ) = D ( x , t ) (cid:0) ∆ (cid:15) S ( x , t ) − x ∆ χ S ( x , t ) (cid:1) . Introducing relation (7) in (5), beam section inelastic constitutive equation reads∆ q = K S ∆ e S with tangent stiffness matrix(8) K S ( x , t ) = (cid:20)(cid:90) S (cid:18) − x (cid:19) D ( x , t ) (cid:0) − x (cid:1) d S (cid:21) . Numerical implementation (structural level).
The finite element methodis used to approximate the displacement fields. Classically, we have for each ele-ment u ( x , t ) = N ( x ) d ( t ) and e S ( x , t ) = B ( x ) d ( t ), where the vector d gathers thedisplacements u S , u S , θ S at the element nodes, and matrices N and B gather theclassical shape functions for Euler-Bernoulli kinematics. Choosing, in equation (6), δ u = N δ d and δ e S = B δ d , and then linearizing the resulting relation, we have(9) K ( k ) n ∆ d ( k ) n = r ( k ) n , where K = (cid:82) L B T K S B dx and r = f − (cid:82) L B T q dx . f is the vector of nodalforces calculated from b and from any concentrated force applied at a beam ends.Subscript n refers to any time step in the loading history; superscript k refers toany Newton-Raphson iteration.For any integrable function g , line integrals are numerically approximated as (cid:82) L g ( x ) dx ≈ (cid:80) N l l =1 g ( x l ) W l where subscript l refers to a quadrature point and W l denotes quadrature weight and length. Section integrals are estimated as (cid:82) S l g ( x l ) d S l ≈ (cid:80) N F F =1 A F g ( x Fl ), where A F is the section area of the so-called fiber F and x Fl is the position of the fiber centroid in the control section S l at quadraturepoint l . 3. Multi-scale uniaxial cyclic model for concrete
In this section, we present the core objective of the paper, which is a nonlin-ear uniaxial constitutive model capable of representing salient features of concretecompressive response in cyclic loading. It is based on a simple local (meso-scale)elasto-plastic constitutive relation, for which the yield stress is modeled as a randomfield. The spatial fluctuations of the yield stress induce at macro-scale constitutiverelation Σ( E ) which resembles that encountered experimentally. Again, as alreadystated in the introduction, the idea of considering an elasto-plastic constitutive re-lation with fluctuating yield stress is not novel per se [29, 5, 21, 30]. It has beenproposed to assess effects of uncertain parameters on model outputs of interest forengineering practice while our aim here is to stress on the fact that this can be seenas a way of parameterizing material nonlinear constitutive relations. In particu-lar, it will be shown that, in some circumstances, even if randomness is present atmeso-scale, our model can predict non-random outputs.3.1. Meso- and macro-scale modeling of concrete.
Concrete is a heteroge-neous material (see figure 3). Two scales are classically considered for its modeling:(i) a micro-scale at which each phase (aggregate, concrete, cement paste) is clearlyidentified and modeled with its own constitutive relation; and (ii) a macro-scale atwhich concrete is considered as homogeneous. The macro-scale is the relevant scalefor structural engineering applications but the behavior at that scale is stronglyinfluenced by phenomena occurring at the micro-scale. In particular, the geome-try of the phases is important as it controls in a large part local concentrationsof stresses. As pointed out in the introduction to this paper, formulating and im-plementing concrete constitutive laws at the macro-scale can then turn out to bechallenging, even in the uniaxial case, and especially when it comes to accountingfor material energy dissipation sources. We follow here another path, consideringa meso-scale at which the parameters of the constitutive relation are assumed tovary continuously. This scale is intermediary between the macro-scale, at which the parameters are homogeneous, and the micro-scale, at which the parameters arediscontinuous. Figure 3.
Polished concrete section where the two phases humaneyes can see are represented: aggregates (crushed gravel and sand)and cement paste in-between.and stress at macro-scale Σ is computed as the spatial mean stress over R :(11) Σ( E ; θ ) = 1 |R| ! R σ ( ϵ ; θ ) d R , where |R| denotes the area of the fiber section R , σ and ϵ the stress and strain atmeso-scale in the uniaxial case. That is, analogously to E and Σ, ϵ = ϵ · i ⊗ i and σ = σ · i ⊗ i , where ϵ and σ are the strain and stress tensors at meso-scale suchthat the constitutive model is considered in a 1D setting. Besides, we enhance thefact that Σ is computed as the spatial mean of σ ( x ) and not as the sample mean of σ ( θ ). Local constitutive relation presented in section 3.2 below governs the relationbetween heterogeneous stress field σ ( x , t ) and strain field ϵ ( x , t ).Then, we introduce the tangent modulus at meso-scale D as ∆ σ = D × ∆ ϵ .According to equations (7), (10) and (11), we have the tangent modulus at macro-scale(12) D ( x , t ; θ ) = 1 |R| ! R D ( x , t ; θ ) d R . We will see in the examples below that, for a wide range of relative correlationlength and size of the section R , Σ and D do not depend on θ . In that case,even though the meso-scale model of concrete is stochastic, the resulting macro-scale model is deterministic, and independent of the actual realization of the localparameters that is being considered.Finally, we recall that concrete specimens exhibit quasi-brittle behavior in ten-sion with tensile strength generally 10 times smaller than compressive strength.Besides, we point out here that it is at macro-scale that these notions of compres-sion (Σ ≤
0) and tension (Σ >
0) are relevant for the modeling presented in thiswork.3.2.
Model of the uniaxial cyclic behavior at meso-scale.
We now concen-trate on the local uniaxial cyclic inelastic constitutive relation that will be consid-ered in this paper at meso-scale. We define it to be a simple elasto-plastic modelwith linear kinematic hardening, as illustrated in Figure 4. We provide here, inthe setting of computational inelasticity [41, 22], the assumptions and resultingequations corresponding to this relation:
Figure 3.
Polished concrete section where the two phases humaneyes can see are represented: aggregates (crushed gravel and sand)and cement paste in-between.For practical implementation, the heterogeneity will be conveyed in our modelby the fluctuations of a random field p ( x , θ ), where θ represents randomness. Con-sistently with the fiber beam formulation presented in the previous section, weconsider a mesh of fibers F spanning beam cross-sections S . These fibers have acentroid located at position x Fl and a cross-section denoted by R . In the spirit ofstrain-controlled tests on concrete specimens (see figure 1 along with its caption),strain E is assumed homogeneous over R :(10) (cid:15) ( x , t ) = E ( x Fl , t ) ∀ x ∈ R , and stress at macro-scale Σ is computed as the spatial mean stress over R :(11) Σ( E ; θ ) = 1 |R| (cid:90) R σ ( (cid:15) ; θ ) d R , where |R| denotes the area of the fiber section R , σ and (cid:15) the stress and strain atmeso-scale in the uniaxial case. That is, analogously to E and Σ, (cid:15) = (cid:15) · i ⊗ i and σ = σ · i ⊗ i , where (cid:15) and σ are the strain and stress tensors at meso-scale suchthat the constitutive model is considered in a 1D setting. Besides, we enhance thefact that Σ is computed as the spatial mean of σ ( x ) and not as the sample mean of σ ( θ ). Local constitutive relation presented in section 3.2 below governs the relationbetween heterogeneous stress field σ ( x , t ) and strain field (cid:15) ( x , t ).Then, we introduce the tangent modulus at meso-scale D as ∆ σ = D × ∆ (cid:15) .According to equations (7), (10) and (11), we have the tangent modulus at macro-scale(12) D ( x , t ; θ ) = 1 |R| (cid:90) R D ( x , t ; θ ) d R . We will see in the examples below that, for a wide range of relative correlationlength and size of the section R , Σ and D do not depend on θ . In that case, even though the meso-scale model of concrete is stochastic, the resulting macro-scale model is deterministic, and independent of the actual realization of the localparameters that is being considered.Finally, we recall that concrete specimens exhibit quasi-brittle behavior in ten-sion with tensile strength generally 10 times smaller than compressive strength.Besides, we point out here that it is at macro-scale that these notions of compres-sion (Σ ≤
0) and tension (Σ >
0) are relevant for the modeling presented in thiswork.3.2.
Model of the uniaxial cyclic behavior at meso-scale.
We now concen-trate on the local uniaxial cyclic inelastic constitutive relation that will be consid-ered in this paper at meso-scale. We define it to be a simple elasto-plastic modelwith linear kinematic hardening, as illustrated in Figure 4. We provide here, inthe setting of computational inelasticity [41, 22], the assumptions and resultingequations corresponding to this relation:(i) The total deformation (cid:15) is split into elastic ( (cid:15) e ) and plastic ( (cid:15) p ) parts:(13) (cid:15) = (cid:15) e + (cid:15) p . (ii) The following state equation holds (upper dot denotes derivative with re-spect to time):(14) ˙ σ = C ˙ (cid:15) e , where C is the elastic modulus.(iii) We impose that the stress σ corrected by α , the so-called back stress dueto kinematic hardening, satisfies yielding criterion(15) φ p = | σ + α | − σ y ≤ , where σ y ≥ φ p ( σ, α ) is negative,the material is elastic; otherwise, plasticity is activated and the materialstate evolves such that the condition φ p ( σ, α ) = 0 is satisfied.(iv) A change in (cid:15) p can only take place if φ p = 0 and yielding occurs in thedirection of σ + α , with a constant rate ˙ γ p ≥ (cid:15) p = (cid:26) ˙ γ p sign( σ + α ) if φ p ( σ, α ) = 00 otherwise . ˙ γ p is the so-called plastic multiplier.(v) With H the kinematic hardening modulus, the evolution of α is defined as:(17) ˙ α = − H ˙ (cid:15) p = − ˙ γ p H sign( σ + α ) . Accordingly, admissible stresses σ and α remain in the set K e = { ( σ, α ) | φ p ≤ } and two kinds of evolutions are possible:(i) If ( σ, α ) ∈ ¯ K e = { ( σ, α ) | φ p < } , the response is elastic:(18) ˙ (cid:15) p = 0 ⇒ ˙ σ = C ˙ (cid:15) . (ii) If ( σ, α ) ∈ ∂ K e = { ( σ, α ) | φ p = 0 } , any evolution is possible only if ˙ φ p = 0:(19) ∂φ p ∂σ ˙ σ + ∂φ p ∂α ˙ α = 0 ⇒ ˙ γ p = C sign( σ + α ) ˙ (cid:15)C + H ⇒ ˙ σ = CHC + H ˙ (cid:15) . (i) The total deformation ϵ is split into elastic ( ϵ e ) and plastic ( ϵ p ) parts:(13) ϵ = ϵ e + ϵ p . (ii) The following state equation holds (upper dot denotes derivative with re-spect to time):(14) ˙ σ = C ˙ ϵ e , where C is the elastic modulus.(iii) We impose that the stress σ corrected by α , the so-called back stress dueto kinematic hardening, satisfies yielding criterion(15) φ p = | σ + α | − σ y ≤ , where σ y ≥ φ p ( σ, α ) is negative,the material is elastic; otherwise, plasticity is activated and the materialstate evolves such that the condition φ p ( σ, α ) = 0 is satisfied.(iv) A change in ϵ p can only take place if φ p = 0 and yielding occurs in thedirection of σ + α , with a constant rate ˙ γ p ≥ ϵ p = ! ˙ γ p sign( σ + α ) if φ p ( σ, α ) = 00 otherwise . ˙ γ p is the so-called plastic multiplier.(v) With H the kinematic hardening modulus, the evolution of α is defined as:(17) ˙ α = − H ˙ ϵ p = − ˙ γ p H sign( σ + α ) . Figure 4.
Compressive cyclic behavior at meso-scale. The yieldstress σ y ( x ) fluctuates over the concrete section R : two local re-sponses at two distinct material points x [left] and x [right] arerepresented in the figure. During elastic loading/unloading, theslope is C ; in yielding phases, the slope is CHC + H < C . With thistype of model first yielding occurs once σ = σ y and the elastic do-main keeps constant amplitude 2 σ y . For low values of σ y , yieldingcan be observed in compression both during loading and unloading[right].Accordingly, admissible stresses σ and α remain in the set K e = { ( σ, α ) | φ p ≤ } and two kinds of evolutions are possible:(i) If ( σ, α ) ∈ ¯ K e = { ( σ, α ) | φ p < } , the response is elastic:(18) ˙ ϵ p = 0 ⇒ ˙ σ = C ˙ ϵ . (ii) If ( σ, α ) ∈ ∂ K e = { ( σ, α ) | φ p = 0 } , any evolution is possible only if ˙ φ p = 0:(19) ∂φ p ∂σ ˙ σ + ∂φ p ∂α ˙ α = 0 ⇒ ˙ γ p = C sign( σ + α ) ˙ ϵC + H ⇒ ˙ σ = C HC + H ˙ ϵ . Figure 4.
Compressive cyclic behavior at meso-scale. The yieldstress σ y ( x ) fluctuates over the concrete section R : two local re-sponses at two distinct material points x [left] and x [right] arerepresented in the figure. During elastic loading/unloading, theslope is C ; in yielding phases, the slope is CHC + H < C . With thistype of model first yielding occurs once σ = σ y and the elastic do-main keeps constant amplitude 2 σ y . For low values of σ y , yieldingcan be observed in compression both during loading and unloading[right].It is then possible, from equations (18) and (19), to give the expression of thetangent modulus D :(20) ˙ σ = D ˙ (cid:15) with D = (cid:26) C if ( σ, α ) ∈ ¯ K eCHC + H if ( σ, α ) ∈ ∂ K e . It should be reminded at this point that, as mentioned in the introduction andillustrated in Figure 4, the yield stress is assumed heterogeneous. The relationpresented here is therefore defined between stress and strain in each point in spacewith a different yield stress.3.3.
Description of the yield stress random field.
In this section, we describethe choice that is made for the modeling of the heterogeneous yield stress: the yieldstress is represented by a 2D log-normal homogeneous random field over the con-crete area R . We note here that, to the best of our knowledge, there currently existsno experimental dataset of local stress-strain uniaxial concrete responses recordedat many points over a concrete area. Here, the choice of using random fields to con-vey heterogeneity of the yield stress is mainly motivated by the effectiveness of themethod. We hope that this proposed interpretation of concrete meso-structure willfoster interaction between numerical and material scientists and help designing ex-perimental investigations that would eventually support or invalidate the numericalmodel we propose in this paper.Let us then consider a probability space (Θ , Ω , Pr), where Ω is a σ -algebra ofelements of Θ and Pr is a probability measure. The 2D random field of yield stressis constructed as a nonlinear point-wise transformation [17] S y ( x , θ ) = f ( G ( x ; θ ))of a homogeneous unit centered Gaussian random field G ( x ; θ ) with given powerspectral density (PSD) S GG ( κ ). The PSD is chosen here as the product of trianglefunctions with identical properties in the two orthogonal directions of the 2D plane, denoted by the subscript 1 and 2 throughout sections 3.3 and 3.4:(21) S GG ( κ ) = 1 κ u Λ (cid:18) κ κ u (cid:19) Λ (cid:18) κ κ u (cid:19) , where Λ( κ ) = 1 − | κ | if | κ | ≤ κ u, = κ u, = κ u . Forwave numbers above the cut-off κ u , the spectral density vanishes. In the spatialdomain, this PSD corresponds to the following autocorrelation function (the Fouriertransform of S GG ( κ )):(22) R GG ( ζ ) = sinc (cid:16) κ u π ζ (cid:17) sinc (cid:16) κ u π ζ (cid:17) , where sinc( x ) = sin( πx ) / ( πx ). The random field G ( x ; θ ) therefore fluctuates overtypical lengths (cid:96) c, = (cid:96) c, = (cid:96) c = 2 π/κ u , the so-called correlation length.The nonlinear point-wise transformation f controls the first-order marginal dis-tribution of σ y . In particular, it controls the desired expectation m and variance s of the yield stress homogeneous random field. In this paper, we choose to considera log-normal first-order marginal density, to ensure that the realizations are almost-surely and almost everywhere positive, as expected. The nonlinear transformationis then given by:(23) S y ( x , θ ) = exp( m G + s G × G ( x , θ )) > , where(24) m G = − ln (cid:32) m (cid:114) s m (cid:33) and s G = (cid:115) ln (cid:18) s m (cid:19) . Other first-order marginal densities could be considered, for example using themaximum entropy principle [37, 47, 42, 10, 11] or Bayesian identification [4, 19].Also, the PSD function is translated by the nonlinear transformation f so that thePSD of the yield stress and of the underlying Gaussian field are different, withpossible incompatibilities with the chosen first-order marginal density [17, 31, 38].These important but technical issues go beyond the scope of this paper and willnot be further discussed here. where(24) m G = − ln ! m " s m and s G = $ ln % s m & . Other first-order marginal densities could be considered, for example using themaximum entropy principle [37, 47, 42, 10, 11] or Bayesian identification [4, 19].Also, the PSD function is translated by the nonlinear transformation f so that thePSD of the yield stress and of the underlying Gaussian field are different, withpossible incompatibilities with the chosen first-order marginal density [17, 31, 38].These important but technical issues go beyond the scope of this paper and willnot be further discussed here. x /dx /dσ y /m Figure 5.
Realizations of log-normal random fields over a squareof size d for different correlation lengths: [left] ℓ c /d = 1, [center] ℓ c /d = 0 .
1, [right] ℓ c /d = 0 (white noise). Coordinates x and x , previously introduced in the description of the beam elementin section 2, are reused here to recall that the random fields aregenerated to parameter heterogeneous yield stress over beam cross-section areas.3.4. Numerical implementation (at each quadrature point in each beamfiber).
At each material point x Fl (quadrature point l , beam fiber F ), relations (10),(11) and (12) have to be calculated.On the one hand, Gaussian random field G is digitized using the spectral rep-resentation method, in its FFT implementation [39, 13, 40]. As an illustration,considering a 2D random field with identical properties in orthogonal directions 1and 2 (see [40] for more details):(25) G ( p ∆ x, p ∆ x ; θ ) = Re M − ’ n =0 M − ’ n =0 ( B n n ( θ ) exp ( iπ ( n p M + n p M )) + ˜ B n n ( θ ) exp ( iπ ( n p M − n p M ))) where i = −
1, ( p , p ) ∈ [0 , . . . , M − , B n n = 2∆ κ * S GG ( n ∆ κ, n ∆ κ ) exp( iφ n n ( θ ))˜ B n n = 2∆ κ * S GG ( n ∆ κ, − n ∆ κ ) exp( iψ n n ( θ ))(26)and ∆ κ = κ u /N ( N ∈ N ⋆ → ∞ ), M ≥ N , φ n n ( θ ) and ψ n n ( θ ) are independentrandom phase angles uniformly distributed in [0 , π ]. The resulting 2D random Figure 5.
Realizations of log-normal random fields over a squareof size d for different correlation lengths: [left] (cid:96) c /d = 1, [center] (cid:96) c /d = 0 .
1, [right] (cid:96) c /d = 0 (white noise). Coordinates x and x , previously introduced in the description of the beam elementin section 2, are reused here to recall that the random fields aregenerated to parameter heterogeneous yield stress over beam cross-section areas. Numerical implementation (at each quadrature point in each beamfiber).
At each material point x Fl (quadrature point l , beam fiber F ), relations (10),(11) and (12) have to be calculated.On the one hand, Gaussian random field G is digitized using the spectral rep-resentation method, in its FFT implementation [39, 13, 40]. As an illustration,considering a 2D random field with identical properties in orthogonal directions 1and 2 (see [40] for more details):(25) G ( p ∆ x, p ∆ x ; θ ) = Re M − (cid:88) n =0 M − (cid:88) n =0 (cid:16) B n n ( θ ) exp (cid:16) iπ (cid:16) n p M + n p M (cid:17)(cid:17) + ˜ B n n ( θ ) exp (cid:16) iπ (cid:16) n p M − n p M (cid:17)(cid:17)(cid:17) where i = −
1, ( p , p ) ∈ [0 , . . . , M − , B n n = 2∆ κ (cid:112) S GG ( n ∆ κ, n ∆ κ ) exp( iφ n n ( θ ))˜ B n n = 2∆ κ (cid:112) S GG ( n ∆ κ, − n ∆ κ ) exp( iψ n n ( θ ))(26)and ∆ κ = κ u /N ( N ∈ N (cid:63) → ∞ ), M ≥ N , φ n n ( θ ) and ψ n n ( θ ) are independentrandom phase angles uniformly distributed in [0 , π ]. The resulting 2D randomfield is periodic with a two-dimensional period L × L with L = M ∆ x = 2 π/ ∆ κ . Realizations of log-normal random fields with different correlation lengths areshown in figure 5.On the other hand, a generic concrete section R is built as a square with edgeof length d and R is meshed by a square grid of N f identical squares. Then, themesh size is d/N f and, for any integrable function g , (cid:82) R g ( x ) d R ≈ d N f (cid:80) N f f =1 g ( x f ),where x f is the position of the centroid of the f -th mesh over R .Then, digitized random field S y is mapped onto the x f ’s over R . To this purpose,we impose L ≥ d , that is |R| is smaller or equal to a period of the random field, andmapping is performed according to the following method. First, N f is calculatedas:(27) d = Int (cid:18) d ∆ x (cid:19) ∆ x + Res ⇒ N f = (cid:26) Int( d/ ∆ x ) if Res = 0Int( d/ ∆ x ) + 1 otherwise . Then, at the N f points x f ∈ R , S y ( x f ; θ ) is calculated as the linear interpolationof the four digitized values of S y ( x ; θ ) in ] x f − ∆ x, x f + ∆ x ] , as illustrated infigure 6.With the spatially variable yield stress now known at each point x f , f ∈ [1 , .., N f ]in R , the equations presented in section 3.2 can be solved to update the variables atmeso-scale. This is done numerically at each of the N f positions following classicalreturn-mapping computational procedure [41, 22].Finally, as a transition from compression to tension is detected during globalNewton-Raphson iterative process to solve structural equilibrium equations, that isΣ ( k +1) n > ( k ) n ≤
0, a local Newton-Raphson precess is implemented to findthe strain E c for which Σ ( k +1) n ( E c ) = 0, to update the meso-structure accordingly,and to set Σ ( k +1) n = 0 and D ( k +1) n = 0.Before observing on numerical tests the shape of the stress-strain curves obtainedwith this model, we turn to the simple case of vanishing correlation length ( (cid:96) c → field is periodic with a two-dimensional period L × L with L = M ∆ x = 2 π/ ∆ κ . Realizations of log-normal random fields with different correlation lengths areshown in figure 5.On the other hand, a generic concrete section R is built as a square with edgeof length d and R is meshed by a square grid of N f identical squares. Then, themesh size is d/N f and, for any integrable function g , ! R g ( x ) d R ≈ d N f " N f f =1 g ( x f ),where x f is the position of the centroid of the f -th mesh over R .Then, digitized random field S y is mapped onto the x f ’s over R . To this purpose,we impose L ≥ d , that is |R| is smaller or equal to a period of the random field, andmapping is performed according to the following method. First, N f is calculatedas:(27) d = Int d ∆ x $ ∆ x + Res ⇒ N f = % Int( d/ ∆ x ) if Res = 0Int( d/ ∆ x ) + 1 otherwise . Then, at the N f points x f ∈ R , S y ( x f ; θ ) is calculated as the linear interpolationof the four digitized values of S y ( x ; θ ) in ] x f − ∆ x, x f + ∆ x ] , as illustrated infigure 6. Figure 6.
On the one hand, the random field is digitized on asquare grid (dashed lines) of ( M + 1) points — here M = 6 —spanning an area of size L × L , with L = M × ∆ x . On the otherhand, generic concrete square section R has an area of size d × d ( d ≤ L ) that is divided into N f identical meshes (plain lines) —here N f = 5. The values of the random field at the positions x f occupied by the centroids of the N f meshes ( × ) is calculated asthe linear interpolation of the four surrounding digitized values ofthe random field ( ◦ ).With the spatially variable yield stress now known at each point x f , f ∈ [1 , .., N f ]in R , the equations presented in section 3.2 can be solved to update the variables atmeso-scale. This is done numerically at each of the N f positions following classicalreturn-mapping computational procedure [41, 22].Finally, as a transition from compression to tension is detected during globalNewton-Raphson iterative process to solve structural equilibrium equations, that isΣ ( k +1) n > ( k ) n ≤
0, a local Newton-Raphson precess is implemented to findthe strain E c for which Σ ( k +1) n ( E c ) = 0, to update the meso-structure accordingly,and to set Σ ( k +1) n = 0 and D ( k +1) n = 0. Figure 6.
On the one hand, the random field is digitized on asquare grid (dashed lines) of ( M + 1) points — here M = 6 —spanning an area of size L × L , with L = M × ∆ x . On the otherhand, generic concrete square section R has an area of size d × d ( d ≤ L ) that is divided into N f identical meshes (plain lines) —here N f = 5. The values of the random field at the positions x f occupied by the centroids of the N f meshes ( × ) is calculated asthe linear interpolation of the four surrounding digitized values ofthe random field ( ◦ ).derived, so that discussion is more straightforward. The more general case withfinite correlation length will be considered later in Section 5.1.4. A particular case: vanishing correlation length and monotonicloading
Preliminaries.
The case of vanishing correlation length along with uniaxialcyclic loading has been treated in a general setting. Indeed in [24], the stress-strain uniaxial response is given as a probability density function ( pdf ) of stresswith respect to the time-dependent strain and a second-order exact expression ofthe pdf evolution is computed solving the Fokker-Planck-Kolmogorov equation thatgoverns the problem. This latter method is valid for monotonic as well as cyclicloading. Hereafter, the validity of the results is limited to monotonic loading, butthe problem is cast in a different and simpler mathematical setting that can besolved analytically. These analytical developments shed light on some capabilitiesof the model introduced in the previous section and that will be retrieved in themore general case of non-zero correlation in Section 5.1.4.2.
Constitutive response at macro-scale.
Let respectively denote R e and R p the shares of a fiber cross-section that remain elastic and yield. Accordingto the developments in section 3.2: R e ( t ; θ ) = { x ∈ R | D ( x , t ; θ ) = C } and R p ( t ; θ ) = { x ∈ R | D ( x , t ; θ ) = CH/ ( C + H ) } . Also, R e ∩R p = ∅ and R = R e ∪R p .Note that R e and R p are time-dependent because D depends on the loading history.We denote by |•| the area of • . Then, considering a subset A of R , we have, ∀ x ∈ R ,the probability measure Pr[ x ∈ A ] = |A| / |R| . Using the fact that |R| = |R e | + |R p | , we first rewrite the tangent modulus atmacro-scale in equation (12) as:(28) D = 1 |R| (cid:18) |R e | C + |R p | CHC + H (cid:19) = CC + H (cid:18) |R e ||R| C + H (cid:19) . We now seek an explicit expression for |R e | / |R| .First, suppose the state of the material is known at time t , then we define thetrial stresses(29) σ tr ( x , t ) = σ ( x ) + C ( (cid:15) ( t ) − (cid:15) ) and α tr ( x , t ) = α ( x ) , where subscript 0 refers to time t . In the particular case of monotonic loading, anecessary and sufficient condition for x to be in R p at time t > t is φ p,tr ( x , t ) ≥ σ y ( x ) ≤ | σ tr ( x , t ) + α ( x ) | (see equation 15). We then have:(30) |R e | / |R| = Pr[ x ∈ R e ] = 1 − Pr[ S y ( x ) ≤ | σ tr ( x , t ) + α ( x ) | ] . Then, in the particular case of vanishing correlation length, the random variables S y ( x ) are independent and identically distributed over R . For the log-normaldistribution assumption made throughout this work, it means that the cumulativedensity function of S y ( x ) is, ∀ x ∈ R :(31) F S y ( x ) ( σ y ) = Pr[ S y ( x ) ≤ σ y ] = 12 (cid:18) (cid:18) ln σ y − m G √ s G (cid:19)(cid:19) , where erf is the so-called error function.Finally, for the sake of simplicity and without any loss of generality, we assume σ = α = (cid:15) = 0. Accordingly, and using equations (29), along with (24) toreplace m G and s G by the mean m and standard deviation s of the homogeneouslog-normal random field S y , it comes:(32) |R e ||R| = 1 − F S y ( x ) ( | C(cid:15) ( t ) | ) = 12 − erf ln (cid:18) | C(cid:15) ( t ) | m (cid:113) s m (cid:19)(cid:113) (cid:0) s m (cid:1) . Equations (28) and (32) are used to plot figure 7 where the response of the modelat macro-scale is shown for different sets of mean and variance parameters for thelog-normal random yield stress field S y .4.3. Asymptotic response of the model at macro-scale.
The following as-ymptotic behaviors can be observed at macro-scale:(i) Suppose s /m approaches 0. Then, according to equation (32), |R e | / |R| approaches the Heaviside’s function H ( m − | CE ( t ) | ), that is |R e | / |R| = 0 if | CE ( t ) | > m and |R e | / |R| = 1 if | CE ( t ) | ≤ m . According to equation (28),the model response at macro-scale is then as follows:(33) ˙Σ( t ) = D ( t ) ˙ E ( t ) where D = (cid:26) C if | CE ( t ) | ≤ m CHC + H if | CE ( t ) | > m . (ii) If s /m → ∞ , then |R e | / |R| → t ) → CH/ ( C + H ) ˙ E ( t ). Finally, for the sake of simplicity and without any loss of generality, we assume σ = α = ϵ = 0. Accordingly, and using equations (29), along with (24) toreplace m G and s G by the mean m and standard deviation s of the homogeneouslog-normal random field S y , it comes:(32) |R e ||R| = 1 − F S y ( x ) ( | Cϵ ( t ) | ) = 12 ⎛⎜⎜⎝ − erf ⎛⎜⎜⎝ ln $ | Cϵ ( t ) | m % s m &% ’ s m ( ⎞⎟⎟⎠⎞⎟⎟⎠ . Equations (28) and (32) are used to plot figure 7 where the response of the modelat macro-scale is shown for different sets of mean and variance parameters for thelog-normal random yield stress field S y . − − − − − − − − − − − − − − − − − − E Σ [ M P a ] increasing s − − − − − − − − − − − − − − − − − − − E Σ [ M P a ] increasing m Figure 7.
Monotonic macro-scale response of the model withvanishing correlation length. Parameters C = 30 GPa and H = 10 GPa are used. [top] m = 30 MPa and s/m = 10 − , . , , , ; [bottom] s = 30 MPa and s/m =10 − , . , , . , . Plain curves show asymptotic behaviors as s/m → ∞ .4.3. Asymptotic response of the model at macro-scale.
The following as-ymptotic behaviors can be observed at macro-scale:(i) Suppose s /m approaches 0. Then, according to equation (32), |R e | / |R| approaches the Heaviside’s function H ( m − | CE ( t ) | ), that is |R e | / |R| = 0 if Finally, for the sake of simplicity and without any loss of generality, we assume σ = α = ϵ = 0. Accordingly, and using equations (29), along with (24) toreplace m G and s G by the mean m and standard deviation s of the homogeneouslog-normal random field S y , it comes:(32) |R e ||R| = 1 − F S y ( x ) ( | Cϵ ( t ) | ) = 12 ⎛⎜⎜⎝ − erf ⎛⎜⎜⎝ ln $ | Cϵ ( t ) | m % s m &% ’ s m ( ⎞⎟⎟⎠⎞⎟⎟⎠ . Equations (28) and (32) are used to plot figure 7 where the response of the modelat macro-scale is shown for different sets of mean and variance parameters for thelog-normal random yield stress field S y . − − − − − − − − − − − − − − − − − − E Σ [ M P a ] increasing s − − − − − − − − − − − − − − − − − − − E Σ [ M P a ] increasing m Figure 7.
Monotonic macro-scale response of the model withvanishing correlation length. Parameters C = 30 GPa and H = 10 GPa are used. [top] m = 30 MPa and s/m = 10 − , . , , , ; [bottom] s = 30 MPa and s/m =10 − , . , , . , . Plain curves show asymptotic behaviors as s/m → ∞ .4.3. Asymptotic response of the model at macro-scale.
The following as-ymptotic behaviors can be observed at macro-scale:(i) Suppose s /m approaches 0. Then, according to equation (32), |R e | / |R| approaches the Heaviside’s function H ( m − | CE ( t ) | ), that is |R e | / |R| = 0 if Figure 7.
Monotonic macro-scale response of the model withvanishing correlation length. Parameters C = 30 GPa and H = 10 GPa are used. [top] m = 30 MPa and s/m = 10 − , . , , , ; [bottom] s = 30 MPa and s/m =10 − , . , , . , . Plain curves show asymptotic behaviors as s/m → ∞ .(iii) Now with finite and non-zero s /m :(34) ˙Σ( t ) = D ( t ) ˙ E ( t ) where D → (cid:26) C if E ( t ) → CHC + H if E ( t ) → ∞ . These asymptotic responses at macro-scale are illustrated in figure 7 (plain lines).5.
Numerical applications
Concrete uniaxial compressive cyclic response at macro-scale.
Firstnumerical applications aim at demonstrating the capability of the model introducedabove in section 3 to represent the response of concrete in uniaxial compressivecyclic loading. Five model parameters need to be considered: elastic and harden-ing moduli C and H , along with mean m , standard deviation s and correlationlength (cid:96) c used to build realizations of a homogeneous log-normal random field thatparameterizes the fluctuations of the yield stress σ y over beam sections.The effects of m , s and (cid:96) c on the material response at macro-scale will be furtherinvestigated below. Right now however, we set: • C = 27 . • H = 0 according to both (i) the fact that H controls tangent modulus atmacro scale as strain becomes large (see equation (34)), and (ii) that weseek a numerical response that ultimately exhibits null tangent modulusin monotonic loading at macro-scale. We anticipate here stressing thatthe model developed in previous sections is not capable of representingthe softening phase as strain increases while stress decreases (non-positivetangent modulus). Figure 8.
Sample mean (thick plain line) plus/minus standarddeviation (boundaries of the shaded areas) monotonic response atmacro-scale computed from a sample of 100 realizations of thematerial structure at meso-scale with m = 30 MPa and s/m = 1for the log-normal marginal law. Meso-structures are generatedwith different correlation lengths: [left] ℓ c /d = 0 .
1, [center] ℓ c /d = 0 .
2, [right] ℓ c /d = 0 .
4. Cyclic response for one particularrealization of the meso-structure is also shown (thin plain line).increases. The area R defined with ℓ c /d = 0 . N = 10 and M = 64, is statisticallyrepresentative in the sense that there is almost independence between the randomrealization of the meso-structure and the response at macro-scale. We remark herethat the model is capable of representing variability from one concrete sample toanother and that this variability can bring information on the correlations in themeso-structure. Suppose indeed that we had 100 concrete samples and a variabilityof the responses at macro-scale close to that shown by the grey area in figure 8[center] for instance. Then, the meso-structure of the tested concrete would be bestrepresented by the ratio ℓ c /d = 0 .
2. Consequently, the sample standard deviationcan bring information on the actual correlation length.Finally, we can notice that the shape of the cyclic response (thin line) representsmost of the salient features exhibited experimentally in uniaxial compression testfor concrete (remember figure 1). We point out here that strength degradation(softening) along with stiffness degradation (observed experimentally one cycle afteranother) are not represented by this model. However, a key point for representingmaterial damping is the capability of the model to generate local hysteresis loopsin unloading-loading cycles.5.1.2.
Influence of m and s on the macroscopic response. We use here, beside N =10 and M = 64 ( N f = 64), ℓ c /d = 0 . m and standard deviation s of the log-normal randomfield that conveys spatial variability in the material structure at meso-scale.It can be observed that for a small value of s , response approaches bi-linear elasto-plastic behavior (actually perfectly plastic because H is set to zero here) where thereis almost no hysteresis observed during unloading-loading cycle (figure 9 [left] and[right]). This comes from the fact that, if s approaches 0, there is almost no spatialvariability of the yield stress because it is almost homogeneous over R and takes Figure 8.
Sample mean (thick plain line) plus/minus standarddeviation (boundaries of the shaded areas) monotonic response atmacro-scale computed from a sample of 100 realizations of thematerial structure at meso-scale with m = 30 MPa and s/m = 1for the log-normal marginal law. Meso-structures are generatedwith different correlation lengths: [left] (cid:96) c /d = 0 .
1, [center] (cid:96) c /d = 0 .
2, [right] (cid:96) c /d = 0 .
4. Cyclic response for one particularrealization of the meso-structure is also shown (thin plain line).5.1.1.
Influence of (cid:96) c on the macroscopic response. We first illustrate how (cid:96) c in-fluences the macroscopic response by considering the three following cases: (i) (cid:96) c /d = 0 .
1, (ii) (cid:96) c /d = 0 . (cid:96) c /d = 0 .
4. For each of these three cases, wetake N = 10 and M = 64, that is ∆ x/d = 0 .
016 and N f = 64 (see equation (27)).We recall that the random field characteristics are taken as identical in both spacedirections ( (cid:96) c = (cid:96) c, = (cid:96) c, , N = N = N , . . . ). Besides, a sample of 100 indepen-dent homogeneous log-normal random fields with targeted mean m = 30 MPa andcoefficient of variation s/m = 1 for the marginal log-normal law is generated foreach case.Resulting material responses at macro-scale are shown in figure 8. A first ob-vious observation is that model response at macro-scale is much richer than atmeso-scale (see figure 4). We can then observe that sample mean response (thickline) is not sensitive to the correlation length. However the variability of the macro-scopic response from one realization of the meso-structure to another depends onthe correlation length: it is almost null for (cid:96) c /d = 0 . (cid:96) c increases. The area R defined with (cid:96) c /d = 0 . N = 10 and M = 64, is statisticallyrepresentative in the sense that there is almost independence between the randomrealization of the meso-structure and the response at macro-scale. We remark herethat the model is capable of representing variability from one concrete sample toanother and that this variability can bring information on the correlations in themeso-structure. Suppose indeed that we had 100 concrete samples and a variabilityof the responses at macro-scale close to that shown by the grey area in figure 8[center] for instance. Then, the meso-structure of the tested concrete would be bestrepresented by the ratio (cid:96) c /d = 0 .
2. Consequently, the sample standard deviationcan bring information on the actual correlation length.Finally, we can notice that the shape of the cyclic response (thin line) representsmost of the salient features exhibited experimentally in uniaxial compression testfor concrete (remember figure 1). We point out here that strength degradation(softening) along with stiffness degradation (observed experimentally one cycle after another) are not represented by this model. However, a key point for representingmaterial damping is the capability of the model to generate local hysteresis loopsin unloading-loading cycles. − − − − − − − − − − − − − − − − − − E Σ / m increasing s − − − − − − − − − − − − − − − − − − E Σ / m increasing m − − − − − − − − − − − − − − − − − − E Σ / m increasing m = s Figure 9.
Sample mean (plain line) plus/minus standard devia-tion (dashed lines) response at macro-scale computed from a sam-ple of 100 different realizations of the material structure at meso-scale. Meso-structures are generated with different targeted mean m and coefficients of variation s/m for the log-normal marginallaw: [left] m = 30 MPa and s/m = 0 . , ,
3; [center] s = 30 MPaand s/m = 0 . , ,
3; [right] s/m = 1 and m = 10 , ,
50 MPa.values close to its mean m ; then the response at macro-scale coincides with thatat meso-scale (elasto-plasticity with H = 0 here). This is also in accordance withwhat was shown already in figure 7.Responses shown in figure 9 [left] lie in-between this latter extreme case andthe other extreme case of s approaching infinity. In this situation, log-normaldistribution approaches 0 all over the positive real semi-line and, consequently,plasticity is activated almost everywhere over R resulting in a macro-scale responsethat is perfectly plastic without elastic phase (that is here Σ = 0 for any E because H = 0).Also, it is shown in figure 9 that the value of the strain E at which stress Σreaches zero when unloading (residual plastic deformation) is much more sensitiveto parameter m than s . Furthermore, the thickness of the hysteresis loops obviouslydepends on the so-called coefficient of variation s/m but it is not clear whether itis more sensitive to either of the two parameters. Finally, the variability in thesample of responses at macro-scale increases with s/m and is more sensitive to m than s , at least as far as the range of values chosen here for both parameters isconcerned.5.2. Damping in a reinforced concrete column in free vibration.
We nowshow how the material model developed in the previous sections can be used torepresent the experimental backbone curve of a concrete test specimen in uniaxialloading. Then, we implement this material law in the fiber beam element presentedin section 2 and show how damping is generated in a reinforced concrete (RC)column in free vibration. The observed damping does not result from the additionof damping forces in the balance equation of the RC column but from the hysteresisloops in the concrete material law at macro-scale.5.2.1.
Geometry of the column and loading.
The column considered here corre-sponds to the 1st-floor external column of the ductile ( R = 4) RC frame testedin [15, 25]. The loading is however here different: it consists of a mass M = 500 kgimposed step by step and kept constant while the column oscillates in free vibration Figure 9.
Sample mean (plain line) plus/minus standard devia-tion (dashed lines) response at macro-scale computed from a sam-ple of 100 different realizations of the material structure at meso-scale. Meso-structures are generated with different targeted mean m and coefficients of variation s/m for the log-normal marginallaw: [left] m = 30 MPa and s/m = 0 . , ,
3; [center] s = 30 MPaand s/m = 0 . , ,
3; [right] s/m = 1 and m = 10 , ,
50 MPa.5.1.2.
Influence of m and s on the macroscopic response. We use here, beside N =10 and M = 64 ( N f = 64), (cid:96) c /d = 0 . m and standard deviation s of the log-normal randomfield that conveys spatial variability in the material structure at meso-scale.It can be observed that for a small value of s , response approaches bi-linear elasto-plastic behavior (actually perfectly plastic because H is set to zero here) where thereis almost no hysteresis observed during unloading-loading cycle (figure 9 [left] and[right]). This comes from the fact that, if s approaches 0, there is almost no spatialvariability of the yield stress because it is almost homogeneous over R and takesvalues close to its mean m ; then the response at macro-scale coincides with thatat meso-scale (elasto-plasticity with H = 0 here). This is also in accordance withwhat was shown already in figure 7.Responses shown in figure 9 [left] lie in-between this latter extreme case andthe other extreme case of s approaching infinity. In this situation, log-normaldistribution approaches 0 all over the positive real semi-line and, consequently,plasticity is activated almost everywhere over R resulting in a macro-scale responsethat is perfectly plastic without elastic phase (that is here Σ = 0 for any E because H = 0).Also, it is shown in figure 9 that the value of the strain E at which stress Σreaches zero when unloading (residual plastic deformation) is much more sensitiveto parameter m than s . Furthermore, the thickness of the hysteresis loops obviouslydepends on the so-called coefficient of variation s/m but it is not clear whether itis more sensitive to either of the two parameters. Finally, the variability in thesample of responses at macro-scale increases with s/m and is more sensitive to m than s , at least as far as the range of values chosen here for both parameters isconcerned.5.2. Damping in a reinforced concrete column in free vibration.
We nowshow how the material model developed in the previous sections can be used torepresent the experimental backbone curve of a concrete test specimen in uniaxialloading. Then, we implement this material law in the fiber beam element presentedin section 2 and show how damping is generated in a reinforced concrete (RC)column in free vibration. The observed damping does not result from the additionof damping forces in the balance equation of the RC column but from the hysteresisloops in the concrete material law at macro-scale.5.2.1.
Geometry of the column and loading.
The column considered here corre-sponds to the 1st-floor external column of the ductile ( R = 4) RC frame testedin [15, 25]. The loading is however here different: it consists of a mass M = 500 kgimposed step by step and kept constant while the column oscillates in free vibrationconsequently to a horizontal force F ( t ). The geometrical and loading characteristicsof the column are depicted in figure 10. consequently to a horizontal force F ( t ). The geometrical and loading characteristicsof the column are depicted in figure 10. Figure 10.
Geometry and loading of the column.5.2.2.
Concrete constitutive model.
In [25], the monotonic uniaxial response (back-bone curve) of the concrete cast to build the RC column is detailed. It is usedhere as the baseline to identify the parameters of the concrete model. Accordingto this report, we set C = 27 . H = 0; then we use N = 10, M = 64( N f = 64) and ℓ c /d = 0 . m and s are identified.Figure 11 shows model response macro-scale (plain line) with m = 30 . s/m = 0 . − − − − − − − − − − − − − − E Σ [ M P a ] Figure 11.
Sample mean (—) response at macro-scale obtainednumerically from a sample of 2,000 meso-structures, along withbackbone curve (- -) recorded during uniaxial test on a specimen ofthe concrete used to build the RC column of interest here. Targetedmean and standard deviation of the marginal log-normal law are m = 30 . s/m = 0 . Figure 10.
Geometry and loading of the column.5.2.2.
Concrete constitutive model.
In [25], the monotonic uniaxial response (back-bone curve) of the concrete cast to build the RC column is detailed. It is usedhere as the baseline to identify the parameters of the concrete model. Accordingto this report, we set C = 27 . H = 0; then we use N = 10, M = 64( N f = 64) and (cid:96) c /d = 0 . m and s are identified.Figure 11 shows model response macro-scale (plain line) with m = 30 . s/m = 0 . Steel cyclic model.
Young modulus C s = 224 . y = 438MPa and ultimate stress Σ u = 601 MPa have been experimentally measured during consequently to a horizontal force F ( t ). The geometrical and loading characteristicsof the column are depicted in figure 10. Figure 10.
Geometry and loading of the column.5.2.2.
Concrete constitutive model.
In [25], the monotonic uniaxial response (back-bone curve) of the concrete cast to build the RC column is detailed. It is usedhere as the baseline to identify the parameters of the concrete model. Accordingto this report, we set C = 27 . H = 0; then we use N = 10, M = 64( N f = 64) and ℓ c /d = 0 . m and s are identified.Figure 11 shows model response macro-scale (plain line) with m = 30 . s/m = 0 . − − − − − − − − − − − − − − E Σ [ M P a ] Figure 11.
Sample mean (—) response at macro-scale obtainednumerically from a sample of 2,000 meso-structures, along withbackbone curve (- -) recorded during uniaxial test on a specimen ofthe concrete used to build the RC column of interest here. Targetedmean and standard deviation of the marginal log-normal law are m = 30 . s/m = 0 . Figure 11.
Sample mean (—) response at macro-scale obtainednumerically from a sample of 2,000 meso-structures, along withbackbone curve (- -) recorded during uniaxial test on a specimen ofthe concrete used to build the RC column of interest here. Targetedmean and standard deviation of the marginal log-normal law are m = 30 . s/m = 0 . Steel cyclic model.
Young modulus C s = 224 . y = 438MPa and ultimate stress Σ u = 601 MPa have been experimentally measured duringuniaxial tests on longitudinal steel rebars [25]. An elasto-plastic model with kine-matic hardening is used to represent steel response in cyclic loading. The modelimplemented with these latter measured parameters is shown in figure 12. − − − − − − − − − − E Σ [ M P a ] Figure 12.
Numerical cyclic response of a steel longitudinal rebarused to build the frame. Cyclic behavior has not been observedexperimentally.5.2.4.
Free vibration – Structural damping.
Those concrete and steel uniaxial con-stitutive models are implemented in the fiber frame element presented in section 2.The column is modeled with one frame element with N l = 2 control sections and N F = 6 fibers (actually layers here in the case of a 2D problem). As already men-tioned in section 5.2.1, the mass M = 500 kg is imposed step by step and keptconstant while the column oscillates in free vibration consequently to the horizon-tal force F ( t ). The column possibly exhibits nonlinear response while the mass M is applied and while F ( t ) increases from 0 to F . Figure 13 shows typical columntop-displacement time histories for two different values of F . It can be observedthat damping depends on the amplitude of the oscillations: the column is clearlydamped for F = 15 kN (grey curve) while damping is much lower for F = 5 kN(black curve). One can also notice the different vibration periods for both hori-zontal forces; this is due to the fact that the larger force activates some nonlinearmechanisms in the structure, which leads to an elongation of the structural vibra-tion period. We finally stress again here that there is no damping force addedin the dynamic balance equations, such as for instance Rayleigh damping forces:the damping effect shown in figure 13 only comes from the hysteresis loops in theconcrete response during unloading-reloading cycles.We now define what we will refer to as “viscous-like damping ratio” and hereafterdenote by ξ v . Considering the column top-displacement time history X top ( t ) in freevibration, we appeal to the so-called log-decrement method to evaluate the modaldamping ratio ξ v (see e.g. [8, § ξ v = 12 πN c ln X peaktop ( t N ) X peaktop ( t N ) X peaktop ( t N ) and X peaktop ( t N ) are the amplitudes of any two peaks separated by N c = N − N cycles. It is worth recalling here that this is only valid in case damping is Figure 12.
Numerical cyclic response of a steel longitudinal rebarused to build the frame. Cyclic behavior has not been observedexperimentally.5.2.4.
Free vibration – Structural damping.
Those concrete and steel uniaxial con-stitutive models are implemented in the fiber frame element presented in section 2.The column is modeled with one frame element with N l = 2 control sections and N F = 6 fibers (actually layers here in the case of a 2D problem). As already men-tioned in section 5.2.1, the mass M = 500 kg is imposed step by step and keptconstant while the column oscillates in free vibration consequently to the horizon-tal force F ( t ). The column possibly exhibits nonlinear response while the mass M is applied and while F ( t ) increases from 0 to F . Figure 13 shows typical column top-displacement time histories for two different values of F . It can be observedthat damping depends on the amplitude of the oscillations: the column is clearlydamped for F = 15 kN (grey curve) while damping is much lower for F = 5 kN(black curve). One can also notice the different vibration periods for both hori-zontal forces; this is due to the fact that the larger force activates some nonlinearmechanisms in the structure, which leads to an elongation of the structural vibra-tion period. We finally stress again here that there is no damping force addedin the dynamic balance equations, such as for instance Rayleigh damping forces:the damping effect shown in figure 13 only comes from the hysteresis loops in theconcrete response during unloading-reloading cycles. − − − time [s] t o p - d i s p l a ce m e n t [ m ] Figure 13.
Top displacement time history in free vibration for F = 5 kN (black) and F = 15 kN (grey). Mass M and horizon-tal forces F are applied step-by-step during the first and secondseconds, then horizontal force F abruptly drops to zero and thecolumn oscillates in free vibration.linear viscous, which in our case is not necessarily the case. Indeed, equations (35)comes from the assumption that the envelope of the decaying top-displacement isdescribed as X top ( t ) = X e − ξπf t with f the modal frequency. Hence the terms“viscous-like” to characterize the calculated damping ratios. Number of cycles N in free vibration V i s c o u s - li k e d a m p i n g r a t i o ξ v [ % ] Figure 14.
Viscous-like damping ratio time history ξ v for F =5 kN (black) and F = 15 kN (grey).Based on the top-displacement time histories in figure 13, figure 14 shows how ξ v decreases throughout free vibration time history for both values of F . Viscous-likedamping ratios ξ v ( t N ), are computed according to equations (35), with N c = 5.Note that such damping ratios depend on the parameters m and s/m of the randomfield along with the hardening parameter H at meso-scale. For the sake of illus-tration, figure 15 shows other results for another set of material parameters that isnot optimal for representing the monotonic response in compression of the concreteused to build the tested column. The capability of the proposed material modelfor generating structural damping has been demonstrated and the development ofan automatic procedure for identifying the full set of parameters targeting accuraterepresentation of both cyclic concrete response and damping is left for future work.6. Conclusions
In this paper, a multi-scale stochastic uniaxial cyclic model suitable for rep-resenting most of the salient features of concrete nonlinear response observed in
Figure 13.
Top displacement time history in free vibration for F = 5 kN (black) and F = 15 kN (grey). Mass M and horizon-tal forces F are applied step-by-step during the first and secondseconds, then horizontal force F abruptly drops to zero and thecolumn oscillates in free vibration.We now define what we will refer to as “viscous-like damping ratio” and hereafterdenote by ξ v . Considering the column top-displacement time history X top ( t ) in freevibration, we appeal to the so-called log-decrement method to evaluate the modaldamping ratio ξ v (see e.g. [8, § ξ v = 12 πN c ln X peaktop ( t N ) X peaktop ( t N ) X peaktop ( t N ) and X peaktop ( t N ) are the amplitudes of any two peaks separated by N c = N − N cycles. It is worth recalling here that this is only valid in case damping islinear viscous, which in our case is not necessarily the case. Indeed, equations (35)comes from the assumption that the envelope of the decaying top-displacement isdescribed as X top ( t ) = X e − ξπft with f the modal frequency. Hence the terms“viscous-like” to characterize the calculated damping ratios.Based on the top-displacement time histories in figure 13, figure 14 shows how ξ v decreases throughout free vibration time history for both values of F . Viscous-likedamping ratios ξ v ( t N ), are computed according to equations (35), with N c = 5.Note that such damping ratios depend on the parameters m and s/m of the randomfield along with the hardening parameter H at meso-scale. For the sake of illus-tration, figure 15 shows other results for another set of material parameters that is − − − time [s] t o p - d i s p l a ce m e n t [ m ] Figure 13.
Top displacement time history in free vibration for F = 5 kN (black) and F = 15 kN (grey). Mass M and horizon-tal forces F are applied step-by-step during the first and secondseconds, then horizontal force F abruptly drops to zero and thecolumn oscillates in free vibration.linear viscous, which in our case is not necessarily the case. Indeed, equations (35)comes from the assumption that the envelope of the decaying top-displacement isdescribed as X top ( t ) = X e − ξπft with f the modal frequency. Hence the terms“viscous-like” to characterize the calculated damping ratios. Number of cycles N in free vibration V i s c o u s - li k e d a m p i n g r a t i o ξ v [ % ] Figure 14.
Viscous-like damping ratio time history ξ v for F =5 kN (black) and F = 15 kN (grey).Based on the top-displacement time histories in figure 13, figure 14 shows how ξ v decreases throughout free vibration time history for both values of F . Viscous-likedamping ratios ξ v ( t N ), are computed according to equations (35), with N c = 5.Note that such damping ratios depend on the parameters m and s/m of the randomfield along with the hardening parameter H at meso-scale. For the sake of illus-tration, figure 15 shows other results for another set of material parameters that isnot optimal for representing the monotonic response in compression of the concreteused to build the tested column. The capability of the proposed material modelfor generating structural damping has been demonstrated and the development ofan automatic procedure for identifying the full set of parameters targeting accuraterepresentation of both cyclic concrete response and damping is left for future work.6. Conclusions
In this paper, a multi-scale stochastic uniaxial cyclic model suitable for rep-resenting most of the salient features of concrete nonlinear response observed in
Figure 14.
Viscous-like damping ratio time history ξ v for F =5 kN (black) and F = 15 kN (grey).not optimal for representing the monotonic response in compression of the concreteused to build the tested column. The capability of the proposed material modelfor generating structural damping has been demonstrated and the development ofan automatic procedure for identifying the full set of parameters targeting accuraterepresentation of both cyclic concrete response and damping is left for future work. − − − time [s] t o p - d i s p l a ce m e n t [ m ] Number of cycles N in free vibration V i s c o u s - li k e d a m p i n g r a t i o ξ v [ % ] Figure 15. [top] Top-displacement time history in free vibrationand [bottom] viscous-like damping ratio time history ξ v for a set ofmaterial parameters that is not optimal for the column consideredabove: m = 20 MPa, s/m = 10, C = 27 . H = 10 GPa.compressive experimental tests has been developed. It is based on the constructionof a meso-scale where the response at each material point is elasto-plastic withkinematic hardening and heterogeneous yield stress. This implies that the transi-tion from elastic to plastic regime occurs at a loading level that is different in eachmaterial point. Heterogeneity is parameterized by a 2D homogeneous log-normalrandom field. As a first illustration of the capabilities of the model, some analyt-ical results are derived in the particular case of monotonic loading and vanishingcorrelation length for the random field. Then, numerical simulations are performedand the effects of the parameters of the random field – that is the mean m , coef-ficient of variation s/m and correlation length ℓ c – are investigated. It is shownthat for small values of the correlation length, material response at macro-scaledoes not depend on the realization of the random field, showing that the devel-oped model is suitable for an objective representation of the material behavior.Besides, it is shown that the mean m and standard deviation s can be identifiedso that the monotonic compressive response of an actual concrete test specimencan be accurately represented by the developed model. The developed model how-ever lacks the ingredients for representing both strength and stiffness degradationmechanisms. Finally, the developed material model is implemented in a frame ele-ment in the purpose of representing the dynamic response of an actual reinforcedconcrete column. The numerical analysis of the column in free vibration shows thecapability of the developed material model to create patterns classically associated − − − time [s] t o p - d i s p l a ce m e n t [ m ] Number of cycles N in free vibration V i s c o u s - li k e d a m p i n g r a t i o ξ v [ % ] Figure 15. [top] Top-displacement time history in free vibrationand [bottom] viscous-like damping ratio time history ξ v for a set ofmaterial parameters that is not optimal for the column consideredabove: m = 20 MPa, s/m = 10, C = 27 . H = 10 GPa.compressive experimental tests has been developed. It is based on the constructionof a meso-scale where the response at each material point is elasto-plastic withkinematic hardening and heterogeneous yield stress. This implies that the transi-tion from elastic to plastic regime occurs at a loading level that is different in eachmaterial point. Heterogeneity is parameterized by a 2D homogeneous log-normalrandom field. As a first illustration of the capabilities of the model, some analyt-ical results are derived in the particular case of monotonic loading and vanishingcorrelation length for the random field. Then, numerical simulations are performedand the effects of the parameters of the random field – that is the mean m , coef-ficient of variation s/m and correlation length ℓ c – are investigated. It is shownthat for small values of the correlation length, material response at macro-scaledoes not depend on the realization of the random field, showing that the devel-oped model is suitable for an objective representation of the material behavior.Besides, it is shown that the mean m and standard deviation s can be identifiedso that the monotonic compressive response of an actual concrete test specimencan be accurately represented by the developed model. The developed model how-ever lacks the ingredients for representing both strength and stiffness degradationmechanisms. Finally, the developed material model is implemented in a frame ele-ment in the purpose of representing the dynamic response of an actual reinforcedconcrete column. The numerical analysis of the column in free vibration shows thecapability of the developed material model to create patterns classically associated Figure 15. [top] Top-displacement time history in free vibrationand [bottom] viscous-like damping ratio time history ξ v for a set ofmaterial parameters that is not optimal for the column consideredabove: m = 20 MPa, s/m = 10, C = 27 . H = 10 GPa.6. Conclusions
In this paper, a multi-scale stochastic uniaxial cyclic model suitable for rep-resenting most of the salient features of concrete nonlinear response observed incompressive experimental tests has been developed. It is based on the constructionof a meso-scale where the response at each material point is elasto-plastic withkinematic hardening and heterogeneous yield stress. This implies that the transi-tion from elastic to plastic regime occurs at a loading level that is different in eachmaterial point. Heterogeneity is parameterized by a 2D homogeneous log-normalrandom field. As a first illustration of the capabilities of the model, some analyt-ical results are derived in the particular case of monotonic loading and vanishing correlation length for the random field. Then, numerical simulations are performedand the effects of the parameters of the random field – that is the mean m , coef-ficient of variation s/m and correlation length (cid:96) c – are investigated. It is shownthat for small values of the correlation length, material response at macro-scaledoes not depend on the realization of the random field, showing that the devel-oped model is suitable for an objective representation of the material behavior.Besides, it is shown that the mean m and standard deviation s can be identifiedso that the monotonic compressive response of an actual concrete test specimencan be accurately represented by the developed model. The developed model how-ever lacks the ingredients for representing both strength and stiffness degradationmechanisms. Finally, the developed material model is implemented in a frame ele-ment in the purpose of representing the dynamic response of an actual reinforcedconcrete column. The numerical analysis of the column in free vibration shows thecapability of the developed material model to create patterns classically associatedto damping effects. In this simulation, damping does no come from some dampingforces added in the dynamic balance equation (e.g. Rayleigh damping) but fromthe multi-scale stochastic nonlinear model. Although the underlying model is sto-chastic, the simulations and results shown are the same for any realization of thestochastic model.The main research prospects lie (i) in the enhancement of the model at meso-scale so that it can represent stiffness and strength degradation mechanisms atmacro-scale; (ii) in the precise characterization of the stochastic model based oninformation from lower scales. This will consist in choosing, based on rationalarguments, the type of first-order marginal law and correlation model, as well asthe value of the corresponding parameters (mean, variance and correlation length).Although in another context, such an interaction between structural and materialscientists has already been appealed for in [6]. Also, these issues could be consideredin the context of stochastic micro-meso scale transition [43, 3, 12]. Acknowledgement
The first author is supported by a Marie Curie International Outgoing Fellow-ship within the 7th European Community Framework Programme (proposal No.275928). The second author, working within the SINAPS@ project, benefited fromFrench state funding managed by the National Research Agency under programRNSR Future Investments bearing reference No. ANR-11-RSNR-0022-04.
References [1] M. Anders and M. Hori. Three-dimensional stochastic finite element method for elasto-plasticbodies.
Int. J. Numer. Meth. Engr. , 51(4):449–478, 2001.[2] Applied Technology Council. Modeling and acceptance criteria for seismic design and anal-ysis of tall buildings. Technical Report PEER/ATC-72-1, Pacific Earthquake EngineeringResearch Center, Richmond (CA), October 2010.[3] M Arnst and R Ghanem. Probabilistic equivalence and stochastic model reduction in mul-tiscale analysis.
Computational Methods in Applied Mechanics and Engineering , 197:3584–3592, 2008.[4] J. L. Beck and L. S. Katafygiotis. Updating models and their uncertainties. Part I: bayesianstatistical framework.
J. Engr. Mech. ASCE , 124(4):455–461, 1998.[5] C. E. Brenner and C. Bucher. A contribution to the SFE-based reliability assessment ofnonlinear structures under dynamic loading.
Prob. Engr. Mech. , 10(4):265–273, 1995. [6] D C Charmpis, G I Schu¨eller, and M F Pellissetti. The need for linking micromechanics ofmaterials with stochastic finite elements: A challenge for materials science. ComputationalMaterials Science , 41:27–37, 2007.[7] Finley A Charney. Unintended consequences of modeling damping in structures.
Journal ofStructural Engineering , 134(4):581–592, 2008.[8] Ray W Clough and Joseph Penzien.
Dynamics of structures . McGraw-Hill, Inc., 1975.[9] Computers & Structures Inc. (CSI).
Perform3D User’s manual . Berkeley, CA, USA, 2007.[10] R. Cottereau, D. Clouteau, and C. Soize. Construction of a probabilistic model for impedancematrices.
Comp. Meth. Appl. Mech. Engr. , 196(17-20):2252–2268, 2007.[11] R. Cottereau, D. Clouteau, and C. Soize. Probabilistic impedance of foundation: impact onthe seismic design on uncertain soils.
Earth. Engr. Struct. Dyn. , 37(6):899–918, 2008.[12] R´egis Cottereau, Didier Clouteau, Hachmi Ben Dhia, and C´edric Zaccardi. A stochastic-deterministic coupling method for continuum mechanics.
Computer Methods in Applied Me-chanics and Engineering , 200:3280–3288, 2011.[13] G. Deodatis. Non-stationary stochastic vector processes: seismic ground motion applications.
Prob. Engr. Mech. , 11:149–168, 1996.[14] Norberto Dom´ınguez and Adnan Ibrahimbegovic. A non-linear thermodynamical model forsteel–concrete bonding.
Computers and Structures , 106–107:29–45, 2012.[15] Andr´e Filiatrault, ´Eric Lachapelle, and Patrick Lamontagne. Seismic performance of ductileand nominally ductile reinforced concrete moment resisting frames. I. Experimental study.
Canadian Journal of Civil Engineering , 25:331–341, 1998.[16] G. N. Frantziskonis. Stochastic modeling of heterogeneous materials – a process for the anal-ysis and evaluation of alternative formulations.
Mech. Mater. , 27(3):165–175, 1998.[17] M. Grigoriu. Simulation of stationary non-gaussian translation processes.
J. Engr. Mech.ASCE , 124(2):121–126, 1998.[18] J F Hall. Problems encountered from the use (or misuse) of Rayleigh damping.
EarthquakeEngineering and Structural Dynamics , 35:525–545, 2006.[19] C. Howson and P. Urbach.
Scientific reasoning. The Bayesian approach . Open Court Pub-lishing Company, 3rd edition, 2005.[20] C. Huet. An integrated micromechanics and statistical continuum thermodynamics approachfor studying the fracture behaviour of microcracked heterogeneous materials with delayedresponse.
Engr. Fracture Mech. , 58(5-6):459–463 465–556, 1997.[21] J. Huh and A. Haldar. Stochastic finite-element-based seismic risk of nonlinear structures.
J.Struct. Engr. ASCE , 127(3):323–329, 2001.[22] Adnan Ibrahimbegovic.
Nonlinear solid mechanics: Theoretical formulations and finite ele-ment solution methods . Springer, 2009.[23] Pierre Jehel, Luc Davenne, Adnan Ibrahimbegovic, and Pierre L´eger. Towards robustviscoelastic-plastic-damage material model with different hardenings / softenings capableof representing salient phenomena in seismic loading applications.
Computers and Concrete ,7(4):365–386, 2010.[24] B. Jeremi´c, K. Sett, and M. L. Kavvas. Probabilistic elasto-plasticity: formulation in 1D.
Acta Geotechnica , 2(3):197–210, 2007.[25] Patrick Lamontagne. Comportement sismique d’une ossature ductile en b´eton arm´e (R=4).M.A.Sc. thesis (in French), ´Ecole Polytechnique de Montr´eal, Montreal, QC, Canada, April1997.[26] Tae-Hyung Lee and Khalid M Mosalam. Probabilistic fiber element modeling of reinforcedconcrete structures.
Computers and Structures , 82:2285–2299, 2004.[27] Pierino Lestuzzi and Hugo Bachmann. Displacement ductility and energy assessment fromshaking table tests on RC structural walls.
Engineering Structures , 29:1708–1721, 2007.[28] J. Li and J. B. Chen. Probability density evolution method for dynamic response analysis ofstructures with uncertain parameters.
Comp. Mech. , 34(5):400–409, 2004.[29] W. K. Liu, T. Belytschko, and A. Mani. Probabilistic finite elements for nonlinear structuraldynamics.
Comp. Meth. Appl. Mech. Engr. , 56(1):61–81, 1986.[30] T. Namikawa and J. Koseki. Effects of spatial correlation on the compression behavior ofcement-treated column.
J. Geotech. Geoenviron. Engr. , 13(8):1346–1359, 2013.[31] B. Puig and J.-L. Akian. Non-gaussian simulation using Hermite polynomials expansion andmaximum entropy principle.
Prob. Engr. Mech. , 19(4):293–305, 2004. [32] Fr´ed´eric Ragueneau, Christian La Borderie, and Jacky Mazars. Damage model for concrete-like materials coupling cracking and friction, contribution towards structural damping: firstuniaxial applications. Mechanics of cohesive-frictional materials , 5:607–625, 2000.[33] S Ramtani.
Contribution to the modeling of the multi-axial behavior of damaged concretewith description of the unilateral characteristics . PhD Thesis (in French), Paris 6 University,1990.[34] Mauro de Vasconcellos Real, Am´erico Campos Filho, and S´ergio Roberto Maestrini. Re-sponse variability in reinforced concrete structures with uncertain geometrical and materialproperties.
Nuclear Engineering and Design , 226:205–220, 2003.[35] Bojana Rosi´c and Hermann G Matthies. Computational approaches to inelastic media withuncertain parameters.
Journal of the Serbian Society for Computational Mechanics , 2(1):28–43, 2008.[36] G. I. Schu¨eller and H. J. Pradlwarter. On the stochastic response of nonlinear FE models.
Arch. Appl. Mech. , 69(9-10):765–784, 1999.[37] C E Shannon. A mathematical theory of communication.
The Bell System Technical Journal ,27(3):379—423, 1948.[38] M D Shields and George Deodatis. A simple and efficient methodology to approximate ageneral non-Gaussian stationary stochastic vector process by a translation process with ap-plications in wind velocity simulation.
Probabilistic Engineering Mechanics , 31:19–29, 2013.[39] Masanobu Shinozuka and George Deodatis. Simulation of stochastic processes by spectralrepresentation.
Applied Mechanics Reviews , 44(4):191—203, 1991.[40] Masanobu Shinozuka and George Deodatis. Simulation of multidimensional gaussian stochas-tic fields by spectral representation.
Appl. Mech. Rev. , 49(1):29–53, 1996.[41] J C Simo and T J R Hughes.
Computational Inelasticity . Springer, Berlin, 1998.[42] C. Soize. A nonparametric model of random uncertainties for reduced matrix models instructural dynamics.
Prob. Engr. Mech. , 15:277–294, 2000.[43] Christian Soize. Tensor-valued random fields for meso-scale stochastic model of anisotropicelastic microstructure and probabilistic analysis of representative volume element size.
Prob-abilistic Engineering Mechanics , 23:307–323, 2008.[44] George Stefanou and Michalis Fragiadakis. Nonlinear dynamic analysis of frames with sto-chastic non-gaussian material properties.
Engineering Structures , 31:1841–1850, 2009.[45] P. Stroeven, J. Hu, and H. Chen. Stochastic heterogeneity as fundamental basis for the designand evaluation of experiments.
Cement Concrete Composites , 30(6):506–514, 2008.[46] R L Taylor, F C Filippou, A Saritas, and F Auricchio. A mixed finite element method forbeam and frame problems.
Computational Mechanics , 31:192–203, 2003.[47] F. E. Udwadia. Some results on maximum entropy distributions for parameters known to liein finite intervals.
SIAM Rev. , 31(1):103–109, 1989.[48] P S Wong and F J Vecchio.
VecTor2 & Formworks User’s Manuals . University of Toronto,Department of Civil Engineering, Toronto, ON, Canada, 2002.[49] P. Wriggers and S. O. Moftah. Mesoscale models for concrete: homogenisation and damage.
Finite Element Anal. Design , 42(7):623–636, 2006.[50] Z. J. Yang, X. T. Su, J. F. Chen, and G. H. Liu. Monte Carlo simulation of complex cohesivefracture in random heterogeneous quasi-brittle materials.