FiniteStrainCrystalPlasticity is deprecated. Use the CrystalPlasticityUOBase based crystal plasticity system.

The crystal plasticity material class in TensorMechanics is based on plastic flow on individual slip systems to obtain the inelastic deformation in materials. The effect of grain-scale anisotropy on the overall plastic behavior can be incorporated through this model. The present formulation considers large deformation and is based on a stress update algorithm. Some of the important functions associated with this class are outlined below.

Function - computeQpStress(): User can override this function.

Input:

  • 2$$$^{nd}$$$ Piola-Kirchhoff (PK2) stress in the intermediate configuration at previous increment (variable _pk2_old, symbol $$$T_n$$$);

  • Plastic deformation gradient at previous increment (variable _fp_old, symbol $$$F^p_n$$$);

  • Slip system resistances at previous increment (variable _gss_old, symbol $$$g^\alpha_n$$$);

  • Current deformation gradient (variable _dfgrd, symbol $$$F$$$).

Output:

  • Current PK2 stress (variable _pk2, symbol $$$T$$$);

  • Current plastic deformation gradient (variable _fp, symbol $$$F^p$$$);

  • Current slip system resistances (variable _gss, symbol $$$g^\alpha$$$).

  • Slip system resistances $$$g^\alpha$$$ are solved iteratively using a predictor corrector algorithm. $$$T$$$ is solved using Newton Raphson.

Flowchart:

i: At iteration i=0, $$$g^\alpha_i=g^\alpha_n$$$

ii: Calculate the residual r and Jacobian J=$$$\partial r/\partial T$$$ from the function calc_resid_jacob

iii: Update T as: $$$T_{i+1}=T_{i}-J^{-1}r$$$

iv: Check $$$||r||_2 <$$$ _rtol (Optional user defined parameter)

If False then go to Step ii

Else

v: Calculate $$$g^\alpha_{i+1}$$$ from function update_gss

vi: Check $$$max|g^\alpha_{i+1}-g^\alpha_{i}| <$$$ _gtol (Optional user defined parameter)

If False then go to Step ii

Else

vii: Obtain rotation tensor R from function getmatrot

viii: Update rotation as $$$\tilde{R}=RQ$$$ where Q is the rotation tensor from Euler angles

Exit Function

Function - calc_resid_jacob: User can override this function.

Input:

  • PK2 stress at previous iteration (variable _pk2, symbol $$$T_i$$$, i - iteration);

  • Inverse of plastic deformation gradient at previous increment (variable _fp_old_inv, symbol $$$F^{p-1}_n$$$);

  • Current deformation gradient (variable _dfgrd, symbol $$$F$$$).

Output:

  • Residual (variable resid, symbol r);

  • Jacobian (variable jac, symbol J);

  • Updated inverse of plastic deformation gradient (variable _fp_inv, symbol $$$F^{p-1}$$$)

Residual r as implemented:

i. Elastic component of deformation gradient (variable fe) $$$F^e=FF^{p-1}_n$$$.

ii. Elastic right Cauchy-Green Tensor (variable ce) $$$C^e=F^{eT}F^e$$$.

iii. Elastic Green-Lagrangian strain tensor (variable ee) $$$E^e=\frac{1}{2}\left(C^e -I\right)$$$.

iv. Shear stresses resolved on slip systems (variable tau) $$$\tau^\alpha=C^eT_i:\left( m_0^\alpha \otimes n_0^\alpha \right)/J^e$$$, where $$$m_0^\alpha$$$ and $$$n_0^\alpha$$$ are the slip direction and slip system normal in the intermediate configuration respectively. $$$J^e=det(F^e)$$$.

v. Slip increment on slip systems $$$\Delta \gamma^\alpha$$$ obtained from function get_slip_incr.

vi. Resultant slip increment (variable eqv_slip_incr) $$$\Delta \gamma_{res}=I-\sum \limits_\alpha \left( \Delta \gamma^\alpha S_0^\alpha\right)$$$ where $$$S_0^\alpha=m_0^\alpha \otimes n_0^\alpha$$$ is the Schmid tensor.

vii. Current plastic component of deformation gradient $$$F^{p-1}=F^{p-1}_n\Delta \gamma_{res}$$$.

viii. Elastic component of deformation gradient in iteration i+1, $$$F^e=FF^{p-1}$$$.

ix: PK2 stress due to $$$T_i$$$ and associated slip increment in current iteration i+1 (variable pk2_new), $$$T^\prime=C:\frac{1}{2}\left(F^{eT}F^e-I\right)$$$.

x. Residual $$$r=T_i-T^\prime$$$.

Jacobian J formulation:

i. $$$r=T_i-T^\prime$$$ and $$$J=\frac{\partial r}{\partial T_i}=I-\frac{\partial T^\prime}{\partial T_i}$$$.

ii. $$$T^\prime=C:E^e$$$ which gives $$$\frac{\partial T^\prime}{\partial T_i}=C: \frac{\partial E^e}{\partial T_i}$$$.

iii. $$$E^e=\frac{1}{2}\left(F^{eT}F^e -I\right)$$$ which gives in indicial notation $$$\frac{\partial E^e_{ij}}{\partial F^e_{kl}}=\frac{1}{2}\left( \delta_{il}F^e_{kj}+\delta_{jl}F^e_{ki}\right)$$$ where $$$\delta_{ij}$$$ is the Kroneckar delta.

iv. $$$F^e=FF^{p-1}$$$ which gives in indicial notation $$$\frac{\partial F^e_{ij}}{\partial F^{p-1}_{kl}}=F_{ik}\delta_{lj}$$$.

v. $$$F^{p-1}=F^{p-1}_n\Delta \gamma_{res}$$$ which gives $$$\frac{\partial F^{p-1}}{\partial \Delta \gamma_{res}}=F^{p-1}_n$$$.

vi. $$$\Delta \gamma_{res}=I-\sum \limits_\alpha \left( \Delta \gamma^\alpha S_0^\alpha\right)$ which gives $\frac{\partial \Delta \gamma_{res}}{\partial \Delta \gamma^\alpha}=-\sum\limits_\alpha S_0^\alpha$$$.

vii. $$$\Delta \gamma^\alpha=f(\tau^\alpha, g^\alpha)$$$ which gives $$$\frac{\partial \Delta \gamma^\alpha}{\partial \tau^\alpha}$$$.

viii. $$$\tau^\alpha=C^eT_i:S_0^\alpha/J^e$$$ which assuming $$$C^e \approx 1$$$ and $$$J^e \approx 1$$$ gives $$$\frac{\partial \tau^\alpha}{\partial T_i}=S_0^\alpha$$$.

ix. Hence $$$\frac{\partial T^\prime}{\partial T_i}=C: \frac{\partial E^e}{\partial F^e} \frac{\partial F^e}{\partial F^{p-1}} \frac{\partial F^{p-1}}{\partial \Delta \gamma_{res}}\frac{\partial \Delta \gamma_{res}}{\partial \Delta \gamma^\alpha}\frac{\partial \Delta \gamma^\alpha}{\partial \tau^\alpha}\frac{\partial \tau^\alpha}{\partial T_i}$$$.

Jacobian J as implemented:

i. Variable dtaudpk2, $$$\frac{\partial \tau^\alpha}{\partial T_i}$$$.

ii. Variable dfpinvdslip, $$$\frac{\partial F^{p-1}}{\partial \Delta \gamma^\alpha}$$$.

iii. Variable dfedfpinv, $$$\frac{\partial F^e}{\partial F^{p-1}}$$$.

iv: Variable deedfe, $$$\frac{\partial E^e}{\partial F^e}$$$.

v. Variable dfpinvdpk2, $$$\frac{\partial F^{p-1}}{\partial T_i}=\frac{\partial F^{p-1}}{\partial \Delta \gamma^\alpha} \frac{\partial \Delta \gamma^\alpha}{\partial \tau^\alpha} \frac{\partial \tau^\alpha}{\partial T_i}$$$.

vi. Variable dslipdtau, $$$\frac{\partial \Delta \gamma^\alpha}{\partial \tau^\alpha}$$$ obtained from function get_slip_incr.

vii. Variable jac, $$$J=C:\frac{\partial E^e}{\partial F^e} \frac{\partial F^e}{\partial F^{p-1}} \frac{\partial F^{p-1}}{\partial T_i}$$$, followed by $$$J=\mathfrak{I}-J$$$ where $$$\mathfrak{I}$$$ is the fourth order identity tensor.

Update Cauchy stress (variable sig) $$$\sigma=F^e T_i F^{eT}/J^e$$$.

Function - update_gss: User can override this function.

Input:

  • Slip increment (variable slip_incr, symbol $$$\Delta \gamma^\alpha$$$);

  • Slip system resistance from previous increment (variable _gss_old, symbol $$$g^\alpha_n$$$);

  • Current slip system resistance (variable _gss, symbol $$$g^\alpha$$$);

  • Hardness properties read from .i file (variables _hprops, _h0; symbols h,$$$h_0$$$);

  • Saturation and initial slip system resistances read from .i file (variables _tau_sat, _tau_init, symbols $$$\tau_s$$$, $$$\tau_0$$$);

  • Latent hardening variable read from .i file (variable _r, symbol q).

Output:

  • Current slip system resistance (variable _gss, symbol $$$g^\alpha$$$).

Implementation:

i. Variable hb $$$h_\beta=h_0 \left( 1-\frac{g^\beta}{\tau_s} \right)^a$$$

ii. Update _gss as $$$g^\alpha=g^\alpha_n+\sum\limits_\beta q_{\alpha \beta}h_\beta$$$ where $$$q_{\alpha \beta}=q$$$ for out of plane slip systems and $$$q_{\alpha \beta}=1$$$ for coplanar systems.

The hardening model provided in [Kalidindi et. al., J. Mech. Phys. Solids, Vol. 40, pp. 537-569, 1992] has been implemented.

Function - get_slip_incr: User can override this function.

Input:

  • Resolved shear stresses (variable tau, symbol $$$\tau^\alpha$$$);

  • Current slip system resistance (variable _gss, symbol $$$g^\alpha$$$);

  • Time increment (variable _dt, symbol $$$\Delta t$$$);

  • Flow properties read from .i file (variable _a0, _xm; symbols $$$a_0, m$$$).

Output:

  • Slip increment (variable slip_incr, symbol $$$\Delta \gamma^\alpha$$$);

  • Derivative of slip increment w.r.t. resolved shear stress (variable dslipdtau, symbol $$$\frac{\partial \Delta \gamma^\alpha}{\partial \tau^\alpha}$$$).

Implementation:

i. Slip increment $$$\Delta \gamma^\alpha=a_0\left|\frac{\tau^\alpha}{g^\alpha}\right|^\frac{1}{m}sign(\tau^\alpha)\Delta t$$$.

ii. Derivative $$$\frac{\partial \Delta \gamma^\alpha}{\partial \tau^\alpha} =\frac{a_0}{m g^\alpha}\left|\frac{\tau^\alpha}{g^\alpha}\right|^{\frac{1}{m}-1} \Delta t$$$.

Function - get_mat_rot: Performs RU decomposition of F where R is orthogonal and U is symmetric postive definite.

Input:

  • Second order positive definite tensor (variable a, symbol F).

Output:

  • Orthogonal tensor (variable rot, symbol R).

Implementation:

i. Variable C, $$$C=F^TF=U^TU=U^2$$$. C has the same eigen vectors as U but the eigen values are $$$\lambda_C=\lambda_U^2$$$.

ii. Perform eigen value analysis of C using function dsyev_ and obtain diagonal matrix $$$D_C$$$ of eigen values and eigen vector P.

iii. Construct diagonal matrix $$$D_U=D_C^{1/2}$$$.

iv. Obtain $$$U=PD_UP^T$$$.

v. Obtain $$$R=FU^{-1}$$$.

Function - get_euler_ang:

Input:

  • Euler angle input file name defined in .i file. Euler angle for every element must be prescribed as element no., euler 1, euler 2, euler 3.

  • Optional user input: Euler angles can also be read from _euler_angle_1, _euler_angle_2 and _euler_angle_3 variables in FiniteStrainCrystalPlasticity material for every block. Default is 0.0, 0.0, 0.0.

Output:

  • _euler_angle_1, _euler_angle_2 and _euler_angle_3.

Function - computeQpElasticityTensor: Defines tangent moduli K and can be used as preconditioner for JFNK. User can override this function.

Input:

  • Current plastic deformation gradient (variable _fp, symbol $$$F^p$$$);

  • Current deformation gradient (variable _dfgrd, symbol $$$F$$$)

Output:

  • Tangent moduli (variable _Jacobian_mult,symbol K).

Formulation:

i. $$$K=\frac{\partial \sigma}{\partial F}$$$.

ii. $$$\sigma=F^eTF^{eT}/J^e$$$ which gives $$$\frac{\partial \sigma}{\partial F}=\frac{1}{J^e}\left(\frac{\partial F^e}{\partial F}TF^{eT}+F^e\frac{\partial T}{\partial F}F^{eT}+F^eT\frac{\partial F^{eT}}{\partial F}\right)$$$.

iii. In indicial notation $$$\frac{\partial F^e_{ij}}{\partial F_{kl}}=\delta_{ki}F^{p-1}_{lj}$$$.

iv. $$$T=C:E^e$$$ which gives $$$\frac{\partial T}{\partial F}=\frac{\partial T}{\partial E^e}\frac{\partial E^e}{\partial F^e}\frac{\partial F^e}{\partial F}=C:\frac{\partial E^e}{\partial F^e}\frac{\partial F^e}{\partial F}$$$.

v. In indicial notation $$$\frac{\partial E^e_{ij}}{\partial F^e_{kl}}=\frac{1}{2}\left(\delta_{il}F^e_{kj}+\delta_{jl}F^e_{ki}\right)$$$.

Implementation:

i. Variable dfedf calculates $$$\frac{\partial F^e}{\partial F}$$$.

ii. Variable deedfe calculates $$$\frac{\partial E^e}{\partial F^e}$$$.

iii. Variable dsigspk2dfe calculates $$$F^{eT} C:\frac{\partial E^e}{\partial F^e} F^e$$$.