From 93862cea249e4441beb714e719120da3bd2b7d0a Mon Sep 17 00:00:00 2001 From: Sam Shaner Date: Tue, 23 Aug 2016 11:33:23 -0400 Subject: [PATCH] changed current tallying so in/out partial currents are tallies separately --- src/cmfd_data.F90 | 122 ++++++++++++++++---------------------------- src/cmfd_header.F90 | 2 +- src/cmfd_input.F90 | 7 ++- src/constants.F90 | 18 ++++--- src/input_xml.F90 | 11 ++-- src/output.F90 | 85 +++++++++++++++++++++++++----- src/tally.F90 | 72 +++++++++++++------------- src/trigger.F90 | 12 ++--- 8 files changed, 185 insertions(+), 144 deletions(-) diff --git a/src/cmfd_data.F90 b/src/cmfd_data.F90 index 0a1490bb9..befdc9c46 100644 --- a/src/cmfd_data.F90 +++ b/src/cmfd_data.F90 @@ -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) - & diff --git a/src/cmfd_header.F90 b/src/cmfd_header.F90 index e743cc928..7d253c308 100644 --- a/src/cmfd_header.F90 +++ b/src/cmfd_header.F90 @@ -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)) diff --git a/src/cmfd_input.F90 b/src/cmfd_input.F90 index 77e8f4348..50f078a18 100644 --- a/src/cmfd_input.F90 +++ b/src/cmfd_input.F90 @@ -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 diff --git a/src/constants.F90 b/src/constants.F90 index ceb68b9e2..481ca8e03 100644 --- a/src/constants.F90 +++ b/src/constants.F90 @@ -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 :: & diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 8de198386..bcdf310e8 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -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) diff --git a/src/output.F90 b/src/output.F90 index 7aacd6713..eae654787 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -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 diff --git a/src/tally.F90 b/src/tally.F90 index 6bb36d9ce..7f743d134 100644 --- a/src/tally.F90 +++ b/src/tally.F90 @@ -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 diff --git a/src/trigger.F90 b/src/trigger.F90 index 6f886cd8e..4af03a2de 100644 --- a/src/trigger.F90 +++ b/src/trigger.F90 @@ -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)