OpenMCCellAverageProblem.h

Go to the documentation of this file.

/********************************************************************/
/*                  SOFTWARE COPYRIGHT NOTIFICATION                 */
/*                             Cardinal                             */
/*                                                                  */
/*                  (c) 2021 UChicago Argonne, LLC                  */
/*                        ALL RIGHTS RESERVED                       */
/*                                                                  */
/*                 Prepared by UChicago Argonne, LLC                */
/*               Under Contract No. DE-AC02-06CH11357               */
/*                With the U. S. Department of Energy               */
/*                                                                  */
/*             Prepared by Battelle Energy Alliance, LLC            */
/*               Under Contract No. DE-AC07-05ID14517               */
/*                With the U. S. Department of Energy               */
/*                                                                  */
/*                 See LICENSE for full restrictions                */
/********************************************************************/
 
#pragma once
 
#include "OpenMCProblemBase.h"
#include "SymmetryPointGenerator.h"
 
#include "TallyBase.h"
#include "FilterBase.h"
 
#ifdef ENABLE_DAGMC
#include "MoabSkinner.h"
#include "DagMC.hpp"
#endif
 
class OpenMCVolumeCalculation;
 
class OpenMCCellAverageProblem : public OpenMCProblemBase
{
public:
  OpenMCCellAverageProblem(const InputParameters & params);
  static InputParameters validParams();
 
  virtual void initialSetup() override;
  virtual void externalSolve() override;
  virtual void syncSolutions(ExternalProblem::Direction direction) override;
  virtual bool converged(unsigned int) override { return true; }
 
  void read2DBlockParameters(const std::string name,
                             std::vector<std::vector<SubdomainName>> & names,
                             std::vector<SubdomainID> & flattened_ids);
 
  void checkBlocksInMesh(const std::string name,
                         const std::vector<SubdomainID> & ids,
                         const std::vector<SubdomainName> & names) const;
 
  void setupProblem();
 
  virtual void addExternalVariables() override;
 
  virtual Real cellVolume(const cellInfo & cell_info) const;
 
  virtual const OpenMCVolumeCalculation * volumeCalculation() const { return _volume_calc; }
 
  virtual const std::map<cellInfo, std::vector<unsigned int>> & cellToElem() const
  {
    return _cell_to_elem;
  }
  virtual const std::map<cellInfo, std::vector<unsigned int>> & cellToElem() const {…}
 
  virtual std::unordered_set<SubdomainID> getCellToElementSub(const cellInfo & info)
  {
    return _cell_to_elem_subdomain.at(info);
  }
  virtual std::unordered_set<SubdomainID> getCellToElementSub(const cellInfo & info) {…}
 
  virtual bool hasPointTransformations() const { return _symmetry != nullptr; }
 
  virtual const std::vector<std::string> & getTallyScores() const { return _all_tally_scores; }
 
  bool hasScore(const std::string & score)
  {
    return std::find(_all_tally_scores.begin(), _all_tally_scores.end(), score) !=
           _all_tally_scores.end();
  }
  bool hasScore(const std::string & score) {…}
 
  std::vector<const TallyBase *> getTalliesByScore(const std::string & score);
 
  std::vector<const MooseVariableFE<Real> *> getTallyScoreVariables(const std::string & score,
                                                                    THREAD_ID tid);
 
  std::vector<const VariableValue *> getTallyScoreVariableValues(const std::string & score,
                                                                 THREAD_ID tid);
 
  std::vector<const VariableValue *> getTallyScoreNeighborVariableValues(const std::string & score,
                                                                         THREAD_ID tid);
 
  virtual Point transformPoint(const Point & pt) const
  {
    return this->hasPointTransformations() ? _symmetry->transformPoint(pt) : pt;
  }
  virtual Point transformPoint(const Point & pt) const {…}
 
  virtual void createQRules(QuadratureType type,
                            Order order,
                            Order volume_order,
                            Order face_order,
                            SubdomainID block,
                            bool allow_negative_weights = true) override;
 
  typedef std::unordered_map<int32_t, std::vector<int32_t>> containedCells;
 
  int32_t elemToCellIndex(const int & elem_id) const { return elemToCellInfo(elem_id).first; }
 
  int32_t elemToCellID(const int & elem_id) const { return cellID(elemToCellIndex(elem_id)); }
 
  int32_t elemToCellInstance(const int & elem_id) const { return elemToCellInfo(elem_id).second; }
 
  cellInfo elemToCellInfo(const int & elem_id) const { return _elem_to_cell[elem_id]; }
 
  int32_t cellToMaterialIndex(const cellInfo & cell_info) const
  {
    return _cell_to_material.at(cell_info);
  }
  int32_t cellToMaterialIndex(const cellInfo & cell_info) const {…}
 
  coupling::CouplingFields cellFeedback(const cellInfo & cell_info) const;
 
  bool hasDensityFeedback(const cellInfo & cell_info) const
  {
    std::vector<coupling::CouplingFields> phase = {coupling::density,
                                                   coupling::density_and_temperature};
    return std::find(phase.begin(), phase.end(), cellFeedback(cell_info)) != phase.end();
  }
  bool hasDensityFeedback(const cellInfo & cell_info) const {…}
 
  bool hasTemperatureFeedback(const cellInfo & cell_info) const
  {
    std::vector<coupling::CouplingFields> phase = {coupling::temperature,
                                                   coupling::density_and_temperature};
    return std::find(phase.begin(), phase.end(), cellFeedback(cell_info)) != phase.end();
  }
  bool hasTemperatureFeedback(const cellInfo & cell_info) const {…}
 
  bool hasFilter(const std::string & filter_name) const { return _filters.count(filter_name) > 0; }
 
  std::shared_ptr<FilterBase> & getFilter(const std::string & filter_name)
  {
    return _filters.at(filter_name);
  }
  std::shared_ptr<FilterBase> & getFilter(const std::string & filter_name) {…}
 
  const std::vector<std::shared_ptr<TallyBase>> & getLocalTally() const { return _local_tallies; }
 
  double cellTemperature(const cellInfo & cell_info) const;
 
  double cellMappedVolume(const cellInfo & cell_info) const;
 
  void reloadDAGMC();
 
  void addFilter(const std::string & type,
                 const std::string & name,
                 InputParameters & moose_object_pars);
 
  void
  addTally(const std::string & type, const std::string & name, InputParameters & moose_object_pars);
 
  Real tallyMultiplier(unsigned int global_score) const;
 
  template <typename T>
  void checkEmptyVector(const std::vector<T> & vector, const std::string & name) const
  {
    if (vector.empty())
      mooseError(name + " cannot be empty!");
  }
  void checkEmptyVector(const std::vector<T> & vector, const std::string & name) const {…}
 
  int fixedPointIteration() const { return _fixed_point_iteration; }
 
  bool hasAdaptivity() const { return _has_adaptivity; }
 
  static constexpr int32_t UNMAPPED{-1};
 
  static constexpr int DIMENSION{3};
 
protected:
  unsigned int getCellLevel(const Point & c) const;
 
  void
  readBlockVariables(const std::string & param,
                     const std::string & default_name,
                     std::map<std::string, std::vector<SubdomainName>> & vars_to_specified_blocks,
                     std::vector<SubdomainID> & specified_blocks);
 
  bool cellHasIdenticalFill(const cellInfo & cell_info) const;
 
  containedCells shiftCellInstances(const cellInfo & cell_info) const;
 
  bool cellMapsToSubdomain(const cellInfo & cell_info,
                           const std::unordered_set<SubdomainID> & id) const;
 
  cellInfo firstContainedMaterialCell(const cellInfo & cell_info) const;
 
  containedCells containedMaterialCells(const cellInfo & cell_info) const;
 
  void updateOpenMCGeometry();
 
  virtual void updateMaterials();
 
  std::vector<std::string> getMaterialInEachSubdomain() const;
 
  Point transformPointToOpenMC(const Point & pt) const;
 
  void checkNormalization(const Real & sum, unsigned int global_score) const;
 
  void
  printTrisoHelp(const std::chrono::time_point<std::chrono::high_resolution_clock> & start) const;
 
  void printAuxVariableIO();
 
  std::vector<int32_t> materialsInCells(const containedCells & contained_cells) const;
 
  void subdomainsToMaterials();
 
  std::set<SubdomainID> coupledSubdomains() const;
 
  template <typename T>
  void gatherCellSum(std::vector<T> & local, std::map<cellInfo, T> & global) const;
 
  template <typename T>
  void gatherCellVector(std::vector<T> & local, std::vector<unsigned int> & n_local, std::map<cellInfo, std::vector<T>> & global);
 
  coupling::CouplingFields elemFeedback(const Elem * elem) const;
 
  void getTallyTriggerParameters(const InputParameters & params);
 
  void readBlockParameters(const std::string name, std::unordered_set<SubdomainID> & blocks);
 
  void cacheContainedCells();
 
  void setContainedCells(const cellInfo & cell_info,
                         const Point & hint,
                         std::map<cellInfo, containedCells> & map);
 
  void checkContainedCellsStructure(const cellInfo & cell_info,
                                    containedCells & reference,
                                    containedCells & compare) const;
 
  void setMinimumVolumeQRules(Order & volume_order, const std::string & type);
 
  std::string printNewline() const
  {
    if (_verbose)
      return "\n";
    else
      return "";
  }
  std::string printNewline() const {…}
 
  void storeElementPhase();
 
  void relaxAndNormalizeTally(unsigned int global_score,
                              unsigned int local_score,
                              std::shared_ptr<TallyBase> local_tally);
 
  void getCellMappedPhase();
 
  void checkCellMappedPhase();
 
  void getCellMappedSubdomains();
 
  void computeCellMappedVolumes();
 
  void initializeElementToCellMapping();
 
  void mapElemsToCells();
 
  void validateLocalTallies();
 
  void initializeTallies();
 
  void resetTallies();
 
  void getMaterialFills();
 
  void getPointInCell();
 
  std::map<cellInfo, Real> computeVolumeWeightedCellInput(
      const std::map<SubdomainID, std::pair<unsigned int, std::string>> & var_num,
      const std::vector<coupling::CouplingFields> * phase) const;
 
  void sendTemperatureToOpenMC() const;
 
  void sendDensityToOpenMC() const;
 
  Real tallyNormalization(unsigned int global_score) const;
 
  void checkTallySum(const unsigned int & score) const;
 
  void latticeOuterCheck(const Point & c, int level) const;
 
  void latticeOuterError(const Point & c, int level) const;
 
  bool findCell(const Point & point);
 
  void compareContainedCells(std::map<cellInfo, containedCells> & reference,
                             std::map<cellInfo, containedCells> & compare) const;
 
  NumericVector<Number> & _serialized_solution;
 
  virtual std::vector<int32_t> getMappedTallyIDs() const override;
 
  const bool & _output_cell_mapping;
 
  const coupling::OpenMCInitialCondition _initial_condition;
 
  const relaxation::RelaxationEnum _relaxation;
 
  const trigger::TallyTriggerTypeEnum _k_trigger;
 
  unsigned int _cell_level;
 
  const bool & _export_properties;
 
  const bool _normalize_by_global;
 
  const bool _using_skinner;
 
  bool _need_to_reinit_coupling;
 
  const bool _check_tally_sum;
 
  const Real & _relaxation_factor;
 
  const bool _has_identical_cell_fills;
 
  const bool & _check_identical_cell_fills;
 
  const bool & _assume_separate_tallies;
 
  bool _map_density_by_cell;
 
  const bool _specified_density_feedback;
 
  const bool _specified_temperature_feedback;
 
  bool _has_cell_tallies = false;
 
  bool _needs_to_map_cells;
 
  const bool _needs_global_tally;
 
  std::map<std::string, std::shared_ptr<FilterBase>> _filters;
 
  std::vector<std::shared_ptr<TallyBase>> _local_tallies;
 
  std::vector<std::string> _all_tally_scores;
 
  std::vector<std::map<std::string, int>> _local_tally_score_map;
 
  std::vector<std::vector<unsigned int>> _tally_var_ids;
 
  std::vector<std::vector<std::vector<unsigned int>>> _tally_ext_var_ids;
 
  std::vector<SubdomainID> _density_blocks;
 
  std::vector<SubdomainID> _temp_blocks;
 
  std::unordered_set<SubdomainID> _identical_cell_fill_blocks;
 
  std::vector<cellInfo> _elem_to_cell{};
 
  std::map<cellInfo, coupling::CouplingFields> _cell_phase;
 
  int _n_moose_density_elems;
 
  int _n_moose_temp_elems;
 
  int _n_moose_temp_density_elems;
 
  int _n_moose_none_elems;
 
  int _n_mapped_temp_elems;
 
  int _n_mapped_density_elems;
 
  int _n_mapped_temp_density_elems;
 
  int _n_mapped_none_elems;
 
  Real _uncoupled_volume;
 
  bool _material_cells_only{true};
 
  std::map<cellInfo, std::vector<unsigned int>> _cell_to_elem;
 
  std::map<cellInfo, std::vector<unsigned int>> _local_cell_to_elem;
 
  std::map<cellInfo, std::unordered_set<SubdomainID>> _cell_to_elem_subdomain;
 
  std::map<SubdomainID, std::set<int32_t>> _subdomain_to_material;
 
  std::map<cellInfo, Point> _cell_to_point;
 
  std::map<cellInfo, Real> _cell_to_elem_volume;
 
  std::map<cellInfo, Real> _cell_volume;
 
  std::map<cellInfo, int32_t> _cell_to_material;
 
  std::map<cellInfo, containedCells> _cell_to_contained_material_cells;
 
  std::map<cellInfo, int32_t> _cell_to_n_contained;
 
  std::vector<openmc::Tally *> _global_tallies;
 
  std::vector<std::vector<std::string>> _global_tally_scores;
 
  std::vector<openmc::TallyEstimator> _global_tally_estimators;
 
  std::vector<Real> _global_sum_tally;
 
  std::vector<Real> _local_sum_tally;
 
  std::vector<Real> _local_mean_tally;
 
  static bool _first_transfer;
 
  static bool _printed_initial;
 
  static bool _printed_triso_warning;
 
  openmc::Particle _particle;
 
  unsigned int _n_particles_1;
 
  std::map<std::string, std::vector<SubdomainName>> _temp_vars_to_blocks;
 
  std::map<std::string, std::vector<SubdomainName>> _density_vars_to_blocks;
 
  OpenMCVolumeCalculation * _volume_calc;
 
  const SymmetryPointGenerator * _symmetry;
 
  std::map<cellInfo, int> _n_temp;
 
  std::map<cellInfo, int> _n_rho;
 
  std::map<cellInfo, int> _n_temp_rho;
 
  std::map<cellInfo, int> _n_none;
 
  unsigned int _global_tally_index = 0;
 
  unsigned int _source_rate_index;
 
#ifdef ENABLE_DAGMC
  MoabSkinner * _skinner = nullptr;
 
  std::shared_ptr<moab::DagMC> _dagmc = nullptr;
#endif
 
  long unsigned int _n_openmc_cells;
 
  int32_t _dagmc_universe_id;
 
  bool _dagmc_root_universe = true;
 
  int32_t _cell_using_dagmc_universe_id;
 
  const int32_t _initial_num_openmc_surfaces;
 
  static constexpr Real EV_TO_JOULE = 1.6022e-19;
 
  static constexpr Real ZERO_TALLY_THRESHOLD = 1e-12;
 
private:
  void dufekGudowskiParticleUpdate();
 
  std::vector<int32_t> _flattened_ids;
 
  std::vector<int32_t> _flattened_instances;
 
  containedCells _instance_offsets;
 
  std::map<cellInfo, int32_t> _n_offset;
 
  cellInfo _first_identical_cell;
 
  std::vector<int32_t> _first_identical_cell_materials;
 
  std::map<SubdomainID, std::pair<unsigned int, std::string>> _subdomain_to_temp_vars;
 
  std::map<SubdomainID, std::pair<unsigned int, std::string>> _subdomain_to_density_vars;
};
class OpenMCCellAverageProblem : public OpenMCProblemBase {…};