分享以及讨论 欧拉拉格朗日模型中 ParcelCollector



  • 下面是我修改的ParcelCollector.C 及.H 基于OF4.1
    目的是通过在domain中建立一个面(同心圆或多边形),对碰到该面的parcel进行一些统计。
    原code中只计算了mass 以及massflowrate。
    我在做喷雾,在validation中需要做一个sauter mean diameter(D32)VS axial diameter 的图。
    在试过其他方法并没有特别好的效果的情况下,试着修改了一下parcelcollector这个cloudfunction。

    #include "ParticleCollector.H"
    #include "Pstream.H"
    #include "surfaceWriter.H"
    #include "unitConversion.H"
    #include "Random.H"
    #include "triangle.H"
    #include "cloud.H"
    
    // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
    
    template<class CloudType>
    void Foam::ParticleCollector<CloudType>::makeLogFile
    (
        const faceList& faces,
        const Field<point>& points,
        const Field<scalar>& area
    )
    {
        // Create the output file if not already created
        if (log_)
        {
            if (debug)
            {
                Info<< "Creating output file" << endl;
            }
    
            if (Pstream::master())
            {
                // Create directory if does not exist
                mkDir(this->writeTimeDir());
    
                // Open new file at start up
                outputFilePtr_.reset
                (
                    new OFstream(this->writeTimeDir()/(type() + ".dat"))
                );
    
                outputFilePtr_()
                    << "# Source     : " << type() << nl
                    << "# Bins       : " << faces.size() << nl
                    << "# Total area : " << sum(area) << nl;
    
                outputFilePtr_()
                    << "# Geometry   :" << nl
                    << '#'
                    << tab << "Bin"
                    << tab << "(Centre_x Centre_y Centre_z)"
                    << tab << "Area"
                    << nl;
    
                forAll(faces, i)
                {
                    outputFilePtr_()
                        << '#'
                        << tab << i
                        << tab << faces[i].centre(points)
                        << tab << area[i]
                        << nl;
                }
    
                outputFilePtr_()
                    << '#' << nl
                    << "# Output format:" << nl;
    
                forAll(faces, i)
                {
                    word id = Foam::name(i);
                    word binId = "bin_" + id;
    
                    outputFilePtr_()
                        << '#'
                        << tab << "Time"
                        << tab << binId
                        << tab << "mass[" << id << "]"
                        << tab << "massFlowRate[" << id << "]"
                        << tab << "D3[" << id << "]"
                        << tab << "D2[" << id << "]"
                        << endl;
                }
            }
        }
    }
    
    
    template<class CloudType>
    void Foam::ParticleCollector<CloudType>::initPolygons
    (
        const List<Field<point>>& polygons
    )
    {
        mode_ = mtPolygon;
    
        label nPoints = 0;
        forAll(polygons, polyI)
        {
            label np = polygons[polyI].size();
            if (np < 3)
            {
                FatalIOErrorInFunction(this->coeffDict())
                    << "polygons must consist of at least 3 points"
                    << exit(FatalIOError);
            }
    
            nPoints += np;
        }
    
        label pointOffset = 0;
        points_.setSize(nPoints);
        faces_.setSize(polygons.size());
        faceTris_.setSize(polygons.size());
        area_.setSize(polygons.size());
        forAll(faces_, facei)
        {
            const Field<point>& polyPoints = polygons[facei];
            face f(identity(polyPoints.size()) + pointOffset);
            UIndirectList<point>(points_, f) = polyPoints;
            area_[facei] = f.mag(points_);
    
            DynamicList<face> tris;
            f.triangles(points_, tris);
            faceTris_[facei].transfer(tris);
    
            faces_[facei].transfer(f);
    
            pointOffset += polyPoints.size();
        }
    }
    
    
    template<class CloudType>
    void Foam::ParticleCollector<CloudType>::initConcentricCircles()
    {
        mode_ = mtConcentricCircle;
    
        vector origin(this->coeffDict().lookup("origin"));
    
        radius_ = this->coeffDict().lookup("radius");
        nSector_ = readLabel(this->coeffDict().lookup("nSector"));
    
        label nS = nSector_;
    
        vector refDir;
        if (nSector_ > 1)
        {
            refDir = this->coeffDict().lookup("refDir");
            refDir -= normal_[0]*(normal_[0] & refDir);
            refDir /= mag(refDir);
        }
        else
        {
            // set 4 quadrants for single sector cases
            nS = 4;
    
            vector tangent = Zero;
            scalar magTangent = 0.0;
    
            Random rnd(1234);
            while (magTangent < SMALL)
            {
                vector v = rnd.vector01();
    
                tangent = v - (v & normal_[0])*normal_[0];
                magTangent = mag(tangent);
            }
    
            refDir = tangent/magTangent;
        }
    
        scalar dTheta = 5.0;
        scalar dThetaSector = 360.0/scalar(nS);
        label intervalPerSector = max(1, ceil(dThetaSector/dTheta));
        dTheta = dThetaSector/scalar(intervalPerSector);
    
        label nPointPerSector = intervalPerSector + 1;
    
        label nPointPerRadius = nS*(nPointPerSector - 1);
        label nPoint = radius_.size()*nPointPerRadius;
        label nFace = radius_.size()*nS;
    
    
    
        // add origin
        nPoint++;
    
        points_.setSize(nPoint);
        faces_.setSize(nFace);
        area_.setSize(nFace);
    
        coordSys_ = cylindricalCS("coordSys", origin, normal_[0], refDir, false);
    
        List<label> ptIDs(identity(nPointPerRadius));
    
        points_[0] = origin;
    
        // points
        forAll(radius_, radI)
        {
            label pointOffset = radI*nPointPerRadius + 1;
    
            for (label i = 0; i < nPointPerRadius; i++)
            {
                label pI = i + pointOffset;
                point pCyl(radius_[radI], degToRad(i*dTheta), 0.0);
                points_[pI] = coordSys_.globalPosition(pCyl);
            }
        }
    
        // faces
        DynamicList<label> facePts(2*nPointPerSector);
        forAll(radius_, radI)
        {
            if (radI == 0)
            {
                for (label secI = 0; secI < nS; secI++)
                {
                    facePts.clear();
    
                    // append origin point
                    facePts.append(0);
    
                    for (label ptI = 0; ptI < nPointPerSector; ptI++)
                    {
                        label i = ptI + secI*(nPointPerSector - 1);
                        label id = ptIDs.fcIndex(i - 1) + 1;
                        facePts.append(id);
                    }
    
                    label facei = secI + radI*nS;
    
                    faces_[facei] = face(facePts);
                    area_[facei] = faces_[facei].mag(points_);
                }
            }
            else
            {
                for (label secI = 0; secI < nS; secI++)
                {
                    facePts.clear();
    
                    label offset = (radI - 1)*nPointPerRadius + 1;
    
                    for (label ptI = 0; ptI < nPointPerSector; ptI++)
                    {
                        label i = ptI + secI*(nPointPerSector - 1);
                        label id = offset + ptIDs.fcIndex(i - 1);
                        facePts.append(id);
                    }
                    for (label ptI = nPointPerSector-1; ptI >= 0; ptI--)
                    {
                        label i = ptI + secI*(nPointPerSector - 1);
                        label id = offset + nPointPerRadius + ptIDs.fcIndex(i - 1);
                        facePts.append(id);
                    }
    
                    label facei = secI + radI*nS;
    
                    faces_[facei] = face(facePts);
                    area_[facei] = faces_[facei].mag(points_);
                }
            }
        }
    }
    
    
    template<class CloudType>
    void Foam::ParticleCollector<CloudType>::collectParcelPolygon
    (
        const point& p1,
        const point& p2
    ) const
    {
        label dummyNearType = -1;
        label dummyNearLabel = -1;
    
        forAll(faces_, facei)
        {
            const label facePoint0 = faces_[facei][0];
    
            const point& pf = points_[facePoint0];
    
            const scalar d1 = normal_[facei] & (p1 - pf);
            const scalar d2 = normal_[facei] & (p2 - pf);
    
            if (sign(d1) == sign(d2))
            {
                // did not cross polygon plane
                continue;
            }
    
            // intersection point
            const point pIntersect = p1 + (d1/(d1 - d2))*(p2 - p1);
    
            const List<face>& tris = faceTris_[facei];
    
            // identify if point is within poly bounds
            forAll(tris, triI)
            {
                const face& tri = tris[triI];
                triPointRef t
                (
                    points_[tri[0]],
                    points_[tri[1]],
                    points_[tri[2]]
                );
    
                if (t.classify(pIntersect, dummyNearType, dummyNearLabel))
                {
                    hitFaceIDs_.append(facei);
                }
            }
        }
    }
    
    
    template<class CloudType>
    void Foam::ParticleCollector<CloudType>::collectParcelConcentricCircles
    (
        const point& p1,
        const point& p2
    ) const
    {
        label secI = -1;
    
        const scalar d1 = normal_[0] & (p1 - coordSys_.origin());
        const scalar d2 = normal_[0] & (p2 - coordSys_.origin());
    
        if (sign(d1) == sign(d2))
        {
            // did not cross plane
            return;
        }
    
        // intersection point in cylindrical co-ordinate system
        const point pCyl = coordSys_.localPosition(p1 + (d1/(d1 - d2))*(p2 - p1));
    
        scalar r = pCyl[0];
    
        if (r < radius_.last())
        {
            label radI = 0;
            while (r > radius_[radI])
            {
                radI++;
            }
    
            if (nSector_ == 1)
            {
                secI = 4*radI;
            }
            else
            {
                scalar theta = pCyl[1] + constant::mathematical::pi;
    
                secI =
                    nSector_*radI
                  + floor
                    (
                        scalar(nSector_)*theta/constant::mathematical::twoPi
                    );
            }
        }
    
        hitFaceIDs_.append(secI);
    }
    
    
    template<class CloudType>
    void Foam::ParticleCollector<CloudType>::write()
    {
        const fvMesh& mesh = this->owner().mesh();
        const Time& time = mesh.time();
        scalar timeNew = time.value();
        scalar timeElapsed = timeNew - timeOld_;
    
        totalTime_ += timeElapsed;
    
        const scalar alpha = (totalTime_ - timeElapsed)/totalTime_;
        const scalar beta = timeElapsed/totalTime_;
    
        forAll(faces_, facei)
        {
            massFlowRate_[facei] =
                alpha*massFlowRate_[facei] + beta*mass_[facei]/timeElapsed;
            massTotal_[facei] += mass_[facei];
            d3Total_[facei] += d3_[facei];
            d2Total_[facei] += d2_[facei];
        }
    
        const label proci = Pstream::myProcNo();
    
        Field<scalar> faceMassTotal(mass_.size(), 0.0);
        this->getModelProperty("massTotal", faceMassTotal);
    
        Field<scalar> faceMassFlowRate(massFlowRate_.size(), 0.0);
        this->getModelProperty("massFlowRate", faceMassFlowRate);
    
        Field<scalar> faceD3Total(d3_.size(), 0.0);
        this->getModelProperty("d3Total", faceD3Total);
    
        Field<scalar> faceD2Total(d2_.size(), 0.0);
        this->getModelProperty("d2Total", faceD2Total);
    
    
        scalar sumTotalMass = 0.0;
        scalar sumAverageMFR = 0.0;
        scalar sumTotalD3 = 0.0;
        scalar sumTotalD2 = 0.0;
        forAll(faces_, facei)
        {
    
            scalarList allProcMass(Pstream::nProcs());
            allProcMass[proci] = massTotal_[facei];
            Pstream::gatherList(allProcMass);
            faceMassTotal[facei] += sum(allProcMass);
    
            scalarList allProcMassFlowRate(Pstream::nProcs());
            allProcMassFlowRate[proci] = massFlowRate_[facei];
            Pstream::gatherList(allProcMassFlowRate);
            faceMassFlowRate[facei] += sum(allProcMassFlowRate);
    
            scalarList allProcD3(Pstream::nProcs());
            allProcD3[proci] = d3Total_[facei];
            Pstream::gatherList(allProcD3);
            faceD3Total[facei] += sum(allProcD3);
    
            scalarList allProcD2(Pstream::nProcs());
            allProcD2[proci] = d2Total_[facei];
            Pstream::gatherList(allProcD2);
            faceD2Total[facei] += sum(allProcD2);
    
            sumTotalMass += faceMassTotal[facei];
            sumAverageMFR += faceMassFlowRate[facei];
            sumTotalD3 += faceD3Total[facei];
            sumTotalD2 += faceD2Total[facei];
    
            //- Display values for each face
            if (outputFilePtr_.valid())
            {
                outputFilePtr_()
                    << time.timeName()
                    << tab << facei
                    << tab << faceMassTotal[facei]
                    << tab << faceMassFlowRate[facei]
                    << tab << faceD3Total[facei]
                    << tab << faceD2Total[facei]
                    << endl;
            }
        }
    
        Info<< "    sum(total mass) = "             << sumTotalMass     << nl
            << "    sum(average mass flow rate) = " << sumAverageMFR    << nl
            << "    sum(total d3) = "               << sumTotalD3       << nl
            << "    sum(total d2) = "               << sumTotalD2       << nl
            << "    radius.size = "                 << radius_.size()   << nl
    
            << endl;
    
        if (surfaceFormat_ != "none")
        {
            if (Pstream::master())
            {
                autoPtr<surfaceWriter> writer(surfaceWriter::New(surfaceFormat_));
    
                writer->write
                (
                    this->writeTimeDir(),
                    "collector",
                    points_,
                    faces_,
                    "massTotal",
                    faceMassTotal,
                    false
                );
    
                writer->write
                (
                    this->writeTimeDir(),
                    "collector",
                    points_,
                    faces_,
                    "massFlowRate",
                    faceMassFlowRate,
                    false
                );
    
                 writer->write
                (
                    this->writeTimeDir(),
                    "collector",
                    points_,
                    faces_,
                    "d3Total",
                    faceD3Total,
                    false
                );
    
                 writer->write
                (
                    this->writeTimeDir(),
                    "collector",
                    points_,
                    faces_,
                    "d2Total",
                    faceD2Total,
                    false
                );
            }
        }
    
        if (resetOnWrite_)
        {
            Field<scalar> dummy(faceMassTotal.size(), 0.0);
            this->setModelProperty("massTotal", dummy);
            this->setModelProperty("massFlowRate", dummy);
            this->setModelProperty("d3Total", dummy);
            this->setModelProperty("d2Total", dummy);
    
            timeOld_ = timeNew;
            totalTime_ = 0.0;
        }
        else
        {
            this->setModelProperty("massTotal", faceMassTotal);
            this->setModelProperty("massFlowRate", faceMassFlowRate);
            this->setModelProperty("d3Total", faceD3Total);
            this->setModelProperty("d2Total", faceD2Total);
        }
    
        forAll(faces_, facei)
        {
            mass_[facei] = 0.0;
            massTotal_[facei] = 0.0;
            massFlowRate_[facei] = 0.0;
    
            d3_[facei] = 0.0;
            d3Total_[facei] = 0.0;
            d2_[facei] = 0.0;
            d2Total_[facei] = 0.0;
        }
    
    }
    
    
    // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
    
    template<class CloudType>
    Foam::ParticleCollector<CloudType>::ParticleCollector
    (
        const dictionary& dict,
        CloudType& owner,
        const word& modelName
    )
    :
        CloudFunctionObject<CloudType>(dict, owner, modelName, typeName),
        mode_(mtUnknown),
        parcelType_(this->coeffDict().lookupOrDefault("parcelType", -1)),
        removeCollected_(this->coeffDict().lookup("removeCollected")),
        points_(),
        faces_(),
        faceTris_(),
        nSector_(0),
        radius_(),
        coordSys_(false),
        normal_(),
        negateParcelsOppositeNormal_
        (
            readBool(this->coeffDict().lookup("negateParcelsOppositeNormal"))
        ),
        surfaceFormat_(this->coeffDict().lookup("surfaceFormat")),
        resetOnWrite_(this->coeffDict().lookup("resetOnWrite")),
        totalTime_(0.0),
        mass_(),
        massTotal_(),
        massFlowRate_(),
        d3_(),
        d2_(),
        d3Total_(),
        d2Total_(),
        log_(this->coeffDict().lookup("log")),
        outputFilePtr_(),
        timeOld_(owner.mesh().time().value()),
        hitFaceIDs_()
    {
        normal_ /= mag(normal_);
    
        word mode(this->coeffDict().lookup("mode"));
        if (mode == "polygon")
        {
            List<Field<point>> polygons(this->coeffDict().lookup("polygons"));
    
            initPolygons(polygons);
    
            vector n0(this->coeffDict().lookup("normal"));
            normal_ = vectorField(faces_.size(), n0);
        }
        else if (mode == "polygonWithNormal")
        {
            List<Tuple2<Field<point>, vector>> polygonAndNormal
            (
                this->coeffDict().lookup("polygons")
            );
    
            List<Field<point>> polygons(polygonAndNormal.size());
            normal_.setSize(polygonAndNormal.size());
    
            forAll(polygons, polyI)
            {
                polygons[polyI] = polygonAndNormal[polyI].first();
                normal_[polyI] = polygonAndNormal[polyI].second();
                normal_[polyI] /= mag(normal_[polyI]) + ROOTVSMALL;
            }
    
            initPolygons(polygons);
        }
        else if (mode == "concentricCircle")
        {
            vector n0(this->coeffDict().lookup("normal"));
            normal_ = vectorField(1, n0);
    
            initConcentricCircles();
        }
        else
        {
            FatalIOErrorInFunction(this->coeffDict())
                << "Unknown mode " << mode << ".  Available options are "
                << "polygon, polygonWithNormal and concentricCircle"
                << exit(FatalIOError);
        }
    
        mass_.setSize(faces_.size(), 0.0);
        massTotal_.setSize(faces_.size(), 0.0);
        massFlowRate_.setSize(faces_.size(), 0.0);
    
        d3_.setSize(faces_.size(), 0.0);
        d2_.setSize(faces_.size(), 0.0);
        d3Total_.setSize(faces_.size(), 0.0);
        d2Total_.setSize(faces_.size(), 0.0);
    
        makeLogFile(faces_, points_, area_);
    }
    
    
    template<class CloudType>
    Foam::ParticleCollector<CloudType>::ParticleCollector
    (
        const ParticleCollector<CloudType>& pc
    )
    :
        CloudFunctionObject<CloudType>(pc),
        mode_(pc.mode_),
        parcelType_(pc.parcelType_),
        removeCollected_(pc.removeCollected_),
        points_(pc.points_),
        faces_(pc.faces_),
        faceTris_(pc.faceTris_),
        nSector_(pc.nSector_),
        radius_(pc.radius_),
        coordSys_(pc.coordSys_),
        normal_(pc.normal_),
        negateParcelsOppositeNormal_(pc.negateParcelsOppositeNormal_),
        surfaceFormat_(pc.surfaceFormat_),
        resetOnWrite_(pc.resetOnWrite_),
        totalTime_(pc.totalTime_),
        mass_(pc.mass_),
        massTotal_(pc.massTotal_),
        massFlowRate_(pc.massFlowRate_),
        d3_(pc.d3_),
        d2_(pc.d2_),
        d3Total_(pc.d3Total_),
        d2Total_(pc.d2Total_),
        log_(pc.log_),
        outputFilePtr_(),
        timeOld_(0.0),
        hitFaceIDs_()
    {}
    
    
    // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
    
    template<class CloudType>
    Foam::ParticleCollector<CloudType>::~ParticleCollector()
    {}
    
    
    // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
    
    template<class CloudType>
    void Foam::ParticleCollector<CloudType>::postMove
    (
        parcelType& p,
        const label celli,
        const scalar dt,
        const point& position0,
        bool& keepParticle
    )
    {
        if ((parcelType_ != -1) && (parcelType_ != p.typeId()))
        {
            return;
        }
    
        // slightly extend end position to avoid falling within tracking tolerances
        const point position1 = position0 + 1.0001*(p.position() - position0);
    
        hitFaceIDs_.clear();
    
        switch (mode_)
        {
            case mtPolygon:
            {
                collectParcelPolygon(position0, position1);
                break;
            }
            case mtConcentricCircle:
            {
                collectParcelConcentricCircles(position0, position1);
                break;
            }
            default:
            {
            }
        }
    
    
        forAll(hitFaceIDs_, i)
        {
            label facei = hitFaceIDs_[i];
            scalar m = p.nParticle() * p.mass();
            scalar d3p = p.nParticle() * p.d() * p.d()* p.d();
            scalar d2p = p.nParticle() * p.d() * p.d();
            if (negateParcelsOppositeNormal_)
            {
                vector Uhat = p.U();
                Uhat /= mag(Uhat) + ROOTVSMALL;
                if ((Uhat & normal_[facei]) < 0)
                {
                    m *= -1.0;
                }
            }
    
            // add mass contribution
            mass_[facei] += m;
            d2_[facei] += d2p;
            d3_[facei] += d3p;
    
            if (nSector_ == 1)
            {
                mass_[facei + 1] += m;
                mass_[facei + 2] += m;
                mass_[facei + 3] += m;
    
                d2_[facei + 1] += d2p;
                d2_[facei + 2] += d2p;
                d2_[facei + 3] += d2p;
    
                d3_[facei + 1] += d3p;
                d3_[facei + 2] += d3p;
                d3_[facei + 3] += d3p;
    
            }
    
            if (removeCollected_)
            {
                keepParticle = false;
            }
    
        }
    
    }
    
    
    #ifndef ParticleCollector_H
    #define ParticleCollector_H
    
    #include "CloudFunctionObject.H"
    #include "cylindricalCS.H"
    #include "face.H"
    #include "Switch.H"
    #include "OFstream.H"
    
    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    
    namespace Foam
    {
    
    /*---------------------------------------------------------------------------*\
                          Class ParticleCollector Declaration
    \*---------------------------------------------------------------------------*/
    
    template<class CloudType>
    class ParticleCollector
    :
        public CloudFunctionObject<CloudType>
    {
    public:
    
        enum modeType
        {
            mtPolygon,
            mtConcentricCircle,
            mtUnknown
        };
    
    
    private:
    
        // Private Data
    
            // Typedefs
    
                //- Convenience typedef for parcel type
                typedef typename CloudType::parcelType parcelType;
    
            //- Collector mode type
            modeType mode_;
    
            //- Index of parcel types to collect (-1 by default = all particles)
            const label parcelType_;
    
            //- Flag to remove collected particles
            Switch removeCollected_;
    
            //- List of points
            Field<point> points_;
    
            //- List of faces
            List<face> faces_;
    
    
            // Polygon collector
    
                //- Triangulation of faces
                List<List<face>> faceTris_;
    
            // Concentric circles collector
    
                //- Number of sectors per circle
                label nSector_;
    
                //- List of radii
                List<scalar> radius_;
    
                //- Cylindrical co-ordinate system
                cylindricalCS coordSys_;
    
    
            //- Face areas
            Field<scalar> area_;
    
            //- Polygon normal vector per face
            Field<vector> normal_;
    
            //- Remove mass of parcel travelling in opposite direction to normal_
            bool negateParcelsOppositeNormal_;
    
            //- Surface output format
            const word surfaceFormat_;
    
            //- Flag to indicate whether data should be reset/cleared on writing
            Switch resetOnWrite_;
    
            //- Total time
            scalar totalTime_;
    
            //- Mass storage
            List<scalar> mass_;
    
            //- Mass total storage
            List<scalar> massTotal_;
    
            //- Mass flow rate storage
            List<scalar> massFlowRate_;
    
            List<scalar> d3_;
            List<scalar> d2_;
            List<scalar> d3Total_;
            List<scalar> d2Total_;
    
            //- Flag to indicate whether data should be written to file
            Switch log_;
    
            //- Output file pointer
            autoPtr<OFstream> outputFilePtr_;
    
            //- Last calculation time
            scalar timeOld_;
    
            //- Work list to store which faces are hit
            mutable DynamicList<label> hitFaceIDs_;
    
    
        // Private Member Functions
    
            //- Helper function to create log files
            void makeLogFile
            (
                const faceList& faces,
                const Field<point>& points,
                const Field<scalar>& area
            );
    
            //- Initialise polygon collectors
            void initPolygons(const List<Field<point>>& polygons);
    
            //- Initialise concentric circle collectors
            void initConcentricCircles();
    
            //- Collect parcels in polygon collectors
            void collectParcelPolygon
            (
                const point& p1,
                const point& p2
            ) const;
    
            //- Collect parcels in concentric circle collectors
            void collectParcelConcentricCircles
            (
                const point& p1,
                const point& p2
            ) const;
    
    
    protected:
    
        // Protected Member Functions
    
            //- Write post-processing info
            void write();
    
    
    public:
    
        //- Runtime type information
        TypeName("particleCollector");
    
    
        // Constructors
    
            //- Construct from dictionary
            ParticleCollector
            (
                const dictionary& dict,
                CloudType& owner,
                const word& modelName
            );
    
            //- Construct copy
            ParticleCollector(const ParticleCollector<CloudType>& pc);
    
            //- Construct and return a clone
            virtual autoPtr<CloudFunctionObject<CloudType>> clone() const
            {
                return autoPtr<CloudFunctionObject<CloudType>>
                (
                    new ParticleCollector<CloudType>(*this)
                );
            }
    
    
        //- Destructor
        virtual ~ParticleCollector();
    
    
        // Member Functions
    
            // Access
    
                //- Return const access to the reset on write flag
                inline const Switch& resetOnWrite() const;
    
    
            // Evaluation
    
                //- Post-move hook
                virtual void postMove
                (
                    parcelType& p,
                    const label celli,
                    const scalar dt,
                    const point& position0,
                    bool& keepParticle
                );
    };
    
    
    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    
    } // End namespace Foam
    
    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    
    #include "ParticleCollectorI.H"
    
    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    
    #ifdef NoRepository
        #include "ParticleCollector.C"
    #endif
    
    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    
    #endif
    
    // ************************************************************************* //
    


  • 现在的情况是这样的,在相同位置有不同的particlecollector,但是得到的数据就完全不同。

     particleCollector1
        {
            type            particleCollector;
            mode            concentricCircle;
            origin          (0.0 0.0 0.975);
            radius          (0.19999);
            nSector         1;
            refDir          (1 0 0);
            normal          (0 0 1);
            negateParcelsOppositeNormal no;
            removeCollected no;
            surfaceFormat   vtk;
            resetOnWrite    no;
            log             yes;
       }
    
     particleCollector2
        {
            type            particleCollector;
            mode            concentricCircle;
            origin          (0.0 0.0 0.975);
            radius          (0.0005 0.0015 0.0025 0.0035 
    0.0045 0.0055 0.0065 0.0075 
    0.0085 0.0095 0.0105 0.0115 
    0.0125 0.0135 0.0145 0.0155 
    0.0165 0.0175 0.0185 0.0195 
    0.0205 0.0215 0.0225 0.0235 
    0.0245 0.0255 0.0265 0.0275 
    0.0285 0.0295 0.0305 0.05 
    0.075 0.1 0.15 0.19999);
            nSector         1;
            refDir          (1 0 0);
            normal          (0 0 1);
            negateParcelsOppositeNormal no;
            removeCollected no;
            surfaceFormat   vtk;
            resetOnWrite    no;
            log             yes;
        }
    
     particleCollector3
        {
            type            particleCollector;
            mode            concentricCircle;
            origin          (0.0 0.0 0.975);
            radius          (0.0005 0.0015 0.0025 0.0035 0.05 0.075 0.1 0.15 0.19999);
            nSector         1;
            refDir          (1 0 0);
            normal          (0 0 1);
            negateParcelsOppositeNormal no;
            removeCollected no;
            surfaceFormat   vtk;
            resetOnWrite    no;
            log             yes;
        }
    
    

    输出结果如下。

        sum(total mass) = 1.75085e-05
        sum(average mass flow rate) = 0.0875421
        sum(total d3) = 3.34387e-08
        sum(total d2) = 0.000615982
        radius.size = 1
        nSector_ = 1
    
        sum(total mass) = 3.51273e-05
        sum(average mass flow rate) = 0.175637
        sum(total d3) = 6.70882e-08
        sum(total d2) = 0.00122567
        radius.size = 36
        nSector_ = 1
    
        sum(total mass) = 1.51653e-05
        sum(average mass flow rate) = 0.0758266
        sum(total d3) = 2.89636e-08
        sum(total d2) = 0.000605451
        radius.size = 9
        nSector_ = 1
    
    

    如果说 mass和 massflow rate是相同的话,我可以理解,因为这个是原生的部分,我并没有修改什么,但是现在就输出成这个样子,我就很费解了。。希望大神能帮忙看一下


Log in to reply
 

CFD中文网 2016 - 2020 | 京ICP备15017992号-2