From 15c433b3379703f0aa4700f7968d0b5f7df942af Mon Sep 17 00:00:00 2001 From: mattijs Date: Tue, 4 Jun 2013 16:50:33 +0100 Subject: [PATCH 01/12] ENH: wmake: 'wmake dep' to build dependencies and softlinks --- Allwmake | 6 +++--- applications/Allwmake | 4 ++-- src/Allwmake | 2 +- src/OSspecific/POSIX/Allwmake | 3 ++- src/postProcessing/Allwmake | 2 +- wmake/wmake | 12 ++++++++---- 6 files changed, 17 insertions(+), 12 deletions(-) diff --git a/Allwmake b/Allwmake index 16edfd1f6f..9187532553 100755 --- a/Allwmake +++ b/Allwmake @@ -26,12 +26,12 @@ else fi # build OpenFOAM libraries and applications -src/Allwmake -applications/Allwmake +src/Allwmake $* +applications/Allwmake $* if [ "$1" = doc ] then - doc/Allwmake + doc/Allwmake $* fi # ----------------------------------------------------------------- end-of-file diff --git a/applications/Allwmake b/applications/Allwmake index 00bc326641..7f7ac6718d 100755 --- a/applications/Allwmake +++ b/applications/Allwmake @@ -16,7 +16,7 @@ wmakeCheckPwd "$WM_PROJECT_DIR/applications" || { set -x -wmake all utilities -wmake all solvers +wmake all utilities $* +wmake all solvers $* # ----------------------------------------------------------------- end-of-file diff --git a/src/Allwmake b/src/Allwmake index 55b8981ce2..d5177cfafc 100755 --- a/src/Allwmake +++ b/src/Allwmake @@ -24,7 +24,7 @@ wmakeLnInclude OpenFOAM wmakeLnInclude OSspecific/${WM_OSTYPE:-POSIX} Pstream/Allwmake $* -OSspecific/${WM_OSTYPE:-POSIX}/Allwmake +OSspecific/${WM_OSTYPE:-POSIX}/Allwmake $* wmake $makeType OpenFOAM wmake $makeType fileFormats diff --git a/src/OSspecific/POSIX/Allwmake b/src/OSspecific/POSIX/Allwmake index c68f089436..52a54646f0 100755 --- a/src/OSspecific/POSIX/Allwmake +++ b/src/OSspecific/POSIX/Allwmake @@ -1,5 +1,6 @@ #!/bin/sh cd ${0%/*} || exit 1 # run from this directory +makeType=${1:-libo} unset COMP_FLAGS LINK_FLAGS @@ -19,6 +20,6 @@ fi # make (non-shared) object -wmake libo +wmake $makeType # ----------------------------------------------------------------- end-of-file diff --git a/src/postProcessing/Allwmake b/src/postProcessing/Allwmake index 24b764c420..3522721c4d 100755 --- a/src/postProcessing/Allwmake +++ b/src/postProcessing/Allwmake @@ -3,7 +3,7 @@ cd ${0%/*} || exit 1 # run from this directory makeType=${1:-libso} set -x -wmake libo postCalc +wmake ${1:-libo} postCalc wmake $makeType foamCalcFunctions functionObjects/Allwmake $* diff --git a/wmake/wmake b/wmake/wmake index ebcd38cb6e..f6de0bd0c4 100755 --- a/wmake/wmake +++ b/wmake/wmake @@ -56,6 +56,7 @@ or a special target: libo build statically linked lib (.o) libso build dynamically linked lib (.so) jar build Java jar + dep build lnInclude and dependencies only USAGE exit 1 @@ -244,7 +245,7 @@ OBJECTS_DIR=$MakeDir/$WM_OPTIONS touch $OBJECTS_DIR/dontIncludeDeps case "$makeType" in -lib | libo | libso ) +lib | libo | libso | dep ) $make -s -f $WM_DIR/Makefile MAKE_DIR=$MakeDir INCLUDE_DEPS=$OBJECTS_DIR/dontIncludeDeps lnInclude/uptodate ;; esac @@ -258,8 +259,11 @@ rc=$? # make the object files and link #------------------------------------------------------------------------------ -cmd="$make -f $WM_DIR/Makefile MAKE_DIR=$MakeDir INCLUDE_DEPS=$OBJECTS_DIR/includeDeps $makeType" -# echo "cmd=$cmd" -exec $cmd +if [ "$makeType" != dep ] +then + cmd="$make -f $WM_DIR/Makefile MAKE_DIR=$MakeDir INCLUDE_DEPS=$OBJECTS_DIR/includeDeps $makeType" + # echo "cmd=$cmd" + exec $cmd +fi #------------------------------------------------------------------------------ From 54b212e769369ad1a0748a6c6936d6fafbadf004 Mon Sep 17 00:00:00 2001 From: mattijs Date: Tue, 11 Jun 2013 15:42:29 +0100 Subject: [PATCH 02/12] ENH: triSurface: improved error message --- src/triSurface/triSurface/triSurface.C | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/triSurface/triSurface/triSurface.C b/src/triSurface/triSurface/triSurface.C index e0e6c265ff..258958ca69 100644 --- a/src/triSurface/triSurface/triSurface.C +++ b/src/triSurface/triSurface/triSurface.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 @@ -494,6 +494,7 @@ void Foam::triSurface::write ( "triSurface::write(const fileName&, const word&, const bool)" ) << "unknown file extension " << ext + << " for file " << name << ". Supported extensions are '.ftr', '.stl', '.stlb', " << "'.gts', '.obj', '.vtk'" << ", '.off', '.dx', '.smesh', '.ac' and '.tri'" From dee88bceddb8c23f28b122d9326e0dbe1d8052ea Mon Sep 17 00:00:00 2001 From: mattijs Date: Wed, 12 Jun 2013 10:33:36 +0100 Subject: [PATCH 03/12] BUG: displacementLayeredMotionSolver: too early evaluation of boundary conditions --- .../displacementInterpolationMotionSolver.H | 5 +- .../displacementLayeredMotionMotionSolver.C | 104 +++++++++--------- .../displacementLayeredMotionMotionSolver.H | 10 +- 3 files changed, 57 insertions(+), 62 deletions(-) diff --git a/src/fvMotionSolver/fvMotionSolvers/displacement/interpolation/displacementInterpolationMotionSolver.H b/src/fvMotionSolver/fvMotionSolvers/displacement/interpolation/displacementInterpolationMotionSolver.H index ec49b1e223..b5da0b5b3b 100644 --- a/src/fvMotionSolver/fvMotionSolvers/displacement/interpolation/displacementInterpolationMotionSolver.H +++ b/src/fvMotionSolver/fvMotionSolvers/displacement/interpolation/displacementInterpolationMotionSolver.H @@ -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 @@ -36,9 +36,6 @@ Description Tables are in the \a constant/tables directory. -Note - could be a motionSolver - does not use any fvMesh structure. - SourceFiles displacementInterpolationMotionSolver.C diff --git a/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.C b/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.C index f491ee586d..ace954062b 100644 --- a/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.C +++ b/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.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 @@ -253,6 +253,16 @@ Foam::displacementLayeredMotionMotionSolver::faceZoneEvaluate // Only on boundary faces - follow boundary conditions fld = vectorField(pointDisplacement_, meshPoints); } + else if (type == "uniformFollow") + { + // Reads name of name of patch. Then get average point dislacement on + // patch. That becomes the value of fld. + const word patchName(dict.lookup("patch")); + label patchID = mesh().boundaryMesh().findPatchID(patchName); + pointField pdf = + pointDisplacement_.boundaryField()[patchID].patchInternalField(); + fld = gAverage(pdf); + } else { FatalIOErrorIn @@ -399,22 +409,6 @@ Info<< "For cellZone:" << cellZoneI // Implement real bc. patchDisp[patchI].correctBoundaryConditions(); - -//Info<< "Writing displacement for faceZone " << fz.name() -// << " to " << patchDisp[patchI].name() << endl; -//patchDisp[patchI].write(); - -// // Copy into pointDisplacement for other fields to use -// forAll(isZonePoint, pointI) -// { -// if (isZonePoint[pointI]) -// { -// pointDisplacement_[pointI] = patchDisp[patchI][pointI]; -// } -// } -// pointDisplacement_.correctBoundaryConditions(); - - patchI++; } @@ -423,37 +417,40 @@ Info<< "For cellZone:" << cellZoneI // ~~~~~ // solving the interior is just interpolating -// // Get normalised distance -// pointScalarField distance -// ( -// IOobject -// ( -// "distance", -// mesh().time().timeName(), -// mesh(), -// IOobject::NO_READ, -// IOobject::NO_WRITE, -// false -// ), -// pointMesh::New(mesh()), -// dimensionedScalar("distance", dimLength, 0.0) -// ); -// forAll(distance, pointI) -// { -// if (isZonePoint[pointI]) -// { -// scalar d1 = patchDist[0][pointI]; -// scalar d2 = patchDist[1][pointI]; -// if (d1+d2 > SMALL) -// { -// scalar s = d1/(d1+d2); -// distance[pointI] = s; -// } -// } -// } -// Info<< "Writing distance pointScalarField to " << mesh().time().timeName() -// << endl; -// distance.write(); + if (debug) + { + // Get normalised distance + pointScalarField distance + ( + IOobject + ( + "distance", + mesh().time().timeName(), + mesh(), + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + pointMesh::New(mesh()), + dimensionedScalar("distance", dimLength, 0.0) + ); + forAll(distance, pointI) + { + if (isZonePoint[pointI]) + { + scalar d1 = patchDist[0][pointI]; + scalar d2 = patchDist[1][pointI]; + if (d1+d2 > SMALL) + { + scalar s = d1/(d1+d2); + distance[pointI] = s; + } + } + } + Info<< "Writing distance pointScalarField to " + << mesh().time().timeName() << endl; + distance.write(); + } // Average forAll(pointDisplacement_, pointI) @@ -470,7 +467,6 @@ Info<< "For cellZone:" << cellZoneI + s*patchDisp[1][pointI]; } } - pointDisplacement_.correctBoundaryConditions(); } @@ -484,9 +480,7 @@ displacementLayeredMotionMotionSolver ) : displacementMotionSolver(mesh, dict, typeName) -{ - pointDisplacement_.correctBoundaryConditions(); -} +{} // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // @@ -518,6 +512,9 @@ void Foam::displacementLayeredMotionMotionSolver::solve() // the motionSolver accordingly movePoints(mesh().points()); + // Apply boundary conditions + pointDisplacement_.boundaryField().updateCoeffs(); + // Apply all regions (=cellZones) const dictionary& regionDicts = coeffDict().subDict("regions"); @@ -544,6 +541,9 @@ void Foam::displacementLayeredMotionMotionSolver::solve() cellZoneSolve(zoneI, regionDict); } + + // Update pointDisplacement for solved values + pointDisplacement_.correctBoundaryConditions(); } diff --git a/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.H b/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.H index 0f29c4f88f..1aa7884f52 100644 --- a/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.H +++ b/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.H @@ -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 @@ -28,11 +28,6 @@ Description Mesh motion solver for an (multi-block) extruded fvMesh. Gets given the structure of the mesh blocks and boundary conditions on these blocks. - Note: should not be an fvMotionSolver but just a motionSolver. Only here - so we can reuse displacementFvMotionSolver functionality (e.g. surface - following boundary conditions) - - The displacementLayeredMotionCoeffs subdict of dynamicMeshDict specifies per region (=cellZone) the boundary conditions on two opposing patches (=faceZones). It then interpolates the boundary values using topological @@ -44,6 +39,9 @@ Description Use this for faceZones on boundary faces (so it uses the proper boundary conditions on the pointDisplacement). + uniformFollow: like 'follow' but takes the average value of + a specified 'patch' (which is not necessarily colocated) + fixedValue: fixed value. timeVaryingUniformFixedValue: table-driven fixed value. From e9932a0be73ac2ecf8cf1cb73113f0d5112215eb Mon Sep 17 00:00:00 2001 From: Henry Date: Wed, 12 Jun 2013 12:51:56 +0100 Subject: [PATCH 04/12] probes: write probe data in user-time rather than s --- src/sampling/probes/patchProbesTemplates.C | 8 ++++++-- src/sampling/probes/probesTemplates.C | 14 +++++++------- 2 files changed, 13 insertions(+), 9 deletions(-) diff --git a/src/sampling/probes/patchProbesTemplates.C b/src/sampling/probes/patchProbesTemplates.C index 8f37e2f07c..7e8ddd74b4 100644 --- a/src/sampling/probes/patchProbesTemplates.C +++ b/src/sampling/probes/patchProbesTemplates.C @@ -43,7 +43,9 @@ void Foam::patchProbes::sampleAndWrite unsigned int w = IOstream::defaultPrecision() + 7; OFstream& probeStream = *probeFilePtrs_[vField.name()]; - probeStream << setw(w) << vField.time().value(); + probeStream + << setw(w) + << vField.time().timeToUserTime(vField.time().value()); forAll(values, probeI) { @@ -67,7 +69,9 @@ void Foam::patchProbes::sampleAndWrite unsigned int w = IOstream::defaultPrecision() + 7; OFstream& probeStream = *probeFilePtrs_[sField.name()]; - probeStream << setw(w) << sField.time().value(); + probeStream + << setw(w) + << sField.time().timeToUserTime(sField.time().value()); forAll(values, probeI) { diff --git a/src/sampling/probes/probesTemplates.C b/src/sampling/probes/probesTemplates.C index 3a34c1c264..b64bbb577f 100644 --- a/src/sampling/probes/probesTemplates.C +++ b/src/sampling/probes/probesTemplates.C @@ -76,7 +76,7 @@ void Foam::probes::sampleAndWrite unsigned int w = IOstream::defaultPrecision() + 7; OFstream& os = *probeFilePtrs_[vField.name()]; - os << setw(w) << vField.time().value(); + os << setw(w) << vField.time().timeToUserTime(vField.time().value()); forAll(values, probeI) { @@ -90,17 +90,17 @@ void Foam::probes::sampleAndWrite template void Foam::probes::sampleAndWrite ( - const GeometricField& vField + const GeometricField& sField ) { - Field values(sample(vField)); + Field values(sample(sField)); if (Pstream::master()) { unsigned int w = IOstream::defaultPrecision() + 7; - OFstream& os = *probeFilePtrs_[vField.name()]; + OFstream& os = *probeFilePtrs_[sField.name()]; - os << setw(w) << vField.time().value(); + os << setw(w) << sField.time().timeToUserTime(sField.time().value()); forAll(values, probeI) { @@ -259,7 +259,7 @@ template Foam::tmp > Foam::probes::sample ( - const GeometricField& vField + const GeometricField& sField ) const { const Type unsetVal(-VGREAT*pTraits::one); @@ -275,7 +275,7 @@ Foam::probes::sample { if (faceList_[probeI] >= 0) { - values[probeI] = vField[faceList_[probeI]]; + values[probeI] = sField[faceList_[probeI]]; } } From 43ac6073864eac444310a2e52a31d35c4ed1b2c8 Mon Sep 17 00:00:00 2001 From: Henry Date: Wed, 12 Jun 2013 12:53:47 +0100 Subject: [PATCH 05/12] GAMGSolverInterpolate: Correct interpolate for conservation such that the sum of the values in the fines cells is equal to the value in the coarse cell --- .../lduMatrix/solvers/GAMG/GAMGSolver.H | 3 ++- .../solvers/GAMG/GAMGSolverInterpolate.C | 27 ++++++++++++++++++- .../lduMatrix/solvers/GAMG/GAMGSolverSolve.C | 8 +++--- 3 files changed, 33 insertions(+), 5 deletions(-) diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.H b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.H index ed49d00827..7f462840f5 100644 --- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.H +++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.H @@ -230,7 +230,8 @@ class GAMGSolver const lduMatrix& m, const FieldField& interfaceBouCoeffs, const lduInterfaceFieldPtrsList& interfaces, - const scalarField& source, + const labelList& restrictAddressing, + const scalarField& psiC, const direction cmpt ) const; diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverInterpolate.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverInterpolate.C index a24612645b..6c18b43d0a 100644 --- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverInterpolate.C +++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverInterpolate.C @@ -34,7 +34,8 @@ void Foam::GAMGSolver::interpolate const lduMatrix& m, const FieldField& interfaceBouCoeffs, const lduInterfaceFieldPtrsList& interfaces, - const scalarField& source, + const labelList& restrictAddressing, + const scalarField& psiC, const direction cmpt ) const { @@ -80,6 +81,30 @@ void Foam::GAMGSolver::interpolate { psiPtr[celli] = -ApsiPtr[celli]/(diagPtr[celli]); } + + register const label nCCells = psiC.size(); + scalarField corrC(nCCells, 0); + scalarField diagC(nCCells, 0); + + for (register label celli=0; celli(ACf.operator const scalarField&()); - if (interpolateCorrection_) + if (interpolateCorrection_ && leveli < coarsestLevel - 2) { interpolate ( @@ -322,7 +322,8 @@ void Foam::GAMGSolver::Vcycle matrixLevels_[leveli], interfaceLevelsBouCoeffs_[leveli], interfaceLevels_[leveli], - coarseSources[leveli], + agglomeration_.restrictAddressing(leveli + 1), + coarseCorrFields[leveli + 1], cmpt ); } @@ -382,7 +383,8 @@ void Foam::GAMGSolver::Vcycle matrix_, interfaceBouCoeffs_, interfaces_, - finestResidual, + agglomeration_.restrictAddressing(0), + coarseCorrFields[0], cmpt ); } From 6b054eab14026feac01dbee74497411f0596bf6c Mon Sep 17 00:00:00 2001 From: Henry Date: Wed, 12 Jun 2013 12:55:37 +0100 Subject: [PATCH 06/12] MGridGenGAMGAgglomerate: Use sum magSf rather than mag(sum Sf) for face-weighting This is more consistent with the faceArea agglomeration and generates more reliable agglomerations. --- .../MGridGenGamgAgglomeration/MGridGenGAMGAgglomerate.C | 5 ----- .../MGridGenGamgAgglomeration/MGridGenGAMGAgglomeration.C | 6 +++++- 2 files changed, 5 insertions(+), 6 deletions(-) diff --git a/src/fvAgglomerationMethods/MGridGenGamgAgglomeration/MGridGenGAMGAgglomerate.C b/src/fvAgglomerationMethods/MGridGenGamgAgglomeration/MGridGenGAMGAgglomerate.C index e3663315aa..4d32bfa89f 100644 --- a/src/fvAgglomerationMethods/MGridGenGamgAgglomeration/MGridGenGAMGAgglomerate.C +++ b/src/fvAgglomerationMethods/MGridGenGamgAgglomeration/MGridGenGAMGAgglomerate.C @@ -30,11 +30,6 @@ Description #include "fvMesh.H" #include "syncTools.H" -//extern "C" -//{ -//# include "mgridgen.h" -//} - // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // void Foam::MGridGenGAMGAgglomeration:: diff --git a/src/fvAgglomerationMethods/MGridGenGamgAgglomeration/MGridGenGAMGAgglomeration.C b/src/fvAgglomerationMethods/MGridGenGamgAgglomeration/MGridGenGAMGAgglomeration.C index 59c45bc11b..cc259c5038 100644 --- a/src/fvAgglomerationMethods/MGridGenGamgAgglomeration/MGridGenGAMGAgglomeration.C +++ b/src/fvAgglomerationMethods/MGridGenGamgAgglomeration/MGridGenGAMGAgglomeration.C @@ -59,7 +59,11 @@ Foam::MGridGenGAMGAgglomeration::MGridGenGAMGAgglomeration // Start geometric agglomeration from the cell volumes and areas of the mesh scalarField* VPtr = const_cast(&fvMesh_.cellVolumes()); - SubField Sf(fvMesh_.faceAreas(), fvMesh_.nInternalFaces()); + + vectorField magFaceAreas(vector::one*mag(fvMesh_.faceAreas())); + SubField Sf(magFaceAreas, fvMesh_.nInternalFaces()); + //SubField Sf(fvMesh_.faceAreas(), fvMesh_.nInternalFaces()); + vectorField* SfPtr = const_cast ( &Sf.operator const vectorField&() From 818f6a549fefcbc70ce512143b43ea249668504c Mon Sep 17 00:00:00 2001 From: Henry Date: Wed, 12 Jun 2013 12:57:50 +0100 Subject: [PATCH 07/12] pairGAMGAgglomeration: Support reverse cell-scanning without reversing the order of the generated agglomeration Gives better cache coherency --- .../pairGAMGAgglomerate.C | 355 +++++++++--------- .../pairGAMGAgglomeration.C | 5 +- .../pairGAMGAgglomeration.H | 4 +- 3 files changed, 190 insertions(+), 174 deletions(-) diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/pairGAMGAgglomeration/pairGAMGAgglomerate.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/pairGAMGAgglomeration/pairGAMGAgglomerate.C index b0f5e757da..508c0c3f94 100644 --- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/pairGAMGAgglomeration/pairGAMGAgglomerate.C +++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/pairGAMGAgglomeration/pairGAMGAgglomerate.C @@ -26,177 +26,7 @@ License #include "pairGAMGAgglomeration.H" #include "lduAddressing.H" -// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // - -Foam::tmp Foam::pairGAMGAgglomeration::agglomerate -( - label& nCoarseCells, - const lduAddressing& fineMatrixAddressing, - const scalarField& faceWeights -) -{ - const label nFineCells = fineMatrixAddressing.size(); - - const labelUList& upperAddr = fineMatrixAddressing.upperAddr(); - const labelUList& lowerAddr = fineMatrixAddressing.lowerAddr(); - - // For each cell calculate faces - labelList cellFaces(upperAddr.size() + lowerAddr.size()); - labelList cellFaceOffsets(nFineCells + 1); - - // memory management - { - labelList nNbrs(nFineCells, 0); - - forAll(upperAddr, facei) - { - nNbrs[upperAddr[facei]]++; - } - - forAll(lowerAddr, facei) - { - nNbrs[lowerAddr[facei]]++; - } - - cellFaceOffsets[0] = 0; - forAll(nNbrs, celli) - { - cellFaceOffsets[celli+1] = cellFaceOffsets[celli] + nNbrs[celli]; - } - - // reset the whole list to use as counter - nNbrs = 0; - - forAll(upperAddr, facei) - { - cellFaces - [ - cellFaceOffsets[upperAddr[facei]] + nNbrs[upperAddr[facei]] - ] = facei; - - nNbrs[upperAddr[facei]]++; - } - - forAll(lowerAddr, facei) - { - cellFaces - [ - cellFaceOffsets[lowerAddr[facei]] + nNbrs[lowerAddr[facei]] - ] = facei; - - nNbrs[lowerAddr[facei]]++; - } - } - - - // go through the faces and create clusters - - tmp tcoarseCellMap(new labelField(nFineCells, -1)); - labelField& coarseCellMap = tcoarseCellMap(); - - nCoarseCells = 0; - - for (label celli=0; celli maxFaceWeight - ) - { - // Match found. Pick up all the necessary data - matchFaceNo = facei; - maxFaceWeight = faceWeights[facei]; - } - } - - if (matchFaceNo >= 0) - { - // Make a new group - coarseCellMap[upperAddr[matchFaceNo]] = nCoarseCells; - coarseCellMap[lowerAddr[matchFaceNo]] = nCoarseCells; - nCoarseCells++; - } - else - { - // No match. Find the best neighbouring cluster and - // put the cell there - label clusterMatchFaceNo = -1; - scalar clusterMaxFaceCoeff = -GREAT; - - for - ( - label faceOs=cellFaceOffsets[celli]; - faceOs clusterMaxFaceCoeff) - { - clusterMatchFaceNo = facei; - clusterMaxFaceCoeff = faceWeights[facei]; - } - } - - if (clusterMatchFaceNo >= 0) - { - // Add the cell to the best cluster - coarseCellMap[celli] = max - ( - coarseCellMap[upperAddr[clusterMatchFaceNo]], - coarseCellMap[lowerAddr[clusterMatchFaceNo]] - ); - } - } - } - } - - - // Check that all cells are part of clusters, - // if not create single-cell "clusters" for each - for (label celli=0; celli Foam::pairGAMGAgglomeration::agglomerate +( + label& nCoarseCells, + const lduAddressing& fineMatrixAddressing, + const scalarField& faceWeights +) +{ + const label nFineCells = fineMatrixAddressing.size(); + + const labelUList& upperAddr = fineMatrixAddressing.upperAddr(); + const labelUList& lowerAddr = fineMatrixAddressing.lowerAddr(); + + // For each cell calculate faces + labelList cellFaces(upperAddr.size() + lowerAddr.size()); + labelList cellFaceOffsets(nFineCells + 1); + + // memory management + { + labelList nNbrs(nFineCells, 0); + + forAll(upperAddr, facei) + { + nNbrs[upperAddr[facei]]++; + } + + forAll(lowerAddr, facei) + { + nNbrs[lowerAddr[facei]]++; + } + + cellFaceOffsets[0] = 0; + forAll(nNbrs, celli) + { + cellFaceOffsets[celli+1] = cellFaceOffsets[celli] + nNbrs[celli]; + } + + // reset the whole list to use as counter + nNbrs = 0; + + forAll(upperAddr, facei) + { + cellFaces + [ + cellFaceOffsets[upperAddr[facei]] + nNbrs[upperAddr[facei]] + ] = facei; + + nNbrs[upperAddr[facei]]++; + } + + forAll(lowerAddr, facei) + { + cellFaces + [ + cellFaceOffsets[lowerAddr[facei]] + nNbrs[lowerAddr[facei]] + ] = facei; + + nNbrs[lowerAddr[facei]]++; + } + } + + + // go through the faces and create clusters + + tmp tcoarseCellMap(new labelField(nFineCells, -1)); + labelField& coarseCellMap = tcoarseCellMap(); + + nCoarseCells = 0; + label celli; + + for (label cellfi=0; cellfi maxFaceWeight + ) + { + // Match found. Pick up all the necessary data + matchFaceNo = facei; + maxFaceWeight = faceWeights[facei]; + } + } + + if (matchFaceNo >= 0) + { + // Make a new group + coarseCellMap[upperAddr[matchFaceNo]] = nCoarseCells; + coarseCellMap[lowerAddr[matchFaceNo]] = nCoarseCells; + nCoarseCells++; + } + else + { + // No match. Find the best neighbouring cluster and + // put the cell there + label clusterMatchFaceNo = -1; + scalar clusterMaxFaceCoeff = -GREAT; + + for + ( + label faceOs=cellFaceOffsets[celli]; + faceOs clusterMaxFaceCoeff) + { + clusterMatchFaceNo = facei; + clusterMaxFaceCoeff = faceWeights[facei]; + } + } + + if (clusterMatchFaceNo >= 0) + { + // Add the cell to the best cluster + coarseCellMap[celli] = max + ( + coarseCellMap[upperAddr[clusterMatchFaceNo]], + coarseCellMap[lowerAddr[clusterMatchFaceNo]] + ); + } + } + } + } + + // Check that all cells are part of clusters, + // if not create single-cell "clusters" for each + for (label cellfi=0; cellfi Date: Wed, 12 Jun 2013 12:58:39 +0100 Subject: [PATCH 08/12] twoPhaseEulerFoam/bed tutorial: update for new location of alpha solver settings --- .../multiphase/twoPhaseEulerFoam/bed/system/fvSolution | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/tutorials/multiphase/twoPhaseEulerFoam/bed/system/fvSolution b/tutorials/multiphase/twoPhaseEulerFoam/bed/system/fvSolution index 5800e2a4b4..3a47e5857c 100644 --- a/tutorials/multiphase/twoPhaseEulerFoam/bed/system/fvSolution +++ b/tutorials/multiphase/twoPhaseEulerFoam/bed/system/fvSolution @@ -53,6 +53,12 @@ solvers relTol 0; } + alpha1 + { + nAlphaCorr 2; + nAlphaSubCycles 3; + } + alpha { solver PCG; @@ -72,12 +78,9 @@ PIMPLE { nCorrectors 2; nNonOrthogonalCorrectors 0; - nAlphaCorr 2; - nAlphaSubCycles 3; correctAlpha yes; pRefCell 0; pRefValue 0; } - // ************************************************************************* // From a38ff1757c99160469c66158b6d7c441dfc61a60 Mon Sep 17 00:00:00 2001 From: Henry Date: Wed, 12 Jun 2013 13:07:58 +0100 Subject: [PATCH 09/12] displacementLayeredMotionMotionSolver: correct constructor for clang --- .../layeredSolver/displacementLayeredMotionMotionSolver.C | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.C b/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.C index ace954062b..3ea1d3ed39 100644 --- a/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.C +++ b/src/fvMotionSolver/fvMotionSolvers/displacement/layeredSolver/displacementLayeredMotionMotionSolver.C @@ -259,8 +259,10 @@ Foam::displacementLayeredMotionMotionSolver::faceZoneEvaluate // patch. That becomes the value of fld. const word patchName(dict.lookup("patch")); label patchID = mesh().boundaryMesh().findPatchID(patchName); - pointField pdf = - pointDisplacement_.boundaryField()[patchID].patchInternalField(); + pointField pdf + ( + pointDisplacement_.boundaryField()[patchID].patchInternalField() + ); fld = gAverage(pdf); } else From 838639cd5dc138c3d98fe3db76fc2d3dc067c476 Mon Sep 17 00:00:00 2001 From: laurence Date: Wed, 12 Jun 2013 13:08:15 +0100 Subject: [PATCH 10/12] ENH: Snapping to surfaces if point is far outside --- .../conformalVoronoiMeshCalcDualMesh.C | 87 +++++++++++++++++++ 1 file changed, 87 insertions(+) diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C index 9cfa9e316c..c771256f4b 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C @@ -1787,6 +1787,11 @@ void Foam::conformalVoronoiMesh::indexDualVertices } } + OBJstream snapping1("snapToSurface1.obj"); + OBJstream snapping2("snapToSurface2.obj"); + OFstream tetToSnapTo("tetsToSnapTo.obj"); + label offset = 0; + for ( Delaunay::Finite_cells_iterator cit = finite_cells_begin(); @@ -1892,6 +1897,88 @@ void Foam::conformalVoronoiMesh::indexDualVertices } } + { + // Snapping points far outside + if (cit->boundaryDualVertex()) + { + pointFromPoint dual = cit->dual(); + + pointIndexHit hitInfo; + label surfHit; + + // Find nearest surface point + geometryToConformTo_.findSurfaceNearest + ( + dual, + sqr(targetCellSize(dual)), + hitInfo, + surfHit + ); + + if (!hitInfo.hit()) + { + // Project dual to nearest point on tet + + tetPointRef tet + ( + topoint(cit->vertex(0)->point()), + topoint(cit->vertex(1)->point()), + topoint(cit->vertex(2)->point()), + topoint(cit->vertex(3)->point()) + ); + + pointFromPoint nearestPointOnTet = + tet.nearestPoint(dual).rawPoint(); + + // Get nearest point on surface from tet. + geometryToConformTo_.findSurfaceNearest + ( + nearestPointOnTet, + sqr(targetCellSize(nearestPointOnTet)), + hitInfo, + surfHit + ); + + vector snapDir = nearestPointOnTet - dual; + snapDir /= mag(snapDir) + SMALL; + + drawDelaunayCell(tetToSnapTo, cit, offset); + offset += 1; + + vectorField norm(1); + allGeometry_[surfHit].getNormal + ( + List(1, hitInfo), + norm + ); + norm[0] /= mag(norm[0]) + SMALL; + + if + ( + hitInfo.hit() + && (mag(snapDir & norm[0]) > 0.5) + ) + { + snapping1.write + ( + linePointRef(dual, nearestPointOnTet) + ); + + snapping2.write + ( + linePointRef + ( + nearestPointOnTet, + hitInfo.hitPoint() + ) + ); + + pts[cit->cellIndex()] = hitInfo.hitPoint(); + } + } + } + } + if (cit->boundaryDualVertex()) { if (cit->featureEdgeDualVertex()) From 6cf69b164426dabdb5d84270a243960c5493a42d Mon Sep 17 00:00:00 2001 From: laurence Date: Wed, 12 Jun 2013 13:16:04 +0100 Subject: [PATCH 11/12] ENH: foamyHexMesh: Enable only the conformation step to be run if internal points are all fixed --- .../conformalVoronoiMesh.C | 65 +++++++++++++++++-- .../conformalVoronoiMesh.H | 4 ++ .../initialPointsMethod/pointFile/pointFile.C | 2 +- .../generation/foamyHexMesh/foamyHexMesh.C | 39 +++++++++-- 4 files changed, 98 insertions(+), 12 deletions(-) diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C index 3e91db48cb..577d81cdcb 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C @@ -1046,6 +1046,18 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh ) ), decomposition_() +{} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::conformalVoronoiMesh::~conformalVoronoiMesh() +{} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + +void Foam::conformalVoronoiMesh::initialiseForMotion() { if (foamyHexMeshControls().objOutput()) { @@ -1061,7 +1073,10 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh runTime_, rndGen_, geometryToConformTo_, - foamyHexMeshDict.subDict("backgroundMeshDecomposition") + foamyHexMeshControls().foamyHexMeshDict().subDict + ( + "backgroundMeshDecomposition" + ) ) ); } @@ -1125,13 +1140,53 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh } -// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // +void Foam::conformalVoronoiMesh::initialiseForConformation() +{ + if (Pstream::parRun()) + { + decomposition_.reset + ( + new backgroundMeshDecomposition + ( + runTime_, + rndGen_, + geometryToConformTo_, + foamyHexMeshControls().foamyHexMeshDict().subDict + ( + "backgroundMeshDecomposition" + ) + ) + ); + } -Foam::conformalVoronoiMesh::~conformalVoronoiMesh() -{} + insertInitialPoints(); + insertFeaturePoints(); + + // Improve the guess that the backgroundMeshDecomposition makes with the + // initial positions. Use before building the surface conformation to + // better balance the surface conformation load. + distributeBackground(*this); + + buildSurfaceConformation(); + + // The introduction of the surface conformation may have distorted the + // balance of vertices, distribute if necessary. + distributeBackground(*this); + + if (Pstream::parRun()) + { + sync(decomposition_().procBounds()); + } + + cellSizeMeshOverlapsBackground(); + + if (foamyHexMeshControls().printVertexInfo()) + { + printVertexInfo(Info); + } +} -// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // void Foam::conformalVoronoiMesh::move() { diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H index ab9b6d955b..6b514dfdf5 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H @@ -1038,6 +1038,10 @@ public: // Member Functions + void initialiseForMotion(); + + void initialiseForConformation(); + //- Move the vertices according to the controller, re-conforming to the // surface as required void move(); diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/pointFile/pointFile.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/pointFile/pointFile.C index 32ce53d376..a44103923b 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/pointFile/pointFile.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/pointFile/pointFile.C @@ -58,7 +58,7 @@ List pointFile::initialPoints() const IOobject ( pointFileName_.name(), - foamyHexMesh_.time().constant(), + foamyHexMesh_.time().timeName(), foamyHexMesh_.time(), IOobject::MUST_READ, IOobject::NO_WRITE diff --git a/applications/utilities/mesh/generation/foamyHexMesh/foamyHexMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/foamyHexMesh.C index 8133723819..87c1dbe22d 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/foamyHexMesh.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/foamyHexMesh.C @@ -45,12 +45,22 @@ int main(int argc, char *argv[]) "check all surface geometry for quality" ); + Foam::argList::addBoolOption + ( + "conformationOnly", + "conform to the initial points without any point motion" + ); + + #include "addOverwriteOption.H" + #include "setRootCase.H" #include "createTime.H" runTime.functionObjects().off(); const bool checkGeometry = args.optionFound("checkGeometry"); + const bool conformationOnly = args.optionFound("conformationOnly"); + const bool overwrite = args.optionFound("overwrite"); IOdictionary foamyHexMeshDict ( @@ -104,16 +114,33 @@ int main(int argc, char *argv[]) conformalVoronoiMesh mesh(runTime, foamyHexMeshDict); - while (runTime.loop()) + if (conformationOnly) { - Info<< nl << "Time = " << runTime.timeName() << endl; + mesh.initialiseForConformation(); - mesh.move(); + if (!overwrite) + { + runTime++; + } - Info<< nl << "ExecutionTime = " << runTime.elapsedCpuTime() << " s" - << " ClockTime = " << runTime.elapsedClockTime() << " s" - << endl; + mesh.writeMesh(runTime.timeName()); } + else + { + mesh.initialiseForMotion(); + + while (runTime.loop()) + { + Info<< nl << "Time = " << runTime.timeName() << endl; + + mesh.move(); + + Info<< nl << "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << " ClockTime = " << runTime.elapsedClockTime() << " s" + << endl; + } + } + Info<< nl << "End" << nl << endl; From 11d22dbf47de4f0349328b666bd16712c9921e80 Mon Sep 17 00:00:00 2001 From: laurence Date: Thu, 13 Jun 2013 09:30:05 +0100 Subject: [PATCH 12/12] REVERT: Change parallel referral algorithm back to original --- .../DelaunayMesh/DistributedDelaunayMesh.C | 106 ++++++++++++------ .../conformalVoronoiMeshCalcDualMesh.C | 2 +- 2 files changed, 73 insertions(+), 35 deletions(-) diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DistributedDelaunayMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DistributedDelaunayMesh.C index 2370c16ba3..7b316e5c44 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DistributedDelaunayMesh.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DistributedDelaunayMesh.C @@ -254,47 +254,85 @@ void Foam::DistributedDelaunayMesh::findProcessorBoundaryCells /Pstream::nProcs() ); - std::list infinite_cells; - Triangulation::incident_cells - ( - Triangulation::infinite_vertex(), - std::back_inserter(infinite_cells) - ); - - for - ( - typename std::list::iterator vcit = infinite_cells.begin(); - vcit != infinite_cells.end(); - ++vcit - ) - { - Cell_handle cit = *vcit; - - // Index of infinite vertex in this cell. - int i = cit->index(Triangulation::infinite_vertex()); - - Cell_handle c = cit->neighbor(i); - - if (c->unassigned()) - { - c->cellIndex() = this->getNewCellIndex(); - - if (checkProcBoundaryCell(c, circumsphereOverlaps)) - { - cellToCheck.insert(c->cellIndex()); - } - } - } +// std::list infinite_cells; +// Triangulation::incident_cells +// ( +// Triangulation::infinite_vertex(), +// std::back_inserter(infinite_cells) +// ); +// +// for +// ( +// typename std::list::iterator vcit +// = infinite_cells.begin(); +// vcit != infinite_cells.end(); +// ++vcit +// ) +// { +// Cell_handle cit = *vcit; +// +// // Index of infinite vertex in this cell. +// int i = cit->index(Triangulation::infinite_vertex()); +// +// Cell_handle c = cit->neighbor(i); +// +// if (c->unassigned()) +// { +// c->cellIndex() = this->getNewCellIndex(); +// +// if (checkProcBoundaryCell(c, circumsphereOverlaps)) +// { +// cellToCheck.insert(c->cellIndex()); +// } +// } +// } +// +// +// for +// ( +// Finite_cells_iterator cit = Triangulation::finite_cells_begin(); +// cit != Triangulation::finite_cells_end(); +// ++cit +// ) +// { +// if (cit->parallelDualVertex()) +// { +// if (cit->unassigned()) +// { +// if (checkProcBoundaryCell(cit, circumsphereOverlaps)) +// { +// cellToCheck.insert(cit->cellIndex()); +// } +// } +// } +// } for ( - Finite_cells_iterator cit = Triangulation::finite_cells_begin(); - cit != Triangulation::finite_cells_end(); + All_cells_iterator cit = Triangulation::all_cells_begin(); + cit != Triangulation::all_cells_end(); ++cit ) { - if (cit->parallelDualVertex()) + if (Triangulation::is_infinite(cit)) + { + // Index of infinite vertex in this cell. + int i = cit->index(Triangulation::infinite_vertex()); + + Cell_handle c = cit->neighbor(i); + + if (c->unassigned()) + { + c->cellIndex() = this->getNewCellIndex(); + + if (checkProcBoundaryCell(c, circumsphereOverlaps)) + { + cellToCheck.insert(c->cellIndex()); + } + } + } + else if (cit->parallelDualVertex()) { if (cit->unassigned()) { diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C index c771256f4b..e33ce7c6c7 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C @@ -1899,7 +1899,7 @@ void Foam::conformalVoronoiMesh::indexDualVertices { // Snapping points far outside - if (cit->boundaryDualVertex()) + if (cit->boundaryDualVertex() && !cit->parallelDualVertex()) { pointFromPoint dual = cit->dual();