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  v_ = new Type[len];
42  }
43 }
44 
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 template<class Form, class Type>
50 :
51  mRows_(0),
52  nCols_(0),
53  v_(nullptr)
54 {}
55 
56 
57 template<class Form, class Type>
59 :
60  Matrix<Form, Type>(dims.first(), dims.second())
61 {}
62 
63 
64 template<class Form, class Type>
66 :
67  Matrix<Form, Type>(dims.first(), dims.second(), Zero)
68 {}
69 
70 
71 template<class Form, class Type>
72 inline Foam::Matrix<Form, Type>::Matrix(const labelPair& dims, const Type& val)
73 :
74  Matrix<Form, Type>(dims.first(), dims.second(), val)
75 {}
76 
77 
78 template<class Form, class Type>
81 {
83 }
84 
85 
86 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
87 
88 template<class Form, class Type>
90 {
91  return NullObjectRef<Matrix<Form, Type>>();
92 }
93 
94 
95 template<class Form, class Type>
96 inline Foam::label Foam::Matrix<Form, Type>::m() const noexcept
97 {
98  return mRows_;
99 }
100 
101 
102 template<class Form, class Type>
103 inline Foam::label Foam::Matrix<Form, Type>::n() const noexcept
104 {
105  return nCols_;
106 }
107 
108 
109 template<class Form, class Type>
110 inline Foam::label Foam::Matrix<Form, Type>::size() const
111 {
112  return mRows_ * nCols_;
113 }
114 
115 
116 template<class Form, class Type>
118 {
119  return labelPair(mRows_, nCols_);
120 }
121 
122 
123 template<class Form, class Type>
125 {
126  return !mRows_ || !nCols_;
127 }
128 
129 
130 template<class Form, class Type>
131 inline void Foam::Matrix<Form, Type>::checki(const label i) const
132 {
133  if (!mRows_ || !nCols_)
134  {
136  << "Attempt to access element from empty matrix"
137  << abort(FatalError);
138  }
139  if (i < 0 || mRows_ <= i)
140  {
142  << "Index " << i << " out of range 0 ... " << mRows_-1
143  << abort(FatalError);
144  }
145 }
146 
147 
148 template<class Form, class Type>
149 inline void Foam::Matrix<Form, Type>::checkj(const label j) const
150 {
151  if (!mRows_ || !nCols_)
152  {
154  << "Attempt to access element from empty matrix"
155  << abort(FatalError);
156  }
157  if (j < 0 || nCols_ <= j)
158  {
160  << "index " << j << " out of range 0 ... " << nCols_-1
161  << abort(FatalError);
162  }
163 }
164 
165 
166 template<class Form, class Type>
167 inline void Foam::Matrix<Form, Type>::checkSize() const
168 {
169  if (mRows_ < 0 || nCols_ < 0)
170  {
172  << "Incorrect size (" << mRows_ << ", " << nCols_ << ')' << nl
174  }
175  // Could also check for odd sizes, like (0, N) and make (0,0)
176 }
177 
178 
179 template<class Form, class Type>
180 inline bool Foam::Matrix<Form, Type>::uniform() const
181 {
182  const label len = size();
183 
184  if (len == 0)
185  {
186  return false;
187  }
188 
189  for (label idx = 1; idx < len; ++idx)
190  {
191  if (v_[0] != v_[idx])
192  {
193  return false;
194  }
195  }
196 
197  return true;
198 }
199 
200 
201 template<class Form, class Type>
202 inline const Type* Foam::Matrix<Form, Type>::cdata() const noexcept
203 {
204  return v_;
205 }
206 
207 
208 template<class Form, class Type>
210 {
211  return v_;
212 }
213 
214 
215 template<class Form, class Type>
217 {
218  return reinterpret_cast<const char*>(v_);
219 }
220 
221 
222 template<class Form, class Type>
224 {
225  return reinterpret_cast<char*>(v_);
226 }
227 
228 
229 template<class Form, class Type>
230 inline std::streamsize Foam::Matrix<Form, Type>::size_bytes() const noexcept
231 {
232  return mRows_*nCols_*sizeof(Type);
233 }
234 
235 
236 template<class Form, class Type>
237 inline const Type* Foam::Matrix<Form, Type>::rowData(const label irow) const
238 {
239  #ifdef FULLDEBUG
240  checki(irow);
241  #endif
242  return (v_ + irow*nCols_);
243 }
244 
245 
246 template<class Form, class Type>
247 inline Type* Foam::Matrix<Form, Type>::rowData(const label irow)
248 {
249  #ifdef FULLDEBUG
250  checki(irow);
251  #endif
252  return (v_ + irow*nCols_);
253 }
254 
255 
256 template<class Form, class Type>
257 inline const Type& Foam::Matrix<Form, Type>::at(const label idx) const
258 {
259  #ifdef FULLDEBUG
260  if (idx < 0 || this->size() <= idx)
261  {
263  << "index " << idx << " out of range 0 ... " << this->size()
264  << abort(FatalError);
265  }
266  #endif
267  return *(v_ + idx);
268 }
269 
270 
271 template<class Form, class Type>
272 inline Type& Foam::Matrix<Form, Type>::at(const label idx)
273 {
274  #ifdef FULLDEBUG
275  if (idx < 0 || this->size() <= idx)
276  {
278  << "index " << idx << " out of range 0 ... " << this->size()
279  << abort(FatalError);
280  }
281  #endif
282  return *(v_ + idx);
283 }
284 
285 
286 template<class Form, class Type>
289 (
290  const label colIndex,
291  const label rowIndex,
292  label len
293 ) const
294 {
295  if (len < 0)
296  {
297  len = mRows_ - rowIndex;
298  }
299 
300  return ConstMatrixBlock<mType>
301  (
302  *this,
303  len, // rows
304  1,
305  rowIndex,
306  colIndex
307  );
308 }
309 
310 
311 template<class Form, class Type>
314 (
315  const label rowIndex,
316  const label colIndex,
317  label len
318 ) const
319 {
320  if (len < 0)
321  {
322  len = nCols_ - colIndex;
323  }
324 
325  return ConstMatrixBlock<mType>
326  (
327  *this,
328  1,
329  len, // columns
330  rowIndex,
331  colIndex
332  );
333 }
334 
335 
336 template<class Form, class Type>
339 (
340  const label rowIndex,
341  const label colIndex,
342  label szRows,
343  label szCols
344 ) const
345 {
346  if (szRows < 0) szRows = mRows_ - rowIndex;
347  if (szCols < 0) szCols = nCols_ - colIndex;
348 
350  (
351  *this,
352  szRows,
353  szCols,
354  rowIndex,
355  colIndex
356  );
357 }
358 
359 
360 template<class Form, class Type>
361 template<class VectorSpace>
364 (
365  const label rowIndex,
366  const label colIndex
367 ) const
368 {
370  (
371  *this,
372  VectorSpace::mRows,
373  VectorSpace::nCols,
374  rowIndex,
375  colIndex
376  );
377 }
378 
379 
380 template<class Form, class Type>
383 (
384  const label colIndex,
385  const label rowIndex,
386  label len
387 )
388 {
389  if (len < 0)
390  {
391  len = mRows_ - rowIndex;
392  }
393 
394  return MatrixBlock<mType>
395  (
396  *this,
397  len, // rows
398  1,
399  rowIndex,
400  colIndex
401  );
402 }
403 
404 
405 template<class Form, class Type>
408 (
409  const label rowIndex,
410  const label colIndex,
411  label len
412 )
413 {
414  if (len < 0)
415  {
416  len = nCols_ - colIndex;
417  }
418 
419  return MatrixBlock<mType>
420  (
421  *this,
422  1,
423  len, // columns
424  rowIndex,
425  colIndex
426  );
427 }
428 
429 
430 template<class Form, class Type>
433 (
434  const label rowIndex,
435  const label colIndex,
436  label szRows,
437  label szCols
438 )
439 {
440  if (szRows < 0) szRows = mRows_ - rowIndex;
441  if (szCols < 0) szCols = nCols_ - colIndex;
442 
443  return MatrixBlock<mType>
444  (
445  *this,
446  szRows,
447  szCols,
448  rowIndex,
449  colIndex
450  );
451 }
452 
453 
454 template<class Form, class Type>
455 template<class VectorSpace>
458 (
459  const label rowIndex,
460  const label colIndex
461 )
462 {
463  return MatrixBlock<mType>
464  (
465  *this,
466  VectorSpace::mRows,
467  VectorSpace::nCols,
468  rowIndex,
469  colIndex
470  );
471 }
472 
473 
474 template<class Form, class Type>
475 inline void Foam::Matrix<Form, Type>::setSize(const label m, const label n)
476 {
477  resize(m, n);
478 }
479 
480 
481 template<class Form, class Type>
482 void Foam::Matrix<Form, Type>::shallowResize(const label m, const label n)
483 {
484  mRows_ = m;
485  nCols_ = n;
486 }
487 
488 
489 template<class Form, class Type>
491 (
492  const UList<Type>& x
493 ) const
494 {
495  return this->AmulImpl(x);
496 }
497 
498 
499 template<class Form, class Type>
500 template<class Addr>
502 (
504 ) const
505 {
506  return this->AmulImpl(x);
507 }
508 
509 
510 template<class Form, class Type>
512 (
513  const UList<Type>& x
514 ) const
515 {
516  return this->TmulImpl(x);
517 }
518 
519 
520 template<class Form, class Type>
521 template<class Addr>
523 (
525 ) const
526 {
527  return this->TmulImpl(x);
528 }
529 
530 
531 // * * * * * * * * * * * * * * * * * Iterators * * * * * * * * * * * * * * * //
532 
533 template<class Form, class Type>
536 {
537  return v_;
538 }
539 
540 
541 template<class Form, class Type>
544 {
545  return v_ + (mRows_ * nCols_);
546 }
547 
548 
549 template<class Form, class Type>
552 {
553  return v_;
554 }
555 
556 
557 template<class Form, class Type>
560 {
561  return v_ + (mRows_ * nCols_);
562 }
563 
564 
565 template<class Form, class Type>
568 {
569  return v_;
570 }
571 
572 
573 template<class Form, class Type>
576 {
577  return v_ + (mRows_ * nCols_);
578 }
579 
580 
581 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
582 
583 template<class Form, class Type>
584 inline const Type& Foam::Matrix<Form, Type>::operator()
585 (
586  const label irow,
587  const label jcol
588 ) const
589 {
590  #ifdef FULLDEBUG
591  checki(irow);
592  checkj(jcol);
593  #endif
594  return v_[irow*nCols_ + jcol];
595 }
596 
597 
598 template<class Form, class Type>
600 (
601  const label irow,
602  const label jcol
603 )
604 {
605  #ifdef FULLDEBUG
606  checki(irow);
607  checkj(jcol);
608  #endif
609  return v_[irow*nCols_ + jcol];
610 }
611 
612 
613 template<class Form, class Type>
614 inline const Type* Foam::Matrix<Form, Type>::operator[](const label irow) const
615 {
616  #ifdef FULLDEBUG
617  checki(irow);
618  #endif
619  return v_ + irow*nCols_;
620 }
621 
622 
623 template<class Form, class Type>
624 inline Type* Foam::Matrix<Form, Type>::operator[](const label irow)
625 {
626  #ifdef FULLDEBUG
627  checki(irow);
628  #endif
629  return v_ + irow*nCols_;
630 }
631 
632 
633 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
634 
635 namespace Foam
636 {
637 
638 // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
639 
640 //- Matrix-vector multiplication (A * x), where x is a column vector
641 template<class Form, class Type>
642 inline tmp<Field<Type>> operator*
643 (
644  const Matrix<Form, Type>& mat,
645  const UList<Type>& x
646 )
647 {
648  return mat.Amul(x);
649 }
650 
651 
652 //- Matrix-vector multiplication (A * x), where x is a column vector
653 template<class Form, class Type, class Addr>
654 inline tmp<Field<Type>> operator*
655 (
656  const Matrix<Form, Type>& mat,
658 )
659 {
660  return mat.Amul(x);
661 }
662 
663 
664 //- Vector-Matrix multiplication (x * A), where x is a row vector
665 template<class Form, class Type>
666 inline tmp<Field<Type>> operator*
667 (
668  const UList<Type>& x,
669  const Matrix<Form, Type>& mat
670 )
671 {
672  return mat.Tmul(x);
673 }
674 
675 
676 //- Vector-Matrix multiplication (x * A), where x is a row vector
677 template<class Form, class Type, class Addr>
678 inline tmp<Field<Type>> operator*
679 (
681  const Matrix<Form, Type>& mat
682 )
683 {
684  return mat.Tmul(x);
685 }
686 
687 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
688 
689 } // End namespace Foam
690 
691 // ************************************************************************* //
void shallowResize(const label m, const label n)
Resize Matrix without reallocating storage (unsafe)
Definition: MatrixI.H:475
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:282
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:307
autoPtr< mType > clone() const
Clone.
Definition: MatrixI.H:73
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:505
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:578
constexpr char nl
The newline &#39;\n&#39; character (0x0a)
Definition: Ostream.H:49
tmp< Field< Type > > Amul(const UList< Type > &x) const
Right-multiply Matrix by a column vector (A * x)
Definition: MatrixI.H:484
label n() const noexcept
The number of columns.
Definition: MatrixI.H:96
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:468
const Type * operator[](const label irow) const
Return const pointer to data in the specified row - rowData().
Definition: MatrixI.H:607
iterator begin() noexcept
Return an iterator to begin traversing a Matrix.
Definition: MatrixI.H:528
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:117
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:332
char * data_bytes() noexcept
Return pointer to the underlying array serving as data storage, reinterpreted as byte data...
Definition: MatrixI.H:216
static const Matrix< Form, Type > & null()
Return a null Matrix.
Definition: MatrixI.H:82
label size() const
The number of elements in Matrix (m*n)
Definition: MatrixI.H:103
const Type & at(const label idx) const
Linear addressing const element access.
Definition: MatrixI.H:250
Type * iterator
Random access iterator for traversing a Matrix.
Definition: Matrix.H:137
labelPair sizes() const
Return row/column sizes.
Definition: MatrixI.H:110
bool uniform() const
True if all entries have identical values, and Matrix is non-empty.
Definition: MatrixI.H:173
iterator end() noexcept
Return an iterator to end traversing a Matrix.
Definition: MatrixI.H:536
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:195
label m() const noexcept
The number of rows.
Definition: MatrixI.H:89
void checki(const label irow) const
Check index i is within valid range [0, m)
Definition: MatrixI.H:124
A templated (m x n) matrix of objects of <T>. The layout is (mRows x nCols) - row-major order: ...
const direction noexcept
Definition: Scalar.H:258
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:223
Pair< label > labelPair
A pair of labels.
Definition: Pair.H:51
const char * cdata_bytes() const noexcept
Return pointer to the underlying array serving as data storage, reinterpreted as byte data...
Definition: MatrixI.H:209
const Type * const_iterator
Random access iterator for traversing a Matrix.
Definition: Matrix.H:142
const_iterator cend() const noexcept
Return const_iterator to end traversing a constant Matrix.
Definition: MatrixI.H:552
const Type * rowData(const label irow) const
Return const pointer to data in the specified row.
Definition: MatrixI.H:230
void checkj(const label jcol) const
Check index j is within valid range [0, n)
Definition: MatrixI.H:142
Matrix() noexcept
Default construct (empty matrix)
Definition: MatrixI.H:42
Type * data() noexcept
Return pointer to the first data element, which can also be used to address into Matrix contents...
Definition: MatrixI.H:202
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:544
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.
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:133
void checkSize() const
Check that dimensions are positive, non-zero.
Definition: MatrixI.H:160