Initial commit from Polytechnique Montreal

1 files changed, 395 insertions, 0 deletions

diff --git a/Dragon/src/SNT.f b/Dragon/src/SNT.f
new file mode 100644
index 0000000..fecbc83
--- /dev/null
+++ b/Dragon/src/SNT.f
@@ -0,0 +1,395 @@
+*DECK SNT
+      SUBROUTINE SNT(NENTRY,HENTRY,IENTRY,JENTRY,KENTRY)
+*
+*-----------------------------------------------------------------------
+*
+*Purpose:
+* SN method tracking operator.
+*
+*Copyright:
+* Copyright (C) 2005 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
+*
+*Parameters: input/output
+* NENTRY  number of LCM objects or files used by the operator.
+* HENTRY  name of each LCM object or file:
+*         HENTRY(1) create or modification type(L_TRACK);
+*         HENTRY(2) read-only type(L_GEOM).
+* IENTRY  type of each LCM object or file:
+*         =1 LCM memory object; =2 XSM file; =3 sequential binary file;
+*         =4 sequential ascii file.
+* JENTRY  access of each LCM object or file:
+*         =0 the LCM object or file is created;
+*         =1 the LCM object or file is open for modifications;
+*         =2 the LCM object or file is open in read-only mode.
+* KENTRY  LCM object address or file unit number.
+*
+*-----------------------------------------------------------------------
+*
+      USE GANLIB
+*----
+*  SUBROUTINE ARGUMENTS
+*----
+      INTEGER      NENTRY,IENTRY(NENTRY),JENTRY(NENTRY)
+      TYPE(C_PTR)  KENTRY(NENTRY)
+      CHARACTER    HENTRY(NENTRY)*12
+*----
+*  LOCAL VARIABLES
+*----
+      PARAMETER(NSTATE=40,IOUT=6)
+      TYPE(C_PTR) IPTRK,IPGEOM
+      CHARACTER TEXT4*4,TEXT12*12,TITLE*72,HSIGN*12
+      DOUBLE PRECISION DFLOTT
+      REAL EPSI
+      LOGICAL LOG,LFIXUP,LDSA,LBIHET,LIVO,LSHOOT
+      INTEGER IGP(NSTATE),ISTATE(NSTATE),NCODE(6),NITMA,EELEM,ESCHM,IGLK
+*----
+*  PARAMETER VALIDATION
+*----
+      IF(NENTRY.LE.1) CALL XABORT('SNT: TWO PARAMETERS EXPECTED.')
+      IF((IENTRY(1).NE.1).AND.(IENTRY(1).NE.2)) CALL XABORT('SNT: L'
+     1 //'CM OBJECT EXPECTED AT LHS.')
+      IF((JENTRY(1).NE.0).AND.(JENTRY(1).NE.1)) CALL XABORT('SNT: E'
+     1 //'NTRY IN CREATE OR MODIFICATION MODE EXPECTED.')
+      IF((JENTRY(2).NE.2).OR.((IENTRY(2).NE.1).AND.(IENTRY(2).NE.2)))
+     1 CALL XABORT('SNT: LCM OBJECT IN READ-ONLY MODE EXPECTED AT R'
+     2 //'HS.')
+      CALL LCMGTC(KENTRY(2),'SIGNATURE',12,HSIGN)
+      IF(HSIGN.NE.'L_GEOM') THEN
+         TEXT12=HENTRY(2)
+         CALL XABORT('SNT: SIGNATURE OF '//TEXT12//' IS '//HSIGN//
+     1   '. L_GEOM EXPECTED.')
+      ENDIF
+      IPTRK=KENTRY(1)
+      IPGEOM=KENTRY(2)
+      CALL LCMGET(IPGEOM,'STATE-VECTOR',ISTATE)
+      ITYPE=ISTATE(1)
+*
+      LBIHET=.FALSE.
+      IMPX=1
+      TITLE=' '
+      IF(JENTRY(1).EQ.0) THEN
+         HSIGN='L_TRACK'
+         CALL LCMPTC(IPTRK,'SIGNATURE',12,HSIGN)
+         HSIGN='SN'
+         CALL LCMPTC(IPTRK,'TRACK-TYPE',12,HSIGN)
+         MAXPTS=ISTATE(6)
+         ISCHM=1
+         ESCHM=1
+         IELEM=1
+         EELEM=1
+         NLF=0
+         ISCAT=0
+         IQUAD=2
+         LFIXUP=.FALSE.
+         LDSA=.FALSE.
+         LIVO=.TRUE.
+         IGLK=0
+         ICL1=3
+         ICL2=3
+         ISPLH=0
+         IF((ITYPE.EQ.8).OR.(ITYPE.EQ.9)) THEN
+            CALL LCMLEN(IPGEOM,'SPLITL',ILEN,ITYLCM)
+            IF(ILEN.GT.0)THEN
+               CALL LCMGET(IPGEOM,'SPLITL',ISPLH)
+               IF(ISPLH.EQ.0) ISPLH=1
+               MAXPTS=MAXPTS*3*ISPLH**2
+            ELSE
+               CALL XABORT('SNT: SPLITL SPECIFIER NEEDED FOR SN '//
+     1          'WITH HEXAGONAL GEOMETRY.')
+            ENDIF
+         ENDIF
+         NSTART=0
+         NSDSA=1000
+         IELEMSA=1
+         ISOLVSA=2
+         MAXIT=100
+         EPSI=1.0E-5
+         CALL LCMGET(IPGEOM,'NCODE',NCODE)
+         LOG=.FALSE.
+         DO 10 I=1,4
+         LOG=LOG.OR.(NCODE(I).EQ.3)
+   10    CONTINUE
+         IF(LOG) MAXPTS=2*MAXPTS
+         IQUA10=0
+         IBIHET=2
+         INSB=0
+         MCELL=0
+         NMPI=0
+         NFOU=0
+         LSHOOT=.TRUE.
+         CALL LCMLEN(IPGEOM,'BIHET',ILONG,ITYLCM)
+         LBIHET=(ILONG.NE.0)
+         IF(LBIHET) IQUA10=5
+         IBFP=0
+      ELSE IF(JENTRY(1).EQ.1) THEN
+         CALL LCMGTC(IPTRK,'SIGNATURE',12,HSIGN)
+         IF(HSIGN.NE.'L_TRACK') THEN
+            TEXT12=HENTRY(1)
+            CALL XABORT('SNT: SIGNATURE OF '//TEXT12//' IS '//HSIGN
+     1      //'. L_TRACK EXPECTED.')
+         ENDIF
+         CALL LCMGTC(IPTRK,'TRACK-TYPE',12,HSIGN)
+         IF(HSIGN.NE.'SN') THEN
+            TEXT12=HENTRY(1)
+            CALL XABORT('SNT: TRACK-TYPE OF '//TEXT12//' IS '//HSIGN
+     1      //'. SN EXPECTED.')
+         ENDIF
+         CALL LCMGET(IPTRK,'STATE-VECTOR',IGP)
+         MAXPTS=IGP(1)
+         IELEM=IGP(8)
+         ISCHM=IGP(10)
+         NLF=IGP(15)
+         ISCAT=IGP(16)
+         IQUAD=IGP(17)
+         LFIXUP=IGP(18).EQ.1
+         LDSA=IGP(19).EQ.1
+         NSTART=IGP(20)
+         NSDSA=IGP(21)
+         MAXIT=IGP(22)
+         LIVO=IGP(23).EQ.1
+         ICL1=IGP(24)
+         ICL2=IGP(25)
+         ISPLH=IGP(26)
+         INSB=IGP(27)
+         MCELL=IGP(28)
+         LSHOOT=IGP(30).EQ.1
+         NMPI=IGP(32)
+         ISOLVSA=IGP(33)
+         NFOU=IGP(34)
+         EELEM=IGP(35)
+         ESCHM=IGP(36)
+         IGLK=IGP(37)
+         LBIHET=IGP(40).GT.0
+         EPSI=1.0E-5
+         IF(LBIHET) THEN
+            CALL LCMSIX(IPTRK,'BIHET',1)
+            CALL LCMGET(IPTRK,'PARAM',IGP)
+            CALL LCMSIX(IPTRK,'BIHET',2)
+            IBIHET=IGP(6)
+            IQUA10=IGP(8)
+         ELSE
+            IBIHET=0
+            IQUA10=0
+         ENDIF
+         CALL LCMLEN(IPTRK,'TITLE',LENGT,ITYLCM)
+         IF(LENGT.GT.0) CALL LCMGTC(IPTRK,'TITLE',72,TITLE)
+      ENDIF
+      FRTM=0.05
+   15 CALL REDGET(INDIC,NITMA,FLOTT,TEXT4,DFLOTT)
+   16 CONTINUE
+      IF(INDIC.EQ.10) GO TO 30
+      IF(INDIC.NE.3) CALL XABORT('SNT: CHARACTER DATA EXPECTED.')
+      IF(TEXT4.EQ.'EDIT') THEN
+         CALL REDGET(INDIC,IMPX,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(1).')
+      ELSE IF(TEXT4.EQ.'TITL') THEN
+         CALL REDGET(INDIC,NITMA,FLOTT,TITLE,DFLOTT)
+         IF(INDIC.NE.3) CALL XABORT('SNT: TITLE EXPECTED.')
+      ELSE IF(TEXT4.EQ.'MAXR') THEN
+         CALL REDGET(INDIC,MAXPTS,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(2).')
+      ELSE IF(TEXT4.EQ.'SCHM') THEN
+         CALL REDGET(INDIC,ISCHM,FLOTT,TEXT4,DFLOTT)
+         IF(ISCHM.EQ.2) IELEM=2
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(2.1).')
+      ELSE IF(TEXT4.EQ.'ESCH') THEN
+         CALL REDGET(INDIC,ESCHM,FLOTT,TEXT4,DFLOTT)
+         IF(ESCHM.EQ.2) EELEM=2
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(2.2).')
+      ELSE IF(TEXT4.EQ.'DIAM') THEN
+         CALL REDGET(INDIC,IELEM,FLOTT,TEXT4,DFLOTT)
+         IELEM=IELEM+1
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(3).')
+      ELSE IF(TEXT4.EQ.'EDIA') THEN
+         CALL REDGET(INDIC,EELEM,FLOTT,TEXT4,DFLOTT)
+         EELEM=EELEM+1
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(3.1).')
+      ELSE IF(TEXT4.EQ.'SN') THEN
+         CALL REDGET(INDIC,NLF,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(4).')
+         IF(MOD(NLF,2).EQ.1) CALL XABORT('SNT: EVEN SN ORDER EXPECTED.')
+         ISCAT=NLF
+      ELSE IF(TEXT4.EQ.'SCAT') THEN
+         IF(NLF.EQ.0) CALL XABORT('SNT: DEFINE SN FIRST.')
+         CALL REDGET(INDIC,ISCAT,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(5).')
+      ELSE IF(TEXT4.EQ.'MAXI') THEN
+         CALL REDGET(INDIC,MAXIT,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(6).')
+      ELSE IF(TEXT4.EQ.'EPSI') THEN
+         CALL REDGET(INDIC,NITMA,EPSI,TEXT4,DFLOTT)
+         IF(INDIC.NE.2) CALL XABORT('SNT: REAL DATA EXPECTED(1).')
+      ELSE IF(TEXT4.EQ.'QUAD') THEN
+         CALL REDGET(INDIC,IQUAD,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(7).')
+         IF((IQUAD.LT.1).OR.(IQUAD.GT.20)) CALL XABORT('SNT: INVALID'//
+     1   ' QUAD VALUE.')
+      ELSE IF(TEXT4.EQ.'NFIX') THEN
+         LFIXUP=.TRUE.
+         CALL REDGET(INDIC,ICL1,FLOTT,TEXT4,DFLOTT)
+      ELSE IF(TEXT4.EQ.'LIVO') THEN
+         LIVO=.TRUE.
+         CALL REDGET(INDIC,ICL1,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(8).')
+         CALL REDGET(INDIC,ICL2,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(9).')
+      ELSE IF(TEXT4.EQ.'NLIV') THEN
+         LIVO=.FALSE.
+      ELSE IF(TEXT4.EQ.'DSA') THEN
+         LDSA=.TRUE.
+         CALL REDGET(INDIC,NSDSA,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(10).')
+         CALL REDGET(INDIC,IELEMSA,FLOTT,TEXT4,DFLOTT)
+         IELEMSA=IELEMSA + 1
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(11).')
+         CALL REDGET(INDIC,ISOLVSA,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(12).')
+      ELSE IF(TEXT4.EQ.'NDSA') THEN
+         LDSA=.FALSE.
+      ELSE IF(TEXT4.EQ.'GMRE') THEN
+         CALL REDGET(INDIC,NSTART,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(13).')
+         IF(NSTART.LT.0) CALL XABORT('SNT: POSITIVE VALUE EXPECTED.')
+      ELSE IF(TEXT4.EQ.'QUAB') THEN
+         CALL REDGET(INDIC,IQUA10,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(14).')
+      ELSE IF(TEXT4.EQ.'SAPO') THEN
+         IBIHET=1
+      ELSE IF(TEXT4.EQ.'HEBE') THEN
+         IBIHET=2
+      ELSE IF(TEXT4.EQ.'SLSI') THEN
+         IBIHET=3
+         CALL REDGET(INDIC,NITMA,FLOTT,TEXT4,DFLOTT)
+         IF (INDIC.NE.2) GOTO 16
+         FRTM=FLOTT
+      ELSE IF(TEXT4.EQ.'SLSS') THEN
+         IBIHET=4
+         CALL REDGET(INDIC,NITMA,FLOTT,TEXT4,DFLOTT)
+         IF (INDIC.NE.2) GOTO 16
+         FRTM=FLOTT
+      ELSE IF(TEXT4.EQ.'ONEG') THEN
+         INSB=0
+      ELSE IF(TEXT4.EQ.'ALLG') THEN
+         INSB=1
+      ELSE IF(TEXT4.EQ.'KBA') THEN
+         ! Enable parallelisation over macrocells in wavefront (and 
+         ! energies depending ONEG/ALLG)
+         CALL REDGET(INDIC,MCELL,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(15).')
+         IF(MCELL.LT.0) CALL XABORT('SNT: POSITIVE INTEGER EXPECTED.')
+      ELSE IF(TEXT4.EQ.'MPIM') THEN
+         ! Compute graph and mpi rank/process allocation depending on 
+         ! the problem for Wyvern
+         CALL REDGET(INDIC,NMPI,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(16).')
+         IF(NMPI.LE.0) CALL XABORT('SNT: POSITIVE INTEGER EXPECTED.')
+      ELSE IF(TEXT4.EQ.'NSHT') THEN
+         LSHOOT=.FALSE.
+      ELSE IF(TEXT4.EQ.'BTE') THEN
+         ! Boltzmann transport equation
+         IBFP=0
+      ELSE IF(TEXT4.EQ.'BFPG') THEN
+         ! Boltzmann Fokker-Planck equation with Galerkin energy
+         ! propagation factors
+         IBFP=1
+      ELSE IF(TEXT4.EQ.'BFPL') THEN
+         ! Boltzmann Fokker-Planck equation with Przybylski and Ligou
+         ! energy propagation factors
+         IBFP=2
+      ELSE IF(TEXT4.EQ.'FOUR') THEN
+         ! Perform Fourier Analysis
+         CALL REDGET(INDIC,NFOU,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(17).')
+      ELSE IF(TEXT4.EQ.'GQ') THEN
+         ! Galerkin quadrature
+         CALL REDGET(INDIC,IGLK,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('SNT: INTEGER DATA EXPECTED(14).')
+      ELSE IF(TEXT4.EQ.';') THEN
+         GO TO 30
+      ELSE
+         CALL XABORT('SNT: '//TEXT4//' IS AN INVALID KEY WORD.')
+      ENDIF
+      GO TO 15
+*
+   30 IF(TITLE.NE.' ') CALL LCMPTC(IPTRK,'TITLE',72,TITLE)
+      TEXT12=HENTRY(2)
+      CALL LCMPTC(IPTRK,'LINK.GEOM',12,TEXT12)
+      IF(IMPX.GT.1) WRITE(IOUT,100) TITLE
+*
+      IF(MAXPTS.EQ.0) CALL XABORT('SNT: MAXPTS NOT DEFINED.')
+      IF(IGLK.EQ.0) ISCAT=MIN(ISCAT,NLF)
+      CALL SNTRK(MAXPTS,IPTRK,IPGEOM,IMPX,ISCHM,IELEM,ISPLH,INSB,
+     1 NLF,MAXIT,EPSI,ISCAT,IQUAD,LFIXUP,LIVO,ICL1,ICL2,LDSA,NSTART,
+     2 NSDSA,IELEMSA,ISOLVSA,LBIHET,LSHOOT,IBFP,MCELL,NMPI,NFOU,
+     3 EELEM,ESCHM,IGLK)
+*
+      IF(IMPX.GT.1) THEN
+         CALL LCMGET(IPTRK,'STATE-VECTOR',IGP)
+         WRITE(IOUT,110) (IGP(I),I=1,15)
+         WRITE(IOUT,120) (IGP(I),I=16,31),IGP(35),IGP(36),IGP(37),
+     1   IGP(40),EPSI
+      ENDIF
+*----
+*  PROCESS DOUBLE HETEROGENEITY (BIHET) DATA (IF AVAILABLE)
+*----
+      IF(LBIHET) CALL XDRTBH(IPGEOM,IPTRK,IQUA10,IBIHET,IMPX,FRTM)
+      RETURN
+*
+  100 FORMAT(1H1,19H SSSSS  NN      NN ,95(1H*)/
+     1 20H SSSSSSS NNN     NN ,57(1H*),
+     2 38H MULTIGROUP VERSION.  A. HEBERT (2005)/
+     3 19H SS   SS NNNN    NN/19H  SS     NN NN   NN/
+     4 19H    SS   NN  NN  NN/19H SS   SS NN   NN NN/
+     5 19H SSSSSSS NN    NNNN/19H  SSSSS  NN     NNN//1X,A72//)
+  110 FORMAT(/14H STATE VECTOR:/
+     1 7H NREG  ,I8,22H   (NUMBER OF REGIONS)/
+     2 7H NUN   ,I8,23H   (NUMBER OF UNKNOWNS)/
+     3 7H ILK   ,I8,39H   (0=LEAKAGE PRESENT/1=LEAKAGE ABSENT)/
+     4 7H NBMIX ,I8,36H   (MAXIMUM NUMBER OF MIXTURES USED)/
+     5 7H NSURF ,I8,29H   (NUMBER OF OUTER SURFACES)/
+     6 7H ITYPE ,I8,21H   (TYPE OF GEOMETRY)/
+     7 7H NFUNL ,I8,45H   (NUMBER OF SPHERICAL HARMONICS COMPONENTS)/
+     8 7H IELEM ,I8,46H   (ORDER OF POLYNOMIAL USED IN SPATIAL APPROX,
+     9 38H./1=CONSTANT/2=LINEAR/3=PARABOLIC/...)/
+     1 7H NDIM  ,I8,35H   (NUMBER OF GEOMETRIC DIMENSIONS)/
+     2 7H ISCHM ,I8,46H   (METHOD OF SPATIAL DISCRETISATION/1=HODD/2=,
+     3 9HDG/3=AWD)/
+     3 7H LL4   ,I8,36H   (ORDER OF THE MATRICES PER GROUP)/
+     4 7H LX    ,I8,38H   (NUMBER OF MESHES ALONG THE X AXIS)/
+     5 7H LY    ,I8,38H   (NUMBER OF MESHES ALONG THE Y AXIS)/
+     6 7H LZ    ,I8,38H   (NUMBER OF MESHES ALONG THE Z AXIS)/
+     7 7H NLF   ,I8,13H   (SN ORDER))
+  120 FORMAT(
+     1 7H ISCAT ,I8,47H   (1=ISOTROPIC SOURCE/2=LINEARLY ANISOTROPIC S,
+     2 6HOURCE)/
+     3 7H IQUAD ,I8,47H   (<4=LEVEL-SYMMETRIC OF TYPE IQUAD/4=LEGENDRE,
+     4 58H-CHEBYSHEV/5=SYMMETRIC LEGENDRE-CHEBYSHEV/6=QR/>9=PRODUCT)/
+     5 7H IFIX  ,I8,29H   (0/1: NEGATIVE FLUX FIXUP)/
+     6 7H IDSA  ,I8,32H   (0/1: SYNTHETIC ACCELERATION)/
+     7 7H NSTART,I8,32H   (NUMBER OF RESTARTS IN GMRES)/
+     8 7H NSDSA ,I8,45H   (NUMBER OF ITERATIONS BEFORE ENABLING DSA)/
+     9 7H MAXIT ,I8,39H   (MAXIMUM NUMBER OF INNER ITERATIONS)/
+     1 7H LIVO  ,I8,38H   (0/1: LIVOLANT ACCELERATION OFF/ON)/
+     2 7H ICL1  ,I8,39H   (NUMBER OF FREE ITERATIONS IN LIVO.)/
+     3 7H ICL2  ,I8,46H   (NUMBER OF ACCELERATED ITERATIONS IN LIVO.)/
+     4 7H ISPLH ,I8,46H   (DEGREE OF LOZENGE SPLITTING FOR HEXAGONAL ,
+     5 9HGEOMETRY)/
+     6 7H INSB  ,I8,36H   (0/1: GROUP VECTORIZATION OFF/ON)/
+     7 7H MCELL ,I8,37H   (0/>0: KBA WAVEFRONT SWEEP OFF/ON)/
+     8 7H IGAV  ,I8,45H   (CONDITION AT AXIAL AXIS FOR CYL./SPH. 1D)/
+     9 7H LSHOOT,I8,38H   (0/1: SHOOTING METHOD IN 1D OFF/ON)/
+     1 7H IBFP  ,I8,43H   (0/1/2: BFP SOLUTION OFF/GALERKIN/LIGOU)/
+     2 7H EELEM ,I8,45H   (ORDER OF POLYNOMIAL USED IN ENERGY APPROX,
+     3 38H./1=CONSTANT/2=LINEAR/3=PARABOLIC/...)/
+     4 7H ESCHM ,I8,45H   (METHOD OF ENERGY DISCRETISATION/1=HODD/2=,
+     5 9HDG/3=AWD)/
+     6 7H IGLK  ,I8,42H   (0=CLASSICAL SN/>0=GALERKIN QUADRATURE)/
+     7 7H IBIHET,I8,47H   (0/1: DOUBLE HETEROGENEITY IS NOT/IS ACTIVE)/
+     8 7H EPSI  ,E11.1,45H  (CONVERGENCE CRITERION ON INNER ITERATIONS))
+      END