Elasto‑plastic constitutive law for solid elements at constant temperature (non-associated) with linear anisotropic elasticity. Isotropic hardening/softening of friction angle and cohesion is possible.
This law is a combination of ORTHO3D (for the elastic part of the law) and PLA3D (for the plastic behaviour).
This law is only used for mechanical analysis of elasto-plastic anisotropic porous media undergoing large strains.
Prepro: LORTHOPLA.F
See PLASOL law
See PLASOL law (above)
The material cohesion depends on the angle $\alpha_{\sigma_1}$ between the major compressive principal stress $\vec{\sigma_1}$ and the normal to the bedding plane $\vec{e_3}$ : \[ \alpha_{\sigma_1} = \arccos\left( \frac{\vec{e_3}\vec{\sigma_{1}}'}{||\vec{e_3}||\ ||\vec{\sigma_1}'||} \right) \]
Fig. 1: Schematic view of the angle between the normal to bedding plane and the
direction of major principal stress
Three cohesion values are defined ($c_{0^{\circ}}, c_{min}, c_{90^{\circ}}$), for major principal stress parallel $\alpha_{\sigma_1} = 0^{\circ}$ (perpendicular), perpendicular $\alpha_{\sigma_1} = 90^{\circ}$ (parallel) and with an angle of $\alpha_{\sigma_1, min}$ with respect to the normal to bedding plane (with respect to the bedding plane) (Salehnia, 2015)1). Between those values, cohesion varies linearly with $\alpha_{\sigma_1}$. The mathematical expression of the cohesion is as follows: \[ c = \max \left[\left( \frac{c_{min} - c_{0^{\circ}}}{\alpha_{\sigma_1, min}} \right)\alpha_{\sigma_1} + c_{0^{\circ}} ; \left( \frac{c_{90^{\circ}} - c_{min}}{90^{\circ} - \alpha_{\sigma_1, min}} \right)\left( \alpha_{\sigma_1} - \alpha_{\sigma_1, min} \right)+ c_{0^{\circ}} \right] \]
Fig. 2: Schematic view of the cohesion evolution as a function of the angle between the
normal to bedding plane and the direction of major principal stress
The material cohesion anisotropy is defined with a microstructure tensor $a_{ij}$, which is a measure of the material fabric. The eigenvectors of this tensor correspond to the preferred or principal material axes (orthotropic axes) $e_1, e_2, e_3$. The cohesion corresponds to the projection of this tensor on a generalized loading vector l , therefore the cohesion specifies the effect of load orientation relative to the material axes. \[c = a_{ij}l_il_j\] \[l_i = \sqrt{\frac{\sigma_{l1}^2 + \sigma_{l2}^2 +\sigma_{l3}^2}{\sigma_{ij}\sigma_{ij}}}\]
Where $\sigma_{ij}$ expressed in reference to the material axes. The cohesion can be expressed as : \[c = c_0\left( 1 + A_{ij}l_il_j\right)\] \[A_{ij} = \frac{dev(a_{ij})}{c_0} =\frac{a_{ij}}{c_0} - \delta_{ij}\] \[c_0 = \frac{a_{kk}}{3}\]
where $A_{ij}$ is a symmetric traceless operator, $A_{kk}=0$. The above expression can be generalized by considering higher order tensors :
\[c= c_0 \left( 1+A_{ij}l_il_j + b_1(A_{ij}l_il_j)^2 + b_2(A_{ij}l_il_j)^3 + … \right)\]
where $B_1,b_2,…$ are constants.
Considering cross-anisotropy, i.e. transverse isotropy, and refering the problem to the principal material axes implies $A_{ij} = 0$ for $i \neq j$, $A_{ii} = A_{11}+A_{22}+A_{33} = 0$, $A_{11} = A_{33}$ if the bedding plane is in ($e_1, e_3$) anisotropic plane, $A_{22} = -2A_{11}$, implying : \[A_{ij}l_il_j = A_{l1}(1-3l_2^2)\] where $A_{11}$ is the component of the microstructure operator $A_{ij}$ in the isotropic (bedding) plane. The late expression for cohesion becomes (Pardoen, 2015)2): \[c= c_0 \left( 1+A_{l1}(1-3l_2^2) + b_1A_{l1}^2(1-3l_2^2)^2 + b_2A_{l1}^3(1-3l_2^2)^3 + … \right)\]
The constants $c_0, A_{11}, b_1, b_2, …$ can be obtained from experimental data and laboratory tests. For uniaxial compression : $l_2 = \cos(\alpha)$ , $\alpha$ being the angle between the compression direction and the normal to the bedding plane (figure 1), $\alpha=0^{\circ}$ if the loading is perpendicular to the bedding plane and $\alpha = 90^{\circ}$ if parallel. The cohesion and the uniaxial compressive straight $f_c$ are linked by $f_c = 2c\cos\phi /(1-\sin\phi)$. The constants can be obtained from results coming from tests performed for differents orientation $\alpha$ (figure 3). For axial compression $\sigma_1$ with confinement $p_0$ (triaxial) : $l_2^2 = \frac{p_0^2\sin^2(\alpha)+\sigma_1^2\cos^2(\alpha)}{2p_0^2+\sigma_1^2}$.
Fig. 3: Schematic view of the cohesion and uniaxial compressive strength evolution as a
function of the angle between the normal to bedding plane and the direction of axial loading
See PLASOL
Remark : For anisotropic Biot’s coeffcient, the deviatoric stress is calculated from the effective stresses (more details about this anisotropy are available in the definition of element CSOL2 and ISOL=9 in Appendix 7).
Plane stress state | NO |
Plane strain state | YES |
Axisymmetric state | YES |
3D state | YES |
Generalized plane state | NO |
Line 1 (2I5, 60A1) | |
---|---|
IL | Law number |
ITYPE | 608 |
COMMENT | Any comment (up to 60 characters) that will be reproduced on the output listing |
Line 1 (14I5) | |
---|---|
NINTV | = number of sub-steps used to integrate numerically the constitutive equation in a time step. |
If NINTV = 0 : number of sub-steps is based on the norm of the deformation increment and on DIV | |
ISOL | = 0 : use of total stresses in the constitutive law |
$\neq$ 0 : use of effective stresses in the constitutive law. See Appendix 7 | |
ICBIF | = 0 : nothing |
1 : Rice bifurcation criterion is computed (only for 2D plane strain analysis) | |
ILODEF | Shape of the yield surface in the deviatoric plane : |
= 1 : circle in the deviatoric plane | |
= 2 : smoothed irregular hexagon in the deviatoric plane | |
ILODEG | Shape of the flow surface in the deviator plane : |
= 1 : circle in the deviatoric plane | |
= 2 : smoothed irregular hexagon in the deviatoric plane | |
IECPS | = 0 : $\psi$ is defined with PSIC and PSIE |
= 1 : $\psi$ is defined with PHMPS | |
KMETH | = 2 : actualised VGRAD integration |
= 3 : Mean VGRAD integration (Default value) | |
IREDUC | = 0 : nothing |
1 : Phi-C reduction method | |
ICOCA | = 0 : nothing |
1 : Capillary cohesion formulation ($c = c_0 + AK1.s + AK2.s^2$) | |
2 : Capillary cohesion formulation ($c = c_0 + AK1.log (s + 0.0001) + AK2$) | |
3 : Young's modulus is a function of suction | |
IBEDDING | Bedding orientation (normal to the bedding plane) |
1 : bedding plane in $e_1e_2$ anisotropic plane (normal $e_3$) | |
2 : bedding plane in $e_1e_3$ anisotropic plane (normal $e_2$) | |
3 : bedding plane in $e_2e_3$ anisotropic plane (normal $e_1$) | |
IANISO | = 0 : anisotropy of cohesion with major principal stress orientation relative to bedding |
= 1 : anisotropy of cohesion by microstructure fabric tensor | |
IVISCO | = 0 : nothing |
1 to 3 : viscoplastic model | |
IDAM | = 0 : nothing |
1 : Damage on elastic properties | |
IHSS | = 0 : nothing |
1 : Coupling stiffness-deformation |
Line 1 (3G10.0) | |
---|---|
ALPHA | Angle of rotation of the anisotropic axis around Z axis (see figure above) |
THETA | Angle of rotation of the anisotropic axis around e1 axis (see figure above) |
PHI | Angle of rotation of the anisotropic axis around e2 axis (see figure above) |
Line 2 (9G10.0) | |
E1 | Elastic Young modulus E($e_1$) |
E2 | Elastic Young modulus E($e_2$) |
E3 | Elastic Young modulus E($e_3$) |
G12 | Elastic shear modulus G($e_1e_2$) |
G13 | Elastic shear modulus G($e_1e_3$) |
G23 | Elastic shear modulus G($e_1e_3$) |
AE1 | E1 = E1 + AE1 * (-SIG(M)), SIG(M) is confinement pressure |
AE2 | E2 = E2 + AE2 * (-SIG(M)) |
AE3 | E3 = E3 + AE3 * (-SIG(M)) |
Line 3 (5G10.0) | |
ANU12 | Poisson ratio NU($e_1e_2$) |
ANU13 | Poisson ratio NU($e_1e_3$) |
ANU23 | Poisson ratio NU($e_2e_3$) |
RHO | Specific mass |
DIV | Size of sub-steps for computation of NINTV (only if NINTV=0; Default value=5.D-3) |
Line 4 (7G10.0) | |
PSIC | Dilatancy angle (in degrees) for compressive paths |
PSIE | Dilatancy angle (in degrees) for extensive paths (ssi ILODEG=2) |
PHMPS | Constant value for definition of |
BIOPT | Bifurcation computation parameter |
AK1 | Capillary cohesion first parameter |
AK2 | Capillary cohesion second parameter |
DECCOH | Cohesion hardening shifting |
Line 5 (7G10.0) | |
PHICF | Final Coulomb's angle (in degrees) for compressive paths |
PHIEF | Final Coulomb’s angle (in degrees) for extensive paths (ssi ILODEF = 2) |
RAYPHIC | Ratio between initial and residual friction angle for compressive paths |
BPHI | Only if there is hardening/softening |
AN | Van Eekelen exponent (default value=-0.229) |
DECPHI | Coulomb’s angle hardening shifting |
RAYPHIE | Ratio between initial and residual friction angle for extensive paths (ssi ILODEF = 2) |
Line 6 (6G10.0) | |
COHF0 | Residual value of cohesion for major principal stress 1 perpendicular to the bedding plane (parallel to the normal to the bedding plane) (if IANISO = 0) |
COHC0 = c0 (if IANISO = 1), see Cohesion anisotropy by microstructure fabric tensor | |
COHFMIN | Minimal residual value of cohesion (if IANISO = 0) |
COHAISO = A11 (if IANISO = 1), see Cohesion anisotropy by microstructure fabric tensor | |
COHF90 | Residual value of cohesion for major principal stress 1 parallel to the bedding plane (perpendicular to the normal to the bedding plane) (if IANISO = 0) |
COHB1 = b1 (if IANISO = 1), see Cohesion anisotropy by microstructure fabric tensor | |
ANGLEMIN | Angle between the normal to the bedding plane and the major principal stress 1 for which the cohesion is minimum (if IANISO = 0) |
COHB2 = b2 (if IANISO = 1), see Cohesion anisotropy by microstructure fabric tensor | |
RAYCOH | Ratio between initial and residual cohesion |
BCOH | Only if there is hardening/softening |
Line 7 (4G10.0) (Only if IHSS = 1) | |
E1F | Final elastic Young modulus E($e_{1f}$) |
E2F | Final elastic Young modulus E($e_{2f}$) |
E3F | Final elastic Young modulus E($e_{3f}$) |
Gamma7 | equivalent strain at which the Young's modulus has reduced to 0.7 times |
Aa | Fitting parameter |
Line 8 (7G10.0) (Only if IECPS = 2 or 3) | |
PSICPEAK | Peak of dilatancy angle for compressive paths (If IECPS=2 then PSICPEAK is the initial value of dilatancy angle |
PSICLIM | Limit value of dilatancy angle for compressive paths |
RATPSI | Ratio between initial and peak of dilatancy angle |
BPSI | Value of EEQU for which PSIC=0.5 (PSICPEAK - PSICLIM) |
PSIEPEAK | Peak of dilatancy angle for extensive paths (If IECPS=2 then PSIEPEAK is the initial value of dilatancy angle) |
PSIELIM | Limit value of dilatancy angle for extensive paths |
DECPSI | Value of EEQU when the dilatancy angle has been half decreased between its initial and final values |
Line 9 (2G10.0) (Only if IDAM = 1) | |
P | Parameter controlling the damage evolution rate |
YD0 | Initial threshold |
= 4 : for 2D analysis
= 6 : for 3D analysis
The stresses are the components of CAUCHY stress tensor in global (X,Y,Z) coordinates.
For 2D analysis :
SIG(1) | $\sigma_{xx}$ |
SIG(2) | $\sigma_{yy}$ |
SIG(3) | $\sigma_{xy}$ |
SIG(4) | $\sigma_{zz}$ |
For the 3-D analysis :
SIG(1) | $\sigma_{xx}$ |
SIG(2) | $\sigma_{yy}$ |
SIG(3) | $\sigma_{zz}$ |
SIG(4) | $\sigma_{xy}$ |
SIG(5) | $\sigma_{xz}$ |
SIG(6) | $\sigma_{yz}$ |
= 48 : for 2D plane strain analysis with bifurcation criterion (ICBIF=1)
= 36 : in all the other cases
Q(1) | = 1 in plane strain state |
= circumferential strain rate ($\dot{\varepsilon_{\theta}}$) in axisymmetrical state | |
Q(2) | = actualised specific mass |
Q(3) | = Reduced deviatoric stress (varies from 0 to 1) |
Q(4) | = 0 if the current state is elastic |
= 1 if the current state is elasto-plastic | |
Q(5) | = equivalent viscoplastic shear strain, i.e. the generalized plastic distorsion, which increment is $\dot{\gamma_{vp}} = \sqrt{\frac{2}{3}\dot{e}_{ij}^{vp}\dot{e}_{ij}^{vp}}$ (see PLASOL) |
Q(6) | = equivalent strain $n^o$1 $\varepsilon_{eq1} = \int \Delta \dot{\varepsilon}_{eq}\Delta t$ |
Q(7) | = equivalent strain indicator $n^o 1$ (Villote $n^o 1$) $\alpha_1 = (\Delta\dot{\varepsilon}_{eq}\Delta t ) / \varepsilon_{eq1}$ |
Q(8) | = $\varepsilon_{xx}$ |
Q(9) | = $\varepsilon_{yy}$ |
Q(10) | = $\varepsilon_{zz}$ |
Q(11) | = $\gamma_{xy} = 2.\varepsilon_{xy}$ |
Q(12) | = equivalent strain $n^o 2$ $\varepsilon_{eq2} = \int \Delta \varepsilon_{eq}$ |
Q(13) | = equivalent strain indicator $n^o 2$ (Villote $n^o 2$) $\alpha_2 = \Delta\varepsilon_{eq} / \varepsilon_{eq2}$ |
Q(14) | = actualised value of equivalent plastic strain $\varepsilon_{ep}^{p}$ |
Q(15) | = actualised value of cohesion for bedding perpendicular to 1st principal stress (COHF0, if IANISO = 0) or COHCO ( if IANISO=1) |
Q(16) | = actualized value of cohesion $c$ |
Q(17) | = actualised value of Coulomb’s friction angle for compr. paths $\phi_C$ |
Q(18) | = actualised value of Coulomb’s friction angle for ext. paths $\phi_E$ |
Q(19) | = 0 : if the stress state is not at the criterion apex |
= 1 : if the stress state is at the criterion apex | |
Q(20) | = number of sub-intervals used for the integration |
Q(21) | = memory of localisation calculated during the re-meshing |
Q(22) | = ? |
Q(23) | = ? |
Q(24) | = ORIENTBED |
Q(25) | = dilatancy angle in compression |
Q(26) | = dilatancy angle in extension |
Q(27) | = damage variable |
Q(28) | = x plastic deformation |
Q(29) | = y plastic deformation |
Q(30) | = z plastic deformation |
Q(31) | = xy plastic deformation |
Q(32) | = ? |
Q(33) | = ? |
Q(34)$\rightarrow$ Q(36) | = reserved for small strain stiffness (E1, E2, E3) |
Q(37)$\rightarrow$ Q(48) | = reserved for bifurcation |