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 #if defined(Epetra_SHOW_DEPRECATED_WARNINGS)
48 #ifdef __GNUC__
49 #warning "The Epetra package is deprecated"
50 #endif
51 #endif
52 
53 
54 
55 #include <Epetra_ConfigDefs.h>
56 #include <Epetra_CrsMatrix.h>
57 #include <Epetra_CombineMode.h>
58 
59 #include <vector>
60 
61 class Epetra_Map;
63 
64 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
66 #endif
67 
69 class Epetra_FECrsGraph;
70 
128 class EPETRA_LIB_DLL_EXPORT Epetra_FECrsMatrix : public Epetra_CrsMatrix {
129  public:
132  const Epetra_Map& RowMap,
133  int* NumEntriesPerRow,
134  bool ignoreNonLocalEntries=false);
135 
138  const Epetra_Map& RowMap,
139  int NumEntriesPerRow,
140  bool ignoreNonLocalEntries=false);
141 
144  const Epetra_Map& RowMap,
145  const Epetra_Map& ColMap,
146  int* NumEntriesPerRow,
147  bool ignoreNonLocalEntries=false);
148 
151  const Epetra_Map& RowMap,
152  const Epetra_Map& ColMap,
153  int NumEntriesPerRow,
154  bool ignoreNonLocalEntries=false);
155 
158  const Epetra_CrsGraph& Graph,
159  bool ignoreNonLocalEntries=false);
160 
163  const Epetra_FECrsGraph& Graph,
164  bool ignoreNonLocalEntries=false);
165 
168 
170  virtual ~Epetra_FECrsMatrix();
171 
174 
175  enum { ROW_MAJOR = 0, COLUMN_MAJOR = 3 };
176 
177  void Print(std::ostream& os) const;
178 #if !defined(EPETRA_NO_32BIT_GLOBAL_INDICES) || !defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
182 #endif
183 
185 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
186  int SumIntoGlobalValues(int GlobalRow, int NumEntries,
187  const double* Values, const int* Indices);
188 #endif
189 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
190  int SumIntoGlobalValues(long long GlobalRow, int NumEntries,
191  const double* Values, const long long* Indices);
192 #endif
193 
195 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
196  int InsertGlobalValues(int GlobalRow, int NumEntries,
197  const double* Values, const int* Indices);
198 #endif
199 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
200  int InsertGlobalValues(long long GlobalRow, int NumEntries,
201  const double* Values, const long long* Indices);
202 #endif
203 
205 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
206  int InsertGlobalValues(int GlobalRow, int NumEntries,
207  double* Values, int* Indices);
208 #endif
209 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
210  int InsertGlobalValues(long long GlobalRow, int NumEntries,
211  double* Values, long long* Indices);
212 #endif
213 
215 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
216  int ReplaceGlobalValues(int GlobalRow, int NumEntries,
217  const double* Values, const int* Indices);
218 #endif
219 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
220  int ReplaceGlobalValues(long long GlobalRow, int NumEntries,
221  const double* Values, const long long* Indices);
222 #endif
223 
239 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
240  int SumIntoGlobalValues(int numIndices, const int* indices,
241  const double* values,
243 #endif
244 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
245  int SumIntoGlobalValues(int numIndices, const long long* indices,
246  const double* values,
248 #endif
249 
266 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
267  int SumIntoGlobalValues(int numRows, const int* rows,
268  int numCols, const int* cols,
269  const double* values,
271 #endif
272 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
273  int SumIntoGlobalValues(int numRows, const long long* rows,
274  int numCols, const long long* cols,
275  const double* values,
277 #endif
278 
293 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
294  int SumIntoGlobalValues(int numIndices, const int* indices,
295  const double* const* values,
296  int format=Epetra_FECrsMatrix::ROW_MAJOR);
297 #endif
298 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
299  int SumIntoGlobalValues(int numIndices, const long long* indices,
300  const double* const* values,
301  int format=Epetra_FECrsMatrix::ROW_MAJOR);
302 #endif
303 
319 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
320  int SumIntoGlobalValues(int numRows, const int* rows,
321  int numCols, const int* cols,
322  const double* const* values,
323  int format=Epetra_FECrsMatrix::ROW_MAJOR);
324 #endif
325 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
326  int SumIntoGlobalValues(int numRows, const long long* rows,
327  int numCols, const long long* cols,
328  const double* const* values,
329  int format=Epetra_FECrsMatrix::ROW_MAJOR);
330 #endif
331 
346 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
347  int InsertGlobalValues(int numIndices, const int* indices,
348  const double* values,
350 #endif
351 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
352  int InsertGlobalValues(int numIndices, const long long* indices,
353  const double* values,
355 #endif
356 
372 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
373  int InsertGlobalValues(int numRows, const int* rows,
374  int numCols, const int* cols,
375  const double* values,
377 #endif
378 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
379  int InsertGlobalValues(int numRows, const long long* rows,
380  int numCols, const long long* cols,
381  const double* values,
383 #endif
384 
398 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
399  int InsertGlobalValues(int numIndices, const int* indices,
400  const double* const* values,
401  int format=Epetra_FECrsMatrix::ROW_MAJOR);
402 #endif
403 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
404  int InsertGlobalValues(int numIndices, const long long* indices,
405  const double* const* values,
406  int format=Epetra_FECrsMatrix::ROW_MAJOR);
407 #endif
408 
423 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
424  int InsertGlobalValues(int numRows, const int* rows,
425  int numCols, const int* cols,
426  const double* const* values,
427  int format=Epetra_FECrsMatrix::ROW_MAJOR);
428 #endif
429 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
430  int InsertGlobalValues(int numRows, const long long* rows,
431  int numCols, const long long* cols,
432  const double* const* values,
433  int format=Epetra_FECrsMatrix::ROW_MAJOR);
434 #endif
435 
451 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
452  int ReplaceGlobalValues(int numIndices, const int* indices,
453  const double* values,
455 #endif
456 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
457  int ReplaceGlobalValues(int numIndices, const long long* indices,
458  const double* values,
460 #endif
461 
479 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
480  int ReplaceGlobalValues(int numRows, const int* rows,
481  int numCols, const int* cols,
482  const double* values,
484 #endif
485 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
486  int ReplaceGlobalValues(int numRows, const long long* rows,
487  int numCols, const long long* cols,
488  const double* values,
490 #endif
491 
506 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
507  int ReplaceGlobalValues(int numIndices, const int* indices,
508  const double* const* values,
509  int format=Epetra_FECrsMatrix::ROW_MAJOR);
510 #endif
511 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
512  int ReplaceGlobalValues(int numIndices, const long long* indices,
513  const double* const* values,
514  int format=Epetra_FECrsMatrix::ROW_MAJOR);
515 #endif
516 
532 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
533  int ReplaceGlobalValues(int numRows, const int* rows,
534  int numCols, const int* cols,
535  const double* const* values,
536  int format=Epetra_FECrsMatrix::ROW_MAJOR);
537 #endif
538 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
539  int ReplaceGlobalValues(int numRows, const long long* rows,
540  int numCols, const long long* cols,
541  const double* const* values,
542  int format=Epetra_FECrsMatrix::ROW_MAJOR);
543 #endif
544 
555 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
557  const Epetra_SerialDenseMatrix& values,
559 #endif
560 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
562  const Epetra_SerialDenseMatrix& values,
564 #endif
565 
580 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
582  const Epetra_IntSerialDenseVector& cols,
583  const Epetra_SerialDenseMatrix& values,
585 #endif
586 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
589  const Epetra_SerialDenseMatrix& values,
591 #endif
592 
603 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
605  const Epetra_SerialDenseMatrix& values,
607 #endif
608 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
610  const Epetra_SerialDenseMatrix& values,
612 #endif
613 
628 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
630  const Epetra_IntSerialDenseVector& cols,
631  const Epetra_SerialDenseMatrix& values,
633 #endif
634 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
637  const Epetra_SerialDenseMatrix& values,
639 #endif
640 
652 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
654  const Epetra_SerialDenseMatrix& values,
656 #endif
657 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
659  const Epetra_SerialDenseMatrix& values,
661 #endif
662 
677 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
679  const Epetra_IntSerialDenseVector& cols,
680  const Epetra_SerialDenseMatrix& values,
682 #endif
683 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
686  const Epetra_SerialDenseMatrix& values,
688 #endif
689 
711  int GlobalAssemble(bool callFillComplete=true,
712  Epetra_CombineMode combineMode=Add,
713  bool save_off_and_reuse_map_exporter=false);
714 
740  int GlobalAssemble(const Epetra_Map& domain_map,
741  const Epetra_Map& range_map,
742  bool callFillComplete=true,
743  Epetra_CombineMode combineMode=Add,
744  bool save_off_and_reuse_map_exporter=false);
745 
749  void setIgnoreNonLocalEntries(bool flag) {
750  ignoreNonLocalEntries_ = flag;
751  }
752 
753  private:
754  void DeleteMemory();
755 
756  enum {SUMINTO = 0, REPLACE = 1, INSERT = 2};
757 
758  template<typename int_type>
759  int InputGlobalValues(int numRows, const int_type* rows,
760  int numCols, const int_type* cols,
761  const double* const* values,
762  int format,
763  int mode);
764 
765  template<typename int_type>
766  int InputGlobalValues(int numRows, const int_type* rows,
767  int numCols, const int_type* cols,
768  const double* values,
769  int format,
770  int mode);
771 
772  template<typename int_type>
773  int InputNonlocalGlobalValues(int_type row,
774  int numCols, const int_type* cols,
775  const double* values,
776  int mode);
777 
778  template<typename int_type>
779  int InputGlobalValues_RowMajor(
780  int numRows, const int_type* rows,
781  int numCols, const int_type* cols,
782  const double* values,
783  int mode);
784 
785  template<typename int_type>
786  int InsertNonlocalRow(int_type row, typename std::vector<int_type>::iterator offset);
787 
788  template<typename int_type>
789  int InputNonlocalValue(int rowoffset,
790  int_type col, double value,
791  int mode);
792 
793  long long myFirstRow_;
795 
797 
798 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
799  std::vector<int> nonlocalRows_int_;
800  std::vector<std::vector<int> > nonlocalCols_int_;
801 #endif
802 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
803  std::vector<long long> nonlocalRows_LL_;
804  std::vector<std::vector<long long> > nonlocalCols_LL_;
805 #endif
806 
807  template<typename int_type> std::vector<int_type>& nonlocalRows();
808  template<typename int_type> std::vector<std::vector<int_type> >& nonlocalCols();
809 
810  std::vector<std::vector<double> > nonlocalCoefs_;
811 
812  //IMPORTANT NOTE: The use of class-member work-data arrays is
813  //**NOT** thread-safe.
814  std::vector<double> workData_;
815  std::vector<const double*> workData2d_;
817 
820 
825 
826  template<typename int_type>
827  int SumIntoGlobalValues(int_type GlobalRow, int NumEntries, const double* values, const int_type* Indices);
828 
829  template<typename int_type>
830  int GlobalAssemble(const Epetra_Map& domain_map,
831  const Epetra_Map& range_map,
832  bool callFillComplete=true,
833  Epetra_CombineMode combineMode=Add,
834  bool save_off_and_reuse_map_exporter=false);
835 
836 };//class Epetra_FECrsMatrix
837 
838 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
839 template<> inline std::vector<int>& Epetra_FECrsMatrix::nonlocalRows<int>()
840 {
841  return nonlocalRows_int_;
842 }
843 #endif
844 
845 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
846 template<> inline std::vector<long long>& Epetra_FECrsMatrix::nonlocalRows<long long>()
847 {
848  return nonlocalRows_LL_;
849 }
850 #endif
851 
852 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
853 template<> inline std::vector<std::vector<int> >& Epetra_FECrsMatrix::nonlocalCols<int>()
854 {
855  return nonlocalCols_int_;
856 }
857 #endif
858 
859 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
860 template<> inline std::vector<std::vector<long long> >& Epetra_FECrsMatrix::nonlocalCols<long long>()
861 {
862  return nonlocalCols_LL_;
863 }
864 #endif
865 
866 
867 #endif /* EPETRA_FECRSMATRIX_H */
std::vector< const double * > workData2d_
Epetra_Map: A class for partitioning vectors and matrices.
Definition: Epetra_Map.h:127
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:70
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_