diff options
Diffstat (limited to 'Utilib/src/PLLEMK.f')
| -rw-r--r-- | Utilib/src/PLLEMK.f | 165 |
1 files changed, 165 insertions, 0 deletions
diff --git a/Utilib/src/PLLEMK.f b/Utilib/src/PLLEMK.f new file mode 100644 index 0000000..d50db6f --- /dev/null +++ b/Utilib/src/PLLEMK.f @@ -0,0 +1,165 @@ +*DECK PLLEMK
+ SUBROUTINE PLLEMK(N,M,EPS,IMPR,P,IROW,ICOL,IERR)
+*
+*-----------------------------------------------------------------------
+*
+*Purpose:
+* This subroutine solve the parametric linear complementary problem.
+* PLLEMK = Linear Programming LEMKe
+*
+*Copyright:
+* Copyright (C) 2002 Ecole Polytechnique de Montreal
+* This library is free software; you can redistribute it and/or
+* modify it under the terms of the GNU Lesser General Public
+* License as published by the Free Software Foundation; either
+* version 2.1 of the License, or (at your option) any later version
+*
+*Author(s):
+* A. Hebert and T. Falcon
+*
+*Parameters: input
+* N number of control variables.
+* M number of constraints.
+* EPS tolerence used for pivoting.
+* IMPR print flag.
+* P coefficient matrix.
+* IROW permutation vector for row elements.
+* ICOL permutation vector for column elements.
+*
+*Parameters: ouput
+* P coefficient matrix.
+* IROW permutation vector for row elements.
+* ICOL permutation vector for column elements.
+* IERR return code (=0: normal completion).
+*
+*-----------------------------------------------------------------------
+*
+ IMPLICIT NONE
+*----
+* SUBROUTINE ARGUMENTS
+*----
+ INTEGER M,N,IERR,IMPR,IROW(N),ICOL(N+1)
+ DOUBLE PRECISION EPS,P(N,M)
+*----
+* LOCAL VARIABLES
+*----
+ DOUBLE PRECISION WWW,S
+ INTEGER MNOP,NP1,IC,L,I,J,K,JJ,J1,LGAR
+ DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: PLJ
+*----
+* SCRATCH STORAGE ALLOCATION
+*----
+ ALLOCATE(PLJ(M))
+*
+ MNOP = N + M
+ NP1 = N + 1
+ IC = 0
+*
+ S = 0.0D0
+ L = 0
+*
+ DO I=1,N
+ P(I,NP1) = 1.0D0
+ IF(P(I,M).GE.0.0D0) CYCLE
+ WWW = P(I,M)
+ IF(WWW.LE.S) THEN
+ L = I
+ S = WWW
+ ENDIF
+ ENDDO
+*
+ IF(L.EQ.0) GO TO 30
+ LGAR = L
+*
+ K = NP1
+*
+ 10 IC = IC + 1
+ IF(IC.GT.MNOP) THEN
+ IERR = 1
+ GOTO 40
+ ENDIF
+*
+ IF(IMPR.GE.5) WRITE (6,1000) L,K
+ IF(IMPR.GE.6) THEN
+ DO J=1,M,12
+ JJ = MIN0(J+11,M)
+ WRITE (6,3000) (J1,J1=J,JJ)
+ DO I=1,N
+ WRITE (6,4000) I,(P(I,J1),J1=J,JJ)
+ ENDDO
+ IF(JJ.NE.M) WRITE (6,5000)
+ ENDDO
+ ENDIF
+*
+ CALL PLPIVT(N,M,L,K,P,IROW,ICOL)
+*
+ IF(ICOL(K).EQ.-NP1) GOTO 30
+ DO J=1,NP1
+ IF(J.EQ.K) CYCLE
+ IF(IABS(ICOL(K)).EQ.IABS(ICOL(J))) THEN
+ K = J
+ GOTO 20
+ ENDIF
+ ENDDO
+*
+ IERR = 2
+ GOTO 40
+*
+ 20 DO I=1,N
+ IF(P(I,K)/P(I,M).LT.-EPS) THEN
+ PLJ(I) = -P(I,M)/P(I,K)
+ ELSE
+ PLJ(I) = 1.0D50
+ ENDIF
+ ENDDO
+*
+ S = PLJ(LGAR)
+ L = LGAR
+*
+ IF(IMPR.GE.7) THEN
+ WRITE (6,*) 'K=',K
+ WRITE (6,6000) L,(PLJ(I),I=1,N)
+ ENDIF
+*
+ DO I=1,N
+ WWW = PLJ(I)
+ IF((ABS(WWW-S).GT.ABS(S)*EPS).AND.(WWW.LT.S)) THEN
+ S = WWW
+ L = I
+ ENDIF
+ ENDDO
+*
+ IF(S.EQ.1.0D50) THEN
+ IERR = 3
+ GOTO 40
+ ENDIF
+*
+ GOTO 10
+ 30 IERR = 0
+*
+ 40 IF(IMPR.GE.6) THEN
+ WRITE (6,2000) IC
+ DO J=1,M,12
+ JJ = MIN0(J+11,M)
+ WRITE (6,3000) (J1,J1=J,JJ)
+ DO I=1,N
+ WRITE (6,4000) I,(P(I,J1),J1=J,JJ)
+ ENDDO
+ IF(JJ.NE.M) WRITE (6,5000)
+ ENDDO
+ ENDIF
+*----
+* SCRATCH STORAGE DEALLOCATION
+*----
+ DEALLOCATE(PLJ)
+ RETURN
+*
+1000 FORMAT (//,10X,'+ + + + + MATRIX P(N,N+3) :',5X,'PIVOT = (',
+ > I4,' ,',I4,' )'/)
+2000 FORMAT (//,10X,'NUMBER OF PIVOTS =',I5,
+ > //,10X,'+ + + + + FINAL MATRIX P(N,N+3) :'/)
+3000 FORMAT (5X,12I12)
+4000 FORMAT (1X,I4,1P,10E12.4/(5X,1P,10E12.4))
+5000 FORMAT (//)
+6000 FORMAT (//,10X,'LGAR =',I4/(1X,1P,8E12.4))
+ END
|