[SPH] Hutchison 2018 Multigrain

[tdavidcl]

No description provided.

[github-actions]

Thanks @tdavidcl for opening this PR!

You can do multiple things directly here:
1 - Comment pre-commit.ci run to run pre-commit checks.
2 - Comment pre-commit.ci autofix to apply fixes.
3 - Add label autofix.ci to fix authorship & pre-commit for every commit made.
4 - Add label light-ci to only trigger a reduced & faster version of the CI (need the full one before merge).
5 - Add label trigger-ci to create an empty commit to trigger the CI.

Once the workflow completes a message will appear displaying informations related to the run.

Also the PR gets automatically reviewed by gemini, you can:
1 - Comment /gemini review to trigger a review
2 - Comment /gemini summary for a summary
3 - Tag it using @gemini-code-assist either in the PR or in review comments on files

[gemini-code-assist]

Summary of Changes

Hello @tdavidcl, I'm Gemini Code Assist¹! I'm currently reviewing this pull request and will post my feedback shortly. In the meantime, here's a summary to help you and other reviewers quickly get up to speed!

This pull request integrates a new multigrain dust model into the SPH solver, specifically focusing on the MonofluidTVI configuration. It introduces new particle fields and implements the numerical calculation for the time evolution of a key dust-related variable, S_j, following the methodology outlined in the Hutchison 2018 paper. Additionally, it enhances the solver's configuration safety by restricting experimental dust features.

Highlights

• New Dust Model Field: Introduced new particle data fields, 's_j' and 'ds_j_dt', for the MonofluidTVI dust model, representing a key variable and its time derivative.
• Derivative Calculation Implementation: Implemented the update_derivs_dust_monofluid_tvi_Sj function to compute the time derivative of the 's_j' field, following Equation 51 from the Hutchison 2018 paper.
• Experimental Feature Safeguard: Added a check in the solver configuration to prevent the use of experimental dust configurations unless explicitly allowed, throwing a runtime error otherwise.
• Configuration Helper Methods: Added utility methods is_monofluid_tvi, is_monofluid_complete, and has_Sj_field to SolverConfig.hpp to query the active dust model configuration.

Changelog

• src/shammodels/sph/include/shammodels/sph/SolverConfig.hpp
  • Added is_monofluid_tvi(), is_monofluid_complete(), and has_Sj_field() methods to query the current dust configuration.
  • Modified check_config() to throw a std::runtime_error if experimental dust configurations are used without shamrock::are_experimental_features_allowed() being true.
• src/shammodels/sph/include/shammodels/sph/modules/UpdateDerivs.hpp
  • Declared the update_derivs_dust_monofluid_tvi_Sj(DustConfig cfg) function for calculating dust derivatives.
• src/shammodels/sph/src/SolverConfig.cpp
  • Added 's_j' and 'ds_j_dt' fields to the PatchDataLayer and 's_j' to the ghost layout fields when dust_config.has_Sj_field() is true.
• src/shammodels/sph/src/modules/UpdateDerivs.cpp
  • Integrated the call to update_derivs_dust_monofluid_tvi_Sj into the main update_derivs() function if the dust configuration has the 'Sj' field.
  • Implemented the update_derivs_dust_monofluid_tvi_Sj function, which computes ds_j_dt based on Equation 51 from Hutchison 2018, utilizing various particle properties and neighbor interactions.

Activity

• New functionality for a specific dust model has been introduced.
• Modifications span across header and source files related to solver configuration and derivative updates.
• The changes include the addition of new data fields, helper methods, and a core computational function for dust dynamics.

Using Gemini Code Assist

The full guide for Gemini Code Assist can be found on our documentation page, here are some quick tips.

Invoking Gemini

You can request assistance from Gemini at any point by creating a comment using either /gemini <command> or @gemini-code-assist <command>. Below is a summary of the supported commands on the current page.

Feature	Command	Description
Code Review	/gemini review	Performs a code review for the current pull request in its current state.
Pull Request Summary	/gemini summary	Provides a summary of the current pull request in its current state.
Comment	@gemini-code-assist	Responds in comments when explicitly tagged, both in pull request comments and review comments.
Help	/gemini help	Displays a list of available commands.

Customization

To customize Gemini Code Assist for GitHub experience, repository maintainers can create a configuration file and/or provide a custom code review style guide (such as PEP-8 for Python) by creating and adding files to a .gemini/ folder in the base of the repository. Detailed instructions can be found here.

Limitations & Feedback

Gemini Code Assist may make mistakes. Please leave feedback on any instances where its feedback is incorrect or counter productive. You can react with 👍 and 👎 on @gemini-code-assist comments. If you're interested in giving your feedback about your experience with Gemini Code Assist for Github and other Google products, sign up here.

You can also get AI-powered code generation, chat, as well as code reviews directly in the IDE at no cost with the Gemini Code Assist IDE Extension.

Footnotes

1. Review the Privacy Notices, Generative AI Prohibited Use Policy, Terms of Service, and learn how to configure Gemini Code Assist in GitHub here. Gemini can make mistakes, so double check it and use code with caution. ↩

[gemini-code-assist]

Code Review

This pull request introduces functionality for multigrain dust based on the Hutchison (2018) model. The changes span configuration and solver logic, adding new fields and a dedicated derivative update function. My review has identified a critical issue in the new derivative calculation where a key term is hardcoded to zero, which will lead to incorrect results. Additionally, there are several instances of unused variables and redundant code, likely from copy-pasting, which impact readability and performance. I've provided suggestions to clean this up and improve maintainability.

[github-actions]

Workflow report

workflow report corresponding to commit 4214b55
Commiter email is timothee.davidcleris@proton.me

Light CI is enabled. This will only run the basic tests and not the full tests.
Merging a PR require the job "on PR / all" to pass which is disabled in this case.

Pre-commit check report

Some failures were detected in base source checks checks.
Check the On PR / Linting / Base source checks (pull_request) job in the tests for more detailled output

❌ ruff-format

1 file reformatted, 200 files left unchanged

Suggested changes

Detailed changes :

diff --git a/examples/sph/run_dusty_settle.py b/examples/sph/run_dusty_settle.py
index 5e595f71..3cb3b0ad 100644
--- a/examples/sph/run_dusty_settle.py
+++ b/examples/sph/run_dusty_settle.py
@@ -173,6 +173,7 @@ model.set_cfl_force(0.1)
 
 model.timestep()
 
+
 def compute_sj_new(patchdata):
     hpart = patchdata["hpart"]
     rho = pmass * (model.get_hfact() / np.array(hpart)) ** 3
@@ -183,10 +184,6 @@ def compute_sj_new(patchdata):
     return s
 
 
-
-
-
-
 # TODO: add function to modify fields e.g. get rho and do stuff according to it
 
 tnext = 0
@@ -195,9 +192,9 @@ for j in range(1000):
         tnext += 0.02
         model.evolve_until(tnext)
 
-        if(j == 70):
+        if j == 70:
             for k in range(ndust):
-                model.overwrite_field_value_f64("s_j",compute_sj_new,k)
+                model.overwrite_field_value_f64("s_j", compute_sj_new, k)
 
     dic = ctx.collect_data()
     print(dic["s_j"])
diff --git a/src/shammodels/sph/include/shammodels/sph/Model.hpp b/src/shammodels/sph/include/shammodels/sph/Model.hpp
index 6e43a7bb..334eff75 100644
--- a/src/shammodels/sph/include/shammodels/sph/Model.hpp
+++ b/src/shammodels/sph/include/shammodels/sph/Model.hpp
@@ -981,9 +981,9 @@ namespace shammodels::sph {
                         * (1.
                            + sycl::sin(shambase::constants::pi<Tscal> / (2. * Hwarp) * (R - Rwarp)))
                         * sycl::sin(incl_rad));
-                    psi = shambase::constants::pi<Tscal>
-                          * Rwarp / (4. * Hwarp) * sycl::sin(incl_rad)
-                          / sycl::sqrt(1. - (0.5 * sycl::pow(sycl::sin(incl_rad), 2)));
+                    psi          = shambase::constants::pi<Tscal>
+                                   * Rwarp / (4. * Hwarp) * sycl::sin(incl_rad)
+                                   / sycl::sqrt(1. - (0.5 * sycl::pow(sycl::sin(incl_rad), 2)));
                     Tscal psimax = sycl::max(psimax, psi);
                     Tscal x      = pos[i].x();
                     Tscal y      = pos[i].y();