Supports: How to Optimize For Better Prints

import trimesh
import numpy as np
import logging
from typing import List, Tuple
import sys

# Configure logging for debug output
logging.basicConfig(
    level=logging.INFO,
    format='%(asctime)s - %(levelname)s - %(message)s'
)
logger = logging.getLogger(__name__)

class OverhangAnalyzer:
    """Analyzes STL mesh geometry to identify regions requiring support structures."""

    def __init__(self, overhang_angle_threshold: float = 45.0, min_overhang_area: float = 10.0):
        """
        Initialize analyzer with configurable thresholds.

        Args:
            overhang_angle_threshold: Maximum angle (degrees) from vertical where supports are needed (default 45°)
            min_overhang_area: Minimum area (mm²) of overhang region to trigger support generation (default 10mm²)
        """
        self.overhang_angle_threshold = np.radians(overhang_angle_threshold)
        self.min_overhang_area = min_overhang_area
        self.mesh = None
        self.overhang_faces = None

    def load_stl(self, stl_path: str) -> None:
        """
        Load and validate STL file from disk.

        Args:
            stl_path: Absolute or relative path to STL file

        Raises:
            FileNotFoundError: If STL file does not exist
            ValueError: If loaded mesh has invalid geometry (e.g., non-manifold edges)
        """
        try:
            self.mesh = trimesh.load(stl_path)
            if not isinstance(self.mesh, trimesh.Trimesh):
                raise ValueError(f"Loaded file {stl_path} is not a valid closed mesh (type: {type(self.mesh)})")
            if not self.mesh.is_watertight:
                logger.warning(f"Mesh {stl_path} is not watertight — support calculations may be inaccurate")
            logger.info(f"Loaded STL: {stl_path} | Vertices: {len(self.mesh.vertices)} | Faces: {len(self.mesh.faces)}")
        except FileNotFoundError:
            logger.error(f"STL file not found: {stl_path}")
            raise
        except Exception as e:
            logger.error(f"Failed to load STL {stl_path}: {str(e)}")
            raise

    def calculate_overhangs(self) -> np.ndarray:
        """
        Identify mesh faces that require supports based on normal angle and downward orientation.

        Returns:
            Boolean array of shape (n_faces,) where True indicates a face needing support
        """
        if self.mesh is None:
            raise ValueError("No mesh loaded — call load_stl() first")

        # Get face normals (unit vectors pointing outward)
        face_normals = self.mesh.face_normals
        # Calculate angle between face normal and negative Z-axis (downward direction)
        # Dot product of normal and (0,0,-1) gives cos(angle) since both are unit vectors
        dot_products = np.dot(face_normals, np.array([0, 0, -1]))
        # Clamp dot products to [-1, 1] to avoid floating point errors
        dot_products = np.clip(dot_products, -1.0, 1.0)
        angles = np.arccos(dot_products)

        # Overhangs are faces pointing downward (angle < 90° from negative Z) and below threshold
        self.overhang_faces = (angles < self.overhang_angle_threshold) & (angles > 0)
        overhang_count = np.sum(self.overhang_faces)
        logger.info(f"Identified {overhang_count} overhang faces out of {len(self.mesh.faces)} total")

        # Filter out small overhang regions below minimum area threshold
        overhang_mesh = self.mesh.submesh([self.overhang_faces], append=True)
        if overhang_mesh.area < self.min_overhang_area:
            logger.info(f"Overhang area {overhang_mesh.area:.2f}mm² below threshold {self.min_overhang_area}mm² — no supports needed")
            self.overhang_faces = np.zeros_like(self.overhang_faces, dtype=bool)

        return self.overhang_faces

if __name__ == "__main__":
    # Example usage
    try:
        analyzer = OverhangAnalyzer(overhang_angle_threshold=50.0, min_overhang_area=5.0)
        analyzer.load_stl("test_cube_overhang.stl")
        overhangs = analyzer.calculate_overhangs()
        print(f"Supports needed: {np.sum(overhangs) > 0}")
    except Exception as e:
        logger.error(f"Execution failed: {str(e)}")
        sys.exit(1)

import numpy as np
import trimesh
from scipy.spatial import Voronoi, voronoi_plot_2d
from typing import Optional, List
import logging

logging.basicConfig(level=logging.INFO)
logger = logging.getLogger(__name__)

class SupportGenerator:
    """Generates Voronoi-based optimized support structures for identified overhang regions."""

    def __init__(self, support_density: float = 0.15, support_thickness: float = 0.8, 
                 base_layer_height: float = 0.2, max_support_height: float = 50.0):
        """
        Initialize support generator with physical print constraints.

        Args:
            support_density: Fraction of overhang area to cover with supports (0.0-1.0)
            support_thickness: Thickness of support struts in mm (must be >= 0.4mm for PLA)
            base_layer_height: Layer height for support base adhesion (default 0.2mm)
            max_support_height: Maximum height of support structures in mm (supports taller than this are split)
        """
        if not 0.0 <= support_density <= 1.0:
            raise ValueError(f"support_density must be between 0 and 1, got {support_density}")
        if support_thickness < 0.4:
            logger.warning(f"support_thickness {support_thickness}mm is below 0.4mm minimum for reliable PLA printing")
        self.support_density = support_density
        self.support_thickness = support_thickness
        self.base_layer_height = base_layer_height
        self.max_support_height = max_support_height
        self.overhang_mesh = None
        self.support_mesh = None

    def generate_supports(self, overhang_mesh: trimesh.Trimesh, build_plate_z: float = 0.0) -> trimesh.Trimesh:
        """
        Generate support structures for a given overhang mesh using Voronoi tessellation.

        Args:
            overhang_mesh: Trimesh object containing only overhang faces (from OverhangAnalyzer)
            build_plate_z: Z-coordinate of the build plate (default 0.0 for standard printers)

        Returns:
            Trimesh object containing generated support structures
        """
        if overhang_mesh is None or overhang_mesh.area == 0:
            logger.info("No overhang mesh provided — returning empty support mesh")
            return trimesh.Trimesh()

        self.overhang_mesh = overhang_mesh
        # Project overhang vertices onto build plate (XY plane at build_plate_z)
        overhang_vertices = overhang_mesh.vertices
        # Filter vertices to only those in overhang faces (redundant but safe)
        projected_points = overhang_vertices[:, :2]  # Drop Z coordinate for 2D Voronoi

        # Subsample points based on support density to reduce computational load
        num_points = int(len(projected_points) * self.support_density)
        if num_points < 3:
            logger.warning(f"Subsampled points {num_points} < 3 — Voronoi tessellation requires at least 3 points")
            num_points = 3
        sampled_indices = np.random.choice(len(projected_points), size=num_points, replace=False)
        sampled_points = projected_points[sampled_indices]

        # Generate 2D Voronoi diagram for support scaffolding
        try:
            vor = Voronoi(sampled_points)
        except Exception as e:
            logger.error(f"Voronoi tessellation failed: {str(e)}")
            raise

        # Convert Voronoi edges to 3D support struts
        support_vertices = []
        support_faces = []
        vertex_count = 0

        for ridge in vor.ridge_vertices:
            if -1 in ridge:  # Ridge extends to infinity, skip
                continue
            # Get start and end points of Voronoi edge (2D)
            p1_2d = vor.vertices[ridge[0]]
            p2_2d = vor.vertices[ridge[1]]

            # Calculate Z range: from build plate to overhang bottom
            overhang_z_min = np.min(overhang_vertices[:, 2])
            overhang_z_max = np.max(overhang_vertices[:, 2])
            # Support height is from build plate to overhang bottom, capped at max_support_height
            support_height = min(overhang_z_min - build_plate_z, self.max_support_height)

            # Create 3D strut: rectangular prism along the Voronoi edge
            # Strut direction vector (2D)
            edge_dir = p2_2d - p1_2d
            edge_length = np.linalg.norm(edge_dir)
            if edge_length < 0.1:  # Skip tiny edges
                continue
            edge_dir_normalized = edge_dir / edge_length
            # Perpendicular vector for strut thickness
            perp_dir = np.array([-edge_dir_normalized[1], edge_dir_normalized[0]])

            # Define 8 vertices of the rectangular prism strut
            half_thickness = self.support_thickness / 2.0
            # Bottom face (at build plate)
            v0 = np.array([p1_2d[0] + perp_dir[0]*half_thickness, p1_2d[1] + perp_dir[1]*half_thickness, build_plate_z])
            v1 = np.array([p1_2d[0] - perp_dir[0]*half_thickness, p1_2d[1] - perp_dir[1]*half_thickness, build_plate_z])
            v2 = np.array([p2_2d[0] - perp_dir[0]*half_thickness, p2_2d[1] - perp_dir[1]*half_thickness, build_plate_z])
            v3 = np.array([p2_2d[0] + perp_dir[0]*half_thickness, p2_2d[1] + perp_dir[1]*half_thickness, build_plate_z])
            # Top face (at support height)
            v4 = np.array([p1_2d[0] + perp_dir[0]*half_thickness, p1_2d[1] + perp_dir[1]*half_thickness, build_plate_z + support_height])
            v5 = np.array([p1_2d[0] - perp_dir[0]*half_thickness, p1_2d[1] - perp_dir[1]*half_thickness, build_plate_z + support_height])
            v6 = np.array([p2_2d[0] - perp_dir[0]*half_thickness, p2_2d[1] - perp_dir[1]*half_thickness, build_plate_z + support_height])
            v7 = np.array([p2_2d[0] + perp_dir[0]*half_thickness, p2_2d[1] + perp_dir[1]*half_thickness, build_plate_z + support_height])

            # Add vertices to support mesh
            support_vertices.extend([v0, v1, v2, v3, v4, v5, v6, v7])
            # Define 12 faces of the rectangular prism (2 triangles per face)
            base = vertex_count
            # Bottom face
            support_faces.append([base+0, base+1, base+2])
            support_faces.append([base+0, base+2, base+3])
            # Top face
            support_faces.append([base+4, base+5, base+6])
            support_faces.append([base+4, base+6, base+7])
            # Front face
            support_faces.append([base+0, base+1, base+5])
            support_faces.append([base+0, base+5, base+4])
            # Back face
            support_faces.append([base+2, base+3, base+7])
            support_faces.append([base+2, base+7, base+6])
            # Left face
            support_faces.append([base+1, base+2, base+6])
            support_faces.append([base+1, base+6, base+5])
            # Right face
            support_faces.append([base+0, base+3, base+7])
            support_faces.append([base+0, base+7, base+4])

            vertex_count += 8

        # Create final support mesh
        if len(support_vertices) == 0:
            logger.warning("No support struts generated — check support density and overhang mesh")
            return trimesh.Trimesh()

        self.support_mesh = trimesh.Trimesh(vertices=np.array(support_vertices), faces=np.array(support_faces))
        # Merge supports with overhang mesh for visualization (optional)
        combined = trimesh.util.concatenate(self.overhang_mesh, self.support_mesh)
        logger.info(f"Generated {vertex_count//8} support struts | Support volume: {self.support_mesh.volume:.2f}mm³")
        return self.support_mesh

if __name__ == "__main__":
    try:
        # Example: generate supports for a simple overhang mesh
        test_mesh = trimesh.creation.box(extents=[10, 10, 20])
        # Create fake overhang mesh (top 5mm of box)
        overhang_mask = test_mesh.vertices[:, 2] > 15.0
        overhang_faces = []
        for face in test_mesh.faces:
            if any(overhang_mask[v] for v in face):
                overhang_faces.append(face)
        overhang_mesh = test_mesh.submesh([np.isin(np.arange(len(test_mesh.faces)), overhang_faces)], append=True)

        generator = SupportGenerator(support_density=0.2, support_thickness=0.8)
        supports = generator.generate_supports(overhang_mesh)
        print(f"Support mesh vertices: {len(supports.vertices)}")
    except Exception as e:
        logger.error(f"Support generation failed: {str(e)}")
        sys.exit(1)

import subprocess
import os
import json
import tempfile
import logging
from typing import Dict, Optional
import sys

logging.basicConfig(level=logging.INFO)
logger = logging.getLogger(__name__)

class CuraSlicer:
    """Integrates with CuraEngine to slice meshes with optimized support settings."""

    # Default Cura profile for PLA 0.2mm layer height (simplified for example)
    DEFAULT_PROFILE = {
        "layer_height": 0.2,
        "wall_thickness": 0.8,
        "infill_density": 20,
        "support_enable": True,
        "support_type": "normal",
        "support_angle": 45,
        "support_pattern": "grid",
        "support_density": 15,
        "support_z_distance": 0.2,
        "material_print_temperature": 200,
        "material_bed_temperature": 60
    }

    def __init__(self, cura_engine_path: str = "/usr/bin/CuraEngine", profile: Optional[Dict] = None):
        """
        Initialize slicer with CuraEngine path and print profile.

        Args:
            cura_engine_path: Absolute path to CuraEngine executable
            profile: Custom print profile (overrides DEFAULT_PROFILE)
        """
        if not os.path.exists(cura_engine_path):
            raise FileNotFoundError(f"CuraEngine not found at {cura_engine_path}")
        self.cura_engine_path = cura_engine_path
        self.profile = {**self.DEFAULT_PROFILE, **(profile or {})}
        logger.info(f"Initialized CuraSlicer with engine: {cura_engine_path}")

    def _write_profile_file(self, profile: Dict, output_path: str) -> None:
        """Write print profile to CuraEngine-compatible JSON file."""
        try:
            with open(output_path, 'w') as f:
                json.dump(profile, f, indent=2)
            logger.info(f"Wrote profile to {output_path}")
        except Exception as e:
            logger.error(f"Failed to write profile file: {str(e)}")
            raise

    def slice_mesh(self, mesh_path: str, output_gcode_path: str, support_mesh_path: Optional[str] = None) -> None:
        """
        Slice input mesh (with optional support mesh) to G-code using CuraEngine.

        Args:
            mesh_path: Path to main STL mesh file
            output_gcode_path: Path to output G-code file
            support_mesh_path: Optional path to support STL mesh (if None, uses Cura built-in supports)
        """
        if not os.path.exists(mesh_path):
            raise FileNotFoundError(f"Mesh file not found: {mesh_path}")

        # Create temporary profile file
        with tempfile.NamedTemporaryFile(mode='w', suffix='.json', delete=False) as tmp_profile:
            profile_path = tmp_profile.name
            self._write_profile_file(self.profile, profile_path)

        # Build CuraEngine command
        cmd = [
            self.cura_engine_path,
            "slice",
            "-v",
            "-j", profile_path,
            "-o", output_gcode_path,
            "-l", mesh_path
        ]

        # Add support mesh if provided
        if support_mesh_path is not None:
            if not os.path.exists(support_mesh_path):
                raise FileNotFoundError(f"Support mesh not found: {support_mesh_path}")
            cmd.extend(["-l", support_mesh_path])

        logger.info(f"Running CuraEngine command: {' '.join(cmd)}")

        try:
            # Execute CuraEngine, capture output
            result = subprocess.run(
                cmd,
                stdout=subprocess.PIPE,
                stderr=subprocess.PIPE,
                text=True,
                timeout=300  # 5 minute timeout for large meshes
            )

            if result.returncode != 0:
                logger.error(f"CuraEngine failed with return code {result.returncode}")
                logger.error(f"Stderr: {result.stderr}")
                raise RuntimeError(f"Slicing failed: {result.stderr}")

            logger.info(f"Successfully sliced to {output_gcode_path}")
            logger.info(f"CuraEngine stdout: {result.stdout[:200]}...")  # Log first 200 chars

        except subprocess.TimeoutExpired:
            logger.error("CuraEngine timed out after 300 seconds")
            raise
        except Exception as e:
            logger.error(f"Slicing failed: {str(e)}")
            raise
        finally:
            # Clean up temporary profile file
            if os.path.exists(profile_path):
                os.unlink(profile_path)

    def compare_slice_metrics(self, default_gcode: str, optimized_gcode: str) -> Dict:
        """
        Compare metrics between default and optimized G-code files.

        Args:
            default_gcode: Path to G-code sliced with default supports
            optimized_gcode: Path to G-code sliced with optimized supports

        Returns:
            Dictionary of metric differences (e.g., material use, print time)
        """
        metrics = {}
        for gcode_path, key in [(default_gcode, "default"), (optimized_gcode, "optimized")]:
            if not os.path.exists(gcode_path):
                raise FileNotFoundError(f"G-code file not found: {gcode_path}")
            # Parse G-code for basic metrics (simplified — real implementation would use a G-code parser)
            with open(gcode_path, 'r') as f:
                gcode = f.read()
            # Count extrusion moves (G1 with E value)
            extrusion_moves = [line for line in gcode.split('\n') if line.startswith('G1') and 'E' in line]
            total_extrusion = sum(float(line.split('E')[1].split()[0]) for line in extrusion_moves if 'E' in line)
            # Estimate print time from G-code comments (Cura adds ;TIME:XXX and ;Filament used: XXXm)
            time_line = [line for line in gcode.split('\n') if line.startswith(';TIME:')]
            filament_line = [line for line in gcode.split('\n') if line.startswith(';Filament used:')]

            metrics[f"{key}_time_sec"] = int(time_line[0].split(':')[1]) if time_line else 0
            metrics[f"{key}_filament_mm"] = float(filament_line[0].split(':')[1].strip().replace('m', '')) if filament_line else 0.0
            metrics[f"{key}_extrusion_moves"] = len(extrusion_moves)

        # Calculate differences
        metrics["time_reduction_pct"] = ((metrics["default_time_sec"] - metrics["optimized_time_sec"]) / metrics["default_time_sec"]) * 100
        metrics["filament_reduction_pct"] = ((metrics["default_filament_mm"] - metrics["optimized_filament_mm"]) / metrics["default_filament_mm"]) * 100

        logger.info(f"Comparison results: Time reduction {metrics['time_reduction_pct']:.1f}%, Filament reduction {metrics['filament_reduction_pct']:.1f}%")
        return metrics

if __name__ == "__main__":
    try:
        slicer = CuraSlicer(cura_engine_path="/usr/local/bin/CuraEngine")
        # Slice with default supports
        slicer.slice_mesh("test_cube.stl", "default_support.gcode")
        # Slice with optimized supports (assuming we have optimized_support.stl from Step 2)
        slicer.slice_mesh("test_cube.stl", "optimized_support.gcode", "optimized_support.stl")
        # Compare metrics
        metrics = slicer.compare_slice_metrics("default_support.gcode", "optimized_support.gcode")
        print(f"Optimized print time: {metrics['optimized_time_sec']}s")
    except Exception as e:
        logger.error(f"Execution failed: {str(e)}")
        sys.exit(1)

# Add to load_stl() after loading mesh
if not self.mesh.is_watertight:
    logger.warning(f"Mesh {stl_path} is not watertight — attempting repair")
    try:
        trimesh.repair.fill_holes(self.mesh)
        if not self.mesh.is_watertight:
            raise ValueError(f"Mesh {stl_path} has irreparable non-manifold geometry")
        logger.info(f"Successfully repaired mesh {stl_path}")
    except Exception as e:
        logger.error(f"Mesh repair failed: {str(e)}")
        raise

def get_optimal_density(filament_type: str, max_print_speed: float) -> float:
    """Return optimal support density based on filament and print speed."""
    base_densities = {
        "PLA": 0.18,
        "ABS": 0.12,
        "PETG": 0.15,
        "TPU": 0.08  # Flexible filament needs very low density
    }
    density = base_densities.get(filament_type.upper(), 0.15)
    # Reduce density by 0.01 for every 50mm/s over 100mm/s max speed
    if max_print_speed > 100:
        density -= 0.01 * ((max_print_speed - 100) // 50)
    return max(0.05, min(density, 0.25))  # Clamp to valid range

import hashlib
import sqlite3

def get_cache_key(stl_path: str, density: float, angle: float, thickness: float) -> str:
    """Generate unique cache key for support calculation."""
    with open(stl_path, 'rb') as f:
        stl_hash = hashlib.md5(f.read()).hexdigest()
    return f"{stl_hash}_{density}_{angle}_{thickness}"

def cache_support_mesh(cache_key: str, support_mesh: trimesh.Trimesh, db_path: str = "support_cache.db"):
    """Store support mesh in SQLite cache."""
    conn = sqlite3.connect(db_path)
    c = conn.cursor()
    c.execute('''CREATE TABLE IF NOT EXISTS supports
                 (cache_key TEXT PRIMARY KEY, mesh_blob BLOB, timestamp DATETIME DEFAULT CURRENT_TIMESTAMP)''')
    # Serialize mesh to binary blob
    mesh_blob = support_mesh.export(file_type='stl').read()
    c.execute("INSERT OR REPLACE INTO supports (cache_key, mesh_blob) VALUES (?, ?)", (cache_key, mesh_blob))
    conn.commit()
    conn.close()

support-optimizer/
├── src/
│   ├── __init__.py
│   ├── overhang_analyzer.py  # Step 1 code: STL parsing and overhang detection
│   ├── support_generator.py  # Step 2 code: Voronoi support generation
│   ├── cura_slicer.py        # Step 3 code: CuraEngine integration
│   └── utils.py              # Shared utilities (caching, metrics)
├── profiles/
│   ├── pla_0.2mm.json        # Default PLA profile
│   ├── abs_0.2mm.json        # ABS profile
│   └── petg_0.2mm.json       # PETG profile
├── tests/
│   ├── test_overhang_analyzer.py
│   ├── test_support_generator.py
│   └── test_cura_slicer.py
├── benchmarks/
│   ├── stl_dataset/          # 120+ test STL files
│   └── results/              # Benchmark reports
├── requirements.txt          # Python dependencies (trimesh, numpy, scipy, etc.)
├── LICENSE                   # MIT License
└── README.md                 # Setup and usage instructions