diff --git a/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMesh.C b/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMesh.C index bafe30550467b75f30c33effec03008170addf4e..cdd6f52d15fce2954c01764e5f51700260adb9f3 100644 --- a/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMesh.C +++ b/applications/utilities/mesh/generation/extrude/extrudeMesh/extrudeMesh.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -383,12 +383,7 @@ int main(int argc, char *argv[]) // Determine extrudePatch normal pointField extrudePatchPointNormals ( - PatchTools::pointNormals //calcNormals - ( - mesh, - extrudePatch, - meshFaces - ) + PatchTools::pointNormals(mesh, extrudePatch) ); @@ -629,12 +624,13 @@ int main(int argc, char *argv[]) const labelListList& layerFaces = layerExtrude.layerFaces(); backPatchFaces.setSize(layerFaces.size()); frontPatchFaces.setSize(layerFaces.size()); - forAll(backPatchFaces, i) + forAll(backPatchFaces, patchFaceI) { - backPatchFaces[i] = layerFaces[i].first(); - frontPatchFaces[i] = layerFaces[i].last(); + backPatchFaces[patchFaceI] = layerFaces[patchFaceI].first(); + frontPatchFaces[patchFaceI] = layerFaces[patchFaceI].last(); } + // Create dummy fvSchemes, fvSolution createDummyFvMeshFiles(mesh, regionDir); @@ -654,6 +650,13 @@ int main(int argc, char *argv[]) mesh ); + layerExtrude.updateMesh + ( + map(), + identity(extrudePatch.size()), + identity(extrudePatch.nPoints()) + ); + // Calculate face labels for front and back. frontPatchFaces = renumber ( diff --git a/applications/utilities/mesh/generation/extrude/extrudeToRegionMesh/Make/options b/applications/utilities/mesh/generation/extrude/extrudeToRegionMesh/Make/options index 9a1bd9942c97cc352b8d14e5a92e62de35537a7a..49cad25f210163dffa9ff87fe916b6e6dffc4b10 100644 --- a/applications/utilities/mesh/generation/extrude/extrudeToRegionMesh/Make/options +++ b/applications/utilities/mesh/generation/extrude/extrudeToRegionMesh/Make/options @@ -1,10 +1,12 @@ EXE_INC = \ + -I$(LIB_SRC)/surfMesh/lnInclude \ -I$(LIB_SRC)/finiteVolume/lnInclude \ -I$(LIB_SRC)/meshTools/lnInclude \ -I$(LIB_SRC)/dynamicMesh/lnInclude \ -I$(LIB_SRC)/mesh/extrudeModel/lnInclude EXE_LIBS = \ + -lsurfMesh \ -lfiniteVolume \ -lmeshTools \ -ldynamicMesh \ diff --git a/applications/utilities/mesh/generation/extrude/extrudeToRegionMesh/extrudeToRegionMesh.C b/applications/utilities/mesh/generation/extrude/extrudeToRegionMesh/extrudeToRegionMesh.C index b159f9c096dd5e8d5813568a2fe52a9d3abfeee3..112de47faf92be168361d9d626ee8a027baea00f 100644 --- a/applications/utilities/mesh/generation/extrude/extrudeToRegionMesh/extrudeToRegionMesh.C +++ b/applications/utilities/mesh/generation/extrude/extrudeToRegionMesh/extrudeToRegionMesh.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -133,6 +133,8 @@ Notes: #include "pointFields.H" //#include "ReadFields.H" #include "fvMeshTools.H" +#include "OBJstream.H" +#include "PatchTools.H" using namespace Foam; @@ -1233,6 +1235,252 @@ void setCouplingInfo } +// Extrude and write geometric properties +void extrudeGeometricProperties +( + const polyMesh& mesh, + const primitiveFacePatch& extrudePatch, + const createShellMesh& extruder, + const polyMesh& regionMesh, + const extrudeModel& model +) +{ + const pointIOField patchFaceCentres + ( + IOobject + ( + "patchFaceCentres", + mesh.pointsInstance(), + mesh.meshSubDir, + mesh, + IOobject::MUST_READ + ) + ); + + const pointIOField patchEdgeCentres + ( + IOobject + ( + "patchEdgeCentres", + mesh.pointsInstance(), + mesh.meshSubDir, + mesh, + IOobject::MUST_READ + ) + ); + + //forAll(extrudePatch.edges(), edgeI) + //{ + // const edge& e = extrudePatch.edges()[edgeI]; + // Pout<< "Edge:" << e.centre(extrudePatch.localPoints()) << nl + // << "read:" << patchEdgeCentres[edgeI] + // << endl; + //} + + + // Determine edge normals on original patch + labelList patchEdges; + labelList coupledEdges; + PackedBoolList sameEdgeOrientation; + PatchTools::matchEdges + ( + extrudePatch, + mesh.globalData().coupledPatch(), + patchEdges, + coupledEdges, + sameEdgeOrientation + ); + + pointField patchEdgeNormals + ( + PatchTools::edgeNormals + ( + mesh, + extrudePatch, + patchEdges, + coupledEdges + ) + ); + + + pointIOField faceCentres + ( + IOobject + ( + "faceCentres", + regionMesh.pointsInstance(), + regionMesh.meshSubDir, + regionMesh, + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + regionMesh.nFaces() + ); + + + // Work out layers. Guaranteed in columns so no fancy parallel bits. + + + forAll(extruder.faceToFaceMap(), faceI) + { + if (extruder.faceToFaceMap()[faceI] != 0) + { + // 'horizontal' face + label patchFaceI = mag(extruder.faceToFaceMap()[faceI])-1; + + label cellI = regionMesh.faceOwner()[faceI]; + if (regionMesh.isInternalFace(faceI)) + { + cellI = max(cellI, regionMesh.faceNeighbour()[faceI]); + } + + // Calculate layer from cell numbering (see createShellMesh) + label layerI = (cellI % model.nLayers()); + + if + ( + !regionMesh.isInternalFace(faceI) + && extruder.faceToFaceMap()[faceI] > 0 + ) + { + // Top face + layerI++; + } + + + // Recalculate based on extrusion model + faceCentres[faceI] = model + ( + patchFaceCentres[patchFaceI], + extrudePatch.faceNormals()[patchFaceI], + layerI + ); + } + else + { + // 'vertical face + label patchEdgeI = extruder.faceToEdgeMap()[faceI]; + label layerI = + ( + regionMesh.faceOwner()[faceI] + % model.nLayers() + ); + + // Extrude patch edge centre to this layer + point pt0 = model + ( + patchEdgeCentres[patchEdgeI], + patchEdgeNormals[patchEdgeI], + layerI + ); + // Extrude patch edge centre to next layer + point pt1 = model + ( + patchEdgeCentres[patchEdgeI], + patchEdgeNormals[patchEdgeI], + layerI+1 + ); + + // Interpolate + faceCentres[faceI] = 0.5*(pt0+pt1); + } + } + + pointIOField cellCentres + ( + IOobject + ( + "cellCentres", + regionMesh.pointsInstance(), + regionMesh.meshSubDir, + regionMesh, + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + regionMesh.nCells() + ); + + forAll(extruder.cellToFaceMap(), cellI) + { + label patchFaceI = extruder.cellToFaceMap()[cellI]; + + // Calculate layer from cell numbering (see createShellMesh) + label layerI = (cellI % model.nLayers()); + + // Recalculate based on extrusion model + point pt0 = model + ( + patchFaceCentres[patchFaceI], + extrudePatch.faceNormals()[patchFaceI], + layerI + ); + point pt1 = model + ( + patchFaceCentres[patchFaceI], + extrudePatch.faceNormals()[patchFaceI], + layerI+1 + ); + + // Interpolate + cellCentres[cellI] = 0.5*(pt0+pt1); + } + + + // Bit of checking + if (false) + { + OBJstream faceStr(regionMesh.time().path()/"faceCentres.obj"); + OBJstream cellStr(regionMesh.time().path()/"cellCentres.obj"); + + forAll(faceCentres, faceI) + { + Pout<< "Model :" << faceCentres[faceI] << endl + << "regionMesh:" << regionMesh.faceCentres()[faceI] << endl; + faceStr.write + ( + linePointRef + ( + faceCentres[faceI], + regionMesh.faceCentres()[faceI] + ) + ); + } + forAll(cellCentres, cellI) + { + Pout<< "Model :" << cellCentres[cellI] << endl + << "regionMesh:" << regionMesh.cellCentres()[cellI] << endl; + cellStr.write + ( + linePointRef + ( + cellCentres[cellI], + regionMesh.cellCentres()[cellI] + ) + ); + } + } + + + + Info<< "Writing geometric properties for mesh " << regionMesh.name() + << " to " << regionMesh.pointsInstance() << nl + << endl; + + bool ok = faceCentres.write() && cellCentres.write(); + + if (!ok) + { + FatalErrorIn("extrudeGeometricProperties(..)") + << "Failed writing " << faceCentres.objectPath() + << " and " << cellCentres.objectPath() + << exit(FatalError); + } +} + + + // Main program: int main(int argc, char *argv[]) @@ -2393,6 +2641,36 @@ int main(int argc, char *argv[]) } + // See if we need to extrude coordinates as well + { + autoPtr<pointIOField> patchFaceCentresPtr; + + IOobject io + ( + "patchFaceCentres", + mesh.pointsInstance(), + mesh.meshSubDir, + mesh, + IOobject::MUST_READ + ); + if (io.headerOk()) + { + // Read patchFaceCentres and patchEdgeCentres + Info<< "Reading patch face,edge centres : " + << io.name() << " and patchEdgeCentres" << endl; + + extrudeGeometricProperties + ( + mesh, + extrudePatch, + extruder, + regionMesh, + model() + ); + } + } + + // Insert baffles into original mesh diff --git a/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C b/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C index 12e4ebb6eb99b59da4dcd6f82c1fcbbddca44e6f..49ee1a445e70abd754853161e3d791d9aabe4786 100644 --- a/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C +++ b/src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -2052,6 +2052,9 @@ void Foam::polyTopoChange::reorderCoupledFaces if (anyChanged) { + // Reorder faces according to oldToNew. + reorderCompactFaces(oldToNew.size(), oldToNew); + // Rotate faces (rotation is already in new face indices). forAll(rotation, faceI) { @@ -2060,9 +2063,6 @@ void Foam::polyTopoChange::reorderCoupledFaces inplaceRotateList<List, label>(faces_[faceI], rotation[faceI]); } } - - // Reorder faces according to oldToNew. - reorderCompactFaces(oldToNew.size(), oldToNew); } }