ERROR ESTIMATES FOR THE FINITE ELEMENT DISCRETIZATION OF SEMI-INFINITE ELLIPTIC
OPTIMAL CONTROL PROBLEMS
Pedro Merino †∗ , Ira Neitzel † , Fredi Tr¨ oltzsch †
† Technische Universit¨ at Berlin Institut f¨ ur Mathematik, Germany
∗ Escuela Polit´ ecnica Nacional Departamento de Matem´ atica, Ecuador e-mails: {merino,neitzel,troeltz}@math.tu-berlin.de
Abstract
In this paper we derive a priori error estimates for linear-quadratic elliptic optimal control problems with finite dimensional control space and state constraints in the whole domain, which can be written as semi-infinite optimization problems. Numerical experiments are con- ducted to ilustrate our theory.
Keywords: elliptic optimal control problem, state constraints, error estimates, finite element discretization.
2000 Mathematics Subject Classification: 49J20, 80M10, 49N05, 41A25, 90C34.
1. Introduction
In this paper we consider elliptic optimal control problems with finite di- mensional controls and pointwise state constraints in a compact subset Ω 0 of the spatial domain Ω of the form
(P)
u∈U minad
J(y, u) = 1 2
Z
Ω
(y − y d ) 2 dx + κ 2
M
X
i=1
u 2 i
subject to y(x) ≤ b, ∀x ∈ Ω 0 ,
where y is the solution to the state equation
(1) A y(x) =
M
X
i=1
u i e i (x) in Ω, y(x) = 0 on Γ,
with a uniformly elliptic second order differential operator A and fixed func- tions e i , i = 1, . . . , M . There is a wide range of literature on a priori error analysis for elliptic optimal control problems governed by partial differential equations where the controls are given as functions. We mention for exam- ple [1, 3, 18, 10, 14, 4] or [6] for state constrained-problems. However, there are not many published results on problems with finite-dimensional control space, although they are very common for applications. In this paper, we aim at extending the optimal error estimates from [12], where a semilinear el- liptic control problem with finite dimensional control space as well as finitely many state constraints has been considered. There, error estimates of order h 2 | log(h)| for the control have been derived. For our model, the situation is more difficult, since the presence of pointwise state constraints in the domain Ω 0 rather than in finitely many points does not allow to reduce the prob- lem to a finite dimensional nonlinear programming problem. Instead of, we obtain a semi-infinite programming problem formulation. Well-established theory for semi-infinite optimization is available, we refer for example to [17, 21, 2] and the references therein for an overview, as well as to [9, 11, 20]
for numerical aspects. We also point out [19], where a discretization ap- proach is considered and a rate of convergence for the discrete solution is shown depending not only on the mesh size but also on the choice of the mesh. Yet, we are looking at additional challenges not usually found in semi- infinite programming. In contrast to the majority of papers on semi-infinite programming problems, our objective function and the constraint function are not given in explicit terms. Both are implicitly defined by the solu- tion of a PDE, such that aspects of finite-element discretization have to be considered in the numerical analysis, and the smoothness assumptions with respect to perturbations, which are standard in semi-infinite optimization, cannot be expected.
The main focus of the paper is on estimating the error in the optimal control due to a finite element discretization of the problem. We are able to prove an order of hp| log h|. Then, we improve this order to h 2 | log h|
under certain conditions, and also construct examples where this higher or-
der cannot be expected. We conclude the paper with a section on numerical
experiments.
2. Analysis of the optimal control problem
For the analysis of Problem (P ) we make the following general assumptions:
Assumption 1. Throughout the paper, let Ω ⊂ R 2 be a convex polygonal spatial domain and denote by Ω 0 ⊂ Ω a compact interior subset. The differential equation is characterized by a uniformly elliptic and symmetric differential operator A of order two. For simplicity, we choose A := −∆.
Furthermore, let κ ∈ R + be a regularization parameter, and consider bounds for the control and state, respectively, that are given by real numbers u a <
u b , and b. Moreover, y d is a given function from L 2 (Ω). For a given positive number M ∈ N consider fixed basis functions e i ∈ C 0,β (Ω), i = 1, . . . , M , for the control, with some 0 < β < 1, that are linearly independent on each open set. For convenience, we define the set of admissible controls U ad = {u ∈ R M : u a ≤ u ≤ u b }, where the inequality is to be understood component-wise. Alternatively, U ad = R M may be considered due to the presence of the regularization parameter κ > 0. By k · k, we denote the natural norm in L 2 (Ω), and (·, ·) will denote the associated inner product.
The Euclidean norm in R M will be denoted by | · |, and the inner product in R M will be denoted by h·, ·i. Last, let B r (x) denote the open ball in R 2 centered in x and with radius r.
We point out that for each basis function e i , i = 1, . . . , M , there exists a unique solution y i ∈ C 2,β (Ω) of the equation
−∆y i = e i in Ω, y i = 0 on Γ.
This follows from the regularity results in [7, Theorem 6.13] by the convexity of Ω. Moreover, by H 2 -regularity according to [8] we obtain:
Theorem 1. For each u ∈ U ad , there exists a unique solution y(u) ∈ H 2 (Ω) ∩ C 2,β (Ω) of the state equation (1) and the mapping u 7→ y is con- tinuous from R M to H 2 (Ω).
Due to the linearity of the state equation we can use the superposition principle and deduce that the solution y(u) associated with u ∈ U ad takes the form y(u)(x) = P M
i=1 u i y i (x). This allows to rewrite Problem (P ) as a
semi-infinite programming problem of the form
(P)
u∈U minad
f (u) := 1 2 k
M
X
i=1
u i y i − y d k 2 + κ 2 |u| 2 subject to
M
P
i=1
u i y i (x) ≤ b, ∀x ∈ Ω 0 .
The existence of a unique solution ¯ u ∈ U ad to (P ) with associated optimal state ¯ y follows by standard arguments, if the feasible set
U f eas := {u ∈ U ad : y(u)(x) ≤ b ∀x ∈ Ω 0 } is not empty. Next, we assume the Slater condition.
Assumption 2. There exist ˜ u ∈ U ad and ε > 0 such that (2) y(˜ u)(x) ≤ b − ε ∀x ∈ Ω 0 .
Assumption 2 guarantees the existence of a regular Borel measure as La- grange multiplier such that the first order optimality conditions can be for- mulated as a Karush-Kuhn-Tucker system. However, for the moment let us handle the state and control constraints in the set of feasible controls U f eas rather than by means of a Lagrange multiplier. The convex nature of the problem yields the standard variational inequality f 0 (¯ u)(u − ¯ u) ≥ 0 for all u ∈ U f eas , and hence we obtain:
Theorem 2. Let ¯ u be the optimal control for Problem (P ). Then, the fol- lowing inequality is fulfilled:
(3)
M
X
j=1
[κ¯ u j + (¯ y − y d , y j )] (v j − ¯ u j ) ≥ 0 ∀v ∈ U f eas .
3. Discretization and error estimates 3.1. The discrete problem formulation
In order to solve the problem numerically, we apply a standard finite-element
discretization of the problem on regular meshes with piecewise linear and
continuous ansatz functions. In this section, we are interested in the error
in the optimal control between the solution of Problem (P ) and the solution of Problem (P h ) defined below. We assume that the reader is familiar with the concept of the FEM and we do not explain the discretization in detail.
Let us denote by y i h the finite element approximates of y i , i = 1, . . . , M, and let a discretization parameter h measure the mesh size of the discretization.
Assumption 3. We assume the following accuracy of the approximation:
ky i h − y i k ≤ Ch 2 , ky i h − y i k L∞(Ω
0) ≤ Ch 2 | log h|
with a constant C > 0 not depending on h.
Remark 1. For later use, we introduce the notation α(h) := ch 2 | log h| with a generic constant c > 0, i.e., α(h) stands for O(h 2 | log h|).
The first estimate is known to hold without restrictive assumptions. We refer for example to [5] for a discussion of the discretization of an elliptic problem. The second estimate was shown in [16] under assumptions that are met in our problem setting. Note that these estimates are also valid for any linear combination y(u) = P M
i=1 u i y i and y h (u) = P M
i=1 u i y i h with u ∈ U ad . By the FE discretization, we obtain the discretized problem formulation
(P h )
u∈U minad
f h (u) := 1 2
M
X
i=1
u i y h i − y d
2 + κ 2 |u| 2 subject to P M
i=1
u i y h i (x) ≤ b, ∀x ∈ Ω 0 .
Thanks to Assumption 2, it can be shown that the feasible set of (P h ) is not empty, since the convex combination ˆ u = ¯ u + t(˜ u − ¯ u) ∈ U ad for 0 < t < 1 is feasible for (P h ) for h small enough. As a direct consequence we obtain the existence of a unique optimal solution ¯ u h of Problem (P h ) for all sufficiently small h > 0. Expressing the constraints by the feasible set
U f eas h :=
(
u ∈ U ad :
M
X
i=1
u i y h i (x) ≤ b ∀x ∈ Ω 0
)
we obtain the following first-order optimality condition:
(4)
M
X
j=1
h
κ¯ u j,h + (¯ y h − y d , y h j ) i
(v j,h − ¯ u j,h ) ≥ 0 ∀v h ∈ U f eas h .
3.2. Convergence analysis
We will prove at first the convergence result for the controls of order hp| log h| and then improve it under certain conditions. Let us mention the standard result that the Slater point ˜ u from Assumption 2 is also a Slater point for the discretized problem (P h ) with associated constant ε h = ε 2 , if h > 0 is small enough. This is not difficult to show. By means of this Slater point, we now construct an auxiliary sequence of controls feasible for (P h ), but converging to ¯ u, as well as another auxiliary sequence feasible for (P ) but converging to ¯ u h : Define u t := ¯ u + t(h)(˜ u − ¯ u) with t(h) tending to zero as h tends to zero, to be defined below. Obviously, {u t } t(h) converges to ¯ u as h tends to zero. By considering
M
X
i=1
u t,i y i h = (1 − t(h))
M
X
i=1
¯
u i y i + (1 − t(h))
M
X
i=1
¯
u i (y i h − y i ) + t(h)
M
X
i=1
˜ u i y h i
≤ (1 − t(h))b + (1 − t(h))Ch 2 | log h| + t(h)b − t(h) ε 2 ≤ b
in Ω 0 , we obtain the feasibility of u t for (P h ), where the last inequality follows by choosing t(h) = Ch Ch2 2| log h|
| log h|+ε/2 . With τ (h) = Ch Ch
22| log h|
| log h|+ε and u τ,h :=
¯
u h +τ (h)(˜ u − ¯ u h ), we obtain the existence of the second sequence converging to ¯ u h , but feasible for (P ), and an associated estimate.
Lemma 1. Let ¯ u be the optimal solution of Problem (P ) and let ¯ u h be the optimal control for (P h ). Then, with a C > 0, there holds
|¯ u − ¯ u h | ≤ Chp| log h|.
This result follows by inserting u t(h) ∈ U f eas h as a test function in the variational inequality (4) for ¯ u h , and u τ,h as a test function in the varia- tional inequality (3) for ¯ u and adding both variational inequalities. Then
|¯ u − ¯ u h | 2 ≤ Ch 2 | log h| is obtained with some C > 0 depending on the Tikhonov parameter κ > 0, and the assertion follows immediately. We omit the details, since this is a common technique used in many papers. For instance, we refer to [13, Theorem 5.1].
This error estimate is true without any other assumption than the Slater
condition (2). In order to improve the error estimate in some cases, we
impose an additional assumption on the structure of the active set of (P ).
Assumption 4. Let Ω 0 have a nonempty interior. The optimal state y ¯ is active in exactly N points ¯ x 1 , . . . , ¯ x N ∈ int Ω 0 , i.e., y(¯ ¯ x i ) = b, and all associated Lagrange multipliers are not smaller than some µ 0 > 0. We say that all x ¯ i are strongly active. Moreover, there exists σ > 0 such that (5) − hξ, ∇ 2 y(¯ ¯ x j )ξi ≥ σ|ξ| 2 ∀ξ ∈ R n , ∀j = 1, . . . , N.
Notice that ∇¯ y(¯ x i ) = 0 holds for all i = 1, . . . , N , since the ¯ x i are local maxima of ¯ y. In general, the structure of the active set can be quite di- verse as associated examples show, see [12]. For instance, ¯ y can be active on a nonempty open set or pieces of curves. However, these cases are, in some sense pathological. For instance, if ¯ y is active on a nonempty open set, then, since b is constant, −∆y vanishes. This can only happen, if ¯ u = 0 or the functions e i are linearly dependent on this set. Therefore, we will consider the situation of finitely many active points discussed in the follow- ing. Moreover, we will assume that U ad = R M . In the case of constraints, strong activity of the active constraints is usually required for convergence results. Then, however, we would readily obtain that for h small enough the associated discrete controls are also active, and hence known. Then the analysis could be restricted to the inactive constraints.
As a consequence of Assumption 4 we obtain, by Taylor expansion and the H¨older-continuity of ∇ 2 y, the existence of a real number R ¯ 1 > 0 such that
¯
y(x) ≤ ¯ y(¯ x j ) − σ
4 |x − ¯ x j | 2 = b − σ
4 |x − ¯ x j | 2 ∀x ∈ Ω 0 with |x − ¯ x j | ≤ R 1 . We know by assumption and continuity of ¯ y that for x ∈ Ω 0 \ S N
j=1 B R
1(¯ x j ) there exists δ > 0 such that ¯ y(x) ≤ b − δ. Moreover, from the convergence of ¯ u h to ¯ u we can conclude that ¯ y h converges uniformly to ¯ y in Ω 0 . Then we obtain the existence of an h 0 > 0 such that ¯ y h (x) ≤ b − δ/2 for all x ∈ Ω 0 \ S N
j=1 B R
1(¯ x j ) for all h ≤ h 0 . This implies that the discrete state can only be active in a neighborhood of the continuous active points ¯ x j , j = 1, . . . , n. Indeed, if ¯ x h j ∈ B R1(¯ x j ) is an active point of (P h ), we obtain
b = ¯ y h (¯ x h j ) = ¯ y(¯ x h j ) + pα(h) ≤ ¯y(¯x j ) − σ
4 |¯ x j − ¯ x h j | 2 + pα(h)
by Assumption 4 and Lemma 1. From ¯ y(¯ x j ) = b it follows that |¯ x j − ¯ x h j | ≤
α(h)
14. Hence, ¯ x h j ∈ B r(h) (¯ x j ) where r(h) tends to zero with order α(h)
14.
Moreover, we can show the existence of at least one associated active point ¯ x h j of Problem (P h ) in such a ball B r(h) (¯ x j ) assuming the contrary. If there were no associated discrete active point, all approximated Lagrange multipliers would vanish in all node points in B r(h) (¯ x h j ) with r(h) tending to zero with order α(h)
14and we finally would arrive at vanishing Lagrange multipliers in ¯ x j for the continuous problem (P ), which is a contradiction to the strong activity required by Assumption 4. For brevity, we leave the details to the reader.
After these considerations, let us now point out that the control ¯ u h is optimal for (P h ) if and only if it is optimal for
( ˆ P h )
u∈U minad
f h (u) := 1 2
M
X
i=1
u i y h i − y d
2
+ κ 2 |u| 2 subject to P M
i=1
u i y h i (x j ) ≤ b, ∀x j ∈ ˆ C h ,
where ˆ C h is the set of nodes of the given triangulation T h of Ω in Ω 0 . Note that ( ˆ P h ) is a completely finite-dimensional problem. To see this, we argue as follows: Let T h ∈ T h denote a triangular element. In any T h ⊂ Ω 0 , we have
¯
y h (x) ≤ b for all x ∈ T h if and only if ¯ y h (x j ) ≤ b for all corners x j of T h , since
¯
y h is linear in T h . The triangles in Ω \ Ω 0 need not be considered. All other triangles T h intersect ∂Ω 0 and we can assume T h ∩ Ω 0 ⊂ Ω 0 \ S M
j=1 B R
1(¯ x j ) for small h. Therefore, we have ¯ y h (x) ≤ b − δ/2 for all x ∈ T h ∩ Ω 0 . By continuity of ¯ y and the uniform convergence of ¯ y h towards ¯ y we find that
¯
y h (x) ≤ b − δ/4 for x ∈ T h \ Ω 0 if T h is a triangle intersecting ∂Ω. Hence, even if constraints are imposed in these triangles lying outside Ω 0 they will remain inactive if h is sufficiently small. Therefore, it is not relevant for optimality of ¯ u h whether the constraints are considered in x j ∈ ˆ C h or in all x ∈ Ω 0 . Hence, for simplicity we will denote ( ˆ P h ) by (P h ).
Lemma 2. For any j = 1, . . . , N, there exists some C > 0 and at least one grid point x ¯ h j ∈ B r(h) (¯ x j ) of Problem (P h ) where ¯ y h is active, with
|¯ x j − ¯ x h j | ≤ Chp| log h|.
P roof. We present only the key ideas of the proof. Let ¯ x j be an active
point of Problem (P ) and let ¯ x h j ∈ B r(h) (¯ x j ) with |r(h)| ≤ α(h)
14be an
associated active point for Problem (P h ), whose existence has already been argued.
• We consider first the auxiliary state ˜ y h := P M
i=1 u ¯ i,h y i and define F (x, u) := P M
i=1 u i ∇y i (x). By Assumption 4 we know that F (¯ x j , ¯ u) = 0 and
∂F
∂x (¯ x j , ¯ u) is not singular. Hence, by applying the implicit function theo- rem, we obtain the existence of ρ, τ, c > 0 such that for all u ∈ R M with
|u − ¯ u| ≤ ρ, there exists a unique ˜ x j (u) ∈ B ρ (¯ x j ) with F (˜ x j (u), u) = 0 and
|˜ x j (u) − ¯ x j | ≤ c|u − ¯ u|. Applying this to u := ¯ u h yields the existence of
˜
x h j := ˜ x j (¯ u h ) with |˜ x h j − ¯ x j | ≤ pα(h) by Lemma 1. By Assumption 4 and the H¨older continuity of ∇ 2 y(¯ ¯ x j ) we obtain coercivity of −∇ 2 y ˜ h (˜ x h j ) so that
˜
y h has a strict local maximum in ˜ x h j . Note, however, that ˜ y h may violate the constraints.
• We obtain ˜ y h (x) = ¯ y h (x) + P M
i=1 u ¯ i,h (y i (x) − y h i (x)) ≤ b + α(h) and it is clear that ˜ y h converges uniformly towards ¯ y as h tends to zero. By taking δ > 0 sufficiently small, we can therefore assume w.l.o.g. that in addition to
¯
y h (x) ≤ b − δ/2, it also holds ˜ y h (x) ≤ b − δ/2 for all x ∈ Ω 0 \ S M
i=1 B R
1(¯ x j ).
Moreover, ¯ y h (x) = ˜ y h (x) + P M
i=1 u ¯ i,h (y i h − y i ) ≤ ˜ y h + α(h) for all x ∈ Ω 0 , from which we deduce that ¯ y h can only be active where ˜ y h (x) ≥ b − α(h) holds. By the uniform estimate ˜ y h (x) ≤ b − δ/2 stated above, this can only hold inside the balls B R
1(¯ x j ). By Taylor expansion we obtain
˜
y h (x) = ˜ y h (˜ x h j ) + 1
2 hx − ˜ x h j , ∇ 2 y ˜ h (x θ j )(x − ˜ x h j )i
≤ ˜ y h (˜ x h j ) + 1
2 hx − ˜ x h j , ∇ 2 y ˜ h (x h j )(x − ˜ x h j )i + L
2 |x − ˜ x h j | β |x − ˜ x h j | 2 , with some x θ j = x + θ(˜ x h j − x), θ ∈ (0, 1), and β ∈ (0, 1) by H¨older continuity of ∇ 2 y ˜ h and ∇˜ y h (˜ x h j ) = 0. Hence, by coercivity of −∇ 2 y(˜ ˜ x h j ), we obtain the existence of R 2 > 0 not depending on h such that ˜ y h (x) ≤ ˜ y h (˜ x h j )− σ 8 |x− ˜ x h j | 2 if h is small enough and |x − ˜ x h j | ≤ R 2 . Note that by |¯ x h j − ˜ x h j | ≤ |¯ x h j − ¯ x j | +
|¯ x j − ˜ x j | ≤ α(h)
14, for h small enough we have ¯ x h j ∈ B R2(˜ x h j ). Knowing that
˜
y h (˜ x h j ) = ¯ y h (˜ x h j ) + P M
i=1 u ¯ i,h (y i − y i h ) ≤ b + α(h), we obtain from the last inequality for ˜ y h , that ˜ y h (x) ≤ b+α(h)− σ 8 |x− ˜ x h j | 2 for all x ∈ B R
2(˜ x h j ) if h is sufficiently small. Hence, collecting all estimates, we find that ¯ y h (x) can only be active if b + α(h) − σ 8 |x − ˜ x h j | 2 ≥ b − α(h), which implies |x − ˜ x h j | ≤ pα(h) for an x ∈ S N
j=1 B R
2(˜ x h j ), where ¯ y h is active, i.e., |¯ x h j − ˜ x h j | ≤ pα(h). The
assertion then follows from |¯ x h j − ¯ x j | ≤ |¯ x h j − ˜ x h j | + |˜ x h j − ¯ x j | ≤ pα(h).
Still, even the structural Assumption 4 will not necessarily guarantee a bet- ter error estimate, as an example in Section 4 will show. It is, however, possible to improve the estimate under yet an additional assumption:
Theorem 3. Let ¯ u be the optimal solution of Problem (P ), let ¯ u h be optimal for (P h ), and let the Assumptions 1–4 be satisfied. Assume further that the number of active points N is equal to the number of control variables M . Moreover, let the (M, M )-matrix Y with entries y ij = (y i (¯ x j )), i, j = 1, . . . , M , be regular. Then the following estimate is true:
|¯ u − ¯ u h | ≤ Ch 2 | log h|.
P roof. For each active point ¯ x j , j = 1, . . . , M , choose one associated dis- crete active point ¯ x h j with |¯ x j −¯ x h j | ≤ pα(h). We obtain for all j = 1, . . . , M:
M
X
i=1
¯
u i y i (¯ x j ) = b =
M
X
i=1
¯
u i,h y h i (¯ x h j ) =
M
X
i=1
¯
u i,h (y i h (¯ x h j ) − y i (¯ x h j )) + ¯ u i,h y i (¯ x h j ).
Since ¯ u h is bounded and |y i h (¯ x h j ) − y i (¯ x h j )| ≤ α(h) by Assumption 3,
α(h) ≥
M
X
i=1
¯
u i y i (¯ x j ) − ¯ u i,h y i (¯ x h j )
=
M
X
i=1
(¯ u i − ¯ u i,h )y i (¯ x j ) − ¯ u i,h ∇y i (¯ x j )(¯ x h j − ¯ x j ) + O(|¯ x h j − ¯ x j | 2 )
≥
M
X
i=1
(¯ u i − ¯ u i,h )(y i (¯ x j ) − ∇y i (¯ x j )(¯ x h j − ¯ x j )) + O(|¯ x h j − ¯ x j | 2 ) is obtained by Taylor expansion and Assumption 4, which implies
∇¯ y(¯ x j ) = 0. With Lemmas 1 and 2 we obtain
P M
i=1 (¯ u i −¯ u i,h )y i (¯ x j )
≤ α(h) for all j = 1, . . . , M. The last inequality is equivalent to |Y (¯ u − ¯ u h )| ≤ α(h).
By regularity of Y , the assertion is obtained.
4. Examples and numerical experiments
We show by means of a simple example that the result of Theorem 3, which
was proven under quite strong assumptions, cannot generally be expected
under Assumption 4 only. Therefore, consider the optimization problem without relation to PDEs,
(P 1 )
u
1min ,u2∈R J(u) := (u 1 − 1) 2 + 1
4 u 2 2 − u 2 + 1 4 u 1 subject to − u 1 x 2 + u 2 x ≤ 1 4 ∀x ∈ (−1, 1).
The reader may verify that the unique optimal solution to this problem is given by ¯ u 1 = ¯ u 2 = 1, with exactly one active point ¯ x = 1 2 . Prescribing the constraints only in the grid points x i , i = 0, . . . , n, of a discretiza- tion of (−1, 1) yields an associated discrete problem formulation. Note that this can be interpreted as approximation of the functions y 1 = x and y 2 = −x 2 by their piecewise linear nodal interpolants y h 1 , y 2 h . We choose a grid with discretization parameter h = 2 n and grid points x i = −1 + (i + 1 3 )h, i = 1, . . . , n − 1, as well as x 0 = −1 and x n = 1. Note the special struc- ture of the grid, where the active point of the continuous solution, ¯ x = 1 2 , has the distance 2 3 h to its nearest neighboring grid point to the left, and distance 1 3 h to its nearest neighboring grid point to the right for each h.
For a given h > 0, the unique optimal solution to the discretized problem is given by
¯ u h 1 =
17 4 + 2 3 h
4 + 1 2 + 1 3 h 2 , u ¯ h 2 = 2 − 2
1 2 + 1 3 h
2 ¯ u h 1 − 1 + 1 4
with one active grid point ¯ x h = 1 2 + 1 3 h. Obviously, (¯ u h 1 , ¯ u h 2 ) converges to (¯ u 1 , ¯ u 2 ) with order h only, even though ky h i − y i k C( ¯ Ω) ≤ ch 2 , i = 1, 2.
Note here that the p| log h|-term in Lemma 1 originates in the FEM-error
of the state equation. Since here we use the piecewise linear interpolants
instead of finite element approximations, the logarithmic term does not ap-
pear in the error estimate. For completeness, we have solved the problem in
Matlab, using the solution routine quadprog, and show in Figure 1 the
experimental error in the control in logarithmic scale. Clearly, linear con-
vergence is observed. This example suggests that also in the case of PDEs
the convergence result of Lemma 1 cannot generally be improved except for
special cases such as discussed in Theorem 3. To illustrate this, we con-
sider the following one-dimensional elliptic PDE example motivated by the
problem above.
(P 2 )
u
1min ,u2∈R
1
2 ku 1 y 1 + u 2 y 2 − y d k 2 + (u 1 − 1) 2 + 1 4 u 1 + 1
4 u 2 2 − 235 228 u 2 subject to the constraints:
−∆y 1 (x) + y 1 (x) = u 1 (2 − x 2 ) −∆y 2 (x) + y 2 (x) = u 2 x
y 1 (0) = 0 y 2 (0) = 0
y 1 (1) = −1 y 2 (1) = 1
y(x) ≤ 1 4 , ∀x ∈ (0, 1).
The example is constructed such that y 1 = −x 2 and y 2 = x are solutions of the PDEs and the optimal solution of Problem (P 2 ) is again given by
¯
u 1 = ¯ u 2 = 1 with one active point at ¯ x = 1 2 . We point out that, strictly speaking, this example does not fit into our theoretical setting, since y 1 and y 2 do not admit homogeneous Dirichlet boundary conditions on the given boundary. Nevertheless, we compute linear finite element approximations y 1 h and y h 2 of y 1 and y 2 , respectively. We choose the same grid as in example (P 1 h ) and solve the associated discrete problem with Matlab’s optimization routine quadprog. Figure 2 shows the error in the control in logarithmic scale, which indicates linear convergence, only. Note that we do not expect to see the influence of the | log |-term numerically.
10−4 10−3 10−2 10−1 100
10−5 10−4 10−3 10−2 10−1 100
h error in uh
Computed error linear order of convergence
Figure 1. Example (P 1 ).
10−4 10−3 10−2 10−1 100
10−5 10−4 10−3 10−2 10−1 100
h error in uh
Computed error Linear order of convergence
Figure 2. Example (P 2 ).
Let us also show an example where the higher order of convergence from
Theorem 3 is to be expected. The following example is the semi-infinite
linear-quadratic version of Example 1 in [12], which was motivated by [15],
with a linear partial differential equation in two dimensions and one single
strongly active point located at the origin.
(P 3 )
min u∈R
1
2 ky(u) − y d k 2 + 1
2 |u − u d | 2 subject to the constraints:
−∆y(x) + y(x) = u(1 − 1/5(x 2 1 + x 2 2 )) in Ω = B(0, 1)
y(x) = 0. on Γ = S(0, 1)
y(x) ≤ 1, ∀x ∈ B 0.9 (0).
The desired control and state are chosen to be u d = 5 + 19/80 and y d = 1−(x 2 1 +x 2 2 )+ 2π 1 log(|x|), respectively. We compare the numerical solution of this problem computed using MATLAB’s optimization routine quadprog with the known solution ¯ u = 5 with associated optimal state ¯ y = 1 − (x 2 1 + x 2 2 ). We choose an initial grid such that the continuous active point
¯
x = 0 is neither a grid point nor exactly in the center of the containing triangle. We observe that the distance of the active points, |¯ x − ¯ x h |, is not decreasing uniformly, which seems to influence the convergence process.
Nevertheless, the numerical results indicate quadratic convergence, as can be seen in Figure 3, where we show the computed error in the control in logarithmic scale, compared to a quadratic error bound, and also include the distance |¯ x − ¯ x h |. To explain the nonuniform decrease in the error, we consider a simple, non-PDE-related example given by
(P 4 )
min u∈R J(u) := (u − 1) 2 − 1 4 u
subject to u(−x 2 + x) ≤ 1 4 ∀x ∈ (−1, 1).
The unique optimal solution to this problem is given by ¯ u = 1, admitting exactly one active point ¯ x = 1 2 , such that the number of controls equals the number of active points. From the convergence result of Theorem 3 we expect an order of h 2 for the error in the control variable, without the
| log |-term if the piecewise linear interpolants of the state are used instead
of a finite-element discretization. We use this example to show how the
grid influences not the order but the constants of the error estimate. For
the numerical approximation, we choose two different grids with principal
mesh-size h, x i = (i + 1 3 )h, and x i = (i + 1 5 )h for i = 1, . . . , n − 1, with
h = n 2 . The optimal solution on the first grid is given by ¯ u h/3 = 9−4h 9 2,
the optimal solution on the second grid is ¯ u h/5 = 25−4h 25 2, which obviously
, which obviously
converge to ¯ u = 1 with order h 2 . In Figure 4 we show the experimental order of convergence for the control alternating between both grids as h decreases.
For comparison, we show lines indicating quadratic order of convergence with two computed constants associated with the two grids. It becomes clear that the error shows quadratic convergence behavior with grid-dependent constants. Here, the constants depend on the distance of the active points,
|¯ x − ¯ x h |. We expect to see this behavior in PDE examples, but point out that due to the FEM-discretization this effect is blurred by an additional error.
10−2 10−1 100
10−5 10−4 10−3 10−2 10−1 100
h
Quadratic order of convergence Distance of active points |x− xh| Computed error |u−uh|
Figure 3. Example (P 3 ).
10−4 10−3 10−2 10−1 100
10−7 10−6 10−5 10−4 10−3 10−2 10−1 100
h error in uh
Computed error
quadratic order of convergence for h/3 grid quadratic order of convergence for h/5 grid
