diff --git a/bin/tools/pre-commit-hook b/bin/tools/pre-commit-hook index 4265cf1bf5..95ca2278be 100755 --- a/bin/tools/pre-commit-hook +++ b/bin/tools/pre-commit-hook @@ -273,7 +273,7 @@ checkCopyright() for f in $fileList do - sYear=`grep "Copyright.*OpenFOAM" $f | sed 's/[^0-9]//g' | cut -c 5-9` + sYear=`grep "Copyright.*OpenFOAM" $f | sed 's/[^0-9]//g'` if [ "$year" != "" ] && [ "$year" != "$sYear" ]; then echo "Updated copyright for: $f" sed -i "s/$sYear OpenFOAM/$year OpenFOAM/g" $f diff --git a/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointPatchFieldNew.C b/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointPatchFieldNew.C index 01a707a8ba..84d26c9998 100644 --- a/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointPatchFieldNew.C +++ b/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointPatchFieldNew.C @@ -1,25 +1,25 @@ /*---------------------------------------------------------------------------*\ ========= | - \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / F ield |2011 OpenFOAM: The Open Source CFD Toolbox \\ / O peration | \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License - This file is part of OpenFOAM. + This file is part of2011 OpenFOAM. - OpenFOAM is free software: you can redistribute it and/or modify it + 2011 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 + 2011 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 . + along with2011 OpenFOAM. If not, see . \*---------------------------------------------------------------------------*/ @@ -169,7 +169,7 @@ Foam::autoPtr > Foam::pointPatchField::New typename dictionaryConstructorTable::iterator patchTypeCstrIter = dictionaryConstructorTablePtr_->find(p.type()); - if (patchTypeCstrIter == pointPatchConstructorTablePtr_->end()) + if (patchTypeCstrIter == dictionaryConstructorTablePtr_->end()) { FatalIOErrorIn ( diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.C b/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.C index cf433a9cb9..07108b67c4 100644 --- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.C +++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.C @@ -375,7 +375,7 @@ Foam::primitiveMesh::cellTree() const overallBb, 8, // maxLevel 10, // leafsize - 3.0 // duplicity + 5.0 // duplicity ); } diff --git a/src/lagrangian/basic/InteractionLists/InteractionLists.C b/src/lagrangian/basic/InteractionLists/InteractionLists.C index f8a36f656d..572c66522e 100644 --- a/src/lagrangian/basic/InteractionLists/InteractionLists.C +++ b/src/lagrangian/basic/InteractionLists/InteractionLists.C @@ -167,7 +167,7 @@ void Foam::InteractionLists::buildInteractionLists() procBbRndExt, 8, // maxLevel, 10, // leafSize, - 100.0 + 100.0 // duplicity ); ril_.setSize(cellBbsToExchange.size()); @@ -386,7 +386,7 @@ void Foam::InteractionLists::buildInteractionLists() procBbRndExt, 8, // maxLevel, 10, // leafSize, - 100.0 + 100.0 // duplicity ); rwfil_.setSize(wallFaceBbsToExchange.size()); diff --git a/src/meshTools/meshSearch/meshSearch.C b/src/meshTools/meshSearch/meshSearch.C index cb166a452d..115fc3272a 100644 --- a/src/meshTools/meshSearch/meshSearch.C +++ b/src/meshTools/meshSearch/meshSearch.C @@ -30,7 +30,6 @@ License #include "demandDrivenData.H" #include "treeDataCell.H" #include "treeDataFace.H" -#include "treeDataPoint.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -508,6 +507,21 @@ Foam::meshSearch::meshSearch(const polyMesh& mesh, const bool faceDecomp) {} +// Construct with a custom bounding box +Foam::meshSearch::meshSearch +( + const polyMesh& mesh, + const treeBoundBox& bb, + const bool faceDecomp +) +: + mesh_(mesh), + faceDecomp_(faceDecomp) +{ + overallBbPtr_.reset(new treeBoundBox(bb)); +} + + // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // Foam::meshSearch::~meshSearch() @@ -527,6 +541,21 @@ const Foam::indexedOctree& Foam::meshSearch::boundaryTree() // Construct tree // + if (!overallBbPtr_.valid()) + { + Random rndGen(261782); + overallBbPtr_.reset + ( + new treeBoundBox(mesh_.points()) + ); + + treeBoundBox& overallBb = overallBbPtr_(); + // Extend slightly and make 3D + overallBb = overallBb.extend(rndGen, 1E-4); + overallBb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL); + overallBb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL); + } + // all boundary faces (not just walls) labelList bndFaces(mesh_.nFaces()-mesh_.nInternalFaces()); forAll(bndFaces, i) @@ -534,12 +563,6 @@ const Foam::indexedOctree& Foam::meshSearch::boundaryTree() bndFaces[i] = mesh_.nInternalFaces() + i; } - treeBoundBox overallBb(mesh_.points()); - Random rndGen(123456); - overallBb = overallBb.extend(rndGen, 1E-4); - overallBb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL); - overallBb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL); - boundaryTreePtr_.reset ( new indexedOctree @@ -550,7 +573,7 @@ const Foam::indexedOctree& Foam::meshSearch::boundaryTree() mesh_, bndFaces // boundary faces only ), - overallBb, // overall search domain + overallBbPtr_(), // overall search domain 8, // maxLevel 10, // leafsize 3.0 // duplicity @@ -567,13 +590,24 @@ const { if (!cellTreePtr_.valid()) { - treeBoundBox overallBb(mesh_.points()); + // + // Construct tree + // - Random rndGen(261782); + if (!overallBbPtr_.valid()) + { + Random rndGen(261782); + overallBbPtr_.reset + ( + new treeBoundBox(mesh_.points()) + ); - overallBb = overallBb.extend(rndGen, 1E-4); - overallBb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL); - overallBb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL); + treeBoundBox& overallBb = overallBbPtr_(); + // Extend slightly and make 3D + overallBb = overallBb.extend(rndGen, 1E-4); + overallBb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL); + overallBb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL); + } cellTreePtr_.reset ( @@ -584,10 +618,10 @@ const false, // not cache bb mesh_ ), - overallBb, + overallBbPtr_(), 8, // maxLevel 10, // leafsize - 3.0 // duplicity + 6.0 // duplicity ) ); } @@ -904,6 +938,7 @@ void Foam::meshSearch::clearOut() { boundaryTreePtr_.clear(); cellTreePtr_.clear(); + overallBbPtr_.clear(); } diff --git a/src/meshTools/meshSearch/meshSearch.H b/src/meshTools/meshSearch/meshSearch.H index 0488b64684..bf39bca8ff 100644 --- a/src/meshTools/meshSearch/meshSearch.H +++ b/src/meshTools/meshSearch/meshSearch.H @@ -49,6 +49,7 @@ class polyMesh; class treeDataCell; class treeDataFace; template class indexedOctree; +class treeBoundBox; /*---------------------------------------------------------------------------*\ Class meshSearch Declaration @@ -64,8 +65,10 @@ class meshSearch //- Whether to use face decomposition for all geometric tests const bool faceDecomp_; - //- demand driven octrees + //- data bounding box + mutable autoPtr overallBbPtr_; + //- demand driven octrees mutable autoPtr > boundaryTreePtr_; mutable autoPtr > cellTreePtr_; @@ -163,9 +166,19 @@ public: // Constructors - //- Construct from components + //- Construct from components. Constructs bb slightly bigger than + // mesh points bb. meshSearch(const polyMesh& mesh, const bool faceDecomp = true); + //- Construct with a custom bounding box. Any mesh element outside + // bb will not be found. Up to user to make sure bb + // extends slightly beyond wanted elements. + meshSearch + ( + const polyMesh& mesh, + const treeBoundBox& bb, + const bool faceDecomp = true + ); //- Destructor ~meshSearch(); diff --git a/src/postProcessing/functionObjects/field/fieldAverage/controlDict b/src/postProcessing/functionObjects/field/fieldAverage/controlDict index 2c8d3cac79..7d59d161f4 100644 --- a/src/postProcessing/functionObjects/field/fieldAverage/controlDict +++ b/src/postProcessing/functionObjects/field/fieldAverage/controlDict @@ -68,6 +68,7 @@ functions mean on; prime2Mean on; base time; + window 0.01; // optional averaging window } p @@ -75,6 +76,7 @@ functions mean on; prime2Mean on; base time; + window 0.01; // optional averaging window } ); } diff --git a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.C b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.C index 46e64afc3e..be8b03a5b5 100644 --- a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.C +++ b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.C @@ -377,8 +377,7 @@ void Foam::fieldAverage::execute() void Foam::fieldAverage::end() -{ -} +{} void Foam::fieldAverage::write() diff --git a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverageTemplates.C b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverageTemplates.C index e067eb2079..0fdf344531 100644 --- a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverageTemplates.C +++ b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverageTemplates.C @@ -149,7 +149,7 @@ const { typedef GeometricField fieldType; - const scalar dt = obr_.time().deltaTValue(); + scalar dt = obr_.time().deltaTValue(); forAll(faItems_, i) { @@ -163,17 +163,24 @@ const obr_.lookupObject(meanFieldList[i]) ); - scalar alpha = 0.0; - scalar beta = 0.0; - if (faItems_[i].timeBase()) + scalar Dt = totalTime_[i]; + if (faItems_[i].iterBase()) { - alpha = (totalTime_[i] - dt)/totalTime_[i]; - beta = dt/totalTime_[i]; + dt = 1.0; + Dt = scalar(totalIter_[i]); } - else + + scalar alpha = (Dt - dt)/Dt; + scalar beta = dt/Dt; + if (faItems_[i].window() > 0) { - alpha = scalar(totalIter_[i] - 1)/scalar(totalIter_[i]); - beta = 1.0/scalar(totalIter_[i]); + const scalar w = faItems_[i].window(); + + if (Dt - dt >= w) + { + alpha = (w - dt)/w; + beta = dt/w; + } } meanField = alpha*meanField + beta*baseField; @@ -192,7 +199,7 @@ void Foam::fieldAverage::calculatePrime2MeanFields typedef GeometricField fieldType1; typedef GeometricField fieldType2; - const scalar dt = obr_.time().deltaTValue(); + scalar dt = obr_.time().deltaTValue(); forAll(faItems_, i) { @@ -213,17 +220,24 @@ void Foam::fieldAverage::calculatePrime2MeanFields obr_.lookupObject(prime2MeanFieldList[i]) ); - scalar alpha = 0.0; - scalar beta = 0.0; - if (faItems_[i].timeBase()) + scalar Dt = totalTime_[i]; + if (faItems_[i].iterBase()) { - alpha = (totalTime_[i] - dt)/totalTime_[i]; - beta = dt/totalTime_[i]; + dt = 1.0; + Dt = scalar(totalIter_[i]); } - else + + scalar alpha = (Dt - dt)/Dt; + scalar beta = dt/Dt; + if (faItems_[i].window() > 0) { - alpha = scalar(totalIter_[i] - 1)/scalar(totalIter_[i]); - beta = 1.0/scalar(totalIter_[i]); + const scalar w = faItems_[i].window(); + + if (Dt - dt >= w) + { + alpha = (w - dt)/w; + beta = dt/w; + } } prime2MeanField = diff --git a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItem.C b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItem.C index 3c04fc64d4..39bbfb4fa9 100644 --- a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItem.C +++ b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItem.C @@ -53,7 +53,8 @@ Foam::fieldAverageItem::fieldAverageItem() fieldName_("unknown"), mean_(0), prime2Mean_(0), - base_(ITER) + base_(ITER), + window_(-1.0) {} @@ -62,7 +63,8 @@ Foam::fieldAverageItem::fieldAverageItem(const fieldAverageItem& faItem) fieldName_(faItem.fieldName_), mean_(faItem.mean_), prime2Mean_(faItem.prime2Mean_), - base_(faItem.base_) + base_(faItem.base_), + window_(faItem.window_) {} @@ -91,6 +93,7 @@ void Foam::fieldAverageItem::operator=(const fieldAverageItem& rhs) mean_ = rhs.mean_; prime2Mean_ = rhs.prime2Mean_; base_ = rhs.base_; + window_ = rhs.window_; } diff --git a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItem.H b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItem.H index f08328b3d7..10d1e08ad3 100644 --- a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItem.H +++ b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItem.H @@ -33,9 +33,13 @@ Description mean on; prime2Mean on; base time; // iteration + window 200; // optional averaging window } \endverbatim + The averaging window corresponds to the averaging interval (iters or time) + If not specified, the averaging is over 'all iters/time' + SourceFiles fieldAverageItem.C fieldAverageItemIO.C @@ -100,6 +104,9 @@ private: //- Averaging base type baseType base_; + //- Averaging window - defaults to -1 for 'all iters/time' + scalar window_; + public: @@ -148,7 +155,7 @@ public: } //- Return true if base is ITER - Switch ITERBase() const + Switch iterBase() const { return base_ == ITER; } @@ -159,6 +166,11 @@ public: return base_ == TIME; } + scalar window() const + { + return window_; + } + // Member Operators @@ -177,7 +189,8 @@ public: a.fieldName_ == b.fieldName_ && a.mean_ == b.mean_ && a.prime2Mean_ == b.prime2Mean_ - && a.base_ == b.base_; + && a.base_ == b.base_ + && a.window_ == b.window_; } friend bool operator!= diff --git a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItemIO.C b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItemIO.C index 945e438867..4ebad485a7 100644 --- a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItemIO.C +++ b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverageItem/fieldAverageItemIO.C @@ -33,7 +33,9 @@ Foam::fieldAverageItem::fieldAverageItem(Istream& is) : fieldName_("unknown"), mean_(0), - prime2Mean_(0) + prime2Mean_(0), + base_(ITER), + window_(-1.0) { is.check("Foam::fieldAverageItem::fieldAverageItem(Foam::Istream&)"); @@ -43,6 +45,7 @@ Foam::fieldAverageItem::fieldAverageItem(Istream& is) entry.lookup("mean") >> mean_; entry.lookup("prime2Mean") >> prime2Mean_; base_ = baseTypeNames_[entry.lookup("base")]; + window_ = entry.lookupOrDefault("window", -1.0); } @@ -62,6 +65,7 @@ Foam::Istream& Foam::operator>>(Istream& is, fieldAverageItem& faItem) entry.lookup("mean") >> faItem.mean_; entry.lookup("prime2Mean") >> faItem.prime2Mean_; faItem.base_ = faItem.baseTypeNames_[entry.lookup("base")]; + faItem.window_ = entry.lookupOrDefault("window", -1.0); return is; } @@ -80,7 +84,15 @@ Foam::Ostream& Foam::operator<<(Ostream& os, const fieldAverageItem& faItem) os.writeKeyword("prime2Mean") << faItem.mean_ << token::END_STATEMENT << nl; os.writeKeyword("base") << faItem.baseTypeNames_[faItem.base_] - << token::END_STATEMENT << nl << token::END_BLOCK << nl; + << token::END_STATEMENT << nl; + + if (faItem.window_ > 0) + { + os.writeKeyword("window") << faItem.window_ + << token::END_STATEMENT << nl; + } + + os << token::END_BLOCK << nl; os.check ( diff --git a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.C b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.C index a5021f6aa1..ebf14d47f0 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.C +++ b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.C @@ -145,22 +145,7 @@ void Foam::fieldValues::cellSource::initialise(const dictionary& dict) if (operation_ == opWeightedAverage) { dict.lookup("weightField") >> weightFieldName_; - if - ( - obr().foundObject(weightFieldName_) - ) - { - Info<< " weight field = " << weightFieldName_; - } - else - { - FatalErrorIn("cellSource::initialise()") - << type() << " " << name_ << ": " - << sourceTypeNames_[source_] << "(" << sourceName_ << "):" - << nl << " Weight field " << weightFieldName_ - << " must be a " << volScalarField::typeName - << nl << exit(FatalError); - } + Info<< " weight field = " << weightFieldName_; } Info<< nl << endl; diff --git a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.H b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.H index a485310391..212d1730cb 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.H +++ b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSource.H @@ -157,7 +157,8 @@ protected: template tmp > setFieldValues ( - const word& fieldName + const word& fieldName, + const bool mustGet = false ) const; //- Apply the 'operation' to the values diff --git a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSourceTemplates.C b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSourceTemplates.C index 5051c50ca9..6c3fe94438 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSourceTemplates.C +++ b/src/postProcessing/functionObjects/field/fieldValues/cellSource/cellSourceTemplates.C @@ -45,7 +45,8 @@ bool Foam::fieldValues::cellSource::validField(const word& fieldName) const template Foam::tmp > Foam::fieldValues::cellSource::setFieldValues ( - const word& fieldName + const word& fieldName, + const bool mustGet ) const { typedef GeometricField vf; @@ -55,6 +56,20 @@ Foam::tmp > Foam::fieldValues::cellSource::setFieldValues return filterField(obr_.lookupObject(fieldName)); } + if (mustGet) + { + FatalErrorIn + ( + "Foam::tmp > " + "Foam::fieldValues::cellSource::setFieldValues" + "(" + "const word&, " + "const bool" + ") const" + ) << "Field " << fieldName << " not found in database" + << abort(FatalError); + } + return tmp >(new Field(0.0)); } @@ -125,7 +140,13 @@ bool Foam::fieldValues::cellSource::writeValues(const word& fieldName) scalarField V(filterField(mesh().V())); combineFields(V); - scalarField weightField(setFieldValues(weightFieldName_)); + scalarField weightField; + + if (operation_ == opWeightedAverage) + { + weightField = setFieldValues(weightFieldName_, true); + } + combineFields(weightField); if (Pstream::master()) diff --git a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.C b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.C index 10fb1966d7..5141604eed 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.C +++ b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.C @@ -284,23 +284,7 @@ void Foam::fieldValues::faceSource::initialise(const dictionary& dict) if (operation_ == opWeightedAverage) { dict.lookup("weightField") >> weightFieldName_; - if - ( - obr().foundObject(weightFieldName_) - || obr().foundObject(weightFieldName_) - ) - { - Info<< " weight field = " << weightFieldName_; - } - else - { - FatalErrorIn("faceSource::initialise()") - << type() << " " << name_ << ": " - << sourceTypeNames_[source_] << "(" << sourceName_ << "):" - << nl << " Weight field " << weightFieldName_ - << " must be either a " << volScalarField::typeName << " or " - << surfaceScalarField::typeName << nl << exit(FatalError); - } + Info<< " weight field = " << weightFieldName_; } Info<< nl << endl; diff --git a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.H b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.H index d9c2bca621..8ad835e2b2 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.H +++ b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.H @@ -200,7 +200,11 @@ protected: //- Return field values by looking up field name template - tmp > getFieldValues(const word& fieldName) const; + tmp > getFieldValues + ( + const word& fieldName, + const bool mustGet = false + ) const; //- Apply the 'operation' to the values template diff --git a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSourceTemplates.C b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSourceTemplates.C index fabab09707..47fb3a61e1 100644 --- a/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSourceTemplates.C +++ b/src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSourceTemplates.C @@ -52,7 +52,8 @@ bool Foam::fieldValues::faceSource::validField(const word& fieldName) const template Foam::tmp > Foam::fieldValues::faceSource::getFieldValues ( - const word& fieldName + const word& fieldName, + const bool mustGet ) const { typedef GeometricField sf; @@ -74,6 +75,20 @@ Foam::tmp > Foam::fieldValues::faceSource::getFieldValues } } + if (mustGet) + { + FatalErrorIn + ( + "Foam::tmp > " + "Foam::fieldValues::faceSource::getFieldValues" + "(" + "const word&, " + "const bool" + ") const" + ) << "Field " << fieldName << " not found in database" + << abort(FatalError); + } + return tmp >(new Field(0)); } @@ -139,7 +154,13 @@ bool Foam::fieldValues::faceSource::writeValues(const word& fieldName) if (ok) { Field values(getFieldValues(fieldName)); - scalarField weightField(getFieldValues(weightFieldName_)); + scalarField weightField; + + if (operation_ == opWeightedAverage) + { + weightField = getFieldValues(weightFieldName_, true); + } + scalarField magSf; if (surfacePtr_.valid()) diff --git a/src/sampling/meshToMeshInterpolation/meshToMesh/calculateMeshToMeshAddressing.C b/src/sampling/meshToMeshInterpolation/meshToMesh/calculateMeshToMeshAddressing.C index 3b600120ed..a1566c28b3 100644 --- a/src/sampling/meshToMeshInterpolation/meshToMesh/calculateMeshToMeshAddressing.C +++ b/src/sampling/meshToMeshInterpolation/meshToMesh/calculateMeshToMeshAddressing.C @@ -106,7 +106,7 @@ void Foam::meshToMesh::calcAddressing() shiftedBb, // overall bounding box 8, // maxLevel 10, // leafsize - 3.0 // duplicity + 6.0 // duplicity ); if (debug) diff --git a/src/sampling/sampledSet/sampledSet/sampledSet.C b/src/sampling/sampledSet/sampledSet/sampledSet.C index ec53a756fe..36704bbbef 100644 --- a/src/sampling/sampledSet/sampledSet/sampledSet.C +++ b/src/sampling/sampledSet/sampledSet/sampledSet.C @@ -28,6 +28,7 @@ License #include "primitiveMesh.H" #include "meshSearch.H" #include "writer.H" +#include "particle.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -183,26 +184,58 @@ Foam::point Foam::sampledSet::pushIn ) const { label cellI = mesh().faceOwner()[faceI]; + const point& cC = mesh().cellCentres()[cellI]; - const point& cellCtr = mesh().cellCentres()[cellI]; + point newPosition = facePt; - point newSample = - facePt + tol*(cellCtr - facePt); + // Taken from particle::initCellFacePt() + label tetFaceI; + label tetPtI; + mesh().findTetFacePt(cellI, facePt, tetFaceI, tetPtI); - if (!searchEngine().pointInCell(newSample, cellI)) + if (tetFaceI == -1 || tetPtI == -1) + { + newPosition = facePt; + + label trap(1.0/particle::trackingCorrectionTol + 1); + + label iterNo = 0; + + do + { + newPosition += particle::trackingCorrectionTol*(cC - facePt); + + mesh().findTetFacePt + ( + cellI, + newPosition, + tetFaceI, + tetPtI + ); + + iterNo++; + + } while (tetFaceI < 0 && iterNo <= trap); + } + + if (tetFaceI == -1) { FatalErrorIn ( "sampledSet::pushIn(const point&, const label)" - ) << "After pushing " << facePt << " to " << newSample - << " it is still outside faceI " << faceI << endl + ) << "After pushing " << facePt << " to " << newPosition + << " it is still outside face " << faceI + << " at " << mesh().faceCentres()[faceI] + << " of cell " << cellI + << " at " << cC << endl << "Please change your starting point" << abort(FatalError); } - //Info<< "pushIn : moved " << facePt << " to " << newSample + + //Info<< "pushIn : moved " << facePt << " to " << newPosition // << endl; - return newSample; + return newPosition; }