summaryrefslogtreecommitdiff
path: root/Donjon/src/RESGEO.f
diff options
context:
space:
mode:
Diffstat (limited to 'Donjon/src/RESGEO.f')
-rw-r--r--Donjon/src/RESGEO.f304
1 files changed, 304 insertions, 0 deletions
diff --git a/Donjon/src/RESGEO.f b/Donjon/src/RESGEO.f
new file mode 100644
index 0000000..0c168f7
--- /dev/null
+++ b/Donjon/src/RESGEO.f
@@ -0,0 +1,304 @@
+*DECK RESGEO
+ SUBROUTINE RESGEO(IPMAP,IPMTX,LX,LY,LZ,NFUEL,IMPX,IGEO,NX,NY,NZ,
+ 1 NCH,NB,NTOT,LNAP,IPCPO)
+*
+*-----------------------------------------------------------------------
+*
+*Purpose:
+* Create and check the fuel-map geometry.
+*
+*Copyright:
+* Copyright (C) 2007 Ecole Polytechnique de Montreal.
+*
+*Author(s):
+* E. Varin, D. Sekki and V. Descotes
+*
+*Update(s):
+* R. Chambon 2014
+*
+*Parameters: input
+* IPMAP pointer to fuel-map information.
+* IPMTX pointer to matex information.
+* LX number of elements along x-axis in geometry.
+* LY number of elements along y-axis in geometry.
+* LZ number of elements along z-axis in geometry.
+* NFUEL number of fuel types.
+* IMPX printing index (=0 for no print).
+* IGEO type of geometry (=7 or =9)
+*
+*Parameters: output
+* NX number of elements along x-axis in fuel map.
+* NY number of elements along y-axis in fuel map.
+* NZ number of elements along z-axis in fuel map.
+* NCH number of reactor channels.
+* NB number of fuel bundles per channel.
+* NTOT total number of fuel bundles.
+* LNAP Flag to call NAP: module to unfold geometry at assembly level
+* IPCPO pointer to multicompo information
+*
+*-----------------------------------------------------------------------
+*
+ USE GANLIB
+*----
+* SUBROUTINE ARGUMENTS
+*----
+ TYPE(C_PTR) IPMAP,IPMTX,IPCPO,IPGNW
+ INTEGER LX,LY,LZ,NFUEL,IGEO,NX,NY,NZ,NCH,NB,NTOT
+ LOGICAL LNAP
+*----
+* LOCAL VARIABLES
+*----
+ PARAMETER(NSTATE=40,IOUT=6,EPSI=1.0E-4)
+ INTEGER ISTATE(NSTATE),JENT(1),IENT(1),JENT2(3),IENT2(3),NCODE(6),
+ 1 ICODE(6)
+ TYPE(C_PTR) KENT(1),KENT2(3)
+ REAL GEOXX(LX+1),GEOYY(LY+1),GEOZZ(LZ+1),GEOSI,GMAPSI,ZCODE(6)
+ CHARACTER HENT(1)*12,HENT2(3)*12,TEXT*12
+ DOUBLE PRECISION DFLOT
+*----
+* ALLOCATABLE STATEMENTS
+*----
+ INTEGER, ALLOCATABLE, DIMENSION(:) :: ISPLX,ISPLY,ISPLZ,MAT
+ REAL, ALLOCATABLE, DIMENSION(:) :: GMAPX,GMAPY,GMAPZ
+*----
+* FUEL-MAP GEOMETRY
+*----
+ IF(IMPX.GT.1)WRITE(IOUT,*)'** CREATING FUEL-MAP GEOMETRY **'
+ CALL LCMSIX(IPMAP,'GEOMAP',1)
+ NENT=1
+ JENT(1)=0
+ HENT(1)='GEOMAP'
+ IENT(1)=1
+ KENT(1)=IPMAP
+ CALL GEOD(NENT,HENT,IENT,JENT,KENT)
+ IF(IMPX.GT.3)CALL LCMLIB(IPMAP)
+*
+ IF(LNAP) THEN
+*----
+* FUEL-MAP GEOMETRY UNFOLDING WITH NAP:
+*----
+ IF(.NOT.C_ASSOCIATED(IPCPO)) THEN
+ CALL XABORT('RESGEO: COMPO LCM OBJECT MISSING AT RHS.')
+ ENDIF
+ CALL REDGET(ITYP,NITMA,FLOT,TEXT,DFLOT)
+ IF(ITYP.NE.3)CALL XABORT('@RESGEO: CHARACTER DATA EXPECTED.')
+ IF(TEXT.NE.':::') CALL XABORT('@RESGEO: ::: keyword EXPECTED.')
+ CALL REDGET(ITYP,NITMA,FLOT,TEXT,DFLOT)
+ IF(ITYP.NE.3)CALL XABORT('@RESGEO: CHARACTER DATA EXPECTED.')
+ IF(TEXT.NE.'NAP:') CALL XABORT('@RESGEO: NAP: keyword '
+ 1 //'EXPECTED.')
+ CALL LCMOP(IPGNW,'GEONEW',0,1,0)
+ CALL LCMSIX(IPMAP,' ',0)
+ CALL LCMSIX(IPMAP,'GEOMAP',1)
+ NENT2=3
+ JENT2(1)=0
+ JENT2(2)=2
+ JENT2(3)=2
+ HENT2(1)='GEONEW'
+ HENT2(1)='GEOOLD'
+ HENT2(1)='COMPO'
+ IENT2(1)=1
+ IENT2(2)=1
+ IENT2(3)=1
+ KENT2(1)=IPGNW
+ KENT2(2)=IPMAP
+ KENT2(3)=IPCPO
+ CALL NAP(NENT2,HENT2,IENT2,JENT2,KENT2)
+ CALL LCMSIX(IPMAP,' ',0)
+ IF(IMPX.GT.3)CALL LCMLIB(IPMAP)
+ CALL LCMDEL(IPMAP,'GEOMAP')
+ IF(IMPX.GT.3)CALL LCMLIB(IPMAP)
+ CALL LCMSIX(IPMAP,'GEOMAP',1)
+ IF(IMPX.GT.3)CALL LCMLIB(IPMAP)
+ CALL LCMEQU(IPGNW,IPMAP)
+ IF(IMPX.GT.3)CALL LCMLIB(IPMAP)
+ CALL LCMCL(IPGNW,1)
+ ENDIF
+****
+ CALL LCMSIX(IPMAP,' ',0)
+ IF(IMPX.GT.3)CALL LCMLIB(IPMAP)
+ CALL LCMSIX(IPMAP,'GEOMAP',1)
+ IF(IMPX.GT.3)CALL LCMLIB(IPMAP)
+****
+ ISTATE(:NSTATE)=0
+ CALL LCMGET(IPMAP,'STATE-VECTOR',ISTATE)
+ IF(ISTATE(1).NE.IGEO) CALL XABORT('@RESGEO: THE GEOMETRY '
+ 1 // 'IN FUEL-MAP MUST HAVE THE SAME TYPE AS IN THE MATEX-OBJECT')
+ IGEO=ISTATE(1)
+ NX=ISTATE(3)
+ NY=ISTATE(4)
+ NZ=ISTATE(5)
+*----
+* READ FUEL-MAP GEOMETRY AND PERFORM MESH-SPLITTING
+*----
+ IMPX0=MAX(0,IMPX-1)
+ NX2=NX
+ NY2=NY
+ IF(IGEO.GE.8) NY2=1
+ NZ2=NZ
+ ALLOCATE(ISPLX(NX2),ISPLY(NY2),ISPLZ(NZ2))
+ ISPLTL=0
+ ISPLTH=0
+ IHEX=0
+ CALL LCMLEN(IPMAP,'SPLITL',ILEN,ITYLCM)
+ IF(ILEN.GT.0) CALL LCMGET(IPMAP,'SPLITL',ISPLTL)
+ CALL LCMLEN(IPMAP,'SPLITH',ILEN,ITYLCM)
+ IF(ILEN.GT.0) CALL LCMGET(IPMAP,'SPLITH',ISPLTH)
+ IF(IGEO.LT.8) THEN
+ CALL LCMLEN(IPMAP,'SPLITX',ILONG,ITYLCM)
+ IF(ILONG.GT.0) THEN
+ CALL LCMGET(IPMAP,'SPLITX',ISPLX)
+ NX2=0
+ DO IOLD=1,NX
+ NX2=NX2+ISPLX(IOLD)
+ ENDDO
+ ENDIF
+ CALL LCMLEN(IPMAP,'SPLITY',ILONG,ITYLCM)
+ IF(ILONG.GT.0) THEN
+ CALL LCMGET(IPMAP,'SPLITY',ISPLY)
+ NY2=0
+ DO IOLD=1,NY
+ NY2=NY2+ISPLY(IOLD)
+ ENDDO
+ ENDIF
+ ELSEIF((ISPLTH.NE.0).AND.((IGEO.EQ.8).OR.(IGEO.EQ.9))) THEN
+ NX2=NX*6*(ISPLTH**2)
+ CALL LCMGET(IPMAP,'IHEX',IHEX)
+ ELSEIF((ISPLTL.NE.0).AND.((IGEO.EQ.8).OR.(IGEO.EQ.9))) THEN
+ NX2=NX*3*(ISPLTL**2)
+ CALL LCMGET(IPMAP,'IHEX',IHEX)
+ ENDIF
+ CALL LCMLEN(IPMAP,'SPLITZ',ILONG,ITYLCM)
+ IF(ILONG.GT.0) THEN
+ CALL LCMGET(IPMAP,'SPLITZ',ISPLZ)
+ NZ2=0
+ DO IOLD=1,NZ
+ NZ2=NZ2+ISPLZ(IOLD)
+ ENDDO
+ ENDIF
+ MAXPTS=NX2*NY2*NZ2
+ MAXX=NX2
+ IF(IHEX.EQ.1) THEN
+ MAXPTS=12*MAXPTS
+ MAXX=12*MAXX
+ ELSE IF((IHEX.EQ.2).OR.(IHEX.EQ.3)) THEN
+ MAXPTS=6*MAXPTS
+ MAXX=6*MAXX
+ ELSE IF(IHEX.EQ.4) THEN
+ MAXPTS=4*MAXPTS
+ MAXX=4*MAXX
+ ELSE IF(IHEX.EQ.5) THEN
+ MAXPTS=3*MAXPTS
+ MAXX=3*MAXX
+ ELSE IF((IHEX.GE.6).AND.(IHEX.LE.8)) THEN
+ MAXPTS=2*MAXPTS
+ MAXX=2*MAXX
+ ENDIF
+ ALLOCATE(MAT(MAXPTS),GMAPX(MAXX+1),GMAPY(NY2+1),GMAPZ(NZ2+1))
+ CALL READ3D(MAXX,NY2,NZ2,MAXPTS,IPMAP,IHEX,IR,ILK,SIDE,GMAPX,
+ 1 GMAPY,GMAPZ,IMPX0,NX2,NY2,NZ2,MAT,NEL,NCODE,ICODE,ZCODE,ISPLX,
+ 2 ISPLY,ISPLZ,ISPLH,ISPLL)
+ IF((NEL.NE.NX2*NY2*NZ2).AND.(IHEX.EQ.0))CALL XABORT('@RESGEO: WR'
+ 1 // 'ONG GEOMETRY.')
+ IF((NEL.NE.NX2*NZ2).AND.(IHEX.NE.0))CALL XABORT('@RESGEO: WRONG '
+ 1 // 'HEXAGONAL GEOMETRY, WRONG NUMBER OF ELEMENTS.')
+ DEALLOCATE(MAT,ISPLZ,ISPLY,ISPLX)
+ IF(IMPX.GT.2)WRITE(IOUT,*)'CHECKING FUEL-MAP GEOMETRY'
+ IF((IGEO.NE.7).AND.(IGEO.NE.9))CALL XABORT('@RESGEO: ONLY '
+ 1 //'3D-CARTESIAN OR 3D-HEXAGONAL GEOMETRY ALLOWED.')
+ IF(IHEX.EQ.0) THEN
+ IF((LX.LT.NX).OR.(LY.LT.NY).OR.(LZ.LT.NZ)) THEN
+ WRITE(IOUT,*) 'Geometry LX=',LX,', LY=',LY,' and LZ=',LZ,
+ 1 ' must be greater or equal to map ',
+ 2 'NX=',NX,' NY=',NY,' and NZ=',NZ
+ CALL XABORT('@RESGEO: WRONG GEOMETRY DEFINITION.')
+ ENDIF
+ ELSE
+ IF((LX.LT.NX).OR.(LZ.LT.NZ)) THEN
+ WRITE(IOUT,*) 'Geometry LX=',LX,' and LZ=',LZ,
+ 1 ' must be greater or equal to map ',
+ 2 'NX=',NX,' and NZ=',NZ
+ CALL XABORT('@RESGEO: WRONG GEOMETRY DEFINITION.')
+ ENDIF
+ ENDIF
+ IF(NZ.LT.NB)THEN
+ WRITE(IOUT,*)'@RESGEO: FOUND NZ=',NZ,' LESS THAN NB=',NB
+ CALL XABORT('@RESGEO: WRONG FUEL-MAP GEOMETRY DEFINITION.')
+ ENDIF
+*----
+* CHECK MESHX OR SIDE
+*----
+ IF(IGEO.EQ.7) THEN
+ GEOXX(:LX+1)=0.0
+ CALL LCMGET(IPMTX,'MESHX',GEOXX)
+ DO 10 IMP=1,NX+1
+ DO IGM=1,LX+1
+ IF(ABS(GMAPX(IMP)-GEOXX(IGM)).LT.EPSI)THEN
+ GEOXX(IGM)=GMAPX(IMP)
+ GOTO 10
+ ENDIF
+ ENDDO
+ WRITE(IOUT,*)'@RESGEO: MESHX IN L_MAP ',GMAPX(IMP)
+ CALL XABORT('@RESGEO: UNABLE TO FIND THIS MESHX IN L_GEOM.')
+ 10 CONTINUE
+ CALL LCMPUT(IPMTX,'MESHX',LX+1,2,GEOXX)
+ ELSE IF(IGEO.EQ.9) THEN
+ ISPLTL=0
+ NY=1
+ CALL LCMGET(IPMAP,'SIDE',GMAPSI)
+ CALL LCMLEN(IPMAP,'SPLITL',ILONG,ITYLCM)
+ IF(ILONG.GT.0) CALL LCMGET(IPMAP,'SPLITL',ISPLTL)
+ IF(ISPLTL.EQ.0) ISPLTL=1
+ GMAPSI=GMAPSI/REAL(ISPLTL)
+ CALL LCMGET(IPMTX,'SIDE',GEOSI)
+ IF(ABS(GMAPSI-GEOSI).LT.EPSI)THEN
+ GEOSI=GMAPSI
+ GOTO 20
+ ENDIF
+ WRITE(IOUT,*)'@RESGEO: SIDE IN L_MAP ',GMAPSI, GEOSI
+ CALL XABORT('@RESGEO: UNABLE TO FIND THIS SIDE IN L_GEOM.')
+ 20 CONTINUE
+ CALL LCMPUT(IPMTX,'SIDE',1,2,GEOSI)
+ ENDIF
+*----
+* CHECK MESHY (ONLY IF 3D-CARTESIAN GEOMETRY)
+*----
+ IF(IGEO.EQ.7) THEN
+ GEOYY(:LY+1)=0.0
+ CALL LCMGET(IPMTX,'MESHY',GEOYY)
+ DO 30 IMP=1,NY+1
+ DO IGM=1,LY+1
+ IF(ABS(GMAPY(IMP)-GEOYY(IGM)).LT.EPSI)THEN
+ GEOYY(IGM)=GMAPY(IMP)
+ GOTO 30
+ ENDIF
+ ENDDO
+ WRITE(IOUT,*)'@RESGEO: MESHY IN FUEL MAP ',GMAPY(IMP)
+ CALL XABORT('@RESGEO: UNABLE TO FIND THIS MESHY IN L_GEOM.')
+ 30 CONTINUE
+ CALL LCMPUT(IPMTX,'MESHY',LY+1,2,GEOYY)
+ ENDIF
+*----
+* CHECK MESHZ
+*----
+ GEOZZ(:LZ+1)=0.0
+ CALL LCMGET(IPMTX,'MESHZ',GEOZZ)
+ DO 50 IMP=1,NZ+1
+ DO IGM=1,LZ+1
+ IF(ABS(GMAPZ(IMP)-GEOZZ(IGM)).LT.EPSI)THEN
+ GEOZZ(IGM)=GMAPZ(IMP)
+ GOTO 50
+ ENDIF
+ ENDDO
+ WRITE(IOUT,*)'@RESGEO: MESHZ IN FUEL MAP ',GMAPZ(IMP)
+ CALL XABORT('@RESGEO: UNABLE TO FIND THIS MESHZ IN L_GEOM.')
+ 50 CONTINUE
+ CALL LCMPUT(IPMTX,'MESHZ',LZ+1,2,GEOZZ)
+ DEALLOCATE(GMAPZ,GMAPY,GMAPX)
+*----
+* CHECK FUEL MIXTURES
+*----
+ CALL RESPFM(IPMAP,IPMTX,NX,NY,NZ,LX,LY,LZ,NFUEL,IMPX,IGEO,NCH,NB,
+ 1 NTOT)
+ RETURN
+ END
generated by cgit v1.2.3 (git 2.46.0) at 2026-02-16 14:55:00 +0000