CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > OpenFOAM Running, Solving & CFD

thermophysicalProperties with temperature dependance

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

Like Tree1Likes
  • 1 Post By olivierG

Reply
 
LinkBack Thread Tools Display Modes
Old   December 14, 2015, 05:06
Default thermophysicalProperties with temperature dependance
  #1
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 3
hcl734 is on a distinguished road
Hi all,

is it possible to set upper and lower limits for let's say density when using a polynomial approach?
I want to establish something like that:

Density as a function of temperature
rho_min: 1000
T < 240: 1100
T>240 : C1* T^2 + C2 *T + C3

It would also be possible to use step functions, if available.


Best
hcl734 is offline   Reply With Quote

Old   December 14, 2015, 05:34
Default
  #2
Senior Member
 
Olivier
Join Date: Jun 2009
Location: France, grenoble
Posts: 266
Rep Power: 10
olivierG is on a distinguished road
hello,
This is not a feature in OF, but it's easy to add this : check how to modify the polynomial density, and add your limiter.
Take a loot at :http://www.tfd.chalmers.se/~hani/kur...nFoam%20v2.pdf

This is for viscosity, but do the same for density...
The basic step are:
1) copy $FOAM_SRC/thermophysicalModels/specie/equationOfState/icoPolynomial/ where you want.
2) rename all "icoPolynomial" by what you want like "icoLimitPolynomail"(file name AND inside !)
3) modify in file to add what you want
4) add a make file (make sur lib will be in FOAM_USER_LIBBIN)
4) compile (wmake) and use it ! (do not forget the libname in you controlDict)


regards,
olivier
olivierG is offline   Reply With Quote

Old   December 14, 2015, 06:00
Default
  #3
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 3
hcl734 is on a distinguished road
Ok found the solution

Use icoPolynomial for density as explained here http://cfd.direct/openfoam/user-guide/thermophysical/

and set the upper and lower limits for rho in the fvSolution
hcl734 is offline   Reply With Quote

Old   December 14, 2015, 07:23
Default
  #4
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 3
hcl734 is on a distinguished road
Thanks for showing me that paper, but to be honest I'm not very familiar with OF programming and I found this simple solution for density which works good for me.

BUT I want to do the same thing for heat capacit in chtMultiRegionSimpleFoam-solver.
For density the feature is implemented in the pressure equation as

Quote:
rho = thermo.rho();
rho = max(rho, rhoMin[i]);
rho = min(rho, rhoMax[i]);
rho.relax();
And an entry in createFluidsField.H
Code:
    PtrList<dimensionedScalar> CpMax(fluidRegions.size());
    PtrList<dimensionedScalar> CpMin(fluidRegions.size());


        rhoMax.set(i, new dimensionedScalar(simpleDict.lookup("rhoMax")));
        rhoMin.set(i, new dimensionedScalar(simpleDict.lookup("rhoMin")));
And added in createFluidsField.H


Code:
        CpMax.set(i, new dimensionedScalar(simpleDict.lookup("CpMax")));
        CpMin.set(i, new dimensionedScalar(simpleDict.lookup("CpMin")));
But when I take a look at the EEqn.H

Code:
{
    volScalarField& he = thermo.he();

    fvScalarMatrix EEqn
    (
        fvm::div(phi, he)
      + (
            he.name() == "e"
          ? fvc::div(phi, volScalarField("Ekp", 0.5*magSqr(U) + p/rho))
          : fvc::div(phi, volScalarField("K", 0.5*magSqr(U)))
        )
      - fvm::laplacian(turb.alphaEff(), he)
     ==
        rad.Sh(thermo)
      + fvOptions(rho, he)
    );

    EEqn.relax();

    fvOptions.constrain(EEqn);

    EEqn.solve();

    fvOptions.correct(he);

    thermo.correct();
    rad.correct();

    Info<< "Min/max T:" << min(thermo.T()).value() << ' '
        << max(thermo.T()).value() << endl;
}
I can't even see where they are using specific heat capacity?
Could you help me with this one?
It would also be nice to write the Cp-field while solving to check whether it's correct or not.
hcl734 is offline   Reply With Quote

Old   December 14, 2015, 08:42
Default Viscosity depending on temperature exponentially
  #5
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 3
hcl734 is on a distinguished road
2) Question

Is it also possible to set an exponential dependance between viscosity and temperature?
Something like:

(C1 + C2*exp(-(T-Tref)/C5) + C3*exp(-(T-Tref)/C6) + C4*exp(-(T-Tref)/C7)
hcl734 is offline   Reply With Quote

Old   December 14, 2015, 09:15
Default
  #6
Senior Member
 
Olivier
Join Date: Jun 2009
Location: France, grenoble
Posts: 266
Rep Power: 10
olivierG is on a distinguished road
hello,

You are wrong about the limit in fvSolution => this bound rho, but NOT in a conservative way. This is only for stability purpose. You should use a correct model for density.

About 1): Cp is inside alphaEff() ( ~ something like k/Cp + kt/Cp)
You can use polynomal for Cp => use janaf model

About 2) by defaut, OF has only Sutherland viscosity(T) for gas. Take a look at the link i give you: this is how to add a custom viscosity.

BTW, you should not be afraid about coding with OF: this is the way do do.

regards,
olivier
olivierG is offline   Reply With Quote

Old   December 14, 2015, 12:44
Default
  #7
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 3
hcl734 is on a distinguished road
Thanks for your advice, I followed the tutorial you showed me for implementing a temperature dependent viscosity, but this tutorial is only for incompressible flow using a viscosityProperties dictionary. But I want to calculate compressible using thermophysicalProperties.

Code:
thermoType
{
    type            heRhoThermo;
    mixture         pureMixture;
    transport       polynomial;
    thermo          hPolynomial;
    equationOfState icoPolynomial;
    specie          specie;
    energy          sensibleEnthalpy;
}
with chtMultiRegionSimpleFoam



When I look at heRhoThermo.C

Code:
#include "heRhoThermo.H"

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

template<class BasicPsiThermo, class MixtureType>
void Foam::heRhoThermo<BasicPsiThermo, MixtureType>::calculate()
{
    const scalarField& hCells = this->he().internalField();
    const scalarField& pCells = this->p_.internalField();

    scalarField& TCells = this->T_.internalField();
    scalarField& psiCells = this->psi_.internalField();
    scalarField& rhoCells = this->rho_.internalField();
    scalarField& muCells = this->mu_.internalField();
    scalarField& alphaCells = this->alpha_.internalField();

    forAll(TCells, celli)
    {
        const typename MixtureType::thermoType& mixture_ =
            this->cellMixture(celli);

        TCells[celli] = mixture_.THE
        (
            hCells[celli],
            pCells[celli],
            TCells[celli]
        );

        psiCells[celli] = mixture_.psi(pCells[celli], TCells[celli]);
        rhoCells[celli] = mixture_.rho(pCells[celli], TCells[celli]);

        muCells[celli] = mixture_.mu(pCells[celli], TCells[celli]);
        alphaCells[celli] = mixture_.alphah(pCells[celli], TCells[celli]);
    }

    forAll(this->T_.boundaryField(), patchi)
    {
        fvPatchScalarField& pp = this->p_.boundaryField()[patchi];
        fvPatchScalarField& pT = this->T_.boundaryField()[patchi];
        fvPatchScalarField& ppsi = this->psi_.boundaryField()[patchi];
        fvPatchScalarField& prho = this->rho_.boundaryField()[patchi];

        fvPatchScalarField& ph = this->he().boundaryField()[patchi];

        fvPatchScalarField& pmu = this->mu_.boundaryField()[patchi];
        fvPatchScalarField& palpha = this->alpha_.boundaryField()[patchi];

        if (pT.fixesValue())
        {
            forAll(pT, facei)
            {
                const typename MixtureType::thermoType& mixture_ =
                    this->patchFaceMixture(patchi, facei);

                ph[facei] = mixture_.HE(pp[facei], pT[facei]);

                ppsi[facei] = mixture_.psi(pp[facei], pT[facei]);
                prho[facei] = mixture_.rho(pp[facei], pT[facei]);
                pmu[facei] = mixture_.mu(pp[facei], pT[facei]);
                palpha[facei] = mixture_.alphah(pp[facei], pT[facei]);
            }
        }
        else
        {
            forAll(pT, facei)
            {
                const typename MixtureType::thermoType& mixture_ =
                    this->patchFaceMixture(patchi, facei);

                pT[facei] = mixture_.THE(ph[facei], pp[facei], pT[facei]);

                ppsi[facei] = mixture_.psi(pp[facei], pT[facei]);
                prho[facei] = mixture_.rho(pp[facei], pT[facei]);
                pmu[facei] = mixture_.mu(pp[facei], pT[facei]);
                palpha[facei] = mixture_.alphah(pp[facei], pT[facei]);
            }
        }
    }
}


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

template<class BasicPsiThermo, class MixtureType>
Foam::heRhoThermo<BasicPsiThermo, MixtureType>::heRhoThermo
(
    const fvMesh& mesh,
    const word& phaseName
)
:
    heThermo<BasicPsiThermo, MixtureType>(mesh, phaseName)
{
    calculate();
}


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

template<class BasicPsiThermo, class MixtureType>
Foam::heRhoThermo<BasicPsiThermo, MixtureType>::~heRhoThermo()
{}


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

template<class BasicPsiThermo, class MixtureType>
void Foam::heRhoThermo<BasicPsiThermo, MixtureType>::correct()
{
    if (debug)
    {
        Info<< "entering heRhoThermo<MixtureType>::correct()" << endl;
    }

    calculate();

    if (debug)
    {
        Info<< "exiting heRhoThermo<MixtureType>::correct()" << endl;
    }
}


// ************************************************************************* //
I can't see much resemblance between this and the viscosityModel files.
Where should I start when modifying the compressible case?
I mean basically I only need an exponential function for defining mu and I need to know how to write the transportProperties as fields.
Thx again for your help!
hcl734 is offline   Reply With Quote

Old   December 15, 2015, 04:21
Default
  #8
Senior Member
 
Olivier
Join Date: Jun 2009
Location: France, grenoble
Posts: 266
Rep Power: 10
olivierG is on a distinguished road
hello,

You should add a custom thermo model.
Check this : How to add a new Equation of State

regards,
olivier
wayne14 likes this.
olivierG is offline   Reply With Quote

Old   December 19, 2015, 14:14
Default
  #9
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 3
hcl734 is on a distinguished road
Ok I checked the post you recommended to me and I am trying to get into programming with OpenFoam.
But since this is quite tedious I am wondering whether it is possible to define a kind of lookup table for a thermophysical propertie in OpenFoam where I can put some data points (calculated externally) and OpenFoam does linear interpolation between those.
I hope you get my idea.
hcl734 is offline   Reply With Quote

Old   February 12, 2016, 08:13
Default
  #10
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 3
hcl734 is on a distinguished road
Didn't checked it but for everybody coming after me

This gives you the possibility to assign exponential terms
https://github.com/OpenFOAM/OpenFOAM...1ec9af6f399210
hcl734 is offline   Reply With Quote

Reply

Thread Tools
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 On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Shadow Wall and temperature norger FLUENT 8 April 22, 2016 10:51
UDF for Back-flow Temperature G340 Fluent UDF and Scheme Programming 3 August 21, 2013 04:56
Inlet won't apply UDF and has temperature at 0K! tccruise Fluent UDF and Scheme Programming 2 September 14, 2012 06:08
where is the calculation of the temperature field Tobi OpenFOAM 1 July 30, 2012 10:40
monitoring point of total temperature rogbrito FLUENT 0 June 21, 2009 17:31


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