A decomposition result for Kirchhoff plate bending problems and a new discretization approach

TL;DR

This paper introduces a novel mixed variational formulation and discretization method for Kirchhoff plate bending problems, enabling the use of standard second-order problem solvers within a modular framework.

Contribution

It develops a new regular decomposition-based approach for Kirchhoff plates, extending previous biharmonic problem results to more complex boundary conditions.

Findings

New discretization methods based on second-order problem solvers

Flexible approach adaptable to various boundary conditions

Theoretical validation of the regular decomposition and Brezzi's conditions

Abstract

A new approach is introduced for deriving a mixed variational formulation for Kirchhoff plate bending problems with mixed boundary conditions involving clamped, simply supported, and free boundary parts. Based on a regular decomposition of an appropriate nonstandard Sobolev space for the bending moments, the fourth-order problem can be equivalently written as a system of three (consecutively to solve) second-order problems in standard Sobolev spaces. This leads to new discretization methods, which are flexible in the sense, that any existing and well-working discretization method and solution strategy for standard second-order problems can be used as a modular building block of the new method. Similar results for the first biharmonic problem have been obtained in our previous work [W. Krendl, K. Rafetseder and W. Zulehner, A decomposition result for biharmonic problems and the Hellan-Herrmann-Johnson method, ETNA, 2016]. The extension to more general boundary conditions encounters several difficulties including the construction of an appropriate nonstandard Sobolev space, the verification of Brezzi's conditions, and the adaptation of the regular decomposition.

Tables

Table 1. Table 1: Discretization errors, k = 1 𝑘 1 k=1 .

$L$	${\Vert w-w_h\Vert }_1$	order	${\Vert 𝑴-{𝑴}_h\Vert }_0$	order	${\Vert p-p_h\Vert }_0$	order	${\Vert \varphi -{\varphi }_h\Vert }_1$	order
4	$1.09\cdot {10}^{-1}$	$0.992$	$1.24\cdot {10}^{-1}$	$0.974$	$1.43\cdot {10}^{-2}$	$1.985$	$1.04\cdot {10}^{-1}$	$1.050$
5	$5.47\cdot {10}^{-2}$	$0.998$	$6.26\cdot {10}^{-2}$	$0.993$	$3.59\cdot {10}^{-3}$	$1.999$	$5.17\cdot {10}^{-2}$	$1.017$
6	$2.73\cdot {10}^{-2}$	$0.999$	$3.13\cdot {10}^{-2}$	$0.998$	$8.90\cdot {10}^{-4}$	$2.011$	$2.56\cdot {10}^{-2}$	$1.012$
7	$1.36\cdot {10}^{-2}$	$0.999$	$1.56\cdot {10}^{-2}$	$0.999$	$2.14\cdot {10}^{-4}$	$2.052$	$1.24\cdot {10}^{-2}$	$1.036$

Table 2. Table 2: Discretization errors, k = 2 𝑘 2 k=2 .

$L$	${\Vert w-w_h\Vert }_1$	order	${\Vert 𝑴-{𝑴}_h\Vert }_0$	order	${\Vert p-p_h\Vert }_0$	order	${\Vert \varphi -{\varphi }_h\Vert }_1$	order
4	$4.33\cdot {10}^{-2}$	$2.071$	$1.75\cdot {10}^{-1}$	$2.104$	$1.22\cdot {10}^{-2}$	$3.160$	$2.04\cdot {10}^{-1}$	$2.055$
5	$1.06\cdot {10}^{-2}$	$2.018$	$4.29\cdot {10}^{-2}$	$2.030$	$1.49\cdot {10}^{-3}$	$3.042$	$5.05\cdot {10}^{-2}$	$2.017$
6	$2.66\cdot {10}^{-3}$	$2.004$	$1.06\cdot {10}^{-2}$	$2.008$	$1.85\cdot {10}^{-4}$	$3.010$	$1.25\cdot {10}^{-2}$	$2.005$
7	$6.65\cdot {10}^{-4}$	$2.001$	$2.66\cdot {10}^{-3}$	$2.002$	$2.30\cdot {10}^{-5}$	$3.002$	$3.13\cdot {10}^{-3}$	$2.003$

Table 3. Table 3: Discretization errors, k = 3 𝑘 3 k=3 .

$L$	${\Vert w-w_h\Vert }_1$	order	${\Vert 𝑴-{𝑴}_h\Vert }_0$	order	${\Vert p-p_h\Vert }_0$	order	${\Vert \varphi -{\varphi }_h\Vert }_1$	order
4	$2.75\cdot {10}^{-3}$	$3.084$	$1.10\cdot {10}^{-2}$	$3.105$	$7.69\cdot {10}^{-4}$	$4.219$	$1.27\cdot {10}^{-2}$	$3.019$
5	$3.46\cdot {10}^{-4}$	$2.994$	$1.38\cdot {10}^{-3}$	$2.989$	$4.63\cdot {10}^{-5}$	$4.054$	$1.62\cdot {10}^{-3}$	$2.966$
6	$4.37\cdot {10}^{-5}$	$2.985$	$1.75\cdot {10}^{-4}$	$2.978$	$2.86\cdot {10}^{-6}$	$4.012$	$2.07\cdot {10}^{-4}$	$2.972$
7	$5.50\cdot {10}^{-6}$	$2.989$	$2.22\cdot {10}^{-5}$	$2.985$	$1.78\cdot {10}^{-7}$	$4.005$	$2.60\cdot {10}^{-5}$	$2.995$

Equations

$$\operatorname{div}\operatorname{Div}\big{(}\mathcal{C}\nabla^{2}w\big{)}=f\quad\text{in }\Omega$$

$$\Gamma ={\mathcal{V}}_{\Gamma ,}\cup {\mathcal{E}}_{\Gamma ,}\text{with}{\mathcal{E}}_{\Gamma ,}=\bigcup _{k=1}^K{E}_k,$$

$$\mathbf{M}=-\mathcal{C}{\nabla }^{2}w,$$

$${\mathbf{M}}_{nn}=\mathbf{M}n\cdot n,{\mathbf{M}}_{nt}=\mathbf{M}n\cdot t$$

$$-\mathrm{div}\mathrm{Div}\mathbf{M}=f\text{in }\Omega ,\text{with}\mathbf{M}=-\mathcal{C}{\nabla }^{2}w$$

$$w=0,$$

$${\mathbf{M}}_{nn}=0,$$

$$[[{\mathbf{M}}_{nt}{]]}_x={\mathbf{M}}_{n_1{t}_1}(x)-{\mathbf{M}}_{n_2{t}_2}(x)=0\text{for all}x\in {\mathcal{V}}_{\Gamma ,f},$$

$${\partial }_{n}w=0,{\partial }_t{\mathbf{M}}_{nt}+\mathrm{Div}\mathbf{M}\cdot n=0$$

$${\int }_{\Omega }\mathcal{C}{\nabla }^{2}w:{\nabla }^{2}vdx=⟨F,v⟩\text{for all}v\in W,$$

$$W=\{v\in H^2(\Omega ):v=0,{\partial }_{n}v=0\text{on}{\Gamma }_c,v=0\text{on}{\Gamma }_s\}$$

$$\Vert w{\Vert }_W\le c\Vert F{\Vert }_{W^{\ast }},$$

$$\nabla v,\nabla \psi ,\mathrm{curl}v,\mathrm{Curl}\psi ,\mathrm{div}\psi ,\mathrm{Div}\mathbf{N},\mathrm{rot}\psi ,\mathrm{Rot}\mathbf{N}$$

$$H^1(\Omega ),(H^1(\Omega ){)}^2,H(\mathrm{curl},\Omega ),H(\mathrm{Curl},\Omega ),\dots .$$

$$\nabla v$$

$$\nabla \psi$$

$$\mathrm{div}\psi$$

$$\mathrm{Div}\mathbf{N}$$

$$\epsilon (\psi )=\frac{1}{2}(\nabla \psi +(\nabla \psi {)}^T),\mathrm{symCurl}\psi =\frac{1}{2}(\mathrm{Curl}\psi +(\mathrm{Curl}\psi {)}^T).$$

$$⟨Bx,y⟩=⟨G,x⟩\text{for all}x\in D(B).$$

$$⟨B^{\ast }y,x⟩=⟨Bx,y⟩\text{for all}x\in D(B),y\in D(B^{\ast }).$$

$$⟨{\nabla }^{2}v,\mathbf{N}⟩={\int }_{\Omega }{\nabla }^{2}v:\mathbf{N}dx\text{for}v\in D(B)=W,\mathbf{N}\in Y={\mathbf{L}}^2(\Omega {)}_{\mathrm{sym}},$$

$${\mathbf{L}}^2(\Omega {)}_{\mathrm{sym}}=\{\mathbf{N}:{\mathbf{N}}_{ij}={\mathbf{N}}_{ji}\in L^2(\Omega ),i,j=1,2\}$$

$$Q=H_{0,{\Gamma }_c\cup {\Gamma }_s}^1(\Omega )=\{v\in H^1(\Omega ):v=0\text{on}{\Gamma }_c\cup {\Gamma }_s\},$$

$$\mathrm{div}\mathrm{Div}\text{and}\mathbf{H}(\mathrm{div}\mathrm{Div},\Omega ;Q^{\ast }{)}_{\mathrm{sym}}$$

$$\mathbf{H}(\mathrm{div}\mathrm{Div},\Omega ;Q^{\ast }{)}_{\mathrm{sym}}=\{\mathbf{N}\in {\mathbf{L}}^2(\Omega {)}_{\mathrm{sym}}:\text{the functional}$$

$$G:v\mapsto {\int }_{\Omega }{\nabla }^{2}v:\mathbf{N}dx,v\in W,\text{is bounded w.r.t. the }Q\text{-norm}\},$$

$$(\mathbf{M},\mathbf{N}{)}_{{\mathcal{C}}^{-1}}$$

$$⟨\mathrm{div}\mathrm{Div}\mathbf{M},v⟩$$

$$\mathbf{V}=\mathbf{H}(\mathrm{div}\mathrm{Div},\Omega ;Q^{\ast }{)}_{\mathrm{sym}},Q=H_{0,{\Gamma }_c\cup {\Gamma }_s}^1(\Omega ),$$

$$a(\mathbf{M},\mathbf{N})$$

$$b(\mathbf{M},v)$$

$$⟨\mathcal{A}(\mathbf{M},w),(\mathbf{N},v)⟩=a(\mathbf{M},\mathbf{N})+b(\mathbf{N},w)+b(\mathbf{M},v).$$

Full text

Similar results for the first biharmonic problem have been obtained in our previous work [W. Krendl, K. Rafetseder and W. Zulehner, A decomposition result for biharmonic problems and the Hellan-Herrmann-Johnson method, ETNA, 2016]. The extension to more general boundary conditions encounters several difficulties including the construction of an appropriate nonstandard Sobolev space, the verification of Brezzi’s conditions, and the adaptation of the regular decomposition.

Key words. Kirchhoff plates, mixed boundary conditions, free boundary, mixed methods, regular decomposition

AMS subject classifications. 65N30, 65N22, 74K20

1 Introduction

We consider the Kirchhoff plate bending problem: For a given load $f$, find the deflection $w$ such that

[TABLE]

with appropriate boundary conditions. Here $\Omega$ is a bounded domain in $\mathbb{R}^{2}$ with a polygonal Lipschitz boundary $\Gamma$, $\operatorname{div}$ denotes the standard divergence of a vector-valued function, $\operatorname{Div}$ the row-wise divergence of a matrix-valued function, $\nabla^{2}$ the Hessian and $\mathcal{C}$ the material tensor. Note that (1) reduces to the biharmonic equation, if $\mathcal{C}$ is the identity. In this paper we focus on mixed methods for the original unknown $w$ and the bending moments $\bm{M}=-\mathcal{C}\nabla^{2}w$ as additional unknowns, which are often quantities of interest on their own.

In our previous work [20] the first biharmonic boundary value problem is studied, which corresponds to the situation of a purely clamped plate. For this model problem a new mixed variational formulation is derived, which satisfies Brezzi’s conditions and is equivalent to the original problem. However, these important properties come at the cost of an appropriate nonstandard Sobolev space ${\bm{V}}$ for $\bm{M}$. Based on a regular decomposition of ${\bm{V}}$, the fourth-order problem can be rewritten as a sequence of three (consecutively to solve) second-order elliptic problems in standard Sobolev spaces. This leads to a new interpretation of known mixed finite element methods as well as to the construction of new discretization methods, see [20] for details. This approach fits into an abstract framework recently presented in [8] for the decomposition of higher-order problems. However, more general boundary conditions (beyond a purely clamped plate) are not considered there as well.

The aim of this paper is to extend the ideas of [20] to the more general situation of a Kirchhoff plate with mixed boundary conditions involving clamped, simply supported, and free boundary parts. This is by far not straight-forward.

The first difficulty arises in the derivation of the new mixed formulation. Contrary to clamped boundary parts, appropriate boundary conditions for $\bm{M}$ have to be incorporated into the definition of the nonstandard Sobolev space ${\bm{V}}$ for simply supported and free boundary parts. In this paper we do this indirectly using the framework of (possibly unbounded) densely defined operators in Hilbert spaces. This approach avoids the direct use of trace operators in nonstandard Sobolev spaces, which would be technically rather involved. A second difficulty is the verification of Brezzi’s conditions. In [20] the main ingredient for the proof of an inf-sup condition is the property that matrix-valued functions of the form $v\bm{I}$ are contained in ${\bm{V}}$, where $\bm{I}$ denotes the identity matrix and $v$ satisfies homogeneous Dirichlet boundary conditions induced by the boundary conditions for $w$. This inclusion is no longer true for problems with free boundary parts. So, a new technique is required for proving the inf-sup condition. A third difficulty arises in the regular decomposition for a similar reason. The first component of the decomposition in [20] is of the form $p\bm{I}$ and is not contained in ${\bm{V}}$ for problems with free boundary parts. So, a new approach is required to pursue the decomposition for problems with free boundary parts. It is shown in this paper how to overcome all these difficulties and how to achieve again a decomposition of the fourth-order problem into three (consecutively to solve) second-order elliptic problems in standard Sobolev spaces.

So far in literature, mixed methods for (1) in $w$ and $\bm{M}$ have been formulated as linear operator equations in function spaces for which either the associated linear operator is not an isomorphism or the involved norms contain a mesh-dependent variant of the $H^{2}$-norm for $w$, see, e.g., [7, 12, 2, 5, 6]. This lack of easy-to-access knowledge on the mapping properties of the involved operators on the continuous level makes it hard to design efficient preconditioners on the discrete level. Our new mixed formulation satisfies Brezzi’s conditions. Therefore, the associated linear operator is an isomorphism. Additionally, the operator is of triangular structure. This naturally leads to the construction of efficient solvers. Moreover, the new mixed formulation is equivalent to the original problem without additional convexity assumptions on $\Omega$, while most of the papers in literature (except for [5]) require $\Omega$ to be convex.

For alternative discretization methods such as conforming, non-conforming, and discontinuous Galerkin methods for the primal formulation as well as alternative mixed methods we refer to the short discussion in [20] and the references cited there. Our approach leads to new discretization methods, which are flexible in the sense, that any existing and well-working discretization method and solution strategy for second-order problems can be used as a modular building block of the new method. One option would be to choose standard $C^{0}$ finite elements for each of the three second-order elliptic problems (for two scalar fields and one vector field) resulting in approximate solutions to $w$ and $\bm{M}$. In [3, 4] a different method was proposed that also uses only standard $C^{0}$ finite element spaces for second-order problems for a formulation in the kinematic variables $w$ and $\nabla w$. For reaching approximate solutions of comparable accuracy the method in [3, 4] requires the approximation of one scalar field less than the approach presented here. However, the linear system resulting from the method in [3, 4] is a coupled system of all degrees of freedoms of one scalar and one vector field, while the method presented here requires to solve linear systems for the degrees of freedom separately for each of the two scalar and the vector field. This reduces the computational costs for direct solvers. For the use of iterative solvers efficent methods for standard second-order problems like multigrid methods can be directly used for each of the three linear systems. Preconditioning is not addressed in [3, 4]. So we feel that our method is competitive with respect to the overall computational efficiency.

The paper is organized as follows. In Section 2 the Kirchhoff plate bending problem is introduced. Section 3 contains a new mixed formulation, for which well-posedness and equivalence to the original problem is shown. A regular decomposition of the nonstandard Sobolev space for $\bm{M}$ is derived in Section 4 and the resulting decoupled formulation is presented. The decoupled formulation leads in a natural way to the construction of a new discretization method, which is introduced in Section 5, and for which a priori error estimates are derived. The paper closes with numerical experiments in Section 6.

2 The Kirchhoff plate bending problem

We consider the Kirchhoff plate bending problem of a linearly elastic plate where the undeformed mid-surface is described by a domain $\Omega\subset\mathbb{R}^{2}$ with a polygonal Lipschitz boundary $\Gamma$. In what follows, let the boundary $\Gamma$ be written in the form

[TABLE]

where $E_{k}$, $k=1,2,\ldots,K$, are the edges of $\Gamma$, considered as open line segments and $\mathcal{V}_{\Gamma,}$ denotes the set of corner points in $\Gamma$. Furthermore, $n=(n_{1},n_{2})^{T}$ and $t=(-n_{2},n_{1})^{T}$ represent the unit outer normal vector and the unit counterclockwise tangent vector to $\Gamma$, respectively.

The plate is considered to be clamped on a part $\Gamma_{c}\subset\Gamma$, simply supported on $\Gamma_{s}\subset\Gamma$, free on $\Gamma_{f}\subset\Gamma$ with $\Gamma=\Gamma_{c}\cup\Gamma_{s}\cup\Gamma_{f}$. We assume that each edge $E\in\mathcal{E}_{\Gamma,}$ is contained in exactly one of the sets $\Gamma_{c}$, $\Gamma_{s}$, $\Gamma_{f}$, and the edges are maximal in the sense that two edges with the same boundary condition do not meet at an angle of $\pi$. Recall the definition of the bending moments $\bm{M}$ by the Hessian of the deflection $w$

[TABLE]

where $\mathcal{C}$ is the fourth-order material tensor. The tensor $\mathcal{C}$ is assumed to be symmetric and positive definite on symmetric matrices, $\lambda_{min}(\mathcal{C})$ and $\lambda_{max}(\mathcal{C})$ denote the minimal and maximal eigenvalue of $\mathcal{C}$, respectively. For example, for isotropic materials with Poisson ratio $\nu$, the material tensor $\mathcal{C}$ is given by \mathcal{C}\bm{N}=D\bigl{(}(1-\nu)\bm{N}+\nu\ \mathrm{tr}(\bm{N})\bm{I}\bigr{)} for matrices $\bm{N}$, where $D>0$ depends on material constants, $\bm{I}$ is the identity matrix and $\mathrm{tr}$ is the trace operator for matrices (cf. [27]). We introduce the following notations

[TABLE]

for the normal-normal component and the normal-tangential component of $\bm{M}$, where the symbol $\operatorname{\cdot}$ denotes the Euclidean inner product. The classical Kirchhoff plate bending problem reads as follows (cf. [27]): For given load $f$, find a deflection $w$ such that

[TABLE]

and the boundary conditions

[TABLE]

where $\partial_{t}$ denotes the tangential derivative, and the corner conditions

[TABLE]

where $\mathcal{V}_{\Gamma,f}$ denotes the set of corner points whose two adjacent edges (with corresponding normal and tangent vectors $n_{1}$, $t_{1}$ and $n_{2}$, $t_{2}$) belong to $\Gamma_{f}$.

Remark 2.1*.*

There is a fourth type of boundary condition given by

[TABLE]

with corner conditions of the form (4), which appears, e.g., in boundary value problems of the Cahn-Hilliard equation. The theory presented in the following can easily be extended to mixed boundary conditions including also this fourth type.

A standard (primal) variational formulation of (3) is given as follows: find $w\in W$ such that

[TABLE]

with the Frobenius inner product $\bm{A}:\bm{B}=\sum_{i,j}\bm{A}_{ij}\,\bm{B}_{ij}$ for matrices $\bm{A},\bm{B}$, the right-hand side $\langle F,v\rangle=\int_{\Omega}f\,v\ dx$, and the function space

[TABLE]

with associated norm $\|v\|_{W}=\|v\|_{2}$. Following, e.g., [1, 22], here and throughout the paper $L^{2}(\Omega)$ and $H^{m}(\Omega)$ denote the standard Lebesgue and Sobolev spaces of functions on $\Omega$ with corresponding norms $\|.\|_{0}$ and $\|.\|_{m}$ for positive integers $m$. For functions on $\Gamma$ we use $L^{2}(\Gamma)$ and $H^{\frac{1}{2}}(\Gamma)$ to denote the Lebesgue space and the trace space of $H^{1}(\Omega)$ with corresponding norms $\|.\|_{0,\Gamma}$ and $\|.\|_{\frac{1}{2},\Gamma}$. Moreover, $H^{1}_{0,\Gamma^{\prime}}(\Omega)$ denotes the set of functions in $H^{1}(\Omega)$ which vanish on a part $\Gamma^{\prime}$ of $\Gamma$. The $L^{2}$-inner product on $\Omega$ and $\Gamma^{\prime}$ are always denoted by $(.,.)$ and $(.,.)_{\Gamma^{\prime}}$, respectively, no matter whether it is used for scalar, vector-valued, or matrix-valued functions.

In order to avoid technicalities related to rigid body motions, we assume throughout the paper that $\Gamma_{c}$ contains at least one non-trivial edge $E\in\mathcal{E}_{\Gamma,}$. Then existence and uniqueness of a solution $w$ to (5) are guaranteed by the theorem of Lax-Milgram (see, e.g., [21, 23]) for even more general right-hand sides $\langle F,v\rangle$, where $F\in W^{*}$. Here we use $H^{*}$ to denote the dual of a Hilbert space $H$ and $\langle.,.\rangle$ for the duality product on $H^{*}\times H$. Moreover, the solution $w$ depends continuously on $F$

[TABLE]

with $c=c^{\prime}/\lambda_{min}(\mathcal{C})$, where $c^{\prime}$ depends only on the constant $c_{F}$ of Friedrichs’ inequality. All results of this paper can easily be extended to the case $\Gamma_{c}=\emptyset$ under appropriate compatibility conditions for the right-hand side $F$.

For scalar functions $v$, vector-valued functions $\psi$, and matrix-valued functions $\bm{N}$ the first order differential expressions

[TABLE]

are defined in the weak sense on the corresponding domains of definition

[TABLE]

In case that all components are in $H^{1}(\Omega)$ they take on their classical form given as follows:

[TABLE]

Moreover, the symmetric gradient and the symmetric $\operatorname{Curl}$ are introduced by

[TABLE]

3 A new mixed variational formulation

For the new mixed variational formulation we introduce the bending moments $\bm{M}$, as defined in (2), as auxiliary variable. Then the Kirchhoff plate bending problem reads in terms of $\bm{M}$ as stated in (3). Note, the involved operators are the second order differential operators $\nabla^{2}$ and $\operatorname{div}\operatorname{Div}$. In the following we give a formally precise definition of them.

Throughout the paper, the differential expression $\nabla^{2}v$ is only considered for functions $v\in W\subset H^{2}(\Omega)$. Therefore, we define $\nabla^{2}v$ in the standard way as the matrix consisting of all second order partial derivatives. In order to introduce the operator $\operatorname{div}\operatorname{Div}$ we use the classical concept of (possibly unbounded) densely defined linear operators $B$. Later on we consider instead of a general operator $B$ the Hessian $\nabla^{2}$ and define $\operatorname{div}\operatorname{Div}$ as its adjoint.

We consider an operator $B\colon D(B)\subset X\rightarrow Y^{*}$, where $X$ and $Y$ are Hilbert spaces and $D(B)$, the domain of definition of $B$, is dense in $X$. The adjoint $B^{*}\colon D(B^{*})\subset Y\rightarrow X^{*}$ is then defined as follows: $y\in D(B^{*})$ if and only if $y\in Y$ and there is a linear functional $G\in X^{*}$ such that

[TABLE]

In this case we define $B^{*}y=G$. Note that $\langle B^{*}y,x\rangle$ is well-defined for $x\in X$ and $y\in D(B^{*})$ and we have in particular

[TABLE]

The domain $D(B^{*})$ is a Hilbert space w.r.t. the graph norm $\|y\|_{D(B^{*})}=(\|y\|_{Y}^{2}+\|B^{*}y\|_{X^{*}}^{2})^{\frac{1}{2}}$.

As already indicated above, it is quite natural to choose $D(B)=W$ and to define $B=\nabla^{2}$ as an operator mapping to $Y=\bm{L}^{2}(\Omega)_{\mathrm{sym}}$, (or, more precisely, to the dual of $Y$,) given by

[TABLE]

where $\bm{L}^{2}(\Omega)_{\mathrm{sym}}$ denotes the space of symmetric matrix-valued functions given by

[TABLE]

and equipped with the standard $L^{2}$-norm $\|\bm{N}\|_{0}$ for a matrix-valued function $\bm{N}$.

Keep in mind, in the following we always fix $D(B)=W$ and obtain for the adjoint $B^{*}$ different domains of definition $D(B^{*})$, which strongly depend on the choice of $X$. There are several options how to choose $X$. However, according to the discussion from above, there is a restriction to meet: $D(B)$ is a dense subset of $X$. Now we discuss three possible choices for $X$. A first and trivial option would be $X=W$. Then it is easy to see that $D(B^{*})=\bm{L}^{2}(\Omega)_{\mathrm{sym}}$. Note that for this choice we have $X\subset H^{2}(\Omega)$, so the disadvantage for the mixed method is to work with a second-order Sobolev space for $w$. A second option would be $X=L^{2}(\Omega)$. Then it turns out that $D(B^{*})\subset\bm{H}(\operatorname{div}\operatorname{Div},\Omega)_{\text{sym}}=\{\bm{N}\in\bm{L}^{2}(\Omega)_{\mathrm{sym}}:\operatorname{div}\operatorname{Div}\bm{N}\in L^{2}(\Omega)\}$, where here $\operatorname{div}\operatorname{Div}$ is defined in the distributional sense. This time the disadvantage for the mixed method is to work with a second-order Sobolev space for $\bm{M}$.

The idea for the new mixed formulation is to distribute the smoothness requirements evenly among $w$ and $\bm{M}$ by choosing the space $X$ as an intermediate space between $W$ and $L^{2}(\Omega)$. In particular, we propose to set $X$ equal to $Q$, given by

[TABLE]

equipped with the norm $\|v\|_{Q}=\|v\|_{1}$.

Remark 3.1*.*

The space ${Q}$ is the interpolation space between $W$ and $L^{2}(\Omega)$. Note that we only use interpolation as motivation, but do not rely in the following on results from interpolation theory.

This choice for $X$ meets the required condition:

Lemma 3.2**.**

The subspace $W$ is dense in ${Q}$.

Proof.

We follow the lines of the proof of [13, Theorem 1.6.1]. In order to verify the density, we have to check that the trace space $\gamma(W)$ is dense in the trace space $\gamma({Q})$, where $\gamma$ is the standard $H^{1}$-trace representing the value on the boundary. By considering the situation locally near each corner in $\mathcal{V}_{\Gamma,}$, the required density follows from the density of $C^{\infty}_{0}(\mathbb{R}^{+})$ in $H^{\frac{1}{2}}(\mathbb{R}^{+})$ and $\tilde{H}^{\frac{1}{2}}(\mathbb{R}^{+})$; see [14]. ∎

Now we leave the abstract framework and use, from now on, the notations

[TABLE]

instead of $B^{*}$ and $D(B^{*})$, respectively, with $X={Q}$ and unchanged $D(B)=W$ and $Y=\bm{L}^{2}(\Omega)_{\mathrm{sym}}$. In consistence with the abstract framework, the Hilbert space $\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$ is explicitly given by

[TABLE]

equipped with the norm $\|\bm{N}\|_{\operatorname{div}\operatorname{Div};{Q}^{*}}=(\|\bm{N}\|^{2}_{0}+\|\operatorname{div}\operatorname{Div}\bm{N}\|^{2}_{{Q}^{*}})^{1/2}$.

This motivates the new mixed formulation as follows: For $F\in{Q}^{*}$, find $\bm{M}\in{\bm{V}}$ and $w\in{Q}$ such that

[TABLE]

with the function spaces

[TABLE]

equipped with the norms $\|\bm{N}\|_{\bm{V}}=\|\bm{N}\|_{\operatorname{div}\operatorname{Div};{Q}^{*}}$ and $\|v\|_{Q}=\|v\|_{1}$. Here, we use the notation $(\bm{M},\bm{N})_{\mathcal{C}^{-1}}=(\mathcal{C}^{-1}\bm{M},\bm{N})$.

The first line in (12) comes from the relation $\mathcal{C}^{-1}\bm{M}+\nabla^{2}w=0$ between bending moment $\bm{M}$ and deflection $w$, the second line originates from (3). Note that we require additional regularity of $F$, namely $F\in{Q}^{*}$. In contrast, for the primal problem (5) we need only $F\in W^{*}$.

Remark 3.3*.*

The operator $\operatorname{div}\operatorname{Div}$ as defined in (10) and the composition of the first order operators $\operatorname{div}$ and $\operatorname{Div}$ introduced at the end of Section 2 differ in two ways. First of all, their domains of definition are different. While $\operatorname{div}(\operatorname{Div}\bm{N})$ is only well-defined for functions $\bm{N}\in{\bm{L}^{2}(\Omega)}$, where $\operatorname{div}\bm{N}$ and $\operatorname{div}\operatorname{Div}\bm{N}$ are $L^{2}$-functions as well, the domain of definition of $\operatorname{div}\operatorname{Div}$ in (11) contains functions $\bm{N}\in\bm{L}^{2}(\Omega)_{\mathrm{sym}}$ with the less restrictive requirement $\operatorname{div}\operatorname{Div}\bm{N}\in{Q}^{*}$. The domain of definition of $\operatorname{div}\operatorname{Div}$ in (11) includes the boundary conditions $\bm{N}_{nn}=0$ on $\Gamma_{s}\cup\Gamma_{f}$ and $\llbracket\bm{N}_{nt}\rrbracket_{x}=0$ for all $x\in\mathcal{V}_{\Gamma,f}$ in a weak sense, as we show later in Theorem 3.9.

But even on the intersection of the domains of definition $\operatorname{div}(\operatorname{Div}\bm{N})$ coincides with $\operatorname{div}\operatorname{Div}\bm{N}$, only if $\bm{N}$ satisfies the boundary condition $\partial_{t}\bm{N}_{nt}+\operatorname{Div}\bm{N}\operatorname{\cdot}n=0$ on $\Gamma_{f}$.

Remark 3.4*.*

In [24, 25, 26] a similar nonstandard Sobolev space is introduced. Note, our way of definition is different and well-suited for the further considerations.

Problem (12) has the typical structure of a saddle point problem

[TABLE]

whose associated linear operator $\mathcal{A}\colon{\bm{V}}\times{Q}\longrightarrow({\bm{V}}\times{Q})^{*}$ is given by

[TABLE]

If the bilinear form $a$ is symmetric, i.e., $a(\bm{M},\bm{N})=a(\bm{N},\bm{M})$, and non-negative, i.e., $a(\bm{N},\bm{N})\geq 0$, which is fulfilled for (12), it is well-known that $\mathcal{A}$ is an isomorphism from ${\bm{V}}\times{Q}$ onto $({\bm{V}}\times{Q})^{*}$, if and only if the following conditions are satisfied; see, e.g., [6]:

$a$ is bounded: There is a constant $\|a\|>0$ such that

[TABLE] 2. 2.

$b$ is bounded: There is a constant $\|b\|>0$ such that

[TABLE] 3. 3.

$a$ is coercive on the kernel of $b$: There is a constant $\alpha>0$ such that

[TABLE]

with $\mathrm{Ker}B=\{\bm{N}\in{\bm{V}}\colon b(\bm{N},v)=0\quad\text{for all}\ v\in{Q}\}$. 4. 4.

$b$ satisfies the inf-sup condition: There is a constant $\beta>0$ such that

[TABLE]

We will refer to these conditions as Brezzi’s conditions with constants $\|a\|$, $\|b\|$, $\alpha$, and $\beta$.

In order to verify Brezzi’s conditions for (12), we need the following result on the relation between the primal problem (5) and the new mixed problem (12).

Theorem 3.5**.**

Let $w$ be the solution of the primal problem (5) for $F\in{Q}^{*}$. Then we have $\bm{M}=-\mathcal{C}\nabla^{2}w\ \in{\bm{V}}$ and $(\bm{M},w)$ solves the mixed problem (12).

Proof.

Since $w\in W$ solves (5), it follows that {\color[rgb]{0,0,0}{\bm{M}}}\in\bm{L}^{2}(\Omega)_{\mathrm{sym}} and

[TABLE]

From the definition of the domain of the adjoint operator in (8) we obtain that $\operatorname{div}\operatorname{Div}\bm{M}=-F\in{Q}^{*}$, which shows that $\bm{M}\in{\bm{V}}$ and that the second equation of (12) is satisfied. Using (9), we receive

[TABLE]

for all $\bm{N}\in{\bm{V}}$, which proves the first equation. ∎

Theorem 3.6**.**

The mixed problem defined by (12) and (13) satisfies Brezzi’s conditions with the constants $\|a\|=1/\lambda_{min}(\mathcal{C})$, $\|b\|=1$, $\alpha=1/\lambda_{max}(\mathcal{C})$ and $\beta=(1+c)^{-1/2}$, where $c=c^{\prime}\,\lambda_{max}(\mathcal{C})/\lambda_{min}(\mathcal{C})$ and $c^{\prime}$ as in (7).

Proof.

The verification of the first three parts of Brezzi’s conditions is simple and, therefore, omitted. For showing the inf-sup condition, let $w^{v}$ be the solution of the primal problem (5) with the right-hand side $F^{v}=-(v,.)_{Q}\in Q^{*}$ for a fixed but arbitrary $v\in{Q}$. From Theorem 3.5 it follows that $\bm{M}^{v}=-\mathcal{C}\nabla^{2}w^{v}\ \in{\bm{V}}$, and $(\bm{M}^{v},w^{v})$ is solution of the corresponding mixed problem (12). From the second line of the mixed formulation (12) we obtain

[TABLE]

and

[TABLE]

Using the stability estimate (7) we obtain

[TABLE]

with $c=c^{\prime}\,\lambda_{max}(\mathcal{C})/\lambda_{min}(\mathcal{C})$. Hence,

[TABLE]

Therefore,

[TABLE]

which completes the proof. ∎

Remark 3.7*.*

The choice of $\bm{M}^{v}$ for proving the inf-sup condition is different to the choice $\bm{M}^{v}=v\,\bm{I}$ as it is used, e.g., in [7, 20], which would not work here, since $v\,\bm{I}\notin{\bm{V}}$ for problems with a free boundary part.

Corollary 3.8**.**

For $F\in{Q}^{*}$, the primal problem (5) and the mixed problem (12) are equivalent in the following sense: If $w$ solves (5), then $\bm{M}=-\mathcal{C}\nabla^{2}w\ \in{\bm{V}}$ and $(\bm{M},w)$ solves (12). Vice versa, if $(\bm{M},w)$ solves (12), then $w\in W$ and solves (5).

Proof.

The first part has already been shown in Theorem 3.5. Since, both problems are uniquely solvable the reverse direction is true as well. ∎

In the case of a purely clamped plate, we have $Q=H_{0}^{1}(\Omega)$. So, in this mixed setting, $w=0$ on $\Gamma$ is treated as an essential boundary condition, while $\partial_{n}w=0$ on $\Gamma$ becomes a natural boundary condition incorporated in the variational formulation. No boundary conditions are prescribed for $\bm{M}$, which makes the definition of an appropriate space ${\bm{V}}$ for $\bm{M}$ much easier. The space ${\bm{V}}$ can be introduced directly as $\{\bm{N}\in\bm{L}^{2}(\Omega)_{\mathrm{sym}}\colon\operatorname{div}\operatorname{Div}\in H^{-1}(\Omega)\}$, where here $\operatorname{div}\operatorname{Div}$ is to be interpreted in the distributional sense, see [20].

The situation is more involved for mixed boundary conditions. Here we have ${Q}=H^{1}_{0,\Gamma_{c}\cup\Gamma_{s}}(\Omega)$. So $w=0$ on $\Gamma_{c}\cup\Gamma_{s}$ is treated as an essential boundary condition, while $\partial_{n}w=0$ on $\Gamma_{c}$ and $\partial_{t}\bm{M}_{nt}+\operatorname{Div}\bm{M}\operatorname{\cdot}n$ on $\Gamma_{f}$ become natural boundary conditions incorporated in the variational formulation. The remaining boundary conditions $\bm{M}_{nn}=0$ on $\Gamma_{s}\cup\Gamma_{f}$ and the corner conditions (4) are treated as essential boundary conditions: They are not explicitly visible but are hidden in the definition of the space $\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$. We doubt that there is an easy and direct way of formulating the corner conditions (4) with the help of pointwise trace operators in $\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$. However, for sufficiently smooth functions $\bm{N}$, the corner conditions can be explicitly extracted using the definition of $\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$, as we will see in the next theorem.

Theorem 3.9**.**

Let $\bm{N}\in\bm{L}^{2}(\Omega)_{\mathrm{sym}}\cap\bm{C}^{1}(\overline{\Omega})$. Then $\bm{N}\in\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$ if and only if

[TABLE]

Proof.

Recall the representation of $\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$ in (11). For $v\in W$, we obtain by integration by parts

[TABLE]

Assume now that $\bm{N}$ satisfies the boundary conditions (14). Then we have

[TABLE]

which is obviously bounded w.r.t. the $H^{1}$-norm. Hence $\bm{N}\in\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$.

On the other hand, if $\bm{N}\in\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$, then the functional $G$ given by (15) is bounded w.r.t. the $H^{1}$-norm. For $v\in H^{2}_{N,\mathcal{V}}(\Omega)\cap W$ we obtain

[TABLE]

where $H^{2}_{N,\mathcal{V}}(\Omega)=\{v\in H^{2}(\Omega):\partial_{n}v=0\ \text{on}\ \Gamma,\quad v(x)=0\ \text{for all}\ x\in\mathcal{V}_{\Gamma,}\}$. Note that all expressions in (16) are continuous in $v$ w.r.t. the $H^{1}$-norm. Analoguously to the proof of Lemma 3.2, one can show that $H^{2}_{N,\mathcal{V}}(\Omega)\cap W$ is dense in $W$ w.r.t. the $H^{1}$-norm. Then it follows that (16) is valid for all $v\in W$. This implies together with (15) that

[TABLE]

From this the boundary conditions in (14) follow by standard arguments. ∎

Remark 3.10*.*

A similar characterization can be derived for piecewise smoothfunctions from $\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$, e.g., for functions from finite element spaces.

4 Regular decomposition

The rather simple proof of the well-posedness of the new mixed formulation (12) comes at the cost of the nonstandard Sobolev space $\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$. The next theorem provides a regular decomposition of this space, which makes $\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$ computationally accessible. In order to derive our decomposition, we need a characterization of the kernel of the distributional $\operatorname{div}\operatorname{Div}$, given in the next lemma, see [17, 20] for a proof.

Lemma 4.1**.**

Let $\Omega$ be simply connected. For $\bm{N}\in\bm{L}^{2}(\Omega)_{\mathrm{sym}}$, we have $\operatorname{div}\operatorname{Div}\bm{N}=0$ in the distributional sense iff $\bm{N}=\operatorname{symCurl}\psi$ for some function $\psi\in(H^{1}(\Omega))^{2}$. The function $\psi$ is unique up to an element from $\text{RT}_{0}=\{ax+b\colon a\in\mathbb{R},\ b\in\mathbb{R}^{2}\}$.

In preparation for the next theorem, observe that $\operatorname{Div}\operatorname{Curl}\psi=0$. Therefore, $\operatorname{Curl}\psi\in\bm{H}(\operatorname{Div},\Omega)=\{\bm{N}\in{\bm{L}^{2}(\Omega)}:\operatorname{Div}\bm{N}\in(L^{2}(\Omega))^{2}\}$, the normal component of $\operatorname{Curl}\psi$ is well-defined on the boundary

[TABLE]

and we have by integration by parts

[TABLE]

for all $\psi\in(H^{1}(\Omega))^{2}$ and $\xi\in(H^{1}(\Omega))^{2}$. For smooth functions $(\operatorname{Curl}\psi)n$ coincides with the tangential derivative of $\psi$, so we use the notation

[TABLE]

Theorem 4.2**.**

Let $\Omega$ be simply connected. For each $\bm{N}\in\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$ there exists a decomposition

[TABLE]

with $\rm\in{Q}=H^{1}_{0,\Gamma_{c}\cup\Gamma_{s}}(\Omega)$ and $\psi\in(H^{1}(\Omega))^{2}$ satisfying the coupling condition

[TABLE]

The function $\rm\in Q$ is the unique solution of the Poisson problem

[TABLE]

and $\psi\in(H^{1}(\Omega))^{2}$ is unique up to an element from $\text{RT}_{0}$. Vice versa, for each $\bm{N}$ given by (17) with $\rm\in Q$ and $\psi\in(H^{1}(\Omega))^{2}$ satisfying (18), it follows that $\bm{N}\in\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$ with $\langle\operatorname{div}\operatorname{Div}\bm{N},v\rangle=-\int_{\Omega}\nabla\rm\operatorname{\cdot}\nabla v\ dx$ for all $v\in{Q}$. Moreover,

[TABLE]

with positive constants $\underline{c}$ and $\overline{c}$, which depend only on the constant $c_{F}$ of Friedrichs’ inequality.

Proof.

Let $\rm\in{Q}$ be the unique solution of the variational problem

[TABLE]

For $v\in C^{\infty}_{0}(\Omega)$ we receive from integration by parts

[TABLE]

This implies $\operatorname{div}\operatorname{Div}(\bm{N}-\rm\bm{I})=0$ in the distributional sense. According to Lemma 4.1, there exists a function $\psi\in(H^{1}(\Omega))^{2}$ such that

[TABLE]

For $v\in W$ integration by parts provides

[TABLE]

With (21) it follows that

[TABLE]

For the reverse direction, assume that (17) and (18) hold. By using the integration by parts formula from above we obtain

[TABLE]

where (18) makes the boundary contributions vanish. This immediately implies that $G\colon v\mapsto\langle\nabla^{2}v,\bm{N}\rangle$ is bounded w.r.t. the ${Q}$-norm, i.e. $\bm{N}\in\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$ with $\operatorname{div}\operatorname{Div}\bm{N}=G\in{Q}^{*}$.

In order to show (20), note that (19) implies

[TABLE]

Hence,

[TABLE]

using Friedrichs’ inequality $\|\rm\|_{0}\leq c_{F}\,|\rm|_{1}$, where $|\rm|_{1}$ denotes the $H^{1}$-semi-norm.

With these inequalities we obtain for the estimate from above

[TABLE]

and for the estimate from below

[TABLE]

Therefore, (20) holds with $1/\underline{c}=\max(2,(1+5c_{F}^{2})(1+c_{F}^{2}))$ and $\overline{c}=4$. ∎

Remark 4.3*.*

By applying Korn’s inequality to $\psi^{\perp}=(-\psi_{2},\psi_{1})^{T}$ we obtain

[TABLE]

where $H/RT_{0}$ denotes the $L^{2}$-orthogonal complement of $RT_{0}$ in $H$ for spaces $H\subset(H^{1}(\Omega))^{2}$. Then it follows that

[TABLE]

So, provided the unique element $\psi\in(H^{1}(\Omega))^{2}/RT_{0}$ is chosen for the decomposition, stability follows from (20) in standard $H^{1}$-norms.

4.1 The coupling condition

Theorem 4.2 shows that each function $\bm{N}\in{\bm{V}}=\bm{H}(\operatorname{div}\operatorname{Div},\Omega;{Q}^{*})_{\mathrm{sym}}$ can be represented by a pair of functions $(\rm,\psi)\in V$ with

[TABLE]

Observe that (18) involves only traces of and $\psi$ on $\Gamma$. We obviously have

[TABLE]

where $\Psi$ is given by the following definition:

Definition 4.4**.**

For fixed $\rm\in{Q}$ we define

[TABLE]

So $\Psi$ consists of those functions $\psi\in(H^{1}(\Omega))^{2}$ that fulfill (18), which becomes in this context a boundary condition for $\psi$ with given . In particular, functions from $\Psi_{0}$, the space associated to $\rm=0$, satisfy the corresponding homogeneous boundary conditions.

It is essential for deriving a decoupled formulation that $V$ is a direct sum of two subspaces, which we will show next. For this we construct, for each $\rm\in{Q}$, a specific function $\psi[\rm]\in\Psi$ as follows: Let $E$ be a fixed edge on $\Gamma_{c}$ with boundary $\partial E=\{x_{A},x_{B}\}$, where $x_{A}$, $x_{B}$ are ordered in counterclockwise direction. We set

[TABLE]

Here $(\sigma,\gamma(\sigma))$ denotes the arc length parametrization of $\Gamma$ in counterclockwise direction with $x_{B}=\gamma(0)$ and $x_{A}=\gamma(\sigma_{E})$ for some $\sigma_{E}>0$. Next we extend $\psi_{\Gamma}[\rm]$ on the whole boundary $\Gamma$ by connecting its values at $\partial E$ linearly on $E$. By this $\psi_{\Gamma}[\rm]$ becomes a continuous function on $\Gamma$ with a weak tangential derivative in $(L^{2}(\Gamma))^{2}$ satisfying

[TABLE]

Finally, let $\psi[\rm]\in(H^{1}(\Omega))^{2}$ be the harmonic extension of $\psi_{\Gamma}[\rm]$, which is well-defined, since we obviously have $\psi_{\Gamma}[\rm]\in(H^{\frac{1}{2}}(\Gamma))^{2}$.

Lemma 4.5**.**

For each $\rm\in Q$, let \psi[q]\in{\color[rgb]{0,0,0}{(H^{1}(\Omega))^{2}}} be given as described above. Then $\psi[\rm]\in\Psi$ and we have

[TABLE]

Moreover, there is a constant $c>0$ such that

[TABLE]

Proof.

For all $v\in W$, we have

[TABLE]

since $\nabla v=0$ on $E$, which shows that $\psi[\rm]\in\Psi$. The decomposition (23) can easily be derived from the representation

[TABLE]

The sum is direct, since $\psi[q]=0$ for $q=0$. Finally, for showing (24), observe that

[TABLE]

where in the counterclockwise parametrization of $E$ we have $x_{A}=\gamma(\sigma_{E})$ and $x_{B}=\gamma(\sigma_{E}+|E|)$ and $|E|$ denotes the length of $E$. From this and (22) it easily follows that

[TABLE]

Here $H^{1}(\Gamma)$ denotes the space of functions with weak tangential derivative with corresponding norm $\|.\|_{1,\Gamma}$, see, e.g., [22]. The rest follows from standard trace and inverse trace inequalities. ∎

Remark 4.6*.*

The above construction relies on the assumption that we have at least one clamped edge, which we made at the beginning of this paper. However, this construction can be extended for general mixed boundary conditions using the freedom we have in the choice of $\psi_{\Gamma}[\rm]$ on clamped and simply supported edges and the compatibility conditions on $F$.

Lemma 4.7**.**

The space $\Psi_{0}$ is equal to the set of all functions $\psi\in(H^{1}(\Omega))^{2}$ that satisfy the following boundary conditions:

[TABLE]

with some constants $c_{E}\in\mathbb{R}$ for each edge $E\subset\Gamma_{s}$ and functions $r_{C}\in RT_{0}$ for each connected component $C\subset\Gamma_{f}$ satisfying the compatibility conditions $c_{E}=r_{C}(x)\operatorname{\cdot}n_{E}$, if $E\subset\Gamma_{s}$ is an adjacent edge to a component $C$ of $\Gamma_{f}$, where $x$ is the enclosed corner point and $n_{E}$ is the normal on $E$. Moreover, for each set of constants $c_{E}\in\mathbb{R}$ and functions $r_{C}\in RT_{0}$ that satisfy the compatibility conditions, there is a function $\psi\in(H^{1}(\Omega))^{2}$, for which the boundary conditions (25) hold.

Proof.

For $\psi\in(H^{1}(\Omega))^{2}$ and $v\in C^{\infty}(\overline{\Omega})\cap W$ we have

[TABLE]

using integration by parts and the boundary conditions for $v$. Now let $\psi\in\Psi_{0}$ and let $E$ be an edge from $\Gamma_{s}\cup\Gamma_{f}$. For each function $\varphi\in C_{0}^{\infty}(E)$, one can easily construct a function $v\in C_{0}^{\infty}(\Omega\cup E)$ with $v=0$ and $\partial_{n}v=\varphi$ on $E$. Since $\psi\in\Psi_{0}$, we have $\langle\partial_{t}\psi,\nabla v\rangle_{\Gamma}=0$, which reduces to

[TABLE]

by using (26). Hence, $\psi\operatorname{\cdot}n$ is equal to a constant $c_{E}$ on $E$. By a similar argument it follows that $\partial_{t}(\psi\operatorname{\cdot}t)$ is equal to a constant $a_{E}$ on each edge $E\in\Gamma_{f}$.

Using that $\psi\operatorname{\cdot}n$ is edgewise constant on {\color[rgb]{0,0,0}{\Gamma_{s}\cup\Gamma_{f}}} and $\partial_{t}(\psi\operatorname{\cdot}t)$ is edgewise constant on $\Gamma_{f}$, we obtain from (26) by edgewise integration by parts:

[TABLE]

Now let $E$ and $E^{\prime}$ be two adjacent edges from $\Gamma_{f}$ with the common corner point $x$ and let $v\in C_{0}^{\infty}(\Omega\cup E\cup E^{\prime}\cup\{x\})$. Then the condition $\langle\partial_{t}\psi,\nabla v\rangle_{\Gamma}=0$ reduces to

[TABLE]

Observe that $\psi\in P_{1}$ on each edge of $\Gamma_{f}$ and is the trace of an $H^{1}$-function. Therefore, $\psi$ must be continuous on $\Gamma_{f}$, which implies $\llbracket\psi\rrbracket_{x}=0$. Since $v(x)$ can be chosen arbitrarily, it follows that $\llbracket\partial_{t}(\psi\operatorname{\cdot}t)\rrbracket_{x}=0$. So $\partial_{t}(\psi\operatorname{\cdot}t)$ is not only constant on each edge $E$ of $\Gamma_{f}$ but it is equal to the same constant $a_{C}$ on each edge of a connected component $C$ of $\Gamma_{f}$. Since we additionally know from above that $\psi\operatorname{\cdot}n$ is constant on such an edge, it easily follows that $\partial_{t}\widetilde{\psi}=0$ on each edge of $C$ with $\widetilde{\psi}(x)=\psi(x)-a_{C}x$. Since $\psi$ is continuous on $\Gamma$, $\widetilde{\psi}$ is continuous, too. Then it follows that $\widetilde{\psi}$ is equal to a common constant $b_{C}\in\mathbb{R}^{2}$ on $C$. Hence, $\psi(x)=a_{C}x+b_{C}$ on $C$.

Let $E\subset\Gamma_{s}$ be an edge adjacent to a connected component $C$ of $\Gamma_{f}$ with enclosed corner point $x$. By a similar argument as above one can deduce from $\langle\partial_{t}\psi,\nabla v\rangle_{\Gamma}=0$ and (27) the compatibility condition $c_{E}=r_{C}(x)\operatorname{\cdot}n_{E}$.

For proving the reverse direction, let $\psi\in(H^{1}(\Omega))^{2}$ satisfy the boundary conditions (25). Using (27) we have

[TABLE]

for all $v\in C^{\infty}(\overline{\Omega})\cap W$ for the following reasons: $\partial_{n}v(x)=0$ on all interior corner points of $\Gamma_{s}$, $\psi\operatorname{\cdot}\nabla v-\partial_{t}(\psi\operatorname{\cdot}t)v$ is continuous on $\Gamma_{f}$, $v=0$ and $\nabla v=0$ on corner points on the interface of $\Gamma_{s}$ or $\Gamma_{f}$ with $\Gamma_{c}$, and the compatibility conditions on corner points on the interface of $\Gamma_{s}$ and $\Gamma_{f}$. So it follows that $\psi\in\Psi_{0}$.

Note, $\partial_{n}v(x)=0$ on all interior corner points of $\Gamma_{s}$, since $v=0$ on $\Gamma_{s}$. This makes the tangential derivative vanish at this corner points for two linear independent tangential directions. Therefore, we obtain $\nabla v(x)=0$ and $\partial_{n}v(x)=0$

For the last part, let a set of data $c_{E}\in\mathbb{R}$ for each edge $E\subset\Gamma_{s}$ and $r_{C}\in RT_{0}$ for each connected component $C\subset\Gamma_{f}$ be given, which satisfy the compatibility conditions. We will first construct a continuous function $\hat{\psi}_{\Gamma}$ on $\Gamma$ with $\psi_{\Gamma}|_{E}\in P_{1}$ for each edge $E\in\mathcal{E}_{\Gamma,}$, where $P_{1}$ is the set of polynomials of degree $\leq 1$, by prescribing its values on all corner points $x\in\mathcal{V}_{\Gamma,}$ as follows:

[TABLE]

Here $\psi_{s}(x)$ and $\psi_{cs}(x)$ are defined as follows: For a common corner point $x$ of $E,E^{\prime}\subset\Gamma_{s}$, $\psi_{s}(x)$ is uniquely given by

[TABLE]

For a common corner point $x$ of $E\subset\Gamma_{s}$ and $E^{\prime}\subset\Gamma_{c}$, $\psi_{sc}(x)$ is uniquely given by

[TABLE]

Let $\hat{\psi}\in(H^{1}(\Omega))^{2}$ be the harmonic extension of $\hat{\psi}_{\Gamma}$. It is easy to verify that the boundary conditions (25) are satisfied. ∎

In order to eliminate $c_{E}$ and $r_{C}$ in (25) we introduce the following Clément-type projection operator on $\Gamma$.

Definition 4.8**.**

Let $\psi_{\Gamma}\in L^{2}(\Gamma)$. Then the projection $\Pi_{\Gamma}\colon L^{2}(\Gamma)\rightarrow L^{2}(\Gamma)$ is given by $\hat{\psi}_{\Gamma}=\Pi_{\Gamma}\psi_{\Gamma}$, where $\hat{\psi}_{\Gamma}$ is constructed as in the last part of the proof of Lemma 4.7 from the data $c_{E}(\psi_{\Gamma})$ and $r_{C}(\psi_{\Gamma})$, which are the $L^{2}$-projection of $\psi_{\Gamma}\operatorname{\cdot}n_{E}$ onto the set $P_{0}$ of constant functions and the $L^{2}$-projection of $\psi_{\Gamma}|_{C}$ onto $RT_{0}$, respectively. Except if $E\subset\Gamma_{s}$ is an adjacent edge to a component $C$ of $\Gamma_{f}$, then $c_{E}(\psi_{\Gamma})=r_{C}(\psi_{\Gamma})(x)\operatorname{\cdot}n_{E}$, where $x$ is the enclosed corner point and $n_{E}$ is the normal on $E$, in order to enforce the compatibility conditions.

With this notation the boundary conditions (25) can be rewritten as

[TABLE]

with $\operatorname{P}{\psi}=(I-\Pi_{\Gamma})\psi$.

4.2 Decoupled formulation

Using the representations

[TABLE]

together with (23) leads to the following equivalent formulation of (12): Find $p\in{Q}$, $\phi\in\Psi_{p}=\psi[p]+\Psi_{0}$ and $w\in{Q}$ such that

[TABLE]

for all $\rm\in{Q}$, $\psi_{0}\in\Psi_{0}$ and $v\in{Q}$. For the representations of $\operatorname{div}\operatorname{Div}\bm{M}$ and $\operatorname{div}\operatorname{Div}\bm{N}$ recall Theorem 4.2, in particular for $\bm{M}$ we use identity (19) and for $\bm{N}$ we rely on the reverse direction. Therefore, the mixed formulation of the Kirchhoff plate bending problem is equivalent to three (consecutively to solve) elliptic second-order problems:

The $p$-problem: Find $p\in{Q}$ such that

[TABLE] 2. 2.

The $\phi$-problem: For given $p\in{Q}$, find $\phi\in\Psi_{p}=\psi[p]+\Psi_{0}$ such that

[TABLE] 3. 3.

The $w$-problem: For given $\bm{M}=p\bm{I}+\operatorname{symCurl}\phi$, find $w\in{Q}$ such that

[TABLE]

Remark 4.9*.*

For the rotated function $\phi^{\perp}=(-\phi_{2},\phi_{1})^{T}$ the $\phi$-problem becomes a linear elasticity problem

[TABLE]

with a appropriately rotated material tensor $\hat{\mathcal{C}}^{-1}$.

Remark 4.10*.*

If we only have clamped and simply supported boundary parts, we know from Section 4.1 that $\psi[\rm]=0$ on $\Omega$. Therefore, the terms involving $\psi[\rm]$ just vanish. Furthermore, in the purely clamped case $\Psi_{p}$ and $\Psi_{0}$ become $(H^{1}(\Omega))^{2}$.

Remark 4.11*.*

The $p$-problem and the $w$-problem are standard Poisson problems with mixed boundary conditions on $\Gamma_{c}\cup\Gamma_{s}$ and $\Gamma_{f}$. 2. 2.

From (30) we obtain

[TABLE]

therefore, $\mathcal{C}^{-1}\bm{M}\in\bm{H}(\operatorname{Rot},\Omega)=\{\bm{N}\in{\bm{L}^{2}(\Omega)}:\operatorname{Rot}\bm{N}\in(L^{2}(\Omega))^{2}\}$. With the representation $\bm{M}=p\bm{I}+\operatorname{symCurl}\phi$, using $p\in{Q}\subset H^{1}(\Omega)$, we receive that the solution of the $\phi$-problem satisfies the second-order differential equation

[TABLE]

the (essential) boundary conditions (see (28))

[TABLE]

and the following condition for the flux $\chi=\left(\mathcal{C}^{-1}\operatorname{symCurl}\phi\right)t$:

[TABLE]

If $\chi\in L^{2}(\Gamma)$, then (35) is equivalent to the following (natural) boundary conditions

[TABLE]

and

[TABLE]

where $p=0$ on $\Gamma_{c}\cup\Gamma_{s}$ is used (see the $p$-problem). 3. 3.

For the rotated function $\phi^{\perp}$, (33) is a linear elasticity problem

[TABLE]

with the corresponding boundary conditions for $\phi^{\perp}$ and \big{(}\hat{\mathcal{C}}^{-1}\varepsilon(\phi^{\perp})\big{)}n, which can be interpreted as displacement and traction.

Remark 4.12*.*

So far we have only considered homogeneous boundary conditions. In the following we indicate how to adapt the decoupled formulation introduced above to inhomogeneous boundary conditions of the form (cf. [27, 18]):

[TABLE]

and the corner forces

[TABLE]

In the $p$-problem (29) the additional contribution on the right-hand side is given as

[TABLE] 2. 2.

In the $\phi$-problem (30) the boundary value $\psi_{\Gamma}[p]$ needed in the construction of $\psi[p]$ has to be adapted. With the same notations as used for (22) we set

[TABLE]

where $c$ is edgewise constant, given by

[TABLE]

Here, the edges $E_{k}$ for $k=1,\dots,K$ are numbered consecutively in counterclockwise direction with $E_{1}$ starting from $x_{B}=\gamma(0)$. Furthermore, we denote the vertex at the end point of $\overline{E}_{k}$ by $x_{k}$ and use the convention $\hat{R}_{x_{k}}=0$ for corner points $x_{k}\not\in\mathcal{V}_{\Gamma,f}$.

As additional contribution on the right-hand side of the $\phi$-problem (30) we obtain

[TABLE]

provided $\partial_{t}\psi\in(L^{2}(\Omega))^{2}$. 3. 3.

In the $w$-problem (31) we get as additional contribution $-(\nabla\overline{w},\nabla\rm)$ on the right-hand side, where $\overline{w}\in H^{1}(\Omega)$ is any extension of the Dirichlet data $\hat{w}$.

5 The discretization method

Let $(\mathcal{T}_{h})_{h\in\mathcal{H}}$ be a shape-regular family of subdivisions of the domain $\Omega$ into polygonal elements. The diameter of an element $T\in\mathcal{T}_{h}$ is denoted by $h_{T}$ and we define $h=\max\{h_{T}:T\in\mathcal{T}_{h}\}$. We denote the set of all edges $e$ of elements $T\in\mathcal{T}_{h}$ with $e\subset\Gamma_{s}$ by $\mathcal{E}_{h,s}$ and with $e\subset\Gamma_{f}$ by $\mathcal{E}_{h,f}$. The length of an edge $e$ is denoted by $h_{e}$. Moreover, we assume that the total number of edges in $\mathcal{E}_{h,s}\cup\mathcal{E}_{h,f}$ is bounded by $c\,h^{-1}$. Here and in the sequel $c$ denotes a generic constant independent of $h$, possibly different at each occurrence. Note, this assumption is weaker than requiring a quasiuniform family of subdivisions. It can be viewed as a quasiuniformity condition in a neighborhood of the boundary.

Remark 5.1*.*

Observe that the symbol $e$ is used to indicate element edges, while the symbol $E$ used in the preceding sections is reserved for edges of the domain as introduced at the beginning of Section 2.

Let $\mathcal{S}_{h}$ be a finite dimensional subspace of $H^{1}(\Omega)$ of piecewise polynomials with degree $k$ associated with $\mathcal{T}_{h}$ and we set $\mathcal{S}_{h,0}=\mathcal{S}_{h}\cap H^{1}_{0,\Gamma_{c}\cup\Gamma_{s}}(\Omega)$. For $v\in H^{s}(\Omega)$, with $1\leq s\leq k+1$, we assume the standard approximation property

[TABLE]

for all $l\in\{0,\dots,s\}$. Moreover, we require the discrete trace inequality

[TABLE]

and the continuous trace inequality

[TABLE]

Here and in the following the $L^{2}$-norm on an element $T$ and an edge $e$ are denoted by $\|.\|_{0,T}$ and $\|.\|_{0,e}$, respectively. The properties (A1), (A2), (A3) are satisfied for standard finite element spaces or isogeometric B-spline discretization spaces under standard assumptions.

The method we propose consists of three consecutive steps:

The discrete $p$-problem: Find $p_{h}\in\mathcal{S}_{h,0}$ such that

[TABLE] 2. 2.

The discrete $\phi$-problem: For given $p_{h}\in\mathcal{S}_{h,0}$, find $\phi_{h}\in(\mathcal{S}_{h})^{2}/RT_{0}$ such that

[TABLE]

with

[TABLE]

where

[TABLE]

for some penalty parameter $\eta>0$. 3. 3.

The discrete $w$-problem: For given $\bm{M}_{h}=p_{h}\bm{I}+\operatorname{symCurl}\phi_{h}$, find $w_{h}\in\mathcal{S}_{h,0}$ such that

[TABLE]

with

[TABLE]

Remark 5.2*.*

1. (a)

The discrete $p$-problem is the standard Galerkin method applied to (29). 2. (b)

The discrete $\phi$-problem is a Nitsche method applied to (30), which is derived as follows. We start from the identity

[TABLE]

which follows from (32) by multiplying with a test function and using integration by parts.

Plugging in the second term of (42) the representation of $\bm{M}$ leads to

[TABLE]

with $\chi=(\mathcal{C}^{-1}\operatorname{symCurl}\phi)t$. For the next step note that $\psi-\operatorname{P}{\psi}=\Pi_{\Gamma}\psi$ on $\Gamma$ and according to Lemma 4.7 there exists a $\psi_{0}\in\Psi_{0}$ such that $\psi_{0}=\Pi_{\Gamma}\psi$ on $\Gamma$. Then the natural boundary conditions (36) and (37) (with $\psi$ replaced by $\psi_{0}$) lead to

[TABLE]

provided $\chi\in L^{2}(\Gamma)$. Then the method is obtained by first extending (42) by the terms $s(\psi,\phi-\psi_{\Gamma}[p])$ and $r_{h}(\phi-\psi_{\Gamma}[p],\psi)$, which vanish for the exact solution (see (34)), and then by replacing $p$, $\phi$, and $\psi$ with $p_{h}$, $\phi_{h}$, and $\psi_{h}$. Following [11] we call $s(\phi,\psi)+c(p,\psi)$, $s(\psi,\phi)$, and $r_{h}(\phi,\psi)$ the consistency, symmetry and penalty terms, respectively. 3. (c)

The discrete $w$-problem is the standard Galerkin method applied to (31), where the right-hand side is reformulated as follows. By using (42) and (44) with $\psi=\psi[\rm]$ we obtain for the right-hand side of (31):

[TABLE]

where we additionally extend by the term $r_{h}(\phi-\psi_{\Gamma}[p],\psi_{\Gamma}[\rm])$, which vanishes for the exact solution. Then (41) is obtained by replacing $p$, $\phi$, and with $p_{h}$, $\phi_{h}$, and h.

Remark 5.3*.*

Since only $\psi_{\Gamma}[p_{h}]$ and $\psi_{\Gamma}[\rm_{h}]$ appear in the numerical method, the extensions of $\psi[p_{h}]$ and $\psi[\rm_{h}]$ to the interior are not needed. Although the functions $\psi_{\Gamma}[\rm_{h}]$ do not have local support, the linear systems can still be assembled with optimal complexity.

Remark 5.4*.*

Our decomposition of the continuous problem also leads to a new interpretation of the well-known Hellan-Herrmann-Johnson (HHJ) method; see [15, 16, 19]. We can proceed similar as in Theorem 4.2 and derive a discrete regular decomposition for the approximation space of the auxiliary variable, leading to a reformulation of the HHJ method in form of three consecutively to solve discretized second-order problems, as it has already been worked out in details in the purely clamped situation in [20]. The HHJ method is mainly restricted to triangular meshes. In [28] a HHJ-type method on rectangular meshes is considered. The new method introduced above is more flexible in the sense that triangular and general quadrilateral meshes can be handled and also isogeometric B-spline discretization spaces can be used.

The main result of this section is the following a priori discretization error estimate for the proposed method.

Theorem 5.5**.**

Assume that $\Omega$ is convex. Let $(w,\bm{M})$, with $\bm{M}=p\bm{I}+\operatorname{symCurl}\phi$, be the solution of the mixed formulation $\eqref{eq:newMixed_formulation}$ and $(w_{h},\bm{M}_{h})$, with $\bm{M}_{h}=p_{h}\bm{I}+\operatorname{symCurl}\phi_{h}$, be the approximate solution, given by (38), (39), (40). For $w\in H^{s_{w}}(\Omega)$, $p\in H^{s_{p}}(\Omega)$ and $\phi\in(H^{s_{\phi}}(\Omega))^{2}$, with $1\leq s_{w},s_{p}\leq k+1$ and $\frac{3}{2}<s_{\phi}\leq k+1$, we have the estimate

[TABLE]

Especially, for $s_{p}=k+1$, $s_{\phi}=k+1$ and $s_{w}=k+1$ we obtain

[TABLE]

We will derive these estimates by discussing the discretization errors of the $p$-problem, the $\phi$-problem, and the $w$-problem consecutively.

5.1 Error estimates for the $p$-problem

We start with the $p$-problem (29) and its discretization (38). It is well-known that the following error estimates hold:

Lemma 5.6** ($p$-problem).**

Under the assumptions of Theorem 5.5 we have

[TABLE]

We refer to standard literature for the proof. Note that the convexity of $\Omega$ is used only for the $L^{2}$-estimates.

5.2 Error estimates for the $\phi$-problem

We follow the standard approach as outlined, e.g., in [11] and introduce two mesh-dependent semi-norms: the jump semi-norm $|\psi|_{\frac{1}{2},h}$ and the average semi-norm $|\chi|_{-\frac{1}{2},h}$, which here are given by

[TABLE]

The analysis relies on the discrete coercivity and the boundedness of the bilinear form $a_{\phi,h}$ in appropriate norms.

Lemma 5.7** (Coercivity).**

There is a constant $c>0$ such that

[TABLE]

provided $\eta$ is sufficiently large, where the mesh-dependent norm $\|\phi\|_{h}$ is given by

[TABLE]

Lemma 5.8** (Boundedness).**

There is a constant $c>0$ such that

[TABLE]

for all $\phi\in(H^{s}(\Omega))^{2}+(\mathcal{S}_{h})^{2}$, with $s>\frac{3}{2}$ and $\psi_{h}\in(\mathcal{S}_{h})^{2}$, where the mesh-dependent norm $\|\phi\|_{h,*}$ is given by

[TABLE]

The proofs of these two lemmas are analogous to the proofs of similar results in [11] and are, therefore, omitted. However, for later use, we explicitly mention here the fundamental estimates for the consistency, symmetry, and penalty terms which are used for proving coercivity and boundedness: For all $\psi\in(H^{1}(\Omega))^{2}$, $\xi\in(L^{2}(\Gamma))^{2}$, $q\in L^{2}(\Gamma)$ we have

[TABLE]

Obviously, (46) simplifies to

[TABLE]

Additionally we need an estimate of the consistency error.

Lemma 5.9** (Consistency error).**

Under the assumptions of Theorem 5.5 we have

[TABLE]

Proof.

From Remark 5.2 it follows that $\phi$ satisfies

[TABLE]

for all $\psi_{h}\in(\mathcal{S}_{h})^{2}/RT_{0}$. Therefore, the difference to the right-hand side in (39) is given as

[TABLE]

From the Cauchy inequality on $\Omega$ and (47), (49), (48) we obtain

[TABLE]

From the continuous trace inequality (A3) it follows that

[TABLE]

From the definition of $\Pi_{\Gamma}\xi$ and $\psi_{\Gamma}[\rm]$ one obtains

[TABLE]

for all $e\in\mathcal{E}_{h,s}\cup\mathcal{E}_{h,f}$, $\xi\in(L^{2}(\Gamma))^{2}$, and $\rm\in L^{2}(\Gamma)$.

For the first inequality note that

[TABLE]

where $c_{E}(\xi)\in\mathbb{R}$ and $r_{C}(\xi)(x)=a_{C}(\xi)+b_{C}(\xi)\,x$ with $a_{C}(\xi)\in\mathbb{R}^{2}$ and $b_{C}(\xi)\in\mathbb{R}$ are the data used in the construction of $\Pi_{\Gamma}\xi$ in Definition 4.8. On an edge $E\subset\Gamma_{s}$ we have

[TABLE]

By similar arguments analogous results hold for $|a_{C}(\xi)|$ and $|b_{C}(\xi)|$. Combining these estimates with

[TABLE]

provides the first inequality. The second inequality holds since

[TABLE]

Using that the total number of edges in $\mathcal{E}_{h,s}\cup\mathcal{E}_{h,f}$ is bounded by $c\,h^{-1}$ it follows that

[TABLE]

and, therefore,

[TABLE]

Then the estimate immediately follows from the error estimates in Lemma 5.6. ∎

From the last three lemmas we obtain the following error estimate for the $\phi$-problem:

Lemma 5.10** ($\phi$-problem).**

Under the assumptions of Theorem 5.5 we have

[TABLE]

Proof.

From coercivity and boundedness we obtain by standard arguments:

[TABLE]

Since $\phi\in(H^{s}(\Omega))^{2}$, with $\frac{3}{2}<s\leq k+1$, it follows from assumptions (A1) and (A3) that

[TABLE]

This approximation property and Lemma 5.9 directly imply the error estimate. ∎

5.3 Error estimates for the $w$-problem

Lemma 5.11** ($w$-problem).**

Under the assumptions of Theorem 5.5 we have

[TABLE]

Proof.

The first lemma of Strang provides

[TABLE]

The first term can be estimated by the approximation property (A1). It remains to estimate the consistency error. Using (31) with (45) and (41) we have

[TABLE]

Subtracting (42) and (39) leads to

[TABLE]

By subtracting the last two equations we obtain

[TABLE]

The five terms on the right-hand side, denoted by $T_{1}$, $T_{2}$, …$T_{5}$ in consecutive order of their appearance, are estimated as follows: From the Cauchy inequality on $\Omega$ and (46), (47), (49), (48) we obtain

[TABLE]

for all $\psi_{h}\in(\mathcal{S}_{h})^{2}$. In particular, we choose $\psi_{h}=\Pi_{h}(\psi[\rm_{h}])$, where $\Pi_{h}$ denotes the $L^{2}$-orthogonal projection onto $(\mathcal{S}_{h})^{2}$. Then, following , e.g., [11], it can be shown that

[TABLE]

With the stability estimate (24) it follows that

[TABLE]

Observe that

[TABLE]

and

[TABLE]

see the estimates in the proof of Lemma 5.10. Then the result follows directly from the estimates in Lemma 5.6, Lemma 5.10, and (A3*). ∎

Finally, we obtain the proof of the main result:

Proof of Theorem 5.5.

By combining the results of Lemma 5.6 and Lemma 5.10 we obtain

[TABLE]

Together with Lemma 5.11 this completes the proof. ∎

6 Numerical experiments

We consider a square plate $\Omega=(-1,1)^{2}$ with simply supported north and south boundary, clamped west boundary and free east boundary. The material tensor $\mathcal{C}$ is the identity, and the load is given by

[TABLE]

The exact solution is of the form

[TABLE]

which satisfies the boundary conditions on the simply supported boundary parts anyway. The constants $a,b,c$ and $d$ are chosen such that the four remaining boundary conditions (on the clamped and free boundary parts) are fulfilled, for details, see [27].

In order to illustrate the flexibility of our discretization method we use for $\mathcal{S}_{h}$ isogeometric B-spline discretization spaces of degree $k$ with maximum smoothness; see, e.g, [9, 10] for information on isogeometric analysis. For $k=1$, the discretization space $\mathcal{S}_{h}$ coincides with the standard finite element space of continuous and piecewise bilinear elements. In all experiments a sparse direct solver is used for each of the three sub-problems. The implementation is done in the framework of G+Smo ("Geometry + Simulation Modules"), an object-oriented C++ library, see https://ricamsvn.ricam.oeaw.ac.at/trac/gismo/wiki/WikiStart.

In Tables 1, 2, 3 the discretization errors for $k=1,2,3$ are shown. The first column shows the refinement level $L$, i.e. the number of uniform $h$-refinements of $\Omega$. The column "order" contains the error reduction relative to the previous level. The experiments show optimal convergence rates for $w$ and $\bm{M}$ as predicted by the analysis. For the columns containing the errors for $p$ and $\phi$, the (analytically not available) exact solutions $p$ and $\phi$ are replaced by their numerical solutions on level $L=9$. Note that also for $p$ and $\phi$ optimal convergence rates are observed.