diff --git a/doc/content/bib/raccoon.bib b/doc/content/bib/raccoon.bib index b7acd1f4c8a..b7769834f0d 100644 --- a/doc/content/bib/raccoon.bib +++ b/doc/content/bib/raccoon.bib @@ -8,7 +8,7 @@ @article{bourdin2007numerical year={2007} } -@Article{Kumar2022, +@Article{Kumar2022, author={Kumar, A. and Ravi-Chandar, K. and Lopez-Pamies, O.}, @@ -53,3 +53,11 @@ @article{Drucker-Prager volume = {10}, year = {1952} } + +@article{larsen2024, + title={A variational formulation of Griffith phase-field fracture with material strength}, + author={C. J. Larsen and J. E. Dolbow and O. Lopez-Pamies}, + journal = {Arxiv}, + year={2024}, + url={https://arxiv.org/abs/2401.13938}, +} diff --git a/doc/content/source/materials/KLBFNucleationMicroForce.md b/doc/content/source/materials/nucleation_models/KLBFNucleationMicroForce.md similarity index 95% rename from doc/content/source/materials/KLBFNucleationMicroForce.md rename to doc/content/source/materials/nucleation_models/KLBFNucleationMicroForce.md index f1092fc3e2c..80a1c516cf2 100644 --- a/doc/content/source/materials/KLBFNucleationMicroForce.md +++ b/doc/content/source/materials/nucleation_models/KLBFNucleationMicroForce.md @@ -8,7 +8,7 @@ In general, model KLR published in [!cite](Kumar2022) is recommended over model ### KLR (Kumar, Lopez, Ravi-Chandar) Model 2022 -Please see [KLRNucleationMicroForce](source/materials/KLRNucleationMicroForce.md) +Please see [KLRNucleationMicroForce](nucleation_models/KLRNucleationMicroForce.md) ### KLBF (Kumar, Lopez, Bourdin, Francfort) Model 2020 diff --git a/doc/content/source/materials/KLRNucleationMicroForce.md b/doc/content/source/materials/nucleation_models/KLRNucleationMicroForce.md similarity index 94% rename from doc/content/source/materials/KLRNucleationMicroForce.md rename to doc/content/source/materials/nucleation_models/KLRNucleationMicroForce.md index 51c35e2e572..064df7621fe 100644 --- a/doc/content/source/materials/KLRNucleationMicroForce.md +++ b/doc/content/source/materials/nucleation_models/KLRNucleationMicroForce.md @@ -74,13 +74,9 @@ l_e = \frac{l_0}{\sqrt{1+\delta}},\quad (\delta > -1) could deviate from $\ell$. The mesh should be able to resolve $l_e$ as well. -Here is an example input deck - -!listing /tutorials/surfing_boundary_problem/fracture.i block=Materials/kumar_material - ### KLBF (Kumar, Lopez, Bourdin, Francfort) Model 2020 -See [KLBFNucleationMicroForce](source/materials/KLBFNucleationMicroForce.md) +See [KLBFNucleationMicroForce](nucleation_models/KLBFNucleationMicroForce.md) ## Example Input File Syntax diff --git a/doc/content/source/materials/nucleation_models/LDLNucleationMicroForce.md b/doc/content/source/materials/nucleation_models/LDLNucleationMicroForce.md new file mode 100644 index 00000000000..dbcd098c5ad --- /dev/null +++ b/doc/content/source/materials/nucleation_models/LDLNucleationMicroForce.md @@ -0,0 +1,82 @@ +# LDLNucleationMicroForce + +!syntax description /Materials/LDLNucleationMicroForce + +See *A variational formulation of Griffith phase-field fracture with material strength* by [!cite](larsen2024) + +This is a variational recast of the phase-field fracture formulation put forth by Kumar, Francfort and Lopez-Pamies (2018), where $\delta$ is explicitly evaluated. + +## Overview + +### LDL (Larsen, Dolbow, Lopez) Model 2024 + +The fracture functional reads + +\begin{equation} +\mathcal{E}_f^l(\boldsymbol{u}, d):= \int_\Omega g(d) \psi_e(\boldsymbol{\varepsilon}(\boldsymbol{u})) \;\mathrm{dV} + \frac{\delta^l G_c}{c_0}\int_\Omega\left(\frac{\alpha(d)}{l} + l\nabla d\cdot \nabla d\right) \;\mathrm{dV} + \int_\Omega\left(\int_0^1 g(d) \;\mathrm{dd}\right)\widehat{c_e}(\boldsymbol{\varepsilon}(\boldsymbol{u})) \;\mathrm{dV}. +\end{equation} + +The governing equation for fracture is +\begin{equation} +-\nabla\cdot \dfrac{2 G_c l \delta^l}{c_0}\nabla d + \dfrac{\partial}{\partial d}\left( \alpha \dfrac{\delta^l G_c}{c_0 l} - g(d)\psi_e\right) - g(d)\widehat{c_e} \ge 0. +\end{equation} + +Note that the following derivation follows settings of `AT-1`, i.e., $\alpha(d)=d$ and $c_0=\dfrac{8}{3}$. In this model, the strict irreversibility condition for $d$ is adopted. + +The "undamaged" nucleation driving force is defined as +\begin{equation} +\widehat{c_e}(\boldsymbol{\varepsilon}(\boldsymbol{u}))= \alpha_2 \sqrt{J_2} + \alpha_1 I_1 - (1-d)\left(1 - \dfrac{\sqrt{I_1^2}}{I_1}\right)\psi_e. +\end{equation} + +In this recast, $\delta^l$ is evaluated explicitly as +\begin{equation} +\text{without h correction:}\quad\delta^l=\left(\frac{\sigma_{ts}+ (1+2\sqrt{3})\sigma_{hs}}{(8+3\sqrt{3})\sigma_{hs}}\right)\frac{3G_c}{16 \psi_{ts} l} + \frac{3}{8}, +\end{equation} + +\begin{equation} +\text{with h correction:}\quad\delta^l=\left(1+\frac{3 h}{8 l}\right)^{-2}\left(\frac{\sigma_{ts}+(1+2\sqrt{3})\sigma_{hs}}{(8+3\sqrt{3})\sigma_{hs}}\right)\frac{3G_c}{16 \psi_{ts} l} + \left(1 + \frac{3h}{8l}\right)^{-1}\frac{2}{5}, +\end{equation} + +where +\begin{equation} +\alpha_1 = - \frac{1}{\sigma_{hs}}\delta^l \frac{G_c}{8l} + \frac{2\psi_{hs}}{3\sigma_{hs}}, \qquad \alpha_2= - \frac{\sqrt{3}(3\sigma_{hs} - \sigma_{ts})}{\sigma_{hs}\sigma_{ts}}\delta^l \frac{G_c}{8l} - \frac{2\psi_{hs}}{\sqrt{3}\sigma_{hs}} + \frac{2\sqrt{3}\psi_{ts}}{\sigma_{ts}}, +\end{equation} + +with +\begin{equation} +\psi_{ts}=\frac{\sigma^2_{ts}}{2E}, \quad \psi_{hs}=\frac{\sigma^2_{hs}}{2\kappa}. +\end{equation} + +The material properties used are the bulk modulus $\kappa$ and Young's modulus $E$ and the shear modulus $\mu$. + +The strength surface is presented with uniaxial tensile $\sigma_{ts}$ and hydrostatic strength $\sigma_{hs}$, where $\sigma_{hs}=\dfrac{2}{3} \dfrac{\sigma_{ts}\sigma_{cs}}{\sigma_{cs} - \sigma_{ts}}$. + +### KLR (Kumar, Lopez, Ravi-Chandar) Model 2022 + +Please see [KLRNucleationMicroForce](nucleation_models/KLRNucleationMicroForce.md) + +### KLBF (Kumar, Lopez, Bourdin, Francfort) Model 2020 + +Please see [KLBFNucleationMicroForce](nucleation_models/KLBFNucleationMicroForce.md) + +## Example Input File Syntax + +``` +[Materials] + [micro_force] + type = LDLNucleationMicroForce + regularization_length = l + normalization_constant = c0 + tensile_strength = sigma_ts + hydrostatic_strength = sigma_hs + external_driving_force_name = ce + h_correction = true + [] +[] +``` + +!syntax parameters /Materials/LDLNucleationMicroForce + +!syntax inputs /Materials/LDLNucleationMicroForce + +!syntax children /Materials/LDLNucleationMicroForce \ No newline at end of file diff --git a/include/materials/KLRNucleationMicroForce.h b/include/materials/KLRNucleationMicroForce.h deleted file mode 100644 index 93fce169dc3..00000000000 --- a/include/materials/KLRNucleationMicroForce.h +++ /dev/null @@ -1,77 +0,0 @@ -//* This file is part of the RACCOON application -//* being developed at Dolbow lab at Duke University -//* http://dolbow.pratt.duke.edu - -#pragma once - -#include "Material.h" -#include "BaseNameInterface.h" -#include "DerivativeMaterialPropertyNameInterface.h" - -/** - * The class implements the external driving force to recover a Drucker-Prager - * strength envelope developed by Kumar et al. in 2022. See Kumar, A., Ravi-Chandar, K. & - * Lopez-Pamies, O. The revisited phase-field approach to brittle fracture: application to - * indentation and notch problems. Int J Fract 237, 83–100 (2022). - * https://doi.org/10.1007/s10704-022-00653-z. all parameters are required to be Material type, not - * double type - */ -class KLRNucleationMicroForce : public Material, - public BaseNameInterface, - public DerivativeMaterialPropertyNameInterface -{ -public: - static InputParameters validParams(); - - KLRNucleationMicroForce(const InputParameters & parameters); - -protected: - virtual void computeQpProperties() override; - - /// Name of the external driving force - const MaterialPropertyName _ex_driving_name; - - /// The external driving force - ADMaterialProperty & _ex_driving; - - ///@{ Phase field properties - /// The fracture toughness - const ADMaterialProperty & _Gc; - /// The normalization constant - const ADMaterialProperty & _c0; - /// phase field regularization length - const ADMaterialProperty & _L; - ///@} - - /// Lame's first parameter - const ADMaterialProperty & _lambda; - /// The shear modulus - const ADMaterialProperty & _mu; - - /// The critical tensile strength - const ADMaterialProperty & _sigma_ts; - - /// The critical compressive strength - const ADMaterialProperty & _sigma_cs; - - /// The materiel and model dependent parameter - const ADMaterialProperty & _delta; - - /// The stress tensor - const ADMaterialProperty & _stress; - - /// Name of the stress space balance - const MaterialPropertyName _stress_balance_name; - /// Quantifying how far is the stress state from stress surface - ADMaterialProperty & _stress_balance; - - ADMaterialProperty & _druck_prager_balance; - - /// Name of the phase-field variable - const VariableName _d_name; - // @{ The degradation function and its derivative w/r/t damage - const MaterialPropertyName _g_name; - const ADMaterialProperty & _g; - const ADMaterialProperty & _dg_dd; - // @} -}; diff --git a/include/materials/nucleation_models/KLBFNucleationMicroForce.h b/include/materials/nucleation_models/KLBFNucleationMicroForce.h new file mode 100644 index 00000000000..fd547326d3e --- /dev/null +++ b/include/materials/nucleation_models/KLBFNucleationMicroForce.h @@ -0,0 +1,31 @@ +//* This file is part of the RACCOON application +//* being developed at Dolbow lab at Duke University +//* http://dolbow.pratt.duke.edu + +#pragma once + +#include "NucleationMicroForceBase.h" + +/** + * The class implements the external driving force to recover a Drucker-Prager + * strength envelope. See Kumar et. al. https://doi.org/10.1016/j.jmps.2020.104027 for model 2020. + */ +class KLBFNucleationMicroForce : public NucleationMicroForceBase +{ +public: + static InputParameters validParams(); + + KLBFNucleationMicroForce(const InputParameters & parameters); + +protected: + virtual void computeQpProperties() override; + + /// The critical tensile strength + const ADMaterialProperty & _sigma_ts; + + /// The critical compressive strength + const ADMaterialProperty & _sigma_cs; + + /// The materiel and model dependent parameter + const ADMaterialProperty & _delta; +}; diff --git a/include/materials/nucleation_models/KLRNucleationMicroForce.h b/include/materials/nucleation_models/KLRNucleationMicroForce.h new file mode 100644 index 00000000000..88cdd60be0a --- /dev/null +++ b/include/materials/nucleation_models/KLRNucleationMicroForce.h @@ -0,0 +1,37 @@ +//* This file is part of the RACCOON application +//* being developed at Dolbow lab at Duke University +//* http://dolbow.pratt.duke.edu + +#pragma once + +#include "NucleationMicroForceBase.h" + +/** + * The class implements the external driving force to recover a Drucker-Prager + * strength envelope developed by Kumar et al. in 2022. See Kumar, A., Ravi-Chandar, K. & + * Lopez-Pamies, O. The revisited phase-field approach to brittle fracture: application to + * indentation and notch problems. Int J Fract 237, 83–100 (2022). + * https://doi.org/10.1007/s10704-022-00653-z. + */ +class KLRNucleationMicroForce : public NucleationMicroForceBase +{ +public: + static InputParameters validParams(); + + KLRNucleationMicroForce(const InputParameters & parameters); + +protected: + virtual void computeQpProperties() override; + + /// The critical tensile strength + const ADMaterialProperty & _sigma_ts; + + /// The critical compressive strength + const ADMaterialProperty & _sigma_cs; + + /// The materiel and model dependent parameter + const ADMaterialProperty & _delta; + + /// Quantifying how far is the stress state from Drucker Prager + ADMaterialProperty & _druck_prager_balance; +}; diff --git a/include/materials/nucleation_models/LDLNucleationMicroForce.h b/include/materials/nucleation_models/LDLNucleationMicroForce.h new file mode 100644 index 00000000000..e8fb4688c5c --- /dev/null +++ b/include/materials/nucleation_models/LDLNucleationMicroForce.h @@ -0,0 +1,35 @@ +//* This file is part of the RACCOON application +//* being developed at Dolbow lab at Duke University +//* http://dolbow.pratt.duke.edu + +#pragma once + +#include "NucleationMicroForceBase.h" + +/** + * The class implements the external driving force to recover a Drucker-Prager + * strength envelope. See Larsen et. al. https://doi.org/10.48550/arXiv.2401.13938 for model 2024. + */ + +class LDLNucleationMicroForce : public NucleationMicroForceBase +{ +public: + static InputParameters validParams(); + + LDLNucleationMicroForce(const InputParameters & parameters); + +protected: + virtual void computeQpProperties() override; + + /// The critical tensile strength + const ADMaterialProperty & _sigma_ts; + + /// The critical hydrostatic strength + const ADMaterialProperty & _sigma_hs; + + /// The materiel and model dependent parameter + ADMaterialProperty & _delta; + + /// Whether to use h correction formula for delta + bool _h_correction; +}; diff --git a/include/materials/KLBFNucleationMicroForce.h b/include/materials/nucleation_models/NucleationMicroForceBase.h similarity index 62% rename from include/materials/KLBFNucleationMicroForce.h rename to include/materials/nucleation_models/NucleationMicroForceBase.h index 21f87c71435..05f2a63b032 100644 --- a/include/materials/KLBFNucleationMicroForce.h +++ b/include/materials/nucleation_models/NucleationMicroForceBase.h @@ -8,30 +8,16 @@ #include "BaseNameInterface.h" #include "DerivativeMaterialPropertyNameInterface.h" -/** - * The class implements the external driving force to recover a Drucker-Prager - * strength envelope. See Kumar et. al. https://doi.org/10.1016/j.jmps.2020.104027 for model 2020. - * all parameters are required to be Material type, not double type - */ -class KLBFNucleationMicroForce : public Material, +class NucleationMicroForceBase : public Material, public BaseNameInterface, public DerivativeMaterialPropertyNameInterface { public: static InputParameters validParams(); - - KLBFNucleationMicroForce(const InputParameters & parameters); + NucleationMicroForceBase(const InputParameters & parameters); protected: - virtual void computeQpProperties() override; - - /// Name of the external driving force - const MaterialPropertyName _ex_driving_name; - - /// The external driving force - ADMaterialProperty & _ex_driving; - - ///@{ Phase field properties + ///@{ fracture properties /// The fracture toughness const ADMaterialProperty & _Gc; /// The normalization constant @@ -45,20 +31,22 @@ class KLBFNucleationMicroForce : public Material, /// The shear modulus const ADMaterialProperty & _mu; - /// The critical tensile strength - const ADMaterialProperty & _sigma_ts; - - /// The critical compressive strength - const ADMaterialProperty & _sigma_cs; - - /// The materiel and model dependent parameter - const ADMaterialProperty & _delta; - /// The stress tensor const ADMaterialProperty & _stress; - /// Name of the stress space balance const MaterialPropertyName _stress_balance_name; /// Quantifying how far is the stress state from stress surface ADMaterialProperty & _stress_balance; + + /// Name of the external driving force + const MaterialPropertyName _ex_driving_name; + /// The external nucleation driving force + ADMaterialProperty & _ex_driving; + + /// @{ The degradation function and its derivative w/r/t damage + const VariableName _d_name; + const MaterialPropertyName _g_name; + const ADMaterialProperty & _g; + const ADMaterialProperty & _dg_dd; + /// @} }; diff --git a/src/materials/KLBFNucleationMicroForce.C b/src/materials/nucleation_models/KLBFNucleationMicroForce.C similarity index 58% rename from src/materials/KLBFNucleationMicroForce.C rename to src/materials/nucleation_models/KLBFNucleationMicroForce.C index ee769ed7ff3..dc68886a711 100644 --- a/src/materials/KLBFNucleationMicroForce.C +++ b/src/materials/nucleation_models/KLBFNucleationMicroForce.C @@ -2,66 +2,34 @@ //* being developed at Dolbow lab at Duke University //* http://dolbow.pratt.duke.edu -#include "Function.h" #include "KLBFNucleationMicroForce.h" -registerADMooseObject("raccoonApp", KLBFNucleationMicroForce); +registerMooseObjectReplaced("raccoonApp", + KLBFNucleationMicroForce, + "12/31/2024 23:59", + LDLNucleationMicroForce); InputParameters KLBFNucleationMicroForce::validParams() { - InputParameters params = Material::validParams(); - params += BaseNameInterface::validParams(); + InputParameters params = NucleationMicroForceBase::validParams(); params.addClassDescription("This class computes the external driving force for nucleation given " "a Drucker-Prager strength envelope developed by Kumar et al. (2020)"); - - params.addParam( - "fracture_toughness", "Gc", "energy release rate or fracture toughness"); - params.addParam( - "normalization_constant", "c0", "The normalization constant $c_0$"); - params.addParam( - "regularization_length", "l", "the phase field regularization length"); - - params.addParam("lambda", "lambda", "Lame's first parameter lambda"); - params.addParam("shear_modulus", "G", "shear modulus mu or G"); - params.addRequiredParam( "tensile_strength", "The tensile strength of the material beyond which the material fails."); - params.addRequiredParam( "compressive_strength", "The compressive strength of the material beyond which the material fails."); - params.addRequiredParam("delta", "delta"); - params.addParam( - "external_driving_force_name", - "ex_driving", - "Name of the material that holds the external_driving_force"); - params.addParam( - "stress_balance_name", - "stress_balance", - "Name of the stress balance function $F= \\dfrac{J_2}{\\mu} + \\dfrac{I_1^2}{9\\kappa} - c_e " - "-\\dfrac{3\\Gc}{8\\delta}=0 $. This value tells how close the material is to stress " - "surface."); - params.addParam("stress_name", "stress", "Name of the stress tensor"); return params; } KLBFNucleationMicroForce::KLBFNucleationMicroForce(const InputParameters & parameters) - : Material(parameters), - BaseNameInterface(parameters), - _ex_driving(declareADProperty(prependBaseName("external_driving_force_name", true))), - _Gc(getADMaterialProperty(prependBaseName("fracture_toughness", true))), - _c0(getADMaterialProperty(prependBaseName("normalization_constant", true))), - _L(getADMaterialProperty(prependBaseName("regularization_length", true))), - _lambda(getADMaterialProperty(prependBaseName("lambda", true))), - _mu(getADMaterialProperty(prependBaseName("shear_modulus", true))), + : NucleationMicroForceBase(parameters), _sigma_ts(getADMaterialProperty(prependBaseName("tensile_strength", true))), _sigma_cs(getADMaterialProperty(prependBaseName("compressive_strength", true))), - _delta(getADMaterialProperty(prependBaseName("delta", true))), - _stress(getADMaterialProperty(prependBaseName("stress_name", true))), - _stress_balance(declareADProperty(prependBaseName("stress_balance_name", true))) + _delta(getADMaterialProperty(prependBaseName("delta", true))) { } diff --git a/src/materials/KLRNucleationMicroForce.C b/src/materials/nucleation_models/KLRNucleationMicroForce.C similarity index 56% rename from src/materials/KLRNucleationMicroForce.C rename to src/materials/nucleation_models/KLRNucleationMicroForce.C index 8287e403901..92f825fb2c2 100644 --- a/src/materials/KLRNucleationMicroForce.C +++ b/src/materials/nucleation_models/KLRNucleationMicroForce.C @@ -2,74 +2,36 @@ //* being developed at Dolbow lab at Duke University //* http://dolbow.pratt.duke.edu -#include "Function.h" #include "KLRNucleationMicroForce.h" -registerADMooseObject("raccoonApp", KLRNucleationMicroForce); +registerMooseObjectReplaced("raccoonApp", + KLRNucleationMicroForce, + "12/31/2024 23:59", + LDLNucleationMicroForce); InputParameters KLRNucleationMicroForce::validParams() { - InputParameters params = Material::validParams(); - params += BaseNameInterface::validParams(); + InputParameters params = NucleationMicroForceBase::validParams(); params.addClassDescription("This class computes the external driving force for nucleation given " "a Drucker-Prager strength envelope developed by Kumar et al. (2022)"); - params.addParam( - "fracture_toughness", "Gc", "energy release rate or fracture toughness"); - params.addParam( - "normalization_constant", "c0", "The normalization constant $c_0$"); - params.addParam( - "regularization_length", "l", "the phase field regularization length"); - - params.addParam("lambda", "lambda", "Lame's first parameter lambda"); - params.addParam("shear_modulus", "G", "shear modulus mu or G"); - params.addRequiredParam( "tensile_strength", "The tensile strength of the material beyond which the material fails."); - params.addRequiredParam( "compressive_strength", "The compressive strength of the material beyond which the material fails."); - params.addRequiredParam("delta", "delta"); - params.addParam( - "external_driving_force_name", - "ex_driving", - "Name of the material that holds the external_driving_force"); - params.addParam( - "stress_balance_name", - "stress_balance", - "Name of the stress balance function $F= \\dfrac{J_2}{\\mu} + \\dfrac{I_1^2}{9\\kappa} - c_e " - "-\\dfrac{3\\Gc}{8\\delta}=0 $. This value tells how close the material is to stress " - "surface."); - params.addParam("stress_name", "stress", "Name of the stress tensor"); - params.addRequiredCoupledVar("phase_field", "Name of the phase-field (damage) variable"); - - params.addParam("degradation_function", "g", "The degradation function"); return params; } KLRNucleationMicroForce::KLRNucleationMicroForce(const InputParameters & parameters) - : Material(parameters), - BaseNameInterface(parameters), - _ex_driving(declareADProperty(prependBaseName("external_driving_force_name", true))), - _Gc(getADMaterialProperty(prependBaseName("fracture_toughness", true))), - _c0(getADMaterialProperty(prependBaseName("normalization_constant", true))), - _L(getADMaterialProperty(prependBaseName("regularization_length", true))), - _lambda(getADMaterialProperty(prependBaseName("lambda", true))), - _mu(getADMaterialProperty(prependBaseName("shear_modulus", true))), + : NucleationMicroForceBase(parameters), _sigma_ts(getADMaterialProperty(prependBaseName("tensile_strength", true))), _sigma_cs(getADMaterialProperty(prependBaseName("compressive_strength", true))), _delta(getADMaterialProperty(prependBaseName("delta", true))), - _stress(getADMaterialProperty(prependBaseName("stress_name", true))), - _stress_balance(declareADProperty(prependBaseName("stress_balance_name", true))), - _druck_prager_balance(declareADProperty("druck_prager_balance")), - _d_name(getVar("phase_field", 0)->name()), - _g_name(prependBaseName("degradation_function", true)), - _g(getADMaterialProperty(_g_name)), - _dg_dd(getADMaterialProperty(derivativePropertyName(_g_name, {_d_name}))) + _druck_prager_balance(declareADProperty("druck_prager_balance")) { } diff --git a/src/materials/nucleation_models/LDLNucleationMicroForce.C b/src/materials/nucleation_models/LDLNucleationMicroForce.C new file mode 100644 index 00000000000..4460be65cfa --- /dev/null +++ b/src/materials/nucleation_models/LDLNucleationMicroForce.C @@ -0,0 +1,104 @@ +//* This file is part of the RACCOON application +//* being developed at Dolbow lab at Duke University +//* http://dolbow.pratt.duke.edu + +#include "LDLNucleationMicroForce.h" + +registerADMooseObject("raccoonApp", LDLNucleationMicroForce); + +InputParameters +LDLNucleationMicroForce::validParams() +{ + InputParameters params = NucleationMicroForceBase::validParams(); + + params.addClassDescription( + "This class computes the external driving force for nucleation given " + "a Drucker-Prager strength envelope developed by Larsen et al. (2024)"); + params.addRequiredParam( + "tensile_strength", "The tensile strength of the material beyond which the material fails."); + params.addRequiredParam( + "hydrostatic_strength", + "The hydrostatic strength of the material beyond which the material fails."); + params.addParam("delta", "delta", "Name of the unitless coefficient delta"); + params.addParam("h_correction", false, "Whether to use h correction formula for delta"); + return params; +} + +LDLNucleationMicroForce::LDLNucleationMicroForce(const InputParameters & parameters) + : NucleationMicroForceBase(parameters), + _sigma_ts(getADMaterialProperty(prependBaseName("tensile_strength", true))), + _sigma_hs(getADMaterialProperty(prependBaseName("hydrostatic_strength", true))), + _delta(declareADProperty(prependBaseName("delta", true))), + _h_correction(getParam("h_correction")) +{ +} + +void +LDLNucleationMicroForce::computeQpProperties() +{ + // The bulk modulus + ADReal K = _lambda[_qp] + 2.0 * _mu[_qp] / 3.0; + + // The Young's modulus + ADReal E = 9.0 * _mu[_qp] * K / (_mu[_qp] + 3.0 * K); + + // The mobility + ADReal M = _Gc[_qp] / _L[_qp] / _c0[_qp]; + + // Invariants of the stress + ADReal I1 = _stress[_qp].trace(); + ADRankTwoTensor stress_dev = _stress[_qp].deviatoric(); + ADReal J2 = 0.5 * stress_dev.doubleContraction(stress_dev); + + // Just to be extra careful... J2 is for sure non-negative but discretization and interpolation + // might bring surprise + if (J2 < 0) + mooseException("Negative J2"); + + // define zero J2's derivative + if (MooseUtils::absoluteFuzzyEqual(J2, 0)) + J2.value() = libMesh::TOLERANCE * libMesh::TOLERANCE; + + // Compute critical energy + ADReal W_ts = _sigma_ts[_qp] * _sigma_ts[_qp] / 2.0 / E; + ADReal W_hs = _sigma_hs[_qp] * _sigma_hs[_qp] / 2.0 / K; + + // Compute delta + if (!_h_correction) + { + // Use formula without h correction + _delta[_qp] = (_sigma_ts[_qp] + (1 + 2 * std::sqrt(3)) * _sigma_hs[_qp]) / + (8 + 3 * std::sqrt(3)) / _sigma_hs[_qp] * 3.0 / 16.0 * + (_Gc[_qp] / W_ts / _L[_qp]) + + 3.0 / 8.0; + } + else + { + // Get mesh size of current element + ADReal h = _current_elem->hmin(); + + // Use formula with h correction + _delta[_qp] = std::pow(1 + 3.0 / 8.0 * h / _L[_qp], -2) * + (_sigma_ts[_qp] + (1 + 2 * std::sqrt(3.0)) * _sigma_hs[_qp]) / + (8 + 3 * std::sqrt(3.0)) / _sigma_hs[_qp] * 3 / 16 * + (_Gc[_qp] / W_ts / _L[_qp]) + + std::pow(1 + 3.0 / 8.0 * h / _L[_qp], -1) * 2 / 5; + } + + // Parameters in the strength surface + ADReal alpha_1 = + -_delta[_qp] * _Gc[_qp] / 8.0 / _sigma_hs[_qp] / _L[_qp] + 2.0 / 3.0 * W_hs / _sigma_hs[_qp]; + ADReal alpha_2 = -(std::sqrt(3.0) / 8.0 * _delta[_qp] * (3.0 * _sigma_hs[_qp] - _sigma_ts[_qp]) / + (_sigma_hs[_qp] * _sigma_ts[_qp]) * _Gc[_qp] / _L[_qp] + + 2.0 / std::sqrt(3.0) * W_hs / _sigma_hs[_qp] - + 2.0 * std::sqrt(3.0) * W_ts / _sigma_ts[_qp]); + + // Compute the external driving force required to recover the desired strength envelope. + _ex_driving[_qp] = + alpha_2 * std::sqrt(J2) + alpha_1 * I1 + + (1.0 - std::sqrt(I1 * I1) / I1) / std::pow(_g[_qp], 1.5) * + (J2 / 2.0 / _mu[_qp] + I1 * I1 / 6.0 / (3.0 * _lambda[_qp] + 2.0 * _mu[_qp])); + + _stress_balance[_qp] = + J2 / _mu[_qp] + std::pow(I1, 2) / 9.0 / K - _ex_driving[_qp] - M * _delta[_qp]; +} diff --git a/src/materials/nucleation_models/NucleationMicroForceBase.C b/src/materials/nucleation_models/NucleationMicroForceBase.C new file mode 100644 index 00000000000..53b4a556791 --- /dev/null +++ b/src/materials/nucleation_models/NucleationMicroForceBase.C @@ -0,0 +1,51 @@ +//* This file is part of the RACCOON application +//* being developed at Dolbow lab at Duke University +//* http://dolbow.pratt.duke.edu + +#include "InputParameters.h" +#include "NucleationMicroForceBase.h" + +InputParameters +NucleationMicroForceBase::validParams() +{ + InputParameters params = Material::validParams(); + params += BaseNameInterface::validParams(); + + // common parameters + params.addParam( + "fracture_toughness", "Gc", "energy release rate or fracture toughness"); + params.addParam( + "normalization_constant", "c0", "The normalization constant $c_0$"); + params.addParam( + "regularization_length", "l", "the phase field regularization length"); + params.addParam("lambda", "lambda", "Lame's first parameter lambda"); + params.addParam("shear_modulus", "G", "shear modulus mu or G"); + params.addParam( + "stress_balance_name", "stress_balance", "Name of the stress balance function $F(\\sigma)$"); + params.addParam("stress_name", "stress", "Name of the stress tensor"); + params.addParam( + "external_driving_force_name", + "ex_driving", + "Name of the material that holds the external_driving_force"); + params.addRequiredCoupledVar("phase_field", "Name of the phase-field (damage) variable"); + params.addParam("degradation_function", "g", "The degradation function"); + return params; +} + +NucleationMicroForceBase::NucleationMicroForceBase(const InputParameters & parameters) + : Material(parameters), + BaseNameInterface(parameters), + _Gc(getADMaterialProperty(prependBaseName("fracture_toughness", true))), + _c0(getADMaterialProperty(prependBaseName("normalization_constant", true))), + _L(getADMaterialProperty(prependBaseName("regularization_length", true))), + _lambda(getADMaterialProperty(prependBaseName("lambda", true))), + _mu(getADMaterialProperty(prependBaseName("shear_modulus", true))), + _stress(getADMaterialProperty(prependBaseName("stress_name", true))), + _stress_balance(declareADProperty(prependBaseName("stress_balance_name", true))), + _ex_driving(declareADProperty(prependBaseName("external_driving_force_name", true))), + _d_name(getVar("phase_field", 0)->name()), + _g_name(prependBaseName("degradation_function", true)), + _g(getADMaterialProperty(_g_name)), + _dg_dd(getADMaterialProperty(derivativePropertyName(_g_name, {_d_name}))) +{ +} diff --git a/test/tests/nucleation_models/klbf.i b/test/tests/nucleation_models/klbf.i index 047006a35b8..77f847cdb0f 100644 --- a/test/tests/nucleation_models/klbf.i +++ b/test/tests/nucleation_models/klbf.i @@ -130,6 +130,7 @@ [] [kumar_material] type = KLBFNucleationMicroForce + phase_field = d normalization_constant = c0 tensile_strength = sigma_ts compressive_strength = sigma_cs diff --git a/tutorials/surfing_boundary_problem/elasticity.i b/tutorials/surfing_boundary_problem/elasticity.i index d698ccaae2f..a64a614c447 100644 --- a/tutorials/surfing_boundary_problem/elasticity.i +++ b/tutorials/surfing_boundary_problem/elasticity.i @@ -15,9 +15,8 @@ Gc = 9.1e-2 # 91N/m l = 0.35 sigma_ts = 27 sigma_cs = 77 -# for model 2022 (KLR), select delta=0.2 -# for model 2020 (KLBF), select delta=1.16 -delta = 0.2 #1.16 +sigma_hs = '${fparse 2/3*sigma_ts*sigma_cs/(sigma_cs - sigma_ts)}' + c1 = '${fparse (1+nu)*sqrt(Gc)/sqrt(2*pi*E)}' c2 = '${fparse (3-nu)/(1+nu)}' ahead = 2 @@ -40,7 +39,7 @@ refine = 3 [fracture] type = TransientMultiApp input_files = fracture.i - cli_args = 'E=${E};K=${K};G=${G};Lambda=${Lambda};Gc=${Gc};l=${l};nx=${nx};ny=${ny};refine=${refine};sigma_ts=${sigma_ts};sigma_cs=${sigma_cs};delta=${delta}' + cli_args = 'E=${E};K=${K};G=${G};Lambda=${Lambda};Gc=${Gc};l=${l};nx=${nx};ny=${ny};refine=${refine};sigma_ts=${sigma_ts};sigma_hs=${sigma_hs}' execute_on = 'TIMESTEP_END' [] [] @@ -201,9 +200,9 @@ refine = 3 boundary = 'left top right bottom' [] [Jint_over_Gc] - type = ScalePostprocessor - value = 'Jint' - scaling_factor = '${fparse 1.0/Gc}' + type = ParsedPostprocessor + pp_names = 'Jint' + expression = 'Jint/${Gc}' [] [] diff --git a/tutorials/surfing_boundary_problem/fracture.i b/tutorials/surfing_boundary_problem/fracture.i index b9b360550fa..6eae16b1ba3 100644 --- a/tutorials/surfing_boundary_problem/fracture.i +++ b/tutorials/surfing_boundary_problem/fracture.i @@ -53,11 +53,10 @@ [Bounds] [conditional] - type = ConditionalBoundsAux + type = VariableOldValueBounds variable = bounds_dummy bounded_variable = d - fixed_bound_value = 0 - threshold_value = 0.95 + bound_type = lower [] [upper] type = ConstantBounds @@ -85,14 +84,15 @@ type = ADCoefMatSource variable = d prop_names = 'ce' + coefficient = 1.0 [] [] [Materials] [fracture_properties] type = ADGenericConstantMaterial - prop_names = 'E K G lambda Gc l sigma_ts sigma_cs delta' - prop_values = '${E} ${K} ${G} ${Lambda} ${Gc} ${l} ${sigma_ts} ${sigma_cs} ${delta}' + prop_names = 'E K G lambda Gc l sigma_ts sigma_hs' + prop_values = '${E} ${K} ${G} ${Lambda} ${Gc} ${l} ${sigma_ts} ${sigma_hs}' [] [degradation] type = PowerDegradationFunction @@ -111,21 +111,23 @@ [psi] type = ADDerivativeParsedMaterial property_name = psi - expression = 'g*psie_active+(Gc/c0/l)*alpha' + expression = 'g*psie_active+(Gc*delta/c0/l)*alpha' coupled_variables = 'd psie_active' - material_property_names = 'alpha(d) g(d) Gc c0 l' + material_property_names = 'delta alpha(d) g(d) Gc c0 l' derivative_order = 1 [] - [kumar_material] - type = KLRNucleationMicroForce + [nucleation_micro_force] + type = LDLNucleationMicroForce phase_field = d - # type = KLBFNucleationMicroForce + degradation_function = g + regularization_length = l normalization_constant = c0 + fracture_toughness = Gc tensile_strength = sigma_ts - compressive_strength = sigma_cs + hydrostatic_strength = sigma_hs delta = delta + h_correction = true external_driving_force_name = ce - output_properties = 'ce' [] [strain] type = ADComputePlaneSmallStrain diff --git a/tutorials/surfing_boundary_problem/gold/elasticity_out.e b/tutorials/surfing_boundary_problem/gold/elasticity_out.e index 77e4620ab8b..5bce570a529 100644 Binary files a/tutorials/surfing_boundary_problem/gold/elasticity_out.e and b/tutorials/surfing_boundary_problem/gold/elasticity_out.e differ diff --git a/tutorials/surfing_boundary_problem/tests b/tutorials/surfing_boundary_problem/tests index d4fafd81b8a..9c8be9e7aab 100644 --- a/tutorials/surfing_boundary_problem/tests +++ b/tutorials/surfing_boundary_problem/tests @@ -1,5 +1,5 @@ [Tests] - [kumar_fracture] + [nucforce_fracture] type = Exodiff input = 'elasticity.i' exodiff = 'elasticity_out.e'