Skip to content

Use pyamrex in spacecraft charging example #4817

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 18 commits into from
Apr 18, 2025

Conversation

RemiLehe
Copy link
Member

@RemiLehe RemiLehe commented Apr 2, 2024

Overview

This PR modernizes the spacecraft charging example to use the newer pyamrex interface for multifabs, instead of the fields wrapper.

The pyamrex interface is in principle more flexible (since it gives access to the AMReX functions that operate on MultiFabs and should have all the functionalities that the old field wrappers had. (See also AMReX-Codes/pyamrex#370). This is also the interface that is documented in the WarpX documentation

Resetting the checksum

I had to reset the checksum, as they differed by ~1e-2, which is not negligible!

Manually looking at the fields seem to indicate that the old interface and new interface differ in how they were handling the guard cells (at least in the way I use them in this test). While I don't fully understand the differences, this PR seems to produce more correct results.
(In hindsight, it is striking that we missed the unphysical Ez field when we implemented this test. The reason might be that we were focusing phi, which does not exhibit this unphysical feature in such an obvious way.)

Screenshot 2025-04-18 at 8 55 31 AM

@RemiLehe RemiLehe marked this pull request as draft April 2, 2024 21:10
@ax3l ax3l added component: tests Tests and CI component: Python Python layer labels Apr 3, 2024
@RemiLehe RemiLehe changed the title [WIP] Use pyamrex in spacecraft charging example Use pyamrex in spacecraft charging example Sep 23, 2024
@RemiLehe RemiLehe marked this pull request as ready for review September 23, 2024 23:06
@ax3l ax3l self-requested a review September 23, 2024 23:09
@ax3l ax3l self-assigned this Sep 23, 2024
Copy link
Member

@ax3l ax3l left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

LGTM! Ran locally for me, minus checksum stability.

This still uses the global API from Dave (but now from pyAMReX); if we want to parallelize it further we could do another follow-up :)

Use the one that does not rely on naming conventions but has
named arguments. (The other one worked as well, but is deprecated.)
@ax3l ax3l enabled auto-merge (squash) January 26, 2025 01:15
@ax3l ax3l disabled auto-merge January 26, 2025 17:54
ax3l added a commit to AMReX-Codes/pyamrex that referenced this pull request Jan 29, 2025
We already have `.size` and we can iterate the `FabArray`, so we should
support `len(mf)`, too.

Seen in BLAST-WarpX/warpx#4817
2 * len(phi)
nr = phi.shape[0]
r = np.linspace(rmin, rmax, nr, endpoint=False) + (rmax - rmin) / (
2 * nr
) # shift of the r points because the derivaties are calculated in the middle
face_z0 = (
2 * np.pi * 1.0 / dz * ((phi[:, 0] - phi[:, 1]) * r).sum() * dr
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks correct/as before: the : is without ghost cells, () is with ghost cells.
https://pyamrex.readthedocs.io/en/latest/usage/compute.html

Comment on lines +81 to +84
self.normalized_Er = warpx.multifab("Efield_fp", dir=self.dir_r, level=0).copy()
self.normalized_Er.mult(1 / q_v, 0)
self.normalized_Ez = warpx.multifab("Efield_fp", dir=self.dir_z, level=0).copy()
self.normalized_Ez.mult(1 / q_v, 0)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is now perfectly parallel in MPI contexts, no global comms.

phi = PhiFPWrapper(include_ghosts=True)
phi[...] += (q - q_v) * self.normalized_phi
Er = warpx.multifab("Efield_fp", dir=self.dir_r, level=0)
Er.saxpy(Er, (q - q_v), self.normalized_Er, 0, 0, 1, 0)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is now perfectly parallel in MPI contexts, no global comms.

@RemiLehe RemiLehe assigned EZoni and unassigned ax3l Apr 18, 2025
r = np.linspace(rmin, rmax, len(phi), endpoint=False) + (rmax - rmin) / (
2 * len(phi)
nr = phi.shape[0]
r = np.linspace(rmin, rmax, nr, endpoint=False) + (rmax - rmin) / (
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This works, but I had to stare at it for a moment to figure out what it is doing. I'll often do

nr = phi.shape[0] - 1
dr = (rmax - rmin)/nr
r = np.linspace(rmin + dr/2., rmax - dr/2., nr)

Copy link
Member

@dpgrote dpgrote left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This looks good to me!

@RemiLehe RemiLehe merged commit 428d0e9 into BLAST-WarpX:development Apr 18, 2025
37 checks passed
@RemiLehe RemiLehe deleted the use_pyamrex branch April 18, 2025 16:38
atmyers pushed a commit to atmyers/WarpX that referenced this pull request Jul 3, 2025
# Overview

This PR modernizes the spacecraft charging example to use the newer
`pyamrex` interface for multifabs, instead of the fields wrapper.

The `pyamrex` interface is in principle more flexible (since it gives
access to the AMReX functions that operate on `MultiFabs` and should
have all the functionalities that the old field wrappers had. (See also
AMReX-Codes/pyamrex#370). This is also the
interface that is documented in the [WarpX
documentation](https://warpx.readthedocs.io/en/latest/usage/workflows/python_extend.html#fields)

# Resetting the checksum

I had to reset the checksum, as they differed by `~1e-2`, which is not
negligible!

Manually looking at the fields seem to indicate that the old interface
and new interface differ in how they were handling the guard cells (at
least in the way I use them in this test). While I don't fully
understand the differences, this PR seems to produce more correct
results.
(In hindsight, it is striking that we missed the unphysical `Ez` field
when we implemented this test. The reason might be that we were focusing
`phi`, which does not exhibit this unphysical feature in such an obvious
way.)

<img width="1251" alt="Screenshot 2025-04-18 at 8 55 31 AM"
src="https://github.com/user-attachments/assets/683235aa-c03b-425c-8b55-d98016d77bdb"
/>

---------

Co-authored-by: Axel Huebl <[email protected]>
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
component: Python Python layer component: tests Tests and CI
Projects
None yet
Development

Successfully merging this pull request may close these issues.

4 participants