Optimal Control Algorithm for Complex Heat Transfer Model Alexander Chebotarev, Gleb Grenkin, and Andrey Kovtanyuk ? 1 Far Eastern Federal University, Sukhanova st. 8, 690950 Vladivostok, Russia, 2 Institute for Applied Mathematics, Radio st. 7, 690041 Vladivostok, Russia {cheb@iam.dvo.ru,glebgrenkin@gmail.com,kovtanyuk.ae@dvfu.ru} Abstract. An optimal control problem for a nonlinear steady-state heat trans- fer model accounting for heat radiation effects is considered. The problem con- sists in minimization of a given cost functional by controlling the sources in the heat equation. The solvability of this control problem is proved, optimality con- ditions are derived, and an iterative algorithm for solving the optimal control problem is constructed. Keywords: optimal control; radiative heat transfer; conductive heat transfer; gradient descent method 1 Introduction The interest in studying problems of complex heat transfer (where the radiative, convec- tive, and conductive contributions are simultaneously taken into account) is motivated by their importance for many engineering applications. Here, the following examples can be mentioned: modeling and predicting the heat transfer in molten glass [1–3], nanofluids [4, 5], etc. A considerable number of works devoted to optimal control problems of complex heat transfer models consider the evolutionary systems (see, e.g., [1–3, 6–10]). In the mentioned works, the radiation transfer is described by steady-state radiative transfer equation. The temperature field is simulated by the conventional evolutionary heat transfer equation with additional source terms describing the contribution of the ra- diative heat transfer. Theoretical analysis of optimal control problems for steady-state systems of com- plex heat transfer with source terms in the heat equation is an open question. It is worth to mention the work [11], where the problem of optimal boundary multiplicative control for a steady-state complex heat transfer model was considered. The problem was formulated as the maximization of the energy outflow from the model domain by ? The research was supported by the Ministry of education and science of Russian Federation (project 14.Y26.31.0003) Copyright c by the paper’s authors. Copying permitted for private and academic purposes. In: A. Kononov et al. (eds.): DOOR 2016, Vladivostok, Russia, published at http://ceur-ws.org 166 A. Chebotarev, G. Grenkin, A. Kovtanyuk controlling reflection properties of the boundary. On the basis of new a priori estimates of solutions of the control system, the solvability of the optimal control problem was proved. The main result there was the proof of an analogue of the bang-bang principle arising in control theory for ordinary differential equations. In this paper, an optimal control problem of obtaining a desired temperature and(or) radiative intensity distributions in a part of the model domain by control- ling the sources in the heat equation is considered. Analogous problems appear in many engineering applications and draw attention of many researchers. For example, similar optimal control problems for non-stationary complex heat transfer models were studied in [1–3, 8] in context of glass manufacturing. In the current work, the optimal control problem for a steady-state model is studied. The solvability of this problem is proved, and an optimality system is derived. Moreover, an iterative algorithm based on the gradient descent method is constructed, and results of numerical experiments are presented. 2 Problem formulation The following steady-state normalized diffusion (P1 ) model (see [12–15]) describing radiative and conductive heat transfer in a bounded domain Ω ⊂ R3 is under consid- eration: −a∆θ + bκa (|θ|θ3 − ϕ) = u, −α∆ϕ + κa (ϕ − |θ|θ3 ) = 0, (1) a∂n θ + β(θ − θb )|Γ = 0, α∂n ϕ + γ(ϕ − |θb |θb3 )|Γ = 0. (2) Here, θ is the normalized temperature, ϕ the normalized radiation intensity averaged over all directions, and κa the absorption coefficient. The physical sense of the param- eters a, b, α, β, γ can be found in [13–15]. The control function u describes the internal sources of heat. The symbol ∂n denotes the derivative in the outward normal direction n on the boundary Γ := ∂Ω. The problem of optimal control consists in the determination of functions u, θ, and ϕ which satisfy the conditions (1), (2) and minimize a cost functional Jµ (θ, ϕ, u), i.e. µ Jµ (θ, ϕ, u) = J(θ, ϕ) + kuk2L2 (Ω) → inf, u ∈ Uad . (3) 2 Here, Uad ⊂ L2 (Ω) is the set of admissible controls, µ ≥ 0 is a given cost parameter. In particular, the functional J can describe the L2 -deviation of the temperature and radiation fields from prescribed distributions, say θd and ϕd . Thus, e.g. J(θ, ϕ) = aθ kθ − θd k2L2 (Ω) + aϕ kϕ − ϕd k2L2 (Ω) , where aθ and aϕ are nonnegative weights. 3 Formalization of the optimal control problem Suppose that the model data satisfy the following conditions: (i) β, γ ∈ L∞ (Γ ), β ≥ β0 > 0, γ ≥ γ0 > 0, β0 , γ0 = const, θb ∈ L∞ (Γ ); Optimal Control Algorithm 167 (ii) Uad is a closed convex set; Uad is a bounded set, if µ = 0; (iii) The cost functional J : H 1 (Ω) × H 1 (Ω) → R is weakly lower semicontinuous and bounded from below. Here and further, the Sobolev space W2s (Ω) is denoted by H s (Ω), s ≥ 0, and (f, g) and kf k denote respectively the inner product and the norm of the space L2 (Ω). Denote H = L2 (Ω), V = H 1 (Ω), Y = V × V . Identifying H with the dual space H yields the Gelfand triple V ⊂ H = H 0 ⊂ V 0 . Let the value of a functional f ∈ V 0 0 on an element v ∈ V be denoted by (f, v). Notice that (f, v) is the inner product in H if f and v are elements of H. Assuming that θ, ϕ, v are arbitrary elements of V , define operators and functionals A1 , A2 : V → V 0 , f, g ∈ V 0 by the following relations: Z Z (A1 θ, v) = a(∇θ, ∇v) + βθvdΓ, (A2 ϕ, v) = α(∇ϕ, ∇v) + γϕvdΓ, Γ Γ Z Z (f, v) = βθb vdΓ, (g, v) = γ|θb |θb3 vdΓ. Γ Γ A pair {θ, ϕ} ∈ V is called weak solution of the problem (1), (2) if A1 θ + bκa (|θ|θ3 − ϕ) = f + u, A2 ϕ + κa (ϕ − |θ|θ3 ) = g. (4) The optimal control problem consists in the minimization of a functional Jµ defined on solutions of system (4) provided that u ∈ Uad . That is, Jµ (θ, ϕ, u) → inf, {θ, ϕ} are solutions of (4) yielded by u ∈ Uad . (5) A pair {θ, b ϕ} b minimizing Jµ and corresponding to a function u b is called optimal state, and u b is called optimal control. 4 Solvability of the optimal control problem To prove the solvability of the problem (5), establish some properties of the boundary value problem (1), (2). Lemma 1. If the conditions (i) hold and u ∈ Uad , then for a weak solution, {θ, ϕ}, of the problem (1), (2) the following estimate is true: kθk2V + kϕk2V ≤ C. (6) Here, a positive constant C depends only on a, b, α, κa , β, γ, kuk, and Ω. Proof. Let hp (s) := |s|p signs, p > 0, s ∈ R. Denote ϕ1 = h1/4 (ϕ) and, for ε > 0, define  ϕ1 − ε,  ϕ1 > ε, wε = 0, |ϕ1 | ≤ ε,  ϕ1 + ε, ϕ1 < −ε.  168 A. Chebotarev, G. Grenkin, A. Kovtanyuk Notice that if ϕ ∈ V , then ϕ1 ∈ L24 (Ω), ϕ1 |Γ ∈ L16 (Γ ), wε ∈ V , and ( 1 |ϕ|−3/4 ∇ϕ, |ϕ1 | > ε, ∇wε = 4 0, otherwise. It is important that Z γϕwε dΓ − κa (h4 (θ) − ϕ, wε ) − (g, wε ) = Γ Z = γ(ϕ − h4 (θb ))ϕ1 dΓ − κa (h4 (θ) − ϕ, ϕ1 ) + cε . (7) Γ In this expression |cε | ≤ Cε, where C > 0 does not depend on ε. Multiply, in the sense of the inner product of H, the first equation of (4) by θ, the second equation by bwε , and add the equalities. Then, taking into account monotonicity of (h4 (θ) − ϕ)(θ − h1/4 (ϕ)) ≥ 0, we obtain the inequality Z Z Z 16 ak∇θk2 + βθ2 dΓ + αb |∇ψ|2 dx + b γψ 2 dΓ Γ 25 |ψ|>ε5/2 Γ Z ≤ (f + u, θ) + b γh4 (θb )ϕ1 dΓ − bcε . (8) Γ Here, ψ = h5/8 (ϕ), ϕ1 = h2/5 (ψ). Passing to the limit in inequality (8) as ε → +0, we obtain ψ ∈ V and Z k1 kθk2V + k2 kψk2V ≤ |(f + u, θ)| + b γ|h4 (θb )h2/5 (ψ)|dΓ. (9) Γ 16 Here, k1 = min{a, β0 }, k2 = b min{ 25 α, γ0 }. The norm in the space V is defined by the following equality: Z kvk2V = k∇vk2 + v 2 dΓ. Γ Taking into account the continuity of the trace operator from V into L2 (Γ ), we obtain from (9):   kθk2V + kψk2V ≤ K1 kf + uk2V 0 + kθb k4L5 (Γ ) . (10) Here, K1 depends only on a, α, b, β0 , γ0 , kγkL∞ (Γ ) , and the domain Ω. The estimate of kθkV allows to obtain the estimate of kϕkV . Multiplying the second equation of (4) by ϕ in the sense of the inner product of H, and denoting k3 = min{α, γ0 }, we obtain the inequality Z k3 kϕk2V + κa kϕk2 ≤ κa |(h4 (θ), ϕ)| + γ|h4 (θb )ϕ|dΓ. Γ Using Hölder and Young inequalities with parameter δ > 0, we estimate: δ 1 |(h4 (θ), ϕ)| ≤ kϕk2L3 (Ω) + kθk8L6 (Ω) , 2 2δ Optimal Control Algorithm 169 Z   δ 1 γ|h4 (θb )ϕ|dΓ ≤ kγkL∞ (Γ ) kϕk2L4 (Γ ) + kθb k8L16/3 (Γ ) . Γ 2 2δ Taking into account the continuity of the embedding of V into L6 (Ω), the continuity of the trace operator from V into L4 (Γ ), and a sufficiently small δ, we obtain the estimate of kϕkV :   kϕk2V ≤ K2 kθb k8L16/3 (Γ ) + kθk8V . (11) Here, K2 depends only on α, γ0 , κa kγkL∞ (Γ ) , and the domain Ω. The estimates (10) and (11) prove the lemma. t u On the base of the estimate (6), similarly as in [11], the solvability of the problem (5) is proved. Theorem 1. If the conditions (i)–(iii) hold, then there exists at least one solution of the problem (5). 5 The necessary conditions of optimality To derive optimality relations, add to conditions (i)-(iii) the following assumption: (iv) J : Y → R is Fréchet differentiable. Introduce a constraint operator F : Y × H → Y 0 as follows: F (y, u) = {A1 θ + bκa (|θ|θ3 − ϕ) − f − u, A2 ϕ + κa (ϕ − |θ|θ3 ) − g}, where y = {θ, ϕ} ∈ Y, u ∈ H. Lemma 2. For any y ∈ Y the map Fy0 : Y → Y 0 is epimorphic, Im Fy0 = Y 0 . Proof. Equation Fy0 q = z = {z1 , z2 } ∈ Y 0 is equivalent to the following boundary value problem: A1 q1 + bκa (4|θ|3 q1 − q2 ) = z1 , A2 q2 + κa (q2 − 4|θ|3 q1 ) = z2 , q = {q1 , q2 } ∈ Y. (12) To prove the solvability of a Fredholm problem (12), it suffices to prove the uniqueness of its solutions. Let sign s = s/|s|, if s 6= 0, sign 0 = [−1, 1]. Let us consider the function µδ which is regularization of multivalued function sign, µδ (s) = s/|s|, if |s| ≥ δ, µδ (s) = s/δ, if |s| < δ. Let z = 0, h = 4|θ|3 . Multiplying, in the sense of the inner product of H, the first equation of (12) by µδ (q1 ), the second equation by bµδ (q2 ), and adding these equalities, we obtain (A1 q1 , µδ (q1 )) + b(A2 q2 , µδ (q2 )) + bκa (hq1 − q2 , µδ (q1 ) − µδ (q2 )) = 0. Notice that Z Z (A1 q1 , µδ (q1 )) = a(∇q1 , µ0δ (q1 )∇q1 ) + βq1 µδ (q1 )dΓ ≥ βq1 µδ (q1 )dΓ. Γ Γ 170 A. Chebotarev, G. Grenkin, A. Kovtanyuk and Z Z (A2 q2 , µδ (q2 )) = α(∇q2 , µ0δ (q2 )∇q2 ) + γq2 µδ (q2 )dΓ ≥ γq2 µδ (q2 )dΓ. Γ Γ Therefore, Z Z βq1 µδ (q1 )dΓ + b γq2 µδ (q2 )dΓ + bκa (hq1 − q2 , µδ (q1 ) − µδ (q2 )) ≤ 0. (13) Γ Γ Passing to the limit as δ → 0, from inequality (13), we obtain Z Z β|q1 |dΓ + b γ|q2 |dΓ + bκa (hq1 − q2 , sign q1 − sign q2 ) ≤ 0. Γ Γ From monotonicity of the function sign, it follows the conditions q1 |Γ = q2 |Γ = 0. Further, notice that A1 q1 + bA2 q2 = 0. Scalarly multiplying this equation by aq1 + αbq2 , and taking into account zero boundary values of q1 , q2 , we obtain k∇(aq1 + αbq2 )k2 = 0. Hence, aq1 + αbq2 = 0. Therefore, a(∇q1 , ∇v) + bκa ((h + a/αb)q1 , v) = 0 ∀v ∈ V. (14) Assuming v = q1 in (14), we obtain q1 = 0, and therefore q2 = 0. t u Applying the principle of Lagrange for smooth convex extremal problems [16], we can prove the following result. Theorem 2. Let yb = {θ, b ∈ Y, u b ϕ} b ∈ Uad be a solution of the control problem (5). Then there exists an adjoint state p = {p1 , p2 } ∈ Y such that the triple {b y, u b, p} satisfies the conditions b 3 κa (bp1 − p2 ) = −J 0 (b 0 A1 p1 + 4|θ| θ y ), A2 p2 + κa (p2 − bp1 ) = −Jϕ (b y ), (15) u − p1 , v − u (µb b) ≥ 0, ∀v ∈ Uad . (16) 6 Example of optimality system Let G1,2 ⊂ Ω be subdomains of Ω. Consider the following cost functional: Z Z 1 1 J(θ, ϕ) = (θ − θd )2 dx + (ϕ − ϕd )2 dx, (17) 2 2 G1 G2 where θd ∈ L2 (G1 ) and ϕd ∈ L2 (G2 ) are given functions. It is easy to see that Z Z (Jθ0 (θ, ϕ), η) = (θ − θd )ηdx, (Jϕ0 (θ, ϕ), η) = (ϕ − ϕd )ηdx, η ∈ V. G1 G2 Optimal Control Algorithm 171 In this case, if Uad = L2 (Ω) and µ > 0, the optimality system assumes the form b θb3 − ϕ) − a∆θb + bκa (|θ| b =u b, −α∆ϕ b + κa (ϕ b θb3 ) = 0, b − |θ| a∂n θb + β(θb − θb )|Γ = 0, α∂n ϕ b − |θb |θb3 )|Γ = 0, (18) b + γ(ϕ b 3 (bp1 − p2 ) = −χG (θb − θd ), − a∆p1 + 4κa |θ| 1 − α∆p2 + κa (p2 − bp1 ) = −χG2 (ϕ b − ϕd ), a∂n p1 + βp1 |Γ = 0, α∂n p2 + γp2 |Γ = 0, (19) and u b = p1 /µ. Here, χG1,2 are the characteristic functions of the subdomains G1,2 , respectively. 7 Iterative algorithm For the numerical solution of the optimality system (18), (19), we can apply the method of gradient descent:   (k) uk+1 = uk − λk µuk − p1 , k = 0, 1, 2, . . . , where u0 ∈ H is given. (k) (k) Here λk > 0 is a step size, p(k) = {p1 , p2 } is a pair satisfying the system (18), (19), where ub := uk . If Uad 6= H we can apply the gradient projection method (see, e.g, [17]):    (k) uk+1 = PUad uk − λk µuk − p1 , k = 0, 1, 2, . . . , where PUad : H → Uad is the projection operator. The method of choosing the step size λk is adjusted as required for decreasing the cost functional. Unlike methods based on the Armijo rule, this method does not need an inner loop for adjustment of λk that requires computing the value of J and, therefore, solving the problem (18). Define J(u) b = Jµ (θ(u), ϕ(u), u), where {θ(u), ϕ(u)} is a solution of system (4). The b k+1 ) ≥ J(u method is as follows. If J(u b k ), then return back to the control uk and reduce λk by a factor of 2. Additionally, if J(u b k ), k = s, s + 1, . . . , s + m0 − 1, b k+1 ) < J(u then λk is increased by a factor of 2. Here, m0 ≥ 1 is a prescribed integer parameter of a quantity of decreases of the cost functional, which is enough for increasing the step size λk . The pseudocode of the algorithm is presented below. 172 A. Chebotarev, G. Grenkin, A. Kovtanyuk Algorithm 1: Gradient descent method with a variable step size Choose the parameters λ0 and m0 . Choose the initial guess u0 . cost func decreases ← 0; for k ← 0, 1, 2, . . . do For the given uk , find {θk , ϕk } from (18). Compute J(ub k ). if k ≥ 1 and J(u b k ) ≥ J(u b k−1 ) then uk ← uk−1 ; λk ← λk−1 /2; p(k) ← p(k−1) ; cost func decreases ← 0; else if k ≥ 1 then cost func decreases ← cost func decreases + 1; Find p(k) from (19). if cost func decreases = m0 then λk ← 2λk−1 ; cost func decreases ← 0; else if k ≥ 1 then λk ← λk−1 ;    (k) uk+1 ← PUad uk − λk µuk − p1 ; 8 Numerical example Consider an example for the two-dimensional domain Ω = {(x, y) : 0 ≤ x, y ≤ L} which can be interpreted as a long rectangular channel in the three-dimensional space. The parameters values are taken as follow: L = 10 [cm], α = 3.3 . . . [cm], κa = 0.01 [cm−1 ], β = 1.5 [cm/s], and γ = ε/2(2 − ε), where ε = 0.7 is the emissivity coefficient of the boundary. The thermodynamical characteristics of the medium correspond to air at the normal atmospheric pressure and the temperature of 400 ◦ C. The maximum temperature is chosen as Tmax = 773 K. This yields a = 0.92 [cm2 /s] and b = 18.7 [cm/s]. Notice that the absolute temperature is T = Tmax θ. The boundary temperature θb = 0.5. The cost functional is defined by (3) and (17), where G1 = Ω \ S, G2 = ∅, θd = 0.7, and µ = 0.01. Let the set of admissible controls be Uad = {u ∈ L2 (Ω) : u1 ≤ u ≤ u2 }, where u1 = 0, u2 = 1 in S = [x1 , x2 ] × [y1 , y2 ], and u1 = u2 = 0 in Ω \ S. Assume x1 = 0.65L, x2 = 0.85L, y1 = L/12, y2 = 5L/12. For the numerical solution we use the software FreeFem++ [18]. The boundary- value problem (18) is solved by Newton’s method. The initial guess for the optimal Optimal Control Algorithm 173 control is chosen as u0 = 0, and parameters of the optimization algorithm are λ0 = 5, m0 = 3. The computed optimal control is presented in Fig. 1. The graph of the optimal temperature is depicted in Fig. 2. The values of J(u b k ) and λk for different k are indicated in Figs. 3, 4. As it is seen in Fig. 4, the most frequent value of λk is 10, and the step size is adjusted as needed. The optimal controls for µ = 0.1 and µ = 0.001 are presented in Figs. 5, 6 for comparison. It can be easily proved from (16) that in the case of µ = 0 the optimal control satisfies an analog of the bang-bang principle, that is u b(x) = u1 (x) or u2 (x) for a.e. x ∈ Ω where p1 (x) 6= 0. Notice that the optimal control comes near to a bang-bang control as µ → 0. The optimal control for µ = 0 is depicted in Fig. 7. This bang-bang control was computed by an optimization algorithm of the gradient descent type, see [9, 10]. 4.0 u=1 3.5 3.0 1 0.9 0.7 0.8 u=0 2.5 0.3 0.3 0.5 2.0 0.5 0 0.6 0.4 1.5 1 0.7 0.5 0.6 1.0 0.9 0.3 0.4 6.5 7.0 7.5 8.0 8.5 Fig. 1. Optimal control (µ = 0.01) 174 A. Chebotarev, G. Grenkin, A. Kovtanyuk 10 0.51 0.51 8 0.55 6 0.7 0.53 0.6 4 1 0.9 2 1.1 0.8 0.51 0.55 0 0 2 4 6 8 10 Fig. 2. Optimal temperature (µ = 0.01) 0.9755 1.8 0.9754 0.9753 1.6 Cost functional 0.9752 0.9751 1.4 0.9750 0.9749 1.2 6 8 10 12 14 16 18 20 22 24 1.0 0 1 2 3 4 5 6 7 8 Iterations Fig. 3. Cost functional J(u b k ) at different iterations (µ = 0.01) Optimal Control Algorithm 175 40 30 Step size 20 10 0 0 5 10 15 20 25 Iterations Fig. 4. Step size λk at different iterations (µ = 0.01) 4.0 1 0.9 0.7 3.5 0.8 0.5 3.0 0.6 2.5 0.4 2.0 0.3 1.5 1.0 6.5 7.0 7.5 8.0 8.5 Fig. 5. Optimal control (µ = 0.1) 176 A. Chebotarev, G. Grenkin, A. Kovtanyuk 4.0 3.5 u=1 3.0 2.5 1 u=0 0 2.0 0 1.5 0.5 1 0.5 1 1.0 0 6.5 7.0 7.5 8.0 8.5 Fig. 6. Optimal control (µ = 0.001) 4.0 3.5 3.0 u=1 2.5 u=0 u=0 2.0 1.5 1.0 u=0 6.5 7.0 7.5 8.0 8.5 Fig. 7. Bang-bang optimal control (µ = 0) Optimal Control Algorithm 177 References 1. Pinnau, R., Thömmes, G.: Optimal boundary control of glass cooling processes. Math. Methods Appl. Sci. 27(11), 1261–1281 (2004) 2. Frank, M., Klar, A., Pinnau, R.: Optimal control of glass cooling using simplified PN theory, Transport Theory Statist. Phys. 39(2–4), 282–311 (2010) 3. Clever, D., Lang, J.: Optimal control of radiative heat transfer in glass cooling with re- strictions on the temperature gradient. Optimal Control Appl. Methods. 33(2), 157–175 (2012) 4. Mabood, F., Shateyi, S., Rashidi, M.M., Momoniat, E., Freidoonimehr, N.: MHD stag- nation point flow heat and mass transfer of nanofluids in porous medium with radiation, viscous dissipation and chemical reaction. Advanced Powder Technology 27, 742–749 (2016) 5. Rashidi, M.M., Ganesh, N.V., Hakeem, A.K.A., Ganga, B.: Buoyancy effect on MHD flow of nanofluid over a stretching sheet in the presence of thermal radiation. J. Molecular Liquids 198, 234-238 (2014) 6. Pinnau, R.: Analysis of optimal boundary control for radiative heat transfer modeled by the SP1 system. Commun. Math. Sci. 5(4), 951–969 (2007) 7. Tse, O., Pinnau, R.: Optimal control of a simplified natural convection-radiation model. Commun. Math. Sci. 11(3), 679–707 (2013) 8. Clever, D., Lang, J., Schröder, D.: Model hierarchy-based optimal control of radiative heat transfer. Int. J. Comput. Sci. Eng. 9(5–6), 509–525 (2014) 9. Grenkin, G.V., Chebotarev, A.Yu., Kovtanyuk, A.E., Botkin, N.D., Hoffmann, K.-H.: Boundary optimal control problem of complex heat transfer model. J. Math. Anal. Appl. 433(2), 1243–1260 (2016) 10. Grenkin, G.V.: An algorithm for solving the problem of boundary optimal control in a complex heat transfer model. Dal’nevost. Mat. Zh. 16(1), 24–38 (2016) 11. Kovtanyuk, A.E., Chebotarev, A.Yu., Botkin, N.D., Hoffmann, K.-H.: Theoretical analysis of an optimal control problem of conductive-convective-radiative heat transfer. J. Math. Anal. Appl. 412(1), 520–528 (2014) 12. Modest, M.F.: Radiative Heat Transfer, Academic Press (2003) 13. Kovtanyuk, A.E., Chebotarev, A.Yu., Botkin, N.D., Hoffmann, K.-H.: The unique solv- ability of a complex 3D heat transfer problem. J. Math. Anal. Appl. 409(2), 808–815 (2014) 14. Kovtanyuk, A.E., Chebotarev, A.Yu.: Steady-state problem of complex heat transfer. Comput. Math. Math. Phys. 54(4), 719–726 (2014) 15. Kovtanyuk, A.E., Chebotarev, A.Yu., Botkin, N.D., Hoffmann, K.-H.: Unique solvability of a steady-state complex heat transfer model. Commun. Nonlinear Sci. Numer. Simul. 20(3), 776–784 (2015) 16. Ioffe, A.D., Tikhomirov, V.M.: Theory of Extremal Problems. North-Holland, Amsterdam (1979) 17. Hinze, M., Pinnau, R., Ulbrich, M., Ulbrich, S.: Optimization with PDE Constraints, Springer (2009) 18. Hecht, F.: New development in FreeFem++. J. Numer. Math. 20(3–4), 251–266 (2012)