Merge branch 'master' of ssh://dm/home/dm4/OpenFOAM/OpenFOAM-dev

This commit is contained in:
Henry
2013-09-19 15:11:27 +01:00
19 changed files with 908 additions and 759 deletions

View File

@ -8,6 +8,7 @@ conformalVoronoiMesh/indexedCell/indexedCellEnum.C
conformalVoronoiMesh/conformalVoronoiMesh.C conformalVoronoiMesh/conformalVoronoiMesh.C
conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
conformalVoronoiMesh/conformalVoronoiMeshZones.C
conformalVoronoiMesh/conformalVoronoiMeshIO.C conformalVoronoiMesh/conformalVoronoiMeshIO.C
conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C

View File

@ -29,10 +29,10 @@ Description
SourceFiles SourceFiles
conformalVoronoiMeshI.H conformalVoronoiMeshI.H
conformalVoronoiMesh.C conformalVoronoiMesh.C
conformalVoronoiMeshZones.C
conformalVoronoiMeshIO.C conformalVoronoiMeshIO.C
conformalVoronoiMeshConformToSurface.C conformalVoronoiMeshConformToSurface.C
conformalVoronoiMeshFeaturePoints.C conformalVoronoiMeshFeaturePoints.C
conformalVoronoiMeshFeaturePointSpecialisations.C
conformalVoronoiMeshCalcDualMesh.C conformalVoronoiMeshCalcDualMesh.C
conformalVoronoiMeshTemplates.C conformalVoronoiMeshTemplates.C
@ -636,6 +636,9 @@ private:
boolList& flipMap boolList& flipMap
) const; ) const;
//- Add zones to the polyMesh
void addZones(polyMesh& mesh, const pointField& cellCentres) const;
//- Tet mesh calculation //- Tet mesh calculation
void calcTetMesh void calcTetMesh
( (

View File

@ -35,12 +35,7 @@ License
#include "pointMesh.H" #include "pointMesh.H"
#include "indexedVertexOps.H" #include "indexedVertexOps.H"
#include "DelaunayMeshTools.H" #include "DelaunayMeshTools.H"
#include "surfaceZonesInfo.H"
#include "polyModifyCell.H"
#include "polyModifyFace.H"
#include "syncTools.H" #include "syncTools.H"
#include "regionSplit.H"
#include "OBJstream.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@ -408,543 +403,6 @@ void Foam::conformalVoronoiMesh::writeMesh(const fileName& instance)
} }
void Foam::conformalVoronoiMesh::calcNeighbourCellCentres
(
const polyMesh& mesh,
const pointField& cellCentres,
pointField& neiCc
) const
{
label nBoundaryFaces = mesh.nFaces() - mesh.nInternalFaces();
if (neiCc.size() != nBoundaryFaces)
{
FatalErrorIn("conformalVoronoiMesh::calcNeighbourCellCentres(..)")
<< "nBoundaries:" << nBoundaryFaces
<< " neiCc:" << neiCc.size()
<< abort(FatalError);
}
const polyBoundaryMesh& patches = mesh.boundaryMesh();
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
const labelUList& faceCells = pp.faceCells();
label bFaceI = pp.start() - mesh.nInternalFaces();
if (pp.coupled())
{
forAll(faceCells, i)
{
neiCc[bFaceI] = cellCentres[faceCells[i]];
bFaceI++;
}
}
}
// Swap coupled boundaries. Apply separation to cc since is coordinate.
syncTools::swapBoundaryFacePositions(mesh, neiCc);
}
void Foam::conformalVoronoiMesh::selectSeparatedCoupledFaces
(
const polyMesh& mesh,
boolList& selected
) const
{
const polyBoundaryMesh& patches = mesh.boundaryMesh();
forAll(patches, patchI)
{
// Check all coupled. Avoid using .coupled() so we also pick up AMI.
if (isA<coupledPolyPatch>(patches[patchI]))
{
const coupledPolyPatch& cpp = refCast<const coupledPolyPatch>
(
patches[patchI]
);
if (cpp.separated() || !cpp.parallel())
{
forAll(cpp, i)
{
selected[cpp.start()+i] = true;
}
}
}
}
}
void Foam::conformalVoronoiMesh::findCellZoneInsideWalk
(
const polyMesh& mesh,
const labelList& locationSurfaces, // indices of surfaces with inside point
const labelList& faceToSurface, // per face index of named surface
labelList& cellToSurface
) const
{
// Analyse regions. Reuse regionsplit
boolList blockedFace(mesh.nFaces());
selectSeparatedCoupledFaces(mesh, blockedFace);
forAll(faceToSurface, faceI)
{
if (faceToSurface[faceI] == -1)
{
blockedFace[faceI] = false;
}
else
{
blockedFace[faceI] = true;
}
}
// No need to sync since namedSurfaceIndex already is synced
// Set region per cell based on walking
regionSplit cellRegion(mesh, blockedFace);
blockedFace.clear();
// Force calculation of face decomposition (used in findCell)
(void)mesh.tetBasePtIs();
const PtrList<surfaceZonesInfo>& surfZones =
geometryToConformTo().surfZones();
// For all locationSurface find the cell
forAll(locationSurfaces, i)
{
label surfI = locationSurfaces[i];
const Foam::point& insidePoint = surfZones[surfI].zoneInsidePoint();
const word& surfName = geometryToConformTo().geometry().names()[surfI];
Info<< " For surface " << surfName
<< " finding inside point " << insidePoint
<< endl;
// Find the region containing the insidePoint
label keepRegionI = -1;
label cellI = mesh.findCell(insidePoint);
if (cellI != -1)
{
keepRegionI = cellRegion[cellI];
}
reduce(keepRegionI, maxOp<label>());
Info<< " For surface " << surfName
<< " found point " << insidePoint << " in cell " << cellI
<< " in global region " << keepRegionI
<< " out of " << cellRegion.nRegions() << " regions." << endl;
if (keepRegionI == -1)
{
FatalErrorIn
(
"conformalVoronoiMesh::findCellZoneInsideWalk"
"(const polyMesh&, const labelList&"
", const labelList&, labelList&)"
) << "Point " << insidePoint
<< " is not inside the mesh." << nl
<< "Bounding box of the mesh:" << mesh.bounds()
<< exit(FatalError);
}
// Set all cells with this region
forAll(cellRegion, cellI)
{
if (cellRegion[cellI] == keepRegionI)
{
if (cellToSurface[cellI] == -2)
{
cellToSurface[cellI] = surfI;
}
else if (cellToSurface[cellI] != surfI)
{
WarningIn
(
"conformalVoronoiMesh::findCellZoneInsideWalk"
"(const labelList&, const labelList&"
", const labelList&, const labelList&)"
) << "Cell " << cellI
<< " at " << mesh.cellCentres()[cellI]
<< " is inside surface " << surfName
<< " but already marked as being in zone "
<< cellToSurface[cellI] << endl
<< "This can happen if your surfaces are not"
<< " (sufficiently) closed."
<< endl;
}
}
}
}
}
Foam::labelList Foam::conformalVoronoiMesh::calcCellZones
(
const pointField& cellCentres
) const
{
labelList cellToSurface(cellCentres.size(), -1);
const PtrList<surfaceZonesInfo>& surfZones =
geometryToConformTo().surfZones();
// Get list of closed surfaces
labelList closedNamedSurfaces
(
surfaceZonesInfo::getAllClosedNamedSurfaces
(
surfZones,
geometryToConformTo().geometry(),
geometryToConformTo().surfaces()
)
);
forAll(closedNamedSurfaces, i)
{
label surfI = closedNamedSurfaces[i];
const searchableSurface& surface =
allGeometry()[geometryToConformTo().surfaces()[surfI]];
const surfaceZonesInfo::areaSelectionAlgo selectionMethod =
surfZones[surfI].zoneInside();
if
(
selectionMethod != surfaceZonesInfo::INSIDE
&& selectionMethod != surfaceZonesInfo::OUTSIDE
&& selectionMethod != surfaceZonesInfo::INSIDEPOINT
)
{
FatalErrorIn("conformalVoronoiMesh::calcCellZones(..)")
<< "Trying to use surface "
<< surface.name()
<< " which has non-geometric inside selection method "
<< surfaceZonesInfo::areaSelectionAlgoNames[selectionMethod]
<< exit(FatalError);
}
if (surface.hasVolumeType())
{
List<volumeType> volType;
surface.getVolumeType(cellCentres, volType);
bool selectInside = true;
if (selectionMethod == surfaceZonesInfo::INSIDEPOINT)
{
List<volumeType> volTypeInsidePoint;
surface.getVolumeType
(
pointField(1, surfZones[surfI].zoneInsidePoint()),
volTypeInsidePoint
);
if (volTypeInsidePoint[0] == volumeType::OUTSIDE)
{
selectInside = false;
}
}
else if (selectionMethod == surfaceZonesInfo::OUTSIDE)
{
selectInside = false;
}
forAll(volType, pointI)
{
if (cellToSurface[pointI] == -1)
{
if
(
(
volType[pointI] == volumeType::INSIDE
&& selectInside
)
|| (
volType[pointI] == volumeType::OUTSIDE
&& !selectInside
)
)
{
cellToSurface[pointI] = surfI;
}
}
}
}
}
return cellToSurface;
}
void Foam::conformalVoronoiMesh::calcFaceZones
(
const polyMesh& mesh,
const pointField& cellCentres,
const labelList& cellToSurface,
labelList& faceToSurface,
boolList& flipMap
) const
{
faceToSurface.setSize(mesh.nFaces(), -1);
flipMap.setSize(mesh.nFaces(), false);
const faceList& faces = mesh.faces();
const labelList& faceOwner = mesh.faceOwner();
const labelList& faceNeighbour = mesh.faceNeighbour();
labelList neiFaceOwner(mesh.nFaces() - mesh.nInternalFaces(), -1);
const polyBoundaryMesh& patches = mesh.boundaryMesh();
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
const labelUList& faceCells = pp.faceCells();
label bFaceI = pp.start() - mesh.nInternalFaces();
if (pp.coupled())
{
forAll(faceCells, i)
{
neiFaceOwner[bFaceI] = cellToSurface[faceCells[i]];
bFaceI++;
}
}
}
syncTools::swapBoundaryFaceList(mesh, neiFaceOwner);
forAll(faces, faceI)
{
const label ownerSurfaceI = cellToSurface[faceOwner[faceI]];
if (faceToSurface[faceI] >= 0)
{
continue;
}
if (mesh.isInternalFace(faceI))
{
const label neiSurfaceI = cellToSurface[faceNeighbour[faceI]];
if
(
(ownerSurfaceI >= 0 || neiSurfaceI >= 0)
&& ownerSurfaceI != neiSurfaceI
)
{
flipMap[faceI] =
(
ownerSurfaceI == max(ownerSurfaceI, neiSurfaceI)
? false
: true
);
faceToSurface[faceI] = max(ownerSurfaceI, neiSurfaceI);
}
}
else
{
label patchID = mesh.boundaryMesh().whichPatch(faceI);
if (mesh.boundaryMesh()[patchID].coupled())
{
const label neiSurfaceI =
neiFaceOwner[faceI - mesh.nInternalFaces()];
if
(
(ownerSurfaceI >= 0 || neiSurfaceI >= 0)
&& ownerSurfaceI != neiSurfaceI
)
{
flipMap[faceI] =
(
ownerSurfaceI == max(ownerSurfaceI, neiSurfaceI)
? false
: true
);
faceToSurface[faceI] = max(ownerSurfaceI, neiSurfaceI);
}
}
else
{
if (ownerSurfaceI >= 0)
{
faceToSurface[faceI] = ownerSurfaceI;
}
}
}
}
const PtrList<surfaceZonesInfo>& surfZones =
geometryToConformTo().surfZones();
labelList unclosedSurfaces
(
surfaceZonesInfo::getUnclosedNamedSurfaces
(
surfZones,
geometryToConformTo().geometry(),
geometryToConformTo().surfaces()
)
);
pointField neiCc(mesh.nFaces() - mesh.nInternalFaces());
calcNeighbourCellCentres
(
mesh,
cellCentres,
neiCc
);
OBJstream intersections(time().path()/"ints.obj");
OBJstream intersectionFaces(time().path()/"intFaces.obj");
// Use intersection of cellCentre connections
forAll(faces, faceI)
{
if (faceToSurface[faceI] >= 0)
{
continue;
}
label patchID = mesh.boundaryMesh().whichPatch(faceI);
const label own = faceOwner[faceI];
List<pointIndexHit> surfHit;
labelList hitSurface;
if (mesh.isInternalFace(faceI))
{
const label nei = faceNeighbour[faceI];
geometryToConformTo().findSurfaceAllIntersections
(
cellCentres[own],
cellCentres[nei],
surfHit,
hitSurface
);
}
else if (patchID != -1 && mesh.boundaryMesh()[patchID].coupled())
{
geometryToConformTo().findSurfaceAllIntersections
(
cellCentres[own],
neiCc[faceI - mesh.nInternalFaces()],
surfHit,
hitSurface
);
if (surfHit.size() == 1 && surfHit[0].hit())
{
intersections.write
(
linePointRef
(
cellCentres[own],
neiCc[faceI - mesh.nInternalFaces()]
)
);
}
}
// If there are multiple intersections then do not add to
// a faceZone
if (surfHit.size() == 1 && surfHit[0].hit())
{
if (findIndex(unclosedSurfaces, hitSurface[0]) != -1)
{
vectorField norm;
geometryToConformTo().getNormal
(
hitSurface[0],
List<pointIndexHit>(1, surfHit[0]),
norm
);
vector fN = faces[faceI].normal(mesh.points());
fN /= mag(fN) + SMALL;
if ((norm[0] & fN) < 0)
{
flipMap[faceI] = true;
}
else
{
flipMap[faceI] = false;
}
faceToSurface[faceI] = hitSurface[0];
intersectionFaces.write(faces[faceI], mesh.points());
}
}
}
// labelList neiCellSurface(mesh.nFaces()-mesh.nInternalFaces());
//
// forAll(patches, patchI)
// {
// const polyPatch& pp = patches[patchI];
//
// if (pp.coupled())
// {
// forAll(pp, i)
// {
// label faceI = pp.start()+i;
// label ownSurface = cellToSurface[faceOwner[faceI]];
// neiCellSurface[faceI - mesh.nInternalFaces()] = ownSurface;
// }
// }
// }
// syncTools::swapBoundaryFaceList(mesh, neiCellSurface);
//
// forAll(patches, patchI)
// {
// const polyPatch& pp = patches[patchI];
//
// if (pp.coupled())
// {
// forAll(pp, i)
// {
// label faceI = pp.start()+i;
// label ownSurface = cellToSurface[faceOwner[faceI]];
// label neiSurface =
// neiCellSurface[faceI-mesh.nInternalFaces()];
//
// if (faceToSurface[faceI] == -1 && (ownSurface != neiSurface))
// {
// // Give face the max cell zone
// faceToSurface[faceI] = max(ownSurface, neiSurface);
// }
// }
// }
// }
//
// Sync
syncTools::syncFaceList(mesh, faceToSurface, maxEqOp<label>());
}
Foam::autoPtr<Foam::fvMesh> Foam::conformalVoronoiMesh::createDummyMesh Foam::autoPtr<Foam::fvMesh> Foam::conformalVoronoiMesh::createDummyMesh
( (
const IOobject& io, const IOobject& io,
@ -1304,7 +762,7 @@ void Foam::conformalVoronoiMesh::reorderProcessorPatches
{ {
if (rotation[faceI] != 0) if (rotation[faceI] != 0)
{ {
inplaceRotateList<List, label>(faces[faceI], rotation[faceI]); faces[faceI] = rotateList(faces[faceI], rotation[faceI]);
nRotated++; nRotated++;
} }
} }
@ -1448,159 +906,9 @@ void Foam::conformalVoronoiMesh::writeMesh
mesh.addFvPatches(patches); mesh.addFvPatches(patches);
// Add zones to mesh
{
Info<< " Adding zones to mesh" << endl;
const PtrList<surfaceZonesInfo>& surfZones =
geometryToConformTo().surfZones();
labelList cellToSurface(calcCellZones(cellCentres));
labelList faceToSurface;
boolList flipMap;
calcFaceZones
(
mesh,
cellCentres,
cellToSurface,
faceToSurface,
flipMap
);
labelList insidePointNamedSurfaces
(
surfaceZonesInfo::getInsidePointNamedSurfaces(surfZones)
);
findCellZoneInsideWalk
(
mesh,
insidePointNamedSurfaces,
faceToSurface,
cellToSurface
);
labelList namedSurfaces(surfaceZonesInfo::getNamedSurfaces(surfZones));
forAll(namedSurfaces, i)
{
label surfI = namedSurfaces[i];
Info<< incrIndent << indent << "Surface : "
<< geometryToConformTo().geometry().names()[surfI] << nl
<< indent << " faceZone : "
<< surfZones[surfI].faceZoneName() << nl
<< indent << " cellZone : "
<< surfZones[surfI].cellZoneName()
<< decrIndent << endl;
}
// Add zones to mesh
labelList surfaceToFaceZone =
surfaceZonesInfo::addFaceZonesToMesh
(
surfZones,
namedSurfaces,
mesh
);
labelList surfaceToCellZone =
surfaceZonesInfo::addCellZonesToMesh
(
surfZones,
namedSurfaces,
mesh
);
// Topochange container
polyTopoChange meshMod(mesh);
forAll(cellToSurface, cellI)
{
label surfaceI = cellToSurface[cellI];
if (surfaceI >= 0)
{
label zoneI = surfaceToCellZone[surfaceI];
if (zoneI >= 0)
{
meshMod.setAction
(
polyModifyCell
(
cellI,
false, // removeFromZone
zoneI
)
);
}
}
}
const labelList& faceOwner = mesh.faceOwner();
const labelList& faceNeighbour = mesh.faceNeighbour();
forAll(faceToSurface, faceI)
{
label surfaceI = faceToSurface[faceI];
if (surfaceI < 0)
{
continue;
}
label patchID = mesh.boundaryMesh().whichPatch(faceI);
if (mesh.isInternalFace(faceI))
{
label own = faceOwner[faceI];
label nei = faceNeighbour[faceI];
meshMod.setAction
(
polyModifyFace
(
mesh.faces()[faceI], // modified face
faceI, // label of face
own, // owner
nei, // neighbour
false, // face flip
-1, // patch for face
false, // remove from zone
surfaceToFaceZone[surfaceI], // zone for face
flipMap[faceI] // face flip in zone
)
);
}
else if (patchID != -1 && mesh.boundaryMesh()[patchID].coupled())
{
label own = faceOwner[faceI];
meshMod.setAction
(
polyModifyFace
(
mesh.faces()[faceI], // modified face
faceI, // label of face
own, // owner
-1, // neighbour
false, // face flip
patchID, // patch for face
false, // remove from zone
surfaceToFaceZone[surfaceI], // zone for face
flipMap[faceI] // face flip in zone
)
);
}
}
// Change the mesh (no inflation, parallel sync)
autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false, true);
}
// Add zones to the mesh
addZones(mesh, cellCentres);
@ -1611,12 +919,18 @@ void Foam::conformalVoronoiMesh::writeMesh
forAll(boundaryFacesToRemove, faceI) forAll(boundaryFacesToRemove, faceI)
{ {
if (boundaryFacesToRemove[faceI]) if
(
boundaryFacesToRemove[faceI]
&& mesh.faceZones().whichZone(faceI) == -1
)
{ {
addr[count++] = faceI; addr[count++] = faceI;
} }
} }
addr.setSize(count);
label sz = mesh.faceZones().size(); label sz = mesh.faceZones().size();
boolList flip(addr.size(), false); boolList flip(addr.size(), false);
mesh.faceZones().setSize(sz + 1); mesh.faceZones().setSize(sz + 1);

View File

@ -0,0 +1,713 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013 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 "conformalVoronoiMesh.H"
#include "polyModifyFace.H"
#include "polyModifyCell.H"
#include "syncTools.H"
#include "regionSplit.H"
#include "surfaceZonesInfo.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::conformalVoronoiMesh::calcNeighbourCellCentres
(
const polyMesh& mesh,
const pointField& cellCentres,
pointField& neiCc
) const
{
label nBoundaryFaces = mesh.nFaces() - mesh.nInternalFaces();
if (neiCc.size() != nBoundaryFaces)
{
FatalErrorIn("conformalVoronoiMesh::calcNeighbourCellCentres(..)")
<< "nBoundaries:" << nBoundaryFaces
<< " neiCc:" << neiCc.size()
<< abort(FatalError);
}
const polyBoundaryMesh& patches = mesh.boundaryMesh();
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
const labelUList& faceCells = pp.faceCells();
label bFaceI = pp.start() - mesh.nInternalFaces();
if (pp.coupled())
{
forAll(faceCells, i)
{
neiCc[bFaceI] = cellCentres[faceCells[i]];
bFaceI++;
}
}
}
// Swap coupled boundaries. Apply separation to cc since is coordinate.
syncTools::swapBoundaryFacePositions(mesh, neiCc);
}
void Foam::conformalVoronoiMesh::selectSeparatedCoupledFaces
(
const polyMesh& mesh,
boolList& selected
) const
{
const polyBoundaryMesh& patches = mesh.boundaryMesh();
forAll(patches, patchI)
{
// Check all coupled. Avoid using .coupled() so we also pick up AMI.
if (isA<coupledPolyPatch>(patches[patchI]))
{
const coupledPolyPatch& cpp = refCast<const coupledPolyPatch>
(
patches[patchI]
);
if (cpp.separated() || !cpp.parallel())
{
forAll(cpp, i)
{
selected[cpp.start()+i] = true;
}
}
}
}
}
void Foam::conformalVoronoiMesh::findCellZoneInsideWalk
(
const polyMesh& mesh,
const labelList& locationSurfaces, // indices of surfaces with inside point
const labelList& faceToSurface, // per face index of named surface
labelList& cellToSurface
) const
{
// Analyse regions. Reuse regionsplit
boolList blockedFace(mesh.nFaces());
selectSeparatedCoupledFaces(mesh, blockedFace);
forAll(faceToSurface, faceI)
{
if (faceToSurface[faceI] == -1)
{
blockedFace[faceI] = false;
}
else
{
blockedFace[faceI] = true;
}
}
// No need to sync since namedSurfaceIndex already is synced
// Set region per cell based on walking
regionSplit cellRegion(mesh, blockedFace);
blockedFace.clear();
// Force calculation of face decomposition (used in findCell)
(void)mesh.tetBasePtIs();
const PtrList<surfaceZonesInfo>& surfZones =
geometryToConformTo().surfZones();
// For all locationSurface find the cell
forAll(locationSurfaces, i)
{
label surfI = locationSurfaces[i];
const Foam::point& insidePoint = surfZones[surfI].zoneInsidePoint();
const word& surfName = geometryToConformTo().geometry().names()[surfI];
Info<< " For surface " << surfName
<< " finding inside point " << insidePoint
<< endl;
// Find the region containing the insidePoint
label keepRegionI = -1;
label cellI = mesh.findCell(insidePoint);
if (cellI != -1)
{
keepRegionI = cellRegion[cellI];
}
reduce(keepRegionI, maxOp<label>());
Info<< " For surface " << surfName
<< " found point " << insidePoint << " in cell " << cellI
<< " in global region " << keepRegionI
<< " out of " << cellRegion.nRegions() << " regions." << endl;
if (keepRegionI == -1)
{
FatalErrorIn
(
"conformalVoronoiMesh::findCellZoneInsideWalk"
"(const polyMesh&, const labelList&"
", const labelList&, labelList&)"
) << "Point " << insidePoint
<< " is not inside the mesh." << nl
<< "Bounding box of the mesh:" << mesh.bounds()
<< exit(FatalError);
}
// Set all cells with this region
forAll(cellRegion, cellI)
{
if (cellRegion[cellI] == keepRegionI)
{
if (cellToSurface[cellI] == -2)
{
cellToSurface[cellI] = surfI;
}
else if (cellToSurface[cellI] != surfI)
{
WarningIn
(
"conformalVoronoiMesh::findCellZoneInsideWalk"
"(const labelList&, const labelList&"
", const labelList&, const labelList&)"
) << "Cell " << cellI
<< " at " << mesh.cellCentres()[cellI]
<< " is inside surface " << surfName
<< " but already marked as being in zone "
<< cellToSurface[cellI] << endl
<< "This can happen if your surfaces are not"
<< " (sufficiently) closed."
<< endl;
}
}
}
}
}
Foam::labelList Foam::conformalVoronoiMesh::calcCellZones
(
const pointField& cellCentres
) const
{
labelList cellToSurface(cellCentres.size(), -1);
const PtrList<surfaceZonesInfo>& surfZones =
geometryToConformTo().surfZones();
// Get list of closed surfaces
labelList closedNamedSurfaces
(
surfaceZonesInfo::getAllClosedNamedSurfaces
(
surfZones,
geometryToConformTo().geometry(),
geometryToConformTo().surfaces()
)
);
forAll(closedNamedSurfaces, i)
{
label surfI = closedNamedSurfaces[i];
const searchableSurface& surface =
allGeometry()[geometryToConformTo().surfaces()[surfI]];
const surfaceZonesInfo::areaSelectionAlgo selectionMethod =
surfZones[surfI].zoneInside();
if
(
selectionMethod != surfaceZonesInfo::INSIDE
&& selectionMethod != surfaceZonesInfo::OUTSIDE
&& selectionMethod != surfaceZonesInfo::INSIDEPOINT
)
{
FatalErrorIn("conformalVoronoiMesh::calcCellZones(..)")
<< "Trying to use surface "
<< surface.name()
<< " which has non-geometric inside selection method "
<< surfaceZonesInfo::areaSelectionAlgoNames[selectionMethod]
<< exit(FatalError);
}
if (surface.hasVolumeType())
{
List<volumeType> volType;
surface.getVolumeType(cellCentres, volType);
bool selectInside = true;
if (selectionMethod == surfaceZonesInfo::INSIDEPOINT)
{
List<volumeType> volTypeInsidePoint;
surface.getVolumeType
(
pointField(1, surfZones[surfI].zoneInsidePoint()),
volTypeInsidePoint
);
if (volTypeInsidePoint[0] == volumeType::OUTSIDE)
{
selectInside = false;
}
}
else if (selectionMethod == surfaceZonesInfo::OUTSIDE)
{
selectInside = false;
}
forAll(volType, pointI)
{
if (cellToSurface[pointI] == -1)
{
if
(
(
volType[pointI] == volumeType::INSIDE
&& selectInside
)
|| (
volType[pointI] == volumeType::OUTSIDE
&& !selectInside
)
)
{
cellToSurface[pointI] = surfI;
}
}
}
}
}
return cellToSurface;
}
void Foam::conformalVoronoiMesh::calcFaceZones
(
const polyMesh& mesh,
const pointField& cellCentres,
const labelList& cellToSurface,
labelList& faceToSurface,
boolList& flipMap
) const
{
faceToSurface.setSize(mesh.nFaces(), -1);
flipMap.setSize(mesh.nFaces(), false);
const faceList& faces = mesh.faces();
const labelList& faceOwner = mesh.faceOwner();
const labelList& faceNeighbour = mesh.faceNeighbour();
labelList neiFaceOwner(mesh.nFaces() - mesh.nInternalFaces(), -1);
const polyBoundaryMesh& patches = mesh.boundaryMesh();
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
const labelUList& faceCells = pp.faceCells();
label bFaceI = pp.start() - mesh.nInternalFaces();
if (pp.coupled())
{
forAll(faceCells, i)
{
neiFaceOwner[bFaceI] = cellToSurface[faceCells[i]];
bFaceI++;
}
}
}
syncTools::swapBoundaryFaceList(mesh, neiFaceOwner);
forAll(faces, faceI)
{
const label ownerSurfaceI = cellToSurface[faceOwner[faceI]];
if (faceToSurface[faceI] >= 0)
{
continue;
}
if (mesh.isInternalFace(faceI))
{
const label neiSurfaceI = cellToSurface[faceNeighbour[faceI]];
if
(
(ownerSurfaceI >= 0 || neiSurfaceI >= 0)
&& ownerSurfaceI != neiSurfaceI
)
{
flipMap[faceI] =
(
ownerSurfaceI == max(ownerSurfaceI, neiSurfaceI)
? false
: true
);
faceToSurface[faceI] = max(ownerSurfaceI, neiSurfaceI);
}
}
else
{
label patchID = mesh.boundaryMesh().whichPatch(faceI);
if (mesh.boundaryMesh()[patchID].coupled())
{
const label neiSurfaceI =
neiFaceOwner[faceI - mesh.nInternalFaces()];
if
(
(ownerSurfaceI >= 0 || neiSurfaceI >= 0)
&& ownerSurfaceI != neiSurfaceI
)
{
flipMap[faceI] =
(
ownerSurfaceI == max(ownerSurfaceI, neiSurfaceI)
? false
: true
);
faceToSurface[faceI] = max(ownerSurfaceI, neiSurfaceI);
}
}
else
{
if (ownerSurfaceI >= 0)
{
faceToSurface[faceI] = ownerSurfaceI;
}
}
}
}
const PtrList<surfaceZonesInfo>& surfZones =
geometryToConformTo().surfZones();
labelList unclosedSurfaces
(
surfaceZonesInfo::getUnclosedNamedSurfaces
(
surfZones,
geometryToConformTo().geometry(),
geometryToConformTo().surfaces()
)
);
pointField neiCc(mesh.nFaces() - mesh.nInternalFaces());
calcNeighbourCellCentres
(
mesh,
cellCentres,
neiCc
);
// Use intersection of cellCentre connections
forAll(faces, faceI)
{
if (faceToSurface[faceI] >= 0)
{
continue;
}
label patchID = mesh.boundaryMesh().whichPatch(faceI);
const label own = faceOwner[faceI];
List<pointIndexHit> surfHit;
labelList hitSurface;
if (mesh.isInternalFace(faceI))
{
const label nei = faceNeighbour[faceI];
geometryToConformTo().findSurfaceAllIntersections
(
cellCentres[own],
cellCentres[nei],
surfHit,
hitSurface
);
}
else if (patchID != -1 && mesh.boundaryMesh()[patchID].coupled())
{
geometryToConformTo().findSurfaceAllIntersections
(
cellCentres[own],
neiCc[faceI - mesh.nInternalFaces()],
surfHit,
hitSurface
);
}
// If there are multiple intersections then do not add to
// a faceZone
if (surfHit.size() == 1 && surfHit[0].hit())
{
if (findIndex(unclosedSurfaces, hitSurface[0]) != -1)
{
vectorField norm;
geometryToConformTo().getNormal
(
hitSurface[0],
List<pointIndexHit>(1, surfHit[0]),
norm
);
vector fN = faces[faceI].normal(mesh.points());
fN /= mag(fN) + SMALL;
if ((norm[0] & fN) < 0)
{
flipMap[faceI] = true;
}
else
{
flipMap[faceI] = false;
}
faceToSurface[faceI] = hitSurface[0];
}
}
}
// labelList neiCellSurface(mesh.nFaces()-mesh.nInternalFaces());
//
// forAll(patches, patchI)
// {
// const polyPatch& pp = patches[patchI];
//
// if (pp.coupled())
// {
// forAll(pp, i)
// {
// label faceI = pp.start()+i;
// label ownSurface = cellToSurface[faceOwner[faceI]];
// neiCellSurface[faceI - mesh.nInternalFaces()] = ownSurface;
// }
// }
// }
// syncTools::swapBoundaryFaceList(mesh, neiCellSurface);
//
// forAll(patches, patchI)
// {
// const polyPatch& pp = patches[patchI];
//
// if (pp.coupled())
// {
// forAll(pp, i)
// {
// label faceI = pp.start()+i;
// label ownSurface = cellToSurface[faceOwner[faceI]];
// label neiSurface =
// neiCellSurface[faceI-mesh.nInternalFaces()];
//
// if (faceToSurface[faceI] == -1 && (ownSurface != neiSurface))
// {
// // Give face the max cell zone
// faceToSurface[faceI] = max(ownSurface, neiSurface);
// }
// }
// }
// }
// Sync
syncTools::syncFaceList(mesh, faceToSurface, maxEqOp<label>());
}
void Foam::conformalVoronoiMesh::addZones
(
polyMesh& mesh,
const pointField& cellCentres
) const
{
Info<< " Adding zones to mesh" << endl;
const PtrList<surfaceZonesInfo>& surfZones =
geometryToConformTo().surfZones();
labelList cellToSurface(calcCellZones(cellCentres));
labelList faceToSurface;
boolList flipMap;
calcFaceZones
(
mesh,
cellCentres,
cellToSurface,
faceToSurface,
flipMap
);
labelList insidePointNamedSurfaces
(
surfaceZonesInfo::getInsidePointNamedSurfaces(surfZones)
);
findCellZoneInsideWalk
(
mesh,
insidePointNamedSurfaces,
faceToSurface,
cellToSurface
);
labelList namedSurfaces(surfaceZonesInfo::getNamedSurfaces(surfZones));
forAll(namedSurfaces, i)
{
label surfI = namedSurfaces[i];
Info<< incrIndent << indent << "Surface : "
<< geometryToConformTo().geometry().names()[surfI] << nl
<< indent << " faceZone : "
<< surfZones[surfI].faceZoneName() << nl
<< indent << " cellZone : "
<< surfZones[surfI].cellZoneName()
<< decrIndent << endl;
}
// Add zones to mesh
labelList surfaceToFaceZone =
surfaceZonesInfo::addFaceZonesToMesh
(
surfZones,
namedSurfaces,
mesh
);
labelList surfaceToCellZone =
surfaceZonesInfo::addCellZonesToMesh
(
surfZones,
namedSurfaces,
mesh
);
// Topochange container
polyTopoChange meshMod(mesh);
forAll(cellToSurface, cellI)
{
label surfaceI = cellToSurface[cellI];
if (surfaceI >= 0)
{
label zoneI = surfaceToCellZone[surfaceI];
if (zoneI >= 0)
{
meshMod.setAction
(
polyModifyCell
(
cellI,
false, // removeFromZone
zoneI
)
);
}
}
}
const labelList& faceOwner = mesh.faceOwner();
const labelList& faceNeighbour = mesh.faceNeighbour();
forAll(faceToSurface, faceI)
{
label surfaceI = faceToSurface[faceI];
if (surfaceI < 0)
{
continue;
}
label patchID = mesh.boundaryMesh().whichPatch(faceI);
if (mesh.isInternalFace(faceI))
{
label own = faceOwner[faceI];
label nei = faceNeighbour[faceI];
meshMod.setAction
(
polyModifyFace
(
mesh.faces()[faceI], // modified face
faceI, // label of face
own, // owner
nei, // neighbour
false, // face flip
-1, // patch for face
false, // remove from zone
surfaceToFaceZone[surfaceI], // zone for face
flipMap[faceI] // face flip in zone
)
);
}
else if (patchID != -1 && mesh.boundaryMesh()[patchID].coupled())
{
label own = faceOwner[faceI];
meshMod.setAction
(
polyModifyFace
(
mesh.faces()[faceI], // modified face
faceI, // label of face
own, // owner
-1, // neighbour
false, // face flip
patchID, // patch for face
false, // remove from zone
surfaceToFaceZone[surfaceI], // zone for face
flipMap[faceI] // face flip in zone
)
);
}
}
// Change the mesh (no inflation, parallel sync)
autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false, true);
}
// ************************************************************************* //

View File

@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
@ -36,7 +36,7 @@ Description
Storage is allocated on free-store during construction. Storage is allocated on free-store during construction.
As a special case a null-contructed CompactListList has an empty As a special case a null-constructed CompactListList has an empty
offsets_ (instead of size 1). offsets_ (instead of size 1).
SourceFiles SourceFiles

View File

@ -63,23 +63,23 @@ bool Foam::IOobject::IOobject::fileNameComponents
// called with directory // called with directory
if (isDir(path)) if (isDir(path))
{ {
WarningIn("IOobject::fileNameComponents(const fileName&, ...)") WarningIn
<< " called with directory: " << path << "\n"; (
"IOobject::fileNameComponents"
"("
"const fileName&, "
"fileName&, "
"fileName&, "
"word&"
")"
)
<< " called with directory: " << path << endl;
return false; return false;
} }
string::size_type first = path.find('/'); if (path.isAbsolute())
if (first == string::npos)
{ {
// no '/' found - no instance or local
// check afterwards
name.string::operator=(path);
}
else if (first == 0)
{
// Leading '/'. Absolute fileName
string::size_type last = path.rfind('/'); string::size_type last = path.rfind('/');
instance = path.substr(0, last); instance = path.substr(0, last);
@ -88,26 +88,48 @@ bool Foam::IOobject::IOobject::fileNameComponents
} }
else else
{ {
instance = path.substr(0, first); string::size_type first = path.find('/');
string::size_type last = path.rfind('/'); if (first == string::npos)
if (last > first)
{ {
// with local // no '/' found - no instance or local
local = path.substr(first+1, last-first-1);
}
// check afterwards // check afterwards
name.string::operator=(path.substr(last+1)); name.string::operator=(path);
}
else
{
instance = path.substr(0, first);
string::size_type last = path.rfind('/');
if (last > first)
{
// with local
local = path.substr(first+1, last-first-1);
}
// check afterwards
name.string::operator=(path.substr(last+1));
}
} }
// check for valid (and stripped) name, regardless of the debug level // check for valid (and stripped) name, regardless of the debug level
if (name.empty() || string::stripInvalid<word>(name)) if (name.empty() || string::stripInvalid<word>(name))
{ {
WarningIn("IOobject::fileNameComponents(const fileName&, ...)") WarningIn
(
"IOobject::fileNameComponents"
"("
"const fileName&, "
"fileName&, "
"fileName&, "
"word&"
")"
)
<< "has invalid word for name: \"" << name << "has invalid word for name: \"" << name
<< "\"\nwhile processing path: " << path << "\n"; << "\"\nwhile processing path: " << path << endl;
return false; return false;
} }
@ -202,9 +224,16 @@ Foam::IOobject::IOobject
{ {
FatalErrorIn FatalErrorIn
( (
"IOobject::IOobject" "(const fileName&, const objectRegistry&, ...)" "IOobject::IOobject"
"("
"const fileName&, "
"const objectRegistry&, "
"readOption, "
"writeOption, "
"bool"
")"
) )
<< " invalid path specification\n" << " invalid path specification"
<< exit(FatalError); << exit(FatalError);
} }

View File

@ -294,15 +294,26 @@ void Foam::cyclicACMIFvPatchField<Type>::updateCoeffs()
} }
template<class Type>
void Foam::cyclicACMIFvPatchField<Type>::initEvaluate
(
const Pstream::commsTypes comms
)
{
// update non-overlap patch
const fvPatchField<Type>& npf = nonOverlapPatchField();
const_cast<fvPatchField<Type>&>(npf).evaluate(comms);
}
template<class Type> template<class Type>
void Foam::cyclicACMIFvPatchField<Type>::evaluate void Foam::cyclicACMIFvPatchField<Type>::evaluate
( (
const Pstream::commsTypes comms const Pstream::commsTypes comms
) )
{ {
// blend contrubutions from the coupled and non-overlap patches // blend contributions from the coupled and non-overlap patches
const fvPatchField<Type>& npf = nonOverlapPatchField(); const fvPatchField<Type>& npf = nonOverlapPatchField();
const_cast<fvPatchField<Type>&>(npf).evaluate();
coupledFvPatchField<Type>::evaluate(comms); coupledFvPatchField<Type>::evaluate(comms);
const Field<Type>& cpf = *this; const Field<Type>& cpf = *this;
@ -392,7 +403,7 @@ void Foam::cyclicACMIFvPatchField<Type>::manipulateMatrix
fvMatrix<Type>& matrix fvMatrix<Type>& matrix
) )
{ {
// blend contrubutions from the coupled and non-overlap patches // blend contributions from the coupled and non-overlap patches
const fvPatchField<Type>& npf = nonOverlapPatchField(); const fvPatchField<Type>& npf = nonOverlapPatchField();
const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask(); const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask();

View File

@ -206,6 +206,12 @@ public:
//- Update the coefficients associated with the patch field //- Update the coefficients associated with the patch field
void updateCoeffs(); void updateCoeffs();
//- Initialise the evaluation of the patch field
virtual void initEvaluate
(
const Pstream::commsTypes commsType
);
//- Evaluate the patch field //- Evaluate the patch field
virtual void evaluate virtual void evaluate
( (

View File

@ -397,6 +397,9 @@ public:
// inline const typename parcelType::forceType& forces() const; // inline const typename parcelType::forceType& forces() const;
inline const forceType& forces() const; inline const forceType& forces() const;
//- Return the optional particle forces
inline forceType& forces();
//- Optional cloud function objects //- Optional cloud function objects
inline functionType& functions(); inline functionType& functions();
@ -487,6 +490,9 @@ public:
//- Total linear kinetic energy in the system //- Total linear kinetic energy in the system
inline scalar linearKineticEnergyOfSystem() const; inline scalar linearKineticEnergyOfSystem() const;
//- Total rotational kinetic energy in the system
inline scalar rotationalKineticEnergyOfSystem() const;
//- Penetration for fraction [0-1] of the current total mass //- Penetration for fraction [0-1] of the current total mass
inline scalar penetration(const scalar fraction) const; inline scalar penetration(const scalar fraction) const;

View File

@ -156,6 +156,14 @@ Foam::KinematicCloud<CloudType>::forces() const
} }
template<class CloudType>
inline typename Foam::KinematicCloud<CloudType>::forceType&
Foam::KinematicCloud<CloudType>::forces()
{
return forces_;
}
template<class CloudType> template<class CloudType>
inline typename Foam::KinematicCloud<CloudType>::functionType& inline typename Foam::KinematicCloud<CloudType>::functionType&
Foam::KinematicCloud<CloudType>::functions() Foam::KinematicCloud<CloudType>::functions()

View File

@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
@ -38,7 +38,9 @@ Foam::ParticleForceList<CloudType>::ParticleForceList
PtrList<ParticleForce<CloudType> >(), PtrList<ParticleForce<CloudType> >(),
owner_(owner), owner_(owner),
mesh_(mesh), mesh_(mesh),
dict_(dictionary::null) dict_(dictionary::null),
calcCoupled_(true),
calcNonCoupled_(true)
{} {}
@ -54,7 +56,9 @@ Foam::ParticleForceList<CloudType>::ParticleForceList
PtrList<ParticleForce<CloudType> >(), PtrList<ParticleForce<CloudType> >(),
owner_(owner), owner_(owner),
mesh_(mesh), mesh_(mesh),
dict_(dict) dict_(dict),
calcCoupled_(true),
calcNonCoupled_(true)
{ {
if (readFields) if (readFields)
{ {
@ -151,9 +155,13 @@ Foam::forceSuSp Foam::ParticleForceList<CloudType>::calcCoupled
) const ) const
{ {
forceSuSp value(vector::zero, 0.0); forceSuSp value(vector::zero, 0.0);
forAll(*this, i)
if (calcCoupled_)
{ {
value += this->operator[](i).calcCoupled(p, dt, mass, Re, muc); forAll(*this, i)
{
value += this->operator[](i).calcCoupled(p, dt, mass, Re, muc);
}
} }
return value; return value;
@ -171,9 +179,13 @@ Foam::forceSuSp Foam::ParticleForceList<CloudType>::calcNonCoupled
) const ) const
{ {
forceSuSp value(vector::zero, 0.0); forceSuSp value(vector::zero, 0.0);
forAll(*this, i)
if (calcNonCoupled_)
{ {
value += this->operator[](i).calcNonCoupled(p, dt, mass, Re, muc); forAll(*this, i)
{
value += this->operator[](i).calcNonCoupled(p, dt, mass, Re, muc);
}
} }
return value; return value;

View File

@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
@ -64,6 +64,12 @@ class ParticleForceList
//- Forces dictionary //- Forces dictionary
const dictionary dict_; const dictionary dict_;
//- Calculate coupled forces flag
bool calcCoupled_;
//- Calculate non-coupled forces flag
bool calcNonCoupled_;
public: public:
@ -105,6 +111,12 @@ public:
//- Return the forces dictionary //- Return the forces dictionary
inline const dictionary& dict() const; inline const dictionary& dict() const;
//- Set the calcCoupled flag
inline void setCalcCoupled(bool flag);
//- Set the calcNonCoupled flag
inline void setCalcNonCoupled(bool flag);
// Evaluation // Evaluation

View File

@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
@ -53,4 +53,18 @@ inline const Foam::dictionary& Foam::ParticleForceList<CloudType>::dict() const
} }
template<class CloudType>
inline void Foam::ParticleForceList<CloudType>::setCalcCoupled(bool flag)
{
calcCoupled_ = flag;
}
template<class CloudType>
inline void Foam::ParticleForceList<CloudType>::setCalcNonCoupled(bool flag)
{
calcNonCoupled_ = flag;
}
// ************************************************************************* // // ************************************************************************* //

View File

@ -68,13 +68,11 @@ void Foam::SprayParcel<ParcelType>::calc
const CompositionModel<reactingCloudType>& composition = const CompositionModel<reactingCloudType>& composition =
td.cloud().composition(); td.cloud().composition();
bool coupled = td.cloud().solution().coupled();
// check if parcel belongs to liquid core // check if parcel belongs to liquid core
if (liquidCore() > 0.5) if (liquidCore() > 0.5)
{ {
// liquid core parcels should not interact with the gas // liquid core parcels should not experience coupled forces
td.cloud().solution().coupled() = false; td.cloud().forces().setCalcCoupled(false);
} }
// get old mixture composition // get old mixture composition
@ -138,8 +136,8 @@ void Foam::SprayParcel<ParcelType>::calc
} }
} }
// restore coupled // restore coupled forces
td.cloud().solution().coupled() = coupled; td.cloud().forces().setCalcCoupled(true);
} }
@ -246,7 +244,7 @@ void Foam::SprayParcel<ParcelType>::calcBreakup
const scalar mass = p.mass(); const scalar mass = p.mass();
const forceSuSp Fcp = forces.calcCoupled(p, dt, mass, Re, muAv); const forceSuSp Fcp = forces.calcCoupled(p, dt, mass, Re, muAv);
const forceSuSp Fncp = forces.calcNonCoupled(p, dt, mass, Re, muAv); const forceSuSp Fncp = forces.calcNonCoupled(p, dt, mass, Re, muAv);
scalar tMom = mass/(Fcp.Sp() + Fncp.Sp()); this->tMom() = mass/(Fcp.Sp() + Fncp.Sp());
const vector g = td.cloud().g().value(); const vector g = td.cloud().g().value();
@ -274,7 +272,7 @@ void Foam::SprayParcel<ParcelType>::calcBreakup
muAv, muAv,
Urel, Urel,
Urmag, Urmag,
tMom, this->tMom(),
dChild, dChild,
massChild massChild
) )

View File

@ -48,7 +48,7 @@ Foam::BlobsSheetAtomization<CloudType>::BlobsSheetAtomization
: :
AtomizationModel<CloudType>(am), AtomizationModel<CloudType>(am),
B_(am.B_), B_(am.B_),
angle_(am.B_) angle_(am.angle_)
{} {}

View File

@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
@ -112,7 +112,7 @@ public:
const scalar rho, const scalar rho,
const scalar mu, const scalar mu,
const scalar sigma, const scalar sigma,
const scalar volHlowRate, const scalar volFlowRate,
const scalar rhoAv, const scalar rhoAv,
const scalar Urel, const scalar Urel,
const vector& pos, const vector& pos,

View File

@ -52,7 +52,8 @@ namespace Foam
Foam::localAxesRotation::localAxesRotation Foam::localAxesRotation::localAxesRotation
( (
const dictionary& dict, const objectRegistry& orb const dictionary& dict,
const objectRegistry& orb
) )
: :
Rptr_(), Rptr_(),
@ -88,10 +89,10 @@ Foam::localAxesRotation::localAxesRotation
e3_() e3_()
{ {
FatalErrorIn("localAxesRotation(const dictionary&)") FatalErrorIn("localAxesRotation(const dictionary&)")
<< " localAxesRotation can not be contructed from dictionary " << " localAxesRotation can not be constructed from dictionary "
<< " use the construtctor : " << " use the construtctor : "
"(" "("
" const dictionary& dict, const objectRegistry& orb" " const dictionary&, const objectRegistry&"
")" ")"
<< exit(FatalIOError); << exit(FatalIOError);
} }
@ -112,7 +113,7 @@ Foam::vector Foam::localAxesRotation::transform(const vector& st) const
{ {
notImplemented notImplemented
( (
"vector Foam::localAxesRotation::transform(const vector&) const" "vector localAxesRotation::transform(const vector&) const"
); );
return vector::zero; return vector::zero;
} }
@ -122,7 +123,7 @@ Foam::vector Foam::localAxesRotation::invTransform(const vector& st) const
{ {
notImplemented notImplemented
( (
"vector Foam::localAxesRotation::invTransform(const vector&) const" "vector localAxesRotation::invTransform(const vector&) const"
); );
return vector::zero; return vector::zero;
} }
@ -135,7 +136,10 @@ Foam::tmp<Foam::vectorField> Foam::localAxesRotation::transform
{ {
if (Rptr_->size() != st.size()) if (Rptr_->size() != st.size())
{ {
FatalErrorIn("localAxesRotation::transform(const vectorField&)") FatalErrorIn
(
"tmp<vectorField> localAxesRotation::transform(const vectorField&)"
)
<< "vectorField st has different size to tensorField " << "vectorField st has different size to tensorField "
<< abort(FatalError); << abort(FatalError);
} }
@ -160,7 +164,13 @@ Foam::tmp<Foam::tensorField> Foam::localAxesRotation::transformTensor
{ {
if (Rptr_->size() != st.size()) if (Rptr_->size() != st.size())
{ {
FatalErrorIn("localAxesRotation::transformTensor(const tensorField&)") FatalErrorIn
(
"tmp<tensorField> localAxesRotation::transformTensor"
"("
"const tensorField&"
")"
)
<< "tensorField st has different size to tensorField Tr" << "tensorField st has different size to tensorField Tr"
<< abort(FatalError); << abort(FatalError);
} }
@ -173,7 +183,11 @@ Foam::tensor Foam::localAxesRotation::transformTensor
const tensor& st const tensor& st
) const ) const
{ {
notImplemented("tensor localAxesRotation::transformTensor() const"); notImplemented
(
"tensor localAxesRotation::transformTensor(const tensor&) const"
);
return tensor::zero; return tensor::zero;
} }
@ -186,7 +200,14 @@ Foam::tmp<Foam::tensorField> Foam::localAxesRotation::transformTensor
{ {
if (cellMap.size() != st.size()) if (cellMap.size() != st.size())
{ {
FatalErrorIn("localAxesRotation::transformTensor(const tensorField&)") FatalErrorIn
(
"tmp<tensorField> localAxesRotation::transformTensor"
"("
"const tensorField&"
"const labelList&"
")"
)
<< "tensorField st has different size to tensorField Tr" << "tensorField st has different size to tensorField Tr"
<< abort(FatalError); << abort(FatalError);
} }

View File

@ -85,6 +85,7 @@ private:
//- Is surface closed //- Is surface closed
mutable label surfaceClosed_; mutable label surfaceClosed_;
// Private Member Functions // Private Member Functions
////- Helper: find instance of files without header ////- Helper: find instance of files without header
@ -166,7 +167,7 @@ public:
//- Move points //- Move points
virtual void movePoints(const pointField&); virtual void movePoints(const pointField&);
//- Demand driven contruction of octree for boundary edges //- Demand driven construction of octree for boundary edges
const indexedOctree<treeDataEdge>& edgeTree() const; const indexedOctree<treeDataEdge>& edgeTree() const;
@ -260,6 +261,7 @@ public:
List<volumeType>& List<volumeType>&
) const; ) const;
// Other // Other
//- WIP. Store element-wise field. //- WIP. Store element-wise field.
@ -285,7 +287,6 @@ public:
IOstream::versionNumber ver, IOstream::versionNumber ver,
IOstream::compressionType cmp IOstream::compressionType cmp
) const; ) const;
}; };

View File

@ -155,7 +155,7 @@ createSetsAndZone floor -$tol $floorMax -$tol $floorMax -$tol $tol
echo "cellSet floorCells new faceToCell floorFaces owner" >> $tmpSetSet echo "cellSet floorCells new faceToCell floorFaces owner" >> $tmpSetSet
echo "faceZoneSet floorFaces new setsToFaceZone floorFaces floorCells" >> $tmpSetSet echo "faceZoneSet floorFaces new setsToFaceZone floorFaces floorCells" >> $tmpSetSet
setSet -batch $tmpSetSet > log.setSet.patchifyObstacles >/dev/null 2>&1 setSet -batch $tmpSetSet > log.setSet.patchifyObstacles 2>&1
# ************************************************************************* # *************************************************************************