From 7dfcc480ba1e19bd3232349fc733caef94034292 Mon Sep 17 00:00:00 2001
From: stainer_t <thomas.stainer@oecd-nea.org>
Date: Mon, 8 Sep 2025 13:48:49 +0200
Subject: Initial commit from Polytechnique Montreal

---
 Dragon/src/SPHCPO.f | 431 ++++++++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 431 insertions(+)
 create mode 100644 Dragon/src/SPHCPO.f

(limited to 'Dragon/src/SPHCPO.f')

diff --git a/Dragon/src/SPHCPO.f b/Dragon/src/SPHCPO.f
new file mode 100644
index 0000000..fa474a3
--- /dev/null
+++ b/Dragon/src/SPHCPO.f
@@ -0,0 +1,431 @@
+*DECK SPHCPO
+      SUBROUTINE SPHCPO(MAXISO,IPLIB,IPCPO,NMIL,NGRP,IMPX,ICAL,ILUPS,B2)
+*
+*-----------------------------------------------------------------------
+*
+*Purpose:
+* Extract a Microlib corresponding to an elementary calculation in a
+* Multicompo.
+*
+*Copyright:
+* Copyright (C) 2012 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
+* MAXISO  maximum allocated space for output microlib TOC information.
+* IPLIB   address of the output microlib LCM object.
+* IPCPO   address of the multicompo object.
+* NMIL    number of mixtures in the elementary calculation.
+* NGRP    number of energy groups.
+* IMPX    print parameter (equal to zero for no print).
+* ICAL    index of the elementary calculation being considered.
+* ILUPS   up-scattering removing flag (=1 to remove up-scattering from
+*         output cross-sections).
+*
+*Parameters: output
+* B2      buckling recovered from the Multicompo.
+*
+*-----------------------------------------------------------------------
+*
+      USE GANLIB
+*----
+*  SUBROUTINE ARGUMENTS
+*----
+      TYPE(C_PTR) IPLIB,IPCPO
+      INTEGER MAXISO,NMIL,NGRP,IMPX,ICAL,ILUPS
+      REAL B2
+*----
+*  LOCAL VARIABLES
+*----
+      PARAMETER (NSTATE=40,MAXED=50,IOUT=6)
+      CHARACTER TEXT12*12,HSMG*131,HVECT1(MAXED)*8,HVECT2(MAXED)*8
+      INTEGER ISTATE(NSTATE)
+      TYPE(C_PTR) JPLIB,KPLIB,JPCPO,KPCPO,LPCPO,MPCPO,NPCPO,OPCPO
+      INTEGER, ALLOCATABLE, DIMENSION(:) :: ITYP1,ITOD1,IMIX2,ITYP2,
+     1 ITOD2,MILVO,MUP
+      INTEGER, ALLOCATABLE, DIMENSION(:,:) :: HUSE1,HNAM1,HUSE2,HNAM2
+      REAL, ALLOCATABLE, DIMENSION(:) :: DENS1,TEMP1,VOL1,DENS2,TEMP2,
+     1 VOL2,ENER,DELT,VOLMI2,GAR4
+      REAL, ALLOCATABLE, DIMENSION(:,:) :: GAR6
+      REAL, ALLOCATABLE, DIMENSION(:,:,:) :: ADF2
+      CHARACTER(LEN=8), ALLOCATABLE, DIMENSION(:) :: HADF,HADF2
+*----
+*  SCRATCH STORAGE ALLOCATION
+*----
+      ALLOCATE(HUSE1(3,MAXISO),HNAM1(3,MAXISO),ITYP1(MAXISO),
+     1 ITOD1(MAXISO),IMIX2(MAXISO),ITYP2(MAXISO),ITOD2(MAXISO),
+     2 HUSE2(3,MAXISO),HNAM2(3,MAXISO),MILVO(NMIL))
+      ALLOCATE(DENS1(MAXISO),TEMP1(MAXISO),VOL1(MAXISO),DENS2(MAXISO),
+     1 TEMP2(MAXISO),VOL2(MAXISO),ENER(NGRP+1),DELT(NGRP),VOLMI2(NMIL))
+*----
+*  MICROLIB INITIALIZATION
+*----
+      IF(ILUPS.EQ.1) CALL XABORT('SPHCPO: UPS OPTION NOT IMPLEMENTED.')
+      NBISO2=0
+      NCOMB2=0
+      NED2=0
+      TEXT12='L_LIBRARY'
+      CALL LCMPTC(IPLIB,'SIGNATURE',12,TEXT12)
+      DENS2(:MAXISO)=0.0
+      TEMP2(:MAXISO)=0.0
+      VOL2(:MAXISO)=0.0
+      VOLMI2(:NMIL)=0.0
+      IMIX2(:MAXISO)=0
+      ITYP2(:MAXISO)=0
+      ITOD2(:MAXISO)=0
+*----
+*  RECOVER NDEPL
+*----
+      NDEPL=0
+      CALL LCMLEN(IPCPO,'DEPL-CHAIN',ILONG,ITYLCM)
+      IF(ILONG.NE.0) THEN
+         CALL LCMSIX(IPCPO,'DEPL-CHAIN',1)
+         CALL LCMGET(IPCPO,'STATE-VECTOR',ISTATE)
+         NDEPL=ISTATE(1)
+         CALL LCMSIX(IPCPO,' ',2)
+      ENDIF
+*----
+*  DETECT DISCONTINUITY FACTORS
+*----
+      JPCPO=LCMGID(IPCPO,'MIXTURES')
+      KPCPO=LCMGIL(JPCPO,1)
+      LPCPO=LCMGID(KPCPO,'CALCULATIONS')
+      MPCPO=LCMGIL(LPCPO,1)
+      CALL LCMGET(MPCPO,'STATE-VECTOR',ISTATE)
+      IDF=ISTATE(24)
+      NTYPE=0
+      IF(IDF.GT.0) THEN
+        IF(IDF.EQ.1) THEN
+          NTYPE=2
+        ELSE IF(IDF.EQ.2) THEN
+          CALL LCMSIX(MPCPO,'MACROLIB',1)
+          CALL LCMLEN(MPCPO,'ADF',ILONG,ITYLCM)
+          IF(ILONG.EQ.0) CALL XABORT('SPHCPO: MISSING ADF DIRECTORY '
+     1    //'IN MULTICOMPO OBJECT.')
+          CALL LCMSIX(MPCPO,'ADF',1)
+          CALL LCMGET(MPCPO,'NTYPE',NTYPE)
+          CALL LCMSIX(MPCPO,' ',2)
+          CALL LCMSIX(MPCPO,' ',2)
+        ELSE
+          CALL XABORT('SPHCPO: MATRIX ADF IS NOT IMPLEMENTED.')
+        ENDIF
+      ENDIF
+      ALLOCATE(HADF(NTYPE),ADF2(NMIL,NGRP,NTYPE),GAR4(NGRP),
+     1 GAR6(NGRP,2))
+      IF((IDF.GE.2).AND.(NTYPE.GT.0)) THEN
+        CALL LCMSIX(MPCPO,'MACROLIB',1)
+        CALL LCMSIX(MPCPO,'ADF',1)
+        CALL LCMGTC(MPCPO,'HADF',8,NTYPE,HADF)
+        CALL LCMSIX(MPCPO,' ',2)
+        CALL LCMSIX(MPCPO,' ',2)
+      ENDIF
+*----
+*  LOOP OVER MICROLIB MIXTURES
+*----
+      MILVO(:NMIL)=0
+      NCOMB=0
+      JPCPO=LCMGID(IPCPO,'MIXTURES')
+      ITRANC=0
+      NDEL=0
+      NDFI=0
+      NL=0
+      NW=0
+      DO 190 IBM=1,NMIL
+      KPCPO=LCMGIL(JPCPO,IBM)
+      LPCPO=LCMGID(KPCPO,'CALCULATIONS')
+*----
+*  SELECT ICAL-TH ELEMENTARY CALCULATION
+*----
+      IF(IMPX.GT.0) THEN
+         WRITE(IOUT,'(33H SPHCPO: COMPO ACCESS FOR MIXTURE,I6,6H AND C,
+     1   10HALCULATION,I5)') IBM,ICAL
+      ENDIF
+      MPCPO=LCMGIL(LPCPO,ICAL)
+      IF(IMPX.GT.50) CALL LCMLIB(MPCPO)
+      CALL LCMGET(MPCPO,'STATE-VECTOR',ISTATE)
+      NL=ISTATE(4)
+      ITRANC=ISTATE(5)
+      NDEPL=MAX(ISTATE(11),NDEPL)
+      NDEL=ISTATE(19)
+      NDFI=ISTATE(20)
+      NW=MAX(NW,ISTATE(25))
+      IF(ISTATE(1).NE.1) CALL XABORT('SPHCPO: INVALID NUMBER OF MATERI'
+     1 //'AL MIXTURES IN THE COMPO.')
+      IF(ISTATE(3).NE.NGRP) CALL XABORT('SPHCPO: INVALID NUMBER OF ENE'
+     1 //'RGY GROUPS IN THE COMPO.')
+      IF(ISTATE(24).NE.IDF) CALL XABORT('SPHCPO: INVALID NUMBER OF DIS'
+     1 //'CONTINUITY FACTORS IN THE COMPO.')
+      NBISO1=ISTATE(2)
+      IF(NBISO1.GT.MAXISO) CALL XABORT('SPHCPO: MAXISO OVERFLOW(1).')
+      NED1=ISTATE(13)
+      IF(NED1.GT.MAXED) CALL XABORT('SPHCPO: MAXED OVERFLOW(1).')
+      CALL LCMLEN(MPCPO,'MIXTURESVOL',ILONG,ITYLCM)
+      IF(ILONG.GT.0) CALL LCMGET(MPCPO,'MIXTURESVOL',VOLMI2(IBM))
+      CALL LCMGET(MPCPO,'ISOTOPESUSED',HUSE1)
+      CALL LCMGET(MPCPO,'ISOTOPERNAME',HNAM1)
+      CALL LCMGET(MPCPO,'ISOTOPESDENS',DENS1)
+      CALL LCMGET(MPCPO,'ISOTOPESTYPE',ITYP1)
+      CALL LCMGET(MPCPO,'ISOTOPESTODO',ITOD1)
+      CALL LCMGET(MPCPO,'ISOTOPESVOL',VOL1)
+      CALL LCMGET(MPCPO,'ISOTOPESTEMP',TEMP1)
+      B2=0.0
+      CALL LCMLEN(MPCPO,'B2  B1HOM',ILONG,ITYLCM)
+      IF(ILONG.GT.0) CALL LCMGET(MPCPO,'B2  B1HOM',B2)
+      IF(NED1.GT.0) CALL LCMGTC(MPCPO,'ADDXSNAME-P0',8,NED1,HVECT1)
+      CALL LCMGET(MPCPO,'ENERGY',ENER)
+      CALL LCMGET(MPCPO,'DELTAU',DELT)
+      DO 30 IED1=1,NED1
+      DO 20 IED2=1,NED2
+      IF(HVECT1(IED1).EQ.HVECT2(IED2)) GO TO 30
+   20 CONTINUE
+      NED2=NED2+1
+      IF(NED2.GT.MAXED) CALL XABORT('SPHCPO: MAXED OVERFLOW(2).')
+      HVECT2(NED2)=HVECT1(IED1)
+   30 CONTINUE
+      IF(IBM.GT.9999) CALL XABORT('SPHCPO: IBM OVERFLOW.')
+      DO 100 ISO=1,NBISO1 ! compo isotope
+      WRITE(TEXT12,'(2A4,I4.4)') (HUSE1(I,ISO),I=1,2),IBM
+      DO 60 JSO=1,NBISO2 ! microlib isotope
+      IF((HUSE1(1,ISO).EQ.HUSE2(1,JSO)).AND.(HUSE1(2,ISO).EQ.
+     1 HUSE2(2,JSO)).AND.(IMIX2(JSO).EQ.IBM)) THEN
+         IF(ITYP1(ISO).NE.ITYP2(JSO)) THEN
+            WRITE(HSMG,500) 'ITYP',ISO,ITYP1(ISO),ITYP2(JSO)
+            CALL XABORT(HSMG)
+         ENDIF
+         JSO1=JSO
+         GO TO 90
+      ENDIF
+   60 CONTINUE
+      NBISO2=NBISO2+1
+      IF(NBISO2.GT.MAXISO) THEN
+         WRITE(IOUT,'(/16H SPHCPO: NBISO2=,I6,8H MAXISO=,I6)') NBISO2,
+     1   MAXISO
+         CALL XABORT('SPHCPO: MAXISO OVERFLOW(2).')
+      ENDIF
+      READ(TEXT12,'(3A4)') (HUSE2(I0,NBISO2),I0=1,3)
+      DO 70 I0=1,3
+      HNAM2(I0,NBISO2)=HNAM1(I0,ISO)
+   70 CONTINUE
+      IMIX2(NBISO2)=IBM
+      ITYP2(NBISO2)=ITYP1(ISO)
+      ITOD2(NBISO2)=ITOD1(ISO)
+      IF(ITYP2(NBISO2).EQ.1) ITOD2(NBISO2)=1
+      DENS2(NBISO2)=0.0
+      JSO1=NBISO2
+      IF(ITOD2(NBISO2).NE.1) THEN
+         DO 80 J=1,NCOMB
+         IF(IBM.EQ.MILVO(J)) GO TO 90
+   80    CONTINUE
+         NCOMB=NCOMB+1
+         IF(NCOMB.GT.NMIL) CALL XABORT('SPHCPO: MILVO OVERFLOW.')
+         MILVO(NCOMB)=IBM
+      ENDIF
+   90 DENS2(JSO1)=DENS1(ISO)
+      VOL2(JSO1)=VOL1(ISO)
+      TEMP2(JSO1)=TEMP1(ISO)
+  100 CONTINUE
+*----
+*  PROCESS ISOTOPE DIRECTORIES FOR MICROLIB MIXTURE IBM
+*----
+      JPLIB=LCMLID(IPLIB,'ISOTOPESLIST',NBISO2)
+      DO 180 ISO=1,NBISO2 ! microlib isotope
+      IF(IMIX2(ISO).NE.IBM) GO TO 180
+      DO 120 JSO=1,NBISO1 ! compo isotope
+      IF((HUSE1(1,JSO).EQ.HUSE2(1,ISO)).AND.(HUSE1(2,JSO).EQ.
+     1 HUSE2(2,ISO))) THEN
+        JSO1=JSO
+        GO TO 130
+      ENDIF
+  120 CONTINUE
+      WRITE(TEXT12,'(3A4)') (HUSE2(I0,ISO),I0=1,3)
+      CALL XABORT('SPHCPO: UNABLE TO FIND '//TEXT12//'.')
+  130 KPLIB=LCMDIL(JPLIB,ISO) ! set ISO-th isotope
+      MPCPO=LCMGIL(LPCPO,ICAL)
+      NPCPO=LCMGID(MPCPO,'ISOTOPESLIST')
+      CALL LCMLEL(NPCPO,JSO1,ILENG,ITYLCM)
+      IF(ILENG.NE.0) THEN
+         OPCPO=LCMGIL(NPCPO,JSO1) ! set JSO1-th isotope
+         CALL LCMEQU(OPCPO,KPLIB)
+      ENDIF
+  180 CONTINUE
+*----
+*  PROCESS ADF INFORMATION
+*----
+      IF(IDF.GT.0) THEN
+        CALL LCMSIX(MPCPO,'MACROLIB',1)
+        CALL LCMLEN(MPCPO,'ADF',ILONG,ITYLCM)
+        IF(ILONG.NE.0) THEN
+          CALL LCMSIX(MPCPO,'ADF',1)
+          IF(IDF.EQ.1) THEN
+            GAR6(:NGRP,:2)=0.0
+            CALL LCMGET(MPCPO,'ALBS00',GAR6)
+            DO IGR=1,NGRP
+              ADF2(IBM,IGR,:2)=GAR6(IGR,:2)
+            ENDDO
+          ELSE IF((IDF.EQ.2).OR.(IDF.EQ.3)) THEN
+            CALL LCMGET(MPCPO,'NTYPE',NTYPE2)
+            ALLOCATE(HADF2(NTYPE2))
+            CALL LCMGTC(MPCPO,'HADF',8,NTYPE2,HADF2)
+            IF(NTYPE2.EQ.1) THEN
+*             assign the same ADF to all sides.
+              CALL LCMLEN(MPCPO,HADF2(1),ILONG,ITYLCM)
+              IF(ILONG.NE.NGRP) CALL XABORT('SPHCPO: INVALID ADF LENGT'
+     1        //'H(1).')
+              CALL LCMGET(MPCPO,HADF2(1),GAR4)
+              DO ITY=1,NTYPE
+                DO IGR=1,NGRP
+                  ADF2(IBM,IGR,ITY)=GAR4(IGR)
+                ENDDO
+              ENDDO
+            ELSE
+              IF(NTYPE2.GT.NTYPE) CALL XABORT('SPHCPO: NTYPE OVERFLOW.')
+              DO ITY2=1,NTYPE2
+                ITY=0
+                DO JTY=1,NTYPE
+                  IF(HADF2(ITY2).EQ.HADF(JTY)) THEN
+                    ITY=JTY
+                    GO TO 185
+                  ENDIF
+                ENDDO
+                WRITE(HSMG,'(18HSPHCPO: ADF NAMED ,A,11H NOT FOUND.)')
+     1          TRIM(HADF2(ITY2))
+                CALL XABORT(HSMG)
+  185           CALL LCMLEN(MPCPO,HADF2(ITY2),ILONG,ITYLCM)
+                IF(ILONG.NE.NGRP) CALL XABORT('SPHCPO: INVALID ADF LEN'
+     1          //'GTH(2).')
+                CALL LCMGET(MPCPO,HADF2(ITY2),GAR4)
+                DO IGR=1,NGRP
+                  ADF2(IBM,IGR,ITY)=GAR4(IGR)
+                ENDDO
+              ENDDO
+            ENDIF
+            DEALLOCATE(HADF2)
+          ENDIF
+          CALL LCMSIX(MPCPO,' ',2)
+        ENDIF
+        CALL LCMSIX(MPCPO,' ',2)
+      ENDIF
+  190 CONTINUE
+*----
+*  MICROLIB FINALIZATION
+*----
+      ISTATE(:NSTATE)=0
+      ISTATE(1)=NMIL
+      ISTATE(2)=NBISO2
+      ISTATE(3)=NGRP
+      ISTATE(4)=NL
+      ISTATE(5)=ITRANC
+      ISTATE(7)=1
+      ISTATE(11)=NDEPL
+      ISTATE(12)=NCOMB+NCOMB2
+      ISTATE(13)=NED2
+      ISTATE(14)=NMIL
+      ISTATE(18)=1
+      ISTATE(19)=NDEL
+      ISTATE(20)=NDFI
+      ISTATE(22)=MAXISO/NMIL
+      ISTATE(25)=NW
+      IF(NBISO2.EQ.0) CALL XABORT('SPHCPO: NBISO2=0.')
+      CALL LCMPUT(IPLIB,'STATE-VECTOR',NSTATE,1,ISTATE)
+      CALL LCMPUT(IPLIB,'MIXTURESVOL',NMIL,2,VOLMI2)
+      CALL LCMPUT(IPLIB,'ISOTOPESUSED',3*NBISO2,3,HUSE2)
+      CALL LCMPUT(IPLIB,'ISOTOPERNAME',3*NBISO2,3,HNAM2)
+      CALL LCMPUT(IPLIB,'ISOTOPESDENS',NBISO2,2,DENS2)
+      CALL LCMPUT(IPLIB,'ISOTOPESMIX',NBISO2,1,IMIX2)
+      CALL LCMPUT(IPLIB,'ISOTOPESTYPE',NBISO2,1,ITYP2)
+      CALL LCMPUT(IPLIB,'ISOTOPESTODO',NBISO2,1,ITOD2)
+      CALL LCMPUT(IPLIB,'ISOTOPESVOL',NBISO2,2,VOL2)
+      CALL LCMPUT(IPLIB,'ISOTOPESTEMP',NBISO2,2,TEMP2)
+      IF(NED2.GT.0) CALL LCMPTC(IPLIB,'ADDXSNAME-P0',8,NED2,HVECT2)
+      CALL LCMPUT(IPLIB,'ENERGY',NGRP+1,2,ENER)
+      CALL LCMPUT(IPLIB,'DELTAU',NGRP,2,DELT)
+      IF(IMPX.GT.5) CALL LCMLIB(IPLIB)
+*----
+*  RECOVER GENERAL INFORMATION FROM MIXTURE 1
+*----
+      B2=0.0
+      KPCPO=LCMGIL(JPCPO,1)
+      LPCPO=LCMGID(KPCPO,'CALCULATIONS')
+      MPCPO=LCMGIL(LPCPO,ICAL)
+      CALL LCMLEN(MPCPO,'K-EFFECTIVE',ILENG,ITYLCM)
+      IF(ILENG.EQ.1) THEN
+        CALL LCMGET(MPCPO,'K-EFFECTIVE',FLOTT)
+        CALL LCMPUT(IPLIB,'K-EFFECTIVE',1,2,FLOTT)
+        IF(IMPX.GT.1) THEN
+          WRITE(6,'(22H SPHCPO: K-EFFECTIVE =,1P,E13.6)') FLOTT
+        ENDIF
+      ENDIF
+      CALL LCMLEN(MPCPO,'K-INFINITY',ILENG,ITYLCM)
+      IF(ILENG.EQ.1) THEN
+        CALL LCMGET(MPCPO,'K-INFINITY',FLOTT)
+        CALL LCMPUT(IPLIB,'K-INFINITY',1,2,FLOTT)
+        IF(IMPX.GT.1) THEN
+          WRITE(6,'(21H SPHCPO: K-INFINITY =,1P,E13.6)') FLOTT
+        ENDIF
+      ENDIF
+      CALL LCMLEN(MPCPO,'B2  B1HOM',ILENG,ITYLCM)
+      IF(ILENG.EQ.1) THEN
+        CALL LCMGET(MPCPO,'B2  B1HOM',B2)
+        CALL LCMPUT(IPLIB,'B2  B1HOM',1,2,B2)
+        IF(IMPX.GT.1) THEN
+          WRITE(6,'(13H SPHCPO: B2 =,1P,E14.6)') B2
+        ENDIF
+      ENDIF
+*----
+*  BUILD EMBEDDED MACROLIB
+*----
+      ALLOCATE(MUP(NMIL))
+      MUP(:NMIL)=1
+      CALL SPHEMB(IPLIB,IPCPO,NGRP,NMIL,MUP,IMPX)
+      DEALLOCATE(MUP)
+*----
+*  WRITE ADF INFORMATION
+*----
+      IF(IDF.EQ.1) THEN
+        CALL LCMSIX(IPLIB,'MACROLIB',1)
+        CALL LCMSIX(IPLIB,'ADF',1)
+        CALL LCMPUT(IPLIB,'ALBS00',NMIL*NGRP*2,2,ADF2(1,1,1))
+        CALL LCMSIX(IPLIB,' ',2)
+        CALL LCMGET(IPLIB,'STATE-VECTOR',ISTATE)
+        ISTATE(12)=IDF
+        CALL LCMPUT(IPLIB,'STATE-VECTOR',NSTATE,1,ISTATE)
+        CALL LCMSIX(IPLIB,' ',2)
+      ELSE IF((IDF.EQ.2).OR.(IDF.EQ.3)) THEN
+        CALL LCMSIX(IPLIB,'MACROLIB',1)
+        CALL LCMSIX(IPLIB,'ADF',1)
+        DO ITYPE=1,NTYPE
+          CALL LCMPUT(IPLIB,HADF(ITYPE),NMIL*NGRP,2,
+     1    ADF2(1,1,ITYPE))
+        ENDDO
+        CALL LCMSIX(IPLIB,' ',2)
+        CALL LCMGET(IPLIB,'STATE-VECTOR',ISTATE)
+        ISTATE(12)=IDF
+        CALL LCMPUT(IPLIB,'STATE-VECTOR',NSTATE,1,ISTATE)
+        CALL LCMSIX(IPLIB,' ',2)
+        IF(IMPX.GT.1) THEN
+          DO IBM=1,NMIL
+            WRITE(6,'(/40H SPHCPO: DISCONTINUITY FACTORS - MIXTURE,I5)')
+     1      IBM
+            DO ITYPE=1,NTYPE
+              WRITE(6,'(1X,A,1H:,1P,(5X,10E12.4))') TRIM(HADF(ITYPE)),
+     1        (ADF2(IBM,IGR,ITYPE),IGR=1,NGRP)
+            ENDDO
+          ENDDO
+        ENDIF
+      ENDIF
+*----
+*  SCRATCH STORAGE DEALLOCATION
+*----
+      DEALLOCATE(GAR6,GAR4,ADF2,HADF)
+      DEALLOCATE(VOLMI2,DELT,ENER,VOL2,TEMP2,DENS2,VOL1,TEMP1,DENS1)
+      DEALLOCATE(MILVO,HNAM2,HUSE2,ITOD2,ITYP2,IMIX2,ITOD1,ITYP1,HNAM1,
+     1 HUSE1)
+      RETURN
+*
+  500 FORMAT(8HSPHCPO: ,A,1H(,I4,2H)=,2I5)
+      END
-- 
cgit v1.2.3