EpetraExt Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
test/Copy/cxx_main.cpp
Go to the documentation of this file.
1//@HEADER
2// ***********************************************************************
3//
4// EpetraExt: Epetra Extended - Linear Algebra Services Package
5// Copyright (2011) Sandia Corporation
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8// the U.S. Government retains certain rights in this software.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ***********************************************************************
40//@HEADER
41
42#include "Epetra_Map.h"
43#include "Epetra_Time.h"
44#include "Epetra_CrsMatrix.h"
45#include "Epetra_Vector.h"
47#include "Epetra_Flops.h"
48#ifdef EPETRA_MPI
49#include "Epetra_MpiComm.h"
50#include "mpi.h"
51#else
52#include "Epetra_SerialComm.h"
53#endif
54#include "Epetra_Version.h"
55#include "Trilinos_Util_CrsMatrixGallery.h"
57using namespace Trilinos_Util;
58
59// prototypes
60
61
62int main(int argc, char *argv[])
63{
64 int ierr = 0, i, forierr = 0;
65 bool debug = false;
66
67#ifdef EPETRA_MPI
68
69 // Initialize MPI
70
71 MPI_Init(&argc,&argv);
72 int size, rank; // Number of MPI processes, My process ID
73
74 MPI_Comm_size(MPI_COMM_WORLD, &size);
75 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
76 Epetra_MpiComm Comm( MPI_COMM_WORLD );
77
78#else
79
80 int size = 1; // Serial case (not using MPI)
81 int rank = 0;
83
84#endif
85
86 bool verbose = false;
87
88 // Check if we should print results to standard out
89 if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
90
91 // char tmp;
92 // if (rank==0) cout << "Press any key to continue..."<< endl;
93 // if (rank==0) cin >> tmp;
94 // Comm.Barrier();
95
96 int MyPID = Comm.MyPID();
97 int NumProc = Comm.NumProc();
98
99 if(verbose && MyPID==0)
100 cout << Epetra_Version() << endl << endl;
101
102 if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
103 << " is alive."<<endl;
104
105// Generate Laplace matrix with 4 rows on each processor
106 CrsMatrixGallery laplace_2d("laplace_2d", Comm);
107 laplace_2d.Set("problem_size",Comm.NumProc()*Comm.NumProc()*4);
108 Epetra_CrsMatrix * laplace_2d_matrix = laplace_2d.GetMatrix();
109 if (verbose) cout << "Orig matrix = " << *laplace_2d_matrix << endl;
110 const Epetra_Map * origRowMap = laplace_2d.GetMap();
111 int * origGids = origRowMap->MyGlobalElements();
112
113 // Next generate a sub map of the original problem, taking every other GID.
114 Epetra_IntSerialDenseVector newRowMapGids(origRowMap->NumMyElements()/2 + 1);
115 int numNewRowGids = 0;
116 for (int i=0; i<origRowMap->NumMyElements(); i=i+2)
117 newRowMapGids[numNewRowGids++] = origGids[i];
118 Epetra_Map newRowMap(-1, numNewRowGids, newRowMapGids.Values(), 0, origRowMap->Comm());
119
120
121 // Use this new map to create a SubCopy transform that we can use on-demand to create submatrices
122 EpetraExt::CrsMatrix_SubCopy subMatrixTransform(newRowMap);
123
124 // Use the transform we just created to create our submatrix
125 Epetra_CrsMatrix & subA = subMatrixTransform(*laplace_2d_matrix);
126 if (verbose) cout << "Sub matrix (every other row/column) = " << subA << endl;
127
128
129 // Now test it using the fwd() method to see if changes in the original matrix are carried to the submatrix
130 (*laplace_2d_matrix)[0][0] = 12.0;
131 if (verbose) cout << "Orig matrix = " << *laplace_2d_matrix << endl;
132 subMatrixTransform.fwd();
133 assert(subA[0][0]==12.0);
134 if (verbose) cout << "Sub matrix (every other row/column) = " << subA << endl;
135
136
137 // Now change the submatrix and use the rvs() method to see if changes are carried to the original matrix
138 subA[0][0] = 24.0;
139 if (verbose) cout << "Sub matrix (every other row/column) = " << subA << endl;
140 subMatrixTransform.rvs();
141 assert((*laplace_2d_matrix)[0][0]==24.0);
142 if (verbose) cout << "Orig matrix = " << *laplace_2d_matrix << endl;
143
144 if (Comm.MyPID()==0) cout << "EpetraExt::CrsMatrix_SubCopy tests passed." << endl;
145#ifdef EPETRA_MPI
146 MPI_Finalize() ;
147#endif
148
149return 0;
150}
151
152
std::string Epetra_Version()
Generates a sub-block view of a Epetra_CrsMatrix.
bool fwd()
Forward transfer of data from orig object input in the operator() method call to the new object creat...
bool rvs()
Reverse transfer of data from new object created in the operator() method call to the orig object inp...
int MyGlobalElements(int *MyGlobalElementList) const
const Epetra_Comm & Comm() const
int NumMyElements() const
int NumProc() const
int MyPID() const
int main(int argc, char *argv[])