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

TEqn in compressibleInterFoam

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

Like Tree1Likes
  • 1 Post By sabrahams

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   August 9, 2023, 10:24
Question TEqn in compressibleInterFoam
  #1
New Member
 
S Abrahams
Join Date: Mar 2022
Location: UK
Posts: 14
Rep Power: 4
sabrahams is on a distinguished road
Hi everyone!

I wonder if anyone can help me to understand the formulation of the TEqn in compressibleInterFoam in OF9. The TEqn is:

Code:
fvm::ddt(rho, T) + fvm::div(rhoPhi, T) - fvm::Sp(contErr, T)
      - fvm::laplacian(turbulence.alphaEff(), T)
      + (
            fvc::div(fvc::absolute(phi, U), p)()() // - contErr/rho*p
          + (fvc::ddt(rho, K) + fvc::div(rhoPhi, K))()()
          - (U()&(fvModels.source(rho, U)&U)()) - contErr*K 
        )
       *(
           alpha1()/mixture.thermo1().Cv()()
         + alpha2()/mixture.thermo2().Cv()()
        )
        
     ==
           fvModels.source(rho, T)
I believe C_v should be inside the derivatives since C_v is averaged over the two phases and can therefore vary in time and space.

I believe the TEqn should be equivalent to
\frac{\partial}{\partial t}(\rho C_v T)+\nabla\cdot(\rho C_v U  T) + \frac{\partial}{\partial t}(\rho k)+\nabla\cdot(\rho U k)+\nabla  \cdot (pU)-\nabla\cdot (K \nabla T)=S,

but the formulation above seems to be equivalent to
\frac{C_{1}C_{2}}{{\alpha_1 C_{2}+\alpha_2 C_{1}}}  \left[\frac{\partial}{\partial t}(\rho T)+\nabla\cdot(\rho  U T)  \right]+  \frac{\partial}{\partial t}(\rho k)+\nabla\cdot(\rho U  k)+\nabla \cdot  (pU)-\nabla\cdot (K \nabla T)=S.

Can anyone help me to understand how C_v has been taken out in this way?

Many thanks,
SAbrahams

Last edited by sabrahams; August 9, 2023 at 13:09. Reason: typos
sabrahams is offline   Reply With Quote

Old   August 15, 2023, 09:34
Default
  #2
New Member
 
S Abrahams
Join Date: Mar 2022
Location: UK
Posts: 14
Rep Power: 4
sabrahams is on a distinguished road
Does no-one know the answer to this?
sabrahams is offline   Reply With Quote

Old   August 15, 2023, 10:55
Default
  #3
Senior Member
 
Join Date: Apr 2020
Location: UK
Posts: 670
Rep Power: 14
Tobermory will become famous soon enough
I think the problem lies with your first Teqn equation ... Start from basics - the underlying energy equation is stated in terms of the internal energy, i.e. in conservative form:

\frac{\partial \rho e}{\partial t} + \frac{\partial \rho U_i e}{\partial x_i} + ...

which we can expand out using continuity to the non-conservative form:

\rho \frac{\partial e}{\partial t} + \rho U_i \frac{\partial  e}{\partial x_i} + ...
Now using the definition of the specific heat, we write de = C_v dT, so substitute and you have:

\rho C_v \frac{\partial T}{\partial t} + \rho C_v U_i \frac{\partial  T}{\partial x_i} + ...

or in terms of the material derivative:
\rho C_v \frac{DT}{Dt} + ...

with the specific heat capacity on the outside. Hopefully that helps.
Tobermory is offline   Reply With Quote

Old   August 15, 2023, 11:22
Default
  #4
New Member
 
S Abrahams
Join Date: Mar 2022
Location: UK
Posts: 14
Rep Power: 4
sabrahams is on a distinguished road
Hi Tobermory,

Thanks for your reply. However, since C_v is not a constant in this case, I don't believe we can say de = C_v dT. Instead, it should be de = C_v dT + T dC_v.

In a two phase system, as we have in compressibleInterFoam, C_v is averaged over the two phases, that is
C_v = \alpha_1 C_{v,1}+\alpha_2 C_{v,2}.

While C_{v,1} and C_{v,2} may be constant, C_v is not since \alpha_1 and \alpha_2 vary in time and space.

Furthermore, TEqn in compressibleInterFoam (as shown in my original post) is not
\rho C_v \frac{\partial T}{\partial t}+ ...
it seems to be (if you divide through by \alpha_1/C_{v,1}+\alpha_2/C_{v,2}, which is multiplying the second term in TEqn)
\frac{C_{v,1} C_{v,2}}{\alpha_1 C_{v,2}+\alpha_2 C_{v,1}} \frac{\partial \rho T}{\partial t}+ ...

So I'm trying to understand how the energy equation has been rearranged to get the form in the TEqn as in my original post. Do you know how this is derived?


Thanks again,
SAbrahams

Last edited by sabrahams; August 15, 2023 at 11:28. Reason: correcting mistake in equation
sabrahams is offline   Reply With Quote

Old   August 16, 2023, 11:32
Default
  #5
Senior Member
 
Join Date: Apr 2020
Location: UK
Posts: 670
Rep Power: 14
Tobermory will become famous soon enough
Okay - let me break this down. First up, the following is not correct I am afraid:

Quote:
Originally Posted by sabrahams View Post
However, since C_v is not a constant in this case, I don't believe we can say de = C_v dT. Instead, it should be de = C_v dT + T dC_v.
You need to start with the definition of the specific heat capacity, which is C_v=(\frac{\partial e}{\partial T})_V ... and now can you see where I got my original de = C_v dT from? We have said nothing so far about the variability of Cv in space, just that the process is assumed to be a constant volume process (which makes sense for two immiscible fluids). If you're still not convinced, then think of the chain rule:

\frac{\partial e}{\partial t} = \frac{\partial e}{\partial T} \frac{\partial T}{\partial t} = C_v \frac{\partial T}{\partial t}

etc. Now consider the single phase energy equation, which in OpenFOAM is (refer to buoyantPimpleFoam for example):

\frac{\partial \rho e}{\partial t} + \frac{\partial \rho U_j e}{\partial x_j} +
\frac{\partial \rho K}{\partial t} + \frac{\partial \rho U_j K}{\partial x_j} +
\frac{\partial U_j p}{\partial x_j} 
- \frac{\partial}{\partial x_j}(\rho\alpha_{eff} \frac{\partial e}{\partial x_j})

expand out the first two terms into the nonconservative form using continuity, substitute for \partial e and convert back to conservative form and you get:

\frac{\partial \rho T}{\partial t} + \frac{\partial \rho U_j T}{\partial x_j} +
\frac{1}{C_v}[\frac{\partial \rho K}{\partial t} + \frac{\partial \rho U_j K}{\partial x_j} + \frac{\partial U_j p}{\partial x_j}]
- \frac{1}{C_v}\frac{\partial}{\partial x_j}(\rho C_v \alpha_{eff} \frac{\partial T}{\partial x_j})

NOW a minor approximation is made in the laplacian, where it assumed that local variations in Cv are small, so that the terms cancel leaving:

\frac{\partial \rho T}{\partial t} + \frac{\partial \rho U_j T}{\partial x_j} +
\frac{1}{C_v}[\frac{\partial \rho K}{\partial t} + \frac{\partial \rho U_j K}{\partial x_j} + \frac{\partial U_j p}{\partial x_j}]
- \frac{\partial}{\partial x_j}(\rho \alpha_{eff} \frac{\partial T}{\partial x_j})

and we are almost there. All that is left is to realise that the finite volume discretised form of the energy (or temperature) budget equation is the sum of alpha1 times the equation for phase 1 plus alpha2 times that for phase two, which results in:

(\alpha_1 + \alpha_2)[\frac{\partial \rho T}{\partial t} + \frac{\partial \rho U_j T}{\partial x_j}]+ (\frac{\alpha_1}{C_{v1}}+\frac{\alpha_2}{C_{v2}})[\frac{\partial \rho K}{\partial t} + \frac{\partial \rho U_j K}{\partial x_j} + \frac{\partial U_j p}{\partial x_j}]
- (\alpha_1 + \alpha_2)\frac{\partial}{\partial x_j}(\rho \alpha_{eff} \frac{\partial T}{\partial x_j})

or on noting that \alpha_1 + \alpha_2=1, we have:

\frac{\partial \rho T}{\partial t} + \frac{\partial \rho U_j T}{\partial x_j} +
(\frac{\alpha_1}{C_{v1}}+\frac{\alpha_2}{C_{v2}})[\frac{\partial \rho K}{\partial t} + \frac{\partial \rho U_j K}{\partial x_j} + \frac{\partial U_j p}{\partial x_j}]
- \frac{\partial}{\partial x_j}(\rho \alpha_{eff} \frac{\partial T}{\partial x_j})

noting that these should now not be differentials etc., but should instead be volume integrals of the differentials yielding volume average values or surface fluxes (I didn't have the time to tidy up the terminology). This is the form that appears in TEqn.h for compressibleInterFoam.

There are clearly a bunch of other implicit assumptions here ... but this is my reverse engineering of the code. I couldn't find any explicit reference to the source of the model equations, and suspect it was authored by one of the OpenFOAM team. If you do find a reference - share it with the forum please.

Last edited by Tobermory; August 16, 2023 at 12:36.
Tobermory is offline   Reply With Quote

Old   August 17, 2023, 10:38
Default
  #6
New Member
 
S Abrahams
Join Date: Mar 2022
Location: UK
Posts: 14
Rep Power: 4
sabrahams is on a distinguished road
That makes perfect sense. Thank you for the clear explanation!


I haven't found any references for the equations used in compressibleInterFoam in OF9. If I find anything I'll certainly share. UEqn and alphaEqn appear to be in the same form as equations (7) and (25) of Shi et al.. Although, this paper was published later than OF9 so isn't a reference but just a useful resource to understand the equations (for anyone else reading this).


Thanks again Tobermory.
Tobermory likes this.
sabrahams is offline   Reply With Quote

Old   August 23, 2023, 10:36
Default
  #7
New Member
 
S Abrahams
Join Date: Mar 2022
Location: UK
Posts: 14
Rep Power: 4
sabrahams is on a distinguished road
Hi again Tobermory,

I've just been working through this and I'm stuck on this step:
Quote:
Originally Posted by Tobermory View Post
All that is left is to realise that the finite volume discretised form of the energy (or temperature) budget equation is the sum of alpha1 times the equation for phase 1 plus alpha2 times that for phase two, which results in:

(\alpha_1 + \alpha_2)[\frac{\partial \rho T}{\partial t} + \frac{\partial \rho U_j T}{\partial x_j}]+ (\frac{\alpha_1}{C_{v1}}+\frac{\alpha_2}{C_{v2}})[\frac{\partial \rho K}{\partial t} + \frac{\partial \rho U_j K}{\partial x_j} + \frac{\partial U_j p}{\partial x_j}]
- (\alpha_1 + \alpha_2)\frac{\partial}{\partial x_j}(\rho \alpha_{eff} \frac{\partial T}{\partial x_j})
I believe that the equation for phase 1 would have \rho_1 and the equation for phase 2 would have \rho_2 so I'm struggling to see how the two equations could be summed in the way that you've shown above. I'm getting an equation of the form:


\alpha_1 [\frac{\partial \rho_1 T}{\partial t} +  \frac{\partial \rho_1 U_j T}{\partial x_j}]+
\alpha_2[\frac{\partial \rho_2 T}{\partial t} +  \frac{\partial \rho_2 U_j T}{\partial x_j}] + ...

Am I missing a step or an assumption somewhere?
sabrahams is offline   Reply With Quote

Old   August 24, 2023, 10:34
Default
  #8
Senior Member
 
Join Date: Apr 2020
Location: UK
Posts: 670
Rep Power: 14
Tobermory will become famous soon enough
The way to think of it, I reckon, is as follows: we need to integrate the differential equations over the cell volume to get the finite volume expressions. Focus just on the first two terms, for simplicity:

\frac{\partial \rho T}{\partial t} + \frac{\partial \rho U_j T}{\partial x_j} + ...

integrate these to get:

\int_V( \frac{\partial \rho T}{\partial t}).dV +
\int_V(\frac{\partial \rho U_j T}{\partial x_j}).dV + ...

\frac{\partial }{\partial t}(\int_V \rho T . dV) +
\int_S(\rho U_j T).dA + ...

\frac{\partial }{\partial t}(\rho_P T_P V_P) +
\sum(\phi_f T_f) + ...

where it's assumed that the mesh is not moving or morphing, \phi_f are the face fluxes, and the P subscript denotes a volume-averaged value, e.g.

\int_V T . dV = T_P V_p

Note though that we are integrating over volume \alpha_1 V_P for phase 1, and \alpha_2 V_P for phase 2. This is the origin of the \alpha_1 and \alpha_2 coefficients in each of the terms. This is what I meant rather laconically in my earlier post by "noting that these should now not be differentials etc.". Hope this helps.
Tobermory is offline   Reply With Quote

Reply

Tags
heat capacity, temperature, teqn

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
Correction of T after TEqn in compressibleInterFoam krikre OpenFOAM Programming & Development 3 June 26, 2022 22:07
TEqn (passive scalar) diverges all of a sudden backscatter OpenFOAM Running, Solving & CFD 11 October 21, 2018 07:56
TEqn. in PorousSimpleFoam svramana Main CFD Forum 0 February 14, 2018 08:53
EEqn to TEqn jishnuhari25 OpenFOAM Programming & Development 1 June 4, 2015 11:06
Replace H equation with Teqn in chtMultiregionFoam OF 2.3 Moncef OpenFOAM Running, Solving & CFD 0 May 22, 2014 20:37


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