www.mooseframework.org
LineMaterialSamplerBase.h
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #ifndef LINEMATERIALSAMPLERBASE_H
16 #define LINEMATERIALSAMPLERBASE_H
17 
18 // MOOSE includes
20 #include "SamplerBase.h"
21 #include "BlockRestrictable.h"
22 #include "LineSegment.h"
23 #include "RayTracing.h" // Moose::elementsIntersectedByLine()
24 #include "Assembly.h" // Assembly::qRule()
25 #include "MooseMesh.h" // MooseMesh::getMesh()
26 #include "SwapBackSentinel.h"
27 #include "FEProblem.h"
28 
29 #include "libmesh/quadrature.h" // _qrule->n_points()
30 
31 // Forward Declarations
32 class MooseMesh;
33 template <typename T>
35 
36 template <>
37 InputParameters validParams<LineMaterialSamplerBase<Real>>();
38 
47 template <typename T>
49  public SamplerBase,
50  public BlockRestrictable
51 {
52 public:
59 
64  virtual void initialize() override;
65 
70  virtual void execute() override;
71 
76  virtual void finalize() override;
77 
84  virtual Real getScalarFromProperty(const T & property, const Point & curr_point) = 0;
85 
86 protected:
88  Point _start;
89 
91  Point _end;
92 
94  std::vector<const MaterialProperty<T> *> _material_properties;
95 
98 
100  QBase *& _qrule;
101 
104 };
105 
106 template <typename T>
108  : GeneralVectorPostprocessor(parameters),
109  SamplerBase(parameters, this, _communicator),
110  BlockRestrictable(this),
111  _start(getParam<Point>("start")),
112  _end(getParam<Point>("end")),
113  _mesh(_subproblem.mesh()),
114  _qrule(_subproblem.assembly(_tid).qRule()),
115  _q_point(_subproblem.assembly(_tid).qPoints())
116 {
117  std::vector<std::string> material_property_names = getParam<std::vector<std::string>>("property");
118  for (unsigned int i = 0; i < material_property_names.size(); ++i)
119  {
120  if (!hasMaterialProperty<T>(material_property_names[i]))
121  mooseError("In LineMaterialSamplerBase material property: " + material_property_names[i] +
122  " does not exist.");
123  _material_properties.push_back(&getMaterialProperty<T>(material_property_names[i]));
124  }
125 
126  SamplerBase::setupVariables(material_property_names);
127 }
128 
129 template <typename T>
130 void
132 {
134 }
135 
136 template <typename T>
137 void
139 {
140  std::vector<Elem *> intersected_elems;
141  std::vector<LineSegment> segments;
142 
143  std::unique_ptr<PointLocatorBase> pl = _mesh.getPointLocator();
144  Moose::elementsIntersectedByLine(_start, _end, _mesh, *pl, intersected_elems, segments);
145 
146  const RealVectorValue line_vec = _end - _start;
147  const Real line_length(line_vec.norm());
148  const RealVectorValue line_unit_vec = line_vec / line_length;
149  std::vector<Real> values(_material_properties.size());
150 
151  std::set<unsigned int> needed_mat_props;
152  const std::set<unsigned int> & mp_deps = getMatPropDependencies();
153  needed_mat_props.insert(mp_deps.begin(), mp_deps.end());
154  _fe_problem.setActiveMaterialProperties(needed_mat_props, _tid);
155 
156  for (const auto & elem : intersected_elems)
157  {
158  if (elem->processor_id() != processor_id())
159  continue;
160 
161  if (!hasBlocks(elem->subdomain_id()))
162  continue;
163 
164  _subproblem.prepare(elem, _tid);
165  _subproblem.reinitElem(elem, _tid);
166 
167  // Set up Sentinel class so that, even if reinitMaterials() throws, we
168  // still remember to swap back during stack unwinding.
170  _fe_problem.reinitMaterials(elem->subdomain_id(), _tid);
171 
172  for (unsigned int qp = 0; qp < _qrule->n_points(); ++qp)
173  {
174  const RealVectorValue qp_pos(_q_point[qp]);
175 
176  const RealVectorValue start_to_qp(qp_pos - _start);
177  const Real qp_proj_dist_along_line = start_to_qp * line_unit_vec;
178 
179  if (qp_proj_dist_along_line < 0 || qp_proj_dist_along_line > line_length)
180  continue;
181 
182  for (unsigned int j = 0; j < _material_properties.size(); ++j)
183  values[j] = getScalarFromProperty((*_material_properties[j])[qp], _q_point[qp]);
184 
185  addSample(_q_point[qp], qp_proj_dist_along_line, values);
186  }
187  }
189 }
190 
191 template <typename T>
192 void
194 {
196 }
197 
198 #endif
Base class for VectorPostprocessors that need to do "sampling" of values in the domain.
Definition: SamplerBase.h:46
virtual void initialize()
Initialize the datastructures.
Definition: SamplerBase.C:75
VectorValue< Real > RealVectorValue
Definition: Assembly.h:40
QBase *& _qrule
The quadrature rule.
This class is here to combine the VectorPostprocessor interface and the base class VectorPostprocesso...
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()
THREAD_ID _tid
Thread ID of this postprocessor.
Definition: UserObject.h:152
Point _start
The beginning of the line.
virtual void reinitMaterials(SubdomainID blk_id, THREAD_ID tid, bool swap_stateful=true)
virtual void setActiveMaterialProperties(const std::set< unsigned int > &mat_prop_ids, THREAD_ID tid) override
Record and set the material properties required by the current computing thread.
Point _end
The end of the line.
SubProblem & _subproblem
Reference to the Subproblem for this user object.
Definition: UserObject.h:146
virtual void finalize() override
Finalize Calls through to base class&#39;s finalize()
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.
virtual void clearActiveMaterialProperties(THREAD_ID tid) override
Clear the active material properties.
virtual void prepare(const Elem *elem, THREAD_ID tid)=0
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
bool hasBlocks(const SubdomainName &name) const
Test if the supplied block name is valid for this object.
virtual void swapBackMaterials(THREAD_ID tid)
const T & getParam(const std::string &name) const
Retrieve a parameter for the object.
Definition: MooseObject.h:122
virtual std::unique_ptr< PointLocatorBase > getPointLocator() const
Proxy function to get a (sub)PointLocator from either the underlying libmesh mesh (default)...
Definition: MooseMesh.C:2601
virtual void reinitElem(const Elem *elem, THREAD_ID tid)=0
virtual Real getScalarFromProperty(const T &property, const Point &curr_point)=0
Reduce the material property to a scalar for output.
virtual void finalize()
Finalize the values.
Definition: SamplerBase.C:87
virtual void addSample(const Point &p, const Real &id, const std::vector< Real > &values)
Call this with the value of every variable at each point you want to sample at.
Definition: SamplerBase.C:61
const MooseArray< Point > & _q_point
The quadrature points.
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:197
FEProblemBase & _fe_problem
Reference to the FEProblemBase for this user object.
Definition: UserObject.h:149
const std::set< unsigned int > & getMatPropDependencies() const
Retrieve the set of material properties that this object depends on.
An interface that restricts an object to subdomains via the &#39;blocks&#39; input parameter.
const InputParameters & parameters() const
Get the parameters of the object.
Definition: MooseObject.h:53
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
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}.