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