On a discrete-to-continuum convergence result for a two dimensional brittle material in the small displacement regime

We consider a two-dimensional atomic mass spring system and show that in the small displacement regime the corresponding discrete energies can be related to a continuum Griffith energy functional in the sense of Gamma-convergence. We also analyze the continuum problem for a rectangular bar under tensile boundary conditions and find that depending on the boundary loading the minimizers are either homogeneous elastic deformations or configurations that are completely cracked generically along a crystallographic line. As applications we discuss cleavage properties of strained crystals and an effective continuum fracture energy for magnets.


Introduction
A fundamental problem in static fracture mechanics is to determine the behavior of a brittle material which is subject to certain displacements imposed at its boundary. Of particular interest is the identification of critical loads at which failure occurs. A natural framework to treat such free discontinuity problems with variational methods is given by Griffith energy functionals introduced by Francfort and Marigo [21] comprising elastic bulk contributions and surface terms comparable to the size of the crack (see also [17]). Often these models contain anisotropic surface terms (see e.g. [2,19,26]) modeling the fact that due to the crystalline structure of the materials certain directions for the formation of cracks are energetically favored. Indeed, fracture typically occurs in the form of cleavage along crystallographic planes. Ultimately, such a continuum model should be identified as an effective theory derived from atomistic interactions.
Specifying the set-up even further, a basic experiment to infer material properties of brittle materials is to probe the specimen by applying a uniaxial tensile strain which allows to determine its Poisson ratio in the elastic regime and a critical load beyond which the body fails due to fracture. From a theoretical point of view this problem has been studied recently by Mora-Corral in [25], where he investigates a rectangular bar of brittle, incompressible, homogeneous and isotropic material subject to uniaxial extension and shows that, depending on the loading, the minimizers are either given by purely elastic configurations or deformations with horizontal fracture.
An atomistic model problem with surface contributions sensitive to the crack geometry has been studied by the authors in [22] leading to a complete analysis of the asymptotically optimal configurations under uniaxial extension in the discrete-to-continuum limit: The body shows pure elastic behavior in the subcritical case and for supercritical boundary values generically cleavage occurs along a specific crystallographic line. However, for a certain symmetric orientation of the lattice cleavage may fail more complicated crack geometries are possible.
The goal of this work is to show that in the small displacement regime the energies associated to such a discrete system can be related to a continuum Griffith energy functional with anisotropic surface contributions in the sense of Γconvergence. Moreover, we analyze the continuum problem under tensile boundary conditions. In this way we (1) obtain convergence scheme which in certain applications to be discussed below allows to identify effective continuum fracture energies, (2) extend the results of [25] to anisotropic and compressible materials and (3) re-derive in part the aforementioned convergence results of [22].
In the theory of fracture mechanics the passage from discrete systems to continuum models via Γ-convergence is by now well understood for one-dimensional chains, see e.g. [7,8,9]. In the higher dimensional setting there are results for scalar valued models (see [10]) and approximations of vector valued free discontinuity problems where the elastic bulk part of the energy is characterized by linearized terms (see [2]) or by a quasiconvex stored energy density (see [19]). However, in more than one dimension the energy density of discrete systems such as well-known mass spring models is in general not given in terms of a discretized continuum quasiconvex function. For large strains these lattices typically become even unstable, see e.g. the basic model discussed in [24]. Consequently, in the regime of finite elasticity it is a subtle question if minimizers for given boundary data exist at all. On the other hand, for sufficiently small strains one may expect the Cauchy-Born rule to apply so that individual atoms do in fact follow a macroscopic deformation gradient, see [24,14]. In particular this applies to the regime of infinitesimal elastic strains. For purely elastic interactions this relation has also been obtained in the sense of Γ-convergence for a simultaneous passage from discrete to continuum and linearization process in [12,27].
The model considered in [22] as well as the one-dimesional seminal paper [11] suggest that the most interesting regime for the elastic strains is given by √ ε (ε denotes the typical interatomic distance) as in this particular regime the elastic and the crack energy are of the same order. This is in accordance to the observation that brittle materials develop cracks already at moderately large strains. Moreover, it shows that a discrete-to-continuum Γ-limit for the discrete energies under consideration naturally involves a linearization process.
Identifying all possible limiting continuum configurations and energies is a challenging task as necessary smallness assumptions on the discrete gradient can not be inferred from suitable energy bounds and deriving rigidity estimates being essential in the passage from nonlinear to linearized theory (see [12,27]) is a subtle problem. Partial results have been obtain in [22] for almost minimizers of a boundary value problem describing uniaxial extension. A general analysis in two dimensions is deferred to a subsequent work. In the present context we make the simplifying assumption that we consider deformations lying √ ε-close to the identity mapping. However, we will also see that there are physically interesting applications e.g. to magnetic materials where such an assumption can be justified rigorously.
It then turns out that the derivation of the continuum limit is an issue similar to those considered in [2,10,19]. Nevertheless, we believe that the present Γ-convergence result is interesting as (1) it gives rise to a limiting Griffith functional in the realm of linearized elasticity which can be explicitly investigated for cleavage, (2) there are applications to systems with small displacements for small energies and (3) to the best of our knowledge our approach to the problem differs from techniques which are predominantly used when treating discrete systems in the framework of fracture mechanics.
The reduction to one-dimensional sections using slicing properties for (special) functions of bounded variation turned out to be a useful tool not only to derive general properties of these function spaces but also to study discrete systems and variational approximation of free discontinuity problems. E.g., the original proofs of the main compactness and closure theorems in SBV (see [3]) as well as the Γ-convergence results in [10,19] make use of this integralgeometric approach. Similar to the fact that there are simplified proofs of these compactness theorems being derived without the slicing technique (see [1]), we show that in our framework the lower bound of the Γ-limit can be achieved in a different way. In fact, we carefully construct the crack shapes of discrete configurations in an explicit way which allows us to directly appeal to lower semicontinuity results for SBV functions.
The paper is organized as follows. We first introduce our discrete model and state our main results in Section 2. Here we also briefly discuss how these results shed new light on our findings in [22] on crystal cleavage and study an application to fractured magnets in an external field. Section 3 is devoted to the derivation of the continuum energy functional via Γ-convergence. The main idea for the lower bound relies on a separation of the energy into elastic and surface contributions by introducing an interpolation with discontinuities on triangles where large expansion occurs. By constructing the set of discontinuity points in a suitable way the surface energy can be estimated using lower semicontinuity results for SBV functions. The elastic part can be treated similarly as in [23,27].
Finally, in Section 4 we analyze the continuum problem under tensile boundary values and extend the results obtained in [25] to anisotropic and compressible materials. A careful analysis of the anisotropic surface contribution shows that in the generic case there is a unique optimal direction for the formation of fracture, while in a symmetrically degenerate case cleavage fails and all energetically optimal crack geometries can be characterized by specific Lipschitz curves. As in [25] the proof makes use of a qualitative rigidity result for SBV functions (see [13]) and of the structure theorem on the boundary of sets of finite perimeter by Federer [18].
, respectively. Without loss of generality we can assume φ ∈ [0, π 3 ). We collect the basic lattice vectors in the The macroscopic region Ω ⊂ R 2 occupied by the body is supposed to be a bounded domain with Lipschitz boundary. In the reference configuration the positions of the specimen's atoms are given by the points of the scaled lattice εL that lie within Ω. Here ε is a small parameter defining the length scale of the typical interatomic distances.
The deformations of our system are mappings y : εL ∩ Ω → R 2 . The energy associated to such a deformation y is assumed to be given by nearest neighbor interactions as Note that the scaling factor 1 ε in the argument of W takes account of the scaling of the interatomic distances with ε. The pair interaction potential W : [0, ∞) → [0, ∞] is supposed to be of 'Lennard-Jones-type': (i) W ≥ 0 and W (r) = 0 if and only if r = 1.
In order to analyze the passage to the limit as ε → 0 it will be useful to interpolate and rewrite the energy as an integral functional. Let C ε be the set of equilateral triangles △ ⊂ Ω of sidelength ε with vertices in εL and define Ω ε = △∈Cε △. Byỹ : Ω ε → R 2 we denote the interpolation of y, which is affine on each △ ∈ C ε . The derivative ofỹ is denoted by ∇ỹ, whereas we write (y) △ for the (constant) value of the derivative on a triangle △ ∈ C ε . Then (1) can be rewritten as where Here we used that |△| = √ 3ε 2 /4. The boundary term is the sum of pair interaction energies 1 4 ) over nearest neighbor pairs which form the side of only one or no triangle in C ε , respectively.
Due to the discreteness of the underlying atomic lattice, Dirichlet boundary conditions have to be imposed in a small neighborhood of the boundary as otherwise cracks near the boundary may become energetically more favorable. Assume thatΩ ⊃ Ω is a bounded, open domain in R 2 with Lipschitz boundary defining the Dirichlet boundary ∂ D Ω = ∂Ω ∩Ω of Ω. For g ∈ W 1,∞ (Ω) we define the class of discrete displacements assuming the boundary value g on ∂ D Ω as where Ω D,ε := {x ∈Ω : dist(x,Ω \ Ω) ≤ ε}. For the corresponding deformations y = id + u this amounts to requiring y(x) = x + g(x) for x ∈ εL ∩ Ω D,ε . Similar as before, we letC ε be the set of equilateral triangles △ ⊂Ω with vertices in εL and defineΩ ε = △∈Cε △. Byỹ :Ω ε → R 2 we again denote the piecewise affine interpolation of y.
It is easy to see that the formation of a crack of finite length resulting from a number of largely deformed triangles scaling with 1 ε leads to an energy contribution to E ε scaling with ε. The most interesting regime is when the elastic energy contributions to E ε and the energy cost of a cracked configurations are of the same order. We are thus particularly interested in boundary displacements g ε scaling with √ ε. For then there are also completely elastic deformations for which E ε scales with ε, e.g. E ε (id + g ε ) = O(ε).
Moreover, we will assume u ∈ A gε for some g ε ∈ W 1,∞ (Ω). We also introduce the functionals E χ ε which arise from E ε by replacing W ∆ by W ∆,χ = W ∆ + χ, where χ : R 2×2 → [0, ∞] is a frame indifferent penalty term with χ ≥ c χ > 0 in a neighborhood of O(2) \ SO(2) and χ ≡ 0 in a neighborhood of SO(2) ∪ {∞}. This term is a mild extra assumption to assure that the orientation of the triangles is preserved in the elastic regime and unphysical effects are avoided.

Convergence of the variational problems
Our convergence analysis applies to discrete deformations which may elongate a number scaling with 1 ε of springs very largely, leading to cracks of finite length in the continuum limit. On triangles not adjacent to such essentially broken springs, the defomations are √ ε-close to the identity mapping, so that the accordingly rescaled displacements are of bounded L 2 -norm. Note that the first of these assumptions can be inferred from suitable energy bounds. By way of example, however, we will see that this cannot be true for the displacement estimates in the bulk: The sequence of functionals (E ε ) ε is not equicoercive. Nevertheless, it is interesting to investigate this regime in order to identify a corresponding continuum functional which describes the system in the realm of Griffith models with linearized elasticity. In fact, below we will discuss two specific models where external fields or boundary conditions break the rotational symmetry whence the sequence (E χ ε ) ε satisfies suitable equicoercivity conditions. Recall that the space SBV (Ω; R 2 ), abbreviated as SBV (Ω) hereafter, of special functions of bounded variation consists of functions u ∈ L 1 (Ω; R 2 ) whose distributional derivative Du is a finite Radon measure, which splits into an absolutely continuous part with density ∇u with respect to Lebesgue measure and a singular part D j u whose Cantor part vanishes and thus is of the form where H 1 denotes the one-dimensional Hausdorff measure, J u (the 'crack path') is an H 1 -rectifiable set in Ω, ν u is a normal of J u and [u] = u + − u − (the 'crack opening') with u ± being the one-sided limits of u at J u . If in addition ∇u ∈ L 2 (Ω) and H 1 (J u ) < ∞, we write u ∈ SBV 2 (Ω). See [5] for the basic properties of these function spaces.
The sense in which discrete displacements are considered convergent to a limiting displacement in SBV is made precise in the following definition.
and there exists a sequence C * ε ⊂C ε with #C * ε ≤ C ε for a constant C independent of ε such that The main idea will be to separate the energy into elastic and crack surface contributions by introducing a threshold such that triangles △ with (y) △ beyond that threshold are considered as cracked andỹ is modified there to a discontinuous function. The treatment of the elastic part draws ideas from [27] and [23]. To derive the crack energy, one could use a slicing technique, see, e.g., [10]. Although also possible in our framework, we follow a different approach here: We carefully construct crack shapes of discrete configurations in an explicit way which allows us to directly appeal to lower semicontinuity results for SBV functions in order to derive the main energy estimates.
Consider the limiting functional , where e(u) = 1 2 ∇u T + ∇u denotes the symmetric part of the gradient. Q is the linearization of W △ about the identity matrix Id (see Lemma 3.2 for its explicit form). Note that for a displacement field u, which is the limit of a sequence (u ε ) ⊂ A gε converging in the sense of Definition 2.1, we get u = g onΩ \ Ω, where g = L 1 -lim ε→0 g ε . Therefore, if u| Ω does not attain the boundary condition g on the Dirichlet boundary ∂ D Ω (in the sense of traces), this will be penalized in the energy E(u). In Section 3 we prove the following Γ-convergence result (see [16] for an exhaustive treatment of Γ-convergence): Note that the recovery sequence is obtained for the energy E χ ε which includes the frame indifferent penalty term. Due to the frame indifference of W , (E ε ) and (E χ ε ) are not equicoercive as the following example shows.
Example. Assume that the specimen satisfying the boundary conditions is broken into three parts by two even cracks where the middle part is subject to a rotation R = Id so that In particular, the energy of the configuration is of order 1. But for p ≤ x 1 ≤ q Thus, ∇ũ ε is not bounded in L 1 and so u ε does not converge. We now add a term to E ε such that the sequence becomes equicoercive. Let m : R 2×2 → S 1 be a function satisfyinĝ Moreover, assume thatm is C 2 in a neighborhood of SO(2) and R 2×2 sym ⊂ ker(Dm(Id)).

Theorem 2.3
The assertions of Theorem 2.2 remain true when E ε , E χ ε and E are replaced by F ε , F χ ε and F , respectively.

Analysis of a limiting cleavage problem
We now analyze the limiting functional E for a rectangular slab Ω = (0, l) × (0, 1) with l ≥ 1 √ 3 under uniaxial extension in e 1 direction. We determine the minimizers and prove uniqueness up to translation of the specimen and the crack line for the boundary conditions (More precisely: u ∈ SBV 2 ((−η, l + η) × (0, 1)) with u 1 (x) = 0 for x ≤ 0 and u 1 (x) = al for x ≥ l.) Note that we can investigate the limiting problem without any assumption on the second component of the boundary displacement. Let and v γ is unique if and only if φ = 0. It turns out that the specimen shows perfect elastic behavior up to the critical boundary displacement Beyond critical loading the body fails by breaking into two pieces. Then All minimizers of E subject to (6) are of the following form: (i) If a < a crit , then x for some s ∈ R.
(iii) If a > a crit and φ = 0 then a.e. and constants s, t ∈ R. This theorem will be addressed in Section 3. An analogous result for isotropic, incompressible materials has been obtained recently by Mora-Corral [25]. Theorem 2.4 is an extension of this result to anisotropic, compressible brittle materials in the framework of linearized elasticity.
In particular, as mentioned above we see that all the optimal configurations show purely elastic behavior in the subcritical case and complete fracture in the supercritical regime. The crack minimizer in (ii) for φ = 0 is broken parallel to Rv γ which proves that cleavage occurs along crystallographic lines, while in the symmetric case φ = 0 cleavage in general fails.

Applications: Cleaved crystals and fractured magnets
As applications of the converging results for the energy functionals E χ ε and F χ ε we consider cleaved crystals and fractured magnets, respectively. In the first model a mild equicoercivity of the sequence (E χ ε ) ε is guaranteed by investigating a specific boundary value problem, in the latter model an external field provides an even stronger equicoercivity condition.

Uniaxially strained crystals
where G(a) := {g ∈ W 1,∞ (Ω) : g 1 (x) = 0 for x 1 ≤ ε, g 1 (x) = al for x 1 ≥ l − ε}. In [22, Theorem 2.1] we proved that the limiting minimal energy leads to a universal cleavage law of the form independent of the particular shape of the interatomic potential W . Optimal configurations are given by the constant sequences u ε = u el in the subcritical case a ≤ a crit and u ε = u cr in the supercritical case a ≥ a crit , respectively, with u el and u cr as in Theorem 2.4. In fact, the above given configurations provide a characterization of all minimizing sequences in the sense that, all low energy sequences (u ε ) ε satisfying and sup ε u ε ∞ < ∞ converge-up to subsequences-in the sense of Definition 2.1 to u el if a < a crit or u cr if a > a crit for suitable s, t, p and g, respectively. This is a direct consequence of [ ] is that, under the tensile boundary conditions u ε ∈ A(a), the requirement that u ε be an almost energy minimizer satisfying (7), guarantees the existence of a subsequence converging in the sense of Definition 2.1. In particular, the sequence (E χ ε ) is mildly equicoercive. A fundamental theorem of Γ-convergence (see, e.g., [6,Theorem 1.21]) implies that such low energy sequences converge to limiting configurations u el , respectively, u cr , in the sense of Definition 2.1. Consequently, in this way we have re-derived the convergence result [22,Corollary 2.4] (in the sense of Definition 2.1).

Permanent magnets in an external field
Assume that the material is a permanent magnet and let e 1 be the magnetization direction. We suppose that there is a constitutive relation between ∇ỹ(x) and the local magnetization directionm(ỹ, x) ∈ S 1 of the deformed configurationỹ at some point x ∈ Ω, which is of the formm(ỹ, x) =m(∇ỹ(x)) withm as defined in Section 2.2. Let H ext : R 2 → R 2 be an external magnetic field. The magnetic energy corresponding to the deformation y = id + √ εu is then given by i.e. alignment of the magnetization direction with the external field is energetically favored. The total energy of the system is given by We now suppose that the external field is homogeneous and satisfies without restriction H ext = κe 1 for κ > 0. We then see that with f κ as in (5)  We consider a boundary value problem min u∈Ag E tot (u) for g ∈ W 1,∞ (Ω). Since the sequence (E tot ε ) ε is equicoercive as discussed in Section 2.2, the theory of Γ-convergence implies lim ε→0 min u∈Ag F ε (u) = min u∈Ag F tot (u) and also convergence of the corresponding (almost) minimizers of E tot ε in the sense of Definition 2.1 is guaranteed. In this context, note that by a truncation argument taking g ∈ W 1,∞ (Ω) into account, we may indeed assume that a low energy sequence satisfies sup ε u ε ε < +∞.

Preparations
The goal of this section is the derivation of the Γ-convergence result for E ε . We first collect some properties of the cell energy W △ proven in [22,Section 3] provided that W satisfies the assumptions (i), (ii) and (iii).
In particular, Q(G) only depends on the symmetric part G T + G /2 of G. Q is positive semidefinite and thus convex on R 2×2 and positive definite and strictly convex on the subspace R 2×2 sym of symmetric matrices. The following lemma provides useful lower bounds for the energy W △ and the pair interaction potential W . (ii) We defineψ (r) = ηr for 0 ≤ r ≤ β η , β for r ≥ β η , for some η > 0 (depending on ρ) such thatψ − ρ ≤ W . Then we set ψ ρ (r) = ψ(r + 1). As ψ ρ is a concave function with ψ ρ (0) > 0, it is subadditive.
Lemma 3.4 Let T > √ 2. Then there are constants C 1 , C 2 > 0 such that for all F ∈ R 2×2 with |F | ≤ T we obtain Proof. (i) Without restriction we may assume that |F − R(F )| is small as otherwise the assertion is clear. So in particular, R(F ) is uniquely determined. Moreover, it suffices to consider F ∈ R 2×2 sym and R(F ) = Id. Indeed, once this is proved, we find |m sym ⊂ ker(Dm(Id)). (2)) for c 0 > 0 and C = C(c 0 , T ) sufficiently small. Otherwise, we again apply Lemma 3.3(i) to obtain for c 0 small enough
As a further preparation we modify the interpolationỹ on triangles with large deformation: We fix a threshold explicitly as R = 7 and letC ε ⊂C ε be the set of those triangles where |(ỹ) △ | > R. By definition of the boundary values in (4) we findC ε ⊂ C ε for ε small enough. We introduce another interpolation y ′ which leavesỹ unchanged on △ ∈C ε \C ε and replacesỹ on △ ∈C ε by a discontinuous function with constant derivative satisfying |(y ′ ) △ | ≤ R. In fact, by introducing jumps we achieve a release of the elastic energy. Note that y ′ ∈ SBV (Ω ε ).
More precisely, note that on △ ∈C ε we have |(ỹ) △ v| ≥ 2 for at least two springs v ∈ V. Indeed, using the elementary identity for any H ∈ R 2×2 sym , we find that |F | > 7 implies and so max v∈V |F v| 4 > 7 4 8 > 4 4 . Hence, |F v| > 4 for at least one v ∈ V and at least two springs are elongated by a factor larger than 2. For m = 2, 3 let C ε,m ⊂C ε be the set of triangles where |(ỹ) △ v| ≥ 2 holds for exactly m springs v ∈ V. For i, j, k = 1, 2, 3 pairwise distinct let h i denote the segment between the centers of the sides in v j and v k direction and define the set We now construct y ′ ∈ SBV 2 (Ω ε ). On △ ∈C ε \C ε we simply set y ′ =ỹ. On Moreover, we ask that y ′ =ỹ at the three vertices and on the side orientated in v i direction. This can and will be done in such a way that y ′ is continuous on int(△) \ h i . We note that the definition of (y ′ ) △ is unique up to a reflection, unless (ỹ) △ v i = 0. We may and will assume that For △ ∈C ε,3 we set (y ′ ) △ = Id and y ′ =ỹ at the three vertices such that y ′ is continuous on int(△) \ V i for some i ∈ {1, 2, 3}. Here, the index i can be taken arbitrarily at first. However, in what follows it will also be necessary to use the following unambiguously defined 'variants' of y ′ : If on every △ ∈C ε,3 the set V i is chosen as the jump set of y ′ we denote this interplation explicitly as y ′ V i . We define the interpolation u ′ for the rescaled displacement field by u ′ = 1 √ ε (y ′ − id). We note that by construction also on an edge [p, q] ⊂ ∂△ for △ ∈C ε jumps may occur. There, however, the jump height |[u ε ]| can be bounded by for a constant c > 0 independent of ε and x ∈ [p, q]. This holds since the interpolations are continuous at the vertices.
The following lemma shows that we may pass fromũ ε to u ′ ε without changing the limit.
Proof. We first note that there is some M > 0 such that for all ε > 0. To see this, we just recall that every triangle △ ∈C ε provides at least the energy ε inf {W (r) : r ≥ 2}. In fact we may assume that C * ε =C ε in Definition 2.1 as for ∆ It follows that χΩ ε ∇u ′ ε is bounded uniformly in L 2 and, in particular, equiintegrable. Finally, the jump lengths H 1 (J u ′ ε ) are readlily seen to be bounded by Cε#C ε ≤ C. But then Ambrosio's compactness Theorem for GSBV [4,Theorem 2.2] shows that indeed χΩ ε ∇u ′ ε ⇀ ∇u in L 2 (Ω).

The Γ-lim inf-inequality
With the above preparations at hand, we may now prove the Γ-lim inf-inequality in Theorem 2.2.
Proof of Theorem 2.2(i). Let (g ε ) ε ∈ W 1,∞ (Ω) with sup ε g ε W 1,∞ (Ω) < +∞ be given. Let u ∈ SBV 2 (Ω) and consider a sequence u ε ⊂ SBV 2 (Ω ε ) with u ε ∈ A gε converging to u in SBV 2 in the sense of Definition 2.1. We split up the energy into bulk and crack parts neglecting the contribution εE boundary ε from the boundary layers: We note that by contruction of the interpolation u ′ ε we may take the integral over Ω ε . As both parts separate completely in the limit, we discuss them individually.
Elastic energy. We first concern ourselves with the elastic part of the energy. We Note that for F ∈ R 2×2 , r > 0 one has Q(rF ) = r 2 Q(F ). We compute The second term of the integral can be bounded by Since ∇u ′ ε is bounded in L 2 and χ ε By assumption χ Ωε ∇u ′ ε ⇀ ∇u weakly in L 2 . As χ ε → 1 boundedly in measure on Ω, it follows χ Ωε χ ε ∇u ′ ε ⇀ u weakly in L 2 (Ω). By lower semicontinuity (Q is convex by Lemma 3.2) we conclude recalling that Q only depends on the symmetric part of the gradient: Crack energy. By construction the functions u ′ ε have jumps on destroyed triangles △ ∈C ε . We now write the energy of such a triangle in terms of the jump height [u] = u + − u − . We first concern ourselves with a triangle △ ∈C ε,3 . For the variant u ′ ε,V i , i = 1, 2, 3 we consider the springs in v j , v k direction for j, k = i. Thus, we compute where [u ′ ε,V i ] hv k denotes the jump height on the set h v k . Here and in the following equations, the same holds true if we interchange the roles of j and k. We claim that Otherwise, applying (13) we compute for ε small enough: Let ρ > 0 sufficiently small. Applying Lemma 3.3(ii) there is an increasing subadditive function ψ ρ 0 with ψ ρ (r − 1) − ρ ≤ W (r) for r ≥ 1. We definẽ ψ ρ = ψ ρ − ρ. The monotonicity of ψ ρ and (14) yield Now for △ ∈C ε,3 we may estimate the energy as follows: where i, j, k = 1, 2, 3 are pairwise distinct. With ν The factors in front occur since H 1 (h v j ) = ε 2 and, letting ν j be a normal of h v j , one has |ν j · v j | = 0 and |ν j · v k | = √ 3 2 . Consequently, defining φ ρ i (r, ν) = ψ ρ (r) (|v j · ν| + |v k · ν|) andφ ρ i (r, ν) =ψ ρ (r) (|v j · ν| + |v k · ν|), respectively, we get on every △ ∈C ε,3 . For △ ∈C ε,2 we proceed analogously. Assuming |(ỹ ε ) △ v i | ≤ 2 we compute for the springs in v j , v k direction (abbreviated by v j,k ) as in (13) Note that in this case we do not have to take a special variant of u ′ ε into account. Repeating the steps (14) and (15) we find Noting that |v j · ν i | = |v k · ν i | = √ 3 2 , |v i · ν i | = 0 and that every of these terms occurs twice in the sum of the right hand side of the following formula, it is not hard to see that this energy satisfies the same integral representation formula as W △,3 : (Recall that the interpolation variant u ′ ε,V i and its crack normal ν (i) u do not depend on i on △ ∈C ε,2 .) Let σ > 0. Note thatC ε ⊂ C ε for ε sufficiently small as sup ε g ε W 1,∞ (Ω) < +∞. Thus, the crack energy can be estimated by where E ρ ε,∪∂△ (ỹ ε ) compensates for the extra contribution provided by jumps lying on the boundary of some △ ∈C ε . We will show that this term vanishes in the limit. Now by construction the φ ρ i (r, ν), i = 1, 2, 3, are products of a positive, increasing and concave function in r and a norm in ν. Moreover, u ′ ε and its variants converge to u in L 1 with ∇u ′ ε bounded in L 2 and thus equiintegrable. By Ambrosio's lower semicontinuity Theorem [4, Theorem 3.7] we obtain where we used that sup ε H 1 (J u ′ ε ) ≤ CM for a constant C > 0 by (11). We recall that ψ ρ (r) → β for r → ∞. In the limit σ → 0 this yields Taking (10) and (11) This proves lim sup ε |E ρ ε,∪∂△ (ỹ ε ) | ≤CMρ for someC > 0. We finally let ρ → 0 in (17). This finishes the proof of (i).

Recovery sequences
It remains to construct recovery sequences in order to complete the proof of Theorem 2.2.
Proof of Theorem 2.2(ii). The basic tool for the proof of the Γ-limsup-inequality is a density result for SBV functions due to Cortesani and Toader [15]. Moreover, a proof very similar to that of Proposition 2.5 in [20] shows that we may also impose suitable boundary conditions on the approximating sequence. We suppose W(Ω, R 2 ) is the space of all SBV functions u ∈ SBV (Ω, R 2 ) such that J u is a finite union of (disjoint) segments and u ∈ W k,∞ (Ω \ J u , R 2 ) for all k. Then W(Ω, R 2 ) is dense in SBV 2 (Ω, R 2 ) ∩ L ∞ (Ω, R 2 ) in the following way: For every u ∈ SBV 2 (Ω, R 2 ) ∩ L ∞ (Ω, R 2 ) with u = g onΩ \ Ω, there exists a sequence u n and a sequence of neighborhoods U n ⊂Ω ofΩ \ Ω such that u n = g on Ω D, 1 n (recall (4)), u n ∈ W 1,∞ (U n ) and u n | Vn ∈ W(V n , R 2 ), where V n ⊂ Ω is some neighborhood of Ω \ U n , such that u n ∞ ≤ u ∞ and (i) u n → u strongly in L 1 (Ω, R 2 ), ∇u n → ∇u strongly in L 2 (Ω, R 2 ), (ii) lim sup n→∞ Ju n φ(ν un )dH 1 ≤ Ju φ(ν u )dH 1 for every upper semicontinuous function φ : Let u ∈ SBV 2 (Ω, R 2 ) with u = g onΩ\Ω. Without restriction we can assume u ∈ ∩L ∞ (Ω, R 2 ) as this hypothesis my be dropped by applying a truncation argument and taking Q(F ) ≤ C|F | 2 into account. In fact, it suffices to provide a recovery sequence for an approximation u n defined above. Although our notion of convergence in Definition 2.1 is not given in terms of a specific metric, similarly to a general density result in the theory of Γ-convergence this can be seen by a diagonal sequence argument. The crucial point is that due to (20) below we may assume that for ε sufficiently small (depending on n) where C is independent of n and ε. If (u n,ε ) ε is a recovery for u n , one may therefore pass to a diagonal sequence which is a recovery sequence for u, in particular converging to u the sense of Definition 2.1. For simplicity write u instead of u n in what follows. Let δ > 0 and define J δ | is Lipschitz continuous on J u , it cannot oscillate infinitely often between values ≤ δ and values ≥ 2δ on a single segment. Consequently, there is a finite number N δ u of disjoint subsegments S 1 , . . . , S N δ u in J u such that |[u]| < 2δ on every S j and |[u]| > δ on J u \(S 1 ∪. . .∪S N δ u ). Note that H 1 ( N δ u i=1 S i ) ≤ H 1 (J u \J 2δ u ) =: ρ(δ) → 0 for δ → 0. We cover S 1 , . . . , S N δ u by pairwise disjoint rectangles Q 1 , . . . Q N δ u which satisfy j H 1 (∂Q i ) + |Q i | ≤ Cρ(δ). It is not hard to see that |u(x) − u(y)| ≤ CH 1 (∂Q i ) + 2δ for x, y ∈ Q j as ∇u ∈ L ∞ (Ω).
We modify u on the rectangles Q i : Let u δ = u onΩ \ N δ u i=1 Q j and define u δ = c j on Q j for c j ∈ R 2 in such a way that J u δ = J δ u δ up to an H 1 -negligible set. As u ∈ L ∞ (Ω), ∇u ∈ L ∞ (Ω) we find u δ → u in L 1 (Ω) and ∇u δ → ∇u in L 2 (Ω). Moreover, we have Consequently, it suffices to establish a recovery sequence for a function u ∈ W(Ω) with u = g in a neighborhood ofΩ \ Ω and J u = J δ u for some δ > 0. Note after the above modification the segments of J u might not be pairwise disjoint.
We define u ε (x) = u(x) for x ∈ L ε ∩Ω and let y ε (x) = id + √ εu ε (x). Clearly we have u ε ∈ A gε for all ε. Byũ ε , u ′ ε we again denote the interpolations onΩ ε . Up to considering a translation of u of order ε, we may assume that J u ∩ L ε = ∅. Let D ε be the sets of triangles where J u crosses at least one side of the triangle. Then for a constant C > 0 independent of u ∈ W(Ω, R 2 ) and ε, where N u denotes the (smallest) number of segments whose union gives J u . From now on for the local nature of the arguments we may assume that J u consists of one segment only. Indeed, if J u consists of segments S 1 , . . . , S Nu , which are possibly not disjoint, the number of triangles ∆ ∈C ε with △ ∩ S i 1 ∩ S i 2 = ∅ for 1 ≤ i i < i 2 ≤ N u scales like N u and therefore their energy contribution is negligible in the limit. We show for ε small enough. Let △ ∈ D ε . We see that, if J u = J δ u crosses a spring v at point x * , say, then a computation similar as in (16) together with ∇u ∈ L ∞ shows Thus, △ ∈C ε for ε small enough. On the other hand, if we assume △ / ∈ D ε , then for at least two springs This is clear for △ / ∈ D ε =C ε as ∇u ∈ L ∞ . For △ ∈C ε,3 it follows by construction.  (9) we get a rotation R ε ∈ SO(2) such that . We note that χΩ εũ ε → u in L 1 as u and thus everyũ ε is bounded uniformly in L ∞ and, u being Lipschitz away from J u ,ũ ε → u uniformly onΩ ε \ △∈Dε △, where | △∈Dε △| ≤ Cε. Letting C * ε = D ε this shows that u ε → u in the sense of Definition 2.1 recalling (20) and the fact that |(ũ ε ) △ | = O(1) for △ / ∈ D ε . We next establish an even stronger convergence of the derivatives. Consider ∇ũ ε on triangles in C ε \ D ε . As u is Lipschitz there, the oscillation on such a triangle, x, x ′ ∈ △}, tends to zero uniformly (i.e., not depending on the choice of the triangle). We thus obtain for ε → 0, so that even χΩ ε\∪△∈D ε △ ∇ũ ε → ∇u strongly in L 2 (Ω). Note that in fact χΩ ε ∇u ′ ε → ∇u in L 2 (Ω). Indeed, recall #D ε ≤ Cε −1 by (20). Using (22) on the set of broken triangles we then get for ε → 0. We now split up the energy in bulk and surface parts as defined in (12). Note that indeed the contribution εE boundary ε is of order O(ε) as ∇u ∈ L ∞ (Ω) and J u ⊂ Ω since u = g in a neighborhood ofΩ \ Ω. We first observe that 1 ε Ωε χ(∇ỹ ε ) = 0 for ε small enough. Indeed, for ∆ ∈C ε this follows from (21). For ∆ / ∈ D ε it suffices to recall |(ũ ε ) △ | = O(1) which implies that (ũ ε ) △ is near SO (2). Repeating the steps in the elastic energy estimate in (i), applying χ Ωε ∇u ′ ε → ∇u strongly in L 2 (Ω), (22) and Q(F ) ≤ C|F | 2 for a constant C > 0 we conclude that It is elementary to see that J u crosses springs in v-direction for v ∈ V, where ν u is a normal to the segment J u . Recalling (21), the crack energy may be estimated by lim sup This together with (23) and (24) shows that u ε is a recovery sequence for u. Finally, we construct recovery sequences for the functionals F χ ε to conclude the proof of Theorem 2.3. Proof of Theorem 2.3, second part. Following the proof of Theorem 2.2(ii) it suffices to show lim First, by (21) and the definition of f κ we get ∆∈Dε ∆ f κ (∇ỹ ε ) = 0 for ε small enough. For ∆ / ∈ D ε we have (∇ỹ ε ) ∆ = (∇y ′ ε ) ∆ and thus we find f κ ((∇ỹ ε ) ∆ ) = −ε κ 2Q ((∇u ′ ε ) ∆ ) − κω( √ ε∇(u ′ ε ) ∆ ) by (19). We obtain 1 ε Ωε f κ (∇ỹ ε ) = 1 ε Ωε\ ∆∈Dε ∆ f κ (∇y ′ ε ) By Lemma 3.2 we obtain min{Q(F ) : e T 1 F e 1 = r} = α 2 r 2 . Then using the slicing method (see, e.g., [5, Section 3.11]) we get where #S x 2 denotes the number of jumps on a slice (0, l) × {x 2 } and In case #S x 2 (u) ≥ 1, the inner integral in (26) is obviously bounded from below by 2β γ . If #S x 2 (u) = 0, by applyig Jensen's inequality we find that this term is bounded from below by αla 2 due to the boundary conditions. We thus obtain inf E ≥ min αla 2 √ 3 , 2β γ . On the other hand, it is straighforward to check that E(u el ) = αla 2 and E(u cr ) = 2β γ , which shows that u el is a minimizer for a < a crit and u cr is a minimizer for a > a crit . It remains to prove uniqueness: (i) Let a < a crit and u be a minimizer of E. Since E(u) = E(u el ) we infer from (26) that u has no jump on a.e. slice (0, l) × {x 2 } and satisfies e T 1 ∇u e 1 = a a.e. by the imposed boundary values and strict convexity of the mapping t → t 2 on [0, ∞). Thus, if J u = ∅, a crack normal must satisfy ν u = ±e 2 H 1 -a.e. Taking E γ (u) and the fact that P (γ, e 2 ) > 0 for γ ∈ [ √ 3 2 , 1] into account, we then may assume J u = ∅ up to an H 1 negligible set, i.e., u ∈ H 1 (Ω). We find u 1 (x 1 , x 2 ) = ax 1 +f (x 2 ) a.e. for a suitable function f , and the boundary condition u 1 (0, x 2 ) = 0 yields f = 0 a.e. In particular, e T 1 ∇u e 2 = 0 a.e. Applying strict convexity of Q on symmetric matrices (Lemma 3.2) we now observe e T 2 ∇u e 2 = − a 3 and e T 1 ∇u e 2 + e T 2 ∇u e 1 = 0 a.e. So the derivative has the form for a.e. x.
(ii) Let a > a crit , φ = 0 and u be a minimizer of E. We again consider the lower bound (26) for the energy E and now obtain that on a.e. slice (0, l) × {x 2 } a minimizer u has precisely one jump and that e T 1 ∇u e 1 = 0 a.e. Now Lemma 3.2 shows that ∇u is antisymmetric a.e. As a consequence, the linearized rigidity estimate for SBD functions of Chambolle, Giacomini and Ponsiglione [13] yields that there is a Caccioppoli partition (E i ) of Ω such that where A T i = −A i ∈ R 2×2 and b i ∈ R 2 . (See [5] for the definition and basic properties of Caccioppoli partitions.) As E γ (u) = 0, we also note that ν u ⊥ v γ