src/OpenFOAM/primitives/VectorSpace/VectorSpaceI.H

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright held by original author
   6      \\/     M anipulation  |
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of OpenFOAM.
  10
  11     OpenFOAM is free software; you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     option) any later version.
  15
  16     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
  17     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  18     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  19     for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with OpenFOAM; if not, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 \*---------------------------------------------------------------------------*/
  26
  27 #include "error.H"
  28 #include "products.H"
  29 #include "VectorSpaceM.H"
  30 #include "ops.H"
  31 #include "pTraits.H"
  32
  33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  34
  35 namespace Foam
  36 {
  37
  38 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  39
  40 template<class Form, class Cmpt, int nCmpt>
  41 inline VectorSpace<Form, Cmpt, nCmpt>::VectorSpace()
  42 {}
  43
  44
  45 template<class Form, class Cmpt, int nCmpt>
  46 inline VectorSpace<Form, Cmpt, nCmpt>::VectorSpace
  47 (
  48     const VectorSpace<Form, Cmpt, nCmpt>& vs
  49 )
  50 {
  51     VectorSpaceOps<nCmpt,0>::eqOp(*this, vs, eqOp<Cmpt>());
  52 }
  53
  54
  55 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  56
  57 template<class Form, class Cmpt, int nCmpt>
  58 inline label VectorSpace<Form, Cmpt, nCmpt>::size() const
  59 {
  60     return nCmpt;
  61 }
  62
  63
  64 template<class Form, class Cmpt, int nCmpt>
  65 inline const Cmpt& VectorSpace<Form, Cmpt, nCmpt>::component
  66 (
  67     const direction d
  68 ) const
  69 {
  70 #   ifdef FULLDEBUG
  71     if (d >= nCmpt)
  72     {
  73         FatalErrorIn
  74         (
  75             "VectorSpace<Form, Cmpt, nCmpt>::component(direction) const"
  76         )   << "index out of range"
  77             << abort(FatalError);
  78     }
  79 #   endif
  80
  81     return v_[d];
  82 }
  83
  84
  85 template<class Form, class Cmpt, int nCmpt>
  86 inline Cmpt& VectorSpace<Form, Cmpt, nCmpt>::component
  87 (
  88     const direction d
  89 )
  90 {
  91 #   ifdef FULLDEBUG
  92     if (d >= nCmpt)
  93     {
  94         FatalErrorIn("VectorSpace<Form, Cmpt, nCmpt>::component(direction)")
  95             << "index out of range"
  96             << abort(FatalError);
  97     }
  98 #   endif
  99
 100     return v_[d];
 101 }
 102
 103
 104 template<class Form, class Cmpt, int nCmpt>
 105 inline void VectorSpace<Form, Cmpt, nCmpt>::component
 106 (
 107     Cmpt& c,
 108     const direction d
 109 ) const
 110 {
 111 #   ifdef FULLDEBUG
 112     if (d >= nCmpt)
 113     {
 114         FatalErrorIn
 115         (
 116             "VectorSpace<Form, Cmpt, nCmpt>::component(Cmpt&, direction) const"
 117         )   << "index out of range"
 118             << abort(FatalError);
 119     }
 120 #   endif
 121
 122     c = v_[d];
 123 }
 124
 125
 126 template<class Form, class Cmpt, int nCmpt>
 127 inline void VectorSpace<Form, Cmpt, nCmpt>::replace
 128 (
 129     const direction d,
 130     const Cmpt& c
 131 )
 132 {
 133 #   ifdef FULLDEBUG
 134     if (d >= nCmpt)
 135     {
 136         FatalErrorIn
 137         (
 138             "VectorSpace<Form, Cmpt, nCmpt>::"
 139             "replace(direction, const Cmpt&) const"
 140         )   << "index out of range"
 141             << abort(FatalError);
 142     }
 143 #   endif
 144
 145     v_[d] = c;
 146 }
 147
 148
 149 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 150
 151 template<class Form, class Cmpt, int nCmpt>
 152 inline const Cmpt& VectorSpace<Form, Cmpt, nCmpt>::operator[]
 153 (
 154     const direction d
 155 ) const
 156 {
 157 #   ifdef FULLDEBUG
 158     if (d >= nCmpt)
 159     {
 160         FatalErrorIn
 161         (
 162             "VectorSpace<Form, Cmpt, nCmpt>::operator[](direction d) const"
 163         )   << "index out of range"
 164             << abort(FatalError);
 165     }
 166 #   endif
 167
 168     return v_[d];
 169 }
 170
 171
 172 template<class Form, class Cmpt, int nCmpt>
 173 inline Cmpt& VectorSpace<Form, Cmpt, nCmpt>::operator[]
 174 (
 175     const direction d
 176 )
 177 {
 178 #   ifdef FULLDEBUG
 179     if (d >= nCmpt)
 180     {
 181         FatalErrorIn("VectorSpace<Form, Cmpt, nCmpt>::operator[](direction d)")
 182             << "index out of range"
 183             << abort(FatalError);
 184     }
 185 #   endif
 186
 187     return v_[d];
 188 }
 189
 190
 191 template<class Form, class Cmpt, int nCmpt>
 192 inline void VectorSpace<Form, Cmpt, nCmpt>::operator=
 193 (
 194     const VectorSpace<Form, Cmpt, nCmpt>& vs
 195 )
 196 {
 197     VectorSpaceOps<nCmpt,0>::eqOp(*this, vs, eqOp<Cmpt>());
 198 }
 199
 200
 201 template<class Form, class Cmpt, int nCmpt>
 202 inline void VectorSpace<Form, Cmpt, nCmpt>::operator+=
 203 (
 204     const VectorSpace<Form, Cmpt, nCmpt>& vs
 205 )
 206 {
 207     VectorSpaceOps<nCmpt,0>::eqOp(*this, vs, plusEqOp<Cmpt>());
 208 }
 209
 210
 211 template<class Form, class Cmpt, int nCmpt>
 212 inline void VectorSpace<Form, Cmpt, nCmpt>::operator-=
 213 (
 214     const VectorSpace<Form, Cmpt, nCmpt>& vs
 215 )
 216 {
 217     VectorSpaceOps<nCmpt,0>::eqOp(*this, vs, minusEqOp<Cmpt>());
 218 }
 219
 220
 221 template<class Form, class Cmpt, int nCmpt>
 222 inline void VectorSpace<Form, Cmpt, nCmpt>::operator*=
 223 (
 224     const scalar s
 225 )
 226 {
 227     VectorSpaceOps<nCmpt,0>::eqOpS(*this, s, multiplyEqOp2<Cmpt, scalar>());
 228 }
 229
 230
 231 template<class Form, class Cmpt, int nCmpt>
 232 inline void VectorSpace<Form, Cmpt, nCmpt>::operator/=
 233 (
 234     const scalar s
 235 )
 236 {
 237     VectorSpaceOps<nCmpt,0>::eqOpS(*this, s, divideEqOp2<Cmpt, scalar>());
 238 }
 239
 240
 241 // * * * * * * * * * * * * * * * Global Functions  * * * * * * * * * * * * * //
 242
 243 template<class Form, class Cmpt, int nCmpt>
 244 inline Cmpt& setComponent
 245 (
 246     VectorSpace<Form, Cmpt, nCmpt>& vs,
 247     const direction d
 248 )
 249 {
 250     return vs.component(d);
 251 }
 252
 253
 254 template<class Form, class Cmpt, int nCmpt>
 255 inline const Cmpt& component
 256 (
 257     const VectorSpace<Form, Cmpt, nCmpt>& vs,
 258     const direction d
 259 )
 260 {
 261     return vs.component(d);
 262 }
 263
 264
 265 // Powers of a Form
 266 // Equivalent to outer-products between the Form and itself
 267 // Form^0 = 1.0
 268 template<class Form, class Cmpt, int nCmpt>
 269 inline typename powProduct<Form, 0>::type pow
 270 (
 271     const VectorSpace<Form, Cmpt, nCmpt>&,
 272     typename powProduct<Form, 0>::type
 273     = pTraits<typename powProduct<Form, 0>::type>::zero
 274 )
 275 {
 276     return 1.0;
 277 }
 278
 279 // Form^1 = Form
 280 template<class Form, class Cmpt, int nCmpt>
 281 inline typename powProduct<Form, 1>::type pow
 282 (
 283     const VectorSpace<Form, Cmpt, nCmpt>& v,
 284     typename powProduct<Form, 1>::type
 285   = pTraits<typename powProduct<Form, 1>::type>::zero
 286 )
 287 {
 288     return static_cast<const Form&>(v);
 289 }
 290
 291
 292 // Form^2 = sqr(Form)
 293 template<class Form, class Cmpt, int nCmpt>
 294 inline typename powProduct<Form, 2>::type pow
 295 (
 296     const VectorSpace<Form, Cmpt, nCmpt>& v,
 297     typename powProduct<Form, 2>::type
 298   = pTraits<typename powProduct<Form, 2>::type>::zero
 299 )
 300 {
 301     return sqr(static_cast<const Form&>(v));
 302 }
 303
 304
 305 template<class Form, class Cmpt, int nCmpt>
 306 inline scalar magSqr
 307 (
 308     const VectorSpace<Form, Cmpt, nCmpt>& vs
 309 )
 310 {
 311     scalar ms = magSqr(vs.v_[0]);
 312     VectorSpaceOps<nCmpt,1>::SeqOp(ms, vs, plusEqMagSqrOp2<scalar, Cmpt>());
 313     return ms;
 314 }
 315
 316
 317 template<class Form, class Cmpt, int nCmpt>
 318 inline scalar mag
 319 (
 320     const VectorSpace<Form, Cmpt, nCmpt>& vs
 321 )
 322 {
 323     return ::sqrt(magSqr(static_cast<const Form&>(vs)));
 324 }
 325
 326
 327 template<class Form, class Cmpt, int nCmpt>
 328 inline VectorSpace<Form, Cmpt, nCmpt> cmptMultiply
 329 (
 330     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 331     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 332 )
 333 {
 334     Form v;
 335     VectorSpaceOps<nCmpt,0>::op(v, vs1, vs2, cmptMultiplyOp<Cmpt>());
 336     return v;
 337 }
 338
 339
 340 template<class Form, class Cmpt, int nCmpt>
 341 inline VectorSpace<Form, Cmpt, nCmpt> cmptDivide
 342 (
 343     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 344     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 345 )
 346 {
 347     Form v;
 348     VectorSpaceOps<nCmpt,0>::op(v, vs1, vs2, cmptDivideOp<Cmpt>());
 349     return v;
 350 }
 351
 352
 353 template<class Form, class Cmpt, int nCmpt>
 354 inline Cmpt cmptSumMultiply
 355 (
 356     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 357     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 358 )
 359 {
 360     return cmptSum(cmptMultiply(vs1, vs2));
 361 }
 362
 363
 364 template<class Form, class Cmpt, int nCmpt>
 365 inline VectorSpace<Form, Cmpt, nCmpt> stabilise
 366 (
 367     const VectorSpace<Form, Cmpt, nCmpt>& vs,
 368     const Cmpt& small
 369 )
 370 {
 371     Form v;
 372     VectorSpaceOps<nCmpt,0>::opVS(v, vs, small, stabiliseOp<Cmpt>());
 373     return v;
 374 }
 375
 376
 377 template<class Form, class Cmpt, int nCmpt>
 378 inline Cmpt cmptMax
 379 (
 380     const VectorSpace<Form, Cmpt, nCmpt>& vs
 381 )
 382 {
 383     Cmpt cMax = vs.v_[0];
 384     VectorSpaceOps<nCmpt,1>::SeqOp(cMax, vs, maxEqOp<Cmpt>());
 385     return cMax;
 386 }
 387
 388
 389 template<class Form, class Cmpt, int nCmpt>
 390 inline Cmpt cmptMin
 391 (
 392     const VectorSpace<Form, Cmpt, nCmpt>& vs
 393 )
 394 {
 395     Cmpt cMin = vs.v_[0];
 396     VectorSpaceOps<nCmpt,1>::SeqOp(cMin, vs, minEqOp<Cmpt>());
 397     return cMin;
 398 }
 399
 400
 401 template<class Form, class Cmpt, int nCmpt>
 402 inline Cmpt cmptSum
 403 (
 404     const VectorSpace<Form, Cmpt, nCmpt>& vs
 405 )
 406 {
 407     Cmpt sum = vs.v_[0];
 408     VectorSpaceOps<nCmpt,1>::SeqOp(sum, vs, plusEqOp<Cmpt>());
 409     return sum;
 410 }
 411
 412
 413 template<class Form, class Cmpt, int nCmpt>
 414 inline Cmpt cmptAv
 415 (
 416     const VectorSpace<Form, Cmpt, nCmpt>& vs
 417 )
 418 {
 419     return cmptSum(vs)/nCmpt;
 420 }
 421
 422
 423 template<class Form, class Cmpt, int nCmpt>
 424 inline Form cmptMag
 425 (
 426     const VectorSpace<Form, Cmpt, nCmpt>& vs
 427 )
 428 {
 429     Form v;
 430     VectorSpaceOps<nCmpt,0>::eqOp(v, vs, eqMagOp<Cmpt>());
 431     return v;
 432 }
 433
 434
 435 template<class Form, class Cmpt, int nCmpt>
 436 inline Form cmptStabilise
 437 (
 438     const VectorSpace<Form, Cmpt, nCmpt>& vs,
 439     const scalar small,
 440     const scalar value
 441 )
 442 {
 443     Form v = vs;
 444     for (int i=0; i<nCmpt; i++)
 445     {
 446         if (mag(vs.v_[i]) < small)
 447         {
 448             v.v_[i] = sign(vs.v_[i])*value;
 449         }
 450     }
 451     return v;
 452 }
 453
 454
 455 template<class Form, class Cmpt, int nCmpt>
 456 inline Form max
 457 (
 458     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 459     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 460 )
 461 {
 462     Form v;
 463     VectorSpaceOps<nCmpt,0>::op(v, vs1, vs2, maxOp<Cmpt>());
 464     return v;
 465 }
 466
 467
 468 template<class Form, class Cmpt, int nCmpt>
 469 inline Form min
 470 (
 471     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 472     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 473 )
 474 {
 475     Form v;
 476     VectorSpaceOps<nCmpt,0>::op(v, vs1, vs2, minOp<Cmpt>());
 477     return v;
 478 }
 479
 480
 481 template<class Form, class Cmpt, int nCmpt>
 482 inline Form minMod
 483 (
 484     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 485     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 486 )
 487 {
 488     Form v;
 489     VectorSpaceOps<nCmpt,0>::op(v, vs1, vs2, minModOp<Cmpt>());
 490     return v;
 491 }
 492
 493
 494 template<class Type>
 495 inline Type dot(const scalar s, const Type& t)
 496 {
 497     return s * t;
 498 }
 499
 500
 501 template<class Type>
 502 inline Type dot(const Type& t, const scalar s)
 503 {
 504     return t * s;
 505 }
 506
 507
 508 template
 509 <
 510     class Form1, class Cmpt1, int nCmpt1,
 511     class Form2, class Cmpt2, int nCmpt2
 512 >
 513 inline typename innerProduct<Form1, Form2>::type dot
 514 (
 515     const VectorSpace<Form1, Cmpt1, nCmpt1>& t1,
 516     const VectorSpace<Form2, Cmpt2, nCmpt2>& t2
 517 )
 518 {
 519     return static_cast<const Form1&>(t1) & static_cast<const Form2&>(t2);
 520 }
 521
 522
 523 // * * * * * * * * * * * * * * * Global Operators  * * * * * * * * * * * * * //
 524
 525 template<class Form, class Cmpt, int nCmpt>
 526 inline Form operator-
 527 (
 528     const VectorSpace<Form, Cmpt, nCmpt>& vs
 529 )
 530 {
 531     Form v;
 532     VectorSpaceOps<nCmpt,0>::eqOp(v, vs, eqMinusOp<Cmpt>());
 533     return v;
 534 }
 535
 536
 537 template<class Form, class Cmpt, int nCmpt>
 538 inline Form operator+
 539 (
 540     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 541     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 542 )
 543 {
 544     Form v;
 545     VectorSpaceOps<nCmpt,0>::op(v, vs1, vs2, plusOp<Cmpt>());
 546     return v;
 547 }
 548
 549
 550 template<class Form, class Cmpt, int nCmpt>
 551 inline Form operator-
 552 (
 553     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 554     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 555 )
 556 {
 557     Form v;
 558     VectorSpaceOps<nCmpt,0>::op(v, vs1, vs2, minusOp<Cmpt>());
 559     return v;
 560 }
 561
 562
 563 template<class Form, class Cmpt, int nCmpt>
 564 inline Form operator*
 565 (
 566     scalar s,
 567     const VectorSpace<Form, Cmpt, nCmpt>& vs
 568 )
 569 {
 570     Form v;
 571     VectorSpaceOps<nCmpt,0>::opSV(v, s, vs, multiplyOp3<Cmpt, scalar, Cmpt>());
 572     return v;
 573 }
 574
 575
 576 template<class Form, class Cmpt, int nCmpt>
 577 inline Form operator*
 578 (
 579     const VectorSpace<Form, Cmpt, nCmpt>& vs,
 580     scalar s
 581 )
 582 {
 583     Form v;
 584     VectorSpaceOps<nCmpt,0>::opVS(v, vs, s, multiplyOp3<Cmpt, Cmpt, scalar>());
 585     return v;
 586 }
 587
 588
 589 template<class Form, class Cmpt, int nCmpt>
 590 inline Form operator/
 591 (
 592     const VectorSpace<Form, Cmpt, nCmpt>& vs,
 593     scalar s
 594 )
 595 {
 596     Form v;
 597     VectorSpaceOps<nCmpt,0>::opVS(v, vs, s, divideOp3<Cmpt, Cmpt, scalar>());
 598     return v;
 599 }
 600
 601 /*
 602 template<class Form, class Cmpt, int nCmpt>
 603 inline Form operator/
 604 (
 605     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 606     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 607 )
 608 {
 609     Form v;
 610     VectorSpaceOps<nCmpt,0>::op(v, vs1, vs2, divideOp<Cmpt>());
 611     return v;
 612 }
 613
 614
 615 template<class Form, class Cmpt, int nCmpt>
 616 inline Form operator/
 617 (
 618     scalar s,
 619     const VectorSpace<Form, Cmpt, nCmpt>& vs
 620 )
 621 {
 622     Form v;
 623     VectorSpaceOps<nCmpt,0>::opSV(v, s, vs, divideOp2<scalar, Cmpt>());
 624     return v;
 625 }
 626 */
 627
 628
 629 template<class Form, class Cmpt, int nCmpt>
 630 inline Cmpt operator&&
 631 (
 632     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 633     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 634 )
 635 {
 636     Cmpt ddProd = vs1.v_[0]*vs2.v_[0];
 637     for (int i=1; i<nCmpt; ++i)
 638     {
 639         ddProd += vs1.v_[i]*vs2.v_[i];
 640     }
 641     return ddProd;
 642 }
 643
 644
 645 template<class Form, class Cmpt, int nCmpt>
 646 inline bool operator==
 647 (
 648     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 649     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 650 )
 651 {
 652     bool eq = true;
 653     for (int i=0; i<nCmpt; ++i)
 654     {
 655         if (!(eq &= (equal(vs1.v_[i], vs2.v_[i])))) break;
 656     }
 657     return eq;
 658 }
 659
 660
 661 template<class Form, class Cmpt, int nCmpt>
 662 inline bool operator!=
 663 (
 664     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 665     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 666 )
 667 {
 668     return !(vs1 == vs2);
 669 }
 670
 671
 672 template<class Form, class Cmpt, int nCmpt>
 673 inline bool operator>
 674 (
 675     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 676     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 677 )
 678 {
 679     bool gt = true;
 680     for (int i=0; i<nCmpt; ++i)
 681     {
 682         if (!(gt &= vs1.v_[i] > vs2.v_[i])) break;
 683     }
 684     return gt;
 685 }
 686
 687
 688 template<class Form, class Cmpt, int nCmpt>
 689 inline bool operator<
 690 (
 691     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 692     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 693 )
 694 {
 695     bool lt = true;
 696     for (int i=0; i<nCmpt; ++i)
 697     {
 698         if (!(lt &= vs1.v_[i] < vs2.v_[i])) break;
 699     }
 700     return lt;
 701 }
 702
 703
 704 template<class Form, class Cmpt, int nCmpt>
 705 inline bool operator>=
 706 (
 707     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 708     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 709 )
 710 {
 711     return !(vs1 < vs2);
 712 }
 713
 714
 715 template<class Form, class Cmpt, int nCmpt>
 716 inline bool operator<=
 717 (
 718     const VectorSpace<Form, Cmpt, nCmpt>& vs1,
 719     const VectorSpace<Form, Cmpt, nCmpt>& vs2
 720 )
 721 {
 722     return !(vs1 > vs2);
 723 }
 724
 725
 726 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 727
 728 } // End namespace Foam
 729
 730 // ************************************************************************* //