From ea1591128670b9fc9a68085f11689b5c5961421e Mon Sep 17 00:00:00 2001 From: mattijs Date: Thu, 14 Jul 2011 22:19:26 +0100 Subject: [PATCH] ENH: foamToTetDualMesh: map data onto tetDualMesh mesh --- .../foamToTetDualMesh/Make/files | 4 + .../foamToTetDualMesh/Make/options | 5 + .../foamToTetDualMesh/foamToTetDualMesh.C | 307 ++++++++++++++++++ 3 files changed, 316 insertions(+) create mode 100644 applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/Make/files create mode 100644 applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/Make/options create mode 100644 applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/foamToTetDualMesh.C diff --git a/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/Make/files b/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/Make/files new file mode 100644 index 0000000000..54d6d3e589 --- /dev/null +++ b/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/Make/files @@ -0,0 +1,4 @@ + +foamToTetDualMesh.C + +EXE = $(FOAM_USER_APPBIN)/foamToTetDualMesh diff --git a/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/Make/options b/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/Make/options new file mode 100644 index 0000000000..fa15f12452 --- /dev/null +++ b/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/Make/options @@ -0,0 +1,5 @@ +EXE_INC = \ + -I$(LIB_SRC)/finiteVolume/lnInclude + +EXE_LIBS = \ + -lfiniteVolume diff --git a/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/foamToTetDualMesh.C b/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/foamToTetDualMesh.C new file mode 100644 index 0000000000..716f419440 --- /dev/null +++ b/applications/utilities/postProcessing/dataConversion/foamToTetDualMesh/foamToTetDualMesh.C @@ -0,0 +1,307 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 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 . + +Application + foamToTetDualMesh + +Description + Converts polyMesh results to tetDualMesh. + +\*---------------------------------------------------------------------------*/ + + +#include "argList.H" +#include "fvMesh.H" +#include "volFields.H" +#include "pointFields.H" +#include "Time.H" +#include "IOobjectList.H" + +using namespace Foam; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +template +void ReadAndMapFields +( + const fvMesh& mesh, + const IOobjectList& objects, + const fvMesh& tetDualMesh, + const labelList& map, + const typename MappedGeoField::value_type& nullValue, + PtrList& tetFields +) +{ + typedef typename MappedGeoField::value_type Type; + + // Search list of objects for wanted type + IOobjectList fieldObjects(objects.lookupClass(ReadGeoField::typeName)); + + tetFields.setSize(fieldObjects.size()); + + label i = 0; + forAllConstIter(IOobjectList, fieldObjects, iter) + { + Info<< "Converting " << ReadGeoField::typeName << ' ' << iter.key() + << endl; + + ReadGeoField readField(*iter(), mesh); + + tetFields.set + ( + i, + new MappedGeoField + ( + IOobject + ( + readField.name(), + readField.instance(), + readField.local(), + tetDualMesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE, + readField.registerObject() + ), + pointMesh::New(tetDualMesh), + dimensioned + ( + "zero", + readField.dimensions(), + pTraits::zero + ) + ) + ); + + Field& fld = tetFields[i].internalField(); + + // Map from read field. Set unmapped entries to nullValue. + fld.setSize(map.size(), nullValue); + forAll(map, pointI) + { + label index = map[pointI]; + + if (index > 0) + { + label cellI = index-1; + fld[pointI] = readField[cellI]; + } + else if (index < 0) + { + label faceI = -index-1; + label bFaceI = faceI - mesh.nInternalFaces(); + if (bFaceI < 0) + { + FatalErrorIn("ReadAndMapFields(..)") + << "Face " << faceI << " from index " << index + << " is not a boundary face." << abort(FatalError); + } + + label patchI = mesh.boundaryMesh().patchID()[bFaceI]; + label localFaceI = mesh.boundaryMesh()[patchI].whichFace(faceI); + fld[pointI] = readField.boundaryField()[patchI][localFaceI]; + } + //else + //{ + // WarningIn("ReadAndMapFields(..)") + // << "Point " << pointI << " at " + // << tetDualMesh.points()[pointI] + // << " has no dual correspondence." << endl; + //} + } + tetFields[i].correctBoundaryConditions(); + + i++; + } +} + + + +// Main program: + +int main(int argc, char *argv[]) +{ +# include "addOverwriteOption.H" +# include "addTimeOptions.H" + +# include "setRootCase.H" +# include "createTime.H" + // Get times list + instantList Times = runTime.times(); +# include "checkTimeOptions.H" + runTime.setTime(Times[startTime], startTime); + + + // Read the mesh +# include "createMesh.H" + + // Read the tetDualMesh + Info<< "Create tetDualMesh for time = " + << runTime.timeName() << nl << endl; + + fvMesh tetDualMesh + ( + IOobject + ( + "tetDualMesh", + runTime.timeName(), + runTime, + IOobject::MUST_READ + ) + ); + // From tet vertices to poly cells/faces + const labelIOList pointDualAddressing + ( + IOobject + ( + "pointDualAddressing", + tetDualMesh.facesInstance(), + tetDualMesh.meshSubDir, + tetDualMesh, + IOobject::MUST_READ + ) + ); + + if (pointDualAddressing.size() != tetDualMesh.nPoints()) + { + FatalErrorIn(args.executable()) + << "Size " << pointDualAddressing.size() + << " of addressing map " << pointDualAddressing.objectPath() + << " differs from number of points in mesh " + << tetDualMesh.nPoints() + << exit(FatalError); + } + + + // Some stats on addressing + label nCells = 0; + label nPatchFaces = 0; + label nUnmapped = 0; + forAll(pointDualAddressing, pointI) + { + label index = pointDualAddressing[pointI]; + + if (index > 0) + { + nCells++; + } + else if (index == 0) + { + nUnmapped++; + } + else + { + label faceI = -index-1; + if (faceI < tetDualMesh.nInternalFaces()) + { + FatalErrorIn(args.executable()) + << "Face " << faceI << " from index " << index + << " is not a boundary face." + << " nInternalFaces:" << tetDualMesh.nInternalFaces() + << exit(FatalError); + } + else + { + nPatchFaces++; + } + } + } + + reduce(nCells, sumOp