usable MHD solver #942
usable MHD solver #942
Conversation
updated + tested mhd_state, updated mhd_model revert back my typo on hybrid_tagger fixed typo in test_mhd_state preparing mhd usablevecfield new field initialisation logic Added end of file newline in field_user_initializer.hpp added mhd state fixture wrote test for mhd_state_fixtures added mhd yee grid, modified gridlayout to handle both hybrid and mhd up init_functions MHD solver solver mhd update more solver updates preparing rebase formatting riemann reconstruction quick commit fix solved conflict more solver updates initial branch commit: new architecture for reconstruction + flux computation improved test, fixed some bugs in godunov-fluxes better timing directory creation (PHAREHUB#914) safer nu (PHAREHUB#911) try fallback for dl on keyerror + ruffage (PHAREHUB#916) convenient utils for mpi/hierarchies keep pyattrs on compute_from_hier nan min/max to handle possible nan ghosts (PHAREHUB#923) * nan min/max to handle possible nan ghosts * flatten rm atefact file fixed missing template keyword for macos-12 build more explicit unwrapping in godunov fluxes with variadic arguments of unclear size fixed lambda capture problem fixed lambda capture problem added time integration added projection functions for centering, primitive/conservetive converter added projection functions for centering, primitive/conservetive converter solver ref binding in a lambda issue added template deduction guide on for macOS compile on a struct of MHDSolver New approach for template deduction for the struct in MHDSolver added usable ghost cells syntaxe fixes for Werror full boundary conditions in solver changed ghost cells function names to be more consistant with hybrid corrected dummyhierarchy in testmhdsolver (operation on nullptr) MHDMock simulator setup added pybind wrappers and first tests successful orszag-tang test fixed typo fixed typo fixed typo bug fixes (wrong usage of prevIndex nextIndex)
corrected van leer slope limiter second order corrected van leer slope limiter
hll solver
📝 WalkthroughWalkthroughThis pull request introduces a comprehensive enhancement of the PHARE framework's MHD (Magnetohydrodynamics) simulation capabilities. The changes span multiple components, including core data structures, numerical methods, solvers, and Python bindings. The modifications introduce new classes, methods, and utilities for handling MHD quantities, grid layouts, and simulation processes across different dimensions. The implementation focuses on providing a flexible and extensible framework for magnetohydrodynamic simulations. Changes
Sequence DiagramsequenceDiagram
participant Simulation
participant MHDModel
participant Simulator
participant Solver
Simulation->>MHDModel: Initialize with parameters
MHDModel-->>Simulation: Model configured
Simulation->>Simulator: Create simulator instance
Simulator->>Solver: Set up solver
Simulator->>Simulator: Configure simulation parameters
Simulator->>Solver: Advance simulation
Solver-->>Simulator: Update simulation state
Simulator->>Simulator: Write output
This sequence diagram illustrates the high-level interaction between key components during an MHD simulation, showing how the simulation is initialized, configured, and advanced through the solver. Finishing Touches
Thank you for using CodeRabbit. We offer it for free to the OSS community and would appreciate your support in helping us grow. If you find it useful, would you consider giving us a shout-out on your favorite social media? 🪧 TipsChatThere are 3 ways to chat with CodeRabbit:
Note: Be mindful of the bot's finite context window. It's strongly recommended to break down tasks such as reading entire modules into smaller chunks. For a focused discussion, use review comments to chat about specific files and their changes, instead of using the PR comments. CodeRabbit Commands (Invoked using PR comments)
Other keywords and placeholders
Documentation and Community
|
![github-advanced-security[bot]](https://avatars.githubusercontent.com/in/57789?s=60&v=4){kind=small}
| from . import global_vars | ||
|
|
||
| class MHDModel(object): | ||
| def defaulter(self, input, value): |
Check notice
Code scanning / CodeQL
Explicit returns mixed with implicit (fall through) returns Note
![coderabbitai[bot]](https://avatars.githubusercontent.com/in/347564?s=60&v=4){kind=small}
There was a problem hiding this comment.
Actionable comments posted: 23
🔭 Outside diff range comments (1)
tests/initializer/init_functions.hpp (1)
Confirmed: Files are identical and should be consolidated
The review comment is accurate. Both files contain identical code:
tests/core/data/mhd_state/init_functions.hpptests/initializer/init_functions.hppRecommendation:
- Move the shared code to a common test utilities location (e.g.,
tests/common/init_functions.hpp)- Update both test suites to include from the new location
- Delete the duplicate files
🔗 Analysis chain
Line range hint
1-83: Eliminate code duplication with tests/core/data/mhd_state/init_functions.hpp.This file appears to be identical to
tests/core/data/mhd_state/init_functions.hpp. Consider:
- Moving the shared code to a common location
- Having both test files include the common code
🏁 Scripts executed
The following scripts were executed for the analysis:
Script:
#!/bin/bash # Verify the duplication diff <(cat tests/initializer/init_functions.hpp) <(cat tests/core/data/mhd_state/init_functions.hpp)Length of output: 212
Script:
#!/bin/bash # First, find the files echo "Searching for files:" fd "init_functions.hpp" --type f # If files exist, show their contents for comparison echo -e "\nFile contents:" for f in $(fd "init_functions.hpp" --type f); do echo -e "\n=== $f ===" cat "$f" doneLength of output: 6930
🧹 Nitpick comments (54)
pyphare/pyphare/mock_mhd_simulator/mhd_model.py (2)
15-19: Ensure consistent handling of spatial variables in default functionsIn the
defaultermethod, the default function fordim == 1includes the spatial variablex, while for higher dimensions, the spatial variables are not utilized. This inconsistency might lead to confusion.Consider making the default functions consistent:
if self.dim == 1: - return lambda x: value + x * 0 + return lambda x: valueAlternatively, if the intention is to use spatial variables, adjust the functions to include them appropriately.
3-44: Add docstrings toMHDModelclass and its methodsProviding docstrings enhances code readability and helps others understand the purpose and usage of your class and methods.
Consider adding descriptive docstrings following PEP 257 conventions for the class and each method.
🧰 Tools
🪛 GitHub Check: CodeQL
[notice] 4-4: Explicit returns mixed with implicit (fall through) returns
Mixing implicit and explicit returns may indicate an error as implicit returns always return None.src/amr/physical_models/mhd_model.hpp (2)
93-93: Avoid unnecessary variable assignment when discarding return valueIn line 93, instead of assigning the result to
auto _, you can call the function without assigning its return value or explicitly cast it tovoidto indicate that the return value is intentionally discarded.Apply this diff to simplify the code:
-auto _ = this->resourcesManager->setOnPatch(*patch, state); +this->resourcesManager->setOnPatch(*patch, state);
59-64: Remove unnecessarystd::moveon const referenceIn the constructor,
_resourcesManageris aconstreference. Applyingstd::moveto aconstreference does not have any effect and may cause confusion. It's unnecessary here.Apply this diff to remove the unnecessary
std::move:explicit MHDModel(PHARE::initializer::PHAREDict const& dict, std::shared_ptr<resources_manager_type> const& _resourcesManager) : IPhysicalModel<AMR_Types>{model_name} , state{dict} - , resourcesManager{std::move(_resourcesManager)} + , resourcesManager{_resourcesManager} { }src/core/numerics/finite_volume_euler/finite_volume_euler.hpp (2)
41-45: Avoid unnecessary copying ofGridLayoutIn the constructor of
FiniteVolumeEuler_ref,GridLayoutis copied. IfGridLayoutis a large object, consider storing aconstreference to avoid unnecessary copying.Modify the code as follows:
class FiniteVolumeEuler_ref { public: - FiniteVolumeEuler_ref(GridLayout const& layout, double const dt) + FiniteVolumeEuler_ref(GridLayout const& layout, double const dt) : layout_{layout} , dt_{dt} { } private: - GridLayout layout_; + GridLayout const& layout_; double dt_; // ... };Ensure that the lifetime of the
layout_reference exceeds that of theFiniteVolumeEuler_refinstance to prevent dangling references.
59-111: Reduce code duplication by generalizing finite volume computationThe
finite_volume_euler_method contains duplicated code for different dimensions. Consider refactoring to reduce duplication and improve maintainability.You can generalize the finite volume computation across dimensions by iterating over the directions and fluxes dynamically. Here's an example:
void finite_volume_euler_(Field const& U, Field& Unew, MeshIndex<Field::dimension> index, const Fluxes&... fluxes) const { constexpr std::array<Direction, dimension> directions = {Direction::X, Direction::Y, Direction::Z}; const std::array<double, dimension> inverseMeshSizes = { layout_.inverseMeshSize(Direction::X), layout_.inverseMeshSize(Direction::Y), layout_.inverseMeshSize(Direction::Z) }; auto fluxes_array = {fluxes...}; auto fluxCenterings = {layout_.centering(flux.physicalQuantity())...}; Unew(index) = U(index); size_t dir = 0; for (const auto& flux : fluxes_array) { auto centering = fluxCenterings[dir]; auto nextIdx = index; nextIdx[dir] = layout_.nextIndex(centering[dir], index[dir]); Unew(index) -= dt_ * inverseMeshSizes[dir] * (flux(nextIdx) - flux(index)); ++dir; } }This approach eliminates the separate cases for each dimension and makes the code more flexible.
src/core/numerics/constrained_transport/constrained_transport.hpp (2)
25-27: Use a more specific exception type for layout validationInstead of throwing a general
std::runtime_errorwhen the layout is not set, consider using a more specific exception likestd::logic_errorto clearly indicate that this is a logical error related to the program's state.Apply this diff to change the exception type:
if (!this->hasLayout()) - throw std::runtime_error( + throw std::logic_error( "Error - ConstrainedTransport - GridLayout not set, cannot proceed to computation");
59-62: Explicitly capture variables in lambda functions for clarityIn the lambda functions passed to
layout_.evalOnBox, you are capturing all variables by reference using[&]. To improve code clarity and prevent unintended captures, consider explicitly listing the variables you need to capture.Modify the lambdas to capture only the necessary variables:
- layout_.evalOnBox(Ey, [&](auto&... args) mutable { Ey_(Ey, {args...}, Bz_x); }); + layout_.evalOnBox(Ey, [this, &Ey, &Bz_x](auto&... args) mutable { Ey_(Ey, {args...}, Bz_x); }); - layout_.evalOnBox(Ez, [&](auto&... args) mutable { Ez_(Ez, {args...}, By_x); }); + layout_.evalOnBox(Ez, [this, &Ez, &By_x](auto&... args) mutable { Ez_(Ez, {args...}, By_x); });Repeat similar changes for other lambda functions in the code.
src/core/data/grid/gridlayoutdefs.hpp (3)
91-161: Refactor to eliminate code duplication betweengridDataTandgridDataT_mhdBoth structs share several identical definitions, such as directions and centerings. To improve maintainability and reduce redundancy, consider creating a base struct that encapsulates the common elements, and have both
gridDataTandgridDataT_mhdinherit from it.Example:
struct gridDataBase { static constexpr Direction dirX = Direction::X; static constexpr Direction dirY = Direction::Y; static constexpr Direction dirZ = Direction::Z; static constexpr QtyCentering primal = QtyCentering::primal; static constexpr QtyCentering dual = QtyCentering::dual; static constexpr std::uint32_t idirX = static_cast<std::uint32_t>(Direction::X); static constexpr std::uint32_t idirY = static_cast<std::uint32_t>(Direction::Y); static constexpr std::uint32_t idirZ = static_cast<std::uint32_t>(Direction::Z); }; struct gridDataT : public gridDataBase { // Existing hybrid quantity constants... }; struct gridDataT_mhd : public gridDataBase { // MHD-specific quantity constants... };
97-98: Avoid redefining already defined constantsThe constants
primalanddualare defined in bothgridDataTandgridDataT_mhd. Since they represent universal concepts, consider defining them once in a shared scope or base struct to prevent redundancy.
104-161: Consistency in constant naming conventionsThe naming of constants varies between abbreviated (e.g.,
iVx,iBx) and more descriptive forms (e.g.,iScalarFlux_x). For better readability and maintainability, ensure that constant names follow a consistent naming convention throughout the codebase.Consider standardizing the prefixing and casing, such as using
iVx,iVy,iVzconsistently for vector components, and similarly for other quantities.src/amr/solvers/solver_mhd_model_view.hpp (1)
70-72: Consider exception handling around core function callsWhile calling core functions like
constrained_transport_, any exceptions thrown may propagate upwards. To enhance robustness, consider wrapping these calls in try-catch blocks to handle exceptions gracefully and provide meaningful error messages.Example:
try { constrained_transport_(*E[i], *fluxes[i]...); } catch (const std::exception& e) { // Handle exception, possibly logging the error }src/core/numerics/time_integrator/time_integrator.hpp (2)
29-69: Consider refactoring to reduce code duplication betweentvdrk2andtvdrk3The methods
tvdrk2andtvdrk3share similar code structures and repetitive blocks. Refactoring the common code into helper functions or using templates can improve maintainability and reduce redundancy.
139-139: Unused member variableintegrator_The member variable
std::string const integrator_is declared but not used within theTimeIntegratorclass. Consider removing it if it's unnecessary or utilizing it appropriately.src/amr/solvers/solver_mhd.hpp (1)
102-104: Empty implementation offillMessengerInfoThe
fillMessengerInfomethod is empty. If no action is required, consider adding a comment explaining why. Otherwise, implement the necessary functionality.tests/core/numerics/mock_mhd_simulator/test_mhd_solver.hpp (1)
635-635: Avoid usingstd::coutfor logging in simulationsUsing
std::coutfor logging can clutter console output and is not suitable for simulations. Consider using a logging framework or control output verbosity through configuration options.src/core/numerics/godunov_fluxes/godunov_fluxes.hpp (4)
42-113: Refactor to eliminate code duplication across dimensions.The code blocks handling different dimensions (1D, 2D, 3D) have similar structures with repetitive code. Consider refactoring these sections to abstract the common functionality into shared methods or templates. This will enhance maintainability and reduce potential errors due to duplicated logic.
173-190: Remove or document commented-out code.There is a block of commented-out code related to resistive contributions and Laplacians. If this code is no longer needed, consider removing it to keep the codebase clean. If it's intended for future use, add comments explaining its purpose or consider using feature flags.
676-677: Enhance error messages for unknown reconstruction methods.When an unknown reconstruction method is encountered, a generic runtime error is thrown. Include the invalid method name in the error message to assist in debugging and quickly identify the issue.
121-124: Use enumerations instead of strings for method identifiers.The member variables
terms_,reconstruction_,limiter_, andriemann_are stored as strings and compared using string literals. Consider using enumerations or string constants to represent these identifiers. This approach improves type safety, reduces the risk of typos, and enhances performance by avoiding costly string comparisons.src/core/data/grid/gridlayoutimplyee_mhd.hpp (1)
142-190: Refactorcenteringmethods to reduce code duplication.The
centeringmethods for scalar quantities contain extensive switch-case statements that are largely duplicated across different dimensions (dim == 1,dim == 2,dim == 3). This duplication increases maintenance effort and the chance of inconsistencies. Consider refactoring these methods using templates, data-driven approaches, or mapping tables to eliminate redundancy and improve code clarity.Also applies to: 194-276, 280-405
src/core/data/grid/gridlayout.hpp (1)
837-838: Address the TODO comment regarding_QtyCenteringdefinition.The code contains a TODO comment suggesting that
_QtyCenteringshould be defined per dimension to simplify computations. Implementing this suggestion can enhance code maintainability and readability by allowing direct usage of_QtyCentering[iQty]without additional loops.Do you want me to help refactor the
_QtyCenteringdefinition to be dimension-specific for improved clarity?pyphare/pyphare/mock_mhd_simulator/simulation.py (2)
10-11: Add attribute validationDirect use of
object.__setattr__without validation could lead to setting invalid or unexpected attributes.Consider adding attribute validation:
VALID_ATTRIBUTES = {'model', 'dimension', 'time_step'} # add your valid attributes def __init__(self, **kwargs): for k, v in kwargs.items(): if k not in VALID_ATTRIBUTES: raise ValueError(f"Invalid attribute: {k}") object.__setattr__(self, k, v)
15-16: Add type checking for modelThe
set_modelmethod should validate the model type to prevent setting invalid objects.from .mhd_model import MHDModel # adjust import path as needed def set_model(self, model): if not isinstance(model, MHDModel): raise TypeError("Model must be an instance of MHDModel") self.model = modeltests/core/data/field/test_usable_field_fixtures_mhd.hpp (1)
12-20: Add class documentation to explain the purpose and usage.Consider adding Doxygen-style documentation to explain:
- Purpose of the UsableFieldMHD class
- Template parameter requirements
- Relationship with the base Field class
- Usage examples
tests/core/numerics/mock_mhd_simulator/alfven1dY.py (2)
6-24: Add docstring to explain simulation configuration and parameters.The
configfunction should include documentation explaining:
- Physical setup of the 1D Alfven wave
- Choice of parameters (timestep, cells, etc.)
- Units of measurement
25-51: Document physical constants and their significance.Add comments explaining:
- The choice of wave number
kx- The amplitude values (1e-6)
- The background field values
tests/core/numerics/mock_mhd_simulator/alfven1d.py (1)
55-57: Make output filename consistent with simulation configuration.The output filename "alfven1d.h5" should reflect the wave direction to distinguish it from alfven1dY.py output.
- MHDMockSimulator(config()).run("alfven1d.h5") + MHDMockSimulator(config()).run("alfven1d_x.h5")pyphare/pyphare/mock_mhd_simulator/simulator.py (1)
15-48: Add type hints and improve documentation.Consider adding:
- Type hints for parameters and return values
- Class docstring explaining the simulator's purpose
- Method documentation
Example:
from typing import Optional from .simulation import Simulation class MHDMockSimulator: """Mock simulator for MHD equations testing. Manages the lifecycle of C++ simulator instance and provides a Python interface for simulation execution. """ def __init__(self, simulation: Simulation) -> None: """Initialize simulator with given configuration."""tests/core/numerics/mock_mhd_simulator/reversed_orszag_tang.py (2)
25-25: Document the physical significance of reversed fieldsThe velocity and magnetic field components are intentionally reversed from the standard Orszag-Tang setup. Please add a comment explaining the physical motivation and expected differences in the solution.
Also applies to: 30-44
55-57: Add error handling for file operationsThe simulation output might fail if the directory is not writable or if there's insufficient disk space.
def main(): - MHDMockSimulator(config()).run("reversed_orszag_tang.h5", dumpfrequency = 1) + try: + MHDMockSimulator(config()).run("reversed_orszag_tang.h5", dumpfrequency = 1) + except IOError as e: + print(f"Failed to write output file: {e}") + raisetests/core/numerics/mock_mhd_simulator/read.py (2)
1-4: Remove unused importThe
numpymodule is imported but never used in the code.import h5py import matplotlib.pyplot as plt -import numpy as np🧰 Tools
🪛 Ruff (0.8.2)
3-3:
numpyimported but unusedRemove unused import:
numpy(F401)
32-38: Enhance plot customization optionsThe plotting parameters (figsize, colormap, etc.) should be configurable through function parameters.
-def read_and_plot_h5(file_path, timestep, quantity): +def read_and_plot_h5(file_path, timestep, quantity, figsize=(8, 6), cmap='coolwarm'): # ... - plt.figure(figsize=(8, 6)) - plt.imshow(data.T, origin="lower", cmap="coolwarm", aspect="auto") + plt.figure(figsize=figsize) + plt.imshow(data.T, origin="lower", cmap=cmap, aspect="auto")tests/core/numerics/mock_mhd_simulator/orszag_tang.py (1)
58-58: Consider reducing output frequencyA dump frequency of 80 with timestep 0.001 will generate output every 0.08 time units. This might create large output files.
Consider reducing the frequency if storage is a concern:
- MHDMockSimulator(config()).run("orszag_tang.h5", dumpfrequency=80) + MHDMockSimulator(config()).run("orszag_tang.h5", dumpfrequency=200)tests/core/numerics/mock_mhd_simulator/alfven2d.py (2)
30-35: Document wave parametersThe wave angle and wave numbers calculation should be documented for clarity.
+ # Set wave propagation angle to arctan(0.5) ≈ 26.57° alpha = np.arctan(0.5) cosalpha = np.cos(alpha) sinalpha = np.sin(alpha) + # Calculate wave numbers ensuring one complete wavelength in the domain kx = 2 * np.pi / (cosalpha * lx) ky = 2 * np.pi / (sinalpha * ly)
39-53: Consider extracting wave amplitude parameterThe perturbation amplitude (1e-4) is used in multiple places. Consider defining it as a constant.
+ # Small amplitude for linear regime + PERTURBATION_AMPLITUDE = 1e-4 + def vx(x, y): - return (-1e-4 * np.sin(kx * x * cosalpha + ky * y * sinalpha)) / (2.0 * sinalpha) + return (-PERTURBATION_AMPLITUDE * np.sin(kx * x * cosalpha + ky * y * sinalpha)) / (2.0 * sinalpha)tests/core/numerics/mock_mhd_simulator/whistler1d.py (2)
24-24: Consider making the timestep calculation dynamicThe hardcoded timestep of 0.0006 might not be optimal for all configurations. Consider calculating it based on CFL conditions for numerical stability.
37-38: Consider documenting the significance of mode numbersThe choice of modes [1, 2, 4, 8] seems deliberate. Adding a comment explaining the physical significance of these modes would improve code maintainability.
tests/core/numerics/mock_mhd_simulator/rotor.py (2)
25-26: Document the physical units and normalizationThe magnetic field strength
B0and velocityv0constants should include comments about their physical units and normalization scheme.
41-42: Consider vectorizing the conditional calculationsThe repeated use of
np.wherefor conditional calculations could be optimized by pre-computing the condition masks.+ mask_inner = r_ <= r0 + mask_transition = (r_ > r0) & (r_ < r1) + def density(x, y): r_ = r(x, y) f_ = f(r_) - rho_values = np.where(r_ <= r0, 10.0, np.where(r_ < r1, 1.0 + 9.0 * f_, 1.0)) + rho_values = np.ones_like(r_) + rho_values[mask_transition] = 1.0 + 9.0 * f_[mask_transition] + rho_values[mask_inner] = 10.0 return rho_valuesAlso applies to: 48-49, 55-56
tests/core/data/tensorfield/test_tensorfield_fixtures_mhd.hpp (2)
14-18: Enhance documentation with usage examplesWhile the current documentation explains the purpose, adding a brief usage example would make it more developer-friendly.
38-43: Consider adding validation in set_on methodThe
set_onmethod should validate that the input tensorfield has the correct number of components before setting buffers.void set_on(Super& tensorfield) { + if (tensorfield.size() != N_elements) + throw std::invalid_argument("Tensorfield has incorrect number of components"); // used for setting on normal model tensorfields for (std::size_t i = 0; i < N_elements; ++i) tensorfield[i].setBuffer(&xyz[i]); }tests/core/numerics/mock_mhd_simulator/harris.py (1)
32-33: Rename parameter for better clarityThe parameter
lis ambiguous. Consider renaming it tothicknessorscale_lengthfor better readability.-def S(y, y0, l): - return 0.5 * (1.0 + np.tanh((y - y0) / l)) +def S(y, y0, thickness): + return 0.5 * (1.0 + np.tanh((y - y0) / thickness))🧰 Tools
🪛 Ruff (0.8.2)
32-32: Ambiguous variable name:
l(E741)
tests/core/data/mhd_state/test_mhd_state_fixtures.hpp (2)
54-66: Consider adding noexcept to move constructor.The move constructor should be marked
noexceptas it doesn't throw exceptions and this helps optimize move operations.- UsableMHDState(UsableMHDState&& that) + UsableMHDState(UsableMHDState&& that) noexcept
73-82: Consider grouping related member variables.The member variables could be organized better by grouping related quantities together with comments explaining their physical meaning.
+ // Primary MHD quantities Grid_t rho; UsableVecFieldMHD<dim> V; UsableVecFieldMHD<dim> B; Grid_t P; + // Derived quantities UsableVecFieldMHD<dim> rhoV; Grid_t Etot; + // Electromagnetic quantities UsableVecFieldMHD<dim> J, E;tests/core/data/mhd_state/test_mhd_state.cpp (2)
20-38: Consider parameterizing the test data.The
getDict()function hardcodes test data. Consider using a test data provider pattern for better test maintainability.+struct MHDTestData { + using initfunc = InitFunction<1>; + static PHAREDict createDict(initfunc density_fn, initfunc vx_fn, /*...*/) { + PHAREDict dict; + dict["density"]["initializer"] = static_cast<initfunc>(density_fn); + // ... other initializations + return dict; + } +};
40-46: Empty destructor can be defaulted.The empty destructor implementation can be simplified using
= default.- virtual ~AnMHDState(); + virtual ~AnMHDState() = default;And remove the separate implementation at line 46.
src/core/numerics/primite_conservative_converter/to_conservative_converter.hpp (1)
24-42: Consider vectorization opportunities in vToRhoV.The component-wise multiplication in
vToRhoVcould benefit from SIMD optimization.Consider using compiler vectorization hints or explicit SIMD instructions for better performance.
src/core/mhd/mhd_quantities.hpp (2)
13-58: Consider enhancing documentation with physical units and constraints.The enums provide good type safety, but would benefit from additional documentation:
- Physical units for each quantity
- Valid ranges or constraints
- Relationships between quantities (e.g., how E, B, and J are related)
70-99: Consider using switch statement for better performance.The if-else chain could be replaced with a switch statement for potentially better performance. Also, consider making the error message more specific by including the invalid vector value.
- if (qty == Vector::V) - return {{Scalar::Vx, Scalar::Vy, Scalar::Vz}}; - if (qty == Vector::B) - return {{Scalar::Bx, Scalar::By, Scalar::Bz}}; + switch(qty) { + case Vector::V: + return {{Scalar::Vx, Scalar::Vy, Scalar::Vz}}; + case Vector::B: + return {{Scalar::Bx, Scalar::By, Scalar::Bz}}; + // ... other cases ... + default: + throw std::runtime_error("Invalid Vector value: " + std::to_string(static_cast<int>(qty))); + }pyphare/pyphare/mock_mhd_simulator/__init__.py (3)
14-26: Add input validation for better error handling.Consider adding validation for:
- Empty arrays
- NaN/Inf values
- Array shape consistency
def __call__(self, *xyz): if not xyz: raise ValueError("No input arguments provided") args = [np.asarray(arg) for arg in xyz] if any(np.any(~np.isfinite(arg)) for arg in args): raise ValueError("Input contains NaN or Inf values") if any(arg.size == 0 for arg in args): raise ValueError("Empty array in input") ret = self.fn(*args) if isinstance(ret, list): ret = np.asarray(ret) if is_scalar(ret): ret = np.full(len(args[-1]), ret) return ret
30-40: Consider moving import to class level.The cpp_etc_lib import could be moved to class level to avoid repeated imports on each call.
class fn_wrapper(py_fn_wrapper): + from pyphare.cpp import cpp_etc_lib + def __init__(self, fn): super().__init__(fn) def __call__(self, *xyz): - from pyphare.cpp import cpp_etc_lib - # convert numpy array to C++ SubSpan # couples vector init functions to C++ - return cpp_etc_lib().makePyArrayWrapper(super().__call__(*xyz)) + return self.cpp_etc_lib().makePyArrayWrapper(super().__call__(*xyz))
65-102: Consider using constants for dictionary paths.The dictionary paths could be defined as constants to improve maintainability and reduce typos.
# At module level PATHS = { 'TIME_STEP': 'time_step', 'FINAL_TIME': 'final_time', 'DENSITY_INIT': 'density/initializer', # ... other paths } # In populateDict add_double(PATHS['TIME_STEP'], simulation.timestep) add_double(PATHS['FINAL_TIME'], simulation.final_time) addInitFunction(PATHS['DENSITY_INIT'], fn_wrapper(d["density"]))src/python3/CMakeLists.txt (1)
16-16: Consider adding optimization and link flags properties.For consistency with other modules (cpp, cpp_etc), consider adding the following properties:
set_target_properties(pyMHD PROPERTIES LIBRARY_OUTPUT_DIRECTORY "${CMAKE_BINARY_DIR}/pybindlibs" ) + set_property(TARGET pyMHD PROPERTY INTERPROCEDURAL_OPTIMIZATION ${PHARE_INTERPROCEDURAL_OPTIMIZATION}) + set_property(TARGET pyMHD APPEND_STRING PROPERTY LINK_FLAGS " ${PHARE_LINK_FLAGS}")
📜 Review details
Configuration used: .coderabbit.yaml
Review profile: CHILL
Plan: Pro
📒 Files selected for processing (52)
.gitignore(1 hunks)pyphare/pyphare/mock_mhd_simulator/__init__.py(1 hunks)pyphare/pyphare/mock_mhd_simulator/global_vars.py(1 hunks)pyphare/pyphare/mock_mhd_simulator/mhd_model.py(1 hunks)pyphare/pyphare/mock_mhd_simulator/simulation.py(1 hunks)pyphare/pyphare/mock_mhd_simulator/simulator.py(1 hunks)res/cmake/test.cmake(1 hunks)src/amr/messengers/mhd_messenger_info.hpp(2 hunks)src/amr/physical_models/mhd_model.hpp(1 hunks)src/amr/solvers/solver_mhd.hpp(1 hunks)src/amr/solvers/solver_mhd_model_view.hpp(1 hunks)src/core/data/field/field.hpp(0 hunks)src/core/data/field/initializers/field_user_initializer.hpp(2 hunks)src/core/data/grid/gridlayout.hpp(12 hunks)src/core/data/grid/gridlayout_impl.hpp(1 hunks)src/core/data/grid/gridlayoutdefs.hpp(2 hunks)src/core/data/grid/gridlayoutimplyee.hpp(5 hunks)src/core/data/grid/gridlayoutimplyee_mhd.hpp(1 hunks)src/core/data/vecfield/vecfield_initializer.hpp(1 hunks)src/core/mhd/mhd_quantities.hpp(1 hunks)src/core/models/mhd_state.hpp(1 hunks)src/core/numerics/constrained_transport/constrained_transport.hpp(1 hunks)src/core/numerics/faraday/faraday.hpp(6 hunks)src/core/numerics/finite_volume_euler/finite_volume_euler.hpp(1 hunks)src/core/numerics/godunov_fluxes/godunov_fluxes.hpp(1 hunks)src/core/numerics/primite_conservative_converter/to_conservative_converter.hpp(1 hunks)src/core/numerics/primite_conservative_converter/to_primitive_converter.hpp(1 hunks)src/core/numerics/time_integrator/time_integrator.hpp(1 hunks)src/python3/CMakeLists.txt(1 hunks)src/python3/pyMHD.cpp(1 hunks)tests/amr/messengers/test_messengers.cpp(1 hunks)tests/amr/multiphysics_integrator/test_multiphysics_integrator.cpp(1 hunks)tests/core/data/electromag/test_electromag_fixtures.hpp(1 hunks)tests/core/data/field/test_field_fixtures_mhd.hpp(1 hunks)tests/core/data/field/test_usable_field_fixtures_mhd.hpp(1 hunks)tests/core/data/mhd_state/CMakeLists.txt(1 hunks)tests/core/data/mhd_state/init_functions.hpp(1 hunks)tests/core/data/mhd_state/test_mhd_state.cpp(1 hunks)tests/core/data/mhd_state/test_mhd_state_fixtures.hpp(1 hunks)tests/core/data/tensorfield/test_tensorfield_fixtures_mhd.hpp(1 hunks)tests/core/data/vecfield/test_vecfield_fixtures_mhd.hpp(1 hunks)tests/core/numerics/mock_mhd_simulator/alfven1d.py(1 hunks)tests/core/numerics/mock_mhd_simulator/alfven1dY.py(1 hunks)tests/core/numerics/mock_mhd_simulator/alfven2d.py(1 hunks)tests/core/numerics/mock_mhd_simulator/harris.py(1 hunks)tests/core/numerics/mock_mhd_simulator/orszag_tang.py(1 hunks)tests/core/numerics/mock_mhd_simulator/read.py(1 hunks)tests/core/numerics/mock_mhd_simulator/reversed_orszag_tang.py(1 hunks)tests/core/numerics/mock_mhd_simulator/rotor.py(1 hunks)tests/core/numerics/mock_mhd_simulator/test_mhd_solver.hpp(1 hunks)tests/core/numerics/mock_mhd_simulator/whistler1d.py(1 hunks)tests/initializer/init_functions.hpp(1 hunks)
💤 Files with no reviewable changes (1)
- src/core/data/field/field.hpp
✅ Files skipped from review due to trivial changes (7)
- pyphare/pyphare/mock_mhd_simulator/global_vars.py
- src/core/data/grid/gridlayout_impl.hpp
- src/core/data/field/initializers/field_user_initializer.hpp
- src/core/data/vecfield/vecfield_initializer.hpp
- .gitignore
- tests/core/data/electromag/test_electromag_fixtures.hpp
- tests/amr/messengers/test_messengers.cpp
🧰 Additional context used
📓 Path-based instructions (25)
tests/core/data/field/test_field_fixtures_mhd.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/data/grid/gridlayoutdefs.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
tests/core/data/mhd_state/test_mhd_state_fixtures.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
tests/core/data/vecfield/test_vecfield_fixtures_mhd.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/numerics/constrained_transport/constrained_transport.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
tests/core/data/tensorfield/test_tensorfield_fixtures_mhd.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
tests/initializer/init_functions.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/numerics/finite_volume_euler/finite_volume_euler.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/amr/messengers/mhd_messenger_info.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/models/mhd_state.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/numerics/primite_conservative_converter/to_conservative_converter.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/numerics/godunov_fluxes/godunov_fluxes.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
tests/core/numerics/mock_mhd_simulator/test_mhd_solver.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/numerics/faraday/faraday.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/data/grid/gridlayoutimplyee.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/numerics/primite_conservative_converter/to_primitive_converter.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/numerics/time_integrator/time_integrator.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/data/grid/gridlayoutimplyee_mhd.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/amr/physical_models/mhd_model.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
tests/core/data/field/test_usable_field_fixtures_mhd.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/amr/solvers/solver_mhd.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/mhd/mhd_quantities.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
tests/core/data/mhd_state/init_functions.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/amr/solvers/solver_mhd_model_view.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
src/core/data/grid/gridlayout.hpp (1)
Pattern **/*.hpp: Review the C++ code, point out issues relative to principles of clean code, expressiveness, and performance.
🪛 Ruff (0.8.2)
tests/core/numerics/mock_mhd_simulator/harris.py
32-32: Ambiguous variable name: l
(E741)
pyphare/pyphare/mock_mhd_simulator/simulator.py
10-10: pybindlibs.pyMHD imported but unused
Remove unused import: pybindlibs.pyMHD
(F401)
tests/core/numerics/mock_mhd_simulator/read.py
3-3: numpy imported but unused
Remove unused import: numpy
(F401)
🪛 GitHub Check: CodeQL
pyphare/pyphare/mock_mhd_simulator/mhd_model.py
[notice] 4-4: Explicit returns mixed with implicit (fall through) returns
Mixing implicit and explicit returns may indicate an error as implicit returns always return None.
⏰ Context from checks skipped due to timeout of 90000ms (6)
- GitHub Check: build (macos-latest)
- GitHub Check: build (macos-14)
- GitHub Check: Analyze (cpp)
- GitHub Check: build (macos-13)
- GitHub Check: build (ubuntu-latest)
- GitHub Check: Analyze (python)
🔇 Additional comments (34)
src/python3/pyMHD.cpp (1)
1-69: LGTM!The
pyMHDmodule is well-implemented usingpybind11. The use of templates and modern C++ features effectively exposes C++ classes to Python, and the code adheres to best practices.src/core/numerics/primite_conservative_converter/to_primitive_converter.hpp (1)
1-77: Overall implementation aligns with best practicesThe
ToPrimitiveConverterclass is well-structured and follows clean coding principles.src/amr/physical_models/mhd_model.hpp (1)
88-98: Initialization logic is clear and well-structuredThe
initializemethod correctly iterates over the level patches, sets resources on each patch, and initializes the state with the corresponding layout.src/core/models/mhd_state.hpp (1)
28-39:isUsableandisSettablemethods are comprehensiveThe methods now include checks for all member variables, ensuring that the state is fully validated before use.
src/core/numerics/finite_volume_euler/finite_volume_euler.hpp (1)
17-32: Implementation ofFiniteVolumeEuleris correctThe
FiniteVolumeEulerclass effectively delegates computations toFiniteVolumeEuler_refafter validating the layout.src/core/numerics/constrained_transport/constrained_transport.hpp (2)
149-154: Inconsistent handling of dimensions inEy_methodIn the
Ey_method, fordimension <= 2, you assignEy(index) = Bz_x(index);, whereas fordimension == 3, a different computation is performed. Verify that this behavior is intentional and consistent with the physical equations you are modeling.Ensure that the simplified assignment for lower dimensions is appropriate and does not omit necessary calculations.
116-119:⚠️ Potential issuePotential out-of-bounds access due to index subtraction
In the
Ex_method for 3D, you subtract 1 fromindex[2]andindex[1]when accessing array elements. Ensure thatindex[2]andindex[1]are greater than zero to prevent accessing invalid memory locations, which can lead to undefined behavior or crashes.Consider adding checks or adjusting the loop boundaries to ensure that
index[2] - 1andindex[1] - 1are valid indices.src/amr/solvers/solver_mhd_model_view.hpp (2)
114-119: Typographical error in variable nameto_primtitve_The variable
to_primtitve_appears to be misspelled. The correct spelling should beto_primitive_. Correcting this enhances code clarity and prevents confusion during future maintenance.[typographical_and_grammatical_issue]
Apply this diff to correct the typo:
- core_type to_primtitve_; + core_type to_primitive_;Ensure that all instances of
to_primtitve_in the code are updated accordingly.
135-137: Verify state consistency after time integrationAfter performing time integration steps in the
eulermethod, ensure that the updated states (rho,rhoV,B,Etot,E) are consistent and valid. Any unintended side effects or errors in the integration process could affect simulation accuracy.Consider adding validation checks or tests to confirm that the states remain within expected physical bounds.
src/amr/solvers/solver_mhd.hpp (1)
106-110:allocatemethod contains TODO but lacks implementationThe
allocatemethod includes a TODO comment but no implementation. Ensure this method is implemented before merging to avoid potential issues during allocation.tests/core/data/field/test_field_fixtures_mhd.hpp (1)
10-11: LGTM! Clean type alias implementationThe template-based type alias for
FieldMHDis well-designed and follows modern C++ practices.src/amr/messengers/mhd_messenger_info.hpp (2)
18-18: Address TODO commentThe comment "What is needed here ?" suggests incomplete implementation or missing requirements.
Please clarify the requirements for this class. Consider:
- What additional methods or data members are needed?
- Are there any pure virtual methods from
IMessengerInfothat need implementation?
16-16: LGTM! Clean type aliasThe type alias for
VecFieldNamesis well-implemented.tests/core/data/vecfield/test_vecfield_fixtures_mhd.hpp (1)
12-17: LGTM! Well-designed type aliasesThe template-based type aliases for
VecFieldMHDandUsableVecFieldMHDare well-implemented with clear documentation for the rank parameter.tests/core/data/field/test_usable_field_fixtures_mhd.hpp (1)
21-28: Add parameter validation in the constructor.The constructor should validate that the layout dimensions match the template parameter
dimto prevent runtime errors.Consider adding assertions:
template<typename GridLayout> UsableFieldMHD(std::string const& name, GridLayout const& layout, scalar_t qty) : Super{name, qty} , xyz{name, layout, qty} { + static_assert(GridLayout::dimension == dim, + "Layout dimension must match template parameter"); super().setBuffer(&xyz); }✅ Verification successful
The static_assert suggestion is appropriate and follows established patterns.
The suggested validation aligns with similar dimension checks throughout the codebase, particularly in
vecfield_initializer.hppandgridlayout.hpp. This addition will maintain consistency and prevent potential runtime errors through compile-time validation.🏁 Scripts executed
The following scripts were executed for the analysis:
Error: No shell scripts found
Script:
#!/bin/bash # Find the complete class definition ast-grep --pattern 'template$$$class UsableFieldMHD$$$' # Find GridLayout related files fd -e hpp -e h -e cpp GridLayout # Look for similar dimension validation patterns rg "static_assert.*dimension" -A 2Length of output: 3660
Script:
#!/bin/bash # Find UsableFieldMHD class with more context rg -B 5 -A 10 "class UsableFieldMHD" # Alternative search using ast-grep with a different pattern ast-grep --pattern 'template<$$$> class UsableFieldMHD$$$'Length of output: 1555
pyphare/pyphare/mock_mhd_simulator/simulator.py (1)
45-48: Add validation for simulation parameters.The
runmethod should validate:
- Positive dumpfrequency
- Valid filename format
def run(self, filename: str, dumpfrequency: int = 1) -> 'MHDMockSimulator': + if dumpfrequency < 1: + raise ValueError("dumpfrequency must be positive") + if not filename.endswith('.h5'): + raise ValueError("filename must have .h5 extension") self._check_init() self.cpp_sim.advance(filename, dumpfrequency) return self✅ Verification successful
Parameter validation is necessary and consistent with codebase practices
The suggested validations are appropriate:
- The .h5 extension requirement aligns with the codebase's data serialization format
- Positive dumpfrequency validation prevents invalid simulation configurations
🏁 Scripts executed
The following scripts were executed for the analysis:
Error: No shell scripts found
Script:
#!/bin/bash # Search for .h5 file usage and validation patterns rg -l "\.h5" # Look for existing parameter validations in simulator files rg "frequency|dump.*frequency" pyphare/pyphare/mock_mhd_simulator/ # Search for validation patterns ast-grep --pattern 'if not $_.endswith(".h5")' ast-grep --pattern 'if $_ < 1:' # Check cpp_sim advance implementation ast-grep --pattern 'def advance($$$)'Length of output: 1541
Script:
#!/bin/bash # Check error handling in pharesee modules rg "raise|error|exception" pyphare/pyphare/pharesee/ # Look for file validation in run module cat pyphare/pyphare/pharesee/run/run.py # Check for any validation in fromh5 cat pyphare/pyphare/pharesee/hierarchy/fromh5.pyLength of output: 26561
tests/core/numerics/mock_mhd_simulator/reversed_orszag_tang.py (1)
11-12: Consider reducing the timestep for better stabilityThe current timestep of 0.014 might be too large for the MHD solver. For MHD simulations, the timestep should satisfy both the CFL condition and the magnetic diffusion condition.
Consider reducing it to 0.001 initially and validate stability.
tests/core/numerics/mock_mhd_simulator/orszag_tang.py (1)
19-23: LGTM! Well-configured solver optionsThe choice of reconstruction ("linear"), limiter ("minmod"), and Riemann solver ("hll") with TVDRK2 integration is appropriate for this test case.
tests/core/numerics/mock_mhd_simulator/whistler1d.py (1)
46-50: Verify numerical stability of velocity calculationsThe velocity components use summation of trigonometric functions with increasing mode numbers. Ensure this doesn't lead to numerical instabilities at high mode numbers.
✅ Verification successful
Implementation is numerically stable
The velocity calculations are safe from numerical instabilities due to:
- Limited mode numbers [1,2,4,8]
- Small amplitude scaling (1e-2)
- Controlled wavenumber based on domain length
🏁 Scripts executed
The following scripts were executed for the analysis:
Script:
#!/bin/bash # Search for similar implementations in the codebase to verify the approach ast-grep --pattern 'sum(np.$_(k * x * m + phi) * $_ for m, phi in zip($_, $_))'Length of output: 221
Script:
#!/bin/bash # Find definition of modes and k variables rg -A 5 "modes\s*=|k\s*=" tests/core/numerics/mock_mhd_simulator/whistler1d.py # Look for test cases using these velocity functions rg -B 5 "vy\(|vz\(" tests/core/numerics/mock_mhd_simulator/whistler1d.pyLength of output: 553
res/cmake/test.cmake (1)
20-20: LGTM!The addition of the MHD state tests directory follows the established pattern and is correctly placed within the core data test group.
tests/core/data/mhd_state/test_mhd_state_fixtures.hpp (2)
15-17: LGTM! Well-structured class definition.The template parameter for dimension and inheritance from
MHDStateprovide good flexibility while maintaining type safety.
20-31: Consider adding error handling in_set().The
_set()method performs critical buffer assignments but lacks error handling for potential null buffers or invalid states.Consider adding validation:
void _set() { auto&& [_rho, _V, _B, _P, _rhoV, _Etot, _J, _E] = Super::getCompileTimeResourcesViewList(); + if (!_rho.isValid() || !_V.isValid() || !_B.isValid() || !_P.isValid() || + !_rhoV.isValid() || !_Etot.isValid() || !_J.isValid() || !_E.isValid()) { + throw std::runtime_error("Invalid resource view in MHD state"); + } _rho.setBuffer(&rho); V.set_on(_V); B.set_on(_B); _P.setBuffer(&P); rhoV.set_on(_rhoV); _Etot.setBuffer(&Etot); J.set_on(_J); E.set_on(_E); }tests/core/data/mhd_state/test_mhd_state.cpp (1)
48-78: Add more test coverage for edge cases.The current tests cover basic functionality but miss important edge cases:
- Invalid initialization values
- Boundary conditions
- State transitions
Consider adding these test cases to improve coverage.
src/core/mhd/mhd_quantities.hpp (1)
59-69: LGTM! Well-designed component accessors.Good practices observed:
- Use of NO_DISCARD prevents accidental value drops
- constexpr enables compile-time evaluation
- Consistent naming convention
pyphare/pyphare/mock_mhd_simulator/__init__.py (1)
4-10: LGTM! Clear and focused utility functions.The type checking functions are well-designed with clear responsibilities.
tests/core/data/mhd_state/init_functions.hpp (1)
1-83: Duplicate file - see previous review comments.src/core/numerics/faraday/faraday.hpp (3)
Line range hint
13-35: Well-structured implementation of Faraday_ref class!The class follows good C++ practices with proper const correctness, clean error handling, and a modular design that separates concerns between Faraday and Faraday_ref.
62-64: Efficient implementation using evalOnGhostBox!The consistent use of lambda functions with evalOnGhostBox provides a clean and efficient way to handle component calculations.
Line range hint
99-129: Mathematically sound implementation of field equations!The dimension-specific implementations correctly handle the mathematical formulations for 1D, 2D, and 3D cases using compile-time conditionals.
tests/amr/multiphysics_integrator/test_multiphysics_integrator.cpp (1)
86-102: Verify test coverage for MHD model data allocation.The removal of MHD model allocation tests could lead to reduced test coverage. Please ensure that these tests are either:
- Replaced with equivalent tests in another location
- No longer needed due to architectural changes
✅ Verification successful
Commented out MHD tests are appropriately handled
The removal of these MHD allocation tests is justified as:
- Similar resource allocation testing is now handled by the hybrid model tests in the same file
- Comprehensive MHD state testing exists in dedicated test files under
tests/core/data/mhd_state/andtests/core/numerics/mock_mhd_simulator/🏁 Scripts executed
The following scripts were executed for the analysis:
Script:
#!/bin/bash # Search for any other MHD model allocation tests rg -l "MHDModel.*allocation" tests/Length of output: 37
Script:
#!/bin/bash # Search for MHD-related tests with broader patterns echo "=== Searching for MHD tests ===" rg -l "MHD" tests/ echo -e "\n=== Searching for B and V state variable tests ===" rg -A 5 "state\.(B|V)" tests/ echo -e "\n=== Searching for resource allocation patterns ===" rg -A 5 "resourcesManager.*getIDs" tests/ echo -e "\n=== Looking for test files with 'multiphysics' in name ===" fd "test.*multiphysics.*\.(cpp|hpp)" tests/Length of output: 9504
src/core/data/grid/gridlayoutimplyee.hpp (2)
39-39: Good type alias addition and naming consistency!The addition of quantity_type and consistent naming with underscore prefix improves code clarity and maintainability.
Also applies to: 108-111
41-55: Verify the status of LinearCombination implementation.A significant block of LinearCombination-related code is commented out. Please verify:
- If this functionality has been moved elsewhere
- If this is still needed for the MHD solver implementation
✅ Verification successful
LinearCombination implementation is properly maintained
The commented-out LinearCombination declarations are fully implemented in the codebase:
- Type definition exists in
gridlayoutdefs.hpp- Complete implementations are present in
gridlayoutimplyee.hppandgridlayoutimplyee_mhd.hpp- All functionality is actively used for grid layout calculations and field transformations
🏁 Scripts executed
The following scripts were executed for the analysis:
Script:
#!/bin/bash # Search for LinearCombination implementations rg -l "class.*LinearCombination|struct.*LinearCombination" src/Length of output: 65
Script:
#!/bin/bash # Search for any mention of LinearCombination rg "LinearCombination" --type cpp # Search for similar naming patterns that might indicate a rename rg "Combination|Linear.*Combination" --type cpp # Check recent changes in the specific file git log -p -- src/core/data/grid/gridlayoutimplyee.hpp | head -n 100Length of output: 5383
Script:
#!/bin/bash # Check WeightPoint definition and usage rg "WeightPoint" --type cpp -A 5 # Look for actual implementations of the linear combination methods ast-grep --pattern 'GridLayoutImpl<$_>::$_ToMoment() const { $$$ }'Length of output: 54686
tests/core/data/mhd_state/CMakeLists.txt (1)
1-21: Well-structured CMake configuration!The CMake configuration properly sets up the test-mhdstate project with:
- Correct minimum CMake version
- Proper include directories
- Required library dependencies
- Appropriate test registration
src/python3/CMakeLists.txt (1)
5-16: LGTM! Module configuration follows project standards.The new pyMHD module is properly configured with appropriate include directories, library linkage, and compile options, maintaining consistency with other modules in the project.
| if self.dim == 1: | ||
| return lambda x: value + x * 0 | ||
| if self.dim == 2: | ||
| return lambda x, y: value | ||
| if self.dim == 3: | ||
| return lambda x, y, z: value |
There was a problem hiding this comment.
Ensure defaulter method handles all dimensions to avoid unintended None returns
The defaulter method does not handle cases where self.dim is not 1, 2, or 3, which could lead to an implicit None return. This might cause unexpected errors if the method's output is used elsewhere.
Add an else clause to handle unsupported dimensions explicitly:
if self.dim == 1:
return lambda x: value + x * 0
if self.dim == 2:
return lambda x, y: value
if self.dim == 3:
return lambda x, y, z: value
+ else:
+ raise ValueError(f"Unsupported dimension: {self.dim}")📝 Committable suggestion
‼️ IMPORTANT
Carefully review the code before committing. Ensure that it accurately replaces the highlighted code, contains no missing lines, and has no issues with indentation. Thoroughly test & benchmark the code to ensure it meets the requirements.
| if self.dim == 1: | |
| return lambda x: value + x * 0 | |
| if self.dim == 2: | |
| return lambda x, y: value | |
| if self.dim == 3: | |
| return lambda x, y, z: value | |
| if self.dim == 1: | |
| return lambda x: value + x * 0 | |
| if self.dim == 2: | |
| return lambda x, y: value | |
| if self.dim == 3: | |
| return lambda x, y, z: value | |
| else: | |
| raise ValueError(f"Unsupported dimension: {self.dim}") |
| assert len(params) | ||
| param_per_dim = len(params) == self.dim | ||
| has_vargs = params[0].kind == inspect.Parameter.VAR_POSITIONAL | ||
| assert param_per_dim or has_vargs |
There was a problem hiding this comment.
Replace assert statements with explicit exceptions for parameter validation
Using assert statements for input validation is discouraged in production code because they can be disabled with optimization flags (e.g., -O), potentially skipping essential checks.
Replace assert statements with explicit exception handling to ensure validations are always enforced:
params = list(inspect.signature(input).parameters.values())
- assert len(params)
+ if not len(params):
+ raise ValueError("Input function must have at least one parameter")
param_per_dim = len(params) == self.dim
has_vargs = params[0].kind == inspect.Parameter.VAR_POSITIONAL
- assert param_per_dim or has_vargs
+ if not (param_per_dim or has_vargs):
+ raise ValueError("Function parameters do not match the model's dimensions")📝 Committable suggestion
‼️ IMPORTANT
Carefully review the code before committing. Ensure that it accurately replaces the highlighted code, contains no missing lines, and has no issues with indentation. Thoroughly test & benchmark the code to ensure it meets the requirements.
| assert len(params) | |
| param_per_dim = len(params) == self.dim | |
| has_vargs = params[0].kind == inspect.Parameter.VAR_POSITIONAL | |
| assert param_per_dim or has_vargs | |
| params = list(inspect.signature(input).parameters.values()) | |
| if not len(params): | |
| raise ValueError("Input function must have at least one parameter") | |
| param_per_dim = len(params) == self.dim | |
| has_vargs = params[0].kind == inspect.Parameter.VAR_POSITIONAL | |
| if not (param_per_dim or has_vargs): | |
| raise ValueError("Function parameters do not match the model's dimensions") |
| auto const& rhoVx = rhoV(Component::X); | ||
| auto const& rhoVy = rhoV(Component::Y); | ||
| auto const& rhoVz = rhoV(Component::Z); | ||
|
|
||
| auto& Vx = V(Component::X); | ||
| auto& Vy = V(Component::Y); | ||
| auto& Vz = V(Component::Z); | ||
|
|
||
| auto rhov_to_v = [&](const Field& rhoV, Field& V, MeshIndex<Field::dimension> index) { | ||
| V(index) = rhoV(index) / rho(index); | ||
| }; | ||
|
|
||
| layout_->evalOnBox(Vx, [&](auto&... args) mutable { rhov_to_v(rhoVx, Vx, {args...}); }); | ||
| layout_->evalOnBox(Vy, [&](auto&... args) mutable { rhov_to_v(rhoVy, Vy, {args...}); }); | ||
| layout_->evalOnBox(Vz, [&](auto&... args) mutable { rhov_to_v(rhoVz, Vz, {args...}); }); | ||
| } |
There was a problem hiding this comment.
Prevent division by zero in rhoVToV method
When computing velocities, dividing by zero may occur if rho(index) is zero, leading to runtime errors.
Add checks to handle zero density values:
auto rhov_to_v = [&](const Field& rhoV, Field& V, MeshIndex<Field::dimension> index) {
+ auto rho_val = rho(index);
+ if (rho_val != 0.0) {
V(index) = rhoV(index) / rho_val;
+ } else {
+ // Handle zero density case appropriately
+ V(index) = 0.0; // or another suitable default
+ }
};📝 Committable suggestion
‼️ IMPORTANT
Carefully review the code before committing. Ensure that it accurately replaces the highlighted code, contains no missing lines, and has no issues with indentation. Thoroughly test & benchmark the code to ensure it meets the requirements.
| auto const& rhoVx = rhoV(Component::X); | |
| auto const& rhoVy = rhoV(Component::Y); | |
| auto const& rhoVz = rhoV(Component::Z); | |
| auto& Vx = V(Component::X); | |
| auto& Vy = V(Component::Y); | |
| auto& Vz = V(Component::Z); | |
| auto rhov_to_v = [&](const Field& rhoV, Field& V, MeshIndex<Field::dimension> index) { | |
| V(index) = rhoV(index) / rho(index); | |
| }; | |
| layout_->evalOnBox(Vx, [&](auto&... args) mutable { rhov_to_v(rhoVx, Vx, {args...}); }); | |
| layout_->evalOnBox(Vy, [&](auto&... args) mutable { rhov_to_v(rhoVy, Vy, {args...}); }); | |
| layout_->evalOnBox(Vz, [&](auto&... args) mutable { rhov_to_v(rhoVz, Vz, {args...}); }); | |
| } | |
| auto const& rhoVx = rhoV(Component::X); | |
| auto const& rhoVy = rhoV(Component::Y); | |
| auto const& rhoVz = rhoV(Component::Z); | |
| auto& Vx = V(Component::X); | |
| auto& Vy = V(Component::Y); | |
| auto& Vz = V(Component::Z); | |
| auto rhov_to_v = [&](const Field& rhoV, Field& V, MeshIndex<Field::dimension> index) { | |
| auto rho_val = rho(index); | |
| if (rho_val != 0.0) { | |
| V(index) = rhoV(index) / rho_val; | |
| } else { | |
| // Handle zero density case appropriately | |
| V(index) = 0.0; // or another suitable default | |
| } | |
| }; | |
| layout_->evalOnBox(Vx, [&](auto&... args) mutable { rhov_to_v(rhoVx, Vx, {args...}); }); | |
| layout_->evalOnBox(Vy, [&](auto&... args) mutable { rhov_to_v(rhoVy, Vy, {args...}); }); | |
| layout_->evalOnBox(Vz, [&](auto&... args) mutable { rhov_to_v(rhoVz, Vz, {args...}); }); | |
| } |
| auto vx = rhoVx(index) / rho(index); | ||
| auto vy = rhoVy(index) / rho(index); | ||
| auto vz = rhoVz(index) / rho(index); | ||
| auto bx = GridLayout::project(Bx, index, GridLayout::faceXToCellCenter()); | ||
| auto by = GridLayout::project(By, index, GridLayout::faceYToCellCenter()); | ||
| auto bz = GridLayout::project(Bz, index, GridLayout::faceZToCellCenter()); | ||
| auto v2 = vx * vx + vy * vy + vz * vz; | ||
| auto b2 = bx * bx + by * by + bz * bz; | ||
| P(index) = (gamma_ - 1.0) * (Etot(index) - 0.5 * rho(index) * v2 - 0.5 * b2); | ||
| }; |
There was a problem hiding this comment.
Handle potential division by zero in eosEtotToP method
Divisions by rho(index) without checks can cause runtime errors if rho(index) is zero.
Include zero-density checks to prevent exceptions:
auto eos_etot_to_p = [&](MeshIndex<Field::dimension> index) {
+ auto rho_val = rho(index);
+ if (rho_val != 0.0) {
auto vx = rhoVx(index) / rho_val;
auto vy = rhoVy(index) / rho_val;
auto vz = rhoVz(index) / rho_val;
+ } else {
+ // Handle zero density case appropriately
+ auto vx = 0.0;
+ auto vy = 0.0;
+ auto vz = 0.0;
+ }
auto bx = GridLayout::project(Bx, index, GridLayout::faceXToCellCenter());
auto by = GridLayout::project(By, index, GridLayout::faceYToCellCenter());
auto bz = GridLayout::project(Bz, index, GridLayout::faceZToCellCenter());
auto v2 = vx * vx + vy * vy + vz * vz;
auto b2 = bx * bx + by * by + bz * bz;
P(index) = (gamma_ - 1.0) * (Etot(index) - 0.5 * rho(index) * v2 - 0.5 * b2);
};Committable suggestion skipped: line range outside the PR's diff.
| void initialize(GridLayout const& layout) | ||
| { | ||
| FieldUserFunctionInitializer::initialize(rho, layout, rhoinit_); | ||
| Vinit_.initialize(V, layout); | ||
| Binit_.initialize(B, layout); | ||
| FieldUserFunctionInitializer::initialize(P, layout, Pinit_); | ||
| } |
There was a problem hiding this comment.
🛠️ Refactor suggestion
Ensure all state variables are initialized
In the initialize method, only rho, V, B, and P are initialized. However, rhoV, Etot, E, and J are also member variables and should be initialized to ensure the state is fully set up.
Consider adding initialization for the additional state variables. Apply this diff to include their initialization:
void initialize(GridLayout const& layout)
{
FieldUserFunctionInitializer::initialize(rho, layout, rhoinit_);
Vinit_.initialize(V, layout);
Binit_.initialize(B, layout);
FieldUserFunctionInitializer::initialize(P, layout, Pinit_);
+ // Initialize rhoV, Etot, E, and J
+ rhoV.allocate(layout);
+ Etot.allocate(layout);
+ E.allocate(layout);
+ J.allocate(layout);
}Committable suggestion skipped: line range outside the PR's diff.
| return sim | ||
|
|
||
| def main(): | ||
| MHDMockSimulator(config()).run("alfven2d.h5") |
There was a problem hiding this comment.
🛠️ Refactor suggestion
Add dump frequency parameter
Unlike other test cases, the dump frequency is not specified. This might lead to excessive output.
- MHDMockSimulator(config()).run("alfven2d.h5")
+ MHDMockSimulator(config()).run("alfven2d.h5", dumpfrequency=100)📝 Committable suggestion
‼️ IMPORTANT
Carefully review the code before committing. Ensure that it accurately replaces the highlighted code, contains no missing lines, and has no issues with indentation. Thoroughly test & benchmark the code to ensure it meets the requirements.
| MHDMockSimulator(config()).run("alfven2d.h5") | |
| MHDMockSimulator(config()).run("alfven2d.h5", dumpfrequency=100) |
| def bx(x, y): | ||
| w1 = 0.2 | ||
| w2 = 1.0 | ||
| x0 = x - 0.5 * Lx | ||
| y1 = y - 0.3 * Ly | ||
| y2 = y - 0.7 * Ly | ||
| w3 = np.exp(-(x0 * x0 + y1 * y1) / (w2 * w2)) | ||
| w4 = np.exp(-(x0 * x0 + y2 * y2) / (w2 * w2)) | ||
| w5 = 2.0 * w1 / w2 | ||
| v1 = -1 | ||
| v2 = 1.0 | ||
| return ( | ||
| v1 | ||
| + (v2 - v1) * (S(y, Ly * 0.3, 0.5) - S(y, Ly * 0.7, 0.5)) | ||
| + (-w5 * y1 * w3) | ||
| + (+w5 * y2 * w4) | ||
| ) | ||
|
|
There was a problem hiding this comment.
🛠️ Refactor suggestion
Consider extracting common magnetic field calculations
The bx and by functions share significant calculation logic. Consider extracting common calculations to reduce code duplication.
def _compute_magnetic_params(x, y, Lx, Ly):
w1, w2 = 0.2, 1.0
x0 = x - 0.5 * Lx
y1 = y - 0.3 * Ly
y2 = y - 0.7 * Ly
w3 = np.exp(-(x0 * x0 + y1 * y1) / (w2 * w2))
w4 = np.exp(-(x0 * x0 + y2 * y2) / (w2 * w2))
w5 = 2.0 * w1 / w2
return x0, y1, y2, w3, w4, w5Also applies to: 69-79
| def p(x, y): | ||
| return 1.0 - (bx(x, y) ** 2 + by(x, y) ** 2) / 2.0 |
There was a problem hiding this comment.
💡 Codebase verification
Add assertion to prevent negative pressure
Add validation to ensure pressure remains positive, similar to the temperature validation in the Harris sheet implementations:
def p(x, y):
pressure = 1.0 - (bx(x, y) ** 2 + by(x, y) ** 2) / 2.0
assert pressure > 0, "Pressure must be positive"
return pressure🔗 Analysis chain
Add validation for pressure calculation
The pressure calculation should ensure it doesn't become negative due to strong magnetic fields.
🏁 Scripts executed
The following scripts were executed for the analysis:
Script:
#!/bin/bash
# Search for similar pressure calculations in the codebase
ast-grep --pattern 'p($_, $_) { return $_ - ($_ ** 2 + $_ ** 2) / $_; }'
Length of output: 74
Script:
#!/bin/bash
# Search for pressure-related functions and their implementations
ast-grep --pattern 'def p($$$) { $$$ }'
# Look for test files related to pressure calculations
fd test.*mhd.*\.py
# Search for validation patterns in numeric calculations
rg -l "assert.*>.*0|if.*<.*0" tests/
Length of output: 769
Script:
#!/bin/bash
# Check Harris test files for pressure-related content
cat tests/functional/harris/harris_2d.py tests/functional/harris/harris_2d_lb.py
# Search for pressure-related assertions or validations
rg "assert.*p.*>|pressure.*>" tests/functional/harris/
rg "assert.*p.*>|pressure.*>" tests/core/numerics/mock_mhd_simulator/
Length of output: 10138
| ToConservativeConverter(PHARE::initializer::PHAREDict const& dict) | ||
| : gamma_{dict["heat_capacity_ratio"].template to<double>()} | ||
| { | ||
| } |
There was a problem hiding this comment.
Add validation for gamma value.
The heat capacity ratio (gamma) should be validated to ensure it's physically meaningful (typically > 1).
ToConservativeConverter(PHARE::initializer::PHAREDict const& dict)
: gamma_{dict["heat_capacity_ratio"].template to<double>()}
{
+ if (gamma_ <= 1.0) {
+ throw std::invalid_argument("Heat capacity ratio must be greater than 1");
+ }
}📝 Committable suggestion
‼️ IMPORTANT
Carefully review the code before committing. Ensure that it accurately replaces the highlighted code, contains no missing lines, and has no issues with indentation. Thoroughly test & benchmark the code to ensure it meets the requirements.
| ToConservativeConverter(PHARE::initializer::PHAREDict const& dict) | |
| : gamma_{dict["heat_capacity_ratio"].template to<double>()} | |
| { | |
| } | |
| ToConservativeConverter(PHARE::initializer::PHAREDict const& dict) | |
| : gamma_{dict["heat_capacity_ratio"].template to<double>()} | |
| { | |
| if (gamma_ <= 1.0) { | |
| throw std::invalid_argument("Heat capacity ratio must be greater than 1"); | |
| } | |
| } |
| template<typename Field, typename VecField> | ||
| void eosPToEtot(const Field& rho, const VecField& V, const VecField& B, const Field& P, | ||
| Field& Etot) const | ||
| { | ||
| auto const& Vx = V(Component::X); | ||
| auto const& Vy = V(Component::Y); | ||
| auto const& Vz = V(Component::Z); | ||
|
|
||
| auto const& Bx = B(Component::X); | ||
| auto const& By = B(Component::Y); | ||
| auto const& Bz = B(Component::Z); | ||
|
|
||
| auto eos_p_to_etot = [&](MeshIndex<Field::dimension> index) { | ||
| auto bx = GridLayout::project(Bx, index, GridLayout::faceXToCellCenter()); | ||
| auto by = GridLayout::project(By, index, GridLayout::faceYToCellCenter()); | ||
| auto bz = GridLayout::project(Bz, index, GridLayout::faceZToCellCenter()); | ||
| auto v2 = Vx(index) * Vx(index) + Vy(index) * Vy(index) + Vz(index) * Vz(index); | ||
| auto b2 = bx * bx + by * by + bz * bz; | ||
| Etot(index) = P(index) / (gamma_ - 1.0) + 0.5 * rho(index) * v2 + 0.5 * b2; | ||
| }; | ||
|
|
||
| layout_->evalOnBox(Etot, [&](auto&... args) mutable { eos_p_to_etot({args...}); }); | ||
| } |
There was a problem hiding this comment.
Add numerical stability checks in eosPToEtot.
The energy calculation could potentially overflow or produce NaN results.
auto eos_p_to_etot = [&](MeshIndex<Field::dimension> index) {
auto bx = GridLayout::project(Bx, index, GridLayout::faceXToCellCenter());
auto by = GridLayout::project(By, index, GridLayout::faceYToCellCenter());
auto bz = GridLayout::project(Bz, index, GridLayout::faceZToCellCenter());
auto v2 = Vx(index) * Vx(index) + Vy(index) * Vy(index) + Vz(index) * Vz(index);
auto b2 = bx * bx + by * by + bz * bz;
+ // Check for potential overflow or NaN
+ if (!std::isfinite(v2) || !std::isfinite(b2) || !std::isfinite(P(index)) ||
+ !std::isfinite(rho(index))) {
+ throw std::runtime_error("Non-finite values detected in energy calculation");
+ }
Etot(index) = P(index) / (gamma_ - 1.0) + 0.5 * rho(index) * v2 + 0.5 * b2;
};Committable suggestion skipped: line range outside the PR's diff.
Usable MHD solver with:
Summary by CodeRabbit
Based on the comprehensive summary of changes, here are the concise release notes:
New Features
Improvements
Test Coverage
Documentation