U.S. Department of Transportation
Federal Highway Administration
1200 New Jersey Avenue, SE
Washington, DC 20590
2023664000
Federal Highway Administration Research and Technology
Coordinating, Developing, and Delivering Highway Transportation Innovations
This report is an archived publication and may contain dated technical, contact, and link information 

Publication Number: FHWAHRT05062
Date: May 2007 

Users Manual for LSDYNA Concrete Material Model 159PDF Version (1.49 KB)
PDF files can be viewed with the Acrobat® Reader® Chapter 3. Users ManualThis chapter is intended to be a brief Users Manual for those who want to run the model with a cursory, rather than indepth, understanding of the underlying theory and equations. This chapter begins with a description of the LSDYNA concrete model input, followed by a brief theoretical description of the model theory. All information in this chapter is included in the LSDYNA Keyword Users Manual.^{(6)} LSDYNA Input*MAT_CSCM {OPTION} This is material type 159. This is a smooth or continuous surface cap model and is available for solid elements in LSDYNA. The user has the option of inputting his or her own material properties (<BLANK> option), or requesting default material properties for normal strength concrete (CONCRETE). Options include: <BLANK> CONCRETE In selecting an option the keyword cards appear: *MAT_CSCM *MAT_CSCM _CONCRETE Define the following card for all options. Card Format
Define the following card for the CONCRETE option. Do not define for the <BLANK> option.
Define the following cards for the <BLANK> option. Do not define for CONCRETE.
Define for all options.
Define for CONCRETE option.
Remarks: Default concrete input parameters are for normal strength concrete with unconfined compression strengths between about 20 and 58 MPa (2,901 and 8,412 psi). Define for <BLANK> option only.
Model Formulation and Input ParametersThis is a cap model with a smooth intersection between the shear yield surface and hardening cap, as shown in Figure 84. The initial damage surface coincides with the yield surface. Rate effects are modeled with viscoplasticity. Figure 84. Illustration. General shape of the concrete model yield surface in two dimensions. Stress Invariants. The yield surface is formulated in terms of three stress invariants: J_{1}, the first invariant of the stress tensor; J¢_{2}, the second invariant of the deviatoric stress tensor; and J¢_{3}, the third invariant of the deviatoric stress tensor. The invariants are defined in terms of the deviatoric stress tensor, S_{ij} and pressure, P, as shown in Figure 85: Figure 85. Equation. Three stress invariants, J_{1}, J′_{2}, J′_{3}. Plasticity Surface. The three invariant yield function is based on these three invariants, and the cap hardening parameter, κ, as shown in Figure 86: Figure 86. Equation. Plasticity yield function f. Here F_{f} is the shear failure surface, F_{c} is the hardening cap, and ℜ is the Rubin threeinvariant reduction factor. The cap hardening parameter κ is the value of the pressure invariant at the intersection of the cap and shear surfaces. Trial elastic stress invariants are temporarily updated via the trial elastic stress tensor, σ^{T}. These are denoted J_{1}^{T}, J¢_{2}^{T}, and J¢_{3}^{T}. Elastic stress states are modeled when f (J_{1}^{T}, J¢_{2}^{T}, J¢_{3}^{T}, κ^{T} ) < 0. Elasticplastic stress states are modeled when f (J_{1}^{T}, J¢_{2}^{T}, J¢_{3}^{T}, κ^{T} ) > 0. In this case, the plasticity algorithm returns the stress state to the yield surface in a manner that f (J_{1}^{P}, J¢_{2}^{P}, J¢_{3}^{P}, κ^{P}) = 0. Shear Failure Surface. The strength of concrete is modeled by the shear surface in the tensile and low confining pressure regimes (see Figure 87): Figure 87. Equation. Shear surface function F_{f}. Here the values of are selected by fitting the model surface to strength measurements from TXC tests conducted on plain concrete cylinders. Rubin Scaling Function. Concrete fails at lower values of J¢_{2}(principal stress difference) for TXE and TOR tests than it does for TXC tests conducted at the same pressure. The Rubin scaling function ℜ determines the strength of concrete for any state of stress relative to the strength for TXE, via ℜF_{f}. Strength in TOR is modeled as Q_{1}F_{f} . Strength in TXE is modeled as Q_{2}F_{f}, in Figure 88, where: Figure 88. Equation. Most general form for scaling functions Q_{1}, Q_{2}. Cap Hardening Surface. The strength of concrete is modeled by a combination of the cap and shear surfaces in the low to high confining pressure regimes. The cap is used to model plastic volume change related to pore collapse (although the pores are not explicitly modeled). The isotropic hardening cap is a twopart function that is either unity or an ellipse, as shown in Figure 89: Figure 89. Equation. Cap surface function, F_{c}. where L(κ) is defined in Figure 90: Figure 90. Equation. Definition of L of kappa. The equation for F_{c} is equal to unity for J_{1} £ L(κ). It describes the ellipse for J_{1} > L(κ). The intersection of the shear surface and the cap is at J_{1} = κ. κ_{0} is the value of J_{1} at the initial intersection of the cap and shear surfaces before hardening is engaged (before the cap moves). The equation for L(κ) restrains the cap from retracting past its initial location at κ_{0}. The intersection of the cap with the J_{1} axis is at J_{1} = X(κ). This intersection depends on the cap ellipticity ratio R, where R is the ratio of its major to minor axes (see Figure 91): Figure 91. Equation. Pressure invariant X as a function of kappa. The cap moves to simulate plastic volume change. The cap expands (X(κ) and κincrease) to simulate plastic volume compaction. The cap contracts (X(κ) and κ decrease) to simulate plastic volume expansion, called dilation. The motion (expansion and contraction) of the cap is based on the hardening rule, as shown in Figure 92: Figure 92. Equation. Plastic volume strain hardening rule, ε ^{p}_{v}. Here ε ^{ρ}_{v} is the plastic volume strain, W is the maximum plastic volume strain, and D_{1} and D_{2} are model input parameters. X_{0} is the initial location of the cap when κ= κ_{0}. The five input parameters (X_{0}, W, D_{1}, D_{2}, and R) are obtained from fits to the pressurevolumetric strain curves in isotropic compression and uniaxial strain. X_{0} determines the pressure at which compaction initiates in isotropic compression. R, combined with X_{0}, determines the pressure at which compaction initiates in uniaxial strain. D_{1} and D_{2} determine the shape of the pressurevolumetric strain curves. W determines the maximum plastic volume compaction. Shear Hardening Surface. In unconfined compression, the stressstrain behavior of concrete exhibits nonlinearity and dilation before the peak. This type of behavior is modeled with an initial shear yield surface, N_{H}F_{f}, which hardens until it coincides with the ultimate shear yield surface, F_{f}. Two input parameters are required. One parameter, N_{H}, initiates hardening by setting the location of the initial yield surface. A second parameter, C_{H}, determines the rate of hardening (amount of nonlinearity). Damage. Concrete exhibits softening in the tensile and low to moderate compressive regimes (see Figure 93): Figure 93. Equation. Transformation of viscoplastic stress to damaged stress, σ^{d}_{ij}. A scalar damage parameter, d, transforms the viscoplastic stress tensor without damage, denoted s ^{vp}, into the stress tensor with damage, denoted s ^{d}. Damage accumulation is based on two distinct formulations, which are called brittle damage and ductile damage. The initial damage threshold is coincident with the shear plasticity surface; thus, the user does not have to specify the threshold. Ductile Damage. Ductile damage accumulates when the pressure (P) is compressive and an energytype term, τ_{d}, exceeds the damage threshold, r_{0d}. Ductile damage accumulation depends on the total strain components, ε_{i}_{j}, as shown in Figure 94: Figure 94. Equation. Ductile damage accumulation, τ _{d}. The stress components s_{ij} are the elastoplastic stresses (with kinematic hardening) calculated before application of damage and rate effects. Brittle Damage. Brittle damage accumulates when the pressure is tensile and an energytype term, τ_{b}, exceeds the damage threshold, r_{0b}. Brittle damage accumulation depends on the maximum principal strain, ε_{max}, as shown in Figure 95: Figure 95. Equation. Brittle damage accumulation, τ_{b}. Softening Function. As damage accumulates, the damage parameter d increases from an initial value of zero, towards a maximum value of 1, via Figures 96 and 97: Figure 96. Equation. Brittle damage, d of t_{b}. Figure 97 Equation. Ductile damage, d of τ_{d}. The damage parameter applied to the six stresses is equal to the current maximum of the brittle or ductile damage parameter. The parameters A and B or C and D set the shape of the softening curve plotted as stressdisplacement or stressstrain. The parameter dmax is the maximum damage level that can be attained. It is calculated internally and is less than 1 at moderate confining pressures. The compressive softening parameter, A, may also be reduced with confinement, using the input parameter pmod, as shown in Figure 98: Figure 98. Equation. Reduction of A with confinement. Regulating Mesh Size Sensitivity. The concrete model maintains constant fracture energy, regardless of element size. The fracture energy is defined here as the area under the stressdisplacement curve from peak strength to zero strength. Constant fracture energy is achieved by internally formulating the softening parameters A and C in terms of the element length, L (cube root of the element volume), the fracture energy, G_{f}, the initial damage threshold, τ_{0t} or τ_{0c} , and the softening shape parameters, D or B. The fracture energy is calculated from up to five userspecified input parameters (G_{fc}, G_{ft}, G_{fs}, pwrc, pwrc). The user specifies three distinct fracture energy values. These are the fracture energy in uniaxial tensile stress, G_{ft}, pure shear stress, G_{fs}, and uniaxial compressive stress, G_{fc.} The model internally selects the fracture energy from equations, which interpolate between the three fracture energy values as a function of the stress state (expressed via two stress invariants). The interpolation equations depend on the userspecified input powers pwrc and pwrt, as shown in Figure 99: Figure 99. Equation. Brittle and ductile damage thresholds, G_{f}. The internal parameter trans is limited to range between 0 and 1. Element Erosion. An element loses all strength and stiffness as d®1. To prevent computational difficulties with very low stiffness, element erosion is available as a user option. An element erodes when d > 0.99 and the maximum principal strain is greater than a usersupplied input value, 1ERODE. Viscoplastic Rate Effects. At each time step, the viscoplastic algorithm interpolates between the elastic trial stress, , and the inviscid stress (without rate effects), , to set the viscoplastic stress (with rate effects), : as shown in Figure 100: Figure 100. Equation. Viscoplastic stress, σ^{vp}_{ij}. This interpolation depends on the effective fluidity coefficient, η, and the time step, Δt. The effective fluidity coefficient is internally calculated from five usersupplied input parameters and interpolation equations shown in Figure 101: Figure 101. Equation. Variation of the fluidity parameter η in tension and compression. The input parameters are η_{0t} and N_{t} for fitting uniaxial tensile stress data, η_{0c} and N_{c} for fitting the uniaxial compressive stress data, and Srate for fitting shear stress data. The effective strain rate is shown in Figure 102: Figure 102. Definition of effective strain rate. This viscoplastic model may predict substantial rate effects at high strain rates (>100). To limit rate effects at high strain rates, the user may input overstress limits in tension (overt) and compression (overc). These input parameters limit calculation of the fluidity parameter, as shown in Figure 103: Figure 103. Equation. Overstress limit of η. where over = overt when the pressure is tensile, and over = overc when the pressure is compressive. The user has the option of increasing the fracture energy as a function of effective strain rate via the repow input parameter, as shown in Figure 104: Figure 104. Equation. Fracture energy with rate effects, . Here is the fracture energy enhanced by rate effects, and f¢ is the yield strength before application of rate effects (which is calculated internally by the model). The term in brackets is greater than, or equal to 1, and is the approximate ratio of the dynamic to static strength. 