diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml
index 048406a..94544a8 100644
--- a/.github/workflows/CI.yml
+++ b/.github/workflows/CI.yml
@@ -1,28 +1,48 @@
 name: CI
 on:
   push:
-    branches: [main]
+    branches:
+      - main
     tags: ['*']
   pull_request:
-permissions:
-  contents: read
+  workflow_dispatch:
 concurrency:
+  # Skip intermediate builds: always.
+  # Cancel intermediate builds: only if it is a pull request build.
   group: ${{ github.workflow }}-${{ github.ref }}
-  cancel-in-progress: true
+  cancel-in-progress: ${{ startsWith(github.ref, 'refs/pull/') }}
 jobs:
   test:
-    name: Julia ${{ matrix.version }} - ${{ matrix.os }}
+    name: Julia ${{ matrix.version }} - ${{ matrix.os }} - ${{ matrix.arch }} - ${{ github.event_name }}
     runs-on: ${{ matrix.os }}
+    timeout-minutes: 60
+    permissions: # needed to allow julia-actions/cache to proactively delete old caches that it has created
+      actions: write
+      contents: read
     strategy:
       fail-fast: false
       matrix:
-        version: ['1.10', '1.11']
-        os: [ubuntu-latest]
+        version:
+          - '1.12'
+          - '1.11'
+          - '1.10'
+        os:
+          - ubuntu-latest
+        arch:
+          - x64
+          # - x86
     steps:
-      - uses: actions/checkout@v4
-      - uses: julia-actions/setup-julia@v2
+      - uses: actions/checkout@v7
+      - uses: julia-actions/setup-julia@v3
         with:
           version: ${{ matrix.version }}
-      - uses: julia-actions/cache@v2
+          arch: ${{ matrix.arch }}
+      - uses: julia-actions/cache@v3
       - uses: julia-actions/julia-buildpkg@v1
       - uses: julia-actions/julia-runtest@v1
+      - uses: julia-actions/julia-processcoverage@v1
+      - uses: codecov/codecov-action@v7
+        with:
+          files: lcov.info
+          token: ${{ secrets.CODECOV_TOKEN }}
+          fail_ci_if_error: false
diff --git a/.github/workflows/CompatHelper.yml b/.github/workflows/CompatHelper.yml
new file mode 100644
index 0000000..43d2889
--- /dev/null
+++ b/.github/workflows/CompatHelper.yml
@@ -0,0 +1,18 @@
+name: CompatHelper
+on:
+  schedule:
+    - cron: 0 0 * * *
+  workflow_dispatch:
+permissions:
+  contents: read
+jobs:
+  CompatHelper:
+    runs-on: ubuntu-latest
+    steps:
+      - name: Pkg.add("CompatHelper")
+        run: julia -e 'using Pkg; Pkg.add("CompatHelper")'
+      - name: CompatHelper.main()
+        env:
+          GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }}
+          COMPATHELPER_PRIV: ${{ secrets.DOCUMENTER_KEY }}
+        run: julia -e 'using CompatHelper; CompatHelper.main()'
diff --git a/.github/workflows/Formatter.yml b/.github/workflows/Formatter.yml
new file mode 100644
index 0000000..950e235
--- /dev/null
+++ b/.github/workflows/Formatter.yml
@@ -0,0 +1,43 @@
+name: Formatter
+on:
+  pull_request:
+  push:
+    branches:
+      - main
+  workflow_dispatch:
+permissions:
+  contents: read
+jobs:
+  julia-format:
+    name: Formatter
+    runs-on: ubuntu-latest
+    steps:
+      - uses: julia-actions/setup-julia@latest
+      - uses: actions/checkout@v6
+        with:
+          token: ${{ secrets.GITHUB_TOKEN }}
+      - name: Install JuliaFormatter and format
+        run: |
+          julia --color=yes -e 'using Pkg; Pkg.add(PackageSpec(name="JuliaFormatter"))'
+          julia --color=yes -e 'using JuliaFormatter; format(".", verbose=true)'
+      - name: Commit formatted code
+        if: github.event_name == 'workflow_dispatch'
+        run: |
+          git config user.name "github-actions[bot]"
+          git config user.email "github-actions[bot]@users.noreply.github.com"
+          git diff --quiet && exit 0
+          git add -A
+          git commit -m "chore: autoformat"
+          git push
+      - name: Format check
+        if: github.event_name != 'workflow_dispatch'
+        run: |
+          julia --color=yes -e '
+          out = Cmd(`git diff --name-only`) |> read |> String
+          if out == ""
+            exit(0)
+          else
+            @error "Some files have not been formatted !!!"
+            write(stdout, out)
+            exit(1)
+          end'
diff --git a/.github/workflows/docs.yml b/.github/workflows/docs.yml
new file mode 100644
index 0000000..b7e4800
--- /dev/null
+++ b/.github/workflows/docs.yml
@@ -0,0 +1,75 @@
+name: Documentation
+on:
+  pull_request:
+  push:
+    branches:
+      - main
+    tags: ["*"]
+
+  repository_dispatch:
+    types: [tagbot-release-created]
+
+concurrency:
+  # Skip intermediate builds: always.
+  # Cancel intermediate builds: only if it is a pull request build.
+  group: ${{ github.workflow }}-${{ github.ref }}
+  cancel-in-progress: ${{ startsWith(github.ref, 'refs/pull/') }}
+jobs:
+  docs:
+    name: Documentation
+    runs-on: ubuntu-latest
+    timeout-minutes: 30
+    permissions:
+      contents: write
+      statuses: write
+    env:
+      DOC_TEMPLATE_VERSION: "v0.7.0" # Change this to the specific tag version you want
+    steps:
+      - uses: actions/checkout@v7
+      - uses: julia-actions/setup-julia@v3
+      - uses: julia-actions/cache@v3
+      - name: Use Documentation Template
+        run: |
+          ./docs/get_docs_utils.sh ${{ env.DOC_TEMPLATE_VERSION }}
+      - uses: julia-actions/julia-buildpkg@v1
+      - uses: julia-actions/julia-docdeploy@v1
+        env:
+          GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }}
+      - name: Upload documentation artifacts
+        uses: actions/upload-artifact@v6
+        with:
+          name: documentation-build
+          path: docs/build/
+          retention-days: 1
+      - run: |
+          julia --project=docs -e '
+            using Documenter: DocMeta, doctest
+            using Intonato
+            DocMeta.setdocmeta!(Intonato, :DocTestSetup, :(using Intonato); recursive=true)
+            doctest(Intonato)'
+
+  tagbot-dispatch:
+    name: Dispatch on TagBot Release
+    runs-on: ubuntu-latest
+    needs: docs
+    if: github.event_name == 'repository_dispatch'
+    steps:
+      - name: Dispatch PiccoloMultiDocs workflow
+        uses: peter-evans/repository-dispatch@v4
+        with:
+          token: ${{ secrets.REPO_ACCESS_TOKEN }}
+          repository: harmoniqs/PiccoloMultiDocs.jl
+          event-type: rebuild-docs
+
+  tag-push-dispatch:
+    name: Dispatch on Tag Push
+    runs-on: ubuntu-latest
+    needs: docs
+    if: github.ref_type == 'tag'
+    steps:
+      - name: Dispatch PiccoloMultiDocs workflow
+        uses: peter-evans/repository-dispatch@v4
+        with:
+          token: ${{ secrets.REPO_ACCESS_TOKEN }}
+          repository: harmoniqs/PiccoloMultiDocs.jl
+          event-type: rebuild-docs
diff --git a/.github/workflows/nightly.yml b/.github/workflows/nightly.yml
new file mode 100644
index 0000000..cc4da68
--- /dev/null
+++ b/.github/workflows/nightly.yml
@@ -0,0 +1,34 @@
+name: Nightly
+on:
+  schedule:
+    - cron: '00 00 * * *'
+  workflow_dispatch:
+concurrency:
+  group: ${{ github.workflow }}-${{ github.ref }}
+  cancel-in-progress: ${{ startsWith(github.ref, 'refs/pull/') }}
+jobs:
+  test:
+    name: Julia nightly - ${{ matrix.os }} - ${{ matrix.arch }} - ${{ github.event_name }}
+    runs-on: ${{ matrix.os }}
+    timeout-minutes: 60
+    permissions:
+      actions: write
+      contents: read
+    strategy:
+      fail-fast: false
+      matrix:
+        version:
+          - 'pre'
+        os:
+          - ubuntu-latest
+        arch:
+          - x64
+    steps:
+      - uses: actions/checkout@v7
+      - uses: julia-actions/setup-julia@v3
+        with:
+          version: ${{ matrix.version }}
+          arch: ${{ matrix.arch }}
+      - uses: julia-actions/cache@v3
+      - uses: julia-actions/julia-buildpkg@v1
+      - uses: julia-actions/julia-runtest@v1
diff --git a/.gitignore b/.gitignore
index 5bddb7c..9d92b7a 100644
--- a/.gitignore
+++ b/.gitignore
@@ -1,4 +1,59 @@
+# Files generated by invoking Julia with --code-coverage
+*.jl.cov
+*.jl.*.cov
+
+# Files generated by invoking Julia with --track-allocation
+*.jl.mem
+
+# System-specific files and directories generated by the BinaryProvider and BinDeps packages
+# They contain absolute paths specific to the host computer, and so should not be committed
+deps/deps.jl
+deps/build.log
+deps/downloads/
+deps/usr/
+deps/src/
+
+# Build artifacts for creating documentation generated by the Documenter package
+docs/build/
+docs/site/
+docs/clones/
+docs/out/
+docs/src/generated/
+
+docs/src/assets/
+
+# Build artifact for using README as index (probably shouldn't have this checked into version control)
+docs/src/index.md
+
+# File generated by Pkg, the package manager, based on a corresponding Project.toml
+# It records a fixed state of all packages used by the project. As such, it should not be
+# committed for packages, but should be committed for applications that require a static
+# environment.
 Manifest.toml
+docs/Manifest.toml
+
+# Project specific ignores below
+# generated example artifacts
+/examples/**/plots/
+/examples/**/trajectories/
+
+# external pkgs and configs
+pardiso.lic
+/.CondaPkg/
+*.code-workspace
+
+# generated build folder
+build/
+
+# VS code
+.vscode/settings.json
+
+# doc_template stuff
+# Temporary directory for doc_template cloning
+doc_template_temp/
+
+# This file is updated via script
+docs/utils.jl
 
 # Agent working memory (not package content)
 .claude/
diff --git a/docs/Project.toml b/docs/Project.toml
new file mode 100644
index 0000000..abe1536
--- /dev/null
+++ b/docs/Project.toml
@@ -0,0 +1,7 @@
+[deps]
+Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
+Intonato = "4c8581c7-0eaf-45eb-b8fa-a3474c50779c"
+Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306"
+LiveServer = "16fef848-5104-11e9-1b77-fb7a48bbb589"
+PiccoloDocsTemplate = "a90a139f-c522-4b23-980b-4210ddb8d065"
+Revise = "295af30f-e4ad-537b-8983-00126c2a3abe"
diff --git a/docs/get_docs_utils.sh b/docs/get_docs_utils.sh
new file mode 100755
index 0000000..9723a50
--- /dev/null
+++ b/docs/get_docs_utils.sh
@@ -0,0 +1,38 @@
+#!/bin/bash
+
+set -euo pipefail
+
+# Get the directory where this script is located
+SCRIPT_DIR="$(cd "$(dirname "${BASH_SOURCE[0]}")" && pwd)"
+PROJECT_ROOT="$SCRIPT_DIR/.."
+
+# if argument is provided, use it as the DOC_TEMPLATE_VERSION
+if [[ $# -gt 0 ]]; then
+    DOC_TEMPLATE_VERSION="$1"
+else
+    WORKFLOW_FILE="$PROJECT_ROOT/.github/workflows/docs.yml"
+
+    # Check if workflow file exists
+    if [[ ! -f "$WORKFLOW_FILE" ]]; then
+        echo "GitHub workflow file not found at: $WORKFLOW_FILE"
+        exit 1
+    fi
+
+    DOC_TEMPLATE_VERSION=$(grep -E '^\s*DOC_TEMPLATE_VERSION:' "$WORKFLOW_FILE" | sed -E 's/.*DOC_TEMPLATE_VERSION:\s*"([^"]+)".*/\1/')
+fi
+
+if [[ -z "$DOC_TEMPLATE_VERSION" ]]; then
+    echo "DOC_TEMPLATE_VERSION is not set"
+    echo "Please provide a version tag as an arg or ensure it is set in $WORKFLOW_FILE"
+    echo "Could not extract DOC_TEMPLATE_VERSION from $WORKFLOW_FILE"
+    echo "Expected format: DOC_TEMPLATE_VERSION: \"<version tag here>\""
+    exit 1
+fi
+
+# Clone the repository
+echo "Grabbing PiccoloDocsTemplate at version $DOC_TEMPLATE_VERSION"
+julia --project="$PROJECT_ROOT/docs" -e "
+using Pkg; Pkg.add(url=\"https://github.com/harmoniqs/PiccoloDocsTemplate.jl\", rev=\"$DOC_TEMPLATE_VERSION\")
+"
+
+echo "Successfully updated PiccoloDocsTemplate with version $DOC_TEMPLATE_VERSION"
\ No newline at end of file
diff --git a/docs/make.jl b/docs/make.jl
new file mode 100644
index 0000000..cc74057
--- /dev/null
+++ b/docs/make.jl
@@ -0,0 +1,15 @@
+using Intonato
+using PiccoloDocsTemplate
+
+pages = ["Home" => "index.md", "Library" => "lib.md"]
+
+generate_docs(
+    @__DIR__,
+    "Intonato",
+    [Intonato],
+    pages;
+    make_literate = false,
+    make_assets = false,
+    format_kwargs = (canonical = "https://docs.harmoniqs.co/Intonato.jl",),
+    versions = ["dev" => "dev", "stable" => "v^", "v#.#"],
+)
diff --git a/docs/src/lib.md b/docs/src/lib.md
new file mode 100644
index 0000000..b705e42
--- /dev/null
+++ b/docs/src/lib.md
@@ -0,0 +1,22 @@
+# API
+
+```@meta
+CollapsedDocStrings = true
+```
+
+The central entry point. `solve!(::PulseTuningProblem)` is defined as a
+`Piccolo.solve!` method (the function-name binding is owned by the `Piccolo` /
+`DirectTrajOpt` module, in scope here via Intonato's `@reexport using Piccolo`).
+It is documented explicitly so it is guaranteed to appear even if a future
+refactor changes which module owns the binding; the `@autodocs` block below is
+filtered to skip this method so it is not documented twice.
+
+```@docs
+solve!(::Intonato.PulseTuningProblem)
+```
+
+```@autodocs
+Modules = [Intonato]
+Order = [:type, :function]
+Filter = b -> b !== solve!
+```
diff --git a/src/Intonato.jl b/src/Intonato.jl
index bbfe017..b5d57f9 100644
--- a/src/Intonato.jl
+++ b/src/Intonato.jl
@@ -49,7 +49,8 @@ include("problems/test.jl")
 # Types
 export Measurement, MeasurementModel
 export AbstractExperiment, SimulatedExperiment, HardwareExperiment
-export AbstractMeasurement, DeterministicMeasurement, ShotNoiseMeasurement, KnownCovarianceMeasurement
+export AbstractMeasurement,
+    DeterministicMeasurement, ShotNoiseMeasurement, KnownCovarianceMeasurement
 export pauli_covariance, population_covariance, wigner_covariance
 export pauli, pop
 export AbstractHardwareBackend
@@ -81,7 +82,8 @@ export AbstractDeviceModel, NominalModel, predict, adapt!
 # IdentityStrategy stand-in. Concrete tuning strategies plug in via this
 # interface.
 export AbstractPulseTuningProblem, AbstractTuningStrategy, IdentityStrategy, step
-export prepare_strategy, tuning_goal, candidate_trajectory, last_timings, accepts_global_data
+export prepare_strategy,
+    tuning_goal, candidate_trajectory, last_timings, accepts_global_data
 
 # Closed-loop tuning chassis. The chassis is strategy-generic; its
 # result/record types are public.
diff --git a/src/device_models/abstract.jl b/src/device_models/abstract.jl
index 0e3bfed..08f288f 100644
--- a/src/device_models/abstract.jl
+++ b/src/device_models/abstract.jl
@@ -26,7 +26,7 @@ dispatched away at compile time — `NominalModel` carries no learned parameters
 - `ψ_init::K`: initial state
 - `ψ_goal::G`: goal state
 """
-struct NominalModel{S, K, G} <: AbstractDeviceModel
+struct NominalModel{S,K,G} <: AbstractDeviceModel
     system::S            # QuantumSystem / OpenQuantumSystem
     ψ_init::K
     ψ_goal::G
@@ -71,12 +71,13 @@ adapt!(m::NominalModel, ::AbstractPulse, ::Any) = m   # no-op via dispatch
     σz = ComplexF64[1 0; 0 -1]
     sys = QuantumSystem(1.0 * σz, [σx], [1.0])
 
-    N = 11; T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    N = 11;
+    T = 5.0
+    times = range(0.0, T, length = N) |> collect
     ψ0 = ComplexF64[1.0, 0.0]
     ψg = ComplexF64[0.0, 1.0]
     pulse = LinearSplinePulse(0.1 * ones(1, N), times)
-    model = MeasurementModel(:ψ̃, [full_state for _ in 1:N], collect(1:N))
+    model = MeasurementModel(:ψ̃, [full_state for _ = 1:N], collect(1:N))
 
     # A NominalModel wraps a QuantumSystem; predict rolls the pulse through it
     # and evaluates the measurement model. adapt! is a no-op.
diff --git a/src/measurement_functions/displaced_parity.jl b/src/measurement_functions/displaced_parity.jl
index 9a8272f..5eb47be 100644
--- a/src/measurement_functions/displaced_parity.jl
+++ b/src/measurement_functions/displaced_parity.jl
@@ -32,7 +32,7 @@ end
 Parity operator Π = (-1)^(a†a) = diag(1, -1, 1, -1, …).
 """
 function _parity_operator(n_max::Int)
-    return Diagonal([iseven(n) ? 1.0 : -1.0 for n in 0:n_max-1])
+    return Diagonal([iseven(n) ? 1.0 : -1.0 for n = 0:(n_max-1)])
 end
 
 """
@@ -42,7 +42,7 @@ Bosonic annihilation operator a in the Fock basis.
 """
 function _annihilation_operator(n_max::Int)
     a = zeros(ComplexF64, n_max, n_max)
-    for n in 1:n_max-1
+    for n = 1:(n_max-1)
         a[n, n+1] = sqrt(n)
     end
     return a
@@ -54,6 +54,6 @@ end
 Create a closure g(x_iso) → [displaced_parity(x_iso, α)] for MeasurementModel.
 """
 displaced_parity_at(α::Complex; n_max::Int) =
-    x_iso -> [displaced_parity(x_iso, α; n_max=n_max)]
+    x_iso -> [displaced_parity(x_iso, α; n_max = n_max)]
 
 # ──── Tests ──────────────────────────────────────────────────────────────────
diff --git a/src/measurement_functions/partial_trace.jl b/src/measurement_functions/partial_trace.jl
index acfefef..d8d8730 100644
--- a/src/measurement_functions/partial_trace.jl
+++ b/src/measurement_functions/partial_trace.jl
@@ -16,7 +16,7 @@ function partial_trace_B(ρ_iso::AbstractVector, dims::Tuple{Int,Int})
 
     # Trace over B: ρ_A[i_A, j_A] = Σ_b ρ_tensor[b, i_A, b, j_A]
     ρ_A = zeros(eltype(ρ), d_A, d_A)
-    for b in 1:d_B
+    for b = 1:d_B
         ρ_A .+= ρ_tensor[b, :, b, :]
     end
 
diff --git a/src/measurement_functions/state_measurements.jl b/src/measurement_functions/state_measurements.jl
index 555ff0f..8f11a4f 100644
--- a/src/measurement_functions/state_measurements.jl
+++ b/src/measurement_functions/state_measurements.jl
@@ -8,7 +8,7 @@ Compute level populations |ψ_j|² from iso-vec ket (Re(ψ), Im(ψ)).
 function populations(x::AbstractVector)
     n = length(x) ÷ 2
     x_re = @view x[1:n]
-    x_im = @view x[n+1:2n]
+    x_im = @view x[(n+1):2n]
     return x_re .^ 2 .+ x_im .^ 2
 end
 
@@ -69,10 +69,7 @@ end
 
 Compute expectations for multiple observables.
 """
-function observable_expectations(
-    x_iso::AbstractVector,
-    Os_iso::Vector{<:AbstractMatrix},
-)
+function observable_expectations(x_iso::AbstractVector, Os_iso::Vector{<:AbstractMatrix})
     return [dot(x_iso, O * x_iso) for O in Os_iso]
 end
 
diff --git a/src/measurement_functions/test.jl b/src/measurement_functions/test.jl
index 619bb48..c8cb0fc 100644
--- a/src/measurement_functions/test.jl
+++ b/src/measurement_functions/test.jl
@@ -13,7 +13,7 @@
     x0 = ket_to_iso(ψ0)
 
     # Parity at origin for vacuum: ⟨0|Π|0⟩ = 1
-    P0 = displaced_parity(x0, 0.0 + 0.0im; n_max=n_max)
+    P0 = displaced_parity(x0, 0.0 + 0.0im; n_max = n_max)
     @test P0 ≈ 1.0 atol = 1e-10
 
     # Displaced parity is related to Wigner: P(α) = (π/2) W(α)
@@ -29,7 +29,7 @@ end
     x = ket_to_iso(ψ)
 
     α = 0.3 + 0.2im
-    g = x -> displaced_parity(x, α; n_max=n_max)
+    g = x -> displaced_parity(x, α; n_max = n_max)
     grad = ForwardDiff.gradient(g, x)
     @test length(grad) == length(x)
     @test all(isfinite, grad)
@@ -100,9 +100,11 @@ end
     # Finite-difference check
     ε = 1e-7
     J_fd = zeros(2, 4)
-    for i in 1:4
-        x_p = copy(x); x_p[i] += ε
-        x_m = copy(x); x_m[i] -= ε
+    for i = 1:4
+        x_p = copy(x);
+        x_p[i] += ε
+        x_m = copy(x);
+        x_m[i] -= ε
         J_fd[:, i] = (populations(x_p) - populations(x_m)) / (2ε)
     end
     @test J ≈ J_fd atol = 1e-5
@@ -154,9 +156,11 @@ end
     # Finite-difference check
     ε = 1e-7
     J_fd = zeros(2, 4)
-    for i in 1:4
-        x_p = copy(ρ̃); x_p[i] += ε
-        x_m = copy(ρ̃); x_m[i] -= ε
+    for i = 1:4
+        x_p = copy(ρ̃);
+        x_p[i] += ε
+        x_m = copy(ρ̃);
+        x_m[i] -= ε
         J_fd[:, i] = (populations_density(x_p) - populations_density(x_m)) / (2ε)
     end
     @test J ≈ J_fd atol = 1e-5
@@ -176,12 +180,12 @@ end
     ρ_iso = density_to_iso_vec(ρ)
 
     # Wigner at origin should be 2/π for vacuum
-    W0 = wigner(ρ_iso, 0.0 + 0.0im; n_max=n_max)
+    W0 = wigner(ρ_iso, 0.0 + 0.0im; n_max = n_max)
     @test W0 ≈ 2 / π atol = 1e-10
 
     # Wigner should be positive everywhere for vacuum (Gaussian)
-    for r in 0.0:0.5:2.0
-        W = wigner(ρ_iso, r + 0.0im; n_max=n_max)
+    for r = 0.0:0.5:2.0
+        W = wigner(ρ_iso, r + 0.0im; n_max = n_max)
         @test W > -1e-10  # non-negative (up to numerics)
     end
 end
@@ -196,7 +200,7 @@ end
     ρ_iso = density_to_iso_vec(ρ)
 
     α = 0.5 + 0.3im
-    g = ρ_iso -> wigner(ρ_iso, α; n_max=n_max)
+    g = ρ_iso -> wigner(ρ_iso, α; n_max = n_max)
     grad = ForwardDiff.gradient(g, ρ_iso)
     @test length(grad) == length(ρ_iso)
     @test all(isfinite, grad)
diff --git a/src/measurement_functions/wigner.jl b/src/measurement_functions/wigner.jl
index b8551f7..3ff0f63 100644
--- a/src/measurement_functions/wigner.jl
+++ b/src/measurement_functions/wigner.jl
@@ -24,8 +24,8 @@ function _wigner_from_density(ρ::AbstractMatrix, α::Complex, n_max::Int)
     x = 4 * abs2(α)
     W = 0.0
 
-    for m in 0:n_max-1
-        for n in 0:n_max-1
+    for m = 0:(n_max-1)
+        for n = 0:(n_max-1)
             if m >= n
                 # L_n^{m-n}(x) via recurrence
                 k = m - n  # order of associated Laguerre
@@ -66,7 +66,7 @@ function _laguerre(n::Int, k::Int, x)
     L_prev2 = one(x)
     L_prev1 = 1 + k - x
 
-    for i in 2:n
+    for i = 2:n
         L_curr = ((2i - 1 + k - x) * L_prev1 - (i - 1 + k) * L_prev2) / i
         L_prev2 = L_prev1
         L_prev1 = L_curr
@@ -82,6 +82,6 @@ using SpecialFunctions: loggamma
 
 Create a closure g(ρ_iso) → [W(α)] for use in MeasurementModel.
 """
-wigner_at(α::Complex; n_max::Int) = ρ_iso -> [wigner(ρ_iso, α; n_max=n_max)]
+wigner_at(α::Complex; n_max::Int) = ρ_iso -> [wigner(ρ_iso, α; n_max = n_max)]
 
 # ──── Tests ──────────────────────────────────────────────────────────────────
diff --git a/src/optimizers/test.jl b/src/optimizers/test.jl
index 80bef04..5d5e1c6 100644
--- a/src/optimizers/test.jl
+++ b/src/optimizers/test.jl
@@ -41,7 +41,7 @@ end
 
     N = 11
     T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    times = range(0.0, T, length = N) |> collect
 
     ψ_init = ComplexF64[1.0, 0.0]
     ψ_goal = ComplexF64[0.0, 1.0]
@@ -82,7 +82,7 @@ end
 
     N = 11
     T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    times = range(0.0, T, length = N) |> collect
 
     ψ_init = ComplexF64[1.0, 0.0]
     ψ_goal = ComplexF64[0.0, 1.0]
diff --git a/src/problems/pulse_tuning_problem.jl b/src/problems/pulse_tuning_problem.jl
index 602294d..5a72121 100644
--- a/src/problems/pulse_tuning_problem.jl
+++ b/src/problems/pulse_tuning_problem.jl
@@ -23,10 +23,7 @@ struct IterationRecord
     # Per-phase wall-clock timings (seconds) for this outer iteration.
     # `nlp` is the inner subproblem solve; `total` excludes post-iter
     # rollout/record overhead. Phases that didn't run for a given iter stay 0.0.
-    timing::NamedTuple{
-        (:experiment, :sysid, :nlp, :armijo, :total),
-        NTuple{5, Float64},
-    }
+    timing::NamedTuple{(:experiment, :sysid, :nlp, :armijo, :total),NTuple{5,Float64}}
 end
 
 """
@@ -63,8 +60,12 @@ poorly-converged warm-start is wasteful. Silently passes if `fidelity` isn't
 defined for the QCP type. Pass `threshold = 0` to skip. `path_label` is used
 only in log messages.
 """
-function _check_nominal_fidelity(qcp, threshold::Real;
-                                 verbose::Bool, path_label::AbstractString)
+function _check_nominal_fidelity(
+    qcp,
+    threshold::Real;
+    verbose::Bool,
+    path_label::AbstractString,
+)
     threshold > 0 || return nothing
     F_nom = try
         Piccolo.fidelity(qcp)
@@ -78,7 +79,7 @@ function _check_nominal_fidelity(qcp, threshold::Real;
             "min_nominal_fidelity ($threshold). The QCP likely did not " *
             "converge — running $path_label on this pulse is wasteful. " *
             "Either solve the QCP more carefully (more iterations, better " *
-            "init, or curriculum) or pass min_nominal_fidelity=0 to skip."
+            "init, or curriculum) or pass min_nominal_fidelity=0 to skip.",
         )
     end
     verbose && @info "$path_label: nominal fidelity check passed" F_nom threshold
@@ -115,16 +116,17 @@ tuned result.
   (a no-op for `NominalModel`).
 - `result::Union{Nothing, TuningResult}`: populated by `solve!`
 """
-mutable struct PulseTuningProblem{S<:AbstractTuningStrategy, M<:AbstractDeviceModel} <: AbstractPulseTuningProblem
+mutable struct PulseTuningProblem{S<:AbstractTuningStrategy,M<:AbstractDeviceModel} <:
+               AbstractPulseTuningProblem
     qcp::QuantumControlProblem
     experiment::AbstractExperiment
     measurement_model::MeasurementModel
     R_tr::NamedTuple
-    Q_meas::Union{Float64, Vector{Float64}}
+    Q_meas::Union{Float64,Vector{Float64}}
     # Chassis/strategy split: the inner step and the predictive device model.
     strategy::S
     device_model::M
-    result::Union{Nothing, TuningResult}
+    result::Union{Nothing,TuningResult}
 end
 
 """
@@ -147,9 +149,9 @@ function PulseTuningProblem(
     experiment::AbstractExperiment,
     measurement_model::MeasurementModel;
     R_tr::NamedTuple = (;),
-    Q_meas::Union{Float64, Vector{Float64}} = 1.0,
-    strategy::Union{Nothing, AbstractTuningStrategy} = nothing,
-    device_model::Union{Nothing, AbstractDeviceModel} = nothing,
+    Q_meas::Union{Float64,Vector{Float64}} = 1.0,
+    strategy::Union{Nothing,AbstractTuningStrategy} = nothing,
+    device_model::Union{Nothing,AbstractDeviceModel} = nothing,
 )
     # Default strategy: a lightweight no-op placeholder. A tuning strategy is
     # provided by passing `strategy=`; the strategy carries its own config and
@@ -164,7 +166,14 @@ function PulseTuningProblem(
     end
 
     return PulseTuningProblem(
-        qcp, experiment, measurement_model, R_tr, Q_meas, strat, devmodel, nothing
+        qcp,
+        experiment,
+        measurement_model,
+        R_tr,
+        Q_meas,
+        strat,
+        devmodel,
+        nothing,
     )
 end
 
@@ -205,13 +214,12 @@ function Piccolo.solve!(
     line_search::Bool = true,
     verbose::Bool = true,
     γ::Float64 = 0.8,
-    max_rejections::Union{Nothing, Int} = nothing,
+    max_rejections::Union{Nothing,Int} = nothing,
     min_nominal_fidelity::Float64 = 0.8,
     polyak_avg::Int = 0,
 )
     # Sanity check: refuse to run on a poorly-converged nominal QCP.
-    _check_nominal_fidelity(ptp.qcp, min_nominal_fidelity;
-                            verbose, path_label="QILC")
+    _check_nominal_fidelity(ptp.qcp, min_nominal_fidelity; verbose, path_label = "QILC")
 
     qcp = ptp.qcp
     z_ref = qcp.prob.trajectory
@@ -238,15 +246,18 @@ function Piccolo.solve!(
     # at the end of solve!. Reduces variance of the final iterate vs picking
     # the literal final or "best-J" iterate, which can chase shot noise on
     # hardware. polyak_avg=0 disables (current behavior preserved).
-    polyak_data_acc   = polyak_avg > 0 ? zero(z_ref.data) : nothing
-    polyak_global_acc = polyak_avg > 0 && z_ref.global_dim > 0 ?
-                        zero(z_ref.global_data) : nothing
+    polyak_data_acc = polyak_avg > 0 ? zero(z_ref.data) : nothing
+    polyak_global_acc =
+        polyak_avg > 0 && z_ref.global_dim > 0 ? zero(z_ref.global_data) : nothing
     polyak_count = 0
 
-    for i in 1:max_iter
+    for i = 1:max_iter
         # Per-iter timing accumulator. Phases that don't run for an iter
         # stay at 0.0.
-        t_experiment = 0.0; t_sysid = 0.0; t_nlp = 0.0; t_armijo = 0.0
+        t_experiment = 0.0;
+        t_sysid = 0.0;
+        t_nlp = 0.0;
+        t_armijo = 0.0
         t_iter_start = time()
 
         # 1. Extract pulse from current trajectory and run experiment
@@ -262,9 +273,23 @@ function Piccolo.solve!(
 
         if J_exp ≤ tol
             t_total = time() - t_iter_start
-            push!(history, IterationRecord(y_exp, J_exp, 1.0, tr_scale, pulse,
-                (experiment = t_experiment, sysid = 0.0, nlp = 0.0,
-                 armijo = 0.0, total = t_total)))
+            push!(
+                history,
+                IterationRecord(
+                    y_exp,
+                    J_exp,
+                    1.0,
+                    tr_scale,
+                    pulse,
+                    (
+                        experiment = t_experiment,
+                        sysid = 0.0,
+                        nlp = 0.0,
+                        armijo = 0.0,
+                        total = t_total,
+                    ),
+                ),
+            )
             converged = true
             break
         end
@@ -274,8 +299,16 @@ function Piccolo.solve!(
         # candidate trajectory populated for acceptance; it stashes the
         # per-phase sysid/NLP timings for this iteration's IterationRecord.
         ctx = (;
-            pulse, y_exp, J_exp, z_ref, iter = i, tr_scale,
-            ipopt_options, verbose, qcp, y_goal,
+            pulse,
+            y_exp,
+            J_exp,
+            z_ref,
+            iter = i,
+            tr_scale,
+            ipopt_options,
+            verbose,
+            qcp,
+            y_goal,
             device_model = ptp.device_model,
         )
         pulse_cand = step(strategy, ctx)
@@ -296,9 +329,8 @@ function Piccolo.solve!(
         if line_search
             local α_value
             t_armijo = @elapsed begin
-                α_value, ls_evals = armijo_line_search(
-                    ptp.experiment, pulse, pulse_cand, y_goal, J_exp
-                )
+                α_value, ls_evals =
+                    armijo_line_search(ptp.experiment, pulse, pulse_cand, y_goal, J_exp)
             end
             α = α_value
             n_experiments += ls_evals
@@ -318,9 +350,8 @@ function Piccolo.solve!(
             # don't overwrite the owned value; a strategy that co-optimizes
             # globals with the controls reports `true`.
             if z_ref.global_dim > 0 && accepts_global_data(strategy)
-                z_ref.global_data .= (
-                    (1 - α) .* z_ref.global_data .+ α .* cand_traj.global_data
-                )
+                z_ref.global_data .=
+                    ((1 - α) .* z_ref.global_data .+ α .* cand_traj.global_data)
             end
             consecutive_rejections = 0
         elseif α > 0.0
@@ -331,10 +362,25 @@ function Piccolo.solve!(
             consecutive_rejections += 1
             if !isnothing(max_rejections) && consecutive_rejections > max_rejections
                 t_total = time() - t_iter_start
-                push!(history, IterationRecord(y_exp, J_exp, α, tr_scale, pulse,
-                    (experiment = t_experiment, sysid = t_sysid, nlp = t_nlp,
-                     armijo = t_armijo, total = t_total)))
-                verbose && @info "QILC: stopping after $consecutive_rejections consecutive rejections"
+                push!(
+                    history,
+                    IterationRecord(
+                        y_exp,
+                        J_exp,
+                        α,
+                        tr_scale,
+                        pulse,
+                        (
+                            experiment = t_experiment,
+                            sysid = t_sysid,
+                            nlp = t_nlp,
+                            armijo = t_armijo,
+                            total = t_total,
+                        ),
+                    ),
+                )
+                verbose &&
+                    @info "QILC: stopping after $consecutive_rejections consecutive rejections"
                 break
             end
         end
@@ -352,11 +398,26 @@ function Piccolo.solve!(
         end
 
         t_total = time() - t_iter_start
-        push!(history, IterationRecord(y_exp, J_exp, α, tr_scale, pulse,
-            (experiment = t_experiment, sysid = t_sysid, nlp = t_nlp,
-             armijo = t_armijo, total = t_total)))
+        push!(
+            history,
+            IterationRecord(
+                y_exp,
+                J_exp,
+                α,
+                tr_scale,
+                pulse,
+                (
+                    experiment = t_experiment,
+                    sysid = t_sysid,
+                    nlp = t_nlp,
+                    armijo = t_armijo,
+                    total = t_total,
+                ),
+            ),
+        )
 
-        verbose && @info "QILC iter $i timing (s)" experiment=t_experiment sysid=t_sysid nlp=t_nlp armijo=t_armijo total=t_total
+        verbose &&
+            @info "QILC iter $i timing (s)" experiment=t_experiment sysid=t_sysid nlp=t_nlp armijo=t_armijo total=t_total
 
         # Polyak-Ruppert averaging: accumulate the last polyak_avg iters'
         # trajectory data into a running mean. Skips iters that line search
diff --git a/src/problems/test.jl b/src/problems/test.jl
index 591a8bb..c9d49b6 100644
--- a/src/problems/test.jl
+++ b/src/problems/test.jl
@@ -15,24 +15,24 @@
 
     N = 11
     T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    times = range(0.0, T, length = N) |> collect
     pulse = LinearSplinePulse(0.1 * randn(1, N), times)
 
     ψ_init = ComplexF64[1.0, 0.0]
     ψ_goal = ComplexF64[0.0, 1.0]
 
     qtraj_nom = KetTrajectory(sys_nom, pulse, ψ_init, ψ_goal)
-    qcp = SplinePulseProblem(qtraj_nom, N; Q=100.0, R=1e-2)
+    qcp = SplinePulseProblem(qtraj_nom, N; Q = 100.0, R = 1e-2)
 
     model = MeasurementModel(:ψ̃, [populations], [N])
     qtraj_true = KetTrajectory(sys_true, pulse, ψ_init, ψ_goal)
     experiment = SimulatedExperiment(qtraj_true, model)
 
-    ptp = PulseTuningProblem(qcp, experiment, model; R_tr=(u=1e-2,))
+    ptp = PulseTuningProblem(qcp, experiment, model; R_tr = (u = 1e-2,))
 
     @test ptp.qcp === qcp
     @test ptp.experiment === experiment
-    @test ptp.R_tr == (u=1e-2,)
+    @test ptp.R_tr == (u = 1e-2,)
     @test ptp.Q_meas == 1.0
     @test isnothing(ptp.result)
     # Default strategy is the no-op IdentityStrategy until a concrete strategy
diff --git a/src/pulse_ops/interpolation.jl b/src/pulse_ops/interpolation.jl
index b60b06d..5e903b3 100644
--- a/src/pulse_ops/interpolation.jl
+++ b/src/pulse_ops/interpolation.jl
@@ -8,18 +8,18 @@ function interpolate_pulse end
 
 function interpolate_pulse(p1::ZeroOrderPulse, p2::ZeroOrderPulse, α::Real)
     u = (1 - α) * p1.controls.u + α * p2.controls.u
-    return ZeroOrderPulse(u, get_knot_times(p1); drive_name=p1.drive_name)
+    return ZeroOrderPulse(u, get_knot_times(p1); drive_name = p1.drive_name)
 end
 
 function interpolate_pulse(p1::LinearSplinePulse, p2::LinearSplinePulse, α::Real)
     u = (1 - α) * get_knot_values(p1) + α * get_knot_values(p2)
-    return LinearSplinePulse(u, get_knot_times(p1); drive_name=p1.drive_name)
+    return LinearSplinePulse(u, get_knot_times(p1); drive_name = p1.drive_name)
 end
 
 function interpolate_pulse(p1::CubicSplinePulse, p2::CubicSplinePulse, α::Real)
-    u  = (1 - α) * get_knot_values(p1)      + α * get_knot_values(p2)
-    du = (1 - α) * get_knot_derivatives(p1)  + α * get_knot_derivatives(p2)
-    return CubicSplinePulse(u, du, get_knot_times(p1); drive_name=p1.drive_name)
+    u = (1 - α) * get_knot_values(p1) + α * get_knot_values(p2)
+    du = (1 - α) * get_knot_derivatives(p1) + α * get_knot_derivatives(p2)
+    return CubicSplinePulse(u, du, get_knot_times(p1); drive_name = p1.drive_name)
 end
 
 # ============================================================================ #
diff --git a/src/pulse_ops/test.jl b/src/pulse_ops/test.jl
index a1baace..9d0cdb2 100644
--- a/src/pulse_ops/test.jl
+++ b/src/pulse_ops/test.jl
@@ -7,12 +7,13 @@
 @testitem "interpolate_pulse LinearSplinePulse" begin
     using Intonato
 
-    N = 11; T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    N = 11;
+    T = 5.0
+    times = range(0.0, T, length = N) |> collect
     u1 = randn(1, N)
     u2 = randn(1, N)
-    p1 = LinearSplinePulse(u1, times; drive_name=:Ω)
-    p2 = LinearSplinePulse(u2, times; drive_name=:Ω)
+    p1 = LinearSplinePulse(u1, times; drive_name = :Ω)
+    p2 = LinearSplinePulse(u2, times; drive_name = :Ω)
 
     # α = 0 recovers p1
     p0 = interpolate_pulse(p1, p2, 0.0)
@@ -34,12 +35,13 @@ end
 @testitem "interpolate_pulse ZeroOrderPulse" begin
     using Intonato
 
-    N = 11; T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    N = 11;
+    T = 5.0
+    times = range(0.0, T, length = N) |> collect
     u1 = randn(1, N)
     u2 = randn(1, N)
-    p1 = ZeroOrderPulse(u1, times; drive_name=:Ω)
-    p2 = ZeroOrderPulse(u2, times; drive_name=:Ω)
+    p1 = ZeroOrderPulse(u1, times; drive_name = :Ω)
+    p2 = ZeroOrderPulse(u2, times; drive_name = :Ω)
 
     # α = 0 recovers p1
     p0 = interpolate_pulse(p1, p2, 0.0)
@@ -61,14 +63,15 @@ end
 @testitem "interpolate_pulse CubicSplinePulse" begin
     using Intonato
 
-    N = 11; T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    N = 11;
+    T = 5.0
+    times = range(0.0, T, length = N) |> collect
     u1 = randn(1, N)
     u2 = randn(1, N)
     du1 = randn(1, N)
     du2 = randn(1, N)
-    p1 = CubicSplinePulse(u1, du1, times; drive_name=:Ω)
-    p2 = CubicSplinePulse(u2, du2, times; drive_name=:Ω)
+    p1 = CubicSplinePulse(u1, du1, times; drive_name = :Ω)
+    p2 = CubicSplinePulse(u2, du2, times; drive_name = :Ω)
 
     # α = 0 recovers p1
     p0 = interpolate_pulse(p1, p2, 0.0)
@@ -97,10 +100,11 @@ end
 @testitem "truncate_pulse LinearSplinePulse at knot" begin
     using Intonato
 
-    N = 11; T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    N = 11;
+    T = 5.0
+    times = range(0.0, T, length = N) |> collect
     u = randn(1, N)
-    pulse = LinearSplinePulse(u, times; drive_name=:Ω)
+    pulse = LinearSplinePulse(u, times; drive_name = :Ω)
 
     # Truncate at the 6th knot (index 6, time = times[6])
     t_end = times[6]
@@ -118,10 +122,11 @@ end
 @testitem "truncate_pulse LinearSplinePulse between knots" begin
     using Intonato
 
-    N = 11; T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    N = 11;
+    T = 5.0
+    times = range(0.0, T, length = N) |> collect
     u = randn(1, N)
-    pulse = LinearSplinePulse(u, times; drive_name=:Ω)
+    pulse = LinearSplinePulse(u, times; drive_name = :Ω)
 
     # Truncate between knots 6 and 7
     t_end = (times[6] + times[7]) / 2.0
@@ -143,11 +148,12 @@ end
     using Intonato
     using ForwardDiff
 
-    N = 11; T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    N = 11;
+    T = 5.0
+    times = range(0.0, T, length = N) |> collect
     u = randn(1, N)
     du = randn(1, N)
-    pulse = CubicSplinePulse(u, du, times; drive_name=:Ω)
+    pulse = CubicSplinePulse(u, du, times; drive_name = :Ω)
 
     # Truncate between knots 6 and 7
     t_end = (times[6] + times[7]) / 2.0
diff --git a/src/pulse_ops/truncation.jl b/src/pulse_ops/truncation.jl
index bee22ea..3061d8e 100644
--- a/src/pulse_ops/truncation.jl
+++ b/src/pulse_ops/truncation.jl
@@ -12,7 +12,7 @@ function truncate_pulse(pulse::ZeroOrderPulse, t_end::Real)
     return ZeroOrderPulse(
         get_knot_values(pulse)[:, mask],
         times[mask];
-        drive_name=pulse.drive_name,
+        drive_name = pulse.drive_name,
     )
 end
 
@@ -28,7 +28,7 @@ function truncate_pulse(pulse::LinearSplinePulse, t_end::Real)
         t = vcat(t, Float64(t_end))
     end
 
-    return LinearSplinePulse(u, t; drive_name=pulse.drive_name)
+    return LinearSplinePulse(u, t; drive_name = pulse.drive_name)
 end
 
 function truncate_pulse(pulse::CubicSplinePulse, t_end::Real)
@@ -46,7 +46,7 @@ function truncate_pulse(pulse::CubicSplinePulse, t_end::Real)
         t = vcat(t, Float64(t_end))
     end
 
-    return CubicSplinePulse(u, du, t; drive_name=pulse.drive_name)
+    return CubicSplinePulse(u, du, t; drive_name = pulse.drive_name)
 end
 
 # ============================================================================ #
diff --git a/src/types/experiment_logger.jl b/src/types/experiment_logger.jl
index bb70e1b..910e9ab 100644
--- a/src/types/experiment_logger.jl
+++ b/src/types/experiment_logger.jl
@@ -70,8 +70,9 @@ record!(lg::InMemoryExperimentLogger, r::ExperimentRecord) = push!(lg.records, r
     σz = ComplexF64[1 0; 0 -1]
     sys = QuantumSystem(1.0 * σz, [σx], [1.0])
 
-    N = 11; T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    N = 11;
+    T = 5.0
+    times = range(0.0, T, length = N) |> collect
     ψ0 = ComplexF64[1.0, 0.0]
     ψg = ComplexF64[0.0, 1.0]
     pulse = LinearSplinePulse(0.1 * ones(1, N), times)
diff --git a/src/types/experiment_record.jl b/src/types/experiment_record.jl
index 47a1a42..7b95d1e 100644
--- a/src/types/experiment_record.jl
+++ b/src/types/experiment_record.jl
@@ -49,7 +49,7 @@ end
     # Single-qubit Rabi fixture (mirrors src/device_models/abstract.jl).
     N = 11
     T = 5.0
-    times = range(0.0, T, length=N) |> collect
+    times = range(0.0, T, length = N) |> collect
     somepulse = LinearSplinePulse(0.1 * ones(1, N), times)
     model = MeasurementModel(:ψ̃, [populations], [N])
     ms = [Measurement([0.7, 0.3], N)]
@@ -59,10 +59,17 @@ end
         measurement_model = model,
         measurements = ms,
         raw = nothing,
-        metadata = (; device="sim", shots=[1000, 1000], basis=:Z,
-                      timestamp="2026-06-14T00:00:00", config_hash="abc",
-                      calibration="nominal", seed=42, versions="Manifest:def",
-                      pulse_hash="ghi"),
+        metadata = (;
+            device = "sim",
+            shots = [1000, 1000],
+            basis = :Z,
+            timestamp = "2026-06-14T00:00:00",
+            config_hash = "abc",
+            calibration = "nominal",
+            seed = 42,
+            versions = "Manifest:def",
+            pulse_hash = "ghi",
+        ),
     )
 
     @test rec.pulse === somepulse
@@ -74,11 +81,7 @@ end
     @test rec.metadata.shots == [1000, 1000]
 
     # `raw` and `metadata` are optional (kwdef defaults).
-    rec2 = ExperimentRecord(
-        pulse = somepulse,
-        measurement_model = model,
-        measurements = ms,
-    )
+    rec2 = ExperimentRecord(pulse = somepulse, measurement_model = model, measurements = ms)
     @test rec2.raw === nothing
     @test rec2.metadata == (;)
 end
diff --git a/src/types/experiments.jl b/src/types/experiments.jl
index db1f76e..573a2b3 100644
--- a/src/types/experiments.jl
+++ b/src/types/experiments.jl
@@ -75,7 +75,7 @@ struct HardwareExperiment <: AbstractExperiment
     run::Function   # (pulse::AbstractPulse) → Vector{Measurement}
     # Optional measurement model used only to annotate logged ExperimentRecords.
     # `nothing` when the user's `run` closure owns the measurement structure.
-    measurement_model::Union{Nothing, MeasurementModel}
+    measurement_model::Union{Nothing,MeasurementModel}
 end
 
 HardwareExperiment(run::Function) = HardwareExperiment(run, nothing)
@@ -97,7 +97,8 @@ function run_experiment(
 )
     measurements = exp.run(pulse)
     if !(logger isa NullExperimentLogger)
-        model = isnothing(exp.measurement_model) ?
+        model =
+            isnothing(exp.measurement_model) ?
             _placeholder_measurement_model(measurements) : exp.measurement_model
         _maybe_record!(logger, exp, pulse, model, measurements; device = "hardware")
     end
@@ -194,23 +195,18 @@ fidelity using a single number-operator generator on the full space.
   authoritative value via a per-subsystem free-phase fidelity that grid-searches
   each independent phase generator separately.
 """
-function phase_max_fidelity(
-    ψ_T::AbstractVector,
-    ψ_goal::AbstractVector;
-    n_grid::Int = 128,
-)
+function phase_max_fidelity(ψ_T::AbstractVector, ψ_goal::AbstractVector; n_grid::Int = 128)
     d = length(ψ_T)
     @assert length(ψ_goal) == d "ψ_T and ψ_goal must have the same length"
     F_max = 0.0
-    @inbounds for i in 0:(n_grid - 1)
+    @inbounds for i = 0:(n_grid-1)
         φ = 2π * i / n_grid
         s = zero(ComplexF64)
-        for k in 0:(d - 1)
-            s += conj(cis(k * φ) * ψ_goal[k + 1]) * ψ_T[k + 1]
+        for k = 0:(d-1)
+            s += conj(cis(k * φ) * ψ_goal[k+1]) * ψ_T[k+1]
         end
         F = abs2(s)
         F > F_max && (F_max = F)
     end
     return F_max
 end
-
diff --git a/src/types/measurement_models.jl b/src/types/measurement_models.jl
index 8a4639a..8408e5a 100644
--- a/src/types/measurement_models.jl
+++ b/src/types/measurement_models.jl
@@ -48,10 +48,7 @@ Evaluate measurement functions on the model-side trajectory at each knot index.
 """
 function model_predict(traj::NamedTrajectory, model::MeasurementModel)
     return [
-        Measurement(
-            model.measurements[j](traj[k][model.state_name]),
-            k
-        )
-        for (j, k) in enumerate(model.indices)
+        Measurement(model.measurements[j](traj[k][model.state_name]), k) for
+        (j, k) in enumerate(model.indices)
     ]
 end
diff --git a/src/types/measurements.jl b/src/types/measurements.jl
index 5449d4f..683a407 100644
--- a/src/types/measurements.jl
+++ b/src/types/measurements.jl
@@ -20,5 +20,5 @@ Base.length(m::Measurement) = length(m.data)
 Sum of squared differences between two measurement collections.
 """
 function measurement_error(y1::Vector{Measurement}, y2::Vector{Measurement})
-    return sum(sum((m1.data .- m2.data).^2) for (m1, m2) in zip(y1, y2))
+    return sum(sum((m1.data .- m2.data) .^ 2) for (m1, m2) in zip(y1, y2))
 end
diff --git a/src/types/noise_models.jl b/src/types/noise_models.jl
index 8a80cb2..d328149 100644
--- a/src/types/noise_models.jl
+++ b/src/types/noise_models.jl
@@ -64,7 +64,7 @@ end
 Covariance for Pauli expectation measurements `⟨σ_j⟩`.
 `Var[⟨σ_j⟩] = (1 - ⟨σ_j⟩²) / n_shots` (since σ_j² = I for Pauli matrices).
 """
-pauli_covariance(y, n) = Diagonal((1 .- y.^2) ./ n)
+pauli_covariance(y, n) = Diagonal((1 .- y .^ 2) ./ n)
 
 """
     population_covariance(y, n_shots) → Matrix
@@ -81,7 +81,7 @@ Covariance for Wigner / displaced parity measurements.
 `Var[W(α)] = (1 - W(α)²) / n_shots` (Rademacher mean).
 Same formula as Pauli, different physics.
 """
-wigner_covariance(y, n) = Diagonal((1 .- y.^2) ./ n)
+wigner_covariance(y, n) = Diagonal((1 .- y .^ 2) ./ n)
 
 # ============================================================================ #
 #                      Measurement presets
@@ -93,9 +93,10 @@ wigner_covariance(y, n) = Diagonal((1 .- y.^2) ./ n)
 Pauli expectation measurement preset. Returns `DeterministicMeasurement` when
 `n_shots` is omitted, `ShotNoiseMeasurement` with `pauli_covariance` otherwise.
 """
-function pauli(ops; n_shots::Union{Nothing, Int} = nothing)
+function pauli(ops; n_shots::Union{Nothing,Int} = nothing)
     g = expect(ops)
-    isnothing(n_shots) ? DeterministicMeasurement(g) : ShotNoiseMeasurement(g, n_shots, pauli_covariance)
+    isnothing(n_shots) ? DeterministicMeasurement(g) :
+    ShotNoiseMeasurement(g, n_shots, pauli_covariance)
 end
 
 """
@@ -104,8 +105,9 @@ end
 Population measurement preset. Returns `DeterministicMeasurement` when
 `n_shots` is omitted, `ShotNoiseMeasurement` with `population_covariance` otherwise.
 """
-function pop(; n_shots::Union{Nothing, Int} = nothing)
-    isnothing(n_shots) ? DeterministicMeasurement(populations) : ShotNoiseMeasurement(populations, n_shots, population_covariance)
+function pop(; n_shots::Union{Nothing,Int} = nothing)
+    isnothing(n_shots) ? DeterministicMeasurement(populations) :
+    ShotNoiseMeasurement(populations, n_shots, population_covariance)
 end
 
 # ============================================================================ #
diff --git a/src/types/test.jl b/src/types/test.jl
index e2b32ba..54b876a 100644
--- a/src/types/test.jl
+++ b/src/types/test.jl
@@ -15,7 +15,7 @@
     N = 11
 
     # New API: pass AbstractMeasurement directly
-    m = pauli([σx_iso, σy_iso, σz_iso]; n_shots=1000)
+    m = pauli([σx_iso, σy_iso, σz_iso]; n_shots = 1000)
     model = MeasurementModel(:ψ̃, [m], [N])
     @test model.measurements[1] isa ShotNoiseMeasurement
     @test length(model.measurements) == 1
@@ -29,7 +29,8 @@
     ψ_iso = Intonato.ket_to_iso(ComplexF64[1.0, 0.0])
     traj = NamedTrajectory(
         (ψ̃ = repeat(ψ_iso, 1, N), u = zeros(1, N), Δt = fill(0.1, 1, N));
-        timestep=:Δt, controls=(:u,)
+        timestep = :Δt,
+        controls = (:u,),
     )
     y = model_predict(traj, model_compat)
     @test length(y) == 1
@@ -51,15 +52,15 @@ end
 
     @test Σ isa Diagonal
     @test size(Σ) == (3, 3)
-    @test Σ[1,1] ≈ (1 - 0.25) / 1000   # (1 - 0.5²) / 1000
-    @test Σ[2,2] ≈ (1 - 0.09) / 1000   # (1 - 0.3²) / 1000
-    @test Σ[3,3] ≈ (1 - 0.64) / 1000   # (1 - 0.8²) / 1000
+    @test Σ[1, 1] ≈ (1 - 0.25) / 1000   # (1 - 0.5²) / 1000
+    @test Σ[2, 2] ≈ (1 - 0.09) / 1000   # (1 - 0.3²) / 1000
+    @test Σ[3, 3] ≈ (1 - 0.64) / 1000   # (1 - 0.8²) / 1000
 
     # Edge: y = ±1 → variance → 0
     Σ_edge = pauli_covariance([1.0, -1.0, 0.0], 100)
-    @test Σ_edge[1,1] ≈ 0.0 atol=1e-15
-    @test Σ_edge[2,2] ≈ 0.0 atol=1e-15
-    @test Σ_edge[3,3] ≈ 1/100  # maximum variance at equator
+    @test Σ_edge[1, 1] ≈ 0.0 atol=1e-15
+    @test Σ_edge[2, 2] ≈ 0.0 atol=1e-15
+    @test Σ_edge[3, 3] ≈ 1/100  # maximum variance at equator
 end
 
 @testitem "population_covariance formula" begin
@@ -72,9 +73,9 @@ end
     Σ = population_covariance(p, n)
 
     @test size(Σ) == (2, 2)
-    @test Σ[1,1] ≈ 0.7 * 0.3 / 500     # p(1-p)/N
-    @test Σ[2,2] ≈ 0.3 * 0.7 / 500
-    @test Σ[1,2] ≈ -0.7 * 0.3 / 500     # -p_j*p_k/N
+    @test Σ[1, 1] ≈ 0.7 * 0.3 / 500     # p(1-p)/N
+    @test Σ[2, 2] ≈ 0.3 * 0.7 / 500
+    @test Σ[1, 2] ≈ -0.7 * 0.3 / 500     # -p_j*p_k/N
     @test issymmetric(Σ)
 end
 
@@ -92,7 +93,7 @@ end
     sn = ShotNoiseMeasurement(g, 100, pauli_covariance)
     @test sn(x) == [1.0, 4.0, 9.0]
 
-    kc = KnownCovarianceMeasurement(g, zeros(3,3))
+    kc = KnownCovarianceMeasurement(g, zeros(3, 3))
     @test kc(x) == [1.0, 4.0, 9.0]
 end
 
@@ -110,7 +111,7 @@ end
     @test m_det isa DeterministicMeasurement
 
     # pauli with n_shots → ShotNoise
-    m_sn = pauli(ops; n_shots=1000)
+    m_sn = pauli(ops; n_shots = 1000)
     @test m_sn isa ShotNoiseMeasurement
     @test m_sn.n_shots == 1000
 
@@ -121,7 +122,7 @@ end
     # pop presets
     p_det = pop()
     @test p_det isa DeterministicMeasurement
-    p_sn = pop(; n_shots=500)
+    p_sn = pop(; n_shots = 500)
     @test p_sn isa ShotNoiseMeasurement
     @test p_sn.n_shots == 500
 end
@@ -138,7 +139,7 @@ end
     # number states with DIFFERENT n. Use a superposition so a *relative* phase
     # on |1⟩ degrades the raw overlap but is exactly removable by the free phase.
     ψ_goal = ComplexF64[1.0, 1.0] / sqrt(2)             # (|0⟩+|1⟩)/√2
-    ψ_T    = ComplexF64[1.0, cis(0.7)] / sqrt(2)        # (|0⟩+e^{i·0.7}|1⟩)/√2
+    ψ_T = ComplexF64[1.0, cis(0.7)] / sqrt(2)        # (|0⟩+e^{i·0.7}|1⟩)/√2
 
     raw = abs2(dot(ψ_goal, ψ_T))
     F = phase_max_fidelity(ψ_T, ψ_goal)
@@ -152,4 +153,3 @@ end
     # A genuine population mismatch can't be fixed by a phase: F < 1.
     @test phase_max_fidelity(ComplexF64[0.0, 1.0], ψ_goal) < 0.75
 end
-