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 "FunctorExtremaPositions.h"
      11             : #include "libmesh/parallel.h"
      12             : #include "libmesh/parallel_algebra.h"
      13             : #include "libmesh/vector_value.h"
      14             : 
      15             : registerMooseObject("MooseApp", FunctorExtremaPositions);
      16             : 
      17             : InputParameters
      18       14291 : FunctorExtremaPositions::validParams()
      19             : {
      20       14291 :   InputParameters params = Positions::validParams();
      21       14291 :   params += NonADFunctorInterface::validParams();
      22       14291 :   params += BlockRestrictable::validParams();
      23             : 
      24       14291 :   params.addRequiredParam<MooseFunctorName>(
      25             :       "functor", "Functor of which the extremas are going to form the list of positions");
      26       42873 :   params.addParam<MooseEnum>("extrema_type",
      27       28582 :                              MooseEnum("max=0 min=1 max_abs=2", "max"),
      28             :                              "Type of extreme value to return. 'max' "
      29             :                              "returns the location(s) of the maximum value(s). 'min' returns "
      30             :                              "the location(s) of the minimum value(s). 'max_abs' returns the "
      31             :                              "locations of the maximum(a) of the absolute value.");
      32       14291 :   params.addParam<unsigned int>("num_extrema", 1, "Number of extrema to look for");
      33             : 
      34             :   // Future arguments we could conceive:
      35             :   // - min_dist to force a certain distance between extrema
      36             :   // - functor argument type. We use elements, we could do nodes, qps, sides, ...
      37             : 
      38             :   // Use user-provided ordering
      39       14291 :   params.set<bool>("auto_sort") = false;
      40       14291 :   params.suppressParameter<bool>("auto_sort");
      41             :   // The broadcasting has to be custom to preserve the order of the extremas
      42       14291 :   params.set<bool>("auto_broadcast") = false;
      43       14291 :   params.suppressParameter<bool>("auto_broadcast");
      44             : 
      45             :   // Keep as up-to-date as possible given the generality of functors
      46       42873 :   params.set<ExecFlagEnum>("execute_on") = {EXEC_LINEAR, EXEC_TIMESTEP_BEGIN};
      47             : 
      48       14291 :   params.addClassDescription("Searches the geometry for the extrema of the functors");
      49       14291 :   return params;
      50       14291 : }
      51             : 
      52          13 : FunctorExtremaPositions::FunctorExtremaPositions(const InputParameters & parameters)
      53             :   : Positions(parameters),
      54             :     NonADFunctorInterface(this),
      55             :     BlockRestrictable(this),
      56          13 :     _functor(getFunctor<Real>("functor")),
      57          13 :     _n_extrema(getParam<unsigned int>("num_extrema")),
      58          13 :     _type(getParam<MooseEnum>("extrema_type").getEnum<ExtremeType>()),
      59          13 :     _positions_values(declareValueByName<std::vector<Real>, ReporterVectorContext<Real>>(
      60          13 :         "functor_extrema", REPORTER_MODE_REPLICATED))
      61             : {
      62             :   // Constants and postprocessors are not interesting here
      63          13 :   const auto & functor_name = getParam<MooseFunctorName>("functor");
      64          26 :   if (_fe_problem.hasPostprocessorValueByName(functor_name) ||
      65          13 :       MooseUtils::parsesToReal(functor_name))
      66           0 :     paramError(
      67             :         "functor",
      68             :         "Postprocessors and constants do not have extrema, they are constant over the domain.");
      69             : 
      70          13 :   if (getParam<ExecFlagEnum>("execute_on").contains(EXEC_NONE))
      71           0 :     paramError("execute_on",
      72             :                "NONE execution flag not supported. Most functors (functions, variables, spatial "
      73             :                "user objects for example) are not initialized at construction.");
      74          13 : }
      75             : 
      76             : void
      77          48 : FunctorExtremaPositions::initialize()
      78             : {
      79          48 :   clearPositions();
      80          48 :   _positions.resize(_n_extrema);
      81          48 :   _positions_values.resize(_n_extrema);
      82             : 
      83          48 :   std::list<Real> extrema;
      84          48 :   std::list<Point> extrema_locs;
      85          48 :   if (_type == ExtremeType::MAX || _type == ExtremeType::MAX_ABS)
      86          48 :     extrema.resize(_n_extrema, -std::numeric_limits<Real>::max());
      87             :   else
      88           0 :     extrema.resize(_n_extrema, std::numeric_limits<Real>::max());
      89          48 :   extrema_locs.resize(_n_extrema);
      90             : 
      91          48 :   const auto time_arg = determineState();
      92             : 
      93             :   // Loop over the local mesh, keep track of all extrema
      94          48 :   for (const auto & elem :
      95        7296 :        _fe_problem.mesh().getMesh().active_local_subdomain_set_element_ptr_range(blockIDs()))
      96             :   {
      97        3600 :     const Moose::ElemArg elem_arg = {elem, false};
      98        3600 :     auto value = _functor(elem_arg, time_arg);
      99             : 
     100        3600 :     auto extremum = extrema.begin();
     101        3600 :     auto extremum_loc = extrema_locs.begin();
     102             : 
     103       10400 :     for (const auto i : make_range(_n_extrema))
     104             :     {
     105       10400 :       libmesh_ignore(i);
     106       10400 :       bool extrema_found = false;
     107       10400 :       switch (_type)
     108             :       {
     109       10400 :         case ExtremeType::MAX_ABS:
     110       10400 :           value = std::abs(value);
     111             :           // fallthrough
     112       10400 :         case ExtremeType::MAX:
     113             :         {
     114       10400 :           if (value > *extremum)
     115        3600 :             extrema_found = true;
     116       10400 :           break;
     117             :         }
     118           0 :         case ExtremeType::MIN:
     119             :         {
     120           0 :           if (value < *extremum)
     121           0 :             extrema_found = true;
     122           0 :           break;
     123             :         }
     124             :       }
     125       10400 :       if (extrema_found)
     126             :       {
     127             :         // Move the extrema down the list of extrema
     128        3600 :         extrema.insert(extremum, value);
     129        3600 :         extrema_locs.insert(extremum_loc, elem->true_centroid());
     130        3600 :         extrema.pop_back();
     131        3600 :         extrema_locs.pop_back();
     132             :         // Found an extremum, filled the extrema vector, done
     133        3600 :         break;
     134             :       }
     135             :       // Consider the next extremum
     136        6800 :       std::advance(extremum, 1);
     137        6800 :       std::advance(extremum_loc, 1);
     138             :     }
     139          48 :   }
     140             : 
     141             :   // Synchronize across all ranks
     142          48 :   auto current_candidate = extrema.begin();
     143          48 :   auto current_candidate_loc = extrema_locs.begin();
     144         528 :   for (const auto i : make_range(_n_extrema))
     145             :   {
     146             :     unsigned int rank;
     147             :     // dont modify the extremum, it might be the global n-th extremum in a later call
     148         480 :     auto copy = *current_candidate;
     149         480 :     if (_type == ExtremeType::MAX || _type == ExtremeType::MAX_ABS)
     150         480 :       comm().maxloc(copy, rank);
     151             :     else
     152           0 :       comm().minloc(copy, rank);
     153             : 
     154         480 :     RealVectorValue extreme_point(0., 0., 0.);
     155         480 :     if (rank == processor_id())
     156             :     {
     157         360 :       extreme_point = *current_candidate_loc;
     158             :       // Our candidate for ith max got accepted, move on to offering the next extremal one
     159         360 :       std::advance(current_candidate, 1);
     160         360 :       std::advance(current_candidate_loc, 1);
     161             :     }
     162             : 
     163             :     // Send the position
     164             :     // The broadcast and scatter are for root 0 sending to all. The sum works just as well
     165         480 :     comm().sum(extreme_point);
     166         480 :     _positions[i] = extreme_point;
     167         480 :     _positions_values[i] = copy;
     168             :   }
     169             :   mooseAssert(comm().verify(_positions_values),
     170             :               "Positions extrema should be the same across all MPI processes.");
     171          48 :   _initialized = true;
     172          48 : }