crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-27 13:45:36 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

changed partial current tallies to net current tallies

Browse source
This commit is contained in:
Sam Shaner 2016-08-23 09:41:53 -04:00
parent 83869d2ab6
commit bc387186ba
10 changed files with 356 additions and 118 deletions
Download patch file Download diff file Expand all files Collapse all files

14
examples/xml/pincell/geometry.xml
View file

@ -8,16 +8,16 @@
-->
<!-- Surfaces for fuel, gap, cladding. Dimensions from Figure 2 in BEAVRS -->
<surface id="1" type="z-cylinder" coeffs="0. 0. 0.39218" /> <!-- Fuel OR -->
<surface id="2" type="z-cylinder" coeffs="0. 0. 0.40005" /> <!-- Clad IR -->
<surface id="3" type="z-cylinder" coeffs="0. 0. 0.45720" /> <!-- Clad OR -->
<surface id="1" type="z-cylinder" coeffs="0. 0. 3.39218" /> <!-- Fuel OR -->
<surface id="2" type="z-cylinder" coeffs="0. 0. 3.40005" /> <!-- Clad IR -->
<surface id="3" type="z-cylinder" coeffs="0. 0. 3.45720" /> <!-- Clad OR -->
<!-- Reflective surfaces on outside of pin-cell. The lattice pitch is 1.25984
cm (taken from Table 2 in BEAVRS). -->
<surface id="4" type="x-plane" coeffs="-0.62992" boundary="reflective" />
<surface id="5" type="x-plane" coeffs=" 0.62992" boundary="reflective" />
<surface id="6" type="y-plane" coeffs="-0.62992" boundary="reflective" />
<surface id="7" type="y-plane" coeffs=" 0.62992" boundary="reflective" />
<surface id="4" type="x-plane" coeffs="-8.62992" boundary="reflective" />
<surface id="5" type="x-plane" coeffs=" 8.62992" boundary="vacuum" />
<surface id="6" type="y-plane" coeffs="-8.62992" boundary="reflective" />
<surface id="7" type="y-plane" coeffs=" 8.62992" boundary="reflective" />
<cell id="1" material="1" region=" -1" /> <!-- UO2 Fuel -->
<cell id="2" material="2" region="1 -2" /> <!-- Helium gap -->

12
examples/xml/pincell/settings.xml
View file

@ -3,7 +3,7 @@
<!-- Define how many particles to run and for how many batches -->
<eigenvalue>
<batches>100</batches>
<batches>20</batches>
<inactive>10</inactive>
<particles>1000</particles>
</eigenvalue>
@ -14,8 +14,8 @@
<source>
<space type="box">
<parameters>
-0.62992 -0.62992 -1.
0.62992 0.62992 1.
-8.62992 -8.62992 -1.
8.62992 8.62992 1.
</parameters>
</space>
</source>
@ -24,9 +24,9 @@
bounds for a mesh over which the Shannon entropy should be
calculated. The extent in the z direction is made arbitrarily large. -->
<entropy>
<lower_left>-0.39218 -0.39218 -1.e50</lower_left>
<upper_right>0.39218 0.39218 1.e50</upper_right>
<lower_left>-3.39218 -3.39218 -1.e50</lower_left>
<upper_right>3.39218 3.39218 1.e50</upper_right>
<dimension>10 10 1</dimension>
</entropy>
</settings>
</settings>

11
examples/xml/pincell/tallies.xml
View file

@ -2,9 +2,9 @@
<tallies>
<mesh id="1" type="regular">
<dimension>100 100 1</dimension>
<lower_left>-0.62992 -0.62992 -1.e50</lower_left>
<upper_right>0.62992 0.62992 1.e50</upper_right>
<dimension>2 2 1</dimension>
<lower_left>-8.62992 -8.62992 -1.e50</lower_left>
<upper_right>8.62992 8.62992 1.e50</upper_right>
</mesh>
<tally id="1">
@ -13,4 +13,9 @@
<scores>flux fission nu-fission</scores>
</tally>
<tally id="2">
<filter type="mesh" bins="1" />
<scores>current</scores>
</tally>
</tallies>

28
src/cmfd_data.F90
View file

@ -51,8 +51,8 @@ contains
subroutine compute_xs()
use constants, only: FILTER_MESH, FILTER_ENERGYIN, FILTER_ENERGYOUT, &
FILTER_SURFACE, OUT_LEFT, OUT_RIGHT, OUT_BACK, &
OUT_FRONT, OUT_BOTTOM, OUT_TOP, CMFD_NOACCEL, &
FILTER_SURFACE, LEFT, RIGHT, BACK, &
FRONT, BOTTOM, TOP, CMFD_NOACCEL, &
ZERO, ONE, TINY_BIT
use error, only: fatal_error
use global, only: cmfd, n_cmfd_tallies, cmfd_tallies, meshes,&
@ -237,7 +237,7 @@ contains
! Left surface
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /))
matching_bins(i_filter_surf) = OUT_LEFT
matching_bins(i_filter_surf) = LEFT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t%stride) + 1 ! outgoing
cmfd % current(1,h,i,j,k) = t % results(1,score_index) % sum
@ -245,7 +245,7 @@ contains
if (i > 1) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i-1, j, k /))
matching_bins(i_filter_surf) = OUT_RIGHT
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
@ -255,7 +255,7 @@ contains
if (i < nx) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i+1, j, k /) )
matching_bins(i_filter_surf) = OUT_LEFT
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
@ -263,7 +263,7 @@ contains
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /) )
matching_bins(i_filter_surf) = OUT_RIGHT
matching_bins(i_filter_surf) = RIGHT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! outgoing
cmfd % current(4,h,i,j,k) = t % results(1,score_index) % sum
@ -271,7 +271,7 @@ contains
! Back surface
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /))
matching_bins(i_filter_surf) = OUT_BACK
matching_bins(i_filter_surf) = BACK
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! outgoing
cmfd % current(5,h,i,j,k) = t % results(1,score_index) % sum
@ -279,7 +279,7 @@ contains
if (j > 1) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j-1, k /))
matching_bins(i_filter_surf) = OUT_FRONT
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
@ -289,7 +289,7 @@ contains
if (j < ny) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j+1, k /))
matching_bins(i_filter_surf) = OUT_BACK
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
@ -297,7 +297,7 @@ contains
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /))
matching_bins(i_filter_surf) = OUT_FRONT
matching_bins(i_filter_surf) = FRONT
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! outgoing
cmfd % current(8,h,i,j,k) = t % results(1,score_index) % sum
@ -305,7 +305,7 @@ contains
! Bottom surface
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /))
matching_bins(i_filter_surf) = OUT_BOTTOM
matching_bins(i_filter_surf) = BOTTOM
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! outgoing
cmfd % current(9,h,i,j,k) = t % results(1,score_index) % sum
@ -313,7 +313,7 @@ contains
if (k > 1) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k-1 /))
matching_bins(i_filter_surf) = OUT_TOP
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
@ -323,7 +323,7 @@ contains
if (k < nz) then
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k+1 /))
matching_bins(i_filter_surf) = OUT_BOTTOM
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
@ -331,7 +331,7 @@ contains
matching_bins(i_filter_mesh) = mesh_indices_to_bin(m, &
(/ i, j, k /))
matching_bins(i_filter_surf) = OUT_TOP
matching_bins(i_filter_surf) = TOP
score_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1 ! outgoing
cmfd % current(12,h,i,j,k) = t % results(1,score_index) % sum

5
src/cmfd_input.F90
View file

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

12
src/constants.F90
View file

@ -359,12 +359,12 @@ module constants
! Tally surface current directions
integer, parameter :: &
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
LEFT = 1, & ! x min
RIGHT = 2, & ! x max
BACK = 3, & ! y min
FRONT = 4, & ! y max
BOTTOM = 5, & ! z min
TOP = 6 ! z max
! Tally trigger types and threshold
integer, parameter :: &

5
src/input_xml.F90
View file

@ -3701,10 +3701,9 @@ contains
filt % n_bins = 2 * m % n_dimension
allocate(filt % surfaces(2 * m % n_dimension))
if (m % n_dimension == 2) then
filt % surfaces = (/ OUT_LEFT, OUT_RIGHT, OUT_BACK, OUT_FRONT /)
filt % surfaces = (/ LEFT, RIGHT, BACK, FRONT /)
elseif (m % n_dimension == 3) then
filt % surfaces = (/ OUT_LEFT, OUT_RIGHT, OUT_BACK, OUT_FRONT,&
OUT_BOTTOM, OUT_TOP /)
filt % surfaces = (/ LEFT, RIGHT, BACK, FRONT, BOTTOM, TOP /)
end if
end select
t % find_filter(FILTER_SURFACE) = size(t % filters)

24
src/output.F90
View file

@ -1037,66 +1037,66 @@ contains
! Left Surface
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = OUT_LEFT
matching_bins(i_filter_surf) = LEFT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current to Left", &
"Net 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) = OUT_RIGHT
matching_bins(i_filter_surf) = RIGHT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current to Right", &
"Net 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) = OUT_BACK
matching_bins(i_filter_surf) = BACK
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current to Back", &
"Net 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) = OUT_FRONT
matching_bins(i_filter_surf) = FRONT
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current to Front", &
"Net Current on Front", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
! Bottom Surface
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, (/ i, j, k /))
matching_bins(i_filter_surf) = OUT_BOTTOM
matching_bins(i_filter_surf) = BOTTOM
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current to Bottom", &
"Net 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) = OUT_TOP
matching_bins(i_filter_surf) = TOP
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
write(UNIT=unit_tally, FMT='(5X,A,T35,A,"+/- ",A)') &
"Outgoing Current to Top", &
"Net Current on Top", &
to_str(t % results(1,filter_index) % sum), &
trim(to_str(t % results(1,filter_index) % sum_sq))
end do

351
src/tally.F90
View file

@ -2367,7 +2367,7 @@ contains
call get_mesh_indices(m, xyz0, ijk0(:m % n_dimension), start_in_mesh)
call get_mesh_indices(m, xyz1, ijk1(:m % n_dimension), end_in_mesh)
! Check to if start or end is in mesh -- if not, check if track still
! Check to see if start or end is in mesh -- if not, check if track still
! intersects with mesh
if ((.not. start_in_mesh) .and. (.not. end_in_mesh)) then
if (m % n_dimension == 2) then
@ -2407,8 +2407,10 @@ contains
if (uvw(3) > 0) then
do j = ijk0(3), ijk1(3) - 1
ijk0(3) = j
! OUT_TOP
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_TOP
matching_bins(i_filter_surf) = TOP
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2417,12 +2419,14 @@ contains
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
end if
end do
else
do j = ijk0(3), ijk1(3) + 1, -1
ijk0(3) = j
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_BOTTOM
! IN_BOTTOM
if (ijk0(1) >= 1 .and. ijk0(2) >= 1 .and. ijk0(3) >= 0 .and. &
ijk0(1) <= m % dimension(1) .and. &
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_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2430,6 +2434,40 @@ contains
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
ijk0(3) = ijk0(3) - 1
end if
end do
else
do j = ijk0(3), ijk1(3) + 1, -1
ijk0(3) = j
! OUT_BOTTOM
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = 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
end if
! IN_TOP
if (ijk0(1) >= 1 .and. ijk0(2) >= 1 .and. ijk0(3) > 1 .and. &
ijk0(1) <= m % dimension(1) .and. &
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_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
ijk0(3) = ijk0(3) + 1
end if
end do
end if
@ -2439,8 +2477,10 @@ contains
if (uvw(2) > 0) then
do j = ijk0(2), ijk1(2) - 1
ijk0(2) = j
! OUT_FRONT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_FRONT
matching_bins(i_filter_surf) = FRONT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2449,12 +2489,14 @@ contains
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
end if
end do
else
do j = ijk0(2), ijk1(2) + 1, -1
ijk0(2) = j
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_BACK
! IN_BACK
if (ijk0(1) >= 1 .and. ijk0(2) >= 0 .and. ijk0(3) >= 1 .and. &
ijk0(1) <= m % dimension(1) .and. &
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_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2462,6 +2504,40 @@ contains
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
ijk0(2) = ijk0(2) - 1
end if
end do
else
do j = ijk0(2), ijk1(2) + 1, -1
ijk0(2) = j
! OUT_BACK
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = 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
end if
! IN_FRONT
if (ijk0(1) >= 1 .and. ijk0(2) > 1 .and. ijk0(3) >= 1 .and. &
ijk0(1) <= m % dimension(1) .and. &
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_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
ijk0(2) = ijk0(2) + 1
end if
end do
end if
@ -2471,8 +2547,10 @@ contains
if (uvw(1) > 0) then
do j = ijk0(1), ijk1(1) - 1
ijk0(1) = j
! OUT_RIGHT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_RIGHT
matching_bins(i_filter_surf) = RIGHT
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2481,12 +2559,14 @@ contains
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
end if
end do
else
do j = ijk0(1), ijk1(1) + 1, -1
ijk0(1) = j
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_LEFT
! IN_LEFT
if (ijk0(1) >= 0 .and. ijk0(2) >= 1 .and. ijk0(3) >= 1 .and. &
ijk0(1) < m % dimension(1) .and. &
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_mesh) = &
mesh_indices_to_bin(m, ijk0)
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
@ -2494,6 +2574,40 @@ contains
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
ijk0(1) = ijk0(1) - 1
end if
end do
else
do j = ijk0(1), ijk1(1) + 1, -1
ijk0(1) = j
! OUT_LEFT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = 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
end if
! IN_RIGHT
if (ijk0(1) > 1 .and. ijk0(2) >= 1 .and. ijk0(3) >= 1 .and. &
ijk0(1) <= m % dimension(1) + 1 .and. &
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_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
ijk0(1) = ijk0(1) + 1
end if
end do
end if
@ -2533,94 +2647,215 @@ contains
distance = minval(d)
! Now use the minimum distance and diretion of the particle to
! Now use the minimum distance and direction of the particle to
! determine which surface was crossed
if (distance == d(1)) then
if (uvw(1) > 0) then
! Crossing into right mesh cell -- this is treated as outgoing
! current from (i,j,k)
! OUT_RIGHT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_RIGHT
matching_bins(i_filter_surf) = 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
end if
! IN_LEFT
if (ijk0(1) >= 0 .and. ijk0(2) >= 1 .and. ijk0(3) >= 1 .and. &
ijk0(1) < m % dimension(1) .and. &
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_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
ijk0(1) = ijk0(1) - 1
end if
ijk0(1) = ijk0(1) + 1
xyz_cross(1) = xyz_cross(1) + m % width(1)
else
! Crossing into left mesh cell -- this is treated as outgoing
! current in (i,j,k)
! OUT_LEFT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_LEFT
matching_bins(i_filter_surf) = 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
end if
! IN_RIGHT
if (ijk0(1) > 1 .and. ijk0(2) >= 1 .and. ijk0(3) >= 1 .and. &
ijk0(1) <= m % dimension(1) + 1 .and. &
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_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
ijk0(1) = ijk0(1) + 1
end if
ijk0(1) = ijk0(1) - 1
xyz_cross(1) = xyz_cross(1) - m % width(1)
end if
elseif (distance == d(2)) then
if (uvw(2) > 0) then
! Crossing into front mesh cell -- this is treated as outgoing
! current in (i,j,k)
! OUT_FRONT
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_FRONT
matching_bins(i_filter_surf) = 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
end if
! IN_BACK
if (ijk0(1) >= 1 .and. ijk0(2) >= 0 .and. ijk0(3) >= 1 .and. &
ijk0(1) <= m % dimension(1) .and. &
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_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
ijk0(2) = ijk0(2) - 1
end if
ijk0(2) = ijk0(2) + 1
xyz_cross(2) = xyz_cross(2) + m % width(2)
else
! Crossing into back mesh cell -- this is treated as outgoing
! current in (i,j,k)
! OUT_BACK
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_BACK
matching_bins(i_filter_surf) = 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
end if
! IN_FRONT
if (ijk0(1) >= 1 .and. ijk0(2) > 1 .and. ijk0(3) >= 1 .and. &
ijk0(1) <= m % dimension(1) .and. &
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_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
ijk0(2) = ijk0(2) + 1
end if
ijk0(2) = ijk0(2) - 1
xyz_cross(2) = xyz_cross(2) - m % width(2)
end if
else if (distance == d(3)) then
if (uvw(3) > 0) then
! Crossing into top mesh cell -- this is treated as outgoing
! current in (i,j,k)
! OUT_TOP
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_TOP
matching_bins(i_filter_surf) = 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
end if
! IN_BOTTOM
if (ijk0(1) >= 1 .and. ijk0(2) >= 1 .and. ijk0(3) >= 0 .and. &
ijk0(1) <= m % dimension(1) .and. &
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_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
ijk0(3) = ijk0(3) - 1
end if
ijk0(3) = ijk0(3) + 1
xyz_cross(3) = xyz_cross(3) + m % width(3)
else
! Crossing into bottom mesh cell -- this is treated as outgoing
! current in (i,j,k)
! OUT_BOTTOM
if (all(ijk0 >= 1) .and. all(ijk0 <= m % dimension)) then
matching_bins(i_filter_surf) = OUT_BOTTOM
matching_bins(i_filter_surf) = 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
end if
! IN_TOP
if (ijk0(1) >= 1 .and. ijk0(2) >= 1 .and. ijk0(3) > 1 .and. &
ijk0(1) <= m % dimension(1) .and. &
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_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
ijk0(3) = ijk0(3) + 1
end if
ijk0(3) = ijk0(3) - 1
xyz_cross(3) = xyz_cross(3) - m % width(3)
end if
end if
! Determine scoring index
if (matching_bins(i_filter_surf) > 0) then
filter_index = sum((matching_bins(1:size(t % filters)) - 1) &
* t % stride) + 1
! Check for errors
if (filter_index <= 0 .or. filter_index > &
t % total_filter_bins) then
call fatal_error("Score index outside range.")
end if
! Add to surface current tally
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
end if
! Calculate new coordinates
xyz0 = xyz0 + distance * uvw
end do

12
src/trigger.F90
View file

@ -332,7 +332,7 @@ contains
mesh_indices_to_bin(m, (/ i, j, k /))
! Left Surface
matching_bins(i_filter_surf) = OUT_LEFT
matching_bins(i_filter_surf) = 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) = OUT_RIGHT
matching_bins(i_filter_surf) = 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) = OUT_BACK
matching_bins(i_filter_surf) = 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) = OUT_FRONT
matching_bins(i_filter_surf) = 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) = OUT_BOTTOM
matching_bins(i_filter_surf) = 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) = OUT_TOP
matching_bins(i_filter_surf) = 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)