Cuda pr 9 kernel folders - #663
Draft
max-models wants to merge 84 commits into
Draft
max-models wants to merge 84 commits into
max-models wants to merge 84 commits into
Conversation
The SPH linear smoothing kernels returned a non-zero gradient at zero separation, so every particle pushed on itself. This PR makes all SPH kernel gradients vanish at r=0, which matches the symmetric value the trigonometric and gaussian kernels already give. Closes #437 ### What was wrong - `grad_linear_uni(0, h)` returned `+1/h**2`, the left derivative. That value carries into `grad_linear_1d`, `grad_linear_2d_{1,2}` and `grad_linear_3d_{1,2,3}`, and it applies whenever the component's own coordinate is 0, not only at the origin. - `grad_linear_isotropic_3d_{1,2,3}` returned `-1/h / (C h^3)` at `r == 0`. - As a result the SPH gradient of a lattice-loaded constant density was not zero: it came out at about 3.5 for 1d linear, 3.8 for 3d tensor linear and -0.58 for isotropic in the check below. ### What changed - `src/struphy/pic/sph_smoothing_kernels.py`: at the cusp the gradient now returns 0 (`x == 0` in `grad_linear_uni`, `r == 0` in the isotropic variants). Values away from 0 are unchanged. The trigonometric and gaussian gradients were checked too and already give 0 there. - `src/struphy/pic/tests/test_kernel_setup.py`: - New test `test_sph_kernel_gradients_vanish_at_zero` covers every gradient type in `smoothing_kernel`, at the origin and on the coordinate planes. - `test_sph_tensor_destinations` relied on the old self-gradient `1/h**2`. It now uses a second neighbouring particle so the viscosity-tensor destinations still get non-trivial values. ### Verification - A standalone check calls all 21 gradient kernel types at the origin and on the coordinate planes. Before the fix, 9 types were non-zero at the origin and all 4 plane cases were non-zero. After the fix all of them are 0. At random points away from 0 the values are bit-identical before and after. - The gradient of a lattice-loaded constant density is now ~1e-16 in 1d linear, 3d tensor linear and 3d isotropic. - `pytest src/struphy/pic/tests/test_kernel_setup.py`: 74 passed. The new test fails without the fix. - A subset of `test_sph.py::test_sph_evaluation_1d` (linear_1d, periodic) passes. ### Risks / follow-ups - SPH results that use the linear kernels change slightly, because the self-contribution is gone. That removal is the intended correction. - #436 (viscosity strain without DF^-T) is in the same area and is not addressed here. - This branch includes the commit bumping `feectools` to 0.1.11 so that CI can run. 🤖 Generated with [Claude Code](https://claude.com/claude-code) --------- Co-authored-by: Claude Opus 5.5 <noreply@anthropic.com>
…ld just be a proof of concept
CudaKernel no longer flattens argument objects or casts scalars at call time: the arguments must already be in the format of a cupy.RawKernel (flat tuple of CuPy arrays and NumPy scalars, built once at setup), and the number of threads is passed as n_threads. Add docstrings to the kernel and argument classes. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Converting the scalars and joining the .values of the CUDA argument classes costs about 1 us per call (the kernel launch alone about 70 us), prevents silently misaligned arguments when a Python int is passed as a 64-bit value, and keeps the calls the same on both backends. Arrays are still never converted or copied; n_threads stays an explicit argument. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
…sses CUDA reads each kernel argument with the size declared in the signature, so Python int/float already arrive correctly in int/double parameters; the casts did not change the result, and did not prevent the actual failure cases (e.g. an integer passed to a double parameter). Checking scalars against the kernel signature is noted as a follow-up in CUDA_STRATEGY.md. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Add CudaKernel.from_file, which reads the CUDA source from a <name>_cuda.cu file and takes the kernel name from the file name, and ship .cu/.cuh files as package data. Step 2 of CUDA_STRATEGY.md. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
KernelCatalog.from_package pairs <name>/<name>_kernels.py (pyccel) with <name>/<name>_cuda.cu (CUDA, optional). The CUDA kernel of a Kernel is now optional: on the CuPy backend, a kernel without a CUDA version raises NotImplementedError naming the expected .cu file. Step 3 of CUDA_STRATEGY.md. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Pusher now takes a Kernel or, as before, a PyccelKernel (wrapped into a Kernel without CUDA version) and calls get_kernel() in its constructor. On the CuPy backend, a pusher whose kernel has no CUDA version thus fails when it is created instead of in the time loop. No changes to the propagators, no behaviour change on the CPU. Step 4 of CUDA_STRATEGY.md. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
- Test only analytic mappings on the CuPy backend: spline mappings such as IGAPolarCylinder cannot be created there yet (interp_mapping passes CuPy arrays to scipy.sparse.csc_matrix); correct CUDA_STRATEGY.md. - Build the CUDA domain arguments with cupy instead of xp, so they can be built whichever backend is active (the arrays are on the device); test this. - Move the new tests above the __main__ block of test_domain.py. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
# Conflicts: # CUDA_STRATEGY.md
…o cuda-pr-7-derham
The kernel for the active backend is chosen once, when the setup is created, as in Pusher; on the CuPy backend a kernel without CUDA version fails there. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
The 43 kernels of pusher_kernels.py, pusher_kernels_gc.py, pusher_kernels_sph.py, eval_kernels_gc.py and eval_kernels_sph.py are copied, with unchanged function bodies, to kernels/<name>/<name>_kernels.py; kernels/__init__.py collects them in a KernelCatalog. Each module imports only what its kernel uses, plus the pusher_args_kernels module that dependencies.py needs. push_gc_bxEstar_discrete_gradient_1st_order_newton is renamed to push_gc_bxEstar_dg_1st_order_newton: pyccel's Fortran wrapper module bind_c_<name>_kernels must fit Fortran's 63-character limit for names. test_pushing_catalog checks the layout and the name length. The old modules are still used; the call sites switch in the next commit. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
The propagators, Particles and the tests take their kernels from struphy.pic.pushing.kernels.catalog; the old modules pusher_kernels.py, pusher_kernels_gc.py, pusher_kernels_sph.py, eval_kernels_gc.py and eval_kernels_sph.py are removed. No re-export modules: their names would have to contain "kernels", so struphy compile would try to compile them. CurrentCoupling5DGradB calls its kernels itself and selects them once with get_kernel(); Particles runs the SPH evaluation kernels on its host bundle with catalog[...].pyccel_kernel, as before. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
9 of 20 tasks
…uphy into cuda-pr-9-kernel-folders
This branch has not been deployed
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
No description provided.