skip to main content

Title: Stabilized formulation for phase‐field fracture in nearly incompressible hyperelasticity

This work presents a stabilized formulation for phase‐field fracture of hyperelastic materials near the limit of incompressibility. At this limit, traditional mixed displacement and pressure formulations must satisfy the inf‐sup condition for solution stability. The mixed formulation coupled with the damage field can lead to an inhibition of crack opening as volumetric changes are severely penalized effectively creating a pressure‐bubble. To overcome this bottleneck, we utilize a mixed formulation with a perturbed Lagrangian formulation which enforces the incompressibility constraint in the undamaged material and reduces the pressure effect in the damaged material. A mesh‐dependent stabilization technique based on the residuals of the Euler–Lagrange equations multiplied with a differential operator acting on the weight space is used, allowing for linear interpolation of all field variables of the elastic subproblem. This formulation was validated with three examples at finite deformations: a plane‐stress pure‐shear test, a two‐dimensional geometry in plane‐stress, and a three‐dimensional notched sample. In the last example, we incorporate a hybrid formulation with an additive strain energy decomposition to account for different behaviors in tension and compression. The results show close agreement with analytical solutions for crack tip opening displacements and performs well at the limit of incompressibility.

 ;  ;  
Award ID(s):
Publication Date:
Journal Name:
International Journal for Numerical Methods in Engineering
Page Range or eLocation-ID:
p. 4655-4673
Wiley Blackwell (John Wiley & Sons)
Sponsoring Org:
National Science Foundation
More Like this
  1. Abstract

    Brittle fracture propagation in rocks is a complex process due to significant grain‐scale heterogeneity and evolving stress states under dynamic loading conditions. In this work, we use digital image correlation and linear elastic fracture mechanics to make instantaneous measurements of the opening (mode I) and in plane shear (mode II) components of the stress intensity field during dynamic mixed mode crack initiation and propagation in crystalline and granular rocks. Both rock types display some similar fracture behaviors as observed in engineered materials, including rate dependent fracture initiation toughness and a direct relationship between propagation toughness and crack velocity; however, measured propagation toughness is higher than quasi‐static values at crack velocities well below the branching velocity in both rocks. Additionally, due to grain scale controls on the fracture process, mixed mode crack propagation is fundamentally different between these two rock types. Mixed mode propagation is energetically more favorable than pure opening mode propagation in sandstone, while the opposite is true in granite. Furthermore, following initiation, propagation in granite occurs so as to minimize the mode II contribution, irrespective of the initiation conditions, while fractures in sandstone maintain a non‐negligible mode II contribution during propagation across the sample.


    We introduce a new finite-element (FE) based computational framework to solve forward and inverse elastic deformation problems for earthquake faulting via the adjoint method. Based on two advanced computational libraries, FEniCS and hIPPYlib for the forward and inverse problems, respectively, this framework is flexible, transparent and easily extensible. We represent a fault discontinuity through a mixed FE elasticity formulation, which approximates the stress with higher order accuracy and exposes the prescribed slip explicitly in the variational form without using conventional split node and decomposition discrete approaches. This also allows the first order optimality condition, that is the vanishing of the gradient, to be expressed in continuous form, which leads to consistent discretizations of all field variables, including the slip. We show comparisons with the standard, pure displacement formulation and a model containing an in-plane mode II crack, whose slip is prescribed via the split node technique. We demonstrate the potential of this new computational framework by performing a linear coseismic slip inversion through adjoint-based optimization methods, without requiring computation of elastic Green’s functions. Specifically, we consider a penalized least squares formulation, which in a Bayesian setting—under the assumption of Gaussian noise and prior—reflects the negative log of the posteriormore »distribution. The comparison of the inversion results with a standard, linear inverse theory approach based on Okada’s solutions shows analogous results. Preliminary uncertainties are estimated via eigenvalue analysis of the Hessian of the penalized least squares objective function. Our implementation is fully open-source and Jupyter notebooks to reproduce our results are provided. The extension to a fully Bayesian framework for detailed uncertainty quantification and non-linear inversions, including for heterogeneous media earthquake problems, will be analysed in a forthcoming paper.

    « less
  3. Abstract In the standard fracture test specimens, the crack-parallel normal stress is negligible. However, its effect can be strong, as revealed by a new type of experiment, briefly named the gap test. It consists of a simple modification of the standard three-point-bend test whose main idea is to use plastic pads with a near-perfect yield plateau to generate a constant crack-parallel compression and install the end supports with a gap that closes only when the pads yield. This way, the test beam transits from one statically determinate loading configuration to another, making evaluation unambiguous. For concrete, the gap test showed that moderate crack-parallel compressive stress can increase up to 1.8 times the Mode I (opening) fracture energy of concrete, and reduce it to almost zero on approach to the compressive stress limit. To model it, the fracture process zone must be characterized tensorially. We use computer simulations with crack-band microplane model, considering both in-plane and out-of-plane crack-parallel stresses for plain and fiber-reinforced concretes, and anisotropic shale. The results have broad implications for all quasibrittle materials, including shale, fiber composites, coarse ceramics, sea ice, foams, and fone. Except for negligible crack-parallel stress, the line crack models are shown to be inapplicable.more »Nevertheless, as an approximation ignoring stress tensor history, the crack-parallel stress effect may be introduced parametrically, by a formula. Finally we show that the standard tensorial strength models such as Drucker–Prager cannot reproduce these effects realistically.« less
  4. Summary

    Methods to compute the stress intensity factors along a three‐dimensional (3D) crack front often display a tenuous rate of convergence under mesh refinement or, worse, do not converge, particularly when applied on unstructured meshes. In this work, we propose an alternative formulation of the interaction integral functional and a method to compute stress intensity factors along the crack front which can be shown to converge. The novelty of our method is the decoupling of the two discretizations: the bulk mesh for the finite element solution and the mesh along the crack front for the numerical stress intensity factors, and hence we term it the multiple mesh interaction integral (MMII) method. Through analysis of the convergence of the functional and method, we find scalings of these two mesh sizes to guarantee convergence of the computed stress intensity factors in a variety of norms, including maximum pointwise error and total variation. We demonstrate the MMII on four examples: a semiinfinite straight crack with the asymptotic displacement fields, the same geometry with a nonuniform stress intensity factor along the crack front, a spherical cap crack in a cylinder under tension, and the elliptical crack under far‐field tension and shear.

  5. The discrete damage model presented in this paper accounts for 42 non-interacting crack microplanes directions. At the scale of the representative volume element, the free enthalpy is the sum of the elastic energy stored in the non-damaged bulk material and in the displacement jumps at crack faces. Closed cracks propagate in the pure mode II, whereas open cracks propagate in the mixed mode (I/II). The elastic domain is at the intersection of the yield surfaces of the activated crack families, and thus describes a non-smooth surface. In order to solve for the 42 crack densities, a Closest Point Projection algorithm is adopted locally. The representative volume element inelastic strain is calculated iteratively using the Newton–Raphson method. The proposed damage model was rigorously calibrated for both compressive and tensile stress paths. Finite element method simulations of triaxial compression tests showed that the transition between brittle and ductile behavior at increasing confining pressure can be captured. The cracks’ density, orientation, and location predicted in the simulations are in agreement with experimental observations made during compression and tension tests, and accurately show the difference between tensile and compressive strength. Plane stress tension tests simulated for a fiber-reinforced brittle material also demonstrated that themore »model can be used to interpret crack patterns, design composite structures and recommend reparation techniques for structural elements subjected to multiple damage mechanisms.« less