[andrea93]

Dear all,

I am trying to simulate the Nusselt condensation problem. I used both interPhaseChangeFoam (incompressible) and compressibleInterFoam. In both cases I implemented by myself what I needed to reproduce the correct physics. My fluid is water.

I realized that at 219 bar and 646K everything works fine, while at 1bar and 373.15K the pressure equation gets lost (max number of iterations is reached) and, consequently, the Courant number explodes (>10000). This means that the solver can't solve in the right way the sharp interface (at 1bar water and liquid properties are very different each other, while at 219bar and 646K they are very similar, so the discontinuity is smoother than before - this helps the solver a lot).

To have a further confirmation, I tried to simulate a single phase flow of water (pure fluid-dynamic simulation) in my domain (a rectangle 10*0.05mm^2) : it crashes. Therefore, I realized that MAYBE the problem is not related to the thermal model which I implemented, but to the Vof method, that is used both in interPhaseChangeFoam and bin compressibleInterFoam.

My questions therefore are:
1) how is the U-p coupling threated in OpenFoam's VoF method?
2) how is the interface threated ?
2) and, of course.. do you know how to solve my problem?

Thank you for your availability.
Andrea Iovenitti

[Taataa]

1- OF uses a segregation approach for dealing with the coupling. In the case of transient solver it is PIMPLE = SIMPLE + PISO
2- a VOF approach that solves a transport equation for volume fraction with an additional compression term and using MULES.
3- I am afraid no one would be able to help effectively until you provide more information about the way you setup your case like schemes, solutions, boundary conditions and so on.

[AliE]

Hi,

If you want to learn more about the vof in openfoam look at this paper:
http://iopscience.iop.org/article/10.../1/014016/meta

There you will find a good explaination of the method.

[andrea93]

Thank you very much for the answer. I am sorry to reply so late.

First, I want to specify my problem. As already said, when I try to simulate the Nusselt condensation problem with water at 219 bar and 646 K everything works fine. Conversely, when I use water at 1 bar and 373.15 K the simulation crashes. I am almost sure the reason is related to the way the interface is solved.

I tried to modify my solver in order to solve the problem. As you can see in the files "TEqn.H", "pEqn.H" and "alphaEqn.H", I have implemented the smearing proposed by Kunkelmann (see pag. 35 and 36 of the following link: http://tuprints.ulb.tu-darmstadt.de/...Kunkelmann.pdf). I commented those lines because they were not useful (I guess I implemented Kulkelmann smearing in the wrong way). This also explains why the files started with "cond" and "evap", which are in each time step folder, are not really important).

Since this approach did not help, I tried to follow another path. I have tried to initialize the test case at 1 bar and 373.15 K as follows.
First, I set up a simulation where all the parameters in "transportProperties" were those at 1 bar and 373.15 K. The only exceptions have been the superficial tension "sigma" and kinematic viscosities of water liquid (phase1) and vapour (phase2). I ran the simulation until the final time (1 s) was reached.
After that, I initialized a new simulation with the results I got at time 1 s in the previous one. I changed progressively the values of the kinematic viscosities until the real values at 1 bar and 373.15 K, leaving the superficial tension unchanged. In brief, I have tried to provide a better initialization for an increasingly complex problem (assuming that the complessity for the solverr is associated to the values of sigma and of kinematic viscosities at 1 bar and 373.15 K).
Every time a new simulation was started (all ending at 1 s) the iniziatilazion was given in such a way T, p_rgh, U and alpha.phase1 were those at 1 s in the previous simulation. The simulations did not crash.
At this point, I repeated the procedure increasing sigma from 0.07e-3 to 0.0582 Pa*s, which is the value at 1 bar and 373.15 K.

I attach the solver I use (it is based on interPhaseChangeFoam) and the test case, which is "the last one"; the files in "0" are those computed at 1 s with sigma = 0.01e-1 Pa*s.
The simulation runs if "maxDeltaT 0.5E-06;", but the results are denefintely unphysical.

Do you have any suggestions to solve the issue?

Best regards,
Andrea