Skip to content

Pass Dirichlet boundary conditions as std::span in the assembly API #4436

Description

@jhale

AGENTS.md asks for std::span for contiguous read-only array views, and reserves std::vector for parameters that are stored, mutated in place, or are container element types. The bcs parameter of the assembly API is none of those, but is spelled as a vector everywhere:

  • fem::assemble_matrixcpp/dolfinx/fem/assembler.h:543, :590
  • fem::set_diagonalcpp/dolfinx/fem/assembler.h:672
  • fem::petsc::set_bccpp/dolfinx/fem/petsc.h:577
  • fem::petsc::apply_liftingcpp/dolfinx/fem/petsc.h:407, :500 (nested, vector<vector<...>>)

all taking const std::vector<std::reference_wrapper<const DirichletBC<T, U>>>&.

Proposal

std::span<const std::reference_wrapper<const DirichletBC<T, U>>> bcs

The outer const makes the elements non-assignable, the inner one keeps the boundary conditions read-only. It binds implicitly from the vectors callers already build, so most call sites are unaffected.

Why it matters beyond style

A function that accepts boundary conditions and forwards them has to spell the same vector type, or allocate on every call. fem::petsc::assemble_residual and assemble_jacobian (added in #4433, cpp/dolfinx/fem/petsc.h:713, :789) take const std::vector<std::reference_wrapper<...>>& purely to forward to assemble_matrix, set_diagonal, set_bc and apply_lifting. They are called once per residual or Jacobian evaluation, so once per Newton step and once per line search trial point. With spans throughout, nothing converts anywhere.

Caveats

  • std::span is not constructible from an initializer_list until C++26 (P2447), so braced call sites need a named local. fem::apply_lifting(b.array(), {a}, {{bc}}, {}, T(1)) in the poisson, biharmonic and mixed_poisson demos is the pattern affected. An empty {}, as in assemble_matrix(A.mat_add_values(), a, {}), still works.
  • The nested vector<vector<...>> of apply_lifting is more disruptive than the flat case, since a span<const span<...>> needs the caller to own the inner sequences. It could keep its current type, or be handled separately.
  • The nanobind wrappers build vectors from Python lists and would bind unchanged, provided the vector outlives the call.

There are around 60 call sites of these functions in the repository.

AI assistance: this issue was drafted with Claude Code (Opus 5) while reviewing #4433. I reviewed and take responsibility for its contents.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions