CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM (http://www.cfd-online.com/Forums/openfoam/)
-   -   Two-phase Eulerian model for nucleate subcooled boiling (http://www.cfd-online.com/Forums/openfoam/81010-two-phase-eulerian-model-nucleate-subcooled-boiling.html)

 Edy October 13, 2010 10:56

Two-phase Eulerian model for nucleate subcooled boiling

Hi Foamers,

I want to model nucleate boiling using a Eulerian-Eulerian approach (no interface tracking or level set method but derivation of mass, momentum and energy conservation equations for both liquid and vapor phase).

Basically my geometry is a 2D semi-channel in the plane (x,y)
* x=0-axis corresponds to the centerline of the real channel (bulk coolant)
* x=xmax corresponds to the heated wall
* water is flowing along the y-axis

Physically speaking, part of the water will boil in the vicinity of the heated wall, and formed bubbles will subsequently move towards the center of the real channel, i.e towards x=0 and condensate due to the contact with the surrounding subcooled water. In this model, one assumes that the whole heat given by the wall is given to the liquid phase. The vapor phase can only gain energy through evaporation of water.

I started from the twoPhaseEulerFoam solver, and I implemented interfacial forces for drag, lift, virtual mass, turbulent dispersion, wall lubrication. In a first time, I just modeled the condensation of bubbles, assuming that at the inlet, there is some void and these bubbles are condensing while flowing along the channel. BUT FOR THE MOMENT NO HEAT WAS SUPPLY TO THE WALL, so no evaporation was considered. The results were quite satisfying.

Now i want to go one step further and add in my model the evaporation of water close to the wall. Before developing a more complicated boiling model, I tried something very simple. I just took a constant evaporation rate that I applied to the cells in contact with the wall. Then in my momentum and energy equations, I took into account both the variations due to condensation and evaporation, as well as the heat supplied by the wall to the liquid phase in contact with the wall.

However, my simulations crashed very quickly and I get this error message:

Code:

```#0  Foam::error::printStack(Foam::Ostream&) in "/home/michta/OpenFOAM/OpenFOAM-1.7.x/lib/linux64GccDPOpt/libOpenFOAM.so" #1  Foam::sigFpe::sigFpeHandler(int) in "/home/michta/OpenFOAM/OpenFOAM-1.7.x/lib/linux64GccDPOpt/libOpenFOAM.so" #2  in "/lib/libc.so.6" #3  Foam::LimitedScheme<double, Foam::limitedLinearLimiter<Foam::NVDTVD>, Foam::limitFuncs::magSqr>::limiter(Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> const&) const in "/home/michta/OpenFOAM/OpenFOAM-1.7.x/lib/linux64GccDPOpt/libfiniteVolume.so" #4  Foam::limitedSurfaceInterpolationScheme<double>::weights(Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> const&) const in "/home/michta/OpenFOAM/OpenFOAM-1.7.x/lib/linux64GccDPOpt/libfiniteVolume.so" #5  Foam::fv::gaussConvectionScheme<double>::fvmDiv(Foam::GeometricField<double, Foam::fvsPatchField, Foam::surfaceMesh> const&, Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>&) const in "/home/michta/OpenFOAM/OpenFOAM-1.7.x/lib/linux64GccDPOpt/libfiniteVolume.so" #6   in "/home/michta/OpenFOAM/michta-1.7.x/applications/bin/linux64GccDPOpt/EulerEulerBoilingFoam" #7   in "/home/michta/OpenFOAM/michta-1.7.x/applications/bin/linux64GccDPOpt/EulerEulerBoilingFoam" #8   in "/home/michta/OpenFOAM/michta-1.7.x/applications/bin/linux64GccDPOpt/EulerEulerBoilingFoam" #9  __libc_start_main in "/lib/libc.so.6" #10   in "/home/michta/OpenFOAM/michta-1.7.x/applications/bin/linux64GccDPOpt/EulerEulerBoilingFoam" Floating point exception```
It seems that my velocity fields for both phases is diverging. My UEqn for the gas phase looks like that (similar for liquid phase):

Code:

```        UaEqn =         (             (scalar(1) + CVM*rhob*alpha*beta/rhoa)*             (                 fvm::ddt(Ua)               + fvm::div(phia, Ua, "div(phia,Ua)")               - fvm::Sp(fvc::div(phia), Ua)             )           - fvm::laplacian(nuEffa, Ua)           + fvc::div(Rca)           + fvm::div(phiRa, Ua, "div(phia,Ua)")           - fvm::Sp(fvc::div(phiRa), Ua)           + (fvc::grad(alpha)/(fvc::average(alpha) + scalar(0.001)) & Rca)         ==         //  g                          // Buoyancy term transfered to p-equation           - fvm::Sp(beta/rhoa*AD, Ua)         //+ beta/rhoa*AD*Ub            // Explicit drag transfered to p-equation           - beta/rhoa*(liftForce - virtualMassForcea + turbulentDispersionForce + wallLubricationForce)           + (GammaGL-GammaLG)/(alpha+scalar(0.001))/rhoa*Ub           - fvm::Sp((GammaGL-GammaLG)/(alpha+scalar(0.001))/rhoa, Ua)         );```
Do you have any idea what could be wrong?
If you need more details, don't hesitate to ask.

/Ed

 chiven October 14, 2010 20:51

Hi Edy, I am also doing the simulation on two-phase flow with phase change using Euler-Euler model. But you are way ahead of me. This time I am struggling at the stage without considering the mass and heat transfer.

Back to your question, I think you need to add a source term to the pressure equation considering the condensation or evaporation. Hope this can help you.

Best regards,
Chiven

 Edy October 15, 2010 07:19

Hi Chiven,

Your answer is actually very helpful, i'll dig deeper in this p-Eqn and try to see which term to add. With condensation only, I kept the pEqn unchanged and it worked. But my evaporation rate is much bigger and it might cause some problems.

I'll check this issue next week cause I found out mistakes in my interfacial forces implementation that I need to solve first.

If you have questions while developing your model, dont hesitate to ask, perhaps I will have encountered the same problem before.

Thanks again.
Best

/Ed

 yasin.ramzani September 26, 2013 11:31

Quote:
 Originally Posted by Edy (Post 279342) Hi Chiven, Thanks for your reply. Your answer is actually very helpful, i'll dig deeper in this p-Eqn and try to see which term to add. With condensation only, I kept the pEqn unchanged and it worked. But my evaporation rate is much bigger and it might cause some problems. I'll check this issue next week cause I found out mistakes in my interfacial forces implementation that I need to solve first. If you have questions while developing your model, dont hesitate to ask, perhaps I will have encountered the same problem before. Thanks again. Best /Ed
Hi
maybe its better to see master thesis of Edouard Michta (Modeling of Subcooled Nucleate Boiling with OpenFOAM)

 All times are GMT -4. The time now is 05:40.