CC-BY-NC-SA 4.0 arxiv 2211.12246v4

###### Acknowledgements.

This research was supported by the German Research Foundation (DFG) under
grant numbers WA 3626/3-2 within the priority program “Non-smooth and
Complementarity-based Distributed Parameter Systems: Simulation and
Hierarchical Optimization” (SPP 1962).

# A topological derivative-based algorithm to solve optimal control problems
with $L^{0}(\Omega)$ control cost

Daniel Wachsmuth Institut für Mathematik, Universität Würzburg, 97074
Würzburg, Germany, , 0000-0001-7828-5614<EMAIL_ADDRESS>wuerzburg.de

(2023-10-04)

###### Abstract

In this paper, we consider optimization problems with $L^{0}$-cost of the
controls. Here, we take the support of the control as independent optimization
variable. Topological derivatives of the corresponding value function with
respect to variations of the support are derived. These topological
derivatives are used in a novel gradient descent algorithm with Armijo line-
search. Under suitable assumptions, the algorithm produces a minimizing
sequence.

###### keywords:

Topological derivative, control support optimization, sparse optimal control,
$L^{0}$ optimization.

> MSC (2020) 49M05, 49K40, 65K10

##  1 Introduction

In this paper we are interested in the following optimal control problem:
Minimize

(1.1)
$\min\frac{1}{2}\|y-y_{d}\|_{L^{2}(\Omega)}^{2}+\frac{\alpha}{2}\|u\|_{L^{2}(\Omega)}^{2}+\beta\|u\|_{0}$

over all $(y,u)$ satisfying

(1.2) $\displaystyle-\Delta y$ $\displaystyle=u\quad\text{ on }\Omega$
$\displaystyle y$ $\displaystyle=0\quad\text{ on }\partial\Omega$

and

(1.3) $u_{a}\leq u\leq u_{b}.$

Here, $\|u\|_{0}$ is the measure of the support of $u$. This optimal control
problem can be interpreted in the context of optimal actuator placement: Find
a (possibly small) measurable set $A\subseteq\Omega$ such that controls
supported on $A$ can still minimize a certain objective functional. We remark
that the elliptic equation in (1.2) can be replaced by other types of partial
differential equations, for example parabolic or hyperbolic equations. A
control support optimization subject to the wave equation with terminal
contraints is performed in [30, 31].

In this work, we will take the support of the control $u$ as own optimization
variable $A\subseteq\Omega$. In addition, we will allow for a more general
control problem as above. The abstract problem we are interested in is:
Minimize with respect to $u\in L^{2}(\Omega)$ and measurable
$A\subseteq\Omega$ the functional

(1.4)
$J(u,A):=\frac{1}{2}\|S(\chi_{A}u)-y_{d}\|_{Y}^{2}+\int_{\Omega}g(u(x))+\chi_{A}(x)\beta(x)\,\text{{}d}x,$

where $S:L^{2}(\Omega)\to Y$ is a solution operator of a linear partial
differential equation, $Y$ is a Hilbert space, $y_{d}\in Y$ is given,
$g:\mathbb{R}\to\bar{\mathbb{R}}$ is a strongly convex function, and $\beta\in
L^{1}(\Omega)$. Then for fixed measurable $A\subseteq\Omega$ the map $u\mapsto
J(u,A)$ is convex, which is crucial for our analysis. As already mentioned, we
can allow for time-dependent differential equations as well, see Examples 2.2
and 2.3 below.

The main contribution of this paper is to derive an algorithm to solve (1.4).
In contrast to existing approaches [27, 28, 39], the algorithm proposed in
this paper produces a minimizing sequence for (1.4). Such an algorithm is not
available in the literature.

Given $A\subseteq\Omega$ measurable, the functional $u\mapsto J(u,A)$ admits
minimizers, and we can study the value function

(1.5) $J(A):=\min_{u\in L^{2}(\Omega)}J(u,A),$

where the minimization is carried out over measurable sets $A\subseteq\Omega$.
We will investigate topological derivatives of the value function. In
additions, we are interested in the shape optimization problem

(1.6) $\min_{A\subseteq\Omega}J(A).$

The topological derivative $DJ(A)$ is the main result of Theorem 4.2. It can
be extended to non-strongly convex $g$, see Theorem 5.9. These results
generalize available results in the literature [3, 28, 30, 31], as we allow
for non-smooth $g$ and incorporate control constraints. In comparison to [3,
28], we will use less smoothness assumptions, in particular no continuity of
controls and adjoints is required. An optimality condition for (1.6) can be
given in terms of the topological derivatives as follows: If $B$ is a solution
of (1.6) then $DJ(B)\geq 0$, see Theorem 4.4.

The concept of topological derivatives goes back to the seminal work [36]. It
was applied to an optimal control problem in [37], where the topological
derivative with respect to changes of the domain but not of the control domain
was computed. In these works, asymptotic analysis with respect to radius of
small inclusions/exclusions was performed. For an introduction and overview of
available results regarding topological derivatives, we refer to the
monographs [32, 33] and the recent introductory exposition [4].

While topological perturbations of source terms is a well understood topic,
this is not true for topological perturbations of the control domain in
control problems like (1.6). The topological derivative of a value function of
an optimal control problem subject to the wave equation with terminal
constraints was given in [30, 31] for $A=\Omega$ without proof. Topological
derivatives of value functions of optimal control problems were derived in
[28] for problems without control constraints, however the result and its
proof are wrong. In addition, the underlying abstract theory only allows to
compute the topological derivative at one fixed point, which necessitates
continuity assumptions on that point $x\in\Omega$. In our proof, we get the
topological derivative at almost all $x\in\Omega$ at once using the Lebesgue
differentiation theorem. Moreover, we can allow for control constraints and
non-smooth functions $g$.

In addition to the development of the topological derivative, we also
investigate an algorithm to solve the problem at hand. In the algorithm,
variations of a given set $A_{k}\subseteq\Omega$ are performed at points,
where the topological derivative $DJ(A_{k})$ has the wrong sign. Let
$\rho_{k}:=DJ(A_{k})^{-}$ be the residual in the optimality condition for
(1.6), and set $R_{k}:=\\{x:\ \rho_{k}\neq 0\\}$. Then the new iterate
$A_{k+1}$ is defined as $A_{k+1}:=A_{k}\triangle D_{k,t}$, where $\triangle$
denotes the symmetric difference of sets. The set $D_{k,t}\subseteq R_{k}$ is
chosen such that

$\|\rho_{k}\|_{L^{1}(D_{k,t})}\geq
t\|\rho_{k}\|_{L^{1}(\Omega)},\qquad|D_{k,t}|\leq t|R_{k}|,$

where $t\in(0,1]$ is determined by a linesearch strategy to guarantee a
sufficient decrease of $J(A_{k+1})$. These choices enable a satisfying
convergence theory: the method produces a minimizing sequence, see Theorem
6.12. A related algorithm can be found in [13, 14]. We comment on the
differences to our method in Remark 6.1.

Our choice of $D_{k,t}$ is different to the classical topological optimization
algorithm, where one sets

$D_{k,t}:=\\{x:\ |\rho_{k}(x)|\geq t\|\rho_{k}\|_{L^{\infty}(\Omega)}\\}.$

The parameter $t$ is determined to enforce a volume constraint, see, e.g.,
[12, 20, 21], or to ensure decrease of the functional, e.g., [10, 26, 23, 24].
An alternative algorithmic idea is to introduce topological changes near local
maxima of $|\rho_{k}|$, see e.g., [1, 7]. No convergence results are given in
these works.

Another method was introduced in [5]: there a (simplified) level-set method
was suggested, where the evolution of the level-set function is done using the
topological derivative. This method is applied to a wide variety of problems,
see e.g., [3, 5, 19, 34]. An convergence proof can be found in [3]. However,
for the proof the functional has to be replaced by its $H^{s}$-lower
semicontinuous envelope.

Topological optimization problems are related to binary control or
0-1-optimization problems. This connection is exploited in [2, 6]. An trust-
region method to solve such problems is analyzed in the recent contributions
[22, 29], where it is proven that a certain criticality measure converges to
zero during the iteration, which is comparable to our result.

Let us emphasize that the convergence analysis in this paper is enabled by the
particular structure of the problem. It is expected that the analysis carries
over to related problems (e.g., source identification problems). However, it
is unclear how to transfer these results to harder problems, where the
optimization variables appear in the main part of the operator as in, e.g.,
material or topology optimization problems. In fact, the question of
convergence of topological derivative-based methods is mentioned as an open
problem in [34, Section 5].

#### Notation

We will denote the Lebesgue measure of a measurable set
$A\subseteq\mathbb{R}^{d}$ by $|A|$. For $r>0$ and $x\in\mathbb{R}^{d}$, let
$B_{r}(x)$ be the open ball with radius $r$ centered at $x$. Its Lebesgue
measure will be denoted by $|B_{r}|$. We set
$\bar{\mathbb{R}}:=\mathbb{R}\cup\\{+\infty\\}$. For a function
$g:\mathbb{R}\to\bar{\mathbb{R}}$, we set
$\operatorname{dom}g:=\\{u\in\mathbb{R}:\ g(u)<+\infty\\}$. The
subdifferential of a convex function $g$ at $u$ will be denoted by $\partial
g(u)$. We will write $x^{+}:=\max(x,0)$ and $x^{-}:=\min(x,0)$ for
$x\in\mathbb{R}$.

#### Convention

In the following, we will always take Lebesgue measurable subsets of $\Omega$
only, without explicitly mentioning.

##  2 Assumptions and preliminary results

Throughout this paper, we will work with the following assumptions concerning
the problem (1.4)

  1. (A1)

$\Omega\subseteq\mathbb{R}^{d}$ is Lebesgue measurable with $|\Omega|<\infty$.

  2. (A2)

$Y$ is a real Hilbert space, $S\in\mathcal{L}(L^{2}(\Omega),Y)$, $y_{d}\in Y$.

  3. (A3)

$g:\mathbb{R}\to\bar{\mathbb{R}}$ is proper, convex, lower semi-continuous. In
addition, $g(u)\geq 0$ for all $u\in\mathbb{R}$, and $g(u)=0$ if and only if
$u=0$.

  4. (A4)

There is $\mu>0$ such that

$\frac{\mu}{2}\lambda(1-\lambda)|u-v|^{2}+g(\lambda u+(1-\lambda)v)\leq\lambda
g(u)+(1-\lambda)g(v)$

for all $u,v\in\operatorname{dom}g$, $\lambda\in(0,1)$.

  5. (A5)

There is $q>6$ such that $S^{*}S\in\mathcal{L}(L^{2}(\Omega),L^{q}(\Omega))$
and $S^{*}y_{d}\in L^{q}(\Omega)$, where
$S^{*}\in\mathcal{L}(Y,L^{2}(\Omega))$ denotes the Hilbert space-adjoint of
$S$.

  6. (A6)

$\beta\in L^{1}(\Omega)$.

Let us comment on these assumptions. As we plan to use the Lebesgue
differentiation theorem, we assume that the underlying measure space is
induced by the Lebesgue measure of $\mathbb{R}^{d}$ in (A1). Conditions (A2),
(A3), (A4) imply the well-posedness of the problem $\min_{u\in
L^{2}(\Omega)}J(u,A)$ for fixed $A$. Assumption (A4) is strong convexity of
the function $g$. The results of the paper are still valid in the non-strongly
convex case ($\mu=0$) under slightly strengthened assumptions on $g$ and $S$,
we will comment on this in Section 5. Condition (A5) is a mild assumption on
the smoothing properties of the operators $S$ and $S^{*}$. It implies that
certain remainder terms in the expansion of topological derivatives are of
higher order, see Theorem 4.2.

We will explicitly mention in upcoming, important results (theorems and
propositions), which of these assumptions are used. If the strong convexity
assumption is not mentioned then $\mu$ can be taken equal to zero.

######  Example 2.1.

Let us comment on the fulfillment of these assumptions for the introductory
example (1.1)–(1.3). Let $\Omega\subseteq\mathbb{R}^{d}$ be a bounded domain.
The partial differential equation (1.2) is uniquely solvable in the weak
sense, so that the mapping $S:u\mapsto y$ is linear and continuous from
$L^{2}(\Omega)$ to $H^{1}_{0}(\Omega)$. The functional (1.1) requires the
choice $Y:=L^{2}(\Omega)$, so that $S$ will be considered as operator on
$L^{2}(\Omega)$, which makes $S$ self-adjoint. Due to the classical result
[38], $S$ is continuous from $L^{2}(\Omega)$ to $L^{\infty}(\Omega)$ if $d\leq
3$. By [9, Theorem 18], $S$ and $S^{*}$ are in
$\mathcal{L}(L^{2}(\Omega),L^{10/3}(\Omega))$ and
$\mathcal{L}(L^{10/3}(\Omega)),L^{10}(\Omega))$ for all $d\leq 10$. And (A2)
and (A5) are satisfied for this example provided $y_{d}\in L^{10/3}(\Omega)$
and $d\leq 10$.

######  Example 2.2.

The following distributed control problem subject to a parabolic equation can
also be put into the framework above: Minimize
$\int_{Q}\frac{1}{2}(y-y_{d})^{2}+g(u)\,\text{{}d}x\,\text{{}d}t$ subject to
the parabolic equation $y_{t}-\Delta y=u$ in $Q$, $y=0$ on
$(0,T)\times\partial D$, and $y(\cdot,0)=0$ on $D$, where
$D\subseteq\mathbb{R}^{d}$ is a bounded domain, $T>0$, and $Q:=(0,T)\times D$.
We set $Y:=L^{2}(Q)$. The corresponding solution operator $S$ is continuous
from $L^{2}(Q)$ to $L^{2}(0,T;H^{1}_{0}(D))\cap H^{1}(0,T;H^{-1}(D))$. Its
adjoint operator $S^{*}$ is given as the solution operator of the adjoint
equation, i.e., $p=S^{*}z$, where $p$ solves $-p_{t}-\Delta p=z$, $p(T)=0$.
According to [11, Theorem 2.3], $S$ and $S^{*}$ are in
$\mathcal{L}(L^{2}(Q),L^{10/3}(Q))$ and $\mathcal{L}(L^{10/3}(Q)),L^{10}(Q))$
for all $d\leq 8$. If $y_{d}\in L^{10/3}(Q)$ then (A2) and (A5) are fulfilled
for all $d\leq 8$. Note that in this example the control domain $\Omega$ has
to be set to $\Omega:=Q=(0,T)\times D$.

######  Example 2.3.

One can also consider distributed control of the wave equation, where we are
interested in minimizing the same functional as in Example 2.2 subject to the
wave equation $y_{tt}-\Delta y=u$ in $Q$, $y=0$ on $(0,T)\times\partial D$,
and $y(\cdot,0)=y_{t}(\cdot,0)=0$ on $D$, where $D\subseteq\mathbb{R}^{d}$ is
a bounded domain, $T>0$, and $Q:=(0,T)\times D$. Given $u\in L^{2}(Q)$, there
is a unique weak solution $y\in L^{\infty}(0,T;H^{1}_{0}(D))$, where
$L^{\infty}(0,T;H^{1}_{0}(D))$ is continuously embedded into $L^{q}(Q)$ for
all $q<\infty$ if $d\leq 2$. And (A2) and (A5) are fulfilled for $d\leq 2$.
Using improved regularity results and Strichartz estimates, see, e.g., [18,
Section 7.2], it should be possible to relax the requirement on $d$.

###  2.1 Existence of minimizers of $J$ for fixed $A$

Let $A\subseteq\Omega$ be given. Here, we consider the problem

(PA) $\min_{u\in L^{2}(\Omega)}J(u,A).$

where $J$ is given by (1.4). Note that due to the construction of $J$ and
(A3), we have

(2.1) $J(\chi_{A}u,A)\leq J(u,A)$

for all $u\in L^{2}(\Omega)$. Due to strong convexity of $g$ and $g(0)=0$ by
(A3), (A4), we have

(2.2) $g(u)\geq\mu|u|^{2}\quad\forall u\in\operatorname{dom}g.$

######  Proposition 2.4.

Assume (A1), (A2), (A3), (A4). Let $A\subseteq\Omega$ be given. Then there is
a uniquely determined minimizer $u_{A}$ of (PA). Moreover,

(2.3) $\chi_{A}u_{A}=u_{A}.$

######  Proof 2.5.

Due to (2.2), minimizing sequences of $J(\cdot,A)$ are bounded in
$L^{2}(\Omega)$. In addition, $u\mapsto J(u,A)$ is weakly lower semi-
continuous from $L^{2}(\Omega)$ to $\mathbb{R}$ because of (A2) and (A3). The
existence of solutions follows now by standard arguments. Uniqueness of
solutions is a consequence of strong convexity of g (A4). The last claim
follows from (2.1).

In all what follows, we will not make use of the unique solvability of (PA).
We will just use that $u_{A}$ is any solution of (PA).

###  2.2 Optimality conditions for (PA)

Let $A\subseteq\Omega$ be given, and let $u_{A}$ be a solution of (PA). Let
us denote the associated state by

(2.4) $y_{A}:=S(\chi_{A}u)$

and adjoint state by

(2.5) $p_{A}:=S^{*}(y_{A}-y_{d})=S^{*}(S(\chi_{A}u)-y_{d}).$

Let $u\in L^{2}(\Omega)$ and $B\subseteq\Omega$ be given. Let
$y:=S(\chi_{B}u)$. Then by elementary calculations, we find

(2.6)
$\displaystyle\frac{1}{2}\|y-y_{d}\|_{Y}^{2}-\frac{1}{2}\|y_{A}-y_{d}\|_{Y}^{2}$
$\displaystyle=(y_{A}-y_{d},y-y_{A})_{Y}+\frac{1}{2}\|y-y_{A}\|_{Y}^{2}$
$\displaystyle=(p_{A},\chi_{B}u-\chi_{A}u_{A})+\frac{1}{2}\|y-y_{A}\|_{Y}^{2}.$

For $B=A$, we get

$\frac{1}{2}\|y-y_{d}\|_{Y}^{2}-\frac{1}{2}\|y_{A}-y_{d}\|_{Y}^{2}=(p_{A},\chi_{A}(u-u_{A}))+\frac{1}{2}\|y-y_{A}\|_{Y}^{2}.$

Hence, $\chi_{A}p_{A}\in L^{2}(\Omega)$ is the Fréchet derivative of
$u\mapsto\frac{1}{2}\|S(\chi_{A}u)-y_{d}\|_{Y}^{2}$ at $u_{A}$.

######  Proposition 2.6.

Assume (A1), (A2), (A3). Let $A\subseteq\Omega$ and let $u_{A}$ be a solution
of (PA). Let $p_{A}$ be given by (2.5). Then it holds

(2.7) $-\chi_{A}(x)p_{A}(x)\in\partial g(u(x))\quad\text{ for almost all
}x\in\Omega$

and

(2.8)
$u_{A}(x)=\operatorname*{arg\,min}_{u\in\mathbb{R}}\chi_{A}(x)p_{A}(x)\cdot
u+g(u)\quad\text{ for almost all }x\in\Omega.$

######  Proof 2.7.

Let us denote $G(u):=\int_{\Omega}g(u(x))\,\text{{}d}x$. As argued above,
$\chi_{A}p_{A}\in L^{2}(\Omega)$ is the Fréchet derivative of
$u\mapsto\frac{1}{2}\|S(\chi_{A}u)-y_{d}\|_{Y}^{2}$ at $u_{A}$. Then by well-
known results, see, e.g., [17, Proposition II.2.2], we get $-p_{A}\in\partial
G(u_{A})$. Using [35, Theorem 3A], this is equivalent to the pointwise a.e.
inclusion (2.7), which in turn is equivalent to (2.8).

Condition (2.8) can be interpreted as Pontryagin’s maximum principle for
(PA).

###  2.3 Boundedness results for solutions of (PA)

In this section, we will derive bounds on $(u_{A},y_{A},p_{A})$ that are
uniform with respect to $A\subseteq\Omega$.

######  Lemma 2.8.

There is $M>0$ such that

$\|y_{A}-y_{d}\|_{Y}+\|u_{A}\|_{L^{2}(\Omega)}\leq M$

for all $A\subseteq\Omega$.

######  Proof 2.9.

This follows directly from $J(A,u_{A})\leq J(A,0)$ and (2.2).

######  Corollary 2.10.

There is $P>0$ such that

$\|u_{A}\|_{L^{q}(\Omega)}\leq P,\quad\|p_{A}\|_{L^{q}(\Omega)}\leq P$

for all $A\subseteq\Omega$, where $q$ is from (A5).

######  Proof 2.11.

First, we have

$\displaystyle\|p_{A}\|_{L^{q}(\Omega)}$
$\displaystyle\leq\|S^{*}S\|_{\mathcal{L}(L^{2}(\Omega),L^{q}(\Omega))}\|u_{A}\|_{L^{2}(\Omega)}+\|S^{*}y_{d}\|_{L^{q}(\Omega)}$
$\displaystyle\leq\|S^{*}S\|_{\mathcal{L}(L^{2}(\Omega),L^{q}(\Omega))}M+\|S^{*}y_{d}\|_{L^{q}(\Omega)}$

by (A5) with $M$ as in Lemma 2.8. Using (2.8) with $u=0$, (A3), and (2.2), we
have for almost all $x\in\Omega$

$\mu|u_{A}(x)|^{2}\leq g(u_{A}(x))\leq-\chi_{A}(x)p_{A}(x)u_{A}(x)$

which implies $\mu|u_{A}(x)|\leq|p_{A}(x)|$ and
$\|u_{A}\|_{L^{q}(\Omega)}\leq\mu^{-1}\|p_{A}\|_{L^{q}(\Omega)}$.

##  3 Analysis of the value function

In this section, we will investigate stability properties of
$A\mapsto(u_{A},y_{A},p_{A})$, where $y_{A}$ and $p_{A}$ solve (2.4) and
(2.5). The goal is to derive formulas for the topological derivative of
$A\mapsto J(A)$, where $J(A)$ is the value function defined in (1.5) by

$J(A)=\min_{u\in L^{2}(\Omega)}J(u,A).$

For brevity, we refer to tuples $(u_{A},y_{A},p_{A})$, where $u_{A}$ solves
(PA) and $y_{A},p_{A}$ are given by (2.4) and (2.5) as solutions of (PA).

###  3.1 Sensitivity analysis of (PA) with respect to $A$

Let us start with the following preliminary expansion.

######  Lemma 3.1.

Let $A,B\subseteq\Omega$, and let $(u_{A},y_{A},p_{A})$ and
$(u_{B},y_{B},p_{B})$ be solutions of (PA) and (PB). Then it holds

$J(A,u_{A})-J(B,u_{B})+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\
=\int_{\Omega}g(u_{A})-g(u_{B})+\chi_{A}p_{A}(u_{A}-u_{B})+(\chi_{A}-\chi_{B})(\beta+p_{A}u_{B})\,\text{{}d}x.$

######  Proof 3.2.

Doing the expansion of $y\mapsto\frac{1}{2}\|y-y_{d}\|_{Y}^{2}$ similarly as
in (2.6), we have

(3.1) $J(A,u_{A})-J(B,u_{B})+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\
=\int_{\Omega}g(u_{A})-g(u_{B})+p_{A}(\chi_{A}u_{A}-\chi_{B}u_{B})+(\chi_{A}-\chi_{B})\beta\,\text{{}d}x.$

In addition, we have

$\int_{\Omega}p_{A}(\chi_{A}u_{A}-\chi_{B}u_{B})\,\text{{}d}x=\int_{\Omega}\chi_{A}p_{A}(u_{A}-u_{B})-(\chi_{B}-\chi_{A})p_{A}u_{B}\,\text{{}d}x,$

which is the claim.

######  Lemma 3.3.

Let $A,B\subseteq\Omega$, and let $(u_{A},y_{A},p_{A})$ and
$(u_{B},y_{B},p_{B})$ be solutions of (PA) and (PB). Then it holds

$\mu\|u_{B}-u_{A}\|_{L^{2}(\Omega)}^{2}+\|y_{B}-y_{A}\|_{Y}^{2}\leq\int_{\Omega}(\chi_{A}-\chi_{B})(p_{A}u_{B}-p_{B}u_{A})\,\text{{}d}x$

with $\mu\geq 0$ as in (A4).

######  Proof 3.4.

Due to Lemma 3.1, we have

$J(A,u_{A})-J(B,u_{B})+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\
=\int_{\Omega}g(u_{A})-g(u_{B})+\chi_{A}p_{A}(u_{A}-u_{B})+(\chi_{A}-\chi_{B})(\beta+p_{A}u_{B})\,\text{{}d}x$

as well as

$J(B,u_{B})-J(A,u_{A})+\frac{1}{2}\|y_{A}-y_{B}\|_{Y}^{2}\\\
=\int_{\Omega}g(u_{B})-g(u_{A})+\chi_{B}p_{B}(u_{B}-u_{A})+(\chi_{B}-\chi_{A})(\beta+p_{B}u_{A})\,\text{{}d}x.$

Adding both equations gives

$\|y_{A}-y_{B}\|_{Y}^{2}=\int_{\Omega}(\chi_{A}p_{A}-\chi_{B}p_{B})(u_{A}-u_{B})+(\chi_{A}-\chi_{B})(p_{A}u_{B}-p_{B}u_{A})\,\text{{}d}x.$

Due the inequality in (A4) and the optimality condition (2.7), we have

(3.2)
$\int_{\Omega}(\chi_{A}p_{A}-\chi_{B}p_{B})(u_{A}-u_{B})\leq-\mu\|u_{A}-u_{B}\|_{L^{2}(\Omega)}^{2},$

and the claim is proven.

Note that the previous result remains true with $\mu=0$ in the non-strongly
convex case. Now we can prove the main result of this section, which is a
stability estimate of solutions of (PA) with respect to variations of $A$ (or
$\chi_{A}$). In the proof, we will use the fact that for characteristic
functions

$\|\chi_{A}-\chi_{B}\|_{L^{s}(\Omega)}=\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{\frac{1}{s}}\quad\forall
s\in(1,\infty).$

######  Theorem 3.5.

Assume (A1), (A2), (A3), (A4), (A5). Then there is a constant $K>0$ such that
for all $A,B\subseteq\Omega$

$\|p_{A}-p_{B}\|_{L^{q}(\Omega)}+\|u_{A}-u_{B}\|_{L^{2}(\Omega)}+\|y_{B}-y_{A}\|_{L^{2}(\Omega)}\leq
K\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{\frac{1}{2}-\frac{1}{q}},$

where $(u_{A},y_{A},p_{A})$ and $(u_{B},y_{B},p_{B})$ are solutions of (PA)
and (PB), and $q$ is from (A5).

######  Proof 3.6.

From (A5), we find

$\|p_{A}-p_{B}\|_{L^{q}(\Omega)}\leq\|S^{*}S\|_{\mathcal{L}(L^{2}(\Omega),L^{q}(\Omega))}\|u_{A}-u_{B}\|_{L^{2}(\Omega)}.$

Define
$\mu^{\prime}:=\mu/\|S^{*}S\|_{\mathcal{L}(L^{2}(\Omega),L^{q}(\Omega))}^{2}$.
Let $s$ be such that $\frac{1}{s}+\frac{1}{q}+\frac{1}{2}=1$. From the
inequality of Lemma 3.3, we obtain with Hölder’s inequality

$\frac{\mu^{\prime}}{2}\|p_{A}-p_{B}\|_{L^{q}(\Omega)}^{2}+\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(\Omega)}^{2}+\|y_{B}-y_{A}\|_{Y}^{2}\\\
\begin{aligned}
&\leq\mu\|u_{B}-u_{A}\|_{L^{2}(\Omega)}^{2}+\|y_{B}-y_{A}\|_{Y}^{2}\\\
&\leq\int_{\Omega}(\chi_{A}-\chi_{B})(p_{A}u_{B}-p_{B}u_{A})\,\text{{}d}x\\\
&\leq\|\chi_{A}-\chi_{B}\|_{L^{s}(\Omega)}(\|p_{A}\|_{L^{q}(\Omega)}\|u_{B}-u_{A}\|_{L^{2}(\Omega)}+\|p_{A}-p_{B}\|_{L^{q}(\Omega)}\|u_{A}\|_{L^{2}(\Omega)})\\\
&\leq(P+M)\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{\frac{1}{s}}(\|u_{B}-u_{A}\|_{L^{2}(\Omega)}+\|p_{A}-p_{B}\|_{L^{q}(\Omega)}),\end{aligned}$

where $P$ and $M$ are from Corollary 2.10 and Lemma 2.8, and the claim is
proven.

###  3.2 Expansions of the value function

Let us define $H:\mathbb{R}\times\mathbb{R}\to\bar{\mathbb{R}}$ by

$H(u,p):=p\cdot u+g(u).$

This function reminds of the Hamiltonian of optimal control problems. In the
sequel, we need its infimum with respect to $u$,

$\min_{u\in\mathbb{R}}H(u,p)=\min_{u\in\mathbb{R}}\left(p\cdot
u+g(u)\right)=-\sup_{u\in\mathbb{R}}\left(-p\cdot u-g(u)\right)=-g^{*}(-p),$

where $g^{*}$ is the convex conjugate to $g$. The existence of this minimum
follows from the properties of $g$ in (A3) and the coercivity estimate (2.2).
Let us denote this function by $\bar{H}$, i.e.,

$\bar{H}(p):=\min_{u\in\mathbb{R}}H(u,p)=-g^{*}(-p).$

If $(u_{A},y_{A},p_{A})$ is a solution of (PA) then we have
$\bar{H}(p_{A}(x))=H(u_{A}(x),p_{A}(x))$ for almost all $x\in A$ by (2.8). We
will need some Lipschitz estimates of $\bar{H}$.

######  Lemma 3.7.

Let $A,B\subseteq\Omega$, and let $(u_{A},y_{A},p_{A})$ and
$(u_{B},y_{B},p_{B})$ be solutions of (PA) and (PB). Then we have

$\|\bar{H}(p_{A})-\bar{H}(p_{B})\|_{L^{q/2}(\Omega)}\leq
P\|p_{A}-p_{B}\|_{L^{q}(\Omega)}\leq
PK\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{\frac{1}{2}-\frac{1}{q}},$

where $P$ and $K$ are from Corollary 2.10 and Theorem 3.5, respectively.

######  Proof 3.8.

Let $p_{1},p_{2}\in\mathbb{R}$ be given. Let
$u_{i}=\operatorname*{arg\,min}_{v\in[u_{a},u_{b}]}H(p_{i},v)$ for $i=1,2$.
Then we get by the properties of $\bar{H}$

$\bar{H}(p_{1})\leq
H(p_{1},u_{2})=(p_{1}-p_{2})u_{2}+H(p_{2},u_{2})=(p_{1}-p_{2})u_{2}+\bar{H}(p_{2}).$

This implies

$\bar{H}(p_{2})\leq-(p_{1}-p_{2})u_{1}+\bar{H}(p_{1})$

by exchanging $(p_{1},u_{1}$) and $(p_{2},u_{2})$ in the above estimate.
Summarizing, we obtain

$|\bar{H}(p_{1})-\bar{H}(p_{2})|\leq|p_{1}-p_{2}|\max(|u_{1}|,|u_{2}|).$

Using Corollary 2.10 yields the claim.

We will proceed with the following expansion of the value function. Note that
in the non-strongly convex case the claim is valid with $\mu=0$.

######  Lemma 3.9.

Let $A,B\subseteq\Omega$, and let $(u_{A},y_{A},p_{A})$ and
$(u_{B},y_{B},p_{B})$ be solutions of (PA) and (PB). Then it holds

$J(A,u_{A})-J(B,u_{B})+\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(A\cap
B)}^{2}+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\
\leq\int_{\Omega}(\chi_{A}-\chi_{B})(\beta+\bar{H}(p_{A}))\,\text{{}d}x,$

where $\mu\geq 0$ as in (A4).

######  Proof 3.10.

From Lemma 3.1 we get

(3.3) $J(A,u_{A})-J(B,u_{B})+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\
=\int_{\Omega}g(u_{A})-g(u_{B})+\chi_{A}p_{A}(u_{A}-u_{B})+(\chi_{A}-\chi_{B})(\beta+p_{A}u_{B})\,\text{{}d}x.$

We will now split the integral on the right-hand side into integrals on $A\cap
B$, $A\setminus B$, and $B\setminus A$. This is sufficient as the integrand
vanishes outside of $A\cup B$. For the integral on $A\cap B$, we can use the
optimality condition (2.7) and the (possibly strong) convexity of $g$ to
obtain

$\int_{A\cap
B}g(u_{A})-g(u_{B})+\chi_{A}p_{A}(u_{A}-u_{B})\,\text{{}d}x\leq-\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(A\cap
B)}^{2}.$

Moreover, $u_{B}$ vanishes on $A\setminus B$, while $u_{A}$ vanishes on
$B\setminus A$. This allows to simplify

(3.4)
$\int_{\Omega}g(u_{A})-g(u_{B})+\chi_{A}p_{A}(u_{A}-u_{B})+(\chi_{A}-\chi_{B})p_{A}u_{B}\,\text{{}d}x\\\
\leq-\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(A\cap B)}^{2}+\int_{A\setminus
B}g(u_{A})+p_{A}u_{A}\,\text{{}d}x-\int_{B\setminus
A}g(u_{B})+p_{A}u_{B}\,\text{{}d}x$

Using $H$ and $\bar{H}$, we can write

(3.5) $\int_{A\setminus B}g(u_{A})+p_{A}u_{A}\,\text{{}d}x-\int_{B\setminus
A}g(u_{B})+p_{A}u_{B}\,\text{{}d}x\\\ =\int_{A\setminus
B}\bar{H}(p_{A})\,\text{{}d}x-\int_{B\setminus
A}H(u_{B},p_{A})\,\text{{}d}x\leq\int_{A\setminus
B}\bar{H}(p_{A})\,\text{{}d}x-\int_{B\setminus A}\bar{H}(p_{A})\,\text{{}d}x.$

Applying (3.4) and (3.5), in (3.3), results in the upper bound

$J(A,u_{A})-J(B,u_{B})+\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(A\cap
B)}^{2}+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\
\leq\int_{\Omega}(\chi_{A}-\chi_{B})(\beta+\bar{H}(p_{A}))\,\text{{}d}x,$

which is the claim.

The next result gives an expansion of the value function $J(A)$ together with
an remainder term that is of higher order in
$\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}$.

######  Theorem 3.11.

Assume (A1), (A2), (A3), (A4), (A5), (A6). Let $A,B\subseteq\Omega$, and let
$(u_{A},y_{A},p_{A})$ and $(u_{B},y_{B},p_{B})$ be solutions of (PA) and
(PB). Then it holds

$\left|J(A,u_{A})-J(B,u_{B})-\int_{\Omega}(\chi_{A}-\chi_{B})(\beta+\bar{H}(p_{B}))\,\text{{}d}x\right|\\\
+\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(A\cap
B)}^{2}+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\ \leq
PK\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{3(\frac{1}{2}-\frac{1}{q})},$

where $P$, $K$, $q$ are from Corollary 2.10, Theorem 3.5, and (A5),
respectively.

######  Proof 3.12.

Using the result of Lemma 3.9, we get

$J(A,u_{A})-J(B,u_{B})+\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(A\cap
B)}^{2}+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\
\leq\int_{\Omega}(\chi_{A}-\chi_{B})(\beta+\bar{H}(p_{B})-\bar{H}(p_{B})+\bar{H}(p_{A}))\,\text{{}d}x.$

Using Lemma 3.7 and Theorem 3.5, we can estimate the integral involving
$\bar{H}(p_{A})-\bar{H}(p_{B})$ as

(3.6)
$\int_{\Omega}(\chi_{A}-\chi_{B})(\bar{H}(p_{A})-\bar{H}(p_{B}))\,\text{{}d}x\\\
\leq\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{1-\frac{2}{q}}\cdot
P\|p_{A}-p_{B}\|_{L^{q}(\Omega)}\leq
PK\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{3(\frac{1}{2}-\frac{1}{q})}.$

This results in the upper bound

$J(A,u_{A})-J(B,u_{B})+\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(A\cap
B)}^{2}+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\
\leq\int_{\Omega}(\chi_{A}-\chi_{B})(\beta+\bar{H}(p_{B}))\,\text{{}d}x+PK\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{3(\frac{1}{2}-\frac{1}{q})}.$

To obtain a lower bound, we use the result of Lemma 3.9 but with the roles of
$A$ and $B$ reversed (and multiplying the resulting inequality by $-1$), which
yields

$J(A,u_{A})-J(B,u_{B})-\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(A\cap
B)}^{2}-\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\
\geq\int_{\Omega}(\chi_{A}-\chi_{B})(\beta+\bar{H}(p_{B}))\,\text{{}d}x.$

Both inequalities together prove the claim.

As a by-product of the previous proof, we get the strengthened stability
estimate

$\|u_{B}-u_{A}\|_{L^{2}(A\cap B)}+\|y_{B}-y_{A}\|_{Y}\leq
K^{\prime}\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{\frac{3}{2}(\frac{1}{2}-\frac{1}{q})},$

which improves the exponent from Theorem 3.5 by a factor $\frac{3}{2}$.

######  Remark 3.13.

If $S^{*}\in\mathcal{L}(Y,L^{q}(\Omega))$ then the estimate can improved to

$\|u_{B}-u_{A}\|_{L^{2}(A\cap B)}+\|y_{B}-y_{A}\|_{Y}\leq
K^{\prime}\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{2(\frac{1}{2}-\frac{1}{q})}$

by estimating $\|p_{A}-p_{B}\|_{L^{q}(\Omega)}$ against $\|y_{A}-y_{B}\|_{Y}$
in the estimate (3.6).

##  4 Topological derivatives

######  Definition 4.1.

Let $B\subseteq\Omega$. Then the topological derivative of $J$ at $B$ at the
point $x$ is defined by

$DJ(B)(x)=\begin{cases}\displaystyle\lim_{r\searrow 0}\frac{J(B\cup
B_{r}(x))-J(B)}{|B_{r}|}&\text{ if }x\not\in B\\\ \\\
\displaystyle\lim_{r\searrow 0}\frac{J(B\setminus
B_{r}(x))-J(B)}{|B_{r}|}&\text{ if }x\in B.\end{cases}$

The existence of the topological derivative is now a consequence of the
expansion in Theorem 3.11 and the Lebesgue differentiation theorem.

######  Theorem 4.2.

Assume (A1), (A2), (A3), (A4), (A5), (A6). Let $B\subseteq\Omega$, and let
$(u_{B},y_{B},p_{B})$ be a solution of (PB).

Then for almost all $x\in\Omega$ the topological derivative $DJ(B)(x)$ exists,
and is given by

$DJ(B)(x)=\sigma(B,x)(\beta(x)+\bar{H}(p_{B}(x)))$

with

$\sigma(B,x):=\begin{cases}+1&\text{ if }x\not\in B\\\ -1&\text{ if }x\in
B.\end{cases}$

######  Proof 4.3.

Let $x_{0}\in B$. Let $r>0$. Define $A(x_{0},r):=B\setminus B_{r}(x_{0})$.
Then it follows $\chi_{A(x_{0},r)}-\chi_{B}=-\chi_{B\cap B_{r}(x_{0})}$, which
implies $\|\chi_{A(x_{0},r)}-\chi_{B}\|_{L^{1}(\Omega)}\leq|B_{r}|$. Using
this in the result of Theorem 3.11, we find

(4.1) $\left|J(A(x_{0},r))-J(B)+\int_{B\cap
B_{r}(x_{0})}\beta+\bar{H}(p_{B})\,\text{{}d}x\right|\leq
PK|B_{r}|^{3(\frac{1}{2}-\frac{1}{q})}.$

Let us now define

$v(x_{0},r):=\frac{1}{|B_{r}|}\int_{B_{r}(x_{0})}\chi_{B}\cdot(\beta+\bar{H}(p_{B}))\,\text{{}d}x.$

By the Lebesgue differentiation theorem, we have

$\lim_{r\searrow 0}v(x,r)=\chi_{B}(x)\cdot(\beta(x)+\bar{H}(p_{B}(x)))$

for almost all $x\in\Omega$. This implies together with (4.1)

$\lim_{r\searrow
0}\frac{J(A(x,r))-J(B)}{|B_{r}|}=-(\beta(x)+\bar{H}(p_{B}(x)))$

for almost all $x\in B$. Here we used that $3(\frac{1}{2}-\frac{1}{q})>1$ by
(A5). This proves the claim for $x\in B$.

The claim for $x\not\in B$ can be proven completely analogously: this time we
set $A(x_{0},r):=B\cup B_{r}(x_{0})$ for $x_{0}\not\in B$, which implies
$\chi_{A(x_{0},r)}-\chi_{B}=\chi_{B_{r}(x_{0})\setminus B}$, resulting in the
different sign of the topological derivative.

Note that in contrast to other works, we do not need to impose continuity of
$u_{B}$ near $x_{0}$ as in [28, Corollary 4.1], nor do we need to argue by
Hölder continuity of the adjoint as in [3, Corollary 3.2].

We can now formulate a necessary optimality condition for (1.6) using the
topological derivative.

######  Theorem 4.4.

Assume (A1), (A2), (A3), (A4), (A5), (A6). Let $B$ be a solution of (1.6).
Then

$DJ(B)(x)\geq 0\quad\text{ for a.a. }x\in\Omega.$

######  Proof 4.5.

The result follows immediately from Theorem 4.2.

######  Remark 4.6.

Using the celebrated Ekeland’s variational principle [16], the following
result can be proven for $\epsilon$-solutions: There is an
$\epsilon$-solution, such that optimality conditions are satisfied up to
$\epsilon$. We briefly sketch the proof.

Let $V$ be the metric space of characteristic functions $\chi_{B}$,
$B\subseteq\Omega$ measurable, supplied with the $L^{1}(\Omega)$-metric, which
makes it a complete space. Applying [16, Theorem 1.1] with $\epsilon>0$ and
$\lambda=1$ there is $B_{\epsilon}\subseteq\Omega$ such that

(4.2) $J(B_{\epsilon})\leq\inf_{B\subseteq\Omega}J(B)+\epsilon$

and

(4.3) $J(A)\geq
J(B_{\epsilon})-\epsilon\|\chi_{A}-\chi_{B_{\epsilon}}\|_{L^{1}(\Omega)}$

for all $A\subseteq\Omega$. Owing to (4.2) the set $B_{\epsilon}$ is then an
$\epsilon$-solution of (1.6). Due to inequality (4.3), we can consider
variations of $J(B_{\epsilon})$ to obtain estimates of the topological
derivative:

For $x_{0}\in\Omega$ and $r>0$, define $A(x_{0},r)$ as in the proof of Theorem
4.2. Then $\frac{1}{|B_{r}|}(J(A(x_{0},r))-J(B_{\epsilon}))\geq-\epsilon$ by
(4.3), which results in $DJ(B_{\epsilon})(x_{0})\geq-\epsilon$ for almost all
$x_{0}$. This proves the existence of an $\epsilon$-solution that satisfies
the optimality condition up to an $\epsilon$.

In addition, the defect in the optimality condition of Theorem 4.4 can be used
to get an error estimate as follows.

######  Corollary 4.7.

Assume (A1), (A2), (A3), (A6). Let $A\subseteq\Omega$, let
$(u_{A},y_{A},p_{A})$ be a solution of (PA). Let the defect $\delta_{A}$ be
defined by

$\delta_{A}:=\int_{A}\left(\beta+\bar{H}(p_{A})\right)^{+}\,\text{{}d}x-\int_{\Omega\setminus
A}\left(\beta+\bar{H}(p_{A})\right)^{-}\,\text{{}d}x=-\int_{\Omega}(DJ(B))^{-}\,\text{{}d}x.$

Then we have

$J(A)-\inf_{B\subseteq\Omega}J(B)\leq\delta_{A},$

and $A$ is a $\delta_{A}$-solution. If $B$ is a solution of (1.6) then we have
the error estimate

$J(A,u_{A})-J(B,u_{B})+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}+\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(A\cap
B)}^{2}\leq\delta_{A}.$

######  Proof 4.8.

Let $B\subseteq\Omega$ and $(u_{B},y_{B},p_{B})$ be a solution of (PB). By
Lemma 3.9, we have

$J(A,u_{A})-J(B,u_{B})+\frac{\mu}{2}\|u_{B}-u_{A}\|_{L^{2}(A\cap
B)}^{2}+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\ \begin{aligned}
&\leq\int_{\Omega}(\chi_{A}-\chi_{B})(\beta+\bar{H}(p_{A}))\,\text{{}d}x,\\\
&=\int_{A\setminus B}\beta+\bar{H}(p_{A})\,\text{{}d}x-\int_{B\setminus
A}\beta+\bar{H}(p_{A})\,\text{{}d}x\\\
&\leq\int_{A}\left(\beta+\bar{H}(p_{A})\right)^{+}\,\text{{}d}x-\int_{\Omega\setminus
A}\left(\beta+\bar{H}(p_{A})\right)^{-}\,\text{{}d}x=\delta_{A}.\end{aligned}$

If $B$ is a solution of (1.6) then the claim follows. Otherwise, we take the
supremum of $-J(B,u_{B})$ on the left-hand side.

##  5 The non-strongly convex case

Let us briefly comment on the non-strongly convex case. That is, we no longer
assume the strong convexity of $g$ as in (A4). We will replace (A4) and (A5)
by the following two assumptions.

  1. (A4’)

$\operatorname{dom}g$ is a bounded subset of $\mathbb{R}$,

  2. (A5’)

There is $q>3$ such that $S^{*}\in\mathcal{L}(Y,L^{q}(\Omega))$, where
$S^{*}\in\mathcal{L}(Y,L^{2}(\Omega))$ denotes the Hilbert space-adjoint of
$S$.

(A4’) implies the solvability of (PA). In addition, solutions $u_{A}$ of
(PA) will be in $L^{\infty}(\Omega)$. Due to the missing strong convexity, we
have to replace the assumption on $S^{*}S$ in (A5) by an assumption on
$S^{*}$. The $L^{\infty}(\Omega)$-regularity of optimal controls will allow us
to work with a smaller exponent $q$ in (A5’) when compared to (A5). Condition
(A5’) is fulfilled for Examples 2.1 and 2.2.

Note that we do not add assumptions that imply unique solvability of (PA).

######  Proposition 5.1.

Let $A\subseteq\Omega$ be given. Then there is a minimizer $u_{A}$ of (PA).
Moreover, $\chi_{A}u_{A}$ is also a minimizer of (PA).

######  Proof 5.2.

Due to (A4’) minimizing sequences of $u\mapsto J(A,u)$ are bounded in
$L^{\infty}(\Omega)$. Then the proof of existence follows as in Proposition
2.4. The last claim is a consequence of (2.1).

In the sequel, we will assume that a solution $u_{A}$ of (PA) satisfies
$\chi_{A}u_{A}=u_{A}$. Due to the previous result, this is not restriction at
all, as for every minimizer $u_{A}$ also $\chi_{A}u_{A}$ is a minimizer. Let
us start with a replacement of Lemma 2.8 and Corollary 2.10.

######  Lemma 5.3.

There is $M>0$ and $P^{\prime}>0$ such that

$\|y_{A}-y_{d}\|_{Y}\leq M$

and

$\|u_{A}\|_{L^{\infty}(\Omega)}\leq
P^{\prime},\quad\|p_{A}\|_{L^{q}(\Omega)}\leq P^{\prime}$

for all $A\subseteq\Omega$ and all solutions $(u_{A},y_{A},p_{A})$ of (PA).
Here, $q$ is as in (A5’).

######  Proof 5.4.

The bound of $y_{A}$ can be obtained as in Lemma 2.8, the bounds of $u_{A}$
and $p_{A}$ are consequences of (A4’) and (A5’).

Due to the missing strong convexity of $g$, we cannot expect stability of
controls as in Theorem 3.5. Here, we have the following replacement.

######  Theorem 5.5.

Assume (A1), (A2), (A3), (A4’), (A5’). Then there is a constant $K^{\prime}>0$
such that for all $A,B\subseteq\Omega$

$\|p_{A}-p_{B}\|_{L^{q}(\Omega)}+\|y_{B}-y_{A}\|_{L^{2}(\Omega)}\leq
K^{\prime}\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{\frac{1}{2}(1-\frac{1}{q})},$

where $(u_{A},y_{A},p_{A})$ and $(u_{B},y_{B},p_{B})$ are solutions of (PA)
and (PB), and $q$ is from (A5’).

######  Proof 5.6.

From (A5’), we get
$\|p_{A}-p_{B}\|_{L^{q}(\Omega)}\leq\|S^{*}\|_{\mathcal{L}(Y,L^{q}(\Omega))}\|y_{A}-y_{B}\|_{L^{2}(\Omega)}$.
Define $\mu^{\prime}:=1/\|S^{*}\|_{\mathcal{L}(Y,L^{q}(\Omega))}^{2}$. Let
$q^{\prime}$ be such that $\frac{1}{q^{\prime}}+\frac{1}{q}=1$. From the
inequality of Lemma 3.3, we obtain with Hölder’s inequality

$\displaystyle\frac{\mu^{\prime}}{2}\|p_{A}-p_{B}\|_{L^{q}(\Omega)}^{2}+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}$
$\displaystyle\leq\int_{\Omega}(\chi_{A}-\chi_{B})(p_{A}u_{B}-p_{B}u_{A})\,\text{{}d}x$
$\displaystyle\leq
2(P^{\prime})^{2}\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{1-\frac{1}{q}},$

where $P^{\prime}$ is from Lemma 5.3, and the claim is proven.

This stability result has to replace Theorem 3.5 in the proof of Theorem 3.11.
The result corresponding to the latter theorem now reads as follows. Note that
due to (A4’), $\bar{H}(p)$ is well-defined and finite for all
$p\in\mathbb{R}$.

######  Theorem 5.7.

Assume (A1), (A2), (A3), (A4’), (A5’), (A6). Let $A,B\subseteq\Omega$, and let
$(u_{A},y_{A},p_{A})$ and $(u_{B},y_{B},p_{B})$ be solutions of (PA) and
(PB). Then it holds

$\left|J(A,u_{A})-J(B,u_{B})-\int_{\Omega}(\chi_{A}-\chi_{B})(\beta+\bar{H}(p_{B}))\,\text{{}d}x\right|+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}\\\
\leq
P^{\prime}K^{\prime}\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{\frac{3}{2}(1-\frac{1}{q})},$

where $P^{\prime}$, $K^{\prime}$, $q$ are from Lemma 5.3 and (A5’),
respectively.

######  Proof 5.8.

We can proceed exactly as in the proof of Theorem 3.11 but now with $\mu=0$.
Only the estimate (3.6) has to be modified. The estimate of Lemma 3.7 has to
be changed to

(5.1) $\|\bar{H}(p_{A})-\bar{H}(p_{B})\|_{L^{q}(\Omega)}\leq
P^{\prime}\|p_{A}-p_{B}\|_{L^{q}(\Omega)}$

using the $L^{\infty}(\Omega)$-bound of optimal controls in Lemma 5.3, as well
as the estimate of $\bar{H}$ from the proof of Lemma 3.7. Then the error term
of (3.6) can be estimated using Eq. 5.1 and Theorem 5.5 as

$\int_{\Omega}(\chi_{A}-\chi_{B})(\bar{H}(p_{A})-\bar{H}(p_{B}))\,\text{{}d}x\\\
\leq\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{1-\frac{1}{q}}\cdot
P^{\prime}\|p_{A}-p_{B}\|_{L^{q}(\Omega)}\leq
P^{\prime}K^{\prime}\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{\frac{3}{2}(1-\frac{1}{q})}.$

The claimed estimate can now be obtained with the same arguments as in the
proof of Theorem 3.11.

######  Theorem 5.9.

Assume (A1), (A2), (A3), (A4’), (A5’), (A6). Let $B\subseteq\Omega$.

Then for almost all $x\in\Omega$ the topological derivative $DJ(B)(x)$ exists,
and it is given by the expression in Theorem 4.2.

##  6 Optimization method based on the topological derivative

In this section, we introduce an optimization algorithm that is motivated by
the work on the topological derivative. Here, we work under the set of
assumptions of Theorem 4.2 or Theorem 5.9.

Let $A_{k}\subseteq\Omega$ be a given iterate together with solutions of
(P${}_{A_{k}}$). Let us define the residual in the optimality condition of
Theorem 4.4 as

(6.1) $\rho_{k}:=(DJ(A_{k}))^{-}.$

Let us denote the support of $\rho_{k}$ by

$R_{k}:=\\{x:\ \rho_{k}(x)\neq 0\\}.$

New iterates $A_{k+1}$ will now be defined by adding/removing points to/from
$A_{k}$, where $\rho_{k}$ is non-zero. That is, we will choose
$D_{k,t}\subseteq R_{k}$ and denote the candidate for a new iterate by
$A_{k}\triangle D_{k,t}$. Given $t\in(0,1]$, we select $D_{k,t}$ with the
properties:

(6.2) $\displaystyle D_{k,t}$ $\displaystyle\subseteq R_{k}:$
$\displaystyle\|\rho_{k}\|_{L^{1}(D_{k,t})}$ $\displaystyle\geq
t\|\rho_{k}\|_{L^{1}(\Omega)},$ $\displaystyle|D_{k,t}|$ $\displaystyle\leq
t|R_{k}|.$

That is, we are looking for a set $D_{k,t}$ with prescribed bound on its
measure that captures a certain part of the mass of the residual $\rho_{k}$.
Sets satisfying the conditions of (6.2) exist, and can be found by solving

(6.3) $\max_{D\subseteq R_{k}:\ |D|\leq tR_{k}}\|\rho_{k}\|_{L^{1}(D)}.$

In [13], such a problem is used to generate search directions. A related
procedure to compute such sets is given in [22, Procedure 1, Lemma 9]. Note
that we need both conditions on $D_{k,t}$: the condition on
$\|\rho_{k}\|_{L^{1}(D_{k,t})}$ will give descent of values of the functional
$J$, while the condition on $|D_{k,t}|$ will help to control the error in the
expansion of the functional $J$.

######  Remark 6.1.

Let us comment on related algorithms based on topological derivatives.

In the seminal work [14] a similar idea was developed. In our notation, their
algorithm reads: find $t>0$ and $D_{k,t}\subseteq R_{k}$ such that
$|D_{k,t}|\leq t$ and $\|\rho_{k}\|_{L^{1}(D)}\geq\frac{\tilde{M}}{2}t^{2}$,
where $\tilde{M}$ is larger than the Lipschitz constant of the Fréchet
derivative of $J$, when considered as a function from $L^{1}(\Omega)$ to
$\mathbb{R}$. This choice of $D_{k,t}$ guarantees a sufficient decrease of
$J(A_{k})$. However, the knowledge of this Lipschitz constant is necessary to
implement this condition. In addition, this Lipschitz condition is not
satisfied for our problem, it would imply that the remainder term in Theorem
3.11 is of order $\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{2}$, while our
analysis works under a weaker estimate of the remainder.

Another popular algorithm is the choice $D_{k,t}=\\{x:\ |\rho_{k}(x)|\geq
t\\}$, where $t\in(0,\|\rho_{k}\|_{L^{\infty}(\Omega)})$, see, e.g., [12, 20,
23, 24], where $t$ is chosen such that the new iterate respects a volume
constraint or a descent condition. This approach is successfully used in
practice. From the viewpoint of optimization, it has the following theoretical
drawback: if $\rho_{k}$ is a constant function, then $D_{k,t}$ is either empty
or equal to $R_{k}$, and a linesearch in $t$ is not guaranteed to succeed.

An algorithm based on trust-region ideas can be found in [22, 29]. In these
works binary control problems are considered. The method proposed there can
also be applied to our problema, and would lead to similar convergence
results.

In [5] a simplified level-set method was introduced, where the level-set
function is updated using the topological derivative. The corresponding
linesearch method was analyzed in [3]. There, small values of the linesearch
parameter may only lead to boundary variations.

Let us first prove that (6.3) can be used to find sets satisfying the
condition (6.2) for fixed $t\in(0,1]$.

######  Lemma 6.2.

Let $t\in(0,1]$. There exist sets $D_{k,t}$ satisfying conditions (6.2).
Problem (6.3) is solvable, and every solution $D_{k,t}$ of (6.3) satisfies
(6.2).

######  Proof 6.3.

The existence of sets $D_{k,t}$ satisfying conditions (6.2) with equality is a
consequence of the Lyapunov convexity theorem for vector measures [15,
Corollary IX.5]. Consequently, solutions of (6.3) satisfy (6.2). Let us proof
the solvability of (6.3). Given $s\geq 0$ define

$B_{>s}:=\\{x:\ |\rho_{k}(x)|>s\\},\ B_{\geq s}:=\\{x:\ |\rho_{k}(x)|\geq
s\\}.$

Then the monotonically decreasing functions $s\mapsto|B_{>s}|$ and
$s\mapsto|B_{\geq s}|$ are continuous from the right and from the left,
respectively. In addition, we have
$\lim_{s\to+\infty}|B_{>s}|=\lim_{s\to+\infty}|B_{\geq s}|=0$.

Assume that there is $s\geq 0$ and a set $B$ with $B_{>s}\subseteq B\subseteq
B_{\geq s}$ and $|B|=tR_{k}$. We will now argue that $B$ is a solution of
(6.3). Let us take $D\subseteq R_{k}$ with $|D|\leq tR_{k}$. Then we get using
the definitions of $B$, $B_{>s}$, $B_{\geq s}$

$\displaystyle\int_{D}|\rho_{k}|\,\text{{}d}x-\int_{B}|\rho_{k}|\,\text{{}d}x$
$\displaystyle=\int_{D\setminus B}|\rho_{k}|\,\text{{}d}x-\int_{B\setminus
D}|\rho_{k}|\,\text{{}d}x$ $\displaystyle\leq s|D\setminus B|-s|B\setminus
D|=s(|D|-|D\cap B|-|B)$ $\displaystyle=s(|D|-tR_{k})-s|D\cap B|\leq 0,$

and $B_{>s}$ solves (6.3).

Hence, if there is $s\geq 0$ such that $|B_{>s}|=tR_{k}$ or $|B_{\geq
s}|=tR_{k}$, then $B_{>s}$ or $B_{\geq s}$ is a solution of (6.3). If this is
not the case, then there is $s>0$ such that $|B_{>s}|<tR_{k}<|B_{\geq s}|$.
Since the measure space is atom-free, we can choose $C\subseteq B_{\geq
s}\setminus B_{>s}$ with $|C|=tR_{k}-|B_{>s}|$ by the Sierpinski theorem [8,
Corollary 1.2.10]. And $B_{>s}\cup C$ is a solution of (6.3) as argued above.

The symmetric difference of $A_{k}$ and $D_{k,t}$ is a candidate for the new
iterate. Let us set

(6.4) $A_{k,t}:=A_{k}\triangle D_{k,t}.$

The step-size $t$ in (6.3) will be determined by an Armijo-like line-search
incorporating the descent condition

(6.5) $J(A_{k,t})\leq
J(A_{k})+\sigma\int_{\Omega}|\chi_{A_{k,t}}-\chi_{A_{k}}|\cdot\rho_{k}\,\text{{}d}x$

where $\sigma\in(0,1)$ is a fixed parameter.

The validity of this approach is enabled by the following observation, which
shows that (6.5) ensures descent of $J$.

######  Lemma 6.4.

Let $D_{k,t}$ satisfy (6.2), let $A_{k,t}$ as in (6.4). Then it holds

$\int_{\Omega}|\chi_{A_{k,t}}-\chi_{A_{k}}|\cdot\rho_{k}\,\text{{}d}x=-\|\rho_{k}\|_{L^{1}(D_{k,t})}\leq-t\|\rho_{k}\|_{L^{1}(\Omega)}.$

######  Proof 6.5.

This follows directly from the properties of $\rho_{k}$ and $D_{k,t}$, see
(6.1).

Let us prove that there are step-sizes $t$ such that the descent condition is
satisfied. Here, the following estimate of $J(A_{k,t})-J(A_{k})$ will prove
useful, which is a consequence of the results of the previous sections.

######  Lemma 6.6.

There is $C>0$ and $\nu>0$ such that for all $A_{k}\subseteq\Omega$,
$t\in(0,1]$, $A_{k,t}$ as in (6.4), it holds

$J(A_{k,t})-J(A_{k})\leq\int_{\Omega}|\chi_{A_{k,t}}-\chi_{A_{k}}|\cdot\rho_{k}\,\text{{}d}x+C\|\chi_{A_{k,t}}-\chi_{A_{k}}\|_{L^{1}(\Omega)}^{1+\nu}.$

######  Proof 6.7.

As a consequence of the expansion of $J$ in Theorems 3.11 and 5.7 and of the
definition of $DJ$ in Theorem 4.2, we get

$J(A_{k,t})-J(A_{k})\leq\int_{\Omega}(\chi_{A_{k,t}}-\chi_{A_{k}})\sigma(A_{k})DJ(A_{k})\,\text{{}d}x+C\|\chi_{A_{k,t}}-\chi_{A_{k}}\|_{L^{1}(\Omega)}^{1+\nu},$

where $\sigma(A_{k})$ is as in Theorem 4.2. Note that
$|\chi_{A_{k,t}}-\chi_{A_{k}}|=(\chi_{A_{k,t}}-\chi_{A_{k}})\sigma(A_{k})$,
which proves the claim.

Now, we are in the position to prove the existence of step-sizes such that the
descent condition (6.5) is satisfied.

######  Lemma 6.8.

Let $\sigma\in(0,1)$. Assume $\rho_{k}\neq 0$. Then there is $\tilde{t}_{k}>0$
such that for all $t\in(0,\tilde{t}_{k}]$ and all sets $D_{k,t}$ satisfying
(6.2) the set $A_{k,t}:=A_{k}\triangle D_{k,t}$ satisfies the descent
condition (6.5).

######  Proof 6.9.

Due to (6.2), we have

$\|\chi_{A_{k,t}}-\chi_{A_{k}}\|_{L^{1}(\Omega)}=|A_{k,t}\triangle
A_{k}|=|D_{k,t}|\leq t|R_{k}|.$

With the results of Lemmas 6.4 and 6.6, we get

(6.6)
$J(A_{k,t})-J(A_{k})-\sigma\int_{\Omega}|\chi_{A_{k,t}}-\chi_{A_{k}}|\cdot\rho_{k}\,\text{{}d}x\\\
\begin{aligned}
&\leq(1-\sigma)\int_{\Omega}|\chi_{A_{k,t}}-\chi_{A_{k}}|\cdot\rho_{k}\,\text{{}d}x+C\|\chi_{A_{k,t}}-\chi_{A_{k}}\|_{L^{1}(\Omega)}^{1+\nu}\\\
&\leq-(1-\sigma)\|\rho_{k}\|_{L^{1}(\Omega)}t+C|R_{k}|^{1+\nu}t^{1+\nu}.\end{aligned}$

Clearly, the right-hand side is negative for $t$ small enough.

Choose $\tau\in(0,1)$, $\sigma\in(0,1)$, $A_{0}\subseteq\Omega$,
$\delta_{\mathrm{tol}}\geq 0$. Set $k:=0$.

loop$\triangleright$ Gradient descent

Compute a solution $(u_{k},y_{k},p_{k})$ of (P${}_{A_{k}}$).

Compute $\rho_{k}$ as in (6.1).

if $\|\rho_{k}\|_{L^{1}(\Omega)}\leq\delta_{\mathrm{tol}}$ then
$\triangleright$ Termination criterion

return $A_{k}$

end if

$t:=1$.

loop$\triangleright$ Armijo line-search

Compute $D_{k,t}$ as a solution of (6.3).

Compute $J(A_{k,t})$ for $A_{k,t}=A_{k}\triangle D_{k,t}$.

if $A_{k,t}$ satisfies (6.5) then

break

end if

$t:=\tau\cdot t$.

end loop

$t_{k}:=t$. $\triangleright$ Update

$A_{k+1}:=A_{k}\triangle D_{k,t}$.

$k:=k+1$.

end loop

Algorithm 1 Topological gradient descent algorithm

The resulting algorithm is Algorithm 1. Let us comment on it in detail. The
algorithm stops if $\|\rho_{k}\|_{L^{1}(\Omega)}\leq\delta_{\mathrm{tol}}$ for
some prescribed tolerance $\delta_{\mathrm{tol}}\geq 0$. This is motivated by
Corollary 4.7: if $\|\rho_{k}\|_{L^{1}(\Omega)}$ is less than some tolerance
$\delta_{\mathrm{tol}}>0$ then $A_{k}$ is a $\delta_{\mathrm{tol}}$-solution
of (1.6). Due to Lemma 6.8, the Armijo line-search will terminate in finitely
many steps, and the algorithm is well-defined. If the algorithm does not stop
after finitely many iterations, then it will produce an infinite sequence of
sets $(A_{k})$, such that $(J(A_{k}))$ is monotonically decreasing. We have
the following basic convergence result.

######  Lemma 6.10.

Let $(A_{k})$ be an infinite sequence of iterates of Algorithm 1. Then it
holds

$\sum_{k=0}^{\infty}t_{k}\|\rho_{k}\|_{L^{1}(\Omega)}<+\infty.$

######  Proof 6.11.

Due to the descent condition (6.5) and the result of Lemma 6.4, we have the
chain of inequalities

$J(A_{k+1})-J(A_{k})\leq\sigma\int_{\Omega}|\chi_{A_{k+1}}-\chi_{A_{k}}|\cdot\rho_{k}\,\text{{}d}x\leq-\sigma
t_{k}\|\rho_{k}\|_{L^{1}(\Omega)}<0.$

Since $J$ is bounded from below, we can sum the above inequalities for
$k=0,\dots$, which proves the claim.

######  Theorem 6.12.

Let $(A_{k})$ be the iterates of Algorithm 1. Then exactly one of the
following statements is true:

  1. (i)

The algorithm returns after $m\in\mathbb{N}$ iterations with a
$\delta_{\mathrm{tol}}$-solution $A_{m}$ of (1.6).

  2. (ii)

The algorithm produces an infinite sequence $(A_{k})$ with
$\lim_{k\to\infty}\|\rho_{k}\|_{L^{1}(\Omega)}=0$, and $(A_{k})$ is a
minimizing sequence of (1.6).

In particular, if $\delta_{\mathrm{tol}}>0$ then the algorithm terminates
after finitely many iterations.

######  Proof 6.13.

Suppose the algorithm returns after $m$ iterations. Then
$\|\rho_{m}\|_{L^{1}(\Omega)}\leq\delta_{\mathrm{tol}}$, and by Corollary 4.7
the set $A_{m}$ is a $\delta_{\mathrm{tol}}$-solution $A_{m}$ of (1.6).

Now suppose that the algorithm produces an infinite sequence $(A_{k})$ of
iterates. Let $k$ be such that $t_{k}<1$. Due to the line-search procedure of
Algorithm 1, it follows that $t:=t_{k}/\tau\leq 1$ violates the descent
condition (6.5), that is

$0\leq
J(A_{k,t})-J(A_{k})-\sigma\int_{\Omega}|\chi_{A_{k,t}}-\chi_{A_{k}}|\cdot\rho_{k}.$

Using estimate (6.6), this implies

$0\leq-(1-\sigma)\|\rho_{k}\|_{L^{1}(\Omega)}t+C|R_{k}|^{1+\nu}t^{1+\nu},$

or equivalently

$(1-\sigma)\|\rho_{k}\|_{L^{1}(\Omega)}\leq
C|R_{k}|^{1+\nu}t^{\nu}=C|R_{k}|^{1+\nu}\tau^{-\nu}\ t_{k}^{\nu}.$

Note that $|R_{k}|\leq|\Omega|$. Together with the result of Lemma 6.10, we
obtain

$\sum_{k:\,t_{k}=1}\|\rho_{k}\|_{L^{1}(\Omega)}+\sum_{k:\,t_{k}<1}\|\rho_{k}\|_{L^{1}(\Omega)}^{1+\frac{1}{\nu}}<+\infty,$

which results in $\lim_{k\to\infty}\|\rho_{k}\|_{L^{1}(\Omega)}=0$. Hence, the
algorithm stops after finitely many iterations if $\delta_{\mathrm{tol}}>0$.

The sequence $(\chi_{A_{k}})$ of characteristic functions of the iterates
admits weak-star converging subsequences in
$L^{\infty}(\Omega)=L^{1}(\Omega)^{*}$. However, the corresponding limits are
not guaranteed to be characteristic functions. If a subsequence converges
weak-star to a characteristic function, i.e.,
$\chi_{A_{k^{\prime}}}\rightharpoonup^{*}\chi_{A}$ in $L^{\infty}(\Omega)$,
then the convergence is strong in every $L^{q}(\Omega)$, $q<\infty$, and the
limit $A$ is a solution of (1.6). If (1.6) is unsolvable then weak-star
sequential limit points of $(\chi_{A_{k}})$ cannot be characteristic
functions.

A similar convergence result for gradient descent without linesearch can be
found in [13, Theorem 2.1], whereas convergence of a trust-region-type
algorithm can be found in [29, Theorem 4.1].

######  Remark 6.14.

Algorithm 1 can be generalized to optimization problems of the type
$\min_{A\subseteq\Omega}J(A)$, where the minimum is taken over measurable
subsets $A$. In order to apply the above analysis, we would need the following
“differentiability” condition on $J$: there is $\eta:[0,\infty)\to[0,\infty)$
with $\lim_{t\searrow 0}\frac{\eta(t)}{t}=0$ such that for each
$A\subseteq\Omega$ there is $DJ(A)\in L^{1}(\Omega)$ such that

(6.7)
$\left|J(B)-J(A)-\int_{\Omega}(\chi_{B}-\chi_{A})DJ(A)\,\text{{}d}x\right|\leq\eta(|B\triangle
A|)\quad\forall B\subseteq\Omega.$

Then the result of Lemma 6.8 is valid, where assumption (6.7) would act as a
substitute for Lemma 6.6. The statement of Theorem 6.12 has to be changed to:
either $\rho_{m}=0$ for some finite $m$ or
$\lim_{k\to\infty}\|\rho_{k}\|_{L^{1}(\Omega)}=0$.

Note that condition (6.7) is weaker than the assumptions of [13, 14] and [29,
Assumption 1.1(a)–(c)]. In the latter reference that stronger assumption was
used to prove a statement analogous to Theorem 6.12.

##  7 Numerical experiments

###  7.1 Optimal control problem with $L^{0}$-control cost

Let us report about numerical results of the application of Algorithm 1 to the
following problem: Minimize

$\min\frac{1}{2}\|y-y_{d}\|_{L^{2}(\Omega)}^{2}+\frac{\alpha}{2}\|u\|_{L^{2}(\Omega)}^{2}+\beta\|u\|_{0}$

over all $(y,u)\in H^{1}_{0}(\Omega)\times L^{2}(\Omega)$ satisfying

$-\Delta y=u\quad\text{ a.e.\@ on }\Omega$

and

$u_{a}\leq u\leq u_{b}\quad\text{ a.e.\@ on }\Omega.$

This corresponds to the abstract setting with the choices
$S:=(-\Delta)^{-1}:L^{2}(\Omega)\to H^{1}_{0}(\Omega)\hookrightarrow
L^{2}(\Omega)$, $Y:=L^{2}(\Omega)$,
$g(u):=\frac{\alpha}{2}u^{2}+I_{[u_{a},u_{b}]}(u)$, $\beta(x):=\beta$. Here,
$I_{C}$ denotes the indicator function of the convex set $C$, defined by
$I_{C}(x)=0$ for $x\in C$ and $I_{C}(x)=+\infty$ for $x\not\in C$. The
assumptions are all satisfied. In particular $g$ is strongly convex with
modulus $\mu:=\alpha$. In the numerical experiment, we used
$\Omega\subseteq\mathbb{R}^{2}$ bounded, so that assumptions (A5) and (A5’)
are satisfied with $q=\infty$ due to Stampacchia’s result [38].

####  7.1.1 Example from [39]

We choose $\Omega=(0,1)^{2}$. We used a standard finite-element discretization
on a shape-regular mesh on $\Omega$. State and adjoint variables (i.e., $y$,
$p$) were discretized using continuous piecewise linear functions, while the
control variable was discretized using piecewise constant functions. Let us
remark that for the finest discretization, the control functions have
$2,000,000$ degrees of freedom. The subproblems (PA) were solved by a
semismooth Newton implementation. The parameters in the line-search of
Algorithm 1 were chosen to be $\tau=0.5$ and $\sigma=0.1$. The algorithm was
stopped if one of the following conditions was fulfilled:
$\|\rho_{k}\|_{L^{\infty}(\Omega)}\leq 10^{-12}$, the support of $\rho_{k}$
contained $\leq 3$ elements, or the line-search failed to find a valid step-
size. Termination due to the latter condition can happen if the relevant
quantities in (6.5), are very small so that errors in the inexact solve of the
sub-problem (PA) are of the same order.

In addition, we used the following data

$y_{d}(x_{1},x_{2})=10x_{1}\sin(5x_{1})\cos(7x_{2}),\ \alpha=0.01,\
\beta=0.01,\ u_{a}=-4,\ u_{b}=+4,$

which was also used in [27, 39]. The computed optimal control, which is
obtained by the last iterate of Algorithm 1 on the finest mesh, can be seen in
Fig. 1. Due to the presence of the $L^{0}$-term in the objective, the control
is zero on a relatively large part of $\Omega$.

Figure 1: Solution $u$ for $h=1.41\cdot 10^{-3}$, Section 7.1.1

The results of the computations for different meshes can be seen in Table 1.
There, $h$ denotes the mesh-size of the triangulation, $J$ denotes the value
of the functional $J$ at the final iterate, similarly
$\|\chi\|_{L^{1}(\Omega)}$ is the size of the support of the optimal control,
and $\|\rho\|_{L^{1}(\Omega)}$ is the error estimate from the topological
derivative at the final iteration. The values corresponding to the mesh-size
$h=2.83\cdot 10^{-3}$ are in agreement with those from [39]. For this example,
all computations stopped due to the support of $\rho_{k}$ containing less than
3 elements. In addition, for this example, the step-size $t=1$ was always
accepted. Algorithm 1 was started with the initial choice $A_{0}=\Omega$. As
can be seen from Table 1, the optimal values of $J$ and
$\|\chi\|_{L^{1}(\Omega)}$ converge for $h\searrow 0$, and
$\|\rho\|_{L^{1}(\Omega)}\to 0$ for $h\searrow 0$. According to Theorem 6.12,
this strongly suggests that the iterates are a minimizing sequence of (1.6).

$h$ | $J$ | $\|\chi\|_{L^{1}(\Omega)}$ | $\|\rho\|_{L^{1}(\Omega)}$
---|---|---|---
$4.42\cdot 10^{-2}$ | $4.712$ | $0.43896$ | $4.33\cdot 10^{-3}$
$2.21\cdot 10^{-2}$ | $5.054$ | $0.44299$ | $2.12\cdot 10^{-8}$
$1.13\cdot 10^{-2}$ | $5.216$ | $0.44352$ | $2.09\cdot 10^{-8}$
$5.66\cdot 10^{-3}$ | $5.299$ | $0.44432$ | $2.04\cdot 10^{-8}$
$2.83\cdot 10^{-3}$ | $5.340$ | $0.44455$ | $2.11\cdot 10^{-11}$
$1.41\cdot 10^{-3}$ | $5.360$ | $0.44460$ | $4.05\cdot 10^{-11}$
Table 1: Results of optimization, Section 7.1.1

Let us report about the influence of the choice of the initial guess
$A_{0}\subseteq\Omega$. Here we chose the following set of parameters: $y_{d}$
was as above, and

$\alpha=0.001,\ \beta=0.1,\ u_{a}=-40,\ u_{b}=+40.$

For this example, the method returned the same solution independent of the
initial guess. We depicted the iteration history for different choices of
$A_{0}$ in Fig. 2. In general, the method was faster when starting from
$A_{0}=\Omega$ than from $A_{0}=\emptyset$. As one can see from Fig. 2, the
convergence of $\|\rho_{k}\|_{L^{1}(\Omega)}$ is stable with respect to mesh
refinement.

Figure 2: Comparison of iteration history of $\|\rho_{k}\|_{L^{1}(\Omega)}$
for different choice of $A_{0}$: $A_{0}=\emptyset$ (left), $A_{0}=\Omega$
(right), Section 7.1.1

####  7.1.2 An unsolvable problem

Let us report about observations when applying our algorithm to an unsolvable
problem. This problem is taken from [39, Section 4.5]. It is very similar to
the above problem. The partial differential equation is chosen in such a way
that a constant control $u$ leads to a constant solution $y$. The problem
reads as follows: Minimize the functional

$\frac{1}{2}\|y-y_{d}\|_{L^{2}(\Omega)}^{2}+\frac{\alpha}{2}\|u\|_{L^{2}(\Omega)}^{2}+\beta\|u\|_{0},$

where $y$ denotes the weak solution of the elliptic partial differential
equation with Neumann boundary conditions

$-\Delta y+y=u\quad\text{ in }\Omega,\quad\frac{\partial y}{\partial
n}=0\quad\text{ on }\partial\Omega.$

Let $\alpha>0$ and $\beta>0$ be given, and set

$y_{d}(x)=-\sqrt{\frac{\beta}{\alpha}}-\sqrt{2\alpha\beta}.$

As argued in [39, Section 4.5] this optimal control problem is unsolvable.
This implies that there is also no minimizer $A$ of the value function $J(A)$.

$h$ | $J$ | $\|\chi_{A}\|_{L^{1}(\Omega)}$ | $\|\rho\|_{L^{1}(\Omega)}$
---|---|---|---
$4.42\cdot 10^{-2}$ | $10.014$ | $0.70703$ | $6.15\cdot 10^{-9}$
$2.21\cdot 10^{-2}$ | $10.014$ | $0.70703$ | $8.63\cdot 10^{-10}$
$1.13\cdot 10^{-2}$ | $10.014$ | $0.70701$ | $9.29\cdot 10^{-10}$
$5.66\cdot 10^{-3}$ | $10.014$ | $0.70670$ | $2.75\cdot 10^{-9}$
$2.83\cdot 10^{-3}$ | $10.014$ | $0.70674$ | $2.49\cdot 10^{-9}$
$1.41\cdot 10^{-3}$ | $10.014$ | $0.70693$ | $1.37\cdot 10^{-9}$
Table 2: Results of optimization, Section 7.1.2

In the computations, we used

$\beta=0.01,\quad\alpha=1000,\ \Omega=(0,1)^{2}.$

The results are shown in Table 2. There, the values of $J$,
$\|\chi_{A}\|_{L^{1}(\Omega)}$, and $\|\rho\|_{L^{1}(\Omega)}$ are shown for
the final iterate on meshes with different mesh-sizes $h$. In contrast to the
results of [39], where resulting controls had always support equal to
$\Omega$, the measure of the sets $A_{k}$ was in the order of $0.7$ after a
few steps of Algorithm 1. In our experiments, the line-search took much more
steps than in the previous example, and only small modifications of the
$A_{k}$ were accepted, resulting in very slow convergence. While
$\|\rho_{k}\|_{L^{1}(\Omega)}$ was very small after a few steps of the
algorithm, the support of $\rho_{k}$ never gets as small as for the previous
example. For this example, we stopped the algorithm after 100 iterations.
Still, according to Theorem 6.12, the algorithm produces a minimizing sequence
for the unsolvable example, which is not the case for the thresholding method
of [39].

###  7.2 Binary control problems

Following the ideas of [2, 3, 22, 29], we will apply our algorithm to a binary
control problem, where controls can only take values in $\\{0,+1\\}$. We will
use a problem considered in [2, 3], which reads: Minimize

$\min\frac{1}{2}\|y-y_{d}\|_{L^{2}(\Omega)}^{2}+\nu\|u\|_{L^{1}(\Omega)}$

over all $(y,u)\in H^{1}_{0}(\Omega)\times L^{2}(\Omega)$ satisfying

$-\Delta y=u\quad\text{ a.e.\@ on }\Omega$

and

$u(x)\in\\{0,1\\}\text{ f.a.a.\@ }x\in\Omega.$

Hence, $u$ itself is a characteristic function. And the above problem can be
written in our setting as: Minimize

$J(A,u):=\frac{1}{2}\|y-y_{d}\|_{L^{2}(\Omega)}^{2}+\nu\int_{A}\,\text{{}d}x$

over all $(y,u)\in H^{1}_{0}(\Omega)\times L^{2}(\Omega)$ satisfying

$-\Delta y=\chi_{A}u\quad\text{ a.e.\@ on }\Omega$

and the (trivial) constraint

$u=1\text{ a.e.\@ on }\Omega.$

This setting does not directly fit into our framework. Still we can compute
the topological derivative as follows. The solution of $u\mapsto J(A,u)$ is
given by $u_{A}\equiv 1$, which greatly simplifies the computations of Section
3. And we have the following result concering the topological derivative of
the value function.

######  Theorem 7.1.

The topological derivative $DJ(B)(x)$ of the value function of the binary
control problem exists for almost all $x\in\Omega$, and is given by

$DJ(B)(x)=\sigma(B,x)(\beta(x)+p_{B}(x))$

with $\sigma(B,x)$ as in Theorem 4.2.

######  Proof 7.2.

The result of Lemma 3.1 in this situation has to be modified to

$J(A,u_{A})-J(B,u_{B})+\frac{1}{2}\|y_{B}-y_{A}\|_{Y}^{2}=\int_{\Omega}(\chi_{A}-\chi_{B})(\beta+p_{A})\,\text{{}d}x.$

where we have used set $g=0$ and
$\chi_{A}u_{A}-\chi_{B}u_{B}=\chi_{A}-\chi_{B}$ in (3.1). Since
$p_{A}-p_{B}=S^{*}S(\chi_{A}-\chi_{B})$, we have the estimate
$\|p_{A}-p_{B}\|_{L^{q}(\Omega)}\leq
c\|\chi_{A}-\chi_{B}\|_{L^{2}(\Omega)}=c\|\chi_{A}-\chi_{B}\|_{L^{1}(\Omega)}^{1/2}$,
which replaces the result of Theorem 3.5. Now the claim can be proven as in
the proof of Theorem 4.2.

The topological derivative coincides with the result [2, Corollary 3.2]. The
computation of the topological derivative does not involve the solution of any
optimization problem: given $A$, only $y_{A}$ and $p_{A}$ have to be computed.

Let us report about the results for the following choice of parameters,
corresponding to Case 3 in [2, Section 9]:

$y_{d}=0.05,\quad\nu=0.002.$

The computed control on the finest discretization can be seen in Fig. 3, which
agrees with [2, Figure 4]. The results of the optimization runs for different
discretizations can be seen in Table 3. In all cases, the algorithm stopped
due to a failed line-search. Nevertheless, the error quantity
$\|\rho\|_{L^{1}(\Omega)}$ is very small, and is decreasing with decreasing
mesh-size. According to Theorem 6.12 this indicates that the algorithm
produces a minimizing sequence.

$h$ | $J$ | $\|\chi\|_{L^{1}(\Omega)}$ | $\|\rho\|_{L^{1}(\Omega)}$
---|---|---|---
$6.99\cdot 10^{-2}$ | $1.799\cdot 10^{-3}$ | $1.63770$ | $9.07\cdot 10^{-8}$
$3.49\cdot 10^{-2}$ | $1.872\cdot 10^{-3}$ | $1.63794$ | $9.15\cdot 10^{-8}$
$1.79\cdot 10^{-2}$ | $1.909\cdot 10^{-3}$ | $1.63808$ | $3.97\cdot 10^{-9}$
$8.94\cdot 10^{-3}$ | $1.928\cdot 10^{-3}$ | $1.63806$ | $2.12\cdot 10^{-11}$
$4.47\cdot 10^{-3}$ | $1.938\cdot 10^{-3}$ | $1.63802$ | $1.51\cdot 10^{-12}$
$2.24\cdot 10^{-3}$ | $1.943\cdot 10^{-3}$ | $1.63802$ | $3.85\cdot 10^{-14}$
Table 3: Results of optimization, Section 7.2 Figure 3: Solution for
$h=2.24\cdot 10^{-3}$, Section 7.2

We also implemented the level-set method from [3, 5]. All parameters were
chosen as in [3]. We took the same example as above and performed the
computations on the grid with $h=2.24\cdot 10^{-3}$. The resulting iteration
history of $\|\rho_{k}\|_{L^{1}(\Omega)}$ for our method and the level-set
method can be seen in Fig. 4. In the level-set method, the domain $A$ is
determined as $A:=\\{x:\ \psi(x)\geq 0\\}$, where $\psi$ is the level-set
function. Here, it makes a difference, whether $\psi$ is discretized by
piecewise constant (P0) or piecewise linear (P1) finite elements. We
implemented both variants. Our method was implemented using piecewise constant
functions for the control. The computations in all three methods were stopped
as soon as $\|\rho_{k}\|_{L^{1}(\Omega)}<10^{-11}$. The results can be seen in
Fig. 4: the black plus-signs correspond to the piecewise constant (P0) case,
while the blue x’s refer to the piecewise linear (P1) case of the level-set
method, where all integrals are computed exactly following the ideas of [25].
While the level-set method in the piecewise linear case seems to be the
fastest of these three methods, the comparison to our method with piecewise
constant discretization is a bit unfair, as the piecewise linear method can
resolve the interface much finer. Still our method is faster than the
piecewise constant variant of the level-set method.

Figure 4: Comparison of Algorithm 1 (red circles) and level-set method [3]
(black +’s: P0, blue x’s: P1): Iteration history of
$\|\rho_{k}\|_{L^{1}(\Omega)}$, Section 7.2

##  8 Conclusion

In this paper, we developed an algorithm to solve control problems with
$L^{0}$-cost. As showed in Theorem 6.12, the algorithm produces a minimizing
sequence, which is remarkable as the $L^{0}$-optimization problem is not
convex. The algorithm was motivated by the concept of topological derivatives
applied to a suitable chosen sub-problem.

## Acknowledgements

The author thanks two anonymous referees, Bastian Dittrich, Gerd Wachsmuth,
and Daniel Walther for remarks on an earlier version of the manuscript.

## References

  * [1] S.Amstutz, The topological asymptotic for the Navier-Stokes equations, _ESAIM Control Optim. Calc. Var._ 11 (2005), 401–425, doi:10.1051/cocv:2005012.
  * [2] S.Amstutz, A semismooth Newton method for topology optimization, _Nonlinear Anal._ 73 (2010), 1585–1595, doi:10.1016/j.na.2010.04.065.
  * [3] S.Amstutz, Analysis of a level set method for topology optimization, _Optim. Methods Softw._ 26 (2011), 555–573, doi:10.1080/10556788.2010.521557.
  * [4] S.Amstutz, An introduction to the topological derivative, _Engineering Computations_ (2021), https://hal.science/hal-03311517.
  * [5] S.Amstutz and H.Andrä, A new algorithm for topology optimization using a level-set method, _J. Comput. Phys._ 216 (2006), 573–588, doi:10.1016/j.jcp.2005.12.015.
  * [6] S.Amstutz and A.Laurain, A semismooth Newton method for a class of semilinear optimal control problems with box and volume constraints, _Comput. Optim. Appl._ 56 (2013), 369–403, doi:10.1007/s10589-013-9555-6.
  * [7] A.Ben Abda, M.Hassine, M.Jaoua, and M.Masmoudi, Topological sensitivity analysis for the location of small cavities in Stokes flow, _SIAM J. Control Optim._ 48 (2009/10), 2871–2900, doi:10.1137/070704332.
  * [8] V. I.Bogachev, _Measure theory. Vol. I_ , Springer-Verlag, Berlin, 2007, doi:10.1007/978-3-540-34514-5.
  * [9] H.Brézis and W. A.Strauss, Semi-linear second-order elliptic equations in $L^{1}$, _J. Math. Soc. Japan_ 25 (1973), 565–590, doi:10.2969/jmsj/02540565.
  * [10] A.Carpio and M. L.Rapún, Hybrid topological derivative and gradient-based methods for electrical impedance tomography, _Inverse Problems_ 28 (2012), 095010, 22, doi:10.1088/0266-5611/28/9/095010.
  * [11] E.Casas and D.Wachsmuth, A Note on Existence of Solutions to Control Problems of Semilinear Partial Differential Equations, _SIAM J. Control Optim._ 61 (2023), 1095–1112, doi:10.1137/22m1486418.
  * [12] J.Céa, S.Garreau, P.Guillaume, and M.Masmoudi, The shape and topological optimizations connection, volume 188, 2000, 713–726, doi:10.1016/s0045-7825(99)00357-6.  IV WCCM, Part II (Buenos Aires, 1998).
  * [13] J.Céa, A.Gioan, and J.Michel, Quelques résultats sur l’identification de domaines, _Calcolo_ 10 (1973), 207–232, doi:10.1007/bf02575843.
  * [14] J.Céa, A.Gioan, and J.Michel, Adaptation de la méthode du gradient à un problème d’identification de domaine, in _Computing methods in applied sciences and engineering (Proc. Internat. Sympos., Versailles, 1973), Part 2_ , Lecture Notes in Comput. Sci., Vol. 11, Springer, Berlin, 1974, 391–402.
  * [15] J.Diestel and J. J.Uhl, Jr., _Vector measures_ , Mathematical Surveys, No. 15, American Mathematical Society, Providence, R.I., 1977.  With a foreword by B. J. Pettis.
  * [16] I.Ekeland, On the variational principle, _J. Math. Anal. Appl._ 47 (1974), 324–353, doi:10.1016/0022-247x(74)90025-0.
  * [17] I.Ekeland and R.Témam, _Convex analysis and variational problems_ , volume 28 of Classics in Applied Mathematics, SIAM, Philadelphia, PA, 1999, doi:10.1137/1.9781611971088.
  * [18] L. C.Evans, _Partial differential equations_ , volume 19 of Graduate Studies in Mathematics, American Mathematical Society, Providence, RI, second edition, 2010, doi:10.1090/gsm/019.
  * [19] P.Fulmański, A.Laurain, J. F.Scheid, and J.Sokołowski, A level set method in shape and topology optimization for variational inequalities, _Int. J. Appl. Math. Comput. Sci._ 17 (2007), 413–430, doi:10.2478/v10006-007-0034-z.
  * [20] S.Garreau, P.Guillaume, and M.Masmoudi, The topological asymptotic for PDE systems: the elasticity case, _SIAM J. Control Optim._ 39 (2001), 1756–1778, doi:10.1137/s0363012900369538.
  * [21] P.Guillaume and K.Sid Idris, The topological asymptotic expansion for the Dirichlet problem, _SIAM J. Control Optim._ 41 (2002), 1042–1072, doi:10.1137/s0363012901384193.
  * [22] M.Hahn, S.Leyffer, and S.Sager, Binary optimal control by trust-region steepest descent, _Math. Program._ 197 (2023), 147–190, doi:10.1007/s10107-021-01733-z.
  * [23] M.Hassine, S.Jan, and M.Masmoudi, From differential calculus to 0-1 topological optimization, _SIAM J. Control Optim._ 45 (2007), 1965–1987, doi:10.1137/050631720.
  * [24] M.Hintermüller and A.Laurain, A shape and topology optimization technique for solving a class of linear complementarity problems in function space, _Comput. Optim. Appl._ 46 (2010), 535–569, doi:10.1007/s10589-008-9201-x.
  * [25] M.Hinze, A variational discretization concept in control constrained optimization: the linear-quadratic case, _Comput. Optim. Appl._ 30 (2005), 45–61, doi:10.1007/s10589-005-4559-5.
  * [26] M.Hrizi and M.Hassine, One-iteration reconstruction algorithm for geometric inverse source problem, _J. Elliptic Parabol. Equ._ 4 (2018), 177–205, doi:10.1007/s41808-018-0015-4.
  * [27] K.Ito and K.Kunisch, Optimal control with $L^{p}(\Omega)$, $p\in[0,1)$, control cost, _SIAM J. Control Optim._ 52 (2014), 1251–1275, doi:10.1137/120896529.
  * [28] D.Kalise, K.Kunisch, and K.Sturm, Optimal actuator design based on shape calculus, _Math. Models Methods Appl. Sci._ 28 (2018), 2667–2717, doi:10.1142/s0218202518500586.
  * [29] P.Manns, M.Hahn, C.Kirches, S.Leyffer, and S.Sager, On convergence of binary trust-region steepest descent, _J. Nonsmooth Anal. Optim._ 4 (2023), doi:10.46298/jnsao-2023-10164.
  * [30] A.Münch, Optimal design of the support of the control for the 2-D wave equation: a numerical method, _Int. J. Numer. Anal. Model._ 5 (2008), 331–351.
  * [31] A.Münch, Optimal location of the support of the control for the 1-D wave equation: numerical investigations, _Comput. Optim. Appl._ 42 (2009), 443–470, doi:10.1007/s10589-007-9133-x.
  * [32] A. A.Novotny and J.Sokołowski, _Topological derivatives in shape optimization_ , Interaction of Mechanics and Mathematics, Springer, Heidelberg, 2013, doi:10.1007/978-3-642-35245-4.
  * [33] A. A.Novotny, J.Sokołowski, and A.Żochowski, _Applications of the topological derivative method_ , volume 188 of Studies in Systems, Decision and Control, Springer, Cham, 2019, doi:10.1007/978-3-030-05432-8.  With a foreword by Michel Delfour.
  * [34] A. A.Novotny, J.Sokołowski, and A.Żochowski, Topological derivatives of shape functionals. Part II: First-order method and applications, _J. Optim. Theory Appl._ 180 (2019), 683–710, doi:10.1007/s10957-018-1419-x.
  * [35] R. T.Rockafellar, Integral functionals, normal integrands and measurable selections, in _Nonlinear operators and the calculus of variations (Summer School, Univ. Libre Bruxelles, Brussels, 1975)_ , Lecture Notes in Math., Vol. 543, Springer, Berlin, 1976, 157–207.
  * [36] J.Sokołowski and A.Żochowski, On the topological derivative in shape optimization, _SIAM J. Control Optim._ 37 (1999), 1251–1272, doi:10.1137/s0363012997323230.
  * [37] J.Sokołowski and A.Żochowski, Topological derivative for optimal control problems, _Control Cybernet._ 28 (1999), 611–625.  Recent advances in control of PDEs.
  * [38] G.Stampacchia, Le problème de Dirichlet pour les équations elliptiques du second ordre à coefficients discontinus, _Ann. Inst. Fourier (Grenoble)_ 15 (1965), 189–258, http://www.numdam.org/item?id=AIF_1965__15_1_189_0.
  * [39] D.Wachsmuth, Iterative hard-thresholding applied to optimal control problems with $L^{0}(\Omega)$ control cost, _SIAM J. Control Optim._ 57 (2019), 854–879, doi:10.1137/18m1194602.