Zoltan2
GraphModel2ndAdjsFromAdjs.cpp
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 //
46 // Basic testing of Zoltan2::PamgenMeshAdapter
47 
48 #include <Zoltan2_GraphModel.hpp>
49 #include <Zoltan2_ModelHelpers.hpp>
51 
52 // Teuchos includes
53 #include "Teuchos_XMLParameterListHelpers.hpp"
54 
55 // Pamgen includes
56 #include "create_inline_mesh.h"
57 
58 using namespace std;
59 using Teuchos::RCP;
60 
61 /*********************************************************/
62 /* Typedefs */
63 /*********************************************************/
64 //Tpetra typedefs
65 typedef Tpetra::DefaultPlatform::DefaultPlatformType Platform;
66 typedef Tpetra::MultiVector<double> tMVector_t;
67 
68 
69 
70 /*****************************************************************************/
71 /******************************** MAIN ***************************************/
72 /*****************************************************************************/
73 
74 int main(int narg, char *arg[]) {
75 
76  Teuchos::GlobalMPISession mpiSession(&narg, &arg,0);
77  Platform &platform = Tpetra::DefaultPlatform::getDefaultPlatform();
78  RCP<const Teuchos::Comm<int> > CommT = platform.getComm();
79 
80  int me = CommT->getRank();
81  int numProcs = CommT->getSize();
82 
83  /***************************************************************************/
84  /*************************** GET XML INPUTS ********************************/
85  /***************************************************************************/
86 
87  // default values for command-line arguments
88  std::string xmlMeshInFileName("Poisson.xml");
89 
90  // Read run-time options.
91  Teuchos::CommandLineProcessor cmdp (false, false);
92  cmdp.setOption("xmlfile", &xmlMeshInFileName,
93  "XML file with PamGen specifications");
94  cmdp.parse(narg, arg);
95 
96  // Read xml file into parameter list
97  ParameterList inputMeshList;
98 
99  if(xmlMeshInFileName.length()) {
100  if (me == 0) {
101  cout << "\nReading parameter list from the XML file \""
102  <<xmlMeshInFileName<<"\" ...\n\n";
103  }
104  Teuchos::updateParametersFromXmlFile(xmlMeshInFileName,
105  Teuchos::inoutArg(inputMeshList));
106  if (me == 0) {
107  inputMeshList.print(cout,2,true,true);
108  cout << "\n";
109  }
110  }
111  else {
112  cout << "Cannot read input file: " << xmlMeshInFileName << "\n";
113  return 5;
114  }
115 
116  // Get pamgen mesh definition
117  std::string meshInput = Teuchos::getParameter<std::string>(inputMeshList,
118  "meshInput");
119 
120  /***************************************************************************/
121  /********************** GET CELL TOPOLOGY **********************************/
122  /***************************************************************************/
123 
124  // Get dimensions
125  int dim = 3;
126 
127  /***************************************************************************/
128  /***************************** GENERATE MESH *******************************/
129  /***************************************************************************/
130 
131  if (me == 0) cout << "Generating mesh ... \n\n";
132 
133  // Generate mesh with Pamgen
134  long long maxInt = std::numeric_limits<long long>::max();
135  Create_Pamgen_Mesh(meshInput.c_str(), dim, me, numProcs, maxInt);
136 
137  // Creating mesh adapter
138  if (me == 0) cout << "Creating mesh adapter ... \n\n";
139 
140  typedef Zoltan2::PamgenMeshAdapter<tMVector_t> inputAdapter_t;
142 
143  inputAdapter_t ia(*CommT, "region");
144  inputAdapter_t ia2(*CommT, "vertex");
145  inputAdapter_t::gno_t const *adjacencyIds=NULL;
146  inputAdapter_t::lno_t const *offsets=NULL;
147  Teuchos::ArrayRCP<inputAdapter_t::lno_t> moffsets;
148  Teuchos::ArrayRCP<inputAdapter_t::gno_t> madjacencyIds;
149  ia.print(me);
150 
151  if (me == 0) std::cout << "REGION-BASED TEST" << std::endl;
152  Zoltan2::MeshEntityType primaryEType = ia.getPrimaryEntityType();
153  Zoltan2::MeshEntityType adjEType = ia.getAdjacencyEntityType();
154  Zoltan2::MeshEntityType secondAdjEType = ia.getSecondAdjacencyEntityType();
155  RCP<const base_adapter_t> baseInputAdapter;
156  RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment);
157  std::bitset<Zoltan2::NUM_MODEL_FLAGS> modelFlags;
158 
159  if (ia.availAdjs(primaryEType, adjEType)) {
160  // Create a GraphModel based on this input data.
161 
162  if (me == 0) std::cout << " Creating GraphModel" << std::endl;
163 
164  baseInputAdapter = (rcp(dynamic_cast<const base_adapter_t *>(&ia), false));
165 
166  Zoltan2::GraphModel<base_adapter_t> graphModel(baseInputAdapter, env,
167  CommT, modelFlags);
168 
169  if (me == 0)
170  std::cout << " Calling get2ndAdjsViewFromAdjs" << std::endl;
171  Zoltan2::get2ndAdjsViewFromAdjs(baseInputAdapter, graphModel.getComm(),
172  primaryEType,
173  secondAdjEType, moffsets, madjacencyIds);
174 
175  if (me == 0) std::cout << " Checking results" << std::endl;
176  if (ia.avail2ndAdjs(primaryEType, secondAdjEType)) {
177  ia.get2ndAdjsView(primaryEType, secondAdjEType, offsets, adjacencyIds);
178  }
179  else{
180  std::cout << "2nd adjacencies not available; cannot check results"
181  << std::endl;
182  return 2;
183  }
184 
185  for (size_t telct = 0; telct < ia.getLocalNumOf(primaryEType); telct++) {
186  if (offsets[telct+1]-offsets[telct]!=moffsets[telct+1]-moffsets[telct]) {
187  std::cout << "Number of adjacencies do not match" << std::endl;
188  return 3;
189  }
190 
191  for (inputAdapter_t::lno_t j=moffsets[telct]; j<moffsets[telct+1]; j++) {
192  ssize_t in_list = -1;
193 
194  for (inputAdapter_t::lno_t k=offsets[telct]; k<offsets[telct+1]; k++) {
195  if (adjacencyIds[k] == madjacencyIds[j]) {
196  in_list = k;
197  break;
198  }
199  }
200 
201  if (in_list < 0) {
202  std::cout << "Adjacency missing" << std::endl;
203  return 4;
204  }
205  }
206  }
207  }
208  else{
209  std::cout << "Adjacencies not available" << std::endl;
210  return 1;
211  }
212 
213  if (me == 0) std::cout << "VERTEX-BASED TEST" << std::endl;
214  primaryEType = ia2.getPrimaryEntityType();
215  adjEType = ia2.getAdjacencyEntityType();
216  secondAdjEType = ia2.getSecondAdjacencyEntityType();
217 
218  if (ia2.availAdjs(primaryEType, adjEType)) {
219  if (ia2.avail2ndAdjs(primaryEType, secondAdjEType)) {
220  ia2.get2ndAdjsView(primaryEType, secondAdjEType, offsets, adjacencyIds);
221  }
222  else{
223  std::cout << "2nd adjacencies not available" << std::endl;
224  return 2;
225  }
226 
227  // Create a GraphModel based on this input data.
228 
229  if (me == 0) std::cout << " Creating GraphModel" << std::endl;
230 
231  baseInputAdapter = (rcp(dynamic_cast<const base_adapter_t *>(&ia2),false));
232 
233  Zoltan2::GraphModel<base_adapter_t> graphModel2(baseInputAdapter, env,
234  CommT, modelFlags);
235 
236  if (me == 0)
237  std::cout << " Calling get2ndAdjsViewFromAdjs" << std::endl;
238  Zoltan2::get2ndAdjsViewFromAdjs(baseInputAdapter, graphModel2.getComm(),
239  primaryEType,
240  secondAdjEType, moffsets, madjacencyIds);
241 
242  if (me == 0) std::cout << " Checking results" << std::endl;
243 
244  for (size_t tnoct = 0; tnoct < ia2.getLocalNumOf(primaryEType); tnoct++) {
245  if (offsets[tnoct+1]-offsets[tnoct]!=moffsets[tnoct+1]-moffsets[tnoct]) {
246  std::cout << "Number of adjacencies do not match" << std::endl;
247  return 3;
248  }
249 
250  for (inputAdapter_t::lno_t j=moffsets[tnoct]; j<moffsets[tnoct+1]; j++) {
251  ssize_t in_list = -1;
252 
253  for (inputAdapter_t::lno_t k=offsets[tnoct]; k<offsets[tnoct+1]; k++) {
254  if (adjacencyIds[k] == madjacencyIds[j]) {
255  in_list = k;
256  break;
257  }
258  }
259 
260  if (in_list < 0) {
261  std::cout << "Adjacency missing" << std::endl;
262  return 4;
263  }
264  }
265  }
266  }
267  else{
268  std::cout << "Adjacencies not available" << std::endl;
269  return 1;
270  }
271 
272  // delete mesh
273  if (me == 0) cout << "Deleting the mesh ... \n\n";
274 
275  Delete_Pamgen_Mesh();
276 
277  if (me == 0)
278  std::cout << "PASS" << std::endl;
279 
280  return 0;
281 }
282 /*****************************************************************************/
283 /********************************* END MAIN **********************************/
284 /*****************************************************************************/
Zoltan2::BaseAdapter< userTypes_t > base_adapter_t
void get2ndAdjsViewFromAdjs(const Teuchos::RCP< const MeshAdapter< User > > &ia, const RCP< const Comm< int > > comm, Zoltan2::MeshEntityType sourcetarget, Zoltan2::MeshEntityType through, Teuchos::ArrayRCP< typename MeshAdapter< User >::lno_t > &offsets, Teuchos::ArrayRCP< typename MeshAdapter< User >::gno_t > &adjacencyIds)
Defines the PamgenMeshAdapter class.
Defines helper functions for use in the models.
Tpetra::DefaultPlatform::DefaultPlatformType Platform
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
int main(int narg, char *arg[])
GraphModel defines the interface required for graph models.
Tpetra::DefaultPlatform::DefaultPlatformType Platform
MeshEntityType
Enumerate entity types for meshes: Regions, Faces, Edges, or Vertices.
Defines the GraphModel interface.
const RCP< const Comm< int > > getComm()
Return the communicator used by the model.
This class represents a mesh.
Tpetra::MultiVector< double > tMVector_t