mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-24 03:55:38 -04:00
Move FILL_LATTICE chunk of find_cell to C++
This commit is contained in:
Sterling Harper
parent
cddf39b201
commit
918735732a
5 changed files with 73 additions and 85 deletions
|
|
@ -65,78 +65,15 @@ contains
|
|||
!===============================================================================
|
||||
|
||||
recursive subroutine find_cell(p, found, search_cells)
|
||||
|
||||
type(Particle), intent(inout) :: p
|
||||
logical, intent(inout) :: found
|
||||
integer, optional :: search_cells(:)
|
||||
integer :: j ! coordinate level index
|
||||
integer :: i_xyz(3) ! indices in lattice
|
||||
integer :: i_cell ! index in cells array
|
||||
|
||||
if (present(search_cells)) then
|
||||
found = find_cell_c(p, size(search_cells), search_cells-1)
|
||||
else
|
||||
found = find_cell_c(p, 0)
|
||||
end if
|
||||
j = p % n_coord
|
||||
i_cell = p % coord(j) % cell
|
||||
|
||||
if (found) then
|
||||
associate(c => cells(i_cell))
|
||||
CELL_TYPE: if (c % type() == FILL_UNIVERSE) then
|
||||
j = j + 1
|
||||
p % n_coord = j
|
||||
call find_cell(p, found)
|
||||
j = p % n_coord
|
||||
|
||||
elseif (c % type() == FILL_LATTICE) then CELL_TYPE
|
||||
! ======================================================================
|
||||
! CELL CONTAINS LATTICE, RECURSIVELY FIND CELL
|
||||
|
||||
associate (lat => lattices(c % fill() + 1) % obj)
|
||||
! Determine lattice indices
|
||||
i_xyz = lat % get_indices(p % coord(j) % xyz + TINY_BIT * p % coord(j) % uvw)
|
||||
|
||||
! Store lower level coordinates
|
||||
p % coord(j + 1) % xyz = lat % get_local_xyz(p % coord(j) % xyz, i_xyz)
|
||||
p % coord(j + 1) % uvw = p % coord(j) % uvw
|
||||
|
||||
! set particle lattice indices
|
||||
p % coord(j + 1) % lattice = c % fill() + 1
|
||||
p % coord(j + 1) % lattice_x = i_xyz(1)
|
||||
p % coord(j + 1) % lattice_y = i_xyz(2)
|
||||
p % coord(j + 1) % lattice_z = i_xyz(3)
|
||||
|
||||
! Set the next lowest coordinate level.
|
||||
if (lat % are_valid_indices(i_xyz)) then
|
||||
! Particle is inside the lattice.
|
||||
p % coord(j + 1) % universe = &
|
||||
lat % get([i_xyz(1), i_xyz(2), i_xyz(3)]) + 1
|
||||
|
||||
else
|
||||
! Particle is outside the lattice.
|
||||
if (lat % outer() == NO_OUTER_UNIVERSE) then
|
||||
call warning("Particle " // trim(to_str(p %id)) &
|
||||
// " is outside lattice " // trim(to_str(lat % id())) &
|
||||
// " but the lattice has no defined outer universe.")
|
||||
found = .false.
|
||||
return
|
||||
else
|
||||
p % coord(j + 1) % universe = lat % outer() + 1
|
||||
end if
|
||||
end if
|
||||
end associate
|
||||
|
||||
! Move particle to next level and search for the lower cells.
|
||||
j = j + 1
|
||||
p % n_coord = j
|
||||
|
||||
call find_cell(p, found)
|
||||
j = p % n_coord
|
||||
|
||||
end if CELL_TYPE
|
||||
end associate
|
||||
end if
|
||||
|
||||
end subroutine find_cell
|
||||
|
||||
|
|
|
|||
|
|
@ -8,9 +8,6 @@
|
|||
#include "lattice.h"
|
||||
#include "settings.h"
|
||||
|
||||
//TODO: remove this include
|
||||
#include <iostream>
|
||||
|
||||
|
||||
namespace openmc {
|
||||
|
||||
|
|
@ -153,8 +150,7 @@ find_cell(Particle* p, int n_search_cells, int* search_cells) {
|
|||
//========================================================================
|
||||
//! Found a lower universe, update this coord level then search the next.
|
||||
|
||||
// Add another coordinate level.
|
||||
//++p->n_coord;
|
||||
// Set the lower coordinate level universe.
|
||||
p->coord[p->n_coord].universe = c.fill + 1;
|
||||
|
||||
// Set the position and direction.
|
||||
|
|
@ -190,6 +186,57 @@ find_cell(Particle* p, int n_search_cells, int* search_cells) {
|
|||
+ w*c.rotation[11];
|
||||
p->coord[p->n_coord].rotated = true;
|
||||
}
|
||||
|
||||
// Update the coordinate level and recurse.
|
||||
++p->n_coord;
|
||||
find_cell(p, 0, nullptr);
|
||||
|
||||
} else if (c.type == FILL_LATTICE) {
|
||||
//========================================================================
|
||||
//! Found a lower lattice, update this coord level then search the next.
|
||||
|
||||
Lattice& lat {*lattices_c[c.fill]};
|
||||
|
||||
// Determine lattice indices.
|
||||
Position r {p->coord[p->n_coord-1].xyz};
|
||||
Direction u {p->coord[p->n_coord-1].uvw};
|
||||
r += TINY_BIT * u;
|
||||
auto i_xyz = lat.get_indices(r);
|
||||
|
||||
// Store lower level coordinates.
|
||||
r = lat.get_local_position(p->coord[p->n_coord-1].xyz, i_xyz);
|
||||
p->coord[p->n_coord].xyz[0] = r.x;
|
||||
p->coord[p->n_coord].xyz[1] = r.y;
|
||||
p->coord[p->n_coord].xyz[2] = r.z;
|
||||
p->coord[p->n_coord].uvw[0] = u.x;
|
||||
p->coord[p->n_coord].uvw[1] = u.y;
|
||||
p->coord[p->n_coord].uvw[2] = u.z;
|
||||
|
||||
// Set lattice indices.
|
||||
p->coord[p->n_coord].lattice = c.fill + 1;
|
||||
p->coord[p->n_coord].lattice_x = i_xyz[0];
|
||||
p->coord[p->n_coord].lattice_y = i_xyz[1];
|
||||
p->coord[p->n_coord].lattice_z = i_xyz[2];
|
||||
|
||||
// Set the lower coordinate level universe.
|
||||
if (lat.are_valid_indices(i_xyz)) {
|
||||
p->coord[p->n_coord].universe = lat[i_xyz] + 1;
|
||||
} else {
|
||||
if (lat.outer != NO_OUTER_UNIVERSE) {
|
||||
p->coord[p->n_coord].universe = lat.outer + 1;
|
||||
} else {
|
||||
std::stringstream err_msg;
|
||||
err_msg << "Particle " << p->id << " is outside lattice "
|
||||
<< lat.id << " but the lattice has no defined outer "
|
||||
"universe.";
|
||||
warning(err_msg);
|
||||
found = false;
|
||||
}
|
||||
}
|
||||
|
||||
// Update the coordinate level and recurse.
|
||||
++p->n_coord;
|
||||
find_cell(p, 0, nullptr);
|
||||
}
|
||||
}
|
||||
|
||||
|
|
|
|||
|
|
@ -191,13 +191,6 @@ module geometry_header
|
|||
integer(C_INT32_T) :: offset
|
||||
end function lattice_offset_c
|
||||
|
||||
function lattice_outer_c(lat_ptr) bind(C, name='lattice_outer') &
|
||||
result(outer)
|
||||
import C_PTR, C_INT32_T
|
||||
type(C_PTR), intent(in), value :: lat_ptr
|
||||
integer(C_INT32_T) :: outer
|
||||
end function lattice_outer_c
|
||||
|
||||
function lattice_universe_c(lat_ptr, i_xyz) &
|
||||
bind(C, name='lattice_universe') result(univ)
|
||||
import C_PTR, C_INT32_t, C_INT
|
||||
|
|
@ -238,7 +231,6 @@ module geometry_header
|
|||
procedure :: get_indices => lattice_get_indices
|
||||
procedure :: get_local_xyz => lattice_get_local_xyz
|
||||
procedure :: offset => lattice_offset
|
||||
procedure :: outer => lattice_outer
|
||||
procedure :: to_hdf5 => lattice_to_hdf5
|
||||
end type Lattice
|
||||
|
||||
|
|
@ -275,7 +267,6 @@ module geometry_header
|
|||
! Boolean expression of half-spaces
|
||||
|
||||
! Rotation matrix and translation vector
|
||||
real(8), allocatable :: translation(:)
|
||||
real(8), allocatable :: rotation(:)
|
||||
real(8), allocatable :: rotation_matrix(:,:)
|
||||
|
||||
|
|
@ -370,12 +361,6 @@ contains
|
|||
offset = lattice_offset_c(this % ptr, map, i_xyz)
|
||||
end function lattice_offset
|
||||
|
||||
function lattice_outer(this) result(outer)
|
||||
class(Lattice), intent(in) :: this
|
||||
integer(C_INT32_T) :: outer
|
||||
outer = lattice_outer_c(this % ptr)
|
||||
end function lattice_outer
|
||||
|
||||
subroutine lattice_to_hdf5(this, group)
|
||||
class(Lattice), intent(in) :: this
|
||||
integer(HID_T), intent(in) :: group
|
||||
|
|
|
|||
|
|
@ -923,8 +923,6 @@ extern "C" {
|
|||
int32_t lattice_offset(Lattice *lat, int map, const int i_xyz[3])
|
||||
{return lat->offset(map, i_xyz);}
|
||||
|
||||
int32_t lattice_outer(Lattice *lat) {return lat->outer;}
|
||||
|
||||
void lattice_to_hdf5(Lattice *lat, hid_t group) {lat->to_hdf5(group);}
|
||||
|
||||
int32_t lattice_universe(Lattice *lat, const int i_xyz[3])
|
||||
|
|
|
|||
|
|
@ -54,6 +54,13 @@ public:
|
|||
|
||||
virtual int32_t& operator[](const int i_xyz[3]) = 0;
|
||||
|
||||
int32_t&
|
||||
operator[](std::array<int, 3> i_xyz)
|
||||
{
|
||||
int i_xyz_[3] {i_xyz[0], i_xyz[1], i_xyz[2]};
|
||||
return operator[](i_xyz_);
|
||||
}
|
||||
|
||||
virtual LatticeIter begin();
|
||||
LatticeIter end();
|
||||
|
||||
|
|
@ -76,6 +83,13 @@ public:
|
|||
//! otherwise.
|
||||
virtual bool are_valid_indices(const int i_xyz[3]) const = 0;
|
||||
|
||||
bool
|
||||
are_valid_indices(std::array<int, 3> i_xyz) const
|
||||
{
|
||||
int i_xyz_[3] {i_xyz[0], i_xyz[1], i_xyz[2]};
|
||||
return are_valid_indices(i_xyz_);
|
||||
}
|
||||
|
||||
//! \brief Find the next lattice surface crossing
|
||||
//! \param r A 3D Cartesian coordinate.
|
||||
//! \param u A 3D Cartesian direction.
|
||||
|
|
@ -98,6 +112,13 @@ public:
|
|||
virtual Position
|
||||
get_local_position(Position r, const int i_xyz[3]) const = 0;
|
||||
|
||||
Position
|
||||
get_local_position(Position r, std::array<int, 3> i_xyz) const
|
||||
{
|
||||
int i_xyz_[3] {i_xyz[0], i_xyz[1], i_xyz[2]};
|
||||
return get_local_position(r, i_xyz_);
|
||||
}
|
||||
|
||||
//! \brief Check flattened lattice index.
|
||||
//! \param indx The index for a lattice tile.
|
||||
//! \return true if the given index fit within the lattice bounds. False
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue