diff --git a/applications/utilities/mesh/generation/cvMesh/cvMeshDict b/applications/utilities/mesh/generation/cvMesh/cvMeshDict
index 1a07724361..6b8703c799 100644
--- a/applications/utilities/mesh/generation/cvMesh/cvMeshDict
+++ b/applications/utilities/mesh/generation/cvMesh/cvMeshDict
@@ -210,16 +210,15 @@ initialPoints
// to the surface. Is fraction of local target cell size (see below)
minimumSurfaceDistanceCoeff 0.55;
- initialPointsMethod hierarchicalDensityWeightedStochastic;
+ initialPointsMethod autoDensity;
// initialPointsMethod uniformGrid;
// initialPointsMethod bodyCentredCubic;
// initialPointsMethod pointFile;
- // initialPointsMethod densityWeightedStochastic;
// Take boundbox of all geometry. Sample with this box. If too much
// samples in box (due to target cell size) split box.
- hierarchicalDensityWeightedStochasticDetails
+ autoDensityDetails
{
// Initial number of refinement levels. Needs to be enough to pick
// up features due to size ratio. If not enough it will take longer
@@ -233,12 +232,6 @@ initialPoints
surfaceSampleResolution 5;
}
- densityWeightedStochasticDetails
- {
- // Volume of space to be meshed. Use e.g. surfaceInertia to find this.
- totalVolume 1.56e-05;
- }
-
uniformGridDetails
{
// Absolute cell size.
diff --git a/applications/utilities/mesh/manipulation/topoSet/Make/files b/applications/utilities/mesh/manipulation/topoSet/Make/files
index 872e324100..5514b34344 100644
--- a/applications/utilities/mesh/manipulation/topoSet/Make/files
+++ b/applications/utilities/mesh/manipulation/topoSet/Make/files
@@ -1,3 +1,4 @@
+timeSelector.C
topoSet.C
EXE = $(FOAM_APPBIN)/topoSet
diff --git a/applications/utilities/mesh/manipulation/topoSet/timeSelector.C b/applications/utilities/mesh/manipulation/topoSet/timeSelector.C
new file mode 100644
index 0000000000..6fa85f58d5
--- /dev/null
+++ b/applications/utilities/mesh/manipulation/topoSet/timeSelector.C
@@ -0,0 +1,289 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see .
+
+\*---------------------------------------------------------------------------*/
+
+#include "timeSelector.H"
+#include "ListOps.H"
+#include "argList.H"
+#include "Time.H"
+#include "IStringStream.H"
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::timeSelector::timeSelector()
+:
+ scalarRanges()
+{}
+
+
+Foam::timeSelector::timeSelector(Istream& is)
+:
+ scalarRanges(is)
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+bool Foam::timeSelector::selected(const instant& value) const
+{
+ return scalarRanges::selected(value.value());
+}
+
+
+Foam::List Foam::timeSelector::selected(const instantList& Times) const
+{
+ List lst(Times.size(), false);
+
+ // check ranges, avoid false positive on constant/
+ forAll(Times, timeI)
+ {
+ if (Times[timeI].name() != "constant" && selected(Times[timeI]))
+ {
+ lst[timeI] = true;
+ }
+ }
+
+ // check specific values
+ forAll(*this, rangeI)
+ {
+ if (operator[](rangeI).isExact())
+ {
+ scalar target = operator[](rangeI).value();
+
+ int nearestIndex = -1;
+ scalar nearestDiff = Foam::GREAT;
+
+ forAll(Times, timeI)
+ {
+ if (Times[timeI].name() == "constant") continue;
+
+ scalar diff = fabs(Times[timeI].value() - target);
+ if (diff < nearestDiff)
+ {
+ nearestDiff = diff;
+ nearestIndex = timeI;
+ }
+ }
+
+ if (nearestIndex >= 0)
+ {
+ lst[nearestIndex] = true;
+ }
+ }
+ }
+
+ return lst;
+}
+
+
+Foam::List Foam::timeSelector::select
+(
+ const instantList& Times
+) const
+{
+ return subset(selected(Times), Times);
+}
+
+
+void Foam::timeSelector::inplaceSelect(instantList& Times) const
+{
+ inplaceSubset(selected(Times), Times);
+}
+
+
+void Foam::timeSelector::addOptions
+(
+ const bool constant,
+ const bool zeroTime
+)
+{
+ if (constant)
+ {
+ argList::addBoolOption
+ (
+ "constant",
+ "include the 'constant/' dir in the times list"
+ );
+ }
+ if (zeroTime)
+ {
+ argList::addBoolOption
+ (
+ "zeroTime",
+ "include the '0/' dir in the times list"
+ );
+ }
+ argList::addBoolOption
+ (
+ "noZero",
+ "exclude the '0/' dir from the times list, "
+ "has precedence over the -zeroTime option"
+ );
+ argList::addBoolOption
+ (
+ "latestTime",
+ "select the latest time"
+ );
+ argList::addOption
+ (
+ "time",
+ "ranges",
+ "comma-separated time ranges - eg, ':10,20,40-70,1000:'"
+ );
+}
+
+
+Foam::List Foam::timeSelector::select
+(
+ const instantList& timeDirs,
+ const argList& args
+)
+{
+ if (timeDirs.size())
+ {
+ List selectTimes(timeDirs.size(), true);
+
+ // determine locations of constant/ and 0/ directories
+ label constantIdx = -1;
+ label zeroIdx = -1;
+
+ forAll(timeDirs, timeI)
+ {
+ if (timeDirs[timeI].name() == "constant")
+ {
+ constantIdx = timeI;
+ }
+ else if (timeDirs[timeI].value() == 0)
+ {
+ zeroIdx = timeI;
+ }
+
+ if (constantIdx >= 0 && zeroIdx >= 0)
+ {
+ break;
+ }
+ }
+
+ // determine latestTime selection (if any)
+ // this must appear before the -time option processing
+ label latestIdx = -1;
+ if (args.optionFound("latestTime"))
+ {
+ selectTimes = false;
+ latestIdx = timeDirs.size() - 1;
+
+ // avoid false match on constant/
+ if (latestIdx == constantIdx)
+ {
+ latestIdx = -1;
+ }
+ }
+
+ if (args.optionFound("time"))
+ {
+ // can match 0/, but can never match constant/
+ selectTimes = timeSelector
+ (
+ args.optionLookup("time")()
+ ).selected(timeDirs);
+ }
+
+
+ // add in latestTime (if selected)
+ if (latestIdx >= 0)
+ {
+ selectTimes[latestIdx] = true;
+ }
+
+ if (constantIdx >= 0)
+ {
+ // only add constant/ if specifically requested
+ selectTimes[constantIdx] = args.optionFound("constant");
+ }
+
+ // special treatment for 0/
+ if (zeroIdx >= 0)
+ {
+ if (args.optionFound("noZero"))
+ {
+ // exclude 0/ if specifically requested
+ selectTimes[zeroIdx] = false;
+ }
+ else if (argList::validOptions.found("zeroTime"))
+ {
+ // with -zeroTime enabled, drop 0/ unless specifically requested
+ selectTimes[zeroIdx] = args.optionFound("zeroTime");
+ }
+ }
+
+ return subset(selectTimes, timeDirs);
+ }
+ else
+ {
+ return timeDirs;
+ }
+}
+
+
+Foam::List Foam::timeSelector::select0
+(
+ Time& runTime,
+ const argList& args,
+ const bool useOptionsOnly
+)
+{
+ if
+ (
+ useOptionsOnly
+ || (
+ args.optionFound("latestTime")
+ || args.optionFound("time")
+ || args.optionFound("constant")
+ || args.optionFound("noZero")
+ || args.optionFound("zeroTime")
+ )
+ )
+ {
+ instantList timeDirs = timeSelector::select(runTime.times(), args);
+
+ if (timeDirs.empty())
+ {
+ FatalErrorIn(args.executable())
+ << "No times selected"
+ << exit(FatalError);
+ }
+
+ runTime.setTime(timeDirs[0], 0);
+
+ return timeDirs;
+ }
+ else
+ {
+ // No timeSelector option specified. Do not change runTime.
+ return instantList(1, instant(runTime.value(), runTime.timeName()));
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/manipulation/topoSet/timeSelector.H b/applications/utilities/mesh/manipulation/topoSet/timeSelector.H
new file mode 100644
index 0000000000..daf821c905
--- /dev/null
+++ b/applications/utilities/mesh/manipulation/topoSet/timeSelector.H
@@ -0,0 +1,163 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see .
+
+Class
+ Foam::timeSelector
+
+Description
+ A List of scalarRange for selecting times.
+
+ The timeSelector provides a convenient means of selecting multiple
+ times. A typical use would be the following:
+
+ \verbatim
+ timeSelector::addOptions();
+ // add other options
+ #include "setRootCase.H"
+ #include "createTime.H"
+ instantList timeDirs = timeSelector::select0(runTime, args);
+ ...
+ forAll(timeDirs, timeI)
+ {
+ ...
+ }
+ \endverbatim
+
+ The result program would receive \b -time, @b -latestTime, @b -constant
+ and \b -noZero options. The @b -constant option explicitly includes the
+ \c constant/ directory in the time list and the \b -noZero option
+ explicitly excludes the \c 0/ directory from the time list.
+
+ There may however also be many cases in which neither the \c constant/
+ directory nor the \c 0/ directory contain particularly relevant
+ information. This might occur, for example, when post-processing
+ results. In this case, addOptions is called with optional boolean
+ arguments.
+
+ \verbatim
+ timeSelector::addOptions(false, true);
+ \endverbatim
+
+ The first argument avoids adding the \b -constant option. The second
+ argument adds an additional \b -zeroTime option and also prevents the
+ \c 0/ directory from being included in the default time range and in the
+ \b -latestTime selection.
+
+SourceFiles
+ timeSelector.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef timeSelector_H
+#define timeSelector_H
+
+#include "scalarRanges.H"
+#include "instantList.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward declaration of classes
+class argList;
+class Time;
+
+/*---------------------------------------------------------------------------*\
+ Class timeSelector Declaration
+\*---------------------------------------------------------------------------*/
+
+class timeSelector
+:
+ public scalarRanges
+{
+public:
+
+ // Constructors
+
+ //- Construct Null
+ timeSelector();
+
+ //- Construct from Istream
+ timeSelector(Istream&);
+
+
+ // Member Functions
+
+ //- Return true if the given instant is within the ranges
+ bool selected(const instant&) const;
+
+ //- Return the set of selected instants in the given list that are
+ // within the ranges
+ List selected(const List&) const;
+
+ //- Select a list of Time values that are within the ranges
+ instantList select(const List&) const;
+
+ //- Select a list of Time values that are within the ranges
+ void inplaceSelect(List&) const;
+
+ //- Add the options handled by timeSelector to argList::validOptions
+ //
+ // \param constant
+ // Add the \b -constant option to include the \c constant/ directory
+ //
+ // \param zeroTime
+ // Enable the \b -zeroTime option and alter the normal time selection
+ // behaviour (and \b -latestTime behaviour) to exclude the \c 0/
+ // directory. The \c 0/ directory will only be included when
+ // \b -zeroTime is specified.
+ // The \b -noZero option has precedence over the @b -zeroTime option.
+ static void addOptions
+ (
+ const bool constant=true,
+ const bool zeroTime=false
+ );
+
+ //- Return the set of times selected based on the argList options
+ static instantList select
+ (
+ const instantList&,
+ const argList& args
+ );
+
+ //- Return the set of times selected based on the argList options
+ // also set the runTime to the first instance
+ static instantList select0
+ (
+ Time& runTime,
+ const argList& args,
+ const bool useOptionsOnly
+ );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/manipulation/topoSet/topoSet.C b/applications/utilities/mesh/manipulation/topoSet/topoSet.C
index a377f6133e..d34cc58257 100644
--- a/applications/utilities/mesh/manipulation/topoSet/topoSet.C
+++ b/applications/utilities/mesh/manipulation/topoSet/topoSet.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
+ \\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -26,6 +26,7 @@ Description
\*---------------------------------------------------------------------------*/
+#include "timeSelector.H"
#include "argList.H"
#include "Time.H"
#include "polyMesh.H"
@@ -33,15 +34,71 @@ Description
#include "cellSet.H"
#include "faceSet.H"
#include "pointSet.H"
+#include "globalMeshData.H"
using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+void printMesh(const Time& runTime, const polyMesh& mesh)
+{
+ Info<< "Time:" << runTime.timeName()
+ << " cells:" << mesh.globalData().nTotalCells()
+ << " faces:" << mesh.globalData().nTotalFaces()
+ << " points:" << mesh.globalData().nTotalPoints()
+ << " patches:" << mesh.boundaryMesh().size()
+ << " bb:" << mesh.bounds() << nl;
+}
+
+
+polyMesh::readUpdateState meshReadUpdate(polyMesh& mesh)
+{
+ polyMesh::readUpdateState stat = mesh.readUpdate();
+
+ switch(stat)
+ {
+ case polyMesh::UNCHANGED:
+ {
+ Info<< " mesh not changed." << endl;
+ break;
+ }
+ case polyMesh::POINTS_MOVED:
+ {
+ Info<< " points moved; topology unchanged." << endl;
+ break;
+ }
+ case polyMesh::TOPO_CHANGE:
+ {
+ Info<< " topology changed; patches unchanged." << nl
+ << " ";
+ printMesh(mesh.time(), mesh);
+ break;
+ }
+ case polyMesh::TOPO_PATCH_CHANGE:
+ {
+ Info<< " topology changed and patches changed." << nl
+ << " ";
+ printMesh(mesh.time(), mesh);
+
+ break;
+ }
+ default:
+ {
+ FatalErrorIn("meshReadUpdate(polyMesh&)")
+ << "Illegal mesh update state "
+ << stat << abort(FatalError);
+ break;
+ }
+ }
+ return stat;
+}
+
+
// Main program:
int main(int argc, char *argv[])
{
+ timeSelector::addOptions(true, false);
argList::addOption
(
"dict",
@@ -57,6 +114,14 @@ int main(int argc, char *argv[])
# include "setRootCase.H"
# include "createTime.H"
+ instantList timeDirs = timeSelector::select0
+ (
+ runTime,
+ args,
+ false // not override runTime if no time options
+ );
+
+
# include "createNamedPolyMesh.H"
const bool noSync = args.optionFound("noSync");
@@ -99,128 +164,143 @@ int main(int argc, char *argv[])
// Read set construct info from dictionary
- PtrList patchSources(topoSetDict.lookup("actions"));
+ PtrList actions(topoSetDict.lookup("actions"));
- forAll(patchSources, i)
+
+ forAll(timeDirs, timeI)
{
- const dictionary& dict = patchSources[i];
+ runTime.setTime(timeDirs[timeI], timeI);
+ Info<< "Time = " << runTime.timeName() << endl;
- const word setName(dict.lookup("name"));
- const word actionName(dict.lookup("action"));
- const word setType(dict.lookup("type"));
+ // Optionally re-read mesh
+ meshReadUpdate(mesh);
-
- topoSetSource::setAction action = topoSetSource::toAction(actionName);
-
- autoPtr currentSet;
- if
- (
- (action == topoSetSource::NEW)
- || (action == topoSetSource::CLEAR)
- )
+ // Execute all actions
+ forAll(actions, i)
{
- currentSet = topoSet::New(setType, mesh, setName, 10000);
- Info<< "Created set " << setName << endl;
- }
- else if (action == topoSetSource::REMOVE)
- {
- //?
- }
- else
- {
- currentSet = topoSet::New
+ const dictionary& dict = actions[i];
+
+ const word setName(dict.lookup("name"));
+ const word actionName(dict.lookup("action"));
+ const word setType(dict.lookup("type"));
+
+
+ topoSetSource::setAction action = topoSetSource::toAction
(
- setType,
- mesh,
- setName,
- IOobject::MUST_READ
+ actionName
);
- Info<< "Read set " << setName << " with size "
- << currentSet().size() << endl;
- }
-
-
- // Handle special actions (clear, invert) locally, rest through sources.
- switch (action)
- {
- case topoSetSource::NEW:
- case topoSetSource::ADD:
- case topoSetSource::DELETE:
+ autoPtr currentSet;
+ if
+ (
+ (action == topoSetSource::NEW)
+ || (action == topoSetSource::CLEAR)
+ )
{
- Info<< " Applying source " << word(dict.lookup("source"))
- << endl;
- autoPtr source = topoSetSource::New
- (
- dict.lookup("source"),
- mesh,
- dict.subDict("sourceInfo")
- );
-
- source().applyToSet(action, currentSet());
- // Synchronize for coupled patches.
- if (!noSync) currentSet().sync(mesh);
- currentSet().write();
+ currentSet = topoSet::New(setType, mesh, setName, 10000);
+ Info<< "Created set " << setName << endl;
}
- break;
-
- case topoSetSource::SUBSET:
+ else if (action == topoSetSource::REMOVE)
{
- Info<< " Applying source " << word(dict.lookup("source"))
- << endl;
- autoPtr source = topoSetSource::New
+ //?
+ }
+ else
+ {
+ currentSet = topoSet::New
(
- dict.lookup("source"),
+ setType,
mesh,
- dict.subDict("sourceInfo")
+ setName,
+ IOobject::MUST_READ
);
+ Info<< "Read set " << setName << " with size "
+ << currentSet().size() << endl;
+ }
- // Backup current set.
- autoPtr oldSet
- (
- topoSet::New
+
+
+ // Handle special actions (clear, invert) locally, rest through
+ // sources.
+ switch (action)
+ {
+ case topoSetSource::NEW:
+ case topoSetSource::ADD:
+ case topoSetSource::DELETE:
+ {
+ Info<< " Applying source " << word(dict.lookup("source"))
+ << endl;
+ autoPtr source = topoSetSource::New
(
- setType,
+ dict.lookup("source"),
mesh,
- currentSet().name() + "_old2",
- currentSet()
- )
- );
+ dict.subDict("sourceInfo")
+ );
- currentSet().clear();
- source().applyToSet(topoSetSource::NEW, currentSet());
+ source().applyToSet(action, currentSet());
+ // Synchronize for coupled patches.
+ if (!noSync) currentSet().sync(mesh);
+ currentSet().write();
+ }
+ break;
- // Combine new value of currentSet with old one.
- currentSet().subset(oldSet());
- // Synchronize for coupled patches.
- if (!noSync) currentSet().sync(mesh);
- currentSet().write();
+ case topoSetSource::SUBSET:
+ {
+ Info<< " Applying source " << word(dict.lookup("source"))
+ << endl;
+ autoPtr source = topoSetSource::New
+ (
+ dict.lookup("source"),
+ mesh,
+ dict.subDict("sourceInfo")
+ );
+
+ // Backup current set.
+ autoPtr oldSet
+ (
+ topoSet::New
+ (
+ setType,
+ mesh,
+ currentSet().name() + "_old2",
+ currentSet()
+ )
+ );
+
+ currentSet().clear();
+ source().applyToSet(topoSetSource::NEW, currentSet());
+
+ // Combine new value of currentSet with old one.
+ currentSet().subset(oldSet());
+ // Synchronize for coupled patches.
+ if (!noSync) currentSet().sync(mesh);
+ currentSet().write();
+ }
+ break;
+
+ case topoSetSource::CLEAR:
+ Info<< " Clearing set" << endl;
+ currentSet().clear();
+ currentSet().write();
+ break;
+
+ case topoSetSource::INVERT:
+ Info<< " Inverting set" << endl;
+ currentSet().invert(currentSet().maxSize(mesh));
+ currentSet().write();
+ break;
+
+ default:
+ WarningIn(args.executable())
+ << "Unhandled action " << action << endl;
+ break;
}
- break;
- case topoSetSource::CLEAR:
- Info<< " Clearing set" << endl;
- currentSet().clear();
- currentSet().write();
- break;
-
- case topoSetSource::INVERT:
- Info<< " Inverting set" << endl;
- currentSet().invert(currentSet().maxSize(mesh));
- currentSet().write();
- break;
-
- default:
- WarningIn(args.executable())
- << "Unhandled action " << action << endl;
- break;
- }
-
- if (currentSet.valid())
- {
- Info<< " Set " << currentSet().name()
- << " now size " << currentSet().size()
- << endl;
+ if (currentSet.valid())
+ {
+ Info<< " Set " << currentSet().name()
+ << " now size " << currentSet().size()
+ << endl;
+ }
}
}
diff --git a/applications/utilities/surface/surfaceClean/collapseBase.C b/applications/utilities/surface/surfaceClean/collapseBase.C
index df38af390b..ff982c571f 100644
--- a/applications/utilities/surface/surfaceClean/collapseBase.C
+++ b/applications/utilities/surface/surfaceClean/collapseBase.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
+ \\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -66,8 +66,8 @@ static void writeRegionOBJ
triSurface regionSurf(surf.subsetMesh(include, pointMap, faceMap));
- //Pout<< "Region " << regionI << " surface:" << nl;
- //regionSurf.writeStats(Pout);
+ Pout<< "Region " << regionI << " surface:" << nl;
+ regionSurf.writeStats(Pout);
regionSurf.write(regionName);
@@ -97,7 +97,7 @@ static void splitTri
DynamicList& tris
)
{
- label oldNTris = tris.size();
+ //label oldNTris = tris.size();
label fp = findIndex(f, e[0]);
label fp1 = f.fcIndex(fp);
@@ -175,7 +175,7 @@ static void splitTri
}
else
{
- FatalErrorIn("splitTri")
+ FatalErrorIn("splitTri(..)")
<< "Edge " << e << " not part of triangle " << f
<< " fp:" << fp
<< " fp1:" << fp1
@@ -183,13 +183,13 @@ static void splitTri
<< abort(FatalError);
}
- Pout<< "Split face " << f << " along edge " << e
- << " into triangles:" << endl;
-
- for (label i = oldNTris; i < tris.size(); i++)
- {
- Pout<< " " << tris[i] << nl;
- }
+ //Pout<< "Split face " << f << " along edge " << e
+ // << " into triangles:" << endl;
+ //
+ //for (label i = oldNTris; i < tris.size(); i++)
+ //{
+ // Pout<< " " << tris[i] << nl;
+ //}
}
@@ -206,14 +206,14 @@ static bool insertSorted
{
if (findIndex(sortedVerts, vertI) != -1)
{
- FatalErrorIn("insertSorted") << "Trying to insert vertex " << vertI
+ FatalErrorIn("insertSorted(..)") << "Trying to insert vertex " << vertI
<< " which is already in list of sorted vertices "
<< sortedVerts << abort(FatalError);
}
if (weight <= 0 || weight >= 1)
{
- FatalErrorIn("insertSorted") << "Trying to insert vertex " << vertI
+ FatalErrorIn("insertSorted(..)") << "Trying to insert vertex " << vertI
<< " with illegal weight " << weight
<< " into list of sorted vertices "
<< sortedVerts << abort(FatalError);
@@ -228,7 +228,7 @@ static bool insertSorted
if (mag(w - weight) < SMALL)
{
- WarningIn("insertSorted")
+ WarningIn("insertSorted(..)")
<< "Trying to insert weight " << weight << " which is close to"
<< " existing weight " << w << " in " << sortedWeights
<< endl;
@@ -263,64 +263,103 @@ static bool insertSorted
}
+// Is triangle candidate for collapse? Small height or small quality
+bool isSliver
+(
+ const triSurface& surf,
+ const scalar minLen,
+ const scalar minQuality,
+ const label faceI,
+ const label edgeI
+)
+{
+ const pointField& localPoints = surf.localPoints();
+
+ // Check
+ // - opposite vertex projects onto base edge
+ // - normal distance is small
+ // - or triangle quality is small
+
+ label opposite0 =
+ triSurfaceTools::oppositeVertex
+ (
+ surf,
+ faceI,
+ edgeI
+ );
+
+ const edge& e = surf.edges()[edgeI];
+ const labelledTri& f = surf[faceI];
+
+ pointHit pHit =
+ e.line(localPoints).nearestDist
+ (
+ localPoints[opposite0]
+ );
+
+ if
+ (
+ pHit.hit()
+ && (
+ pHit.distance() < minLen
+ || f.tri(surf.points()).quality() < minQuality
+ )
+ )
+ {
+ // Remove faceI and split all other faces using this
+ // edge. This is done by 'replacing' the edgeI with the
+ // opposite0 vertex
+ //Pout<< "Splitting face " << faceI << " since distance "
+ // << pHit.distance()
+ // << " from vertex " << opposite0
+ // << " to edge " << edgeI
+ // << " points "
+ // << localPoints[e[0]]
+ // << localPoints[e[1]]
+ // << " is too small or triangle quality "
+ // << f.tri(surf.points()).quality()
+ // << " too small." << endl;
+
+ return true;
+ }
+ else
+ {
+ return false;
+ }
+}
+
+
// Mark all faces that are going to be collapsed.
// faceToEdge: per face -1 or the base edge of the face.
static void markCollapsedFaces
(
const triSurface& surf,
const scalar minLen,
+ const scalar minQuality,
labelList& faceToEdge
)
{
faceToEdge.setSize(surf.size());
faceToEdge = -1;
- const pointField& localPoints = surf.localPoints();
const labelListList& edgeFaces = surf.edgeFaces();
forAll(edgeFaces, edgeI)
{
- const edge& e = surf.edges()[edgeI];
-
const labelList& eFaces = surf.edgeFaces()[edgeI];
forAll(eFaces, i)
{
label faceI = eFaces[i];
- // Check distance of vertex to edge.
- label opposite0 =
- triSurfaceTools::oppositeVertex
- (
- surf,
- faceI,
- edgeI
- );
+ bool isCandidate = isSliver(surf, minLen, minQuality, faceI, edgeI);
- pointHit pHit =
- e.line(localPoints).nearestDist
- (
- localPoints[opposite0]
- );
-
- if (pHit.hit() && pHit.distance() < minLen)
+ if (isCandidate)
{
- // Remove faceI and split all other faces using this
- // edge. This is done by 'replacing' the edgeI with the
- // opposite0 vertex
- Pout<< "Splitting face " << faceI << " since distance "
- << pHit.distance()
- << " from vertex " << opposite0
- << " to edge " << edgeI
- << " points "
- << localPoints[e[0]]
- << localPoints[e[1]]
- << " is too small" << endl;
-
// Mark face as being collapsed
if (faceToEdge[faceI] != -1)
{
- FatalErrorIn("markCollapsedFaces")
+ FatalErrorIn("markCollapsedFaces(..)")
<< "Cannot collapse face " << faceI << " since "
<< " is marked to be collapsed both to edge "
<< faceToEdge[faceI] << " and " << edgeI
@@ -347,7 +386,7 @@ static void markRegion
{
if (faceToEdge[faceI] == -1 || collapseRegion[faceI] != -1)
{
- FatalErrorIn("markRegion")
+ FatalErrorIn("markRegion(..)")
<< "Problem : crossed into uncollapsed/regionized face"
<< abort(FatalError);
}
@@ -383,7 +422,7 @@ static void markRegion
}
else if (collapseRegion[nbrFaceI] != regionI)
{
- FatalErrorIn("markRegion")
+ FatalErrorIn("markRegion(..)")
<< "Edge:" << edgeI << " between face " << faceI
<< " with region " << regionI
<< " and face " << nbrFaceI
@@ -411,8 +450,8 @@ static label markRegions
{
if (collapseRegion[faceI] == -1 && faceToEdge[faceI] != -1)
{
- Pout<< "markRegions : Marking region:" << regionI
- << " starting from face " << faceI << endl;
+ //Pout<< "markRegions : Marking region:" << regionI
+ // << " starting from face " << faceI << endl;
// Collapsed face. Mark connected region with current region number
markRegion(surf, faceToEdge, regionI++, faceI, collapseRegion);
@@ -728,7 +767,12 @@ static void getSplitVerts
}
-label collapseBase(triSurface& surf, const scalar minLen)
+label collapseBase
+(
+ triSurface& surf,
+ const scalar minLen,
+ const scalar minQuality
+)
{
label nTotalSplit = 0;
@@ -743,7 +787,7 @@ label collapseBase(triSurface& surf, const scalar minLen)
labelList faceToEdge(surf.size(), -1);
// Calculate faceToEdge (face collapses)
- markCollapsedFaces(surf, minLen, faceToEdge);
+ markCollapsedFaces(surf, minLen, minQuality, faceToEdge);
// Find regions of connected collapsed faces
@@ -754,8 +798,8 @@ label collapseBase(triSurface& surf, const scalar minLen)
label nRegions = markRegions(surf, faceToEdge, collapseRegion);
- Pout<< "Detected " << nRegions << " regions of faces to be collapsed"
- << nl << endl;
+ //Pout<< "Detected " << nRegions << " regions of faces to be collapsed"
+ // << nl << endl;
// Pick up all vertices on outside of region
labelListList outsideVerts
@@ -772,10 +816,10 @@ label collapseBase(triSurface& surf, const scalar minLen)
{
spanPoints[regionI] = getSpanPoints(surf, outsideVerts[regionI]);
- Pout<< "For region " << regionI << " found extrema at points "
- << surf.localPoints()[spanPoints[regionI][0]]
- << surf.localPoints()[spanPoints[regionI][1]]
- << endl;
+ //Pout<< "For region " << regionI << " found extrema at points "
+ // << surf.localPoints()[spanPoints[regionI][0]]
+ // << surf.localPoints()[spanPoints[regionI][1]]
+ // << endl;
// Project all non-span points onto the span edge.
projectNonSpanPoints
@@ -787,21 +831,21 @@ label collapseBase(triSurface& surf, const scalar minLen)
orderedWeights[regionI]
);
- Pout<< "For region:" << regionI
- << " span:" << spanPoints[regionI]
- << " orderedVerts:" << orderedVertices[regionI]
- << " orderedWeights:" << orderedWeights[regionI]
- << endl;
+ //Pout<< "For region:" << regionI
+ // << " span:" << spanPoints[regionI]
+ // << " orderedVerts:" << orderedVertices[regionI]
+ // << " orderedWeights:" << orderedWeights[regionI]
+ // << endl;
- writeRegionOBJ
- (
- surf,
- regionI,
- collapseRegion,
- outsideVerts[regionI]
- );
+ //writeRegionOBJ
+ //(
+ // surf,
+ // regionI,
+ // collapseRegion,
+ // outsideVerts[regionI]
+ //);
- Pout<< endl;
+ //Pout<< endl;
}
@@ -864,20 +908,19 @@ label collapseBase(triSurface& surf, const scalar minLen)
// Split edge using splitVerts. All non-collapsed triangles
// using edge will get split.
-
- {
- const pointField& localPoints = surf.localPoints();
- Pout<< "edge " << edgeI << ' ' << e
- << " points "
- << localPoints[e[0]] << ' ' << localPoints[e[1]]
- << " split into edges with extra points:"
- << endl;
- forAll(splitVerts, i)
- {
- Pout<< " " << splitVerts[i] << " weight "
- << splitWeights[i] << nl;
- }
- }
+ //{
+ // const pointField& localPoints = surf.localPoints();
+ // Pout<< "edge " << edgeI << ' ' << e
+ // << " points "
+ // << localPoints[e[0]] << ' ' << localPoints[e[1]]
+ // << " split into edges with extra points:"
+ // << endl;
+ // forAll(splitVerts, i)
+ // {
+ // Pout<< " " << splitVerts[i] << " weight "
+ // << splitWeights[i] << nl;
+ // }
+ //}
const labelList& eFaces = surf.edgeFaces()[edgeI];
@@ -914,7 +957,8 @@ label collapseBase(triSurface& surf, const scalar minLen)
}
}
- Pout<< "collapseBase : splitting " << nSplit << " triangles"
+ Info<< "collapseBase : collapsing " << nSplit
+ << " triangles by splitting their base edge."
<< endl;
nTotalSplit += nSplit;
@@ -927,15 +971,15 @@ label collapseBase(triSurface& surf, const scalar minLen)
// Pack the triangles
newTris.shrink();
- Pout<< "Resetting surface from " << surf.size() << " to "
- << newTris.size() << " triangles" << endl;
+ //Pout<< "Resetting surface from " << surf.size() << " to "
+ // << newTris.size() << " triangles" << endl;
surf = triSurface(newTris, surf.patches(), surf.localPoints());
- {
- fileName fName("bla" + name(iter) + ".obj");
- Pout<< "Writing surf to " << fName << endl;
- surf.write(fName);
- }
+ //{
+ // fileName fName("bla" + name(iter) + ".obj");
+ // Pout<< "Writing surf to " << fName << endl;
+ // surf.write(fName);
+ //}
iter++;
}
diff --git a/applications/utilities/surface/surfaceClean/collapseBase.H b/applications/utilities/surface/surfaceClean/collapseBase.H
index 5909c4f846..5baf102d4f 100644
--- a/applications/utilities/surface/surfaceClean/collapseBase.H
+++ b/applications/utilities/surface/surfaceClean/collapseBase.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
+ \\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -41,9 +41,14 @@ using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-//- Keep collapsing all triangles whose height is < minLen.
+//- Keep collapsing all triangles whose height is < minLen or quality < minQ.
// Returns number of triangles collapsed.
-label collapseBase(triSurface& surf, const scalar minLen);
+label collapseBase
+(
+ triSurface& surf,
+ const scalar minLen,
+ const scalar minQuality
+);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/applications/utilities/surface/surfaceClean/collapseEdge.C b/applications/utilities/surface/surfaceClean/collapseEdge.C
index b0cd2778f7..addef9c1d5 100644
--- a/applications/utilities/surface/surfaceClean/collapseEdge.C
+++ b/applications/utilities/surface/surfaceClean/collapseEdge.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
+ \\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -121,8 +121,8 @@ label collapseEdge(triSurface& surf, const scalar minLen)
pointMap[v1] = v;
newPoints[v] = 0.5*(localPoints[v1] + localPoints[v]);
- Pout<< "Collapsing triange " << faceI << " to edge mid "
- << newPoints[v] << endl;
+ //Pout<< "Collapsing triange " << faceI
+ // << " to edge mid " << newPoints[v] << endl;
nCollapsed++;
okToCollapse[faceI] = false;
@@ -136,7 +136,8 @@ label collapseEdge(triSurface& surf, const scalar minLen)
}
}
- Pout<< "collapseEdge : collapsing " << nCollapsed << " triangles"
+ Info<< "collapseEdge : collapsing " << nCollapsed
+ << " triangles to a single edge."
<< endl;
nTotalCollapsed += nCollapsed;
diff --git a/applications/utilities/surface/surfaceClean/surfaceClean.C b/applications/utilities/surface/surfaceClean/surfaceClean.C
index 5d25ea4859..c534b13970 100644
--- a/applications/utilities/surface/surfaceClean/surfaceClean.C
+++ b/applications/utilities/surface/surfaceClean/surfaceClean.C
@@ -50,6 +50,7 @@ int main(int argc, char *argv[])
argList::noParallel();
argList::validArgs.append("surfaceFile");
argList::validArgs.append("min length");
+ argList::validArgs.append("min quality");
argList::validArgs.append("output surfaceFile");
argList::addBoolOption
(
@@ -60,10 +61,13 @@ int main(int argc, char *argv[])
const fileName inFileName = args[1];
const scalar minLen = args.argRead(2);
- const fileName outFileName = args[3];
+ const scalar minQuality = args.argRead(3);
+ const fileName outFileName = args[4];
Info<< "Reading surface " << inFileName << nl
- << "Collapsing all triangles with edges or heights < " << minLen << nl
+ << "Collapsing all triangles with" << nl
+ << " edges or heights < " << minLen << nl
+ << " quality < " << minQuality << nl
<< "Writing result to " << outFileName << nl << endl;
@@ -90,7 +94,7 @@ int main(int argc, char *argv[])
}
while (true)
{
- label nSplitEdge = collapseBase(surf, minLen);
+ label nSplitEdge = collapseBase(surf, minLen, minQuality);
if (nSplitEdge == 0)
{
diff --git a/applications/utilities/surface/surfaceFeatureExtract/CGALPolyhedron/buildCGALPolyhedron.C b/applications/utilities/surface/surfaceFeatureExtract/CGALPolyhedron/buildCGALPolyhedron.C
index bc1c46a571..369d712ab1 100644
--- a/applications/utilities/surface/surfaceFeatureExtract/CGALPolyhedron/buildCGALPolyhedron.C
+++ b/applications/utilities/surface/surfaceFeatureExtract/CGALPolyhedron/buildCGALPolyhedron.C
@@ -51,8 +51,8 @@ void Foam::buildCGALPolyhedron::operator()
HalfedgeDS& hds
)
{
- typedef typename HalfedgeDS::Traits Traits;
- typedef typename Traits::Point_3 Point;
+ typedef HalfedgeDS::Traits Traits;
+ typedef Traits::Point_3 Point;
// Postcondition: `hds' is a valid polyhedral surface.
CGAL::Polyhedron_incremental_builder_3 B(hds, false);
diff --git a/differencesWrtDev.txt b/differencesWrtDev.txt
index 3d87e6d01a..f3b0568ca2 100644
--- a/differencesWrtDev.txt
+++ b/differencesWrtDev.txt
@@ -1,7 +1,7 @@
src/
mesh/conformalVoronoiMesh
- - all the meshing.
+ - all the meshing. See separate section below.
polyMeshGeometry:
- disabled tetquality check on face-center decomp tet.
@@ -93,5 +93,43 @@ surfaceSplitByTopology:
+conformalVoronoiMesh/
+--------------------
+- Make sure the surface does not have any sliver triangles. These are
+hard to get the surface normal correct for so might cause bleeding.
+- Use surfaceCheck to find out triangle quality and size of smallest edges
+- Use surfaceClean .. 1e-5 .. to get rid of any edges < 1e-5m.
+
+- If you get bleeding you might see in face filtering:
+
+ ...
+ cells with with zero or one non-boundary face : 1
+ ...
+ Initial check before face collapse, found 48 bad quality faces
+
+but this was real - the cell that got created inside the cone and sphere by the
+bad triangle was actually attached to a lot of faces. This screwed up the
+subsequent filtering as it stopped too early.
+
+I ran:
+
+ surfaceClean coneAndSphere.obj 1e-5 coneAndSphere_clean.obj
+
+and re-ran with that surface and got
+
+ ...
+ cells with with zero or one non-boundary face : 0
+ ...
+ Initial check before face collapse, found 0 bad quality faces
+
+and the bad cells inside are gone.
+
+That group of cells would be picked up at the end of the meshing as the
+
+ cvMesh_remainingProtrusions
+
+cellSet which can be deleted.
+
+
diff --git a/doc/cvMesh.org b/doc/cvMesh.org
new file mode 100644
index 0000000000..ffc168c06e
--- /dev/null
+++ b/doc/cvMesh.org
@@ -0,0 +1,166 @@
+# -*- mode: org; -*-
+#
+#+TITLE: cvMesh usage
+#+AUTHOR: OpenCFD Ltd.
+#+DATE: June 2011
+#+LINK: http://www.OpenFOAM.com
+#+OPTIONS: author:nil ^:{}
+#+STARTUP: hidestars
+#+STARTUP: odd
+
+* Method
+ + initial point seeding. These are the background points to
+ - start the triangulation from (a subset)
+ - do cellsize and alignment queries on
+ + mesh movement based on cell size and alignment
+ + surface conformation by introducing duplicate points
+ + do final simplification by face filtering
+
+
+* Typical work flow
+ + clean surfaces (surfaceClean, surfaceOrient)
+ + set in system/controlDict
+ - startTime to latestTime
+ - endTime to number of meshing iterations
+ - writeInterval to large value
+ + select defaultCellSize [m]
+ + select refinement patterns based on e.g. surface
+ (cellSizeControlGeometry)
+ + select autoDensity initial point seeding
+ + force output of unfiltered mesh (to save time) by specifying unobtainable
+ mesh quality:
+ continueFilteringOnBadInitialPolyMesh false;
+ maxNonOrtho 1;
+ + check visually
+ - features resolved?
+ - cell size ok? Load in targetCellSize volScalarField.
+ - if any surface protusions that are unacceptable
+ - lower targetCellSize
+ - lower maxSurfaceProtrusionCoeff (usually 0.1 of targetCellSize)
+ + restart. Repeat until surface good enough.
+ + produce final mesh
+ - set initialPointsMethod to pointFile so it reuses the Delaunay vertices
+ - re-enable mesh filtering
+ continueFilteringOnBadInitialPolyMesh true;
+ maxNonOrtho 65; //or whatever;
+
+
+* Debugging
+ + run surfaceFeatureExtract with -writeObj to dump all features
+ + set initialPointsMethod to pointFile
+ + set writeInterval to 1
+ + run e.g. topoSet to create slices
+
+
+* Annotated output
+#+BEGIN_EXAMPLE
+Create time
+
+Reading geometryToConformTo
+
+Reading cellSizeControlGeometry
+
+Selecting initialPointsMethod autoDensity
+
+Selecting relaxationModel adaptiveLinear
+
+Selecting faceAreaWeightModel piecewiseLinearRamp
+
+Conforming to feature points
+ Inserted 32 feature vertices
+
+#+END_EXAMPLE
+In this particular case there are 8 straight corners so these can be recreated
+with a triplet each.
+
+#+BEGIN_EXAMPLE
+Inserting initial points
+
+ autoDensity
+ 6922 points placed
+ 9139 locations queried
+ 0.7574132837 success rate
+#+END_EXAMPLE
+The initial points get created by the =autoDensity= method. It will accept only
+part (75% in this cae) of the initial points to create randomly distributed locations.
+#+BEGIN_EXAMPLE
+
+ 6922 points to insert...
+ 6922 points inserted
+
+Build coarse surface conformation
+
+Initial conformation
+ Number of vertices 6962
+ Number of surface hits 2245
+ Number of edge hits 238
+#+END_EXAMPLE
+The initial points get inserted into the triangulation and it detects the number
+of surface pierces (surface hits) and the number of XXX (edge hits).
+#+BEGIN_EXAMPLE
+
+Stopping iterations when:
+ total number of hits drops below 0.001 of initial hits (2)
+ or
+ maximum number of iterations (15) is reached
+#+END_EXAMPLE
+It will introduce additional duplicate vertices on either side of the surface
+pierce to conform to the geometry. This needs to repeated until the number of
+pierces is low enough.
+#+BEGIN_EXAMPLE
+
+Conformation iteration 0
+ Number of vertices 12166
+ Number of surface hits 234
+ Number of edge hits 38
+
+..
+
+Conformation iteration 4
+ Number of vertices 12908
+ Number of surface hits 0
+ Number of edge hits 0
+
+--- [ cpuTime 5.46 s, delta 0.8 s, "Conformation iteration 4" ] ---
+--- [ 5.46 s, mem size 114228 kB, mem peak 114228 kB, mem rss 28712 kB ] ---
+
+Total hits (0) less than limit (2), stopping iterations
+
+Storing surface conformation
+ Stored 5946 vertices
+Store size and alignment
+
+#+END_EXAMPLE
+The vertices needed to create the conformation are stored together with the size
+and alignment.
+#+BEGIN_EXAMPLE
+Time = 1
+
+Relaxation = 1
+
+Looking up target cell alignment and size
+
+Determining vertex displacements
+
+Reinserting stored feature points
+ Reinserted 32 vertices
+
+Inserting displaced tessellation
+ 4328 points to insert...
+ 4328 points inserted
+
+Reinserting stored surface conformation
+ Reinserted 5946 vertices
+
+Total displacement = (0.03629314022 0.1324815449 -0.2774017193)
+Total distance = 159.3482907
+
+ExecutionTime = 8.16 s ClockTime = 8 s
+
+
+Time = 2
+
+Relaxation = 0.9875
+
+#+END_EXAMPLE
+
diff --git a/src/mesh/conformalVoronoiMesh/cellSizeControlSurfaces/cellSizeControlSurfaces.H b/src/mesh/conformalVoronoiMesh/cellSizeControlSurfaces/cellSizeControlSurfaces.H
index f12547e7e2..206e374ab2 100644
--- a/src/mesh/conformalVoronoiMesh/cellSizeControlSurfaces/cellSizeControlSurfaces.H
+++ b/src/mesh/conformalVoronoiMesh/cellSizeControlSurfaces/cellSizeControlSurfaces.H
@@ -37,7 +37,6 @@ SourceFiles
#include "searchableSurfaces.H"
#include "searchableSurfacesQueries.H"
-#include "triSurface.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
index 0be7efa5f9..ce69696d6a 100644
--- a/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
+++ b/src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
@@ -52,7 +52,6 @@ SourceFiles
#include "plane.H"
#include "SortableList.H"
#include "meshTools.H"
-#include "triSurfaceTools.H"
#include "indexedOctree.H"
#include "treeDataPoint.H"
#include "unitConversion.H"
@@ -683,17 +682,18 @@ private:
bool filterFaces
);
- // //- Tet mesh calculation
- // void calcTetMesh
- // (
- // pointField& points,
- // faceList& faces,
- // labelList& owner,
- // labelList& neighbour,
- // wordList& patchNames,
- // labelList& patchSizes,
- // labelList& patchStarts
- // );
+ //- Tet mesh calculation
+ void calcTetMesh
+ (
+ pointField& points,
+ faceList& faces,
+ labelList& owner,
+ labelList& neighbour,
+ wordList& patchTypes,
+ wordList& patchNames,
+ labelList& patchSizes,
+ labelList& patchStarts
+ );
//- Determines if the dual face constructed by the Delaunay
// edge is a boundary face
diff --git a/src/mesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.H b/src/mesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.H
index 29e1a6d9b7..066e9f3228 100644
--- a/src/mesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.H
+++ b/src/mesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.H
@@ -38,7 +38,7 @@ SourceFiles
#include "searchableSurfaces.H"
#include "searchableSurfacesQueries.H"
#include "extendedFeatureEdgeMesh.H"
-#include "triSurface.H"
+#include "boolList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.C b/src/meshTools/searchableSurface/triSurfaceMesh.C
index f91b78df47..4abec4ae3d 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.C
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.C
@@ -331,6 +331,7 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
),
triSurface(s),
tolerance_(indexedOctree::perturbTol()),
+ minQuality_(-1),
maxTreeDepth_(10),
surfaceClosed_(-1)
{
@@ -377,6 +378,7 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io)
)
),
tolerance_(indexedOctree::perturbTol()),
+ minQuality_(-1),
maxTreeDepth_(10),
surfaceClosed_(-1)
{
@@ -426,6 +428,7 @@ Foam::triSurfaceMesh::triSurfaceMesh
)
),
tolerance_(indexedOctree::perturbTol()),
+ minQuality_(-1),
maxTreeDepth_(10),
surfaceClosed_(-1)
{
@@ -449,6 +452,13 @@ Foam::triSurfaceMesh::triSurfaceMesh
<< tolerance_ << endl;
}
+ // Have optional minimum quality for normal calculation
+ if (dict.readIfPresent("minQuality", minQuality_) && minQuality_ > 0)
+ {
+ Info<< searchableSurface::name()
+ << " : ignoring triangles with quality < "
+ << minQuality_ << " for normals calculation." << endl;
+ }
// Have optional non-standard tree-depth to limit storage.
if (dict.readIfPresent("maxTreeDepth", maxTreeDepth_) && maxTreeDepth_ > 0)
@@ -820,22 +830,70 @@ void Foam::triSurfaceMesh::getNormal
{
normal.setSize(info.size());
- forAll(info, i)
+ if (minQuality_ >= 0)
{
- if (info[i].hit())
- {
- label faceI = info[i].index();
- //- Cached:
- //normal[i] = faceNormals()[faceI];
+ // Make sure we don't use triangles with low quality since
+ // normal is not reliable.
- //- Uncached
- normal[i] = triSurface::operator[](faceI).normal(points());
- normal[i] /= mag(normal[i]) + VSMALL;
- }
- else
+ const triSurface& s = static_cast(*this);
+ const labelListList& faceFaces = s.faceFaces();
+
+ forAll(info, i)
{
- // Set to what?
- normal[i] = vector::zero;
+ if (info[i].hit())
+ {
+ label faceI = info[i].index();
+
+ scalar qual = s[faceI].tri(points()).quality();
+
+ if (qual < minQuality_)
+ {
+ // Search neighbouring triangles
+ const labelList& fFaces = faceFaces[faceI];
+
+ forAll(fFaces, j)
+ {
+ label nbrI = fFaces[j];
+ scalar nbrQual = s[nbrI].tri(points()).quality();
+ if (nbrQual > qual)
+ {
+ qual = nbrQual;
+ normal[i] = s[nbrI].normal(points());
+ }
+ }
+ }
+ else
+ {
+ normal[i] = s[faceI].normal(points());
+ }
+ normal[i] /= mag(normal[i]);
+ }
+ else
+ {
+ // Set to what?
+ normal[i] = vector::zero;
+ }
+ }
+ }
+ else
+ {
+ forAll(info, i)
+ {
+ if (info[i].hit())
+ {
+ label faceI = info[i].index();
+ //- Cached:
+ //normal[i] = faceNormals()[faceI];
+
+ //- Uncached
+ normal[i] = triSurface::operator[](faceI).normal(points());
+ normal[i] /= mag(normal[i]) + VSMALL;
+ }
+ else
+ {
+ // Set to what?
+ normal[i] = vector::zero;
+ }
}
}
}
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.H b/src/meshTools/searchableSurface/triSurfaceMesh.H
index 911d12cf73..2b4acaf153 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.H
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
+ \\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -31,6 +31,7 @@ Description
- scale : scaling factor.
- tolerance : relative tolerance for doing intersections
(see triangle::intersection)
+ - minQuality: discard triangles with low quality when getting normal
SourceFiles
triSurfaceMesh.C
@@ -70,6 +71,10 @@ private:
//- Optional tolerance to use in searches
scalar tolerance_;
+ //- Optional min triangle quality. Triangles below this get
+ // ignored for normal calculation
+ scalar minQuality_;
+
//- Optional max tree depth of octree
label maxTreeDepth_;