# An augmented Lagrangian method for the Signorini boundary value problem with BEM

Shougui Zhang, Xiaolin Li

► To cite this version:

Shougui Zhang, Xiaolin Li. An augmented Lagrangian method for the Signorini boundary value problem with BEM. Boundary Value Problems, 2016, 2016, pp.62. 10.1186/s13661-016-0570-2. hal-01308846

HAL Id: hal-01308846

https://hal.science/hal-01308846

Submitted on 28 Apr 2016

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# An augmented Lagrangian method for the Signorini boundary value problem with BEM

Shougui Zhang* and Xiaolin Li

**Abstract**

We analyze augmented Lagrangian and boundary element methods for the Signorini boundary value problem of Laplacian. The boundary variational formulation is presented by the boundary integral operators, and the Signorini boundary conditions are formulated as a fixed point problem. Semismooth Newton methods are applied for the numerical solution of the problem. We prove the convergence of the method and confirm the theory by some numerical experiments.

**Keywords:** Signorini boundary conditions; fixed point; augmented Lagrangian; boundary element; semismooth Newton methods

## 1 Introduction

Signorini boundary value problems are of great importance in many applications, such as the electropaint process and contact problem [1–3], etc. These problems involve inequality constraints on a part of the boundary that make them nonlinear. Therefore, development of efficient numerical algorithms for Signorini boundary value problems is one of the most important branches of modern computational mathematics and mechanics. Usually Signorini boundary value problems have been considered mathematically and numerically with variational inequalities, especially by the finite element method (FEM) [1–7] and the boundary element method (BEM) [8–13]. Recently, sound and efficient algorithms for the solution of these problems were still a very active field of research (see [14–20]).

As we know, the fixed point method based on projection theory is a powerful tool to deal with complementary problems and variational inequalities in finite dimensional space [21]. The main idea of this method is to establish the equivalence between the original problem and the fixed point problem by using projection. This equivalent formulation plays a significant role in developing various iterative methods for solving an original problem. During the last 20 years, a number of projection methods have been studied extensively [4, 5, 22, 23], which are perfectly efficient for solving the problem. In these methods, the problem has been formulated by a projection algorithm, and no other inequality constraint is needed. Consequently, the method is easy to implement for the numerical solution and the proof of convergence is very simple.

In the case of Signorini boundary conditions, the unknown boundary values are the potential and its derivative on the boundary, which are considered primary variables in BEM. They can be obtained directly using boundary integral equations [24, 25]. Therefore, the method combining the projection method with the BEM is more appropriate for

<footnote>*Correspondence:
shgzhang@cqnu.edu.cn
School of Mathematical Sciences,
Chongqing Normal University,
Chongqing, 401331, P.R. China</footnote>
Signorini boundary value problems [18–20]. However, few investigations have been done on the Signorini boundary value problem by the augmented Lagrangian and fixed point methods with the BEM.

In this paper, we focus on the boundary augmented Lagrangian methods (BALM) for the solution of Signorini boundary value problems, which is inspired by classical augmented Lagrangian methods (ALM) [7] and the BEM [11]. First of all, we deduce the boundary weak formulation with Steklov-Poincaré operator. Second, we use the projection technique to deal with the Signorini boundary conditions by an equality and the projection. Although the new problem is still nonlinear on the boundary, this problem no longer has the inequality constraint and can be solved by the semismooth Newton method with local superlinear convergence rate [7, 14]. Using these transformations, we propose a BALM for the Signorini boundary value problem, which only needs the iteration for boundary function and the computation of the boundary variational equation. Then we can use properties of the projection and the boundary integral operator to prove the convergence of the method. Numerical results show that our method is accurate and efficient.

The structure of the paper is as follows. In Section 2, we first describe the classical Signorini boundary value problem and use the Steklov-Poincaré operator to introduce the variational formulation. Then we establish equivalent formulations between the non-linear boundary conditions and the fixed point problem, and propose a new ALM for the problem. Section 3 is devoted to the convergence analysis of the method, which shows monotone convergence properties of the numerical solution toward the solution of the original problem. In Section 4, we present some numerical experiments to investigate the performance of our method, and finally a brief conclusion is given in Section 5.

## 2 The weak formulations for the Signorini boundary-value problem

Let $\Omega \subset \mathbb{R}^2$ be a bounded domain with a Lipschitz boundary $\Gamma$ and outward unit normal $n$. This boundary consists of three disjoint parts, the Dirichlet boundary $\Gamma_D$, the Neumann boundary $\Gamma_N$ and the Signorini boundary $\Gamma_S \neq \emptyset$. For simplicity, we consider the Signorini boundary value problem for the Laplace equation: find $u \in H^1(\Omega)$ such that

$$ \Delta u = 0 \quad \text{in } \Omega, \tag{2.1} $$

$$ u = g \quad \text{on } \Gamma_D, \tag{2.2} $$

$$ \lambda = f \quad \text{on } \Gamma_N, \tag{2.3} $$

$$ u \geq g, \quad \lambda \geq f, \quad (u-g)(\lambda-f) = 0 \quad \text{on } \Gamma_S, \tag{2.4} $$

where $\lambda := \frac{\partial u}{\partial n}$, $g \in H^{1/2}(\Gamma \setminus \bar{\Gamma}_N)$, and $f \in H^{-1/2}(\Gamma \setminus \bar{\Gamma}_D)$ are given. We note that $H^{-1/2}(\Gamma \setminus \bar{\Gamma}_D)$ is defined by the dual of $H_0^{1/2}(\Gamma \setminus \bar{\Gamma}_D) := \{u \in H^{1/2}(\Gamma \setminus \bar{\Gamma}_D) | \text{supp } u \subset \bar{\Omega}\}$ with the norm

$$ \|f\|_{H^{-1/2}(\Gamma \setminus \bar{\Gamma}_D)} := \sup_{0 \neq \nu \in H^{1/2}(\Gamma \setminus \bar{\Gamma}_D)} \frac{|\langle f, \nu \rangle_{\Gamma \setminus \bar{\Gamma}_D}|}{\|\nu\|_{H^{1/2}(\Gamma \setminus \bar{\Gamma}_D)}}. $$

Here $\langle \cdot, \cdot \rangle_{\Gamma \setminus \bar{\Gamma}_D}$ denotes the extension of the usual $L^2(\Gamma \setminus \bar{\Gamma}_D)$ scalar product to $H^{-1/2}(\Gamma \setminus \bar{\Gamma}_D) \times H^{1/2}(\Gamma \setminus \bar{\Gamma}_D)$. It is well known in the theory of variational inequalities that this problem has a unique solution if $\Gamma_S \neq \emptyset$ or $\int_{(\Gamma \setminus \bar{\Gamma}_D)} f \, ds < 0$ [1, 8, 9].
We define the space of functions as

$$H_D^1(\Omega) := \{v \in H^1(\Omega), v = g \text{ on } \Gamma_D\}.$$

Applying the Green's formula, we obtain the variational formulation from (2.1)-(2.4) as follows:

$$\int_{\Omega} \nabla u \nabla v dx = \int_{\Gamma_N \cup \Gamma_S} \lambda v ds, \quad \forall v \in H_D^1(\Omega). \qquad (2.5)$$

Further, we introduce the Dirichlet-to-Neumann map on $\Gamma$ [11, 25, 26]

$$S: H^{1/2}(\Gamma) \to H^{-1/2}(\Gamma) \\ u|_{\Gamma} \mapsto \lambda|_{\Gamma}.$$

Inserting $S(u|_{\Gamma}) = \lambda|_{\Gamma}$ into (2.5) yields

$$\int_{\Omega} \nabla u \nabla v dx = \int_{\Gamma} S(u|_{\Gamma}) v ds, \quad \forall v \in H_D^1(\Omega),$$

where the Steklov-Poincaré operator $S$ is defined by

$$(Su)(x) = \left[ D + \left( \frac{1}{2}I + K' \right) V^{-1} \left( \frac{1}{2}I + K \right) \right] u(x) \qquad (2.6)$$

with the boundary integral operators

$$(V\lambda)(x) = \int_{\Gamma} U(x, y)\lambda(y) ds_y, \quad V: H^{-1/2}(\Gamma) \to H^{1/2}(\Gamma),$$

$$(Ku)(x) = \int_{\Gamma} \frac{\partial}{\partial n_y} U(x, y)u(y) ds_y, \quad K: H^{1/2}(\Gamma) \to H^{1/2}(\Gamma),$$

$$(K'\lambda)(x) = \int_{\Gamma} \frac{\partial}{\partial n_x} U(x, y)\lambda(y) ds_y, \quad K': H^{-1/2}(\Gamma) \to H^{-1/2}(\Gamma),$$

$$(Du)(x) = -\frac{\partial}{\partial n_x} \int_{\Gamma} \frac{\partial}{\partial n_y} U(x, y)u(y) ds_y, \quad D: H^{1/2}(\Gamma) \to H^{-1/2}(\Gamma).$$

Here the function $U$ is the fundamental solution of the two-dimensional Laplace equation:

$$U(x, y) = -\frac{1}{2\pi} \ln |x - y|.$$

Let us define

$$H_D^{1/2}(\Gamma) := \{v \in H^{1/2}(\Gamma), v = g \text{ on } \Gamma_D\},$$

$$H_0^{1/2}(\Gamma, \Gamma_D) := \{v \in H^{1/2}(\Gamma), v = 0 \text{ on } \Gamma_D\},$$

$$\langle Su, v \rangle_{\Gamma} := \int_{\Gamma} Su(x)v(x) ds_x,$$

$$\langle \lambda, v \rangle_{\Gamma_S} := \int_{\Gamma_S} \lambda(x)v(x) ds_x,$$
$$L(v) := \int_{\Gamma_N} f(x)v(x) ds_x.$$

We then can obtain the boundary weak formulation of the problem (2.1)-(2.4) as follows:

$$\langle Su, v \rangle_{\Gamma} - \langle \lambda, v \rangle_{\Gamma_S} = L(v), \quad \forall v \in H_D^{1/2}(\Gamma). \qquad (2.7)$$

It can be proved that the Steklov-Poincaré operator $S$ is bounded and symmetric, and semielliptic on $H^{1/2}(\Gamma)$. We will use the following property for the operator $S$ [25].

**Lemma 2.1** The Steklov-Poincaré operator $S$ defined by (2.6) is elliptic on $H_0^{1/2}(\Gamma, \Gamma_D)$,

$$\langle Sv, v \rangle_{L^2(\Gamma)} \geq \alpha \|v\|_{H^{1/2}(\Gamma)}^2, \qquad (2.8)$$

where $\alpha$ is a positive constant.

On the other hand, we can transfer the nonlinear boundary conditions (2.4) to a fixed point problem [4, 5, 18, 19, 22, 23]. Let us introduce the projection notation for $v \in \mathbb{R}$:

$$[v]_+ = \begin{cases} v & \text{if } v > 0, \\ 0 & \text{otherwise.} \end{cases}$$

Consequently, we obtain the following result [7, 20].

**Lemma 2.2** For all $\rho > 0$, the boundary conditions (2.4) on $\Gamma_S$ are equivalent to:

$$u - g - [u - g - \rho(\lambda - f)]_+ = 0 \quad \text{on } \Gamma_S. \qquad (2.9)$$

*Proof* Let $u$ be such that (2.4) holds. Considering the condition $u \ge g$ we have either $u > g$ or $u = g$. First, suppose that $u > g$. Then the condition $(u - g)(\lambda - f) = 0$ implies that $\lambda = f$. Therefore, we have

$$[u - g - \rho(\lambda - f)]_+ = [u - g]_+ = u - g.$$

Then suppose that $u = g$. The condition $\lambda \ge f$ can also be rewritten as $[-\rho(\lambda - f)]_+ = 0$, so it follows that

$$[u - g - \rho(\lambda - f)]_+ = [-\rho(\lambda - f)]_+ = 0 = u - g.$$

On the other hand, let $u$ and $\lambda$ be such that (2.9) holds. Note first that it implies $u \ge g$. Consider first the case $u = g$. Equation (2.9) can be rewritten as $[-\rho(\lambda - f)]_+ = 0$, which is equivalent to the condition

$$\lambda \ge f.$$

Since $u = g$, the following condition also holds:

$$(u - g)(\lambda - f) = 0.$$
We now consider the case $u > g$. From (2.9), $[u - g - \rho(\lambda - f)]_+ > 0$, so that in this case
$u - g = [u - g - \rho(\lambda - f)]_+ = u - g - \rho(\lambda - f)$, from which we have

$\lambda = f,$

so that all conditions in (2.4) hold.

<figure><img src="image_3.png" /></figure>

Now, we obtain the boundary weak formulation (2.7) by the Steklov-Poincaré operator $S$
and the fixed point problem (2.9) for the problem (2.1)-(2.4), which only involves a bound-
ary integral operator and no inequality constraint. This alternative equivalent formulation
is also useful from the numerical and theoretical analysis point of view. With the above
preparations, we can now list our boundary augmented Lagrangian method (BALM) for
the Signorini boundary value problem as below.

**Algorithm BALM**

Step 0: Choose $u^{(0)} \in L^2(\Gamma_S)$, $\rho \in R^+$, and set $k := 0$.

Step 1: Solve

$$
\langle S u^{(k+1)}, v \rangle_{\Gamma} + \langle \lambda^{(k+1)}, v \rangle_{\Gamma_S} = L(v), \quad \forall v \in H_D^{1/2}(\Gamma), \qquad (2.10)
$$

with boundary condition

$$
u^{(k+1)} - g - [u^{(k)} - g - \rho(\lambda^{(k+1)} - f)]_+ = 0 \quad \text{on } \Gamma_S, \qquad (2.11)
$$

and obtain $u^{(k+1)}$ and $\lambda^{(k+1)}$ on $\Gamma_S$.

Step 2: Update (2.11) by $k := k + 1$ and return to Step 1.

For the nonlinear problem (2.11), we can apply the semismooth Newton method for its
solution [7, 14]. In the next section the convergence of the algorithm is analyzed.

**3 Convergence of the algorithm**

Let $u^*$ and $\lambda^*$ denote the solution of Signorini boundary value problem and the corre-
sponding derivative on the boundary $\Gamma$, respectively. In order to analyze the convergence
of the method in Section 2, we suppose that the sequences $u^{(k)} \in H^{1/2}(\Gamma)$ and $\lambda^{(k)} \in L^2(\Gamma_S)$.
Now, we define

$$
B_{\rho}(u, \lambda) := u - g - \rho(\lambda - f),
$$

and we have the following projection property on $\Gamma_S$ [7, 22, 23].

**Lemma 3.1** For all $u^{(k)} \in H^{1/2}(\Gamma)$ and $\lambda^{(k)} \in L^2(\Gamma_S)$ generated by the BALM,

$$
\langle u^{(k+1)} - u^*, B_\rho(u^{(k)}, \lambda^{(k+1)}) - B_\rho(u^*, \lambda^*) \rangle_{\Gamma_S} \ge \|u^{(k+1)} - u^*\|_{\Gamma_S}^2. \quad (3.1)
$$

*Proof* We separate $\Gamma_S$ into the four subparts $\Gamma_{S1}$, $\Gamma_{S2}$, $\Gamma_{S3}$, and $\Gamma_{S4}$, where

$$
B_{\rho}(u^{(k)}, \lambda^{(k+1)}) \geq 0, B_{\rho}(u^{*}, \lambda^{*}) \geq 0 \quad \text{on } \Gamma_{S1},
$$
$$B_{\rho}(u^{(k)}, \lambda^{(k+1)}) \geq 0, B_{\rho}(u^{*}, \lambda^{*}) < 0 \quad \text{on } \Gamma_{S2},$$

$$B_{\rho}(u^{(k)}, \lambda^{(k+1)}) < 0, B_{\rho}(u^{*}, \lambda^{*}) \geq 0 \quad \text{on } \Gamma_{S3},$$

$$B_{\rho}(u^{(k)}, \lambda^{(k+1)}) < 0, B_{\rho}(u^{*}, \lambda^{*}) < 0 \quad \text{on } \Gamma_{S4}.$$

Using (2.9) and (2.11) we then get

$$u^{(k+1)} - u^{*} = B_{\rho}(u^{(k)}, \lambda^{(k+1)}) - B_{\rho}(u^{*}, \lambda^{*}) \quad \text{on } \Gamma_{S1},$$

$$0 \leq u^{(k+1)} - u^{*} = B_{\rho}(u^{(k)}, \lambda^{(k+1)}) - 0 < B_{\rho}(u^{(k)}, \lambda^{(k+1)}) - B_{\rho}(u^{*}, \lambda^{*}) \quad \text{on } \Gamma_{S2},$$

$$0 \geq u^{(k+1)} - u^{*} = 0 - B_{\rho}(u^{*}, \lambda^{*}) > B_{\rho}(u^{(k)}, \lambda^{(k+1)}) - B_{\rho}(u^{*}, \lambda^{*}) \quad \text{on } \Gamma_{S3},$$

$$u^{(k+1)} - u^{*} = 0 - 0 = 0 \quad \text{on } \Gamma_{S4}.$$

Consequently,

$$\langle u^{(k+1)} - u^{*}, B_{\rho}(u^{(k)}, \lambda^{(k+1)}) - B_{\rho}(u^{*}, \lambda^{*}) \rangle_{\Gamma_S} \geq \|u^{(k+1)} - u^{*}\|_{\Gamma_S}^2. \quad \square$$

**Theorem 3.1** Let $\{(u^{(k)}, \lambda^{(k)})\}$ be the sequence generated by the BALM, for all $k$ then $u^{(k)}$ converges to $u^{*}$ in $H^{1/2}(\Gamma)$ and $\lambda^{(k)}$ converges to $\lambda^{*}$ in $L^2(\Gamma_S)$ as $k \to \infty$.

*Proof* Note $\delta_u^{(k)} := u^{(k)} - u^* \in H_0^{1/2}(\Gamma, \Gamma_D)$ and $\delta_\lambda^{(k)} := \lambda^{(k)} - \lambda^* \in L^2(\Gamma_S)$. Considering that $(u^*, \lambda^*)$ satisfies (2.7) we obtain

$$\langle Su^*, \delta_u^{(k+1)} \rangle_{\Gamma} - \langle \lambda^*, \delta_u^{(k+1)} \rangle_{\Gamma_S} = L(\delta_u^{(k+1)}). \quad (3.2)$$

From (2.10) of BALM we have

$$\langle Su^{(k+1)}, \delta_u^{(k+1)} \rangle_{\Gamma} - \langle \lambda^{(k+1)}, \delta_u^{(k+1)} \rangle_{\Gamma_S} = L(\delta_u^{(k+1)}). \quad (3.3)$$

Subtracting (3.2) from (3.3) yields

$$\langle S\delta_u^{(k+1)}, \delta_u^{(k+1)} \rangle_{\Gamma} = \langle \delta_\lambda^{(k+1)}, \delta_u^{(k+1)} \rangle_{\Gamma_S}. \quad (3.4)$$

Using Lemma 3.1 and Young's inequality

$$\langle \delta_\lambda^{(k)}, \delta_u^{(k+1)} \rangle_{\Gamma_S} \leq \frac{1}{2} (\| \delta_\lambda^{(k)} \|_{\Gamma_S}^2 + \| \delta_u^{(k+1)} \|_{\Gamma_S}^2),$$

we obtain

$$\begin{aligned} \langle \delta_\lambda^{(k+1)}, \delta_u^{(k+1)} \rangle_{\Gamma_S} &= \rho^{-1} \langle \delta_u^{(k)}, \delta_u^{(k+1)} \rangle_{\Gamma_S} - \rho^{-1} \langle B_\rho(u^{(k)}, \lambda^{(k+1)}), B_\rho(u^{*}, \lambda^{*}) \rangle_{\Gamma_S} \\ &\leq \rho^{-1} \langle \delta_u^{(k)}, \delta_u^{(k+1)} \rangle_{\Gamma_S} - \rho^{-1} \| \delta_u^{(k+1)} \|_{\Gamma_S}^2 \\ &\leq (2\rho)^{-1} \| \delta_u^{(k)} \|_{\Gamma_S}^2 - (2\rho)^{-1} \| \delta_u^{(k+1)} \|_{\Gamma_S}^2. \end{aligned}$$

From (3.4) and Lemma 2.1 we then have

$$\langle S\delta_u^{(k+1)}, \delta_u^{(k+1)} \rangle_{\Gamma} \leq (2\rho)^{-1} \|\delta_u^{(k)}\|_{\Gamma_S}^2 - (2\rho)^{-1} \|\delta_u^{(k+1)}\|_{\Gamma_S}^2 \quad (3.5)$$
and

$$ \langle S\delta_u^{(k+1)}, \delta_u^{(k+1)} \rangle_{\Gamma} \geq \alpha \|\delta_u^{(k+1)}\|_{H^{1/2}(\Gamma)}^2 \quad (3.6) $$

It follows from (3.5) and (3.6) that

$$ \alpha \|\delta_u^{(k+1)}\|_{H^{1/2}(\Gamma)}^2 \leq (2\rho)^{-1} \|\delta_u^{(k)}\|_{\Gamma_S}^2 - (2\rho)^{-1} \|\delta_u^{(k+1)}\|_{\Gamma_S}^2 \quad (3.7) $$

Then

$$ \sum_{k=0}^{\infty} \alpha \|\delta_u^{(k+1)}\|_{H^{1/2}(\Gamma)}^2 \leq (2\rho)^{-1} \|\delta_u^{(0)}\|_{\Gamma_S}^2 < \infty. $$

Consequently,

$$ \lim_{k \to \infty} \|\delta_u^{(k+1)}\|_{H^{1/2}(\Gamma)}^2 = 0. $$

Thus $u^{(k)}$ converges to $u^*$ in $H^{1/2}(\Gamma)$ and from (2.11) of BALM $\lambda^{(k)}$ converges to $\lambda^*$ in $L^2(\Gamma_S)$ as $k \to \infty$. $\square$

According to (3.7), the sequence {$u^{(k)}$} is bounded and the sequence {$\|\delta_u^{(k)}\|_{\Gamma_S}$} is monotonically decreasing. Furthermore, we can see that the larger value of the parameter $\rho$ results in faster convergence to the algorithm.

# 4 Numerical experiments

To test the numerical verification of the theory, the algorithm above has been implemented and applied to some examples of Signorini problems in this section. An analytic solution is known for the first example, and the exact solution for the two other examples is not known. For the sake of simplicity, we apply a constant BEM to the problem (2.10) with iterations. Here, we choose $\|u^{(k+1)} - u^{(k)}\|_{\infty, \Gamma_S} \leq 10^{-10} \|u^{(k+1)}\|_{\infty, \Gamma_S}$ as a stopping criterion.

## 4.1 Dirichlet-Signorini problem

First we consider a Signorini boundary value problem in an annular domain $\Omega = \{(x, y) : a < \sqrt{x^2 + y^2} < b\}$ with a Dirichlet boundary condition on the boundary $\Gamma_D = \{(x, y) : \sqrt{x^2 + y^2} = b\} \cup \{(x, y) : \sqrt{x^2 + y^2} = a, y \geq 0\}$ and a Signorini boundary condition on the $\Gamma_S = \Gamma \setminus \bar{\Gamma}_D$. The analytic solution in the domain $\Omega$ is given by the following complex function:

$$ u(x, y) = \operatorname{Im} \omega^3(x + iy), $$

where

$$ \omega(x + iy) = \sqrt{ \frac{1}{2} \left( \frac{x^2 - y^2}{r^2} \right)^2 + \frac{1}{4} \left( \frac{r^2}{a^2} - \frac{a^2}{r^2} \right)^2 + \frac{1}{4} \frac{x^2 - y^2}{r^2} \left( \frac{r^2}{a^2} + \frac{a^2}{r^2} \right)^2 } \operatorname{sign} x \\ + i \sqrt{ \frac{1}{2} \left( \frac{x^2 - y^2}{r^2} \right)^2 + \frac{1}{4} \left( \frac{r^2}{a^2} - \frac{a^2}{r^2} \right)^2 - \frac{1}{4} \frac{x^2 - y^2}{r^2} \left( \frac{r^2}{a^2} + \frac{a^2}{r^2} \right)^2 } \operatorname{sign} y, $$
and $r = \sqrt{x^2 + y^2}, \sqrt{x^2 + y^2} \ge a$. For this problem, we can easily obtain the Dirichlet boundary condition on $\Gamma_D$ from the analytic solution, and the Signorini boundary condition is given by

$$u \ge 0, \quad \lambda \ge 0, \quad u\lambda = 0 \quad \text{on } \Gamma_S.$$

The analytic solution of this problem and its normal derivative on the Signorini boundary $\Gamma_S$ are

$$u(x,y) = \sqrt{\max\left(0, \frac{x^2 - y^2}{a^2}\right)^3} \operatorname{sign} y, \qquad (4.1)$$

$$\lambda(x,y) = -\frac{6}{a^3} \sqrt{\max\left(0, \frac{x^2 - y^2}{a^2}\right)} |x| y. \qquad (4.2)$$

Using the BEM, this problem has been considered by the decomposition-coordination method [9], the projection iterative algorithm [19], and the linear complementary method [20], etc.

We consider the cases $a = 0.1$ and $b = 0.25$ and introduce the parameterizations $t \to (a \cos \pi t, -a \sin \pi t)$ and $t \to (b \cos \pi t, b \sin \pi t)$. First we choose $\rho = 100$, and apply our method to this problem on a uniform grid for $t$. The discretization includes 40 boundary elements on $\Gamma_S$ and 120 boundary elements on $\Gamma_D$, so $N = 160$. Figure 1 plots the approximations and exact solutions for the potential $u$ on $\Gamma_S$. The results for the normal derivative $\lambda$ are shown in Figure 2. These figures show that our results are in excellent agreement with the exact solution (4.1) and (4.2).

To verify the convergence of our method, we solve the problem by choosing different numbers of boundary elements $N$ and various values of the parameter $\rho$. Table 1 shows the number of iterations with $N = 40, 80, 160, 320$, and 640 on $\Gamma$ and $\rho = 10^0, 10^1, 10^2, 10^3, 10^4$, and $10^5$. We can observe that the numerical results show good convergence as the parameter $\rho$ increases. Moreover, the number of iterations increases

<figure><img src="image_4.png" /><figcaption>Figure 1 Analytic and approximate solutions for $u$ on $\Gamma_S$.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>Table 1 Number of iterations with different values of N and $\rho$</figcaption></figure>

<table><caption>Figure 3 Log-log plot of the convergence for approximate solutions u on $\Gamma_5$.</caption><thead><tr><th></th><th>N = 40</th><th>N = 80</th><th>N = 160</th><th>N = 320</th><th>N = 640</th></tr></thead><tbody><tr><td>ρ = 10<sup>0</sup></td><td>9</td><td>9</td><td>10</td><td>10</td><td>12</td></tr><tr><td>ρ = 10<sup>1</sup></td><td>6</td><td>6</td><td>7</td><td>8</td><td>8</td></tr><tr><td>ρ = 10<sup>2</sup></td><td>5</td><td>5</td><td>6</td><td>7</td><td>8</td></tr><tr><td>ρ = 10<sup>3</sup></td><td>4</td><td>4</td><td>5</td><td>6</td><td>7</td></tr><tr><td>ρ = 10<sup>4</sup></td><td>3</td><td>4</td><td>5</td><td>6</td><td>7</td></tr><tr><td>ρ = 10<sup>5</sup></td><td>3</td><td>4</td><td>5</td><td>6</td><td>7</td></tr></tbody></table>

<figure><img src="image_6.png" /></figure>

slowly as N increases. We define the error

$$e(u) = \frac{1}{N_S} \sum_{i=1}^{N_S} |u(x_i) - u_h(x_i)|,$$

where $N_S$ denotes the number of boundary elements on $\Gamma_5$, and $u(x_i)$ and $u_h(x_i)$ denote exact and numerical solutions, respectively. Here we draw the error with the logarithmic scale depending on the step $h$. Figures 3 and 4 display the change trend of the error for $u$ and $\lambda$. We can see that our method yields very accurate results and converges superlinearly.
<figure><img src="image_7.png" /><figcaption>Figure 4 Log-log plot of the convergence for approximate solutions $\lambda$ on $\Gamma_S$.</figcaption></figure>

## 4.2 Dirichlet-Signorini problem

Second, we consider the following Signorini boundary value problem known as the steady-state electropainting problem:

$$ 
\begin{array}{ll}
\Delta u = 0 & \text{in } \Omega = (-0.5, 0.5) \times (0, 1), \\
u = 1 & \text{on } \Gamma_D = \{(x, y) : -0.5 \le x \le 0.5, y = 0\},
\end{array}
$$

with the following Signorini boundary conditions

$$ u \ge 0, \quad \lambda \ge -\varepsilon, \quad u(\lambda + \varepsilon) = 0 \quad \text{on } \Gamma_S, $$

where $\Gamma_S = \Gamma \setminus \bar{\Gamma}_D$. This problem has been extensively considered by different methods [2, 13, 16, 20].

In this problem the Signorini boundary conditions describe the location of these painted and unpainted parts on $\Gamma_S$, and the solution of the problem depends on the value of $\varepsilon$. We now apply our method to this problem, and four cases with different values of $\varepsilon$ are considered. Here, we choose $N = 160$ and $\rho = 100$ again, and we let $\frac{u}{10\varepsilon}$ denote the paint thickness. We only consider the paint distribution over half the boundary because the problem is symmetric, and the numerical results corresponding to $\varepsilon = 0.4$, $\varepsilon = 0.5$, $\varepsilon = 0.55$, and $\varepsilon = 0.7$ are presented in Figures 5-8, respectively. It can be observed that our results are in good agreement with the corresponding results in [2, 16, 20].

We also investigate the convergence behavior of our method for this example. Tables 2-5 display the number of iterations for the four cases with different $N$ and $\rho$. As can be seen from our tests, our method converges quickly when $\rho$ is sufficiently large and the number of iterations depends only weakly on $N$.

## 4.3 Signorini problem

Finally, we consider a Signorini boundary value problem in a domain defined by an ellipse [9],

$$ \Delta u = 0 \quad \text{in } \Omega = E(a, b), $$
<figure><img src="image_8.png" /><figcaption>Figure 6 Paint distributions for $\epsilon = 0.5$.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Figure 7 Paint distributions for $\epsilon = 0.55$.</figcaption></figure>

<figure><img src="image_10.png" /></figure>
<figure><img src="image_11.png" /><figcaption>Table 2 Number of iterations for $\epsilon = 0.4$ with different $N$ and $\rho$</figcaption></figure>

<table><caption>Table 3 Number of iterations for $\epsilon = 0.5$ with different $N$ and $\rho$</caption><thead><tr><th></th><th>N = 40</th><th>N = 80</th><th>N = 160</th><th>N = 320</th><th>N = 640</th></tr></thead><tbody><tr><td>ρ = 10<sup>0</sup></td><td>48</td><td>48</td><td>47</td><td>47</td><td>46</td></tr><tr><td>ρ = 10<sup>1</sup></td><td>15</td><td>15</td><td>16</td><td>16</td><td>16</td></tr><tr><td>ρ = 10<sup>2</sup></td><td>9</td><td>9</td><td>10</td><td>11</td><td>11</td></tr><tr><td>ρ = 10<sup>3</sup></td><td>8</td><td>8</td><td>9</td><td>9</td><td>9</td></tr><tr><td>ρ = 10<sup>4</sup></td><td>7</td><td>7</td><td>8</td><td>9</td><td>9</td></tr><tr><td>ρ = 10<sup>5</sup></td><td>6</td><td>6</td><td>7</td><td>8</td><td>8</td></tr></tbody></table>

<table><caption>Table 4 Number of iterations for $\epsilon = 0.55$ with different $N$ and $\rho$</caption><thead><tr><th></th><th>N = 40</th><th>N = 80</th><th>N = 160</th><th>N = 320</th><th>N = 640</th></tr></thead><tbody><tr><td>ρ = 10<sup>0</sup></td><td>22</td><td>25</td><td>24</td><td>25</td><td>24</td></tr><tr><td>ρ = 10<sup>1</sup></td><td>11</td><td>11</td><td>12</td><td>13</td><td>13</td></tr><tr><td>ρ = 10<sup>2</sup></td><td>8</td><td>9</td><td>9</td><td>11</td><td>11</td></tr><tr><td>ρ = 10<sup>3</sup></td><td>7</td><td>8</td><td>9</td><td>10</td><td>10</td></tr><tr><td>ρ = 10<sup>4</sup></td><td>6</td><td>7</td><td>8</td><td>10</td><td>10</td></tr><tr><td>ρ = 10<sup>5</sup></td><td>6</td><td>7</td><td>8</td><td>9</td><td>9</td></tr></tbody></table>

<table><caption>Table 5 Number of iterations for $\epsilon = 0.7$ with different $N$ and $\rho$</caption><thead><tr><th></th><th>N = 40</th><th>N = 80</th><th>N = 160</th><th>N = 320</th><th>N = 640</th></tr></thead><tbody><tr><td>ρ = 10<sup>0</sup></td><td>22</td><td>20</td><td>21</td><td>21</td><td>21</td></tr><tr><td>ρ = 10<sup>1</sup></td><td>10</td><td>11</td><td>11</td><td>12</td><td>12</td></tr><tr><td>ρ = 10<sup>2</sup></td><td>8</td><td>8</td><td>9</td><td>10</td><td>10</td></tr><tr><td>ρ = 10<sup>3</sup></td><td>7</td><td>7</td><td>8</td><td>9</td><td>10</td></tr><tr><td>ρ = 10<sup>4</sup></td><td>6</td><td>6</td><td>8</td><td>9</td><td>10</td></tr><tr><td>ρ = 10<sup>5</sup></td><td>6</td><td>6</td><td>8</td><td>9</td><td>9</td></tr></tbody></table>

<table><thead><tr><th></th><th>N = 40</th><th>N = 80</th><th>N = 160</th><th>N = 320</th><th>N = 640</th></tr></thead><tbody><tr><td>ρ = 10<sup>0</sup></td><td>17</td><td>17</td><td>16</td><td>16</td><td>16</td></tr><tr><td>ρ = 10<sup>1</sup></td><td>9</td><td>10</td><td>10</td><td>10</td><td>11</td></tr><tr><td>ρ = 10<sup>2</sup></td><td>7</td><td>8</td><td>8</td><td>9</td><td>10</td></tr><tr><td>ρ = 10<sup>3</sup></td><td>6</td><td>7</td><td>7</td><td>8</td><td>9</td></tr><tr><td>ρ = 10<sup>4</sup></td><td>5</td><td>7</td><td>7</td><td>8</td><td>9</td></tr><tr><td>ρ = 10<sup>5</sup></td><td>5</td><td>7</td><td>7</td><td>8</td><td>9</td></tr></tbody></table>
<figure><img src="image_12.png" /><figcaption>Figure 10 Approximate solutions for $\lambda$ on $\Gamma$.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Table 6 Number of iterations with different $N$ and $\rho$</figcaption></figure>

<table><thead><tr><th></th><th>N = 40</th><th>N = 80</th><th>N = 160</th><th>N = 320</th><th>N = 640</th></tr></thead><tbody><tr><td>ρ = 10<sup>0</sup></td><td>43</td><td>46</td><td>46</td><td>45</td><td>46</td></tr><tr><td>ρ = 10<sup>1</sup></td><td>15</td><td>16</td><td>17</td><td>17</td><td>18</td></tr><tr><td>ρ = 10<sup>2</sup></td><td>10</td><td>11</td><td>13</td><td>13</td><td>14</td></tr><tr><td>ρ = 10<sup>3</sup></td><td>8</td><td>10</td><td>12</td><td>12</td><td>13</td></tr><tr><td>ρ = 10<sup>4</sup></td><td>8</td><td>9</td><td>11</td><td>11</td><td>13</td></tr><tr><td>ρ = 10<sup>5</sup></td><td>7</td><td>9</td><td>11</td><td>11</td><td>13</td></tr></tbody></table>

with the Signorini boundary conditions

$$u \ge x + y^3, \quad \lambda \ge -x^2, \quad (u - x - y^3)(\lambda + x^2) = 0 \quad \text{on } \Gamma = \partial E(a, b),$$

where $E(a, b)$ denotes the ellipse $((x, y) : (x/a)^2 + (y/b)^2 < 1)$. The presented algorithm is applied with the parameterizations $t \to (a \cos 2\pi t, -b \sin 2\pi t)$ and $a = 0.25, b = 0.125$. For $N = 160$ and parameter $\rho = 100$, the numerical results are given in Figures 9 and 10. It can be seen that our results are again in excellent agreement with the corresponding results in [9]. Table 6 shows the numbers of iterations for different values of $\rho$ and various $N$. We observe
that the algorithm also converges quickly and the numbers of iterations decrease as $\rho$
increases, and a different $N$ has little effect on the numbers of iterations.

**5 Conclusion**

In this paper, we have proposed a new ALM for the solution of Signorini boundary value problems and proven its convergence. Using the BEM and the fixed point method, we can easily apply this algorithm to the Signorini boundary value problems defined in domains of arbitrary shape. Each iteration only needs to solve an elliptic variational problem and the semismooth Newton method is used to find the solution. The examples tested show the feasibility and effectiveness of the algorithm.

**Competing interests**

The authors declare that they have no competing interests.

**Authors' contributions**

All authors equally have made contributions. All authors read and approved the final manuscript.

**Acknowledgements**

This work was funded by China Scholarship Council, the Natural Science Foundation Project of CQ CSTC of China (Grant Nos. cstc2013jcyjA30001 and cstc2013jcyjA10049) and Fundamental Research Funds of Chongqing Normal University of China (Grant No. 13XLB001), the National Natural Science Foundation of China (Grant Nos. 11471063 and 11301575).

Received: 9 October 2015 Accepted: 5 March 2016 Published online: 10 March 2016

**References**

1. Glowinski, R: Numerical Methods for Nonlinear Variational Problems. Springer, Berlin (2008)

2. Aitchison, JM, Lacey, AA, Shillor, M: A model for an electropaint process. IMA J. Appl. Math. **33**, 269-287 (1984)

3. Chouly, F, Hild, P: A Nitsche-based method for unilateral contact problems: numerical analysis. SIAM J. Numer. Anal. **51**, 1295-1307 (2013)

4. Ryoo, CS: Numerical verification of solutions for Signorini problems using Newton-like method. Int. J. Numer. Methods Eng. **73**, 1181-1196 (2008)

5. Ryoo, CS: An approach to the numerical verification of solutions for variational inequalities using Schauder fixed point theory. Bound. Value Probl. **2014**, Article ID 235 (2014)

6. Hild, P, Renard, Y: An improved *a priori* error analysis for finite element approximations of Signorini's problem. SIAM J. Numer. Anal. **50**, 2400-2419 (2012)

7. Stadler, G: Path-following and augmented Lagrangian methods for contact problems in linear elasticity. J. Comput. Appl. Math. **203**, 533-547 (2007)

8. Han, HD: A direct boundary element method for Signorini problems. Math. Comput. **55**, 115-128 (1990)

9. Spann, P: On the boundary element method for the Signorini problem of the Laplacian. Numer. Math. **65**, 337-356 (1993)

10. Steinbach, O: Boundary element methods for variational inequalities. Numer. Math. **126**, 173-197 (2014)

11. Maischak, M, Stephan, EP: Adaptive hp-versions of BEM for Signorini problems. Appl. Numer. Math. **54**, 425-449 (2005)

12. Mel'nyk, TM: Homogenization of the Signorini boundary-value problem in a thick junction and boundary integral equations for the homogenized problem. Math. Methods Appl. Sci. **34**, 758-775 (2011)

13. Aitchison, JM, Poole, MW: A numerical algorithm for the solution of Signorini problems. J. Comput. Appl. Math. **94**, 55-67 (1998)

14. Ito, K, Kunisch, K: Semi-smooth Newton methods for the Signorini problem. Appl. Math. **53**, 455-468 (2008)

15. Bustinza, R, Sayas, FJ: Error estimates for an LDG method applied to Signorini type problems. J. Sci. Comput. **52**, 322-339 (2012)

16. Karageorghis, A, Lesnic, D, Marin, L: The method of fundamental solutions for solving direct and inverse Signorini problems. Comput. Struct. **151**, 11-19 (2015)

17. Wang, F, Han, WM, Cheng, XL: Discontinuous Galerkin methods for solving the Signorini problem. IMA J. Numer. Anal. **31**, 1754-1772 (2011)

18. Li, F, Li, XL: The interpolating boundary element-free method for unilateral problems arising in variational inequalities. Math. Probl. Eng. **2014**, Article ID 518727 (2014)

19. Zhang, SG, Zhu, JL: A projection iterative algorithm boundary element method for the Signorini problem. Eng. Anal. Bound. Elem. **37**, 176-181 (2013)

20. Zhang, SG: A projection iterative algorithm for the Signorini problem using the boundary element method. Eng. Anal. Bound. Elem. **50**, 313-319 (2015)

21. Cottle, RW, Pang, JS, Stone, RE: The Linear Complementarity Problem. Academic Press, Boston (1992)

22. He, BS, Liao, LZ: Improvements of some projection methods for monotone nonlinear variational inequalities. J. Optim. Theory Appl. **11**, 111-128 (2002)

23. Noor, MA: Some developments in general variational inequalities. Appl. Math. Comput. **152**, 199-277 (2004)

24. Hsiao, GC: Boundary Integral Equations. Springer, Berlin (2008)

25. Steinbach, O: Numerical Approximation Methods for Elliptic Boundary Value Problems. Springer, New York (2008)

26. Steinbach, O: Stability Estimates for Hybrid Coupled Domain Decomposition Methods. Springer, Berlin (2003)