crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-28 22:26:08 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Move cross_lattice to C++

Browse source
This commit is contained in:
Sterling Harper 2018-08-20 22:45:28 -04:00
parent 918735732a
commit c98142dd23
4 changed files with 89 additions and 78 deletions
Download patch file Download diff file Expand all files Collapse all files

82
src/geometry.F90
View file

@ -41,11 +41,18 @@ module geometry
function find_cell_c(p, n_search_cells, search_cells) &
bind(C, name="find_cell") result(found)
import Particle, C_INT, C_BOOL
type(Particle), intent(in) :: p
type(Particle), intent(inout) :: p
integer(C_INT), intent(in), value :: n_search_cells
integer(C_INT), intent(in), optional :: search_cells(n_search_cells)
logical(C_BOOL) :: found
end function find_cell_c
subroutine cross_lattice(p, lattice_translation) &
bind(C, name="cross_lattice")
import Particle, C_INT
type(Particle), intent(inout) :: p
integer(C_INT), intent(in) :: lattice_translation(3)
end subroutine cross_lattice
end interface
contains
@ -77,79 +84,6 @@ contains
end subroutine find_cell
!===============================================================================
! CROSS_LATTICE moves a particle into a new lattice element
!===============================================================================
subroutine cross_lattice(p, lattice_translation)
type(Particle), intent(inout) :: p
integer, intent(in) :: lattice_translation(3)
integer :: j
integer :: i_xyz(3) ! indices in lattice
logical :: found ! particle found in cell?
class(Lattice), pointer :: lat
j = p % n_coord
lat => lattices(p % coord(j) % lattice) % obj
if (verbosity >= 10 .or. trace) then
call write_message(" Crossing lattice " // trim(to_str(lat % id())) &
&// ". Current position (" // trim(to_str(p % coord(j) % lattice_x)) &
&// "," // trim(to_str(p % coord(j) % lattice_y)) // "," &
&// trim(to_str(p % coord(j) % lattice_z)) // ")")
end if
! Set the lattice indices.
p % coord(j) % lattice_x = p % coord(j) % lattice_x + lattice_translation(1)
p % coord(j) % lattice_y = p % coord(j) % lattice_y + lattice_translation(2)
p % coord(j) % lattice_z = p % coord(j) % lattice_z + lattice_translation(3)
i_xyz(1) = p % coord(j) % lattice_x
i_xyz(2) = p % coord(j) % lattice_y
i_xyz(3) = p % coord(j) % lattice_z
! Set the new coordinate position.
p % coord(j) % xyz = lat % get_local_xyz(p % coord(j - 1) % xyz, i_xyz)
OUTSIDE_LAT: if (.not. lat % are_valid_indices(i_xyz)) then
! The particle is outside the lattice. Search for it from base coord
p % n_coord = 1
call find_cell(p, found)
if (.not. found) then
if (p % alive) then ! Particle may have been killed in find_cell
call particle_mark_as_lost(p, "Could not locate particle " &
// trim(to_str(p % id)) // " after crossing a lattice boundary.")
return
end if
end if
else OUTSIDE_LAT
! Find cell in next lattice element
p % coord(j) % universe = &
lat % get([i_xyz(1), i_xyz(2), i_xyz(3)]) + 1
call find_cell(p, found)
if (.not. found) then
! In some circumstances, a particle crossing the corner of a cell may
! not be able to be found in the next universe. In this scenario we cut
! off all lower-level coordinates and search from universe zero
! Remove lower coordinates
p % n_coord = 1
! Search for particle
call find_cell(p, found)
if (.not. found) then
call particle_mark_as_lost(p, "Could not locate particle " // &
trim(to_str(p % id)) // " after crossing a lattice boundary.")
return
end if
end if
end if OUTSIDE_LAT
end subroutine cross_lattice
!===============================================================================
! DISTANCE_TO_BOUNDARY calculates the distance to the nearest boundary for a
! particle 'p' traveling in a certain direction. For a cell in a subuniverse

64
src/geometry.cpp
View file

@ -1,5 +1,6 @@
#include "geometry.h"
#include <array>
#include <sstream>
#include "cell.h"
@ -230,7 +231,7 @@ find_cell(Particle* p, int n_search_cells, int* search_cells) {
<< lat.id << " but the lattice has no defined outer "
"universe.";
warning(err_msg);
found = false;
return false;
}
}
@ -243,4 +244,65 @@ find_cell(Particle* p, int n_search_cells, int* search_cells) {
return found;
}
//==============================================================================
extern "C" void
cross_lattice(Particle* p, int lattice_translation[3])
{
Lattice& lat {*lattices_c[p->coord[p->n_coord-1].lattice-1]};
if (openmc_verbosity >= 10 || openmc_trace) {
std::stringstream msg;
msg << " Crossing lattice " << lat.id << ". Current position ("
<< p->coord[p->n_coord-1].lattice_x << ","
<< p->coord[p->n_coord-1].lattice_y << ","
<< p->coord[p->n_coord-1].lattice_z << ")";
write_message(msg, 1);
}
// Set the lattice indices.
p->coord[p->n_coord-1].lattice_x += lattice_translation[0];
p->coord[p->n_coord-1].lattice_y += lattice_translation[1];
p->coord[p->n_coord-1].lattice_z += lattice_translation[2];
std::array<int, 3> i_xyz {p->coord[p->n_coord-1].lattice_x,
p->coord[p->n_coord-1].lattice_y,
p->coord[p->n_coord-1].lattice_z};
// Set the new coordinate position.
auto r = lat.get_local_position(p->coord[p->n_coord-2].xyz, i_xyz);
p->coord[p->n_coord-1].xyz[0] = r.x;
p->coord[p->n_coord-1].xyz[1] = r.y;
p->coord[p->n_coord-1].xyz[2] = r.z;
if (!lat.are_valid_indices(i_xyz)) {
// The particle is outside the lattice. Search for it from the base coords.
p->n_coord = 1;
bool found = find_cell(p, 0, nullptr);
if (!found && p->alive) {
std::stringstream err_msg;
err_msg << "Could not locate particle " << p->id
<< " after crossing a lattice boundary";
p->mark_as_lost(err_msg);
}
} else {
// Find cell in next lattice element.
p->coord[p->n_coord-1].universe = lat[i_xyz] + 1;
bool found = find_cell(p, 0, nullptr);
if (!found) {
// A particle crossing the corner of a lattice tile may not be found. In
// this case, search for it from the base coords.
p->n_coord = 1;
bool found = find_cell(p, 0, nullptr);
if (!found && p->alive) {
std::stringstream err_msg;
err_msg << "Could not locate particle " << p->id
<< " after crossing a lattice boundary";
p->mark_as_lost(err_msg);
}
}
}
}
} // namespace openmc

11
src/geometry.h
View file

@ -13,19 +13,26 @@ extern "C" int openmc_root_universe;
extern std::vector<int64_t> overlap_check_count;
//==============================================================================
//! Check for overlapping cells at the particle's position.
//! Check for overlapping cells at a particle's position.
//==============================================================================
extern "C" bool
check_cell_overlap(Particle* p);
//==============================================================================
//! Locate the particle in the geometry tree and set its geometry data fields.
//! Locate a particle in the geometry tree and set its geometry data fields.
//==============================================================================
extern "C" bool
find_cell(Particle* p, int n_search_cells, int* search_cells);
//==============================================================================
//! Move a particle into a new lattice tile.
//==============================================================================
extern "C" void
cross_lattice(Particle* p, int lattice_translation[3]);
} // namespace openmc
#endif // OPENMC_GEOMETRY_H

10
src/particle.h
View file

@ -4,8 +4,10 @@
//! \file particle.h
//! \brief Particle type
#include <cstdint>
#include <array>
#include <cstdint>
#include <sstream>
#include <string>
#include "openmc.h"
@ -152,6 +154,12 @@ extern "C" {
//! \param message A warning message to display
void mark_as_lost(const char* message);
void mark_as_lost(const std::string& message)
{mark_as_lost(message.c_str());}
void mark_as_lost(const std::stringstream& message)
{mark_as_lost(message.str());}
//! create a particle restart HDF5 file
void write_restart();
};