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 32904 : DerivativeParsedMaterialHelperTempl<is_ad>::validParams()
20 : {
21 32904 : InputParameters params = ParsedMaterialHelper<is_ad>::validParams();
22 65808 : params.addClassDescription("Parsed Function Material with automatic derivatives.");
23 131616 : params.addParam<unsigned int>("derivative_order", 3, "Maximum order of derivatives taken");
24 98712 : 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 32904 : return params;
30 0 : }
31 :
32 : template <bool is_ad>
33 2619 : 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 2619 : _derivative_order(this->template getParam<unsigned int>("derivative_order")),
39 5238 : _dmatvar_base("matpropautoderiv"),
40 7857 : _dmatvar_index(0)
41 : {
42 2619 : }
43 :
44 : template <bool is_ad>
45 : void
46 2619 : DerivativeParsedMaterialHelperTempl<is_ad>::functionsPostParse()
47 : {
48 : // set up the list of all FParser symbols to derive w.r.t.
49 5238 : 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 2619 : _derivative_symbol_table.reserve(_nargs + additional_derivative_symbols.size());
55 4824 : for (unsigned i = 0; i < _nargs; ++i)
56 2205 : _derivative_symbol_table.push_back(i);
57 :
58 : // then we add the additional derivative symbols
59 2703 : for (auto & ads : additional_derivative_symbols)
60 : {
61 : // find the index of the symbol if it exists
62 96 : auto it = std::find(_symbol_names.begin(), _symbol_names.end(), ads);
63 96 : if (it == _symbol_names.end())
64 8 : this->paramError("additional_derivative_symbols", "Invalid symbol name '", ads, "'.");
65 92 : 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 92 : if (idx < _nargs)
69 8 : 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 88 : if (std::count(
79 88 : additional_derivative_symbols.begin(), additional_derivative_symbols.end(), ads) > 1)
80 8 : 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 84 : _derivative_symbol_table.push_back(idx);
85 : }
86 :
87 : // optimize base function
88 2607 : ParsedMaterialHelper<is_ad>::functionsOptimize(_func_F);
89 :
90 : // generate derivatives
91 2607 : assembleDerivatives();
92 :
93 : // force a value update to get the property at least once and register it for the dependencies
94 5739 : for (auto & mpd : _mat_prop_descriptors)
95 3132 : mpd.value();
96 :
97 : // batch compilation
98 2607 : }
99 :
100 : template <bool is_ad>
101 : void
102 9957 : 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 9957 : if (order > _derivative_order)
107 3303 : return;
108 :
109 : // FParser symbol we are deriving w.r.t.
110 6654 : auto derivative_var = _derivative_symbol_table[var];
111 6654 : auto derivative_symbol = _symbol_names[derivative_var];
112 :
113 : // generate parent material property descriptors derivatives
114 6654 : MatPropDescriptorList mpd_list;
115 11109 : 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 2637 : if (!parent_mpd.dependsOn(derivative_symbol))
119 819 : continue;
120 :
121 : // otherwise add it to _mat_prop_descriptors
122 1818 : FunctionMaterialPropertyDescriptor<is_ad> mpd(parent_mpd);
123 1818 : mpd.addDerivative(derivative_symbol);
124 :
125 : // create a new symbol name for it
126 1818 : std::string newvarname = _dmatvar_base + Moose::stringify(_dmatvar_index++);
127 1818 : mpd.setSymbolName(newvarname);
128 :
129 : // add the new mpd and register it as the current variable derivative of the parent mpd
130 1818 : _bulk_rules.emplace_back(parent_mpd.getSymbolName(), derivative_symbol, newvarname);
131 :
132 : // append to list
133 1818 : mpd_list.push_back(mpd);
134 : }
135 :
136 : // append material property descriptors
137 8472 : for (const auto & mpd : mpd_list)
138 1818 : _mat_prop_descriptors.push_back(mpd);
139 :
140 : // go one order deeper
141 14517 : for (unsigned int i = var; i < _derivative_symbol_table.size(); ++i)
142 7863 : recurseMatProps(i, order + 1, mpd_list);
143 6654 : }
144 :
145 : template <bool is_ad>
146 : void
147 9672 : 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 9672 : if (order > _derivative_order)
153 3096 : return;
154 :
155 : // FParser symbol we are deriving w.r.t.
156 6576 : auto derivative_var = _derivative_symbol_table[var];
157 6576 : auto derivative_symbol = _symbol_names[derivative_var];
158 :
159 : // current derivative starts off of the parent function
160 6576 : Derivative current;
161 6576 : current._darg_names = parent_derivative._darg_names;
162 :
163 : // the moose variable name goes into the derivative property name
164 6576 : current._darg_names.push_back(derivative_var < _nargs ? _arg_names[derivative_var]
165 : : derivative_symbol);
166 :
167 6576 : current._F = std::make_shared<SymFunction>(*parent_derivative._F);
168 :
169 : // execute derivative
170 6576 : if (current._F->AutoDiff(derivative_symbol) != -1)
171 0 : mooseError("Failed to take order ", order, " derivative in material ", _name);
172 :
173 : // optimize
174 6576 : if (!_disable_fpoptimizer)
175 6576 : current._F->Optimize();
176 :
177 : // proceed only if the derivative is not zero
178 6576 : if (!current._F->isZero() || is_ad)
179 : {
180 : // compile
181 6369 : 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 13947 : for (unsigned int i = var; i < _derivative_symbol_table.size(); ++i)
186 7578 : recurseDerivative(i, order + 1, current);
187 :
188 : // set up a material property for the derivative
189 6369 : current._mat_prop =
190 6369 : &this->template declarePropertyDerivative<Real, is_ad>(_F_name, current._darg_names);
191 :
192 : // save off current derivative
193 6369 : _derivatives.push_back(current);
194 : }
195 6576 : }
196 :
197 : /**
198 : * Perform construction of all requested derivatives.
199 : */
200 : template <bool is_ad>
201 : void
202 2607 : DerivativeParsedMaterialHelperTempl<is_ad>::assembleDerivatives()
203 : {
204 : // need to check for zero derivatives here, otherwise at least one order is generated
205 2607 : if (_derivative_order < 1)
206 222 : 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 2607 : if (_tid > 0)
211 : {
212 : // get the master object from thread 0
213 222 : const MaterialWarehouse & material_warehouse = this->_fe_problem.getMaterialWarehouse();
214 222 : const MooseObjectWarehouse<MaterialBase> & warehouse = material_warehouse[_material_data_type];
215 :
216 222 : MooseSharedPointer<DerivativeParsedMaterialHelperTempl> master =
217 : MooseSharedNamespace::dynamic_pointer_cast<DerivativeParsedMaterialHelperTempl>(
218 : warehouse.getActiveObject(name()));
219 :
220 : // copy parsers and declare properties
221 759 : for (const auto & D : master->_derivatives)
222 : {
223 537 : Derivative newderivative;
224 537 : newderivative._mat_prop =
225 537 : &this->template declarePropertyDerivative<Real, is_ad>(_F_name, D._darg_names);
226 537 : newderivative._F = std::make_shared<SymFunction>(*D._F);
227 537 : _derivatives.push_back(newderivative);
228 : }
229 :
230 : // copy coupled material properties
231 222 : auto start = _mat_prop_descriptors.size();
232 384 : for (MooseIndex(master->_mat_prop_descriptors) i = start;
233 384 : 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 222 : _func_params.resize(master->_func_params.size());
243 222 : return;
244 222 : }
245 :
246 : // generate all coupled material property mappings
247 4479 : for (unsigned int i = 0; i < _derivative_symbol_table.size(); ++i)
248 2094 : 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 2385 : if (!_bulk_rules.empty())
253 : {
254 606 : std::string vars = _bulk_rules[0]._child;
255 1818 : for (MooseIndex(_bulk_rules) i = 1; i < _bulk_rules.size(); ++i)
256 1212 : vars += ',' + _bulk_rules[i]._child;
257 606 : _func_F->AddVariable(vars);
258 :
259 2424 : for (const auto & rule : _bulk_rules)
260 1818 : _func_F->RegisterDerivative(rule._parent, rule._var, rule._child);
261 606 : }
262 :
263 : // generate all derivatives
264 2385 : Derivative root;
265 2385 : root._F = _func_F;
266 2385 : root._mat_prop = nullptr;
267 4479 : for (unsigned int i = 0; i < _derivative_symbol_table.size(); ++i)
268 2094 : recurseDerivative(i, 1, root);
269 :
270 : // increase the parameter buffer to provide storage for the material property derivatives
271 4770 : _func_params.resize(_nargs + _mat_prop_descriptors.size() + _postprocessor_values.size() +
272 2385 : _extra_symbols.size() + _functors.size());
273 2385 : }
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 2418277 : DerivativeParsedMaterialHelperTempl<is_ad>::computeQpProperties()
285 : {
286 2418277 : ParsedMaterialHelper<is_ad>::computeQpProperties();
287 :
288 : // set derivatives
289 9750314 : for (auto & D : _derivatives)
290 7332070 : (*D._mat_prop)[_qp] = evaluate(D._F, _name);
291 2418244 : }
292 :
293 : // explicit instantiation
294 : template class DerivativeParsedMaterialHelperTempl<false>;
295 : template class DerivativeParsedMaterialHelperTempl<true>;
|
