0
# Helper Functions and Utilities
1

2
Essential utility functions for mesh processing, geometric transformations, file I/O, and mesh optimization that complement the core geometry operations. These functions provide crucial functionality for working with meshes, performing geometric calculations, and optimizing mesh quality.
3

4
## Capabilities
5

6
### File I/O Operations
7

8
Write geometry and mesh data to various file formats supported by Gmsh, enabling integration with external analysis tools and mesh processing workflows.
9

10
```python { .api }
11
def write(filename: str):
12
    """
13
    Write current Gmsh model to file.
14
    
15
    Supports multiple file formats based on extension:
16
    - .msh: Gmsh mesh format (versions 2.2, 4.0, 4.1)
17
    - .vtk: VTK format for ParaView/VisIt visualization
18
    - .vtu: VTK unstructured grid format
19
    - .med: MED format for Code_Aster/Salome
20
    - .stl: STL format for 3D printing
21
    - .ply: PLY format for 3D graphics
22
    - .obj: Wavefront OBJ format
23
    - .off: Object File Format
24
    
25
    Parameters:
26
    - filename: str, output filename with extension determining format
27
    
28
    Returns:
29
    None (writes file to disk)
30
    """
31
```
32

33
### Geometric Transformations
34

35
Perform essential geometric calculations including rotation matrices and line orientation for complex geometric modeling operations.
36

37
```python { .api }
38
def rotation_matrix(u, theta):
39
    """
40
    Create rotation matrix for rotation around vector u by angle theta.
41
    
42
    Based on Rodrigues' rotation formula:
43
    https://en.wikipedia.org/wiki/Rotation_matrix#Rotation_matrix_from_axis_and_angle
44
    
45
    Parameters:
46
    - u: array-like, rotation axis vector (must be unit vector)
47
    - theta: float, rotation angle in radians
48
    
49
    Returns:
50
    numpy.ndarray: 3x3 rotation matrix
51
    
52
    Raises:
53
    AssertionError: if u is not a unit vector (|u| != 1)
54
    """
55

56
def orient_lines(lines):
57
    """
58
    Reorder and reorient unordered lines to form a closed polygon.
59
    
60
    Takes a sequence of unordered and potentially mis-oriented line segments
61
    that define a closed polygon and returns them in proper order with 
62
    correct orientations to form a valid curve loop.
63
    
64
    Parameters:
65
    - lines: sequence of Line objects defining a closed polygon
66
    
67
    Returns:
68
    list: Reordered and reoriented Line objects forming closed loop
69
    
70
    Algorithm:
71
    1. Analyzes vertex connectivity of all line segments
72
    2. Reorders lines to form continuous path
73
    3. Reverses line orientation where needed for consistency
74
    4. Ensures resulting loop is closed and properly oriented
75
    """
76
```
77

78
### Mesh Optimization
79

80
Optimize mesh quality using various algorithms to improve element shapes, reduce distortion, and enhance numerical accuracy.
81

82
```python { .api }
83
def optimize(mesh, method="", verbose=False):
84
    """
85
    Optimize mesh quality using Gmsh optimization methods.
86
    
87
    Applies mesh optimization algorithms to improve element quality metrics
88
    such as aspect ratio, skewness, and Jacobian determinant. Optimization
89
    can significantly improve numerical accuracy in finite element analysis.
90
    
91
    Parameters:
92
    - mesh: meshio.Mesh, input mesh object to optimize
93
    - method: str, optimization method name (empty string uses default)
94
      Available methods include:
95
      - "": Default optimization (Gmsh automatic selection)
96
      - "Netgen": Netgen-based optimization
97
      - "HighOrder": High-order node positioning optimization
98
      - "HighOrderElastic": Elastic high-order optimization
99
      - "HighOrderFastCurving": Fast curving for high-order elements
100
      - "Laplace2D": Laplacian smoothing for 2D meshes
101
      - "Relocate2D": Node relocation for 2D meshes
102
    - verbose: bool, enable verbose output during optimization
103
    
104
    Returns:
105
    meshio.Mesh: Optimized mesh object with improved element quality
106
    
107
    Note:
108
    Optimization modifies node positions and may add/remove elements
109
    to improve overall mesh quality. Physical groups and mesh regions
110
    are preserved during the optimization process.
111
    """
112
```
113

114
### Mesh Extraction
115

116
Extract mesh data from Gmsh into meshio format for further processing, analysis, or export to other tools.
117

118
```python { .api }
119
def extract_to_meshio():
120
    """
121
    Extract current Gmsh mesh to meshio format.
122
    
123
    Converts the active Gmsh mesh model into a meshio.Mesh object,
124
    including all mesh entities (points, lines, surfaces, volumes),
125
    physical groups, and material properties.
126
    
127
    Returns:
128
    meshio.Mesh: Mesh object containing:
129
    - points: numpy.ndarray, node coordinates (N×3)
130
    - cells: list of CellBlock objects with connectivity data
131
    - cell_sets: dict, physical groups mapped to cell indices
132
    - point_data: dict, node-based data fields (if any)
133
    - cell_data: dict, element-based data fields (if any)
134
    
135
    Supported element types:
136
    - 1D: Line elements (line, line3)
137
    - 2D: Triangle, quadrilateral elements (triangle, triangle6, quad, quad9)
138
    - 3D: Tetrahedron, hexahedron, prism, pyramid elements
139
    
140
    Physical groups are converted to cell_sets with original names,
141
    enabling material assignment and boundary condition specification
142
    in external analysis codes.
143
    """
144
```
145

146
## Usage Examples
147

148
### Optimizing Mesh Quality
149

150
```python
151
import pygmsh
152

153
# Create geometry and generate initial mesh
154
with pygmsh.occ.Geometry() as geom:
155
    # Create complex geometry
156
    sphere = geom.add_ball([0, 0, 0], 1.0, mesh_size=0.2)
157
    box = geom.add_box([-0.5, -0.5, -0.5], [1, 1, 1], mesh_size=0.15)
158
    complex_shape = geom.boolean_union([sphere, box])
159
    
160
    # Generate initial mesh
161
    initial_mesh = geom.generate_mesh(dim=3)
162

163
# Optimize mesh quality
164
optimized_mesh = pygmsh.optimize(initial_mesh, method="HighOrder", verbose=True)
165

166
# Write optimized mesh to file
167
pygmsh.write("optimized_mesh.msh")
168

169
# Compare mesh statistics
170
print(f"Initial elements: {len(initial_mesh.cells[0].data)}")
171
print(f"Optimized elements: {len(optimized_mesh.cells[0].data)}")
172
```
173

174
### Creating Proper Line Loops from Unordered Lines
175

176
```python
177
import pygmsh
178
import numpy as np
179

180
with pygmsh.geo.Geometry() as geom:
181
    # Create unordered line segments for a polygon
182
    points = [
183
        geom.add_point([0, 0, 0]),
184
        geom.add_point([1, 0, 0]),
185
        geom.add_point([1, 1, 0]),
186
        geom.add_point([0.5, 1.5, 0]),
187
        geom.add_point([0, 1, 0])
188
    ]
189
    
190
    # Create lines in random order (simulating imported data)
191
    unordered_lines = [
192
        geom.add_line(points[2], points[3]),  # line 2-3
193
        geom.add_line(points[0], points[1]),  # line 0-1
194
        geom.add_line(points[4], points[0]),  # line 4-0
195
        geom.add_line(points[1], points[2]),  # line 1-2
196
        geom.add_line(points[3], points[4])   # line 3-4
197
    ]
198
    
199
    # Orient lines to form proper closed loop
200
    oriented_lines = pygmsh.orient_lines(unordered_lines)
201
    
202
    # Create surface from properly oriented lines
203
    curve_loop = geom.add_curve_loop(oriented_lines)
204
    surface = geom.add_plane_surface(curve_loop)
205
    
206
    mesh = geom.generate_mesh(dim=2)
207
```
208

209
### Using Rotation Matrices for Geometric Transformations
210

211
```python
212
import pygmsh
213
import numpy as np
214

215
with pygmsh.occ.Geometry() as geom:
216
    # Create base geometry
217
    cylinder = geom.add_cylinder([0, 0, 0], [0, 0, 1], 0.5)
218
    
219
    # Create rotation matrix for 45-degree rotation around [1,1,1] axis
220
    axis = np.array([1, 1, 1]) / np.sqrt(3)  # Normalize to unit vector
221
    angle = np.pi / 4  # 45 degrees
222
    R = pygmsh.rotation_matrix(axis, angle)
223
    
224
    # Apply rotation to create multiple rotated copies
225
    rotated_cylinders = [cylinder]
226
    for i in range(1, 8):
227
        # Accumulate rotations
228
        total_angle = i * angle
229
        R_total = pygmsh.rotation_matrix(axis, total_angle)
230
        
231
        # Create rotated copy (conceptual - actual rotation done via Gmsh)
232
        rotated_cyl = geom.copy(cylinder)
233
        geom.rotate(rotated_cyl, [0, 0, 0], total_angle, axis)
234
        rotated_cylinders.append(rotated_cyl)
235
    
236
    # Union all cylinders
237
    final_shape = geom.boolean_union(rotated_cylinders)
238
    
239
    mesh = geom.generate_mesh(dim=3)
240
```
241

242
### Mesh Format Conversion and Analysis
243

244
```python
245
import pygmsh
246
import meshio
247

248
# Create and optimize mesh
249
with pygmsh.geo.Geometry() as geom:
250
    circle = geom.add_circle([0, 0, 0], 1.0, mesh_size=0.1)
251
    mesh = geom.generate_mesh(dim=2)
252

253
# Optimize mesh quality
254
optimized_mesh = pygmsh.optimize(mesh, verbose=True)
255

256
# Write to multiple formats for different tools
257
pygmsh.write("mesh.msh")        # Gmsh format
258
pygmsh.write("mesh.vtk")        # VTK for ParaView  
259
pygmsh.write("mesh.med")        # MED for Code_Aster
260
pygmsh.write("mesh.stl")        # STL for 3D printing
261

262
# Extract to meshio for analysis
263
mesh_data = pygmsh.extract_to_meshio()
264

265
# Analyze mesh statistics
266
print(f"Number of nodes: {len(mesh_data.points)}")
267
print(f"Number of elements: {sum(len(cell.data) for cell in mesh_data.cells)}")
268
print(f"Physical groups: {list(mesh_data.cell_sets.keys())}")
269

270
# Compute mesh quality metrics
271
element_areas = []
272
for cell_block in mesh_data.cells:
273
    if cell_block.type == "triangle":
274
        # Compute triangle areas
275
        triangles = cell_block.data
276
        points = mesh_data.points
277
        for tri in triangles:
278
            p0, p1, p2 = points[tri]
279
            area = 0.5 * np.linalg.norm(np.cross(p1 - p0, p2 - p0))
280
            element_areas.append(area)
281

282
print(f"Average element area: {np.mean(element_areas):.6f}")
283
print(f"Element area std dev: {np.std(element_areas):.6f}")
284
```
285

286
### Advanced Mesh Processing Pipeline
287

288
```python
289
import pygmsh
290
import numpy as np
291

292
def create_and_process_mesh(geometry_func, optimization_method="HighOrder"):
293
    """
294
    Complete mesh processing pipeline with geometry creation, 
295
    mesh generation, optimization, and quality analysis.
296
    """
297
    
298
    # Create geometry
299
    with pygmsh.occ.Geometry() as geom:
300
        # Apply user-defined geometry function
301
        geometry_entities = geometry_func(geom)
302
        
303
        # Generate initial mesh
304
        initial_mesh = geom.generate_mesh(dim=3)
305
    
306
    # Optimize mesh
307
    optimized_mesh = pygmsh.optimize(
308
        initial_mesh, 
309
        method=optimization_method, 
310
        verbose=True
311
    )
312
    
313
    # Extract mesh data for analysis
314
    mesh_data = pygmsh.extract_to_meshio()
315
    
316
    # Write mesh in multiple formats
317
    formats = ["mesh.msh", "mesh.vtk", "mesh.med"]
318
    for fmt in formats:
319
        pygmsh.write(fmt)
320
    
321
    return {
322
        'initial_mesh': initial_mesh,
323
        'optimized_mesh': optimized_mesh,
324
        'mesh_data': mesh_data,
325
        'num_nodes': len(mesh_data.points),
326
        'num_elements': sum(len(cell.data) for cell in mesh_data.cells)
327
    }
328

329
# Example usage
330
def complex_geometry(geom):
331
    """Create complex geometry for testing."""
332
    main_body = geom.add_ball([0, 0, 0], 2.0)
333
    
334
    holes = [
335
        geom.add_cylinder([x, y, -3], [0, 0, 1], 0.3)
336
        for x, y in [(-1, -1), (1, -1), (1, 1), (-1, 1)]
337
    ]
338
    
339
    return geom.boolean_difference(main_body, holes)
340

341
# Process mesh with pipeline
342
results = create_and_process_mesh(complex_geometry, "HighOrderElastic")
343
print(f"Generated mesh with {results['num_elements']} elements")
344
```
345

346
## Command Line Interface
347

348
### Mesh Optimization CLI
349

350
Command-line interface for mesh optimization operations, enabling batch processing and integration with external workflows.
351

352
```python { .api }
353
def optimize_cli(argv=None):
354
    """
355
    Command-line interface for mesh optimization.
356
    
357
    Provides a command-line tool for optimizing mesh files using
358
    various Gmsh optimization algorithms. Supports all meshio-compatible
359
    formats for input and output.
360
    
361
    Command line usage:
362
    pygmsh-optimize input.vtk output.xdmf [--method METHOD] [--verbose]
363
    
364
    Parameters:
365
    - argv: list[str], optional command line arguments (uses sys.argv if None)
366
    
367
    Supported formats:
368
    - Input: .msh, .vtk, .vtu, .stl, .ply, .obj, .off, .med
369
    - Output: .msh, .xdmf, .vtk, .vtu, .stl, .ply, .obj, .off, .med
370
    
371
    Returns:
372
    int: Exit code (0 for success, non-zero for errors)
373
    
374
    Examples:
375
    # Optimize VTK mesh and save as XDMF
376
    pygmsh-optimize mesh.vtk optimized.xdmf
377
    
378
    # Use specific optimization method with verbose output
379
    pygmsh-optimize input.msh output.msh --method "HighOrderOptimize" --verbose
380
    """
381
```