**2. Elliptic variational inequalities**

The theory of elliptic variational inequalities has become a rich source of inspiration in both mathematical and engineering sciences. Elliptic variational inequalities are an effective tool for studying the existence of solutions of constrained problems arising in mechanics, optimization and control, operation research, engineering science, etc. [47–52]. It is the aim of this chapter to introduce a numerical technique to verify the solutions for elliptic variational inequalities. The basic approach of this technique consists of the fixed point formulation of elliptic variational inequalities and construction of the function set, on computer, satisfying the validation condition of a certain infinite dimensional fixed point theorem. For fixed point formulation, we consider a candidate set which possibly contains a solution. In order to get such a candidate set, we divide the verification procedure into two phases: one is the computation of a projection into a closed convex subset of some finite dimensional subspace (rounding); the other is the estimation of the error for the projection (rounding error). Combining these methods with some iterative technique, the exact solution can be enclosed by sum of rounding parts, which is a subset of finite dimensional space, and the rounding error, which is indicated by a nonnegative real number. These two procedures enable us to treat infinite dimensional problems as finite procedures, thta is, by computer.

### **Notations**


A bilinear form *a*ð Þ �, � is said to be *V*-*elliptic* if there exists a positive constant *α* such that *a v*ð Þ , *<sup>v</sup>* <sup>≥</sup>*α*∥*v*∥<sup>2</sup> , ∀*v*∈*V:*.

In general we do not assume *a*ð Þ �, � to be symmetric, since in some applications nonsymmetric bilinear forms may occur naturally.


*Numerical Verification Method of Solutions for Elliptic Variational Inequalities DOI: http://dx.doi.org/10.5772/intechopen.101357*

• *j*ð Þ� : *V* ! **R**∪f g ∞ is a convex lower semicontinuous (l.s.c) and proper functional (*j*ð Þ� is proper if *j v*ð Þ> � ∞, ∀*v*∈*V* and *j* 6¼ þ∞Þ.

#### **The two types of elliptic variational inequalities.**

We consider two classes of elliptic variational inequalities.

• Elliptic variational inequalities of the first kind: *Find u*∈*V such that u is a solution of the problem*

$$a(u, v - u) \ge L(v - u), \forall v \in K, u \in K.$$

• Elliptic variational inequalities of the second kind: *Find u* ∈*Vsuch that u is a solution of the problem*

$$a(u, -u) + j(v) - j(u) \ge L(v - u), \forall v \in V, u \in V.$$

#### **2.1 The problem and the fixed point formulation**

Let us first set a few notations [1, 47, 49, 50, 53–61]. In what follows we shall make use of the Sobolev spaces *<sup>W</sup><sup>k</sup>*,*<sup>p</sup>*ð Þ <sup>Ω</sup> of functions which possess generalized derivatives integrable with the *p*th power up to and including the *k*th order. For *<sup>p</sup>* <sup>¼</sup> 2, we shall write *<sup>W</sup><sup>k</sup>*,*<sup>p</sup>*ð Þ¼ <sup>Ω</sup> *<sup>H</sup><sup>k</sup>* ð Þ <sup>Ω</sup> , *<sup>H</sup>*<sup>0</sup>ð Þ¼ <sup>Ω</sup> *<sup>L</sup>*<sup>2</sup> ð Þ Ω *:* Further, we introduce the scalar product in *L*<sup>2</sup> ð Þ Ω by

$$(f, \mathbf{g}) = \int\_{\tilde{\Omega}} f(\mathbf{x}) \mathbf{g}(\mathbf{x}) d\mathbf{x}.$$

The norm in *H<sup>k</sup>* ð Þ <sup>Ω</sup> will be denoted by <sup>∥</sup> � <sup>∥</sup>*Hk*ð Þ <sup>Ω</sup> *:* The symbol j j� *<sup>H</sup>k*ð Þ <sup>Ω</sup> will stand for the seminorm,

$$|u|\_{H^k(\mathfrak{Q})} = \left(\sum\_{|\alpha|=k} \|D^\alpha u\|\_{L^2(\mathfrak{Q})}^2\right)^{\frac{1}{2}}, \quad \|u\|\_{H^k(\mathfrak{Q})} = \left(\sum\_{j=0}^k |u|\_{H^j(\mathfrak{Q})}^2\right)^{\frac{1}{2}}.$$

Let *V* be a real Hilbert space with a scalar product ð Þ �, � *<sup>V</sup>* and an associated norm <sup>∥</sup> � <sup>∥</sup>*V*, *<sup>V</sup>* <sup>∗</sup> its dual space. *<sup>K</sup>* denotes a nonempty closed convex subset of *V*, *a*ð Þ �, � : *V* � *V* ! **R** is a bilinear, symmetric, continuous and elliptic form of *V*, *a*ð Þ �, � : *V* � *V* ! **R** is a bilinear, symmetric, continuous and elliptic form of *V* � *V*; that is, there exist constants *α* >0, and *β* >0 such that*a u*ð Þ , *v* ≤*α*∥*u*∥*V*∥*v*∥*V*, ∀*u*, *v*∈*V* and *a v*ð Þ , *<sup>v</sup>* <sup>≥</sup>*β*∥*v*∥<sup>2</sup> *<sup>V</sup>*, <sup>∀</sup>*v*∈*V*. The pairing between *<sup>V</sup>* and *<sup>V</sup>* <sup>∗</sup> is denoted by <sup>&</sup>lt; � , � <sup>&</sup>gt;. Let <sup>Λ</sup> be a canonical isomorphism from *<sup>V</sup>* <sup>∗</sup> onto *<sup>V</sup>* defined, for *<sup>g</sup>* <sup>∈</sup>*<sup>V</sup>* <sup>∗</sup> , by <sup>&</sup>lt;*g*, *<sup>v</sup>*<sup>&</sup>gt; <sup>¼</sup> ð Þ <sup>Λ</sup>*g*, *<sup>v</sup> <sup>V</sup>*, <sup>∀</sup>*v*∈*V:* We can easily see that <sup>∥</sup>Λ∥*<sup>V</sup>* <sup>∗</sup> <sup>¼</sup> <sup>∥</sup>Λ�<sup>1</sup> ∥*<sup>V</sup>* ¼ 1. Now, let us consider the following variational inequality:

$$\text{Find } \ u \in K \text{ such that } a(u, v - u) \ge < f(u), v - u > , \forall v \in K,\tag{1}$$

where *<sup>f</sup>* is a nonlinear operator such that *f u*ð Þ∈*<sup>V</sup>* <sup>∗</sup> *:*.

In order to obtain a fixed point formulation of variational inequality (1) we need the following standard result.

**Lemma 1.** *Let K be a closed convex subset of V. Then u* ¼ *PKω, the projection of ω on K, if and only if*

$$u \in K : (u - a, v - u)\_V \ge 0, \forall v \in K. \tag{2}$$

For some constant *ρ* >0, let us define a mapping *G* : *V* ! *V* by

$$G(u) = P\_K \Lambda \Phi(u),\tag{3}$$

where *<sup>u</sup>* <sup>∈</sup>*V*, <sup>Φ</sup>ð Þ *<sup>u</sup>* <sup>∈</sup>*<sup>V</sup>* <sup>∗</sup> is defined by

$$<\Phi(u), v> = (u, v)\_V - \rho a(u, v) + \rho < f(u), v>, \forall v \in V. \tag{4}$$

For some constant *ρ* >0, problem (1) can be written as

$$\left\{ (u, v - u)\_V - \left\{ (u, v - u)\_V - \rho a(u, v - u) + \rho \lhd f(u), v - u > \right\} \ge 0, \forall v \in K. \right\}$$

Using (4) in the above inequality, problem (1) is equivalent to that of finding *u*∈ *K* such that

$$(u - \Lambda \Phi(u), v - u)\_V \ge 0, \forall v \in K. \tag{5}$$

By (2) and (5), we now have the following fixed point problem for the operator *G*:

$$
\mu = P\_K \Lambda \Phi(\mu) = \mathcal{G}(\mu). \tag{6}
$$

Under appropriate conditions on the space *V* and the operator *G* : *V* ! *V*(e.g., continuity, compactness), which usually have to be verified by theoretical means, fixed point theorem yields the existence of a solution *u* of the problem (1) in some suitable set *U* ⊂*V*, provided that

$$\mathcal{G}(U) \subset U. \tag{7}$$

In order to compute an explicit inclusion, we must therefore construct *U* explicitly. For the numerical verification of condition (7), we have to use interval analysis on many levels between basic interval arithmetic and functional analysis. For the appropriate and suitable choice of the operator *f*, the form *a*ð Þ �, � , and the convex set *K*; one encounters problems governed by the elliptic variational inequality as special cases from the problem (1) [48–52]. Inbrief, it is clear that the problem (1) is the most common. Up to now, devising a verification technique for the problem (1) is still an open problem. It is an important and interesting area of future research to find the numerical inclusion methods for the problem (1) by using (6). In this paper, we suppose that *V* ⊂*L*<sup>2</sup> ð Þ <sup>Ω</sup> and the nonlinear map *<sup>f</sup>*ð Þ� : *<sup>V</sup>* ! *<sup>L</sup>*<sup>2</sup> ð Þ Ω satisfies the following assumptions.

**A1.** *f* is a continuous map from *V* to *L*<sup>2</sup> ð Þ Ω *:*

**A2.** For each bounded subset *<sup>W</sup>* <sup>⊂</sup>*V*, *f W*ð Þ is also bounded in *<sup>L</sup>*<sup>2</sup> ð Þ Ω *:*

If we restrict the nonlinear map *f* as above, then it can be shown that the problem (1) can be characterized by a class of variational inequality of the type,

$$\text{find } \ u \in K \text{ such that } \ a(u, v - u) \ge (\ f(u), v - u), \forall v \in K. \tag{8}$$

The problem (8) has the restricted condition; even so (8) is an important and very useful class of nonlinear problems arising in mathematical physics, mechanics, engineering sciences, etc. In Section 3, we briefly consider a particular example of interest in applications. Another example is given in [13, 24–46]. In the special case in which *K* � *V*, (8) yields the variational theory of the boundary value problems

*Numerical Verification Method of Solutions for Elliptic Variational Inequalities DOI: http://dx.doi.org/10.5772/intechopen.101357*

for partial differential equations. We will discuss existence and inclusion methods for problem (8). These are methods providing the existence of a solution of the problem (8) within explicitly computable bounds. As we have seen before, the transformation of problem (8) into some fixed point formulation (6) can be carried out in the same way. In a conclusion problem (8) is equivalent to the fixed point problem of finding *u* ∈*K* such that

$$
u = \mathbb{S}(u),\tag{9}$$

where *S* denotes a specific operator, not necessarily the same as in (6). In particular for a given problem, we reduced the problem (8) to the fixed point formulation (9) and the continuity and compactness of *S* is discussed. For this reason, we shall say nothing about this problem for which we refer to [13, 24–46]. In order to simplify argument we assume that *S* is a continuous and compact operator. Since *S* is continuous and compact, as a result of Schauder's fixed point theorem, if there exists a nonempty, bounded, convex, and closed subset *U* such that *S U*ð Þ⊂ *U*, then there exists a solution of *u* ¼ *S u*ð Þ in *U*. In Sections 2.2 and 2.3, we describe how to construct *U* explicitly.

#### **2.2 Verification by a simple iteration method**

In this subsection, we describe a simple iteration method for numerical verification of solutions for elliptic variational inequalities. In order to treat functions and variational inequalities in the infinite dimensional space *V* by computer, we introduce two concepts, rounding and rounding error. Now, let *Vh* be a finite dimensional subspace of *V* dependent on *h*ð Þ 0<*h*<1 and let *Kh* be a nonempty closed convex subset of *Vh*. Usually, *Vh* is taken to be a finite element subspace with mesh size *h*. For the sake of simplicity, we shall define *Kh*, an approximate subset of *K*, by *Kh* ¼ *Vh* ∩*K: Kh* is a closed convex subset of *Vh*. In practical applications, the construction of *Kh* is one of the difficulties presented by variational inequalities. For a given problem, several approximations are available. For a general study of the approximation of convex sets, we refer the reader to the work of Mosco [51]. We define the projection *PKh* from *V* into *Kh* [49, 50]. That is, *vh* ¼ *PKh* ð Þ *u* , the projection of *u* into *Kh*, is defined as follows:

$$u = \mathcal{S}(u), v\_h \in \mathcal{K}\_h : \quad (v\_h, \zeta - v\_h)\_V \ge (u, \zeta - v\_h)\_V, \quad \forall \zeta \in \mathcal{K}\_h. \tag{10}$$

To verify the existence of a solution of (9), we determine a set *W* for a bounded, convex, and closed subset *U* ⊂*V* as

$$W = \{ v \in V \, : \, v = S(u), u \in U \}.$$

From Schauder's fixed point theorem, if *W* ⊂ *U* holds, then there exists a solution of (8) in the set *U*. Our goal is to find a set *U* which includes *W:* For any subset *W* ⊂*V*, we define *R W*ð Þ⊂*Kh* by the projection of *V* to *Kh*, which is called the rounding of *W*. Additionally, we define RE(W), the rounding error of *W*, as a subset of *V* so that *W* ⊂ *R W*ð Þþ *RE W*ð Þ holds. Using *R W*ð Þþ *RE W*ð Þ instead of *W*, the verification condition becomes

$$R(W) + RE(W) \subset U. \tag{11}$$

Let us describe the procedure more concretely. First, we consider the auxiliary problem: given *g* ∈*L*<sup>2</sup> ð Þ Ω ,

$$\text{find } u \in K \text{ such that } a(u, v - u) \ge (\mathfrak{g}, v - u), \forall v \in K. \tag{12}$$

We note that, by well known result [49], there is a unique element *u* which satisfies (12).

Secondly, we define the approximate problem corresponding to (12) as

$$a(u\_h, v\_h - u\_h) \ge (\mathbf{g}, v\_h - u\_h), \forall v\_h \in K\_h, u\_h \in K\_h \tag{13}$$

and (13) admit one and only one solution [49]. Error estimates for the variational inequalities can be found in [48, 49, 52], etc. Now, using (10), (12), (13) and error estimates, we make the following assumption.

**A3.** For each *u* ∈*V*, there exists a positive constant *C*, independent of *u* and *h*, such that

$$\|\|u - P\_{K\_h}u\|\|\_{V} \le \mathcal{C}h \|\|\mathbf{g}\|\|\_{L^2(\Omega)}.\tag{14}$$

In order to verify the solutions numerically, it is necessary to determine the constant *C* that appears in a priori error estimations; this constant will be discussed later.

In order to construct the set *U* satisfying the verification condition (11) in a computer, we use an iterative procedure, that is, the sequential iteration. We propose a computer algorithm to obtain the set *U* which satisfies the condition (11).

(1) First, we obtain an approximate solution *v* ð Þ 0 *<sup>h</sup>* ∈*Kh* to (8) by an appropriate method. Set *U*ð Þ <sup>0</sup> *<sup>h</sup>* ¼ *v* ð Þ 0 *h* n o and *<sup>α</sup>*<sup>0</sup> <sup>¼</sup> <sup>0</sup>*:*.

(2) Next we will define *R W*ð Þ*<sup>i</sup>* � � and *RE W*ð Þ*<sup>i</sup>* � � for *i* ≥0, where *W*ð Þ*<sup>i</sup>* is the set defined as follows:

$$\mathcal{W}^{(i)} = \left\{ v^{(i)} \in V : v^{(i)} = \mathcal{S} \left( u^{(i)} \right), \quad u^{(i)} \in U^{(i)} \right\}.$$

*R W*ð Þ*<sup>i</sup>* � � is defined by the subset of *Kh* which consists of all the elements *v* ð Þ*i <sup>h</sup>* ∈*Kh* such that

$$a\left(v\_h^{(i)}, \boldsymbol{\psi} - v\_h^{(i)}\right) \ge \left(f\left(u^{(i)}\right), \boldsymbol{\psi} - v\_h^{(i)}\right), \quad \forall \boldsymbol{\psi} \in K\_h,\tag{15}$$

holds for some *u*ð Þ*<sup>i</sup>* ∈ *U*ð Þ*<sup>i</sup> :* Note that *R W*ð Þ*<sup>i</sup>* � � can be enclosed by *R W*ð Þ*<sup>i</sup>* � �⊂P*<sup>M</sup> <sup>j</sup>*¼<sup>1</sup>*<sup>A</sup> <sup>j</sup>ϕj*, where *<sup>A</sup> <sup>j</sup>* <sup>¼</sup> *<sup>A</sup> <sup>j</sup>*, *<sup>A</sup> <sup>j</sup>* h i are intervals, *<sup>ϕ</sup><sup>j</sup>* n o*<sup>M</sup> j*¼1 is a basis of *Vh*, and *M* ¼ dim*Vh*. For details of the interval calculation, we refer the reader to Nakao [6, 7, 12]. Next *RE W*ð Þ*<sup>i</sup>* � � is defined as

$$RE\left(\mathcal{W}^{(i)}\right) = \left\{ v \in V : \quad \|v\|\_{V} \le Ch \sup\_{u^{(i)} \in U^{(i)}} \|f\left(u^{(i)}\right)\|\_{L^2(\Omega)} \right\}.\tag{16}$$

Here, *<sup>C</sup>* is the same constant as in (14). Hence, *<sup>W</sup>*ð Þ*<sup>i</sup>* <sup>⊂</sup>*R W*ð Þ*<sup>i</sup>* � � <sup>þ</sup> *RE W*ð Þ*<sup>i</sup>* � � holds. (3) Check the verification condition:

$$R\left(\boldsymbol{W}^{(i)}\right) + RE\left(\boldsymbol{W}^{(i)}\right) \subset U^{(i)}.\tag{17}$$

If the condition is satisfied, then *U*ð Þ*<sup>i</sup>* is the desired set, and a solution to (8) exists in *W*ð Þ*<sup>i</sup>* , and hence in *U*ð Þ*<sup>i</sup> :*

*Numerical Verification Method of Solutions for Elliptic Variational Inequalities DOI: http://dx.doi.org/10.5772/intechopen.101357*

(4) If the condition is not satisfied, we continue the simple iteration by using *δ* � *inflation*; that is, let *δ* be a certain positive constant given beforehand, and take

$$\begin{aligned} a\_{i+1} &= \operatorname{Chsup}\_{\mathfrak{u}^{(i)} \in U^{(i)}} \| f \left( \mathfrak{u}^{(i)} \right) \|\_{L^2(\Omega)} + \delta, \\ [a\_{i+1}] &= \{ v \in V : \| \boldsymbol{v} \|\|\_{V} \le a\_{i+1} \}, \\ U\_h^{(i+1)} &= \sum\_{j=1}^{M} \left[ \underline{A}\_j - \delta, \ \overline{A\_j} + \delta \right] \phi\_j, \\ U^{(i+1)} &= U\_h^{(i+1)} + [a\_{i+1}], \end{aligned}$$

and then go back to the second step. The reader may refer to [26–46] for the details. If the condition (17) is satisfied, in our inclusion method of solutions for (9), the solution *u* is enclosed in the set *U*ð Þ*<sup>i</sup>* , which we call 'a candidate set' of the form *<sup>U</sup>*ð Þ*<sup>i</sup>* <sup>¼</sup> *<sup>U</sup>*ð Þ*<sup>i</sup> <sup>h</sup>* þ *α<sup>i</sup>* ½ �.

#### **2.3 Verification by a Newton-like method**

The significance of a Newton-like operator was already pointed out in [29, 43]. Hence we will not discuss it in detail here. In Subsection 2.1, numerical verification of solutions for elliptic variational inequalities using a finite element method have been discussed only for simple iteration method. The method proposed in Subsection 2.2 is such that *<sup>U</sup>*ð Þ*<sup>i</sup> <sup>h</sup>* , *<sup>α</sup><sup>i</sup>* n o � � always converges to the limit value f g ð Þ *Uh*, *<sup>α</sup>* from an arbitrary initial value *<sup>U</sup>*ð Þ <sup>0</sup> *<sup>h</sup>* , *<sup>α</sup>*<sup>0</sup> n o � � if *S* in (9) is retractive operator (we refer to Zeidler [59–61] for the definition of retraction), while no convergence can generally be expected if *S* is not retractive operator. Briefly, for not retractive operator in the neighborhood of the solution, it is difficult to use the previous scheme proposed in Subsection 2.2. To overcome such a difficulty, in this section, we newly formulate a verification method using the Newton-like method. This approach enables us to remove the restriction in Subsection 2.2 to the retraction property of the operator in the neighborhood of the solution. Namely, this technique can be applied to general variational inequalities without any retraction property of the associated operator *S*. We refer to [29, 43] for a detailed study of the properties of the Newton-like Method.

In this subsection, we use the notation of Section 2.2. We assume that *Kh* ¼ *Vh* ∩*K* is a closed convex cone with vertex at 0 and *K*<sup>∗</sup> *<sup>h</sup>* its dual. We note that *K* <sup>∗</sup> *<sup>h</sup>* is also a closed convex cone with vertex at 0, which is the only point common to *Kh* and *K*<sup>∗</sup> *<sup>h</sup>* . From (10) it follows that *K* <sup>∗</sup> *<sup>h</sup>* is the set of points whose projections into *Kh* is 0. We need some additional lemma.

**Lemma 2.** *Any u*∈*V can be uniquely decomposed into the sum of two orthogonal elements. That is,*

$$
\mu = P\_{K\_h} \mathfrak{u} \oplus \left( I - P\_{K\_h} \right) \mathfrak{u} = P\_{K\_h} \mathfrak{u} \oplus P\_{K\_h^\*} \mathfrak{u} \dots
$$

*Here,* ⊕ *denotes the sum of two orthogonal elements in the sense of V:.* Note that (9) can be rewritten as the following decomposed form in *Kh* and *K* <sup>∗</sup> *h* :

$$\begin{cases} P\_{K\_k} \mathfrak{u} = P\_{K\_k} \mathcal{S}(\mathfrak{u}), \\ \left( I - P\_{K\_k} \right) \mathfrak{u} = \left( I - P\_{K\_k} \right) \mathcal{S}(\mathfrak{u}). \end{cases} \tag{18}$$

In order to formulate a Newton-like verification condition for (18), we need a Fréchet derivative of the operator *S*. For most of the variational inequalities, the *S* in (9) is not Fréchet differentiable at all. Therefore, in order to use a Newton-like type method, a major difficulty in numerically solving the fixed point formulation *u* ¼ *S u*ð Þ is the treatment of the non-differentiable operator *S*. We need a suitable modification of the Fréchet derivative of *S*. Using some techniques, we can devise the approximate Fréchet derivative of *S*. Hence we shall assume that *DS u* e ð Þ is the approximate Fréchet derivative of the *S u*ð Þ at *u* as the linear operator. Let *DS u* e ð Þ be designated as the Fréchet-like derivative of *S* at *u*.

To consider the Newton-like operator for (18), we define the nonlinear operator *Nh* : *V* ! *Vh* as

$$N\_h(\boldsymbol{\mu}) \equiv P\_{K\_h} \boldsymbol{\mu} - \left[\boldsymbol{I} - \widetilde{\boldsymbol{D}} \mathbf{S}(\boldsymbol{\mu}\_h)\right]\_h^{-1} \left(P\_{K\_h} - P\_{K\_h} \mathbf{S}\right)(\boldsymbol{\mu}).$$

Here *I* is the identity operator and *I* � *DS u* e ð Þ*<sup>h</sup>* h i�<sup>1</sup> *h* denotes the inverse on *Vh* of the restriction operator *I* � *DS u* e ð Þ*<sup>h</sup>* h i� � � *Vh :* Note that we will verify the existence of the inverse operator *I* � *DS u* e ð Þ*<sup>h</sup>* h i�<sup>1</sup> *<sup>h</sup>* from the nonsingularity of the matrix corresponding to *I* � *DS u* e ð Þ*<sup>h</sup>* h i� � � *Vh* in actual calculations.

Next we define the operator *T* : *V* ! *V* as follows:

$$T(u) \equiv N\_h(u) + \left(I - P\_{K\_h}\right) \mathcal{S}(u). \tag{19}$$

Then *T* is considered as the Newton-like operator for the former part of (18), but as the simple iterative operator for the latter part. *T* becomes a compact and continuous map on *V* by properties of *S*. Using some techniques, for a given problem we can not only define the Newton-like operator, but also devise a Newton-like Method. Furthermore, we obtain the following proposition and theorem.

**Proposition 3.** *Given the assumption that Nh*ð Þ *u* ∈*Kh*,

$$
\mathfrak{u} = \mathfrak{S}(\mathfrak{u}) \Leftrightarrow \mathfrak{u} = T(\mathfrak{u}).\tag{20}
$$

**Theorem 4.** *If there exists a nonempty, bounded, convex, and closed subset U* ⊂*K such that T U*ð Þ¼ f g *T u*ð Þj*u*∈ *U* ⊂ *U, then by the Schauder fixed point theorem, there exists a solution u*∈ *U of u* ¼ *S u*ð Þ*.*

When we decompose the set *U* as *U* ¼ *Uh* ⊕ *U*<sup>⊥</sup> in Theorem 8.1, where *Uh* ⊂ *Kh* and *U*<sup>⊥</sup> ⊂*K* <sup>∗</sup> *<sup>h</sup>* , the verification condition can be written by

$$\begin{cases} N\_h(U) \subset U\_h, \\ (I - P\_{K\_\mathbb{L}}) S(U) \subset U\_\perp. \end{cases} \tag{21}$$

Here, *Uh* is represented as the linear combination of the base functions of *Vh* with interval coefficients, whereas *U*<sup>⊥</sup> is the intersection of *K*<sup>∗</sup> *<sup>h</sup>* with a ball in *V*. That is,

$$U\_h = \left\{ \boldsymbol{\varrho}\_h \in K\_h : \boldsymbol{\varrho}\_h = \sum\_{j=1}^M A\_j \boldsymbol{\phi}\_j \text{ with } \boldsymbol{a}\_j \in \left[ \underline{A\_j}, \overline{A\_j} \right] \right\},$$

$$U\_\perp = \left\{ \boldsymbol{\varrho} \in K\_h^\* \, : \quad \|\boldsymbol{\varrho}\|\, \boldsymbol{\varrho} \le \boldsymbol{a} \right\},$$

respectively.

Note that *Nh*ð Þ *U* can be directly computed from *Uh* and *U*<sup>⊥</sup> with additional information on the a priori error estimates. On the other hand, *I* � *PKh* � �*S U*ð Þ is evaluated using (14), by the following constructive error estimates for the finite approximate solution of variational inequality (8):

$$\|\|\left(I - P\_{K\_h}\right)S(U)\|\|\_{V} \le Ch \sup\_{\mathfrak{u} \in U} \|\|f(\mathfrak{u})\|\|\_{L^2(\Omega)}.$$

Therefore, the former condition in (21) is validated as the inclusion relations of corresponding coefficient intervals; the latter part can be checked by comparing two nonnegative real numbers.

Next we show a computer algorithm to construct the set *U* which satisfies the verification condition (21). In order to realize it, we use the iteration method described in Subsection 2.2. Similarly to that in Subsection 2.2, we now generate the following iteration sequence *U*ð Þ *<sup>n</sup> <sup>h</sup>* , *α<sup>n</sup>* n o � � for *<sup>n</sup>* <sup>¼</sup> 0, 1, 2, <sup>⋯</sup>*:* For *<sup>n</sup>*≥1, the *δ*-inflation of *U*ð Þ *<sup>n</sup>*�<sup>1</sup> *<sup>h</sup>* , *α<sup>n</sup>*�<sup>1</sup> � � is denoted by *<sup>U</sup>*eð Þ *<sup>n</sup>*�<sup>1</sup> *<sup>h</sup>* , <sup>e</sup>*α<sup>n</sup>*�<sup>1</sup> � �. Next, for the set *<sup>U</sup>*eð Þ *<sup>n</sup>*�<sup>1</sup> <sup>¼</sup> *<sup>U</sup>*eð Þ *<sup>n</sup>*�<sup>1</sup> *<sup>h</sup>* <sup>⊕</sup> ½ � <sup>e</sup>*α<sup>n</sup>*�<sup>1</sup> , define *<sup>U</sup>*ð Þ *<sup>n</sup> <sup>h</sup>* , *α<sup>n</sup>* � � by

$$\begin{cases} \boldsymbol{U}\_h^{(n)} \supset \boldsymbol{N}\_h \left( \boldsymbol{\widetilde{U}}^{(n-1)} \right), \\ \boldsymbol{a}\_n = \boldsymbol{C}h \sup\_{\boldsymbol{u} \in \boldsymbol{(n-1)}} \|\boldsymbol{f}(\boldsymbol{u})\|\_{L^2(\Omega)}. \end{cases} \tag{22}$$

Finally, the verification condition in a computer is given by the following theorem. The proof of Theorem 4 will be given here for the sake of completeness; it is based on Proposition 3 and Schauder's fixed point theorem.

**Theorem 5.** *For an integer N, if two relationships*

$$U\_h^{(N)} \subset \tilde{U}\_h^{(N-1)} \quad \text{and} \quad a\_N < \tilde{a}\_{N-1} \tag{23}$$

hold, then there exists a solution *u* of (8) in *U*ð Þ *<sup>N</sup> <sup>h</sup>* ⊕ ½ � *α<sup>N</sup>* . Here, the first term of (21) means the strict inclusion in the sense of each coefficient interval of *U*ð Þ *<sup>N</sup> <sup>h</sup>* and *<sup>U</sup>*eð Þ *<sup>N</sup>*�<sup>1</sup> *<sup>h</sup>* .

### **3. Applications**

The study for the numerical verification method for elliptic variational inequalities has been still made less progress than for the differential equation case. The author's method in the present chapter can be also applied, in principal, to the verification of solutions of the practical problems. Namely, in Section 3.1, we first give, a slightly detailed descroption of the basic principle and formulation of our numerical verification method for the solution of obstacle problems with a homogeneous condition. This should be an appropriate introduction to another applications of our idea. The basic approach of the method consisits of the fixed point formulation of the problems and construction of the function set, in a computer, satisying the validation condition of a certain infinite dimensional fixed point theorem. We also mention that it is possible to extend the method to more general problems with non-homogenerous obstacles. Moreover, in order to apply our method to the problem whose associated operator is not retractive in a neighborhood of the solution, a Newton-like method is introduced. Next, in Section 3.2, we apply our method to another type of free boundary problem with appears in the

elasto-plastic deformation theory. This problem causes some properties of nonsmoothness in tha associated finite dimensional equations. But, we can also overcome such a difficulty by appling the solution method for non-smooth problems developed by [29, 32, 33]. In the Section 3.3, we briefly remark that our enclosure method can also be applied to the so-called simplified Signorini problem which is a simplified version of a problem accurring in the elasticity theory [43]. Finally, in Section 3.4, we show the way to apply our approach to elliptic variational inequalities of the second kind appearing in the flow problems of a viscos-plastic fulied in a pipe.

### **3.1 Obstacle problems**

We introduce the verification method for solutions of the obstacle problem which is known as a free boundary problem to cahracterize the contacted zone by an obstacle *ψ* in an elastic membrane region.

#### *3.1.1 Homogeneous case*

Here, 'homogeneous'stands for the case that obstacle *ψ* � 0 in the whole domain.

### *3.1.1.1 Basic formulation of verification*

Though the basic idea of verification is given in other places [26–28], in order to keep the paper as self-contained as possible, we describe rather detailed formulation and verification procedure for the present case.

Let Ω be a bounded convex domain in *<sup>n</sup>*, 1≤*n*≤ 2, with piecewise smooth boundary <sup>∂</sup>Ω. We set *<sup>V</sup>* � *<sup>H</sup>*<sup>1</sup> <sup>0</sup>ð Þ¼f <sup>Ω</sup> *<sup>v</sup>*<sup>∈</sup> *<sup>H</sup>*<sup>1</sup> ð Þ Ω : *v* � � <sup>∂</sup><sup>Ω</sup> ¼ 0g and

$$a(u, v) = (\nabla u, \nabla v)$$

which is adopted as the inner product on *V*, where ð Þ �, � stands for the inner product on *L*<sup>2</sup> ð Þ Ω . We define *K* ≔ f g *v*∈*V* : *v*≥0 *a:e: on* Ω *:*.

First, we note that, by well-known result [49], for any *g* ∈*L*<sup>2</sup> ð Þ Ω , the problem:

$$a(u, v - u) \ge (\mathbf{g}, v - u), \quad \forall v \in K, \quad u \in K,\tag{24}$$

has a unique solution *u*∈*V* ∩ *H*<sup>2</sup> ð Þ Ω , and the estimate

$$\|\mathfrak{u}\|\_{H^{2}(\mathfrak{Q})} \leq \|\mathfrak{g}\|\_{L^{2}(\mathfrak{Q})} \tag{25}$$

holds [49], where j j *<sup>w</sup> <sup>H</sup>*<sup>2</sup> implies the semi-norm of *<sup>w</sup>* in *<sup>H</sup>*<sup>2</sup> ð Þ Ω defined by

$$\|\boldsymbol{w}\|\_{H^{2}(\Omega)}^{2} \equiv \sum\_{i,j=1}^{n} \|\frac{\partial^{2}\boldsymbol{w}}{\partial\boldsymbol{\kappa}\_{i}\partial\boldsymbol{\kappa}\_{j}}\|\_{L^{2}(\Omega)}^{2}.$$

Now consider the following elliptic variational inequalities with nonlinear righthand side;

$$\begin{cases} \text{Find } \ w \in K \quad \text{such that} \\ a(w, v - w) \ge (f(w), v - w), \quad \forall v \in K. \end{cases} \tag{26}$$

*Numerical Verification Method of Solutions for Elliptic Variational Inequalities DOI: http://dx.doi.org/10.5772/intechopen.101357*

We take an appropriate finite dimensional subspace *Vh* of *V* for 0 <*h*<1. Usually, *Vh* is taken to be a finite element subspace with mesh size *h*. We then define *Kh*, an approximation of *K*, by

$$K\_h = V\_h \cap K = \left\{ v\_h \, | \, v\_h \in V\_h, \quad v\_h \ge 0 \quad \text{on} \quad \overline{\Omega} \right\}.$$

We also define the projection *PK* from *V* onto *K*. That is, *v* ¼ *PK*ð Þ *w* , the projection of *w* ∈*V* into *K*, is defined as the unique solution of the following problem:

$$v \in K: \quad a(v, \zeta - v) \ge a(w, \zeta - v), \quad \forall \zeta \in K. \tag{27}$$

And define the projection *PKh* from *V* onto *Kh*. That is, *vh* ¼ *PKh* ð Þ *w* , the projection of *w* into *Kh*, is defined as follows:

$$w\_h \in K\_h : \quad a(v\_h, \zeta - v\_h) \ge a(w, \zeta - v\_h), \quad \forall \zeta \in K\_h. \tag{28}$$

Now, as one of the approximation properties of *Kh*, assume that.

For each *w* ∈*K* ∩ *H*<sup>2</sup> ð Þ Ω , there exists a positive constant *C*1, independent of *h*, such that

$$\|\boldsymbol{w} - P\_{K\_h}\boldsymbol{w}\|\_{V} \leq C\_1 h |\boldsymbol{w}|\_{H^2(\Omega)}.\tag{29}$$

Here, *C*<sup>1</sup> has to be numerically determined. For example, it is known that we may take *<sup>C</sup>*<sup>1</sup> <sup>¼</sup> ffiffi 5 p *<sup>π</sup>* for the linear element in one dimensional case [27]. Furthermore, it will be readily seen that the same constant can be taken for the two dimensional bilinear element from the consideration on the proof of Theorem 5.1 in [27]. To verify the existence of a solution of (26) in a computer, we use the fixed point formulation.

First, note that, for each *w* ∈*V*, there exists a unique *F w*ð Þ∈*V* such that

$$(\nabla F(w), \nabla v) = (\ulcorner f(w), v), \quad \forall v \in V,\tag{30}$$

which also implies that

$$\begin{cases} -\Delta F(w) = f(w) \quad \text{in } \Omega, \\\ F(w) = 0 \quad \text{on } \partial\Omega. \end{cases} \tag{31}$$

Then the map *F* : *V* ! *V* is compact. By (30), the problem (26) is equivalent to finding *w* ∈*V* such that

$$a(w, v - w) \ge a(F(w), v - w), \qquad \forall v \in K. \tag{32}$$

Using the definition (27) and (32), we now have the following fixed point problem for the compact operator *PKF*.

$$\text{Find } \exists w \in V \text{ such that } w = P\_K F(w). \tag{33}$$

#### *3.1.1.2 Verification condition*

We introduce two concepts, rounding and rounding error, which enable us to deal with the infinite dimensional problem by finite procedures, that is, in a computer.

Now we define the dual cone of *Kh* by

$$K\_h^\* = \{ w \in V : a(w, v) \le 0, \quad \forall v \in K\_h \},$$

and note that *K* <sup>∗</sup> *<sup>h</sup>* is also closed convex cone in *V* with vertex at 0 which is the only point common to *Kh* and *K* <sup>∗</sup> *<sup>h</sup>* . From (28) it follows that *K* <sup>∗</sup> *<sup>h</sup>* is the set of points whose projections into *Kh* is 0.

**Lemma 6.** *Any w* ∈*V can be uniquely decomposed into the sum of two orthogonal elements. That is,*

$$\boldsymbol{w} = \boldsymbol{P\_{K\_{h}}} \boldsymbol{w} \oplus \left(\boldsymbol{I} - \boldsymbol{P\_{K\_{h}}}\right) \boldsymbol{w} = \boldsymbol{P\_{K\_{h}}} \boldsymbol{w} \oplus \boldsymbol{P\_{K\_{h}^\*}} \boldsymbol{w}.$$

*Here,* ⊕ *denotes the sum of two orthogonal elements in the sense of V.*

For any *w* ∈*V*, we now define the rounding *R P*ð Þ *KF w*ð Þ ∈*Kh* by the solution of the following problem:

$$(a(R(P\_KF(w)), \upsilon\_h - R(P\_KF(w)))) \ge (\ f(w), \upsilon\_h - R(P\_KF(w))), \quad \forall \upsilon\_h \in K\_h.$$

Next, for any subset *W* ⊂*V*, we define the rounding *R P*ð Þ *KFW* ⊂*Kh* by

$$R(P\_K W W) = \{ w\_h \in K\_h : w\_h = R(P\_K F(w)), \quad w \in W \}.$$

Usually, *R P*ð Þ *KFW* is enclosed and represented as a linear conbination of the base functions in *Vh* with interval coefficients.

Moreover, for *W* ⊂*V*, we define *RE P*ð Þ *KFW* , the rounding error of *PKFW*, as a subset of *K*<sup>∗</sup> *<sup>h</sup>* , that is,

$$RE(P\_KFW) = \left\{ v \in K\_h^\* \; : \; \|v\|\_V \le C\_0 h \|f(W)\|\_{L^2} \right\},\tag{34}$$

where

$$\|f(\mathcal{W})\|\_{L^2} \equiv \sup\_{w \in \mathcal{W}} \|f(w)\|\_{L^2}.$$

Here, *C*<sup>0</sup> � *C*1*C*2, where *C*<sup>1</sup> is the same positive constant as in (29), and *C*<sup>2</sup> is determined by the following regularity estimate for the solution to (24) of the form

$$\|u\|\_{H^2} \le \mathcal{C}\_2 \|\underline{\mathbf{g}}\|\_{L^2}.\tag{35}$$

Thus we may take as *C*<sup>2</sup> ¼ 1 for the present case from (25). Then, we have

$$P\_K F(w) - R(P\_K F(w)) \in RE(P\_K F(w)), \quad \forall w \in W.$$

Therefore, the following verification condition is obtained by Schauder's fixed point theorem.

**Lemma 7.** *If there exists a nonempty, bounded, convex, and closed subset W* ⊂*K such that*

$$R(P\_K \text{FW}) \oplus RE(P\_K \text{FW}) \subset \mathcal{W},\tag{36}$$

then there exists a solution of *w* ¼ *PKF w*ð Þ in *W*.

We sometimes refer the above set *W* as *a candidate set*, which we generate in computer so that it satisfies the condition (36).

*Numerical Verification Method of Solutions for Elliptic Variational Inequalities DOI: http://dx.doi.org/10.5772/intechopen.101357*

#### *3.1.1.3 Verification procedures*

We describe the method to find a set *W* satisfying (36) in th ebelow. Consider the following approximate solution *wh* ∈*Kh* of (24):

$$a(w\_h, v\_h - w\_h) \ge (\mathbf{g}, v\_h - w\_h), \quad \forall v\_h \in \mathcal{K}\_h, \quad w\_h \in \mathcal{K}\_h. \tag{37}$$

Since the bilinear form *a*ð Þ �, � is symmetric, (37) is reduced to the quadratic programming problem:

$$\min\_{\nu} \left[ \frac{1}{2} a(\nu, \nu) - (\mathbf{g}, \nu) \right]. \tag{38}$$

Let *ϕ<sup>j</sup>* n o *j*¼1⋯*M* be a basis of *Vh* with usual linear functions such that *ϕj* ð Þ *x* ≥ 0, ∀*x*∈ Ω and satisfying

$$\phi\_j(\mathbf{x}\_i) = \begin{cases} \mathbf{1}, & i = j, \\ \mathbf{0}, & i \neq j, \end{cases}$$

where *xi* is an interior node of the finite element mesh. Then (38) reduces to the following vector form:

$$\min\_{w \ge 0} \left[ \frac{1}{2} w' Dw - P'w \right],\tag{39}$$

where *w* ≥ 0 means the componentwise relation. Here, *D* ≔ *dij* � � <sup>1</sup>≤*i*,*<sup>j</sup>* <sup>≤</sup> *<sup>M</sup>* with *dij* ¼ ∇*ϕi*, ∇*ϕ<sup>j</sup>* � �, and *<sup>w</sup>* is the coefficient vector with *<sup>ϕ</sup><sup>j</sup>* n o of the function *<sup>v</sup>* in (38). Also, *P* ≔ *g*, *ϕ<sup>j</sup>* � � � � <sup>1</sup>≤*<sup>j</sup>* <sup>≤</sup> *<sup>M</sup>*.

Furthermore, we define for any *α* ∈*R*þ, nonnegative real number, we set

$$[a] \equiv \{ \phi \in K\_h^\* \; ; \quad \|\phi\|\_V \le a \}.$$

Then, for a given candidate set *W* ¼ *Wh* ⊕ ½ � *α* with *Wh* ⊂*Kh*, the computation of the rounding *R P*ð Þ *KFW* reduces to enclose an interval vector *Z* ¼ *Z <sup>j</sup>* � � and *<sup>Y</sup>* <sup>¼</sup> *<sup>Y</sup> <sup>j</sup>* � � satisfying the following nonlinear system of equations [27]:

$$\begin{cases} Y - DZ = -\left( f(W), \phi\_j \right), \quad 1 \le j \le M, \\ Y\_j Z\_j = 0, \quad 1 \le j \le M. \end{cases} \tag{40}$$

Here, *f W*ð Þ, *ϕ<sup>j</sup>* � � is evaluated as an interval *<sup>B</sup> <sup>j</sup>* such that *f w*ð Þ, *<sup>ϕ</sup> <sup>j</sup>* � �j*<sup>w</sup>* <sup>∈</sup>*<sup>W</sup>* n o⊂*<sup>B</sup> <sup>j</sup>:* In order to solve (40) with guaranteed accuracy, we use some interval approaches for nonlinear system of equations [19, 20]. Thus, using the solution of (40), we can enclose the set *R P*ð Þ *KFW* in (36). Combining this with (34), we can successfully compute the left-hand side of (36) for any candidate set *W* ¼ *Wh* ⊕ ½ � *α* .

Thus we can present a computational verification condition. In the actual computation, we use an iterative procedure with *δ*-*inflation* technique to find the set *W* satisfying (36). Several numerical examples for verification are presented in [27] for one dimensional problem using linear finite element.

#### *3.1.2 Non-homogeneous case*

In this subsection, we consider the two-dimensional case. In order to verify solutions numerically, it is necessary to determine some constants that appear in the a priori error estimates. For the non-homogeneous case, we define *K* ≔ f g *v*∈*V* : *v*≥*ψ* a*:*e*:* on Ω , where *ψ* is a given *H*<sup>2</sup> ð Þ Ω function such that *ψ* ≤0 on ∂Ω and is not identically equal to 0. Let Ω be a square with side 1 and let T *<sup>h</sup>* be the uniform triangulation of Ω. We introduce Σ*<sup>h</sup>* ¼ *p*; *p*∈ Ω, *p* is a vertex of *T* ∈T *<sup>h</sup>* � � and define the approximate *Vh* of *H*<sup>1</sup> <sup>0</sup>ð Þ <sup>Ω</sup> by *Vh* ¼ f*vh*; *vh* <sup>∈</sup> *<sup>H</sup>*<sup>1</sup> <sup>0</sup>ð Þ <sup>Ω</sup> <sup>∩</sup>*C*<sup>0</sup> <sup>Ω</sup> � �, *vh* � � *<sup>T</sup>* ∈ *P*1, ∀*T* ∈T *<sup>h</sup>*g*:* Here, *vh*j *T* denotes the restriction of *vh* to *T* and *P*<sup>1</sup> representing the space of polynomials in two variables of degree ≤1. It is then quite natural to approximate *K* by

$$K\_h = \{ v\_h \in V\_h \colon v\_h(p) \ge \varphi(p), \forall p \in \Sigma\_h \} \dots$$

Note that, in general, *Kh* 6¼ *Vh* ∩ *K*. Then, *PK* and *PKh* are similarly defined as before, and we also have the constructive error estimates of the form, ∀*vh* ∈*Kh* and ∀*v*∈*K*,

$$\|\mathfrak{u}\_h - \mathfrak{u}\|\_{H^1\_0(\Omega)} \le \mathcal{C}(\mathfrak{g}, \boldsymbol{\varphi}, h), \tag{41}$$

where,

$$\mathcal{C}(\mathfrak{g},\boldsymbol{\nu},h) \leq \sup\_{\boldsymbol{\mathcal{g}} \in L^2(\mathfrak{Q})} \sqrt{(\mathbb{0}.494)^2 h^2 |\boldsymbol{u}|\_{\boldsymbol{H}^2}^2 + 2\left(\|\boldsymbol{\mathcal{g}}\|\_{\boldsymbol{L}^2} + \|\boldsymbol{A}\boldsymbol{u}\|\_{\boldsymbol{L}^2}\right) \left((\mathbb{0}.494)^2 h^2 |\boldsymbol{u}|\_{\boldsymbol{H}^2} + \mathsf{G} h^2 |\boldsymbol{\mathcal{y}}|\_{\boldsymbol{H}^2}\right)}.$$

We provide a numerical example of verification in the two-dimensional case according to the procedures described in the previous section. Let Ω ¼ ð Þ� 0, 1 ð Þ 0, 1 . We consider the case *f u*ð Þ¼ *Ku* þ sin *πx* sin 2*πy* and *ψ* ¼ sin *πx* sin *πy*. For simplicity, we only consider the uniform mesh here. First, we divide the domain into small triangles with a uniform mesh size *h* and choose the basis of *Vh* as the pyramid functions.

The execution conditions are as follows (**Figures 1**–**3**):

*K* ¼ 0*:*1, dim*Vh* ¼ 10 Obstacle function *ψ* ¼ sin *π x* sin *π y* the outline of *ψ* is shown in Figure 1*:* Initial value : *u*ð Þ <sup>0</sup> *<sup>h</sup>* ¼ Galerkin approximation, *α*<sup>0</sup> ¼ 0 the outline of *u*ð Þ <sup>0</sup> *<sup>h</sup>* is shown in Figure 2*:* Illustration of contact zone between obstacle and approximate solution is shown in Figure 3*:* Extension parameters : *<sup>δ</sup>* <sup>¼</sup> <sup>10</sup>�<sup>5</sup> *:*

Results are as follows:

$$\begin{aligned} & \text{Iteration numbers for verification} \quad : 2\\ & H\_0^1(\Omega) - \text{error bound} \quad : 0.15437\\ & \text{Maximum width of coefficient intervals in } \left\{ A\_j^{(N)} \right\} = 0.00001. \end{aligned}$$

Detailed arguments and with numerical examples are presented in [42].

#### *Numerical Verification Method of Solutions for Elliptic Variational Inequalities DOI: http://dx.doi.org/10.5772/intechopen.101357*

**Figure 1.** *Obstacle function ψ.*

**Figure 2.** *Approximate solution u*ð Þ <sup>0</sup> *h .*

**Figure 3.** *Illustration of the contact zone.*

#### *3.1.3 A Newton-type verification method*

The idea of the enclosure method for solutions of obstacle problems is based upon simply sequential iterations for the original fixed point operator *PKF*. Therefore, it is difficult to apply the method to the problem of which associated operator is not retractive in a neighborhood of the solution. In order to overcome such a difficulty, we introduce an another formulation using a Newton-like operator. The essential point is the way to devise the Newton-like operator for a kind of nondifferentiable map which defines the original problem.

To formulate a Newton-type verification condition, we need a Fréchet derivative of the operator *PKF*. However, *PKF* is not Fréchet differentiable at all. Therefore, we define the approximate Fréchet-like derivative *D*e*KF u*ð Þ*<sup>h</sup>* on *Vh* for some *uh* ∈ *Kh* instead of the Fréchet derivative. Assume that *ϕ<sup>j</sup>* n o *j*¼1⋯*M* is a basis of *Vh*, where *M* ¼ dim*Vh*, such that *ϕ<sup>j</sup>* ð Þ *x* ≥ 0 on Ω and satisfying

$$\phi\_j(\mathbf{x}\_i) = \begin{cases} \mathbf{1}, & i = j, \\ \mathbf{0}, & i \neq j, \end{cases}$$

where *xi* is an interior node of the finite element mesh. And, for *vh* ∈*Vh*, we represent it such as

$$v\_h = \sum\_{j=1}^{M} v\_{hj} \phi\_j.$$

Here, *vhj* � � *<sup>j</sup>*¼1,⋯,*<sup>M</sup>* is called as the coefficient vector of *vh*. Now we take a fixed subset *N*<sup>0</sup> ⊂ f g 1, 2, ⋯, *M* , define *Vh*,*N*<sup>0</sup> , the closed subspace of *Vh*, by

$$V\_{h,N\_0} = \left\{ v\_h | v\_h \in V\_h, v\_{hj} = 0 \mid \text{for } \ j \notin N\_0 \right\}.$$

*Numerical Verification Method of Solutions for Elliptic Variational Inequalities DOI: http://dx.doi.org/10.5772/intechopen.101357*

And let *Ph*,*N*<sup>0</sup> be a *H*<sup>1</sup> 0-projection from *V* onto *Vh*,*N*<sup>0</sup> defined by

*a u* � *Ph*,*N*<sup>0</sup> ð *u*, *v*Þ ¼ 0, ∀*v*∈*Vh*,*N*<sup>0</sup> , *Ph*,*N*<sup>0</sup> *u*∈*Vh*,*N*<sup>0</sup> *:*

In order to define *D*e*KF u*ð Þ*<sup>h</sup>* : *Vh* ! *Vh*,*N*<sup>0</sup> , we differentiate the first equation of (40) in *W* at *W* ¼ *uh* to get, for arbitrary *δ*∈*Vh*,

$$
\partial Y^\* - D \partial Z^\* = -\left\{ \left( f'(u\_h)\delta, \phi\_j \right) \right\}\_{1 \le j \le M}.\tag{42}
$$

Here, *<sup>∂</sup><sup>Y</sup>* <sup>∗</sup> <sup>¼</sup> *<sup>Y</sup>*<sup>e</sup> <sup>∗</sup> *j* � � 1≤*j* ≤ *M* and *<sup>∂</sup><sup>Z</sup>* <sup>∗</sup> <sup>¼</sup> *<sup>Z</sup>*<sup>e</sup> <sup>∗</sup> *j* � � 1≤*j*≤ *M* , where *<sup>Y</sup>*<sup>e</sup> <sup>∗</sup> *<sup>j</sup>* ¼ 0 for *j*∈ *N*<sup>0</sup> and *Z*e ∗ *<sup>j</sup>* ¼ 0 for *j* ∉ *N*0, respectively.

Then we define the approximate Fréchet-like derivative of *PKF u*ð Þ at *u* ¼ *uh*, as the linear map *D*e*KF u*ð Þ*<sup>h</sup>* : *Vh* ! *Vh*,*N*<sup>0</sup> such that, for each *δ*∈*Vh*,

$$
\widetilde{D}\_K F(\mathfrak{u}\_h)(\delta) \coloneqq \sum\_{j=1}^M \widetilde{Z\_j}^\* \phi\_j.
$$

We now assume that.

**A4.** The restriction to *Vh*,*N*<sup>0</sup> of the operator *Ph*,*N*<sup>0</sup> *I* � *D*e*KF u*ð Þ*<sup>h</sup>* h i : *Vh* ! *Vh*,*N*<sup>0</sup> has the inverse operator

$$\left[P\_{h,N\_0} - \widetilde{D}\_K F(\mathfrak{u}\_h)\right]\_h^{-1} : \mathcal{V}\_{h,N\_0} \to \mathcal{V}\_{h,N\_0}.$$

Here, *I* means the identity map on *Vh*.

By using the above approximate Fréchet-like derivative, we define the Newtonlike operator *Nh* : *V* ! *Vh* by

$$N\_h(w) \equiv P\_{K\_h} w - \left[P\_{h,N\_0} - \tilde{D}\_K F(u\_h)\right]\_h^{-1} P\_{h,N\_0} \left(P\_{K\_h} - P\_{K\_h} P\_K F\right)(w) ./\ .$$

Next we define the operator *T* : *V* ! *V* as follows:

$$T(\boldsymbol{w}) \equiv \boldsymbol{N}\_h(\boldsymbol{w}) + (\boldsymbol{I} - \boldsymbol{P}\_{\boldsymbol{K}\_h})\boldsymbol{P}\_K \boldsymbol{F}(\boldsymbol{w}).$$

Then *T* becomes a compact map on *V* and it follows the fixed point problem *w* ¼ *PKFw* is equivalent to *w* ¼ *T w*ð Þ. Detailed arguments and with numerical examples are presented in [35].

#### **3.2 Elasto-plastic torsion problems**

In this subsection, we consider an enclosure metnod of solutions for elastoplastic torsion problems governed by an elliptic variational inequalities [25, 32, 33]. The nonlinear elasto-plastic torsion problem is defined as the same type elliptic variational inequalities as (26) with

$$K \coloneqq \left\{ v \in H\_0^1(\Omega); \left| \nabla v \right| \; \le 1 \right\} \quad \text{a.e. on } \left\| \Omega \right\|. \tag{43}$$

As is well known [56, 58], two sub-domains Ω*<sup>p</sup>* and Ω*<sup>e</sup>* defined by

$$\mathfrak{Q}\_p = \{ \mathfrak{x}; \mathfrak{x} \in \mathfrak{Q}, |\nabla u| = 1 \},$$

and

$$
\mathfrak{Q}\_{\mathfrak{e}} = \mathfrak{Q} \backslash \mathfrak{Q}\_p = \{ \mathfrak{x}; \mathfrak{x} \in \mathfrak{Q}, \left| \nabla u \right|\_{} < 1 \},
$$

correspond to the plastic and elastic regions, respectively. The elastic region Ω*<sup>e</sup>* and the plastic region Ω*<sup>p</sup>* are not known beforehand and should be determined, therefore *∂*Ω*<sup>e</sup>* ∩ *∂*Ω*<sup>p</sup>* is actually the free boundary of the problem (26). The problem (26) has been formulated as the problem of finding *u* satisfying

$$\begin{cases} -\Delta u = f(u) \quad \text{in} \quad \Omega\_{\varepsilon}, \\ \quad |\nabla u| \quad = 1 \quad \text{in} \quad \Omega\_{p}, \\ \quad u = 0 \quad \text{on} \quad \partial\Omega. \end{cases} \tag{44}$$

The finite dimensional convex subset *Kh* is also defined similarly as before:

$$K\_h \coloneqq V\_h \cap K = \{ \upsilon\_h \mid \upsilon\_h \in V\_h, \ |\nabla \upsilon\_h| \le 1 \text{ a.e. on } \Omega \}. \tag{45}$$

In order to formulate the verification procedure, we need a verified computational method for solving the finite dimensional part (rounding) and a constructive estimates for infinite dimensional part (rounding error) as in the previous subsection.

Following [49, 56], we define the Lagrangian functional L associated with (1) by

$$\mathcal{L}(\boldsymbol{\nu},\mu) = \frac{1}{2} \int\_{\Omega} \left| \nabla \boldsymbol{\nu} \right|^2 d\boldsymbol{x} - (\boldsymbol{g}, \boldsymbol{\nu}) + \frac{1}{2} \int\_{\Omega} \mu \left( \left| \nabla \boldsymbol{\nu} \right|^2 - \mathbf{1} \right) d\boldsymbol{x} \dots$$

It follows, from [49, 56], that if *<sup>L</sup>* has a saddle point f g *<sup>u</sup>*, *<sup>λ</sup>* <sup>∈</sup> *<sup>H</sup>*<sup>1</sup> <sup>0</sup>ð Þ� <sup>Ω</sup> *<sup>L</sup>*<sup>∞</sup> þð Þ Ω , then *u* is a solution of (1), where *L*<sup>∞</sup> þð Þ¼ <sup>Ω</sup> *<sup>q</sup>*∈*L*<sup>∞</sup> f g ð Þ <sup>Ω</sup> ; *<sup>q</sup>*<sup>≥</sup> <sup>0</sup> *<sup>a</sup>:e: in* <sup>Ω</sup> *:* We use the Uzawa algorithm to solve (1). Thus we can claculate the rounding *R P*ð Þ *KF W*ð Þ , for a candidate set *W*, by solving the following problem with guaranteed error bounds:

$$\begin{cases} \text{Find } \left\{ \boldsymbol{u}\_{h}, \boldsymbol{\lambda}\_{h} \right\} \in \boldsymbol{K}\_{h} \times \boldsymbol{\Lambda}\_{h} & \text{such that} \\ \boldsymbol{\lambda}\_{h} = \max \left[ \boldsymbol{\lambda}\_{h} + \rho \left( |\nabla \boldsymbol{u}\_{h}|^{2} - 1 \right), \mathbf{0} \right] & \text{with } \rho > 0. \\ \int\_{\Omega} (\boldsymbol{1} + \boldsymbol{\lambda}\_{h}) \nabla \boldsymbol{u}\_{h} \cdot \nabla \boldsymbol{v}\_{h} d\mathbf{x} = (\boldsymbol{f}(\mathcal{W}), \boldsymbol{v}\_{h}), \forall \boldsymbol{v}\_{h} \in \boldsymbol{V}\_{h}, \quad \boldsymbol{u}\_{h} \in \boldsymbol{V}\_{h}, \end{cases} \tag{46}$$

The problem (46) can be formulated as a system of nonlinear and nonsmooth (nondifferentiable) equations. A verification method for nonsmooth equations by a generalized Krawczyk operator is studied in [1, 55]. We briefly describe the method presented by [55] in the below.

We consider the following equivalent system of nonlinear(and nondifferentiable) equation to (46) for a fixed *w* ∈*W*

$$H(\mathbf{x}) = \mathbf{0}.\tag{47}$$

Here, we assume that *<sup>H</sup>* : *<sup>R</sup><sup>n</sup>* ! *<sup>R</sup><sup>n</sup>* is locally Lipschitz continuous. The equivalence means that *x*<sup>∗</sup> solves (46) if and only if *x*<sup>∗</sup> solves (47). The method is based on the mean value theorem for local Lipschitz functions of the form

$$H(\mathfrak{x}) - H(\mathfrak{y}) \in co\partial H([\mathfrak{x}]) (\mathfrak{x} - \mathfrak{y}), \text{ for all } \mathfrak{x}, \mathfrak{y} \in [\mathfrak{x}], \mathfrak{y}$$

*Numerical Verification Method of Solutions for Elliptic Variational Inequalities DOI: http://dx.doi.org/10.5772/intechopen.101357*

where [x] stands for an interval vector, "co" denotes the convex hull, and *∂H* the generalized Jacobian in Clarke's sense [57], which is also considered as a slope function, and

$$co\partial H([\varkappa]) \coloneqq co\{V \in \partial H(\varkappa);\ \varkappa \in [\varkappa]\}.$$

Let *L*½ � *<sup>x</sup>* � � be an interval matrix such that *co∂H x*ð Þ ½ � <sup>⊆</sup> *<sup>L</sup>*½ � *<sup>x</sup>* � �*:* Then for any *<sup>x</sup>*, *<sup>y</sup>*∈½ � *<sup>x</sup>* <sup>⊆</sup>*R<sup>n</sup>* it holds that *H x*ð Þ� *H y*ð Þ<sup>∈</sup> *<sup>L</sup>*½ � *<sup>x</sup>* � �ð Þ *<sup>x</sup>* � *<sup>y</sup>* .

Then an interval operator for nonsmooth equations is defined by

$$G(\mathbf{x}, A, [\mathbf{x}]) \coloneqq \mathbf{x} - A^{-1} H(\mathbf{x}) + \left(I - A^{-1} [L\_{[\mathbf{x}]}]\right) ([\mathbf{x}] - \mathbf{x}).\tag{48}$$

The mapping *G x*ð Þ , *A*, ½ � *x* is called a generalized Krawczyk operator. Therefore, the verification condition of solutions for (46) in ½ � *x* is given by

$$G(x, A, [\chi]) \subseteq [\chi] \subset D.$$

Thus, we can compute the solution of (46) with guaranteed accuracy. That is, we can enclose the rounding *R P*ð Þ *KF U*ð Þ . On the other hand, in order for the calculation of the rounding error *RE P*ð Þ *KF U*ð Þ , the similar arguements can also be applied for one dimensional problem. Actually, we can prove that the same constant *<sup>C</sup>*<sup>0</sup> <sup>¼</sup> ffiffi 5 p *<sup>π</sup>* is also valid for the present problem in one dimensinal case, which implies that we can give a verification procedure besed on the same principle as before [25, 32, 33]. In [33], we extended the approach to the numerical proof of existence of solutions for elasto-plastic torsion problems as well as gave a numerical example for one dimensional case. The verification method in [33] is based on the generalized Krawczyk operator for solving a system of nonsmooth (nondifferentiable) equations. In order to use the generalized Krawczyk operator, we need to calculate the Jacobian. In that case, we need some complicated techniques. However, in many cases, calculating the generalized Jacobian is very difficult. To overcome such difficulties, we proposed a numerical verification method without using the generalized Krawczyk operator. This method is attractive, since calculating the generalized Jacobian is not required in the computational performance. Furthermore, up to know, our verification methods are mainly based on the enclosure of solutions in the sense of *L*<sup>2</sup> or *H*<sup>1</sup> norms. We considered a numerical verification method with guaranteed *L*<sup>∞</sup> error bounds for the solution of elasto-plastic torsion problem.

#### **3.3 Simplified Signorini problems**

A simplified Signorini problem is also given by the elliptic variational inequalities of the form (26) with

$$K \coloneqq \left\{ v \in H\_0^1(\Omega) ; \quad v \ge 0 \quad \text{on} \quad \partial \Omega \Big| \right\} \tag{49}$$

and

$$\mathcal{A}(u,v) = \int\_{\Omega} \nabla u \cdot \nabla v d\mathfrak{x} + \int\_{\Omega} u v d\mathfrak{x}.\tag{50}$$

where

$$
\nabla u \cdot \nabla v = \frac{\partial u}{\partial \mathbf{x}\_1} \frac{\partial v}{\partial \mathbf{x}\_1} + \frac{\partial u}{\partial \mathbf{x}\_2} \frac{\partial v}{\partial \mathbf{x}\_2} \dots
$$

As well known, the solution *u* of this elliptic variational inequalities can be characterized as a solution of the following free boundary problem finding *u* and two subsets Γ<sup>0</sup> and Γ<sup>þ</sup> such that Γ0∪Γ<sup>þ</sup> ¼ ∂Ω and Γ<sup>0</sup> ∩ Γ<sup>þ</sup> ¼ Ø

$$\begin{cases} -\Delta u + u = f(u) \quad \text{in} \quad \Omega, \\\quad u = 0 \quad \text{on} \quad \Gamma\_0, \frac{\partial u}{\partial n} \ge 0 \quad \text{on} \quad \Gamma\_0, \\\quad u > 0 \quad \text{on} \quad \Gamma\_+, \frac{\partial u}{\partial n} = 0 \quad \text{on} \quad \Gamma\_+, \end{cases} \tag{51}$$

where *<sup>∂</sup> <sup>∂</sup><sup>n</sup>* the outer normal derivative on ∂Ω. In the present case, the approximation subspace *Kh* is taken as

$$K\_h \coloneqq V\_h \cap K = \{v\_h \mid \quad v\_h \in V\_h, \quad v\_h \ge 0 \quad \text{on} \quad \partial\Omega\}.\tag{52}$$

For a candidate set *W*, the computation of rounding *R P*ð Þ *KF W*ð Þ is also reduced to the quadratic programming problem as in the Section 3.1 [56].

Since the constant *C*<sup>2</sup> in (25) is easily estimated as *C*<sup>2</sup> ¼ 1, the standard approximation property of the interpolation by *Kh* gives a constructive error estimates to compute the rounding error *RE P*ð Þ *KF W*ð Þ . For a simplified Signorini problem [43], we constructed a computing algorithm which automatically encloses the solution within guaranteed error bounds. In particular, the method proposed in [43] enables us to verify the free boundary of a simplified Signorini problem, which has been impossible so far. Concerning the numerical verification of solutions for elliptic variational inequalities, we would like to mention that the inclusion method described in this article can be applied to the solution of the elliptic variational inequalities on large space domains.

#### **3.4 Some other problems**

In this subsection, we show that our idea of verification method can also be applied to the elliptic variational inequalities of the second kind.

Now, we define the functional *j v*ð Þ¼ <sup>Ð</sup> <sup>Ω</sup>∣∇*v*∣*dx:* We consider the following problem of the flow of a viscous plastic fluid in a pipe:

$$\begin{cases} \text{Find } \ u \in H\_0^1(\Omega) \quad \text{such that} \\ a(u, v - u) + j(v) - j(u) \ge (\ f(u), v - u), \quad \forall v \in H\_0^1(\Omega). \end{cases} \tag{53}$$

As in the previous section, we consider the following auxiliary problem associated with (53) for a given *g* ∈ *L*<sup>2</sup> ð Þ Ω :

$$a(u, v - u) + j(v) - j(u) \ge (\mathbf{g}, v - u), \forall v \in H\_0^1(\Omega), u \in H\_0^1(\Omega). \tag{54}$$

By the well known result, we have the following lemma.

**Lemma 8.** There exists a unique solution *u*∈ *H*<sup>1</sup> <sup>0</sup>ð Þ <sup>Ω</sup> <sup>∩</sup> *<sup>H</sup>*<sup>2</sup> ð Þ Ω of (54) for any *g* ∈*L*<sup>2</sup> , such that

$$\|\mathfrak{u}\|\_{H^{2}(\mathfrak{\Omega})} \leq \hat{C} \|\mathfrak{g}\|\_{L^{2}(\mathfrak{\Omega})}.$$

When we denote the solution *u* of (54) by *u* ¼ *Ag* and define the composite map *F* on *H*<sup>1</sup> <sup>0</sup>ð Þ Ω by *F u*ð Þ� *Af u*ð Þ, which is a little bit of different from the previously appeared symbol *F* in Section 2, we have.

*Numerical Verification Method of Solutions for Elliptic Variational Inequalities DOI: http://dx.doi.org/10.5772/intechopen.101357*

**Theorem 9.** *F* is compact on *H*<sup>1</sup> <sup>0</sup>ð Þ Ω and the problem (53) is equivalent to the fixed point problem

$$u = F(u).$$

*Proof.* First, for a bounded subset *U* ⊂*L*<sup>2</sup> ð Þ <sup>Ω</sup> , we show that *AU* <sup>⊂</sup> *<sup>H</sup>*<sup>1</sup> <sup>0</sup>ð Þ Ω is relatively compact. Secondly, prove that *A* : *L*<sup>2</sup> ð Þ! <sup>Ω</sup> *<sup>H</sup>*<sup>1</sup> <sup>0</sup>ð Þ Ω is continuous. By Lemma 3, *AU* ⊂ *H*<sup>2</sup> ð Þ <sup>Ω</sup> <sup>∩</sup> *<sup>H</sup>*<sup>1</sup> <sup>0</sup>ð Þ <sup>Ω</sup> and *AU* is bounded in *<sup>H</sup>*<sup>2</sup> ð Þ Ω . Since *U* is bounded in *L*2 ð Þ Ω , by the Sobolev imbedding theorem, we have *AU* is relatively compact in *H*1 <sup>0</sup>ð Þ <sup>Ω</sup> *:* Next, for arbitrary *<sup>f</sup>* <sup>1</sup>, *<sup>f</sup>* <sup>2</sup> <sup>∈</sup> *<sup>L</sup>*<sup>2</sup> ð Þ Ω , setting *u*<sup>1</sup> ¼ *Af* <sup>1</sup> and *u*<sup>2</sup> ¼ *Af* <sup>2</sup>, by using (54), we obtain

$$a(u\_1, u\_2 - u\_1) + j(u\_2) - j(u\_1) \ge \left(\begin{array}{c} f\_1, u\_2 - u\_1 \end{array} \right),$$

$$a(u\_2, u\_1 - u\_2) + j(u\_1) - j(u\_2) \ge \left(\begin{array}{c} f\_2, u\_1 - u\_2 \end{array} \right).$$

With the above inequalities, we obtain *a u*ð <sup>2</sup> � *u*1, *u*<sup>2</sup> � *u*1޼�*a u*ð Þþ 1, *u*<sup>2</sup> � *u*<sup>1</sup> *a u*ð Þ 2, *u*<sup>2</sup> � *u*<sup>1</sup> ≤*j u*ð Þ�<sup>2</sup> T *<sup>h</sup>* Hence, by the Poincaré inequality, we have

$$\left\|\|\boldsymbol{u}\_{2}-\boldsymbol{u}\_{1}\|\right\|\_{H^{1}\_{0}(\Omega)}^{2} \leq \left\|\boldsymbol{f}\_{2}-\boldsymbol{f}\_{1}\right\|\_{L^{2}(\Omega)}\left\|\boldsymbol{u}\_{2}-\boldsymbol{u}\_{1}\right\|\_{L^{2}(\Omega)} \leq \overline{C}\|\boldsymbol{f}\_{2}-\boldsymbol{f}\_{1}\|\_{L^{2}}\left\|\boldsymbol{u}\_{2}-\boldsymbol{u}\_{1}\right\|\_{H^{1}\_{0}(\Omega)}.$$

Therefore, we obtain

$$\|\mu\_2 - \mu\_1\|\_{H^1\_0(\Omega)} \le C \|f\_2 - f\_1\|\_{L^2(\Omega)}.$$

That is, *A* is Lipschitz continuous as a map *L*<sup>2</sup> ð Þ! <sup>Ω</sup> *<sup>H</sup>*<sup>1</sup> <sup>0</sup>ð Þ Ω *:* Hence *A* is compact. The latter half in the theorem is straightforward from the definition of *F*.

We now define the approximate problem corresponding to (54) as

$$a(u\_h, v\_h - u\_h) + j(v\_h) - j(u\_h) \ge (\mathbf{g}, v\_h - u\_h), \forall v\_h \in V\_h, u\_h \in V\_h. \tag{55}$$

In order to apply our verification method to enclose the solutions of (53), we need a guaranteed computation of the exact solution of the problem (55), *a rounding procedure*, as well as the constructive error estimates between the solution of (54) and (55), *rounding error estimates*.

A major difficulty in solving the problem (55) numerically is the processing of the nondifferentiable term *j u*ð Þ¼ <sup>Ð</sup> <sup>Ω</sup>∣∇*u*∣*dx:* One approach is the method of Lagrange multiplier on that term, whose continuous version is as follows [56].

Let us define <sup>Λ</sup> <sup>¼</sup> *<sup>q</sup>* <sup>j</sup> *<sup>q</sup>*<sup>∈</sup> *<sup>L</sup>*<sup>2</sup> ð Þ� <sup>Ω</sup> *<sup>L</sup>*<sup>2</sup> ð Þ <sup>Ω</sup> , <sup>j</sup>*q x*ð Þj <sup>≤</sup><sup>1</sup> *<sup>a</sup>:e: <sup>x</sup>*<sup>∈</sup> <sup>Ω</sup> � � with ∣*q x*ð Þ∣ ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi *q*1ð Þ *x* <sup>2</sup> <sup>þ</sup> *<sup>q</sup>*2ð Þ *<sup>x</sup>* 2 q *:* Then the solution *u* of (54) is equivalent to the existence of *q* satisfying

$$\begin{cases} a(u, v) + \int\_{\Omega} q \cdot \nabla v = (\mathsf{g}, v), \forall v \in H\_0^1(\Omega), u \in H\_0^1(\Omega), \\ q \cdot \nabla u = |\nabla u| \quad a.e. \ , q \in \Lambda. \end{cases} \tag{56}$$

Moreover, it is known that (56) is equivalent to the following problem:

$$\begin{cases} a(u,v) + \int\_{\Omega} q \cdot \nabla v = (\mathsf{g}, \nu), \forall \nu \in H\_0^1(\Omega), u \in H\_0^1(\Omega), \\\\ q = \frac{q + \rho \nabla u}{\sup(\mathsf{1}, |q + \rho \nabla u|)}. \end{cases} \tag{57}$$

Here *ρ* is a positive constant. Let T *<sup>h</sup>* be a triangulation of Ω, and let define *Lh* and <sup>Λ</sup>*<sup>h</sup>* (approximation of *<sup>L</sup>*<sup>∞</sup>ð Þ� <sup>Ω</sup> *<sup>L</sup>*<sup>∞</sup>ð Þ <sup>Ω</sup> and <sup>Λ</sup>, respectively) by

$$L\_h = \left\{ q\_h \vert q\_h = \sum\_{\tau \in \mathcal{T}\_h} q\_\tau \chi\_\tau, q\_\tau \in \mathbb{R}^2 \right\} \text{ and } \; \Lambda\_h = \Lambda \cap L\_h, \; \text{respectively},$$

where *χτ* is the characteristic function of *τ*.

Then our first purpose, computing the rounding *RF U*ð Þ, is to enclose the solution of the following approximation problem of (57):

$$\begin{cases} \boldsymbol{a}(\boldsymbol{u}\_{h}, \boldsymbol{v}\_{h}) + \int\_{\Omega} \boldsymbol{q}\_{h} \cdot \nabla \boldsymbol{v}\_{h} = (\boldsymbol{g}, \boldsymbol{v}\_{h}), \forall \boldsymbol{v}\_{h} \in \boldsymbol{V}\_{h}, \boldsymbol{u}\_{h} \in \boldsymbol{V}\_{h}, \\\\ \boldsymbol{q}\_{h} = \frac{\boldsymbol{q}\_{h} + \rho \nabla \boldsymbol{u}\_{h}}{\sup\left(\mathbf{1}, |\boldsymbol{q}\_{h} + \rho \nabla \boldsymbol{u}\_{h}|\right)}. \end{cases} \tag{58}$$

The Eq. (58) leads to a kind of finite dimensional, nonlinear but nondifferentiable problem. We use a slope function method proposed by Rump [18–20] to enclose the solutions of (58) with *g* ¼ *f W*ð Þ for a candidate set *W*. On the other hand, the rounding error *RE F U* ð Þ ð Þ can be computed by using the following constructive error estimates:

**Theorem 10.** Let *u* and *uh* be solutions of (54) and (55), respectively. If *g* ∈*L*<sup>2</sup> ð Þ Ω , then there exists a constant *C h*ð Þ such that

$$\|\mathfrak{u}\_h - \mathfrak{u}\|\_{H^1\_0(\Omega)} \le C(h) \|\mathfrak{g}\|\_{L^2(\Omega)}.$$

Here, we may take *C h*ð Þ¼ ffiffi 5 p *<sup>π</sup> h* for the linear element in one dimensional case, and *<sup>C</sup>* is also numerically estimated such that *C h*ð Þ<sup>≈</sup> *O h*<sup>1</sup> 2 � � for the two dimensional linear element. A proof of this theorem is described in Ryoo and Nakao [34]. Thus we can also implement the verification algorithm for the solution of (53) as in the previous section. For details on this subsection, please refer to Ref. [47].

#### **4. Conclusions**

We have surveyed numerical verification methods for differential equations, especially around partial differential equations, variational inequalities and the author's works. But the period of this research is shorter than the history of the numerical methods for differential equations by computer and we can say it is still in the stage of case studies. Indeed, recently, this kind of studies have been referred little by little for practical applications in PDEs and variational inequalities but there are many open problems to be resolve. Therefore, we can make no safe prediction that these approaches will grow into really useful methods for various kinds of equations and variational inequalities in mathematical analysis. Also, since the program description of the verification algorithm is very complicated in general, there is another problem like software technology associated with assurance for the correctness of the verification program itself. Actually, some of the mathematician would not give credit the computer assisted proof in analysis as correct as they believe the theoretical proof, which might cause a kind of seriously emotional problem in the methodology of mathematical sciences. And there is another difficulty from the huge scale of numerical computations which often exceed the capacity of the concurrent computing facilities.

*Numerical Verification Method of Solutions for Elliptic Variational Inequalities DOI: http://dx.doi.org/10.5772/intechopen.101357*

However, in the twenty-first century, the computing environment would make more and more rapid progress, which should be beyond conception in the present state. In any case, a realistic study for partial differential equations and variational inequalities should be the future subject of the numerical computations with guaranteed accuracy. The authors believe that numerical methods with guaranteed accuracy for differential equations and variational inequalities would highly improve the reliability in the numerical simulation of the complicated phenomena in both mathematical and engineering sciences.
