Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Epetra_SerialDenseMatrix.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_SERIALDENSEMATRIX_H
45#define EPETRA_SERIALDENSEMATRIX_H
46
47#include "Epetra_ConfigDefs.h"
48#include "Epetra_Object.h"
49#include "Epetra_CompObject.h"
50#include "Epetra_BLAS.h"
54
56
106//=========================================================================
107class EPETRA_LIB_DLL_EXPORT Epetra_SerialDenseMatrix : public Epetra_CompObject, public Epetra_Object, public Epetra_SerialDenseOperator, public Epetra_BLAS {
108
109 public:
110
112
113
119 Epetra_SerialDenseMatrix(bool set_object_label=true);
120
122
133 Epetra_SerialDenseMatrix(int NumRows, int NumCols, bool set_object_label=true);
134
136
151 Epetra_SerialDenseMatrix(Epetra_DataAccess CV, double* A_in, int LDA_in, int NumRows, int NumCols,
152 bool set_object_label=true);
153
155
157
159 virtual ~Epetra_SerialDenseMatrix ();
161
163
164
177 int Shape(int NumRows, int NumCols);
178
180
193 int Reshape(int NumRows, int NumCols);
195
197
198
200
219 int Multiply(char TransA, char TransB, double ScalarAB,
222 double ScalarThis);
223
225 /* This method is intended to imitate the semantics of the matrix-vector
226 multiplication provided by Epetra's sparse matrices. The 'vector' arguments
227 are actually matrices; this method will return an error if the
228 dimensions of 'x' are not compatible. 'y' will be reshaped if necessary.
229 */
230 int Multiply(bool transA,
233
235
256 int Multiply(char SideA, double ScalarAB,
259 double ScalarThis);
260
262
270 int Scale(double ScalarA);
271
273
276 virtual double NormOne() const;
277
279 virtual double NormInf() const;
280
282
284
285
287
294
296
299 bool operator==(const Epetra_SerialDenseMatrix& rhs) const;
300
302
305 { return !(*this == rhs); }
306
308
315 Epetra_SerialDenseMatrix & operator += (const Epetra_SerialDenseMatrix& Source);
316
318
327 double& operator () (int RowIndex, int ColIndex);
328
330
339 const double& operator () (int RowIndex, int ColIndex) const;
340
342
352 double* operator [] (int ColIndex);
353
355
365 const double* operator [] (int ColIndex) const;
366
368
374 int Random();
375
377 int M() const {return(M_);};
378
380 int N() const {return(N_);};
381
383 double* A() const {return(A_);};
384
386 double* A() {return(A_);};
387
389 int LDA() const {return(LDA_);};
390
392 Epetra_DataAccess CV() const {return(CV_);};
394
396
397
398 virtual void Print(std::ostream& os) const;
400
402
403
405
408 virtual double OneNorm() const {return(NormOne());};
409
411 virtual double InfNorm() const {return(NormInf());};
413
415
416
418
427 virtual int SetUseTranspose(bool UseTranspose_in) { UseTranspose_ = UseTranspose_in; return (0); }
428
430
439
441
451 {
452 (void)X;//prevents unused variable compiler warning
453 (void)Y;
454 return (-1);
455 }
456
458 virtual const char * Label() const { return Epetra_Object::Label(); }
459
461 virtual bool UseTranspose() const { return UseTranspose_; }
462
464 virtual bool HasNormInf() const { return true; }
465
467 virtual int RowDim() const { return M(); }
468
470 virtual int ColDim() const { return N(); }
472
473 protected:
474
475 void CopyMat(const double* Source, int Source_LDA, int NumRows, int NumCols,
476 double* Target, int Target_LDA, bool add=false);
477 void CleanupData();
478
479 int M_;
480 int N_;
483
484 //For performance reasons, it's better if Epetra_VbrMatrix can access the
485 //LDA_ and A_ members of this class directly without going through an
486 //accessor method. Rather than making them public members, we'll make
487 //Epetra_VbrMatrix a friend class.
488
489 friend class Epetra_VbrMatrix;
490
491 int LDA_;
492 double* A_;
493
495};
496
497// inlined definitions of op() and op[]
498//=========================================================================
499inline double& Epetra_SerialDenseMatrix::operator () (int RowIndex, int ColIndex) {
500#ifdef HAVE_EPETRA_ARRAY_BOUNDS_CHECK
501 if (RowIndex >= M_ || RowIndex < 0)
502 throw ReportError("Row index = " +toString(RowIndex) +
503 " Out of Range 0 - " + toString(M_-1),-1);
504 if (ColIndex >= N_ || ColIndex < 0)
505 throw ReportError("Column index = " +toString(ColIndex) +
506 " Out of Range 0 - " + toString(N_-1),-2);
507#endif
508 return(A_[ColIndex*LDA_ + RowIndex]);
509}
510//=========================================================================
511inline const double& Epetra_SerialDenseMatrix::operator () (int RowIndex, int ColIndex) const {
512#ifdef HAVE_EPETRA_ARRAY_BOUNDS_CHECK
513 if (RowIndex >= M_ || RowIndex < 0)
514 throw ReportError("Row index = " +toString(RowIndex) +
515 " Out of Range 0 - " + toString(M_-1),-1);
516 if (ColIndex >= N_ || ColIndex < 0)
517 throw ReportError("Column index = " +toString(ColIndex) +
518 " Out of Range 0 - " + toString(N_-1),-2);
519#endif
520 return(A_[ColIndex*LDA_ + RowIndex]);
521}
522//=========================================================================
523inline double* Epetra_SerialDenseMatrix::operator [] (int ColIndex) {
524#ifdef HAVE_EPETRA_ARRAY_BOUNDS_CHECK
525 if (ColIndex >= N_ || ColIndex < 0)
526 throw ReportError("Column index = " +toString(ColIndex) +
527 " Out of Range 0 - " + toString(N_-1),-2);
528#endif
529 return(A_ + ColIndex*LDA_);
530}
531//=========================================================================
532inline const double* Epetra_SerialDenseMatrix::operator [] (int ColIndex) const {
533#ifdef HAVE_EPETRA_ARRAY_BOUNDS_CHECK
534 if (ColIndex >= N_ || ColIndex < 0)
535 throw ReportError("Column index = " +toString(ColIndex) +
536 " Out of Range 0 - " + toString(N_-1),-2);
537#endif
538 return(A_+ ColIndex*LDA_);
539}
540//=========================================================================
541
542#endif /* EPETRA_SERIALDENSEMATRIX_H */
Epetra_DataAccess
Epetra_BLAS: The Epetra BLAS Wrapper Class.
Definition Epetra_BLAS.h:70
Epetra_CompObject: Functionality and data that is common to all computational classes.
Epetra_CompObject & operator=(const Epetra_CompObject &src)
Epetra_Object: The base Epetra class.
virtual void Print(std::ostream &os) const
Print object to an output stream Print method.
virtual int ReportError(const std::string Message, int ErrorCode) const
Error reporting method.
virtual const char * Label() const
Epetra_Object Label access funtion.
std::string toString(const int &x) const
Epetra_SerialDenseMatrix: A class for constructing and using real double precision general dense matr...
virtual bool UseTranspose() const
Returns the current UseTranspose setting.
double * A() const
Returns pointer to the this matrix.
double * operator[](int ColIndex)
Column access function.
int LDA() const
Returns the leading dimension of the this matrix.
bool operator!=(const Epetra_SerialDenseMatrix &rhs) const
Inequality operator.
virtual double OneNorm() const
Computes the 1-Norm of the this matrix (identical to NormOne() method).
virtual int ApplyInverse(const Epetra_SerialDenseMatrix &X, Epetra_SerialDenseMatrix &Y)
Returns the result of a Epetra_SerialDenseOperator inverse applied to an Epetra_SerialDenseMatrix X i...
virtual bool HasNormInf() const
Returns true if the this object can provide an approximate Inf-norm, false otherwise.
virtual double InfNorm() const
Computes the Infinity-Norm of the this matrix (identical to NormInf() method).
virtual int ColDim() const
Returns the column dimension of operator.
virtual int SetUseTranspose(bool UseTranspose_in)
If set true, transpose of this operator will be applied.
int M() const
Returns row dimension of system.
double * A()
Returns pointer to the this matrix.
double & operator()(int RowIndex, int ColIndex)
Element access function.
int N() const
Returns column dimension of system.
Epetra_DataAccess CV() const
Returns the data access mode of the this matrix.
virtual const char * Label() const
Returns a character string describing the operator.
virtual int RowDim() const
Returns the row dimension of operator.
Epetra_SerialDenseOperator: A pure virtual class for using real-valued double-precision operators.
virtual int Apply(const Epetra_SerialDenseMatrix &X, Epetra_SerialDenseMatrix &Y)=0
Returns the result of a Epetra_SerialDenseOperator applied to a Epetra_SerialDenseMatrix X in Y.
virtual double NormInf() const =0
Returns the infinity norm of the global matrix.
Epetra_SerialSymDenseMatrix: A class for constructing and using symmetric positive definite dense mat...
Epetra_VbrMatrix: A class for the construction and use of real-valued double-precision variable block...