changed current tallying so in/out partial currents are tallies separately

This commit is contained in:
Sam Shaner 2016-08-23 11:33:23 -04:00
parent bc387186ba
commit 93862cea24
8 changed files with 185 additions and 144 deletions

View file

@ -51,8 +51,9 @@ contains
subroutine compute_xs()
use constants, only: FILTER_MESH, FILTER_ENERGYIN, FILTER_ENERGYOUT, &
FILTER_SURFACE, LEFT, RIGHT, BACK, &
FRONT, BOTTOM, TOP, CMFD_NOACCEL, &
FILTER_SURFACE, OUT_LEFT, OUT_RIGHT, OUT_BACK, &
OUT_FRONT, OUT_BOTTOM, OUT_TOP, IN_LEFT, IN_RIGHT, &
IN_BACK, IN_FRONT, IN_BOTTOM, IN_TOP, CMFD_NOACCEL, &
ZERO, ONE, TINY_BIT
use error, only: fatal_error
use global, only: cmfd, n_cmfd_tallies, cmfd_tallies, meshes,&
@ -234,106 +235,74 @@ contains
matching_bins(i_filter_ein) = ng - h + 1
end if
! Left surface
! Get the bin for this mesh cell
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /))
matching_bins(i_filter_surf) = LEFT
! Left surface
matching_bins(i_filter_surf) = OUT_LEFT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t%stride) + 1 ! outgoing
* t % stride) + 1
cmfd % current(1,h,i,j,k) = t % results(1,score_index) % sum
if (i > 1) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i-1, j, k /))
matching_bins(i_filter_surf) = RIGHT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! incoming
cmfd % current(2,h,i,j,k) = t % results(1,score_index) % sum
end if
matching_bins(i_filter_surf) = IN_LEFT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(2,h,i,j,k) = t % results(1,score_index) % sum
! Right surface
if (i < nx) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i+1, j, k /) )
matching_bins(i_filter_surf) = LEFT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! incoming
cmfd % current(3,h,i,j,k) = t % results(1,score_index) % sum
end if
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /) )
matching_bins(i_filter_surf) = RIGHT
matching_bins(i_filter_surf) = OUT_RIGHT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! outgoing
* t % stride) + 1
cmfd % current(3,h,i,j,k) = t % results(1,score_index) % sum
matching_bins(i_filter_surf) = IN_RIGHT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(4,h,i,j,k) = t % results(1,score_index) % sum
! Back surface
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /))
matching_bins(i_filter_surf) = BACK
matching_bins(i_filter_surf) = OUT_BACK
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! outgoing
* t % stride) + 1
cmfd % current(5,h,i,j,k) = t % results(1,score_index) % sum
if (j > 1) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j-1, k /))
matching_bins(i_filter_surf) = FRONT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! incoming
cmfd % current(6,h,i,j,k) = t % results(1,score_index) % sum
end if
matching_bins(i_filter_surf) = IN_BACK
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(6,h,i,j,k) = t % results(1,score_index) % sum
! Front surface
if (j < ny) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j+1, k /))
matching_bins(i_filter_surf) = BACK
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! incoming
cmfd % current(7,h,i,j,k) = t % results(1,score_index) % sum
end if
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /))
matching_bins(i_filter_surf) = FRONT
matching_bins(i_filter_surf) = OUT_FRONT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! outgoing
* t % stride) + 1
cmfd % current(7,h,i,j,k) = t % results(1,score_index) % sum
matching_bins(i_filter_surf) = IN_FRONT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(8,h,i,j,k) = t % results(1,score_index) % sum
! Bottom surface
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /))
matching_bins(i_filter_surf) = BOTTOM
matching_bins(i_filter_surf) = OUT_BOTTOM
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! outgoing
* t % stride) + 1
cmfd % current(9,h,i,j,k) = t % results(1,score_index) % sum
if (k > 1) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k-1 /))
matching_bins(i_filter_surf) = TOP
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! incoming
cmfd % current(10,h,i,j,k) = t % results(1,score_index) % sum
end if
matching_bins(i_filter_surf) = IN_BOTTOM
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(10,h,i,j,k) = t % results(1,score_index) % sum
! Top surface
if (k < nz) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k+1 /))
matching_bins(i_filter_surf) = BOTTOM
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! incoming
cmfd % current(11,h,i,j,k) = t % results(1,score_index) % sum
end if
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /))
matching_bins(i_filter_surf) = TOP
matching_bins(i_filter_surf) = OUT_TOP
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! outgoing
* t % stride) + 1
cmfd % current(11,h,i,j,k) = t % results(1,score_index) % sum
matching_bins(i_filter_surf) = IN_TOP
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
cmfd % current(12,h,i,j,k) = t % results(1,score_index) % sum
end if TALLY
@ -478,7 +447,6 @@ contains
! Get leakage
leakage = ZERO
LEAK: do l = 1, 3
leakage = leakage + ((cmfd % current(4*l,g,i,j,k) - &
cmfd % current(4*l-1,g,i,j,k))) - &
((cmfd % current(4*l-2,g,i,j,k) - &

View file

@ -126,7 +126,7 @@ contains
if (.not. allocated(this % hxyz)) allocate(this % hxyz(3,nx,ny,nz))
! Allocate surface currents
if (.not. allocated(this % current)) allocate(this % current(12,ng,nx,ny,nz))
if (.not. allocated(this % current)) allocate(this % current(6,ng,nx,ny,nz))
! Allocate source distributions
if (.not. allocated(this % cmfd_src)) allocate(this % cmfd_src(ng,nx,ny,nz))

View file

@ -534,9 +534,12 @@ contains
filt % n_bins = 2 * m % n_dimension
allocate(filt % surfaces(2 * m % n_dimension))
if (m % n_dimension == 2) then
filt % surfaces = (/ LEFT, RIGHT, BACK, FRONT /)
filt % surfaces = (/ OUT_LEFT, OUT_RIGHT, OUT_BACK, OUT_FRONT, &
IN_LEFT, IN_RIGHT, IN_BACK, IN_FRONT /)
elseif (m % n_dimension == 3) then
filt % surfaces = (/ LEFT, RIGHT, BACK, FRONT, BOTTOM, TOP /)
filt % surfaces = (/ OUT_LEFT, OUT_RIGHT, OUT_BACK, OUT_FRONT, &
OUT_BOTTOM, OUT_TOP, IN_LEFT, IN_RIGHT, IN_BACK, IN_FRONT, &
IN_BOTTOM, IN_TOP /)
end if
end select
t % find_filter(FILTER_SURFACE) = n_filters

View file

@ -359,12 +359,18 @@ module constants
! Tally surface current directions
integer, parameter :: &
LEFT = 1, & ! x min
RIGHT = 2, & ! x max
BACK = 3, & ! y min
FRONT = 4, & ! y max
BOTTOM = 5, & ! z min
TOP = 6 ! z max
OUT_LEFT = 1, & ! x min
OUT_RIGHT = 2, & ! x max
OUT_BACK = 3, & ! y min
OUT_FRONT = 4, & ! y max
OUT_BOTTOM = 5, & ! z min
OUT_TOP = 6, & ! z max
IN_LEFT = 7, & ! x min
IN_RIGHT = 8, & ! x max
IN_BACK = 9, & ! y min
IN_FRONT = 10, & ! y max
IN_BOTTOM = 11, & ! z min
IN_TOP = 12 ! z max
! Tally trigger types and threshold
integer, parameter :: &

View file

@ -3698,12 +3698,15 @@ contains
allocate(SurfaceFilter :: t % filters(n_filters) % obj)
select type (filt => t % filters(size(t % filters)) % obj)
type is (SurfaceFilter)
filt % n_bins = 2 * m % n_dimension
allocate(filt % surfaces(2 * m % n_dimension))
filt % n_bins = 4 * m % n_dimension
allocate(filt % surfaces(4 * m % n_dimension))
if (m % n_dimension == 2) then
filt % surfaces = (/ LEFT, RIGHT, BACK, FRONT /)
filt % surfaces = (/ OUT_LEFT, OUT_RIGHT, OUT_BACK, OUT_FRONT, &
IN_LEFT, IN_RIGHT, IN_BACK, IN_FRONT /)
elseif (m % n_dimension == 3) then
filt % surfaces = (/ LEFT, RIGHT, BACK, FRONT, BOTTOM, TOP /)
filt % surfaces = (/ OUT_LEFT, OUT_RIGHT, OUT_BACK, OUT_FRONT, &
OUT_BOTTOM, OUT_TOP, IN_LEFT, IN_RIGHT, IN_BACK, &
IN_FRONT, IN_BOTTOM, IN_TOP /)
end if
end select
t % find_filter(FILTER_SURFACE) = size(t % filters)

View file

@ -1034,43 +1034,85 @@ contains
matching_bins(i_filter_ein)))
end if
! Get the bin for this mesh cell
! Left Surface
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = LEFT
matching_bins(i_filter_surf) = OUT_LEFT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Net Current on Left", &
"Outgoing Current on Left", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = IN_LEFT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Incoming Current on Left", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
! Right Surface
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = RIGHT
matching_bins(i_filter_surf) = OUT_RIGHT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Net Current on Right", &
"Outgoing Current on Right", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = IN_RIGHT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Incoming Current on Right", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
! Back Surface
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = BACK
matching_bins(i_filter_surf) = OUT_BACK
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Net Current on Back", &
"Outgoing Current on Back", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = IN_BACK
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Incoming Current on Back", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
! Front Surface
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = FRONT
matching_bins(i_filter_surf) = OUT_FRONT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Net Current on Front", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = IN_FRONT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
@ -1081,26 +1123,45 @@ contains
! Bottom Surface
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = BOTTOM
matching_bins(i_filter_surf) = OUT_BOTTOM
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Net Current on Bottom", &
"Outgoing Current on Bottom", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = IN_BOTTOM
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Incoming Current on Bottom", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
! Top Surface
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = TOP
matching_bins(i_filter_surf) = OUT_TOP
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Net Current on Top", &
"Outgoing Current on Top", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = IN_TOP
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Incoming Current on Top", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
end do
end do
end do
end do

View file

@ -2410,7 +2410,7 @@ contains
! OUT_TOP
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = TOP
matching_bins(i_filter_surf) = OUT_TOP
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2426,7 +2426,7 @@ contains
ijk0(2) <= m % dimension(2) .and. &
ijk0(3) < m % dimension(3)) then
ijk0(3) = ijk0(3) + 1
matching_bins(i_filter_surf) = BOTTOM
matching_bins(i_filter_surf) = IN_BOTTOM
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2443,14 +2443,14 @@ contains
! OUT_BOTTOM
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = BOTTOM
matching_bins(i_filter_surf) = OUT_BOTTOM
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
end if
! IN_TOP
@ -2459,14 +2459,14 @@ contains
ijk0(2) <= m % dimension(2) .and. &
ijk0(3) <= m % dimension(3) + 1) then
ijk0(3) = ijk0(3) - 1
matching_bins(i_filter_surf) = TOP
matching_bins(i_filter_surf) = IN_TOP
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
ijk0(3) = ijk0(3) + 1
end if
end do
@ -2480,7 +2480,7 @@ contains
! OUT_FRONT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = FRONT
matching_bins(i_filter_surf) = OUT_FRONT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2496,7 +2496,7 @@ contains
ijk0(2) < m % dimension(2) .and. &
ijk0(3) <= m % dimension(3)) then
ijk0(2) = ijk0(2) + 1
matching_bins(i_filter_surf) = BACK
matching_bins(i_filter_surf) = IN_BACK
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2513,14 +2513,14 @@ contains
! OUT_BACK
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = BACK
matching_bins(i_filter_surf) = OUT_BACK
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
end if
! IN_FRONT
@ -2529,14 +2529,14 @@ contains
ijk0(2) <= m % dimension(2) + 1 .and. &
ijk0(3) <= m % dimension(3)) then
ijk0(2) = ijk0(2) - 1
matching_bins(i_filter_surf) = FRONT
matching_bins(i_filter_surf) = IN_FRONT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
ijk0(2) = ijk0(2) + 1
end if
end do
@ -2550,7 +2550,7 @@ contains
! OUT_RIGHT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = RIGHT
matching_bins(i_filter_surf) = OUT_RIGHT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2566,7 +2566,7 @@ contains
ijk0(2) <= m % dimension(2) .and. &
ijk0(3) <= m % dimension(3)) then
ijk0(1) = ijk0(1) + 1
matching_bins(i_filter_surf) = LEFT
matching_bins(i_filter_surf) = IN_LEFT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2583,14 +2583,14 @@ contains
! OUT_LEFT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = LEFT
matching_bins(i_filter_surf) = OUT_LEFT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
end if
! IN_RIGHT
@ -2599,14 +2599,14 @@ contains
ijk0(2) <= m % dimension(2) .and. &
ijk0(3) <= m % dimension(3)) then
ijk0(1) = ijk0(1) - 1
matching_bins(i_filter_surf) = RIGHT
matching_bins(i_filter_surf) = IN_RIGHT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
ijk0(1) = ijk0(1) + 1
end if
end do
@ -2655,7 +2655,7 @@ contains
! OUT_RIGHT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = RIGHT
matching_bins(i_filter_surf) = OUT_RIGHT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2671,7 +2671,7 @@ contains
ijk0(2) <= m % dimension(2) .and. &
ijk0(3) <= m % dimension(3)) then
ijk0(1) = ijk0(1) + 1
matching_bins(i_filter_surf) = LEFT
matching_bins(i_filter_surf) = IN_LEFT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2688,14 +2688,14 @@ contains
! OUT_LEFT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = LEFT
matching_bins(i_filter_surf) = OUT_LEFT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
end if
! IN_RIGHT
@ -2704,14 +2704,14 @@ contains
ijk0(2) <= m % dimension(2) .and. &
ijk0(3) <= m % dimension(3)) then
ijk0(1) = ijk0(1) - 1
matching_bins(i_filter_surf) = RIGHT
matching_bins(i_filter_surf) = IN_RIGHT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
ijk0(1) = ijk0(1) + 1
end if
@ -2723,7 +2723,7 @@ contains
! OUT_FRONT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = FRONT
matching_bins(i_filter_surf) = OUT_FRONT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2739,7 +2739,7 @@ contains
ijk0(2) < m % dimension(2) .and. &
ijk0(3) <= m % dimension(3)) then
ijk0(2) = ijk0(2) + 1
matching_bins(i_filter_surf) = BACK
matching_bins(i_filter_surf) = IN_BACK
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2756,14 +2756,14 @@ contains
! OUT_BACK
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = BACK
matching_bins(i_filter_surf) = OUT_BACK
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
end if
! IN_FRONT
@ -2772,14 +2772,14 @@ contains
ijk0(2) <= m % dimension(2) + 1 .and. &
ijk0(3) <= m % dimension(3)) then
ijk0(2) = ijk0(2) - 1
matching_bins(i_filter_surf) = FRONT
matching_bins(i_filter_surf) = IN_FRONT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
ijk0(2) = ijk0(2) + 1
end if
@ -2791,7 +2791,7 @@ contains
! OUT_TOP
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = TOP
matching_bins(i_filter_surf) = OUT_TOP
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2807,7 +2807,7 @@ contains
ijk0(2) <= m % dimension(2) .and. &
ijk0(3) < m % dimension(3)) then
ijk0(3) = ijk0(3) + 1
matching_bins(i_filter_surf) = BOTTOM
matching_bins(i_filter_surf) = IN_BOTTOM
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2824,14 +2824,14 @@ contains
! OUT_BOTTOM
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = BOTTOM
matching_bins(i_filter_surf) = OUT_BOTTOM
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
end if
! IN_TOP
@ -2840,14 +2840,14 @@ contains
ijk0(2) <= m % dimension(2) .and. &
ijk0(3) <= m % dimension(3) + 1) then
ijk0(3) = ijk0(3) - 1
matching_bins(i_filter_surf) = TOP
matching_bins(i_filter_surf) = IN_TOP
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value - p % wgt
t % results(1, filter_index) % value + p % wgt
ijk0(3) = ijk0(3) + 1
end if

View file

@ -332,7 +332,7 @@ contains
mesh_indices_to_bin(m, (/ i, j, k /))
! Left Surface
matching_bins(i_filter_surf) = LEFT
matching_bins(i_filter_surf) = OUT_LEFT
filter_index = &
sum((matching_bins(1:size(t % filters)) - 1) * t % stride) + 1
call get_trigger_uncertainty(std_dev, rel_err, 1, filter_index, t)
@ -345,7 +345,7 @@ contains
trigger % variance = std_dev**2
! Right Surface
matching_bins(i_filter_surf) = RIGHT
matching_bins(i_filter_surf) = OUT_RIGHT
filter_index = &
sum((matching_bins(1:size(t % filters)) - 1) * t % stride) + 1
call get_trigger_uncertainty(std_dev, rel_err, 1, filter_index, t)
@ -358,7 +358,7 @@ contains
trigger % variance = trigger % std_dev**2
! Back Surface
matching_bins(i_filter_surf) = BACK
matching_bins(i_filter_surf) = OUT_BACK
filter_index = &
sum((matching_bins(1:size(t % filters)) - 1) * t % stride) + 1
call get_trigger_uncertainty(std_dev, rel_err, 1, filter_index, t)
@ -371,7 +371,7 @@ contains
trigger % variance = trigger % std_dev**2
! Front Surface
matching_bins(i_filter_surf) = FRONT
matching_bins(i_filter_surf) = OUT_FRONT
filter_index = &
sum((matching_bins(1:size(t % filters)) - 1) * t % stride) + 1
call get_trigger_uncertainty(std_dev, rel_err, 1, filter_index, t)
@ -384,7 +384,7 @@ contains
trigger % variance = trigger % std_dev**2
! Bottom Surface
matching_bins(i_filter_surf) = BOTTOM
matching_bins(i_filter_surf) = OUT_BOTTOM
filter_index = &
sum((matching_bins(1:size(t % filters)) - 1) * t % stride) + 1
call get_trigger_uncertainty(std_dev, rel_err, 1, filter_index, t)
@ -397,7 +397,7 @@ contains
trigger % variance = trigger % std_dev**2
! Top Surface
matching_bins(i_filter_surf) = TOP
matching_bins(i_filter_surf) = OUT_TOP
filter_index = &
sum((matching_bins(1:size(t % filters)) - 1) * t % stride) + 1
call get_trigger_uncertainty(std_dev, rel_err, 1, filter_index, t)