LARGE-SCALE
LINEAR PROGRAMMING
Proceedings of a I IASA Workshop, 2-6 June 1980
Volume 1
George B. Dantzig, M.A.H. Dernpster, and Markku Kallio
Editors
INTERNATIONAL INSTITUTE FOR APPLIED SYSTEMS ANALYSIS Laxenburg, Austria
1981
lnternational Standard Book Number 3-7045-0006-2
Volumes in the IIASA Collaborative Proceedings Series contain papers offered a t llASA professional meetings, and are designed to be issued promptly, with a minimum of editing and review.
The views or opinions expressed in this volume do not necessarily represent those of the Institute or the National Member Organizations that support it.
Copyright
O 1981.
International Institute for Applied Systems Analysis A-2361 Laxenburg, AustriaAll rights reserved. No part of this publication may be reproduced or transmitted in any form or by any means, electronic or mechanical, including photocopy, recording, or any information storage or retrieval system, without permission in writing from the publisher.
FOREWORD
The lnternational lnstitute for Applied Systems Analysis i s a nongovernmental, multi- disciplinary, international research institution whose goal i s t o bring together scientists from around the world t o work on problems of common interest.
IlASA pursues this goal, not only by pursuing a research program at the lnstitute i n col- laboration with many other institutions, but also by holding a wide variety of scientific and technical meetings. Often the interest in these meetings extends beyond the concerns of the participants, and proceedings are issued. Carefully edited and reviewed proceedings occasionally appear in the International Series on Applied Systems Analysis (published by John Wiley and Sons Limited, Chichester, England); edited proceedings appear in the IIASA Proceedings Series (published by Pergamon Press Limited, Oxford, England).
When relatively quick publication is desired, unedited and only lightly reviewed proceed- ings reproduced from manuscripts provided by the authors o f the papers appear in this new IIASA Collaborative Proceedings Series. Volumes in this series are available from the lnstitute at moderate cost.
PREFACE
During the week of June 2-6, 1980, the System and Decision Sciences Area of the Inter- national Institute for Applied Systems Analysis organized a workshop on large-scale linear programming in collaboration with the Systems Optimization Laboratory (SOL) of Stan- ford University, and cosponsored by the Mathematical Programming Society (MPS). The participants in the meeting were invited from amongst those who actively contribute t o research in large-scale linear programming methodology (including development of algo- rithms and software). Although primarily methodologically oriented scientists attended the workshop, i t s theme was the improvement of the long range applicability of linear pro- gramming (LP) techniques. Besides the exchange of ideas and experience
-
and sugges- tions for future research directions and international cooperation-
fostered by the meet- ing, it wasa general feeling of the participants that a proceedings would reflect the current state of large-scale linear programming in both East and West.To this end, it was considered important t o produce the proceedings volumes in a lecture note format as quickly as possible, so as t o secure a complete record of the papers presented at the workshop
-
including those destined for publication elsewhere-
together with several papers solicited by the editors in order t o extend coverage. I n some cases, papers presented at llASA have been revised by their authors i n the two months following the meeting; in others, no revisions have been made. Although a standard title page format has been used, the papers have been largely reproduced from camera-ready copy supplied bytheirauthors.Most have not been refereed, edited or proofread for typographical errors.Papers are grouped together in chapters by topic and are listed in alphabetical order by author in each cha~ter.
The first volume of these Proceedings contains five chapters. The first i s an historical review by George 6. Dantzig of his own and related research in time-staged linear program- ming problems. Chapter
2
contains five papers which address various techniques for exploit- ing sparsity and degeneracy in the now standard L U decomposition of the basis used with the simplex algorithm for standard (unstructured) problems. The six papers of Chapter 3 concern aspects of variants of the simplex method which take into account through basis factorization the specific block-angular structure of constraint matrices generated by dynamic and/or stochastic linear programs. By means of these techniques it i s hoped t o extend the size of solvable LP's beyond the range of current commercial codes for specific problems in the fields of energy, resource and macro/economic modeling (including eco- nomic planning models). I n Chapter 4, five papers address extensions of the original Dantzig-Wolfe procedure for utilizing the structure of planning problems by decomposing the original LP into LP subproblems coordinated by a relatively simple LP master problem of a certain type. Two of these papers concern the recent idea of applying this approach re- cursively t o the subproblems themselves. Chapter5
contains four papers which constitute a mini-symposium on the now famous Shor-Khachian ellipsoidal method applied t o both real and integer linear programs. This completes the description of the contents of Volume1.
The first chapter o f Volume
2
contains three papers on non-simplex methods for linear programming. This chapter concludes reports in the mainstream of current research on solution algorithms in large-scale linear programming. The remaining chapters of Volume2
concern more peripheral-
but no less important -topics of present interest i n the field.Techniques for exploiting network structure in LP problems are the topic of the three papers of Chapter
7.
In the next chapter, the emphasis turns t o the practically crucial and inter-related issues of automatic LP model generation and structure identification. The seven papers of this chapter discuss software both for model and matrix generation and for model reduction through detection of imbedded special constraint structure. The final chapter,9,
contains a number of applications of large-scale LP techniques t o practical problems in industrial and agricultural production and economic planning. Some of these involve multi-criteria optimization, and two of the eight papers deal explicitly with imple- mentations of new approaches t o the multi-criteria problem. A bibliography of large-scale linear programming research completes Volume2.
The editors wish t o take this opportunity on behalf of the participants t o thank
I
IASA, SOL and MPS for their cooperation and t o thank IlASA as well as various Academies of Sciences and governmental agencies of several countries for making the resources available t o hold the Large-scale Linear Programming Workshop and t o publish these Proceedings, I n particular, we are grateful t o the Communications Department at I IASA for their cheer- ful cooperation in expediting publication of this record of an important and memorable international meeting.George B. Dan tzig M.A. H. Dernpster Markku Kallio Stanford, California August
1980
CONTENTS VOLUME 1
Page
1. TIME-STAGED METHODS IN LINEAR PROGRAMMING: COMMENTS ANDEARLY HISTORY G. B. Dantzig.
2.
THE SIMPLEX METHOD FOR NONSTRUCTURED LINEAR PROGRAMSP. Huard. Solving large scale linear programs without structure.
A.F. Perold. Exploiting degeneracy in the simplex method.
A.F. Perold. A degeneracy exploiting LU factorization for the simplex method.
R.D. McBride. Controlling the size o f minikernels.
E. Toczylowski. Algorithms for block triangularization of basis matrices and exploi- tation of dual degeneracy in the dual simplex method.
3. THE SIMPLEX METHOD FOR DYNAMIC AND BLOCK-ANGULAR LINEAR PROGRAMS
T. Aonuma. A resourcedirective basis decomposition algorithm for weakly coupled dynamic linear programs.
M. Bastian. Aspects of basis factorization for block-angular systems with coupling rows.
R. Fourer. Solving staircase linear programs by the simplex method: 1. Inversion.
P. Gille and E. Loute. A basis factorization technique for staircase linear programs.
P. Kall. Computational methods for solving two-stage stochastic linear programming problems.
A. Propoi and V. Krivonozhko. The simplex method for dynamic linear programs.
4. DECOMPOSITION ALGORITHMS
P. Abrahamson. A nested decomposition approach for solving staircase-structured linear programs.
D. Ament, J. Ho, E. Loute, and M. Remmelswaal. L I FT: A nested decomposition algorithm for solving lower block triangular linear programs.
K. Beer. Solving linear programming problems by resource allocation methods.
J.K. Ho and E. Loute. An advanced implementation of the Dantzig-Wolfe decom- position algorithm for linear programming.
J. Stahl. On a decomposition procedure for doubly coupled LP's.
5. THE ELLIPSOIDAL ALGORITHM
M.L. Fisher. An equivalence between the subgradient and ellipsoidal algorithms. 479 P.E. Gill, W. Murray, M.A. Saunders, and M.H. Wright. A numerical investigation of
ellipsoid algorithms for large-scale linear programming. 48 7 M. Grotschel, L. Lovasz, and A. Schrijver. The ellipsoid method and its conse-
quences in combinatorial optimization. 51 1
S. Walukiewicz. Ellipsoidal algorithms for linear programming. 547
VOLUME 2
Page
6. OTHER ITERATIVE ALGORITHMSM.L. Fisher. The Lagrangian relaxation method for solving integer programming
problems. 58 1
E.G. Gol'shtein. An iterative linear programming algorithm based on the modified
Lagrangian. 61 7
M. Kallio and W. Orchard-Hays. Experiments with the reduced gradient method for
general and dynamic linear programming. 631
7. NETWORK PROBLEMS
J.K. Ho. A successive linear optimization approach to the dynamic traffic assign-
ment problem. 665
K. V. Kim, B.R. Frenkin, and B. V. Cherkassky. An efficient algorithm for updating
the basis in bicomponent linear problems. 687
U.H. Malkov, G.G. Padchin, and N.A. Sokolov. Some techniques to improve the effi-
ciency of solving linear programming problems. 699
8 . MODEL GENERATION AND STRUCTURE IDENTIFICATION
J. Bisschop and A. Meeraus. Toward successful modeling applications in a strategic planning environment.
G.G. Brown and D.S. Thomen. Automatic identification of generalized upper bounds in large scale optimization models.
G.G. Brown and W.G. Wright. Automatic identification of embedded structure in large-scale optimization models.
V.A. Bulavskiy. Hierarchical block-structure and factorization methods.
J.J.M. Evers and T.L. Knol. Constructing large linear input-output systems with recursively generated matrices.
G. Knolmayer. Computational experiments in the formulation of large-scale linear programs.
W. Orchard-Hays. Problems of symbology and recent experience.
9. APPLICATIONS AND MULTlCRlTERlA OPTIMIZATION
T. Aonuma. A decomposition-coordination approach to large-scale linear program- ming models: An aspect of applications of Aonuma's decomposition method.
I. Oeak, J. Ho ffer, J. Mayer, A. Nemeth, B. Potecz, A. Prekopa, and B. Strazicky.
Optimal daily scheduling of electricity production in Hungary. 923
2.
Harnos. A mathematical model for the determination of optimal crop productionstructures. 961
J.K. Ho. Holistic preference evaluation i n multiple criteria optimization. 97 7 M. Kallio, A. Lewandowski, and W. Orchard-Hays. An implementation of the refer-
ence point approach for multiobjective optimization. 1025
M. Kallio, A. Propoi, and R. Seppala. A model for the forest sector. 1055 A. T. Langeveld. Operational use of multiperiod LP models for planning and schedul-
ing. 1101
Y. Smeers. A column generationlnested decomposition algorithm for dynamic in-
putloutput models. 1113
BIBLIOGRAPHY ON LARGE-SCALE SYSTEMS Part
I.
General 1949-1 966Part
l I.
Classified 1949-1 980List of Participants
TIME-STAGED METHODS IN LINEAR PROGRAMMING: COMMENTS AND EARLY
HISTORY
TlhrlESTAGED METHODS I N LINEAR PROGRAMMING: COMMENTS AND EARLY HISTORY
George B. Dantzig
Department of Operations Research Stanford University
The Workshop on Large-scale Linear Programming reflects the active research taking place in many parts of the world along a very broad front, namely on:
the theory of solution, software development,
experiments on representative problems, application t o real problems,
matrix input generators, matrix analyzers, output report generators,
alternative methods of formulation.
This paper i s a historical review of the author's interest in one important facet of this field
-
the solution of time-staged programs. Indeed it was dynamic LP that initiated the linear programming field back in 1947. Over the years, many good ideas have been pro- posed, some that still merit serious consideration. This Workshop may provide the answer t o the question whether or not we have begun at last t o achieve the efficiency of solution necessary for successful application.This paper is a more polished version of the talk which I delivered opening the International Institute for Applied Systems Analysis Workshop on Large-Scale Linear Programming at Laxenburg Austria, June 2-6, 1 9 8 0 . Except for a short review of large-
scale methods also presented, but omitted here, my perspective is historical.
TIME-STAGED STAIRCASE SYSTEMS
The first formal papers about the new field of linear pro- gramming (that started in 19U7) appeared in Econometrica July
-
October 19U9. At the very beginning, the emphasis was on solving time-staged (dynamic) linear programs. That this is so, is clear from the following quote from [ I ] :
T h i s p a p e r i s c o n c e r n e d w i t h improved t e c h n i q u e s o f program p l a n n i n g , p a r t i c u l a r l y a s t h e y a p p l y t o t h e s c h e d u l i n g of a c t i v i t i e s o v e r t i m e w i t h i n a n o r g a n i z a t i o n o r economy i n which t h e a c t i v i t i e s must s h a r e i n t h e u s e of l i m i t e d amounts of v a r i o u s c o m o d i t i e s . The c o n t e m p l a t e d u s e o f e l e c t r o n i c c o m p u t e r s f o r r a p i d l y computing programs and t h e a s s u m p t i o n s u n d e r l y i n g t h e m a t h e m a t i c a l model a r e d i s c u s s e d . The p a p e r i s conclu{ed by a n i l l u s t r a t i v e example, [ B e r l i n A i r l i f t , -4 Time-Staged Dynamic L i n e a r Program].
The X a t h e m a t i c a l Xodel d i s c u s s e d h e r e i s a g e n e r a l i z a t i o n of t h e L e o n t i e f I n t e r - I n d u s t r y Model. I t i s c l o s e l y r e l a t r d t o t h e o n e found i n von Neumann's p a p e r "A Yodel of G e n e r a l Economic E q u i l i b r i u m " . I t s c h i e f p o i n t s of d i f f e r e n c e l : e i n i t s emphasis on dynamic, r a t h e r t h a n e q u i l i b r i u m o r s t e a d y s t a t e s . I t s p u r p o s e i s c l o s e c o n t r o l of a n o r g a n i z a t i o n - -
hence i t must be q u i t e d e t a i l e d ; i t i s d e s i g n e d t o h a n d l e h i g h l y dynamic problems--hence g r e a t e r s m p h a s i s on t i m e l a g s and c a p i t a l equipment; i t t a k e s i n t o c o n s i d e r a t i o n t h e many d i f f e r e n t ways o f d o i n g things--hence i t e x p l i c i t l y
i n t r o d u c e s a l t e r n a t i v e a c t i v i t i e s ; and i t r e c o g n i z e s t h a t any p a r t i c u l a r c h o i c e o f a dynamic program d e p e n d s on t h e
" o b j e c t i v e s " of t h e "economy", --hence t h e s e l e c t i o n and t y p e s of a c t i v i t i e s a r e made t o depend o n t h e m a x i m i z a t i o n of a n o b j e c t i v e f u n c t i o n .
In the companion paper [ 2 \ , the time staged staircase model is displayed and its relationship t o Leontief Input-Output model and continuous-time models is discussed:
where t h e x ( ~ ) a r e v e c t o r s of n o n n e g a t i v e e l e m e n t s .
.
When t h e m a t r i c e s and 9 ( t ) ( t = 1 , 2 , . . ,T) a r e s q u a r e and n o n s i n g u l a r , a d i r e c t s o l u t i o n i s p o s s i b l e t h a t may l e a d , however, t o n e g a t i v e and n o n n e g a t i v e a c t i v i t y l e v e l s ( i n which c a s e no f e a s i b l e s o l u t i o n e x i s t s ) .
I t s h o u l d b e n o t e d t h a t t h e g e n e r a l m a t h e m a t i c a l problem r e d u c e s i n t h e l i n e a r programming c a s e t o c o n s i d e r a t i o n of a s y s t e m of e q u a t i o n s of n o n n e g a t i v e v a r i a b l e s v h o s e m a t r i x of c o e f f i c i e n t s i s composed m o s t l y of b l o c k s of z e r o s e x c e p t f o r s u b m a t r i c e s a l o n g and j u s t o f f t h e " d i a g o n a l " . Thus any good c o m p u c a t i o n a l t e c h n i q u e f o r s o l v i n g programs would prob- a b l y t a k e a d v a n t a g e of t h i s f a c t .
Having f o m u l a t e d the time-staged model, it soon became clear that the techniques at hand at the time were inadequate. In a companion paper [ 3 1 , first presented in 1 9 4 9 , appeared the follow- ing statement:
Computing t e c h n i q u e s a r e nov a v a i l a b l e f o r s o l u t i o n of s m a l l l i n e a r programming problems. However, f o r a c c u r a t e o v e r - a l l A i r F o r c e p l a n n i n g , t h e s i z e o f t h e r e q u i r e d model i s s u c h
t h a t c o n v e n t i o n a l punched c a r d computing equipment, o r e v e n t h e i n t e r i m e l e c t r o n i c computer b e i n g b u i l t f o r t h e A i r F o r c e by t h e N a t i o n a l Bureau o f S t a n d a r d s , is n o t s u f f i c i e n t l y p o w e r f u l t o c o p e s a t i s f a c t o r i l y v i t h t h e problem of c h o o s i n g t h e optimum a c t i v i t i e s and a c t i v i t y l e v e l s o v e r t i m e .
I n o r d e r t o o b t a i n a programming p r o c e d u r e which would be i m m e d i a t e l y u s e f u l w i t h p r e s e n t l y a v a i l a b l e computing e q u i p - ment, we h a v e been forced t o u s e a d e t e r m i n a t e , and h e n c e l e s s g e n e r a l f o r m u l a t i o n of t h e programming problem t h a t p a r a l l e l s c l o s e l y t h e s c a f f p r o c e d u r e .
( ? j
Exogenous r _(t) _(t) -(t)
1
-
7-
3-
4t - l t = 2 c=3 t 3 4
I I n i t i a l 1 I
I
We have c a l l e d t h i s a zrianguZar ~ o d e i b e c a u s e i n i t t h e m a t r i x of d e t a c h e d c o e f f i c i e n t s , when a r r a g n e d a s i n t h e T a b l e . and o m i t t i n g t h e " i n i t i a l " p a r t , a s s u m e s a t r i a n - g u l a r f o r m , w i t h a l l c o e f f i c i e n t s a b o v e and t o t h e r i g h t of t h e p r i n c i p a l d i a g o n a l b e i n g z e r o . Thus t h e a c t i v i t i e s and i t e m s a r e s o o r d e r e d t h a t t h e l e v e l s o f a n y o n e a c t i v - i t y o v e r t i m e depend o n l y o n t h e l e v e l s o f t h e a c t i v i t i e s which p r e c e d e i t i n t h e h i e r a r c h y . T h i s means t h a t i n t h e c o m p u t a t i o n of t h e program we s u c c e s s i v e l y work down t h e h i e r a r c h y , a t e a c h s t e p s o l v i n g c o m p l e t e l y f o r t h e l e v e l s of e a c h a c t i v i t y i n e a c h of t h e t i m e p e r i o d s b e f o r e pro- c e e d i n g t o t h e n e x t a c t i v i t y ( s e e f i g u r e a b o v e ) .
The triangular model technique is a powerful empirical method when there is a natural hierarchy of activities and output items.
Certain energy models, for example, currently in vogue use such an approach.
BLOCK TRIANGULARITY
My p a p e r [ 4 1 , i s my f i r s t on methods f o r s o l v i n g l a r g e s y s - tems :
With t h e growing a w a r e n e s s of t h e p o t e n t i a l i t i e s of t h e l i n e a r programming a p p r o a c h t o b o t h dynamic and s t a t i c p r o b l e m s of i n d u s t r y , of t h e economy, and of t h e m i l i t a r y , t h e main o b s t a c l e toward f u l l a p p l i c a t i o n i s t h e i n a b i l i t y of c u r r e n t c o m p u t a t i o n a l methods t o c o p e w i t h t h e magnlt- ude of t h e t e c h n o l o g i c a l m a t r i c e s f o r even t h e s i m p l e s t s i t u a t i o n s . However, i n c e r t a i n c a s e s , s u c h a s t h e now c l a s s i c a l Hitchcock-Koopmans t r a n s p o r t a t i o n model, i t h a s been p o s s i b l e t o s o l v e t h e l i n e a r i n e q u a l i t y s y s t e m i n s p i t e of s i z e b e c a u s e of s i m p l e p r o p e r t i e s of t h e system.
T h i s s u g g e s t s t h a t c o n s i d e r a b l e r e s e a r c h be u n d e r t a k e n t o e x p l o i t c e r t a i n s p e c i a l matrix s t r u c t u r e s i n o r d e r t o f a c - i l i t a t e r e a d y s o l u t i o n of l a r g e r s y s t e m s .
I n d e e d , r e c e n t c o m p u t a t i o n a l e x p e r i e n c e h a s made i t c l e a r t h a t s t a n d a r d t e c h n i q u e s s u c h a s t h e s i m p l e x a l g o r i t h m , which h a v e been used t o s o l v e s u c c e s s f u l l y g e n e r a l s y s t e m s i n v o l v i n g o n e hundred e q u a t i o n s ( i n any r e a s o n a b l e number of n o n n e g a t i v e unknowns), a r e t o o t e d i o u s and l e n g t h y t o be p r a c t i c a l f o r e x t e n s i o n s much beyond t h i s f i g u r e . Our p u r p o s e h e r e w i l l b e t o d e v e l o p s h o r t - c u t c o m p u t a t i o n a l methods Eor s o l v i n g a n i m p o r t a n t c l a s s of s y s t e m s whose m a t r i c e s may b e g e n e r a l l y d e s c r i b e d a s " b l o c k t r i a n g u l a r " . By "block" t r i a n g u l a r we mean t h a t i f o n e p a r t i t i o n s t h e m a t r i x of c o e f f i c i e n t s of t h e t e c h n o l o g y m a t r i x i n t o sub- m a t r i c e s , t h e s u b m a t r i c e s ( o r b l o c k s ) c o n s i d e r e d a s e l e - ments form a t r i m g u l a r syscsm,
7 9
For example, von Neumann, i n c o n s i d e r i n g a c o n s t a n t l y ex- p a n d i n g economy, d e v e l o p e d a l i n e a r dynamic model whose m a t r i x of c o e f f i c i e n t s may b e w r i t t e n i n t h e form,
where A i s t h e s u b m a t r i x of c o e f f i c i e n t s of a c t i v i t i e s i n - i t i a t e d i n p e r i o d t , and B is t h e s u b m a t r i x of o u t p u t co- e f f i c i e n t s of t h e s e a c t i v i t i e s i n t h e f o l l o w i n g p e r i o d .
Now t h e main o b s t a c l e toward the f u l l a p p l i c a t i o n of s t a n - dard l i n e a r p r o g r a m i n g t e c h n i q u e s t o dynamic systems i s t h e magnitude of t h e m a t r i x f o r even t h e s i m p l e s t s i t u a t i o n s . For example, a t r i v i a l 15-activity--7-item s t a t i c model, when s e t up a s a 12-period dynamic model, would become a 1 8 0 - a c t i v i t y by 84-item system, which i s c o n s i d e r e d a l a r g e problem f o r a p p l i c a t i o n of t h e s t a n d a r d simplex method. A fancy model i n v o l v i n g , s a y , 200 a c t i v i t i e s and 100 items f o r a s t a t i c c a s e would become a 2000 x 1000 m a t r i x i f r e - c a s t a s a 10-period model. I t i s c l e a r t h a t dynamic models must be t r e a t e d w i t h s p e c i a l t o o l s i f any p r o g r e s s is t o be made toward s o l u t i o n s of t h e s e systems.
From a computational p o i n t of view, t h e r e a r e a number of observed c h a r a c t e r i s t i c s of t h e dynamic models which a r e o f t e n t r u e f o r s t a t i c models a s w e l l .
These a r e :
(1) The m a t r i x ( o r i t s t r a n s p o s e ) can be arranged i n t r i - a n g u l a r form
( 2 ) Most s u b m a t r i c e s A i j a r e e i t h e r z e r o m a t r i c e s o r com- posed of elements, most of which a r e zero.
( 3 ) A b a s i s f o r t h e simplex method i s o f t e n block t r i a n - g u l a r w i t h i t s d i a g o n a l submatrices s q u a r e and non- s i n g u l a r ( r e f e r r e d t o a s a "square b l o c k t r i a n g u l a r "
b a s i s ) .
( 4 ) For dynamic models s i m i l a r type a c t i v i t i e s a r e l i k e l y t o p e r s i s t i n t h e b a s i s f o r s e v e r a l p e r i o d s .
To i l l u s t r a t e , c o n s i d e r a dynamic v e r s i o n of t h e Leontief model i n which ( a ) a l t e r n a t i v e a c t i v i t i e s a r e p e r m i t t e d
(a s i m p l e c a s e would be where s t e e l can be o b t a i n e d from d i r e c t p r o d u c t i o n o r s t o r a g e ) ; ( b ) i n p u t s t o an a c t i v i t y f o r p r o d u c t i o n i n t h e t t h time period may occur i n t h e same o r e a r l i e r time p e r i o d s . It can be shown i n t h i s model t h a t
( a ) a b a s i c s o l u t i o n w i l l have e x a c t l y m a c t i v i t i e s i n each time period (where n
-
number of time dependent e q u a t i o n s ) . ( b ) each s h i f t i n b a s i s w i l l b r i n g i n a s u b s t i t u t e a c t i v i t y i n t h e same time p e r i o d , and ( c ) o p t i m i z a t i o n can be c a r r i e d o u t a s a sequence of one-period o p t i m i z a t i o n problems; i . e . , t h e optimum c h o i c e of a c t i v i t i e s ( b u t n o t t h e i r amounts) can be determined f o r t h e f i r s t time period (independent of t h e l a t e r p e r i o d s ) t h i s permits a d e t e r m i n a t i o n f o r t h e second time period (independent of t h e l a t e r p e r i o d s ) , e t c e t e r a . When flow models a r e r e p l a c e d with more complex models whichi n c l u d e i n i t i a l i n v e n t o r i e s , c a p a c i t i e s , and t h e b u i l d i n g of new c a p a c i t i e s , t h e i d e a l s t r u c t u r e of a b a s i s ( s e e t h i r d c h a r a c t e r i s t i c above) no longer holds. However, t e s t s ( c a r - r i e d on s i n c e 1950) on a number of c a s e s i n d i c a t e t h a t bases.
w h i l e o f t e n not square btock t r C a p l a r i n the s m s e above, could be nude so by changing r e l a t i ~ ~ e l y few n o l m s i n the b c s i s (e.g.. one o r two a c t i v i t i e s i n s m a l l models). T h i s c h a r a c t e r i s t i c of near-square t Z o c i t r i a n g u l a r i t y of t h e b a s i s , i . e . . w i t h n o n s i n g u l a r s q u a r e submatrices down t h e d i a g o n a l , i s , of c o u r s e , c o m p u t a t i o n a l l y convenient and t h i s paper w i l l be concerned with ways t o e x p l o i t i t .
Towards the end of the above paper can be found the following:
F i n a l l y , may I make a s h o r t p l e a t h a t l i n e a r programmers pay g r e a t e r a t t e n t i o n t o s p e c i a l methods t o r s o l v i n g t h e l a r g e r m a t r i c e s t h a t a r e encountered i n p r a c t i c e . The es- c e l l e n t work of Jacobs on t h e c a t e r e r problem and t h e work of J a c o b s , Hoffman, Johnson on t h e production smoothing problem a r e examples of what may be done with c e r t a i n dyn- amic models w i t h a simple r e p e t i t i v e s t r u c t u r e . Cooper and Charnes have employed i n t h e i r work a number of s h o r t c u t s t h a t have p e r m i t t e d r e s o l u t i o n of c e r t a i n l a r g e s c a l e sys- tems. A t RAND we have found e f f i c i e n t ways t o hand compute g e n e r a l i z e d t r a n s p o r t a t i o n problems, and Markowicz has pro- posed a g e n e r a l procedure i n t h i s a r e a t h a t i s promising.
Many models e x h i b i t a block t r i a n g u l a r s t r u c t u r e and c e r - t a i n p a r t i t i o n i n g methods have been proposed which t a k e advantage of t h i s type of s t r u c t u r e . There i s need f o r t h o s e of you who a r e f o r e s i g h t e d t o do s e r i o u s r e s e a r c h i n t h i s a r e a .
A t t h e p r e s e n t time (1955), i t i s p o s s i b l e t o s o l v e r a p i d l y problems i n t h e o r d e r of a hundred e q u a t i o n s . The Orchard- Hays 701 Simplex Code has solved many problems of t h i s s i z e w i t h a s high a s 1,500 unknowns and machine times of f i v e t o e i g h t hours a s a r u l e - - a l l w i t h e x c e l l e n t s t a n d a r d s of ac- c u r a c y . However, i t is s e l f - e v i d e n t t h a t no m a t t e r how much t h e g e n e r a l purpose codes a r e p e r f e c t e d they w i l l be unable t o cope w i t h t h e n e x t g e n e r a t i o n of problems which w i l l be l a r g e r i n s i z e . I t i s a l s o e v i d e n t t h a t t h e models c u r r e n t l y being run could have been handled more e f f e c t i v e l y by t h e proposed s p e c i a l methods.
There a r e c e r t a i n c h a r a c t e r i s t i c s common t o many models which I b e l i e v e should be emphasized:
(1) Host f a c t o r s i n t h e c o e f f i c i e n t m a t r i x a r e zero.
( 2 ) I n dynamic s t r u c t u r e s t h e c o e f f i c i e n t s a r e o f t e n t h e same from one time period t o t h e n e x t .
( 3 ) In dynamic s o l u t i o n s t h e a c t i v i t i e s employed o f t e n p e r s i s t from one period t o t h e next.
(O) T r a n s p o r t a t i o n t y p e s u b m a t r i c e s a r e common.
( 5 ) Block t r i a n g u l a r submatrices a r e common.
P a r t of t h e r e s e a r c h i n t h i s a r e a should c e r t a i n l y be de- voted t o a b e t t e r understanding of t h e p o t e n t i a l i t i e s of t e c h n i q u e s o t h e r than t h e simplex method.
UNCERTAINTY
In a related paper [ 5 1 , published in 1 9 5 6 , appears the following I n t h e p a s t few months t h e r e have been important developments t h a t p o i n t t o t h e s p p i i c a t i s n of l i n e a r progruwning neshods xnder uncertainty. By way of background l e t us r e c a l l t h a t t h e r e a r e i n connnon use two e s s e n t i a l l y d i f f e r e n t t y p e s of sched- u l i n g applications--one designed f o r t h e s h o r t run and t h o s e
f o r t h e l o n g r u n . For t h e l a t t e r t h e e f f e c t of p r o b a b i l i s t i c o r c h a n c e e v e n t s i s r e d u c e d t o a minimum, by t h e u s u a l t e c h - n i q u e of p r o v i d i n g p l e n t y of fat i n t h e system. For example, conswrrprion razes, a t t r t z i o n r a i e s , wecrr-~ut r a t e s a r e a l l planned o n t h e h i g h s i d e . 7 h e s t o ship, time t o t r a v e l , tz;7les zo produce a r e a l w a y s made w e l l above a c t u a l needs.
I n d e e d , t h e e n t i r e s y s t e m i s p u t t o g e t h e r w i t h p l e n t y of slack and fat w i t h t h e hope t h a t t h e y w i l l be t h e shock d s o r b e r s which w i l l p e r m i t t h e g e n e r a l o b j e c t i v e s and tim- i n g of t h e p l a n t o be e x e c u t e d i n s p i t e of u n f o r e s e e n e v e n t s . I n t h e g e n e r a l c o u r s e of t h i n g s , long-range p l a n s a r e r e - v i s e d f r e q u e n t l y b e c a u s e t h e s t o c h a s t i c s e l e m e n t s of t h e problem have a n a s t y way of i n t r u d i n g . For t h i s r e a s o n a l s o t h e c h i e f c o n t r i b u t i o n , i f a n y , of t h e long-range p l a n , is t o e f f e c t a n i n m e d i a t e d e c i s i o n - s u c h a s t h e a p p r o p i a t i o n of f u n d s o r t h e i n i t i a t i o n of a n i m p o r t a n t development con- t r a c t .
For s h o r t - r u n s c h e d u l i n g , many of t h e s l u k and fat t e c h - n i q u e s of i t s l o n g - r a n g e b r o t h e r a r e employed. The p r i n c i - p l e d i f f e r e n c e s a r e a t t e n t i o n t o d e t a i l and t h e s h o r t time- h o r i z o n . As l o n g a s c a p u b i i i t i e s a r e w e l l a b o v e require- ments ( o r demands) o r i f t h e demands c a n b e s h i f t e d i n t i m e , t h i s a p p r o a c h p r e s e n t s no problems s i n c e i t i s f e a s i b l e t o implement t h e s c h e d u l e i n d e t a i l . However, where t h e r e a r e s h o r t a g e s , t h e p r o j e c t e d p l a n based on s u c h t e c h n i q u e s may l e a d t o a c t i o n s f a r from o p t i m a l , whereas t h e s e new methods, where a p p l i c a b l e , may r e s u l t i n c o n s i d e r a b l e s a v i n g s . I s h a l l s u b s t a n t i a t e t h i s l a t e r by r e f e r e n c e t o a problem of A. F e r g u s o n o n t h e r o u t i n g of a i r c r a f t .
With r e g a r d t o t h e p o s s i b i l i t i e s of s o l v i n g l a r g e s c a l e l i n - e a r p r o g r a m i n g problems, one c a n sound b o t h a n o p t i m i s c i c and a p e s s i r n i a t i c ~ n o t e . The p e s s i m i s t i c n o t e c o n c e r n s t h e a b i l i t y of t h e problem f o r m u l a c o r , e i t h e r a m a t e u r o r p r o f e s - s i o n a l , t o d e v e l o p models t h a t a r e l a r g e s c a l e . The p e s s i - m i s t i c n o t e a l s o c o n c e r n s t h e i n a b i l i t y of t h e problem s o v l e r
t o compute models by general t e c h i p e s when t h e y a r e l a r g e s c a l e . I f t h i s i s s o , i s n o t t h e g r e a t p r o m i s e t h a t t h e l i n - e a r programming a p p r o a c h w i l l s o l v e s c h e d u l i n g and l o n g r a n g e p l a n n i n g problems w i t h s u b s t a n t i a l s a v i n g s t o t h e o r g a n i z a t i o n s a d o p t i n g t h e s e methods b u t a n i l l u s i o n and a s n a r e ? Are t h e b i g problems g o i n g t o be s o l v e d a s t h e y have a l w a y s been solved--by a d e t a i l e d s y s t e m of on-the-spot somewhat n a t u r a l s e t of p r i o r i t i e s t h a t r e s o l v e e v e r y p o s s i b l e a l t e r n a t i v e a s i t a r i s e s ?
The status of problems involving uncertainty as far as prac- tical solutions are concerned, has not changed much since 1956.
The following, sums up the 1965 situation:
When o n e c o n s i d e r s i n s t e a d , a d i r e c t a t t a c k on u n c e r t a i n t y v i a m a t h e m a t i c a l programming, i t i n e v i t a b l y l e a d s t o t h e con-
s i d e r a t i o n of large-scale s y s t m s . Problems w i t h t h e i r s t r u c - t u r e , have p r o v e n d i f f i c u l t of s o l u t i o n s o f a r . I b e l i e v e t h a t t h e y w i l l b e t h e s u b j e c t of i n t e n s i v e i n v e s t i g a t i o n i n t h e f u t u r e .
DECOMPOSITION PRINCIPLE
The Decomposition Principle [ 61 arose in 1958 in connection with a military tactical problem which was too large to handle by conventional linear programming problem. A good summary o f the approach c a n be found in my 1965 survey article:
R e c e n t l y t h e a u t h o r , j o i n t l y w i t h P h i l i p Wolfe, developed a new p r o c e d u r e t h a t i s p a r t i c u l a r l y a p p l i c a b l e t o a n g u l a r s y s t e m s and m u l t i s t a g e s y s t e m s of t h e s t a i r c a s e t y p e T h i s i s r e p o r t e d i n p r e l i m i n a r y form i n RAND P-1544 (Nov.10, 1958) under t h e t i t l e , "A Decomposition P r i n c i p l e f o r L i n e a r Programs?. The s y s t e m c o n s i s t s of c e r t a i n goods s h a r e d i n common among s e v e r a l p a r t s and c e r t a i n goods ( i n c l u d i n g f a c - i l i t i e s , raw m a t e r i a l s ) p e c u l i a r t o e a c h p a r t . I n s h o r t t h e s y s t e m i s a n g u l a r i n s t r u c t u r e .
Although t h e e n t i r e p r o c e d u r e i s one i n t e n d e d t o b e c a r r i e d o u t i n t e r n a l l y i n a n e l e c t r o n i c computer i t may a l s o b e viewed a s a decentratized decision making process. Each indepen- d e n t p a r t i n i t i a l l y o f f e r s a p o s s i b l e b i l l o f goods ( a vec- t o r o f t h e comon o u t p u t s and s u p p o r t i n g i n p u t s i n c l u d i n g o u t s i d e c o s t s ) t o a c e n t r a l c o o r d i n a t i n g agency. A s a s e t t h e s e a r e m u t u a l l y f e a s i b l e v i t h e a c h o t h e r and t h e g i v e n common r e s o u r c e s and demands from o u t s i d e t h e system. The c o o r d i n a t o r works o u t a s y s t e m of " p r i c e s " f o r paying f o r e a c h component of t h e v e c t o r p l u s a s p e c i a l s u b s i d y f o r e a c h p a r t t h a t j u s t b a l a n c e s t h e c o s t .
The management of each p a r t t h e n o f f e r s , based on t h e s e p r i c e s , a nev f e a s i b l e program f o r h i s p a r t w i t h lower c o s t ' r ) i ~ h o u t regard t o w;zetner it i s ;^easibZe for t h e sgszem a s a w k t e . The c o o r d i n a t o r , however, combines t h e s e new o f f e r s v i t h t h e s e t of e a r l i e r o f f e r s s o a s t o p r e s e r v e m u t u a l f e a - s i b i l i t y and c o n s i s t e n c y v i t h exogeneous demand and s u p p l y and t o minimize c o s t . Using t h e improved o v e r - a l l s o l u t i o n h e g e n e r a t e s a r e v i s e d s e t of p r i c e s , s u b s i d i e s , and r e c e i v e s new o f f e r s . The e s s e n t i a l i d e a i s t h a t o l d o f f e r s a r e n e v e r f o r g o t t e n by t h e c e n t r a l agency ( u n l e s s u s i n g " c u r r e n t "
p r i c e s t h e y a r e u n p r o f i t a b l e ) ; t h e f o w e r a r e mixed w i t h t h e new o f f e r s t o form new p r i c e s .
In the original paper [61 appears this abstract:
A t e c h n i q u e i s p r e s e n t e d f o r t h e d e c o m p o s i t i o n of a l i n e a r program t h a t p e r m i t s t h e problem t o be s o l v e d by a l t e r n a t e s o l u t i o n s of l i n e a r sub-programs r e p r e s e n t i n g i t s s e v e r a l p a r t s and a c o o r d i n a t i n g program t h a t i s o b t a i n e d from t h e p a r t s by l i n e a r t r a n s f o r m a t i o n s . The c o o r d i n a t i n g program g e n e r a t e s a t e a c h c y c l e new o b j e c t i v e forms f o r e a c h p a r t , and e a c h p a r t g e n e r a t e s i n t u r n (from i t s o p t i m a l b a s i c f e a - s i b l e s o l u t i o n s ) new a c t i v i t i e s (columns) f o r t h e i n t e r c o n - n e c t i n g program. Viewed a s a n i n s t a n c e of a ' g e n e r a l i z e d programming problem' whose columns a r e drawn f r e e l y from g i v e n convex s e t s . s u c h a problem c a n b e s t u d i e d by a n ap- p r o p r i a t e g e n e r a l i z a t i o n of t h e d u a L i t y theorem f o r l i n e a r
p r o g r a m i n g , which p e r m i t s a s h a r p d i s t i n c t i o n t o b e made b e t v e e n t h o s e c o n s t r a i n t s t h a t p e r t a i n o n l y t o a p a r t of t h e problem and t h o s e t h a t c o n n e c t i t s p a r t s . T h i s l e a d s t o a g e n e r a l i z a t i o n of t h e Simplex Algorithm, f o r v h i c h t h e de- c o m p o s i t i o n p r o c e d u r e becomes a s p e c i a l c a s e .
The reported experience with solving structured linear pro- grams by means of the decomposition principle varies from very good to poor, In general it appears that if the decomposition between master and sub is a "naturaln one, it can perform very well. Like the simplex method, there is rapid improvement for the early iterations followed by a long tail except here the tail is much longer.
COMPACT BASIS INVERSES
From 1 9 6 2 onwards there has been growing interest in schemes for compactly representing the inverse of the basis for the simplex method. This effort goes under various names: compact basis tri- angularization, LU basis factorization. One must worry not only about the compactness but also about the stability of the solution to small changes in the original data. My 1 9 6 2 paper [ 7 1 was dir- ected to finding a compact representation of a basis for staircase systems.
Alex Orden was t h e f i r s t t o p o i n t o u t t h a t t h e i n v e r s e of t h e b a s i s i n t h e s i m p l e x method s e r v e s no f u n c t i o n e x c e p t a s a means f o r o b t a i n i n g t h e r e p r e s e n t a t i o n of t h e v e c t o r e r t e r i n g t h e b a s i s and f o r d e t e r m i n i n g t h e new p r i c e v e c t o r . For t h i s p u r p o s e one of t h e many forms of " s u b s t i t u t e i n - v e r s e s " ( s u c h a s t h e w e l l known p r o d u c t form of t h e i n v e r s e ) v o u l d do j u s t a s w e l l and i n f a c t may have c e r t a i n advan-
t a g e s i n c o m p u t a t i o n .
Harry Xarkowitz was i n t e r e s t e d i n d e v e l o p i n g , f o r a s p a r s e m a t r i x , a s u b s t i t u t e i n v e r s e w i t h a s few nonzero e n t r i e s a s p o s s i b l e . He s u g g e s t e d s e v e r a l ways t o do t h i s a p p r o x i m a t e l y . For example, t h e b a s i s c o u l d b e reduced t o t r i a n g u l a r form by s u c c e s s i v e l y s e l e c t i n g f o r p i v o t p o s i t i o n t h a t row and column whose p r o d u c t of nonzero e n t r i e s ( e x c l u d i n g t h e p i v o t ) is minimum. He a l s o p o i n t e d o u t t h a t , f o r b a s e s whose non- z e r o s a p p e a r i n a band on a s t a i r c a s e a b o u t t h e d i a g n o n a l . p r o p e r s e l e c t i o n of p i v o t s c o u l d r e s u l t i n a compact sub- s t i t u t e w i t h no more n o n z e r o s t h a n t h e o r i g i n a l b a s i s . We s h a l l a d o p t M a r k m i t z ' s s u g g e s t i o n . However, i n s t e a d of r e c o r d i n g t h e s u c c e s s i v e t r a n s f o r m a t i o n s of one b a s i s t o t h e n e x t i n p r o d u c t form. we s h a l l show t h a t i t i s e f f i c i e n t t o g e n e r a t e e a c h s u b s t i t u t e i n v e r s e i n t u r n from i t s predeces- s o r . The s u b s t i t u t e i n v e r s e remains compact. of f i x e d s i z e . Thus " r e i n v e r s i o n s " a r e u n n e c e s s a r y ( e x c e p t i n s o f a r a s t h e y a r e needed t o r e s t o r e l o s s of a c c u r a c y d u e t o cumula- t i v e round-off e r r o r ) .
The p r o c e d u r e which we s h a l l g i v e c a n b e a p p l i e d t o a gen- e r a l s x m b a s i s w i t h o u t s p e c i a l s t r u c t u r e . A s s u c h , i t i s
p r o b a b l y c o m p e t i t i v e w i t h t h e s t a n d a r d p r o d u c t form, f o r i t may have a l l of i t s a d v a n t a g e s and none of i t s d i s a d v a n t a g e s . With c e r t a i n m a t r i x s t r u c t u r e s , moreover, i t a p p e a r s t o b e p a r t i c u l a r l y a t t r a c t i v e .
We s h a l l f o c u s o u r remarks o n szaircase strucedes. The r e a d e r w i l l f i n d no d i f f i c u l t y i n f i n d i n g a n e q u a l l y e f f i - c i e n t way t o compact b l o c k - a n g u l a r s t r u c t u r e s .
STATUS AS OF 1967
A s u m m a r y o f t h e s t a t u s o f s o l v i n g l a r g e - a c a l e p r o b l e m s c a n be f o u n d i n my 1967 p a p e r [ 8 ]
.
From i t s v e r y i n c e p t i o n , i t was e n v i s i o n e d t h a t l i n e a r pro- gramming would b e a p p l i e d t o v e r y l a r g e . d e t a i l e d models of economic and m i l i t a r y s y s t e m s . K a n t o r o v i t c h ' s 1939 propos- a l s , which w e r e b e f o r e t h e a d v e n t o f t h e e l e c t r o n i c computer, mentioned s u c h p o s s i b i l i t i e s . L i n e a r programming e v o l v e d o u t of t h e U.S. A i r F o r c e i n t e r e s t i n 1947 i n f i n d i n g o p t i m a l t i m e - s t a g e d deployment p l a n s i n c a s e o f war; a problem whose m a t h e m a t i c a l s t r u c t u r e is s i m i l a r t o t h a t o f f i n d i n g a n op- timal g r o w t h p a t t e r n o f a d e v e l o p i n g economy and s i m i l a r t o o t h e r c o n t r o l problems. S t r u c t u r a l l y t h e dynamic p r o b l e m s a r e c h a r a c t e r i z e d i n d i s c r e t e form by s t a i r c a s e m a t r i c e s r e p r e s e n t i n g t h e i c p u t s and o u t p u t s from o n e t i m e p e r i o d t o t h e n e x t . T r e a t e d a s a n o r d i n a r y l i n e a r program, t h e number of rows and columns grows i n p r o p o r t i o n t o t h e number of t i m e p e r i o d s T and t h e c o m p u t a t i o n a l e f f o r t grows by T~ and p o s s i b l y h i g h e r . T h i s f a c t h a s l i m i t e d t h e u s e of l i n e a r programming a s a t o o l f o r p l a n n i n g o v e r many t i m e p e r i o d s . At t h e p r e s e n t 1967 s t a g e of t h e computer r e v o l u t i o n , t h e r e is growing i n t e r e s t on t h e p a r t of p r a c t i c a l u s e r s o f l i n e a r programming models t o s o l v e l a r g e r and l a r g e r s y s t e m s . Such a p p l i c a t i o n s imply t h a t e v e n t u a l l y automated s y s t e m s w i l l o b t a i n information from c o u n t e r s and s e n s i n g d e v i c e s , pro- c e s s d a t a i n t o t h e p r o p e r form f o r o p t i m i z a t i o n and f i n a l l y implement t h e r e s u l t s by c o n t r o l d e v i c e s . T h e r e h a s been s t e a d y p r o g r e s s i n t h i s m e c h a n i z a t i o n of f l o w t o and from t h e computer. H i t h e r t o , t h i s h a s been one of t h e o b s t a c l e s e n c o u n t e r e d i n s e t t i n g - u p and s o l v i n g l a r g e - s c a l e s y s t e m s . The s e c o n d o b s t a c l e h a s been t h e c o s t and t h e t i m e r e q u i r e d t o s u c c e s s f u l l y s o l v e l a r g e problems.
It is d i f f i c u l t t o measure t h e p o t e n t i a l o f l a r g e - s c a l e p l a n n i n g . C e r t a i n d e v e l o p i n g c o u n t r i e s a p p e a r , a c c o r d i n g t o o p t i m a l c a l c u l a t i o n s on s i m p l i f i e d models t o be a b l e t o grow a t t h e r a t e o f 15% p e r y e a r i m p l y i n g a d o u b l i n g of t h e i r i n d u s t r i a l b a s e i n f i v e y e a r s . But a d m i n i s t r a t o r s a p p a r e n t l y i g n o r e p l a n s and make d e c i s i o n s based on p o l i t - i c a l e x p e d i e n c y which r e s t r i c t growth t o 2 o r 3% o r some- t i m e s -2%. It i s t h e b e l i e f o f t h e a u t h o r t h a t t h e mech- a n i z a t i o n of d a t a f l o w ( a t l e a s t i n advanced c o u n t r i e s ) i n t h e n e x t d e c a d e w i l l p r o v i d e pathways f o r c o n s t r u c t i n g
l a r g e models and t h e e f f e c t i v e implementation of t h e r e s u l t s of o p t i m i z a t i o n . T h i s a p p l i c a t i o n of mathematics t o d e c i s i o n p r o c e s s e s w i l l e v e n t u a l l y become a s i m p o r t a n t a s t h e c l a s s i c a l a p p l i c a t i o n s t o p h y s i c s and w i l l , i n time, change t h e emph- a s i s i n p u r e mathematics.
In this paper the following unsolved problem was posed:
I t h a s been d i s c o v e r e d r e c e n t l y t h a t t h e s i z e of t h e i n v e r s e r e p r e s e n t a t i o n of t h e b a s i s i n t h e simplex method c o u l d have a n i m p o r t a n t e f f e c t on r u n n i n g time. T h e r e f o r e , compact- i n v e r s e schemes a l o n g t h e l i n e s f i r s t proposed by Harry Markovitz of RAND have become i n c r e a s i n g l y i m p o r t a n t . Re- c e n t l y , two groups working i n d e p e n d e n t l y , developed t h i s a p p r o a c h w i t h a s t o u n d i n g r e s u l t s . For example, t h e S t a n d a r d O i l Company of C a l i f o r n i a group r e p o r t s r u n n i n g - t h e on some of t h e i r t y p i c a l l a r g e problems c u t t o 114.
How t o f i n d t h e most compact i n v e r s e r e p r e s e n t a t i o n of a s p a r s e m a t r i x i s s t i l l a n unsolved problem:
CONJECTURE: I f a non-singuLzr ,mtri.z has K non-zero zLe.rrmts, it i s always possible t o represent them as G pro-
&cC of e Z z . ~ m t a r g ma$rices such chat she f c r a l nunber of non-zero e n s r i z s ;ezcLuding t h e f r ii- aooral u n i t elements) i s a s mosr i(. [Inc*Lden~ally, -he -ml;ir.:zni schemes just mentioned o f f e n have no more than K+IC"ii non-zeros in t h e %verse rs- presentation.
1
STATUS TO THE PRESENT ( 1 9 8 0 )
From 1 9 6 7 onwards there has been an increasing interest in techniques for solving large-scale linear programs. A number of conferences have been exclusively concerned with the topic. Most general operations research and management science meetings have at least one session devoted to it. A selected reference list which I use in my seminars (mostly published during the period 1 9 7 0 - 7 8 ) contain 2 3 7 titles which I have arranged by sub area.
General Books
(10 e x c l u s i v e l y l a r g e s c a l e , 2 s p a r s e methods, 8 o t h e r ) Survey articles
GUB, G-GUB and the decompositioh principle Variants of above
Block Triangularity
Linear optimal control and dynamic systems Nested decomposition
Column generation, convex and nonlinear programs Sparse matrix techniques
Large networks and related problems Applications
Software
Total
Some i d e a o f t h e r e c e n t r e s e a r c h o f t h e S y s t e m s O p t i m i z a t i o n L a b o r a t o r y o f t h e O p e r a t i o n s R e s e a r c h D e p a r t m e n t a t S t a n f o r d c a n b e g l e a n e d f r o m t h e t i t l e s t h a t f o l l o w :
A n d r e P e r o l d : " F u n d a m e n t a l s o f a C o n t i n u o u s Time S i m p l e x M e t h o d " .
A n d r e P e r o l d a n d G e o r g e B. D a n t z i g : "A B a s i s F a c t o r - i z a t i o n Method f o r B l o c k T r i a n g u l a r L i n e a r P r o g r a m s " . Bob F o u r e r : " S o l v i n g S t a i r c a s e - s t r u c t u r e d L i n e a r
P r o g r a m s by A d a p t a t i o n o f t h e S i m p l e x M e t h o d " . Ron D a v i s : "New Jump C o n d i t i o n s f o r S t a t e C o n s t r a i n e d
O p t i m a l C o n t r o l P r o b l e m s " .
P h i l i p Abrahamson a n d G e o r g e B. D a n t z i g : "Imbedded D u a l D e c o m p o s i t i o n A p p r o a c h t o S t a i r c a s e S y s t e m s " . J o h n B i r g e : " S o l v i n g S t a i r c a s e S y s t e m s u n d e r U n c e r t a i n t y " . T h i s Workshop may w e l l mark t h e p o i n t i n t i m e when e f f i c i e n t m e t h o d s f o r s o l v i n g l a r g e d y n a m i c s y s t e m s may b e more t h a n j u s t a p r o m i s e . T h i r t y t h r e e y e a r s f r o m t h e t i m e t h e h o p e was f i r s t e x - p r e s s e d t h a t s u c h m e t h o d s b e f o u n d , t h e y may s o o n become a r e a l i t y :
REFERENCES
[I] Wood, Marshall K., and Dantzig, George B., Programming of Interdependent Activities, I, Z c o n o r n e t r i c a , July-October 1949.
[2] Dantzig, George B., Programming of Interdependent Activities, 11, Z c o n o r n e t r i c a , July-October 1949.
[3] Wood, M.K., and Geister, M.A., Development of Dynamic Modehs for Program Planning, d c ; ; v C z ~ Aca Zy s d a 3 f F r r 2 : i ~ ' i ~ q n J Z Z .-I l ' o c z s i o r , , T. Koopmans, Ed., \Jlley and Co. 1951.
[4] Dantzig, George B., Upper Bounds, Secondary Constraints and Block Triangularity in Linear Programming, Econometrics, 23, April 1955, pp 174-183.
[5] Dantzig, George B., Recent Advances in Linear Programming, Management S c i e n c e Vol. 2, No. 2, January 1956, pp 131- 144.
[6] Dantzig, George B., and wolfe,Philip, Decomposition Principle for Linear Programs, O ? e r a t i o n s E'esearch 8, Jan-Feb.1960, pp 101-111.
Another version of this paper appeared in Z c o n o n e t r i c - 29, Oct. 1960, pp 767-778.
[7] Dantzig, George B., Compact Basis Triangularization for the Simplex Method, ,?ecens A d v a n c e s i n 'Yo:hernaticai Frogram- m i n g , Graves, R.L., and Wolfe, P., Eds., McGraw Hill, 1963, N.Y. pp 125-133.
[81 Dantzig, George B., Large-Scale Systems and the Computer Revolution, Proceedings of the Princeton Symposium on Mathematical Programming, H.W. Kuhn, Ed. Princeton University Press, Princeton N.J. Aug. 1967 pp 51-72.
T H E SIMPLEX METHOD FOR
NONSTRUCTURED LINEAR PROGRAMS
SOLVING LARGE SCALE LINEAR PROGRAMS WITHOUT STRUCTURE
P.
HuardDirection des Etudes et Recherches Electricite de France
A variant of the simplex method is adapted for the solution of large-size linear program- ming problems with a very sparse constraint matrix. Instead of using the inverse of the basis, three sparse linear systems are directly solved at each step, using a suitable pivoting method. Two advantages of this variant compared t o standard procedure are:
Memory volume requirements are proportional to the number of constraints (and not t o i t s square).
Calculation may be faster; the appropriate numerical tests are described in the paper.
1 .
-
IBTRODUCTIONWith regard to the resolution of large linear programs, the basis of a variant of the Simplex method, using only a snall amount of memory, has already been briefly described C31.
The aim of the present paper is to give a detailed study of this method and of the numerical experiments that validate it.
In its classiral form, the Simplex method uses a square matrix, the inverse of the basic matrix. whose value is updated at each iteration.
The number of nonzero elements of this matrix increases rapidly as the iterations go along and it is necessary in practice, when using the explicit form of the inverse, to have on hand a number of memories equal to the square of its dimension, say m2 for a linear program with m constraints. Thus it becomes difficult to handle problems having several hundred constraints, without using disks or tapes; then the overhead time may becomes prohibitive, because of their repetitive use and the large number of iterations.
Some special structures of the matrix of the linear program
-
likefor example the block-angular one
-
allow for various interesting decompositions of the inverse of the basic matrix, which is similar to the solving of smaller linear programming problems. Then the amount of ncccssary mr.mory varies only lincari;. with t!ie size of the program, if the dirnensic>n of t h e biocks is a co:i-,~snt. Fortuaately, such a block-angular structure is rather often encountered (dynamic problems, regionalization problems) anti various decomposition methods have been proposed (see for example [ 5 ] ) .
However, many linear programs do not have any structures suitable for decomposition. This is the case for problems related to a graph
-
e.g. flow-problems
-
which contain the problem of electrical dispatching, as far as its structure is concerned.Large linear programs, issued from "real life", have a very sparse matrix : only a few percent of the elements are nonzero. Of course, this sparsity appears in each basic matrix, but it disappears from the inverse matrix. The variant of the Simplex method, which follows, uses the basic matrix itself, instead of its irmerse, and then eliminates the need for mL memory positions. However, in the calculations, products of a matrix by a vector are replaced by
resolutions of linear systems of the same dimensionality. The complexity of these two operations would be of m2 and m3 order respectively, if the matrices were full, which would rule out the proposed variant. But, as will be seen below, two factors may make it competitive. One is the difference in sparsity betveen the basic matrix and its inverse. The other is the fact that generally, the basic aatrix is almost triangular, or more precisely
"triangular-band-wise". In other words : after having performed a suitable permutation of rows and col,mns, nonzero elements lie bclow an extra- diagonal line, located at a small distance p above the diagonal. Such a linear system is easily solved through a specialized pivoting method that we call below the method of parameters. The amount of calculations is proportional to p p m2, where p is the proportion of nonzero elements, p the width of the band located above the diagonal, and m the dimension of the matrix-(a large number, by hypothesis). In large problems, of real origin, that we have known of, p is often between p m and 2 p m. If p '
is the proportion of nonzero elements (density) of the inverse matrix is normally much larger than p), the respective amounts of computatjon for one iteration of the Simplex method are roughly in the ratio 4(p/p') m. 2 For p '
-
60 p and m =lo3,
this is practically 1. In actual fact, numerical comparisons of Section 7, involving linear programs of up to 900 constraints, exhibit a very good speed for the proposed variant. In Section 8 the detailed costs for one iteration of the Simplex method are given with a comparison between the two variants.- 2 2 -
2.
-
TIE RREQUIRED CALCULATIONS D U R I N G ONE ITERATION OF THE SIW'LEX METHODThe l i n e a r program t o b e s o l v e d i s g i v e n i n s t a n d a r d form
-
Maximise f x s u b j e c t t o
A x P a
X ' P
where A i s a f u l l - r a n k m a t r i x ; i t s rows a r e i n d e x e d by M = ( I , 2,
. . .,
rn]and i t s columns by N = ( 1 , 2 ,
...,
n l .A t e a c h i t e r a t i o n , a b a s i s I i s c o n s i d e r e d , i . e . a s u b s e t I s u c h t h a t : I c N
Ill =
A' i n v e r t i b l e
where ',A t h e b a s i c m a t r i x r e l a t i v e t o I , i s composed o f t h e columns
A', j c I.
To t h e b a s i s I i s a s s o c i a t e d t h e s o - c a l l e d b a s i c s o l u t i o n of t h e b a s i s I , d e f i n e d by
x- = 0 I
where
7
i s t h e complement o f I i n N .The s l ~ c c e s s i v e b a s e s g e n e r a t e d by t h e Simplex method, a r e s u c h t h a t xI Z 0 ; h e n c e , t h e c o n s i d e r e d b a s i c s o l u t i o n s a r e a l l f e a s i b l e ( t h e y s a t i s f y c o n d i t i o n s ( I ) and ( 2 ) ) .
An iteration consists of changing the basis I inco a neighboring basis 1', that is a basis obtained by exchanging an index r € I with an index s E
i
:To determine r and s, one can compute, in order :
-
Iwhere u, fl, d are row-vectors. This allows the candidate s E
y,
to bechosen with the condition dS > 0. Then :
where xI, a, TS and AS are column vectors. This gives r c I by the condition
Once r and s are determined, it remains to update the inverse of the basic matrix, i.e. to compute (A'')-'. This is clasiicaly done from
(A1)-' through the relation :
where E is an elementary matrix, explicitly known (see figure 1)
Figure 1
Thus t h e n e c e s s a r y c a l c u l a t i o n s a r e r e p r e s e n t e d by r e l a t i o n s ( 4 ) t o (9). and t h e i n v e r s e of t h e b a s i c m a t r i x i s u s e d i n ( 4 ) . ( 6 ) . ( 7 ) . These l a s t r e l a t i o n s c a n b e r e p l a c e d by
i . e . t h r e e l i n e a r s y s t e m s t o s o l v e . I n t h e f i r s t one, t h e m a t r i x i s t h e t r a n s p o s e of t h e b a s i c m a t r i x , i n t h e l a s t two, i t i s t h e b a s i c m a t r i x i t s e l f : t h e s e s y s t e m s e n j o y t h e s p a r s i t y o f t h e A m a t r i x , and s o l v i n g them c a n be done w i t h o u t s t o r i n g and u s i n g t h e i n v e r s e .
3.
-
DIRECT RESOLUTION OF THE LINEAR SYSTEMThe systems (4'), (6'). (7') have long been successfully solved directly in the case of classical transportation problems. These very special linear programs can be stated :
. .
Minimize Z cLJ xij subject to i j
Z
.
x.. 1 = b. 1,
i = 1.2,...,
qJ
x 2 0
,
V i j i jHere the A matrix has no more than 2 nonzero elements per column.
which are equal to I , and the basic matrices are triangular. Thus solving the three linear systems is particularly easy and fast (it is not even necessary, here, to solve (7')).
An extension to problems of flow with gains was proposed by MAURR4S [4]
in 1972. In this type of linear programs, the A matrix still has no more than 2 nonzero elements per column, but of any real value. Systems (4').
(6') or (7') are almost as simple as a triangular system. The method of solution consists of particularizing one unknown as a parameter, and in expressing one after the other the (ml) remaining unknowns as functions of this parameter, using (wl) equations. Eliminating these (m-I) unknowns from the last equation
-
not yet used-
gives the value of the parameter.Plugging this value in the expression of the (m-I) unknowns completes the solution. The choice of the particularized unknown is guided by an interpretation of the structure of the A matrix, as incidence matrix of a graph. Of course, it is not possible to extend this theory to matrices with more than 2 nonzero elements per column. However, a study of many square matrices, very large and very sparse, issued from real problems.
shows that they often have a triangular-band-wise structure (after suitable permutations of rows and columns); their band-width has the samc order of magnitude as the average number of nonzero elements per column or pcr row. Elorc .precisely, these square matrices are such that
