Initial commit from Polytechnique Montreal

1 files changed, 785 insertions, 0 deletions

diff --git a/Dragon/src/LELCSY.f b/Dragon/src/LELCSY.f
new file mode 100644
index 0000000..1aedfa0
--- /dev/null
+++ b/Dragon/src/LELCSY.f
@@ -0,0 +1,785 @@
+*DECK LELCSY
+      FUNCTION LELCSY(IPGEOM,IPRT   ,GEONAM,GEOCV ,NXYZ  ,IGT   ,
+     >                GMESH ,ISPLG ,TMESH ,ISPLT )
+*
+*-----------------------------------------------------------------------
+*
+*Purpose:
+* Verify if geometry satisfies an intrinsic symmetry.
+*
+*Copyright:
+* Copyright (C) 2003 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): R. Roy and G. Marleau
+*
+*Parameters: input
+* IPGEOM  pointer to the reference geometry data structure.
+* IPRT    intermediate printing level for output.
+* GEONAM  name of the reference geometry.
+* GEOCV   name of the geometry to analyze if different from 
+*         reference geometry.
+* NXYZ    maximum mesh size in directions $x$, $y$ and $z$.
+* IGT     geometry turn number.
+*
+*Parameters: scratch
+* GMESH   general mesh description for geometry to analyze.
+* ISPLG   general split desctiption for geometry to analyze.
+* TMESH   temporary mesh description for geometry comparison.
+* ISPLT   temporary split desctiption for geometry comparison.
+* LELCSY  result of geometry testing (true or false).
+*
+*-----------------------------------------------------------------------
+*
+      USE                GANLIB
+      IMPLICIT           NONE
+*----
+*  SUBROUTINE ARGUMENTS
+*----
+      LOGICAL            LELCSY
+      TYPE(C_PTR)        IPGEOM
+      INTEGER            IPRT,NXYZ
+      CHARACTER*12       GEONAM,GEOCV 
+      INTEGER            IGT(2)
+      INTEGER            ISPLG(0:NXYZ-1)
+      INTEGER            ISPLT(3,3,0:NXYZ-1)
+      REAL               GMESH(0:NXYZ)
+      DOUBLE PRECISION   TMESH(3,3,0:NXYZ)
+*----
+*  LOCAL VARIABLES
+*----
+      INTEGER            IOUT,NSTATE
+      CHARACTER          NAMSBR*6
+      PARAMETER         (IOUT=6,NSTATE=40,NAMSBR='LELCSY')
+*----
+*  ALLOCATABLE ARRAYS
+*----
+      INTEGER, ALLOCATABLE, TARGET, DIMENSION(:) :: IMIX,IMIX1,IM1,IM2
+      INTEGER, POINTER, DIMENSION(:) :: IMIX2
+*----
+*  LOCAL PARAMETERS
+*----
+      INTEGER            ISTATE(NSTATE),IKT(2)
+      INTEGER            NR,NX,NY,NZ,NK,
+     >                   NTM(2,3)
+      REAL               OFFCEN(3),OFFTR(3,2)
+      INTEGER            ILCMLN,ILCMTY
+      INTEGER            IG,IDIR,ITMI,IDMI,
+     >                   IR,IX,IY,IZ
+      DOUBLE PRECISION   DDM(3),DDD 
+      INTEGER            ITYPG,IMTN,IMTO
+      LELCSY=.TRUE. 
+      IF(IPRT .GE. 10) THEN 
+        IF(GEOCV .EQ. ' ') THEN
+          WRITE(IOUT,6000) NAMSBR,GEONAM,IGT(1),IGT(2)
+        ELSE
+          WRITE(IOUT,6000) NAMSBR,GEOCV,IGT(1),IGT(2)
+        ENDIF
+      ENDIF
+*----
+*  READ GEOMETRY INFORMATION
+*----
+      IF(GEOCV .NE. '            ') THEN  
+        CALL LCMSIX(IPGEOM,GEOCV,1) 
+      ENDIF
+*----
+*  STATE-VECTOR
+*----
+      ISTATE(:NSTATE)=0
+      CALL LCMGET(IPGEOM,'STATE-VECTOR',ISTATE)
+      ITYPG= ISTATE(1)
+      NR=ISTATE(2)
+      NX=ISTATE(3)
+      NY=ISTATE(4)
+      NZ=ISTATE(5)
+      NK=ISTATE(6)
+      TMESH(:3,:3,0:NXYZ)=0.0D0
+      ISPLT(:3,:3,0:NXYZ-1)=1
+*----
+*  MESHX AND SPLITX
+*----
+      IDIR=1
+      GMESH(0:NXYZ)=0.0
+      ISPLG(0:NXYZ-1)=1
+      CALL LCMLEN(IPGEOM,'MESHX',ILCMLN,ILCMTY)
+      IF(ILCMLN .EQ. NX+1)
+     >  CALL LCMGET(IPGEOM,'MESHX',GMESH)
+      CALL LCMLEN(IPGEOM,'SPLITX',ILCMLN,ILCMTY)
+      IF(ILCMLN .EQ. NX)
+     >  CALL LCMGET(IPGEOM,'SPLITX',ISPLG)
+      DDM(IDIR)=0.0D0 
+      DO 10 IX=0,NX-1
+        TMESH(3,IDIR,IX)=DBLE(GMESH(IX))
+        DDM(IDIR)=MAX(DDM(IDIR),ABS(TMESH(3,IDIR,IX)))
+        ISPLT(3,IDIR,IX)=ISPLG(IX)
+ 10   CONTINUE 
+      TMESH(3,IDIR,NX)=DBLE(GMESH(NX))
+      DDM(IDIR)=MAX(DDM(IDIR),ABS(TMESH(3,IDIR,NX)))
+*----
+*  MESHY AND SPLITY
+*----
+      IDIR=2
+      GMESH(0:NXYZ)=0.0
+      ISPLG(0:NXYZ-1)=1
+      CALL LCMLEN(IPGEOM,'MESHY',ILCMLN,ILCMTY)
+      IF(ILCMLN .EQ. NY+1)
+     >  CALL LCMGET(IPGEOM,'MESHY',GMESH) 
+      CALL LCMLEN(IPGEOM,'SPLITY',ILCMLN,ILCMTY)
+      IF(ILCMLN .EQ. NY)
+     >  CALL LCMGET(IPGEOM,'SPLITY',ISPLG) 
+      DDM(IDIR)=0.0D0 
+      DO 11 IY=0,NY-1
+        TMESH(3,IDIR,IY)=DBLE(GMESH(IY))
+        DDM(IDIR)=MAX(DDM(IDIR),ABS(TMESH(3,IDIR,IX)))
+        ISPLT(3,IDIR,IY)=ISPLG(IY)
+ 11   CONTINUE 
+      TMESH(3,IDIR,NY)=DBLE(GMESH(NY))
+      DDM(IDIR)=MAX(DDM(IDIR),ABS(TMESH(3,IDIR,NX)))
+*----
+*  MESHZ AND SPLITZ
+*----
+      IDIR=3
+      DDM(IDIR)=0.0D0 
+      IF(NZ .GT. 0) THEN
+        GMESH(0:NXYZ)=0.0
+        ISPLG(0:NXYZ-1)=1
+        CALL LCMLEN(IPGEOM,'MESHZ',ILCMLN,ILCMTY)
+        IF(ILCMLN .EQ. NZ+1) 
+     >    CALL LCMGET(IPGEOM,'MESHZ',GMESH)
+        CALL LCMLEN(IPGEOM,'SPLITZ',ILCMLN,ILCMTY)
+        IF(ILCMLN .EQ. NZ) 
+     >    CALL LCMGET(IPGEOM,'SPLITZ',ISPLG)
+        DO 12 IZ=0,NZ-1
+          TMESH(3,IDIR,IZ)=DBLE(GMESH(IZ))
+          DDM(IDIR)=MAX(DDM(IDIR),ABS(TMESH(3,IDIR,IX)))
+          ISPLT(3,IDIR,IZ)=ISPLG(IZ)
+ 12     CONTINUE 
+        TMESH(3,IDIR,NZ)=DBLE(GMESH(NZ))
+        DDM(IDIR)=MAX(DDM(IDIR),ABS(TMESH(3,IDIR,NX)))
+      ELSE
+        NZ=1
+      ENDIF
+C-- MIX
+      ALLOCATE(IMIX(NK))    
+      IMIX(:NK)=0
+      CALL LCMLEN(IPGEOM,'MIX',ILCMLN,ILCMTY)
+      IF(ITYPG .EQ. 21) THEN
+C----
+C  FOR CARCELX REORDER MIXTURE
+C----
+        ALLOCATE(IMIX1(NK))    
+        IF(ILCMLN .GT. 0 .AND. ILCMLN .LE. NK)
+     >    CALL LCMGET(IPGEOM,'MIX',IMIX1)
+        IMTN=0
+        DO 20 IZ=1,NZ
+          DO 21 IY=1,NY
+            DO 22 IX=1,NX
+              DO 23 IR=1,NR+1
+                IMTN=IMTN+1 
+                IMTO=(IX-1)*NZ*NY*(NR+1)+(IZ-1)*NY*(NR+1)+
+     >               (IY-1)*(NR+1)+IR
+                IMIX(IMTN)=IMIX1(IMTO)
+ 23           CONTINUE
+ 22         CONTINUE
+ 21       CONTINUE
+ 20     CONTINUE
+        DEALLOCATE(IMIX1)    
+      ELSE IF(ITYPG .EQ. 22) THEN
+C----
+C  FOR CARCELY REORDER MIXTURE
+C----
+        ALLOCATE(IMIX1(NK))    
+        IF(ILCMLN .GT. 0 .AND. ILCMLN .LE. NK)
+     >    CALL LCMGET(IPGEOM,'MIX',IMIX1)
+        IMTN=0
+        DO 30 IZ=1,NZ
+          DO 31 IY=1,NY
+            DO 32 IX=1,NX
+              DO 33 IR=1,NR+1
+                IMTN=IMTN+1 
+                IMTO=(IY-1)*NZ*NX*(NR+1)+(IX-1)*NZ*(NR+1)+
+     >               (IZ-1)*(NR+1)+IR
+                IMIX(IMTN)=IMIX1(IMTO)
+ 33           CONTINUE
+ 32         CONTINUE
+ 31       CONTINUE
+ 30     CONTINUE
+        DEALLOCATE(IMIX1)    
+      ELSE
+        IF(ILCMLN .GT. 0 .AND. ILCMLN .LE. NK)
+     >    CALL LCMGET(IPGEOM,'MIX',IMIX) 
+      ENDIF
+C- OFF CENTER 
+      OFFCEN(:3)=0.0
+      CALL LCMLEN(IPGEOM,'OFFCENTER',ILCMLN,ILCMTY)
+      IF(ILCMLN .EQ. 3)
+     >  CALL LCMGET(IPGEOM,'OFFCENTER',OFFCEN)
+C----
+C  PRINT INITIAL MESH
+C----
+      IF(IPRT .GE. 20) THEN 
+        WRITE(IOUT,6001) NX,NY,NZ,NR,NK
+        WRITE(IOUT,6002) 'MESHX ='
+        WRITE(IOUT,6003) (TMESH(3,1,IX),IX=0,NX)
+        WRITE(IOUT,6002) 'SPLTX ='
+        WRITE(IOUT,6004) (ISPLT(3,1,IX),IX=0,NX-1)
+        WRITE(IOUT,6002) 'MESHY ='
+        WRITE(IOUT,6003) (TMESH(3,2,IY),IY=0,NY)
+        WRITE(IOUT,6002) 'SPLTY ='
+        WRITE(IOUT,6004) (ISPLT(3,2,IY),IY=0,NY-1)
+        WRITE(IOUT,6002) 'MESHZ ='
+        WRITE(IOUT,6003) (TMESH(3,3,IZ),IZ=0,NZ)
+        WRITE(IOUT,6002) 'SPLTZ ='
+        WRITE(IOUT,6004) (ISPLT(3,3,IZ),IZ=0,NZ-1)
+        WRITE(IOUT,6002) 'MIXT  ='
+        WRITE(IOUT,6004) (IMIX(IX),IX=1,NK)
+        WRITE(IOUT,6002) 'OFFC  ='
+        WRITE(IOUT,6003) (OFFCEN(IX),IX=1,3)
+      ENDIF
+C----
+C  TURN GEOMETRY WITH IGT
+C---- 
+      DO 1000 IG=1,2
+        IF(IGT(IG) .GT. 12 ) THEN
+          IKT(IG)=12-IGT(IG)
+        ELSE
+          IKT(IG)=IGT(IG)
+        ENDIF
+        IF(IG.EQ.1) THEN
+          ALLOCATE(IM1(2*NK))
+          IMIX2=>IM1
+        ELSE
+          ALLOCATE(IM2(2*NK))
+          IMIX2=>IM2
+        ENDIF
+        IMIX2(:NK)=0
+        IF(IKT(IG) .LT. 0) THEN 
+          OFFTR(3,IG)=-OFFCEN(3)
+        ELSE
+          OFFTR(3,IG)=OFFCEN(3)
+        ENDIF
+        IF(ABS(IKT(IG)) .EQ. 1) THEN
+          NTM(IG,1)=NX
+          NTM(IG,2)=NY
+          NTM(IG,3)=NZ
+          DO 100 IX=0,NX-1 
+            TMESH(IG,1,IX)=TMESH(3,1,IX+1)-TMESH(3,1,IX)
+            ISPLT(IG,1,IX)=ISPLT(3,1,IX)
+ 100      CONTINUE
+          DO 110 IY=0,NY-1
+            TMESH(IG,2,IY)=TMESH(3,2,IY+1)-TMESH(3,2,IY)
+            ISPLT(IG,2,IY)=ISPLT(3,2,IY)
+ 110      CONTINUE
+          OFFTR(1,IG)=OFFCEN(1)
+          OFFTR(2,IG)=OFFCEN(2)
+          IF(IKT(IG) .LT. 0) THEN
+            DO 120 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,NZ-IZ)-TMESH(3,3,NZ-IZ-1)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,NZ-IZ-1)
+              ITMI=IZ*NX*NY*(NR+1)
+              IDMI=(NZ-IZ-1)*NX*NY*(NR+1)
+              DO 121 IY=0,NY-1
+                DO 122 IX=0,NX-1
+                  DO 123 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 123              CONTINUE
+ 122            CONTINUE
+ 121          CONTINUE
+ 120        CONTINUE
+          ELSE
+            DO 130 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,IZ+1)-TMESH(3,3,IZ)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,IZ)
+              ITMI=IZ*NX*NY*(NR+1)
+              IDMI=ITMI
+              DO 131 IY=0,NY-1
+                DO 132 IX=0,NX-1
+                  DO 133 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 133              CONTINUE
+ 132            CONTINUE
+ 131          CONTINUE
+ 130        CONTINUE
+          ENDIF
+        ELSE IF(ABS(IKT(IG)) .EQ. 2) THEN 
+C----
+C  ROTATION OF PI/2 
+C----
+          NTM(IG,1)=NY
+          NTM(IG,2)=NX
+          NTM(IG,3)=NZ
+          DDD=DDM(2)
+          DDM(2)=DDM(1)
+          DDM(1)=DDD
+          DO 200 IX=0,NY-1 
+            TMESH(IG,1,IX)=TMESH(3,2,IX+1)-TMESH(3,2,IX)
+            ISPLT(IG,1,IX)=ISPLT(3,2,IX)
+ 200      CONTINUE
+          DO 210 IY=0,NX-1
+            TMESH(IG,2,IY)=TMESH(3,1,NX-IY)-TMESH(3,1,NX-IY-1)
+            ISPLT(IG,2,IY)=ISPLT(3,1,NX-IY-1)
+ 210      CONTINUE
+          OFFTR(1,IG)=OFFCEN(2)
+          OFFTR(2,IG)=-OFFCEN(1)
+          IF(IKT(IG) .LT. 0) THEN
+            DO 220 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,NZ-IZ)-TMESH(3,3,NZ-IZ-1)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,NZ-IZ-1)
+              DO 221 IY=0,NX-1
+                DO 222 IX=0,NY-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NY*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=(NZ-IZ-1)*NX*NY*(NR+1)+(NY-IX-1)*NX*(NR+1)+
+     >                 IY*(NR+1)
+                  DO 223 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 223              CONTINUE
+ 222            CONTINUE
+ 221          CONTINUE
+ 220        CONTINUE
+          ELSE
+            DO 230 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,IZ+1)-TMESH(3,3,IZ)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,IZ)
+              DO 231 IY=0,NX-1
+                DO 232 IX=0,NY-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NY*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=IZ*NX*NY*(NR+1)+IX*NX*(NR+1)+
+     >                 (NX-IY-1)*(NR+1)
+                  DO 233 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 233              CONTINUE
+ 232            CONTINUE
+ 231          CONTINUE
+ 230        CONTINUE 
+          ENDIF
+        ELSE IF(ABS(IKT(IG)) .EQ. 3) THEN
+C----
+C  ROTATION OF PI 
+C----
+          NTM(IG,1)=NX
+          NTM(IG,2)=NY
+          NTM(IG,3)=NZ
+          DO 300 IX=0,NX-1 
+            TMESH(IG,1,IX)=TMESH(3,1,NX-IX)-TMESH(3,1,NX-IX-1)
+            ISPLT(IG,1,IX)=ISPLT(3,1,NX-IX-1)
+ 300      CONTINUE
+          DO 310 IY=0,NY-1
+            TMESH(IG,2,IY)=TMESH(3,2,NY-IY)-TMESH(3,2,NY-IY-1)
+            ISPLT(IG,2,IY)=ISPLT(3,2,NY-IY-1)
+ 310      CONTINUE
+          OFFTR(1,IG)=-OFFCEN(1)
+          OFFTR(2,IG)=-OFFCEN(2)
+          IF(IKT(IG) .LT. 0) THEN
+            DO 320 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,NZ-IZ)-TMESH(3,3,NZ-IZ-1)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,NZ-IZ-1)
+              DO 321 IY=0,NY-1
+                DO 322 IX=0,NX-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NX*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=(NZ-IZ-1)*NX*NY*(NR+1)+(NY-IY-1)*NX*(NR+1)+
+     >                 (NX-IX-1)*(NR+1)
+                  DO 323 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 323              CONTINUE
+ 322            CONTINUE
+ 321          CONTINUE
+ 320        CONTINUE
+          ELSE
+            DO 330 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,IZ+1)-TMESH(3,3,IZ)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,IZ)
+              DO 331 IY=0,NY-1
+                DO 332 IX=0,NX-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NX*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=IZ*NX*NY*(NR+1)+(NY-IY-1)*NX*(NR+1)+
+     >                 (NX-IX-1)*(NR+1)
+                  DO 333 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 333              CONTINUE
+ 332            CONTINUE
+ 331          CONTINUE
+ 330        CONTINUE
+          ENDIF
+        ELSE IF(ABS(IKT(IG)) .EQ. 4) THEN
+C----
+C  ROTATION OF 3*PI/2 
+C----
+          NTM(IG,1)=NY
+          NTM(IG,2)=NX
+          NTM(IG,3)=NZ
+          DDD=DDM(2)
+          DDM(2)=DDM(1)
+          DDM(1)=DDD
+          DO 400 IX=0,NY-1 
+            TMESH(IG,1,IX)=TMESH(3,2,NY-IX)-TMESH(3,2,NY-IX-1)
+            ISPLT(IG,1,IX)=ISPLT(3,2,NY-IX-1)
+ 400      CONTINUE
+          DO 410 IY=0,NX-1
+            TMESH(IG,2,IY)=TMESH(3,1,IY+1)-TMESH(3,1,IY)
+            ISPLT(IG,2,IY)=ISPLT(3,1,IY)
+ 410      CONTINUE
+          OFFTR(1,IG)=-OFFCEN(2)
+          OFFTR(2,IG)=OFFCEN(1)
+          IF(IKT(IG) .LT. 0) THEN
+            DO 420 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,NZ-IZ)-TMESH(3,3,NZ-IZ-1)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,NZ-IZ-1)
+              DO 421 IY=0,NX-1
+                DO 422 IX=0,NY-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NY*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=(NZ-IZ-1)*NX*NY*(NR+1)+(NY-IX-1)*NX*(NR+1)+
+     >                 IY*(NR+1)
+                  DO 423 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 423              CONTINUE
+ 422            CONTINUE
+ 421          CONTINUE
+ 420        CONTINUE
+          ELSE
+            DO 430 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,IZ+1)-TMESH(3,3,IZ)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,IZ)
+              DO 431 IY=0,NX-1
+                DO 432 IX=0,NY-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NY*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=IZ*NX*NY*(NR+1)+(NY-IX-1)*NX*(NR+1)+
+     >                 IY*(NR+1)
+                  DO 433 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 433              CONTINUE
+ 432            CONTINUE
+ 431          CONTINUE
+ 430        CONTINUE
+          ENDIF
+        ELSE IF(ABS(IKT(IG)) .EQ. 5) THEN
+C----
+C  REFLECTION WITH RESPECT TO AXIS  // TO Y
+C----
+          NTM(IG,1)=NX
+          NTM(IG,2)=NY
+          NTM(IG,3)=NZ
+          DO 500 IX=0,NX-1 
+            TMESH(IG,1,IX)=TMESH(3,1,NX-IX)-TMESH(3,1,NX-IX-1)
+            ISPLT(IG,1,IX)=ISPLT(3,1,NX-IX-1)
+ 500      CONTINUE
+          DO 510 IY=0,NY-1
+            TMESH(IG,2,IY)=TMESH(3,2,IY+1)-TMESH(3,2,IY)
+            ISPLT(IG,2,IY)=ISPLT(3,2,IY)
+ 510      CONTINUE
+          OFFTR(1,IG)=-OFFCEN(1)
+          OFFTR(2,IG)=OFFCEN(2)
+          IF(IKT(IG) .LT. 0) THEN
+            DO 520 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,NZ-IZ)-TMESH(3,3,NZ-IZ-1)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,NZ-IZ-1)
+              DO 521 IY=0,NY-1
+                DO 522 IX=0,NX-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NX*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=(NZ-IZ-1)*NX*NY*(NR+1)+IY*NX*(NR+1)+
+     >                 (NX-IX-1)*(NR+1)
+                  DO 523 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 523              CONTINUE
+ 522            CONTINUE
+ 521          CONTINUE
+ 520        CONTINUE
+          ELSE
+            DO 530 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,IZ+1)-TMESH(3,3,IZ)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,IZ)
+              DO 531 IY=0,NY-1
+                DO 532 IX=0,NX-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NX*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=IZ*NX*NY*(NR+1)+IY*NX*(NR+1)+
+     >                 (NX-IX-1)*(NR+1)
+                  DO 533 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 533              CONTINUE
+ 532            CONTINUE
+ 531          CONTINUE
+ 530        CONTINUE
+          ENDIF
+        ELSE IF(ABS(IKT(IG)) .EQ. 6) THEN
+C----
+C  ROTATION OF PI/2 FOLLOWED BY 
+C  REFLECTION WITH RESPECT TO AXIS  // TO Y
+C----
+          NTM(IG,1)=NY
+          NTM(IG,2)=NX
+          NTM(IG,3)=NZ
+          DDD=DDM(2)
+          DDM(2)=DDM(1)
+          DDM(1)=DDD
+          DO 600 IX=0,NY-1 
+            TMESH(IG,1,IX)=TMESH(3,2,IX+1)-TMESH(3,2,IX)
+            ISPLT(IG,1,IX)=ISPLT(3,2,IX)
+ 600      CONTINUE
+          DO 610 IY=0,NX-1
+            TMESH(IG,2,IY)=TMESH(3,1,IY+1)-TMESH(3,1,IY)
+            ISPLT(IG,2,IY)=ISPLT(3,1,IY)
+ 610      CONTINUE
+          OFFTR(1,IG)=OFFCEN(2)
+          OFFTR(2,IG)=OFFCEN(1)
+          IF(IKT(IG) .LT. 0) THEN
+            DO 620 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,NZ-IZ)-TMESH(3,3,NZ-IZ-1)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,NZ-IZ-1)
+              DO 621 IY=0,NX-1
+                DO 622 IX=0,NY-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NY*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=(NZ-IZ-1)*NX*NY*(NR+1)+IX*NX*(NR+1)+
+     >                 IY*(NR+1)
+                  DO 623 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 623              CONTINUE
+ 622            CONTINUE
+ 621          CONTINUE
+ 620        CONTINUE
+          ELSE
+            DO 630 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,IZ+1)-TMESH(3,3,IZ)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,IZ)
+              DO 631 IY=0,NX-1
+                DO 632 IX=0,NY-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NY*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=IZ*NX*NY*(NR+1)+IX*NX*(NR+1)+
+     >                 IY*(NR+1)
+                  DO 633 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 633              CONTINUE
+ 632            CONTINUE
+ 631          CONTINUE
+ 630        CONTINUE 
+          ENDIF
+        ELSE IF(ABS(IKT(IG)) .EQ. 7) THEN
+C----
+C  REFLECTION WITH RESPECT TO AXIS // TO X
+C----
+          NTM(IG,1)=NX
+          NTM(IG,2)=NY
+          NTM(IG,3)=NZ
+          DO 700 IX=0,NX-1 
+            TMESH(IG,1,IX)=TMESH(3,1,IX+1)-TMESH(3,1,IX)
+            ISPLT(IG,1,IX)=ISPLT(3,1,IX)
+ 700      CONTINUE
+          DO 710 IY=0,NY-1
+            TMESH(IG,2,IY)=TMESH(3,2,NY-IY)-TMESH(3,2,NY-IY-1)
+            ISPLT(IG,2,IY)=ISPLT(3,2,NY-IY-1)
+ 710      CONTINUE
+          OFFTR(1,IG)=OFFCEN(1)
+          OFFTR(2,IG)=-OFFCEN(2)
+          IF(IKT(IG) .LT. 0) THEN
+            DO 720 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,NZ-IZ)-TMESH(3,3,NZ-IZ-1)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,NZ-IZ-1)
+              DO 721 IY=0,NY-1
+                ITMI=IZ*NX*NY*(NR+1)+IY*NX*(NR+1)
+                IDMI=(NZ-IZ-1)*NX*NY*(NR+1)+(NY-IY-1)*NX*(NR+1)
+                DO 722 IX=0,NX-1
+                  DO 723 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 723              CONTINUE
+ 722            CONTINUE
+ 721          CONTINUE
+ 720        CONTINUE
+          ELSE
+            DO 730 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,IZ+1)-TMESH(3,3,IZ)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,IZ)
+              DO 731 IY=0,NY-1
+                ITMI=IZ*NX*NY*(NR+1)+IY*NX*(NR+1)
+                IDMI=IZ*NX*NY*(NR+1)+(NY-IY-1)*NX*(NR+1)
+                DO 732 IX=0,NX-1
+                  DO 733 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 733              CONTINUE
+ 732            CONTINUE
+ 731          CONTINUE
+ 730        CONTINUE
+          ENDIF
+        ELSE IF(ABS(IKT(IG)) .EQ. 8) THEN
+C----
+C  ROTATION OF PI/2 FOLLOWED BY
+C  REFLECTION WITH RESPECT TO AXIS // TO X
+C----
+          NTM(IG,1)=NY
+          NTM(IG,2)=NX
+          NTM(IG,3)=NZ
+          DDD=DDM(2)
+          DDM(2)=DDM(1)
+          DDM(1)=DDD
+          DO 800 IX=0,NY-1 
+            TMESH(IG,1,IX)=TMESH(3,2,NY-IX)-TMESH(3,2,NY-IX-1)
+            ISPLT(IG,1,IX)=ISPLT(3,2,NY-IX-1)
+ 800      CONTINUE
+          DO 810 IY=0,NX-1
+            TMESH(IG,2,IY)=TMESH(3,1,NX-IY)-TMESH(3,1,NX-IY-1)
+            ISPLT(IG,2,IY)=ISPLT(3,1,NX-IY-1)
+ 810      CONTINUE
+          OFFTR(1,IG)=-OFFCEN(2)
+          OFFTR(2,IG)=-OFFCEN(1)
+          IF(IKT(IG) .LT. 0) THEN
+            DO 820 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,NZ-IZ)-TMESH(3,3,NZ-IZ-1)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,NZ-IZ-1)
+              DO 821 IY=0,NX-1
+                DO 822 IX=0,NY-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NY*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=(NZ-IZ-1)*NX*NY*(NR+1)+(NY-IX-1)*NX*(NR+1)+
+     >                 (NX-IY-1)*(NR+1)
+                  DO 823 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 823              CONTINUE
+ 822            CONTINUE
+ 821          CONTINUE
+ 820        CONTINUE
+          ELSE
+            DO 830 IZ=0,NZ-1
+              TMESH(IG,3,IZ)=TMESH(3,3,IZ+1)-TMESH(3,3,IZ)
+              ISPLT(IG,3,IZ)=ISPLT(3,3,IZ)
+              DO 831 IY=0,NX-1
+                DO 832 IX=0,NY-1
+                  ITMI=IZ*NX*NY*(NR+1)+IY*NY*(NR+1)+
+     >                 IX*(NR+1)
+                  IDMI=IZ*NX*NY*(NR+1)+(NY-IX-1)*NX*(NR+1)+
+     >                 (NX-IY-1)*(NR+1)
+                  DO 833 IR=0,NR
+                    IMIX2(ITMI+1)=IMIX(IDMI+1)
+                    ITMI=ITMI+1
+                    IDMI=IDMI+1
+ 833              CONTINUE
+ 832            CONTINUE
+ 831          CONTINUE
+ 830        CONTINUE 
+          ENDIF
+        ENDIF
+C----
+C  PRINT TURNED MESH
+C----
+        IF(IPRT .GE. 20) THEN 
+          WRITE(IOUT,6010) IKT(IG)
+          WRITE(IOUT,6002) 'MESHX ='
+          WRITE(IOUT,6003) (TMESH(IG,1,IX),IX=0,NX)
+          WRITE(IOUT,6002) 'SPLTX ='
+          WRITE(IOUT,6004) (ISPLT(IG,1,IX),IX=0,NX-1)
+          WRITE(IOUT,6002) 'MESHY ='
+          WRITE(IOUT,6003) (TMESH(IG,2,IY),IY=0,NY)
+          WRITE(IOUT,6002) 'SPLTY ='
+          WRITE(IOUT,6004) (ISPLT(IG,2,IY),IY=0,NY-1)
+          WRITE(IOUT,6002) 'MESHZ ='
+          WRITE(IOUT,6003) (TMESH(IG,3,IZ),IZ=0,NZ)
+          WRITE(IOUT,6002) 'SPLTZ ='
+          WRITE(IOUT,6004) (ISPLT(IG,3,IZ),IZ=0,NZ-1)
+          WRITE(IOUT,6002) 'MIXT  ='
+          WRITE(IOUT,6004) (IMIX2(IX),IX=1,NK)
+          WRITE(IOUT,6002) 'OFFC  ='
+          WRITE(IOUT,6003) (OFFTR(IX,IG),IX=1,3) 
+        ENDIF
+ 1000 CONTINUE
+C----
+C  COMPARE GEOMETRY
+C  1- MESH AND SPLIT IN X, Y AND Z
+C  2- MIXTURES
+C  3- OFFCENTER
+C----
+      DO 900 IDIR=1,3 
+        IF(NTM(1,IDIR) .EQ. NTM(2,IDIR)) THEN
+          DO 910 IX=0,NTM(1,IDIR)-1
+            DDD=ABS(TMESH(2,IDIR,IX)-TMESH(1,IDIR,IX))
+            IF(DDD .GT. 1.0D-6*ABS(DDM(IDIR))       .OR.
+     >         ISPLT(2,IDIR,IX) .NE. ISPLT(1,IDIR,IX) ) THEN
+               WRITE(IOUT,6020) IDIR,IX,
+     >           ISPLT(1,IDIR,IX),ISPLT(2,IDIR,IX),
+     >           TMESH(1,IDIR,IX),TMESH(2,IDIR,IX),DDD 
+              LELCSY=.FALSE. 
+              GO TO 995
+            ENDIF
+ 910      CONTINUE
+        ELSE
+          LELCSY=.FALSE.
+          GO TO 995
+        ENDIF
+ 900  CONTINUE
+      DO 920 IX=1,NK
+        IF(IM1(IX) .NE. IM2(IX) ) THEN
+          LELCSY=.FALSE.
+          WRITE(IOUT,6021) IX,IM1(IX),IM2(IX)
+          GO TO 995
+        ENDIF
+ 920  CONTINUE
+      IF(OFFTR(1,1) .NE. OFFTR(1,2) .OR.
+     >   OFFTR(2,1) .NE. OFFTR(2,2) .OR.   
+     >   OFFTR(3,1) .NE. OFFTR(3,2) ) THEN
+        LELCSY=.FALSE.
+        GO TO 995
+      ENDIF
+ 995  CONTINUE
+C----
+C  RELEASE MEMORY
+C----
+      DEALLOCATE(IM2,IM1,IMIX)
+      IF(GEOCV .NE. ' ') THEN  
+        CALL LCMSIX(IPGEOM,GEOCV,2) 
+      ENDIF
+C----
+C  RETURN 
+C----
+      RETURN
+C----
+C  FORMATS
+C----
+ 6000 FORMAT(1X,A6,'-- ANALYZING :',A12,2I10)
+ 6001 FORMAT(1X,' DIMENSIONS =',5I10/1X,' ORIGINAL MESH ')
+ 6002 FORMAT(1X,A7)
+ 6003 FORMAT(5F15.9)
+ 6004 FORMAT(5I15) 
+ 6010 FORMAT(1X,' GEOMETRY AFTER TURN = ',I10)
+ 6020 FORMAT(1X,'ERROR FOR DIRECTION ',2I10/
+     >       1X,'SPLIT = ',2I10/
+     >       1X,'MESH  = ',3F15.9)
+ 6021 FORMAT(1X,'ERROR FOR MIXTURE ',I10,2I10)
+      END