MatrixI.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2019-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "MatrixBlock.H"
30 
31 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
32 
33 template<class Form, class Type>
35 {
36  const label len = size();
37 
38  if (len > 0)
39  {
40  // With sign-check to avoid spurious -Walloc-size-larger-than
41  this->v_ = ListPolicy::allocate<Type>(len);
42 
43  }
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
49 template<class Form, class Type>
51 :
52  mRows_(0),
53  nCols_(0),
54  v_(nullptr)
55 {}
56 
57 
58 template<class Form, class Type>
60 :
61  Matrix<Form, Type>(dims.first(), dims.second())
62 {}
63 
64 
65 template<class Form, class Type>
67 :
68  Matrix<Form, Type>(dims.first(), dims.second(), Foam::zero{})
69 {}
70 
71 
72 template<class Form, class Type>
73 inline Foam::Matrix<Form, Type>::Matrix(const labelPair& dims, const Type& val)
74 :
75  Matrix<Form, Type>(dims.first(), dims.second(), val)
76 {}
77 
78 
79 template<class Form, class Type>
82 {
84 }
85 
86 
87 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
88 
89 template<class Form, class Type>
90 inline Foam::label Foam::Matrix<Form, Type>::size() const noexcept
91 {
92  return (mRows_ * nCols_);
93 }
94 
95 
96 template<class Form, class Type>
98 {
99  return labelPair(mRows_, nCols_);
100 }
101 
102 
103 template<class Form, class Type>
105 {
106  return !mRows_ || !nCols_;
107 }
108 
109 
110 template<class Form, class Type>
111 inline void Foam::Matrix<Form, Type>::checki(const label i) const
112 {
113  if (!mRows_ || !nCols_)
114  {
116  << "Attempt to access element from empty matrix"
117  << abort(FatalError);
118  }
119  if (i < 0 || mRows_ <= i)
120  {
122  << "Index " << i << " out of range 0 ... " << mRows_-1
123  << abort(FatalError);
124  }
125 }
126 
127 
128 template<class Form, class Type>
129 inline void Foam::Matrix<Form, Type>::checkj(const label j) const
130 {
131  if (!mRows_ || !nCols_)
132  {
134  << "Attempt to access element from empty matrix"
135  << abort(FatalError);
136  }
137  if (j < 0 || nCols_ <= j)
138  {
140  << "index " << j << " out of range 0 ... " << nCols_-1
141  << abort(FatalError);
142  }
143 }
144 
145 
146 template<class Form, class Type>
147 inline void Foam::Matrix<Form, Type>::checkSize() const
148 {
149  if (mRows_ < 0 || nCols_ < 0)
150  {
152  << "Incorrect size (" << mRows_ << ", " << nCols_ << ')' << nl
154  }
155  // Could also check for odd sizes, like (0, N) and make (0,0)
156 }
157 
158 
159 template<class Form, class Type>
160 inline bool Foam::Matrix<Form, Type>::uniform() const
161 {
162  // Assume Type has 'operator!=' defined
163  return
164  (
166  == Foam::ListPolicy::check_uniformity(this->v_, this->v_ + size())
167  );
168 }
169 
170 
171 template<class Form, class Type>
172 inline const Type* Foam::Matrix<Form, Type>::cdata() const noexcept
173 {
174  return v_;
175 }
176 
177 
178 template<class Form, class Type>
180 {
181  return v_;
182 }
183 
184 
185 template<class Form, class Type>
187 {
188  return reinterpret_cast<const char*>(v_);
189 }
190 
191 
192 template<class Form, class Type>
194 {
195  return reinterpret_cast<char*>(v_);
196 }
197 
198 
199 template<class Form, class Type>
200 inline std::streamsize Foam::Matrix<Form, Type>::size_bytes() const noexcept
201 {
202  return std::streamsize(mRows_*nCols_)*sizeof(Type);
203 }
204 
205 
206 template<class Form, class Type>
207 inline const Type* Foam::Matrix<Form, Type>::rowData(const label irow) const
208 {
209  #ifdef FULLDEBUG
210  checki(irow);
211  #endif
212  return (v_ + irow*nCols_);
213 }
214 
215 
216 template<class Form, class Type>
217 inline Type* Foam::Matrix<Form, Type>::rowData(const label irow)
218 {
219  #ifdef FULLDEBUG
220  checki(irow);
221  #endif
222  return (v_ + irow*nCols_);
223 }
224 
225 
226 template<class Form, class Type>
227 inline const Type& Foam::Matrix<Form, Type>::at(const label idx) const
228 {
229  #ifdef FULLDEBUG
230  if (idx < 0 || this->size() <= idx)
231  {
233  << "index " << idx << " out of range 0 ... " << this->size()
234  << abort(FatalError);
235  }
236  #endif
237  return *(v_ + idx);
238 }
239 
240 
241 template<class Form, class Type>
242 inline Type& Foam::Matrix<Form, Type>::at(const label idx)
243 {
244  #ifdef FULLDEBUG
245  if (idx < 0 || this->size() <= idx)
246  {
248  << "index " << idx << " out of range 0 ... " << this->size()
249  << abort(FatalError);
250  }
251  #endif
252  return *(v_ + idx);
253 }
254 
255 
256 template<class Form, class Type>
259 (
260  const label colIndex,
261  const label rowIndex,
262  label len
263 ) const
264 {
265  if (len < 0)
266  {
267  len = mRows_ - rowIndex;
268  }
269 
270  return ConstMatrixBlock<mType>
271  (
272  *this,
273  len, // rows
274  1,
275  rowIndex,
276  colIndex
277  );
278 }
279 
280 
281 template<class Form, class Type>
284 (
285  const label rowIndex,
286  const label colIndex,
287  label len
288 ) const
289 {
290  if (len < 0)
291  {
292  len = nCols_ - colIndex;
293  }
294 
295  return ConstMatrixBlock<mType>
296  (
297  *this,
298  1,
299  len, // columns
300  rowIndex,
301  colIndex
302  );
303 }
304 
305 
306 template<class Form, class Type>
309 (
310  const label rowIndex,
311  const label colIndex,
312  label szRows,
313  label szCols
314 ) const
315 {
316  if (szRows < 0) szRows = mRows_ - rowIndex;
317  if (szCols < 0) szCols = nCols_ - colIndex;
318 
320  (
321  *this,
322  szRows,
323  szCols,
324  rowIndex,
325  colIndex
326  );
327 }
328 
329 
330 template<class Form, class Type>
331 template<class VectorSpace>
334 (
335  const label rowIndex,
336  const label colIndex
337 ) const
338 {
340  (
341  *this,
342  VectorSpace::mRows,
343  VectorSpace::nCols,
344  rowIndex,
345  colIndex
346  );
347 }
348 
349 
350 template<class Form, class Type>
353 (
354  const label colIndex,
355  const label rowIndex,
356  label len
357 )
358 {
359  if (len < 0)
360  {
361  len = mRows_ - rowIndex;
362  }
363 
364  return MatrixBlock<mType>
365  (
366  *this,
367  len, // rows
368  1,
369  rowIndex,
370  colIndex
371  );
372 }
373 
374 
375 template<class Form, class Type>
378 (
379  const label rowIndex,
380  const label colIndex,
381  label len
382 )
383 {
384  if (len < 0)
385  {
386  len = nCols_ - colIndex;
387  }
388 
389  return MatrixBlock<mType>
390  (
391  *this,
392  1,
393  len, // columns
394  rowIndex,
395  colIndex
396  );
397 }
398 
399 
400 template<class Form, class Type>
403 (
404  const label rowIndex,
405  const label colIndex,
406  label szRows,
407  label szCols
408 )
409 {
410  if (szRows < 0) szRows = mRows_ - rowIndex;
411  if (szCols < 0) szCols = nCols_ - colIndex;
412 
413  return MatrixBlock<mType>
414  (
415  *this,
416  szRows,
417  szCols,
418  rowIndex,
419  colIndex
420  );
421 }
422 
423 
424 template<class Form, class Type>
425 template<class VectorSpace>
428 (
429  const label rowIndex,
430  const label colIndex
431 )
432 {
433  return MatrixBlock<mType>
434  (
435  *this,
436  VectorSpace::mRows,
437  VectorSpace::nCols,
438  rowIndex,
439  colIndex
440  );
441 }
442 
443 
444 template<class Form, class Type>
445 inline void Foam::Matrix<Form, Type>::setSize(const label m, const label n)
446 {
447  resize(m, n);
448 }
449 
450 
451 template<class Form, class Type>
452 void Foam::Matrix<Form, Type>::shallowResize(const label m, const label n)
453 {
454  mRows_ = m;
455  nCols_ = n;
456 }
457 
458 
459 template<class Form, class Type>
461 (
462  const UList<Type>& x
463 ) const
464 {
465  return this->AmulImpl(x);
466 }
467 
468 
469 template<class Form, class Type>
470 template<class Addr>
472 (
474 ) const
475 {
476  return this->AmulImpl(x);
477 }
478 
479 
480 template<class Form, class Type>
482 (
483  const UList<Type>& x
484 ) const
485 {
486  return this->TmulImpl(x);
487 }
488 
489 
490 template<class Form, class Type>
491 template<class Addr>
493 (
495 ) const
496 {
497  return this->TmulImpl(x);
498 }
499 
500 
501 // * * * * * * * * * * * * * * * * * Iterators * * * * * * * * * * * * * * * //
502 
503 template<class Form, class Type>
506 {
507  return v_;
508 }
509 
510 
511 template<class Form, class Type>
514 {
515  return v_ + (mRows_ * nCols_);
516 }
517 
518 
519 template<class Form, class Type>
522 {
523  return v_;
524 }
525 
526 
527 template<class Form, class Type>
530 {
531  return v_ + (mRows_ * nCols_);
532 }
533 
534 
535 template<class Form, class Type>
538 {
539  return v_;
540 }
541 
542 
543 template<class Form, class Type>
546 {
547  return v_ + (mRows_ * nCols_);
548 }
549 
550 
551 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
552 
553 template<class Form, class Type>
554 inline const Type& Foam::Matrix<Form, Type>::operator()
555 (
556  const label irow,
557  const label jcol
558 ) const
559 {
560  #ifdef FULLDEBUG
561  checki(irow);
562  checkj(jcol);
563  #endif
564  return v_[irow*nCols_ + jcol];
565 }
566 
567 
568 template<class Form, class Type>
570 (
571  const label irow,
572  const label jcol
573 )
574 {
575  #ifdef FULLDEBUG
576  checki(irow);
577  checkj(jcol);
578  #endif
579  return v_[irow*nCols_ + jcol];
580 }
581 
582 
583 template<class Form, class Type>
584 inline const Type* Foam::Matrix<Form, Type>::operator[](const label irow) const
585 {
586  #ifdef FULLDEBUG
587  checki(irow);
588  #endif
589  return v_ + irow*nCols_;
590 }
591 
592 
593 template<class Form, class Type>
594 inline Type* Foam::Matrix<Form, Type>::operator[](const label irow)
595 {
596  #ifdef FULLDEBUG
597  checki(irow);
598  #endif
599  return v_ + irow*nCols_;
600 }
601 
602 
603 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
604 
605 namespace Foam
606 {
607 
608 // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
609 
610 //- Matrix-vector multiplication (A * x), where x is a column vector
611 template<class Form, class Type>
612 inline tmp<Field<Type>> operator*
613 (
614  const Matrix<Form, Type>& mat,
615  const UList<Type>& x
616 )
617 {
618  return mat.Amul(x);
619 }
620 
621 
622 //- Matrix-vector multiplication (A * x), where x is a column vector
623 template<class Form, class Type, class Addr>
624 inline tmp<Field<Type>> operator*
625 (
626  const Matrix<Form, Type>& mat,
628 )
629 {
630  return mat.Amul(x);
631 }
632 
633 
634 //- Vector-Matrix multiplication (x * A), where x is a row vector
635 template<class Form, class Type>
636 inline tmp<Field<Type>> operator*
637 (
638  const UList<Type>& x,
639  const Matrix<Form, Type>& mat
640 )
641 {
642  return mat.Tmul(x);
643 }
644 
645 
646 //- Vector-Matrix multiplication (x * A), where x is a row vector
647 template<class Form, class Type, class Addr>
648 inline tmp<Field<Type>> operator*
649 (
651  const Matrix<Form, Type>& mat
652 )
653 {
654  return mat.Tmul(x);
655 }
656 
657 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
658 
659 } // End namespace Foam
660 
661 // ************************************************************************* //
void shallowResize(const label m, const label n)
Resize Matrix without reallocating storage (unsafe)
Definition: MatrixI.H:445
ConstMatrixBlock< mType > subColumn(const label colIndex, const label rowIndex=0, label len=-1) const
Return const column or column&#39;s subset of Matrix.
Definition: MatrixI.H:252
ConstMatrixBlock< mType > subRow(const label rowIndex, const label colIndex=0, label len=-1) const
Return const row or const row&#39;s subset of Matrix.
Definition: MatrixI.H:277
autoPtr< mType > clone() const
Clone.
Definition: MatrixI.H:74
constexpr Matrix() noexcept
Default construct (empty matrix)
Definition: MatrixI.H:43
patchWriters resize(patchIds.size())
tmp< Field< Type > > Tmul(const UList< Type > &x) const
Left-multiply Matrix by a row vector (x * A)
Definition: MatrixI.H:475
error FatalError
Error stream (stdout output on all processes), with additional &#39;FOAM FATAL ERROR&#39; header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:652
constexpr char nl
The newline &#39;\n&#39; character (0x0a)
Definition: Ostream.H:50
tmp< Field< Type > > Amul(const UList< Type > &x) const
Right-multiply Matrix by a column vector (A * x)
Definition: MatrixI.H:454
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
A templated block of an (m x n) matrix of type <MatrixType>.
Definition: Matrix.H:62
void setSize(const label m, const label n)
Change Matrix dimensions, preserving the elements.
Definition: MatrixI.H:438
const Type * operator[](const label irow) const
Return const pointer to data in the specified row - rowData().
Definition: MatrixI.H:577
enum uniformity check_uniformity(InputIt first, InputIt last)
Algorithm to determine list/container uniformity.
Definition: ListPolicy.H:372
iterator begin() noexcept
Return an iterator to begin traversing a Matrix.
Definition: MatrixI.H:498
Base for lists with indirect addressing, templated on the list contents type and the addressing type...
bool empty() const noexcept
Return true if Matrix is empty (i.e., size() is zero)
Definition: MatrixI.H:97
ConstMatrixBlock< mType > subMatrix(const label rowIndex, const label colIndex, label szRows=-1, label szCols=-1) const
Return const sub-block of Matrix.
Definition: MatrixI.H:302
char * data_bytes() noexcept
Return pointer to the underlying array serving as data storage, reinterpreted as byte data...
Definition: MatrixI.H:186
labelPair sizes() const noexcept
Return row/column sizes.
Definition: MatrixI.H:90
const Type & at(const label idx) const
Linear addressing const element access.
Definition: MatrixI.H:220
Type * iterator
Random access iterator for traversing a Matrix.
Definition: Matrix.H:138
bool uniform() const
True if all entries have identical values, and Matrix is non-empty.
Definition: MatrixI.H:153
iterator end() noexcept
Return an iterator to end traversing a Matrix.
Definition: MatrixI.H:506
Container (non-empty) with identical values.
Definition: ListPolicy.H:363
label size() const noexcept
The number of elements in Matrix (m*n)
Definition: MatrixI.H:83
errorManip< error > abort(error &err)
Definition: errorManip.H:139
const Type * cdata() const noexcept
Return const pointer to the first data element, which can also be used to address into Matrix content...
Definition: MatrixI.H:165
void checki(const label irow) const
Check index i is within valid range [0, m)
Definition: MatrixI.H:104
A templated (m x n) matrix of objects of <T>. The layout is (mRows x nCols) - row-major order: ...
const direction noexcept
Definition: scalarImpl.H:265
std::streamsize size_bytes() const noexcept
Number of contiguous bytes for the Matrix data, no runtime check that the type is actually contiguous...
Definition: MatrixI.H:193
Pair< label > labelPair
A pair of labels.
Definition: Pair.H:54
const char * cdata_bytes() const noexcept
Return pointer to the underlying array serving as data storage, reinterpreted as byte data...
Definition: MatrixI.H:179
const Type * const_iterator
Random access iterator for traversing a Matrix.
Definition: Matrix.H:143
const_iterator cend() const noexcept
Return const_iterator to end traversing a constant Matrix.
Definition: MatrixI.H:522
const Type * rowData(const label irow) const
Return const pointer to data in the specified row.
Definition: MatrixI.H:200
void checkj(const label jcol) const
Check index j is within valid range [0, n)
Definition: MatrixI.H:122
Type * data() noexcept
Return pointer to the first data element, which can also be used to address into Matrix contents...
Definition: MatrixI.H:172
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
Definition: zero.H:57
label n
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
Definition: HashPtrTable.H:48
const_iterator cbegin() const noexcept
Return const_iterator to begin traversing a constant Matrix.
Definition: MatrixI.H:514
A class for managing temporary objects.
Definition: HashPtrTable.H:50
ConstMatrixBlock< mType > block(const label rowIndex, const label colIndex) const
Access Field as a ConstMatrixBlock.
Namespace for OpenFOAM.
void checkSize() const
Check that dimensions are positive, non-zero.
Definition: MatrixI.H:140