Ifpack Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Ifpack_ex_Reordering_LL.cpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // IFPACK
5 // Copyright (2004) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // This library is free software; you can redistribute it and/or modify
11 // it under the terms of the GNU Lesser General Public License as
12 // published by the Free Software Foundation; either version 2.1 of the
13 // License, or (at your option) any later version.
14 //
15 // This library is distributed in the hope that it will be useful, but
16 // WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 // Lesser General Public License for more details.
19 //
20 // You should have received a copy of the GNU Lesser General Public
21 // License along with this library; if not, write to the Free Software
22 // Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23 // USA
24 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
25 //
26 // ***********************************************************************
27 // @HEADER
28 
29 #include "Ifpack_ConfigDefs.h"
30 #ifdef HAVE_MPI
31 #include "Epetra_MpiComm.h"
32 #else
33 #include "Epetra_SerialComm.h"
34 #endif
35 #include "Epetra_Map.h"
36 #include "Epetra_CrsMatrix.h"
37 #include "Ifpack_Reordering.h"
38 #include "Ifpack_RCMReordering.h"
39 #include "Ifpack_ReorderFilter.h"
40 #include "Ifpack_Utils.h"
41 #include "Teuchos_RefCountPtr.hpp"
42 
43 //==============================================================================
44 int main(int argc, char *argv[])
45 {
46  using std::cout;
47  using std::endl;
48 
49 #ifdef HAVE_MPI
50  MPI_Init(&argc,&argv);
51  Epetra_MpiComm Comm(MPI_COMM_WORLD);
52 #else
53  Epetra_SerialComm Comm;
54 #endif
55 
56  // only one process
57  if (Comm.NumProc() != 1) {
58 #ifdef HAVE_MPI
59  MPI_Finalize();
60 #endif
61  if (Comm.MyPID() == 0)
62  cout << "Please run this test with one process only" << endl;
63  // return success not to break the tests
64  exit(EXIT_SUCCESS);
65  }
66 
67  // ======================================================== //
68  // now create the famous "upper arrow" matrix, which //
69  // should be reordered as a "lower arrow". Sparsity pattern //
70  // will be printed on screen. //
71  // ======================================================== //
72 
73  int NumPoints = 16;
74 
75 #if !defined(EPETRA_NO_32BIT_GLOBAL_INDICES) || !defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
76  Epetra_Map Map(-1LL,NumPoints,0LL,Comm);
77 #else
78  Epetra_Map Map;
79 #endif
80 
81  std::vector<long long> Indices(NumPoints);
82  std::vector<double> Values(NumPoints);
83 
84  Teuchos::RefCountPtr<Epetra_CrsMatrix> A = Teuchos::rcp( new Epetra_CrsMatrix(Copy,Map,0) );
85  for (int i = 0 ; i < NumPoints ; ++i) {
86 
87  int NumEntries;
88  if (i == 0) {
89  NumEntries = NumPoints;
90  for (int j = 0 ; j < NumPoints ; ++j) {
91  Indices[j] = j;
92  Values[j] = 1.0;
93  }
94  }
95  else {
96  NumEntries = 2;
97  Indices[0] = 0;
98  Indices[1] = i;
99  Values[0] = 1.0;
100  Values[1] = 1.0;
101  }
102 
103 #if !defined(EPETRA_NO_32BIT_GLOBAL_INDICES) || !defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
104  A->InsertGlobalValues(i, NumEntries, &Values[0], &Indices[0]);
105 #endif
106  }
107 
108  A->FillComplete();
109 
110  // print the sparsity to file, postscript format
112 
113  // create the reordering...
114  Teuchos::RefCountPtr<Ifpack_RCMReordering> Reorder = Teuchos::rcp( new Ifpack_RCMReordering() );
115  // and compute is on A
116  IFPACK_CHK_ERR(Reorder->Compute(*A));
117 
118  // cout information
119  cout << *Reorder;
120 
121  // create a reordered matrix
122  Ifpack_ReorderFilter ReordA(A, Reorder);
123 
124  // print the sparsity to file, postscript format
126 
127 #ifdef HAVE_MPI
128  MPI_Finalize();
129 #endif
130  return(EXIT_SUCCESS);
131 
132 }
int MyPID() const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Ifpack_RCMReordering: reverse Cuthill-McKee reordering.
int main(int argc, char *argv[])
int NumProc() const
#define IFPACK_CHK_ERR(ifpack_err)
Ifpack_ReorderFilter: a class for light-weight reorder of local rows and columns of an Epetra_RowMatr...