C      ALGORITHM 431, COLLECTED ALGORITHMS FROM ACM.
C      THIS WORK PUBLISHED IN COMMUNICATIONS OF THE ACM
C      VOL. 15, NO. 9, September, 1972, PP.818--820.
#! /bin/sh
# This is a shell archive, meaning:
# 1. Remove everything above the #! /bin/sh line.
# 2. Save the resulting text in a file.
# 3. Execute the file with /bin/sh (not csh) to create the files:
#	Fortran/
#	Fortran/Sp/
#	Fortran/Sp/Drivers/
#	Fortran/Sp/Drivers/Makefile
#	Fortran/Sp/Drivers/data
#	Fortran/Sp/Drivers/res
#	Fortran/Sp/Src/
#	Fortran/Sp/Src/src.f
# This archive created: Fri Mar 10 09:34:27 2006
export PATH; PATH=/bin:$PATH
if test ! -d 'Fortran'
then
	mkdir 'Fortran'
fi
cd 'Fortran'
if test ! -d 'Sp'
then
	mkdir 'Sp'
fi
cd 'Sp'
if test ! -d 'Drivers'
then
	mkdir 'Drivers'
fi
cd 'Drivers'
if test -f 'Makefile'
then
	echo shar: will not over-write existing file "'Makefile'"
else
cat << "SHAR_EOF" > 'Makefile'
all: Res

src.o: src.f
	$(F77) $(F77OPTS) -c src.f

DRIVERS= driver
RESULTS= Res

Objs1= src.o
driver: $(Objs1)
	$(F77) $(F77OPTS) -o driver $(Objs1) $(SRCLIBS)
Res: driver  data
	./driver < data >Res

diffres:Res res
	echo "Differences in results from driver"
	$(DIFF) Res res

clean: 
	rm -rf *.o $(DRIVERS) $(CLEANUP) $(RESULTS)
SHAR_EOF
fi # end of overwriting check
if test -f 'data'
then
	echo shar: will not over-write existing file "'data'"
else
cat << "SHAR_EOF" > 'data'
 1
 7
   0.00000   0.00000   0.00000   1.00000  -2.00000   0.00000   1.00000
   0.00000   0.00000   0.00000  -1.00000   0.00000   1.00000   4.00000
   0.00000   0.00000   0.00000   1.00000   3.00000   0.00000  -1.00000
  -1.00000   1.00000  -1.00000   0.00000   0.00000   0.00000   0.00000
   2.00000   0.00000  -3.00000   0.00000   0.00000   0.00000   0.00000
   0.00000  -1.00000   0.00000   0.00000   0.00000   0.00000   0.00000
  -1.00000  -4.00000   1.00000   0.00000   0.00000   0.00000   0.00000
   3.00000   2.00000   4.00000  -5.00000   1.00000   2.00000  -4.00000

SHAR_EOF
fi # end of overwriting check
if test -f 'res'
then
	echo shar: will not over-write existing file "'res'"
else
cat << "SHAR_EOF" > 'res'
1          PROBLEM NO. 1



     INITIAL ALMOST COMPLEMENTARY SOLUTION
          W(   1)=        8.00000
          W(   2)=        7.00000
          W(   3)=        9.00000
          W(   4)=        0.00000
          W(   5)=        6.00000
          W(   6)=        7.00000
          W(   7)=        1.00000
          Z0=        5.00000
     COMPLEMENTARY SOLUTION
          ITERATION NO.   7
          Z(   4)=        3.14286
          Z(   7)=        1.28571
          Z(   5)=        0.71429
          Z(   2)=        0.61905
          Z(   1)=        3.57143
          W(   6)=        2.61905
          Z(   3)=        2.04762
SHAR_EOF
fi # end of overwriting check
cd ..
if test ! -d 'Src'
then
	mkdir 'Src'
fi
cd 'Src'
if test -f 'src.f'
then
	echo shar: will not over-write existing file "'src.f'"
else
cat << "SHAR_EOF" > 'src.f'
C  REMARKS
C  SINCE THIS PROGRAM IS COMPLETE IN ALL RESPECTS, IT CAN BE
C  RUN AS IT IS WITHOUT ANY ADDITIONAL MODIFICATIONS OR
C  INSTRUCTIONS.  IN SUCH CASE, FOLLOW THE INPUT FORMAT AS GIVEN.
C
C  PROGRAM FOR SOLVING LINEAR AND QUADRATIC PROGRAMMING
C  PROBLEMS IN THE FORM W = M * Z + Q, W * Z = 0, W AND Z NONNEGATIVE
C  BY LEMKE'S ALGORITHM.
C
C  MAIN PROGRAM WHICH CALLS THE SIX SUBROUTINES - MATRIX,
C  INITIA, NEWBAS, SORT, PIVOT AND PPRINT IN PROPER ORDER.
C
      IMPLICIT NONE

      REAL A(50)
      REAL AM(50,50)
      REAL B(50,50)
      INTEGER IP
      INTEGER IR
      INTEGER L1
      INTEGER MBASIS(100)
      INTEGER N 
      INTEGER NE1
      INTEGER NE2
      INTEGER NL1
      INTEGER NL2
      INTEGER NO
      REAL Q(50)
      REAL W(50)
      REAL Z(50)

      COMMON AM,Q,L1,B,NL1,NL2,A,NE1,NE2,IR,MBASIS,W,Z
C
C  DESCRIPTION OF PARAMETERS IN COMMON.
C
C    AM     A TWO DIMENSIONAL ARRAY CONTAINING THE
C           ELEMENTS OF MATRIX M.
C
C    Q      A SINGLY SUBSCRIPTED ARRAY CONTAINING THE
C           ELEMENTS OF VECTOR Q.
C
C    L1     AN INTEGER VARIABLE INDICATING THE NUMBER OF
C           ITERATIONS TAKEN FOR EACH PROBLEM.
C
C    A      A TWO DIMENSIONAL ARRAY CONTAINING THE
C           ELEMENTS OF THE INVERSE OF THE CURRENT BASIS.
C
C    W      A SINGLY SUBSCRIPTED ARRAY CONTAINING THE VALUES
C           OF W VARIABLES IN EACH SOLUTION.
C
C    Z      A SINGLY SUBSCRIPTED ARRAY CONTAINING THE VALUES
C           OF Z VARIABLES IN EACH SOLUTION.
C
C    NL1    AN INTEGER VARIABLE TAKING VALUE 1 OR 2 DEPEND-
C           ING ON WHETHER VARIABLE W OR Z LEAVES THE BASIS.
C
C    NE1    SIMILAR TO NL1 BUT INDICATES VARIABLE ENTERING.
C
C    NL2    AN INTEGER VARIABLE INDICATING WHAT COMPONENT
C           OF W OR Z VARIABLE LEAVES THE BASIS.
C
C    NE2    SIMILAR TO NL2 BUT INDICATES VARIABLE ENTERING.
C
C    A      A SINGLY SUBSCRIPTED ARRAY CONTAINING THE
C           ELEMENTS OF THE TRANSFORMED COLUMN THAT IS
C           ENTERING THE BASIS.
C
C    IR     AN INTEGER VARIABLE DENOTING THE PIVOT ROW AT
C           EACH ITERATION.  ALSO USED TO INDICATE TERMINA-
C           TION OF A PROBLEM BY GIVING IT A VALUE OF 1000.
C
C    MBASIS A SINGLE SUBSCRIPTED ARRAY-INDICATOR FOR THE
C           BASIS VARIABLES.  TWO INDICATORS ARE USED FOR
C           EACH BASIC VARIABLE - ONE INDICATING WHETHER
C           IT IS A W OR Z AND ANOTHER INDICATING WHAT
C           COMPONENT OF W OR Z.
C
C  READ IN THE VALUE OF VARIABLE IP INDICATING THE
C  NUMBER OF PROBLEMS TO BE SOLVED.
C
      READ ( 5, 3 ) IP
C
C  VARIABLE NO INDICATES THE CURRENT PROBLEM BEING SOLVED.
C
      NO = 0
1     NO = NO + 1
      IF ( NO .GT. IP ) GO TO 5
      WRITE ( 6, 2 ) NO
C
C  READ IN THE SIZE OF THE MATRIX M.
C
2     FORMAT ( "1", 10X, "PROBLEM NO.", I2 )
      READ ( 5, 3 ) N
C
C  PROGRAM CALLING SEQUENCE.
C
3     FORMAT ( I2 )
      CALL MATRIX ( N )
C
C  PARAMETER N INDICATES THE PROBLEM SIZE.
C
      CALL INITIA ( N )
C
C  SINCE FOR ANY PROBLEM TERMINATION CAN OCCUR IN INITIA,
C  NEWBAS OR SORT SUBROUTINE, THE VALUE OF IR IS MATCHED WITH
C  1000 TO CHECK WHETHER TO CONTINUE OR TO TO NEXT PROBLEM.
C
      IF ( IR .EQ. 1000 ) GO TO 1
4     CALL NEWBAS ( N )
      IF ( IR .EQ. 1000 ) GO TO 1
      CALL SORT ( N )
      IF ( IR .EQ. 1000 ) GO TO 1
      CALL PIVOT ( N )
      GO TO 4
5     STOP
      END
      SUBROUTINE MATRIX ( N )
C
C  PURPOSE - TO INITIALIZE AND READ IN THE VARIOUS INPUT DATA.
C
      IMPLICIT NONE

      REAL A(50)
      REAL AM(50,50)
      REAL B(50,50)
      INTEGER I
      INTEGER IR
      INTEGER J
      INTEGER L1
      INTEGER MBASIS(100)
      INTEGER N
      INTEGER NE1
      INTEGER NE2
      INTEGER NL1
      INTEGER NL2
      REAL Q(50)
      REAL W(50)
      REAL Z(50)

      COMMON AM,Q,L1,B,NL1,NL2,A,NE1,NE2,IR,MBASIS,W,Z
C
C  READ THE ELEMENTS OF M MATRIX COLUMN BY COLUMN.
C
      DO 1 J = 1, N
        READ ( 5, 2 ) ( AM(I,J), I = 1, N )
1     CONTINUE
C
C  READ THE ELEMENTS OF Q VECTOR.
C
2     FORMAT ( 7F10.5 )
      READ ( 5, 2 ) ( Q(I), I = 1, N )
C
C  IN ITERATION 1, BASIS INVERSE IS AN IDENTITY MATRIX.
C
      DO 5 J = 1, N
        DO 4 I = 1, N
          IF ( I .EQ. J ) GO TO 3
          B(I,J) = 0.0E+00
          GO TO 4
3         B(I,J) = 1.0E+00
4       CONTINUE
5     CONTINUE

      RETURN
      END
      SUBROUTINE INITIA ( N )
C
C  PURPOSE - TO FIND THE INITIAL ALMOST COMPLEMENTARY SOLUTION
C  BY ADDING AN ARTIFICIAL VARIABLE Z0.
C
      IMPLICIT NONE

      REAL A(50)
      REAL AM(50,50)
      REAL B(50,50)
      INTEGER I
      INTEGER IR, IZR
      INTEGER J
      INTEGER L
      INTEGER L1
      INTEGER MBASIS(100)
      INTEGER N
      INTEGER NE1
      INTEGER NE2
      INTEGER NL1
      INTEGER NL2
      REAL Q(50)
      REAL T1
      REAL W(50)
      REAL Z(50)
      REAL Z0

      COMMON AM,Q,L1,B,NL1,NL2,A,NE1,NE2,IR,MBASIS,W,Z, IZR
C
C  SET Z0 EQUAL TO THE MOST NEGATIVE Q(I).
C
      I = 1
      J = 2
1     IF ( Q(I) .LE. Q(J) ) GO TO 2
      I = J
2     J = J + 1
      IF ( J .LE. N ) GO TO 1
C
C  UPDATE Q VECTOR
C
      IR = I
      T1 = -Q(IR)
      IF ( T1 .LE. 0.0E+00 ) GO TO 9

      DO 3 I = 1, N
        Q(I) = Q(I) + T1
3     CONTINUE

      Q(IR) = T1
C
C  UPDATE BASIS INVERSE AND INDICATOR VECTOR
C  OF BASIC VARIABLES.
C
      DO 4 J = 1, N
        B(J,IR) = -1.0E+00
        W(J) = Q(J)
        Z(J) = 0.0E+00
        MBASIS(J) = 1
        L = N + J
        MBASIS(L) = J
4     CONTINUE

      
      NL1 = 1
      L = N + IR
      NL2 = IR
      MBASIS(IR) = 3
      MBASIS(L) = 0
      W(IR) = 0.0E+00
      Z0 = Q(IR)
      L1 = 1
C
C  PRINT THE INITIAL ALMOST COMPLEMENTARY SOLUTION.
C
      WRITE ( 6, 5 )
5     FORMAT ( 3(/), 5X, "INITIAL ALMOST COMPLEMENTARY ",  "SOLUTION" )
      DO 7 I = 1, N
        WRITE ( 6, 6 ) I, W(I)
6       FORMAT ( 10X, "W(", I4, ")=", F15.5 )
7     CONTINUE
      WRITE ( 6, 8 ) Z0
8     FORMAT ( 10X, "Z0=", F15.5 )
      RETURN
9     WRITE ( 6, 10 )
10    FORMAT ( 5X, "PROBLEM HAS A TRIVIAL COMPLEMENTARY ",  "SOLUTION WI
     +TH W=Q, Z=0." )
      IR = 1000
      RETURN
      END
      SUBROUTINE NEWBAS ( N )
C
C  PURPOSE - TO FIND THE NEW BASIS COLUMN TO ENTER IN
C  TERMS OF THE CURRENT BASIS.
C
      IMPLICIT NONE

      REAL A(50)
      REAL AM(50,50)
      REAL B(50,50)
      INTEGER I
      INTEGER IR
      INTEGER J
      INTEGER L1
      INTEGER MBASIS(100)
      INTEGER N
      INTEGER NE1
      INTEGER NE2
      INTEGER NL1
      INTEGER NL2
      REAL Q(50)
      REAL TI
      REAL W(50)
      REAL Z(50)

      COMMON AM,Q,L1,B,NL1,NL2,A,NE1,NE2,IR,MBASIS,W,Z
C
C  IF NL1 IS NEITHER 1 NOR 2 THEN THE VARIABLE ZO LEAVES THE
C  BASIS INDICATING TERMINATION WITH A COMPLEMENTARY SOLUTION.
C
      IF ( NL1 .EQ. 1 ) GO TO 2
      IF ( NL1 .EQ. 2 ) GO TO 5
      WRITE ( 6, 1 )
1     FORMAT ( 5X, "COMPLEMENTARY SOLUTION" )
      CALL PPRINT ( N )
      IR = 1000
      RETURN
2     NE1 = 2
      NE2 = NL2
C
C  UPDATE NEW BASIS COLUMN BY MULTIPLYING BY BASIS INVERSE.
C
      DO 4 I = 1, N
        TI = 0.0E+00
        DO 3 J = 1, N
          TI = TI - B(I,J) * AM(J,NE2)
3     CONTINUE
        A(I) = TI
4     CONTINUE

      RETURN

5     NE1 = 1
      NE2 = NL2
      DO 6 I = 1, N
        A(I) = B(I,NE2)
6     CONTINUE

      RETURN
      END
      SUBROUTINE SORT ( N )
C
C  PURPOSE - TO FIND THE PIVOT ROW FOR THE NEXT ITERATION BY THE
C  USE OF (SIMPLEX-TYPE) MINIMUM RATIO RULE.
C
      IMPLICIT NONE

      INTEGER, PARAMETER :: NB = digits(0.0e0)-11
      REAL A(50)
      REAL AM(50,50)
      REAL B(50,50)
      INTEGER I
      INTEGER IR, IZR
      INTEGER L1
      INTEGER MBASIS(100)
      INTEGER N
      INTEGER NE1
      INTEGER NE2
      INTEGER NL1
      INTEGER NL2
      REAL T1, AMX, TOL
      REAL T2
      REAL Q(50)
      REAL W(50)
      REAL Z(50)

      COMMON AM,Q,L1,B,NL1,NL2,A,NE1,NE2,IR,MBASIS,W,Z, IZR

      AMX = ABS(A(1))
      DO 10 I = 2, N
	AMX = MAX(AMX, ABS(A(I)))
10    CONTINUE
      TOL = AMX*2.0**(-NB)

      I = 1
1     IF ( A(I) .GT. TOL ) GO TO 2
      I = I + 1
      IF ( I .GT. N ) GO TO 9
      GO TO 1

2     T1 = Q(I) / A(I)
      IR = I

3     I = I + 1
      IF ( I .GT. N ) GO TO 5
      IF ( A(I) .GT. TOL ) GO TO 4
      GO TO 3

4     T2 = Q(I) / A(I)
      IF ( T2 .GE. T1 ) GO TO 3
      IR = I
      T1 = T2
      GO TO 3

5     RETURN

9     IF(Q(IZR).GT.TOL) GO TO 6
      WRITE(6,11)
11    FORMAT(5X,'COMPLEMENTARY SOLUTION')
      CALL PPRINT(N)
      IR = 1000
      RETURN
C
C  FAILURE OF THE RATIO RULE INDICATES TERMINATION WITH
C  NO COMPLEMENTARY SOLUTION.
C
6     WRITE ( 6, 7 )
7     FORMAT ( 5X, "PROBLEM HAS NO COMPLEMENTARY SOLUTION" )
      WRITE ( 6, 8 ) L1
8     FORMAT ( 10X, "ITERATION NO.", I4 )
      IR = 1000
      RETURN
      END
      SUBROUTINE PIVOT ( N )
C
C  PURPOSE - TO PERFORM THE PIVOT OPERATION BY UPDATING THE
C  INVERSE OF THE BASIS AND Q VECTOR.
C
      IMPLICIT NONE

      REAL A(50)
      REAL AM(50,50)
      REAL B(50,50)
      INTEGER I
      INTEGER IR
      INTEGER J
      INTEGER L
      INTEGER L1
      INTEGER MBASIS(100)
      INTEGER N
      INTEGER NE1
      INTEGER NE2
      INTEGER NL1
      INTEGER NL2
      REAL Q(50)
      REAL W(50)
      REAL Z(50)

      COMMON AM,Q,L1,B,NL1,NL2,A,NE1,NE2,IR,MBASIS,W,Z

      DO 1 I = 1, N
        B(IR,I) = B(IR,I) / A(IR)
1     CONTINUE

      Q(IR) = Q(IR) / A(IR)

      DO 3 I = 1, N
        IF ( I .EQ. IR ) GO TO 3
        Q(I) = Q(I) - Q(IR) * A(I)
        DO 2 J = 1, N
          B(I,J) = B(I,J) - B(IR,J) * A(I)
2       CONTINUE
3     CONTINUE
C
C  UPDATE THE INDICATOR VECTOR OF BASIC VARIABLES.
C
      NL1 = MBASIS(IR)
      L = N + IR
      NL2 = MBASIS(L)
      MBASIS(IR) = NE1
      MBASIS(L) = NE2
      L1 = L1 + 1

      RETURN
      END
      SUBROUTINE PPRINT ( N )
C
C  PURPOSE - TO PRINT THE CURRENT SOLUTION TO COMPLEMENTARY
C  PROBLEM AND THE ITERATION NUMBER.
C
      IMPLICIT NONE

      REAL A(50)
      REAL AM(50,50)
      REAL B(50,50)
      INTEGER I
      INTEGER IR
      INTEGER J
      INTEGER K1
      INTEGER K2
      INTEGER L1
      INTEGER MBASIS(100)
      INTEGER N
      INTEGER NE1
      INTEGER NE2
      INTEGER NL1
      INTEGER NL2
      REAL Q(50)
      REAL W(50)
      REAL Z(50)

      COMMON AM,Q,L1,B,NL1,NL2,A,NE1,NE2,IR,MBASIS,W,Z

      WRITE ( 6, 1 ) L1
1     FORMAT ( 10X, "ITERATION NO.", I4 )
      I = N + 1
      J = 1

2     K1 = MBASIS(I)
      K2 = MBASIS(J)
      IF ( Q(J) .GE. 0.0E+00 ) GO TO 3
      Q(J) = 0.0E+00

3     IF ( K2 .EQ. 1 ) GO TO 5
      WRITE ( 6, 4 ) K1, Q(J)
4     FORMAT ( 10X, "Z(", I4, ")=", F15.5 )
      GO TO 7

5     WRITE ( 6, 6 ) K1, Q(J)
6     FORMAT ( 10X, "W(", I4, ")=", F15.5 )

7     I = I + 1
      J = J + 1
      IF ( J .LE. N ) GO TO 2

      RETURN
      END
SHAR_EOF
fi # end of overwriting check
cd ..
cd ..
cd ..
#       End of shell archive
exit 0