Numerical solutions for optimal control problem governed by elliptic system on Lipschitz domains

ABSTRACT In this paper, an optimal boundary control problem for a distributed elliptic system on Lipschitz domains with boundary homogeneous Dirichlet conditions and independently with Neumann conditions is analysed. The necessary and sufficient optimality conditions for such problems with the quadratic cost functionals are obtained. A Jacobi spectral Galerkin method is introduced to develop a direct solution technique for the numerical solution of elliptic problems subject to Dirichlet and Neumann conditions in one and two dimensions. The numerical examples are included to demonstrate the validity and applicability of the techniques and comparison is made with the existing results. The method is easy to implement and yields very accurate results.


Introduction
The boundary value problems on Lipschitz domains have attracted attention in both pure and applied mathematics (see for example [1][2][3]). On the other hand, there are a lot of engineering problems such as point sensor placements at vertex points of the boundary in the mechanical engineering, corrosive engineering and spacecraft, point controllers at vertex points of the boundary in the stabilization of structural dynamics, etc., which require very careful local analysis around vertex points of the boundary (see [4] and the references therein). The advantage of the Sobolev spaces on Lipschitz domains is connected with the fact that the study of the elliptic partial differential equations becomes very simple and elegant within these spaces. Thus, Sobolev spaces on Lipschitz domains play a very important role in those studies. Most properties of Sobolev spaces on Lipschitz domains are rigorously proved (see e.g. [1,2]).
Optimal control problem for a distributed elliptic system has been a research area of its own (see for example [5][6][7][8]). The optimization problem associated with the optimal control of second-order distributed parameter systems defined on smooth domains of R n has been studied by Lions [5] and recently in [9]. For other cases, we refer for example to [10][11][12].
Recently, spectral methods have received considerable attention and become increasingly popular for solving different kinds of differential equations. This is definitively due to their superior accuracy for approximating the solutions of problems whose analytical solutions are smooth functions. In practice, this means that highly accurate solutions can be achieved with a small number of discretization nodes. Compared with algebraic convergence rates of the more traditional, finite element, volume, and difference methods spectral methods converge exponentially very fast [13][14][15][16][17].
The spectral methods are divided into Galerkin method, tau method, and collocation method [18][19][20]. For the convenience of different conditions, three classes of functions are always chosen as the basis for the spectral method: a power basis, a polynomial basis and a trigonometric basis. The important characteristic of the Galerkin method is that the expansion solution must individually satisfy the boundary conditions, whereas in collocation and tau methods, they are not necessary. Physical problems almost always be constrained by known boundary conditions which can be fully incorporated in a Galerkin method [21][22][23]. Finding a fast and precise numerical solutions for elliptic equations are often a demand in the process of solving problems of fluid dynamics [13,24]. The Galerkin method provides the most suitable numerical approximation for such problems. In his context, Doha and Bhrawy [25][26][27] introduced a new bases of some orthogonal polynomials for solving the two-dimensional Helmholtz and biharmonic problems by the diagonalization process.
In the present article, we consider a distributed elliptic system of second-order partial differential equation on an open set with a Lipschitz boundary in R n . The existence and uniqueness theorems for the Dirichlet and Neumann problems are proved. The optimal control is characterized by the adjoint problem, and then particular properties of optimal control are also obtained. Furthermore, we aim to develop some efficient spectral algorithms based on the Jacobi Galerkin method (JGM) to approximate the linear second-order elliptic problems such that they can be implemented efficiently. Finally, the accuracy of the proposed algorithms is demonstrated by test problems.
The paper is organized as follows. In Section 2, we recall some basic notation and definitions of the Sobolev space on Lipschitz domains. Section 3 is devoted to formulate the Dirichlet and Neumann problems for elliptic systems on Lipschitz domains. Consequently, the distributed control problems for the Dirichlet problem are studied. We also give the necessary and sufficient conditions for the control to be an optimal for the two problems. Finally, we study the boundary control problem for the Neumann problem. In Section 4, we introduce some algorithms for solving second-order elliptic problems in one and two dimensions by using Jacobi Galerkin spectral methods. Also, the proposed methods are applied to several examples.

Sobolev spaces on Lipschitz domains
We first give brief discussions on the Sobolev spaces on the whole domains of R n . As usual R n denotes the n-dimensional real Euclidian space of points x = (x 1 , . . . , x n ). By a multi-index α, we understand a vector α = (α 1 , . . . , α n ) of non-negative integer components α i . We set |α| = n i=1 α i . Let be a non-empty subset of R n . If f (x) is a sufficiently smooth function on , the partial derivative of f (x) is defined by In the following, we shall recall some basic function spaces: For any r ≥ 0, we set Obviously, C ∞ 0 ( ) is a real vector space and can be turned into a topological vector space by a proper topology. The space C ∞ 0 ( ) equipped with the following topology is denoted by D( ): a sequence of functions Evidently T f , thus defined, is a linear functional on D( ).
If f ∈ C |α| ( ), then integrating by parts |α|-times implies that there are no boundary terms, this is because φ has a compact support in and thus φ and its derivatives vanish near ∂ . This identity motivates the following definition of the derivative D α T of a distribution T ∈ D ( ): We are now able to define the concept of a function being the weak derivative of another function. Let α be a multi-index and let Then the function v is called the weak or distributional α th derivative of u and is denoted by With the notion of distributions, it is easy now to give a definition for the Sobolev spaces in a form needed in our purpose.
with a finite norm u m, given by where || · || is the usual L 2 -norm.

Lemma 2.1: The Sobolev space H m ( ) a Hilbert space together with the inner product
where (·, ·) is the scalar product in L 2 ( ).
Negative order Sobolev spaces: The closure of where the duality pairing on

Proposition 2.2 (Embedding of dual Sobolev spaces):
For any bounded domain in R n , one has the following chain of embeddings: Proof: The first embedding is trivial (by definition). The equality follows from the Riesz representation theorem and L 2 ( ) being a Hilbert space. For the last embedding observe that which implies The proposition is proved.
and an open neighbourhood U of p such that for some Lipschitz function χ. i.e. χ : R n−1 → R and satisfying the Lipschitz condition The smallest M in which (2.1) holds is called the bound of the Lipschitz constant. By choosing a finite open covering {U i } of ∂ , the Lipschitz constant for a Lipschitz domain is the smallest M such that the Lipschitz constant is bounded by M in every U i . We recall also that a Lipschitz function is almost everywhere differentiable (see e.g. [11]).

Example 2.3:
All the smooth domains are Lipschitz. In particular, a domain with smooth C 1 -boundary is Lipschitz. A very significant non-smooth domain which is Lipschitz is that every polygonal domain in R 2 or polyhedron in R 3 . A more interesting example is that every convex domain in R n is Lipschitz. A simple example of non-Lipschitz domain is two polygons touching at one vertex only.

Dirichlet problem for the elliptic system
Let be an open set in R n with a Lipschitz boundary . Consider the following distributed parameter elliptic system with the homogeneous Dirichlet boundary condition: where = n i=1 (2/x 2 i ) is a second-order self-adjoint elliptic partial differential operator maps H 1 0 ( ) onto H −1 ( ). We associate with the bilinear: and then Proof: It is well known that the ellipticity of is sufficient for the coerciveness of π(y, φ) on H 1 0 ( ). By virtue of (4) and (5), we have π(y, y) = n i,j=1 where C is a positive constant.
One can easily show the symmetric property of π(y, ·), π(y, φ) = π(φ, y), for all y, φ ∈ H 1 0 ( ). Under the above assumptions, we have the following existence and uniqueness theorem. Proof: By the assumptions, the bilinear form π(y, ·) and the linear functional F(φ) = (f , φ) satisfy the hypotheses of the Lax-Milgram theorem for H 1 0 ( ). There exists thus a unique element u ∈ H 1 0 ( ) such that This proves the theorem, since y is the weak solution to the Dirichlet problem (2), (3).

The control problem and the optimization theorem
For a control u ∈ L 2 ( ), the state y(u) of the system (2), (3) is given by The observation equation is given by Here, M is a given function in L ∞ ( ), n is the exterior normal to at and cos(n, x i ) is an ith direction cosine of n.
The cost functional J(v) is given by where z d is a given element in L 2 ( ) and S : L 2 ( ) → L 2 ( ) is a linear bounded self-adjoint operator (the solution operator for the problem (6), (7)) satisfying If U ad (set of admissible controls) is a closed convex subset of L 2 ( ), then making use of Lions's techniques we shall derive the necessary and sufficient conditions of minimizing the cost functional (8). The solving of the formulated optimal control problem is equivalent to seeking a u ∈ U ad such that

Theorem 3.3:
Under the above assumptions, the problem (9) admits a unique solution given by (6), (7) if, and only if, where p(u) is the adjoint state defined as the solution of the adjoint equation: Proof: It follows from Theorem 3.1 of Lions [5] that the control u ∈ U ad is optimal if and only if This condition, when explicitly calculated, gives Multiplying (10) by (y(v) − y(u)) and applying Green's formula, we obtain In view of (6) and (7), we have y(v) − y(u) − (y(v) − y(u))) = v − u and y(u) = 0 for each u ∈ . Then, the second boundary integral on the right-hand side of (13) vanishes and then it becomes: (14) Substituting (11) into (14), we get and hence (12) is equivalent to which completes the proof.

Neumann problem for the elliptic system
Let us consider the distributed parameter elliptic system with Neumann boundary condition: The classical result of Meyers and Serrin (cf., e.g. Theorem 3.17 in [1]) enables one to have the following chain of embeddings: Thus, the coerciveness condition (5) holds on H 1 ( ).
Repeating the same arguments as in Section 3.3 and using the trace theorem for the Sobolev spaces, we can prove the following existence and uniqueness theorem for the Neumann problem.

Optimization theorem for the boundary control problem
We consider the space U = L 2 ( ) (the space of controls), for every control u ∈ U, the state y(u) of the system is given as the solution of equation which may be interpreted as and the observation is given by Finally, the cost function is given by where z d is a given element in L 2 ( ) and S : L 2 ( ) → L 2 ( ) is a linear, bounded, self-adjoint operator satisfying Now, we have the following optimal control problem: where U ad (set of admissible controls) is closed convex subset of U. Repeating arguments in the proof of Theorem 3.3, we have the following theorem.

Theorem 3.5:
Under the above hypotheses, the optimal control u of (19) is characterized by (17), (18) together with where p(u) is the adjoint state.

Example 3.6:
In the case of no constraints on the control (U ad = U), we obtain Then, the optimal control is determined by simultaneously solving (17), (18), (20), (21) (where we eliminate u with the aid of (23)) and then utilizing (23).

Numerical results and comparisons
We are interested in using the JGM to solve the secondorder elliptic equation where I = (−1, 1), d = 1,2, β 1 is a nonnegative constant, and f is a given source function.
These polynomials are eigenfunctions of the following singular Sturm-Liouville equation: Let us first introduce some basic notation that will be used in the sequel. We set then the standard Jacobi-Galerkin approximation to (24) is to find u N ∈ V d N such that where ω(x) = d i=1 (1 − x i ) α (1 + x i ) β and (u, ν) w = uνω dx is the scalar product in the weighted space L 2 ω ( ). The norm in L 2 ω ( ) will be denoted by . ω . Let us denote H s ω ( ) as the weighted Sobolev spaces with the norm · s, ,ω . It is well known (see [13]) that for β 1 ≥ 0, s ≥ 1, and u ∈ H s ω ( ), the following optimal error estimate holds: (27) In practice, the proper choice of a basis for V d N is the key to achieve the optimal convergence rate (26) and to obtain a simple linear system as possible. In this section, we introduce some numerical results based on some recent works [25,26,31]. We consider the following examples.
This problem is solved in [32] using an LGM in which the integrals on the right-hand side in the resulting linear system are approximated using (N + 1)-point Legendre Gauss quadrature. The maximum pointwise errors of the proposed JGM using a tensor product process are presented in Table 2. We should note that for all values of N, the present method is always more accurate than the result of LGM [32], which shows the spectral accuracy of our method.
It should be noted that the Dirichlet boundary conditions in this example are nonhomogeneous. In [32], the nonhomogeneous Dirichlet boundary conditions at the corners are dealt with using collocation. Then, the nonhomogeneous Dirichlet boundary conditions on each side, excluding the endpoints, are treated using the L 2 projection onto the space of polynomials vanishing at the endpoints. In [31], the nonhomogeneous Dirichlet boundary conditions are handled by determining two functions u 1 and u 2 , defined on ∂ , such that u = u 1 + u 2 on ∂ . In Table 3, we present the maximum pointwise errors of u − u N using the JGM with various choices of α, β, N for the four cases of boundary conditions. Auteri et al. [18] proposed a direct spectral Galerkin-Legendre solver for the Neumann problem (Case 3) based on a double diagonalization process in which the Gauss-Legendre quadrature formula with N+1 points is used to evaluate the integrals defining the perturbed Neumann condition on the four sides and the term involving the integral of f. The corresponding integrals in the proposed method are evaluated exactly using the Jacobi polynomials.