CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM > OpenFOAM Programming & Development

Creating a new field from terms of the turbulence model

Register Blogs Members List Search Today's Posts Mark Forums Read

Like Tree2Likes
  • 2 Post By HaZe

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   September 23, 2013, 08:09
Default Creating a new field from terms of the turbulence model
  #1
New Member
 
Join Date: Jan 2012
Location: Germay - Stuttgart
Posts: 21
Rep Power: 0
HaZe is on a distinguished road
Hi Foamers,

I'm trying to create new fields from the single terms of the k and Epsilon equations of the turbulence model for every time step of writeInterval.

Unfortunately this is my first attempt to do some OF progamming and C++ programming at all. I read a lot in the forum, user and programmer's guide and tried different things recommanded in [1]-[3].
But I did't make any progress until now.

To get an idea of what to do generally I'm trying to create a new field of the diffusion term of k
Code:
- fvm::laplacian(DkEff(), k_)
hoping this will work for all other terms like convection, time derivate and so on accordingly.

My first problem is where to implement the code. In the solvers createFields.H I tried to add

Code:
    volSymmTensorField kDiffusion
    (
        IOobject
        (
            "kDiffusion",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ
            IOobject::AUTO_WRITE
        ),
        - fvm::laplacian(DkEff(), k_)
    );
but there are DkEff and k_ not declared. So I tried to do it in kEpsilon.C in a way like

Code:
tmp<???> kEpsilon::kDiffusion() const
{
    return tmp<???>
    (
        new ???
        (
            IOobject
            (
                "kDiffusion",
                runTime_.timeName(),
                mesh_,
                IOobject::NO_READ,
                IOobject::NO_WRITE
            ),
           - fvm::laplacian(DkEff(), k_)
        )
    );
}
Here are two problems: at first I don't understand of which class ??? is and second would this be the right way to do it anyway?

I already did read about swak4foam which could maybe help me but I'm bound to OF1.5-dev so I can't use it.

And just for my understanding: is it necessary to compile the solver after I made a change in a adjeced lib (libcompressibleRASModels in my special case)? Or do I just have to compile the lib and everything will be ok?

[1] http://openfoamwiki.net/index.php/In...IOobject_class
[2] http://openfoamwiki.net/index.php/Op...objectRegistry
[3] http://www.cfd-online.com/Forums/ope...ta-output.html
HaZe is offline   Reply With Quote

Old   September 24, 2013, 15:48
Default
  #2
New Member
 
Matt Laurita
Join Date: Jan 2012
Location: New York
Posts: 8
Rep Power: 14
mlaurita is on a distinguished road
I think that since both k_ and DfEff() are volScalarFields, the laplacian should also return a volScalarField (i.e. this is what the ??? should be)

Second, using fvm::laplacian is an implicit operation that actually returns an fvScalarMatrix, not what you want here. If you use fvc::laplacian you perform an explicit calculation on the two existing volScalarFields, returning your kDiffusion also as a volScalarField.

-MJL
mlaurita is offline   Reply With Quote

Old   September 25, 2013, 03:59
Default
  #3
New Member
 
Join Date: Jan 2012
Location: Germay - Stuttgart
Posts: 21
Rep Power: 0
HaZe is on a distinguished road
Yeah thanks! Now it's at least compiling without errors and I did understand the difference of fvc/fvm a little bit more even if I had read the Programmer's guide several times :-)

Next problem: I expected a new file in my time directories called "kDiffusion". But there is none.

My code so far (green for my additions/changes):

kEpsilon.H
Code:
 //- Return the turbulence kinetic energy dissipation rate
        tmp<volScalarField> epsilon() const
        {
            return epsilon_;
        }
//k Epsilon Details by HaZe 25.9.2013

        //- Return the kinetic energy diffusion
       tmp<volScalarField> kDiffusion() const;
            
//HaZe end
        //- Return the Reynolds stress tensor
        tmp<volSymmTensorField> R() const;
kEpsilon.C
Code:
tmp<volSymmTensorField> kEpsilon::devRhoReff() const
{
    return tmp<volSymmTensorField>
    (
        new volSymmTensorField
        (
            IOobject
            (
                "devRhoReff",
                runTime_.timeName(),
                mesh_,
                IOobject::NO_READ,
                IOobject::NO_WRITE
            ),
           -muEff()*dev(twoSymm(fvc::grad(U_)))
        )
    );
}
// HaZe 25.9.2013
tmp<volScalarField> kEpsilon::kDiffusion() const
{
    return tmp<volScalarField>
    (
        new volScalarField
        (
            IOobject
            (
                "kDiffusion",
                runTime_.timeName(),
                mesh_,
                IOobject::NO_READ,
                IOobject::AUTO_WRITE
            ),
            -fvc::laplacian(DkEff(), k_)
        )
    );
}
//haze ende

tmp<fvVectorMatrix> kEpsilon::divDevRhoReff(volVectorField& U) const
{
    return
    (
      - fvm::laplacian(muEff(), U)
      - fvc::div(muEff()*dev2(fvc::grad(U)().T()))
    );
}
After changing the code I compiled the lib and the solver. The solver runs but I'm missing the output of "kDiffusion". Is there anything else to change?
HaZe is offline   Reply With Quote

Old   September 25, 2013, 11:12
Default
  #4
New Member
 
Matt Laurita
Join Date: Jan 2012
Location: New York
Posts: 8
Rep Power: 14
mlaurita is on a distinguished road
Well, you're part of the way there. You've successfully added a function to your turbulence model that calculates kDiffusion, but you haven't actually called it from anywhere.

For example, somewhere in your kEpsilon.C file (probably inside of the update() function), you could declare a variable, e.g.:

volScalarField kDiffusion
(
IOobject
(
"kDiffusion",
runTime_.timeName(),
mesh_,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
kDiffusion()
);

Here, the value returned by your kDiffusion() function is assigned to a variable also called kDiffusion. If you do this after the k equation has been solved you'll get kDiffusion based on the latest values. Then you can write the field with:

kDiffusion.write();

There may also be a way to do it by adding kDiffusion_ as a member of the kEpsilon class and assigning it in the constructor. This could lead to it auto writing with the rest of the variables, but not having tried this myself I'm not sure of the details.

MJL
mlaurita is offline   Reply With Quote

Old   September 26, 2013, 08:56
Default
  #5
New Member
 
Join Date: Jan 2012
Location: Germay - Stuttgart
Posts: 21
Rep Power: 0
HaZe is on a distinguished road
Quote:
Originally Posted by mlaurita View Post
Well, you're part of the way there. You've successfully added a function to your turbulence model that calculates kDiffusion, but you haven't actually called it from anywhere.
This sounds logical to me.

I followed your instructions but it's still not working. There's no update() function as you mentioned so I'm not sure where to place it. I tried several locations inside kEpsilon.C but it was always wrong. This is my actual code for kEpsilon.C:

Code:
   
#include "kEpsilon.H"
#include "addToRunTimeSelectionTable.H"
#include "wallFvPatch.H"
            
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
            
namespace Foam
{
namespace compressible
{
namespace RASModels
{
            
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
               
defineTypeNameAndDebug(kEpsilon, 0);
addToRunTimeSelectionTable(RASModel, kEpsilon, dictionary);
                
// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
     
kEpsilon::kEpsilon
(
    const volScalarField& rho,
    const volVectorField& U,
    const surfaceScalarField& phi,
    basicThermo& thermophysicalModel
)
:
    RASModel(typeName, rho, U, phi, thermophysicalModel),
      
    Cmu_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "Cmu",
            coeffDict_,
            0.09
        )
    ),
    C1_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "C1",
            coeffDict_,
            1.44
        )
    ),   
    C2_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "C2",
            coeffDict_,
            1.92
        )
    ),          
    C3_ 
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "C3",
            coeffDict_,
            -0.33 
        )
    ),
    alphak_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "alphak",
            coeffDict_,
            1.0
        )
    ),
    alphaEps_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "alphaEps",
            coeffDict_,
            0.76923
        )
    ),
    alphah_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "alphah",
            coeffDict_,
            1.0
        )
    ),
                
    k_
    (
        IOobject
        (
            "k",
            runTime_.timeName(),
            mesh_,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh_
    ),
     
    epsilon_
    (
        IOobject
        (
            "epsilon",
            runTime_.timeName(),
            mesh_,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh_
    ),   
      
    mut_
    (
        IOobject
        (
            "mut",
            runTime_.timeName(),
            mesh_,
            IOobject::NO_READ,
            IOobject::NO_WRITE
        ),
        Cmu_*rho_*sqr(k_)/(epsilon_ + epsilonSmall_)
    )                

{
#   include "wallViscosityI.H"   

    volScalarField kDiffusion
    (
        IOobject
        (
            "kDiffusion",
            runTime_.timeName(),
            mesh_,
            IOobject::NO_READ,
            IOobject::AUTO_WRITE
        ),
        kDiffusion()
    );
    
    kDiffusion.write();

    printCoeffs();
}
    
     
// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
      
tmp<volSymmTensorField> kEpsilon::R() const
{
    return tmp<volSymmTensorField>
    (
        new volSymmTensorField
        (
            IOobject
            (
                "R",
                runTime_.timeName(),  
                mesh_,
                IOobject::NO_READ,
                IOobject::AUTO_WRITE
            ),
            ((2.0/3.0)*I)*k_ - (mut_/rho_)*dev(twoSymm(fvc::grad(U_))),
            k_.boundaryField().types()
        )
    );  
}
        
tmp<volSymmTensorField> kEpsilon::devRhoReff() const
{
    return tmp<volSymmTensorField>
    (
        new volSymmTensorField
        (
            IOobject
            (
                "devRhoReff",
                runTime_.timeName(),
                mesh_,
                IOobject::NO_READ,
                IOobject::NO_WRITE
            ),
           -muEff()*dev(twoSymm(fvc::grad(U_)))
        )
    );
}

//haze
tmp<volScalarField> kEpsilon::kDiffusion() const
{
    return tmp<volScalarField>
    (
        new volScalarField
        (
            IOobject
            (
                "kDiffusion",
                runTime_.timeName(),
                mesh_,        
                IOobject::NO_READ,
                IOobject::AUTO_WRITE
            ),
            -fvc::laplacian(DkEff(), k_)
        )
    );
}
//haze ende

tmp<fvVectorMatrix> kEpsilon::divDevRhoReff(volVectorField& U) const
{
    return
    (
      - fvm::laplacian(muEff(), U)
      - fvc::div(muEff()*dev2(fvc::grad(U)().T()))
    );
}
        
bool kEpsilon::read()
{
    if (RASModel::read())
    {
        Cmu_.readIfPresent(coeffDict_);   
        C1_.readIfPresent(coeffDict_);
        C2_.readIfPresent(coeffDict_);
        C3_.readIfPresent(coeffDict_);
        alphak_.readIfPresent(coeffDict_);
        alphaEps_.readIfPresent(coeffDict_);
        alphah_.readIfPresent(coeffDict_);

        return true;
    }
    else
    {
        return false;
    }
}


void kEpsilon::correct()
{
    // Bound in case of topological change.  Signalling currently not available
    // HJ, 22/Aug/2007
    bound(k_, k0_);
    bound(epsilon_, epsilon0_);

    if (!turbulence_)
    {
        // Re-calculate viscosity
        mut_ = rho_*Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
#       include "wallViscosityI.H"
        return;
    }

    RASModel::correct();

    volScalarField divU = fvc::div(phi_/fvc::interpolate(rho_));

    if (mesh_.moving())
    {
        divU += fvc::div(mesh_.phi());
    }

    tmp<volTensorField> tgradU = fvc::grad(U_);
    volScalarField G = mut_*(tgradU() && dev(twoSymm(tgradU())));
    tgradU.clear();

#   include "wallFunctionsI.H"

    // Dissipation equation
    tmp<fvScalarMatrix> epsEqn
    (
        fvm::ddt(rho_, epsilon_)
      + fvm::div(phi_, epsilon_)
      - fvm::laplacian(DepsilonEff(), epsilon_)
     ==
        C1_*G*epsilon_/k_
      - fvm::SuSp(((2.0/3.0)*C1_ + C3_)*rho_*divU, epsilon_)
      - fvm::Sp(C2_*rho_*epsilon_/k_, epsilon_)
    );

#   include "wallDissipationI.H"

    epsEqn().relax();

    solve(epsEqn);
    bound(epsilon_, epsilon0_);


    // Turbulent kinetic energy equation

    tmp<fvScalarMatrix> kEqn
    (
        fvm::ddt(rho_, k_)
      + fvm::div(phi_, k_)
      - fvm::laplacian(DkEff(), k_)
     ==
        G
      - fvm::SuSp((2.0/3.0)*rho_*divU, k_)
      - fvm::Sp(rho_*epsilon_/k_, k_)
    );

    kEqn().relax();
    solve(kEqn);
    bound(k_, k0_);


    // Re-calculate viscosity
    mut_ = rho_*Cmu_*sqr(k_)/epsilon_;

#   include "wallViscosityI.H"

}


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

} // End namespace RASModels
} // End namespace compressible
} // End namespace Foam
I get the compiler error:

Code:
kEpsilon/kEpsilon.C: In constructor ‘Foam::compressible::RASModels::kEpsilon::kEpsilon(const Foam::volScalarField&, const Foam::volVectorField&, const Foam::surfaceScalarField&, Foam::basicThermo&)’:
kEpsilon/kEpsilon.C:172: error: no match for call to ‘(Foam::volScalarField) ()’
Among a lot of other variants where to place the code I also tried according to this post
and the following post.
Code:
  
#include "kEpsilon.H"
#include "addToRunTimeSelectionTable.H"
#include "wallFvPatch.H"
            
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
            
namespace Foam
{
namespace compressible
{
namespace RASModels
{
            
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
               
defineTypeNameAndDebug(kEpsilon, 0);
addToRunTimeSelectionTable(RASModel, kEpsilon, dictionary);
                
// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
     
kEpsilon::kEpsilon
(
    const volScalarField& rho,
    const volVectorField& U,
    const surfaceScalarField& phi,
    basicThermo& thermophysicalModel
)
:
    RASModel(typeName, rho, U, phi, thermophysicalModel),
      
    Cmu_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "Cmu",
            coeffDict_,
            0.09
        )
    ),
    C1_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "C1",
            coeffDict_,
            1.44
        )
    ),   
    C2_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "C2",
            coeffDict_,
            1.92
        )
    ),          
    C3_ 
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "C3",
            coeffDict_,
            -0.33 
        )
    ),
    alphak_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "alphak",
            coeffDict_,
            1.0
        )
    ),
    alphaEps_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "alphaEps",
            coeffDict_,
            0.76923
        )
    ),
    alphah_
    (
        dimensioned<scalar>::lookupOrAddToDict
        (
            "alphah",
            coeffDict_,
            1.0
        )
    ),
                
    k_
    (
        IOobject
        (
            "k",
            runTime_.timeName(),
            mesh_,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh_
    ),
     
    epsilon_
    (
        IOobject
        (
            "epsilon",
            runTime_.timeName(),
            mesh_,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh_
    ),   
      
    mut_
    (
        IOobject
        (
            "mut",
            runTime_.timeName(),
            mesh_,
            IOobject::NO_READ,
            IOobject::NO_WRITE
        ),
        Cmu_*rho_*sqr(k_)/(epsilon_ + epsilonSmall_)
    ),
                
    kDiffusion_ 
    (
        IOobject
        (
            "kDiffusion",
            runTime_.timeName(),
            mesh_,
            IOobject::NO_READ,
            IOobject::AUTO_WRITE
        ),
        kDiffusion()
    )

{
#   include "wallViscosityI.H"   

if(runTime_.outputTime())
    {
        kDiffusion_.write();
    }
    printCoeffs();
}
    
     
// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
      
tmp<volSymmTensorField> kEpsilon::R() const
{
    return tmp<volSymmTensorField>
    (
        new volSymmTensorField
        (
            IOobject
            (
                "R",
                runTime_.timeName(),  
                mesh_,
                IOobject::NO_READ,
                IOobject::AUTO_WRITE
            ),
            ((2.0/3.0)*I)*k_ - (mut_/rho_)*dev(twoSymm(fvc::grad(U_))),
            k_.boundaryField().types()
        )
    );  
}
        
tmp<volSymmTensorField> kEpsilon::devRhoReff() const
{
    return tmp<volSymmTensorField>
    (
        new volSymmTensorField
        (
            IOobject
            (
                "devRhoReff",
                runTime_.timeName(),
                mesh_,
                IOobject::NO_READ,
                IOobject::NO_WRITE
            ),
           -muEff()*dev(twoSymm(fvc::grad(U_)))
        )
    );
}

//haze
tmp<volScalarField> kEpsilon::kDiffusion() const
{
    return tmp<volScalarField>
    (
        new volScalarField
        (
            IOobject
            (
                "kDiffusion",
                runTime_.timeName(),
                mesh_,        
                IOobject::NO_READ,
                IOobject::AUTO_WRITE
            ),
            -fvc::laplacian(DkEff(), k_)
        )
    );
}
//haze ende

tmp<fvVectorMatrix> kEpsilon::divDevRhoReff(volVectorField& U) const
{
    return
    (
      - fvm::laplacian(muEff(), U)
      - fvc::div(muEff()*dev2(fvc::grad(U)().T()))
    );
}
        
bool kEpsilon::read()
{
    if (RASModel::read())
    {
        Cmu_.readIfPresent(coeffDict_);   
        C1_.readIfPresent(coeffDict_);
        C2_.readIfPresent(coeffDict_);
        C3_.readIfPresent(coeffDict_);
        alphak_.readIfPresent(coeffDict_);
        alphaEps_.readIfPresent(coeffDict_);
        alphah_.readIfPresent(coeffDict_);

        return true;
    }
    else
    {
        return false;
    }
}


void kEpsilon::correct()
{
    // Bound in case of topological change.  Signalling currently not available
    // HJ, 22/Aug/2007
    bound(k_, k0_);
    bound(epsilon_, epsilon0_);

    if (!turbulence_)
    {
        // Re-calculate viscosity
        mut_ = rho_*Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
#       include "wallViscosityI.H"
        return;
    }

    RASModel::correct();

    volScalarField divU = fvc::div(phi_/fvc::interpolate(rho_));

    if (mesh_.moving())
    {
        divU += fvc::div(mesh_.phi());
    }

    tmp<volTensorField> tgradU = fvc::grad(U_);
    volScalarField G = mut_*(tgradU() && dev(twoSymm(tgradU())));
    tgradU.clear();

#   include "wallFunctionsI.H"

    // Dissipation equation
    tmp<fvScalarMatrix> epsEqn
    (
        fvm::ddt(rho_, epsilon_)
      + fvm::div(phi_, epsilon_)
      - fvm::laplacian(DepsilonEff(), epsilon_)
     ==
        C1_*G*epsilon_/k_
      - fvm::SuSp(((2.0/3.0)*C1_ + C3_)*rho_*divU, epsilon_)
      - fvm::Sp(C2_*rho_*epsilon_/k_, epsilon_)
    );

#   include "wallDissipationI.H"

    epsEqn().relax();

    solve(epsEqn);
    bound(epsilon_, epsilon0_);


    // Turbulent kinetic energy equation

    tmp<fvScalarMatrix> kEqn
    (
        fvm::ddt(rho_, k_)
      + fvm::div(phi_, k_)
      - fvm::laplacian(DkEff(), k_)
     ==
        G
      - fvm::SuSp((2.0/3.0)*rho_*divU, k_)
      - fvm::Sp(rho_*epsilon_/k_, k_)
    );

    kEqn().relax();
    solve(kEqn);
    bound(k_, k0_);


    // Re-calculate viscosity
    mut_ = rho_*Cmu_*sqr(k_)/epsilon_;

#   include "wallViscosityI.H"

}


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

} // End namespace RASModels
} // End namespace compressible
} // End namespace Foam

// ************************************************************************* //
with kEpsilon.H:

Code:
\*---------------------------------------------------------------------------*/
         
#ifndef compressiblekEpsilon_H
#define compressiblekEpsilon_H

#include "RASModel.H"
        
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
            
namespace Foam
{
namespace compressible
{
namespace RASModels
{
        
/*---------------------------------------------------------------------------*\
                           Class kEpsilon Declaration
\*---------------------------------------------------------------------------*/
                
class kEpsilon
:
    public RASModel
{
    // Private data
         
        // Model coefficients
             
//            dimensionedScalar Cmu;
            dimensionedScalar Cmu_;
            dimensionedScalar C1_;
            dimensionedScalar C2_;
            dimensionedScalar C3_;
            dimensionedScalar alphak_;
            dimensionedScalar alphaEps_;
            dimensionedScalar alphah_;

        // Fields

            volScalarField k_;
            volScalarField epsilon_;
            volScalarField mut_;
//Haze add
            volScalarField kDiffusion_;
//End Haze

public:

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

    // Constructors

        //- from components
        kEpsilon
        (
            const volScalarField& rho,
            const volVectorField& U,
            const surfaceScalarField& phi,
            basicThermo& thermophysicalModel
        );


    // Destructor

        ~kEpsilon()
        {}


    // Member Functions

      //- Return the turbulence viscosity
        tmp<volScalarField> mut() const
        {
            return mut_;
        }

        //- Return the effective diffusivity for k
        tmp<volScalarField> DkEff() const
        {
            return tmp<volScalarField>
            (
                new volScalarField("DkEff", alphak_*mut_ + mu())
            );
        }   

        //- Return the effective diffusivity for epsilon
        tmp<volScalarField> DepsilonEff() const
        {
            return tmp<volScalarField>
            (
                new volScalarField("DepsilonEff", alphaEps_*mut_ + mu())
            );
        }

        //- Return the effective turbulent thermal diffusivity
        tmp<volScalarField> alphaEff() const
        {
            return tmp<volScalarField>
            (
                new volScalarField("alphaEff", alphah_*mut_ + alpha())
            );
        }
             
        //- Return the turbulence kinetic energy
        tmp<volScalarField> k() const
        {
            return k_;
        }
            
        //- Return the turbulence kinetic energy dissipation rate
        tmp<volScalarField> epsilon() const
        {
            return epsilon_;
        }
//k Epsilon Details by HaZE 23.9.2013
            
        //- Return the kinetic energy diffusion
        tmp<volScalarField> kDiffusion() const;
            
//Details end
        //- Return the Reynolds stress tensor
        tmp<volSymmTensorField> R() const;

        //- Return the effective stress tensor including the laminar stress
        tmp<volSymmTensorField> devRhoReff() const;

        //- Return the source term for the momentum equation
        tmp<fvVectorMatrix> divDevRhoReff(volVectorField& U) const;
        
        //- Solve the turbulence equations and correct the turbulence viscosity
        void correct();
            
        //- Read turbulenceProperties dictionary
        bool read();
};
          

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
    
} // End namespace RASModels
} // End namespace compressible
} // End namespace Foam
Produces:

Code:
kEpsilon/kEpsilon.C: In constructor ‘Foam::compressible::RASModels::kEpsilon::kEpsilon(const Foam::volScalarField&, const Foam::volVectorField&, const Foam::surfaceScalarField&, Foam::basicThermo&)’:
kEpsilon/kEpsilon.C:176: error: ‘((Foam::compressible::RASModels::kEpsilon*)this)->Foam::compressible::RASModels::kEpsilon::kDiffusion’ does not have class type
Now I'm really stuck and don't know what else to try.

Could you please help me another time. I also don't understand why I have to create a second IOobject for the variable "kDiffusion". Isn't this already done in my member function?
HaZe is offline   Reply With Quote

Old   September 27, 2013, 10:14
Default
  #6
New Member
 
Matt Laurita
Join Date: Jan 2012
Location: New York
Posts: 8
Rep Power: 14
mlaurita is on a distinguished road
My mistake, I meant correct() (not update()).

But after thinking about this a little more I think you're almost there...
in the assignment section of the constructor try:

mut_ ( IOobject ( "mut", runTime_.timeName(), mesh_, IOobject::NO_READ, IOobject::NO_WRITE ), Cmu_*rho_*sqr(k_)/(epsilon_ + epsilonSmall_) ), kDiffusion_ ( IOobject ( "kDiffusion", runTime_.timeName(), mesh_, IOobject::MUST_READ, IOobject::AUTO_WRITE ), mesh_ )
{
kDiffusion_ = kDiffusion();
kDiffusion_.updateBoundaryConditions();
}

Then, your kDiffusion() function can simply be:
tmp<volScalarField> kEpsilon::kDiffusion() const { return (
-fvc::laplacian(DkEff(), k_) );
}


And finally, at the end of the correct() function you can add:
mut_ = rho_*Cmu_*sqr(k_)/epsilon_; # include "wallViscosityI.H"

kDiffusion_ = kDiffusion();
kDiffusion_.updateBoundaryConditions();

}


This creates kDiffusion_ as a variable that is updated with the latest value of k_ and then automatically written to the disk whenever the solution is written. I'm not 100% sure you need the call to updateBoundaryConditions(), but I was basing this off of the way that nut is handled in the incompressible version. One catch is that you need to have a "kDiffusion" file in the 0/ directory, but you can specify an internal field of 0 and boundary fields as "calculated".

I tried this with the incompressible version of kEpsilon and it worked fine, though I don't know if the resulting "kDiffusion" values that are output are accurate, but that'll be up to you.

One last thing, it's always good practice to copy the turbulence model to a new folder and give the files a new name ("myKEpsilon.H", "myKEpsilon.C", for example). Also be sure to change the typename in both the .H and .C files. That way you still have the original one available and you can select your model by specifying "RASModel myKEpsilon;" in the RASProperties dictionary.

Good luck!

MJL
mlaurita is offline   Reply With Quote

Old   September 27, 2013, 10:18
Default
  #7
New Member
 
Matt Laurita
Join Date: Jan 2012
Location: New York
Posts: 8
Rep Power: 14
mlaurita is on a distinguished road
Sorry, let me try again with the code tags this time!

Code:
    mut_
    (
        IOobject
        (
            "mut",              
            runTime_.timeName(),             
            mesh_,             
            IOobject::NO_READ,              
            IOobject::NO_WRITE         
        ), 
        Cmu_*rho_*sqr(k_)/(epsilon_ + epsilonSmall_)     
    ),                     
    kDiffusion_     
    (         
        IOobject
        (              
            "kDiffusion", 
            runTime_.timeName(),
            mesh_,              
            IOobject::MUST_READ,             
            IOobject::AUTO_WRITE         
        ), 
        mesh_     
    )
{
         kDiffusion_ = kDiffusion();
         kDiffusion_.updateBoundaryConditions();
}
Code:
tmp<volScalarField> kEpsilon::kDiffusion() const 
{     
    return     
    (
                -fvc::laplacian(DkEff(), k_)     
    ); 
}
Code:
mut_ = rho_*Cmu_*sqr(k_)/epsilon_;  

#   include "wallViscosityI.H"

    kDiffusion_ = kDiffusion();
    kDiffusion_.updateBoundaryConditions();
mlaurita is offline   Reply With Quote

Old   September 30, 2013, 13:38
Default
  #8
Senior Member
 
Huynh Phong Thanh
Join Date: Aug 2013
Location: Ho Chi Minh City
Posts: 105
Rep Power: 12
hiuluom is on a distinguished road
Hi all.

I have just used openfoam. I has LES model and it runs good. And now, I would like to calculate some terms in LES as first term1 = nu{p*[d2(ui/xi)]} and second term2 = p*[d((uiuj)/xj)]. After that, get average term1 - average term2.
d2 is second derivative and d is first derivative.

Anybody could help me how to setup all math in openfoam?

Thanks.

Thanh.
hiuluom is offline   Reply With Quote

Old   October 1, 2013, 03:33
Default
  #9
New Member
 
Join Date: Jan 2012
Location: Germay - Stuttgart
Posts: 21
Rep Power: 0
HaZe is on a distinguished road
@hiuluom: Would you please start a new thread for your question. Otherwise it's getting complicated with more than one discussion per thread.

@mlaurita: Thank you so much! Now it's working!

At first the compiler complained about "volScalarField" has no member kDiffusion_.updateBoundaryConditions(). But I changed kDiffusion_.updateBoundaryConditions() to kDiffusion_.write and now it's working.

I will do some testing and then post the complete working code.

Thanks for everything and keep FOAMing.
HaZe is offline   Reply With Quote

Old   November 21, 2013, 09:54
Default The working model so far
  #10
New Member
 
Join Date: Jan 2012
Location: Germay - Stuttgart
Posts: 21
Rep Power: 0
HaZe is on a distinguished road
Thank you mlaurita! Here's the extended kEpsilon Model so far. It's working with OF-1.5-dev.
Attached Files
File Type: zip kEpsilonDetailed.zip (10.7 KB, 93 views)
benqing and kevinyou like this.
HaZe is offline   Reply With Quote

Old   December 3, 2013, 23:18
Default how to specify the Boundary Condition for kDiffusion
  #11
New Member
 
zdong
Join Date: Aug 2011
Posts: 15
Rep Power: 14
abbott.hn is on a distinguished road
Quote:
Originally Posted by mlaurita View Post

kDiffusion_
(
IOobject
(
"kDiffusion",
runTime_.timeName(),
mesh_,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh_
)
{
kDiffusion_ = kDiffusion();
kDiffusion_.updateBoundaryConditions();
}
when define the member kDiffusion_ like this ,how should I specify the BC for it? thanks for your reply
abbott.hn is offline   Reply With Quote

Old   December 4, 2013, 03:24
Default
  #12
New Member
 
Join Date: Jan 2012
Location: Germay - Stuttgart
Posts: 21
Rep Power: 0
HaZe is on a distinguished road
Sorry for not mentioning you have to make a new file in your time directory for every new k** and epsilon** field.
The boundary condition is
Code:
type calculated;
for every patch. The dimension for k** is

Code:
dimensions      [1 -1 -3 0 0 0 0];
and for epsilon**
Code:
dimensions      [1 -1 -4 0 0 0 0];
And at last I have to say I didn't take the "updateCoefficients" line into my code. It didn't compile with it, so I don't know yet if it has any effects on the results.

It would be great if anyone using this code would share his/her experiences with it.
HaZe is offline   Reply With Quote

Old   March 6, 2014, 02:40
Default
  #13
New Member
 
yafuji aki
Join Date: Jul 2010
Location: Japan
Posts: 14
Rep Power: 15
aki_yafuji is on a distinguished road
Hi Haze, and everyone!

Thank you for sharing the extended kEpsilon code! It's working with OF-2.2 by modifying a little, and I can output new fields such as kProduction!
I would like to use these new fields in the interFoam solver, how should I do??

In interFoam.C, I add
/////////
volScalarField A1 = turbulence->nut(); // OK
volScalarField A2 = turbulence->k(); // OK
volSymmTensorField A3 = turbulence->R(); // OK
volScalarField A4 = turbulence->kProduction(); // did not work ...
/////////

but, I recieved the following error message. ..

interFoam.C:106:40: error: ‘class Foam::incompressible::turbulenceModel’ has no member named ‘kProduction’

If you wouldn't mind, would you please give me some advice?

aki
aki_yafuji is offline   Reply With Quote

Old   March 6, 2014, 03:28
Default
  #14
New Member
 
Join Date: Jan 2012
Location: Germay - Stuttgart
Posts: 21
Rep Power: 0
HaZe is on a distinguished road
That should be easy. My code works for the compressible Turbulence Models. If you want to use it with incompressible solvers you have to adjust it for the incompressible k and epsilon function and add it to the incompressible Turbulence Models, too.

Meanwhile I figured out the dimensions of the k** and Epsilon** terms are scaled by density (look at the dimensions). I don't know why and where this happens. I already completed my studies and it was possible for me to fix this bug post processing.

If you want to use the fields directly in the solver, don't call the new functions from the solver. I tried this and didn't make it this way. For me it was more elegant to get the values from the new calculated fields.
This works something like:
Code:
volScalarField& kProductionPtr = const_cast<volScalarField&>(mesh.lookupObject<volScalarField>("kProduction"));
HaZe is offline   Reply With Quote

Old   November 24, 2014, 12:50
Default create turbulence fields for additional terms of k and epsilon
  #15
New Member
 
Nitin
Join Date: Mar 2012
Location: Bombay
Posts: 16
Rep Power: 14
Nitin Minocha is on a distinguished road
Hello
First of all, i would like to thank you for this post. I tried to use the modified kEpsilonDetailed model, in order to solve additional terms of turbulent kinetic energy and dissipation rate. I am using OpenFOAM 2.2.0. I am using buoyantBoussinesqPisoFoam solver for incompressible fluid. I got the following error message while compiling the code.
KEmake: Warning: File `options' has modification time 1.8e+03 s in the future
make: warning: Clock skew detected. Your build may be incomplete.
make: Warning: File `linux64Gcc46DPOpt/options' has modification time 15 s in the future
make: warning: Clock skew detected. Your build may be incomplete.
make: Warning: File `Make/linux64Gcc46DPOpt/dontIncludeDeps' has modification time 15 s in the future
make: *** No rule to make target `kEpsilonDetailed/kEpsilonDetailed.dep', needed by `Make/linux64Gcc46DPOpt/dependencies'. Stop.
I am sending you my .c and .h files alongwith make directory. Please have a look and help me solve this problem.

Last edited by Nitin Minocha; August 13, 2016 at 08:10.
Nitin Minocha is offline   Reply With Quote

Old   November 24, 2014, 14:51
Default
  #16
New Member
 
Join Date: Jan 2012
Location: Germay - Stuttgart
Posts: 21
Rep Power: 0
HaZe is on a distinguished road
Quote:
Originally Posted by Nitin Minocha View Post
Hello
First of all, i would like to thank you for this post. I tried to use the modified kEpsilonDetailed model, in order to solve additional terms of turbulent kinetic energy and dissipation rate. I am using OpenFOAM 2.2.0. I am using buoyantBoussinesqPisoFoam solver for incompressible fluid. I got the following error message while compiling the code.
KEmake: Warning: File `options' has modification time 1.8e+03 s in the future
make: warning: Clock skew detected. Your build may be incomplete.
make: Warning: File `linux64Gcc46DPOpt/options' has modification time 15 s in the future
make: warning: Clock skew detected. Your build may be incomplete.
make: Warning: File `Make/linux64Gcc46DPOpt/dontIncludeDeps' has modification time 15 s in the future
make: *** No rule to make target `kEpsilonDetailed/kEpsilonDetailed.dep', needed by `Make/linux64Gcc46DPOpt/dependencies'. Stop.
I am sending you my .c and .h files alongwith make directory. Please have a look and help me solve this problem.
Hi there,
it's great to hear somebody is using my work. Since I didn't work with OF for the last months I'm not sure if I can provide any helpful hints.

But there are some points you can check before further problems occur:
- Did you recognise my model is for compressible fluids?
- Did you recognise my model is made for OF-1.5-dev (I don't know what changes are made between OF-1.5-dev and your OF-version but I'm sure it's a lot!)
- Were you able to compile my model without problems and use it?
- Did you read my last post? In my special case I was able to correct my calculations post processing - but that's no nice way to do so.
- At last I didn't understand why you are trying to extend my model and not extending the existing one from OF-2.2.0. What I did is only to write some new fields to visualize all terms of the k-Epsilon-Model - nothing more.

And last of all would you please start a new thread with your questions.

Regards,
HaZe
HaZe is offline   Reply With Quote

Reply

Thread Tools Search this Thread
Search this Thread:

Advanced Search
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Implement an les turbulence model pante OpenFOAM Programming & Development 19 December 5, 2014 17:16
Compressible turbulence model issues 351Cleveland OpenFOAM 5 October 24, 2013 16:41
Water subcooled boiling Attesz CFX 7 January 5, 2013 04:32
Low Reynolds k-epsilon model YJZ ANSYS 1 August 20, 2010 14:57
build your own turbulence model with buoyancy Thomas Baumann OpenFOAM 11 November 23, 2009 09:53


All times are GMT -4. The time now is 18:33.