Aurelien Doitrand
email: aurelien.doitrand@insa-lyon.fr

Gergely Molnár

Rafael Estevez

Anthony Gravouil

Strength-based regularization length in phase field fracture

Keywords: Phase field, coupled criterion, internal length

We propose a length-free (LF) implementation of phase field (PF) approach to fracture. The input parameters of the LF-PF model are the same as those of the coupled criterion (CC), namely the critical energy release rate G c and the material tensile strength 𝜎 c . This formulation relies on the previously determined correlation between the internal length in classical PF fracture models and the material tensile strength. Similar failure loads are obtained using either AT1 and AT2 versions of the LF-PF, a larger diffused damaged zone being observed using AT1. The LF-PF formulation allows retrieving the size effect for a crack under remote tensile stress in infinite medium and shows square hole specimen failure stresses close to the predictions obtained using the CC, which allows using both models in a complementary manner. Finally, the proposed formulation overcomes the problem of identifying different internal lengths for different specimen configurations.

Introduction

Crack propagation has been studied for many years in the fracture mechanics framework proposed by [START_REF] Griffith | The Phenomena of Rupture and Flow in Solids[END_REF]; [START_REF] Griffith | The theory of rupture[END_REF]. However, two main limitations of this approach can be highlighted. First, it is based on the assumption of an existing crack and enables studying its propagation but fails to predict its nucleation. Second, its implementation in a finite element (FE) code to study crack propagation either requires a priori the knowledge of the crack path or robust remeshing methods to update the crack configuration at each propagation step [START_REF] Chiaruttini | Advanced remeshing techniques for complex 3D crack propagation[END_REF][START_REF] Vattré | Singularity-free theory and adaptive finite element computations of arbitrarily-shaped dislocation loop dynamics in 3D heterogeneous material structures[END_REF]. The first limitation was overcome by [START_REF] Leguillon | Strength or toughness? A criterion for crack onset at a notch[END_REF] who developed the coupled criterion in the finite fracture mechanics (FFM) framework [START_REF] Nairn | A critical evaluation of theories for predicting microcracking in composite laminates[END_REF][START_REF] Hashin | Finite thermoelastic fracture criterion with application to laminate cracking analysis[END_REF], dedicated to study the nucleation of a crack. This approach proved to be a robust and efficient way to predict crack initiation in many configurations such as, e.g., weak or strong singularities, non-singular stress raisers [START_REF] Weißgraeber | A review of Finite Fracture Mechanics: crack initiation at singular and non-singular stress raisers[END_REF]Doitrand et al. 2020b). Initially developed under small deformation assumption and linear elastic 2D framework, it has since then been extended to 3D [START_REF] Leguillon | Cracks Pattern Formation and Spalling in Functionalized Thin Films[END_REF][START_REF] Yosibash | A 3-D failure initiation criterion from a sharp V-notch edge in elastic brittle structures[END_REF][START_REF] García | A numerical study of transverse cracking in cross-ply laminates by 3D finite fracture mechanics[END_REF]Doitrand et al. 2018a;Doitrand et al. 2018b), to consider material or geometry nonlinearities [START_REF] Leguillon | Failure initiation at V-notch tips in quasi-brittle materials[END_REF][START_REF] Li | Numerical implementation of the coupled criterion for damaged materials[END_REF][START_REF] Rosendahl | Nonlinear elastic finite fracture mechanics: modeling mixed-mode crack nucleation in structural glazing silicone sealants[END_REF]Doitrand et al. 2020c;[START_REF] Leite | Crack onset in stretched open hole PMMA plates considering linear and non-linear elastic behaviours[END_REF], as well as dynamic crack initiation [START_REF] Doitrand | Dynamic crack initiation assessment with the coupled criterion[END_REF][START_REF] Chao Correas | Finite Fracture Mechanics extension to dynamic loading scenarios[END_REF]. It also revealed efficient for small scale fracture assessment (Doitrand et al. 2020a;[START_REF] Gallo | Brittle failure of nanoscale notched silicon cantilevers: a finite fracture mechanics approach[END_REF][START_REF] Jimenez Alfaro | Finite Fracture Mechanics at the micro-scale. Application to bending tests of micro cantilever beams[END_REF]. The second limitation of Griffith's approach was solved by first reconsidering the local energy criterion as a global minimization problem [START_REF] Francfort | Revisiting brittle fracture as an energy minimization problem[END_REF], which was the first step towards a variational formulation of fracture problems. Indeed, this approach suffered from not knowing the crack path described as a surface in the volume. The introduction of a crack density function inspired from [START_REF] Mumford | Optimal approximations by piecewise smooth functions and associated variational problems[END_REF]; [START_REF] Tortorelli | Approximation of functional depending on jumps by elliptic functional via t-convergence[END_REF] combined with an internal length scale made it possible to obtain a variational formulation and thus its FE implementation (Bourdin et al. 2000b;Miehe et al. 2010a). Phase-field (PF) fracture approach has become a common method to assess fracture problems and has been improved to consider plasticity [START_REF] Ambati | Phase-field modeling of ductile fracture[END_REF][START_REF] Fang | Phase field fracture in elasto-plastic solids: Abaqus implementation and case studies[END_REF]), dynamic effects [START_REF] Borden | A phase-field description of dynamic brittle fracture[END_REF]Molnár et al. 2020a), fatigue [START_REF] Lo | A phase-field model for fatigue crack growth[END_REF][START_REF] Mesgarnejad | Phase-field models for fatigue crack growth[END_REF]), interfacial damages [START_REF] Nguyen | Phase field modeling of interfacial damage in heterogeneous media with stiff and soft interphases[END_REF], hydrogen assisted cracking [START_REF] Martínez-Pañeda | A phase field formulation for hydrogen assisted cracking[END_REF], hydraulic fractures [START_REF] Wilson | Phase-field modeling of hydraulic fracture[END_REF][START_REF] Xia | Phase field modeling of hydraulic fracturing with interfacial damage in highly heterogeneous fluid-saturated porous media[END_REF][START_REF] Zhou | Phase field method for quasi-static hydro-fracture in porous media under stress boundary condition considering the effect of initial stress field[END_REF] or phase transformation-induced fracture [START_REF] Wang | Capturing the stress evolution in electrode materials that undergo phase transformations during electrochemical cycling[END_REF][START_REF] Djeumen | Modeling diffusive phase transformation and fracture in viscoplastic materials[END_REF].

The main idea of PF approach consists in approximating the sharp crack discontinuity by a smeared damage field description through the use of an internal length (ℓ c ). This length controls the magnitude of the damage diffusion. Initially, ℓ c was used to facilitate the numerical solution and avoid mesh dependence of the crack path, with the idea of reducing ℓ c to 0 to retrieve the original Griffith theory, thus simply considering ℓ c as a numerical parameter without physical meaning. Some authors considered this parameter as a material internal length that must be identified (Pham et al. 2010a;Pham et al. 2010b;[START_REF] Freddi | Regularized variational theories of fracture: A unified approach[END_REF][START_REF] Pham | Gradient damage models and their use to approximate brittle fracture[END_REF][START_REF] Del Piero | A variational approach to fracture and other inelastic phenomena[END_REF], for instance based on the failure loading measured experimentally (Nguyen et al. 2016a;[START_REF] Tanné | Crack nucleation in variational phase-field models of brittle fracture[END_REF]Kumar et al. 2020). Some authors also tried to directly measure the PF internal length based on microtomograph (Nguyen et al. 2016a;Nguyen et al. 2016b) or fractography [START_REF] Darban | Experimental Determination of the Length-Scale Parameter for the Phase-Field Modeling of Macroscale Fracture in Cr-Al2O3 Composites Fabricated by Powder Metallurgy[END_REF]) observations of the damage region extent.

The CC and PF approaches may somehow be related since they can be used in order to model the same fracture problems. Both approaches were for instance used for Hertzian indentation induced fracture [START_REF] Strobl | Analysis of Hertzian indentation fracture in the framework of finite fracture mechanics[END_REF][START_REF] Strobl | Phase field modeling of Hertzian indentation fracture[END_REF], notched thin ply laminates fracture [START_REF] Reinoso | Strength prediction of notched thin ply laminates using finite fracture mechanics and the phase field approach[END_REF] or notched specimens under bending [START_REF] Abaza | Prediction of crack nucleation and propagation in porous ceramics using the phase-field approach[END_REF]. Based on the main idea proposed by [START_REF] Leguillon | Strength or toughness? A criterion for crack onset at a notch[END_REF], i.e., considering not only the material critical energy release rate but also its strength, Kumar et al. (2020) recently proposed to revisit nucleation in the PF approach. Both the CC and PF approaches involve the material critical energy release rate, therefore, the main difference in the input parameters of both methods is the tensile strength in the CC, replaced by the internal length in PF approach. An analytical relation between the internal length and the material tensile strength can actually be derived based on the PF homogeneous solution under uniaxial loading, thus neglecting the damage gradient [START_REF] Benallal | Bifurcation and stability issues in gradient theories with softening[END_REF][START_REF] Amor | Regularized formulation of the variational brittle fracture with unilateral contact: Numerical experiments[END_REF][START_REF] Pham | Gradient damage models and their use to approximate brittle fracture[END_REF][START_REF] Borden | A phase-field description of dynamic brittle fracture[END_REF][START_REF] Kuhn | On degradation functions in phase field fracture models[END_REF][START_REF] Mesgarnejad | Validation simulations for the variational approach to fracture[END_REF]Nguyen et al. 2016a;[START_REF] Pham | Experimental validation of a phase-field model for fracture[END_REF][START_REF] Zhang | Numerical evaluation of the phase-field model for brittle fracture with emphasis on the length scale[END_REF][START_REF] Tanné | Crack nucleation in variational phase-field models of brittle fracture[END_REF]Kumar et al. 2020;Kumar et al. 2020b). This relation can actually be extended to any multiaxial loadings, which results in a surface describing the correlation between 𝜎 c and ℓ c , depending on the local principal stress ratio and the Poisson's ratio (Molnár et al. 2020b). The correlation between 𝜎 c and ℓ c was further investigated by confronting the CC and PF approaches to mode I crack propagation, shear fracture and crack arrest configuration (Molnár et al. 2020b). Interestingly, it was shown that similarly to the CC, the PF actually followed the principle of fulfilling both stress and energy criteria. Both approaches provided similar quantitative and qualitative results, which enabled deriving a unique correlation between the internal length, the tensile strength and the initiation length.

It was thus highlighted that for a given geometry and loading configuration, it is possible to identify ℓ c that corresponds to a given value of 𝜎 c (and vice-versa). Neverthless, [START_REF] Abaza | Prediction of crack nucleation and propagation in porous ceramics using the phase-field approach[END_REF] recently observed that two significanly different internal lengths had to be identified to capture the initiation loading level predicted by the CC for two different ceramic specimens containing either a V-notch (identified ℓ c = 0.02 mm) or a crack blunted by a hole (identified ℓ c = 0.007 mm). Therefore, it is not straithforward to determine a unique internal length for a given material and it raised the question of which internal length should be chosen if the configuration under investigation contains both a V-notch and a crack blunted by a hole. The objective of this work is to overcome this problem and propose a PF formulation in which the input parameters are the critical energy release rate G c and the tensile strength 𝜎 c based on the previously derived correlation between ℓ c and 𝜎 c , which is recalled in Section 2. The PF formulation is then presented in Section 3 and illustrated on different examples in Section 4.

Correlation for the homogeneous solution

We first recall the CC and PF approaches which enabled deriving a correlation between the PF internal length ℓ c and the material tensile strength 𝜎 𝑐 .

The coupled criterion

The main idea of the CC relies on two necessary conditions that must be simultaneously fulfilled in order to make the nucleation of a crack possible. The first condition invokes a sufficiently high stress just before initiation over the whole area on which a crack initiates. It thus consists in comparing the tensile stress along the crack path before initiation to the material tensile strength. The second condition results from a balance of the energies between and after crack initiation, which results in comparing the incremental energy release rate (defined as the opposite of the variation of potential energy divided by the finite initiation crack increment) to the material critical energy release rate. Solving the CC reverts to determining the minimum imposed loading for which both conditions are fulfilled, which also allows determining the initiation crack length. The input parameters of the CC are the elastic properties of the material, the critical energy release rate G c and the material tensile strength 𝜎 c . The outputs of the model are the critical loading at crack initiation and the initiation crack length. More details about the CC and its implementation can be found in [START_REF] Leguillon | Strength or toughness? A criterion for crack onset at a notch[END_REF]Doitrand et al. 2020b;Molnár et al. 2020b).

Correlation between PF and CC

The PF approach shares some similar input parameters with the CC, namely the elastic properties of the material and the critical energy release rate G c . However, the material tensile strength does not intervene in the classical PF formulation, it is replaced by an internal length ℓ c . More details about the PF formulation and its FE implementation can for instance be found in (Molnár et al. 2020b;[START_REF] Molnár | Thermodynamically consistent linear-gradient damage model in Abaqus[END_REF] as well as in Section 3. Similarly to the CC, the PF approach also enables determining the critical loading at crack initiation as well as crack geometry (Molnár et al. 2020b), which was the starting point to find a correlation between 𝜎 c and ℓ c . Based on the homogeneous solution and different benchmark examples, we previously showed (Molnár et al. 2020b) that the correlation between 𝜎 c and ℓ c depends on the Poisson's ratio 𝜈 and principal stress ratio:

𝜎 c = 𝜂 (𝜈, 𝜎 𝐼 𝐼 𝜎 𝐼 ) √ 𝐸G c ℓ c , (1)
where 𝜂 can somehow be understood as a normalized tensile strength which varies depending on the local stress state and on the Poisson's ratio. This relation equivalently writes:

ℓ c = 𝜂 (𝜈, 𝜎 𝐼 𝐼 𝜎 𝐼 ) 2 𝐸G c 𝜎 2 c = 𝜂 (𝜈, 𝜎 𝐼 𝐼 𝜎 𝐼 ) 2 ℓ mat . (2)
The methodology to compute 𝜂 and its numerical implementation are provided in (Molnár et al. 2020b). As a matter of example, the particular values of 𝜂 under uniaxial tension (𝜎 𝐼 𝐼 /𝜎 𝐼 = 0) for AT1 and AT2 implementations are respectively 1 shows the variation of the normalized tensile strength 𝜂 as a function of 𝜎 𝐼 𝐼 /𝜎 𝐼 for several values of 𝜈 obtained using the AT1 or AT2 implementation of the PF approach. For a given stress state and tensile strength, AT1 implementation results in a larger regularization length than for AT2. This correlation is of main interest when implementing a PF approach for a material whose tensile strength is known, or equivalently determine the tensile strength in a configuration for which ℓ c is known. For configurations in which only one crack nucleates and propagates under a non-evolutive loading, it enables making a clear link between the CC and the PF approaches. It also leads to some issues when considering a structure made of the same material (described by fracture properties G c and 𝜎 c ), in which several cracks may develop from different stress raisers or singularities with different local stress states. From the correlation determined above, it would require the use of several internal lengths to be adapted at each possible crack nucleation locations. These observations thus motivates the implementation of a strength-based PF formulation in which the inputs are no longer ℓ c and G c but 𝜎 c and G c , the introduction of 𝜎 c in the PF formulation being based on the above-mentioned correlation to ℓ c . Figure 2 shows the failure enveloppe corresponding to the homogeneous solution obtained using either AT1 PF implementation (ℓ 𝑐 = 0.2 mm) or the strength-based PF implementation. Under homogeneous stress state, the strength-based PF implementation actually reverts to a Rankine-like fracture enveloppe, which means that for given Poisson's ratio and local stress state, the regularization length is adapted so that failure occurs for 𝜎 = 𝜎 c .

𝜂 𝐴𝑇 1 (𝜈, 0) = √ 3 8(1-𝜈 2 ) and 𝜂 𝐴𝑇 2 (𝜈, 0) = √ 27 256(1-𝜈 2 ) . Figure
-10 -8 -6 -4 -2 0 0.2 0.4 0.6 0.8 1 1.2 𝜎 𝐼 𝐼 /𝜎 𝐼 𝜂 𝐴𝑇 1 𝜈 = 0 𝜈 = 0.1 𝜈 = 0.2 𝜈 = 0.3 𝜈 = 0.4 𝜈 = 0.45 (a) -10 -8 -6 -4 -2 0 0.1 0.2 0.3 0.4 0.5 0.6 𝜎 𝐼 𝐼 /𝜎 𝐼 𝜂 𝐴𝑇 2 𝜈 = 0 𝜈 = 0.1 𝜈 = 0.2 𝜈 = 0.3 𝜈 = 0.4 𝜈 = 0.45 (b)
𝜎 𝐼 𝐼 (MPa) PF-ℓ 𝑐 , 𝜈 = 0 PF-ℓ 𝑐 , 𝜈 = 0.3 PF-𝜎 𝑐 , 𝜈 = 0 PF-𝜎 𝑐 , 𝜈 = 0.3

Phase field formulation

The proposed length free phase field for fracture (LF-PF) formulation is expected to be reliable in homogeneous stress configurations, since it will be constructed based on a 𝜎 cℓ c correlation constructed from the PF homogeneous solution. It is also expected to be reliable in configurations driven by energy (for instance crack length large compared to the internal length/material characteristic length, (Molnár et al. 2020b)) since in these configurations the failure load and crack configuration will mainly depend on G c . In intermediate configurations, where fracture is driven both by stress and energy, it should incorporate both aspects as shown in (Molnár et al. 2020b).

Internal length free formulation

We consider 2D phase field implementation for fracture under quasi-static conditions and small deformation assumption. The energy functional that must be minimized to determine the displacements (𝒖) and the phase-field (𝑑) writes:

L = Π int (𝒖, 𝑑) + 𝑊 (𝑑) -Π ext , (3)
with the internal and external strain energies defined as follows:

Π int (𝒖, 𝑑) = ∫ Ω 𝜓 (𝜀, 𝑑) 𝑑 Ω, Π ext = ∫ Ω 𝜸 • 𝒖𝑑𝑉 + ∫ 𝜕Ω 𝒕 • 𝒖𝑑𝐴.
(4)

version of November 10, 2022

where 𝜸 and 𝒕 are external volumetric and boundary forces. The strain energy density (𝜓 ) depends on the phase field variable 𝑑 varying between 0 (undamaged state) and 1 (crack represented by a total material stiffness and resistance loss) so that:

𝜓 (𝜺, 𝑑) = 𝑔 (𝑑) 𝜓 + 0 (𝜺) + 𝜓 - 0 (𝜺) . (5)
The elastic strain energy is split into positive (tensile) and negative (compression) components and only the positive part of the elastic strain energy is damaged, which ensures avoiding damage due to compression. Therefore, the material stiffness and the elastic stress are reduced due to increasing damage variable only when the material is under tension. The degradation function is defined as: 𝑔 (𝑑) = (1 -𝑑) 2 + 𝜉, where 𝜉 = 10 -12 ensures a good solution conditioning. The decomposition is based on the principal strains (𝜀 𝑖 ) and on the trace of the strain tensor:

       𝜓 + 0 ( ε) = 𝑖 𝜇⟨𝜀 𝑖 ⟩ 2 + + 𝜆 2 ⟨𝑡𝑟 (𝜺)⟩ 2 + , 𝜓 - 0 ( ε) = 𝑖 𝜇⟨𝜀 𝑖 ⟩ 2 -+ 𝜆 2 ⟨𝑡𝑟 (𝜺)⟩ 2 -. (6)
with 𝜆 and 𝜇 the Lamé parameters. 𝑊 (𝑑) represents the energy dissipated by opening the crack, described by a diffuse representation of the localized discontinuity. In classical implementation of phase-field approach for fracture, it writes:

𝑊 (𝑑) = G c Γ (𝑑) = ∫ Ω G c 𝛾 (𝑑, ∇𝑑) 𝑑 Ω = ∫ Ω G c 𝑐 𝜔 ℓ c 𝜔 (𝑑) + ℓ 2 c |∇𝑑 | 2 𝑑 Ω, ( 7
)
where Γ is the overall crack surface, 𝛾 the crack surface density and G c is the material critical energy release rate. The crack surface density induces a smeared representation of a sharp crack.

It is usually calculated using the internal length ℓ c , which controls the extent of the damaged zone. We recalled in Section 2 the correlation between the internal length and the material tensile strength that was derived based on a confrontation between PF and CC. We now exploit this correlation (Equation ( 2)) in order to reformulate the crack surface density so that it no longer depends on an input internal length:

𝛾 (𝑑, ∇𝑑) = 𝜎 2 c 𝑐 𝜔 𝜂 2 𝐸G c 𝜔 (𝑑) + 𝜂 2 𝐸G c 𝜎 2 c 2 |∇𝑑 | 2 . ( 8
)
Based on this formulation, the crack surface density only depends on the material properties (Young's modulus, Poisson's ratio, critical energy release rate and tensile strength) and the local stress state through 𝜂. As a consequence, the total energy dissipated by the crack opening writes:

𝑊 (𝑑) = ∫ Ω 𝜎 2 c 𝑐 𝜔 𝜂 2 𝐸 𝜔 (𝑑) + 𝜂 2 𝐸G c 𝜎 2 c 2 |∇𝑑 | 2 𝑑 Ω, (9)
where 𝜔 (𝑑) is a geometric function that controls the topology of the phase field distribution and 𝑐 𝜔 is a normalization constant so that the total crack surface density function is consistent with the theoretical value of the corresponding discrete crack representation:

𝑐 𝜔 = 4 ∫ 1 0 √ 𝜔 (𝑠)𝑑𝑠. ( 10
)
We consider two models classically used, namely AT1 [START_REF] Pham | Gradient damage models and their use to approximate brittle fracture[END_REF]) and AT2 (Bourdin et al. 2000a;Miehe et al. 2010a), which differs in the crack surface density through the definition of 𝜔 (𝑑):

𝜔 AT1 (𝑑) = 𝑑, 𝜔 AT2 (𝑑) = 𝑑 2 , (11)
with the corresponding normalization constants 𝑐 AT1 𝜔 = 8/3 and 𝑐 AT2 𝜔 = 2. Damage irreversibility is ensured so that 𝑑 > 0.

Numerical implementation

Solving the fully coupled problem in a monolithic way results in convergence issues, which can be avoided by decoupling both problems and successively solving the uncoupled mechanical (considering fixed phase field) and phase-field (considering fixed displacements) problems, respectively:

𝒖 𝑛+1 = Arg inf 𝒖 ∫ Ω 𝜓 (𝒖, 𝑑 𝑛 ) -γ • 𝒖 𝑑𝑉 - ∫ 𝜕Ω t • 𝒖𝑑𝐴 , ( 12
)
𝑑 𝑛+1 = Arg inf 𝑑 ∫ Ω G 𝑐 𝛾 (𝑑, ∇𝑑) + (1 -𝑑) 2 𝜓 + 0 (0, 𝑛) 𝑑𝑉 , (13)
Damage irreversibility is enforced by using Lagrange multipliers corresponding to the phase field variable increment 𝑓 (𝑑 n+1 ) = 𝑑 n -𝑑 n+1 ⩽ 0 and thus modifying the energy functional

L = Π int (𝒖, 𝑑) + 𝑊 (𝑑) -Π ext + ∑ 𝑗={𝑑 n >𝑑 n+1 } 𝜆 𝑗 𝑓 𝑗 (𝑑) , (14)
where 𝑗 denotes the active constraints and 𝜆 𝑗 are the Lagrange multipliers. More details about the implementation of AT1 and AT2 models can be found in previous papers [START_REF] Molnár | 2D and 3D Abaqus implementation of a robust staggered phase-field solution for modeling brittle fracture[END_REF]Molnár et al. 2020b;[START_REF] Molnár | Thermodynamically consistent linear-gradient damage model in Abaqus[END_REF].

The numerical implementation of LF-PF formulation is actually based on these previous implementations. The main difference of the proposed implementation is due to the change in input parameters since the internal length ℓ c is replaced by the material tensile strength 𝜎 c . Therefore, the internal length is no longer fixed but calculated at each integration point using Equation ( 2), which requires the knowledge of 𝜂. For the sake of simplicity, 𝜂 is computed a priori as a function of the principal stress ratio 𝜎 𝐼 𝐼 /𝜎 𝐼 and the Poisson' ratio. Then, at each iteration, the value of 𝜂 is adapted at each integration points to calculate the crack surface density function (Equation ( 8)) and thus the total energy dissipated by the crack opening (Equation ( 9)). A constant value of 𝜂 is finally set as soon as the phase field at the integration point becomes larger than 0.99. The proposed LF-PF formulation finally reverts to using the conventional phase-field scheme in which the regularization length at each integration point is locally defined based on the 𝜂 function in (Equation (2)).

Examples

Crack in infinite media under tension

We first evaluate the ability of the LF-PF implementation to reproduce the size effect for a crack in infinite media under remote tensile stress (Figure 3(a)). We consider isotropic homogeneous material (Young's modulus: 𝐸 = 3 GPa, Poisson's ratio 𝜈 = 0.37 and critical energy release rate G c = 300 J/m 2 , typical properties of brittle polymers). The size effect can thus be studied by varying the normalized initial crack length 𝑎/ℓ mat , i.e., by modifying either the initial crack length or the material tensile strength. Figure 3(b) shows the variation of the normalized critical stress as a function of the normalized initial crack length obtained with the LF-PF implementation and compared to the one obtained either with classical AT2 implementation or with the CC (Molnár et al. 2020b). The LF-PF implementation is able to capture the transition between a regime following Griffith's solution for sufficiently large normalized initial crack lengths to a regime governed by the material tensile strength for small normalized initial crack lengths. Indeed, when the initial crack length becomes small enough compared to the material characteristic length ℓ mat , the material behaves like if there was no defect and the critical stress tends towards the material tensile strength, hence retrieving the homogeneous solution. Moreover, when the initial crack length is sufficiently large compared to ℓ mat , we retrieve the Griffith's configuration of a semi-infinite crack in infinite media. In this configuration, the critical stress is mainly controlled by the material critical energy release rate. It can be observed that the critical stress variation starts to deviate from Griffith's solution when the initial crack length becomes smaller than ℓ mat .

Square hole specimens under tension

The second example concerns mixed mode crack initiation at square hole corner under tensile loading (Figure 4) [START_REF] Doitrand | Experimental and theoretical characterization of mixed mode brittle failure from square holes[END_REF]). The material under investigation is PMMA, considered as isotropic and homogeneous (Young's modulus: 𝐸 = 2.7 GPa, Poisson's ratio 𝜈 = 0.39 and critical energy release rate G c = 92 J/m 2 ). More details about experiments and CC implementation can be found in [START_REF] Doitrand | Experimental and theoretical characterization of mixed mode brittle failure from square holes[END_REF]. The mesh size near the crack initiation location is set to at most 𝑐/250, where 𝑐 is the hole side, differences on the crack initiation loading smaller than 1% are obtained using finer meshes. Figure 5 shows the principal stress ratio distribution before crack initiation and after a certain propagation for 𝑐 = 15.34 mm square hole side. This quantity is of primary importance for the LF-PF formulation since it enables calculating the coefficient 𝜂 required for the crack surface density calculation (Equation ( 8)). The principal stress ratio varies between 0 and 0.7 before crack initiation (Figure 5a-b) and between 0 and 0.55 after propagation (Figure 5c-e). As a consequence, 𝜂 varies between 0.433 and 0.457 before initiation (respectively 0.458 after crack propagation). It results in variations of the corresponding internal length calculated based on the homogeneous solution through Equation (2) between 25.4 µm and 28.2 µm before crack initiation (respectively 28.3 µm after crack propagation). Therefore, the LF-PF formulation can be understood as adapting the internal length locally depending on the principal stress ratio so that the tensile strength of the material is homogeneous and equal to 𝜎 c . Figure 6 shows the phase field distribution after a certain crack propagation and at final failure of the specimens obtained using either AT1 or AT2 versions of the LF-PF implementation. The final crack path predicted by both model versions are similar and can almost be superimposed. However, a noticeable difference concerns the extent of the damage zone, which is larger for AT1 than for AT2 version. This can be explained by the fact that both model versions have the same input parameters, namely G c and 𝜎 c . Therefore, based on Equation (1) and Figure 1, for a given principal stress ratio, the coefficient 𝜂 is larger for AT1 than for AT2, which results in a larger corresponding internal length (Equation ( 2)). Damage widening is observed near the specimen free edge due to the interaction between the crack tip and the free edge, which provides a local stress state corresponding to a larger local 𝜂 value and thus a wider damage zone. The LF-PF formulation is now compared to the results obtained using the coupled criterion [START_REF] Doitrand | Experimental and theoretical characterization of mixed mode brittle failure from square holes[END_REF]. Figure 7 shows the failure stress obtained using the CC, AT1 or AT2 versions of the LF-PF formulation. Results are presented for the LF-PF formulation considering or not a 𝑑 = 1 Dirichlet damage boundary condition on the node located at the V-notch tip. Both the CC and LF-PF are able to capture the failure force decrease with increasing hole size observed experimentally [START_REF] Doitrand | Experimental and theoretical characterization of mixed mode brittle failure from square holes[END_REF]. Moreover, similar failure stress decreases are obtained using either AT1 or AT2 version of LF-PF. This property results from the fact that classical AT1 or AT2 versions give two different correlations between ℓ c and 𝜎 c . The LF-PF formulation uses the coefficient 𝜂 AT1 or 𝜂 AT2 derived from these correlations so that the tensile strength of the material is 𝜎 𝑐 for both models. Considering a Dirichlet boundary condition at the V-notch tip results in a decrease in the failure stress and results in a slightly better agreement with the CC, which might not be surprising since Dirichlet boundary conditions on the initial crack nodes were also used in order to determine the correlation between 𝜎 c and ℓ c (Molnár et al. 2020b).

𝐿 SH 𝑐 𝑤 SH 𝜎 0

Abaza's benchmark

The last example highlights the interest of using a PF formulation with the same input parameters to assess failure in specimens made of the same material and containing either a singularity (V-notch) or a stress concentration (crack blunted by a hole) (Figure 8) [START_REF] Abaza | Prediction of crack nucleation and propagation in porous ceramics using the phase-field approach[END_REF]. Two 2022) highlighted that using a classical PF implementation, two significantly different internal lengths had to be chosen in order to capture the crack initiation loading levels for these two configurations, which thus raises the question of which internal lengths should be chosen if we considered a specimen containing both a V-notch and a crack blunted by a hole. A possibility would be to adapt locally the internal length depending on the local singularity or stress concentration. This is somehow what is intrinsically done when using the LF-PF model through varying the coefficient 𝜂 based on the variation of principal stress ratio. The mesh size near the crack initiation location is set to 𝑎 𝑛 /400, where 𝑎 𝑛 is the notch depth (Figure 8), differences on crack initiation loading smaller than 1% are obtained using finer meshes. The load-displacement curves corresponding to the studied cases obtained with the LF-PF model are linear up to failure that corresponds to a sudden force drop to zero. The apparent stress intensity factor is calculated based on the maximum force before failure. good agreement between the apparent generalized stress intensity factors obtained using the CC and the LF-PF model, especially when if the Dirichlet damage boundary condition is used. Contrary to the classical PF implementation, the LF-PF model enables simulating crack initiation and propagation in two different configurations using the same input parameters. Moreover, the relevance of the LF-PF formulation with respect to the CC is highlighted on two different materials. The principal stress ratio and damage distribution obtained before initiation and after a certain crack propagation are shown in Figure 11 for 3YZ specimen containing a 𝛽 = 90 deg V-notch angle. Before crack initiation, the principal stress ratio varies between 0 and 0.76 for the specimen containing the blunted hole and between 0 and 0.68 for the specimen containing the V-notch. The corresponding range of values for 𝜂 are between 0.36 and 0.38 (specimen with blunted crack) and between 0.37 and 0.38 (specimen with V-notch). The corresponding range of internal lengths calculated based on the homogeneous solution through Equation ( 2) is between 1.67 mm and 1.87 mm (specimen with blunted crack) and between 1.75 mm and 1.87 mm (specimen with V-notch). After propagation to a certain crack length, the principal stress ratio are relatively similar between both cases and varies between 0 and 0.99. The corresponding range of values for 𝜂 are between 0.33 and 0.38. The corresponding range of internal lengths calculated based on the homogeneous solution through Equation ( 2) is between 1.37 mm and 1.87 mm.

Discussion

PF approaches for fracture were originally built in order to retrieve Griffith's approach [START_REF] Griffith | The Phenomena of Rupture and Flow in Solids[END_REF][START_REF] Griffith | The theory of rupture[END_REF] when the regularization length tends towards zero. Using an internal length was thus a convenient way to define crack surface densities and propose a variational formulation and associated FE implementation of fracture problems. Griffith's approach is dedicated to study the propagation of an existing crack, but actually fails to predicts its initiation at singular point or stress raiser for instance, generally because the energy release rate tends towards zero when the crack length tends towards zero. The CC overcomes this problem by combining a stress criterion to the energy criterion, which enables assessing both crack initiation and propagation. Indeed, the CC reverts to Griffith's approach when applied to assess crack propagation. It requires as input parameters G c and 𝜎 c , the material tensile strength. We previously showed that there exists a correlation between 𝜎 c and ℓ c (Molnár et al. 2020b). Basically, increasing 𝜎 c in the CC reverted to decreasing ℓ c in PF, following a relation involving both the material properties and the local stress state. Therefore, considering ℓ c that tends towards zero reverts to consider 𝜎 c that tends towards infinity (see Equation ( 2)) thus retrieving the discrete solution of Griffith. From the CC point of view, of course it still enables dealing with crack propagation in Griffith-like configurations, i.e., configurations only driven by energy such as a semi-infinite crack in infinite medium. In such configurations, any finite values of 𝜎 c may be used, and similarly the influence of the internal length magnitude is negligible. For Griffith-like configurations, the internal length thus becomes a numerical parameter that is convenient for the numerical implementation but that does not play a significant role from a physical point of view. However, considering ℓ c that tends towards zero prevents assessing Griffith-unlike configurations for which a stress criterion is also needed, such as crack initiation.

Instead of trying to only recover Griffith-like configurations of crack propagation, PF models are actually able to retrieve both crack propagation and initiation provided ℓ c is chosen to obtain the maximum stress value of 𝜎 c under uniaxial homogeneous tension. As a consequence, ℓ c is not only convenient from a numerical point but also has a clear physical meaning from a macroscopic point of view: It describes the relation between the material properties (Irwin length and Poisson's ratio) and the local stress state. This relation motivated the proposition of the LF-PF model with physically based input parameters in which ℓ c is locally adapted based on this relation. A conclusion that naturally arises from this reasoning is that the influence of the regularization length in PF models will be highlighted only in Griffith-unlike configurations, such as crack initiation, short initial crack or several cracks interacting.

A major consequence is that the internal length is not related to the material microstructure or elastic representative volume element (RVE) and can even be very large compared to the RVE size. With the proposed LF-PF approach, ℓ c does not appear explicitely as an input parameter. However, a local value at each material point can be calculated based on Equation (2), resulting in ℓ c values and diffused damage zone that are neither spatially nor temporally constant depending on the stress state time and space variation.

Conclusion

In complement of the material critical energy release rate G c , the LF-PF formulation requires as input the material tensile strength 𝜎 c instead of the internal length classically used in PF models. The LF-PF model overcomes the problem of identifying different internal lengths for different testing configurations for the same material, since there is no need to choose or identify an internal length in this formulation. The LF-PF implementation ressembles the classical PF implementation except that the internal length is replaced by an expression involving the material characteristic length and the normalized tensile strength 𝜂 that depends on the Poisson's ratio and the local principal stress ratio. The normalized tensile strength is determined based on the homogeneous solution to provide a correlation between the internal length and the material tensile strength. Similar results are obtained using either the AT1 or AT2 version of the LF-PF. The LF-PF formulation gives results close to that obtained with the coupled criterion, which means that it is able to assess crack initiation in a similar way as the coupled criterion, still conserving the possibility to deal with subsequent crack propagation. In some sense, the LF-PF can thus be considered as a phase field implementation of the coupled criterion, both models can be used in a complementary manner since their input parameters are the same and that they provide similar results. The proposed approach could be generalized to the 3D case, the main difference being the calculation of the correlation between ℓ c and 𝜎 c that must be based on the 3D homogeneous solution (provided in [START_REF] Molnár | Thermodynamically consistent linear-gradient damage model in Abaqus[END_REF]) instead of the 2D homogeneous solution.

Figure 1 :

 1 Figure 1: Variation of the normalized tensile strength 𝜂 as a function of the principal stress ratio 𝜎 𝐼 𝐼 /𝜎 𝐼 obtained for a) AT1 and b) AT2 implementation of the classical PF approach.

Figure 2 :

 2 Figure 2: Homogeneous solution failure surface in the space of principal stresses obtained for 𝐸 = 3 GPa, G 𝑐 = 300 J/m 2 , 𝜈 = 0 or 𝜈 = 0.3, using either AT1 implementation (ℓ 𝑐 = 0.2 mm) or the proposed strength-based PF formulation (𝜎 𝑐 = 41 MPa for 𝜈 = 0 and 𝜎 𝑐 = 46 MPa for 𝜈 = 0.3).

Figure 3 :

 3 Figure 3: (a) Infinite media containing a crack of length 2𝑎 subjected to a remote stress 𝜎 0 and (b) variation of the critical stress at crack propagation normalized by the material tensile strength as a function of the initial crack length normalized by the material characteristic length obtained using either the CC, the classical PF or the LF-PF model. The blue dashed line indicates Griffith's solution.

Figure 4 :

 4 Figure 4: Square hole specimen under tensile loading in which crack initiation occurs at the square corners.

Figure 5 :

 5 Figure 5: Principal stress ratio 𝜎 𝐼 𝐼 /𝜎 𝐼 distribution around the square hole corner (a-b) before crack initiation and (c,e) after crack propagation on a certain distance and (d) corresponding damage distribution (the superimposed black zone in (e) corresponds to an area where damage is larger than 0.8. A 𝑑 = 1 Dirichlet boundary condition is prescribed on the V-notch tip node in this example.

Figure 6 :Figure 7 :

 67 Figure 6: Damage distribution around the square hole in the specimen obtained using (left) AT1 or (right) AT2 versions of the LF-PF model. The final crack path obtained with AT2 version is compared to an experimental observation of the specimen after failure (Doitrand et al. 2021).

Figure 8 :

 8 Figure 8: Abaza's benchmark geometries under investigation: (a) three-point bending specimen containing a V-notch and (b) three point bending specimen containing a crack blunted by a hole.

Figure 9 :Figure 10 :

 910 Figure 9: Apparent (a) stress intensity factor variation as a function of the hole radius and (b) generalized stress intensity factor variation as a function of the V-notch angle obtained using either the CC or the LF-PF model with or without Dirichlet damage boundary condition for 3YSZ ceramic.

Figure 11 :

 11 Figure 11: Principal stress ratio 𝜎 𝐼 𝐼 /𝜎 𝐼 distribution around (a,c) the hole and (d,f) the V-notch (a,d) before crack initiation and (c,f) after crack propagation on a certain distance and (b,e) corresponding damage distribution (the black zone in (c) and (e) corresponds to an area where damage is larger than 0.8.

version of November 10, 2022

Authors' contributions Aurélien Doitrand performed numerical simulations and drafted the manuscript. Gergely Molnár helped with the numerical implementation. All authors developed the methodology, conceived of the study, and participated in its design, coordination, and critical review of the manuscript. All authors read and approved the final manuscript.