From 743dea87d257bcc30bba5f45a2efb58de8ad69b1 Mon Sep 17 00:00:00 2001 From: Will Bainbridge Date: Fri, 28 Apr 2017 08:03:44 +0100 Subject: [PATCH 001/535] Lagrangian: Rewrite of the particle tracking algorithm to function in terms of the local barycentric coordinates of the current tetrahedron, rather than the global coordinate system. Barycentric tracking works on any mesh, irrespective of mesh quality. Particles do not get "lost", and tracking does not require ad-hoc "corrections" or "rescues" to function robustly, because the calculation of particle-face intersections is unambiguous and reproducible, even at small angles of incidence. Each particle position is defined by topology (i.e. the decomposed tet cell it is in) and geometry (i.e. where it is in the cell). No search operations are needed on restart or reconstruct, unlike when particle positions are stored in the global coordinate system. The particle positions file now contains particles' local coordinates and topology, rather than the global coordinates and cell. This change to the output format is not backwards compatible. Existing cases with Lagrangian data will not restart, but they will still run from time zero without any modification. This change was necessary in order to guarantee that the loaded particle is valid, and therefore fundamentally prevent "loss" and "search-failure" type bugs (e.g., 2517, 2442, 2286, 1836, 1461, 1341, 1097). The tracking functions have also been converted to function in terms of displacement, rather than end position. This helps remove floating point error issues, particularly towards the end of a tracking step. Wall bounded streamlines have been removed. The implementation proved incompatible with the new tracking algorithm. ParaView has a surface LIC plugin which provides equivalent, or better, functionality. Additionally, bug report is resolved by this change. --- .../DPMFoam/DPMDyMFoam/DPMDyMFoam.C | 3 + .../icoUncoupledKinematicParcelDyMFoam.C | 4 +- .../sprayFoam/sprayDyMFoam/sprayDyMFoam.C | 5 +- .../Make/files | 3 + .../Make/options | 36 + .../uncoupledKinematicParcelDyMFoam.C | 90 ++ .../decomposePar/lagrangianFieldDecomposer.C | 40 +- .../preProcessing/mapFields/mapLagrangian.C | 22 +- .../mapFieldsPar/mapLagrangian.C | 23 +- etc/caseDicts/meshQualityDict | 3 +- src/OpenFOAM/Make/files | 2 + src/OpenFOAM/meshes/polyMesh/polyMesh.C | 51 +- src/OpenFOAM/meshes/polyMesh/polyMesh.H | 13 +- src/OpenFOAM/meshes/polyMesh/polyMeshClear.C | 29 +- src/OpenFOAM/meshes/polyMesh/polyMeshIO.C | 18 +- .../primitives/Barycentric/Barycentric.H | 117 ++ .../primitives/Barycentric/BarycentricI.H | 132 +++ .../Barycentric/BarycentricTensor.H | 139 +++ .../Barycentric/BarycentricTensorI.H | 196 ++++ .../Barycentric/barycentric/barycentric.C | 95 ++ .../Barycentric/barycentric/barycentric.H | 74 ++ .../barycentricTensor/barycentricTensor.H | 64 ++ .../polynomialEqns/cubicEqn/cubicEqn.H | 3 + .../polynomialEqns/cubicEqn/cubicEqnI.H | 6 + .../polynomialEqns/linearEqn/linearEqn.H | 3 + .../polynomialEqns/linearEqn/linearEqnI.H | 6 + .../quadraticEqn/quadraticEqn.H | 3 + .../quadraticEqn/quadraticEqnI.H | 6 + .../motionSmoother/motionSmootherAlgo.C | 2 +- src/functionObjects/field/Make/files | 5 - .../field/nearWallFields/findCellParticle.C | 50 +- .../field/nearWallFields/findCellParticle.H | 42 +- .../field/nearWallFields/nearWallFields.C | 39 +- .../field/streamLine/streamLineParticle.C | 104 +- .../field/streamLine/streamLineParticle.H | 60 +- .../clouds/Templates/DSMCCloud/DSMCCloud.C | 35 +- .../clouds/Templates/DSMCCloud/DSMCCloud.H | 7 +- .../parcels/Templates/DSMCParcel/DSMCParcel.C | 21 +- .../parcels/Templates/DSMCParcel/DSMCParcel.H | 20 +- .../Templates/DSMCParcel/DSMCParcelI.H | 28 +- .../FreeStream/FreeStream.C | 13 +- src/lagrangian/basic/Cloud/Cloud.C | 91 +- src/lagrangian/basic/Cloud/Cloud.H | 14 +- src/lagrangian/basic/Cloud/CloudIO.C | 9 +- .../basic/InteractionLists/InteractionLists.C | 20 +- .../basic/indexedParticle/indexedParticle.H | 19 +- src/lagrangian/basic/particle/particle.C | 998 +++++++++++++++++- src/lagrangian/basic/particle/particle.H | 384 ++++--- src/lagrangian/basic/particle/particleI.H | 753 +------------ src/lagrangian/basic/particle/particleIO.C | 42 +- .../basic/particle/particleTemplates.C | 810 ++------------ .../basic/passiveParticle/passiveParticle.H | 16 +- .../Templates/KinematicCloud/KinematicCloud.C | 12 +- .../Templates/ReactingCloud/ReactingCloud.C | 8 +- .../ReactingMultiphaseCloud.C | 9 +- .../Templates/ThermoCloud/ThermoCloud.C | 8 +- .../CollidingParcel/CollidingParcel.H | 17 +- .../CollidingParcel/CollidingParcelI.H | 24 +- .../KinematicParcel/KinematicParcel.C | 83 +- .../KinematicParcel/KinematicParcel.H | 20 +- .../KinematicParcel/KinematicParcelI.H | 48 +- .../Templates/MPPICParcel/MPPICParcel.H | 17 +- .../Templates/MPPICParcel/MPPICParcelI.H | 21 +- .../ReactingMultiphaseParcel.H | 16 +- .../ReactingMultiphaseParcelI.H | 24 +- .../Templates/ReactingParcel/ReactingParcel.H | 17 +- .../ReactingParcel/ReactingParcelI.H | 23 +- .../Templates/ThermoParcel/ThermoParcel.H | 17 +- .../Templates/ThermoParcel/ThermoParcelI.H | 24 +- .../InjectionModel/InjectionModel.C | 8 +- .../SurfaceFilmModel/SurfaceFilmModel.C | 21 +- .../ThermoSurfaceFilm/ThermoSurfaceFilm.C | 2 +- .../molecule/molecule/molecule.C | 20 +- .../molecule/molecule/molecule.H | 22 +- .../molecule/molecule/moleculeI.H | 41 +- .../molecule/moleculeCloud/moleculeCloud.C | 67 +- .../molecule/moleculeCloud/moleculeCloud.H | 4 +- src/lagrangian/solidParticle/solidParticle.C | 31 +- src/lagrangian/solidParticle/solidParticle.H | 4 +- src/lagrangian/solidParticle/solidParticleI.H | 6 +- .../Templates/SprayParcel/SprayParcel.H | 17 +- .../Templates/SprayParcel/SprayParcelI.H | 37 +- .../meshRefinement/meshRefinementRefine.C | 37 +- .../trackedParticle/trackedParticle.C | 54 +- .../trackedParticle/trackedParticle.H | 27 +- .../reconstructLagrangianPositions.C | 21 +- src/sampling/sampledSet/face/faceOnlySet.C | 8 +- .../sampledSet/polyLine/polyLineSet.C | 23 +- .../sampledSet/sampledSet/sampledSet.C | 12 +- src/sampling/sampledSet/uniform/uniformSet.C | 33 +- .../simpleFoam/motorBike/system/controlDict | 1 - 91 files changed, 3209 insertions(+), 2443 deletions(-) create mode 100644 applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/Make/files create mode 100644 applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/Make/options create mode 100644 applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/uncoupledKinematicParcelDyMFoam.C create mode 100644 src/OpenFOAM/primitives/Barycentric/Barycentric.H create mode 100644 src/OpenFOAM/primitives/Barycentric/BarycentricI.H create mode 100644 src/OpenFOAM/primitives/Barycentric/BarycentricTensor.H create mode 100644 src/OpenFOAM/primitives/Barycentric/BarycentricTensorI.H create mode 100644 src/OpenFOAM/primitives/Barycentric/barycentric/barycentric.C create mode 100644 src/OpenFOAM/primitives/Barycentric/barycentric/barycentric.H create mode 100644 src/OpenFOAM/primitives/Barycentric/barycentricTensor/barycentricTensor.H diff --git a/applications/solvers/lagrangian/DPMFoam/DPMDyMFoam/DPMDyMFoam.C b/applications/solvers/lagrangian/DPMFoam/DPMDyMFoam/DPMDyMFoam.C index aec4c64304..ded2af139b 100644 --- a/applications/solvers/lagrangian/DPMFoam/DPMDyMFoam/DPMDyMFoam.C +++ b/applications/solvers/lagrangian/DPMFoam/DPMDyMFoam/DPMDyMFoam.C @@ -77,6 +77,9 @@ int main(int argc, char *argv[]) Info<< "Time = " << runTime.timeName() << nl << endl; + // Store the particle positions + kinematicCloud.storeGlobalPositions(); + mesh.update(); // Calculate absolute flux from the mapped surface velocity diff --git a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelDyMFoam/icoUncoupledKinematicParcelDyMFoam.C b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelDyMFoam/icoUncoupledKinematicParcelDyMFoam.C index be0905796b..00c66823ea 100644 --- a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelDyMFoam/icoUncoupledKinematicParcelDyMFoam.C +++ b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelDyMFoam/icoUncoupledKinematicParcelDyMFoam.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -68,6 +68,8 @@ int main(int argc, char *argv[]) { Info<< "Time = " << runTime.timeName() << nl << endl; + kinematicCloud.storeGlobalPositions(); + mesh.update(); U.correctBoundaryConditions(); diff --git a/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/sprayDyMFoam.C b/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/sprayDyMFoam.C index 16d2d524c5..40acccae9e 100644 --- a/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/sprayDyMFoam.C +++ b/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/sprayDyMFoam.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2015-2016 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2015-2017 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -91,6 +91,9 @@ int main(int argc, char *argv[]) // Store momentum to set rhoUf for introduced faces. volVectorField rhoU("rhoU", rho*U); + // Store the particle positions + parcels.storeGlobalPositions(); + // Do any mesh changes mesh.update(); diff --git a/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/Make/files b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/Make/files new file mode 100644 index 0000000000..430b1df24f --- /dev/null +++ b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/Make/files @@ -0,0 +1,3 @@ +uncoupledKinematicParcelDyMFoam.C + +EXE = $(FOAM_APPBIN)/uncoupledKinematicParcelDyMFoam diff --git a/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/Make/options b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/Make/options new file mode 100644 index 0000000000..2d385f9853 --- /dev/null +++ b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/Make/options @@ -0,0 +1,36 @@ +EXE_INC = \ + -I.. \ + -I$(LIB_SRC)/lagrangian/basic/lnInclude \ + -I$(LIB_SRC)/lagrangian/intermediate/lnInclude \ + -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \ + -I$(LIB_SRC)/transportModels/compressible/lnInclude \ + -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ + -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \ + -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + -I$(LIB_SRC)/meshTools/lnInclude \ + -I$(LIB_SRC)/regionModels/regionModel/lnInclude \ + -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \ + -I$(LIB_SRC)/dynamicMesh/lnInclude \ + -I$(LIB_SRC)/dynamicFvMesh/lnInclude + +EXE_LIBS = \ + -llagrangian \ + -llagrangianIntermediate \ + -llagrangianTurbulence \ + -lcompressibleTransportModels \ + -lfluidThermophysicalModels \ + -lspecie \ + -lradiationModels \ + -lturbulenceModels \ + -lcompressibleTurbulenceModels \ + -lfiniteVolume \ + -lfvOptions \ + -lmeshTools \ + -lregionModels \ + -lsurfaceFilmModels \ + -ldynamicMesh \ + -ldynamicFvMesh \ + -ltopoChangerFvMesh diff --git a/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/uncoupledKinematicParcelDyMFoam.C b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/uncoupledKinematicParcelDyMFoam.C new file mode 100644 index 0000000000..18f8e63869 --- /dev/null +++ b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/uncoupledKinematicParcelDyMFoam.C @@ -0,0 +1,90 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2017 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 . + +Application + uncoupledKinematicParcelDyMFoam + +Description + Transient solver for the passive transport of a particle cloud. + + Uses a pre- calculated velocity field to evolve the cloud. + +\*---------------------------------------------------------------------------*/ + +#include "fvCFD.H" +#include "dynamicFvMesh.H" +#include "psiThermo.H" +#include "turbulentFluidThermoModel.H" +#include "basicKinematicCloud.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + argList::addOption + ( + "cloudName", + "name", + "specify alternative cloud name. default is 'kinematicCloud'" + ); + + #define NO_CONTROL + #include "postProcess.H" + + #include "setRootCase.H" + #include "createTime.H" + #include "createDynamicFvMesh.H" + #include "createFields.H" + #include "compressibleCourantNo.H" + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + Info<< "\nStarting time loop\n" << endl; + + while (runTime.loop()) + { + Info<< "Time = " << runTime.timeName() << nl << endl; + + kinematicCloud.storeGlobalPositions(); + + mesh.update(); + + U.correctBoundaryConditions(); + + Info<< "Evolving " << kinematicCloud.name() << endl; + kinematicCloud.evolve(); + + runTime.write(); + + Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << " ClockTime = " << runTime.elapsedClockTime() << " s" + << nl << endl; + } + + Info<< "End\n" << endl; + + return 0; +} + + +// ************************************************************************* // diff --git a/applications/utilities/parallelProcessing/decomposePar/lagrangianFieldDecomposer.C b/applications/utilities/parallelProcessing/decomposePar/lagrangianFieldDecomposer.C index efb7dd2a23..98457dc2b4 100644 --- a/applications/utilities/parallelProcessing/decomposePar/lagrangianFieldDecomposer.C +++ b/applications/utilities/parallelProcessing/decomposePar/lagrangianFieldDecomposer.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -47,13 +47,12 @@ Foam::lagrangianFieldDecomposer::lagrangianFieldDecomposer { label pi = 0; - // faceProcAddressing not required currently - // labelList decodedProcFaceAddressing(faceProcAddressing.size()); + labelList decodedProcFaceAddressing(faceProcAddressing.size()); - // forAll(faceProcAddressing, i) - // { - // decodedProcFaceAddressing[i] = mag(faceProcAddressing[i]) - 1; - // } + forAll(faceProcAddressing, i) + { + decodedProcFaceAddressing[i] = mag(faceProcAddressing[i]) - 1; + } forAll(cellProcAddressing, procCelli) { @@ -68,27 +67,28 @@ Foam::lagrangianFieldDecomposer::lagrangianFieldDecomposer const indexedParticle& ppi = *iter(); particleIndices_[pi++] = ppi.index(); - // label mappedTetFace = findIndex - // ( - // decodedProcFaceAddressing, - // ppi.tetFace() - // ); + label mappedTetFace = findIndex + ( + decodedProcFaceAddressing, + ppi.tetFace() + ); - // if (mappedTetFace == -1) - // { - // FatalErrorInFunction - // << "Face lookup failure." << nl - // << abort(FatalError); - // } + if (mappedTetFace == -1) + { + FatalErrorInFunction + << "Face lookup failure." << nl + << abort(FatalError); + } positions_.append ( new passiveParticle ( procMesh, - ppi.position(), + ppi.coordinates(), procCelli, - false + mappedTetFace, + ppi.procTetPt(procMesh, procCelli, mappedTetFace) ) ); } diff --git a/applications/utilities/preProcessing/mapFields/mapLagrangian.C b/applications/utilities/preProcessing/mapFields/mapLagrangian.C index ad80ba3d07..45cd8e4bb5 100644 --- a/applications/utilities/preProcessing/mapFields/mapLagrangian.C +++ b/applications/utilities/preProcessing/mapFields/mapLagrangian.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -98,7 +98,6 @@ void mapLagrangian(const meshToMesh0& meshToMesh0Interp) const fvMesh& meshSource = meshToMesh0Interp.fromMesh(); const fvMesh& meshTarget = meshToMesh0Interp.toMesh(); - const pointField& targetCc = meshTarget.cellCentres(); fileNameList cloudDirs @@ -182,15 +181,17 @@ void mapLagrangian(const meshToMesh0& meshToMesh0Interp) new passiveParticle ( meshTarget, - targetCc[targetCells[i]], - targetCells[i] + barycentric(1, 0, 0, 0), + targetCells[i], + meshTarget.cells()[targetCells[i]][0], + 1 ) ); passiveParticle& newP = newPtr(); - label facei = newP.track(iter().position(), td); + newP.track(iter().position() - newP.position(), 0); - if (facei < 0 && newP.cell() >= 0) + if (!newP.onFace()) { // Hit position. foundCell = true; @@ -233,11 +234,16 @@ void mapLagrangian(const meshToMesh0& meshToMesh0Interp) { unmappedSource.erase(sourceParticleI); addParticles.append(sourceParticleI); - iter().cell() = targetCell; targetParcels.addParticle ( - sourceParcels.remove(&iter()) + new passiveParticle + ( + meshTarget, + iter().position(), + targetCell + ) ); + sourceParcels.remove(&iter()); } } sourceParticleI++; diff --git a/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C b/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C index 06f7edc910..2f41a6eee6 100644 --- a/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C +++ b/applications/utilities/preProcessing/mapFieldsPar/mapLagrangian.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -89,8 +89,6 @@ void mapLagrangian(const meshToMesh& interp) const polyMesh& meshTarget = interp.tgtRegion(); const labelListList& sourceToTarget = interp.srcToTgtCellAddr(); - const pointField& targetCc = meshTarget.cellCentres(); - fileNameList cloudDirs ( readDir @@ -173,15 +171,17 @@ void mapLagrangian(const meshToMesh& interp) new passiveParticle ( meshTarget, - targetCc[targetCells[i]], - targetCells[i] + barycentric(1, 0, 0, 0), + targetCells[i], + meshTarget.cells()[targetCells[i]][0], + 1 ) ); passiveParticle& newP = newPtr(); - label facei = newP.track(iter().position(), td); + newP.track(iter().position() - newP.position(), 0); - if (facei < 0 && newP.cell() >= 0) + if (!newP.onFace()) { // Hit position. foundCell = true; @@ -224,11 +224,16 @@ void mapLagrangian(const meshToMesh& interp) { unmappedSource.erase(sourceParticleI); addParticles.append(sourceParticleI); - iter().cell() = targetCell; targetParcels.addParticle ( - sourceParcels.remove(&iter()) + new passiveParticle + ( + meshTarget, + iter().position(), + targetCell + ) ); + sourceParcels.remove(&iter()); } } sourceParticleI++; diff --git a/etc/caseDicts/meshQualityDict b/etc/caseDicts/meshQualityDict index e69207563a..fe5a1f321a 100644 --- a/etc/caseDicts/meshQualityDict +++ b/etc/caseDicts/meshQualityDict @@ -33,8 +33,7 @@ minVol 1e-13; //- Minimum quality of the tet formed by the face-centre // and variable base point minimum decomposition triangles and -// the cell centre. This has to be a positive number for tracking -// to work. Set to very negative number (e.g. -1E30) to +// the cell centre. Set to very negative number (e.g. -1E30) to // disable. // <0 = inside out tet, // 0 = flat tet diff --git a/src/OpenFOAM/Make/files b/src/OpenFOAM/Make/files index d9f2104c18..f79161fa78 100644 --- a/src/OpenFOAM/Make/files +++ b/src/OpenFOAM/Make/files @@ -133,6 +133,8 @@ $(spatialVectorAlgebra)/CompactSpatialTensor/compactSpatialTensor/compactSpatial primitives/polynomialEqns/cubicEqn/cubicEqn.C primitives/polynomialEqns/quadraticEqn/quadraticEqn.C +primitives/Barycentric/barycentric/barycentric.C + containers/HashTables/HashTable/HashTableCore.C containers/HashTables/StaticHashTable/StaticHashTableCore.C containers/Lists/SortableList/ParSortableListName.C diff --git a/src/OpenFOAM/meshes/polyMesh/polyMesh.C b/src/OpenFOAM/meshes/polyMesh/polyMesh.C index 8437bcde5c..52b8c4c54c 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyMesh.C +++ b/src/OpenFOAM/meshes/polyMesh/polyMesh.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation \\/ M anipulation | Copyright (C) 2016 OpenCFD Ltd. ------------------------------------------------------------------------------- License @@ -135,6 +135,29 @@ void Foam::polyMesh::calcDirections() const } +Foam::autoPtr Foam::polyMesh::readTetBasePtIs() const +{ + IOobject io + ( + "tetBasePtIs", + instance(), + meshSubDir, + *this, + IOobject::READ_IF_PRESENT, + IOobject::NO_WRITE + ); + + if (io.headerOk()) + { + return autoPtr(new labelIOList(io)); + } + else + { + return autoPtr(nullptr); + } +} + + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // Foam::polyMesh::polyMesh(const IOobject& io) @@ -207,7 +230,7 @@ Foam::polyMesh::polyMesh(const IOobject& io) comm_(UPstream::worldComm), geometricD_(Zero), solutionD_(Zero), - tetBasePtIsPtr_(nullptr), + tetBasePtIsPtr_(readTetBasePtIs()), cellTreePtr_(nullptr), pointZones_ ( @@ -401,7 +424,7 @@ Foam::polyMesh::polyMesh comm_(UPstream::worldComm), geometricD_(Zero), solutionD_(Zero), - tetBasePtIsPtr_(nullptr), + tetBasePtIsPtr_(readTetBasePtIs()), cellTreePtr_(nullptr), pointZones_ ( @@ -552,7 +575,7 @@ Foam::polyMesh::polyMesh comm_(UPstream::worldComm), geometricD_(Zero), solutionD_(Zero), - tetBasePtIsPtr_(nullptr), + tetBasePtIsPtr_(readTetBasePtIs()), cellTreePtr_(nullptr), pointZones_ ( @@ -815,7 +838,7 @@ Foam::label Foam::polyMesh::nSolutionD() const } -const Foam::labelList& Foam::polyMesh::tetBasePtIs() const +const Foam::labelIOList& Foam::polyMesh::tetBasePtIs() const { if (tetBasePtIsPtr_.empty()) { @@ -828,8 +851,17 @@ const Foam::labelList& Foam::polyMesh::tetBasePtIs() const tetBasePtIsPtr_.reset ( - new labelList + new labelIOList ( + IOobject + ( + "tetBasePtIs", + instance(), + meshSubDir, + *this, + IOobject::READ_IF_PRESENT, + IOobject::NO_WRITE + ), polyMeshTetDecomposition::findFaceBasePts(*this) ) ); @@ -1105,6 +1137,13 @@ Foam::tmp Foam::polyMesh::movePoints points_.instance() = time().timeName(); points_.eventNo() = getEvent(); + if (tetBasePtIsPtr_.valid()) + { + tetBasePtIsPtr_().writeOpt() = IOobject::AUTO_WRITE; + tetBasePtIsPtr_().instance() = time().timeName(); + tetBasePtIsPtr_().eventNo() = getEvent(); + } + tmp sweptVols = primitiveMesh::movePoints ( points_, diff --git a/src/OpenFOAM/meshes/polyMesh/polyMesh.H b/src/OpenFOAM/meshes/polyMesh/polyMesh.H index 1f61f73828..00ff5aa2f9 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyMesh.H +++ b/src/OpenFOAM/meshes/polyMesh/polyMesh.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -150,7 +150,7 @@ private: mutable Vector