API Reference

Domain

Macchiato.Domain Type

Domain{M, C}

Central container that ties together a point cloud, boundary conditions, and physics models.

Fields

• cloud::PointCloud{M, C}: The discretized geometry (boundary + interior points)

• boundaries::Dict{Symbol, Tuple{UnitRange, AbstractBoundaryCondition}}: Mapping from surface names to (index_range, bc) pairs, where index_range gives the global indices of that surface's points in the assembled system

• models::AbstractVector{<:AbstractModel}: Physics models (e.g., SolidEnergy, LinearElasticity)

• name::Symbol: Domain identifier (defaults to :domain1)

Constructors

Domain(cloud, boundaries, model)    # cloud + BCs + model(s)

Boundary conditions are required: a physics model on a domain with no boundary conditions is ill-posed. At construction, the Domain validates that: 2. Every boundary condition key matches a surface name in the point cloud

3. Every surface in the point cloud has a corresponding boundary condition entry

source

Macchiato.add! Function

add!(domain::Domain, model::AbstractModel)

Append a physics model to the domain's model list.

source

add!(domain::Domain, boundary::AbstractBoundaryCondition, name::Symbol)

Replace the boundary condition on the named surface, preserving that surface's index range. The surface must already exist in the domain (boundary conditions are assigned at construction).

source

Base.delete! Method

delete!(domain::Domain, model::AbstractModel)

Remove a physics model from the domain.

source

Models

Energy

Macchiato.SolidEnergy Type

SolidEnergy(; k, ρ, cₚ, source=nothing)

Solid-body energy (heat) transport model.

Solves the heat equation in a solid medium:

• Steady-state: k ∇²T = -f (Poisson equation)

• Transient: ρ cₚ ∂T/∂t = k ∇²T + f

Fields

• k: Thermal conductivity

• ρ: Density

• cₚ: Specific heat capacity

• source: Optional volumetric source term f(x, t) -> value (default: nothing)

Example

model = SolidEnergy(k=50.0, ρ=7800.0, cₚ=500.0)
model = SolidEnergy(k=1.0, ρ=1.0, cₚ=1.0, source=(x, t) -> -4.0)

source

Mechanics

Macchiato.LinearElasticity Type

LinearElasticity{E, Nu, Rho, F} <: Solid

Linear isotropic elasticity model for solid mechanics (Navier-Cauchy equations).

Supports 2D plane stress formulation. The governing equations in displacement form:

(λ*+2μ) ∂²u/∂x² + μ ∂²u/∂y² + (λ*+μ) ∂²v/∂x∂y + fₓ = 0
(λ*+μ) ∂²u/∂x∂y + μ ∂²v/∂x² + (λ*+2μ) ∂²v/∂y² + fᵧ = 0

Fields

• E: Young's modulus

• ν: Poisson's ratio

• ρ: Density (optional, for body forces / dynamics)

• body_force: Body force function f(x) -> (fx, fy) (optional)

Example

model = LinearElasticity(E=200e3, ν=0.3)
model = LinearElasticity(E=200e3, ν=0.3, body_force=x -> (0.0, -9.81))

source

Macchiato.lame_parameters Function

lame_parameters(model::LinearElasticity)

Compute Lamé parameters for plane stress:

• μ = E / (2(1+ν))

• λ = Eν / ((1+ν)(1-2ν))

• λ* = 2μλ / (λ+2μ) (plane stress modification)

Returns (μ, λstar).

source

Fluids

Macchiato.IncompressibleNavierStokes Type

IncompressibleNavierStokes(; μ, ρ)
IncompressibleNavierStokes(μ::Real, ρ)

Incompressible Navier-Stokes fluid model.

When μ is a plain number, it is automatically wrapped in NewtonianViscosity. Pass an AbstractViscosity subtype (e.g., CarreauYasudaViscosity) for non-Newtonian behavior.

Warning

The fluid solver is under active development and not yet fully functional.

Fields

• μ::AbstractViscosity: Dynamic viscosity model

• ρ: Fluid density

Example

model = IncompressibleNavierStokes(μ=0.001, ρ=1000.0)
model = IncompressibleNavierStokes(μ=CarreauYasudaViscosity(μ_inf=0.0035, μ_0=0.056), ρ=1060.0)

source

Macchiato.AbstractViscosity Type

AbstractViscosity

Abstract supertype for viscosity models used with IncompressibleNavierStokes. Concrete subtypes must be callable as μ(γ̇) returning the dynamic viscosity at shear rate γ̇.

source

Macchiato.NewtonianViscosity Type

NewtonianViscosity(μ)

Constant (Newtonian) viscosity model. Returns the same viscosity μ regardless of shear rate.

Example

visc = NewtonianViscosity(0.001)  # water-like viscosity
visc(100.0)  # => 0.001

source

Macchiato.CarreauYasudaViscosity Type

CarreauYasudaViscosity(; μ_inf, μ_0, n=0.333, λ=0.31, a=2)

Generalized Newtonian viscosity model for shear-thinning (or shear-thickening) fluids.

The Carreau-Yasuda model:

μ(γ̇) = μ_inf + (μ_0 - μ_inf) * (1 + (λ γ̇)^a)^((n - 1) / a)

Fields

• μ_inf: Infinite-shear-rate viscosity

• μ_0: Zero-shear-rate viscosity

• n: Power-law index (< 1 for shear-thinning)

• λ: Relaxation time

• a: Yasuda parameter (a = 2 recovers the Carreau model)

Example

visc = CarreauYasudaViscosity(μ_inf=0.0035, μ_0=0.056, n=0.333, λ=0.31)

source

Simulation Modes

Macchiato.AbstractSimulationMode Type

AbstractSimulationMode

Abstract supertype for simulation modes. See Steady and Transient.

source

Macchiato.Steady Type

Steady()

Mode for steady-state simulations solved via LinearSolve.

source

Macchiato.Transient Type

Transient(; Δt, stop_time, solver=FBDF())

Mode for transient (time-dependent) simulations solved via OrdinaryDiffEq.solve.

Arguments

• Δt: Time step size

• stop_time: End time for the simulation

• solver: ODE solver. Defaults to FBDF(), an implicit stiff solver — semidiscrete diffusion is stiff, so explicit methods need prohibitively small steps. Pass an explicit solver (e.g. Tsit5()) only if you accept the Δt ≲ Δx²/α CFL limit.

Examples

Transient(Δt=0.001, stop_time=1.0)
Transient(Δt=0.001, stop_time=1.0, solver=Tsit5())

source

Model Interface

These are the functions to implement when defining a custom PDE.

Macchiato.AbstractModel Type

AbstractModel

Abstract supertype for all PDE models. Subtype this to define a custom PDE.

See Custom PDEs for a complete walkthrough.

source

Macchiato.num_vars Function

num_vars(model::AbstractModel, dim) -> Int

Return the number of solution variables per point for model in dim dimensions.

Examples: 1 for scalar PDEs, dim for vector PDEs, dim + 1 for velocity + pressure.

source

Macchiato.node_coordinates Function

node_coordinates(domain::Domain; strip_units=true)
node_coordinates(cloud; strip_units=true)

Coordinates of every node in the domain (boundary first, then interior) as a Vector of SVectors — with units stripped by default, which is the form the RadialBasisFunctions.jl operator constructors (laplacian, partial, …) accept.

Pass strip_units=false to keep the Unitful quantities.

source

Macchiato.make_system Function

make_system(model::AbstractModel, domain; kwargs...) -> (A, b)

Assemble the system matrix A and right-hand side b for steady-state solving.

Required for steady-state simulations. Macchiato applies boundary conditions and solves Ax = b.

source

Macchiato.make_f Function

make_f(model::AbstractModel, domain; kwargs...) -> f

Return an in-place ODE function f(du, u, p, t) for transient integration.

Required for transient simulations. Macchiato passes the returned function to OrdinaryDiffEq.jl.

source

Boundary Conditions

Core Types

Macchiato.AbstractBoundaryCondition Type

AbstractBoundaryCondition

Base abstract type for all boundary conditions. All BCs must subtype one of: Dirichlet, Neumann, or Robin.

source

Macchiato.Dirichlet Type

Dirichlet <: AbstractBoundaryCondition

Essential boundary conditions that prescribe values at the boundary. Examples: Temperature, VelocityInlet, Displacement.

source

Macchiato.Neumann Type

Neumann <: DerivativeBoundaryCondition

Natural boundary conditions that prescribe normal derivatives: ∂u/∂n = g. Examples: HeatFlux, Adiabatic, VelocityOutlet.

source

Macchiato.Robin Type

Robin <: DerivativeBoundaryCondition

Mixed boundary conditions combining value and derivative: β·∂u/∂n + α·u = g. Example: Convection.

source

Macchiato.DerivativeBoundaryCondition Type

DerivativeBoundaryCondition <: AbstractBoundaryCondition

Abstract type for BCs involving derivatives (Neumann and Robin).

source

Generic BC Types

Macchiato.PrescribedValue Type

PrescribedValue{F<:Function} <: Dirichlet

Generic Dirichlet BC that prescribes a value via a function.

The function has signature f(x, t) -> value where:

• x: spatial coordinate of the boundary point

• t: time

• Returns the prescribed value at that location and time

Fields

• f: Function with signature (x, t) -> value

• name: Display name (e.g., :Temperature, :PrescribedValue)

For built-in physics, use named constructors (e.g., Temperature, VelocityInlet). For custom PDEs, use the unparameterized constructors directly: PrescribedValue(0.0).

source

Macchiato.PrescribedFlux Type

PrescribedFlux{F<:Function} <: Neumann

Generic Neumann BC that prescribes a flux (normal derivative) via a function.

The flux condition is: ∂u/∂n = f(x, t)

The function has signature f(x, t) -> flux_value where:

• x: spatial coordinate of the boundary point

• t: time

• Returns the prescribed flux value at that location and time

Fields

• f: Function with signature (x, t) -> flux

• name: Display name (e.g., :HeatFlux, :PrescribedFlux)

For built-in physics, use named constructors (e.g., HeatFlux, Traction). For custom PDEs, use the unparameterized constructors directly: PrescribedFlux(1.0).

source

Macchiato.ZeroFlux Type

ZeroFlux <: Neumann

Generic Neumann BC with zero flux: ∂u/∂n = 0.

Represents symmetry, insulation, or fully-developed flow depending on context.

Fields

• name: Display name (e.g., :Adiabatic, :VelocityOutlet, :ZeroFlux)

For built-in physics, use named constructors (e.g., Adiabatic, VelocityOutlet). For custom PDEs, use the unparameterized constructor directly: ZeroFlux().

source

Note

These generic types work with any physics model. For custom PDEs, use them directly: PrescribedValue(0.0), PrescribedFlux(1.0), ZeroFlux(). See Custom PDEs for a complete example.

Energy BCs

Macchiato.Temperature Function

Temperature(value)

Prescribed temperature BC. Value can be a Number or Function (x, t) -> value.

source

Macchiato.HeatFlux Function

HeatFlux(flux)

Prescribed heat flux BC: ∂T/∂n = q. Flux can be a Number or Function (x, t) -> flux.

source

Macchiato.Adiabatic Function

Adiabatic()

Thermally insulated boundary: ∂T/∂n = 0.

source

Macchiato.Convection Type

Convection(h, k, T∞)

Convective heat transfer: h·T + k·∂T/∂n = h·T∞(x,t). T∞ can be a Number or Function (x, t) -> ambient_temp.

source

Mechanics BCs

Macchiato.Displacement Type

Displacement{F<:Function} <: Dirichlet

Prescribed displacement BC for solid mechanics. The function returns a tuple of displacement components: f(x, t) -> (ux, uy) for 2D.

Constructors

Displacement((x, t) -> (0.0, 0.0))          # Function returning tuple
Displacement(0.0, 0.0)                       # Constant displacement
Displacement(ux::Function, uy::Function)     # Per-component functions

source

Macchiato.Traction Type

Traction{F<:Function} <: Neumann

Prescribed traction BC for solid mechanics. The function returns a tuple of traction components: f(x, t) -> (tx, ty) for 2D.

In terms of stress: t = σ·n where n is the outward normal.

Constructors

Traction((x, t) -> (0.0, -1000.0))          # Function returning tuple
Traction(tx::Number, ty::Number)             # Constant traction
Traction(tx::Function, ty::Function)         # Per-component functions

source

Macchiato.TractionFree Function

TractionFree()

Zero-traction (free surface) BC: σ·n = 0. Convenience for Traction(0.0, 0.0).

source

Fluid BCs

Macchiato.VelocityInlet Function

VelocityInlet(velocity)

Prescribed velocity at inlet. Value can be a Number or Function (x, t) -> velocity.

source

Macchiato.PressureOutlet Function

PressureOutlet(pressure)

Prescribed pressure at outlet. Value can be a Number or Function (x, t) -> pressure.

source

Macchiato.VelocityOutlet Function

VelocityOutlet()

Zero-gradient velocity outlet: ∂v/∂n = 0. Used for fully developed outflow.

source

Wall BC

Macchiato.Wall Function

Wall(velocity)
Wall()

No-slip wall or moving wall BC. Value can be a Number or Function (x, t) -> velocity. No arguments creates stationary wall (v=0).

source

Simulation

Macchiato.Simulation Type

Simulation{M, C, Mode}

High-level simulation container that manages solving and solution storage.

Fields

• domain::Domain{M, C}: The computational domain with models and boundary conditions

• mode::Mode: Simulation mode (Steady or Transient)

• u0: Initial condition vector (transient only)

• time: Current simulation time

• running: Whether simulation is currently running

• _solution: Solution vector

Constructors

Steady-state simulation

sim = Simulation(domain)
sim = Simulation(domain, Steady())

Transient simulation

sim = Simulation(domain, Transient(Δt=0.001, stop_time=1.0))

source

Macchiato.run! Function

run!(sim::Simulation; kwargs...)

Execute the simulation. Dispatches to transient or steady-state path based on mode. Extra kwargs are forwarded to OrdinaryDiffEq.solve (transient) or LinearSolve.LinearProblem (steady).

Returns the simulation object.

source

Macchiato.set! Function

set!(sim::Simulation; kwargs...)

Set initial conditions for simulation fields.

Arguments

• sim: Simulation to set initial conditions for

• kwargs: Field name/value pairs

Supported value types

• Number: Uniform value for entire field

• Function: Called as f(x) where x is coordinate vector [x, y] or [x, y, z]

• Vector: Direct assignment (must match field length)

Examples

set!(sim, T=300.0)                           # Uniform temperature
set!(sim, T=x -> 300 + 10*x[1])              # Temperature function of position
set!(sim, u=0.0, v=0.0, p=0.0)               # Multiple fields

source

Field Extraction

Macchiato.solution Function

solution(sim) -> Vector{Float64}

Return the raw solution vector for the simulation.

For built-in models, prefer the typed accessors temperature, velocity, pressure, and displacement. For custom PDEs this is the supported way to read results after run!. Falls back to the initial condition if run! has not been called; throws an ArgumentError if no solution or initial condition is available.

source

Macchiato.temperature Function

temperature(sim) -> Vector{Float64}

Extract temperature field from simulation.

source

Macchiato.velocity Function

velocity(sim) -> Tuple{Vector{Float64}, ...}

Extract velocity components from simulation. Returns (u, v) for 2D or (u, v, w) for 3D.

source

Macchiato.pressure Function

pressure(sim) -> Vector{Float64}

Extract pressure field from simulation.

source

Macchiato.displacement Function

displacement(sim) -> Tuple{Vector{Float64}, ...}

Extract displacement components from simulation. Returns (ux, uy) for 2D or (ux, uy, uz) for 3D.

source

Solvers

SciMLBase.LinearProblem Type

LinearSolve.LinearProblem(domain::Domain; scheme=nothing, verbose=false, kwargs...)

Construct a LinearProblem for steady-state simulation from a Domain.

Assembles the system matrix A and right-hand side b from the physics model, then applies boundary conditions by modifying the appropriate rows of A and b. The resulting system Ax = b is solved with LinearSolve.solve. Pass verbose=true to print assembly progress.

Steady-state solving currently supports a single physics model per domain.

source

I/O

Macchiato.exportvtk Function

exportvtk(filename, points, data, names)

Export point-based simulation results to a VTK file.

Arguments

• filename::String: Output file path (without .vtu extension)

• points::AbstractVector: Point cloud or coordinate vector

• data::AbstractVector{<:AbstractVector}: Field data arrays to export

• names::AbstractVector: Corresponding field names (e.g., ["T", "u"])

Example

exportvtk("results/temperature", points(cloud), [T_values], ["T"])

source

Macchiato.savevtk! Function

savevtk!(vtkfile)

Write the VTK file to disk.

source

Operators

Macchiato.upwind Function

upwind(data, eval_points, dim[, basis]; Δ=nothing, k=autoselect_k(data, basis))
upwind(data, dim[, basis]; Δ=nothing, k=autoselect_k(data, basis))

Build an upwind finite-difference-style operator using RBF interpolation.

Computes backward, forward, and centered partial derivatives with respect to dimension dim, then returns a function (ϕ, v, θ) that blends them based on flow direction v and upwind parameter θ ∈ [0, 1] (1 = full upwind, 0 = centered). The returned callable computes the advective term v .* ∂ϕ (the derivative scaled by the local velocity v), not the bare derivative.

The single-argument form upwind(data, dim) evaluates at the data points themselves.

Arguments

• data: Stencil points for RBF approximation

• eval_points: Points where the derivative is evaluated

• dim: Spatial dimension (1 = x, 2 = y, …)

• basis: Radial basis function (default: PHS(3; poly_deg=2))

• Δ: Virtual node offset distance (auto-detected if nothing)

• k: Number of nearest neighbors for stencil

source