From 31c08efeeec82aab8fd4417cffbe84e68c06c760 Mon Sep 17 00:00:00 2001 From: Bryan Herman Date: Wed, 17 Oct 2012 06:24:31 -0700 Subject: [PATCH] added parameter CMFD_NOACCEL that is associated with non-accelerated regions on coarse mesh overlay --- src/DEPENDENCIES | 5 ++ src/Makefile | 10 +-- src/cmfd_data.F90 | 137 ++++++++++--------------------------- src/cmfd_execute.F90 | 5 +- src/cmfd_header.F90 | 6 +- src/cmfd_loss_operator.F90 | 10 +-- src/cmfd_prod_operator.F90 | 5 +- src/constants.F90 | 6 ++ src/global.F90 | 3 + 9 files changed, 69 insertions(+), 118 deletions(-) diff --git a/src/DEPENDENCIES b/src/DEPENDENCIES index 03f739bd6..23832c6f6 100644 --- a/src/DEPENDENCIES +++ b/src/DEPENDENCIES @@ -13,6 +13,7 @@ ace.o: string.o ace_header.o: constants.o ace_header.o: endf_header.o +cmfd_data.o: constants.o cmfd_data.o: global.o cmfd_data.o: mesh.o cmfd_data.o: mesh_header.o @@ -25,6 +26,8 @@ cmfd_execute.o: cmfd_snes_solver.o cmfd_execute.o: tally.o cmfd_execute.o: timing.o +cmfd_header.o: constants.o + cmfd_input.o: datatypes.o cmfd_input.o: error.o cmfd_input.o: global.o @@ -38,6 +41,7 @@ cmfd_jacobian_operator.o: global.o cmfd_jacobian_operator.o: cmfd_loss_operator.o cmfd_jacobian_operator.o: cmfd_prod_operator.o +cmfd_loss_operator.o: constants.o cmfd_loss_operator.o: global.o cmfd_message_passing.o: cmfd_header.o @@ -53,6 +57,7 @@ cmfd_output.o: timing.o cmfd_power_solver.o: cmfd_loss_operator.o cmfd_power_solver.o: cmfd_prod_operator.o +cmfd_prod_operator.o: constants.o cmfd_prod_operator.o: global.o cmfd_snes_solver.o: cmfd_jacobian_operator.o diff --git a/src/Makefile b/src/Makefile index 72694ad7f..b61e5c855 100644 --- a/src/Makefile +++ b/src/Makefile @@ -15,13 +15,13 @@ include OBJECTS # User Options #=============================================================================== -COMPILER = gnu +COMPILER = intel DEBUG = no PROFILE = no -OPTIMIZE = no -MPI = no -HDF5 = no -PETSC = no +OPTIMIZE = yes +MPI = yes +HDF5 = yes +PETSC = yes #=============================================================================== # External Library Paths diff --git a/src/cmfd_data.F90 b/src/cmfd_data.F90 index b179f57fd..2b38db761 100644 --- a/src/cmfd_data.F90 +++ b/src/cmfd_data.F90 @@ -20,20 +20,21 @@ contains subroutine set_up_cmfd() - use global, only: cmfd, cmfd_coremap, cmfd_run_2grp use cmfd_header, only: allocate_cmfd + use constants, only: CMFD_NOACCEL + use global, only: cmfd, cmfd_coremap, cmfd_run_2grp ! initialize cmfd object if (.not.allocated(cmfd%flux)) call allocate_cmfd(cmfd) ! check for core map and set it up - if ((cmfd_coremap) .and. (cmfd%mat_dim == 9999)) call set_coremap() + if ((cmfd_coremap) .and. (cmfd%mat_dim == CMFD_NOACCEL)) call set_coremap() ! calculate all cross sections based on reaction rates from last batch if (.not. cmfd_run_2grp) call compute_xs() ! check neutron balance - call neutron_balance(670) +! call neutron_balance(670) ! fix 2 grp cross sections if (cmfd_run_2grp) call fix_2_grp() @@ -55,7 +56,7 @@ contains use constants, only: FILTER_MESH, FILTER_ENERGYIN, FILTER_ENERGYOUT& , SURF_FILTER_ENERGYIN, SURF_FILTER_SURFACE, & IN_RIGHT, OUT_RIGHT, IN_FRONT, OUT_FRONT, & - IN_TOP, OUT_TOP, N_FILTER_TYPES + IN_TOP, OUT_TOP, N_FILTER_TYPES, CMFD_NOACCEL use global, only: cmfd, message, n_user_tallies, n_tallies, & tallies, meshes use error, only: fatal_error @@ -117,7 +118,7 @@ contains ! check for active mesh cell if (allocated(cmfd%coremap)) then - if (cmfd%coremap(i,j,k) == 99999) then + if (cmfd%coremap(i,j,k) == CMFD_NOACCEL) then cycle end if end if @@ -293,7 +294,8 @@ contains subroutine set_coremap() - use global, only: cmfd + use constants, only: CMFD_NOACCEL + use global, only: cmfd integer :: kount=1 ! counter for unique fuel assemblies integer :: nx ! number of mesh cells in x direction @@ -326,7 +328,7 @@ contains if (cmfd % coremap(i,j,k) == 1) then ! reset value to 99999 - cmfd % coremap(i,j,k) = 99999 + cmfd % coremap(i,j,k) = CMFD_NOACCEL else @@ -353,7 +355,7 @@ contains subroutine neutron_balance(uid) - use constants, only: ONE + use constants, only: ONE, CMFD_NOACCEL use global, only: cmfd, keff integer :: nx ! number of mesh cells in x direction @@ -396,7 +398,7 @@ contains ! check for active mesh if (allocated(cmfd%coremap)) then - if (cmfd%coremap(i,j,k) == 99999) then + if (cmfd%coremap(i,j,k) == CMFD_NOACCEL) then cmfd%resnb(g,i,j,k) = 99999.0 cycle end if @@ -439,14 +441,15 @@ contains cmfd%resnb(g,i,j,k) = res ! write out info to file - write(uid,*) 'Location',i,j,k,' Group:',g,'Balance:',res - write(uid,*) 'Leakage:',leakage - write(uid,*) 'Interactions:',interactions - write(uid,*) 'Scattering:',scattering - write(uid,*) 'Fission:',fission - write(uid,*) 'k-eff:',keff - write(uid,*) 'k-eff balanced:',fission/(leakage+interactions-scattering) - write(uid,*) 'Balance:',keff - fission/(leakage+interactions-scattering) + write(uid,'(A,1X,I0,1X,I0,1X,I0,1X,A,1X,I0,1X,A,1PE11.4)') & + 'Location',i,j,k,' Group:',g,'Balance:',res + write(uid,100) 'Leakage:',leakage + write(uid,100) 'Interactions:',interactions + write(uid,100) 'Scattering:',scattering + write(uid,100) 'Fission:',fission + write(uid,100) 'k-eff:',keff + write(uid,100) 'k-eff balanced:',fission/(leakage+interactions-scattering) + write(uid,100) 'Balance:',keff - fission/(leakage+interactions-scattering) end do GROUPG @@ -456,6 +459,8 @@ contains end do ZLOOP + 100 FORMAT(A,1X,1PE11.4) + end subroutine neutron_balance !=============================================================================== @@ -464,7 +469,8 @@ contains subroutine compute_dtilde() - use global, only: cmfd, cmfd_coremap + use constants, only: CMFD_NOACCEL + use global, only: cmfd, cmfd_coremap integer :: nx ! maximum number of cells in x direction integer :: ny ! maximum number of cells in y direction @@ -514,7 +520,7 @@ contains ! check for active mesh cell if (allocated(cmfd%coremap)) then - if (cmfd%coremap(i,j,k) == 99999) cycle + if (cmfd%coremap(i,j,k) == CMFD_NOACCEL) cycle end if ! get cell data @@ -604,7 +610,8 @@ contains subroutine compute_dhat() - use global, only: cmfd, cmfd_coremap + use constants, only: CMFD_NOACCEL + use global, only: cmfd, cmfd_coremap integer :: nx ! maximum number of cells in x direction integer :: ny ! maximum number of cells in y direction @@ -650,7 +657,7 @@ contains ! check for active mesh cell if (allocated(cmfd%coremap)) then - if (cmfd%coremap(i,j,k) == 99999) then + if (cmfd%coremap(i,j,k) == CMFD_NOACCEL) then cycle end if end if @@ -696,7 +703,7 @@ contains if (cmfd_coremap) then if (cmfd % coremap(neig_idx(1),neig_idx(2),neig_idx(3)) == & - 99999 .and. cmfd % coremap(i,j,k) /= 99999) then + CMFD_NOACCEL .and. cmfd % coremap(i,j,k) /= CMFD_NOACCEL) then ! compute dhat dhat = (net_current - shift_idx*cell_dtilde(l)*cell_flux) /& @@ -813,15 +820,15 @@ contains dhat_reset = .true. end subroutine fix_2_grp -# ifdef HIDE + !=============================================================================== ! FIX_NEUTRON_BALANCE !=============================================================================== subroutine fix_neutron_balance() - use constants, only: ONE, ZERO - use global, only: cmfd, cmfd_balance, keff + use constants, only: ONE, ZERO, CMFD_NOACCEL + use global, only: cmfd, cmfd_balance, keff integer :: nx ! number of mesh cells in x direction integer :: ny ! number of mesh cells in y direction @@ -870,7 +877,7 @@ contains ! check for active mesh if (allocated(cmfd%coremap)) then - if (cmfd%coremap(i,j,k) == 99999) cycle + if (cmfd%coremap(i,j,k) == CMFD_NOACCEL) cycle end if ! compute leakage in groups 1 and 2 @@ -937,26 +944,6 @@ contains ! zero out upscatter cross section cmfd % scattxs(2,1,i,j,k) = ZERO -! print *,'Leakage g=1:',leak1 -! print *,'Leakage g=2:',leak2 -! print *,'Flux g=1:',flux1 -! print *,'Flux g=2:',flux2 -! print *,'Total g=1:',sigt1 -! print *,'Total g=2:',sigt2 -! print *,'Scattering 1-->1:',sigs11 -! print *,'Scattering 2-->1:',sigs21 -! print *,'Scattering 1-->2:',sigs12 -! print *,'Scattering 2-->2:',sigs22 -! print *,'Fission 1-->1:',nsigf11 -! print *,'Fission 2-->1:',nsigf21 -! print *,'Fission 1-->2:',nsigf12 -! print *,'Fission 2-->2:',nsigf22 -! print *,'keff:',keff -! print *,'Effective Downscatter:',sigs12_eff -! -! print *,'Group 1 balance:',leak1+sigt1*flux1-sigs11*flux1-ONE/keff*(nsigf11*flux1+nsigf21*flux2) -! print *,'Group 2 balance:',leak2+sigt2*flux2-sigs12_eff*flux1-sigs22*flux2 -! stop end do XLOOP end do YLOOP @@ -971,8 +958,8 @@ contains subroutine compute_effective_downscatter() - use constants, only: ZERO - use global, only: cmfd, cmfd_downscatter, keff + use constants, only: ZERO, CMFD_NOACCEL + use global, only: cmfd, cmfd_downscatter, keff integer :: nx ! number of mesh cells in x direction integer :: ny ! number of mesh cells in y direction @@ -1014,7 +1001,7 @@ contains ! check for active mesh if (allocated(cmfd%coremap)) then - if (cmfd%coremap(i,j,k) == 99999) cycle + if (cmfd%coremap(i,j,k) == CMFD_NOACCEL) cycle end if ! extract cross sections and flux from object @@ -1056,56 +1043,4 @@ contains end subroutine compute_effective_downscatter -!============================================================================== -! EXTRACT_ACCUM_FSRC -!============================================================================== -! -! subroutine extract_accum_fsrc(i,j,k) -! -! use global -! use mesh, only: mesh_indices_to_bin -! use mesh_header, only: StructuredMesh -! use tally_header, only: TallyObject, TallyScore -! -! integer :: nx ! number of mesh cells in x direction -! integer :: ny ! number of mesh cells in y direction -! integer :: nz ! number of mesh cells in z direction -! integer :: ng ! number of energy groups -! integer :: i ! iteration counter for x -! integer :: j ! iteration counter for y -! integer :: k ! iteration counter for z -! integer :: g ! iteration counter for g -! integer :: h ! iteration counter for outgoing groups -! integer :: ital ! tally object index -! integer :: ijk(3) ! indices for mesh cell -! integer :: score_index ! index to pull from tally object -! integer :: bins(N_FILTER_TYPES) ! bins for filters -! -! real(8) :: flux ! temp variable for flux -! -! type(TallyObject), pointer :: t ! pointer for tally object -! type(StructuredMesh), pointer :: m ! pointer for mesh object -! -! ! associate tallies and mesh -! t => tallies(1) -! m => meshes(t % mesh) -! -! ! reset all bins to 1 -! bins = 1 -! -! ! set ijk as mesh indices -! ijk = (/ i, j, k /) -! -! ! get bin number for mesh indices -! bins(FILTER_MESH) = mesh_indices_to_bin(m,ijk) -! -! ! calculate score index from bins -! score_index = sum((bins - 1) * t%stride) + 1 -! -! ! bank source -! cmfd % openmc_accum_src(i,j,k) = t % scores(1,score_index) % sum -! -! end subroutine extract_accum_fsrc -# endif - end module cmfd_data diff --git a/src/cmfd_execute.F90 b/src/cmfd_execute.F90 index 844fb003c..10bb12bea 100644 --- a/src/cmfd_execute.F90 +++ b/src/cmfd_execute.F90 @@ -198,7 +198,8 @@ contains subroutine calc_fission_source() - use global, only: cmfd, cmfd_coremap, master, mpi_err, entropy_on + use constants, only: CMFD_NOACCEL + use global, only: cmfd, cmfd_coremap, master, mpi_err, entropy_on integer :: nx ! maximum number of cells in x direction integer :: ny ! maximum number of cells in y direction @@ -241,7 +242,7 @@ contains ! check for core map if (cmfd_coremap) then - if (cmfd%coremap(i,j,k) == 99999) then + if (cmfd%coremap(i,j,k) == CMFD_NOACCEL) then cycle end if end if diff --git a/src/cmfd_header.F90 b/src/cmfd_header.F90 index be6dc830f..ed2253077 100644 --- a/src/cmfd_header.F90 +++ b/src/cmfd_header.F90 @@ -1,8 +1,6 @@ module cmfd_header -#ifdef HDF5 - use hdf5 -#endif + use constants, only: CMFD_NOACCEL implicit none private @@ -19,7 +17,7 @@ module cmfd_header ! core overlay map integer, allocatable :: coremap(:,:,:) integer, allocatable :: indexmap(:,:) - integer :: mat_dim = 9999 + integer :: mat_dim = CMFD_NOACCEL ! energy grid real(8), allocatable :: egrid(:) diff --git a/src/cmfd_loss_operator.F90 b/src/cmfd_loss_operator.F90 index 6f6c8b642..74f392c8f 100644 --- a/src/cmfd_loss_operator.F90 +++ b/src/cmfd_loss_operator.F90 @@ -84,7 +84,8 @@ contains subroutine preallocate_loss_matrix(this) - use global, only: cmfd, cmfd_coremap + use constants, only: CMFD_NOACCEL + use global, only: cmfd, cmfd_coremap type(loss_operator) :: this @@ -182,7 +183,7 @@ contains ! check for neighbor that is non-acceleartred if (cmfd % coremap(neig_idx(1),neig_idx(2),neig_idx(3)) /= & - 99999) then + CMFD_NOACCEL) then ! get neighbor matrix index call indices_to_matrix(g,neig_idx(1), neig_idx(2), & @@ -247,7 +248,8 @@ contains subroutine build_loss_matrix(this, adjoint) - use global, only: cmfd, cmfd_coremap, cmfd_write_matrices + use constants, only: CMFD_NOACCEL + use global, only: cmfd, cmfd_coremap, cmfd_write_matrices type(loss_operator) :: this logical, optional :: adjoint ! set up the adjoint @@ -338,7 +340,7 @@ contains ! check that neighbor is not reflector if (cmfd % coremap(neig_idx(1),neig_idx(2),neig_idx(3)) /= & - 99999) then + CMFD_NOACCEL) then ! compute leakage coefficient for neighbor jn = -dtilde(l) + shift_idx*dhat(l) diff --git a/src/cmfd_prod_operator.F90 b/src/cmfd_prod_operator.F90 index 734781efd..d350b8aed 100644 --- a/src/cmfd_prod_operator.F90 +++ b/src/cmfd_prod_operator.F90 @@ -179,7 +179,8 @@ contains subroutine build_prod_matrix(this, adjoint) - use global, only: cmfd, cmfd_coremap, cmfd_write_matrices + use constants, only: CMFD_NOACCEL + use global, only: cmfd, cmfd_coremap, cmfd_write_matrices type(prod_operator) :: this logical, optional :: adjoint @@ -217,7 +218,7 @@ contains if (cmfd_coremap) then ! check if at a reflector - if (cmfd % coremap(i,j,k) == 99999) then + if (cmfd % coremap(i,j,k) == CMFD_NOACCEL) then cycle end if diff --git a/src/constants.F90 b/src/constants.F90 index 3af95db7a..92022dad5 100644 --- a/src/constants.F90 +++ b/src/constants.F90 @@ -365,4 +365,10 @@ module constants integer, parameter :: UNIT_STATE = 16 ! unit # for writing state point integer, parameter :: CMFD_BALANCE = 17 ! unit # for writing cmfd balance file + !============================================================================= + ! CMFD CONSTANTS + + ! for non-accelerated regions on coarse mesh overlay + integer, parameter :: CMFD_NOACCEL = 99999 + end module constants diff --git a/src/global.F90 b/src/global.F90 index 13962383f..58c6cbee3 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -342,6 +342,9 @@ module global ! batch to last flush before active batches integer :: cmfd_act_flush = 0 + ! compute effective downscatter cross section + logical :: cmfd_downscatter = .false. + ! convergence monitoring logical :: cmfd_snes_monitor = .false. logical :: cmfd_ksp_monitor = .false.