[Alczem]

Hey!


After testing a lot of parameters, I wanted to report what worked for me (numerically speaking ) for an overset simulation of a body falling/moving in a fluid domain linked to the movement of the overset mesh, to avoid meshing the whole background, with OpenFoam v2306. This was done to study fluttering and tumbling motion of bodies in a pool.

• I tested overPimpleDyMFoam, overBuoyantPimpleDyMFoam and overRhoPimpleDyMFoam. overRhoPimpleDyMFoam was best suited because it can use the moveMeshOuterCorrectors setting to compute mesh movement for each outer corrector instead of once per timestep (overPimpleDyMFoam ignores this switch, and compute the mesh movement only once per timestep). A constant density was used for the fluid. overBuoyantPimpleDyMFoam was given me trouble with velocity exploding near the falling body, despite underrelaxing or reducing the maxCo.

• nOuterCorrectors = 2 was enough to have a stable simulation, without any spurious velocities, with a maxCo = 3.

• No underrelaxation in fvSolution or dynamicMeshDict was necessary. In particular, accelerationRelax is useful only when using moveMeshOuterCorrectors, and the simulation remains time accurate, but you have to have enough outer correctors. accelerationDamping should not be used to capture accurately the transient evolution of the body or the flow, but only to reach a steady state.

• Newmark or symplectic solver gave similar results with marginal differences. Newmark seemed more robust when using larger CFLs.

• inverseDistance and cellVolumeWeight gave similar results, with inverseDistance being waaaay faster (as reported in several papers and threads) so I would recommend it over cellVolumeWeight except if nothing else works.

• Set nCorrectors at least to 2. The solver would diverge if nCorrectors = 1

• Setting the relaxation factors of fvSolution to 0.8 or 0.7 (even for the final iteration) did not affect significatively the results when using residualControl and the Pimple algorithm to have a converged timestep and allowed for a more stable run, but differences were minimal compared to nOuterCorrectors = 2 and no relaxation.

• The sides of the moving domain were treated as one moving boundary where fluid could go through freely, so inletOutlet was used for turbulence quantities, totalPressure for pressure, and zeroGradient for velocity (pressureInletOutletVelocity was not tested)

Here are the settings:


fvSolution

solvers
{
    cellDisplacement
    {
        solver          PCG;
        preconditioner  DIC;

        tolerance       1e-06;
        relTol          0;
        maxIter         100;
    }
    
    rho
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-8;
        relTol          0.1;
    }

    rhoFinal
    {
        $rho;
        tolerance       1e-8;
        relTol          0;
    }

    "(p|p_rgh)"
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
    }

    "(p|p_rgh)Final"
    {
        $p;
    }

    pcorr
    {
        $p;
        //solver          PCG;
        //preconditioner  DIC;
    }

    pcorrFinal
    {
        $pcorr;
        relTol          0;
    }
    
    "h.*"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-6;
        relTol          0;
        maxIter            0;
    }

    "(U|k|epsilon|omega)"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-6;
        relTol          0;
    }

    "(U|k|epsilon|omega)Final"
    {
        $U;
        tolerance       1e-6;
        relTol          0;
    }

    "yPsi.*"
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-5;
        relTol          0.0;
        minIter            1;
    }

}

PIMPLE
{
    momentumPredictor   no;
    nOuterCorrectors    2;
    nCorrectors         3;
    nNonOrthogonalCorrectors 1;
    moveMeshOuterCorrectors yes;
    turbOnFinalIterOnly no;
    
    consistent            no;

    oversetAdjustPhi    yes; //replaces ddtCorr, correctPhi and massFluxInterpolation
    
     residualControl
    {
        "p.*"
        {
                tolerance  1e-5;
                relTol      0;
        }

     }
}


relaxationFactors
{
    fields
    {
        "p_rgh.*"       1;
    }
    equations
    {
        "U.*"                1;
        "k.*"                 1;
        "epsilon.*"         1;
        "omega.*"             1;
    }
}

cache
{
    grad(U);
}

fvSchemes

ddtSchemes
{
    default         Euler;
    //default         backward;
}

gradSchemes
{
    default         Gauss linear;
    limited            cellLimited Gauss linear 0.5;
}

divSchemes
{
    default         none;
    div(phi,U)      Gauss linearUpwindV limited; //Gauss upwind;

    div(phi,epsilon) Gauss linearUpwind limited; //Gauss upwind;
    div(phi,omega) Gauss linearUpwind limited; //Gauss upwind;
    div(phi,k)       Gauss linearUpwind limited; //Gauss upwind;
    div(phi,h)       Gauss linearUpwind limited;
    
    div(phi,K)      Gauss linear;

    div((nuEff*dev2(T(grad(U))))) Gauss linear;
    div(((rho*nuEff)*dev2(T(grad(U)))))    Gauss linear;
    
    div(meshPhi,p)  Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear corrected;
    laplacian(diffusivity,cellDisplacement)  Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected ;
}

oversetInterpolation
{
    //method          cellVolumeWeight;
    method          inverseDistance;
    //method          leastSquares;

    // The inverseDistance method uses a 'voxel' like search structure.
    // Optionally specify the extent and number of divisions n.
    // Note that it will allocate an array of nx*ny*nz. If not specified:
    // - searchBox          : local mesh bounding box
    // - searchBoxDivisions : root (2D) or cube-root(3D) of number of cells
    //searchBox           (-0.1 -0.175 -0.175) (0.55 0.175 0.175);
    //searchBoxDivisions  (10 10 10);

    //layerRelax 0.5; //control the number of interpolated cells
    
    //allowInterpolatedDonors    false;
    //layerRelax          0.3;
    //holeLayers            3;
    //useLayer            2;// Layer used
    
    
    
    // Faster but less accurate
    //method              trackingInverseDistance;
    //searchBox           (0 0 0)(0.02 0.01 0.01);
    //searchBoxDivisions  2{(64 64 64)};
    //allowInterpolatedDonors     false;
    
}


oversetInterpolationSuppressed
{
    grad(p_rgh);
    surfaceIntegrate(phiHbyA);
}


fluxRequired
{
    default         no;
    pcorr           ;
    p               ;
}


wallDist
{
    method        Poisson;
    //correctWalls    true;
}

dynamicMeshDict

motionSolverLibs    (sixDoFRigidBodyMotion);

dynamicFvMesh       dynamicOversetFvMesh;

solver          sixDoFRigidBodyMotion;



solvers
{
    
     domain
     {
        motionSolverLibs (fvMotionSolvers);

        motionSolver            solidBody;
        solidBodyMotionFunction drivenLinearMotion;

        cellSet                 c0;
        cOfGdisplacement        CofG;
    }
    



    plate
    {
        motionSolverLibs   (sixDoFRigidBodyMotion);
        motionSolver        sixDoFRigidBodyMotion;

        cellSet         c1;
        cOfGdisplacement    CofG;

        patches         (plate);
        innerDistance   1e3;
        outerDistance   1.5e3;

        centreOfMass    (0 0 20);

        //rho rhoInf;

        //rhoInf 1;
        

        // Cuboid mass
        mass            21.25;  //mass minus the buoyancy of the body
        
        g (0 0 -9.81);
        
        velocity        (0 0 -5);
        
        momentOfInertia (1.40625 0.71015 0.71015);

        report          on;
        accelerationRelaxation 1; //0.5; // for newmark
        accelerationDamping    1; //0.8;

        solver
        {
            type Newmark;
        }

    
        constraints
        {
            /*
            line_restraint
            {
               sixDoFRigidBodyMotionConstraint   line;
               direction                              (0 0 1);
            }
            
            orientation_restraint
            {
                sixDoFRigidBodyMotionConstraint   orientation;
            }
            */
            
        }
    
    }
}

Sorry for the messy post but I wanted to have some written record as much for me as for others. This is all based on my recent experience with overset simulations and whatever I could find online, so take it with a grain of salt!

I will be happy to hear from others about recommended practices overset computations and moving domains. There are a lot of settings I have not yet tested, such as the ddtScheme. Hopefully I will update this thread with more detailed information as I learn more about it.

[NotOverUnderated]

Hi,

Thank you for sharing this. This is certainly very helpful.

I have some questions about some dictionaries in your fvSolution files:

1) Why you are using a small number for nOuterCorrectors? I read somewhere here in the forum (here: https://openfoamwiki.net/index.php/O...hm_in_OpenFOAM) that you must use a large number.

    nOuterCorrectors    2; // I think this should be a large number like 30

2) This is a followup question from the question 1):

    nOuterCorrectors    2;

    ...
    ...
    
     residualControl
    {
        "p.*"
        {
                tolerance  1e-5;
                relTol      0;
        }

     }

- What is the point of using residualControl in this case when you are using only 2 nOuterCorrectors? as far as I understand the residualControl dictionary controls the number of outer correctors: The solver will do nOuterCorrector iterations unless the initial residuals of the fields satisfy the tolerance specified in residualControl.

Kind regards

[xiangxiang]

Originally Posted by Alczem

Hey!


After testing a lot of parameters, I wanted to report what worked for me (numerically speaking ) for an overset simulation of a body falling/moving in a fluid domain linked to the movement of the overset mesh, to avoid meshing the whole background, with OpenFoam v2306. This was done to study fluttering and tumbling motion of bodies in a pool.

• I tested overPimpleDyMFoam, overBuoyantPimpleDyMFoam and overRhoPimpleDyMFoam. overRhoPimpleDyMFoam was best suited because it can use the moveMeshOuterCorrectors setting to compute mesh movement for each outer corrector instead of once per timestep (overPimpleDyMFoam ignores this switch, and compute the mesh movement only once per timestep). A constant density was used for the fluid. overBuoyantPimpleDyMFoam was given me trouble with velocity exploding near the falling body, despite underrelaxing or reducing the maxCo.

• nOuterCorrectors = 2 was enough to have a stable simulation, without any spurious velocities, with a maxCo = 3.

• No underrelaxation in fvSolution or dynamicMeshDict was necessary. In particular, accelerationRelax is useful only when using moveMeshOuterCorrectors, and the simulation remains time accurate, but you have to have enough outer correctors. accelerationDamping should not be used to capture accurately the transient evolution of the body or the flow, but only to reach a steady state.

• Newmark or symplectic solver gave similar results with marginal differences. Newmark seemed more robust when using larger CFLs.

• inverseDistance and cellVolumeWeight gave similar results, with inverseDistance being waaaay faster (as reported in several papers and threads) so I would recommend it over cellVolumeWeight except if nothing else works.

• Set nCorrectors at least to 2. The solver would diverge if nCorrectors = 1

• Setting the relaxation factors of fvSolution to 0.8 or 0.7 (even for the final iteration) did not affect significatively the results when using residualControl and the Pimple algorithm to have a converged timestep and allowed for a more stable run, but differences were minimal compared to nOuterCorrectors = 2 and no relaxation.

• The sides of the moving domain were treated as one moving boundary where fluid could go through freely, so inletOutlet was used for turbulence quantities, totalPressure for pressure, and zeroGradient for velocity (pressureInletOutletVelocity was not tested)

Here are the settings:


fvSolution

solvers
{
    cellDisplacement
    {
        solver          PCG;
        preconditioner  DIC;

        tolerance       1e-06;
        relTol          0;
        maxIter         100;
    }
    
    rho
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-8;
        relTol          0.1;
    }

    rhoFinal
    {
        $rho;
        tolerance       1e-8;
        relTol          0;
    }

    "(p|p_rgh)"
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
    }

    "(p|p_rgh)Final"
    {
        $p;
    }

    pcorr
    {
        $p;
        //solver          PCG;
        //preconditioner  DIC;
    }

    pcorrFinal
    {
        $pcorr;
        relTol          0;
    }
    
    "h.*"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-6;
        relTol          0;
        maxIter            0;
    }

    "(U|k|epsilon|omega)"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-6;
        relTol          0;
    }

    "(U|k|epsilon|omega)Final"
    {
        $U;
        tolerance       1e-6;
        relTol          0;
    }

    "yPsi.*"
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-5;
        relTol          0.0;
        minIter            1;
    }

}

PIMPLE
{
    momentumPredictor   no;
    nOuterCorrectors    2;
    nCorrectors         3;
    nNonOrthogonalCorrectors 1;
    moveMeshOuterCorrectors yes;
    turbOnFinalIterOnly no;
    
    consistent            no;

    oversetAdjustPhi    yes; //replaces ddtCorr, correctPhi and massFluxInterpolation
    
     residualControl
    {
        "p.*"
        {
                tolerance  1e-5;
                relTol      0;
        }

     }
}


relaxationFactors
{
    fields
    {
        "p_rgh.*"       1;
    }
    equations
    {
        "U.*"                1;
        "k.*"                 1;
        "epsilon.*"         1;
        "omega.*"             1;
    }
}

cache
{
    grad(U);
}

fvSchemes

ddtSchemes
{
    default         Euler;
    //default         backward;
}

gradSchemes
{
    default         Gauss linear;
    limited            cellLimited Gauss linear 0.5;
}

divSchemes
{
    default         none;
    div(phi,U)      Gauss linearUpwindV limited; //Gauss upwind;

    div(phi,epsilon) Gauss linearUpwind limited; //Gauss upwind;
    div(phi,omega) Gauss linearUpwind limited; //Gauss upwind;
    div(phi,k)       Gauss linearUpwind limited; //Gauss upwind;
    div(phi,h)       Gauss linearUpwind limited;
    
    div(phi,K)      Gauss linear;

    div((nuEff*dev2(T(grad(U))))) Gauss linear;
    div(((rho*nuEff)*dev2(T(grad(U)))))    Gauss linear;
    
    div(meshPhi,p)  Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear corrected;
    laplacian(diffusivity,cellDisplacement)  Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected ;
}

oversetInterpolation
{
    //method          cellVolumeWeight;
    method          inverseDistance;
    //method          leastSquares;

    // The inverseDistance method uses a 'voxel' like search structure.
    // Optionally specify the extent and number of divisions n.
    // Note that it will allocate an array of nx*ny*nz. If not specified:
    // - searchBox          : local mesh bounding box
    // - searchBoxDivisions : root (2D) or cube-root(3D) of number of cells
    //searchBox           (-0.1 -0.175 -0.175) (0.55 0.175 0.175);
    //searchBoxDivisions  (10 10 10);

    //layerRelax 0.5; //control the number of interpolated cells
    
    //allowInterpolatedDonors    false;
    //layerRelax          0.3;
    //holeLayers            3;
    //useLayer            2;// Layer used
    
    
    
    // Faster but less accurate
    //method              trackingInverseDistance;
    //searchBox           (0 0 0)(0.02 0.01 0.01);
    //searchBoxDivisions  2{(64 64 64)};
    //allowInterpolatedDonors     false;
    
}


oversetInterpolationSuppressed
{
    grad(p_rgh);
    surfaceIntegrate(phiHbyA);
}


fluxRequired
{
    default         no;
    pcorr           ;
    p               ;
}


wallDist
{
    method        Poisson;
    //correctWalls    true;
}

dynamicMeshDict

motionSolverLibs    (sixDoFRigidBodyMotion);

dynamicFvMesh       dynamicOversetFvMesh;

solver          sixDoFRigidBodyMotion;



solvers
{
    
     domain
     {
        motionSolverLibs (fvMotionSolvers);

        motionSolver            solidBody;
        solidBodyMotionFunction drivenLinearMotion;

        cellSet                 c0;
        cOfGdisplacement        CofG;
    }
    



    plate
    {
        motionSolverLibs   (sixDoFRigidBodyMotion);
        motionSolver        sixDoFRigidBodyMotion;

        cellSet         c1;
        cOfGdisplacement    CofG;

        patches         (plate);
        innerDistance   1e3;
        outerDistance   1.5e3;

        centreOfMass    (0 0 20);

        //rho rhoInf;

        //rhoInf 1;
        

        // Cuboid mass
        mass            21.25;  //mass minus the buoyancy of the body
        
        g (0 0 -9.81);
        
        velocity        (0 0 -5);
        
        momentOfInertia (1.40625 0.71015 0.71015);

        report          on;
        accelerationRelaxation 1; //0.5; // for newmark
        accelerationDamping    1; //0.8;

        solver
        {
            type Newmark;
        }

    
        constraints
        {
            /*
            line_restraint
            {
               sixDoFRigidBodyMotionConstraint   line;
               direction                              (0 0 1);
            }
            
            orientation_restraint
            {
                sixDoFRigidBodyMotionConstraint   orientation;
            }
            */
            
        }
    
    }
}

Sorry for the messy post but I wanted to have some written record as much for me as for others. This is all based on my recent experience with overset simulations and whatever I could find online, so take it with a grain of salt!

I will be happy to hear from others about recommended practices overset computations and moving domains. There are a lot of settings I have not yet tested, such as the ddtScheme. Hopefully I will update this thread with more detailed information as I learn more about it.

Hi, Alczem
Here's my personal opinion:
1. Euler in ddtSchemes is of first-order precision, and usually we prefer to get a result with second-order precision, such as using backward
2. The number of external loops of nOuterCorrectors is usually set to be relatively large, and I usually set it to nOuterCorrectors =50~100, which does not affect the calculation speed, because the external loops will jump out after the convergence condition is reached.

Kind regards
xiangxiang

[Alczem]

Hey!


About the nOuterCorrectors set to 2, I found out that letting the solver converge with nOuterCorrectors set to 50 and residualControl versus having a fixed number of loops (2 or 4 in my case) made almost no difference in the computed position and velocity after 1 second of simulated time. I should add that I found a a tutorial by Alletto, investigating the Newmark and symplectic solvers for oversets, and it also uses a fixed number of loops for Newmark, which gave me confidence in using 2 loops. Here is his collection of tutorials:


https://wiki.openfoam.com/Collection...ichael_Alletto



You are right that using PIMPLE and residualControl would usually result in a few iterations to converge (more are needed with a tighter tolerance of course). But since I used moveMeshOuterCorrectors, additional loops are more expansive than usual. That is why I stuck with 2 loops, but I agree that if you need the best accuracy regardless of runtime, having PIMPLE manage the number of loops is the way! The residualControl subdictionary in my fvSolution is just there so I didn't have to re-edit the file when I was trying things out. It does not prevent the simulation from running even if it does not converge each timestep.



I have not properly tested the ddtSchemes, I hope I will have some time to compare the different schemes. But similar to the nOuterCorrectors discussion, the settings I shared are meant to be a compromise between accuracy and cost, since overset simulations are so heavy


(As additional info, my test case was a plate falling with an initial vertical velocity in water, and I monitored the position and velocity after 1 second. This may not be chaotic or long enough to see the differences in accuracy when using more loops or different schemes though)

[Alczem]

Quick edit about the velocity boundary condition for the domain boundaries:


I used zeroGradient, but running the simulation for a longer time eventually created some kind of succion effect at the "inlet" (the bottom boundary since my body is falling). My understanding is that the fluid getting pushed downward by the falling body reaches the zeroGradient boundary, which registers this downward current and keep applying this downward velocity, resulting in a succion instead of an immobile fluid. As a a fix, I tried to enlarge my background mesh, with better results, but whats seems to be best is to use an inletOutlet condition with inletValue set to (0 0 0) for all the boundaries that allow the fluid to enter or leave freely. That way, the fluid is introduced into the domain with a null velocity (duh, I don't know how I missed that, it is also similar to the rigidBodyHull tutorial for the interFoam overset solver).


With this change, results are similar between a cube falling with a fixed background mesh and a cube falling with a moving domain linked to it (terminal velocity and forces are equal), though the background mesh for the moving domain still needs to be large enough to mitigate any boundary effect on the moving body.

[Karim ahmed]

Hello Alczem ,

Thank for this useful informations. can you please provide the 0 folder , just to check the BC , i am trying to solve something near this case and i want to check my BCs

thanks in advance

[Alczem]

Originally Posted by Karim ahmed

Hello Alczem ,

Thank for this useful informations. can you please provide the 0 folder , just to check the BC , i am trying to solve something near this case and i want to check my BCs

thanks in advance

Hey,


Here you go! Don't hesitate to report back if you find issues with this setup

[Redrakham]

Originally Posted by Alczem

Hey!


After testing a lot of parameters, I wanted to report what worked for me (numerically speaking ) for an overset simulation of a body falling/moving in a fluid domain linked to the movement of the overset mesh, to avoid meshing the whole background, with OpenFoam v2306. This was done to study fluttering and tumbling motion of bodies in a pool.

• I tested overPimpleDyMFoam, overBuoyantPimpleDyMFoam and overRhoPimpleDyMFoam. overRhoPimpleDyMFoam was best suited because it can use the moveMeshOuterCorrectors setting to compute mesh movement for each outer corrector instead of once per timestep (overPimpleDyMFoam ignores this switch, and compute the mesh movement only once per timestep). A constant density was used for the fluid. overBuoyantPimpleDyMFoam was given me trouble with velocity exploding near the falling body, despite underrelaxing or reducing the maxCo.

• nOuterCorrectors = 2 was enough to have a stable simulation, without any spurious velocities, with a maxCo = 3.

• No underrelaxation in fvSolution or dynamicMeshDict was necessary. In particular, accelerationRelax is useful only when using moveMeshOuterCorrectors, and the simulation remains time accurate, but you have to have enough outer correctors. accelerationDamping should not be used to capture accurately the transient evolution of the body or the flow, but only to reach a steady state.

• Newmark or symplectic solver gave similar results with marginal differences. Newmark seemed more robust when using larger CFLs.

• inverseDistance and cellVolumeWeight gave similar results, with inverseDistance being waaaay faster (as reported in several papers and threads) so I would recommend it over cellVolumeWeight except if nothing else works.

• Set nCorrectors at least to 2. The solver would diverge if nCorrectors = 1

• Setting the relaxation factors of fvSolution to 0.8 or 0.7 (even for the final iteration) did not affect significatively the results when using residualControl and the Pimple algorithm to have a converged timestep and allowed for a more stable run, but differences were minimal compared to nOuterCorrectors = 2 and no relaxation.

• The sides of the moving domain were treated as one moving boundary where fluid could go through freely, so inletOutlet was used for turbulence quantities, totalPressure for pressure, and zeroGradient for velocity (pressureInletOutletVelocity was not tested)

Here are the settings:


fvSolution

solvers
{
    cellDisplacement
    {
        solver          PCG;
        preconditioner  DIC;

        tolerance       1e-06;
        relTol          0;
        maxIter         100;
    }
    
    rho
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-8;
        relTol          0.1;
    }

    rhoFinal
    {
        $rho;
        tolerance       1e-8;
        relTol          0;
    }

    "(p|p_rgh)"
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0;
    }

    "(p|p_rgh)Final"
    {
        $p;
    }

    pcorr
    {
        $p;
        //solver          PCG;
        //preconditioner  DIC;
    }

    pcorrFinal
    {
        $pcorr;
        relTol          0;
    }
    
    "h.*"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-6;
        relTol          0;
        maxIter            0;
    }

    "(U|k|epsilon|omega)"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-6;
        relTol          0;
    }

    "(U|k|epsilon|omega)Final"
    {
        $U;
        tolerance       1e-6;
        relTol          0;
    }

    "yPsi.*"
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-5;
        relTol          0.0;
        minIter            1;
    }

}

PIMPLE
{
    momentumPredictor   no;
    nOuterCorrectors    2;
    nCorrectors         3;
    nNonOrthogonalCorrectors 1;
    moveMeshOuterCorrectors yes;
    turbOnFinalIterOnly no;
    
    consistent            no;

    oversetAdjustPhi    yes; //replaces ddtCorr, correctPhi and massFluxInterpolation
    
     residualControl
    {
        "p.*"
        {
                tolerance  1e-5;
                relTol      0;
        }

     }
}


relaxationFactors
{
    fields
    {
        "p_rgh.*"       1;
    }
    equations
    {
        "U.*"                1;
        "k.*"                 1;
        "epsilon.*"         1;
        "omega.*"             1;
    }
}

cache
{
    grad(U);
}

fvSchemes

ddtSchemes
{
    default         Euler;
    //default         backward;
}

gradSchemes
{
    default         Gauss linear;
    limited            cellLimited Gauss linear 0.5;
}

divSchemes
{
    default         none;
    div(phi,U)      Gauss linearUpwindV limited; //Gauss upwind;

    div(phi,epsilon) Gauss linearUpwind limited; //Gauss upwind;
    div(phi,omega) Gauss linearUpwind limited; //Gauss upwind;
    div(phi,k)       Gauss linearUpwind limited; //Gauss upwind;
    div(phi,h)       Gauss linearUpwind limited;
    
    div(phi,K)      Gauss linear;

    div((nuEff*dev2(T(grad(U))))) Gauss linear;
    div(((rho*nuEff)*dev2(T(grad(U)))))    Gauss linear;
    
    div(meshPhi,p)  Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear corrected;
    laplacian(diffusivity,cellDisplacement)  Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected ;
}

oversetInterpolation
{
    //method          cellVolumeWeight;
    method          inverseDistance;
    //method          leastSquares;

    // The inverseDistance method uses a 'voxel' like search structure.
    // Optionally specify the extent and number of divisions n.
    // Note that it will allocate an array of nx*ny*nz. If not specified:
    // - searchBox          : local mesh bounding box
    // - searchBoxDivisions : root (2D) or cube-root(3D) of number of cells
    //searchBox           (-0.1 -0.175 -0.175) (0.55 0.175 0.175);
    //searchBoxDivisions  (10 10 10);

    //layerRelax 0.5; //control the number of interpolated cells
    
    //allowInterpolatedDonors    false;
    //layerRelax          0.3;
    //holeLayers            3;
    //useLayer            2;// Layer used
    
    
    
    // Faster but less accurate
    //method              trackingInverseDistance;
    //searchBox           (0 0 0)(0.02 0.01 0.01);
    //searchBoxDivisions  2{(64 64 64)};
    //allowInterpolatedDonors     false;
    
}


oversetInterpolationSuppressed
{
    grad(p_rgh);
    surfaceIntegrate(phiHbyA);
}


fluxRequired
{
    default         no;
    pcorr           ;
    p               ;
}


wallDist
{
    method        Poisson;
    //correctWalls    true;
}

dynamicMeshDict

motionSolverLibs    (sixDoFRigidBodyMotion);

dynamicFvMesh       dynamicOversetFvMesh;

solver          sixDoFRigidBodyMotion;



solvers
{
    
     domain
     {
        motionSolverLibs (fvMotionSolvers);

        motionSolver            solidBody;
        solidBodyMotionFunction drivenLinearMotion;

        cellSet                 c0;
        cOfGdisplacement        CofG;
    }
    



    plate
    {
        motionSolverLibs   (sixDoFRigidBodyMotion);
        motionSolver        sixDoFRigidBodyMotion;

        cellSet         c1;
        cOfGdisplacement    CofG;

        patches         (plate);
        innerDistance   1e3;
        outerDistance   1.5e3;

        centreOfMass    (0 0 20);

        //rho rhoInf;

        //rhoInf 1;
        

        // Cuboid mass
        mass            21.25;  //mass minus the buoyancy of the body
        
        g (0 0 -9.81);
        
        velocity        (0 0 -5);
        
        momentOfInertia (1.40625 0.71015 0.71015);

        report          on;
        accelerationRelaxation 1; //0.5; // for newmark
        accelerationDamping    1; //0.8;

        solver
        {
            type Newmark;
        }

    
        constraints
        {
            /*
            line_restraint
            {
               sixDoFRigidBodyMotionConstraint   line;
               direction                              (0 0 1);
            }
            
            orientation_restraint
            {
                sixDoFRigidBodyMotionConstraint   orientation;
            }
            */
            
        }
    
    }
}

Sorry for the messy post but I wanted to have some written record as much for me as for others. This is all based on my recent experience with overset simulations and whatever I could find online, so take it with a grain of salt!

I will be happy to hear from others about recommended practices overset computations and moving domains. There are a lot of settings I have not yet tested, such as the ddtScheme. Hopefully I will update this thread with more detailed information as I learn more about it.

I am doing similar overset problem of Vortex induced vibration of cylinder in 3D. But I am using following part:
PIMPLE
{
momentumPredictor true;
correctPhi true; //true;
//oversetAdjustPhi false;
nOuterCorrectors 50;
nCorrectors 5;
nNonOrthogonalCorrectors 2;

//ddtCorr false;

//checkMeshCourantNo yes;
pRefCell 0;
pRefValue 0;
moveMeshOuterCorrector yes;
residualControl
{
p
{
relTol 0;
tolerance 1e-5;
}

U
{
relTol 0;
tolerance 1e-5;
}
}
}
After every 50 outer iterations it is not converged. Is it something to do with momentumPredictor? Thanks in advance

[Alczem]

Hey!


Could you share the case? It may have to do with your max Courant number, the mesh or maybe the schemes. People can have a look and maybe find out the issue.

[Redrakham]

The case is attached.

[Jim1310]

Hello Alczem!


Very interesting what you have done with the moving background. Id like to ask you a few questions:
1) Is it possible to do this with a deforming background mesh (using for e.g. displacementLaplacian)? I Have tried it but i got an error that was saying:

"--> FOAM FATAL ERROR: (openfoam-2306)
cannot be called for a calculatedFvPatchField
on patch floatingObject of field cellDisplacement in file "/home/dimitris/OpenFOAM/dimitris-v2106/run/tutorials/multiphase/overInterDyMFoam/floatingBody2D_flap/background/0/cellDisplacement"
You are probably trying to solve for a field with a default boundary condition."


2) Could you also provide information about the boundary conditions of the fields especially those of poindDisplacement?


Thank you very much!!

[Alczem]

Hey!


1) My guess is that since you are solving a 2D case, the background patch should have an empty condition in the pointDisplacement file. But you also need the calculated condition. To achieve that, you should have two entries for background, like this:

background
{
type calculated;
patchType empty;

}

Not sure about this, but that is what I understand from the error message.



2) So overset boundaries should be set to type zeroGradient and patchType overset. Moving boundaries should be set to calculated, and fixed boundaries should be set to fixedValue.

[Jim1310]

Hello Alczem and thank you for your quick answer!


I am not sure that i understand to which patch is the background patch referred to. I have 4 active patches for background (inlet, outlet, top, bottom) and two empty ones (i do solve in 2D). For the overset mesh i have a wall for the floating object which is of "type calculated" and an overset patchType, with type zeroGradient.


What i undertand is that the laplacian solver is trying to solve the laplace equation for the entire mesh (including the overset mesh) and doesn't know what to do with the boundary condition of the floatingObject. I have tried to restrict the action of the plaplacian solver by providing a cellSet (corresponding to the cells of the bacground mesh) or by applying a frozenPpointsZone, but to no avail...

[Alczem]

Hey!


My mistake about the background patch, I read the error message a little bit too fast, the path name got me confused


At this point, share your files if you can. I still think the error comes from a misused boundary condition, but it could also be the fact that you cannot use a deforming mesh with overset meshes (I neved did, maybe someone else can chime in).

[Jim1310]

Hello!



Yes of course here is the case attached. I run it with v2306, and using a modified kind of overset mesh which allows you to use various solvers for each cellZone. You can find the code from this post: Ahmed et al. - 2024 - A Modified Overset Method in OpenFOAM for Simultaneous Motion..


In the attached case i have tried various things and it is a little messy. I have reached the following conclusions:


-I am using the displacementLaplacian solver for a cellZone corresponding to all the cells of the background mesh and applying a uniform diffusivity for the mesh (at some point i was using inverseDistance diffusivity which lead to various problems. To partially solve them i was using a part of the background which includes all the cells from the left up to x=-0.2m. That led the to crushes, since there are some cells (the ones in location x=-0.2m) which are getting squished. A way to avoid this squishing is to include all the cells of the background mesh to the cellZone used for the displacementLaplacian solver, which would again lead to problems because of the diffusivity around the body).


-The main problem is the boundary condition of the pointDisplacement at the boundary of the body. If we use the sixDoFRigidBodyMotion for the deformation of the overset mesh we have to use the B.C. 'calculated' there. But this B.C. is not working (see also this post for 'calculated' how to use the 'calculated' boundary condition). By inspection displacementLaplacian code we see that the solver creates a field named cellDisplacement which takes exactly the same boundary conditions as pointDisplacement. This is where we think is the problem. Although the displacement is applied to the background mesh it is previously calculated for the entire mesh. Therefore, during the solution of cellDisplacement there is no information for what to do with the 'type calculated' boundary. We have tested a different (prescribed) boundary condition instead of calculated and it seems to work. Specifically, we use an oscillating rotation for the moving body. We also use displacementLaplacian again for the deformation of the overset mesh. Maybe a different solver which would apply rigid motion would be better but i dont know at the moment. A remark here is that the boundaries of the overset mesh are not moving even though they are set to zeroGradient in pointDisplacement.
-Maybe a solution would be to create a modified version of claculated B.C. where the functions valueInernalCoeffs, valueBoundaryCoeffs, gradientInternalCoeffs, gradientBoundaryCoeffs are implemented (instead of throwing the error: "You are probably trying to solve for a field with a default boundary condition") in something like zeroGradient, so that the solver for cellDisplacement knows what to do (even though it doesnt apply the deformation ultimately).

[Jim1310]

It seems that modification of the 'calculated' boundary condition is the way to go! If you look at /src/fvMotionSolver/displacement/laplacian/displacementLaplacianfvMotionSolver.C you will se the definition of cellDisplacement:

cellDisplacement_
    (
        IOobject
        (
            "cellDisplacement",
            mesh.time().timeName(),
            mesh,
            IOobject::READ_IF_PRESENT,
            IOobject::AUTO_WRITE
        ),
        fvMesh_,
        dimensionedVector(pointDisplacement_.dimensions(), Zero),
        cellMotionBoundaryTypes<vector>(pointDisplacement_.boundaryField())
    )

We also see that cellDisplacement (cD) is a volVectorField whereas pointDisplacement (pD) is a pointVectorField. Therefore, the boundary conditions of pD (which are of type pointPatchField since we have a point<Type>Field) are redefined for cD using their fvPatchField counterpart (fvPatchField is for vol<Type>Field since cD is a volVectorField). Therefore, we had to reimplement the calculated boundary condition for the one located in /src/finiteVolume/fields/fvPatchFields/basic/calculated (which is of type fvPatchField) and the one located in /src/OpenFOAM/fields/pointPatchFields/basic/calculated (which is of type fvPointPatchField). The boundary condition for cD was named (TypeName) myCalculated, inherits from claculated and the only change is that the functions valueInernalCoeffs, valueBoundaryCoeffs, gradientInternalCoeffs, gradientBoundaryCoeffs are implemented. We have chosen the implementation of zeroGradient boundary condition. For the calculated of pD we didn't change anything except from the name which was set to myCalculated (it is important to have the same name). We include the library of the boundary conditions in the controlDict and that's it! This way the simulation run successfully and i am attaching the case below.