Commit 466b95af authored by mattijs's avatar mattijs
Browse files

use extrapolated cell-centre for testing;moved added patch information into meshRefinement class

parent a2535471
......@@ -170,6 +170,13 @@ int main(int argc, char *argv[])
const dictionary& layerDict = meshDict.subDict("addLayersControls");
const scalar mergeDist = getMergeDistance
(
mesh,
readScalar(meshDict.lookup("mergeTolerance"))
);
// Debug
// ~~~~~
......@@ -192,8 +199,9 @@ int main(int argc, char *argv[])
IOobject
(
"abc", // dummy name
mesh.time().constant(), // directory
"triSurface", // instance
//mesh.time().constant(), // instance
mesh.time().findInstance("triSurface", word::null),// instance
"triSurface", // local
mesh.time(), // registry
IOobject::MUST_READ,
IOobject::NO_WRITE
......@@ -235,6 +243,33 @@ int main(int argc, char *argv[])
<< mesh.time().cpuTimeIncrement() << " s" << nl << endl;
// Refinement engine
// ~~~~~~~~~~~~~~~~~
Info<< nl
<< "Determining initial surface intersections" << nl
<< "-----------------------------------------" << nl
<< endl;
// Main refinement engine
meshRefinement meshRefiner
(
mesh,
mergeDist, // tolerance used in sorting coordinates
surfaces, // for surface intersection refinement
shells // for volume (inside/outside) refinement
);
Info<< "Calculated surface intersections in = "
<< mesh.time().cpuTimeIncrement() << " s" << nl << endl;
// Some stats
meshRefiner.printMeshInfo(debug, "Initial mesh");
meshRefiner.write
(
debug&meshRefinement::OBJINTERSECTIONS,
mesh.time().path()/mesh.time().timeName()
);
// Add all the surface regions as patches
......@@ -265,9 +300,8 @@ int main(int argc, char *argv[])
forAll(regNames, i)
{
label patchI = meshRefinement::addPatch
label patchI = meshRefiner.addMeshedPatch
(
mesh,
regNames[i],
wallPolyPatch::typeName
);
......@@ -308,45 +342,10 @@ int main(int argc, char *argv[])
<< exit(FatalError);
}
const scalar mergeDist = getMergeDistance
(
mesh,
readScalar(meshDict.lookup("mergeTolerance"))
);
// Mesh distribution engine (uses tolerance to reconstruct meshes)
fvMeshDistribute distributor(mesh, mergeDist);
// Refinement engine
// ~~~~~~~~~~~~~~~~~
Info<< nl
<< "Determining initial surface intersections" << nl
<< "-----------------------------------------" << nl
<< endl;
// Main refinement engine
meshRefinement meshRefiner
(
mesh,
mergeDist, // tolerance used in sorting coordinates
surfaces, // for surface intersection refinement
shells // for volume (inside/outside) refinement
);
Info<< "Calculated surface intersections in = "
<< mesh.time().cpuTimeIncrement() << " s" << nl << endl;
// Some stats
meshRefiner.printMeshInfo(debug, "Initial mesh");
meshRefiner.write
(
debug&meshRefinement::OBJINTERSECTIONS,
mesh.time().path()/mesh.time().timeName()
);
......@@ -403,11 +402,7 @@ int main(int argc, char *argv[])
if (wantLayers)
{
autoLayerDriver layerDriver
(
meshRefiner,
globalToPatch
);
autoLayerDriver layerDriver(meshRefiner);
// Layer addition parameters
layerParameters layerParams(layerDict, mesh.boundaryMesh());
......@@ -435,7 +430,7 @@ int main(int argc, char *argv[])
Info<< "End\n" << endl;
return 0;
return(0);
}
......
......@@ -292,6 +292,40 @@ Foam::autoHexMeshDriver::autoHexMeshDriver
meshRefinement::checkCoupledFaceZones(mesh_);
// Refinement engine
// ~~~~~~~~~~~~~~~~~
{
Info<< nl
<< "Determining initial surface intersections" << nl
<< "-----------------------------------------" << nl
<< endl;
// Main refinement engine
meshRefinerPtr_.reset
(
new meshRefinement
(
mesh,
mergeDist_, // tolerance used in sorting coordinates
surfaces(),
shells()
)
);
Info<< "Calculated surface intersections in = "
<< mesh_.time().cpuTimeIncrement() << " s" << endl;
// Some stats
meshRefinerPtr_().printMeshInfo(debug_, "Initial mesh");
meshRefinerPtr_().write
(
debug_&meshRefinement::OBJINTERSECTIONS,
mesh_.time().path()/mesh_.time().timeName()
);
}
// Add all the surface regions as patches
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
......@@ -319,9 +353,8 @@ Foam::autoHexMeshDriver::autoHexMeshDriver
forAll(regNames, i)
{
label patchI = meshRefinement::addPatch
label patchI = meshRefinerPtr_().addMeshedPatch
(
mesh,
regNames[i],
wallPolyPatch::typeName
);
......@@ -404,40 +437,6 @@ Foam::autoHexMeshDriver::autoHexMeshDriver
// Mesh distribution engine (uses tolerance to reconstruct meshes)
distributorPtr_.reset(new fvMeshDistribute(mesh_, mergeDist_));
}
// Refinement engine
// ~~~~~~~~~~~~~~~~~
{
Info<< nl
<< "Determining initial surface intersections" << nl
<< "-----------------------------------------" << nl
<< endl;
// Main refinement engine
meshRefinerPtr_.reset
(
new meshRefinement
(
mesh,
mergeDist_, // tolerance used in sorting coordinates
surfaces(),
shells()
)
);
Info<< "Calculated surface intersections in = "
<< mesh_.time().cpuTimeIncrement() << " s" << endl;
// Some stats
meshRefinerPtr_().printMeshInfo(debug_, "Initial mesh");
meshRefinerPtr_().write
(
debug_&meshRefinement::OBJINTERSECTIONS,
mesh_.time().path()/mesh_.time().timeName()
);
}
}
......@@ -522,11 +521,7 @@ void Foam::autoHexMeshDriver::doMesh()
const dictionary& shrinkDict = dict_.subDict("shrinkDict");
PtrList<dictionary> surfaceDicts(dict_.lookup("surfaces"));
autoLayerDriver layerDriver
(
meshRefinerPtr_(),
globalToPatch_
);
autoLayerDriver layerDriver(meshRefinerPtr_());
// Get all the layer specific params
layerParameters layerParams
......
......@@ -75,7 +75,7 @@ Foam::label Foam::autoLayerDriver::mergePatchFacesUndo
labelHashSet boundaryCells(mesh.nFaces()-mesh.nInternalFaces());
{
labelList patchIDs(meshRefinement::addedPatches(globalToPatch_));
labelList patchIDs(meshRefiner_.meshedPatches());
const polyBoundaryMesh& patches = mesh.boundaryMesh();
......@@ -2446,14 +2446,9 @@ void Foam::autoLayerDriver::getLayerCellsFaces
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::autoLayerDriver::autoLayerDriver
(
meshRefinement& meshRefiner,
const labelList& globalToPatch
)
Foam::autoLayerDriver::autoLayerDriver(meshRefinement& meshRefiner)
:
meshRefiner_(meshRefiner),
globalToPatch_(globalToPatch)
meshRefiner_(meshRefiner)
{}
......
......@@ -101,9 +101,6 @@ class autoLayerDriver
//- Mesh+surface
meshRefinement& meshRefiner_;
//- From surface region to patch
const labelList globalToPatch_;
// Private Member Functions
......@@ -509,11 +506,7 @@ public:
// Constructors
//- Construct from components
autoLayerDriver
(
meshRefinement& meshRefiner,
const labelList& globalToPatch
);
autoLayerDriver(meshRefinement& meshRefiner);
// Member Functions
......
......@@ -680,7 +680,7 @@ void Foam::autoRefineDriver::mergePatchFaces
(
Foam::cos(45*mathematicalConstant::pi/180.0),
Foam::cos(45*mathematicalConstant::pi/180.0),
meshRefinement::addedPatches(globalToPatch_)
meshRefiner_.meshedPatches()
);
if (debug)
......
......@@ -830,35 +830,6 @@ Foam::scalarField Foam::autoSnapDriver::calcSnapDistance
}
//// Invert globalToPatch_ to get the patches related to surfaces.
//Foam::labelList Foam::autoSnapDriver::getSurfacePatches() const
//{
// // Set of patches originating from surface
// labelHashSet surfacePatchSet(globalToPatch_.size());
//
// forAll(globalToPatch_, i)
// {
// if (globalToPatch_[i] != -1)
// {
// surfacePatchSet.insert(globalToPatch_[i]);
// }
// }
//
// const fvMesh& mesh = meshRefiner_.mesh();
//
// DynamicList<label> surfacePatches(surfacePatchSet.size());
//
// for (label patchI = 0; patchI < mesh.boundaryMesh().size(); patchI++)
// {
// if (surfacePatchSet.found(patchI))
// {
// surfacePatches.append(patchI);
// }
// }
// return surfacePatches.shrink();
//}
void Foam::autoSnapDriver::preSmoothPatch
(
const snapParameters& snapParams,
......@@ -1479,7 +1450,7 @@ void Foam::autoSnapDriver::doSnap
const_cast<Time&>(mesh.time())++;
// Get the labels of added patches.
labelList adaptPatchIDs(meshRefinement::addedPatches(globalToPatch_));
labelList adaptPatchIDs(meshRefiner_.meshedPatches());
// Create baffles (pairs of faces that share the same points)
// Baffles stored as owner and neighbour face that have been created.
......
......@@ -170,9 +170,6 @@ public:
const indirectPrimitivePatch&
) const;
////- Get patches generated for surfaces.
//labelList getSurfacePatches() const;
//- Smooth the mesh (patch and internal) to increase visibility
// of surface points (on castellated mesh) w.r.t. surface.
void preSmoothPatch
......
......@@ -84,12 +84,15 @@ void Foam::meshRefinement::calcNeighbourData
const polyBoundaryMesh& patches = mesh_.boundaryMesh();
labelHashSet addedPatchIDSet(meshedPatches());
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
const unallocLabelList& faceCells = pp.faceCells();
const vectorField::subField faceCentres = pp.faceCentres();
const vectorField::subField faceAreas = pp.faceAreas();
label bFaceI = pp.start()-mesh_.nInternalFaces();
......@@ -102,6 +105,36 @@ void Foam::meshRefinement::calcNeighbourData
bFaceI++;
}
}
else if (addedPatchIDSet.found(patchI))
{
// Face was introduced from cell-cell intersection. Try to
// reconstruct other side cell(centre). Three possibilities:
// - cells same size.
// - preserved cell smaller. Not handled.
// - preserved cell larger.
forAll(faceCells, i)
{
// Extrapolate the face centre.
vector fn = faceAreas[i];
fn /= mag(fn)+VSMALL;
// Normal distance from face centre to cell centre
scalar d = ((faceCentres[i] - cellCentres[i]) & fn);
label own = faceCells[i];
label ownLevel = cellLevel[own];
label faceLevel = meshCutter_.getAnchorLevel(pp.start()+i);
if (faceLevel > ownLevel)
{
// Other cell more refined. Adjust normal distance
d *= 0.5;
}
neiLevel[bFaceI] = cellLevel[ownLevel];
// Calculate other cell centre by extrapolation
neiCc[bFaceI] = faceCentres[i] + d*fn;
bFaceI++;
}
}
else
{
forAll(faceCells, i)
......@@ -1195,7 +1228,6 @@ Foam::labelList Foam::meshRefinement::intersectedFaces() const
// Helper function to get points used by faces
Foam::labelList Foam::meshRefinement::intersectedPoints
(
// const labelList& globalToPatch
) const
{
const faceList& faces = mesh_.faces();
......@@ -1221,9 +1253,10 @@ Foam::labelList Foam::meshRefinement::intersectedPoints
}
//// Insert all meshed patches.
//forAll(globalToPatch, i)
//labelList adaptPatchIDs(meshedPatches());
//forAll(adaptPatchIDs, i)
//{
// label patchI = globalToPatch[i];
// label patchI = adaptPatchIDs[i];
//
// if (patchI != -1)
// {
......@@ -1262,27 +1295,6 @@ Foam::labelList Foam::meshRefinement::intersectedPoints
}
Foam::labelList Foam::meshRefinement::addedPatches
(
const labelList& globalToPatch
)
{
labelList patchIDs(globalToPatch.size());
label addedI = 0;
forAll(globalToPatch, i)
{
if (globalToPatch[i] != -1)
{
patchIDs[addedI++] = globalToPatch[i];
}
}
patchIDs.setSize(addedI);
return patchIDs;
}
//- Create patch from set of patches
Foam::autoPtr<Foam::indirectPrimitivePatch> Foam::meshRefinement::makePatch
(
......@@ -1653,6 +1665,53 @@ Foam::label Foam::meshRefinement::addPatch
}
Foam::label Foam::meshRefinement::addMeshedPatch
(
const word& name,
const word& type
)
{
label meshedI = findIndex(meshedPatches_, name);
if (meshedI != -1)
{
// Already there. Get corresponding polypatch
return mesh_.boundaryMesh().findPatchID(name);
}
else
{
// Add patch
label patchI = addPatch(mesh_, name, type);
// Store
label sz = meshedPatches_.size();
meshedPatches_.setSize(sz+1);
meshedPatches_[sz] = name;
return patchI;
}
}
Foam::labelList Foam::meshRefinement::meshedPatches() const
{
labelList patchIDs(meshedPatches_.size());
forAll(meshedPatches_, i)
{
patchIDs[i] = mesh_.boundaryMesh().findPatchID(meshedPatches_[i]);
if (patchIDs[i] == -1)
{
FatalErrorIn("meshRefinement::meshedPatches() const")
<< "Problem : did not find patch " << meshedPatches_[i]
<< abort(FatalError);
}
}
return patchIDs;
}
Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::splitMeshRegions
(
const point& keepPoint
......
......@@ -125,6 +125,10 @@ private:
//- user supplied face based data.
List<Tuple2<mapType, labelList> > userFaceData_;
//- Meshed patches - are treated differently. Stored as wordList since
// order changes.
wordList meshedPatches_;
// Private Member Functions
......@@ -400,12 +404,11 @@ private:
const labelList& globalToPatch
) const;
//- Initial test of marking faces using geometric information.
labelList markFacesOnProblemCellsGeometric
(
const dictionary& motionDict,
const labelList& globalToPatch
) const;
////- Initial test of marking faces using geometric information.
//labelList markFacesOnProblemCellsGeometric
//(
// const dictionary& motionDict
//) const;
// Baffle merging
......@@ -578,9 +581,6 @@ public:
//- Get points on surfaces with intersection and boundary faces.
labelList intersectedPoints() const;
//- Get added patches (inverse of globalToPatch)
static labelList addedPatches(const labelList& globalToPatch);
//- Create patch from set of patches
static autoPtr<indirectPrimitivePatch> makePatch
(
......@@ -688,9 +688,16 @@ public:
// Other topo changes
//- Helper function to add patch to mesh
//- Helper:add patch to mesh. Update all registered fields.
// Use addMeshedPatch to add patches originating from surfaces.
static label addPatch(fvMesh&, const word& name, const word& type);
//- Add patch originating from meshing. Update meshedPatches_.
label addMeshedPatch(const word& name, const word& type);
//- Get patchIDs for patches added in addMeshedPatch.
labelList meshedPatches() const;
//- Split mesh. Keep part containing point.
autoPtr<mapPolyMesh> splitMeshRegions(const point& keepPoint);
......@@ -699,7 +706,11 @@ public:
//- Update for external change to mesh. changedFaces are in new mesh
// face labels.
void updateMesh(const mapPolyMesh&, const labelList& changedFaces);
void updateMesh
(
const mapPolyMesh&,
const labelList& changedFaces
);
// Restoring : is where other processes delete and reinsert data.
......
......@@ -1511,11 +1511,7 @@ void Foam::meshRefinement::baffleAndSplitMesh
perpendicularAngle,
globalToPatch
)
//markFacesOnProblemCellsGeometric
//(
// motionDict,
// globalToPatch
//)
//markFacesOnProblemCellsGeometric(motionDict)
);
Info<< "Analyzed problem cells in = "
<< runTime.cpuTimeIncrement() << " s\n" << nl << endl;
......@@ -1665,7 +1661,7 @@ void Foam::meshRefinement::baffleAndSplitMesh
// Split off (with optional buffer layers) unreachable areas of mesh.
Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::splitMesh
Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::splitMesh
(
const label nBufferLayers,
const labelList& globalToPatch,
......
......@@ -136,15 +136,13 @@ Foam::Map<Foam::label> Foam::meshRefinement::findEdgeConnectedProblemCells
const labelList& globalToPatch
) const
{