diff options
Diffstat (limited to 'Trivac/src/TRIVAA.f')
| -rwxr-xr-x | Trivac/src/TRIVAA.f | 303 |
1 files changed, 303 insertions, 0 deletions
diff --git a/Trivac/src/TRIVAA.f b/Trivac/src/TRIVAA.f new file mode 100755 index 0000000..acee907 --- /dev/null +++ b/Trivac/src/TRIVAA.f @@ -0,0 +1,303 @@ +*DECK TRIVAA + SUBROUTINE TRIVAA(NENTRY,HENTRY,IENTRY,JENTRY,KENTRY) +* +*----------------------------------------------------------------------- +* +*Purpose: +* TRIVAC type (3-D and ADI) system matrix assembly operator. +* +*Copyright: +* Copyright (C) 2002 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/output +* NENTRY number of LCM objects or files used by the operator. +* HENTRY name of each LCM object or file: +* HENTRY(1): create or modification type(L_SYSTEM); +* HENTRY(2): read-only type(L_MACROLIB) (unperturbed); +* HENTRY(3): read-only type(L_TRACK); +* HENTRY(4): optional read-only type(L_MACROLIB) (perturbed). +* IENTRY type of each LCM object or file: +* =1 LCM memory object; =2 XSM file; =3 sequential binary file; +* =4 sequential ascii file. +* JENTRY access of each LCM object or file: +* =0 the LCM object or file is created; +* =1 the LCM object or file is open for modifications; +* =2 the LCM object or file is open in read-only mode. +* KENTRY LCM object address or file unit number. +* +*Comments: +* The TRIVAA: calling specifications are: +* SYST := TRIVAA: [ SYST ] MACRO TRACK [ DMACRO ] :: (trivaa\_data) ; +* where +* SYST : name of the \emph{lcm} object (type L\_SYSTEM) containing the +* system matrices. If SYST appears on the RHS, the system matrices +* previously stored in SYST are kept. +* MACRO : name of the \emph{lcm} object (type L\_MACROLIB) containing the +* macroscopic cross sections and diffusion coefficients. +* TRACK : name of the \emph{lcm} object (type L\_TRIVAC) containing the +* TRIVAC \emph{tracking}. +* DMACRO : name of the \emph{lcm} object (type L\_MACROLIB) containing +* derivatives or perturbations of the macroscopic cross sections and +* diffusion coefficients. If DMACRO is given, only the derivatives or +* perturbations of the system matrices are computed. +* trivaa\_data : structure containing the data to module TRIVAA: +* +*----------------------------------------------------------------------- +* + 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 TEXT4*4,TEXT11*12,TEXT12*12,HSMG*131,TITLE*72,CNAM*12 + DOUBLE PRECISION DFLOTT + INTEGER IGP(NSTATE),IPAR(NSTATE),ITR(NSTATE) + LOGICAL LDIFF + TYPE(C_PTR) IPSYS,JPSYS,KPSYS,IPMACR,JPMACR,KPMACR,IPTRK,IPMACP + INTEGER, DIMENSION(:), ALLOCATABLE :: MAT + REAL, DIMENSION(:), ALLOCATABLE :: VOL,UN,VII +*---- +* PARAMETER VALIDATION. +*---- + IF(NENTRY.LE.2) CALL XABORT('TRIVAA: THREE PARAMETERS EXPECTED.') + IF((IENTRY(1).NE.1).AND.(IENTRY(1).NE.2)) CALL XABORT('TRIVAA: L' + 1 //'CM OBJECT EXPECTED AT LHS.') + IF((JENTRY(1).NE.0).AND.(JENTRY(1).NE.1)) CALL XABORT('TRIVAA: E' + 1 //'NTRY IN CREATE OR MODIFICATION MODE EXPECTED.') + IF((JENTRY(3).NE.2).OR.((IENTRY(3).NE.1).AND.(IENTRY(3).NE.2))) + 1 CALL XABORT('TRIVAA: LCM OBJECT IN READ-ONLY MODE EXPECTED AT F' + 2 //'IRST RHS.') + IF((JENTRY(2).NE.2).OR.((IENTRY(2).NE.1).AND.(IENTRY(2).NE.2))) + 1 CALL XABORT('TRIVAA: LCM OBJECT IN READ-ONLY MODE EXPECTED AT S' + 2 //'ECOND RHS.') + CALL LCMGTC(KENTRY(3),'SIGNATURE',12,TEXT11) + IF(TEXT11.NE.'L_TRACK') THEN + TEXT12=HENTRY(3) + CALL XABORT('TRIVAA: SIGNATURE OF '//TEXT12//' IS '//TEXT11// + 1 '. L_TRACK EXPECTED.') + ENDIF + CALL LCMGTC(KENTRY(3),'TRACK-TYPE',12,TEXT11) + IF(TEXT11.NE.'TRIVAC') THEN + TEXT12=HENTRY(3) + CALL XABORT('TRIVAA: TRACK-TYPE OF '//TEXT12//' IS '//TEXT11 + 1 //'. TRIVAC EXPECTED.') + ENDIF + TEXT11='L_SYSTEM' + IPSYS=KENTRY(1) + CALL LCMPTC(IPSYS,'SIGNATURE',12,TEXT11) + IPMACR=KENTRY(2) + IPTRK=KENTRY(3) + TEXT12=HENTRY(2) + CALL LCMPTC(IPSYS,'LINK.MACRO',12,TEXT12) + TEXT12=HENTRY(3) + CALL LCMPTC(IPSYS,'LINK.TRACK',12,TEXT12) +*---- +* RECOVER GENERAL TRACKING INFORMATION. +*---- + CALL LCMGET(IPTRK,'STATE-VECTOR',IGP) + NEL=IGP(1) + NLF=IGP(30) + ISCAT=IGP(32) + LDIFF=(ISCAT.LT.0) + ISCAT=ABS(ISCAT) + IF((NLF.NE.0).AND.(IGP(31).NE.1)) CALL XABORT('TRIVAA: ONLY SPN ' + 1 //'DISCRETIZATIONS ARE ALLOWED.') + ITY=2 + IF(IGP(12).EQ.2) ITY=3 + ALLOCATE(MAT(NEL),VOL(NEL)) + CALL LCMGET(IPTRK,'MATCOD',MAT) + CALL LCMGET(IPTRK,'VOLUME',VOL) + CALL LCMLEN(IPTRK,'TITLE',LENGT,ITYLCM) + IF(LENGT.GT.0) THEN + CALL LCMGTC(IPTRK,'TITLE',72,TITLE) + ELSE + TITLE='*** NO TITLE PROVIDED ***' + ENDIF +*---- +* RECOVER MACROLIB PARAMETERS. +*---- + CALL LCMGTC(IPMACR,'SIGNATURE',12,TEXT11) + IF(TEXT11.NE.'L_MACROLIB') THEN + TEXT12=HENTRY(2) + CALL XABORT('TRIVAA: SIGNATURE OF '//TEXT12//' IS '//TEXT11// + 1 '. L_MACROLIB EXPECTED.') + ENDIF + CALL LCMGET(IPMACR,'STATE-VECTOR',IPAR) + NGRP=IPAR(1) + NBMIX=IPAR(2) + NANI=IPAR(3) + NBFIS=IPAR(4) + NALBP=IPAR(8) + IF(IGP(4).GT.NBMIX) THEN + WRITE(HSMG,'(46HTRIVAA: THE NUMBER OF MIXTURES IN THE TRACKING, + 1 2H (,I5,51H) IS GREATER THAN THE NUMBER OF MIXTURES IN THE MAC, + 2 7HROLIB (,I5,2H).)') IGP(4),NBMIX + CALL XABORT(HSMG) + ENDIF +* + IMPX=1 + IASM=0 + IPR=0 + IUNIT=0 + IOVEL=0 + NSTEP=0 + 10 CALL REDGET(INDIC,NITMA,FLOTT,TEXT4,DFLOTT) + IF(INDIC.EQ.10) GO TO 30 + IF(INDIC.NE.3) CALL XABORT('TRIVAA: CHARACTER DATA EXPECTED.') + IF(TEXT4.EQ.'EDIT') THEN + CALL REDGET(INDIC,IMPX,FLOTT,TEXT4,DFLOTT) + IF(INDIC.NE.1) CALL XABORT('TRIVAA: INTEGER DATA EXPECTED(1).') + ELSE IF(TEXT4.EQ.'SKIP') THEN +* OPTION TO SKIP THE SYSTEM MATRIX ASSEMBLY (DO NOT SKIP THE +* LDLT FACTORIZATION). + IASM=1 + ELSE IF(TEXT4.EQ.'DERI') THEN + IPR=1 + WRITE(6,'(/43H TRIVAA: USE DERIVATIVE OF SYSTEM MATRICES.)') + ELSE IF(TEXT4.EQ.'PERT') THEN + IPR=2 + WRITE(6,'(/41H TRIVAA: PERTURBATION OF SYSTEM MATRICES.)') + ELSE IF(TEXT4.EQ.'UNIT') THEN +* COMPUTE THE UNITARY WEIGHTING MATRIX. + IUNIT=1 + ALLOCATE(UN(NBMIX)) + UN(:NBMIX)=1.0 + CALL TRIDIG('RM',IPTRK,IPSYS,IMPX,NBMIX,NEL,0,MAT,VOL,UN) + DEALLOCATE(UN) + ELSE IF(TEXT4.EQ.'OVEL') THEN +* COMPUTE THE RECIPROCAL NEUTRON VELOCITIES MATRIX. + IOVEL=1 + JPMACR=LCMGID(IPMACR,'GROUP') + ALLOCATE(VII(NBMIX)) + DO 25 IGR=1,NGRP + KPMACR=LCMGIL(JPMACR,IGR) + CALL LCMLEN(KPMACR,'OVERV',LENGT,ITYLCM) + IF(LENGT.EQ.0) THEN + CALL XABORT('TRIVAA: NO ''VELOCITY'' INFORMATION.') + ELSE IF(LENGT.GT.NBMIX) THEN + CALL XABORT('TRIVAA: INVALID LENGTH FOR ''VELOCITY'' IN' + 1 //'FORMATION.') + ENDIF + CALL LCMGET(KPMACR,'OVERV',VII) + WRITE (CNAM,'(1HV,2I3.3)') IGR,IGR + CALL TRIDIG(CNAM,IPTRK,IPSYS,IMPX,NBMIX,NEL,0,MAT,VOL,VII) + 25 CONTINUE + DEALLOCATE(VII) + ELSE IF(TEXT4.EQ.';') THEN + GO TO 30 + ELSE + CALL XABORT('TRIVAA: '//TEXT4//' IS AN INVALID KEY WORD.') + ENDIF + GO TO 10 +*---- +* 2-MACROLIBS PERTURBATION CALCULATION. +*---- + 30 IF(IPR.GT.0) THEN + IF(NENTRY.LE.3) CALL XABORT('TRIVAA: 4 PARAMETERS EXPECTED WIT' + 1 //'H DERI OR PERT OPTIONS.') + IF((JENTRY(4).NE.2).OR.((IENTRY(4).NE.1).AND.(IENTRY(4).NE.2))) + 1 CALL XABORT('TRIVAA: LINKED LIST OR XSM FILE IN READ-ONLY MODE' + 2 //' EXPECTED AT THIRD RHS.') + IPMACP=KENTRY(4) + CALL LCMGTC(IPMACP,'SIGNATURE',12,TEXT11) + IF(TEXT11.NE.'L_MACROLIB') THEN + TEXT12=HENTRY(4) + CALL XABORT('TRIVAA: SIGNATURE OF '//TEXT12//' IS ' + 1 //TEXT11//'. L_MACROLIB EXPECTED.') + ENDIF + CALL LCMGET(IPMACP,'STATE-VECTOR',IPAR) + NSTEP=IPAR(11) + IF((IPAR(1).NE.NGRP).OR.(IPAR(2).GT.NBMIX)) THEN + WRITE(HSMG,'(43HTRIVAA: INCONSISTENT PERTURBATION MACROLIB , + 1 1H'',A12,8H''. NGRP=,2I5,7H NBMIX=,2I9)') HENTRY(4),IPAR(1), + 2 NGRP,IPAR(2),NBMIX + CALL XABORT(HSMG) + ENDIF + ENDIF +*---- +* SET THE STATE VECTOR FOR THE L_SYSTEM OBJECT +*---- + ITR(:NSTATE)=0 + ITR(1)=NGRP + ITR(2)=IGP(11) + ITR(3)=0 + ITR(4)=ITY + IF((NLF.GT.0).AND.(ITY.GE.3)) ITR(4)=10+ITR(4) + IF(IUNIT.EQ.1) ITR(5)=1 + ITR(6)=NSTEP + ITR(7)=NBMIX + NAN=MIN(ISCAT,NANI) + ITR(8)=NLF + ITR(9)=IPR + CALL LCMPUT(IPSYS,'STATE-VECTOR',NSTATE,1,ITR) +*---- +* SYSTEM MATRIX ASSEMBLY. +*---- + IF((IASM.EQ.0).AND.(IPR.EQ.0)) THEN + IF(NLF.EQ.0) THEN +* DIFFUSION THEORY. + CALL TRISYS(IPTRK,IPMACR,IPMACR,IPSYS,IMPX,NGRP,NEL,NBFIS, + 1 NALBP,IPR,MAT,VOL,NBMIX) + ELSE +* SIMPLIFIED PN THEORY. + CALL TRISPS(IPTRK,IPMACR,IPMACR,IPSYS,IMPX,NGRP,NEL,NLF, + 1 NAN,NBFIS,NALBP,LDIFF,IPR,MAT,VOL,NBMIX) + ENDIF + ELSE IF((IASM.EQ.1).AND.(IPR.EQ.0)) THEN +* PERFORM FACTORIZATION WITHOUT ASSEMBLY. + DO 40 I=1,NGRP + WRITE(TEXT11,'(1HA,2I3.3)') I,I + CALL MTLDLF(TEXT11,IPTRK,IPSYS,ITY,IMPX) + 40 CONTINUE + ELSE IF((IPR.GT.0).AND.(NSTEP.EQ.0)) THEN +* ASSEMBLY OF PERTURBED SYSTEM MATRICES (NO STEP DIRECTORIES). + IF(NLF.EQ.0) THEN +* DIFFUSION THEORY. + CALL TRISYS(IPTRK,IPMACR,IPMACP,IPSYS,IMPX,NGRP,NEL,NBFIS, + 1 NALBP,IPR,MAT,VOL,NBMIX) + ELSE +* SIMPLIFIED PN THEORY. + CALL TRISPS(IPTRK,IPMACR,IPMACP,IPSYS,IMPX,NGRP,NEL,NLF, + 1 NAN,NBFIS,NALBP,LDIFF,IPR,MAT,VOL,NBMIX) + ENDIF + ELSE IF(NSTEP.GT.0) THEN +* ASSEMBLY OF PERTURBED SYSTEM MATRICES (WITH STEP DIRECTORIES). + JPMACR=LCMGID(IPMACP,'STEP') + JPSYS=LCMLID(IPSYS,'STEP',NSTEP) + DO 50 ISTEP=1,NSTEP + KPMACR=LCMGIL(JPMACR,ISTEP) + KPSYS=LCMDIL(JPSYS,ISTEP) + CALL LCMPUT(KPSYS,'STATE-VECTOR',NSTATE,1,ITR) + IF(NLF.EQ.0) THEN +* DIFFUSION THEORY. + CALL TRISYS(IPTRK,IPMACR,KPMACR,KPSYS,IMPX,NGRP,NEL,NBFIS, + 1 NALBP,IPR,MAT,VOL,NBMIX) + ELSE +* SIMPLIFIED PN THEORY. + CALL TRISPS(IPTRK,IPMACR,KPMACR,KPSYS,IMPX,NGRP,NEL,NLF, + 1 NAN,NBFIS,NALBP,LDIFF,IPR,MAT,VOL,NBMIX) + ENDIF + 50 CONTINUE + ELSE + CALL XABORT('TRIVAA: INVALID REQUEST.') + ENDIF +* + IF(IMPX.GE.3) CALL LCMLIB(IPSYS) +*---- +* RELEASE GENERAL TRACKING INFORMATION. +*---- + DEALLOCATE(VOL,MAT) + RETURN + END |