www.mooseframework.org
LineMaterialSamplerBase.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
12 // MOOSE includes
14 #include "SamplerBase.h"
15 #include "BlockRestrictable.h"
16 #include "LineSegment.h"
17 #include "RayTracing.h" // Moose::elementsIntersectedByLine()
18 #include "Assembly.h" // Assembly::qRule()
19 #include "MooseMesh.h" // MooseMesh::getMesh()
20 #include "SwapBackSentinel.h"
21 #include "FEProblem.h"
22 
23 #include "libmesh/quadrature.h" // _qrule->n_points()
24 
25 // Forward Declarations
26 class MooseMesh;
27 
36 template <typename T>
38  public SamplerBase,
39  public BlockRestrictable
40 {
41 public:
43 
50 
55  virtual void initialize() override;
56 
61  virtual void execute() override;
62 
67  virtual void finalize() override;
68 
75  virtual Real getScalarFromProperty(const T & property, const Point & curr_point) = 0;
76 
77 protected:
79  Point _start;
80 
82  Point _end;
83 
85  std::vector<const MaterialProperty<T> *> _material_properties;
86 
89 
91  const QBase * const & _qrule;
92 
95 };
96 
97 template <typename T>
100 {
102  params += SamplerBase::validParams();
103  params += BlockRestrictable::validParams();
104  params.addRequiredParam<Point>("start", "The beginning of the line");
105  params.addRequiredParam<Point>("end", "The end of the line");
106  params.addRequiredParam<std::vector<std::string>>(
107  "property", "Name of the material property to be output along a line");
108 
109  return params;
110 }
111 
112 template <typename T>
114  : GeneralVectorPostprocessor(parameters),
115  SamplerBase(parameters, this, _communicator),
116  BlockRestrictable(this),
117  _start(getParam<Point>("start")),
118  _end(getParam<Point>("end")),
119  _mesh(_subproblem.mesh()),
120  _qrule(_subproblem.assembly(_tid, 0).qRule()),
121  _q_point(_subproblem.assembly(_tid, 0).qPoints())
122 {
123  // See https://github.com/idaholab/moose/issues/21865
124  _mesh.errorIfDistributedMesh("LineMaterialSamplerBase");
125 
126  std::vector<std::string> material_property_names = getParam<std::vector<std::string>>("property");
127  for (unsigned int i = 0; i < material_property_names.size(); ++i)
128  {
129  if (!hasMaterialProperty<T>(material_property_names[i]))
130  mooseError("In LineMaterialSamplerBase material property: " + material_property_names[i] +
131  " does not exist.");
132  _material_properties.push_back(&getMaterialProperty<T>(material_property_names[i]));
133  }
134 
135  SamplerBase::setupVariables(material_property_names);
136 }
137 
138 template <typename T>
139 void
141 {
143 }
144 
145 template <typename T>
146 void
148 {
149  std::vector<Elem *> intersected_elems;
150  std::vector<LineSegment> segments;
151 
152  std::unique_ptr<PointLocatorBase> pl = _mesh.getPointLocator();
153  Moose::elementsIntersectedByLine(_start, _end, _mesh, *pl, intersected_elems, segments);
154 
155  const RealVectorValue line_vec = _end - _start;
156  const Real line_length(line_vec.norm());
157  const RealVectorValue line_unit_vec = line_vec / line_length;
158  std::vector<Real> values(_material_properties.size());
159 
160  std::unordered_set<unsigned int> needed_mat_props;
161  const auto & mp_deps = getMatPropDependencies();
162  needed_mat_props.insert(mp_deps.begin(), mp_deps.end());
163  _fe_problem.setActiveMaterialProperties(needed_mat_props, _tid);
164 
165  for (const auto & elem : intersected_elems)
166  {
167  if (elem->processor_id() != processor_id())
168  continue;
169 
170  if (!hasBlocks(elem->subdomain_id()))
171  continue;
172 
173  _subproblem.setCurrentSubdomainID(elem, _tid);
174  _subproblem.prepare(elem, _tid);
175  _subproblem.reinitElem(elem, _tid);
176 
177  // Set up Sentinel class so that, even if reinitMaterials() throws, we
178  // still remember to swap back during stack unwinding.
179  SwapBackSentinel sentinel(_fe_problem, &FEProblem::swapBackMaterials, _tid);
180  _fe_problem.reinitMaterials(elem->subdomain_id(), _tid);
181 
182  for (unsigned int qp = 0; qp < _qrule->n_points(); ++qp)
183  {
184  const RealVectorValue qp_pos(_q_point[qp]);
185 
186  const RealVectorValue start_to_qp(qp_pos - _start);
187  const Real qp_proj_dist_along_line = start_to_qp * line_unit_vec;
188 
189  if (qp_proj_dist_along_line < 0 || qp_proj_dist_along_line > line_length)
190  continue;
191 
192  for (unsigned int j = 0; j < _material_properties.size(); ++j)
193  values[j] = getScalarFromProperty((*_material_properties[j])[qp], _q_point[qp]);
194 
195  addSample(_q_point[qp], qp_proj_dist_along_line, values);
196  }
197  }
198  _fe_problem.clearActiveMaterialProperties(_tid);
199 }
200 
201 template <typename T>
202 void
204 {
206 }
Base class for VectorPostprocessors that need to do "sampling" of values in the domain.
Definition: SamplerBase.h:36
virtual void initialize()
Initialize the datastructures.
Definition: SamplerBase.C:75
auto norm() const -> decltype(std::norm(Real()))
static InputParameters validParams()
This class is here to combine the VectorPostprocessor interface and the base class VectorPostprocesso...
MeshBase & mesh
LineMaterialSamplerBase(const InputParameters &parameters)
Class constructor Sets up variables for output based on the properties to be output.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
std::vector< const MaterialProperty< T > * > _material_properties
The material properties to be output.
virtual void initialize() override
Initialize Calls through to base class&#39;s initialize()
Point _start
The beginning of the line.
Point _end
The end of the line.
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
static InputParameters validParams()
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
Definition: MooseMesh.C:3367
virtual void finalize() override
Finalize Calls through to base class&#39;s finalize()
const QBase *const & _qrule
The quadrature rule.
void setupVariables(const std::vector< std::string > &variable_names)
You MUST call this in the constructor of the child class and pass down the name of the variables...
Definition: SamplerBase.C:51
virtual void execute() override
Finds all elements along the user-defined line, loops through them, and samples their material proper...
This is a base class for sampling material properties for the integration points in all elements that...
MooseMesh & _mesh
The mesh.
static InputParameters validParams()
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:88
virtual Real getScalarFromProperty(const T &property, const Point &curr_point)=0
Reduce the material property to a scalar for output.
static InputParameters validParams()
Definition: SamplerBase.C:22
virtual void finalize()
Finalize the values.
Definition: SamplerBase.C:91
virtual void swapBackMaterials(const THREAD_ID tid)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const MooseArray< Point > & _q_point
The quadrature points.
An interface that restricts an object to subdomains via the &#39;blocks&#39; input parameter.
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
const InputParameters & parameters() const
Get the parameters of the object.
void elementsIntersectedByLine(const Point &p0, const Point &p1, const MeshBase &mesh, const PointLocatorBase &point_locator, std::vector< Elem *> &intersected_elems, std::vector< LineSegment > &segments)
Find all of the elements intersected by a line.
Definition: RayTracing.C:192
The "SwapBackSentinel" class&#39;s destructor guarantees that FEProblemBase::swapBackMaterials{Face,Neighbor}() is called even when an exception is thrown from FEProblemBase::reinitMaterials{Face,Neighbor}.