[SPH] Add MHD block to Phantom I/O#1889
Conversation
…licate sinks being added once per range chunk
|
Thanks @JackNarvaez for opening this PR! You can do multiple things directly here: Once the workflow completes a message will appear displaying informations related to the run. Also the PR gets automatically reviewed by gemini, you can: |
![gemini-code-assist[bot]](https://avatars.githubusercontent.com/in/956858?s=60&v=4){kind=small}
There was a problem hiding this comment.
Code Review
This pull request introduces Magnetohydrodynamics (MHD) support to the SPH model, enabling the loading, processing, and writing of MHD fields (such as magnetic fields and psi/ch) from and to Phantom dump files. Key changes include adding MHD configuration helper functions, loading MHD fields during model initialization, and appending MHD data blocks when generating Phantom dumps. The review feedback highlights two main improvements: validating the sizes of the loaded MHD vector components to prevent potential out-of-bounds memory access, and removing an unnecessary empty third block (block_3rd) when writing the Phantom dump to avoid compatibility issues with external tools.
Important
The consumer version of Gemini Code Assist on GitHub is being sunset. Starting June 18, 2026, new organization installations will be blocked, and all code review activity will officially cease on July 17, 2026.
For more details on the timeline and next steps, please review the Help Documentation.
| for (u32 i = 0; i < Brhox.size(); i++) { | ||
| Brhoxyz.push_back({Brhox[i], Brhoy[i], Brhoz[i]}); | ||
| } |
There was a problem hiding this comment.
Mismatched vector sizes can lead to out-of-bounds memory access and undefined behavior/crashes. We should validate that:
- The individual components of the magnetic field (
Brhox,Brhoy,Brhoz) have the same size before constructingBrhoxyz. - The constructed
Brhoxyzandpsichvectors have the same size as the coordinate vectorxyz(if they are populated), ensuring that subsequent indexing usingsel_indexis safe.
if (Brhox.size() != Brhoy.size() || Brhox.size() != Brhoz.size()) {
shambase::throw_with_loc<std::runtime_error>(
"MHD B/rho component vector sizes do not match");
}
for (u32 i = 0; i < Brhox.size(); i++) {
Brhoxyz.push_back({Brhox[i], Brhoy[i], Brhoz[i]});
}
if (!Brhoxyz.empty() && Brhoxyz.size() != xyz.size()) {
shambase::throw_with_loc<std::runtime_error>(
"MHD B/rho vector size does not match particle count");
}
if (!psich.empty() && psich.size() != xyz.size()) {
shambase::throw_with_loc<std::runtime_error>(
"MHD psi/ch vector size does not match particle count");
}| PhantomDumpBlock block_3rd; | ||
| dump.blocks.push_back(std::move(block_3rd)); | ||
| dump.blocks.push_back(std::move(block_mhd)); |
There was a problem hiding this comment.
An empty block_3rd is written to the Phantom dump file. Since ndust is 0, there is no dust data, and writing an empty block is unnecessary. Furthermore, init_from_phantom_dump searches for MHD fields by looping over all blocks, so the exact block index of block_mhd does not need to be padded. Writing an empty block may cause compatibility issues with other tools (like Phantom or Splash) that expect a specific block structure based on the header parameters.
| PhantomDumpBlock block_3rd; | |
| dump.blocks.push_back(std::move(block_3rd)); | |
| dump.blocks.push_back(std::move(block_mhd)); | |
| dump.blocks.push_back(std::move(block_mhd)); |
y-lapeyre
changed the title
for more information, see https://pre-commit.ci
|
Hi @JackNarvaez thanks for the PR ! |
|
Hi @tdavidcl |
y-lapeyre
left a comment
y-lapeyre
left a comment
There was a problem hiding this comment.
Hi! Thanks for the PR. Currently does not work on a simple load test
filename = "mhdwave_00000"
import shamrock
dump = shamrock.load_phantom_dump(filename)
dump.print_state()libc++abi: terminating due to uncaught exception of type pybind11::error_already_set: RuntimeError: the byte count is not correctI'll investigate now :)
Codecov Report❌ Patch coverage is 📢 Thoughts on this report? Let us know! |
y-lapeyre
left a comment
y-lapeyre
left a comment
There was a problem hiding this comment.
Looks good to me. @JackNarvaez Could you add a file named run_mhd_open_phnatom_dump.py in examples/sph
"""
Open a phantom dump
============================
This simple example shows how to open a phantom dump in shamrock,
including the MHD fields
"""
# %%
# Download a phantom dump
dump_folder = "_to_trash"
import os
os.system("mkdir -p " + dump_folder)
url = "https://raw.githubusercontent.com/Shamrock-code/reference-files/refs/heads/main/briowu_00000"
filename = dump_folder + "/briowu_00000"
from urllib.request import urlretrieve
urlretrieve(url, filename)
# %%
# Open the phantom dump
import shamrock
dump = shamrock.load_phantom_dump(filename)
# %%
# Print the data
dump.print_state()
Workflow reportworkflow report corresponding to commit c957664 Light CI is enabled. This will only run the basic tests and not the full tests. Pre-commit check reportPre-commit check: ✅ Test pipeline can run. Doxygen diff with
|
3 participants
Adds support for MHD fields (B/rho and psi/ch) in dump files.
Notes:
set_artif_viscosity_None(), dump files should not containalpha.