From 7dfcc480ba1e19bd3232349fc733caef94034292 Mon Sep 17 00:00:00 2001 From: stainer_t Date: Mon, 8 Sep 2025 13:48:49 +0200 Subject: Initial commit from Polytechnique Montreal --- Dragon/src/MOCDS2.f | 207 ++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 207 insertions(+) create mode 100644 Dragon/src/MOCDS2.f (limited to 'Dragon/src/MOCDS2.f') diff --git a/Dragon/src/MOCDS2.f b/Dragon/src/MOCDS2.f new file mode 100644 index 0000000..b6e72c8 --- /dev/null +++ b/Dragon/src/MOCDS2.f @@ -0,0 +1,207 @@ +*DECK MOCDS2 + SUBROUTINE MOCDS2(SUBDSC,LC,M,N,H,NOM,NZON,TR,SC,W,NFI,DIAGF, + 1 DIAGQ,CA,CQ,PREV,NEXT,DINV2,A2,B2) +* +*----------------------------------------------------------------------- +* +*Purpose: +* Calculation of contribution in second-order ACA coefficients on one +* track. +* +*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): I. Suslov and R. Le Tellier +* +*Parameters: input +* SUBDSC ACA coefficients calculation subroutine. +* LC dimension of vector MCU. +* M number of material mixtures. +* N number of elements for this track. +* H tracking widths. +* NOM integer tracking elements. +* NZON index-number of the mixture type assigned to each volume. +* TR macroscopic total cross section. +* SC macroscopic P0 scattering cross section. +* W weight associated with this track. +* NFI total number of volumes for which specific values +* of the neutron flux and reactions rates are required. +* +*Parameters: input/output +* CA undefined. +* CQ undefined. +* DIAGQ undefined. +* DIAGF undefined. +* +*Parameters: scratch +* PREV undefined. +* NEXT undefined. +* DINV2 undefined. +* A2 undefined. +* B2 undefined. +* +*----------------------------------------------------------------------- +* + IMPLICIT NONE +* +*--- +* SUBROUTINE ARGUMENTS +*--- + INTEGER LC,M,N,NFI,NZON(NFI),NOM(N),PREV(N),NEXT(N) + DOUBLE PRECISION W,H(N),CA(LC),DIAGF(NFI),DINV2(N),A2(N),B2(N) + REAL TR(0:M),SC(0:M),DIAGQ(NFI),CQ(LC) + EXTERNAL SUBDSC +*--- +* LOCAL VARIABLES +*--- + INTEGER IBCV + DOUBLE PRECISION DMINV,DMAX + PARAMETER(DMINV=2.D-2,IBCV=-7) + DOUBLE PRECISION WW,AAW,CAW,AQW,CQW,CAWP,CQWP,B,A,DINV,DH,DHP, + 1 BN,AN,DINVN,BP,AP,DINVP + INTEGER I,I1,NOMI,NZI,NOMIN,NZIN,NOMIP,NZIP,ICN,ICP +* + DMAX=1.D0/DMINV + WW=W +*--- +* CALCULATE COEFFICIENTS OF THIS TRACK +*--- +* MCGDS2A: Tabulated Exponentials +* MCGDS2E: Exact Exponentials + CALL SUBDSC(N,M,NFI,NOM,NZON,H,TR,SC,DINV2,B2,A2) +*---- +* CONSTRUCTION OF ACA MATRICES +*--- + NOMIN=0 + NZIN=0 + NOMI=0 + NZI=0 + DHP=0.0D0 + DH=0.0D0 + AN=0.0D0 + BN=0.0D0 + AP=0.0D0 + BP=0.0D0 + DO I=1,N + ICN=NEXT(I) + ICP=PREV(I) + IF (I.GT.1) THEN + NOMIP=NOMI + NZIP=NZI + NOMI=NOMIN + NZI=NZIN + ELSE + NOMIP=NOM(N) + NZIP=NZON(NOMIP) + NOMI=NOM(1) + NZI=NZON(NOMI) + ENDIF + IF (I.LT.N) THEN + NOMIN=NOM(I+1) + NZIN=NZON(NOMIN) + ELSE + NOMIN=NOM(1) + NZIN=NZON(NOMIN) + ENDIF + IF (NZI.GE.0) THEN +* ----------- +* Volume Cell +* ----------- + DINV=DINV2(I) + B=B2(I) + A=A2(I) + IF (NZIN.GE.0) THEN +* next cell is a volume + IF (I.EQ.N) THEN + I1=1 + ELSE + I1=I+1 + ENDIF + DINVN=DINV2(I1) + BN=B2(I1) + AN=A2(I1) + IF (ABS(DINV+DINVN).LT.DMINV) THEN + DH=DMAX + ELSE + DH=1.D0/(DINV+DINVN) + ENDIF + ELSEIF (NZIN.EQ.IBCV) THEN +* next cell is a fixed boundary condition + DH=1.D0/(1.D0+DINV) + ENDIF + IF (NZIP.GE.0) THEN +* previous cell is a volume + IF (I.EQ.1) THEN + I1=N + ELSE + I1=I-1 + ENDIF + DINVP=DINV2(I1) + BP=B2(I1) + AP=A2(I1) + IF (ABS(DINV+DINVP).LT.DMINV) THEN + DHP=DMAX + ELSE + DHP=1.D0/(DINV+DINVP) + ENDIF + ELSEIF (NZIP.EQ.IBCV) THEN +* previous cell is a fixed boundary condition + DHP=1.D0/(1.D0+DINV) + ENDIF +* +* assembling coefficients + AAW=DH*A + AQW=DH*B + IF(NZIN.GE.0) THEN +* next cell is a volume + CAW=DH*AN + CQW=DH*BN + ELSE +* next cell is a voided boundary condition + CAW=0.D0 + CQW=0.D0 + ENDIF +* + AAW=AAW+DHP*A + AQW=AQW+DHP*B + IF(NZIP.GE.0) THEN +* previous cell is a volume + CAWP=DHP*AP + CQWP=DHP*BP + ELSE +* previous cell is a voided boundary condition + CAWP=0.D0 + CQWP=0.D0 + ENDIF +* +* assembling matrices + DIAGF(NOMI)=DIAGF(NOMI)+AAW*WW + DIAGQ(NOMI)=DIAGQ(NOMI)-REAL(W*AQW) + IF(ICN.GT.0) THEN +* next cell is a volume different from this one + CA(ICN)=CA(ICN)-CAW*WW + CQ(ICN)=CQ(ICN)+REAL(W*CQW) + ELSE +* next cell is a voided boundary or a volume identical to this one + DIAGF(NOMI)=DIAGF(NOMI)-CAW*WW + DIAGQ(NOMI)=DIAGQ(NOMI)+REAL(W*CQW) + ENDIF + IF(ICP.GT.0) THEN +* previous cell is a volume different from this one + CA(ICP)=CA(ICP)-CAWP*WW + CQ(ICP)=CQ(ICP)+REAL(W*CQWP) + ELSE +* previous cell is a voided boundary or a volume identical to this one + DIAGF(NOMI)=DIAGF(NOMI)-CAWP*WW + DIAGQ(NOMI)=DIAGQ(NOMI)+REAL(W*CQWP) + ENDIF + ENDIF +* ----------- + ENDDO +* + RETURN + END -- cgit v1.2.3