EpetraExt Development
Loading...
Searching...
No Matches
EpetraExt_SolverMap_CrsMatrix.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
43
44#include <Epetra_CrsGraph.h>
45#include <Epetra_CrsMatrix.h>
46#include <Epetra_Map.h>
47#include <Epetra_Comm.h>
48
49#include <vector>
50
51namespace EpetraExt {
52
55{
56 if( newObj_ && newObj_ != origObj_ ) delete newObj_;
57 if( NewGraph_ ) delete NewGraph_;
58 if( NewColMap_ ) delete NewColMap_;
59}
60
61 template<typename int_type>
64construct( OriginalTypeRef orig )
65{
66 origObj_ = &orig;
67
68 assert( !orig.IndicesAreGlobal() );
69
70 //test if matrix has missing local columns in its col std::map
71 const Epetra_Map & RowMap = orig.RowMap();
72 const Epetra_Map & DomainMap = orig.DomainMap();
73 const Epetra_Map & ColMap = orig.ColMap();
74 const Epetra_Comm & Comm = RowMap.Comm();
75 int NumMyRows = RowMap.NumMyElements();
76 int NumCols = DomainMap.NumMyElements();
77 int Match = 0;
78 for( int i = 0; i < NumCols; ++i )
79 if( DomainMap.GID64(i) != ColMap.GID64(i) )
80 {
81 Match = 1;
82 break;
83 }
84
85 int MatchAll = 0;
86 Comm.SumAll( &Match, &MatchAll, 1 );
87
88 if( !MatchAll )
89 {
91 }
92 else
93 {
94 //create ColMap with all local rows included
95 std::vector<int_type> Cols(NumCols);
96 //fill Cols list with GIDs of all local columns
97 for( int i = 0; i < NumCols; ++i )
98 Cols[i] = (int_type) DomainMap.GID64(i);
99
100 //now append to Cols any ghost column entries
101 int NumMyCols = ColMap.NumMyElements();
102 for( int i = 0; i < NumMyCols; ++i )
103 if( !DomainMap.MyGID( ColMap.GID64(i) ) ) Cols.push_back( (int_type) ColMap.GID64(i) );
104
105 int NewNumMyCols = Cols.size();
106 int NewNumGlobalCols;
107 Comm.SumAll( &NewNumMyCols, &NewNumGlobalCols, 1 );
108 //create new column std::map
109 NewColMap_ = new Epetra_Map( NewNumGlobalCols, NewNumMyCols,&Cols[0], DomainMap.IndexBase64(), Comm );
110
111 //New Graph
112 std::vector<int> NumIndicesPerRow( NumMyRows );
113 for( int i = 0; i < NumMyRows; ++i )
114 NumIndicesPerRow[i] = orig.NumMyEntries(i);
115 NewGraph_ = new Epetra_CrsGraph( Copy, RowMap, *NewColMap_, &NumIndicesPerRow[0] );
116
117 int MaxNumEntries = orig.MaxNumEntries();
118 int NumEntries;
119 std::vector<int_type> Indices( MaxNumEntries );
120 for( int i = 0; i < NumMyRows; ++i )
121 {
122 int_type RowGID = (int_type) RowMap.GID64(i);
123 orig.Graph().ExtractGlobalRowCopy( RowGID, MaxNumEntries, NumEntries, &Indices[0] );
124 NewGraph_->InsertGlobalIndices( RowGID, NumEntries, &Indices[0] );
125 }
126 const Epetra_Map & RangeMap = orig.RangeMap();
127 NewGraph_->FillComplete(DomainMap,RangeMap);
128
129 //intial construction of matrix
130 Epetra_CrsMatrix * NewMatrix = new Epetra_CrsMatrix( View, *NewGraph_ );
131
132 //insert views of row values
133 int * myIndices;
134 double * myValues;
135 int indicesCnt;
136 int numMyRows = NewMatrix->NumMyRows();
137 for( int i = 0; i < numMyRows; ++i )
138 {
139 orig.ExtractMyRowView( i, indicesCnt, myValues, myIndices );
140 NewGraph_->ExtractMyRowView( i, indicesCnt, myIndices );
141
142 NewMatrix->InsertMyValues( i, indicesCnt, myValues, myIndices );
143 }
144
145 NewMatrix->FillComplete(DomainMap,RangeMap);
146
147 newObj_ = NewMatrix;
148 }
149
150 return *newObj_;
151}
152
155operator()( OriginalTypeRef orig )
156{
157#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
158 if(orig.RowMap().GlobalIndicesInt()) {
159 return construct<int>(orig);
160 }
161 else
162#endif
163#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
164 if(orig.RowMap().GlobalIndicesLongLong()) {
165 return construct<long long>(orig);
166 }
167 else
168#endif
169 throw "CrsMatrix_SolverMap::operator(): GlobalIndices type unknown";
170}
171} // namespace EpetraExt
172
NewTypeRef operator()(OriginalTypeRef orig)
Constructs fixed view of Epetra_CrsMatrix as necessary.
virtual NewTypeRef construct()
Construction of new object as a result of the transform.
const Epetra_Comm & Comm() const
int NumMyElements() const
bool MyGID(int GID_in) const
virtual int SumAll(double *PartialSums, double *GlobalSums, int Count) const=0
int InsertGlobalIndices(int GlobalRow, int NumIndices, int *Indices)
int ExtractMyRowView(int LocalRow, int &NumIndices, int *&Indices) const
EpetraExt::BlockCrsMatrix: A class for constructing a distributed block matrix.