[ovie]

That was my initial approach but of course the result didnt converge. Now I dont know if the lack of convergence resulted from the incorrect approximation of the divergence term or the unsteady term or the diffusion term in the heat equation. I would have to do one at a time to see which actually caused the problem.

Again, reading through Edin's comments, he seems to argue that its not correct to apply the same linear interpolation to kappa as is done to density. So there are two issues here, one is the use of face centered values as against the cell centered values and the method of computing the "mixture" average: linear or harmonic interpolation.

[feijooos]

I don't see why it's not possible to do a simple linear interpolation, since the calculation is based on gamma, which has already been calculated and corrected.

Edin, would you mind explaining why it is not possible to do the same interpolation for kappa as is done with rhoCp?

Thanks,

Eelco

[chiven]

Hi, I meet a problem. That is I want to update the rho1 according to the liquid temperature, thus rho1 should be changed to volScalarField from dimensionedScalar. However, my C++ knowledge is limited, and now I have no idea about this. Would you help me? thanks.
Chiven

part of the files of twoPhaseMixture.C and twoPhaseMixture.H are attached.

for the file of twoPhaseMixture.C
rho1_(nuModel1_->viscosityProperties().lookup("rho")),
rho2_(nuModel2_->viscosityProperties().lookup("rho")),

bool twoPhaseMixture::read()
{
if (transportModel::read())
{
if
(
nuModel1_().read(subDict(phase1Name_))
&& nuModel2_().read(subDict(phase2Name_))
)
{
nuModel1_->viscosityProperties().lookup("rho") >> rho1_;
nuModel2_->viscosityProperties().lookup("rho") >> rho2_;
return true;
}
else
{
return false;
}
}
else
{
return false;
}
}

for the file of twoPhaseMixture.C

dimensionedScalar rho1_;
dimensionedScalar rho2_;
//- Return const-access to phase1 density
const dimensionedScalar& rho1() const
{
return rho1_;
}
//- Return const-access to phase2 density
const dimensionedScalar& rho2() const
{
return rho2_;
};

[ovie]

Hi Chiven,

Have you been able to implement the temperature dependent density in your code? Please let me know if you still need any help with that.

Regards

[chiven]

Hi Ovie, thank you, this time I am struggling for the interPhaseChangeFoam with condensation (vapour to water), I have finished the code implementation, but the calculated temperature fields are not right. I am thinking how to improve the codes.
How about your progress?

Best regards,
Chiven

[ovie]

Hi Chiven,

I am also working on phase change but evaporation. In my case the properties of my fluid are all dependent on temperature and another scalar variable (dry solids mass fraction) so I have another scalar transport equation in addition to the energy transport equation.

Right now, I am having two challenges:

1.) I have implemented the temperature dependence of the various fluid transport properties but when I tried to validate the code on a single phase 2-D channel flow the results dont converge. The temperature field looks reasonable but not "smooth" and same goes for velocity. I have included a thread on this at :
http://www.cfd-online.com/Forums/ope...-pisofoam.html

I am still not quite sure how to fix this.

2.) The second issue is how to implement the production and annihilation source terms for the alpha1 equation for a temperature dependent phase change phenomena. I have also included a thread on this at:
http://www.cfd-online.com/Forums/ope...hangefoam.html

These are the tow challenges I am facing.

I am sure though that you have gotten around implementing your viscosity models. If you have any ideas on how to proceed with my two questions, I would very much appreciate it.
Regards and thanks...
ovie

[romant]

Hej,

I am working on the same problem. I have a model for the phase change from vapor to water, already implemented and at that point i am struggling with the implementation for the temperature equation. The first idea was to implement the energy equation only on one side, since I am only interested on the water side at the moment. I am using saturated steam on subcooled water.

Have zou guys come further with implementation for the energy equation in interPhaseChangeFoam?

[ovie]

Hey Romant,

If you are using the single field representation which the interFoam family of solvers is based on, then the energy equation must apply to the entire flow field which is both the liquid and vapour phases of the flow. So I am not quite sure what you meant by implementing the energy equation for only the liquid side as you indicated in your comment.

The implementation of the energy equation is pretty straight forward just like any scalar associated with the flow. However, it must be expressed in conservative form. i.e.
d(rho*cp*T)/dt + div(cp*T*rho*V) = div(kappa*grad T) + source terms

The first term (energy stored) on the LHS of the equation can be rewritten as:
d(rho*cp*T)/dt = d(rhoCp*T)

where rhoCp = "rho*cp".
rhoCp is for the mixture. To evaluate rhoCp, you would go to the alphaEqnSubCycle.H code. At the end, you find the line where rho is updated for the solver. Include a line that computes rhoCp by doing:

rhoCp = alpha1*rho1*cp1 + (scalar(1) - alpha1)*rho2*cp2

The second term (energy flux) on the LHS is rewritten (based on OF terminology):
div(cp*T*rho*V) = div(cp*T*rhophi).

rhophi is the face mass flux. obtained from the alphaEqn.H code.

Again, this can be written as:

div(cp*T*rhophi) = div(rhophiCp*T)

rhophi is updated for the solver in the alphaEqn.H code. To update rhophiCp, go the alphaEqn.H code and just after the line that computes rhophi, include the following line:

rhoPhiCp = phiAlpha*(rho1*cp1 - rho2*cp2) + phi*rho2*cp2;

For the energy diffusion term i.e. first term on the RHS, you have to include a line that computes kappa in your code or read it off the transportProperties dictionary.

The last term on the RHS are the source terms. These are the models which should describe energy consumption when liquid is converted to vapour (depending on the temperature and saturation pressure) and energy released when vapour reverts to liquid. So you would have something like
source terms = (mass_of_vapour_to_liquid - mass_of_liquid_to_vapour) * latent heat of vapourization.

This is how I have implemented my solver so far. However, I am yet to come up with appropriate models for the evaporation and condensation processes.

Have you been able to get past that?

I hope I have been able to help somewhat.

[romant]

Hej,

I think i am gonna try that as soon as possible as you described. at the moment I have considered saturated steam only and direct contact condensation, so there is no evaporation. the model that i am using is something like

m = h * A (Tsat - Tbulk) / latentHeat

where A is modelled based on the alpha per cell with

Amax = sqrt(3) * Vcell^(2/3)

and then the interface becomes Al = Amax * (2*alpha) for 0<alpha<0.5
and Al = Amax * 2 * (1-alpha) for 0.5<alpha<1

for the heat transfer coefficient, h, I used a model based on steam falling onto a flat plate for the moment. of course this model is exchangeable, as well as the interface area. evaporation is set to 0.

this is also the reason I also wanted to model the temperature on only the water side, because the other is supposed to stay constant.

[metro]

Hey guys,

I am new to openfoam and c++ programming and I am trying to implement Ovie's temperature field code. As I am trying to compile it I keep receiving the following error message

"metro@ubuntu:~/OpenFOAM/OpenFOAM-1.6.x/applications/solvers/multiphase/interTempFoam$ wmake
SOURCE=interTempFoam.C ; g++ -m64 -Dlinux64 -DWM_DP -Wall -Wno-strict-aliasing -Wextra -Wno-unused-parameter -Wold-style-cast -Wnon-virtual-dtor -O3 -DNoRepository -ftemplate-depth-40 -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/transportModels -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/transportModels/incompressibleTemp/lnInclude -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/transportModels/interfaceProperties/lnInclude -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/turbulenceModels/incompressible/turbulenceModel -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude -IlnInclude -I. -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/OpenFOAM/lnInclude -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/OSspecific/POSIX/lnInclude -fPIC -c $SOURCE -o Make/linux64GccDPOpt/interTempFoam.o
In file included from interTempFoam.C:58:
createFields.H: In function ‘int main(int, char**)’:
createFields.H:26: error: ‘alpha1’ was not declared in this scope
createFields.H:26: error: ‘rho1’ was not declared in this scope
createFields.H:26: error: ‘cp1’ was not declared in this scope
createFields.H:26: error: ‘rho2’ was not declared in this scope
createFields.H:26: error: ‘cp2’ was not declared in this scope
createFields.H:42: error: ‘rhoPhi’ was not declared in this scope
In file included from interTempFoam.C:85:
TEqn.H:5: error: no matching function for call to ‘ddt(Foam::volScalarField&, <unresolved overloaded function type>)’
TEqn.H:6: error: no matching function for call to ‘div(Foam::surfaceScalarField&, <unresolved overloaded function type>)’
TEqn.H:7: error: no matching function for call to ‘laplacian(Foam::surfaceScalarField&, <unresolved overloaded function type>)’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:11: warning: unused variable ‘transonic’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:14: warning: unused variable ‘nOuterCorr’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:3: warning: unused variable ‘nCorr’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:8: warning: unused variable ‘momentumPredictor’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:11: warning: unused variable ‘transonic’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:14: warning: unused variable ‘nOuterCorr’
make: *** [Make/linux64GccDPOpt/interTempFoam.o] Error 1"

Can anyone help me? Thanks in advance

Regards

Metro

[ovie]

Originally Posted by metro

Hey guys,

I am new to openfoam and c++ programming and I am trying to implement Ovie's temperature field code. As I am trying to compile it I keep receiving the following error message

"metro@ubuntu:~/OpenFOAM/OpenFOAM-1.6.x/applications/solvers/multiphase/interTempFoam$ wmake
SOURCE=interTempFoam.C ; g++ -m64 -Dlinux64 -DWM_DP -Wall -Wno-strict-aliasing -Wextra -Wno-unused-parameter -Wold-style-cast -Wnon-virtual-dtor -O3 -DNoRepository -ftemplate-depth-40 -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/transportModels -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/transportModels/incompressibleTemp/lnInclude -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/transportModels/interfaceProperties/lnInclude -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/turbulenceModels/incompressible/turbulenceModel -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude -IlnInclude -I. -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/OpenFOAM/lnInclude -I/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/OSspecific/POSIX/lnInclude -fPIC -c $SOURCE -o Make/linux64GccDPOpt/interTempFoam.o
In file included from interTempFoam.C:58:
createFields.H: In function ‘int main(int, char**)’:
createFields.H:26: error: ‘alpha1’ was not declared in this scope
createFields.H:26: error: ‘rho1’ was not declared in this scope
createFields.H:26: error: ‘cp1’ was not declared in this scope
createFields.H:26: error: ‘rho2’ was not declared in this scope
createFields.H:26: error: ‘cp2’ was not declared in this scope
createFields.H:42: error: ‘rhoPhi’ was not declared in this scope
In file included from interTempFoam.C:85:
TEqn.H:5: error: no matching function for call to ‘ddt(Foam::volScalarField&, <unresolved overloaded function type>)’
TEqn.H:6: error: no matching function for call to ‘div(Foam::surfaceScalarField&, <unresolved overloaded function type>)’
TEqn.H:7: error: no matching function for call to ‘laplacian(Foam::surfaceScalarField&, <unresolved overloaded function type>)’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:11: warning: unused variable ‘transonic’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:14: warning: unused variable ‘nOuterCorr’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:3: warning: unused variable ‘nCorr’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:8: warning: unused variable ‘momentumPredictor’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:11: warning: unused variable ‘transonic’
/home/metro/OpenFOAM/OpenFOAM-1.6.x/src/finiteVolume/lnInclude/readPISOControls.H:14: warning: unused variable ‘nOuterCorr’
make: *** [Make/linux64GccDPOpt/interTempFoam.o] Error 1"

Can anyone help me? Thanks in advance

Regards

Metro

Hi Metro,

You still need to make additional modifications to the codes for the interTempFoam solver like creating variables in the createFields.H. I am attaching a copy of my codes for your review to see if this would help.

ovie

[ovie]

Originally Posted by ovie

Hi Metro,

You still need to make additional modifications to the codes for the interTempFoam solver like creating variables in the createFields.H. I am attaching a copy of my codes for your review to see if this would help.

ovie

I missed something:

In addition, rather than making hard coded changes to the twoPhaseMixture class, I implemented a derived class for the thermal properties which inherits the functions in the incompressibleTwoPhaseMixture class. This I called incompressibleTwoPhaseThermalMixture class. This directory should be in the $FOAM_SRC/transportModels/incompressible directory

Dont forget to recompile this as well.

ovie

[metro]

Ovie,

This is great thank you for your response and sharing your files. Very much appreciated. I am just in the progress of reviewing it and finding what is wrong with mine.

I noticed in the previous responses that you were considering adding mass transport due to evaporation and condensation to your interFoam solver. Did you have any success?

[ovie]

Originally Posted by metro

Ovie,

I noticed in the previous responses that you were considering adding mass transport due to evaporation and condensation to your interFoam solver. Did you have any success?

metro,

I have not been able to implement that yet. I am still trying to benchmark the results for my falling film flow with experimental results. Are you working on something similar?

ovie

[metro]

Ovie,

Thank you again for attaching the additional files. As I said before I just need to review and compare it to my code. I am currently looking at condensation/evaporation in two phase horizontal pipes so in affect one can roughly look at it at a radial falling film problem. How accurate is our model compared to its benchmarks?

Regards

Metro

[ovie]

metro

Originally Posted by metro

Thank you again for attaching the additional files.

No worries. We all need the help and the forum has truly been invaluable.

Originally Posted by metro

I am currently looking at condensation/evaporation in two phase horizontal pipes so in affect one can roughly look at it at a radial falling film problem.

Metro

Not quite sure you can make that sort of direct comparison i.e. confined pipe flow vs free surface falling films. For free surface flows, the deformation of the free surface significantly enhances the transport process (some literature claim up to 30% increase in local heat transfer coefficients) when compared to flat surfaces without deformation. This is something you dont observe with confined flows.

The part relating to phase change however is definitely similar so long as temperature is the driving force for evaporation/condensation.

Originally Posted by metro

Ovie,

How accurate is our model compared to its benchmarks?

Regards

Metro

I have run some computations for a non heated flat plate. I am comparing the wave characteristics (instabillities that form on the free surface) for different inlet perturbations with results from the following papers:
1. D. Gao et al :Numerical simulation of wavy falling film flow using VOF method (Journal of Computational Physics)
2. Z. F Xu et al: Mass transfer across the falling film: Simulations and Experiments (Chemical Engineering Science)

The results are pretty encouraging so far.

[metro]

Hey Ovie,

I have implemented your code and I am now having stability issues. When I run my case, the temperature diverges and ends up reaching inf. I tested the case/boundary conditions with interFoam and no issues. When it came to implementing the div(rho*phi*cp,T) in fvScheme, I used a standard linear Gaussian discretisation. Is this correct? Or did you utilise another methodology? Can you think of any other reason why this is happening?

Regards

Metro

[ovie]

Originally Posted by metro

Hey Ovie,
When I run my case, the temperature diverges and ends up reaching inf

Try a higher grid resolution. I had similar problems earlier. The instability is typically around the interface region and simply grows unboundedly with time. Somewhere in the forum someone said this could be a consequence of the large difference in thermal properties for both fluids. A better resolved grid could help.

Originally Posted by metro

Hey Ovie,
When it came to implementing the div(rho*phi*cp,T) in fvScheme, I used a standard linear Gaussian discretisation. Is this correct?
Metro

I think that is a fair scheme to employ for the discretization.

Goodluck!

[ovie]

Metro,

please note that this solver is specifically for laminar flows. You have to make some modifications to extend to turbulent flows. Just needed to clarify this so u are using it for the right flow problem.

rgds..

[ovie]

Originally Posted by metro

Hey Ovie,

I have implemented your code and I am now having stability issues. When I run my case, the temperature diverges and ends up reaching inf. I tested the case/boundary conditions with interFoam and no issues. When it came to implementing the div(rho*phi*cp,T) in fvScheme, I used a standard linear Gaussian discretisation. Is this correct? Or did you utilise another methodology? Can you think of any other reason why this is happening?

Regards

Metro

Hi Metro,

Like I said in an earlier reply, if the instability you speak of is at the interfacial region then its almost a consequence of the sharp jump in material properties at the interface. In other words, if the difference in specific heat (cp) between the two fluids is significant then grid resolution might not solve the problem. Besides, since cp is such a huge value compared to the other fluid properties, it only serves to blow up the convective term in the discretized governing equation relative to the other terms. So the instability is persistent.

However, you may modify the equations if certain conditions apply: if the specific heat of the fluids is not dependent on temperature or still weakly dependent or if the temperature differences within the flow (eg between the heated or cooled wall and the incoming fluid) is not "large" (note that you would have to evaluate the fluid properties at the mean temperature), you can move cp from the convective term to the diffusive term so that the governing equation to solve is

d(rho*T)/dt + div(rho*V*T) = div((kappa/cp)*grad T) + source terms/cp

This is a more stable problem. For your TEqn.H, you can do:

volScalarField kappa = twoPhaseProperties.kappa();
volScalarField Cp = alpha1*cp1 + (scalar(1) - alpha1)*cp2;

fvScalarMatrix TEqn
(
fvm::ddt(rho, T)
+ fvm::div(rhoPhi, T)
- fvm::laplacian((kappa/Cp), T)
);

TEqn.solve();

rho and rhoPhi already come from alphaEqnSubCycle.H and alphaEqn.H respectively. In this case, you would no longer require rhoCp and rhoPhiCp.

Hope this helps.