libMesh
sfc_partitioner.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2024 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/libmesh_config.h"
22 #include "libmesh/elem.h"
23 #include "libmesh/enum_partitioner_type.h"
24 #include "libmesh/libmesh_logging.h"
25 #include "libmesh/mesh_base.h"
26 #include "libmesh/sfc_partitioner.h"
27 
28 #ifdef LIBMESH_HAVE_SFCURVES
29 namespace Sfc {
30 extern "C" {
31 # include "sfcurves.h"
32 }
33 }
34 #else
35 # include "libmesh/linear_partitioner.h"
36 #endif
37 
38 namespace libMesh
39 {
40 
41 
43 {
44  return SFC_PARTITIONER;
45 }
46 
47 
51  unsigned int n)
52 {
53  // Check for easy returns
54  if (beg == end)
55  return;
56 
57  if (n == 1)
58  {
59  this->single_partition_range (beg, end);
60  return;
61  }
62 
63  libmesh_assert_greater (n, 0);
64 
65  // What to do if the sfcurves library IS NOT present
66 #ifndef LIBMESH_HAVE_SFCURVES
67 
68  libmesh_do_once(
69  libMesh::out << "ERROR: The library has been built without" << std::endl
70  << "Space Filling Curve support. Using a linear" << std::endl
71  << "partitioner instead!" << std::endl;);
72 
74  lp.partition_range (mesh, beg, end, n);
75 
76  // What to do if the sfcurves library IS present
77 #else
78 
79  LOG_SCOPE("partition_range()", "SFCPartitioner");
80 
81  // We don't yet support distributed meshes with this Partitioner
82  if (!mesh.is_serial())
83  libmesh_not_implemented();
84 
85  const dof_id_type n_range_elem = std::distance(beg, end);
86  const dof_id_type max_elem_id = mesh.max_elem_id();
87 
88  // The forward_map maps the range's element ids into a contiguous
89  // block of indices.
90  std::vector<dof_id_type> forward_map (max_elem_id, DofObject::invalid_id);
91 
92  // the reverse_map maps the contiguous ids back
93  // to active elements
94  std::vector<Elem *> reverse_map (n_range_elem, nullptr);
95 
96  std::vector<double> x (n_range_elem);
97  std::vector<double> y (n_range_elem);
98  std::vector<double> z (n_range_elem);
99  std::vector<int> table (n_range_elem);
100 
101  // Map the range's element ids into a contiguous range.
102  dof_id_type el_num = 0;
103 
104  for (auto & elem : as_range(beg, end))
105  {
106  libmesh_assert_less (elem->id(), forward_map.size());
107  libmesh_assert_less (el_num, reverse_map.size());
108 
109  forward_map[elem->id()] = el_num;
110  reverse_map[el_num] = elem;
111  el_num++;
112  }
113  libmesh_assert_equal_to (el_num, n_range_elem);
114 
115  // Get the vertex average for each range element.
116  for (const auto & elem : as_range(beg, end))
117  {
118  libmesh_assert_less (elem->id(), forward_map.size());
119 
120  const Point p = elem->vertex_average();
121 
122  x[forward_map[elem->id()]] = double(p(0));
123  y[forward_map[elem->id()]] = double(p(1));
124  z[forward_map[elem->id()]] = double(p(2));
125  }
126 
127  // We need an integer reference to pass to the Sfc interface.
128  int size = static_cast<int>(n_range_elem);
129 
130  // build the space-filling curve
131  if (_sfc_type == "Hilbert")
132  Sfc::hilbert (x.data(),
133  y.data(),
134  z.data(),
135  &size,
136  table.data());
137 
138  else if (_sfc_type == "Morton")
139  Sfc::morton (x.data(),
140  y.data(),
141  z.data(),
142  &size,
143  table.data());
144 
145  else
146  {
147  libMesh::out << "ERROR: Unknown type: " << _sfc_type << std::endl
148  << " Valid types are" << std::endl
149  << " \"Hilbert\"" << std::endl
150  << " \"Morton\"" << std::endl
151  << " " << std::endl
152  << "Partitioning with a Hilbert curve." << std::endl;
153 
154  Sfc::hilbert (x.data(),
155  y.data(),
156  z.data(),
157  &size,
158  table.data());
159  }
160 
161 
162  // Assign the partitioning to the range elements
163  {
164  // {
165  // std::ofstream out ("sfc.dat");
166  // out << "variables=x,y,z" << std::endl;
167  // out << "zone f=point" << std::endl;
168  // for (unsigned int i=0; i<n_range_elem; i++)
169  // out << x[i] << " " << y[i] << " " << z[i] << std::endl;
170  // }
171 
172  const dof_id_type blksize = (n_range_elem + n - 1) / n;
173 
174  for (dof_id_type i=0; i<n_range_elem; i++)
175  {
176  libmesh_assert_less (table[i] - 1, reverse_map.size());
177 
178  Elem * elem = reverse_map[table[i] - 1];
179 
180  elem->processor_id() = cast_int<processor_id_type>(i/blksize);
181  }
182  }
183 
184 #endif
185 }
186 
187 
188 
190  const unsigned int n)
191 {
192  this->partition_range(mesh,
193  mesh.active_elements_begin(),
194  mesh.active_elements_end(),
195  n);
196 }
197 
198 } // namespace libMesh
virtual void partition_range(MeshBase &mesh, MeshBase::element_iterator it, MeshBase::element_iterator end, const unsigned int n) override
Called by the SubdomainPartitioner to partition elements in the range (it, end).
The definition of the element_iterator struct.
Definition: mesh_base.h:2103
virtual void _do_partition(MeshBase &mesh, const unsigned int n) override
Partition the MeshBase into n subdomains.
bool single_partition_range(MeshBase::element_iterator it, MeshBase::element_iterator end)
Slightly generalized version of single_partition which acts on a range of elements defined by the pai...
Definition: partitioner.C:319
The LinearPartitioner simply takes the element list and splits it into equal-sized chunks assigned to...
std::string _sfc_type
The type of space-filling curve to use.
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
MeshBase & mesh
The libMesh namespace provides an interface to certain functionality in the library.
Real distance(const Point &p)
This is the MeshBase class.
Definition: mesh_base.h:74
PartitionerType
Defines an enum for mesh partitioner types.
SimpleRange< IndexType > as_range(const std::pair< IndexType, IndexType > &p)
Helper function that allows us to treat a homogenous pair as a range.
Definition: simple_range.h:57
virtual PartitionerType type() const override
static const dof_id_type invalid_id
An invalid id to distinguish an uninitialized DofObject.
Definition: dof_object.h:477
OStreamProxy out
virtual void partition_range(MeshBase &mesh, MeshBase::element_iterator it, MeshBase::element_iterator end, const unsigned int n) override
Called by the SubdomainPartitioner to partition elements in the range (it, end).
processor_id_type processor_id() const
Definition: dof_object.h:898
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
uint8_t dof_id_type
Definition: id_types.h:67