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 "EigenproblemEquationSystem.h"
      13             : #include "MFEMEigensolverBase.h"
      14             : #include "libmesh/int_range.h"
      15             : 
      16             : namespace Moose::MFEM
      17             : {
      18             : 
      19             : void
      20          52 : EigenproblemEquationSystem::ApplyEssentialBCs()
      21             : {
      22          52 :   _ess_tdof_lists.resize(1);
      23          52 :   mfem::ParGridFunction & trial_gf = *(_var_ess_constraints.at(0));
      24          52 :   _global_ess_markers.SetSize(trial_gf.ParFESpace()->GetParMesh()->bdr_attributes.Max());
      25          52 :   _global_ess_markers = 0;
      26          52 :   trial_gf.Update();
      27          52 :   trial_gf = _gfuncs->GetRef(_trial_var_names.at(0));
      28          52 :   trial_gf.ParFESpace()->GetParMesh()->MarkExternalBoundaries(_global_ess_markers);
      29          52 :   trial_gf.ParFESpace()->GetEssentialTrueDofs(_global_ess_markers, _ess_tdof_lists.at(0));
      30          52 : }
      31             : 
      32             : void
      33          26 : EigenproblemEquationSystem::FormEigenproblemMatrix(mfem::OperatorHandle & op)
      34             : {
      35          26 :   auto & test_var_name = _test_var_names.at(0);
      36          26 :   auto blf = _blfs.Get(test_var_name);
      37             : 
      38          26 :   blf->EliminateEssentialBCDiag(_global_ess_markers, 1.0);
      39          26 :   blf->Finalize();
      40          26 :   op.Reset(blf->ParallelAssemble());
      41          26 : }
      42             : 
      43             : void
      44          26 : EigenproblemEquationSystem::FormMassMatrix(mfem::OperatorHandle & op)
      45             : {
      46          26 :   mfem::ConstantCoefficient one(1.0);
      47          26 :   mfem::ParFiniteElementSpace * fespace = _test_pfespaces.at(0);
      48          26 :   std::unique_ptr<mfem::ParBilinearForm> m = std::make_unique<mfem::ParBilinearForm>(fespace);
      49             : 
      50          26 :   if (fespace->GetTypicalFE()->GetRangeType() == mfem::FiniteElement::SCALAR)
      51          13 :     m->AddDomainIntegrator(new mfem::MassIntegrator(one));
      52             :   else
      53          13 :     m->AddDomainIntegrator(new mfem::VectorFEMassIntegrator(one));
      54             : 
      55          26 :   m->Assemble();
      56             :   // Shift the eigenvalue corresponding to eliminated dofs to a large value. The BC DoFs on the
      57             :   // stiffness matrix are set to 1 and the mass matrix BC DoFs are set to a small value eps, such
      58             :   // that the eigenvaluesd associate with these DOFs are ~1/eps.
      59          26 :   m->EliminateEssentialBCDiag(_global_ess_markers, std::numeric_limits<mfem::real_t>::min());
      60          26 :   m->Finalize();
      61          26 :   op.Reset(m->ParallelAssemble());
      62          26 : }
      63             : 
      64             : void
      65          26 : EigenproblemEquationSystem::BuildEigenproblemJacobian(mfem::BlockVector & trueX,
      66             :                                                       mfem::OperatorHandle & massRHS)
      67             : {
      68             :   mooseAssert(_test_var_names.size() == 1 && (_test_var_names.size() == _trial_var_names.size()) &&
      69             :                   (_test_var_names.at(0) == _trial_var_names.at(0)),
      70             :               "Eigensolve is only supported for single-variable, square systems");
      71             : 
      72          26 :   height = trueX.Size();
      73          26 :   width = trueX.Size();
      74          26 :   ApplyEssentialBCs();
      75          26 :   FormEigenproblemMatrix(_jacobian);
      76          26 :   FormMassMatrix(massRHS);
      77          26 : }
      78             : 
      79             : } // namespace Moose::MFEM
      80             : 
      81             : #endif