Iterative Solution of Mesh Constrained Optimal Control Problems with Two-Level Mesh Approximations of Parabolic State Equation ()
1. Introduction
Optimal control of time-dependent production processes plays an important role in many real world applications. State constraints in optimal control of systems governed by partial differential equations have to be often included in the mathematical models. For instance, in continuous casting process a need to prevent the cracks in a slab and the solidification at a wrong place leads to the bounds on the temperature variable. Similar demands arise in the processes of crystal growth and cooling of glass melts (see articles [1] [2] [3] [4] [5] and bibliography therein). The introduction of pointwise state constraints yields adjoint variables and multipliers which only admit low regularity complicating both theoretical analysis and the constructing appropriate numerical methods.
There are few results in the area of numerical solution methods for the constrained parabolic optimal control problems. One of the approaches to overcome the difficulty connected with low regularity of the solutions is using Lavrentiev regularization. This approach has been used in [6] [7]. In [8] [9] [10] [11] a priori error estimates for space-time discretizations of linear-quadratic parabolic optimal control problems have been obtained for problems.
The implementations of the iterative methods for parabolic optimal control problems include the solution of the parabolic equation and corresponding adjoint parabolic equation at each iteration and this is the most time consuming part of the algorithms if applying the implicit (backward Euler) approximation of parabolic state equation. On the other hand, easily implementable explicit (forward Euler) approximation of a parabolic equation with a constant step in time requests extremely restrictive constraint for this step. Using explicit approximations of the parabolic equations with special series of non-uniform time steps allows partially avoid this deficiency. Such kind approximation is well-known for the differential equations [12] and they demonstrate the advantage in time of calculations in relation to the implicit schemes. Similar approximation have been used for the continuous casting problem [13] and recently applied to a parabolic state constrained problem [14]. One more approach to construct effective algorithm for parabolic optimal control problem is to use parareal approximation of the state equation [15] [16].
In this article we continue the investigations of [14] [16] [17] [18] [19] [20] on the iterative solution methods for the constrained saddle point problems with applications to optimal control problems. In the cited papers parabolic optimal control problems have been solved by using either backward or forward Euler approximating schemes for the state equation.
The main purpose of this article is to generalize these results for the case when any two-level scheme is used for the approximation of the parabolic state equation, including different splitting (locally one-dimensional) schemes.
To simplify the exposition we restrict ourselves to consider a problem in unit square with distributed control and observation, and to use finite difference schemes to approximate the state equation, while the most of the results can be extended for the case of lowest order finite element method, for a control in Neumann boundary condition, for final observation etc.
2. Formulation and Approximations of the Problem
Consider homogeneous Dirichlet initial-boundary value problem
(1)
in the cylinder
,
, with lateral surface
. We call
and
as state and control functions. It is well-known that for any
there exists a unique weak solution of problem (1) such that
, and the following inequality takes place ([21] [22] [23]):
(2)
Define objective function
with a given functions
, and the sets of constraints for state and control:
Above
and
. We solve the following optimal control problem:
(3)
Problem (3) has a unique solution (cf., e.g. [14]).
We construct the finite-difference approximations of problem (3) using a uniform in x and t mesh
in
.
Let
be the uniform mesh of the meshsize
on
,
be the set of the boundary nodes while
. By
we denote the space of mesh functions which are defined on
and vanish in the boundary nodes
, let
. The mesh on the time segment we denote by
.
We will use the notations
for the mesh functions from
and for the vectors of their nodal values as well. We also don’t distinguish the linear operators from
to
and corresponding them
matrices acting on the vectors of nodal values of mesh functions from
.
By
we denote the values on a time level
of a mesh function of x and t, and by
-euclidian norm in the space
. We use the following notations:
,
and
are the inner product and the norm in
. Let
be
unit matrix while
be the unit matrix in
.
We denote by
the linear operator such that
and similarly for
Obviously,
is the mesh Laplasian with homogeneous Dirichlet boundary conditions. Recall that we use the same notations for corresponding matrices. Matrices
are symmetric, commute and positive definite with spectrum in a segment
, where
is of order
, while
is bounded below by a constant which doesn’t depend on
.
Let us introduce
matrices U, D and R, where U is a symmetric and positive matrix, and approximate state problem (1) by two-level finite difference scheme
(4)
Below we give several examples of two-level finite difference scheme (4).
1) Finite difference scheme with weights:
(5)
with
can be written in the form (4) with
,
and
. Recall that (5) contains forward Euler (
), backward Euler (
) and Crank-Nicolson (
) schemes.
Scheme (5) is unconditionally stable for
and stable in the case
if
, where
is the maximal eigenvalue of
. The stability estimate is:
(6)
2) Two-level scheme with factorized preconditioner:
(7)
Mesh scheme (7) is unconditionally stable for
with stability estimate (6).
3) Fractional steps scheme:
(8)
System (8) can be written in form (4) with
Scheme (8) is unconditionally stable with stability estimate (6).
Let us further use notation
for an
-approximation of the function
. Define a mesh goal function and the sets of the constraints:
(9)
The mesh optimal control problem reads as follows:
(10)
Problem (10) has a unique solution because the set
is a convex compact and the quadratical function
is continuous and strictly convex on
.
3. Saddle Point Problem and Iterative Solution Method
Let us rewrite problem (10) in a “vector-matrix” form. Let the matrices
be defined by the equalities:
Denote by
and
the indicator functions of the sets
and
, respectively. We obtain the following algebraic form of mesh optimal control problem (10):
(11)
Below we suppose that for a two-level finite difference scheme (4) the following assumption is valid:
(12)
The stability estimate (6) approve that for all cited above particular cases of (4) assumption (12) is true (with constant
).
Remark 1 More well-known finite difference schemes can be written in form (4) and satisfy assumption (12): different kinds of ADI schemes proposed in [24], [25] and “sequential” variant of fractional steps scheme [26] [27] etc. Also a more general variant of scheme (8) (see [28]) with positive weights
, such that
instead of
as in (8) can be considered.
We construct Lagrange function for problem (11):
A saddle point of this Lagrangian satisfies the following system:
(13)
where
and
are the subdifferentials of the corresponding functions.
With the notations
,
,
and
(11) becomes a particular case of minimization problem
(14)
while (13)-a particular case of saddle point problem
(15)
We will use the following results from the article [17]:
Proposition 1 Let
(16)
(17)
(18)
(19)
Then
1) Problem (15) has a non-empty set of the solutions
, where z is unique solution of (14).
2) Uzawa-type iterative method
(20)
with a symmetric preconditioner
satisfying the inequality
(21)
converges for any initial guess
:
for
.
Theorem 1 Problem (13) has a solution
with unique pair
, which coincides with the solution of problem (11).
Proof. Obviously, all assumptions (16)-(18) and (19) of proposition 1 are satisfied for problem (13). In particular, vector
satisfies assumption (19).
Below we construct for problem (13) the easily implementable preconditioner D which is spectrally equivalent to
with constants independent on
and h. As this preconditioner we take
. Then method (20) for problem (13) with this preconditioner reads as follows:
(22)
Theorem 2 Method (22) converges for any initial guess
if
(23)
Proof. The matrix
is spectrally equivalent to
, namely,
(24)
where constant
is defined in (12). In fact, by direct calculation we find
. Further, since
, then
. This inequality is equivalent to
(25)
whence the result.
In virtue of inequality (24) the convergence condition (21) of proposition 1 is true for the parameter
from (23).
In the case of optimal control problem without state constraints we can estimate a rate of convergence for iterative method (22) and give an optimal iterative parameter
. Namely, the following statement holds:
Theorem 3 Let
. Then there exists a unique solution
of saddle point problem (13), and for theoretically optimal iterative parameter
the following estimate for the rate of convergence of method (22) is
valid:
(26)
Proof. The uniqueness of
is proved in theorem 1. The uniqueness of
in the case
follows from the equation
.
Vector
is the solution of the equation
(27)
while iterative method (22) can be written in the form
(28)
It is well-known (cf., e.g. [29]) that the operator
is co-coercive:
Because of this and (25) the operator
satisfies the following properties (strong monotonicity and Lipshitz-continuity):
The rest of the prove is quite standard (cf. [18]). Namely, with the notations
,
and
we have the equation
We multiply this equation by
and obtain the equality
Due to the properties of
the following estimate holds:
Substituting this estimate in the previous equality we get
whence
if
and rate of convergence (26) is true for optimal parameter
.
Remark 2 Due to the equalities
,
, we have the following estimate for the rate of convergence of
:
For the sequence of control vectors
the estimate is as follows:
Remark 3 For the state constrained problems there are no estimates for rate of convergence for iterative method (22). In this case instead of (27) we have the equation
which operator is only co-coercive. The convergence of the iterative methods for such kind of equations have been investigated in [18].
On the other hand, numerous calculations show that in this general case the choice of the iterative parameter
as in theorem 3 is practical and seems to be close to optimal one.
When implementing method (22) one has to solve the inclusions with respect
and
, and a system of linear equations with matrix
. Concerning solving the inclusions we underline that the matrices and the operators in them have diagonal form, so, their solving reduces to easy pointwise projection.
In turn, solving equation with the matrix
is equivalent to solving direct (with
) and adjoint (with
) parabolic mesh schemes. In the case of mesh schemes with factorized preconditioner (7) or fractional steps scheme (8) approximating state equation their solving reduces to solving set of non-coupled “one-dimensional” mesh problems-systems of linear algebraic equations with tridiagonal matrices
,
. Obviously, these systems can be solved by a well-known direct method and parallel.
4. Variants and Generalizations
The effectiveness of the implementation of iterative method (22) is based on two main properties:
・ Preconditioner has factorized form
and is spectrally equivalent to “main” matrix of the problem;
・ Equations with the matrices
and
as in (7) and (8) can be easily implementable.
The first property is ensured by the inequality (12):
. Just this inequality allows us to prove spectral equivalency of the matrix
and
with the constants independent on mesh parameters
and
. In turn, this inequality is nothing but a stability estimate for a corresponding two-level approximation of a parabolic state equation. Numerous classes of stable two-level finite difference schemes for the parabolic equations can be found in [27].
The second property-easy solution of the equations with matrices
,
in (7) and (8)-is the consequence of their “local one-dimensional” structure. This imposes several limitations to the domains, boundary conditions and using orthogonal meshes. Nevertheless, a lot of different mesh schemes with factorized preconditioner of the form (7), satisfying stability property (12) is known (cf., e.g. [24] [27] and the bibliography therein).
The results of this paper can be extended to the parabolic optimal control problems with other state and control constraints, such as, for example, in [14] and [30].
Acknowledgements
The research was supported by the RFBR grant No. 16-01-00408 and by grants No. 296348 and No. 296928 from Academy of Finland.