Nonlinear elasticity with vanishing nonlocal self-repulsion

Abstract

We prove that for nonlinear elastic energies with strong enough energetic control of the outer distortion of admissible deformations, almost everywhere global invertibility as constraint can be obtained in the $\Gamma$-limit of the elastic energy with an added nonlocal self-repulsion term with asymptocially vanishing coefficient. The self-repulsion term considered here formally coincides with a Sobolev–Slobodeckiĭ seminorm of the inverse deformation. Variants near the boundary or on the surface of the domain are also studied.

1. Introduction

A basic requirement for many models of deformable solids is that they should prevent interpenetration of mass. In context of hyperelasticity, i.e., nonlinear elasticity fully determined by a stored elastic energy function (see e.g. [3, 10] for an introduction), this is ensured by a strong local resistance to compression built into the energy density, which in particular prevents local change of orientation, combined with a constraint preventing global self-penetration, usually the Ciarlet–Nečas condition [11], see (CN).

In this article, we study the approximation of the latter by augmenting the local elastic energy with a nonlocal functional with self-repulsive properties, formally corresponding to suitable Sobolev–Slobodeckiĭ seminorms of the inverse deformation. While all results presented here are purely analytical, our motivation is mainly numerical, related to the fact that the Ciarlet–Nečas condition is hard to handle numerically in such a way that the algorithm maintains an acceptable computational cost while still provably converging. In particular, there is still no known projection onto the Ciarlet–Nečas condition which is rigorous with acceptable computational cost, see [1] for some partial results. There is a well-known straightforward penalty term that rigorously reproduces the Ciarlet–Nečas condition in the limit (see e.g. [24]), but it is hard to implement, non-smooth and computationally very expensive as a double integral on the full domain. Recent results on more practical rigorous approximation of the Ciarlet–Nečas condition via nonlocal penalty terms added to the elastic energy were obtained in [21, 22], but these require additional regularity of elastic deformations which possibly interferes with the Lavrentiev phenomenon which is known to appear at least in particular nonlinear elastic models [14].

Using the language of $\Gamma$-convergence (see e.g. [8]), we show that in combination with local nonlinear elastic energies, the self-repulsive terms studied here also provide a rigorous approximation of the Ciarlet–Nečas condition without requiring regularity of deformations beyond what is naturally provided by the nonlinear elastic energy (theorem 3.3). In addition, these admit natural variants near or on the boundary (theorems 3.8 and 3.12), which are significantly cheaper to compute in practice. The latter crucially rely on a global invertibility property of orientation preserving maps exploiting topological information on the boundary [20] (for related results see also [2, 17]).

Our results here still do not cover the full range of hyperelastic energies which are known to be variationally well-posed, though. In fact, we require lower bounds on the energy density which are strong enough so that deformation maps with finite elastic energy are automatically continuous, open and discrete, the latter two by the theory of functions of bounded distortion [18]. In our proofs, this is essential so that all local regularity is controlled by the elastic energy, while the nonlocal self-repulsive term asymptotically only controls global self-contact.

For related results concerning self-avoiding curves and surfaces in more geometrical context with higher regularity, we refer to [5–7, 29] and references therein.

General assumptions

Let $\Omega \subset {\mathbb {R}}^{d}$ be a bounded Lipschitz domain, $d\ge 2$, $p\in (d,\,\infty )$, $r>0$, $ {\varepsilon }\ge 0$, $q\ge 1$ and $s\geq 0$. By $W_+^{1,p}(\Omega,\, {\mathbb {R}}^{d})$ we denote the set of all functions $y\in W^{1,p}(\Omega,\, {\mathbb {R}}^{d})$ with $\det \nabla y(x)>0$ for a.e. $x\in \Omega$.

We consider an integral functional modelling the internal elastic energy of a deformation $y\in W_+^{1,p}(\Omega,\, {\mathbb {R}}^{d})$ of a nonlinear hyperelastic solid. For simplicity, we restrict ourselves to the following generalized ‘neo-Hookean’ form given by

(1.1)\begin{equation} {\mathcal{E}}(y) = \int_{\Omega}{\left\lvert \nabla y(x) \right\rvert}^{p}{\,\mathrm{d}} x + \int_{\Omega}\frac{{\,\mathrm{d}} x}{{\left( \det \nabla y(x) \right)}^{r}}. \end{equation}

We are interested in precluding deformations corresponding to self-interpenetration of matter, i.e., non-injective $y$. Classically, the latter is imposed by adding the Ciarlet–Nečas condition [11] as a constraint:

(CN)\begin{align} \int_{\Omega}{\left\lvert \det \nabla y(x) \right\rvert}{\,\mathrm{d}} x= {\left\lvert y(\Omega) \right\rvert}. \end{align}

Remark 1.1 By the area formula, the inequality ‘$\ge$’ in (CN) always holds true and (CN) is equivalent to a.e. injectivity of $y$ provided that $\det \nabla y>0$ a.e. As the latter is usually given, (CN) can also be expressed in the more standard form

\[ \int_{\Omega}\det \nabla y(x){\,\mathrm{d}} x\le {\left\lvert y(\Omega) \right\rvert}. \]

As a step towards a possible (numerical) approximation of (CN), we regularize $ {\mathcal {E}}$ by adding a singular nonlocal contribution $ {\mathcal {D}}$. Below, we will show that (CN) automatically holds whenever $ {\mathcal {E}}(y)<\infty$ and $ {\mathcal {D}}(y)<\infty$ with various examples for $ {\mathcal {D}}$, see propositions 3.4, 3.9 and 3.14.

The first such example for $ {\mathcal {D}}$ is given by

(1.2)\begin{equation} {\mathcal{D}}_{U}(y) = \iint_{U\times U}\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{{\left\lvert y(x)-y(\tilde x) \right\rvert}^{d +sq}}{\left\lvert \det \nabla y(x) \right\rvert}{\left\lvert \det \nabla y(\tilde x) \right\rvert} {\,\mathrm{d}} x{\,\mathrm{d}}\tilde x \end{equation}

where $U\subset \Omega$ is some open neighbourhood of $\partial \Omega$ in $\Omega$ and suitable parameters $q\in [1,\,\infty )$, $s\in [0,\,1)$. In particular, we can choose $U=\Omega$. Transforming the integral and invoking [9, Prop. 2] reveals that the integral is singular if $s\ge 1$.

Formally, after a change of variables, $ {\mathcal {D}}_{U}$ is the Sobolev–Slobodeckiĭsemi-norm of $y^{-1}$ in the space $W^{s,q}(y(U),\, {\mathbb {R}}^{d})$. As long as $sq\ge 0$, the functional $ {\mathcal {D}}_{U}$ effectively prevents self-interpenetration, i.e., a loss of injectivity of $y$, as shown in proposition 3.4. To the best of our knowledge, variants of $ {\mathcal {D}}_{U}$ for curves, with such a purpose in mind, first appeared in a master thesis [27] supervised by Dziuk, and have subsequently been studied in another master thesis [19]. The functional $ {\mathcal {D}}_{U}$ can be interpreted to be a sort of ‘relaxation’ of the bi-Lipschitz constant. In this sense it is a rather weak quantity in comparison with similar concepts that have been introduced earlier [15, 25].

We also study the boundary variant

\[ {\widetilde{\mathcal{D}}}_{\partial\Omega}(y) = \iint_{\partial\Omega\times\partial\Omega}\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{{\left\lvert y(x)-y(\tilde x) \right\rvert}^{d-1 +sq}} {\,\mathrm{d}} A(x){\,\mathrm{d}} A(\tilde x) \]

where $A$ denotes the $(d-1)$-dimensional Hausdorff measure.

We give a rigorous statement of this approximation by establishing $\Gamma$-convergence which is the main result of this paper. To this end, we consider for $y\in W^{1,p}(\Omega,\, {\mathbb {R}}^{d})$, $ {\varepsilon }>0$,

\[ E_{{\varepsilon}} (y) = \begin{cases} {\mathcal{E}}(y) + {\varepsilon}{\mathcal{D}}(y) & \text{if }y\in W_+^{1,p}(\Omega,{\mathbb{R}}^{d}), \\ + \infty & \text{else}. \end{cases} \]

We reserve the symbol $E_{0}$ for the $\Gamma$-limit which will turn out to be

(1.3)\begin{equation} E_{0}(y) := \begin{cases} {\mathcal{E}}(y) & \text{if $y\in W_+^{1,p}(\Omega,\,{\mathbb{R}}^{d})$ satisfies (CN),} \\ + \infty & \text{else.} \end{cases} \end{equation}

2. Preliminary results

Change of variables in quite general form is important for us throughout, in form of the following special case of the area formula due to Marcus and Mizel. We use the convention that $f(y(x)){\left \lvert \det \nabla y(x) \right \rvert }=0$ whenever ${\left \lvert \det \nabla y(x) \right \rvert }=0$ for some $x\in E$ and abbreviate $N_{y}(z,\,E)=\#(y^{-1}(z)\cap E)$ for any $z\in {\mathbb {R}}^{d}$, where $\#$ denotes the counting measure.

Lemma 2.1 cf. [23, Theorem 2]

Let $y\in W^{1,p}(\Lambda,\, {\mathbb {R}}^{d})$ with $p>d$, where $\Lambda$ is a bounded domain in $ {\mathbb {R}}^{d}$. Moreover, assume that $f: {\mathbb {R}}^{d}\to {\mathbb {R}}$ is measurable and $E\subset \Lambda$ is measurable. Then, if one of the functions $x\mapsto f(y(x)){\left \lvert \det \nabla y(x) \right \rvert }$ and $z\mapsto f(z)N_{f}(z,\,E)$ is integrable, so is the other one and the identity

\[ \int_{E} f(y(x)){\left\lvert \det\nabla y(x) \right\rvert}{\,\mathrm{d}} x = \int_{{\mathbb{R}}^{m}} f(z)N_{f}(z,E) {\,\mathrm{d}} z \]

holds.

Proposition 2.2 For $p>d$ and $r>0$, the functional $ {\mathcal {E}}:W^{1,p}(\Omega ; {\mathbb {R}}^d)\to [0,\,\infty ]$ is lower semicontinuous with respect to weak convergence in $W^{1,p}$.

Proof. The integrand of $ {\mathcal {E}}$ is polyconvex, since $(F,\,J)\mapsto |F|^p+J^{-r}$, $ {\mathbb {R}}^{d\times d}\times (0,\,\infty )\to [0,\,\infty ],$ is convex. As shown in detail by Ball [3], sequential weak lower semicontinuity of $ {\mathcal {E}}$ therefore follows from the weak continuity of the determinant, i.e., $y\mapsto \det \nabla y$ as a map between $W^{1,p}(\Omega ; {\mathbb {R}}^d)$ and $L^{p/d}(\Omega )$, where both spaces are endowed with their weak topologies.

The Ciarlet–Nečas condition is a viable constraint for direct methods:

Lemma 2.3 Weak stability of (CN) [11, p. 185]

Let $y_{k}\rightharpoonup y_{\infty }$ in $W_+^{1,p}(\Omega,\, {\mathbb {R}}^{n})$, $p>d$, and assume that (CN) holds for all $y_{k}$, $k\in {\mathbb {N}}$. Then (CN) applies to $y_{\infty }$ as well.

Using the theory of maps of bounded distortion, we can obtain even more. For sufficiently large $p$ and $r$, deformations with finite energy are open and discrete (see below) due to a result of Villamor and Manfredi [28]. With added global topological information, say, in form of (CN), finite energy maps are even necessarily homeomorphisms [16, Section 3] (see also [20] for related results). In summary, we have the following.

Proposition 2.4 Let $d\geq 2$,

\[ p>d(d-1)\quad(\ge d),\qquad r>\frac{p(d-1)}{p-d(d-1)}\quad({>}d-1), \]

and let $y\in W_+^{1,p}(\Omega,\, {\mathbb {R}}^{d})$ such that $ {\mathcal {E}}(y)<\infty$. Then the continuous representative of $y$ is open (i.e., $y$ maps open subset of $\Omega$ to open sets in $ {\mathbb {R}}^{d}$) and discrete (i.e., for each $z\in {\mathbb {R}}^d$, $y^{-1}(\{z\})$ does not have accumulation points in $\Omega$). In particular, $y(\Omega )$ is open in $ {\mathbb {R}}^d$. If, in addition, (CN) holds, then $y$ is a homeomorphism on $\Omega$ and $y^{-1}\in W_+^{1,\sigma }(y(\Omega ),\,\Omega )$, where

\[ \sigma:=\frac{(r+1)p}{r(d-1)+p}>d. \]

Remark 2.5 Possible self-contact on $\partial \Omega$ is not ruled out, and so $y$ is not necessarily a homeomorphism on $\overline \Omega$.

Proof Proof of proposition 2.4

With $F=\nabla y(x)$,

\[ K^O(F):=\frac{{\left\lvert F \right\rvert}^d}{\det F} \]

is the outer distortion of $y$ at $x$ (or dilatation in the terminology of [28]). We infer from Young's inequality for some $\kappa,\,\rho \in (1,\,\infty )$ that

\[ {\left\lvert K^O(F) \right\rvert}^{\kappa} \leq C_{\rho}\left({\left\lvert F \right\rvert}^{d\kappa\rho} + {{\left( \frac1{\det F} \right)}^{\frac{\kappa\rho}{\rho-1}}}\right). \]

Given $y\in W_+^{1,p}(\Omega,\, {\mathbb {R}}^{d})$, we find that $K^O(\nabla y)\in L^{\kappa }(\Omega )$ provided $ {\mathcal {E}}(y)<\infty$ as well as ${\rm d}\kappa \rho =p$ and $\kappa =r(1-\frac 1\rho )$. Since $\frac 1\kappa =\frac 1r+\frac {\rm d}p<\frac {p-d(d-1)}{p(d-1)} + \frac {\rm d}p = \frac 1{d-1}$ and $\rho =1+\frac p{{\rm d}r}>1$ we may conclude that $y$ is open and discrete as shown by Villamor and Manfredi [28, Theorem 1].

Finally, (CN) implies that $y\in W_+^{1,p}$ is a map of (Brouwer's) degree $1$ for values of its image (e.g. $y\in \rm {DEG}1$ by [20, Remark 2.19(b)]). By [20, Thm. 6.8], it now follows that $y:\Omega \to y(\Omega )$ is a homeomorphism with weakly differentiable inverse, and $\nabla (y^{-1})\in L^d(y(\Omega ); {\mathbb {R}}^{d\times d})$.

Now that $y$ is invertible with weakly differentiable inverse, we may improve the last conclusion to $L^{\sigma }$. To this end, we first apply a change of variables and use $F^{-1}=\frac {\operatorname {cof} F}{\det F}$, whence $|F^{-1}|\leq c{\left \lvert F \right \rvert }^{d-1}|\det F|^{-1}$. Then the assertion follows again by invoking Young's inequality to bound ${\left \lvert F \right \rvert }^{(d-1)\sigma }|\det F|^{-(\sigma -1)}$.

The final two lemmas of this section will be crucial ingredients in the construction of a recovery sequence in the proof of theorem 3.3. They are also used in [21, 22] in similar fashion. For a closely related result and further references we refer to [4, Theorem 5.1].

Lemma 2.6 domain shrinking

Let $\Omega \subset {\mathbb {R}}^d$ be a bounded Lipschitz domain. Then there exists a sequence of $C^\infty$-diffeomorphisms

\[ \Psi_j:\overline{\Omega}\to \Psi_j(\overline{\Omega})\subset{\subset} \Omega \]

such that as $j\to \infty$, $\Psi _j\to \operatorname {id}$ in $C^m(\overline {\Omega }; {\mathbb {R}}^d)$ for all $m\in {\mathbb {N}}$.

Lemma 2.7 composition with domain shrinking is continuous

Let $\Omega \subset {\mathbb {R}}^d$ be a bounded Lipschitz domain, $k\in {\mathbb {N}}_0$, $1\leq p <\infty$ and $f\in W^{k,p}(\Omega ; {\mathbb {R}}^m)$, $m\in {\mathbb {N}}$. Let $\Psi _j$ be a sequence of maps as in lemma 2.6, we then have that $f\circ \Psi _j\to f$ in $W^{k,p}(\Omega ; {\mathbb {R}}^n)$.

Proof Proof of lemma 2.6

If $\Omega$ is strictly star-shaped with respect to a point $x_0\in \Omega$, one may take $\Psi _j(x):=x_0+\frac {j-1}{j}(x-x_0)$. For a general Lipschitz domain, one can combine local constructions near the boundary using a smooth decomposition of unity: if, locally in some open cube $Q$, the set $\Omega$ is given as a Lipschitz subgraph, i.e., $\Omega \cap Q=\{x\in Q\mid x\cdot e\leq f(x')\}$ and $\partial \Omega \cap Q=\{x'+ef(x')\mid x\in Q\}$, where $e$ is a unit vector orthogonal to one of the faces of $Q$, $x':=x-(x\cdot e) e$ and $f$ is a Lipschitz function, we define

\[ \hat{\Psi}_j(x;Q):=x'+\alpha_0 e+\frac{j-1}{j}(e\cdot x-\alpha_0)e\quad\text{for} \ x\in Q, \ \text{with} \ \alpha_0:=\inf_{x\in Q}e\cdot x. \]

Notice that $\hat {\Psi }_j(\cdot ;Q)$ pulls the local boundary piece $\partial \Omega \cap Q$ ‘down’ (in direction $-e$) into the original domain while leaving the ‘lower’ face of $Q$ fixed. Since $\partial \Omega$ can be covered by finitely many such cubes, we can write $\overline {\Omega }\subset Q_0\cup \bigcup _{k=1}^{n} Q_k$ with some open interior set $Q_0\subset \subset \Omega$. For a smooth decomposition of unity $1=\sum _{k=0}^n \varphi _k$ subordinate to this covering of $\Omega$ (i.e., $\varphi _k$ smooth, non-negative and compactly supported in $Q_k$),

\[ \Psi_j(x):=\varphi_0(x)x+\sum_{k=1}^n \varphi_j(x) \hat{\Psi}_j(x;Q_k) \]

now has the asserted properties.

Proof Proof of lemma 2.7

We only provide a proof for the case $k=1$, which will include the argument for $k=0$. For $k\geq 2$, the assertion follows inductively. It suffices to show that as $j\to \infty$, $\partial _n [f\circ \Psi _j-f] \to 0$ in $L^p$, for each partial derivative $\partial _n$, $n=1,\,\ldots,\, d$. By the chain rule,

(2.1)\begin{align} & \partial_n [f\circ \Psi_j-f]= [(\nabla f)\circ \Psi_j] \partial_n \Psi_j -\partial_n f\nonumber\\ & \quad =\Big([(\nabla f)\circ \Psi_j] \partial_n \Psi_j-(\partial_n f)\circ \Psi_j \Big) +\Big((\partial_n f)\circ \Psi_j -\partial_n f\Big). \end{align}

The first term above converges to zero in $L^p$ since $(\nabla f) e_n=\partial _n f$ for the $n$th unit vector $e_n$, and $\partial _n\Psi _j \to \partial _n\operatorname {id}=e_n$ uniformly. The convergence of the second term corresponds to our assertion for the case $k=0$, with $\tilde {f}:=\partial _n f\in L^p$. It can be proved in the same way as the well-known continuity of the shift in $L^p$: if $\tilde {f}$ is smooth and can be extended to a smooth function on $ {\mathbb {R}}^d$, we have

(2.2)\begin{equation} {\left\lVert \tilde{f}\circ \Psi_j -\tilde{f} \right\rVert}_{L^p(\Omega;{\mathbb{R}}^m)}\leq {\left\lVert \nabla\tilde{f} \right\rVert}_{L^\infty({\mathbb{R}}^d;{\mathbb{R}}^{m\times d})} {\left\lVert \Psi_j-\operatorname{id} \right\rVert}_{L^p(\Omega;{\mathbb{R}}^d)}\underset{j\to\infty}{\longrightarrow} 0. \end{equation}

The general case follows by approximation of $\tilde {f}$ in $L^p$ with such smooth functions, by first extending $\tilde {f}$ by zero to all of $ {\mathbb {R}}^d$, and then mollifying. Here, notice that for the mollified function, $\|\nabla \tilde {f}\|_{L^\infty }$ in (2.2) is unbounded in general as a function of the mollification parameter, but one can always choose the latter to converge slow enough with respect to $j$ so that (2.2) still holds.

3. Elasticity with vanishing nonlocal self-repulsion

In this section, we will study energies of the form

\[ E_{{\varepsilon}} (y) = \begin{cases} {\mathcal{E}}(y) + {\varepsilon} {\mathcal{D}}(y) & \text{if }y\in W_+^{1,p}(\Omega,{\mathbb{R}}^{d}), \\ + \infty & \text{else}, \end{cases} \]

where $ {\mathcal {E}}$ is defined in (1.1), in the limit $ {\varepsilon }\to 0^+$, in the sense of $\Gamma$-convergence with respect to the weak topology of $W^{1,p}$. Here, we say that $E_{ {\varepsilon }}$ $\Gamma$-converges to a functional $E_0$ if the following two properties hold for every sequence $ {\varepsilon }(k)\to 0^+$ and every $y\in W^{1,p}(\Omega ; {\mathbb {R}}^d)$:

1. (i) (lower bound) For all sequences $y_k\rightharpoonup y$ weakly in $W^{1,p}$,

   \[ \liminf_k E_{{\varepsilon}(k)}(y_k)\geq E_0(y). \]

2. (ii) (recovery sequence) There exists a sequence $y_k\rightharpoonup y$ weakly in $W^{1,p}$ s.t.

   \[ \limsup_k E_{{\varepsilon}(k)}(y_k)\leq E_0(y). \]

Remark 3.1 Notice that we do not require compactness here, i.e., that any sequence $(y_k)$ with bounded $E_{ {\varepsilon }(k)}(y_k)$ has a subsequence weakly converging in $W^{1,p}$. This is automatic as soon as bounded energy implies a bound in the norm of $W^{1,p}$. However, in the most basic form, the energies we study are translation invariant and only control $\nabla y$ but not $y$. Of course, this would change as soon as a Poincaré inequality can be used due to a suitable boundary condition or other controls on $y$ or its average added via constraint or additional energy terms of lower order.

We will discuss three different examples for $ {\mathcal {D}}$, each preventing self-interpenetration, i.e., a loss of injectivity of $y$, in a different way. Recall that, according to (1.3),

\[ E_{0}(y) = \begin{cases} {\mathcal{E}}(y) & \text{if} \ y\in W_+^{1,p}(\Omega,{\mathbb{R}}^{d})\ \text{satisfies}\; ({\textrm CN}), \\ + \infty & \text{else.} \end{cases} \]

Throughout this section, we assume that $d,\,p,\,r,\,\sigma$ satisfy the assumptions of proposition 2.4, namely,

(3.1)\begin{equation} d\geq 2,\quad p>d(d-1),\quad r>\frac{p(d-1)}{p-d(d-1)}, \quad \sigma=\frac{(r+1)p}{r(d-1)+p}. \end{equation}

3.1 Bulk self-repulsion

Here we consider the energy $E_ {\varepsilon }$ with $ {\mathcal {D}} := {\mathcal {D}}_{\Omega }$ introduced in (1.2), i.e.,

\[ {\mathcal{D}}(y) = {\mathcal{D}}_{\Omega}(y) = \iint_{\Omega\times\Omega}\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{{\left\lvert y(x)-y(\tilde x) \right\rvert}^{d +sq}}{\left\lvert \det \nabla y(x) \right\rvert}{\left\lvert \det \nabla y(\tilde x) \right\rvert} {\,\mathrm{d}} x{\,\mathrm{d}}\tilde x \]

for $y\in W^{1,p}(\Omega,\, {\mathbb {R}}^{d})$, $q\in [1,\,\infty )$, $s\in [0,\,1)$.

The following statement is actually not required for our main result. However, together with its counterpart for the elastic energy (cf. proposition 2.2) it ensures well-posedness of the variational model.

Proposition 3.2 For $p>d$, the functional $ {\mathcal {D}}_\Omega :W^{1,p}(\Omega ; {\mathbb {R}}^d)\to [0,\,\infty ]$ is lower semicontinuous with respect to the weak convergence in $W^{1,p}$.

Proof. For $\delta >0$ define

\[ {\mathcal{D}}^{[\delta]}(y):=\iint_{\Omega\times\Omega}\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{\max\{\delta,{\left\lvert y(x)-y(\tilde x) \right\rvert}^{d+sq}\}}{\left\lvert \det \nabla y(x) \right\rvert}{\left\lvert \det \nabla y(\tilde x) \right\rvert}\,{\rm d}x\,{\rm d}\tilde x\leq {\mathcal{D}}(y). \]

Let $(y_k)\subset W^{1,p}(\Omega ; {\mathbb {R}}^d)$ with $y_k\rightharpoonup y$ in $W^{1,p}$, for some $y\in W^{1,p}(\Omega ; {\mathbb {R}}^d)$. In particular, $y_k\to y$ in $C(\overline \Omega ; {\mathbb {R}}^d)$ by embedding, and consequently,

(3.2)\begin{equation} \liminf_{k\to\infty}{\mathcal{D}}^{[\delta]}(y_k)=\liminf_{k\to\infty} \iint_{\Omega\times\Omega} W_{\delta,y}(x,\tilde x,\det \nabla y_k(x),\det \nabla y_k(\tilde x))\,{\rm d}x\,{\rm d}\tilde x \end{equation}

where

\[ W_{\delta,y}(x,\tilde x,J,\tilde{J}):=\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{\max\{\delta,{\left\lvert y(x)-y(\tilde x) \right\rvert}^{d+sq}\}}{\left\lvert J \right\rvert}|\tilde{J}| \quad \text{for}\ x,\tilde x\in \Omega, \ J,\tilde{J}\in {\mathbb{R}}. \]

Clearly, $W_{\delta,y}$ is symmetric in $(x,\,\tilde x)$ and $(J,\,\tilde {J})$, as well as separately convex in $(J,\,\tilde {J})$, i.e., convex in $J$ with $x,\,\tilde x,\,\tilde {J}$ fixed and convex in $\tilde {J}$ with $x,\,\tilde x,\,J$ fixed. By [26, Theorem 2.5] (see also the related earlier result [12, Theorem 11]), this implies weak lower semicontinuity of $J\mapsto \iint _{\Omega \times \Omega } W_{\delta,y}(x,\,\tilde x,\,J(x),\,J(\tilde {x}))\,{\rm d}x\,{\rm d}\tilde x$ in $L^\alpha (\Omega )$, in particular for $\alpha :=\frac {p}{d}$. Again exploiting the weak continuity of the determinant, i.e., $J_k:=\det \nabla y_k\rightharpoonup J:=\det \nabla y$ weakly in $L^{p/d}$, we thus get that

(3.3)\begin{align} & \liminf_{k\to\infty}\iint_{\Omega\times\Omega} W_{\delta,y}(x,\tilde x,\det \nabla y_k(x),\det \nabla y_k(\tilde x))\,{\rm d}x\,{\rm d}\tilde x \nonumber\\ & \quad \geq \iint_{\Omega\times\Omega} W_{\delta,y}(x,\tilde x,\det \nabla y(x),\det \nabla y(\tilde x))\,{\rm d}x\,{\rm d}\tilde x & = {\mathcal{D}}^{[\delta]}(y). \end{align}

Combining (3.2) and (3.3), we see that $ {\mathcal {D}}^{[\delta ]}$ is weakly lower semicontinuous for each $\delta >0$. Since $ {\mathcal {D}}_\Omega (y)=\sup _{\delta >0} {\mathcal {D}}^{[\delta ]}(y)$, this implies weakly lower semicontinuity of $ {\mathcal {D}}_\Omega$.

Theorem 3.3 Let $\Omega \subset {\mathbb {R}}^d$, be a bounded Lipschitz domain, and suppose that $q\geq 1$ and $s\in [0,\,1)$. In addition, suppose that (3.1) holds together with

(3.4)\begin{equation} s-\frac {\rm d}q\le 1-\frac {\rm d}\sigma. \end{equation}

Then the functionals $E_{ {\varepsilon }}$ $\Gamma$-converge to $E_{0}$ as $ {\varepsilon }\searrow 0$, with respect to the weak topology of $W^{1,p}(\Omega,\, {\mathbb {R}}^{d})$.

For the proof, we additionally need the following two propositions.

Proposition 3.4 Finite $E_ {\varepsilon }(y)$ implies (CN)

Suppose that $s,\,q\geq 0$ and (3.1) holds, and let $y\in W_+^{1,p}(\Omega,\, {\mathbb {R}}^{d})$ such that $ {\mathcal {E}}(y)<\infty$ and $ {\mathcal {D}}(y)<\infty$. Then $y$ satisfies (CN).

Recall that a function is said to satisfy Lusin's property (N) if it maps sets of measure zero to maps of measure zero. In particular, it applies to Sobolev functions $W^{1,p}(\Omega,\, {\mathbb {R}}^{n})$ where $p>n$; cf. [14].

Proof. The proof is indirect. Suppose that (CN) does not hold. By the area formula (cf. Lemma 2.1), this means that $Z_2:=\{z\in {\mathbb {R}}^d\mid N_y(z,\,\Omega )\geq 2\}$ has positive measure. As a consequence, $X_2:=y^{-1}(Z_2)$ also has positive measure, because $y$ satisfies Lusin's property (N) as a map in $W^{1,p}$ with $p>d$. In addition, we claim that $X_2$ is open. For a proof, take any $x\in X_2\neq \emptyset$. By definition of $X_2$, there exists another point $\tilde {x}\in X_2\setminus \{x\}$ such that $y(x)=y(\tilde {x})$. If we choose disjoint open neighbourhoods $U,\,\tilde {U}\subset \Omega$ of $x,\,\tilde {x}$, respectively, then $y(U)$ and $y(\tilde {U})$ are open sets because $y$ is open by proposition 2.4. Hence, their intersection $y(U)\cap y(\tilde {U})\subset Z_2$ is also open, and it contains $y(x)=y(\tilde {x})$. By continuity of $y$, we conclude that $y^{-1}(y(U)\cap y(\tilde {U}))$ is now an open subset of $X_2$ containing $x$ (and $\tilde {x}$).

The above construction in particular shows that we can have two open, nonempty sets $V,\,W\subset X_2\subset \Omega$ with $x\in V\subset U\cap y^{-1}(y(U)\cap y(\tilde {U}))$ and $W:=\tilde {U}\cap y^{-1}(y(V))$ such that $\overline V\cap \overline W=\emptyset$ and $y(W)\subset y(V)$ are open. Hence, with

\[ \delta:=\min\{ {\lvert x-\tilde{x} \rvert} \mid x\in \overline V,~\tilde{x}\in \overline W\}>0, \]

we have that

\[ {\mathcal{D}}(y)\geq \iint_{V\times W}\frac{{\delta}^{q}}{{{\lvert y(x)-y(\tilde{x}) \rvert}}^{d+sq}} {\lvert \det\nabla y(x) \rvert}{\lvert \det\nabla y(\tilde{x}) \rvert} \,{\,\mathrm{d}} x{\,\mathrm{d}}\tilde x. \]

Changing variables in both integrals using lemma 2.1, also using that $N_y\geq 1$ on the image of $y$, we infer that

\[ {\mathcal{D}}(y)\geq {\delta}^{q}\int_{y(W)}\int_{y(V)} \frac{1}{{\lvert \xi-\tilde{\xi} \rvert}^{d+sq}} \,{\,\mathrm{d}}\xi{\,\mathrm{d}}\tilde\xi. \]

The double integral is infinite since $y(W)\subset y(V)$ and $sq\geq 0$. This implies that $ {\mathcal {D}}(y)=+\infty$, contradicting our assumption.

Proposition 3.5 Let $y\in W_+^{1,p}(\Omega,\, {\mathbb {R}}^{d})$ be a homeomorphism $\Omega \to y(\Omega )$ with $y^{-1}\in W^{1,\sigma }(y(\Omega ),\,\Omega )$, $q\in [1,\,\infty )$, $s\in [0,\,1)$. In addition, suppose that (3.1) and (3.4) hold. If $\Omega '$ is open and $\Omega '\subset \subset \Omega$ then $ {\mathcal {D}}_{\Omega '}(y|_{\Omega '})<\infty$.

Proof. We apply lemma 2.1 twice to change variables in each of the two integrals in $\mathcal {D}_{\Omega '}$, with $E=\Omega '$. First, change variables in the inner integral, say, over $x$, using $z=\xi$ and $f(z)={{\lvert y^{-1}(z)-\tilde {x} \rvert }^{q}}{{\lvert z-y(\tilde {x}) \rvert }^{-(d+sq)}}$ for any fixed $\tilde {x}\in \Omega '$. Afterwards, use Fubini's theorem to change the order of integration and change variables in the integral over $\tilde {x}$, now for any fixed $\xi$ using $z=\tilde {\xi }$ and $f(z)={{\lvert y^{-1}(\xi )-y^{-1}(z) \rvert }^{q}}{{\lvert \xi -z \rvert }^{-(d+sq)}}$. We thus obtain that

\[ \mathcal{D}_{\Omega'}(y|_{\Omega'}) =\iint_{y(\Omega')\times y(\Omega')}\frac{{\lvert y^{{-}1}(\xi)-y^{{-}1}(\tilde \xi) \rvert}^{q}}{{\lvert \xi-\tilde\xi \rvert}^{d+sq}} {\,\mathrm{d}}\xi{\,\mathrm{d}}\tilde\xi. \]

The right-hand side is just the $q$th power of the seminorm belonging to the Sobolev–Slobodeckiĭ space $W^{s,q}(y(\Omega '),\, {\mathbb {R}}^{d})$. As $\overline {y(\Omega ')}\subset y(\Omega )$ we may find $\psi \in C^{\infty }( {\mathbb {R}}^{d})$ supported in $y(\Omega )$ with $\psi =1$ on $y(\Omega ')$. Choosing any regular value $r\in (0,\,1)$, the set $\psi ^{-1}((r,\,1])$ has a smooth boundary. We denote its component containing $y(\Omega ')$ by $\Upsilon$. In case $s\in (0,\,1)$, using $y(\Omega ')\subset \Upsilon \subset y(\Omega )$ and applying the embedding theorem, we infer

\begin{align*} \mathcal{D}_{\Omega'}(y|_{\Omega'}) & \le {\left[y^{{-}1}|_{y(\Omega')} \right]}_{W^{s,q}(y(\Omega'),{\mathbb{R}}^{d})}^{q} \le {\left[y^{{-}1}|_{\Upsilon} \right]}_{W^{s,q}(\Upsilon,{\mathbb{R}}^{d})}^{q} \\ & \le C_{d,q,\sigma,\Upsilon}{\left\lVert y^{{-}1}|_{\Upsilon} \right\rVert}_{W^{1,\sigma}(\Upsilon,{\mathbb{R}}^{d})}^{q} \le C_{d,q,\sigma,\Upsilon}{\left\lVert y^{{-}1}|_{y(\Omega)} \right\rVert}_{W^{1,\sigma}(\Omega,{\mathbb{R}}^{d})}^{q}. \end{align*}

The case $s=0$ is similar; in the intermediate step above, we now use $W^{\tilde {s},\,q}$ with some $\tilde {s}>0$ small enough so that $W^{1,\sigma }$ still embeds into $W^{\tilde {s},\,q}$ since $\sigma > \frac {{\rm d}q} {q+d}$.

Proof Proof of theorem 3.3

Lower bound ( $\mathit \Gamma$-lim inf-inequality): Assume that $y_{k}\rightharpoonup y$ in $W^{1,p}$ and $ {\varepsilon }_{k}\searrow 0$. If $\liminf _k E_{ {\varepsilon }_k}(y_k)=+\infty$, there is nothing to show. Hence, passing to a suitable subsequence (not relabelled), we may assume that the $\liminf$ is a limit and $E_{ {\varepsilon }_k}(y_k)$ is bounded. Since $ {\mathcal {D}}\ge 0$ and $ {\mathcal {E}}$ is weakly lower semicontinuous, we get that

(3.5)\begin{equation} \lim_k E_{{\varepsilon}_k}(y_k)\geq \liminf_k {\mathcal{E}}(y_k)\geq {\mathcal{E}}(y). \end{equation}

Moreover, by proposition 3.4, we see that $y_k$ satisfies (CN) for all $k$, and by lemma 2.3, this implies that $y$ satisfies (CN). Hence, $ {\mathcal {E}}(y)=E_0(y)$, and (3.5) thus implies the asserted lower bound.

Upper bound (construction of a recovery sequence): Let $y\in W^{1,p}(\Omega ; {\mathbb {R}}^d)$ be given. We may assume that $E_0(y)<\infty$, because otherwise there is nothing to show. We therefore have that $y\in W_+^{1,p}(\Omega ; {\mathbb {R}}^d)$, $y$ satisfies (CN) and $ {\mathcal {E}}(y)<\infty$.

By proposition 2.4, the map $y:\Omega \to y(\Omega )$ is a homeomorphism. We choose $j\in {\mathbb {N}}$ and shrink the domain $\Omega$ to $\Omega _{j}=\Psi _j(\Omega )$, using lemma 2.6. Now define

\[ y_{j}:\Omega\to y(\Omega_{j}), \qquad y_{j}:=y|_{\Omega_{j}}\circ \Psi_j. \]

As $j\to \infty$, $y_{j}\to y$ in $W^{1,p}(\Omega ; {\mathbb {R}}^d)$ and $ {\mathcal {E}}(y_{j})\to {\mathcal {E}}(y)$ by lemma 2.7. Here, concerning the term $(\det \nabla y)^{-r}$ in $ {\mathcal {E}}$, notice that by the chain rule and the multiplicativity of the determinant, we have that

\[ \frac{1}{(\det\nabla y_j)^{r}} =(f\circ\Psi_j) \frac{1}{(\det\nabla \Psi_j)^{r}},\quad\text{with} \ f:=\frac{1}{(\det\nabla y)^{r}}\in L^1(\Omega). \]

Combined with the fact that $\det \nabla \Psi _j\to 1$ uniformly, lemma 2.7 can therefore indeed be applied with $k=0$ and $p=1$ to obtain convergence of this singular term in $ {\mathcal {E}}$.

Since $\Omega _j\subset \subset \Omega$ for each $j$, we also have that

\[ {\mathcal{D}}_{\Omega}(y_j)={\mathcal{D}}_{\Omega_j}(y) <\infty \]

by change of variables and proposition 3.5.

Now let $ {\varepsilon }_{k}\searrow 0$ be given. We choose $j_{k}\to \infty$ such that $ {\varepsilon }_{k} {\mathcal {D}}(y_{j_{k}})\xrightarrow {k\to \infty }0$. Here, notice that this can be always achieved by choosing $j_{k}\to \infty$ slow enough (depending on $ {\varepsilon }_{k}$), even if $ {\mathcal {D}}(y_{j})\to +\infty$. So $E_{ {\varepsilon }_{k}}(y_{j_{k}}) = {\mathcal {E}}(y_{j_{k}}) + {\varepsilon }_{k} {\mathcal {D}}(y_{j_{k}})\xrightarrow {k\to \infty } {\mathcal {E}}(y) = E_{0}(y)$.

3.2 Bulk self-repulsion near the boundary

Here, we consider $ {\mathcal {D}} := {\mathcal {D}}_{U_{\delta }}$ where, for any $\delta >0$, the set $U_{\delta }\subset \Omega$ can be chosen as any open neighbourhood of $\partial \Omega$ which is at least $\delta$-thick in the sense that

(3.6)\begin{equation} U_{\delta}\supset {\left( \partial\Omega \right)}^{(\delta)} = {\left\{x\in\Omega \,\middle|\, \operatorname{dist}(x,\partial\Omega)<\delta \right\}}. \end{equation}

For the ease of notation, we abbreviate $ {\mathcal {D}}_{\delta } := {\mathcal {D}}_{U_{\delta }}$ so that

\[ {\mathcal{D}}(y) = {\mathcal{D}}_{\delta}(y) = \iint_{U_{\delta}\times U_{\delta}}\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{{\left\lvert y(x)-y(\tilde x) \right\rvert}^{d +sq}}{\left\lvert \det \nabla y(x) \right\rvert}{\left\lvert \det \nabla y(\tilde x) \right\rvert} {\,\mathrm{d}} x{\,\mathrm{d}}\tilde x \]

for $y\in W^{1,p}(\Omega,\, {\mathbb {R}}^{d})$, $q\in [1,\,\infty )$, $s\in [0,\,1)$. The combined energy $E_ {\varepsilon }$ now also depends on the choice of $U_\delta$, and to make this more visible, we will now write

\[ E_{{\varepsilon},\delta} (y) = \begin{cases} {\mathcal{E}}(y) + {\varepsilon} {\mathcal{D}}_\delta(y) & \text{if }y\in W_+^{1,p}(\Omega,{\mathbb{R}}^{d}), \\ + \infty & \text{else}. \end{cases} \]

We will see that $E_0$ is still the correct limit functional independently of $\delta$. In fact, we can even allow the simultaneous limit as $( {\varepsilon },\,\delta )\to (0,\,0)$.

Remark 3.6 The fact that the limit as $\delta \to 0^+$ is admissible offers an attractive choice of $U_\delta$ for numerical purposes: a single boundary layer of the triangulation, which requires $\delta$ of the order of the grid size $h$. In that case, the cost of a numerical evaluation of $ {\mathcal {D}}_\delta$ scales like $h^{-2(d-1)}$ (like a double integral on the surface), which is much cheaper than for $ {\mathcal {D}}_\Omega$ which scales like $h^{-2d}$.

As before, for fixed $ {\varepsilon },\,\delta >0$, the functional $E_{ {\varepsilon },\delta }$ is well suited for minimization by the direct method:

Proposition 3.7 For $p>d$, the functional $ {\mathcal {D}}_\delta :W^{1,p}(\Omega ; {\mathbb {R}}^d)\to [0,\,\infty ]$ is lower semicontinuous with respect to weak convergence in $W^{1,p}$.

Proof. This is proposition 3.2 with $\Omega$ replaced by $U_\delta$. Here, notice that no boundary regularity of $U_\delta$ is required: if needed, we can cover $U_\delta$ from inside with open domains with smooth boundary, and since the integrand of $ {\mathcal {D}}_\delta$ is nonnegative, we can therefore write $ {\mathcal {D}}_\delta$ as a supremum of weakly lower semicontinuous functionals using the smooth smaller domains as domain of integration.

To prove the main result in this and the following subsection, we will employ results of [20] which require that $\Omega$ does not have ‘holes’ as made precise in the following statement.

Theorem 3.8 Let $\Omega \subset {\mathbb {R}}^d$ be a bounded Lipschitz domain such that $ {\mathbb {R}}^d\setminus \partial \Omega$ has exactly two connected components, $q\in [1,\,\infty )$, $s\in [0,\,1)$. In addition, suppose that (3.1) and (3.4) hold. Then for any $\delta _0\in [0,\,\infty ]$, the functionals $E_{ {\varepsilon },\delta }$ $\Gamma$-converge to $E_{0}$ as $( {\varepsilon },\,\delta )\to (0,\,\delta _0)$ ($ {\varepsilon },\,\delta >0$), with respect to the weak topology of $W^{1,p}(\Omega,\, {\mathbb {R}}^{d})$.

For the proof, we additionally need the following modification of proposition 3.4.

Proposition 3.9 Finite $E_{ {\varepsilon },\delta }(y)$ implies (CN)

Let $\Omega \subset {\mathbb {R}}^d$ be a bounded Lipschitz domain such that $ {\mathbb {R}}^d\setminus \partial \Omega$ has exactly two connected components, $q\in [1,\,\infty )$, $s\in [0,\,1)$. In addition, suppose that (3.1) holds. Moreover, let $y\in W_+^{1,p}(\Omega,\, {\mathbb {R}}^{d})$ such that $ {\mathcal {E}}(y)<\infty$ and $ {\mathcal {D}}_\delta (y)<\infty$. Then $y$ satisfies (CN) on $\Omega$.

Proof. Analogously to proposition 3.4, we infer that $y$ satisfies (CN) on $U_\delta$. By proposition 2.4, we obtain that $y:U_\delta \to y(U_\delta )$ is a homeomorphism. Here, notice that for this conclusion, we do not need any regularity of the boundary of $U_\delta$, since it is enough to apply proposition 2.4 with $\Omega$ replaced by subdomains of $U_\delta$ with smooth boundary, and a sequence of such subdomains covers $U_\delta$ from the inside.

Such an inner covering can also be used for $\Omega$: choose open $\Omega _{j}\nearrow \Omega$ such that $\partial \Omega _{j}$ is smooth, say, Lipschitz. In addition, using that $\partial \Omega$ itself is also Lipschitz, we can make sure that as for $\Omega$, we have that $ {\mathbb {R}}^d\setminus \partial \Omega _j$ has exactly two connected components. For any fixed $\delta$, there exists a sufficiently large $j$ such that $\partial \Omega _{j}$ is contained in the open $\delta$-neighbourhood $(\partial \Omega )^{(\delta )}$ of $\partial \Omega$, and therefore $\partial \Omega _{j}\subset U_\delta$. Consequently, $y|_{\partial \Omega _j}$ is injective, which implies that $y\in \operatorname {AIB}(\Omega _j)$ in the sense of [20, Def. 2.1 and 2.2]. In addition, we know that $y\in W_+^{1,p}(\Omega _j; {\mathbb {R}}^d)$. By [20, Thm. 6.1 and Rem. 6.3], we infer that $y$ satisfies (CN) on $\Omega _j$. As the latter holds for all $j$, we conclude that $y$ satisfies (CN) on $\Omega$ by monotone convergence.

Remark 3.10 The proof of proposition 3.9 exploits that $y$ is a homeomorphism near the boundary, which we obtain from proposition 2.4. This forces the relatively restrictive assumptions on $p$ and $r$. While this may be technical to some degree, some restrictions are definitely needed. In fact, by itself, (CN) on a boundary strip like $U_\delta$ is not strong enough to provide the necessary global topological information: if one can squeeze surfaces to points with a deformation of finite elastic energy (possible if $p$ and $r$ are small enough), then self-penetration on $U_\delta$ is indeed possible for a $y\in W_+^{1,p}$ which is injective on $U_\delta$ outside a set of dimension $d-1$. Such a set of $d$-dimensional measure zero is invisible to (CN).

Proof Proof of theorem 3.8

Lower bound ( $\mathit \Gamma$-lim inf-inequality): This is completely analogous to the proof of theorem 3.3, using proposition 3.9 instead of proposition 3.4.

Upper bound (construction of a recovery sequence): Again, we can follow the proof of theorem 3.3 step by step, using the domain shrinking maps $\Psi _j$ to define $y_j:=y\circ \Psi _j$ as before. In particular, changing variables we now observe that

\[ {\mathcal{D}}_{\delta}(y_j)={\mathcal{D}}_{\Psi_j(U_\delta)}(y) \leq {\mathcal{D}}_{\Omega_j}(y) <\infty \]

by proposition 3.5, for any fixed $j$. Given $( {\varepsilon }_k,\,\delta _k)\to (0,\,\delta _0)$, we thus again get a suitable recovery sequence given by $(y_{j(k)})$ as long as $j(k)\to \infty$ slow enough so that $ {\varepsilon }_k {\mathcal {D}}_{\delta _k}(y_{j(k)})\to 0$.

3.3 Surface self-repulsion

Here we look at $ {\mathcal {D}}:= {\widetilde {\mathcal {D}}}_{\partial \Omega }$ where

\[ {\widetilde{\mathcal{D}}}_{\partial\Omega}(y) = \iint_{\partial\Omega\times\partial\Omega}\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{{\left\lvert y(x)-y(\tilde x) \right\rvert}^{d-1 +sq}} {\,\mathrm{d}} A(x){\,\mathrm{d}} A(\tilde x), \quad q\in[1,\infty), s\in[0,1), \]

and $A(\cdot )$ denotes the $(d-1)$-dimensional Hausdorff measure. Again, this is a term well-suited for minimization via the direct method:

Proposition 3.11 For $p>d$, the functional $ {\widetilde {\mathcal {D}}}_{\partial \Omega }:W^{1,p}(\Omega ; {\mathbb {R}}^d)\to [0,\,\infty ]$ is lower semicontinuous with respect to weak convergence in $W^{1,p}$.

Proof. As the trace operator from $W^{1,p}(\Omega,\, {\mathbb {R}}^{d})$ to $L^1(\partial \Omega,\, {\mathbb {R}}^{d})$ is compact, we obtain pointwise a.e. convergence of the integrand. So the claim immediately follows from Fatou's lemma.

Theorem 3.12 Let $\Omega \subset {\mathbb {R}}^d$ be a bounded Lipschitz domain such that $ {\mathbb {R}}^{d}\setminus \partial \Omega$ has exactly two connected components. Given (3.1), we require $q\geq 1$ and $s\in [0,\,1]$ to be chosen such that

(3.7)\begin{equation} q {\left( (1-s)\sigma-d \right)}>d^2-\sigma \end{equation}

and

(3.8)\begin{equation} sq\ge (d-1)\frac{p+d}{p-d}. \end{equation}

Then the functionals $E_{ {\varepsilon }}$ $\Gamma$-converge to $E_{0}$ as $ {\varepsilon }\searrow 0$, with respect to the weak topology of $W^{1,p}(\Omega,\, {\mathbb {R}}^{d})$.

Remark 3.13 With $\sigma =\sigma (r,\,p,\,d)>d$ as defined in (3.1), the conditions (3.7) and (3.8) are met if $0< s<1-\frac {\rm d}\sigma$ and $q>\max \left \{\frac {d^2-\sigma }{(1-s)\sigma -d},\,\frac {d-1}s\cdot \frac {p+d}{p-d}\right \}$

For the proof, we additionally need the following two propositions.

Proposition 3.14 Finite $E_ {\varepsilon }(y)$ implies (CN)

Let $\Omega \subset {\mathbb {R}}^d$ be a bounded Lipschitz domain such that $ {\mathbb {R}}^{d}\setminus \partial \Omega$ has exactly two connected components. Suppose that (3.1) and (3.8) hold. Moreover, let $y\in W_+^{1,p}(\Omega,\, {\mathbb {R}}^{d})$ such that $ {\mathcal {E}}(y)<\infty$ and $ {\widetilde {\mathcal {D}}}_{\partial \Omega }(y)<\infty$. Then $y$ satisfies (CN).

Proof. According to [20, Cor. 6.5] it is enough to show injectivity of $y|_{\partial \Omega }$. If the latter is not the case, we may choose $x_{0},\,\tilde x_{0}\in \partial \Omega$, $x_{0}\ne \tilde x_{0}$, such that $y(x_{0})=y(\tilde x_{0})$. Recalling that $y\in C^{0,\alpha }(\overline \Omega,\, {\mathbb {R}}^{d})$, $\alpha =1-\frac {\rm d}p$, and abbreviating $ {\varepsilon }=\tfrac 13{\left \lvert x_{0}-\tilde x_{0} \right \rvert }$, $t=d-1+sq$, we infer

\begin{align*} {\widetilde{\mathcal{D}}}_{\partial\Omega}(y) & =\iint_{\partial\Omega\times\partial\Omega}\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{{\left\lvert y(x)-y(\tilde x) \right\rvert}^{t}} {\,\mathrm{d}} A(x){\,\mathrm{d}} A(\tilde x) \\ & \ge\iint_{(\partial\Omega\cap B_{{\varepsilon}}(x_{0}))\times(\partial\Omega\cap B_{{\varepsilon}}(\tilde x_{0}))}\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{{\left( {\left\lvert y(x)-y(x_{0}) \right\rvert}+{\left\lvert y(\tilde x_{0})-y(\tilde x) \right\rvert} \right)}^{t}} {\,\mathrm{d}} A(x){\,\mathrm{d}} A(\tilde x) \\ & \ge c_{\alpha}\iint_{(\partial\Omega\cap B_{{\varepsilon}}(x_{0}))\times(\partial\Omega\cap B_{{\varepsilon}}(\tilde x_{0}))}\frac{{\left( \frac{\varepsilon}3 \right)}^{q}}{{\left( {\left\lvert x-x_{0} \right\rvert}^{\alpha}+{\left\lvert \tilde x_{0}-\tilde x \right\rvert}^{\alpha} \right)}^{t}} {\,\mathrm{d}} A(x){\,\mathrm{d}} A(\tilde x). \end{align*}

Introducing local bi-Lipschitz charts $\Phi :V\to \partial \Omega \cap B_{ {\varepsilon }}(x_{0})$, $\tilde \Phi :\tilde V\to \partial \Omega \cap B_{ {\varepsilon }}(\tilde x_{0})$ where $V,\,\tilde V\subset {\mathbb {R}}^{d-1}$ are open sets and $\Phi (0)=x_{0}$, $\tilde \Phi (0)=\tilde x_{0}$, we arrive at

\begin{align*} {\widetilde{\mathcal{D}}}_{\partial\Omega}(y) & \ge c_{\alpha,{\varepsilon},q,t}\iint_{V\times\tilde V}\frac{\sqrt{\det D\Phi(\xi)^{\top}D\Phi(\xi)}\sqrt{\det D\tilde\Phi(\tilde\xi)^{\top}D\tilde\Phi(\tilde\xi)}}{{\left( {\left\lvert \Phi(\xi)-\Phi(0) \right\rvert}^{2}+{\left\lvert \tilde\Phi(\tilde\xi)-\tilde\Phi(0) \right\rvert}^{2} \right)}^{\alpha t/2}} {\,\mathrm{d}}\xi{\,\mathrm{d}}\tilde\xi \\ & \ge c_{\alpha,{\varepsilon},q,t,\Phi}\iint_{V\times\tilde V}\frac{{\,\mathrm{d}}\xi{\,\mathrm{d}}\tilde\xi}{{\left( {\left\lvert \xi \right\rvert}^{2}+{\left\lvert \tilde\xi \right\rvert}^{2} \right)}^{\alpha t/2}}. \end{align*}

By assumption, both $V$ and $\tilde V$ contain $B_{\delta }(0)\subset {\mathbb {R}}^{d-1}$ for some $\delta >0$. Decomposing $(\xi ^{\top },\,\tilde \xi ^{\top })=\rho \eta ^{\top }\in {\mathbb {R}}^{2d-2}$ where $\rho >0$, $\eta \in \mathbb {S}^{2d-3}$, yields

\[ {\widetilde{\mathcal{D}}}_{\partial\Omega}(y) \ge c_{\alpha,{\varepsilon},q,t,\Phi,d}\int_{0}^{\delta}\frac{\rho^{2d-3}}{\rho^{\alpha t}}{\,\mathrm{d}}\rho. \]

This term is infinite provided $2d-3-\alpha t\le -1$ which is equivalent to (3.8). This contradicts our assumption that $ {\widetilde {\mathcal {D}}}_{\partial \Omega }(y)$ is finite.

Proposition 3.15 Let $\Omega \subset {\mathbb {R}}^d$ be a bounded Lipschitz domain. Assume that $y\in W_+^{1,p}(\Omega,\, {\mathbb {R}}^{d})$ is a homeomorphism $\Omega \to y(\Omega )$ with $y^{-1}\in W^{1,\sigma }(y(\Omega ),\,\Omega )$ for which (3.7) applies. If $\Omega '$ is open and $\Omega '\subset \subset \Omega$ and $\partial \Omega '$ is Lipschitz, then $ {\widetilde {\mathcal {D}}}_{\partial \Omega '}(y)<\infty$.

Proof. First notice that $y(\Omega )$ is open and bounded in $ {\mathbb {R}}^{d}$, the former by invariance of domain (see e.g. [13, Theorem 3.30]) and the latter due to the fact that $y\in C(\overline {\Omega }; {\mathbb {R}}^{d})$ by embedding. Hence, $y(\overline {\Omega '})$ is a compact and connected subset of $y(\Omega )$ with positive distance to $\partial [y(\Omega )]$. We choose a domain $\Lambda \subset {\mathbb {R}}^{d}$ with smooth boundary such that $y(\Omega ')\subset \subset \Lambda \subset \subset y(\Omega )$. By embedding, $y^{-1}\in C^{0,\beta }(\Lambda,\, {\mathbb {R}}^{d})$, $\beta =1-\frac {\rm d}\sigma$. Abbreviating $t=d-1+sq$, we arrive at

\begin{align*} {\widetilde{\mathcal{D}}}_{\partial\Omega'}(y) & =\iint_{\partial\Omega'\times\partial\Omega'}\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{{\left\lvert y(x)-y(\tilde x) \right\rvert}^{t}} {\,\mathrm{d}} A(x){\,\mathrm{d}} A(\tilde x) \\ & \le C_{\beta}\iint_{\partial\Omega'\times\partial\Omega'}{\left\lvert x-\tilde x \right\rvert}^{q-t/\beta} {\,\mathrm{d}} A(x){\,\mathrm{d}} A(\tilde x). \end{align*}

The term ${\left \lvert x-\tilde x \right \rvert }$ is bounded above since $\Omega$ is bounded. It approaches zero only in a neighbourhood of the diagonal. In order to show that $ {\widetilde {\mathcal {D}}}_{\partial \Omega '}(y)$ is finite we only have to consider $\iint _{\Phi (V)\times \Phi (V)}{\left \lvert x-\tilde x \right \rvert }^{q-t/\beta } {\,\mathrm {d}} A(x) {\,\mathrm {d}} A(\tilde x)$ where $\Phi :V\to U\subset \partial \Omega '$ is a chart and $V\subset B_{R}(0)\subset {\mathbb {R}}^{d-1}$ is an open set. Decomposing $\xi =\rho \eta \in {\mathbb {R}}^{d-1}$ where $\rho >0$, $\eta \in \mathbb {S}^{d-2}$, yields

\begin{align*} & \iint_{\Phi(V)\times\Phi(V)}{\left\lvert x-\tilde x \right\rvert}^{q-t/\beta} {\,\mathrm{d}} A(x){\,\mathrm{d}} A(\tilde x) \\ & =\iint_{V\times V}{\left\lvert \Phi(\xi)-\Phi(\tilde\xi) \right\rvert}^{q-t/\beta} {\sqrt{\det D\Phi(\xi)^{\top}D\Phi(\xi)}\sqrt{\det D\Phi(\tilde\xi)^{\top}D\Phi(\tilde\xi)}} {\,\mathrm{d}}\xi{\,\mathrm{d}}\tilde\xi \\ & \le C_{\Phi}\iint_{V\times V}{\left\lvert \xi-\tilde\xi \right\rvert}^{q-t/\beta}{\,\mathrm{d}}\xi{\,\mathrm{d}}\tilde\xi \\ & \le C_{\Phi}\int_{B_{R}(0)}\int_{B_{R}(0)}{\left\lvert \xi-\tilde\xi \right\rvert}^{q-t/\beta} {\,\mathrm{d}}\xi{\,\mathrm{d}}\tilde\xi \\ & \le C_{\Phi}\int_{B_{3R}(0)}\int_{B_{R}(\tilde\xi)}{\left\lvert \xi \right\rvert}^{q-t/\beta} {\,\mathrm{d}}\xi{\,\mathrm{d}}\tilde\xi \\ & \le C_{\Phi,d,R}\int_{0}^{R}\rho^{d-1+q-t/\beta}{\,\mathrm{d}}\rho. \end{align*}

The right-hand side is finite if $d-1+q-t/\beta >-1$ which is equivalent to (3.7).

Proof Proof of theorem 3.12

We proceed as in the proof of theorem 3.3. For the lower bound we use proposition 3.14 in place of proposition 3.4. To see that the recovery sequence also works for $ {\widetilde {\mathcal {D}}}_{\partial \Omega }$, we compute

\begin{align*} {\widetilde{\mathcal{D}}}_{\partial\Omega}(y_{j}) & ={\widetilde{\mathcal{D}}}_{\partial\Omega}(y|_{\Omega_{j}}\circ \Psi_j) \\ & =\iint_{\partial\Omega\times\partial\Omega}\frac{{\left\lvert \xi-\tilde\xi \right\rvert}^{q}}{{\left\lvert y(\Psi_{j}(\xi))-y(\Psi_{j}(\tilde\xi)) \right\rvert}^{d-1+sq}} {\,\mathrm{d}} A(\xi){\,\mathrm{d}} A(\tilde\xi) \\ & \le C_{\Psi_{j}}\iint_{\partial\Omega_j\times\partial\Omega_j}\frac{{\left\lvert \Psi_{j}^{{-}1}(x)-\Psi_{j}^{{-}1}(\tilde x) \right\rvert}^{q}}{{\left\lvert y(x)-y(\tilde x) \right\rvert}^{d-1+sq}} {\,\mathrm{d}} A(x){\,\mathrm{d}} A(\tilde x) \\ & \le C_{\Psi_{j}}{\left\lVert \nabla\Psi_{j}^{{-}1} \right\rVert}_{L^{\infty}}^{q}\iint_{\partial\Omega_j\times\partial\Omega_j}\frac{{\left\lvert x-\tilde x \right\rvert}^{q}}{{\left\lvert y(x)-y(\tilde x) \right\rvert}^{d-1+sq}} {\,\mathrm{d}} A(x){\,\mathrm{d}} A(\tilde x) \end{align*}

where $C_{\Psi _{j}}$ denotes a factor that bounds the terms arising from the change of variables. Now we deduce from proposition 3.15 (instead of proposition 3.5) that the right-hand side is finite.

3.4 Further generalizations and remarks

Remark 3.16 More general elastic energies

It is easy to see that throughout, the integrand of $ {\mathcal {E}}$ can be replaced by any polyconvex function admitting the original integrand as a lower bound (up to multiplicative and additive constants). Moreover, the latter is only exploited for the application of the theory of functions of finite distortion in proposition 2.4. More precisely, theorems 3.3, 3.8 and 3.12 also hold for any elastic energy of the form

\[ {\mathcal{E}}(y) = \int_{\Omega} W(\nabla y(x))\,{\,\mathrm{d}} x \]

such that

1. (i) $W: {\mathbb {R}}^{d\times d}\to (-\infty,\,+\infty ]$ is continuous and polyconvex, and $W(F)<\infty$ if and only if $\det F>0$,

2. (ii) $W(F)\geq c |F|^p-C$ for all $F\in {\mathbb {R}}^{d\times d}$, where $p>d$,

3. (iii) $W(F)\geq c (\frac {|F|^{d}}{\det F})^\beta -C$ for all $F\in {\mathbb {R}}^{d\times d}$ with $\det F>0$, where $\beta >d-1$

Here, $p>d,\,\beta >d-1$, $c>0$, and $C\in {\mathbb {R}}$ are constants. Notice that (iii) directly provides the bound on the outer distortion we need to generalize proposition 2.4.

Remark 3.17 Boundary conditions and force terms

Due to the stability of $\Gamma$-convergence with respect to addition of continuous functionals, our main results continue to hold if $ {\mathcal {E}}$ is modified by adding a term which is continuous in the weak topology of $W^{1,p}$ (typically either linear or lower order, exploiting a compact embedding). This includes many classical force potentials for body forces and surface tractions. Additional boundary conditions, say, a Dirichlet condition of the form $y=y_0$ on a part $\Lambda$ of $\partial \Omega$, are in principle also possible but not trivial to add, as they require modified recovery sequences in the proof of the theorems. In particular, we would need a suitable modification of lemma 2.6 which keeps the Dirichlet part of the boundary fixed, as well as additional assumptions on $y_0$ which at the very least should map $\overline {\Lambda }$ to a reasonably smooth set out of self-contact. The easiest way to set up a meaningful model with full coercivity in $W^{1,p}$ which is compatible with our theorems is to confine the deformed material to a box by constraint ($y(\Omega )\subset \mathcal {B}$ for a given compact $\mathcal {B}\subset {\mathbb {R}}^d$ with non-empty interior).

Remark 3.18 More general nonlocal self-repulsive terms

It is clear that our general proof strategy can also be applied to other nonlocal terms $ {\mathcal {D}}$. The only key features of such a term $ {\mathcal {D}}$ are the following:

1. (i) for any deformation $y$ with finite elastic energy $ {\mathcal {E}}(y)$, finite $ {\mathcal {D}}(y)$ implies (CN) (cf. propositions 3.4, 3.9 and 3.14);

2. (ii) for any homeomorphisms $y\in W^{1,p}_+$ ($p>d$) whose inverse has the Sobolev regularity $W^{1,\sigma }$ ($\sigma >d$) obtained from the control of its distortion through the elastic energy (see proposition 2.4), we obtain $ {\mathcal {D}}(y)<\infty$, at least if we move to a slightly smaller domain $\Omega '\subset \subset \Omega$ (cf. propositions 3.5 and 3.15).

Moreover, it is in principle possible to work with added penalty terms with properties closer to ones used in [21, 22], where, unlike here, the added terms are finite on all a-priori admissible deformations, even those that exhibit self-interpenetration. For instance, we could introduce an everywhere finite $ {\mathcal {D}}_ {\varepsilon }$ instead of $ {\varepsilon } {\mathcal {D}}$ as, say, a suitable truncation of the latter. For sufficiently small positive $ {\varepsilon }$, the self-repulsive property (i) in such a scenario can still hold for deformations satisfying a fixed energy bound. Here, the basic idea is to find at least one deformation $y_0$ so that $e_0:= {\mathcal {E}}(y_0)+\sup _{0< {\varepsilon }\leq 1} {\mathcal {D}}_ {\varepsilon }(y_0)<\infty$, for instance the identity or another map far from self-contact. Then check if (i) still holds in such a case if we replace the assumption $ {\mathcal {D}}(y)<\infty$ by $ {\mathcal {E}}(y)+ {\mathcal {D}}_ {\varepsilon }(y)\leq e_0$ (for sufficiently small $ {\varepsilon }$ independently of $y$).

Remark 3.19 Mosco-convergence and recovery by homeomorphisms

Our proofs of theorems 3.3, 3.8 and 3.12 actually provide more than $\Gamma$-convergence: the recovery sequence we construct always converges strongly in $W^{1,p}$, which means that we actually proved so called Mosco-convergence. Moreover, as constructed, each member of the recovery sequence is a homeomorphism on $\overline {\Omega }$. In particular, any admissible $y$ with finite $E_0(y)$ is always contained in the $C^0$-closure of these homeomorphisms, i.e., $y\in AI(\overline {\Omega })$ in the notation of [20]. Our results here therefore also show that within $W^{1,p}_+(\Omega ; {\mathbb {R}}^d)$ with $p>d$, $AI(\overline {\Omega })$ coincides the class of maps satisfying (CN) if we also impose strong enough a-priori bounds on the outer distortion to apply the result of Villamor and Manfredi as in proposition 2.4. The general case is still not clear, cf. [20, Remark 2.19].