Bilinear mixed mode traction separation law

This class implements the bilinear mixed mode traction separation law described in Camanho and Dávila (2002).

Softening onset prediction

The initiation of the softening process is predicted using the quadratic failure criterion given below,

(τ3N)2+(τ2S)2+(τ1T)2=1\left(\frac{\langle\tau_3\rangle}{N}\right)^2 + \left(\frac{\tau_2}{S} \right)^2 + \left(\frac{\tau_1}{T} \right)^2 = 1

The total mixed-mode relative displacement δm\delta_m is defined as δm=δ12+δ22+δ32=δshear2+δ32\delta_m = \sqrt{\delta_1^2+\delta_2^2+\langle\delta_3\rangle^2} = \sqrt{\delta_{shear}^2+\langle\delta_3\rangle^2} where δshear\delta_{shear} represents the norm of the vector defining the tangential relative displacements of the element.

Using the same penalty stiffness in Modes I, II and III, the tractions before softening onset are: τi=Kδi,   i=1,2,3\tau_i = K\delta_i, ~~~\text{i=1,2,3}

Assuming S=TS=T, the single mode relative displacements at softening onset are: δ30=NK\delta_3^0=\frac{N}{K}

δ10=δ20=δshear0=SK\delta_1^0=\delta_2^0=\delta_{shear}^0=\frac{S}{K}

For an opening displacement δ3\delta_3 greater than zero, the mode mixity ratio β\beta is defined as: β=δshearδ3\beta = \frac{\delta_{shear}}{\delta_3}

The mixed-mode relative displacement corresponding to the onset of softening δm0\delta_m^0 is given as δm0={δ30δ101+β2(δ10)2+(βδ30)2,δ3>0δshear0,δ30\delta_m^0= \begin{cases} \delta_3^0\delta_1^0\sqrt{\frac{1+\beta^2}{(\delta_1^0)^2+(\beta\delta_3^0)^2}} , & \delta_3> 0\\ \delta_{shear}^0, & \delta_3\leq 0 \end{cases}

Delamination propagation prediction

Power law criterion

The power law criterion is given as (GIGIC)α+(GIIGIIC)α=1\left(\frac{G_{I}}{G_{IC}}\right)^{\alpha} + \left(\frac{G_{II}}{G_{IIC}}\right)^{\alpha} = 1

The mixed-mode displacements corresponding to total decohesion is given as: δmf={ 2(1+β2)Kδm0[(1GIC)α+(β2GIIC)α]1/α,δ3>0 (δ1f)2+(δ2f)2,δ30\delta_m^f= \begin{cases} \ \frac{2(1+\beta^2)}{K\delta_m^0}\left[\left(\frac{1}{G_{IC}}\right)^{\alpha} + \left(\frac{\beta^2}{G_{IIC}}\right)^{\alpha}\right]^{-1/\alpha} , & \delta_3> 0\\ \ \sqrt{(\delta_1^f)^2+(\delta_2^f)^2}, & \delta_3\leq 0 \end{cases}

B-K criterion

The mixed-mode criterion proposed by Benzeggagh and Kenane is given as (B-K criterion): GIC+(GIICGIC)(GshearGT)η=GC with GT=GI+GshearG_{IC}+(G_{IIC}-G_{IC})\left(\frac{G_{shear}}{G_T}\right)^{\eta} = G_C~\text{with}~G_T=G_I+G_{shear}

The mixed-mode displacements corresponding to total decohesion is given as: δmf={ 2Kδm0[GIC+(GIICGIC(β21+β2)η)],δ3>0 (δ1f)2+(δ2f)2,δ30\delta_m^f= \begin{cases} \ \frac{2}{K\delta_m^0}\left[G_{IC}+(G_{IIC}-G_{IC}\left(\frac{\beta^2}{1+\beta^2}\right)^\eta)\right] , & \delta_3> 0\\ \ \sqrt{(\delta_1^f)^2+(\delta_2^f)^2}, & \delta_3\leq 0 \end{cases}

Constitutive equation for mixed-mode loading

The constitutive equation for mixed-mode loading is given as τs=Dsrδr\tau_s = D_{sr}\delta_r

Dsr={δˉsrK,if δmmaxδm0δˉsr[(1d)K+Kdδˉs3δ3δ3],if δm0<δmmaxδmfδˉs3δˉ3rδ3δ3K,if δmmaxδmfD_{sr} = \begin{cases} \bar{\delta}_{sr}K, \text{if}~\delta_m^{max} \leq \delta_m^0 \\ \bar{\delta}_{sr}\left[(1-d)K+Kd\bar{\delta}_{s3}\frac{\langle-\delta_3\rangle}{-\delta_3}\right], \text{if}~\delta_m^0 < \delta_m^{max} \leq \delta_m^f \\ \bar{\delta}_{s3}\bar{\delta}_{3r}\frac{\langle-\delta_3\rangle}{-\delta_3}K, \text{if}~\delta_m^{max} \geq \delta_m^f \end{cases} (1)

d=δmf(δmmaxδm0)δmmax(δmfδm0),d[0,1]d= \frac{\delta_m^f(\delta_m^{max}-\delta_m^0)}{\delta_m^{max}(\delta_m^f-\delta_m^0)}, d\in[0,1] (2)

Solver options

Viscous regularization

Cohesive zone models exhibiting softening behavior and stiffness degradation often lead to convergence difficulties in an implicit solver. The traction-separation laws can be regularized using viscosity. The viscous damage variable dvd_v is defined by dv˙=1μ(ddv)\dot{d_v}=\frac{1}{\mu}(d-d_v) where μ\mu is the viscosity parameter representing the relaxation time of the viscous system. An analytical expression of dvd_v can be obtained by using the backward Euler method. With viscous regularization, the dd will be replaced by dvd_v in Eq. (1) to compute traction.

Lag separation state

It is typically useful to improve convergence by lagging the separation state. When lag_separation_state = true, the δ3\delta_3, δmmax\delta_m^{max}, δm0\delta_m^0 and δmf\delta_m^f will be replaced by their old values from previous time step.

Use Regularized Heaviside Function

The step (Heaviside) function δ3δ3\frac{\langle-\delta_3\rangle}{-\delta_3} in Eq. (1) usually makes convergence bad. In the code, we replaced it with the regularized Heaviside function which provides a C0 continuity. The regularization parameter can be set by alpha parameter.

Examples

[Materials]
  [czm]
    type = BiLinearMixedModeTraction
    boundary = 'interface'
    penalty_stiffness = 1e6
    GI_c = 1e3
    GII_c = 1e2
    normal_strength = 1e4
    shear_strength = 1e3
    displacements = 'disp_x disp_y'
    eta = 2.2
    viscosity = 1e-3
  []
[]
(contrib/moose/modules/solid_mechanics/test/tests/cohesive_zone_model/bilinear_mixed.i)

References

  1. Pedro P. Camanho and Carlos G. Dávila. Mixed-mode decohesion finite elements for the simulation of delamination in composite materials. Technical Report NASA/TM-2002-211737, National Aeronautics and Space Administration, 2002.[BibTeX]