Initial commit from Polytechnique Montreal

1 files changed, 346 insertions, 0 deletions

diff --git a/Dragon/src/USS.f b/Dragon/src/USS.f
new file mode 100644
index 0000000..33ea918
--- /dev/null
+++ b/Dragon/src/USS.f
@@ -0,0 +1,346 @@
+*DECK USS
+      SUBROUTINE USS(NENTRY,HENTRY,IENTRY,JENTRY,KENTRY)
+*
+*-----------------------------------------------------------------------
+*
+*Purpose:
+* Universal self-shielding operator.
+*
+*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
+*
+*Parameters: input/output
+* NENTRY  number of LCM objects or files used by the operator.
+* HENTRY  name of each LCM object or file:
+*         HENTRY(1) creation or modification type(L_LIBRARY) (no
+*         subgroups);
+*         HENTRY(2) read-only type(L_LIBRARY) (with subgroups);
+*         HENTRY(3) read-only type(L_TRACK);
+*         HENTRY(4) optional read-only sequential binary tracking file.
+* 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,MAXRSS=300,IOUT=6)
+      TYPE(C_PTR) IPLI0,IPLIB,IPTRK
+      CHARACTER TEXT4*4,HSIGN*12,TEXT12*12,HSMG*131,CDOOR*12,TITR*72,
+     1 HISOT*12
+      DOUBLE PRECISION DFLOTT
+      LOGICAL LEAKSW,LFLAT
+      INTEGER IGP(NSTATE),IPAR(NSTATE),IPAS(NSTATE),IRSS(MAXRSS)
+*----
+*  PARAMETER VALIDATION.
+*----
+      IF(NENTRY.LE.2) CALL XABORT('USS: THREE PARAMETERS EXPECTED.')
+      IF((IENTRY(1).NE.1).AND.(IENTRY(1).NE.2)) CALL XABORT('USS: LI'
+     1 //'NKED LIST OR XSM FILE EXPECTED AT LHS.')
+      IF((JENTRY(1).NE.0).AND.(JENTRY(1).NE.1)) CALL XABORT('USS: EN'
+     1 //'TRY IN CREATE OR MODIFICATION MODE EXPECTED.')
+      IF((JENTRY(2).NE.2).OR.((IENTRY(2).NE.1).AND.(IENTRY(2).NE.2)))
+     1 CALL XABORT('USS: LINKED LIST OR XSM FILE IN READ-ONLY MODE E'
+     2 //'XPECTED AT FIRST RHS.')
+      IF((JENTRY(3).NE.2).OR.((IENTRY(3).NE.1).AND.(IENTRY(3).NE.2)))
+     1 CALL XABORT('USS: LINKED LIST OR XSM FILE IN READ-ONLY MODE E'
+     2 //'XPECTED AT SECOND RHS.')
+      CALL LCMGTC(KENTRY(3),'SIGNATURE',12,HSIGN)
+      IF(HSIGN.NE.'L_TRACK') THEN
+         TEXT12=HENTRY(3)
+         CALL XABORT('USS: SIGNATURE OF '//TEXT12//' IS '//HSIGN//
+     1   '. L_TRACK EXPECTED.')
+      ENDIF
+      CALL LCMGTC(KENTRY(3),'TRACK-TYPE',12,CDOOR)
+      IPLI0=KENTRY(1)
+      IPLIB=KENTRY(2)
+      IPTRK=KENTRY(3)
+      INDREC=0
+      IF(JENTRY(1).EQ.0) THEN
+         INDREC=1
+         HSIGN='L_LIBRARY'
+         CALL LCMPTC(IPLI0,'SIGNATURE',12,HSIGN)
+      ELSE IF(JENTRY(1).EQ.1) THEN
+         CALL LCMGTC(IPLI0,'SIGNATURE',12,HSIGN)
+         IF(HSIGN.NE.'L_LIBRARY') THEN
+            TEXT12=HENTRY(1)
+            CALL XABORT('USS: SIGNATURE OF '//TEXT12//' IS '//HSIGN//
+     1      '. L_LIBRARY EXPECTED.')
+         ENDIF
+         INDREC=2
+      ENDIF
+*----
+*  RECOVER TABULATED FUNCTIONS.
+*----
+      CALL XDRTA2
+*----
+*  RECOVER TRACKING FILE INFORMATION.
+*----
+      IFTRAK=0
+      IF(NENTRY.GE.4) THEN
+        IF(IENTRY(4).EQ.3) THEN
+          IF(JENTRY(4).NE.2) CALL XABORT('USS: BINARY TRACKING FILE NA'
+     1    //'MED '//HENTRY(4)//' IS NOT IN REAL-ONLY MODE.')
+          IFTRAK=FILUNIT(KENTRY(4))
+        ENDIF
+      ENDIF
+*----
+*  RECOVER GENERAL TRACKING INFORMATION.
+*----
+      CALL LCMGET(IPTRK,'STATE-VECTOR',IGP)
+      NREG=IGP(1)
+      NUN=IGP(2)
+      LEAKSW=IGP(3).EQ.0
+      CALL LCMLEN(IPTRK,'TITLE',LENGT,ITYLCM)
+      IF(LENGT.GT.0) THEN
+         CALL LCMGTC(IPTRK,'TITLE',72,TITR)
+      ELSE
+         TITR='*** NO TITLE PROVIDED ***'
+      ENDIF
+*----
+*  RECOVER INTERNAL LIBRARY PARAMETERS.
+*----
+      CALL LCMGTC(IPLIB,'SIGNATURE',12,HSIGN)
+      IF(HSIGN.NE.'L_LIBRARY') THEN
+         TEXT12=HENTRY(2)
+         CALL XABORT('USS: SIGNATURE OF '//TEXT12//' IS '//HSIGN//
+     1   '. L_LIBRARY EXPECTED.')
+      ENDIF
+      CALL LCMGET(IPLIB,'STATE-VECTOR',IPAR)
+      NBISO=IPAR(2)
+      NGRP=IPAR(3)
+      NL=IPAR(4)
+      ITRANC=IPAR(5)
+      IGRMIN=IPAR(9)+1
+      IGRMAX=IPAR(10)
+      NED=IPAR(13)
+      NBMIX=IPAR(14)
+      NRES=IPAR(15)
+      ISUBG=IPAR(17)
+      NDEL=IPAR(19)
+      IF(IGP(4).GT.NBMIX) THEN
+         WRITE(HSMG,'(45HUSS: THE NUMBER OF MIXTURES IN THE TRACKING (,
+     1   I5,55H) IS GREATER THAN THE NUMBER OF MIXTURES IN THE INTERNA,
+     2   11HL LIBRARY (,I5,2H).)') IGP(4),NBMIX
+         CALL XABORT(HSMG)
+      ELSE IF((ISUBG.LE.0).OR.(ISUBG.EQ.2)) THEN
+         CALL XABORT('USS: THE INPUT INTERNAL LIBRARY HAS NO SUBGROUP'
+     1   //'S.')
+      ENDIF
+      IF(INDREC.EQ.2) THEN
+         CALL LCMGET(IPLI0,'STATE-VECTOR',IPAR)
+         IF(IPAR(2).NE.NBISO) CALL XABORT('USS: INVALID LIBRARY.')
+      ENDIF
+*
+      IMPX=1
+      LFLAT=.FALSE.
+      CALL LCMLEN(IPLI0,'SHIBA_SG',LENLCM,ITYLCM)
+      IF(LENLCM.NE.0) THEN
+         CALL LCMSIX(IPLI0,'SHIBA_SG',1)
+         CALL LCMGET(IPLI0,'STATE-VECTOR',IPAS)
+         CALL LCMSIX(IPLI0,' ',2)
+         IGRMIN=IPAS(1)
+         IGRMAX=IPAS(2)
+         KSPH=IPAS(3)
+         ITRANZ=IPAS(4)
+         NPASS=IPAS(5)
+         IPHASE=IPAS(6)
+         ICALC=IPAS(8)
+         ICORR=IPAS(9)
+         MAXST=IPAS(10)
+      ELSE
+         KSPH=1
+         ITRANZ=ITRANC
+         NPASS=2
+         IF(CDOOR.EQ.'SYBIL') THEN
+            IPHASE=2
+         ELSE IF(CDOOR.EQ.'EXCELL') THEN
+            IPHASE=2
+         ELSE
+            IPHASE=1
+         ENDIF
+         ICALC=0
+         ICORR=0
+         MAXST=50
+         IF(ISUBG.EQ.6) MAXST=20
+      ENDIF
+*
+   30 CALL REDGET(INDIC,NITMA,FLOTT,TEXT4,DFLOTT)
+      IF(INDIC.EQ.10) GO TO 70
+      IF(INDIC.NE.3) CALL XABORT('USS: CHARACTER DATA EXPECTED(1).')
+      IF(TEXT4.EQ.'EDIT') THEN
+         CALL REDGET(INDIC,IMPX,FLOTT,TEXT4,DFLOTT)
+         IF(INDIC.NE.1) CALL XABORT('USS: INTEGER DATA EXPECTED(1).')
+      ELSE IF(TEXT4.EQ.'GRMI') THEN
+         CALL REDGET(ITYPLU,IGRMIN,FLOTT,TEXT4,DFLOTT)
+         IF(ITYPLU.NE.1) CALL XABORT('USS: INTEGER DATA EXPECTED(2).')
+      ELSE IF(TEXT4.EQ.'GRMA') THEN
+         CALL REDGET(ITYPLU,IGRMAX,FLOTT,TEXT4,DFLOTT)
+         IF(ITYPLU.NE.1) CALL XABORT('USS: INTEGER DATA EXPECTED(3).')
+         IF(IGRMAX.GT.NGRP) THEN
+            CALL XABORT('USS: ILLEGAL NUMBER OF GROUP IN LIBRARY.')
+         ENDIF
+      ELSE IF(TEXT4.EQ.'NOSP') THEN
+         KSPH=0
+      ELSE IF(TEXT4.EQ.'NOTR') THEN
+         ITRANZ=0
+      ELSE IF(TEXT4.EQ.'TRAN') THEN
+         ITRANZ=1
+      ELSE IF(TEXT4.EQ.'PASS') THEN
+         CALL REDGET(ITYPLU,NPASS,FLOTT,TEXT4,DFLOTT)
+         IF(ITYPLU.NE.1) CALL XABORT('USS: INTEGER DATA EXPECTED(5).')
+         IF(NPASS.LE.0) CALL XABORT('USS: POSITIVE PASS EXPECTED.')
+      ELSE IF(TEXT4.EQ.'ARM') THEN
+         IPHASE=1
+      ELSE IF(TEXT4.EQ.'PIJ') THEN
+         IPHASE=2
+      ELSE IF(TEXT4.EQ.'CALC') THEN
+         ICALC=1
+         CALL LCMSIX(IPLI0,'SHIBA_SG',1)
+         CALL LCMSIX(IPLI0,'-DATA-CALC-',1)
+         CALL REDGET(ITYPLU,NITMA,FLOTT,TEXT12,DFLOTT)
+         IF(ITYPLU.NE.3) CALL XABORT('USS: CHARACTER DATA EXPECTED(2).')
+   40    IF(TEXT12.EQ.'ENDC') THEN
+            CALL LCMSIX(IPLI0,' ',2)
+            CALL LCMSIX(IPLI0,' ',2)
+            GO TO 30
+         ENDIF
+         IF(TEXT12.NE.'REGI') CALL XABORT('USS: REGI KEYWORD EXPECTED.')
+         CALL REDGET(ITYPLU,NITMA,FLOTT,TEXT12,DFLOTT)
+         IF(ITYPLU.NE.3) CALL XABORT('USS: CHARACTER DATA EXPECTED(3).')
+         IF(TEXT12(5:).NE.' ') CALL XABORT('USS: 4-CHARACTER NAME EXPE'
+     1   //'CTED.')
+         CALL LCMSIX(IPLI0,TEXT12(:4),1)
+         CALL REDGET(ITYPLU,NITMA,FLOTT,TEXT12,DFLOTT)
+         IF(ITYPLU.NE.3) CALL XABORT('USS: CHARACTER DATA EXPECTED(4).')
+   50    IF((TEXT12.EQ.'ENDC').OR.(TEXT12.EQ.'REGI')) THEN
+            CALL LCMSIX(IPLI0,' ',2)
+            GO TO 40
+         ENDIF
+         HISOT=TEXT12
+         NRSS=0
+   60    CALL REDGET(ITYPLU,NITMA,FLOTT,TEXT12,DFLOTT)
+         IF(ITYPLU.EQ.3) THEN
+            IF(TEXT12.EQ.'ALL') THEN
+               NRSS=1
+               IRSS(1)=-999
+               GO TO 60
+            ENDIF
+            IF(NRSS.EQ.0) CALL XABORT('USS: INTEGER ARRAY EXPECTED.')
+            CALL LCMPUT(IPLI0,HISOT,NRSS,1,IRSS)
+            GO TO 50
+         ENDIF
+         IF(ITYPLU.NE.1) CALL XABORT('USS: INTEGER DATA EXPECTED(6).')
+         NRSS=NRSS+1
+         IF(NRSS.GT.MAXRSS) CALL XABORT('USS: MAXRSS OVERFLOW.')
+         IF((NITMA.LE.0).OR.(NITMA.GT.NBMIX)) THEN
+            WRITE(HSMG,'(43HUSS: REGI KEYWORD -- INVALID MIXTURE INDEX=,
+     1      I5,1H.)') NITMA
+            CALL XABORT(HSMG)
+         ENDIF
+         IRSS(NRSS)=NITMA
+         GO TO 60
+      ELSE IF(TEXT4.EQ.'NOCO') THEN
+         ICORR=1
+      ELSE IF(TEXT4.EQ.'MAXS') THEN
+         CALL REDGET(ITYPLU,MAXST,FLOTT,TEXT4,DFLOTT)
+         IF(ITYPLU.NE.1) CALL XABORT('USS: INTEGER DATA EXPECTED(7).')
+      ELSE IF(TEXT4.EQ.'FLAT') THEN
+         IF(INDREC.EQ.1) CALL XABORT('USS: OUTPUT MICROLIB IN MODIFICA'
+     1   //'TION MODE EXPECTED.')
+         LFLAT=.TRUE.
+      ELSE IF(TEXT4.EQ.';') THEN
+         GO TO 70
+      ELSE
+         CALL XABORT('USS: '//TEXT4//' IS AN INVALID KEY WORD.')
+      ENDIF
+      GO TO 30
+*----
+*  THE NUMBER OF UNKNOWNS OF A CURRENT-BASED MULTICELL ITERATION IS
+*  INCREASED TO HOLD INTERFACE CURRENT COMPONENTS.
+*----
+   70 IF(IPHASE.EQ.1) THEN
+        IF(CDOOR.EQ.'SYBIL') NUN=NUN+IGP(9)
+        IF((CDOOR.EQ.'EXCELL').AND.(IGP(7).EQ.5)) NUN=NUN+IGP(28)
+      ENDIF
+*----
+*  CALL USS: DRIVER.
+*----
+      IF(IMPX.GT.0) THEN
+        IF(INDREC.EQ.1) WRITE(IOUT,100)
+        WRITE(IOUT,110) TITR,CDOOR(:9),ISUBG,IGRMIN,IGRMAX,KSPH,ITRANZ,
+     1  NPASS,IPHASE,ICALC,ICORR,MAXST
+      ENDIF
+*----
+*  PERFORM SELF-SHIELDING.
+*----
+      IF(CDOOR.EQ.'MCCG') THEN
+         CALL LCMLEN(IPTRK,'KEYFLX',LKFL,ITYLCM)
+         NFUNL=LKFL/NREG
+      ELSE
+         NFUNL=1
+      ENDIF
+      CALL USSDRV(IPLI0,IPTRK,IPLIB,IFTRAK,INDREC,CDOOR,IMPX,IGRMIN,
+     1 IGRMAX,NGRP,NBMIX,NREG,NUN,NBISO,NL,NED,NDEL,LEAKSW,ITRANZ,
+     2 IPHASE,TITR,KSPH,NRES,NPASS,ICALC,ICORR,ISUBG,MAXST,LFLAT)
+*----
+*  STORE THE GENERAL SHELF-SHIELDING PARAMETERS.
+*----
+      IPAS(:NSTATE)=0
+      IPAS(1)=IGRMIN
+      IPAS(2)=IGRMAX
+      IPAS(3)=KSPH
+      IPAS(4)=ITRANZ
+      IPAS(5)=NPASS
+      IPAS(6)=IPHASE
+      IPAS(8)=ICALC
+      IPAS(9)=ICORR
+      IPAS(10)=MAXST
+      CALL LCMSIX(IPLI0,'SHIBA_SG',1)
+      CALL LCMPUT(IPLI0,'STATE-VECTOR',NSTATE,1,IPAS)
+      CALL LCMSIX(IPLI0,' ',2)
+      RETURN
+*
+  100 FORMAT(1H1,23HUU   UU  SSSSS   SSSSS ,107(1H*)/
+     1 25H UU   UU SSSSSSS SSSSSSS ,63(1H*),
+     2 43H UNIVERSAL SELF-SHIELDING MODEL.  A. HEBERT/
+     3 24H UU   UU SS   SS SS   SS/21H UU   UU  SSS     SSS/
+     4 23H UU   UU    SSS     SSS/24H UU   UU SS   SS SS   SS/
+     5 24H UUUUUUU SSSSSSS SSSSSSS/23H  UUUUU   SSSSS   SSSSS/)
+  110 FORMAT(/1X,A72//8H OPTIONS/8H -------/
+     1 7H CDOOR ,A8,30H   (NAME OF THE SOLUTION DOOR)/
+     2 7H ISUBG ,I8,47H   (=1: SUBG; =3: PT; =4: PTSL; =5: PTMC; =6: R,
+     3 3HSE)/
+     4 7H IGRMIN,I8,27H   (FIRST GROUP TO PROCESS)/
+     5 7H IGRMAX,I8,34H   (MOST THERMAL GROUP TO PROCESS)/
+     6 7H KSPH  ,I8,47H   (=0: NO SPH CORRECTION; =1: SPH CORRECTION I,
+     7 19HN RESONANT REGIONS)/
+     8 7H ITRANZ,I8,45H   (0=NO TRANSPORT CORRECTION/1=APOLLO TYPE/2,
+     9 57H=RECOVER FROM LIBRARY/3=WIMS-D TYPE/4=LEAKAGE CORRECTION)/
+     1 7H NPASS ,I8,31H   (NUMBER OF OUTER ITERATIONS)/
+     2 7H IPHASE,I8,37H   (=1: NATIVE ASSEMBLY; =2: USE PIJ)/
+     3 7H ICALC ,I8,48H   (=0: NO &CALC DATA; =1: &CALC DATA AVAILABLE)/
+     4 7H ICORR ,I8,47H   (=1: SUPPRESS MUTUAL RESONANCE SHIELDING EFF,
+     5 4HECT)/
+     6 7H MAXST ,I8,36H   (MAXIMUM NUMBER OF ST ITERATIONS))
+      END