Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
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 : #include "ElementMaterialSampler.h"
11 : #include "Material.h"
12 : #include "IndirectSort.h"
13 : #include "MooseMesh.h"
14 :
15 : #include "libmesh/quadrature.h"
16 :
17 : #include <numeric>
18 :
19 : registerMooseObject("MooseApp", ElementMaterialSampler);
20 : registerMooseObjectRenamed("MooseApp",
21 : MaterialVectorPostprocessor,
22 : "06/30/2025 24:00",
23 : ElementMaterialSampler);
24 :
25 : InputParameters
26 28701 : ElementMaterialSampler::validParams()
27 : {
28 28701 : InputParameters params = ElementVectorPostprocessor::validParams();
29 28701 : params.addClassDescription("Records all Real-valued material properties of a material object, "
30 : "or Real-valued material properties of the supplied property names "
31 : "on quadrature points on elements at the indicated execution points.");
32 28701 : params.addParam<MaterialName>("material", "Material for which all properties will be recorded.");
33 28701 : params.addParam<std::vector<MaterialPropertyName>>(
34 : "property", "Material property names that will be recorded.");
35 28701 : params.addParam<std::vector<dof_id_type>>(
36 : "elem_ids",
37 : "Subset of element IDs to print data for. If omitted, all elements will be printed.");
38 28701 : return params;
39 0 : }
40 :
41 95 : ElementMaterialSampler::ElementMaterialSampler(const InputParameters & parameters)
42 : : ElementVectorPostprocessor(parameters),
43 95 : _elem_ids(declareVector("elem_id")),
44 95 : _qp_ids(declareVector("qp_id")),
45 95 : _x_coords(declareVector("x")),
46 95 : _y_coords(declareVector("y")),
47 190 : _z_coords(declareVector("z"))
48 : {
49 : // Check either "material" or "property" was set but not both
50 95 : if (parameters.isParamValid("material") && parameters.isParamValid("property"))
51 0 : mooseError("Setting both 'material' and 'property' is not allowed. Use one or the other.");
52 95 : if (!parameters.isParamValid("material") && !parameters.isParamValid("property"))
53 0 : mooseError("Either 'material' and 'property' needs to be set.");
54 :
55 : // List of property names to collect
56 95 : std::vector<MaterialName> prop_names;
57 :
58 : // Get list of elements from user
59 95 : if (parameters.isParamValid("elem_ids"))
60 : {
61 34 : const auto & ids = getParam<std::vector<dof_id_type>>("elem_ids");
62 34 : _elem_filter = std::set<dof_id_type>(ids.begin(), ids.end());
63 : }
64 :
65 : // If Material is used, get all properties.
66 95 : if (parameters.isParamValid("material"))
67 : {
68 51 : auto & mat = getMaterialByName(getParam<MaterialName>("material"), true);
69 47 : if (mat.isBoundaryMaterial())
70 4 : mooseError(name(), ": boundary materials (i.e. ", mat.name(), ") cannot be used");
71 :
72 : // Get property names from the Material
73 43 : auto & props = mat.getSuppliedItems(); // returns std::set
74 43 : prop_names = std::vector<MaterialName>(props.begin(), props.end());
75 :
76 : // Check requested materials are available
77 43 : if (_elem_filter)
78 : {
79 82 : for (const auto & id : *_elem_filter)
80 : {
81 69 : auto el = _mesh.getMesh().query_elem_ptr(id);
82 :
83 : // We'd better have found the requested element on *some*
84 : // processor.
85 69 : bool found_elem = (el != nullptr);
86 69 : this->comm().max(found_elem);
87 :
88 : // We might not have el on this processor in a distributed mesh,
89 : // but it should be somewhere and it ought to have a material
90 : // defined for its subdomain
91 69 : if (!found_elem || (el && !mat.hasBlocks(el->subdomain_id())))
92 4 : mooseError(name(), ": material ", mat.name(), " is not defined on element ", id);
93 : }
94 : }
95 : }
96 : else
97 : {
98 :
99 : // Properties supplied by user
100 44 : auto & props = getParam<std::vector<MaterialPropertyName>>("property");
101 44 : prop_names = std::vector<MaterialName>(props.begin(), props.end());
102 : }
103 :
104 : // Check properties are valid and store references
105 348 : for (auto & prop : prop_names)
106 : {
107 265 : if (hasMaterialProperty<Real>(prop))
108 265 : _prop_refs.push_back(&getMaterialProperty<Real>(prop));
109 0 : else if (hasMaterialProperty<unsigned int>(prop))
110 0 : _prop_refs.push_back(&getMaterialProperty<unsigned int>(prop));
111 0 : else if (hasMaterialProperty<int>(prop))
112 0 : _prop_refs.push_back(&getMaterialProperty<int>(prop));
113 : else
114 : {
115 0 : mooseWarning("property " + prop +
116 : " is of unsupported type and skipped by ElementMaterialSampler");
117 0 : continue;
118 : }
119 265 : _prop_vecs.push_back(&declareVector(prop));
120 265 : _prop_names.push_back(prop);
121 : }
122 83 : }
123 :
124 : void
125 180 : ElementMaterialSampler::initialize()
126 : {
127 180 : if (!containsCompleteHistory())
128 : {
129 180 : _elem_ids.clear();
130 180 : _qp_ids.clear();
131 180 : _x_coords.clear();
132 180 : _y_coords.clear();
133 180 : _z_coords.clear();
134 792 : for (auto vec : _prop_vecs)
135 612 : vec->clear();
136 : }
137 180 : }
138 :
139 : void
140 2032 : ElementMaterialSampler::execute()
141 : {
142 : // skip execution if element not in filter, assuming filter was used
143 2032 : const auto elem_id = _current_elem->id();
144 2032 : if (_elem_filter && !_elem_filter->count(elem_id))
145 1568 : return;
146 :
147 464 : unsigned int nqp = _qrule->n_points();
148 2320 : for (unsigned int qp = 0; qp < nqp; qp++)
149 : {
150 1856 : _elem_ids.push_back(elem_id);
151 1856 : _qp_ids.push_back(qp);
152 1856 : _x_coords.push_back(_q_point[qp](0));
153 1856 : _y_coords.push_back(_q_point[qp](1));
154 1856 : _z_coords.push_back(_q_point[qp](2));
155 : }
156 :
157 2000 : for (unsigned int i = 0; i < _prop_names.size(); i++)
158 : {
159 1536 : auto prop_name = _prop_names[i];
160 1536 : auto prop = _prop_vecs[i];
161 1536 : std::vector<Real> vals;
162 1536 : if (hasMaterialProperty<Real>(prop_name))
163 : {
164 1536 : auto vals = dynamic_cast<const MaterialProperty<Real> *>(_prop_refs[i]);
165 7680 : for (unsigned int qp = 0; qp < nqp; qp++)
166 6144 : prop->push_back((*vals)[qp]);
167 : }
168 0 : else if (hasMaterialProperty<unsigned int>(prop_name))
169 : {
170 0 : auto vals = dynamic_cast<const MaterialProperty<unsigned int> *>(_prop_refs[i]);
171 0 : for (unsigned int qp = 0; qp < nqp; qp++)
172 0 : prop->push_back((*vals)[qp]);
173 : }
174 0 : else if (hasMaterialProperty<int>(prop_name))
175 : {
176 0 : auto vals = dynamic_cast<const MaterialProperty<int> *>(_prop_refs[i]);
177 0 : for (unsigned int qp = 0; qp < nqp; qp++)
178 0 : prop->push_back((*vals)[qp]);
179 : }
180 1536 : }
181 : }
182 :
183 : void
184 165 : ElementMaterialSampler::finalize()
185 : {
186 : // collect all processor data
187 165 : comm().gather(0, _elem_ids);
188 165 : comm().gather(0, _qp_ids);
189 165 : comm().gather(0, _x_coords);
190 165 : comm().gather(0, _y_coords);
191 165 : comm().gather(0, _z_coords);
192 726 : for (auto vec : _prop_vecs)
193 561 : comm().gather(0, *vec);
194 165 : sortVecs();
195 165 : }
196 :
197 : void
198 15 : ElementMaterialSampler::threadJoin(const UserObject & y)
199 : {
200 15 : const auto & vpp = static_cast<const ElementMaterialSampler &>(y);
201 15 : _elem_ids.insert(_elem_ids.end(), vpp._elem_ids.begin(), vpp._elem_ids.end());
202 15 : _qp_ids.insert(_qp_ids.end(), vpp._qp_ids.begin(), vpp._qp_ids.end());
203 15 : _x_coords.insert(_x_coords.end(), vpp._x_coords.begin(), vpp._x_coords.end());
204 15 : _y_coords.insert(_y_coords.end(), vpp._y_coords.begin(), vpp._y_coords.end());
205 15 : _z_coords.insert(_z_coords.end(), vpp._z_coords.begin(), vpp._z_coords.end());
206 :
207 66 : for (unsigned int i = 0; i < _prop_vecs.size(); i++)
208 : {
209 51 : auto & vec = *_prop_vecs[i];
210 51 : auto & othervec = *vpp._prop_vecs[i];
211 51 : vec.insert(vec.end(), othervec.begin(), othervec.end());
212 : }
213 15 : sortVecs();
214 15 : }
215 :
216 : void
217 180 : ElementMaterialSampler::sortVecs()
218 : {
219 180 : std::vector<size_t> ind;
220 180 : ind.resize(_elem_ids.size());
221 180 : std::iota(ind.begin(), ind.end(), 0);
222 180 : std::sort(ind.begin(),
223 : ind.end(),
224 6363 : [&](size_t a, size_t b) -> bool
225 : {
226 6363 : if (_elem_ids[a] == _elem_ids[b])
227 : {
228 2584 : return _qp_ids[a] < _qp_ids[b];
229 : }
230 3779 : return _elem_ids[a] < _elem_ids[b];
231 : });
232 :
233 180 : Moose::applyIndices(_elem_ids, ind);
234 180 : Moose::applyIndices(_qp_ids, ind);
235 180 : Moose::applyIndices(_x_coords, ind);
236 180 : Moose::applyIndices(_y_coords, ind);
237 180 : Moose::applyIndices(_z_coords, ind);
238 792 : for (auto vec : _prop_vecs)
239 612 : Moose::applyIndices(*vec, ind);
240 180 : }
|