[dplthuy]

Hi,

I want to use adaptive mesh refinement for the simulation of droplets in OpenFOAM v2106. I am using the interIsoFoam solver, combined with dynamicRefineFvMesh, which should refine based on the volume fraction field of the VOF simulations.
When performing a simple case of advection of droplets, I found that the adaptive mesh refinement causes breakup of the droplets near the interface. This is problematic since I am ultimately interested in the droplet size distribution of a spray.

So far, I have found three parameters that seem to influence the shattering of the droplets: the refineInterval, maxRefinement in the dynamicMeshDict, and the reconstructionScheme for alpha in fvSolution. In short, they have the following effect:

- refineInterval: Droplets are more likely to breakup when the number of timesteps between remeshing is larger. I have varied this parameter between 1 and 5.

- maxRefinement: Droplets are more likely to breakup when the maximum refinement is larger. I have varied this parameter between 1 and 3.

- reconstructionScheme: Droplets are more likely to breakup when the isoAlpha scheme is used, as compared to plicRDF. Higher values of the previous two parameters are allowed when using plicRDF. A comparison of the performance of these schemes can be seen in the figures, with maxRefinement = 2 and refineInterval = 3.

The comparison of the reconstructionSchemes also shows worrysome behavior of the surface after the first mesh refinement. For both schemes, the droplets do not have a smooth surface anymore after refinement. This is more pronounced for isoAlpha than for plicRDF. It is my suspicion that this causes the issue of the breakup of the droplets. In the final image, the differences can clearly be seen.

I will be glad to take any suggestions as to how I can avoid this behavior. In my opinion, the surface of the droplets should stay smooth when applying the adaptive mesh refinement, but I don't know how to achieve this. I have attached my complete case for those interested.

Thank you for your insights!

[Jun_93]

Hi Dennis,


Have you found the reason/solution for this problem. I recently also tried running interIsoFoam with ARM and encountered the same issue you mentioned. I applied interIsoFoam with the same AMR configuration to test a rising bubble case and a droplets' collision case. The combination runs well for the rising bubble case, but badly for the droplets' collision case. In the collision case, I observed the same unphysical scattered small droplets as you displayed in your results. Hm, I don't know how to insert images into the reply. If you want to know more about my cases, PM me plz.

BR,
Jun

[dplthuy]

Hi Jun,

Thanks for your reply. Have you tried running your case with both the iso-Alpha and the plicRDF reconstruction scheme? Do you observe the same difference that I do?

After reading the publication by Scheufler & Roenby on the different surface reconstruction schemes (https://www.sciencedirect.com/scienc...269?via%3Dihub), I have a hypothesis as to why iso-Alpha shows the shattering and plicRDF does not.

In section 3.3 you can read that plicRDF uses the surface normal obtained in a previous time step as an initial guess in the surface reconstruction of the current time step. Iso-Alpha does not do this, and calculates the surface normals independently every time step.

When applying the first refinement to the spheres, it seems that the smoothness of the surface is not exactly preserved (see my previous post). However, since plicRDF uses the surface normals of the previous time step (where the surface was still smooth), it is able to suppress this. The iso-Alpha scheme is not able to do so, since it does not use this information.

This hypothesis leads to a question, which I am currently trying to answer:
Does the shattering also occur when the VOF interface is not initiated by a setFieldsDict entry?

In other words, if we have an inlet of alpha=1 into a domain with initially alpha=0 everywhere, will we also observe the shattering? I would expect not, since in this case the cells at the interface are already refined before the interface enters the domain. I am investigating this with a case of cylindrical jet breakup (see attached image). So far, I have not observed the shattering for iso-Alpha, nor for plicRDF. However, results are still preliminary.

Please let me know if my train of thought makes any sense to you!

[Zand]

Dear Dennis,

Have you found a solution for this problem?

I am trying to execute AMR on a free rising bubble and I am facing the same problem. I am using interIsoFoam so it is set on plicRDF and it does not solve my problem.

When the mesh gets refined/unrefined it generates unrealistic alpha values inside (like liquid pockets inside the bubble) or outside the bubble (smaller satellite bubbles).

I have solved the issue for the inside by not letting the mesh to unrefine for values higher than 0.5 and there are no more liquid pockets. But I cannot solve the issue of smaller bubbles generated downstream.

P.S An interesting thing to mention. Increasing the background liquid viscosity (from water up to 5 or 10 times) seems to resolve the issue!! But that’s not gonna be the way for me.

[dplthuy]

Hi,

I think I did manage a solution. It is my suspicion that these small 'inclusions' of the opposite phase are related to the recombination of cells in the AMR. This recombination can cause the value of alpha to deviate slightly from 0 or 1, and this seems to happen cumulatively, until at some point the surfCellTol is exceeded.
When this happens, isoAdvector will start to see these cells as interface cells, thus creating these inclusions.
I believe there are two ways to avoid this behavior:

1. Refine the entire bubble/droplet, instead of only the interface. This avoids the recombination of cells within, and thus the creation of these inclusions. But of course it increases the computing cost massively.
2. In fvSolution, set equal values for surfCellTol and snapTol. This means that if a cell is not considered an interface cell, its alpha value will always be either 0 or 1.

Hope this helps you.

Kind regards,
Dennis

[Zand]

Hi Dennis,


Thank you very much for your reply.

Yes I agree with the first solution that you mentioned. That solves the problem for inside the bubble. Now I tried the second solution, and here is how it went:


Setting the snapTol to surfCellTol's default value (1e-10) does not seem to make much changes.
Setting the surfCellTol to snapTol's default value (1e-12) aggravates the problem.


Would you mind taking a look at my files to see if I am doing something wrong here?


DynamicMeshDict:

dynamicFvMesh dynamicRefineFvMesh;

// How often to refine
refineInterval 5;

// Field to be refinement on
field alpha.air;

// Refine field inbetween lower..upper
lowerRefineLevel 0.01;
upperRefineLevel 0.99;

// If value < unrefineLevel unrefine
unrefineLevel 0.5;

// Have slower than 2:1 refinement
nBufferLayers 2;

// Refine cells only up to maxRefinement levels
maxRefinement 2;

// Stop refinement if maxCells reached
maxCells 20000000;

// Flux field and corresponding velocity field. Fluxes on changed
// faces get recalculated by interpolating the velocity. Use 'none'
// on surfaceScalarFields that do not need to be reinterpolated.
correctFluxes
(
(phi none)
(nHatf none)
(rhoPhi none)
(alphaPhi_ none)
(ghf none)
(phi0 none)
(dVf_ none)
);

// Write the refinement level as a volScalarField
dumpLevel true;

fvSolution

solvers
{
"alpha.air.*"
{
interfaceMethod "isoAdvector";
isoFaceTol 1e-10;
surfCellTol 1e-12;
snapTol 1e-12;
nAlphaBounds 3;
clip true;
reconstructionScheme plicRDF;




nAlphaCorr 1;
nAlphaSubCycles 1;
cAlpha 1;
MULESCorr no;
nLimiterIter 3;
solver smoothSolver;
smoother symGaussSeidel;
tolerance 1e-8;
relTol 0;

}

"pcorr.*"
{
solver PCG;
preconditioner
{
preconditioner GAMG;
tolerance 1e-05;
relTol 0;
smoother DICGaussSeidel;
}
tolerance 1e-05;
relTol 0;
maxIter 100;
}

p_rgh
{
solver GAMG;
tolerance 1e-08;
relTol 0.001;
smoother DIC;
}

p_rghFinal
{
solver PCG;
preconditioner
{
preconditioner GAMG;
tolerance 1e-07;
relTol 0;
nVcycles 2;
smoother DICGaussSeidel;
nPreSweeps 2;
}
tolerance 1e-07;
relTol 0;
maxIter 20;
}

"(U|k|C)"
{
solver smoothSolver;
smoother GaussSeidel;
tolerance 1e-06;
relTol 0.1;
nSweeps 1;
}

"(U|k|C)Final"
{
solver smoothSolver;
smoother symGaussSeidel;
tolerance 1e-08;
relTol 0;
}
}

PIMPLE
{
momentumPredictor no;
nCorrectors 2;
nOuterCorrectors 1;
nNonOrthogonalCorrectors 1;

pRefCell 0;
pRefValue 0;
correctPhi yes;

}