diff options
| author | stainer_t <thomas.stainer@oecd-nea.org> | 2025-09-08 13:48:49 +0200 |
|---|---|---|
| committer | stainer_t <thomas.stainer@oecd-nea.org> | 2025-09-08 13:48:49 +0200 |
| commit | 7dfcc480ba1e19bd3232349fc733caef94034292 (patch) | |
| tree | 03ee104eb8846d5cc1a981d267687a729185d3f3 /Donjon/src/DSPH.f | |
Initial commit from Polytechnique Montreal
Diffstat (limited to 'Donjon/src/DSPH.f')
| -rw-r--r-- | Donjon/src/DSPH.f | 544 |
1 files changed, 544 insertions, 0 deletions
diff --git a/Donjon/src/DSPH.f b/Donjon/src/DSPH.f new file mode 100644 index 0000000..f52b811 --- /dev/null +++ b/Donjon/src/DSPH.f @@ -0,0 +1,544 @@ +*DECK DSPH
+ SUBROUTINE DSPH(NENTRY,HENTRY,IENTRY,JENTRY,KENTRY)
+*
+*-----------------------------------------------------------------------
+*
+*Purpose:
+* Create a delta Macrolib with respect to a SPH correction.
+*
+*Copyright:
+* Copyright (C) 2017 Ecole Polytechnique de Montreal.
+*
+*Author(s):
+* A. Hebert
+*
+*Parameters: input
+* NENTRY number of data structures transfered to this module.
+* HENTRY name of the data structures.
+* IENTRY data structure type where:
+* IENTRY=1 for LCM memory object;
+* IENTRY=2 for XSM file;
+* IENTRY=3 for sequential binary file;
+* IENTRY=4 for sequential ASCII file.
+* JENTRY access permission for the data structure where:
+* JENTRY=0 for a data structure in creation mode;
+* JENTRY=1 for a data structure in modifications mode;
+* JENTRY=2 for a data structure in read-only mode.
+* KENTRY data structure pointer.
+*
+*-----------------------------------------------------------------------
+*
+ USE GANLIB
+*----
+* SUBROUTINE ARGUMENTS
+*----
+ INTEGER NENTRY,IENTRY(NENTRY),JENTRY(NENTRY)
+ CHARACTER HENTRY(NENTRY)*12
+ TYPE(C_PTR) KENTRY(NENTRY)
+*----
+* LOCAL VARIABLES
+*----
+ PARAMETER(NSTATE=40)
+ CHARACTER HSIGN*12,TEXT2*2,TEXT8*8,TEXT12*12
+ DOUBLE PRECISION DFLOTT
+ INTEGER ISTATE(NSTATE)
+ DOUBLE PRECISION OPTPRR(NSTATE)
+ TYPE(C_PTR) IPOPT,IPNEW,IPOLD,JPNEW,JPOLD,KPNEW,KPOLD,LPNEW,MPNEW
+ INTEGER, ALLOCATABLE, DIMENSION(:) :: IJJ,NJJ,IPOS
+ INTEGER, ALLOCATABLE, DIMENSION(:,:) :: IHEDIT
+ REAL, ALLOCATABLE, DIMENSION(:) :: DIFHOM,GAR,PER,GAR1,PER1
+ REAL, ALLOCATABLE, DIMENSION(:,:) :: SPH,GAR2,PER2,ALBP,PALBP
+ REAL, ALLOCATABLE, DIMENSION(:,:,:) :: PSIGS
+ DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: VARV
+*----
+* PARAMETER VALIDATION
+*----
+ IF(NENTRY.NE.3)CALL XABORT('DSPH: THREE PARAMETERS EXPECTED.')
+ IF((IENTRY(1).NE.1).AND.(IENTRY(1).NE.2))CALL XABORT('@DSPH'
+ 1 //': LCM OBJECT EXPECTED AT LHS.')
+ IF(JENTRY(1).EQ.0)THEN
+ HSIGN='L_MACROLIB'
+ CALL LCMPTC(KENTRY(1),'SIGNATURE',12,HSIGN)
+ ELSE
+ CALL XABORT('DSPH: EMPTY DELTA MACROLIB EXPECTED AT LHS.')
+ ENDIF
+ IPNEW=KENTRY(1)
+ IF((IENTRY(2).NE.1).AND.(IENTRY(2).NE.2))CALL XABORT('DSPH: LCM '
+ 1 //'OBJECT EXPECTED AT LHS.')
+ IF(JENTRY(2).EQ.0)THEN
+ HSIGN='L_OPTIMIZE'
+ CALL LCMPTC(KENTRY(2),'SIGNATURE',12,HSIGN)
+ ELSE IF(JENTRY(2).EQ.1)THEN
+ CALL LCMGTC(KENTRY(2),'SIGNATURE',12,HSIGN)
+ IF(HSIGN.NE.'L_OPTIMIZE')THEN
+ CALL XABORT('DSPH: SIGNATURE OF '//HENTRY(2)//' IS '//HSIGN//
+ 1 '. L_OPTIMIZE EXPECTED.')
+ ENDIF
+ ELSE IF(JENTRY(2).EQ.2)THEN
+ CALL XABORT('DSPH: OPTIMIZE OBJECT IN CREATION OR MODIFICATION'
+ 1 //' MODE EXPECTED.')
+ ENDIF
+ IPOPT=KENTRY(2)
+ IF((IENTRY(3).NE.1).AND.(IENTRY(3).NE.2))CALL XABORT('DSPH: LCM '
+ 1 //'OBJECT EXPECTED AT RHS.')
+ IF(JENTRY(3).NE.2)CALL XABORT('DSPH: MACROLIB IN READ-ONLY MODE '
+ 1 //'EXPECTED AT RHS.')
+ CALL LCMGTC(KENTRY(3),'SIGNATURE',12,HSIGN)
+ IF(HSIGN.NE.'L_MACROLIB')THEN
+ CALL XABORT('DSPH: SIGNATURE OF '//HENTRY(3)//' IS '//HSIGN//
+ 1 '. L_MACROLIB EXPECTED.')
+ ENDIF
+ IPOLD=KENTRY(3)
+ CALL LCMGET(IPOLD,'STATE-VECTOR',ISTATE)
+ NGRP=ISTATE(1)
+ NMIX=ISTATE(2)
+ NL=ISTATE(3)
+ NIFISS=ISTATE(4)
+ NED=ISTATE(5)
+ NDEL=ISTATE(7)
+ NALBP=ISTATE(8)
+ ILEAKS=ISTATE(9)
+*----
+* READ THE INPUT DATA
+*----
+ IMPX=1
+ IMC=2
+ NGR1=1
+ NGR2=NGRP
+ NMIXP=NMIX
+ 20 CALL REDGET(INDIC,NITMA,FLOTT,TEXT12,DFLOTT)
+ IF(INDIC.EQ.10) GO TO 30
+ IF(INDIC.NE.3) CALL XABORT('DSPH: CHARACTER DATA EXPECTED(1).')
+ IF(TEXT12.EQ.'EDIT') THEN
+* READ THE PRINT INDEX.
+ CALL REDGET(INDIC,IMPX,FLOTT,TEXT12,DFLOTT)
+ IF(INDIC.NE.1) CALL XABORT('DSPH: INTEGER DATA EXPECTED(1).')
+ ELSE IF(TEXT12.EQ.'SPH') THEN
+* READ THE TYPE OF SPH CORRECTION.
+ CALL REDGET(INDIC,NITMA,FLOTT,TEXT12,DFLOTT)
+ IF(INDIC.NE.3) CALL XABORT('DSPH: CHARACTER DATA EXPECTED(2).')
+ IF(TEXT12.EQ.'PN') THEN
+ IMC=1
+ ELSE IF(TEXT12.EQ.'SN') THEN
+ IMC=2
+ ELSE IF(TEXT12.EQ.'ALBEDO') THEN
+ IMC=3
+ ELSE
+ CALL XABORT('DSPH: INVALID TYPE OF SPH CORRECTION.')
+ ENDIF
+ ELSE IF(TEXT12.EQ.'GRPMIN') THEN
+* READ THE MINIMUM GROUP INDEX.
+ CALL REDGET(INDIC,NGR1,FLOTT,TEXT12,DFLOTT)
+ IF(INDIC.NE.1) CALL XABORT('DSPH: INTEGER DATA EXPECTED(4).')
+ IF((NGR1.LE.0).OR.(NGR1.GT.NGRP)) CALL XABORT('DSPH: INVALID '
+ 1 //'VALUE OF GRPMIN.')
+ ELSE IF(TEXT12.EQ.'GRPMAX') THEN
+* READ THE MAXIMUM GROUP INDEX.
+ CALL REDGET(INDIC,NGR2,FLOTT,TEXT12,DFLOTT)
+ IF(INDIC.NE.1) CALL XABORT('DSPH: INTEGER DATA EXPECTED(5).')
+ IF((NGR2.LT.NGR1).OR.(NGR2.GT.NGRP)) CALL XABORT('DSPH: INVAL'
+ 1 //'ID VALUE OF GRPMAX.')
+ ELSE IF(TEXT12.EQ.';') THEN
+ GO TO 30
+ ELSE
+ CALL XABORT('DSPH: '//TEXT12//' IS AN INVALID KEYWORD.')
+ ENDIF
+ GO TO 20
+ 30 IF(NGR2.LT.NGR1) CALL XABORT('DSPH: INVALID GROUP INDICES.')
+ NMIXP=NMIX
+ IF(IMC.EQ.3) NMIXP=0
+ NPERT=(NMIXP+NALBP)*(NGR2-NGR1+1)
+ IF(IMPX.GT.0) WRITE(6,'(/36H DSPH: NUMBER OF CROSS-SECTION PERTU,
+ 1 9HRBATIONS=,I5)') NPERT
+*----
+* SET THE PERTURBED MACROLIB
+*----
+ JPNEW=LCMLID(IPNEW,'STEP',NPERT)
+ JPOLD=LCMGID(IPOLD,'GROUP')
+ IPERT=0
+ ALLOCATE(SPH(NMIXP+NALBP,NGRP),VARV(NPERT),ALBP(NALBP,NGRP),
+ 1 PALBP(NALBP,NGRP))
+ ALLOCATE(IHEDIT(2,NED+1),IJJ(NMIX),NJJ(NMIX),IPOS(NMIX))
+ ALLOCATE(DIFHOM(NGRP),GAR(NMIX),PER(NMIX),GAR1(NMIX*NGRP),
+ 1 PER1(NMIX*NGRP),GAR2(NMIX,NIFISS),PER2(NMIX,NIFISS),
+ 2 PSIGS(NMIX,NGRP,NL))
+*----
+* RECOVER SPH FACTORS
+*----
+ IF(NALBP.GT.0) CALL LCMGET(IPOLD,'ALBEDO',ALBP)
+ SPH(:NMIXP+NALBP,:NGRP)=1.0
+ DO 40 IGRP=NGR1,NGR2
+ KPOLD=LCMGIL(JPOLD,IGRP)
+ CALL LCMLEN(KPOLD,'NSPH',ILCMLN,ITYLCM)
+ IF(ILCMLN.EQ.NMIX) THEN
+ CALL LCMGET(KPOLD,'NSPH',SPH(1,IGRP))
+ IF(NALBP.GT.0) SPH(NMIXP+1:NMIXP+NALBP,IGRP)=1.0
+ ELSE
+ SPH(:NMIXP+NALBP,IGRP)=1.0
+ ENDIF
+ 40 CONTINUE
+*----
+* MACROSCOPIC TOTAL CROSS SECTIONS
+*----
+ DO 190 IGRP=NGR1,NGR2
+ DO 130 IBMP=1,NMIXP
+ PSIGS(:NMIX,:NGRP,:NL)=0.0
+ IPERT=IPERT+1
+ IF(IPERT.GT.NPERT) CALL XABORT('DSPH: NPERT OVERFLOW(1).')
+ VARV(IPERT)=SPH(IBMP,IGRP)
+ KPNEW=LCMDIL(JPNEW,IPERT)
+ IF(NALBP.GT.0) THEN
+ PALBP(:NALBP,:NGRP)=1.0
+ CALL LCMPUT(KPNEW,'ALBEDO',NALBP*NGRP,2,PALBP)
+ ENDIF
+ LPNEW=LCMLID(KPNEW,'GROUP',NGRP)
+ DO 110 IGR=1,NGRP
+ MPNEW=LCMDIL(LPNEW,IGR)
+ KPOLD=LCMGIL(JPOLD,IGR)
+ GAR(:NMIX)=0.0
+ NJJ(:NMIX)=1
+ IJJ(:NMIX)=IGR
+ CALL LCMPUT(MPNEW,'NTOT0',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'SIGS00',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'SIGW00',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'SCAT00',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'NJJS00',NMIX,1,NJJ)
+ CALL LCMPUT(MPNEW,'IJJS00',NMIX,1,IJJ)
+ CALL LCMPUT(MPNEW,'IPOS00',NMIX,1,NJJ)
+*----
+* MACROSCOPIC TOTAL CROSS SECTIONS
+*----
+ PER(:NMIX)=0.0
+ CALL LCMLEN(KPOLD,'NTOT0',ILCMLN,ITYLCM)
+ IF(ILCMLN.EQ.0) CALL XABORT('DSPH: MISSING NTOT0 INFO')
+ CALL LCMGET(KPOLD,'NTOT0',GAR)
+ IF(IMC.EQ.1) THEN
+ IF(IGR.EQ.IGRP) PER(IBMP)=GAR(IBMP)/SPH(IBMP,IGR)
+ CALL LCMPUT(MPNEW,'NTOT0',NMIX,2,PER)
+ ENDIF
+ PER(:NMIX)=0.0
+ CALL LCMLEN(KPOLD,'NTOT1',ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) CALL LCMGET(KPOLD,'NTOT1',GAR)
+ IF(IMC.EQ.1) THEN
+ IF(IGR.EQ.IGRP) PER(IBMP)=-GAR(IBMP)/SPH(IBMP,IGR)
+ CALL LCMPUT(MPNEW,'NTOT1',NMIX,2,PER)
+ ENDIF
+*----
+* MACROSCOPIC NU*FISSION CROSS SECTIONS (STEADY-STATE AND DELAYED)
+*----
+ IF(NIFISS.GT.0) THEN
+ PER2(:NMIX,:NIFISS)=0.0
+ CALL LCMGET(KPOLD,'NUSIGF',GAR2)
+ IF(IGR.EQ.IGRP) THEN
+ DO 50 IFIS=1,NIFISS
+ PER2(IBMP,IFIS)=GAR2(IBMP,IFIS)/SPH(IBMP,IGR)
+ 50 CONTINUE
+ ENDIF
+ CALL LCMPUT(MPNEW,'NUSIGF',NMIX*NIFISS,2,PER2)
+ DO 70 IDEL=1,NDEL
+ WRITE(TEXT12,'(6HNUSIGF,I2.2)') IDEL
+ PER2(:NMIX,:NIFISS)=0.0
+ CALL LCMGET(KPOLD,TEXT12,GAR2)
+ IF(IGR.EQ.IGRP) THEN
+ DO 60 IFIS=1,NIFISS
+ PER2(IBMP,IFIS)=GAR2(IBMP,IFIS)/SPH(IBMP,IGR)
+ 60 CONTINUE
+ ENDIF
+ CALL LCMPUT(MPNEW,TEXT12,NMIX*NIFISS,2,PER2)
+ 70 CONTINUE
+ ENDIF
+*----
+* MACROSCOPIC SCATTERING CROSS SECTIONS
+*----
+ DO 90 IL=1,NL
+ WRITE(TEXT2,'(I2.2)') IL-1
+ CALL LCMLEN(KPOLD,'NJJS'//TEXT2,ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) THEN
+ CALL LCMGET(KPOLD,'NJJS'//TEXT2,NJJ)
+ CALL LCMGET(KPOLD,'IJJS'//TEXT2,IJJ)
+ CALL LCMGET(KPOLD,'IPOS'//TEXT2,IPOS)
+ CALL LCMGET(KPOLD,'SCAT'//TEXT2,GAR1)
+ CALL LCMLEN(KPOLD,'NTOT1',ILCMLN,ITYLCM)
+ IF((MOD(IL-1,2).EQ.1).AND.(ILCMLN.GT.0)) THEN
+ CALL LCMGET(KPOLD,'NTOT1',GAR)
+ ELSE
+ CALL LCMGET(KPOLD,'NTOT0',GAR)
+ ENDIF
+ PER1(:NMIX*NGRP)=0.0
+ IPO=IPOS(IBMP)
+ DO 80 JGR=IJJ(IBMP),IJJ(IBMP)-NJJ(IBMP)+1,-1
+ IF(MOD(IL-1,2).EQ.0) THEN
+ IF((IGR.EQ.JGR).AND.(IMC.GT.1)) THEN
+ IF(IGR.EQ.IGRP) THEN
+ PER1(IPO)=GAR1(IPO)/SPH(IBMP,IGR)-GAR(IBMP)/SPH(IBMP,IGR)
+ ENDIF
+ ELSE
+ IF(JGR.EQ.IGRP) THEN
+ PER1(IPO)=GAR1(IPO)/SPH(IBMP,JGR) ! IGR <- JGR
+ ENDIF
+ ENDIF
+ ELSE
+ IF((IGR.EQ.JGR).AND.(IMC.GT.1)) THEN
+ IF(IGR.EQ.IGRP) THEN
+ PER1(IPO)=-GAR1(IPO)/SPH(IBMP,IGR)+GAR(IBMP)/SPH(IBMP,IGR)
+ ENDIF
+ ELSE
+ IF(IGR.EQ.IGRP) THEN
+ PER1(IPO)=-GAR1(IPO)/SPH(IBMP,IGR)
+ ENDIF
+ ENDIF
+ ENDIF
+ PSIGS(IBMP,IGR,IL)=PSIGS(IBMP,IGR,IL)+PER1(IPO)
+ IPO=IPO+1
+ 80 CONTINUE
+ CALL LCMPUT(MPNEW,'NJJS'//TEXT2,NMIX,1,NJJ)
+ CALL LCMPUT(MPNEW,'IJJS'//TEXT2,NMIX,1,IJJ)
+ CALL LCMPUT(MPNEW,'IPOS'//TEXT2,NMIX,1,IPOS)
+ CALL LCMPUT(MPNEW,'SCAT'//TEXT2,IPOS(NMIX)+NJJ(NMIX)-1,2,PER1)
+ ENDIF
+ CALL LCMLEN(KPOLD,'SIGW'//TEXT2,ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) THEN
+ PER(:NMIX)=0.0
+ CALL LCMGET(KPOLD,'SIGW'//TEXT2,GAR1)
+ CALL LCMLEN(KPOLD,'NTOT1',ILCMLN,ITYLCM)
+ IF((MOD(IL-1,2).EQ.1).AND.(ILCMLN.GT.0)) THEN
+ CALL LCMGET(KPOLD,'NTOT1',GAR)
+ ELSE
+ CALL LCMGET(KPOLD,'NTOT0',GAR)
+ ENDIF
+ IF(MOD(IL-1,2).EQ.0) THEN
+ IF(IMC.EQ.1) THEN
+ IF(IGR.EQ.IGRP) THEN
+ PER(IBMP)=GAR1(IBMP)/SPH(IBMP,IGR)
+ ENDIF
+ ELSE
+ IF(IGR.EQ.IGRP) THEN
+ PER(IBMP)=-GAR1(IBMP)/SPH(IBMP,IGR)+GAR(IBMP)/SPH(IBMP,IGR)
+ ENDIF
+ ENDIF
+ ELSE
+ IF(IMC.EQ.1) THEN
+ IF(IGR.EQ.IGRP) THEN
+ PER(IBMP)=-GAR1(IBMP)/SPH(IBMP,IGR)
+ ENDIF
+ ELSE
+ IF(IGR.EQ.IGRP) THEN
+ PER(IBMP)=-GAR1(IBMP)/SPH(IBMP,IGR)+GAR(IBMP)/SPH(IBMP,IGR)
+ ENDIF
+ ENDIF
+ ENDIF
+ CALL LCMPUT(MPNEW,'SIGW'//TEXT2,NMIX,2,PER)
+ ENDIF
+ 90 CONTINUE
+*----
+* DIFFUSION COEFFICIENTS
+*----
+ IF(ILEAKS.EQ.1) THEN
+ CALL LCMLEN(KPOLD,'DIFF',ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) THEN
+ PER(:NMIX)=0.0
+ CALL LCMGET(KPOLD,'DIFF',GAR)
+ IF(IGR.EQ.IGRP) PER(IBMP)=GAR(IBMP)/SPH(IBMP,IGR)
+ ELSE
+ PER(:NMIX)=0.0
+ CALL LCMGET(IPOLD,'DIFHOMB1HOM',DIFHOM)
+ IF(IGR.EQ.IGRP) PER(IBMP)=DIFHOM(IGR)/SPH(IBMP,IGR)
+ ENDIF
+ CALL LCMPUT(MPNEW,'DIFF',NMIX,2,PER)
+ ELSE IF(ILEAKS.EQ.2) THEN
+ CALL LCMLEN(KPOLD,'DIFFX',ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) THEN
+ PER(:NMIX)=0.0
+ CALL LCMGET(KPOLD,'DIFFX',GAR)
+ IF(IGR.EQ.IGRP) PER(IBMP)=GAR(IBMP)/SPH(IBMP,IGR)
+ CALL LCMPUT(MPNEW,'DIFFX',NMIX,2,PER)
+ ENDIF
+ CALL LCMLEN(KPOLD,'DIFFY',ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) THEN
+ PER(:NMIX)=0.0
+ CALL LCMGET(KPOLD,'DIFFY',GAR)
+ IF(IGR.EQ.IGRP) PER(IBMP)=GAR(IBMP)/SPH(IBMP,IGR)
+ CALL LCMPUT(MPNEW,'DIFFY',NMIX,2,PER)
+ ENDIF
+ CALL LCMLEN(KPOLD,'DIFFZ',ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) THEN
+ PER(:NMIX)=0.0
+ CALL LCMGET(KPOLD,'DIFFZ',GAR)
+ IF(IGR.EQ.IGRP) PER(IBMP)=GAR(IBMP)/SPH(IBMP,IGR)
+ CALL LCMPUT(MPNEW,'DIFFZ',NMIX,2,PER)
+ ENDIF
+ ENDIF
+*----
+* SPECIFIC REACTIONS
+*----
+ CALL LCMLEN(KPOLD,'TRANC',ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) THEN
+ PER(:NMIX)=0.0
+ CALL LCMGET(KPOLD,'TRANC',GAR)
+ IF(IGR.EQ.IGRP) PER(IBMP)=-GAR(IBMP)/SPH(IBMP,IGR)
+ CALL LCMPUT(MPNEW,'TRANC',NMIX,2,PER)
+ ENDIF
+*----
+* ADDITIONAL PHI-WEIGHTED EDITS
+*----
+ DO 100 IED=1,NED
+ WRITE(TEXT8,'(2A4)') (IHEDIT(I0,IED),I0=1,2)
+ IF(TEXT8(:5).EQ.'TRANC') GO TO 100
+ CALL LCMLEN(KPOLD,TEXT8,ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) THEN
+ PER(:NMIX)=0.0
+ CALL LCMGET(KPOLD,TEXT8,GAR)
+ IF(TEXT8(:4).EQ.'STRD') THEN
+ IF(IGR.EQ.IGRP) PER(IBMP)=GAR(IBMP)/SPH(IBMP,IGR)
+ ELSE
+ IF(IGR.EQ.IGRP) PER(IBMP)=GAR(IBMP)*SPH(IBMP,IGR)
+ ENDIF
+ CALL LCMPUT(MPNEW,TEXT8,NMIX,2,PER)
+ ENDIF
+ 100 CONTINUE
+ 110 CONTINUE
+*----
+* STORE SCATTERING CROSS SECTIONS
+*----
+ DO 125 IGR=1,NGRP
+ MPNEW=LCMDIL(LPNEW,IGR)
+ KPOLD=LCMGIL(JPOLD,IGR)
+ CALL LCMLEN(KPOLD,'SIGS00',ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) THEN
+ PSIGS(:NMIX,IGR,1)=0.0
+ CALL LCMGET(KPOLD,'SIGS00',GAR1)
+ CALL LCMGET(KPOLD,'NTOT0',GAR)
+ IF(IMC.EQ.1) THEN
+ IF(IGR.EQ.IGRP) PSIGS(IBMP,IGR,1)=GAR1(IBMP)/SPH(IBMP,IGR)
+ ELSE
+ IF(IGR.EQ.IGRP) PSIGS(IBMP,IGR,1)=GAR1(IBMP)/SPH(IBMP,IGR)-
+ > GAR(IBMP)/SPH(IBMP,IGR)
+ ENDIF
+ ENDIF
+ DO 120 IL=1,NL
+ WRITE(TEXT2,'(I2.2)') IL-1
+ CALL LCMLEN(KPOLD,'SIGS'//TEXT2,ILCMLN,ITYLCM)
+ IF(ILCMLN.GT.0) THEN
+ CALL LCMPUT(MPNEW,'SIGS'//TEXT2,NMIX,2,PSIGS(1,IGR,IL))
+ ENDIF
+ 120 CONTINUE
+ 125 CONTINUE
+ 130 CONTINUE
+*----
+* DERIVATIVE RELATIVE TO PHYSICAL ALBEDOS
+*----
+ DO 180 IALP=1,NALBP
+ IPERT=IPERT+1
+ IF(IPERT.GT.NPERT) CALL XABORT('DSPH: NPERT OVERFLOW(2).')
+ VARV(IPERT)=SPH(NMIXP+IALP,IGRP)
+ KPNEW=LCMDIL(JPNEW,IPERT)
+ PALBP(:NALBP,:NGRP)=1.0
+ FAT=0.5*(1.0-ALBP(IALP,IGRP))/(1.0+ALBP(IALP,IGRP))/
+ 1 REAL(VARV(IPERT))
+ PALBP(IALP,IGRP)=(1.0-2.0*FAT)/(1.0+2.0*FAT)
+ LPNEW=LCMLID(KPNEW,'GROUP',NGRP)
+ CALL LCMPUT(KPNEW,'ALBEDO',NALBP*NGRP,2,PALBP)
+ DO 170 IGR=1,NGRP
+ MPNEW=LCMDIL(LPNEW,IGR)
+ GAR(:NMIX)=0.0
+ NJJ(:NMIX)=1
+ DO 140 IMIX=1,NMIX
+ IJJ(IMIX)=IGR
+ 140 CONTINUE
+ CALL LCMPUT(MPNEW,'NTOT0',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'SIGS00',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'SIGW00',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'SCAT00',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'NJJS00',NMIX,1,NJJ)
+ CALL LCMPUT(MPNEW,'IJJS00',NMIX,1,IJJ)
+ CALL LCMPUT(MPNEW,'IPOS00',NMIX,1,NJJ)
+ CALL LCMPUT(MPNEW,'NTOT0',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'NTOT1',NMIX,2,GAR)
+ IF(NIFISS.GT.0) THEN
+ PER2(:NMIX,:NIFISS)=0.0
+ CALL LCMPUT(MPNEW,'NUSIGF',NMIX*NIFISS,2,PER2)
+ DO 150 IDEL=1,NDEL
+ WRITE(TEXT12,'(6HNUSIGF,I2.2)') IDEL
+ PER2(:NMIX,:NIFISS)=0.0
+ CALL LCMPUT(MPNEW,TEXT12,NMIX*NIFISS,2,PER2)
+ 150 CONTINUE
+ ENDIF
+ IF(ILEAKS.EQ.1) THEN
+ CALL LCMPUT(MPNEW,'DIFF',NMIX,2,GAR)
+ ELSE IF(ILEAKS.EQ.2) THEN
+ CALL LCMPUT(MPNEW,'DIFFX',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'DIFFY',NMIX,2,GAR)
+ CALL LCMPUT(MPNEW,'DIFFZ',NMIX,2,GAR)
+ ENDIF
+ DO 160 IED=1,NED
+ WRITE(TEXT8,'(2A4)') (IHEDIT(I0,IED),I0=1,2)
+ IF(TEXT8(:5).EQ.'TRANC') GO TO 160
+ CALL LCMPUT(MPNEW,TEXT8,NMIX,2,GAR)
+ 160 CONTINUE
+*----
+* END OF LOOP OVER PERTURBED MACROLIBS
+*----
+ 170 CONTINUE
+ 180 CONTINUE
+ 190 CONTINUE
+ DEALLOCATE(PSIGS,PER2,GAR2,PER1,GAR1,PER,GAR,DIFHOM)
+ DEALLOCATE(IPOS,NJJ,IJJ,IHEDIT)
+ DEALLOCATE(PALBP,ALBP,SPH)
+*----
+* SET THE PERTURBED MACROLIB STATE-VECTOR
+*----
+ ISTATE(:NSTATE)=0
+ ISTATE(1)=NGRP
+ ISTATE(2)=NMIX
+ ISTATE(3)=1
+ ISTATE(7)=NALBP
+ ISTATE(9)=ILEAKS
+ ISTATE(11)=NPERT
+ CALL LCMPUT(IPNEW,'STATE-VECTOR',NSTATE,1,ISTATE)
+ IF(IMPX.GT.1) CALL LCMLIB(IPNEW)
+*----
+* PUT OPTIMIZE OBJECT INFORMATION
+*----
+ CALL LCMPUT(IPOPT,'VAR-VALUE',NPERT,4,VARV)
+ DEALLOCATE(VARV)
+ IF(JENTRY(2).EQ.0)THEN
+ ISTATE(:NSTATE)=0
+ ISTATE(1)=NGRP
+ ISTATE(2)=NMIX
+ ISTATE(3)=1
+ ISTATE(4)=2+IMC
+ ISTATE(5)=NGR1
+ ISTATE(6)=NGR2
+ ISTATE(7)=1
+ ISTATE(8)=NMIX
+ ISTATE(9)=NALBP
+ IF(IMPX.GT.0) WRITE(6,200) (ISTATE(I),I=1,6)
+ CALL LCMPUT(IPOPT,'DEL-STATE',NSTATE,1,ISTATE)
+ ISTATE(:NSTATE)=0
+ ISTATE(1)=NPERT
+ ISTATE(2)=0
+ ISTATE(3)=1
+ ISTATE(4)=0
+ ISTATE(5)=0
+ ISTATE(6)=2
+ ISTATE(9)=2
+ ISTATE(10)=0
+ CALL LCMPUT(IPOPT,'STATE-VECTOR',NSTATE,1,ISTATE)
+ OPTPRR(:NSTATE)=0.0D0
+ OPTPRR(1)=1.0D0
+ OPTPRR(2)=0.1D0
+ OPTPRR(3)=1.0D-4
+ OPTPRR(4)=1.0D-4
+ OPTPRR(5)=1.0D-4
+ CALL LCMPUT(IPOPT,'OPT-PARAM-R',NSTATE,4,OPTPRR)
+ ENDIF
+ RETURN
+*
+ 200 FORMAT(/18H DEL-STATE OPTIONS/18H -----------------/
+ 1 7H NGRP ,I8,28H (NUMBER OF ENERGY GROUPS)/
+ 2 7H NMIX ,I8,32H (NUMBER OF MATERIAL MIXTURES)/
+ 3 7H ITYPE ,I8,13H (NOT USED)/
+ 4 7H IDELTA,I8,43H (=3/4/5: USE PN-TYPE/USE SN-TYPE/ALBEDO)/
+ 5 7H NGR1 ,I8,24H (MINIMUM GROUP INDEX)/
+ 6 7H NGR2 ,I8,24H (MAXIMUM GROUP INDEX))
+ END
|