merged Paul's changes on tally handling, changed cmfd set tally and read tally routines appropriately

This commit is contained in:
Bryan Herman 2012-01-29 13:26:20 -05:00
parent 254aa57ef9
commit 6fb307a803
4 changed files with 90 additions and 43 deletions

View file

@ -27,6 +27,7 @@ cmfd_execute.o: cmfd_snes_solver.o
cmfd_input.o: datatypes.o
cmfd_input.o: error.o
cmfd_input.o: global.o
cmfd_input.o: hdf5_interface.o
cmfd_input.o: mesh_header.o
cmfd_input.o: string.o
cmfd_input.o: tally_header.o

View file

@ -28,7 +28,7 @@ USE_HDF5 = yes
ifeq ($(COMPILER),gnu)
F90 = gfortran
F90FLAGS := -cpp -fbacktrace
F90FLAGS := -cpp -fbacktrace -DN2008
LDFLAGS =
# Debugging

View file

@ -65,7 +65,7 @@ contains
integer :: h ! iteration counter for outgoing groups
integer :: ijk(3) ! indices for mesh cell
integer :: score_index ! index to pull from tally object
integer :: bins(TALLY_TYPES) ! bins for filters
integer :: bins(N_FILTER_TYPES) ! bins for filters
real(8) :: flux ! temp variable for flux
@ -103,10 +103,10 @@ contains
ijk = (/ i, j, k /)
! get bin number for mesh indices
bins(T_MESH) = mesh_indices_to_bin(m,ijk)
bins(FILTER_MESH) = mesh_indices_to_bin(m,ijk)
! apply energy in filter
bins(T_ENERGYIN) = ng - h + 1
bins(FILTER_ENERGYIN) = ng - h + 1
! calculate score index from bins
score_index = sum((bins - 1) * t%stride) + 1
@ -148,7 +148,7 @@ contains
t => tallies(2)
! set energy out bin
bins(T_ENERGYOUT) = ng - g + 1
bins(FILTER_ENERGYOUT) = ng - g + 1
! calculate score index from bins
score_index = sum((bins - 1) * t%stride) + 1
@ -166,59 +166,59 @@ contains
! initialize and filter for energy
bins = 1
bins(TS_ENERGYIN) = ng - h + 1
bins(SURF_FILTER_ENERGYIN) = ng - h + 1
! left surface
bins(1:3) = (/ i-1, j, k /) + 1
bins(TS_SURFACE) = IN_RIGHT
bins(SURF_FILTER_SURFACE) = IN_RIGHT
score_index = sum((bins - 1) * t % stride) + 1 ! outgoing
cmfd % current(1,h,i,j,k) = t % scores(score_index,1) % val
bins(TS_SURFACE) = OUT_RIGHT
bins(SURF_FILTER_SURFACE) = OUT_RIGHT
score_index = sum((bins - 1) * t % stride) + 1 ! incoming
cmfd % current(2,h,i,j,k) = t % scores(score_index,1) % val
! right surface
bins(1:3) = (/ i, j, k /) + 1
bins(TS_SURFACE) = IN_RIGHT
bins(SURF_FILTER_SURFACE) = IN_RIGHT
score_index = sum((bins - 1) * t % stride) + 1 ! incoming
cmfd % current(3,h,i,j,k) = t % scores(score_index,1) % val
bins(TS_SURFACE) = OUT_RIGHT
bins(SURF_FILTER_SURFACE) = OUT_RIGHT
score_index = sum((bins - 1) * t % stride) + 1 ! outgoing
cmfd % current(4,h,i,j,k) = t % scores(score_index,1) % val
! back surface
bins(1:3) = (/ i, j-1, k /) + 1
bins(TS_SURFACE) = IN_FRONT
bins(SURF_FILTER_SURFACE) = IN_FRONT
score_index = sum((bins - 1) * t % stride) + 1 ! outgoing
cmfd % current(5,h,i,j,k) = t % scores(score_index,1) % val
bins(TS_SURFACE) = OUT_FRONT
bins(SURF_FILTER_SURFACE) = OUT_FRONT
score_index = sum((bins - 1) * t % stride) + 1 ! incoming
cmfd % current(6,h,i,j,k) = t % scores(score_index,1) % val
! front surface
bins(1:3) = (/ i, j, k /) + 1
bins(TS_SURFACE) = IN_FRONT
bins(SURF_FILTER_SURFACE) = IN_FRONT
score_index = sum((bins - 1) * t % stride) + 1 ! incoming
cmfd % current(7,h,i,j,k) = t % scores(score_index,1) % val
bins(TS_SURFACE) = OUT_FRONT
bins(SURF_FILTER_SURFACE) = OUT_FRONT
score_index = sum((bins - 1) * t % stride) + 1 ! outgoing
cmfd % current(8,h,i,j,k) = t % scores(score_index,1) % val
! bottom surface
bins(1:3) = (/ i, j, k-1 /) + 1
bins(TS_SURFACE) = IN_TOP
bins(SURF_FILTER_SURFACE) = IN_TOP
score_index = sum((bins - 1) * t % stride) + 1 ! outgoing
cmfd % current(9,h,i,j,k) = t % scores(score_index,1) % val
bins(TS_SURFACE) = OUT_TOP
bins(SURF_FILTER_SURFACE) = OUT_TOP
score_index = sum((bins - 1) * t % stride) + 1 ! incoming
cmfd % current(10,h,i,j,k) = t % scores(score_index,1) % val
! top surface
bins(1:3) = (/ i, j, k /) + 1
bins(TS_SURFACE) = IN_TOP
bins(SURF_FILTER_SURFACE) = IN_TOP
score_index = sum((bins - 1) * t % stride) + 1 ! incoming
cmfd % current(11,h,i,j,k) = t % scores(score_index,1) % val
bins(TS_SURFACE) = OUT_TOP
bins(SURF_FILTER_SURFACE) = OUT_TOP
score_index = sum((bins - 1) * t % stride) + 1 ! outgoing
cmfd % current(12,h,i,j,k) = t % scores(score_index,1) % val

View file

@ -87,6 +87,8 @@ contains
integer :: n ! size of arrays in mesh specification
integer :: ng=1 ! number of energy groups (default 1)
integer :: n_words ! number of words read
integer :: n_filters ! number of filters
integer :: filters(N_FILTER_TYPES) ! temp list of filters
character(MAX_LINE_LEN) :: filename
character(MAX_WORD_LEN) :: words(MAX_WORDS)
type(TallyObject), pointer :: t => null()
@ -152,14 +154,20 @@ contains
! begin loop around tallies
do i = 1,n_tallies
! set n filters to 0
n_filters = 0
filters = 0
! point t to tally variable
t => tallies(i)
! allocate arrays for number of bins and stride in scores array
allocate(t % n_bins(TALLY_TYPES))
allocate(t % stride(TALLY_TYPES))
allocate(t % n_filter_bins(N_FILTER_TYPES))
allocate(t % stride(N_FILTER_TYPES))
! initialize number of bins and stride
t % n_bins = 0
t % n_filter_bins = 0
t % stride = 0
! record tally id which is equivalent to loop number
@ -168,7 +176,10 @@ contains
! set mesh filter mesh id = 1
t % mesh = 1
m => meshes(1)
t % n_bins(T_MESH) = t % n_bins(T_MESH) + product(m % dimension)
t % n_filter_bins(FILTER_MESH) = t % n_filter_bins(FILTER_MESH) + &
& product(m % dimension)
n_filters = n_filters + 1
filters(n_filters) = FILTER_MESH
! read and set incoming energy mesh filter
if (len_trim(mesh_ % energy) > 0) then
@ -178,22 +189,41 @@ contains
do j = 1,n_words
t % energy_in(j) = str_to_real(words(j))
end do
t % n_bins(T_ENERGYIN) = n_words - 1
t % n_filter_bins(FILTER_ENERGYIN) = n_words - 1
n_filters = n_filters + 1
filters(n_filters) = FILTER_ENERGYIN
end if
if (i == 1) then
! set tally estimator to analog
t % estimator = ESTIMATOR_ANALOG
! set tally type to volume
t % type = TALLY_VOLUME
! allocate and set filters
t % n_filters = n_filters
allocate(t % filters(n_filters))
t % filters = filters(1:n_filters)
! allocate macro reactions
allocate(t % macro_bins(3))
t % n_macro_bins = 3
allocate(t % score_bins(3))
t % n_score_bins = 3
! set macro_bins
t % macro_bins(1) % scalar = MACRO_FLUX
t % macro_bins(2) % scalar = MACRO_TOTAL
t % macro_bins(3) % scalar = MACRO_SCATTER_1
t % score_bins(1) % scalar = SCORE_FLUX
t % score_bins(2) % scalar = SCORE_TOTAL
t % score_bins(3) % scalar = SCORE_SCATTER_1
else if (i == 2) then
! set tally estimator to analog
t % estimator = ESTIMATOR_ANALOG
! set tally type to volume
t % type = TALLY_VOLUME
! read and set outgoing energy mesh filter
if (len_trim(mesh_ % energy) > 0) then
call split_string(mesh_ % energy, words, n_words)
@ -201,31 +231,47 @@ contains
do j = 1, n_words
t % energy_out(j) = str_to_real(words(j))
end do
t % n_bins(T_ENERGYOUT) = n_words - 1
t % n_filter_bins(FILTER_ENERGYOUT) = n_words - 1
n_filters = n_filters + 1
filters(n_filters) = FILTER_ENERGYOUT
end if
! allocate and set filters
t % n_filters = n_filters
allocate(t % filters(n_filters))
t % filters = filters(1:n_filters)
! allocate macro reactions
allocate(t % macro_bins(2))
t % n_macro_bins = 2
allocate(t % score_bins(2))
t % n_score_bins = 2
! set macro_bins
t % macro_bins(1) % scalar = MACRO_NU_SCATTER
t % macro_bins(2) % scalar = MACRO_NU_FISSION
t % score_bins(1) % scalar = SCORE_NU_SCATTER
t % score_bins(2) % scalar = SCORE_NU_FISSION
else if (i == 3) then
! set tally estimator to analog
t % estimator = ESTIMATOR_ANALOG
! allocate and set filters
t % n_filters = n_filters
allocate(t % filters(n_filters))
t % filters = filters(1:n_filters)
! allocate macro reactions
allocate(t % macro_bins(1))
t % n_macro_bins = 1
allocate(t % score_bins(1))
t % n_score_bins = 1
! set macro bins
t % macro_bins(1) % scalar = MACRO_CURRENT
t % surface_current = .true.
t % score_bins(1) % scalar = SCORE_CURRENT
t % type = TALLY_SURFACE_CURRENT
! since the number of bins for the mesh filter was already set
! assuming it was a flux tally, we need to adjust the number of
! bins
t % n_bins(T_MESH) = t % n_bins(T_MESH) - product(m % dimension)
t % n_filter_bins(FILTER_MESH) = t % n_filter_bins(FILTER_MESH) &
& - product(m % dimension)
! get pointer to mesh
id = t % mesh
@ -234,12 +280,12 @@ contains
! we need to increase the dimension by one since we also need
! currents coming into and out of the boundary mesh cells.
if (size(m % dimension) == 2) then
t % n_bins(T_MESH) = t % n_bins(T_MESH) + &
& product(m % dimension + 1) * 4
if (size(m % dimension) == 2) then ! these lines need changing
t % n_filter_bins(FILTER_MESH) = t % n_filter_bins(FILTER_MESH) + &
& product(m % dimension + 1) * 4
elseif (size(m % dimension) == 3) then
t % n_bins(T_MESH) = t % n_bins(T_MESH) + &
product(m % dimension + 1) * 6
t % n_filter_bins(FILTER_MESH) = t % n_filter_bins(FILTER_MESH) + &
& product(m % dimension + 1) * 6
end if
end if