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

Adding temperature to pisoFoam

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

Like Tree8Likes
  • 2 Post By Red Devil
  • 2 Post By alexeym
  • 1 Post By alexeym
  • 2 Post By alexeym
  • 1 Post By jherb

Reply
 
LinkBack Thread Tools Display Modes
Old   April 17, 2014, 14:31
Default Adding temperature to pisoFoam
  #1
New Member
 
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3
Red Devil is on a distinguished road
Hi All,

If someone has successfully added temperature to pisoFoam solver, I request you to share your Solver folder (one with *.C and createFields.H files) and any one of your case folders (with '0', 'Constant' and 'System' folders) which you have successfully run.

I have successfully added temperature to icoFoam using the instructions given here:
HTML Code:
http://openfoamwiki.net/index.php/How_to_add_temperature_to_icoFoam
.
Using the same procedure, I added temperature to pisoFoam also and the solver runs fine. I just want to confirm the correctness of my solver. Please provide some help.

Thanks in advance.
ashvinc9 and y_jiang like this.
Red Devil is offline   Reply With Quote

Old   April 17, 2014, 23:41
Default
  #2
Senior Member
 
Dongyue Li
Join Date: Jun 2012
Location: Torino, Italy
Posts: 667
Rep Power: 8
sharonyue is on a distinguished road
Createfields
Code:

    Info<< "Reading field p\n" << endl;
    volScalarField p
    (
        IOobject
        (
            "p",
            runTime.timeName(),
            mesh,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh
    );

    Info<< "Reading field U\n" << endl;
    volVectorField U
    (
        IOobject
        (
            "U",
            runTime.timeName(),
            mesh,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh
    );
Info<< "Reading field T\n" <<endl;
volScalarField T
(
IOobject
(
"T",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
   
IOdictionary transportProperties
(
IOobject
(
"transportProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
);

dimensionedScalar DT
(
transportProperties.lookup("DT")
);
#   include "createPhi.H"


    label pRefCell = 0;
    scalar pRefValue = 0.0;
    setRefCell(p, mesh.solutionDict().subDict("PISO"), pRefCell, pRefValue);


    singlePhaseTransportModel laminarTransport(U, phi);

    autoPtr<incompressible::turbulenceModel> turbulence
    (
        incompressible::turbulenceModel::New(U, phi, laminarTransport)
    );
solver:
Code:
int main(int argc, char *argv[])
{
    #include "setRootCase.H"

    #include "createTime.H"
    #include "createMesh.H"
    #include "createFields.H"
    #include "initContinuityErrs.H"

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

    Info<< "\nStarting time loop\n" << endl;

    while (runTime.loop())
    {
        Info<< "Time = " << runTime.timeName() << nl << endl;

        #include "readPISOControls.H"
        #include "CourantNo.H"

        // Pressure-velocity PISO corrector
        {
            // Momentum predictor

            fvVectorMatrix UEqn
            (
                fvm::ddt(U)
              + fvm::div(phi, U)
              + turbulence->divDevReff(U)
            );

            UEqn.relax();

            if (momentumPredictor)
            {
                solve(UEqn == -fvc::grad(p));
            }

            // --- PISO loop

            for (int corr=0; corr<nCorr; corr++)
            {
                volScalarField rAU(1.0/UEqn.A());

                volVectorField HbyA("HbyA", U);
                HbyA = rAU*UEqn.H();
                surfaceScalarField phiHbyA
                (
                    "phiHbyA",
                    (fvc::interpolate(HbyA) & mesh.Sf())
                  + fvc::ddtPhiCorr(rAU, U, phi)
                );

                adjustPhi(phiHbyA, U, p);

                // Non-orthogonal pressure corrector loop
                for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
                {
                    // Pressure corrector

                    fvScalarMatrix pEqn
                    (
                        fvm::laplacian(rAU, p) == fvc::div(phiHbyA)
                    );

                    pEqn.setReference(pRefCell, pRefValue);

                    if
                    (
                        corr == nCorr-1
                     && nonOrth == nNonOrthCorr
                    )
                    {
                        pEqn.solve(mesh.solver("pFinal"));
                    }
                    else
                    {
                        pEqn.solve();
                    }

                    if (nonOrth == nNonOrthCorr)
                    {
                        phi = phiHbyA - pEqn.flux();
                    }
                }

                #include "continuityErrs.H"

                U = HbyA - rAU*fvc::grad(p);
                U.correctBoundaryConditions();
		fvScalarMatrix TEqn
		(
fvm::ddt(T)
+ fvm::div(phi, T)
- fvm::laplacian(DT, T)
);


TEqn.solve();
            }
        }

        turbulence->correct();

        runTime.write();

        Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
            << "  ClockTime = " << runTime.elapsedClockTime() << " s"
            << nl << endl;
    }

    Info<< "End\n" << endl;

    return 0;
}
Is this right?
sharonyue is offline   Reply With Quote

Old   April 18, 2014, 14:31
Default Adding Temperature to pisoFoam
  #3
New Member
 
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3
Red Devil is on a distinguished road
Thanks for posting your files sharonyue.

I have one change in my solver file. I kept the TEqn outside the piso loop. Which one is correct?

Also, I've seen many people using kappat while adding temperature to pisoFoam. Please help me to understand kappat and it's significance
Red Devil is offline   Reply With Quote

Old   April 18, 2014, 15:48
Default
  #4
Senior Member
 
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19
alexeym will become famous soon enoughalexeym will become famous soon enough
Hi,

would you like to reimplement buoyantBoussinesqPimpleFoam (with thermal expansion coefficient equal to zero, as you did not modify momentum equations)? In that solver temperature equation is solver just after velocity equation.

UPD. kappat is turbulent heat conductivity. If you're planning to simulate heat transfer in turbulent flow, then you should account for it, otherwise, you can forget about it.
alexeym is offline   Reply With Quote

Old   April 18, 2014, 21:05
Default
  #5
Senior Member
 
Dongyue Li
Join Date: Jun 2012
Location: Torino, Italy
Posts: 667
Rep Power: 8
sharonyue is on a distinguished road
Quote:
Originally Posted by alexeym View Post
Hi,

would you like to reimplement buoyantBoussinesqPimpleFoam (with thermal expansion coefficient equal to zero, as you did not modify momentum equations)? In that solver temperature equation is solver just after velocity equation.

UPD. kappat is turbulent heat conductivity. If you're planning to simulate heat transfer in turbulent flow, then you should account for it, otherwise, you can forget about it.
Absolutely rite! I once made a comparison with this TpisoFoam and buoyantBoussinesqPimpleFoam with CFX, the latter is more accurate. Maybe we should consider this kappat into a turbulent solver.

And whether to put it in or out piso loop(after turbulence correction). I just made a comparison too. It is totally the same without any difference... See the picture. So I think put it out of the loop is better to save more calculating time.

Regards,
Attached Images
File Type: jpg QQ??20140419090314.jpg (65.3 KB, 14 views)
sharonyue is offline   Reply With Quote

Old   April 19, 2014, 01:15
Default
  #6
New Member
 
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3
Red Devil is on a distinguished road
Thanks alexeym and sharonyue for your inputs.
I'm trying to create a turbulent solver to solve forced convection heat transfer problem. What all changes I have to make for that?
Red Devil is offline   Reply With Quote

Old   April 19, 2014, 01:41
Default
  #7
Senior Member
 
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19
alexeym will become famous soon enoughalexeym will become famous soon enough
Hi,

if you'd like to create this solver for the sake of creation of the solver, then you should also add buoyancy force to momentum equation.

If you'd like to have this solver for practical reasons (running simulations, obtaining results) then you can just use buoyantBousinesqPimpleFoam. And if you like it to behave like "pisoFoam with temperature and buoyancy" set number of outer correctors (nOuterCorrectors in PIMPLE dictionary) to 1.
alexeym is offline   Reply With Quote

Old   April 19, 2014, 01:53
Default
  #8
New Member
 
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3
Red Devil is on a distinguished road
Hi alexeym,
I'm not considering buoyancy in my problem.
Red Devil is offline   Reply With Quote

Old   April 19, 2014, 02:00
Default
  #9
Senior Member
 
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19
alexeym will become famous soon enoughalexeym will become famous soon enough
OK. Set beta (thermal expansion coefficient) to zero and there will be no buoyancy in buoyantBoussinesqPimpleFoam.

Or if you'd like to go with "your own" solver, here is how turbulent heat diffusivity is calculated:

Code:
    alphat = turbulence->nut()/Prt;
    alphat.correctBoundaryConditions();

    volScalarField alphaEff("alphaEff", turbulence->nu()/Pr + alphat);
then alphaEff is used in laplacian of temperature equation.
sharonyue and Red Devil like this.
alexeym is offline   Reply With Quote

Old   April 19, 2014, 02:05
Default
  #10
New Member
 
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3
Red Devil is on a distinguished road
Thanks a lot alexeym. Will try that.
Red Devil is offline   Reply With Quote

Old   April 19, 2014, 15:32
Default
  #11
New Member
 
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3
Red Devil is on a distinguished road
Hi alexeym,
I've successfully added temperature to pisoFoam by following your inputs. Thanks a lot.

I've another question. I want to set a constant heat flux (q w/m2) boundary condition on a wall. I did it this way:
Code:
    wall
    {
        type      fixedGradient;
        gradient      uniform 100;
    }
where the uniform value of temperature gradient is calculated as q/k (k-thermal conductivity).
Is this correct?
Red Devil is offline   Reply With Quote

Old   April 19, 2014, 16:11
Default
  #12
Senior Member
 
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19
alexeym will become famous soon enoughalexeym will become famous soon enough
Hi,

no, it won't be correct, you need heat transfer coefficient (http://en.wikipedia.org/wiki/Heat_transfer_coefficient).
alexeym is offline   Reply With Quote

Old   April 20, 2014, 02:38
Default
  #13
New Member
 
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3
Red Devil is on a distinguished road
Oh! Then how to implement this constant heat flux bc? I'm totally confused.
Red Devil is offline   Reply With Quote

Old   April 20, 2014, 05:15
Default
  #14
Senior Member
 
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19
alexeym will become famous soon enoughalexeym will become famous soon enough
Well...

start with description of the problem, as "...forced convection heat transfer problem..." is not enough (you'd like to use q = hc*(T - Tinf) relation for calculation of heat flux, it OK, but to calculate that hc one needs to know more about problem).

Then

Code:
$ cd $FOAM_SRC
$ find . -iname '*wallHeat*'
this will give you possible starting points for implementation of your BC.
Red Devil likes this.
alexeym is offline   Reply With Quote

Old   April 20, 2014, 06:58
Default
  #15
New Member
 
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3
Red Devil is on a distinguished road
Problem is flow of water through a square duct. At one wall of the duct, I want to give constant heat flux. Want to run both laminar and turbulent cases. For laminar case i'm using modified icoFoam and for turbulent, modified pisoFoam.

So I have to find h for my problem, then from q=h*(T-Tinf) I can get T.
Then the temperature bc is a uniform fixedValue? Is that the way?
Red Devil is offline   Reply With Quote

Old   April 20, 2014, 12:28
Default
  #16
Senior Member
 
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19
alexeym will become famous soon enoughalexeym will become famous soon enough
Quote:
Originally Posted by Red Devil View Post
Problem is flow of water through a square duct. At one wall of the duct, I want to give constant heat flux. Want to run both laminar and turbulent cases. For laminar case i'm using modified icoFoam and for turbulent, modified pisoFoam.
If you just go with buoyantBoussinesqPimpleFoam (with beta equal zero cause you don't need buoyancy), you can avoid this zoo of solvers, as it is capable of simulating momentum and heat transfer in laminar and turbulent regimes.

Quote:
So I have to find h for my problem, then from q=h*(T-Tinf) I can get T.
Then the temperature bc is a uniform fixedValue? Is that the way?
Actually everything is described in the wikipedia article. But if you'd like to see the math...

1. You've got heat transfer from fluid into a wall and them heat transfer inside the wall. And three temperatures T_fluid, T_wall (temperature of a wall near the fluid interface), T_inf (temperature of the wall at the infinity). You'd like to estimate hc of the system, so you can write q = hc*(T_fluid - T_inf)

2. 1/hc = 1/hc_flud->wall + 1/hc_wall. As we have got fluid, we can assume ideal contact with the wall, so 1/hc_fluid->wall will be equal 0 (so T_wall = T_fluid).

3. Heat transfer coefficient inside the wall can be estimated as kappa_wall/L_wall, where kappa_wall is heat conductivity of the wall material, L_wall is the thickness of the wall.

And finally if your wall thickness is zero, then yes, BC should be fixed value with T_inf.
sharonyue and Red Devil like this.
alexeym is offline   Reply With Quote

Old   April 22, 2014, 17:50
Default
  #17
Senior Member
 
Joachim Herb
Join Date: Sep 2010
Posts: 289
Rep Power: 9
jherb is on a distinguished road
I think you are looking for the turbulentHeatFluxTemperature boundary condition:
https://github.com/OpenFOAM/OpenFOAM...hScalarField.H
Red Devil likes this.
jherb is offline   Reply With Quote

Reply

Tags
pisofoam, temperature.add

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
Problems adding temperature to icoFoam sur4j OpenFOAM Programming & Development 8 January 11, 2015 02:14
chtMultiRegionSimpleFoam: strange error samiam1000 OpenFOAM Running, Solving & CFD 25 August 29, 2014 17:28
pisoFoam LES turbulent model + temperature Byxon OpenFOAM Programming & Development 0 October 4, 2013 06:55
Adding temperature equation to porousSimpleFoam David_010 OpenFOAM Programming & Development 7 February 18, 2013 09:42
Adding Temperature in channelFoam dhruv OpenFOAM 10 May 5, 2012 05:40


All times are GMT -4. The time now is 07:22.