From 7dfcc480ba1e19bd3232349fc733caef94034292 Mon Sep 17 00:00:00 2001 From: stainer_t Date: Mon, 8 Sep 2025 13:48:49 +0200 Subject: Initial commit from Polytechnique Montreal --- Utilib/src/PLLEMK.f | 165 ++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 165 insertions(+) create mode 100644 Utilib/src/PLLEMK.f (limited to 'Utilib/src/PLLEMK.f') 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 -- cgit v1.2.3