 # Two-sided bounds for cost functionals of time-periodic parabolic optimal control problems

In this paper, a new technique is applied on deriving computable, guaranteed lower bounds of functional type (minorants) for two different cost functionals subject to a parabolic time-periodic boundary value problem. Together with previous results on upper bounds (majorants) for one of the cost functionals, both minorants and majorants lead to two-sided estimates of functional type for the optimal control problem. Moreover, a different cost functional is discussed with the same PDE-constraints. Both upper and lower bounds are derived. The time-periodic optimal control problems are discretized by the multiharmonic finite element method leading to large systems of linear equations having a saddle point structure. The derivation of preconditioners for the minimal residual method for the new optimal control problem is discussed in more detail. Finally, several numerical experiments for both optimal control problems are presented confirming the theoretical results obtained.

## Authors

##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1. Introduction

The optimal control of partial differential equations is the subject matter of lots of different works (see, e.g.,

[39, 16, 48, 6] and the references therein) due to their importance in research and application. Optimal control problems with time-periodic parabolic state equations often arise in different practical applications such as in electromagnetics, chemistry, biology, or heat transfer. We refer as examples to the recent works  and  considering problems in electromagnetics and biochemistry, respectively. The multiharmonic finite element method (MhFEM) – also called harmonic-balanced finite element method – is a natural approach for the solution of parabolic time-periodic problems, where the state functions are expanded into truncated Fourier series in time. The space-dependent Fourier coefficients are then approximated by the finite element method (FEM). The MhFEM was successfully used for the simulation of electromagnetic devices described by nonlinear eddy current problems with harmonic excitations, see, e.g., [50, 3, 4, 5, 11]. Later, the MhFEM has been applied to linear time-periodic parabolic boundary value and optimal control problems [19, 20, 25, 31] and to linear time-periodic eddy current problems and the corresponding optimal control problems [21, 22]. Recent works on robust preconditioners for time-periodic parabolic and eddy-current optimal control problems are discussed in  and , respectively.

A posteriori estimates of functional type provide a useful machinery to derive computable and guaranteed quantities for the desired unknown solution. These estimation techniques were earlier introduced by S. Repin, see, e.g., the papers on parabolic problems [42, 14]. Recent works on new estimates for parabolic problems and parabolic optimal control problems can be found in  and , respectively. A posteriori estimates of functional type for elliptic optimal control problems can be found in [12, 13, 43, 34]. First functional type estimates for inverse problems can be found in [44, 10]. Moreover, recent results on guaranteed computable estimates for convection-dominated diffusion problems are presented in .

In , majorants (fully computable upper bounds) for the optimality system and the cost functional of a time-periodic parabolic optimal control were presented. In this work, we want to extend the analysis by deriving the corresponding minorants (fully computable lower bounds) for this cost functional using the technique presented in  applying also earlier results derived by Mikhlin  (see also  for a more recent discussion). We mention here that Mali  presents a different approach for the derivation of a lower bound for a class of elliptic optimal control problems. Moreover, we extend the analysis and consider a second cost functional with respect to the same parabolic time-periodic boundary value problem. The results on computable lower bounds together with the upper bounds lead to two-sided estimates which can be used to derive majorants for the discretization error in state and control. These majorants and minorants provide a new formulation of the optimization problems since they can, in principle, be used as objects of direct minimization on their difference. In this work, robust preconditioners for the preconditioned minimal residual (MinRes) method (see ) are discussed for the second minimization problem, which are new for this case.

The paper is organized as follows: In Section 2, we discuss both time-periodic parabolic optimal control problems including some preliminary results. The majorants and new minorants for minimization problem I are presented in Section 3 followed by the results for minimization problem II in Section 4. The multiharmonic finite element discretization of the optimization problems is considered in Section 5. Section 6 presents robust preconditioners for applying the preconditioned MinRes method on the problems discretized by the MhFEM. Finally, Section 7 discusses detailed a set of various numerical experiments for both problems and presents their results.

## 2. Two Cost Functionals with Respect to Time-Periodic Parabolic PDE-Constraints

Let denote the space-time cylinder and its lateral surface, where is a given time interval. The spatial domain , , is assumed to be a bounded Lipschitz domain with boundary . The minimization problems are both subject to the following parabolic time-periodic boundary value problem:

 (1) ⎧⎪⎨⎪⎩σ(x)∂ty(x,t)−div(ν(x)∇y(x,t))=u(x,t)(x,t)∈QT,y(x,t)=0(x,t)∈ΣT,y(x,0)=y(x,T)x∈¯¯¯¯Ω,

where and denote the state and control, respectively. The coefficients and are uniformly bounded and satisfy the conditions

 (2) 0<σ––≤σ(x)≤¯¯¯σ,and0<ν––≤ν(x)≤¯¯¯ν,x∈Ω,

where , , and are positive constants. In the following, we present a proper functional space setting for time-periodic problems. The notation is similar to that was used in [27, 28]. We define the Hilbert spaces

 H1,0(QT) ={v∈L2(QT):∇v∈[L2(QT)]d}, H0,1(QT) ={v∈L2(QT):∂tv∈L2(QT)}, H1,1(QT) ={v∈L2(QT):∇v∈[L2(QT)]d,∂tv∈L2(QT)},

which are endowed with the norms

 ∥v∥1,0 :=(∫QT(v(x,t)2+|∇v(x,t)|2)dxdt)1/2, ∥v∥0,1 :=(∫QT(v(x,t)2+|∂tv(x,t)|2)dxdt)1/2, ∥v∥1,1 :=(∫QT(v(x,t)2+|∇v(x,t)|2+|∂tv(x,t)|2)dxdt)1/2,

respectively. Here, is the spatial gradient and denotes the weak derivative with respect to time. Also we define subspaces of the above introduced spaces by putting subindex zero if the functions satisfy the homogeneous Dirichlet condition on and subindex per if they satisfy the time-periodical condition . All inner products and norms in related to the whole space-time domain are denoted by and , respectively. If they are associated with the spatial domain , then we write and . The symbols and denote the standard inner products and norms of the space .

The functions used in our analysis are presented in terms of Fourier series, e.g., the Fourier representation of the function is given by the series

 (3) v(x,t)=vc0(x)+∞∑k=1(vck(x)cos(kωt)+vsk(x)sin(kωt)),

where

are the Fourier coefficients, denotes the period, and is the frequency. The choice of this type of representation is natural because the problem has time-periodical conditions. We also introduce the spaces

 H0,12per(QT) ={v∈L2(QT):∥∥∂1/2tv∥∥<∞}, H1,12per(QT) ={v∈H1,0(QT):∥∥∂1/2tv∥∥<∞}, H1,120,per(QT) ={v∈H1,12per(QT):v=0 on ΣT},

where is defined in the Fourier space by the relation

 (4) ∥∥∂1/2tv∥∥2:=|v|20,12:=T2∞∑k=1kω∥vk∥2Ω.

Here, for all , see also . These spaces can be considered as Hilbert spaces if we introduce the following (equivalent) products:

 ⟨∂1/2ty,∂1/2tv⟩:=T2∞∑k=1kω⟨yk,vk⟩Ω,⟨σ∂1/2ty,∂1/2tv⟩:=T2∞∑k=1kω⟨σyk,vk⟩Ω.

The above introduced spaces allow us to operate with “symmetrized” formulations of the time-periodic problems. The seminorm and norm in are defined by the relations

 |v|21,12 :=∥∇v∥2+∥∂1/2tv∥2=T∥∇vc0∥2Ω+T2∞∑k=1(kω∥vk∥2Ω+∥∇vk∥2Ω) and ∥v∥21,12 :=∥v∥2+|v|21,12=T(∥vc0∥2Ω+∥∇vc0∥2Ω)+T2∞∑k=1((1+kω)∥vk∥2Ω+∥∇vk∥2Ω),

respectively. Using Fourier type series, it is easy to define the function “orthogonal” to :

 v⊥(x,t):=∞∑k=1(−vck(x)sin(kωt)+vsk(x)cos(kωt))=∞∑k=1(vsk(x),−vck(x))=:(−v⊥k)T⋅(cos(kωt)sin(kωt)).

Obviously, and we find that

 ∥∥∂1/2tv⊥∥∥2=T2∞∑k=1kω∥v⊥k∥2Ω=T2∞∑k=1kω∥vk∥2Ω=∥∥∂1/2tv∥∥2∀v∈H0,12per(QT).

The identities

 (5) ⟨σ∂1/2ty,∂1/2tv⟩=⟨σ∂ty,v⊥⟩ and ⟨σ∂1/2ty,∂1/2tv⊥⟩=⟨σ∂ty,v⟩

are valid for all and , see also . Also, we define the product

 (6) ⟨ξ,∂1/2tv⟩:=T2∞∑k=1(kω)1/2⟨ξk,vk⟩Ω.

We recall the orthogonality relations

 (7)

and the identity (in the sense of Fourier series)

 (8) ∫QTξ∂1/2tv⊥dxdt=−∫QT∂1/2tξ⊥vdxdt∀ξ,v∈H0,12per(QT),

see . We note that for functions presented in terms of Fourier series the standard Friedrichs inequality holds in :

 (9) ∥∇u∥2=T∥∇uc0∥2Ω+T2∞∑k=1∥∇uk∥2Ω≥1C2F(T∥uc0∥2Ω+T2∞∑k=1∥uk∥2Ω)=1C2F∥u∥2.

In the following, the parameter denotes the regularization or cost parameter.

### 2.1. Minimization problem I

In the first case, we want to minimize the following cost functional with respect to the unknown state and control :

 (10) J(y,u):=12∥y−yd∥2+λ2∥u∥2

subject to the time-periodic boundary value problem (1), where is the given desired state. The cost functional defined in (10) can be written as

 (11) J(y,u)=TJ0(yc0,uc0)+T2∞∑k=1Jk(yk,uk),

where and . In , the corresponding optimality system is derived, which is given in weak formulation as follows: Given , find such that

 (12) ⎧⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪⎩∫QT(yz−ν(x)∇p⋅∇z+σ(x)∂1/2tp∂1/2tz⊥)dxdt=∫QTydzdxdt,∫QT(ν(x)∇y⋅∇q+σ(x)∂1/2ty∂1/2tq⊥+λ−1pq)dxdt=0,

for all test functions . Using the Fourier series ansatz (3) in (12) and exploiting the orthogonality of and , we arrive at the following problem: Find such that

 (13)

for all test functions . The system (13) must be solved for every mode . For , we obtain a reduced problem: Find such that

 (14) ⎧⎪ ⎪ ⎪⎨⎪ ⎪ ⎪⎩∫Ω(yc0⋅zc0−ν(x)∇pc0⋅∇zc0)dx=∫Ωycd0⋅zc0dx,∫Ω(ν(x)∇yc0⋅∇qc0+λ−1pc0⋅qc0)dx=0,

for all test functions . The problems (13) and (14) have unique solutions, see . In , the a posteriori error analysis for the optimality system as well as majorants for the cost functional are presented. This work extends and deepens these results by deriving minorants for this cost functional as well as for a different one, which is introduced in the next subsection.

### 2.2. Minimization problem II

In the second case, we want to minimize the following cost functional with respect to the unknown state and control :

 (15) ~J(y,u):=12∥∇y−gd∥2+λ2∥u∥2

subject to the time-periodic boundary value problem (1), where is the given desired gradient. The optimality system can analogously be derived as for minimization problem I using the Lagrange functional

 ~L(y,u,p)=~J(y,u)−∫QT(σ(x)∂ty−div(ν(x)∇y)−u)pdxdt.

Only the equation in the optimality conditions is different for problem II. The optimality conditions are given in weak form as follows: Given , find such that

 (16) ⎧⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪⎩∫QT(∇y⋅∇z−ν(x)∇p⋅∇z+σ(x)∂1/2tp∂1/2tz⊥)dxdt=∫QTgd⋅∇zdxdt,∫QT(ν(x)∇y⋅∇q+σ(x)∂1/2ty∂1/2tq⊥+λ−1pq)dxdt=0,

for all test functions . The optimality systems corresponding to every mode are analogously derived as for minimization problem I (similar to (13) and (14)). In Section 4, we will derive the two-sided bounds for minimization problem II.

## 3. Two-sided Bounds for Minimization Problem I

### 3.1. Majorant for cost functional (10)

First, the results of  on upper bounds for minimization problem I are presented, since they are needed later to derive the two-sided estimate. Let be the corresponding state to an arbitrary control . The following upper bound can be proved:

 J(y(v),v)≤J⊕(α,β;η,τ,v)∀v∈L2(QT),

for arbitrary , and

 τ∈H(div,QT):={τ∈[L2(QT)]d:divxτ(⋅,t)∈L2(Ω) for a.e. t∈(0,T)},

where, for any , the identity

 (17) ∫Ωdivτwdx=−∫Ωτ⋅∇wdx∀w∈H10(Ω)

is valid. The guaranteed and fully computable majorant is given by

 (18) ⎧⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎩J⊕(α,β;η,τ,v):=1+α2∥η−yd∥2+(1+α)(1+β)C2F2αμ1–––2∥R2(η,τ)∥2+(1+α)(1+β)C4F2αβμ1–––2∥R1(η,τ,v)∥2+λ2∥v∥2,

where and is the constant coming from the Friedrichs inequality. The parameters have been introduced in order to obtain a quadratic functional by applying Young’s inequality. The arbitrary functions and can be taken as the approximate solutions of the optimal control problem (10) subject to (1) and represents the image of the exact flux . For the derivation of (18), the following estimate for the approximation error has been used:

 (19)

where

 (20) R1(η,τ,v):=σ∂tη−divτ−v and R2(η,τ):=τ−ν∇η.

The derivation of estimate (19) can be found in . The function is a sharp upper bound on for arbitrary but fixed as well as on the optimal value , i.e.,

 (21) infη∈H1,10,per(QT),τ∈H(div,QT)v∈L2(QT),α,β>0J⊕(α,β;η,τ,v)=J(y(u),u),

since the infimum is attained for the optimal control , its corresponding state and its exact flux , and for going to zero. Therefore, we have the estimate

 (22) J(y(u),u)≤J⊕(α,β;η,τ,v)∀η∈H1,10,per(QT),τ∈H(div,QT),v∈L2(QT),α,β>0.

### 3.2. Minorant for cost functional (10)

In this work, we complement the guaranteed upper bounds for the discretization error in state and control of minimizing cost functional defined in (10) subject to (1). This is done by obtaining fully computable lower bounds for in the spirit of  leading to two-sided bounds for the cost functional (10). Let be the optimal state corresponding to the optimal control function , which is connected with the adjoint state by the identity . Then is the solution of the variational formulation

 (23) ∫QT(ν∇y⋅∇q+σ∂1/2ty∂1/2tq⊥+λ−1pq)dxdt=0∀q∈H1,120,per(QT)

of the boundary value problem (1) (see also (12)). For any , we have that

 J(y(u),u)=12∥y−η∥2+∫QT(y−η)(η−yd)dxdt+12∥η−yd∥2+λ2∥u∥2.

Since and using the identity , we can estimate from below by

 (24) J(y(u),u)=J(y(u),p(u))≥12∥η−yd∥2+12λ∥p∥2+∫QT(y−η)(η−yd)dxdt.

For , let be the solutions to the equations

 (25) ∫QT(ν∇pη⋅∇z−σ∂1/2tpη∂1/2tz⊥)dxdt=∫QT(η−yd)zdxdt, (26) ∫QT(ν∇η⋅∇q+σ∂1/2tη∂1/2tq⊥+λ−1~pηq)dxdt=0,

for all test functions , respectively.

###### Remark 1.

Note that we assumed that according to the derivation of the majorant, but so far the assumption would be enough.

Adding and subtracting in (24) together with yields the estimate

 J(y(u),u) =J(y(u),p(u)) ≥12∥η−yd∥2+12λ∥pη∥2+∫QT(y−η)(η−yd)dxdt+∫QTλ−1(p−pη)pηdxdt.

By using (25) and identity (8), it follows that

 (27)

Using the equations (23) and (26) leads to the estimate

 (28) J(y(u),u)≥ 12∥η−yd∥2+12λ∥pη∥2+∫QTλ−1(~pη−pη)pηdxdt.

We introduce now the arbitrary function

. Note that at the moment

would be enough, but the higher regularity in time will be needed in another step. This goes along with the higher regularity assumption on (see Remark 1). Since , we have that

 J(y(u),u)≥ 12∥η−yd∥2+12λ∥ζ∥2+∫QTλ−1(pηζ−ζ2+~pηpη−p2η)dxdt.

Now we add and subtract in the last integral as well as use equation (26) again. Moreover, we exploit the fact that we have assumed that , hence, we can apply the identities (5). Altogether these steps yield the estimate

 (29) ⎧⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪⎩J(y(u),u)≥12∥η−yd∥2+12λ∥ζ∥2−∫QT(ν∇η⋅∇ζ+σ∂tηζ+λ−1ζ2)dxdt+∫QTλ−1(ζ−pη)(pη−~pη)dxdt.

In the following, we need to estimate the last integral of this expression in order to formulate a computable lower bound for the cost functional. For that let us first prove a computable upper bound for the error in the adjoint state, which is presented in the following theorem. Note that here we will need the higher regularity assumption (in time) on .

###### Theorem 1.

Let be given and let meet equation (25) with satisfying assumptions (2). For any , we have that

 (30) ∥∇(pη−ζ)∥≤1μ1–––(CF∥R3(ζ,ρ,η)∥+∥R4(ζ,ρ)∥),

where , and with and is the constant coming from the Friedrichs inequality.

###### Proof.

Let us consider the adjoint equation (25). Adding and subtracting in the left-hand side of the equation leads to

 (31) ⎧⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪⎩∫QT(ν∇(pη−ζ)⋅∇z−σ∂1/2t(pη−ζ)∂1/2tz⊥)dxdt=∫QT(η−yd)zdxdt−∫QTν∇ζ⋅∇zdxdt+∫QTσ∂1/2tζ∂1/2tz⊥dxdt.

Next we introduce the auxiliary variable . Together with using that as well as applying the Cauchy-Schwarz and Friedrichs inequalities, the following estimate for the right-hand side of (31) can be obtained:

 sup0≠z∈H1,120(QT)∫QT(η−yd+divρ+σ∂tζ)zdxdt+∫QT(ρ−ν∇ζ)⋅∇zdxdt|z|1,12 ≤sup0≠z∈H1,120(QT)∥η−yd+divρ+σ∂tζ∥∥z∥+∥ρ−ν∇ζ∥∥∇z∥|z|1,12 ≤sup0≠z∈H1,120(QT)CF∥η−yd+divρ+σ∂tζ∥∥∇z∥+∥ρ−ν∇ζ∥∥∇z∥∥∇z∥ ≤CF∥η−yd+divρ+σ∂tζ∥+∥ρ−ν∇ζ∥.

Using (2) and the orthogonality relations (7), we can prove the following estimate from below for the left-hand side of (31):

 sup0≠z∈H1,120(QT)∫QT(ν∇(pη−ζ)⋅∇z−σ∂1/2t(pη−ζ)∂1/2tz⊥)dxdt|z|1,12≥μ1–––|pη−ζ|1,12.

Combining now both estimates together with we finally derive estimate (30). ∎

Now we have all the tools in order to estimate the last term of (29) as follows

 (32) ⎧⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎩∫QTλ−1(ζ−pη)(pη−ζ+ζ−~pη)dxdt=∫QT(λ−1(ζ−pη)(pη−ζ)+λ−1(ζ−pη)(ζ−~pη))dxdt=−λ−1∥ζ−pη∥2+∫QT(λ−1(ζ−pη)(ζ−~pη))dxdt=−λ−1∥ζ−pη∥2+∫