Skip to content

More DMStag: local_to_local!, with_product_coordinates, set_matrix_preallocate_only! - #296

Merged
filoferra merged 2 commits into
JuliaParallel:mainfrom
filoferra:dmstag-halos
Sep 28, 2026
Merged

filoferra merged 2 commits into
JuliaParallel:mainfrom
filoferra:dmstag-halos

Conversation

@filoferra

Copy link
Copy Markdown
Collaborator

Stacked on #295 : review that one first. Until it merges, this diff includes its commit.

The rest of what DMStag usersstill might do through LibPETSc:

  • local_to_local!(dst, dm, src) and the in-place local_to_local!(v, dm) refresh the ghost points of a local vector from its neighbours. Tested on 4 ranks.
  • with_product_coordinates(f, dm) gives a block read-only access to a DMStag's per-axis coordinates (x[i, 1] is the lower face of element i, x[i, 2] its centre) and hands them back when the block returns or throws.
  • set_matrix_preallocate_only!(dm, true) makes a matrix from dm get its nonzero pattern from the first assembly, instead of storing every coupling the stencil allows as a zero. It's a setter rather than a PetscMat keyword because the flag stays on the DM.

Fixes:

  • LibPETSc.DMStagRestoreProductCoordinateArraysRead called the writable restore, which doesn't match the read-only get.
  • The KSP and SNES docstrings now say when the constructor's options are applied: once at construction for KSP, at every solve! for SNES (and TS, whose docstring already said so).

@filoferra

Copy link
Copy Markdown
Collaborator Author

Added a third piece to this PR: with_field_views!(f, dm, vecs...; fields, read, write).

It checks out DMStag local vectors and gives the block one view per field (location => dof), indexed by element like stencil, ghosts included:

flow = (face_location(dm, 1) => 0, face_location(dm, 2) => 0, element_location(dm) => 0)
with_field_views!(dm, x, r; fields = flow, write = (false, true)) do (Vx, Vy, P), (Rx, Ry, Rp)
    ...
end

Without fields, the block gets the whole [I..., slot] array. read/write work as in with_local_array!.

  • The views are concretely typed, so user code no longer has to spell out PETSc.jl's internal array types.
  • A checkout costs 9 allocations per vector, the same at 8² and 64².
  • An exception in the block hands every vector back. A global vector, or one from another DM, throws ArgumentError.

Tests: test/dmstag_views.jl.

local_to_local!(dst, dm, src) and the in-place local_to_local!(v, dm) refresh ghost points. with_product_coordinates(f, dm) hands out a DMStag's per-axis coordinates read-only. set_matrix_preallocate_only! lets the first assembly define a DM matrix's pattern. The Read restore of the product coordinates called the writable restore. The KSP and SNES docstrings say when constructor options are applied.
with_field_views!(f, dm, vecs...; fields, read, write) checks out local vectors and hands f one view per location => dof pair, indexed like stencil, ghosts included, or the whole [I..., slot] array without fields.

The views are concretely typed; a checkout costs 9 allocations per vector (independent of the grid size)
@filoferra
filoferra marked this pull request as ready for review September 28, 2026 08:04
@filoferra
filoferra merged commit 8ebdb4b into JuliaParallel:main Sep 28, 2026
6 checks passed
@filoferra
filoferra deleted the dmstag-halos branch September 28, 2026 08:06
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant