diff --git a/src/postProcessing/functionObjects/field/nearWallFields/findCellParticle.C b/src/postProcessing/functionObjects/field/nearWallFields/findCellParticle.C new file mode 100644 index 0000000000..5fe71cac6e --- /dev/null +++ b/src/postProcessing/functionObjects/field/nearWallFields/findCellParticle.C @@ -0,0 +1,226 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 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 . + +\*---------------------------------------------------------------------------*/ + +#include "findCellParticle.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::findCellParticle::findCellParticle +( + const polyMesh& mesh, + const vector& position, + const label cellI, + const label tetFaceI, + const label tetPtI, + const point& end, + const label data +) +: + particle(mesh, position, cellI, tetFaceI, tetPtI), + end_(end), + data_(data) +{} + + +Foam::findCellParticle::findCellParticle +( + const polyMesh& mesh, + Istream& is, + bool readFields +) +: + particle(mesh, is, readFields) +{ + if (readFields) + { + if (is.format() == IOstream::ASCII) + { + is >> end_; + data_ = readLabel(is); + } + else + { + is.read + ( + reinterpret_cast(&end_), + sizeof(end_) + sizeof(data_) + ); + } + } + + // Check state of Istream + is.check + ( + "findCellParticle::findCellParticle" + "(const Cloud&, Istream&, bool)" + ); +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +bool Foam::findCellParticle::move +( + trackingData& td, + const scalar maxTrackLen +) +{ + td.switchProcessor = false; + td.keepParticle = true; + + scalar tEnd = (1.0 - stepFraction())*maxTrackLen; + scalar dtMax = tEnd; + + while (td.keepParticle && !td.switchProcessor && tEnd > SMALL) + { + // set the lagrangian time-step + scalar dt = min(dtMax, tEnd); + + dt *= trackToFace(end_, td); + + tEnd -= dt; + stepFraction() = 1.0 - tEnd/maxTrackLen; + } + + if (tEnd < SMALL || !td.keepParticle) + { + // Hit endpoint or patch. If patch hit could do fancy stuff but just + // to use the patch point is good enough for now. + td.cellToData()[cell()].append(data()); + td.cellToEnd()[cell()].append(position()); + } + + return td.keepParticle; +} + + +bool Foam::findCellParticle::hitPatch +( + const polyPatch&, + trackingData& td, + const label patchI, + const scalar trackFraction, + const tetIndices& tetIs +) +{ + return false; +} + + +void Foam::findCellParticle::hitWedgePatch +( + const wedgePolyPatch&, + trackingData& td +) +{ + // Remove particle + td.keepParticle = false; +} + + +void Foam::findCellParticle::hitSymmetryPatch +( + const symmetryPolyPatch&, + trackingData& td +) +{ + // Remove particle + td.keepParticle = false; +} + + +void Foam::findCellParticle::hitCyclicPatch +( + const cyclicPolyPatch&, + trackingData& td +) +{ + // Remove particle + td.keepParticle = false; +} + + +void Foam::findCellParticle::hitProcessorPatch +( + const processorPolyPatch&, + trackingData& td +) +{ + // Remove particle + td.switchProcessor = true; +} + + +void Foam::findCellParticle::hitWallPatch +( + const wallPolyPatch& wpp, + trackingData& td, + const tetIndices& +) +{ + // Remove particle + td.keepParticle = false; +} + + +void Foam::findCellParticle::hitPatch +( + const polyPatch& wpp, + trackingData& td +) +{ + // Remove particle + td.keepParticle = false; +} + + +// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * // + +Foam::Ostream& Foam::operator<<(Ostream& os, const findCellParticle& p) +{ + if (os.format() == IOstream::ASCII) + { + os << static_cast(p) + << token::SPACE << p.end_ + << token::SPACE << p.data_; + } + else + { + os << static_cast(p); + os.write + ( + reinterpret_cast(&p.end_), + sizeof(p.end_) + sizeof(p.data_) + ); + } + + // Check state of Ostream + os.check("Ostream& operator<<(Ostream&, const findCellParticle&)"); + + return os; +} + + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/nearWallFields/findCellParticle.H b/src/postProcessing/functionObjects/field/nearWallFields/findCellParticle.H new file mode 100644 index 0000000000..7f63671c4d --- /dev/null +++ b/src/postProcessing/functionObjects/field/nearWallFields/findCellParticle.H @@ -0,0 +1,259 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 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 . + +Class + Foam::findCellParticle + +Description + Particle class that finds cells by tracking + +SourceFiles + findCellParticle.C + +\*---------------------------------------------------------------------------*/ + +#ifndef findCellParticle_H +#define findCellParticle_H + +#include "particle.H" +#include "autoPtr.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +class findCellParticleCloud; + +/*---------------------------------------------------------------------------*\ + Class findCellParticle Declaration +\*---------------------------------------------------------------------------*/ + +class findCellParticle +: + public particle +{ + // Private data + + //- end point to track to + point end_; + + //- passive data + label data_; + +public: + + friend class Cloud; + + //- Class used to pass tracking data to the trackToFace function + class trackingData + : + public particle::TrackingData > + { + labelListList& cellToData_; + List >& cellToEnd_; + + public: + + + // Constructors + + trackingData + ( + Cloud& cloud, + labelListList& cellToData, + List >& cellToEnd + ) + : + particle::TrackingData >(cloud), + cellToData_(cellToData), + cellToEnd_(cellToEnd) + {} + + + // Member functions + + labelListList& cellToData() + { + return cellToData_; + } + List >& cellToEnd() + { + return cellToEnd_; + } + }; + + + + // Constructors + + //- Construct from components + findCellParticle + ( + const polyMesh& mesh, + const vector& position, + const label cellI, + const label tetFaceI, + const label tetPtI, + const point& end, + const label data + ); + + //- Construct from Istream + findCellParticle + ( + const polyMesh& mesh, + Istream& is, + bool readFields = true + ); + + //- Construct and return a clone + autoPtr clone() const + { + return autoPtr(new findCellParticle(*this)); + } + + //- Factory class to read-construct particles used for + // parallel transfer + class iNew + { + const polyMesh& mesh_; + + public: + + iNew(const polyMesh& mesh) + : + mesh_(mesh) + {} + + autoPtr operator()(Istream& is) const + { + return autoPtr + ( + new findCellParticle(mesh_, is, true) + ); + } + }; + + + // Member Functions + + //- point to track to + const point& end() const + { + return end_; + } + + //- transported label + label data() const + { + return data_; + } + + + + // Tracking + + //- Track all particles to their end point + bool move(trackingData&, const scalar); + + + //- Overridable function to handle the particle hitting a patch + // Executed before other patch-hitting functions + bool hitPatch + ( + const polyPatch&, + trackingData& td, + const label patchI, + const scalar trackFraction, + const tetIndices& tetIs + ); + + //- Overridable function to handle the particle hitting a wedge + void hitWedgePatch + ( + const wedgePolyPatch&, + trackingData& td + ); + + //- Overridable function to handle the particle hitting a + // symmetryPlane + void hitSymmetryPatch + ( + const symmetryPolyPatch&, + trackingData& td + ); + + //- Overridable function to handle the particle hitting a cyclic + void hitCyclicPatch + ( + const cyclicPolyPatch&, + trackingData& td + ); + + //- Overridable function to handle the particle hitting a + //- processorPatch + void hitProcessorPatch + ( + const processorPolyPatch&, + trackingData& td + ); + + //- Overridable function to handle the particle hitting a wallPatch + void hitWallPatch + ( + const wallPolyPatch&, + trackingData& td, + const tetIndices& + ); + + //- Overridable function to handle the particle hitting a polyPatch + void hitPatch + ( + const polyPatch&, + trackingData& td + ); + + + // Ostream Operator + + friend Ostream& operator<<(Ostream&, const findCellParticle&); +}; + + +template<> +inline bool contiguous() +{ + return true; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/nearWallFields/findCellParticleCloud.C b/src/postProcessing/functionObjects/field/nearWallFields/findCellParticleCloud.C new file mode 100644 index 0000000000..0d5fdf5ea4 --- /dev/null +++ b/src/postProcessing/functionObjects/field/nearWallFields/findCellParticleCloud.C @@ -0,0 +1,42 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 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 . + +\*---------------------------------------------------------------------------*/ + +#include "findCellParticle.H" +#include "Cloud.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +defineTemplateTypeNameAndDebug(Cloud, 0); + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.C b/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.C index 71258bff20..67271ff7cb 100644 --- a/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.C +++ b/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.C @@ -25,6 +25,9 @@ License #include "nearWallFields.H" #include "wordReList.H" +#include "findCellParticle.H" +#include "mappedPatchBase.H" +#include "OBJstream.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -34,6 +37,189 @@ defineTypeNameAndDebug(nearWallFields, 0); } +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +void Foam::nearWallFields::calcAddressing() +{ + const fvMesh& mesh = refCast(obr_); + + // Count number of faces + label nPatchFaces = 0; + forAllConstIter(labelHashSet, patchSet_, iter) + { + label patchI = iter.key(); + nPatchFaces += mesh.boundary()[patchI].size(); + } + + // Global indexing + globalIndex globalCells(mesh.nCells()); + globalIndex globalWalls(nPatchFaces); + + if (debug) + { + Info<< "nearWallFields::calcAddressing() :" + << " nPatchFaces:" << globalWalls.size() << endl; + } + + // Construct cloud + Cloud cloud(mesh, IDLList()); + + // Add particles to track to sample locations + nPatchFaces = 0; + + forAllConstIter(labelHashSet, patchSet_, iter) + { + label patchI = iter.key(); + const fvPatch& patch = mesh.boundary()[patchI]; + + vectorField nf = patch.nf(); + vectorField faceCellCentres = patch.patch().faceCellCentres(); + + forAll(patch, patchFaceI) + { + label meshFaceI = patch.start()+patchFaceI; + + // Find starting point on face (since faceCentre might not + // be on face-diagonal decomposition) + pointIndexHit startInfo + ( + mappedPatchBase::facePoint + ( + mesh, + meshFaceI, + polyMesh::FACEDIAGTETS + ) + ); + + + point start; + if (startInfo.hit()) + { + start = startInfo.hitPoint(); + } + else + { + // Fallback: start tracking from neighbouring cell centre + start = faceCellCentres[patchFaceI]; + } + + const point end = start-distance_*nf[patchFaceI]; + + // Find tet for starting location + label cellI = -1; + label tetFaceI = -1; + label tetPtI = -1; + mesh.findCellFacePt(start, cellI, tetFaceI, tetPtI); + + // Add to cloud. Add originating face as passive data + cloud.addParticle + ( + new findCellParticle + ( + mesh, + start, + cellI, + tetFaceI, + tetPtI, + end, + globalWalls.toGlobal(nPatchFaces) // passive data + ) + ); + + nPatchFaces++; + } + } + + + + if (debug) + { + // Dump particles + OBJstream str + ( + mesh.time().path() + /"wantedTracks_" + mesh.time().timeName() + ".obj" + ); + Info<< "nearWallFields::calcAddressing() :" + << "Dumping tracks to " << str.name() << endl; + + forAllConstIter(Cloud, cloud, iter) + { + const findCellParticle& tp = iter(); + str.write(linePointRef(tp.position(), tp.end())); + } + } + + + + // Per cell: empty or global wall index and end location + cellToWalls_.setSize(mesh.nCells()); + cellToSamples_.setSize(mesh.nCells()); + + // Database to pass into findCellParticle::move + findCellParticle::trackingData td(cloud, cellToWalls_, cellToSamples_); + + // Track all particles to their end position. + scalar maxTrackLen = 2.0*mesh.bounds().mag(); + + + //Debug: collect start points + pointField start; + if (debug) + { + start.setSize(nPatchFaces); + nPatchFaces = 0; + forAllConstIter(Cloud, cloud, iter) + { + const findCellParticle& tp = iter(); + start[nPatchFaces++] = tp.position(); + } + } + + + cloud.move(td, maxTrackLen); + + + // Rework cell-to-globalpatchface into a map + List > compactMap; + getPatchDataMapPtr_.reset + ( + new mapDistribute + ( + globalWalls, + cellToWalls_, + compactMap + ) + ); + + + // Debug: dump resulting tracks + if (debug) + { + getPatchDataMapPtr_().distribute(start); + { + OBJstream str + ( + mesh.time().path() + /"obtainedTracks_" + mesh.time().timeName() + ".obj" + ); + Info<< "nearWallFields::calcAddressing() :" + << "Dumping obtained to " << str.name() << endl; + + forAll(cellToWalls_, cellI) + { + const List& ends = cellToSamples_[cellI]; + const labelList& cData = cellToWalls_[cellI]; + forAll(cData, i) + { + str.write(linePointRef(ends[i], start[cData[i]])); + } + } + } + } +} + + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // Foam::nearWallFields::nearWallFields @@ -127,16 +313,11 @@ void Foam::nearWallFields::read(const dictionary& dict) reverseFieldMap_.insert(sampleFldName, fldName); } - Info<< type() << " " << name_ << ": Creating " << fieldMap_.size() + Info<< type() << " " << name_ << ": Sampling " << fieldMap_.size() << " fields" << endl; - createFields(vsf_); - createFields(vvf_); - createFields(vSpheretf_); - createFields(vSymmtf_); - createFields(vtf_); - - Info<< endl; + // Do analysis + calcAddressing(); } } @@ -147,15 +328,6 @@ void Foam::nearWallFields::execute() { Info<< "nearWallFields:execute()" << endl; } - - //if (active_) - //{ - // sampleFields(vsf_); - // sampleFields(vvf_); - // sampleFields(vSpheretf_); - // sampleFields(vSymmtf_); - // sampleFields(vtf_); - //} } @@ -165,8 +337,6 @@ void Foam::nearWallFields::end() { Info<< "nearWallFields:end()" << endl; } - // Update fields - execute(); } @@ -186,6 +356,28 @@ void Foam::nearWallFields::write() // Do nothing if (active_) { + if + ( + fieldMap_.size() + && vsf_.empty() + && vvf_.empty() + && vSpheretf_.empty() + && vSymmtf_.empty() + && vtf_.empty() + ) + { + Info<< type() << " " << name_ << ": Creating " << fieldMap_.size() + << " fields" << endl; + + createFields(vsf_); + createFields(vvf_); + createFields(vSpheretf_); + createFields(vSymmtf_); + createFields(vtf_); + + Info<< endl; + } + Info<< type() << " " << name_ << " output:" << nl; Info<< " Writing sampled fields to " << obr_.time().timeName() diff --git a/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.H b/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.H index cea4e42139..014313b2aa 100644 --- a/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.H +++ b/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.H @@ -76,6 +76,7 @@ SourceFiles #include "OFstream.H" #include "volFields.H" #include "Tuple2.H" +#include "interpolationCellPoint.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -122,21 +123,32 @@ protected: //- From resulting back to original field HashTable reverseFieldMap_; - //- Locally constructed fields - PtrList vsf_; - PtrList vvf_; - PtrList vSpheretf_; - PtrList vSymmtf_; - PtrList vtf_; + + // Calculated addressing + + //- From cell to seed patch faces + labelListList cellToWalls_; + + //- From cell to tracked end point + List > cellToSamples_; + + //- Map from cell based data back to patch based data + autoPtr getPatchDataMapPtr_; + + + // Locally constructed fields + + PtrList vsf_; + PtrList vvf_; + PtrList vSpheretf_; + PtrList vSymmtf_; + PtrList vtf_; // Protected Member Functions - //- Disallow default bitwise copy construct - nearWallFields(const nearWallFields&); - - //- Disallow default bitwise assignment - void operator=(const nearWallFields&); + //- Calculate addressing from cells back to patch faces + void calcAddressing(); template void createFields @@ -144,6 +156,14 @@ protected: PtrList >& ) const; + //- Override boundary fields with sampled values + template + void sampleBoundaryField + ( + const interpolationCellPoint& interpolator, + GeometricField& fld + ) const; + template void sampleFields ( @@ -151,6 +171,14 @@ protected: ) const; +private: + + //- Disallow default bitwise copy construct + nearWallFields(const nearWallFields&); + + //- Disallow default bitwise assignment + void operator=(const nearWallFields&); + public: //- Runtime type information diff --git a/src/postProcessing/functionObjects/field/nearWallFields/nearWallFieldsTemplates.C b/src/postProcessing/functionObjects/field/nearWallFields/nearWallFieldsTemplates.C index 2541b6fefc..42b1b2abef 100644 --- a/src/postProcessing/functionObjects/field/nearWallFields/nearWallFieldsTemplates.C +++ b/src/postProcessing/functionObjects/field/nearWallFields/nearWallFieldsTemplates.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -24,8 +24,6 @@ License \*---------------------------------------------------------------------------*/ #include "nearWallFields.H" -#include "mappedFieldFvPatchFields.H" -#include "interpolationCellPoint.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -64,35 +62,8 @@ void Foam::nearWallFields::createFields io.rename(sampleFldName); sflds.set(sz, new vfType(io, fld)); - vfType& sampleFld = sflds[sz]; - // Reset the bcs to be mapped - forAllConstIter(labelHashSet, patchSet_, iter) - { - label patchI = iter.key(); - - sampleFld.boundaryField().set - ( - patchI, - new mappedFieldFvPatchField - ( - sampleFld.mesh().boundary()[patchI], - sampleFld.dimensionedInternalField(), - - sampleFld.mesh().name(), - mappedPatchBase::NEARESTCELL, - word::null, // samplePatch - -distance_, - - sampleFld.name(), // fieldName - false, // setAverage - pTraits::zero, // average - interpolationCellPoint::typeName - ) - ); - } - - Info<< " created " << sampleFld.name() << " to sample " + Info<< " created " << sflds[sz].name() << " to sample " << fld.name() << endl; } } @@ -100,6 +71,53 @@ void Foam::nearWallFields::createFields } +template +void Foam::nearWallFields::sampleBoundaryField +( + const interpolationCellPoint& interpolator, + GeometricField& fld +) const +{ + // Construct flat fields for all patch faces to be sampled + Field sampledValues(getPatchDataMapPtr_().constructSize()); + + forAll(cellToWalls_, cellI) + { + const labelList& cData = cellToWalls_[cellI]; + + forAll(cData, i) + { + const point& samplePt = cellToSamples_[cellI][i]; + sampledValues[cData[i]] = interpolator.interpolate(samplePt, cellI); + } + } + + // Send back sampled values to patch faces + getPatchDataMapPtr_().reverseDistribute + ( + getPatchDataMapPtr_().constructSize(), + sampledValues + ); + + // Pick up data + label nPatchFaces = 0; + forAllConstIter(labelHashSet, patchSet_, iter) + { + label patchI = iter.key(); + + fvPatchField& pfld = fld.boundaryField()[patchI]; + + Field newFld(pfld.size()); + forAll(pfld, i) + { + newFld[i] = sampledValues[nPatchFaces++]; + } + + pfld == newFld; + } +} + + template void Foam::nearWallFields::sampleFields ( @@ -115,8 +133,12 @@ void Foam::nearWallFields::sampleFields // Take over internal and boundary values sflds[i] == fld; - // Evaluate to update the mapped - sflds[i].correctBoundaryConditions(); + + // Construct interpolation method + interpolationCellPoint interpolator(fld); + + // Override sampled values + sampleBoundaryField(interpolator, sflds[i]); } }