Use scratch arrays in CG, coupled CG, and GMRES on the CPU - #2754
Use scratch arrays in CG, coupled CG, and GMRES on the CPU#2754timofeymukha wants to merge 3 commits into
Conversation
There was a problem hiding this comment.
🟡 Changes recommended
The new scratch-array requests pass uninitialised pointers/indices to request_host_array dummy arguments declared intent(inout), which is nonconforming and can fail under runtime/debug checks.
Once you've addressed the issues Copilot identified, you can request another Copilot review.
Pull request overview
This PR reduces CPU-side memory overhead for Krylov solvers by switching GMRES, CG, and coupled CG to borrow their large work arrays from the global scratch registry during solve() calls, instead of retaining dedicated solver-owned storage.
Changes:
- Replaced persistent allocatable work arrays in CPU
cg_t,cg_cpld_t, andgmres_twith pointers that are associated with scratch-registryhost_array_tstorage for the duration of each solve. - Updated build dependencies to include
registries/scratch_registryanddata_types/host_arrayfor the affected CPU Krylov backends. - Documented the change in
CHANGELOG.mdand improved/expanded solver docstrings.
File summaries
| File | Description |
|---|---|
| src/krylov/bcknd/cpu/gmres.f90 | Switch GMRES work vectors/bases to scratch-registry host arrays during solve(). |
| src/krylov/bcknd/cpu/cg.f90 | Switch CG workspace (w/r/p/z/alpha) to scratch-registry host arrays during solve(). |
| src/krylov/bcknd/cpu/cg_coupled.f90 | Switch coupled CG component work vectors to scratch-registry host arrays during solve_coupled(). |
| src/.depends | Add scratch-registry and host-array dependencies for the updated CPU Krylov objects. |
| CHANGELOG.md | Note the CPU Krylov workspace change in the Develop section. |
Review details
- Files reviewed: 5/5 changed files
- Comments generated: 3
- Review effort level: Lite
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
njansson
left a comment
There was a problem hiding this comment.
Same as for bicgstab, can we keep the pointers with the same name in the type, it eats some bytes but greatly improves readability for linear solver people looking at the code
We do! |
Co-authored-by: Niclas Jansson <njansson@kth.se>
| call neko_scratch_registry%request_host_array(w_tmp, temp_indices(1), & | ||
| n, .false.) | ||
| call neko_scratch_registry%request_host_array(r_tmp, temp_indices(2), & | ||
| n, .false.) | ||
| call neko_scratch_registry%request_host_array(p_tmp, temp_indices(3), & | ||
| n * CG_P_SPACE, .false.) | ||
| call neko_scratch_registry%request_host_array(z_tmp, temp_indices(4), & | ||
| n, .false.) | ||
| call neko_scratch_registry%request_host_array(alpha_tmp, & | ||
| temp_indices(5), CG_P_SPACE, .false.) | ||
|
|
||
| this%w => w_tmp%x | ||
| this%r => r_tmp%x | ||
| this%p(1:n, 1:CG_P_SPACE) => p_tmp%x | ||
| this%z => z_tmp%x | ||
| this%alpha => alpha_tmp%x | ||
|
|
||
| associate(w => this%w, r => this%r, p => this%p, & | ||
| z => this%z, alpha => this%alpha) |
There was a problem hiding this comment.
| call neko_scratch_registry%request_host_array(w_tmp, temp_indices(1), & | |
| n, .false.) | |
| call neko_scratch_registry%request_host_array(r_tmp, temp_indices(2), & | |
| n, .false.) | |
| call neko_scratch_registry%request_host_array(p_tmp, temp_indices(3), & | |
| n * CG_P_SPACE, .false.) | |
| call neko_scratch_registry%request_host_array(z_tmp, temp_indices(4), & | |
| n, .false.) | |
| call neko_scratch_registry%request_host_array(alpha_tmp, & | |
| temp_indices(5), CG_P_SPACE, .false.) | |
| this%w => w_tmp%x | |
| this%r => r_tmp%x | |
| this%p(1:n, 1:CG_P_SPACE) => p_tmp%x | |
| this%z => z_tmp%x | |
| this%alpha => alpha_tmp%x | |
| associate(w => this%w, r => this%r, p => this%p, & | |
| z => this%z, alpha => this%alpha) | |
| call neko_scratch_registry%request(this%w, temp_indices(1), n, .false.) | |
| call neko_scratch_registry%request(this%r, temp_indices(2), n, .false.) | |
| call neko_scratch_registry%request(this%z, temp_indices(4), n, .false.) | |
| call neko_scratch_registry%request(this%alpha, temp_indices(5), CG_P_SPACE, .false.) | |
| call neko_scratch_registry%request(p_tmp, temp_indices(3), & | |
| n * CG_P_SPACE, .false.) | |
| this%p(1:n, 1:CG_P_SPACE) => p_tmp%x | |
| associate(w => this%w%x, r => this%r%x, p => this%p%x, & | |
| z => this%z%x, alpha => this%alpha%x) |
timfelle
left a comment
There was a problem hiding this comment.
This can be done much cleaner with out having the subroutine temporaries if you use the host_array_t in the class.
|
|
||
| this%w => w_tmp%x | ||
| this%r => r_tmp%x | ||
| this%p(1:n, 1:CG_P_SPACE) => p_tmp%x |
There was a problem hiding this comment.
This makes me a bit nervous. Can't we add a rank-2 registry as well (@timfelle )?
There was a problem hiding this comment.
Array reshaping like this seem to be completely standard in the fortran community, But sure we can easily add another host_array_t to the system. But explicitly for a host implementation the only difference between host_array_t and vector_t is a %x_d=C_NULL_PTR, so I actually don't think the host array thing makes a difference at all since it should never be active in a device build.
There was a problem hiding this comment.
So my point is, you should probably just use a matrix_t here.
There was a problem hiding this comment.
Array reshaping like this seem to be completely standard in the fortran community, But sure we can easily add another host_array_t to the system. But explicitly for a host implementation the only difference between host_array_t and vector_t is a %x_d=C_NULL_PTR, so I actually don't think the host array thing makes a difference at all since it should never be active in a device build.
standard, and working on compilers are different things, as we know ;) But sure let's give it a try
There was a problem hiding this comment.
Sure, but we do it all over the place already, we just do it at subroutine interfaces, with out of bound indexing and such.
Should remove quite a bit of memory overhead on the CPU, since CG and GMRES are typically used in tandem.
The types retain pointers (instead of persistent storage) that are redirected to scratch arrays in the solve routine.