NPS55-82-012 NAVAL POSTGRADUATE SCHOOL Monterey, California SOLVING GENERALIZED NETWORKS by Gerald G. Brown and Richard D. McBride March 1982 Approved for public release; distribution unlimited Prepared for: M — 1 Postgraduate School FEDDOCS )re *> California 93943 D 208.14/2 NPS-55-82-012 P< ZZ-OIZ NAVAL POSTGRADUATE SCHUUL MONTEREY, CALIFORNIA Commodore K. H. Shumaker uavid A - ^chrady Superintendent Provost Reproduction of all or part of this report is authorized. UNCLASSIFIED SECURITY CLASSIFICATION OF THIS PAGE (Whtn Data Entered) SkR£V CA 93943-5101 REPORT DOCUMENTATION PAGE READ INSTRUCTIONS BEFORE COMPLETING FORM 1. REPORT NUMBER NPS55-82-012 2. GOVT ACCESSION NO 3. RECIPIENT'S CATALOG NUMBER 4. TITLE (and Subtitle) SOLVING GENERALIZED NETWORKS 5. TYPE OF REPORT & PERIOD COVERED Technical Report March 19 3 2 6. PERFORMING ORG. REPORT NUMBER 7. AUTHORS Gerald G. Brown Richard D. McBride 8. CONTRACT OR GRANT NUMBERS 9. PERFORMING ORGANIZATION NAME AND ADDRESS Naval Postgraduate School Monterey, CA 93943 10. PROGRAM ELEMENT, PROJECT, TASK AREA & WORK UNIT NUMBERS II. CONTROLLING OFFICE NAME AND ADDRESS Naval Postgraduate School Monterey, CA 93943 12. REPORT DATE March 1982 13. NUMBER OF PAGES 60 14. MONITORING AGENCY NAME & ADDRESS^' different from Controlling Office) 15. SECURITY CLASS, (of thlt report) Unclassified 15«. DECLASSIFI CATION/ DOWN GRADING SCHEDULE 16. DISTRIBUTION STATEMENT (ol this Report) Approved for public release; distribution unlimited. 17. DISTRIBUTION STATEMENT (of the mbatract entered In Block 20, II different from Report) 18. SUPPLEMENTARY NOTES 19. KEY WORDS (Continue on reverse elde if neceeamry and Identity by block number) Generalized Networks, Flow Networks, Minimum Cost Flow Networks, Generalized Transshipment Problem, Primal Simplex Network Optimization 20. ABSTRACT (Continue on reverae aide It neceeamry and Identity by block number) A complete, unified description is given of the design, implemen- tation and use of a family of very fast and efficient large scale minimum-cost (primal simplex) network programs. The class of capacitated generalized transshipment problems solved includes the capacitated and uncapacitated generalized transportation problems and the continuous generalized assignment problem, as well as the pure network flow models which are specializations of these DD , ^ N RM 73 1473 EDITION OF 1 NOV 65 IS OBSOLETE S/N 0102- LF-014-6601 UNCLASSIFIED SECURITY CLASSIFICATION OF THIS PAGE (Whan Data Bnfred) UNCLASSIFIED SECURITY CLASSIFICATION OF THIS PAGE (Whan Dmtm Bnftud) problems. These formulations are used for a large number of diverse applications to determine how (or at what rate) flows through the arcs of a network can minimize total shipment costs, A generalized network problem can also be viewed as a linear program with at most two non-zero entries in each column of the constraint matrix; this property is exploited in the mathematical presentation with special emphasis on data structures for basis representation, basis manipulation, and pricing mechanisms. A literature review accompanies computational testing of promising ideas, and extensive experimentation is reported which has pro- duced GENNET, an extremely efficient family of generalized network systems. S'N 0102- LF- 014-6601 UNCLASSIFIED SECURITY CLASSIFICATION OF THIS PAGEfWhen Datm Entafd) SOLVING GENERALIZED NETWORKS by Gerald G. brown Naval Postgraduate School Monterey, CA 93943 USA and Richard D. McBride University of Southern California Los Angeles, CA 9UUU" 7 USA August 1981 (Rev 84.U2.U1) (Copyright ly84) DUDLEY KM ■ ABSTRACT A complete, unified description is given ot the design, implementation and use of a family ot very fast and etticient larye scale minimum-cost (primal simplex) network proyrams. The class ot capacitated generalized transshipment problems solved includes the capacitated and uncapaci tated generalized transportation problems and the continuous generalized assign- ment problem, as well as the pure network flow models which are specializations ot these problems. These formulations are used for a large number of diverse applications to determine how (or at what rate) flows through the arcs of a network can minimize total shipment costs. A generalized network problem can also be viewed as a linear program with at most two non-zero entries in each column ot the constraint matrix; this property is ex- ploited in the mathematical presentation with special emphasis on data structures for basis representation, basis manipulation, and pricing mechanisms. a literature review accompanies compu- tational testing ot promising ideas, and extensive experimenta- tion is reported which has produced GfciMNhiT, an extremely effi- cient family of generalized network systems. 1 . INTRODUCTION This paper reports the development ot a larye-scaie primal network code tor solving capacitated generalized transshipment problems. The capacitated generalized transshipment problem is the most general of the minimum cost tlow models in continuous variables, whicn include the capacitated and uncapaci tated trans- portation problems and the continuous generalized assignment problem as well as the pure network specializations ot these problems. These models are used tor a large number ot diverse applications that include transportation ot goods, desiyn ot reservoir, com- munications, and pipeline systems, assignment ot personnel ana machinery to jobs, bid evaluation; currency exchange and cash management; production, sales and inventory planning; and many others. For further discussion of these applications see survey articles such as Bradley [bj, or Glover et al., [17, 18], and textbooks (as well as their cited references) such as Dantzig [13], Jensen and Barnes [23], and Kennington and Helgason [27], The capacitated generalized transshipment model and its specializations are minimum-cost network flow problems. The goal is to determine how (or at what rate) flows through the arcs ot a network can minimize shipment costs. The network is a directed graph, G , defined by a set ot nodes, N , and a set ot arcs, A , with ordered pairs ot nodes (tail, head) as elements indexed by k: ( t , t, ) . for each arc there is a shipping cost per unit tlow, c , a minimum allowable tlow (or lower bound), I, , and a maximum allowable tlow (or upper bound, or capacity), u k In addition, there are coefficients (or multipliers, or gains, or losses), a and b, which can change the magnitude of (or — K K amplify, or attenuate) each unit of flow respectively entering and leaving arc k . Each node is either a supply node where units of the good enter the network, a demand node where units leave, or a transshipment node. The problem is to minimize total costs with flows, x, , that satisfy the associated lower and upper bounds and preserve the conservation of flow at each node : (GN) MIN I c, x keA k k s . t . i keA with f, =i k a. x. — k k keA with t. =i x, = b i e N \ < x R < u R , keA, where : b i ■ supply if i is a supply node; - demand it i is a demand node; U otherwise. Note that any linear program with at most two nonzero coefficients associated with each variable is a generalized network (GN). for- mulation (Gin) is a generalized transshipment model (GT) if all a. are +1, in which case the corresponding coefficient is called the multiplier m, = -b. . For purposes of exposition, we will address (GT) since assuming the existence of a finite upper bound on each variable it is possible to transform the (GN) coet t icients--by scaling or by reflecting variables with respect to their upper bounds--so that one coefficient is +1 tor each variable. The other coefficient for each variable then becomes the multiplier value, m. , and the generalized transshipment (GT) network flow interpre- tation results with a node tor each constraint and a directed arc tor each variable. If a variable has two nonzero coet t icients , its arc is directed away from the node corresponding to the constraint with the +1 coefficient; a variable with just one nonzero coefficient (±1) corresponds to an arc forming a self-cycle, leaving and returning to the same node. Some generalized networks (GN), such as those obtained by relaxing integer restrictions on flow variables, cannot be scaled conveniently to (GT). These (GN) are accommodated by obvious minor modifications in the following (GT) presentation. We have developed codes tor solving (GN) and specialized them tor ( GT ) problems. Generalized networks can be solved as linear programming problems, but contemporary commercial linear programming sys- tems consume much more computer time and data storage region than special purpose network codes. Indeed, the advantage of network codes is so pronounced that it is even worthwhile to develop special linear programming procedures to exploit intrin- sic network structure found embedded within more general models (e.g., Kennington and Helgason [27J, McBride [3U], Brown and Graves [8], and Brown and McBride [9]). Brown and Wright [11] and Brown, McBride and Wood [10] show that many real-life linear programs contain a large embedded generalized network structure. 3 These models are widely used because they accurately de- scribe a large variety of important applications. Generalized networks not only directly represent gains with m, > 1 (e.g. interest return on investment, heat gain, etc.) and losses with < m, < 1 (e.g., evaporation, voltage drop, attrition, etc.) rK in the flows, but also admit conversions of units for these flows (e.g., machine time to output pieces, lira to yen, etc.). In addition, m, < represents situations without obvious physical K flow interpretation (e.g., flows which "enter the head" of arc k in proportion m. of those entering the tail), but which nonethe- less provide valuable modelling tools. There has been continuing growth of interest in network models because efficient computer programs have made possible the reliable, economic solution of problems with more variables than virtually any other optimiza- tion technique (e.g., the pure network system, GNET [6], has been installed at hundreds of sites worldwide and is now cited as a routine research tool). Perhaps most important, networks are readily accepted by nonanalysts and are consequently extremely popular operations research models. Although several papers have been written in this general area, and significant computational breakthroughs have been reported, there has not previously been a single, unified de- scription of a complete implementation, nor have "new generation" computer programs been made generally available to the academic community. Here we report the research and computational experi- ments which have produced GfcJNNET, an extremely efficient family of network optimization systems. GENNET exploits pure network structure embedded in generalized networks, and specializes in- trinsically to GNET lb], when both systems are applied to pure networks, (the floating point arithmetic in) GENNET requires about lb percent more time than GNET. An important objective ot this paper is to make these new approaches easily accessible to a wide audience via a clear exposition and concrete examples of efficient FORTRAN programs. Further, the availability ot the ( GT ) and (GN) computer programs will now make it possible tor other investigators to reproduce and extend our experimental results . Bradley, Brown and Graves lb] trace the historical develop- ments leading to contemporary primal simplex pure network algo- rithms, their supporting data structures, and etticient imple- mentations such as GNET. For other sub-classes ot generalized networks, algorithms have been reported by Jewell [2b], Eisemann [14], Maurras [2y], Glover, Hultz, Klingman and Stutz [17,la,ly], Balachandran [2], and Jensen and Bhaumik 124]. An efficient algorithm for large generalized network problems has been de- veloped by Glover, Klingman, Hultz, Stutz, Karney and Elam [15,17,18,21,22]. However, their contributions are scattered among the papers referenced; Kennington and Helgason [27] and Jensen and Barnes [23] provide textbook descriptions of compu- tations apparently gleaned from these papers, providing an ade- quate treatment of the graphical algorithm with more computational advice than the seminal presentation by Dantzig [13] but tew details ot etticient basis updating. 2. THE APPROACH Our approach continues with generalized networks in pre- cisely the philosophical vein of the pure network exposition of Bradley, Brown and Graves: we seek data structures and algo- rithms that yield efficient implementations without abandoning the flexibility of a general large-scale mathematical program- ming perspective [6, p. 3 ff]. We introduce few of the details of the general bounded-variable simplex algorithm, and we repeat little of the underlying pure network material; the assiduous reader might well review the prior paper for which this is an intimate companion. We continue with a brief description of the algebraic specialization of the simplex method for generalized networks. Specific design decisions and experiments carried out with GENNET are described, including computational tests of alternate ap- proaches. Some extensions of GENNET are presented to further exploit special problem structure. J. GENERALIZED NETWORK SPECIALIZATION Efficient primal simplex specialization to the generalized network case depends upon the well-known result that any generalized network basis can be put in neirly (upper) tri- angular form by simple permutation of rows and columns. This inherent near triangularity can be exploited by direct solution of the simplex equations with modif ied forward, or back substi- tution. Fortuitously, this basis structure also leads to ex- tremely fast solution updates orchestrated in concert with efficient dynamic reorganization of each new basis. Theorem (e.g., Uantzig [13]) Any basis B extracted from a generalized network problem can be put in the form (1) by rearranging rows and columns . ,P (1) where each square submatrix component B is either upper triangular or nearly upper triangular with only one element below the diagonal . Proof. Our proof is constructive. We know that the basis has at least one nonzero entry in each row, and one, or two, entries in each column. Define a pairing as the associ- ation of a row with a column sharing an entry in b , a deferral as a temporary deletion of a pair from consider- ation, and an ass ignment as the fixing of a pair on the diagonal in the ordering of the rearranged sequence, followed by the assignment of all deterred pairs which have a column with an entry in an assigned row. Step 1) Defer singleton rows . Locate a row with one entry, pair the row with the column of the entry and defer the pair. Repeat Step 1 until no row remains with one entry. Comment: Step 1 reduces B to a submatrix with exactly two entries in each row and in each column . To see this, note that each column can have at most two entries, or 2m entries for m columns. Each row has at least two entries. Suppose that some row has more than two entries; then at least 2m + 1 entries exist, leading to a contradiction. Thus, exactly 2m entries remain; consequently, each column has exactly two entries, as does each row. Comment: Step 1 deters an upper triangular set of pairs. To see this, diagonalize the pairs in reverse of the order of their deferral, and place this sequence at the end of the rearrangement of B . Comment: A sequence of assignments in the rearrangement begins. step 2) Assign a Component . 2.1) If a row remains which is neither assigned nor deterred, begin a near triangular component: pair the row and a column and ass ign this pair next in the rearranged sequence. utherwise, go to Step 2.3. 2.2) Apply Step 1 and ass ign each newly deterred pair next in the sequence. Go to Step 2.4. 2.3) begin a triangular component: assign the most recently deferred pair next in the sequence. 2.4) Repeat Step 2 until all rows and columns are assigned. Comment: It Step 2.1 assigns a pair to a component, then the component is nearly triangular with one element below the diagonal in the first column. utherwise, the component is triangular. Consider the generalized network problem given in Figure 1, a generalized transshipment problem with 2 sources, 1U sinks, 15 nodes, and 3U arcs. toe will use this problem to illustrate concepts and etticient solution methods. In this problem the multipliers are all positive. tor arc k directed trom node i to node j , if tlow x, leaves node i , then ^i^k arrives at node j . When m > 1 the amount arriving at node j is greater than the amount leaving node i . This would be the case in cash flow problems when the arc corresponds to an investment ana m, = (1 + r. ) with r_, the rate of return. When m, < 1 then the amount arriving at node ] is less than — k the amount leaving node i . here the loss could correspond to evaporation, taxation, transmission loss, brokerage fees, seepage or deterioration. Pure network arcs are indicated by m. = 1 — k (and may be even more abundant in real problems than in our example,). for clarity, minimum allowable flows are zero, and maximum allowable flows are not specified. Arc from To Cost Multiplier k f k fc k C k ^k 1 4 3 33.84 .99 Node Supply 2 3 2 1 3 5 15.47 53.54 1.00 .74 1 22.86 4 2 5 26.76 .74 2 177.14 5 3 5 73.49 1.00 6 5 5 52.52 1.00 Node Demand 7 3 6 35.12 .91 8 5 6 11 .12 1 .00 6 19.39 9 4 7 59.56 1.17 7 3.64 10 2 7 88.38 1.06 8 24.92 11 4 8 84.12 1.00 9 9.38 12 2 8 21 .86 .92 10 14.07 13 4 9 3.46 1.00 11 56.91 14 3 9 29.72 1 .00 12 2.45 15 4 10 6.12 1.00 13 30.93 16 2 10 31.08 .96 14 21.76 17 3 10 1.07 1.07 15 16.55 18 5 10 44 .44 1 .00 19 1 11 67.15 .91 20 2 11 59.83 .79 21 3 11 50.46 1.17 22 5 11 71.42 1 .00 23 2 12 8.88 1 . 18 24 1 13 28. 22 .83 25 4 13 77.34 1.00 26 3 13 45.60 1 .00 27 5 13 20.67 . 88 28 4 14 37.76 1.13 29 2 14 18. 16 .98 3U 3 15 67.62 1 .00 Figure 1. A Single Commodity Generalized Transshipment Problem ( GT ) 10 figure 2 shows a basis for the problem introduced in Figure 1. (In this simple case, there is only one component: p = 1.) A unique subgraph partition ot G denoted G corre- sponds to b . Let A = {a-^la- 1 is an arc associated with b J a column ot b} , then G^ = [N,A. J denotes the directed graph £ associated with the basis B . To each submatrix b ot b £ £ £ there corresponds a component ot G denoted by G = [w ,A ] . £ £ £ N is the set ot nodes corresponding to the rows ot b and A £ is the set ot arcs corresponding to the columns ot B . It is £ known (e.g., [13]) that G is either a rooted tree or a one- l tree (a tree with an additional arc forming one cycle). If G £ £ is a rooted tree, then B is upper triangular. If G is a £ is a one-tree, then B is upper triangular except tor one element below the diagonal. Each component can be viewed as having one cycle if we assert that each rooted tree has a self- cycle corresponding to its root node. In our example basis, the subgraph G has only one com - pone nt , (one-tree) shown in Figure 2: G is a one-tree, with nodes 2, 3 and 11 composing its cycle. As in the pure network case, this near tr iangulat ion and associated sub-graph are naturally represented by a pre decessor t_uncticm p( ), and a predecessor graph (which does not preserve the orientation of arcs in the original network). The predeces- sor function can be used iteratively to construct the unique backpath from any node to the root (or cycle); the backpath includes all nodes on the cycle. The immedia te successors ot a 11 Column 2 4 12 29 23 10 20 21 17 26 24 7 30 14 13 Row 2 5 14 12 7 11 3 10 13 1 6 15 9 4 11111 -.74 -.92 -.98 -1.18 ! 1 1 | -1 1-1.06 -.79 -1.17 111 111 -1.07 -1.0 j -.83 i 1 -.91 -1.0 -1.0 -1.0 A Preorder Near-Triangulation (The underlined coefficient is below the diagonal.) Join 00000 0©0©0 t t One-Tree Figure 2. A Generalized Transshipment Basis (For the problem in Figure 1.) 12 node, if any, are the first nodes encountered on all paths ex- cept the backpath to the root, and all the nodes on these paths are called successors . Note that each basis may have many near tr ianyulat ions . However, all such near trianyulat ions yield the same predecessor function and graph (where the right to left ordering of succes- sors of any node is immaterial). Thus, the predecessor graph does not completely represent a near tr iangulat ion without additional information: an ordering of the rows (nodes). For algebraic reasons, we restrict such partial ordering to preorder [6], in which a node i always precedes its successors, if any and in which all its successors, if any, precede any node which does not precede node i . 13 4. IMPLEMENTATION For didactic reasons, we bey in by introducing a complete primal generalized network algorithm using a preorder traversal method. Controversial alternatives are deferred until this paradigm is presented. Hereafter, notation with upper case roman letters followed by parentheses indicates a program data array. For instance, the predecessor array is referred to as P( ). Static arc storage is used for tails T( ), heads H( ), costs C( ), multipliers MUL( ), and capacities CP( ). Contigu- ous storage by tail, or by head node reduces T( ) , or H ( ) to an hierarchical node-length entry point array. (GENNET uses contiguous storage by head node, as does GNET.) Lower bounds on arc flow are translated out prior to solution, with appropri- ate adjustment of the initial right-hand side of (GT), and of CP( ). The sign bit of CP( ) is available to indicate arcs nonbas ic at their upper bounds ( ref lected with flow -CP( ) ). The predecessor function and its array P( ) are defined so that the basic arcs in a cycle are oriented uniformly in a directed cycle. All basic arcs not on a cycle are oriented so that a backpath is created to a cycle. To obtain this orienta- tion the direction of some arcs must be reversed, and the sign bit of the predecessor array is used to indicate: if P(I) < 0, then the orientation of arc (I, -P(I)) is reversea from its original orientation. This is the complement of the discipline used in GNET [6] 14 A depth array D( ) reveals for each node the number of nodes on the backpath before encountering a cycle. Nodes on a cycle have depth zero. Number of successors , or preorder distance are acceptable substitutes tor depth [6], but are not discussed here. A preorder traversal array IT( ) is maintained so that all preorder successors of a cycle node are encountered before another cycle node. It is convenient to make this a circular list for each near triangular basis component by setting IT( ) of the last preorder node in the component equal to the first preorder node in the component. The components of G are not inter-connected, or equiva- i lently, the sub-matrices B in (1) do not have common rows or columns. Consequently, the p components of a basis may be represented in a single set of node-length arrays. The array X( ) contains the values of basic variables, values of dual variables (or simplex multipliers, or node poten- tials) are stored in U( ), and IVAR( ) gives the location of basic variables in the arc arrays. The array FAC ( ) contains the cycle factors , defined later. figure 3 shows these arrays for the basis given in figure 2. Generalized networks do not exhibit totally unimodular bases. Consequently, floating point representation is required for x( ), u( ) and FAC ( ), and is desirable for arc-length arrays C( ) , MUL( ) , and CP( ) . 15 Node Predecessor Depth Traversal basic Variable Dual Cycle Factor P( ) D( ) IT( ) IVAR( ) X( ) U( ) FAC ( ) 1 13 2 6 24 22.860 16.665 * 2 3 5 2 118.169 47.148 -.48101 3 11 1U 21 45.826 31.678 -.48101 4 9 2 2 13 0.0 5.418 * 5 -2 1 8 4 0.0 27.552 * 6 -3 1 15 7 21.308 -3.782 * 7 -2 1 11 10 3.434 -38.898 * 8 -2 1 14 12 27.U87 27.487 * y ^3 1 4 14 9.380 1.958 * 10 -3 1 13 17 13.150 28.606 * n -2 3 20 4.170 -16.053 -.48101 12 -2 1 7 23 2.076 32.431 * 13 -3 1 1 26 11.956 -13.922 * 14 -2 1 12 29 22.204 29.580 * 15 -3 1 y 30 16.550 -35.942 * Figure 3. GENNtT Basis Representation Arrays (for Basis in Figure 2) 16 Step Si, Priceout The reduced cost for nonbasic arc k , oriented from f. k to t, is (given the current dual solution u and column N ): r. = c, - uN k k = c k " u f + ^kV k k = C(k) - U(f ) + MUL(k)*U(t, ) . K K. (If CP(k) < 0, arc k is reflected and the siyn of r, is reversed.) At most, one multiplication, addition and subtrac- tion are required. Note that the multiplication is unneccessary if | m, | = 1; further specialization is possible for sets of priced arcs with common attributes. If arc k is a logical arc (slack, artificial, or surplus variable) then C(k) can be logically generated, rather than explicitly stored, and r k = C(k) ± U(f R depending upon the sign of P(f ) . K i From the example, r 2l =? C(Z1) - U(b) + MUL( zl )*U( li) = 2U.57 - (27.552) + ,8b(-13.922) = -19.133 . 17 This variable will be used as the entering non-basic variable for further illustration. Step S2, Ratio Test For the determination of the arc to leave the basis the system of equations BZ k = w k , k k must be solved for the transformed column Z . ( -N is used if arc k is reflected.) Due to the near triangularity of B , this incoming column transformation can be combined with the ratio test in a single integrated process. Suppose that N has two nonzero coefficients representing an arc oriented from f, to t, , and that the arc is not reflected. An apparent complication arises if f and t, are s t in separate components of B in (1), say B and B , respec- tively. In this case two disjoint subsystems must be solved and the results added to determine the nonzero elements in Z The subsystems are: B (J = e f , and k B fc g k = -m R e , k 18 f k t k (e is the t, th unit vector; (J and Q are disjoint k k components of Z . ) This complication is inconsequential. In order to see this/ consider solving one of the systems: ut„ fc k b g = -m.e k The only nonzero elements of Q will be those that correspond to the nodes in the backpath from node t . As we shall see, this follows from the manner in which the coefficient -m, in — k row t propagates during substitution solution. K. Suppose that G is a rooted tree. The backpath from node t, can be denoted by iteration of the predecessor function: t , K. K p(t, ), p(p(t )),..., root; this sequence is shown below as a, K K d-l,...,0, analogous to the backpath length remaining to the root Root P(p(t k )) l d-2 d-1 P(t k ) l d-l 19 The values of a and b for each basic arc depend upon the original orientation of the arc, given by the sign of P( ). For original orientation, P( ) > implies that a = +1 and b is the multiplier value, P( ) < U reverses these definitions. The triangular system corresponding to this backpath is: a Q b x a d-2 b d-l a d-l b d Q U • u -m, By back substitution, its solution is: m, — k — d g 6 = " b 6+l q 6+l = 3 6 Now suppose that G is a one-tree. Let node t, be on the cycle. The backpath is t , p(t. ) , p(p(t, )), ..., c, with p(c) t and length p ; this sequence is shown below as -1 ,-2 , . . . , -p , analogous to the backpath length beyond the cycle start . 20 The corresponding near triangular system is a b ,, -p -p+1 '-P + 1 a -3 b -2 ~P a -2 b -l l -l -m, By modified back substitution, its solution is _m -k 1 " b 6+l q 6+l for 6 = -2, . . . ,-p , where -p -b. f = l - n 6=-l a 6 21 The cycle factor (or loop factor [13]), f , is common to all nodes in the cycle and can be computed when the cycle is created and stored in FAC ( ) for all nodes on the cycle so that it is immediately available at this step. Suppose that the entering arc is ( f , t. ) with coefficient entries (a,,b ) and that the f, and t backpaths converye. Using nomenclature for the backpath sequences introduced above, the corresponding system is: join a u b i 3 £-l b £ w w d-1 d a d-X b d 22 by back substitution, its solution is: b k q, = — , ( d begins t. backpath sequence), Q a . k — d " b 6+l q 6+1^6+1 for 6 = d, d-1 , . . . , r , l d a , (d bey ins f, backpath sequence) , " b 6+l q 6+ tor 6 = d, d-1 , . . . ,w , <±o- i-i - (b w^w + b rV a £-l b 6+l q 6+i q x = tor 6 = 1-2 , 1- i , . . . , U . Suppose that entering arc (f, , t, ) with coetticients (a, ,b. ) enters and that the backpaths converge on a cycle. The corresponding system is: 23 a -p b -p + l a -p+l a -*-l h -l w ■-* a b -s-1 -s -s a -2 b -l ~P -1 a b ^i w w+1 a w+l a d-l b d a r b r+l a r+l \ a d-l b d 24 By modified back substitution, its soiution is: q. = — — , (d bey ins t, backpath sequence) , l d a " b 6+i q 6+ for 6 = d, d-i , . . . , r, "b q 3-8-1 r^r X a , t -s-1 " b 6 + i y 6 + l for 6 = -s-2 ,...,-£, -£-1 ,... ,-s , q-, = — , (d begins t, backpath sequence) , q a , k d b 6+i q 6+l for 6 = d , . . . ,w , <J»-i * -b q *- a -*-i t D 6+l q &+l . q,. = for 6 = X- £ f % • • f P ( ~I-r*«*fO/«*«/ JC / q. = q„ + q, tor 6 = -l,...,-p . 25 For the cycle in Figure 2, i -(-1) -(1) -(-1.17) f = X J - x ^779 x 1 = -.481U1 . By now it should be apparent that one composite back sub- stitution scheme will suffice for all cases. The cycle factor is applied once when, and if, a cycle is encountered on a backpath . If the backpaths of f, and t. converge, let join be the first node on the t backpath that is also on the f backpath. If the backpaths converge on a cycle and it the leaving arc pre- ceeds the join on the f, backpath, define join as the first cycle node encountered on the f backpath. If the backpaths do not converger join = <f> . Several schemes are available for identifying the join efficiently [6J. The depth (or number of successors, or pre- order distance) of nodes on the backpaths can be used to avoid iterating either backpath past the join. Depth, the number of nodes on the backpath until a cycle is encountered, can be used to indicate which backpath node is deeper and should be iterated. When both backpath nodes have matching depths, zero depths indi- cate that each backpath is on a root cycle. By remembering the first root cycle node on each backpath, further iteration will either reveal the root cycles to be distinct with join = $ , or coincident with join defined as above. When both backpath nodes have matching depths greater than zero, the nodes are compared for 26 equality. A match indicates the join, and a mismatch indicates that both backpaths should be iterated tor another comparison. If a join is encountered, the backpaths have merged and either one can be used to complete the ratio test (GENNET continues the t backpath ) . When a cycle is encountered the q values are computed on the first pass around the cycle. On the second pass around the cycle the q values are computed and the ratio test is completed. As the backpaths are iterated, the column transformation is applied and the resulting terms of transformed column, Z , are used in the rat io test , seeking the minimum ratio: MIN i CP(k) XU) '9. CP( IVAKU ) )-XU ) - z the capacity of the incoming arc, for z > 0, node I, for z < . If a zero ratio is encountered during this process, the ratio test may be preemptively terminated. Step S3, Pivot IF CP(k) is selected as the minimum ratio, then the enter- ing variable remains nonbasic and is reflected to its opposite bound. Only the flows X( ) need be updated. To do this, the backpaths are iterated again ana for each node I encountered, X( ) is reduced by z x CP(k) . 27 If a basis exchange is required, an efficient update of the basis representation must preserve the rooted cycle orien- tation, updating some entries in P( ) and D( ). Also, some flows X( ), and some dual variables U( ), must be changed. FAC ( ) must be established for nodes on newly created cycles. Some bookkeeping in IVAR( ) and CP( ) may also be required. The apparent intricacy of our task is deceiving. Careful analysis yields an elegant solution. However, the supporting arguments require close attention. To simplify the explanation, reorient the incoming arc (f, ,t.) to (i,j) or (j,i) (if necessary) so that the minimum ratio is on the j backpath . Let the entering arc (i,j) have the outgoing arc (c,d) on its j-backpath. Also, reorient the outgoing arc (if necessary) so that the first node encountered on the j backpath is c . Figure 2 shows a case for which both reorientations are necessary. In our example, arc 20, oriented from node 2 to node 11, leaves the basis and arc 27, oriented from node 5 to node 13 enters. We call the backpath segment from j to arc (c,d) the j - stem . In Figure 2, i, j, c, and d are shown. The j-stem is composed of nodes 5, 2, 3, and 11. The skeletal update: a. Reverses the orientation of arcs within the j-stem; and b. Urients the entering arc so that it precedes this redirected path with the same orientation. If the i and j backpaths merge, the node where they merge is called the join. The join is node 3 in Figure 2. 28 If the leaving arc lies beyond the join on the backpaths a new cycle is created in the basis exchange. In this case the portion of the i backpath from node i to the node with the join as its predecessor is called the i-stem. When node i is on the j backpath the i-stem is null. In Figure 2, the i-stem is node 13. Figure 4 displays the pivot logic to be applied. The algorithm visits each node affected by the basis update exactly once . It proceeds up the j-stem one node at a time visiting the preorder successors of each stem node via IT( ). If a join is encountered it switches to proceed up the i-stem one node at a time and then returns to the j-stem. At each j-stem node, the successors of the next lower stem node have already been visited. The unvisited successors of the current stem node can be divided into two groups: the left successors are the nodes visited in preorder by iterating IT( ) from the current stem node until the next lower stem node is encountered, and the right successors of the stem node are the remaining unvisited nodes reached by further iteration of IT( ). In Figure 2, nodes 8, 14, 12, and 7 are right successors of node 2 and nodes 10, 13, 1, 6, 15, 9, and 4 are left successors of node 3. As we climb the i-stem the traversal IT( ) is moditied so that the last of the left successors (if any) points to the first of the right successors and the last of the right successors (it any) points to the previous stem node. As we climb the j-stem, the traversal is modified so that the last of the left succes- sors (if any) points to the first of the right successors and 29 f BITEH J Set first J-stem node to visit. Visit left successor Figure 4. Pivot Traversal Scheme 30 the last of the right successors (if any) points to the next node up the stem (because the update reverses predecessor orientation for the j-stem). The immediate switch to the i-stem upon encountering the join during the j-stem iteration is motivated by a subtle complication: while visiting the left and right successors of the j-stem nodes, the nodes on the i-stem and their succes- ors must be skipped if encountered. Because the i-stem nodes are successors of the join, visiting the i-stem as soon as the join is encountered (if one exists) on the j-stem leaves us with the preorder successor of the last i-stem node visited. This valuable artifact enables the subsequent j-stem iteration to immediately skip all i-stem nodes and their successors should they be encountered. This is the key step preserving an effi- cient one-pass basis update. In Figure 2, i-stem node 13 and its successor node 1 are successors of the join, node 3. The preorder successor of the last i-stem node (called the preorder link in Figure 4) is node b. A stem node may have a right successor which is on the root cycle (with depth zero). The preorder traversal array is organized so that all successors of a cycle node are encountered before the next cycle node. This implies that a cycle being broken by a leaving arc will always be encountered as a right successor of a stem node. In Figure 2, the broken cycle is encountered as a right successor of node 2; if arc (2,11) were not the leaving arc, then node 11 would be a preorder successor of node 2. 31 The basic arc flows, X( ), are changed as each arc is visited on the j-stem, and (if a new cycle is created) on the i-stem. If no cycle is created, X( ) is changed only it the minimum ratio is nonzero, and then only on the arcs visited on the backpaths. During the update, the dual variables must be changed so that c k = u - m R u , k k for every basic arc k oriented from node t to node f , K K With incoming arc (i,j) (reoriented as in Figure 2 so that the outgoing arc is on the j-stem), this relationship is retained for all nodes except for those which the update changes to be successors of i: (i.e. the nodes on the j-stem and their successors ) . If a cycle is not created, this update proceeds for each j-stem node and its successors as these nodes are visited in preorder. for node s , and associated basic arc £ oriented from s to p (s ) , U(s) = CU) + MULU) * U(P(s) ) , while the reverse orientation -P(s) to s yields U(s) = (CU) - U(-P(s)))/(-MULU)) . 32 As in pure networks, the preorder traversal assures that a value of the dual variable of a predecessor node is always determined prior to its use by any immediate successor nodes. When a cycle is created, the dual variable must be deter- mined for one of the cycle nodes. Then, the dual variables of the remaining cycle nodes and their descendents can be found one at a time in preorder traversal. This key dual variable is computed immediately after the ratio test predicts creation of a new cycle (the leaving arc lies beyond the join on the ratio backpaths). Consider the current context for a new cycle: rr / * P(i) P(j) from which the new cycle will be formed (with modified pre- decessor function p ( ). ) : /a p<p(j>> rv rx p(j) ^ — ^d , ., a „ v b ^ a ,, ^ b , -p+1 d -3 b -2 a -2 - P -1 33 The corresponding near-trianyular system is: (u -p' u - P+ i / u_ 1 ) a b , -p -p + 1 -p + 1 -p a_3 b_ 2 a -2 b -l = (c ,c -p - p + i . 1 r.../C_ i ) where the indexing of c is understood to yield basis arc costs (accomplished by indexing with IVAR( ) ). The determination of one term (say, u_, ) of the solution of this system is induced by modified forward substitution to be (e .g . , [ 11] ) : - U' , X — , -1 ""I where (the cycle factor) f is defined (again) as -b f = 1 6=-l and u' is the corresponding term of the solution of the strictly upper triangular system (omitting the cyclic coefficient b , -p and using forward substitution): 34 c - p a -p i & 6 6-1 r , , , u; = for 6 = -p+l,...,l . o a „ Note that computation of u' requires that we traverse the new cycle in a direction opposite to its predecessor orientation. However, before the update creating the new cycle, the j-stem exhibits proper orientation for at least part of our work. Thus, we can complete the first portion of the forward substitution for u\ and accumulate the associated partial product component of f while iterating the j-stem before the update. The remainder of the new cycle is accessed by iterating the i-stem. As we proceed up the i-stem the remaining product terms of f are accumulated, and the reverse i-^tem path is stored (e.g., using U( ) locations, which contain obsolete dual values to be replaced during the imminent update). Reaching the join, this stored reverse i-stem path is then accessed to complete computation of u' . The newly created cycle in Figure 5 is composed of j-stem nodes 5, 2, and 3, and i-stem node 13. The dual system becomes ~ u 13 + U 3 = 4b,6(J u + u = lb .47 u - -74u 5 = 26.76 88u + u = 20.67 , 13 5 35 which yields uj 3 = -19.0142, f = 0,3488 (the new cycle factor), and u,~ - -54.5132 (the new dual for node 13) Thus, when a new cycle is to be formed, the new cycle nodes must be visited once (after the ratio test and prior to the pivotal update). The one-pass preorder traversal update can now proceed as presented in Figure 4. The basis representation arrays are all modified on-the-fly during this traversal. The update of nodes on a newly created cycle (if any) includes establishing the new cycle factor PAC ( ). Changes for P( ), D( ), IT( ), IVAR( ), X( ), and U ( ) proceed analogously to the pure network case (e.g., GNET [6]) with simple modification of generators for X( ) and U( ) to accommodate generalized network coefficients for basic arcs. Figure 5 shows the new basis (derived from the example in Figure 2) before restoring near-triangulat ion with the update. A new cycle is to be created and the new cycle factor and a dual solution (for node 13 on the i-stem) are found at this stage. Figure 6 displays the new basis restored to near-triangular form. At this point, all basis representation arrays are updated with the values shown in Figure 7 (data in italics has been changed by the update ) . 36 Column 2 4 12 29 23 10 21 21 17 26 24 7 30 14 13 Row 2 5 8 14 12 7 11 3 10 13 1 6 15 11111 -.74 -.92 -.98 -1.18 i 1 ! i 1 ! ! -1 -1.06 -1.17 111 111 -1.07 -:ti -1.0 J -.83 1 -.91 -1.0 -1.0 -1.0 1 (The entering column is italicized. ) 11 •< *//< •- /// One-Tree Figure 5. New Basis Before Restoring Near-Triangulation (entering column 27, arc (5,13)) 37 Column 4 12 29 23 10 2 17 7 30 14 13 21 26 24 27 Bow 2 8 14 12 7 3 10 6 IS 9 4 11 13 1 5 11111 -.92 -.98 -1.18 -1.06 1 1 -1.0 1 1 1 1 -1.07 -.91 -1.0 -1.0 j 1 1 -1.0 V 1 „ -1.17 -1.0 -.83 -.88 1 1 A Preorder Near-Triangulation Figure 6. New Basis Near-Triangulated 38 GENNET is designed to exploit intrinsic pure network struc- ture commonly found embedded in generalized network problems. Note that when this algorithm is applied to pure network prob- lems, _i_t automatically adapts to a minor variant ot GNET. of course/ GENNET uses floating point arithmetic operations which are intrinsically slower on most computers than the pure additive integer arithmetic of GNET (also, floating point arithmetic requires some extra editing for mantissa truncation errors). To mitigate this disadvantage, GENNET can test logically for pure network arcs (with unit multipliers) and avoid unnecessary floating point multiplication and division operations. GENNET also employs an automatic dual basis aggregation refinement ([6], p. 26 f f ) . Explicit values for D( ), U( ) and IT( ) are maintained only for nodes with successors. An array, A( ), records for each node the current number of its aggregated successors . When an aggregated node is encountered in the priceout, its dual is generated from that of the immediate predecessor of the node. When a backpath ot an entering arc begins with an aggregated node, it is disaggregated, and when the leaving arc isolates nodes with no successors, they are aggregated . Figure 8 shows the arrays affected by the dual aggregation scheme for the basis in Figures 2 and 3. IT( ) indicates aggre- gated nodes with (J entries, and these nodes have broken outlines in the one-tree depiction. The priceout of arc (5,1J) requires that U(5) be generated using the predecessor dual U(2). Node b subsequently starts the j-backpath and is disaggregated. The update leaves node 11 with no successors, and thus aggregated. 39 Node Predecessor Depth Traversal Basic Variable Dual Cyc l e Fa ctor P( ) D( ) IT( ) IVAR( ) X( ) U( ) FAC ( ) 1 13 1 5 24 22.86U -17.026 * 2 5 8 4 3.883 6.557 0.3488 3 -2 10 2 118.456 -8.913 0.3488 4 9 2 lj, 13 O.U -35.173 * 5 U 2 .27 2.872 -27.302 0.3488 6 -3 1 11 7 21.308 -48.388 * 7 -2 1 3 10 3.434 -77.192 * 8 -2 1 14 12 27.087 -16.634 * 9 -3 1 4 14 9.380 -38.633 * 10 -3 1 6^ 17 13.150 - 9.330 * 11 -3 2, 13 21 48.641 -50.746 * 12 -2 i 7 23 2.076 -1.969 * 13 -3 1 26 9.428 -54.513 0.3488 14 -2 i 12 29 22.204 -11.834 * 15 -3 i 9 20 16.550 -76.533 * Figure 7. GENNET Basis Representation Arrays (for Basis in Figure 5) 40 Aggregated Node Depth Traversal Dual Successors D( ) IT( ) U( ) A( ) 1 2 3 47.148 5 3 U 13 31 .678 3 4 5 6 u U 7 8 9 1 11 9.380 1 10 11 2 4.170 12 13 1 9 11.956 1 14 15 Aggregated One-Tree \ ^^^ Figure 8. Aggregated Basis Representation (for Figures 2 and 3) 41 b. COMPUTATIONAL EXPERIENCE Significant design alternatives tor GENNET have each been evaluated by extensive experimentation at large scale. Illus- trative computational experience is abstracted in this section for some of the prototype systems tested. Departing somewhat from the style of the paper documenting such work for GNET lb], relative performance is reported even lor some competitive design features subsequently rejected tor adoption (some readers of [bj have concluded, quite incorrectly, that only those teatures reported tor GNET were tested). This should help other re- searchers avoid our mistakes, and may even change some widely held misconceptions and correct a few translation errors in textbooks . Among the key issues to be resolved are: a) Static Storage of Arcs . The arc lists can be stored in arbitrary sequence, or, to save space, arcs can be stored contiguously by tail node, or by head node, thereby replacing an arc-length index array by a (presumably much shorter) hierarchical node-length entry point array. b) Preorder Manipulation of Basis . The triple label (aug- mented predecessor index) method [zu,22\ will be presented and compared with the preorder traversal method . c) basis Aggregation . An aggregated basis representation will compete with an explicit representation. 42 d ) i^lci n y Sch emes . Candidate list schemes and explicit arc pricing mechanisms widely used in general linear programming systems will vie with dynamic candidate queue disciplines. e) Pure Network Speci alization . Generalized network algorithms would ideally adapt to pure networks with efficiency comparable to pure network codes. t ■) Star tin g , Tunin g and Tai lor ing . Which algorithm parameters and settings lead to high efficiency tor interesting classes of problems? Are heuristics tor advanced starting solutions worthwhile? ^ ) Generali za tions. Advanced teatures and generalizations will be suggested. Computational tests have been made with many problems, including a benchmark suite ot pure network problems generated with NETGEN [28J, and generalized network problems generated by NETGENG [17,1»J. figure y gives some problem characteristics Static arc storage has been implemented in three ways: Contiguous by head node with hierarchical node- length entry point array, Contiguous by tail node, with hierarchical node- length entry point array, and Explicit arcs in arbitrary order. In competition with our basis manipulation using preorder traversal, a triple label representation originally suggested by h'. Johnson [2bJ has been implemented tor pure networks and called the augmented predecessor indexing method by Glover, Karney, and Klingman [2UJ. Glover, Klingman ana stutz 122J 43 Percent Problem Nodes Sources S inks Arcs Capacitated NETGEN (Pure ) NG15 4UU 200 200 4,500 NG18 400 8 60 1,306 20 NGly 400 8 bO 2,443 20 NG2 2 400 8 60 1,416 40 NG2 3 400 8 60 2,836 40 NG26 400 4 12 1,382 80 NG27* 400 4 12 2,676 80 NG28 1,000 50 50 2,yoo NG29 1,000 bO 50 3,400 NG30 1,000 50 50 4,400 NG31 1,000 50 bO 4,800 NG3 2 1,500 75 75 4,342 NG33 l,b00 75 75 4,385 NG34 1,500 75 75 5,107 NG3b 1,500 75 75 5,730 NETGENG (Generalized ) GTU1 200 100 100 1,500 GTU2 200 100 100 2,000 100 GTU7 300 135 115 4,000 GT12 400 20 100 5,000 GT 1 5 1,000 5 yy5 4,000 100 GTlb 1 ,000 20 100 6,000 100 GT 1 8 1,000 30 400 7,000 Figure y. bpme Benchmark Problems (Problem NG27 is extensively studied in [6J 44 report that the method has been extended to generalized networks, but reveal no details. We have implemented our own elticient version of the triple label scheme. The triple label representation uses predecessor, successor, and brother pointers lor each node. Figure 1U shows these arrays lor the basis in Figure 2. To briefly illustrate the triple label scheme tor gener- alized networks, let our situation be the same as tor the pre- order example (the j-backpath of the incoming arc is arranged to include the outgoing arc and to encounter node c first on that backpath). Let V d+1 = X ' and the backpath from j to c v d' V d-1 C The skeletal update scheme is 45 Node Predecessor successor brother P( ) S( ) b( ) 1 13 u 2 3 11 10 3 11 2 4 9 U 5 -2 8 6 -3 o lb 7 -2 u 8 -2 u 14 9 -3 4 10 -3 U 13 11 -2 3 b 12 -2 U 7 13 -3 1 b 14 -2 u 12 15 -3 u 9 Figure lu. Triple Label Basis Representation ( for t' igure 2 ) 46 1. Set 6 = d 2. It s ( v 6 _^) = v 6 ' ^oto Step J Otherwise, it possible, tind a node v* such that P(v*) = v 6-1 and B( v* ) = v r , — — > — o and set B(v*) «■ b(v ), Go to Step 4 3. Set b (v (5 _ 1 ) + B(v 6 ) , 4. Set P(v 6 ) ♦ v 6+1 , B <V * S(v 6 + 1 } ' b. Set 6 <■ 6-1, Then, it v x + i * c, Go to ^tep_ 2, Otherwise, Stop . (Step 2 exhibits the key extension tor generalized networks of the pure network triple label scheme.) These tive steps reverse the orientation ot all arcs on the j-stem and orient ( v ^ + _i' v h) so that it beyins this redirected path. All other triple label operations are obvious alterations of the preorder traversal procedure. A static cand idate list pricing strateyy (e.g., [I7,la,jl])) an explic it arc pricing method reminiscent ot yeneral linear proyraioming systems, and a dynamic candi date queue lb] have been tested. The (L,,L.J candidate list procedure is a simple strategy. Nodes with leaving arcs (or entering arcs for contiguous head node arc lists) are sequentially priced, placing the most nega- tive candidate (it any) trom each node on the candidate list 47 until L,. candidates have been located (arbitrarily organized arc lists require explicit arc pricing). Entering arcs are chosen trom the candidate list by most negative reduced cost until L, iterations have been carried out/ or until all list entries have non-negative reduced costs. Arcs with non-negative reduced costs are dropped from the list and replaced by continuing to scan the nodes until L.. candidates have been tound. It an exhaustive pass through the nodes results in less than L candidates, then an optional closing gambit sets L equal to the actual candidates found/ and reduces L, by halt unless L, equals one. Glover/ Hultz, Klingman and Stutz [ 1 7 , 1 8 J report (L,,L y ) ot (b,lU) to be best in their work. The candidate queue is a dynamic list ot interesting arcs and nodes, scanned in a cyclic manner. The entering arc is selected from the queue by pricing nne entries; if an interesting node is encountered it is replaced by its best-priced entering arc (or leaving arc tor contiguous tail node arc lists). Arcs pricing favorably are retained in the queue. When the end of queue is encountered/ the queue is refreshed by pricing 1PG nodes in a cyclic general arc scan. During an opening gambit ot nns pivots, the nodes incident to the entering basic arc are added to the queue. There is no closing gambit, since the queue automa- tically shrinks and finally collapses at optimality. Bradley, brown, and Graves 16J suggest WNhi = 62, NNh = 3m/4 and 1PG = m/lU + 1 for pure network problems with m nodes. A rule to break ties in the ratio test which guarantees finite convergence for pure transshipment problems has been 48 developed by Cunningham [12]. Bradley, Brown and Graves 16] show that the conditions necessary for finite convergence are naturally satisfied by GNET on over yu% of its degenerate pivots. Elam, Glover and Klingman lib] have observed that the results of Cunningham can only be extended to the generalized network case when the multipliers are positive. We have not used Cunningham's modif ication . Bland [4] presents a class ot restrictions of pricing and ratio tests for general linear programs which relies exclusively on primal simplex representation and guarantees finite conver- gence. These rules are easily modified to produce an efticient finite simplex algorithm. The modification interferes with ef- fective pricing strategies only during degenerate pivot sequences, and the restrictions increase in severity only with the number of pivots in that sequence. However, during a degenerate pivot se- quence restriction records must be accumulated (e.g., a list with each incoming variable in one ot our schemes ) . This record is naturally accommodated by the dynamic candidate queue, but not by a static candidate list or explicit pricing. No purpose is served by reporting such unbalanced competition. A starting strategy has also been tested in conjunction with pricing alternatives. A straightforward starting method ex- amines each node with supply (or demand tor contiguous arc stor- age by head node) and assigns as much flow as possible to its least cost leaving (entering) arc. The procedure stops when an exhaustive pass of the nodes makes no additional flow assignments. 49 The starting solution achieved is not necessarily leasible. (F.g., "exhaustive pass sequential source minimum start" [Its].) Artificial arcs are driven from the basis in all experiments using a big-M method (e.g., [b]). This choice is principally motivated by the comparability of competitive tests between pure and generalized network codes on pure and generalized network problems. (A two-phase method is employed in production use.) Choosing the best big-M value is a bit tricky. The smallest big-M value which yields a feasible optimal solution (if one exists) is best in our experience. Small big-M values may tail to produce feasibility, and large values inflict numerical difficulties. In practice, a default value is used and an automatic restart recovery is applied it an mteasible solution persists. If a restart with a higher big-M value tails to reduce the total inf eas ibility , a terminally inteasible solution is declared. Figures 11 and 12 indicate the multiple ot maximum absolute arc cost used for big-M in each problem. Computational tests have been performed on various computer systems. The times reported here are accurate to the precision displayed tor IbM 37u/lbb-3 using the FORTRAN H compiler with OPT(z). Solution times exclude input/output overhead. Ch'NNtT uses double precision (IbM RhJAL*b) arithmetic for floating point operations and storage. solution times are given in Figure 11 for the pure network test problems. Performance is given tor three pure network codes (two versions ot GNfcT and SUPFkk) as well as tor several repre- sentative generalized network prototype systems. we can thus 50 compare the best generalized network scheme with the best pure network code exhibiting equivalent features (GNfcT, or its aggregated successor variant [6J). The times for SUPEKK also provide the only available objec- tive means lor comparison of our implementations to that of Glover, Hultz, Klingman and btutz [14]; they report that their generalized network code, N£TG, is about as fast as stiPfcKK (a fast out-ot-kilter code for pure networks lij) when both are used to solve pure networks. Using their version of SUPtKK on our computer we have shown that our implementation of NETG (called TLA in our nomenclature) is at least as efficient as their claims for NfcTG. 51 ^ un 3 X> f ^H 00 r- oo CN 00 r» CN n ^r CN r- * —i oo oo JO .— 1 OO rH in ^r JO ^ 00 oo 00 r> 00 a a< ."N r-H (N r-i Ox) .—1 ^r ^r in -n r» X> X oo ^r ^ £> co a CD ■o O :; 00 r- oo vO 30 oo o P~ xO 3 jo oo oo OO ■=T 3 oo H II >< 1— 1 oo JO 00 ■^ 00 "* m <=r *=}< JO cn JO xj> 3 r- -* aJ a a u Z z .-> —1 i— 1 r-H — i f— 1 cn N CN oo 00 o 'J z x: rH <D Z £> <-< r- om r^ 00 00 VO 00 r-H ^H r- oo X> —1 =* r- r- CD H II oo ^ KD <tf «tf 00 ^ >r JO ■JS sD •JO JO 00 r- oo u a a a, 3 2 z m r-H r— t f— 1 — i H cn •N r*i N -H X 3 z 33 CN ■M Li fO to ■30 II a z z CN OO Ox) 00 II a z z .-M 00 II a z z g I Oil X - < in a g <L> z r-H a .Q :o H Li a Oi z in jo o o r-\ CN OM OM mouOLDO-n-non H CN H M ,M H H M M UO JO r- ^ 3 x> ^H JO 00 r^ OO o OO CN 00 oo o — 1 o •<* r^ t r» ■xf JO xO 00 oo JO i— I =tf JO 00 00 N r-H -H rx) Ox] CO oO OO oo r- CN o 00 00 ^ OO CN r-H 00 v£) r~- r- OS) JO r- 00 X m v£> r^ <£> xO m r- cH 00 OM OJ ■n <—\ ■* JO J-l CN OJ rH Osl .N 00 ~n 00 00 oo OJ On) o 00 o r~ Nf ^r 00 r-\ OJ o 3 og x> 00 r~ r- rH 00 oo m 00 JO r- r- f—l o ^r ^r tt 00 v£> 00 00 OM 00 oo 00 ^r ^r ^r n oo .^0 r- O O rH UO oo ^r 3 OnI 30 CO OnI 00 ^ ,00 r-^ rH VO r— 1 sO r- ^r r~ -NI ,00 JO JO 00 r-{ On) 30 r- .N rH OM Osl og OJ 00 ^r ^r .00 i-H .00 00 oo oo f—t r- o oo rH CN xO N 'X> 00 oo CN oo OO 1— 1 oo r~ oo vD 3 rH >^ 00 in ^3" r~ 00 CN r-\ 30 JO Ov! r-< r-t rH rH "* JO ^r JO 00 oo 00 3 00 • • to c JO 00 oo CN 00 liO r- OO oo 3 rH og OO ^ JO c CO r-t r-t rH Osl OS] CN CN CN CN .00 OO OO OO .00 00 3h ■3 '■J 'J) 'J ^J a '■J ■J) 'J 'J :d ■~J :d :•) CJ5 o Z z z Z z z z z z z z z, z z, Z to u <D > (0 g u U) ■iJ CO c •H g c u CD •H fO o rH 4J x: -o CO A fO o Li l N <D u •r-l e CD CO % c r^ CO <T3 ai 3-1 0) 4-> 31 c N. o CO Vj •f-i r-t o <D u Q) 3 H •H A CO U tl (0 ^: iJ a r-{ TD Li fD CD CD CD 4J 5 rH 3 rH <T3 jJ •H CD 04 Ol CD <n 3. •H CD Z Hi o< u U \ \ +J Ji CD ■a JJ rji Li rO ID M fC 3 fl) •H rV :u XI J| <C| XI • • • "O —1 c Q) CD 3 Li CD 3 J 01 s I en ■h r> in in j~) in H x> ■^ <T\ x> r- *& ro rM lO a2 rH X H <T\ lD r~ r-H 30 [» rH Z II Ol 2 u o cn ro rsi ro n CTi 30 u z X cm o z CN o r- CN CN -Y r 00 ro 3 i— i 00 3 <J\ r- sf r-H II a< a 3 CM «* ro r in 3 ^r 2 X i— 1 ro CN ro II W z m 33 4J CN Ll CO CN fO II < •_n ■U a 3 w z X r-H CN 00 II < M O Z H Z VD ro 3 x> . ■0 7i .•-o r-> Li": =* Xl ro =* LD r~ rH P~ CT> -r «* 00 3 VO CN o ro 50 Ti CO r-H rH ro o IT) CN r-H J-'. -? - tf r "=1" i-H < x> .—1 CTi in H 00 rH ^r ». ij • LO EH rH -f N ^r t N rH CT\ rH <3< ro X -r ro rH 23 r-H x> 00 t, kO r- •^r r~ 3 =* 3 00 CO CN x> --;■ -O N • rH 3 •~n X> 3 o =-r 37> r o ^t* 3 00 3 ro .--. -J* X ro CN ro ro CO CN r rH CN X> N ro j X) c -H c ro ■H ■n •U U ■" OJ N > •H . c e u re CO 4-J rj) • H Li c Li ro CD £ a Cfl 'J U Li Li Li 03 O CO E 03 w Li C/3 +J cu| 03 •H c L) " ■ H — U •H ; 3 H •H jQ W a L : X a r- 1 a Ml '. \ 03 03 jj -J. 3 rH rO ■H 03 3, •j] ro 3 Hi Ol Li S. X \ -I-) L^ -o iJ rji • M ro m T xl <l XI • e 73 . ' rji 03 J 1) H ~ L a< .0 03 H '3 '-' C V c Z N 03 : - 3 4J Li ro 4J W3 uO 00 lO o rN rH ro rH X> ^r m ro .>..' ^r J> r- <* rH CO N o ro Ti N rH r ^ lO c- 'J 03 Z r 1 M .a 'J) H M M Q, z CM r r-H -N f) X 3 3 rH r-i rH H rH H r-H H H ^1 H Eh C3 CJ ; ^ r J J 'jo .--• 'TO c c ro r: o 03 j" 53 However, in these tests TLA is substantially outperformed by the alternate systems. For the pure network problems, best performance is achieved by (Figure 11): Contiguous arcs by head node, Candidate queue pricing, Preorder traversal, and Aggregated successors. TLA does not incorporate any of these features, and is generally less than half as fast as competitors. Although GENNET (HgPX) should in theory rival GNET (HyPX) with pure network problems, the overhead of testing in GENNET for more general basis structure and the additional computa- tional burden of floating point arithmetic exact a performance penalty of about 15 percent. Figure 12 shows solution times for the generalized trans- shipment network problems. Note that the starting strategy helps candidate list performance and hinders the candidate queue. Arranging arcs contiguously by head node dominates both tail node and explicit arc list designs. The candidate queue pro- vides good performance _if_ accompanied with contiguous arcs by head node. Preorder traversal continues to provide better performance than triple label representation in all design contexts. Aggregated successors offer a pronounced advantage. GENNET (HyPX) provides best overall performance. It offers a decided advantage on problems with many more sinks than sources, a situation common in real life. 54 Figure 13 displays performance of GENNET(GN, HPQX) applied to a set of (GN) problems extracted from a collection of real- world LP/MIP models [1UJ. Despite the slight additional floating point arithmetic required to solve (GN), GENNET solves these problems much faster than would be predicted by experience with randomly generated GT problems. Tuning of the pricing mechanism greatly enhances this difference. Problem Node-; Arcs Seconds Pivots AIRLP 17U 3,040 2.62 420 CUAL 170 3,923 1.80 471 STEEL 422 1,279 .39 499 FOAM 951 4,953 3.74 1,258 ODSAS(GN) 1,431 4,615 3.22 1,427 ALUMINUM(GN) 2,178 7,216 3.57 2,794 REFINE (GN) 3,110 6,617 4.72 3,322 FOOD(GN) 3,716 13,907 12.11 7,004 (Big-M = 10 x largest cost coefficient) Figure 13. (GN) Test Problems (GN rows extracted from real-world LP/MIP models [10 with null columns deleted and slack arcs added.) 55 Close scrutiny of solution trajectories lends some insight into GENNET's good performance. GENNET has a one-pass itera- tion unless a cycle is formed; a cycle is formed on only b-to-24 percent of all iterations for these problem sets — 5-to-lU percent for most problems. Also, the explicit (non-aggregated) subset of the nodes is remarkably small, seldom numbering much more than the number of source nodes. Finally, the length of backpaths is quite short, averaging about the number of echelons (path length from sources to sinks) in the model, or just more than 2 in these problems . 56 CONCLUSION The generalized network, system GENNh'T is small, fast and easy to modify. Adaptations have already included using GENNfcT in a system to solve generalized networks with complicating side constraints an/or complicating variables (Mcbride [3UJ). GENNET has also been incorporated in a powerful microcomputer-based network optimization system by Brown, Duff and t'inley [7,lbJ using an APPLE-II host and PASCAL implementation language. Modifications for mixed integer generalized networks have also been tested (though not with care sufficient to warrant publication at this time). GENNET has proven to be a worthy successor of GNET [6]. Preorder traversal is appealing for its mathematical and implementation elegance, and has proven to be efficient and flexible for generalized networks (as it was tor pure networks). (Adolphson and Heum [1] have also suspected this and have in- dependently pursued this avenue.) Experience shows that the GENNET design pertorms much more efficiently on real models than on randomly generated test problems of nominally equivalent size; this design is also tech- nically and philosophically compatible with the various systems we have devised tor solviny other more general classes ot optimization models. b7 The EQRTRAN programs GENNET — (GT) and (GN) versions — (Copyright 1984) are licensed to researchers lor a nominal charge on an exclusive use basis. For further information write the authors via P.O. Box 1832, Alexandria, Virginia, 22313, USA. 5b 7. ACKNOWLEDGMENTS Gordon Bradley and Glenn Graves have always olfered us complete support in our work: we gratefully acknowledge their help and encouragement. Rick Rosenthal contributed many valu- able suggestions to clarify the presentation. John Tomlin has discovered some subtle imprecisions in our description. Kevin Wood has taught from early manuscripts and suffered their shortcomings nobly. by REFERENCES 1. Adolphson, D. and Heum, L., "Computational Experiments on a Threaded Index Generalized Network Code," presented at the Houston URSA/TIMS Meeting, October 14, 1981. 2. Balachandran , V., "An Integer Generalized Transportation Model for Optimal Job Assignment in Computer Networks," Operations Research , Vol, 24, No. 4 (1976), pp. 742-759. 3. Barr, R., Glover, F., and Klingman, D., "An Improved Version of the Out-of-Kilter Method and a Comparative Study of Computer Codes," Mathematical Programming , Vol. 7, No. 1 (1974), pp. 60-87. 4. Bland, R., "New Finite Pivoting Rules for the Simplex Method," Mathematics of Operations Research , Vol. 2, No. 2 (1977), pp. 1U3-1U7. 5. Bradley, G., "Survey of Deterministic Networks," AIIE Transactions , Vol. 7, No. 3 (1975), pp. 222-234. 6. , Brown, G., and Graves, G., "Design and Implemen- tation of Large Scale Primal Transshipment Algorithms," Management Science , Vol. 24 (1977), No. 1, pp. 1-34. 7. Brown, G., Duff, R., and Finley, M., "Design and Demonstra- tion of a Microcomputer-Based Network Optimization System," presented and demonstrated at the Detroit ORSA/TIMS Meeting, April 19, 1982. 8. , and Graves, G., "Design and Implementation of a Large Scale (Mixed Integer, Nonlinear) Optimization System, " paper presented at the Las Vegas ORSA/TIMS Meeting, November 1975. 9. , and McBride, R., "Exploiting Large-Scale Networks with Gains," paper presented at the Houston ORSA/TIMS Meeting, October 14, 1981. 10. , , and Wood, K., "Extracting Embedded Generalized Networks from Linear Programming Problems," Technical Report, University of Southern California, August 1983. 11. , and Wright, W. , "Automatic Identification of Embedded Network Rows in Large-scale Optimization Models," Mathematical Programming ( to appear ) . 12. Cunningham, W., "A Network simplex Method," Mathematical Programming , Vol. II, No. 2 (1976), pp. 1U5-116. 13. Dantzig, G., Linear Programming and Extensions , Princeton University Press, Princeton, New Jersey, 1963. 60 14. Eisemann, D., "The Generalized Stepping Stone Method tor the Machine Loading Model," Man agement S c ienc e, Vol. 11, No. 1 (1964), pp. 154-177. 15. Elam, J., Glover, F., and Klingman, D., "A Strongly Con- vergent Primal Simplex Algorithm for Generalized Networks,' Mat h ematic s of Operations Resear ch, Vol. 4, No. 1 (1979), pp." 39-59. 16. b'inley, M., "An Extended Microcomputer-based Network Optimization Package," MS Thesis, Naval Postgraduate School, Monterey, September 1982. 17. Glover, F., Hultz, J., Klingman, D., and Stutz, J., "A New Computer-based Planning Tool," Research Report CCS 289, Center for Cybernetic Studies, University of Texas at Austin, 1977. 18. , , , and , "Generalized Networks: A Fundamental Computer-based Planning Tool," M anag e ment Science , Vol. 24, No. 12 (1978), pp. 12U9-122U. 19. ^, Klingman, D. and Stutz, J., "Augmented Threaded Index Method tor Network optimization," IN FO R, Vol. 12, No. 3 (1974), p. 293-298. 20. , Karney, D., and Klingman, U., "The Augmented Predecessor Index Method for Locating Stepping Stone Paths and Assigning Dual Prices in Distribution Problems," Transportation Science , Vol. 6, No. 2 (1972), pp. 171-179. 21. , and Klingman, D., "A Note on Computational Simplifications in Solving Generalized Transportation Problems," Tran sport ation Science , Vol. ~! , No. 4 (1973), pp. 351-361. 22. , , and Stutz, J., "Extensions of the Aug- mented Predecessor Index Method to Generalized Network Problems," Transportation Science , Vol. 7, No. 4 (19 7 3), pp. 377-384. 23. Jensen, P. and barnes, J., Network Flo w Pr ogramming , John Wiley and Sons, New York, 198U. 24. , and bhaumik, G., "A Flow Augmentation Approach to the Network With Gains Minimum Cost Flow Problem," Man a gement Science , Vol. 23, No. 6 (1977), pp. 631-64J. 25. Jewell, w., "uptimal Flow Through Networks with Gains," Operation s Research , Vol. 1U (1962), pp. 4' 7 6-499. 26. Johnson, E., "Networks and Basic Solutions," Operations Research, Vol. 14, No. 4 (1966), pp. 619-623. 61 27. Kenninyton, J., and Helgason, R., Algorithms t or Ne twor k Progr amming , John Wiley & sons, New York, 1980. 28. Klingman, D., Napier, A., and Stutz, J., "Nh'TGFN — A Program for Generating Large Scale (Un) Capacitated Assignment, Transportation, and Minimum Cost Flow Network Problems," Manag e me nt Science , Vol. 20, No. 5 (1974), pp. 814-821. 29. Maurras, J., "Optimization of the Flow Through Networks with Gains," Mathematical Programmin g, Vol. 3 (1972), pp. 13S-144. 30. McBride, R. , "Solving Embedded Generalized Network Problems," Working Paper, School of Business Administration, University of Southern California, April 1981. 31. Mulvey, J., "Pivot Strategies for Primal-Simplex Network Codes," Journal of the- Association t or C o mp ut ing Machinery, Vol. 25 (1978), p. 266-270. 62 DISTRIBUTION LIST NO. OF COPIES Library, Code 0142 2 Naval Postgraduate School Monterey, CA 93940 Dean of Research 1 Code 012A Naval Postgraduate Schoo Monterey, CA 93940 Library, Code 5b 1 Naval Postgraduate School Monterey, CA 93940 Professor G. G. Brown 60 Code b5Bw Naval Postgraduate School Monterey, CA 9394U Defense Technical Information Center 2 Cameron Station Alexandria, VA 22314 63 DUDLEY KNOX LIBRARY 3 2768 00347374 5