diff options
Diffstat (limited to 'Dragon/src/XELVOL.f')
| -rw-r--r-- | Dragon/src/XELVOL.f | 228 |
1 files changed, 228 insertions, 0 deletions
diff --git a/Dragon/src/XELVOL.f b/Dragon/src/XELVOL.f new file mode 100644 index 0000000..1b1e0f4 --- /dev/null +++ b/Dragon/src/XELVOL.f @@ -0,0 +1,228 @@ +*DECK XELVOL + SUBROUTINE XELVOL( IPRT, NDIM, NEXTGE, NCPC, MINDO, MAXDO, + > ICORDO, NSURO, NVOLO, IDLGEO, INDEXO, + > MAXC, REMESH, MATGEO, VOLSO ) +* +*----------------------------------------------------------------------- +* +*Purpose: +* Compute volumes and surfaces. +* +*Copyright: +* Copyright (C) 1987 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): R. Roy +* +*Parameters: input +* IPRT intermediate printing level for output. +* NDIM number of dimensions (2 or 3). +* NEXTGE rectangular(0)/circular(1) boundary. +* NCPC dimension for MINDO. +* MINDO min index values for all axes (rect/cyl). +* MAXDO max index values for all axes (rect/cyl). +* ICORDO principal axis directions (X/Y/Z) for meshes. +* NSURO number of surface of the geometry. +* NVOLO max. number of track segments in a single track. +* IDLGEO relative index of geometry in VOLSO. +* INDEXO to retrieve zones in geometry. +* MAXC dimension of REMESH. +* REMESH real meshes values (rect/cyl). +* MATGEO material numbers used in the geometry. +* +*Parameters: output +* VOLSO volumes and surfaces for each geometry. +* +*----------------------------------------------------------------------- +* + IMPLICIT NONE +* + INTEGER IPRT,NDIM,NEXTGE,NCPC,NSURO,NVOLO,IDLGEO,MAXC + REAL VOLSO(*),REMESH(MAXC) + INTEGER MINDO(NCPC), MAXDO(NCPC), ICORDO(NCPC), + > INDEXO(4,*), MATGEO(*) +* + DOUBLE PRECISION PI, PIO2 + PARAMETER ( PI=3.14159265358979323846D0,PIO2= 0.5D0*PI) + INTEGER IOUT + PARAMETER ( IOUT=6 ) + INTEGER ICUR(4), IND, I, KSUR, JX, JY, JZ, JCY, JRAY, + > IVO, IXYZ, ICT, ICX, ICY, IVO1, IVO2, IP, IR, + > NESUR, NSUNN, NSURC, NSURM, + > NEVOL, NVONN, NVOLC, NVOLM + INTEGER NVOL0 + DOUBLE PRECISION CENTEC(2),RAYONC(2),DXX(2),DYY(2),DZZ, + > XYPOS(2,2),DELTA(3),PRODUC,VOLARE,SURF2,AREAI + LOGICAL LELCRN, SWZCYL + CHARACTER*4 CORIEN(-6:0) + SAVE CORIEN +* + DATA CORIEN + > / ' Z+ ',' Z- ',' Y+ ',' Y- ',' X+ ',' X- ',' ' / +* + IND(I)= IDLGEO + I +* +* VOL & SURF CALCULATION (CARTESIAN MESHES) + RAYONC(1)= 0.0D0 + KSUR= MOD(NDIM+1,3) + DO 50 JX = MINDO(1)-1, MAXDO(1) + ICUR(1)= JX + IF( JX.EQ.MINDO(1)-1 .OR. JX.EQ.MAXDO(1) )THEN + DELTA(1) = 0.25D0 + ELSE + DELTA(1)= DBLE(REMESH(JX+1))- DBLE(REMESH(JX)) + ENDIF + DO 40 JY = MINDO(2)-1, MAXDO(2) + ICUR(2)= JY + IF( JY.EQ.MINDO(2)-1 .OR. JY.EQ.MAXDO(2) )THEN + DELTA(2) = 0.25D0 + ELSE + DELTA(2)= DBLE(REMESH(JY+1))-DBLE(REMESH(JY)) + ENDIF + DO 30 JZ = MINDO(3)-KSUR, MAXDO(3)+KSUR-1 + ICUR(3)= JZ + IF( JZ.EQ.MINDO(3)-1 .OR. JZ.EQ.MAXDO(3) )THEN + DELTA(3) = 0.25D0 + ELSE + DELTA(3)= DBLE(REMESH(JZ+1))-DBLE(REMESH(JZ)) + ENDIF + PRODUC= DELTA(1)*DELTA(2)*DELTA(3) + DO 20 IVO = NSURO, NVOLO + DO 10 IXYZ= 1, 3 + IF(INDEXO(IXYZ,IND(IVO)).NE.ICUR(IXYZ))GO TO 20 + 10 CONTINUE + VOLSO(IND(IVO))= REAL(PRODUC) + 20 CONTINUE + 30 CONTINUE + 40 CONTINUE + 50 CONTINUE + NVOL0=0 + DO 130 JCY= 4, NCPC + ICT= ICORDO(JCY) + ICX= MOD(ICT , 3) + 1 + ICY= MOD(ICT+1, 3) + 1 + CENTEC(1)= DBLE(REMESH(MINDO(JCY)-2)) + CENTEC(2)= DBLE(REMESH(MINDO(JCY)-1)) + DO 120 JRAY= MAXDO(JCY), MINDO(JCY), -1 + RAYONC(2)= DBLE(REMESH(JRAY)) + DO 110 JX = MINDO(ICX), MAXDO(ICX)-1 + ICUR(ICX)= JX + DXX(1) = DBLE(REMESH(JX)) + DXX(2) = DBLE(REMESH(JX+1)) + XYPOS(1,1) = DXX(1)-CENTEC(1) + XYPOS(2,1) = DXX(2)-CENTEC(1) + DO 100 JY = MINDO(ICY), MAXDO(ICY)-1 + ICUR(ICY)= JY + DYY(1) = DBLE(REMESH(JY)) + DYY(2) = DBLE(REMESH(JY+1)) + XYPOS(1,2) = DYY(1)-CENTEC(2) + XYPOS(2,2) = DYY(2)-CENTEC(2) + IF( .NOT.LELCRN(CENTEC,RAYONC,DXX,DYY )) + > GO TO 100 + CALL XELCRN(IPRT,RAYONC(2),1,1,XYPOS,AREAI) + DO 90 JZ = MINDO(ICT)-KSUR, MAXDO(ICT)+KSUR-1 + ICUR(ICT)= JZ + IF( JZ.EQ.MINDO(ICT)-1 .OR. JZ.EQ.MAXDO(ICT) )THEN + DZZ = 0.25D0 + ELSE + DZZ = DBLE(REMESH(JZ+1)) - DBLE(REMESH(JZ)) + ENDIF + VOLARE= AREAI * DZZ + DO 85 IVO1= NSURO, NVOLO + ICUR(4)= JRAY-1 + DO 70 IXYZ= 1, 4 + IF(INDEXO(IXYZ,IND(IVO1)).NE.ICUR(IXYZ)) GO TO 85 + 70 CONTINUE + ICUR(4)= JRAY + DO 80 IVO2= NSURO, NVOLO + DO 75 IXYZ= 1, 4 + IF(INDEXO(IXYZ,IND(IVO2)).NE.ICUR(IXYZ)) GO TO 80 + 75 CONTINUE + IF(VOLARE .LE. 0.0D0) THEN + IF(NVOL0 .EQ. 0) WRITE(IOUT,8000) + NVOL0=NVOL0+1 + VOLARE=0.0D0 + ELSE IF(VOLSO(IND(IVO2))-VOLARE .LE. 0.0) THEN + IF(NVOL0 .EQ. 0) WRITE(IOUT,8000) + NVOL0=NVOL0+1 + VOLARE=DBLE(VOLSO(IND(IVO2))) + ENDIF + VOLSO(IND(IVO1))= REAL(VOLARE) + VOLSO(IND(IVO2))= VOLSO(IND(IVO2)) - REAL(VOLARE) + GO TO 85 + 80 CONTINUE + 85 CONTINUE + 90 CONTINUE + 100 CONTINUE + 110 CONTINUE + 120 CONTINUE + 130 CONTINUE + NESUR= 0 + NEVOL= 0 + SWZCYL=.TRUE. + IF( NEXTGE.EQ.1 )THEN + DO 700 IVO1= NSURO, NVOLO + IF( IVO1.LT.0.AND.MATGEO(IND(IVO1)).EQ.0 ) NESUR= NESUR-1 + IF( IVO1.GT.0.AND.MATGEO(IND(IVO1)).LT.0 ) NEVOL= NEVOL+1 + 700 CONTINUE + JRAY = MAXDO(4) + SURF2= PIO2 * SQRT(DBLE(REMESH(JRAY))) + ICT= ICORDO(4) + DO 720 JZ = MINDO(ICT), MAXDO(ICT)-1 + DZZ = DBLE(REMESH(JZ+1)) - DBLE(REMESH(JZ)) + DO 710 IVO1= 1, NVOLO + IF( INDEXO(ICT,IND(IVO1)).NE.JZ ) GO TO 710 + IF( INDEXO( 4,IND(IVO1)).NE.JRAY ) GO TO 710 + SWZCYL= SWZCYL.AND.(MATGEO(IND(IVO1)).LT.0) + VOLSO(IND(MATGEO(IND(IVO1))))= REAL(DZZ*SURF2) + 710 CONTINUE + 720 CONTINUE + ENDIF +* + IF( IPRT.GT.1 )THEN + NSUNN = NSURO-NESUR-NEVOL + NSURC = -1 + WRITE(IOUT,'(1H )') + WRITE(IOUT,'(/35H CELL SURFACES BEFORE ASSEMBLING )') + DO 180 IP = 1, (9 - NSUNN) / 10 + NSURM= MAX( NSUNN, NSURC-9 ) + WRITE(IOUT,'(10X,10(A5,I7))') + > (' SUR ',-IR,IR= NSURC, NSURM, -1) + WRITE(IOUT,'(8H VALUE ,2X,1P,10E12.4)') + > (4.*VOLSO(IND(IR)),IR=NSURC,NSURM,-1) + WRITE(IOUT,'(8H ORIENT ,2X,10A12)') + > (CORIEN(MATGEO(IND(IR))),IR=NSURC,NSURM,-1) + WRITE(IOUT,'(1H )') + NSURC = NSURC - 10 + 180 CONTINUE + NVONN= NVOLO-NEVOL + NVOLC= 1 + WRITE(IOUT,'(1H )') + WRITE(IOUT,'( 35H CELL VOLUMES BEFORE ASSEMBLING )') + DO 190 IP = 1, (9 + NVONN) / 10 + NVOLM= MIN( NVONN, NVOLC+9 ) + WRITE(IOUT,'(10X,10(A5,I7))') + > (' VOL ',IR,IR=NVOLC,NVOLM, 1) + WRITE(IOUT,'(8H VALUE ,2X,1P,10E12.4)') + > (VOLSO(IND(IR)),IR=NVOLC,NVOLM, 1) + WRITE(IOUT,'(8H MERGE ,2X,10I12)') + > (MATGEO(IND(IR)),IR=NVOLC,NVOLM, 1) + WRITE(IOUT,'(1H )') + NVOLC = NVOLC + 10 + 190 CONTINUE + IF( .NOT.SWZCYL ) + > CALL XABORT( 'XELVOL: '// + > 'CIRCULAR ZONES INCORRECTLY DEFINED' ) + ENDIF +* + RETURN +*---- +* Formats +*---- + 8000 FORMAT(1X,'****** WARNING IN XELVOL ******'/ + > 1X,'AT LEAST ONE REGION WITH NEGATIVE VOLUME'/ + > 1X,'THIS VOLUME IS RESET TO 0.0') + END |
