MatrixI.H 13.4 KB
Newer Older
1 2 3 4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
OpenFOAM bot's avatar
OpenFOAM bot committed
5
    \\  /    A nd           | www.openfoam.com
6
     \\/     M anipulation  |
OpenFOAM bot's avatar
OpenFOAM bot committed
7
-------------------------------------------------------------------------------
OpenFOAM bot's avatar
OpenFOAM bot committed
8 9
    Copyright (C) 2011-2016 OpenFOAM Foundation
    Copyright (C) 2019 OpenCFD Ltd.
10 11 12 13
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

14 15 16 17
    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.
18 19 20 21 22 23 24

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    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 30
#include "MatrixBlock.H"

31 32 33 34 35 36 37 38 39 40 41 42 43 44
// * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * * //

template<class Form, class Type>
inline void Foam::Matrix<Form, Type>::doAlloc()
{
    const label len = size();

    if (len)
    {
        v_ = new Type[len];
    }
}


45 46
// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

47 48
template<class Form, class Type>
inline Foam::Matrix<Form, Type>::Matrix()
49
:
50
    mRows_(0),
51
    nCols_(0),
52
    v_(nullptr)
53 54 55
{}


56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76
template<class Form, class Type>
inline Foam::Matrix<Form, Type>::Matrix(const labelPair& dims)
:
    Matrix<Form, Type>(dims.first(), dims.second())
{}


template<class Form, class Type>
inline Foam::Matrix<Form, Type>::Matrix(const labelPair& dims, const zero)
:
    Matrix<Form, Type>(dims.first(), dims.second(), Zero)
{}


template<class Form, class Type>
inline Foam::Matrix<Form, Type>::Matrix(const labelPair& dims, const Type& val)
:
    Matrix<Form, Type>(dims.first(), dims.second(), val)
{}


77
template<class Form, class Type>
78 79
inline Foam::autoPtr<Foam::Matrix<Form, Type>>
Foam::Matrix<Form, Type>::clone() const
80
{
81
    return autoPtr<Matrix<Form, Type>>::New(*this);
82 83 84 85 86
}


// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //

87 88 89
template<class Form, class Type>
inline const Foam::Matrix<Form, Type>& Foam::Matrix<Form, Type>::null()
{
90
    return NullObjectRef<Matrix<Form, Type>>();
91 92 93
}


94
template<class Form, class Type>
95
inline Foam::label Foam::Matrix<Form, Type>::m() const noexcept
96
{
97
    return mRows_;
98 99 100
}


101
template<class Form, class Type>
102
inline Foam::label Foam::Matrix<Form, Type>::n() const noexcept
103
{
104
    return nCols_;
105 106 107
}


108 109
template<class Form, class Type>
inline Foam::label Foam::Matrix<Form, Type>::size() const
110
{
111 112 113 114
    return mRows_ * nCols_;
}


115 116 117 118 119 120 121
template<class Form, class Type>
inline Foam::labelPair Foam::Matrix<Form, Type>::sizes() const
{
    return labelPair(mRows_, nCols_);
}


122 123 124 125
template<class Form, class Type>
inline bool Foam::Matrix<Form, Type>::empty() const noexcept
{
    return !mRows_ || !nCols_;
126 127 128
}


129 130
template<class Form, class Type>
inline void Foam::Matrix<Form, Type>::checki(const label i) const
131
{
132
    if (!mRows_ || !nCols_)
133
    {
134
        FatalErrorInFunction
135
            << "Attempt to access element from empty matrix"
136 137
            << abort(FatalError);
    }
138
    if (i < 0 || mRows_ <= i)
139
    {
140
        FatalErrorInFunction
141
            << "Index " << i << " out of range 0 ... " << mRows_-1
142 143 144 145 146
            << abort(FatalError);
    }
}


147 148
template<class Form, class Type>
inline void Foam::Matrix<Form, Type>::checkj(const label j) const
149
{
150
    if (!mRows_ || !nCols_)
151
    {
152
        FatalErrorInFunction
153
            << "Attempt to access element from empty matrix"
154 155
            << abort(FatalError);
    }
156
    if (j < 0 || nCols_ <= j)
157
    {
158
        FatalErrorInFunction
159
            << "index " << j << " out of range 0 ... " << nCols_-1
160 161 162 163 164
            << abort(FatalError);
    }
}


165
template<class Form, class Type>
166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187
inline void Foam::Matrix<Form, Type>::checkSize() const
{
    if (mRows_ < 0 || nCols_ < 0)
    {
        FatalErrorInFunction
            << "Incorrect size (" << mRows_ << ", " << nCols_ << ')' << nl
            << abort(FatalError);
    }
    // Could also check for odd sizes, like (0, N) and make (0,0)
}


template<class Form, class Type>
inline bool Foam::Matrix<Form, Type>::uniform() const
{
    const label len = size();

    if (len == 0)
    {
        return false;
    }

188
    for (label idx = 1; idx < len; ++idx)
189 190 191 192 193 194 195 196 197 198 199 200 201
    {
        if (v_[0] != v_[idx])
        {
            return false;
        }
    }

    return true;
}


template<class Form, class Type>
inline const Type* Foam::Matrix<Form, Type>::cdata() const noexcept
202 203 204 205 206 207
{
    return v_;
}


template<class Form, class Type>
208
inline Type* Foam::Matrix<Form, Type>::data() noexcept
209 210 211 212 213
{
    return v_;
}


214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237
template<class Form, class Type>
inline const Type* Foam::Matrix<Form, Type>::rowData(const label irow) const
{
    #ifdef FULLDEBUG
    checki(irow);
    #endif
    return (v_ + irow*nCols_);
}


template<class Form, class Type>
inline Type* Foam::Matrix<Form, Type>::rowData(const label irow)
{
    #ifdef FULLDEBUG
    checki(irow);
    #endif
    return (v_ + irow*nCols_);
}


template<class Form, class Type>
inline const Type& Foam::Matrix<Form, Type>::at(const label idx) const
{
    #ifdef FULLDEBUG
238
    if (idx < 0 || this->size() <= idx)
239 240 241 242 243 244
    {
        FatalErrorInFunction
            << "index " << idx << " out of range 0 ... " << this->size()
            << abort(FatalError);
    }
    #endif
245
    return *(v_ + idx);
246 247 248 249 250 251 252
}


template<class Form, class Type>
inline Type& Foam::Matrix<Form, Type>::at(const label idx)
{
    #ifdef FULLDEBUG
253
    if (idx < 0 || this->size() <= idx)
254 255 256 257 258 259
    {
        FatalErrorInFunction
            << "index " << idx << " out of range 0 ... " << this->size()
            << abort(FatalError);
    }
    #endif
260
    return *(v_ + idx);
261 262 263
}


264
template<class Form, class Type>
265
inline Foam::ConstMatrixBlock<Foam::Matrix<Form, Type>>
266
Foam::Matrix<Form, Type>::subColumn
267
(
268 269 270
    const label colIndex,
    const label rowIndex,
    label len
271 272
) const
{
273 274 275 276 277
    if (len < 0)
    {
        len = mRows_ - rowIndex;
    }

278
    return ConstMatrixBlock<mType>
279 280
    (
        *this,
281 282 283 284
        len, // rows
        1,
        rowIndex,
        colIndex
285 286 287 288 289
    );
}


template<class Form, class Type>
290
inline Foam::ConstMatrixBlock<Foam::Matrix<Form, Type>>
291
Foam::Matrix<Form, Type>::subRow
292
(
293 294 295
    const label rowIndex,
    const label colIndex,
    label len
296 297
) const
{
298 299 300 301 302
    if (len < 0)
    {
        len = nCols_ - colIndex;
    }

303
    return ConstMatrixBlock<mType>
304 305
    (
        *this,
306 307 308 309
        1,
        len, // columns
        rowIndex,
        colIndex
310 311 312 313 314
    );
}


template<class Form, class Type>
315
inline Foam::ConstMatrixBlock<Foam::Matrix<Form, Type>>
316
Foam::Matrix<Form, Type>::subMatrix
317
(
318 319 320 321
    const label rowIndex,
    const label colIndex,
    label szRows,
    label szCols
322 323
) const
{
324 325 326
    if (szRows < 0) szRows = mRows_ - rowIndex;
    if (szCols < 0) szCols = nCols_ - colIndex;

327
    return ConstMatrixBlock<mType>
328 329
    (
        *this,
330 331 332 333
        szRows,
        szCols,
        rowIndex,
        colIndex
334 335 336 337 338
    );
}


template<class Form, class Type>
339
template<class VectorSpace>
340
inline Foam::ConstMatrixBlock<Foam::Matrix<Form, Type>>
341
Foam::Matrix<Form, Type>::block
342
(
343 344
    const label rowIndex,
    const label colIndex
345 346
) const
{
347
    return ConstMatrixBlock<mType>
348 349
    (
        *this,
350 351 352 353
        VectorSpace::mRows,
        VectorSpace::nCols,
        rowIndex,
        colIndex
354 355 356 357 358
    );
}


template<class Form, class Type>
359
inline Foam::MatrixBlock<Foam::Matrix<Form, Type>>
360
Foam::Matrix<Form, Type>::subColumn
361
(
362 363 364
    const label colIndex,
    const label rowIndex,
    label len
365 366
)
{
367 368 369 370 371
    if (len < 0)
    {
        len = mRows_ - rowIndex;
    }

372
    return MatrixBlock<mType>
373 374
    (
        *this,
375 376 377 378
        len, // rows
        1,
        rowIndex,
        colIndex
379 380 381 382 383
    );
}


template<class Form, class Type>
384
inline Foam::MatrixBlock<Foam::Matrix<Form, Type>>
385 386 387 388 389 390
Foam::Matrix<Form, Type>::subRow
(
    const label rowIndex,
    const label colIndex,
    label len
)
391
{
392 393 394 395 396
    if (len < 0)
    {
        len = nCols_ - colIndex;
    }

397
    return MatrixBlock<mType>
398 399
    (
        *this,
400 401 402 403
        1,
        len, // columns
        rowIndex,
        colIndex
404 405 406 407 408
    );
}


template<class Form, class Type>
409
inline Foam::MatrixBlock<Foam::Matrix<Form, Type>>
410 411 412 413 414 415 416
Foam::Matrix<Form, Type>::subMatrix
(
    const label rowIndex,
    const label colIndex,
    label szRows,
    label szCols
)
417
{
418 419 420
    if (szRows < 0) szRows = mRows_ - rowIndex;
    if (szCols < 0) szCols = nCols_ - colIndex;

421
    return MatrixBlock<mType>
422 423
    (
        *this,
424 425 426 427
        szRows,
        szCols,
        rowIndex,
        colIndex
428 429 430 431 432
    );
}


template<class Form, class Type>
433
template<class VectorSpace>
434
inline Foam::MatrixBlock<Foam::Matrix<Form, Type>>
435
Foam::Matrix<Form, Type>::block
436
(
437 438
    const label rowIndex,
    const label colIndex
439 440
)
{
441
    return MatrixBlock<mType>
442 443
    (
        *this,
444 445 446 447
        VectorSpace::mRows,
        VectorSpace::nCols,
        rowIndex,
        colIndex
448 449 450 451
    );
}


452 453 454 455 456 457 458
template<class Form, class Type>
inline void Foam::Matrix<Form, Type>::setSize(const label m, const label n)
{
    resize(m, n);
}


459 460 461 462 463 464 465 466
template<class Form, class Type>
void Foam::Matrix<Form, Type>::shallowResize(const label m, const label n)
{
    mRows_ = m;
    nCols_ = n;
}


467
template<class Form, class Type>
468
inline Foam::tmp<Foam::Field<Type>> Foam::Matrix<Form, Type>::Amul
469
(
470
    const UList<Type>& x
471 472
) const
{
473
    return this->AmulImpl(x);
474 475 476 477 478
}


template<class Form, class Type>
template<class Addr>
479
inline Foam::tmp<Foam::Field<Type>> Foam::Matrix<Form, Type>::Amul
480
(
481
    const IndirectListBase<Type, Addr>& x
482 483
) const
{
484
    return this->AmulImpl(x);
485 486 487 488
}


template<class Form, class Type>
489
inline Foam::tmp<Foam::Field<Type>> Foam::Matrix<Form, Type>::Tmul
490
(
491
    const UList<Type>& x
492 493
) const
{
494
    return this->TmulImpl(x);
495 496 497 498 499
}


template<class Form, class Type>
template<class Addr>
500
inline Foam::tmp<Foam::Field<Type>> Foam::Matrix<Form, Type>::Tmul
501
(
502
    const IndirectListBase<Type, Addr>& x
503 504
) const
{
505
    return this->TmulImpl(x);
506 507 508
}


509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556 557 558
// * * * * * * * * * * * * * * * * * Iterators * * * * * * * * * * * * * * * //

template<class Form, class Type>
inline typename Foam::Matrix<Form, Type>::iterator
Foam::Matrix<Form, Type>::begin()
{
    return v_;
}


template<class Form, class Type>
inline typename Foam::Matrix<Form, Type>::iterator
Foam::Matrix<Form, Type>::end()
{
    return v_ + (mRows_ * nCols_);
}


template<class Form, class Type>
inline typename Foam::Matrix<Form, Type>::const_iterator
Foam::Matrix<Form, Type>::cbegin() const
{
    return v_;
}


template<class Form, class Type>
inline typename Foam::Matrix<Form, Type>::const_iterator
Foam::Matrix<Form, Type>::cend() const
{
    return v_ + (mRows_ * nCols_);
}


template<class Form, class Type>
inline typename Foam::Matrix<Form, Type>::const_iterator
Foam::Matrix<Form, Type>::begin() const
{
    return v_;
}


template<class Form, class Type>
inline typename Foam::Matrix<Form, Type>::const_iterator
Foam::Matrix<Form, Type>::end() const
{
    return v_ + (mRows_ * nCols_);
}


559 560
// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //

561 562 563
template<class Form, class Type>
inline const Type& Foam::Matrix<Form, Type>::operator()
(
564 565
    const label irow,
    const label jcol
566
) const
567
{
568
    #ifdef FULLDEBUG
569 570
    checki(irow);
    checkj(jcol);
571
    #endif
572
    return v_[irow*nCols_ + jcol];
573 574 575 576 577 578
}


template<class Form, class Type>
inline Type& Foam::Matrix<Form, Type>::operator()
(
579 580
    const label irow,
    const label jcol
581 582
)
{
583
    #ifdef FULLDEBUG
584 585
    checki(irow);
    checkj(jcol);
586
    #endif
587
    return v_[irow*nCols_ + jcol];
588 589 590
}


591
template<class Form, class Type>
592
inline const Type* Foam::Matrix<Form, Type>::operator[](const label irow) const
593
{
594
    #ifdef FULLDEBUG
595
    checki(irow);
596
    #endif
597
    return v_ + irow*nCols_;
598 599 600 601
}


template<class Form, class Type>
602
inline Type* Foam::Matrix<Form, Type>::operator[](const label irow)
603
{
604
    #ifdef FULLDEBUG
605
    checki(irow);
606
    #endif
607
    return v_ + irow*nCols_;
608 609 610
}


611 612 613 614 615 616 617 618
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace Foam
{

// * * * * * * * * * * * * * * * Global Operators  * * * * * * * * * * * * * //

//- Matrix-vector multiplication (A * x), where x is a column vector
619
template<class Form, class Type>
620
inline tmp<Field<Type>> operator*
621 622 623 624 625
(
    const Matrix<Form, Type>& mat,
    const UList<Type>& x
)
{
626
    return mat.Amul(x);
627 628 629
}


630
//- Matrix-vector multiplication (A * x), where x is a column vector
631
template<class Form, class Type, class Addr>
632
inline tmp<Field<Type>> operator*
633 634 635 636 637
(
    const Matrix<Form, Type>& mat,
    const IndirectListBase<Type, Addr>& x
)
{
638
    return mat.Amul(x);
639 640 641
}


642
//- Vector-Matrix multiplication (x * A), where x is a row vector
643
template<class Form, class Type>
644
inline tmp<Field<Type>> operator*
645 646 647 648 649
(
    const UList<Type>& x,
    const Matrix<Form, Type>& mat
)
{
650
    return mat.Tmul(x);
651 652 653
}


654
//- Vector-Matrix multiplication (x * A), where x is a row vector
655
template<class Form, class Type, class Addr>
656
inline tmp<Field<Type>> operator*
657 658 659 660 661
(
    const IndirectListBase<Type, Addr>& x,
    const Matrix<Form, Type>& mat
)
{
662
    return mat.Tmul(x);
663 664
}

665 666 667
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace Foam
668

669
// ************************************************************************* //