Accepted Manuscript Heterogeneous multiscale method for optimal control problem governed by parabolic equation with highly oscillatory coefficients Liang Ge, Yanzhen Chang, Danping Yang
PII: DOI: Reference:
S0377-0427(17)30329-1 http://dx.doi.org/10.1016/j.cam.2017.06.029 CAM 11207
To appear in:
Journal of Computational and Applied Mathematics
Received date : 7 September 2016 Revised date : 14 April 2017 Please cite this article as: L. Ge, Y. Chang, D. Yang, Heterogeneous multiscale method for optimal control problem governed by parabolic equation with highly oscillatory coefficients, Journal of Computational and Applied Mathematics (2017), http://dx.doi.org/10.1016/j.cam.2017.06.029 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Heterogeneous Multiscale Method for Optimal Control Problem Governed by Parabolic Equation with Highly Oscillatory CoefficientsI Liang Gea , Yanzhen Changb,∗, Danping Yangc a Shandong
Provincial Key Laboratory of Computer Networks,Shandong Computer Science Center (National Supercomputer Center in Jinan), Jinan 250014, China b Department
c Department
of Mathematics, Beijing University of Chemical Technology, Beijing, 100029, China
of Mathematics, Shanghai Key Laboratory of Pure Mathematics and Mathematical Practice, East China Normal University, Shanghai, 200062, China
Abstract In this paper, we investigate the heterogeneous multiscale method (HMM) for the optimal control problem governed by the parabolic equation with highly oscillatory coefficients. The state variable and adjoint state variable are approximated by the multiscale discretization scheme that relies on coupled macro and micro finite elements, while the control variable is discretized by the piecewise constants. By applying the well-known Lions’ Lemma to the discretized optimal control problem, we obtain the necessary and sufficient optimality conditions. A priori error estimates are derived for the state, co-state and the control with uniform bounded constants. Finally, numerical results are presented to illustrate our theoretical findings. Keywords: optimal control problem, heterogeneous multiscale finite element, a priori error estimate 2010 MSC: 49J20, 65N30
1. Introduction In modern scientific and engineering computation, optimal control problems have been so frequently met in all kinds of practical problems. The research on developing fast numerical algorithms for optimal control problems is continuously and rapidly expanding, which is simply impossible to give even a very brief review 5
here. Moreover, because the practical procedures are very complicated, many problems show the multiscale character, such as composite materials with thermal/electrical conductivity, flow through the heterogeneous porous media, and time scale of the chemical reactions, etc. Therefore, optimal control problems with multiscale character become an especially important research branch. I The work is supported by the Natural Science Foundation of Shandong Province, China(Grant No. ZR2011FQ030, ZR2013FQ001, ZR2013FM025), Natural Science Foundation of China(Grant No. 11501326), and the Shandong Academy of Sciences Youth Fund Project(Grant No. 2013QN007), BUCT Fund for Disciplines Construction and Development(Grant NO. XK1523). ∗ Corresponding author Email address:
[email protected](Liang Ge),
[email protected](Yanzhen Chang),
[email protected](Danping Yang) (Danping Yang)
Preprint submitted to Journal of Computational and Applied Mathematics
July 10, 2017
As well known, the multiscale problems are often described by partial differential equations(PDEs) with 10
highly oscillatory coefficients [1]. Thus there come up with many kinds of multiscale computation methods to deal with the intrinsic complexity of the multiscale problems, such as multigrid method [2],[3], adaptive method [4], domain decomposition method [5], homogenization method [6], wavelet-based numerical homogenization method [7], partition of unity method [8], multiscale finite element method [9], heterogeneous multiscale method [10], variational multiscale method [11], etc.
15
Similarly, the multiscale optimal control problems can be described by optimal control problems governed by PDEs with highly oscillatory coefficients. Compared with the classical optimal control problems, the multiscale optimal control problems are more difficult due to the intrinsic complexity of the state equations. In [12], Lions studied the multiscale asymptotic methods which give the asymptotic expansions for the optimal control with small parameter ε. Cao gave the multiscale asymptotic expansions to the boundary control [13], optimal
20
control for elliptic systems with constraints [14] and optimal control for linear parabolic equations [15]. In [16, 17, 18, 19], the authors studied the homogenization methods which could give the overall behavior by incorporating the fluctuations due to the heterogeneities [14]. But the homogenized methods are specifically designed for periodic homogenization problems. From the point of computing, in contrast with high performance that multiscale finite element methods have shown on other multiscale problems, there exists little work about
25
multiscale optimal control problems. To our best knowledge, there are only some results for the mulitscale optimal control problems in [1], which gave the result for elliptic equation optimal control problem and obtained just the first order convergence for the control, state, and co-state. Besides, the computation cost of multiscale finite element methods in general is equivalent to solve the full fine scale problems, see [20]. Compared with the above methods, the heterogeneous multiscale method in [21, 22, 23, 10, 24] which has more widely applications
30
and higher computational efficiency, is more suitable for the constrained multiscale optimal control problems. Therefore, we investigate the heterogeneous multiscale methods for the constrained optimal control problems with highly oscillatory coefficients. As to the elliptic optimal control problem with highly oscillatory coefficients, the related results have been given in our another work [25]. However, it is much more complicated to implement computational schemes for evolutional control problems. So in this paper, the heterogeneous multiscale finite
35
element method is applied to solve the optimal control problem governed by parabolic equation with rapidly oscillating coefficients ε. We derive the continuous and discrete first-order optimality conditions by Lagrange multiplier method. Furthermore, we give several a priori error estimates in L2 (0, T ; L2 (Ω)) norm for the state variable, the adjoint state variable, and the control variable with uniform bounded constants as well as a priori error estimates in L∞ (0, T ; L2 (Ω)) norm for the state variable and the adjoint state variable with uniform
40
bounded constants. Numerical experiments are given to illustrate the validity of the estimates. Compared with the multiscale finite element methods for optimal control problem in [1], this work gives improved results in convergence order. The remainder of the paper is organized as follows. In the next section, we introduce the mathematical
2
model, and derive the optimal conditions for the parabolic optimal control problem with rapidly oscillating 45
coefficients. In section 3, the heterogeneous multiscale finite element scheme for the optimal control problem is given and the discretized optimal conditions are obtained. In Section 4, we prove a priori error estimates for the state variable, the co-state variable and the control variable, respectively. In Section 5, we carry out the numerical experiments to confirm our theoretical findings. Finally, we draw some concluding remarks in Section 6.
50
In what follows, C > 0 or c > 0 denotes a generic constant, independent of ε, h and hU . We denote Ω ⊂ Rd (d = 2, 3) a bounded domain with Lipschitz boundary ∂Ω and ΩU ⊂ Rd another bounded domain with Lipschitz boundary ∂ΩU . Generally, ΩU can be a subset of Ω. In the special case, ΩU = Ω. We adopt the standard notation W m,q (Ω) for Sobolev space on Ω with a norm ∥ · ∥m,q,Ω and a seminorm | · |m,q,Ω . We set
W0m,q (Ω) = {v ∈ W m,q (Ω) : v|∂Ω = 0}. We denote H m (Ω) = W m,2 (Ω) and ∥ · ∥m,Ω = ∥ · ∥m,2,Ω . We denote 55
by Ls (0, T ; W m,q (Ω)) the Banach space of all Ls integrable functions from (0, T ) into W m,q (Ω) with norm ∫T 1 ∥v∥Ls (0,T ;W m,q (Ω)) = ( 0 ∥v∥sW m,q (Ω) dt) s for s ∈ [1, ∞) and the standard modification for s = ∞. Similarly,
we define the spaces H 1 (0, T ; W m,q (Ω)) and C l (0, T ; W m,q (Ω)).
2. Optimal control problem governed by parabolic equation with oscillatory coefficients In this paper, we are interested in the following constrained optimal control problem: ∫ T ε min J(u ) = min { g(y ε (uε )) + j(uε )dt}, ε ε u ∈K
60
u ∈K
(2.1)
0
subject to the parabolic equation: d ∑ ∂y ε ∂ ∂y ε − (aε (x, t) ) = ∂t ∂xj ∂xi i,j=1
yε
=
y(x, 0) =
f + Buε ,
x ∈ Ω, t ∈ (0, T ],
0, on ∂Ω, t ∈ (0, T ], y0 (x),
x ∈ Ω,
where B is a linear operator from L2 (0, T ; L2 (ΩU )) to L2 (0, T ; H −1 (Ω)). Here uε (x) is a control function, g(·) and j(·) are two convex differentiable functions on L2 (0, T ; H01 (Ω)) and L2 (0, T ; L2 (ΩU )) respectively. We set f ∈ L2 (0, T ; H −1 (Ω)), y0 ∈ H01 (Ω) and ε > 0 is a small periodic parameter which represents the relative size
of a periodic cell. Here, K is a closed convex subset in the control space L2 (0, T ; L2 (ΩU )). In the following
65
context, we will give the choice of K. Similarly to [26], we make the following assumptions. (A1 ) : Let ξ =
x ε
and aε (x, t) = aε (x, xε , t) = a(x, ξ, t) is smooth and periodic for the variable ξ with period
Y = (− 21 , 12 )d . (A2 ) : aε (x, t) = (aij (x, xε , t))d×d is a symmetric positive definite matrix, which satisfies the uniform condi3
70
tion i.e., λ|ξ|2 ≤ aε (x, t)ξ · ξ ≤ Λ|ξ|2 , ∀ ξ ∈ Rd , where λ, Λ are positive constants independent of ε. (A3 ) : Let the Gˆateaux derivative g ′ (·) of the functional g(·) be Lipschitz continuous, i.e., ∥g ′ (u) − g ′ (v)∥ ≤ M ∥u − v∥L2 (Ω) , ∀ u, v ∈ L2 (0, T ; L2 (ΩU )), and (j ′ (u) − j ′ (v), u − v)U ≥ γ∥u − v∥2L2 (ΩU ) , ∀ u, v ∈ L2 (0, T ; L2 (ΩU )), where M, γ are positive constants independent of ε.
75
(A4 ) : The operator B is bounded. To consider the heterogeneous multiscale approximation of the above optimal control problem, we firstly derive a weak formulation for the state equation. Then the weak formulation of the state equation is to find ∩ y ε ∈ L2 (0, T ; H01 (Ω)) H 1 (0, T ; L2 (Ω)) such that (∂t y ε , w) + A(y ε , w) ∫ ∫ ∂y ε = wdx + aε (x, t)∇y ε · ∇wdx ∂t Ω ∫Ω ε = (f + Bu )wdx = (f + Buε , w), ∀ w ∈ L2 (0, T ; H01 (Ω)).
(2.2)
Ω
Therefore the variational formulation corresponding to (2.1) can be rewritten as the following, which we 80
shall label (OCP):
(OCP ) :
where y ε ∈ L2 (0, T ; H01 (Ω))
∩
∫T min J(uε ) = min { 0 g(y ε ) + j(uε )dt}, ε ∈K ε ∈K u u (∂t y ε , w) + A(y ε , w) = (f + Buε , w), ∀ w ∈ H01 (Ω), y ε (0) = y (x), 0
H 1 (0, T ; L2 (Ω)), uε ∈ L2 (0, T ; L2 (ΩU ))
∩
(2.3)
K.
According to Lions’ theorem ([27]), there exists a unique minimizer, which satisfies the following variational inequality: J ′ (u)(w − u) ≥ 0, ∀ w ∈ K. Here, the directional derivative of functional J at u ∈ K along the direction w ∈ K is defined by J ′ (u)(w) = lim+ s→0
J(u + sw) − J(u) . s
4
85
From the above assumptions we can derive (see e.g.[27], or [14] Lemma 1.1) that the control problem (2.3) has a unique solution (y ε , uε ), and that a pair (y ε , uε ) is the solution of (2.3) if and only if there is a co∩ state pε ∈ L2 (0, T ; H01 (Ω)) H 1 (0, T ; L2 (Ω)), such that (y ε , pε , uε ) satisfies the following optimality conditions, which we shall label (OCP-OPT): (∂t y ε , w) + A(y ε , w) = (f + Buε , w), ∀ w ∈ L2 (0, T ; H01 (Ω)), y ε (0) = y0 , (OCP − OP T ) : −(∂t pε , q) + A(q, pε ) = (g ′ (y ε ), q), ∀ q ∈ L2 (0, T ; H01 (Ω)), pε (T ) = 0, ∫T (j ′ (uε ) + B ∗ pε , v − uε )U dt ≥ 0, ∀ v ∈ K, 0
(2.4)
where B ∗ is the adjoint operator of B. 90
The explicit solution of the variational inequality of (2.4) depends heavily on the choice of the convex set ∫T ∫ ∫ K. For example, if J(uε ) = 0 ( 12 Ω (y ε − yd )2 dx + α2 Ω (uε )2 dx)dt, we can have the following explicit solutions for some cases (see, e.g. [27, 28]). Now, let K be given by
K = {uε ∈ L2 (0, T ; L2 (Ω)) : uε ≥ 0, a.e. in Ω},
(2.5)
1 uε (x, t) = max{0, − pε (x, t)}, a.e. x ∈ Ω. α
(2.6)
then the solution is
3. Heterogeneous multiscale method for optimal control problem
95
In this section we will give the heterogeneous multiscale finite element approximation of (2.1). The numerical method is based on a macro and a micro FEM. It has been firstly introduced by E and Engquist [10] and developed by [21, 22, 23, 24, 29]. At first, we describe the spatial discretization in Section 3.1 based on the HMM framework, using macro and micro FEM. Then, we discuss the time-discretization in Section 3.2. 3.1. The semi-discrete approximation scheme
100
Let Ωh be a polygonal approximation to Ω with a boundary ∂Ωh . For simplicity, we assume that Ωh = Ω ¯h = ∪ in this paper. Let T h be a partitioning of Ωh into disjoint regular n-simplices τ , so that Ω ¯ . Each τ ∈T h τ element has at most one face on ∂Ωh , and τ¯ and τ¯′ have either only one common vertex or a whole edge or face if τ and τ ′ ∈ T h . We further require that Pi ∈ ∂Ωh ⇒ Pi ∈ ∂Ω, where {Pi }(i = 1 . . . J) is the vertex set associated with the triangulation T h .
105
¯ h ), such that χ|τ are polynomials of m-order Associated with T h is a finite dimensional subspace S h of C(Ω (m ≥ 1) ∀χ ∈ S h and τ ∈ T h . Let Y h = S h ∩ H01 (Ω). It is easy to see that Y h ⊂ V = H01 (Ω).
¯h = ∪ Let TUh be a partitioning of ΩhU into disjoint regular n-simplices τU , so that Ω ¯U . Assume that U τU ∈T h τ U
ΩhU = ΩU , and τ¯U and τ¯U′ have either only one common vertex or a whole face or are disjoint if τU and τU′ ∈ TUh . 5
Figure 1: Illustration of HMM for solving y ε and pε . The dots are the quadrature points.
Associated with TUh is another finite dimensional subspace U h of L2 (ΩhU ), such that χ|τU are polynomials 110
of m-order (m ≥ 0) ∀χ ∈ U h and τU ∈ TUh . Here there is no requirement for the continuity. It is easy to see
that U h ⊂ L2 (ΩU ).
In this paper, we will only consider the simplest finite element spaces, i.e., m = 1 for Y h and m = 0 for U h . Let hτ (hτU ) denote the maximum diameter of the element τ (τU ) in T h (TUh ), let h = maxτ ∈T h {hτ }, hU = maxτU ∈TUh {hτU }. Moreover, let K h = U h ∩ K. Then we have that K h ⊂ K.
115
Next, we consider the HMM scheme for the parabolic partial differential equation (see, e.g., [29]): d ∑ ∂y ∂ x ∂y − (aij (x, , t) ) = f, ∂t i,j=1 ∂xj ε ∂xi
in Ω,
(3.1)
y = 0, on ∂Ω.
Let vh ∈ Y h be the solution of the problem (
∂vh , wh ) + Ah (vh , wh ) = (f, wh ), ∂t
∀ wh ∈ Y h ,
(3.2)
and the bilinear form Ah (v, w) for any v, w ∈ Y h can be written as ∫ ∑ ∑ ∫ |τ |∇w · A(xτ , t)∇v, Ah (v, w) = ∇w · A(x, t)∇vdx = ∇w · A(x, t)∇vdx ≃ Ω
τ ∈T h
τ
τ ∈T h
where xτ is the barycenter of τ . The approximation of A(xτ , t) can be given by solving the following CauchyDirichlet problem:
∂v ε ∂t
−
∑d
∂ x ∂v ε i,j=1 ∂xj (aij (x, ε , t) ∂xi )
v ε = v, v ε (t) = v(t),
= 0,
in (xτ + Iδ ) × (t, t + τt ), on ∂Iδ × (t, t + τt ), in xτ + Iδ .
6
(3.3)
120
where τt denotes the microsimulation time that evolves in time t as well as Iδ = δ · Y with the unit cell Y := (− 21 , 12 )d and small size δ. Let Qt = (xτ + Iδ ) × (t, t + τt ), we have Ah (v, w) =
d ∑ |τ | ∫ ∑ x ∂v ε ∂wε aij (x, , t) dxdt. |Qt | Qt i,j=1 ε ∂xi ∂xj h
τ ∈T
Then, HMM shceme of problem (3.1) is to find yh ∈ Y h such that (∂t yh , wh ) + Ah (yh , wh ) = (f, wh ), ∀ wh ∈ Y h . Therefore, a possible semi-discrete finite element heterogeneous multiscale approximation of optimal control
125
problem (2.1) is as follows, which we shall label (OCP )h : ∫T min Jh (uh ) = min { 0 g(yh ) + j(uh )dt}, h h u ∈K u ∈K h h h h h (OCP ) : ( ∂y ∂t , wh ) + Ah (yh , wh ) = (f + Buh , wh ), ∀ wh ∈ Y , yh (x, 0) = yh0 (x),
(3.4)
where Ah (·, ·) is defined as above and yh0 (x) ∈ Y h is an approximation of y 0 (x).
Next, the following theorem gives the existence and uniqueness of the solution for the above system. Theorem 3.1. Assume that the condition (A1 ) − (A4 ) hold, then there exists an unique solution (uh , yh ) for
the minimization problem (OCP )h such that uh ∈ L2 (0, T ; U h ), yh ∈ L2 (0, T ; Y h ).
h ∞ Proof. Let {(uhn , yhn )}∞ n=1 be a minimization sequence for the system (OCP ) , then it is clear that {uhn }n=1 130
∞ are bounded in L2 (0, T ; U h ). Thus there is a subsequence of {uhn }∞ n=1 (still denote by {uhn }n=1 ) such that
uhn converges to uh weakly in L2 (0, T ; U h ). For the subsequence {uhn }∞ n=1 , we have (
∂yhn , wh ) + Ah (yhn , wh ) = (f + Buhn , wh ), ∀ wh ∈ Y h . ∂t
(3.5)
Noticing that if taking wh = yhn , we get (
∂yhn , yhn ) + Ah (yhn , yhn ) = (f + Buhn , yhn ), ∂t
i.e. 1 d ∥yhn ∥2L2 (Ω) + Ah (yhn , yhn ) = (f + Buhn , yhn ). 2 dt
(3.6)
From [24], it is obvious that Ah (yhn , yhn ) ≥ c∥yhn ∥2H 1 (Ω) . Integrating (3.6) from 0 to t, we obtain ∥yhn (t)∥2L2 (Ω) + c 135
∫
t 0
∥yhn ∥2H 1 (Ω) dτ ≤ ∥yh0 ∥2L2 (Ω) +
∫
t
∥f ∥2H −1 (Ω) + ∥uhn ∥2L2 (ΩU ) dτ.
0
Applying Gronwall’s inequality to (3.7) yields ∥yhn ∥2L∞ (0,T ;L2 (Ω))
+
∥yhn ∥2L2 (0,T ;H 1 (Ω))
≤
C(∥yh0 ∥2L2 (Ω) 7
+
∫
0
T
(∥f ∥2H −1 (Ω) + ∥uhn ∥2L2 (ΩU ) )dτ ).
(3.7)
∩ 2 ∞ 2 Thus, we have that {yhn }∞ L (0, T ; Y h ), and n=1 is a bounded set in L (0, T ; L (Ω)) u −→ u , weakly in L2 (0, T ; U h ), hn h y −→ y , weakly in L2 (0, T ; Y h ). hn
h
So, we can derive
(
∂yh , wh ) + Ah (yh , wh ) = (f + Buh , wh ), ∀ wh ∈ Y h . ∂t
Since g(·) is a convex functional on space L2 (Ω) and j(·) is a strictly convex functional on L2 (ΩU ), we know ∫
0
T
∫ g(yh ) + j(uh )dt ≤ lim{
g(yhn ) + j(uhn )dt}.
0
This implies (yh , uh ) is the unique solution of (OCP )h . 140
T
Using the Lagrange multiplier method, we drive the optimality conditions of the problem (OCP )h in the following theorem. Theorem 3.2. Assume that the condition (A1 ) − (A4 ) hold, we have that a pair (yh , uh ) is the solution of
(OCP )h , if and only if there is a co-state ph ∈ L2 (0, T ; Y h ), such that (yh , ph , uh ) satisfies the following
145
optimality conditions, which we shall label (OCP − OP T )h : h h ( ∂y ∂t , wh ) + Ah (yh , wh ) = (f + Buh , wh ), ∀ wh ∈ Y , y (x, 0) = y 0 (x), h h h ′ h h (OCP − OP T ) : −( ∂p ∂t , qh ) + Ah (qh , ph ) = (g (yh ), qh ), ∀ qh ∈ Y , ph (x, T ) = 0, ∫ T (j ′ (u ) + B ∗ p , v − u ) dt ≥ 0, ∀ v ∈ K h . h h h h U h 0
(3.8)
Proof. By the same argument of [27], the optimal condition reads
Jh′ (uh )(vh − uh ) ≥ 0, ∀ vh ∈ K h . At first, from the definition of directional derivative we have that
= = =
Jh′ (uh )(vh − uh ) ∫ ∫ T } 1{ T g(yh (uh + s(vh − uh ))) + j(uh + s(vh − uh ))dt − g(yh (uh )) + j(uh )dt lim s→0 s 0 0 ∫ ∫ 1 T 1 T lim g(yh (uh + s(vh − uh ))) − g(yh (uh ))dt + lim j(uh + s(vh − uh )) − j(uh )dt s→0 s 0 s→0 s 0 ∫ T (g ′ (yh (uh )), yh′ (uh )(vh − uh )) + (j ′ (uh ), vh − uh )dt. 0
8
(3.9)
Next, let us differentiate the state equation of (OCP )h at uh in the direction vh − uh , } 1 { ∂(yh (uh + s(vh − uh )) − yh (uh )) ( , wh ) + Ah (yh (uh + s(vh − uh )) − yh (uh ), wh ) s ∂t 1 { ∂(yh (uh + s(vh − uh )) − yh (uh )) = ( , wh ) s ∂t ∑ |τ | } x + (A(x, )∇(yhε (uh + s(vh − uh )) − yhε (uh )), ∇whε )Qt |Qt | ε τ ∈Th
= (B(vh − uh ), wh ), ∀ wh ∈ Y h . Taking the limit as s → 0, we obtain (
∂yh′ (uh ) , wh ) + Ah (yh′ (uh )(vh − uh ), wh ) = (B(vh − uh ), wh ), ∀ v h ∈ K h , ∀ wh ∈ Y h . ∂t
(3.10)
Define the co-state ph ∈ Y h satisfying ∂ph , qh ) + Ah (qh , ph ) = (g ′ (yh ), qh ), ∀ qh ∈ Y h , ∂t ph (x, T ) = 0.
−(
150
(3.11)
Adopting wh = ph in (3.10) and integrating by parts, it implies from (3.11) ∫ T ∫ T (B(vh − uh ), ph )dt = (vh − uh , B ∗ ph )U dt 0
=
∫
0
T
0
=
∫
0
T
∂y ′ (uh ) ( h , ph ) + Ah (yh′ (uh )(vh − uh ), ph )dt ∂t
(g ′ (yh ), yh′ (uh )(vh − uh ))dt.
(3.12)
So by (3.9) and (3.12), the optimality condition reads ∫ T Jh′ (uh )(vh − uh ) = (j ′ (uh ) + B ∗ ph , vh − uh )U dt ≥ 0, 0
where ph is defined in (3.11). This completes the proof of Theorem 3.2.
∀ vh ∈ K h ,
3.2. The fully discrete approximation scheme Let 0 = t0 < t1 < . . . < tN −1 < tN = T , ki = ti − ti−1 , i = 1, 2, . . . , N , k = maxi∈[1,N ] {ki }. For i =
155
1, 2, . . . , N , construct the finite element spaces Y h,i ∈ H01 (Ω) (similar to Y h ) with the mesh T h,i , and construct
the finite element spaces U h,i ∈ L2 (ΩU ) (similarly as U h ) with the mesh (TUh )i . Let hτ i (hτUi ) denote the
maximum diameter of the element τ i (τUi ) in Tih ((TUh )i ). Define mesh functions τ (·), τU (·) and mesh functions
hτ (·), hτU (·) such that τ (t)|t∈(ti−1 ,ti ] = τ i , τU (t)|t∈(ti−1 ,ti ] = τUi , hτ (t)|t∈(ti−1 ,ti ] = hτ i , hτU (t)|t∈(ti−1 ,ti ] = hτUi .
160
For ease of exposition, in the following we will denote τ (t), τU (t), hτ (t) and hτU (t) by τ , τU , hτ and hτU , ∩ respectively. Let Kih ⊂ Uih K. The fully discrete approximation scheme (QCP )hk is to find (yhi , uih ) ∈ Yih × Kih , i = 1, 2, . . . , N , such that min {
uih ∈Kih
N ∑
ki (g(yhi ) + j(uih ))}
i=1
9
yhi − yhi−1 , wh ) + Ah (yhi , wh ) = (f (x, ti ) + Buih , wh ), ∀ wh ∈ Yih ⊂ H01 (Ω), i = 1, 2, . . . , N, ki yh0 (x) = y0h (x), x ∈ Ω. (
Similarly to the semi-discrete case, we know that the control problem (QCP )hk has a unique solution (Yhi , Uhi ), i = 1, 2, . . . , N , and that a pair (Yhi , Uhi ) ∈ Yih × Kih , i = 1, 2, . . . , N , is the solution of (QCP )hk if and
only if there is a co-state Phi−1 ∈ Yih , i = 1, 2, . . . , N , such that the triplet (Yhi , Phi−1 , Uhi ) ∈ Yih × Yih × Kih ,
165
i = 1, 2, . . . , N , satisfies the following optimality conditions: (QCP − OP T )hk Y i −Y i−1 ( h kih , wh ) + Ah (Yhi , wh ) = (f (x, ti ) + BUhi , wh ), ∀ wh ∈ Yih , i = 1, 2, . . . , N, Yh0 (x) = y0h (x), x ∈ Ω, P i−1 −P i ( h ki h , qh ) + Ah (qh , Phi−1 ) = (g ′ (Yhi ), qh ), ∀ qh ∈ Yih , i = N, N − 1, . . . , 1, P N (x) = 0, x ∈ Ω, h (j ′ (U i ) + B ∗ P i−1 , v − U i ) ≥ 0, U i ∈ K h , ∀ v ∈ K h ⊂ U h ∩ K, i = 1, 2, . . . , N. h
h
h
h U
h
i
h
i
(3.13)
i
For i = 1, 2, . . . , N , let
(ti − t)Yhi−1 + (t − ti−1 )Yhi , ki (ti − t)Phi−1 + (t − ti−1 )Phi Ph |(ti−1 ,ti ] = , ki Uh |(ti−1 ,ti ] = Uhi . Yh |(ti−1 ,ti ] =
For any function ω ∈ C(0, T ; L2 (Ω)), let ω ˆ (x, t)|t∈(ti−1 ,ti ] = ω(x, ti ), ω ˜ (x, t)|t∈(ti−1 ,ti ] = ω(x, ti−1 ). Then the optimality conditions (3.13) can be restated as h h ˆ ˆ ( ∂Y ∂t , wh ) + Ah (Yh , wh ) = (f + BUh , wh ), ∀ wh ∈ Yi , i = 1, 2, . . . , N, Yh (x, 0) = y0h (x), x ∈ Ω, ′ ˆ h h ˜ −( ∂P ∂t , qh ) + Ah (qh , Ph ) = (g (Yh ), qh ), ∀ qh ∈ Yi , i = N, N − 1, . . . , 1, Ph (x, T ) = 0, x ∈ Ω, (j ′ (U ) + B ∗ P˜ , v − U ) ≥ 0, U ∈ K h , ∀ v ∈ K h ⊂ U h ∩ K, i = 1, 2, . . . , N. h h h h U h h i i i
10
(3.14)
4. A priori error estimate
170
In this section, we consider the a priori estimates for the full-discrete heterogeneous multiscale finite element approximation (OCP − OP T )hk . Similarly to the continuous case, we need an auxiliary problem: Jh′ (uε )(w − uε ) =
N ∑ i=1
ki (j ′ (uεi ) + B ∗ Phi−1 (uε ), wi − uεi ), ∀ w ∈ K ⊂ U,
where (Yhi (uε ), Phi (uε )) ∈ Yih × Yih is the solution of the system: Y i (uε )−Yhi−1 (uε ) , wh ) + Ah (Yhi (uε ), wh ) = (f (x, ti ) + Buεi , wh ), ∀ wh ∈ Yih , i = 1, 2, . . . , N, ( h ki Y 0 (uε ) = y0h (x), x ∈ Ω, h
P i−1 (uε )−Phi (uε ) ( h , qh ) + Ah (qh , Phi−1 (uε )) = (g ′ (Yhi (uε )), qh ), ∀ qh ∈ Yih , i = N, N − 1, . . . , 1, ki P N (uε ) = 0, x ∈ Ω. h
(4.1)
(4.2)
For y ∈ Y h , we can define the norm about discretization time t, ∥y∥Lp (0,T ;H 1 (Ω)) = (
N ∑ i=1
1
ki ∥y(ti , x)∥pH 1 (Ω) ) p ,
for 1 ≤ p < ∞, and if p = ∞, let ∥y∥L∞ (0,T ;H 1 (Ω)) = max ∥y(ti , x)∥H 1 (Ω) . 1≤i≤N
175
Similarly, we can define ∥u∥Lp (0,T ;L2 (ΩU )) and ∥y∥Lp (0,T ;L2 (Ω)) . Lemma 4.1. Under the definition of (4.1), we have the following estimate: Jh′ (w)(w − uε ) − Jh′ (uε )(w − uε ) ≥ α∥w − uε ∥2L2 (0,T ;L2 (ΩU )) .
(4.3)
Proof. From (4.1), we have Jh′ (w)(w − uε ) − Jh′ (uε )(w − uε )
=
N ∑ i=1
(4.4)
ki ((B ∗ (Phi−1 (w) − Phi−1 (uε )), wi − uεi ) + (j ′ (wi ) − j ′ (uεi ), wi − uεi )).
Note that g(·) and j(·) are two convex differentiable functions, then from (4.2) we have N ∑ i=1
=
N ∑ i=1
ki (B ∗ (Phi−1 (w) − Phi−1 (uε )), wi − uεi ) ki ((Phi−1 (w) − Phi−1 (uε )), B(wi − uεi )) 11
(4.5)
=
N ∑ {ki Ah (Yhi (w) − Yhi (uε ), Phi−1 (w) − Phi−1 (uε )) i=1
+((Yhi (w) − Yhi (uε )) − (Yhi−1 (w) − Yhi−1 (uε )), Phi−1 (w) − Phi−1 (uε ))}
=
N ∑ {ki (g ′ (Yhi (w)) − g ′ (Yhi (uε )), Yhi (w) − Yhi (uε )) + (Phi (w) − Phi (uε ), Yhi (w) − Yhi (uε )) i=1
−(Phi−1 (w) − Phi−1 (uε ), Yhi−1 (w) − Yhi−1 (uε ))} =
N ∑ i=1
ki (g ′ (Yhi (w)) − g ′ (Yhi (uε )), Yhi (w) − Yhi (uε )) ≥ 0.
Furthermore, by (4.4) and assumption A3 , we can obtain (4.3).
180
Lemma 4.2. Let (Yhi , Phi ) and (Yhi (uε ), Phi (uε )) be the solutions of (3.13) and (4.2) respectively, then ∥Yh − Yh (uε )∥L∞ (0,T ;H 1 (Ω)) ≤ C∥uε − U h ∥L2 (0,T ;L2 (ΩU )) ,
(4.6)
∥Ph − Ph (uε )∥L∞ (0,T ;H 1 (Ω)) ≤ C∥uε − U h ∥L2 (0,T ;L2 (ΩU )) .
(4.7)
Proof. Let β i = Yhi − Yhi (uε ) and θi = Phi − Phi (uε ). From (3.13) and (4.2), we can get β i − β i−1 , wh ) + Ah (β i , wh ) = (B(Uhi − uεi ), wh ), ki
(4.8)
θi−1 − θi , qh ) + Ah (qh , θi−1 ) = (g ′ (Yhi ) − g ′ (Yhi (uε ), qh ). ki
(4.9)
( ( In equation (4.8), let wh = ki ( 185
β i −β i−1 , ki
we have
β i − β i−1 β i − β i−1 β i − β i−1 , ) + Ah (β i , β i ) = Ah (β i , β i−1 ) + ki (B(Uhi − uεi ), ). ki ki ki
(4.10)
Then, from triangle inequality and Cauchy-Schwartz inequality, we deduce ki ( ≤ So we can get
β i − β i−1 β i − β i−1 , ) + Ah (β i , β i ) ki ki
(4.11)
1{ β i − β i−1 β i − β i−1 } Ah (β i , β i ) + Ah (β i−1 , β i−1 ) + ki ∥Uhi − uεi ∥2L2 (ΩU ) + ki ( , ) . 2 ki ki Ah (β i , β i ) ≤ Ah (β i−1 , β i−1 ) + ki ∥Uhi − uεi ∥2L2 (ΩU ) .
(4.12)
Then for 1 ≤ i ≤ N , it follows Ah (β i , β i ) ≤ Ah (β 0 , β 0 ) +
i ∑ l=1
kl ∥Uhl − uεl ∥2L2 (ΩU ) ≤ ∥uε − U h ∥L2 (0,T ;L2 (ΩU )) .
Combining the last inequality with Poincare inequality, it follows (4.6). Similarly, from (4.9), we can obtain the inequality (4.7).
12
(4.13)
190
Now, we introduce the local averaging operator πhτ i given by U
∫
τi
(πhτ i w)|τUi = ∫ U U
i τU
w 1
, ∀ τUi ∈ (TUh )i .
(4.14)
It is obvious that if w ∈ K then πhτ i w ∈ Kih . Finally, we can give the main result of this paper. U
Theorem 4.1. Let (uε , y ε , pε ) and (Uh , Yh , Ph ) be the solutions of (2.4) and (3.14), respectively. Then the following estimates hold: 1
∥uε − Uh ∥L2 (0,T ;L2 (ΩU )) + ∥y ε − Yh ∥L∞ (0,T ;L2 (Ω)) + ∥pε − Ph ∥L∞ (0,T ;L2 (Ω)) ≤ C(k + hU + h2 + ε 2 ), 1
∥uε − Uh ∥L2 (0,T ;L2 (ΩU )) + ∥y ε − Yh ∥L2 (0,T ;L2 (Ω)) + ∥pε − Ph ∥L2 (0,T ;L2 (Ω)) ≤ C(k + hU + h2 + ε 2 ). 195
(4.15) (4.16)
Proof. From (2.4), (3.14), (4.1) and (4.3), we know c∥uε − Uh ∥2L2 (0,T ;L2 (ΩU ))
≤ Jh′ (uε )(uε − Uh ) − Jh′ (Uh )(uε − Uh ) =
N ∑ i=1
=
N ∑ i=1
+
ki (j ′ (uεi ) + B ∗ Phi−1 (uε ), uεi − Uhi ) − ki (j ′ (uεi ) + B ∗ pεi , uεi − Uhi ) +
N ∑ i=1
≤ For the term
∑N
i=1
N ∑ i=1
N ∑ i=1
U
ki (j ′ (Uhi ) + B ∗ Phi−1 , πhτ i uεi − uεi ) + U
N ∑ i=1
N ∑ i=1
ki (B ∗ Phi−1 (uε ) − B ∗ pεi , uεi − Uhi )
ki (B ∗ Phi−1 (uε ) − B ∗ pεi , uεi − Uhi ).
ki (B ∗ Phi−1 , πhτ i uεi − uεi ), we have U
=
U
N ∑ i=1
ki (B ∗ pεi−1 , πhτ i uεi − uεi )
+
N ∑
N ∑ i=1
Note that
N ∑ i=1
(4.18)
U
i=1
i=1
i=1
i=1
ki (j ′ (Uhi ) + B ∗ Phi−1 , uεi − Uhi )
ki (j ′ (Uhi ) + B ∗ Phi−1 , Uhi − πhτ i uεi )
U
ki (B ∗ Phi−1 , πhτ i uεi − uεi )
N ∑
N ∑
ki (j ′ (Uhi ) + B ∗ Phi−1 , πhτ i uεi − uεi ) +
+
and
(4.17) N ∑
ki (B ∗ (pεi−1 − Phi−1 (uε )), uεi − πhτ i uεi ) U
ki (B ∗ (Phi−1 (uε ) − Phi−1 ), uεi − πhτ i uεi ). U
ki (j ′ (Uhi ), πhτ i uεi − uεi ) = 0
ki (B ∗ pεi−1 , πhτ i uεi − uεi ) = U
(4.19)
U
N ∑ i=1
ki (B ∗ pεi−1 − πhτ i B ∗ pεi−1 , πhτ i uεi − uεi ), U
(4.20)
U
combining (4.17)-(4.20) with triangle inequality and Cauchy-Schwartz inequality, we have c∥uε − Uh ∥2L2 (0,T ;L2 (ΩU ))
(4.21) 13
≤
C(δ)
N ∑ i=1
+C(δ)
ki ∥uεi − πhτ i uεi ∥2L2 (ΩU ) + C(δ) U
N ∑ i=1
+Cδ
N ∑ i=1
200
N ∑ i=1
ki ∥pεi−1 − Phi−1 (uε )∥2L2 (Ω) + C(δ)
ki ∥Phi−1 (uε ) − Phi−1 ∥2L2 (Ω) + Cδ
N ∑ i=1
ki ∥pεi−1 − πhτ i pεi−1 ∥2L2 (Ω) U
N ∑ i=1
ki ∥pεi−1 − pεi ∥2L2 (Ω)
ki ∥uεi − Uhi ∥2L2 (Ω) .
It follows from the standard interpolation error estimate (see, e.g., [30]) and (4.7) that c∥uε − Uh ∥2L2 (0,T ;L2 (ΩU )) ≤ Ch2U + Ck 2 + C(δ)
N ∑ i=1
(4.22)
ki ∥pεi−1 − Phi−1 (uε )∥2L2 (Ω) + Cδ∥uε − Uh ∥2L2 (0,T ;L2 (ΩU )) .
Note that Ph (uε ) and Yh (uε ) are the standard HMM approximations of y ε and pε , respectively. It has been proved in Theorem 1.1 and Theorem 1.2 of [29] that N ∑ i=1
N ∑ i=1
205
i=1
1
ki ∥y εi−1 − Yhi−1 (uε )∥H 1 (Ω) ≤ C(k + h2 + ε 2 ), 1
1≤i≤N
i=1
i=1
N ∑
ki ∥y εi−1 − Yhi−1 (uε )∥L2 (Ω) ≤ CT max ∥y εi−1 − Yhi−1 (uε )∥L2 (Ω) ≤ C(k + h2 + ε 2 ),
N ∑ N ∑
ki ∥y εi−1 − Yhi−1 (uε )∥L2 (Ω) ≤
ki ∥pεi−1 − Phi−1 (uε )∥L2 (Ω) ≤
N ∑ i=1
1
ki ∥pεi−1 − Phi−1 (uε )∥H 1 (Ω) ≤ C(k + h2 + ε 2 ), 1
ki ∥pεi−1 − Phi−1 (uε )∥L2 (Ω) ≤ CT max ∥pεi−1 − Phi−1 (uε )∥L2 (Ω) ≤ C(k + h2 + ε 2 ). 1≤i≤N
Combing (4.22)-(4.26), let δ =
c 2C
and from Lemma 4.2, we can obtain (4.15) and (4.16).
(4.23)
(4.24)
(4.25)
(4.26)
Remark 4.1. The conclusion of Theorem 4.1 is given for the case Ω ̸= ΩU . So from [31] Theorem 6.1, we 1
know that if hU ≤ h, the L2 error of uε − uh could be O(k + hU + h2 + ε 2 ) and the L2 error of y ε − yh or 1
1
pε − ph could be O(k + h2 + ε 2 ). If Ω = ΩU and Th = ThU , the L2 error of uε − uh could be O(k + h + ε 2 ) and 1
the L2 error of y ε − yh or pε − ph could be O(k + h2 + ε 2 ).
210
5. Numerical Example In this section, we present our numerical results to validate the theoretical error estimates derived for the control, state and co-state variables. For this purpose, we consider the following optimal control problem. For the simplicity, we only carry out it under the case ΩU = Ω and use the same mesh for the control, state and adjoint state. Since it is difficult to give an exact solution for the parabolic equation with oscillatory coefficients,
215
we follow the example as in [22, 32] and give the model of our optimal control problem.
14
Example 1 This example is to solve the following problem: ∫ ∫ ∫ ∫ 1 T α T ε 2 ε min J(u ) = min (y − y ) dxdt + (uε )2 dxdt d uε ∈K u∈K 2 0 2 0 Ω Ω ∂y ε ∂t
s.t.
− ∇(A∇y ε ) = uε + f, y ε (x) = 0, ε
on ∂Ω × (0, T ],
y (x, 0) = 0, Here, the matrix A=
(5.1)
inΩ × (0, T ],
on Ω.
2 , x 2+cos(2π ε1 )
1 8π 2 0.0,
0.0 1 + 0.5 cos(2π xε1 )
,
the regularization parameter α = 1, space domain Ω = (0, 1)2 , time T = 1 and the source term f and yd will be of the form f (x, t) =
ε x1 cos 2πx1 sin 2πx2 sin 2π ) 2 ε 1 x1 +t(sin 2πx1 sin 2πx2 + sin 2πx1 sin 2πx2 cos 2π 4 ε 1 + 2 cos 2π xε1 1 ε x1 + sin 2πx1 sin 2πx2 sin 2πx2 cos 2πx1 sin 2π x1 2 + 2 (2 + cos 2π ε ) 2 ε
(sin 2πx1 sin 2πx2 +
+
ε sin 2πx2 cos 2πx1 sin 2π xε1 (cos 2π xε1 )2 ) − uε (x, t), 8 2 + cos 2π xε1
220
yd (x, t) =
ε x1 cos 2πx1 sin 2πx2 sin 2π ) 2 ε 1 x1 −(1 − t)(sin 2πx1 sin 2πx2 + sin 2πx1 sin 2πx2 cos 2π 4 ε 1 + 2 cos 2π xε1 1 ε x1 + sin 2πx1 sin 2πx2 + sin 2πx2 cos 2πx1 sin 2π 2 (2 + cos 2π xε1 )2 2 ε
y ε (x, t) − (sin 2πx1 sin 2πx2 +
+
ε sin 2πx2 cos 2πx1 sin 2π xε1 (cos 2π xε1 )2 ). 8 2 + cos 2π xε1
To show the order of convergence, we would require the exact solution of the above problem. So we adopt K = {v ∈ L2 (0, T ; L2 (Ω))|v ≥ 0}. Moreover, with the choice of the source term f and desired state yd , the
exact state y ε and co-state pε will be given in the following manner:
ε x1 cos 2πx1 sin 2πx2 sin 2π ), 2 ε ε x1 pε (x, t) = (1 − t)(sin 2πx1 sin 2πx2 + cos 2πx1 sin 2πx2 sin 2π ). 2 ε y ε (x, t) = t(sin 2πx1 sin 2πx2 +
225
Therefore, the control variable is defined as 1 uε (x, t) = max(− pε (x, t), 0). α Firstly, in order to show the convergence about h, we use the uniform time step with k = 0.1 and ε = 1.0e−6. The control variable uε is approximated by piecewise constant functions, and the state y ε and adjoint state pε are approximated by the linear FE-HMM method. The numerical results are shown in Table 1. 15
mesh
Uh − uε
Yh − y ε
Ph − pε
h
L2 (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
0.0625
0.0540
0.0457
0.0489
0.0562
0.0944
0.03125
0.0253
0.0101
0.0124
0.0131
0.0223
0.015625
0.0128
0.0024
0.0031
0.0031
0.0053
0.0078125
0.0070
5.9448e-004
7.6552e-004
7.5981e-004
0.0013
Table 1: k = 0.1 and ε = 1.0e − 6.
From Table 1, we can see that the error is reduced when the mesh decrease, and also the convergence rate 230
O(h) for control is consistent with our theoretical result given in Theorem 4.1 for hU = h. However, better than the theoretical result O(h), the convergence rates for both the state and co-state are about O(h2 ) which achieve optimal case even when hU = h. The possible explanation is the better regularity for the state and co-state than the control. Secondly, for testing the convergence of k, we take h = 0.0078125 and ε = 1.0e − 6. The related results are
235
shown in Table 2. From Table 2, we can see that the error is reduced following the decrease of the time step. The convergence rates for control, state and co-state are consistent with our theoretical results for k. Time step
Uh − uε
Yh − y ε
Ph − pε
k
L2 (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
0.2
0.0153
1.2043e-003
1.5712e-003
1.52163e-003
0.0025
0.1
0.0070
5.9448e-004
7.6552e-004
7.5981e-004
0.0013
0.05
0.0038
3.1057e-004
4.0341e-004
3.9925e-004
6.9357e-004
0.025
0.0021
1.6417e-004
2.1255e-004
2.0915e-004
3.8241e-004
Table 2: h = 0.0078125 and ε = 1.0e − 6.
Finally, we also estimate the convergence of ε. Here, we adopt k = 0.01 and h = 0.0078125 for this computation. Then the results are shown in Table 3. Uh − uε
Yh − y ε
Ph − pε
ε
L2 (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
0.5
0.0147
1.1491e-003
1.4878e-003
1.4641e-003
0.0026
0.25
0.0107
8.2131e-004
1.0534e-003
1.0311e-003
0.0019
0.125
0.0079
5.9241e-004
7.5131e-004
7.3124e-004
0.0013
0.0625
0.0057
4.2901e-004
5.3053e-004
5.2003e-004
9.5780e-004
Table 3: h = 0.0078125 and k = 0.01.
16
From Table 3, we can see that the error is reduced as the ε decrease. However, the convergence rates for control, state and co-state with ε are not obvious from Table 3. For this purpose, now we define the convergence rates 240
for control, state and co-state as following: RL2 (L2 ) (v) =
RL∞ (L2 ) (v) =
log(∥v ε − Vh ∥L2 (0,T ;L2 (Ω)) /∥v ε − Vh ∥L2 (0,T ;L2 (Ω)) ) , log(ε/¯ ε)
log(∥v ε − Vh ∥L∞ (0,T ;L2 (Ω)) /∥v ε − Vh ∥L∞ (0,T ;L2 (Ω)) ) . log(ε/¯ ε)
Then the results for ε are shown in Table 4, and we can see that the convergence rates for control, state and co-state are all consistent with our theoretical results for given ε in Theorem 4.1. ε
RL2 (L2 ) (u)
RL2 (L2 ) (y)
RL∞ (L2 ) (y)
RL2 (L2 ) (p)
RL∞ (L2 ) (p)
0.5
-
-
-
-
-
0.25
0.4582
0.4845
0.4981
0.5058
0.4525
0.125
0.4377
0.4713
0.4876
0.4958
0.5475
0.0625
0.4709
0.4656
0.5020
0.4918
0.4407
Table 4: h = 0.0078125 and k = 0.01.
Example 2 This example is to solve the problem as Example 1 and the matrix is: x x sin(2π ε1 ) sin(2π ε2 ) + + sin(4x x ) + 1, 0.0 x2 x1 1 2 , A = 0.1 sin(2π ε ) cos(2π ε ) x x sin(2π ε1 ) sin(2π ε2 ) 0.0, + + sin(4x x ) + 1 x2 x1 1 2 sin(2π cos(2π ) ) ε
ε
the regularization parameter α = 0.1, space domain Ω = (0, 1)2 , time T = 1, K = {v ∈ L2 (0, T ; L2 (Ω))|v ≥ 0} 245
and the source term f and yd will be of the form f (x, t)
(
( ) ( ) ( )) π sin 2 πεx1 sin 2 πεx2 π cos 2 πεx1 2 x2 cos(4 x1 x2 ) ( )+ = sin(2 π t) + ( )2 5 5 ε sin 2 πεx2 5 ε cos 2 πεx1 ( ( ) ) 2 π x1 x1 sin(2 π x2 ) + sin(2 π x2 ) (x1 − 1) + 2 π x2 cos (x2 − 1) ε ( ) ( ) 4 π 2 x2 sin 2 πεx1 (x2 − 1) + sin(2 π t) 2 sin(2 π x2 ) − ε ( ) ( 2 π x2 ) ( 2 π x1 ) sin ε sin ε sin(4 x1 x2 ) 1 ( ( )+ )+ + 10 10 10 cos 2 πεx1 10 sin 2 πεx2 ( ( ) ) 2 π x1 −2 π cos(2 π t) x1 sin(2 π x2 ) (x1 − 1) + ε x2 sin (x2 − 1) ε ( ( ( ) ) ( )) π cos 2 πεx2 sin 2 πεx1 π cos 2 πεx2 2 x1 cos(4 x1 x2 ) ( )− + sin(2 π t) + ( )2 5 5 ε cos 2 πεx1 5 ε sin 2 π x2 ε
17
( ( ) ( ) ) 2 π x1 2 π x1 ε x2 sin + ε sin (x2 − 1) + 2 π x1 cos(2 π x2 ) (x1 − 1) ε ε ( ( ) ) 2 π x1 + sin(2 π t) 2 ε sin − 4 π 2 x1 sin(2 π x2 ) (x1 − 1) ε ( ) (2π x ) ( ) sin ε 2 sin 2 πεx1 sin(4 x1 x2 ) 1 ( )+ ( )+ + − uε (x, t), 10 10 10 cos 2 πεx1 10 sin 2 πεx2
yd (x, t)
( ) ( ) ( )) π sin 2 πεx1 sin 2 πεx2 π cos 2 πεx1 2 x2 cos(4 x1 x2 ) )+ ( = y (x, t) + sin(2 π t) + ( )2 5 5 ε sin 2 πεx2 5 ε cos 2 πεx1 ( ( ) ) 2 π x1 x1 sin(2 π x2 ) + sin(2 π x2 ) (x1 − 1) + 2 π x2 cos (x2 − 1) ε ( ) ( ) 4 π 2 x2 sin 2 πεx1 (x2 − 1) + sin(2 π t) 2 sin(2 π x2 ) − ε ( ) ( 2 π x2 ) ( 2 π x1 ) sin ε sin ε sin(4 x1 x2 ) 1 ( )+ ( )+ + 10 10 10 cos 2 πεx1 10 sin 2 πεx2 ( ( ) ) 2 π x1 +2 π cos(2 π t) x1 sin(2 π x2 ) (x1 − 1) + ε x2 sin (x2 − 1) ε ( ( 2 π x2 ) ( ) ( )) π cos ε π cos 2 πεx2 sin 2 πεx1 2 x1 cos(4 x1 x2 ) ( )− + sin(2 π t) + ( )2 5 5 ε cos 2 πεx1 5 ε sin 2 πεx2 ( ( ) ( ) ) 2 π x1 2 π x1 ε x2 sin + ε sin (x2 − 1) + 2 π x1 cos(2 π x2 ) (x1 − 1) ε ε ( ( ) ) 2 π x1 + sin(2 π t) 2 ε sin − 4 π 2 x1 sin(2 π x2 ) (x1 − 1) ε ) ( ( 2 π x2 ) ( ) sin ε sin 2 πεx1 1 sin(4 x1 x2 ) )+ )+ ( ( + . 10 10 10 cos 2 πεx1 10 sin 2 πεx2 ε
(
The exact solutions of the above problem are as following:
x1 )(1 − x2 )x2 ), ε x1 pε (x, t) = − sin(2πt)((1 − x1 )x1 sin(2πx2 ) + ε sin(2π )(1 − x2 )x2 ), ε 1 ε ε u (x, t) = max(− p (x, t), 0). α y ε (x, t) = sin(2πt)((1 − x1 )x1 sin(2πx2 ) + ε sin(2π
250
In this example, the control variable uε is approximated by piecewise constant functions, and the state y ε and adjoint state pε are approximated by the linear FE-HMM method. The numerical results for this example are shown in Table 5 - Table 8. The solutions for control, state and co-state are shown in Figure 2 - Figure 5 when ε = 0.0625 and t = 0.125, 0.25, 0.625, 0.75. From Table 5 - Table 8, we can see all the numerical results show the similar convergent behaviors as
255
the first example. Moreover, from the Figure 2 - Figure 5, it is obvious the HMM can catch the multiscale character well. Over all, from the computational load and the operation with the high oscillation, our numerical experiments indicate that the proposed algorithm can handle efficiently the optimal control problems with the highly oscillatory coefficients. 18
mesh
Uh − uε
Yh − y ε
Ph − pε
h
L2 (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
0.0625
3.7264e-002
3.1072e-003
3.6296e-003
3.1072e-003
3.6294e-003
0.03125
1.7941e-002
7.5786e-004
8.7423e-004
7.5786e-004
8.7431e-004
0.015625
9.0821e-003
1.8760e-004
2.0984e-004
1.8760e-004
2.0988e-004
0.0078125
4.4521e-003
4.4631e-005
5.0729e-005
4.4631e-005
5.0741e-005
Table 5: k = 0.1 and ε = 1.0e − 6.
Time step
Uh − uε
Yh − y ε
Ph − pε
k
L2 (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
0.2
9.1104e-003
9.4279e-005
1.1479e-4
9.4279e-005
1.1476e-4
0.1
4.4521e-003
4.4631e-005
5.0729e-5
4.4631e-005
5.0735e-5
0.05
2.2736e-003
2.2617e-005
2.4476e-005
2.2617e-005
2.4468e-005
0.025
1.0769e-003
1.1174e-005
1.2038e-005
1.1174e-005
1.2027e-005
Table 6: h = 0.0078125 and ε = 1.0e − 6.
Uh − uε
Yh − y ε
Ph − pε
ε
L2 (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
L2 (0, T ; L2 (Ω))
L∞ (0, T ; L2 (Ω))
0.5
7.7860e-003
8.1555e-005
8.4196e-005
8.1578e-005
8.4184e-005
0.25
5.5476e-003
5.6643e-005
5.9358e-005
5.6649e-005
5.9346e-005
0.125
4.0914e-003
4.1435e-005
4.2786e-005
4.1437e-005
4.2782e-005
0.0625
3.0133e-003
3.0093e-005
3.0241e-005
3.0084e-005
3.0249e-005
Table 7: h = 0.0078125 and k = 0.01.
ε
RL2 (L2 ) (u)
RL2 (L2 ) (y)
RL∞ (L2 ) (y)
RL2 (L2 ) (p)
RL∞ (L2 ) (p)
0.5
-
-
-
-
-
0.25
0.4890
0.5259
0.5043
0.5261
0.5044
0.125
0.4393
0.4510
0.4723
0.4511
0.4721
0.0625
0.4413
0.4614
0.5006
0.4619
0.5001
Table 8: h = 0.0078125 and k = 0.01.
6. Summary and conclusion
260
In this article, we study the heterogeneous multiscale finite element method for the constrained optimal control problem governed by parabolic equation with highly oscillatory coefficients, develop the fully discrete
19
Figure 2: The control(left), state(middle), co-state(right) in t = 0.125 with ε = 0.0625.
Figure 3: The control(left), state(middle), co-state(right) in t = 0.25 with ε = 0.0625.
Figure 4: The control(left), state(middle), co-state(right) in t = 0.625 with ε = 0.0625.
Figure 5: The control(left), state(middle), co-state(right) in t = 0.75 with ε = 0.0625.
20
scheme and derive a priori error estimates. From our numerical results, we can see that our HMM shows high performances in catching the multiscale character caused by the highly oscillatory coefficients and win much computational efficiency. Both from the theoretical analysis and numerical experiments, the HMM 265
is particularly suitable for such control problems. Whether the HMM can also work well for other more complicated multiscale optimal control problems, we will consider it in the further research.
References [1] J. Li, A multiscale finite element method for optimal control problems governed by the elliptic homogenization equations, Comput. Math. Appl. 60 (2010) 390–398. 270
[2] A. Brandt, Multi-level adaptive solutions to boundary value problems, Math. Comput. 31 (1977) 333–390. [3] J. Xu, L. Zikatanov, On multigrid methods for generalized finite element methods, Comput. Sci. Eng. 26 (2003) 401–418. [4] C. Duarte, J. Oden, An h-p adaptive method using clouds, Comput. Methods Anal. Mech. Engrg. 139 (1996) 237–262.
275
[5] M. Sarkis, Nonstandard coase spaces and schwarz methods for elliptic problems with discontinous coefficients using non-conforming elements, Numer. Math. 77 (1997) 383–406. [6] G. Allaire, Homogenization and two-scale convergence, SIAM J. Math. Anal. 23 (1992) 1482–1518. [7] M. Dorobantu, B. Engquist, Wavelet-based numerical homogenization, SIAM J. Numer. Anal. 35 (1998) 540–559.
280
[8] J. Melenk, I. Babuska, The partition of unity finite element method: basic theory and applications, Comput. Methods Appl. Mech. Engrg. 139 (1996) 289–314. [9] T. Hou, X. Wu, A multiscale finite element method for elliptic problems in composite materials and porous media, J. Comput. Phys. 134 (1997) 169–189. [10] W. E, B. Engquist, The heterogeneous multi-scale methods, Commun. Math. Sci. 1 (2003) 87–132.
285
[11] T. Hughes, G. Feijo, L. Mazzei, J. Quincy, The variational multiscale method a paradigm for computational mechanics, Comput. Methods Appl. Mech. Engrg. 166 (1998) 3–24. [12] J. Lions, Some methods in the mathematical analysis of systems and their control, Science Press, Beijing, 1981. [13] L. Cao, Multiscale asymptotic method of optimal control on the boundary for heat equations of composite
290
materials, J. Math. Anal. Appl. 343 (2008) 1103–1118.
21
[14] J. Liu, L. Cao, N. Yan, Multiscale asymptotic analysis and computation of optimal control for elliptic systems with constraints, SIAM J. Numer. Anal. 31 (2013) 1978–2004. [15] L. Cao, J. Liu, W. Allegretto, N. Yan, A multiscale approach for optimal control problems of linear parabolic equations, SIAM J. Control Optim. 50 (2012) 3269–3291. 295
[16] C. Fabre, J. Puel, E. Zuazua, Approximate controllability of the semilinear heat equation, Proc. Roy. Soc. Edinburgh Sect. A 125 (1995) 31–62. [17] S. Kesavan, J. Saint, Homogenization of an optimal control problem, SIAM J. Control Optim. 35 (1997) 1557–1573. [18] J. Lions, Remarks on the theory of optimal control of distributed systems, Control Theory of Systems
300
Governed by Partial Differential Equations, Academic Press, New York (1997) 1–103. [19] H. Zoubairi, Optimal control and homogenization in a mixture of fluids separated by a rapidly oscillating interface, Electron. J. Diff. Equa. 27 (2002) 1–32. [20] P. Ming, X. Yue, Numerical methods for multiscale elliptic problems, J. Comput. Phy. 214 (2006) 421–445. [21] A. Abdulle, On a priori error analysis of fully discrete heterogeneous multiscale fem, SIAM Multiscale
305
Model. Simu. 4 (2005) 447–459. [22] A. Abdulle, The finite element heterogeneous multiscale method: a computational strategy for multiscale pdes, GAKUTO Int. Ser. Math. Sci. Appl. 31 (2009) 135–184. [23] A. Abdulle, A. Nonnenmacher, A short and versatile finite element multiscale code for homogenization problems, Comput. Methods Appl. Mech. Engrg. 198 (2009) 2839–2859.
310
[24] W. E, P. Ming, P. Zhang, Analysis of the hetegeous multiscale method for elliptic homogenization problems, J. Amer. Math. Soc. 18 (2005) 121–156. [25] L. Ge, N. Yan, H. Wang, W. Liu, D. Yang, Heterogeneous multiscale method for optimal control problem governed by elliptic equations with highly oscillatory coefficients, J. Comput. Math. accepted. [26] V. Jikov, S. Kozlov, O. Oleinik, Homogenization of differential operators and integral functionals, Springer-
315
Verlag, Berlin, Heidelberg, 1994. [27] J. Lions, Optimal control of systems governed by partial differential equations, Springer-Verlag, Berlin, 1971. [28] W. Liu, N. Yan, Adaptive finite element methods for optimal control governed by pdes, Series in Information and Computational Science 41, Science Press, Beijing, 2008.
22
320
[29] P. Ming, P. Zhang, Analysic of the heterogeneous multiscale method for parabolic homogenization problems, Math. Comput. 76 (2007) 153–177. [30] P. Ciarlet, The finite element method for elliptic problems, North-Holland, 1978. [31] D. Yang, W. Liu, Y. Chang, A priori error estimate and superconvergence analysis for an optimal control problem of bilinear type, J. Comput. Math. 26 (2008) 471–487.
325
[32] P. Henning, M. Ohlberger, B. Schweizer, An adaptive multiscale finite element method, SIAM J. Multiscale Model. Simu. 12 (2014) 1078–1107.
23