applications/utilities/mesh/advanced/refineWallLayer/refineWallLayer.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
   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 Description
  25     Utility to refine cells next to patches.
  26
  27     Takes a patchName and number of layers to refine. Works out cells within
  28     these layers and refines those in the wall-normal direction.
  29
  30 \*---------------------------------------------------------------------------*/
  31
  32 #include "argList.H"
  33 #include "Time.H"
  34 #include "polyTopoChange.H"
  35 #include "polyTopoChanger.H"
  36 #include "mapPolyMesh.H"
  37 #include "polyMesh.H"
  38 #include "cellCuts.H"
  39 #include "cellSet.H"
  40 #include "meshCutter.H"
  41
  42 using namespace Foam;
  43
  44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  45 // Main program:
  46
  47 int main(int argc, char *argv[])
  48 {
  49     #include "addOverwriteOption.H"
  50     argList::noParallel();
  51     argList::validArgs.append("patchName");
  52     argList::validArgs.append("edgeWeight");
  53
  54     argList::addOption
  55     (
  56         "useSet",
  57         "name",
  58         "restrict cells to refine based on specified cellSet name"
  59     );
  60
  61
  62 #   include "setRootCase.H"
  63 #   include "createTime.H"
  64     runTime.functionObjects().off();
  65 #   include "createPolyMesh.H"
  66     const word oldInstance = mesh.pointsInstance();
  67
  68     const word patchName = args[1];
  69     const scalar weight  = args.argRead<scalar>(2);
  70     const bool overwrite = args.optionFound("overwrite");
  71
  72     label patchID = mesh.boundaryMesh().findPatchID(patchName);
  73
  74     if (patchID == -1)
  75     {
  76         FatalErrorIn(args.executable())
  77             << "Cannot find patch " << patchName << endl
  78             << "Valid patches are " << mesh.boundaryMesh().names()
  79             << exit(FatalError);
  80     }
  81     const polyPatch& pp = mesh.boundaryMesh()[patchID];
  82
  83
  84     // Cells cut
  85
  86     labelHashSet cutCells(4*pp.size());
  87
  88     const labelList& meshPoints = pp.meshPoints();
  89
  90     forAll(meshPoints, pointI)
  91     {
  92         label meshPointI = meshPoints[pointI];
  93
  94         const labelList& pCells = mesh.pointCells()[meshPointI];
  95
  96         forAll(pCells, pCellI)
  97         {
  98             cutCells.insert(pCells[pCellI]);
  99         }
 100     }
 101
 102     Info<< "Selected " << cutCells.size()
 103         << " cells connected to patch " << pp.name() << endl << endl;
 104
 105     //
 106     // List of cells to refine
 107     //
 108
 109     word setName;
 110     if (args.optionReadIfPresent("useSet", setName))
 111     {
 112         Info<< "Subsetting cells to cut based on cellSet"
 113             << setName << nl << endl;
 114
 115         cellSet cells(mesh, setName);
 116
 117         Info<< "Read " << cells.size() << " cells from cellSet "
 118             << cells.instance()/cells.local()/cells.name()
 119             << nl << endl;
 120
 121         forAllConstIter(cellSet, cells, iter)
 122         {
 123             cutCells.erase(iter.key());
 124         }
 125         Info<< "Removed from cells to cut all the ones not in set "
 126             << setName << nl << endl;
 127     }
 128
 129     // Mark all meshpoints on patch
 130
 131     boolList vertOnPatch(mesh.nPoints(), false);
 132
 133     forAll(meshPoints, pointI)
 134     {
 135         const label meshPointI = meshPoints[pointI];
 136
 137         vertOnPatch[meshPointI] = true;
 138     }
 139
 140
 141     // Mark cut edges.
 142
 143     DynamicList<label> allCutEdges(pp.nEdges());
 144
 145     DynamicList<scalar> allCutEdgeWeights(pp.nEdges());
 146
 147     forAll(meshPoints, pointI)
 148     {
 149         label meshPointI = meshPoints[pointI];
 150
 151         const labelList& pEdges = mesh.pointEdges()[meshPointI];
 152
 153         forAll(pEdges, pEdgeI)
 154         {
 155             const label edgeI = pEdges[pEdgeI];
 156             const edge& e = mesh.edges()[edgeI];
 157
 158             label otherPointI = e.otherVertex(meshPointI);
 159
 160             if (!vertOnPatch[otherPointI])
 161             {
 162                 allCutEdges.append(edgeI);
 163
 164                 if (e.start() == meshPointI)
 165                 {
 166                     allCutEdgeWeights.append(weight);
 167                 }
 168                 else
 169                 {
 170                     allCutEdgeWeights.append(1 - weight);
 171                 }
 172             }
 173         }
 174     }
 175
 176     allCutEdges.shrink();
 177     allCutEdgeWeights.shrink();
 178
 179     Info<< "Cutting:" << nl
 180         << "    cells:" << cutCells.size() << nl
 181         << "    edges:" << allCutEdges.size() << nl
 182         << endl;
 183
 184     // Transfer DynamicLists to straight ones.
 185     scalarField cutEdgeWeights;
 186     cutEdgeWeights.transfer(allCutEdgeWeights);
 187     allCutEdgeWeights.clear();
 188
 189
 190     // Gets cuts across cells from cuts through edges.
 191     cellCuts cuts
 192     (
 193         mesh,
 194         cutCells.toc(),     // cells candidate for cutting
 195         labelList(0),       // cut vertices
 196         allCutEdges,        // cut edges
 197         cutEdgeWeights      // weight on cut edges
 198     );
 199
 200     polyTopoChange meshMod(mesh);
 201
 202     // Cutting engine
 203     meshCutter cutter(mesh);
 204
 205     // Insert mesh refinement into polyTopoChange.
 206     cutter.setRefinement(cuts, meshMod);
 207
 208     // Do all changes
 209     Info<< "Morphing ..." << endl;
 210
 211     if (!overwrite)
 212     {
 213         runTime++;
 214     }
 215
 216     autoPtr<mapPolyMesh> morphMap = meshMod.changeMesh(mesh, false);
 217
 218     if (morphMap().hasMotionPoints())
 219     {
 220         mesh.movePoints(morphMap().preMotionPoints());
 221     }
 222
 223     // Update stored labels on meshCutter.
 224     cutter.updateMesh(morphMap());
 225
 226     if (overwrite)
 227     {
 228         mesh.setInstance(oldInstance);
 229     }
 230
 231     // Write resulting mesh
 232     Info<< "Writing refined morphMesh to time " << runTime.timeName() << endl;
 233
 234     mesh.write();
 235
 236     Info<< "End\n" << endl;
 237
 238     return 0;
 239 }
 240
 241
 242 // ************************************************************************* //