Coupled Dislocations and Fracture dynamics at finite deformation: model derivation, and physical questions

A continuum mechanical model of coupled dislocation based plasticity and fracture at finite deformation is proposed. Motivating questions and target applications of the model are sketched.


Introduction
Based on experience and insights gathered from the partial differential equation (PDE) based modeling of dislocation dynamics in Zhang et al. (2015); Arora and Acharya (2020); Arora et al. (2020); Garg et al. (2015) and fracture (Acharya 2018(Acharya , 2020;;Morin and Acharya 2021), a coupled model of fracture and dislocation based plasticity at finite deformation is explored.Even though plasticity, whether fundamentally rooted in the mechanics of dislocations or in the phenomenology of slip, and fracture are much studied subjects, e.g.(Freund (1998); Hutchinson (1979); Hirth and Lothe (1982); Bulatov and Cai (2006); Asaro (1983); Havner (1992), and the literature reviews in the papers mentioned above), to our knowledge, a full-blown continuum PDE model for their coupled mechanics does not exist and can be useful in the understanding of the deformation, flow, and fracture of solids (e.g., metals or glaciers), and the mutual interactions of these phenomena as, e.g., addressed in the seminal work (Rice and Thomson 1974).While a full thermodynamically consistent model is presented, it is recognized that this is merely a beginning that sets the stage for future computation and analysis of a simply stated, but intricate, nonlinear model which is expected to have some bearing on its target applications.
An outline of the paper is as follows: in "Governing equations: mechanics" section the mechanical equations of the model are presented.In "Guidance for constitutive assumptions from the second law of thermodynamics" section a possible set of thermodynamically consistent constitutive equations are proposed.In "Motivating questions for the development of the model" section some target problems motivating the development of the theory are sketched.It is understood that most of the questions posed are beyond the reach of rigorous methods of PDE analysis, but it is felt that demonstrating a dynamical theoretical setup where such questions can at least begin to be clearly posed and, consequently, at least be approached with finite-dimensional methods of approximation and rigorous mathematical guidance, even if short of the 'proven-theorem' variety, can be helpful in advancing the science of deformation, flow, and fracture of solids.
A few words on notation: tensor components (when invoked) are written with respect to the basis of a fixed Rectangular Cartesian coordinate system.All spatial differential operators are w.r.t.position on the current configuration.A superposed dot represents a material time derivative.X will be the alternating tensor and the curl operator acting on tensor fields may simply be thought of as row-wise curls of the corresponding matrix field of components.

Governing equations: mechanics
Based on the detailed kinematic motivations presented in Acharya (2011Acharya ( , 2018Acharya ( , 2020)), the governing equations of the model are given by where T = T T ensures balance of angular momentum, and is the velocity gradient.In the above, T is the Cauchy stress, ρ is the mass density, v is the material velocity, f is the prescribed body force density, W is the inverse elastic distortion (a 2-point tensor field), V α is the dislocation velocity (vector field), L p is a meso-macro- scale construct not used in the fundamental microscale theory, the plastic distortion rate of dislocations (tensor field) that are 'averaged out' in terms of their charge (the meaning of this can be made precise in terms of microscopic quantities), c is the crack (vector field), and V t is the crack-tip velocity (vector field).The magnitude of the crack vector field encapsulates the degree of damage at a material point, while its orientation reflects that of the crack face normal at that point.An independent vector-valued field representing the crack-face normal as a fundamental kinematic ingredient in a PDE model of fracture was introduced in Acharya (2018); Steinke et al. (2019), and is beginning to be used (Morin and Acharya 2021;Hakimzadeh et al. 2022;Steinke et al. 2022).The dislocation and crack-tip velocity fields are relative velocities of the motion of the dislocation density field α and the crack-tip field t, respectively, w.r.t. the material velocity.Defining the dislocation and crack-tip line density fields Equations (1c) and (1d) imply (Evolution of the crack field) Physically, Eqs.(1c, 1d) are motivated from the conservation of topological charge statements Eqs.(3a, 3b) on assuming that a 'free' gradient that arises in the process in each case vanishes.
In the above, W and c are considered to be dimensionless physical quantities.

Guidance for constitutive assumptions from the second law of thermodynamics
We consider the free-energy density per unit mass and require that the power supplied by external agents be greater than or equal to the rate of change of the sum of the free energy and kinetic energy of the body: for any process in which the mechanical equations hold, where C is the (time-varying) current configuration of the body; and this is ensured by choosing constitutive assumptions for T , V α , V c , L p that guarantee Eq. ( 5).Now, so that Eq. ( 5) can be expressed as (3a) Since the response function for ψ is invariant under superposed rigid body motions, it can be shown (see Appendix) that the term highlighted in blue in Eq. ( 6) is symmetric.Since the stress is symmetric due to balance of angular momentum, this also implies that the spin ( L skw ) does not appear in the dissipation for the model making the latter invariant under superposed rigid motions.This is an important consistency check on the kinematic structure of the model.
Viewing L p , V α , and V t (in the bulk and at the boundary) as the sole dissipative mecha- nisms of the model, one recovers the stress-relation of the model from the consideration of energetically reversible, purely elastic processes: Finally, a sufficient condition for non-negative dissipation is obtained by choosing the constitutive assumptions for the dissipative mechanisms to be in the 'direction' of their respective 'driving forces, ' as exposed in Eq. ( 6), in the bulk and at the boundary.

Define (and we caution that F is not, in general, the gradient of a deformation w.r.t.a fixed global reference configuration).
With I the fourth-order identity tensor on the space of symmetric second order tensors, the intact elastic modulus given by C , the damaged elastic modulus by C , > 0, µ > 0 the intact Lamé parameters, and and µ the Lamé parameters for the damaged material, define where , µ are positive, monotone decreasing functions of |r| from the intact values of the parameters to some small (positive) residual values.Let ρ 0 > 0 be the mass density of the intact, unstretched elastic material, and H(x) = 0 for x ≤ 0 and H (x) = 1 for x > 0 be the Heaviside function (an appropriately smoothed rep- resentation will also suffice).We now define the strain energy density of the material accounting for damage due to cracking as (cf.Morin and Acharya (2021)) ψ e has physical dimensions of energy per unit mass, and A • 4 B := A ijkl B ijkl for fourth- order tensors A, B. (7) It can be shown (using arguments given in Acharya (2011); Morin and Acharya (2021)) that for a crack-only model with ψ = ψ e (W , c) = ψ e F T F , c frame-indifferent, the 'elastic-distortion driven' part of the Cauchy stress, say T (F ) , is given by T and the normal stress component to the crack, |c| −2 c • T (F ) c , at a damaged point where H(|r|) = 1 is given, up to a factor of ρ|r| 2 ρ 0 |c| 2 , by (tr(E) − E r ) + E r + 2µE r if the material point experiences compressive strain charac- terized by H (E r ) = 0 , and by tr(E) + 2 µE r if the point experiences tensile normal strain perpendicular to the crack.
We also introduce a crack-resistance energy density function η(|c|) with physical dimensions of energy per unit mass.A typical example reflecting no residual energy stored in damaged regions is where a ≥ 0 (with physical dimensions of stress) and c sat > 0 (dimensionless) are material constants.Another example, modeling Griffith type 'surface energy' (but not dependent on crack-length) is With these constructs a physically reasonable constitutive assumption for the free energy density of our material is where 0 < µ α , µ t = O(µ) are material constants with dimensions of stress, and l α , l t > 0 are material constants with dimensions of length.The lengths involved are expected to be much smaller than typical macroscopic dimensions.
Turning to the constitutive equation for the dislocation velocity in the bulk, motivated by the 'driving force' for dislocation motion in Eq. ( 6), define and a dislocation mobility tensor of the form where f = 0 or 1, and m gl , m cl ≥ 0 are material constants with physical dimensions of length 2 stress.timefor f = 0 and length stress.timefor f = 1 .Then we propose the constitutive assumption and note that when ψ is independent of α, ( which is the natural generalization of the form of the Peach-Koehler force of classical dislocation theory to finite deformation.
For the crack velocity we assume a simple isotropic mobility: where m cr > 0 (with same physical dimensions as m cl or m gl ) is a material constant reflecting crack mobility and f = 0 or 1. Ignoring the contribution of the second term in the crack driving force, We note that the crack velocity is not restricted to be in the direction of t × c , allowing crack-tip motions off of the local crack-plane (defined by c ⊥ ).Equations ( 7), ( 11), ( 13), and ( 14) form the constitutive assumptions of a specific model.

Motivating questions for the development of the model
Here, we outline some fundamental physical problems that served as the motivation for the development of the theory, and which can be used to evaluate its predictive capability through analysis and computation in the future.

'Stokes flow' from nonlinear elasticity with defects
Consider the quasi-static approximation for balance of linear momentum without body force: Since this holds for all times, it can be shown that this is equivalent to with div T = 0 initially.Assuming, for simplicity, that ψ = ψ(F , c) , we have that T = T (F , c) so that Ṫ = ∂ F T : Ḟ + ∂ c T • ċ , and combining with Eq. (1c) written in the form and Eq. ( 1d) we obtain Combining Eqs. ( 17) and ( 15) and defining we have Evidently, the evolution of cracks and dislocations play the role of a 'body-force forcing' in the evolution of the deformation of the body.
In the presence of dislocations and cracks in general, but in the absence of their motion relative to the material, we have with standard combinations of Dirichlet b.c. on the velocity on the boundary and Neumann conditions related to the First Piola-Kirchhoff traction (w.r.t the current configuration as the reference) rate.
Evidently, Eq. ( 19) does not reduce to an isotropic 4 th -order tensor acting on the stretch- ing tensor L sym = D , but this system is energetically and mechanically (in terms of applied loads) reversible, whereas 'viscous Stokes flow' is only mechanically reversible.

Relation between defects in elastic solids and viscous fluids, and the mechanical load induced solid-fluid transition
In an elastic solid, (dislocation) defects can be said to arise when the inverse elastic distortion is no longer curl-free, i.e.
In a fluid one might say that defects arise when the velocity gradient develops a 'singular part, ' thinking, roughly, that the velocity field is discontinuous across 2-d surfaces.
What might be the connection between these two ideas?Can such a connection, in the context of a specific constitutive model, be used to study the transition of a solid to a fluid due to a proliferation of defects?
Noting Eq. (1c) rewritten in the suggestive form of Eq. ( 16) and assuming L p = 0 (a coarse-scale 'homogenized' effect), when the α field is a distribution of superposed core fields moving with velocity V α , α × V α very much looks like a singular distribution (when viewed macroscopically), and then Eq. ( 16) suggests that ḞF −1 -the elastic part of the velocity gradient (of the solid (fluid?)) -is its 'regular part' (the absolutely continuous part), with F (α × V α ) being its 'singular' part.

Dislocation nucleation
Here we consider a model with no cracks and L p = 0 .Assume ψ = ψ(F ); T = T (F ).

Quasi-static balance of forces
The governing equations for v, W , α are: (although the fields v, W suffice, nucleation related questions are best dealt with the α equation).
• Question: Do perturbations in α from a dislocation-free state α = 0 grow?Charac- terize the instability in terms of the class of elastic distortion fields F and energy densities ψ(F ) .The constitutive choices for V α can be as in Eq. ( 13) and further simpli- fied as necessary, e.g.assume isotropic mobility.The initial state satisfies div T = 0 and loading is required.

Dynamic balance of forces
In Eq. (20a) replace Eq. ( 20a) with balance of linear momentum and balance of mass and ask the same question as in "Quasi-static balance of forces" section.

Crack nucleation
Here we consider a model with no dislocation or plasticity, α, L p = 0 , and ψ = ψ(F , c); T = T (F , c) .Here, W is a gradient on the current configuration and F is as well, on the reference defined by the inverse deformation which is a potential for W since curl W = 0.

Quasi-static balance of forces
The governing equations for v, W , c, t are: In the above Ẇ + WL = 0 ⇒ ḞF −1 = L = ∇v .As in the dislocation case, one of Eqs. ( 20e) and (20f ) suffices, but can be used as necessary.
• Question: Do perturbations in t from a crack-free state c = 0 =⇒ t = 0 grow?Char- acterize the instability in terms of the class of elastic distortion fields F and energy densities ψ(F ) .The constitutive choices for V t can be as in Eq. ( 14).The initial state satisfies div T = 0 and and loads are required.

Dynamic balance of forces
In Eq. (20d) replace Eq. ( 20d) with balance of linear momentum and balance of mass and ask the same question as in "Quasi-static balance of forces" section. (20c Remark 1 For dislocation or crack nucleation from an undislocated or uncracked state, respectively, the corresponding evolution equations for the perturbation in dislocation density ( α ) and crack-tip density ( t ) are given by where, for simplicity, we assume that i = 1, 2 and only straight dislocation/crack-tips in the 3-direction are allowed.Based on the above, it seems that the distinction between crack and dislocation nucleation in this ansatz is a matter of nonlinear stability.We note that for the purposes of linear stability, B = 0 at the crack-free state.

Brittle-ductile transition
This is a coupled crack-dislocation problem.The initial condition is that of an unloaded body with an edge crack as shown in Fig. 1.
• Question: Under, say Mode I, loading, i.e., Dirichlet conditions on velocity v 2 = 0 on top boundary with bottom fixed as shown, does the stress field of the crack with a concentration at the notch-tip nucleate a dislocation (or a dipole) in the body which then moves (expands) causing plasticity (ductile behavior), or does the crack propagate without any dislocation nucleation and propagation (brittle behavior)?Characterize based on material parameters of the model (for a large body).

Macroscopic model of elasto-viscoplasticity
Consider the dislocation-only model and define ε := l α H , recall Eq. ( 11), where H is a representative dimension of the body and we will be interested in ε → 0 with l α fixed.
Fig. 1 Schematic of a body with an edge crack under load.In the brittle-ductile transition, the question is whether under load the crack propagates or a dislocation (dipole) nucleates and moves Consider the system subjected to a constraint on the initial condition a boundary condition where v is a given function, and • Question: (assuming existence of solutions for ε > 0 , or plausible demonstration of such in finite-dimensional computational settings), what is the limit model that arises as ε → 0 ?Here, the 'limit model' is the question of what model the evolution of the weak limits (roughly, all smoothly weighted space-time averages) of the fundamental fields of the model, say α ε , v ε here, satisfy.Since the microscopic equations are nonlinear and averages of nonlinear functions of field quantities ( say, e.g., α ε × V ε ), do not equal the same functions evaluated on the averages of the said quantities, this requires the determination of the evolution of the weak limits of a further set of quantities, like say |α ε | 2 , defined on a sequence of solutions of the microscopic model.Moreover, to be useful, such evolution of the limits must be 'closed' in the sense that it must need information only on the state of only the limits of these quantities at any given time.In particular, there is good physical intuition behind the expectation that lim ε→0 α ε × (V α ) ε pro- duces an extra term in the limit, related to the plastic strain rate produced by the expansion of 'sub-grid' loops, the latter not sensed by lim ε→0 α ε .In fact, this is the reason for the phenomenological introduction of the term L p (and only this term as representative of the plastic strain rate) in macroscopic models.Is the limit parametrized by constant α ?

Macroscopic model of damage
Consider the crack-only model and define ε := l t H , recall Eq. ( 11), where H is a representative dimension of the body and we will be interested in ε → 0 with l t fixed.
Consider the system In the above, the equivalence is not strict since the bottom equation implies the top one up to the gradient of a scalar field which is assumed to vanish based on the assumption that microscopically the crack-tip flux can occur only in the presence of a crack-tip at a point (much like microscopic plastic strain rate/slipping rate at a point can arise only if a dislocation is present at a point).
Let the system be subjected to a constraint on the initial condition a boundary condition where v is a given function, and • Question: As in "Macroscopic model of elasto-viscoplasticity" section, what is the limit model as ε → 0 ?Does a natural connection arise with the type of coupled brit- tle-ductile model of fracture proposed in Acharya (2020)?
Classical elasto-viscoplasticity, viscoplasticity (a non-Newtonian viscous fluid), as limit models Consider the dislocation-only model from Eq. (1a).The classical, phenomenological model of elasto-viscoplasticity is given by the system Eq.(1a) with V α = 0 and L p and T(F) specified by constitutive assumptions.The strain-rate decomposition Eq. ( 16), that follows from the conservation of Burgers vector Eq.(3a), then takes the form Recall that the model does not involve a reference configuration of any sort and F is not a deformation gradient in general (customarily it is written as F e , but a 'multiplica- tive decomposition' of a deformation gradient from any reference plays no role in our development).Define The constitutive equation for L p specifies D p and ω p ; in describing the elastoviscoplas- ticity of polycrystals without texture, it is customary to assume (but not for single crystals or strongly textured polycrystals).