mirror of
https://develop.openfoam.com/Development/openfoam.git
synced 2025-11-28 03:28:01 +00:00
Merge branch 'master' of /home/dm4/OpenFOAM/OpenFOAM-dev
This commit is contained in:
andy
9
applications/solvers/compressible/rhoSimpleFoam/Allwclean
Executable file
9
applications/solvers/compressible/rhoSimpleFoam/Allwclean
Executable file
@ -0,0 +1,9 @@
|
||||
#!/bin/sh
|
||||
cd ${0%/*} || exit 1 # run from this directory
|
||||
set -x
|
||||
|
||||
wclean
|
||||
wclean rhoPorousMRFSimpleFoam
|
||||
wclean rhoSimplecFoam
|
||||
|
||||
# ----------------------------------------------------------------- end-of-file
|
||||
@ -4,5 +4,6 @@ set -x
|
||||
|
||||
wmake
|
||||
wmake rhoPorousMRFSimpleFoam
|
||||
wmake rhoSimplecFoam
|
||||
|
||||
# ----------------------------------------------------------------- end-of-file
|
||||
|
||||
@ -3,6 +3,7 @@
|
||||
tmp<fvVectorMatrix> UEqn
|
||||
(
|
||||
fvm::div(phi, U)
|
||||
- fvm::Sp(fvc::div(phi), U)
|
||||
+ turbulence->divDevRhoReff(U)
|
||||
);
|
||||
|
||||
|
||||
@ -1,11 +1,13 @@
|
||||
{
|
||||
volScalarField K("K", 0.5*magSqr(U));
|
||||
|
||||
fvScalarMatrix hEqn
|
||||
(
|
||||
fvm::div(phi, h)
|
||||
- fvm::Sp(fvc::div(phi), h)
|
||||
- fvm::laplacian(turbulence->alphaEff(), h)
|
||||
==
|
||||
- fvc::div(phi, 0.5*magSqr(U), "div(phi,K)")
|
||||
fvc::div(phi)*K - fvc::div(phi, K)
|
||||
);
|
||||
|
||||
hEqn.relax();
|
||||
|
||||
@ -28,7 +28,7 @@ if (simple.transonic())
|
||||
);
|
||||
|
||||
// Relax the pressure equation to ensure diagonal-dominance
|
||||
pEqn.relax(mesh.equationRelaxationFactor("pEqn"));
|
||||
pEqn.relax();
|
||||
|
||||
pEqn.setReference(pRefCell, pRefValue);
|
||||
|
||||
|
||||
@ -1,5 +1,5 @@
|
||||
EXE_INC = \
|
||||
-I../rhoSimpleFoam \
|
||||
-I.. \
|
||||
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
|
||||
-I$(LIB_SRC)/turbulenceModels \
|
||||
-I$(LIB_SRC)/turbulenceModels/compressible/RAS/RASModel \
|
||||
@ -29,21 +29,17 @@ if (simple.transonic())
|
||||
(
|
||||
"phic",
|
||||
fvc::interpolate(rho/AtU - rho/AU)*fvc::snGrad(p)*mesh.magSf()
|
||||
+ phid*(fvc::interpolate(p) - fvc::interpolate(p, "UD"))
|
||||
);
|
||||
|
||||
//refCast<mixedFvPatchScalarField>(p.boundaryField()[1]).refValue()
|
||||
// = p.boundaryField()[1];
|
||||
|
||||
fvScalarMatrix pEqn
|
||||
(
|
||||
fvm::div(phid, p)
|
||||
+ fvc::div(phic)
|
||||
- fvm::Sp(fvc::div(phid), p)
|
||||
+ fvc::div(phid)*p
|
||||
- fvm::laplacian(rho/AtU, p)
|
||||
);
|
||||
//pEqn.relax();
|
||||
|
||||
// Relax the pressure equation to ensure diagonal-dominance
|
||||
pEqn.relax();
|
||||
|
||||
pEqn.setReference(pRefCell, pRefValue);
|
||||
|
||||
@ -71,7 +67,6 @@ else
|
||||
fvScalarMatrix pEqn
|
||||
(
|
||||
fvc::div(phiHbyA)
|
||||
//- fvm::laplacian(rho/AU, p)
|
||||
- fvm::laplacian(rho/AtU, p)
|
||||
);
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -310,14 +310,21 @@
|
||||
|
||||
Info << "dispersedPhase is " << dispersedPhase << endl;
|
||||
|
||||
scalar minInterfaceAlpha
|
||||
scalar residualPhaseFraction
|
||||
(
|
||||
readScalar
|
||||
(
|
||||
interfacialProperties.lookup("minInterfaceAlpha")
|
||||
interfacialProperties.lookup("residualPhaseFraction")
|
||||
)
|
||||
);
|
||||
|
||||
dimensionedScalar residualSlip
|
||||
(
|
||||
"residualSlip",
|
||||
dimVelocity,
|
||||
interfacialProperties.lookup("residualSlip")
|
||||
);
|
||||
|
||||
kineticTheoryModel kineticTheory
|
||||
(
|
||||
phase1,
|
||||
|
||||
@ -36,7 +36,7 @@ volScalarField heatTransferCoeff
|
||||
|
||||
{
|
||||
volVectorField Ur(U1 - U2);
|
||||
volScalarField magUr(mag(Ur));
|
||||
volScalarField magUr(mag(Ur) + residualSlip);
|
||||
|
||||
if (dispersedPhase == "1")
|
||||
{
|
||||
@ -69,12 +69,9 @@ volScalarField heatTransferCoeff
|
||||
<< exit(FatalError);
|
||||
}
|
||||
|
||||
volScalarField alpha1Coeff
|
||||
(
|
||||
(alpha1 + minInterfaceAlpha)*(alpha2 + minInterfaceAlpha)
|
||||
);
|
||||
dragCoeff *= alpha1Coeff;
|
||||
heatTransferCoeff *= alpha1Coeff;
|
||||
volScalarField alphaCoeff(max(alpha1*alpha2, residualPhaseFraction));
|
||||
dragCoeff *= alphaCoeff;
|
||||
heatTransferCoeff *= alphaCoeff;
|
||||
|
||||
liftForce = Cl*(alpha1*rho1 + alpha2*rho2)*(Ur ^ fvc::curl(U));
|
||||
|
||||
|
||||
@ -32,7 +32,7 @@
|
||||
|
||||
dimensionedScalar rAUf("(1|A(U))", dimTime/rho.dimensions(), 1.0);
|
||||
|
||||
adjustPhi(phi, U, pcorr);
|
||||
//adjustPhi(phi, U, pcorr);
|
||||
|
||||
while (pimple.correctNonOrthogonal())
|
||||
{
|
||||
|
||||
3
applications/test/tetTetOverlap/Make/files
Normal file
3
applications/test/tetTetOverlap/Make/files
Normal file
@ -0,0 +1,3 @@
|
||||
Test-tetTetOverlap.C
|
||||
|
||||
EXE = $(FOAM_USER_APPBIN)/Test-tetTetOverlap
|
||||
5
applications/test/tetTetOverlap/Make/options
Normal file
5
applications/test/tetTetOverlap/Make/options
Normal file
@ -0,0 +1,5 @@
|
||||
EXE_INC = \
|
||||
-I$(LIB_SRC)/meshTools/lnInclude
|
||||
|
||||
EXE_LIBS = \
|
||||
-lmeshTools
|
||||
164
applications/test/tetTetOverlap/Test-tetTetOverlap.C
Normal file
164
applications/test/tetTetOverlap/Test-tetTetOverlap.C
Normal file
@ -0,0 +1,164 @@
|
||||
/*---------------------------------------------------------------------------*\
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
This file is part of OpenFOAM.
|
||||
|
||||
OpenFOAM is free software: you can redistribute it and/or modify it
|
||||
under the terms of the GNU General Public License as published by
|
||||
the Free Software Foundation, either version 3 of the License, or
|
||||
(at your option) any later version.
|
||||
|
||||
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
|
||||
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
|
||||
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
|
||||
for more details.
|
||||
|
||||
You should have received a copy of the GNU General Public License
|
||||
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
Application
|
||||
Test-tetTetOverlap
|
||||
|
||||
Description
|
||||
Overlap volume of two tets
|
||||
|
||||
\*---------------------------------------------------------------------------*/
|
||||
|
||||
#include "tetrahedron.H"
|
||||
#include "OFstream.H"
|
||||
#include "meshTools.H"
|
||||
|
||||
using namespace Foam;
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
void writeOBJ
|
||||
(
|
||||
Ostream& os,
|
||||
label& vertI,
|
||||
const tetPoints& tet
|
||||
)
|
||||
{
|
||||
forAll(tet, fp)
|
||||
{
|
||||
meshTools::writeOBJ(os, tet[fp]);
|
||||
}
|
||||
os << "l " << vertI+1 << ' ' << vertI+2 << nl
|
||||
<< "l " << vertI+1 << ' ' << vertI+3 << nl
|
||||
<< "l " << vertI+1 << ' ' << vertI+4 << nl
|
||||
<< "l " << vertI+2 << ' ' << vertI+3 << nl
|
||||
<< "l " << vertI+2 << ' ' << vertI+4 << nl
|
||||
<< "l " << vertI+3 << ' ' << vertI+4 << nl;
|
||||
vertI += 4;
|
||||
}
|
||||
|
||||
|
||||
int main(int argc, char *argv[])
|
||||
{
|
||||
tetPoints A
|
||||
(
|
||||
point(0, 0, 0),
|
||||
point(1, 0, 0),
|
||||
point(1, 1, 0),
|
||||
point(1, 1, 1)
|
||||
);
|
||||
const tetPointRef tetA = A.tet();
|
||||
|
||||
tetPoints B
|
||||
(
|
||||
point(0.1, 0.1, 0.1),
|
||||
point(1.1, 0.1, 0.1),
|
||||
point(1.1, 1.1, 0.1),
|
||||
point(1.1, 1.1, 1.1)
|
||||
);
|
||||
const tetPointRef tetB = B.tet();
|
||||
|
||||
|
||||
tetPointRef::tetIntersectionList insideTets;
|
||||
label nInside = 0;
|
||||
tetPointRef::tetIntersectionList outsideTets;
|
||||
label nOutside = 0;
|
||||
|
||||
tetA.tetOverlap
|
||||
(
|
||||
tetB,
|
||||
insideTets,
|
||||
nInside,
|
||||
outsideTets,
|
||||
nOutside
|
||||
);
|
||||
|
||||
|
||||
// Dump to file
|
||||
// ~~~~~~~~~~~~
|
||||
|
||||
{
|
||||
OFstream str("tetA.obj");
|
||||
Info<< "Writing A to " << str.name() << endl;
|
||||
label vertI = 0;
|
||||
writeOBJ(str, vertI, A);
|
||||
}
|
||||
{
|
||||
OFstream str("tetB.obj");
|
||||
Info<< "Writing B to " << str.name() << endl;
|
||||
label vertI = 0;
|
||||
writeOBJ(str, vertI, B);
|
||||
}
|
||||
{
|
||||
OFstream str("inside.obj");
|
||||
Info<< "Writing parts of A inside B to " << str.name() << endl;
|
||||
label vertI = 0;
|
||||
for (label i = 0; i < nInside; ++i)
|
||||
{
|
||||
writeOBJ(str, vertI, insideTets[i]);
|
||||
}
|
||||
}
|
||||
{
|
||||
OFstream str("outside.obj");
|
||||
Info<< "Writing parts of A outside B to " << str.name() << endl;
|
||||
label vertI = 0;
|
||||
for (label i = 0; i < nOutside; ++i)
|
||||
{
|
||||
writeOBJ(str, vertI, outsideTets[i]);
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
// Check
|
||||
// ~~~~~
|
||||
|
||||
Info<< "Vol A:" << tetA.mag() << endl;
|
||||
|
||||
scalar volInside = 0;
|
||||
for (label i = 0; i < nInside; ++i)
|
||||
{
|
||||
volInside += insideTets[i].tet().mag();
|
||||
}
|
||||
Info<< "Vol A inside B:" << volInside << endl;
|
||||
|
||||
scalar volOutside = 0;
|
||||
for (label i = 0; i < nOutside; ++i)
|
||||
{
|
||||
volOutside += outsideTets[i].tet().mag();
|
||||
}
|
||||
Info<< "Vol A outside B:" << volOutside << endl;
|
||||
|
||||
Info<< "Sum inside and outside:" << volInside+volOutside << endl;
|
||||
|
||||
if (mag(volInside+volOutside-tetA.mag()) > SMALL)
|
||||
{
|
||||
FatalErrorIn("Test-tetetOverlap")
|
||||
<< "Tet volumes do not sum up to input tet."
|
||||
<< exit(FatalError);
|
||||
}
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
||||
|
||||
// ************************************************************************* //
|
||||
@ -4,5 +4,7 @@ set -x
|
||||
|
||||
wclean libso conformalVoronoiMesh
|
||||
wclean
|
||||
wclean cvMeshSurfaceSimplify
|
||||
wclean cvMeshBackgroundMesh
|
||||
|
||||
# ----------------------------------------------------------------- end-of-file
|
||||
|
||||
@ -4,6 +4,7 @@ set -x
|
||||
|
||||
wmake libso conformalVoronoiMesh
|
||||
wmake
|
||||
wmake cvMeshSurfaceSimplify
|
||||
wmake cvMeshBackgroundMesh
|
||||
(cd cvMeshSurfaceSimplify && ./Allwmake)
|
||||
|
||||
# ----------------------------------------------------------------- end-of-file
|
||||
|
||||
@ -2,11 +2,15 @@ EXE_DEBUG = -DFULLDEBUG -g -O0
|
||||
EXE_FROUNDING_MATH = -frounding-math
|
||||
EXE_NDEBUG = -DNDEBUG
|
||||
|
||||
CGAL_EXACT =
|
||||
CGAL_INEXACT = -DCGAL_INEXACT
|
||||
|
||||
include $(GENERAL_RULES)/CGAL
|
||||
|
||||
EXE_INC = \
|
||||
${EXE_FROUNDING_MATH} \
|
||||
${EXE_NDEBUG} \
|
||||
${CGAL_INEXACT} \
|
||||
${CGAL_INC} \
|
||||
-IconformalVoronoiMesh/lnInclude \
|
||||
-I$(LIB_SRC)/finiteVolume/lnInclude \
|
||||
|
||||
@ -2,12 +2,16 @@ EXE_DEBUG = -DFULLDEBUG -g -O0
|
||||
EXE_FROUNDING_MATH = -frounding-math
|
||||
EXE_NDEBUG = -DNDEBUG
|
||||
|
||||
CGAL_EXACT =
|
||||
CGAL_INEXACT = -DCGAL_INEXACT
|
||||
|
||||
include $(GENERAL_RULES)/CGAL
|
||||
FFLAGS = -DCGAL_FILES='"${CGAL_ARCH_PATH}/share/files"'
|
||||
|
||||
EXE_INC = \
|
||||
${EXE_FROUNDING_MATH} \
|
||||
${EXE_NDEBUG} \
|
||||
${CGAL_INEXACT} \
|
||||
${CGAL_INC} \
|
||||
-I$(LIB_SRC)/finiteVolume/lnInclude \
|
||||
-I$(LIB_SRC)/meshTools/lnInclude \
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -142,7 +142,7 @@ void Foam::backgroundMeshDecomposition::initialRefinement()
|
||||
zeroGradientFvPatchScalarField::typeName
|
||||
);
|
||||
|
||||
const conformationSurfaces& geometry = cvMesh_.geometryToConformTo();
|
||||
const conformationSurfaces& geometry = geometryToConformTo_;
|
||||
|
||||
decompositionMethod& decomposer = decomposerPtr_();
|
||||
|
||||
@ -512,7 +512,7 @@ bool Foam::backgroundMeshDecomposition::refineCell
|
||||
// Sample the box to find an estimate of the min size, and a volume
|
||||
// estimate when overlapping == true.
|
||||
|
||||
const conformationSurfaces& geometry = cvMesh_.geometryToConformTo();
|
||||
const conformationSurfaces& geometry = geometryToConformTo_;
|
||||
|
||||
treeBoundBox cellBb
|
||||
(
|
||||
@ -578,7 +578,7 @@ bool Foam::backgroundMeshDecomposition::refineCell
|
||||
|
||||
forAll(samplePoints, i)
|
||||
{
|
||||
scalar s = cvMesh_.cellSizeControl().cellSize
|
||||
scalar s = cellSizeControl_.cellSize
|
||||
(
|
||||
hitInfo[i].hitPoint()
|
||||
);
|
||||
@ -693,7 +693,7 @@ void Foam::backgroundMeshDecomposition::buildPatchAndTree()
|
||||
// Overall bb
|
||||
treeBoundBox overallBb(boundaryFacesPtr_().localPoints());
|
||||
|
||||
Random& rnd = cvMesh_.rndGen();
|
||||
Random& rnd = rndGen_;
|
||||
|
||||
bFTreePtr_.reset
|
||||
(
|
||||
@ -726,11 +726,11 @@ void Foam::backgroundMeshDecomposition::buildPatchAndTree()
|
||||
|
||||
octreeNearestDistances_ = bFTreePtr_().calcNearestDistance();
|
||||
|
||||
if (cvMesh_.cvMeshControls().objOutput())
|
||||
if (cvMeshControls_.objOutput())
|
||||
{
|
||||
OFstream fStr
|
||||
(
|
||||
cvMesh_.time().path()
|
||||
mesh_.time().path()
|
||||
/"backgroundMeshDecomposition_proc_"
|
||||
+ name(Pstream::myProcNo())
|
||||
+ "_boundaryFaces.obj"
|
||||
@ -777,15 +777,18 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
|
||||
const conformalVoronoiMesh& cvMesh
|
||||
)
|
||||
:
|
||||
coeffsDict_(coeffsDict),
|
||||
cvMesh_(cvMesh),
|
||||
runTime_(cvMesh.time()),
|
||||
geometryToConformTo_(cvMesh.geometryToConformTo()),
|
||||
cellSizeControl_(cvMesh.cellSizeControl()),
|
||||
rndGen_(cvMesh.rndGen()),
|
||||
cvMeshControls_(cvMesh.cvMeshControls()),
|
||||
mesh_
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
fvMesh::defaultRegion,
|
||||
cvMesh_.time().timeName(),
|
||||
cvMesh_.time(),
|
||||
runTime_.timeName(),
|
||||
runTime_,
|
||||
IOobject::MUST_READ
|
||||
)
|
||||
),
|
||||
@ -805,22 +808,22 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
|
||||
IOobject
|
||||
(
|
||||
"decomposeParDict",
|
||||
cvMesh_.time().system(),
|
||||
cvMesh_.time(),
|
||||
runTime_.system(),
|
||||
runTime_,
|
||||
IOobject::MUST_READ_IF_MODIFIED,
|
||||
IOobject::NO_WRITE
|
||||
)
|
||||
),
|
||||
decomposerPtr_(decompositionMethod::New(decomposeDict_)),
|
||||
mergeDist_(1e-6*mesh_.bounds().mag()),
|
||||
spanScale_(readScalar(coeffsDict_.lookup("spanScale"))),
|
||||
spanScale_(readScalar(coeffsDict.lookup("spanScale"))),
|
||||
minCellSizeLimit_
|
||||
(
|
||||
coeffsDict_.lookupOrDefault<scalar>("minCellSizeLimit", 0.0)
|
||||
coeffsDict.lookupOrDefault<scalar>("minCellSizeLimit", 0.0)
|
||||
),
|
||||
minLevels_(readLabel(coeffsDict_.lookup("minLevels"))),
|
||||
volRes_(readLabel(coeffsDict_.lookup("sampleResolution"))),
|
||||
maxCellWeightCoeff_(readScalar(coeffsDict_.lookup("maxCellWeightCoeff")))
|
||||
minLevels_(readLabel(coeffsDict.lookup("minLevels"))),
|
||||
volRes_(readLabel(coeffsDict.lookup("sampleResolution"))),
|
||||
maxCellWeightCoeff_(readScalar(coeffsDict.lookup("maxCellWeightCoeff")))
|
||||
{
|
||||
if (!Pstream::parRun())
|
||||
{
|
||||
@ -854,6 +857,74 @@ Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
|
||||
}
|
||||
|
||||
|
||||
Foam::backgroundMeshDecomposition::backgroundMeshDecomposition
|
||||
(
|
||||
const scalar spanScale,
|
||||
const scalar minCellSizeLimit,
|
||||
const label minLevels,
|
||||
const label volRes,
|
||||
const scalar maxCellWeightCoeff,
|
||||
|
||||
const Time& runTime,
|
||||
const conformationSurfaces& geometryToConformTo,
|
||||
const cellSizeControlSurfaces& cellSizeControl,
|
||||
Random& rndGen,
|
||||
const cvControls& cvMeshControls
|
||||
)
|
||||
:
|
||||
runTime_(runTime),
|
||||
geometryToConformTo_(geometryToConformTo),
|
||||
cellSizeControl_(cellSizeControl),
|
||||
rndGen_(rndGen),
|
||||
cvMeshControls_(cvMeshControls),
|
||||
mesh_
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
fvMesh::defaultRegion,
|
||||
runTime_.timeName(),
|
||||
runTime_,
|
||||
IOobject::MUST_READ
|
||||
)
|
||||
),
|
||||
meshCutter_
|
||||
(
|
||||
mesh_,
|
||||
labelList(mesh_.nCells(), 0),
|
||||
labelList(mesh_.nPoints(), 0)
|
||||
),
|
||||
boundaryFacesPtr_(),
|
||||
bFTreePtr_(),
|
||||
octreeNearestDistances_(),
|
||||
allBackgroundMeshBounds_(Pstream::nProcs()),
|
||||
globalBackgroundBounds_(),
|
||||
decomposeDict_
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"decomposeParDict",
|
||||
runTime_.system(),
|
||||
runTime_,
|
||||
IOobject::MUST_READ_IF_MODIFIED,
|
||||
IOobject::NO_WRITE
|
||||
)
|
||||
),
|
||||
decomposerPtr_(decompositionMethod::New(decomposeDict_)),
|
||||
mergeDist_(1e-6*mesh_.bounds().mag()),
|
||||
spanScale_(spanScale),
|
||||
minCellSizeLimit_(minCellSizeLimit),
|
||||
minLevels_(minLevels),
|
||||
volRes_(volRes),
|
||||
maxCellWeightCoeff_(maxCellWeightCoeff)
|
||||
{
|
||||
// Stand-alone operation
|
||||
|
||||
Info<< nl << "Building initial background mesh decomposition" << endl;
|
||||
|
||||
initialRefinement();
|
||||
}
|
||||
|
||||
|
||||
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
|
||||
|
||||
Foam::backgroundMeshDecomposition::~backgroundMeshDecomposition()
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -92,10 +92,22 @@ class backgroundMeshDecomposition
|
||||
// Private data
|
||||
|
||||
//- Method details dictionary
|
||||
dictionary coeffsDict_;
|
||||
//dictionary coeffsDict_;
|
||||
|
||||
//- Reference to the conformalVoronoiMesh holding this object
|
||||
const conformalVoronoiMesh& cvMesh_;
|
||||
//- Reference to runtime
|
||||
const Time& runTime_;
|
||||
|
||||
//- Reference to surface
|
||||
const conformationSurfaces& geometryToConformTo_;
|
||||
|
||||
//- The cell size control object
|
||||
const cellSizeControlSurfaces& cellSizeControl_;
|
||||
|
||||
//- Random number generator
|
||||
Random& rndGen_;
|
||||
|
||||
//- Controls
|
||||
const cvControls& cvMeshControls_;
|
||||
|
||||
//- Mesh stored on for this processor, specifiying the domain that it
|
||||
// is responsible for.
|
||||
@ -191,13 +203,29 @@ public:
|
||||
|
||||
// Constructors
|
||||
|
||||
//- Construct from components
|
||||
//- Construct from components in cvMesh operation
|
||||
backgroundMeshDecomposition
|
||||
(
|
||||
const dictionary& coeffsDict,
|
||||
const conformalVoronoiMesh& cvMesh
|
||||
);
|
||||
|
||||
//- Construct from components for standalone operation
|
||||
backgroundMeshDecomposition
|
||||
(
|
||||
const scalar spanScale,
|
||||
const scalar minCellSizeLimit,
|
||||
const label minLevels,
|
||||
const label volRes,
|
||||
const scalar maxCellWeightCoeff,
|
||||
|
||||
const Time& runTime,
|
||||
const conformationSurfaces& geometryToConformTo,
|
||||
const cellSizeControlSurfaces& cellSizeControl,
|
||||
Random& rndGen,
|
||||
const cvControls& cvMeshControls
|
||||
);
|
||||
|
||||
|
||||
//- Destructor
|
||||
~backgroundMeshDecomposition();
|
||||
@ -299,6 +327,15 @@ public:
|
||||
|
||||
//- Return the boundBox of this processor
|
||||
inline const treeBoundBox& procBounds() const;
|
||||
|
||||
//- Return the cell level of the underlying mesh
|
||||
inline const labelList& cellLevel() const;
|
||||
|
||||
//- Return the point level of the underlying mesh
|
||||
inline const labelList& pointLevel() const;
|
||||
|
||||
//- Return the current decomposition method
|
||||
inline const decompositionMethod& decomposer() const;
|
||||
};
|
||||
|
||||
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -43,9 +43,30 @@ Foam::backgroundMeshDecomposition::octreeNearestDistances() const
|
||||
}
|
||||
|
||||
|
||||
const Foam::treeBoundBox& Foam::backgroundMeshDecomposition::procBounds() const
|
||||
const Foam::treeBoundBox&
|
||||
Foam::backgroundMeshDecomposition::procBounds() const
|
||||
{
|
||||
return allBackgroundMeshBounds_[Pstream::myProcNo()];
|
||||
}
|
||||
|
||||
|
||||
const Foam::labelList& Foam::backgroundMeshDecomposition::cellLevel() const
|
||||
{
|
||||
return meshCutter_.cellLevel();
|
||||
}
|
||||
|
||||
|
||||
const Foam::labelList& Foam::backgroundMeshDecomposition::pointLevel() const
|
||||
{
|
||||
return meshCutter_.pointLevel();
|
||||
}
|
||||
|
||||
|
||||
const Foam::decompositionMethod&
|
||||
Foam::backgroundMeshDecomposition::decomposer() const
|
||||
{
|
||||
return decomposerPtr_();
|
||||
}
|
||||
|
||||
|
||||
// ************************************************************************* //
|
||||
|
||||
@ -249,6 +249,7 @@ Foam::conformationSurfaces::conformationSurfaces
|
||||
|
||||
Info<< endl
|
||||
<< "Testing for locationInMesh " << locationInMesh_ << endl;
|
||||
|
||||
forAll(surfaces_, s)
|
||||
{
|
||||
const searchableSurface& surface(allGeometry_[surfaces_[s]]);
|
||||
@ -267,7 +268,6 @@ Foam::conformationSurfaces::conformationSurfaces
|
||||
|
||||
referenceVolumeTypes_[s] = vTypes[0];
|
||||
|
||||
|
||||
Info<< " is "
|
||||
<< searchableSurface::volumeTypeNames[referenceVolumeTypes_[s]]
|
||||
<< " surface " << surface.name()
|
||||
|
||||
@ -303,6 +303,9 @@ public:
|
||||
//- Find which patch is closest to the point
|
||||
label findPatch(const point& pt) const;
|
||||
|
||||
//- Is the surface a baffle
|
||||
inline bool isBaffle(const label index) const;
|
||||
|
||||
|
||||
// Write
|
||||
|
||||
|
||||
@ -62,4 +62,10 @@ const Foam::treeBoundBox& Foam::conformationSurfaces::globalBounds() const
|
||||
}
|
||||
|
||||
|
||||
bool Foam::conformationSurfaces::isBaffle(const label index) const
|
||||
{
|
||||
return baffleSurfaces_[index];
|
||||
}
|
||||
|
||||
|
||||
// ************************************************************************* //
|
||||
|
||||
@ -61,7 +61,7 @@ Foam::cvControls::cvControls
|
||||
|
||||
specialiseFeaturePoints_ = Switch
|
||||
(
|
||||
surfDict.lookupOrDefault<Switch>("specialiseFeaturePoints", false)
|
||||
surfDict.lookup("specialiseFeaturePoints")
|
||||
);
|
||||
|
||||
surfaceSearchDistanceCoeff_ = readScalar
|
||||
|
||||
@ -0,0 +1,3 @@
|
||||
cvMeshBackgroundMesh.C
|
||||
|
||||
EXE = $(FOAM_APPBIN)/cvMeshBackgroundMesh
|
||||
@ -0,0 +1,29 @@
|
||||
EXE_DEBUG = -DFULLDEBUG -g -O0
|
||||
EXE_FROUNDING_MATH = -frounding-math
|
||||
EXE_NDEBUG = -DNDEBUG
|
||||
|
||||
include $(GENERAL_RULES)/CGAL
|
||||
|
||||
EXE_INC = \
|
||||
${EXE_FROUNDING_MATH} \
|
||||
${EXE_NDEBUG} \
|
||||
${CGAL_INC} \
|
||||
-I../conformalVoronoiMesh/lnInclude \
|
||||
-I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
|
||||
-I$(LIB_SRC)/edgeMesh/lnInclude \
|
||||
-I$(LIB_SRC)/sampling/lnInclude \
|
||||
-I$(LIB_SRC)/triSurface/lnInclude \
|
||||
-I$(LIB_SRC)/finiteVolume/lnInclude \
|
||||
-I$(LIB_SRC)/meshTools/lnInclude \
|
||||
-I$(LIB_SRC)/dynamicMesh/lnInclude
|
||||
|
||||
EXE_LIBS = \
|
||||
$(CGAL_LIBS) \
|
||||
-lconformalVoronoiMesh \
|
||||
-ldecompositionMethods /* -L$(FOAM_LIBBIN)/dummy -lscotchDecomp */ \
|
||||
-ledgeMesh \
|
||||
-lsampling \
|
||||
-ltriSurface \
|
||||
-lmeshTools \
|
||||
-ldynamicMesh \
|
||||
-lfiniteVolume
|
||||
@ -0,0 +1,770 @@
|
||||
/*---------------------------------------------------------------------------*\
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
This file is part of OpenFOAM.
|
||||
|
||||
OpenFOAM is free software: you can redistribute it and/or modify it
|
||||
under the terms of the GNU General Public License as published by
|
||||
the Free Software Foundation, either version 3 of the License, or
|
||||
(at your option) any later version.
|
||||
|
||||
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
|
||||
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
|
||||
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
|
||||
for more details.
|
||||
|
||||
You should have received a copy of the GNU General Public License
|
||||
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
Application
|
||||
cvMeshBackGroundMesh
|
||||
|
||||
Description
|
||||
Writes out background mesh as constructed by cvMesh and constructs
|
||||
distanceSurface.
|
||||
|
||||
\*---------------------------------------------------------------------------*/
|
||||
|
||||
#include "PatchTools.H"
|
||||
#include "argList.H"
|
||||
#include "Time.H"
|
||||
#include "triSurface.H"
|
||||
#include "searchableSurfaces.H"
|
||||
#include "conformationSurfaces.H"
|
||||
#include "cellSizeControlSurfaces.H"
|
||||
#include "backgroundMeshDecomposition.H"
|
||||
#include "cellShape.H"
|
||||
#include "cellModeller.H"
|
||||
#include "DynamicField.H"
|
||||
#include "isoSurfaceCell.H"
|
||||
#include "vtkSurfaceWriter.H"
|
||||
#include "syncTools.H"
|
||||
|
||||
using namespace Foam;
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
// Tolerance (as fraction of the bounding box). Needs to be fairly lax since
|
||||
// usually meshes get written with limited precision (6 digits)
|
||||
static const scalar defaultMergeTol = 1E-6;
|
||||
|
||||
// Get merging distance when matching face centres
|
||||
scalar getMergeDistance
|
||||
(
|
||||
const argList& args,
|
||||
const Time& runTime,
|
||||
const boundBox& bb
|
||||
)
|
||||
{
|
||||
scalar mergeTol = defaultMergeTol;
|
||||
args.optionReadIfPresent("mergeTol", mergeTol);
|
||||
|
||||
scalar writeTol =
|
||||
Foam::pow(scalar(10.0), -scalar(IOstream::defaultPrecision()));
|
||||
|
||||
Info<< "Merge tolerance : " << mergeTol << nl
|
||||
<< "Write tolerance : " << writeTol << endl;
|
||||
|
||||
if (runTime.writeFormat() == IOstream::ASCII && mergeTol < writeTol)
|
||||
{
|
||||
FatalErrorIn("getMergeDistance")
|
||||
<< "Your current settings specify ASCII writing with "
|
||||
<< IOstream::defaultPrecision() << " digits precision." << endl
|
||||
<< "Your merging tolerance (" << mergeTol << ") is finer than this."
|
||||
<< endl
|
||||
<< "Please change your writeFormat to binary"
|
||||
<< " or increase the writePrecision" << endl
|
||||
<< "or adjust the merge tolerance (-mergeTol)."
|
||||
<< exit(FatalError);
|
||||
}
|
||||
|
||||
scalar mergeDist = mergeTol * bb.mag();
|
||||
|
||||
Info<< "Overall meshes bounding box : " << bb << nl
|
||||
<< "Relative tolerance : " << mergeTol << nl
|
||||
<< "Absolute matching distance : " << mergeDist << nl
|
||||
<< endl;
|
||||
|
||||
return mergeDist;
|
||||
}
|
||||
|
||||
|
||||
void printMeshData(const polyMesh& mesh)
|
||||
{
|
||||
// Collect all data on master
|
||||
|
||||
globalIndex globalCells(mesh.nCells());
|
||||
labelListList patchNeiProcNo(Pstream::nProcs());
|
||||
labelListList patchSize(Pstream::nProcs());
|
||||
const labelList& pPatches = mesh.globalData().processorPatches();
|
||||
patchNeiProcNo[Pstream::myProcNo()].setSize(pPatches.size());
|
||||
patchSize[Pstream::myProcNo()].setSize(pPatches.size());
|
||||
forAll(pPatches, i)
|
||||
{
|
||||
const processorPolyPatch& ppp = refCast<const processorPolyPatch>
|
||||
(
|
||||
mesh.boundaryMesh()[pPatches[i]]
|
||||
);
|
||||
patchNeiProcNo[Pstream::myProcNo()][i] = ppp.neighbProcNo();
|
||||
patchSize[Pstream::myProcNo()][i] = ppp.size();
|
||||
}
|
||||
Pstream::gatherList(patchNeiProcNo);
|
||||
Pstream::gatherList(patchSize);
|
||||
|
||||
|
||||
// Print stats
|
||||
|
||||
globalIndex globalBoundaryFaces(mesh.nFaces()-mesh.nInternalFaces());
|
||||
|
||||
label maxProcCells = 0;
|
||||
label totProcFaces = 0;
|
||||
label maxProcPatches = 0;
|
||||
label totProcPatches = 0;
|
||||
label maxProcFaces = 0;
|
||||
|
||||
for (label procI = 0; procI < Pstream::nProcs(); procI++)
|
||||
{
|
||||
Info<< endl
|
||||
<< "Processor " << procI << nl
|
||||
<< " Number of cells = " << globalCells.localSize(procI)
|
||||
<< endl;
|
||||
|
||||
label nProcFaces = 0;
|
||||
|
||||
const labelList& nei = patchNeiProcNo[procI];
|
||||
|
||||
forAll(patchNeiProcNo[procI], i)
|
||||
{
|
||||
Info<< " Number of faces shared with processor "
|
||||
<< patchNeiProcNo[procI][i] << " = " << patchSize[procI][i]
|
||||
<< endl;
|
||||
|
||||
nProcFaces += patchSize[procI][i];
|
||||
}
|
||||
|
||||
Info<< " Number of processor patches = " << nei.size() << nl
|
||||
<< " Number of processor faces = " << nProcFaces << nl
|
||||
<< " Number of boundary faces = "
|
||||
<< globalBoundaryFaces.localSize(procI) << endl;
|
||||
|
||||
maxProcCells = max(maxProcCells, globalCells.localSize(procI));
|
||||
totProcFaces += nProcFaces;
|
||||
totProcPatches += nei.size();
|
||||
maxProcPatches = max(maxProcPatches, nei.size());
|
||||
maxProcFaces = max(maxProcFaces, nProcFaces);
|
||||
}
|
||||
|
||||
// Stats
|
||||
|
||||
scalar avgProcCells = scalar(globalCells.size())/Pstream::nProcs();
|
||||
scalar avgProcPatches = scalar(totProcPatches)/Pstream::nProcs();
|
||||
scalar avgProcFaces = scalar(totProcFaces)/Pstream::nProcs();
|
||||
|
||||
// In case of all faces on one processor. Just to avoid division by 0.
|
||||
if (totProcPatches == 0)
|
||||
{
|
||||
avgProcPatches = 1;
|
||||
}
|
||||
if (totProcFaces == 0)
|
||||
{
|
||||
avgProcFaces = 1;
|
||||
}
|
||||
|
||||
Info<< nl
|
||||
<< "Number of processor faces = " << totProcFaces/2 << nl
|
||||
<< "Max number of cells = " << maxProcCells
|
||||
<< " (" << 100.0*(maxProcCells-avgProcCells)/avgProcCells
|
||||
<< "% above average " << avgProcCells << ")" << nl
|
||||
<< "Max number of processor patches = " << maxProcPatches
|
||||
<< " (" << 100.0*(maxProcPatches-avgProcPatches)/avgProcPatches
|
||||
<< "% above average " << avgProcPatches << ")" << nl
|
||||
<< "Max number of faces between processors = " << maxProcFaces
|
||||
<< " (" << 100.0*(maxProcFaces-avgProcFaces)/avgProcFaces
|
||||
<< "% above average " << avgProcFaces << ")" << nl
|
||||
<< endl;
|
||||
}
|
||||
|
||||
|
||||
// Return cellID
|
||||
label cellLabel
|
||||
(
|
||||
const Vector<label>& nCells,
|
||||
const label i,
|
||||
const label j,
|
||||
const label k
|
||||
)
|
||||
{
|
||||
return i*nCells[1]*nCells[2]+j*nCells[2]+k;
|
||||
}
|
||||
|
||||
label vtxLabel
|
||||
(
|
||||
const Vector<label>& nCells,
|
||||
const label i,
|
||||
const label j,
|
||||
const label k
|
||||
)
|
||||
{
|
||||
Vector<label> nPoints
|
||||
(
|
||||
nCells[0]+1,
|
||||
nCells[1]+1,
|
||||
nCells[2]+1
|
||||
);
|
||||
return i*nPoints[1]*nPoints[2]+j*nPoints[2]+k;
|
||||
}
|
||||
|
||||
|
||||
autoPtr<polyMesh> generateHexMesh
|
||||
(
|
||||
const IOobject& io,
|
||||
const point& origin,
|
||||
const vector& cellSize,
|
||||
const Vector<label>& nCells
|
||||
)
|
||||
{
|
||||
pointField points;
|
||||
if (nCells[0]+nCells[1]+nCells[2] > 0)
|
||||
{
|
||||
points.setSize((nCells[0]+1)*(nCells[1]+1)*(nCells[2]+1));
|
||||
|
||||
// Generate points
|
||||
for (label i = 0; i <= nCells[0]; i++)
|
||||
{
|
||||
for (label j = 0; j <= nCells[1]; j++)
|
||||
{
|
||||
for (label k = 0; k <= nCells[2]; k++)
|
||||
{
|
||||
point pt = origin;
|
||||
pt.x() += i*cellSize[0];
|
||||
pt.y() += j*cellSize[1];
|
||||
pt.z() += k*cellSize[2];
|
||||
points[vtxLabel(nCells, i, j, k)] = pt;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
const cellModel& hex = *(cellModeller::lookup("hex"));
|
||||
cellShapeList cellShapes(nCells[0]*nCells[1]*nCells[2]);
|
||||
|
||||
labelList hexPoints(8);
|
||||
label cellI = 0;
|
||||
for (label i = 0; i < nCells[0]; i++)
|
||||
{
|
||||
for (label j = 0; j < nCells[1]; j++)
|
||||
{
|
||||
for (label k = 0; k < nCells[2]; k++)
|
||||
{
|
||||
hexPoints[0] = vtxLabel(nCells, i, j, k);
|
||||
hexPoints[1] = vtxLabel(nCells, i+1, j, k);
|
||||
hexPoints[2] = vtxLabel(nCells, i+1, j+1, k);
|
||||
hexPoints[3] = vtxLabel(nCells, i, j+1, k);
|
||||
hexPoints[4] = vtxLabel(nCells, i, j, k+1);
|
||||
hexPoints[5] = vtxLabel(nCells, i+1, j, k+1);
|
||||
hexPoints[6] = vtxLabel(nCells, i+1, j+1, k+1);
|
||||
hexPoints[7] = vtxLabel(nCells, i, j+1, k+1);
|
||||
cellShapes[cellI++] = cellShape(hex, hexPoints);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
faceListList boundary(0);
|
||||
wordList patchNames(0);
|
||||
wordList patchTypes(0);
|
||||
word defaultFacesName = "defaultFaces";
|
||||
word defaultFacesType = polyPatch::typeName;
|
||||
wordList patchPhysicalTypes(0);
|
||||
|
||||
return autoPtr<polyMesh>
|
||||
(
|
||||
new polyMesh
|
||||
(
|
||||
io,
|
||||
xferMoveTo<pointField>(points),
|
||||
cellShapes,
|
||||
boundary,
|
||||
patchNames,
|
||||
patchTypes,
|
||||
defaultFacesName,
|
||||
defaultFacesType,
|
||||
patchPhysicalTypes
|
||||
)
|
||||
);
|
||||
}
|
||||
|
||||
|
||||
// Determine for every point a signed distance to the nearest surface
|
||||
// (outside is positive)
|
||||
tmp<scalarField> signedDistance
|
||||
(
|
||||
const scalarField& distSqr,
|
||||
const pointField& points,
|
||||
const searchableSurfaces& geometry,
|
||||
const labelList& surfaces
|
||||
)
|
||||
{
|
||||
tmp<scalarField> tfld(new scalarField(points.size(), Foam::sqr(GREAT)));
|
||||
scalarField& fld = tfld();
|
||||
|
||||
// Find nearest
|
||||
List<pointIndexHit> nearest;
|
||||
labelList nearestSurfaces;
|
||||
searchableSurfacesQueries::findNearest
|
||||
(
|
||||
geometry,
|
||||
surfaces,
|
||||
points,
|
||||
scalarField(points.size(), Foam::sqr(GREAT)),//distSqr
|
||||
nearestSurfaces,
|
||||
nearest
|
||||
);
|
||||
|
||||
// Determine sign of nearest. Sort by surface to do this.
|
||||
DynamicField<point> surfPoints(points.size());
|
||||
DynamicList<label> surfIndices(points.size());
|
||||
|
||||
forAll(surfaces, surfI)
|
||||
{
|
||||
// Extract points on this surface
|
||||
surfPoints.clear();
|
||||
surfIndices.clear();
|
||||
forAll(nearestSurfaces, i)
|
||||
{
|
||||
if (nearestSurfaces[i] == surfI)
|
||||
{
|
||||
surfPoints.append(points[i]);
|
||||
surfIndices.append(i);
|
||||
}
|
||||
}
|
||||
|
||||
// Calculate sideness of these surface points
|
||||
label geomI = surfaces[surfI];
|
||||
List<searchableSurface::volumeType> volType;
|
||||
geometry[geomI].getVolumeType(surfPoints, volType);
|
||||
|
||||
// Push back to original
|
||||
forAll(volType, i)
|
||||
{
|
||||
label pointI = surfIndices[i];
|
||||
scalar dist = mag(points[pointI] - nearest[pointI].hitPoint());
|
||||
|
||||
searchableSurface::volumeType vT = volType[i];
|
||||
|
||||
if (vT == searchableSurface::OUTSIDE)
|
||||
{
|
||||
fld[pointI] = dist;
|
||||
}
|
||||
else if (vT == searchableSurface::INSIDE)
|
||||
{
|
||||
fld[i] = -dist;
|
||||
}
|
||||
else
|
||||
{
|
||||
FatalErrorIn("signedDistance()")
|
||||
<< "getVolumeType failure, neither INSIDE or OUTSIDE"
|
||||
<< exit(FatalError);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
return tfld;
|
||||
}
|
||||
|
||||
|
||||
|
||||
// Main program:
|
||||
|
||||
int main(int argc, char *argv[])
|
||||
{
|
||||
argList::addNote
|
||||
(
|
||||
"Generate cvMesh-consistent representation of surfaces"
|
||||
);
|
||||
argList::addBoolOption
|
||||
(
|
||||
"writeMesh",
|
||||
"write the resulting mesh and distance fields"
|
||||
);
|
||||
argList::addOption
|
||||
(
|
||||
"mergeTol",
|
||||
"scalar",
|
||||
"specify the merge distance relative to the bounding box size "
|
||||
"(default 1E-6)"
|
||||
);
|
||||
|
||||
#include "setRootCase.H"
|
||||
#include "createTime.H"
|
||||
runTime.functionObjects().off();
|
||||
|
||||
const bool writeMesh = args.optionFound("writeMesh");
|
||||
|
||||
if (writeMesh)
|
||||
{
|
||||
Info<< "Writing resulting mesh and cellDistance, pointDistance fields."
|
||||
<< nl << endl;
|
||||
}
|
||||
|
||||
|
||||
IOdictionary cvMeshDict
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"cvMeshDict",
|
||||
runTime.system(),
|
||||
runTime,
|
||||
IOobject::MUST_READ_IF_MODIFIED,
|
||||
IOobject::NO_WRITE
|
||||
)
|
||||
);
|
||||
|
||||
// Define/load all geometry
|
||||
searchableSurfaces allGeometry
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"cvSearchableSurfaces",
|
||||
runTime.constant(),
|
||||
"triSurface",
|
||||
runTime,
|
||||
IOobject::MUST_READ,
|
||||
IOobject::NO_WRITE
|
||||
),
|
||||
cvMeshDict.subDict("geometry")
|
||||
);
|
||||
|
||||
Random rndGen(64293*Pstream::myProcNo());
|
||||
|
||||
conformationSurfaces geometryToConformTo
|
||||
(
|
||||
runTime,
|
||||
rndGen,
|
||||
allGeometry,
|
||||
cvMeshDict.subDict("surfaceConformation")
|
||||
);
|
||||
|
||||
cellSizeControlSurfaces cellSizeControl
|
||||
(
|
||||
allGeometry,
|
||||
cvMeshDict.subDict("motionControl")
|
||||
);
|
||||
|
||||
|
||||
// Generate starting block mesh
|
||||
vector cellSize;
|
||||
{
|
||||
const treeBoundBox& bb = geometryToConformTo.globalBounds();
|
||||
|
||||
// Determine the number of cells in each direction.
|
||||
const vector span = bb.span();
|
||||
vector nScalarCells = span/cellSizeControl.defaultCellSize();
|
||||
|
||||
// Calculate initial cell size to be a little bit smaller than the
|
||||
// defaultCellSize to avoid initial refinement triggering.
|
||||
Vector<label> nCells = Vector<label>
|
||||
(
|
||||
label(nScalarCells.x())+2,
|
||||
label(nScalarCells.y())+2,
|
||||
label(nScalarCells.z())+2
|
||||
);
|
||||
cellSize = vector
|
||||
(
|
||||
span[0]/nCells[0],
|
||||
span[1]/nCells[1],
|
||||
span[2]/nCells[2]
|
||||
);
|
||||
|
||||
Info<< "Generating initial hex mesh with" << nl
|
||||
<< " bounding box : " << bb << nl
|
||||
<< " nCells : " << nCells << nl
|
||||
<< " cellSize : " << cellSize << nl
|
||||
<< endl;
|
||||
|
||||
autoPtr<polyMesh> meshPtr
|
||||
(
|
||||
generateHexMesh
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
polyMesh::defaultRegion,
|
||||
runTime.constant(),
|
||||
runTime
|
||||
),
|
||||
bb.min(),
|
||||
cellSize,
|
||||
(
|
||||
Pstream::master()
|
||||
? nCells
|
||||
: Vector<label>(0, 0, 0)
|
||||
)
|
||||
)
|
||||
);
|
||||
Info<< "Writing initial hex mesh to " << meshPtr().instance() << nl
|
||||
<< endl;
|
||||
meshPtr().write();
|
||||
}
|
||||
|
||||
// Distribute the initial mesh
|
||||
if (Pstream::parRun())
|
||||
{
|
||||
# include "createMesh.H"
|
||||
Info<< "Loaded mesh:" << endl;
|
||||
printMeshData(mesh);
|
||||
|
||||
// Allocate a decomposer
|
||||
IOdictionary decompositionDict
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"decomposeParDict",
|
||||
runTime.system(),
|
||||
mesh,
|
||||
IOobject::MUST_READ_IF_MODIFIED,
|
||||
IOobject::NO_WRITE
|
||||
)
|
||||
);
|
||||
|
||||
autoPtr<decompositionMethod> decomposer
|
||||
(
|
||||
decompositionMethod::New
|
||||
(
|
||||
decompositionDict
|
||||
)
|
||||
);
|
||||
|
||||
labelList decomp = decomposer().decompose(mesh, mesh.cellCentres());
|
||||
|
||||
// Global matching tolerance
|
||||
const scalar tolDim = getMergeDistance
|
||||
(
|
||||
args,
|
||||
runTime,
|
||||
mesh.bounds()
|
||||
);
|
||||
|
||||
// Mesh distribution engine
|
||||
fvMeshDistribute distributor(mesh, tolDim);
|
||||
|
||||
Info<< "Wanted distribution:"
|
||||
<< distributor.countCells(decomp) << nl << endl;
|
||||
|
||||
// Do actual sending/receiving of mesh
|
||||
autoPtr<mapDistributePolyMesh> map = distributor.distribute(decomp);
|
||||
|
||||
// Print some statistics
|
||||
//Info<< "After distribution:" << endl;
|
||||
//printMeshData(mesh);
|
||||
|
||||
mesh.setInstance(runTime.constant());
|
||||
Info<< "Writing redistributed mesh" << nl << endl;
|
||||
mesh.write();
|
||||
}
|
||||
|
||||
|
||||
Info<< "Refining backgroud mesh according to cell size specification" << nl
|
||||
<< endl;
|
||||
|
||||
backgroundMeshDecomposition backgroundMesh
|
||||
(
|
||||
1.0, //spanScale,ratio of poly cell size v.s. hex cell size
|
||||
0.0, //minCellSizeLimit
|
||||
0, //minLevels
|
||||
4, //volRes, check multiple points per cell
|
||||
20.0, //maxCellWeightCoeff
|
||||
runTime,
|
||||
geometryToConformTo,
|
||||
cellSizeControl,
|
||||
rndGen,
|
||||
cvMeshDict
|
||||
);
|
||||
|
||||
if (writeMesh)
|
||||
{
|
||||
runTime++;
|
||||
Info<< "Writing mesh to " << runTime.timeName() << endl;
|
||||
backgroundMesh.mesh().write();
|
||||
}
|
||||
|
||||
const scalar tolDim = getMergeDistance
|
||||
(
|
||||
args,
|
||||
runTime,
|
||||
backgroundMesh.mesh().bounds()
|
||||
);
|
||||
|
||||
|
||||
faceList isoFaces;
|
||||
pointField isoPoints;
|
||||
|
||||
{
|
||||
// Apply a distanceSurface to it.
|
||||
const fvMesh& fvm = backgroundMesh.mesh();
|
||||
|
||||
volScalarField cellDistance
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"cellDistance",
|
||||
fvm.time().timeName(),
|
||||
fvm.time(),
|
||||
IOobject::NO_READ,
|
||||
IOobject::NO_WRITE,
|
||||
false
|
||||
),
|
||||
fvm,
|
||||
dimensionedScalar("zero", dimLength, 0)
|
||||
);
|
||||
|
||||
const searchableSurfaces& geometry = geometryToConformTo.geometry();
|
||||
const labelList& surfaces = geometryToConformTo.surfaces();
|
||||
|
||||
|
||||
// Get maximum search size per cell
|
||||
scalarField distSqr(cellDistance.size());
|
||||
|
||||
const labelList& cellLevel = backgroundMesh.cellLevel();
|
||||
forAll(cellLevel, cellI)
|
||||
{
|
||||
// The largest edge of the cell will always be less than the
|
||||
// span of the bounding box of the cell.
|
||||
distSqr[cellI] = magSqr(cellSize)/pow(2, cellLevel[cellI]);
|
||||
}
|
||||
|
||||
{
|
||||
// Internal field
|
||||
cellDistance.internalField() = signedDistance
|
||||
(
|
||||
distSqr,
|
||||
fvm.C(),
|
||||
geometry,
|
||||
surfaces
|
||||
);
|
||||
// Patch fields
|
||||
forAll(fvm.C().boundaryField(), patchI)
|
||||
{
|
||||
const pointField& cc = fvm.C().boundaryField()[patchI];
|
||||
fvPatchScalarField& fld = cellDistance.boundaryField()[patchI];
|
||||
scalarField patchDistSqr
|
||||
(
|
||||
fld.patch().patchInternalField(distSqr)
|
||||
);
|
||||
fld = signedDistance(patchDistSqr, cc, geometry, surfaces);
|
||||
}
|
||||
|
||||
// On processor patches the fvm.C() will already be the cell centre
|
||||
// on the opposite side so no need to swap cellDistance.
|
||||
|
||||
if (writeMesh)
|
||||
{
|
||||
cellDistance.write();
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
// Distance to points
|
||||
pointScalarField pointDistance
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"pointDistance",
|
||||
fvm.time().timeName(),
|
||||
fvm.time(),
|
||||
IOobject::NO_READ,
|
||||
IOobject::NO_WRITE,
|
||||
false
|
||||
),
|
||||
pointMesh::New(fvm),
|
||||
dimensionedScalar("zero", dimLength, 0)
|
||||
);
|
||||
{
|
||||
scalarField pointDistSqr(fvm.nPoints(), -sqr(GREAT));
|
||||
for (label faceI = 0; faceI < fvm.nInternalFaces(); faceI++)
|
||||
{
|
||||
label own = fvm.faceOwner()[faceI];
|
||||
label ownDistSqr = distSqr[own];
|
||||
|
||||
const face& f = fvm.faces()[faceI];
|
||||
forAll(f, fp)
|
||||
{
|
||||
pointDistSqr[f[fp]] = max(pointDistSqr[f[fp]], ownDistSqr);
|
||||
}
|
||||
}
|
||||
syncTools::syncPointList
|
||||
(
|
||||
fvm,
|
||||
pointDistSqr,
|
||||
maxEqOp<scalar>(),
|
||||
-sqr(GREAT) // null value
|
||||
);
|
||||
|
||||
pointDistance.internalField() = signedDistance
|
||||
(
|
||||
pointDistSqr,
|
||||
fvm.points(),
|
||||
geometry,
|
||||
surfaces
|
||||
);
|
||||
|
||||
if (writeMesh)
|
||||
{
|
||||
pointDistance.write();
|
||||
}
|
||||
}
|
||||
|
||||
isoSurfaceCell iso
|
||||
(
|
||||
fvm,
|
||||
cellDistance,
|
||||
pointDistance,
|
||||
0, //distance,
|
||||
false //regularise
|
||||
);
|
||||
|
||||
isoFaces.setSize(iso.size());
|
||||
forAll(isoFaces, i)
|
||||
{
|
||||
isoFaces[i] = iso[i].triFaceFace();
|
||||
}
|
||||
isoPoints = iso.points();
|
||||
}
|
||||
|
||||
|
||||
pointField mergedPoints;
|
||||
faceList mergedFaces;
|
||||
labelList pointMergeMap;
|
||||
PatchTools::gatherAndMerge
|
||||
(
|
||||
tolDim,
|
||||
primitivePatch
|
||||
(
|
||||
SubList<face>(isoFaces, isoFaces.size()),
|
||||
isoPoints
|
||||
),
|
||||
mergedPoints,
|
||||
mergedFaces,
|
||||
pointMergeMap
|
||||
);
|
||||
|
||||
vtkSurfaceWriter writer;
|
||||
writer.write
|
||||
(
|
||||
runTime.path(),
|
||||
"iso",
|
||||
mergedPoints,
|
||||
mergedFaces
|
||||
);
|
||||
|
||||
Info<< "End\n" << endl;
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
||||
|
||||
// ************************************************************************* //
|
||||
10
applications/utilities/mesh/generation/cvMesh/cvMeshSurfaceSimplify/Allwmake
Executable file
10
applications/utilities/mesh/generation/cvMesh/cvMeshSurfaceSimplify/Allwmake
Executable file
@ -0,0 +1,10 @@
|
||||
#!/bin/sh
|
||||
cd ${0%/*} || exit 1 # run from this directory
|
||||
set -x
|
||||
|
||||
if [ -d "${FASTDUALOCTREE_SRC_PATH}" ]
|
||||
then
|
||||
wmake
|
||||
fi
|
||||
|
||||
# ----------------------------------------------------------------- end-of-file
|
||||
@ -1,6 +1,22 @@
|
||||
cvMeshSurfaceSimplify.C
|
||||
|
||||
/*
|
||||
cvMeshSurfaceSimplify_non_octree.C
|
||||
MarchingCubes/MarchingCubes.cpp
|
||||
MarchingCubes/ply.c
|
||||
*/
|
||||
|
||||
/*
|
||||
MarchingCubes = fastdualoctree_sgp
|
||||
|
||||
$(MarchingCubes)/data_access.cpp
|
||||
$(MarchingCubes)/fparser.cpp
|
||||
$(MarchingCubes)/fpoptimizer.cpp
|
||||
$(MarchingCubes)/MarchingCubes.cpp
|
||||
$(MarchingCubes)/mc_draw.cpp
|
||||
$(MarchingCubes)/morton.cpp
|
||||
$(MarchingCubes)/opt_octree.cpp
|
||||
$(MarchingCubes)/hash_octree.cpp
|
||||
*/
|
||||
|
||||
cvMeshSurfaceSimplify.C
|
||||
|
||||
EXE = $(FOAM_APPBIN)/cvMeshSurfaceSimplify
|
||||
|
||||
@ -1,7 +1,12 @@
|
||||
MarchingCubes = fastdualoctree_sgp
|
||||
|
||||
include $(GENERAL_RULES)/CGAL
|
||||
|
||||
EXE_INC = \
|
||||
-IMarchingCubes \
|
||||
-DUNIX \
|
||||
-Wno-old-style-cast \
|
||||
/* -IMarchingCubes */ \
|
||||
-I$(FASTDUALOCTREE_SRC_PATH) \
|
||||
-I../conformalVoronoiMesh/lnInclude \
|
||||
-I$(LIB_SRC)/edgeMesh/lnInclude \
|
||||
-I$(LIB_SRC)/triSurface/lnInclude \
|
||||
@ -9,6 +14,8 @@ EXE_INC = \
|
||||
|
||||
EXE_LIBS = \
|
||||
$(CGAL_LIBS) \
|
||||
-L$(FASTDUALOCTREE_SRC_PATH) -lperf_main \
|
||||
-lGL \
|
||||
-lconformalVoronoiMesh \
|
||||
-ldecompositionMethods -L$(FOAM_LIBBIN)/dummy -lscotchDecomp \
|
||||
-ledgeMesh \
|
||||
|
||||
File diff suppressed because it is too large
Load Diff
File diff suppressed because it is too large
Load Diff
@ -1,343 +0,0 @@
|
||||
/**
|
||||
* @file MarchingCubes.h
|
||||
* @author Thomas Lewiner <thomas.lewiner@polytechnique.org>
|
||||
* @author Math Dept, PUC-Rio
|
||||
* @version 0.2
|
||||
* @date 12/08/2002
|
||||
*
|
||||
* @brief MarchingCubes Algorithm
|
||||
*/
|
||||
//________________________________________________
|
||||
|
||||
|
||||
#ifndef _MARCHINGCUBES_H_
|
||||
#define _MARCHINGCUBES_H_
|
||||
|
||||
#if !defined(WIN32) || defined(__CYGWIN__)
|
||||
#pragma interface
|
||||
#endif // WIN32
|
||||
|
||||
|
||||
//_____________________________________________________________________________
|
||||
// types
|
||||
/** unsigned char alias */
|
||||
typedef unsigned char uchar ;
|
||||
/** signed char alias */
|
||||
typedef signed char schar ;
|
||||
/** isovalue alias */
|
||||
typedef float real ;
|
||||
|
||||
//-----------------------------------------------------------------------------
|
||||
// Vertex structure
|
||||
/** \struct Vertex "MarchingCubes.h" MarchingCubes
|
||||
* Position and normal of a vertex
|
||||
* \brief vertex structure
|
||||
* \param x X coordinate
|
||||
* \param y Y coordinate
|
||||
* \param z Z coordinate
|
||||
* \param nx X component of the normal
|
||||
* \param ny Y component of the normal
|
||||
* \param nz Z component of the normal
|
||||
*/
|
||||
typedef struct
|
||||
{
|
||||
real x, y, z ; /**< Vertex coordinates */
|
||||
real nx, ny, nz ; /**< Vertex normal */
|
||||
} Vertex ;
|
||||
|
||||
//-----------------------------------------------------------------------------
|
||||
// Triangle structure
|
||||
/** \struct Triangle "MarchingCubes.h" MarchingCubes
|
||||
* Indices of the oriented triange vertices
|
||||
* \brief triangle structure
|
||||
* \param v1 First vertex index
|
||||
* \param v2 Second vertex index
|
||||
* \param v3 Third vertex index
|
||||
*/
|
||||
typedef struct
|
||||
{
|
||||
int v1,v2,v3 ; /**< Triangle vertices */
|
||||
} Triangle ;
|
||||
//_____________________________________________________________________________
|
||||
|
||||
|
||||
|
||||
//_____________________________________________________________________________
|
||||
/** Marching Cubes algorithm wrapper */
|
||||
/** \class MarchingCubes
|
||||
* \brief Marching Cubes algorithm.
|
||||
*/
|
||||
class MarchingCubes
|
||||
//-----------------------------------------------------------------------------
|
||||
{
|
||||
// Constructors
|
||||
public :
|
||||
/**
|
||||
* Main and default constructor
|
||||
* \brief constructor
|
||||
* \param size_x width of the grid
|
||||
* \param size_y depth of the grid
|
||||
* \param size_z height of the grid
|
||||
*/
|
||||
MarchingCubes ( const int size_x = -1, const int size_y = -1, const int size_z = -1 ) ;
|
||||
/** Destructor */
|
||||
~MarchingCubes() ;
|
||||
|
||||
//-----------------------------------------------------------------------------
|
||||
// Accessors
|
||||
public :
|
||||
/** accesses the number of vertices of the generated mesh */
|
||||
inline const int nverts() const { return _nverts ; }
|
||||
/** accesses the number of triangles of the generated mesh */
|
||||
inline const int ntrigs() const { return _ntrigs ; }
|
||||
/** accesses a specific vertex of the generated mesh */
|
||||
inline Vertex * vert( const int i ) const { if( i < 0 || i >= _nverts ) return ( Vertex *)NULL ; return _vertices + i ; }
|
||||
/** accesses a specific triangle of the generated mesh */
|
||||
inline Triangle * trig( const int i ) const { if( i < 0 || i >= _ntrigs ) return (Triangle*)NULL ; return _triangles + i ; }
|
||||
|
||||
/** accesses the vertex buffer of the generated mesh */
|
||||
inline Vertex *vertices () { return _vertices ; }
|
||||
/** accesses the triangle buffer of the generated mesh */
|
||||
inline Triangle *triangles() { return _triangles ; }
|
||||
|
||||
/** accesses the width of the grid */
|
||||
inline const int size_x() const { return _size_x ; }
|
||||
/** accesses the depth of the grid */
|
||||
inline const int size_y() const { return _size_y ; }
|
||||
/** accesses the height of the grid */
|
||||
inline const int size_z() const { return _size_z ; }
|
||||
|
||||
/**
|
||||
* changes the size of the grid
|
||||
* \param size_x width of the grid
|
||||
* \param size_y depth of the grid
|
||||
* \param size_z height of the grid
|
||||
*/
|
||||
inline void set_resolution( const int size_x, const int size_y, const int size_z ) { _size_x = size_x ; _size_y = size_y ; _size_z = size_z ; }
|
||||
/**
|
||||
* selects wether the algorithm will use the enhanced topologically controlled lookup table or the original MarchingCubes
|
||||
* \param originalMC true for the original Marching Cubes
|
||||
*/
|
||||
inline void set_method ( const bool originalMC = false ) { _originalMC = originalMC ; }
|
||||
/**
|
||||
* selects to use data from another class
|
||||
* \param data is the pointer to the external data, allocated as a size_x*size_y*size_z vector running in x first
|
||||
*/
|
||||
inline void set_ext_data ( real *data )
|
||||
{ if( !_ext_data ) delete [] _data ; _ext_data = data != NULL ; if( _ext_data ) _data = data ; }
|
||||
/**
|
||||
* selects to allocate data
|
||||
*/
|
||||
inline void set_int_data () { _ext_data = false ; _data = NULL ; }
|
||||
|
||||
// Data access
|
||||
/**
|
||||
* accesses a specific cube of the grid
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
inline const real get_data ( const int i, const int j, const int k ) const { return _data[ i + j*_size_x + k*_size_x*_size_y] ; }
|
||||
/**
|
||||
* sets a specific cube of the grid
|
||||
* \param val new value for the cube
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
inline void set_data ( const real val, const int i, const int j, const int k ) { _data[ i + j*_size_x + k*_size_x*_size_y] = val ; }
|
||||
|
||||
// Data initialization
|
||||
/** inits temporary structures (must set sizes before call) : the grid and the vertex index per cube */
|
||||
void init_temps () ;
|
||||
/** inits all structures (must set sizes before call) : the temporary structures and the mesh buffers */
|
||||
void init_all () ;
|
||||
/** clears temporary structures : the grid and the main */
|
||||
void clean_temps() ;
|
||||
/** clears all structures : the temporary structures and the mesh buffers */
|
||||
void clean_all () ;
|
||||
|
||||
|
||||
//-----------------------------------------------------------------------------
|
||||
// Exportation
|
||||
public :
|
||||
/**
|
||||
* PLY exportation of the generated mesh
|
||||
* \param fn name of the PLY file to create
|
||||
* \param bin if true, the PLY will be written in binary mode
|
||||
*/
|
||||
void writePLY( const char *fn, bool bin = false ) ;
|
||||
|
||||
/**
|
||||
* PLY importation of a mesh
|
||||
* \param fn name of the PLY file to read from
|
||||
*/
|
||||
void readPLY( const char *fn ) ;
|
||||
|
||||
/**
|
||||
* VRML / Open Inventor exportation of the generated mesh
|
||||
* \param fn name of the IV file to create
|
||||
*/
|
||||
void writeIV ( const char *fn ) ;
|
||||
|
||||
/**
|
||||
* ISO exportation of the input grid
|
||||
* \param fn name of the ISO file to create
|
||||
*/
|
||||
void writeISO( const char *fn ) ;
|
||||
|
||||
|
||||
//-----------------------------------------------------------------------------
|
||||
// Algorithm
|
||||
public :
|
||||
/**
|
||||
* Main algorithm : must be called after init_all
|
||||
* \param iso isovalue
|
||||
*/
|
||||
void run( real iso = (real)0.0 ) ;
|
||||
|
||||
protected :
|
||||
/** tesselates one cube */
|
||||
void process_cube () ;
|
||||
/** tests if the components of the tesselation of the cube should be connected by the interior of an ambiguous face */
|
||||
bool test_face ( schar face ) ;
|
||||
/** tests if the components of the tesselation of the cube should be connected through the interior of the cube */
|
||||
bool test_interior( schar s ) ;
|
||||
|
||||
|
||||
//-----------------------------------------------------------------------------
|
||||
// Operations
|
||||
protected :
|
||||
/**
|
||||
* computes almost all the vertices of the mesh by interpolation along the cubes edges
|
||||
* \param iso isovalue
|
||||
*/
|
||||
void compute_intersection_points( real iso ) ;
|
||||
|
||||
/**
|
||||
* routine to add a triangle to the mesh
|
||||
* \param trig the code for the triangle as a sequence of edges index
|
||||
* \param n the number of triangles to produce
|
||||
* \param v12 the index of the interior vertex to use, if necessary
|
||||
*/
|
||||
void add_triangle ( const char* trig, char n, int v12 = -1 ) ;
|
||||
|
||||
/** tests and eventually doubles the vertex buffer capacity for a new vertex insertion */
|
||||
void test_vertex_addition() ;
|
||||
/** adds a vertex on the current horizontal edge */
|
||||
int add_x_vertex() ;
|
||||
/** adds a vertex on the current longitudinal edge */
|
||||
int add_y_vertex() ;
|
||||
/** adds a vertex on the current vertical edge */
|
||||
int add_z_vertex() ;
|
||||
/** adds a vertex inside the current cube */
|
||||
int add_c_vertex() ;
|
||||
|
||||
/**
|
||||
* interpolates the horizontal gradient of the implicit function at the lower vertex of the specified cube
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
real get_x_grad( const int i, const int j, const int k ) const ;
|
||||
/**
|
||||
* interpolates the longitudinal gradient of the implicit function at the lower vertex of the specified cube
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
real get_y_grad( const int i, const int j, const int k ) const ;
|
||||
/**
|
||||
* interpolates the vertical gradient of the implicit function at the lower vertex of the specified cube
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
real get_z_grad( const int i, const int j, const int k ) const ;
|
||||
|
||||
/**
|
||||
* accesses the pre-computed vertex index on the lower horizontal edge of a specific cube
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
inline int get_x_vert( const int i, const int j, const int k ) const { return _x_verts[ i + j*_size_x + k*_size_x*_size_y] ; }
|
||||
/**
|
||||
* accesses the pre-computed vertex index on the lower longitudinal edge of a specific cube
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
inline int get_y_vert( const int i, const int j, const int k ) const { return _y_verts[ i + j*_size_x + k*_size_x*_size_y] ; }
|
||||
/**
|
||||
* accesses the pre-computed vertex index on the lower vertical edge of a specific cube
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
inline int get_z_vert( const int i, const int j, const int k ) const { return _z_verts[ i + j*_size_x + k*_size_x*_size_y] ; }
|
||||
|
||||
/**
|
||||
* sets the pre-computed vertex index on the lower horizontal edge of a specific cube
|
||||
* \param val the index of the new vertex
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
inline void set_x_vert( const int val, const int i, const int j, const int k ) { _x_verts[ i + j*_size_x + k*_size_x*_size_y] = val ; }
|
||||
/**
|
||||
* sets the pre-computed vertex index on the lower longitudinal edge of a specific cube
|
||||
* \param val the index of the new vertex
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
inline void set_y_vert( const int val, const int i, const int j, const int k ) { _y_verts[ i + j*_size_x + k*_size_x*_size_y] = val ; }
|
||||
/**
|
||||
* sets the pre-computed vertex index on the lower vertical edge of a specific cube
|
||||
* \param val the index of the new vertex
|
||||
* \param i abscisse of the cube
|
||||
* \param j ordinate of the cube
|
||||
* \param k height of the cube
|
||||
*/
|
||||
inline void set_z_vert( const int val, const int i, const int j, const int k ) { _z_verts[ i + j*_size_x + k*_size_x*_size_y] = val ; }
|
||||
|
||||
/** prints cube for debug */
|
||||
void print_cube() ;
|
||||
|
||||
//-----------------------------------------------------------------------------
|
||||
// Elements
|
||||
protected :
|
||||
bool _originalMC ; /**< selects wether the algorithm will use the enhanced topologically controlled lookup table or the original MarchingCubes */
|
||||
bool _ext_data ; /**< selects wether to allocate data or use data from another class */
|
||||
|
||||
int _size_x ; /**< width of the grid */
|
||||
int _size_y ; /**< depth of the grid */
|
||||
int _size_z ; /**< height of the grid */
|
||||
real *_data ; /**< implicit function values sampled on the grid */
|
||||
|
||||
int *_x_verts ; /**< pre-computed vertex indices on the lower horizontal edge of each cube */
|
||||
int *_y_verts ; /**< pre-computed vertex indices on the lower longitudinal edge of each cube */
|
||||
int *_z_verts ; /**< pre-computed vertex indices on the lower vertical edge of each cube */
|
||||
|
||||
int _nverts ; /**< number of allocated vertices in the vertex buffer */
|
||||
int _ntrigs ; /**< number of allocated triangles in the triangle buffer */
|
||||
int _Nverts ; /**< size of the vertex buffer */
|
||||
int _Ntrigs ; /**< size of the triangle buffer */
|
||||
Vertex *_vertices ; /**< vertex buffer */
|
||||
Triangle *_triangles ; /**< triangle buffer */
|
||||
|
||||
int _i ; /**< abscisse of the active cube */
|
||||
int _j ; /**< height of the active cube */
|
||||
int _k ; /**< ordinate of the active cube */
|
||||
|
||||
real _cube[8] ; /**< values of the implicit function on the active cube */
|
||||
uchar _lut_entry ; /**< cube sign representation in [0..255] */
|
||||
uchar _case ; /**< case of the active cube in [0..15] */
|
||||
uchar _config ; /**< configuration of the active cube */
|
||||
uchar _subconfig ; /**< subconfiguration of the active cube */
|
||||
};
|
||||
//_____________________________________________________________________________
|
||||
|
||||
|
||||
#endif // _MARCHINGCUBES_H_
|
||||
File diff suppressed because it is too large
Load Diff
@ -1,233 +0,0 @@
|
||||
/*
|
||||
|
||||
Header for PLY polygon files.
|
||||
|
||||
- Greg Turk
|
||||
|
||||
A PLY file contains a single polygonal _object_.
|
||||
|
||||
An object is composed of lists of _elements_. Typical elements are
|
||||
vertices, faces, edges and materials.
|
||||
|
||||
Each type of element for a given object has one or more _properties_
|
||||
associated with the element type. For instance, a vertex element may
|
||||
have as properties three floating-point values x,y,z and three unsigned
|
||||
chars for red, green and blue.
|
||||
|
||||
-----------------------------------------------------------------------
|
||||
|
||||
Copyright (c) 1998 Georgia Institute of Technology. All rights reserved.
|
||||
|
||||
Permission to use, copy, modify and distribute this software and its
|
||||
documentation for any purpose is hereby granted without fee, provided
|
||||
that the above copyright notice and this permission notice appear in
|
||||
all copies of this software and that you do not sell the software.
|
||||
|
||||
THE SOFTWARE IS PROVIDED "AS IS" AND WITHOUT WARRANTY OF ANY KIND,
|
||||
EXPRESS, IMPLIED OR OTHERWISE, INCLUDING WITHOUT LIMITATION, ANY
|
||||
WARRANTY OF MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE.
|
||||
|
||||
*/
|
||||
|
||||
#ifndef __PLY_H__
|
||||
#define __PLY_H__
|
||||
|
||||
#ifdef __cplusplus
|
||||
extern "C" {
|
||||
#endif
|
||||
|
||||
#include <stdio.h>
|
||||
#include <stddef.h>
|
||||
|
||||
#define PLY_ASCII 1 /* ascii PLY file */
|
||||
#define PLY_BINARY_BE 2 /* binary PLY file, big endian */
|
||||
#define PLY_BINARY_LE 3 /* binary PLY file, little endian */
|
||||
|
||||
#define PLY_OKAY 0 /* ply routine worked okay */
|
||||
#define PLY_ERROR -1 /* error in ply routine */
|
||||
|
||||
/* scalar data types supported by PLY format */
|
||||
|
||||
#define StartType 0
|
||||
#define Int8 1
|
||||
#define Int16 2
|
||||
#define Int32 3
|
||||
#define Uint8 4
|
||||
#define Uint16 5
|
||||
#define Uint32 6
|
||||
#define Float32 7
|
||||
#define Float64 8
|
||||
#define EndType 9
|
||||
|
||||
#define PLY_SCALAR 0
|
||||
#define PLY_LIST 1
|
||||
#define PLY_STRING 2
|
||||
|
||||
|
||||
typedef struct PlyProperty { /* description of a property */
|
||||
|
||||
char *name; /* property name */
|
||||
int external_type; /* file's data type */
|
||||
int internal_type; /* program's data type */
|
||||
int offset; /* offset bytes of prop in a struct */
|
||||
|
||||
int is_list; /* 0 = scalar, 1 = list, 2 = char string */
|
||||
int count_external; /* file's count type */
|
||||
int count_internal; /* program's count type */
|
||||
int count_offset; /* offset byte for list count */
|
||||
|
||||
} PlyProperty;
|
||||
|
||||
typedef struct PlyElement { /* description of an element */
|
||||
char *name; /* element name */
|
||||
int num; /* number of elements in this object */
|
||||
int size; /* size of element (bytes) or -1 if variable */
|
||||
int nprops; /* number of properties for this element */
|
||||
PlyProperty **props; /* list of properties in the file */
|
||||
char *store_prop; /* flags: property wanted by user? */
|
||||
int other_offset; /* offset to un-asked-for props, or -1 if none*/
|
||||
int other_size; /* size of other_props structure */
|
||||
} PlyElement;
|
||||
|
||||
typedef struct PlyOtherProp { /* describes other properties in an element */
|
||||
char *name; /* element name */
|
||||
int size; /* size of other_props */
|
||||
int nprops; /* number of properties in other_props */
|
||||
PlyProperty **props; /* list of properties in other_props */
|
||||
} PlyOtherProp;
|
||||
|
||||
typedef struct OtherData { /* for storing other_props for an other element */
|
||||
void *other_props;
|
||||
} OtherData;
|
||||
|
||||
typedef struct OtherElem { /* data for one "other" element */
|
||||
char *elem_name; /* names of other elements */
|
||||
int elem_count; /* count of instances of each element */
|
||||
OtherData **other_data; /* actual property data for the elements */
|
||||
PlyOtherProp *other_props; /* description of the property data */
|
||||
} OtherElem;
|
||||
|
||||
typedef struct PlyOtherElems { /* "other" elements, not interpreted by user */
|
||||
int num_elems; /* number of other elements */
|
||||
OtherElem *other_list; /* list of data for other elements */
|
||||
} PlyOtherElems;
|
||||
|
||||
#define AVERAGE_RULE 1
|
||||
#define MAJORITY_RULE 2
|
||||
#define MINIMUM_RULE 3
|
||||
#define MAXIMUM_RULE 4
|
||||
#define SAME_RULE 5
|
||||
#define RANDOM_RULE 6
|
||||
|
||||
typedef struct PlyPropRules { /* rules for combining "other" properties */
|
||||
PlyElement *elem; /* element whose rules we are making */
|
||||
int *rule_list; /* types of rules (AVERAGE_PLY, MAJORITY_PLY, etc.) */
|
||||
int nprops; /* number of properties we're combining so far */
|
||||
int max_props; /* maximum number of properties we have room for now */
|
||||
void **props; /* list of properties we're combining */
|
||||
float *weights; /* list of weights of the properties */
|
||||
} PlyPropRules;
|
||||
|
||||
typedef struct PlyRuleList {
|
||||
char *name; /* name of the rule */
|
||||
char *element; /* name of element that rule applies to */
|
||||
char *property; /* name of property that rule applies to */
|
||||
struct PlyRuleList *next; /* pointer for linked list of rules */
|
||||
} PlyRuleList;
|
||||
|
||||
typedef struct PlyFile { /* description of PLY file */
|
||||
FILE *fp; /* file pointer */
|
||||
int file_type; /* ascii or binary */
|
||||
float version; /* version number of file */
|
||||
int num_elem_types; /* number of element types of object */
|
||||
PlyElement **elems; /* list of elements */
|
||||
int num_comments; /* number of comments */
|
||||
char **comments; /* list of comments */
|
||||
int num_obj_info; /* number of items of object information */
|
||||
char **obj_info; /* list of object info items */
|
||||
PlyElement *which_elem; /* element we're currently reading or writing */
|
||||
PlyOtherElems *other_elems; /* "other" elements from a PLY file */
|
||||
PlyPropRules *current_rules; /* current propagation rules */
|
||||
PlyRuleList *rule_list; /* rule list from user */
|
||||
} PlyFile;
|
||||
|
||||
/* memory allocation */
|
||||
/*
|
||||
extern char *my_alloc();
|
||||
*/
|
||||
#define myalloc(mem_size) my_alloc((mem_size), __LINE__, __FILE__)
|
||||
|
||||
|
||||
/* old routines */
|
||||
|
||||
#if 0
|
||||
extern PlyFile *ply_write(FILE *, int, char **, int);
|
||||
extern PlyFile *ply_read(FILE *, int *, char ***);
|
||||
extern PlyFile *ply_open_for_reading( char *, int *, char ***, int *, float *);
|
||||
extern void ply_close(PlyFile *);
|
||||
extern PlyOtherProp *ply_get_other_properties(PlyFile *, char *, int);
|
||||
#endif
|
||||
|
||||
extern void ply_describe_property( PlyFile * , char * , PlyProperty * );
|
||||
extern void ply_get_property( PlyFile * , char * , PlyProperty * );
|
||||
extern void ply_get_element( PlyFile * , void * );
|
||||
|
||||
|
||||
/*** delcaration of routines ***/
|
||||
|
||||
PlyOtherElems *get_other_element_ply( PlyFile * );
|
||||
|
||||
PlyFile *read_ply( FILE * );
|
||||
PlyFile *write_ply( FILE * , int, char ** , int );
|
||||
extern PlyFile *open_for_writing_ply( char * , int, char ** , int );
|
||||
void close_ply( PlyFile * );
|
||||
void free_ply( PlyFile * );
|
||||
|
||||
void get_info_ply( PlyFile * , float * , int * );
|
||||
void free_other_elements_ply( PlyOtherElems * );
|
||||
|
||||
void append_comment_ply( PlyFile * , char * );
|
||||
void append_obj_info_ply( PlyFile * , char * );
|
||||
void copy_comments_ply( PlyFile * , PlyFile * );
|
||||
void copy_obj_info_ply( PlyFile * , PlyFile * );
|
||||
char* *get_comments_ply( PlyFile * , int * );
|
||||
char* *get_obj_info_ply( PlyFile * , int * );
|
||||
|
||||
char* *get_element_list_ply( PlyFile * , int * );
|
||||
int setup_property_ply( PlyFile * , PlyProperty * );
|
||||
void get_element_ply( PlyFile * , void * );
|
||||
char *setup_element_read_ply( PlyFile * , int, int * );
|
||||
PlyOtherProp *get_other_properties_ply( PlyFile * , int );
|
||||
|
||||
void element_count_ply( PlyFile * , char * , int );
|
||||
void describe_element_ply( PlyFile * , char * , int );
|
||||
void describe_property_ply( PlyFile * , PlyProperty * );
|
||||
void describe_other_properties_ply( PlyFile * , PlyOtherProp * , int );
|
||||
void describe_other_elements_ply( PlyFile * , PlyOtherElems * );
|
||||
void get_element_setup_ply( PlyFile * , char * , int, PlyProperty * );
|
||||
PlyProperty* *get_element_description_ply( PlyFile * , char * , int * , int * );
|
||||
void element_layout_ply( PlyFile * , char * , int, int, PlyProperty * );
|
||||
|
||||
void header_complete_ply( PlyFile * );
|
||||
void put_element_setup_ply( PlyFile * , char * );
|
||||
void put_element_ply( PlyFile * , void * );
|
||||
void put_other_elements_ply( PlyFile * );
|
||||
|
||||
PlyPropRules *init_rule_ply( PlyFile * , char * );
|
||||
void modify_rule_ply( PlyPropRules * , char * , int );
|
||||
void start_props_ply( PlyFile * , PlyPropRules * );
|
||||
void weight_props_ply( PlyFile * , float, void * );
|
||||
void *get_new_props_ply( PlyFile * );
|
||||
void set_prop_rules_ply( PlyFile * , PlyRuleList * );
|
||||
PlyRuleList *append_prop_rule( PlyRuleList * , char * , char * );
|
||||
int matches_rule_name( char * );
|
||||
|
||||
int equal_strings( char * , char * );
|
||||
char *recreate_command_line( int, char *argv[] );
|
||||
|
||||
|
||||
#ifdef __cplusplus
|
||||
}
|
||||
#endif
|
||||
#endif /* !__PLY_H__ */
|
||||
|
||||
@ -28,6 +28,7 @@ Description
|
||||
Simplifies surfaces by resampling.
|
||||
|
||||
Uses Thomas Lewiner's topology preserving MarchingCubes.
|
||||
(http://zeus.mat.puc-rio.br/tomlew/tomlew_uk.php)
|
||||
|
||||
\*---------------------------------------------------------------------------*/
|
||||
|
||||
@ -37,13 +38,328 @@ Description
|
||||
#include "conformationSurfaces.H"
|
||||
#include "triSurfaceMesh.H"
|
||||
|
||||
#include "MarchingCubes.h"
|
||||
|
||||
#include "opt_octree.h"
|
||||
#include "cube.h"
|
||||
|
||||
using namespace Foam;
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
class pointConversion
|
||||
{
|
||||
const vector scale_;
|
||||
|
||||
const vector offset_;
|
||||
|
||||
public:
|
||||
|
||||
//- Construct from components
|
||||
pointConversion
|
||||
(
|
||||
const vector scale,
|
||||
const vector offset
|
||||
)
|
||||
:
|
||||
scale_(scale),
|
||||
offset_(offset)
|
||||
{}
|
||||
|
||||
inline Point toLocal(const Foam::point& pt) const
|
||||
{
|
||||
Foam::point p = cmptMultiply(scale_, (pt + offset_));
|
||||
return Point(p.x(), p.y(), p.z());
|
||||
}
|
||||
|
||||
inline Foam::point toGlobal(const Point& pt) const
|
||||
{
|
||||
point p(pt.x(), pt.y(), pt.z());
|
||||
return cmptDivide(p, scale_) - offset_;
|
||||
}
|
||||
};
|
||||
|
||||
|
||||
|
||||
|
||||
// For use in Fast-Dual Octree from Thomas Lewiner
|
||||
class distanceCalc
|
||||
:
|
||||
public ::data_access
|
||||
{
|
||||
|
||||
const Level min_level_;
|
||||
|
||||
const conformationSurfaces& geometryToConformTo_;
|
||||
|
||||
const pointConversion& converter_;
|
||||
|
||||
|
||||
// Private Member Functions
|
||||
|
||||
scalar signedDistance(const Foam::point& pt) const
|
||||
{
|
||||
const searchableSurfaces& geometry = geometryToConformTo_.geometry();
|
||||
const labelList& surfaces = geometryToConformTo_.surfaces();
|
||||
|
||||
static labelList nearestSurfaces;
|
||||
static scalarField distance;
|
||||
|
||||
static pointField samples(1);
|
||||
samples[0] = pt;
|
||||
|
||||
searchableSurfacesQueries::signedDistance
|
||||
(
|
||||
geometry,
|
||||
surfaces,
|
||||
samples,
|
||||
scalarField(1, GREAT),
|
||||
searchableSurface::OUTSIDE,
|
||||
nearestSurfaces,
|
||||
distance
|
||||
);
|
||||
|
||||
return distance[0];
|
||||
}
|
||||
|
||||
|
||||
public:
|
||||
|
||||
// Constructors
|
||||
|
||||
//- Construct from components
|
||||
distanceCalc
|
||||
(
|
||||
Level max_level_,
|
||||
real iso_val_,
|
||||
Level min_level,
|
||||
const conformationSurfaces& geometryToConformTo,
|
||||
const pointConversion& converter
|
||||
)
|
||||
:
|
||||
data_access(max_level_,iso_val_),
|
||||
min_level_(min_level),
|
||||
geometryToConformTo_(geometryToConformTo),
|
||||
converter_(converter)
|
||||
{}
|
||||
|
||||
|
||||
//- Destructor
|
||||
virtual ~distanceCalc()
|
||||
{}
|
||||
|
||||
|
||||
// Member Functions
|
||||
|
||||
//- test function
|
||||
virtual bool need_refine( const Cube &c )
|
||||
{
|
||||
int l = c.lv() ;
|
||||
|
||||
if( l >= _max_level ) return false;
|
||||
if( l < min_level_ ) return true;
|
||||
|
||||
treeBoundBox bb
|
||||
(
|
||||
converter_.toGlobal
|
||||
(
|
||||
Point
|
||||
(
|
||||
c.xmin(),
|
||||
c.ymin(),
|
||||
c.zmin()
|
||||
)
|
||||
),
|
||||
converter_.toGlobal
|
||||
(
|
||||
Point
|
||||
(
|
||||
c.xmax(),
|
||||
c.ymax(),
|
||||
c.zmax()
|
||||
)
|
||||
)
|
||||
);
|
||||
|
||||
const searchableSurfaces& geometry =
|
||||
geometryToConformTo_.geometry();
|
||||
const labelList& surfaces =
|
||||
geometryToConformTo_.surfaces();
|
||||
|
||||
|
||||
//- Uniform refinement around surface
|
||||
{
|
||||
forAll(surfaces, i)
|
||||
{
|
||||
if (geometry[surfaces[i]].overlaps(bb))
|
||||
{
|
||||
return true;
|
||||
}
|
||||
}
|
||||
return false;
|
||||
}
|
||||
|
||||
|
||||
////- Surface intersects bb (but not using intersection test)
|
||||
//scalar ccDist = signedDistance(bb.midpoint());
|
||||
//scalar ccVal = ccDist - _iso_val;
|
||||
//if (mag(ccVal) < SMALL)
|
||||
//{
|
||||
// return true;
|
||||
//}
|
||||
//const pointField points(bb.points());
|
||||
//forAll(points, pointI)
|
||||
//{
|
||||
// scalar pointVal = signedDistance(points[pointI]) - _iso_val;
|
||||
// if (ccVal*pointVal < 0)
|
||||
// {
|
||||
// return true;
|
||||
// }
|
||||
//}
|
||||
//return false;
|
||||
|
||||
|
||||
////- Refinement based on intersection with multiple planes.
|
||||
//// Does not work well - too high a ratio between
|
||||
//// neighbouring cubes.
|
||||
//const pointField points(bb.points());
|
||||
//const edgeList& edges = treeBoundBox::edges;
|
||||
//pointField start(edges.size());
|
||||
//pointField end(edges.size());
|
||||
//forAll(edges, i)
|
||||
//{
|
||||
// start[i] = points[edges[i][0]];
|
||||
// end[i] = points[edges[i][1]];
|
||||
//}
|
||||
//Foam::List<Foam::List<pointIndexHit> > hitInfo;
|
||||
//labelListList hitSurfaces;
|
||||
//searchableSurfacesQueries::findAllIntersections
|
||||
//(
|
||||
// geometry,
|
||||
// surfaces,
|
||||
// start,
|
||||
// end,
|
||||
// hitSurfaces,
|
||||
// hitInfo
|
||||
//);
|
||||
//
|
||||
//// Count number of intersections
|
||||
//label nInt = 0;
|
||||
//forAll(hitSurfaces, edgeI)
|
||||
//{
|
||||
// nInt += hitSurfaces[edgeI].size();
|
||||
//}
|
||||
//
|
||||
//if (nInt == 0)
|
||||
//{
|
||||
// // No surface intersected. See if there is one inside
|
||||
// forAll(surfaces, i)
|
||||
// {
|
||||
// if (geometry[surfaces[i]].overlaps(bb))
|
||||
// {
|
||||
// return true;
|
||||
// }
|
||||
// }
|
||||
// return false;
|
||||
//}
|
||||
//
|
||||
//// Check multiple surfaces
|
||||
//label baseSurfI = -1;
|
||||
//forAll(hitSurfaces, edgeI)
|
||||
//{
|
||||
// const labelList& hSurfs = hitSurfaces[edgeI];
|
||||
// forAll(hSurfs, i)
|
||||
// {
|
||||
// if (baseSurfI == -1)
|
||||
// {
|
||||
// baseSurfI = hSurfs[i];
|
||||
// }
|
||||
// else if (baseSurfI != hSurfs[i])
|
||||
// {
|
||||
// // Multiple surfaces
|
||||
// return true;
|
||||
// }
|
||||
// }
|
||||
//}
|
||||
//
|
||||
//// Get normals
|
||||
//DynamicList<pointIndexHit> baseInfo(nInt);
|
||||
//forAll(hitInfo, edgeI)
|
||||
//{
|
||||
// const Foam::List<pointIndexHit>& hits = hitInfo[edgeI];
|
||||
// forAll(hits, i)
|
||||
// {
|
||||
// (void)hits[i].hitPoint();
|
||||
// baseInfo.append(hits[i]);
|
||||
// }
|
||||
//}
|
||||
//vectorField normals;
|
||||
//geometry[surfaces[baseSurfI]].getNormal(baseInfo, normals);
|
||||
//for (label i = 1; i < normals.size(); ++i)
|
||||
//{
|
||||
// if ((normals[0] & normals[i]) < 0.9)
|
||||
// {
|
||||
// return true;
|
||||
// }
|
||||
//}
|
||||
//labelList regions;
|
||||
//geometry[surfaces[baseSurfI]].getRegion(baseInfo, regions);
|
||||
//for (label i = 1; i < regions.size(); ++i)
|
||||
//{
|
||||
// if (regions[0] != regions[i])
|
||||
// {
|
||||
// return true;
|
||||
// }
|
||||
//}
|
||||
//return false;
|
||||
|
||||
|
||||
|
||||
//samples[0] = point(c.xmin(), c.ymin(), c.zmin());
|
||||
//samples[1] = point(c.xmax(), c.ymin(), c.zmin());
|
||||
//samples[2] = point(c.xmax(), c.ymax(), c.zmin());
|
||||
//samples[3] = point(c.xmin(), c.ymax(), c.zmin());
|
||||
//
|
||||
//samples[4] = point(c.xmin(), c.ymin(), c.zmax());
|
||||
//samples[5] = point(c.xmax(), c.ymin(), c.zmax());
|
||||
//samples[6] = point(c.xmax(), c.ymax(), c.zmax());
|
||||
//samples[7] = point(c.xmin(), c.ymax(), c.zmax());
|
||||
|
||||
//scalarField nearestDistSqr(8, GREAT);
|
||||
//
|
||||
//Foam::List<pointIndexHit> nearestInfo;
|
||||
//surf_.findNearest(samples, nearestDistSqr, nearestInfo);
|
||||
//vectorField normals;
|
||||
//surf_.getNormal(nearestInfo, normals);
|
||||
//
|
||||
//for (label i = 1; i < normals.size(); ++i)
|
||||
//{
|
||||
// if ((normals[0] & normals[i]) < 0.5)
|
||||
// {
|
||||
// return true;
|
||||
// }
|
||||
//}
|
||||
//return false;
|
||||
|
||||
//// Check if surface octree same level
|
||||
//const labelList elems(surf_.tree().findBox(bb));
|
||||
//
|
||||
//if (elems.size() > 1)
|
||||
//{
|
||||
// return true;
|
||||
//}
|
||||
//else
|
||||
//{
|
||||
// return false;
|
||||
//}
|
||||
}
|
||||
|
||||
//- data function
|
||||
virtual real value_at( const Cube &c )
|
||||
{
|
||||
return signedDistance(converter_.toGlobal(c)) - _iso_val;
|
||||
}
|
||||
};
|
||||
|
||||
|
||||
// Main program:
|
||||
|
||||
int main(int argc, char *argv[])
|
||||
@ -52,19 +368,16 @@ int main(int argc, char *argv[])
|
||||
(
|
||||
"Re-sample surfaces used in cvMesh operation"
|
||||
);
|
||||
//argList::validArgs.append("inputFile");
|
||||
argList::validArgs.append("(nx ny nz)");
|
||||
argList::validArgs.append("outputName");
|
||||
|
||||
#include "setRootCase.H"
|
||||
#include "createTime.H"
|
||||
runTime.functionObjects().off();
|
||||
|
||||
const Vector<label> n(IStringStream(args.args()[1])());
|
||||
const fileName exportName = args.args()[2];
|
||||
const fileName exportName = args.args()[1];
|
||||
|
||||
Info<< "Reading surfaces as specified in the cvMeshDict and"
|
||||
<< " writing re-sampled " << n << " to " << exportName
|
||||
<< " writing a re-sampled surface to " << exportName
|
||||
<< nl << endl;
|
||||
|
||||
cpuTime timer;
|
||||
@ -114,124 +427,100 @@ int main(int argc, char *argv[])
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
|
||||
|
||||
// Extend
|
||||
treeBoundBox bb = geometryToConformTo.globalBounds();
|
||||
{
|
||||
const vector smallVec = 0.1*bb.span();
|
||||
bb.min() -= smallVec;
|
||||
bb.max() += smallVec;
|
||||
}
|
||||
|
||||
Info<< "Meshing to bounding box " << bb << nl << endl;
|
||||
|
||||
const vector span(bb.span());
|
||||
const vector d
|
||||
(
|
||||
span.x()/(n.x()-1),
|
||||
span.y()/(n.y()-1),
|
||||
span.z()/(n.z()-1)
|
||||
);
|
||||
|
||||
MarchingCubes mc(span.x(), span.y(), span.z() ) ;
|
||||
mc.set_resolution(n.x(), n.y(), n.z());
|
||||
mc.init_all() ;
|
||||
|
||||
|
||||
// Generate points
|
||||
pointField points(mc.size_x()*mc.size_y()*mc.size_z());
|
||||
label pointI = 0;
|
||||
|
||||
point pt;
|
||||
for( int k = 0 ; k < mc.size_z() ; k++ )
|
||||
{
|
||||
pt.z() = bb.min().z() + k*d.z();
|
||||
for( int j = 0 ; j < mc.size_y() ; j++ )
|
||||
{
|
||||
pt.y() = bb.min().y() + j*d.y();
|
||||
for( int i = 0 ; i < mc.size_x() ; i++ )
|
||||
{
|
||||
pt.x() = bb.min().x() + i*d.x();
|
||||
points[pointI++] = pt;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
Info<< "Generated " << points.size() << " sampling points in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
|
||||
// Compute data
|
||||
const searchableSurfaces& geometry = geometryToConformTo.geometry();
|
||||
const labelList& surfaces = geometryToConformTo.surfaces();
|
||||
|
||||
scalarField signedDist;
|
||||
labelList nearestSurfaces;
|
||||
searchableSurfacesQueries::signedDistance
|
||||
|
||||
const label minLevel = 2;
|
||||
|
||||
// The max cube size follows from the minLevel and the default cube size
|
||||
// (1)
|
||||
const scalar maxSize = 1.0 / (1 << minLevel);
|
||||
const scalar halfMaxSize = 0.5*maxSize;
|
||||
|
||||
|
||||
// Scale the geometry to fit within
|
||||
// halfMaxSize .. 1-halfMaxSize
|
||||
|
||||
scalar wantedRange = 1.0-maxSize;
|
||||
|
||||
const treeBoundBox bb = geometryToConformTo.globalBounds();
|
||||
|
||||
const vector scale = cmptDivide
|
||||
(
|
||||
geometry,
|
||||
surfaces,
|
||||
points,
|
||||
scalarField(points.size(), sqr(GREAT)),
|
||||
nearestSurfaces,
|
||||
signedDist
|
||||
point(wantedRange, wantedRange, wantedRange),
|
||||
bb.span()
|
||||
);
|
||||
const vector offset =
|
||||
cmptDivide
|
||||
(
|
||||
point(halfMaxSize, halfMaxSize, halfMaxSize),
|
||||
scale
|
||||
)
|
||||
-bb.min();
|
||||
|
||||
|
||||
const pointConversion converter(scale, offset);
|
||||
|
||||
|
||||
// Marching cubes
|
||||
|
||||
OptOctree octree;
|
||||
|
||||
distanceCalc ref
|
||||
(
|
||||
8, //maxLevel
|
||||
0.0, //distance
|
||||
minLevel, //minLevel
|
||||
geometryToConformTo,
|
||||
converter
|
||||
);
|
||||
|
||||
// Fill elements
|
||||
pointI = 0;
|
||||
for( int k = 0 ; k < mc.size_z() ; k++ )
|
||||
{
|
||||
for( int j = 0 ; j < mc.size_y() ; j++ )
|
||||
{
|
||||
for( int i = 0 ; i < mc.size_x() ; i++ )
|
||||
{
|
||||
mc.set_data(float(signedDist[pointI++]), i, j, k);
|
||||
}
|
||||
}
|
||||
}
|
||||
octree.refine(&ref);
|
||||
octree.set_impl(&ref);
|
||||
|
||||
Info<< "Determined signed distance in = "
|
||||
Info<< "Calculated octree in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
MarchingCubes& mc = octree.mc();
|
||||
|
||||
mc.run() ;
|
||||
mc.clean_all() ;
|
||||
octree.build_isosurface(&ref) ;
|
||||
|
||||
Info<< "Constructed iso surface in = "
|
||||
Info<< "Constructed iso surface of distance in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
|
||||
mc.clean_temps() ;
|
||||
|
||||
|
||||
|
||||
// Write output file
|
||||
if (mc.ntrigs() > 0)
|
||||
{
|
||||
Triangle* triangles = mc.triangles();
|
||||
List<labelledTri> tris(mc.ntrigs());
|
||||
forAll(tris, triI)
|
||||
label nTris = mc.ntrigs();
|
||||
Foam::DynamicList<labelledTri> tris(mc.ntrigs());
|
||||
for (label triI = 0; triI < nTris; ++triI)
|
||||
{
|
||||
tris[triI] = labelledTri
|
||||
(
|
||||
triangles[triI].v1,
|
||||
triangles[triI].v2,
|
||||
triangles[triI].v3,
|
||||
0 // region
|
||||
);
|
||||
const Triangle& t = triangles[triI];
|
||||
if (t.v1 != t.v2 && t.v1 != t.v3 && t.v2 != t.v3)
|
||||
{
|
||||
tris.append
|
||||
(
|
||||
labelledTri
|
||||
(
|
||||
triangles[triI].v1,
|
||||
triangles[triI].v2,
|
||||
triangles[triI].v3,
|
||||
0 // region
|
||||
)
|
||||
);
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
Vertex* vertices = mc.vertices();
|
||||
Point* vertices = mc.vertices();
|
||||
pointField points(mc.nverts());
|
||||
forAll(points, pointI)
|
||||
{
|
||||
Vertex& v = vertices[pointI];
|
||||
points[pointI] = point
|
||||
(
|
||||
bb.min().x() + v.x*span.x()/mc.size_x(),
|
||||
bb.min().y() + v.y*span.y()/mc.size_y(),
|
||||
bb.min().z() + v.z*span.z()/mc.size_z()
|
||||
);
|
||||
const Point& v = vertices[pointI];
|
||||
points[pointI] = converter.toGlobal(v);
|
||||
}
|
||||
|
||||
|
||||
@ -265,6 +554,7 @@ int main(int argc, char *argv[])
|
||||
}
|
||||
|
||||
triSurface s(tris, patches, points, true);
|
||||
tris.clearStorage();
|
||||
|
||||
Info<< "Extracted triSurface in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
@ -272,7 +562,7 @@ int main(int argc, char *argv[])
|
||||
|
||||
// Find out region on local surface of nearest point
|
||||
{
|
||||
List<pointIndexHit> hitInfo;
|
||||
Foam::List<pointIndexHit> hitInfo;
|
||||
labelList hitSurfaces;
|
||||
geometryToConformTo.findSurfaceNearest
|
||||
(
|
||||
@ -340,7 +630,6 @@ int main(int argc, char *argv[])
|
||||
|
||||
mc.clean_all() ;
|
||||
|
||||
|
||||
Info<< "End\n" << endl;
|
||||
|
||||
return 0;
|
||||
|
||||
@ -0,0 +1,352 @@
|
||||
/*---------------------------------------------------------------------------*\
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
This file is part of OpenFOAM.
|
||||
|
||||
OpenFOAM is free software: you can redistribute it and/or modify it
|
||||
under the terms of the GNU General Public License as published by
|
||||
the Free Software Foundation, either version 3 of the License, or
|
||||
(at your option) any later version.
|
||||
|
||||
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
|
||||
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
|
||||
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
|
||||
for more details.
|
||||
|
||||
You should have received a copy of the GNU General Public License
|
||||
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
Application
|
||||
cvMeshSurfaceSimplify
|
||||
|
||||
Description
|
||||
Simplifies surfaces by resampling.
|
||||
|
||||
Uses Thomas Lewiner's topology preserving MarchingCubes.
|
||||
|
||||
\*---------------------------------------------------------------------------*/
|
||||
|
||||
#include "argList.H"
|
||||
#include "Time.H"
|
||||
#include "searchableSurfaces.H"
|
||||
#include "conformationSurfaces.H"
|
||||
#include "triSurfaceMesh.H"
|
||||
|
||||
#include "MarchingCubes.h"
|
||||
|
||||
|
||||
using namespace Foam;
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
// Main program:
|
||||
|
||||
int main(int argc, char *argv[])
|
||||
{
|
||||
argList::addNote
|
||||
(
|
||||
"Re-sample surfaces used in cvMesh operation"
|
||||
);
|
||||
//argList::validArgs.append("inputFile");
|
||||
argList::validArgs.append("(nx ny nz)");
|
||||
argList::validArgs.append("outputName");
|
||||
|
||||
#include "setRootCase.H"
|
||||
#include "createTime.H"
|
||||
runTime.functionObjects().off();
|
||||
|
||||
const Vector<label> n(IStringStream(args.args()[1])());
|
||||
const fileName exportName = args.args()[2];
|
||||
|
||||
Info<< "Reading surfaces as specified in the cvMeshDict and"
|
||||
<< " writing re-sampled " << n << " to " << exportName
|
||||
<< nl << endl;
|
||||
|
||||
cpuTime timer;
|
||||
|
||||
IOdictionary cvMeshDict
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"cvMeshDict",
|
||||
runTime.system(),
|
||||
runTime,
|
||||
IOobject::MUST_READ_IF_MODIFIED,
|
||||
IOobject::NO_WRITE
|
||||
)
|
||||
);
|
||||
|
||||
// Define/load all geometry
|
||||
searchableSurfaces allGeometry
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
"cvSearchableSurfaces",
|
||||
runTime.constant(),
|
||||
"triSurface",
|
||||
runTime,
|
||||
IOobject::MUST_READ,
|
||||
IOobject::NO_WRITE
|
||||
),
|
||||
cvMeshDict.subDict("geometry")
|
||||
);
|
||||
|
||||
Info<< "Geometry read in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
|
||||
Random rndGen(64293*Pstream::myProcNo());
|
||||
|
||||
conformationSurfaces geometryToConformTo
|
||||
(
|
||||
runTime,
|
||||
rndGen,
|
||||
allGeometry,
|
||||
cvMeshDict.subDict("surfaceConformation")
|
||||
);
|
||||
|
||||
Info<< "Set up geometry in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
|
||||
|
||||
// Extend
|
||||
treeBoundBox bb = geometryToConformTo.globalBounds();
|
||||
{
|
||||
const vector smallVec = 0.1*bb.span();
|
||||
bb.min() -= smallVec;
|
||||
bb.max() += smallVec;
|
||||
}
|
||||
|
||||
Info<< "Meshing to bounding box " << bb << nl << endl;
|
||||
|
||||
const vector span(bb.span());
|
||||
const vector d
|
||||
(
|
||||
span.x()/(n.x()-1),
|
||||
span.y()/(n.y()-1),
|
||||
span.z()/(n.z()-1)
|
||||
);
|
||||
|
||||
MarchingCubes mc(span.x(), span.y(), span.z() ) ;
|
||||
mc.set_resolution(n.x(), n.y(), n.z());
|
||||
mc.init_all() ;
|
||||
|
||||
|
||||
// Generate points
|
||||
pointField points(mc.size_x()*mc.size_y()*mc.size_z());
|
||||
label pointI = 0;
|
||||
|
||||
point pt;
|
||||
for( int k = 0 ; k < mc.size_z() ; k++ )
|
||||
{
|
||||
pt.z() = bb.min().z() + k*d.z();
|
||||
for( int j = 0 ; j < mc.size_y() ; j++ )
|
||||
{
|
||||
pt.y() = bb.min().y() + j*d.y();
|
||||
for( int i = 0 ; i < mc.size_x() ; i++ )
|
||||
{
|
||||
pt.x() = bb.min().x() + i*d.x();
|
||||
points[pointI++] = pt;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
Info<< "Generated " << points.size() << " sampling points in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
|
||||
// Compute data
|
||||
const searchableSurfaces& geometry = geometryToConformTo.geometry();
|
||||
const labelList& surfaces = geometryToConformTo.surfaces();
|
||||
|
||||
scalarField signedDist;
|
||||
labelList nearestSurfaces;
|
||||
searchableSurfacesQueries::signedDistance
|
||||
(
|
||||
geometry,
|
||||
surfaces,
|
||||
points,
|
||||
scalarField(points.size(), sqr(GREAT)),
|
||||
searchableSurface::OUTSIDE, // for non-closed surfaces treat as
|
||||
// outside
|
||||
nearestSurfaces,
|
||||
signedDist
|
||||
);
|
||||
|
||||
// Fill elements
|
||||
pointI = 0;
|
||||
for( int k = 0 ; k < mc.size_z() ; k++ )
|
||||
{
|
||||
for( int j = 0 ; j < mc.size_y() ; j++ )
|
||||
{
|
||||
for( int i = 0 ; i < mc.size_x() ; i++ )
|
||||
{
|
||||
mc.set_data(float(signedDist[pointI++]), i, j, k);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
Info<< "Determined signed distance in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
|
||||
mc.run() ;
|
||||
|
||||
Info<< "Constructed iso surface in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
|
||||
mc.clean_temps() ;
|
||||
|
||||
|
||||
|
||||
// Write output file
|
||||
if (mc.ntrigs() > 0)
|
||||
{
|
||||
Triangle* triangles = mc.triangles();
|
||||
List<labelledTri> tris(mc.ntrigs());
|
||||
forAll(tris, triI)
|
||||
{
|
||||
tris[triI] = labelledTri
|
||||
(
|
||||
triangles[triI].v1,
|
||||
triangles[triI].v2,
|
||||
triangles[triI].v3,
|
||||
0 // region
|
||||
);
|
||||
}
|
||||
|
||||
|
||||
Vertex* vertices = mc.vertices();
|
||||
pointField points(mc.nverts());
|
||||
forAll(points, pointI)
|
||||
{
|
||||
Vertex& v = vertices[pointI];
|
||||
points[pointI] = point
|
||||
(
|
||||
bb.min().x() + v.x*span.x()/mc.size_x(),
|
||||
bb.min().y() + v.y*span.y()/mc.size_y(),
|
||||
bb.min().z() + v.z*span.z()/mc.size_z()
|
||||
);
|
||||
}
|
||||
|
||||
|
||||
// Find correspondence to original surfaces
|
||||
labelList regionOffsets(surfaces.size());
|
||||
label nRegions = 0;
|
||||
forAll(surfaces, i)
|
||||
{
|
||||
const wordList& regions = geometry[surfaces[i]].regions();
|
||||
regionOffsets[i] = nRegions;
|
||||
nRegions += regions.size();
|
||||
}
|
||||
|
||||
|
||||
geometricSurfacePatchList patches(nRegions);
|
||||
nRegions = 0;
|
||||
forAll(surfaces, i)
|
||||
{
|
||||
const wordList& regions = geometry[surfaces[i]].regions();
|
||||
|
||||
forAll(regions, regionI)
|
||||
{
|
||||
patches[nRegions] = geometricSurfacePatch
|
||||
(
|
||||
"patch",
|
||||
geometry[surfaces[i]].name() + "_" + regions[regionI],
|
||||
nRegions
|
||||
);
|
||||
nRegions++;
|
||||
}
|
||||
}
|
||||
|
||||
triSurface s(tris, patches, points, true);
|
||||
|
||||
Info<< "Extracted triSurface in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
|
||||
// Find out region on local surface of nearest point
|
||||
{
|
||||
List<pointIndexHit> hitInfo;
|
||||
labelList hitSurfaces;
|
||||
geometryToConformTo.findSurfaceNearest
|
||||
(
|
||||
s.faceCentres(),
|
||||
scalarField(s.size(), sqr(GREAT)),
|
||||
hitInfo,
|
||||
hitSurfaces
|
||||
);
|
||||
|
||||
// Get region
|
||||
DynamicList<pointIndexHit> surfInfo(hitSurfaces.size());
|
||||
DynamicList<label> surfIndices(hitSurfaces.size());
|
||||
|
||||
forAll(surfaces, surfI)
|
||||
{
|
||||
// Extract info on this surface
|
||||
surfInfo.clear();
|
||||
surfIndices.clear();
|
||||
forAll(hitSurfaces, triI)
|
||||
{
|
||||
if (hitSurfaces[triI] == surfI)
|
||||
{
|
||||
surfInfo.append(hitInfo[triI]);
|
||||
surfIndices.append(triI);
|
||||
}
|
||||
}
|
||||
|
||||
// Calculate sideness of these surface points
|
||||
labelList region;
|
||||
geometry[surfaces[surfI]].getRegion(surfInfo, region);
|
||||
|
||||
forAll(region, i)
|
||||
{
|
||||
label triI = surfIndices[i];
|
||||
s[triI].region() = regionOffsets[surfI]+region[i];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
Info<< "Re-patched surface in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
|
||||
triSurfaceMesh smesh
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
exportName,
|
||||
runTime.constant(), // instance
|
||||
"triSurface",
|
||||
runTime, // registry
|
||||
IOobject::NO_READ,
|
||||
IOobject::AUTO_WRITE,
|
||||
false
|
||||
),
|
||||
s
|
||||
);
|
||||
|
||||
Info<< "writing surfMesh:\n "
|
||||
<< smesh.searchableSurface::objectPath() << nl << endl;
|
||||
smesh.searchableSurface::write();
|
||||
|
||||
Info<< "Written surface in = "
|
||||
<< timer.cpuTimeIncrement() << " s." << nl << endl;
|
||||
}
|
||||
|
||||
mc.clean_all() ;
|
||||
|
||||
|
||||
Info<< "End\n" << endl;
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
||||
|
||||
// ************************************************************************* //
|
||||
Binary file not shown.
@ -154,6 +154,13 @@ castellatedMeshControls
|
||||
//faceZone sphere;
|
||||
//cellZone sphere;
|
||||
//cellZoneInside inside; //outside/insidePoint
|
||||
//- Optional specification of what to do with faceZone faces:
|
||||
// internal : keep them as internal faces (default)
|
||||
// baffle : create baffles from them. This gives more
|
||||
// freedom in mesh motion
|
||||
// boundary : create free-standing boundary faces (baffles
|
||||
// but without the shared points)
|
||||
//faceType internal;
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
@ -641,9 +641,12 @@ int main(int argc, char *argv[])
|
||||
Info<< "Mesh size: " << mesh.globalData().nTotalCells() << nl
|
||||
<< "Before renumbering :" << nl
|
||||
<< " band : " << band << nl
|
||||
<< " profile : " << profile << nl
|
||||
<< " rms frontwidth : " << rmsFrontwidth << nl
|
||||
<< endl;
|
||||
<< " profile : " << profile << nl;
|
||||
if (doFrontWidth)
|
||||
{
|
||||
Info<< " rms frontwidth : " << rmsFrontwidth << nl;
|
||||
}
|
||||
Info<< endl;
|
||||
|
||||
bool sortCoupledFaceCells = false;
|
||||
bool writeMaps = false;
|
||||
|
||||
@ -1,3 +1,4 @@
|
||||
loadOrCreateMesh.C
|
||||
redistributePar.C
|
||||
|
||||
EXE = $(FOAM_APPBIN)/redistributePar
|
||||
|
||||
@ -0,0 +1,323 @@
|
||||
/*---------------------------------------------------------------------------*\
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
This file is part of OpenFOAM.
|
||||
|
||||
OpenFOAM is free software: you can redistribute it and/or modify it
|
||||
under the terms of the GNU General Public License as published by
|
||||
the Free Software Foundation, either version 3 of the License, or
|
||||
(at your option) any later version.
|
||||
|
||||
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
|
||||
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
|
||||
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
|
||||
for more details.
|
||||
|
||||
You should have received a copy of the GNU General Public License
|
||||
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
\*---------------------------------------------------------------------------*/
|
||||
|
||||
#include "loadOrCreateMesh.H"
|
||||
#include "processorPolyPatch.H"
|
||||
#include "Time.H"
|
||||
|
||||
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
|
||||
|
||||
// Read mesh if available. Otherwise create empty mesh with same non-proc
|
||||
// patches as proc0 mesh. Requires all processors to have all patches
|
||||
// (and in same order).
|
||||
Foam::autoPtr<Foam::fvMesh> Foam::loadOrCreateMesh
|
||||
(
|
||||
const IOobject& io
|
||||
)
|
||||
{
|
||||
fileName meshSubDir;
|
||||
|
||||
if (io.name() == polyMesh::defaultRegion)
|
||||
{
|
||||
meshSubDir = polyMesh::meshSubDir;
|
||||
}
|
||||
else
|
||||
{
|
||||
meshSubDir = io.name()/polyMesh::meshSubDir;
|
||||
}
|
||||
|
||||
// Check who has a mesh
|
||||
const bool haveMesh = isDir(io.time().path()/io.instance()/meshSubDir);
|
||||
|
||||
Pout<< "meshpath:" << io.time().path()/io.instance()/meshSubDir << endl;
|
||||
Pout<< "haveMesh:" << haveMesh << endl;
|
||||
|
||||
if (!haveMesh)
|
||||
{
|
||||
// Create dummy mesh. Only used on procs that don't have mesh.
|
||||
IOobject noReadIO(io);
|
||||
noReadIO.readOpt() = IOobject::NO_READ;
|
||||
fvMesh dummyMesh
|
||||
(
|
||||
noReadIO,
|
||||
xferCopy(pointField()),
|
||||
xferCopy(faceList()),
|
||||
xferCopy(labelList()),
|
||||
xferCopy(labelList()),
|
||||
false
|
||||
);
|
||||
// Add some dummy zones so upon reading it does not read them
|
||||
// from the undecomposed case. Should be done as extra argument to
|
||||
// regIOobject::readStream?
|
||||
List<pointZone*> pz
|
||||
(
|
||||
1,
|
||||
new pointZone
|
||||
(
|
||||
"dummyPointZone",
|
||||
labelList(0),
|
||||
0,
|
||||
dummyMesh.pointZones()
|
||||
)
|
||||
);
|
||||
List<faceZone*> fz
|
||||
(
|
||||
1,
|
||||
new faceZone
|
||||
(
|
||||
"dummyFaceZone",
|
||||
labelList(0),
|
||||
boolList(0),
|
||||
0,
|
||||
dummyMesh.faceZones()
|
||||
)
|
||||
);
|
||||
List<cellZone*> cz
|
||||
(
|
||||
1,
|
||||
new cellZone
|
||||
(
|
||||
"dummyCellZone",
|
||||
labelList(0),
|
||||
0,
|
||||
dummyMesh.cellZones()
|
||||
)
|
||||
);
|
||||
dummyMesh.addZones(pz, fz, cz);
|
||||
//Pout<< "Writing dummy mesh to " << dummyMesh.polyMesh::objectPath()
|
||||
// << endl;
|
||||
dummyMesh.write();
|
||||
}
|
||||
|
||||
//Pout<< "Reading mesh from " << io.objectPath() << endl;
|
||||
autoPtr<fvMesh> meshPtr(new fvMesh(io));
|
||||
fvMesh& mesh = meshPtr();
|
||||
|
||||
|
||||
// Sync patches
|
||||
// ~~~~~~~~~~~~
|
||||
|
||||
if (Pstream::master())
|
||||
{
|
||||
// Send patches
|
||||
for
|
||||
(
|
||||
int slave=Pstream::firstSlave();
|
||||
slave<=Pstream::lastSlave();
|
||||
slave++
|
||||
)
|
||||
{
|
||||
OPstream toSlave(Pstream::scheduled, slave);
|
||||
toSlave << mesh.boundaryMesh();
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
// Receive patches
|
||||
IPstream fromMaster(Pstream::scheduled, Pstream::masterNo());
|
||||
PtrList<entry> patchEntries(fromMaster);
|
||||
|
||||
if (haveMesh)
|
||||
{
|
||||
// Check master names against mine
|
||||
|
||||
const polyBoundaryMesh& patches = mesh.boundaryMesh();
|
||||
|
||||
forAll(patchEntries, patchI)
|
||||
{
|
||||
const entry& e = patchEntries[patchI];
|
||||
const word type(e.dict().lookup("type"));
|
||||
const word& name = e.keyword();
|
||||
|
||||
if (type == processorPolyPatch::typeName)
|
||||
{
|
||||
break;
|
||||
}
|
||||
|
||||
if (patchI >= patches.size())
|
||||
{
|
||||
FatalErrorIn
|
||||
(
|
||||
"createMesh(const Time&, const fileName&, const bool)"
|
||||
) << "Non-processor patches not synchronised."
|
||||
<< endl
|
||||
<< "Processor " << Pstream::myProcNo()
|
||||
<< " has only " << patches.size()
|
||||
<< " patches, master has "
|
||||
<< patchI
|
||||
<< exit(FatalError);
|
||||
}
|
||||
|
||||
if
|
||||
(
|
||||
type != patches[patchI].type()
|
||||
|| name != patches[patchI].name()
|
||||
)
|
||||
{
|
||||
FatalErrorIn
|
||||
(
|
||||
"createMesh(const Time&, const fileName&, const bool)"
|
||||
) << "Non-processor patches not synchronised."
|
||||
<< endl
|
||||
<< "Master patch " << patchI
|
||||
<< " name:" << type
|
||||
<< " type:" << type << endl
|
||||
<< "Processor " << Pstream::myProcNo()
|
||||
<< " patch " << patchI
|
||||
<< " has name:" << patches[patchI].name()
|
||||
<< " type:" << patches[patchI].type()
|
||||
<< exit(FatalError);
|
||||
}
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
// Add patch
|
||||
List<polyPatch*> patches(patchEntries.size());
|
||||
label nPatches = 0;
|
||||
|
||||
forAll(patchEntries, patchI)
|
||||
{
|
||||
const entry& e = patchEntries[patchI];
|
||||
const word type(e.dict().lookup("type"));
|
||||
const word& name = e.keyword();
|
||||
|
||||
if (type == processorPolyPatch::typeName)
|
||||
{
|
||||
break;
|
||||
}
|
||||
|
||||
//Pout<< "Adding patch:" << nPatches
|
||||
// << " name:" << name << " type:" << type << endl;
|
||||
|
||||
dictionary patchDict(e.dict());
|
||||
patchDict.remove("nFaces");
|
||||
patchDict.add("nFaces", 0);
|
||||
patchDict.remove("startFace");
|
||||
patchDict.add("startFace", 0);
|
||||
|
||||
patches[patchI] = polyPatch::New
|
||||
(
|
||||
name,
|
||||
patchDict,
|
||||
nPatches++,
|
||||
mesh.boundaryMesh()
|
||||
).ptr();
|
||||
}
|
||||
patches.setSize(nPatches);
|
||||
mesh.addFvPatches(patches, false); // no parallel comms
|
||||
|
||||
//// Write empty mesh now we have correct patches
|
||||
//meshPtr().write();
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
// Determine zones
|
||||
// ~~~~~~~~~~~~~~~
|
||||
|
||||
wordList pointZoneNames(mesh.pointZones().names());
|
||||
Pstream::scatter(pointZoneNames);
|
||||
wordList faceZoneNames(mesh.faceZones().names());
|
||||
Pstream::scatter(faceZoneNames);
|
||||
wordList cellZoneNames(mesh.cellZones().names());
|
||||
Pstream::scatter(cellZoneNames);
|
||||
|
||||
if (!haveMesh)
|
||||
{
|
||||
// Add the zones. Make sure to remove the old dummy ones first
|
||||
mesh.pointZones().clear();
|
||||
mesh.faceZones().clear();
|
||||
mesh.cellZones().clear();
|
||||
|
||||
List<pointZone*> pz(pointZoneNames.size());
|
||||
forAll(pointZoneNames, i)
|
||||
{
|
||||
pz[i] = new pointZone
|
||||
(
|
||||
pointZoneNames[i],
|
||||
labelList(0),
|
||||
i,
|
||||
mesh.pointZones()
|
||||
);
|
||||
}
|
||||
List<faceZone*> fz(faceZoneNames.size());
|
||||
forAll(faceZoneNames, i)
|
||||
{
|
||||
fz[i] = new faceZone
|
||||
(
|
||||
faceZoneNames[i],
|
||||
labelList(0),
|
||||
boolList(0),
|
||||
i,
|
||||
mesh.faceZones()
|
||||
);
|
||||
}
|
||||
List<cellZone*> cz(cellZoneNames.size());
|
||||
forAll(cellZoneNames, i)
|
||||
{
|
||||
cz[i] = new cellZone
|
||||
(
|
||||
cellZoneNames[i],
|
||||
labelList(0),
|
||||
i,
|
||||
mesh.cellZones()
|
||||
);
|
||||
}
|
||||
mesh.addZones(pz, fz, cz);
|
||||
}
|
||||
|
||||
|
||||
if (!haveMesh)
|
||||
{
|
||||
// We created a dummy mesh file above. Delete it.
|
||||
//Pout<< "Removing dummy mesh " << io.objectPath() << endl;
|
||||
rmDir(io.objectPath());
|
||||
}
|
||||
|
||||
// Force recreation of globalMeshData.
|
||||
mesh.clearOut();
|
||||
mesh.globalData();
|
||||
|
||||
|
||||
// Do some checks.
|
||||
|
||||
// Check if the boundary definition is unique
|
||||
mesh.boundaryMesh().checkDefinition(true);
|
||||
// Check if the boundary processor patches are correct
|
||||
mesh.boundaryMesh().checkParallelSync(true);
|
||||
// Check names of zones are equal
|
||||
mesh.cellZones().checkDefinition(true);
|
||||
mesh.cellZones().checkParallelSync(true);
|
||||
mesh.faceZones().checkDefinition(true);
|
||||
mesh.faceZones().checkParallelSync(true);
|
||||
mesh.pointZones().checkDefinition(true);
|
||||
mesh.pointZones().checkParallelSync(true);
|
||||
|
||||
return meshPtr;
|
||||
}
|
||||
|
||||
|
||||
// ************************************************************************* //
|
||||
@ -0,0 +1,58 @@
|
||||
/*---------------------------------------------------------------------------*\
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
This file is part of OpenFOAM.
|
||||
|
||||
OpenFOAM is free software: you can redistribute it and/or modify it
|
||||
under the terms of the GNU General Public License as published by
|
||||
the Free Software Foundation, either version 3 of the License, or
|
||||
(at your option) any later version.
|
||||
|
||||
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
|
||||
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
|
||||
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
|
||||
for more details.
|
||||
|
||||
You should have received a copy of the GNU General Public License
|
||||
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
InNamespace
|
||||
Foam
|
||||
|
||||
Description
|
||||
Load or create (0 size) a mesh. Used in distributing meshes to a
|
||||
larger number of processors
|
||||
|
||||
SourceFiles
|
||||
loadOrCreateMesh.C
|
||||
|
||||
\*---------------------------------------------------------------------------*/
|
||||
|
||||
#ifndef loadOrCreateMesh_H
|
||||
#define loadOrCreateMesh_H
|
||||
|
||||
#include "fvMesh.H"
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
namespace Foam
|
||||
{
|
||||
|
||||
//- Load (if it exists) or create zero cell mesh given an IOobject:
|
||||
// name : regionName
|
||||
// instance : exact directory where to find mesh (i.e. does not
|
||||
// do a findInstance
|
||||
autoPtr<fvMesh> loadOrCreateMesh(const IOobject& io);
|
||||
|
||||
}
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
#endif
|
||||
|
||||
// ************************************************************************* //
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -54,6 +54,7 @@ Description
|
||||
#include "mapDistributePolyMesh.H"
|
||||
#include "IOobjectList.H"
|
||||
#include "globalIndex.H"
|
||||
#include "loadOrCreateMesh.H"
|
||||
|
||||
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
||||
|
||||
@ -62,292 +63,292 @@ Description
|
||||
static const scalar defaultMergeTol = 1E-6;
|
||||
|
||||
|
||||
// Read mesh if available. Otherwise create empty mesh with same non-proc
|
||||
// patches as proc0 mesh. Requires all processors to have all patches
|
||||
// (and in same order).
|
||||
autoPtr<fvMesh> createMesh
|
||||
(
|
||||
const Time& runTime,
|
||||
const word& regionName,
|
||||
const fileName& instDir,
|
||||
const bool haveMesh
|
||||
)
|
||||
{
|
||||
//Pout<< "Create mesh for time = "
|
||||
// << runTime.timeName() << nl << endl;
|
||||
|
||||
IOobject io
|
||||
(
|
||||
regionName,
|
||||
instDir,
|
||||
runTime,
|
||||
IOobject::MUST_READ
|
||||
);
|
||||
|
||||
if (!haveMesh)
|
||||
{
|
||||
// Create dummy mesh. Only used on procs that don't have mesh.
|
||||
IOobject noReadIO(io);
|
||||
noReadIO.readOpt() = IOobject::NO_READ;
|
||||
fvMesh dummyMesh
|
||||
(
|
||||
noReadIO,
|
||||
xferCopy(pointField()),
|
||||
xferCopy(faceList()),
|
||||
xferCopy(labelList()),
|
||||
xferCopy(labelList()),
|
||||
false
|
||||
);
|
||||
// Add some dummy zones so upon reading it does not read them
|
||||
// from the undecomposed case. Should be done as extra argument to
|
||||
// regIOobject::readStream?
|
||||
List<pointZone*> pz
|
||||
(
|
||||
1,
|
||||
new pointZone
|
||||
(
|
||||
"dummyPointZone",
|
||||
labelList(0),
|
||||
0,
|
||||
dummyMesh.pointZones()
|
||||
)
|
||||
);
|
||||
List<faceZone*> fz
|
||||
(
|
||||
1,
|
||||
new faceZone
|
||||
(
|
||||
"dummyFaceZone",
|
||||
labelList(0),
|
||||
boolList(0),
|
||||
0,
|
||||
dummyMesh.faceZones()
|
||||
)
|
||||
);
|
||||
List<cellZone*> cz
|
||||
(
|
||||
1,
|
||||
new cellZone
|
||||
(
|
||||
"dummyCellZone",
|
||||
labelList(0),
|
||||
0,
|
||||
dummyMesh.cellZones()
|
||||
)
|
||||
);
|
||||
dummyMesh.addZones(pz, fz, cz);
|
||||
//Pout<< "Writing dummy mesh to " << dummyMesh.polyMesh::objectPath()
|
||||
// << endl;
|
||||
dummyMesh.write();
|
||||
}
|
||||
|
||||
//Pout<< "Reading mesh from " << io.objectPath() << endl;
|
||||
autoPtr<fvMesh> meshPtr(new fvMesh(io));
|
||||
fvMesh& mesh = meshPtr();
|
||||
|
||||
|
||||
// Sync patches
|
||||
// ~~~~~~~~~~~~
|
||||
|
||||
if (Pstream::master())
|
||||
{
|
||||
// Send patches
|
||||
for
|
||||
(
|
||||
int slave=Pstream::firstSlave();
|
||||
slave<=Pstream::lastSlave();
|
||||
slave++
|
||||
)
|
||||
{
|
||||
OPstream toSlave(Pstream::scheduled, slave);
|
||||
toSlave << mesh.boundaryMesh();
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
// Receive patches
|
||||
IPstream fromMaster(Pstream::scheduled, Pstream::masterNo());
|
||||
PtrList<entry> patchEntries(fromMaster);
|
||||
|
||||
if (haveMesh)
|
||||
{
|
||||
// Check master names against mine
|
||||
|
||||
const polyBoundaryMesh& patches = mesh.boundaryMesh();
|
||||
|
||||
forAll(patchEntries, patchI)
|
||||
{
|
||||
const entry& e = patchEntries[patchI];
|
||||
const word type(e.dict().lookup("type"));
|
||||
const word& name = e.keyword();
|
||||
|
||||
if (type == processorPolyPatch::typeName)
|
||||
{
|
||||
break;
|
||||
}
|
||||
|
||||
if (patchI >= patches.size())
|
||||
{
|
||||
FatalErrorIn
|
||||
(
|
||||
"createMesh(const Time&, const fileName&, const bool)"
|
||||
) << "Non-processor patches not synchronised."
|
||||
<< endl
|
||||
<< "Processor " << Pstream::myProcNo()
|
||||
<< " has only " << patches.size()
|
||||
<< " patches, master has "
|
||||
<< patchI
|
||||
<< exit(FatalError);
|
||||
}
|
||||
|
||||
if
|
||||
(
|
||||
type != patches[patchI].type()
|
||||
|| name != patches[patchI].name()
|
||||
)
|
||||
{
|
||||
FatalErrorIn
|
||||
(
|
||||
"createMesh(const Time&, const fileName&, const bool)"
|
||||
) << "Non-processor patches not synchronised."
|
||||
<< endl
|
||||
<< "Master patch " << patchI
|
||||
<< " name:" << type
|
||||
<< " type:" << type << endl
|
||||
<< "Processor " << Pstream::myProcNo()
|
||||
<< " patch " << patchI
|
||||
<< " has name:" << patches[patchI].name()
|
||||
<< " type:" << patches[patchI].type()
|
||||
<< exit(FatalError);
|
||||
}
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
// Add patch
|
||||
List<polyPatch*> patches(patchEntries.size());
|
||||
label nPatches = 0;
|
||||
|
||||
forAll(patchEntries, patchI)
|
||||
{
|
||||
const entry& e = patchEntries[patchI];
|
||||
const word type(e.dict().lookup("type"));
|
||||
const word& name = e.keyword();
|
||||
|
||||
if (type == processorPolyPatch::typeName)
|
||||
{
|
||||
break;
|
||||
}
|
||||
|
||||
//Pout<< "Adding patch:" << nPatches
|
||||
// << " name:" << name << " type:" << type << endl;
|
||||
|
||||
dictionary patchDict(e.dict());
|
||||
patchDict.remove("nFaces");
|
||||
patchDict.add("nFaces", 0);
|
||||
patchDict.remove("startFace");
|
||||
patchDict.add("startFace", 0);
|
||||
|
||||
patches[patchI] = polyPatch::New
|
||||
(
|
||||
name,
|
||||
patchDict,
|
||||
nPatches++,
|
||||
mesh.boundaryMesh()
|
||||
).ptr();
|
||||
}
|
||||
patches.setSize(nPatches);
|
||||
mesh.addFvPatches(patches, false); // no parallel comms
|
||||
|
||||
//// Write empty mesh now we have correct patches
|
||||
//meshPtr().write();
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
// Determine zones
|
||||
// ~~~~~~~~~~~~~~~
|
||||
|
||||
wordList pointZoneNames(mesh.pointZones().names());
|
||||
Pstream::scatter(pointZoneNames);
|
||||
wordList faceZoneNames(mesh.faceZones().names());
|
||||
Pstream::scatter(faceZoneNames);
|
||||
wordList cellZoneNames(mesh.cellZones().names());
|
||||
Pstream::scatter(cellZoneNames);
|
||||
|
||||
if (!haveMesh)
|
||||
{
|
||||
// Add the zones. Make sure to remove the old dummy ones first
|
||||
mesh.pointZones().clear();
|
||||
mesh.faceZones().clear();
|
||||
mesh.cellZones().clear();
|
||||
|
||||
List<pointZone*> pz(pointZoneNames.size());
|
||||
forAll(pointZoneNames, i)
|
||||
{
|
||||
pz[i] = new pointZone
|
||||
(
|
||||
pointZoneNames[i],
|
||||
labelList(0),
|
||||
i,
|
||||
mesh.pointZones()
|
||||
);
|
||||
}
|
||||
List<faceZone*> fz(faceZoneNames.size());
|
||||
forAll(faceZoneNames, i)
|
||||
{
|
||||
fz[i] = new faceZone
|
||||
(
|
||||
faceZoneNames[i],
|
||||
labelList(0),
|
||||
boolList(0),
|
||||
i,
|
||||
mesh.faceZones()
|
||||
);
|
||||
}
|
||||
List<cellZone*> cz(cellZoneNames.size());
|
||||
forAll(cellZoneNames, i)
|
||||
{
|
||||
cz[i] = new cellZone
|
||||
(
|
||||
cellZoneNames[i],
|
||||
labelList(0),
|
||||
i,
|
||||
mesh.cellZones()
|
||||
);
|
||||
}
|
||||
mesh.addZones(pz, fz, cz);
|
||||
}
|
||||
|
||||
|
||||
if (!haveMesh)
|
||||
{
|
||||
// We created a dummy mesh file above. Delete it.
|
||||
//Pout<< "Removing dummy mesh " << io.objectPath() << endl;
|
||||
rmDir(io.objectPath());
|
||||
}
|
||||
|
||||
// Force recreation of globalMeshData.
|
||||
mesh.clearOut();
|
||||
mesh.globalData();
|
||||
|
||||
|
||||
// Do some checks.
|
||||
|
||||
// Check if the boundary definition is unique
|
||||
mesh.boundaryMesh().checkDefinition(true);
|
||||
// Check if the boundary processor patches are correct
|
||||
mesh.boundaryMesh().checkParallelSync(true);
|
||||
// Check names of zones are equal
|
||||
mesh.cellZones().checkDefinition(true);
|
||||
mesh.cellZones().checkParallelSync(true);
|
||||
mesh.faceZones().checkDefinition(true);
|
||||
mesh.faceZones().checkParallelSync(true);
|
||||
mesh.pointZones().checkDefinition(true);
|
||||
mesh.pointZones().checkParallelSync(true);
|
||||
|
||||
return meshPtr;
|
||||
}
|
||||
//// Read mesh if available. Otherwise create empty mesh with same non-proc
|
||||
//// patches as proc0 mesh. Requires all processors to have all patches
|
||||
//// (and in same order).
|
||||
//autoPtr<fvMesh> createMesh
|
||||
//(
|
||||
// const Time& runTime,
|
||||
// const word& regionName,
|
||||
// const fileName& instDir,
|
||||
// const bool haveMesh
|
||||
//)
|
||||
//{
|
||||
// //Pout<< "Create mesh for time = "
|
||||
// // << runTime.timeName() << nl << endl;
|
||||
//
|
||||
// IOobject io
|
||||
// (
|
||||
// regionName,
|
||||
// instDir,
|
||||
// runTime,
|
||||
// IOobject::MUST_READ
|
||||
// );
|
||||
//
|
||||
// if (!haveMesh)
|
||||
// {
|
||||
// // Create dummy mesh. Only used on procs that don't have mesh.
|
||||
// IOobject noReadIO(io);
|
||||
// noReadIO.readOpt() = IOobject::NO_READ;
|
||||
// fvMesh dummyMesh
|
||||
// (
|
||||
// noReadIO,
|
||||
// xferCopy(pointField()),
|
||||
// xferCopy(faceList()),
|
||||
// xferCopy(labelList()),
|
||||
// xferCopy(labelList()),
|
||||
// false
|
||||
// );
|
||||
// // Add some dummy zones so upon reading it does not read them
|
||||
// // from the undecomposed case. Should be done as extra argument to
|
||||
// // regIOobject::readStream?
|
||||
// List<pointZone*> pz
|
||||
// (
|
||||
// 1,
|
||||
// new pointZone
|
||||
// (
|
||||
// "dummyPointZone",
|
||||
// labelList(0),
|
||||
// 0,
|
||||
// dummyMesh.pointZones()
|
||||
// )
|
||||
// );
|
||||
// List<faceZone*> fz
|
||||
// (
|
||||
// 1,
|
||||
// new faceZone
|
||||
// (
|
||||
// "dummyFaceZone",
|
||||
// labelList(0),
|
||||
// boolList(0),
|
||||
// 0,
|
||||
// dummyMesh.faceZones()
|
||||
// )
|
||||
// );
|
||||
// List<cellZone*> cz
|
||||
// (
|
||||
// 1,
|
||||
// new cellZone
|
||||
// (
|
||||
// "dummyCellZone",
|
||||
// labelList(0),
|
||||
// 0,
|
||||
// dummyMesh.cellZones()
|
||||
// )
|
||||
// );
|
||||
// dummyMesh.addZones(pz, fz, cz);
|
||||
// //Pout<< "Writing dummy mesh to " << dummyMesh.polyMesh::objectPath()
|
||||
// // << endl;
|
||||
// dummyMesh.write();
|
||||
// }
|
||||
//
|
||||
// //Pout<< "Reading mesh from " << io.objectPath() << endl;
|
||||
// autoPtr<fvMesh> meshPtr(new fvMesh(io));
|
||||
// fvMesh& mesh = meshPtr();
|
||||
//
|
||||
//
|
||||
// // Sync patches
|
||||
// // ~~~~~~~~~~~~
|
||||
//
|
||||
// if (Pstream::master())
|
||||
// {
|
||||
// // Send patches
|
||||
// for
|
||||
// (
|
||||
// int slave=Pstream::firstSlave();
|
||||
// slave<=Pstream::lastSlave();
|
||||
// slave++
|
||||
// )
|
||||
// {
|
||||
// OPstream toSlave(Pstream::scheduled, slave);
|
||||
// toSlave << mesh.boundaryMesh();
|
||||
// }
|
||||
// }
|
||||
// else
|
||||
// {
|
||||
// // Receive patches
|
||||
// IPstream fromMaster(Pstream::scheduled, Pstream::masterNo());
|
||||
// PtrList<entry> patchEntries(fromMaster);
|
||||
//
|
||||
// if (haveMesh)
|
||||
// {
|
||||
// // Check master names against mine
|
||||
//
|
||||
// const polyBoundaryMesh& patches = mesh.boundaryMesh();
|
||||
//
|
||||
// forAll(patchEntries, patchI)
|
||||
// {
|
||||
// const entry& e = patchEntries[patchI];
|
||||
// const word type(e.dict().lookup("type"));
|
||||
// const word& name = e.keyword();
|
||||
//
|
||||
// if (type == processorPolyPatch::typeName)
|
||||
// {
|
||||
// break;
|
||||
// }
|
||||
//
|
||||
// if (patchI >= patches.size())
|
||||
// {
|
||||
// FatalErrorIn
|
||||
// (
|
||||
// "createMesh(const Time&, const fileName&, const bool)"
|
||||
// ) << "Non-processor patches not synchronised."
|
||||
// << endl
|
||||
// << "Processor " << Pstream::myProcNo()
|
||||
// << " has only " << patches.size()
|
||||
// << " patches, master has "
|
||||
// << patchI
|
||||
// << exit(FatalError);
|
||||
// }
|
||||
//
|
||||
// if
|
||||
// (
|
||||
// type != patches[patchI].type()
|
||||
// || name != patches[patchI].name()
|
||||
// )
|
||||
// {
|
||||
// FatalErrorIn
|
||||
// (
|
||||
// "createMesh(const Time&, const fileName&, const bool)"
|
||||
// ) << "Non-processor patches not synchronised."
|
||||
// << endl
|
||||
// << "Master patch " << patchI
|
||||
// << " name:" << type
|
||||
// << " type:" << type << endl
|
||||
// << "Processor " << Pstream::myProcNo()
|
||||
// << " patch " << patchI
|
||||
// << " has name:" << patches[patchI].name()
|
||||
// << " type:" << patches[patchI].type()
|
||||
// << exit(FatalError);
|
||||
// }
|
||||
// }
|
||||
// }
|
||||
// else
|
||||
// {
|
||||
// // Add patch
|
||||
// List<polyPatch*> patches(patchEntries.size());
|
||||
// label nPatches = 0;
|
||||
//
|
||||
// forAll(patchEntries, patchI)
|
||||
// {
|
||||
// const entry& e = patchEntries[patchI];
|
||||
// const word type(e.dict().lookup("type"));
|
||||
// const word& name = e.keyword();
|
||||
//
|
||||
// if (type == processorPolyPatch::typeName)
|
||||
// {
|
||||
// break;
|
||||
// }
|
||||
//
|
||||
// //Pout<< "Adding patch:" << nPatches
|
||||
// // << " name:" << name << " type:" << type << endl;
|
||||
//
|
||||
// dictionary patchDict(e.dict());
|
||||
// patchDict.remove("nFaces");
|
||||
// patchDict.add("nFaces", 0);
|
||||
// patchDict.remove("startFace");
|
||||
// patchDict.add("startFace", 0);
|
||||
//
|
||||
// patches[patchI] = polyPatch::New
|
||||
// (
|
||||
// name,
|
||||
// patchDict,
|
||||
// nPatches++,
|
||||
// mesh.boundaryMesh()
|
||||
// ).ptr();
|
||||
// }
|
||||
// patches.setSize(nPatches);
|
||||
// mesh.addFvPatches(patches, false); // no parallel comms
|
||||
//
|
||||
// //// Write empty mesh now we have correct patches
|
||||
// //meshPtr().write();
|
||||
// }
|
||||
// }
|
||||
//
|
||||
//
|
||||
// // Determine zones
|
||||
// // ~~~~~~~~~~~~~~~
|
||||
//
|
||||
// wordList pointZoneNames(mesh.pointZones().names());
|
||||
// Pstream::scatter(pointZoneNames);
|
||||
// wordList faceZoneNames(mesh.faceZones().names());
|
||||
// Pstream::scatter(faceZoneNames);
|
||||
// wordList cellZoneNames(mesh.cellZones().names());
|
||||
// Pstream::scatter(cellZoneNames);
|
||||
//
|
||||
// if (!haveMesh)
|
||||
// {
|
||||
// // Add the zones. Make sure to remove the old dummy ones first
|
||||
// mesh.pointZones().clear();
|
||||
// mesh.faceZones().clear();
|
||||
// mesh.cellZones().clear();
|
||||
//
|
||||
// List<pointZone*> pz(pointZoneNames.size());
|
||||
// forAll(pointZoneNames, i)
|
||||
// {
|
||||
// pz[i] = new pointZone
|
||||
// (
|
||||
// pointZoneNames[i],
|
||||
// labelList(0),
|
||||
// i,
|
||||
// mesh.pointZones()
|
||||
// );
|
||||
// }
|
||||
// List<faceZone*> fz(faceZoneNames.size());
|
||||
// forAll(faceZoneNames, i)
|
||||
// {
|
||||
// fz[i] = new faceZone
|
||||
// (
|
||||
// faceZoneNames[i],
|
||||
// labelList(0),
|
||||
// boolList(0),
|
||||
// i,
|
||||
// mesh.faceZones()
|
||||
// );
|
||||
// }
|
||||
// List<cellZone*> cz(cellZoneNames.size());
|
||||
// forAll(cellZoneNames, i)
|
||||
// {
|
||||
// cz[i] = new cellZone
|
||||
// (
|
||||
// cellZoneNames[i],
|
||||
// labelList(0),
|
||||
// i,
|
||||
// mesh.cellZones()
|
||||
// );
|
||||
// }
|
||||
// mesh.addZones(pz, fz, cz);
|
||||
// }
|
||||
//
|
||||
//
|
||||
// if (!haveMesh)
|
||||
// {
|
||||
// // We created a dummy mesh file above. Delete it.
|
||||
// //Pout<< "Removing dummy mesh " << io.objectPath() << endl;
|
||||
// rmDir(io.objectPath());
|
||||
// }
|
||||
//
|
||||
// // Force recreation of globalMeshData.
|
||||
// mesh.clearOut();
|
||||
// mesh.globalData();
|
||||
//
|
||||
//
|
||||
// // Do some checks.
|
||||
//
|
||||
// // Check if the boundary definition is unique
|
||||
// mesh.boundaryMesh().checkDefinition(true);
|
||||
// // Check if the boundary processor patches are correct
|
||||
// mesh.boundaryMesh().checkParallelSync(true);
|
||||
// // Check names of zones are equal
|
||||
// mesh.cellZones().checkDefinition(true);
|
||||
// mesh.cellZones().checkParallelSync(true);
|
||||
// mesh.faceZones().checkDefinition(true);
|
||||
// mesh.faceZones().checkParallelSync(true);
|
||||
// mesh.pointZones().checkDefinition(true);
|
||||
// mesh.pointZones().checkParallelSync(true);
|
||||
//
|
||||
// return meshPtr;
|
||||
//}
|
||||
|
||||
|
||||
// Get merging distance when matching face centres
|
||||
@ -786,13 +787,24 @@ int main(int argc, char *argv[])
|
||||
const bool allHaveMesh = (findIndex(haveMesh, false) == -1);
|
||||
|
||||
// Create mesh
|
||||
autoPtr<fvMesh> meshPtr = createMesh
|
||||
//autoPtr<fvMesh> meshPtr = createMesh
|
||||
//(
|
||||
// runTime,
|
||||
// regionName,
|
||||
// masterInstDir,
|
||||
// haveMesh[Pstream::myProcNo()]
|
||||
//);
|
||||
autoPtr<fvMesh> meshPtr = loadOrCreateMesh
|
||||
(
|
||||
runTime,
|
||||
regionName,
|
||||
masterInstDir,
|
||||
haveMesh[Pstream::myProcNo()]
|
||||
IOobject
|
||||
(
|
||||
regionName,
|
||||
masterInstDir,
|
||||
runTime,
|
||||
Foam::IOobject::MUST_READ
|
||||
)
|
||||
);
|
||||
|
||||
fvMesh& mesh = meshPtr();
|
||||
|
||||
// Print some statistics
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -216,7 +216,6 @@ bool intersectSurfaces
|
||||
|
||||
int main(int argc, char *argv[])
|
||||
{
|
||||
argList::validArgs.clear();
|
||||
argList::noParallel();
|
||||
argList::validArgs.append("action");
|
||||
argList::validArgs.append("surface file");
|
||||
|
||||
@ -525,12 +525,6 @@ int main(int argc, char *argv[])
|
||||
surfaceDict.lookupOrAddDefault<Switch>("writeVTK", "off");
|
||||
const Switch writeObj =
|
||||
surfaceDict.lookupOrAddDefault<Switch>("writeObj", "off");
|
||||
const Switch writeFeatureEdgeMesh =
|
||||
surfaceDict.lookupOrAddDefault<Switch>
|
||||
(
|
||||
"writeFeatureEdgeMesh",
|
||||
"off"
|
||||
);
|
||||
|
||||
const Switch curvature =
|
||||
surfaceDict.lookupOrAddDefault<Switch>("curvature", "off");
|
||||
@ -592,7 +586,7 @@ int main(int argc, char *argv[])
|
||||
if (extractionMethod == "extractFromFile")
|
||||
{
|
||||
const fileName featureEdgeFile =
|
||||
surfaceDict.subDict("extractFromFile").lookup
|
||||
surfaceDict.subDict("extractFromFileCoeffs").lookup
|
||||
(
|
||||
"featureEdgeFile"
|
||||
);
|
||||
@ -612,7 +606,7 @@ int main(int argc, char *argv[])
|
||||
const scalar includedAngle =
|
||||
readScalar
|
||||
(
|
||||
surfaceDict.subDict("extractFromSurface").lookup
|
||||
surfaceDict.subDict("extractFromSurfaceCoeffs").lookup
|
||||
(
|
||||
"includedAngle"
|
||||
)
|
||||
@ -747,10 +741,10 @@ int main(int argc, char *argv[])
|
||||
surfaceFeatures newSet(surf);
|
||||
newSet.setFromStatus(edgeStat);
|
||||
|
||||
if (writeObj)
|
||||
{
|
||||
newSet.writeObj("final");
|
||||
}
|
||||
//if (writeObj)
|
||||
//{
|
||||
// newSet.writeObj("final");
|
||||
//}
|
||||
|
||||
Info<< nl
|
||||
<< "Final feature set after trimming and subsetting:" << nl
|
||||
@ -773,37 +767,36 @@ int main(int argc, char *argv[])
|
||||
Info<< nl << "Writing extendedFeatureEdgeMesh to "
|
||||
<< feMesh.objectPath() << endl;
|
||||
|
||||
mkDir(feMesh.path());
|
||||
|
||||
if (writeObj)
|
||||
{
|
||||
feMesh.writeObj(surfFileName.lessExt().name());
|
||||
feMesh.writeObj(feMesh.path()/surfFileName.lessExt().name());
|
||||
}
|
||||
|
||||
feMesh.write();
|
||||
|
||||
// Write a featureEdgeMesh for backwards compatibility
|
||||
if (writeFeatureEdgeMesh)
|
||||
{
|
||||
featureEdgeMesh bfeMesh
|
||||
featureEdgeMesh bfeMesh
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
IOobject
|
||||
(
|
||||
surfFileName.lessExt().name() + ".eMesh", // name
|
||||
runTime.constant(), // instance
|
||||
"triSurface",
|
||||
runTime, // registry
|
||||
IOobject::NO_READ,
|
||||
IOobject::AUTO_WRITE,
|
||||
false
|
||||
),
|
||||
feMesh.points(),
|
||||
feMesh.edges()
|
||||
);
|
||||
surfFileName.lessExt().name() + ".eMesh", // name
|
||||
runTime.constant(), // instance
|
||||
"triSurface",
|
||||
runTime, // registry
|
||||
IOobject::NO_READ,
|
||||
IOobject::AUTO_WRITE,
|
||||
false
|
||||
),
|
||||
feMesh.points(),
|
||||
feMesh.edges()
|
||||
);
|
||||
|
||||
Info<< nl << "Writing featureEdgeMesh to "
|
||||
<< bfeMesh.objectPath() << endl;
|
||||
Info<< nl << "Writing featureEdgeMesh to "
|
||||
<< bfeMesh.objectPath() << endl;
|
||||
|
||||
bfeMesh.regIOobject::write();
|
||||
}
|
||||
bfeMesh.regIOobject::write();
|
||||
|
||||
triSurfaceMesh searchSurf
|
||||
(
|
||||
|
||||
@ -16,10 +16,34 @@ FoamFile
|
||||
|
||||
surface1.stl
|
||||
{
|
||||
// extractFromFile || extractFromSurface
|
||||
// How to obtain raw features (extractFromFile || extractFromSurface)
|
||||
extractionMethod extractFromSurface;
|
||||
|
||||
extractFromSurfaceCoeffs
|
||||
{
|
||||
// Mark edges whose adjacent surface normals are at an angle less
|
||||
// than includedAngle as features
|
||||
// - 0 : selects no edges
|
||||
// - 180: selects all edges
|
||||
includedAngle 120;
|
||||
}
|
||||
|
||||
// Write options
|
||||
|
||||
// Write .eMesh file (for snappyHexMesh)
|
||||
writeFeatureEdgeMesh yes;
|
||||
|
||||
// Write features to obj format for postprocessing
|
||||
writeObj yes;
|
||||
}
|
||||
|
||||
|
||||
surface2.nas
|
||||
{
|
||||
// How to obtain raw features (extractFromFile || extractFromSurface)
|
||||
extractionMethod extractFromFile;
|
||||
|
||||
extractFromFile
|
||||
extractFromFileCoeffs
|
||||
{
|
||||
// Load from an existing feature edge file
|
||||
featureEdgeFile "constant/triSurface/featureEdges.nas";
|
||||
@ -62,24 +86,13 @@ surface1.stl
|
||||
closeness no;
|
||||
|
||||
// Write options
|
||||
writeVTK no;
|
||||
writeObj yes;
|
||||
writeFeatureEdgeMesh no;
|
||||
}
|
||||
|
||||
// Write features to obj format for postprocessing
|
||||
writeObj yes;
|
||||
|
||||
surface2.nas
|
||||
{
|
||||
extractionMethod extractFromSurface;
|
||||
|
||||
extractFromSurface
|
||||
{
|
||||
// Mark edges whose adjacent surface normals are at an angle less
|
||||
// than includedAngle as features
|
||||
// - 0 : selects no edges
|
||||
// - 180: selects all edges
|
||||
includedAngle 120;
|
||||
}
|
||||
// Write surface proximity and curvature fields to vtk format
|
||||
// for postprocessing
|
||||
writeVTK no;
|
||||
}
|
||||
|
||||
|
||||
|
||||
@ -2,7 +2,7 @@
|
||||
========= |
|
||||
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
||||
\\ / O peration |
|
||||
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
|
||||
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
|
||||
\\/ M anipulation |
|
||||
-------------------------------------------------------------------------------
|
||||
License
|
||||
@ -312,7 +312,7 @@ int main(int argc, char *argv[])
|
||||
|
||||
triSurface allSurf(allFaces, allPatches, allPoints);
|
||||
|
||||
// Cleanup (which does point merge as well
|
||||
// Cleanup (which does point merge as well)
|
||||
allSurf.cleanup(false);
|
||||
|
||||
// Write surface mesh
|
||||
|
||||
Reference in New Issue
Block a user