diff --git a/applications/utilities/preProcessing/mapFieldsNew/Make/files b/applications/utilities/preProcessing/mapFieldsNew/Make/files new file mode 100644 index 0000000000000000000000000000000000000000..d4d439083d4a643251a51731916b639e36bde601 --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/Make/files @@ -0,0 +1,4 @@ +mapLagrangian.C +mapFields.C + +EXE = $(FOAM_APPBIN)/mapFieldsNew diff --git a/applications/utilities/preProcessing/mapFieldsNew/Make/options b/applications/utilities/preProcessing/mapFieldsNew/Make/options new file mode 100644 index 0000000000000000000000000000000000000000..b4ed87bf74c96bb804916b668df8d38ec702d335 --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/Make/options @@ -0,0 +1,13 @@ +EXE_INC = \ + -DFULLDEBUG -g -O0 \ + -I$(LIB_SRC)/meshTools/lnInclude \ + -I$(LIB_SRC)/lagrangian/basic/lnInclude \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + -I$(LIB_SRC)/sampling/lnInclude + +EXE_LIBS = \ + -lsampling \ + -lmeshTools \ + -llagrangian \ + -lfiniteVolume \ + -lgenericPatchFields diff --git a/applications/utilities/preProcessing/mapFieldsNew/MapLagrangianFields.H b/applications/utilities/preProcessing/mapFieldsNew/MapLagrangianFields.H new file mode 100644 index 0000000000000000000000000000000000000000..5f3c4230d5e78a132b25f6131c918a3c1725927c --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/MapLagrangianFields.H @@ -0,0 +1,184 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-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/>. + +InNamespace + Foam + +Description + Gets the indices of (source)particles that have been appended to the + target cloud and maps the lagrangian fields accordingly. + +\*---------------------------------------------------------------------------*/ + +#ifndef MapLagrangianFields_H +#define MapLagrangianFields_H + +#include "cloud.H" +#include "GeometricField.H" +#include "meshToMeshNew.H" +#include "IOobjectList.H" +#include "CompactIOField.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +//- Gets the indices of (source)particles that have been appended to the +// target cloud and maps the lagrangian fields accordingly. +template<class Type> +void MapLagrangianFields +( + const string& cloudName, + const IOobjectList& objects, + const polyMesh& meshTarget, + const labelList& addParticles +) +{ + { + IOobjectList fields = objects.lookupClass(IOField<Type>::typeName); + + forAllIter(IOobjectList, fields, fieldIter) + { + const word& fieldName = fieldIter()->name(); + + Info<< " mapping lagrangian field " << fieldName << endl; + + // Read field (does not need mesh) + IOField<Type> fieldSource(*fieldIter()); + + // Map + IOField<Type> fieldTarget + ( + IOobject + ( + fieldName, + meshTarget.time().timeName(), + cloud::prefix/cloudName, + meshTarget, + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + addParticles.size() + ); + + forAll(addParticles, i) + { + fieldTarget[i] = fieldSource[addParticles[i]]; + } + + // Write field + fieldTarget.write(); + } + } + + { + IOobjectList fieldFields = + objects.lookupClass(IOField<Field<Type> >::typeName); + + forAllIter(IOobjectList, fieldFields, fieldIter) + { + const word& fieldName = fieldIter()->name(); + + Info<< " mapping lagrangian fieldField " << fieldName << endl; + + // Read field (does not need mesh) + IOField<Field<Type> > fieldSource(*fieldIter()); + + // Map - use CompactIOField to automatically write in + // compact form for binary format. + CompactIOField<Field<Type>, Type> fieldTarget + ( + IOobject + ( + fieldName, + meshTarget.time().timeName(), + cloud::prefix/cloudName, + meshTarget, + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + addParticles.size() + ); + + forAll(addParticles, i) + { + fieldTarget[i] = fieldSource[addParticles[i]]; + } + + // Write field + fieldTarget.write(); + } + } + + { + IOobjectList fieldFields = + objects.lookupClass(CompactIOField<Field<Type>, Type>::typeName); + + forAllIter(IOobjectList, fieldFields, fieldIter) + { + Info<< " mapping lagrangian fieldField " + << fieldIter()->name() << endl; + + // Read field (does not need mesh) + CompactIOField<Field<Type>, Type> fieldSource(*fieldIter()); + + // Map + CompactIOField<Field<Type>, Type> fieldTarget + ( + IOobject + ( + fieldIter()->name(), + meshTarget.time().timeName(), + cloud::prefix/cloudName, + meshTarget, + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + addParticles.size() + ); + + forAll(addParticles, i) + { + fieldTarget[i] = fieldSource[addParticles[i]]; + } + + // Write field + fieldTarget.write(); + } + } +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/utilities/preProcessing/mapFieldsNew/MapMeshes.H b/applications/utilities/preProcessing/mapFieldsNew/MapMeshes.H new file mode 100644 index 0000000000000000000000000000000000000000..c03c1a76eab00e498d98cc6ce228c74a402ba05c --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/MapMeshes.H @@ -0,0 +1,131 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-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/>. + +\*---------------------------------------------------------------------------*/ + +#ifndef MapMeshes_H +#define MapMeshes_H + +#include "MapVolFields.H" +#include "mapLagrangian.H" +#include "UnMapped.H" +#include "pointMesh.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +template<template<class> class CombineOp> +void MapMesh +( + const meshToMeshNew& interp, + const HashSet<word>& selectedFields, + const bool noLagrangian +) +{ + { + const polyMesh& meshSource = interp.srcRegion(); + + // Search for list of objects for this time + IOobjectList objects(meshSource, meshSource.time().timeName()); + + // Map volFields + // ~~~~~~~~~~~~~ + MapVolFields<scalar> + ( + objects, + selectedFields, + interp, + CombineOp<scalar>() + ); + + MapVolFields<vector> + ( + objects, + selectedFields, + interp, + CombineOp<vector>() + ); + MapVolFields<sphericalTensor> + ( + objects, + selectedFields, + interp, + CombineOp<sphericalTensor>() + ); + MapVolFields<symmTensor> + ( + objects, + selectedFields, + interp, + CombineOp<symmTensor>() + ); + MapVolFields<tensor> + ( + objects, + selectedFields, + interp, + CombineOp<tensor>() + ); + } + + { + const polyMesh& meshTarget = interp.tgtRegion(); + + // Search for list of target objects for this time + IOobjectList objects(meshTarget, meshTarget.time().timeName()); + + // Mark surfaceFields as unmapped + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + UnMapped<surfaceScalarField>(objects); + UnMapped<surfaceVectorField>(objects); + UnMapped<surfaceSphericalTensorField>(objects); + UnMapped<surfaceSymmTensorField>(objects); + UnMapped<surfaceTensorField>(objects); + + // Mark pointFields as unmapped + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + UnMapped<pointScalarField>(objects); + UnMapped<pointVectorField>(objects); + UnMapped<pointSphericalTensorField>(objects); + UnMapped<pointSymmTensorField>(objects); + UnMapped<pointTensorField>(objects); + } + + if (!noLagrangian) + { + mapLagrangian(interp); + } +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/utilities/preProcessing/mapFieldsNew/MapVolFields.H b/applications/utilities/preProcessing/mapFieldsNew/MapVolFields.H new file mode 100644 index 0000000000000000000000000000000000000000..fc59c5733656d71fe6936699934a8b18ce6575ab --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/MapVolFields.H @@ -0,0 +1,104 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-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/>. + +\*---------------------------------------------------------------------------*/ + +#ifndef MapConsistentVolFields_H +#define MapConsistentVolFields_H + +#include "GeometricField.H" +#include "meshToMeshNew.H" +#include "IOobjectList.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +template<class Type, class CombineOp> +void MapVolFields +( + const IOobjectList& objects, + const HashSet<word>& selectedFields, + const meshToMeshNew& interp, + const CombineOp& cop +) +{ + typedef GeometricField<Type, fvPatchField, volMesh> fieldType; + + const fvMesh& meshSource = static_cast<const fvMesh&>(interp.srcRegion()); + const fvMesh& meshTarget = static_cast<const fvMesh&>(interp.tgtRegion()); + + IOobjectList fields = objects.lookupClass(fieldType::typeName); + + forAllIter(IOobjectList, fields, fieldIter) + { + const word& fieldName = fieldIter()->name(); + + if (selectedFields.empty() || selectedFields.found(fieldName)) + { + Info<< " interpolating " << fieldName << endl; + + const fieldType fieldSource(*fieldIter(), meshSource); + + IOobject targetIO + ( + fieldName, + meshTarget.time().timeName(), + meshTarget, + IOobject::MUST_READ + ); + + if (targetIO.headerOk()) + { + fieldType fieldTarget(targetIO, meshTarget); + + interp.mapSrcToTgt(fieldSource, cop, fieldTarget); + + fieldTarget.write(); + } + else + { + targetIO.readOpt() = IOobject::NO_READ; + + tmp<fieldType> + tfieldTarget(interp.mapSrcToTgt(fieldSource, cop)); + + fieldType fieldTarget(targetIO, tfieldTarget); + + fieldTarget.write(); + } + } + } +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/utilities/preProcessing/mapFieldsNew/UnMapped.H b/applications/utilities/preProcessing/mapFieldsNew/UnMapped.H new file mode 100644 index 0000000000000000000000000000000000000000..aa67e3c666c5bc6494a0284361466220dafd9dfd --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/UnMapped.H @@ -0,0 +1,57 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-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/>. + +\*---------------------------------------------------------------------------*/ + +#ifndef UnMapped_H +#define UnMapped_H + +#include "IOobjectList.H" +#include "OSspecific.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +template<class Type> +void UnMapped(const IOobjectList& objects) +{ + IOobjectList fields = objects.lookupClass(Type::typeName); + + forAllConstIter(IOobjectList, fields, fieldIter) + { + mvBak(fieldIter()->objectPath(), "unmapped"); + } +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/utilities/preProcessing/mapFieldsNew/createTimes.H b/applications/utilities/preProcessing/mapFieldsNew/createTimes.H new file mode 100644 index 0000000000000000000000000000000000000000..f5fd92b607544939d01f0b67f0109c5098c7c07b --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/createTimes.H @@ -0,0 +1,11 @@ + Info<< "\nCreate databases as time" << endl; + + HashTable<string> srcOptions(args.options()); + srcOptions.erase("case"); + srcOptions.insert("case", fileName(rootDirSource/caseDirSource)); + + argList argsSrc(args, srcOptions, false, false, false); + + Time runTimeSource(Time::controlDictName, argsSrc); + + Time runTimeTarget(Time::controlDictName, args); diff --git a/applications/utilities/preProcessing/mapFieldsNew/mapFields.C b/applications/utilities/preProcessing/mapFieldsNew/mapFields.C new file mode 100644 index 0000000000000000000000000000000000000000..2c8757fbf0821da1034b42d1ab0f1f8bb7180e23 --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/mapFields.C @@ -0,0 +1,343 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-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/>. + +Application + mapFields + +Description + Maps volume fields from one mesh to another, reading and + interpolating all fields present in the time directory of both cases. + +\*---------------------------------------------------------------------------*/ + +#include "fvCFD.H" +#include "meshToMeshNew.H" +#include "processorPolyPatch.H" +#include "MapMeshes.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +void mapConsistentMesh +( + const fvMesh& meshSource, + const fvMesh& meshTarget, + const meshToMeshNew::interpolationMethod& mapMethod, + const bool subtract, + const HashSet<word>& selectedFields, + const bool noLagrangian +) +{ + Info<< nl << "Consistently creating and mapping fields for time " + << meshSource.time().timeName() << nl << endl; + + meshToMeshNew interp(meshSource, meshTarget, mapMethod); + + if (subtract) + { + MapMesh<minusEqOp> + ( + interp, + selectedFields, + noLagrangian + ); + } + else + { + MapMesh<plusEqOp> + ( + interp, + selectedFields, + noLagrangian + ); + } +} + + +void mapSubMesh +( + const fvMesh& meshSource, + const fvMesh& meshTarget, + const HashTable<word>& patchMap, + const wordList& cuttingPatches, + const meshToMeshNew::interpolationMethod& mapMethod, + const bool subtract, + const HashSet<word>& selectedFields, + const bool noLagrangian +) +{ + Info<< nl << "Creating and mapping fields for time " + << meshSource.time().timeName() << nl << endl; + + meshToMeshNew interp + ( + meshSource, + meshTarget, + mapMethod, + patchMap, + cuttingPatches + ); + + if (subtract) + { + MapMesh<minusEqOp> + ( + interp, + selectedFields, + noLagrangian + ); + } + else + { + MapMesh<plusEqOp> + ( + interp, + selectedFields, + noLagrangian + ); + } +} + + +wordList addProcessorPatches +( + const fvMesh& meshTarget, + const wordList& cuttingPatches +) +{ + // Add the processor patches to the cutting list + HashSet<word> cuttingPatchTable; + forAll(cuttingPatches, i) + { + cuttingPatchTable.insert(cuttingPatches[i]); + } + + const polyBoundaryMesh& pbm = meshTarget.boundaryMesh(); + + forAll(pbm, patchI) + { + if (isA<processorPolyPatch>(pbm[patchI])) + { + const word& patchName = pbm[patchI].name(); + cuttingPatchTable.insert(patchName); + } + } + + return cuttingPatchTable.toc(); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + argList::addNote + ( + "map volume fields from one mesh to another" + ); + + argList::validArgs.append("sourceCase"); + + argList::addOption + ( + "sourceTime", + "scalar|'latestTime'", + "specify the source time" + ); + argList::addOption + ( + "sourceRegion", + "word", + "specify the source region" + ); + argList::addOption + ( + "targetRegion", + "word", + "specify the target region" + ); + argList::addBoolOption + ( + "consistent", + "source and target geometry and boundary conditions identical" + ); + argList::addOption + ( + "mapMethod", + "word", + "specify the mapping method" + ); + argList::addBoolOption + ( + "subtract", + "subtract mapped source from target" + ); + argList::addOption + ( + "fields", + "list", + "specify a list of fields to be mapped. Eg, '(U T p)' - " + "regular expressions not currently supported" + ); + argList::addBoolOption + ( + "noLagrangian", + "skip mapping lagrangian positions and fields" + ); + + argList args(argc, argv); + + fileName rootDirTarget(args.rootPath()); + fileName caseDirTarget(args.globalCaseName()); + + const fileName casePath = args[1]; + const fileName rootDirSource = casePath.path(); + const fileName caseDirSource = casePath.name(); + + Info<< "Source: " << rootDirSource << " " << caseDirSource << endl; + word sourceRegion = fvMesh::defaultRegion; + if (args.optionFound("sourceRegion")) + { + sourceRegion = args["sourceRegion"]; + Info<< "Source region: " << sourceRegion << endl; + } + + Info<< "Target: " << rootDirTarget << " " << caseDirTarget << endl; + word targetRegion = fvMesh::defaultRegion; + if (args.optionFound("targetRegion")) + { + targetRegion = args["targetRegion"]; + Info<< "Target region: " << targetRegion << endl; + } + + const bool consistent = args.optionFound("consistent"); + + meshToMeshNew::interpolationMethod mapMethod = + meshToMeshNew::imCellVolumeWeight; + + if (args.optionFound("mapMethod")) + { + mapMethod = meshToMeshNew::interpolationMethodNames_[args["mapMethod"]]; + + Info<< "Mapping method: " + << meshToMeshNew::interpolationMethodNames_[mapMethod] << endl; + } + + const bool subtract = args.optionFound("subtract"); + if (subtract) + { + Info<< "Subtracting mapped source field from target" << endl; + } + + HashSet<word> selectedFields; + if (args.optionFound("fields")) + { + args.optionLookup("fields")() >> selectedFields; + } + + const bool noLagrangian = args.optionFound("noLagrangian"); + + #include "createTimes.H" + + HashTable<word> patchMap; + wordList cuttingPatches; + + if (!consistent) + { + IOdictionary mapFieldsDict + ( + IOobject + ( + "mapFieldsDict", + runTimeTarget.system(), + runTimeTarget, + IOobject::MUST_READ_IF_MODIFIED, + IOobject::NO_WRITE, + false + ) + ); + + mapFieldsDict.lookup("patchMap") >> patchMap; + mapFieldsDict.lookup("cuttingPatches") >> cuttingPatches; + } + + #include "setTimeIndex.H" + + Info<< "\nCreate meshes\n" << endl; + + fvMesh meshSource + ( + IOobject + ( + sourceRegion, + runTimeSource.timeName(), + runTimeSource + ) + ); + + fvMesh meshTarget + ( + IOobject + ( + targetRegion, + runTimeTarget.timeName(), + runTimeTarget + ) + ); + + Info<< "Source mesh size: " << meshSource.nCells() << tab + << "Target mesh size: " << meshTarget.nCells() << nl << endl; + + if (consistent) + { + mapConsistentMesh + ( + meshSource, + meshTarget, + mapMethod, + subtract, + selectedFields, + noLagrangian + ); + } + else + { + mapSubMesh + ( + meshSource, + meshTarget, + patchMap, + addProcessorPatches(meshTarget, cuttingPatches), + mapMethod, + subtract, + selectedFields, + noLagrangian + ); + } + + Info<< "\nEnd\n" << endl; + + return 0; +} + + +// ************************************************************************* // diff --git a/applications/utilities/preProcessing/mapFieldsNew/mapFieldsDict b/applications/utilities/preProcessing/mapFieldsNew/mapFieldsDict new file mode 100644 index 0000000000000000000000000000000000000000..00a190cc77e4ee16f157198cee7af4e5ce19e8f1 --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/mapFieldsDict @@ -0,0 +1,30 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: dev | +| \\ / A nd | Web: www.OpenFOAM.org | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object mapFieldsDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +// List of pairs of source/target patches for mapping +patchMap +( + lid movingWall +); + +// List of target patches cutting the source domain (these need to be +// handled specially e.g. interpolated from internal values) +cuttingPatches +( + fixedWalls +); + +// ************************************************************************* // diff --git a/applications/utilities/preProcessing/mapFieldsNew/mapLagrangian.C b/applications/utilities/preProcessing/mapFieldsNew/mapLagrangian.C new file mode 100644 index 0000000000000000000000000000000000000000..c868336e2f9652d39295de7607d59966fd457d8c --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/mapLagrangian.C @@ -0,0 +1,303 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-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 "MapLagrangianFields.H" +#include "passiveParticleCloud.H" +#include "meshSearch.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +static const scalar perturbFactor = 1e-6; + + +// Special version of findCell that generates a cell guaranteed to be +// compatible with tracking. +static label findCell(const Cloud<passiveParticle>& cloud, const point& pt) +{ + label cellI = -1; + label tetFaceI = -1; + label tetPtI = -1; + + const polyMesh& mesh = cloud.pMesh(); + + mesh.findCellFacePt(pt, cellI, tetFaceI, tetPtI); + + if (cellI >= 0) + { + return cellI; + } + else + { + // See if particle on face by finding nearest face and shifting + // particle. + + meshSearch meshSearcher + ( + mesh, + polyMesh::FACEPLANES // no decomposition needed + ); + + label faceI = meshSearcher.findNearestBoundaryFace(pt); + + if (faceI >= 0) + { + const point& cc = mesh.cellCentres()[mesh.faceOwner()[faceI]]; + + const point perturbPt = (1-perturbFactor)*pt+perturbFactor*cc; + + mesh.findCellFacePt(perturbPt, cellI, tetFaceI, tetPtI); + + return cellI; + } + } + + return -1; +} + + +void mapLagrangian(const meshToMeshNew& interp) +{ + // Determine which particles are in meshTarget + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + const polyMesh& meshSource = interp.srcRegion(); + const polyMesh& meshTarget = interp.tgtRegion(); + const labelListList& sourceToTarget = interp.srcToTgtCellAddr(); + + const pointField& targetCc = meshTarget.cellCentres(); + + fileNameList cloudDirs + ( + readDir + ( + meshSource.time().timePath()/cloud::prefix, + fileName::DIRECTORY + ) + ); + + forAll(cloudDirs, cloudI) + { + // Search for list of lagrangian objects for this time + IOobjectList objects + ( + meshSource, + meshSource.time().timeName(), + cloud::prefix/cloudDirs[cloudI] + ); + + IOobject* positionsPtr = objects.lookup(word("positions")); + + if (positionsPtr) + { + Info<< nl << " processing cloud " << cloudDirs[cloudI] << endl; + + // Read positions & cell + passiveParticleCloud sourceParcels + ( + meshSource, + cloudDirs[cloudI], + false + ); + Info<< " read " << sourceParcels.size() + << " parcels from source mesh." << endl; + + // Construct empty target cloud + passiveParticleCloud targetParcels + ( + meshTarget, + cloudDirs[cloudI], + IDLList<passiveParticle>() + ); + + particle::TrackingData<passiveParticleCloud> td(targetParcels); + + label sourceParticleI = 0; + + // Indices of source particles that get added to targetParcels + DynamicList<label> addParticles(sourceParcels.size()); + + // Unmapped particles + labelHashSet unmappedSource(sourceParcels.size()); + + + // Initial: track from fine-mesh cell centre to particle position + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + // This requires there to be no boundary in the way. + + + forAllConstIter(Cloud<passiveParticle>, sourceParcels, iter) + { + bool foundCell = false; + + // Assume that cell from read parcel is the correct one... + if (iter().cell() >= 0) + { + const labelList& targetCells = + sourceToTarget[iter().cell()]; + + // Particle probably in one of the targetcells. Try + // all by tracking from their cell centre to the parcel + // position. + + forAll(targetCells, i) + { + // Track from its cellcentre to position to make sure. + autoPtr<passiveParticle> newPtr + ( + new passiveParticle + ( + meshTarget, + targetCc[targetCells[i]], + targetCells[i] + ) + ); + passiveParticle& newP = newPtr(); + + label faceI = newP.track(iter().position(), td); + + if (faceI < 0 && newP.cell() >= 0) + { + // Hit position. + foundCell = true; + addParticles.append(sourceParticleI); + targetParcels.addParticle(newPtr.ptr()); + break; + } + } + } + + if (!foundCell) + { + // Store for closer analysis + unmappedSource.insert(sourceParticleI); + } + + sourceParticleI++; + } + + Info<< " after meshToMesh addressing found " + << targetParcels.size() + << " parcels in target mesh." << endl; + + + // Do closer inspection for unmapped particles + // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + + if (unmappedSource.size()) + { + sourceParticleI = 0; + + forAllIter(Cloud<passiveParticle>, sourceParcels, iter) + { + if (unmappedSource.found(sourceParticleI)) + { + label targetCell = + findCell(targetParcels, iter().position()); + + if (targetCell >= 0) + { + unmappedSource.erase(sourceParticleI); + addParticles.append(sourceParticleI); + iter().cell() = targetCell; + targetParcels.addParticle + ( + sourceParcels.remove(&iter()) + ); + } + } + sourceParticleI++; + } + } + addParticles.shrink(); + + Info<< " after additional mesh searching found " + << targetParcels.size() << " parcels in target mesh." << endl; + + if (addParticles.size()) + { + IOPosition<passiveParticleCloud>(targetParcels).write(); + + // addParticles now contains the indices of the sourceMesh + // particles that were appended to the target mesh. + + // Map lagrangian fields + // ~~~~~~~~~~~~~~~~~~~~~ + + MapLagrangianFields<label> + ( + cloudDirs[cloudI], + objects, + meshTarget, + addParticles + ); + MapLagrangianFields<scalar> + ( + cloudDirs[cloudI], + objects, + meshTarget, + addParticles + ); + MapLagrangianFields<vector> + ( + cloudDirs[cloudI], + objects, + meshTarget, + addParticles + ); + MapLagrangianFields<sphericalTensor> + ( + cloudDirs[cloudI], + objects, + meshTarget, + addParticles + ); + MapLagrangianFields<symmTensor> + ( + cloudDirs[cloudI], + objects, + meshTarget, + addParticles + ); + MapLagrangianFields<tensor> + ( + cloudDirs[cloudI], + objects, + meshTarget, + addParticles + ); + } + } + } +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// ************************************************************************* // diff --git a/applications/utilities/preProcessing/mapFieldsNew/mapLagrangian.H b/applications/utilities/preProcessing/mapFieldsNew/mapLagrangian.H new file mode 100644 index 0000000000000000000000000000000000000000..47904fbabdfa43704ccb3809482655017888d931 --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/mapLagrangian.H @@ -0,0 +1,56 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-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/>. + +InNamespace + Foam + +Description + Maps lagrangian positions and fields + +SourceFiles + mapLagrangian.C + +\*---------------------------------------------------------------------------*/ + +#ifndef mapLagrangian_H +#define mapLagrangian_H + +#include "meshToMeshNew.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +//- Maps lagrangian positions and fields +void mapLagrangian(const meshToMeshNew& interp); + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/utilities/preProcessing/mapFieldsNew/setTimeIndex.H b/applications/utilities/preProcessing/mapFieldsNew/setTimeIndex.H new file mode 100644 index 0000000000000000000000000000000000000000..a468d50648ebd176c14fbf3716319c6b0ee5e2fb --- /dev/null +++ b/applications/utilities/preProcessing/mapFieldsNew/setTimeIndex.H @@ -0,0 +1,33 @@ +{ + instantList sourceTimes = runTimeSource.times(); + label sourceTimeIndex = runTimeSource.timeIndex(); + if (args.optionFound("sourceTime")) + { + if (args["sourceTime"] == "latestTime") + { + sourceTimeIndex = sourceTimes.size() - 1; + } + else + { + sourceTimeIndex = Time::findClosestTimeIndex + ( + sourceTimes, + args.optionRead<scalar>("sourceTime") + ); + } + } + else + { + sourceTimeIndex = Time::findClosestTimeIndex + ( + sourceTimes, + runTimeTarget.time().value() + ); + } + + runTimeSource.setTime(sourceTimes[sourceTimeIndex], sourceTimeIndex); + + Info<< "\nSource time: " << runTimeSource.value() + << "\nTarget time: " << runTimeTarget.value() + << endl; +}