[roenby]

Could you post a zipped version of your setup files (0, constant and system folders)?

[kal1943335]

Please find the file attached herewith.
Kalpana

[roenby]

Hi again Kalpana

I now tried to run your case with my interDyMFlow solver, and I get no warnings. It runs fine. However, to be able to use your Allrun, I had to "mv 0 0.orig" and "mv 0.orig/alpha.water.orig 0.orig/alpha.water". Else I get errors already in setFields.

The warning you get basically says that the face cutting procedure of isoAdvector finds more cutpoints (3) along a face than it is supposed to (2). Plotting in the cutpoints as I did in the attached figure, you see that it is at a point where isoAdvector should not be messing around, i.e. deep down into the bulk of the water column. So something is definitely wrong with your setup/solver. I don't know what.

It should be said that I have been running with my OpenFOAM-v1706 installation, but I don't think this is important.

Best,
Johan

[kal1943335]

Johan,
I'll check it over this weekend.

Thank you again for checking it.

Kalpana

[Fauster]

Dear Johan Roenby,

What a fantastic job you did !
Could you tell me if isoAdvector is robust for high CFL number ? For example for ship flow ?
I would like also to use also HRIC and CICSAM schemes. Are they available with the interIsoFoam implemented in ofv1712 (it seems to not to be the case but I can be wrong) ?
Do you think it's hard to code HRIC, mHRIC or CICSAM in interFoam ? I was thinking to use your work available in https://github.com/isoAdvector as a begin.

Thanks for your help
Paul

[roenby]

Hi Paul

Thanks for the nice words.

isoAdvector is based on a model of how the local interface moves across faces into neighbour cells. In this model, the interface in a cell is treated as approximately planar and moving with constant speed during the time step. These model assumptions will not be satisfied for too large time steps and so accuracy will be lost.

Further, in the current implementation the "isoAdvection" is only done to downstream face neighbours of interface cells. Thus, if the time step is so large that the interface moves into these face neighbours and out again "on the other side" during the time step then the advection through the downstream faces of these face neighbours will not be done correctly.

In summary, running isoAdvector in its current state with maxCo and maxAlphaCo > 1 is not a good idea (although if your large velocities are mainly tangential to the interface you might be lucky that it works).

HRIC, mHRIC and CICSAM are not provided with OpenFOAM.They are provided as legacy code at github.com/isoAdvector. I don't use them or keep them updated, so you should expect that modifications are necessary to make them work with the latest OpenFOAM API changes. This should not be too difficult (depending on your OpenFOAM development skills, of course). They were provided by Prof. Hrvoje Jasak.

Kind regards,
Johan

[Fauster]

Originally Posted by roenby

Hi Paul

Thanks for the nice words.

isoAdvector is based on a model of how the local interface moves across faces into neighbour cells. In this model, the interface in a cell is treated as approximately planar and moving with constant speed during the time step. These model assumptions will not be satisfied for too large time steps and so accuracy will be lost.

Further, in the current implementation the "isoAdvection" is only done to downstream face neighbours of interface cells. Thus, if the time step is so large that the interface moves into these face neighbours and out again "on the other side" during the time step then the advection through the downstream faces of these face neighbours will not be done correctly.

In summary, running isoAdvector in its current state with maxCo and maxAlphaCo > 1 is not a good idea (although if your large velocities are mainly tangential to the interface you might be lucky that it works).

HRIC, mHRIC and CICSAM are not provided with OpenFOAM.They are provided as legacy code at github.com/isoAdvector. I don't use them or keep them updated, so you should expect that modifications are necessary to make them work with the latest OpenFOAM API changes. This should not be too difficult (depending on your OpenFOAM development skills, of course). They were provided by Prof. Hrovje Jasak.

Kind regards,
Johan

Thanks for answering me so quickly.
Now I understand the limitations of isoAdvector for large Courant Number. Could it be possible to avoid these limitations ? In most engineering applications having the possibility to not loose too much accuracy or stability with Co >>1 is really important. For me it's one of the problem with OpenFOAM in comparison with commercial code like Fluent or STAR. Layer insertion with SHM is also still a problem ^^!

We will do our best to make them work with latest version of OpenFOAM ! Right now we are learning (doing bibliography) if we need to preserve the compression term or not in the alpha equation (Rusche formulation) with CICSAM or HRIC.

Paul

[saeed.barzegar.v]

Originally Posted by roenby

Hi Kalpana

The isoAdvection construction should not take place in a ".C" file. It should be in the alphaEqn.H file as I wrote. Within the "{}" brackets so it is destroyed properly at each time step.

Also make sure that if your surfCellTol is set to 1e-6 then your p_rghFinal should be lower, e.g. 1e-8 to avoid isoAdvector isocutting cells that are not surface cells but are just cells with bad volume conservation (sum(phi) != 0).

Cheers,
Johan

Hi Johan,

Thanks a lot for all your helps here. I implemented isoAdvector in waveFoam solver in of1712 and ran different cases with your suggested values (https://www.openfoam.com/releases/op...cs-isoadvector).

Here, you mentioned that surfCellTol must be less than tolerance of p_rghFinal, for the sake of mass conservation, and I considered in my caese.

My results (~20 cases with different setup on very fine mesh) showed bubbly water when I used your suggested setup, and I could fix the issue only when I used snapTol less than tolerance of p_rghFinal. Do you have any explanation or idea for this?

At first, I thought this bubbly water might have been because of the gap between surfCellTol (e.g. 1e-6) and snapTol (1e-12), for example alpha=1e-9. But, my results showed this deference is not the source of bubbly water.

Thank you,
Saeed

[arjun]

Originally Posted by Fauster

Thanks for answering me so quickly.
Now I understand the limitations of isoAdvector for large Courant Number. Could it be possible to avoid these limitations ? In most engineering applications having the possibility to not loose too much accuracy or stability with Co >>1 is really important. For me it's one of the problem with OpenFOAM in comparison with commercial code like Fluent or STAR. Layer insertion with SHM is also still a problem ^^!

We will do our best to make them work with latest version of OpenFOAM ! Right now we are learning (doing bibliography) if we need to preserve the compression term or not in the alpha equation (Rusche formulation) with CICSAM or HRIC.

Paul

I do not think that Fluent or StarCCM able to run explicit version of VOF with Co > 1.



I think that their implicit version may run but it is at the cost of too much smearing of the interface. If interface smearing is acceptable to you then you can perhaps tune HRIC and CICSAM to behave like that too.

[roenby]

Originally Posted by saeed.barzegar.v

Hi Johan,

Thanks a lot for all your helps here. I implemented isoAdvector in waveFoam solver in of1712 and ran different cases with your suggested values (https://www.openfoam.com/releases/op...cs-isoadvector).

Here, you mentioned that surfCellTol must be less than tolerance of p_rghFinal, for the sake of mass conservation, and I considered in my caese.

My results (~20 cases with different setup on very fine mesh) showed bubbly water when I used your suggested setup, and I could fix the issue only when I used snapTol less than tolerance of p_rghFinal. Do you have any explanation or idea for this?

At first, I thought this bubbly water might have been because of the gap between surfCellTol (e.g. 1e-6) and snapTol (1e-12), for example alpha=1e-9. But, my results showed this deference is not the source of bubbly water.

Thank you,
Saeed

Hi Saeed

Sorry for the late reply.

Your problem is most likely related to the bug described here.

In short, you have to use the latest version of the isoAdvector code either from openfoam.com or github.com/isoadvector.

Please let us know if this solves your problem.

Best,
Johan

[saeed.barzegar.v]

Originally Posted by roenby

Hi Saeed

Sorry for the late reply.

Your problem is most likely related to the bug described here.

In short, you have to use the latest version of the isoAdvector code either from openfoam.com or github.com/isoadvector.

Please let us know if this solves your problem.

Best,
Johan

Thanks Johan,

I am still getting that bubbly water error.

I downloaded the last version of isoAdvector (which is debugged for those pos and neg functions), renamed everything in all codes from iso(something) to iso222(something) to avoid any duplicative name issues in of-v1712; then compiled (interFlow solver was also checked); and coupled with my solver. But, I am still getting same Warning as I had before. This Warning indicates that I have bubble in my water phase!
I have attached my case set up, as well as that Warning in log file and bubbles in my results.

I really appreciate if you could give me some advise to fix this problem.


Best,
Saeed

[dzordz]

Hi Johan,

hat off for some great work with development of IsoAdvector.
I am wondering are there any plans to implement it in compressibleInterFoam?

Best,
G

[hamsteri15]

I've been testing interIsoFoam from the OF v1806 (development version) in a simple floating box case. When I turn off the box dynamics the solver runs fine but with dynamics (mesh motion) enabled bubbles start to occur in the water region. These bubbles occur only in the region where the mesh is moving (see the attached figures) and so far I haven't been able to get rid of them with different settings. Has anyone else experienced anything like this? I also attached the case which can be run by:
./initCase.sh
./run.sh

[roenby]

Originally Posted by hamsteri15

I've been testing interIsoFoam from the OF v1806 (development version) in a simple floating box case. When I turn off the box dynamics the solver runs fine but with dynamics (mesh motion) enabled bubbles start to occur in the water region. These bubbles occur only in the region where the mesh is moving (see the attached figures) and so far I haven't been able to get rid of them with different settings. Has anyone else experienced anything like this? I also attached the case which can be run by:
./initCase.sh
./run.sh

This is as expected. In the release notes here, we remark that: "Although interIsoFoam can now also be run with other moving mesh classes, the solver has not yet been modified to work consistently with cells of changing shape and volume. The aim is to include this in v1812."


Kind regards,
Johan

[roenby]

Originally Posted by dzordz

Hi Johan,

hat off for some great work with development of IsoAdvector.
I am wondering are there any plans to implement it in compressibleInterFoam?

Best,
G

This is being worked on. But not by me. By all likelihood an isoAdvector based compressibleInterFoam solver will be released in the not too distant future. But unfortunately I cannot say when.



Kind regards,
Johan

[dzordz]

Hi,
I have a thought regarding the spurious currents development.

As I understand the source of spurious currents in interFoam based solvers is the implementation surface tension force term in momentum equation. More specifically the erroneous calculation of the interface normal, the gradient is f.

f=sigma * div(grad(alpha)/abs(grad(alpha))) * grad(alpha)

The issue comes from the fact that alpha field is a step discontinuous function and calculating a gradient of such function will introduce errors and the f term will "overshoot". This overshooting will be balanced by the velocity terms in momentum equation, hence the spurious currents will appear. Duong et. al 2013 have shown that by smoothing the alpha field in the f term, reduces the spurious currents magnitude by an order of magnitude. The advantage comes from solving for gradients on a continuous alpha function and thus reducing the error calculation, reducing the spurious currents.

With the isoAdvector you geometrically reconstruct the isoSurfaces and calculate the unit normal vectors of the actual interface (the normal that we have issues calculating in f term). You then use this to advect the isosurface inside the cell but then this information doesn't get used further. This is valuable info that is not available with MULES algorithm. My thought here is could you somehow use this sub-cell info for the the calculation of the surface tension force (f term) in interfaceProperties.C. Isn't this a superior way to calculate the normals, hence the error should be much smaller and consequently the spurious currents should reduce.

Hope my train of though is clear. Any thought on this?

P.S. It makes sense to me that the spurious currents with your isoAdvector are higher than with MULES as you show in your papers. Since you have a sharper interface, the step function is steeper and the gradient is calculated with bigger errors, hence bigger parasitic currents.

[roenby]

Hi dzordz



Thank you for sharing your idea! Your train of thought is indeed very clear. Lately, I have been collaboration with Henning Scheufler from DRL in Germany on investigating the quality of interfaces reconstructed with the iso-surface based method employed by isoAdvector. It has turned out that it has some of the same problems as the gradient based methods when it comes to mesh convergence, especially for the local interface orientation (to which curvature is of course tightly connected) and for fine unstructured meshes. Henning has worked intensively on developing a new reconstruction method which is accurately converging on general meshes. The initial preprint on this work can be found here. The paper will hopefully soon come out (in a rather extended version) in JCP. The OpenFOAM based codes for the new method will be released together with the paper.



This will then allow the kind of things you talk about: exploiting the improved interface orientation to get improved curvature estimates and hence reduced spurious currents.

[wdx_cfd]

Hi Johan,

Thank you for your execellent work!

I am a new Foamer and still has a longway to go. I have downloaded your code in OF1806, and I want to develop a toolbox with the mass source wave generation method.

Recently, I have finished this kind of work with InterFoam. I found that the alpha value in the source region is always larger than 1.0 because of the added mass. I forced the value in this region equals to 1.0 and it seems better. I wondered that if there will be some solutions in your method?
(Or I would ask:
I found that there is a "subsequent clipping" step in your program to ensure the boundedness of alpha. How will this effect the result of the mass source wave maker?)

Thank you!

Best wishes,
D.X. Wang

[mrishi]

Hi Johan,
It is very delightful to read through your work over the past week or so. I am looking towards such a geometric interface in multiphaseInterFoam setup. Do you have some pointers on how one would proceed, and/or if you know of work on this front?
I suppose the major differences lie in how two secondary phases are addressed by mpIF as compared to iF where this scenario doesn't exist. Based on this convention of phase IDs, we should be able to setup the isoFace normals and calculate dVf_, but problems like triple junction can come into picture.

[roenby]

Originally Posted by mrishi

Hi Johan,
It is very delightful to read through your work over the past week or so. I am looking towards such a geometric interface in multiphaseInterFoam setup. Do you have some pointers on how one would proceed, and/or if you know of work on this front?
I suppose the major differences lie in how two secondary phases are addressed by mpIF as compared to iF where this scenario doesn't exist. Based on this convention of phase IDs, we should be able to setup the isoFace normals and calculate dVf_, but problems like triple junction can come into picture.

I am not working on this subject and don't know if anyone else is. I guess you could include an alpha field for each pair of phases e.g., alpha12, alpha13, alpha23 for three phases. Then if the phases are well-resolved by the mesh there will only be two phases present in most interface cells, so only, say, alpha12 will be between 0 and 1 in a given cell and you can use isoAdvector "out of the box" to propagate the interface in that cell. For cells containing all 3 phases I am not so sure what to do. I want to mention, that actually in the present implementation of isoAdvector the situation with multiple interfaces in a single cell is treated geometrically consistently (e.g. the situation just before a droplets hits a flat water surface so the flat water surface is in the bottom of a cell and the droplet enters the top of the same cell). Treating such situations more consistently is on my todo list, but for now, the solver seems to cope with this situation surprisingly well. You might be lucky that something similar is the case for the triple points in a 3-phase flow, although a complicated splitting of the cell into 3 subcells would be more proper.