First-Order Shape Derivative of the Energy for Elastic Plates with Rigid Inclusions and Interfacial Cracks

Abstract

Within the framework of Kirchhoff–Love plate theory, we analyze a variational model for elastic plates with rigid inclusions and interfacial cracks. The main feature of the model is a fully coupled nonpenetration condition that involves both the normal component of the longitudinal displacements and the normal derivative of the transverse deflection of the crack faces. Without making any artificial assumptions on the crack geometry and shape variation, we prove that the first-order shape derivative of the potential deformation energy is well defined and provide an explicit representation for it. The result is applied to derive the Griffith formula for the energy release rate associated with crack extension.

1 Introduction

Recent and ongoing advances in engineering and material science have increased the need for mathematical tools in order to design and optimize in an efficient way three-dimensional highly inhomogeneous thin structures.

In this paper, we confine our attention to the model of a composite Kirchhoff–Love plate proposed in [17]. More precisely, we consider an elastic plate containing a partially debonded rigid inclusion of regular shape. The displacement field of the inclusion is specified by an equality constraint and determined completely by six unknown constants. In order to eliminate nonphysical interpenetration of the crack faces, an inequality constraint that involves both the normal component of the longitudinal displacements and the normal derivative of the transverse deflection of the crack faces is imposed. This implies a full coupling between the longitudinal displacements and transverse deflection of the plate. Our interest is in investigating the differentiability of the potential deformation energy with respect to the shape variations of the reference configuration and in finding an explicit form of the corresponding first-order shape derivative. The motivation for such a study lies in the fact that this derivative can be used to calculate efficiently the energy release rates associated with variation of defects, which are utmost of importance to predict crack propagation [4, 5, 12, 20].

The principal difficulty in the problem at hand comes from the equality and inequality constraints, which do not preserve their structure under domain variations. To overcome this difficulty, we develop an approach that is purely based on the minimizing properties of variational solutions. The most important idea behind our technique originates from [40], where the Griffith formula for the energy release rate associated with crack extension was deduced in a rigorous way for two-dimensional elastic bodies with curvilinear cracks subjected to the Signorini conditions. Following this idea, in the case of Signorini-type constraints imposed on crack faces, general results on the shape differentiability of the potential deformation energy supported by explicit formulae for the first-order shape derivative were derived in [41] for linear elastic materials, in [32] for a Mindlin–Timoshenko elastic plate model, and in [42] for a Kirchhoff–Love elastic plate model. Another approach to the shape differentiability of the potential deformation energy based on primal-dual Lagrange formulations of the nonlinear crack problems goes back to [25, 26] and the subsequent [15]. Recent developments of the primal-dual shape sensitivity analysis in smooth domains and some applications to fluid mechanics problems with divergence-free constraints were given in [7, 13, 28, 47]. We finally refer to [1, 24] for a derivation of the Griffith formula for elastic materials with cracks by using the Piola transform. We stress here that all of the cited works related to solid mechanics, except [42], tackle inequality constraints on the crack faces that do not involve the gradient of displacement fields, while in the model under consideration the nonpenetration condition depends on the normal derivative of the transverse deflection. This additional difficulty essentially affects the computation of the shape derivative of the potential deformation energy.

In the context of Kirchhoff–Love plate theory, the nonpenetration condition for vertical cracks was put forward in [16]. A dimension reduction from 3D to 2D in the Kirchhoff–Love energy scaling regime for rate-independent fracture processes has been the subject of research in two past decades. To our best knowledge, all of the results are restricted to the use of the notion of global energetic solutions [34]; here we just mention a couple of them. With the help of evolutionary $\Gamma$-convergence [35], a two-dimensional model for the propagation of a single crack in a plate from a three-dimensional linearly elastic model involving a Barenblatt-like cohesive crack surface energy was obtained in [8]. A similar scenario for deriving two-dimensional models in the cases of Griffith-type delamination and an adhesive contact was performed in [9]. The advantages of the energetic formulation are that the initial problems and their two-dimensional counterparts have a derivative-free form and hence do not require any formulae for the energy release rate associated with crack extension; moreover, the crack path is not prescribed in advance.

In the presence of a preexisting interfacial crack between strongly dissimilar materials, it is natural to assume that the interface is the weakest fracture path and to apply the Griffith criterion [12] in order to describe propagation of the crack along the prescribed path. We recall that the Griffith criterion states: a crack is stationary if the energy release rate associated with crack extension is less than the fracture toughness of a material, and otherwise the crack grows. The main result of the present study gives a closed formula for the first-order shape derivative of the potential deformation energy (Theorem 2.2), and its application consists in deducing the Griffith formula for the energy release rate associated with extension of the crack along the interface, under minimal assumptions on the regularity of the crack path (Remark 2.11). Recent results in the same spirit of our work are [29, 39, 43, 45].

Fracture behavior depends on many factors, among them the mutual dislocation and interaction of inclusions and cracks. In a series of papers [23, 30, 31, 44], this issue was treated by an optimization technique based on regular shape variations. A more delicate shape-topological analysis of the nonlinear crack problems aimed at retarding or avoiding the crack propagation process was carried out in [22, 27, 33, 49]. Inverse problems for the nondestructive determination of elastic and rigid inclusions embedded in Kirchhoff–Love bending elastic plates were investigated in [3, 36, 37].

The layout of the paper is as follows. In Sect. 1, we discuss the mechanical model and present the functional framework used to assert the existence of the variational solutions. After the necessary preliminaries regarding the velocity method, which we employ to describe shape variations of the reference configuration, we formulate and prove in Sect. 2 our main result. The paper closes with a short discussion of the case when the elastic properties of the plate are inhomogeneous and an application of the result to fracture mechanics: we derive the Griffith formula for the energy release rate associated with crack extension.

2 Formulation of the Problem

We consider a composite plate, which is a planar thin-walled structure, within the confines of linear Kirchhoff–Love theory. The 3D undamaged plate in the reference configuration is the right cylinder

$$\begin{aligned} \Omega ^\mathrm{3D}=\left\{ \, (x_1,x_2,x_3)\in \mathbb {R}^3:x=(x_1,x_2)\in \Omega ,\;x_3\in (-h,h)\,\right\} , \end{aligned}$$

where $\Omega \subset \mathbb {R}^2$ defines the midsurface of the plate and 2h denotes the thickness, which we assume to be constant. The midsurface $\Omega$ is a bounded simply connected domain with boundary $\partial \Omega$ of class $C^{1,1}$, and let $\omega _0$ be a compact simply connected subdomain of $\Omega$ (i.e., $\omega _0\Subset \Omega$), with boundary $\partial \omega _0$ of class $C^{1,1}$. The boundary $\partial \omega _0$ is the union of two disjoint simple curves $\gamma _0$ and $\partial \omega _0\setminus \gamma _0$, where $\gamma _0$ is a relatively open set. The outward pointing unit normal to $\partial \omega _0$ as well as to $\partial \Omega$ is denoted by $\nu _0=(\nu _{01},\nu _{02})^\mathrm{T}$. Figure 1 provides a schematic visualization of the setup.

Fig. 1

Schematic of the midsurface of the plate

We make also the following assumption on the geometry of the non-Lipschitz domain $\Omega _0=\Omega \setminus \overline{\gamma }_0$:

G1:

There exist two domains $\Omega _1$ and $\Omega _2$ with Lipschitz boundaries such that $\overline{\Omega }_1\cup \overline{\Omega }_2=\overline{\Omega },$ $\gamma _0\subset \partial \Omega _1\cap \partial \Omega _2,$ and $\mathcal H^1(\partial \Omega _i\cap \partial \Omega )>0,$ $i=1,2,$ where $\mathcal H^1$ stands for the one-dimensional Hausdorff measure.

The region $(\Omega \setminus \overline{\omega }_0)^\mathrm{3D}$ is occupied by a linearly elastic homogeneous material whose properties are characterized by a fourth-order tensor $\mathbb {C}^\mathrm{3D}$ that is supposed to be positive, that is, there exists $c_{\mathbb {C}^\mathrm{3D}}>0$ such that

$$\begin{aligned} \mathbb {C}^\mathrm{3D}X:X\ge c_{\mathbb {C}^\mathrm{3D}}|X|^2\quad \text {for all}\quad X\in \mathbb {R}^{3\times 3}_\mathrm{sym}, \end{aligned}$$

symmetric

$$\begin{aligned} \mathbb {C}^\mathrm{3D}_{ijkl}=\mathbb {C}^\mathrm{3D}_{ijlk}=\mathbb {C}^\mathrm{3D}_{klij},\quad i,j,k,l=1,2,3, \end{aligned}$$

and monoclinic symmetric with respect to the $(x_1,x_2)$-plane, which implies

$$\begin{aligned} \mathbb {C}^\mathrm{3D}_{ijk3}=\mathbb {C}^\mathrm{3D}_{i333}=0,\quad i,j,k=1,2. \end{aligned}$$

(1)

We refer the reader to Remark 2.10 for some comments on the choice of the elasticity tensor. The region $\omega _0^\mathrm{3D}$ corresponds to a rigid inclusion, while the cylindrical surface $\gamma _0^\mathrm{3D}$ is an interfacial crack. As usual in the theory of thin-walled structures, we deal with quantities that refer to the midsurface of the undeformed plate $\Omega _0$. Let $W=(w_1,w_2)^\mathrm{T}$ be the longitudinal displacements of the midsurface points in the $(x_1,x_2)$-plane, and w is the vertical deflection of the midsurface points along the $x_3$-axis. In what follows, we use a variational approach to equilibrium and thus look for (W, w) belonging to the Cartesian product of Sobolev spaces $H^1(\Omega _0)^2\times H^2(\Omega _0)$.

Let $\mathbb {C}$ be the fourth-order tensor with components

$$\begin{aligned} \mathbb {C}_{ijkl}=\mathbb {C}^\mathrm{3D}_{ijkl}-\frac{\mathbb {C}^\mathrm{3D}_{ij33}\mathbb {C}^\mathrm{3D}_{kl33}}{\mathbb {C}^\mathrm{3D}_{3333}},\quad i,j,k,l=1,2. \end{aligned}$$

Hence the resultant membrane stress tensor $\sigma$ can be expressed in terms of the membrane strain tensor $\varepsilon (W)=1/2(\nabla _xW+(\nabla _x W)^\mathrm{T})$ via

$$\begin{aligned} \sigma =2h\mathbb {C}\varepsilon (W), \end{aligned}$$

(2)

and the bending moment tensor m is related to the bending strain tensor $-\nabla ^2_x w$ through the constitutive equations

$$\begin{aligned} m=-\frac{2h^3}{3}\mathbb {C}\nabla ^2_x w. \end{aligned}$$

(3)

The restriction of a displacement field (W, w) on the domain $\omega _0$ is an element of the Cartesian product $R_{p}(\omega _0)\times R_{v}(\omega _0),$ where $R_{p}(\omega _0)$ and $R_{v}(\omega _0)$ are the spaces of infinitesimal rigid displacement fields on $\omega _0$ that represent $\ker \varepsilon$ and $\ker \nabla ^2_x,$ respectively:

$$\begin{aligned} R_{p}(\omega _0)= & {} \left\{ \,W\in H^1(\omega _0)^2:\varepsilon (W)=0\;\text {in}\;\omega _0\,\right\} ,\\ R_{v}(\omega _0)= & {} \left\{ \,w\in H^2(\omega _0):\nabla ^2_x w=0\;\text {in}\;\omega _0\,\right\} . \end{aligned}$$

We observe that the following identities hold true [48, Lemma 1.1; Chapter III, Subsection 2.1.1]

$$\begin{aligned} R_{p}(\omega _0)= & {} \left\{ \,W(x)=B\mathrm{Id}_x(x)+C:B\in \mathbb {R}^{2\times 2},\;B=-B^\mathrm{T},\; C\in \mathbb {R}^{2}, \;x\in \omega _0\,\right\} ,\\ R_{v}(\omega _0)= & {} \left\{ \,w(x)=A\cdot \mathrm{Id}_x(x)+a:A\in \mathbb {R}^2,\;a\in \mathbb {R},\;x\in \omega _0\,\right\} . \end{aligned}$$

The regularity assumptions on the domains $\omega _0$ and $\Omega \setminus \overline{\omega }_0$ provide that the trace operators

$$\begin{aligned}&\mathrm{Tr}_{k,\omega _0}:H^k(\omega _0)^{k_*}\rightarrow \prod \limits _{i=0}^{k-1} H^{k-i-1/2}(\partial \omega _0)^{k_*},\\&\quad \mathrm{Tr}_{k,\Omega \setminus \overline{\omega }_0}:H^k(\Omega \setminus \overline{\omega }_0)^{k_*}\rightarrow \prod \limits _{i=0}^{k-1} H^{k-i-1/2}(\partial \omega _0;\partial \Omega )^{k_*}, \end{aligned}$$

where the superscript $k_*$ is equal to 2 if $k=1$ and 1 if $k=2,$ are well defined and continuous. Moreover, these operators have right continuous inverses (lifting operators) [14, Theorem 1.5.1.2], which we denote by $\mathrm{Tr}^{-1}_{k,\omega _0}$ and $\mathrm{Tr}^{-1}_{k,\Omega \setminus \overline{\omega }_0}$, respectively. According to the positive and negative directions of the unit normal $\nu _0$ on $\partial \omega _0$, there is a positive face $\partial \omega ^+_0$ and a negative face $\partial \omega ^-_0$. The double brackets $\llbracket v\rrbracket =v^{+}-v^{-}$ stand for the jump of a function v across $\partial \omega _0$. For convenience, the same notation serves as well for the jump of functions across another smooth curves, it always clear from the context. In order to eliminate interpenetration of the crack faces, we impose the inequality restriction

(4)

which is understood in the sense of traces. We emphasize that, from a mechanical point of view, the nonpenetration condition (4) describes the interaction of the opposite crack faces, admitting their frictionless contact (the equality case) or the absence of contact (the inequality case), the contact surface being unknown a priori.

We finally suppose that the plate is clamped along its outer edge, which leads to the homogeneous boundary conditions

$$\begin{aligned} W=(0,0)^\mathrm{T},\quad w= \nabla _x w\cdot \nu _0=0\quad \text {on} \quad \partial \Omega . \end{aligned}$$

The potential deformation energy of the system associated with a displacement field (W, w) is the stored elastic energy minus the work of external loadings $(F,f)\in C^1(\overline{\Omega })^2\times C^1(\overline{\Omega })$:

(5)

The variational principle of minimum potential deformation energy forces that the composite plate $\Omega ^\mathrm{3D}_0$ clamped along the outer edge, with elastic part $(\Omega \setminus \overline{\omega }_0)^\mathrm{3D}$, rigid part $\omega ^\mathrm{3D}_0$, interfacial crack $\gamma ^\mathrm{3D}_0$, and potential deformation energy given by (5), is in equilibrium if the displacement field of its midsurface $(W_0,w_0)$ is a solution of the minimization problem

$$\begin{aligned} \inf \limits _{({W},{w})\in K_0(\Omega _0)}\mathcal {E}(\Omega _0;{W},{w}), \end{aligned}$$

(6)

where the set of kinematically admissible displacement fields reads as follows

The convex set $K_0(\Omega _0)$ is closed in the Cartesian product $H^1_{\partial \Omega }(\Omega _0)^2\times H^2_{\partial \Omega }(\Omega _0)$ of the Sobolev spaces

$$\begin{aligned} H^1_{\partial \Omega }(\Omega _0)= & {} \left\{ \,v\in H^1(\Omega _0):v=0\; \text {on}\;\partial \Omega \,\right\} ,\\ H^2_{\partial \Omega }(\Omega _0)= & {} \left\{ \,v\in H^2(\Omega _0):v=\nabla _x v\cdot \nu _0=0\; \text {on}\;\partial \Omega \,\right\} . \end{aligned}$$

Since $\mathcal {E}$ is convex and Gâteux-differentiable, the minimization problem (6) is equivalent to the variational inequality

$$\begin{aligned}&(W_0,w_0)\in K_0(\Omega _0),\quad \int \limits _{\Omega \setminus \overline{\omega }_0} \sigma (W_0):\varepsilon (\overline{W}-W_0)\,\mathrm{d}x-\int \limits _{\Omega \setminus \overline{\omega }_0} m (w_0): \nabla ^2_x (\overline{w}-w_0)\,\mathrm{d}x\nonumber \\&\qquad -\int \limits _{\Omega _0} F\cdot (\overline{W}-W_0)\,\mathrm{d}x-\int \limits _{\Omega _0} f(\overline{w}-w_0) \,\mathrm{d}x\ge 0\quad \forall (\overline{W},\overline{w})\in K_0(\Omega _0). \end{aligned}$$

(7)

Thanks to G1, the Korn inequality and the Friedrichs inequality are valid in $\Omega _i,$ $i=1,2$, which implies the coercivity of $\mathcal {E}.$ Moreover, the energy functional $\mathcal {E}$ is weakly lower semicontinuous. The minimization problem (6) (and hence the variational inequality (7)) therefore admits a solution $(W_0,w_0)\in K_0(\Omega _0).$ It is straightforward to check that the solution to (7) is unique.

We are now ready to give the classical formulation of the equilibrium problem. We seek $(W_0,w_0)$ in the domain $\Omega _0,$ $(W^R_0,w^R_0)\in R_p(\omega _0)\times R_v(\omega _0)$, and $(\sigma ,m)$ in the domain $\Omega \setminus \overline{\omega }_0$ satysfying the following equations and boundary conditions

(8)

(9)

(10)

(11)

(12)

(13)

(14)

(15)

(16)

(17)

(18)

(19)

Here the traction $\sigma ^+ \nu _0$, transverse force $t^+_{\nu _0}$, and bending moment $m^+_{\nu _0}$ acting on $\partial \omega _0$ are

$$\begin{aligned}&\sigma ^+ \nu _0=(\sigma _{1j}^+\nu _{0j},\sigma _{2j}^+\nu _{0j})^\mathrm{T},\quad t^+_{\nu _0}=-(\mathrm{div}_x m)^+\cdot \nu _0-\partial _{\tau _0}\left( \left( m^+\nu _0\right) \cdot \tau _0\right) ,\\&\quad m^+_{\nu _0}=-(m^+\nu _0)\cdot \nu _0,\quad \tau _0=(\tau _{01},\tau _{02})^\mathrm{T}=(-\nu _{02},\nu _{01})^\mathrm{T}, \end{aligned}$$

and the traction $\sigma ^+ \nu _0$ is decomposed into normal and tangential components given by

$$\begin{aligned} \sigma ^+_{\nu _0}=\sigma ^+\nu _0\cdot \nu _0,\quad \sigma ^+_{\tau _0}=\sigma ^+\nu _0-\sigma ^+_{\nu _0}\nu _0. \end{aligned}$$

A repeated Roman subscript signifies summation from 1 to 2. Equations (8) and (10) are equilibrium equations for the elastic part of the plate with respect to the $(x_1, x_2)$-plane and $x_3$-axis, respectively. Conditions (13)–(15) are the complementarity conditions for contact of the crack faces. The nonlocal boundary conditions (18) and (19) ensure that the rigid part of the plate in equilibrium, that is, the resultant force and moment acting on it vanish.

It is a routine matter to verify that our minimizer $(W_0,w_0)\in K_0(\Omega _0)$ is indeed a weak solution of the Euler–Lagrange system (8)–(19), see [17].

Remark 1.1

We elucidate briefly in what sense the natural boundary conditions (14)–(19) for the minimizer $(W_0,w_0)$ are satisfied. For each $k=1,2,$ we introduce the Lions–Magenes space

$$\begin{aligned} H^{k-1/2}_{00}(\gamma _0)=\left\{ \,v\in H^{k-1/2}_0(\gamma _0):d^{-1/2}\partial ^{i}_{\tau _0} v\in L^2(\gamma _0),\; 0\le i\le k-1\right\} , \end{aligned}$$

where $\partial ^0_{\tau _0}v=v,$ $\partial ^1_{\tau _0}v=\partial _{\tau _0}v,$ and $d(s)=\mathrm{dist}_{\partial \omega _0}(s,\partial \gamma _0)$. This space equipped with the norm

$$\begin{aligned} \Vert v\Vert _{H^{k-1/2}_{00}(\gamma _0)}=\left( \Vert v\Vert ^2_{H^{k-1/2}(\gamma _0)}+\sum \limits _{i=0}^{k-1}\Vert d^{-1/2}\partial ^{i}_{\tau _0} v\Vert ^2_{L^2(\gamma _0)}\right) ^\frac{1}{2} \end{aligned}$$

is a Banach space. Since $\partial \omega _0\in C^{1,1}$, the following equivalence holds [18, Proposition 1.1]

$$\begin{aligned} v\in H^{k-1/2}_{00}(\gamma _0)\quad \Leftrightarrow \quad \overline{v}={\left\{ \begin{array}{ll} v\quad \text {on}\quad \gamma _0\\ 0\quad \text {on}\quad \partial \omega _0\setminus {\gamma }_0 \end{array}\right. }\in H^{k-1/2}(\partial \omega _0). \end{aligned}$$

Denote by $\langle \cdot ,\cdot \rangle ^{00}_{{k-1/2},\gamma _0}$ the duality pairing between the space $H^{k-1/2}_{00}(\gamma _0)$ and its topological dual space $H^{-k+1/2}_{00}(\gamma _0).$ If $\sigma ,m\in L^2(\Omega \setminus \overline{\omega }_0)^{2\times 2},$ $\mathrm{div}_x\sigma \in L^2(\Omega \setminus \overline{\omega }_0)^{2},$ and $\mathrm{div}_x \mathrm{div}_x m\in L^2(\Omega \setminus \overline{\omega }_0),$ then

$$\begin{aligned} \sigma ^+_{\tau _0}\in H^{-1/2}_{00}(\gamma _0)^2,\quad \sigma ^+_{\nu _0},m^+_{\nu _0}\in H^{-1/2}_{00}(\gamma _0),\quad t^+_{\nu _0}\in H^{-3/2}_{00}(\gamma _0), \end{aligned}$$

and the boundary conditions (14) and (15) read as follows [17]

$$\begin{aligned}&\left\langle \sigma ^+_{\tau _{0}}, v\right\rangle ^{00}_{1/2,\gamma _0}=0\quad \forall \,v\in H^{1/2}_{00}(\gamma _0)^2,\quad {v\cdot \nu _0=0},\\&\qquad \left\langle t^+_{\nu _{0}}, v\right\rangle ^{00}_{3/2,\gamma _0}=0\quad \forall \,v\in H^{3/2}_{00}(\gamma _0),\\&\left\langle \sigma ^+_{\nu _0}\pm h^{-1}m^+_{\nu _0}, v\right\rangle ^{00}_{1/2,\gamma _0}\le 0\quad \forall \,v\in H^{1/2}_{00}(\gamma _0),\quad v\ge 0,\\&\left\langle \sigma ^+_{\nu _{0}}, \llbracket W_0\rrbracket \cdot \nu _0\right\rangle ^{00}_{1/2,\gamma _0}+\left\langle m^+_{\nu _{0}}, \llbracket \nabla _x w_0\rrbracket \cdot \nu _0\right\rangle ^{00}_{1/2,\gamma _0}=0, \end{aligned}$$

whereas the interpretation of the nonlocal conditions (18) and (19) is provided by

$$\begin{aligned}&\left\langle \sigma ^+\nu _0,v\right\rangle _{1/2,\partial \omega _0}=-\int \limits _{\omega _0}F\cdot v\,\mathrm{d} x\quad \forall \,v\in R_p(\omega _0),\\&\left\langle t^+_{\nu _0}, v\right\rangle _{3/2,\partial \omega _0}-\left\langle m^+_{\nu _0},\nabla _x v\cdot \nu _0\right\rangle _{1/2,\partial \omega _0}=\int \limits _{\omega _0}fv\,\mathrm{d} x\quad \forall \,v\in R_v(\omega _0). \end{aligned}$$

Here $\langle \cdot ,\cdot \rangle _{k-1/2,\partial \omega _0}$ is the duality pairing between the spaces $H^{k-1/2}(\partial \omega _0)$ and $H^{-k+1/2}(\partial \omega _0)$.

We exploit some of these conditions to derive an explicit formula for the first-order shape derivative of the potential deformation energy, see Theorem 2.2.

3 Shape Derivative of the Energy

We begin this section by recalling the velocity (speed) method [46, Section 2.9], which we employ to describe shape variations of the reference domain $\Omega _0$. Let $V\in W^{2,\infty }_\mathrm{loc}(\mathbb {R}^2)^2$ be a given shape velocity field. We assume further that V is compactly supported in the domain $\Omega$. This assumption does not limit the generality of our analysis, but it simplifies the presentation. The following proposition is a core of the velocity method and allows us to construct a family of diffeomorphisms of $\mathbb {R}^2$ generated by V.

Proposition 2.1

Let $V\in W^{2,\infty }_\mathrm{loc}(\mathbb {R}^2)^2$ be a given shape velocity field satisfying the property $\mathrm{supp}V\Subset {\Omega }.$ Then there exists $\delta _0>0$ such that:

• (a) the Cauchy problem

  $$\begin{aligned} \frac{d}{d\delta }T_\delta =V\circ T_\delta \quad \text {for}\quad \delta \not =0,\quad T_0=\mathrm{Id}_x, \end{aligned}$$

  has a unique solution $T_{[\cdot ]}\in C^1((-\delta _0,\delta _0); W^{2,\infty }_\mathrm{loc}(\mathbb {R}^2)^2);$

• (b) the Cauchy problem

  $$\begin{aligned} \frac{d}{d\delta }T^{-1}_\delta =-V\circ T_\delta \quad \text {for}\quad \delta \not =0,\quad T^{-1}_0=\mathrm{Id}_y, \end{aligned}$$

  has a unique solution $T^{-1}_{[\cdot ]}\in C^1((-\delta _0,\delta _0); W^{2,\infty }_\mathrm{loc}(\mathbb {R}^2)^2);$

• (c) fixed $\delta \in (-\delta _0,\delta _0),$ the inverse mapping to $T_\delta$ is $T^{-1}_\delta ,$ that is,

  $$\begin{aligned} T_\delta \circ T^{-1}_\delta =\mathrm{Id}_y,\quad T^{-1}_\delta \circ T_\delta =\mathrm{Id}_x,\quad x,y\in \mathbb {R}^2. \end{aligned}$$

Proof

The proof is carried out in [21, Section 2] for $W^{1,\infty }$ vector fields. Without any changes, the arguments are also applicable to $W^{2,\infty }$ vector fields V. $\square$

Before proceeding further, we record here a consequence of [2, Theorem 2.2] that will be useful throughout the paper. For each $k=1,2,$ the space $W^{k,\infty }(\Omega )$ is algebraically and topologically isomorphic to $C^{k-1,1}(\overline{\Omega }).$

Let us fix $\delta \in (-\delta _0,\delta _0)$ and set $\omega _\delta =T_\delta (\omega _0)$ and $\gamma _\delta =T_\delta (\gamma _0).$ Since $T_\delta$ is an orientation-preserving $C^{1,1}$-diffeomorphism between $\Omega$ and $T_\delta (\Omega ),$ coinciding with the identity on $\partial \Omega$, we have $T_\delta (\Omega )=\Omega ,$ see [6, Theorem 5.5-2]. The domain $\omega _\delta \Subset \Omega$ is a simply connected one, and its boundary $\partial \omega _\delta =T_\delta (\partial \omega _0)$ is of class $C^{1,1}$, with the outward pointing unit normal $\nu ^\delta .$ Moreover, the boundary $\partial \omega _\delta$ is the union of two disjoint simple curves $\gamma _\delta$ and $\partial \omega _\delta \setminus \gamma _\delta$, and $\gamma _\delta$ is a relatively open set. So we can put $\Omega _\delta =\Omega \setminus \overline{\gamma }_\delta .$ Against this background, we conclude that the mapping $T_\delta$ induces a $C^{1,1}$-diffeomorphism between the non-Lipschitz domains $\Omega _0$ and $\Omega _\delta .$ To complete our considerations on the geometry of the perturbed domain $\Omega _\delta ,$ let us note that from the condition $\mathrm{supp}V\Subset \Omega$ follows the validity of G1 for $\Omega _\delta .$

As above, the equilibrium of the composite plate $\Omega ^\mathrm{3D}_\delta$ clamped along the outer edge, with elastic part $(\Omega \setminus \overline{\omega }_\delta )^\mathrm{3D}$, rigid part $\omega ^\mathrm{3D}_\delta$, interfacial crack $\gamma ^\mathrm{3D}_\delta$, and potential deformation energy

$$\begin{aligned} \mathcal {E}(\Omega _\delta ;W,w)&=\frac{1}{2}\int \limits _{\Omega \setminus \overline{\omega }_\delta } \sigma (W):\varepsilon (W)\,\mathrm{d}y\\&\quad -\frac{1}{2}\int \limits _{\Omega \setminus \overline{\omega }_\delta } m(w): \nabla ^2_x w\,\mathrm{d}y-\int \limits _{\Omega _\delta } F\cdot W\,\mathrm{d}y-\int \limits _{\Omega _\delta } fw\,\mathrm{d}y, \end{aligned}$$

is achieved when the displacement field of its midsurface $(W^\delta ,w^\delta )$ is a solution of the minimization problem

$$\begin{aligned} \inf \limits _{({W},{w})\in K^\delta (\Omega _\delta )}\mathcal {E}(\Omega _\delta ;{W},{w}), \end{aligned}$$

(20)

the set of kinematically admissible displacement fields being

The minimization problem (20) possesses a unique solution satisfying the variational inequality

$$\begin{aligned}&(W^\delta ,w^\delta )\in K^\delta (\Omega _\delta ),\quad \int \limits _{\Omega \setminus \overline{\omega }_\delta } \sigma (W^\delta ):\varepsilon (\overline{W}-W^\delta )\,\mathrm{d}y-\int \limits _{\Omega \setminus \overline{\omega }^\delta } m (w^\delta ): \nabla ^2_y (\overline{w}-w^\delta )\,\mathrm{d}y\nonumber \\&\quad -\int \limits _{\Omega _\delta } F\cdot (\overline{W}-W^\delta )\,\mathrm{d}y-\int \limits _{\Omega _\delta } f(\overline{w}-w^\delta ) \,\mathrm{d}y\ge 0\quad \forall (\overline{W},\overline{w})\in K^\delta (\Omega _\delta ).\nonumber \\ \end{aligned}$$

(21)

By definition, the first-order shape derivative of the potential deformation energy $\mathcal {E}$ at $\Omega _0$ in the direction of V,  if it exists, is given by the limit

$$\begin{aligned} \lim \limits _{\delta \rightarrow 0}\frac{\mathcal {E}(\Omega _\delta ;W^\delta ,w^\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }. \end{aligned}$$

(22)

The questions we are now interested in are whether (22) does exist and how this shape derivative can be expressed in an explicit form.

To state our main result, let us fix some notation. For $U,W\in H^1(\Omega _0)^2,$ we define the bilinear form

$$\begin{aligned} A_p(V; U,W)=\mathrm{div}_{x}V\sigma (U):\varepsilon (W)-\sigma (U):E\left( \nabla _x V; W\right) -\sigma (W):E\left( \nabla _x V; U\right) , \end{aligned}$$

(23)

where the second-order tensor E reads as

$$\begin{aligned} E\left( X; W\right) =1/2\left( \nabla _x W X+\left( \nabla _x W X\right) ^\mathrm{T}\right) ,\quad X\in \mathbb {R}^{2\times 2}. \end{aligned}$$

Simultaneously, for $u,w\in H^2(\Omega _0),$ we put

$$\begin{aligned} A_v(V; u,w)= & {} -\mathrm{div}_xV m(u):\nabla ^2_x w+m(u):\left( 2\nabla ^2_x w\nabla _xV+N\nabla _xw\right) \nonumber \\&\quad +m(w):\left( 2\nabla ^2_x u\nabla _xV+N\nabla _xu\right) , \end{aligned}$$

(24)

where

$$\begin{aligned} (N\nabla _xw)_{ij}=(\nabla ^2_xV_k)_{ij} (\nabla _x w)_k. \end{aligned}$$

We finally set

$$\begin{aligned} \mathfrak {G}(\Omega _0;V)= & {} \frac{1}{2}\int \limits _{\Omega \setminus \overline{\omega }_0}\left( A_p(V; W_0,W_0)+A_v(V; w_0,w_0)\right) \,\mathrm{d}x\nonumber \\&\quad -\int \limits _{\Omega \setminus \overline{\omega }_0}\left( \mathrm{div}_x(F\otimes V)\cdot W_0+\mathrm{div}_x(fV)w_0\right) \,\mathrm{d}x\nonumber \\&\quad -\left\langle \sigma ^+_{\nu _0}(W_0),\nabla _x V\llbracket W_0\rrbracket \cdot \nu _0\right\rangle ^{00}_{1/2,\gamma _0}\nonumber \\&\quad -\left\langle m^+_{\nu _0}(w_0),\llbracket \nabla _x w_0\rrbracket \cdot \left( \nabla _x V+\left( \nabla _x V\right) ^\mathrm{T}\right) \nu _0\right\rangle ^{00}_{1/2,\gamma _0}\nonumber \\&\quad -\left\langle \sigma ^+(W_0)\nu _0, B_0V\right\rangle _{1/2,\partial \omega _0}+\left\langle t^+_{\nu _0}(w_0),A_0\cdot V\right\rangle _{3/2,\partial \omega _0}\nonumber \\&\quad -\left\langle m^+_{\nu _0}(w_0),\nabla _x(A_0\cdot V)\cdot \nu _0\right\rangle _{1/2,\partial \omega _0}\nonumber \\&\quad -\int \limits _{\omega _0}\left( F\cdot B_0V+fA_0\cdot V\right) \,\mathrm{d}x. \end{aligned}$$

(25)

Theorem 2.2

Under the assumptions of Proposition 2.1, the limit in (22) exists and is equal to (25).

The remainder of the paper is devoted to the quite long and technical proof of Theorem 2.2. The key idea behind the proof consists in applying the mapping ${T}_\delta :\Omega _0\rightarrow \Omega _\delta$ to transform the potential deformation energy $\mathcal {E}(\Omega _\delta ;W^\delta ,w^\delta )$ to the reference domain $\Omega _0$ and calculating the desired limit as $\delta \rightarrow 0$ in $\Omega _0$. The particular steps are presented below.

We begin with the study of properties of the mapping ${T}_\delta :\Omega _0\rightarrow \Omega _\delta$. In view of [38, Chapter 2, Lemma 3.4], it generates an isomorphism between the spaces $H^1_{\partial \Omega }(\Omega _\delta )^2\times H^2_{\partial \Omega }(\Omega _\delta )$ and $H^1_{\partial \Omega }(\Omega _0)^2\times H^2_{\partial \Omega }(\Omega _0)$ via

$$\begin{aligned} \mathbf{T} _{\delta }:H^1_{\partial \Omega }(\Omega _\delta )^2\times H^2_{\partial \Omega }(\Omega _\delta )\rightarrow H^1_{\partial \Omega }(\Omega _0)^2\times H^2_{\partial \Omega }(\Omega _0):(W,w)\mapsto (W,w)\circ T_\delta . \end{aligned}$$

For $(U^\delta ,u^\delta )\in H^1(\Omega _\delta )^2\times H^2(\Omega _\delta ),$ we introduce the notation

$$\begin{aligned} (U_\delta ,u_\delta )(x)=(U^\delta ,u^\delta )(T_\delta (x)),\quad x\in \Omega _0. \end{aligned}$$

Then the derivatives of first and second orders of such transformed functions can be calculated from the relations

$$\begin{aligned}&\qquad \nabla _{x} U_\delta (x)=\nabla _{y} U^\delta (T_\delta (x))\nabla _{x} T_\delta (x),\\&\nabla ^2_{x} u_\delta (x)=(\nabla _{x} T_\delta (x))^\mathrm{T}\nabla ^2_{y} u^\delta (T_\delta (x))\nabla _{x} T_\delta (x)+(L^\delta \nabla _{y} u^\delta )(T_\delta (x)), \end{aligned}$$

where

$$\begin{aligned} (L^\delta \nabla _{y} u^\delta )_{ij}(T_\delta (x))=(\nabla ^2_{x} T_{\delta k}(x))_{ij}(\nabla _y u^\delta )_k(T_\delta (x)). \end{aligned}$$

On the other hand, for $(U_\delta ,u_\delta )\in H^1(\Omega _0)^2\times H^2(\Omega _0),$ with

$$\begin{aligned} (U^\delta ,u^\delta )(y)=(U_\delta ,u_\delta )(T^{-1}(y)),\quad y\in \Omega _\delta , \end{aligned}$$

it holds

$$\begin{aligned}&\nabla _y U^\delta (y)=\nabla _x U_\delta (T^{-1}_\delta (y))(\nabla _x T_\delta )^{-1}(T^{-1}_\delta (y)),\nabla ^2_y u^\delta (y)=(\nabla _x{T}_\delta )^{-\mathrm T}(T^{-1}_\delta (y))\\&\quad \left( \nabla ^2_xu_\delta (T^{-1}_\delta (y))-(M_\delta \nabla _xu_\delta )(T^{-1}_\delta (y))\right) (\nabla _x {T}_\delta )^{-1}(T^{-1}_\delta (y)), \end{aligned}$$

where

$$\begin{aligned} (M_\delta \nabla _xu_\delta )_{ij}(T^{-1}_\delta (y))=(\nabla ^2_{x} {T}_{\delta k})_{ij}(T^{-1}_\delta (y))(\nabla _x{T}_\delta )^{-1}_{lk}(T^{-1}_\delta (y))(\nabla _x u_{\delta })_l(T^{-1}(y)). \end{aligned}$$

Let $K_\delta (\Omega _0)$ be the image of the set of the admissible displacement fields $K^\delta (\Omega _\delta )$ under the mapping $\mathbf{T} _{\delta }$. In order to identify $K_\delta (\Omega _0),$ we first note that any element $(W,w)\in K_\delta (\Omega _0)$ satisfies the inequality constraint

(26)

The transformed outward pointing unit normal $\nu _\delta =\nu ^\delta \circ T_\delta$ is related to $\nu _0$ through the equation [6, Section 1.7]

$$\begin{aligned} \nu _\delta =\frac{(\nabla _{x} T_\delta )^{-\mathrm T}\nu _0}{\left\| (\nabla _{x} T_\delta )^{-\mathrm T}\nu _0\right\| _2}, \end{aligned}$$

where $\Vert \cdot \Vert _2$ stands for the Euclidean norm in $\mathbb {R}^2.$ Inequality (26) then reads

Additionally, the restriction (W, w) on the domain $\omega _0$ belongs to the Cartesian product $R^\delta _{p}(\omega _0)\times R^\delta _{v}(\omega _0)$ of the sets

$$\begin{aligned}&R^\delta _{p}(\omega _0)=\{\,W(x)=BT_\delta (x)+C:B\in \mathbb {R}^{2\times 2}, B=-B^\mathrm{T},\; C\in \mathbb {R}^{2}, \;x\in \omega _0\,\},\\&R^\delta _{v}(\omega _0)=\{\,w(x)=A\cdot T_\delta (x)+ a:A\in \mathbb {R}^2,\;a\in \mathbb {R},\;x\in \omega _0\,\}. \end{aligned}$$

Thus, $K_\delta (\Omega _0)$ may be defined as follows

and so $\mathbf{T} _{\delta }$ is a bijection between the sets $K^\delta (\Omega _\delta )$ and $K_\delta (\Omega _0).$

The change of variables $y=T_\delta (x)$ in the variational inequality (7) implies that $(W_\delta ,w_\delta )\in K_\delta (\Omega _0)$ is a unique solution of the variational inequality

$$\begin{aligned}&2h\int \limits _{\Omega \setminus \overline{\omega }_0}\mathrm{det}(\nabla _x T_\delta )\mathbb {C}E((\nabla _x T_\delta )^{-1};W_\delta ):E((\nabla _x T_\delta )^{-1};\overline{W}-W_\delta )\,\mathrm{d}x\nonumber \\&\qquad -\frac{2h^3}{3}\int \limits _{\Omega \setminus \overline{\omega }_0}\mathrm{det}(\nabla _x T_\delta )\mathbb {C}\left\{ (\nabla _x{T}_\delta )^{-\mathrm T}\left( \nabla ^2_x w_\delta -M_\delta \nabla _x w_\delta \right) (\nabla _x {T}_\delta )^{-1}\right\} \nonumber \\&\qquad :\left\{ (\nabla _x{T}_\delta )^{-\mathrm T}\left( \nabla ^2_x (\overline{w}-w_\delta )-M_\delta \nabla _x (\overline{w}-w_\delta )\right) (\nabla _x {T}_\delta )^{-1}\right\} \,\mathrm{d}x\nonumber \\&\quad \ge \int \limits _{\Omega _0}\mathrm{det}(\nabla _x T_\delta ) F_\delta \cdot (\overline{W}-W_\delta )\,\mathrm{d}x+\int \limits _{\Omega _0}\mathrm{det}(\nabla _x T_\delta ) f_\delta (\overline{w}-w_\delta )\,\mathrm{d}x\nonumber \\&\quad \forall (\overline{W},\overline{w})\in K_\delta (\Omega _0). \end{aligned}$$

(27)

Here, $(F_\delta ,f_\delta )(x)=(F,f)(T_\delta (x)).$

As a preliminary step towards the verification of the uniform boundedness of $(W_\delta ,w_\delta )$ in $H^1_{\partial \Omega }(\Omega _0)^2\times H^2_{\partial \Omega }(\Omega _0)$, we establish the next proposition.

Proposition 2.3

The following Taylor formulae as $\delta \rightarrow 0$ are valid:

$$\begin{aligned} T_\delta =\mathrm{Id}_x+\delta V+r^1_\delta ,\quad \Vert r^1_\delta \Vert _{W^{2,\infty }(\Omega )^2}=o(\delta ), \end{aligned}$$

(28)

$$\begin{aligned} \nabla _x T_\delta =I+\delta \nabla _xV+r^2_\delta ,\quad \Vert r^2_\delta \Vert _{W^{1,\infty }(\Omega )^{2\times 2}}=o(\delta ), \end{aligned}$$

(29)

$$\begin{aligned} \mathrm{det}(\nabla _x T_\delta )=1+\delta \mathrm{div}_xV+r^3_\delta ,\quad \Vert r^3_\delta \Vert _{W^{1,\infty }(\Omega )}=o(\delta ),\end{aligned}$$

(30)

$$\begin{aligned} (\nabla _{x} T_\delta )^{-1}=I-\delta \nabla _xV+r^4_\delta ,\quad \Vert r^4_\delta \Vert _{W^{1,\infty }(\Omega )^{2\times 2}}=o(\delta ), \end{aligned}$$

(31)

Proof

This is an immediate consequence of Proposition 2.1. $\square$

Lemma 2.4

There exists $\delta _1\in (0,\delta _0)$ such that the solution $(W_\delta ,w_\delta )$ of the transformed variational inequality (27) is uniformly bounded with respect to $\delta \in [-\delta _1,\delta _1]$ in $H^1_{\partial \Omega }(\Omega _0)^2\times H^2_{\partial \Omega }(\Omega _0)$ by a constant $c>0:$

$$\begin{aligned} \left\| (W_\delta ,w_\delta )\right\| _{H^1(\Omega _0)^2\times H^2(\Omega _0)}\le c. \end{aligned}$$

(32)

Proof

Since $(W^\delta ,w^\delta )|_{\omega _\delta }\in R_p(\omega _\delta )\times R_v(\omega _\delta ),$ we can replace the integration domain $\Omega \setminus \overline{\omega }_\delta$ by $\Omega _\delta$ in (21) and the integration domain $\Omega \setminus \overline{\omega }_0$ by $\Omega _0$ in (27). The modified variational inequality (27) is referred to as (27)$_{\Omega _0}.$ Let $(W,w), (U,u)\in H^1(\Omega _0)^2\times H^2(\Omega _0).$ Invoking Proposition 2.3, we deduce that the integrals on the left-hand side of (27)$_{\Omega _0}$ admit the following asymptotic expansions as $\delta \rightarrow 0:$

$$\begin{aligned}&2h\int \limits _{\Omega _0}\mathrm{det}(\nabla _x T_\delta )\mathbb {C}E((\nabla _x T_\delta )^{-1};U):E((\nabla _x T_\delta )^{-1};W)\,\mathrm{d}x\nonumber \\&\quad =\int \limits _{\Omega _0}\left( \sigma (U):\varepsilon (W)+\delta A_p(V;U,W)+R_1(\delta ,U,W)\right) \,\mathrm{d}x,\end{aligned}$$

(33)

$$\begin{aligned}&\frac{2h^3}{3}\int \limits _{\Omega _0}\mathrm{det}(\nabla _x T_\delta )\mathbb {C}\left\{ (\nabla _x{T}_\delta )^{-\mathrm T}\left( \nabla ^2_x u-M_\delta \nabla _x u\right) (\nabla _x {T}_\delta )^{-1}\right\} \nonumber \\&\qquad :\left\{ (\nabla _x{T}_\delta )^{-\mathrm T}\left( \nabla ^2_x w-M_\delta \nabla _x w\right) (\nabla _x {T}_\delta )^{-1}\right\} \,\mathrm{d}x\nonumber \\&\quad = \int \limits _{\Omega _0}\left( -m(u):\nabla ^2_xw+\delta A_v(V;u,w)+R_2(\delta ,u,w)\right) \,\mathrm{d}x, \end{aligned}$$

(34)

where

$$\begin{aligned}&\Vert R_1(\delta ,U,W)\Vert _{L^1(\Omega _0)}\le c_1(|\delta |)\Vert U\Vert _{H^1(\Omega _0)^2}\Vert W\Vert _{H^1(\Omega _0)^2},\quad 0\le c_1(|\delta |)=o(\delta ),\nonumber \\ \end{aligned}$$

(35)

$$\begin{aligned}&\Vert R_2(\delta ,u,w)\Vert _{L^1(\Omega _0)}\le c_2(|\delta |)\Vert u\Vert _{H^2(\Omega _0)}\Vert w\Vert _{H^2(\Omega _0)},\quad 0\le c_2(|\delta |)=o(\delta ). \end{aligned}$$

(36)

The terms in the right-hand side of (27)$_{\Omega _0}$ can be rewritten in like manner:

$$\begin{aligned}&\int \limits _{\Omega _0}\mathrm{det}(\nabla _x T_\delta )F_\delta \cdot W\,\mathrm{d}x=\int \limits _{\Omega _0}\left( F\cdot W+\delta \mathrm{div}_x(F\otimes V)\cdot W+R_3(\delta ,W)\right) \,\mathrm{d}x, \end{aligned}$$

(37)

$$\begin{aligned}&\int \limits _{\Omega _0}\mathrm{det}(\nabla _x T_\delta )f_\delta w\,\mathrm{d}x=\int \limits _{\Omega _0}\left( fw+\delta \mathrm{div}_x(fV)w+R_4(\delta ,w)\right) \,\mathrm{d}x, \end{aligned}$$

(38)

where

$$\begin{aligned}&\Vert R_3(\delta ,W)\Vert _{L^1(\Omega _0)}\le c_3(|\delta |)\Vert W\Vert _{L^2(\Omega _0)^2},\quad 0\le c_3(|\delta |)=o(\delta ), \end{aligned}$$

(39)

$$\begin{aligned}&\Vert R_4(\delta ,w)\Vert _{L^1(\Omega _0)}\le c_4(|\delta |)\Vert w\Vert _{L^2(\Omega _0)},\quad 0\le c_4(|\delta |)=o(\delta ). \end{aligned}$$

(40)

Moreover, since $\nabla _x V\in C^{0,1}(\overline{\Omega })^{2\times 2}$ and $(F,f)\in C^1(\overline{\Omega })^2\times C^1(\overline{\Omega }),$ we arrive at the estimates

$$\begin{aligned}&\left\| A_p(V;U,W) \right\| _{L^1(\Omega _0)}\le c\Vert U\Vert _{H^1(\Omega _0)^2}\Vert W\Vert _{H^1(\Omega _0)^2},\nonumber \\&\quad \left\| A_v(V;u,w) \right\| _{L^1(\Omega _0)}\le c\Vert u\Vert _{H^2(\Omega _0)}\Vert w\Vert _{H^2(\Omega _0)}, \end{aligned}$$

(41)

$$\begin{aligned}&\left\| \mathrm{div}_x(F\otimes V)\cdot W\right\| _{L^1(\Omega _0)}\le c \Vert W\Vert _{L^2(\Omega _0)^2},\nonumber \\&\quad \left\| \mathrm{div}_x(fV)w\right\| _{L^1(\Omega _0)}\le c\Vert w\Vert _{L^2(\Omega _0)}. \end{aligned}$$

(42)

We next take $\overline{W}=(0,0)^\mathrm{T}$ and $\overline{w}=0$ into the variational inequality (27)$_{\Omega _0}$ and use the asymptotic expansions (33)–(34) and (37)–(38) to discover

$$\begin{aligned}&\int \limits _{\Omega _0}\sigma (W_\delta ):\varepsilon (W_\delta )\,\mathrm{d}x-\int \limits _{\Omega _0}m(w_\delta ):\nabla ^2_xw_\delta \,\mathrm{d}x\nonumber \\&\quad \le -\delta \int \limits _{\Omega _0}\left( A_p(V;W_\delta ,W_\delta )+A_v(V;w_\delta ,w_\delta )\right) \,\mathrm{d}x\nonumber \\&\qquad +\int \limits _{\Omega _0}\left( F\cdot W_\delta +fw_\delta +\delta \left( \mathrm{div}_x(F\otimes V)\cdot W_\delta +\mathrm{div}_x(fV)w_\delta \right) \right. -R_1(\delta ,W_\delta ,W_\delta )\nonumber \\&\qquad -R_2(\delta ,w_\delta ,w_\delta )\left. +R_3(\delta ,W_\delta )+R_4(\delta ,w_\delta )\right) \,\mathrm{d}x. \end{aligned}$$

(43)

Applying the Korn inequality and the Friedrichs inequality on the left-hand side of (43), the Cauchy inequality and estimates (35)–(36), (39)–(40), and (41)–(42) on the right-hand side of (43) yields that there exists $\delta _1\in (0,\delta _0)$ such that for every $\delta \in [-\delta _1,\delta _1]$ the uniform estimate (32) takes place. This proofs the lemma. $\square$

For convenience purposes, we collect here auxiliary statements that will be employed repeatedly in the sequel.

Proposition 2.5

[19, Lemma 1.14] Let $u\in H^{1/2}(\partial \omega _0)$ and $v\in C^{0,1}(\partial \omega _0).$ Then $uv\in H^{1/2}(\partial \omega _0)$ and

$$\begin{aligned} \Vert uv\Vert _{H^{1/2}(\partial \omega _0)}\le c\Vert u\Vert _{H^{1/2}(\partial \omega _0)}\Vert v\Vert _{C^{0,1}(\partial \omega _0)} \end{aligned}$$

with $c>0$ independent of u and v.

We label the part $\partial \Omega \cup \partial \omega _0\setminus {\gamma }_0$ of the outer boundary of the domain $\Omega _0\setminus \overline{\omega }_0$ as $\Sigma$. For each $k=1,2,$ we introduce the extension operator $\mathrm{Ext}_k: H^k_{\Sigma }(\Omega \setminus \overline{\omega }_0)^{k_*}\rightarrow H^k_{\partial \Omega }(\Omega _0)^{k_*}$ given by

$$\begin{aligned} v\mapsto {\left\{ \begin{array}{ll} (\underbrace{0,\ldots ,0}_{k_*})^\mathrm{T}&{}\quad \text {in}\quad \omega _0,\\ v&{}\quad \text {in}\quad \Omega \setminus \overline{\omega }_0. \end{array}\right. } \end{aligned}$$

It is clearly that $\mathrm{Ext}_k$ are well defined, linear, and continuous operators.

With these preliminaries behind us, we are able to state the key assertion of our analysis that guarantees the existence of certain test elements.

Theorem 2.6

Let $(W_0,w_0)\in K_0(\Omega _0)$ be the solution of the variational inequality (7), and let $(W_\delta ,w_\delta )\in K_\delta (\Omega _0)$ be the solution of the transformed variational inequality (27). Then there exists $\delta _2\in (0,\delta _1]$ and tuples $({W}^i_\delta ,{w}^i_\delta )$ such that for every $\delta \in [-\delta _2,\delta _2]$ the inclusions

$$\begin{aligned}&(\widehat{W}^1_\delta , \widehat{w}^1_\delta )=(W_0,w_0)+\delta ({W}^1_\delta ,{w}^1_\delta )\in K_\delta (\Omega _0),\\&(\widehat{W}^2_\delta , \widehat{w}^2_\delta )=(W_\delta ,w_\delta )-\delta ({W}^2_\delta , {w}^2_\delta )\in K_0(\Omega _0), \end{aligned}$$

and uniform in $\delta$ bounds

$$\begin{aligned} \Vert ({W}^i_\delta , {w}^i_\delta )\Vert _{H^1(\Omega _0)^2\times H^2(\Omega _0)}\le c,\quad i=1,2, \end{aligned}$$

(44)

are valid.

Proof

For each $i=1,2,$ we construct the required tuple $(W^i_\delta ,w^i_\delta )$ as sums $W^i_\delta =P^i_\delta +Q^i_\delta$ and $w^i_\delta =p^i_\delta +q^i_\delta$, where the set $\{P^i_\delta , p^i_\delta \}$ provides the appropriate inequality constraint for $(\widehat{W}^i_\delta , \widehat{w}^i_\delta )$ on $\gamma _0$, and the set $\{Q^i_\delta , q^i_\delta \}$ is responsible for the correct structure of $(\widehat{W}^i_\delta , \widehat{w}^i_\delta )$ in $\omega _0,$ which is specified by the equality constraint. The case $\delta =0$ is trivial, so we assume that $\delta \ne 0.$

We start with the construction of $({W}^1_\delta , {w}^1_\delta )$. Since $\llbracket W_0\rrbracket \in H^{1/2}(\partial \omega _0)^2$ and by Proposition 2.5, the function

$$\begin{aligned} {P}^1_\delta =\mathrm{Ext}_1\circ \mathrm{Tr}^{-1}_{1,\Omega \setminus \overline{\omega }_0}\left( \frac{\nabla _{x} T_\delta -I}{\delta }\llbracket W_0\rrbracket ; (0,0)^\mathrm{T}\right) \end{aligned}$$

(45)

belongs to $H^1_{\partial \Omega }(\Omega _0)^2.$ We define $w^1_\delta$ in the following way. Observe first that $Z_\delta =(\nabla _{x} T_\delta )^{-1}(\nabla _{x} T_\delta )^{-\mathrm T}\in C^{0,1}(\overline{\Omega })^{2\times 2}.$ For any function $v\in H^2(\Omega \setminus \overline{\omega }_0),$ we write

$$\begin{aligned} ((\nabla _xv)^+\cdot \nu _0)\nu _0\cdot Z_\delta \nu _0=(\nabla _x v)^{+}\cdot Z_\delta \nu _0+(\nabla _x v)^{+}\cdot (\nu _0\otimes \nu _0Z_\delta -Z_\delta )\nu _0\quad \text {on}\quad \partial \omega _0. \end{aligned}$$

There exists $\delta _2\in (0,\delta _1]$ and $\alpha >0$ such that $\alpha _\delta (s)=(\nu _0\cdot Z_\delta \nu _0)(s)\ge \alpha$ for every $(\delta ,s)\in [-\delta _2,\delta _2]\times \partial \omega _0$. For $Z_{\delta \tau _0}=(\nu _0\otimes \nu _0Z_\delta -Z_\delta )\nu _0,$ we compute $\nu _0\cdot Z_{\delta {\tau _0}}=0$ on $\partial \omega _0,$ which implies $Z_{\delta {\tau _0}}=\beta _\delta \tau _0$ with $\beta _\delta =\tau _0\cdot Z_{\delta {\tau _0}}.$ By construction, $\beta _\delta \in C^{0,1}(\partial \omega _0).$ With the help of standard arguments in local coordinates (see also [11, Proposition 3.1]), we infer

$$\begin{aligned} (\nabla _x v)^+=(\partial _{\tau _0}v^+)\tau _0+((\nabla _x v)^+\cdot \nu _0)\nu _0\quad \text {on}\quad \partial \omega _0, \end{aligned}$$

and so

$$\begin{aligned} \alpha _\delta (\nabla _x v)^+\cdot \nu _0=(\nabla _x v)^{+}\cdot Z_\delta \nu _0+\beta _\delta \partial _{\tau _0}v^+\quad \text {on}\quad \partial \omega _0. \end{aligned}$$

Exploiting the inclusions $\llbracket w_0\rrbracket \in H^{3/2}(\partial \omega _0)$ and $\llbracket \nabla _x w_0\rrbracket \in H^{1/2}(\partial \omega _0)^2$ and once more Proposition 2.5, we conclude that the function

$$\begin{aligned} p^{1}_\delta =\mathrm{Ext}_2\circ \mathrm{Tr}^{-1}_{2,\Omega \setminus \overline{\omega }_0}\left( \llbracket w_0\rrbracket ; 0; \frac{1}{\alpha _\delta }\left( \llbracket \nabla _x w_0\rrbracket \cdot \frac{I-Z_\delta }{\delta }\nu _0+{\beta _\delta }\partial _{\tau _0}\llbracket w_0\rrbracket \right) ; 0\right) \end{aligned}$$

is from the space $H^2_{\partial \Omega }(\Omega _0).$ Under the circumstances, a straightforward calculation reveals

$$\begin{aligned} \llbracket \nabla _x p^1_\delta \rrbracket \cdot Z_\delta \nu _0&=\left( \nabla _x p^{1}_\delta \right) ^+\cdot Z_\delta \nu _0=\alpha _\delta \left( \nabla _x p^{1}_\delta \right) ^+\cdot \nu _0-\beta _\delta \partial _{\tau _0} p^{1+}_\delta \nonumber \\&=\llbracket \nabla _x w_0\rrbracket \cdot \frac{I-Z_\delta }{\delta }\nu _0\quad \text {on}\quad {\partial \omega _0}. \end{aligned}$$

(46)

Let us set

$$\begin{aligned} Q^1_\delta (x)=B_0\frac{T_\delta -\mathrm{Id}_x}{\delta }(x),\quad q^1_\delta (x)=A_0\cdot \frac{T_\delta -\mathrm{Id}_x}{\delta }(x),\quad x\in \Omega . \end{aligned}$$

The functions $Q^1_\delta$ and $q^1_\delta$ vanish in a neighborhood of $\partial \Omega .$ We now check that $(\widehat{W}^1_\delta , \widehat{w}^1_\delta )\in K_\delta (\Omega _0).$ In fact, the inequality constraint follows from the following chain of equalities and inequalities

Moreover, we have no difficulty in verifying the equality constraint,

$$\begin{aligned} \widehat{W}^1_\delta (x)&=W_0(x)+\delta Q^1_\delta (x)=B_0\mathrm{Id}_x(x)+C_0+B_0(T_\delta (x)-\mathrm{Id}_x(x))\\&=B_0T_\delta (x)+C_0\quad \text {in}\quad \omega _0,\\ \widehat{w}^1_\delta (x)&=w_0(x)+\delta q^1_\delta (x)=A_0\cdot \mathrm{Id}_x(x)+a_0+A_0\cdot (T_\delta (x)-\mathrm{Id}_x(x))\\&= A_0\cdot T_\delta (x)+a_0\quad \text {in}\quad \omega _0. \end{aligned}$$

We next turn our attention to $(W^2_\delta ,w^2_\delta ).$ The construction of $W^2_\delta$ is not much harder than that for $W^1_\delta .$ Indeed, with $\llbracket W_\delta \rrbracket \in H^{1/2}(\partial \omega _0)$ and Proposition 2.5 in mind, the function

$$\begin{aligned} {P}^2_\delta =\mathrm{Ext}_1\circ \mathrm{Tr}^{-1}_{1,\Omega \setminus \overline{\omega }_0}\left( \frac{I-(\nabla _{x} T_\delta )^{-1}}{\delta }\llbracket W_\delta \rrbracket ; (0,0)^\mathrm{T}\right) \end{aligned}$$

(47)

is an element of the space $H^1_{\partial \Omega }(\Omega _0)^2.$ To define $w^2_\delta ,$ we appeal to slightly simpler arguments than for $w^1_\delta$. Taking into account the inclusions $\llbracket w_\delta \rrbracket \in H^{3/2}(\partial \omega _0)$ and $\llbracket \nabla _x w_\delta \rrbracket \in H^{1/2}(\partial \omega _0)^2$ and applying Proposition 2.5 again, we derive that the function

$$\begin{aligned} p^{2}_\delta =\mathrm{Ext}_2\circ \mathrm{Tr}^{-1}_{2,\Omega \setminus \overline{\omega }_0}\left( \llbracket w_\delta \rrbracket ; 0; \llbracket \nabla _x w_\delta \rrbracket \cdot \frac{I-Z_\delta }{\delta }\nu _0; 0\right) \end{aligned}$$

(48)

belongs to $H^{2}_{\partial \Omega }(\Omega _0).$ We finally put

$$\begin{aligned} Q^2_\delta (x)=B_\delta \frac{T_\delta -\mathrm{Id}_x}{\delta }(x),\quad q^2_\delta =A_\delta \cdot \frac{T_\delta -\mathrm{Id}_x}{\delta }(x),\quad x\in \Omega . \end{aligned}$$

The functions $Q^2_\delta$ and $q^2_\delta$ vanish in a neighborhood of $\partial \Omega .$ In view of the above discussion, we obtain

Furthermore,

$$\begin{aligned} \widehat{W}^2_\delta (x)&=W_\delta (x)-\delta Q^2_\delta (x)=B_\delta T_\delta (x)+C_\delta -B_\delta (T_\delta (x)-\mathrm{Id}_x(x))\\&=B_\delta \mathrm{Id}_x(x)+C_\delta \quad \text {in}\quad \omega _0,\\ \widehat{w}^2_\delta (x)&=w_\delta (x)-\delta q^2_\delta (x)=A_\delta \cdot T_\delta (x)+a_\delta -A_\delta \cdot (T_\delta (x)-\mathrm{Id}_x(x))\\&= A_\delta \cdot \mathrm{Id}_x(x)+a_\delta \quad \text {in}\quad \omega _0. \end{aligned}$$

We thus can conclude that $(\widehat{W}^2_\delta , \widehat{w}^2_\delta )\in K_0(\Omega _0),$ as asserted.

In order to complete the proof, we observe that the uniform estimates (44) is an immediate consequence of (32) and continuity of the lifting and extension operators. $\square$

These considerations allow us to establish convergence properties of the transformed solutions $(W_\delta ,w_\delta )$ as $\delta \rightarrow 0.$

Lemma 2.7

There exists a constant $c>0$ that is independent of $\delta$ such that it holds for every $\delta \in [-\delta _2,\delta _2]:$

$$\begin{aligned} \Vert (W_\delta ,w_\delta )-(W_0,w_0)\Vert _{H^1(\Omega _0)^2\times H^2(\Omega _0)}\le c\sqrt{|\delta |}. \end{aligned}$$

Proof

As in the proof of Lemma 2.4, we change the integration domain $\Omega \setminus \overline{\omega }_0$ by $\Omega _0$ in the variational inequality (7). The resulting variational inequality is referred to as (7)$_{\Omega _0}.$ Inserting $(\widehat{W}^1_\delta , \widehat{w}^1_\delta )$ and $(\widehat{W}^2_\delta , \widehat{w}^2_\delta )$ as test functions for (27)$_{\Omega _0}$ and (7)$_{\Omega _0},$ respectively, and adding both variational inequalities, we invoke the asymptotic expansions (33)–(34) and (37)–(38) to obtain

$$\begin{aligned}&\int \limits _{\Omega _0}\sigma (W_\delta -W_0):\varepsilon (W_\delta -W_0)\,\mathrm{d}x-\int \limits _{\Omega _0}m(w_\delta -w_0):\nabla ^2_x(w_\delta -w_0)\,\mathrm{d}x\\&\quad \le \delta \int \limits _{\Omega _0}\left( A_p(V;W_\delta ,W_0+\delta W^1_\delta -W_\delta )+A_v(V;w_\delta ,w_0+\delta w^1_\delta -w_\delta )\right) \,\mathrm{d}x\\&\qquad + \,\delta \int \limits _{\Omega _0}\left( F\cdot (W^2_\delta -W^1_\delta )+f(w^2_\delta -w^1_\delta )-\mathrm{div}_x(F\otimes V)\cdot (W_0+\delta W^1_\delta -W_\delta )\right. \\&\qquad \left. - \,\mathrm{div}_x(fV)(w_0+\delta w^1_\delta -w_\delta )\right) \,\mathrm{d}x\\&\qquad + \,\delta \int \limits _{\Omega _0}\left( \sigma (W_\delta ):\varepsilon (W^1_\delta )-m(w_\delta ):\nabla ^2_x w^1_\delta \right) \,\mathrm{d}x\\&\qquad - \,\delta \int \limits _{\Omega _0}\left( \sigma (W_0):\varepsilon (W^2_\delta )-m(w_0):\nabla ^2_x w^2_\delta \right) \,\mathrm{d}x\\&\qquad +\,\int \limits _{\Omega _0}\left( R_1(\delta ,W_\delta ,W_0+\delta W^1_\delta -W_\delta )+R_2(\delta ,w_\delta ,w_0+\delta w^1_\delta -w_\delta )\right) \,\mathrm{d}x\\&\qquad +\,\int \limits _{\Omega _0}\left( R_3(\delta ,W_0+\delta W^1_\delta -W_\delta )+R_4(\delta ,w_\delta ,w_0+\delta w^1_\delta -w_\delta )\right) \,\mathrm{d}x. \end{aligned}$$

By virtue of estimates (32),(35)–(36),(39)–(42), and (44), it holds for every $\delta \in [-\delta _2,\delta _2]:$

$$\begin{aligned} \int \limits _{\Omega _0}\sigma (W_\delta -W_0):\varepsilon (W_\delta -W_0)\,\mathrm{d}x-\int \limits _{\Omega _0}m(w_\delta -w_0):\nabla ^2_x(w_\delta -w_0)\,\mathrm{d}x\le c|\delta |. \end{aligned}$$

We finally apply the Korn inequality and the Friedrichs inequality to complete the proof. $\square$

Corollary 2.8

The sequence $(B_\delta , C_\delta , A_\delta , a_\delta )$ converges to $(B_0,C_0,A_0,a_0)$ in $\mathbb {R}^{2\times 2}\times \mathbb {R}^2\times \mathbb {R}^2\times \mathbb {R}$ as $\delta \rightarrow 0.$

Proof

By Lemma 2.7, it follows that $(W_\delta ,w_\delta )\rightarrow (W_0,w_0)$ strongly in $H^1(\omega _0)^2\times H^2(\omega _0),$ in other words,

$$\begin{aligned}&B_\delta \mathrm{Id}_x+C_\delta +\delta B_\delta \left( V+\frac{r^1_\delta }{\delta }\right) \rightarrow B_0\mathrm{Id}_x+C_0\quad \text {strongly in}\quad H^1(\omega _0)^2, \end{aligned}$$

(49)

$$\begin{aligned}&A_\delta \cdot \mathrm{Id}_x+a_\delta +\delta A_\delta \cdot \left( V+\frac{r^1_\delta }{\delta }\right) \rightarrow A_0\cdot \mathrm{Id}_x+a_0\quad \text {strongly in}\quad H^2(\omega _0).\nonumber \\ \end{aligned}$$

(50)

Estimate (32) forces that the sequences $B_\delta$ and $A_\delta$ are uniformly bounded with respect to small enough $\delta$. Thus we can rewrite (49)–(50) as

$$\begin{aligned}&B_\delta \mathrm{Id}_x+C_\delta \rightarrow B_0\mathrm{Id}_x+C_0\quad \text {strongly in}\quad H^1(\omega _0)^2\nonumber \\&\quad \text {and (up to a subsequence) a.e. in}\quad \omega _0, \end{aligned}$$

(51)

$$\begin{aligned}&A_\delta \cdot \mathrm{Id}_x+a_\delta \rightarrow A_0\cdot \mathrm{Id}_x+a_0\quad \text {strongly in}\quad H^2(\omega _0)\nonumber \\&\quad \text {and (up to a subsequence) a.e. in}\quad \omega _0. \end{aligned}$$

(52)

The strong convergences $\nabla _x(B_\delta \mathrm{Id}_x+C_\delta )\rightarrow \nabla _x(B_0\mathrm{Id}_x+C_0)$ in $L^2(\omega _0)^{2\times 2}$ and $\nabla _x(A_\delta \cdot \mathrm{Id}_x+a_\delta )\rightarrow \nabla _x(A_0\cdot \mathrm{Id}_x+a_0)$ in $H^1(\omega _0)^2$ mean in particular $B_\delta \rightarrow B_0$ in $\mathbb {R}^{2\times 2}$ and $A_\delta \rightarrow A_0$ in $\mathbb {R}^2.$ Hence (51)–(52) imply the remaining convergences $C_\delta \rightarrow C_0$ in $\mathbb {R}^{2}$ and $a_\delta \rightarrow a_0$ in $\mathbb {R},$ so we are done. $\square$

We intend now to identify the limit of the sequence $(W^i_\delta ,w^i_\delta )$ for each $i=1,2.$

Lemma 2.9

Let $(\widetilde{W}_0, \widetilde{w}_0)=(P_0+Q_0,p_0+q_0),$ where $P_0=\mathrm{Ext}_1\circ \mathrm{Tr}^{-1}_{1,\Omega \setminus \overline{\omega }_0}(\nabla _xV \llbracket W_0\rrbracket ; (0,0)^\mathrm{T}),$ $Q_0=B_0V,$ $p_0=\mathrm{Ext}_2\circ \mathrm{Tr}^{-1}_{2,\Omega \setminus \overline{\omega }_0}\left( \llbracket w_0\rrbracket ; 0;\llbracket \nabla _x w_0\rrbracket \cdot \left( \nabla _x V+\left( \nabla _x V\right) ^\mathrm{T}\right) \nu _0; 0\right) ,$ and $q_0=A_0\cdot V$. Then, for each $i=1,2,$ the sequence $({W}^i_\delta , {w}^i_\delta )$ converges strongly to $(\widetilde{W}_0,\widetilde{w}_0)$ in $H^1_{\partial \Omega }(\Omega _0)^2\times H^2_{\partial \Omega }(\Omega _0)$ as $\delta \rightarrow 0.$

Proof

We first note that, by construction, $(\widetilde{W}_0, \widetilde{w}_0)\in H^1_{\partial \Omega }(\Omega _0)^2\times H^2_{\partial \Omega }(\Omega _0);$ furthermore, $P_0\equiv (0,0)^\mathrm{T}$ and $p_0\equiv 0$ in $\omega _0.$ Since all the functions $P^i_\delta$ and $p^i_\delta$ vanish in $\omega _0,$ it is enough to investigate their convergence in $\Omega \setminus \overline{\omega }_0.$

Thanks to the asymptotic expansion (2930), continuity of the lifting operator $\mathrm{Tr}^{-1}_{1,\Omega \setminus \overline{\omega }_0},$ and Proposition 2.5, we have $P^{1}_\delta \rightarrow P_0$ strongly in $H^1(\Omega \setminus \overline{\omega }_0)^2.$ From the asymptotic expansion (28), it follows that $Q^1_\delta \rightarrow Q_0$ strongly in $H^1(\Omega _0)^2,$ and $q^1_\delta \rightarrow q_0$ strongly in $H^2(\Omega _0).$ The asymptotic expansion

$$\begin{aligned} Z_\delta =I-\delta \left( \nabla _xV+\left( \nabla _xV\right) ^\mathrm{T}\right) +{r^5_\delta },\quad \Vert r^5_\delta \Vert _{W^{1,\infty }(\Omega )^{2\times 2}}=o(\delta ), \end{aligned}$$

(53)

strong convergences $\alpha _\delta \rightarrow 1$ and $\beta _\delta \rightarrow 0$ in $C^{0,1}(\partial \omega _0),$ and Proposition 2.5 imply that $(\nabla _x p^{1}_\delta )^+\cdot \nu _0\rightarrow (\nabla _x p_0)^{+}\cdot \nu _0$ strongly in $H^{1/2}(\partial \omega _0).$ In view of $p^{1+}_\delta \equiv p^{+}_0$ on $\partial \omega _0$ and the continuity of the lifting operator $\mathrm{Tr}^{-1}_{2,\Omega \setminus \overline{\omega }_0},$ we conclude that $p^{1}_\delta \rightarrow p_0$ strongly in $H^2(\Omega \setminus \overline{\omega }_0).$ Consequently, $({W}^1_\delta , {w}^1_\delta )$ converges strongly to $(\widetilde{W}_0, \widetilde{w}_0)$ in $H^1_{\partial \Omega }(\Omega _0)^2\times H^2_{\partial \Omega }(\Omega _0),$ as claimed.

By Lemma 2.7 and the continuity of the trace operators $\mathrm{Tr}_{1,\Omega \setminus \overline{\omega }_0}$ and $\mathrm{Tr}_{1,{\omega }_0}$, we find $\llbracket W_\delta \rrbracket \rightarrow \llbracket W_0\rrbracket$ strongly in $H^{1/2}(\partial \omega _0)^2.$ Therefore, using the asymptotic expansion (31), continuity of the lifting operator $\mathrm{Tr}^{-1}_{1,\Omega \setminus \overline{\omega }_0},$ and Proposition 2.5, we deduce that $P^{2}_\delta \rightarrow P_0$ strongly in $H^1(\Omega \setminus \overline{\omega }_0)^2.$ With the aid of the asymptotic expansion (53) and continuity of the trace operators $\mathrm{Tr}_{2,\Omega \setminus \overline{\omega }_0}$ and $\mathrm{Tr}_{2,{\omega }_0},$ we derive $p^{2+}_\delta \rightarrow p^{+}_0$ strongly in $H^{3/2}(\partial \omega _0)$ and $(\nabla _x p^{2}_\delta )^+\cdot \nu _0\rightarrow (\nabla _x p_0)^+\cdot \nu _0$ strongly in $H^{1/2}(\partial \omega _0).$ Employing the continuity of the lifting operator $\mathrm{Tr}^{-1}_{2,\Omega \setminus \overline{\omega }_0},$ we infer that $p^{2}_\delta \rightarrow p_0$ strongly in $H^2(\Omega \setminus \overline{\omega }_0).$ Additionally, Corollary 2.8 yields $Q^2_\delta \rightarrow Q_0$ strongly in $H^1(\Omega _0)^2$ and $q^2_\delta \rightarrow q_0$ strongly in $H^2(\Omega _0).$ Thus, the strong convergence $({W}^2_\delta , {w}^2_\delta )$ to $(\widetilde{W}_0, \widetilde{w}_0)$ in $H^1_{\partial \Omega }(\Omega _0)^2\times H^2_{\partial \Omega }(\Omega _0)$ is established, and the lemma is entirely proved. $\square$

Now we have enough tools at hand to prove Theorem 2.2.

Proof of Theorem 2.2

For the calculation of limit (22), we first transform the potential deformation energy $\mathcal {E}(\Omega _\delta ;\cdot ,\cdot )$ to the reference configuration $\Omega _0,$ which leads to

$$\begin{aligned} \mathcal {E}_\delta (\Omega _0;W,w)= & {} h\int \limits _{\Omega \setminus {\overline{\omega }_0}}\mathrm{det}(\nabla _x T_\delta )\mathbb {C}E((\nabla _x T_\delta )^{-1};W):E((\nabla _x T_\delta )^{-1};W)\,\mathrm{d}x\\&+\frac{h^3}{3}\int \limits _{\Omega \setminus {\overline{\omega }_0}}\mathrm{det}(\nabla _x T_\delta )\mathbb {C}\left\{ (\nabla _x{T}_\delta )^{-\mathrm T}\left( \nabla ^2_x w-M_\delta \nabla _x w\right) (\nabla _x {T}_\delta )^{-1}\right\} \\&:\left\{ (\nabla _x{T}_\delta )^{-\mathrm T}\left( \nabla ^2_x w-M_\delta \nabla _x w\right) (\nabla _x {T}_\delta )^{-1}\right\} \,\mathrm{d}x\\&-\int \limits _{\Omega _0}\mathrm{det}(\nabla _x T_\delta )F_\delta \cdot W\,\mathrm{d}x-\int \limits _{\Omega _0}\mathrm{det}(\nabla _x T_\delta )f_\delta w\,\mathrm{d}x, \end{aligned}$$

(W, w) belonging to $H^1(\Omega _0)^2\times H^2(\Omega _0).$ Utilizing the same machinery as in the proof of Lemma 2.4 yields

$$\begin{aligned} \mathcal {E}_\delta (\Omega _0;W,w)= & {} \mathcal {E}(\Omega _0;W,w)+\frac{\delta }{2}\int \limits _{\Omega \setminus \overline{\omega }_0}(A_p(V;W,W)+A_v(V;w,w))\,\mathrm{d}x\nonumber \\&-\delta \int \limits _{\Omega _0}\left( \mathrm{div}_x(F\otimes V)\cdot W+\mathrm{div}_x(fV)w\right) \,\mathrm{d}x+\int \limits _{\Omega _0}R_5(\delta ,W,w)\,\mathrm{d}x,\nonumber \\ \end{aligned}$$

(54)

where

$$\begin{aligned}&\Vert R_5(\delta ,W,w)\Vert _{L^1(\Omega _0)}\le c_5(|\delta |)\left( \Vert (W,w)\Vert ^2_{H^1(\Omega _0)^2\times H^2(\Omega _0)}+\Vert (W,w)\Vert _{L^2(\Omega _0)^3}\right) ,\\&\quad \quad 0\le c_5(|\delta |)=o(\delta ). \end{aligned}$$

Since $\mathbf{T} _{\delta }$ is a bijection between the sets $K^\delta (\Omega ^\delta )$ and $K_\delta (\Omega _0),$ it follows that

$$\begin{aligned} \mathcal {E}_\delta (\Omega _0;W_\delta ,w_\delta )&=\mathcal {E}(\Omega _\delta ;W^\delta ,w^\delta )=\inf \limits _{(W,w)\in K^\delta (\Omega _\delta )}\mathcal {E}(\Omega _\delta ;W,w)\nonumber \\&=\inf \limits _{(W,w)\in K_\delta (\Omega _0)}\mathcal {E}_\delta (\Omega _0;W,w). \end{aligned}$$

(55)

Let $\delta >0$. Using $(\widehat{W}^1_\delta , \widehat{w}^1_\delta )\in K_\delta (\Omega _0)$ as a competitor in (55), we discover

$$\begin{aligned} \frac{\mathcal {E}(\Omega _\delta ;W^\delta ,w^\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }&=\frac{\mathcal {E}_\delta (\Omega _0;W_\delta ,w_\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }\nonumber \\&\le \frac{\mathcal {E}_\delta (\Omega _0; \widehat{W}^1_\delta , \widehat{w}^1_\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }, \end{aligned}$$

(56)

which implies

$$\begin{aligned} \limsup \limits _{\delta \searrow 0}\frac{\mathcal {E}(\Omega _\delta ;W^\delta ,w^\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }\le \limsup \limits _{\delta \searrow 0}\frac{\mathcal {E}_\delta (\Omega _0; \widehat{W}^1_\delta , \widehat{w}^1_\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }. \end{aligned}$$

(57)

On the other hand, considering $(\widehat{W}^2_\delta , \widehat{w}^2_\delta )\in K_0(\Omega _0)$ as a competitor in (6), we arrive at

$$\begin{aligned} \frac{\mathcal {E}_\delta (\Omega _0;W_\delta ,w_\delta )-\mathcal {E}(\Omega _0; \widehat{W}^2_\delta , \widehat{w}^2_\delta )}{\delta }\le \frac{\mathcal {E}_\delta (\Omega _0;W_\delta ,w_\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }, \end{aligned}$$

(58)

so that

$$\begin{aligned} \liminf \limits _{\delta \searrow 0}\frac{\mathcal {E}_\delta (\Omega _0;W_\delta ,w_\delta )-\mathcal {E}(\Omega _0; \widehat{W}^2_\delta , \widehat{w}^2_\delta )}{\delta }\le \liminf \limits _{\delta \searrow 0}\frac{\mathcal {E}(\Omega _\delta ;W^\delta ,w^\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }.\end{aligned}$$

(59)

Invoking the asymptotic expansion of the transformed potential deformation energy (54), we can calculate the limit superior on the right-hand side of (57) and the limit inferior on the left-hand side of (59). Direct computations based on Lemmas 2.7 and 2.9 reveal that they are finite and equal. Therefore,

$$\begin{aligned} \lim \limits _{\delta \searrow 0}\frac{\mathcal {E}(\Omega _\delta ;W^\delta ,w^\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }= & {} \frac{1}{2}\int \limits _{\Omega \setminus \overline{\omega }_0}\left( A_p(V; W_0,W_0)+A_v(V; w_0,w_0)\right) \,\mathrm{d}x\nonumber \\&-\int \limits _{\Omega \setminus \overline{\omega }_0}\left( \mathrm{div}_x(F\otimes V)\cdot W+\mathrm{div}_x(fV)w\right) \,\mathrm{d}x\nonumber \\&+\int \limits _{\Omega \setminus \overline{\omega }_0}\sigma (W_0):\varepsilon (\widetilde{W}_0)\,\mathrm{d} x-\int \limits _{\Omega \setminus \overline{\omega }_0} m (w_0):\nabla ^2_x\widetilde{w}_0\,\mathrm{d} x\nonumber \\&-\int \limits _{\Omega _0}(F\cdot \widetilde{W}_0+f\widetilde{w}_0)\,\mathrm{d}x. \end{aligned}$$

(60)

For $\delta <0,$ we just have to flip inequality signs in (56) and (58) and reason as before to prove that

$$\begin{aligned} \lim \limits _{\delta \nearrow 0}\frac{\mathcal {E}(\Omega _\delta ;W^\delta ,w^\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }=\lim \limits _{\delta \searrow 0}\frac{\mathcal {E}(\Omega _\delta ;W^\delta ,w^\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }. \end{aligned}$$

In view of the above discussion, the limit in (22) exists and coincides with the right-hand side of (60). Thus it remains to show that the right-hand side of (60) is equal to $\mathfrak {G}(\Omega _0;V).$ To do so, let us set

$$\begin{aligned} \int \limits _{\Omega \setminus \overline{\omega }_0}\sigma (W_0):\varepsilon (\widetilde{W}_0)\,\mathrm{d} x&- \int \limits _{\Omega \setminus \overline{\omega }_0} m (w_0):\nabla ^2_x\widetilde{w}_0\,\mathrm{d} x\\&-\int \limits _{\Omega _0}(F\cdot \widetilde{W}_0+f\widetilde{w}_0)\,\mathrm{d}x=I_{p}(W_0,\widetilde{W}_0)+I_{v}(w_0,\widetilde{w}_0). \end{aligned}$$

An application of the Green formula [17, Theorem 1.17] gives

$$\begin{aligned} I_{p}(W_0,\widetilde{W}_0)= & {} \int \limits _{\Omega \setminus \overline{\omega }_0}(\mathrm{div}_x\sigma (W_0)-F)\cdot \widetilde{W}_0\,\mathrm{d} x -\int \limits _{\omega _0}F\cdot \widetilde{W}_0\,\mathrm{d}x\\&-\left\langle \sigma ^+(W_0)\nu _0, Q_0\right\rangle _{1/2,\partial \omega _0}\\&-\left\langle \sigma ^+_{\nu _0}(W_0), P^+_0\cdot \nu _0\right\rangle _{1/2,\partial \omega _0}-\left\langle \sigma ^+_{\tau _{0}}(W_0), P^+_{0\tau _0}\right\rangle _{1/2,\partial \omega _0}, \end{aligned}$$

where $P^+_{0\tau _0}=P^+_0-(P^+_0\cdot \nu _0)\nu _0.$ Since $P^+_0=(0,0)^\mathrm{T}$ on $\partial \omega _0\setminus {\gamma }_0,$ we find

$$\begin{aligned}&\left\langle \sigma ^+_{\nu _0}(W_0), P^+_0\cdot \nu _0\right\rangle _{1/2,\partial \omega _0}+\left\langle \sigma ^+_{\tau _{0}}(W_0), P^+_{0\tau _0}\right\rangle _{1/2,\partial \omega _0}= \left\langle \sigma ^+_{\nu _0}(W_0),P^+_0\cdot \nu _0\right\rangle ^{00}_{1/2,\gamma _0}\\&\quad +\left\langle \sigma ^+_{\tau _{0}}(W_0), P^+_{0\tau _0}\right\rangle ^{00}_{1/2,\gamma _0}. \end{aligned}$$

Owing to the equilibrium equations (8) and Remark 1.1, we infer

$$\begin{aligned} I_{p}(W_0,\widetilde{W}_0)= & {} -\int \limits _{\omega _0}F\cdot B_0V\,\mathrm{d}x-\left\langle \sigma ^+(W_0)\nu _0, B_0V\right\rangle _{1/2,\partial \omega _0}\nonumber \\&-\left\langle \sigma ^+_{\nu _0}(W_0),\nabla _x V\llbracket W_0\rrbracket \cdot \nu _0\right\rangle ^{00}_{1/2,\gamma _0}. \end{aligned}$$

(61)

Let us turn now to $I_v(w_0,\widetilde{w}_0)$. We employ the Green formula [17, Theorem 1.19] to deduce

$$\begin{aligned} I_v(w_0,\widetilde{w}_0)= & {} \int \limits _{\Omega \setminus \overline{\omega }_0}(-\mathrm{div}_x\mathrm{div}_x m(w_0)-f)\widetilde{w}_0\,\mathrm{d}x-\int \limits _{\omega _0}f\widetilde{w}_0\,\mathrm{d}x+ \left\langle t^+_{\nu _0}(w_0),{q}_0\right\rangle _{3/2,\partial \omega _0}\\&-\left\langle m^+_{\nu _0}(w_0),\nabla _x{q}_0\cdot \nu _0\right\rangle _{1/2,\partial \omega _0}+ \left\langle t^+_{\nu _0}(w_0),{p}^+_0\right\rangle _{3/2,\partial \omega _0}\\&-\left\langle m^+_{\nu _0}(w_0),(\nabla _x{p}_0)^+\cdot \nu _0\right\rangle _{1/2,\partial \omega _0}. \end{aligned}$$

Taking into account that $p^+_0=0$ and $(\nabla _x{p}_0)^+\cdot \nu _0=0$ on $\partial \omega _0\setminus {\gamma }_0,$ we thereby obtain

$$\begin{aligned}&\left\langle t^+_{\nu _0}(w_0),{p}^+_0\right\rangle _{3/2,\partial \omega _0}-\left\langle m^+_{\nu _0}(w_0),(\nabla _x{p}_0)^+\cdot \nu _0\right\rangle _{1/2,\partial \omega _0}\\&=\left\langle t^+_{\nu _0}(w_0),{p}^+_0\right\rangle ^{00}_{3/2,\gamma _0}-\left\langle m^+_{\nu _0}(w_0),(\nabla _x{p}_0)^+\cdot \nu _0\right\rangle ^{00}_{1/2,\gamma _0}. \end{aligned}$$

Furthermore, the equilibrium equation (10) and Remark 1.1 yield

$$\begin{aligned} I_v(w_0,\widetilde{w}_0)= & {} -\int \limits _{\omega _0}fA_0\cdot V\,\mathrm{d}x+\left\langle t^+_{\nu _0}(w_0),A_0\cdot V\right\rangle _{3/2,\partial \omega _0}\nonumber \\&\nonumber \\&-\left\langle m^+_{\nu _0}(w_0),\nabla _x(A_0\cdot V)\cdot \nu _0\right\rangle _{1/2,\partial \omega _0}\nonumber \\&-\left\langle m^+_{\nu _0}(w_0),\llbracket \nabla _x w_0\rrbracket \cdot \left( \nabla _x V+\left( \nabla _x V\right) ^\mathrm{T}\right) \nu _0\right\rangle ^{00}_{1/2,\gamma _0}. \end{aligned}$$

(62)

Combining (60)–(62) completes the proof. $\square$

Remark 2.10

The structural assumption (1) arises in dimension reductions from $3\mathrm{D}$ to $2\mathrm{D}$ within the Kirchhoff–Love energy scaling regime [8,9,10]. In the case of an undamaged elastic plate, this assumption provides in particular that the limit stored elastic energy is a sum of uncoupled bending and stretching energies.

From a mathematical standpoint, we have no difficulty in treating a more general case in which the elastic properties of the plate depend on the spartial variable x as in [17]. Indeed, let $\mathbb {C}^p, \mathbb {C}^v\in C^1(\overline{\Omega })^{2\times 2\times 2\times 2}$ be fourth-order tensors. For each $n=p,v,$ the tensor $\mathbb {C}^n$ is assumed to be positive, i.e., there exists $c_{\mathbb {C}^n}>0$ such that

$$\begin{aligned} \mathbb {C}^n(x)X:X\ge c_{\mathbb {C}^n}|X|^2\quad \text {for all}\quad X\in \mathbb {R}^{2\times 2}_\mathrm{sym}\quad \text {and all}\quad x\in \Omega , \end{aligned}$$

and symmetric

$$\begin{aligned} \mathbb {C}^n_{ijkl}=\mathbb {C}^n_{ijlk}=\mathbb {C}^n_{klij},\quad i,j,k,l=1,2. \end{aligned}$$

We denote by $\mathbb {C}^{n,V}$ the fourth-order tensor with components $\mathbb {C}^{n,V}_{ijkl}=\nabla _x\mathbb {C}^n_{ijkl}\cdot V,$ $i,j,k,l=1,2,$ and replace the constitutive relations (2),(3) (see also (9),(11)) by

$$\begin{aligned} \sigma =\mathbb {C}^p\varepsilon (W),\quad m=-\mathbb {C}^v\nabla ^2_x w \end{aligned}$$

and the bilinear forms (23),(24) by

$$\begin{aligned} {A}_p(V; U,W)= & {} \mathrm{div}_xV\sigma (U):\varepsilon (W)-\sigma (U):E\left( \nabla _x V; W\right) \\&\quad - \sigma (W):E\left( \nabla _x V; U\right) +\mathbb {C}^{p,V}\varepsilon (U):\varepsilon (W),\\ {A}_v(V; u,w)= & {} - \mathrm{div}_xV m(u):\nabla ^2_x w+m(u):\left( 2\nabla ^2_x w\nabla _xV+N\nabla _xw\right) \\&\quad + m(w):\left( 2\nabla ^2_x u\nabla _xV+N\nabla _xu\right) +\mathbb {C}^{v,V}\nabla ^2_x u:\nabla ^2_x w. \end{aligned}$$

Exploiting the asymptotic representations

$$\begin{aligned} \mathbb {C}^n\circ T_\delta =\mathbb {C}^n+\delta \mathbb {C}^{n,V}+r^{n6}_\delta ,\quad \Vert r^{n6}_\delta \Vert _{L^\infty (\Omega )^{2\times 2\times 2\times 2}}=o(\delta ), \end{aligned}$$

we can easily adapt the arguments above to show that all our results remain valid.

Remark 2.11

(Griffith formula) We now demonstrate how one can apply the result obtained to derive the Griffith formula for the energy release rate associated with extension of the crack along the interface.

We assume that there exists a simple curve $\Lambda$ of class $C^{2,1}$ such that ${\gamma }_0\varsubsetneq \Lambda \subset \partial \omega _0,$ and let $\lambda :(0,\mathcal {H}^1(\Lambda ))\rightarrow \mathbb {R}^2$ be its arc-length parametrization, whereas $\gamma _0=\lambda (0,l)$ with $l=\mathcal {H}^1(\gamma _0).$ For technical reason, we require the validity of G1 for the domain $\Omega \setminus \overline{\Lambda }.$ The cylindrical surface $\Lambda ^\mathrm{3D}$ is a prescribed crack path. For $\delta \in (-l,\mathcal {H}^1(\Lambda \setminus \gamma _0)),$ the admissible crack set is $\gamma ^\mathrm{3D}_\delta =\lambda (0,l+\delta )\times (-h,h),$ and we put $\Omega _\delta =\Omega \setminus \overline{\gamma }_\delta .$

The energy release rate associated with crack extension $\mathcal {ERR}(\Omega _0)$ is formally defined as the opposite of the derivative of the potential deformation energy with respect to the crack length:

$$\begin{aligned} \mathcal {ERR}(\Omega _0)=- \lim \limits _{\delta \rightarrow 0}\frac{\mathcal {E}(\Omega _\delta ;W^\delta ,w^\delta )-\mathcal {E}(\Omega _0;W_0,w_0)}{\delta }, \end{aligned}$$

(63)

where $(W^\delta ,w^\delta )$ and $(W_0,w_0)$ are the displacement fields corresponding to $\Omega _\delta ,$ $\delta \ne 0,$ and $\Omega _0,$ respectively. Thanks to the regularity assumptions on $\Lambda ,$ we can find $\mu >0$ and $\psi _l \in C^{2,1}([\lambda _1(l)-\mu ,\lambda _1(l)+\mu ])$ such that

$$\begin{aligned} \Lambda =\{\,(x_1,\psi _l(x_1)):x_1\in [\lambda _1(l)-\mu ,\lambda _1(l)+\mu ]\,\}. \end{aligned}$$

Choose a cut-off function $\theta \in C^\infty (\mathbb {R}^2)$ so that $\theta =1$ on ${B}_{\mu /3}(0)$ and $\theta =0$ outside ${B}_{\mu /2}(0),$ and define the mapping ${T}_{l,\delta }:\mathbb {R}^2\rightarrow \mathbb {R}^2$ by

$$\begin{aligned} {T}_{l,\delta }(x)=x+ \begin{pmatrix} (\lambda _1(l+\delta )-\lambda _1(l))\theta (\lambda (l)-x)\\ \psi _l(x_1+(\lambda _1(l+\delta )-\lambda _1(l))\theta (\lambda (l)-x))-\psi _l(x_1) \end{pmatrix}. \end{aligned}$$

Hence there exists $\delta _0>0$ such that $T_{l,\cdot }\in C^{2,1}([-\delta _0,\delta _0]\times \mathbb {R}^2)^2$ and, for every $\delta \in [-\delta _0,\delta _0],$ the mapping $T_{l,\delta }$ is a $C^{2,1}$-diffeomorphism on $\mathbb {R}^2$ with $T_{l,\delta }(\Omega _0)=\Omega _\delta ,$ $T_{l,\delta }(\lambda (l))=\lambda (l+\delta ),$ $T_{l,\delta }(\gamma _0)=\gamma _\delta ,$ and $T_{l,\delta }(\omega _0)=\omega _0,$ see [24, Lemma 3.5]. By construction, the vector field

$$\begin{aligned} V_\mathrm{Gr}(x)=\partial _\delta ({T}_{l,\delta }(x))|_{\delta =0}=\lambda ^{\prime }_1(l)\theta (\lambda (l)-x) (1,\psi _{l,1}(x_1))^\mathrm{T} \end{aligned}$$

belongs to $C^{1,1}(\mathbb {R}^2)^2$ and is compactly supported in $\Omega .$ We thus deduce that the energy release rate associated with crack extension (63) is well defined and can be calculated as

$$\begin{aligned} \mathcal {ERR}(\Omega _0)=-\mathfrak {G}(\Omega _0;V_\mathrm{Gr}). \end{aligned}$$

We conclude this remark and paper by noting that, in view of the definition (63), the energy release rate associated with crack extension $\mathcal {ERR}(\Omega _0)$ is independent of the choice of the cut-off function $\theta .$