sparsity_pattern.h

Go to the documentation of this file.

// The libMesh Finite Element Library.
// Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner

// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.

// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.

// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA


#ifndef LIBMESH_SPARSITY_PATTERN_H
#define LIBMESH_SPARSITY_PATTERN_H

// Local Includes
#include "libmesh/elem_range.h"
#include "libmesh/threads_allocators.h"
#include "libmesh/parallel_object.h"

// C++ includes
#include <algorithm> // is_sorted
#include <unordered_set>
#include <vector>

namespace libMesh
{

// Forward declarations
class DofMap;
class CouplingMatrix;
class StaticCondensationDofMap;

namespace SparsityPattern // use a namespace so member classes can be forward-declared.
{
typedef std::vector<dof_id_type, Threads::scalable_allocator<dof_id_type>> Row;
class Graph : public std::vector<Row> {};

class NonlocalGraph : public std::map<dof_id_type, Row> {};

template<typename BidirectionalIterator>
static void sort_row (const BidirectionalIterator begin,
                      BidirectionalIterator       middle,
                      const BidirectionalIterator end);



  class AugmentSparsityPattern
  {
  public:
    virtual ~AugmentSparsityPattern () = default;

    virtual void augment_sparsity_pattern (SparsityPattern::Graph & sparsity,
                                           std::vector<dof_id_type> & n_nz,
                                           std::vector<dof_id_type> & n_oz) = 0;
  };



class Build : public ParallelObject
{
public:
  Build (const DofMap & dof_map_in,
         const CouplingMatrix * dof_coupling_in,
         const std::set<GhostingFunctor *> & coupling_functors_in,
         const bool implicit_neighbor_dofs_in,
         const bool need_full_sparsity_pattern_in,
         const bool calculate_constrained_in = false,
         const StaticCondensationDofMap * sc = nullptr);

  Build (const Build &) = default;
  Build & operator= (const Build &) = delete;
  Build (Build &&) = default;
  Build & operator= (Build &&) = delete;
  ~Build () = default;

  Build (Build & other, Threads::split);

  void operator()(const ConstElemRange & range);

  void join (const Build & other);

  void parallel_sync ();

  const SparsityPattern::Graph & get_sparsity_pattern() const
  { return sparsity_pattern; }

  const SparsityPattern::NonlocalGraph & get_nonlocal_pattern() const
  { return nonlocal_pattern; }

  std::size_t n_nonzeros() const;

  const std::vector<dof_id_type> & get_n_nz() const
  { return n_nz; }

  const std::vector<dof_id_type> & get_n_oz() const
  { return n_oz; }

  void apply_extra_sparsity_object(SparsityPattern::AugmentSparsityPattern & asp);

  void apply_extra_sparsity_function(void (*func)(SparsityPattern::Graph & sparsity,
                                                   std::vector<dof_id_type> & n_nz,
                                                   std::vector<dof_id_type> & n_oz,
                                                   void * context),
                                      void * context)
  { func(sparsity_pattern, n_nz, n_oz, context); }

  void clear_full_sparsity()
  {
    sparsity_pattern.clear();
    nonlocal_pattern.clear();
  }

private:
  const DofMap & dof_map;
  const CouplingMatrix * dof_coupling;
  const std::set<GhostingFunctor *> & coupling_functors;
  const bool implicit_neighbor_dofs;
  const bool need_full_sparsity_pattern;
  const bool calculate_constrained;
  const StaticCondensationDofMap * const sc;

  std::unordered_set<dof_id_type> hashed_dof_sets;

  std::vector<dof_id_type> dummy_vec;

  void handle_vi_vj(const std::vector<dof_id_type> & element_dofs_i,
                    const std::vector<dof_id_type> & element_dofs_j);

  void sorted_connected_dofs(const Elem * elem,
                             std::vector<dof_id_type> & dofs_vi,
                             unsigned int vi);

#ifndef LIBMESH_ENABLE_DEPRECATED
private:
#endif

  SparsityPattern::Graph sparsity_pattern;

  SparsityPattern::NonlocalGraph nonlocal_pattern;

  std::vector<dof_id_type> n_nz;

  std::vector<dof_id_type> n_oz;
};

}



// ------------------------------------------------------------
// SparsityPattern inline member functions
template<typename BidirectionalIterator>
inline
void SparsityPattern::sort_row (const BidirectionalIterator begin,
                                BidirectionalIterator       middle,
                                const BidirectionalIterator end)
{
  // Assure we have the conditions for an inplace_merge
#ifdef DEBUG
  libmesh_assert(std::is_sorted(begin, middle));
  libmesh_assert(std::is_sorted(middle, end));
#endif
  libmesh_assert(std::unique(begin, middle) == middle);
  libmesh_assert(std::unique(middle, end) == end);

  std::inplace_merge(begin, middle, end);

  // Assure the algorithm worked if we are in DEBUG mode
#ifdef DEBUG
  libmesh_assert (std::is_sorted(begin,end));
#endif

  // Make sure the two ranges did not contain any common elements
  libmesh_assert (std::unique (begin, end) == end);
}

} // namespace libMesh

#endif // LIBMESH_SPARSITY_PATTERN_H