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 "DerivativeParsedMaterialHelper.h"
11 : #include "Conversion.h"
12 :
13 : #include <deque>
14 :
15 : #include "libmesh/quadrature.h"
16 :
17 : template <bool is_ad>
18 : InputParameters
19 9492 : DerivativeParsedMaterialHelperTempl<is_ad>::validParams()
20 : {
21 9492 : InputParameters params = ParsedMaterialHelper<is_ad>::validParams();
22 18984 : params.addClassDescription("Parsed Function Material with automatic derivatives.");
23 37968 : params.addParam<unsigned int>("derivative_order", 3, "Maximum order of derivatives taken");
24 28476 : params.addParam<std::vector<SymbolName>>(
25 : "additional_derivative_symbols",
26 : {},
27 : "A list of additional (non-variable) symbols (such as material property or postprocessor "
28 : "names) to take derivatives w.r.t.");
29 9492 : return params;
30 0 : }
31 :
32 : template <bool is_ad>
33 2589 : DerivativeParsedMaterialHelperTempl<is_ad>::DerivativeParsedMaterialHelperTempl(
34 : const InputParameters & parameters,
35 : const VariableNameMappingMode map_mode,
36 : const std::optional<std::string> & function_param_name /* = {} */)
37 : : ParsedMaterialHelper<is_ad>(parameters, map_mode, function_param_name),
38 2589 : _derivative_order(this->template getParam<unsigned int>("derivative_order")),
39 5178 : _dmatvar_base("matpropautoderiv"),
40 7767 : _dmatvar_index(0)
41 : {
42 2589 : }
43 :
44 : template <bool is_ad>
45 : void
46 2589 : DerivativeParsedMaterialHelperTempl<is_ad>::functionsPostParse()
47 : {
48 : // set up the list of all FParser symbols to derive w.r.t.
49 5178 : auto additional_derivative_symbols =
50 : this->template getParam<std::vector<SymbolName>>("additional_derivative_symbols");
51 :
52 : // build a list of indices into _symbol_table for all symbols we want to derive w.r.t.
53 : // we start with all MOOSE variables, which we _always_ take the derivatives w.r.t.
54 2589 : _derivative_symbol_table.reserve(_nargs + additional_derivative_symbols.size());
55 4965 : for (unsigned i = 0; i < _nargs; ++i)
56 2376 : _derivative_symbol_table.push_back(i);
57 :
58 : // then we add the additional derivative symbols
59 2667 : for (auto & ads : additional_derivative_symbols)
60 : {
61 : // find the index of the symbol if it exists
62 87 : auto it = std::find(_symbol_names.begin(), _symbol_names.end(), ads);
63 87 : if (it == _symbol_names.end())
64 6 : this->paramError("additional_derivative_symbols", "Invalid symbol name '", ads, "'.");
65 84 : auto idx = std::distance(_symbol_names.begin(), it);
66 :
67 : // make sure it is not a MOOSE variable (we already added those) - their index is < nargs
68 84 : if (idx < _nargs)
69 6 : this->paramError(
70 : "additional_derivative_symbols",
71 : "Symbol name '",
72 : ads,
73 : "' is a coupled variable. Derivatives w.r.t. coupled variables are always taken "
74 : "and they must not be specified again.");
75 :
76 : // check that the user did not specify any duplicates in the additional_derivative_symbols
77 : // parameter
78 81 : if (std::count(
79 81 : additional_derivative_symbols.begin(), additional_derivative_symbols.end(), ads) > 1)
80 6 : this->paramError(
81 : "additional_derivative_symbols", "Symbol name '", ads, "' was given more than once.");
82 :
83 : // if all checks passed, add the symbol index to the list
84 78 : _derivative_symbol_table.push_back(idx);
85 : }
86 :
87 : // optimize base function
88 2580 : ParsedMaterialHelper<is_ad>::functionsOptimize(_func_F);
89 :
90 : // generate derivatives
91 2580 : assembleDerivatives();
92 :
93 : // force a value update to get the property at least once and register it for the dependencies
94 5448 : for (auto & mpd : _mat_prop_descriptors)
95 2868 : mpd.value();
96 :
97 : // batch compilation
98 2580 : }
99 :
100 : template <bool is_ad>
101 : void
102 9912 : DerivativeParsedMaterialHelperTempl<is_ad>::recurseMatProps(
103 : unsigned int var, unsigned int order, const MatPropDescriptorList & parent_mpd_list)
104 : {
105 : // quit if we have exceeded the requested derivative order
106 9912 : if (order > _derivative_order)
107 3495 : return;
108 :
109 : // FParser symbol we are deriving w.r.t.
110 6417 : auto derivative_var = _derivative_symbol_table[var];
111 6417 : auto derivative_symbol = _symbol_names[derivative_var];
112 :
113 : // generate parent material property descriptors derivatives
114 6417 : MatPropDescriptorList mpd_list;
115 10485 : for (const auto & parent_mpd : parent_mpd_list)
116 : {
117 : // if this material property does not depend on the variable we are deriving w.r.t. skip it
118 2412 : if (!parent_mpd.dependsOn(derivative_symbol))
119 756 : continue;
120 :
121 : // otherwise add it to _mat_prop_descriptors
122 1656 : FunctionMaterialPropertyDescriptor<is_ad> mpd(parent_mpd);
123 1656 : mpd.addDerivative(derivative_symbol);
124 :
125 : // create a new symbol name for it
126 1656 : std::string newvarname = _dmatvar_base + Moose::stringify(_dmatvar_index++);
127 1656 : mpd.setSymbolName(newvarname);
128 :
129 : // add the new mpd and register it as the current variable derivative of the parent mpd
130 1656 : _bulk_rules.emplace_back(parent_mpd.getSymbolName(), derivative_symbol, newvarname);
131 :
132 : // append to list
133 1656 : mpd_list.push_back(mpd);
134 : }
135 :
136 : // append material property descriptors
137 8073 : for (const auto & mpd : mpd_list)
138 1656 : _mat_prop_descriptors.push_back(mpd);
139 :
140 : // go one order deeper
141 14094 : for (unsigned int i = var; i < _derivative_symbol_table.size(); ++i)
142 7677 : recurseMatProps(i, order + 1, mpd_list);
143 6417 : }
144 :
145 : template <bool is_ad>
146 : void
147 9651 : DerivativeParsedMaterialHelperTempl<is_ad>::recurseDerivative(unsigned int var,
148 : unsigned int order,
149 : const Derivative & parent_derivative)
150 : {
151 : // quit if we have exceeded the requested derivative order
152 9651 : if (order > _derivative_order)
153 3306 : return;
154 :
155 : // FParser symbol we are deriving w.r.t.
156 6345 : auto derivative_var = _derivative_symbol_table[var];
157 6345 : auto derivative_symbol = _symbol_names[derivative_var];
158 :
159 : // current derivative starts off of the parent function
160 6345 : Derivative current;
161 6345 : current._darg_names = parent_derivative._darg_names;
162 :
163 : // the moose variable name goes into the derivative property name
164 6345 : current._darg_names.push_back(derivative_var < _nargs ? _arg_names[derivative_var]
165 : : derivative_symbol);
166 :
167 6345 : current._F = std::make_shared<SymFunction>(*parent_derivative._F);
168 :
169 : // execute derivative
170 6345 : if (current._F->AutoDiff(derivative_symbol) != -1)
171 0 : mooseError("Failed to take order ", order, " derivative in material ", _name);
172 :
173 : // optimize
174 6345 : if (!_disable_fpoptimizer)
175 6345 : current._F->Optimize();
176 :
177 : // proceed only if the derivative is not zero
178 6345 : if (!current._F->isZero() || is_ad)
179 : {
180 : // compile
181 6156 : if (_enable_jit && !current._F->JITCompile())
182 0 : mooseInfo("Failed to JIT compile expression, falling back to byte code interpretation.");
183 :
184 : // go one order deeper
185 13572 : for (unsigned int i = var; i < _derivative_symbol_table.size(); ++i)
186 7416 : recurseDerivative(i, order + 1, current);
187 :
188 : // set up a material property for the derivative
189 6156 : current._mat_prop =
190 6156 : &this->template declarePropertyDerivative<Real, is_ad>(_F_name, current._darg_names);
191 :
192 : // save off current derivative
193 6156 : _derivatives.push_back(current);
194 : }
195 6345 : }
196 :
197 : /**
198 : * Perform construction of all requested derivatives.
199 : */
200 : template <bool is_ad>
201 : void
202 2580 : DerivativeParsedMaterialHelperTempl<is_ad>::assembleDerivatives()
203 : {
204 : // need to check for zero derivatives here, otherwise at least one order is generated
205 2580 : if (_derivative_order < 1)
206 237 : return;
207 :
208 : // if we are not on thread 0 we fetch all data from the thread 0 copy that already did all the
209 : // work
210 2580 : if (_tid > 0)
211 : {
212 : // get the master object from thread 0
213 237 : const MaterialWarehouse & material_warehouse = this->_fe_problem.getMaterialWarehouse();
214 237 : const MooseObjectWarehouse<MaterialBase> & warehouse = material_warehouse[_material_data_type];
215 :
216 237 : MooseSharedPointer<DerivativeParsedMaterialHelperTempl> master =
217 : MooseSharedNamespace::dynamic_pointer_cast<DerivativeParsedMaterialHelperTempl>(
218 : warehouse.getActiveObject(name()));
219 :
220 : // copy parsers and declare properties
221 801 : for (const auto & D : master->_derivatives)
222 : {
223 564 : Derivative newderivative;
224 564 : newderivative._mat_prop =
225 564 : &this->template declarePropertyDerivative<Real, is_ad>(_F_name, D._darg_names);
226 564 : newderivative._F = std::make_shared<SymFunction>(*D._F);
227 564 : _derivatives.push_back(newderivative);
228 : }
229 :
230 : // copy coupled material properties
231 237 : auto start = _mat_prop_descriptors.size();
232 399 : for (MooseIndex(master->_mat_prop_descriptors) i = start;
233 399 : i < master->_mat_prop_descriptors.size();
234 : ++i)
235 : {
236 162 : FunctionMaterialPropertyDescriptor<is_ad> newdescriptor(master->_mat_prop_descriptors[i],
237 : this);
238 162 : _mat_prop_descriptors.push_back(newdescriptor);
239 : }
240 :
241 : // size parameter buffer
242 237 : _func_params.resize(master->_func_params.size());
243 237 : return;
244 237 : }
245 :
246 : // generate all coupled material property mappings
247 4578 : for (unsigned int i = 0; i < _derivative_symbol_table.size(); ++i)
248 2235 : recurseMatProps(i, 1, _mat_prop_descriptors);
249 :
250 : // bulk register material property derivative rules to avoid repeated calls
251 : // to the (slow) AddVariable method
252 2343 : if (!_bulk_rules.empty())
253 : {
254 552 : std::string vars = _bulk_rules[0]._child;
255 1656 : for (MooseIndex(_bulk_rules) i = 1; i < _bulk_rules.size(); ++i)
256 1104 : vars += ',' + _bulk_rules[i]._child;
257 552 : _func_F->AddVariable(vars);
258 :
259 2208 : for (const auto & rule : _bulk_rules)
260 1656 : _func_F->RegisterDerivative(rule._parent, rule._var, rule._child);
261 552 : }
262 :
263 : // generate all derivatives
264 2343 : Derivative root;
265 2343 : root._F = _func_F;
266 2343 : root._mat_prop = nullptr;
267 4578 : for (unsigned int i = 0; i < _derivative_symbol_table.size(); ++i)
268 2235 : recurseDerivative(i, 1, root);
269 :
270 : // increase the parameter buffer to provide storage for the material property derivatives
271 4686 : _func_params.resize(_nargs + _mat_prop_descriptors.size() + _postprocessor_values.size() +
272 2343 : _extra_symbols.size() + _functors.size());
273 2343 : }
274 :
275 : template <bool is_ad>
276 : void
277 0 : DerivativeParsedMaterialHelperTempl<is_ad>::initQpStatefulProperties()
278 : {
279 0 : computeQpProperties();
280 0 : }
281 :
282 : template <bool is_ad>
283 : void
284 2195969 : DerivativeParsedMaterialHelperTempl<is_ad>::computeQpProperties()
285 : {
286 2195969 : ParsedMaterialHelper<is_ad>::computeQpProperties();
287 :
288 : // set derivatives
289 8825104 : for (auto & D : _derivatives)
290 6629162 : (*D._mat_prop)[_qp] = evaluate(D._F, _name);
291 2195942 : }
292 :
293 : // explicit instantiation
294 : template class DerivativeParsedMaterialHelperTempl<false>;
295 : template class DerivativeParsedMaterialHelperTempl<true>;
|
