src/dynamicMesh/fvMeshAdder/fvMeshAdderTemplates.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   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
  13     the Free Software Foundation, either version 3 of the License, or
  14     (at your 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, see <http://www.gnu.org/licenses/>.
  23
  24 \*---------------------------------------------------------------------------*/
  25
  26
  27 #include "volFields.H"
  28 #include "surfaceFields.H"
  29 #include "emptyFvPatchField.H"
  30
  31 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  32
  33 template<class Type>
  34 void Foam::fvMeshAdder::map
  35 (
  36     const Field<Type>& oldFld,
  37     const labelList& oldToNew,
  38     Field<Type>& fld
  39 )
  40 {
  41     forAll(oldFld, cellI)
  42     {
  43         label newCellI = oldToNew[cellI];
  44
  45         if (newCellI >= 0 && newCellI < fld.size())
  46         {
  47             fld[newCellI] = oldFld[cellI];
  48         }
  49     }
  50 }
  51
  52
  53 template<class Type>
  54 void Foam::fvMeshAdder::MapVolField
  55 (
  56     const mapAddedPolyMesh& meshMap,
  57
  58     GeometricField<Type, fvPatchField, volMesh>& fld,
  59     const GeometricField<Type, fvPatchField, volMesh>& fldToAdd
  60 )
  61 {
  62     const fvMesh& mesh = fld.mesh();
  63
  64     // Internal field
  65     // ~~~~~~~~~~~~~~
  66
  67     {
  68         // Store old internal field
  69         Field<Type> oldInternalField(fld.internalField());
  70
  71         // Modify internal field
  72         Field<Type>& intFld = fld.internalField();
  73
  74         intFld.setSize(mesh.nCells());
  75
  76         map(oldInternalField, meshMap.oldCellMap(), intFld);
  77         map(fldToAdd.internalField(), meshMap.addedCellMap(), intFld);
  78     }
  79
  80
  81     // Patch fields from old mesh
  82     // ~~~~~~~~~~~~~~~~~~~~~~~~~~
  83
  84     {
  85         const labelList& oldPatchMap = meshMap.oldPatchMap();
  86         const labelList& oldPatchStarts = meshMap.oldPatchStarts();
  87         const labelList& oldPatchSizes = meshMap.oldPatchSizes();
  88
  89         // Reorder old patches in order of new ones. Put removed patches at end.
  90
  91         label unusedPatchI = 0;
  92
  93         forAll(oldPatchMap, patchI)
  94         {
  95             label newPatchI = oldPatchMap[patchI];
  96
  97             if (newPatchI != -1)
  98             {
  99                 unusedPatchI++;
 100             }
 101         }
 102
 103         label nUsedPatches = unusedPatchI;
 104
 105         // Reorder list for patchFields
 106         labelList oldToNew(oldPatchMap.size());
 107
 108         forAll(oldPatchMap, patchI)
 109         {
 110             label newPatchI = oldPatchMap[patchI];
 111
 112             if (newPatchI != -1)
 113             {
 114                 oldToNew[patchI] = newPatchI;
 115             }
 116             else
 117             {
 118                 oldToNew[patchI] = unusedPatchI++;
 119             }
 120         }
 121
 122
 123         // Sort deleted ones last so is now in newPatch ordering
 124         fld.boundaryField().reorder(oldToNew);
 125         // Extend to covers all patches
 126         fld.boundaryField().setSize(mesh.boundaryMesh().size());
 127         // Delete unused patches
 128         for
 129         (
 130             label newPatchI = nUsedPatches;
 131             newPatchI < fld.boundaryField().size();
 132             newPatchI++
 133         )
 134         {
 135             fld.boundaryField().set(newPatchI, NULL);
 136         }
 137
 138
 139         // Map old values
 140         // ~~~~~~~~~~~~~~
 141
 142         forAll(oldPatchMap, patchI)
 143         {
 144             label newPatchI = oldPatchMap[patchI];
 145
 146             if (newPatchI != -1)
 147             {
 148                 labelList newToOld
 149                 (
 150                     calcPatchMap
 151                     (
 152                         oldPatchStarts[patchI],
 153                         oldPatchSizes[patchI],
 154                         meshMap.oldFaceMap(),
 155                         mesh.boundaryMesh()[newPatchI],
 156                         0                   // unmapped value
 157                     )
 158                 );
 159
 160                 directFvPatchFieldMapper patchMapper(newToOld);
 161
 162
 163                 // Create new patchField with same type as existing one.
 164                 // Note:
 165                 // - boundaryField already in new order so access with newPatchI
 166                 // - fld.boundaryField()[newPatchI] both used for type and old
 167                 //   value
 168                 // - hope that field mapping allows aliasing since old and new
 169                 //   are same memory!
 170                 fld.boundaryField().set
 171                 (
 172                     newPatchI,
 173                     fvPatchField<Type>::New
 174                     (
 175                         fld.boundaryField()[newPatchI], // old field
 176                         mesh.boundary()[newPatchI],     // new fvPatch
 177                         fld.dimensionedInternalField(), // new internal field
 178                         patchMapper                     // mapper (new to old)
 179                     )
 180                 );
 181             }
 182         }
 183     }
 184
 185
 186
 187     // Patch fields from added mesh
 188     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 189
 190     {
 191         const labelList& addedPatchMap = meshMap.addedPatchMap();
 192
 193         // Add addedMesh patches
 194         forAll(addedPatchMap, patchI)
 195         {
 196             label newPatchI = addedPatchMap[patchI];
 197
 198             if (newPatchI != -1)
 199             {
 200                 const polyPatch& newPatch = mesh.boundaryMesh()[newPatchI];
 201                 const polyPatch& oldPatch =
 202                     fldToAdd.mesh().boundaryMesh()[patchI];
 203
 204                 if (!fld.boundaryField()(newPatchI))
 205                 {
 206                     // First occurrence of newPatchI. Map from existing
 207                     // patchField
 208
 209                     // From new patch faces to patch faces on added mesh.
 210                     labelList newToAdded
 211                     (
 212                         calcPatchMap
 213                         (
 214                             oldPatch.start(),
 215                             oldPatch.size(),
 216                             meshMap.addedFaceMap(),
 217                             newPatch,
 218                             0                       // unmapped values
 219                         )
 220                     );
 221
 222                     directFvPatchFieldMapper patchMapper(newToAdded);
 223
 224                     fld.boundaryField().set
 225                     (
 226                         newPatchI,
 227                         fvPatchField<Type>::New
 228                         (
 229                             fldToAdd.boundaryField()[patchI], // added field
 230                             mesh.boundary()[newPatchI],       // new fvPatch
 231                             fld.dimensionedInternalField(),   // new int. field
 232                             patchMapper                       // mapper
 233                         )
 234                     );
 235                 }
 236                 else
 237                 {
 238                     // PatchField will have correct size already. Just slot in
 239                     // my elements.
 240
 241                     // From new patch faces to patch faces on added mesh. This
 242                     // time keep unmapped elements -1.
 243                     labelList newToAdded
 244                     (
 245                         calcPatchMap
 246                         (
 247                             oldPatch.start(),
 248                             oldPatch.size(),
 249                             meshMap.addedFaceMap(),
 250                             newPatch,
 251                             -1                      // unmapped values
 252                         )
 253                     );
 254
 255                     const fvPatchField<Type>& addedFld =
 256                         fldToAdd.boundaryField()[patchI];
 257
 258                     fvPatchField<Type>& newFld = fld.boundaryField()[newPatchI];
 259
 260                     forAll(newFld, i)
 261                     {
 262                         label oldFaceI = newToAdded[i];
 263
 264                         if (oldFaceI >= 0 && oldFaceI < addedFld.size())
 265                         {
 266                             newFld[i] = addedFld[oldFaceI];
 267                         }
 268                     }
 269                 }
 270             }
 271         }
 272     }
 273 }
 274
 275
 276 template<class Type>
 277 void Foam::fvMeshAdder::MapVolFields
 278 (
 279     const mapAddedPolyMesh& meshMap,
 280     const fvMesh& mesh,
 281     const fvMesh& meshToAdd
 282 )
 283 {
 284     HashTable<const GeometricField<Type, fvPatchField, volMesh>*> fields
 285     (
 286         mesh.objectRegistry::lookupClass
 287         <GeometricField<Type, fvPatchField, volMesh> >
 288         ()
 289     );
 290
 291     HashTable<const GeometricField<Type, fvPatchField, volMesh>*> fieldsToAdd
 292     (
 293         meshToAdd.objectRegistry::lookupClass
 294         <GeometricField<Type, fvPatchField, volMesh> >
 295         ()
 296     );
 297
 298     // It is necessary to enforce that all old-time fields are stored
 299     // before the mapping is performed.  Otherwise, if the
 300     // old-time-level field is mapped before the field itself, sizes
 301     // will not match.
 302
 303     for
 304     (
 305         typename HashTable<const GeometricField<Type, fvPatchField, volMesh>*>::
 306             iterator fieldIter = fields.begin();
 307         fieldIter != fields.end();
 308         ++fieldIter
 309     )
 310     {
 311         const_cast<GeometricField<Type, fvPatchField, volMesh>*>(fieldIter())
 312             ->storeOldTimes();
 313     }
 314
 315
 316     for
 317     (
 318         typename HashTable<const GeometricField<Type, fvPatchField, volMesh>*>::
 319             iterator fieldIter = fields.begin();
 320         fieldIter != fields.end();
 321         ++fieldIter
 322     )
 323     {
 324         GeometricField<Type, fvPatchField, volMesh>& fld =
 325             const_cast<GeometricField<Type, fvPatchField, volMesh>&>
 326             (
 327                 *fieldIter()
 328             );
 329
 330         if (fieldsToAdd.found(fld.name()))
 331         {
 332             const GeometricField<Type, fvPatchField, volMesh>& fldToAdd =
 333                 *fieldsToAdd[fld.name()];
 334
 335             MapVolField<Type>(meshMap, fld, fldToAdd);
 336         }
 337         else
 338         {
 339             WarningIn("fvMeshAdder::MapVolFields(..)")
 340                 << "Not mapping field " << fld.name()
 341                 << " since not present on mesh to add"
 342                 << endl;
 343         }
 344     }
 345 }
 346
 347
 348 template<class Type>
 349 void Foam::fvMeshAdder::MapSurfaceField
 350 (
 351     const mapAddedPolyMesh& meshMap,
 352
 353     GeometricField<Type, fvsPatchField, surfaceMesh>& fld,
 354     const GeometricField<Type, fvsPatchField, surfaceMesh>& fldToAdd
 355 )
 356 {
 357     const fvMesh& mesh = fld.mesh();
 358     const labelList& oldPatchStarts = meshMap.oldPatchStarts();
 359
 360     // Internal field
 361     // ~~~~~~~~~~~~~~
 362
 363     // Store old internal field
 364     {
 365         Field<Type> oldField(fld);
 366
 367         // Modify internal field
 368         Field<Type>& intFld = fld.internalField();
 369
 370         intFld.setSize(mesh.nInternalFaces());
 371
 372         map(oldField, meshMap.oldFaceMap(), intFld);
 373         map(fldToAdd, meshMap.addedFaceMap(), intFld);
 374
 375         // Faces that were boundary faces but are not anymore.
 376         // Use owner value (so lowest numbered cell, i.e. from 'old' not 'added'
 377         // mesh)
 378         forAll(fld.boundaryField(), patchI)
 379         {
 380             const fvsPatchField<Type>& pf = fld.boundaryField()[patchI];
 381
 382             label start = oldPatchStarts[patchI];
 383
 384             forAll(pf, i)
 385             {
 386                 label newFaceI = meshMap.oldFaceMap()[start + i];
 387
 388                 if (newFaceI >= 0 && newFaceI < mesh.nInternalFaces())
 389                 {
 390                     intFld[newFaceI] = pf[i];
 391                 }
 392             }
 393         }
 394     }
 395
 396
 397     // Patch fields from old mesh
 398     // ~~~~~~~~~~~~~~~~~~~~~~~~~~
 399
 400     {
 401         const labelList& oldPatchMap = meshMap.oldPatchMap();
 402         const labelList& oldPatchSizes = meshMap.oldPatchSizes();
 403
 404         // Reorder old patches in order of new ones. Put removed patches at end.
 405
 406         label unusedPatchI = 0;
 407
 408         forAll(oldPatchMap, patchI)
 409         {
 410             label newPatchI = oldPatchMap[patchI];
 411
 412             if (newPatchI != -1)
 413             {
 414                 unusedPatchI++;
 415             }
 416         }
 417
 418         label nUsedPatches = unusedPatchI;
 419
 420         // Reorder list for patchFields
 421         labelList oldToNew(oldPatchMap.size());
 422
 423         forAll(oldPatchMap, patchI)
 424         {
 425             label newPatchI = oldPatchMap[patchI];
 426
 427             if (newPatchI != -1)
 428             {
 429                 oldToNew[patchI] = newPatchI;
 430             }
 431             else
 432             {
 433                 oldToNew[patchI] = unusedPatchI++;
 434             }
 435         }
 436
 437
 438         // Sort deleted ones last so is now in newPatch ordering
 439         fld.boundaryField().reorder(oldToNew);
 440         // Extend to covers all patches
 441         fld.boundaryField().setSize(mesh.boundaryMesh().size());
 442         // Delete unused patches
 443         for
 444         (
 445             label newPatchI = nUsedPatches;
 446             newPatchI < fld.boundaryField().size();
 447             newPatchI++
 448         )
 449         {
 450             fld.boundaryField().set(newPatchI, NULL);
 451         }
 452
 453
 454         // Map old values
 455         // ~~~~~~~~~~~~~~
 456
 457         forAll(oldPatchMap, patchI)
 458         {
 459             label newPatchI = oldPatchMap[patchI];
 460
 461             if (newPatchI != -1)
 462             {
 463                 labelList newToOld
 464                 (
 465                     calcPatchMap
 466                     (
 467                         oldPatchStarts[patchI],
 468                         oldPatchSizes[patchI],
 469                         meshMap.oldFaceMap(),
 470                         mesh.boundaryMesh()[newPatchI],
 471                         0                   // unmapped value
 472                     )
 473                 );
 474
 475                 directFvPatchFieldMapper patchMapper(newToOld);
 476
 477
 478                 // Create new patchField with same type as existing one.
 479                 // Note:
 480                 // - boundaryField already in new order so access with newPatchI
 481                 // - fld.boundaryField()[newPatchI] both used for type and old
 482                 //   value
 483                 // - hope that field mapping allows aliasing since old and new
 484                 //   are same memory!
 485                 fld.boundaryField().set
 486                 (
 487                     newPatchI,
 488                     fvsPatchField<Type>::New
 489                     (
 490                         fld.boundaryField()[newPatchI], // old field
 491                         mesh.boundary()[newPatchI],     // new fvPatch
 492                         fld.dimensionedInternalField(), // new internal field
 493                         patchMapper                     // mapper (new to old)
 494                     )
 495                 );
 496             }
 497         }
 498     }
 499
 500
 501
 502     // Patch fields from added mesh
 503     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 504
 505     {
 506         const labelList& addedPatchMap = meshMap.addedPatchMap();
 507
 508         // Add addedMesh patches
 509         forAll(addedPatchMap, patchI)
 510         {
 511             label newPatchI = addedPatchMap[patchI];
 512
 513             if (newPatchI != -1)
 514             {
 515                 const polyPatch& newPatch = mesh.boundaryMesh()[newPatchI];
 516                 const polyPatch& oldPatch =
 517                     fldToAdd.mesh().boundaryMesh()[patchI];
 518
 519                 if (!fld.boundaryField()(newPatchI))
 520                 {
 521                     // First occurrence of newPatchI. Map from existing
 522                     // patchField
 523
 524                     // From new patch faces to patch faces on added mesh.
 525                     labelList newToAdded
 526                     (
 527                         calcPatchMap
 528                         (
 529                             oldPatch.start(),
 530                             oldPatch.size(),
 531                             meshMap.addedFaceMap(),
 532                             newPatch,
 533                             0                       // unmapped values
 534                         )
 535                     );
 536
 537                     directFvPatchFieldMapper patchMapper(newToAdded);
 538
 539                     fld.boundaryField().set
 540                     (
 541                         newPatchI,
 542                         fvsPatchField<Type>::New
 543                         (
 544                             fldToAdd.boundaryField()[patchI],// added field
 545                             mesh.boundary()[newPatchI],      // new fvPatch
 546                             fld.dimensionedInternalField(),  // new int. field
 547                             patchMapper                      // mapper
 548                         )
 549                     );
 550                 }
 551                 else
 552                 {
 553                     // PatchField will have correct size already. Just slot in
 554                     // my elements.
 555
 556                     // From new patch faces to patch faces on added mesh. This
 557                     // time keep unmapped elements -1.
 558                     labelList newToAdded
 559                     (
 560                         calcPatchMap
 561                         (
 562                             oldPatch.start(),
 563                             oldPatch.size(),
 564                             meshMap.addedFaceMap(),
 565                             newPatch,
 566                             -1                      // unmapped values
 567                         )
 568                     );
 569
 570                     const fvsPatchField<Type>& addedFld =
 571                         fldToAdd.boundaryField()[patchI];
 572
 573                     fvsPatchField<Type>& newFld =
 574                         fld.boundaryField()[newPatchI];
 575
 576                     forAll(newFld, i)
 577                     {
 578                         label oldFaceI = newToAdded[i];
 579
 580                         if (oldFaceI >= 0 && oldFaceI < addedFld.size())
 581                         {
 582                             newFld[i] = addedFld[oldFaceI];
 583                         }
 584                     }
 585                 }
 586             }
 587         }
 588     }
 589 }
 590
 591
 592 template<class Type>
 593 void Foam::fvMeshAdder::MapSurfaceFields
 594 (
 595     const mapAddedPolyMesh& meshMap,
 596     const fvMesh& mesh,
 597     const fvMesh& meshToAdd
 598 )
 599 {
 600     typedef GeometricField<Type, fvsPatchField, surfaceMesh> fldType;
 601
 602     HashTable<const fldType*> fields
 603     (
 604         mesh.objectRegistry::lookupClass<fldType>()
 605     );
 606
 607     HashTable<const fldType*> fieldsToAdd
 608     (
 609         meshToAdd.objectRegistry::lookupClass<fldType>()
 610     );
 611
 612     // It is necessary to enforce that all old-time fields are stored
 613     // before the mapping is performed.  Otherwise, if the
 614     // old-time-level field is mapped before the field itself, sizes
 615     // will not match.
 616
 617     for
 618     (
 619         typename HashTable<const fldType*>::
 620             iterator fieldIter = fields.begin();
 621         fieldIter != fields.end();
 622         ++fieldIter
 623     )
 624     {
 625         const_cast<fldType*>(fieldIter())
 626             ->storeOldTimes();
 627     }
 628
 629
 630     for
 631     (
 632         typename HashTable<const fldType*>::
 633             iterator fieldIter = fields.begin();
 634         fieldIter != fields.end();
 635         ++fieldIter
 636     )
 637     {
 638         fldType& fld = const_cast<fldType&>(*fieldIter());
 639
 640         if (fieldsToAdd.found(fld.name()))
 641         {
 642             const fldType& fldToAdd = *fieldsToAdd[fld.name()];
 643
 644             MapSurfaceField<Type>(meshMap, fld, fldToAdd);
 645         }
 646         else
 647         {
 648             WarningIn("fvMeshAdder::MapSurfaceFields(..)")
 649                 << "Not mapping field " << fld.name()
 650                 << " since not present on mesh to add"
 651                 << endl;
 652         }
 653     }
 654 }
 655
 656
 657 // ************************************************************************* //