
[Sponsors] 
April 17, 2014, 14:31 
Adding temperature to pisoFoam

#1 
New Member
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3 
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. 

April 17, 2014, 23:41 

#2 
Senior Member
Dongyue Li
Join Date: Jun 2012
Location: Torino, Italy
Posts: 667
Rep Power: 8 
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) ); 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" // Pressurevelocity 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); // Nonorthogonal 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 == nCorr1 && 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; } 

April 18, 2014, 14:31 
Adding Temperature to pisoFoam

#3 
New Member
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3 
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 

April 18, 2014, 15:48 

#4 
Senior Member
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19 
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. 

April 18, 2014, 21:05 

#5  
Senior Member
Dongyue Li
Join Date: Jun 2012
Location: Torino, Italy
Posts: 667
Rep Power: 8 
Quote:
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, 

April 19, 2014, 01:41 

#7 
Senior Member
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19 
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. 

April 19, 2014, 02:00 

#9 
Senior Member
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19 
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); 

April 19, 2014, 15:32 

#11 
New Member
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3 
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; } Is this correct? 

April 19, 2014, 16:11 

#12 
Senior Member
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19 
Hi,
no, it won't be correct, you need heat transfer coefficient (http://en.wikipedia.org/wiki/Heat_transfer_coefficient). 

April 20, 2014, 02:38 

#13 
New Member
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3 
Oh! Then how to implement this constant heat flux bc? I'm totally confused.


April 20, 2014, 05:15 

#14 
Senior Member
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19 
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*' 

April 20, 2014, 06:58 

#15 
New Member
RD
Join Date: Oct 2013
Location: UK
Posts: 18
Rep Power: 3 
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*(TTinf) I can get T. Then the temperature bc is a uniform fixedValue? Is that the way? 

April 20, 2014, 12:28 

#16  
Senior Member
Alexey Matveichev
Join Date: Aug 2011
Location: Nancy, France
Posts: 1,089
Rep Power: 19 
Quote:
Quote:
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. 

April 22, 2014, 17:50 

#17 
Senior Member
Joachim Herb
Join Date: Sep 2010
Posts: 289
Rep Power: 9 
I think you are looking for the turbulentHeatFluxTemperature boundary condition:
https://github.com/OpenFOAM/OpenFOAM...hScalarField.H 

Tags 
pisofoam, temperature.add 
Thread Tools  
Display Modes  


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 