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 MFEM_ENABLED
      11             : 
      12             : #pragma once
      13             : #include <map>
      14             : #include <memory>
      15             : #include <string>
      16             : #include <tuple>
      17             : #include <utility>
      18             : #include <variant>
      19             : #include <vector>
      20             : 
      21             : #include "MooseException.h"
      22             : #include "MooseError.h"
      23             : 
      24             : #include "libmesh/ignore_warnings.h"
      25             : #include <mfem.hpp>
      26             : #include "libmesh/restore_warnings.h"
      27             : #include "libmesh/utility.h"
      28             : 
      29             : namespace Moose::MFEM
      30             : {
      31             : 
      32             : /**
      33             :  * Class to manage MFEM coefficient objects representing material
      34             :  * properties. It can build up piecewise coefficients representing
      35             :  * properties defined across multiple materials.
      36             :  */
      37             : template <class T, class Tpw>
      38             : class CoefficientMap
      39             : {
      40             : public:
      41         964 :   CoefficientMap() = default;
      42             : 
      43             :   /// Make arbitrary coefficients which will be tracked by this object
      44             :   template <class P, class... Args>
      45         879 :   std::shared_ptr<P> make(Args &&... args)
      46             :   {
      47         879 :     auto result = std::make_shared<P>(args...);
      48         879 :     this->_iterable_coefficients.push_back(result);
      49         879 :     return result;
      50           0 :   }
      51             : 
      52             :   /// Add a named global coefficient. It must have been created with
      53             :   /// the `make` method on this object.
      54         752 :   void addCoefficient(const std::string & name, std::shared_ptr<T> coeff)
      55             :   {
      56             :     mooseAssert(std::find(this->_iterable_coefficients.cbegin(),
      57             :                           this->_iterable_coefficients.cend(),
      58             :                           coeff) != this->_iterable_coefficients.cend(),
      59             :                 "Coefficient object was not created by this CoefficientMap.");
      60             : 
      61         752 :     const auto [_, inserted] = this->_coefficients.emplace(name, std::move(coeff));
      62         752 :     if (!inserted)
      63          12 :       mooseError("Coefficient with name '" + name + "' already present in CoefficientMap object");
      64         748 :   }
      65             : 
      66             :   /// Add piecewise material property. The coefficient must have been created with the `make` method on this object.
      67             :   ///
      68             :   /// Note: If you attempt to overwrite an existing block then an exception will be thrown and data
      69             :   /// for that property will be left in an undefined state.
      70         138 :   void addPiecewiseBlocks(const std::string & name,
      71             :                           std::shared_ptr<T> coeff,
      72             :                           const std::vector<std::string> & blocks)
      73             :   {
      74             :     // If list of blocks is empty then treat as a global coefficient
      75         138 :     if (blocks.size() == 0)
      76             :     {
      77           8 :       this->addCoefficient(name, coeff);
      78           8 :       return;
      79             :     }
      80             : 
      81             :     // Initialise property with empty coefficients, if it does not already exist
      82         130 :     if (!this->hasCoefficient(name))
      83          87 :       this->_coefficients.insert({name, this->emptyPWData(coeff)});
      84         130 :     PWData * data = std::get_if<PWData>(&this->_coefficients[name]);
      85             :     // Throw an exception if the data is not piecewise
      86         130 :     if (!data)
      87           3 :       mooseError("Global coefficient with name '" + name + "' already present in CoefficientMap");
      88             :     mooseAssert(std::find(this->_iterable_coefficients.cbegin(),
      89             :                           this->_iterable_coefficients.cend(),
      90             :                           coeff) != this->_iterable_coefficients.cend(),
      91             :                 "Coefficient object was not created by the appropriate coefficient manager.");
      92         129 :     auto & [pw_coeff, coeff_map] = *data;
      93         133 :     this->checkPWData(coeff, pw_coeff, name);
      94             : 
      95         303 :     for (const auto & block : blocks)
      96             :     {
      97         177 :       if (coeff_map.count(block) > 0)
      98           5 :         mooseError("Property with name '" + name + "' already assigned to block " + block +
      99             :                    " in CoefficientMap object");
     100         176 :       coeff_map[block] = coeff;
     101         176 :       pw_coeff->UpdateCoefficient(std::stoi(block), *coeff);
     102             :     }
     103             :   }
     104             : 
     105           7 :   T & getCoefficient(const std::string & name) { return *this->getCoefficientPtr(name); }
     106             : 
     107         736 :   std::shared_ptr<T> getCoefficientPtr(const std::string & name)
     108             :   {
     109             :     try
     110             :     {
     111         736 :       auto & coeff = this->_coefficients.at(name);
     112             :       try
     113             :       {
     114         733 :         return std::get<std::shared_ptr<T>>(coeff);
     115             :       }
     116         380 :       catch (const std::bad_variant_access &)
     117             :       {
     118         190 :         return std::get<0>(std::get<PWData>(coeff));
     119             :       }
     120             :     }
     121           6 :     catch (const std::out_of_range &)
     122             :     {
     123           9 :       mooseError("Property with name '" + name + "' has not been declared.");
     124             :     }
     125             :   }
     126             : 
     127         970 :   bool hasCoefficient(const std::string & name) const
     128             :   {
     129         970 :     return this->_coefficients.count(name) > 0;
     130             :   }
     131             : 
     132         130 :   bool propertyDefinedOnBlock(const std::string & name, const std::string & block) const
     133             :   {
     134         130 :     if (!this->hasCoefficient(name))
     135          78 :       return false;
     136          52 :     auto & coeff = libmesh_map_find(this->_coefficients, name);
     137          52 :     if (std::holds_alternative<std::shared_ptr<T>>(coeff))
     138          21 :       return true;
     139          31 :     auto block_map = std::get<1>(std::get<PWData>(coeff));
     140          31 :     return block_map.count(block) > 0;
     141          31 :   }
     142             : 
     143         204 :   void setTime(const double time)
     144             :   {
     145         600 :     for (auto & coef : this->_iterable_coefficients)
     146         396 :       coef->SetTime(time);
     147         204 :   }
     148             : 
     149             : private:
     150             :   using PWData = std::tuple<std::shared_ptr<Tpw>, std::map<const std::string, std::shared_ptr<T>>>;
     151             :   std::map<const std::string, std::variant<std::shared_ptr<T>, PWData>> _coefficients;
     152             :   std::vector<std::shared_ptr<T>> _iterable_coefficients;
     153             : 
     154          52 :   PWData emptyPWData(std::shared_ptr<T> /*coeff*/)
     155             :   {
     156             :     return std::make_tuple(this->template make<Tpw>(),
     157          52 :                            std::map<const std::string, std::shared_ptr<T>>());
     158             :   }
     159          76 :   void checkPWData(std::shared_ptr<T> /*coeff*/,
     160             :                    std::shared_ptr<Tpw> /* existing_pw*/,
     161             :                    const std::string & /*name*/)
     162             :   {
     163          76 :   }
     164             : };
     165             : 
     166             : using ScalarMap = CoefficientMap<mfem::Coefficient, mfem::PWCoefficient>;
     167             : using VectorMap = CoefficientMap<mfem::VectorCoefficient, mfem::PWVectorCoefficient>;
     168             : using MatrixMap = CoefficientMap<mfem::MatrixCoefficient, mfem::PWMatrixCoefficient>;
     169             : 
     170             : template <>
     171             : VectorMap::PWData VectorMap::emptyPWData(std::shared_ptr<mfem::VectorCoefficient> coeff);
     172             : 
     173             : template <>
     174             : MatrixMap::PWData MatrixMap::emptyPWData(std::shared_ptr<mfem::MatrixCoefficient> coeff);
     175             : 
     176             : template <>
     177             : void VectorMap::checkPWData(std::shared_ptr<mfem::VectorCoefficient> coeff,
     178             :                             std::shared_ptr<mfem::PWVectorCoefficient> existing_pw,
     179             :                             const std::string & name);
     180             : 
     181             : template <>
     182             : void MatrixMap::checkPWData(std::shared_ptr<mfem::MatrixCoefficient> coeff,
     183             :                             std::shared_ptr<mfem::PWMatrixCoefficient> existing_pw,
     184             :                             const std::string & name);
     185             : } // namespace Moose::MFEM
     186             : 
     187             : #endif