Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Epetra_IntSerialDenseMatrix.cpp
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
44#include "Epetra_Util.h"
45
46//=============================================================================
48 : Epetra_Object("Epetra::IntSerialDenseMatrix"),
49 CV_(Copy),
50 A_Copied_(false),
51 M_(0),
52 N_(0),
53 LDA_(0),
54 A_(0)
55{
56}
57
58//=============================================================================
60 : Epetra_Object("Epetra::IntSerialDenseMatrix"),
61 CV_(Copy),
62 A_Copied_(false),
63 M_(0),
64 N_(0),
65 LDA_(0),
66 A_(0)
67{
68 if(NumRows < 0)
69 throw ReportError("NumRows = " + toString(NumRows) + ". Should be >= 0", -1);
70 if(NumCols < 0)
71 throw ReportError("NumCols = " + toString(NumCols) + ". Should be >= 0", -1);
72
73 int errorcode = Shape(NumRows, NumCols);
74 if(errorcode != 0)
75 throw ReportError("Shape returned non-zero (" + toString(errorcode) + ").", -2);
76}
77//=============================================================================
79 int NumRows, int NumCols)
80 : Epetra_Object("Epetra::IntSerialDenseMatrix"),
81 CV_(CV_in),
82 A_Copied_(false),
83 M_(NumRows),
84 N_(NumCols),
85 LDA_(lda),
86 A_(A_in)
87{
88 if(A_in == 0)
89 throw ReportError("Null pointer passed as A_in parameter.", -3);
90 if(NumRows < 0)
91 throw ReportError("NumRows = " + toString(NumRows) + ". Should be >= 0", -1);
92 if(NumCols < 0)
93 throw ReportError("NumCols = " + toString(NumCols) + ". Should be >= 0", -1);
94 if(lda < 0)
95 throw ReportError("LDA = " + toString(lda) + ". Should be >= 0", -1);
96
97 if(CV_in == Copy) {
98 LDA_ = M_;
99 const int newsize = LDA_ * N_;
100 if(newsize > 0) {
101 A_ = new int[newsize];
102 CopyMat(A_in, lda, M_, N_, A_, LDA_);
103 A_Copied_ = true;
104 }
105 else {
106 A_ = 0;
107 }
108 }
109}
110//=============================================================================
112 : Epetra_Object(Source),
113 CV_(Source.CV_),
114 A_Copied_(false),
115 M_(Source.M_),
116 N_(Source.N_),
117 LDA_(Source.LDA_),
118 A_(Source.A_)
119{
120 if(CV_ == Copy) {
121 LDA_ = M_;
122 const int newsize = LDA_ * N_;
123 if(newsize > 0) {
124 A_ = new int[newsize];
125 CopyMat(Source.A_, Source.LDA_, M_, N_, A_, LDA_);
126 A_Copied_ = true;
127 }
128 else {
129 A_ = 0;
130 A_Copied_ = false;
131 }
132 }
133}
134//=============================================================================
135int Epetra_IntSerialDenseMatrix::Reshape(int NumRows, int NumCols) {
136 if(NumRows < 0 || NumCols < 0)
137 return(-1);
138
139 int* A_tmp = 0;
140 const int newsize = NumRows * NumCols;
141
142 if(newsize > 0) {
143 // Allocate space for new matrix
144 A_tmp = new int[newsize];
145 for(int k = 0; k < newsize; k++)
146 A_tmp[k] = 0; // Zero out values
147 int M_tmp = EPETRA_MIN(M_, NumRows);
148 int N_tmp = EPETRA_MIN(N_, NumCols);
149 if(A_ != 0)
150 CopyMat(A_, LDA_, M_tmp, N_tmp, A_tmp, NumRows); // Copy principal submatrix of A to new A
151 }
152
153 CleanupData(); // Get rid of anything that might be already allocated
154 M_ = NumRows;
155 N_ = NumCols;
156 LDA_ = M_;
157 A_ = A_tmp; // Set pointer to new A
158 A_Copied_ = (newsize>0);
159 return(0);
160}
161//=============================================================================
162int Epetra_IntSerialDenseMatrix::Shape(int NumRows, int NumCols) {
163 if(NumRows < 0 || NumCols < 0)
164 return(-1);
165
166 CleanupData(); // Get rid of anything that might be already allocated
167 M_ = NumRows;
168 N_ = NumCols;
169 LDA_ = M_;
170 const int newsize = LDA_ * N_;
171 if(newsize > 0) {
172 A_ = new int[newsize];
173#ifdef EPETRA_HAVE_OMP
174#pragma omp parallel for
175#endif
176 for(int k = 0; k < newsize; k++)
177 A_[k] = 0; // Zero out values
178 A_Copied_ = true;
179 }
180
181 return(0);
182}
183//=============================================================================
188//=============================================================================
190{
191 if(A_Copied_)
192 delete[] A_;
193 A_ = 0;
194 A_Copied_ = false;
195 M_ = 0;
196 N_ = 0;
197 LDA_ = 0;
198}
199//=============================================================================
201 if(this == &Source)
202 return(*this); // Special case of source same as target
203 if((CV_ == View) && (Source.CV_ == View) && (A_ == Source.A_))
204 return(*this); // Special case of both are views to same data.
205
206 if(std::strcmp(Label(), Source.Label()) != 0)
207 throw ReportError("operator= type mismatch (lhs = " + std::string(Label()) +
208 ", rhs = " + std::string(Source.Label()) + ").", -5);
209
210 if(Source.CV_ == View) {
211 if(CV_ == Copy) { // C->V only
212 CleanupData();
213 CV_ = View;
214 }
215 M_ = Source.M_; // C->V and V->V
216 N_ = Source.N_;
217 LDA_ = Source.LDA_;
218 A_ = Source.A_;
219 }
220 else {
221 if(CV_ == View) { // V->C
222 CV_ = Copy;
223 M_ = Source.M_;
224 N_ = Source.N_;
225 LDA_ = Source.M_;
226 const int newsize = LDA_ * N_;
227 if(newsize > 0) {
228 A_ = new int[newsize];
229 A_Copied_ = true;
230 }
231 else {
232 A_ = 0;
233 A_Copied_ = false;
234 }
235 }
236 else { // C->C
237 if((Source.M_ <= LDA_) && (Source.N_ == N_)) { // we don't need to reallocate
238 M_ = Source.M_;
239 N_ = Source.N_;
240 }
241 else { // we need to allocate more space (or less space)
242 CleanupData();
243 M_ = Source.M_;
244 N_ = Source.N_;
245 LDA_ = Source.M_;
246 const int newsize = LDA_ * N_;
247 if(newsize > 0) {
248 A_ = new int[newsize];
249 A_Copied_ = true;
250 }
251 }
252 }
253 CopyMat(Source.A_, Source.LDA_, M_, N_, A_, LDA_); // V->C and C->C
254 }
255
256 return(*this);
257}
258
259//=============================================================================
261{
262 if (M_ != rhs.M_ || N_ != rhs.N_) return(false);
263
264 const int* A_tmp = A_;
265 const int* rhsA = rhs.A_;
266
267 for(int j=0; j<N_; ++j) {
268 int offset = j*LDA_;
269 int rhsOffset = j*rhs.LDA_;
270 for(int i=0; i<M_; ++i) {
271 if (A_tmp[offset+i] != rhsA[rhsOffset+i]) {
272 return(false);
273 }
274 }
275 }
276
277 return(true);
278}
279
280//=============================================================================
282 if(std::strcmp(Label(), Source.Label()) != 0)
283 return(-1);
284
285 if(CV_ == Copy) {
286 CleanupData();
287 CV_ = View;
288 }
289 M_ = Source.M_;
290 N_ = Source.N_;
291 LDA_ = Source.LDA_;
292 A_ = Source.A_;
293
294 return(0);
295}
296
297//=============================================================================
298void Epetra_IntSerialDenseMatrix::CopyMat(int* Source, int Source_LDA, int NumRows, int NumCols,
299 int* Target, int Target_LDA)
300{
301 int i, j;
302 int* targetPtr = Target;
303 int* sourcePtr = 0;
304 for(j = 0; j < NumCols; j++) { // for each column
305 targetPtr = Target + j * Target_LDA; // set pointers to correct stride
306 sourcePtr = Source + j * Source_LDA;
307 for(i = 0; i < NumRows; i++) // for each row
308 *targetPtr++ = *sourcePtr++; // copy element (i,j) and increment pointer to (i,j+1)
309 }
310 return;
311}
312
313//=============================================================================
315 int anorm = 0;
316 int* ptr = 0;
317 for(int j = 0; j < N_; j++) {
318 int sum = 0;
319 ptr = A_ + j*LDA_;
320 for(int i = 0; i < M_; i++)
321 sum += std::abs(*ptr++);
322 anorm = EPETRA_MAX(anorm, sum);
323 }
324 return(anorm);
325}
326
327//=============================================================================
329 int anorm = 0;
330 int* ptr = 0;
331 // Loop across columns in inner loop. Most expensive memory access, but
332 // requires no extra storage.
333 for(int i = 0; i < M_; i++) {
334 int sum = 0;
335 ptr = A_ + i;
336 for(int j = 0; j < N_; j++) {
337 sum += std::abs(*ptr);
338 ptr += LDA_;
339 }
340 anorm = EPETRA_MAX(anorm, sum);
341 }
342 return(anorm);
343}
344
345//=========================================================================
346void Epetra_IntSerialDenseMatrix::Print(std::ostream& os) const {
347 if(CV_ == Copy)
348 os << "Data access mode: Copy" << std::endl;
349 else
350 os << "Data access mode: View" << std::endl;
351 if(A_Copied_)
352 os << "A_Copied: yes" << std::endl;
353 else
354 os << "A_Copied: no" << std::endl;
355 os << "Rows(M): " << M_ << std::endl;
356 os << "Columns(N): " << N_ << std::endl;
357 os << "LDA: " << LDA_ << std::endl;
358 if(M_ == 0 || N_ == 0)
359 os << "(matrix is empty, no values to display)" << std::endl;
360 else
361 for(int i = 0; i < M_; i++) {
362 for(int j = 0; j < N_; j++){
363 os << (*this)(i,j) << " ";
364 }
365 os << std::endl;
366 }
367}
368
369//=========================================================================
371
372 Epetra_Util util;
373
374 for(int j = 0; j < N_; j++) {
375 int* arrayPtr = A_ + (j * LDA_);
376 for(int i = 0; i < M_; i++) {
377 *arrayPtr++ = util.RandomInt();
378 }
379 }
380
381 return(0);
382}
#define EPETRA_MIN(x, y)
#define EPETRA_MAX(x, y)
Epetra_DataAccess
Epetra_IntSerialDenseMatrix: A class for constructing and using general dense integer matrices.
int Shape(int NumRows, int NumCols)
Set dimensions of a Epetra_IntSerialDenseMatrix object; init values to zero.
Epetra_IntSerialDenseMatrix()
Default constructor; defines a zero size object.
int MakeViewOf(const Epetra_IntSerialDenseMatrix &Source)
Reset an existing IntSerialDenseMatrix to point to another Matrix.
int Reshape(int NumRows, int NumCols)
Reshape a Epetra_IntSerialDenseMatrix object.
void CopyMat(int *Source, int Source_LDA, int NumRows, int NumCols, int *Target, int Target_LDA)
bool operator==(const Epetra_IntSerialDenseMatrix &rhs) const
Comparison operator.
virtual void Print(std::ostream &os) const
Print service methods; defines behavior of ostream << operator.
virtual int InfNorm()
Computes the Infinity-Norm of the this matrix.
Epetra_IntSerialDenseMatrix & operator=(const Epetra_IntSerialDenseMatrix &Source)
Copy from one matrix to another.
virtual int OneNorm()
Computes the 1-Norm of the this matrix.
virtual ~Epetra_IntSerialDenseMatrix()
Epetra_IntSerialDenseMatrix destructor.
int Random()
Set matrix values to random numbers.
Epetra_Object: The base Epetra class.
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_Util: The Epetra Util Wrapper Class.
Definition Epetra_Util.h:79
unsigned int RandomInt()
Returns a random integer on the interval (0, 2^31-1)