Time Optimal Control of a Thermoelastic System

This paper considers the numerical approximation for the time optimal control problemof a thermoelastic systemwith some control and state constraints. By the Galerkin finite element method (FEM), the original problem is projected into a semidiscrete optimal control problem governed by a system of ordinary differential equations. Then the optimal time and control parameterization method is applied to reduce the original system to an optimal parameter selection problem, in which both the optimal time and control are taken as decision variables to be optimized.This problem can be solved as a nonlinear optimization problem by a hybrid algorithm consisting of chaotic particle swarm optimization (CPSO) and sequential quadratic programming (SQP) algorithm.The numerical simulations demonstrate the effectiveness of the proposed numerical approximation method.


Introduction
The expansion and contraction of a material exposed to temperature changes are of great importance in many applications such as stabilization of satellite antennas, which is commonly with large temperature variations and characterized by the thermoelastic coupling [1].This phenomenon is modelled by a system consisting of a hyperbolic elasticity equation describing the displacement coupled with a parabolic equation for the temperature.Over the past few decades, there has been much work on control of thermoelastic systems including [1][2][3][4][5][6][7][8] and references therein.For more studies of the related problems of the thermoelasticity, we would like to mention some other related works [9,10].Very little, however, is known about the time optimal control of thermoelastic systems.
As one of the most important optimal control problems, time optimal control problems aim at driving the corresponding trajectory of the system from an initial state to a given target set in the shortest time, through applying constrained controllers.Time optimal control problems in the finite-dimensional case [11,12] were first studied, which were governed by a system of ordinary differential equations.
Then time optimal control problems were extended to the infinite-dimensional case; namely, the controlled equations are partial differential equations (PDEs) or abstract evolution equations in Banach spaces [13][14][15][16].The main theoretical tool to deal with time optimal control problems analytically is the Pontryagin minimum principle [17].However, practical problems are usually becoming too complex to allow analytical techniques alone.Thus numerical methods for solving these problems attracted more and more researchers' attention [18][19][20][21][22].
In this paper, we will study the numerical approximation for the time optimal control problem governed by the onedimensional linear thermoelastic system.As far as we know, little work has been done on numerical methods for solving this problem.Our method follows the discrete-then-optimize approach [23][24][25], whereby the PDE system is first projected onto the finite-dimensional subspace to obtain a system of ordinary differential equations (ODEs), and computational optimal control techniques are applied; then the resulting problem is solved by appropriate optimization methods.Specifically, by making use of the standard finite element method (FEM), the original problem was firstly projected into a time optimal control problem governed by a system of ODEs.Noting that the terminal time of this problem is unspecified, we will, then, introduce a transformation such that the optimal time can be regarded as an unknown system parameter.Moreover, control parametrization method, which involves approximating the control function by a piecewise-constant function with possible discontinuities at a set of preassigned switching points [26,27], will be carried out and the original problem can be reduced to an optimal parameter selection problem.In principle, this problem can be solved as a nonlinear optimization problem by standard mathematical programming algorithms such as SQP.
The SQP method is one of the nonlinear programming methods for constrained optimization.The basic idea of SQP is to make a quadratic approximation of the Lagrange function of the original problem, to form a quadratic programming (QP) subproblem.The solution of QP subproblem is then used to form a search direction for a line search procedure in the next iteration.Although with the strong capability of local searching, the SQP method is dependent on the initial estimation and is easy to be trapped into a local minimum for multimodal problem.Thus a hybrid PSO is introduced to overcome this drawback.
The PSO algorithm is one of the most well-known bioinspired algorithms used in optimization problems, which basically consists of a number of particles that collectively move in search of the global optimum.Meanwhile, similar with other evolution algorithms, the traditional PSO is also prone to get trapped into the local minimum when solving some complex problems.There is a rich literature devoted to the improvement of PSO; for example, refer to [28][29][30] and references therein.In our work, a chaotic method is integrated with PSO to further promote the diversity of PSO [31][32][33].This chaotic PSO algorithm is referred to as CPSO.In this way, the possibility of exploring a global optimum solution in problems is increased.Furthermore, as empirical and theoretical studies show that PSO is sensitive to the control parameters, especially the inertia weight and acceleration coefficients, we also make full use of existing results of parameter controlling strategies to improve the performance of CPSO.
In this paper, a hybrid algorithm consisting of CPSO and SQP is introduced to solve the optimal parameter selection problem, where CPSO can be viewed as the global optimizer while the SQP is employed for the local search.
The rest of this paper is organized as follows.The formulation of the time optimal control problem is introduced in Section 2. In Section 3, we obtain the semidiscrete form of the original optimal control problem with FEM method.In Section 4, an optimal parameter selection problem is obtained by the optimal time and control parameterization method, and a hybrid optimization algorithm CPSO-SQP is employed to solve this problem.In Section 5, the numerical simulations are presented to illustrate the effectiveness of our numerical approximation method.Finally, in Section 6, we conclude the paper with summary comments and suggestions for future research.

Formulation of the Time Optimal Control Problem
Consider a linear one-dimensional thermoelastic rod of length  with unit reference density.The controlled thermoelastic system takes the form where (, ) is the displacement, (, ) represents the temperature deviation from the reference temperature of the bar at (, ),  0 ,  1 , and  0 are the given initial data, (, ⋅), (, ⋅) ∈  2 (0, ) are the external effects acting on this system, and  1 and  2 are the boundary control functions.In the above, the first equation is the momentum equation and the second equation is the energy equation.We assume that the system parameters , , , and  are given constants with ,  > 0, and ,  ̸ = 0, which are dependent on the material properties.
Define  = { ∈  1 (0, ) | (0) = 0} with the norm By the famous Poincare's theorem, we see that  is a Hilbert space.Moreover, let  =  2 (0, ) with the standard norm ‖ ⋅ ‖.Define H =  ×  ×  equipped with the norm Then, H is a Hilbert space and the energy function of the thermoelastic system (1) is where (, ) is the solution of (1).

Mathematical Problems in Engineering 3
Define where   ,   ∈ R,  = 1, 2, are given constants.Let  ∈ R be positive constant and assume that, for given ( 0 ,  1 ,  0 ) ∈ H, there exists at least one control  ∈ U, such that the solution (, ) of ( 1) satisfies that () ≤  for some  > 0. Under this setting, the time optimal control problem of the thermoelastic system, which we shall study, is formulated as follows: where the energy function () = (; ) is defined by ( 4) and (, ) is the solution of ( 1) with respect to ( 0 ,  1 ,  0 ) ∈ H and  ∈ U. Define as the optimal time for the problem (TOCP).Moreover, if there exists  * ∈ U such that ( * ; ) ≤ , we call  * as an optimal control for the problem (TOCP) and the solution ( * ,  * ) of ( 1) corresponding to the optimal control  * as the optimal state.Throughout this paper, we assume that there exists a pair (, ) with  ∈ U such that (; ) ≤ .By the standard arguments [34,35], we can get the existence of the optimal time  * and the optimal control  * .

Finite Element Approximation
In this section, we will apply the Galerkin finite element method [36,37] to discretize the original problem (TOCP), an infinite-dimensional distributed parameter time optimal control problem, into a semidiscrete approximation problem governed by a finite-dimensional lumped parameter system.In order to construct a finite-dimensional subspace  ℎ of the space , we subdivide [0, ] into  subintervals by the grid points: where the length of each subinterval is ℎ  =   −  −1 ,  = 1, 2, . . ., .Then ℎ = max ℎ  is the mesh size of this triangulation and a measure of how fine the partition is.Note that in general the mesh points   need not be equally spaced.
For simplification, we choose a uniform spacing here.Now let  ℎ be the set of functions V ℎ such that V ℎ is linear on each subinterval   and continuous on [0, ]: where  1 (  ) is the space of polynomials of degree at most 1 defined on the subinterval   .We observe that  ℎ is a -dimensional space and  ℎ ⊂ .Let {  ;  = 1, 2, . . ., } be the basis functions of the finitedimensional space  ℎ , defined by for  = 1, 2, . . ., ,  = 1, 2, . . ., , where   is the continuous piecewise linear function that takes value 1 at node point   and the value 0 at other node points.
Substituting ( 14) and ( 15) into ( 13) and taking V ℎ =   ,  = 1, 2, . . ., , yield that Define Thus, it follows from ( 16) that Define Ẋ =  1 and ) . ( Thus, ( 13) can be rewritten as follows: where Note that this system is a set of initial value problems (IVPs) and is equivalent to (13).Define where (⋅, ⋅) R 3 denotes the inner product in R 3 and Now, by the finite element method, we get a new time optimal control problem (TOCP ℎ ) governed by a system of controlled ordinary differential equations from the original problem (TOCP) as follows: where () : [0, ] → R 3 is the solution of (20).

Time Optimal Control Computation
4.1.Optimal Time Parameterization.Obviously, the problem (TOCP ℎ ) is a problem with unspecified terminal time.In this case, we may regard the optimal time  as an unknown system parameter.Write and define the transformation Let Then, ( 20) is equivalent to the following equation: Thus, the problem (TOCP ℎ ) can be transformed to the following form: where  ℎ (; ) = √ ((), Π()) R 3 for  ∈ [0, 1].
By this optimal time parameterization process, the problem (TOCP ℎ ) becomes a combined optimal parameter selection and optimal control problem, which has a fixed time horizon [0, 1].

Control Parameterization. To approximate the problem (TOCP
ℎ ), we will use the control parametrization method [27,38].More specifically, the admissible controls will be restricted to a suitable linear combination of basis functions, where the coefficients in the linear combination are decision variables to be selected optimally.This produces a sequence of finite-dimensional optimization problems for approximating problem (TOCP ℎ ).As piecewise-constant basis functions are normally applied in control parameterization, the optimal control problem (TOCP ℎ ) can be reduced to an optimal parameter selection problem.We subdivide the time horizon [0, 1] into  subintervals [ −1 ,   ),  = 1, 2, . . ., , where  ≥ 1 is the number of approximation subintervals, [ −1 ,   ) is the th approximation subinterval, and   ,  = 0, 1, . . ., , are fixed knot points such that Define The vector-valued control function () can be approximated by a constant vector on each subinterval: Then we restrict the control function to the piecewiseconstant control of the form where  ∈ U  and  [ −1 ,  ) is the indictor function of the interval [ −1 ,   ) defined by Corresponding to optimal time parameterization, the control function is transformed as follows: where By taking the approximation   () = ∑  =1    [ −1 ,  ) () of , (28) takes the following form: and   (⋅ | ) denote the solution of system (37) corresponding to  ∈ U  .Now, we may state the approximate optimal parameter selection problem as follows: where is the solution of (37).

Optimization Implementation.
The finite element method and the control parameterization method have been applied to various optimal control problems.The convergence for approximation schemes of these two methods has been discussed by many authors (see [18,20,39] for the FEM or [26,38,40] for control parameterization method).Inspired by these results, by making use of the classical space discretization scheme with FEM, we first project the original time optimal control problem (TOCP) into a semidiscrete problem (TOCP ℎ ) governed by a system of ODEs.Then, through the time parameterization and control parameterization, the semidiscrete approximation (TOCP ℎ ) of the original problem has been reduced to an optimal parameter selection problem (TOCP ℎ  ) in the canonical form, where  0 ≥ 0 is the time parameter and  ∈ U  is the control parameter.The optimal parameter selection problem is basically nonlinear mathematical programming problem and hence many gradient-based optimization techniques can be applied to optimal parameter selection problem.Combined with the algorithms for computing the gradients of the cost and constraint functions such as the variation method [41,42] and the costate method [27,43], the optimal parameter selection problem (TOCP ℎ  ) can be solved numerically with some mathematical programming algorithms.
In this paper, a hybrid optimization algorithm CPSO-SQP from [31,32], which combines SQP and CPSO, is introduced to solve the problem (TOCP ℎ  ).In this algorithm, the CPSO is taken as the global main optimizer of the algorithm while the SQP is employed for the local search to improve the solution.
In PSO, each solution is called a "particle" flying around in the  dimensional search space with a velocity   = (V 1 , V 2 , . . ., V  ).During flight, each particle of a swarm adjusts its position   = ( 1 ,  2 , . . .,   ) according to the best position encountered by itself (personal best position   = ( 1 ,  2 , . . .,   ) of its previous flights) and its neighbors (global best position   = ( 1 ,  2 , . . .,   ) of the population) as the following formulas: where  = 1, 2, . . ., ,  is the swarm size,  is the inertia weight,  Because of its ergodicity and stochasticity, the piecewise linear chaotic map is used to improve the global search capability of PSO, which is defined as follows: where  +1  denotes the chaotic mapping of th chaotic variable in the th iteration and  ∈ (0, 0.5) ∪ (0.5, 1) is a given parameter.Under this chaotic mapping,  behaves chaotically in (0, 1).Moreover, based on the fact that chaotic search is efficient in small range, we take chaotic mapping within a dynamically shrinking range around the global best position of the swarm by the following formula: where  = 1, 2, . . ., , x+1  denotes the chaotic mapping of the particle  in the th dimension,   is the global best position in the th dimension, and  is a scaling parameter with an empirical value 1.1.After chaos, the particle  should still be clamped to the range of [ min ,  max ] by the following formula: Integrated with chaotic map and SQP, the hybrid optimization algorithm CPSO-SQP is described as follows.
Step 1. Initialize swarm size, maximum function evaluations, other parameters, and the initial position   and velocity   of each particle in the swarm.
Step 2. Evaluate the objective function of each particle, then set the initial personal best solution   , select the global best particle   of CPSO to guide the swarm evolution direction, and set   =   which is the global best particle of CPSO-SQP.
Step 3. Update the velocity and position of each particle according to (40) if the particle is not in stagnation, otherwise, according to (41)- (43).
Step 4. Evaluate the objective function of each particle in the new position.
Step 5.If the current position is better than personal best solution, update   ; otherwise, set this particle in stagnation.
Step 6.If   is better than global best solution of CPSO, update   and set a flag for arousing SQP.Moreover, if   is better than   , update   .
Step 7. If   is just updated (the flag for SQP is true), set   as initial value and run SQP.If the result of SQP is better than   , update   .
Step 8.If maximum function evaluation is reached, output   as the global best solution; otherwise, go back to Step 3 for continuing iterations.
In this algorithm, the powerful global search capability of CPSO to search the solution space and the solution obtained from CPSO is taken as an initial estimate of SQP algorithm to compensate for the weakness of SQP weak global search.With good local convergence performance and strong nonlinear convergence speed, SQP algorithm makes a careful search for the optimization problem to make up shortcomings of the CPSO weak local search.Thus, the possibility of exploring a global minimum in problems with more local optima is increased.Benefiting from the fast globally converging characteristics of CPSO and the effective local search ability of SQP, the proposed method can obtain the global optimal result of our problem through constant iteration.

Numerical Simulation
In this section, we present some numerical simulation results by the proposed method in this paper.Take  =  =  =  = 1,  = 1, and  = 2. where positive constant  is the bound of the control constraint.
For the spatial discretization of the thermoelastic system, we use the basis functions   () for any 1 ≤  ≤  as follows: where  is the number of subintervals in the spatial domain.
Here, we choose  = 20, which means that the finite element triangulation T ℎ divides the space domain [0, ] into 20 subintervals.Through optimal time parameterization, the problem (TOCP ℎ ) has a fixed time horizon [0, 1].Using control parameterization, we subdivide the time interval [0, 1] into  subintervals with  = 20 and approximate the control functions by the piecewise-constant function in each subinterval [ −1 ,   ),  = 1, 2, . . ., .Then, the derived optimal parameter selection problem can be solved by the CPSO-SQP algorithm.
Our numerical simulation study was carried out within the MATLAB programming environment running on a personal computer with the following configuration: Intel Core i5-3470 3.20 GHz CPU, 4 GB RAM, 64-bit Windows 7 Operating System.
Figure 1 shows the relationship of the optimal time  * and the control bound .As the increase of , the shortest time  * , which is needed to drive the solution of the system to the energy constrained range, drops clearly.This means that the control ability of  * with large  is stronger than the that with small .In Figure 2, the computed optimal control is presented for  = 1, 2, 3, 4, which shows that the Bang-Bang

1 /Figure 1 :
Figure 1: The optimal time  * for different control bound .
−  −1 ),  −1 ≤  <   , − ( −  +1 ),   ≤  <  +1 , 0,otherwise, 1and  2 are positive constants called acceleration coefficients for regulating the maximum step size for the particle  to fly towards   and   , respectively, superscript  represents the algorithm iteration count, and rand() is a random value in the range [0, 1].  and   should be clamped to the range of [ min ,  max ] and [ min ,  max ] to control excessive roaming of particles outside the search space, respectively.