Solver Reference

Each solver card shows its numerical scheme, AD strategy, and Tesseract schema; per-(solver, problem) exclusions and explained anomalies come from the problem configs and are listed alongside. Click Details on any solver to see its field tables. The ∂ column marks fields that support automatic differentiation (VJP/JVP).

Legend. Each card shows a Numerics: line — the discretization followed by the numerical method — and three badge categories: language (backend / runtime), AD: strategy (autodiff / adjoint / hybrid / forward-only), and GPU on GPU-capable solvers.

Checks. Where benchmark results are available, the card footer shows how many of the forward, gradient and optimization checks the solver passes on each setup it runs: H heat transfer, S structural mechanics, F2/F3 Navier–Stokes 2D/3D. Checks the method cannot support are left out. See Results for what each check means.

Navier–Stokes (Grid)

Navier–Stokes (Grid)

Exponax

Numerics: Spectral, ETDRK2

jax AD: autodiff GPU

Exponential time-differencing Runge–Kutta (ETDRK) spectral integrator for 3D Navier–Stokes, enforcing incompressibility to machine precision by construction. Gradients via JAX source-transformation AD

F3 periodic · 3/3 pass Details

Upstream docs: https://fkoehler.site/exponax/

Path: mosaic/tesseracts/navier-stokes-grid/exponax

Image: exponax_navier_stokes_grid

Benchmark checks: F3 Triply periodic box (Taylor–Green): forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
v0	GridVectorField	✓	Initial velocity at the cell centers, shape (nx, ny, 1, 2) in 2-D or (nx, ny, nz, 3) in 3-D. Cell (i, j) is centered at ((i + 1/2) dx, (j + 1/2) dx). Default: 64×64 divergence-free random vortex field (seed=42).
viscosity	Array[(1,), Float32]	✓	Kinematic viscosity ν.
dt	Array[(1,), Float32]	✓	Timestep size.
steps	int		Number of timesteps. Total simulated time = steps × dt.
domain_extent	float | list[float]		Physical size of the domain. A scalar L is a square (cube) of side L, and a list gives the length per axis, [Lx, Ly] or [Lx, Ly, Lz]. Periodic boxes take a scalar. The grid spacing must be the same along every axis.
boundary_conditions	GridBC		Per-face boundary conditions. Default: all periodic.
obstacle	GridObstacle | None		Optional solid obstacle embedded in the domain (no-slip walls), in physical coordinates.
inflow_profile	Array[(None,), Float32] | None	✓	Inlet velocity u_x(y) at the x_lo face, sampled at the cell-center heights y_j = (j + 1/2) dx, shape (ny,). u_y is zero on the inlet. Requires x_lo to be ‘dirichlet’ with no value.
drag	Array[(1,), Float32]	✓	Linear drag coefficient
order	int		ETDRK time integration order (0-4)
kolmogorov_forcing	bool		Enable Kolmogorov sinusoidal body forcing for sustained turbulence
injection_mode	int		Wavenumber at which energy is injected (Kolmogorov forcing)
injection_scale	Array[(1,), Float32]	✓	Amplitude of the Kolmogorov body forcing

Outputs

Field	Type	∂	Description
result	GridVectorField	✓	Final velocity at the cell centers, same shape as v0.
pressure	GridScalarField | None		Final kinematic pressure p/ρ at the cell centers, shape (nx, ny, nz). Values inside the obstacle are ignored.
body_force	Array[(None, None), Float32] | None		Force of the fluid on the obstacle per unit depth and density, shape (steps, ndim): row k is (F_x, F_y) at t = (k + 1) dt. Positive along +x and +y, so a flow in +x gives F_x > 0.

Navier–Stokes (Grid)

INS.jl

Numerics: Finite Volume, Projection, RK4

julia AD: autodiff

Julia staggered finite-volume solver with spectral pressure projection and RK4 time integration; run on CPU in this benchmark (the library itself also supports CUDA via KernelAbstractions.jl). Differentiates through the time loop via Zygote.jl reverse-mode AD.

F2 periodic · 2/2 pass F2 cylinder · 0/2 pass F3 periodic · 3/3 pass Details

Upstream docs: https://agdestein.github.io/IncompressibleNavierStokes.jl/dev/

Path: mosaic/tesseracts/navier-stokes-grid/ins-jl

Image: ins_navier_stokes_grid

Benchmark checks: F2 Periodic box (Taylor–Green, multimode): forward pass, gradient pass; F2 Channel flow past a cylinder: forward not run, optimization not run; F3 Triply periodic box (Taylor–Green): forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
v0	GridVectorField	✓	Initial velocity at the cell centers, shape (nx, ny, 1, 2) in 2-D or (nx, ny, nz, 3) in 3-D. Cell (i, j) is centered at ((i + 1/2) dx, (j + 1/2) dx). Default: 64×64 divergence-free random vortex field (seed=42).
viscosity	Array[(1,), Float32]	✓	Kinematic viscosity ν.
dt	Array[(1,), Float32]	✓	Timestep size.
steps	int		Number of timesteps. Total simulated time = steps × dt.
domain_extent	float | list[float]		Physical size of the domain. A scalar L is a square (cube) of side L, and a list gives the length per axis, [Lx, Ly] or [Lx, Ly, Lz]. Periodic boxes take a scalar. The grid spacing must be the same along every axis.
boundary_conditions	GridBC		Per-face boundary conditions. Default: all periodic.
obstacle	GridObstacle | None		Optional solid obstacle embedded in the domain (no-slip walls), in physical coordinates.
inflow_profile	Array[(None,), Float32] | None		Inlet velocity u_x(y) at the x_lo face, sampled at the cell-center heights y_j = (j + 1/2) dx, shape (ny,). u_y is zero on the inlet. Requires x_lo to be ‘dirichlet’ with no value.

Outputs

Field	Type	∂	Description
result	GridVectorField	✓	Final velocity at the cell centers, same shape as v0.
pressure	GridScalarField | None		Final kinematic pressure p/ρ at the cell centers, shape (nx, ny, nz). Values inside the obstacle are ignored.
body_force	Array[(None, None), Float32] | None		Force of the fluid on the obstacle per unit depth and density, shape (steps, ndim): row k is (F_x, F_y) at t = (k + 1) dt. Positive along +x and +y, so a flow in +x gives F_x > 0.

Navier–Stokes (Grid)

jax-cfd

Numerics: Finite Difference, Projection, RK4

jax AD: autodiff GPU

JAX-native incompressible Navier–Stokes solver on a staggered MAC grid with finite-difference advection and spectral FFT pressure projection. Gradients via JAX source-transformation AD through the time loop.

F2 periodic · 2/2 pass F2 cylinder · 0/2 pass Details

Upstream docs: https://github.com/google/jax-cfd

Path: mosaic/tesseracts/navier-stokes-grid/jax-cfd

Image: jax_cfd_navier_stokes_grid

Benchmark checks: F2 Periodic box (Taylor–Green, multimode): forward pass, gradient pass; F2 Channel flow past a cylinder: forward not run, optimization not run — see Results.

Inputs

Field	Type	∂	Description
v0	GridVectorField	✓	Initial velocity at the cell centers, shape (nx, ny, 1, 2) in 2-D or (nx, ny, nz, 3) in 3-D. Cell (i, j) is centered at ((i + 1/2) dx, (j + 1/2) dx). Default: 64×64 divergence-free random vortex field (seed=42).
viscosity	Array[(1,), Float32]	✓	Kinematic viscosity ν.
dt	Array[(1,), Float32]	✓	Timestep size.
steps	int		Number of timesteps. Total simulated time = steps × dt.
domain_extent	float | list[float]		Physical size of the domain. A scalar L is a square (cube) of side L, and a list gives the length per axis, [Lx, Ly] or [Lx, Ly, Lz]. Periodic boxes take a scalar. The grid spacing must be the same along every axis.
boundary_conditions	GridBC		Per-face boundary conditions. Default: all periodic.
obstacle	GridObstacle | None		Optional solid obstacle embedded in the domain (no-slip walls), in physical coordinates.
inflow_profile	Array[(None,), Float32] | None		Inlet velocity u_x(y) at the x_lo face, sampled at the cell-center heights y_j = (j + 1/2) dx, shape (ny,). u_y is zero on the inlet. Requires x_lo to be ‘dirichlet’ with no value.
density	Array[(1,), Float32]	✓	Density of the fluid
inner_steps	int		Solver sub-iterations per timestep (higher = more accurate pressure solve)

Outputs

Field	Type	∂	Description
result	GridVectorField	✓	Final velocity at the cell centers, same shape as v0.
pressure	GridScalarField | None		Final kinematic pressure p/ρ at the cell centers, shape (nx, ny, nz). Values inside the obstacle are ignored.
body_force	Array[(None, None), Float32] | None		Force of the fluid on the obstacle per unit depth and density, shape (steps, ndim): row k is (F_x, F_y) at t = (k + 1) dt. Positive along +x and +y, so a flow in +x gives F_x > 0.

Navier–Stokes (Grid)

OpenFOAM

Numerics: Finite Volume, PISO/PIMPLE, Euler

cpp AD: forward-only

OpenFOAM icoFoam PISO solver run as a forward-only reference baseline (C++ subprocess). No reverse-mode AD; used to anchor cross-solver comparisons.

F2 periodic · 1/1 pass F2 cylinder · 1/1 pass F3 periodic · 1/1 pass Details

Upstream docs: https://www.openfoam.com/documentation/overview

Path: mosaic/tesseracts/navier-stokes-grid/openfoam

Image: openfoam_navier_stokes_grid

Benchmark checks: F2 Periodic box (Taylor–Green, multimode): forward pass, gradient not supported; F2 Channel flow past a cylinder: forward pass, optimization not supported; F3 Triply periodic box (Taylor–Green): forward pass, gradient not supported, optimization not supported — see Results.

Inputs

Field	Type	∂	Description
v0	GridVectorField		Initial velocity at the cell centers, shape (nx, ny, 1, 2) in 2-D or (nx, ny, nz, 3) in 3-D. Cell (i, j) is centered at ((i + 1/2) dx, (j + 1/2) dx). Default: 64×64 divergence-free random vortex field (seed=42).
viscosity	Array[(1,), Float32]		Kinematic viscosity ν.
dt	Array[(1,), Float32]		Timestep size.
steps	int		Number of timesteps. Total simulated time = steps × dt.
domain_extent	float | list[float]		Physical size of the domain. A scalar L is a square (cube) of side L, and a list gives the length per axis, [Lx, Ly] or [Lx, Ly, Lz]. Periodic boxes take a scalar. The grid spacing must be the same along every axis.
boundary_conditions	GridBC		Per-face boundary conditions. Default: all periodic.
obstacle	GridObstacle | None		Optional solid obstacle embedded in the domain (no-slip walls), in physical coordinates.
inflow_profile	Array[(None,), Float32] | None		Inlet velocity u_x(y) at the x_lo face, sampled at the cell-center heights y_j = (j + 1/2) dx, shape (ny,). u_y is zero on the inlet. Requires x_lo to be ‘dirichlet’ with no value.

Outputs

Field	Type	∂	Description
result	GridVectorField		Final velocity at the cell centers, same shape as v0.
pressure	GridScalarField | None		Final kinematic pressure p/ρ at the cell centers, shape (nx, ny, nz). Values inside the obstacle are ignored.
body_force	Array[(None, None), Float32] | None		Force of the fluid on the obstacle per unit depth and density, shape (steps, ndim): row k is (F_x, F_y) at t = (k + 1) dt. Positive along +x and +y, so a flow in +x gives F_x > 0.

Navier–Stokes (Grid)

PhiFlow

Numerics: Finite Difference, 2nd-order central FD (advective form), explicit Euler, CG pressure projection

jax AD: hybrid GPU

2nd-order central-difference advection (explicit Euler) with CG pressure projection. The cylinder channel uses phiflow’s wind-tunnel boundaries and a staircase cylinder. Gradients via JAX autodiff with adjoint CG solves.

F2 periodic · 2/2 pass F2 cylinder · 2/2 pass F3 periodic · 3/3 pass Details

Upstream docs: https://tum-pbs.github.io/PhiFlow/

Path: mosaic/tesseracts/navier-stokes-grid/phiflow

Image: phiflow_navier_stokes_grid

Benchmark checks: F2 Periodic box (Taylor–Green, multimode): forward pass, gradient pass; F2 Channel flow past a cylinder: forward pass, optimization pass; F3 Triply periodic box (Taylor–Green): forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
v0	GridVectorField	✓	Initial velocity at the cell centers, shape (nx, ny, 1, 2) in 2-D or (nx, ny, nz, 3) in 3-D. Cell (i, j) is centered at ((i + 1/2) dx, (j + 1/2) dx). Default: 64×64 divergence-free random vortex field (seed=42).
viscosity	Array[(1,), Float32]	✓	Kinematic viscosity ν.
dt	Array[(1,), Float32]	✓	Timestep size.
steps	int		Number of timesteps. Total simulated time = steps × dt.
domain_extent	float | list[float]		Physical size of the domain. A scalar L is a square (cube) of side L, and a list gives the length per axis, [Lx, Ly] or [Lx, Ly, Lz]. Periodic boxes take a scalar. The grid spacing must be the same along every axis.
boundary_conditions	GridBC		Per-face boundary conditions. Default: all periodic.
obstacle	GridObstacle | None		Optional solid obstacle embedded in the domain (no-slip walls), in physical coordinates.
inflow_profile	Array[(None,), Float32] | None	✓	Inlet velocity u_x(y) at the x_lo face, sampled at the cell-center heights y_j = (j + 1/2) dx, shape (ny,). u_y is zero on the inlet. Requires x_lo to be ‘dirichlet’ with no value.

Outputs

Field	Type	∂	Description
result	GridVectorField	✓	Final velocity at the cell centers, same shape as v0.
pressure	GridScalarField | None		Final kinematic pressure p/ρ at the cell centers, shape (nx, ny, nz). Values inside the obstacle are ignored.
body_force	Array[(None, None), Float32] | None	✓	Force of the fluid on the obstacle per unit depth and density, shape (steps, ndim): row k is (F_x, F_y) at t = (k + 1) dt. Positive along +x and +y, so a flow in +x gives F_x > 0.

Navier–Stokes (Grid)

PICT

Numerics: Finite Volume, PISO, BDF1

pytorch AD: autodiff GPU

GPU-accelerated differentiable PISO solver built on PyTorch with custom CUDA kernels (PISOtorch). VJP w.r.t. v0, inflow_profile, and viscosity flows through the PISO time loop via PISOtorch_diff (tracks VELOCITY, BOUNDARY_VELOCITY, and VISCOSITY gradients). dt is a Python scalar consumed at construction time; zero gradients are returned for it.

F2 periodic · 1/2 pass F2 cylinder · 1/2 pass F3 periodic · 3/3 pass Details

Upstream docs: https://github.com/tum-pbs/PICT

Path: mosaic/tesseracts/navier-stokes-grid/pict

Image: pict_navier_stokes_grid

Benchmark checks: F2 Periodic box (Taylor–Green, multimode): forward borderline, gradient pass; F2 Channel flow past a cylinder: forward pass, optimization not run; F3 Triply periodic box (Taylor–Green): forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
v0	GridVectorField	✓	Initial velocity at the cell centers, shape (nx, ny, 1, 2) in 2-D or (nx, ny, nz, 3) in 3-D. Cell (i, j) is centered at ((i + 1/2) dx, (j + 1/2) dx). Default: 64×64 divergence-free random vortex field (seed=42).
viscosity	Array[(1,), Float32]	✓	Kinematic viscosity ν.
dt	Array[(1,), Float32]		Timestep size.
steps	int		Number of timesteps. Total simulated time = steps × dt.
domain_extent	float | list[float]		Physical size of the domain. A scalar L is a square (cube) of side L, and a list gives the length per axis, [Lx, Ly] or [Lx, Ly, Lz]. Periodic boxes take a scalar. The grid spacing must be the same along every axis.
boundary_conditions	GridBC		Per-face boundary conditions. Default: all periodic.
obstacle	GridObstacle | None		Optional solid obstacle embedded in the domain (no-slip walls), in physical coordinates.
inflow_profile	Array[(None,), Float32] | None	✓	Inlet velocity u_x(y) at the x_lo face, sampled at the cell-center heights y_j = (j + 1/2) dx, shape (ny,). u_y is zero on the inlet. Requires x_lo to be ‘dirichlet’ with no value.

Outputs

Field	Type	∂	Description
result	GridVectorField	✓	Final velocity at the cell centers, same shape as v0.
pressure	GridScalarField | None		Final kinematic pressure p/ρ at the cell centers, shape (nx, ny, nz). Values inside the obstacle are ignored.
body_force	Array[(None, None), Float32] | None	✓	Force of the fluid on the obstacle per unit depth and density, shape (steps, ndim): row k is (F_x, F_y) at t = (k + 1) dt. Positive along +x and +y, so a flow in +x gives F_x > 0.

Navier–Stokes (Grid)

Warp-NS

Numerics: Finite Difference, Projection, Euler

warp AD: autodiff GPU

Chorin projection scheme with explicit Euler time stepping, implemented in NVIDIA Warp CUDA kernels with FFT pressure Poisson on periodic domains. The cylinder channel uses a staggered grid, a staircase cylinder and a conjugate-gradient pressure Poisson solve. Gradients via wp.Tape kernel-level VJP, with adjoint CG solves.

F2 periodic · 2/2 pass F2 cylinder · 2/2 pass F3 periodic · 3/3 pass Details

Upstream docs: https://github.com/NVIDIA/warp

Path: mosaic/tesseracts/navier-stokes-grid/warp-ns

Image: warp_ns_navier_stokes_grid

Benchmark checks: F2 Periodic box (Taylor–Green, multimode): forward pass, gradient pass; F2 Channel flow past a cylinder: forward pass, optimization pass; F3 Triply periodic box (Taylor–Green): forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
v0	GridVectorField	✓	Initial velocity at the cell centers, shape (nx, ny, 1, 2) in 2-D or (nx, ny, nz, 3) in 3-D. Cell (i, j) is centered at ((i + 1/2) dx, (j + 1/2) dx). Default: 64×64 divergence-free random vortex field (seed=42).
viscosity	Array[(1,), Float32]	✓	Kinematic viscosity ν.
dt	Array[(1,), Float32]	✓	Timestep size.
steps	int		Number of timesteps. Total simulated time = steps × dt.
domain_extent	float | list[float]		Physical size of the domain. A scalar L is a square (cube) of side L, and a list gives the length per axis, [Lx, Ly] or [Lx, Ly, Lz]. Periodic boxes take a scalar. The grid spacing must be the same along every axis.
boundary_conditions	GridBC		Per-face boundary conditions. Default: all periodic.
obstacle	GridObstacle | None		Optional solid obstacle embedded in the domain (no-slip walls), in physical coordinates.
inflow_profile	Array[(None,), Float32] | None	✓	Inlet velocity u_x(y) at the x_lo face, sampled at the cell-center heights y_j = (j + 1/2) dx, shape (ny,). u_y is zero on the inlet. Requires x_lo to be ‘dirichlet’ with no value.

Outputs

Field	Type	∂	Description
result	GridVectorField	✓	Final velocity at the cell centers, same shape as v0.
pressure	GridScalarField | None		Final kinematic pressure p/ρ at the cell centers, shape (nx, ny, nz). Values inside the obstacle are ignored.
body_force	Array[(None, None), Float32] | None	✓	Force of the fluid on the obstacle per unit depth and density, shape (steps, ndim): row k is (F_x, F_y) at t = (k + 1) dt. Positive along +x and +y, so a flow in +x gives F_x > 0.

Navier–Stokes (Grid)

XLB

Numerics: Lattice Boltzmann, Streaming

jax AD: autodiff GPU

JAX-accelerated lattice Boltzmann solver (D2Q9 / D3Q27) with BGK or KBC collision, running the streaming/collision loop on GPU. Gradients via source-transformation AD; recovers incompressible Navier–Stokes only to O(Ma²) via Chapman–Enskog expansion.

F2 periodic · 2/2 pass F2 cylinder · 2/2 pass F3 periodic · 3/3 pass Details

Upstream docs: https://github.com/Autodesk/XLB

Path: mosaic/tesseracts/navier-stokes-grid/xlb

Image: xlb_navier_stokes_grid

Benchmark checks: F2 Periodic box (Taylor–Green, multimode): forward pass, gradient pass; F2 Channel flow past a cylinder: forward pass, optimization pass; F3 Triply periodic box (Taylor–Green): forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
v0	GridVectorField	✓	Initial velocity at the cell centers, shape (nx, ny, 1, 2) in 2-D or (nx, ny, nz, 3) in 3-D. Cell (i, j) is centered at ((i + 1/2) dx, (j + 1/2) dx). Default: 64×64 divergence-free random vortex field (seed=42).
viscosity	Array[(1,), Float32]	✓	Kinematic viscosity ν.
dt	Array[(1,), Float32]	✓	Timestep size.
steps	int		Number of timesteps. Total simulated time = steps × dt.
domain_extent	float | list[float]		Physical size of the domain. A scalar L is a square (cube) of side L, and a list gives the length per axis, [Lx, Ly] or [Lx, Ly, Lz]. Periodic boxes take a scalar. The grid spacing must be the same along every axis.
boundary_conditions	GridBC		Per-face boundary conditions. Default: all periodic.
obstacle	GridObstacle | None		Optional solid obstacle embedded in the domain (no-slip walls), in physical coordinates.
inflow_profile	Array[(None,), Float32] | None	✓	Inlet velocity u_x(y) at the x_lo face, sampled at the cell-center heights y_j = (j + 1/2) dx, shape (ny,). u_y is zero on the inlet. Requires x_lo to be ‘dirichlet’ with no value.

Outputs

Field	Type	∂	Description
result	GridVectorField	✓	Final velocity at the cell centers, same shape as v0.
pressure	GridScalarField | None		Final kinematic pressure p/ρ at the cell centers, shape (nx, ny, nz). Values inside the obstacle are ignored.
body_force	Array[(None, None), Float32] | None	✓	Force of the fluid on the obstacle per unit depth and density, shape (steps, ndim): row k is (F_x, F_y) at t = (k + 1) dt. Positive along +x and +y, so a flow in +x gives F_x > 0.

Structural mechanics

Structural mechanics

deal.II

Numerics: Finite Element, Direct (UMFPACK)

cpp AD: forward-only

deal.II Q1 trilinear hexahedral FEM linear-elasticity solver via C++ subprocess (UMFPACK direct). Forward-only — reference implementation for cross-framework validation.

S cantilever · 1/1 pass Details

Upstream docs: https://www.dealii.org/

Path: mosaic/tesseracts/structural-mesh/dealii-structural

Image: dealii_structural_mesh

Benchmark checks: S Linear-elastic cantilever: forward pass, gradient not supported, optimization not supported — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]		Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. Values in [0, 1]: 0 = void, 1 = solid.
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe displacement (values shape (n_groups, 3), None → zero); neumann.mask/values prescribe surface traction (values shape (n_groups, 3)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.
E_max	float		Young’s modulus of the fully solid material [MPa].
nu	float		Poisson’s ratio.
xmin	float		Void stiffness ratio (E_min = xmin * E_max).
penal	float		SIMP penalisation exponent (E(ρ) = E_min + (E_max−E_min)·ρ^penal).

Outputs

Field	Type	∂	Description
compliance	Array[(), Float32]		Structural compliance C = F^T U (scalar). Minimising compliance is equivalent to maximising stiffness.

Structural mechanics

FEniCS

Numerics: Finite Element, Direct (MUMPS)

fenics AD: adjoint

FEniCS/DOLFIN CG1 FEM linear-elasticity solver with SIMP penalisation. Gradient ∂C/∂ρ via dolfin-adjoint tape replay (discrete adjoint derived at the variational level via UFL).

S cantilever · 3/3 pass Details

Upstream docs: https://fenicsproject.org/

Path: mosaic/tesseracts/structural-mesh/fenics-structural

Image: fenics_structural_mesh

Benchmark checks: S Linear-elastic cantilever: forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]	✓	Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. Values in [0, 1]: 0 = void, 1 = solid.
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe displacement (values shape (n_groups, 3), None → zero); neumann.mask/values prescribe surface traction (values shape (n_groups, 3)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.
E_max	float		Young’s modulus of the fully solid material.
nu	float		Poisson’s ratio.
xmin	float		Void stiffness ratio (E_min = xmin * E_max).
penal	float		SIMP penalisation exponent (E(rho) = E_min + (E_max-E_min)*rho^penal).

Outputs

Field	Type	∂	Description
compliance	Array[(), Float32]	✓	Structural compliance C = F^T U (scalar). Minimising compliance is equivalent to maximising stiffness.

Structural mechanics

Firedrake

Numerics: Finite Element, Direct (MUMPS)

firedrake AD: adjoint

Firedrake CG1 FEM linear-elasticity solver with SIMP penalisation. Gradient ∂C/∂ρ via firedrake-adjoint (pyadjoint ReducedFunctional) — an independent tape-based discrete adjoint complement to FEniCS.

S cantilever · 3/3 pass Details

Upstream docs: https://www.firedrakeproject.org/

Path: mosaic/tesseracts/structural-mesh/firedrake-structural

Image: firedrake_structural_mesh

Benchmark checks: S Linear-elastic cantilever: forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]	✓	Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. Values in [0, 1]: 0 = void, 1 = solid.
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe displacement (values shape (n_groups, 3), None → zero); neumann.mask/values prescribe surface traction (values shape (n_groups, 3)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.
E_max	float		Young’s modulus of the fully solid material.
nu	float		Poisson’s ratio.
xmin	float		Void stiffness ratio (E_min = xmin * E_max).
penal	float		SIMP penalisation exponent (E(rho) = E_min + (E_max-E_min)*rho^penal).

Outputs

Field	Type	∂	Description
compliance	Array[(), Float32]	✓	Structural compliance C = F^T U (scalar). Minimising compliance is equivalent to maximising stiffness.

Structural mechanics

JAX-FEM

Numerics: Finite Element, Direct (SuperLU)

jax AD: hybrid GPU

JAX-FEM HEX8 linear-elasticity solver with SIMP stiffness E(ρ) = E_min + (E_max − E_min)·ρ^penal. Gradient ∂C/∂ρ via JAX automatic differentiation through the linear FE solve (hybrid implicit-function-theorem adjoint).

S cantilever · 3/3 pass Details

Upstream docs: https://github.com/deepmodeling/jax-fem

Path: mosaic/tesseracts/structural-mesh/jax-fem

Image: jax_fem_structural_mesh

Benchmark checks: S Linear-elastic cantilever: forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]	✓	Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. Values in [0, 1]: 0 = void, 1 = solid.
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe displacement (values shape (n_groups, 3), None → zero); neumann.mask/values prescribe surface traction (values shape (n_groups, 3)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.

Outputs

Field	Type	∂	Description
compliance	Array[(), Float32]	✓	Structural compliance C = F^T U (scalar). Minimising compliance is equivalent to maximising stiffness.

Structural mechanics

TopOpt.jl

Numerics: Finite Element, Direct (CHOLMOD)

julia AD: adjoint

TopOpt.jl HEX8 FEM linear-elasticity solver with SIMP penalisation, Julia + Ferrite.jl backend. Gradient ∂C/∂ρ via analytical SIMP adjoint sensitivity computed alongside the forward solve.

S cantilever · 3/3 pass Details

Upstream docs: https://github.com/JuliaTopOpt/TopOpt.jl

Path: mosaic/tesseracts/structural-mesh/topopt-jl

Image: topopt_jl_structural_grid

Benchmark checks: S Linear-elastic cantilever: forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]	✓	Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. Values in [0, 1]: 0 = void, 1 = solid.
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe displacement (values shape (n_groups, 3), None → zero); neumann.mask/values prescribe surface traction (values shape (n_groups, 3)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.
E	float		Young’s modulus of the solid material.
nu	float		Poisson’s ratio.
xmin	float		Minimum density (void stiffness) to prevent singular stiffness matrix.

Outputs

Field	Type	∂	Description
compliance	Array[(), Float32]	✓	Structural compliance C = F^T U (scalar). Minimising compliance is equivalent to maximising stiffness.

Structural mechanics

torch-fem

Numerics: Finite Element, CG+Jacobi

pytorch AD: hybrid GPU

torch-fem Solid HEX8 linear-elasticity solver with SIMP stiffness E(ρ) = E_min + (E_max − E_min)·ρ^p — the PyTorch backend in structural-mesh. Gradient ∂C/∂ρ via PyTorch autograd through the FE solve (implicit-function-theorem adjoint); GPU-native.

S cantilever · 3/3 pass Details

Upstream docs: https://github.com/meyer-nils/torch-fem

Path: mosaic/tesseracts/structural-mesh/torch-fem-structural

Image: torch_fem_structural_mesh

Benchmark checks: S Linear-elastic cantilever: forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]	✓	Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. Values in [0, 1]: 0 = void, 1 = solid.
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe displacement (values shape (n_groups, 3), None → zero); neumann.mask/values prescribe surface traction (values shape (n_groups, 3)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.
E_max	float		Young’s modulus of fully solid material (ρ = 1).
nu	float		Poisson’s ratio (density-independent).
xmin	float		Void stiffness ratio E_min / E_max.
p_exp	float		SIMP penalisation exponent p in E(ρ) = E_min + (E_max−E_min)·ρ^p.

Outputs

Field	Type	∂	Description
compliance	Array[(), Float32]	✓	Structural compliance C = F^T U (scalar). Minimising compliance is equivalent to maximising stiffness.

Heat transfer

Heat transfer

deal.II

Numerics: Finite Element, CG+SSOR

cpp AD: forward-only

deal.II Q1 trilinear hexahedral FEM heat-conduction solver with SIMP conductivity (C++ subprocess, SSOR-preconditioned CG). Forward-only — reference implementation for cross-framework validation.

H slab · 1/1 pass Details

Upstream docs: https://www.dealii.org/

Path: mosaic/tesseracts/thermal-mesh/dealii-heat

Image: dealii_heat_thermal_mesh

Benchmark checks: H Steady conduction in a slab: forward pass, gradient not supported, optimization not supported — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]		Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. ρ=1 → fully conducting, ρ=0 → insulating.
source	Array[(None,), Float32]		Per-element volumetric heat source (W/m³), shape (n_cells,). Active slice = source[:hex_mesh.n_faces]. Zero by default (topology-optimisation mode; no body heat source).
target_temperature	Array[(None,), Float32]		Per-node target temperature for the inverse (source-identification) objective, shape (n_nodes,). identification_error =
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe temperature (values shape (n_groups, 1), None → zero); neumann.mask/values prescribe heat flux (values shape (n_groups, 1)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.
k_max	float		Maximum thermal conductivity (fully solid material).
p_exp	float		SIMP penalisation exponent p (k(ρ) = k_min + (k_max−k_min)·ρ^p).

Outputs

Field	Type	∂	Description
thermal_compliance	Array[(), Float32]		Total thermal compliance C = ∮_Γ_N q_n · T dΓ (scalar). Equals the work done by the Neumann heat flux on the temperature field.
identification_error	Array[(), Float32]		Source-identification objective:

Heat transfer

FEniCS

Numerics: Finite Element, Direct (MUMPS)

fenics AD: adjoint

FEniCS P1 FEM heat-conduction solver with SIMP penalisation. Gradients (∂C/∂ρ, ∂C/∂source) via dolfin-adjoint tape replay.

H slab · 3/3 pass Details

Upstream docs: https://fenicsproject.org/

Path: mosaic/tesseracts/thermal-mesh/fenics-heat

Image: fenics_heat_thermal_mesh

Benchmark checks: H Steady conduction in a slab: forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]	✓	Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. ρ=1 → fully conducting, ρ=0 → insulating.
source	Array[(None,), Float32]	✓	Per-element volumetric heat source (W/m³), shape (n_cells,). Active slice = source[:hex_mesh.n_faces]. Zero by default (topology-optimisation mode; no body heat source).
target_temperature	Array[(None,), Float32]		Per-node target temperature for the inverse (source-identification) objective, shape (n_nodes,). identification_error =
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe temperature (values shape (n_groups, 1), None → zero); neumann.mask/values prescribe heat flux (values shape (n_groups, 1)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.
k_max	float		Maximum thermal conductivity (fully solid/conducting material).
p_exp	float		SIMP penalisation exponent p (k(ρ) = k_min + (k_max−k_min)·ρ^p).

Outputs

Field	Type	∂	Description
thermal_compliance	Array[(), Float32]	✓	Total thermal compliance C = ∮_Γ_N q_n · T dΓ (scalar). Equals the work done by the Neumann heat flux on the temperature field.
identification_error	Array[(), Float32]	✓	Source-identification objective:

Heat transfer

Firedrake

Numerics: Finite Element, Direct (MUMPS)

firedrake AD: adjoint

Firedrake P1 FEM heat-conduction solver with SIMP penalisation. Gradients (∂C/∂ρ, ∂C/∂source) via firedrake-adjoint, complementing the FEniCS reference.

H slab · 3/3 pass Details

Upstream docs: https://www.firedrakeproject.org/

Path: mosaic/tesseracts/thermal-mesh/firedrake-heat

Image: firedrake_heat_thermal_mesh

Benchmark checks: H Steady conduction in a slab: forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]	✓	Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. ρ=1 → fully conducting, ρ=0 → insulating.
source	Array[(None,), Float32]	✓	Per-element volumetric heat source (W/m³), shape (n_cells,). Active slice = source[:hex_mesh.n_faces]. Zero by default (topology-optimisation mode; no body heat source).
target_temperature	Array[(None,), Float32]		Per-node target temperature for the inverse (source-identification) objective, shape (n_nodes,). identification_error =
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe temperature (values shape (n_groups, 1), None → zero); neumann.mask/values prescribe heat flux (values shape (n_groups, 1)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.
k_max	float		Maximum thermal conductivity (fully solid/conducting material).
p_exp	float		SIMP penalisation exponent p (k(ρ) = k_min + (k_max−k_min)·ρ^p).

Outputs

Field	Type	∂	Description
thermal_compliance	Array[(), Float32]	✓	Total thermal compliance C = ∮_Γ_N q_n · T dΓ (scalar). Equals the work done by the Neumann heat flux on the temperature field.
identification_error	Array[(), Float32]	✓	Source-identification objective:

Heat transfer

JAX-FEM

Numerics: Finite Element, Direct (SuperLU)

jax AD: hybrid GPU

JAX-FEM HEX8 steady-state heat-conduction solver with SIMP conductivity and optional volumetric source. Gradients ∂C/∂ρ and ∂C/∂source via JAX autodiff through the linear FE solve (hybrid implicit-function-theorem adjoint)

H slab · 3/3 pass Details

Upstream docs: https://github.com/deepmodeling/jax-fem

Path: mosaic/tesseracts/thermal-mesh/jax-fem

Image: jax_fem_thermal_mesh

Benchmark checks: H Steady conduction in a slab: forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]	✓	Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. ρ=1 → fully conducting, ρ=0 → insulating.
source	Array[(None,), Float32]	✓	Per-element volumetric heat source (W/m³), shape (n_cells,). Active slice = source[:hex_mesh.n_faces]. Zero by default (topology-optimisation mode; no body heat source).
target_temperature	Array[(None,), Float32]		Per-node target temperature for the inverse (source-identification) objective, shape (n_nodes,). identification_error =
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe temperature (values shape (n_groups, 1), None → zero); neumann.mask/values prescribe heat flux (values shape (n_groups, 1)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.

Outputs

Field	Type	∂	Description
thermal_compliance	Array[(), Float32]	✓	Total thermal compliance C = ∮_Γ_N q_n · T dΓ (scalar). Equals the work done by the Neumann heat flux on the temperature field.
identification_error	Array[(), Float32]	✓	Source-identification objective:

Heat transfer

torch-fem

Numerics: Finite Element, CG+Jacobi

pytorch AD: hybrid GPU

torch-fem SolidHeat HEX8 heat-conduction solver with IsotropicConductivity3D SIMP material — the first PyTorch backend in thermal-mesh. Gradients via PyTorch autograd through the FE solve; GPU-native.

H slab · 3/3 pass Details

Upstream docs: https://github.com/meyer-nils/torch-fem

Path: mosaic/tesseracts/thermal-mesh/torch-fem-thermal

Image: torch_fem_thermal_mesh

Benchmark checks: H Steady conduction in a slab: forward pass, gradient pass, optimization pass — see Results.

Inputs

Field	Type	∂	Description
rho	Array[(None,), Float32]	✓	Per-cell density field, shape (n_max_cells,). Active slice = rho[:hex_mesh.n_faces]. ρ=1 → fully conducting, ρ=0 → insulating.
source	Array[(None,), Float32]	✓	Per-element volumetric heat source (W/m³), shape (n_cells,). Active slice = source[:hex_mesh.n_faces]. Zero by default (topology-optimisation mode; no body heat source).
target_temperature	Array[(None,), Float32]		Per-node target temperature for the inverse (source-identification) objective, shape (n_nodes,). identification_error =
boundary_conditions	MeshBC		Mesh boundary conditions. dirichlet.mask/values prescribe temperature (values shape (n_groups, 1), None → zero); neumann.mask/values prescribe heat flux (values shape (n_groups, 1)).
hex_mesh	HexMesh		Hexahedral mesh (8-node hex elements). n_points / n_faces give the active slice sizes.
k_max	float		Maximum thermal conductivity (fully solid material).
p_exp	float		SIMP penalisation exponent p (k(ρ) = k_min + (k_max−k_min)·ρ^p).

Outputs

Field	Type	∂	Description
thermal_compliance	Array[(), Float32]	✓	Total thermal compliance C = ∮_Γ_N q_n · T dΓ (scalar). Equals the work done by the Neumann heat flux on the temperature field.
identification_error	Array[(), Float32]	✓	Source-identification objective: