Add finite-volume Euler functionals example#1626
Conversation
Greptile SummaryThis PR introduces a new finite-volume Euler forward-step example under
Important Files Changed
Reviews (1): Last reviewed commit: "Add finite-volume Euler functionals exam..." | Re-trigger Greptile |
![greptile-apps[bot]](https://avatars.githubusercontent.com/in/867647?s=60&v=4){kind=small}
| phi_3: float, | ||
| phi_4: float, | ||
| ): | ||
| rho = wp.max(values[cell, 0], 1.0e-30) |
There was a problem hiding this comment.
_predicted_value is hardcoded to 1.0e-30 instead of using the density_floor parameter that is available in the enclosing _euler_update_cell call. The Warp update kernel passes density_floor to _euler_update_cell and correctly applies it to the reconstructed face values (rho_l = wp.max(rho_l, density_floor)), but the prediction step silently uses a different constant. This means a user who sets density_floor=1e-4 (for example) will still divide by densities as small as 1e-30 during the half-step prediction, producing a divergence with the Torch path and potentially NaN velocities.
| rho = wp.max(values[cell, 0], 1.0e-30) | |
| rho = wp.max(values[cell, 0], density_floor) |
| @wp.func | ||
| def _predicted_value( | ||
| values: wp.array2d(dtype=wp.float32), | ||
| gradients: wp.array3d(dtype=wp.float32), | ||
| cell: int, | ||
| comp: int, | ||
| dt: float, | ||
| gamma: float, | ||
| dims: int, | ||
| phi_0: float, | ||
| phi_1: float, | ||
| phi_2: float, | ||
| phi_3: float, | ||
| phi_4: float, | ||
| ): |
There was a problem hiding this comment.
_predicted_value Warp function signature does not accept density_floor and pressure_floor, so the fix above also requires threading those parameters through. Additionally, _predicted_value does not clamp the pressure contribution — in the Torch path, W_half[:, -1] is clamped with clamp_min(pressure_floor) after the prediction, but _predicted_value has no equivalent guard before computing the energy derivative (gamma * pressure * div_velocity). Consider adding both parameters to the signature so the prediction step is consistent with the rest of the Warp kernel.
| @wp.func | |
| def _predicted_value( | |
| values: wp.array2d(dtype=wp.float32), | |
| gradients: wp.array3d(dtype=wp.float32), | |
| cell: int, | |
| comp: int, | |
| dt: float, | |
| gamma: float, | |
| dims: int, | |
| phi_0: float, | |
| phi_1: float, | |
| phi_2: float, | |
| phi_3: float, | |
| phi_4: float, | |
| ): | |
| @wp.func | |
| def _predicted_value( | |
| values: wp.array2d(dtype=wp.float32), | |
| gradients: wp.array3d(dtype=wp.float32), | |
| cell: int, | |
| comp: int, | |
| dt: float, | |
| gamma: float, | |
| dims: int, | |
| density_floor: float, | |
| pressure_floor: float, | |
| phi_0: float, | |
| phi_1: float, | |
| phi_2: float, | |
| phi_3: float, | |
| phi_4: float, | |
| ): |
|
|
||
| @hydra.main(version_base="1.3", config_path="conf", config_name="config") | ||
| def main(cfg: DictConfig) -> None: | ||
| output_dir = Path("outputs") |
There was a problem hiding this comment.
hydra.job.chdir: True and hydra.run.dir: ./outputs_finite_volume_euler, Hydra changes the working directory to outputs_finite_volume_euler/ before calling main. Path("outputs") is then resolved relative to that directory, so all output files land in outputs_finite_volume_euler/outputs/ rather than directly in outputs_finite_volume_euler/ as the README describes. Using Path(".") writes directly into the Hydra run directory and matches the documented behaviour.
| output_dir = Path("outputs") | |
| output_dir = Path(".") |
|
|
||
| right_cells = cell_neighbors.to(dtype=torch.int64).clamp_min(0) | ||
| right = cell_values[right_cells].clone() | ||
| interior = cell_neighbors >= 0 | ||
| if bool(torch.any(interior).item()): | ||
| neighbor_dx = face_centroids - cell_centroids[right_cells] | ||
| right[interior] = ( | ||
| right[interior] | ||
| + torch.einsum("cfd,cfdk->cfk", neighbor_dx, cell_gradients[right_cells])[ | ||
| interior | ||
| ] | ||
| ) | ||
| return left, right | ||
|
|
||
|
|
||
| def euler_cfl_timestep_torch( | ||
| cell_values: torch.Tensor, | ||
| points: torch.Tensor, | ||
| cells: torch.Tensor, | ||
| cell_neighbors: torch.Tensor, | ||
| boundary_tags: torch.Tensor, | ||
| inflow_state: torch.Tensor, | ||
| gamma: float, | ||
| cfl: float, | ||
| density_floor: float, |
There was a problem hiding this comment.
In _euler_primitive_time_derivative_torch, rho = W[:, 0] is used without any floor clamping when computing the pressure-gradient term grad_pressure[:, dim] / rho. The Warp path explicitly guards against near-zero density in _predicted_value (even if with a hardcoded constant), but the Torch path does not. If a cell's density drops below the numerical noise floor before W_half[:, 0].clamp_min(density_floor) is applied, this division can produce very large or inf velocity increments. Adding rho = W[:, 0].clamp_min(density_floor) (passing density_floor into the function) would make the two paths consistent.
loliverhennigh
force-pushed
the
finite-volume-euler-example
branch
from
98ec53c to
8898696
Compare
|
@peterdsharpe This is something I was playing around with using the Mesh module, basically combining |
|
This is a very cool proof-of-concept, thanks for sharing Oliver! The If you later decide that you want to merge some variant of this (and have the bandwidth to maintain it), let me know and I'm happy to give this PR a review. |
2 participants
Summary
examples/functionals/finite_volume_eulerValidation
uv run python -m ruff format examples/functionals/finite_volume_euler/euler_solver.py examples/functionals/finite_volume_euler/euler_finite_volume.py test/nn/functional/finite_volume/test_euler_example.pyuv run python -m ruff check examples/functionals/finite_volume_euler/euler_solver.py examples/functionals/finite_volume_euler/euler_finite_volume.py test/nn/functional/finite_volume/test_euler_example.pyuv run python -m pytest test/nn/functional/finite_volume -quv run python -m pytest test/nn/functional -q