[bagbagwan]

I need to get yPlus in my particle tracking simulation to make continuous random walk model. I made my own code using "modification of stochastic model in Lagrangian tracking method" by Xu Jundi.

I was able to read yPlus value on my code and was working correctly as I expect right after start. But, program was stopped soon due to lack of memory on computer (my computer has 32 GB).

Whenever yPlus is loaded, it remains on the memory, not to be disappeared.
I think dummy variables are remaining on my memory but I have no idea how to fix cause I am not familiar with C++.

Could you check my code and give me some advice on my code? My code is below. yP is one variable used for yPlus on my code.



% StochasticDispersionCRW.H

/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Copyright (C) 2011-2018 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 <http://www.gnu.org/licenses/>.

Class
    Foam::StochasticDispersionCRW

Description
    The velocity is perturbed in random direction, with a
    Gaussian random number distribution with variance sigma.
    where sigma is defined below

\*---------------------------------------------------------------------------*/

#ifndef StochasticDispersionCRW_H
#define StochasticDispersionCRW_H

#include "DispersionRASModel.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace Foam
{

/*---------------------------------------------------------------------------*\
                   Class StochasticDispersionCRW Declaration
\*---------------------------------------------------------------------------*/

template<class CloudType>
class StochasticDispersionCRW
:
    public DispersionRASModel<CloudType>
{
public:

    //- Runtime type information
    TypeName("stochasticDispersionCRW");

    //yP
    volScalarField yP_;

    scalar f_x;
    scalar f_y;
    scalar f_z;

    // Constructors

        //- Construct from components
        StochasticDispersionCRW(const dictionary& dict, CloudType& owner);

        //- Construct copy
        StochasticDispersionCRW(const StochasticDispersionCRW<CloudType>& dm);

        //- Construct and return a clone
        virtual autoPtr<DispersionModel<CloudType>> clone() const
        {
            return autoPtr<DispersionModel<CloudType>>
            (
                new StochasticDispersionCRW<CloudType>(*this)
            );
        }


    //- Destructor
    virtual ~StochasticDispersionCRW();


    // Member Functions

        //- Update (disperse particles)
        virtual vector update
        (
            const scalar dt,
            const label celli,
            const vector& U,
            const vector& Uc,
            vector& UTurb,
            scalar& tTurb,
            const scalar dxi
        );
};


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace Foam

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#ifdef NoRepository
    #include "StochasticDispersionCRW.C"
#endif

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif

// ************************************************************************* //

% StochasticDispersionCRW.C

/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Copyright (C) 2011-2018 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 <http://www.gnu.org/licenses/>.

\*---------------------------------------------------------------------------*/

#include "StochasticDispersionCRW.H"
#include "constants.H"
#include "nearWallDist.H"
#include "wallFvPatch.H"

using namespace Foam::constant::mathematical;

// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

template<class CloudType>
Foam::StochasticDispersionCRW<CloudType>::StochasticDispersionCRW
(
    const dictionary& dict,
    CloudType& owner
)
:
    DispersionRASModel<CloudType>(dict, owner),
    yP_
       (
         IOobject
         (
          "yP",
          this->owner_.mesh().time().timeName(),
          this->owner_.mesh(),
          IOobject::NO_READ,
          IOobject::AUTO_WRITE
         ),
         this->owner_.mesh(),
         dimensionedScalar("yP", dimless, 0.0)
       ) 
{}


template<class CloudType>
Foam::StochasticDispersionCRW<CloudType>::StochasticDispersionCRW
(
    const StochasticDispllersionCRW<CloudType>& dm
)
:
    DispersionRASModel<CloudType>(dm),
    yP_(dm.yP_)
{}


// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //

template<class CloudType>
Foam::StochasticDispersionCRW<CloudType>::~StochasticDispersionCRW()
{}


// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //

template<class CloudType>
Foam::vector Foam::StochasticDispersionCRW<CloudType>::update
(
    const scalar dt,
    const label celli,
    const vector& U,
    const vector& Uc,
    vector& UTurb,
    scalar& tTurb,
    const scalar dxi
)
{
    Random& rnd = this->owner().rndGen();
    const scalar rndG = rnd.scalar01();
    const symmTensor rooVSmall_R(rootVSmall,rootVSmall,rootVSmall,rootVSmall,rootVSmall,rootVSmall);
    const scalar cps = 0.16432;
    const scalar k = this->kPtr_->primitiveField()[celli];
    const scalar epsilon = this->epsilonPtr_->primitiveField()[celli] + rootVSmall;
    const symmTensor R = this->RPtr_->primitiveField()[celli] + rooVSmall_R;
    const scalar UrelMag = mag(U - Uc - UTurb);
    const scalar tTurbLoc = min(k/epsilon, cps*pow(k, 1.5)/epsilon/(UrelMag + small));

    // Parcel is perturbed by the turbulence
    if (dt < tTurbLoc)
    {
        tTurb += dt;

        if (tTurb > tTurbLoc)
        {
            tTurb = 0 + dxi*0.0  ;

            const fvPatchList& patches = this->owner().mesh().boundary();
            const objectRegistry& obr = this->owner().mesh();
            const word turbName = IOobject::groupName(turbulenceModel::propertiesName,this->owner().U().group());
            const turbulenceModel& turbModel = obr.lookupObject<turbulenceModel>(turbName);
            
            forAll(patches, patchi)
            {
                const fvPatch& currPatch = patches[patchi];
                if (typeid(currPatch) == typeid(wallFvPatch))
                {
                    const scalarField& y_ = turbModel.y()[patchi];
                    const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchi];
                    const tmp<scalarField> tnuw = turbModel.nu(patchi);
                    const scalarField& nuw = tnuw();

                    fvPatchField<scalar>& yPpatch = const_cast<fvPatchField<scalar>&>(yP_.boundaryField()[patchi]);
                    yPpatch = (y_*sqrt(nuw*mag(Uw.snGrad()))/nuw);

                    if (yPpatch[patchi] < 100.0)
                    {
                    scalar f_x = 1+0.285*(yPpatch[patchi]+6)*exp(-0.455*pow(yPpatch[patchi]+6, 0.53))+ small;
                    scalar f_y = 1-exp(-0.02*yPpatch[patchi])+ small;   
                    scalar f_z = sqrt(3-sqr(f_x)*sqr(f_y))+ small;

                    const scalar sigma = sqrt(2*k/3.0+R[0]);
                    const scalar theta = rndG*twoPi;
                    const scalar u = 2*rndG - 1;                  
                    const scalar a = sqrt(1 - sqr(u));
                    const vector dir(f_x*a*cos(theta), f_y*a*sin(theta), f_z*u);;
                    UTurb = sigma*dir;         
                    }
                    else
                    {
                    const scalar sigma = sqrt(2*k/3.0+R[0]);
                    const scalar theta = rndG*twoPi;
                    const scalar u = 2*rndG - 1;
                    const scalar a = sqrt(1 - sqr(u));
                    const vector dir(a*cos(theta), a*sin(theta), u);
                    UTurb = sigma*dir;         
                    }

                }  


     
            }
       
        }
    }

    else
    {
        tTurb = great;
        UTurb = Zero;
    }
    return Uc + UTurb;
}


// ************************************************************************* //

[gxuxg]

Have you got your problem solved? I am having the similar problem when I use the dispersion model created by Jundi. It works fine in serial running model, but breaks down in parallel run. I cannot tell much difference between your version with his, so suspect that the problem for Jundi's code may inherently happen in yours.
I am trying to debug the code now. If you have got your solution, great appreciate your kind update and sharing. Thanks.