Merge pull request #342 from paulromano/atomic

Use OpenMP atomic updates instead of critical sections
This commit is contained in:
Jon Walsh 2014-11-26 08:26:37 -05:00
commit 0f1270e005
6 changed files with 57 additions and 76 deletions

View file

@ -11,7 +11,7 @@ module geometry
use tally, only: score_surface_current
implicit none
contains
!===============================================================================
@ -265,13 +265,13 @@ contains
lattice_edge = .true.
end if
end if
if (lattice_edge) then
! In this case the neutron is leaving the lattice, so we move it
! out, remove all lower coordinate levels and then search from
! universe 0.
p % coord => p % coord0
call deallocate_coord(p % coord % next)
@ -288,7 +288,7 @@ contains
p % last_material = p % material
p % material = c % material
! We'll still make a new coordinate for the particle, as
! We'll still make a new coordinate for the particle, as
! distance_to_boundary will still need to track through lattice
! widths even though there's nothing in them but this material
@ -407,10 +407,9 @@ contains
! Score to global leakage tally
if (tallies_on) then
!$omp critical
!$omp atomic
global_tallies(LEAKAGE) % value = &
global_tallies(LEAKAGE) % value + p % wgt
!$omp end critical
end if
! Display message
@ -645,7 +644,7 @@ contains
return
end if
end if
end subroutine cross_surface
!===============================================================================
@ -907,7 +906,7 @@ contains
if (quad < ZERO) then
! no intersection with cylinder
d = INFINITY
d = INFINITY
elseif (on_surface) then
! particle is on the cylinder, thus one distance is
@ -956,7 +955,7 @@ contains
if (quad < ZERO) then
! no intersection with cylinder
d = INFINITY
d = INFINITY
elseif (on_surface) then
! particle is on the cylinder, thus one distance is
@ -1005,7 +1004,7 @@ contains
if (quad < ZERO) then
! no intersection with cylinder
d = INFINITY
d = INFINITY
elseif (on_surface) then
! particle is on the cylinder, thus one distance is
@ -1052,7 +1051,7 @@ contains
if (quad < ZERO) then
! no intersection with sphere
d = INFINITY
d = INFINITY
elseif (on_surface) then
! particle is on the sphere, thus one distance is
@ -1099,7 +1098,7 @@ contains
if (quad < ZERO) then
! no intersection with cone
d = INFINITY
d = INFINITY
elseif (on_surface) then
! particle is on the cone, thus one distance is positive/negative
@ -1118,7 +1117,7 @@ contains
d = (-k - quad)/a
b = (-k + quad)/a
! determine the smallest positive solution
! determine the smallest positive solution
if (d < ZERO) then
if (b > ZERO) then
d = b
@ -1148,7 +1147,7 @@ contains
if (quad < ZERO) then
! no intersection with cone
d = INFINITY
d = INFINITY
elseif (on_surface) then
! particle is on the cone, thus one distance is positive/negative
@ -1167,7 +1166,7 @@ contains
d = (-k - quad)/a
b = (-k + quad)/a
! determine the smallest positive solution
! determine the smallest positive solution
if (d < ZERO) then
if (b > ZERO) then
d = b
@ -1197,7 +1196,7 @@ contains
if (quad < ZERO) then
! no intersection with cone
d = INFINITY
d = INFINITY
elseif (on_surface) then
! particle is on the cone, thus one distance is positive/negative
@ -1216,7 +1215,7 @@ contains
d = (-k - quad)/a
b = (-k + quad)/a
! determine the smallest positive solution
! determine the smallest positive solution
if (d < ZERO) then
if (b > ZERO) then
d = b
@ -1274,7 +1273,7 @@ contains
! logic here checks whether the relative difference is within floating
! point precision.
if (d < dist) then
if (d < dist) then
if (abs(d - dist)/dist >= FP_REL_PRECISION) then
dist = d
if (u > 0) then
@ -1566,9 +1565,8 @@ contains
! Increment number of lost particles
p % alive = .false.
!$omp critical
!$omp atomic
n_lost_particles = n_lost_particles + 1
!$omp end critical
! Abort the simulation if the maximum number of lost particles has been
! reached

View file

@ -62,7 +62,7 @@ module global
! Cross section arrays
type(Nuclide), allocatable, target :: nuclides(:) ! Nuclide cross-sections
type(SAlphaBeta), allocatable, target :: sab_tables(:) ! S(a,b) tables
type(XsListing), allocatable, target :: xs_listings(:) ! cross_sections.xml listings
type(XsListing), allocatable, target :: xs_listings(:) ! cross_sections.xml listings
! Cross section caches
type(NuclideMicroXS), allocatable :: micro_xs(:) ! Cache for each nuclide
@ -119,7 +119,7 @@ module global
! 2) track-length estimate of k-eff
! 3) leakage fraction
type(TallyResult), target :: global_tallies(N_GLOBAL_TALLIES)
type(TallyResult), allocatable, target :: global_tallies(:)
! Tally map structure
type(TallyMap), allocatable :: tally_maps(:)
@ -300,7 +300,7 @@ module global
logical :: write_initial_source = .false.
! ============================================================================
! CMFD VARIABLES
! CMFD VARIABLES
! Main object
type(cmfd_type) :: cmfd
@ -310,11 +310,11 @@ module global
! CMFD communicator
integer :: cmfd_comm
! Timing objects
type(Timer) :: time_cmfd ! timer for whole cmfd calculation
type(Timer) :: time_cmfdbuild ! timer for matrix build
type(Timer) :: time_cmfdsolve ! timer for solver
type(Timer) :: time_cmfdsolve ! timer for solver
! Flag for active core map
logical :: cmfd_coremap = .false.
@ -390,7 +390,7 @@ module global
! RESONANCE SCATTERING VARIABLES
logical :: treat_res_scat = .false. ! is resonance scattering treated?
integer :: n_res_scatterers_total = 0 ! total number of resonant scatterers
integer :: n_res_scatterers_total = 0 ! total number of resonant scatterers
type(Nuclide0K), allocatable, target :: nuclides_0K(:) ! 0K nuclides info
!$omp threadprivate(micro_xs, material_xs, fission_bank, n_bank, &
@ -399,14 +399,14 @@ module global
contains
!===============================================================================
! FREE_MEMORY deallocates and clears all global allocatable arrays in the
! FREE_MEMORY deallocates and clears all global allocatable arrays in the
! program
!===============================================================================
subroutine free_memory()
integer :: i ! Loop Index
! Deallocate cells, surfaces, materials
if (allocated(cells)) deallocate(cells)
if (allocated(universes)) deallocate(universes)
@ -449,6 +449,7 @@ contains
if (allocated(entropy_p)) deallocate(entropy_p)
! Deallocate tally-related arrays
if (allocated(global_tallies)) deallocate(global_tallies)
if (allocated(meshes)) deallocate(meshes)
if (allocated(tallies)) then
! First call the clear routines
@ -488,7 +489,7 @@ contains
! Deallocate track_identifiers
if (allocated(track_identifiers)) deallocate(track_identifiers)
! Deallocate dictionaries
call cell_dict % clear()
call universe_dict % clear()
@ -525,7 +526,7 @@ contains
if (allocated(ufs_mesh % width)) deallocate(ufs_mesh % width)
deallocate(ufs_mesh)
end if
end subroutine free_memory
end module global

View file

@ -256,22 +256,19 @@ contains
p % last_wgt = p % wgt
! Score implicit absorption estimate of keff
!$omp critical
!$omp atomic
global_tallies(K_ABSORPTION) % value = &
global_tallies(K_ABSORPTION) % value + p % absorb_wgt * &
micro_xs(i_nuclide) % nu_fission / micro_xs(i_nuclide) % absorption
!$omp end critical
else
! See if disappearance reaction happens
if (micro_xs(i_nuclide) % absorption > &
prn() * micro_xs(i_nuclide) % total) then
! Score absorption estimate of keff
!$omp critical
!$omp atomic
global_tallies(K_ABSORPTION) % value = &
global_tallies(K_ABSORPTION) % value + p % wgt * &
micro_xs(i_nuclide) % nu_fission / micro_xs(i_nuclide) % absorption
!$omp end critical
p % alive = .false.
p % event = EVENT_ABSORB
@ -797,7 +794,7 @@ contains
sampling_scheme = 'cxs'
end if
! otherwise, use free gas model
! otherwise, use free gas model
else
if (E >= FREE_GAS_THRESHOLD * kT .and. awr > ONE) then
v_target = ZERO
@ -859,7 +856,7 @@ contains
m = (nuc % elastic_0K(i_E_up + 1) - xs_up) &
& / (nuc % energy_0K(i_E_up + 1) - nuc % energy_0K(i_E_up))
xs_up = xs_up + m * (E_up - nuc % energy_0K(i_E_up))
! get max 0K xs value over range of practical relative energies
xs_max = max(xs_low, &
& maxval(nuc % elastic_0K(i_E_low + 1 : i_E_up - 1)), xs_up)
@ -972,7 +969,7 @@ contains
case default
call fatal_error("Not a recognized resonance scattering treatment!")
end select
end subroutine sample_target_velocity
!===============================================================================

View file

@ -281,10 +281,9 @@ contains
! get the score and tally it
score = last_wgt * calc_pn(n, mu)
!$omp critical
!$omp atomic
t % results(score_index, filter_index) % value = &
t % results(score_index, filter_index) % value + score
!$omp end critical
end do
j = j + t % moment_order(j)
cycle SCORE_LOOP
@ -347,10 +346,9 @@ contains
! get the score and tally it
score = wgt * calc_pn(n, mu)
!$omp critical
!$omp atomic
t % results(score_index, filter_index) % value = &
t % results(score_index, filter_index) % value + score
!$omp end critical
end do
j = j + t % moment_order(j)
cycle SCORE_LOOP
@ -542,10 +540,9 @@ contains
end select
! Add score to tally
!$omp critical
!$omp atomic
t % results(score_index, filter_index) % value = &
t % results(score_index, filter_index) % value + score
!$omp end critical
end do SCORE_LOOP
@ -617,10 +614,9 @@ contains
i_filter = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1
! Add score to tally
!$omp critical
!$omp atomic
t % results(i_score, i_filter) % value = &
t % results(i_score, i_filter) % value + score
!$omp end critical
end do
! reset outgoing energy bin and score index
@ -1015,10 +1011,9 @@ contains
end if
! Add score to tally
!$omp critical
!$omp atomic
t % results(score_index, filter_index) % value = &
t % results(score_index, filter_index) % value + score
!$omp end critical
end do SCORE_LOOP
@ -1214,10 +1209,9 @@ contains
end select
! Add score to tally
!$omp critical
!$omp atomic
t % results(score_index, filter_index) % value = &
t % results(score_index, filter_index) % value + score
!$omp end critical
end do SCORE_LOOP
@ -1366,10 +1360,9 @@ contains
end select
! Add score to tally
!$omp critical
!$omp atomic
t % results(score_index, filter_index) % value = &
t % results(score_index, filter_index) % value + score
!$omp end critical
end do MATERIAL_SCORE_LOOP
@ -1787,10 +1780,9 @@ contains
end if
! Add score to tally
!$omp critical
!$omp atomic
t % results(score_index, filter_index) % value = &
t % results(score_index, filter_index) % value + score
!$omp end critical
end do SCORE_LOOP
@ -2021,10 +2013,9 @@ contains
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0 + 1, .true.)
filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1
!$omp critical
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
!$omp end critical
end if
end do
else
@ -2035,10 +2026,9 @@ contains
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0 + 1, .true.)
filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1
!$omp critical
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
!$omp end critical
end if
end do
end if
@ -2053,10 +2043,9 @@ contains
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0 + 1, .true.)
filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1
!$omp critical
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
!$omp end critical
end if
end do
else
@ -2067,10 +2056,9 @@ contains
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0 + 1, .true.)
filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1
!$omp critical
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
!$omp end critical
end if
end do
end if
@ -2085,10 +2073,9 @@ contains
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0 + 1, .true.)
filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1
!$omp critical
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
!$omp end critical
end if
end do
else
@ -2099,10 +2086,9 @@ contains
matching_bins(i_filter_mesh) = &
mesh_indices_to_bin(m, ijk0 + 1, .true.)
filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1
!$omp critical
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
!$omp end critical
end if
end do
end if
@ -2224,10 +2210,9 @@ contains
end if
! Add to surface current tally
!$omp critical
!$omp atomic
t % results(1, filter_index) % value = &
t % results(1, filter_index) % value + p % wgt
!$omp end critical
end if
! Calculate new coordinates

View file

@ -19,6 +19,9 @@ contains
subroutine configure_tallies()
! Allocate global tallies
allocate(global_tallies(N_GLOBAL_TALLIES))
call setup_tally_arrays()
call setup_tally_maps()

View file

@ -62,9 +62,8 @@ contains
n_event = 0
! Add paricle's starting weight to count for normalizing tallies later
!$omp critical
!$omp atomic
total_weight = total_weight + p % wgt
!$omp end critical
! Force calculation of cross-sections by setting last energy to zero
micro_xs % last_E = ZERO
@ -112,11 +111,10 @@ contains
call score_tracklength_tally(p, distance)
! Score track-length estimate of k-eff
!$omp critical
!$omp atomic
global_tallies(K_TRACKLENGTH) % value = &
global_tallies(K_TRACKLENGTH) % value + p % wgt * distance * &
material_xs % nu_fission
!$omp end critical
if (d_collision > d_boundary) then
! ====================================================================
@ -140,11 +138,10 @@ contains
! PARTICLE HAS COLLISION
! Score collision estimate of keff
!$omp critical
!$omp atomic
global_tallies(K_COLLISION) % value = &
global_tallies(K_COLLISION) % value + p % wgt * &
material_xs % nu_fission / material_xs % total
!$omp end critical
! score surface current tallies -- this has to be done before the collision
! since the direction of the particle will change and we need to use the