diff --git a/.gitignore b/.gitignore
index 96c3093b718b9a8ee6d711583018ae836555158a..150288aae2555094adff1f866bb0e4e950f07c9f 100644
--- a/.gitignore
+++ b/.gitignore
@@ -55,7 +55,7 @@ doc/[Dd]oxygen/man
# ignore .timeStamp in the main directory
/.timeStamp
-# ignore .ebrowse in the main directory
-/.ebrowse
+# ignore .tags in the main directory
+/.tags
# end-of-file
diff --git a/applications/solvers/incompressible/simpleFoam/Make/options b/applications/solvers/incompressible/simpleFoam/Make/options
index 61091811bf04fd0547f08d2b4cd5c3c9385d4d6d..1223bdd06f48178f0b2eee8032d29c936634f328 100644
--- a/applications/solvers/incompressible/simpleFoam/Make/options
+++ b/applications/solvers/incompressible/simpleFoam/Make/options
@@ -8,5 +8,4 @@ EXE_INC = \
EXE_LIBS = \
-lincompressibleRASModels \
-lincompressibleTransportModels \
- -lfiniteVolume \
- -lmeshTools
+ -lfiniteVolume
diff --git a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict
index 43d550f4068b72c6967430b49c0c63ac1516bc99..58f324160e8c66aa5e253cd1a725ff0351b0493d 100644
--- a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict
+++ b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict
@@ -316,6 +316,11 @@ meshQualityControls
//must be >0 for Fluent compatibility
minTriangleTwist -1;
+ //- if >0 : preserve single cells with all points on the surface if the
+ // resulting volume after snapping is larger than minVolFraction times old
+ // volume. If <0 : delete always.
+ minVolFraction 0.1;
+
// Advanced
diff --git a/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C b/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C
index 52f23658a84b947d3e9fbeaf2a51a0bf85bcebff..0adcd309275da4ae41ebb2274a2a29d68e23cb4a 100644
--- a/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C
+++ b/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C
@@ -23,16 +23,13 @@ Foam::label Foam::checkGeometry(const polyMesh& mesh, const bool allGeometry)
scalar minDistSqr = magSqr(1e-6 * globalBb.span());
// Non-empty directions
- const Vector<label> validDirs = (mesh.directions() + Vector<label>::one)/2;
+ const Vector<label> validDirs = (mesh.geometricD() + Vector<label>::one)/2;
+ Info<< " Mesh (non-empty, non-wedge) directions " << validDirs << endl;
- Info<< " Mesh (non-empty) directions " << validDirs << endl;
+ const Vector<label> solDirs = (mesh.solutionD() + Vector<label>::one)/2;
+ Info<< " Mesh (non-empty) directions " << solDirs << endl;
- scalar nGeomDims = mesh.nGeometricD();
-
- Info<< " Mesh (non-empty, non-wedge) dimensions "
- << nGeomDims << endl;
-
- if (nGeomDims < 3)
+ if (mesh.nGeometricD() < 3)
{
pointSet nonAlignedPoints(mesh, "nonAlignedEdges", mesh.nPoints()/100);
diff --git a/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/PV3FoamReader_SM.xml b/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/PV3FoamReader_SM.xml
index 644b794a44637c0e3af84fd2bc62f07252dd28d0..ef25ac03786e43e5dd407cd7ef992d6b1c3b91b4 100644
--- a/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/PV3FoamReader_SM.xml
+++ b/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/PV3FoamReader_SM.xml
@@ -30,16 +30,16 @@
<!-- Global settings -->
- <!-- Extrapolate Walls check-box -->
+ <!-- Extrapolate Patches check-box -->
<IntVectorProperty
- name="ExtrapolateWalls"
- command="SetExtrapolateWalls"
+ name="ExtrapolatePatches"
+ command="SetExtrapolatePatches"
number_of_elements="1"
default_values="0"
animateable="0">
<BooleanDomain name="bool"/>
<Documentation>
- Extrapolate internalField to wall and empty patches
+ Extrapolate internalField to non-constraint patches
</Documentation>
</IntVectorProperty>
diff --git a/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.cxx b/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.cxx
index 8fc518d745edaa19bd0107bcc391c94d6c9bbbb8..c16a43deaef9c5a619a2fb3c06fcd01a3cbce4c5 100644
--- a/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.cxx
+++ b/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.cxx
@@ -64,7 +64,7 @@ vtkPV3FoamReader::vtkPV3FoamReader()
CacheMesh = 1;
- ExtrapolateWalls = 0;
+ ExtrapolatePatches = 0;
IncludeSets = 0;
IncludeZones = 0;
ShowPatchNames = 0;
diff --git a/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.h b/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.h
index fde87527de29c134526ee73c72909d1bd0fed8d9..bc21cd8ce816978129720a4822454c455897780d 100644
--- a/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.h
+++ b/applications/utilities/postProcessing/graphics/PV3FoamReader/PV3FoamReader/vtkPV3FoamReader.h
@@ -65,9 +65,9 @@ public:
vtkGetMacro(CacheMesh, int);
// Description:
- // FOAM extrapolate internal values onto the walls
- vtkSetMacro(ExtrapolateWalls, int);
- vtkGetMacro(ExtrapolateWalls, int);
+ // FOAM extrapolate internal values onto the patches
+ vtkSetMacro(ExtrapolatePatches, int);
+ vtkGetMacro(ExtrapolatePatches, int);
// FOAM read sets control
vtkSetMacro(IncludeSets, int);
@@ -183,7 +183,7 @@ private:
int TimeStepRange[2];
int CacheMesh;
- int ExtrapolateWalls;
+ int ExtrapolatePatches;
int IncludeSets;
int IncludeZones;
int ShowPatchNames;
diff --git a/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.C b/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.C
index 7a8fd89e9ed1529b86acfcf3c82f67a11dd474fd..0ed3d53185000a439edbb5654bc2a300fbf1f3a0 100644
--- a/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.C
+++ b/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3Foam.C
@@ -659,29 +659,55 @@ void Foam::vtkPV3Foam::addPatchNames(vtkRenderer* renderer)
}
}
+ // Count number of zones we're actually going to display. This is truncated
+ // to a max per patch
+
+ const label MAXPATCHZONES = 20;
+
+ label displayZoneI = 0;
+
+ forAll(pbMesh, patchI)
+ {
+ displayZoneI += min(MAXPATCHZONES, nZones[patchI]);
+ }
+
+
zoneCentre.shrink();
if (debug)
{
Info<< "patch zone centres = " << zoneCentre << nl
+ << "displayed zone centres = " << displayZoneI << nl
<< "zones per patch = " << nZones << endl;
}
// Set the size of the patch labels to max number of zones
- patchTextActorsPtrs_.setSize(zoneCentre.size());
+ patchTextActorsPtrs_.setSize(displayZoneI);
if (debug)
{
Info<< "constructing patch labels" << endl;
}
+ // Actor index
+ displayZoneI = 0;
+
+ // Index in zone centres
label globalZoneI = 0;
+
forAll(pbMesh, patchI)
{
const polyPatch& pp = pbMesh[patchI];
// Only selected patches will have a non-zero number of zones
- for (label i=0; i<nZones[patchI]; i++)
+ label nDisplayZones = min(MAXPATCHZONES, nZones[patchI]);
+ label increment = 1;
+ if (nZones[patchI] >= MAXPATCHZONES)
+ {
+ increment = nZones[patchI]/MAXPATCHZONES;
+ }
+
+ for (label i = 0; i < nDisplayZones; i++)
{
if (debug)
{
@@ -719,14 +745,15 @@ void Foam::vtkPV3Foam::addPatchNames(vtkRenderer* renderer)
// Maintain a list of text labels added so that they can be
// removed later
- patchTextActorsPtrs_[globalZoneI] = txt;
+ patchTextActorsPtrs_[displayZoneI] = txt;
- globalZoneI++;
+ globalZoneI += increment;
+ displayZoneI++;
}
}
// Resize the patch names list to the actual number of patch names added
- patchTextActorsPtrs_.setSize(globalZoneI);
+ patchTextActorsPtrs_.setSize(displayZoneI);
if (debug)
{
diff --git a/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3FoamVolFields.H b/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3FoamVolFields.H
index f217d1a8890afa67b7f71d3856d2290282e07753..6c8d32f4702c54810b59c1e63f867513a62365f9 100644
--- a/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3FoamVolFields.H
+++ b/applications/utilities/postProcessing/graphics/PV3FoamReader/vtkPV3Foam/vtkPV3FoamVolFields.H
@@ -132,8 +132,8 @@ void Foam::vtkPV3Foam::convertVolFields
isType<emptyFvPatchField<Type> >(ptf)
||
(
- typeid(patches[patchId]) == typeid(wallPolyPatch)
- && reader_->GetExtrapolateWalls()
+ reader_->GetExtrapolatePatches()
+ && !polyPatch::constraintType(patches[patchId].type())
)
)
{
diff --git a/bin/foamEbrowse b/bin/foamEbrowse
index e303a772643c8b9f74fa5bc60bc565e110a90a20..66fafc1e2b1b46f5384ca1c78f309fb1b785559d 100755
--- a/bin/foamEbrowse
+++ b/bin/foamEbrowse
@@ -27,26 +27,27 @@
# foamEbrowse
#
# Description
-# Build the Ebrowse database for all the .C and .H files
+# Build the Ebrowse database for all the .H and .C files
#
#------------------------------------------------------------------------------
-headersFile=${TMPDIR:-/tmp}/headersFile.$$
sourcesFile=${TMPDIR:-/tmp}/sourcesFile.$$
if [ $# -ne 0 ]; then
echo "Usage : ${0##*/}"
echo ""
- echo "Build the Ebrowse dadbase for all the .C and .H files"
+ echo "Build the Ebrowse dadbase for all the .H and .C files"
echo ""
exit 1
fi
# Clean up on termination and on Ctrl-C
-trap 'rm -f $headersFile $sourcesFile 2>/dev/null; exit 0' EXIT TERM INT
+trap 'rm -f $sourcesFile 2>/dev/null; exit 0' EXIT TERM INT
cd $WM_PROJECT_DIR
-find -H . -name "*.H" | fgrep -v lnInclude > $headersFile
-find -H . -name "*.C" | fgrep -v lnInclude > $sourcesFile
-ebrowse --files=$headersFile --files=$sourcesFile --output-file=.ebrowse
+mkdir .tags 2>/dev/null
+cd .tags
+
+find -H .. \( -name "*.[HC]" -not -name "lnInclude" -not -name "Doxygen" \) -print > $sourcesFile
+ebrowse --files=$sourcesFile --output-file=ebrowse
#------------------------------------------------------------------------------
diff --git a/bin/foamTags b/bin/foamTags
new file mode 100755
index 0000000000000000000000000000000000000000..1cfd186f52b1b32e1d9d8c77fce0d44e5218d507
--- /dev/null
+++ b/bin/foamTags
@@ -0,0 +1,56 @@
+#!/bin/sh
+#------------------------------------------------------------------------------
+# ========= |
+# \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+# \\ / O peration |
+# \\ / A nd | Copyright (C) 1991-2009 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 2 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, write to the Free Software Foundation,
+# Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+#
+# Script
+# foamTags
+#
+# Description
+# Build the tags files for all the .C and .H files
+#
+#------------------------------------------------------------------------------
+
+if [ $# -ne 0 ]; then
+ echo "Usage : ${0##*/}"
+ echo ""
+ echo "Build the tags files for all the .C and .H files"
+ echo ""
+ exit 1
+fi
+
+cd $WM_PROJECT_DIR
+mkdir .tags 2>/dev/null
+
+find -H . \( -name "*.[HC]" -not -name "lnInclude" -not -name "Doxygen" \) | \
+ etags --declarations -l c++ -o .tags/etags -
+find -H . \( -name "*.[HC]" -not -name "lnInclude" -not -name "Doxygen" \) | \
+ etags -l c++ -o .tags/etagsDef -
+find -H . \( -name "*.H" -not -name "lnInclude" -not -name "Doxygen" \) | \
+ etags --declarations -l c++ -o .tags/etagsDec -
+
+gtags -i --gtagsconf bin/tools/gtagsrc .tags
+
+foamEbrowse
+
+#------------------------------------------------------------------------------
diff --git a/bin/tools/gtagsrc b/bin/tools/gtagsrc
new file mode 100644
index 0000000000000000000000000000000000000000..ce66e0c084a8188bf63b707880344ccc18813143
--- /dev/null
+++ b/bin/tools/gtagsrc
@@ -0,0 +1,63 @@
+#------------------------------------------------------------------------------
+# ========= |
+# \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+# \\ / O peration |
+# \\ / A nd | Copyright (C) 1991-2009 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 2 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, write to the Free Software Foundation,
+# Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+#
+# Canfiguration file
+# gtagsrc
+#
+# Description
+# Configuration file for gtags(1).
+#
+#------------------------------------------------------------------------------
+default:\
+ :tc=gtags:tc=htags:
+#------------------------------------------------------------------------------
+# Configuration for gtags(1)
+# See gtags(1).
+#------------------------------------------------------------------------------
+common:\
+ :skip=GPATH,GTAGS,GRTAGS,GSYMS,HTML/,HTML.pub/,html/,tags,TAGS,ID,.ebrowse,.etags,.etagsDef,.etagsDec,y.tab.c,y.tab.h,.notfunction,cscope.out,cscope.po.out,cscope.in.out,.gdbinit,SCCS/,RCS/,CVS/,CVSROOT/,{arch}/,.svn/,.git/,.cvsrc,.cvsignore,.gitignore,.cvspass,.cvswrappers,.deps/,autom4te.cache/,.snprj/:\
+ :langmap=c\:.c.h,yacc\:.y,asm\:.s.S,java\:.java,cpp\:.c++.cc.cpp.cxx.hxx.hpp.C.H,php\:.php.php3.phtml:
+gtags:\
+ :tc=common:\
+ :GTAGS=gtags-parser %s:\
+ :GRTAGS=gtags-parser -r %s:\
+ :GSYMS=gtags-parser -s %s:\
+ :skip=lnInclude/,tutorials/,wmake/,doc/,lib/,etc/:
+#
+#------------------------------------------------------------------------------
+# Configuration for htags(1)
+#------------------------------------------------------------------------------
+htags:\
+ :body_begin=<body text='#191970' bgcolor='#f5f5dc' vlink='gray'>:body_end=</body>:\
+ :table_begin=<table>:table_end=</table>:\
+ :title_begin=<h1><font color='#cc0000'>:title_end=</font></h1>:\
+ :comment_begin=<i><font color='green'>:comment_end=</font></i>:\
+ :sharp_begin=<font color='darkred'>:sharp_end=</font>:\
+ :brace_begin=<font color='red'>:brace_end=</font>:\
+ :warned_line_begin=<span style='background-color\:yellow'>:warned_line_end=</span>:\
+ :reserved_begin=<b>:reserved_end=</b>:script_alias=/cgi-bin/:\
+ :ncol#4:tabs#8:normal_suffix=html:gzipped_suffix=ghtml:\
+ :definition_header=no:
+
+#------------------------------------------------------------------------------
\ No newline at end of file
diff --git a/src/Allwmake b/src/Allwmake
index f22e0729dbe8e50694ee43a373205494b0a2a536..c3ce648ba0aa1462b53dff9a1502db02d5a215aa 100755
--- a/src/Allwmake
+++ b/src/Allwmake
@@ -5,8 +5,8 @@ set -x
# update OpenFOAM version strings if required
wmakePrintBuild -check || /bin/rm -f OpenFOAM/Make/$WM_OPTIONS/global.? 2>/dev/null
-wmakeLnInclude -f OpenFOAM
-wmakeLnInclude -f OSspecific/$WM_OS
+wmakeLnInclude -f OpenFOAM -sf
+wmakeLnInclude -f OSspecific/$WM_OS -sf
Pstream/Allwmake
wmake libo OSspecific/$WM_OS
diff --git a/src/OpenFOAM/db/Time/findInstance.C b/src/OpenFOAM/db/Time/findInstance.C
index f21a524c80919f0a395e36b479c8a30390f9ce75..a19d3bda633975b579013deb8e6149a5006506de 100644
--- a/src/OpenFOAM/db/Time/findInstance.C
+++ b/src/OpenFOAM/db/Time/findInstance.C
@@ -58,8 +58,9 @@ Foam::word Foam::Time::findInstance
{
if (debug)
{
- Info<< "Time::findInstance(const fileName&, const word&) : "
- << "found \"" << name
+ Info<< "Time::findInstance"
+ "(const fileName&, const word&, const IOobject::readOption)"
+ << " : found \"" << name
<< "\" in " << timeName()/dir
<< endl;
}
@@ -98,8 +99,8 @@ Foam::word Foam::Time::findInstance
if (debug)
{
Info<< "Time::findInstance"
- "(const fileName&,const word&) : "
- << "found \"" << name
+ "(const fileName&, const word&, const IOobject::readOption)"
+ << " : found \"" << name
<< "\" in " << ts[instanceI].name()/dir
<< endl;
}
@@ -129,8 +130,8 @@ Foam::word Foam::Time::findInstance
if (debug)
{
Info<< "Time::findInstance"
- "(const fileName&,const word&) : "
- << "found \"" << name
+ "(const fileName&, const word&, const IOobject::readOption)"
+ << " : found \"" << name
<< "\" in " << constant()/dir
<< endl;
}
@@ -141,10 +142,10 @@ Foam::word Foam::Time::findInstance
if (rOpt == IOobject::MUST_READ)
{
FatalErrorIn
- (
- "Time::findInstance(const fileName&,const word&)"
- )
- << "Cannot find file \"" << name << "\" in directory "
+ (
+ "Time::findInstance"
+ "(const fileName&, const word&, const IOobject::readOption)"
+ ) << "Cannot find file \"" << name << "\" in directory "
<< constant()/dir
<< exit(FatalError);
}
diff --git a/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.C b/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.C
index 70f738fcfbd41f9b4d369b19d33e3442d420a312..fa6a2a58d0a2874f56b6ed7db2e7bc15af060f4c 100644
--- a/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.C
+++ b/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.C
@@ -154,7 +154,7 @@ tmp<GeometricField<scalar, PatchField, GeoMesh> > pow
pow
(
gsf1.dimensions(),
- dimensionedScalar("1", 1.0, gsf2.dimensions())
+ dimensionedScalar("1", gsf2.dimensions(), 1.0)
)
)
);
@@ -183,7 +183,7 @@ tmp<GeometricField<scalar, PatchField, GeoMesh> > pow
pow
(
gsf1.dimensions(),
- dimensionedScalar("1", 1.0, gsf2.dimensions())
+ dimensionedScalar("1", gsf2.dimensions(), 1.0)
)
)
);
@@ -214,7 +214,7 @@ tmp<GeometricField<scalar, PatchField, GeoMesh> > pow
pow
(
gsf1.dimensions(),
- dimensionedScalar("1", 1.0, gsf2.dimensions())
+ dimensionedScalar("1", gsf2.dimensions(), 1.0)
)
)
);
@@ -247,7 +247,7 @@ tmp<GeometricField<scalar, PatchField, GeoMesh> > pow
pow
(
gsf1.dimensions(),
- dimensionedScalar("1", 1.0, gsf2.dimensions())
+ dimensionedScalar("1", gsf2.dimensions(), 1.0)
)
)
);
diff --git a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixATmul.C b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixATmul.C
index 666d84a5c1206588bce906d4174fe31a8857cda1..05958edf68ad060e63c51e8f078656fbbb82a5e6 100644
--- a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixATmul.C
+++ b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixATmul.C
@@ -60,7 +60,7 @@ void Foam::lduMatrix::Amul
interfaceBouCoeffs,
interfaces,
psi,
- Apsi,
+ Apsi,
cmpt
);
@@ -79,7 +79,7 @@ void Foam::lduMatrix::Amul
register const label nFaces = upper().size();
#ifdef ICC_IA64_PREFETCH
- #pragma swp
+ #pragma swp
#endif
for (register label face=0; face<nFaces; face++)
{
@@ -227,7 +227,7 @@ void Foam::lduMatrix::sumA
}
#ifdef ICC_IA64_PREFETCH
- #pragma swp
+ #pragma swp
#endif
for (register label face=0; face<nFaces; face++)
@@ -316,7 +316,7 @@ void Foam::lduMatrix::residual
mBouCoeffs,
interfaces,
psi,
- rA,
+ rA,
cmpt
);
@@ -336,7 +336,7 @@ void Foam::lduMatrix::residual
register const label nFaces = upper().size();
#ifdef ICC_IA64_PREFETCH
- #pragma swp
+ #pragma swp
#endif
for (register label face=0; face<nFaces; face++)
{
diff --git a/src/OpenFOAM/meshes/polyMesh/polyMesh.C b/src/OpenFOAM/meshes/polyMesh/polyMesh.C
index 0428b252d7d173713c13b79510a4d5244e64dd39..e0b4c0536c4ff4623475789c7db5589bb9beff7a 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMesh.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyMesh.C
@@ -54,40 +54,79 @@ void Foam::polyMesh::calcDirections() const
{
for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
{
- directions_[cmpt] = 1;
+ solutionD_[cmpt] = 1;
}
+ // Knock out empty and wedge directions. Note:they will be present on all
+ // domains.
+
label nEmptyPatches = 0;
+ label nWedgePatches = 0;
- vector dirVec = vector::zero;
+ vector emptyDirVec = vector::zero;
+ vector wedgeDirVec = vector::zero;
forAll(boundaryMesh(), patchi)
{
- if (isA<emptyPolyPatch>(boundaryMesh()[patchi]))
+ if (boundaryMesh()[patchi].size())
{
- if (boundaryMesh()[patchi].size())
+ if (isA<emptyPolyPatch>(boundaryMesh()[patchi]))
{
nEmptyPatches++;
- dirVec += sum(cmptMag(boundaryMesh()[patchi].faceAreas()));
+ emptyDirVec += sum(cmptMag(boundaryMesh()[patchi].faceAreas()));
+ }
+ else if (isA<wedgePolyPatch>(boundaryMesh()[patchi]))
+ {
+ const wedgePolyPatch& wpp = refCast<const wedgePolyPatch>
+ (
+ boundaryMesh()[patchi]
+ );
+
+ nWedgePatches++;
+ wedgeDirVec += cmptMag(wpp.centreNormal());
}
}
}
if (nEmptyPatches)
{
- reduce(dirVec, sumOp<vector>());
+ reduce(emptyDirVec, sumOp<vector>());
- dirVec /= mag(dirVec);
+ emptyDirVec /= mag(emptyDirVec);
for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
{
- if (dirVec[cmpt] > 1e-6)
+ if (emptyDirVec[cmpt] > 1e-6)
{
- directions_[cmpt] = -1;
+ solutionD_[cmpt] = -1;
}
else
{
- directions_[cmpt] = 1;
+ solutionD_[cmpt] = 1;
+ }
+ }
+ }
+
+
+ // Knock out wedge directions
+
+ geometricD_ = solutionD_;
+
+ if (nWedgePatches)
+ {
+ reduce(wedgeDirVec, sumOp<vector>());
+
+ wedgeDirVec /= mag(wedgeDirVec);
+
+ for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+ {
+ if (wedgeDirVec[cmpt] > 1e-6)
+ {
+ geometricD_[cmpt] = -1;
+ }
+ else
+ {
+ geometricD_[cmpt] = 1;
}
}
}
@@ -163,7 +202,8 @@ Foam::polyMesh::polyMesh(const IOobject& io)
*this
),
bounds_(points_),
- directions_(Vector<label>::zero),
+ geometricD_(Vector<label>::zero),
+ solutionD_(Vector<label>::zero),
pointZones_
(
IOobject
@@ -350,7 +390,8 @@ Foam::polyMesh::polyMesh
0
),
bounds_(points_, syncPar),
- directions_(Vector<label>::zero),
+ geometricD_(Vector<label>::zero),
+ solutionD_(Vector<label>::zero),
pointZones_
(
IOobject
@@ -505,7 +546,8 @@ Foam::polyMesh::polyMesh
0
),
bounds_(points_, syncPar),
- directions_(Vector<label>::zero),
+ geometricD_(Vector<label>::zero),
+ solutionD_(Vector<label>::zero),
pointZones_
(
IOobject
@@ -766,44 +808,37 @@ const Foam::fileName& Foam::polyMesh::facesInstance() const
}
-const Foam::Vector<Foam::label>& Foam::polyMesh::directions() const
+const Foam::Vector<Foam::label>& Foam::polyMesh::geometricD() const
{
- if (directions_.x() == 0)
+ if (geometricD_.x() == 0)
{
calcDirections();
}
- return directions_;
+ return geometricD_;
}
Foam::label Foam::polyMesh::nGeometricD() const
{
- label nWedges = 0;
+ return cmptSum(geometricD() + Vector<label>::one)/2;
+}
- forAll(boundary_, patchi)
- {
- if (isA<wedgePolyPatch>(boundary_[patchi]))
- {
- nWedges++;
- }
- }
- if (nWedges != 0 && nWedges != 2 && nWedges != 4)
+const Foam::Vector<Foam::label>& Foam::polyMesh::solutionD() const
+{
+ if (solutionD_.x() == 0)
{
- FatalErrorIn("label polyMesh::nGeometricD() const")
- << "Number of wedge patches " << nWedges << " is incorrect, "
- "should be 0, 2 or 4"
- << exit(FatalError);
+ calcDirections();
}
- return nSolutionD() - nWedges/2;
+ return solutionD_;
}
Foam::label Foam::polyMesh::nSolutionD() const
{
- return cmptSum(directions() + Vector<label>::one)/2;
+ return cmptSum(solutionD() + Vector<label>::one)/2;
}
@@ -823,6 +858,10 @@ void Foam::polyMesh::addPatches
<< abort(FatalError);
}
+ // Reset valid directions
+ geometricD_ = Vector<label>::zero;
+ solutionD_ = Vector<label>::zero;
+
boundary_.setSize(p.size());
// Copy the patch pointers
@@ -1037,6 +1076,10 @@ Foam::tmp<Foam::scalarField> Foam::polyMesh::movePoints
faceZones_.movePoints(points_);
cellZones_.movePoints(points_);
+ // Reset valid directions (could change with rotation)
+ geometricD_ = Vector<label>::zero;
+ solutionD_ = Vector<label>::zero;
+
// Hack until proper callbacks. Below are all the polyMeh MeshObjects with a
// movePoints function.
diff --git a/src/OpenFOAM/meshes/polyMesh/polyMesh.H b/src/OpenFOAM/meshes/polyMesh/polyMesh.H
index 588394d1c230e771d5155e537dd81673befd4925..ecb805ad7f0d7912a916bbd302cf9d96da219f77 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMesh.H
+++ b/src/OpenFOAM/meshes/polyMesh/polyMesh.H
@@ -120,9 +120,13 @@ private:
// Created from points on construction, updated when the mesh moves
boundBox bounds_;
+ //- vector of non-constrained directions in mesh
+ // defined according to the presence of empty and wedge patches
+ mutable Vector<label> geometricD_;
+
//- vector of valid directions in mesh
// defined according to the presence of empty patches
- mutable Vector<label> directions_;
+ mutable Vector<label> solutionD_;
// Zoning information
@@ -309,17 +313,22 @@ public:
return bounds_;
}
- //- Return the vector of valid directions in mesh.
- // Defined according to the presence of empty patches.
- // 1 indicates valid direction and -1 an invalid direction.
- const Vector<label>& directions() const;
+ //- Return the vector of geometric directions in mesh.
+ // Defined according to the presence of empty and wedge patches.
+ // 1 indicates unconstrained direction and -1 a constrained
+ // direction.
+ const Vector<label>& geometricD() const;
//- Return the number of valid geometric dimensions in the mesh
label nGeometricD() const;
- //- Return the number of valid solution dimensions in the mesh.
- // For wedge cases this includes the circumferential direction
- // in case of swirl.
+ //- Return the vector of solved-for directions in mesh.
+ // Differs from geometricD in that it includes for wedge cases
+ // the circumferential direction in case of swirl.
+ // 1 indicates valid direction and -1 an invalid direction.
+ const Vector<label>& solutionD() const;
+
+ //- Return the number of valid solved-for dimensions in the mesh
label nSolutionD() const;
//- Return point zone mesh
diff --git a/src/OpenFOAM/meshes/polyMesh/polyMeshClear.C b/src/OpenFOAM/meshes/polyMesh/polyMeshClear.C
index 52f44464444b3643044c37fb6ba14563989121da..71831209849efb0106a2adde6d0d1a8522670ecc 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMeshClear.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyMeshClear.C
@@ -68,6 +68,10 @@ void Foam::polyMesh::clearGeom()
boundary_[patchI].clearGeom();
}
+ // Reset valid directions (could change with rotation)
+ geometricD_ = Vector<label>::zero;
+ solutionD_ = Vector<label>::zero;
+
pointMesh::Delete(*this);
}
@@ -87,6 +91,10 @@ void Foam::polyMesh::clearAddressing()
// recalculation
deleteDemandDrivenData(globalMeshDataPtr_);
+ // Reset valid directions
+ geometricD_ = Vector<label>::zero;
+ solutionD_ = Vector<label>::zero;
+
pointMesh::Delete(*this);
}
diff --git a/src/OpenFOAM/meshes/polyMesh/polyMeshFromShapeMesh.C b/src/OpenFOAM/meshes/polyMesh/polyMeshFromShapeMesh.C
index d12dffff62be220a1016e6ac59eedafeba6018ba..f1cc51cb9ec05fd5f3261584f6b92391eda8316b 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMeshFromShapeMesh.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyMeshFromShapeMesh.C
@@ -217,7 +217,8 @@ Foam::polyMesh::polyMesh
boundaryFaces.size() + 1 // add room for a default patch
),
bounds_(points_, syncPar),
- directions_(Vector<label>::zero),
+ geometricD_(Vector<label>::zero),
+ solutionD_(Vector<label>::zero),
pointZones_
(
IOobject
diff --git a/src/OpenFOAM/meshes/polyMesh/polyMeshUpdate.C b/src/OpenFOAM/meshes/polyMesh/polyMeshUpdate.C
index 24d9ebbb136d1f7d751d2f30549d43d34896c533..7dcd3ba6b40eb1fdbf088af60b8a90c5e3f24991 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMeshUpdate.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyMeshUpdate.C
@@ -68,6 +68,11 @@ void Foam::polyMesh::updateMesh(const mapPolyMesh& mpm)
newMotionPoints.map(oldMotionPoints, mpm.pointMap());
}
+ // Reset valid directions (could change by faces put into empty patches)
+ geometricD_ = Vector<label>::zero;
+ solutionD_ = Vector<label>::zero;
+
+
// Hack until proper callbacks. Below are all the polyMesh-MeshObjects.
// pointMesh
diff --git a/src/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C b/src/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C
index ecfff62ac3f9b3bb43cff0b9021f31a71780d08d..ddf162dd0f2910d2aaaf89a0f2465a26b86926b7 100644
--- a/src/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C
+++ b/src/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C
@@ -29,6 +29,7 @@ License
#include "mathematicalConstants.H"
#include "refinementSurfaces.H"
#include "searchableSurfaces.H"
+#include "regExp.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -300,10 +301,44 @@ Foam::layerParameters::layerParameters
// readScalar(layerDict.lookup("minThickness"));
}
}
-}
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+ // Check whether layer specification matches any patches
+ const List<keyType> wildCards = layersDict.keys(true);
+
+ forAll(wildCards, i)
+ {
+ regExp re(wildCards[i]);
+
+ bool hasMatch = false;
+ forAll(boundaryMesh, patchI)
+ {
+ if (re.match(boundaryMesh[patchI].name()))
+ {
+ hasMatch = true;
+ break;
+ }
+ }
+ if (!hasMatch)
+ {
+ IOWarningIn("layerParameters::layerParameters(..)", layersDict)
+ << "Wildcard layer specification for " << wildCards[i]
+ << " does not match any patch." << endl;
+ }
+ }
+
+ const List<keyType> nonWildCards = layersDict.keys(false);
+
+ forAll(nonWildCards, i)
+ {
+ if (boundaryMesh.findPatchID(nonWildCards[i]) == -1)
+ {
+ IOWarningIn("layerParameters::layerParameters(..)", layersDict)
+ << "Layer specification for " << nonWildCards[i]
+ << " does not match any patch." << endl;
+ }
+ }
+}
// ************************************************************************* //
diff --git a/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C b/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
index a626103cd1af105355aee538164e09ed9250a102..6bbb1c8d3987ebf35cc723c407340a3ad86f2962 100644
--- a/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
+++ b/src/autoMesh/autoHexMesh/refinementSurfaces/refinementSurfaces.C
@@ -265,101 +265,112 @@ Foam::refinementSurfaces::refinementSurfaces
zoneInside_(surfacesDict.size()),
regionOffset_(surfacesDict.size())
{
- labelList globalMinLevel(surfacesDict.size(), 0);
- labelList globalMaxLevel(surfacesDict.size(), 0);
- scalarField globalAngle(surfacesDict.size(), -GREAT);
- List<Map<label> > regionMinLevel(surfacesDict.size());
- List<Map<label> > regionMaxLevel(surfacesDict.size());
- List<Map<scalar> > regionAngle(surfacesDict.size());
+ // Wilcard specification : loop over all surface, all regions
+ // and try to find a match.
+ // Count number of surfaces.
label surfI = 0;
- forAllConstIter(dictionary, surfacesDict, iter)
+ forAll(allGeometry.names(), geomI)
{
- names_[surfI] = iter().keyword();
+ const word& geomName = allGeometry_.names()[geomI];
- surfaces_[surfI] = allGeometry_.findSurfaceID(names_[surfI]);
-
- if (surfaces_[surfI] == -1)
+ if (surfacesDict.found(geomName))
{
- FatalErrorIn
- (
- "refinementSurfaces::refinementSurfaces"
- "(const searchableSurfaces&, const dictionary>&"
- ) << "No surface called " << iter().keyword() << endl
- << "Valid surfaces are " << allGeometry_.names()
- << exit(FatalError);
+ surfI++;
}
- const dictionary& dict = surfacesDict.subDict(iter().keyword());
+ }
- const labelPair refLevel(dict.lookup("level"));
- globalMinLevel[surfI] = refLevel[0];
- globalMaxLevel[surfI] = refLevel[1];
+ // Size lists
+ surfaces_.setSize(surfI);
+ names_.setSize(surfI);
+ faceZoneNames_.setSize(surfI);
+ cellZoneNames_.setSize(surfI);
+ zoneInside_.setSize(surfI);
+ regionOffset_.setSize(surfI);
+
+ labelList globalMinLevel(surfI, 0);
+ labelList globalMaxLevel(surfI, 0);
+ scalarField globalAngle(surfI, -GREAT);
+ List<Map<label> > regionMinLevel(surfI);
+ List<Map<label> > regionMaxLevel(surfI);
+ List<Map<scalar> > regionAngle(surfI);
+
+ surfI = 0;
+ forAll(allGeometry.names(), geomI)
+ {
+ const word& geomName = allGeometry_.names()[geomI];
- // Global zone names per surface
- if (dict.found("faceZone"))
+ if (surfacesDict.found(geomName))
{
- dict.lookup("faceZone") >> faceZoneNames_[surfI];
- dict.lookup("cellZone") >> cellZoneNames_[surfI];
- dict.lookup("zoneInside") >> zoneInside_[surfI];
- }
+ const dictionary& dict = surfacesDict.subDict(geomName);
- // Global perpendicular angle
- if (dict.found("perpendicularAngle"))
- {
- globalAngle[surfI] = readScalar(dict.lookup("perpendicularAngle"));
- }
+ names_[surfI] = geomName;
+ surfaces_[surfI] = geomI;
- if (dict.found("regions"))
- {
- const dictionary& regionsDict = dict.subDict("regions");
- const wordList& regionNames =
- allGeometry_[surfaces_[surfI]].regions();
+ const labelPair refLevel(dict.lookup("level"));
+ globalMinLevel[surfI] = refLevel[0];
+ globalMaxLevel[surfI] = refLevel[1];
- forAllConstIter(dictionary, regionsDict, iter)
+ // Global zone names per surface
+ if (dict.found("faceZone"))
{
- const word& key = iter().keyword();
-
- if (regionsDict.isDict(key))
- {
- // Get the dictionary for region iter.keyword()
- const dictionary& regionDict = regionsDict.subDict(key);
-
- label regionI = findIndex(regionNames, key);
- if (regionI == -1)
- {
- FatalErrorIn
- (
- "refinementSurfaces::refinementSurfaces"
- "(const searchableSurfaces&, const dictionary>&"
- ) << "No region called " << key << " on surface "
- << allGeometry_[surfaces_[surfI]].name() << endl
- << "Valid regions are " << regionNames
- << exit(FatalError);
- }
+ dict.lookup("faceZone") >> faceZoneNames_[surfI];
+ dict.lookup("cellZone") >> cellZoneNames_[surfI];
+ dict.lookup("zoneInside") >> zoneInside_[surfI];
+ }
- const labelPair refLevel(regionDict.lookup("level"));
+ // Global perpendicular angle
+ if (dict.found("perpendicularAngle"))
+ {
+ globalAngle[surfI] = readScalar
+ (
+ dict.lookup("perpendicularAngle")
+ );
+ }
- regionMinLevel[surfI].insert(regionI, refLevel[0]);
- regionMaxLevel[surfI].insert(regionI, refLevel[1]);
+ if (dict.found("regions"))
+ {
+ const dictionary& regionsDict = dict.subDict("regions");
+ const wordList& regionNames =
+ allGeometry_[surfaces_[surfI]].regions();
- if (regionDict.found("perpendicularAngle"))
+ forAll(regionNames, regionI)
+ {
+ if (regionsDict.found(regionNames[regionI]))
{
- regionAngle[surfI].insert
+ // Get the dictionary for region
+ const dictionary& regionDict = regionsDict.subDict
(
- regionI,
- readScalar(regionDict.lookup("perpendicularAngle"))
+ regionNames[regionI]
);
+
+ const labelPair refLevel(regionDict.lookup("level"));
+
+ regionMinLevel[surfI].insert(regionI, refLevel[0]);
+ regionMaxLevel[surfI].insert(regionI, refLevel[1]);
+
+ if (regionDict.found("perpendicularAngle"))
+ {
+ regionAngle[surfI].insert
+ (
+ regionI,
+ readScalar
+ (
+ regionDict.lookup("perpendicularAngle")
+ )
+ );
+ }
}
}
}
+ surfI++;
}
- surfI++;
}
// Calculate local to global region offset
label nRegions = 0;
- forAll(surfacesDict, surfI)
+ forAll(surfaces_, surfI)
{
regionOffset_[surfI] = nRegions;
nRegions += allGeometry_[surfaces_[surfI]].regions().size();
diff --git a/src/finiteVolume/cfdTools/general/MRF/MRFZone.C b/src/finiteVolume/cfdTools/general/MRF/MRFZone.C
index 8bc4d3b938109171a57dbb8513f7a4d6ec753537..ab9709e993871b6aeaf16138e3f13408b664a591 100644
--- a/src/finiteVolume/cfdTools/general/MRF/MRFZone.C
+++ b/src/finiteVolume/cfdTools/general/MRF/MRFZone.C
@@ -29,6 +29,10 @@ License
#include "volFields.H"
#include "surfaceFields.H"
#include "fvMatrices.H"
+#include "PackedList.H"
+#include "syncTools.H"
+
+#include "faceSet.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
@@ -39,12 +43,15 @@ Foam::MRFZone::MRFZone(const fvMesh& mesh, Istream& is)
dict_(is),
cellZoneID_(mesh_.cellZones().findZoneID(name_)),
faceZoneID_(mesh_.faceZones().findZoneID(name_)),
+ outsideFaces_(0),
patchNames_(dict_.lookup("patches")),
origin_(dict_.lookup("origin")),
axis_(dict_.lookup("axis")),
omega_(dict_.lookup("omega")),
Omega_("Omega", omega_*axis_)
{
+ const polyBoundaryMesh& patches = mesh_.boundaryMesh();
+
axis_ = axis_/mag(axis_);
Omega_ = omega_*axis_;
@@ -65,18 +72,71 @@ Foam::MRFZone::MRFZone(const fvMesh& mesh, Istream& is)
if (!faceZoneFound)
{
- FatalErrorIn
- (
- "Foam::MRFZone::MRFZone(const fvMesh& , const dictionary&)"
- ) << "cannot find MRF faceZone " << name_
- << exit(FatalError);
+ // Determine faces in cell zone
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+ // (does not construct cells)
+
+ const labelList& own = mesh_.faceOwner();
+ const labelList& nei = mesh_.faceNeighbour();
+
+ // Cells in zone
+ PackedBoolList zoneCell(mesh_.nCells());
+
+ if (cellZoneID_ != -1)
+ {
+ const labelList& cellLabels = mesh_.cellZones()[cellZoneID_];
+ forAll(cellLabels, i)
+ {
+ zoneCell[cellLabels[i]] = 1u;
+ }
+ }
+
+
+ // Faces in zone
+ PackedBoolList zoneFacesSet(mesh_.nFaces());
+
+ for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
+ {
+ if
+ (
+ zoneCell.get(own[faceI]) == 1u
+ || zoneCell.get(nei[faceI]) == 1u
+ )
+ {
+ zoneFacesSet[faceI] = 1u;
+ }
+ }
+ syncTools::syncFaceList(mesh_, zoneFacesSet, orEqOp<unsigned int>());
+
+
+ // Transfer to labelList
+ label n = zoneFacesSet.count();
+ outsideFaces_.setSize(n);
+ n = 0;
+ forAll(zoneFacesSet, faceI)
+ {
+ if (zoneFacesSet.get(faceI) == 1u)
+ {
+ outsideFaces_[n++] = faceI;
+ }
+ }
+
+ Info<< nl
+ << "MRFZone " << name_ << " : did not find a faceZone; using "
+ << returnReduce(outsideFaces_.size(), sumOp<label>())
+ << " faces internal to cellZone instead." << endl;
+
+
+ // Flag use of outsideFaces
+ faceZoneID_ = -2;
}
+
patchLabels_.setSize(patchNames_.size());
forAll(patchNames_, i)
{
- patchLabels_[i] = mesh_.boundaryMesh().findPatchID(patchNames_[i]);
+ patchLabels_[i] = patches.findPatchID(patchNames_[i]);
if (patchLabels_[i] == -1)
{
@@ -125,7 +185,13 @@ void Foam::MRFZone::relativeFlux(surfaceScalarField& phi) const
const vector& origin = origin_.value();
const vector& Omega = Omega_.value();
- const labelList& faces = mesh_.faceZones()[faceZoneID_];
+ // Use either zone faces or outsideFaces_
+ const labelList& faces =
+ (
+ faceZoneID_ == -2
+ ? outsideFaces_
+ : mesh_.faceZones()[faceZoneID_]
+ );
forAll(faces, i)
{
diff --git a/src/finiteVolume/cfdTools/general/MRF/MRFZone.H b/src/finiteVolume/cfdTools/general/MRF/MRFZone.H
index 3ffb4e461a853a7703ee7a5a9038bc99893c5efb..7c489ed25bba06d0e9e44334ff609ced8fbd0f57 100644
--- a/src/finiteVolume/cfdTools/general/MRF/MRFZone.H
+++ b/src/finiteVolume/cfdTools/general/MRF/MRFZone.H
@@ -26,7 +26,7 @@ Class
Foam::MRFZone
Description
- MRF zone definition based on both cell and face zones and parameters
+ MRF zone definition based on cell zone and optional face zone and parameters
obtained from a control dictionary constructed from the given stream.
The rotation of the MRF region is defined by an origin and axis of
@@ -68,18 +68,26 @@ class MRFZone
const fvMesh& mesh_;
- word name_;
+ const word name_;
- dictionary dict_;
+ const dictionary dict_;
label cellZoneID_;
+
+ //- label of face zone with faces on outside of cell zone.
+ // If -2 determines outside faces itself
label faceZoneID_;
- wordList patchNames_;
+
+ //- outside faces (only if faceZoneID = -2)
+ labelList outsideFaces_;
+
+
+ const wordList patchNames_;
labelList patchLabels_;
- dimensionedVector origin_;
+ const dimensionedVector origin_;
dimensionedVector axis_;
- dimensionedScalar omega_;
+ const dimensionedScalar omega_;
dimensionedVector Omega_;
diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/quadraticFitSnGrad/quadraticFitSnGradData.C b/src/finiteVolume/finiteVolume/snGradSchemes/quadraticFitSnGrad/quadraticFitSnGradData.C
index 133f1e79332a307bb6cf73a57a9cdc30077a732f..6bef8de896812095279f28d7a1d066104d81de39 100644
--- a/src/finiteVolume/finiteVolume/snGradSchemes/quadraticFitSnGrad/quadraticFitSnGradData.C
+++ b/src/finiteVolume/finiteVolume/snGradSchemes/quadraticFitSnGrad/quadraticFitSnGradData.C
@@ -140,11 +140,11 @@ void Foam::quadraticFitSnGradData::findFaceDirs
#ifndef SPHERICAL_GEOMETRY
if (mesh.nGeometricD() <= 2) // find the normal direcion
{
- if (mesh.directions()[0] == -1)
+ if (mesh.geometricD()[0] == -1)
{
kdir = vector(1, 0, 0);
}
- else if (mesh.directions()[1] == -1)
+ else if (mesh.geometricD()[1] == -1)
{
kdir = vector(0, 1, 0);
}
@@ -153,7 +153,7 @@ void Foam::quadraticFitSnGradData::findFaceDirs
kdir = vector(0, 0, 1);
}
}
- else // 3D so find a direction in the place of the face
+ else // 3D so find a direction in the plane of the face
{
const face& f = mesh.faces()[faci];
kdir = mesh.points()[f[0]] - mesh.points()[f[1]];
diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
index e1ef6135a6c407e754baa42ea404bcb345aa3c76..2134420eadf769f6e09dd30b5dd3d622757714ec 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
@@ -713,7 +713,7 @@ Foam::fvMatrix<Type>::H() const
(
pow
(
- psi_.mesh().directions(),
+ psi_.mesh().solutionD(),
pTraits<typename powProduct<Vector<label>, Type::rank>::type>::zero
)
);
diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C
index e8f4e6de494b64a43e0f24ffa77ec1838a4ab4ec..21f3738a02692932c687d8d02c07e365b661b0e0 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrixSolve.C
@@ -82,7 +82,7 @@ Foam::lduMatrix::solverPerformance Foam::fvMatrix<Type>::solve
(
pow
(
- psi_.mesh().directions(),
+ psi_.mesh().solutionD(),
pTraits<typename powProduct<Vector<label>, Type::rank>::type>::zero
)
);
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.C b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.C
index 0872cb81a9c3cdc33fd43e20353b745a37d17590..28a3c41e64d68edde533776ab6cb4317eacef0b0 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.C
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.C
@@ -83,11 +83,11 @@ void Foam::FitData<FitDataType, ExtendedStencil, Polynomial>::findFaceDirs
# ifndef SPHERICAL_GEOMETRY
if (mesh.nGeometricD() <= 2) // find the normal direction
{
- if (mesh.directions()[0] == -1)
+ if (mesh.geometricD()[0] == -1)
{
kdir = vector(1, 0, 0);
}
- else if (mesh.directions()[1] == -1)
+ else if (mesh.geometricD()[1] == -1)
{
kdir = vector(0, 1, 0);
}
@@ -115,7 +115,7 @@ void Foam::FitData<FitDataType, ExtendedStencil, Polynomial>::findFaceDirs
if (magk < SMALL)
{
- FatalErrorIn("findFaceDirs") << " calculated kdir = zero"
+ FatalErrorIn("findFaceDirs(..)") << " calculated kdir = zero"
<< exit(FatalError);
}
else
diff --git a/src/meshTools/Make/files b/src/meshTools/Make/files
index a93cd58400f6dfb5c138fa3e1f67b0252df6139a..9f849c58945b11d04d9193da1fe201653c09f01d 100644
--- a/src/meshTools/Make/files
+++ b/src/meshTools/Make/files
@@ -60,6 +60,7 @@ $(searchableSurface)/searchablePlane.C
$(searchableSurface)/searchablePlate.C
$(searchableSurface)/searchableSphere.C
$(searchableSurface)/searchableSurface.C
+$(searchableSurface)/searchableSurfaceCollection.C
$(searchableSurface)/searchableSurfaces.C
$(searchableSurface)/searchableSurfacesQueries.C
$(searchableSurface)/searchableSurfaceWithGaps.C
diff --git a/src/meshTools/cellDist/wallPoint/wallPoint.H b/src/meshTools/cellDist/wallPoint/wallPoint.H
index 0544b6ce361ea8d9f4143987703f83765b3996cf..82d6b385af96d0efe538574416ec70a3d41f579d 100644
--- a/src/meshTools/cellDist/wallPoint/wallPoint.H
+++ b/src/meshTools/cellDist/wallPoint/wallPoint.H
@@ -51,6 +51,12 @@ namespace Foam
// Forward declaration of classes
class polyPatch;
class polyMesh;
+class wallPoint;
+
+// Forward declaration of friend functions and operators
+Ostream& operator<<(Ostream&, const wallPoint&);
+Istream& operator>>(Istream&, wallPoint&);
+
/*---------------------------------------------------------------------------*\
Class wallPoint Declaration
@@ -78,12 +84,15 @@ class wallPoint
const scalar tol
);
+
public:
+
// Static data members
//- initial point far away.
static point greatPoint;
+
// Constructors
//- Construct null
@@ -102,6 +111,7 @@ public:
const wallPoint&
);
+
// Member Functions
// Access
@@ -184,11 +194,11 @@ public:
const scalar tol
);
+
// Member Operators
// Needed for List IO
inline bool operator==(const wallPoint&) const;
-
inline bool operator!=(const wallPoint&) const;
diff --git a/src/meshTools/meshTools/meshTools.C b/src/meshTools/meshTools/meshTools.C
index a4f57a8de844d777f6ba7dd6aa4161a95b3f5dd8..f826675b30a0fd9ed2985c1459cf3d88035cead9 100644
--- a/src/meshTools/meshTools/meshTools.C
+++ b/src/meshTools/meshTools/meshTools.C
@@ -25,7 +25,7 @@ License
\*---------------------------------------------------------------------------*/
#include "meshTools.H"
-#include "primitiveMesh.H"
+#include "polyMesh.H"
#include "hexMatcher.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@@ -635,6 +635,103 @@ Foam::label Foam::meshTools::walkFace
}
+void Foam::meshTools::constrainToMeshCentre(const polyMesh& mesh, point& pt)
+{
+ const Vector<label>& dirs = mesh.geometricD();
+ const point& min = mesh.bounds().min();
+ const point& max = mesh.bounds().max();
+
+ for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+ {
+ if (dirs[cmpt] == -1)
+ {
+ pt[cmpt] = 0.5*(min[cmpt]+max[cmpt]);
+ }
+ }
+}
+
+
+void Foam::meshTools::constrainToMeshCentre
+(
+ const polyMesh& mesh,
+ pointField& pts
+)
+{
+ const Vector<label>& dirs = mesh.geometricD();
+ const point& min = mesh.bounds().min();
+ const point& max = mesh.bounds().max();
+
+ bool isConstrained = false;
+ for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+ {
+ if (dirs[cmpt] == -1)
+ {
+ isConstrained = true;
+ break;
+ }
+ }
+
+ if (isConstrained)
+ {
+ forAll(pts, i)
+ {
+ for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+ {
+ if (dirs[cmpt] == -1)
+ {
+ pts[i][cmpt] = 0.5*(min[cmpt]+max[cmpt]);
+ }
+ }
+ }
+ }
+}
+
+
+//- Set the constrained components of directions/velocity to zero
+void Foam::meshTools::constrainDirection(const polyMesh& mesh, vector& d)
+{
+ const Vector<label>& dirs = mesh.geometricD();
+
+ for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+ {
+ if (dirs[cmpt] == -1)
+ {
+ d[cmpt] = 0.0;
+ }
+ }
+}
+
+
+void Foam::meshTools::constrainDirection(const polyMesh& mesh, vectorField& d)
+{
+ const Vector<label>& dirs = mesh.geometricD();
+
+ bool isConstrained = false;
+ for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+ {
+ if (dirs[cmpt] == -1)
+ {
+ isConstrained = true;
+ break;
+ }
+ }
+
+ if (isConstrained)
+ {
+ forAll(d, i)
+ {
+ for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
+ {
+ if (dirs[cmpt] == -1)
+ {
+ d[i][cmpt] = 0.0;
+ }
+ }
+ }
+ }
+}
+
+
void Foam::meshTools::getParallelEdges
(
const primitiveMesh& mesh,
diff --git a/src/meshTools/meshTools/meshTools.H b/src/meshTools/meshTools/meshTools.H
index 304159e04589cc22006931636049ba4d128f10ef..3710d5af06250dc654465ebcba5bb47697f13465 100644
--- a/src/meshTools/meshTools/meshTools.H
+++ b/src/meshTools/meshTools/meshTools.H
@@ -50,7 +50,7 @@ namespace Foam
{
class primitiveMesh;
-
+class polyMesh;
/*---------------------------------------------------------------------------*\
Namespace meshTools Declaration
@@ -81,196 +81,212 @@ namespace meshTools
static const label pXpYpZMask = 1 << pXpYpZ;
- //- Check if n is in same direction as normals of all faceLabels
- bool visNormal
- (
- const vector& n,
- const vectorField& faceNormals,
- const labelList& faceLabels
- );
-
- //- Calculate point normals on a 'box' mesh (all edges aligned with
- // coordinate axes)
- vectorField calcBoxPointNormals(const primitivePatch& pp);
-
- //- Normalized edge vector
- vector normEdgeVec(const primitiveMesh&, const label edgeI);
-
- //- Write obj representation of point
- void writeOBJ
- (
- Ostream& os,
- const point& pt
- );
-
-
- //- Write obj representation of faces subset
- void writeOBJ
- (
- Ostream& os,
- const faceList&,
- const pointField&,
- const labelList& faceLabels
- );
-
- //- Write obj representation of faces
- void writeOBJ
- (
- Ostream& os,
- const faceList&,
- const pointField&
- );
-
- //- Write obj representation of cell subset
- void writeOBJ
- (
- Ostream& os,
- const cellList&,
- const faceList&,
- const pointField&,
- const labelList& cellLabels
- );
-
- //- Is edge used by cell
- bool edgeOnCell
- (
- const primitiveMesh&,
- const label cellI,
- const label edgeI
- );
-
- //- Is edge used by face
- bool edgeOnFace
- (
- const primitiveMesh&,
- const label faceI,
- const label edgeI
- );
-
- //- Is face used by cell
- bool faceOnCell
- (
- const primitiveMesh&,
- const label cellI,
- const label faceI
- );
-
- //- Return edge among candidates that uses the two vertices.
- label findEdge
- (
- const edgeList& edges,
- const labelList& candidates,
- const label v0,
- const label v1
- );
-
- //- Return edge between two vertices. Returns -1 if no edge.
- label findEdge
- (
- const primitiveMesh&,
- const label v0,
- const label v1
- );
-
- //- Return edge shared by two faces. Throws error if no edge found.
- label getSharedEdge
- (
- const primitiveMesh&,
- const label f0,
- const label f1
- );
-
- //- Return face shared by two cells. Throws error if none found.
- label getSharedFace
- (
- const primitiveMesh&,
- const label cell0,
- const label cell1
- );
-
- //- Get faces on cell using edgeI. Throws error if no two found.
- void getEdgeFaces
- (
- const primitiveMesh&,
- const label cellI,
- const label edgeI,
- label& face0,
- label& face1
- );
-
- //- Return label of other edge (out of candidates edgeLabels)
- // connected to vertex but not edgeI. Throws error if none found.
- label otherEdge
- (
- const primitiveMesh&,
- const labelList& edgeLabels,
- const label edgeI,
- const label vertI
- );
-
- //- Return face on cell using edgeI but not faceI. Throws error
- // if none found.
- label otherFace
- (
- const primitiveMesh&,
- const label cellI,
- const label faceI,
- const label edgeI
- );
-
- //- Return cell on other side of face. Throws error
- // if face not internal.
- label otherCell
- (
- const primitiveMesh&,
- const label cellI,
- const label faceI
- );
-
- //- Returns label of edge nEdges away from startEdge (in the direction
- // of startVertI)
- label walkFace
- (
- const primitiveMesh&,
- const label faceI,
- const label startEdgeI,
- const label startVertI,
- const label nEdges
- );
+ // Normal handling
+
+ //- Check if n is in same direction as normals of all faceLabels
+ bool visNormal
+ (
+ const vector& n,
+ const vectorField& faceNormals,
+ const labelList& faceLabels
+ );
+
+ //- Calculate point normals on a 'box' mesh (all edges aligned with
+ // coordinate axes)
+ vectorField calcBoxPointNormals(const primitivePatch& pp);
+
+ //- Normalized edge vector
+ vector normEdgeVec(const primitiveMesh&, const label edgeI);
+
+
+ // OBJ writing
+
+ //- Write obj representation of point
+ void writeOBJ
+ (
+ Ostream& os,
+ const point& pt
+ );
+
+ //- Write obj representation of faces subset
+ void writeOBJ
+ (
+ Ostream& os,
+ const faceList&,
+ const pointField&,
+ const labelList& faceLabels
+ );
+
+ //- Write obj representation of faces
+ void writeOBJ
+ (
+ Ostream& os,
+ const faceList&,
+ const pointField&
+ );
+
+ //- Write obj representation of cell subset
+ void writeOBJ
+ (
+ Ostream& os,
+ const cellList&,
+ const faceList&,
+ const pointField&,
+ const labelList& cellLabels
+ );
+
+
+ // Cell/face/edge walking
+
+ //- Is edge used by cell
+ bool edgeOnCell
+ (
+ const primitiveMesh&,
+ const label cellI,
+ const label edgeI
+ );
+
+ //- Is edge used by face
+ bool edgeOnFace
+ (
+ const primitiveMesh&,
+ const label faceI,
+ const label edgeI
+ );
+
+ //- Is face used by cell
+ bool faceOnCell
+ (
+ const primitiveMesh&,
+ const label cellI,
+ const label faceI
+ );
+
+ //- Return edge among candidates that uses the two vertices.
+ label findEdge
+ (
+ const edgeList& edges,
+ const labelList& candidates,
+ const label v0,
+ const label v1
+ );
+
+ //- Return edge between two vertices. Returns -1 if no edge.
+ label findEdge
+ (
+ const primitiveMesh&,
+ const label v0,
+ const label v1
+ );
+
+ //- Return edge shared by two faces. Throws error if no edge found.
+ label getSharedEdge
+ (
+ const primitiveMesh&,
+ const label f0,
+ const label f1
+ );
+
+ //- Return face shared by two cells. Throws error if none found.
+ label getSharedFace
+ (
+ const primitiveMesh&,
+ const label cell0,
+ const label cell1
+ );
+
+ //- Get faces on cell using edgeI. Throws error if no two found.
+ void getEdgeFaces
+ (
+ const primitiveMesh&,
+ const label cellI,
+ const label edgeI,
+ label& face0,
+ label& face1
+ );
+
+ //- Return label of other edge (out of candidates edgeLabels)
+ // connected to vertex but not edgeI. Throws error if none found.
+ label otherEdge
+ (
+ const primitiveMesh&,
+ const labelList& edgeLabels,
+ const label edgeI,
+ const label vertI
+ );
+
+ //- Return face on cell using edgeI but not faceI. Throws error
+ // if none found.
+ label otherFace
+ (
+ const primitiveMesh&,
+ const label cellI,
+ const label faceI,
+ const label edgeI
+ );
+
+ //- Return cell on other side of face. Throws error
+ // if face not internal.
+ label otherCell
+ (
+ const primitiveMesh&,
+ const label cellI,
+ const label faceI
+ );
+
+ //- Returns label of edge nEdges away from startEdge (in the direction
+ // of startVertI)
+ label walkFace
+ (
+ const primitiveMesh&,
+ const label faceI,
+ const label startEdgeI,
+ const label startVertI,
+ const label nEdges
+ );
+
+
+ // Constraints on position
+
+ //- Set the constrained components of position to mesh centre
+ void constrainToMeshCentre(const polyMesh&, point&);
+ void constrainToMeshCentre(const polyMesh&, pointField&);
+
+ //- Set the constrained components of directions/velocity to zero
+ void constrainDirection(const polyMesh&, vector&);
+ void constrainDirection(const polyMesh&, vectorField&);
- //
// Hex only functionality.
- //
-
- //- Given edge on hex find other 'parallel', non-connected edges.
- void getParallelEdges
- (
- const primitiveMesh&,
- const label cellI,
- const label e0,
- label&,
- label&,
- label&
- );
-
- //- Given edge on hex find all 'parallel' (i.e. non-connected)
- // edges and average direction of them
- vector edgeToCutDir
- (
- const primitiveMesh&,
- const label cellI,
- const label edgeI
- );
-
- //- Reverse of edgeToCutDir: given direction find edge bundle and
- // return one of them.
- label cutDirToEdge
- (
- const primitiveMesh&,
- const label cellI,
- const vector& cutDir
- );
+
+ //- Given edge on hex find other 'parallel', non-connected edges.
+ void getParallelEdges
+ (
+ const primitiveMesh&,
+ const label cellI,
+ const label e0,
+ label&,
+ label&,
+ label&
+ );
+
+ //- Given edge on hex find all 'parallel' (i.e. non-connected)
+ // edges and average direction of them
+ vector edgeToCutDir
+ (
+ const primitiveMesh&,
+ const label cellI,
+ const label edgeI
+ );
+
+ //- Reverse of edgeToCutDir: given direction find edge bundle and
+ // return one of them.
+ label cutDirToEdge
+ (
+ const primitiveMesh&,
+ const label cellI,
+ const vector& cutDir
+ );
} // End namespace meshTools
diff --git a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C
index 06020467036039319af49714caa501e7d74ae6b5..7f86c0f6581fcd94504b9cb7d39fdc67284a4506 100644
--- a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C
+++ b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.C
@@ -76,10 +76,8 @@ bool Foam::distributedTriSurfaceMesh::read()
// Distribution type
distType_ = distributionTypeNames_.read(dict_.lookup("distributionType"));
- if (dict_.found("mergeDistance"))
- {
- dict_.lookup("mergeDistance") >> mergeDist_;
- }
+ // Merge distance
+ mergeDist_ = readScalar(dict_.lookup("mergeDistance"));
return true;
}
@@ -1341,25 +1339,44 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh
dict_(io, dict)
{
read();
+
+ if (debug)
+ {
+ Info<< "Constructed from triSurface:" << endl;
+ writeStats(Info);
+
+ labelList nTris(Pstream::nProcs());
+ nTris[Pstream::myProcNo()] = triSurface::size();
+ Pstream::gatherList(nTris);
+ Pstream::scatterList(nTris);
+
+ Info<< endl<< "\tproc\ttris\tbb" << endl;
+ forAll(nTris, procI)
+ {
+ Info<< '\t' << procI << '\t' << nTris[procI]
+ << '\t' << procBb_[procI] << endl;
+ }
+ Info<< endl;
+ }
}
Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh(const IOobject& io)
:
- triSurfaceMesh(io),
-// triSurfaceMesh
-// (
-// IOobject
-// (
-// io.name(),
-// io.db().time().findInstanceDir(io.local()),
-// io.local(),
-// io.db(),
-// io.readOpt(),
-// io.writeOpt(),
-// io.registerObject()
-// )
-// ),
+ //triSurfaceMesh(io),
+ triSurfaceMesh
+ (
+ IOobject
+ (
+ io.name(),
+ io.time().findInstance(io.local(), word::null),
+ io.local(),
+ io.db(),
+ io.readOpt(),
+ io.writeOpt(),
+ io.registerObject()
+ )
+ ),
dict_
(
IOobject
@@ -1375,6 +1392,26 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh(const IOobject& io)
)
{
read();
+
+ if (debug)
+ {
+ Info<< "Read distributedTriSurface from " << io.objectPath()
+ << ':' << endl;
+ writeStats(Info);
+
+ labelList nTris(Pstream::nProcs());
+ nTris[Pstream::myProcNo()] = triSurface::size();
+ Pstream::gatherList(nTris);
+ Pstream::scatterList(nTris);
+
+ Info<< endl<< "\tproc\ttris\tbb" << endl;
+ forAll(nTris, procI)
+ {
+ Info<< '\t' << procI << '\t' << nTris[procI]
+ << '\t' << procBb_[procI] << endl;
+ }
+ Info<< endl;
+ }
}
@@ -1384,21 +1421,21 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh
const dictionary& dict
)
:
- triSurfaceMesh(io, dict),
-// triSurfaceMesh
-// (
-// IOobject
-// (
-// io.name(),
-// io.db().time().findInstanceDir(io.local()),
-// io.local(),
-// io.db(),
-// io.readOpt(),
-// io.writeOpt(),
-// io.registerObject()
-// ),
-// dict
-// ),
+ //triSurfaceMesh(io, dict),
+ triSurfaceMesh
+ (
+ IOobject
+ (
+ io.name(),
+ io.time().findInstance(io.local(), word::null),
+ io.local(),
+ io.db(),
+ io.readOpt(),
+ io.writeOpt(),
+ io.registerObject()
+ ),
+ dict
+ ),
dict_
(
IOobject
@@ -1414,6 +1451,26 @@ Foam::distributedTriSurfaceMesh::distributedTriSurfaceMesh
)
{
read();
+
+ if (debug)
+ {
+ Info<< "Read distributedTriSurface from " << io.objectPath()
+ << " and dictionary:" << endl;
+ writeStats(Info);
+
+ labelList nTris(Pstream::nProcs());
+ nTris[Pstream::myProcNo()] = triSurface::size();
+ Pstream::gatherList(nTris);
+ Pstream::scatterList(nTris);
+
+ Info<< endl<< "\tproc\ttris\tbb" << endl;
+ forAll(nTris, procI)
+ {
+ Info<< '\t' << procI << '\t' << nTris[procI]
+ << '\t' << procBb_[procI] << endl;
+ }
+ Info<< endl;
+ }
}
@@ -2416,6 +2473,8 @@ void Foam::distributedTriSurfaceMesh::writeStats(Ostream& os) const
boundBox bb;
label nPoints;
calcBounds(bb, nPoints);
+ reduce(bb.min(), minOp<point>());
+ reduce(bb.max(), maxOp<point>());
os << "Triangles : " << returnReduce(triSurface::size(), sumOp<label>())
<< endl
diff --git a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H
index 733cd24891e21befee9e9318fa301180bc3f8966..4d5f7dfd173edd6ccb3c559221de0564d46edc50 100644
--- a/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H
+++ b/src/meshTools/searchableSurface/distributedTriSurfaceMesh.H
@@ -91,15 +91,14 @@ public:
private:
+ // Private member data
+
//- Merging distance
scalar mergeDist_;
//- Decomposition used when independently decomposing surface.
autoPtr<decompositionMethod> decomposer_;
-
- // Private member data
-
//- Bounding box settings
IOdictionary dict_;
@@ -317,10 +316,11 @@ public:
const dictionary& dict
);
- //- Construct read
+ //- Construct read. Does findInstance to find io.local().
distributedTriSurfaceMesh(const IOobject& io);
- //- Construct from dictionary (used by searchableSurface)
+ //- Construct from dictionary (used by searchableSurface).
+ // Does read. Does findInstance to find io.local().
distributedTriSurfaceMesh
(
const IOobject& io,
diff --git a/src/meshTools/searchableSurface/searchableSphere.C b/src/meshTools/searchableSurface/searchableSphere.C
index ba5a1888cfb343f39605b49b8111083427d7562c..287ae6743c2bb508c0ef64322036fac617123dcd 100644
--- a/src/meshTools/searchableSurface/searchableSphere.C
+++ b/src/meshTools/searchableSurface/searchableSphere.C
@@ -53,7 +53,14 @@ Foam::pointIndexHit Foam::searchableSphere::findNearest
if (nearestDistSqr > sqr(magN-radius_))
{
- info.rawPoint() = centre_ + n/magN*radius_;
+ if (magN < ROOTVSMALL)
+ {
+ info.rawPoint() = centre_ + vector(1,0,0)/magN*radius_;
+ }
+ else
+ {
+ info.rawPoint() = centre_ + n/magN*radius_;
+ }
info.setHit();
info.setIndex(0);
}
diff --git a/src/meshTools/searchableSurface/searchableSurfaceCollection.C b/src/meshTools/searchableSurface/searchableSurfaceCollection.C
new file mode 100644
index 0000000000000000000000000000000000000000..ce6ddd643122703e9214fc6041a5a76265b098f4
--- /dev/null
+++ b/src/meshTools/searchableSurface/searchableSurfaceCollection.C
@@ -0,0 +1,525 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2008 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 2 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, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "searchableSurfaceCollection.H"
+#include "addToRunTimeSelectionTable.H"
+#include "SortableList.H"
+#include "Time.H"
+#include "ListOps.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+defineTypeNameAndDebug(searchableSurfaceCollection, 0);
+addToRunTimeSelectionTable(searchableSurface, searchableSurfaceCollection, dict);
+
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+void Foam::searchableSurfaceCollection::findNearest
+(
+ const pointField& samples,
+ scalarField& minDistSqr,
+ List<pointIndexHit>& nearestInfo,
+ labelList& nearestSurf
+) const
+{
+ // Initialise
+ nearestInfo.setSize(samples.size());
+ nearestInfo = pointIndexHit();
+ nearestSurf.setSize(samples.size());
+ nearestSurf = -1;
+
+ List<pointIndexHit> hitInfo(samples.size());
+
+ const scalarField localMinDistSqr(samples.size(), GREAT);
+
+ forAll(subGeom_, surfI)
+ {
+ // Transform then divide
+ tmp<pointField> localSamples = cmptDivide
+ (
+ transform_[surfI].localPosition(samples),
+ scale_[surfI]
+ );
+
+ subGeom_[surfI].findNearest(localSamples, localMinDistSqr, hitInfo);
+
+ forAll(hitInfo, pointI)
+ {
+ if (hitInfo[pointI].hit())
+ {
+ // Rework back into global coordinate sys. Multiply then
+ // transform
+ point globalPt = transform_[surfI].globalPosition
+ (
+ cmptMultiply
+ (
+ hitInfo[pointI].rawPoint(),
+ scale_[surfI]
+ )
+ );
+
+ scalar distSqr = magSqr(globalPt - samples[pointI]);
+
+ if (distSqr < minDistSqr[pointI])
+ {
+ minDistSqr[pointI] = distSqr;
+ nearestInfo[pointI].setPoint(globalPt);
+ nearestInfo[pointI].setHit();
+ nearestInfo[pointI].setIndex(hitInfo[pointI].index());
+ nearestSurf[pointI] = surfI;
+ }
+ }
+ }
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::searchableSurfaceCollection::searchableSurfaceCollection
+(
+ const IOobject& io,
+ const dictionary& dict
+)
+:
+ searchableSurface(io),
+ instance_(dict.size()),
+ scale_(dict.size()),
+ transform_(dict.size()),
+ subGeom_(dict.size())
+{
+ Info<< "SearchableCollection : " << name() << endl;
+
+ label surfI = 0;
+ forAllConstIter(dictionary, dict, iter)
+ {
+ if (dict.isDict(iter().keyword()))
+ {
+ instance_[surfI] = iter().keyword();
+
+ const dictionary& subDict = dict.subDict(instance_[surfI]);
+
+ scale_[surfI] = subDict.lookup("scale");
+ transform_.set
+ (
+ surfI,
+ coordinateSystem::New
+ (
+ "",
+ subDict.subDict("transform")
+ )
+ );
+
+ const word subGeomName(subDict.lookup("surface"));
+ //Pout<< "Trying to find " << subGeomName << endl;
+
+ const searchableSurface& s =
+ io.db().lookupObject<searchableSurface>(subGeomName);
+
+ subGeom_.set(surfI, &const_cast<searchableSurface&>(s));
+
+ Info<< " instance : " << instance_[surfI] << endl;
+ Info<< " surface : " << s.name() << endl;
+ Info<< " scale : " << scale_[surfI] << endl;
+ Info<< " coordsys : " << transform_[surfI] << endl;
+
+ surfI++;
+ }
+ }
+ instance_.setSize(surfI);
+ scale_.setSize(surfI);
+ transform_.setSize(surfI);
+ subGeom_.setSize(surfI);
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::searchableSurfaceCollection::~searchableSurfaceCollection()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+const Foam::wordList& Foam::searchableSurfaceCollection::regions() const
+{
+ if (regions_.size() == 0)
+ {
+ regionOffset_.setSize(subGeom_.size());
+
+ DynamicList<word> allRegions;
+ forAll(subGeom_, surfI)
+ {
+ regionOffset_[surfI] = allRegions.size();
+
+ const wordList& subRegions = subGeom_[surfI].regions();
+
+ forAll(subRegions, i)
+ {
+ //allRegions.append(subRegions[i] + "_" + Foam::name(surfI));
+ allRegions.append(instance_[surfI] + "_" + subRegions[i]);
+ }
+ }
+ regions_.transfer(allRegions.shrink());
+ }
+ return regions_;
+}
+
+
+Foam::label Foam::searchableSurfaceCollection::size() const
+{
+ label n = 0;
+ forAll(subGeom_, surfI)
+ {
+ n += subGeom_[surfI].size();
+ }
+ return n;
+}
+
+
+void Foam::searchableSurfaceCollection::findNearest
+(
+ const pointField& samples,
+ const scalarField& nearestDistSqr,
+ List<pointIndexHit>& nearestInfo
+) const
+{
+ // How to scale distance?
+ scalarField minDistSqr(nearestDistSqr);
+
+ labelList nearestSurf;
+ findNearest
+ (
+ samples,
+ minDistSqr,
+ nearestInfo,
+ nearestSurf
+ );
+}
+
+
+void Foam::searchableSurfaceCollection::findLine
+(
+ const pointField& start,
+ const pointField& end,
+ List<pointIndexHit>& info
+) const
+{
+ info.setSize(start.size());
+ info = pointIndexHit();
+
+ // Current nearest (to start) intersection
+ pointField nearest(end);
+
+ List<pointIndexHit> hitInfo(start.size());
+
+ forAll(subGeom_, surfI)
+ {
+ // Starting point
+ tmp<pointField> e0 = cmptDivide
+ (
+ transform_[surfI].localPosition
+ (
+ start
+ ),
+ scale_[surfI]
+ );
+
+ // Current best end point
+ tmp<pointField> e1 = cmptDivide
+ (
+ transform_[surfI].localPosition
+ (
+ nearest
+ ),
+ scale_[surfI]
+ );
+
+ subGeom_[surfI].findLine(e0, e1, hitInfo);
+
+ forAll(hitInfo, pointI)
+ {
+ if (hitInfo[pointI].hit())
+ {
+ // Transform back to global coordinate sys.
+ nearest[pointI] = transform_[surfI].globalPosition
+ (
+ cmptMultiply
+ (
+ hitInfo[pointI].rawPoint(),
+ scale_[surfI]
+ )
+ );
+ info[pointI] = hitInfo[pointI];
+ info[pointI].rawPoint() = nearest[pointI];
+ }
+ }
+ }
+
+
+ // Debug check
+ if (false)
+ {
+ forAll(info, pointI)
+ {
+ if (info[pointI].hit())
+ {
+ vector n(end[pointI] - start[pointI]);
+ scalar magN = mag(n);
+
+ if (magN > SMALL)
+ {
+ n /= mag(n);
+
+ scalar s = ((info[pointI].rawPoint()-start[pointI])&n);
+
+ if (s < 0 || s > 1)
+ {
+ FatalErrorIn
+ (
+ "searchableSurfaceCollection::findLine(..)"
+ ) << "point:" << info[pointI]
+ << " s:" << s
+ << " outside vector "
+ << " start:" << start[pointI]
+ << " end:" << end[pointI]
+ << abort(FatalError);
+ }
+ }
+ }
+ }
+ }
+}
+
+
+void Foam::searchableSurfaceCollection::findLineAny
+(
+ const pointField& start,
+ const pointField& end,
+ List<pointIndexHit>& info
+) const
+{
+ // To be done ...
+ findLine(start, end, info);
+}
+
+
+void Foam::searchableSurfaceCollection::findLineAll
+(
+ const pointField& start,
+ const pointField& end,
+ List<List<pointIndexHit> >& info
+) const
+{
+ // To be done. Assume for now only one intersection.
+ List<pointIndexHit> nearestInfo;
+ findLine(start, end, nearestInfo);
+
+ info.setSize(start.size());
+ forAll(info, pointI)
+ {
+ if (nearestInfo[pointI].hit())
+ {
+ info[pointI].setSize(1);
+ info[pointI][0] = nearestInfo[pointI];
+ }
+ else
+ {
+ info[pointI].clear();
+ }
+ }
+}
+
+
+void Foam::searchableSurfaceCollection::getRegion
+(
+ const List<pointIndexHit>& info,
+ labelList& region
+) const
+{
+ if (subGeom_.size() == 0)
+ {}
+ else if (subGeom_.size() == 1)
+ {
+ subGeom_[0].getRegion(info, region);
+ }
+ else
+ {
+ region.setSize(info.size());
+ region = -1;
+
+ // Which region did point come from. Retest for now to see which
+ // surface it originates from - crap solution! Should use global indices
+ // in index inside pointIndexHit to do this better.
+
+ pointField samples(info.size());
+ forAll(info, pointI)
+ {
+ if (info[pointI].hit())
+ {
+ samples[pointI] = info[pointI].hitPoint();
+ }
+ else
+ {
+ samples[pointI] = vector::zero;
+ }
+ }
+ //scalarField minDistSqr(info.size(), SMALL);
+ scalarField minDistSqr(info.size(), GREAT);
+
+ labelList nearestSurf;
+ List<pointIndexHit> nearestInfo;
+ findNearest
+ (
+ samples,
+ minDistSqr,
+ nearestInfo,
+ nearestSurf
+ );
+
+ // Check
+ {
+ forAll(info, pointI)
+ {
+ if (info[pointI].hit() && nearestSurf[pointI] == -1)
+ {
+ FatalErrorIn
+ (
+ "searchableSurfaceCollection::getRegion(..)"
+ ) << "pointI:" << pointI
+ << " sample:" << samples[pointI]
+ << " nearest:" << nearestInfo[pointI]
+ << " nearestsurf:" << nearestSurf[pointI]
+ << abort(FatalError);
+ }
+ }
+ }
+
+ forAll(subGeom_, surfI)
+ {
+ // Collect points from my surface
+ labelList indices(findIndices(nearestSurf, surfI));
+
+ labelList surfRegion;
+ subGeom_[surfI].getRegion
+ (
+ IndirectList<pointIndexHit>(info, indices),
+ surfRegion
+ );
+ forAll(indices, i)
+ {
+ region[indices[i]] = regionOffset_[surfI] + surfRegion[i];
+ }
+ }
+ }
+}
+
+
+void Foam::searchableSurfaceCollection::getNormal
+(
+ const List<pointIndexHit>& info,
+ vectorField& normal
+) const
+{
+ if (subGeom_.size() == 0)
+ {}
+ else if (subGeom_.size() == 1)
+ {
+ subGeom_[0].getNormal(info, normal);
+ }
+ else
+ {
+ normal.setSize(info.size());
+
+ // See above - crap retest to find surface point originates from.
+ pointField samples(info.size());
+ forAll(info, pointI)
+ {
+ if (info[pointI].hit())
+ {
+ samples[pointI] = info[pointI].hitPoint();
+ }
+ else
+ {
+ samples[pointI] = vector::zero;
+ }
+ }
+ //scalarField minDistSqr(info.size(), SMALL);
+ scalarField minDistSqr(info.size(), GREAT);
+
+ labelList nearestSurf;
+ List<pointIndexHit> nearestInfo;
+ findNearest
+ (
+ samples,
+ minDistSqr,
+ nearestInfo,
+ nearestSurf
+ );
+
+
+ forAll(subGeom_, surfI)
+ {
+ // Collect points from my surface
+ labelList indices(findIndices(nearestSurf, surfI));
+
+ vectorField surfNormal;
+ subGeom_[surfI].getNormal
+ (
+ IndirectList<pointIndexHit>(info, indices),
+ surfNormal
+ );
+ forAll(indices, i)
+ {
+ normal[indices[i]] = surfNormal[i];
+ }
+ }
+ }
+}
+
+
+void Foam::searchableSurfaceCollection::getVolumeType
+(
+ const pointField& points,
+ List<volumeType>& volType
+) const
+{
+ FatalErrorIn
+ (
+ "searchableSurfaceCollection::getVolumeType(const pointField&"
+ ", List<volumeType>&) const"
+ ) << "Volume type not supported for collection."
+ << exit(FatalError);
+}
+
+
+// ************************************************************************* //
diff --git a/src/meshTools/searchableSurface/searchableSurfaceCollection.H b/src/meshTools/searchableSurface/searchableSurfaceCollection.H
new file mode 100644
index 0000000000000000000000000000000000000000..c85eaa10264f69939fae9f2c2924522b2cefcb8e
--- /dev/null
+++ b/src/meshTools/searchableSurface/searchableSurfaceCollection.H
@@ -0,0 +1,212 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2008 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 2 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, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::searchableSurfaceCollection
+
+Description
+ Union of transformed searchableSurfaces
+
+SourceFiles
+ searchableSurfaceCollection.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef searchableSurfaceCollection_H
+#define searchableSurfaceCollection_H
+
+#include "searchableSurface.H"
+#include "treeBoundBox.H"
+#include "coordinateSystem.H"
+#include "UPtrList.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward declaration of classes
+
+/*---------------------------------------------------------------------------*\
+ Class searchableSurfaceCollection Declaration
+\*---------------------------------------------------------------------------*/
+
+class searchableSurfaceCollection
+:
+ public searchableSurface
+{
+private:
+
+ // Private Member Data
+
+ // Per instance data
+
+ //- instance name
+ wordList instance_;
+
+ //- scaling vector
+ vectorField scale_;
+
+ //- transformation
+ PtrList<coordinateSystem> transform_;
+
+ UPtrList<searchableSurface> subGeom_;
+
+ //- Region names
+ mutable wordList regions_;
+ //- From individual regions to collection regions
+ mutable labelList regionOffset_;
+
+
+ // Private Member Functions
+
+ //- Find point nearest to sample. Updates minDistSqr. Sets nearestInfo
+ // and surface index
+ void findNearest
+ (
+ const pointField& samples,
+ scalarField& minDistSqr,
+ List<pointIndexHit>& nearestInfo,
+ labelList& nearestSurf
+ ) const;
+
+
+ //- Disallow default bitwise copy construct
+ searchableSurfaceCollection(const searchableSurfaceCollection&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const searchableSurfaceCollection&);
+
+
+public:
+
+ //- Runtime type information
+ TypeName("searchableSurfaceCollection");
+
+
+ // Constructors
+
+ //- Construct from dictionary (used by searchableSurface)
+ searchableSurfaceCollection
+ (
+ const IOobject& io,
+ const dictionary& dict
+ );
+
+ // Destructor
+
+ virtual ~searchableSurfaceCollection();
+
+
+ // Member Functions
+
+ virtual const wordList& regions() const;
+
+ //- Whether supports volume type below
+ virtual bool hasVolumeType() const
+ {
+ return false;
+ }
+
+ //- Range of local indices that can be returned.
+ virtual label size() const;
+
+
+ // Multiple point queries.
+
+ virtual void findNearest
+ (
+ const pointField& sample,
+ const scalarField& nearestDistSqr,
+ List<pointIndexHit>&
+ ) const;
+
+ virtual void findLine
+ (
+ const pointField& start,
+ const pointField& end,
+ List<pointIndexHit>&
+ ) const;
+
+ virtual void findLineAny
+ (
+ const pointField& start,
+ const pointField& end,
+ List<pointIndexHit>&
+ ) const;
+
+ //- Get all intersections in order from start to end.
+ virtual void findLineAll
+ (
+ const pointField& start,
+ const pointField& end,
+ List<List<pointIndexHit> >&
+ ) const;
+
+ //- From a set of points and indices get the region
+ virtual void getRegion
+ (
+ const List<pointIndexHit>&,
+ labelList& region
+ ) const;
+
+ //- From a set of points and indices get the normal
+ virtual void getNormal
+ (
+ const List<pointIndexHit>&,
+ vectorField& normal
+ ) const;
+
+ //- Determine type (inside/outside/mixed) for point. unknown if
+ // cannot be determined (e.g. non-manifold surface)
+ virtual void getVolumeType
+ (
+ const pointField&,
+ List<volumeType>&
+ ) const;
+
+
+ // regIOobject implementation
+
+ bool writeData(Ostream&) const
+ {
+ notImplemented
+ (
+ "searchableSurfaceCollection::writeData(Ostream&) const"
+ );
+ return false;
+ }
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/meshTools/searchableSurface/searchableSurfaces.C b/src/meshTools/searchableSurface/searchableSurfaces.C
index 46c968f15369b62c8f47abe05d5fcad7eb8449cd..cec34bd7e792669ace03bad6f68ea2e0b12ffe9b 100644
--- a/src/meshTools/searchableSurface/searchableSurfaces.C
+++ b/src/meshTools/searchableSurface/searchableSurfaces.C
@@ -39,8 +39,6 @@ defineTypeNameAndDebug(searchableSurfaces, 0);
}
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// Construct with length.
@@ -187,14 +185,6 @@ Foam::searchableSurfaces::searchableSurfaces
// their object name. Maybe have stlTriSurfaceMesh which appends .stl
// when reading/writing?
namedIO().rename(key); // names_[surfI]
- if (namedIO().local() != word::null)
- {
- namedIO().instance() = namedIO().time().findInstance
- (
- namedIO().local(),
- namedIO().name()
- );
- }
// Create and hook surface
set
diff --git a/src/meshTools/searchableSurface/searchableSurfaces.H b/src/meshTools/searchableSurface/searchableSurfaces.H
index ea91fcdee958a6d45a9fdb922b286f8580544477..af02ae188396fd12b255b2e308f6dff3cd970de0 100644
--- a/src/meshTools/searchableSurface/searchableSurfaces.H
+++ b/src/meshTools/searchableSurface/searchableSurfaces.H
@@ -69,6 +69,8 @@ class searchableSurfaces
labelList allSurfaces_;
+ // Private Member Functions
+
//- Disallow default bitwise copy construct
searchableSurfaces(const searchableSurfaces&);
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.C b/src/meshTools/searchableSurface/triSurfaceMesh.C
index 7b9abc1aec6664b9efa541cbb128a7a9e27a5528..9958afd1b155167a1e08bc9f3926f94cad434532 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.C
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.C
@@ -43,6 +43,64 @@ addToRunTimeSelectionTable(searchableSurface, triSurfaceMesh, dict);
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+//// Special version of Time::findInstance that does not check headerOk
+//// to search for instances of raw files
+//Foam::word Foam::triSurfaceMesh::findRawInstance
+//(
+// const Time& runTime,
+// const fileName& dir,
+// const word& name
+//)
+//{
+// // Check current time first
+// if (isFile(runTime.path()/runTime.timeName()/dir/name))
+// {
+// return runTime.timeName();
+// }
+// instantList ts = runTime.times();
+// label instanceI;
+//
+// for (instanceI = ts.size()-1; instanceI >= 0; --instanceI)
+// {
+// if (ts[instanceI].value() <= runTime.timeOutputValue())
+// {
+// break;
+// }
+// }
+//
+// // continue searching from here
+// for (; instanceI >= 0; --instanceI)
+// {
+// if (isFile(runTime.path()/ts[instanceI].name()/dir/name))
+// {
+// return ts[instanceI].name();
+// }
+// }
+//
+//
+// // not in any of the time directories, try constant
+//
+// // Note. This needs to be a hard-coded constant, rather than the
+// // constant function of the time, because the latter points to
+// // the case constant directory in parallel cases
+//
+// if (isFile(runTime.path()/runTime.constant()/dir/name))
+// {
+// return runTime.constant();
+// }
+//
+// FatalErrorIn
+// (
+// "searchableSurfaces::findRawInstance"
+// "(const Time&, const fileName&, const word&)"
+// ) << "Cannot find file \"" << name << "\" in directory "
+// << runTime.constant()/dir
+// << exit(FatalError);
+//
+// return runTime.constant();
+//}
+
+
//- Check file existence
const Foam::fileName& Foam::triSurfaceMesh::checkFile
(
@@ -149,7 +207,19 @@ bool Foam::triSurfaceMesh::isSurfaceClosed() const
Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
:
searchableSurface(io),
- objectRegistry(io),
+ objectRegistry
+ (
+ IOobject
+ (
+ io.name(),
+ io.instance(),
+ io.local(),
+ io.db(),
+ io.readOpt(),
+ io.writeOpt(),
+ false // searchableSurface already registered under name
+ )
+ ),
triSurface(s),
surfaceClosed_(-1)
{}
@@ -157,8 +227,34 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io)
:
+ // Find instance for triSurfaceMesh
searchableSurface(io),
- objectRegistry(io),
+ //(
+ // IOobject
+ // (
+ // io.name(),
+ // io.time().findInstance(io.local(), word::null),
+ // io.local(),
+ // io.db(),
+ // io.readOpt(),
+ // io.writeOpt(),
+ // io.registerObject()
+ // )
+ //),
+ // Reused found instance in objectRegistry
+ objectRegistry
+ (
+ IOobject
+ (
+ io.name(),
+ static_cast<const searchableSurface&>(*this).instance(),
+ io.local(),
+ io.db(),
+ io.readOpt(),
+ io.writeOpt(),
+ false // searchableSurface already registered under name
+ )
+ ),
triSurface
(
checkFile
@@ -178,7 +274,32 @@ Foam::triSurfaceMesh::triSurfaceMesh
)
:
searchableSurface(io),
- objectRegistry(io),
+ //(
+ // IOobject
+ // (
+ // io.name(),
+ // io.time().findInstance(io.local(), word::null),
+ // io.local(),
+ // io.db(),
+ // io.readOpt(),
+ // io.writeOpt(),
+ // io.registerObject()
+ // )
+ //),
+ // Reused found instance in objectRegistry
+ objectRegistry
+ (
+ IOobject
+ (
+ io.name(),
+ static_cast<const searchableSurface&>(*this).instance(),
+ io.local(),
+ io.db(),
+ io.readOpt(),
+ io.writeOpt(),
+ false // searchableSurface already registered under name
+ )
+ ),
triSurface
(
checkFile
diff --git a/src/meshTools/searchableSurface/triSurfaceMesh.H b/src/meshTools/searchableSurface/triSurfaceMesh.H
index aa0af5979c56e582ba96beb8186aec90b4840d6b..a4231457d5183097cff2eceb31ca00e944c6a05d 100644
--- a/src/meshTools/searchableSurface/triSurfaceMesh.H
+++ b/src/meshTools/searchableSurface/triSurfaceMesh.H
@@ -76,6 +76,14 @@ private:
// Private Member Functions
+ ////- Helper: find instance of files without header
+ //static word findRawInstance
+ //(
+ // const Time&,
+ // const fileName&,
+ // const word&
+ //);
+
//- Check file existence
static const fileName& checkFile
(
@@ -105,10 +113,10 @@ public:
//- Construct from triSurface
triSurfaceMesh(const IOobject&, const triSurface&);
- //- Construct read
+ //- Construct read.
triSurfaceMesh(const IOobject& io);
- //- Construct from dictionary (used by searchableSurface)
+ //- Construct from IO and dictionary (used by searchableSurface).
// Dictionary may contain a 'scale' entry (eg, 0.001: mm -> m)
triSurfaceMesh
(
diff --git a/src/sampling/sampledSurface/isoSurface/isoSurface.C b/src/sampling/sampledSurface/isoSurface/isoSurface.C
index c69a51f40731f18bc1c3aaa377aed705f84ce109..78081050ab4e715d3ec6f4378e3f91e1e91dc8da 100644
--- a/src/sampling/sampledSurface/isoSurface/isoSurface.C
+++ b/src/sampling/sampledSurface/isoSurface/isoSurface.C
@@ -30,6 +30,7 @@ License
#include "syncTools.H"
#include "addToRunTimeSelectionTable.H"
#include "slicedVolFields.H"
+#include "surfaceFields.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -84,9 +85,74 @@ bool Foam::isoSurface::isEdgeOfFaceCut
}
+// Get neighbour value and position.
+void Foam::isoSurface::getNeighbour
+(
+ const labelList& boundaryRegion,
+ const volScalarField& cVals,
+ const label cellI,
+ const label faceI,
+ scalar& nbrValue,
+ point& nbrPoint
+) const
+{
+ const labelList& own = mesh_.faceOwner();
+ const labelList& nei = mesh_.faceNeighbour();
+ const surfaceScalarField& weights = mesh_.weights();
+
+ if (mesh_.isInternalFace(faceI))
+ {
+ label nbr = (own[faceI] == cellI ? nei[faceI] : own[faceI]);
+ nbrValue = cVals[nbr];
+ nbrPoint = mesh_.cellCentres()[nbr];
+ }
+ else
+ {
+ label bFaceI = faceI-mesh_.nInternalFaces();
+ label patchI = boundaryRegion[bFaceI];
+ const polyPatch& pp = mesh_.boundaryMesh()[patchI];
+ label patchFaceI = faceI-pp.start();
+
+ if (isA<emptyPolyPatch>(pp))
+ {
+ // Assume zero gradient
+ nbrValue = cVals[own[faceI]];
+ nbrPoint = mesh_.faceCentres()[faceI];
+ }
+ else if (pp.coupled())
+ {
+ if (!refCast<const coupledPolyPatch>(pp).separated())
+ {
+ // Behave as internal face:
+ // other side value
+ nbrValue = cVals.boundaryField()[patchI][patchFaceI];
+ // other side cell centre
+ nbrPoint = mesh_.C().boundaryField()[patchI][patchFaceI];
+ }
+ else
+ {
+ // Do some interpolation for now
+ const scalarField& w = weights.boundaryField()[patchI];
+
+ nbrPoint = mesh_.faceCentres()[faceI];
+ nbrValue =
+ (1-w[patchFaceI])*cVals[own[faceI]]
+ + w[patchFaceI]*cVals.boundaryField()[patchI][patchFaceI];
+ }
+ }
+ else
+ {
+ nbrValue = cVals.boundaryField()[patchI][patchFaceI];
+ nbrPoint = mesh_.faceCentres()[faceI];
+ }
+ }
+}
+
+
// Determine for every face/cell whether it (possibly) generates triangles.
void Foam::isoSurface::calcCutTypes
(
+ const labelList& boundaryRegion,
const volScalarField& cVals,
const scalarField& pVals
)
@@ -102,7 +168,20 @@ void Foam::isoSurface::calcCutTypes
{
// CC edge.
bool ownLower = (cVals[own[faceI]] < iso_);
- bool neiLower = (cVals[nei[faceI]] < iso_);
+
+ scalar nbrValue;
+ point nbrPoint;
+ getNeighbour
+ (
+ boundaryRegion,
+ cVals,
+ own[faceI],
+ faceI,
+ nbrValue,
+ nbrPoint
+ );
+
+ bool neiLower = (nbrValue < iso_);
if (ownLower != neiLower)
{
@@ -120,53 +199,47 @@ void Foam::isoSurface::calcCutTypes
}
}
}
+
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
label faceI = pp.start();
- if (isA<emptyPolyPatch>(pp))
+ forAll(pp, i)
{
- // Assume zero gradient so owner and neighbour/boundary value equal
-
- forAll(pp, i)
- {
- bool ownLower = (cVals[own[faceI]] < iso_);
+ bool ownLower = (cVals[own[faceI]] < iso_);
- const face f = mesh_.faces()[faceI];
+ scalar nbrValue;
+ point nbrPoint;
+ getNeighbour
+ (
+ boundaryRegion,
+ cVals,
+ own[faceI],
+ faceI,
+ nbrValue,
+ nbrPoint
+ );
- if (isEdgeOfFaceCut(pVals, f, ownLower, ownLower))
- {
- faceCutType_[faceI] = CUT;
- }
+ bool neiLower = (nbrValue < iso_);
- faceI++;
+ if (ownLower != neiLower)
+ {
+ faceCutType_[faceI] = CUT;
}
- }
- else
- {
- forAll(pp, i)
+ else
{
- bool ownLower = (cVals[own[faceI]] < iso_);
- bool neiLower = (cVals.boundaryField()[patchI][i] < iso_);
+ // Mesh edge.
+ const face f = mesh_.faces()[faceI];
- if (ownLower != neiLower)
+ if (isEdgeOfFaceCut(pVals, f, ownLower, neiLower))
{
faceCutType_[faceI] = CUT;
}
- else
- {
- // Mesh edge.
- const face f = mesh_.faces()[faceI];
-
- if (isEdgeOfFaceCut(pVals, f, ownLower, neiLower))
- {
- faceCutType_[faceI] = CUT;
- }
- }
- faceI++;
}
+
+ faceI++;
}
}
@@ -330,59 +403,6 @@ Foam::pointIndexHit Foam::isoSurface::collapseSurface
}
-// Get neighbour value and position.
-void Foam::isoSurface::getNeighbour
-(
- const labelList& boundaryRegion,
- const volScalarField& cVals,
- const label cellI,
- const label faceI,
- scalar& nbrValue,
- point& nbrPoint
-) const
-{
- const labelList& own = mesh_.faceOwner();
- const labelList& nei = mesh_.faceNeighbour();
-
- if (mesh_.isInternalFace(faceI))
- {
- label nbr = (own[faceI] == cellI ? nei[faceI] : own[faceI]);
- nbrValue = cVals[nbr];
- nbrPoint = mesh_.cellCentres()[nbr];
- }
- else
- {
- label bFaceI = faceI-mesh_.nInternalFaces();
- label patchI = boundaryRegion[bFaceI];
- const polyPatch& pp = mesh_.boundaryMesh()[patchI];
- label patchFaceI = faceI-pp.start();
-
- if (isA<emptyPolyPatch>(pp))
- {
- // Assume zero gradient
- nbrValue = cVals[own[faceI]];
- nbrPoint = mesh_.faceCentres()[faceI];
- }
- else if
- (
- pp.coupled()
- && !refCast<const coupledPolyPatch>(pp).separated()
- )
- {
- // other side value
- nbrValue = cVals.boundaryField()[patchI][patchFaceI];
- // other side cell centre
- nbrPoint = mesh_.C().boundaryField()[patchI][patchFaceI];
- }
- else
- {
- nbrValue = cVals.boundaryField()[patchI][patchFaceI];
- nbrPoint = mesh_.faceCentres()[faceI];
- }
- }
-}
-
-
// Determine per cell centre whether all the intersections get collapsed
// to a single point
void Foam::isoSurface::calcSnappedCc
@@ -606,13 +626,14 @@ void Foam::isoSurface::calcSnappedPoint
forAll(pFaces, pFaceI)
{
+ // Create points for all intersections close to point
+ // (i.e. from pyramid edges)
+
label faceI = pFaces[pFaceI];
const face& f = mesh_.faces()[faceI];
label own = mesh_.faceOwner()[faceI];
- // Create points for all intersections close to point
- // (i.e. from pyramid edges)
-
+ // Get neighbour value
scalar nbrValue;
point nbrPoint;
getNeighbour
@@ -831,6 +852,7 @@ Foam::triSurface Foam::isoSurface::stitchTriPoints
}
+
// Determine 'flat hole' situation (see RMT paper).
// Two unconnected triangles get connected because (some of) the edges
// separating them get collapsed. Below only checks for duplicate triangles,
@@ -846,22 +868,29 @@ Foam::triSurface Foam::isoSurface::stitchTriPoints
forAll(tris, triI)
{
const labelledTri& tri = tris[triI];
-
const labelList& pFaces = pointFaces[tri[0]];
- // Find the minimum of any duplicates
+ // Find the maximum of any duplicates. Maximum since the tris
+ // below triI
+ // get overwritten so we cannot use them in a comparison.
label dupTriI = -1;
forAll(pFaces, i)
{
- if (pFaces[i] < triI && tris[pFaces[i]] == tri)
+ label nbrTriI = pFaces[i];
+
+ if (nbrTriI > triI && (tris[nbrTriI] == tri))
{
- dupTriI = pFaces[i];
+ //Pout<< "Duplicate : " << triI << " verts:" << tri
+ // << " to " << nbrTriI << " verts:" << tris[nbrTriI]
+ // << endl;
+ dupTriI = nbrTriI;
+ break;
}
}
if (dupTriI == -1)
{
- // There is no lower triangle
+ // There is no (higher numbered) duplicate triangle
label newTriI = newToOldTri.size();
newToOldTri.append(triI);
tris[newTriI] = tris[triI];
@@ -876,6 +905,43 @@ Foam::triSurface Foam::isoSurface::stitchTriPoints
Pout<< "isoSurface : merged from " << nTris
<< " down to " << tris.size() << " unique triangles." << endl;
}
+
+
+ if (debug)
+ {
+ triSurface surf(tris, geometricSurfacePatchList(0), newPoints);
+
+ forAll(surf, faceI)
+ {
+ const labelledTri& f = surf[faceI];
+ const labelList& fFaces = surf.faceFaces()[faceI];
+
+ forAll(fFaces, i)
+ {
+ label nbrFaceI = fFaces[i];
+
+ if (nbrFaceI <= faceI)
+ {
+ // lower numbered faces already checked
+ continue;
+ }
+
+ const labelledTri& nbrF = surf[nbrFaceI];
+
+ if (f == nbrF)
+ {
+ FatalErrorIn("validTri(const triSurface&, const label)")
+ << "Check : "
+ << " triangle " << faceI << " vertices " << f
+ << " fc:" << f.centre(surf.points())
+ << " has the same vertices as triangle " << nbrFaceI
+ << " vertices " << nbrF
+ << " fc:" << nbrF.centre(surf.points())
+ << abort(FatalError);
+ }
+ }
+ }
+ }
}
return triSurface(tris, geometricSurfacePatchList(0), newPoints, true);
@@ -1473,45 +1539,23 @@ Foam::isoSurface::isoSurface
{
if (debug)
{
- Pout<< "isoSurface :" << nl
- << " isoField : " << cVals.name() << nl
- << " isoValue : " << iso << nl
- << " regularise : " << regularise << nl
- << " mergeTol : " << mergeTol << nl
+ Pout<< "isoSurface:" << nl
+ << " isoField : " << cVals.name() << nl
+ << " cell min/max : "
+ << min(cVals.internalField()) << " / "
+ << max(cVals.internalField()) << nl
+ << " point min/max : "
+ << min(pVals) << " / "
+ << max(pVals) << nl
+ << " isoValue : " << iso << nl
+ << " regularise : " << regularise << nl
+ << " mergeTol : " << mergeTol << nl
<< endl;
}
const polyBoundaryMesh& patches = mesh_.boundaryMesh();
- // Check
- forAll(patches, patchI)
- {
- if (isA<emptyPolyPatch>(patches[patchI]))
- {
- FatalErrorIn
- (
- "isoSurface::isoSurface\n"
- "(\n"
- " const volScalarField& cVals,\n"
- " const scalarField& pVals,\n"
- " const scalar iso,\n"
- " const bool regularise,\n"
- " const scalar mergeTol\n"
- ")\n"
- ) << "Iso surfaces not supported on case with empty patches."
- << exit(FatalError);
- }
- }
-
-
- // Determine if any cut through face/cell
- calcCutTypes(cVals, pVals);
-
-
- // Determine if point is on boundary. Points on boundaries are never
- // snapped. Coupled boundaries are handled explicitly so not marked here.
- PackedBoolList isBoundaryPoint(mesh_.nPoints());
-
+ // Pre-calculate patch-per-face to avoid whichPatch call.
labelList boundaryRegion(mesh_.nFaces()-mesh_.nInternalFaces());
forAll(patches, patchI)
@@ -1525,31 +1569,11 @@ Foam::isoSurface::isoSurface
boundaryRegion[faceI-mesh_.nInternalFaces()] = patchI;
faceI++;
}
+ }
- // Mark all points that are not physically coupled (so anything
- // but collocated coupled patches)
- if
- (
- !pp.coupled()
- || refCast<const coupledPolyPatch>(pp).separated()
- )
- {
- label faceI = pp.start();
-
- forAll(pp, i)
- {
- boundaryRegion[faceI-mesh_.nInternalFaces()] = patchI;
-
- const face& f = mesh_.faces()[faceI];
- forAll(f, fp)
- {
- isBoundaryPoint.set(f[fp], 1);
- }
- faceI++;
- }
- }
- }
+ // Determine if any cut through face/cell
+ calcCutTypes(boundaryRegion, cVals, pVals);
DynamicList<point> snappedPoints(nCutCells_);
@@ -1584,6 +1608,39 @@ Foam::isoSurface::isoSurface
label nCellSnaps = snappedPoints.size();
+
+ // Determine if point is on boundary. Points on boundaries are never
+ // snapped. Coupled boundaries are handled explicitly so not marked here.
+ PackedBoolList isBoundaryPoint(mesh_.nPoints());
+
+ forAll(patches, patchI)
+ {
+ // Mark all boundary points that are not physically coupled (so anything
+ // but collocated coupled patches)
+ const polyPatch& pp = patches[patchI];
+
+ if
+ (
+ !pp.coupled()
+ || refCast<const coupledPolyPatch>(pp).separated()
+ )
+ {
+ label faceI = pp.start();
+
+ forAll(pp, i)
+ {
+ const face& f = mesh_.faces()[faceI];
+
+ forAll(f, fp)
+ {
+ isBoundaryPoint.set(f[fp], 1);
+ }
+ faceI++;
+ }
+ }
+ }
+
+
// Per point -1 or a point inside snappedPoints.
labelList snappedPoint;
if (regularise)
@@ -1613,7 +1670,7 @@ Foam::isoSurface::isoSurface
// Generate field to interpolate. This is identical to the mesh.C()
- // except on separated coupled patches.
+ // except on separated coupled patches and on empty patches.
slicedVolVectorField meshC
(
IOobject
@@ -1635,10 +1692,12 @@ Foam::isoSurface::isoSurface
const polyBoundaryMesh& patches = mesh_.boundaryMesh();
forAll(patches, patchI)
{
+ const polyPatch& pp = patches[patchI];
+
if
(
- patches[patchI].coupled()
- && refCast<const coupledPolyPatch>(patches[patchI]).separated()
+ pp.coupled()
+ && refCast<const coupledPolyPatch>(pp).separated()
)
{
fvPatchVectorField& pfld = const_cast<fvPatchVectorField&>
@@ -1647,8 +1706,31 @@ Foam::isoSurface::isoSurface
);
pfld.operator==
(
- patches[patchI].patchSlice(mesh_.faceCentres())
+ pp.patchSlice(mesh_.faceCentres())
+ );
+ }
+ else if (isA<emptyPolyPatch>(pp))
+ {
+ typedef slicedVolVectorField::GeometricBoundaryField bType;
+
+ bType& bfld = const_cast<bType&>(meshC.boundaryField());
+
+ // Clear old value. Cannot resize it since slice.
+ bfld.set(patchI, NULL);
+
+ // Set new value we can change
+ bfld.set
+ (
+ patchI,
+ new calculatedFvPatchField<vector>
+ (
+ mesh_.boundary()[patchI],
+ meshC
+ )
);
+
+ // Change to face centres
+ bfld[patchI] = pp.patchSlice(mesh_.faceCentres());
}
}
}
@@ -1676,7 +1758,8 @@ Foam::isoSurface::isoSurface
if (debug)
{
Pout<< "isoSurface : generated " << triMeshCells.size()
- << " unmerged triangles." << endl;
+ << " unmerged triangles from " << triPoints.size()
+ << " unmerged points." << endl;
}
@@ -1773,6 +1856,14 @@ Foam::isoSurface::isoSurface
orientSurface(*this, faceEdges, edgeFace0, edgeFace1, edgeFacesRest);
//}
+
+
+ if (debug)
+ {
+ fileName stlFile = mesh_.time().path() + ".stl";
+ Pout<< "Dumping surface to " << stlFile << endl;
+ triSurface::write(stlFile);
+ }
}
diff --git a/src/sampling/sampledSurface/isoSurface/isoSurface.H b/src/sampling/sampledSurface/isoSurface/isoSurface.H
index 7dc0588a1c0d0ccb448f5e00e693d09de2dbe721..b859524b2cc40f21aad1351b3dfa2d6f4d798a88 100644
--- a/src/sampling/sampledSurface/isoSurface/isoSurface.H
+++ b/src/sampling/sampledSurface/isoSurface/isoSurface.H
@@ -31,9 +31,6 @@ Description
G.M. Treece, R.W. Prager and A.H. Gee.
Note:
- - not possible on patches of type 'empty'. There are no values on
- 'empty' patch fields so even the api would have to change
- (no volScalarField as argument). Too messy.
- in parallel the regularisation (coarsening) always takes place
and slightly different surfaces will be created compared to non-parallel.
The surface will still be continuous though!
@@ -43,6 +40,7 @@ Description
- uses geometric merge with fraction of bounding box as distance.
- triangles can be between two cell centres so constant sampling
does not make sense.
+ - on empty patches behaves like zero gradient.
- does not do 2D correctly, creates non-flat iso surface.
SourceFiles
@@ -135,9 +133,20 @@ class isoSurface
const bool neiLower
) const;
+ void getNeighbour
+ (
+ const labelList& boundaryRegion,
+ const volScalarField& cVals,
+ const label cellI,
+ const label faceI,
+ scalar& nbrValue,
+ point& nbrPoint
+ ) const;
+
//- Set faceCutType,cellCutType.
void calcCutTypes
(
+ const labelList& boundaryRegion,
const volScalarField& cVals,
const scalarField& pVals
);
@@ -156,16 +165,6 @@ class isoSurface
DynamicList<labelledTri, 64>& localTris
);
- void getNeighbour
- (
- const labelList& boundaryRegion,
- const volScalarField& cVals,
- const label cellI,
- const label faceI,
- scalar& nbrValue,
- point& nbrPoint
- ) const;
-
//- Determine per cc whether all near cuts can be snapped to single
// point.
void calcSnappedCc
@@ -193,37 +192,39 @@ class isoSurface
template<class Type>
Type generatePoint
(
- const DynamicList<Type>& snappedPoints,
-
const scalar s0,
const Type& p0,
- const label p0Index,
+ const bool hasSnap0,
+ const Type& snapP0,
const scalar s1,
const Type& p1,
- const label p1Index
+ const bool hasSnap1,
+ const Type& snapP1
) const;
template<class Type>
void generateTriPoints
(
- const DynamicList<Type>& snapped,
-
const scalar s0,
const Type& p0,
- const label p0Index,
+ const bool hasSnap0,
+ const Type& snapP0,
const scalar s1,
const Type& p1,
- const label p1Index,
+ const bool hasSnap1,
+ const Type& snapP1,
const scalar s2,
const Type& p2,
- const label p2Index,
+ const bool hasSnap2,
+ const Type& snapP2,
const scalar s3,
const Type& p3,
- const label p3Index,
+ const bool hasSnap3,
+ const Type& snapP3,
DynamicList<Type>& points
) const;
@@ -244,7 +245,8 @@ class isoSurface
const scalar neiVal,
const Type& neiPt,
- const label neiSnap,
+ const bool hasNeiSnap,
+ const Type& neiSnapPt,
DynamicList<Type>& triPoints,
DynamicList<label>& triMeshCells
diff --git a/src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C b/src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C
index e2251aad2abdecd824c6d3afc61a93ec7d222432..0760b06a52561de62d50ceb68b00f79b4b016fe7 100644
--- a/src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C
+++ b/src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C
@@ -33,15 +33,15 @@ License
template<class Type>
Type Foam::isoSurface::generatePoint
(
- const DynamicList<Type>& snappedPoints,
-
const scalar s0,
const Type& p0,
- const label p0Index,
+ const bool hasSnap0,
+ const Type& snapP0,
const scalar s1,
const Type& p1,
- const label p1Index
+ const bool hasSnap1,
+ const Type& snapP1
) const
{
scalar d = s1-s0;
@@ -50,13 +50,13 @@ Type Foam::isoSurface::generatePoint
{
scalar s = (iso_-s0)/d;
- if (s >= 0.5 && s <= 1 && p1Index != -1)
+ if (hasSnap1 && s >= 0.5 && s <= 1)
{
- return snappedPoints[p1Index];
+ return snapP1;
}
- else if (s >= 0.0 && s <= 0.5 && p0Index != -1)
+ else if (hasSnap0 && s >= 0.0 && s <= 0.5)
{
- return snappedPoints[p0Index];
+ return snapP0;
}
else
{
@@ -75,23 +75,25 @@ Type Foam::isoSurface::generatePoint
template<class Type>
void Foam::isoSurface::generateTriPoints
(
- const DynamicList<Type>& snapped,
-
const scalar s0,
const Type& p0,
- const label p0Index,
+ const bool hasSnap0,
+ const Type& snapP0,
const scalar s1,
const Type& p1,
- const label p1Index,
+ const bool hasSnap1,
+ const Type& snapP1,
const scalar s2,
const Type& p2,
- const label p2Index,
+ const bool hasSnap2,
+ const Type& snapP2,
const scalar s3,
const Type& p3,
- const label p3Index,
+ const bool hasSnap3,
+ const Type& snapP3,
DynamicList<Type>& points
) const
@@ -123,29 +125,55 @@ void Foam::isoSurface::generateTriPoints
case 0x0E:
case 0x01:
- points.append(generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index));
- points.append(generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index));
- points.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
+ points.append
+ (
+ generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1)
+ );
+ points.append
+ (
+ generatePoint(s0,p0,hasSnap0,snapP0,s2,p2,hasSnap2,snapP2)
+ );
+ points.append
+ (
+ generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
+ );
break;
case 0x0D:
case 0x02:
- points.append(generatePoint(snapped,s1,p1,p1Index,s0,p0,p0Index));
- points.append(generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index));
- points.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
+ points.append
+ (
+ generatePoint(s1,p1,hasSnap1,snapP1,s0,p0,hasSnap0,snapP0)
+ );
+ points.append
+ (
+ generatePoint(s1,p1,hasSnap1,snapP1,s3,p3,hasSnap3,snapP3)
+ );
+ points.append
+ (
+ generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
+ );
break;
case 0x0C:
case 0x03:
{
- Type tp1 = generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index);
- Type tp2 = generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index);
+ Type tp1 =
+ generatePoint(s0,p0,hasSnap0,snapP0,s2,p2,hasSnap2,snapP2);
+ Type tp2 =
+ generatePoint(s1,p1,hasSnap1,snapP1,s3,p3,hasSnap3,snapP3);
- points.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
+ points.append
+ (
+ generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
+ );
points.append(tp1);
points.append(tp2);
points.append(tp2);
- points.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
+ points.append
+ (
+ generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
+ );
points.append(tp1);
}
break;
@@ -153,23 +181,40 @@ void Foam::isoSurface::generateTriPoints
case 0x0B:
case 0x04:
{
- points.append(generatePoint(snapped,s2,p2,p2Index,s0,p0,p0Index));
- points.append(generatePoint(snapped,s2,p2,p2Index,s1,p1,p1Index));
- points.append(generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index));
+ points.append
+ (
+ generatePoint(s2,p2,hasSnap2,snapP2,s0,p0,hasSnap0,snapP0)
+ );
+ points.append
+ (
+ generatePoint(s2,p2,hasSnap2,snapP2,s1,p1,hasSnap1,snapP1)
+ );
+ points.append
+ (
+ generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3)
+ );
}
break;
case 0x0A:
case 0x05:
{
- Type tp0 = generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
- Type tp1 = generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
+ Type tp0 =
+ generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1);
+ Type tp1 =
+ generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3);
points.append(tp0);
points.append(tp1);
- points.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
+ points.append
+ (
+ generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
+ );
points.append(tp0);
- points.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
+ points.append
+ (
+ generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
+ );
points.append(tp1);
}
break;
@@ -177,23 +222,40 @@ void Foam::isoSurface::generateTriPoints
case 0x09:
case 0x06:
{
- Type tp0 = generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
- Type tp1 = generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
+ Type tp0 =
+ generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1);
+ Type tp1 =
+ generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3);
points.append(tp0);
- points.append(generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index));
+ points.append
+ (
+ generatePoint(s1,p1,hasSnap1,snapP1,s3,p3,hasSnap3,snapP3)
+ );
points.append(tp1);
points.append(tp0);
- points.append(generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index));
+ points.append
+ (
+ generatePoint(s0,p0,hasSnap0,snapP0,s2,p2,hasSnap2,snapP2)
+ );
points.append(tp1);
}
break;
case 0x07:
case 0x08:
- points.append(generatePoint(snapped,s3,p3,p3Index,s0,p0,p0Index));
- points.append(generatePoint(snapped,s3,p3,p3Index,s2,p2,p2Index));
- points.append(generatePoint(snapped,s3,p3,p3Index,s1,p1,p1Index));
+ points.append
+ (
+ generatePoint(s3,p3,hasSnap3,snapP3,s0,p0,hasSnap0,snapP0)
+ );
+ points.append
+ (
+ generatePoint(s3,p3,hasSnap3,snapP3,s2,p2,hasSnap2,snapP2)
+ );
+ points.append
+ (
+ generatePoint(s3,p3,hasSnap3,snapP3,s1,p1,hasSnap1,snapP1)
+ );
break;
}
}
@@ -215,7 +277,8 @@ Foam::label Foam::isoSurface::generateFaceTriPoints
const scalar neiVal,
const Type& neiPt,
- const label neiSnap,
+ const bool hasNeiSnap,
+ const Type& neiSnapPt,
DynamicList<Type>& triPoints,
DynamicList<label>& triMeshCells
@@ -234,23 +297,37 @@ Foam::label Foam::isoSurface::generateFaceTriPoints
generateTriPoints
(
- snappedPoints,
-
pVals[pointI],
pCoords[pointI],
- snappedPoint[pointI],
+ snappedPoint[pointI] != -1,
+ (
+ snappedPoint[pointI] != -1
+ ? snappedPoints[snappedPoint[pointI]]
+ : pTraits<Type>::zero
+ ),
pVals[nextPointI],
pCoords[nextPointI],
- snappedPoint[nextPointI],
+ snappedPoint[nextPointI] != -1,
+ (
+ snappedPoint[nextPointI] != -1
+ ? snappedPoints[snappedPoint[nextPointI]]
+ : pTraits<Type>::zero
+ ),
cVals[own],
cCoords[own],
- snappedCc[own],
+ snappedCc[own] != -1,
+ (
+ snappedCc[own] != -1
+ ? snappedPoints[snappedCc[own]]
+ : pTraits<Type>::zero
+ ),
neiVal,
neiPt,
- neiSnap,
+ hasNeiSnap,
+ neiSnapPt,
triPoints
);
@@ -311,25 +388,6 @@ void Foam::isoSurface::generateTriPoints
<< abort(FatalError);
}
- // Determine neighbouring snap status
- labelList neiSnappedCc(mesh_.nFaces()-mesh_.nInternalFaces(), -1);
- forAll(patches, patchI)
- {
- const polyPatch& pp = patches[patchI];
-
- if (pp.coupled())
- {
- label faceI = pp.start();
- forAll(pp, i)
- {
- neiSnappedCc[faceI-mesh_.nInternalFaces()] =
- snappedCc[own[faceI]];
- faceI++;
- }
- }
- }
- syncTools::swapBoundaryFaceList(mesh_, neiSnappedCc, false);
-
// Generate triangle points
@@ -356,7 +414,12 @@ void Foam::isoSurface::generateTriPoints
cVals[nei[faceI]],
cCoords[nei[faceI]],
- snappedCc[nei[faceI]],
+ snappedCc[nei[faceI]] != -1,
+ (
+ snappedCc[nei[faceI]] != -1
+ ? snappedPoints[snappedCc[nei[faceI]]]
+ : pTraits<Type>::zero
+ ),
triPoints,
triMeshCells
@@ -365,52 +428,139 @@ void Foam::isoSurface::generateTriPoints
}
+ // Determine neighbouring snap status
+ boolList neiSnapped(mesh_.nFaces()-mesh_.nInternalFaces(), false);
+ List<Type> neiSnappedPoint(neiSnapped.size(), pTraits<Type>::zero);
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
- if
- (
- isA<processorPolyPatch>(pp)
- && refCast<const processorPolyPatch>(pp).owner()
- && !refCast<const processorPolyPatch>(pp).separated()
- )
+ if (pp.coupled())
{
label faceI = pp.start();
-
forAll(pp, i)
{
- if (faceCutType_[faceI] != NOTCUT)
+ label bFaceI = faceI-mesh_.nInternalFaces();
+ label snappedIndex = snappedCc[own[faceI]];
+
+ if (snappedIndex != -1)
{
- generateFaceTriPoints
- (
- cVals,
- pVals,
+ neiSnapped[bFaceI] = true;
+ neiSnappedPoint[bFaceI] = snappedPoints[snappedIndex];
+ }
+ faceI++;
+ }
+ }
+ }
+ syncTools::swapBoundaryFaceList(mesh_, neiSnapped, false);
+ syncTools::swapBoundaryFaceList(mesh_, neiSnappedPoint, false);
- cCoords,
- pCoords,
- snappedPoints,
- snappedCc,
- snappedPoint,
- faceI,
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
- cVals.boundaryField()[patchI][i],
- cCoords.boundaryField()[patchI][i],
- neiSnappedCc[faceI-mesh_.nInternalFaces()],
+ if
+ (
+ isA<processorPolyPatch>(pp)
+ && !refCast<const processorPolyPatch>(pp).separated()
+ )
+ {
+ //if (refCast<const processorPolyPatch>(pp).owner())
+ {
+ label faceI = pp.start();
- triPoints,
- triMeshCells
- );
+ forAll(pp, i)
+ {
+ if (faceCutType_[faceI] != NOTCUT)
+ {
+ label bFaceI = faceI-mesh_.nInternalFaces();
+ if
+ (
+ neiSnapped[bFaceI]
+ && (neiSnappedPoint[bFaceI]==pTraits<Type>::zero)
+ )
+ {
+ FatalErrorIn("isoSurface::generateTriPoints(..)")
+ << "problem:" << abort(FatalError);
+ }
+
+
+ generateFaceTriPoints
+ (
+ cVals,
+ pVals,
+
+ cCoords,
+ pCoords,
+
+ snappedPoints,
+ snappedCc,
+ snappedPoint,
+ faceI,
+
+ cVals.boundaryField()[patchI][i],
+ cCoords.boundaryField()[patchI][i],
+ neiSnapped[faceI-mesh_.nInternalFaces()],
+ neiSnappedPoint[faceI-mesh_.nInternalFaces()],
+
+ triPoints,
+ triMeshCells
+ );
+ }
+ faceI++;
}
- faceI++;
}
}
else if (isA<emptyPolyPatch>(pp))
{
- // Assume zero-gradient. But what about coordinates?
+ // Check if field is there (when generating geometry the
+ // empty patches have been rewritten to be the face centres),
+ // otherwise use zero-gradient.
+
label faceI = pp.start();
+
+ const fvPatchScalarField& fvp = cVals.boundaryField()[patchI];
+
+ // Owner value of cVals
+ scalarField internalVals;
+ if (fvp.size() == 0)
+ {
+ internalVals.setSize(pp.size());
+ forAll(pp, i)
+ {
+ internalVals[i] = cVals[own[pp.start()+i]];
+ }
+ }
+ const scalarField& bVals =
+ (
+ fvp.size() > 0
+ ? static_cast<const scalarField&>(fvp)
+ : internalVals
+ );
+
+
+ const fvPatchField<Type>& pc = cCoords.boundaryField()[patchI];
+
+ // Owner value of cCoords
+ Field<Type> internalCoords;
+ if (pc.size() == 0)
+ {
+ internalCoords.setSize(pp.size());
+ forAll(pp, i)
+ {
+ internalCoords[i] = cCoords[own[pp.start()+i]];
+ }
+ }
+ const Field<Type>& bCoords =
+ (
+ pc.size() > 0
+ ? static_cast<const Field<Type>&>(pc)
+ : internalCoords
+ );
+
+
forAll(pp, i)
{
if (faceCutType_[faceI] != NOTCUT)
@@ -428,9 +578,10 @@ void Foam::isoSurface::generateTriPoints
snappedPoint,
faceI,
- cVals[own[faceI]],
- cCoords.boundaryField()[patchI][i],
- -1, // fc not snapped
+ bVals[i],
+ bCoords[i],
+ false, // fc not snapped
+ pTraits<Type>::zero,
triPoints,
triMeshCells
@@ -462,7 +613,8 @@ void Foam::isoSurface::generateTriPoints
cVals.boundaryField()[patchI][i],
cCoords.boundaryField()[patchI][i],
- -1, // fc not snapped
+ false, // fc not snapped
+ pTraits<Type>::zero,
triPoints,
triMeshCells
diff --git a/src/turbulenceModels/LES/Allwmake b/src/turbulenceModels/LES/Allwmake
index a80c71ab57d28206f7585707735d49f8be229d95..d83439a9130c64ad25c955e5059ee79e1206c5d8 100755
--- a/src/turbulenceModels/LES/Allwmake
+++ b/src/turbulenceModels/LES/Allwmake
@@ -2,7 +2,7 @@
cd ${0%/*} || exit 1 # run from this directory
set -x
-wmakeLnInclude -f ../incompressible/LES
+wmakeLnInclude -f ../incompressible/LES -sf
wmake libso LESfilters
wmake libso LESdeltas
diff --git a/src/turbulenceModels/LES/LESdeltas/cubeRootVolDelta/cubeRootVolDelta.C b/src/turbulenceModels/LES/LESdeltas/cubeRootVolDelta/cubeRootVolDelta.C
index a3ce3ab3a99b1b94fe2ecc20e69d947dc3d7501d..2002deaa2181d05711c147f2f22bccc104712631 100644
--- a/src/turbulenceModels/LES/LESdeltas/cubeRootVolDelta/cubeRootVolDelta.C
+++ b/src/turbulenceModels/LES/LESdeltas/cubeRootVolDelta/cubeRootVolDelta.C
@@ -42,8 +42,7 @@ addToRunTimeSelectionTable(LESdelta, cubeRootVolDelta, dictionary);
void cubeRootVolDelta::calcDelta()
{
- const Vector<label>& directions = mesh().directions();
- label nD = (directions.nComponents + cmptSum(directions))/2;
+ label nD = mesh().nGeometricD();
if (nD == 3)
{
@@ -55,14 +54,15 @@ void cubeRootVolDelta::calcDelta()
<< "Case is 2D, LES is not strictly applicable\n"
<< endl;
+ const Vector<label>& directions = mesh().geometricD();
+
scalar thickness = 0.0;
for (direction dir=0; dir<directions.nComponents; dir++)
{
if (directions[dir] == -1)
{
- boundBox bb(mesh().points(), false);
-
- thickness = bb.span()[dir];
+ thickness = mesh().bounds().span()[dir];
+ break;
}
}
diff --git a/src/turbulenceModels/incompressible/LES/SpalartAllmarasIDDES/IDDESDelta/IDDESDelta.C b/src/turbulenceModels/incompressible/LES/SpalartAllmarasIDDES/IDDESDelta/IDDESDelta.C
index 3f3aa7699bb335fdebb01684e02f9b32de037ac4..005ee2eb1610b6f43a58db645123f7096f015b44 100644
--- a/src/turbulenceModels/incompressible/LES/SpalartAllmarasIDDES/IDDESDelta/IDDESDelta.C
+++ b/src/turbulenceModels/incompressible/LES/SpalartAllmarasIDDES/IDDESDelta/IDDESDelta.C
@@ -41,8 +41,7 @@ namespace Foam
void Foam::IDDESDelta::calcDelta()
{
- const Vector<label>& directions = mesh().directions();
- label nD = (directions.nComponents + cmptSum(directions))/2;
+ label nD = mesh().nGeometricD();
// initialise hwn as wall distance
volScalarField hwn = wallDist(mesh()).y();
diff --git a/tutorials/incompressible/MRFSimpleFoam/MRFSimpleFoam/MRFSimpleFoam.C b/tutorials/incompressible/MRFSimpleFoam/MRFSimpleFoam/MRFSimpleFoam.C
index 9fc36826fe69d29b6164386dc74e19eb6b5df6ff..b3d7312a064a37cea473bb0de18e0e6687da3133 100644
--- a/tutorials/incompressible/MRFSimpleFoam/MRFSimpleFoam/MRFSimpleFoam.C
+++ b/tutorials/incompressible/MRFSimpleFoam/MRFSimpleFoam/MRFSimpleFoam.C
@@ -32,8 +32,8 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
-#include "incompressible/singlePhaseTransportModel/singlePhaseTransportModel.H"
-#include "incompressible/RASModel/RASModel.H"
+#include "singlePhaseTransportModel.H"
+#include "RASModel.H"
#include "MRFZones.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/tutorials/incompressible/MRFSimpleFoam/MRFSimpleFoam/Make/options b/tutorials/incompressible/MRFSimpleFoam/MRFSimpleFoam/Make/options
index 2a417f954936f537595fba5fa461a73a81650574..1223bdd06f48178f0b2eee8032d29c936634f328 100644
--- a/tutorials/incompressible/MRFSimpleFoam/MRFSimpleFoam/Make/options
+++ b/tutorials/incompressible/MRFSimpleFoam/MRFSimpleFoam/Make/options
@@ -1,7 +1,9 @@
EXE_INC = \
- -I$(LIB_SRC)/finiteVolume/lnInclude \
- -I$(LIB_SRC)/turbulenceModels/RAS \
- -I$(LIB_SRC)/transportModels
+ -I$(LIB_SRC)/turbulenceModels \
+ -I$(LIB_SRC)/turbulenceModels/incompressible/RAS/RASModel \
+ -I$(LIB_SRC)/transportModels \
+ -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
+ -I$(LIB_SRC)/finiteVolume/lnInclude
EXE_LIBS = \
-lincompressibleRASModels \
diff --git a/tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/makeMesh b/tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/makeMesh
index f45c6ff0be3a4fe281484bbfd4789baf4fe00b32..b0be10d8e5f2b0525cec8f23c93585fde73f5879 100755
--- a/tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/makeMesh
+++ b/tutorials/incompressible/MRFSimpleFoam/mixerVessel2D/makeMesh
@@ -3,8 +3,9 @@
m4 < constant/polyMesh/blockMeshDict.m4 > constant/polyMesh/blockMeshDict
blockMesh
cellSet
-cp system/faceSetDict_rotorFaces system/faceSetDict
-faceSet
-cp system/faceSetDict_noBoundaryFaces system/faceSetDict
-faceSet
+#- MRF determines its own faceZone if not supplied
+#cp system/faceSetDict_rotorFaces system/faceSetDict
+#faceSet
+#cp system/faceSetDict_noBoundaryFaces system/faceSetDict
+#faceSet
setsToZones -noFlipMap
diff --git a/tutorials/incompressible/channelFoam/channel395/constant/postChannelDict b/tutorials/incompressible/channelFoam/channel395/constant/postChannelDict
index 89fcd5080961334007bc68fd3ea0921fcdddf1d9..96bbe30d760a249594044f2dfd8637b878fb1eb2 100644
--- a/tutorials/incompressible/channelFoam/channel395/constant/postChannelDict
+++ b/tutorials/incompressible/channelFoam/channel395/constant/postChannelDict
@@ -15,10 +15,14 @@ FoamFile
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+// Seed patches to start layering from
patches ( bottomWall );
+// Direction in which the layers are
component y;
+// Is the mesh symmetric? If so average(symmetric fields) or
+// subtract(asymmetric) contributions from both halves
symmetric true;
diff --git a/tutorials/incompressible/channelFoam/channel395/system/controlDict b/tutorials/incompressible/channelFoam/channel395/system/controlDict
index 2ca28c29a1150383756dc56581f73663dd918fa9..218a76d8617819d086e88f897203bf152dae747f 100644
--- a/tutorials/incompressible/channelFoam/channel395/system/controlDict
+++ b/tutorials/incompressible/channelFoam/channel395/system/controlDict
@@ -43,7 +43,35 @@ timePrecision 6;
runTimeModifiable yes;
-functions ( fieldAverage1 { type fieldAverage ; functionObjectLibs ( "libfieldFunctionObjects.so" ) ; enabled true ; outputControl outputTime ; fields ( U { mean on ; prime2Mean on ; base time ; } p { mean on ; prime2Mean on ; base time ; } ) ; } );
+functions
+(
+ fieldAverage1
+ {
+ type fieldAverage;
+
+ functionObjectLibs ( "libfieldFunctionObjects.so" );
+
+ enabled true;
+
+ outputControl outputTime;
+
+ fields
+ (
+ U
+ {
+ mean on;
+ prime2Mean on;
+ base time;
+ }
+ p
+ {
+ mean on;
+ prime2Mean on;
+ base time;
+ }
+ );
+ }
+);
// ************************************************************************* //
diff --git a/tutorials/mesh/snappyHexMesh/iglooWithFridges/system/snappyHexMeshDict b/tutorials/mesh/snappyHexMesh/iglooWithFridges/system/snappyHexMeshDict
index 2efe6901aa1657a7bb782cd7991cfc944b713132..19ef2afec1d08e169d49fd7907cc518d0eccd94b 100644
--- a/tutorials/mesh/snappyHexMesh/iglooWithFridges/system/snappyHexMeshDict
+++ b/tutorials/mesh/snappyHexMesh/iglooWithFridges/system/snappyHexMeshDict
@@ -2,7 +2,7 @@
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 1.5 |
-| \\ / A nd | Web: www.OpenFOAM.org |
+| \\ / A nd | Web: http://www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
@@ -10,138 +10,373 @@ FoamFile
version 2.0;
format ascii;
class dictionary;
- location "system";
- object snappyHexMeshDict;
+ object autoHexMeshDict;
}
+
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+// Which of the steps to run
castellatedMesh true;
-
snap true;
-
addLayers true;
+
+// Geometry. Definition of all surfaces. All surfaces are of class
+// searchableSurface.
+// Surfaces are used
+// - to specify refinement for any mesh cell intersecting it
+// - to specify refinement for any mesh cell inside/outside/near
+// - to 'snap' the mesh boundary to the surface
geometry
{
- fridgeA
+ igloo
{
- type searchableBox;
- min ( 2 2 0 );
- max ( 3 3 2 );
+ type searchableSphere;
+ centre (3 3 0);
+ radius 4;
}
- fridgeB
+ box1
{
- type searchableBox;
- min ( 3.5 3 0 );
- max ( 4.3 3.8 1.8 );
+ type searchableBox;
+ min (0 0 0);
+ max (1 1 1);
}
-
- igloo
+ fridgeFreezer
{
- type searchableSphere;
- centre ( 3 3 0 );
- radius 4;
+ type searchableSurfaceCollection;
+
+ freezer
+ {
+ surface box1;
+ scale (1 1 1);
+ transform
+ {
+ type cartesian;
+ origin (0 0 0);
+ e1 (1 0 0);
+ e3 (0 0 1);
+ }
+ }
+ fridge
+ {
+ surface box1;
+ scale (1 1 1.1);
+ transform
+ {
+ type cartesian;
+ origin (0 0 1);
+ e1 (1 0 0);
+ e3 (0 0 1);
+ }
+ }
}
-}
+ twoFridgeFreezers
+ {
+ type searchableSurfaceCollection;
+
+ seal
+ {
+ surface fridgeFreezer;
+ scale (1.0 1.0 1.0);
+ transform
+ {
+ type cartesian;
+ origin (2 2 0);
+ e1 (1 0 0);
+ e3 (0 0 1);
+ }
+ }
+ herring
+ {
+ surface fridgeFreezer;
+ scale (1.0 1.0 1.0);
+ transform
+ {
+ type cartesian;
+ origin (3.5 3 0);
+ e1 (1 0 0);
+ e3 (0 0 1);
+ }
+ }
+ }
+};
+
+
+// Settings for the castellatedMesh generation.
castellatedMeshControls
{
- maxLocalCells 1000000;
- maxGlobalCells 2000000;
- minRefinementCells 0;
+
+ // Refinement parameters
+ // ~~~~~~~~~~~~~~~~~~~~~
+
+ // While refining maximum number of cells per processor. This is basically
+ // the number of cells that fit on a processor. If you choose this too small
+ // it will do just more refinement iterations to obtain a similar mesh.
+ maxLocalCells 1000000;
+
+ // Overall cell limit (approximately). Refinement will stop immediately
+ // upon reaching this number so a refinement level might not complete.
+ // Note that this is the number of cells before removing the part which
+ // is not 'visible' from the keepPoint. The final number of cells might
+ // actually be a lot less.
+ maxGlobalCells 2000000;
+
+ // The surface refinement loop might spend lots of iterations refining just a
+ // few cells. This setting will cause refinement to stop if <= minimumRefine
+ // are selected for refinement. Note: it will at least do one iteration
+ // (unless the number of cells to refine is 0)
+ minRefinementCells 100;
+
+ // Number of buffer layers between different levels.
+ // 1 means normal 2:1 refinement restriction, larger means slower
+ // refinement.
nCellsBetweenLevels 1;
- features ( );
+
+
+
+ // Explicit feature edge refinement
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ // Specifies a level for any cell intersected by its edges.
+ // This is a featureEdgeMesh, read from constant/triSurface for now.
+ features
+ (
+// {
+// file "fridgeA.eMesh";
+// level 3;
+// }
+ );
+
+
+
+ // Surface based refinement
+ // ~~~~~~~~~~~~~~~~~~~~~~~~
+
+ // Specifies two levels for every surface. The first is the minimum level,
+ // every cell intersecting a surface gets refined up to the minimum level.
+ // The second level is the maximum level. Cells that 'see' multiple
+ // intersections where the intersections make an
+ // angle > resolveFeatureAngle get refined up to the maximum level.
+
refinementSurfaces
{
- fridgeA
+ twoFridgeFreezers
{
- level ( 2 2 );
- }
+ // Surface-wise min and max refinement level
+ level (2 2);
- fridgeB
- {
- level ( 2 2 );
+ regions
+ {
+ // Region-wise override
+ "cook.*"
+ {
+ level (3 3);
+ }
+ }
}
- igloo
+ "iglo.*"
{
- level ( 1 1 );
+ // Surface-wise min and max refinement level
+ level (1 1);
}
}
+ // Resolve sharp angles on fridges
resolveFeatureAngle 60;
+
+
+ // Region-wise refinement
+ // ~~~~~~~~~~~~~~~~~~~~~~
+
+ // Specifies refinement level for cells in relation to a surface. One of
+ // three modes
+ // - distance. 'levels' specifies per distance to the surface the
+ // wanted refinement level. The distances need to be specified in
+ // descending order.
+ // - inside. 'levels' is only one entry and only the level is used. All
+ // cells inside the surface get refined up to the level. The surface
+ // needs to be closed for this to be possible.
+ // - outside. Same but cells outside.
+
refinementRegions
{
}
- locationInMesh ( 3 0.28 0.43 );
+
+ // Mesh selection
+ // ~~~~~~~~~~~~~~
+
+ // After refinement patches get added for all refinementSurfaces and
+ // all cells intersecting the surfaces get put into these patches. The
+ // section reachable from the locationInMesh is kept.
+ // NOTE: This point should never be on a face, always inside a cell, even
+ // after refinement.
+ locationInMesh (3 0.28 0.43);
}
+
+
+// Settings for the snapping.
snapControls
{
- nSmoothPatch 3;
- tolerance 4;
- nSolveIter 30;
- nRelaxIter 5;
+ //- Number of patch smoothing iterations before finding correspondence
+ // to surface
+ nSmoothPatch 3;
+
+ //- Relative distance for points to be attracted by surface feature point
+ // or edge. True distance is this factor times local
+ // maximum edge length.
+ tolerance 4.0;
+
+ //- Number of mesh displacement relaxation iterations.
+ nSolveIter 30;
+
+ //- Maximum number of snapping relaxation iterations. Should stop
+ // before upon reaching a correct mesh.
+ nRelaxIter 5;
}
+
+
+// Settings for the layer addition.
addLayersControls
{
+ // Per final patch (so not geometry!) the layer information
layers
{
- fridgeA_region0
+ "two.*"
{
- nSurfaceLayers 1;
+ nSurfaceLayers 3;
}
-
- fridgeB_region0
- {
- nSurfaceLayers 1;
- }
-
- igloo_region0
+ "igloo_.*"
{
- nSurfaceLayers 1;
+ nSurfaceLayers 1;
}
}
- expansionRatio 1;
+ // Expansion factor for layer mesh
+ expansionRatio 1.0;
+
+ //- Wanted thickness of final added cell layer. If multiple layers
+ // is the
+ // thickness of the layer furthest away from the wall.
+ // Relative to undistorted size of cell outside layer.
finalLayerRatio 0.5;
- minThickness 0.25;
- nGrow 0;
- featureAngle 60;
- nRelaxIter 5;
+
+ //- Minimum thickness of cell layer. If for any reason layer
+ // cannot be above minThickness do not add layer.
+ // Relative to undistorted size of cell outside layer.
+ minThickness 0.25;
+
+ //- If points get not extruded do nGrow layers of connected faces that are
+ // also not grown. This helps convergence of the layer addition process
+ // close to features.
+ nGrow 0;
+
+
+ // Advanced settings
+
+ //- When not to extrude surface. 0 is flat surface, 90 is when two faces
+ // make straight angle.
+ featureAngle 60;
+
+ //- Maximum number of snapping relaxation iterations. Should stop
+ // before upon reaching a correct mesh.
+ nRelaxIter 5;
+
+ // Number of smoothing iterations of surface normals
nSmoothSurfaceNormals 1;
- nSmoothNormals 3;
+
+ // Number of smoothing iterations of interior mesh movement direction
+ nSmoothNormals 3;
+
+ // Smooth layer thickness over surface patches
nSmoothThickness 10;
+
+ // Stop layer growth on highly warped cells
maxFaceThicknessRatio 0.5;
+
+ // Reduce layer growth where ratio thickness to medial
+ // distance is large
maxThicknessToMedialRatio 0.3;
+
+ // Angle used to pick up medial axis points
minMedianAxisAngle 130;
+
+ // Create buffer region for new layer terminations
nBufferCellsNoExtrude 0;
}
+
+
+// Generic mesh quality settings. At any undoable phase these determine
+// where to undo.
meshQualityControls
{
- maxNonOrtho 65;
+ //- Maximum non-orthogonality allowed. Set to 180 to disable.
+ maxNonOrtho 65;
+
+ //- Max skewness allowed. Set to <0 to disable.
maxBoundarySkewness 20;
maxInternalSkewness 4;
- maxConcave 80;
- minFlatness 0.5;
- minVol 1e-13;
- minArea -1;
- minTwist 0.05;
- minDeterminant 0.001;
- minFaceWeight 0.05;
- minVolRatio 0.01;
+
+ //- Max concaveness allowed. Is angle (in degrees) below which concavity
+ // is allowed. 0 is straight face, <0 would be convex face.
+ // Set to 180 to disable.
+ maxConcave 80;
+
+ //- Minimum projected area v.s. actual area. Set to -1 to disable.
+ minFlatness 0.5;
+
+ //- Minimum pyramid volume. Is absolute volume of cell pyramid.
+ // Set to very negative number (e.g. -1E30) to disable.
+ minVol 1e-13;
+
+ //- Minimum face area. Set to <0 to disable.
+ minArea -1;
+
+ //- Minimum face twist. Set to <-1 to disable. dot product of face normal
+ //- and face centre triangles normal
+ minTwist 0.05;
+
+ //- minimum normalised cell determinant
+ //- 1 = hex, <= 0 = folded or flattened illegal cell
+ minDeterminant 0.001;
+
+ //- minFaceWeight (0 -> 0.5)
+ minFaceWeight 0.05;
+
+ //- minVolRatio (0 -> 1)
+ minVolRatio 0.01;
+
+ //must be >0 for Fluent compatibility
minTriangleTwist -1;
- nSmoothScale 4;
- errorReduction 0.75;
+
+
+ // Advanced
+
+ //- Number of error distribution iterations
+ nSmoothScale 4;
+ //- amount to scale back displacement at error points
+ errorReduction 0.75;
}
-debug 0;
-mergeTolerance 1e-06;
+// Advanced
+
+// Flags for optional output
+// 0 : only write final meshes
+// 1 : write intermediate meshes
+// 2 : write volScalarField with cellLevel for postprocessing
+// 4 : write current intersections as .obj files
+debug 0;
+
+// Merge tolerance. Is fraction of overall bounding box of initial mesh.
+// Note: the write tolerance needs to be higher than this.
+mergeTolerance 1E-6;
// ************************************************************************* //
diff --git a/tutorials/mesh/snappyHexMesh/motorBike/system/snappyHexMeshDict b/tutorials/mesh/snappyHexMesh/motorBike/system/snappyHexMeshDict
index fa094e1b03b4fe1b3d419af45a6b757592cfac4c..3f0ff0417e6b92047878cb285aa71452c8fd3c74 100644
--- a/tutorials/mesh/snappyHexMesh/motorBike/system/snappyHexMeshDict
+++ b/tutorials/mesh/snappyHexMesh/motorBike/system/snappyHexMeshDict
@@ -13,14 +13,21 @@ FoamFile
location "system";
object snappyHexMeshDict;
}
+
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+// Which of the steps to run
castellatedMesh true;
-
snap true;
-
addLayers true;
+
+// Geometry. Definition of all surfaces. All surfaces are of class
+// searchableSurface.
+// Surfaces are used
+// - to specify refinement for any mesh cell intersecting it
+// - to specify refinement for any mesh cell inside/outside/near
+// - to 'snap' the mesh boundary to the surface
geometry
{
motorBike.stl
@@ -37,22 +44,89 @@ geometry
}
}
+
+// Settings for the castellatedMesh generation.
castellatedMeshControls
{
+
+ // Refinement parameters
+ // ~~~~~~~~~~~~~~~~~~~~~
+
+ // While refining maximum number of cells per processor. This is basically
+ // the number of cells that fit on a processor. If you choose this too small
+ // it will do just more refinement iterations to obtain a similar mesh.
maxLocalCells 1000000;
+
+
+ // Overall cell limit (approximately). Refinement will stop immediately
+ // upon reaching this number so a refinement level might not complete.
+ // Note that this is the number of cells before removing the part which
+ // is not 'visible' from the keepPoint. The final number of cells might
+ // actually be a lot less.
maxGlobalCells 2000000;
+
+ // The surface refinement loop might spend lots of iterations refining just a
+ // few cells. This setting will cause refinement to stop if <= minimumRefine
+ // are selected for refinement. Note: it will at least do one iteration
+ // (unless the number of cells to refine is 0)
minRefinementCells 10;
+
+ // Number of buffer layers between different levels.
+ // 1 means normal 2:1 refinement restriction, larger means slower
+ // refinement.
nCellsBetweenLevels 2;
- features ( );
+
+
+
+ // Explicit feature edge refinement
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ // Specifies a level for any cell intersected by its edges.
+ // This is a featureEdgeMesh, read from constant/triSurface for now.
+ features
+ (
+ //{
+ // file "someLine.eMesh";
+ // level 2;
+ //}
+ );
+
+
+
+ // Surface based refinement
+ // ~~~~~~~~~~~~~~~~~~~~~~~~
+
+ // Specifies two levels for every surface. The first is the minimum level,
+ // every cell intersecting a surface gets refined up to the minimum level.
+ // The second level is the maximum level. Cells that 'see' multiple
+ // intersections where the intersections make an
+ // angle > resolveFeatureAngle get refined up to the maximum level.
+
refinementSurfaces
{
motorBike
{
+ // Surface-wise min and max refinement level
level ( 5 6 );
}
}
+ // Resolve sharp angles
resolveFeatureAngle 30;
+
+
+ // Region-wise refinement
+ // ~~~~~~~~~~~~~~~~~~~~~~
+
+ // Specifies refinement level for cells in relation to a surface. One of
+ // three modes
+ // - distance. 'levels' specifies per distance to the surface the
+ // wanted refinement level. The distances need to be specified in
+ // descending order.
+ // - inside. 'levels' is only one entry and only the level is used. All
+ // cells inside the surface get refined up to the level. The surface
+ // needs to be closed for this to be possible.
+ // - outside. Same but cells outside.
refinementRegions
{
refinementBox
@@ -62,19 +136,45 @@ castellatedMeshControls
}
}
+
+
+ // Mesh selection
+ // ~~~~~~~~~~~~~~
+
+ // After refinement patches get added for all refinementSurfaces and
+ // all cells intersecting the surfaces get put into these patches. The
+ // section reachable from the locationInMesh is kept.
+ // NOTE: This point should never be on a face, always inside a cell, even
+ // after refinement.
locationInMesh ( 3 3 0.43 );
}
+
+// Settings for the snapping.
snapControls
{
+ //- Number of patch smoothing iterations before finding correspondence
+ // to surface
nSmoothPatch 3;
+
+ //- Relative distance for points to be attracted by surface feature point
+ // or edge. True distance is this factor times local
+ // maximum edge length.
tolerance 4;
+
+ //- Number of mesh displacement relaxation iterations.
nSolveIter 30;
+
+ //- Maximum number of snapping relaxation iterations. Should stop
+ // before upon reaching a correct mesh.
nRelaxIter 5;
}
+
+// Settings for the layer addition.
addLayersControls
{
+ // Per final patch (so not geometry!) the layer information
layers
{
minZ
@@ -418,41 +518,126 @@ addLayersControls
}
}
+ // Expansion factor for layer mesh
expansionRatio 1;
+
+ //- Wanted thickness of final added cell layer. If multiple layers
+ // is the
+ // thickness of the layer furthest away from the wall.
+ // Relative to undistorted size of cell outside layer.
finalLayerRatio 0.3;
+
+ //- Minimum thickness of cell layer. If for any reason layer
+ // cannot be above minThickness do not add layer.
+ // Relative to undistorted size of cell outside layer.
minThickness 0.1;
+
+ //- If points get not extruded do nGrow layers of connected faces that are
+ // also not grown. This helps convergence of the layer addition process
+ // close to features.
nGrow 1;
+
+
+ // Advanced settings
+
+ //- When not to extrude surface. 0 is flat surface, 90 is when two faces
+ // make straight angle.
featureAngle 30;
+
+ //- Maximum number of snapping relaxation iterations. Should stop
+ // before upon reaching a correct mesh.
nRelaxIter 3;
+
+ // Number of smoothing iterations of surface normals
nSmoothSurfaceNormals 1;
+
+ // Number of smoothing iterations of interior mesh movement direction
nSmoothNormals 3;
+
+ // Smooth layer thickness over surface patches
nSmoothThickness 10;
+
+ // Stop layer growth on highly warped cells
maxFaceThicknessRatio 0.5;
+
+ // Reduce layer growth where ratio thickness to medial
+ // distance is large
maxThicknessToMedialRatio 0.3;
+
+ // Angle used to pick up medial axis points
minMedianAxisAngle 130;
+
+ // Create buffer region for new layer terminations
nBufferCellsNoExtrude 0;
}
+
+// Generic mesh quality settings. At any undoable phase these determine
+// where to undo.
meshQualityControls
{
+ //- Maximum non-orthogonality allowed. Set to 180 to disable.
maxNonOrtho 65;
+
+ //- Max skewness allowed. Set to <0 to disable.
maxBoundarySkewness 20;
maxInternalSkewness 4;
+
+ //- Max concaveness allowed. Is angle (in degrees) below which concavity
+ // is allowed. 0 is straight face, <0 would be convex face.
+ // Set to 180 to disable.
maxConcave 80;
+
+ //- Minimum projected area v.s. actual area. Set to -1 to disable.
minFlatness 0.5;
+
+ //- Minimum pyramid volume. Is absolute volume of cell pyramid.
+ // Set to very negative number (e.g. -1E30) to disable.
minVol 1e-13;
+
+ //- Minimum face area. Set to <0 to disable.
minArea -1;
+
+ //- Minimum face twist. Set to <-1 to disable. dot product of face normal
+ //- and face centre triangles normal
minTwist 0.02;
+
+ //- minimum normalised cell determinant
+ //- 1 = hex, <= 0 = folded or flattened illegal cell
minDeterminant 0.001;
+
+ //- minFaceWeight (0 -> 0.5)
minFaceWeight 0.02;
+
+ //- minVolRatio (0 -> 1)
minVolRatio 0.01;
+
+ //must be >0 for Fluent compatibility
minTriangleTwist -1;
+
+
+ // Advanced
+
+ //- Number of error distribution iterations
nSmoothScale 4;
+ //- amount to scale back displacement at error points
errorReduction 0.75;
}
+
+
+// Advanced
+
+// Flags for optional output
+// 0 : only write final meshes
+// 1 : write intermediate meshes
+// 2 : write volScalarField with cellLevel for postprocessing
+// 4 : write current intersections as .obj files
debug 0;
+
+// Merge tolerance. Is fraction of overall bounding box of initial mesh.
+// Note: the write tolerance needs to be higher than this.
mergeTolerance 1e-06;
diff --git a/wmake/wmakeLnInclude b/wmake/wmakeLnInclude
index bce2cfbb51d3cfcb74c52a54432dd386bbdafc49..1ed1671e7a32e3043027fe2a680eb3e3031048ca 100755
--- a/wmake/wmakeLnInclude
+++ b/wmake/wmakeLnInclude
@@ -145,16 +145,6 @@ find -L . -type l -exec rm \{\} \;
find .. $findOpt \
\( -name lnInclude -o -name Make -o -name config -o -name noLink \) -prune \
-o \( -name '*.[CHh]' -o -name '*.[ch]xx' -o -name '*.[ch]pp' -o -name '*.type' \) \
- -a ! -name ".#*" \
- -print | \
- while read src
- do
- link=$(readlink ${src##*/})
- if [ "$link" != "$src" ]
- then
- rm $link 2>/dev/null
- ln $lnOpt $src .
- fi
- done
+ -exec ln $lnOpt {} . \;
#------------------------------------------------------------------------------