Skip to content
Snippets Groups Projects
Select Git revision
  • gmsh_4_8_4
  • master default protected
  • bl
  • pluginMeshQuality
  • fixBugsAmaury
  • hierarchical-basis
  • alphashapes
  • relaying
  • new_export_boris
  • oras_vs_osm
  • reassign_partitions
  • distributed_fwi
  • rename-classes
  • fix/fortran-api-example-t4
  • robust_partitions
  • reducing_files
  • fix_overlaps
  • 3115-issue-fix
  • 3023-Fillet2D-Update
  • convert_fdivs
  • tmp_jcjc24
  • gmsh_4_14_0
  • gmsh_4_13_1
  • gmsh_4_13_0
  • gmsh_4_12_2
  • gmsh_4_12_1
  • gmsh_4_12_0
  • gmsh_4_11_1
  • gmsh_4_11_0
  • gmsh_4_10_5
  • gmsh_4_10_4
  • gmsh_4_10_3
  • gmsh_4_10_2
  • gmsh_4_10_1
  • gmsh_4_10_0
  • gmsh_4_9_5
  • gmsh_4_9_4
  • gmsh_4_9_3
  • gmsh_4_9_2
  • gmsh_4_9_1
  • gmsh_4_9_0
41 results

t13.geo

Blame
  • t13.geo 2.66 KiB
    // -----------------------------------------------------------------------------
    //
    //  Gmsh GEO tutorial 13
    //
    //  Remeshing an STL file without an underlying CAD model
    //
    // -----------------------------------------------------------------------------
    
    // Let's merge an STL mesh that we would like to remesh.
    Merge "t13_data.stl";
    
    // We first classify ("color") the surfaces by splitting the original surface
    // along sharp geometrical features. This will create new discrete surfaces,
    // curves and points.
    
    DefineConstant[
      // Angle between two triangles above which an edge is considered as sharp
      angle = {40, Min 20, Max 120, Step 1,
        Name "Parameters/Angle for surface detection"},
      // For complex geometries, patches can be too complex, too elongated or too
      // large to be parametrized; setting the following option will force the
      // creation of patches that are amenable to reparametrization:
      forceParametrizablePatches = {0, Choices{0,1},
        Name "Parameters/Create surfaces guaranteed to be parametrizable"},
      // For open surfaces include the boundary edges in the classification process:
      includeBoundary = 1,
      // Force curves to be split on given angle:
      curveAngle = 180
    ];
    ClassifySurfaces{angle * Pi/180, includeBoundary, forceParametrizablePatches,
                     curveAngle * Pi / 180};
    
    // Create a geometry for all the discrete curves and surfaces in the mesh, by
    // computing a parametrization for each one
    CreateGeometry;
    
    // In batch mode the two steps above can be performed with `gmsh t13.stl
    // -reparam 40', which will save `t13.msh' containing the parametrizations, and
    // which can thus subsequently be remeshed.
    
    // Note that if a CAD model (e.g. as a STEP file, see `t20.geo') is available
    // instead of an STL mesh, it is usually better to use that CAD model instead of
    // the geometry created by reparametrizing the mesh. Indeed, CAD geometries will
    // in general be more accurate, with smoother parametrizations, and will lead to
    // more efficient and higher quality meshing. Discrete surface remeshing in Gmsh
    // is optimized to handle dense STL meshes coming from e.g. imaging systems
    // where no CAD is available; it is less well suited for the poor quality STL
    // triangulations (optimized for size, with e.g. very elongated triangles) that
    // are usually generated by CAD tools for e.g. 3D printing.
    
    // Create a volume as usual
    Surface Loop(1) = Surface{:};
    Volume(1) = {1};
    
    // We specify element sizes imposed by a size field, just because we can :-)
    funny = DefineNumber[0, Choices{0,1},
      Name "Parameters/Apply funny mesh size field?" ];
    
    Field[1] = MathEval;
    If(funny)
      Field[1].F = "2*Sin((x+y)/5) + 3";
    Else
      Field[1].F = "4";
    EndIf
    Background Field = 1;