[floquation]

In dynamic mesh cases (e.g. interDyMFoam), the following entry is inside the "dynamicMeshDict" file for a dynamicRefineFvMesh:

${FOAM_TUTORIALS}/multiphase/interDyMFoam/ras/damBreakWithObstacle/constant/dynamicMeshDict

    correctFluxes
    (
        (phi none)
        (nHatf none)
        (rhoPhi none)
        (alphaPhi none)
        (ghf none)
    );

If the right value is not "none", the fluxes are corrected by the dynamic mesh library. In OF40 with a dynamicRefineFvMesh, this happens in the "refine()" method of dynamicRefineFvMesh.C.

However, inside the solver (interDyMFoam), the fluxes seem to be corrected as well, unconditionally:

${FOAM_SOLVERS}/multiphase/interFoam/interDyMFoam/interDyMFoam.C

                mesh.update();

                (...)

                if (mesh.changing() && correctPhi)
                {
                    // Calculate absolute flux from the mapped surface velocity
                    phi = mesh.Sf() & Uf;

                    #include "correctPhi.H"

                    // Make the flux relative to the mesh motion
                    fvc::makeRelative(phi, U);

                    mixture.correct();
                }

I suppose they do the same thing? In fact, the solver will overwite anything the dynamic mesh library did with phi due to the line:

phi = mesh.Sf() & Uf;

I have ran a quick rising bubble case with the phi-flux set to both 'none' and 'U'. The result was identical.

Question:
In other words, am I correct when I say that the values in "correctFluxes" may always be 'none'?
Would the same apply to the other fluxes, e.g. alphaPhi? As they also get computed by the solver (as first occurrence):

"solver"/alphaEqn.H
tmp<surfaceScalarField> talphaPhiUD(alpha1Eqn.flux());
alphaPhi = talphaPhiUD();

[magnushaese]

I come to the same conclusion.

phi = mesh.Sf() & Uf;

seems to be somewhat redundant.

The .update() function of dynamicRefineFvMesh executes .refine() and/or .unrefine() which both use the input of correcFluxes from the dynamicMeshDict. For this mesh class the advantage of using the internal flux correction, instead of the line above, seems to be the exclusive treatment of cells that have been modified by the refinement/unrefinement. Those are being marked with the label "-1" and only those fluxes are corrected. See Foam::dynamicRefineFvMesh::refine(const labelList& cellsToRefine):

const labelList& faceMap = map().faceMap();
const labelList& reverseFaceMap = map().reverseFaceMap();

[...]

            // Recalculate new internal faces.
            for (label faceI = 0; faceI < nInternalFaces(); faceI++)
            {
                label oldFaceI = faceMap[faceI];

                if (oldFaceI == -1)
                {
                    // Inflated/appended
                    phi[faceI] = phiU[faceI];
                }
                else if (reverseFaceMap[oldFaceI] != faceI)
                {
                    // face-from-masterface
                    phi[faceI] = phiU[faceI];
                }
            }

[...]

Nevertheless an issue arises in my case. The usage of flux correction does not seem to work.
I want to extend the simpleFoam solver by the dynamicRefineFvMesh class. In short, bold is the extension:

#include "fvCFD.H"
#include "dynamicFvMesh.H"
#include "singlePhaseTransportModel.H"
#include "turbulentTransportModel.H"
#include "simpleControl.H"
#include "fvIOoptionList.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

int main(int argc, char *argv[])
{
    #include "setRootCase.H"
    #include "createTime.H"
    #include "createDynamicFvMesh.H"
    //#include "createMesh.H"

    simpleControl simple(mesh);

    #include "createFields.H"
    #include "createMRF.H"
    #include "createFvOptions.H"
    #include "initContinuityErrs.H"

    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

    Info<< "\nStarting time loop\n" << endl;

simpleTime.start();

    while (simple.loop())
    {
        Info<< "Time = " << runTime.timeName() << nl << endl;
          mesh.update();
          phi = fvc::interpolate(U) & mesh.Sf();
        // --- Pressure-velocity SIMPLE corrector
        {
             UEqnTime.start();
            #include "UEqn.H"
[...]

In run the code on the PitzDaily tutorial and what i get right after the refinement, based on a scalar k interval, is not what i expect. See attachment k_before and k_refined. The U field has small jumps on the refinement edges which i ascribe to the interpolation scheme of U. However, i cannot explain the unphysical jump in the absolute k values. The mesh can still reach more refinement levels but the solution process crashes sooner or later.
The setting for correctFluxes is in DynamicMeshDict is:

    correctFluxes
    (
    (phi none)
    );

Somehow using the right velocity field e.g.

    correctFluxes
    (
    (phi U)
    );

does result into an error:

   
[...]
Selected 1023 cells for refinement out of 12225.
Refined from 12225 to 19386 cells.
#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::sigSegv::sigHandler(int) at ??:?
#2  ? in "/lib/x86_64-linux-gnu/libc.so.6"
#3  Foam::dynamicRefineFvMesh::refine(Foam::List<long> const&) at ??:?
#4  Foam::dynamicRefineFvMesh::update() at ??:?
#5  ? at ??:?
#6  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
#7  ? at ??:?
Segmentation fault (core dumped)

It shows that the error arises in the .refine() function. Debugging it rudimentarily, i think i have found the error-causing lines, which is the last loop of the .refine() function.

            // Update master faces
            forAllConstIter(labelHashSet, masterFaces, iter)
            {
                label faceI = iter.key();

                if (isInternalFace(faceI))
                {
                    phi[faceI] = phiU[faceI];
                }
                else
                {
                    label patchI = boundaryMesh().whichPatch(faceI);
                    label i = faceI - boundaryMesh()[patchI].start();
//Info<< "Breakpoint" << endl;
                    const fvsPatchScalarField& patchPhiU =
                        phiU.boundaryField()[patchI];

                    fvsPatchScalarField& patchPhi = bphi[patchI];

                    patchPhi[i] = patchPhiU[i];
                }

So, two questions arise.
1. Where do these jumps in k come from? Is the dynamicRefineFvMesh implementation wrong?
2. Why cant the correctFlux dictionary with the correct velocity field not be used?

I am thankful for every hint.