1 files changed, 96 insertions, 0 deletions
diff --git a/Trivac/src/BIVPER.f b/Trivac/src/BIVPER.f
new file mode 100755
index 0000000..f4ba1bc
--- /dev/null
+++ b/Trivac/src/BIVPER.f
@@ -0,0 +1,96 @@
+*DECK BIVPER
+      SUBROUTINE BIVPER (JP,IDIR,LX,LT4,IP,IENV)
+*
+*-----------------------------------------------------------------------
+*
+*Purpose:
+* Compute the permutation vectors in 2-D hexagonal geometry.
+*
+*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. Benaboud
+*
+*Parameters: input
+* JP      first index.
+* IDIR    choice of direction (=1: W axis ; =2: X axis ; =3: Y axis).
+* LX      number of hexagons, including virtual hexagons.
+* LT4     number of non virtual hexagons.
+* IENV    index of non virtual hexagon corresponding to each hexagon.
+*
+*Parameters: output
+* IP      permutation vector.
+*
+*-----------------------------------------------------------------------
+*
+*----
+*  SUBROUTINE ARGUMENTS
+*----
+      INTEGER JP,IDIR,LX,LT4,IP(LX),IENV(LX)
+*----
+*  LOCAL VARIABLES
+*----
+      LOGICAL LPAS
+*
+      LPAS = .TRUE.
+      DO 10 I=1,LT4
+      IP(I)=0
+ 10   CONTINUE
+      NC = INT((SQRT(REAL((4*LX-1)/3))+1.)/2.)
+      IFACE1 = 0
+      IFACE2 = 0
+      IFACE3 = 0
+      IF(IDIR.EQ.1) THEN
+         IFACE1 = 3
+         IFACE2 = 5
+         IFACE3 = 4
+      ELSE IF(IDIR.EQ.2) THEN
+         IFACE1 = 4
+         IFACE2 = 6
+         IFACE3 = 5
+      ELSE IF(IDIR.EQ.3) THEN
+         IFACE1 = 5
+         IFACE2 = 1
+         IFACE3 = 6
+      ELSE IF(IDIR.EQ.5) THEN
+         IFACE1 = 1
+         IFACE2 = 3
+         IFACE3 = 2
+      ELSE
+         CALL XABORT('BIVPER: INVALID DATA')
+      ENDIF
+      JI = JP
+      JS = JP
+      KEL = 0
+      M = JI + 1
+      IF(IENV(JI).GT.0) THEN
+         IP(IENV(JI)) = 1
+         KEL = KEL + 1
+      ENDIF
+      IC = JP + NC - 1
+ 20   IF(KEL.EQ.LT4) RETURN
+      IF(JI.EQ.IC) IFACE2 = IDIR + 3
+ 30   IF(M.LE.LX) THEN
+         IF(IENV(M).GT.0) THEN
+            KEL = KEL + 1
+            IP(IENV(M)) = KEL
+         ENDIF
+         JI = M
+         M = NEIGHB(JI,IFACE1,9,LX,POIDS)
+         GOTO 30
+      ELSE
+ 40      JI = NEIGHB(JS,IFACE2,9,LX,POIDS)
+         IF(JI.GT.LX.AND.LPAS) THEN
+            IFACE2 = IFACE3
+            LPAS = .FALSE.
+            GO TO 40
+         ENDIF
+         M  = JI
+         JS = JI
+         GOTO 20
+      ENDIF
+      END