Pass only retained nodes to MueLu

src/core/fem/src/discretization/4C_fem_discretization.cpp

@@ -563,12 +563,13 @@ void Core::FE::Discretization::replace_dof_set(
 
 /*----------------------------------------------------------------------*
  *----------------------------------------------------------------------*/
-std::map<int, std::vector<int>>* Core::FE::Discretization::get_all_pbc_coupled_col_nodes()
+const std::map<int, std::vector<int>>* Core::FE::Discretization::get_all_pbc_coupled_col_nodes()
+    const
 {
   // check for pbcs
   for (int nds = 0; nds < num_dof_sets(); nds++)
   {
-    std::shared_ptr<Core::DOFSets::PBCDofSet> pbcdofset =
+    std::shared_ptr<const Core::DOFSets::PBCDofSet> pbcdofset =
         std::dynamic_pointer_cast<Core::DOFSets::PBCDofSet>(dofsets_[nds]);
 
     if (pbcdofset != nullptr)
@@ -584,8 +585,8 @@ std::map<int, std::vector<int>>* Core::FE::Discretization::get_all_pbc_coupled_c
 
 /*----------------------------------------------------------------------*
  *----------------------------------------------------------------------*/
-std::shared_ptr<std::map<int, int>>
-Core::FE::Discretization::get_pbc_slave_to_master_node_connectivity()
+std::shared_ptr<const std::map<int, int>>
+Core::FE::Discretization::get_pbc_slave_to_master_node_connectivity() const
 {
   // check for pbcs
   for (int nds = 0; nds < num_dof_sets(); nds++)

src/core/fem/src/discretization/4C_fem_discretization.hpp

@@ -622,12 +622,12 @@ namespace Core::FE
     /*!
     \brief Get master to slave coupling in case of periodic boundary conditions
     */
-    std::map<int, std::vector<int>>* get_all_pbc_coupled_col_nodes();
+    const std::map<int, std::vector<int>>* get_all_pbc_coupled_col_nodes() const;
 
     /*!
     \brief Get slave to master connectivity in case of periodic boundary conditions
     */
-    std::shared_ptr<std::map<int, int>> get_pbc_slave_to_master_node_connectivity();
+    std::shared_ptr<const std::map<int, int>> get_pbc_slave_to_master_node_connectivity() const;
 
     //@}
 

src/core/fem/src/discretization/4C_fem_discretization_faces.cpp

@@ -239,7 +239,7 @@ void Core::FE::DiscretizationFaces::build_faces(const bool verbose)
   std::map<std::vector<int>, std::shared_ptr<Core::Elements::Element>> faces;
 
   // get pbcs
-  std::map<int, std::vector<int>>* col_pbcmapmastertoslave = get_all_pbc_coupled_col_nodes();
+  const std::map<int, std::vector<int>>* col_pbcmapmastertoslave = get_all_pbc_coupled_col_nodes();
 
   std::map<std::vector<int>, InternalFacesData>::iterator face_it;
   for (face_it = surfmapdata.begin(); face_it != surfmapdata.end(); ++face_it)
@@ -383,11 +383,11 @@ void Core::FE::DiscretizationFaces::build_faces(const bool verbose)
             for (std::size_t rr = 0; rr < mymasternodeids.size(); rr++)
             {
               // this master node has one slave node
-              if (((*col_pbcmapmastertoslave)[mymasternodeids[rr]]).size() == 1)
+              if (((*col_pbcmapmastertoslave).at(mymasternodeids[rr])).size() == 1)
               {
-                myslavenodeids.push_back(((*col_pbcmapmastertoslave)[mymasternodeids[rr]])[0]);
+                myslavenodeids.push_back(((*col_pbcmapmastertoslave).at(mymasternodeids[rr]))[0]);
                 local_pbcmapmastertoslave[mymasternodeids[rr]] =
-                    ((*col_pbcmapmastertoslave)[mymasternodeids[rr]])[0];
+                    ((*col_pbcmapmastertoslave).at(mymasternodeids[rr]))[0];
               }
               // this master node has several slave nodes
               // it is a corner or edge node of two or three pbc sets
@@ -443,7 +443,7 @@ void Core::FE::DiscretizationFaces::build_faces(const bool verbose)
                 //  this special case is marked by flag
                 bool three_sets_edge_node = false;
                 if (numpbcpairs == 3 and
-                    ((*col_pbcmapmastertoslave)[mymasternodeids[rr]]).size() == 3)
+                    ((*col_pbcmapmastertoslave).at(mymasternodeids[rr])).size() == 3)
                   three_sets_edge_node = true;
 
                 // pbc id of master face also for the slave
@@ -459,10 +459,10 @@ void Core::FE::DiscretizationFaces::build_faces(const bool verbose)
                 // check which node fulfills above conditions
                 int actslaveid = -999;
                 for (std::size_t islave = 0;
-                    islave < ((*col_pbcmapmastertoslave)[mymasternodeids[rr]]).size(); islave++)
+                    islave < ((*col_pbcmapmastertoslave).at(mymasternodeids[rr])).size(); islave++)
                 {
                   // get id
-                  actslaveid = ((*col_pbcmapmastertoslave)[mymasternodeids[rr]])[islave];
+                  actslaveid = ((*col_pbcmapmastertoslave).at(mymasternodeids[rr]))[islave];
 
                   // check first criterion -> (i)
                   if ((slavetopbcset[slavepbcid]).find(actslaveid) !=
@@ -547,7 +547,7 @@ void Core::FE::DiscretizationFaces::build_faces(const bool verbose)
 
                 // first, look for masters with more than one slave
                 // then, check whether node id (i.e. actnodeid) is contained in slave list
-                std::map<int, std::vector<int>>::iterator master_it;
+                std::map<int, std::vector<int>>::const_iterator master_it;
                 for (master_it = col_pbcmapmastertoslave->begin();
                     master_it != col_pbcmapmastertoslave->end(); master_it++)
                 {

src/core/fem/src/discretization/4C_fem_discretization_partition.cpp

@@ -569,7 +569,7 @@ void Core::FE::Discretization::extended_ghosting(const Core::LinAlg::Map& elecol
     {
       have_pbc = true;
       // fill content of pbcmap int std::map<int, std::set<int> > in preparation for gather_all
-      std::map<int, std::vector<int>>* tmp = pbcdofset->get_coupled_nodes();
+      const std::map<int, std::vector<int>>* tmp = pbcdofset->get_coupled_nodes();
       for (auto& [gid, nodes] : *tmp) pbcmap[gid].insert(nodes.begin(), nodes.end());
 
       // it is assumed that, if one pbc set is available, all other potential dofsets hold the same

src/core/fem/src/discretization/4C_fem_discretization_utils.cpp

@@ -7,6 +7,7 @@
 
 #include "4C_fem_discretization_utils.hpp"
 
+#include "4C_comm_mpi_utils.hpp"
 #include "4C_fem_discretization.hpp"
 #include "4C_fem_general_element.hpp"
 #include "4C_fem_general_node.hpp"
@@ -48,6 +49,57 @@ std::shared_ptr<Core::LinAlg::MultiVector<double>> Core::FE::extract_node_coordi
   return coordinates;
 }
 
+std::unique_ptr<Core::LinAlg::MultiVector<double>> Core::FE::extract_retained_node_coordinates(
+    const Discretization& discretization, const Core::LinAlg::Map& node_row_map)
+{
+  std::set<int> skipped_nodes = {};
+  if (discretization.get_pbc_slave_to_master_node_connectivity())
+  {
+    const auto& pbcconnectivity = *(discretization.get_pbc_slave_to_master_node_connectivity());
+    for (const auto& [slave, _] : pbcconnectivity)
+    {
+      skipped_nodes.insert(slave);
+    }
+  }
+  std::set<int> fully_redundant_skipped_nodes =
+      Core::Communication::all_reduce(skipped_nodes, node_row_map.get_comm());
+
+  std::vector<int> my_filtered_global_node_ids{};
+  my_filtered_global_node_ids.reserve(node_row_map.num_my_elements());
+  for (int lid = 0; lid < node_row_map.num_my_elements(); ++lid)
+  {
+    const int gid = node_row_map.gid(lid);
+
+    // filter pbc slave nodes
+    if (fully_redundant_skipped_nodes.contains(gid)) continue;
+    my_filtered_global_node_ids.emplace_back(gid);
+  }
+  Core::LinAlg::Map filtered_row_map(
+      node_row_map.num_global_elements() - static_cast<int>(fully_redundant_skipped_nodes.size()),
+      static_cast<int>(my_filtered_global_node_ids.size()), my_filtered_global_node_ids.data(), 0,
+      node_row_map.get_comm());
+
+  std::unique_ptr<Core::LinAlg::MultiVector<double>> coordinates =
+      std::make_unique<Core::LinAlg::MultiVector<double>>(filtered_row_map, 3, true);
+
+  for (int lid = 0; lid < filtered_row_map.num_my_elements(); ++lid)
+  {
+    const int gid = filtered_row_map.gid(lid);
+    if (!discretization.have_global_node(gid)) continue;
+
+    auto x = discretization.g_node(gid)->x();
+    for (size_t dim = 0; dim < 3; ++dim)
+    {
+      if (dim >= discretization.n_dim())
+        coordinates->replace_local_value(lid, dim, 0.0);
+      else
+        coordinates->replace_local_value(lid, dim, x[dim]);
+    }
+  }
+
+  return coordinates;
+}
+
 
 /*----------------------------------------------------------------------------*
  *----------------------------------------------------------------------------*/

src/core/fem/src/discretization/4C_fem_discretization_utils.hpp

@@ -60,6 +60,17 @@ namespace Core::FE
   std::shared_ptr<Core::LinAlg::MultiVector<double>> extract_node_coordinates(
       const Core::FE::Discretization& discretization, const Core::LinAlg::Map& node_row_map);
 
+  /**
+   * @brief Extract all retained node coordinates of a discretization.
+   *
+   * @note Allows to specify a node row map of the non-eliminated nodes that is used to obtain the
+   * filtered node row map. The map has to be a submap of the overall full node rowmap of this
+   * discretization.
+   * @note Eliminated nodes (e.g., when using periodic boundary conditions) are skipped.
+   */
+  std::unique_ptr<Core::LinAlg::MultiVector<double>> extract_retained_node_coordinates(
+      const Core::FE::Discretization& discretization, const Core::LinAlg::Map& node_row_map);
+
 
   /** \brief Evaluate the elements of the given discretization and fill the
    *         system matrix and vector

src/core/fem/src/dofset/4C_fem_dofset_pbc.hpp

@@ -70,10 +70,14 @@ namespace Core::DOFSets
     virtual void set_coupled_nodes(std::shared_ptr<std::map<int, std::vector<int>>> couplednodes);
 
     /// Get coupled nodes map (corresponding col format)
-    std::map<int, std::vector<int>>* get_coupled_nodes() { return perbndcouples_.get(); }
+    [[nodiscard]] const std::map<int, std::vector<int>>* get_coupled_nodes() const
+    {
+      return perbndcouples_.get();
+    }
 
     /// Get connectivity map between slave node and its master node
-    virtual std::shared_ptr<std::map<int, int>> get_slave_to_master_node_connectivity()
+    [[nodiscard]] virtual std::shared_ptr<const std::map<int, int>>
+    get_slave_to_master_node_connectivity() const
     {
       return perbnd_slavetomaster_;
     };

src/core/fem/src/general/utils/4C_fem_general_l2_projection.cpp

@@ -135,7 +135,7 @@ std::shared_ptr<Core::LinAlg::MultiVector<double>> Core::FE::compute_nodal_l2_pr
   // build inverse map from slave to master nodes
   std::map<int, int> slavetomastercolnodesmap;
 
-  std::map<int, std::vector<int>>* allcoupledcolnodes = dis.get_all_pbc_coupled_col_nodes();
+  const std::map<int, std::vector<int>>* allcoupledcolnodes = dis.get_all_pbc_coupled_col_nodes();
   if (allcoupledcolnodes)
   {
     for (auto [master_gid, slave_gids] : *allcoupledcolnodes)

src/core/fem/src/general/utils/4C_fem_general_utils_superconvergent_patch_recovery.cpp

@@ -41,7 +41,7 @@ std::shared_ptr<Core::LinAlg::MultiVector<double>> Core::FE::compute_superconver
   // handle pbcs if existing
   // build inverse map from slave to master nodes
   std::map<int, int> slavetomastercolnodesmap;
-  std::map<int, std::vector<int>>* allcoupledcolnodes = dis.get_all_pbc_coupled_col_nodes();
+  const std::map<int, std::vector<int>>* allcoupledcolnodes = dis.get_all_pbc_coupled_col_nodes();
 
   if (allcoupledcolnodes)
   {
@@ -432,8 +432,7 @@ std::shared_ptr<Core::LinAlg::MultiVector<double>> Core::FE::compute_superconver
         // this will result in off processor assembling (boundary node as ghost node)
         if (closestnode->owner() != myrank) continue;
 
-        std::map<int, std::vector<int>>::iterator masternode =
-            allcoupledcolnodes->find(closestnodeid);
+        auto masternode = allcoupledcolnodes->find(closestnodeid);
 
         int numslavenodes = -1;
         if (masternode != allcoupledcolnodes->end())

src/core/linear_solver/src/method/4C_linear_solver_method_parameters.cpp

@@ -86,9 +86,9 @@ void Core::LinearSolver::Parameters::compute_solver_parameters(
   {
     std::shared_ptr<Core::LinAlg::MultiVector<double>> coordinates;
     if (nullspace_node_map == nullptr)
-      coordinates = extract_node_coordinates(dis);
+      coordinates = extract_retained_node_coordinates(dis, *dis.node_row_map());
     else
-      coordinates = extract_node_coordinates(dis, *nullspace_node_map);
+      coordinates = extract_retained_node_coordinates(dis, *nullspace_node_map);
 
     solverlist.set<std::shared_ptr<Core::LinAlg::MultiVector<double>>>("Coordinates", coordinates);
   }

src/fluid_ele/4C_fluid_ele.hpp

@@ -868,7 +868,7 @@ namespace Discret
               lm_masterNodeToPatch,  ///< local map between master nodes and nodes in patch
           std::vector<int>&
               lm_slaveNodeToPatch,  ///< local map between slave nodes and nodes in patch
-          std::shared_ptr<std::map<int, int>>
+          std::shared_ptr<const std::map<int, int>>
               pbcconnectivity  ///< connectivity between slave and PBC's master nodes
       );
 
@@ -892,7 +892,7 @@ namespace Discret
               lm_masterNodeToPatch,  ///< local map between master nodes and nodes in patch
           std::vector<int>&
               lm_slaveNodeToPatch,  ///< local map between slave nodes and nodes in patch
-          std::shared_ptr<std::map<int, int>>
+          std::shared_ptr<const std::map<int, int>>
               pbcconnectivity  ///< connectivity between slave and PBC's master nodes
       );
 

src/fluid_ele/4C_fluid_ele_intfaces.cpp

@@ -122,7 +122,7 @@ void Discret::Elements::FluidIntFace::patch_location_vector(
     std::vector<int>& lm_faceToPatch,          ///< local map between lm_face and lm_patch
     std::vector<int>& lm_masterNodeToPatch,  ///< local map between master nodes and nodes in patch
     std::vector<int>& lm_slaveNodeToPatch,   ///< local map between slave nodes and nodes in patch
-    std::shared_ptr<std::map<int, int>>
+    std::shared_ptr<const std::map<int, int>>
         pbcconnectivity  ///< connectivity between slave and PBC's master nodes
 )
 {
@@ -185,7 +185,7 @@ void Discret::Elements::FluidIntFace::patch_location_vector(
 
     if (has_PBC)  // set the id of the master node if the node is a PBC node
     {
-      std::map<int, int>::iterator slave_it = pbcconnectivity->find(
+      const auto slave_it = pbcconnectivity->find(
           nid);  // find the slave node id, is there a corresponding pbc master node?
 
       if (slave_it != pbcconnectivity->end()) nid = slave_it->second;
@@ -222,7 +222,7 @@ void Discret::Elements::FluidIntFace::patch_location_vector(
 
     if (has_PBC)  // set the id of the master node if the node is a PBC node
     {
-      std::map<int, int>::iterator slave_it = pbcconnectivity->find(
+      const auto slave_it = pbcconnectivity->find(
           nid);  // find the slave node id, is there a corresponding pbc master node?
 
       if (slave_it != pbcconnectivity->end()) nid = slave_it->second;
@@ -287,7 +287,7 @@ void Discret::Elements::FluidIntFace::patch_location_vector(
 
     if (has_PBC)  // set the id of the master node if the node is a PBC node
     {
-      std::map<int, int>::iterator slave_it = pbcconnectivity->find(
+      const auto slave_it = pbcconnectivity->find(
           nid);  // find the slave node id, is there a corresponding pbc master node?
 
       if (slave_it != pbcconnectivity->end()) nid = slave_it->second;
@@ -334,7 +334,7 @@ void Discret::Elements::FluidIntFace::patch_location_vector(
     std::vector<int>& lm_faceToPatch,          ///< local map between lm_face and lm_patch
     std::vector<int>& lm_masterNodeToPatch,  ///< local map between master nodes and nodes in patch
     std::vector<int>& lm_slaveNodeToPatch,   ///< local map between slave nodes and nodes in patch
-    std::shared_ptr<std::map<int, int>>
+    std::shared_ptr<const std::map<int, int>>
         pbcconnectivity  ///< connectivity between slave and PBC's master nodes
 )
 {
@@ -411,7 +411,7 @@ void Discret::Elements::FluidIntFace::patch_location_vector(
 
     if (has_PBC)  // set the id of the master node if the node is a PBC node
     {
-      std::map<int, int>::iterator slave_it = pbcconnectivity->find(
+      const auto slave_it = pbcconnectivity->find(
           nid);  // find the slave node id, is there a corresponding pbc master node?
 
       if (slave_it != pbcconnectivity->end()) nid = slave_it->second;
@@ -453,7 +453,7 @@ void Discret::Elements::FluidIntFace::patch_location_vector(
 
     if (has_PBC)  // set the id of the master node if the node is a PBC node
     {
-      std::map<int, int>::iterator slave_it = pbcconnectivity->find(
+      const auto slave_it = pbcconnectivity->find(
           nid);  // find the slave node id, is there a corresponding pbc master node?
 
       if (slave_it != pbcconnectivity->end()) nid = slave_it->second;
@@ -523,7 +523,7 @@ void Discret::Elements::FluidIntFace::patch_location_vector(
 
     if (has_PBC)  // set the id of the master node if the node is a PBC node
     {
-      std::map<int, int>::iterator slave_it = pbcconnectivity->find(
+      const auto slave_it = pbcconnectivity->find(
           nid);  // find the slave node id, is there a corresponding pbc master node?
 
       if (slave_it != pbcconnectivity->end()) nid = slave_it->second;

src/fluid_ele/4C_fluid_ele_intfaces_calc.cpp

@@ -134,7 +134,7 @@ void Discret::Elements::FluidIntFaceImpl<distype>::assemble_internal_faces_using
 
 
   //---------------------- check for PBCS ------------------
-  std::shared_ptr<std::map<int, int>> pbcconnectivity =
+  std::shared_ptr<const std::map<int, int>> pbcconnectivity =
       discretization.get_pbc_slave_to_master_node_connectivity();
 
   //----------------------- create patchlm -----------------

src/fluid_turbulence/4C_fluid_turbulence_boxfilter.cpp

@@ -433,7 +433,7 @@ void FLD::Boxfilter::apply_box_filter(
     double expression_val;
     double alpha2_val;
 
-    std::map<int, std::vector<int>>* pbcmapmastertoslave =
+    const std::map<int, std::vector<int>>* pbcmapmastertoslave =
         discret_->get_all_pbc_coupled_col_nodes();
     // loop all master nodes on this proc
     if (pbcmapmastertoslave)
@@ -1170,7 +1170,7 @@ void FLD::Boxfilter::apply_box_filter_scatra(
     double phiexpression_val = 0.0;
 
     // loop all master nodes on this proc
-    std::map<int, std::vector<int>>* pbcmapmastertoslave =
+    const std::map<int, std::vector<int>>* pbcmapmastertoslave =
         scatradiscret_->get_all_pbc_coupled_col_nodes();
     // loop all master nodes on this proc
     if (pbcmapmastertoslave)

src/scatra/4C_scatra_timint_implicit_service.cpp

@@ -964,7 +964,7 @@ void ScaTra::ScaTraTimIntImpl::compute_null_space_if_necessary() const
       mllist.set<std::shared_ptr<Core::LinAlg::MultiVector<double>>>("nullspace", nullspace);
 
       std::shared_ptr<Core::LinAlg::MultiVector<double>> coordinates =
-          extract_node_coordinates(*discret_);
+          extract_retained_node_coordinates(*discret_, *discret_->node_row_map());
 
       mllist.set<std::shared_ptr<Core::LinAlg::MultiVector<double>>>("Coordinates", coordinates);
     }

src/scatra/4C_scatra_timint_meshtying_strategy_s2i.cpp

@@ -3661,8 +3661,9 @@ void ScaTra::MeshtyingStrategyS2I::equip_extended_solver_with_null_space_info()
     // build the coordinate object related to the growth dofs
     const auto growth_cond_node_row_map = Core::Conditions::condition_node_row_map(
         *scatratimint_->discretization(), "S2IKineticsGrowth");
-    const auto coordinates =
-        extract_node_coordinates(*scatratimint_->discretization(), *growth_cond_node_row_map);
+    const std::shared_ptr<Core::LinAlg::MultiVector<double>> coordinates =
+        extract_retained_node_coordinates(
+            *scatratimint_->discretization(), *growth_cond_node_row_map);
     mllist.set<std::shared_ptr<Core::LinAlg::MultiVector<double>>>("Coordinates", coordinates);
 
     mllist.set<std::shared_ptr<Core::LinAlg::MultiVector<double>>>("nullspace", nullspace);

src/sti/4C_sti_monolithic.cpp

@@ -1240,7 +1240,8 @@ void STI::Monolithic::compute_null_space_if_necessary(Teuchos::ParameterList& so
     mllist.set<std::shared_ptr<Core::LinAlg::MultiVector<double>>>("nullspace", nullspace);
 
     std::shared_ptr<Core::LinAlg::MultiVector<double>> coordinates =
-        extract_node_coordinates(*scatra_field()->discretization());
+        extract_retained_node_coordinates(
+            *scatra_field()->discretization(), *scatra_field()->discretization()->node_row_map());
 
     mllist.set<std::shared_ptr<Core::LinAlg::MultiVector<double>>>("Coordinates", coordinates);
   }