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