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             : #include "ParsedExtraElementIDGenerator.h"
      11             : 
      12             : #include "Conversion.h"
      13             : #include "MooseMeshUtils.h"
      14             : 
      15             : #include "libmesh/fparser_ad.hh"
      16             : #include "libmesh/elem.h"
      17             : 
      18             : registerMooseObject("MooseApp", ParsedExtraElementIDGenerator);
      19             : 
      20             : InputParameters
      21        3263 : ParsedExtraElementIDGenerator::validParams()
      22             : {
      23        3263 :   InputParameters params = MeshGenerator::validParams();
      24        3263 :   params += FunctionParserUtils<false>::validParams();
      25             : 
      26       13052 :   params.addRequiredParam<MeshGeneratorName>("input", "The mesh we want to modify");
      27       13052 :   params.addRequiredParam<std::string>(
      28             :       "expression", "Function expression to return the extra element ID based on element centroid");
      29       13052 :   params.addRequiredParam<std::string>(
      30             :       "extra_elem_integer_name", "Name of the extra element integer to be added by this generator");
      31       13052 :   params.addParam<std::vector<SubdomainName>>("restricted_subdomains",
      32             :                                               "Only set ids for elements within these restricted "
      33             :                                               "subdomains (empty means the entire domain)");
      34       13052 :   params.addParam<std::vector<std::string>>(
      35             :       "constant_names", {}, "Vector of constants used in the parsed function");
      36       13052 :   params.addParam<std::vector<std::string>>(
      37             :       "constant_expressions",
      38             :       {},
      39             :       "Vector of values for the constants in constant_names (can be an FParser expression)");
      40       13052 :   params.addParam<std::vector<ExtraElementIDName>>(
      41             :       "extra_element_id_names", {}, "Extra element integers used in the parsed expression");
      42        3263 :   params.addClassDescription(
      43             :       "Uses a parsed expression to set an extra element id for elements (via their centroids).");
      44             : 
      45        3263 :   return params;
      46           0 : }
      47             : 
      48         101 : ParsedExtraElementIDGenerator::ParsedExtraElementIDGenerator(const InputParameters & parameters)
      49             :   : MeshGenerator(parameters),
      50             :     FunctionParserUtils<false>(parameters),
      51         101 :     _input(getMesh("input")),
      52         101 :     _function(parameters.get<std::string>("expression")),
      53         202 :     _eeid_names(getParam<std::vector<ExtraElementIDName>>("extra_element_id_names")),
      54         404 :     _elem_id_name(getParam<std::string>("extra_elem_integer_name"))
      55             : {
      56             :   // Form the vectors of constants names (symbols in the expression) and expression (values)
      57         202 :   auto c_names = getParam<std::vector<std::string>>("constant_names");
      58         303 :   auto c_defs = getParam<std::vector<std::string>>("constant_expressions");
      59         101 :   c_names.push_back("invalid_elem_id");
      60         101 :   c_defs.push_back(std::to_string(DofObject::invalid_id));
      61         101 :   c_names.push_back("pi");
      62         101 :   c_defs.push_back(std::to_string(libMesh::pi));
      63         101 :   c_names.push_back("e");
      64         101 :   c_defs.push_back(std::to_string(std::exp(Real(1))));
      65             : 
      66             :   // add the extra element integers
      67         101 :   std::string symbol_str = "x,y,z";
      68         111 :   for (const auto & eeid_name : _eeid_names)
      69          10 :     symbol_str += "," + eeid_name;
      70             : 
      71             :   // base function object
      72         101 :   _func_F = std::make_shared<SymFunction>();
      73         101 :   parsedFunctionSetup(_func_F, _function, symbol_str, c_names, c_defs, comm());
      74             : 
      75         101 :   _func_params.resize(3 + _eeid_names.size());
      76         101 : }
      77             : 
      78             : std::unique_ptr<MeshBase>
      79         100 : ParsedExtraElementIDGenerator::generate()
      80             : {
      81         100 :   std::unique_ptr<MeshBase> mesh = std::move(_input);
      82             : 
      83             :   // Make sure subdomain caches are up to date
      84         100 :   if (!mesh->preparation().has_cached_elem_data)
      85          60 :     mesh->cache_elem_data();
      86             : 
      87             :   // the extra element ids would not have existed at construction so we only do this now
      88         110 :   for (const auto & eeid_name : _eeid_names)
      89          10 :     _eeid_indices.push_back(mesh->get_elem_integer_index(eeid_name));
      90             : 
      91         200 :   bool has_restriction = isParamValid("restricted_subdomains");
      92         100 :   std::set<SubdomainID> restricted_blocks;
      93         100 :   if (has_restriction)
      94             :   {
      95          60 :     auto names = getParam<std::vector<SubdomainName>>("restricted_subdomains");
      96         180 :     for (auto & name : names)
      97             :     {
      98             :       // check that the subdomain exists in the mesh
      99         150 :       if (!MooseMeshUtils::hasSubdomainName(*mesh, name))
     100           0 :         paramError("restricted_subdomains", "The block '", name, "' was not found in the mesh");
     101             : 
     102         150 :       restricted_blocks.insert(MooseMeshUtils::getSubdomainID(name, *mesh));
     103             :     }
     104          30 :   }
     105             : 
     106         100 :   if (!mesh->has_elem_integer(_elem_id_name))
     107          90 :     _extra_elem_id = mesh->add_elem_integer(_elem_id_name);
     108             :   else
     109          10 :     _extra_elem_id = mesh->get_elem_integer_index(_elem_id_name);
     110             : 
     111             :   // Loop over the elements
     112        5975 :   for (const auto & elem : mesh->active_element_ptr_range())
     113             :   {
     114        5875 :     if (has_restriction && restricted_blocks.count(elem->subdomain_id()) == 0)
     115        2048 :       continue;
     116             : 
     117        3827 :     const auto centroid = elem->true_centroid();
     118        3827 :     _func_params[0] = centroid(0);
     119        3827 :     _func_params[1] = centroid(1);
     120        3827 :     _func_params[2] = centroid(2);
     121        4419 :     for (const auto i : index_range(_eeid_indices))
     122         592 :       _func_params[3 + i] = elem->get_extra_integer(_eeid_indices[i]);
     123        7654 :     const auto id_real = evaluate(_func_F);
     124             : 
     125        3827 :     dof_id_type id = id_real;
     126             :     // this is to ensure a more robust conversion between Real and dof_id_type
     127        3827 :     if (id_real == (Real)DofObject::invalid_id)
     128         144 :       id = DofObject::invalid_id;
     129             : 
     130        3827 :     elem->set_extra_integer(_extra_elem_id, id);
     131         100 :   }
     132             : 
     133         200 :   return mesh;
     134         100 : }