From 4bec4063a5b88302b2c6cac1156659bafb3ad64e Mon Sep 17 00:00:00 2001 From: Antoine Date: Mon, 16 Feb 2026 00:01:35 +0000 Subject: [PATCH] =?UTF-8?q?feat:=20add=20adaptive=20isogrid=20tool=20?= =?UTF-8?q?=E2=80=94=20project=20foundations?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit - Python Brain: density field, constrained Delaunay triangulation, pocket profiles, profile assembly, validation modules - NX Hands: skeleton scripts for geometry extraction, AFEM setup, per-iteration solve (require NX environment to develop) - Atomizer integration: 15-param space definition, objective function - Technical spec, README, sample test geometry, requirements.txt - Architecture: Python Brain + NX Hands + Atomizer Manager --- tools/adaptive-isogrid/README.md | 75 ++ tools/adaptive-isogrid/docs/technical-spec.md | 1128 +++++++++++++++++ tools/adaptive-isogrid/requirements.txt | 5 + tools/adaptive-isogrid/src/__init__.py | 0 tools/adaptive-isogrid/src/atomizer_study.py | 48 + tools/adaptive-isogrid/src/brain/__init__.py | 8 + .../src/brain/density_field.py | 146 +++ .../src/brain/pocket_profiles.py | 176 +++ .../src/brain/profile_assembly.py | 100 ++ .../src/brain/triangulation.py | 160 +++ .../adaptive-isogrid/src/brain/validation.py | 106 ++ tools/adaptive-isogrid/src/nx/__init__.py | 6 + .../src/nx/build_interface_model.py | 28 + .../src/nx/extract_geometry.py | 54 + .../src/nx/iteration_solve.py | 35 + .../tests/test_geometries/sample_bracket.json | 49 + 16 files changed, 2124 insertions(+) create mode 100644 tools/adaptive-isogrid/README.md create mode 100644 tools/adaptive-isogrid/docs/technical-spec.md create mode 100644 tools/adaptive-isogrid/requirements.txt create mode 100644 tools/adaptive-isogrid/src/__init__.py create mode 100644 tools/adaptive-isogrid/src/atomizer_study.py create mode 100644 tools/adaptive-isogrid/src/brain/__init__.py create mode 100644 tools/adaptive-isogrid/src/brain/density_field.py create mode 100644 tools/adaptive-isogrid/src/brain/pocket_profiles.py create mode 100644 tools/adaptive-isogrid/src/brain/profile_assembly.py create mode 100644 tools/adaptive-isogrid/src/brain/triangulation.py create mode 100644 tools/adaptive-isogrid/src/brain/validation.py create mode 100644 tools/adaptive-isogrid/src/nx/__init__.py create mode 100644 tools/adaptive-isogrid/src/nx/build_interface_model.py create mode 100644 tools/adaptive-isogrid/src/nx/extract_geometry.py create mode 100644 tools/adaptive-isogrid/src/nx/iteration_solve.py create mode 100644 tools/adaptive-isogrid/tests/test_geometries/sample_bracket.json diff --git a/tools/adaptive-isogrid/README.md b/tools/adaptive-isogrid/README.md new file mode 100644 index 00000000..43eecadb --- /dev/null +++ b/tools/adaptive-isogrid/README.md @@ -0,0 +1,75 @@ +# Adaptive Isogrid — Plate Lightweighting Tool + +**Status:** Foundation / Pre-Implementation +**Architecture:** Python Brain + NX Hands + Atomizer Manager + +## What It Does + +Takes a plate with holes → generates an optimally lightweighted isogrid pattern → produces manufacturing-ready geometry. Isogrid density varies across the plate based on hole importance, edge proximity, and optimization-driven meta-parameters. + +## Architecture + +| Component | Role | Runtime | +|-----------|------|---------| +| **Python Brain** | Density field → Constrained Delaunay → rib profile | ~1-3 sec | +| **NX Hands** | Import profile → mesh → AFEM merge → Nastran solve → extract results | ~60-90 sec | +| **Atomizer Manager** | Optuna TPE sampling → objective evaluation → convergence | 500-2000 trials | + +### Key Insight: Assembly FEM with Superposed Models + +- **Model A** (permanent): Spider elements at holes + edge BC nodes. All loads/BCs applied here. +- **Model B** (variable): 2D shell mesh of ribbed plate. Rebuilt each iteration. +- **Node merge** at fixed interface locations connects them reliably every time. + +Loads and BCs never need re-association. Only the rib pattern changes. + +## Directory Structure + +``` +adaptive-isogrid/ +├── README.md +├── requirements.txt +├── docs/ +│ └── technical-spec.md # Full architecture spec +├── src/ +│ ├── brain/ # Python geometry generator +│ │ ├── __init__.py +│ │ ├── density_field.py # η(x) evaluation +│ │ ├── triangulation.py # Constrained Delaunay + refinement +│ │ ├── pocket_profiles.py # Pocket inset + filleting +│ │ ├── profile_assembly.py # Final plate - pockets - holes +│ │ └── validation.py # Manufacturing constraint checks +│ ├── nx/ # NXOpen journal scripts +│ │ ├── extract_geometry.py # One-time: face → geometry.json +│ │ ├── build_interface_model.py # One-time: Model A + spiders +│ │ └── iteration_solve.py # Per-trial: rebuild Model B + solve +│ └── atomizer_study.py # Atomizer/Optuna integration +└── tests/ + └── test_geometries/ # Sample geometry.json files +``` + +## Implementation Phases + +1. **Python Brain standalone** (1-2 weeks) — geometry generator with matplotlib viz +2. **NX extraction + AFEM setup** (1-2 weeks) — one-time project setup scripts +3. **NX iteration script** (1-2 weeks) — per-trial mesh/solve/extract loop +4. **Atomizer integration** (1 week) — wire objective function + study management +5. **Validation + first real project** (1-2 weeks) — production run on client plate + +## Quick Start (Phase 1) + +```bash +cd tools/adaptive-isogrid +pip install -r requirements.txt +python -m src.brain --geometry tests/test_geometries/sample_bracket.json --params default +``` + +## Parameter Space + +15 continuous parameters optimized by Atomizer (Optuna TPE): +- Density field: η₀, α, R₀, κ, p, β, R_edge +- Spacing: s_min, s_max +- Rib thickness: t_min, t₀, γ +- Manufacturing: w_frame, r_f, d_keep + +See `docs/technical-spec.md` for full formulation. diff --git a/tools/adaptive-isogrid/docs/technical-spec.md b/tools/adaptive-isogrid/docs/technical-spec.md new file mode 100644 index 00000000..7a6599c5 --- /dev/null +++ b/tools/adaptive-isogrid/docs/technical-spec.md @@ -0,0 +1,1128 @@ +# Adaptive Isogrid Plate Lightweighting — Technical Specification + +## System Architecture: "Python Brain + NX Hands + Atomizer Manager" + +**Author:** Atomaste Solution +**Date:** February 2026 +**Status:** Architecture Locked — Ready for Implementation + +--- + +## 1. Project Summary + +### What We're Building + +A semi-automated tool that takes a plate with holes, generates an optimally lightweighted isogrid pattern, and produces a manufacturing-ready geometry. The isogrid density varies across the plate based on hole importance, edge proximity, and optimization-driven meta-parameters. + +### Architecture Decision Record + +After extensive brainstorming, the following decisions are locked: + +| Decision | Choice | Rationale | +|---|---|---| +| Geometry generation | External Python (Constrained Delaunay) | Full access to scipy/triangle/gmsh, debuggable, fast | +| FEA strategy | Assembly FEM with superposed models | Decouples fixed interfaces (loads/BCs) from variable rib pattern | +| FEA solver | NX Simcenter + Nastran (2D shell) | Existing expertise, handles complex BCs, extensible to modal/buckling | +| NX role | Remesh plate model + merge nodes + solve | Geometry is agnostic to rib pattern; loads/BCs never re-associated | +| Optimization | Atomizer (Optuna TPE), pure parametric v1 | One FEA per trial, ~2 min/iteration, stress feedback deferred to v2 | +| Manufacturing | Through-cuts only (waterjet + CNC finish) | Simplifies geometry to 2D profile problem | +| Plate type | Flat, 200–600 mm scale, 6–15 mm thick, 16–30 holes | Shell elements appropriate, fast solves | + +### System Diagram + +``` +ONE-TIME SETUP +══════════════ +User in NX: + ├── Select plate face → tool extracts boundary + holes → geometry.json + ├── Assign hole weights (interactive UI or table) + │ + ├── Build Assembly FEM: + │ ├── Model A — "Interface Model" (PERMANENT, never changes) + │ │ ├── Spider elements (RBE2/RBE3) at each hole center → circumference + │ │ ├── All loads applied to spider center nodes + │ │ ├── All BCs applied to spider center nodes or edge nodes + │ │ └── Edge BC nodes along plate boundary + │ │ + │ ├── Model B — "Plate Model" (VARIABLE, rebuilt each iteration) + │ │ ├── 2D shell mesh of the ribbed plate profile + │ │ ├── Mesh seeds at hole circumference nodes (match Model A spiders) + │ │ └── Mesh seeds at edge BC nodes (match Model A edge nodes) + │ │ + │ └── Assembly FEM + │ ├── Superpose Model A + Model B + │ ├── Merge nodes at hole circumferences + edges + │ └── Solver settings (SOL 101, SOL 103, etc.) + │ + └── Verify: solve dummy Model B (solid plate), confirm results + +OPTIMIZATION LOOP (~2 min/iteration) +════════════════════════════════════ +Atomizer (Optuna TPE) + │ + │ Samples meta-parameters + │ (η₀, α, β, p, R₀, κ, s_min, s_max, t_min, t₀, γ, w_frame, r_f, ...) + │ + ▼ +External Python — "The Brain" (~1-3 sec) + ├── Load plate geometry + hole table (fixed) + ├── Compute density field η(x) + ├── Generate constrained Delaunay triangulation + ├── Compute rib thicknesses from density + ├── Apply manufacturing constraints + ├── Validate geometry (min web, min pocket, no degenerates) + ├── Output: rib_profile.json (curve coordinates) + └── Output: validity_flag + mass_estimate + │ + ▼ +NXOpen Journal — "The Hands" (~60-90 sec) + ├── Delete old Model B geometry + mesh + ├── Import new 2D ribbed profile (from rib_profile.json) + ├── Create sheet body from profile curves + ├── Mesh Model B with seeds at fixed interface nodes + ├── Merge nodes in Assembly FEM (holes + edges) + ├── Solve (Nastran) + └── Extract results (stress, displacement, strain, mass) + │ + ▼ +Result Extraction (~5-10 sec) + ├── Von Mises stress field (nodal) + ├── Displacement magnitude field (nodal) + ├── Strain field (elemental) + ├── Total mass + ├── (Future: modal frequencies, buckling load factors) + └── Report metrics to Atomizer + │ + ▼ +Atomizer updates surrogate, samples next trial +Repeat until convergence (~500-2000 trials) +``` + +### Key Architectural Insight: Why Assembly FEM + +The plate lightweighting problem has a natural separation: + +**What doesn't change:** Hole positions, hole diameters, plate boundary, loads, boundary conditions. These are the structural interfaces — where forces enter and leave the plate. + +**What changes every iteration:** The rib pattern between holes. This is purely the internal load-path topology. + +The Assembly FEM with superposed models exploits this separation directly. Model A captures everything fixed (interfaces, loads, BCs) using spider elements at holes and edge nodes at boundaries. Model B captures everything variable (the ribbed plate mesh). They connect through node merging at the fixed interface locations (hole circumferences and plate edges). + +This means: +- Loads and BCs **never need re-association** — they live permanently on Model A +- Only Model B's mesh is rebuilt each iteration +- Node merging at fixed geometric locations is a reliable, automated operation in NX +- The solver sees one connected model with proper load paths through the spiders into the ribs +Atomizer updates surrogate, samples next trial +Repeat until convergence (~500-2000 trials) +``` + +--- + +## 2. One-Time Setup: Geometry Extraction from NX + +### 2.1 What the User Does + +The user opens their plate model in NX and runs a setup script (NXOpen Python). The interaction is: + +1. **Select the plate face** — click the top (or bottom) face of the plate. The tool reads this face's topology. + +2. **Hole classification** — the tool lists all inner loops (holes) found on the face, showing each hole's center, diameter, and a preview highlight. The user assigns each hole a weight from 0.0 (ignore — just avoid it) to 1.0 (critical — maximum reinforcement). Grouping by class (A/B/C) is optional; raw per-hole weights work fine since they're not optimization variables. + +3. **Review** — the tool displays the extracted boundary and holes overlaid on the model for visual confirmation. + +4. **Export** — writes `geometry.json` containing everything the Python brain needs. + +### 2.2 Geometry Extraction Logic (NXOpen) + +The plate face in NX is a B-rep face bounded by edge loops. Extraction pseudocode: + +```python +# NXOpen extraction script (runs inside NX) +import NXOpen +import json + +def extract_plate_geometry(face, hole_weights): + """ + face: NXOpen.Face — the selected plate face + hole_weights: dict — {loop_index: weight} from user input + Returns: geometry dict for export + """ + # Get all edge loops on this face + loops = face.GetLoops() + + geometry = { + 'outer_boundary': None, + 'holes': [], + 'face_normal': None, + 'thickness': None # can be read from plate body + } + + for loop in loops: + edges = loop.GetEdges() + # Sample each edge as a polyline + points = [] + for edge in edges: + # Get edge curve, sample at intervals + pts = sample_edge(edge, tolerance=0.1) # 0.1 mm chord tol + points.extend(pts) + + if loop.IsOuter(): + geometry['outer_boundary'] = points + else: + # Inner loop = hole + center, diameter = fit_circle(points) # for circular holes + hole_idx = len(geometry['holes']) + geometry['holes'].append({ + 'index': hole_idx, + 'boundary': points, # actual boundary polyline + 'center': center, # [x, y] + 'diameter': diameter, # mm (None if non-circular) + 'is_circular': is_circle(points, tolerance=0.5), + 'weight': hole_weights.get(hole_idx, 0.0) + }) + + # Get plate thickness from body + geometry['thickness'] = measure_plate_thickness(face) + + # Get face normal and establish XY coordinate system + geometry['face_normal'] = get_face_normal(face) + geometry['transform'] = get_face_to_xy_transform(face) + + return geometry + +def export_geometry(geometry, filepath='geometry.json'): + with open(filepath, 'w') as f: + json.dump(geometry, f, indent=2) +``` + +### 2.3 What geometry.json Contains + +```json +{ + "plate_id": "bracket_v3", + "units": "mm", + "thickness": 10.0, + "material": "AL6061-T6", + "outer_boundary": [[0,0], [400,0], [400,300], [0,300]], + "holes": [ + { + "index": 0, + "center": [50, 50], + "diameter": 12.0, + "is_circular": true, + "boundary": [[44,50], [44.2,51.8], ...], + "weight": 1.0 + }, + { + "index": 1, + "center": [200, 150], + "diameter": 8.0, + "is_circular": true, + "boundary": [...], + "weight": 0.3 + } + ] +} +``` + +Non-circular holes (slots, irregular cutouts) carry their full boundary polyline and weight, but `diameter` and `is_circular` are null/false. The density field uses point-to-polygon distance instead of point-to-center distance for these. + +--- + +## 3. The Python Brain: Density Field + Geometry Generation + +### 3.1 Density Field Formulation + +The density field η(x) is the core of the method. It maps every point on the plate to a value between 0 (minimum density — remove maximum material) and 1 (maximum density — retain material). + +**Hole influence term:** + +``` +I(x) = Σᵢ wᵢ · exp( -(dᵢ(x) / Rᵢ)^p ) +``` + +Where: +- `dᵢ(x)` = distance from point x to hole i (center-to-point for circular, boundary-to-point for non-circular) +- `Rᵢ = R₀ · (1 + κ · wᵢ)` = influence radius, scales with hole importance +- `p` = decay exponent (controls transition sharpness) +- `wᵢ` = user-assigned hole weight (fixed, not optimized) + +**Edge reinforcement term:** + +``` +E(x) = exp( -(d_edge(x) / R_edge)^p_edge ) +``` + +Where `d_edge(x)` is the distance from x to the nearest plate boundary edge. + +**Combined density field:** + +``` +η(x) = clamp(0, 1, η₀ + α · I(x) + β · E(x)) +``` + +**Density to local target spacing:** + +``` +s(x) = s_max - (s_max - s_min) · η(x) +``` + +Where s(x) is the target triangle edge length at point x. High density → small spacing (more ribs). Low density → large spacing (fewer ribs). + +**Density to local rib thickness:** + +``` +t(x) = clamp(t_min, t_max, t₀ · (1 + γ · η(x))) +``` + +Where t₀ is the nominal rib thickness and γ controls how much density affects thickness. + +### 3.2 Geometry Generation: Constrained Delaunay Pipeline + +The geometry generation pipeline converts the density field into a manufacturable 2D rib profile. The recommended approach uses Jonathan Shewchuk's Triangle library (Python binding: `triangle`) for constrained Delaunay triangulation with area constraints. + +**Step 1 — Define the Planar Straight Line Graph (PSLG):** + +The PSLG is the input to Triangle. It consists of: +- The outer boundary as a polygon (vertices + segments) +- Each hole boundary as a polygon (vertices + segments) +- Hole markers (points inside each hole, telling Triangle to leave these regions empty) + +```python +import triangle +import numpy as np + +def build_pslg(geometry, keepout_distance): + """ + Build PSLG from plate geometry. + keepout_distance: extra clearance around holes (mm) + """ + vertices = [] + segments = [] + holes_markers = [] + + # Outer boundary + outer = offset_inward(geometry['outer_boundary'], keepout_distance) + v_start = len(vertices) + vertices.extend(outer) + for i in range(len(outer)): + segments.append([v_start + i, v_start + (i+1) % len(outer)]) + + # Each hole boundary (offset outward for keepout) + for hole in geometry['holes']: + hole_boundary = offset_outward(hole['boundary'], keepout_distance) + v_start = len(vertices) + vertices.extend(hole_boundary) + for i in range(len(hole_boundary)): + segments.append([v_start + i, v_start + (i+1) % len(hole_boundary)]) + holes_markers.append(hole['center']) # point inside hole + + return { + 'vertices': np.array(vertices), + 'segments': np.array(segments), + 'holes': np.array(holes_markers) + } +``` + +**Step 2 — Compute area constraints from density field:** + +Triangle supports per-region area constraints via a callback or a maximum area parameter. For spatially varying area, we use an iterative refinement approach: + +```python +def compute_max_area(x, y, params): + """ + Target triangle area at point (x,y) based on density field. + Smaller area = denser triangulation = more ribs. + """ + eta = evaluate_density_field(x, y, params) + s = params['s_max'] - (params['s_max'] - params['s_min']) * eta + # Area of equilateral triangle with side length s + target_area = (np.sqrt(3) / 4) * s**2 + return target_area +``` + +**Step 3 — Run constrained Delaunay triangulation:** + +```python +def generate_triangulation(pslg, params): + """ + Generate adaptive triangulation using Triangle library. + """ + # Initial triangulation with global max area + global_max_area = (np.sqrt(3) / 4) * params['s_max']**2 + + # Triangle options: + # 'p' = triangulate PSLG + # 'q30' = minimum angle 30° (quality mesh) + # 'a' = area constraint + # 'D' = conforming Delaunay + result = triangle.triangulate(pslg, f'pq30Da{global_max_area}') + + # Iterative refinement based on density field + for iteration in range(3): # 2-3 refinement passes + # For each triangle, check if area exceeds local target + triangles = result['triangles'] + vertices = result['vertices'] + + areas = compute_triangle_areas(vertices, triangles) + centroids = compute_centroids(vertices, triangles) + + # Build per-triangle area constraints + max_areas = np.array([ + compute_max_area(cx, cy, params) + for cx, cy in centroids + ]) + + # If all triangles satisfy constraints, done + if np.all(areas <= max_areas * 1.2): # 20% tolerance + break + + # Refine: set area constraint and re-triangulate + # Triangle supports this via the 'r' (refine) flag + result = triangle.triangulate( + result, f'rpq30Da', + # per-triangle area constraints via triangle_max_area + ) + + return result +``` + +**Step 4 — Extract ribs and compute thicknesses:** + +```python +def extract_ribs(triangulation, params, geometry): + """ + Convert triangulation edges to rib definitions. + Each rib = (start_point, end_point, thickness, midpoint_density) + """ + vertices = triangulation['vertices'] + triangles = triangulation['triangles'] + + # Get unique edges from triangle connectivity + edges = set() + for tri in triangles: + for i in range(3): + edge = tuple(sorted([tri[i], tri[(i+1)%3]])) + edges.add(edge) + + ribs = [] + for v1_idx, v2_idx in edges: + p1 = vertices[v1_idx] + p2 = vertices[v2_idx] + midpoint = (p1 + p2) / 2 + + # Skip edges on the boundary (these aren't interior ribs) + if is_boundary_edge(v1_idx, v2_idx, triangulation): + continue + + # Compute local density and rib thickness + eta = evaluate_density_field(midpoint[0], midpoint[1], params) + thickness = compute_rib_thickness(eta, params) + + ribs.append({ + 'start': p1.tolist(), + 'end': p2.tolist(), + 'midpoint': midpoint.tolist(), + 'thickness': thickness, + 'density': eta + }) + + return ribs +``` + +**Step 5 — Generate pocket profiles:** + +Each triangle in the triangulation defines a pocket. The pocket profile is the triangle inset by half the local rib thickness on each edge, with fillet radii at corners. + +```python +def generate_pocket_profiles(triangulation, ribs, params): + """ + For each triangle, compute the pocket outline + (triangle boundary inset by half-rib-width on each edge). + """ + vertices = triangulation['vertices'] + triangles = triangulation['triangles'] + + pockets = [] + for tri_idx, tri in enumerate(triangles): + # Get the three edge thicknesses + edge_thicknesses = get_triangle_edge_thicknesses( + tri, ribs, vertices + ) + + # Inset each edge by half its rib thickness + inset_polygon = inset_triangle( + vertices[tri[0]], vertices[tri[1]], vertices[tri[2]], + edge_thicknesses[0]/2, edge_thicknesses[1]/2, edge_thicknesses[2]/2 + ) + + if inset_polygon is None: + # Triangle too small for pocket — skip (solid region) + continue + + # Check minimum pocket size + inscribed_r = compute_inscribed_radius(inset_polygon) + if inscribed_r < params.get('min_pocket_radius', 1.5): + continue # pocket too small to manufacture + + # Apply fillet to pocket corners + filleted = fillet_polygon(inset_polygon, params['r_f']) + + pockets.append({ + 'triangle_index': tri_idx, + 'vertices': filleted, + 'area': polygon_area(filleted) + }) + + return pockets +``` + +**Step 6 — Assemble the ribbed plate profile:** + +The final output is the plate boundary minus all pocket regions, plus the hole cutouts. This is a 2D profile that NX will mesh as shells. + +```python +def assemble_profile(geometry, pockets, params): + """ + Create the final 2D ribbed plate profile. + Plate boundary - pockets - holes = ribbed plate + """ + from shapely.geometry import Polygon, MultiPolygon + from shapely.ops import unary_union + + # Plate outline (with optional perimeter frame) + plate = Polygon(geometry['outer_boundary']) + + # Inset plate by frame width + if params['w_frame'] > 0: + inner_plate = plate.buffer(-params['w_frame']) + else: + inner_plate = plate + + # Union all pocket polygons + pocket_polys = [Polygon(p['vertices']) for p in pockets] + all_pockets = unary_union(pocket_polys) + + # Clip pockets to inner plate (don't cut into frame) + clipped_pockets = all_pockets.intersection(inner_plate) + + # Subtract pockets from plate + ribbed_plate = plate.difference(clipped_pockets) + + # Subtract holes (with original hole boundaries) + for hole in geometry['holes']: + hole_poly = Polygon(hole['boundary']) + ribbed_plate = ribbed_plate.difference(hole_poly) + + return ribbed_plate +``` + +**Step 7 — Validate and export:** + +```python +def validate_and_export(ribbed_plate, params, output_path): + """ + Check manufacturability and export for NXOpen. + """ + checks = { + 'min_web_width': check_minimum_web(ribbed_plate, params['t_min']), + 'no_islands': check_no_floating_islands(ribbed_plate), + 'no_self_intersections': ribbed_plate.is_valid, + 'mass_estimate': estimate_mass(ribbed_plate, params), + } + + valid = all([ + checks['min_web_width'], + checks['no_islands'], + checks['no_self_intersections'] + ]) + + # Export as JSON (coordinate arrays for NXOpen) + profile_data = { + 'valid': valid, + 'checks': checks, + 'outer_boundary': list(ribbed_plate.exterior.coords), + 'pockets': [list(interior.coords) + for interior in ribbed_plate.interiors + if is_pocket(interior)], # pocket cutouts only + 'hole_boundaries': [list(interior.coords) + for interior in ribbed_plate.interiors + if is_hole(interior)], # original hole cutouts + 'mass_estimate': checks['mass_estimate'], + 'num_pockets': len([i for i in ribbed_plate.interiors if is_pocket(i)]), + 'parameters_used': params + } + + with open(output_path, 'w') as f: + json.dump(profile_data, f) + + return valid, checks +``` + +### 3.3 Manufacturing Constraint Summary + +These constraints are enforced during geometry generation, not as FEA post-checks: + +| Constraint | Value | Enforcement Point | +|---|---|---| +| Minimum rib width | t_min (param, ≥ 2.0 mm) | Rib thickness computation + validation | +| Minimum pocket inscribed radius | 1.5 mm (waterjet pierce requirement) | Pocket generation — skip small pockets | +| Corner fillet radius | r_f (param, ≥ 0.5 mm for waterjet, ≥ tool_radius for CNC) | Pocket profile filleting | +| Hole keepout | d_keep,hole (param, typically 1.5× hole diameter) | PSLG construction | +| Edge keepout / frame | w_frame (param) | Profile assembly | +| Minimum triangle quality | q_min = 30° minimum angle | Triangle library quality flag | +| No floating islands | — | Validation step | +| No self-intersections | — | Shapely validity check | + +--- + +## 4. The NX Hands: Assembly FEM Architecture + +### 4.1 The Two-Model Structure + +The FEA uses an Assembly FEM (AFEM) with two superposed finite element models. This is the standard aerospace approach for decoupling load introduction from structural detail, and it perfectly fits our problem where interfaces are fixed but internal topology varies. + +**Model A — Interface Model (built once, never modified):** + +Model A contains only the structural interface elements. It has no plate geometry — just the load introduction and boundary condition hardware. + +For each hole: +- One node at the hole center (the "master" node) +- N nodes equally spaced on the hole circumference (N = 12–24 depending on hole diameter, typically 1 node per ~2 mm of circumference) +- RBE2 or RBE3 spider element connecting center node to circumference nodes + - RBE2 (rigid): for bolted/pinned connections where the hole is rigidly constrained + - RBE3 (distributing): for bearing loads or distributed force introduction + +For plate boundary (edge BCs): +- Nodes distributed along the plate outer boundary at regular spacing (~2–5 mm) +- These serve as merge targets for Model B's edge mesh nodes + +All loads and boundary conditions are applied to Model A's nodes: +- Bolt forces → hole center nodes +- Fixed constraints → hole center nodes or edge nodes +- Bearing loads → hole center nodes (via RBE3) +- Enforced displacements → relevant center/edge nodes +- Pressures → applied directly to Model B elements (only exception) + +**Model B — Plate Model (rebuilt each iteration):** + +Model B is the 2D shell mesh of the current ribbed plate profile. It is the only model that changes during optimization. Key requirements: + +- Shell elements (CQUAD4/CTRIA3) with PSHELL property at plate thickness +- Mesh must have nodes at exact locations matching Model A's spider circumference nodes and edge BC nodes +- These "interface nodes" are enforced via mesh seed points (hard nodes) in NX's mesher +- The rest of the mesh fills in freely, adapting to the rib pattern geometry + +**Assembly FEM:** +- Superimposes Model A and Model B +- Merges coincident nodes at hole circumferences and plate edges (tolerance ~0.01 mm) +- After merge, loads flow from Model A's spiders through the merged nodes into Model B's shell mesh +- Solver settings, solution sequences, and output requests live at the assembly level + +### 4.2 One-Time Setup Procedure + +This is done once per plate project, before any optimization runs. + +**Step 1 — Extract geometry and build Model A:** + +```python +# NXOpen setup script — build the interface model +import NXOpen +import json +import math + +def build_interface_model(geometry_json_path, fem_part): + """ + Build Model A: spider elements at each hole + edge BC nodes. + """ + with open(geometry_json_path) as f: + geometry = json.load(f) + + interface_nodes = { + 'hole_centers': [], # (hole_idx, node_id, x, y) + 'hole_circumferences': [], # (hole_idx, node_id, x, y) + 'edge_nodes': [] # (node_id, x, y) + } + + # --- Create spider elements for each hole --- + for hole in geometry['holes']: + cx, cy = hole['center'] + radius = hole['diameter'] / 2.0 + + # Create center node + center_node = create_node(fem_part, cx, cy, 0.0) + interface_nodes['hole_centers'].append({ + 'hole_index': hole['index'], + 'node_id': center_node.Label, + 'x': cx, 'y': cy + }) + + # Create circumference nodes + n_circ = max(12, int(math.pi * hole['diameter'] / 2.0)) # ~1 node per 2mm + circ_nodes = [] + for j in range(n_circ): + angle = 2 * math.pi * j / n_circ + nx_ = cx + radius * math.cos(angle) + ny_ = cy + radius * math.sin(angle) + node = create_node(fem_part, nx_, ny_, 0.0) + circ_nodes.append(node) + interface_nodes['hole_circumferences'].append({ + 'hole_index': hole['index'], + 'node_id': node.Label, + 'x': nx_, 'y': ny_ + }) + + # Create RBE2 spider (rigid) — default for mounting holes + # Use RBE3 (distributing) for bearing load holes + spider_type = 'RBE2' # could be per-hole from hole table + create_spider(fem_part, center_node, circ_nodes, spider_type) + + # --- Create edge BC nodes --- + boundary = geometry['outer_boundary'] + edge_spacing = 3.0 # mm between edge nodes + for i in range(len(boundary)): + p1 = boundary[i] + p2 = boundary[(i + 1) % len(boundary)] + edge_length = distance(p1, p2) + n_nodes = max(2, int(edge_length / edge_spacing)) + for j in range(n_nodes): + t = j / n_nodes + ex = p1[0] + t * (p2[0] - p1[0]) + ey = p1[1] + t * (p2[1] - p1[1]) + node = create_node(fem_part, ex, ey, 0.0) + interface_nodes['edge_nodes'].append({ + 'node_id': node.Label, + 'x': ex, 'y': ey + }) + + # Save interface node map for the iteration script + with open('interface_nodes.json', 'w') as f: + json.dump(interface_nodes, f, indent=2) + + return interface_nodes +``` + +**Step 2 — Apply loads and BCs to Model A:** + +This is done manually in NX Simcenter (or scripted for standard load cases). The user applies all structural loads and boundary conditions to Model A's center nodes and edge nodes. These never change. + +Examples: +- Bolt M8 at hole 3: axial force 5000 N on center node of hole 3 +- Fixed constraint at holes 0, 1: fix all DOFs on center nodes of holes 0, 1 +- Simply supported edge: constrain Z-displacement on all edge nodes along one plate edge +- Bearing load at hole 7: 2000 N in X-direction via RBE3 at hole 7 center + +**Step 3 — Build dummy Model B and verify:** + +Create a simple Model B (e.g., the un-lightweighted plate meshed as shells), set up the assembly FEM with merging, and solve to verify the setup is correct. Compare results against a monolithic (non-assembly) FEA to confirm the spider elements and merging behave as expected. + +**Step 4 — Save the template:** + +Save Model A, the assembly FEM structure, and the interface node map. This is the reusable template for all optimization iterations. + +### 4.3 Per-Iteration NXOpen Journal Script + +This is the script that runs inside the optimization loop. It receives the ribbed plate profile from the Python brain and handles Model B rebuild + solve. + +```python +# NXOpen iteration script — rebuild Model B, merge, solve, extract +import NXOpen +import json +import numpy as np + +def iteration_solve(profile_path, interface_nodes_path, afem_part): + """ + Single optimization iteration: + 1. Delete old Model B geometry + mesh + 2. Import new profile + 3. Mesh with interface node seeds + 4. Merge nodes in AFEM + 5. Solve + 6. Extract results + """ + session = NXOpen.Session.GetSession() + + # Load inputs + with open(profile_path) as f: + profile = json.load(f) + with open(interface_nodes_path) as f: + interface_nodes = json.load(f) + + if not profile['valid']: + return {'status': 'invalid_geometry', 'mass': float('inf')} + + # --- Step 1: Clean old Model B --- + model_b = get_model_b_fem(afem_part) + delete_all_mesh(model_b) + delete_all_geometry(model_b) + + # --- Step 2: Import new 2D profile --- + # Create curves from profile coordinate arrays + outer_coords = profile['outer_boundary'] + create_closed_polyline(model_b, outer_coords) + + for pocket_coords in profile['pockets']: + create_closed_polyline(model_b, pocket_coords) + + # Create sheet body from bounded regions + sheet_body = create_sheet_from_curves(model_b) + + # --- Step 3: Mesh with interface seeds --- + mesh_control = get_mesh_collector(model_b) + + # Add hard-point mesh seeds at all interface locations + # These force the mesher to place nodes exactly where + # Model A's spiders attach + for node_info in interface_nodes['hole_circumferences']: + add_mesh_seed_point(mesh_control, node_info['x'], node_info['y']) + + for node_info in interface_nodes['edge_nodes']: + add_mesh_seed_point(mesh_control, node_info['x'], node_info['y']) + + # Set mesh parameters + set_element_size(mesh_control, target=2.0, min_size=0.5) # mm + set_element_type(mesh_control, 'CQUAD4') # prefer quads, allow tris + + # Generate mesh + mesh_control.GenerateMesh() + + # --- Step 4: Merge nodes in Assembly FEM --- + afem = get_assembly_fem(afem_part) + merge_coincident_nodes(afem, tolerance=0.05) # mm + + # Verify merge count matches expected + expected_merges = ( + len(interface_nodes['hole_circumferences']) + + len(interface_nodes['edge_nodes']) + ) + actual_merges = get_merge_count(afem) + if actual_merges < expected_merges * 0.95: + # Merge failed for some nodes — flag as warning + print(f"WARNING: Expected {expected_merges} merges, got {actual_merges}") + + # --- Step 5: Solve --- + solution = get_solution(afem, 'static_analysis') + solve_result = solution.Solve() + + if not solve_result.success: + return {'status': 'solve_failed', 'mass': float('inf')} + + # --- Step 6: Extract results --- + results = extract_results(afem, solution) + + return results + + +def extract_results(afem, solution): + """ + Extract field results from the solved assembly FEM. + Only extracts from Model B elements (the plate mesh), + ignoring Model A's spider elements. + """ + post = get_post_processor(afem) + + # Get results only from Model B's element group + model_b_elements = get_model_b_element_group(afem) + + # Von Mises stress (nodal, averaged at nodes) + stress_data = post.GetNodalResults( + solution, + result_type='Stress', + component='Von Mises', + element_group=model_b_elements + ) + + # Displacement magnitude (nodal) + disp_data = post.GetNodalResults( + solution, + result_type='Displacement', + component='Magnitude', + element_group=model_b_elements + ) + + # Strain (elemental) + strain_data = post.GetElementalResults( + solution, + result_type='Strain', + component='Von Mises', + element_group=model_b_elements + ) + + # Mass from Model B mesh (plate material only, not spiders) + mass = compute_shell_mass(model_b_elements) + + results = { + 'status': 'solved', + 'mass': mass, + 'max_von_mises': float(np.max(stress_data['values'])), + 'max_displacement': float(np.max(disp_data['values'])), + 'mean_von_mises': float(np.mean(stress_data['values'])), + 'stress_field': { + 'nodes_xy': stress_data['coordinates'].tolist(), + 'values': stress_data['values'].tolist() + }, + 'displacement_field': { + 'nodes_xy': disp_data['coordinates'].tolist(), + 'values': disp_data['values'].tolist() + }, + 'strain_field': { + 'elements_xy': strain_data['centroids'].tolist(), + 'values': strain_data['values'].tolist() + } + } + + return results +``` + +### 4.4 Why This Approach Is Robust + +The AFEM node-merge strategy solves the hardest problem in automated FEA iteration: **load and BC persistence across geometry changes.** Here's why it works reliably: + +**Fixed merge locations:** Hole centers, hole circumferences, and plate edges don't move between iterations. The merge is always at the same physical coordinates. NX's node merge by tolerance is a simple, reliable geometric operation. + +**Mesh seed enforcement:** By placing hard-point seeds at all interface locations, the mesher is forced to create nodes at exactly those coordinates. This guarantees that every merge finds its partner node. + +**Rib pattern agnostic:** Model A doesn't know or care what the rib pattern looks like. It only knows where the holes and edges are. Whether there are 50 or 200 pockets, the spiders connect the same way. + +**Easy validation:** After merging, a simple check (did we get the expected number of merged node pairs?) catches any meshing or geometry issues before wasting time on a solve. + +**Extensible:** Adding new load cases, new holes, or new BC types only requires modifying Model A (once). The optimization loop and Model B generation are unaffected. + +### 4.5 Simcenter Configuration Details + +**Shell property (PSHELL):** +- Thickness: from `geometry.json` (e.g., 10.0 mm) +- Material: MAT1 referencing the plate material (e.g., AL6061-T6) +- Applied to all Model B elements + +**Spider elements:** +- RBE2 (rigid): 6 DOF coupling from center to circumference nodes. Use for fixed/bolted connections where the hole acts as a rigid interface. +- RBE3 (weighted average): Distributes loads from center to circumference. Use for bearing loads where the hole deforms under load and you want realistic load distribution. +- Choice is per-hole, set during one-time setup based on connection type. + +**Mesh controls for Model B:** +- Target element size: 1.5–3.0 mm (captures rib geometry adequately) +- Minimum element size: 0.5 mm (allows refinement at narrow rib junctions) +- Element type: CQUAD4 dominant with CTRIA3 fill +- Mesh seeds: hard points at all interface node locations +- Edge mesh control on hole circumferences: N elements matching spider node count + +**Solution sequences:** +- SOL 101: Static analysis (v1) +- SOL 103: Normal modes / modal analysis (v2, for natural frequency constraints) +- SOL 105: Buckling (v2, for thin-rib stability checks) + +--- + +## 5. Atomizer Integration + +### 5.1 Parameter Space Definition + +```python +# Atomizer/Optuna parameter space +PARAM_SPACE = { + # Density field parameters + 'eta_0': {'type': 'float', 'low': 0.0, 'high': 0.4, 'desc': 'Baseline density offset'}, + 'alpha': {'type': 'float', 'low': 0.3, 'high': 2.0, 'desc': 'Hole influence scale'}, + 'R_0': {'type': 'float', 'low': 10.0, 'high': 100.0, 'desc': 'Base influence radius (mm)'}, + 'kappa': {'type': 'float', 'low': 0.0, 'high': 3.0, 'desc': 'Weight-to-radius coupling'}, + 'p': {'type': 'float', 'low': 1.0, 'high': 4.0, 'desc': 'Decay exponent'}, + 'beta': {'type': 'float', 'low': 0.0, 'high': 1.0, 'desc': 'Edge influence scale'}, + 'R_edge': {'type': 'float', 'low': 5.0, 'high': 40.0, 'desc': 'Edge influence radius (mm)'}, + + # Spacing parameters + 's_min': {'type': 'float', 'low': 8.0, 'high': 20.0, 'desc': 'Min cell size (mm)'}, + 's_max': {'type': 'float', 'low': 25.0, 'high': 60.0, 'desc': 'Max cell size (mm)'}, + + # Rib thickness parameters + 't_min': {'type': 'float', 'low': 2.0, 'high': 4.0, 'desc': 'Min rib thickness (mm)'}, + 't_0': {'type': 'float', 'low': 2.0, 'high': 6.0, 'desc': 'Nominal rib thickness (mm)'}, + 'gamma': {'type': 'float', 'low': 0.0, 'high': 3.0, 'desc': 'Density-thickness coupling'}, + + # Manufacturing / frame parameters + 'w_frame': {'type': 'float', 'low': 3.0, 'high': 20.0, 'desc': 'Perimeter frame width (mm)'}, + 'r_f': {'type': 'float', 'low': 0.5, 'high': 3.0, 'desc': 'Pocket fillet radius (mm)'}, + 'd_keep': {'type': 'float', 'low': 1.0, 'high': 3.0, 'desc': 'Hole keepout multiplier (× diameter)'}, +} +``` + +**Total: 15 continuous parameters.** This is a comfortable range for Optuna TPE. Can easily expand to 20-25 if needed (e.g., adding per-class influence overrides, smoothing length, separate edge/hole decay exponents). + +### 5.2 Objective Function + +```python +def objective(trial, geometry, sim_template): + # Sample parameters + params = {} + for name, spec in PARAM_SPACE.items(): + params[name] = trial.suggest_float(name, spec['low'], spec['high']) + + # --- Python Brain: generate geometry --- + profile_path = f'/tmp/isogrid_trial_{trial.number}.json' + valid, checks = generate_isogrid(geometry, params, profile_path) + + if not valid: + # Geometry failed validation — penalize + return float('inf') + + # --- NX Hands: mesh and solve --- + results = nx_import_and_solve(profile_path, sim_template) + + if results['status'] != 'solved': + return float('inf') + + # --- Evaluate --- + mass = results['mass'] + max_stress = results['max_von_mises'] + max_disp = results['max_displacement'] + + # Constraint penalties + penalty = 0.0 + SIGMA_ALLOW = 150.0 # MPa (example for AL6061-T6 with SF) + DELTA_MAX = 0.1 # mm (example) + + if max_stress > SIGMA_ALLOW: + penalty += 1e4 * ((max_stress / SIGMA_ALLOW) - 1.0) ** 2 + if max_disp > DELTA_MAX: + penalty += 1e4 * ((max_disp / DELTA_MAX) - 1.0) ** 2 + + # Store fields for visualization (not used by optimizer) + trial.set_user_attr('stress_field', results.get('stress_field')) + trial.set_user_attr('displacement_field', results.get('displacement_field')) + trial.set_user_attr('mass', mass) + trial.set_user_attr('max_stress', max_stress) + trial.set_user_attr('max_disp', max_disp) + trial.set_user_attr('num_pockets', checks.get('num_pockets', 0)) + + return mass + penalty +``` + +### 5.3 Convergence and Stopping + +With ~2 min/iteration and 15 parameters, expect: + +| Trials | Wall Time | Expected Outcome | +|---|---|---| +| 50 | ~1.5 hours | Random exploration, baseline understanding | +| 200 | ~7 hours | Surrogate learning, good solutions emerging | +| 500 | ~17 hours | Near-optimal, diminishing returns starting | +| 1000 | ~33 hours | Refined optimum, convergence likely | +| 2000 | ~67 hours | Exhaustive, marginal improvement | + +**Recommendation:** Start with 500 trials overnight. Review results. If the Pareto front is still moving, extend to 1000. + +--- + +## 6. V2 Roadmap: Stress-Feedback Enhancement + +Once v1 is running and producing good results, the stress-feedback enhancement adds the FEA stress/displacement fields as inputs to the density field: + +``` +η(x) = clamp(0, 1, η₀ + α·I(x) + β·E(x) + λ·S_prev(x)) +``` + +Where `S_prev(x)` is the normalized, smoothed stress field from the previous FEA of the same trial. This creates a local feedback loop within each Atomizer trial: + +1. Generate isogrid from density field (hole-based only, no stress data) +2. Run FEA → get stress field +3. Regenerate isogrid from updated density field (now including stress) +4. Run FEA → get updated stress field +5. Check convergence (stress field stable?) → if not, repeat from 3 +6. Report final metrics to Atomizer + +Atomizer then optimizes the meta-parameters including λ (stress feedback strength) and a smoothing kernel size for the stress field. + +**New parameters for v2:** λ (stress feedback weight), σ_smooth (stress field smoothing kernel size, mm), n_inner_max (max inner iterations). + +This is more expensive (2-5× FEA per trial) but produces designs that are truly structurally adapted, not just geometrically adapted. + +--- + +## 7. Implementation Sequence + +### Phase 1 — Python Brain Standalone (1-2 weeks) + +Build and test the geometry generator independently of NX: + +- Density field evaluation with exponential kernel +- Constrained Delaunay triangulation (using `triangle` library) +- Rib thickening and pocket profile generation (using `shapely`) +- Manufacturing constraint validation +- Matplotlib visualization (density field heatmap, rib pattern overlay, pocket profiles) +- Test with 3-5 different plate geometries and 50+ parameter sets + +**Deliverable:** A Python module that takes `geometry.json` + parameters → outputs `rib_profile.json` + visualization plots. + +### Phase 2 — NX Geometry Extraction + AFEM Setup Scripts (1-2 weeks) + +Build the NXOpen scripts for one-time project setup: + +- Face selection and hole detection script → exports `geometry.json` +- Hole weight assignment UI (or table import) +- Model A builder: spider elements at all holes + edge BC nodes → exports `interface_nodes.json` +- Assembly FEM creation with Model A + dummy Model B +- Load/BC application to Model A (manual or scripted for standard cases) +- Verification solve on dummy Model B +- Save as reusable template + +**Deliverable:** Complete one-time setup pipeline. Given any plate model, produces the AFEM template ready for optimization. + +### Phase 3 — NX Iteration Script (1-2 weeks) + +Build the NXOpen journal script for the per-iteration loop: + +- Model B geometry cleanup (delete old mesh + geometry) +- Profile import from `rib_profile.json` → NX curves → sheet body +- Mesh with hard-point seeds at interface node locations +- Assembly node merge with verification +- Nastran solve trigger +- Result extraction (stress/displacement/strain fields + scalar metrics) → `results.json` +- End-to-end test: Python Brain → NX journal → results → validate against manual FEA + +**Deliverable:** Complete single-iteration pipeline, verified against known-good manual analysis. + +### Phase 4 — Atomizer Integration (1 week) + +Wire Atomizer to orchestrate the pipeline: + +- Parameter sampling → Python Brain → NX journal trigger → result extraction +- Objective function with constraint penalties +- Study creation, execution, result logging +- Failed iteration handling (geometry validation failures, solve failures, merge warnings) +- Convergence monitoring (plot best mass vs. trial number) + +**Deliverable:** Full automated optimization loop, ready for production runs. + +### Phase 5 — Validation + First Real Project (1-2 weeks) + +Run on an actual client plate: + +- Full optimization campaign (500+ trials) +- Compare optimized mass vs. original solid plate and vs. uniform isogrid +- Manufacturing review (waterjet quote/feasibility from shop) +- Verify optimal design with refined mesh / higher-fidelity analysis +- Iterate on parameter bounds and manufacturing constraints based on feedback + +**Deliverable:** First optimized plate design, manufacturing-ready. + +--- + +## Appendix A: Python Dependencies + +``` +numpy >= 1.24 +scipy >= 1.10 +shapely >= 2.0 +triangle >= 20230923 # Python binding for Shewchuk's Triangle +matplotlib >= 3.7 +``` + +Optional for v2: `gmsh` (alternative mesher), `plotly` (interactive viz). + +## Appendix B: Key Reference Material + +- Shewchuk, J.R. "Triangle: A Two-Dimensional Quality Mesh Generator" — the engine behind the constrained Delaunay step +- NASA CR-124075 "Isogrid Design Handbook" — classical isogrid design equations +- Optuna documentation — TPE sampler configuration and multi-objective support +- NXOpen Python API Reference — for geometry creation and Simcenter automation diff --git a/tools/adaptive-isogrid/requirements.txt b/tools/adaptive-isogrid/requirements.txt new file mode 100644 index 00000000..3706e7da --- /dev/null +++ b/tools/adaptive-isogrid/requirements.txt @@ -0,0 +1,5 @@ +numpy>=1.24 +scipy>=1.10 +shapely>=2.0 +triangle>=20230923 +matplotlib>=3.7 diff --git a/tools/adaptive-isogrid/src/__init__.py b/tools/adaptive-isogrid/src/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/tools/adaptive-isogrid/src/atomizer_study.py b/tools/adaptive-isogrid/src/atomizer_study.py new file mode 100644 index 00000000..354e6395 --- /dev/null +++ b/tools/adaptive-isogrid/src/atomizer_study.py @@ -0,0 +1,48 @@ +""" +Atomizer/Optuna integration for adaptive isogrid optimization. + +Wires: parameter sampling → Python Brain → NX Hands → result extraction. +""" + +# Parameter space definition +PARAM_SPACE = { + # Density field + 'eta_0': {'type': 'float', 'low': 0.0, 'high': 0.4, 'desc': 'Baseline density offset'}, + 'alpha': {'type': 'float', 'low': 0.3, 'high': 2.0, 'desc': 'Hole influence scale'}, + 'R_0': {'type': 'float', 'low': 10.0, 'high': 100.0, 'desc': 'Base influence radius (mm)'}, + 'kappa': {'type': 'float', 'low': 0.0, 'high': 3.0, 'desc': 'Weight-to-radius coupling'}, + 'p': {'type': 'float', 'low': 1.0, 'high': 4.0, 'desc': 'Decay exponent'}, + 'beta': {'type': 'float', 'low': 0.0, 'high': 1.0, 'desc': 'Edge influence scale'}, + 'R_edge': {'type': 'float', 'low': 5.0, 'high': 40.0, 'desc': 'Edge influence radius (mm)'}, + # Spacing + 's_min': {'type': 'float', 'low': 8.0, 'high': 20.0, 'desc': 'Min cell size (mm)'}, + 's_max': {'type': 'float', 'low': 25.0, 'high': 60.0, 'desc': 'Max cell size (mm)'}, + # Rib thickness + 't_min': {'type': 'float', 'low': 2.0, 'high': 4.0, 'desc': 'Min rib thickness (mm)'}, + 't_0': {'type': 'float', 'low': 2.0, 'high': 6.0, 'desc': 'Nominal rib thickness (mm)'}, + 'gamma': {'type': 'float', 'low': 0.0, 'high': 3.0, 'desc': 'Density-thickness coupling'}, + # Manufacturing + 'w_frame': {'type': 'float', 'low': 3.0, 'high': 20.0, 'desc': 'Perimeter frame width (mm)'}, + 'r_f': {'type': 'float', 'low': 0.5, 'high': 3.0, 'desc': 'Pocket fillet radius (mm)'}, + 'd_keep': {'type': 'float', 'low': 1.0, 'high': 3.0, 'desc': 'Hole keepout multiplier (× diameter)'}, +} + +# Default parameters for standalone brain testing +DEFAULT_PARAMS = { + 'eta_0': 0.1, + 'alpha': 1.0, + 'R_0': 30.0, + 'kappa': 1.0, + 'p': 2.0, + 'beta': 0.3, + 'R_edge': 15.0, + 's_min': 12.0, + 's_max': 35.0, + 't_min': 2.5, + 't_0': 3.0, + 'gamma': 1.0, + 'w_frame': 8.0, + 'r_f': 1.5, + 'd_keep': 1.5, + 'min_pocket_radius': 1.5, +} diff --git a/tools/adaptive-isogrid/src/brain/__init__.py b/tools/adaptive-isogrid/src/brain/__init__.py new file mode 100644 index 00000000..6f0960da --- /dev/null +++ b/tools/adaptive-isogrid/src/brain/__init__.py @@ -0,0 +1,8 @@ +""" +Adaptive Isogrid — Python Brain + +Density field evaluation → Constrained Delaunay triangulation → +Rib thickness computation → Pocket profile generation → Validation. + +Takes geometry.json + parameters → outputs rib_profile.json +""" diff --git a/tools/adaptive-isogrid/src/brain/density_field.py b/tools/adaptive-isogrid/src/brain/density_field.py new file mode 100644 index 00000000..af71670e --- /dev/null +++ b/tools/adaptive-isogrid/src/brain/density_field.py @@ -0,0 +1,146 @@ +""" +Density field η(x) — maps every point on the plate to [0, 1]. + +η(x) = clamp(0, 1, η₀ + α·I(x) + β·E(x)) + +Where: + I(x) = Σᵢ wᵢ · exp(-(dᵢ(x)/Rᵢ)^p) — hole influence + E(x) = exp(-(d_edge(x)/R_edge)^p_edge) — edge reinforcement +""" + +import numpy as np +from shapely.geometry import Polygon, Point, LinearRing + + +def compute_hole_influence(x, y, holes, params): + """ + Compute hole influence term I(x) at point (x, y). + + I(x) = Σᵢ wᵢ · exp(-(dᵢ(x) / Rᵢ)^p) + + Parameters + ---------- + x, y : float + Query point coordinates. + holes : list[dict] + Hole definitions from geometry.json. + params : dict + Must contain: R_0, kappa, p + """ + R_0 = params['R_0'] + kappa = params['kappa'] + p = params['p'] + + influence = 0.0 + for hole in holes: + w = hole['weight'] + if w <= 0: + continue + + # Distance to hole + if hole.get('is_circular', True): + cx, cy = hole['center'] + d = np.sqrt((x - cx)**2 + (y - cy)**2) + else: + # Non-circular: distance to boundary polygon + ring = LinearRing(hole['boundary']) + d = Point(x, y).distance(ring) + + # Influence radius scales with hole weight + R_i = R_0 * (1.0 + kappa * w) + + # Exponential decay + influence += w * np.exp(-(d / R_i)**p) + + return influence + + +def compute_edge_influence(x, y, outer_boundary, params): + """ + Compute edge reinforcement term E(x) at point (x, y). + + E(x) = exp(-(d_edge(x) / R_edge)^p) + """ + R_edge = params['R_edge'] + p = params['p'] # reuse same decay exponent (could be separate in v2) + + ring = LinearRing(outer_boundary) + d_edge = Point(x, y).distance(ring) + + return np.exp(-(d_edge / R_edge)**p) + + +def evaluate_density(x, y, geometry, params): + """ + Evaluate the combined density field η(x, y). + + η(x) = clamp(0, 1, η₀ + α·I(x) + β·E(x)) + + Returns + ------- + float : density value in [0, 1] + """ + eta_0 = params['eta_0'] + alpha = params['alpha'] + beta = params['beta'] + + I = compute_hole_influence(x, y, geometry['holes'], params) + E = compute_edge_influence(x, y, geometry['outer_boundary'], params) + + eta = eta_0 + alpha * I + beta * E + return np.clip(eta, 0.0, 1.0) + + +def evaluate_density_grid(geometry, params, resolution=2.0): + """ + Evaluate density field on a regular grid covering the plate. + Useful for visualization (heatmap). + + Parameters + ---------- + geometry : dict + Plate geometry with outer_boundary and holes. + params : dict + Density field parameters. + resolution : float + Grid spacing in mm. + + Returns + ------- + X, Y, eta : ndarray + Meshgrid coordinates and density values. + """ + boundary = np.array(geometry['outer_boundary']) + x_min, y_min = boundary.min(axis=0) + x_max, y_max = boundary.max(axis=0) + + xs = np.arange(x_min, x_max + resolution, resolution) + ys = np.arange(y_min, y_max + resolution, resolution) + X, Y = np.meshgrid(xs, ys) + + plate_poly = Polygon(geometry['outer_boundary']) + eta = np.full_like(X, np.nan) + + for i in range(X.shape[0]): + for j in range(X.shape[1]): + pt = Point(X[i, j], Y[i, j]) + if plate_poly.contains(pt): + eta[i, j] = evaluate_density(X[i, j], Y[i, j], geometry, params) + + return X, Y, eta + + +def density_to_spacing(eta, params): + """Convert density η to local target triangle edge length s(x).""" + s_min = params['s_min'] + s_max = params['s_max'] + return s_max - (s_max - s_min) * eta + + +def density_to_rib_thickness(eta, params): + """Convert density η to local rib thickness t(x).""" + t_min = params['t_min'] + t_0 = params['t_0'] + gamma = params['gamma'] + t_max = t_0 * (1.0 + gamma) # max possible thickness + return np.clip(t_0 * (1.0 + gamma * eta), t_min, t_max) diff --git a/tools/adaptive-isogrid/src/brain/pocket_profiles.py b/tools/adaptive-isogrid/src/brain/pocket_profiles.py new file mode 100644 index 00000000..bce05ba3 --- /dev/null +++ b/tools/adaptive-isogrid/src/brain/pocket_profiles.py @@ -0,0 +1,176 @@ +""" +Pocket profile generation — convert triangulation into manufacturable pocket cutouts. + +Each triangle → inset by half-rib-thickness per edge → fillet corners → pocket polygon. +""" + +import numpy as np +from shapely.geometry import Polygon +from shapely.ops import unary_union + +from .density_field import evaluate_density, density_to_rib_thickness + + +def get_unique_edges(triangles): + """Extract unique edges from triangle connectivity.""" + edges = set() + for tri in triangles: + for i in range(3): + edge = tuple(sorted([tri[i], tri[(i + 1) % 3]])) + edges.add(edge) + return edges + + +def get_triangle_edge_thicknesses(tri_vertices, tri_edges_thickness): + """Get thickness for each edge of a triangle.""" + return [tri_edges_thickness.get(e, 2.0) for e in tri_vertices] + + +def inset_triangle(p0, p1, p2, d0, d1, d2): + """ + Inset a triangle by distances d0, d1, d2 on each edge. + + Edge 0: p0→p1, inset by d0 + Edge 1: p1→p2, inset by d1 + Edge 2: p2→p0, inset by d2 + + Returns inset triangle vertices or None if triangle collapses. + """ + def edge_normal_inward(a, b, opposite): + """Unit inward normal of edge a→b (toward opposite vertex).""" + dx, dy = b[0] - a[0], b[1] - a[1] + length = np.sqrt(dx**2 + dy**2) + if length < 1e-10: + return np.array([0.0, 0.0]) + # Two possible normals + n1 = np.array([-dy / length, dx / length]) + n2 = np.array([dy / length, -dx / length]) + # Pick the one pointing toward the opposite vertex + mid = (np.array(a) + np.array(b)) / 2.0 + if np.dot(n1, np.array(opposite) - mid) > 0: + return n1 + return n2 + + points = [np.array(p0), np.array(p1), np.array(p2)] + offsets = [d0, d1, d2] + + # Compute inset lines for each edge + inset_lines = [] + edges = [(0, 1, 2), (1, 2, 0), (2, 0, 1)] + for (i, j, k), d in zip(edges, offsets): + n = edge_normal_inward(points[i], points[j], points[k]) + # Offset edge inward + a_off = points[i] + n * d + b_off = points[j] + n * d + inset_lines.append((a_off, b_off)) + + # Intersect consecutive inset lines to get new vertices + new_verts = [] + for idx in range(3): + a1, b1 = inset_lines[idx] + a2, b2 = inset_lines[(idx + 1) % 3] + pt = line_intersection(a1, b1, a2, b2) + if pt is None: + return None + new_verts.append(pt) + + # Check if triangle is valid (positive area) + area = triangle_area(new_verts[0], new_verts[1], new_verts[2]) + if area < 0.1: # mm² — too small + return None + + return new_verts + + +def line_intersection(a1, b1, a2, b2): + """Find intersection of two lines (a1→b1) and (a2→b2).""" + d1 = b1 - a1 + d2 = b2 - a2 + cross = d1[0] * d2[1] - d1[1] * d2[0] + if abs(cross) < 1e-12: + return None # parallel + t = ((a2[0] - a1[0]) * d2[1] - (a2[1] - a1[1]) * d2[0]) / cross + return a1 + t * d1 + + +def triangle_area(p0, p1, p2): + """Signed area of triangle.""" + return 0.5 * abs( + (p1[0] - p0[0]) * (p2[1] - p0[1]) - + (p2[0] - p0[0]) * (p1[1] - p0[1]) + ) + + +def generate_pockets(triangulation, geometry, params): + """ + Generate pocket profiles from triangulation. + + Each triangle becomes a pocket (inset + filleted). + Small triangles are left solid (no pocket). + + Returns + ------- + list[dict] : pocket definitions with vertices and metadata. + """ + vertices = triangulation['vertices'] + triangles = triangulation['triangles'] + r_f = params['r_f'] + min_pocket_radius = params.get('min_pocket_radius', 1.5) + + # Precompute edge thicknesses + edge_thickness = {} + for tri in triangles: + for i in range(3): + edge = tuple(sorted([tri[i], tri[(i + 1) % 3]])) + if edge not in edge_thickness: + mid = (vertices[edge[0]] + vertices[edge[1]]) / 2.0 + eta = evaluate_density(mid[0], mid[1], geometry, params) + edge_thickness[edge] = density_to_rib_thickness(eta, params) + + pockets = [] + for tri_idx, tri in enumerate(triangles): + p0 = vertices[tri[0]] + p1 = vertices[tri[1]] + p2 = vertices[tri[2]] + + # Get edge thicknesses (half for inset) + e01 = tuple(sorted([tri[0], tri[1]])) + e12 = tuple(sorted([tri[1], tri[2]])) + e20 = tuple(sorted([tri[2], tri[0]])) + + d0 = edge_thickness[e01] / 2.0 + d1 = edge_thickness[e12] / 2.0 + d2 = edge_thickness[e20] / 2.0 + + inset = inset_triangle(p0, p1, p2, d0, d1, d2) + if inset is None: + continue + + # Check minimum pocket size via inscribed circle approximation + pocket_poly = Polygon(inset) + if not pocket_poly.is_valid or pocket_poly.is_empty: + continue + + # Approximate inscribed radius + inscribed_r = pocket_poly.area / (pocket_poly.length / 2.0) + if inscribed_r < min_pocket_radius: + continue + + # Apply fillet (buffer negative then positive = round inward corners) + fillet_amount = min(r_f, inscribed_r * 0.4) # don't over-fillet + if fillet_amount > 0.1: + filleted = pocket_poly.buffer(-fillet_amount).buffer(fillet_amount) + else: + filleted = pocket_poly + + if filleted.is_empty or not filleted.is_valid: + continue + + coords = list(filleted.exterior.coords) + pockets.append({ + 'triangle_index': tri_idx, + 'vertices': coords, + 'area': filleted.area, + }) + + return pockets diff --git a/tools/adaptive-isogrid/src/brain/profile_assembly.py b/tools/adaptive-isogrid/src/brain/profile_assembly.py new file mode 100644 index 00000000..61d1f573 --- /dev/null +++ b/tools/adaptive-isogrid/src/brain/profile_assembly.py @@ -0,0 +1,100 @@ +""" +Profile assembly — combine plate boundary, pockets, and holes +into the final 2D ribbed plate profile for NX import. + +Output: plate boundary - pockets - holes = ribbed plate (Shapely geometry) +""" + +from shapely.geometry import Polygon +from shapely.ops import unary_union + + +def assemble_profile(geometry, pockets, params): + """ + Create the final 2D ribbed plate profile. + + Plate boundary - pockets - holes = ribbed plate + + Parameters + ---------- + geometry : dict + Plate geometry with outer_boundary and holes. + pockets : list[dict] + Pocket definitions from pocket_profiles.generate_pockets(). + params : dict + Must contain w_frame (perimeter frame width). + + Returns + ------- + Shapely Polygon/MultiPolygon : the ribbed plate profile. + """ + plate = Polygon(geometry['outer_boundary']) + w_frame = params['w_frame'] + + # Inner boundary (frame inset) + if w_frame > 0: + inner_plate = plate.buffer(-w_frame) + else: + inner_plate = plate + + if inner_plate.is_empty: + # Frame too wide for plate — return solid plate minus holes + ribbed = plate + for hole in geometry['holes']: + ribbed = ribbed.difference(Polygon(hole['boundary'])) + return ribbed + + # Union all pocket polygons + pocket_polys = [] + for p in pockets: + poly = Polygon(p['vertices']) + if poly.is_valid and not poly.is_empty: + pocket_polys.append(poly) + + if pocket_polys: + all_pockets = unary_union(pocket_polys) + # Clip pockets to inner plate (don't cut into frame) + clipped_pockets = all_pockets.intersection(inner_plate) + # Subtract pockets from plate + ribbed = plate.difference(clipped_pockets) + else: + ribbed = plate + + # Subtract original hole boundaries + for hole in geometry['holes']: + hole_poly = Polygon(hole['boundary']) + ribbed = ribbed.difference(hole_poly) + + return ribbed + + +def profile_to_json(ribbed_plate, params): + """ + Convert Shapely ribbed plate geometry to JSON-serializable dict + for NXOpen import. + """ + if ribbed_plate.is_empty: + return {'valid': False, 'reason': 'empty_geometry'} + + # Separate exterior and interiors + if ribbed_plate.geom_type == 'MultiPolygon': + # Take the largest polygon (should be the plate) + largest = max(ribbed_plate.geoms, key=lambda g: g.area) + else: + largest = ribbed_plate + + # Classify interiors as pockets or holes + outer_coords = list(largest.exterior.coords) + pocket_coords = [] + hole_coords = [] + + for interior in largest.interiors: + coords = list(interior.coords) + pocket_coords.append(coords) + + return { + 'valid': True, + 'outer_boundary': outer_coords, + 'pockets': pocket_coords, + 'parameters_used': params, + } diff --git a/tools/adaptive-isogrid/src/brain/triangulation.py b/tools/adaptive-isogrid/src/brain/triangulation.py new file mode 100644 index 00000000..9fee6d9b --- /dev/null +++ b/tools/adaptive-isogrid/src/brain/triangulation.py @@ -0,0 +1,160 @@ +""" +Constrained Delaunay triangulation with density-adaptive refinement. + +Uses Shewchuk's Triangle library to generate adaptive mesh that +respects plate boundary, hole keepouts, and density-driven spacing. +""" + +import numpy as np +import triangle as tr +from shapely.geometry import Polygon, LinearRing + +from .density_field import evaluate_density, density_to_spacing + + +def offset_polygon(coords, distance, inward=True): + """ + Offset a polygon boundary by `distance`. + inward=True shrinks, inward=False expands. + + Uses Shapely buffer (negative = inward for exterior ring). + """ + poly = Polygon(coords) + if inward: + buffered = poly.buffer(-distance) + else: + buffered = poly.buffer(distance) + + if buffered.is_empty: + return coords # can't offset that much, return original + + if hasattr(buffered, 'exterior'): + return list(buffered.exterior.coords)[:-1] # remove closing duplicate + return coords + + +def build_pslg(geometry, params): + """ + Build Planar Straight Line Graph for Triangle library. + + Parameters + ---------- + geometry : dict + Plate geometry (outer_boundary, holes). + params : dict + Must contain d_keep (hole keepout multiplier). + + Returns + ------- + dict : Triangle-compatible PSLG with vertices, segments, holes. + """ + d_keep = params['d_keep'] + + vertices = [] + segments = [] + hole_markers = [] + + # Outer boundary (no offset needed — frame is handled in profile assembly) + outer = geometry['outer_boundary'] + v_start = len(vertices) + vertices.extend(outer) + n = len(outer) + for i in range(n): + segments.append([v_start + i, v_start + (i + 1) % n]) + + # Each hole with keepout offset + for hole in geometry['holes']: + keepout_dist = d_keep * (hole.get('diameter', 10.0) or 10.0) / 2.0 + hole_boundary = offset_polygon(hole['boundary'], keepout_dist, inward=False) + + v_start = len(vertices) + vertices.extend(hole_boundary) + n_h = len(hole_boundary) + for i in range(n_h): + segments.append([v_start + i, v_start + (i + 1) % n_h]) + + # Marker inside hole tells Triangle to leave it empty + hole_markers.append(hole['center']) + + result = { + 'vertices': np.array(vertices, dtype=np.float64), + 'segments': np.array(segments, dtype=np.int32), + } + if hole_markers: + result['holes'] = np.array(hole_markers, dtype=np.float64) + + return result + + +def compute_triangle_areas(vertices, triangles): + """Compute area of each triangle.""" + v0 = vertices[triangles[:, 0]] + v1 = vertices[triangles[:, 1]] + v2 = vertices[triangles[:, 2]] + # Cross product / 2 + areas = 0.5 * np.abs( + (v1[:, 0] - v0[:, 0]) * (v2[:, 1] - v0[:, 1]) - + (v2[:, 0] - v0[:, 0]) * (v1[:, 1] - v0[:, 1]) + ) + return areas + + +def compute_centroids(vertices, triangles): + """Compute centroid of each triangle.""" + v0 = vertices[triangles[:, 0]] + v1 = vertices[triangles[:, 1]] + v2 = vertices[triangles[:, 2]] + return (v0 + v1 + v2) / 3.0 + + +def generate_triangulation(geometry, params, max_refinement_passes=3): + """ + Generate density-adaptive constrained Delaunay triangulation. + + Parameters + ---------- + geometry : dict + Plate geometry. + params : dict + Full parameter set (density field + spacing + manufacturing). + max_refinement_passes : int + Number of iterative refinement passes. + + Returns + ------- + dict : Triangle result with 'vertices' and 'triangles'. + """ + pslg = build_pslg(geometry, params) + + # Initial triangulation with global max area + s_max = params['s_max'] + global_max_area = (np.sqrt(3) / 4.0) * s_max**2 + + # Triangle options: p=PSLG, q30=min angle 30°, a=area constraint, D=conforming Delaunay + result = tr.triangulate(pslg, f'pq30Da{global_max_area:.1f}') + + # Iterative refinement based on density field + for iteration in range(max_refinement_passes): + verts = result['vertices'] + tris = result['triangles'] + + areas = compute_triangle_areas(verts, tris) + centroids = compute_centroids(verts, tris) + + # Compute target area for each triangle based on density at centroid + target_areas = np.array([ + (np.sqrt(3) / 4.0) * density_to_spacing( + evaluate_density(cx, cy, geometry, params), params + )**2 + for cx, cy in centroids + ]) + + # Check if all triangles satisfy constraints (20% tolerance) + if np.all(areas <= target_areas * 1.2): + break + + # Set per-triangle max area and refine + result['triangle_max_area'] = target_areas + result = tr.triangulate(result, 'rpq30D') + + return result diff --git a/tools/adaptive-isogrid/src/brain/validation.py b/tools/adaptive-isogrid/src/brain/validation.py new file mode 100644 index 00000000..ed1f4d0e --- /dev/null +++ b/tools/adaptive-isogrid/src/brain/validation.py @@ -0,0 +1,106 @@ +""" +Manufacturing constraint validation for ribbed plate profiles. + +Checks run after profile assembly to catch issues before NX import. +""" + +import numpy as np +from shapely.geometry import Polygon +from shapely.validation import make_valid + + +def check_minimum_web(ribbed_plate, t_min): + """ + Check that no rib section is thinner than t_min. + + Uses negative buffer: if buffering inward by t_min/2 leaves + the geometry connected, minimum web width is satisfied. + + Returns True if minimum web width is OK. + """ + eroded = ribbed_plate.buffer(-t_min / 2.0) + if eroded.is_empty: + return False + # Check connectivity — should still be one piece + if eroded.geom_type == 'MultiPolygon': + # Some thin sections broke off — might be OK if they're tiny + main_area = max(g.area for g in eroded.geoms) + total_area = sum(g.area for g in eroded.geoms) + # If >95% of area is in one piece, it's fine + return main_area / total_area > 0.95 + return True + + +def check_no_islands(ribbed_plate): + """ + Check that the ribbed plate has no floating islands + (disconnected solid regions). + + Returns True if no islands found. + """ + if ribbed_plate.geom_type == 'Polygon': + return True + elif ribbed_plate.geom_type == 'MultiPolygon': + # Multiple polygons = islands. Only OK if one is dominant. + areas = [g.area for g in ribbed_plate.geoms] + return max(areas) / sum(areas) > 0.99 + return False + + +def estimate_mass(ribbed_plate, params, density_kg_m3=2700.0): + """ + Estimate mass of the ribbed plate. + + Parameters + ---------- + ribbed_plate : Shapely geometry + The ribbed plate profile. + params : dict + Must contain thickness info (from geometry). + density_kg_m3 : float + Material density (default: aluminum 6061). + + Returns + ------- + float : estimated mass in grams. + """ + area_mm2 = ribbed_plate.area + thickness_mm = params.get('thickness', 10.0) + volume_mm3 = area_mm2 * thickness_mm + volume_m3 = volume_mm3 * 1e-9 + mass_kg = volume_m3 * density_kg_m3 + return mass_kg * 1000.0 # grams + + +def validate_profile(ribbed_plate, params): + """ + Run all validation checks on the ribbed plate profile. + + Returns + ------- + tuple : (is_valid: bool, checks: dict) + """ + t_min = params['t_min'] + + # Fix any geometry issues + if not ribbed_plate.is_valid: + ribbed_plate = make_valid(ribbed_plate) + + checks = { + 'is_valid_geometry': ribbed_plate.is_valid, + 'min_web_width': check_minimum_web(ribbed_plate, t_min), + 'no_islands': check_no_islands(ribbed_plate), + 'no_self_intersections': ribbed_plate.is_valid, + 'mass_estimate_g': estimate_mass(ribbed_plate, params), + 'area_mm2': ribbed_plate.area, + 'num_interiors': len(list(ribbed_plate.interiors)) if ribbed_plate.geom_type == 'Polygon' else 0, + } + + is_valid = all([ + checks['is_valid_geometry'], + checks['min_web_width'], + checks['no_islands'], + checks['no_self_intersections'], + ]) + + return is_valid, checks diff --git a/tools/adaptive-isogrid/src/nx/__init__.py b/tools/adaptive-isogrid/src/nx/__init__.py new file mode 100644 index 00000000..b1e66c31 --- /dev/null +++ b/tools/adaptive-isogrid/src/nx/__init__.py @@ -0,0 +1,6 @@ +""" +Adaptive Isogrid — NX Hands + +NXOpen journal scripts for geometry extraction, AFEM setup, and per-iteration solve. +These scripts run inside NX Simcenter via the NXOpen Python API. +""" diff --git a/tools/adaptive-isogrid/src/nx/build_interface_model.py b/tools/adaptive-isogrid/src/nx/build_interface_model.py new file mode 100644 index 00000000..2fdea4af --- /dev/null +++ b/tools/adaptive-isogrid/src/nx/build_interface_model.py @@ -0,0 +1,28 @@ +""" +NXOpen script — Build Interface Model (Model A) for Assembly FEM. + +ONE-TIME: Creates spider elements (RBE2/RBE3) at each hole + edge BC nodes. +Exports interface_nodes.json for the iteration script. + +NOTE: Skeleton — requires NX environment for development. +See docs/technical-spec.md Section 4.2 for full pseudocode. +""" + + +def build_interface_model(geometry_json_path, fem_part): + """ + Build Model A: spider elements at each hole + edge BC nodes. + + For each hole: + - Center node (master) + - N circumference nodes (~1 per 2mm of circumference) + - RBE2 (rigid) or RBE3 (distributing) spider + + For plate boundary: + - Edge nodes at ~3mm spacing + + Exports interface_nodes.json with all node IDs and coordinates. + """ + raise NotImplementedError( + "Develop inside NX Simcenter. See docs/technical-spec.md Section 4.2." + ) diff --git a/tools/adaptive-isogrid/src/nx/extract_geometry.py b/tools/adaptive-isogrid/src/nx/extract_geometry.py new file mode 100644 index 00000000..ae1bdd36 --- /dev/null +++ b/tools/adaptive-isogrid/src/nx/extract_geometry.py @@ -0,0 +1,54 @@ +""" +NXOpen script — Extract plate geometry from selected face. + +ONE-TIME: Run inside NX. User selects plate face, assigns hole weights. +Exports geometry.json for the Python brain. + +NOTE: This is pseudocode / skeleton. Actual NXOpen API calls need +NX environment to develop and test. +""" + +# This script runs inside NX — NXOpen is available at runtime +# import NXOpen +# import json +# import math + + +def extract_plate_geometry(face, hole_weights): + """ + Extract plate geometry from an NX face. + + Parameters + ---------- + face : NXOpen.Face + The selected plate face. + hole_weights : dict + {loop_index: weight} from user input. + + Returns + ------- + dict : geometry definition for export. + """ + raise NotImplementedError( + "This script must be developed and tested inside NX Simcenter. " + "See docs/technical-spec.md Section 2 for full pseudocode." + ) + + +def sample_edge(edge, tolerance=0.1): + """Sample edge curve as polyline with given chord tolerance.""" + # NXOpen: edge.GetCurve(), evaluate at intervals + raise NotImplementedError + + +def fit_circle(points): + """Fit a circle to boundary points. Returns (center, diameter).""" + # Least-squares circle fit + raise NotImplementedError + + +def export_geometry(geometry, filepath='geometry.json'): + """Export geometry dict to JSON.""" + import json + with open(filepath, 'w') as f: + json.dump(geometry, f, indent=2) diff --git a/tools/adaptive-isogrid/src/nx/iteration_solve.py b/tools/adaptive-isogrid/src/nx/iteration_solve.py new file mode 100644 index 00000000..c6dc64bf --- /dev/null +++ b/tools/adaptive-isogrid/src/nx/iteration_solve.py @@ -0,0 +1,35 @@ +""" +NXOpen script — Per-iteration Model B rebuild + solve + extract. + +LOOP: Called by Atomizer for each trial. +1. Delete old Model B geometry + mesh +2. Import new 2D ribbed profile from rib_profile.json +3. Mesh with hard-point seeds at interface node locations +4. Merge nodes in Assembly FEM +5. Solve (Nastran) +6. Extract results → results.json + +NOTE: Skeleton — requires NX environment for development. +See docs/technical-spec.md Section 4.3 for full pseudocode. +""" + + +def iteration_solve(profile_path, interface_nodes_path, afem_part): + """ + Single optimization iteration. + + Returns dict with status, mass, stress/displacement fields. + """ + raise NotImplementedError( + "Develop inside NX Simcenter. See docs/technical-spec.md Section 4.3." + ) + + +def extract_results(afem, solution): + """ + Extract field results from solved assembly FEM. + Only from Model B elements (plate mesh), ignoring spiders. + """ + raise NotImplementedError( + "Develop inside NX Simcenter. See docs/technical-spec.md Section 4.3." + ) diff --git a/tools/adaptive-isogrid/tests/test_geometries/sample_bracket.json b/tools/adaptive-isogrid/tests/test_geometries/sample_bracket.json new file mode 100644 index 00000000..4e3d2074 --- /dev/null +++ b/tools/adaptive-isogrid/tests/test_geometries/sample_bracket.json @@ -0,0 +1,49 @@ +{ + "plate_id": "sample_bracket", + "units": "mm", + "thickness": 10.0, + "material": "AL6061-T6", + "outer_boundary": [[0,0], [400,0], [400,300], [0,300]], + "holes": [ + { + "index": 0, + "center": [50, 50], + "diameter": 12.0, + "is_circular": true, + "boundary": [[44,50], [44.2,51.9], [44.8,53.7], [45.8,55.2], [47.1,56.3], [48.6,57.0], [50.2,57.0], [51.7,56.5], [53.0,55.4], [53.9,53.9], [54.4,52.1], [54.4,50.1], [54.0,48.3], [53.1,46.7], [51.8,45.6], [50.3,45.0], [48.6,45.0], [47.1,45.5], [45.8,46.6], [44.9,48.1]], + "weight": 1.0 + }, + { + "index": 1, + "center": [50, 250], + "diameter": 12.0, + "is_circular": true, + "boundary": [[44,250], [44.2,251.9], [44.8,253.7], [45.8,255.2], [47.1,256.3], [48.6,257.0], [50.2,257.0], [51.7,256.5], [53.0,255.4], [53.9,253.9], [54.4,252.1], [54.4,250.1], [54.0,248.3], [53.1,246.7], [51.8,245.6], [50.3,245.0], [48.6,245.0], [47.1,245.5], [45.8,246.6], [44.9,248.1]], + "weight": 1.0 + }, + { + "index": 2, + "center": [200, 150], + "diameter": 8.0, + "is_circular": true, + "boundary": [[196,150], [196.1,151.3], [196.5,152.5], [197.2,153.5], [198.1,154.2], [199.2,154.6], [200.4,154.6], [201.5,154.2], [202.4,153.4], [203.0,152.3], [203.3,151.0], [203.3,149.7], [203.0,148.5], [202.3,147.4], [201.4,146.6], [200.3,146.2], [199.1,146.2], [198.0,146.7], [197.2,147.5], [196.5,148.7]], + "weight": 0.3 + }, + { + "index": 3, + "center": [350, 50], + "diameter": 10.0, + "is_circular": true, + "boundary": [[345,50], [345.1,51.6], [345.6,53.1], [346.5,54.3], [347.6,55.2], [349.0,55.7], [350.5,55.7], [351.9,55.1], [353.0,54.1], [353.8,52.8], [354.2,51.3], [354.2,49.7], [353.8,48.2], [352.9,46.9], [351.7,46.0], [350.3,45.5], [348.8,45.5], [347.5,46.1], [346.5,47.1], [345.6,48.4]], + "weight": 0.7 + }, + { + "index": 4, + "center": [350, 250], + "diameter": 10.0, + "is_circular": true, + "boundary": [[345,250], [345.1,251.6], [345.6,253.1], [346.5,254.3], [347.6,255.2], [349.0,255.7], [350.5,255.7], [351.9,255.1], [353.0,254.1], [353.8,252.8], [354.2,251.3], [354.2,249.7], [353.8,248.2], [352.9,246.9], [351.7,246.0], [350.3,245.5], [348.8,245.5], [347.5,246.1], [346.5,247.1], [345.6,248.4]], + "weight": 0.7 + } + ] +}