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