https://mooseframework.inl.gov
EquationSystem.C
Go to the documentation of this file.
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 "EquationSystem.h"
13 #include "libmesh/int_range.h"
14 
15 namespace Moose::MFEM
16 {
17 
19 
20 void
22 {
23  for (const auto i : make_range(_h_blocks.NumRows()))
24  for (const auto j : make_range(_h_blocks.NumCols()))
25  delete _h_blocks(i, j);
26  _h_blocks.DeleteAll();
27 }
28 
29 bool
30 EquationSystem::VectorContainsName(const std::vector<std::string> & the_vector,
31  const std::string & name) const
32 {
33 
34  auto iter = std::find(the_vector.begin(), the_vector.end(), name);
35 
36  return (iter != the_vector.end());
37 }
38 
39 void
40 EquationSystem::AddCoupledVariableNameIfMissing(const std::string & coupled_var_name)
41 {
42  if (!VectorContainsName(_coupled_var_names, coupled_var_name))
43  _coupled_var_names.push_back(coupled_var_name);
44 }
45 
46 void
47 EquationSystem::AddTestVariableNameIfMissing(const std::string & test_var_name)
48 {
49  if (!VectorContainsName(_test_var_names, test_var_name))
50  _test_var_names.push_back(test_var_name);
51 }
52 
53 void
55 {
56  // If a coupled variable has an equation associated with it,
57  // add it to the set of trial variables.
58  for (const auto & coupled_var_name : _coupled_var_names)
59  {
60  if (VectorContainsName(_test_var_names, coupled_var_name))
61  _trial_var_names.push_back(coupled_var_name);
62  else
63  _eliminated_var_names.push_back(coupled_var_name);
64  }
65 }
66 
67 void
68 EquationSystem::AddKernel(std::shared_ptr<MFEMKernel> kernel)
69 {
70  AddTestVariableNameIfMissing(kernel->getTestVariableName());
71  AddCoupledVariableNameIfMissing(kernel->getTrialVariableName());
72  auto trial_var_name = kernel->getTrialVariableName();
73  auto test_var_name = kernel->getTestVariableName();
74  if (!_kernels_map.Has(test_var_name))
75  {
76  auto kernel_field_map =
77  std::make_shared<Moose::MFEM::NamedFieldsMap<std::vector<std::shared_ptr<MFEMKernel>>>>();
78  _kernels_map.Register(test_var_name, std::move(kernel_field_map));
79  }
80  // Register new kernels map if not present for the test/trial variable
81  // pair
82  if (!_kernels_map.Get(test_var_name)->Has(trial_var_name))
83  {
84  auto kernels = std::make_shared<std::vector<std::shared_ptr<MFEMKernel>>>();
85  _kernels_map.Get(test_var_name)->Register(trial_var_name, std::move(kernels));
86  }
87  _kernels_map.GetRef(test_var_name).Get(trial_var_name)->push_back(std::move(kernel));
88 }
89 
90 void
91 EquationSystem::AddIntegratedBC(std::shared_ptr<MFEMIntegratedBC> bc)
92 {
93  AddTestVariableNameIfMissing(bc->getTestVariableName());
94  AddCoupledVariableNameIfMissing(bc->getTrialVariableName());
95  auto trial_var_name = bc->getTrialVariableName();
96  auto test_var_name = bc->getTestVariableName();
97  if (!_integrated_bc_map.Has(test_var_name))
98  {
99  auto integrated_bc_field_map = std::make_shared<
101  _integrated_bc_map.Register(test_var_name, std::move(integrated_bc_field_map));
102  }
103  // Register new integrated bc map if not present for the test/trial variable
104  // pair
105  if (!_integrated_bc_map.Get(test_var_name)->Has(trial_var_name))
106  {
107  auto bcs = std::make_shared<std::vector<std::shared_ptr<MFEMIntegratedBC>>>();
108  _integrated_bc_map.Get(test_var_name)->Register(trial_var_name, std::move(bcs));
109  }
110  _integrated_bc_map.GetRef(test_var_name).Get(trial_var_name)->push_back(std::move(bc));
111 }
112 
113 void
114 EquationSystem::AddEssentialBC(std::shared_ptr<MFEMEssentialBC> bc)
115 {
116  AddTestVariableNameIfMissing(bc->getTestVariableName());
117  auto test_var_name = bc->getTestVariableName();
118  if (!_essential_bc_map.Has(test_var_name))
119  {
120  auto bcs = std::make_shared<std::vector<std::shared_ptr<MFEMEssentialBC>>>();
121  _essential_bc_map.Register(test_var_name, std::move(bcs));
122  }
123  _essential_bc_map.GetRef(test_var_name).push_back(std::move(bc));
124 }
125 
126 void
128  const Moose::MFEM::FESpaces & /*fespaces*/,
129  mfem::AssemblyLevel assembly_level)
130 {
131  _assembly_level = assembly_level;
132 
133  for (auto & test_var_name : _test_var_names)
134  {
135  if (!gridfunctions.Has(test_var_name))
136  {
137  mooseError("MFEM variable ",
138  test_var_name,
139  " requested by equation system during initialisation was "
140  "not found in gridfunctions");
141  }
142  // Store pointers to test FESpaces
143  _test_pfespaces.push_back(gridfunctions.Get(test_var_name)->ParFESpace());
144  // Create auxiliary gridfunctions for storing essential constraints from Dirichlet conditions
145  _var_ess_constraints.emplace_back(
146  std::make_unique<mfem::ParGridFunction>(gridfunctions.Get(test_var_name)->ParFESpace()));
147  }
148 
149  // Store pointers to FESpaces of all coupled variables
150  for (auto & coupled_var_name : _coupled_var_names)
151  _coupled_pfespaces.push_back(gridfunctions.Get(coupled_var_name)->ParFESpace());
152 
153  // Extract which coupled variables are to be trivially eliminated and which are trial variables
155 
156  // Store pointers to coupled variable GridFunctions that are to be eliminated prior to forming the
157  // jacobian
158  for (auto & eliminated_var_name : _eliminated_var_names)
159  _eliminated_variables.Register(eliminated_var_name,
160  gridfunctions.GetShared(eliminated_var_name));
161 }
162 
163 void
165 {
166  _ess_tdof_lists.resize(_test_var_names.size());
167  for (const auto i : index_range(_test_var_names))
168  {
169  auto test_var_name = _test_var_names.at(i);
170  if (!_essential_bc_map.Has(test_var_name))
171  continue;
172 
173  // Set default value of gridfunction used in essential BC. Values
174  // overwritten in applyEssentialBCs
175  mfem::ParGridFunction & trial_gf(*(_var_ess_constraints.at(i)));
176  auto * const pmesh = _test_pfespaces.at(i)->GetParMesh();
177  mooseAssert(pmesh, "parallel mesh is null");
178 
179  auto bcs = _essential_bc_map.GetRef(test_var_name);
180  mfem::Array<int> global_ess_markers(pmesh->bdr_attributes.Max());
181  global_ess_markers = 0;
182  for (auto & bc : bcs)
183  {
184  bc->ApplyBC(trial_gf);
185 
186  mfem::Array<int> ess_bdrs(bc->getBoundaryMarkers());
187  for (auto it = 0; it != pmesh->bdr_attributes.Max(); ++it)
188  {
189  global_ess_markers[it] = std::max(global_ess_markers[it], ess_bdrs[it]);
190  }
191  }
192  trial_gf.FESpace()->GetEssentialTrueDofs(global_ess_markers, _ess_tdof_lists.at(i));
193  }
194 }
195 
196 void
198 {
199  for (const auto & test_var_name : _test_var_names)
200  {
201  auto lf = _lfs.Get(test_var_name);
202  for (const auto & eliminated_var_name : _eliminated_var_names)
203  {
204  if (_mblfs.Has(test_var_name) && _mblfs.Get(test_var_name)->Has(eliminated_var_name))
205  {
206  auto mblf = _mblfs.Get(test_var_name)->Get(eliminated_var_name);
207  mblf->AddMult(*_eliminated_variables.Get(eliminated_var_name), *lf, -1.0);
208  }
209  }
210  }
211 }
212 
213 void
214 EquationSystem::FormLinearSystem(mfem::OperatorHandle & op,
215  mfem::BlockVector & trueX,
216  mfem::BlockVector & trueRHS)
217 {
218 
219  switch (_assembly_level)
220  {
221  case mfem::AssemblyLevel::LEGACY:
222  FormLegacySystem(op, trueX, trueRHS);
223  break;
224  default:
225  mooseAssert(_test_var_names.size() == 1,
226  "Non-legacy assembly is only supported for single-variable systems");
227  mooseAssert(
228  _test_var_names.size() == _trial_var_names.size(),
229  "Non-legacy assembly is only supported for single test and trial variable systems");
230  FormSystem(op, trueX, trueRHS);
231  }
232 }
233 
234 void
235 EquationSystem::FormSystem(mfem::OperatorHandle & op,
236  mfem::BlockVector & trueX,
237  mfem::BlockVector & trueRHS)
238 {
239  auto & test_var_name = _test_var_names.at(0);
240  auto blf = _blfs.Get(test_var_name);
241  auto lf = _lfs.Get(test_var_name);
242  mfem::BlockVector aux_x, aux_rhs;
243  mfem::OperatorPtr aux_a;
244 
245  blf->FormLinearSystem(
246  _ess_tdof_lists.at(0), *(_var_ess_constraints.at(0)), *lf, aux_a, aux_x, aux_rhs);
247 
248  trueX.GetBlock(0) = aux_x;
249  trueRHS.GetBlock(0) = aux_rhs;
250  trueX.SyncFromBlocks();
251  trueRHS.SyncFromBlocks();
252 
253  op.Reset(aux_a.Ptr());
254  aux_a.SetOperatorOwner(false);
255 }
256 
257 void
258 EquationSystem::FormLegacySystem(mfem::OperatorHandle & op,
259  mfem::BlockVector & trueX,
260  mfem::BlockVector & trueRHS)
261 {
262 
263  // Allocate block operator
264  DeleteAllBlocks();
265  _h_blocks.SetSize(_test_var_names.size(), _test_var_names.size());
266  // Form diagonal blocks.
267  for (const auto i : index_range(_test_var_names))
268  {
269  auto & test_var_name = _test_var_names.at(i);
270  auto blf = _blfs.Get(test_var_name);
271  auto lf = _lfs.Get(test_var_name);
272  mfem::Vector aux_x, aux_rhs;
273  mfem::HypreParMatrix * aux_a = new mfem::HypreParMatrix;
274  blf->FormLinearSystem(
275  _ess_tdof_lists.at(i), *(_var_ess_constraints.at(i)), *lf, *aux_a, aux_x, aux_rhs);
276  _h_blocks(i, i) = aux_a;
277  trueX.GetBlock(i) = aux_x;
278  trueRHS.GetBlock(i) = aux_rhs;
279  }
280 
281  // Form off-diagonal blocks
282  for (const auto i : index_range(_test_var_names))
283  {
284  auto test_var_name = _test_var_names.at(i);
285  for (const auto j : index_range(_trial_var_names))
286  {
287  auto trial_var_name = _trial_var_names.at(j);
288 
289  mfem::Vector aux_x, aux_rhs;
290  mfem::ParLinearForm aux_lf(_test_pfespaces.at(i));
291  aux_lf = 0.0;
292  if (_mblfs.Has(test_var_name) && _mblfs.Get(test_var_name)->Has(trial_var_name))
293  {
294  auto mblf = _mblfs.Get(test_var_name)->Get(trial_var_name);
295  mfem::HypreParMatrix * aux_a = new mfem::HypreParMatrix;
296  mblf->FormRectangularLinearSystem(_ess_tdof_lists.at(j),
297  _ess_tdof_lists.at(i),
298  *(_var_ess_constraints.at(j)),
299  aux_lf,
300  *aux_a,
301  aux_x,
302  aux_rhs);
303  _h_blocks(i, j) = aux_a;
304  trueRHS.GetBlock(i) += aux_rhs;
305  }
306  }
307  }
308  // Sync memory
309  for (const auto i : index_range(_test_var_names))
310  {
311  trueX.GetBlock(i).SyncAliasMemory(trueX);
312  trueRHS.GetBlock(i).SyncAliasMemory(trueRHS);
313  }
314 
315  // Create monolithic matrix
316  op.Reset(mfem::HypreParMatrixFromBlocks(_h_blocks));
317 }
318 
319 void
320 EquationSystem::BuildJacobian(mfem::BlockVector & trueX, mfem::BlockVector & trueRHS)
321 {
322  height = trueX.Size();
323  width = trueRHS.Size();
324  FormLinearSystem(_jacobian, trueX, trueRHS);
325 }
326 
327 void
328 EquationSystem::Mult(const mfem::Vector & x, mfem::Vector & residual) const
329 {
330  _jacobian->Mult(x, residual);
331  x.HostRead();
332  residual.HostRead();
333 }
334 
335 mfem::Operator &
336 EquationSystem::GetGradient(const mfem::Vector &) const
337 {
338  return *_jacobian;
339 }
340 
341 void
342 EquationSystem::RecoverFEMSolution(mfem::BlockVector & trueX,
343  Moose::MFEM::GridFunctions & gridfunctions)
344 {
345  for (const auto i : index_range(_trial_var_names))
346  {
347  auto & trial_var_name = _trial_var_names.at(i);
348  trueX.GetBlock(i).SyncAliasMemory(trueX);
349  gridfunctions.Get(trial_var_name)->Distribute(&(trueX.GetBlock(i)));
350  }
351 }
352 
353 void
355 {
356  // Register linear forms
357  for (const auto i : index_range(_test_var_names))
358  {
359  auto test_var_name = _test_var_names.at(i);
360  _lfs.Register(test_var_name, std::make_shared<mfem::ParLinearForm>(_test_pfespaces.at(i)));
361  _lfs.GetRef(test_var_name) = 0.0;
362  }
363 
364  // Apply boundary conditions
366 
367  for (auto & test_var_name : _test_var_names)
368  {
369  // Apply kernels
370  auto lf = _lfs.GetShared(test_var_name);
371  ApplyDomainLFIntegrators(test_var_name, lf, _kernels_map);
373  lf->Assemble();
374  }
375 
376  // Eliminate trivially eliminated variables by subtracting contributions from linear forms
378 }
379 
380 void
382 {
383  // Register bilinear forms
384  for (const auto i : index_range(_test_var_names))
385  {
386  auto test_var_name = _test_var_names.at(i);
387  _blfs.Register(test_var_name, std::make_shared<mfem::ParBilinearForm>(_test_pfespaces.at(i)));
388 
389  // Apply kernels
390  auto blf = _blfs.GetShared(test_var_name);
391  blf->SetAssemblyLevel(_assembly_level);
392  ApplyBoundaryBLFIntegrators<mfem::ParBilinearForm>(
393  test_var_name, test_var_name, blf, _integrated_bc_map);
394  ApplyDomainBLFIntegrators<mfem::ParBilinearForm>(
395  test_var_name, test_var_name, blf, _kernels_map);
396  // Assemble
397  blf->Assemble();
398  }
399 }
400 
401 void
403 {
404  // Register mixed bilinear forms. Note that not all combinations may
405  // have a kernel
406 
407  // Create mblf for each test/coupled variable pair with an added kernel
408  for (const auto i : index_range(_test_var_names))
409  {
410  auto test_var_name = _test_var_names.at(i);
411  auto test_mblfs = std::make_shared<Moose::MFEM::NamedFieldsMap<mfem::ParMixedBilinearForm>>();
412  for (const auto j : index_range(_coupled_var_names))
413  {
414  const auto & coupled_var_name = _coupled_var_names.at(j);
415  auto mblf = std::make_shared<mfem::ParMixedBilinearForm>(_coupled_pfespaces.at(j),
416  _test_pfespaces.at(i));
417  // Register MixedBilinearForm if kernels exist for it, and assemble
418  // kernels
419  if (_kernels_map.Has(test_var_name) &&
420  _kernels_map.Get(test_var_name)->Has(coupled_var_name) &&
421  test_var_name != coupled_var_name)
422  {
423  mblf->SetAssemblyLevel(_assembly_level);
424  // Apply all mixed kernels with this test/trial pair
425  ApplyDomainBLFIntegrators<mfem::ParMixedBilinearForm>(
426  coupled_var_name, test_var_name, mblf, _kernels_map);
427  // Assemble mixed bilinear forms
428  mblf->Assemble();
429  // Register mixed bilinear forms associated with a single trial variable
430  // for the current test variable
431  test_mblfs->Register(coupled_var_name, mblf);
432  }
433  }
434  // Register all mixed bilinear form sets associated with a single test
435  // variable
436  _mblfs.Register(test_var_name, test_mblfs);
437  }
438 }
439 
440 void
442 {
446 }
447 
449 
450 void
452 {
459 }
460 
461 void
463 {
464  if (fabs(dt - _dt_coef.constant) > 1.0e-12 * dt)
465  {
466  _dt_coef.constant = dt;
467  for (auto test_var_name : _test_var_names)
468  {
469  auto blf = _blfs.Get(test_var_name);
470  blf->Update();
471  blf->Assemble();
472  }
473  }
474 }
475 
476 void
478 {
479  // If a coupled variable has an equation associated with it,
480  // add it to the set of trial variables.
481  for (const auto & coupled_var_name : _coupled_var_names)
482  {
483  for (const auto & test_var_name : _test_var_names)
484  {
485  const auto time_derivative_test_var_name = GetTimeDerivativeName(test_var_name);
486  if (time_derivative_test_var_name == coupled_var_name)
487  {
488  if (!VectorContainsName(_trial_var_names, coupled_var_name))
489  _trial_var_names.push_back(coupled_var_name);
490  }
491  else
492  {
493  if (!VectorContainsName(_eliminated_var_names, coupled_var_name))
494  _eliminated_var_names.push_back(coupled_var_name);
495  }
496  }
497  }
498 }
499 
500 void
501 TimeDependentEquationSystem::AddKernel(std::shared_ptr<MFEMKernel> kernel)
502 {
503  if (kernel->getTrialVariableName() == GetTimeDerivativeName(kernel->getTestVariableName()))
504  {
505  auto trial_var_name = kernel->getTrialVariableName();
506  auto test_var_name = kernel->getTestVariableName();
507  AddTestVariableNameIfMissing(test_var_name);
508  AddCoupledVariableNameIfMissing(test_var_name);
509 
510  if (!_td_kernels_map.Has(test_var_name))
511  {
512  auto kernel_field_map =
513  std::make_shared<Moose::MFEM::NamedFieldsMap<std::vector<std::shared_ptr<MFEMKernel>>>>();
514  _td_kernels_map.Register(test_var_name, std::move(kernel_field_map));
515  }
516  // Register new kernels map if not present for the test variable
517  if (!_td_kernels_map.Get(test_var_name)->Has(test_var_name))
518  {
519  auto kernels = std::make_shared<std::vector<std::shared_ptr<MFEMKernel>>>();
520  _td_kernels_map.Get(test_var_name)->Register(test_var_name, std::move(kernels));
521  }
522  _td_kernels_map.GetRef(test_var_name).Get(test_var_name)->push_back(std::move(kernel));
523  }
524  else
525  {
527  }
528 }
529 
530 void
532 {
534 
535  // Build and assemble bilinear forms acting on time derivatives
536  for (const auto i : index_range(_test_var_names))
537  {
538  auto test_var_name = _test_var_names.at(i);
539 
540  _td_blfs.Register(test_var_name,
541  std::make_shared<mfem::ParBilinearForm>(_test_pfespaces.at(i)));
542 
543  // Apply kernels to td_blf
544  auto td_blf = _td_blfs.GetShared(test_var_name);
545  td_blf->SetAssemblyLevel(_assembly_level);
546  ApplyBoundaryBLFIntegrators<mfem::ParBilinearForm>(
547  test_var_name, test_var_name, td_blf, _integrated_bc_map);
548  ApplyDomainBLFIntegrators<mfem::ParBilinearForm>(
549  test_var_name, test_var_name, td_blf, _td_kernels_map);
550 
551  // Recover and scale integrators from blf. This is to apply the dt*du/dt contributions from the
552  // operator on the trial variable in the implicit integration scheme
553  auto blf = _blfs.Get(test_var_name);
554  auto integs = blf->GetDBFI();
555  auto b_integs = blf->GetBBFI();
556  auto markers = blf->GetBBFI_Marker();
557 
558  // If implicit contributions exist, scale them by timestep and add to td_blf
559  if (integs->Size() || b_integs->Size())
560  {
561  mfem::SumIntegrator * sum = new mfem::SumIntegrator(false);
562  ScaleIntegrator * scaled_sum = new ScaleIntegrator(sum, _dt_coef.constant, true);
563 
564  for (int i = 0; i < integs->Size(); ++i)
565  sum->AddIntegrator(*integs[i]);
566 
567  for (int i = 0; i < b_integs->Size(); ++i)
568  td_blf->AddBoundaryIntegrator(new ScaleIntegrator(*b_integs[i], _dt_coef.constant, false),
569  *(*markers[i]));
570 
571  // scaled_sum is owned by td_blf
572  td_blf->AddDomainIntegrator(scaled_sum);
573  }
574 
575  // Assemble form
576  td_blf->Assemble();
577  }
578 }
579 
580 void
582  mfem::BlockVector & truedXdt,
583  mfem::BlockVector & trueRHS)
584 {
585 
586  // Allocate block operator
587  DeleteAllBlocks();
588  _h_blocks.SetSize(_test_var_names.size(), _test_var_names.size());
589  // Form diagonal blocks.
590  for (const auto i : index_range(_test_var_names))
591  {
592  auto & test_var_name = _test_var_names.at(i);
593  auto td_blf = _td_blfs.Get(test_var_name);
594  auto blf = _blfs.Get(test_var_name);
595  auto lf = _lfs.Get(test_var_name);
596  // if implicit, add contribution to linear form from terms involving state
597  // variable at previous timestep: {
598  blf->AddMult(*_eliminated_variables.Get(test_var_name), *lf, -1.0);
599  // }
600  mfem::Vector aux_x, aux_rhs;
601  // Update solution values on Dirichlet values to be in terms of du/dt instead of u
602  mfem::Vector bc_x = *(_var_ess_constraints.at(i).get());
603  bc_x -= *_eliminated_variables.Get(test_var_name);
604  bc_x /= _dt_coef.constant;
605 
606  // Form linear system for operator acting on vector of du/dt
607  mfem::HypreParMatrix * aux_a = new mfem::HypreParMatrix;
608  td_blf->FormLinearSystem(_ess_tdof_lists.at(i), bc_x, *lf, *aux_a, aux_x, aux_rhs);
609  _h_blocks(i, i) = aux_a;
610  truedXdt.GetBlock(i) = aux_x;
611  trueRHS.GetBlock(i) = aux_rhs;
612  }
613 
614  truedXdt.SyncFromBlocks();
615  trueRHS.SyncFromBlocks();
616 
617  // Create monolithic matrix
618  op.Reset(mfem::HypreParMatrixFromBlocks(_h_blocks));
619 }
620 
621 void
622 TimeDependentEquationSystem::FormSystem(mfem::OperatorHandle & op,
623  mfem::BlockVector & truedXdt,
624  mfem::BlockVector & trueRHS)
625 {
626  auto & test_var_name = _test_var_names.at(0);
627  auto td_blf = _td_blfs.Get(test_var_name);
628  auto blf = _blfs.Get(test_var_name);
629  auto lf = _lfs.Get(test_var_name);
630  // if implicit, add contribution to linear form from terms involving state
631  // variable at previous timestep: {
632 
633  // The AddMult method in mfem::BilinearForm is not defined for non-legacy assembly
634  mfem::Vector lf_prev(lf->Size());
635  blf->Mult(*_eliminated_variables.Get(test_var_name), lf_prev);
636  *lf -= lf_prev;
637  // }
638  mfem::Vector aux_x, aux_rhs;
639  // Update solution values on Dirichlet values to be in terms of du/dt instead of u
640  mfem::Vector bc_x = *(_var_ess_constraints.at(0).get());
641  bc_x -= *_eliminated_variables.Get(test_var_name);
642  bc_x /= _dt_coef.constant;
643 
644  // Form linear system for operator acting on vector of du/dt
645  mfem::OperatorPtr aux_a;
646  td_blf->FormLinearSystem(_ess_tdof_lists.at(0), bc_x, *lf, aux_a, aux_x, aux_rhs);
647 
648  truedXdt.GetBlock(0) = aux_x;
649  trueRHS.GetBlock(0) = aux_rhs;
650  truedXdt.SyncFromBlocks();
651  trueRHS.SyncFromBlocks();
652 
653  // Create monolithic matrix
654  op.Reset(aux_a.Ptr());
655  aux_a.SetOperatorOwner(false);
656 }
657 
658 void
660 {
662 }
663 
664 } // namespace Moose::MFEM
665 
666 #endif
std::string name(const ElemQuality q)
mfem::ConstantCoefficient _dt_coef
Coefficient for timestep scaling.
virtual void EliminateCoupledVariables()
virtual void AddKernel(std::shared_ptr< MFEMKernel > kernel) override
Add kernels.
virtual void AddTestVariableNameIfMissing(const std::string &test_var_name)
Add test variable to EquationSystem.
virtual void AddKernel(std::shared_ptr< MFEMKernel > kernel)
Add kernels.
virtual void RecoverFEMSolution(mfem::BlockVector &trueX, Moose::MFEM::GridFunctions &gridfunctions)
Update variable from solution vector after solve.
virtual void BuildBilinearForms()
Build bilinear forms.
bool Has(const std::string &field_name) const
Predicate to check if a field is registered with name field_name.
std::vector< mfem::ParFiniteElementSpace * > _coupled_pfespaces
Pointers to finite element spaces associated with coupled variables.
Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMKernel > > > > _td_kernels_map
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:333
virtual void FormSystem(mfem::OperatorHandle &op, mfem::BlockVector &truedXdt, mfem::BlockVector &trueRHS) override
void ApplyBoundaryLFIntegrators(const std::string &test_var_name, std::shared_ptr< mfem::ParLinearForm > form, Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMIntegratedBC >>>> &integrated_bc_map)
bool VectorContainsName(const std::vector< std::string > &the_vector, const std::string &name) const
std::vector< std::string > _eliminated_var_names
Names of all coupled variables without a corresponding test variable.
Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< mfem::ParMixedBilinearForm > > _mblfs
Moose::MFEM::GridFunctions _eliminated_variables
Pointers to coupled variables not part of the reduced EquationSystem.
virtual void AddEssentialBC(std::shared_ptr< MFEMEssentialBC > bc)
void AddCoupledVariableNameIfMissing(const std::string &coupled_var_name) override
Add coupled variable to EquationSystem.
std::vector< mfem::Array< int > > _ess_tdof_lists
Lightweight adaptor over an std::map from strings to pointer to T.
mfem::AssemblyLevel _assembly_level
std::string GetTimeDerivativeName(std::string name)
auto max(const L &left, const R &right)
Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMEssentialBC > > > _essential_bc_map
Arrays to store essential BCs to act on each component of weak form.
std::vector< mfem::ParFiniteElementSpace * > _test_pfespaces
Pointers to finite element spaces associated with test variables.
void Mult(const mfem::Vector &u, mfem::Vector &residual) const override
Compute residual y = Mu.
std::vector< std::string > _trial_var_names
Subset of _coupled_var_names of all variables corresponding to gridfunctions with degrees of freedom ...
T * Get(const std::string &field_name) const
Returns a non-owning pointer to the field. This is guaranteed to return a non-null pointer...
virtual void FormLegacySystem(mfem::OperatorHandle &op, mfem::BlockVector &trueX, mfem::BlockVector &trueRHS)
std::shared_ptr< T > GetShared(const std::string &field_name) const
Returns a shared pointer to the field. This is guaranteed to return a non-null shared pointer...
virtual void SetTrialVariableNames() override
Set trial variable names from subset of coupled variables that have an associated test variable...
mfem::OperatorHandle _jacobian
std::vector< std::string > _test_var_names
Names of all test variables corresponding to linear forms in this equation system.
virtual void BuildEquationSystem()
std::vector< std::unique_ptr< mfem::ParGridFunction > > _var_ess_constraints
Gridfunctions holding essential constraints from Dirichlet BCs.
virtual void Init(Moose::MFEM::GridFunctions &gridfunctions, const Moose::MFEM::FESpaces &fespaces, mfem::AssemblyLevel assembly_level)
Initialise.
virtual void FormLegacySystem(mfem::OperatorHandle &op, mfem::BlockVector &truedXdt, mfem::BlockVector &trueRHS) override
virtual void BuildMixedBilinearForms()
Integrator which scales its results by a constant value.
mfem::Array2D< const mfem::HypreParMatrix * > _h_blocks
virtual void SetTrialVariableNames()
Set trial variable names from subset of coupled variables that have an associated test variable...
virtual void AddCoupledVariableNameIfMissing(const std::string &coupled_var_name)
Add coupled variable to EquationSystem.
virtual void FormLinearSystem(mfem::OperatorHandle &op, mfem::BlockVector &trueX, mfem::BlockVector &trueRHS)
Form linear system, with essential boundary conditions accounted for.
std::vector< std::string > _coupled_var_names
Names of all trial variables of kernels and boundary conditions added to this EquationSystem.
void Register(const std::string &field_name, FieldArgs &&... args)
Construct new field with name field_name and register.
Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMIntegratedBC > > > > _integrated_bc_map
Arrays to store integrated BCs to act on each component of weak form.
void DeleteAllBlocks()
Deletes the HypreParMatrix associated with any pointer stored in _h_blocks, and then proceeds to dele...
Moose::MFEM::NamedFieldsMap< mfem::ParBilinearForm > _td_blfs
Container to store contributions to weak form of the form (F du/dt, v)
virtual void BuildLinearForms()
Build linear forms and eliminate constrained DoFs.
IntRange< T > make_range(T beg, T end)
virtual void BuildJacobian(mfem::BlockVector &trueX, mfem::BlockVector &trueRHS)
Build linear system, with essential boundary conditions accounted for.
void ApplyDomainLFIntegrators(const std::string &test_var_name, std::shared_ptr< mfem::ParLinearForm > form, Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMKernel >>>> &kernels_map)
virtual void ApplyEssentialBCs()
Moose::MFEM::NamedFieldsMap< mfem::ParLinearForm > _lfs
virtual void BuildBilinearForms() override
Build bilinear forms.
Moose::MFEM::NamedFieldsMap< Moose::MFEM::NamedFieldsMap< std::vector< std::shared_ptr< MFEMKernel > > > > _kernels_map
Arrays to store kernels to act on each component of weak form.
T & GetRef(const std::string &field_name) const
Returns a reference to a field.
mfem::Operator & GetGradient(const mfem::Vector &u) const override
Compute J = M + grad_H(u)
virtual void FormSystem(mfem::OperatorHandle &op, mfem::BlockVector &trueX, mfem::BlockVector &trueRHS)
Moose::MFEM::NamedFieldsMap< mfem::ParBilinearForm > _blfs
virtual void AddIntegratedBC(std::shared_ptr< MFEMIntegratedBC > kernel)
auto index_range(const T &sizable)