applications/utilities/mesh/manipulation/createBaffles/createBaffles.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     Makes internal faces into boundary faces. Does not duplicate points, unlike
  26     mergeOrSplitBaffles.
  27
  28     Note: if any coupled patch face is selected for baffling the opposite
  29     member has to be selected for baffling as well. Note that this
  30     is the same as repatching. This was added only for convenience so
  31     you don't have to filter coupled boundary out of your set.
  32
  33 \*---------------------------------------------------------------------------*/
  34
  35 #include "syncTools.H"
  36 #include "argList.H"
  37 #include "foamTime.H"
  38 #include "faceSet.H"
  39 #include "directTopoChange.H"
  40 #include "polyModifyFace.H"
  41 #include "polyAddFace.H"
  42 #include "ReadFields.H"
  43 #include "volFields.H"
  44 #include "surfaceFields.H"
  45 #include "ZoneIDs.H"
  46
  47 using namespace Foam;
  48
  49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  50
  51 void modifyOrAddFace
  52 (
  53     directTopoChange& meshMod,
  54     const face& f,
  55     const label faceI,
  56     const label own,
  57     const bool flipFaceFlux,
  58     const label newPatchI,
  59     const label zoneID,
  60     const bool zoneFlip,
  61     PackedBoolList& modifiedFace
  62 )
  63 {
  64     if (!modifiedFace[faceI])
  65     {
  66         // First usage of face. Modify.
  67         meshMod.setAction
  68         (
  69             polyModifyFace
  70             (
  71                 f,                          // modified face
  72                 faceI,                      // label of face
  73                 own,                        // owner
  74                 -1,                         // neighbour
  75                 flipFaceFlux,               // face flip
  76                 newPatchI,                  // patch for face
  77                 false,                      // remove from zone
  78                 zoneID,                     // zone for face
  79                 zoneFlip                    // face flip in zone
  80             )
  81         );
  82         modifiedFace[faceI] = 1;
  83     }
  84     else
  85     {
  86         // Second or more usage of face. Add.
  87         meshMod.setAction
  88         (
  89             polyAddFace
  90             (
  91                 f,                          // modified face
  92                 own,                        // owner
  93                 -1,                         // neighbour
  94                 -1,                         // master point
  95                 -1,                         // master edge
  96                 faceI,                      // master face
  97                 flipFaceFlux,               // face flip
  98                 newPatchI,                  // patch for face
  99                 zoneID,                     // zone for face
 100                 zoneFlip                    // face flip in zone
 101             )
 102         );
 103     }
 104 }
 105
 106
 107 label findPatchID(const polyMesh& mesh, const word& name)
 108 {
 109     label patchI = mesh.boundaryMesh().findPatchID(name);
 110
 111     if (patchI == -1)
 112     {
 113         FatalErrorIn("findPatchID(const polyMesh&, const word&)")
 114             << "Cannot find patch " << name << endl
 115             << "Valid patches are " << mesh.boundaryMesh().names()
 116             << exit(FatalError);
 117     }
 118     return patchI;
 119 }
 120
 121
 122 // Main program:
 123
 124 int main(int argc, char *argv[])
 125 {
 126 #   include "addRegionOption.H"
 127     argList::validArgs.append("faceZone");
 128     argList::validArgs.append("patch");
 129     argList::validOptions.insert("additionalPatches", "(patch2 .. patchN)");
 130     argList::validOptions.insert("internalFacesOnly", "");
 131     argList::validOptions.insert("overwrite", "");
 132
 133 #   include "setRootCase.H"
 134 #   include "createTime.H"
 135     runTime.functionObjects().off();
 136 #   include "createNamedMesh.H"
 137     const word oldInstance = mesh.pointsInstance();
 138
 139     const polyBoundaryMesh& patches = mesh.boundaryMesh();
 140     const faceZoneMesh& faceZones = mesh.faceZones();
 141
 142     // Faces to baffle
 143     faceZoneID zoneID(args.additionalArgs()[0], faceZones);
 144
 145     Info<< "Converting faces on zone " << zoneID.name()
 146         << " into baffles." << nl << endl;
 147
 148     if (zoneID.index() == -1)
 149     {
 150         FatalErrorIn(args.executable()) << "Cannot find faceZone "
 151             << zoneID.name() << endl
 152             << "Valid zones are " << faceZones.names()
 153             << exit(FatalError);
 154     }
 155
 156     const faceZone& fZone = faceZones[zoneID.index()];
 157
 158     Info<< "Found " << returnReduce(fZone.size(), sumOp<label>())
 159         << " faces on zone " << zoneID.name() << nl << endl;
 160
 161     // Make sure patches and zoneFaces are synchronised across couples
 162     patches.checkParallelSync(true);
 163     fZone.checkParallelSync(true);
 164
 165     // Patches to put baffles into
 166     DynamicList<label> newPatches(1);
 167
 168     word patchName(args.additionalArgs()[1]);
 169     newPatches.append(findPatchID(mesh, patchName));
 170     Info<< "Using patch " << patchName
 171         << " at index " << newPatches[0] << endl;
 172
 173
 174     // Additional patches
 175     if (args.optionFound("additionalPatches"))
 176     {
 177         const wordList patchNames
 178         (
 179             args.optionLookup("additionalPatches")()
 180         );
 181
 182         newPatches.reserve(patchNames.size() + 1);
 183         forAll(patchNames, i)
 184         {
 185             newPatches.append(findPatchID(mesh, patchNames[i]));
 186             Info<< "Using additional patch " << patchNames[i]
 187                 << " at index " << newPatches[newPatches.size()-1] << endl;
 188         }
 189     }
 190
 191
 192     bool overwrite = args.optionFound("overwrite");
 193
 194     bool internalFacesOnly = args.optionFound("internalFacesOnly");
 195
 196     if (internalFacesOnly)
 197     {
 198         Info<< "Not converting faces on non-coupled patches." << nl << endl;
 199     }
 200
 201
 202     // Read objects in time directory
 203     IOobjectList objects(mesh, runTime.timeName());
 204
 205     // Read vol fields.
 206
 207     PtrList<volScalarField> vsFlds;
 208     ReadFields(mesh, objects, vsFlds);
 209
 210     PtrList<volVectorField> vvFlds;
 211     ReadFields(mesh, objects, vvFlds);
 212
 213     PtrList<volSphericalTensorField> vstFlds;
 214     ReadFields(mesh, objects, vstFlds);
 215
 216     PtrList<volSymmTensorField> vsymtFlds;
 217     ReadFields(mesh, objects, vsymtFlds);
 218
 219     PtrList<volTensorField> vtFlds;
 220     ReadFields(mesh, objects, vtFlds);
 221
 222     // Read surface fields.
 223
 224     PtrList<surfaceScalarField> ssFlds;
 225     ReadFields(mesh, objects, ssFlds);
 226
 227     PtrList<surfaceVectorField> svFlds;
 228     ReadFields(mesh, objects, svFlds);
 229
 230     PtrList<surfaceSphericalTensorField> sstFlds;
 231     ReadFields(mesh, objects, sstFlds);
 232
 233     PtrList<surfaceSymmTensorField> ssymtFlds;
 234     ReadFields(mesh, objects, ssymtFlds);
 235
 236     PtrList<surfaceTensorField> stFlds;
 237     ReadFields(mesh, objects, stFlds);
 238
 239
 240     // Mesh change container
 241     directTopoChange meshMod(mesh);
 242
 243
 244     // Do the actual changes. Note:
 245     // - loop in incrementing face order (not necessary if faceZone ordered).
 246     //   Preserves any existing ordering on patch faces.
 247     // - two passes, do non-flip faces first and flip faces second. This
 248     //   guarantees that when e.g. creating a cyclic all faces from one
 249     //   side come first and faces from the other side next.
 250
 251     // Whether first use of face (modify) or consecutive (add)
 252     PackedBoolList modifiedFace(mesh.nFaces());
 253     // Never modify coupled faces
 254     forAll(patches, patchI)
 255     {
 256         const polyPatch& pp = patches[patchI];
 257         if (pp.coupled())
 258         {
 259             forAll(pp, i)
 260             {
 261                 modifiedFace[pp.start()+i] = 1;
 262             }
 263         }
 264     }
 265     label nModified = 0;
 266
 267     forAll(newPatches, i)
 268     {
 269         label newPatchI = newPatches[i];
 270
 271         // Pass 1. Do selected side of zone
 272         // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 273
 274         for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
 275         {
 276             label zoneFaceI = fZone.whichFace(faceI);
 277
 278             if (zoneFaceI != -1)
 279             {
 280                 if (!fZone.flipMap()[zoneFaceI])
 281                 {
 282                     // Use owner side of face
 283                     modifyOrAddFace
 284                     (
 285                         meshMod,
 286                         mesh.faces()[faceI],    // modified face
 287                         faceI,                  // label of face
 288                         mesh.faceOwner()[faceI],// owner
 289                         false,                  // face flip
 290                         newPatchI,              // patch for face
 291                         zoneID.index(),         // zone for face
 292                         false,                  // face flip in zone
 293                         modifiedFace            // modify or add status
 294                     );
 295                 }
 296                 else
 297                 {
 298                     // Use neighbour side of face
 299                     modifyOrAddFace
 300                     (
 301                         meshMod,
 302                         mesh.faces()[faceI].reverseFace(),  // modified face
 303                         faceI,                      // label of face
 304                         mesh.faceNeighbour()[faceI],// owner
 305                         true,                       // face flip
 306                         newPatchI,                  // patch for face
 307                         zoneID.index(),             // zone for face
 308                         true,                       // face flip in zone
 309                         modifiedFace                // modify or add status
 310                     );
 311                 }
 312
 313                 nModified++;
 314             }
 315         }
 316
 317
 318         // Pass 2. Do other side of zone
 319         // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 320
 321         for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
 322         {
 323             label zoneFaceI = fZone.whichFace(faceI);
 324
 325             if (zoneFaceI != -1)
 326             {
 327                 if (!fZone.flipMap()[zoneFaceI])
 328                 {
 329                     // Use neighbour side of face
 330                     modifyOrAddFace
 331                     (
 332                         meshMod,
 333                         mesh.faces()[faceI].reverseFace(),  // modified face
 334                         faceI,                              // label of face
 335                         mesh.faceNeighbour()[faceI],        // owner
 336                         true,                               // face flip
 337                         newPatchI,                          // patch for face
 338                         zoneID.index(),                     // zone for face
 339                         true,                               // face flip in zone
 340                         modifiedFace                        // modify or add
 341                     );
 342                 }
 343                 else
 344                 {
 345                     // Use owner side of face
 346                     modifyOrAddFace
 347                     (
 348                         meshMod,
 349                         mesh.faces()[faceI],    // modified face
 350                         faceI,                  // label of face
 351                         mesh.faceOwner()[faceI],// owner
 352                         false,                  // face flip
 353                         newPatchI,              // patch for face
 354                         zoneID.index(),         // zone for face
 355                         false,                  // face flip in zone
 356                         modifiedFace            // modify or add status
 357                     );
 358                 }
 359             }
 360         }
 361
 362
 363         // Modify any boundary faces
 364         // ~~~~~~~~~~~~~~~~~~~~~~~~~
 365
 366         // Normal boundary:
 367         // - move to new patch. Might already be back-to-back baffle
 368         // you want to add cyclic to. Do warn though.
 369         //
 370         // Processor boundary:
 371         // - do not move to cyclic
 372         // - add normal patches though.
 373
 374         // For warning once per patch.
 375         labelHashSet patchWarned;
 376
 377         forAll(patches, patchI)
 378         {
 379             const polyPatch& pp = patches[patchI];
 380
 381             if (pp.coupled() && patches[newPatchI].coupled())
 382             {
 383                 // Do not allow coupled faces to be moved to different coupled
 384                 // patches.
 385             }
 386             else if (pp.coupled() || !internalFacesOnly)
 387             {
 388                 forAll(pp, i)
 389                 {
 390                     label faceI = pp.start()+i;
 391
 392                     label zoneFaceI = fZone.whichFace(faceI);
 393
 394                     if (zoneFaceI != -1)
 395                     {
 396                         if (patchWarned.insert(patchI))
 397                         {
 398                             WarningIn(args.executable())
 399                                 << "Found boundary face (in patch "
 400                                 << pp.name()
 401                                 << ") in faceZone " << fZone.name()
 402                                 << " to convert to baffle patch "
 403                                 << patches[newPatchI].name()
 404                                 << endl
 405                                 << "    Run with -internalFacesOnly option"
 406                                 << " if you don't wish to convert"
 407                                 << " boundary faces." << endl;
 408                         }
 409
 410                         modifyOrAddFace
 411                         (
 412                             meshMod,
 413                             mesh.faces()[faceI],        // modified face
 414                             faceI,                      // label of face
 415                             mesh.faceOwner()[faceI],    // owner
 416                             false,                      // face flip
 417                             newPatchI,                  // patch for face
 418                             zoneID.index(),             // zone for face
 419                             fZone.flipMap()[zoneFaceI], // face flip in zone
 420                             modifiedFace                // modify or add status
 421                         );
 422
 423                         nModified++;
 424                     }
 425                 }
 426             }
 427         }
 428     }
 429
 430
 431     Info<< "Converted " << returnReduce(nModified, sumOp<label>())
 432         << " faces into boundary faces on patch " << patchName << nl << endl;
 433
 434     if (!overwrite)
 435     {
 436         runTime++;
 437     }
 438
 439     // Change the mesh. Change points directly (no inflation).
 440     autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false);
 441
 442     // Update fields
 443     mesh.updateMesh(map);
 444
 445     // Move mesh (since morphing might not do this)
 446     if (map().hasMotionPoints())
 447     {
 448         mesh.movePoints(map().preMotionPoints());
 449     }
 450
 451     if (overwrite)
 452     {
 453         mesh.setInstance(oldInstance);
 454     }
 455     Info<< "Writing mesh to " << runTime.timeName() << endl;
 456
 457     mesh.write();
 458
 459     Info<< "End\n" << endl;
 460
 461     return 0;
 462 }
 463
 464
 465 // ************************************************************************* //