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 : #ifdef MOOSE_MFEM_ENABLED
11 :
12 : #include "MFEMProblem.h"
13 : #include "MFEMVariable.h"
14 : #include "MFEMIndicator.h"
15 : #include "MFEMSubMesh.h"
16 : #include "MFEMFunctorMaterial.h"
17 : #include "MFEMExecutedObject.h"
18 : #include "MFEMVectorUtils.h"
19 : #include "MFEMFESpaceHierarchy.h"
20 : #include "Postprocessor.h"
21 : #include "VectorPostprocessor.h"
22 : #include "MFEMNonlinearSolverBase.h"
23 : #include "DependencyResolver.h"
24 : #include "MooseUtils.h"
25 :
26 : #include "libmesh/string_to_enum.h"
27 :
28 : #include <vector>
29 : #include <algorithm>
30 : #include <map>
31 : #include <deque>
32 : #include <sstream>
33 :
34 : registerMooseObject("MooseApp", MFEMProblem);
35 :
36 : namespace
37 : {
38 : std::vector<MFEMSolverName>
39 2096 : getMFEMSolverDependencies(const InputParameters & parameters)
40 : {
41 2096 : std::vector<MFEMSolverName> dependencies;
42 :
43 52706 : for (const auto & [param_name, _] : parameters)
44 : {
45 50610 : if (parameters.isPrivate(param_name) || !parameters.isParamValid(param_name))
46 32235 : continue;
47 :
48 18375 : if (parameters.have_parameter<MFEMSolverName>(param_name))
49 972 : dependencies.push_back(parameters.get<MFEMSolverName>(param_name));
50 17403 : else if (parameters.have_parameter<std::vector<MFEMSolverName>>(param_name))
51 : {
52 9 : const auto & names = parameters.get<std::vector<MFEMSolverName>>(param_name);
53 9 : dependencies.insert(dependencies.end(), names.begin(), names.end());
54 : }
55 : }
56 :
57 2096 : return dependencies;
58 0 : }
59 : }
60 :
61 : InputParameters
62 8194 : MFEMProblem::validParams()
63 : {
64 8194 : InputParameters params = ExternalProblem::validParams();
65 16388 : params.addClassDescription("Problem type for building and solving the finite element problem "
66 : "using the MFEM finite element library.");
67 32776 : MooseEnum numeric_types("real complex", "real");
68 24582 : params.addParam<MooseEnum>("numeric_type", numeric_types, "Number type used for the problem");
69 :
70 16388 : return params;
71 8194 : }
72 :
73 1971 : MFEMProblem::MFEMProblem(const InputParameters & params)
74 5913 : : ExternalProblem(params), _num_type{static_cast<int>(getParam<MooseEnum>("numeric_type"))}
75 : {
76 : // Initialise Hypre for all MFEM problems.
77 1971 : mfem::Hypre::Init();
78 : // Disable multithreading for all MFEM problems (including any libMesh or MFEM subapps).
79 1971 : libMesh::libMeshPrivateData::_n_threads = 1;
80 : #ifdef LIBMESH_HAVE_OPENMP
81 1971 : omp_set_num_threads(1);
82 : #endif
83 1971 : setMesh();
84 1971 : }
85 :
86 : void
87 1617 : MFEMProblem::initialSetup()
88 : {
89 1617 : ExternalProblem::initialSetup();
90 :
91 : // MFEM indicators create their estimators during addIndicator(); markers still need an explicit
92 : // setup pass because they are no longer initialized through the libMesh/MOOSE user-object path.
93 1617 : std::vector<MFEMRefinementMarker *> markers;
94 1617 : theWarehouse().query().condition<AttribSystem>("Marker").queryInto(markers);
95 1643 : for (auto marker : markers)
96 26 : marker->initialSetup();
97 1617 : }
98 :
99 : void
100 11251 : MFEMProblem::execute(const ExecFlagType & exec_type)
101 : {
102 11251 : setCurrentExecuteOnFlag(exec_type);
103 11251 : executeMFEMObjects(exec_type);
104 :
105 11251 : ExternalProblem::execute(exec_type);
106 11251 : }
107 :
108 : void
109 1971 : MFEMProblem::setMesh()
110 : {
111 1971 : auto pmesh = mesh().getMFEMParMeshPtr();
112 1971 : getProblemData().pmesh = pmesh;
113 1971 : getProblemData().comm = pmesh->GetComm();
114 1971 : getProblemData().num_procs = pmesh->GetNRanks();
115 1971 : getProblemData().myid = pmesh->GetMyRank();
116 1971 : }
117 :
118 : void
119 26 : MFEMProblem::addIndicator(const std::string & indicator_type,
120 : const std::string & name,
121 : InputParameters & parameters)
122 : {
123 52 : auto estimator = addObject<MFEMIndicator>(indicator_type, name, parameters).front();
124 :
125 : // construct the estimator itself
126 26 : estimator->createEstimator();
127 26 : }
128 :
129 : void
130 26 : MFEMProblem::addMarker(const std::string & marker_type,
131 : const std::string & name,
132 : InputParameters & parameters)
133 : {
134 52 : getProblemData().refiner = addObject<MFEMRefinementMarker>(marker_type, name, parameters).front();
135 26 : }
136 :
137 : void
138 2096 : MFEMProblem::addMFEMSolver(const std::string & solver_type,
139 : const std::string & name,
140 : InputParameters & parameters)
141 : {
142 4192 : if (!_mfem_solver_definitions.emplace(name, MFEMSolverDefinition{solver_type, ¶meters})
143 2096 : .second)
144 0 : mooseError("Multiple MFEM solvers named '", name, "' were provided.");
145 2096 : }
146 :
147 : void
148 1617 : MFEMProblem::resolveMFEMSolvers()
149 : {
150 1617 : if (_mfem_solver_definitions.empty())
151 545 : return;
152 :
153 1072 : DependencyResolver<std::string> resolver;
154 :
155 3168 : for (auto & [solver_name, definition] : _mfem_solver_definitions)
156 : {
157 2096 : const auto dependencies = getMFEMSolverDependencies(*definition.parameters);
158 2096 : if (dependencies.empty())
159 1124 : resolver.addNode(solver_name);
160 :
161 3077 : for (const auto & dependency_name : dependencies)
162 : {
163 981 : auto dependency_it = _mfem_solver_definitions.find(dependency_name);
164 981 : if (dependency_it == _mfem_solver_definitions.end())
165 0 : mooseError("MFEM solver '",
166 : solver_name,
167 : "' references MFEM solver '",
168 : dependency_name,
169 : "', but no solver with that name was provided in the [Solvers] block.");
170 :
171 981 : dependency_it->second.referenced = true;
172 981 : resolver.addEdge(dependency_name, solver_name);
173 : }
174 2096 : }
175 :
176 1072 : const std::vector<std::string> * sorted_solver_names = nullptr;
177 : try
178 : {
179 1072 : sorted_solver_names = &resolver.getSortedValues();
180 : }
181 0 : catch (CyclicDependencyException<std::string> & e)
182 : {
183 0 : mooseError("Cyclic MFEM solver dependency detected: ",
184 0 : MooseUtils::join(e.getCyclicDependencies(), " <- "));
185 0 : }
186 :
187 1072 : auto & problem_data = getProblemData();
188 1072 : problem_data.jacobian_solver = nullptr;
189 1072 : problem_data.nonlinear_solver = nullptr;
190 :
191 3168 : for (const auto & solver_name : *sorted_solver_names)
192 : {
193 2096 : auto & definition = libmesh_map_find(_mfem_solver_definitions, solver_name);
194 : auto solver =
195 6288 : addObject<Moose::MFEM::SolverBase>(definition.type, solver_name, *definition.parameters)
196 2096 : .front();
197 :
198 2096 : if (definition.referenced)
199 981 : continue;
200 :
201 1115 : if (auto lin_solver = std::dynamic_pointer_cast<Moose::MFEM::LinearSolverBase>(solver))
202 : {
203 1065 : if (problem_data.jacobian_solver)
204 0 : mooseError("Multiple MFEM linear solver drivers provided. '",
205 0 : problem_data.jacobian_solver->name(),
206 : "' and '",
207 0 : lin_solver->name(),
208 : "' are not referenced by another MFEM solver.");
209 1065 : problem_data.jacobian_solver = lin_solver;
210 : }
211 50 : else if (auto nonlinear_solver =
212 50 : std::dynamic_pointer_cast<Moose::MFEM::NonlinearSolverBase>(solver);
213 50 : nonlinear_solver)
214 : {
215 50 : if (problem_data.nonlinear_solver)
216 0 : mooseError("Multiple MFEM nonlinear solver drivers provided. '",
217 0 : problem_data.nonlinear_solver->name(),
218 : "' and '",
219 0 : nonlinear_solver->name(),
220 : "' are not referenced by another MFEM solver.");
221 50 : problem_data.nonlinear_solver = nonlinear_solver;
222 : }
223 : else
224 0 : mooseError("Unsupported MFEM solver object type '",
225 0 : solver->type(),
226 : "' for solver '",
227 0 : solver->name(),
228 1165 : "'.");
229 2096 : }
230 :
231 1072 : _mfem_solver_definitions.clear();
232 1072 : }
233 :
234 : void
235 1448 : MFEMProblem::addBoundaryCondition(const std::string & bc_name,
236 : const std::string & name,
237 : InputParameters & parameters)
238 : {
239 2896 : auto bc = addObject<MFEMBoundaryCondition>(bc_name, name, parameters).front();
240 1448 : const auto & mfem_bc = *bc;
241 :
242 1448 : if (dynamic_cast<const MFEMIntegratedBC *>(&mfem_bc))
243 : {
244 58 : auto integrated_bc = std::dynamic_pointer_cast<MFEMIntegratedBC>(bc);
245 : auto eqsys =
246 58 : std::dynamic_pointer_cast<Moose::MFEM::EquationSystem>(getProblemData().eqn_system);
247 58 : if (eqsys)
248 58 : eqsys->AddIntegratedBC(std::move(integrated_bc));
249 : else
250 0 : mooseError("Cannot add integrated BC with name '" + name +
251 : "' because there is no corresponding equation system.");
252 58 : }
253 1390 : else if (dynamic_cast<const MFEMComplexIntegratedBC *>(&mfem_bc))
254 : {
255 14 : auto integrated_bc = std::dynamic_pointer_cast<MFEMComplexIntegratedBC>(bc);
256 : auto eqsys =
257 14 : std::dynamic_pointer_cast<Moose::MFEM::ComplexEquationSystem>(getProblemData().eqn_system);
258 14 : if (eqsys)
259 14 : eqsys->AddComplexIntegratedBC(std::move(integrated_bc));
260 : else
261 0 : mooseError("Cannot add complex integrated BC with name '" + name +
262 : "' because there is no corresponding equation system.");
263 14 : }
264 1376 : else if (dynamic_cast<const MFEMComplexEssentialBC *>(&mfem_bc))
265 : {
266 63 : auto essential_bc = std::dynamic_pointer_cast<MFEMComplexEssentialBC>(bc);
267 : auto eqsys =
268 63 : std::dynamic_pointer_cast<Moose::MFEM::ComplexEquationSystem>(getProblemData().eqn_system);
269 63 : if (eqsys)
270 63 : eqsys->AddComplexEssentialBCs(std::move(essential_bc));
271 : else
272 0 : mooseError("Cannot add boundary condition with name '" + name +
273 : "' because there is no corresponding equation system.");
274 63 : }
275 1313 : else if (dynamic_cast<const MFEMEssentialBC *>(&mfem_bc))
276 : {
277 1313 : auto essential_bc = std::dynamic_pointer_cast<MFEMEssentialBC>(bc);
278 : auto eqsys =
279 1313 : std::dynamic_pointer_cast<Moose::MFEM::EquationSystem>(getProblemData().eqn_system);
280 1313 : if (eqsys)
281 1313 : eqsys->AddEssentialBC(std::move(essential_bc));
282 : else
283 0 : mooseError("Cannot add boundary condition with name '" + name +
284 : "' because there is no corresponding equation system.");
285 1313 : }
286 : else
287 : {
288 0 : mooseError("Unsupported bc of type '", bc_name, "' and name '", name, "' detected.");
289 : }
290 1448 : }
291 :
292 : void
293 0 : MFEMProblem::addMaterial(const std::string &, const std::string &, InputParameters &)
294 : {
295 0 : mooseError(
296 : "MFEM materials must be added through the 'FunctorMaterials' block and not 'Materials'");
297 : }
298 :
299 : void
300 339 : MFEMProblem::addFunctorMaterial(const std::string & material_name,
301 : const std::string & name,
302 : InputParameters & parameters)
303 : {
304 694 : addObject<MFEMFunctorMaterial>(material_name, name, parameters);
305 323 : }
306 :
307 : void
308 2749 : MFEMProblem::addFESpace(const std::string & type,
309 : const std::string & name,
310 : InputParameters & parameters)
311 : {
312 2749 : if (getProblemData().fespace_hierarchies.Has(name))
313 0 : mooseError("Cannot add FESpace '",
314 : name,
315 : "': a FESpaceHierarchy with the same name already exists. "
316 : "FESpaces and FESpaceHierarchies share the fespaces namespace.");
317 :
318 5498 : auto & mfem_fespace = *addObject<MFEMFESpace>(type, name, parameters).front();
319 :
320 : // Register fespace and associated fe collection.
321 2749 : getProblemData().fecs.Register(name, mfem_fespace.getFEC());
322 2749 : getProblemData().fespaces.Register(name, mfem_fespace.getFESpace());
323 2749 : }
324 :
325 : void
326 9 : MFEMProblem::addFESpaceHierarchy(const std::string & type,
327 : const std::string & name,
328 : InputParameters & parameters)
329 : {
330 9 : if (getProblemData().fespaces.Has(name))
331 0 : mooseError("Cannot add FESpaceHierarchy '",
332 : name,
333 : "': a FESpace with the same name already exists. "
334 : "FESpaces and FESpaceHierarchies share the fespaces namespace.");
335 :
336 18 : auto hierarchy_obj = addObject<MFEMFESpaceHierarchy>(type, name, parameters).front();
337 9 : auto hierarchy_shared = hierarchy_obj->getHierarchyShared();
338 : // Register the hierarchy for co-ownership by solvers.
339 9 : getProblemData().fespace_hierarchies.Register(name, hierarchy_shared);
340 : // Register the finest-level FESpace in fespaces under the hierarchy name so that
341 : // variables can say `fespace = <hierarchy_name>` without a separate FESpace definition.
342 : // The aliasing shared_ptr keeps the hierarchy alive as long as this entry lives.
343 : auto finest = std::shared_ptr<mfem::ParFiniteElementSpace>(
344 : hierarchy_shared,
345 : &static_cast<mfem::ParFiniteElementSpace &>(
346 9 : hierarchy_obj->getHierarchy().GetFinestFESpace()));
347 9 : getProblemData().fespaces.Register(name, finest);
348 9 : }
349 :
350 : void
351 1892 : MFEMProblem::addVariable(const std::string & var_type,
352 : const std::string & var_name,
353 : InputParameters & parameters)
354 : {
355 1892 : addGridFunction(var_type, var_name, parameters);
356 : // MOOSE variables store DoFs for the trial variable and its time derivatives up to second order;
357 : // MFEM GridFunctions store data for only one set of DoFs each, so we must add additional
358 : // GridFunctions for time derivatives.
359 1892 : if (isTransient())
360 : {
361 : const auto time_derivative_var_name =
362 225 : getMFEMObject<MFEMVariable>("MooseVariableBase", var_name).getTimeDerivativeName();
363 225 : getProblemData().time_derivative_map.addTimeDerivativeAssociation(var_name,
364 : time_derivative_var_name);
365 225 : addGridFunction(var_type, time_derivative_var_name, parameters);
366 225 : }
367 1892 : }
368 :
369 : void
370 3869 : MFEMProblem::addGridFunction(const std::string & var_type,
371 : const std::string & var_name,
372 : InputParameters & parameters)
373 : {
374 :
375 3869 : if (var_type == "MFEMVariable" || var_type == "MFEMComplexVariable")
376 : {
377 : // Add MFEM variable directly.
378 3764 : if (var_type == "MFEMComplexVariable")
379 189 : addObject<MFEMComplexVariable>(var_type, var_name, parameters);
380 : else
381 11103 : addObject<MFEMVariable>(var_type, var_name, parameters);
382 : }
383 : else
384 : {
385 : // Add MOOSE variable.
386 105 : ExternalProblem::addVariable(var_type, var_name, parameters);
387 :
388 : // Add MFEM variable indirectly ("gridfunction").
389 105 : InputParameters mfem_variable_params = addMFEMFESpaceFromMOOSEVariable(parameters);
390 420 : addObject<MFEMVariable>("MFEMVariable", var_name, mfem_variable_params);
391 105 : }
392 :
393 : // Register gridfunction.
394 3869 : if (var_type == "MFEMComplexVariable")
395 : {
396 : MFEMComplexVariable & mfem_variable =
397 63 : getMFEMObject<MFEMComplexVariable>("MooseVariableBase", var_name);
398 63 : getProblemData().cmplx_gridfunctions.Register(var_name, mfem_variable.getComplexGridFunction());
399 63 : mfem_variable.declareCoefficients();
400 : }
401 : else // must be real, but may have been set up indirectly from a MOOSE variable
402 : {
403 3806 : MFEMVariable & mfem_variable = getMFEMObject<MFEMVariable>("MooseVariableBase", var_name);
404 3806 : getProblemData().gridfunctions.Register(var_name, mfem_variable.getGridFunction());
405 3806 : mfem_variable.declareCoefficients();
406 : }
407 3869 : }
408 :
409 : void
410 1596 : MFEMProblem::addAuxVariable(const std::string & var_type,
411 : const std::string & var_name,
412 : InputParameters & parameters)
413 : {
414 : // We handle MFEM AuxVariables just like MFEM Variables, except
415 : // we do not add additional GridFunctions for time derivatives.
416 1596 : addGridFunction(var_type, var_name, parameters);
417 1596 : }
418 :
419 : void
420 518 : MFEMProblem::addAuxKernel(const std::string & kernel_name,
421 : const std::string & name,
422 : InputParameters & parameters)
423 : {
424 1036 : addObject<MFEMExecutedObject>(kernel_name, name, parameters);
425 518 : }
426 :
427 : void
428 1967 : MFEMProblem::addKernel(const std::string & kernel_name,
429 : const std::string & name,
430 : InputParameters & parameters)
431 : {
432 3934 : auto kernel = addObject<MFEMKernel>(kernel_name, name, parameters).front();
433 1967 : const auto & kernel_object = *kernel;
434 :
435 1967 : if (dynamic_cast<const MFEMComplexKernel *>(&kernel_object))
436 : {
437 63 : auto complex_kernel = std::dynamic_pointer_cast<MFEMComplexKernel>(kernel);
438 : auto eqsys =
439 63 : std::dynamic_pointer_cast<Moose::MFEM::ComplexEquationSystem>(getProblemData().eqn_system);
440 63 : if (eqsys)
441 63 : eqsys->AddComplexKernel(std::move(complex_kernel));
442 : else
443 0 : mooseError("Cannot add complex kernel with name '" + name +
444 : "' because there is no corresponding equation system.");
445 63 : }
446 : else
447 : {
448 : auto eqsys =
449 1904 : std::dynamic_pointer_cast<Moose::MFEM::EquationSystem>(getProblemData().eqn_system);
450 1904 : if (eqsys)
451 1904 : eqsys->AddKernel(std::move(kernel));
452 : else
453 0 : mooseError("Cannot add kernel with name '" + name +
454 : "' because there is no corresponding equation system.");
455 1904 : }
456 1967 : }
457 :
458 : void
459 63 : MFEMProblem::addRealComponentToKernel(const std::string & kernel_name,
460 : const std::string & name,
461 : InputParameters & parameters)
462 : {
463 : auto parent_ptr = std::dynamic_pointer_cast<MFEMComplexKernel>(
464 63 : getMFEMObject<MFEMComplexKernel>("Kernel", name).getSharedPtr());
465 252 : parameters.set<VariableName>("variable") = parent_ptr->getParam<VariableName>("variable");
466 126 : auto kernel_ptr = addObject<MFEMKernel>(kernel_name, name + "_real", parameters).front();
467 63 : parent_ptr->setRealKernel(kernel_ptr);
468 63 : }
469 :
470 : void
471 56 : MFEMProblem::addImagComponentToKernel(const std::string & kernel_name,
472 : const std::string & name,
473 : InputParameters & parameters)
474 : {
475 : auto parent_ptr = std::dynamic_pointer_cast<MFEMComplexKernel>(
476 56 : getMFEMObject<MFEMComplexKernel>("Kernel", name).getSharedPtr());
477 224 : parameters.set<VariableName>("variable") = parent_ptr->getParam<VariableName>("variable");
478 112 : auto kernel_ptr = addObject<MFEMKernel>(kernel_name, name + "_imag", parameters).front();
479 56 : parent_ptr->setImagKernel(kernel_ptr);
480 56 : }
481 :
482 : void
483 0 : MFEMProblem::addRealComponentToBC(const std::string & kernel_name,
484 : const std::string & name,
485 : InputParameters & parameters)
486 : {
487 : auto parent_ptr = std::dynamic_pointer_cast<MFEMComplexIntegratedBC>(
488 0 : getMFEMObject<MFEMComplexIntegratedBC>("BoundaryCondition", name).getSharedPtr());
489 0 : parameters.set<VariableName>("variable") = parent_ptr->getParam<VariableName>("variable");
490 0 : parameters.set<std::vector<BoundaryName>>("boundary") =
491 0 : parent_ptr->getParam<std::vector<BoundaryName>>("boundary");
492 : auto bc_ptr = std::dynamic_pointer_cast<MFEMIntegratedBC>(
493 0 : addObject<MFEMBoundaryCondition>(kernel_name, name + "_real", parameters).front());
494 0 : parent_ptr->setRealBC(bc_ptr);
495 0 : }
496 :
497 : void
498 0 : MFEMProblem::addImagComponentToBC(const std::string & kernel_name,
499 : const std::string & name,
500 : InputParameters & parameters)
501 : {
502 : auto parent_ptr = std::dynamic_pointer_cast<MFEMComplexIntegratedBC>(
503 0 : getMFEMObject<MFEMComplexIntegratedBC>("BoundaryCondition", name).getSharedPtr());
504 0 : parameters.set<VariableName>("variable") = parent_ptr->getParam<VariableName>("variable");
505 0 : parameters.set<std::vector<BoundaryName>>("boundary") =
506 0 : parent_ptr->getParam<std::vector<BoundaryName>>("boundary");
507 : auto bc_ptr = std::dynamic_pointer_cast<MFEMIntegratedBC>(
508 0 : addObject<MFEMBoundaryCondition>(kernel_name, name + "_imag", parameters).front());
509 0 : parent_ptr->setImagBC(bc_ptr);
510 0 : }
511 :
512 : int
513 339 : vectorFunctionDim(const std::string & type, const InputParameters & parameters)
514 : {
515 678 : if (parameters.isParamSetByUser("expression_z"))
516 309 : return 3;
517 90 : if (parameters.isParamSetByUser("expression_y") || type == "LevelSetOlssonVortex")
518 28 : return 2;
519 4 : if (parameters.isParamSetByUser("expression_x"))
520 2 : return 1;
521 :
522 0 : return 3;
523 : }
524 :
525 : const std::vector<std::string> SCALAR_FUNCS = {"Axisymmetric2D3DSolutionFunction",
526 : "BicubicSplineFunction",
527 : "CoarsenedPiecewiseLinear",
528 : "CompositeFunction",
529 : "ConstantFunction",
530 : "ImageFunction",
531 : "ParsedFunction",
532 : "ParsedGradFunction",
533 : "PeriodicFunction",
534 : "PiecewiseBilinear",
535 : "PiecewiseConstant",
536 : "PiecewiseConstantFromCSV",
537 : "PiecewiseLinear",
538 : "PiecewiseLinearFromVectorPostprocessor",
539 : "PiecewiseMultiInterpolation",
540 : "PiecewiseMulticonstant",
541 : "SolutionFunction",
542 : "SplineFunction",
543 : "FunctionSeries",
544 : "LevelSetOlssonBubble",
545 : "LevelSetOlssonPlane",
546 : "NearestReporterCoordinatesFunction",
547 : "ParameterMeshFunction",
548 : "ParsedOptimizationFunction",
549 : "FourierNoise",
550 : "MovingPlanarFront",
551 : "MultiControlDrumFunction",
552 : "Grad2ParsedFunction",
553 : "GradParsedFunction",
554 : "ScaledAbsDifferenceDRLRewardFunction",
555 : "CircularAreaHydraulicDiameterFunction",
556 : "CosineHumpFunction",
557 : "CosineTransitionFunction",
558 : "CubicTransitionFunction",
559 : "GeneralizedCircumference",
560 : "PiecewiseFunction",
561 : "TimeRampFunction"},
562 : VECTOR_FUNCS = {"ParsedVectorFunction", "LevelSetOlssonVortex"};
563 :
564 : void
565 1575 : MFEMProblem::addFunction(const std::string & type,
566 : const std::string & name,
567 : InputParameters & parameters)
568 : {
569 1575 : ExternalProblem::addFunction(type, name, parameters);
570 1575 : auto & func = getFunction(name);
571 : // FIXME: Do we want to have optimised versions for when functions
572 : // are only of space or only of time.
573 1575 : if (std::find(SCALAR_FUNCS.begin(), SCALAR_FUNCS.end(), type) != SCALAR_FUNCS.end())
574 : {
575 1106 : getCoefficients().declareScalar<mfem::FunctionCoefficient>(
576 : name,
577 1106 : [&func](const mfem::Vector & p, mfem::real_t t) -> mfem::real_t
578 2302264 : { return func.value(t, Moose::MFEM::libMeshPointFromMFEMVector(p)); });
579 : }
580 469 : else if (std::find(VECTOR_FUNCS.begin(), VECTOR_FUNCS.end(), type) != VECTOR_FUNCS.end())
581 : {
582 339 : int dim = vectorFunctionDim(type, parameters);
583 339 : getCoefficients().declareVector<mfem::VectorFunctionCoefficient>(
584 : name,
585 : dim,
586 339 : [&func, dim](const mfem::Vector & p, mfem::real_t t, mfem::Vector & u)
587 : {
588 : libMesh::RealVectorValue vector_value =
589 1380446 : func.vectorValue(t, Moose::MFEM::libMeshPointFromMFEMVector(p));
590 5203538 : for (int i = 0; i < dim; i++)
591 : {
592 3823092 : u[i] = vector_value(i);
593 : }
594 1380446 : });
595 : }
596 130 : else if ("MFEMParsedFunction" != type)
597 : {
598 2 : mooseWarning("Could not identify whether function ",
599 : type,
600 : " is scalar or vector; no MFEM coefficient object created.");
601 : }
602 1573 : }
603 :
604 : void
605 878 : MFEMProblem::addPostprocessor(const std::string & type,
606 : const std::string & name,
607 : InputParameters & parameters)
608 : {
609 878 : if (parameters.getSystemAttributeName() == "MFEMExecutedObject")
610 : {
611 1564 : checkUserObjectNameCollision(name, "Postprocessor");
612 1564 : addObject<MFEMExecutedObject>(type, name, parameters);
613 782 : const PostprocessorValue & val = getPostprocessorValueByName(name);
614 782 : getCoefficients().declareScalar<mfem::FunctionCoefficient>(
615 784 : name, [&val](const mfem::Vector &) -> mfem::real_t { return val; });
616 : }
617 : else
618 96 : ExternalProblem::addPostprocessor(type, name, parameters);
619 878 : }
620 :
621 : void
622 392 : MFEMProblem::addVectorPostprocessor(const std::string & type,
623 : const std::string & name,
624 : InputParameters & parameters)
625 : {
626 392 : if (parameters.getSystemAttributeName() == "MFEMExecutedObject")
627 : {
628 784 : checkUserObjectNameCollision(name, "VectorPostprocessor");
629 1170 : addObject<MFEMExecutedObject>(type, name, parameters);
630 : }
631 : else
632 0 : ExternalProblem::addVectorPostprocessor(type, name, parameters);
633 386 : }
634 :
635 : InputParameters
636 105 : MFEMProblem::addMFEMFESpaceFromMOOSEVariable(InputParameters & parameters)
637 : {
638 :
639 210 : InputParameters fespace_params = _factory.getValidParams("MFEMGenericFESpace");
640 210 : InputParameters variable_params = _factory.getValidParams("MFEMVariable");
641 :
642 105 : const auto family = Utility::string_to_enum<FEFamily>(parameters.get<MooseEnum>("family"));
643 105 : auto order = static_cast<int>(parameters.get<MooseEnum>("order"));
644 105 : const auto dim = mesh().dimension();
645 :
646 105 : std::string space;
647 105 : int vdim = 1;
648 :
649 105 : switch (family)
650 : {
651 19 : case FEFamily::LAGRANGE:
652 19 : space = "H1";
653 19 : break;
654 19 : case FEFamily::NEDELEC_ONE:
655 19 : space = "ND";
656 19 : break;
657 19 : case FEFamily::RAVIART_THOMAS:
658 19 : space = "RT";
659 19 : --order;
660 19 : break;
661 18 : case FEFamily::MONOMIAL:
662 : case FEFamily::L2_LAGRANGE:
663 18 : space = "L2";
664 18 : break;
665 12 : case FEFamily::LAGRANGE_VEC:
666 12 : space = "H1";
667 12 : vdim = dim;
668 12 : break;
669 18 : case FEFamily::MONOMIAL_VEC:
670 : case FEFamily::L2_LAGRANGE_VEC:
671 18 : space = "L2";
672 18 : vdim = dim;
673 18 : break;
674 0 : default:
675 0 : mooseError("Unable to set MFEM FESpace for MOOSE variable");
676 : break;
677 : }
678 :
679 : // Create fespace name. If this already exists, we will reuse this for
680 : // the mfem variable ("gridfunction"). If using AMR, this implies all
681 : // variables sharing the fespace are affected.
682 105 : const auto fec_name = space + "_" + std::to_string(dim) + "D_P" + std::to_string(order);
683 105 : const auto fes_name = fec_name + "_X" + std::to_string(vdim);
684 :
685 : // Set all fespace parameters.
686 105 : fespace_params.set<std::string>("fec_name") = fec_name;
687 315 : fespace_params.set<int>("vdim") = vdim;
688 :
689 210 : if (!hasMFEMObject("MFEMFESpace", fes_name))
690 138 : addFESpace("MFEMGenericFESpace", fes_name, fespace_params);
691 :
692 315 : variable_params.set<MFEMFESpaceName>("fespace") = fes_name;
693 :
694 210 : return variable_params;
695 105 : }
696 :
697 : void
698 3117 : MFEMProblem::displaceMesh()
699 : {
700 : // Displace mesh
701 3117 : if (mesh().shouldDisplace())
702 : {
703 24 : mesh().displace(static_cast<mfem::GridFunction const &>(*getMeshDisplacementGridFunction()));
704 : // TODO: update FESpaces GridFunctions etc for transient solves
705 : }
706 3117 : }
707 :
708 : std::optional<std::reference_wrapper<mfem::ParGridFunction const>>
709 1641 : MFEMProblem::getMeshDisplacementGridFunction()
710 : {
711 : // If C++23 transform were available this would be easier
712 1641 : auto const displacement_variable = mesh().getMeshDisplacementVariable();
713 1641 : if (displacement_variable)
714 : {
715 48 : return *_problem_data.gridfunctions.Get(displacement_variable.value());
716 : }
717 : else
718 : {
719 1593 : return std::nullopt;
720 : }
721 : }
722 :
723 : void
724 0 : MFEMProblem::rebalanceMesh(mfem::ParMesh & pmesh)
725 : {
726 0 : if (pmesh.Nonconforming())
727 : {
728 0 : pmesh.Rebalance();
729 0 : updateFESpaces();
730 0 : updateGridFunctions();
731 : }
732 0 : }
733 :
734 : void
735 13 : MFEMProblem::updateFESpaces()
736 : {
737 39 : for (const auto & fe_space_pair : _problem_data.fespaces)
738 26 : fe_space_pair.second->Update();
739 13 : }
740 :
741 : void
742 26 : MFEMProblem::updateGridFunctions()
743 : {
744 78 : for (const auto & gridfunction_pair : _problem_data.gridfunctions)
745 52 : gridfunction_pair.second->Update();
746 26 : }
747 :
748 : std::vector<VariableName>
749 0 : MFEMProblem::getAuxVariableNames()
750 : {
751 0 : return systemBaseAuxiliary().getVariableNames();
752 : }
753 :
754 : MFEMMesh &
755 71549 : MFEMProblem::mesh()
756 : {
757 : mooseAssert(ExternalProblem::mesh().type() == "MFEMMesh",
758 : "Please choose the MFEMMesh mesh type for an MFEMProblem\n");
759 71549 : return static_cast<MFEMMesh &>(_mesh);
760 : }
761 :
762 : const MFEMMesh &
763 3535 : MFEMProblem::mesh() const
764 : {
765 3535 : return const_cast<MFEMProblem *>(this)->mesh();
766 : }
767 :
768 : void
769 187 : MFEMProblem::addSubMesh(const std::string & var_type,
770 : const std::string & var_name,
771 : InputParameters & parameters)
772 : {
773 374 : auto & mfem_submesh = *addObject<MFEMSubMesh>(var_type, var_name, parameters).front();
774 : // Register submesh.
775 187 : getProblemData().submeshes.Register(var_name, mfem_submesh.getSubMesh());
776 187 : }
777 :
778 : void
779 34 : MFEMProblem::addQuadratureFunction(const std::string & type,
780 : const std::string & name,
781 : InputParameters & parameters)
782 : {
783 : // The object declares its coefficient with the CoefficientManager on construction.
784 68 : addObject<MFEMObject>(type, name, parameters);
785 34 : }
786 :
787 : void
788 664 : MFEMProblem::addTransfer(const std::string & transfer_name,
789 : const std::string & name,
790 : InputParameters & parameters)
791 : {
792 664 : if (parameters.getBase() == "MFEMSubMeshTransfer")
793 687 : addObject<MFEMExecutedObject>(transfer_name, name, parameters);
794 : else
795 435 : ExternalProblem::addTransfer(transfer_name, name, parameters);
796 664 : }
797 :
798 : void
799 1148 : MFEMProblem::addInitialCondition(const std::string & ic_name,
800 : const std::string & name,
801 : InputParameters & parameters)
802 : {
803 2296 : addObject<MFEMExecutedObject>(ic_name, name, parameters);
804 1148 : }
805 :
806 : void
807 11257 : MFEMProblem::executeMFEMObjects(const ExecFlagType & exec_type)
808 : {
809 11257 : std::vector<MFEMExecutedObject *> objects;
810 11257 : theWarehouse()
811 11257 : .query()
812 11257 : .condition<AttribSystem>("MFEMExecutedObject")
813 11257 : .condition<AttribExecOns>(exec_type)
814 22514 : .condition<AttribThread>(0)
815 11257 : .queryInto(objects);
816 :
817 11257 : std::map<std::string, const MFEMExecutedObject *> suppliers;
818 15596 : for (auto * const object : objects)
819 8680 : for (const auto & item : object->getSuppliedItems())
820 : {
821 4341 : const auto [it, inserted] = suppliers.emplace(item, object);
822 4341 : if (!inserted && it->second != object)
823 2 : mooseError("MFEM executed-object dependency ambiguity on ",
824 : exec_type,
825 : ": both '",
826 2 : it->second->name(),
827 : "' and '",
828 2 : object->name(),
829 : "' supply '",
830 : item,
831 : "'.");
832 : }
833 :
834 15592 : for (auto * const object : objects)
835 : {
836 4337 : object->initialize();
837 4337 : object->execute();
838 4337 : object->finalize();
839 :
840 4337 : if (auto * const pp = dynamic_cast<const Postprocessor *>(object))
841 : {
842 1304 : _reporter_data.finalize(pp->PPName());
843 1304 : setPostprocessorValueByName(pp->PPName(), pp->getValue());
844 : }
845 :
846 4337 : if (auto * const vpp = dynamic_cast<VectorPostprocessor *>(object))
847 998 : _reporter_data.finalize(vpp->PPName());
848 : }
849 11259 : }
850 :
851 : std::string
852 1599 : MFEMProblem::solverTypeString(const unsigned int libmesh_dbg_var(solver_sys_num))
853 : {
854 : mooseAssert(solver_sys_num == 0, "No support for multi-system with MFEM right now");
855 :
856 1599 : std::vector<std::string> solvers;
857 :
858 1599 : if (getProblemData().nonlinear_solver)
859 50 : solvers.push_back(MooseUtils::prettyCppType(getProblemData().nonlinear_solver.get()));
860 :
861 1599 : if (getProblemData().jacobian_solver)
862 : {
863 1065 : solvers.push_back(MooseUtils::prettyCppType(getProblemData().jacobian_solver.get()));
864 1065 : if (const auto * prec = getProblemData().jacobian_solver->GetPreconditioner())
865 954 : solvers.push_back(MooseUtils::prettyCppType(prec));
866 : }
867 :
868 5869 : return solvers.empty() ? "None" : MooseUtils::stringJoin(solvers);
869 1599 : }
870 :
871 : bool
872 105 : MFEMProblem::hasMFEMObject(const std::string & system, const std::string & name) const
873 : {
874 105 : std::vector<MooseObject *> objs;
875 105 : theWarehouse()
876 105 : .query()
877 105 : .condition<AttribSystem>(system)
878 210 : .condition<AttribThread>(0)
879 105 : .condition<AttribName>(name)
880 105 : .queryInto(objs);
881 210 : return !objs.empty();
882 105 : }
883 :
884 : #endif
|
