100' NAME--LINQD1 110' 120' DESCRIPTION--LINEAR OR QUADRATIC PROGRAMMING ALGORITHM 130' 140' SOURCE--M.P. JUHA AND A.O. CONVERSE, JUNE 1968, THAYER SCHOOL 150' OF ENGINEERING 160' 170' INSTRUCTIONS--INSERT DATA AS INDICATED BELOW BETWEEN LINES 3160 AND 3370 180' 190' 200' 1) PROBLEM TYPE: MAX OR MIN 210' 2) NUMBER OF TERMS IN OBJECTIVE FUNCTION 220' E.G., THE EXPRESSION: X1*X2 - 2*X1 - X2 230' CONSISTS OF THREE TERMS 240' 3) NUMBER OF ALLOCATION VARIABLES 250' E.G., THE EXPRESSION: X1*X2 - 2*X1 - X2 260' CONTAINS TWO ALLOCATION VARIABLES: X1 AND X2 270' 4) NUMBER OF CONSTRAINT RELATIONS ON PROBLEM 280' 5) THE TERMS OF THE OBJECTIVE FUNCTION EXPRESSED IN THE 290' FOLLOWING FORMAT: 300' [ENTRY 1],[ENTRY 2],[COEFFICIENT] 310' WHERE: ENTRY 1 IS THE NUMBER OF AN ALLOCATION VARIABLE 320' AND ENTRY 2 IS THE NUMBER OF AN ALLOCATION VARIABLE 330' E.G., THE TERM: X1*X2 WOULD BE EXPRESSED AS 340' 1,2,1 350' AND THE TERM: - 2*X1 WOULD BE 1,0,-2 OR 0,1,-2 360' [THE ZERO SIGNIFIES THAT THE TERM IS LINEAR, THAT IS 370' X0=1 BY DEFINITION] 380' 6) THE COEFFICIENTS OF THE CONSTRAINT RELATIONS WRITTEN 390' AS FOLLOWS: FOR BOTH MINIMIZATION AND MAXIMIZATION 400' PROBLEMS THE CONSTRAINTS MUST ALL BE OF THE FORM: 410' H[X1,X2,...,XN]>=0 420' AS AN EXAMPLE, FOR A MINIMIZATION PROBLEM WHERE ONE 430' OF THE CONSTRAINT RELATIONS IS 440' - 2*X1 - 3*X2 + 6 >= 0 450' THE INPUT DATA WOULD CONSIST OF THE COEFFICIENTS OF 460' X1 THROUGH XN [N=2 FOR THIS CASE] FOLLOWED BY THE 470' CONSTANT TERM. E.G., THE INPUT DATA WOULD BE: 480' -2, -3, 6 490' NOTE FOR NOMENCLATURE CONTINUE LISTING THROUGH LINE 760. 500' AN EXAMPLE IS PROVIDED FOR THE USER STARTING IN LINE 3020. 510' 520' 530' * * * * * * * MAIN PROGRAM * * * * * * * * * * * 540' 550 DIM A(51,51),M(50),Y(50) 560' *****N0MENCLATURE***** 570' T$= PROBLEM TYPE 580' W = NUMBER OF TERMS IN OBJECTIVE FUNCTION 590' Q = NUMBER OF ALLOCATION VARIABLES X(I) 600' U = NUMBER OF CONSTRAINT RELATIONS 610' R = NUMBER OF ROWS IN TABLEAU 620' C = NUMBER OF COLUMNS IN TABLEAU 630' A = WORKING MATRIX 640' M = VECTOR OF VALUES OF LAGRANGE MULTIPLIERS 650' Y = VECTOR OF VALUES OF ALLOCATION VARIABLES AND CONSTRAINTS 660' Z9 = A FLAG. IT IS SET > 1 BY THE USER WHEN HE WISHES TO STOP 670' MANUAL PIVOTING AND PRINT FINAL RESULTS 680' Z8 = A FLAG. IT IS SET > 1 BY THE USER WHEN HE WISHES TO 690' SPECIFY THE PIVOT ELEMENT MANUALLY. 700' Z7 = A FLAG. IT IS SET = 5 BY THE PROGRAM WHEN PIVOTING ON THE 710' TRANSPOSE OF THE ORIGINAL PIVOT ELEMENT. 720' Z6 = A FLAG. IT IS SET = 1 BY THE PROGRAM WHEN THE PROBLEM IS 730' LINEAR. 740' Z5 = A FLAG. IT IS SET > 1 BY THE USER WHEN PRINTOUT OF THE 750' INTERMEDIATE TABLEAUS IS DESIRED. 760' 770' *******DATA INPUT******* 780 LET Z7 = 0 790 READ T$,W,Q,U 800 LET Z6 = 0 810 PRINT"OBJECTIVE IS TO ";T$;"IMIZE THE SUM OF THE FOLLOWING TERMS:" 820 PRINT " " 830 LET R=C=U+Q+1'COMPUTE SIZE OF TABLEAU 840 MAT A = ZER(R,C)'CREATE AN EMPTY TABLEAU 850 LET F=-1 'SET FLAG TO INDICATE QUADRATIC PROGRAMMING 860' 870 FOR I = 1 TO Q 880 LET V$(I) = "X" 890 NEXT I 900 FOR I = Q+1 TO Q+U 910 LET V$(I) = "MU" 920 NEXT I 930 LET V$(Q+U+1) = "1" 940 FOR K=1 TO W 'READ IN OBJECTIVE FUNCTION TERM-BY-TERM 950 READ I,J,A(I,J) 960 IF I=0 THEN 1060 'IS TERM LINEAR? YES, GO TO 400 970 IF J=0 THEN 1060 'IS TERM LINEAR? YES, GO TO 400 980 PRINT " ";A(I,J);"* X(";I;") * X(";J;")" 990 IF I=J THEN 1020 'CHECK THAT TERM IS A SQUARE TERM 1000 LET A(J,I) = A(I,J) 1010 GO TO 1030' THE CURRENT TERM IS NOT IN THE DIAGONAL 1020 LET A(J,I)=A(I,J)+A(J,I) 'ASSURE THAT DIAGONAL TERMS ARE TIMES 2 1030 LET F=0 'INDICATE THAT A QUADRATIC TERM HAS BEEN FOUND 1040 GO TO 1080 'GO ON TO NEXT TERM IN OBJECTIVE FUNCTION 1050'READ IN COEFFICIENT OF TERM LINEAR IN X 1060 LET A(I+J,C)=A(R,I+J)=A(I,J) 1070 PRINT " ";A(I,J);"* X(";I+J;")" 1080 NEXT K 'READ IN NEXT TERM OF OBJECTIVE FUNCTION 1090 PRINT " " 1100 IF T$ = "MIN" THEN 1180 1110 FOR I9 = 1 TO Q 1120 FOR I8 = 1 TO Q+U+1 1130 LET A(I9,I8) = -A(I9,I8) 1140 NEXT I8 1150 NEXT I9 1160' 1170' 1180 PRINT "THE CONSTRAINTS ARE:" 1190'READ IN CONSTRAINT RELATIONS 1200 FOR I=1 TO U 1210 PRINT " " 1220 PRINT " H";I;" =:" 1230 FOR J=1 TO Q 1240 READ A(Q+I,J) 'READ COEFFICIENT TO X(J) 1250 PRINT " ";A(Q+I,J);"* X(";J;")" 1260 LET A(J,Q+I)=-A(Q+I,J) 'CREATE TRANSPOSE TERM 1270 NEXT J 1280 READ A(Q+I,C) 'READ CONSTANT TERM IN CONSTRAINT 1290 PRINT " ";A(Q+I,C) 1300 IF T$ = "MIN" THEN 1330 1310 LET A(R,Q+I) = A(Q+I,C) 1320 GO TO 1340 1330 LET A(R,Q+I)=-A(Q+I,C) 1340 PRINT " >= 0" 1350 NEXT I 1360 PRINT " " 1370 PRINT " " 1380 PRINT " " 1390' 1400 IF F = 0 THEN 1440' SKIP TO 660 IF THIS IS A Q.P. PROBLEM 1410 PRINT "LINEAR PROGRAMMING" 1420 PRINT" " 1430 LET Z6 = 1' SIGNIFIES A LINEAR PROGRAM 1440 PRINT"DO YOU WISH TO SELECT THE PIVOT ELEMENTS MANUALLY,YES OR NO" 1450 INPUT T1$ 1460 IF T1$ = "NO" THEN 1490 1470 LET Z8 = 2 1480 GO TO 1520 1490 PRINT"DO YOU WISH PRINTOUT OF INTERMEDIATE TABLEAUS,YES OR NO" 1500 INPUT T2$ 1510 IF T2$ = "NO" THEN 1530 1520 LET Z5 = 2 1530' 1540'END DATA READ IN PROCEDURE 1550' 1560'INITIALIZE EXTREMUM STORAGE MATRIX 1570 IF T$="MAX" THEN 1610 1580'PROBLEM IS MIN 1590 LET B(R,C)= 1 E 27'SET OBJECTIVE COMPARISON VERY LARGE 1600 GO TO 1660'CONTINUE ON TO SET UP EXCHANGE INDICES 1610'PROBLEM IS MAX 1620 LET B(R,C)=-1 E 27'SET OBJECTIVE COMPARISON VERY SMALL 1630' 1640'SET EXCHANGE VARIABLE INDICES, A(0,I) 1650'INDICES STORE A RECORD OF WHICH VARIABLE IS EXCHANGED 1660 FOR I = 1 TO Q 1670 LET A (0,I) = 100+I' IDENTIFIES THE ALLOCATION VARIABLES 1680 NEXT I 1690 FOR I = Q+1 TO Q+U 1700 LET A(0,I) = 200 + I' IDENTIFIES L.M. ASSOCIATED WITH CONSTRAINTS 1710 NEXT I 1720 FOR J = 1 TO R-1 1730 LET A(J,0)=J 1740 NEXT J 1750' 1760 IF Z5 < 1 THEN 1830' NO INTERMEDIATE PRINTOUT 1770 GO SUB 2090' TO PRINT EXPLANATION OF INITIAL TABLEAU 1780 MAT PRINT V$, 1790 IF Z5 < 1 THEN 1830' NO INTERMEDIATE PRINTOUT 1800' 1810 MAT PRINT A 1820 IF Z8 > 1 THEN 2170' TO SELECT PIVOT MANUALLY 1830 LET P = 1E27' SET FLAG AND SEARCH FOR PIVOT 1840 FOR I = 1 TO R-1 1850 IF A(I,C) >= 0 THEN 1960 1860 IF Z6 = 1 THEN 1880 1870 GO SUB 1990' Q.P. SEARCH FOR PIVOT 1880 FOR J = 1 TO C-1 1890 IF A(I,J) <=0 THEN 1950 1900 LET X = A(R,J)/A(I,J) 1910 IF X >=P THEN 1950 1920 LET P = X 1930 LET K = I 1940 LET L = J 1950 NEXT J 1960 NEXT I 1970 IF P < 1E27 THEN 2210 1980 GO TO 2650' NO PIVOT 1990' *** Q.P. SEARCH FOR PIVOT SUBROUTINE *** 2000' 2010 FOR I = 1 TO R-1 2020 IF A(I,0) > 100 THEN 2060' HAS VARIABLE BEEN EXCHANGED 2030 IF A(I,I) =0 THEN 2060' IS DIAGONAL ELEMENT EMPTY 2040 LET K = L = I' NOTE PIVOT LOCATION 2050 GO TO 2210' GO TO EXCHANGE ALGORITHM 2060 NEXT I 2070 RETURN 2080' 2090' PRINT SUBROUTINE 2100 PRINT" IN THE FIRST TABLEAU, THE FIRST";Q;" COLUMNS ARE ASSOCIATED" 2110 PRINT" WITH THE CORRESPONDING";Q;"ALLOCATION VARIABLES. THE NEXT" 2120 PRINT U;"COLUMNS ARE ASSOCIATED WITH THE ";U;"LAGRANGE MULTTIPLIERS" 2130 PRINT" ASSOCIATED WITH THE CORRESPONDING CONSTRAINTS. " 2140 RETURN 2150' 2160 IF Z8 < 1 THEN 2210' USER WILL NOT CHOOSE PIVOT 2170 PRINT " INPUT PIVOT LOCATION AND Z9 " 2180 INPUT K,L,Z9 2190 IF Z9 > 1 THEN 2650' USER HAS DECIDED TO STOP AND PRINT RESULTS 2200' 2210' ***** EXCHANGE ALGORITHM ***** 2220 PRINT "PIVOT LOCATION";K;L 2230 LET X=A(0,L) 2240 LET A(0,L)=A(K,0) 2250 LET A(K,0)=X 2260 LET X=A(K,L) 2270 REM EXCHANGE ALL ENTRIES NOT IN PIVOT ROW OR COLUMN 2280 FOR I=1 TO R 2290 IF I=K THEN 2340 2300 FOR J=1 TO C 2310 IF J=L THEN 2330 2320 LET A(I,J)=A(I,J)-A(I,L)*A(K,J)/X 2330 NEXT J 2340 NEXT I 2350 REM CALCULATE PIVOT COLUMN VALUES 2360 FOR I=1 TO R 2370 LET A(I,L)=A(I,L)/X 2380 NEXT I 2390 REM CALCULATE PIVOT ROW VALUES 2400 FOR J=1 TO C 2410 LET A(K,J)=-A(K,J)/X 2420 NEXT J 2430 REM CALCULATE PIVOT VALUE 2440 LET A(K,L)=1/X 2450'FINISHED EXCHANGING VARIABLES 2460 IF K = L THEN 2540' SKIP TO 1625 IF PIVOT IS IN DIAGONAL 2470 IF Z7 = 5 THEN 2540' GO TO 1625 IF TRANSPOSED PIVOT HAS BEEN DONE 2480 LET Z7 = K' TEMPORARY STORAGE 2490 LET K = L' TRANSPOSING 2500 LET L = Z7 2510 LET Z7 = 5' FLAG TO SHOW WHEN TRANSPOSED PIVOT HAS BEEN DONE 2520 GO TO 2210' TO PERFORM TRANSPOSE PIVOT 2530' 2540 LET Z7 = 0' FLAG 2550 GO TO 1790 2560' 2570' 2580'PROBLEM END 2590' 2600'OUTPUT RESULTS 2610' 2620'DETERMINE VALUES OF ALLOCATION VARIABLES, LAGRANGE MULTIPLIERS 2630' AND CONSTRAINT RELATIONS. 2640' 2650 MAT M=ZER(Q+U) 'SET UP LAGRANGE MULTIPLIER VECTOR 2660 MAT Y=ZER(Q+U) 'SET UP ALLOCATION VARIABLE/CONSTRAINT VECTOR 2670' 2680 FOR I=1 TO R-1 2690 IF A(I,0) > 200 THEN 2780 2700 IF A(I,0)>100 THEN 2760 2710 IF A(I,0) > Q THEN 2740 2720 LET M(A(I,0)) = A(I,C) 2730 GO TO 2850 'GO ON TO NEXT ROW IN TABLEAU 2740 LET Y(A(I,0)) = A(I,C) 2750 GO TO 2850 2760 LET Y(A(I,0) -100) = A(I,C) 2770 GO TO 2850 'GO ON TO NEXT ROW IN TABLEAU 2780 LET M(A(I,0)-200) = A(I,C) 2790 GO TO 2850 2800'ROW OF TABLEAU IS IN CONSTRAINT RELATION PARTITION 2810 IF A(I,0)<100 THEN 2840 2820 LET M(I)=A(I,C) 'LAGRANGE MULTIPLIER, SAVE IT 2830 GO TO 2850 'GO ON TO NEXT ROW IN TABLEAU 2840 LET Y(I)=A(I,C) 'SAVE CONSTRAINT VALUE 2850 NEXT I 'GO TO NEXT ROW IN TABLEAU 2860' 2870 PRINT"ALLOCATION VARIABLES AND ASSOCIATED LAGRANGE MULTIPLIERS" 2880 PRINT " " 2890 FOR I=1 TO Q 2900 PRINT "X";I;"=";Y(I),"MU";I;"=";M(I) 2910 NEXT I 2920 PRINT " " 2930 PRINT "CONSTRAINT RELATIONS AND ASSOCIATED LAGRANGE MULTIPLIERS" 2940 PRINT " " 2950 FOR I=Q+1 TO Q+U 2960 PRINT "H";I-Q;"=";Y(I),"MU";I;"=";M(I) 2970 NEXT I 2980' 2990 PRINT " " 3000 PRINT "OBJECTIVE FUNCTION =";A(R,C)/2 3010' 3020'*****EXAMPLE PROBLEM***** 3030' 3040' MINIMIZE: 2*X1^2 + 2*X2^2 -X1*X2 -10*X1 -20*X2 3050' 3060' BY SELECTION OF: X1, X2 3070' 3080' SUBJECT TO: 3090' 3100' H1 = 8 - X2 >= 0 3110' 3120' H2 = 10 - X1 - X2 >= 0 3130' 3140' INPUT DATA FOR ALGORITHM 3150' 3160'PROBLEM TYPE 3170 DATA MIN 3180'NUMBER OF TERMS IN OBJECTIVE FUNCTION 3190 DATA 5 3200'NUMBER OF ALLOCATION VARIABLES 3210 DATA 2 3220'NUMBER OF CONSTRAINT RELATIONS 3230 DATA 2 3240'COEFFICIENTS OF OBJECTIVE FUNCTION---TERM-BY-TERM 3250' FIRST TERM 3260 DATA 1,1,2 3270' SECOND TERM 3280 DATA 2,2,2 3290' THIRD TERM 3300 DATA 1,2,-1 3310 DATA 1,0,-10 3320 DATA 2,0,-20 3330' FIRST CONSTRAINT 3340 DATA 0,-1,8 3350' SECOND CONSTRAINT 3360 DATA -1,-1,10 3370'*****END OF DATA ***** 3380 END