Initial commit from Polytechnique Montreal

1 files changed, 174 insertions, 0 deletions

diff --git a/Dragon/src/DUO001.f b/Dragon/src/DUO001.f
new file mode 100644
index 0000000..a85c393
--- /dev/null
+++ b/Dragon/src/DUO001.f
@@ -0,0 +1,174 @@
+*DECK DUO001
+      SUBROUTINE DUO001(IPMAC,IPRINT,NMIX,NGRP,NFIS,IDIV,ZKEFF,RHS,LHS,
+     > FLUX,AFLUX)
+*
+*-----------------------------------------------------------------------
+*
+*Purpose:
+* Processing one of the two macrolibs and return mixture-dependent
+* RHS and LHS matrices.
+*
+*Copyright:
+* Copyright (C) 2013 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   macrolib.
+* IPRINT  print parameter.
+* NMIX    number of mixtures.
+* NGRP    number of energy groups.
+* NFIS    number of fissile isotopes.
+* IDIV    type of divergence term processing (=0: no processing;
+*         =1: direct processing; =2: adjoint processing;
+*         =3: direct-adjoint processing).
+*
+*Parameters: output
+* ZKEFF   effective multiplication factor.
+* RHS     absorption macroscopic cross-section matrix.
+* LHS     production macroscopic cross-section matrix.
+* FLUX    integrated direct flux.
+* AFLUX   integrated adjoint flux.
+*
+*-----------------------------------------------------------------------
+*
+      USE GANLIB
+*----
+*  SUBROUTINE ARGUMENTS
+*----
+      TYPE(C_PTR) IPMAC
+      INTEGER IPRINT,NMIX,NGRP,NFIS,IDIV
+      REAL ZKEFF,RHS(NGRP,NGRP,NMIX),LHS(NGRP,NGRP,NMIX),
+     > FLUX(NGRP,NMIX),AFLUX(NGRP,NMIX)
+*----
+*  LOCAL VARIABLES
+*----
+      TYPE(C_PTR) JPMAC,KPMAC
+      DOUBLE PRECISION SUM
+*----
+*  ALLOCATABLE ARRAYS
+*----
+      INTEGER, ALLOCATABLE, DIMENSION(:) :: IJJ,NJJ,IPOS
+      REAL, ALLOCATABLE, DIMENSION(:) :: VOL,GAR,GAR2,DLK,ALK,V,W
+      REAL, ALLOCATABLE, DIMENSION(:,:) :: NUF
+      REAL, ALLOCATABLE, DIMENSION(:,:,:) :: CHI
+*----
+*  SCRATCH STORAGE ALLOCATION
+*----
+      ALLOCATE(IJJ(NMIX),NJJ(NMIX),IPOS(NMIX))
+      ALLOCATE(VOL(NMIX),GAR(NMIX),CHI(NMIX,NFIS,NGRP),NUF(NMIX,NFIS),
+     > GAR2(NMIX*NGRP))
+*----
+*  COMPUTE THE RHS AND LHS MATRICES
+*----
+      RHS(:NGRP,:NGRP,:NMIX)=0.0
+      LHS(:NGRP,:NGRP,:NMIX)=0.0
+      CALL LCMGET(IPMAC,'K-EFFECTIVE',ZKEFF)
+      IF(IPRINT.GT.1) WRITE(6,'(35H DUO001: EFFECTIVE MULTIPLICATION F,
+     > 6HACTOR=,1P,E12.5)') ZKEFF
+      CALL LCMGET(IPMAC,'VOLUME',VOL)
+      JPMAC=LCMGID(IPMAC,'GROUP')
+      DO IGR=1,NGRP
+        KPMAC=LCMGIL(JPMAC,IGR)
+        CALL LCMGET(KPMAC,'CHI',CHI(1,1,IGR))
+      ENDDO
+      DO IGR=1,NGRP
+        KPMAC=LCMGIL(JPMAC,IGR)
+        CALL LCMGET(KPMAC,'FLUX-INTG',GAR)
+        DO IBM=1,NMIX
+          FLUX(IGR,IBM)=GAR(IBM)/VOL(IBM)
+        ENDDO
+        CALL LCMLEN(KPMAC,'NWAT0',ILONG,ITYLCM)
+        IF(ILONG.EQ.NMIX) THEN
+          CALL LCMGET(KPMAC,'NWAT0',GAR)
+          DO IBM=1,NMIX
+            AFLUX(IGR,IBM)=GAR(IBM)
+          ENDDO
+        ELSE
+          AFLUX(:NMIX,IBM)=1.0
+        ENDIF
+        CALL LCMGET(KPMAC,'NTOT0',GAR)
+        CALL LCMGET(KPMAC,'SCAT00',GAR2)
+        CALL LCMGET(KPMAC,'NJJS00',NJJ)
+        CALL LCMGET(KPMAC,'IJJS00',IJJ)
+        CALL LCMGET(KPMAC,'IPOS00',IPOS)
+        DO IBM=1,NMIX
+          IPOSDE=IPOS(IBM)
+          DO JGR=IJJ(IBM),IJJ(IBM)-NJJ(IBM)+1,-1
+            RHS(IGR,JGR,IBM)=RHS(IGR,JGR,IBM)-GAR2(IPOSDE) ! IGR <-- JGR
+            IPOSDE=IPOSDE+1
+          ENDDO
+          RHS(IGR,IGR,IBM)=RHS(IGR,IGR,IBM)+GAR(IBM)
+        ENDDO
+        CALL LCMGET(KPMAC,'NUSIGF',NUF)
+        DO IBM=1,NMIX
+          DO IFIS=1,NFIS
+            DO JGR=1,NGRP
+              LHS(JGR,IGR,IBM)=LHS(JGR,IGR,IBM)+CHI(IBM,IFIS,JGR)*
+     >        NUF(IBM,IFIS)
+            ENDDO
+          ENDDO
+        ENDDO
+      ENDDO
+*----
+*  INTRODUCE THE DIRECT OR ADJOINT DIVERGENCE COMPONENT IN THE RHS
+*  MATRIX
+*----
+      DO IBM=1,NMIX
+        IF(IDIV.EQ.1) THEN
+          DO JGR=1,NGRP
+            SUM=0.0D0
+            DO IGR=1,NGRP
+              SUM=SUM+(RHS(JGR,IGR,IBM)-LHS(JGR,IGR,IBM)/ZKEFF)*
+     >        FLUX(IGR,IBM)
+            ENDDO
+            RHS(JGR,JGR,IBM)=RHS(JGR,JGR,IBM)-REAL(SUM)/FLUX(JGR,IBM)
+          ENDDO
+        ELSE IF(IDIV.EQ.2) THEN
+          DO IGR=1,NGRP
+            SUM=0.0D0
+            DO JGR=1,NGRP
+              SUM=SUM+(RHS(JGR,IGR,IBM)-LHS(JGR,IGR,IBM)/ZKEFF)*
+     >        AFLUX(JGR,IBM)
+            ENDDO
+            RHS(IGR,IGR,IBM)=RHS(IGR,IGR,IBM)-REAL(SUM)/AFLUX(IGR,IBM)
+          ENDDO
+        ELSE IF(IDIV.EQ.3) THEN
+          ALLOCATE(DLK(NGRP),ALK(NGRP))
+          DO JGR=1,NGRP
+            SUM=0.0D0
+            DO IGR=1,NGRP
+              SUM=SUM+(RHS(JGR,IGR,IBM)-LHS(JGR,IGR,IBM)/ZKEFF)*
+     >        FLUX(IGR,IBM)
+            ENDDO
+            DLK(JGR)=REAL(SUM)
+          ENDDO
+          DO IGR=1,NGRP
+            SUM=0.0D0
+            DO JGR=1,NGRP
+              SUM=SUM+(RHS(JGR,IGR,IBM)-LHS(JGR,IGR,IBM)/ZKEFF)*
+     >        AFLUX(JGR,IBM)
+            ENDDO
+            ALK(IGR)=REAL(SUM)
+          ENDDO
+          ALLOCATE(V(NGRP),W(NGRP))
+          CALL DUO005(NGRP,DLK,ALK,FLUX(1,IBM),AFLUX(1,IBM),V,W)
+          DO IGR=1,NGRP
+            DO JGR=1,NGRP
+              RHS(IGR,JGR,IBM)=RHS(IGR,JGR,IBM)-V(IGR)-W(JGR)
+            ENDDO
+          ENDDO
+          DEALLOCATE(W,V,ALK,DLK)
+        ENDIF
+      ENDDO
+*----
+*  SCRATCH STORAGE DEALLOCATION
+*----
+      DEALLOCATE(GAR2,NUF,CHI,GAR,VOL)
+      DEALLOCATE(IPOS,NJJ,IJJ)
+      RETURN
+      END