mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-21 06:25:30 -04:00
Refactor Ray class into its own file (#3845)
Co-authored-by: Paul Romano <paul.k.romano@gmail.com>
This commit is contained in:
GuySten
GitHub
parent
6f72619729
commit
be4148ad0d
5 changed files with 220 additions and 200 deletions
|
|
@ -412,6 +412,7 @@ list(APPEND libopenmc_SOURCES
|
|||
src/random_ray/linear_source_domain.cpp
|
||||
src/random_ray/moment_matrix.cpp
|
||||
src/random_ray/source_region.cpp
|
||||
src/ray.cpp
|
||||
src/reaction.cpp
|
||||
src/reaction_product.cpp
|
||||
src/scattdata.cpp
|
||||
|
|
|
|||
|
|
@ -17,6 +17,7 @@
|
|||
#include "openmc/particle.h"
|
||||
#include "openmc/position.h"
|
||||
#include "openmc/random_lcg.h"
|
||||
#include "openmc/ray.h"
|
||||
#include "openmc/xml_interface.h"
|
||||
|
||||
namespace openmc {
|
||||
|
|
@ -497,47 +498,6 @@ private:
|
|||
Position light_location_;
|
||||
};
|
||||
|
||||
// Base class that implements ray tracing logic, not necessarily through
|
||||
// defined regions of the geometry but also outside of it.
|
||||
class Ray : public GeometryState {
|
||||
|
||||
public:
|
||||
// Initialize from location and direction
|
||||
Ray(Position r, Direction u) { init_from_r_u(r, u); }
|
||||
|
||||
// Initialize from known geometry state
|
||||
Ray(const GeometryState& p) : GeometryState(p) {}
|
||||
|
||||
// Called at every surface intersection within the model
|
||||
virtual void on_intersection() = 0;
|
||||
|
||||
/*
|
||||
* Traces the ray through the geometry, calling on_intersection
|
||||
* at every surface boundary.
|
||||
*/
|
||||
void trace();
|
||||
|
||||
// Stops the ray and exits tracing when called from on_intersection
|
||||
void stop() { stop_ = true; }
|
||||
|
||||
// Sets the dist_ variable
|
||||
void compute_distance();
|
||||
|
||||
protected:
|
||||
// Records how far the ray has traveled
|
||||
double traversal_distance_ {0.0};
|
||||
|
||||
private:
|
||||
// Max intersections before we assume ray tracing is caught in an infinite
|
||||
// loop:
|
||||
static const int MAX_INTERSECTIONS = 1000000;
|
||||
|
||||
bool hit_something_ {false};
|
||||
bool stop_ {false};
|
||||
|
||||
unsigned event_counter_ {0};
|
||||
};
|
||||
|
||||
class ProjectionRay : public Ray {
|
||||
public:
|
||||
ProjectionRay(Position r, Direction u, const WireframeRayTracePlot& plot,
|
||||
|
|
|
|||
50
include/openmc/ray.h
Normal file
50
include/openmc/ray.h
Normal file
|
|
@ -0,0 +1,50 @@
|
|||
#ifndef OPENMC_RAY_H
|
||||
#define OPENMC_RAY_H
|
||||
|
||||
#include "openmc/particle_data.h"
|
||||
#include "openmc/position.h"
|
||||
|
||||
namespace openmc {
|
||||
|
||||
// Base class that implements ray tracing logic, not necessarily through
|
||||
// defined regions of the geometry but also outside of it.
|
||||
class Ray : public GeometryState {
|
||||
|
||||
public:
|
||||
// Initialize from location and direction
|
||||
Ray(Position r, Direction u) { init_from_r_u(r, u); }
|
||||
|
||||
// Initialize from known geometry state
|
||||
Ray(const GeometryState& p) : GeometryState(p) {}
|
||||
|
||||
// Called at every surface intersection within the model
|
||||
virtual void on_intersection() = 0;
|
||||
|
||||
/*
|
||||
* Traces the ray through the geometry, calling on_intersection
|
||||
* at every surface boundary.
|
||||
*/
|
||||
void trace();
|
||||
|
||||
// Stops the ray and exits tracing when called from on_intersection
|
||||
void stop() { stop_ = true; }
|
||||
|
||||
// Sets the dist_ variable
|
||||
void compute_distance();
|
||||
|
||||
protected:
|
||||
// Records how far the ray has traveled
|
||||
double traversal_distance_ {0.0};
|
||||
|
||||
private:
|
||||
// Max intersections before we assume ray tracing is caught in an infinite
|
||||
// loop:
|
||||
static const int MAX_INTERSECTIONS = 1000000;
|
||||
|
||||
bool stop_ {false};
|
||||
|
||||
unsigned event_counter_ {0};
|
||||
};
|
||||
|
||||
} // namespace openmc
|
||||
#endif // OPENMC_RAY_H
|
||||
159
src/plot.cpp
159
src/plot.cpp
|
|
@ -1649,165 +1649,6 @@ void SolidRayTracePlot::set_diffuse_fraction(pugi::xml_node node)
|
|||
}
|
||||
}
|
||||
|
||||
void Ray::compute_distance()
|
||||
{
|
||||
boundary() = distance_to_boundary(*this);
|
||||
}
|
||||
|
||||
void Ray::trace()
|
||||
{
|
||||
// To trace the ray from its origin all the way through the model, we have
|
||||
// to proceed in two phases. In the first, the ray may or may not be found
|
||||
// inside the model. If the ray is already in the model, phase one can be
|
||||
// skipped. Otherwise, the ray has to be advanced to the boundary of the
|
||||
// model where all the cells are defined. Importantly, this is assuming that
|
||||
// the model is convex, which is a very reasonable assumption for any
|
||||
// radiation transport model.
|
||||
//
|
||||
// After phase one is done, we can starting tracing from cell to cell within
|
||||
// the model. This step can use neighbor lists to accelerate the ray tracing.
|
||||
|
||||
bool inside_cell;
|
||||
// Check for location if the particle is already known
|
||||
if (lowest_coord().cell() == C_NONE) {
|
||||
// The geometry position of the particle is either unknown or outside of the
|
||||
// edge of the model.
|
||||
if (lowest_coord().universe() == C_NONE) {
|
||||
// Attempt to initialize the particle. We may have to
|
||||
// enter a loop to move it up to the edge of the model.
|
||||
inside_cell = exhaustive_find_cell(*this, settings::verbosity >= 10);
|
||||
} else {
|
||||
// It has been already calculated that the current position is outside of
|
||||
// the edge of the model.
|
||||
inside_cell = false;
|
||||
}
|
||||
} else {
|
||||
// Availability of the cell means that the particle is located inside the
|
||||
// edge.
|
||||
inside_cell = true;
|
||||
}
|
||||
|
||||
// Advance to the boundary of the model
|
||||
while (!inside_cell) {
|
||||
advance_to_boundary_from_void();
|
||||
inside_cell = exhaustive_find_cell(*this, settings::verbosity >= 10);
|
||||
|
||||
// If true this means no surface was intersected. See cell.cpp and search
|
||||
// for numeric_limits to see where we return it.
|
||||
if (surface() == std::numeric_limits<int>::max()) {
|
||||
warning(fmt::format("Lost a ray, r = {}, u = {}", r(), u()));
|
||||
return;
|
||||
}
|
||||
|
||||
// Exit this loop and enter into cell-to-cell ray tracing (which uses
|
||||
// neighbor lists)
|
||||
if (inside_cell)
|
||||
break;
|
||||
|
||||
// if there is no intersection with the model, we're done
|
||||
if (boundary().surface() == SURFACE_NONE)
|
||||
return;
|
||||
|
||||
event_counter_++;
|
||||
if (event_counter_ > MAX_INTERSECTIONS) {
|
||||
warning("Likely infinite loop in ray traced plot");
|
||||
return;
|
||||
}
|
||||
}
|
||||
|
||||
// Call the specialized logic for this type of ray. This is for the
|
||||
// intersection for the first intersection if we had one.
|
||||
if (boundary().surface() != SURFACE_NONE) {
|
||||
// set the geometry state's surface attribute to be used for
|
||||
// surface normal computation
|
||||
surface() = boundary().surface();
|
||||
on_intersection();
|
||||
if (stop_)
|
||||
return;
|
||||
}
|
||||
|
||||
// reset surface attribute to zero after the first intersection so that it
|
||||
// doesn't perturb surface crossing logic from here on out
|
||||
surface() = 0;
|
||||
|
||||
// This is the ray tracing loop within the model. It exits after exiting
|
||||
// the model, which is equivalent to assuming that the model is convex.
|
||||
// It would be nice to factor out the on_intersection at the end of this
|
||||
// loop and then do "while (inside_cell)", but we can't guarantee it's
|
||||
// on a surface in that case. There might be some other way to set it
|
||||
// up that is perhaps a little more elegant, but this is what works just
|
||||
// fine.
|
||||
while (true) {
|
||||
|
||||
compute_distance();
|
||||
|
||||
// There are no more intersections to process
|
||||
// if we hit the edge of the model, so stop
|
||||
// the particle in that case. Also, just exit
|
||||
// if a negative distance was somehow computed.
|
||||
if (boundary().distance() == INFTY || boundary().distance() == INFINITY ||
|
||||
boundary().distance() < 0) {
|
||||
return;
|
||||
}
|
||||
|
||||
// See below comment where call_on_intersection is checked in an
|
||||
// if statement for an explanation of this.
|
||||
bool call_on_intersection {true};
|
||||
if (boundary().distance() < 10 * TINY_BIT) {
|
||||
call_on_intersection = false;
|
||||
}
|
||||
|
||||
// DAGMC surfaces expect us to go a little bit further than the advance
|
||||
// distance to properly check cell inclusion.
|
||||
boundary().distance() += TINY_BIT;
|
||||
|
||||
// Advance particle, prepare for next intersection
|
||||
for (int lev = 0; lev < n_coord(); ++lev) {
|
||||
coord(lev).r() += boundary().distance() * coord(lev).u();
|
||||
}
|
||||
surface() = boundary().surface();
|
||||
// Initialize last cells from the current cell, because the cell() variable
|
||||
// does not contain the data for the case of a single-segment ray
|
||||
for (int j = 0; j < n_coord(); ++j) {
|
||||
cell_last(j) = coord(j).cell();
|
||||
}
|
||||
n_coord_last() = n_coord();
|
||||
n_coord() = boundary().coord_level();
|
||||
if (boundary().lattice_translation()[0] != 0 ||
|
||||
boundary().lattice_translation()[1] != 0 ||
|
||||
boundary().lattice_translation()[2] != 0) {
|
||||
cross_lattice(*this, boundary(), settings::verbosity >= 10);
|
||||
}
|
||||
|
||||
// Record how far the ray has traveled
|
||||
traversal_distance_ += boundary().distance();
|
||||
inside_cell = neighbor_list_find_cell(*this, settings::verbosity >= 10);
|
||||
|
||||
// Call the specialized logic for this type of ray. Note that we do not
|
||||
// call this if the advance distance is very small. Unfortunately, it seems
|
||||
// darn near impossible to get the particle advanced to the model boundary
|
||||
// and through it without sometimes accidentally calling on_intersection
|
||||
// twice. This incorrectly shades the region as occluded when it might not
|
||||
// actually be. By screening out intersection distances smaller than a
|
||||
// threshold 10x larger than the scoot distance used to advance up to the
|
||||
// model boundary, we can avoid that situation.
|
||||
if (call_on_intersection) {
|
||||
on_intersection();
|
||||
if (stop_)
|
||||
return;
|
||||
}
|
||||
|
||||
if (!inside_cell)
|
||||
return;
|
||||
|
||||
event_counter_++;
|
||||
if (event_counter_ > MAX_INTERSECTIONS) {
|
||||
warning("Likely infinite loop in ray traced plot");
|
||||
return;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
void ProjectionRay::on_intersection()
|
||||
{
|
||||
// This records a tuple with the following info
|
||||
|
|
|
|||
168
src/ray.cpp
Normal file
168
src/ray.cpp
Normal file
|
|
@ -0,0 +1,168 @@
|
|||
#include "openmc/ray.h"
|
||||
|
||||
#include "openmc/error.h"
|
||||
#include "openmc/geometry.h"
|
||||
#include "openmc/settings.h"
|
||||
|
||||
namespace openmc {
|
||||
|
||||
void Ray::compute_distance()
|
||||
{
|
||||
boundary() = distance_to_boundary(*this);
|
||||
}
|
||||
|
||||
void Ray::trace()
|
||||
{
|
||||
// To trace the ray from its origin all the way through the model, we have
|
||||
// to proceed in two phases. In the first, the ray may or may not be found
|
||||
// inside the model. If the ray is already in the model, phase one can be
|
||||
// skipped. Otherwise, the ray has to be advanced to the boundary of the
|
||||
// model where all the cells are defined. Importantly, this is assuming that
|
||||
// the model is convex, which is a very reasonable assumption for any
|
||||
// radiation transport model.
|
||||
//
|
||||
// After phase one is done, we can starting tracing from cell to cell within
|
||||
// the model. This step can use neighbor lists to accelerate the ray tracing.
|
||||
|
||||
bool inside_cell;
|
||||
// Check for location if the particle is already known
|
||||
if (lowest_coord().cell() == C_NONE) {
|
||||
// The geometry position of the particle is either unknown or outside of the
|
||||
// edge of the model.
|
||||
if (lowest_coord().universe() == C_NONE) {
|
||||
// Attempt to initialize the particle. We may have to
|
||||
// enter a loop to move it up to the edge of the model.
|
||||
inside_cell = exhaustive_find_cell(*this, settings::verbosity >= 10);
|
||||
} else {
|
||||
// It has been already calculated that the current position is outside of
|
||||
// the edge of the model.
|
||||
inside_cell = false;
|
||||
}
|
||||
} else {
|
||||
// Availability of the cell means that the particle is located inside the
|
||||
// edge.
|
||||
inside_cell = true;
|
||||
}
|
||||
|
||||
// Advance to the boundary of the model
|
||||
while (!inside_cell) {
|
||||
advance_to_boundary_from_void();
|
||||
inside_cell = exhaustive_find_cell(*this, settings::verbosity >= 10);
|
||||
|
||||
// If true this means no surface was intersected. See cell.cpp and search
|
||||
// for numeric_limits to see where we return it.
|
||||
if (surface() == std::numeric_limits<int>::max()) {
|
||||
warning(fmt::format("Lost a ray, r = {}, u = {}", r(), u()));
|
||||
return;
|
||||
}
|
||||
|
||||
// Exit this loop and enter into cell-to-cell ray tracing (which uses
|
||||
// neighbor lists)
|
||||
if (inside_cell)
|
||||
break;
|
||||
|
||||
// if there is no intersection with the model, we're done
|
||||
if (boundary().surface() == SURFACE_NONE)
|
||||
return;
|
||||
|
||||
event_counter_++;
|
||||
if (event_counter_ > MAX_INTERSECTIONS) {
|
||||
warning("Likely infinite loop in ray traced plot");
|
||||
return;
|
||||
}
|
||||
}
|
||||
|
||||
// Call the specialized logic for this type of ray. This is for the
|
||||
// intersection for the first intersection if we had one.
|
||||
if (boundary().surface() != SURFACE_NONE) {
|
||||
// set the geometry state's surface attribute to be used for
|
||||
// surface normal computation
|
||||
surface() = boundary().surface();
|
||||
on_intersection();
|
||||
if (stop_)
|
||||
return;
|
||||
}
|
||||
|
||||
// reset surface attribute to zero after the first intersection so that it
|
||||
// doesn't perturb surface crossing logic from here on out
|
||||
surface() = 0;
|
||||
|
||||
// This is the ray tracing loop within the model. It exits after exiting
|
||||
// the model, which is equivalent to assuming that the model is convex.
|
||||
// It would be nice to factor out the on_intersection at the end of this
|
||||
// loop and then do "while (inside_cell)", but we can't guarantee it's
|
||||
// on a surface in that case. There might be some other way to set it
|
||||
// up that is perhaps a little more elegant, but this is what works just
|
||||
// fine.
|
||||
while (true) {
|
||||
|
||||
compute_distance();
|
||||
|
||||
// There are no more intersections to process
|
||||
// if we hit the edge of the model, so stop
|
||||
// the particle in that case. Also, just exit
|
||||
// if a negative distance was somehow computed.
|
||||
if (boundary().distance() == INFTY || boundary().distance() == INFINITY ||
|
||||
boundary().distance() < 0) {
|
||||
return;
|
||||
}
|
||||
|
||||
// See below comment where call_on_intersection is checked in an
|
||||
// if statement for an explanation of this.
|
||||
bool call_on_intersection {true};
|
||||
if (boundary().distance() < 10 * TINY_BIT) {
|
||||
call_on_intersection = false;
|
||||
}
|
||||
|
||||
// DAGMC surfaces expect us to go a little bit further than the advance
|
||||
// distance to properly check cell inclusion.
|
||||
boundary().distance() += TINY_BIT;
|
||||
|
||||
// Advance particle, prepare for next intersection
|
||||
for (int lev = 0; lev < n_coord(); ++lev) {
|
||||
coord(lev).r() += boundary().distance() * coord(lev).u();
|
||||
}
|
||||
surface() = boundary().surface();
|
||||
// Initialize last cells from the current cell, because the cell() variable
|
||||
// does not contain the data for the case of a single-segment ray
|
||||
for (int j = 0; j < n_coord(); ++j) {
|
||||
cell_last(j) = coord(j).cell();
|
||||
}
|
||||
n_coord_last() = n_coord();
|
||||
n_coord() = boundary().coord_level();
|
||||
if (boundary().lattice_translation()[0] != 0 ||
|
||||
boundary().lattice_translation()[1] != 0 ||
|
||||
boundary().lattice_translation()[2] != 0) {
|
||||
cross_lattice(*this, boundary(), settings::verbosity >= 10);
|
||||
}
|
||||
|
||||
// Record how far the ray has traveled
|
||||
traversal_distance_ += boundary().distance();
|
||||
inside_cell = neighbor_list_find_cell(*this, settings::verbosity >= 10);
|
||||
|
||||
// Call the specialized logic for this type of ray. Note that we do not
|
||||
// call this if the advance distance is very small. Unfortunately, it seems
|
||||
// darn near impossible to get the particle advanced to the model boundary
|
||||
// and through it without sometimes accidentally calling on_intersection
|
||||
// twice. This incorrectly shades the region as occluded when it might not
|
||||
// actually be. By screening out intersection distances smaller than a
|
||||
// threshold 10x larger than the scoot distance used to advance up to the
|
||||
// model boundary, we can avoid that situation.
|
||||
if (call_on_intersection) {
|
||||
on_intersection();
|
||||
if (stop_)
|
||||
return;
|
||||
}
|
||||
|
||||
if (!inside_cell)
|
||||
return;
|
||||
|
||||
event_counter_++;
|
||||
if (event_counter_ > MAX_INTERSECTIONS) {
|
||||
warning("Likely infinite loop in ray traced plot");
|
||||
return;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
} // namespace openmc
|
||||
Loading…
Add table
Add a link
Reference in a new issue