Epetra Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Epetra_FECrsMatrix.h
Go to the documentation of this file.
1 /*
2 //@HEADER
3 // ************************************************************************
4 //
5 // Epetra: Linear Algebra Services Package
6 // Copyright 2011 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 Michael A. Heroux (maherou@sandia.gov)
39 //
40 // ************************************************************************
41 //@HEADER
42 */
43 
44 #ifndef EPETRA_FECRSMATRIX_H
45 #define EPETRA_FECRSMATRIX_H
46 
47 #include <Epetra_ConfigDefs.h>
48 #include <Epetra_CrsMatrix.h>
49 #include <Epetra_CombineMode.h>
50 
51 #include <vector>
52 
53 class Epetra_Map;
55 
56 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
58 #endif
59 
61 class Epetra_FECrsGraph;
62 
120 class EPETRA_LIB_DLL_EXPORT Epetra_FECrsMatrix : public Epetra_CrsMatrix {
121  public:
124  const Epetra_Map& RowMap,
125  int* NumEntriesPerRow,
126  bool ignoreNonLocalEntries=false);
127 
130  const Epetra_Map& RowMap,
131  int NumEntriesPerRow,
132  bool ignoreNonLocalEntries=false);
133 
136  const Epetra_Map& RowMap,
137  const Epetra_Map& ColMap,
138  int* NumEntriesPerRow,
139  bool ignoreNonLocalEntries=false);
140 
143  const Epetra_Map& RowMap,
144  const Epetra_Map& ColMap,
145  int NumEntriesPerRow,
146  bool ignoreNonLocalEntries=false);
147 
150  const Epetra_CrsGraph& Graph,
151  bool ignoreNonLocalEntries=false);
152 
155  const Epetra_FECrsGraph& Graph,
156  bool ignoreNonLocalEntries=false);
157 
160 
162  virtual ~Epetra_FECrsMatrix();
163 
166 
167  enum { ROW_MAJOR = 0, COLUMN_MAJOR = 3 };
168 
169  void Print(std::ostream& os) const;
170 #if !defined(EPETRA_NO_32BIT_GLOBAL_INDICES) || !defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
174 #endif
175 
177 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
178  int SumIntoGlobalValues(int GlobalRow, int NumEntries,
179  const double* Values, const int* Indices);
180 #endif
181 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
182  int SumIntoGlobalValues(long long GlobalRow, int NumEntries,
183  const double* Values, const long long* Indices);
184 #endif
185 
187 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
188  int InsertGlobalValues(int GlobalRow, int NumEntries,
189  const double* Values, const int* Indices);
190 #endif
191 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
192  int InsertGlobalValues(long long GlobalRow, int NumEntries,
193  const double* Values, const long long* Indices);
194 #endif
195 
197 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
198  int InsertGlobalValues(int GlobalRow, int NumEntries,
199  double* Values, int* Indices);
200 #endif
201 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
202  int InsertGlobalValues(long long GlobalRow, int NumEntries,
203  double* Values, long long* Indices);
204 #endif
205 
207 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
208  int ReplaceGlobalValues(int GlobalRow, int NumEntries,
209  const double* Values, const int* Indices);
210 #endif
211 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
212  int ReplaceGlobalValues(long long GlobalRow, int NumEntries,
213  const double* Values, const long long* Indices);
214 #endif
215 
231 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
232  int SumIntoGlobalValues(int numIndices, const int* indices,
233  const double* values,
235 #endif
236 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
237  int SumIntoGlobalValues(int numIndices, const long long* indices,
238  const double* values,
240 #endif
241 
258 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
259  int SumIntoGlobalValues(int numRows, const int* rows,
260  int numCols, const int* cols,
261  const double* values,
263 #endif
264 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
265  int SumIntoGlobalValues(int numRows, const long long* rows,
266  int numCols, const long long* cols,
267  const double* values,
269 #endif
270 
285 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
286  int SumIntoGlobalValues(int numIndices, const int* indices,
287  const double* const* values,
288  int format=Epetra_FECrsMatrix::ROW_MAJOR);
289 #endif
290 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
291  int SumIntoGlobalValues(int numIndices, const long long* indices,
292  const double* const* values,
293  int format=Epetra_FECrsMatrix::ROW_MAJOR);
294 #endif
295 
311 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
312  int SumIntoGlobalValues(int numRows, const int* rows,
313  int numCols, const int* cols,
314  const double* const* values,
315  int format=Epetra_FECrsMatrix::ROW_MAJOR);
316 #endif
317 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
318  int SumIntoGlobalValues(int numRows, const long long* rows,
319  int numCols, const long long* cols,
320  const double* const* values,
321  int format=Epetra_FECrsMatrix::ROW_MAJOR);
322 #endif
323 
338 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
339  int InsertGlobalValues(int numIndices, const int* indices,
340  const double* values,
342 #endif
343 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
344  int InsertGlobalValues(int numIndices, const long long* indices,
345  const double* values,
347 #endif
348 
364 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
365  int InsertGlobalValues(int numRows, const int* rows,
366  int numCols, const int* cols,
367  const double* values,
369 #endif
370 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
371  int InsertGlobalValues(int numRows, const long long* rows,
372  int numCols, const long long* cols,
373  const double* values,
375 #endif
376 
390 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
391  int InsertGlobalValues(int numIndices, const int* indices,
392  const double* const* values,
393  int format=Epetra_FECrsMatrix::ROW_MAJOR);
394 #endif
395 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
396  int InsertGlobalValues(int numIndices, const long long* indices,
397  const double* const* values,
398  int format=Epetra_FECrsMatrix::ROW_MAJOR);
399 #endif
400 
415 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
416  int InsertGlobalValues(int numRows, const int* rows,
417  int numCols, const int* cols,
418  const double* const* values,
419  int format=Epetra_FECrsMatrix::ROW_MAJOR);
420 #endif
421 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
422  int InsertGlobalValues(int numRows, const long long* rows,
423  int numCols, const long long* cols,
424  const double* const* values,
425  int format=Epetra_FECrsMatrix::ROW_MAJOR);
426 #endif
427 
443 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
444  int ReplaceGlobalValues(int numIndices, const int* indices,
445  const double* values,
447 #endif
448 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
449  int ReplaceGlobalValues(int numIndices, const long long* indices,
450  const double* values,
452 #endif
453 
471 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
472  int ReplaceGlobalValues(int numRows, const int* rows,
473  int numCols, const int* cols,
474  const double* values,
476 #endif
477 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
478  int ReplaceGlobalValues(int numRows, const long long* rows,
479  int numCols, const long long* cols,
480  const double* values,
482 #endif
483 
498 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
499  int ReplaceGlobalValues(int numIndices, const int* indices,
500  const double* const* values,
501  int format=Epetra_FECrsMatrix::ROW_MAJOR);
502 #endif
503 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
504  int ReplaceGlobalValues(int numIndices, const long long* indices,
505  const double* const* values,
506  int format=Epetra_FECrsMatrix::ROW_MAJOR);
507 #endif
508 
524 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
525  int ReplaceGlobalValues(int numRows, const int* rows,
526  int numCols, const int* cols,
527  const double* const* values,
528  int format=Epetra_FECrsMatrix::ROW_MAJOR);
529 #endif
530 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
531  int ReplaceGlobalValues(int numRows, const long long* rows,
532  int numCols, const long long* cols,
533  const double* const* values,
534  int format=Epetra_FECrsMatrix::ROW_MAJOR);
535 #endif
536 
547 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
549  const Epetra_SerialDenseMatrix& values,
551 #endif
552 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
554  const Epetra_SerialDenseMatrix& values,
556 #endif
557 
572 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
574  const Epetra_IntSerialDenseVector& cols,
575  const Epetra_SerialDenseMatrix& values,
577 #endif
578 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
581  const Epetra_SerialDenseMatrix& values,
583 #endif
584 
595 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
597  const Epetra_SerialDenseMatrix& values,
599 #endif
600 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
602  const Epetra_SerialDenseMatrix& values,
604 #endif
605 
620 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
622  const Epetra_IntSerialDenseVector& cols,
623  const Epetra_SerialDenseMatrix& values,
625 #endif
626 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
629  const Epetra_SerialDenseMatrix& values,
631 #endif
632 
644 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
646  const Epetra_SerialDenseMatrix& values,
648 #endif
649 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
651  const Epetra_SerialDenseMatrix& values,
653 #endif
654 
669 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
671  const Epetra_IntSerialDenseVector& cols,
672  const Epetra_SerialDenseMatrix& values,
674 #endif
675 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
678  const Epetra_SerialDenseMatrix& values,
680 #endif
681 
703  int GlobalAssemble(bool callFillComplete=true,
704  Epetra_CombineMode combineMode=Add,
705  bool save_off_and_reuse_map_exporter=false);
706 
732  int GlobalAssemble(const Epetra_Map& domain_map,
733  const Epetra_Map& range_map,
734  bool callFillComplete=true,
735  Epetra_CombineMode combineMode=Add,
736  bool save_off_and_reuse_map_exporter=false);
737 
741  void setIgnoreNonLocalEntries(bool flag) {
742  ignoreNonLocalEntries_ = flag;
743  }
744 
745  private:
746  void DeleteMemory();
747 
748  enum {SUMINTO = 0, REPLACE = 1, INSERT = 2};
749 
750  template<typename int_type>
751  int InputGlobalValues(int numRows, const int_type* rows,
752  int numCols, const int_type* cols,
753  const double* const* values,
754  int format,
755  int mode);
756 
757  template<typename int_type>
758  int InputGlobalValues(int numRows, const int_type* rows,
759  int numCols, const int_type* cols,
760  const double* values,
761  int format,
762  int mode);
763 
764  template<typename int_type>
765  int InputNonlocalGlobalValues(int_type row,
766  int numCols, const int_type* cols,
767  const double* values,
768  int mode);
769 
770  template<typename int_type>
771  int InputGlobalValues_RowMajor(
772  int numRows, const int_type* rows,
773  int numCols, const int_type* cols,
774  const double* values,
775  int mode);
776 
777  template<typename int_type>
778  int InsertNonlocalRow(int_type row, typename std::vector<int_type>::iterator offset);
779 
780  template<typename int_type>
781  int InputNonlocalValue(int rowoffset,
782  int_type col, double value,
783  int mode);
784 
785  long long myFirstRow_;
787 
789 
790 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
791  std::vector<int> nonlocalRows_int_;
792  std::vector<std::vector<int> > nonlocalCols_int_;
793 #endif
794 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
795  std::vector<long long> nonlocalRows_LL_;
796  std::vector<std::vector<long long> > nonlocalCols_LL_;
797 #endif
798 
799  template<typename int_type> std::vector<int_type>& nonlocalRows();
800  template<typename int_type> std::vector<std::vector<int_type> >& nonlocalCols();
801 
802  std::vector<std::vector<double> > nonlocalCoefs_;
803 
804  //IMPORTANT NOTE: The use of class-member work-data arrays is
805  //**NOT** thread-safe.
806  std::vector<double> workData_;
807  std::vector<const double*> workData2d_;
809 
812 
817 
818  template<typename int_type>
819  int SumIntoGlobalValues(int_type GlobalRow, int NumEntries, const double* values, const int_type* Indices);
820 
821  template<typename int_type>
822  int GlobalAssemble(const Epetra_Map& domain_map,
823  const Epetra_Map& range_map,
824  bool callFillComplete=true,
825  Epetra_CombineMode combineMode=Add,
826  bool save_off_and_reuse_map_exporter=false);
827 
828 };//class Epetra_FECrsMatrix
829 
830 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
831 template<> inline std::vector<int>& Epetra_FECrsMatrix::nonlocalRows<int>()
832 {
833  return nonlocalRows_int_;
834 }
835 #endif
836 
837 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
838 template<> inline std::vector<long long>& Epetra_FECrsMatrix::nonlocalRows<long long>()
839 {
840  return nonlocalRows_LL_;
841 }
842 #endif
843 
844 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
845 template<> inline std::vector<std::vector<int> >& Epetra_FECrsMatrix::nonlocalCols<int>()
846 {
847  return nonlocalCols_int_;
848 }
849 #endif
850 
851 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
852 template<> inline std::vector<std::vector<long long> >& Epetra_FECrsMatrix::nonlocalCols<long long>()
853 {
854  return nonlocalCols_LL_;
855 }
856 #endif
857 
858 
859 #endif /* EPETRA_FECRSMATRIX_H */
std::vector< const double * > workData2d_
Epetra_Map: A class for partitioning vectors and matrices.
Definition: Epetra_Map.h:119
Epetra_IntSerialDenseVector: A class for constructing and using dense vectors.
Epetra_CrsMatrix * tempMat_
std::vector< std::vector< int > > nonlocalCols_int_
virtual int InsertGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
Insert a list of elements in a given global row of the matrix.
Epetra Finite-Element CrsGraph.
virtual int SumIntoGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
Add this list of entries to existing values for a given global row of the matrix. ...
Epetra_Export: This class builds an export object for efficient exporting of off-processor elements...
Definition: Epetra_Export.h:62
Epetra_SerialDenseMatrix: A class for constructing and using real double precision general dense matr...
std::vector< int > nonlocalRows_int_
Epetra_Combine Mode enumerable type.
std::vector< long long > nonlocalRows_LL_
Epetra Finite-Element CrsMatrix.
Epetra_Export * exporter_
Epetra_LongLongSerialDenseVector: A class for constructing and using dense vectors.
virtual void Print(std::ostream &os) const
Print method.
std::vector< std::vector< long long > > nonlocalCols_LL_
std::vector< std::vector< double > > nonlocalCoefs_
void setIgnoreNonLocalEntries(bool flag)
Set whether or not non-local data values should be ignored.
Epetra_CrsMatrix: A class for constructing and using real-valued double-precision sparse compressed r...
Epetra_CrsMatrix & operator=(const Epetra_CrsMatrix &src)
Assignment operator.
Epetra_CombineMode
Epetra_DataAccess
std::vector< double > workData_
virtual int ReplaceGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
Replace specified existing values with this list of entries for a given global row of the matrix...
Epetra_CrsGraph: A class for constructing and using sparse compressed row graphs. ...
Epetra_CrsMatrix * nonlocalMatrix_