23 #include "libmesh/quadrature.h" 61 virtual void execute()
override;
107 "property",
"Name of the material property to be output along a line");
112 template <
typename T>
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())
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)
129 if (!hasMaterialProperty<T>(material_property_names[i]))
130 mooseError(
"In LineMaterialSamplerBase material property: " + material_property_names[i] +
138 template <
typename T>
145 template <
typename T>
149 std::vector<Elem *> intersected_elems;
150 std::vector<LineSegment> segments;
152 std::unique_ptr<PointLocatorBase> pl = _mesh.getPointLocator();
156 const Real line_length(line_vec.
norm());
158 std::vector<Real> values(_material_properties.size());
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);
165 for (
const auto & elem : intersected_elems)
167 if (elem->processor_id() != processor_id())
170 if (!hasBlocks(elem->subdomain_id()))
173 _subproblem.setCurrentSubdomainID(elem, _tid);
174 _subproblem.prepare(elem, _tid);
175 _subproblem.reinitElem(elem, _tid);
180 _fe_problem.reinitMaterials(elem->subdomain_id(), _tid);
182 for (
unsigned int qp = 0; qp < _qrule->n_points(); ++qp)
187 const Real qp_proj_dist_along_line = start_to_qp * line_unit_vec;
189 if (qp_proj_dist_along_line < 0 || qp_proj_dist_along_line > line_length)
192 for (
unsigned int j = 0; j < _material_properties.size(); ++j)
193 values[j] = getScalarFromProperty((*_material_properties[j])[qp], _q_point[qp]);
195 addSample(_q_point[qp], qp_proj_dist_along_line, values);
198 _fe_problem.clearActiveMaterialProperties(_tid);
201 template <
typename T>
Base class for VectorPostprocessors that need to do "sampling" of values in the domain.
virtual void initialize()
Initialize the datastructures.
auto norm() const -> decltype(std::norm(Real()))
static InputParameters validParams()
This class is here to combine the VectorPostprocessor interface and the base class VectorPostprocesso...
LineMaterialSamplerBase(const InputParameters ¶meters)
Class constructor Sets up variables for output based on the properties to be output.
std::vector< const MaterialProperty< T > * > _material_properties
The material properties to be output.
virtual void initialize() override
Initialize Calls through to base class's initialize()
Point _start
The beginning of the line.
Point _end
The end of the line.
static InputParameters validParams()
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
virtual void finalize() override
Finalize Calls through to base class'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...
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...
virtual Real getScalarFromProperty(const T &property, const Point &curr_point)=0
Reduce the material property to a scalar for output.
static InputParameters validParams()
virtual void finalize()
Finalize the values.
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 'blocks' 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.
The "SwapBackSentinel" class's destructor guarantees that FEProblemBase::swapBackMaterials{Face,Neighbor}() is called even when an exception is thrown from FEProblemBase::reinitMaterials{Face,Neighbor}.