Hallo Zhipeng

Before Fabian answers you, i can give you some info. regarding this

DC is a darcy constant which is used to switch off or stiff the velocity in the solid domain of a melting problem during phase change. It is derived from Caman-Kozeny equation used to solve porous flows. The interface between the solid and liquid is called as mushy zone which is treated as porous zone. Moreover
Which makes DC zero in the liquid zone and very high value in solid zone.

This was implemented by many researchers,but however voller and prakash paper is major cited of this kind.


http://ac.els-cdn.com/00179310879031...8aa918c64efc5d

Thanks
Rohith

Hi Rohith , Thank you very much for your help . I will read the paper you recommended for more information . Moreover , I have some difficulties with the following code . Can you give me some explanation ?

Hi Zhipeng


This is however, followed in major buoyancy bousinessq approximations in openfoam. According to my understanding, it is basically balancing the p_rgh ( buoyant pressure) and p(pressure) by settling reference values defined in Peqn. It is a small form of declaration however on the first place.

Go through the OpenFOAM programmers guide for detail on all such definitions.

Regards


Rohith

Hi

this code snipped checks weather a reference cell is needed or not and sets the pRefValue to the reference cell pRefCell. In incompressible flow, pressure is relative (the solution does not depend on the actual value of the pressure but its gradient). Thus, in closed domains you have to specify pressure somehow. This is done by the reference cell.



Regards

Fabian

Hi ,
Thanks ! Both of you has give me great help . But can you explain me some thing about p_rgh ( buoyant pressure) and p(pressure) . Or can you send me a paper about this ?

moreover , Fabian
about this function
The refering value is for p , then what is the "p_rgh" for ?

What's more , rhok is different in z direct, can we use p=p_rgh+gh*rhok
Regards

Zhipeng

Hi,

p_rgh is not the buoyant pressure. It is pressure p minus the static pressure (which is rho*g*h). So p_rgh is the dynamic pressure. As pressure is relative, you set p and p_rgh at the reference cell and calculate the p and p_rgh field relative to this cell.

Cheers

Fabian

Hi Zhipeng

i am sorry for the false nomenclature. It is actually dynamic pressure and static pressure, as fabian said. He is absolutely right.

Actually the Total pressure(p) is the sum of p_rgh and rhogh. The term rhogh actually comes from the formulation of bernoullis theorem, which makes the total pressure the sum of static and dynamic pressures. generally this rho*gh is computed through
rhi*gravity*position of the cell center, this is due to the position of the computational cell. which is in physical sense is due to the virtue of its position.

Thanks
Rohith

big temperature difference convection
 

Hello Fabian ,
your solver about melting uses the Boussinesq Approximation . But as you know ,
boussinesq approximation will bring in large error when the temperature didderent is big .
I am calculating a melting problem about copper contact , and the temperature didderent is very big , can you tell me what solver is suited for my problem ?
Thanks ,
Zhipeng

temperature dependant density for all terms in all conservation equations
 

Hi Zhipeng

You will not find a ready-made solver for that. However, here is what I did for such problems: I switched from the Boussinesq approximation to a temperature dependent density for all terms in all conservation equations. As the density is calculated from temperature only, you can easily change the convMeltFoam solver for your needs (Keeps in mind that you still have an incompressible solver and you do not need the compressible mass conservation equation).
I do not have access to my temperature dependent density melt solver now. So have a try - the changes are of minor kind.

Regards

Fabian

Hello,

i took the solver of Fabian (convMeltfoam, Post 80/81) and implemented it in interFoam. Instead of using rhok and the boussinesq approximation i used a temperature depended rho i read in every iteration (compare Anja's post) everything is working (compiling and calculating) but when it comes to the results of Fabian (see attached pictures, the left one is with the untouched solver of Fabian, the right one is with extended interFoam) I'm wondering about two things:

First the effect of the convection in Fabian's calculation is much bigger (see T and alpha)

Second the values of p_rgh are much different. I didn't touch the pEqn of interFoam. But as p_rgh is the dynamic pressure my results look wrong. Do you have an idea why this is so?

thank you in advance
Thilo

Isothermal melting
 

Hi Foamer's,
Can any one advice me whether the recent method (recently publish by Fabian:http://www.cfd-online.com/Forums/ope...problem-5.html) is appropriate for the isothermal melting ? I'm trying to use the code for isothermal melting taking narrow temperature width (i.e. Tl-Ts <0.5). Please suggest me, whether the method is good or not OR the appropriate method.

Thanks in advance.
#shakil

Modifying the solver to include MultiRegion
 

Hi,

Thanks for this post. I found it very helpful for the foamers working in phase change materials. I need to modify Fabian's solver to solve for a multiregion problem, I know that my first steps should be looking into chtMultiRegionFoam.
Any ideas on how can this be implemented ??

Thanks,
David

running solver in parallel
 

Hi all,

first of all thank you for providing the solvers.

I tried to run the convMeltFoam solver for the OF Versions 2.1 and 2.3. Both are running fine in serial. But I recognized that the simulation stops after some timesteps without crashing when I try to run the tutorial cases in parallel.

I also tried it with another case in 3D with the same result: Running fine in serial and stopping in parallel.

Does anybordy has an idea why the solver has problems in parallel runs?
Did you discover similar problems?

Thanks in advance!

Kind regards

Chrisi

Running in parallel
 

Hi,
Have you checked the space available in your drive? I don't think it will cause any problem for running in parallel. One thing you can check, whether you used gSum instead of sum in the code, because sum will not work for parallel.

Cheers
shakil

Hi Shakil,

thank you for your hint. Indeed there was used "sum" instead of "gSum" in the solver.

But now I removed that function and it still stoping/crashing.

On another machine the parallel run now ends up after 10 or 15 timesteps with the following error:
What I can also observe that the crash occurs when the solver should solve the p_rgh equation.

By the way: Other solvers run fine in parallel. So it should not be a general mpi problem.

What can be the reason therefore?

Kind regards

Chrisi

Hi Chrisi1984,
It seems a bit strange !!! Is it possible to upload your exact solver that you are compiled ?
Cheers
shakil

Hi Shakil,

attached you find the solver that I am using for version 2.1.

It is nealy that solver that was posted in this thread few pages before:

I hope you can help me finding the problem for the parallel simulation.

Kind regards,

Chrisi

Hello Chrisi,

the same which applies for sum and gSum has to be used for max and gMax. If you want to have the sum or the maximum of a field, which could be spread on several nodes when using parallel, the function name is gMax, not max.

Thus, change line 31 in TEqn.H from
to

Then it works, at least on my computer. So I hope, this solves it for you as well. By the way, this should be the same with Version 2.3.

Regards, Anja

Parallel convMeltFoam
 

Hi

attached you find the parallelized convMeltFoam solver. It runs parallel now after decomposing the case. I tested it on a 3D case with 10M cells.

Regards

Fabian

convMeltFoam with conjugate heat transfer
 

Hi fabian_roesler,
Thanks for the parallelized convMeltFoam solver. Have you worked with the melting problem with conjugate heat transfer ? Any help regarding to this will be great for me.

Thanks in advance
shakil

third solid phase for fins in solid liquid phase change
 

Hi shakil

Well, I did not perform real conjugate heat transfer simulations with two different meshes for PCM and solid/gas like in chtMultiRegionFoam. However, I did some simulations to study the influence of fins on the overall melting process.
What id did basically was to add another, third phase to my solver with new thermophysical properties. The energy conservation equation doesn't change and for the momentum conservation equation I introduce a switch off technique that keeps a zero velocity in the fin. I posted the solver in this thread (post #81):

http://www.cfd-online.com/Forums/ope...tml#post467203

Hope this tip pokes you into the right direction. Alternatively you could implement the solid/liquid phase change into chtMultiRegionFoam.

Cheers

Fabian

multi-species solidification
 

Hello everyone,

I wish to model the solidification of water -> ice in the presence of a gas such as air, just in a box as a simple test case. My question is would this be possible using the convMeltFoam solver?

As such I want to create a hybrid case somewhere between damBreak (interFoam) and convMeltFoamOF230.

Having looked at the convMeltFoamOF230 example I'm a little confused as to the significance of the alpha3 scalar field in the initial 0 folder, particularly given that the solver proceeds without writing it out to file. I understand the concept of using a Darcy constant but I'm still not sure if it needs to be there.

Thanks.

compressibleInterFoam and convMeltFoam
 

Hi steph

I did exactly that within by PhD thesis. I defended the thesis one month ago and it will be published soon. Unfortunately it is written in German. I combined the compressibleInterFoam solver with my new melting solver based on the convMeltFoam. The density of the PCM is no longer described by means of a boussinesq approximation but by a temperature dependent density in all terms of the conservation equations. This makes it possible to simulate volume change during melting and solidification in a closed volume/capsule.
As soon as my thesis is published, I will post at least parts of the solver here.

Cheers

Fabian

That's great. When you've finished translating your thesis into English let me know... Natürlich mache ich Spaß. Es gibt die Möglichkeit, dass ich Ihre Doktorarbeit lesen kann. Falls ich Ihre Doktorarbeit/Code nicht verstehe, wegen meines schlechten Deutsches, dann werde ich mit Ihnen oder meinem Professor, der auch ein Deutscher ist, sprechen.

In all seriousness that sounds like a sensible approach. I look forward to seeing the code. :)

Hello Fabian,

thanks for the solver that also runs in parallel.

Just as information for others: The solver that Fabian recently posted for parallel usage was made for OF 2.3.

Kind regards
Chrisi

Hullo All

I am new to Foam but I have been trying to simulate 3D melting phenomena using AVL Fire for my master thesis. AVL works well for hexa mesh but when i go for a tetra mesh it diverges. Now, this has to got me to investigate the possibility of running my simulation in Foam. I set up the case with Foam (thanks to the solver developed by Fabian). Its running well and I can already see some melting.

But the problem I have is : I had to resolve the mesh near heater elements to make sure correct transfer of heat flux occurs. But this has made my calculation very expensive. The reason is I have varying cell size and my max courant number is close to one. But when it calculates at the fine cells, the mean courant number is quite small, of the order 1e-5. I am thinking this is slowing down my calculation. Could anyone kindly tell me -

1) Is it necessary to have a finer/finest resolution near heater elements?

2) If it is necessary, then is there any way to solve the issue of slower convergence by somehow trying to work around the courant number issue?

Cheers!

Increase outer correctors and maxCO
 

Hi

as the solver is based on the PIMPLE algorithm you can increase your max Co above one. This is achieved by outer Simple iterations and thus, the conservation equations can be under relaxed. This could give you a boost, as the number of inner iterations for energy conservation should more or less stay the same for one time step. This is just a guess but try the following:

maxCo 5
nCorrectos 1 or 2
nOuterCorrectors 10
lower the number of maximum iterations for alpha to say 10 to 15
add some under relaxation factors and convergence control

If you play with these options you should get some speedup. Moreover you could use GAMG solver for pressure and smoothSolver for the other equations. And moreover, try to use more CPUs :)

Cheers Fabian

Thanks for the inputs Fabian. I played around a bit and looks like Courant no of 10 gives the fastest time possible. Maybe reducing the nCorrectors to 1 (I have two) might boost my speed up further more. But the overall simulation time is in acceptable range for now. I could do with more CPUs but i am very much limited to 30 cores.

Will post here again if i find something new or have further issues!

Cheers

hullo

I am a little stuck up with boundary conditions for wall. Let me explain. I am solving a 2D melting problem for comparison with AVL Fire. I have a surface named THOT and it is at a temperature of 330K. Rest of the domain has adiabatic wall condition with the IC of the medium being a solid at a temperature below the freezing point. I have run simulations by taking this condition in the T file in 0 folder.

THOT
{
type fixedValue;
value uniform 330;
}

Now I want to change this condition in such a way that the value of 330K is alloted to THOT only at the 0 time step. It should not remian fixed throughout the simulation. Rather the temperature at THOT has to be calculated at every time. How can I do that? Would be glad to hear some suggestions :)

Cheers

Well, do you have a fixed heat flux through the wall or is it an adiabatic wall? Second would be zeroGradient and first depends on what you intend to simulate. For example heat transfer through a wall and an outer/infinite temperature is implemented by

In general: Have a look into the OpenFOAM C++ documentation
http://foam.sourceforge.net/docs/cpp/a10483.html

Cheers

Fabian

Hullo Fabian

I am sorry but there was no problem in my simulation actually. I was making a visualization error in AVL Fire. I took a sample test case(2D) to do some comparison of Foam and Fire. A 2D block with one face being heated at 330K and all other faces at adiabatic conditions. The fluid properties are of AdBlue. All this would form a basis of my master thesis later on and that is the reason I did not take standard test cases available in literature. The dimensions of box is 0.1X0.1X0.003. The region is initially assumed to be solid.

The problem is the melting profiles are not the same for the same BC and IC. I have attached the plots for simulation results at the end of 600 seconds. Can you advise on where I might be going wrong.

Cheers

Attachment 32555

Attachment 32556

And yes, the profile with higher melting is from openfoam....

More info
 

Could you post your setup for OpenFOAM? blockMeshDict, fvSolution, fvSchemes, transportProperties and controlDict. With a litle bit more info the answer could be much more specific.

Cheers Fabian

Why am I getting totally different result using the convMeltFoam (parallel) on OF211?
 

As shown below and the case setup have been attached as well.

Attachment 32764

Attachment 32765

Attachment 32766

Attachment 32767

Vortices during melting of metal at a vertical wall
 

Well, from your transportProperties I can see that you simulate a metal PCM (if I may guess I'd say it is Gallium). Moreover I guess you have a fine mesh.
Have a look into the literature. There you will find lots of articles that describe the melting of metals at a vertical wall. Most authors discovered the formation of several vortices/eddies along the wall. The numerical simulation of the metal melting at vertical walls is very much mesh dependent. So that’s why your results differ from the experimental results from Gau and Viskanta. You dug up a very up-to-date problem that is heavily discussed in literature.

Cheers

Fabian

Thanks for reply. I am following, however, the parameters tabulated in table 1 of your <Heat Mass Transfer> paper on 2011. The only thing I have changed is that I tripled the mesh resolution used in your paper. Why different meshes are giving totally different solutions? Backed up by the experiment conducted by Gau and Viskanta, am I right to say that the converged result on the finer mesh was wrong? But how could it be?

Hi

I do not have access to my literature at the moment. However, what I remember is that the experiments by Gau and Viskanta where stopped at certain times after start of melting; the liquid phase was dumped and the phase front at the remaining solid phase was measured. Then the measurement was restarted and stopped at the next breakpoint. This lead to uncertain results. Campbell and Koster (Visualization of liquid-solid interface morphologies in gallium subject to natural convection) repeated the measurements with a X-Ray analysis and obtained different results. In the younger literature, you will find more and more articles on the phenomena of multiple vortices when melting metal at a vertical wall. Moreover, the results in my paper where conducted with a different solver than the one posted in the forum. I would say your results look reasonable.
Try to make a mesh refinement study. Start with the resolution of Brent and Voller and go up to half a million cells.

Cheers Fabian

Hullo Fabian

I am not at liberty to disclose the properties of PCM. But I feel that Boussinesq approximation is introducing some error here. The approximation is I believe valid for a difference of 2 degrees for water and the PCM i am using is also comparable to water to an extent with respect to density and other properties. But in my setup, I have a temperature difference of 40K between the frozen solid and heat source. Will this introduce some sort of numerical error because fundamentally boussinesq approximation works only for a small range of temperatures?


Cheers

Yes sure, the boussinesq approximation is only valid for small temperature differences. However, the error one faces is not that big and the flow field will not change entirely. I also programed a solver with polynomial temperature dependent density. The simulation domain has a small outlet to account for volume expansion. For the simulations and experiments I compared in my thesis, the only difference was a slightly faster melting of the boussinesq approximation case. However, the PCM I studied was a paraffin with totally different transport and thermal properties.

Cheers

Fabian

Yes, the flow field did not change entirely. When i compare the result with AVL Fire, the difference in profile can be observed. This probably comes with the way the solver is implemented in Fire. I am also noticing a slightly higher melting in openFoam.

Do you have a specific paper/reference for the mathematical formulation of the solver you implemented. Reading that would come in handy for me.

Cheers