User:Cbonneville/sandbox

Computational Methods In Fracture Mechanics
Many methods exist for calculating G with finite elements. Although a direct calculation of the J-integral is possible (using the strains and stresses outputted by finite element analysis), approximate approaches for some types of crack growth exist and provide reasonable accuracy with straightforward calculations.

Nodal Release Method
If the crack is growing straight, the energy release rate can be decomposed as a sum of 3 terms $$G_i$$associated with the energy in each 3 modes. As a result, the Nodal Release method (NR) can be used to determine $$G_i$$from FEA results. The energy release rate is calculated at the nodes of the finite element mesh for the crack at an initial length and extended by a small distance $$\Delta a$$. First, we calculate the displacement variation at the node of interest $$\Delta\vec{u}=\vec{u}^{(t+1)}-\vec{u}^{(t)}$$(before and after the crack tip node is released). Secondly, we keep track of the nodal force $$\vec{F}$$ outputted by FEA. Finally, we can find each components of $$G$$ using the following formulas:

$$G_1^{\text{NR}}=\frac{1}{\Delta a}F_2\frac{\Delta u_2}{2} $$$$G_2^{\text{NR}}=\frac{1}{\Delta a}F_1\frac{\Delta u_1}{2} $$$$G_3^{\text{NR}}=\frac{1}{\Delta a}F_3\frac{\Delta u_3}{2} $$Where $$\Delta a$$ is the width of the element bounding the crack tip. The accuracy of the method highly depends on the mesh refinement, both because the displacement and forces depend on it, and because $$G=\lim_{\Delta a \to 0}G^{\text{NR}}$$. Note that the equations above are derived using the crack closure integral.

If the energy release rate exceed a critical value, the crack will grow. In this case, a new FEA simulation is perfomed (for the next time step) where the node at the crack tip is released. For a bounded substrate, we may simply stop enforcing fixed Dirichlet boundary conditions at the crack tip node of the previous time step (i.e. displacements are no longer restrained). For a symmetric crack, we would need to update the geometry of the domain with a longer crack opening (and therefore generate a new mesh ).

MCCI Method
Terrence is doing this part

Domain Integral Approach for J
The J-integral may be calculated directly using the finite element mesh and shape functions. We consider a domain countour as shown in figure 2 and choose an arbitrary smooth function $$\tilde{q}(x_1,x_2)=\sum_iN_i(x_1,x_2)\tilde{q}_i$$ such that $$\tilde{q}=1$$ on $$\Gamma$$ and  $$\tilde{q}=0$$ on  $$\mathcal{C}_1$$.

For linear elastic cracks growing straight ahead, $$G=J$$. The energy release rate can then be calculated over the area bounded by the contour using an updated formulation:

$$J=\int_\mathcal{A}(\sigma_{ij}u_{i,1}\tilde{q}_{,j}-W\tilde{q}_{,1})d\mathcal{A}$$

The formula above may be applied to any annular area surrounding the crack tip (in particular, a set of neighboring elements can be used). This method is very accurate, even with a coarse mesh arround the crack tip (one may choose an integration domain located far away, with stresses and displacement less sensitive to mesh refinment)


 * {| class="toccolours collapsible collapsed" width="80%" style="text-align:left"

!Derivation of the J-integral for domain integral approach The J-intregral may be expressed over the full contour as follow:

$$G=J=\int_{\Gamma}(Wn_1-t_i\frac{\partial u_i}{\partial x_1})d\Gamma=\int_{\mathcal{C}_1+\mathcal{C}_++\mathcal{C}_--\mathcal{C}}(Wn_1-\sigma_{ij}n_j\frac{\partial u_i}{\partial x_1})\tilde{q}d\Gamma$$

With $$\mathcal{C}=\mathcal{C}_1+\mathcal{C}_++\mathcal{C}_--\Gamma$$. $$\vec{n}=-\vec{m}$$, $$\tilde{q}=0$$ on $$\mathcal{C}_1$$and the work and stresses cancel out on $$\mathcal{C}_+$$and $$\mathcal{C}_-$$, hence by application of the divergence theorem this leads to:

$$J=\int_\mathcal{C}(-Wm_1+\sigma_{ij}m_i\frac{\partial u_i}{\partial x_1})d\mathcal{C}=\int_\mathcal{A}(\frac{\partial\sigma_{ij}u_{i,1}\tilde{q}}{\partial x_j}-\frac{\partial W\tilde{q}}{\partial x_1})d\mathcal{A}$$

Finally, by noting that $$\frac{\partial W}{\partial x_1}=\sigma_{ij}\frac{\partial^2 u_i}{\partial x_j \partial x_1}$$and using the equilibrium equation:

$$J=\int_\mathcal{A}(\sigma_{ij}u_{i,1}\tilde{q}_{,j}-W\tilde{q}_{,1})d\mathcal{A}$$


 * }

Crack Tip Singular Elements
Terrence is doing this part