diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index fc6c6630cc..033ce2a020 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -26,6 +26,7 @@ option(debug "Compile with debug flags" OFF) option(optimize "Turn on all compiler optimization flags" OFF) option(verbose "Create verbose Makefiles" OFF) option(coverage "Compile with flags" OFF) +option(cmfd_debug "Compile with cmfd debug flags" OFF) if (verbose) set(CMAKE_VERBOSE_MAKEFILE on) @@ -214,6 +215,15 @@ if(GIT_SHA1_SUCCESS EQUAL 0) add_definitions(-DGIT_SHA1="${GIT_SHA1}") endif() +#=============================================================================== +# CMFD debug flags +#=============================================================================== + +if (cmfd_debug) + message("-- CMFD Debug flags on") + add_definitions(-DCMFD_DEBUG) +endif() + #=============================================================================== # FoX Fortran XML Library #=============================================================================== diff --git a/src/cmfd_data.F90 b/src/cmfd_data.F90 index b5612dcd05..7962760cbf 100644 --- a/src/cmfd_data.F90 +++ b/src/cmfd_data.F90 @@ -56,7 +56,7 @@ contains ONE, TINY_BIT use error, only: fatal_error use global, only: cmfd, message, n_cmfd_tallies, cmfd_tallies, meshes,& - matching_bins + matching_bins, current_batch use mesh, only: mesh_indices_to_bin use mesh_header, only: StructuredMesh use string, only: to_str @@ -71,6 +71,7 @@ contains integer :: k ! iteration counter for z integer :: g ! iteration counter for g integer :: h ! iteration counter for outgoing groups + integer :: l ! iteration counter for leakage debug integer :: ital ! tally object index integer :: ijk(3) ! indices for mesh cell integer :: score_index ! index to pull from tally object @@ -80,6 +81,8 @@ contains integer :: i_filter_eout ! index for outgoing energy filter integer :: i_filter_surf ! index for surface filter real(8) :: flux ! temp variable for flux + real(8) :: leak1 ! group 1 leakage + real(8) :: leak2 ! group 2 leakage type(TallyObject), pointer :: t => null() ! pointer for tally object type(StructuredMesh), pointer :: m => null() ! pointer for mesh object @@ -301,6 +304,91 @@ contains if (associated(t)) nullify(t) if (associated(m)) nullify(m) +#ifdef CMFD_DEBUG + open(file='cmfd_src_' // trim(to_str(current_batch)) // '.dat', unit=100) + open(file='source_bank_' // trim(to_str(current_batch)) // '.dat', unit=101) + open(file='openmc_src_' // trim(to_str(current_batch)) // '.dat', unit=102) + open(file='totalxs1_' // trim(to_str(current_batch)) // '.dat', unit=103) + open(file='totalxs2_' // trim(to_str(current_batch)) // '.dat', unit=104) + open(file='p1scattxs1_' // trim(to_str(current_batch)) // '.dat', unit=105) + open(file='p1scattxs2_' // trim(to_str(current_batch)) // '.dat', unit=106) + open(file='scattxs11_' // trim(to_str(current_batch)) // '.dat', unit=107) + open(file='scattxs12_' // trim(to_str(current_batch)) // '.dat', unit=108) + open(file='scattxs21_' // trim(to_str(current_batch)) // '.dat', unit=109) + open(file='scattxs22_' // trim(to_str(current_batch)) // '.dat', unit=110) + open(file='nufissxs11_' // trim(to_str(current_batch)) // '.dat', unit=111) + open(file='nufissxs12_' // trim(to_str(current_batch)) // '.dat', unit=112) + open(file='nufissxs21_' // trim(to_str(current_batch)) // '.dat', unit=113) + open(file='nufissxs22_' // trim(to_str(current_batch)) // '.dat', unit=114) + open(file='diff_coef1_' // trim(to_str(current_batch)) // '.dat', unit=115) + open(file='diff_coef2_' // trim(to_str(current_batch)) // '.dat', unit=116) + open(file='flux1_' // trim(to_str(current_batch)) // '.dat', unit=117) + open(file='flux2_' // trim(to_str(current_batch)) // '.dat', unit=118) + open(file='leak1_' // trim(to_str(current_batch)) // '.dat', unit=119) + open(file='leak2_' // trim(to_str(current_batch)) // '.dat', unit=120) + + do i = 1, nx + do j = 1,ny + leak1 = ZERO + leak2 = ZERO + do l = 1, 3 + leak1 = leak1 + ((cmfd % current(4*l,1,i,j,1) - & + cmfd % current(4*l-1,1,i,j,1))) - & + ((cmfd % current(4*l-2,1,i,j,1) - & + cmfd % current(4*l-3,1,i,j,1))) + end do + do l = 1, 3 + leak2 = leak2 + ((cmfd % current(4*l,2,i,j,1) - & + cmfd % current(4*l-1,2,i,j,1))) - & + ((cmfd % current(4*l-2,2,i,j,1) - & + cmfd % current(4*l-3,2,i,j,1))) + end do + write(100,*) cmfd % cmfd_src(1,i,j,1) + write(101,*) cmfd % sourcecounts(1,i,j,1) + write(102,*) cmfd % openmc_src(1,i,j,1) + write(103,*) cmfd % totalxs(1,i,j,1) + write(104,*) cmfd % totalxs(2,i,j,1) + write(105,*) cmfd % p1scattxs(1,i,j,1) + write(106,*) cmfd % p1scattxs(2,i,j,1) + write(107,*) cmfd % scattxs(1,1,i,j,1) + write(108,*) cmfd % scattxs(1,2,i,j,1) + write(109,*) cmfd % scattxs(2,1,i,j,1) + write(110,*) cmfd % scattxs(2,2,i,j,1) + write(111,*) cmfd % nfissxs(1,1,i,j,1) + write(112,*) cmfd % nfissxs(1,2,i,j,1) + write(113,*) cmfd % nfissxs(2,1,i,j,1) + write(114,*) cmfd % nfissxs(2,2,i,j,1) + write(115,*) cmfd % diffcof(1,i,j,1) + write(116,*) cmfd % diffcof(2,i,j,1) + write(117,*) cmfd % flux(1,i,j,1) + write(118,*) cmfd % flux(2,i,j,1) + write(119,*) leak1 + write(120,*) leak2 + end do + end do + close(100) + close(101) + close(102) + close(103) + close(104) + close(105) + close(106) + close(107) + close(108) + close(109) + close(110) + close(111) + close(112) + close(113) + close(114) + close(115) + close(116) + close(117) + close(118) + close(119) + close(120) +#endif + end subroutine compute_xs !=============================================================================== diff --git a/src/cmfd_execute.F90 b/src/cmfd_execute.F90 index e4e849436f..a1b8853042 100644 --- a/src/cmfd_execute.F90 +++ b/src/cmfd_execute.F90 @@ -69,7 +69,7 @@ contains call calc_fission_source() ! calculate weight factors - if (cmfd_feedback) call cmfd_reweight(.true.) + call cmfd_reweight(.true.) ! stop cmfd timer if (master) call time_cmfd % stop() @@ -132,6 +132,7 @@ contains use constants, only: CMFD_NOACCEL, ZERO, TWO use global, only: cmfd, cmfd_coremap, master, entropy_on, current_batch + use string, only: to_str #ifdef MPI use global, only: mpi_err @@ -238,6 +239,15 @@ contains cmfd % src_cmp(current_batch) = sqrt(ONE/cmfd % norm * & sum((cmfd % cmfd_src - cmfd % openmc_src)**2)) +#ifdef CMFD_DEBUG + open(file='cmfd_src_' // trim(to_str(current_batch)) // '.dat', unit=100) + do i = 1, nx + do j = 1,ny + write(100,*) cmfd % sourcecounts(1,i,j,1) + end do + end do + close(100) +#endif end if #ifdef MPI @@ -260,6 +270,7 @@ contains use mesh_header, only: StructuredMesh use mesh, only: count_bank_sites, get_mesh_indices use search, only: binary_search + use string, only: to_str #ifdef MPI use global, only: mpi_err @@ -273,6 +284,7 @@ contains integer :: nz ! maximum number of cells in z direction integer :: ng ! maximum number of energy groups integer :: i ! iteration counter + integer :: j ! iteration counter integer :: ijk(3) ! spatial bin location integer :: e_bin ! energy bin of source particle integer :: n_groups ! number of energy groups @@ -327,8 +339,19 @@ contains cmfd % weightfactors = cmfd % cmfd_src/sum(cmfd % cmfd_src)* & sum(cmfd % sourcecounts) / cmfd % sourcecounts end where +#ifdef CMFD_DEBUG + open(file='source_bank_' // trim(to_str(current_batch)) // '.dat', unit=101) + do i = 1, nx + do j = 1,ny + write(101,*) cmfd % sourcecounts(1,i,j,1) + end do + end do + close(101) +#endif end if +if (.not. cmfd_feedback) return + ! Broadcast weight factors to all procs #ifdef MPI call MPI_BCAST(cmfd % weightfactors, ng*nx*ny*nz, MPI_REAL8, 0, &