7.74 1 976 SLCO1 3374 fooms) DRC- 18 Matrix Modifications and the Solution of Large Linear Discrete-Time Control Problems with Time-Integral Constraints and Upper Bounds by SECrORAL LIBRAI K *1 1,4TERN OA ~,, v* RECOSTRUCTION, ASND DEVEI-OPi' Johannes Bisschop and Alexander Meeraus AUG21 Technicil Note No. 5 Reseaich Project 670-24 October 1976 Development Research Center World Bank 1818 H Street, N.W. Washington, D.C. 20433 57.74 . B475 1976 SLCO13374 Matrix Modifications and the Solution of Large Linear Discrete-Time Control Problems with Time-Integral Constraints and Upper Bounds by Johannes J. Bisschop and Alexander Meeraus 1. Introduction. Discrete-time control models have been formulated and solved in the field of economics. An excellent summary is provided in [ 3 ] . Especially the class of trajectory tracking problems formulated as linear-quadratic control models has found applications in economics. This is partly because of the computational feasibility of solving large problems using the explicit solution technique based on the theory of Ricatti equations. This technique solves problemp with a quadratic objective function and linear dynamics. No additional constraints on either the state or the control vari- ables, or both, are allowed. If in the optimal solution of a tracking pro- blem any undesired values of the state and/or control variables are observed, one usually changes the objective function using weights such that these undesirable values disappear. In this technical note we will first study the discrete-time optimal control problem with linear dynamics and linear objective function allowing for simple upper bounds and time-integral constraints. We will then extend this problem to incorporate additional constraints on the control and state variables for each time period. Time-dynamic economic planning models with resource and/or balance of payments constraints are potential -2- members of these classes of problems. 2. The Problem Structure. The following notation will be used. x t denotes an n-dimensional state vector; ut denotes an m-dimensional control vector; s t denotes an n-dimensional vector of constants; t denotes time and varies from 0 to T ; A , B , C , and D are matrices of sizes nxn , nxm , pxn , and pzn respectively; at ,bt , c , d , and w are vectors of sizes n ,m n , m and p respectively We will consider the problem T-1 Minimize I c x + du t=O s.t. xt+1 = A xt + B ut +z t T-1 I C xtl + Du = w t=0 0 < x < at 0 < u < b (1) -t-t- t- t For notational convenience we will consider only the above problem, but the methodology in this paper can be adapted to handle slightly more general problems. Some of the difference equations may be taplicit rather than explicit which will change I x t+l to G xt+1 . Here G is usually nonsingular and deviates from the identity matrix by only a few rows and columns. Another variant is to make the constant matrices A , B , C , D , c and d time-varying. Usually they obey some simple functional form, and not all matrices need to be stored explicitly. The initial simplex tableau corresponding to the system (1) takes on the form x1 x2 x3 x4 . 0 ul u2 u3 I -B -A I -B -A I -B (2) -A I -B . C C C D D D D C Note that this tableau can be soored compactly by storing A , B , C and D columnwise, and using a set of pointers to these columns to represent (2) . The following theorem will identify a structure that is common to any nonsinglar submatrix of selected columns of (2) . Theorem 1 : Any (Tn+p) x (Tn+p) basis extracted from the above tableau (2) can always be arranged in a bordered block lower triangular structure where each diagonal block is of size nxn , and contains itself two nonsingular diagonal blocks I i and Bi ,i1,2,...,T -4- I B 11 10 I B 12 B2 A, 0 0 B2 0 0 B2 13 3 0 A 0 0 B3 2 03 01 1 C2 {2 C3 f3 01 2 3 Proof: Take any basis. Initially, consider all rows to be in the order as they are found in the tableau (2) , and let all basis columns corresponding to elements of xt and ut-1 be assigned to block t , t-1,2,...,T . The resulting basis has a block angular structure where the number of columns in each diagonal block is greater than or equal to the number of rows. By proper permutation of rows and columns in each block, the upper left portion will contain an identity matrix It of size (n-gt) x (n-gt) where (n-gt) is the number of xt elements that are basic variables, 1 < g < m From here on, consider the permuted matrix. Let t-l . The row rank of the first diagonal block is n as the entire matrix is assumed to be nonsingular. If I has n columns, the remaining columns, if any, corresponding to elements of u0 , are transferred to the very right of the matrix. The resulting nxn diagonal block is nonsingular. If Il has less than n columns, say n-g1, then the nonsingularity of the basis matrix implies that the last g rows of the columns corresponding to the basic elements of u0 are of rank gl , and contain a glxg1 nonsingular submatrix, say B1 . Adjoin the columns for which the last gl rows make up this nonsingular matrix B1 , to the columns corresponding to the basic elements of x1 . The resulting nxn diagonal block is nonsingular. All remaining columns, if any, corresponding to' u0 are transferred to the very right of the matrix. Next assume the theorem is true for t=i-l> 1 . Then the set of n basis rows corresponding to t=E is linearly independent, even if the previously transferred columns corresponding to U0, u1 ..., u-2 are ignored. As all remaining columns have zero elements in the first (f-1)n rows, it must be that the ^th diagonal block contains a nxn nonsingular matrix. The previous arguments used for the case t=l do now apply. This concludes the proof of the Theorem. The above Theorem 1 can be used to find an inverse representation for any basis B . Let B be partitioned as follows. B M P 1 P2 B = [, : P3 P4 where P1 is the TnxTn block lower triangular matrix as identified -1 -1 by the theorem. Then B can be represented using P1 and -1 -1 -1 -1 Q M (P4 - P3 P1 2) Note that P itself can be represented using - 6- -1 only the inverse matrices Bi , and that the size of Q is not dependent on the size of T . The system P1 P2 xl1 bl1 P 3 P4 x2 b 2 can then be solved via the equations -1 - x Q1 (b -P P-1b 2 2 3 1 b1 -1 xl P1 (b I 2 x2 ) 3. Modifications of the Basis Matrix The inverse representation is based on Theorem 1 of the previous section, and can be a mixture of full and sparse matrix technology. Note that these results are closely related to those derived for LP problems with GUB or block-diagonal.structures. This observation will conceptually simplify further extensions of problem (1) , and will add to our previous work [ 2 ] . Using the revised simplex method, one needs to update the inverse basis. This can be done using the process of augmentation as is described in [ 1 ] . After k iterations, the augmented basis matrix i can be partitioned as follows. P1 P2 U1I SP3 4 U 2 V 1 V 2 -I where the diagonal blocks are of sizes TnxTn , pxp , and IxI , respectively, - 7- and I < k . Here, every row of [V1 V2I is a unit vector with a one underneath a leaving column of-B and zeros everywhere else. Every column of consists of the difference between an incoming vector and a leaving U2 vector of B, and is partitioned. There is a series of representations for the inverse of a 3x3 blocked.matrix, assuming the inverses of certain transformed submatrices -1 -1 -1 exist. For the above system construct the matrices P-1 D Q1 v where Q (P P -1 S (P4 -P3 1 P2 and A -1 Q L I- Vl2][1 2 81 3 4 2 Using these matrices, there are several ways to compute the solution of the problem -1 2 UA 1 b1 P3 4 U2 2 ii- b 2 V 1 V2 -1 i3 b 3 (a) If only P1 QI and Q are available, then one may compute X3 -Q (b3 1 2 1 [U] r3 P4 U 2 -8- and is computed via the equations - -- -1 * - - - - Q (b2 3 P b 1) 1 - 1 (b P2 2 Note that the above series of computations requires one to access the matrices Q once, Q- twice, and P-1 four times. -1 1 (b) If in addition to p1 9 Q , and Q one has also stored the matrix P 1 P2 ' then the series of computations can be rearranged such - 1 1 b1 3 4 2 b U 2 -11- andtw acese o P1 b2e 3quired -l-1 (c) If in addition to P1 1 , and Q , one has also stored the matrix Pl 2,teth seisocopttoscnbreragduh [-1 -1] 1: -9- then the series of computations can be rearranged such that one access -1 -1 of Q , and two accesses of P1 are required. (d) If in addition to P 1 , , and Q , one has also stored --the matrices -1 and Pl P , then the series of computations can be rearranged such that one ^1-1 -1 access of Q , one access of Q1 and one access of P are required. -1 -1l- (e) If in addition to P1 Q1 , and Q , one has also stored the marce -1 P-ndP1 matrices P 2 and P1 Ul, then the series of computations can be ^-1 -1 rearranged such that one access of Q , two accesses of Q , and one access of P are required. This is done in the following manner x2 Qb 2P 3 P 1b1 1 3 1 1 and 1 1 b1 Pi '2 Pi 1 Ul x2 x3 where Q . P -P Pm P U - P P U 4 3 1 2 2 3 1 1 v2 1 1 1 2 -1 v1 1 1 - 10 - -1^ and Y = Q b is computed via the equations ^-1 ^- -1 1 Y = (b2- (U2- P3 P U1) Q b ) 2 Q 2 2 P3 P1 1) 1 b1)- -1 6 -1 1 1 (1 2 1 1 2 2 The desirability of any specific representation is determined by storage requirements and computational times involved. - 11 - 4. Added Structure to the Problem Problem (1) in section 2 can be extended to incorporate more .general bounds on the control and/or state vectors. Let the matrices E, -F, H1 and H2 be of sizes exn, fxm, hxn and hxm, respectively. Then the tableau (2) can be extended as follows: xj x_ x3 xU U U3 1 2 x 3 x4 U o U1 U2 _3 I -B -A I -B -A I -B -A I -B E E (3) E E F F F F H1 H 2 H1 H 2 H H 2 HH H2 C C C C D D D D - 12 - In many control problems one may assume that e + f + h < < n. The following theorem can be used to identify a structure that is common to any nonsingular submatrix of selected columns of (3). Theorem 2: Any (T(n + e + f + H) + P) x (T(n + e + f + h) + P) basis extracted from the above table (3) can always be arranged in a bordered block lower triangular structure where each diagonal block Di is of size (n + e + f + h), and takes on the form i i 0 B J. D = i i 3 4 5 J J J i i i 5 where Ji is of size (e + f + h) x (e + f + h). Proof: The proof is essentially that given for Theorem 2. To find the appropriate vectors to make up the third set of columns in the above block, one needs to make sure that Q = (1[ j4 ] - iii i i i 0 B J is a nonsingular matrix. This requires at least (e + f + h) accesses of I -1 -1 and Bi , and the inverted matrix Qi is needed in the overall inverse representation of any diagonal block Di* - 13 - 5. Conclusion By taking advantage of the explicit structure of the-problem provided via the algebraic statement, one can build a basis inverse repre- sentation founded on this explicit structure. By doing so, one can store the overall representation efficiently. Combining this with our updating technique this methodology effectively results in the updating of a reduced working basis. The size of this working basis is PxP and independent of the number of time periods. - 14 - References [1 Bisschop, J.,and A. Meeraus, "Efficient Updating of the Basis Inverse in Linear Programming via Partitioning", DRC Technical Note #2, 1976. [2 Bisschop, J.,and A. Meeraus, "Partitioning, Generalized Upper Bounding, and Angular Structures in Linear Programming", DRC Technical Note # 3, 1976. (3 Kendrick, D., "Applications of Controle Theory to Macroeconomics", to be published in Annals of Economic and Social Measurement, 1976.