diff options
Diffstat (limited to 'Dragon/src/LIBCOR.f')
| -rw-r--r-- | Dragon/src/LIBCOR.f | 177 |
1 files changed, 177 insertions, 0 deletions
diff --git a/Dragon/src/LIBCOR.f b/Dragon/src/LIBCOR.f new file mode 100644 index 0000000..2fe872d --- /dev/null +++ b/Dragon/src/LIBCOR.f @@ -0,0 +1,177 @@ +*DECK LIBCOR + SUBROUTINE LIBCOR (IPLIB,NGRO,ISOT,JSOT,HNAMIS1,HNAMIS2) +* +*----------------------------------------------------------------------- +* +*Purpose: +* Compute the correlation information between a pair of resonant +* isotopes for the CALENDF method. +* +*Copyright: +* Copyright (C) 2003 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 +* IPLIB pointer to the internal library (L_LIBRARY signature). +* NGRO number of energy groups. +* ISOT position in list of the first isotope. +* JSOT position in list of the second isotope. +* HNAMIS1 local name of the first isotope: +* HNAMIS1(1:8) is the local isotope name; +* HNAMIS1(9:12) is a suffix function of the mixture index. +* HNAMIS2 local name of the second isotope +* HNAMIS2(1:8) is the local isotope name; +* HNAMIS2(9:12) is a suffix function of the mixture index. +* +*----------------------------------------------------------------------- +* + USE GANLIB +*---- +* SUBROUTINE ARGUMENTS +*---- + TYPE(C_PTR) IPLIB + INTEGER NGRO,ISOT,JSOT + CHARACTER HNAMIS1*12,HNAMIS2*12 +*---- +* LOCAL VARIABLES +*---- + PARAMETER(MAXNPT=12) + TYPE(C_PTR) JPLIB,IP1,IP2,JP1,JP2,KP1,KP2 + REAL SIGQT1(MAXNPT),SIGQT2(MAXNPT),WSLD1(MAXNPT**2), + 1 WSLD2(MAXNPT**2) + DOUBLE PRECISION SUMA1,SUMB1,SUMA2,SUMB2 + INTEGER, ALLOCATABLE, DIMENSION(:) :: NFS1,NFS2 + REAL, ALLOCATABLE, DIMENSION(:) :: TBIN1,TBIN2,EBIN,DBIN,PROB1, + 1 PROB2 + DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: COMOM +*---- +* SCRATCH STORAGE ALLOCATION +*---- + ALLOCATE(NFS1(NGRO),NFS2(NGRO)) +* + JPLIB=LCMGID(IPLIB,'ISOTOPESLIST') + IP1=LCMGIL(JPLIB,ISOT) ! set ISOT-th isotope + IP2=LCMGIL(JPLIB,JSOT) ! set JSOT-th isotope + CALL LCMLEN(IP1,'BIN-NFS',LENGT1,ITYLCM) + CALL LCMLEN(IP2,'BIN-NFS',LENGT2,ITYLCM) + IF((LENGT1.EQ.0).OR.(LENGT1.NE.LENGT2)) CALL XABORT('LIBCOR: UNA' + 1 //'BLE TO FIND CONSISTENT BIN TYPE INFORMATION.') + CALL LCMGET(IP1,'BIN-NFS',NFS1) + CALL LCMGET(IP2,'BIN-NFS',NFS2) + LBIN=0 + IGRMIN=1 + IGRMAX=NGRO + DO 10 IGRP=NGRO,1,-1 + IF(NFS1(IGRP).NE.NFS2(IGRP)) CALL XABORT('INVALID BIN INFO.') + IF((IGRMAX.EQ.IGRP).AND.(NFS1(IGRP).EQ.0)) IGRMAX=IGRP-1 + LBIN=LBIN+NFS1(IGRP) + 10 CONTINUE + DO 20 IGRP=1,NGRO + IF((IGRMIN.EQ.IGRP).AND.(NFS1(IGRP).EQ.0)) IGRMIN=IGRP+1 + 20 CONTINUE + ALLOCATE(TBIN1(LBIN),TBIN2(LBIN),EBIN(LBIN+1),DBIN(LBIN)) + CALL LCMGET(IP1,'BIN-ENERGY',EBIN) + CALL LCMGET(IP1,'BIN-NTOT0',TBIN1) + CALL LCMGET(IP2,'BIN-NTOT0',TBIN2) + CALL LCMSIX(IP1,'PT-TABLE',1) + CALL LCMSIX(IP2,'PT-TABLE',1) +*--- +* LOOP OVER THE RESONANT ENERGY GROUPS. +*--- + LBIN=0 + JP1=LCMGID(IP1,'GROUP-PT') + JP2=LCMGID(IP2,'GROUP-PT') + DO 130 IGRP=IGRMIN,IGRMAX + SUMA1=0.0D0 + SUMB1=0.0D0 + SUMA2=0.0D0 + SUMB2=0.0D0 + DO 30 IGF=1,NFS1(IGRP) + SIGTA=MAX(0.002,TBIN1(LBIN+IGF)) + SIGTB=MAX(0.002,TBIN2(LBIN+IGF)) + DELM=LOG(EBIN(LBIN+IGF)/EBIN(LBIN+IGF+1)) + SUMA1=SUMA1+TBIN1(LBIN+IGF)*DELM + SUMB1=SUMB1+SIGTA*DELM + SUMA2=SUMA2+TBIN2(LBIN+IGF)*DELM + SUMB2=SUMB2+SIGTB*DELM + TBIN1(LBIN+IGF)=SIGTA + TBIN2(LBIN+IGF)=SIGTB + DBIN(LBIN+IGF)=DELM + 30 CONTINUE + DO 40 IGF=1,NFS1(IGRP) + TBIN1(LBIN+IGF)=TBIN1(LBIN+IGF)*REAL(SUMA1/SUMB1) + TBIN2(LBIN+IGF)=TBIN2(LBIN+IGF)*REAL(SUMA2/SUMB2) + 40 CONTINUE +* + CALL LCMLEL(JP1,IGRP,N1,ITYLCM) + CALL LCMLEL(JP2,IGRP,N2,ITYLCM) + IF((N1.EQ.0).OR.(N2.EQ.0)) GO TO 120 + KP1=LCMGIL(JP1,IGRP) + KP2=LCMGIL(JP2,IGRP) + CALL LCMLEN(KP1,'SIGQT-SIGS',NQT1,ITYLCM) + CALL LCMLEN(KP2,'SIGQT-SIGS',NQT2,ITYLCM) + CALL LCMLEN(KP1,'PROB-TABLE',NQT10,ITYLCM) + CALL LCMLEN(KP2,'PROB-TABLE',NQT20,ITYLCM) + ALLOCATE(PROB1(NQT10),PROB2(NQT20)) + CALL LCMGET(KP1,'PROB-TABLE',PROB1) + CALL LCMGET(KP2,'PROB-TABLE',PROB2) + DO 50 I=1,NQT1 + SIGQT1(I)=PROB1(MAXNPT+I) + 50 CONTINUE + DO 60 I=1,NQT2 + SIGQT2(I)=PROB2(MAXNPT+I) + 60 CONTINUE +* + ALLOCATE(COMOM(NQT1*NQT2)) + CALL LIBCOM(NFS1(IGRP),DBIN(LBIN+1),TBIN1(LBIN+1), + 1 TBIN2(LBIN+1),NQT1,NQT2,COMOM) + CALL LIBOMG(NQT1,1-NQT1/2,SIGQT1,NQT2,1-NQT2/2,SIGQT2, + 1 COMOM,WSLD1) + DEALLOCATE(COMOM) +*--- +* CHECK NORMALIZATION OF THE CORRELATED WEIGHT MATRIX. +*--- + DO 80 I=1,NQT1 + SUM=0.0 + DO 70 J=1,NQT2 + SUM=SUM+WSLD1((J-1)*NQT1+I) + 70 CONTINUE + IF(ABS(SUM-PROB1(I)).GT.1.0E-4) THEN + CALL XABORT('LIBCOR: BAD NORMALIZATION EXCEPTION(1).') + ENDIF + 80 CONTINUE + DO 100 I=1,NQT2 + SUM=0.0 + DO 90 J=1,NQT1 + SUM=SUM+WSLD1((I-1)*NQT1+J) + 90 CONTINUE + IF(ABS(SUM-PROB2(I)).GT.1.0E-4) THEN + CALL XABORT('LIBCOR: BAD NORMALIZATION EXCEPTION(2).') + ENDIF + 100 CONTINUE + DEALLOCATE(PROB2,PROB1) +* + CALL LCMPUT(KP1,HNAMIS2,NQT1*NQT2,2,WSLD1) + DO 115 I=1,NQT1 + DO 110 J=1,NQT2 + WSLD2((I-1)*NQT2+J)=WSLD1((J-1)*NQT1+I) + 110 CONTINUE + 115 CONTINUE + CALL LCMPUT(KP2,HNAMIS1,NQT2*NQT1,2,WSLD2) + 120 LBIN=LBIN+NFS1(IGRP) + 130 CONTINUE +* + CALL LCMSIX(IP2,' ',2) + CALL LCMSIX(IP1,' ',2) + DEALLOCATE(DBIN,EBIN,TBIN2,TBIN1) +*---- +* SCRATCH STORAGE DEALLOCATION +*---- + DEALLOCATE(NFS2,NFS1) + RETURN + END |