[giulio.topazio]

Hi foamers!
I can not extrude the layers in areas where the radius of curvature becomes small! How can I fix this?
Moreover, even in some areas to the radius of curvature wider have the same problem ...
Any suggestion is really appreciated!

[anothr_acc]

Hi Giulio,

Did you fix this? I'm seeing something similar at the ends (90 degree) of a three-dimensional circular cylinder I'm trying to mesh.

Best regards,

Mark.

[giulio.topazio]

I'm sorry Mark but i gave up on using snappy...moreover than the issues with the curvatures i was in very big troubles with importing a well-done .stl geometry so i keep on using ANSA

[Rider]

I have the same problem ...

Originally Posted by giulio.topazio

I'm sorry Mark but i gave up on using snappy...moreover than the issues with the curvatures i was in very big troubles with importing a well-done .stl geometry so i keep on using ANSA

What do you think about ANSA ?

Have you a comparison of the same case ? (between snappyHexMesh and ANSA)

[giulio.topazio]

I have no comparison!but i can tell you that sweated a lot on snappy for 0 results!!Maybe the problem it's me!!
Good luck and good work!

[anothr_acc]

Some good news: I spent more time with snappy and got some results but I had to change to large radius of curvatures. I managed to get a circular cylinder with domed end to run well.

I'm back to sharp edges now and seeing the old problems. I'll persist....

[Rider]

Thanks for the news.

If it is possible, can you post here your new parameters and your conclusions ?

[anothr_acc]

Originally Posted by Rider

Thanks for the news.

If it is possible, can you post here your new parameters and your conclusions ?

Hi Rider, sorry for the slow response. I've split snappyHexMeshDict into
sections below, to ease my pain a bit. I've started performing castellation, snapping and layers in separate runs.

system/snappyHexMesh (header)

FoamFile { version 2.0; format ascii; class dictionary; object autoHexMeshDict; }

// Here, I made changes to bring my file in line with one on the net
// of an airfoil with boundary layer that ran around a sharp detail ok.

// Which of the steps to run {{{
castellatedMesh true; // whether to split cells or not.  Yes!
snap            true; // surface snapping; essential for cylinder.
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

// Specify the import stl file and the user-defined patch name.
// It is possible to use, `regions { one { name sphere; } two { name ellipse; } }'
// if the .stl has regions, `sphere' and `ellipse' defined.  Then, different refinements
// can be applied.  See
// http://www.cfd-online.com/Forums/openfoam-meshing-utilities/61483-snappyhexmesh-local-refinement-one-stlfile.html
// Generate the multi-region stl file by appending multiple exports.

geometry { cc2.stl  { type triSurfaceMesh; name cylinder;  } } ;

    //- Optionally, refine a bit extra around the small centre hole
    // refineHole  { type searchableSphere; centre (0 0 -0.012); radius 0.003; }
    //    };
    // Refine a bit around the cylinder?
    // refineCylinder { type searchableCylinder ; point1 ( 0 0 0 ) ; point2 ( 0 0 0.1 ) ; radius 0.025 ; }
    // no effect.

//}}}

system/snappyHexMesh (castellation section)

// Settings for the castellatedMesh generation. This is the cell splitting. {{{
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 100000;
    // While refining maximum number of cells per processor. This is basically
    // the number of cells that fit on a processor. If you choose this too small
    // it will do just more refinement iterations to obtain a similar mesh.
    maxLocalCells 100000000;


    // 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 200000000;

    // 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.  Keep this high.  From 1 to 5 resulted in factor 700
    // min determinant value and non-orth changing from 36 to 24.
    nCellsBetweenLevels 3;

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

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

    // create cylinder.eMesh
    // Run this from the shell:
    // surfaceFeatureExtract -includedAngle 150 -writeObj constant/triSurface/surface.stl cylinder
    features ( {  file "cc2.eMesh";  level 10;  } );
    // the level seems to make no difference.


    // 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.


    // this makes a big difference to the mesh.  This mesh has no
    // doubly intersecting surfaces so refine all cells by factor 8
    // on the cylinder.  The first number then, is 8.  The second
    // number the same (as the cylinder intersects the boundaries but
    // I don't want additional refinement in that case.

    refinementSurfaces  {  cylinder   {  level (6 6); } }

    // Or, multiple regions in the .stl:
    // refinementSurfaces  {  cylinder   {  level (2 8); regions { one { level (3 3); two { level (2 2); } } } } }
  // 2 2 gave ok mesh but symmetry plane cells pull away from
    // planes.  60 odd steps.
    // 2 4 gave refined mesh on the symmetry plane of
    // the cylinder and smaller cells on the cylinder. 74 steps.
    // 5 or 6 cells across the radius.
    // 4 8 wuth 1 degree gave a mesh that started to look good but
    // took a long time (>20 minutes).  Try 2 6.
    // 2 6 used 1 degree resolveFeatureAngle and this time I made sure the fine
    // gmsh mesh was used.  Previously, large surface elements were used and the
    // snappy mesh converged towards it -- but took too long.  This time at time 22,
    // step 18 the mesh has not been snapped to the surface but also has not peeled
    // away from the symmetry planes.  Nice looking and very regular along the
    // cylinder wall length.  That's with 4 refinements on the cylinder wall.  The
    // refinement level splits the cells once, so four refinements give a cell size
    // 0.5 * 0.5 * 0.5 * 0.5 = 1/16 the original.
    // Wow.  this mesh is stunning but the cells at the symmetry planes have pulled away
    // slightly.  They are small so the error is smaller but still present.
    // Reducing the gmsh mesh size has really helped.
    // Try 150 angles.  Used at command line with surfaceFeatureExtract:  Almost perfect!
    // Change to 30 in this file:  resolveFeatureAngle was 1....
    // This caused excellent behaviour at the interfaces between walls but the cylinder
    // wall, away from interactions with boundaries, is a coarse mesh.  I need to force the
    // wall mesh to be fine.  Then, this is complete.
    // Next:  can I force the cell size on `wall'?
    // or, reduce those angles progressively until the
    // cylinder wall refines but the right angles at interfaces
    // to symmetry planes are not lost




    // Surface-wise min and max refinement level



    //resolveFeatureAngle 30;  // original.
    resolveFeatureAngle 30; // 30 works.  60 stopped the boundary layer being created
    // resolveFeatureAngle 30 seemed ok with 2 4 with 74 steps.  Try 1.
    // 1 gave better repeated behaviour along the cylinder length.
    // 1 is required to refine all over the cylinder but will result
    // in flat surfaces not being refined.  This is the wrong approach.

    // 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.

    // These two work equally well.
    refinementRegions { cylinder { mode inside; levels ((1E15 3)); } }
    // refinementRegions {cylinder { mode distance; levels ( (0.25e-3 2) (0.5e-2 1) (0.01 0)); } }

    // 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 locationInMesh is kept.
    // NOTE: This point should never be on a face, always inside a cell, even
    // after refinement.
    // This is an outside point locationInMesh (-0.033 -0.033 0.0033);
    locationInMesh (0.05 0.05 0.05); // point inside the meshed region

    // 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 true;
} // }}}

system/snappyHexMesh (snapping section)

snapControls { // Settings for snapping {{{
    //- Number of patch smoothing iterations before finding correspondence
    //  to surface
    nSmoothPatch 3; // originally 3. Tried 6; no significant mesh quality difference.

    //- 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.0;  // Originally 1. 2 made no clear change to the mesh.

    //- Number of mesh displacement relaxation iterations.
    nSolveIter 30;
  //- Maximum number of snapping relaxation iterations. Should stop
    //  before upon reaching a correct mesh.
    nRelaxIter 3;  // originally 5. 10 made no difference to checkMesh.

    //- Highly experimental and wip: number of feature edge snapping
    //  iterations. Leave out altogether to disable.
    nFeatureSnapIter 10;
    // nFeatureSnapIter controls the number of iterations of snapping to the
    // feature object.  This, for example should help prevent edges peeling
    // away from the 2D mesh. Changing to 100 made a large difference to the
    // mesh time and little difference to the mesh quality.
}   // }}}

system/snappyHexMesh (layers section)

// 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 true;

    // Per final patch (so not geometry!) the layer information
    // The bike demo used a filename.  The link far above re: different
    // mesh regions referred to the user defined name above.
    layers
    {
        //"surface_.*"
        cylinder_CreatedbyGmsh { nSurfaceLayers 4; }  // changing from 1 to 2 did nothing.
        // If there are multiple regions, also list them here.

    }

    // Expansion factor for layer mesh.  Originally outside `layers'
    expansionRatio 1.3; // above 1.5 and no layer refinement occurs.  1.5 looks good.

    //- Wanted thickness of final added cell layer. If multiple layers
    //  is the thickness of the layer furthest away from the wall.
    //  See relativeSizes parameter.
    finalLayerThickness 0.5;  // this will set the size of the jump
    // between cells in the boundary layer and in the main region.
    // 0.9 sounds good.  I was using 0.3.  0.3 showed better nonOrthogonality but poorer determinants,
    // probably due to the smaller smallest cell size possible.

    //- Minimum thickness of cell layer. If for any reason layer
    //  cannot be above minThickness do not add layer.
    //  See relativeSizes parameter.
    minThickness 0.00001;  // little difference made to the mesh

   //- 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.
    nGrow 1; // was 1.  100 damaged the cylinder wall
   // Advanced settings

    //- When not to extrude surface. 0 is flat surface, 90 is when two faces
    //  make straight angle.
    featureAngle 180; // was 30.  100 seemed to make no difference nor 1.

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

    // Number of smoothing iterations of surface normals
    nSmoothSurfaceNormals 1; // was 1.  This does make for a smoother mesh moving from the b/l into the main region.

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

    // Smooth layer thickness over surface patches
    nSmoothThickness 10;

   // Stop layer growth on highly warped cells
    maxFaceThicknessRatio 0.5; // was 0.5

    // Reduce layer growth where ratio thickness to medial
    // distance is large
    maxThicknessToMedialRatio 0.3; // was 0.3. 0.003 resulted in max aspect ratio increasing from 34 to 62.
    // non-orthog was unaffected but small determinant cells started occurring.  Best keep to 0.3.

    // Angle used to pick up medial axis points
    minMedianAxisAngle 130; // 90.  Change to 1 and 180 no clear difference.

    // Create buffer region for new layer terminations
    nBufferCellsNoExtrude 0; // 0  Changing to 10 removed the b/l layers.


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

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

system/snappyHexMesh (other settings)

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

    //- Maximum non-orthogonality allowed. Set to 180 to disable.
    maxNonOrtho 65;

    //- Max skewness allowed. Set to <0 to disable.
    maxBoundarySkewness 20;
    maxInternalSkewness 4;

   //- Max concaveness allowed. Is angle (in degrees) below which concavity
    //  is allowed. 0 is straight face, <0 would be convex face.
    //  Set to 180 to disable.
    maxConcave 80;

    //- Minimum pyramid volume. Is absolute volume of cell pyramid.
    //  Set to a sensible fraction of the smallest cell volume expected.
    //  Set to very negative number (e.g. -1E30) to disable.
    minVol 1e-13;

    //- Minimum quality of the tet formed by the face-centre
    //  and variable base point minimum decomposition triangles and
    //  the cell centre.  Set to very negative number (e.g. -1E30) to
    //  disable.
    //     <0 = inside out tet,
    //      0 = flat tet
    //      1 = regular tet
    minTetQuality 1e-30;

  //- Minimum face area. Set to <0 to disable.
    minArea -1;

    //- Minimum face twist. Set to <-1 to disable. dot product of face normal
    //- and face centre triangles normal
    minTwist 0.05;

  //- minimum normalised cell determinant
    //- 1 = hex, <= 0 = folded or flattened illegal cell
    minDeterminant 0.001;

    //- minFaceWeight (0 -> 0.5)
    minFaceWeight 0.05;

    //- minVolRatio (0 -> 1)
    minVolRatio 0.01;

    //must be >0 for Fluent compatibility
    minTriangleTwist -1;

  //- if >0 : preserve single cells with all points on the surface if the
    //  resulting volume after snapping (by approximation) is larger than
    //  minVolCollapseRatio times old volume (i.e. not collapsed to flat cell).
    //  If <0 : delete always.
    //minVolCollapseRatio 0.5;


    // Advanced

    //- Number of error distribution iterations
    nSmoothScale 4;
    //- amount to scale back displacement at error points
    errorReduction 0.75;
   // Optional : some meshing phases allow usage of relaxed rules.
    // See e.g. addLayersControls::nRelaxedIter.

    //- Maximum non-orthogonality allowed. Set to 180 to disable.
    relaxed { maxNonOrtho 75; }
    }  // }}}


// Advanced

// Flags for optional output
// 0 : only write final meshes
// 1 : write intermediate meshes
// 2 : write volScalarField with cellLevel for postprocessing
// 4 : write current intersections as .obj files
debug 2;


// Merge tolerance. Is fraction of overall bounding box of initial mesh.
// Note: the write tolerance needs to be higher than this.
mergeTolerance 1E-8;
// might this fix small dodgy cells?  It was 1E-6

I created a circular cylinder with gmsh, with quite a high resolution. The cylinder had curved ends; I'm still having troubles getting boundary layers to wrap around sharp corners.

This mesh solved slowly but well. I initially solved as laminar until converge was clear from the residuals and then turned on Reynolds stress turbulence. Of great importance to me was using wall functions (I omitted them originally and generated poor results) and setting zeroGradients or calculated values on input and output planes for almost all parameters except velocity at the input; instead of say, setting a constant epsilon. When I set zero Reynolds stresses at the input and constant epsilon, the solution blew. I ended up using I used symmetryPlanes on the side walls. Separation was a little late (about 90 degrees) but not too bad. I was pleased.

I've included my 0/* too as this was critical to getting convergence:

[anothr_acc]

FoamFile{version     2.0;  format      ascii;  class       volScalarField;
    location    "0"; object      epsilon;}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -3 0 0 0 0];

internalField   uniform 10;

boundaryField
{
    in { type turbulentMixingLengthDissipationRateInlet; mixingLength 0.04; value uniform 0.04;}
    // in { type         freestream;  freestreamValue uniform 0.94;  value uniform 0.94; }
    out {  type zeroGradient; } // freestream; freestreamValue uniform 0.05;  value  uniform 0.05;  }
    inside { type            symmetryPlane;  }
    outside { type           symmetryPlane;  }
    topbot  { type           symmetryPlane;  }
 // cylinder_CreatedbyGmsh   {type   zeroGradient;  }
    cylinder_CreatedbyGmsh   {type   epsilonWallFunction;  value uniform 0;}
}

FoamFile{version 2.0; format ascii; class volScalarField; object k;}

boundaryField {
   //in{type freestream; freestreamValue uniform 0.375;}
//   in{type turbulentIntensityKineticEnergyInlet; intensity 0.01; value uniform 0.01;}
   in{type turbulentIntensityKineticEnergyInlet; intensity 0.01; value calculated;}
   out{type zeroGradient;} // freestream; freestreamValue uniform 0.0015;}
   outside{type symmetryPlane;}
   inside{type symmetryPlane;}
   topbot{type symmetryPlane;}
   //cylinder_CreatedbyGmsh{type fixedValue; value uniform 1E-6;} //
   cylinder_CreatedbyGmsh{type kqRWallFunction; value uniform 0;} //
}
internalField uniform 0.375;
dimensions [0 2 -2 0 0 0 0];

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    location    "0";
    object      nut;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -1 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    in { type            calculated; value           uniform 0; }
    out {type            calculated; value           uniform 0; }
    inside { type            symmetryPlane; }
    outside { type            symmetryPlane; }
    topbot {  type            symmetryPlane; }
    //cylinder_CreatedbyGmsh { type            calculated;  value           uniform 0; }
    cylinder_CreatedbyGmsh { type    nutkWallFunction; value uniform 0;}
}

FoamFile{version 2.0; format ascii; class volScalarField; object p;}

boundaryField {   in{type freestreamPressure;}  //  zeroGradient
   out{type fixedValue; value uniform 0;} //zeroGradient;} //  value uniform 0;} // zeroGradient; freestreamPressure; fixedValue; value uniform 0
   inside{type symmetryPlane;}
   outside{type symmetryPlane;}
   topbot{type symmetryPlane;}
   cylinder_CreatedbyGmsh{type zeroGradient;}
}
internalField uniform 0;
dimensions [0 2 -2 0 0 0 0];

FoamFile{version 2.0; format ascii; class volSymmTensorField; object R;}

boundaryField{
      //cylinder_CreatedbyGmsh{  type fixedValue; value uniform (0 0 0 0 0 0);}
      cylinder_CreatedbyGmsh{  type kqRWallFunction; value uniform (0 0 0 0 0 0); }
      inside{  type symmetryPlane;}
      outside{  type symmetryPlane;}
      topbot{  type symmetryPlane;}
      in{  type freestream; freestreamValue uniform (0 0 0 0 0 0); } // fixedValue; value uniform (10 0 0);}
      out{  type zeroGradient; } // freestream; freestreamValue uniform (0 0 0 0 0 0);} // zeroGradient;}
   }
internalField uniform (0 0 0 0 0 0);
dimensions [0 2 -2 0 0 0 0];

FoamFile{version 2.0; format ascii; class volVectorField; object U;}

boundaryField{
      cylinder_CreatedbyGmsh{  type fixedValue; value uniform (0 0 0);}
      inside{  type symmetryPlane;}
      outside{  type symmetryPlane;}
      topbot{  type symmetryPlane;}
      in{  type freestream; freestreamValue uniform (10 0  0); } // fixedValue; value uniform (10 0 0);}
      out{  type zeroGradient; } // freestream; freestreamValue uniform (6.42 -7.66 0);} // zeroGradient;}
   }
internalField uniform (10 0  0);
dimensions [0 1 -1 0 0 0 0];

Conclusions: Laminar then turbulent solution was essential. Relaxing the boundary conditions at the input from constant values to calculated ones was essential.

As I say, sharp corners are still upsetting me. Work goes on....

Best regards,

Mark.

[Rider]

Thank you Mark for your very detailed answer!

Best regards,
Rider.

[Eloise]

Originally Posted by anothr_acc

As I say, sharp corners are still upsetting me. Work goes on....

Have you tried in the addLayersControls section of snappyHexMeshDict to use at the same time the two following parameters?
nGrow = 0;
featureAngle = 180;

If I had to explain those parameters with my words, I'd say:
nGrow is the distance from the sharp edge at which the layers will collapse (counted in number of cells).
featureAngle is the maximum edge angle will be surrounded by layers. In your case, 60 or 70 might be enough.

Eloïse