Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Zoltan2_Adapter.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // Zoltan2: A package of combinatorial algorithms for scientific computing
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Karen Devine (kddevin@sandia.gov)
39 // Erik Boman (egboman@sandia.gov)
40 // Siva Rajamanickam (srajama@sandia.gov)
41 //
42 // ***********************************************************************
43 //
44 // @HEADER
45 
50 #ifndef _ZOLTAN2_ADAPTER_HPP_
51 #define _ZOLTAN2_ADAPTER_HPP_
52 
53 #include <Kokkos_Core.hpp>
54 #include <Zoltan2_Standards.hpp>
55 #include <Zoltan2_InputTraits.hpp>
57 
58 namespace Zoltan2 {
59 
70 };
71 
82 public:
83  virtual ~BaseAdapterRoot() {}; // required virtual declaration
84 
90  virtual size_t getLocalNumIDs() const = 0;
91 
97  virtual int getNumWeightsPerID() const { return 0; };
98 };
99 
100 template <typename User>
101  class BaseAdapter : public BaseAdapterRoot {
102 
103 public:
104  typedef typename InputTraits<User>::lno_t lno_t;
105  typedef typename InputTraits<User>::gno_t gno_t;
109 
112  virtual enum BaseAdapterType adapterType() const = 0;
113 
116  virtual ~BaseAdapter() {};
117 
123  virtual void getIDsView(const gno_t *&ids) const {
124  Kokkos::View<gno_t *> kokkosIds;
125  getIDsKokkosView(kokkosIds);
126  ids = kokkosIds.data();
127  }
128 
134  virtual void getIDsKokkosView(Kokkos::View<gno_t *> &/* ids */) const {
136  }
137 
139  // * \param wgt on return a pointer to the weights for this idx
140  // * \param stride on return, the value such that
141  // * the \t nth weight should be found at <tt> wgt[n*stride] </tt>.
142  // * \param idx the weight index, zero or greater
143  // * This function must be implemented in derived adapter if
144  // * getNumWeightsPerID > 0.
145  // * This function should not be called if getNumWeightsPerID is zero.
146  // */
147  virtual void getWeightsView(const scalar_t *&wgt, int &stride,
148  int idx = 0) const {
149  Kokkos::View<scalar_t *> tempWeightsView;
150  getWeightsKokkosView(tempWeightsView, idx);
151  wgt = tempWeightsView.data();
152  stride = 1;
153  }
154 
156  // * \param wgt on return a pointer to the weights for this idx
157  // * \param idx the weight index, zero or greater
158  // * This function must be implemented in derived adapter if
159  // * getNumWeightsPerID > 0.
160  // * This function should not be called if getNumWeightsPerID is zero.
161  // */
162  virtual void getWeightsKokkosView(Kokkos::View<scalar_t *> &/* wgt */,
163  int /* idx */ = 0) const {
165  }
166 
176  void getPartsView(const part_t *&inputPart) const {
177  // Default behavior: return NULL for inputPart array;
178  // assume input part == rank
179  inputPart = NULL;
180  }
181 
199  template <typename Adapter>
200  void applyPartitioningSolution(const User &in, User *&out,
201  const PartitioningSolution<Adapter> &solution) const {
203  }
204 
205 protected:
206 
207  // Write Chaco-formatted graph and assign files echoing adapter input
208  // This routine is serial and may be slow; use it only for debugging
209  // This function does not write edge info to the graph file, as the
210  // BaseAdapter does not know about edge info; it writes
211  // only the Chaco header and vertex weights (if applicable).
212  void generateWeightFileOnly(const char* fileprefix,
213  const Teuchos::Comm<int> &comm) const;
214 
215 };
216 
217 template <typename User>
219  const char *fileprefix,
220  const Teuchos::Comm<int> &comm
221 ) const
222 {
223  int np = comm.getSize();
224  int me = comm.getRank();
225 
226  size_t nLocalIDs = this->getLocalNumIDs();
227 
228  // Write .graph file: header and weights only (no edges)
229  // Adapters with edges have to implement their own generateFiles function
230  // to provide edge info.
231  {
232  // append suffix to filename
233  std::string filenamestr = fileprefix;
234  filenamestr = filenamestr + ".graph";
235  const char *filename = filenamestr.c_str();
236 
237  size_t nGlobalIDs;
238  Teuchos::reduceAll(comm, Teuchos::REDUCE_SUM, 1, &nLocalIDs, &nGlobalIDs);
239 
240  int nWgts = this->getNumWeightsPerID();
241 
242  for (int p = 0; p < np; p++) {
243 
244  // Is it this processor's turn to write to files?
245  if (me == p) {
246 
247  std::ofstream fp;
248 
249  if (me == 0) {
250  // open file for writing
251  fp.open(filename, std::ios::out);
252  // write Chaco header info
253  // this function assumes no edges
254  fp << nGlobalIDs << " " << 0 << " "
255  << (nWgts ? "010" : "000") << " "
256  << (nWgts > 1 ? std::to_string(nWgts) : " ") << std::endl;
257  }
258  else {
259  // open file for appending
260  fp.open(filename, std::ios::app);
261  }
262 
263  if (nWgts) {
264 
265  // get weight data
266  const scalar_t **wgts = new const scalar_t *[nWgts];
267  int *strides = new int[nWgts];
268  for (int n = 0; n < nWgts; n++)
269  getWeightsView(wgts[n], strides[n], n);
270 
271  // write weights to file
272  for (size_t i = 0; i < nLocalIDs; i++) {
273  for (int n = 0; n < nWgts; n++)
274  fp << wgts[n][i*strides[n]] << " ";
275  fp << "\n";
276  }
277 
278  delete [] strides;
279  delete [] wgts;
280  }
281 
282  fp.close();
283  }
284 
285  comm.barrier();
286  }
287  }
288 
289  // write assignments file
290  {
291  std::string filenamestr = fileprefix;
292  filenamestr = filenamestr + ".assign";
293  const char *filename = filenamestr.c_str();
294 
295  for (int p = 0; p < np; p++) {
296 
297  // Is it this processor's turn to write to files?
298  if (me == p) {
299 
300  std::ofstream fp;
301 
302  if (me == 0) {
303  // open file for writing
304  fp.open(filename, std::ios::out);
305  }
306  else {
307  // open file for appending
308  fp.open(filename, std::ios::app);
309  }
310 
311  const part_t *parts;
312  this->getPartsView(parts);
313 
314  for (size_t i = 0; i < nLocalIDs; i++) {
315  fp << (parts != NULL ? parts[i] : me) << "\n";
316  }
317  fp.close();
318  }
319 
320  comm.barrier();
321  }
322  }
323 }
324 
325 } //namespace Zoltan2
326 
327 #endif
InputTraits< User >::scalar_t scalar_t
virtual int getNumWeightsPerID() const
Returns the number of weights per object. Number of weights per object should be zero or greater...
virtual void getIDsView(const gno_t *&ids) const
Provide a pointer to this process&#39; identifiers.
InputTraits< User >::gno_t gno_t
default_part_t part_t
The data type to represent part numbers.
default_offset_t offset_t
The data type to represent offsets.
InputTraits< User >::lno_t lno_t
#define Z2_THROW_NOT_IMPLEMENTED
Defines the PartitioningSolution class.
virtual void getWeightsView(const scalar_t *&wgt, int &stride, int idx=0) const
Provide pointer to a weight array with stride.
virtual ~BaseAdapter()
Destructor.
A PartitioningSolution is a solution to a partitioning problem.
void generateWeightFileOnly(const char *fileprefix, const Teuchos::Comm< int > &comm) const
default_lno_t lno_t
The ordinal type (e.g., int, long, int64_t) that represents local counts and local indices...
BaseAdapterType
An enum to identify general types of adapters.
identifier data, just a list of IDs
BaseAdapter defines methods required by all Adapters.
Traits for application input objects.
InputTraits< User >::part_t part_t
void getPartsView(const part_t *&inputPart) const
Provide pointer to an array containing the input part assignment for each ID. The input part informat...
default_gno_t gno_t
The ordinal type (e.g., int, long, int64_t) that can represent global counts and identifiers.
void applyPartitioningSolution(const User &in, User *&out, const PartitioningSolution< Adapter > &solution) const
Apply a PartitioningSolution to an input.
virtual size_t getLocalNumIDs() const =0
Returns the number of objects on this process.
virtual void getWeightsKokkosView(Kokkos::View< scalar_t * > &, int=0) const
Provide pointer to a weight View.
Gathering definitions used in software development.
virtual void getIDsKokkosView(Kokkos::View< gno_t * > &) const
Provide a pointer to this process&#39; identifiers.
InputTraits< User >::offset_t offset_t
default_scalar_t scalar_t
The data type for weights and coordinates.
virtual enum BaseAdapterType adapterType() const =0
Returns the type of adapter.