Optimal Control of the Laplace-Beltrami operator on compact surfaces - concept and numerical treatment

3 MB·English

by  Michael Hinze

Hamburger Beitr¨age
zur Angewandten Mathematik
1
1
0
2

n
a Optimal Control of the Laplace-Beltrami operator on compact
J
surfaces – concept and numerical treatment
7

]
C
Michael Hinze and Morten Vierling
O
.
h
t
a
m
[

1
v
5
8
3
1
.
1
0
1
1
:
v
i
X
r
a
Nr. 2011-1
January 2011

Abstract: We consider optimal control problems of elliptic PDEs on hypersurfaces Γ in Rn for
n=2,3.TheleadingpartofthePDEisgivenbytheLaplace-Beltramioperator,whichisdiscretizedby
finiteelementsonapolyhedralapproximationofΓ.Thediscreteoptimalcontrolproblemisformulated
on the approximating surface and is solved numerically with a semi-smooth Newton algorithm. We
deriveoptimalapriorierrorestimatesforproblemsincludingcontrolconstraintsandprovidenumerical
examples confirming our analytical findings.
Mathematics Subject Classification (2010): 58J32 , 49J20, 49M15
Keywords:Ellipticoptimalcontrolproblem,Laplace-Beltramioperator,surfaces,controlconstraints,
error estimates,semi-smooth Newton method.
1 Introduction
Weareinterestedinthenumericaltreatmentofthefollowinglinear-quadraticoptimalcontrol
problem on a n-dimensional, sufficiently smooth hypersurface Γ ⊂ Rn+1, n = 1,2.
1 α
min J(u,y) = (cid:107)y−z(cid:107)2 + (cid:107)u(cid:107)2
u∈L2(Γ),y∈H1(Γ) 2 L2(Γ) 2 L2(Γ)
subject to u ∈ U and (1.1)
ad
(cid:90) (cid:90)
∇ y∇ ϕ+cyϕdΓ = uϕdΓ,∀ϕ ∈ H1(Γ)
Γ Γ
Γ Γ
with U = (cid:8)v ∈ L2(Γ) | a ≤ v ≤ b(cid:9), a < b ∈ R . For simplicity we will assume Γ to be
ad
compact and c = 1. In section 4 we briefly investigate the case c = 0, in section 5 we give an
example on a surface with boundary.
Problem (1.1) may serve as a mathematical model for the optimal distribution of surfactants
on a biomembrane Γ with regard to achieving a prescribed desired concentration z of a
quantity y.
It follows by standard arguments that (1.1) admits a unique solution u ∈ U with unique
ad
associated state y = y(u) ∈ H2(Γ).
Our numerical approach uses variational discretization applied to (1.1), see [Hin05] and
[HPUU09], on a discrete surface Γh approximating Γ. The discretization of the state equation
in (1.1) is achieved by the finite element method proposed in [Dzi88], where a priori error
estimates for finite element approximations of the Poisson problem for the Laplace-Beltrami
operator are provided. Let us mention that uniform estimates are presented in [Dem09], and
stepstowardsaposteriorierrorcontrolforellipticPDEsonsurfacesaretakenbyDemlowand
Dziuk in [DD07]. For alternative approaches for the discretization of the state equation by
finite elements see the work of Burger [Bur08]. Finite element methods on moving surfaces
are developed by Dziuk and Elliott in [DE07]. To the best of the authors knowledge, the
present paper contains the first attempt to treat optimal control problems on surfaces.
We assume that Γ is of class C2 with unit normal field ν. As an embedded, compact hy-
persurface in Rn+1 it is orientable and hence the zero level set of a signed distance function
|d(x)| = dist(x,Γ). We assume w.l.o.g. ∇d(x) = ν(x) for x ∈ Γ. Further, there exists an
neighborhood N ⊂ Rn+1 of Γ, such that d is also of class C2 on N and the projection
1
a : N → Γ, a(x) = x−d(x)∇d(x) (1.2)
is unique, see e.g. [GT98, Lemma 14.16]. Note that ∇d(x) = ν(a(x)).
Using a we can extend any function φ : Γ → R to N as φ¯(x) = φ(a(x)). This allows us to
represent the surface gradient in global exterior coordinates ∇ φ = (I −ννT)∇φ¯, with the
Γ
euclidean projection (I −ννT) onto the tangential space of Γ.
We use the Laplace-Beltrami operator ∆ = ∇ ·∇ in its weak form i.e. ∆ : H1(Γ) →
Γ Γ Γ Γ
H1(Γ)∗
(cid:90)
y (cid:55)→ − ∇ y∇ (·)dΓ ∈ H1(Γ)∗.
Γ Γ
Γ
Let S denote the prolongated restricted solution operator of the state equation
S : L2(Γ) → L2(Γ), u (cid:55)→ y −∆ y+cy = u,
Γ
which is compact and constitutes a linear homeomorphism onto H2(Γ), see [Dzi88, 1. Theo-
rem].
By standard arguments we get the following necessary (and here also sufficient) conditions
for optimality of u ∈ U
ad
(cid:104)∇ J(u,y(u)),v−u(cid:105) = (cid:104)αu+S∗(Su−z),v−u(cid:105) ≥ 0 ∀v ∈ U , (1.3)
u L2(Γ) L2(Γ) ad
We rewrite (1.3) as
(cid:18) (cid:19)
1
u = P − S∗(Su−z) , (1.4)
U
ad α
where P denotes the L2-orthogonal projection onto U .
Uad ad
2 Discretization
We now discretize (1.1) using an approximation Γh to Γ which is globally of class C0,1.
FollowingDziuk,weconsiderpolyhedralΓh = (cid:83) Ti consistingoftrianglesTi withcorners
i∈Ih h h
on Γ, whose maximum diameter is denoted by h. With FEM error bounds in mind we assume
the family of triangulations Γh to be regular in the usual sense that the angles of all triangles
are bounded away from zero uniformly in h.
We assume for Γh that a(Γh) = Γ, with a from (1.2). For small h > 0 the projection a also
is injective on Γh. In order to compare functions defined on Γh with functions on Γ we use a
to lift a function y ∈ L2(Γh) to Γ
yl(a(x)) = y(x) ∀x ∈ Γh,
and for y ∈ L2(Γ) and sufficiently small h > 0 we define the inverse lift
y (x) = y(a(x)) ∀x ∈ Γh.
l
For small mesh parameters h the lift operation (·) : L2(Γ) → L2(Γh) defines a linear homeo-
l
morphism with inverse (·)l. Moreover, there exists c > 0 such that
int
1−c h2 ≤ (cid:107)(·) (cid:107)2 ,(cid:107)(·)l(cid:107)2 ≤ 1+c h2, (2.1)
int l L(L2(Γ),L2(Γh)) L(L2(Γh),L2(Γ)) int
as the following lemma shows.
2
Lemma and Definition 2.1. Denote by dΓ the Jacobian of a| : Γh → Γ, i.e. dΓ =
dΓh Γh dΓh
|det(M)| where M ∈ Rn×n represents the Derivative da(x) : T Γh → T Γ with respect to
x a(x)
arbitrary orthonormal bases of the respective tangential space. For small h > 0 there holds
(cid:12) (cid:12)
sup(cid:12)(cid:12)(cid:12)1− ddΓΓh(cid:12)(cid:12)(cid:12) ≤ cinth2,
Γ
Now let dΓh denote |det(M−1)|, so that by the change of variable formula
dΓ
(cid:12)(cid:12)(cid:12)(cid:12)(cid:90) vldΓh−(cid:90) vdΓ(cid:12)(cid:12)(cid:12)(cid:12) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) vddΓΓh −vdΓ(cid:12)(cid:12)(cid:12)(cid:12) ≤ cinth2(cid:107)v(cid:107)L1(Γ).
Γh Γ Γ
Proof. see [DE07, Lemma 5.1]
Problem (1.1) is approximated by the following sequence of optimal control problems
1 α
min J(u,y) = (cid:107)y−z (cid:107)2 + (cid:107)u(cid:107)2
u∈L2(Γh),y∈H1(Γh) 2 l L2(Γh) 2 L2(Γh)
(2.2)
subject to u ∈ Uh and
ad
y = S u,
h
with Uh = (cid:8)v ∈ L2(Γh) | a ≤ v ≤ b(cid:9), i.e. the mesh parameter h enters into U only through
ad ad
Γh . Problem (2.2) may be regarded as the extension of variational discretization introduced
in [Hin05] to optimal control problems on surfaces.
In [Dzi88] it is explained, how to implement a discrete solution operator S : L2(Γh) →
h
L2(Γh), such that
(cid:107)(·)lS (·) −S(cid:107) ≤ C h2, (2.3)
h l L(L2(Γ),L2(Γ)) FE
which we will use throughout this paper. See in partikular [Dzi88, Equation (6)] and [Dzi88,
7. Lemma]. For the convenience of the reader we briefly sketch the method. Consider the
space
(cid:110) (cid:16) (cid:17) (cid:12) (cid:111)
V = ϕ ∈ C0 Γh (cid:12) ∀i ∈ I : ϕ| ∈ P1(Ti) ⊂ H1(Γh)
h (cid:12) h Ti h
h
of piecewise linear, globally continuous functions on Γh. For some u ∈ L2(Γ), to compute
yl = (·)lS (·) u solve
h h l
(cid:90) (cid:90)
∇ y ∇ ϕ +cy ϕ dΓh = u ϕ dΓh, ∀ϕ ∈ V
Γh h Γh i h i l i h
Γh Γh
fory ∈ V .WechooseL2(Γh)ascontrolspace,becauseingeneralwecannotevaluate(cid:82) vdΓ
exacthly, whhereas the expression (cid:82) v dΓh for piecewise polynomials v can be computΓed up
Γh l l
to machine accuracy. Also, the operator S is self-adjoint, while ((·)lS (·) )∗ = (·) ∗S (·)l∗ is
h h l l h
not. The adjoint operators of (·) and (·)l have the shapes
l
dΓh dΓ
∀v ∈ L2(Γh) : ((·) )∗v = vl, ∀v ∈ L2(Γ) : ((·)l)∗v = v , (2.4)
l dΓ dΓh l
hence evaluating (·) ∗ and (·)l∗ requires knowledge of the Jacobians dΓh and dΓ which may
l dΓ dΓh
not be known analytically.
3
Similar to (1.1), problem (2.2) possesses a unique solution u ∈ Uh which satisfies
h ad
(cid:18) (cid:19)
1
u = P − p (u ) . (2.5)
h Uahd α h h
Here P : L2(Γh) → Uh is the L2(Γh)-orthogonal projection onto Uh and for v ∈ L2(Γh)
Uh ad ad
ad
the adjoint state is p (v) = S∗(S v−z ) ∈ H1(Γh).
h h h l
Observe that the projections P and P coincide with the point-wise projection P on
Uad Uahd [a,b]
Γ and Γh, respectively, and hence
(cid:16) (cid:17)l
P (v ) = P (v) (2.6)
Uahd l Uad
for any v ∈ L2(Γ).
Let us now investigate the relation between the optimal control problems (1.1) and (2.2).
Theorem 2.2 (Order of Convergence). Let u ∈ L2(Γ), u ∈ L2(Γh) be the solutions of (1.1)
h
and (2.2), respectively. Then for sufficiently small h > 0 there holds
α(cid:13)(cid:13)ulh−u(cid:13)(cid:13)2L2(Γ)+(cid:13)(cid:13)yhl −y(cid:13)(cid:13)2L2(Γ) ≤ 11+−cciinntthh22(cid:18)α1 (cid:13)(cid:13)(cid:13)(cid:16)(·)lSh∗(·)l −S∗(cid:17)(y−z)(cid:13)(cid:13)(cid:13)2L2(Γ)...
(2.7)
(cid:13)(cid:16) (cid:17) (cid:13)2 (cid:19)
+(cid:13) (·)lS (·) −S u(cid:13) ,
(cid:13) h l (cid:13)
L2(Γ)
with y = Su and y = S u .
h h h
Proof. From (2.6) it follows that the projection of −(cid:0)1p(u)(cid:1) onto Uh is u
α l ad l
(cid:18) (cid:19)
1
u = P − p(u) ,
l Uahd α l
which we insert into the necessary condition of (2.2). This gives
(cid:104)αu +p (u ),u −u (cid:105) ≥ 0.
h h h l h L2(Γh)
On the other hand u is the L2(Γh)-orthogonal projection of −1p(u) , thus
l α l
1
(cid:104)− p(u) −u ,u −u (cid:105) ≤ 0.
α l l h l L2(Γh)
Adding these inequalities yields
α(cid:107)u −u (cid:107)2 ≤(cid:104)(p (u )−p(u) ),u −u (cid:105)
l h L2(Γh) h h l l h L2(Γh)
=(cid:104)p (u )−S∗(y−z) ,u −u (cid:105) +(cid:104)S∗(y−z) −p(u) ,u −u (cid:105) .
h h h l l h L2(Γh) h l l l h L2(Γh)
The first addend is estimated via
(cid:104)p (u )−S∗(y−z) ,u −u (cid:105) = (cid:104)y −y ,S u −y (cid:105)
h h h l l h L2(Γh) h l h l h L2(Γh)
= −(cid:107)y −y (cid:107)2 +(cid:104)y −y ,S u −y (cid:105)
h l L2(Γh) h l h l l L2(Γh)
1 1
≤ − (cid:107)y −y (cid:107)2 + (cid:107)S u −y (cid:107)2 .
2 h l L2(Γh) 2 h l l L2(Γh)
4
The second addend satisfies
α 1
(cid:104)S∗(y−z) −p(u) ,u −u (cid:105) ≤ (cid:107)u −u (cid:107)2 + (cid:107)S∗(y−z) −p(u) (cid:107)2 .
h l l l h L2(Γh) 2 l h L2(Γh) 2α h l l L2(Γh)
Together this yields
1
α(cid:107)u −u (cid:107)2 +(cid:107)y −y (cid:107)2 ≤ (cid:107)S∗(y−z) −p(u) (cid:107)2 +(cid:107)S u −y (cid:107)2
l h L2(Γh) h l L2(Γh) α h l l L2(Γh) h l l L2(Γh)
The claim follows using (2.1) for sufficiently small h > 0.
Because both S and S are self-adjoint, quadratic convergence follows directly from (2.7).
h
For operators that are not self-adjoint one can use
(cid:107)(·) ∗S∗(·)l∗−S∗(cid:107) ≤ C h2. (2.8)
l h L(L2(Γ),L2(Γ)) FE
which is a consequence of (2.3). Equation (2.4) and Lemma 2.1 imply
(cid:107)((·) )∗−(·)l(cid:107) ≤ c h2, (cid:107)((·)l)∗−(·) (cid:107) ≤ c h2. (2.9)
l L(L2(Γh),L2(Γ)) int l L(L2(Γ),L2(Γh)) int
Combine (2.7) with (2.8) and (2.9) to proof quadratic convergence for arbitrary linear elliptic
state equations.
3 Implementation
Inordertosolve(2.5)numerically,weproceedasin[Hin05]usingthefiniteelementtechniques
forPDEsonsurfacesdevelopedin[Dzi88]combinedwiththesemi-smoothNewtontechniques
from [HIK03] and [Ulb03] applied to the equation
(cid:18) (cid:18) (cid:19) (cid:19)
1
G (u ) = u −P − p (u ) . = 0 (3.1)
h h h [a,b] α h h
Since the operator p continuously maps v ∈ L2(Γh) into H1(Γh), Equation (3.1) is semis-
h
mooth and thus is amenable to a semismooth Newton method. The generalized derivative of
G is given by
h
(cid:16) χ (cid:17)
DG (u) = I + S∗S ,
h α h h
where χ : Γh → {0,1} denotes the indicator function of the inactive set I(−1p (u)) =
(cid:8)γ ∈ Γh (cid:12)(cid:12) a < −α1ph(u)[γ] < b(cid:9) α h

 1 on I(−1p (u)) ⊂ Γh
χ = α h ,
 0 elsewhere on Γh
which we use both as a function and as the operator χ : L2(Γh) → L2(Γh) defined as the
point-wise multiplication with the function χ. A step semi-smooth Newton method for (3.1)
then reads
(cid:18) (cid:19)
(cid:16) χ (cid:17) 1 χ
I + S∗S u+ = −G (u)+DG (u)u = P − p (u) + S∗S u.
α h h h h [a,b] α h α h h
Given u the next iterate u+ is computed by performing three steps
5
1. Set ((1−χ)u+)[γ] = (cid:0)(1−χ)P (cid:0)−1p (u)+m(cid:1)(cid:1)[γ], which is either a or b, depend-
[a,b] α h
ing on γ ∈ Γ .
h
2. Solve
(cid:16) χ (cid:17) χ(cid:16) (cid:17)
I + S∗S χu+ = S∗z −S∗S (1−χ)u+
α h h α h l h h
for χu+ by CG iteration over L2(I(−1p (u)).
α h
3. Set u+ = χu++(1−χ)u+.
Details can be found in [HV11] .
4 The case c = 0
Inthissectionweinvestigatethecasec = 0whichcorrespondstoastationary,purelydiffusion
driven process. Since Γ has no boundary, in this case total mass must be conserved, i.e. the
state equation admits a solution only for controls with mean value zero. For such a control
the state is uniquely determined up to a constant. Thus the admissible set U has to be
ad
changed to
(cid:26) (cid:12) (cid:90) (cid:27)
Uad = (cid:8)v ∈ L2(Γ) | a ≤ v ≤ b(cid:9)∩L20(Γ), where L20(Γ) := v ∈ L2(Γ) (cid:12)(cid:12) vdΓ = 0 ,
(cid:12)
Γ
(cid:82)
and a < 0 < b. Problem (1.1) then admits a unique solution (u,y) and there holds ydΓ =
(cid:82) (cid:82) Γ
zdΓ. W.l.o.g we assume zdΓ = 0 and therefore only need to consider states with mean
Γ Γ
value zero. The state equation now reads y = S˜u with the solution operator S˜ : L2(Γ) →
0
L2(Γ) of the equation −∆ y = u, (cid:82) ydΓ = 0.
0 Γ Γ
Using the injection L2(Γ) →ı L2(Γ), S˜ is prolongated as an operator S : L2(Γ) → L2(Γ) by
0
S = ıS˜ı∗. The adjoint ı∗ : L2(Γ) → L2(Γ) of ı is the L2-orthogonal projection onto L2(Γ).
0 0
The unique solution of (1.1) is again characterized by (1.4), where the orthogonal projection
now takes the form
P (v) = P (v+m)
Uad [a,b]
with m ∈ R chosen such that
(cid:90)
P (v+m) dΓ = 0.
[a,b]
Γ
If for v ∈ L2(Γ) the inactive set I(v +m) = {γ ∈ Γ | a < v[γ]+m < b} is non-empty, the
constantm=m(v)isuniquelydeterminedbyv ∈ L2(Γ).Hence,thesolutionu ∈ U satisfies
ad
(cid:18) (cid:18) (cid:19)(cid:19)
1 1
u = P − p(u)+m − p(u) ,
[a,b] α α
withp(u) = S∗(Su−ı∗z) ∈ H2(Γ)denotingtheadjointstateandm(−1p(u)) ∈ Risimplicitly
α
given by (cid:82) udΓ = 0. Note that ı∗ı is the identity on L2(Γ).
In (2.2) weΓnow replace Uh by Uh = (cid:8)v ∈ L2(Γh) | a ≤0v ≤ b(cid:9)∩L2(Γh). Similar as in (2.5),
ad ad 0
the unique solution u then satisfies
h
(cid:18) (cid:19) (cid:18) (cid:18) (cid:19)(cid:19)
1 1 1
u = P − p (u ) = P − p (u )+m − p (u ) , (4.1)
h Uahd α h h [a,b] α h h h α h h
6
with p (v ) = S∗(S v −ı∗z ) ∈ H1(Γh) and m (−1p (u )) ∈ R the unique constant such
h h h h h h l h α h h
that (cid:82) u dΓh = 0. Note that m (cid:0)−1p (u )(cid:1) is semi-smooth with respect to u and thus
Γh h h α h h h
Equation (4.1) is amenable to a semi-smooth Newton method.
The discretization error between the problems (2.2) and (1.1) now decomposes into two
components, one introduced by the discretization of U through the discretization of the
ad
surface, the other by discretization of S.
For the first error we need to investigate the relation between P (u) and P (u), which
Uahd Uad
is now slightly more involved than in (2.6).
Lemma 4.1. Let h > 0 be sufficiently small. There exists a constant C > 0 depending only
m
on Γ, |a| and |b| such that for all v ∈ L2(Γ) with (cid:82) dΓ > 0 there holds
I(v+m(v))
C
|m (v )−m(v)| ≤ m h2.
h l (cid:82)
dΓ
I(v+m(v))
Proof. For v ∈ L2(Γ), (cid:15) > 0 choose δ > 0 and h > 0 so small that the set
(cid:110) (cid:111)
Iδ = γ ∈ Γh | a+δ ≤ v (γ)+m(v) ≤ b−δ .
v l
satisfies (cid:82) dΓh(1+(cid:15)) ≥ (cid:82) dΓ. It is easy to show that hence m (v ) is unique. Set
Iδ I(v+m(v)) h l
v (cid:82)
C = c max(|a|,|b|) dΓ. Decreasing h further if necessary ensures
int Γ
Ch2 Ch2
≤ (1+(cid:15)) ≤ δ.
(cid:82) (cid:82)
dΓh dΓ
Iδ I(v+m(v))
v
For x ∈ R let
(cid:90)
Mh(x) = P (v +x) dΓh.
v [a,b] l
Γh
(cid:82)
Since P (v+m(v)) dΓ = 0, Lemma 2.1 yields
Γ [a,b]
|Mh(m(v))| ≤ c (cid:107)P (v+m(v))(cid:107) h2 ≤ Ch2.
v int [a,b] L1(Γ)
Let us assume w.l.o.g. −Ch2 ≤ Mh(m(v)) ≤ 0. Then
v
(cid:32) (cid:33)
Ch2
Mh m(v)+ ≥ Mh(m(v))+Ch2 ≥ 0
v (cid:82) dΓh v
Iδ
v
implies 0 ≤ m(v)−m (v) ≤ Ch2 ≤ (1+(cid:15))C h2, since Mh(x) is continuous with respect
h (cid:82) dΓh (cid:82) dΓ v
Ivδ I(v+m(v))
to x. This proves the claim.
Because
(cid:16) (cid:17)l
P (v ) −P (v) = P (v+m (v ))−P (v+m(v)) ,
Uahd l Uad [a,b] h l [a,b]
we get the following corollary.
Corollary 4.2. Let h > 0 be sufficiently small and C as in Lemma 4.1. For any fixed
m
v ∈ L2(Γ) with (cid:82) dΓ > 0 we have
I(v+m(v))
(cid:113)
(cid:82)
(cid:13) (cid:13) dΓ
(cid:13)(cid:13)(cid:13)(cid:16)PUahd(vl)(cid:17)l −PUad(v)(cid:13)(cid:13)(cid:13)L2(Γ) ≤ Cm(cid:82)I(v+mΓ(v)) dΓh2.
7
Note that since for u ∈ L2(Γ) the adjoint p(u) is a continuous function on Γ, the corollary is
applicable for v = −1p(u).
α
The following theorem can be proofed along the lines of Theorem 2.2.
Theorem 4.3. Let u ∈ L2(Γ), u ∈ L2(Γh) be the solutions of (1.1) and (2.2), respectively,
h
in the case c = 0. Let u˜ = (cid:16)P (cid:0)−1p(u) (cid:1)(cid:17)l. Then there holds for (cid:15) > 0 and 0 ≤ h < h
h Uh α l (cid:15)
ad
α(cid:107)ulh−u˜h(cid:107)2L2(Γ)+(cid:13)(cid:13)yhl −y(cid:13)(cid:13)2L2(Γ) ≤ (1+(cid:15))(cid:18)α1 (cid:13)(cid:13)(cid:13)(cid:16)(·)lSh∗(·)l −S∗(cid:17)(y−z)(cid:13)(cid:13)(cid:13)2L2(Γ)...
(cid:13) (cid:13)2 (cid:19)
+(cid:13)(·)lS (·) u˜ −y(cid:13) .
(cid:13) h l h (cid:13)
L2(Γ)
Using Corollary 4.2 we conclude from the theorem
(cid:107)ulh−u(cid:107)L2(Γ) ≤C(cid:18)α1 (cid:13)(cid:13)(cid:13)(cid:13)(cid:18)(·)lSh∗(·)l −S∗(cid:19)(y−z)(cid:13)(cid:13)(cid:13)(cid:13)L2(Γ)+ √1α (cid:13)(cid:13)(cid:13)(cid:16)(·)lSh(·)l −S(cid:17)u(cid:13)(cid:13)(cid:13)L2(Γ)...
(cid:113)
(cid:18) (cid:107)S(cid:107)L(L2(Γ),L2(Γ))(cid:19) Cm (cid:82)Γ dΓh2 (cid:19)
+ 1+ √ ,
(cid:82)
α dΓ
I(−1p(u)+m(−1p(u)))
α α
the latter part of which is the error introduced by the discretization of U . Hence one has
ad
h2-convergence of the optimal controls.
5 Numerical Examples
ThefiguresshowsomeselectedNewtonstepsu+.Notethatjumpsofthecolor-codedfunction
values are well observable along the border between active and inactive set. For all examples
Newton’s method is initialized with u ≡ 0.
0
The meshes are generated from a macro triangulation through congruent refinement, new
nodes are projected onto the surface Γ. The maximal edge length h in the triangulation is
not exactly halved in each refinement, but up to an error of order O(h2). Therefore we just
compute our estimated order of convergence (EOC) according to
ln(cid:107)u −u (cid:107) −ln(cid:107)u −u (cid:107)
hi−1 l L2(Γhi−1) hi l L2(Γhi)
EOC = .
i
ln(2)
For different refinement levels, the tables show L2-errors, the corresponding EOC and the
number of Newton iterations before the desired accuracy of 10−6 is reached.
It was shown in [HU04], under certain assumptions on the behaviour of −1p(u), that the
α
undamped Newton Iteration is mesh-intdependent. These assumptions are met by all our
examples, since the surface gradient of −1p(u) is bounded away from zero along the border
α
of the inactive set. Moreover, the displayed number of Newton-Iterations suggests mesh-
independence of the semi-smooth Newton method.
Example 5.1 (Sphere I). We consider the problem
min J(u,y) subject to −∆ y+y = u−r, −1 ≤ u ≤ 1 (5.1)
Γ
u∈L2(Γ),y∈H1(Γ)
8