summaryrefslogtreecommitdiff
path: root/Donjon
diff options
context:
space:
mode:
Diffstat (limited to 'Donjon')
-rw-r--r--Donjon/src/DONDRV.F6
-rw-r--r--Donjon/src/FPSOUT.f150
-rw-r--r--Donjon/src/FPSPH.f472
3 files changed, 0 insertions, 628 deletions
diff --git a/Donjon/src/DONDRV.F b/Donjon/src/DONDRV.F
index ee42c61..1b269d7 100644
--- a/Donjon/src/DONDRV.F
+++ b/Donjon/src/DONDRV.F
@@ -298,12 +298,6 @@
NAM='A. HEBERT'
WRITE(IOUT,1000)HMODUL,COD,DSR,NAM
CALL IDET(NENTRY,HENTRY,IENTRY,JENTRY,KENTRY)
- ELSE IF(HMODUL.EQ.'FPSPH:') THEN
- COD='DONJON'
- DSR='SINGLE SPH FACTOR FIXED POINT ITERATION'
- NAM='A. HEBERT'
- WRITE(IOUT,1000)HMODUL,COD,DSR,NAM
- CALL FPSPH(NENTRY,HENTRY,IENTRY,JENTRY,KENTRY)
ELSE
DONMOD=.FALSE.
DONDRV=KDRDRV(HMODUL,NENTRY,HENTRY,IENTRY,JENTRY,KENTRY)
diff --git a/Donjon/src/FPSOUT.f b/Donjon/src/FPSOUT.f
deleted file mode 100644
index 5d78203..0000000
--- a/Donjon/src/FPSOUT.f
+++ /dev/null
@@ -1,150 +0,0 @@
-*DECK FPSOUT
- SUBROUTINE FPSOUT(IPMAC,IPRINT,NG,NMIL,NFIS,ILEAKS,TEXT9,OUTG)
-*
-*-----------------------------------------------------------------------
-*
-*Purpose:
-* Compute the leakage rate in each energy group
-*
-*Copyright:
-* Copyright (C) 2019 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
-* IPMAC pointer to the macrolib structure.
-* IPRINT print parameter
-* NG number of energy groups.
-* NMIL number of material mixtures.
-* NFIS number of fissile isotopes.
-* ILEAKS type of leakage calculation =0: no leakage; =1: homogeneous
-* leakage (Diffon).
-* TEXT9 type of calculation ('REFERENCE' or 'MACRO').
-*
-*Parameters: output
-* OUTG leakage rates.
-*
-*-----------------------------------------------------------------------
-*
- USE GANLIB
-*----
-* SUBROUTINE ARGUMENTS
-*----
- TYPE(C_PTR) IPMAC
- INTEGER IPRINT,NG,NMIL,NFIS,ILEAKS
- CHARACTER TEXT9*9
- REAL OUTG(NG)
-*----
-* LOCAL VARIABLES
-*----
- TYPE(C_PTR) JPMAC,KPMAC
- CHARACTER HSMG*131
-*----
-* ALLOCATABLE ARRAYS
-*----
- INTEGER, ALLOCATABLE, DIMENSION(:) :: IJJ,NJJ,IPOS
- REAL, ALLOCATABLE, DIMENSION(:) :: GAR,WORK,DIFHOM,DIFF
- REAL, ALLOCATABLE, DIMENSION(:,:) :: PHI,NUF
- REAL, ALLOCATABLE, DIMENSION(:,:,:) :: CHI,RHS,LHS
-*----
-* SCRATCH STORAGE ALLOCATION
-*----
- ALLOCATE(PHI(NMIL,NG),RHS(NMIL,NG,NG),LHS(NMIL,NG,NG))
- ALLOCATE(IJJ(NMIL),NJJ(NMIL),IPOS(NMIL),GAR(NMIL),WORK(NMIL*NG),
- > CHI(NMIL,NFIS,NG),NUF(NMIL,NFIS),DIFHOM(NG),DIFF(NMIL))
-*----
-* COMPUTE THE ACTUAL AND REFERENCE REACTION RATE MATRICES
-*----
- CALL LCMGET(IPMAC,'K-EFFECTIVE',ZKEFF)
- IF(IPRINT.GT.1) WRITE(6,120) TEXT9,ZKEFF
- CALL LCMLEN(IPMAC,'B2 B1HOM',ILCMLN,ITYLCM)
- IF(ILCMLN.EQ.1) THEN
- CALL LCMGET(IPMAC,'B2 B1HOM',B2)
- ELSE
- B2=0.0
- ENDIF
- IF((ILEAKS.EQ.1).AND.(IPRINT.GT.1)) THEN
- WRITE(6,'(/9H FPSOUT: ,A,4H B2=,1P,E12.4)') TEXT9,B2
- ENDIF
- RHS(:NMIL,:NG,:NG)=0.0
- LHS(:NMIL,:NG,:NG)=0.0
- JPMAC=LCMGID(IPMAC,'GROUP')
- DO IG=1,NG
- KPMAC=LCMGIL(JPMAC,IG)
- CALL LCMGET(KPMAC,'CHI',CHI(1,1,IG))
- CALL LCMLEN(KPMAC,'FLUX-INTG',ILG,ITYLCM)
- IF(ILG.NE.NMIL) CALL XABORT('FPSOUT: MISSING REFERENCE FLUX.')
- CALL LCMGET(KPMAC,'FLUX-INTG',PHI(1,IG))
- ENDDO
- DO IG=1,NG
- KPMAC=LCMGIL(JPMAC,IG)
- IF(ILEAKS.EQ.1) THEN
- CALL LCMLEN(KPMAC,'DIFF',ILCMLN,ITYLCM)
- IF(ILCMLN.GT.0) THEN
- CALL LCMGET(KPMAC,'DIFF',DIFF)
- ELSE
- CALL LCMGET(IPMAC,'DIFHOMB1HOM',DIFHOM)
- DO IBM=1,NMIL
- DIFF(IBM)=DIFHOM(IG)
- ENDDO
- ENDIF
- ELSE
- DIFF(:NMIL)=0.0
- ENDIF
- CALL LCMGET(KPMAC,'NTOT0',GAR)
- CALL LCMGET(KPMAC,'SCAT00',WORK)
- CALL LCMGET(KPMAC,'NJJS00',NJJ)
- CALL LCMGET(KPMAC,'IJJS00',IJJ)
- CALL LCMGET(KPMAC,'IPOS00',IPOS)
- DO IBM=1,NMIL
- IPOSDE=IPOS(IBM)
- DO JG=IJJ(IBM),IJJ(IBM)-NJJ(IBM)+1,-1
-* IG <-- JG
- RHS(IBM,IG,JG)=RHS(IBM,IG,JG)-WORK(IPOSDE)*PHI(IBM,JG)
- IPOSDE=IPOSDE+1
- ENDDO
- RHS(IBM,IG,IG)=RHS(IBM,IG,IG)+(GAR(IBM)+B2*DIFF(IBM))*
- > PHI(IBM,IG)
- ENDDO
- CALL LCMGET(KPMAC,'NUSIGF',NUF)
- DO IBM=1,NMIL
- DO IFIS=1,NFIS
- DO JG=1,NG
- LHS(IBM,JG,IG)=LHS(IBM,JG,IG)+CHI(IBM,IFIS,JG)*
- > NUF(IBM,IFIS)*PHI(IBM,IG)
- ENDDO
- ENDDO
- ENDDO
- ENDDO
-*----
-* COMPUTE THE ACTUAL AND REFERENCE ABSORPTION AND FISSION RATES
-*----
- DO IG=1,NG
- OUTG(IG)=0.0
- DO IBM=1,NMIL
- OUTG(IG)=OUTG(IG)+SUM(LHS(IBM,IG,:NG))/ZKEFF-
- 1 SUM(RHS(IBM,IG,:NG))
- ENDDO
- IF(OUTG(IG).LT.-1.0E-6) THEN
- WRITE(HSMG,'(21HFPSOUT: INCONSISTENT ,A,17H LEAKAGE IN GROUP,
- 1 I4,7H. LEAK=,1P,E13.4)') TEXT9,IG,OUTG(IG)
- CALL XABORT(HSMG)
- ENDIF
- IF(IPRINT.GT.1) WRITE(6,130) IG,TEXT9,OUTG(IG)
- ENDDO
-*----
-* SCRATCH STORAGE DEALLOCATION
-*----
- DEALLOCATE(DIFF,DIFHOM,NUF,CHI,WORK,GAR,IPOS,NJJ,IJJ)
- DEALLOCATE(LHS,RHS,PHI)
- RETURN
-*
- 120 FORMAT(/9H FPSOUT: ,A,33H EFFECTIVE MULTIPLICATION FACTOR=,1P,
- 1 E12.4)
- 130 FORMAT(/8H FPSOUT:,5X,6HGROUP=,I4,1X,A,9H LEAKAGE=,1P,E12.4)
- END
diff --git a/Donjon/src/FPSPH.f b/Donjon/src/FPSPH.f
deleted file mode 100644
index 96db43a..0000000
--- a/Donjon/src/FPSPH.f
+++ /dev/null
@@ -1,472 +0,0 @@
-*DECK FPSPH
- SUBROUTINE FPSPH(NENTRY,HENTRY,IENTRY,JENTRY,KENTRY)
-*
-*-----------------------------------------------------------------------
-*
-*Purpose:
-* Perform a single SPH factor fixed point iteration
-*
-*Copyright:
-* Copyright (C) 2019 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
-* 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.
-*
-*Comments:
-* The FPSPH: calling specifications are:
-* OPTIM := FPSPH: [ OPTIM ] MACROLIB MACROREF :: (fpsph\_data) ;
-* where
-* OPTIM : name of the \emph{optimize} object (L\_OPTIMIZE signature)
-* containing the SPH factors. At the first call, object OPTIM must appear on
-* LHS to receive its initial values. On subsequent calls, object OPTIM must
-* appear on both LHS and RHS to be able to update the previous values.
-* MACROLIB : name of the read-only extended \emph{macrolib} object
-* (L\_MACROLIB signature) containing the macroscopic cross sections used by
-* the macro-calculation and fluxes produced by the macro-calculation.
-* MACROREF : name of the read-only extended \emph{macrolib} object
-* (L\_MACROLIB signature) containing the reference macroscopic cross
-* sections and fluxes.
-* (fpsph\_data) : structure containing the data to the module FPSPH:
-*
-*-----------------------------------------------------------------------
-*
- USE GANLIB
-*----
-* SUBROUTINE ARGUMENTS
-*----
- INTEGER NENTRY,IENTRY(NENTRY),JENTRY(NENTRY)
- TYPE(C_PTR) KENTRY(NENTRY)
- CHARACTER HENTRY(NENTRY)*12
-*----
-* LOCAL VARIABLES
-*----
- PARAMETER (NSTATE=40)
- TYPE(C_PTR) IPOPT,IPMAC1,IPMAC2,JPMAC1,JPMAC2,KPMAC1,KPMAC2
- CHARACTER HSIGN*12,TEXT12*12
- INTEGER ISTATE(NSTATE),DNVTST
- DOUBLE PRECISION OPTPRR(NSTATE),DFLOTT,ZNORM1,ZNORM2,EPSPH,ERRT,
- > ERR2,ERROR,SPHMIN,SPHMAX
-*----
-* ALLOCATABLE ARRAYS
-*----
- REAL, ALLOCATABLE, DIMENSION(:) :: SPH,FLUX1,FLUX2,OUTG1,OUTG2
- DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: VARV,VAROLD,XMIN,
- > XMAX,P,FF,UD
- DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: DFF,TDFF
-*----
-* PARAMETER VALIDATION.
-*----
- IF(NENTRY.NE.3) CALL XABORT('FPSPH: THREE PARAMETERS EXPECTED.')
- IF((IENTRY(1).NE.1).AND.(IENTRY(1).NE.2)) CALL XABORT('FPSPH: LCM'
- > //' OBJECT EXPECTED AT LHS.')
- IF(JENTRY(1).EQ.0)THEN
- HSIGN='L_OPTIMIZE'
- CALL LCMPTC(KENTRY(1),'SIGNATURE',12,HSIGN)
- ELSE IF(JENTRY(1).EQ.1)THEN
- CALL LCMGTC(KENTRY(1),'SIGNATURE',12,HSIGN)
- IF(HSIGN.NE.'L_OPTIMIZE')THEN
- CALL XABORT('FPSPH: SIGNATURE OF '//HENTRY(2)//' IS '//HSIGN//
- > '. L_OPTIMIZE EXPECTED.')
- ENDIF
- ELSE IF(JENTRY(1).EQ.2)THEN
- CALL XABORT('FPSPH: OPTIMIZE OBJECT IN CREATION OR MODIFICATIO'
- > //'N MODE EXPECTED.')
- ENDIF
- IPOPT=KENTRY(1)
- IF(JENTRY(1).EQ.1) THEN
- CALL LCMGET(IPOPT,'STATE-VECTOR',ISTATE)
- NVAR=ISTATE(1)
- NFUNC=ISTATE(2)+1
- ITER=ISTATE(5)
- IMETH=ISTATE(8)
- CALL LCMGET(IPOPT,'OPT-PARAM-R',OPTPRR)
- EPSPH=OPTPRR(3)
- CALL LCMGET(IPOPT,'DEL-STATE',ISTATE)
- NGRP=ISTATE(1)
- NMIX=ISTATE(2)
- ICONT=ISTATE(4)
- NGR1=ISTATE(5)
- NGR2=ISTATE(6)
- NALBP=ISTATE(9)
- IF((ICONT.NE.3).AND.(ICONT.NE.4)) CALL XABORT('FPSPH: SPH FACT'
- > //'ORS EXPECTED IN OPTIMIZE OBJECT.')
- IF(NVAR.NE.(NGR2-NGR1+1)*(NMIX+NALBP)) CALL XABORT('FPSPH: INC'
- > //'OHERENT NUMBER OF DECISION VARIABLES.')
- ELSE
- ITER=0
- IMETH=3
- EPSPH=1.0D-4
- NGRP=0
- NMIX=0
- ENDIF
- DO I=2,3
- IF((JENTRY(I).NE.2).OR.((IENTRY(I).NE.1).AND.(IENTRY(I).NE.2)))
- 1 CALL XABORT('FPSPH: LCM OBJECT IN READ-ONLY MODE EXPECTED AT R'
- 2 //'HS.')
- ENDDO
- ITER=ITER+1
-*----
-* RECOVER THE ACTUAL MACROLIB.
-*----
- CALL LCMGTC(KENTRY(2),'SIGNATURE',12,HSIGN)
- IF(HSIGN.EQ.'L_MACROLIB') THEN
- IPMAC1=KENTRY(2)
- ELSE IF(HSIGN.EQ.'L_LIBRARY') THEN
- IPMAC1=LCMGID(KENTRY(5),'MACROLIB')
- ELSE
- TEXT12=HENTRY(2)
- CALL XABORT('FPSPH: SIGNATURE OF '//TEXT12//' IS '//HSIGN//
- > '. ACTUAL L_MACROLIB OR L_LIBRARY EXPECTED.')
- ENDIF
- CALL LCMGET(IPMAC1,'STATE-VECTOR',ISTATE)
- IF(JENTRY(1).EQ.0) THEN
- NGRP=ISTATE(1)
- NMIX=ISTATE(2)
- ELSE IF(ISTATE(1).NE.NGRP) THEN
- CALL XABORT('FPSPH: INVALID NUMBER OF GROUPS.')
- ELSE IF(ISTATE(2).NE.NMIX) THEN
- CALL XABORT('FPSPH: INVALID NUMBER OF MIXTURES.')
- ENDIF
- NFIS1=ISTATE(4)
- ILEAKS=ISTATE(9)
-*----
-* RECOVER THE REFERENCE MACROLIB.
-*----
- CALL LCMGTC(KENTRY(3),'SIGNATURE',12,HSIGN)
- IF(HSIGN.EQ.'L_MACROLIB') THEN
- IPMAC2=KENTRY(3)
- ELSE IF(HSIGN.EQ.'L_LIBRARY') THEN
- IPMAC2=LCMGID(KENTRY(3),'MACROLIB')
- ELSE
- TEXT12=HENTRY(3)
- CALL XABORT('FPSPH: SIGNATURE OF '//TEXT12//' IS '//HSIGN//
- 1 '. REFERENCE L_MACROLIB OR L_LIBRARY EXPECTED.')
- ENDIF
- CALL LCMGET(IPMAC2,'STATE-VECTOR',ISTATE)
- IF(ISTATE(1).NE.NGRP) THEN
- CALL XABORT('FPSPH: INVALID NUMBER OF REFERENCE GROUPS.')
- ELSE IF(ISTATE(2).NE.NMIX) THEN
- CALL XABORT('FPSPH: INVALID NUMBER OF REFERENCE MIXTURES.')
- ELSE IF(ISTATE(9).NE.ILEAKS) THEN
- CALL XABORT('FPSPH: INVALID TYPE OF LEAKAGE.')
- ENDIF
- NFIS2=ISTATE(4)
- NALBP=ISTATE(8)
- IF(NALBP.GT.1) CALL XABORT('FPSPH: NALBP>1 NOT SUPPORTED.')
-*----
-* READ INPUT PARAMETERS
-*----
- IPICK=0
- IPRINT=1
- SPHMIN=0.0D0
- SPHMAX=10.0D0
- IF(JENTRY(1).EQ.0) THEN
- IMC=2
- NGR1=1
- NGR2=NGRP
- ENDIF
- 10 CALL REDGET(INDIC,NITMA,FLOTT,TEXT12,DFLOTT)
- IF(INDIC.EQ.10) GO TO 50
- IF(INDIC.NE.3) CALL XABORT('FPSPH: CHARACTER DATA EXPECTED')
- IF(TEXT12(1:4).EQ.'EDIT') THEN
- CALL REDGET(INDIC,IPRINT,FLOTT,TEXT12,DFLOTT)
- IF(INDIC.NE.1) CALL XABORT('FPSPH: INTEGER DATA EXPECTED FOR I'
- 1 //'PRINT')
- 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('FPSPH: CHARACTER DATA EXPECTED(2).')
- IF(TEXT12.EQ.'PN') THEN
- IMC=1
- ELSE IF(TEXT12.EQ.'SN') THEN
- IMC=2
- ELSE
- CALL XABORT('FPSPH: 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('FPSPH: INTEGER DATA EXPECTED(4).')
- IF((NGR1.LE.0).OR.(NGR1.GT.NGRP)) CALL XABORT('FPSPH: INVALID '
- > //'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('FPSPH: INTEGER DATA EXPECTED(5).')
- IF((NGR2.LT.NGR1).OR.(NGR2.GT.NGRP)) CALL XABORT('FPSPH: INVAL'
- > //'ID VALUE OF GRPMAX.')
- ELSE IF(TEXT12.EQ.'OUT-STEP-EPS') THEN
-* Set the tolerence used for SPH iterations.
- CALL REDGET(INDIC,NITMA,FLOTT,TEXT12,DFLOTT)
- IF(INDIC.EQ.2) THEN
- EPSPH=FLOTT
- ELSE IF(INDIC.EQ.4) THEN
- EPSPH=DFLOTT
- ELSE
- CALL XABORT('FPSPH: REAL OR DOUBLE PRECISION VALUE EXPECTED.')
- ENDIF
- ELSE IF(TEXT12.EQ.'VAR-VAL-MIN') THEN
-* Set the minimum value for SPH dactors.
- CALL REDGET(INDIC,NITMA,FLOTT,TEXT12,DFLOTT)
- IF(INDIC.EQ.2) THEN
- SPHMIN=FLOTT
- ELSE IF(INDIC.EQ.4) THEN
- SPHMIN=DFLOTT
- ELSE
- CALL XABORT('FPSPH: REAL OR DOUBLE PRECISION VALUE EXPECTED.')
- ENDIF
- ELSE IF(TEXT12.EQ.'VAR-VAL-MAX') THEN
-* Set the maximum value for SPH dactors.
- CALL REDGET(INDIC,NITMA,FLOTT,TEXT12,DFLOTT)
- IF(INDIC.EQ.2) THEN
- SPHMAX=FLOTT
- ELSE IF(INDIC.EQ.4) THEN
- SPHMAX=DFLOTT
- ELSE
- CALL XABORT('FPSPH: REAL OR DOUBLE PRECISION VALUE EXPECTED.')
- ENDIF
- ELSE IF(TEXT12.EQ.'OUT-CONV-TST') THEN
-* Convergence test
- IPICK=1
- GO TO 50
- ELSE IF(TEXT12(1:1).EQ.';') THEN
- GO TO 50
- ELSE
- CALL XABORT('FPSPH: '//TEXT12//' IS AN INVALID KEYWORD')
- ENDIF
- GO TO 10
-*----
-* RECOVER SPH FACTORS FROM PREVIOUS ITERATION
-*----
- 50 NPERT=(NGR2-NGR1+1)*(NMIX+NALBP)
- ALLOCATE(VARV(NPERT),VAROLD(NPERT),XMIN(NPERT),XMAX(NPERT))
- CALL LCMLEN(IPOPT,'VAR-VAL-MIN',ILONG,ITYLCM)
- IF(ILONG.EQ.0) THEN
- XMIN(:NPERT)=SPHMIN
- CALL LCMPUT(IPOPT,'VAR-VAL-MIN',NPERT,4,XMIN)
- ELSE
- CALL LCMGET(IPOPT,'VAR-VAL-MIN',XMIN)
- ENDIF
- CALL LCMLEN(IPOPT,'VAR-VAL-MAX',ILONG,ITYLCM)
- IF(ILONG.EQ.0) THEN
- XMAX(:NPERT)=SPHMAX
- CALL LCMPUT(IPOPT,'VAR-VAL-MAX',NPERT,4,XMAX)
- ELSE
- CALL LCMGET(IPOPT,'VAR-VAL-MAX',XMAX)
- ENDIF
- CALL LCMLEN(IPOPT,'VAR-VALUE',ILONG,ITYLCM)
- IF(ILONG.EQ.0) THEN
- VAROLD(:NPERT)=1.0D0
- ELSE
- CALL LCMGET(IPOPT,'VAR-VALUE',VAROLD)
- ENDIF
-*----
-* PERFORM A FIXED POINT SPH ITERATION
-*----
- IF(IPRINT.GT.0) WRITE(6,'(/34H FPSPH: COMPUTE SPH FACTORS AT ITE,
- > 6HRATION,I5,12H WITH METHOD,I2,1H.)') ITER,IMETH
- IF(IMETH.EQ.3) THEN
- IPERT=0
- JPMAC1=LCMGID(IPMAC1,'GROUP')
- JPMAC2=LCMGID(IPMAC2,'GROUP')
- ALLOCATE(SPH(NMIX+NALBP),FLUX1(NMIX),FLUX2(NMIX),OUTG1(NGRP),
- > OUTG2(NGRP))
- IF(IPRINT.GT.4) WRITE(6,'(/32H FPSPH: SPH FACTORS AT ITERATION,
- > I5)') ITER
- IF(NALBP.GT.0) THEN
- CALL FPSOUT(IPMAC1,IPRINT,NGRP,NMIX,NFIS1,ILEAKS,' MACRO',
- > OUTG1)
- CALL FPSOUT(IPMAC2,IPRINT,NGRP,NMIX,NFIS2,ILEAKS,'REFERENCE',
- > OUTG2)
- ENDIF
- DO 120 IGR=NGR1,NGR2
- SPH(:NMIX+NALBP)=1.0
- KPMAC1=LCMGIL(JPMAC1,IGR)
- KPMAC2=LCMGIL(JPMAC2,IGR)
- CALL LCMGET(KPMAC1,'FLUX-INTG',FLUX1)
- CALL LCMGET(KPMAC2,'FLUX-INTG',FLUX2)
- DO 60 IBM=1,NMIX
- SPH(IBM)=FLUX2(IBM)/FLUX1(IBM)
- 60 CONTINUE
- DO 70 IAL=1,NALBP
- IF(OUTG1(IGR).NE.0.0) THEN
- SPH(NMIX+IAL)=REAL(VAROLD(IPERT+NMIX+1))*OUTG2(IGR)/OUTG1(IGR)
- ENDIF
- 70 CONTINUE
- ZNORM1=0.0D0
- ZNORM2=0.0D0
- DO 80 IBM=1,NMIX
- ZNORM1=ZNORM1+FLUX2(IBM)/SPH(IBM)
- ZNORM2=ZNORM2+FLUX2(IBM)
- 80 CONTINUE
- ZNORM1=ZNORM1/ZNORM2
- IF(IPRINT.GT.1) THEN
- WRITE(6,'(/14H FPSPH: GROUP=,I4,22H NORMALIZATION FACTOR=,1P,
- > E12.4)') IGR,ZNORM1
- ENDIF
- DO 90 IBM=1,NMIX+NALBP
- SPH(IBM)=SPH(IBM)*REAL(ZNORM1)
- 90 CONTINUE
- DO 100 IBM=1,NMIX
- IPERT=IPERT+1
- VARV(IPERT)=SPH(IBM)
- 100 CONTINUE
- DO 110 IAL=1,NALBP
- IPERT=IPERT+1
- VARV(IPERT)=SPH(NMIX+IAL)
- 110 CONTINUE
- 120 CONTINUE
- DEALLOCATE(OUTG2,OUTG1,FLUX2,FLUX1,SPH)
-*----
-* PERFORM A NEWTONIAN SPH ITERATION
-*----
- ELSE IF(IMETH.EQ.4) THEN
- ALLOCATE(P(NPERT),FF(NFUNC),DFF(NPERT,NFUNC),TDFF(NFUNC,NPERT),
- > UD(NPERT))
- CALL LCMGET(IPOPT,'FOBJ-CST-VAL',FF)
- CALL LCMGET(IPOPT,'GRADIENT',DFF)
- TDFF=TRANSPOSE(DFF)
- CALL ALST2F(NFUNC,NFUNC,NPERT,TDFF,UD)
- CALL ALST2S(NFUNC,NFUNC,NPERT,TDFF,UD,FF,P)
- DO 130 IPERT=1,NPERT
- VARV(IPERT)=VAROLD(IPERT)-P(IPERT)
- 130 CONTINUE
- DEALLOCATE(UD,TDFF,DFF,FF,P)
- ENDIF
-*----
-* APPLY CONSTRAINTS ON SPH FACTORS
-*----
- DO 135 IPERT=1,NPERT
- VARV(IPERT)=MAX(VARV(IPERT),XMIN(IPERT))
- VARV(IPERT)=MIN(VARV(IPERT),XMAX(IPERT))
- 135 ENDDO
-*----
-* PRINT SPH FACTORS
-*----
- IF(IPRINT.GT.4) THEN
- ALLOCATE(SPH(NMIX+NALBP))
- IPERT=0
- DO 150 IGR=NGR1,NGR2
- DO 140 IBM=1,NMIX+NALBP
- IPERT=IPERT+1
- SPH(IBM)=REAL(VARV(IPERT))
- 140 CONTINUE
- WRITE(6,200) 'NSPH',IGR,(SPH(IBM),IBM=1,NMIX+NALBP)
- 150 CONTINUE
- DEALLOCATE(SPH)
- ENDIF
-*----
-* TEST CONVERGENCE
-*----
- ICONV=0
- IF(JENTRY(1).EQ.1) THEN
- ERROR=0.0
- ERR2=0.0
- DO 160 IPERT=1,NPERT
- ERRT=ABS((VARV(IPERT)-VAROLD(IPERT))/VARV(IPERT))
- ERR2=ERR2+ERRT*ERRT
- ERROR=MAX(ERROR,ERRT)
- 160 CONTINUE
- ERR2=SQRT(ERR2/REAL(NPERT))
- IF(IPRINT.GT.0) WRITE(6,230) ITER,ERROR,ERR2
- IF(ERR2.LT.EPSPH) THEN
- ICONV=1
- IF(IPRINT.GT.0) WRITE(6,220) ITER
- ENDIF
- ELSE
- ERR2=1.0E10
- ENDIF
-*----
-* PUT OPTIMIZE OBJECT INFORMATION
-*----
- CALL LCMPUT(IPOPT,'VAR-VALUE',NPERT,4,VARV)
- DEALLOCATE(XMAX,XMIN,VAROLD,VARV)
- IF(JENTRY(1).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(IPRINT.GT.0) WRITE(6,210) (ISTATE(I),I=1,6)
- CALL LCMPUT(IPOPT,'DEL-STATE',NSTATE,1,ISTATE)
- ISTATE(:NSTATE)=0
- ISTATE(1)=NPERT
- ISTATE(8)=IMETH ! set to fixed point or Newtonian method
- CALL LCMPUT(IPOPT,'STATE-VECTOR',NSTATE,1,ISTATE)
- OPTPRR(:NSTATE)=0.0D0
- OPTPRR(1)=1.0D0
- OPTPRR(2)=0.1D0
- OPTPRR(3)=EPSPH
- OPTPRR(4)=1.0D-4
- OPTPRR(5)=1.0D-4
- CALL LCMPUT(IPOPT,'OPT-PARAM-R',NSTATE,4,OPTPRR)
- ELSE
- CALL LCMGET(IPOPT,'STATE-VECTOR',ISTATE)
- ISTATE(1)=NPERT
- ISTATE(4)=ICONV ! convergence index
- ISTATE(5)=ITER ! number of iterations
- ISTATE(8)=IMETH ! set to fixed point or Newtonian method
- CALL LCMPUT(IPOPT,'STATE-VECTOR',NSTATE,1,ISTATE)
- ENDIF
-*----
-* RECOVER THE CONVERGENCE FLAGS AND SAVE IT IN A CLE-2000 VARIABLE
-*----
- IF(IPICK.EQ.1) THEN
- CALL REDGET(INDIC,NITMA,FLOTT,TEXT12,DFLOTT)
- IF(INDIC.NE.-5) CALL XABORT('FPSPH: OUTPUT LOGICAL EXPECTED.')
- INDIC=5
- IF(ICONV.EQ.0) THEN
- DNVTST=-1 ! not converged
- ELSE IF(ICONV.EQ.1) THEN
- DNVTST=1 ! converged
- ENDIF
- CALL REDPUT(INDIC,DNVTST,FLOTT,TEXT12,DFLOTT)
- CALL REDGET(INDIC,NITMA,FLOTT,TEXT12,DFLOTT)
- IF(INDIC.EQ.-4) THEN
- INDIC=4
- CALL REDPUT(INDIC,NITMA,FLOTT,TEXT12,ERR2)
- CALL REDGET(INDIC,NITMA,FLOTT,TEXT12,DFLOTT)
- ENDIF
- IF((INDIC.NE.3).OR.(TEXT12.NE.';')) THEN
- CALL XABORT('FPSPH: ; CHARACTER EXPECTED.')
- ENDIF
- ENDIF
- RETURN
-*
- 200 FORMAT(/25H FPSPH: VALUES OF VECTOR ,A,9H IN GROUP,I5,4H ARE/
- > (1X,1P,10E13.5))
- 210 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,34H (=3/4: USE PN-TYPE/USE SN-TYPE)/
- 5 7H NGR1 ,I8,24H (MINIMUM GROUP INDEX)/
- 6 7H NGR2 ,I8,24H (MAXIMUM GROUP INDEX))
- 220 FORMAT(/39H FPSPH: CONVERGENCE OF SPH ALGORITHM IN,I5,
- > 12H ITERATIONS.)
- 230 FORMAT(/13H FPSPH: ITER=,I3,4X,6HERROR=,1P,E10.3,1X,5HERR2=,
- > E10.3)
- END