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 "SamplerBase.h"
11 : #include "IndirectSort.h"
12 : #include "InputParameters.h"
13 : #include "MooseEnum.h"
14 : #include "MooseError.h"
15 : #include "VectorPostprocessor.h"
16 : #include "MooseVariableFieldBase.h"
17 : #include "FEProblemBase.h"
18 : #include "MooseApp.h"
19 : #include "TransientBase.h"
20 :
21 : #include "libmesh/point.h"
22 :
23 : InputParameters
24 131893 : SamplerBase::validParams()
25 : {
26 131893 : InputParameters params = emptyInputParameters();
27 :
28 527572 : MooseEnum sort_options("x y z id");
29 527572 : params.addRequiredParam<MooseEnum>("sort_by", sort_options, "What to sort the samples by");
30 :
31 : // The value from this VPP is naturally already on every processor
32 : // TODO: Make this not the case! See #11415
33 131893 : params.set<bool>("_auto_broadcast") = false;
34 :
35 263786 : return params;
36 131893 : }
37 :
38 1778 : SamplerBase::SamplerBase(const InputParameters & parameters,
39 : VectorPostprocessor * vpp,
40 1778 : const libMesh::Parallel::Communicator & comm)
41 1778 : : _sampler_params(parameters),
42 1778 : _vpp(vpp),
43 1778 : _comm(comm),
44 1778 : _sampler_transient(dynamic_cast<TransientBase *>(
45 7112 : parameters.getCheckedPointerParam<FEProblemBase *>("_fe_problem_base")
46 1778 : ->getMooseApp()
47 1778 : .getExecutioner())),
48 1778 : _sort_by(parameters.get<MooseEnum>("sort_by")),
49 3556 : _x(vpp->declareVector("x")),
50 3556 : _y(vpp->declareVector("y")),
51 3556 : _z(vpp->declareVector("z")),
52 8890 : _id(vpp->declareVector("id"))
53 : {
54 1778 : }
55 :
56 : void
57 1754 : SamplerBase::setupVariables(const std::vector<std::string> & variable_names)
58 : {
59 1754 : _variable_names = variable_names;
60 1754 : _values.reserve(variable_names.size());
61 :
62 4175 : for (const auto & variable_name : variable_names)
63 2421 : _values.push_back(&_vpp->declareVector(variable_name));
64 1754 : }
65 :
66 : void
67 233431 : SamplerBase::addSample(const Point & p, const Real & id, const std::vector<Real> & values)
68 : {
69 233431 : _x.push_back(p(0));
70 233431 : _y.push_back(p(1));
71 233431 : _z.push_back(p(2));
72 :
73 233431 : _id.push_back(id);
74 :
75 : mooseAssert(values.size() == _variable_names.size(), "Mismatch of variable names to vector size");
76 473243 : for (MooseIndex(values) i = 0; i < values.size(); ++i)
77 239812 : _values[i]->emplace_back(values[i]);
78 :
79 233431 : _curr_num_samples++;
80 233431 : }
81 :
82 : void
83 9987 : SamplerBase::initialize()
84 : {
85 : // Don't reset the vectors if we want to retain history
86 9987 : if (_vpp->containsCompleteHistory() && _comm.rank() == 0)
87 : {
88 : // If we are repeating the timestep due to an aborted solve, we want to throw away the last
89 : // values
90 228 : if (_sampler_transient && !_sampler_transient->lastSolveConverged())
91 : {
92 : // Convenient to allocate a single variable for all the vpp values
93 6 : std::vector<VectorPostprocessorValue *> vec_ptrs = {{&_x, &_y, &_z, &_id}};
94 3 : vec_ptrs.insert(vec_ptrs.end(), _values.begin(), _values.end());
95 : // Erase the elements from the last execution
96 18 : for (auto vec_ptr : vec_ptrs)
97 : {
98 : // Vector may have already been restored, if so, skip the erasure
99 15 : if (_curr_total_samples > vec_ptr->size())
100 : {
101 : mooseAssert(vec_ptr->size() == (_curr_total_samples - _curr_num_samples),
102 : "Number of samples is not what is expected.");
103 1 : continue;
104 : }
105 28 : for (auto ind : _curr_indices)
106 : {
107 : mooseAssert(ind < vec_ptr->size(), "Trying to remove a sample that doesn't exist.");
108 14 : vec_ptr->erase(vec_ptr->begin() + ind);
109 : }
110 : }
111 3 : }
112 228 : _curr_indices.clear();
113 : }
114 : else
115 : {
116 9759 : _x.clear();
117 9759 : _y.clear();
118 9759 : _z.clear();
119 9759 : _id.clear();
120 :
121 9759 : std::for_each(
122 10641 : _values.begin(), _values.end(), [](VectorPostprocessorValue * vec) { vec->clear(); });
123 : }
124 :
125 9987 : _curr_num_samples = 0;
126 9987 : }
127 :
128 : void
129 1997 : SamplerBase::checkForStandardFieldVariableType(const MooseVariableFieldBase * const var_ptr,
130 : const std::string & var_param_name) const
131 : {
132 : // A pointer to a MooseVariableFieldBase should never be SCALAR
133 : mooseAssert(var_ptr->feType().family != SCALAR,
134 : "Scalar variable '" + var_ptr->name() + "' cannot be sampled.");
135 : mooseAssert(dynamic_cast<const MooseObject *>(_vpp), "Should have succeeded");
136 1997 : if (var_ptr->isVector())
137 16 : dynamic_cast<const MooseObject *>(_vpp)->paramError(
138 : var_param_name,
139 : "The variable '",
140 8 : var_ptr->name(),
141 : "' is a vector variable. Sampling those is not currently supported in the "
142 : "framework. It may be supported using a dedicated object in your application. Use "
143 : "'VectorVariableComponentAux' auxkernel to copy those values into a regular field "
144 : "variable");
145 1989 : if (var_ptr->isArray())
146 16 : dynamic_cast<const MooseObject *>(_vpp)->paramError(
147 : var_param_name,
148 : "The variable '",
149 8 : var_ptr->name(),
150 : "' is an array variable. Sampling those is not currently supported in the "
151 : "framework. It may be supported using a dedicated object in your application. Use "
152 : "'ArrayVariableComponent' auxkernel to copy those values into a regular field variable");
153 1981 : }
154 :
155 : void
156 9567 : SamplerBase::finalize()
157 : {
158 : /**
159 : * We have several vectors that all need to be processed in the same way.
160 : * Rather than enumerate them all, let's just create a vector of pointers
161 : * and work on them that way.
162 : */
163 9567 : constexpr auto NUM_ID_VECTORS = 4;
164 :
165 9567 : std::vector<VectorPostprocessorValue *> vec_ptrs;
166 9567 : vec_ptrs.reserve(_values.size() + NUM_ID_VECTORS);
167 : // Initialize the pointer vector with the position and ID vectors
168 19134 : vec_ptrs = {{&_x, &_y, &_z, &_id}};
169 : // Now extend the vector by all the remaining values vector before processing
170 9567 : vec_ptrs.insert(vec_ptrs.end(), _values.begin(), _values.end());
171 :
172 : // Gather up each of the partial vectors
173 58275 : for (auto vec_ptr : vec_ptrs)
174 48708 : _comm.allgather(*vec_ptr, /* identical buffer lengths = */ false);
175 :
176 : // Now create an index vector by using an indirect sort
177 9567 : std::vector<std::size_t> sorted_indices;
178 9567 : Moose::indirectSort(vec_ptrs[_sort_by]->begin(), vec_ptrs[_sort_by]->end(), sorted_indices);
179 :
180 : /**
181 : * We now have one sorted vector. The remaining vectors need to be sorted according to that
182 : * vector.
183 : * We'll need a temporary vector to hold values as we map them according to the sorted indices.
184 : * After that, we'll swap the vector contents with the sorted vector to get the values
185 : * back into the original vector.
186 : */
187 : // This vector is used as temp storage to sort each of the remaining vectors according to the
188 : // first
189 9567 : auto vector_length = sorted_indices.size();
190 9567 : VectorPostprocessorValue tmp_vector(vector_length);
191 :
192 : #ifndef NDEBUG
193 : for (const auto vec_ptr : vec_ptrs)
194 : if (vec_ptr->size() != vector_length)
195 : mooseError("Vector length mismatch");
196 : #endif
197 :
198 : // Sort each of the vectors using the same sorted indices
199 58275 : for (auto & vec_ptr : vec_ptrs)
200 : {
201 1643639 : for (MooseIndex(sorted_indices) j = 0; j < sorted_indices.size(); ++j)
202 1594931 : tmp_vector[j] = (*vec_ptr)[sorted_indices[j]];
203 :
204 : // Swap vector storage with sorted vector
205 48708 : vec_ptr->swap(tmp_vector);
206 : }
207 :
208 : // Gather the indices of samples from the last execution
209 : // Used to determine which parts of the vector need to be erased if a solve fails
210 9567 : if (_vpp->containsCompleteHistory())
211 : {
212 303 : _comm.sum(_curr_num_samples);
213 303 : if (_comm.rank() == 0)
214 : {
215 228 : _curr_indices.insert(sorted_indices.end() - _curr_num_samples, sorted_indices.end());
216 228 : _curr_total_samples = vec_ptrs[0]->size();
217 : }
218 : }
219 9567 : }
220 :
221 : void
222 416 : SamplerBase::threadJoin(const SamplerBase & y)
223 : {
224 416 : _x.insert(_x.end(), y._x.begin(), y._x.end());
225 416 : _y.insert(_y.end(), y._y.begin(), y._y.end());
226 416 : _z.insert(_z.end(), y._z.begin(), y._z.end());
227 :
228 416 : _id.insert(_id.end(), y._id.begin(), y._id.end());
229 :
230 837 : for (MooseIndex(_variable_names) i = 0; i < _variable_names.size(); i++)
231 421 : _values[i]->insert(_values[i]->end(), y._values[i]->begin(), y._values[i]->end());
232 :
233 416 : _curr_num_samples += y._curr_num_samples;
234 416 : }
|