• Home
  • Features
  • Pricing
  • Docs
  • Announcements
  • Sign In

openmc-dev / openmc / 12954288577

24 Jan 2025 05:13PM UTC coverage: 84.833% (-0.1%) from 84.928%
12954288577

Pull #2671

github

web-flow
Merge d2ca87df5 into 560bd22bc
Pull Request #2671: Adding methods to automatically apply results to existing Tally objects.

53 of 65 new or added lines in 7 files covered. (81.54%)

59 existing lines in 20 files now uncovered.

50089 of 59044 relevant lines covered (84.83%)

34992387.74 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

77.78
/include/openmc/cell.h
1
#ifndef OPENMC_CELL_H
2
#define OPENMC_CELL_H
3

4
#include <cstdint>
5
#include <functional> // for hash
6
#include <limits>
7
#include <string>
8
#include <unordered_map>
9
#include <unordered_set>
10

11
#include "hdf5.h"
12
#include "pugixml.hpp"
13
#include <gsl/gsl-lite.hpp>
14

15
#include "openmc/bounding_box.h"
16
#include "openmc/constants.h"
17
#include "openmc/memory.h" // for unique_ptr
18
#include "openmc/neighbor_list.h"
19
#include "openmc/position.h"
20
#include "openmc/surface.h"
21
#include "openmc/universe.h"
22
#include "openmc/vector.h"
23

24
namespace openmc {
25

26
//==============================================================================
27
// Constants
28
//==============================================================================
29

30
enum class Fill { MATERIAL, UNIVERSE, LATTICE };
31

32
constexpr int32_t OP_LEFT_PAREN {std::numeric_limits<int32_t>::max()};
33
constexpr int32_t OP_RIGHT_PAREN {std::numeric_limits<int32_t>::max() - 1};
34
constexpr int32_t OP_COMPLEMENT {std::numeric_limits<int32_t>::max() - 2};
35
constexpr int32_t OP_INTERSECTION {std::numeric_limits<int32_t>::max() - 3};
36
constexpr int32_t OP_UNION {std::numeric_limits<int32_t>::max() - 4};
37

38
//==============================================================================
39
// Global variables
40
//==============================================================================
41

42
class Cell;
43
class GeometryState;
44
class ParentCell;
45
class CellInstance;
46
class Universe;
47
class UniversePartitioner;
48

49
namespace model {
50
extern std::unordered_map<int32_t, int32_t> cell_map;
51
extern vector<unique_ptr<Cell>> cells;
52

53
} // namespace model
54

55
//==============================================================================
56

57
class Region {
58
public:
59
  //----------------------------------------------------------------------------
60
  // Constructors
61
  Region() {}
31,234✔
62
  explicit Region(std::string region_spec, int32_t cell_id);
63

64
  //----------------------------------------------------------------------------
65
  // Methods
66

67
  //! \brief Determine if a cell contains the particle at a given location.
68
  //!
69
  //! The bounds of the cell are determined by a logical expression involving
70
  //! surface half-spaces. The expression used is given in infix notation
71
  //!
72
  //! The function is split into two cases, one for simple cells (those
73
  //! involving only the intersection of half-spaces) and one for complex cells.
74
  //! Both cases use short circuiting; however, in the case fo complex cells,
75
  //! the complexity increases with the binary operators involved.
76
  //! \param r The 3D Cartesian coordinate to check.
77
  //! \param u A direction used to "break ties" the coordinates are very
78
  //!   close to a surface.
79
  //! \param on_surface The signed index of a surface that the coordinate is
80
  //!   known to be on.  This index takes precedence over surface sense
81
  //!   calculations.
82
  bool contains(Position r, Direction u, int32_t on_surface) const;
83

84
  //! Find the oncoming boundary of this cell.
85
  std::pair<double, int32_t> distance(
86
    Position r, Direction u, int32_t on_surface) const;
87

88
  //! Get the BoundingBox for this cell.
89
  BoundingBox bounding_box(int32_t cell_id) const;
90

91
  //! Get the CSG expression as a string
92
  std::string str() const;
93

94
  //! Get a vector containing all the surfaces in the region expression
95
  vector<int32_t> surfaces() const;
96

97
  //----------------------------------------------------------------------------
98
  // Accessors
99

100
  //! Get Boolean of if the cell is simple or not
101
  bool is_simple() const { return simple_; }
2,147,483,647✔
102

103
private:
104
  //----------------------------------------------------------------------------
105
  // Private Methods
106

107
  //! Get a vector of the region expression in postfix notation
108
  vector<int32_t> generate_postfix(int32_t cell_id) const;
109

110
  //! Determine if a particle is inside the cell for a simple cell (only
111
  //! intersection operators)
112
  bool contains_simple(Position r, Direction u, int32_t on_surface) const;
113

114
  //! Determine if a particle is inside the cell for a complex cell.
115
  //!
116
  //! Uses the comobination of half-spaces and binary operators to determine
117
  //! if short circuiting can be used. Short cicuiting uses the relative and
118
  //! absolute depth of parenthases in the expression.
119
  bool contains_complex(Position r, Direction u, int32_t on_surface) const;
120

121
  //! BoundingBox if the paritcle is in a simple cell.
122
  BoundingBox bounding_box_simple() const;
123

124
  //! BoundingBox if the particle is in a complex cell.
125
  BoundingBox bounding_box_complex(vector<int32_t> postfix) const;
126

127
  //! Enfource precedence: Parenthases, Complement, Intersection, Union
128
  void add_precedence();
129

130
  //! Add parenthesis to enforce precedence
131
  gsl::index add_parentheses(gsl::index start);
132

133
  //! Remove complement operators from the expression
134
  void remove_complement_ops();
135

136
  //! Remove complement operators by using DeMorgan's laws
137
  void apply_demorgan(
138
    vector<int32_t>::iterator start, vector<int32_t>::iterator stop);
139

140
  //----------------------------------------------------------------------------
141
  // Private Data
142

143
  //! Definition of spatial region as Boolean expression of half-spaces
144
  // TODO: Should this be a vector of some other type
145
  vector<int32_t> expression_;
146
  bool simple_; //!< Does the region contain only intersections?
147
};
148

149
//==============================================================================
150

151
class Cell {
152
public:
153
  //----------------------------------------------------------------------------
154
  // Constructors, destructors, factory functions
155

156
  explicit Cell(pugi::xml_node cell_node);
157
  Cell() {};
31,425✔
158
  virtual ~Cell() = default;
31,423✔
UNCOV
159

×
160
  //----------------------------------------------------------------------------
31,423✔
161
  // Methods
162

163
  //! \brief Determine if a cell contains the particle at a given location.
164
  //!
165
  //! The bounds of the cell are detemined by a logical expression involving
166
  //! surface half-spaces. At initialization, the expression was converted
167
  //! to RPN notation.
168
  //!
169
  //! The function is split into two cases, one for simple cells (those
170
  //! involving only the intersection of half-spaces) and one for complex cells.
171
  //! Simple cells can be evaluated with short circuit evaluation, i.e., as soon
172
  //! as we know that one half-space is not satisfied, we can exit. This
173
  //! provides a performance benefit for the common case. In
174
  //! contains_complex, we evaluate the RPN expression using a stack, similar to
175
  //! how a RPN calculator would work.
176
  //! \param r The 3D Cartesian coordinate to check.
177
  //! \param u A direction used to "break ties" the coordinates are very
178
  //!   close to a surface.
179
  //! \param on_surface The signed index of a surface that the coordinate is
180
  //!   known to be on.  This index takes precedence over surface sense
181
  //!   calculations.
182
  virtual bool contains(Position r, Direction u, int32_t on_surface) const = 0;
183

184
  //! Find the oncoming boundary of this cell.
185
  virtual std::pair<double, int32_t> distance(
186
    Position r, Direction u, int32_t on_surface, GeometryState* p) const = 0;
187

188
  //! Write all information needed to reconstruct the cell to an HDF5 group.
189
  //! \param group_id An HDF5 group id.
190
  void to_hdf5(hid_t group_id) const;
191

192
  virtual void to_hdf5_inner(hid_t group_id) const = 0;
193

194
  //! Export physical properties to HDF5
195
  //! \param[in] group  HDF5 group to read from
196
  void export_properties_hdf5(hid_t group) const;
197

198
  //! Import physical properties from HDF5
199
  //! \param[in] group  HDF5 group to write to
200
  void import_properties_hdf5(hid_t group);
201

202
  //! Get the BoundingBox for this cell.
203
  virtual BoundingBox bounding_box() const = 0;
204

205
  //! Get a vector of surfaces in the cell
206
  virtual vector<int32_t> surfaces() const { return vector<int32_t>(); }
207

UNCOV
208
  //! Check if the cell region expression is simple
×
209
  virtual bool is_simple() const { return true; }
210

211
  //----------------------------------------------------------------------------
92,429✔
212
  // Accessors
213

214
  //! Get the temperature of a cell instance
215
  //! \param[in] instance Instance index. If -1 is given, the temperature for
216
  //!   the first instance is returned.
217
  //! \return Temperature in [K]
218
  double temperature(int32_t instance = -1) const;
219

220
  //! Set the temperature of a cell instance
221
  //! \param[in] T Temperature in [K]
222
  //! \param[in] instance Instance index. If -1 is given, the temperature for
223
  //!   all instances is set.
224
  //! \param[in] set_contained If this cell is not filled with a material,
225
  //!   collect all contained cells with material fills and set their
226
  //!   temperatures.
227
  void set_temperature(
228
    double T, int32_t instance = -1, bool set_contained = false);
229

230
  //! Set the rotation matrix of a cell instance
231
  //! \param[in] rot The rotation matrix of length 3 or 9
232
  void set_rotation(const vector<double>& rot);
233

234
  //! Get the name of a cell
235
  //! \return Cell name
236
  const std::string& name() const { return name_; };
237

238
  //! Set the temperature of a cell instance
314✔
239
  //! \param[in] name Cell name
240
  void set_name(const std::string& name) { name_ = name; };
241

242
  //! Get all cell instances contained by this cell
12✔
243
  //! \param[in] instance Instance of the cell for which to get contained cells
244
  //! (default instance is zero)
245
  //! \param[in] hint positional hint for determining the parent cells
246
  //! \return Map with cell indexes as keys and
247
  //! instances as values
248
  std::unordered_map<int32_t, vector<int32_t>> get_contained_cells(
249
    int32_t instance = 0, Position* hint = nullptr) const;
250

251
  //! Determine the material index corresponding to a specific cell instance,
252
  //! taking into account presence of distribcell material
253
  //! \param[in] instance of the cell
254
  //! \return material index
255
  int32_t material(int32_t instance) const
×
256
  {
257
    // If distributed materials are used, then each instance has its own
258
    // material definition. If distributed materials are not used, then
259
    // all instances used the same material stored at material_[0]. The
260
    // presence of distributed materials is inferred from the size of
261
    // the material_ vector being greater than one.
262
    if (material_.size() > 1) {
×
263
      return material_[instance];
×
264
    } else {
265
      return material_[0];
×
266
    }
267
  }
268

269
  //! Determine the temperature index corresponding to a specific cell instance,
270
  //! taking into account presence of distribcell temperature
271
  //! \param[in] instance of the cell
272
  //! \return temperature index
273
  double sqrtkT(int32_t instance) const
274
  {
275
    // If distributed materials are used, then each instance has its own
276
    // temperature definition. If distributed materials are not used, then
277
    // all instances used the same temperature stored at sqrtkT_[0]. The
278
    // presence of distributed materials is inferred from the size of
279
    // the sqrtkT_ vector being greater than one.
280
    if (sqrtkT_.size() > 1) {
281
      return sqrtkT_[instance];
282
    } else {
283
      return sqrtkT_[0];
284
    }
285
  }
286

287
protected:
288
  //! Determine the path to this cell instance in the geometry hierarchy
289
  //! \param[in] instance of the cell to find parent cells for
290
  //! \param[in] r position used to do a fast search for parent cells
291
  //! \return parent cells
292
  vector<ParentCell> find_parent_cells(
293
    int32_t instance, const Position& r) const;
294

295
  //! Determine the path to this cell instance in the geometry hierarchy
296
  //! \param[in] instance of the cell to find parent cells for
297
  //! \param[in] p particle used to do a fast search for parent cells
298
  //! \return parent cells
299
  vector<ParentCell> find_parent_cells(
300
    int32_t instance, GeometryState& p) const;
301

302
  //! Determine the path to this cell instance in the geometry hierarchy
303
  //! \param[in] instance of the cell to find parent cells for
304
  //! \return parent cells
305
  vector<ParentCell> exhaustive_find_parent_cells(int32_t instance) const;
306

307
  //! Inner function for retrieving contained cells
308
  void get_contained_cells_inner(
309
    std::unordered_map<int32_t, vector<int32_t>>& contained_cells,
310
    vector<ParentCell>& parent_cells) const;
311

312
public:
313
  //----------------------------------------------------------------------------
314
  // Data members
315

316
  int32_t id_;              //!< Unique ID
317
  std::string name_;        //!< User-defined name
318
  Fill type_;               //!< Material, universe, or lattice
319
  int32_t universe_;        //!< Universe # this cell is in
320
  int32_t fill_;            //!< Universe # filling this cell
321
  int32_t n_instances_ {0}; //!< Number of instances of this cell
322

323
  //! \brief Index corresponding to this cell in distribcell arrays
324
  int distribcell_index_ {C_NONE};
325

326
  //! \brief Material(s) within this cell.
327
  //!
328
  //! May be multiple materials for distribcell.
329
  vector<int32_t> material_;
330

331
  //! \brief Temperature(s) within this cell.
332
  //!
333
  //! The stored values are actually sqrt(k_Boltzmann * T) for each temperature
334
  //! T. The units are sqrt(eV).
335
  vector<double> sqrtkT_;
336

337
  //! \brief Neighboring cells in the same universe.
338
  NeighborList neighbors_;
339

340
  Position translation_ {0, 0, 0}; //!< Translation vector for filled universe
341

342
  //! \brief Rotational tranfsormation of the filled universe.
343
  //
344
  //! The vector is empty if there is no rotation. Otherwise, the first 9 values
345
  //! give the rotation matrix in row-major order. When the user specifies
346
  //! rotation angles about the x-, y- and z- axes in degrees, these values are
347
  //! also present at the end of the vector, making it of length 12.
348
  vector<double> rotation_;
349

350
  vector<int32_t> offset_; //!< Distribcell offset table
351

352
  // Accessors
353
  const GeometryType& geom_type() const { return geom_type_; }
354
  GeometryType& geom_type() { return geom_type_; }
355

356
private:
31,470✔
357
  GeometryType geom_type_; //!< Geometric representation type (CSG, DAGMC)
358
};
359

360
struct CellInstanceItem {
361
  int32_t index {-1};    //! Index into global cells array
362
  int lattice_indx {-1}; //! Flat index value of the lattice cell
363
};
364

365
//==============================================================================
366

367
class CSGCell : public Cell {
368
public:
369
  //----------------------------------------------------------------------------
370
  // Constructors
371
  CSGCell();
372
  explicit CSGCell(pugi::xml_node cell_node);
373

374
  //----------------------------------------------------------------------------
375
  // Methods
376
  vector<int32_t> surfaces() const override { return region_.surfaces(); }
377

378
  std::pair<double, int32_t> distance(Position r, Direction u,
5,875✔
379
    int32_t on_surface, GeometryState* p) const override
380
  {
2,147,483,647✔
381
    return region_.distance(r, u, on_surface);
382
  }
383

2,147,483,647✔
384
  bool contains(Position r, Direction u, int32_t on_surface) const override
385
  {
386
    return region_.contains(r, u, on_surface);
2,147,483,647✔
387
  }
388

2,147,483,647✔
389
  BoundingBox bounding_box() const override
390
  {
391
    return region_.bounding_box(id_);
96✔
392
  }
393

96✔
394
  void to_hdf5_inner(hid_t group_id) const override;
395

396
  bool is_simple() const override { return region_.is_simple(); }
397

398
protected:
2,147,483,647✔
399
  //! Returns the beginning position of a parenthesis block (immediately before
400
  //! two surface tokens) in the RPN given a starting position at the end of
401
  //! that block (immediately after two surface tokens)
402
  //! \param start Starting position of the search
403
  //! \param rpn The rpn being searched
404
  static vector<int32_t>::iterator find_left_parenthesis(
405
    vector<int32_t>::iterator start, const vector<int32_t>& rpn);
406

407
private:
408
  Region region_;
409
};
410

411
//==============================================================================
412
//! Define an instance of a particular cell
413
//==============================================================================
414

415
//!  Stores information used to identify a unique cell in the model
416
struct CellInstance {
417
  //! Check for equality
418
  bool operator==(const CellInstance& other) const
4,254,520✔
419
  {
420
    return index_cell == other.index_cell && instance == other.instance;
4,254,520✔
421
  }
422

423
  gsl::index index_cell;
424
  gsl::index instance;
425
};
426

427
//! Structure necessary for inserting CellInstance into hashed STL data
428
//! structures
429
struct CellInstanceHash {
430
  std::size_t operator()(const CellInstance& k) const
4,255,210✔
431
  {
432
    return 4096 * k.index_cell + k.instance;
4,255,210✔
433
  }
434
};
435

436
//==============================================================================
437
// Non-member functions
438
//==============================================================================
439

440
void read_cells(pugi::xml_node node);
441

442
//!  Add cells to universes
443
void populate_universes();
444

445
} // namespace openmc
446
#endif // OPENMC_CELL_H
STATUS · Troubleshooting · Open an Issue · Sales · Support · CAREERS · ENTERPRISE · START FREE · SCHEDULE DEMO
ANNOUNCEMENTS · TWITTER · TOS & SLA · Supported CI Services · What's a CI service? · Automated Testing

© 2026 Coveralls, Inc