mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-27 05:35:49 -04:00
somehow changes to lattice methods were lost
This commit is contained in:
Gavin Ridley
parent
41ef034a0b
commit
f0e1110e01
5 changed files with 90 additions and 118 deletions
|
|
@ -1,7 +1,6 @@
|
|||
#ifndef OPENMC_LATTICE_H
|
||||
#define OPENMC_LATTICE_H
|
||||
|
||||
#include <array>
|
||||
#include <cstdint>
|
||||
#include <string>
|
||||
#include <unordered_map>
|
||||
|
|
@ -9,6 +8,7 @@
|
|||
#include "hdf5.h"
|
||||
#include "pugixml.hpp"
|
||||
|
||||
#include "openmc/array.h"
|
||||
#include "openmc/constants.h"
|
||||
#include "openmc/memory.h"
|
||||
#include "openmc/position.h"
|
||||
|
|
@ -26,7 +26,6 @@ enum class LatticeType {
|
|||
rect, hex
|
||||
};
|
||||
|
||||
|
||||
//==============================================================================
|
||||
// Global variables
|
||||
//==============================================================================
|
||||
|
|
@ -60,7 +59,7 @@ public:
|
|||
|
||||
virtual ~Lattice() {}
|
||||
|
||||
virtual int32_t& operator[](std::array<int, 3> i_xyz) = 0;
|
||||
virtual int32_t const& operator[](array<int, 3> const& i_xyz) = 0;
|
||||
|
||||
virtual LatticeIter begin();
|
||||
LatticeIter end();
|
||||
|
|
@ -83,14 +82,7 @@ public:
|
|||
//! \param i_xyz[3] The indices for a lattice tile.
|
||||
//! \return true if the given indices fit within the lattice bounds. False
|
||||
//! 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_);
|
||||
}
|
||||
virtual bool are_valid_indices(array<int, 3> const& i_xyz) const = 0;
|
||||
|
||||
//! \brief Find the next lattice surface crossing
|
||||
//! \param r A 3D Cartesian coordinate.
|
||||
|
|
@ -103,15 +95,17 @@ public:
|
|||
|
||||
//! \brief Find the lattice tile indices for a given point.
|
||||
//! \param r A 3D Cartesian coordinate.
|
||||
//! \return An array containing the indices of a lattice tile.
|
||||
virtual std::array<int, 3> get_indices(Position r, Direction u) const = 0;
|
||||
//! \param u Direction of a particle
|
||||
//! \param result resulting indices to save to
|
||||
virtual void get_indices(
|
||||
Position r, Direction u, array<int, 3>& result) const = 0;
|
||||
|
||||
//! \brief Get coordinates local to a lattice tile.
|
||||
//! \param r A 3D Cartesian coordinate.
|
||||
//! \param i_xyz The indices for a lattice tile.
|
||||
//! \return Local 3D Cartesian coordinates.
|
||||
virtual Position
|
||||
get_local_position(Position r, const std::array<int, 3> i_xyz) const = 0;
|
||||
virtual Position get_local_position(
|
||||
Position r, const array<int, 3>& i_xyz) const = 0;
|
||||
|
||||
//! \brief Check flattened lattice index.
|
||||
//! \param indx The index for a lattice tile.
|
||||
|
|
@ -125,7 +119,7 @@ public:
|
|||
//! \param i_xyz[3] The indices for a lattice tile.
|
||||
//! \return Distribcell offset i.e. the largest instance number for the target
|
||||
//! cell found in the geometry tree under this lattice tile.
|
||||
virtual int32_t& offset(int map, const int i_xyz[3]) = 0;
|
||||
virtual int32_t& offset(int map, array<int, 3> const& i_xyz) = 0;
|
||||
|
||||
//! \brief Get the distribcell offset for a lattice tile.
|
||||
//! \param The map index for the target cell.
|
||||
|
|
@ -211,19 +205,18 @@ class RectLattice : public Lattice
|
|||
public:
|
||||
explicit RectLattice(pugi::xml_node lat_node);
|
||||
|
||||
int32_t& operator[](std::array<int, 3> i_xyz);
|
||||
int32_t const& operator[](array<int, 3> const& i_xyz);
|
||||
|
||||
bool are_valid_indices(const int i_xyz[3]) const;
|
||||
bool are_valid_indices(array<int, 3> const& i_xyz) const;
|
||||
|
||||
std::pair<double, array<int, 3>> distance(
|
||||
Position r, Direction u, const array<int, 3>& i_xyz) const;
|
||||
|
||||
std::array<int, 3> get_indices(Position r, Direction u) const;
|
||||
void get_indices(Position r, Direction u, array<int, 3>& result) const;
|
||||
|
||||
Position
|
||||
get_local_position(Position r, const std::array<int, 3> i_xyz) const;
|
||||
Position get_local_position(Position r, const array<int, 3>& i_xyz) const;
|
||||
|
||||
int32_t& offset(int map, const int i_xyz[3]);
|
||||
int32_t& offset(int map, array<int, 3> const& i_xyz);
|
||||
|
||||
int32_t offset(int map, int indx) const;
|
||||
|
||||
|
|
@ -235,11 +228,6 @@ private:
|
|||
array<int, 3> n_cells_; //!< Number of cells along each axis
|
||||
Position lower_left_; //!< Global lower-left corner of the lattice
|
||||
Position pitch_; //!< Lattice tile width along each axis
|
||||
|
||||
// Convenience aliases
|
||||
int &nx {n_cells_[0]};
|
||||
int &ny {n_cells_[1]};
|
||||
int &nz {n_cells_[2]};
|
||||
};
|
||||
|
||||
//==============================================================================
|
||||
|
|
@ -249,7 +237,7 @@ class HexLattice : public Lattice
|
|||
public:
|
||||
explicit HexLattice(pugi::xml_node lat_node);
|
||||
|
||||
int32_t& operator[](std::array<int, 3> i_xyz);
|
||||
int32_t const& operator[](array<int, 3> const& i_xyz);
|
||||
|
||||
LatticeIter begin();
|
||||
|
||||
|
|
@ -260,14 +248,13 @@ public:
|
|||
std::pair<double, array<int, 3>> distance(
|
||||
Position r, Direction u, const array<int, 3>& i_xyz) const;
|
||||
|
||||
std::array<int, 3> get_indices(Position r, Direction u) const;
|
||||
void get_indices(Position r, Direction u, array<int, 3>& result) const;
|
||||
|
||||
Position
|
||||
get_local_position(Position r, const std::array<int, 3> i_xyz) const;
|
||||
Position get_local_position(Position r, const array<int, 3>& i_xyz) const;
|
||||
|
||||
bool is_valid_index(int indx) const;
|
||||
|
||||
int32_t& offset(int map, const int i_xyz[3]);
|
||||
int32_t& offset(int map, array<int, 3> const& i_xyz);
|
||||
|
||||
int32_t offset(int map, int indx) const;
|
||||
|
||||
|
|
|
|||
|
|
@ -159,9 +159,8 @@ find_cell_inner(Particle& p, const NeighborList* neighbor_list)
|
|||
if (c_i.type_ == Fill::UNIVERSE) {
|
||||
offset += c_i.offset_[c.distribcell_index_];
|
||||
} else if (c_i.type_ == Fill::LATTICE) {
|
||||
auto& lat {*model::lattices[p.coord_[i + 1].lattice]};
|
||||
int i_xyz[3] {p.coord_[i + 1].lattice_x, p.coord_[i + 1].lattice_y,
|
||||
p.coord_[i + 1].lattice_z};
|
||||
auto& lat {*model::lattices[p.coord(i + 1).lattice]};
|
||||
const auto& i_xyz {p.coord(i + 1).lattice_i};
|
||||
if (lat.are_valid_indices(i_xyz)) {
|
||||
offset += lat.offset(c.distribcell_index_, i_xyz);
|
||||
}
|
||||
|
|
@ -226,16 +225,14 @@ find_cell_inner(Particle& p, const NeighborList* neighbor_list)
|
|||
}
|
||||
|
||||
// Determine lattice indices.
|
||||
auto i_xyz = lat.get_indices(coord.r, coord.u);
|
||||
auto& i_xyz {coord.lattice_i};
|
||||
lat.get_indices(coord.r, coord.u, i_xyz);
|
||||
|
||||
// Get local position in appropriate lattice cell
|
||||
coord.r = lat.get_local_position(coord.r, i_xyz);
|
||||
|
||||
// Set lattice indices.
|
||||
coord.lattice = c.fill_;
|
||||
coord.lattice_x = i_xyz[0];
|
||||
coord.lattice_y = i_xyz[1];
|
||||
coord.lattice_z = i_xyz[2];
|
||||
|
||||
// Set the lower coordinate level universe.
|
||||
if (lat.are_valid_indices(i_xyz)) {
|
||||
|
|
@ -318,10 +315,9 @@ cross_lattice(Particle& p, const BoundaryInfo& boundary)
|
|||
}
|
||||
|
||||
// Set the lattice indices.
|
||||
coord.lattice_x += boundary.lattice_translation[0];
|
||||
coord.lattice_y += boundary.lattice_translation[1];
|
||||
coord.lattice_z += boundary.lattice_translation[2];
|
||||
std::array<int, 3> i_xyz {coord.lattice_x, coord.lattice_y, coord.lattice_z};
|
||||
coord.lattice_i[0] += boundary.lattice_translation[0];
|
||||
coord.lattice_i[1] += boundary.lattice_translation[1];
|
||||
coord.lattice_i[2] += boundary.lattice_translation[2];
|
||||
|
||||
// Set the new coordinate position.
|
||||
const auto& upper_coord {p.coord(p.n_coord() - 2)};
|
||||
|
|
@ -331,9 +327,9 @@ cross_lattice(Particle& p, const BoundaryInfo& boundary)
|
|||
if (!cell->rotation_.empty()) {
|
||||
r = r.rotate(cell->rotation_);
|
||||
}
|
||||
p.r_local() = lat.get_local_position(r, i_xyz);
|
||||
p.r_local() = lat.get_local_position(r, coord.lattice_i);
|
||||
|
||||
if (!lat.are_valid_indices(i_xyz)) {
|
||||
if (!lat.are_valid_indices(coord.lattice_i)) {
|
||||
// The particle is outside the lattice. Search for it from the base coords.
|
||||
p.n_coord() = 1;
|
||||
bool found = exhaustive_find_cell(p);
|
||||
|
|
@ -345,7 +341,7 @@ cross_lattice(Particle& p, const BoundaryInfo& boundary)
|
|||
|
||||
} else {
|
||||
// Find cell in next lattice element.
|
||||
p.coord_[p.n_coord_-1].universe = lat[i_xyz];
|
||||
p.coord(p.n_coord() - 1).universe = lat[coord.lattice_i];
|
||||
bool found = exhaustive_find_cell(p);
|
||||
|
||||
if (!found) {
|
||||
|
|
@ -387,13 +383,12 @@ BoundaryInfo distance_to_boundary(Particle& p)
|
|||
// Find the distance to the next lattice tile crossing.
|
||||
if (coord.lattice != C_NONE) {
|
||||
auto& lat {*model::lattices[coord.lattice]};
|
||||
std::array<int, 3> i_xyz {coord.lattice_x, coord.lattice_y, coord.lattice_z};
|
||||
//TODO: refactor so both lattice use the same position argument (which
|
||||
//also means the lat.type attribute can be removed)
|
||||
std::pair<double, array<int, 3>> lattice_distance;
|
||||
switch (lat.type_) {
|
||||
case LatticeType::rect:
|
||||
lattice_distance = lat.distance(r, u, i_xyz);
|
||||
lattice_distance = lat.distance(r, u, coord.lattice_i);
|
||||
break;
|
||||
case LatticeType::hex:
|
||||
auto& cell_above {model::cells[p.coord(i - 1).cell]};
|
||||
|
|
@ -403,7 +398,7 @@ BoundaryInfo distance_to_boundary(Particle& p)
|
|||
r_hex = r_hex.rotate(cell_above->rotation_);
|
||||
}
|
||||
r_hex.z = coord.r.z;
|
||||
lattice_distance = lat.distance(r_hex, u, i_xyz);
|
||||
lattice_distance = lat.distance(r_hex, u, coord.lattice_i);
|
||||
break;
|
||||
}
|
||||
d_lat = lattice_distance.first;
|
||||
|
|
|
|||
118
src/lattice.cpp
118
src/lattice.cpp
|
|
@ -184,16 +184,19 @@ RectLattice::RectLattice(pugi::xml_node lat_node)
|
|||
if (univ_words.size() != n_cells_[0] * n_cells_[1] * n_cells_[2]) {
|
||||
fatal_error(fmt::format(
|
||||
"Expected {} universes for a rectangular lattice of size {}x{}x{} but {} "
|
||||
"were specified.", nx*ny*nz, nx, ny, nz, univ_words.size()));
|
||||
"were specified.",
|
||||
n_cells_[0] * n_cells_[1] * n_cells_[2], n_cells_[0], n_cells_[1],
|
||||
n_cells_[2], univ_words.size()));
|
||||
}
|
||||
|
||||
// Parse the universes.
|
||||
universes_.resize(nx*ny*nz, C_NONE);
|
||||
for (int iz = 0; iz < nz; iz++) {
|
||||
for (int iy = ny-1; iy > -1; iy--) {
|
||||
for (int ix = 0; ix < nx; ix++) {
|
||||
int indx1 = nx*ny*iz + nx*(ny-iy-1) + ix;
|
||||
int indx2 = nx*ny*iz + nx*iy + ix;
|
||||
universes_.resize(n_cells_[0] * n_cells_[1] * n_cells_[2], C_NONE);
|
||||
for (int iz = 0; iz < n_cells_[2]; iz++) {
|
||||
for (int iy = n_cells_[1] - 1; iy > -1; iy--) {
|
||||
for (int ix = 0; ix < n_cells_[0]; ix++) {
|
||||
int indx1 = n_cells_[0] * n_cells_[1] * iz +
|
||||
n_cells_[0] * (n_cells_[1] - iy - 1) + ix;
|
||||
int indx2 = n_cells_[0] * n_cells_[1] * iz + n_cells_[0] * iy + ix;
|
||||
universes_[indx1] = std::stoi(univ_words[indx2]);
|
||||
}
|
||||
}
|
||||
|
|
@ -202,17 +205,16 @@ RectLattice::RectLattice(pugi::xml_node lat_node)
|
|||
|
||||
//==============================================================================
|
||||
|
||||
int32_t&
|
||||
RectLattice::operator[](std::array<int, 3> i_xyz)
|
||||
int32_t const& RectLattice::operator[](array<int, 3> const& i_xyz)
|
||||
{
|
||||
int indx = nx*ny*i_xyz[2] + nx*i_xyz[1] + i_xyz[0];
|
||||
int indx =
|
||||
n_cells_[0] * n_cells_[1] * i_xyz[2] + n_cells_[0] * i_xyz[1] + i_xyz[0];
|
||||
return universes_[indx];
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
|
||||
bool
|
||||
RectLattice::are_valid_indices(const int i_xyz[3]) const
|
||||
bool RectLattice::are_valid_indices(array<int, 3> const& i_xyz) const
|
||||
{
|
||||
return ( (i_xyz[0] >= 0) && (i_xyz[0] < n_cells_[0])
|
||||
&& (i_xyz[1] >= 0) && (i_xyz[1] < n_cells_[1])
|
||||
|
|
@ -279,48 +281,44 @@ std::pair<double, array<int, 3>> RectLattice::distance(
|
|||
|
||||
//==============================================================================
|
||||
|
||||
std::array<int, 3>
|
||||
RectLattice::get_indices(Position r, Direction u) const
|
||||
void RectLattice::get_indices(
|
||||
Position r, Direction u, array<int, 3>& result) const
|
||||
{
|
||||
// Determine x index, accounting for coincidence
|
||||
double ix_ {(r.x - lower_left_.x) / pitch_.x};
|
||||
long ix_close {std::lround(ix_)};
|
||||
int ix;
|
||||
if (coincident(ix_, ix_close)) {
|
||||
ix = (u.x > 0) ? ix_close : ix_close - 1;
|
||||
result[0] = (u.x > 0) ? ix_close : ix_close - 1;
|
||||
} else {
|
||||
ix = std::floor(ix_);
|
||||
result[0] = std::floor(ix_);
|
||||
}
|
||||
|
||||
// Determine y index, accounting for coincidence
|
||||
double iy_ {(r.y - lower_left_.y) / pitch_.y};
|
||||
long iy_close {std::lround(iy_)};
|
||||
int iy;
|
||||
if (coincident(iy_, iy_close)) {
|
||||
iy = (u.y > 0) ? iy_close : iy_close - 1;
|
||||
result[1] = (u.y > 0) ? iy_close : iy_close - 1;
|
||||
} else {
|
||||
iy = std::floor(iy_);
|
||||
result[1] = std::floor(iy_);
|
||||
}
|
||||
|
||||
// Determine z index, accounting for coincidence
|
||||
int iz = 0;
|
||||
result[2] = 0;
|
||||
if (is_3d_) {
|
||||
double iz_ {(r.z - lower_left_.z) / pitch_.z};
|
||||
long iz_close {std::lround(iz_)};
|
||||
if (coincident(iz_, iz_close)) {
|
||||
iz = (u.z > 0) ? iz_close : iz_close - 1;
|
||||
result[2] = (u.z > 0) ? iz_close : iz_close - 1;
|
||||
} else {
|
||||
iz = std::floor(iz_);
|
||||
result[2] = std::floor(iz_);
|
||||
}
|
||||
}
|
||||
return {ix, iy, iz};
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
|
||||
Position
|
||||
RectLattice::get_local_position(Position r, const std::array<int, 3> i_xyz)
|
||||
const
|
||||
Position RectLattice::get_local_position(
|
||||
Position r, const array<int, 3>& i_xyz) const
|
||||
{
|
||||
r.x -= (lower_left_.x + (i_xyz[0] + 0.5)*pitch_.x);
|
||||
r.y -= (lower_left_.y + (i_xyz[1] + 0.5)*pitch_.y);
|
||||
|
|
@ -332,10 +330,11 @@ const
|
|||
|
||||
//==============================================================================
|
||||
|
||||
int32_t&
|
||||
RectLattice::offset(int map, const int i_xyz[3])
|
||||
int32_t& RectLattice::offset(int map, array<int, 3> const& i_xyz)
|
||||
{
|
||||
return offsets_[nx*ny*nz*map + nx*ny*i_xyz[2] + nx*i_xyz[1] + i_xyz[0]];
|
||||
return offsets_[n_cells_[0] * n_cells_[1] * n_cells_[2] * map +
|
||||
n_cells_[0] * n_cells_[1] * i_xyz[2] +
|
||||
n_cells_[0] * i_xyz[1] + i_xyz[0]];
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
|
|
@ -343,7 +342,7 @@ RectLattice::offset(int map, const int i_xyz[3])
|
|||
int32_t
|
||||
RectLattice::offset(int map, int indx) const
|
||||
{
|
||||
return offsets_[nx*ny*nz*map + indx];
|
||||
return offsets_[n_cells_[0] * n_cells_[1] * n_cells_[2] * map + indx];
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
|
|
@ -351,9 +350,9 @@ RectLattice::offset(int map, int indx) const
|
|||
std::string
|
||||
RectLattice::index_to_string(int indx) const
|
||||
{
|
||||
int iz {indx / (nx * ny)};
|
||||
int iy {(indx - nx * ny * iz) / nx};
|
||||
int ix {indx - nx * ny * iz - nx * iy};
|
||||
int iz {indx / (n_cells_[0] * n_cells_[1])};
|
||||
int iy {(indx - n_cells_[0] * n_cells_[1] * iz) / n_cells_[0]};
|
||||
int ix {indx - n_cells_[0] * n_cells_[1] * iz - n_cells_[0] * iy};
|
||||
std::string out {std::to_string(ix)};
|
||||
out += ',';
|
||||
out += std::to_string(iy);
|
||||
|
|
@ -653,8 +652,7 @@ void HexLattice::fill_lattice_y(const vector<std::string>& univ_words)
|
|||
|
||||
//==============================================================================
|
||||
|
||||
int32_t&
|
||||
HexLattice::operator[](std::array<int, 3> i_xyz)
|
||||
int32_t const& HexLattice::operator[](array<int, 3> const& i_xyz)
|
||||
{
|
||||
int indx = (2*n_rings_-1)*(2*n_rings_-1) * i_xyz[2]
|
||||
+ (2*n_rings_-1) * i_xyz[1]
|
||||
|
|
@ -672,8 +670,7 @@ ReverseLatticeIter HexLattice::rbegin()
|
|||
|
||||
//==============================================================================
|
||||
|
||||
bool
|
||||
HexLattice::are_valid_indices(const int i_xyz[3]) const
|
||||
bool HexLattice::are_valid_indices(array<int, 3> const& i_xyz) const
|
||||
{
|
||||
// Check if (x, alpha, z) indices are valid, accounting for number of rings
|
||||
return ((i_xyz[0] >= 0) && (i_xyz[1] >= 0) && (i_xyz[2] >= 0)
|
||||
|
|
@ -826,44 +823,43 @@ std::pair<double, array<int, 3>> HexLattice::distance(
|
|||
|
||||
//==============================================================================
|
||||
|
||||
std::array<int, 3>
|
||||
HexLattice::get_indices(Position r, Direction u) const
|
||||
void HexLattice::get_indices(
|
||||
Position r, Direction u, array<int, 3>& result) const
|
||||
{
|
||||
// Offset the xyz by the lattice center.
|
||||
Position r_o {r.x - center_.x, r.y - center_.y, r.z};
|
||||
if (is_3d_) {r_o.z -= center_.z;}
|
||||
|
||||
// Index the z direction, accounting for coincidence
|
||||
int iz = 0;
|
||||
result[2] = 0;
|
||||
if (is_3d_) {
|
||||
double iz_ {r_o.z / pitch_[1] + 0.5 * n_axial_};
|
||||
long iz_close {std::lround(iz_)};
|
||||
if (coincident(iz_, iz_close)) {
|
||||
iz = (u.z > 0) ? iz_close : iz_close - 1;
|
||||
result[2] = (u.z > 0) ? iz_close : iz_close - 1;
|
||||
} else {
|
||||
iz = std::floor(iz_);
|
||||
result[2] = std::floor(iz_);
|
||||
}
|
||||
}
|
||||
|
||||
int i1, i2;
|
||||
if (orientation_ == Orientation::y) {
|
||||
// Convert coordinates into skewed bases. The (x, alpha) basis is used to
|
||||
// find the index of the global coordinates to within 4 cells.
|
||||
double alpha = r_o.y - r_o.x / std::sqrt(3.0);
|
||||
i1 = std::floor(r_o.x / (0.5*std::sqrt(3.0) * pitch_[0]));
|
||||
i2 = std::floor(alpha / pitch_[0]);
|
||||
result[0] = std::floor(r_o.x / (0.5 * std::sqrt(3.0) * pitch_[0]));
|
||||
result[1] = std::floor(alpha / pitch_[0]);
|
||||
} else {
|
||||
// Convert coordinates into skewed bases. The (alpha, y) basis is used to
|
||||
// find the index of the global coordinates to within 4 cells.
|
||||
double alpha = r_o.y - r_o.x * std::sqrt(3.0);
|
||||
i1 = std::floor(-alpha / (std::sqrt(3.0) * pitch_[0]));
|
||||
i2 = std::floor(r_o.y / (0.5*std::sqrt(3.0) * pitch_[0]));
|
||||
result[0] = std::floor(-alpha / (std::sqrt(3.0) * pitch_[0]));
|
||||
result[1] = std::floor(r_o.y / (0.5 * std::sqrt(3.0) * pitch_[0]));
|
||||
}
|
||||
|
||||
// Add offset to indices (the center cell is (i1, i2) = (0, 0) but
|
||||
// the array is offset so that the indices never go below 0).
|
||||
i1 += n_rings_-1;
|
||||
i2 += n_rings_-1;
|
||||
result[0] += n_rings_ - 1;
|
||||
result[1] += n_rings_ - 1;
|
||||
|
||||
// Calculate the (squared) distance between the particle and the centers of
|
||||
// the four possible cells. Regular hexagonal tiles form a Voronoi
|
||||
|
|
@ -882,14 +878,14 @@ HexLattice::get_indices(Position r, Direction u) const
|
|||
// is kept (i.e. the cell with the lowest dot product as the vectors will be
|
||||
// completely opposed if the particle is moving directly toward the center of
|
||||
// the cell).
|
||||
int i1_chg {};
|
||||
int i2_chg {};
|
||||
int i1_chg;
|
||||
int i2_chg;
|
||||
double d_min {INFTY};
|
||||
double dp_min {INFTY};
|
||||
for (int i = 0; i < 2; i++) {
|
||||
for (int j = 0; j < 2; j++) {
|
||||
// get local coordinates
|
||||
const std::array<int, 3> i_xyz {i1 + j, i2 + i, 0};
|
||||
const array<int, 3> i_xyz {result[0] + j, result[1] + i, 0};
|
||||
Position r_t = get_local_position(r, i_xyz);
|
||||
// calculate distance
|
||||
double d = r_t.x*r_t.x + r_t.y*r_t.y;
|
||||
|
|
@ -916,17 +912,14 @@ HexLattice::get_indices(Position r, Direction u) const
|
|||
}
|
||||
|
||||
// update outgoing indices
|
||||
i1 += i1_chg;
|
||||
i2 += i2_chg;
|
||||
|
||||
return {i1, i2, iz};
|
||||
result[0] += i1_chg;
|
||||
result[1] += i2_chg;
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
|
||||
Position
|
||||
HexLattice::get_local_position(Position r, const std::array<int, 3> i_xyz)
|
||||
const
|
||||
Position HexLattice::get_local_position(
|
||||
Position r, const array<int, 3>& i_xyz) const
|
||||
{
|
||||
if (orientation_ == Orientation::y) {
|
||||
// x_l = x_g - (center + pitch_x*cos(30)*index_x)
|
||||
|
|
@ -961,14 +954,13 @@ HexLattice::is_valid_index(int indx) const
|
|||
int iz = indx / (nx * ny);
|
||||
int iy = (indx - nx*ny*iz) / nx;
|
||||
int ix = indx - nx*ny*iz - nx*iy;
|
||||
int i_xyz[3] {ix, iy, iz};
|
||||
array<int, 3> i_xyz {ix, iy, iz};
|
||||
return are_valid_indices(i_xyz);
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
|
||||
int32_t&
|
||||
HexLattice::offset(int map, const int i_xyz[3])
|
||||
int32_t& HexLattice::offset(int map, array<int, 3> const& i_xyz)
|
||||
{
|
||||
int nx {2*n_rings_ - 1};
|
||||
int ny {2*n_rings_ - 1};
|
||||
|
|
|
|||
|
|
@ -181,8 +181,8 @@ void print_particle(Particle& p)
|
|||
if (p.coord(i).lattice != C_NONE) {
|
||||
const Lattice& lat {*model::lattices[p.coord(i).lattice]};
|
||||
fmt::print(" Lattice = {}\n", lat.id_);
|
||||
fmt::print(" Lattice position = ({},{},{})\n", p.coord_[i].lattice_x,
|
||||
p.coord_[i].lattice_y, p.coord_[i].lattice_z);
|
||||
fmt::print(" Lattice position = ({},{},{})\n", p.coord(i).lattice_i[0],
|
||||
p.coord(i).lattice_i[1], p.coord(i).lattice_i[2]);
|
||||
}
|
||||
|
||||
fmt::print(" r = {}\n", p.coord(i).r);
|
||||
|
|
|
|||
|
|
@ -48,10 +48,8 @@ DistribcellFilter::get_all_bins(const Particle& p, TallyEstimator estimator,
|
|||
if (c.type_ == Fill::UNIVERSE) {
|
||||
offset += c.offset_[distribcell_index];
|
||||
} else if (c.type_ == Fill::LATTICE) {
|
||||
auto& lat {*model::lattices[p.coord_[i+1].lattice]};
|
||||
int i_xyz[3] {p.coord_[i+1].lattice_x,
|
||||
p.coord_[i+1].lattice_y,
|
||||
p.coord_[i+1].lattice_z};
|
||||
auto& lat {*model::lattices[p.coord(i + 1).lattice]};
|
||||
const auto& i_xyz {p.coord(i + 1).lattice_i};
|
||||
if (lat.are_valid_indices(i_xyz)) {
|
||||
offset += lat.offset(distribcell_index, i_xyz);
|
||||
}
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue