Glassvision/JuliaFEM.jl
1
0
Fork 0
mirror of https://github.com/JuliaFEM/JuliaFEM.jl.git synced 2026-09-11 06:08:07 +00:00
Code Issues Packages Projects Releases Wiki Activity

docs(elements): rewrite README for Element{K,P,S,N} and DOFHandler

Browse Source 
Replace outdated NamedTuple dof_indices narrative with flat tuples,
`local_dof_layout`, and `create_elements!` workflows.

- Align examples with `@DOFSet` and multi-field thermo-mechanical setup.
- Document extraction/interpolation entry points and compile-time `field_dof_range`.
- Point readers at dofs/topology/basis/assemblers tests without referencing removed files.
This commit is contained in:
Jukka Aho
2026-05-09 17:02:26 +03:00
parent 00ef826dd3
commit dda48b88fa
1 changed files with 83 additions and 164 deletions
Download Patch File Download Diff File Expand all files Collapse all files
src/elements/README.md
+83 -164
View File
@@ -1,219 +1,138 @@
# Elements Module
# src/elements/
This module implements finite elements following Ciarlet's mathematical definition, adapted for computational efficiency.
The `Element{K, P, S, N}` template plus the DOF-extraction and
field-interpolation utilities used by the assembly kernels.
## Ciarlet's Finite Element Triple (K, P, Σ)
## Mathematical background
### Mathematical Definition
A finite element is, in Ciarlet's sense, a triple `(K, P, Σ)` where
A finite element is a triple **(K, P, Σ)** where:
- `K` is the reference domain (a topology),
- `P` is the local approximation space (a basis family), and
- `Σ` is a set of linear functionals (degrees of freedom).
- **K** ⊂ ℝⁿ: Compact, connected reference domain (geometric set)
- **P**: Finite-dimensional space of functions on K
- **Σ** = {σ₁, ..., σₙ}: Set of linear functionals σᵢ : P → ℝ (dual basis)
In the implementation `Σ` is encoded by a field specification `S` that
together with `K` and `P` uniquely determines the functionals for the
standard Lagrange / Serendipity families. `S` is the DOFSet built by the
`@DOFSet` macro in `src/dofs/`.
### Computational Implementation
We use **(K, P, S)** where:
- **K**: Reference domain type (e.g., `Triangle{3}`, `Tetrahedron{4}`) - **exact match**
- **P**: Polynomial space type (e.g., `Lagrange{1}`, `Lagrange{2}`) - **exact match**
- **S**: Field specification → **uniquely determines Σ** (computational encoding)
### Why S Instead of Σ?
**S does not equal Σ, but S determines Σ uniquely.**
For standard Lagrange elements:
| S specification | Resulting Σ functionals | Example |
|----------------|------------------------|---------|
| `Float64, Vertex` | σᵢ(u) = u(vertex_i) | Point evaluation (nodal values) |
| `Vec{3}, Vertex` | σᵢ(u) = uₐ(vertex_i), α=1,2,3 | Vector point evaluation |
| `Float64, Cell` | σ(u) = (1/\|K\|) ∫_K u dx | Cell-average functional |
| `Float64, Edge` | σ(u) = ∫_edge u ds | Edge integral functional |
**Rationale:**
1. Functionals are never instantiated in computational FEM
2. S contains the essential information: quantity type + entity location
3. Given (K, P, S), the functionals Σ are uniquely determined
4. Type-level encoding = zero runtime cost
## Element Structure
## The element type
```julia
struct Element{K<:AbstractTopology, P<:AbstractBasis, S<:DOFSet, N}
id::UInt # Element identifier (mesh index)
dof_indices::NTuple{N,UInt64} # Flat tuple of global DOF indices
id::UInt
dof_indices::NTuple{N, UInt64}
end
```
### Design Philosophy
- `K` and `P` are types — no runtime fields.
- `S` is a NamedTuple type whose values are `DOF{Quantity, Entity}`
(see `src/dofs/README.md`).
- `N` is the total number of local DOFs (computed by the constructor).
- `dof_indices` is a flat tuple of global DOF indices in the order
defined by `local_dof_layout(::Type{Element{K, P, S, N}})`.
**Everything mathematical lives in the types.** The instance holds only:
Use `create_elements!(mesh, Element{K, P, S})` to build a
`Vector{Element{K, P, S, N}}` together with a `DOFHandler` that already
carries the inverse DOF connectivity.
- Identification (`id`)
- Assignment (`dof_indices`)
## Compile-time DOF layout
No connectivity, no coordinates stored in element! Mesh holds geometric data.
### Type Stability via @generated Constructor
The `dof_indices` field is typed as `NamedTuple` (without parameters), but the `@generated` constructor ensures the concrete type is inferred:
`local_dof_layout(::Type{Element{K, P, S, N}})` is a `@generated`
function returning `NTuple{N, DOFLayoutEntry}`; each entry exposes
`field_idx`, `entity_local`, `component`. The compiler folds the result
into a constant at the call site, so DOF decoding inside hot loops is a
tuple lookup with no arithmetic.
```julia
@generated function Element{K,P,S}(id::UInt, dof_indices::D) where {K,P,S,D<:NamedTuple}
# Julia infers D = @NamedTuple{u::NTuple{12, Int64}} from the argument
# Field access elem.dof_indices.u returns NTuple{12, Int64} (type-stable!)
end
S = @DOFSet{u::DOF{Displacement{3}, Vertex}}
ET = Element{Hex8, Lagrange{1}, S, 24}
local_dof_layout(ET)
```
This achieves zero-allocation performance without adding a 4th type parameter.
## Building elements
## Field Specifications
### Single-Field Elements
### Single-field
```julia
# Heat conduction (scalar field at vertices)
S = @NamedTuple{T::Tuple{Float64, Vertex}}
Element{Triangle{3}, Lagrange{1}, S}(UInt(1), (T=(1, 2, 3),))
S = @DOFSet{u::DOF{Displacement{3}, Vertex}}
ET = Element{Tetrahedron{4}, Lagrange{1}, S}
# 2D elasticity (vector field at vertices)
S = @NamedTuple{u::Tuple{Vec{2}, Vertex}}
Element{Triangle{3}, Lagrange{1}, S}(UInt(1), (u=(1, 2, 3, 4, 5, 6),))
elements, handler = create_elements!(mesh, ET)
```
### Multi-Field Elements
### Multi-field
```julia
# Thermo-mechanical coupling
S = @NamedTuple{
T::Tuple{Float64, Vertex}, # Temperature at vertices
u::Tuple{Vec{3}, Vertex} # Displacement at vertices
}
S = @DOFSet{T::DOF{Temperature, Vertex},
u::DOF{Displacement{3}, Vertex}}
Element{Tetrahedron{4}, Lagrange{1}, S}(
UInt(1),
(T=(1,2,3,4), u=(5,6,7,8,9,10,11,12,13,14,15,16))
)
# Access fields directly
elem.dof_indices.T # (1, 2, 3, 4)
elem.dof_indices.u # (5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16)
elements, handler = create_elements!(mesh, Element{Tetrahedron{4}, Lagrange{1}, S})
```
## DOF Extraction
In both cases `dof_indices` is a flat `NTuple` whose ordering is dictated
by `local_dof_layout`.
Two extraction strategies for getting element DOFs from global solution:
## DOF extraction from a global vector
### Flat Extraction
Two extraction strategies are provided. Both are zero-allocation and
type-stable.
Returns scalars grouped by field:
### Flat extraction
```julia
u_global = [1.0, 2.0, ..., 20.0]
dofs = extract_element_dofs(elem, u_global)
# Returns: (u = (1.0, 2.0, 3.0, ..., 12.0),)
# (u = (1.0, 2.0, ..., 12.0),)
```
### Structured Extraction
### Structured extraction
Returns quantities matching field type (Vec, Tensor, etc.):
Reinterprets the values into the field's quantity type so that they can
be combined with shape-function values directly.
```julia
dofs = extract_element_dofs_structured(elem, u_global)
# Returns: (u = (Vec{3}(1,2,3), Vec{3}(4,5,6), Vec{3}(7,8,9), Vec{3}(10,11,12)),)
# (u = (Vec{3}(1,2,3), Vec{3}(4,5,6), ...),)
u_at_xi = N1 * dofs.u[1] + N2 * dofs.u[2] + N3 * dofs.u[3] + N4 * dofs.u[4]
```
**Use case:** Structured extraction is for interpolation where tuple length must match shape function count:
## Field-block ranges (multi-field)
Per-field local index ranges are computed at compile time and are useful
for picking out coupling sub-blocks of an element matrix:
```julia
u_interp = N1 * u1 + N2 * u2 + N3 * u3 + N4 * u4
T_range = field_dof_range(elem, :T) # 1:4 for Tet4 + Vertex
u_range = field_dof_range(elem, :u) # 5:16
K_Tu = K_local[T_range, u_range]
```
Both are **zero-allocation** (5.5 ns) thanks to type stability and `@generated` functions.
## Local-Global DOF Mapping
For coupled multi-field assembly:
## Type queries
```julia
# Element with 2 fields: T (4 DOFs) + u (12 DOFs) = 16 total
map = local_to_global_map(elem)
# map[1:4] = [1,2,3,4] Temperature DOFs
# map[5:16] = [10,...,21] Displacement DOFs
# Assembly loop
K_local = zeros(16, 16) # Fully coupled local matrix
# ... fill K_local with physics coupling (∂T/∂u, ∂u/∂T, etc.) ...
for i in 1:16, j in 1:16
K_global[map[i], map[j]] += K_local[i, j]
end
topology_type(elem) # K
basis_type(elem) # P
dof_type(elem) # S
n_element_dofs(elem) # N
nnodes(elem) # nnodes(K)
```
### Field-Specific DOF Ranges
## Files
Extract local DOF ranges for field blocks (compile-time computation):
- `elements.jl` — element type, constructors, type queries,
`local_dof_layout`.
- `extract_element_dofs.jl` — flat and structured DOF extraction.
- `interpolate.jl` — field interpolation at points
(`interpolate_field`, `interpolate_fields`,
`interpolate_field_value`, `interpolate_local_fields`).
```julia
T_range = field_dof_range(elem, :T) # 1:4
u_range = field_dof_range(elem, :u) # 5:16
## Related code
# Extract field-field coupling block
K_Tu = K_local[T_range, u_range] # 4×12 temperature-displacement coupling
```
The range is computed at compile time via `@generated` - zero runtime cost.
## Type Queries
```julia
topology_type(elem) # Tetrahedron{4}
basis_type(elem) # Lagrange{1}
dof_type(elem) # @NamedTuple{T::Tuple{Float64,Vertex}, u::Tuple{Vec{3},Vertex}}
nnodes(elem) # 4
```
## Performance Notes
### Type Stability Achievement
The key to zero allocations was ensuring `elem.dof_indices` has a concrete type:
**Before (BAD):**
```julia
dof_indices::NamedTuple # Type instability!
# Field access returns Any → heap allocation
```
**After (GOOD):**
```julia
@generated function Element{K,P,S}(id::UInt, dof_indices::D) where {K,P,S,D<:NamedTuple}
# Julia infers D = @NamedTuple{u::NTuple{12,Int64}}
# Field access returns NTuple{12,Int64} → stack allocation!
end
```
### Benchmark Results
```text
Flat extraction: 5.472 ns (0 allocations: 0 bytes)
Structured extraction: 5.474 ns (0 allocations: 0 bytes)
```
Compared to original implementation: **300× faster**, zero allocations.
## Files in This Module
- `elements.jl` - Element struct, constructors, type queries
- `extract_element_dofs.jl` - DOF extraction (flat and structured)
- `README.md` - This file (module documentation)
## See Also
- `docs/src/developer/dof_extraction.md` - Detailed DOF extraction design
- `test/elements/test_extract_element_dofs.jl` - Comprehensive test suite
- `src/dofs/README.md` — DOFSet, `DOF{Q, E}` and `DOFHandler`.
- `src/topology/` — `Triangle`, `Tetrahedron`, `Hexahedron`,
topological entities.
- `src/basis/README.md` — basis families and interpolation API.
- `src/assemblers/` — assemblers consume `local_dof_layout`
and the element DOF tuples.
- `test/elements/` — public test suite.