[<nico>]

Hello everyone,

I’m still quite new to OpenFoam, so I guess this might be a very basic question regarding pressure values. There are a lot of threads regarding these pressure topics but there is still some confusion on my end.

In my case I want to generate pump curves using the MRF approach. To compare the results I have the experimental curves available.

When searching the forum, a typical combination for pressure & velocity boundaries seems to be the following:
- Pressure inlet: totalPressure
- Pressure outlet: fixedValue
- Velocity inlet: pressureInletVelocity
- Velocity outlet: zeroGradient

But every time I set the outlet fixed value to a higher pressure than the inlet (to simulate the pump curves pressure differences), the solver crashes and I have a back flow.
Find attached also a post processing picture from Paraview (inlet totalPressure = 0; outlet fixedValue 0). This combination converges.

In one thread I read, that for MRF it is not possible to set a higher outlet pressure.
When this is correct, how do I apply the measured pressure differences from the experimental results to my simulation?

Thanks in advance.
Nico
IMG_1448.jpg

[tomf]

Hi Nico,

There may be an issue with the initialization in your setup. When you start, there is no pressure increase from your pump yet, meaning that flow will start from high to low pressure (outlet => inlet). It may work to start at a lower pressure outlet and then ramp up the outlet pressure to your value after some X number of iterations.

Something like using a table type function 1 may work (please modify the names and numbers to your liking):

pressureOutletPatch
{
    
    type          uniformFixedValue;
    uniformValue  table( ( 0 0)
                         ( 100 0)
                         (200 pBackPressure)
                         ( 500000 pBackPressure)
                     );
    value        uniform 0;
}

This means that from iteration 0-100 the outlet pressure is 0 (like in your converged example), between iteration 100 and 200 it will linearly ramp up the outlet pressure to the value you would want and then keep it there.

See also uniformFixedValue.

Regards,
Tom

[<nico>]

Hi Tom,

that definitely makes sense!
I will try your approach, thank you.

In the meantime I tried a different couple of boundary conditions:
- Pressure inlet: totalPressure
- Pressure outlet: zeroGradient
- Velocity inlet: pressureInletVelocity
- Velocity outlet: flowRateInletVelocity (with negative flow rate)

This combination was basically working. I just had a solver crash and my time step continuity error exploded. But I guess this might be due to my mesh quality.
Do you know if this boundary approach is as suitable as the first one for this case?

[tomf]

Hi Nico,

I think it would make more sense to specify the flow rate at the inlet, with zeroGradient pressure and then use fixedValue for pressure at the outlet with zeroGradient on velocity. In this case you are recreating the fan curve by specifying your flow rates and calculating the pressure jump rather than specifying the pressure jump and calculating the flow rate. Your setup may have resulted in incompatible boundary conditions leading to the continuity error.

Regards,
Tom

[<nico>]

Hi Tom,

when setting a fixed inlet flow (random flow from flow curve) and a fixed pressure at the outlet of 0, the solver crashes between 700 and 800 iterations. See also attached flow field state. („Pump_MRF_fixedFlow_0m2/s2_700iter“)
How can the pressure at the outlet raise above my defined outlet condition of 0 m2/s2?

Furthermore, I tried your other idea and ramped the outlet pressure after 100 iterations. But also here the solver crashes. I already raised the ramp to 300 iterations, but it is the same behavior.
See attached also the flow fields for 0 m2/s2 at 100 iterations and 30 m2/s2 at 300 iterations.

Here some more ideas / questions from my side:
- May I try to raise the ramp even further and wait some more iterations because the flow needs to be more converged for the rotor region because it is still the time step continuity error that explodes?
- When seeing pressure regions below 100 m2/s2, is this due to numerical diffusion? Because for water density this would mean below absolute vacuum.
- This is not the final geometry I want to simulate. The real one is a rotor stator geometry and therefore a bit more complex (here it is just the rotor that I simulate). I want to use overset meshes and steady state MRF (or SRF). Is this combination possible in general and does anything change for my inlet / outlet boundaries?

Thank you so much for your support!Pump_MRF_fixedFlow_0m2s2_700iter.jpgPump_MRF_ramp_0m2s2_100iter.jpgPump_MRF_ramp_30m2s2_300iter.jpg

[tomf]

Hi Nico,

For this type of analysis I would not recommend SRF. MRF could work, but there may be issues at the transition from rotating to stationary parts of the mesh.

It could well be that you need more iterations before starting the ramping up and maybe it needs to be slower.

Numerical diffusion would typically reduce extremes of pressure or velocity, so I do not think that this could be the case.

Overset would also not be necessary. as the rotating and stagnant parts cannot overlap, otherwise the pump cannot work. Can you show which part of your mesh is the MRF zone?

Regards,
Tom

[unilord]

Originally Posted by tomf

Hi Nico,

For this type of analysis I would not recommend SRF. MRF could work, but there may be issues at the transition from rotating to stationary parts of the mesh.

It could well be that you need more iterations before starting the ramping up and maybe it needs to be slower.

Numerical diffusion would typically reduce extremes of pressure or velocity, so I do not think that this could be the case.

Overset would also not be necessary. as the rotating and stagnant parts cannot overlap, otherwise the pump cannot work. Can you show which part of your mesh is the MRF zone?

Regards,
Tom

Hey Tom,

I was wondering if you are familiar with the new Non Conformal Coupling "function" of openfoam. It is a replacement to the previous AMI. I am performing a steady state centrifugal pump simulation with kEpsilon turbulence model. I am using NCC in a steady state because I am meshing with salome and adding layers with snappy, and I want to declare the MRF impeller region the whole impeller mesh from salome.

When I perform the command "createNonConformalCouples", I get this output:

patchToPatch: Calculating couplings between 1920 source faces and 2960 target faces
    Source min/average/max coverage = 0.847949/0.998741/1
    Target min/average/max coverage = 0.757311/0.998741/1
    Source average openness/error/depth/angle = 1.18631e-16/0.000519423/0.00256177/0.0791148
    Source max openness/error/depth/angle = 5.33544e-16/0.00632648/0.00492436/0.261164
    10508 couplings calculated in 0.257038s
patchToPatch: Calculating couplings between 6992 source faces and 6992 target faces
    Source min/average/max coverage = 0.999725/0.999999/1
    Target min/average/max coverage = 0.999725/0.999999/1
    Source average openness/error/depth/angle = 1.41142e-16/3.72745e-05/0.000239833/0.00545699
    Source max openness/error/depth/angle = 7.4155e-16/0.00185617/0.00170314/0.0817809
    41086 couplings calculated in 0.336095s
fvMeshStitcher: Connecting
    0/1 small couplings removed/added to nonConformalCouple_23_on_NCC2_3
    Cell min/avg/max openness = 0/3.37538e-17/2.51584e-16

As you can see, the source and target min and max coverages are not 1. There were times that I got 1 (when using tetrahedral mesh with big refinement in the NCC interface).

My question is, how does this affect the simulation, and do we have to change the nonConformalCouple matchTolerance (as seen below) that is created in the constant/polymesh/boundary file after using "createNonConformalCouples" command according to that error? I am attaching some countours of the results (achieved convergence with max residual of 10^-6) so you can see that the coupling is causing me some "errors" or "unphysical" behaviour near the NCC regions.

    nonConformalCouple_12_on_NCC1_2
    {
        type            nonConformalCyclic;
        inGroups        List<word> 1(nonConformalCyclic);
        nFaces          0;
        startFace       2650615;
        matchTolerance  0.001;
        neighbourPatch  nonConformalCouple_12_on_NCC2_1;
        transformType   none;
        originalPatch   NCC1_2;
    }

[tomf]

Hi,

Unfortunately I have not yet worked with NCC, so I cannot comment I am afraid.

What I have read is that it should be more accurate and conservative, but I am not sure about these kind of settings.

Tom

[unilord]

Originally Posted by tomf

Hi,

Unfortunately I have not yet worked with NCC, so I cannot comment I am afraid.

What I have read is that it should be more accurate and conservative, but I am not sure about these kind of settings.

Tom

Thanks, I understand. Then, maybe you could aid me in this. When I switch the numerics to linearUpwind for turbulence (k and epsilon), starting the simulation using the converged results provided from upwind, I obtain residuals like the ones shown in the attached picture.

My fvSchemes are this:

ddtSchemes
{
    default		steadyState;
}

gradSchemes
{   
    //default             leastSquares;
    default		cellMDLimited Gauss linear 0.5;
    grad(U)             cellMDLimited Gauss linear 0.5;
    //grad(p)             cellMDLimited Gauss linear 0.5;
    //grad(epsilon)       cellMDLimited Gauss linear 0.5;
    //grad(k)             cellMDLimited Gauss linear 0.5;
}

divSchemes
{
    default                         none;
    div(phi,k)                      Gauss linearUpwind grad(k);
    //div(phi,omega)                  Gauss linearUpwind limited 1;
    div(phi,omega)                  Gauss linearUpwind grad(omega);
    div(phi,epsilon)                Gauss linearUpwind grad(epsilon);
    div(phi,U)                      bounded Gauss linearUpwind grad(U);
    div((nuEff)*dev2(T(grad(U))))   Gauss linear;
}

laplacianSchemes
{
    default		Gauss linear limited 1;
}

interpolationSchemes
{
    default		linear;
}

snGradSchemes
{
    default		limited 1;
}

wallDist
{
    method Poisson;
}

My fvSolutions are these:

solvers
{
    p
    {
        solver                    GAMG;
        tolerance                 1e-08;
        relTol                    0.0;
        smoother                  GaussSeidel;
        nCellsInCoarsestLevel     1000;
        //nSweeps                   200;
        //minIter                   10;
    }

    U
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-08;
        //nSweeps         200;
        relTol          0.0;
        //type            coupled;
        //solver          PBiCGStab;
        //preconditioner  DILU;
        //tolerance       (1e-08 1e-08 1e-08);
        //tolerance       1e-14
        //relTol          (0 0 0);
        //relTol          0.0;
        //minIter         3;  
    }

    k
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-08;
        relTol          0.0;
        //nSweeps         50;
    }

    omega
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-08;
        //nSweeps         50;
        relTol          0.0;
    }
    
    epsilon
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-08;
        //nSweeps         50;
        relTol          0.0;
    }
    
    Phi
    {
        solver          GAMG;
        smoother        DIC;

        tolerance       1e-06;
        relTol          0.01;
    }
    
    yPsi
    {
        solver          GAMG;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 1000;
        agglomerator    faceAreaPair;
        mergeLevels     1;
        tolerance       1e-7;
        relTol          0;
    }

}

potentialFlow
{
    nNonOrthogonalCorrectors 3;
}

SIMPLE
{
    momentumPredictor         yes;
    nNonOrthogonalCorrectors  2;
    nCorrector                0;
    
    residualControl
    {
    	//p		1e-4;
    	//U		1e-5;
    	//k		1e-7;
    	//epsilon		1e-7;
    }
}

relaxationFactors
{
    fields
    {
        p               0.3;
    }
    equations
    {
        U               0.7;
        k               0.4;
        omega           0.4;
        epsilon		0.4;
    }
}

Thank you,
Pedro Gouveia

[<nico>]

Originally Posted by tomf

Hi Nico,



For this type of analysis I would not recommend SRF. MRF could work, but there may be issues at the transition from rotating to stationary parts of the mesh.



It could well be that you need more iterations before starting the ramping up and maybe it needs to be slower.



Numerical diffusion would typically reduce extremes of pressure or velocity, so I do not think that this could be the case.



Overset would also not be necessary. as the rotating and stagnant parts cannot overlap, otherwise the pump cannot work. Can you show which part of your mesh is the MRF zone?



Regards,

Tom

Hi Tom,

I tried a lot more iterations and changed the outlet pressure from 0 to 10 at about 1000 iterations. The solver was stable and I stopped at 3600 iterations. The only thing is that my residuals don’t go below approximately 1e-3.

I think I am good with my boundaries now, so thank you for your advices already!
The rest might be due to my general setup.

As described, this is a very strong simplification of my real process because first of all I wanted to simulate the basic principle of a pump.
I was familiar with MRF which is why I chose this idea. See attached picture for the actual position of the zone.

In reality the application is a rotor/stator homogenizer. The rotor is far more closed to the pump housing that a small shear gap is given. For me the difficulty is the stator tool because it is kind of merged into the rotor region that it is quite hard to establish a simple cylindrical MRF or even AMI zone.
Thus, I was thinking about overset meshes and put every tool on a separate mesh that I can apply an easy cylindrical zone again. Or is this not possible for steady state or in general?

If not, I guess I need to prepare a very complex rotation zone around the rotor region with a very good mesh resolution for the rotation zone right?

Best regards,
Nico
image001.jpg

[unilord]

Originally Posted by <nico>

Hi Tom,

I tried a lot more iterations and changed the outlet pressure from 0 to 10 at about 1000 iterations. The solver was stable and I stopped at 3600 iterations. The only thing is that my residuals don’t go below approximately 1e-3.

I think I am good with my boundaries now, so thank you for your advices already!
The rest might be due to my general setup.

As described, this is a very strong simplification of my real process because first of all I wanted to simulate the basic principle of a pump.
I was familiar with MRF which is why I chose this idea. See attached picture for the actual position of the zone.

In reality the application is a rotor/stator homogenizer. The rotor is far more closed to the pump housing that a small shear gap is given. For me the difficulty is the stator tool because it is kind of merged into the rotor region that it is quite hard to establish a simple cylindrical MRF or even AMI zone.
Thus, I was thinking about overset meshes and put every tool on a separate mesh that I can apply an easy cylindrical zone again. Or is this not possible for steady state or in general?

If not, I guess I need to prepare a very complex rotation zone around the rotor region with a very good mesh resolution for the rotation zone right?

Best regards,
Nico
Attachment 98267

Hey Nico,

I am with some doubts. The MRF process is basically, using a single mesh, and afterwards using the topoSet function of openfoam in order to indicate where is the MRF region that rotates? If that is so, you are putting the angular velocity inside all the cells of that MRF region right?

[<nico>]

Hi Pedro,

and what if I introduce the cellZone during snappyHexMesh before merging the meshes? Then my zone is already defined and I could use it afterwards for my MRF. Or is this not possible?

Sorry if this is a stupid question, but I am still a beginner and maybe there is a basic understanding issue regarding MRF on my end.

[tomf]

Hi Nico,

Great that you found a solution. Yes you can define cellZones in snappyHexMesh, but using topoSet may work, or commercial tools for meshing may allow it.

The interface between the static and rotating part may not be optimal. If there is some separation in the geometry (there typically is), there can be some problems to get the solution to converge further than a certain limit as it cannot find a fully physical steady state.

Using overset on a steady case could work. I think there is overSimpleFoam? But I have not tried it myself.

Pedro, I am not completely sure about your settings, but like above if the physics is not leading to a steady state then you cannot expect the solver to always get down. You may want to use relTol not equal to 0.0 however. But I cannot tell what settings you should use in your case. It depends on the actual case, mesh quality, etc. But there are some examples to be found on the forum.

Best regards,
Tom

[unilord]

Originally Posted by tomf

Hi Nico,

Great that you found a solution. Yes you can define cellZones in snappyHexMesh, but using topoSet may work, or commercial tools for meshing may allow it.

The interface between the static and rotating part may not be optimal. If there is some separation in the geometry (there typically is), there can be some problems to get the solution to converge further than a certain limit as it cannot find a fully physical steady state.

Using overset on a steady case could work. I think there is overSimpleFoam? But I have not tried it myself.

Pedro, I am not completely sure about your settings, but like above if the physics is not leading to a steady state then you cannot expect the solver to always get down. You may want to use relTol not equal to 0.0 however. But I cannot tell what settings you should use in your case. It depends on the actual case, mesh quality, etc. But there are some examples to be found on the forum.

Best regards,
Tom

Hey Tom,

I didn't know that we should not set relTol equal to 0. I thought that would give the most accurate results for every case. Thank you for the intel.

I found one thing interesting when testing two differrent BC's, and I am wondering why that happens:

• - Fixed Pressure Outlet, Flow Rate Inlet: I get a velocity at inlet of 1.4 m/s at most of the inlet tube length, which is the analytical result for a flow rate of 0.0063 m^3/s and a pipe radius of 38mm
• - Fixed Pressure Inlet, Flow Rate Outlet: I get a velocity of 2m/s. Nevertheless, when comparing the pump head with experimental results provided in a scientific paper from 1989, this condition gives the smallest deviation (around 7.5% against 9%)

Do you have any idea why this happens?

Thank you,
Pedro

[Alczem]

Originally Posted by unilord

I am using NCC in a steady state because I am meshing with salome and adding layers with snappy

Hey!


Hijacking the thread, sorry, but I wanted to ask about adding layers on an already meshed geometry with snappy. Did you still need to feed the STLs to snappy, or is it able to detect the patches of the initial mesh (created by Salome) and create the layers without any additional files? If you could post your sHM dict, or any info, that would be great, thanks!

[tomf]

Originally Posted by unilord

I found one thing interesting when testing two differrent BC's, and I am wondering why that happens:

• - Fixed Pressure Outlet, Flow Rate Inlet: I get a velocity at inlet of 1.4 m/s at most of the inlet tube length, which is the analytical result for a flow rate of 0.0063 m^3/s and a pipe radius of 38mm
• - Fixed Pressure Inlet, Flow Rate Outlet: I get a velocity of 2m/s. Nevertheless, when comparing the pump head with experimental results provided in a scientific paper from 1989, this condition gives the smallest deviation (around 7.5% against 9%)

Do you have any idea why this happens?

Thank you,
Pedro

Hi Pedro,

In both cases is the velocity profile a uniform profile? I would expect it to be like that in the first case, but not in the second. Similarly there are probably going to be differences in the pressure profiles between both cases.

Regards,
Tom

[unilord]

Originally Posted by tomf

Hi Pedro,

In both cases is the velocity profile a uniform profile? I would expect it to be like that in the first case, but not in the second. Similarly there are probably going to be differences in the pressure profiles between both cases.

Regards,
Tom

Hey Tom,

They are uniform in both cases. However, the boundary layer is different. However, I must correct myself. The flow rates are the same, now that I saw the actual flow rate and not the freestream velocity. It is just that the maximum velocity at the tube is different. Please see attached pictures. Compare "InletFlow" with "OutletFlow". I am also attaching a third picture (OutletFlow2), corresponding to wrong (not actually wrong, just not following literature initial condition formulae) initialization variables for k and epsilon. The flow rate is also 0.0063 m^3/s, but the velocity profile is quite different.

Do you know why the change in BC's might affect the velocity profile?

One other question, I am now trying this setup of boundary conditions, since I read in some literature that is gives the most accurate results, despite needing more time steps for convergence.

Velocity:

boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
    "NCC.*"
    {
        type            slip;
    }
    outlet
    {
        type            pressureInletOutletVelocity;
        value           $internalField;
    }
    "wall.*"
    {
        type            noSlip;
    }
    "blade.*"
    {
        type            MRFnoSlip;
        value           uniform (0 0 0);
    }
    inlet
    {
        type            pressureInletOutletVelocity;
        value           uniform (0 0 1.4);
}

Pressure:

boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
    "NCC.*"
    {
        type            zeroGradient;
    }
    outlet
    {
        type            pressure;
        p               uniform 101.325;
    }
    "wall.*"
    {
        type            zeroGradient;
    }
    "blade.*"
    {
        type            zeroGradient;
    }
    inlet
    {
        type            totalPressure;
        p0              uniform 100.0768;
    }
}

However, I do now now what values I should input in totalPressure and pressure (corresponding to static pressure). The ones that I have now, give a flow rate of 0.00938 m^3/s in inlet, which is higher than the desired. I read in another thread that we should adjust the static pressure value in order to obtain the desired flow rate at inlet. Is that actually how we should proceed? For the totalPressure at inlet, following openfoam's formula for incompressible flow:

p_p = p_0 + 0.5 |U_0|^2 - 0.5 |U|^2

where

p_p - pressure at patch
p_0 - external static pressure (i am assuming 101.325 = ambient pressure)
U - velocity (I am using 1.4 m/s)
U_0 - external velocity (I am using 0 m/s)

What do you think about this?

Thank you in advance,
Pedro

[unilord]

Originally Posted by Alczem

Hey!


Hijacking the thread, sorry, but I wanted to ask about adding layers on an already meshed geometry with snappy. Did you still need to feed the STLs to snappy, or is it able to detect the patches of the initial mesh (created by Salome) and create the layers without any additional files? If you could post your sHM dict, or any info, that would be great, thanks!

Hey Alczem,

No problem at all, this is why we are here.

I actually found this snappy only layer thing quite usefull, since it is a powerful tool when salome can't mesh layers correctly (which i find that happens quite a lot of time).

You actually do not need the STL files, snappy can detect the patches, as long as you defined them in SALOME. I assume you know how to do that, but I am explaining here just in case:

- Create geometry, or import from CAD format
- Create groups of faces. You can do this in shaper, but i prefer to do in geometry module.
- In mesh module, just create a new mesh and make sure that the groups of faces also get meshed correctly.
- Export the mesh as .UNV
- Run "ideasUnvToFoam <Name_of_File.unv>" and you should see in terminal the names of the groups.
- SnappyHexMesh Layer should look something like this:

addLayersControls
{
    // Are the thickness parameters below relative to the undistorted
    // size of the refined cell outside layer (true) or absolute sizes (false).
    relativeSizes false;

    // Per final patch (so not geometry!) the layer information
    layers
    {
        wall_3
        {
            nSurfaceLayers      2;
            expansionRatio      1.2;
            firstLayerThickness 0.0012;
            minThickness        0.001;
        }
        
        wall_tongue
        {
            nSurfaceLayers      2;
            expansionRatio      1.1;
            firstLayerThickness 0.0007;
            minThickness        0.0005;
        }
    }

 //MORE LAYERS COMMANDS
}

where wall_3 and wall_tongue are names of groups that are created in SALOME. Below I am leaving the whole snappy directory that I am using. I found out that, having "geometry" and "regions" without anything, in the geometry, castelation and snap parts really help the meshing. Also, do not have the .eMesh surfaceFeatures geometry in the folder, since I also found that it gives some problems for the layer meshing.

/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  10
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
FoamFile
{
    format      ascii;
    class       dictionary;
    object      snappyHexMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

// Which of the steps to run
castellatedMesh false;
snap            false;
addLayers       true;

// Geometry. Definition of all surfaces. All surfaces are of class
// searchableSurface. 
// Surfaces are used
// - to specify refinement for any mesh cell intersecting it
// - to specify refinement for any mesh cell inside/outside/near
// - to 'snap' the mesh boundary to the surface
geometry
{
 
};

// Settings for the castellatedMesh generation.
castellatedMeshControls
{

    // Refinement parameters
    // ~~~~~~~~~~~~~~~~~~~~~

    // If local number of cells is >= maxLocalCells on any processor
    // switches from from refinement followed by balancing
    // (current method) to (weighted) balancing before refinement.
    maxLocalCells 10000000;

    // Overall cell limit (approximately). Refinement will stop immediately
    // upon reaching this number so a refinement level might not complete.
    // Note that this is the number of cells before removing the part which
    // is not 'visible' from the keepPoint. The final number of cells might
    // actually be a lot less.
    maxGlobalCells 20000000;

    // The surface refinement loop might spend lots of iterations refining just a
    // few cells. This setting will cause refinement to stop if <= minimumRefine
    // are selected for refinement. Note: it will at least do one iteration
    // (unless the number of cells to refine is 0)
    minRefinementCells 0;

    // Number of buffer layers between different levels.
    // 1 means normal 2:1 refinement restriction, larger means slower
    // refinement.
    nCellsBetweenLevels 4;



    // Explicit feature edge refinement
    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

    // Specifies a level for any cell intersected by its edges.
    // This is a featureEdgeMesh, read from constant/geometry for now.
    
    features
    (

    );
    


    // Surface based refinement
    // ~~~~~~~~~~~~~~~~~~~~~~~~

    // Specifies two levels for every surface. The first is the minimum level,
    // every cell intersecting a surface gets refined up to the minimum level.
    // The second level is the maximum level. Cells that 'see' multiple
    // intersections where the intersections make an
    // angle > resolveFeatureAngle get refined up to the maximum level.

    refinementSurfaces
    {

    }

    resolveFeatureAngle 10;


    // Region-wise refinement
    // ~~~~~~~~~~~~~~~~~~~~~~

    // Specifies refinement level for cells in relation to a surface. One of
    // three modes
    // - distance. 'levels' specifies per distance to the surface the
    //   wanted refinement level. The distances need to be specified in
    //   descending order.
    // - inside. 'levels' is only one entry and only the level is used. All
    //   cells inside the surface get refined up to the level. The surface
    //   needs to be closed for this to be possible.
    // - outside. Same but cells outside.

    refinementRegions
    {

    }
    
        

    // Mesh selection
    // ~~~~~~~~~~~~~~

    // After refinement patches get added for all refinementSurfaces and
    // all cells intersecting the surfaces get put into these patches. The
    // section reachable from the insidePoint is kept.
    // NOTE: This point should never be on a face, always inside a cell, even
    // after refinement.
    // This is an outside point insidePoint (-0.033 -0.033 0.0033);
    insidePoint (0 0.01 0.06); // Inside point

    // Whether any faceZones (as specified in the refinementSurfaces)
    // are only on the boundary of corresponding cellZones or also allow
    // free-standing zone faces. Not used if there are no faceZones.
    allowFreeStandingZoneFaces false;
}



// Settings for the snapping.
snapControls
{
    //- Number of patch smoothing iterations before finding correspondence
    //  to surface
    nSmoothPatch 1;

    //- Relative distance for points to be attracted by surface feature point
    //  or edge. True distance is this factor times local
    //  maximum edge length.
    tolerance 1.3;

    //- Number of mesh displacement relaxation iterations.
    nSolveIter 150;

    //- Maximum number of snapping relaxation iterations. Should stop
    //  before upon reaching a correct mesh.
    nRelaxIter 30;

    // Feature snapping

        //- Number of feature edge snapping iterations.
        //  Leave out altogether to disable.
        nFeatureSnapIter 30;

        //- Detect (geometric) features by sampling the surface
        implicitFeatureSnap true;

        //- Use castellatedMeshControls::features
        explicitFeatureSnap false;

        //- Detect features between multiple surfaces
        //  (only for explicitFeatureSnap, default = false)
        multiRegionFeatureSnap false;
}



// Settings for the layer addition.
addLayersControls
{
    // Are the thickness parameters below relative to the undistorted
    // size of the refined cell outside layer (true) or absolute sizes (false).
    relativeSizes false;

    // Per final patch (so not geometry!) the layer information
    layers
    {
        wall_3
        {
            nSurfaceLayers      2;
            expansionRatio      1.2;
            firstLayerThickness 0.0012;
            minThickness        0.001;
        }
        
        wall_tongue
        {
            nSurfaceLayers      2;
            expansionRatio      1.1;
            firstLayerThickness 0.0007;
            minThickness        0.0005;
        }
    }

    // Expansion factor for layer mesh
    expansionRatio 1.3;


    // Wanted thickness of final added cell layer. If multiple layers
    // is the thickness of the layer furthest away from the wall.
    // Relative to undistorted size of cell outside layer.
    // See relativeSizes parameter.
    firstLayerThickness 0.0002;

    // Minimum thickness of cell layer. If for any reason layer
    // cannot be above minThickness do not add layer.
    // See relativeSizes parameter.
    minThickness 0.0001;

    // If points get not extruded do nGrow layers of connected faces that are
    // also not grown. This helps convergence of the layer addition process
    // close to features.
    // Note: changed(corrected) w.r.t 17x! (didn't do anything in 17x)
    nGrow 0;

    // Advanced settings


    // Static analysis of starting mesh

        // When not to extrude surface. 0 is flat surface, 90 is when two faces
        // are perpendicular
        featureAngle 70;

        // Stop layer growth on highly warped cells
        maxFaceThicknessRatio 1;


    // Patch displacement

        // Number of smoothing iterations of surface normals
        nSmoothSurfaceNormals 3;

        // Smooth layer thickness over surface patches
        nSmoothThickness 1;



    // Medial axis analysis

        // Angle used to pick up medial axis points
        // Note: changed(corrected) w.r.t 17x! 90 degrees corresponds to 130
        // in 17x.
        minMedialAxisAngle 90;

        // Reduce layer growth where ratio thickness to medial
        // distance is large
        maxThicknessToMedialRatio 1;

        // Number of smoothing iterations of interior mesh movement direction
        nSmoothNormals 6;

        // Optional: limit the number of steps walking away from the surface.
        // Default is unlimited.
        // nMedialAxisIter 10;

        // Optional: smooth displacement after medial axis determination.
        // default is 0.
        // nSmoothDisplacement 90;

        // (wip)Optional: do not extrude a point if none of the surrounding points is
        // not extruded. Default is false.
        // detectExtrusionIsland true;


    // Mesh shrinking

        // Optional: at non-patched sides allow mesh to slip if extrusion
        // direction makes angle larger than slipFeatureAngle. Default is
        // 0.5*featureAngle.
        slipFeatureAngle 35;

        // Maximum number of snapping relaxation iterations. Should stop
        // before upon reaching a correct mesh.
        nRelaxIter 50;

        // Create buffer region for new layer terminations
        nBufferCellsNoExtrude 0;

        // Overall max number of layer addition iterations. The mesher will
        // exit if it reaches this number of iterations; possibly with an
        // illegal mesh.
        nLayerIter 2000;

        // Max number of iterations after which relaxed meshQuality controls
        // get used. Up to nRelaxedIter it uses the settings in
        // meshQualityControls,
        // after nRelaxedIter it uses the values in
        // meshQualityControls::relaxed.
        nRelaxedIter 50;

        // Additional reporting: if there are just a few faces where there
        // are mesh errors (after adding the layers) print their face centres.
        // This helps in tracking down problematic mesh areas.
        // additionalReporting true;
}



// Generic mesh quality settings. At any undoable phase these determine
// where to undo.
meshQualityControls
{
    #include "meshQualityDict"

    // Optional : some meshing phases allow usage of relaxed rules.
    // See e.g. addLayersControls::nRelaxedIter.
    relaxed
    {
        //- Maximum non-orthogonality allowed. Set to 180 to disable.
        maxNonOrtho 70;
    }
}


// Advanced

// Write flags
writeFlags
(
    scalarLevels    // write volScalarField with cellLevel for postprocessing
    layerSets       // write cellSets, faceSets of faces in layer
    layerFields     // write volScalarField for layer coverage
);


// Merge tolerance. Is fraction of overall bounding box of initial mesh.
// Note: the write tolerance needs to be higher than this.
mergeTolerance 1e-6;


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

Good luck,
Pedro

[Alczem]

Originally Posted by unilord

I actually found this snappy only layer thing quite usefull, since it is a powerful tool when salome can't mesh layers correctly (which i find that happens quite a lot of time).

You actually do not need the STL files, snappy can detect the patches, as long as you defined them in SALOME. I assume you know how to do that, but I am explaining here just in case:

- Create geometry, or import from CAD format
- Create groups of faces. You can do this in shaper, but i prefer to do in geometry module.
- In mesh module, just create a new mesh and make sure that the groups of faces also get meshed correctly.
- Export the mesh as .UNV
- Run "ideasUnvToFoam <Name_of_File.unv>" and you should see in terminal the names of the groups.
- SnappyHexMesh Layer should look something like this:

Hey,


Thanks a lot for the detailed answer, I agree with you that Salome can be disappointing when generating layers (my biggest gripe is that you cannot set a relative layer size in Salome). I will give your settings a try, it could definitely by a game changer for mesh generation in my case, thanks again

[unilord]

Originally Posted by Alczem

Hey,


Thanks a lot for the detailed answer, I agree with you that Salome can be disappointing when generating layers (my biggest gripe is that you cannot set a relative layer size in Salome). I will give your settings a try, it could definitely by a game changer for mesh generation in my case, thanks again

No problem, I found out that the less the iterations you put in the layer generation (nSmoothNormals, etc) the better. There is a sweet spot, and that is the difficult part. Another thing, I always try to keep "maxNonOrtho" in meshQualityDict higher than the nonOrtho that came from salome original mesh. nLayerIter is also problematic, if below 200 depending on the mesh. That is why I keep it at 2000. it really speeds up the process, even though having more iterations. I guess it's because it "converges" faster?

Anyway, if you have any question do not hesitate to ask here.

Pedro