src/foam/interpolations/GGIInterpolation/GGIInterpolate.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | foam-extend: Open Source CFD
   4    \\    /   O peration     | Version:     3.2
   5     \\  /    A nd           | Web:         http://www.foam-extend.org
   6      \\/     M anipulation  | For copyright notice see file Copyright
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of foam-extend.
  10
  11     foam-extend 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 3 of the License, or (at your
  14     option) any later version.
  15
  16     foam-extend is distributed in the hope that it will be useful, but
  17     WITHOUT ANY WARRANTY; without even the implied warranty of
  18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19     General Public License for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Description
  25     GGI interpolation functions
  26
  27 Author
  28     Hrvoje Jasak, Wikki Ltd.  All rights reserved
  29
  30 \*---------------------------------------------------------------------------*/
  31
  32 namespace Foam
  33 {
  34
  35 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  36
  37 template<class MasterPatch, class SlavePatch>
  38 template<class Type>
  39 void GGIInterpolation<MasterPatch, SlavePatch>::interpolate
  40 (
  41     const Field<Type>& ff,
  42     Field<Type>& result,
  43     const labelListList& addr,
  44     const scalarListList& weights
  45 )
  46 {
  47     forAll (result, faceI)
  48     {
  49         const labelList& curAddr = addr[faceI];
  50         const scalarList& curWeights = weights[faceI];
  51
  52         result[faceI] = pTraits<Type>::zero;
  53
  54         forAll (curAddr, i)
  55         {
  56             result[faceI] += ff[curAddr[i]]*curWeights[i];
  57         }
  58     }
  59 }
  60
  61
  62 template<class MasterPatch, class SlavePatch>
  63 template<class Type>
  64 void GGIInterpolation<MasterPatch, SlavePatch>::maskedInterpolate
  65 (
  66     const Field<Type>& ff,
  67     Field<Type>& result,
  68     const labelList& mask,
  69     const labelListList& addr,
  70     const scalarListList& weights
  71 )
  72 {
  73     forAll (mask, maskI)
  74     {
  75         // Pick the masked face
  76         const label faceI = mask[maskI];
  77
  78         const labelList& curAddr = addr[faceI];
  79
  80         const scalarList& curWeights = weights[faceI];
  81
  82         // Clear condensed list entry: masked faces only
  83         result[maskI] = pTraits<Type>::zero;
  84
  85         forAll (curAddr, i)
  86         {
  87             // Put the result into condensed list: masked faces only
  88             result[maskI] += ff[curAddr[i]]*curWeights[i];
  89         }
  90     }
  91 }
  92
  93
  94 template<class MasterPatch, class SlavePatch>
  95 template<class Type>
  96 void GGIInterpolation<MasterPatch, SlavePatch>::bridge
  97 (
  98     const Field<Type>& bridgeField,
  99     Field<Type>& ff,
 100     const labelList& addr
 101 )
 102 {
 103     forAll (addr, faceI)
 104     {
 105         ff[addr[faceI]] = bridgeField[addr[faceI]];
 106     }
 107 }
 108
 109
 110 template<class MasterPatch, class SlavePatch>
 111 template<class Type>
 112 void GGIInterpolation<MasterPatch, SlavePatch>::maskedBridge
 113 (
 114     const Field<Type>& bridgeField,
 115     Field<Type>& ff,
 116     const labelList& mask,
 117     const labelList& uncoveredFaces
 118 )
 119 {
 120     // Note: tricky algorithm
 121     // In order for a face to be bridged it needs to be both in the
 122     // mask and in selection of faces that are bridged (addr).
 123     // This implies an n-squared search, but we can use the fact that
 124     // both lists are ordered.
 125
 126     label maskAddrI = 0;
 127
 128     forAll (uncoveredFaces, uncoI)
 129     {
 130         // Pick the uncovered face
 131         const label faceI = uncoveredFaces[uncoI];
 132
 133         // Search through the mask
 134         for (; maskAddrI < mask.size(); maskAddrI++)
 135         {
 136             if (faceI == mask[maskAddrI])
 137             {
 138                 // Found masked bridged face
 139                 // Put the result into condensed list: masked faces only
 140                 ff[maskAddrI] = bridgeField[maskAddrI];
 141
 142                 break;
 143             }
 144             else if (mask[maskAddrI] > faceI)
 145             {
 146                 // Gone beyond my index: my face is not present in the mask
 147                 // Go one back and check for next uncovered face
 148                 if (maskAddrI > 0)
 149                 {
 150                     maskAddrI--;
 151                 }
 152
 153                 break;
 154             }
 155         }
 156     }
 157 }
 158
 159
 160 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 161
 162 template<class MasterPatch, class SlavePatch>
 163 template<class Type>
 164 tmp<Field<Type> >
 165 GGIInterpolation<MasterPatch, SlavePatch>::masterToSlave
 166 (
 167     const Field<Type>& ff
 168 ) const
 169 {
 170     if (ff.size() != masterPatch_.size())
 171     {
 172         FatalErrorIn
 173         (
 174             "GGIInterpolation::masterToSlave(const Field<Type> ff)"
 175         )   << "given field does not correspond to patch. Patch size: "
 176             << masterPatch_.size() << " field size: " << ff.size()
 177             << abort(FatalError);
 178     }
 179
 180     tmp<Field<Type> > tresult
 181     (
 182         new Field<Type>
 183         (
 184             slavePatch_.size(),
 185             pTraits<Type>::zero
 186         )
 187     );
 188
 189     // Do interpolation
 190     Field<Type>& result = tresult();
 191
 192     if (this->doTransform() && pTraits<Type>::rank > 0)
 193     {
 194         // Transform master data to slave
 195         Field<Type> transformFF;
 196
 197         if (reverseT_.size() == 1)
 198         {
 199             // Constant transform
 200             transformFF = transform(reverseT_[0], ff);
 201         }
 202         else
 203         {
 204             // Full patch transform
 205             transformFF = transform(reverseT_, ff);
 206         }
 207
 208         GGIInterpolation<MasterPatch, SlavePatch>::interpolate
 209         (
 210             transformFF,
 211             result,
 212             this->slaveAddr(),
 213             this->slaveWeights()
 214         );
 215     }
 216     else
 217     {
 218         GGIInterpolation<MasterPatch, SlavePatch>::interpolate
 219         (
 220             ff,
 221             result,
 222             this->slaveAddr(),
 223             this->slaveWeights()
 224         );
 225     }
 226
 227     return tresult;
 228 }
 229
 230
 231 template<class MasterPatch, class SlavePatch>
 232 template<class Type>
 233 tmp<Field<Type> >
 234 GGIInterpolation<MasterPatch, SlavePatch>::masterToSlave
 235 (
 236     const tmp<Field<Type> >& tff
 237 ) const
 238 {
 239     tmp<Field<Type> > tint = masterToSlave(tff());
 240     tff.clear();
 241     return tint;
 242 }
 243
 244
 245 template<class MasterPatch, class SlavePatch>
 246 template<class Type>
 247 void GGIInterpolation<MasterPatch, SlavePatch>::maskedMasterToSlave
 248 (
 249     const Field<Type>& ff,
 250     Field<Type>& result,
 251     const labelList& mask
 252 ) const
 253 {
 254     if (ff.size() != masterPatch_.size())
 255     {
 256         FatalErrorIn
 257         (
 258             "void GGIInterpolation::maskedMasterToSlave\n"
 259             "(\n"
 260             "    const Field<Type>& ff,\n"
 261             "    Field<Type>& result,\n"
 262             "    const labelList& mask\n"
 263             ") const"
 264         )   << "given field does not correspond to patch. Patch size: "
 265             << masterPatch_.size() << " field size: " << ff.size()
 266             << abort(FatalError);
 267     }
 268
 269     if (result.size() != mask.size())
 270     {
 271         FatalErrorIn
 272         (
 273             "void GGIInterpolation::maskedMasterToSlave\n"
 274             "(\n"
 275             "    const Field<Type>& ff,\n"
 276             "    Field<Type>& result,\n"
 277             "    const labelList& mask\n"
 278             ") const"
 279         )   << "result field does not correspond to mask. Field size: "
 280             << result.size() << " mask size: " << mask.size()
 281             << abort(FatalError);
 282     }
 283
 284     if (this->doTransform() && pTraits<Type>::rank > 0)
 285     {
 286         // Transform master data to slave
 287         Field<Type> transformFF;
 288
 289         if (reverseT_.size() == 1)
 290         {
 291             // Constant transform
 292             transformFF = transform(reverseT_[0], ff);
 293         }
 294         else
 295         {
 296             // Full patch transform
 297             transformFF = transform(reverseT_, ff);
 298         }
 299
 300         GGIInterpolation<MasterPatch, SlavePatch>::maskedInterpolate
 301         (
 302             transformFF,
 303             result,
 304             mask,
 305             this->slaveAddr(),
 306             this->slaveWeights()
 307         );
 308     }
 309     else
 310     {
 311         GGIInterpolation<MasterPatch, SlavePatch>::maskedInterpolate
 312         (
 313             ff,
 314             result,
 315             mask,
 316             this->slaveAddr(),
 317             this->slaveWeights()
 318         );
 319     }
 320 }
 321
 322
 323 template<class MasterPatch, class SlavePatch>
 324 template<class Type>
 325 tmp<Field<Type> >
 326 GGIInterpolation<MasterPatch, SlavePatch>::slaveToMaster
 327 (
 328     const Field<Type>& ff
 329 ) const
 330 {
 331     if (ff.size() != slavePatch_.size())
 332     {
 333         FatalErrorIn
 334         (
 335             "GGIInterpolation::slaveToMaster"
 336             "(const Field<Type> ff)"
 337         )   << "given field does not correspond to patch. Patch size: "
 338             << slavePatch_.size() << " field size: " << ff.size()
 339             << abort(FatalError);
 340     }
 341
 342     tmp<Field<Type> > tresult
 343     (
 344         new Field<Type>
 345         (
 346             masterPatch_.size(),
 347             pTraits<Type>::zero
 348         )
 349     );
 350
 351     // Do interpolation
 352     Field<Type>& result = tresult();
 353
 354     if (this->doTransform() && pTraits<Type>::rank > 0)
 355     {
 356         // Transform slave data to master
 357         Field<Type> transformFF;
 358         if (forwardT_.size() == 1)
 359         {
 360             // Constant transform
 361             transformFF = transform(forwardT_[0], ff);
 362         }
 363         else
 364         {
 365             // Full patch transform
 366             transformFF = transform(forwardT_, ff);
 367         }
 368
 369         GGIInterpolation<MasterPatch, SlavePatch>::interpolate
 370         (
 371             transformFF,
 372             result,
 373             this->masterAddr(),
 374             this->masterWeights()
 375         );
 376     }
 377     else
 378     {
 379         GGIInterpolation<MasterPatch, SlavePatch>::interpolate
 380         (
 381             ff,
 382             result,
 383             this->masterAddr(),
 384             this->masterWeights()
 385         );
 386     }
 387
 388     return tresult;
 389 }
 390
 391
 392 template<class MasterPatch, class SlavePatch>
 393 template<class Type>
 394 tmp<Field<Type> >
 395 GGIInterpolation<MasterPatch, SlavePatch>::slaveToMaster
 396 (
 397     const tmp<Field<Type> >& tff
 398 ) const
 399 {
 400     tmp<Field<Type> > tint = slaveToMaster(tff());
 401     tff.clear();
 402     return tint;
 403 }
 404
 405
 406 template<class MasterPatch, class SlavePatch>
 407 template<class Type>
 408 void GGIInterpolation<MasterPatch, SlavePatch>::maskedSlaveToMaster
 409 (
 410     const Field<Type>& ff,
 411     Field<Type>& result,
 412     const labelList& mask
 413 ) const
 414 {
 415     if (ff.size() != slavePatch_.size())
 416     {
 417         FatalErrorIn
 418         (
 419             "void GGIInterpolation::maskedSlaveToMaster"
 420             "(\n"
 421             "    const Field<Type>& ff,\n"
 422             "    Field<Type>& result,\n"
 423             "    const labelList& mask\n"
 424             ") const"
 425         )   << "given field does not correspond to patch. Patch size: "
 426             << slavePatch_.size() << " field size: " << ff.size()
 427             << abort(FatalError);
 428     }
 429
 430     if (result.size() != mask.size())
 431     {
 432         FatalErrorIn
 433         (
 434             "void GGIInterpolation::maskedSlaveToMaster\n"
 435             "(\n"
 436             "    const Field<Type>& ff,\n"
 437             "    Field<Type>& result,\n"
 438             "    const labelList& mask\n"
 439             ") const"
 440         )   << "result field does not correspond to mask. Field size: "
 441             << result.size() << " mask size: " << mask.size()
 442             << abort(FatalError);
 443     }
 444
 445     if (this->doTransform() && pTraits<Type>::rank > 0)
 446     {
 447         // Transform slave data to master
 448         Field<Type> transformFF;
 449         if (forwardT_.size() == 1)
 450         {
 451             // Constant transform
 452             transformFF = transform(forwardT_[0], ff);
 453         }
 454         else
 455         {
 456             // Full patch transform
 457             transformFF = transform(forwardT_, ff);
 458         }
 459
 460         GGIInterpolation<MasterPatch, SlavePatch>::maskedInterpolate
 461         (
 462             transformFF,
 463             result,
 464             mask,
 465             this->masterAddr(),
 466             this->masterWeights()
 467         );
 468     }
 469     else
 470     {
 471         GGIInterpolation<MasterPatch, SlavePatch>::maskedInterpolate
 472         (
 473             ff,
 474             result,
 475             mask,
 476             this->masterAddr(),
 477             this->masterWeights()
 478         );
 479     }
 480 }
 481
 482
 483 template<class MasterPatch, class SlavePatch>
 484 template<class Type>
 485 void GGIInterpolation<MasterPatch, SlavePatch>::bridgeMaster
 486 (
 487     const Field<Type>& bridgeField,
 488     Field<Type>& ff
 489 ) const
 490 {
 491     if
 492     (
 493         bridgeField.size() != masterPatch_.size()
 494      || ff.size() != masterPatch_.size())
 495     {
 496         FatalErrorIn
 497         (
 498             "void GGIInterpolation<MasterPatch, SlavePatch>::bridgeMaster\n"
 499             "(\n"
 500             "    const Field<Type>& bridgeField,\n"
 501             "    Field<Type>& ff\n"
 502             ") const"
 503         )   << "given field does not correspond to patch. Patch size: "
 504             << masterPatch_.size()
 505             << " bridge field size: " << bridgeField.size()
 506             << " field size: " << ff.size()
 507             << abort(FatalError);
 508     }
 509
 510     bridge(bridgeField, ff, uncoveredMasterFaces());
 511 }
 512
 513
 514 template<class MasterPatch, class SlavePatch>
 515 template<class Type>
 516 void GGIInterpolation<MasterPatch, SlavePatch>::maskedBridgeMaster
 517 (
 518     const Field<Type>& bridgeField,
 519     Field<Type>& ff,
 520     const labelList& mask
 521 ) const
 522 {
 523     if (ff.size() != mask.size())
 524     {
 525         FatalErrorIn
 526         (
 527             "void GGIInterpolation<MasterPatch, SlavePatch>::"
 528             "maskedBridgeMaster\n"
 529             "(\n"
 530             "    const Field<Type>& bridgeField,\n"
 531             "    Field<Type>& ff,\n"
 532             "    const labelList& mask\n"
 533             ") const"
 534         )   << "given field does not correspond to patch. Patch size: "
 535             << masterPatch_.size()
 536             << " bridge field size: " << bridgeField.size()
 537             << " field size: " << ff.size()
 538             << " mask size: " << mask.size()
 539             << abort(FatalError);
 540     }
 541
 542     maskedBridge(bridgeField, ff, mask, uncoveredMasterFaces());
 543 }
 544
 545
 546 template<class MasterPatch, class SlavePatch>
 547 template<class Type>
 548 void GGIInterpolation<MasterPatch, SlavePatch>::bridgeSlave
 549 (
 550     const Field<Type>& bridgeField,
 551     Field<Type>& ff
 552 ) const
 553 {
 554     if
 555     (
 556         bridgeField.size() != slavePatch_.size()
 557      || ff.size() != slavePatch_.size()
 558     )
 559     {
 560         FatalErrorIn
 561         (
 562             "void GGIInterpolation<MasterPatch, SlavePatch>::"
 563             "bridgeSlave\n"
 564             "(\n"
 565             "    const Field<Type>& bridgeField,\n"
 566             "    Field<Type>& ff"
 567             ") const"
 568         )   << "given field does not correspond to patch. Patch size: "
 569             << slavePatch_.size()
 570             << " bridge field size: " << bridgeField.size()
 571             << " field size: " << ff.size()
 572             << abort(FatalError);
 573     }
 574
 575     bridge(bridgeField, ff, uncoveredSlaveFaces());
 576 }
 577
 578
 579 template<class MasterPatch, class SlavePatch>
 580 template<class Type>
 581 void GGIInterpolation<MasterPatch, SlavePatch>::maskedBridgeSlave
 582 (
 583     const Field<Type>& bridgeField,
 584     Field<Type>& ff,
 585     const labelList& mask
 586 ) const
 587 {
 588     if (ff.size() != mask.size())
 589     {
 590         FatalErrorIn
 591         (
 592             "void GGIInterpolation<MasterPatch, SlavePatch>::"
 593             "maskedBridgeSlave\n"
 594             "(\n"
 595             "    const Field<Type>& bridgeField,\n"
 596             "    Field<Type>& ff\n,"
 597             "    const labelList& mask\n"
 598             ") const"
 599         )   << "given field does not correspond to patch. Patch size: "
 600             << slavePatch_.size()
 601             << " bridge field size: " << bridgeField.size()
 602             << " field size: " << ff.size()
 603             << " mask size: " << mask.size()
 604             << abort(FatalError);
 605     }
 606
 607     maskedBridge(bridgeField, ff, mask, uncoveredSlaveFaces());
 608 }
 609
 610
 611 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 612
 613 } // End namespace Foam
 614
 615 // ************************************************************************* //