Iterative method for solving a problem with mixed boundary conditions for biharmonic equation arising in fracture mechanics

abstract: In this paper we consider a mixed boundary value problem for biharmonic equation of the Airy stress function which models a crack problem of a solid elastic plate. An iterative method for reducing the problem to a sequence of mixed problems for Poisson equations is proposed and investigated. The convergence of the method is established theoretically and illustrated on many numerical experiments.


Introduction
The solution of fourth order differential equations by their reduction to BVP for the second order equations, with the aim of using available efficient algorithms for the latter ones attracts attention from many researchers.Namely, for the biharmonic equation with the Dirichlet boundary condition, there is intensively developed the iterative method, which leads the problem to two problems for the Poisson equation at each iteration (see e.g.[10,12,14]).In 1992, Abramov and Ulijanova [1] proposed an iterative method for the Dirichlet problem for the biharmonic type equation, but the convergence of the method is not proved.In our previous works [3,5,7,8] with the help of boundary or mixed boundary-domain operators appropriately introduced, we constructed iterative methods for biharmonic and biharmonic type equations associated with the Dirichlet, Neumann or simple type of mixed boundary conditions.These iterative methods are originated from our earlier works [2,6].It should be said that the mentioned above problems are 66 Dang Q.A and Mai X. T.
reduced to sequences of second order problems with boundary conditions of only one type on the whole boundary, i.e., the boundary conditions are not mixed.Recently, in [9] we have developed the iterative method for a problem in rectangular domain with rather complicated mixed boundary conditions for biharmonic equation arising in nano physics [15].It leads to the solution of a sequence of problems for the Poisson equation with mixed boundary conditions.But these boundary conditions are weakly mixed in the sense that on each side of the rectangle there is only one type of conditions.This property does not cause difficulties when using the method of complete reduction [17] for solving difference equations for second order differential problems at each iteration.
In this work we develop our technique for a problem with more complicated mixed conditions for biharmonic equation, namely, we consider the following problem where Ω is the rectangle (−1, 1) × (0, 1), and [18] (see also [11], which deals with a two-dimensional solid elastic plate containing a single edge crack, subjected to a uniform inplane load normal to the two edges parallel to the crack, while the remaining edges are stress free.For the problem in general setting (1.1) -(1.3) we propose an iterative method which reduces it to a sequence of problems for the Poisson equation.The convergence of the method is established and performed numerical experiments confirm the efficiency of the method under investigation.

Description of method
First, we assume that the problem (1.1)-(1.3)has a unique solution and it is sufficiently smooth.
As usual, we set Then the problem (1.1)-(1.3) is reduced to the problem ) ( where ϕ as u is unknown function but it is related to u by the second condition in (1.2), i. e., by the relation Now we consider the following iterative process for finding ϕ and simultaneously for finding u: (i) Given ϕ (0) ∈ L 2 (Γ 1 ) , for example, ϕ (0) = 0 on Γ 1 ; (ii) Knowing ϕ (k) on Γ 1 (k = 0, 1, ...) solve consecutively two problems (2.4) (2.5) (iii) Compute the new approximation where τ is an iterative parameter to be chosen later.

Investigation of convergence
In order to investigate the convergence of the iterative process (2.4)-(2.6)firstly we rewrite (2.6) in the canonical form of two-layer iterative scheme [16]: Next, we introduce the operator B defined on boundary functions ϕ by the formula where u is found from the problems: (2.9) (2.10) The properties of the operator B will be investigated in the sequel.Now, let us return to the problem (2.1)-(2.2).We represent their solution in the form where u 1 , v 1 satisfy the problems (2.9)-(2.10)and u 2 , v 2 are the solutions of the problems (2.12) (2.13) According to the definition of the operator B we have where Thus, we have reduced the original problem (1.1)-( 1.3) to the operator equation (2.15), whose right hand side F is completely defined by the data f, g 1 , g 2 , g 3 , g 4 .
Proof.Indeed, if in (2.4), (2.5) we put where u 2 , v 2 are the solutions of Problems (2.14)-(2.15)then we get (2.20) From here we see that Therefore, taking into account the first relation in (2.18) and the above equality, from (2.7) we obtain (2.17).Thus, the proposition is proved.2 Proposition 2.1 enables us to lead the investigation of convergence of the ierative process (2.4)-(2.6) to the study of the iterative scheme (2.17).For this reason we need some properties of the operator B.
Proof.The linearity of B is obvious.To estiblish the other properties of B we consider the inner product (Bϕ, φ) for two arbitrary functions ϕ and φ in L 2 (Γ 1 ).Recall that the operator B acting on ϕ is defined by (2.8)-(2.10).Denote now by v and ū the solutions of (2.9) and (2.10), where instead of ϕ there stands φ.
We have Next, in view of v = ∆u, using the Green formula we have (2.24) It means that the operator B is symmetric.Besides, we have If (Bϕ, ϕ) = 0 then v = 0 almost everywhere in Ω, hence ϕ = v| Γ1 = 0. Thus, B is positive operator.Now, we prove the compactness of the operator B. In order to do this suppose that ϕ ∈ H s (Γ 1 ) with s ≥ 0. Then Problem (2.9) has a unique solution v ∈ H s+1/2 (Ω), and consequently, Problem (2.10) has a unique solution Thus, the proof of the proposition is complete.2 Before stating the result of convergence of the iterative process (2.4)-(2.6)we assume that the data functions f, g 1 , . . ., g 4 have needed smoothness so that the original problem (1.1)-(1.3)has a unique solution u ∈ H 5/2 (Ω).It is guaranteed if . Then, we can see that the function F defined by (2.16) belongs to H 3/2 (Γ 1 ).
We shall consider (2.15) as an operator equation in the space H = L 2 (Γ 1 ).

.25)
Proof.This theorem follows from Lemma A.1 in Appendix A of [9] due to the properties of symmetry, positivity and compactness of the operator B estiblished by Proposition 2.2. 2 It should be said that the determination or estimation of ||B|| is a difficult problem, but in Section 4 by experimental way we can find a interval of τ , for which the iterative process has good convergence.

On numerical realization of the iterative method
From the previous section we see that for realizing the iterative method it is required to solve consecutively two mixed BVPs (2.4) and (2.5).These BVPs are strongly mixed in the sense that the transmission of the Dirichlet and Neumann boundary conditions occurs at a inner point, namely at the middle of the bottom side of the rectangle.For solving these problems we use a domain decomposition method proposed in [4] which reduces the strongly mixed problem to a sequence of weakly mixed problems in subdomains in the sense that on each side of subdomains there is given boundary condition of only one type, either Dirichlet or Neumann type.
Below we briefly describe this domain decomposition method applied to the model problem where 2 (the same as in Figure 1).Denote two parts of the rectangle [−1, 1] × [0, 1] by Ω 1 and Ω 2 and their common boundary by S I .Besides, we denote the outward normal to the boundary of Ω i by ν i and the solution u of Problem (3.1) in Ω i by u i , i.e., u i = u| Ωi (i = 1, 2).
-6 The iterative process for finding u 1 and u 2 is described as follows: (i) Given ϕ (0) ∈ L 2 (S I ) , for example, ϕ (0) = 0 on S I ; (ii) Knowing ϕ (k) on S I (k = 0, 1, ...) solve consecutively two problems Iterative method for mixed BVP for biharmonic equation 73 ∆u (iii) Compute the new approximation where θ is an iterative parameter to be chosen appropriately.
Remark that Problems (3.2) and (3.3) are weakly mixed problems, where the Neumann boundary condition is prescribed on one side of the subdomains and the Dirichlet boundary condition is prescribed on other sides.In order to numerically solve these problems we discretize them on uniform grids by difference schemes of second order approximation obtained by a variational method.After that the system of difference equations are solved following the method of complete reduction with the complexity O(M N ln N ), where M, N are the number of grid nodes on the vertical and horizontal sides of the subdomains.Next, for computing the normal derivative in (3.4) we also use an approximate formula of second order error.We take in the formula (3.4) θ = 0.5 and carry out the iterative process (3.3)-(3.4)until max{ u ∞ } < ε , where ε is a given accuracy taken of the same order as O(h 2 ), h being the stepsize of the grid.
Below we report the results of using the above domain decomposition method for the numerical realization of the iterative process (2.4)- (2.6), where for computing the normal derivative in (2.6) we also use an approximate formula of second order error.

Numerical results
We perform some experiments for testing the convergence of the iterative process (2.4)-(2.6) in both two cases, where the exact solution of the problem (1.1)-(1.3) is known and unknown.
Example 1: Given function u(x, y) as the exact solution of the problem (1.1)-(1.3),calculate corresponding right hand side function f (x, y) and boundary conditions (1.2), (1.3), and then carry out the iterative process (2.4)-(2.6)until ||u (k) − u|| ∞ ≤ ε, where ε is the same given accuracy as in the previous section.The following functions are taken as the exact solutions of the problem The results of computation on the uniform grid of 65 × 65 nodes are given in Table 1, where K i (i = 1, ..., 5) is the number of iterations for achieving the exact solution u i (x, y) with the accuracy ε = 10 −3 .
Looking at Table 1 it appears that the result of computation for the function u = u 5 (x, y) = e x .sin y + sin x.e y is surprising.But this result is completely right because for the function we have ∆u = 0, g 4 = ∂∆u ∂ν = 0 and f = ∆ 2 u = 0.This implies that the solutions of the problems (2.4), (2.5) for k = 0 are v (0) = 0, u (0) = e x .sin y + sin x.e y .So, immediately we achieve the exact solution of the problem, and hence, K 5 =0.
Example 2: Given arbitrary data functions f, g 1 , g 2 , g 3 , g 4 in the problem (1.1)-(1.3),we perform the iterative process (2.4)-( 2 Below we report the results on convergence of the process for two collections of data functions: (i) f = xe −y + y.e x g 1 = sin x. sin y; g 2 = − sin x. sin y + ln(x + y + 2) g 3 = sin y + e y .sin x + x; g 4 = cos y + e x .sin y + x − y The results of computation on the uniform grid of 65×65 nodes are given in Table 2, where K i (i = 1, 2) is the number of iterations for achieving the accuracy ε = 10 −3 for the collections (i) and (ii).
Example 3: Consider the model fracture problem which is depicted in Figure 3 (see [11]).The graph of the obtained approximate solution is shown in Figure 4.

Concluding remarks
In this paper we investigated an iterative method for solving a boundary value problem for biharmonic equation with boundary conditions, which are of different types on different sides of a rectangle and the transmission of boundary conditions occurs not only in vertices but also in the middle point of a side of the rectangle.Due to the latter property we say that the problem is strongly mixed.Our iterative method based on iterative scheme for operator equation reduces the problem to sequence of strongly mixed problems for Poisson equation, and for the latter ones we apply a decomposition method recently developed by ourselves, which in its turn leads the problem to weakly mixed problems for Poisson equation.Finally, the difference method is used for realizing the latter problems.The convergence of the proposed iterative method at continuous level is proved and many different numerical experiments confirmed the efficiency of the method.

Figure 3 :
Figure 3: The model fracture problem

Figure 4 :
Figure 4: The graph of the approximate solution of the model fracture problem

Table 1 :
Convergence of the iterative process in Example 1

Table 2 :
Convergence of the iterative process in Example 2