Line data    Source code

       1             : // The libMesh Finite Element Library.
       2             : // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
       3             : 
       4             : // This library is free software; you can redistribute it and/or
       5             : // modify it under the terms of the GNU Lesser General Public
       6             : // License as published by the Free Software Foundation; either
       7             : // version 2.1 of the License, or (at your option) any later version.
       8             : 
       9             : // This library is distributed in the hope that it will be useful,
      10             : // but WITHOUT ANY WARRANTY; without even the implied warranty of
      11             : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
      12             : // Lesser General Public  License for more details.
      13             : 
      14             : // You should have received a copy of the GNU Lesser General Public
      15             : // License along with this library; if not, write to the Free Software
      16             : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
      17             : 
      18             : 
      19             : 
      20             : // Local Includes
      21             : #include "libmesh/default_coupling.h"
      22             : #include "libmesh/coupling_matrix.h"
      23             : #include "libmesh/elem.h"
      24             : #include "libmesh/periodic_boundaries.h"
      25             : #include "libmesh/remote_elem.h"
      26             : #include "libmesh/int_range.h"
      27             : #include "libmesh/libmesh_logging.h"
      28             : 
      29             : // C++ Includes
      30             : #include <unordered_set>
      31             : 
      32             : namespace libMesh
      33             : {
      34             : 
      35      500636 : void DefaultCoupling::set_dof_coupling(const CouplingMatrix * dof_coupling)
      36             : {
      37             :   // We used to treat an empty 0x0 _dof_coupling matrix as if it
      38             :   // were an NxN all-ones matrix.  We'd like to stop supporting this
      39             :   // behavior, but for now we'll just warn about it, while supporting
      40             :   // it via the preferred mechanism: a nullptr _dof_coupling
      41             :   // matrix pointer is interpreted as a full coupling matrix.
      42      500636 :   if (dof_coupling && dof_coupling->empty())
      43             :     {
      44             :       libmesh_deprecated();
      45           0 :       _dof_coupling = nullptr;
      46             :     }
      47             :   else
      48      500636 :     _dof_coupling = dof_coupling;
      49      500636 : }
      50             : 
      51             : 
      52             : 
      53      435728 : void DefaultCoupling::mesh_reinit()
      54             : {
      55             :   // Unless we have periodic boundary conditions, we don't need
      56             :   // anything precomputed.
      57             : #ifdef LIBMESH_ENABLE_PERIODIC
      58      435728 :   if (!_periodic_bcs || _periodic_bcs->empty())
      59        7548 :     return;
      60             : #endif
      61             : 
      62             :   // If we do have periodic boundary conditions, we'll need a master
      63             :   // point locator, so we'd better have a mesh to build it on.
      64         636 :   libmesh_assert(_mesh);
      65             : 
      66             :   // Make sure an up-to-date master point locator has been
      67             :   // constructed; we'll need to grab sub-locators soon.
      68        2642 :   _mesh->sub_point_locator();
      69             : }
      70             : 
      71             : 
      72             : 
      73     4329778 : void DefaultCoupling::operator()
      74             :   (const MeshBase::const_element_iterator & range_begin,
      75             :    const MeshBase::const_element_iterator & range_end,
      76             :    processor_id_type p,
      77             :    map_type & coupled_elements)
      78             : {
      79      331277 :   LOG_SCOPE("operator()", "DefaultCoupling");
      80             : 
      81             :   // Let us not do assertion at this moment for API upgrade.
      82             :   // There is a functor inside of ElementSideNeighborLayers.
      83             :   // We can not set mesh for that functor because there is no handle
      84             :   // in libmesh. We need to override set_mesh in moose for setting a mesh for the functor.
      85             :   // The set_mesh overridden will not happen until the current change gets in.
      86             :   //libmesh_assert(_mesh);
      87             : 
      88             : #ifdef LIBMESH_ENABLE_PERIODIC
      89             :   bool check_periodic_bcs =
      90     4329778 :     (_periodic_bcs && !_periodic_bcs->empty());
      91             : 
      92     3577232 :   std::unique_ptr<PointLocatorBase> point_locator;
      93      331277 :   if (check_periodic_bcs)
      94             :     {
      95       45725 :       libmesh_assert(_mesh);
      96      271432 :       point_locator = _mesh->sub_point_locator();
      97             :     }
      98             : #endif
      99             : 
     100     4329778 :   if (!this->_n_levels)
     101             :     {
     102    25717574 :       for (const auto & elem : as_range(range_begin, range_end))
     103             :       {
     104             :         //libmesh_assert(_mesh->query_elem_ptr(elem->id()) ==elem);
     105    10585053 :         if (elem->processor_id() != p)
     106     2674302 :           coupled_elements.emplace(elem, _dof_coupling);
     107     2563758 :       }
     108      246675 :       return;
     109             :     }
     110             : 
     111             :   typedef std::unordered_set<const Elem*> set_type;
     112     2629940 :   set_type next_elements_to_check(range_begin, range_end);
     113      169204 :   set_type elements_to_check;
     114      169204 :   set_type elements_checked;
     115             : 
     116     2546914 :   for (unsigned int i=0; i != this->_n_levels; ++i)
     117             :     {
     118       84618 :       elements_to_check.swap(next_elements_to_check);
     119       84618 :       next_elements_to_check.clear();
     120     1189627 :       elements_checked.insert(elements_to_check.begin(), elements_to_check.end());
     121             : 
     122    10995544 :       for (const auto & elem : elements_to_check)
     123             :         {
     124     1282558 :           std::vector<const Elem *> active_neighbors;
     125             : 
     126             :           //libmesh_assert(_mesh->query_elem_ptr(elem->id()) ==elem);
     127             : 
     128     9721299 :           if (elem->processor_id() != p)
     129      678909 :             coupled_elements.emplace(elem, _dof_coupling);
     130             : 
     131    49179081 :           for (auto s : elem->side_index_range())
     132             :             {
     133    38816502 :               const Elem * neigh = elem->neighbor_ptr(s);
     134             : 
     135             : #ifdef LIBMESH_ENABLE_PERIODIC
     136             :               // We might still have a periodic neighbor here
     137    38816502 :               if (!neigh && check_periodic_bcs)
     138             :                 {
     139        5976 :                   libmesh_assert(_mesh);
     140             : 
     141             :                   neigh = elem->topological_neighbor
     142       27394 :                     (s, *_mesh, *point_locator, _periodic_bcs);
     143             :                 }
     144             : #endif
     145             : 
     146             :               // With no regular *or* periodic neighbors we have nothing
     147             :               // to do. *Or* Mesh ghosting might ask us about what we want to
     148             :               // distribute along with non-local elements, and those
     149             :               // non-local elements might have remote neighbors, and
     150             :               // if they do then we can't say anything about them.
     151    38816502 :               if (!neigh || neigh == remote_elem)
     152     2100548 :                 continue;
     153             : 
     154             :               // With any kind of neighbor, we need to couple to all the
     155             :               // active descendants on our side.
     156             : #ifdef LIBMESH_ENABLE_AMR
     157    38955182 :               if (neigh == elem->neighbor_ptr(s))
     158    36568676 :                 neigh->active_family_tree_by_neighbor(active_neighbors,elem);
     159             : #  ifdef LIBMESH_ENABLE_PERIODIC
     160             :               else
     161             :                 neigh->active_family_tree_by_topological_neighbor
     162       15046 :                   (active_neighbors,elem,*_mesh,*point_locator,_periodic_bcs);
     163             : #  endif
     164             : #else
     165             :               active_neighbors.clear();
     166             :               active_neighbors.push_back(neigh);
     167             : #endif
     168             : 
     169    73742397 :               for (const auto & neighbor : active_neighbors)
     170             :                 {
     171     4857996 :                   if (!elements_checked.count(neighbor))
     172      344133 :                     next_elements_to_check.insert(neighbor);
     173             : 
     174    37158675 :                   if (neighbor->processor_id() != p)
     175     5507031 :                     coupled_elements.emplace(neighbor, _dof_coupling);
     176             :                 }
     177             :             }
     178             :         }
     179             :     }
     180     3667223 : }
     181             : 
     182             : 
     183             : } // namespace libMesh