Commit Graph

1117 Commits

Author SHA1 Message Date
Jukka Aho c615012228 feat(demos): Add mock GPU and MPI demonstration (early prototype)
New 327-line mock demonstration (prototype before real hardware version):

MockCUDA module (lines 20-64):
- Mock CuArray type wrapping CPU arrays
- Mock cu() transfer (simulates CPU→GPU)
- Mock @cuda macro (simulates kernel launch)
- Mock thread/block indexing functions
- Demonstrates API without requiring CUDA.jl dependency

GPU kernel example (lines 66-180):
- Type-stable element assembly kernel
- Shows concrete types required (Float64, Matrix{Float64})
- Demonstrates zero-allocation pattern
- Mock execution showing what real CUDA would do

MPI communication examples (lines 182-280):
- Mock MPI module with Send/Recv
- Type-stable data transfer patterns
- Demonstrates fast vs slow paths

Summary (lines 282-327):
- Why type stability matters for GPU/MPI
- GPU: type-unstable code FAILS to compile
- MPI: typed arrays 100× faster than serialization
- Zero allocations required in GPU kernels
- Pattern: typed structures → pre-allocated buffers → type-stable code
- Critical insight: type stability is REQUIREMENT not optimization

Purpose: Educational prototype demonstrating concepts before real hardware.
Superseded by: gpu_mpi_demo.jl (real CUDA and MPI)
2025-11-09 10:51:06 +02:00
Jukka Aho 4daa429760 feat(demos): Add multi-GPU MPI Krylov solver demonstration
New 400-line distributed FEM solver demonstration with 6 parts:

Part 1: Generate test problem (lines 61-100)
- 10×10 SPD system, condition number ~3.45
- Distributed nodal assembly: each rank owns nodes
- Exact solution x=[1,2,...,10], RHS b=A*x

Part 2: Nodal assembly pattern (lines 101-140)
- get_row(i) and get_rhs(i) abstractions
- Row-by-row matrix construction
- Each rank assembles its local rows

Part 3: GPU transfer (lines 141-167)
- Transfer local data to GPU if CUDA available
- Falls back to CPU arrays if no GPU
- Reports bytes transferred per rank

Part 4: Distributed matrix-vector product (lines 168-203)
- matvec_distributed! function
- Each rank computes y_local = A_local * x_global
- GPU acceleration if available, CPU fallback

Part 5: Conjugate Gradient solver (lines 204-318)
- cg_distributed() with MPI collectives
- Allreduce for global dot products
- Allgatherv for vector assembly
- Reports convergence progress per iteration

Part 6: Verification (lines 319-345)
- Compare computed vs exact solution
- Report relative error
- Pass/fail verification (threshold 1e-6)

Summary (lines 346-400):
- Reports what was demonstrated on real hardware
- Nodal assembly, distributed computing, multi-GPU, Krylov CG
- Key insight: type-stable + nodal → scalable
- Relevance to JuliaFEM contact mechanics

Results: Converges in 9 iterations, 7.73×10⁻¹⁴ relative error
Hardware: 2 MPI ranks, NVIDIA RTX A2000 12GB per rank
Run: mpiexec -np 2 julia --project=. demos/krylov_mpi_gpu_demo.jl
2025-11-09 10:49:58 +02:00
Jukka Aho d39a5cd22f feat(demos): Add GPU and MPI real hardware demonstration script
New 322-line demonstration script proving type-stable data flows to GPU and MPI:

Part 1: Type-stable data structures (lines 54-77)
- Creates nodes, connectivity, displacement as typed arrays
- Material properties E, ν as Float64
- All structures explicitly typed (Matrix{Float64}, not Dict)

Part 2: MPI communication (lines 79-123)
- Rank 0 sends 24KB displacement data to rank 1
- Transfers material properties
- Validates data integrity with checksum
- Uses MPI.Send/Recv with typed buffers

Part 3: GPU kernel execution (lines 125-192)
- Defines assemble_element_kernel! for CUDA
- Type-stable kernel: Float64, Int32, no allocations
- Transfers data to GPU (CuArray)
- Launches kernel with thread blocks
- Validates results against expected values

Part 4: Combined GPU+MPI workflow (lines 194-271)
- Rank 0 computes on GPU
- Transfers results via MPI to rank 1
- End-to-end validation

Summary section (lines 273-322):
- Reports hardware used (GPU model, MPI ranks)
- Key insights: type stability required for GPU, enables fast MPI
- Conclusion: type-stable fields are foundation for modern FEM

Requirements: MPI (required), CUDA (optional, detects and uses if available)
Run: mpiexec -np 2 julia --project=. demos/gpu_mpi_demo.jl
2025-11-09 10:48:54 +02:00
Jukka Aho 00771c62d0 docs(demos): Add comprehensive Krylov solver demonstration guide
New 307-line comprehensive guide documenting:
- Overview: type-stable nodal assembly enables distributed solving
- Four key demonstrations: nodal assembly, distributed computing, multi-GPU, Krylov CG
- Running instructions for 2 or 4 MPI processes
- Expected output with all 6 parts (problem generation through verification)
- Technical details: 10×10 SPD system, partitioning, distributed matvec, CG algorithm
- GPU execution: CPU↔GPU transfer, type stability requirement
- MPI communication: Allreduce and Allgatherv patterns
- Performance characteristics: communication cost, computation cost, scaling analysis
- Relevance to JuliaFEM: why nodal assembly, type stability, matrix-free, distributed solving matter
- v0.5.1 vs v1.0 comparison and path forward
- Key insights tables: type stability enables everything, nodal assembly advantages, Krylov vs direct
- Validation results: 9 iterations, 7.73×10⁻¹⁴ error on real hardware
- References: CG method, domain decomposition, GPU computing, MPI
- Conclusion: 5 validated achievements proving the path forward
2025-11-09 10:48:11 +02:00
Jukka Aho 6d583eb30c docs(demos): Add GPU and MPI demonstration guide
New 77-line guide documenting:
- Prerequisites (MPI and CUDA globally installed)
- Running commands for MPI communication test
- Running commands for combined GPU+MPI test
- Single-process GPU test instructions
- What gets demonstrated (type stability requirement, MPI fast transfer, real hardware)
- Success indicators and result interpretation
- Key insight: same patterns enable CPU speedup, GPU execution, and MPI efficiency
2025-11-09 10:47:32 +02:00
Jukka Aho ebf823b5c6 docs(demos): Add README for technology demonstrations directory
New 125-line README documenting:
- Two main demonstrations (GPU+MPI and Krylov solver)
- Requirements (Julia 1.9+, MPI, optional CUDA)
- Key insights: type stability required for GPU/MPI/Krylov
- Nodal assembly pattern explanation
- Architecture validation (v0.5.1 vs v1.0 comparison)
- References to benchmarks and design docs
- Contributing guidelines for new demos
2025-11-09 10:36:51 +02:00
Jukka Aho 0cfe966063 docs(book): Add ADR-002 for topology without hardcoded node counts
Architecture Decision Record documenting topology/basis separation (292 lines):

- Explains decision to remove node counts from topology types
- Documents topology = pure geometry, basis determines node count
- Shows old Code Aster anti-pattern (TRIA3, TRIA6, QUAD4, QUAD8)
- Describes new design: Triangle + Lagrange{Triangle, P}
- Rationale: mathematical correctness, separation of concerns
- Enables edge/face DOFs (Nédélec, Raviart-Thomas)
- Eliminates combinatorial explosion (8 topologies vs hundreds)
- Consequences: extensible, correct, but breaking change
- Implementation strategy and migration plan
- Includes proper YAML front matter for book chapter
2025-11-09 09:33:41 +02:00
Jukka Aho d8fc224c42 docs(book): Add ADR-002 for topology without hardcoded node counts
New Architecture Decision Record for the comprehensive book (292 lines):

- Documents decision to remove node counts from topology types
- Explains topology = pure geometry, basis determines node count
- Shows old Code Aster anti-pattern (TRIA3, TRIA6, QUAD4, QUAD8)
- Describes new design: Triangle + Lagrange{Triangle, P}
- Rationale: mathematical correctness, separation of concerns
- Enables edge/face DOFs (Nédélec, Raviart-Thomas)
- Eliminates combinatorial explosion (8 topologies vs hundreds)
- Consequences: extensible, correct, but breaking change
- Implementation strategy and migration plan
- Backwards compatibility via aliases and shims
2025-11-09 09:31:45 +02:00
Jukka Aho 5e210187d8 docs(architecture): Complete rewrite explaining topology vs basis separation
Major documentation update (380 additions, 167 deletions):

- Explained topology is pure geometry (NO hardcoded node counts)
- Clarified basis determines BOTH polynomial degree AND node count
- Distinguished node count (connectivity) vs DOF count (unknowns)
- Added examples: Nedelec (edge DOFs), Raviart-Thomas (face DOFs)
- Documented Lagrange{Topology, P} parametric architecture
- Showed why Tri3/Quad4/Tet10 names are anti-pattern
- Updated all code examples to use new architecture
- Explained Serendipity vs full Lagrange tensor products
- Added performance implications and trade-offs
- Showed how one Triangle topology works for P1/P2/P3/Nedelec/etc
2025-11-09 09:29:58 +02:00
Jukka Aho 3cf39bb14d refactor(assembly): Comment out Tet10 specialization, fix formatting
- Commented out assemble_mass_matrix! specialization for Element{Tet10}
- Tet10 is now a topology type, not a basis type (name conflict)
- Needs refactoring to use Tet10Basis or new parametric architecture
- Fixed code formatting (spacing around operators, indentation)
- Added TODO comment explaining the issue
2025-11-09 09:29:30 +02:00
Jukka Aho 05febfc938 refactor(exports): Update topology exports for new architecture
- Removed nnodes from topology exports (now in basis module)
- Added exports for new topology names (Triangle, Quadrilateral, etc.)
- Kept old names as exports (they're aliases for backwards compatibility)
- Added explanatory comments about topology vs basis separation
- Added examples showing how node count comes from basis now
- Organized exports by dimension (0D/1D/2D/3D)
2025-11-09 09:29:19 +02:00
Jukka Aho abdbcb37cc refactor(elements): Comment out old topology_to_basis shims
- Commented out topology_to_basis() helper and old Element constructors
- These mapped old names (Tri3→Tri3Basis) which no longer exist
- Need to update for new Lagrange{Topology, P} parametric architecture
- Added TODO comment explaining migration needed
- Temporary measure until new Element constructors are implemented
2025-11-09 09:29:06 +02:00
Jukka Aho 9f42575cf2 feat(basis): Add Lagrange{Topology, P} parametric basis type
- Added Lagrange{T<:AbstractTopology, P} struct for parametric basis
- Implemented nnodes() formulas for all 7 topologies:
  - Segment: P+1 nodes
  - Triangle: (P+1)(P+2)/2 nodes (simplex formula)
  - Quadrilateral: (P+1)² nodes (tensor product)
  - Tetrahedron: (P+1)(P+2)(P+3)/6 nodes (simplex formula)
  - Hexahedron: (P+1)³ nodes (tensor product)
  - Pyramid: hardcoded for P=1,2,3 (no simple formula)
  - Wedge: (P+1)²(P+2)/2 nodes (triangle × segment)
- Added comprehensive documentation with examples
- Exported Lagrange and nnodes
- Node count now comes from basis, not topology
2025-11-09 09:28:45 +02:00
Jukka Aho dcc7f69a75 refactor(topology): Rename Wedge6 to Wedge, remove hardcoded node count
- Changed struct name from Wedge6 to Wedge
- Removed nnodes() method (node count now determined by basis)
- Updated all function signatures to use Wedge
- Added Wedge6 as backwards compatibility alias
- Added note explaining basis determines node count (P1=6, P2=15 nodes)
2025-11-09 09:28:29 +02:00
Jukka Aho ca7c8a4c75 refactor(topology): Rename Pyr5 to Pyramid, remove hardcoded node count
- Changed struct name from Pyr5 to Pyramid
- Removed nnodes() method (node count now determined by basis)
- Updated all function signatures to use Pyramid
- Added Pyr5 as backwards compatibility alias
- Added note explaining basis determines node count (P1=5, P2=13, P3=29)
2025-11-09 09:28:17 +02:00
Jukka Aho 835aac9962 refactor(topology): Rename Hex8 to Hexahedron, remove hardcoded node count
- Changed struct name from Hex8 to Hexahedron
- Removed nnodes() method (node count now determined by basis)
- Updated all function signatures to use Hexahedron
- Added Hex8 as backwards compatibility alias
- Added note explaining basis determines node count (Q1=8, Q2=27 nodes)
2025-11-09 09:28:07 +02:00
Jukka Aho 260ea0b170 refactor(topology): Rename Tet4 to Tetrahedron, remove hardcoded node count
- Changed struct name from Tet4 to Tetrahedron
- Removed nnodes() method (node count now determined by basis)
- Updated all function signatures to use Tetrahedron
- Added Tet4 as backwards compatibility alias
- Added note explaining basis determines node count (P1=4, P2=10 nodes)
2025-11-09 09:27:54 +02:00
Jukka Aho 7d5966c3d3 refactor(topology): Rename Seg2 to Segment, remove hardcoded node count
- Changed struct name from Seg2 to Segment
- Removed nnodes() method (node count now determined by basis)
- Updated all function signatures to use Segment instead of Seg2
- Added Seg2 as backwards compatibility alias
- Added note explaining basis determines node count
2025-11-09 09:27:43 +02:00
Jukka Aho 82758bd7c6 refactor(topology): Rename Quad4 to Quadrilateral, remove hardcoded node count
- Changed struct name from Quad4 to Quadrilateral
- Removed nnodes() method (node count now determined by basis)
- Updated documentation to explain topology vs basis separation
- Added examples showing Lagrange and Serendipity differences
- Added Quad4 as deprecated alias for backwards compatibility
- Clarified that topology defines geometry only, basis determines nodes
2025-11-09 09:27:31 +02:00
Jukka Aho 575136ba24 refactor(topology): Rename Tri3 to Triangle, remove hardcoded node count
- Changed struct name from Tri3 to Triangle
- Removed nnodes() method (node count now determined by basis)
- Updated documentation to explain topology vs basis separation
- Added examples showing Lagrange{Triangle, P} for different degrees
- Added Tri3 as deprecated alias for backwards compatibility
- Clarified that topology defines geometry only, basis determines nodes
2025-11-09 09:27:12 +02:00
Jukka Aho 27ff4b19f0 refactor(basis): Remove individual lagrange basis files
Deleted 7 files:
- src/basis/lagrange_segments.jl (Seg2, Seg3)
- src/basis/lagrange_triangles.jl (Tri3, Tri6)
- src/basis/lagrange_quadrangles.jl (Quad4, Quad8, Quad9)
- src/basis/lagrange_tetrahedrons.jl (Tet4, Tet10)
- src/basis/lagrange_hexahedrons.jl (Hex8, Hex20, Hex27)
- src/basis/lagrange_pyramids.jl (Pyr5)
- src/basis/lagrange_wedges.jl (Wedge6, Wedge15)

Reason: All 15 basis types consolidated into src/basis/lagrange_generated.jl
Generated by: julia --project=. src/basis/lagrange_generator.jl
2025-11-09 08:30:35 +02:00
Jukka Aho ccfec6eea7 refactor(scripts): Remove standalone generation script
Deleted: scripts/generate_lagrange_basis.jl (720 lines)

Reason: Functionality merged into src/basis/lagrange_generator.jl
The generator is now both a library (for inclusion) and a script (for execution).

Run as: julia --project=. src/basis/lagrange_generator.jl
2025-11-09 08:30:17 +02:00
Jukka Aho ca34e8da9d chore: Add auto-generated Lagrange basis functions
New file: src/basis/lagrange_generated.jl (448 lines, machine-generated)

Generated by: julia --project=. src/basis/lagrange_generator.jl
Generated at: 2025-11-09 07:28:58

Contains basis functions for 15 element types:
- 1D: Seg2Basis, Seg3Basis
- 2D triangles: Tri3Basis, Tri6Basis
- 2D quads: Quad4Basis, Quad8Basis, Quad9Basis
- 3D tets: Tet4Basis, Tet10Basis
- 3D hexes: Hex8Basis, Hex20Basis, Hex27Basis
- 3D pyramid: Pyr5Basis
- 3D wedges: Wedge6Basis, Wedge15Basis

All types have "Basis" suffix to avoid conflicts with topology types.

DO NOT EDIT MANUALLY - regenerate with generator script.
2025-11-09 08:29:50 +02:00
Jukka Aho 724eb9cffe refactor(basis): Merge generation script into lagrange_generator.jl
Consolidates scripts/generate_lagrange_basis.jl into src/basis/lagrange_generator.jl

Changes:
- Added Vecish type alias handling for standalone/included execution
- Added vandermonde_matrix() function (~40 lines) for polynomial basis construction
- Added ElementDescription struct with keyword constructor for readability
- Added 15 element definitions with reference coordinates and polynomial ansatz:
  * 1D: Seg2, Seg3
  * 2D triangles: Tri3, Tri6
  * 2D quads: Quad4, Quad8, Quad9
  * 3D tets: Tet4, Tet10
  * 3D hexes: Hex8, Hex20, Hex27
  * 3D pyramid: Pyr5
  * 3D wedges: Wedge6, Wedge15
- Added generation script block (~550 lines) that runs when file executed directly
- Generator now appends "Basis" suffix to all types (Tri3Basis, Quad4Basis, etc.)
- Outputs to src/basis/lagrange_generated.jl with clean formatting
- Includes progress reporting and next steps guidance

Total: 254 → 813 lines (+559 lines)

Run as: julia --project=. src/basis/lagrange_generator.jl
2025-11-09 08:28:50 +02:00
Jukka Aho f2492640e5 refactor(basis): Consolidate Lagrange basis includes
- Removed 7 individual lagrange_*.jl includes (segments, quadrangles, triangles,
  tetrahedrons, hexahedrons, wedges, pyramids)
- Added lagrange_generator.jl (generation infrastructure)
- Added lagrange_generated.jl (auto-generated basis functions for all 15 types)
- Updated comment explaining Basis suffix convention (Tri3Basis vs Tri3 topology)
- Removed TODO about name conflicts (resolved by Basis suffix pattern)
- Comment notes generator script location: scripts/generate_lagrange_basis.jl
2025-11-09 08:27:42 +02:00
Jukka Aho 31ecd6c0dc docs(book): Update Lagrange basis generation references
- Changed generator path: scripts/generate_lagrange_basis.jl → src/basis/lagrange_generator.jl
- Updated execution command: now run directly with julia --project=.
- Consolidated "See Also" section: removed duplicate generator reference
- Clarified generator role: symbolic engine AND generation script in single file
- Updated comment explaining basis functions are pregenerated (not runtime)
2025-11-09 08:27:18 +02:00
Jukka Aho 159e738038 docs(contributor): Reorder sections to emphasize coding standards
- Removed redundant H1 heading (already in YAML frontmatter)
- Moved "Coding Standards" above "Architecture" in What's Here section
- Updated "Before Contributing" list to prioritize standards (now item 2)
- Marked coding standards as REQUIRED for all contributions
- Changed reference from "Code Style" to "Coding Standards" (file renamed)
- Added blank line after "We assume you:" for better formatting
2025-11-09 08:26:49 +02:00
Jukka Aho f8851ef7d6 docs: Add comprehensive coding standards document
New 500-line standards document covering:
- Core principles (readability, type stability, zero allocations, explicit code)
- Variable naming: NO Greek letters in code (critical rule - use u,v,w not ξ,η,ζ)
- Type naming: PascalCase for types, snake_case for functions, Basis suffix pattern
- Performance guidelines: type stability, zero allocations, tuple returns
- Documentation style: docstrings with examples, theory, performance notes
- Testing standards: test organization, floating point comparisons
- Anti-patterns: Dict without types, abstract types in structs, globals, type piracy
- Git commit style: Conventional Commits format with examples
- Editor configuration: .editorconfig and JuliaFormatter.toml settings
- Summary checklist for pre-submission verification

Rationale for no Greek letters: keyboard accessibility, editor compatibility,
copy-paste issues, search/replace problems, terminal rendering, git diffs,
internationalization, and accessibility concerns.
2025-11-09 08:23:52 +02:00
Jukka Aho 241fadc669 docs: Add contributor quick-start guide
New file providing step-by-step onboarding for contributors:
- Quick links to contributor manual, coding standards, and testing philosophy
- 8-step workflow from fork to pull request
- Code of conduct principles (respectful, constructive, welcoming)
- Clear acceptance criteria (type stability, tests, documentation, clean commits)
- Rejection criteria (type instability, no tests, Greek letters, breaking changes)
- Help resources (discussions, issues, PRs)
- MIT license acknowledgment
2025-11-09 08:23:16 +02:00
Jukka Aho 77cb9f6394 docs(readme): Enhance Contributing section with guides and standards
- Expanded contributing text with clearer call to action
- Added links to CONTRIBUTING.md, coding standards, and contributor manual
- Added key requirements list (type stability, tests, Greek letter rule, clean commits)
- Fixed code block syntax highlighting (markdown → text) for Zen and citation
- Improved formatting with proper newlines before code blocks
- Added friendly questions/help invitation
2025-11-09 08:22:51 +02:00
Jukka Aho 6ca17e0569 Integrate topology/integration modules with comprehensive testing
INTEGRATION COMPLETE ✓
=======================

What's New:
-----------
- Integrated 17 topology types into main JuliaFEM module
- Integrated Gauss quadrature integration system
- Added comprehensive standalone test suite (36 tests, all passing)
- Documented topology coordinates for Hex20, Hex27, Pyr5, Quad8, Quad9, Tri7, Wedge6, Wedge15

Changes:
--------
src/JuliaFEM.jl:
  - Added topology module includes (17 topology types)
  - Added integration module includes (integration.jl, gauss.jl)
  - Exported all topology and integration symbols
  - Documented lagrange basis conflict (TODO for Phase 2)

test/test_topology_integration.jl (NEW):
  - Comprehensive test suite for full JuliaFEM integration
  - Tests all 17 topology types (1D, 2D, 3D)
  - Tests integration point generation for all topologies
  - Validates zero-allocation design
  - 370+ lines of test coverage

test/test_topology_standalone.jl (NEW):
  - Standalone validation tests (36/36 passing)
  - Tests topology module independently
  - Tests integration module independently
  - Bypasses name conflicts with old basis system
  - Proves core functionality correct

Topology Fixes:
  - Hex20, Hex27: Added proper node numbering documentation
  - Hex8: Fixed reference coordinates to match standard [-1,1]³
  - Pyr5: Fixed apex coordinate to (0,0,1)
  - Quad8, Quad9: Fixed midpoint coordinates
  - Tri7: Added standard node order
  - Wedge6, Wedge15: Fixed coordinate system

Documentation:
  - Updated book README with integration status
  - Updated contributor test fixes with topology integration notes

Test Results:
-------------
Topology standalone: 23/23 passed
  ✓ Seg2: nnodes, dim, coordinates
  ✓ Tri3: nnodes, dim, coordinates, edges
  ✓ Quad4: nnodes, dim, coordinates, edges
  ✓ Tet4: nnodes, dim, coordinates, edges, faces
  ✓ Hex8: nnodes, dim, coordinates, edges, faces

Integration standalone: 13/13 passed
  ✓ IntegrationPoint structure
  ✓ Gauss{1} + Tri3: 1 point at (1/3, 1/3), weight 0.5
  ✓ Gauss{3} + Tri3: 3 points, weights sum to 0.5
  ✓ Gauss{2} + Quad4: 4 points, weights sum to 4.0
  ✓ Gauss{1} + Tet4: 1 point (3D)
  ✓ Gauss{2} + Hex8: 8 points, weights sum to 8.0

Known Issue:
------------
Name conflict between topology types (Tri3 <: AbstractTopology) and
basis types (Tri3 <: AbstractBasis). Lagrange basis files currently
commented out to allow topology/integration to load. Will be resolved
in Phase 2 by renaming basis types (e.g., Tri3 -> Tri3Basis).

Zero-Allocation Design Verified:
---------------------------------
All topology and integration functions return tuples (immutable, stack-allocated).
No heap allocations in hot paths. Performance-critical design validated.

Next Steps:
-----------
1. Resolve name conflicts (rename basis types with *Basis suffix)
2. Refactor AbstractElement to accept separate topology/basis types
3. Run full test suite with integrated modules
4. Generate code coverage report
2025-11-09 06:13:40 +02:00
Jukka Aho b5fdf61851 feat(integration): Complete integration rule mappings for all topologies
**Added Gauss quadrature mappings for all 17 topology types**

Extended src/integration/gauss.jl to support all element types from 1D to 3D,
both linear and quadratic variants.

## Integration Rule Mappings

### 1D Segments (Seg2, Seg3)
- Tensor product rules: GLSEG1, GLSEG2, GLSEG3, GLSEG4, GLSEG5
- Support for Gauss{1} through Gauss{5}

### 2D Triangles (Tri3, Tri6, Tri7)
- Dedicated triangular rules: GLTRI1, GLTRI3, GLTRI4, GLTRI6, GLTRI7, GLTRI12
- Support for Gauss{1}, Gauss{3}, Gauss{4}, Gauss{6}, Gauss{7}, Gauss{12}
- Same rules used for linear (Tri3) and quadratic (Tri6, Tri7) topologies

### 2D Quadrilaterals (Quad4, Quad8, Quad9)
- Tensor product rules: GLQUAD1, GLQUAD4, GLQUAD9, GLQUAD16, GLQUAD25
- Support for Gauss{1} through Gauss{5}
- Same rules for linear (Quad4) and quadratic (Quad8, Quad9) variants

### 3D Tetrahedra (Tet4, Tet10)
- Dedicated tetrahedral rules: GLTET1, GLTET4, GLTET5, GLTET15
- Support for Gauss{1}, Gauss{4}, Gauss{5}, Gauss{15}

### 3D Hexahedra (Hex8, Hex20, Hex27)
- Tensor product rules: GLHEX1, GLHEX8, GLHEX27, GLHEX64, GLHEX125
- Support for Gauss{1} through Gauss{5}
- Same rules for linear (Hex8) and quadratic (Hex20, Hex27) variants

### 3D Wedges/Prisms (Wedge6, Wedge15)
- Dedicated wedge rules: GLWED6, GLWED21
- Support for Gauss{6}, Gauss{21}

### 3D Pyramids (Pyr5)
- Dedicated pyramid rules: GLPYR5
- Support for Gauss{5}

## Design Notes

**Quadrature rules from src/quadrature/**
All actual integration point data comes from src/quadrature/*.jl files
(consolidated from FEMQuad.jl). This file just maps high-level scheme + topology
to the appropriate low-level rule name.

**Tensor product elements:**
Segments, quads, and hexes use tensor product quadrature generated programmatically
in glquad.jl. Number follows pattern: N_points = N_per_dim^dimension
- GLSEG3 = 3 points in 1D
- GLQUAD9 = 3² = 9 points in 2D
- GLHEX27 = 3³ = 27 points in 3D

**Simplex elements:**
Triangles, tetrahedra use specialized rules (not tensor products) with optimized
point locations. Number roughly indicates integration order capability.

**Quadratic elements use same rules:**
Quadratic variants (Tri6, Quad8, Hex20, etc.) use same quadrature rules as
linear counterparts. User selects integration order via Gauss{N} parameter,
not topology type. Higher order topologies typically need higher N for exact
integration.

**Zero-allocation maintained:**
All functions return tuples, no heap allocation in integration point queries.

## Usage Examples

```julia
# Linear triangle with 1-point rule
ips = integration_points(Gauss{1}(), Tri3())

# Quadratic triangle with 6-point rule (more accurate)
ips = integration_points(Gauss{6}(), Tri6())

# Linear hex with 8-point rule (2³)
ips = integration_points(Gauss{2}(), Hex8())

# Quadratic hex with 27-point rule (3³)
ips = integration_points(Gauss{3}(), Hex27())
```

## Completeness

 All 17 topology types now supported
 Linear and quadratic variants covered
 1D, 2D, and 3D elements complete
 Zero-allocation design maintained

## References

- src/quadrature/glquad.jl (tensor product generation)
- src/quadrature/gltri.jl (triangle rules)
- src/quadrature/gltet.jl (tetrahedron rules)
- src/quadrature/glwed.jl (wedge rules)
- src/quadrature/glpyr.jl (pyramid rules)
- Dunavant, "High degree efficient symmetrical Gaussian quadrature rules for the triangle"
- Abramowitz & Stegun, "Handbook of Mathematical Functions"
2025-11-09 06:01:01 +02:00
Jukka Aho d18622d41c feat(topology): Complete topology library with all element types
**Implemented 14 additional topology types with zero-allocation interfaces**

This completes the topology module with all standard FEM element types from
1D to 3D, both linear and quadratic variants.

## New Topologies

### 1D Elements (Segments)
- Seg2: 2-node linear segment
- Seg3: 3-node quadratic segment

### 2D Elements
**Triangles:**
- Tri6: 6-node quadratic triangle
- Tri7: 7-node quadratic triangle (with center node)

**Quadrilaterals:**
- Quad8: 8-node quadratic quad (Serendipity)
- Quad9: 9-node quadratic quad (with center node)

### 3D Elements
**Tetrahedra:**
- Tet4: 4-node linear tetrahedron
- Tet10: 10-node quadratic tetrahedron

**Hexahedra:**
- Hex8: 8-node linear hexahedron
- Hex20: 20-node biquadratic hexahedron (Serendipity)
- Hex27: 27-node quadratic hexahedron (with face/volume nodes)

**Pyramids:**
- Pyr5: 5-node linear pyramid

**Wedges/Prisms:**
- Wedge6: 6-node linear wedge
- Wedge15: 15-node quadratic wedge

## Design Principles

**Zero-allocation throughout:**
- reference_coordinates() → NTuple{N, NTuple{D, Float64}}
- edges() → NTuple{Ne, Tuple{Int, Int}}
- faces() → NTuple{Nf, NTuple{Nn, Int}} or NTuple{Nf, Tuple{Vararg{Int}}}

All topology data is stack-allocated, compile-time sized tuples. No heap
allocations in hot assembly loops.

**Reference coordinates extracted from existing basis files:**
- src/basis/lagrange_segments.jl
- src/basis/lagrange_triangles.jl
- src/basis/lagrange_quadrangles.jl
- src/basis/lagrange_tetrahedrons.jl
- src/basis/lagrange_hexahedrons.jl
- src/basis/lagrange_pyramids.jl
- src/basis/lagrange_wedges.jl

**Complete topology coverage:**
- 1D: linear and quadratic segments
- 2D: triangles (3,6,7 nodes), quads (4,8,9 nodes)
- 3D: tets (4,10), hexes (8,20,27), pyramids (5), wedges (6,15)

This matches the rich set of elements JuliaFEM supported historically.

## Implementation Notes

**Edge/Face Connectivity:**
- edges(): Corner nodes only (defines element boundary)
- faces(): For 2D elements, all nodes; for 3D elements, corner nodes of each face
- Consistent with standard FEM conventions

**Pyramid Special Case:**
- Pyr5 uses Code Aster convention (from lagrange_pyramids.jl)
- Base at z=-1, apex at z=+1
- Mixed face types: 1 quad base + 4 triangular faces

**Wedge/Prism Special Case:**
- Triangular cross-section extruded along w-axis
- Mixed face types: 2 triangular + 3 quadrilateral faces

## Status

Total topology types: 17 (Seg2, Seg3, Tri3, Tri6, Tri7, Quad4, Quad8, Quad9,
Tet4, Tet10, Hex8, Hex20, Hex27, Pyr5, Wedge6, Wedge15)

**Not yet integrated** into src/JuliaFEM.jl (staged approach).

## Next Steps

1. Update topology.jl to export all types
2. Update src/JuliaFEM.jl to include all topology files
3. Add integration rules for all topologies in src/integration/gauss.jl
4. Generate Lagrange basis functions for all topologies

## References

- Existing basis files in src/basis/ (reference coordinate source)
- Abaqus Theory Manual (standard element definitions)
- Code Aster documentation (pyramid element convention)
- TECHNICAL_VISION.md (zero-allocation design philosophy)
2025-11-09 05:58:34 +02:00
Jukka Aho ee02f9f37a feat: Separation of concerns architecture with zero-allocation foundation
**Architecture Decision: Element = Topology + Interpolation + Integration + Fields**

This commit establishes the architectural foundation for separating orthogonal concerns
in finite element implementation, preventing Abaqus-style combinatorial explosion.

## New Modules (Not Yet Integrated)

### src/topology/
Reference element geometries (pure mathematical objects):
- topology.jl: Abstract interface for reference elements
- tri3.jl: 3-node triangle reference element
- quad4.jl: 4-node quadrilateral reference element

**Zero-allocation design:**
- reference_coordinates() → NTuple{N, NTuple{D, Float64}}
- edges() → NTuple{Ne, Tuple{Int, Int}}
- faces() → NTuple{Nf, NTuple{Nn, Int}}

All topology queries return compile-time sized tuples (stack allocated, no heap).

### src/integration/
High-level integration scheme abstraction:
- integration.jl: Abstract types and IntegrationPoint struct
- gauss.jl: Gauss-Legendre quadrature wrapper around existing src/quadrature/

**Zero-allocation design:**
- integration_points() → Tuple{Vararg{IntegrationPoint{D}}}
- IntegrationPoint.ξ → NTuple{D, Float64}

**Key Insight:** Integration rules already exist in src/quadrature/ (consolidated from
FEMQuad.jl). New code is a thin architectural wrapper, not reimplementation.

## Documentation

### docs/book/element_architecture.md (NEW - 650+ lines)
Complete book chapter explaining:
- What is an Element? (composition of 4 orthogonal concerns)
- The Abaqus anti-pattern (C3D8, C3D8R, C3D8I explosion)
- JuliaFEM approach: Topology + Interpolation + Integration separation
- Type system enforcement
- Performance implications (100× speedup from type stability)
- Extending the system (adding new topologies/bases/quadrature)
- Comparison with Gridap.jl, Ferrite.jl, Deal.II

### llm/ARCHITECTURE.md (UPDATED)
Added "Architectural Decision: Separation of Concerns" section at top:
- Problem statement
- Anti-pattern example
- JuliaFEM solution
- Directory structure rationale
- Type system design
- Migration strategy

### scripts/generate_lagrange_basis.jl (UPDATED)
Added architectural context explaining Lagrange bases are INTERPOLATION SCHEMES
(not topologies, not integration rules).

## Performance: Zero-Allocation Foundation

**Why tuples matter:**
1. **Zero heap allocations** - All data stack-allocated
2. **Compile-time sizes** - Compiler can unroll loops
3. **Cache friendly** - Contiguous memory layout
4. **Type stable** - Concrete tuple types enable optimization
5. **Immutable** - No accidental mutation, thread-safe

**Example impact:**
```julia
# Compiler knows at compile time:
# - Tri3 has exactly 3 edges
# - Each edge has exactly 2 nodes
# → Loop unrolling, no bounds checks, SIMD vectorization

for edge in edges(Tri3())  # Tuple iteration, fully unrolled!
    node1, node2 = edge
    # ... assembly code (zero allocations)
end
```

**Principle from Roadmap to HPC:**
> "Zero allocations in hot paths" - Strategic Decision #2

Topology/integration queries happen billions of times in assembly loops.
Even small Vector allocations accumulate to GC pressure and cache misses.

**Rule:** If size known at compile time → use Tuple, not Vector

## Benefits

 Clear separation of mathematical concepts
 Mix-and-match: Tri3 + Lagrange + Gauss, Tri3 + Hierarchical + Lobatto, etc.
 Type system enforces correctness at compile time
 Compiler generates specialized code for each combination → 100× speedup
 Zero allocations in topology/integration queries
 No code duplication (each concern in one place)
 Educational: teaches proper software engineering

## Status

- **NOT YET INTEGRATED**: New modules not included in src/JuliaFEM.jl
- **SAFE**: Package loads successfully (verified with `using JuliaFEM`)
- **READY**: Architecture documented, zero-alloc foundation established

## Next Steps

1. Create remaining topology files (Tet4, Tet10, Hex8, Hex20, etc.)
2. Update src/JuliaFEM.jl to include new modules
3. Refactor existing Element to use new separation
4. Run generation script with new architecture
5. Integrate with existing codebase

## References

- Abaqus documentation (anti-pattern example)
- Gridap.jl (alternative approach)
- Ferrite.jl (mixed approach)
- Deal.II (C++ template approach)
- llm/ROADMAP_TO_HPC.md (performance philosophy)

See: docs/book/element_architecture.md for complete rationale and examples.
2025-11-09 05:46:34 +02:00
Jukka Aho 91b06b23b6 fix: Re-enable lagrange_generator.jl include for existing basis files
**Problem:** CI documentation build failing with:
  UndefVarError: `create_basis_and_eval` not defined

**Root Cause:**
- Commented out include("basis/lagrange_generator.jl")
- But existing lagrange_*.jl files still use create_basis_and_eval()
- Generator must be loaded at runtime (for now)

**Solution:**
- Re-enable include("basis/lagrange_generator.jl")
- Add TODO comment explaining this is temporary
- Once we generate lagrange_generated.jl, we can:
  1. Remove old lagrange_*.jl includes
  2. Include lagrange_generated.jl instead
  3. Remove lagrange_generator.jl from runtime

**Status:**
- ✓ Package loads successfully
- ✓ Should fix CI documentation build

**Next Steps:**
1. Run scripts/generate_lagrange_basis.jl
2. Create src/basis/lagrange_generated.jl
3. Replace dynamic includes with static generated file
4. Remove __precompile__(false) completely
2025-11-09 05:03:53 +02:00
Jukka Aho c65abfa5cc docs: Add 'Roadmap to HPC' - justifying hard performance choices
**Purpose:** Comprehensive justification for all technical decisions prioritizing
performance over convenience.

**Key Principles:**
- Efficiency > Educativeness (when forced to choose)
- Type stability over everything (100× performance difference)
- No free lunch - Julia doesn't make miracles
- HPC requires discipline and trade-offs

**Core Decisions Justified:**

1. **No Dynamic Field System**
   - field["foo"] = x is 100× slower (Dict{String,Any})
   - Type-stable structs only
   - Sacrifice: Runtime flexibility
   - Gain: Performance

2. **Immutable Data Structures**
   - struct over mutable struct
   - Sacrifice: Convenient mutation
   - Gain: 2-10× speedup, thread-safety, stack allocation

3. **NTuple Over Vector**
   - Compile-time size → SIMD optimization
   - Sacrifice: Dynamic sizing
   - Gain: Zero allocations, type stability

4. **Monolithic Over Multi-Package**
   - Learned from 2015-2019 mistake
   - Sacrifice: Small dependencies
   - Gain: It actually works

5. **Manual Derivatives (hot paths)**
   - 30× faster than AD for Tet10
   - Sacrifice: More code
   - Gain: Assembly loops stay fast

6. **Matrix-Free Methods**
   - Design for 1M+ DOF from day 1
   - Cannot retrofit later

7. **Explicit Over Implicit**
   - No magic, show the steps
   - Debuggable and teachable

**Hierarchy of Values:**
1. Correctness
2. Performance
3. Maintainability
4. Educativeness
5. Convenience

**What We're Giving Up:**
- Runtime flexibility (no element["custom_field"])
- Dynamic problem definition (no runtime topology changes)
- Duck typing convenience
- Small dependencies
- Beginner-friendly magic

**What We're Getting:**
- 10× single-thread speedup target
- 1M DOF contact problems
- Thread/GPU/distributed scalability
- Real HPC capability

**The Hard Truth:**
From Issue #266: "Do like Python, be slow like Python. Know what you do
before compiling, and be fast like C. There's no free lunch."

**Success Metrics:**
-  Zero allocations in assembly
-  Type-stable hot paths
- 🎯 10× faster than v0.5.1
- 🎯 1M DOF in < 1 hour
- 🎯 100+ thread scaling

**Use Cases:**
- "Why can't I use Dict?" → Point here
- "Why immutable?" → Point here
- "Why manual derivatives?" → Point here
- Any "why not convenience?" → Point here

**Status:** Living document, updated as we learn

See: Issue #266, TECHNICAL_VISION.md, benchmark results
2025-11-09 05:01:34 +02:00
Jukka Aho 1636e255fe docs: Add YAML front matter to all documentation files
**Purpose:** Prepare documentation for publishing as blog posts or book

**YAML Headers Include:**
- title: Document title
- subtitle: Optional subtitle for context
- description: Brief summary for SEO/indexing
- date: Creation date
- updated: Last update date (for status docs)
- author: Jukka Aho
- categories: Taxonomic classification
- keywords: Search/indexing keywords
- audience: Target reader (users/contributors/researchers)
- level: Difficulty level (beginner/intermediate/advanced/expert)
- type: Document type (manual/guide/theory/benchmark/status)
- series: Which manual it belongs to
- chapter: Book structure (for The JuliaFEM Book)
- status: Current state (completed/work in progress/active maintenance)
- math: Whether document contains mathematical notation
- prerequisites: Required background knowledge
- tools: Software/packages used (for benchmarks)
- context: Background information

**Files Updated:**
- docs/README.md (main index)
- docs/user/README.md (user manual index)
- docs/contributor/README.md (contributor manual index)
- docs/book/README.md (book index)
- docs/contributor/testing_philosophy.md
- docs/contributor/status.md
- docs/contributor/test_fixes_needed.md
- docs/book/lagrange_basis_functions.md
- docs/book/benchmarks/shape_function_derivatives_ad_vs_manual.md
- scripts/README.md

**Benefits:**
- Ready for static site generators (Jekyll, Hugo, MkDocs)
- Can generate book with proper metadata
- SEO-friendly with descriptions and keywords
- Clear audience/level targeting
- Trackable with dates and status
- Organized by series and chapters

**Compatible With:**
- Jekyll (GitHub Pages)
- Hugo (fast static site generator)
- MkDocs (Python-based documentation)
- Jupyter Book (interactive books)
- Docusaurus (React-based docs)
- Custom publishing scripts
2025-11-09 04:45:12 +02:00
Jukka Aho 626266c990 docs: Reorganize documentation into three-tier structure
**Three Manuals for Three Audiences:**

1. **User Manual** (docs/user/) - "Just Get It Done"
   - For end users, engineers, students
   - Simple, practical, step-by-step
   - Quick start, tutorials, examples, troubleshooting
   - Philosophy: Show me how to solve my problem

2. **Contributor Manual** (docs/contributor/) - "Show Me the Code"
   - For developers, contributors, advanced users
   - Technical, detailed, design rationale
   - Testing, architecture, performance, CI/CD
   - Philosophy: Explain HOW and WHY

3. **The JuliaFEM Book** (docs/book/) - "Let Me Show You How I Think"
   - For researchers, theory nerds, and Jukka
   - Comprehensive, educational, opinionated, personal
   - Math foundations, design philosophy, history, research
   - Philosophy: Mix theory, code, and personal experience

**Reorganization:**
- Moved: TESTING_PHILOSOPHY.md → contributor/testing_philosophy.md
- Moved: STATUS.md → contributor/status.md
- Moved: TEST_FIXES_NEEDED.md → contributor/test_fixes_needed.md
- Moved: lagrange_basis_functions.md → book/lagrange_basis_functions.md
- Moved: benchmarks/ → book/benchmarks/
- Created: docs/README.md (main index explaining structure)
- Created: README.md in each section explaining audience and contents
- Updated: All references in scripts and source files

**Naming:** All docs now lowercase (testing_philosophy not TESTING_PHILOSOPHY)

**Benefits:**
- Clear separation of concerns
- Users don't get overwhelmed with implementation details
- Contributors get technical depth
- Book preserves deep theory and personal insights
- Each manual optimized for its audience

**Next:** Populate each section with appropriate content
2025-11-09 04:38:28 +02:00
Jukka Aho 5141fd6de5 refactor: Move theory docs to src/ with lowercase naming
- Moved docs/theory/lagrange_basis_functions.md → src/lagrange_basis_functions.md
- Updated all references in scripts and source files
- Using lowercase for consistency (no uppercase in filenames)
- Documentation now under src/ for automated doc generation

Rationale: Documentation should be close to implementation and follow
consistent naming conventions (lowercase).
2025-11-09 04:26:26 +02:00
Jukka Aho 63346a0591 style: IDE automatic code formatting
No functional changes - only whitespace and formatting adjustments:
- Removed spaces around = in named tuple syntax (name = → name=)
- Adjusted spacing in array literals
- Standardized spacing around operators
2025-11-09 04:10:37 +02:00
Jukka Aho 31d8463ef0 feat: Pre-generation infrastructure for Lagrange basis functions
**Problem:**
- __precompile__(false) in create_basis.jl causes slow package loading
- Symbolic math evaluated at runtime (100+ ms overhead)
- Dynamic eval() prevents full precompilation
- Difficult to debug generated code

**Solution: Generate Once, Use Forever**
- Renamed: create_basis.jl → lagrange_generator.jl (tool, not runtime code)
- Created: scripts/generate_lagrange_basis.jl (orchestration script)
- Created: scripts/README.md (documentation for generation workflow)
- Created: docs/theory/lagrange_basis_functions.md (mathematical foundation)

**Theory Documentation (400+ lines):**
- Kronecker delta property: N_i(x_j) = δ_ij
- Vandermonde matrix method: Vα_i = e_i
- Worked example: Seg2 linear element (step-by-step derivation)
- Polynomial completeness table (1D/2D/3D orders)
- Complete standard element catalog
- Pre-generation vs runtime comparison
- Numerical stability discussion

**Generation Script:**
- Defines all 15 standard Lagrange element types:
  * 1D: Seg2, Seg3
  * 2D Tri: Tri3, Tri6
  * 2D Quad: Quad4, Quad8, Quad9
  * 3D Tet: Tet4, Tet10
  * 3D Hex: Hex8, Hex20, Hex27
  * 3D Pyr: Pyr5
  * 3D Wedge: Wedge6, Wedge15
- For each: node coordinates + polynomial ansatz
- Calls lagrange_generator symbolic engine
- Writes clean Julia code → src/basis/lagrange_generated.jl (to be created)

**Architecture:**

**Benefits:**
- ~150× faster package loading (150ms → <1ms)
- Full precompilation enabled
- Generated code is readable/debuggable
- Git shows what changed (mathematics visible in diffs)
- Reproducible builds

**Workflow:**
1. Edit element catalog in scripts/generate_lagrange_basis.jl
2. Run: julia --project=. scripts/generate_lagrange_basis.jl
3. Review src/basis/lagrange_generated.jl
4. Test and commit

**Next Steps:**
1. Run generation script → create lagrange_generated.jl
2. Update src/JuliaFEM.jl to include generated file
3. Comment out old lagrange_*.jl includes
4. Remove __precompile__(false)
5. Verify all tests pass
6. Measure package load time improvement

**Also Included:**
- scripts/check_namespace_collisions.jl (consolidation tool)
- scripts/fix_vendor_element_types.py (Element type fixer)

See: docs/theory/lagrange_basis_functions.md for full mathematical explanation
2025-11-09 04:07:28 +02:00
Jukka Aho 6a8f8adc1f docs: Benchmark manual vs AD derivatives for Tet10
RESEARCH QUESTION: Should JuliaFEM use hand-calculated derivatives or AD?

Created comprehensive benchmark comparing:
- Manual: Hand-calculated derivatives (traditional FEM)
- AD: Tensors.jl gradient() (automatic differentiation)

RESULTS (AMD Ryzen 9, Julia 1.12.1):
- Manual: 8.7 ns, 0 allocations
- AD:     268.1 ns, 0 allocations
- AD is 30× SLOWER than manual

KEY FINDINGS:
 Both achieve zero allocations (Tensors.jl is well-optimized)
 AD has 30× compute overhead from dual number arithmetic
⚠️  In assembly loops: millions of calls = 10+ seconds extra per solve

RECOMMENDATION:
- Keep manual derivatives for common elements (Tet10, Hex8, Quad4, etc.)
- Use AD for prototyping and rare elements
- Unit test manual vs AD to catch errors
- Future: Generate derivatives symbolically (Symbolics.jl)

WHY NOT AD EVERYWHERE?
Assembly is hottest path in FEM. 30× overhead = unacceptable for
production code. Users will notice the performance difference.

WHY NOT ABANDON AD?
- Excellent for prototyping
- Required for exotic bases (NURBS)
- Perfect for unit testing manual derivatives
- Zero allocations impressive

Files:
- benchmarks/tet10_derivatives_benchmark.jl (runnable benchmark)
- docs/benchmarks/shape_function_derivatives_ad_vs_manual.md (analysis)

Dependencies added: BenchmarkTools

This answers the research question definitively with data.
2025-11-09 03:41:47 +02:00
Jukka Aho 6b24ed9d76 refactor: Make Point immutable
Changed Point from 'mutable struct' to 'struct'.

The Dict for fields remains a reference type, so field updates via setindex!
and update! still work correctly. This change improves type stability and
enables better compiler optimizations.

Benefits:
- Better compiler optimizations (immutable types)
- Type stability improvements
- Stack allocation when possible
- No breaking changes (Dict fields still mutable)

Tests: All 157 tests passing
2025-11-09 03:30:29 +02:00
Jukka Aho 907ec0b183 refactor: Zero-allocation basis functions and immutable Element
MAJOR PERFORMANCE REFACTORING:

1. Shape functions return tuples instead of allocating vectors:
   - eval_basis!(): Returns NTuple{N,T} directly (zero allocations)
   - eval_dbasis!(): Returns NTuple{N,Vec{D}} directly (zero allocations)
   - API boundary (get_basis/get_dbasis) still returns vectors for compat

2. Element is now immutable with compile-time known structure:
   - connectivity: Vector{UInt} → NTuple{N,UInt}
   - integration_points: Vector{IP} → NTuple{NIP,IP}
   - Element{N,NIP,M,B} parametrized by connectivity/IP count
   - Changed from 'mutable struct' to 'struct'

3. Helper function for immutability:
   - with_integration_points(element, ips) returns new element
   - get_integration_points() returns tuple directly

Benefits:
- Zero allocations in hot paths (basis evaluation)
- Compile-time sizes enable better optimization
- Type stability improvements
- Stack allocation instead of heap

Breaking changes:
- Element.connectivity is now tuple (use collect() for vector)
- Element is immutable (use with_integration_points for updates)

Tests: All 157 tests passing
2025-11-09 03:29:36 +02:00
Jukka Aho 065156b40a style: Format core_types.jl (spacing consistency) 2025-11-09 03:18:41 +02:00
Jukka Aho 06e8276268 fix: Change node and element IDs to UInt (Issue #267)
Gmsh returns node and element IDs as UInt64, so we should use unsigned
integers consistently throughout JuliaFEM to avoid unnecessary conversions.

Changes:
- Point.id: Int → UInt
- Element.id: Int → UInt
- Element.connectivity: Vector{Int} → Vector{UInt}
- Element constructors: Accept Integer (converts to UInt internally)
- Default element_id: -1 → 0 (UInt has no negative values)

Benefits:
- Direct compatibility with Gmsh.jl (no Int/UInt conversions)
- Semantically correct (node/element IDs are never negative)
- Slightly more efficient (no sign checks)

Tests: All 156 tests passing

Closes #267
2025-11-09 03:17:34 +02:00
Jukka Aho 52ebe682e9 fix: Standardize on Tensors.jl Vec type throughout
Major architectural decision: Use Tensors.jl consistently everywhere
for geometric vectors, integration points, and coordinates.

Changes to src/elements/elements.jl:
- get_basis(): Convert ip to Vec, use Vector (not Matrix) for eval_basis!
- get_dbasis(): Convert ip to Vec
- jacobian evaluation: Convert geometry and ip.coords to Vec properly
- Handle both raw coordinates (Tuple) and IP struct transparently

New Tutorial 3: Numerical Integration and Jacobian (49 tests)
- Integration point structure and weights
- Jacobian determinant and matrix evaluation
- Numerical integration (constant, linear, quadratic functions)
- Multiple element types (Quad4, Seg2, Tri3)

Tests: 107 → 156 passing (49 new)
Runtime: ~7 seconds

Closes architectural standardization on Tensors.jl.
Related to Issue #250 (merge conflict resolution).

Why Tensors.jl:
- Type stability (100× performance vs Dict-based)
- Material science compatibility (stress tensors)
- Zero-cost abstractions
- Consistent API across all geometric calculations
2025-11-09 03:10:11 +02:00
Jukka Aho 5a07b3ab21 docs: Document Tutorial 3 API limitations, update test runner
Current state discovery:
- Element basis function evaluation broken (eval_basis! signature mismatch)
- Field interpolation at integration points broken (same root cause)
- Jacobian evaluation at integration points broken
- These are fundamental API issues affecting multiple test paths

Impact:
- Tutorial 3 (basis functions) deferred until API fixed
- Affects any code trying to evaluate fields at integration points
- Related to Quad4 assembly issues discovered in Tutorial 4

Working tutorials (107/107 tests passing):
- Tutorial 1: Element creation (5 tests)
- Tutorial 2: Gmsh mesh reading (72 tests)
- Tutorial 4: 1-element validation (35 tests)

Next: Focus on tutorials using working APIs only
2025-11-09 02:58:09 +02:00
Jukka Aho bafa3af4d0 test: Add Tutorial 4 - 1-element Quad4 validation (35 tests passing)
Educational validation test for Issue #265 use case (JuliaFEM as reference).

Covers:
- Element creation and connectivity
- Field assignment (geometry, material properties)
- Field retrieval with function call syntax
- Hand-calculated constitutive matrix for plane stress
- Geometry validation (dimensions, center, area)
- Material property validation (physical ranges)

Note: Defers stiffness matrix assembly to future work due to current
Quad4 assembly issues. Focus is on element setup validation that
other FEM developers can use as reference.

Tutorial series now: 107/107 tests passing
- Tutorial 1: Creating elements (5 tests)
- Tutorial 2: Gmsh mesh reading (72 tests)
- Tutorial 4: 1-element validation (35 tests - done before Tutorial 3)
2025-11-09 02:45:29 +02:00
Jukka Aho 7571487e86 docs: Update testing philosophy with current progress
Updates based on actual implementation:
- Gmsh chosen over ABAQUS (accessibility, no license needed)
- Co-located mesh files with recipe scripts (reproducible)
- Realistic mesh sizes (~10 elements, not 1-4)
- 1-element validation tests prioritized (Issue #265)
- Progress tracking: 77/77 tests passing (Tutorial 1-2 complete)
- Mesh generation pattern documented (recipe + .msh + test)

New section: 1-Element Validation Tests
- Motivation from Issue #265 (JuliaFEM validated other FEM software)
- Hand-calculable reference solutions
- High priority for Tutorial 4
2025-11-09 02:36:10 +02:00