Merge master into MTHINC: resolve conflicts

- case.md: add muscl_eps row/description alongside int_comp (MTHINC version)
- m_global_parameters.fpp: keep int_comp as integer (0/1/2), add muscl_eps; keep GPU_DECLARE; both GPU_UPDATEs
- m_mpi_proxy.fpp: broadcast both int_comp and collision_model as integers
- case_validator.py: keep recon_type validation from master; drop stale int_comp=boolean check (superseded by check_interface_compression)

Author: sbryngelson

.github/workflows/test.yml

@@ -149,7 +149,7 @@ jobs:
         continue-on-error: true
 
   github:
-    name: ${{ matrix.nvhpc && format('NVHPC {0} ({1})', matrix.nvhpc, matrix.target) || format('Github ({0}, {1}, {2}, intel={3})', matrix.os, matrix.mpi, matrix.debug, matrix.intel) }}
+    name: ${{ matrix.nvhpc && format('NVHPC {0} ({1})', matrix.nvhpc, matrix.target) || format('Github ({0}, {1}, {2}, {3})', matrix.os, matrix.mpi, matrix.debug, matrix.intel && 'intel' || 'GNU') }}
     needs: [lint-gate, file-changes, rebuild-cache]
     if: >-
       !cancelled() &&

CMakeLists.txt

@@ -636,7 +636,7 @@ exit 0
                     target_compile_options(${a_target} PRIVATE -fopenmp)
                     target_link_options(${a_target} PRIVATE -fopenmp)
                 elseif(CMAKE_Fortran_COMPILER_ID STREQUAL "LLVMFlang")
-                    target_compile_options(${a_target} PRIVATE -fopenmp --offload-arch=gfx90a)
+                    target_compile_options(${a_target} PRIVATE -fopenmp --offload-arch=gfx90a -fopenmp-target-fast -fopenmp-assume-threads-oversubscription -fopenmp-assume-teams-oversubscription)
                     target_link_options(${a_target} PRIVATE -fopenmp --offload-arch=gfx90a)
                 endif()
             endif()
@@ -1026,4 +1026,4 @@ site_name(SITE_NAME)
 
 configure_file(
     "${CMAKE_CURRENT_SOURCE_DIR}/toolchain/cmake/configuration.cmake.in"
-    "${CMAKE_CURRENT_BINARY_DIR}/configuration.txt")
+    "${CMAKE_CURRENT_BINARY_DIR}/configuration.txt")

docs/documentation/case.md

@@ -331,6 +331,10 @@ This is enabled by adding ``'elliptic_smoothing': "T",`` and ``'elliptic_smoothi
 | `moving_ibm`         | Integer | Sets the method used for IB movement. |
 | `vel(i)`             | Real    | Initial velocity of the moving IB in the i-th direction. |
 | `angular_vel(i)`     | Real    | Initial angular velocity of the moving IB in the i-th direction. |
+| `coefficient_of_restitution`     | Real    | A number 0 to 1 describing how elastic IB collisions are |
+| `collision_model`     | Integer    | Integer to select the collision model being used for IB collisions. |
+| `collision_time`     | Real    | Amount of simulation time used to resolve collisions |
+| `ib_coefficient_of_friction`     | Real    | Coefficient of friction used in IB collisions |
 
 These parameters should be prepended with `patch_ib(j)%` where $j$ is the patch index.
 
@@ -361,6 +365,14 @@ Additional details on this specification can be found in [NACA airfoil](https://
 
 - `angular_vel(i)` is the initial angular velocity of the IB about the x, y, z axes for i=1, 2, 3 in radians per second. When `moving_ibm` equals 2, this rotation rate is just the starting rate of the object, which will then change due to external torques. If `moving_ibm` equals 1, then this is constant if it is a number, or can be described analytically with an expression. 
 
+- `coefficient_of_restitution` is a number from 0 (exclusive) to 1 (inclusive) describing how elastic IB collisions are. 0 is for perfectly inellastic collisions while 1 is for perfectly ellastic collisions.
+
+- `collision_model` is an integer to select the collision model being used for IB collisions. Using 0 disables collisions and collisiono checking. 1 enables the soft-sphere collision model, where all IBs must be circles or sphere and those IBs can collide with each other as well as walls.
+
+- `collision_time` is approximately the amount of simulation time used to resolve collisions. This is handled by modifying the spring gonstant used to apply collision forces.
+
+- `ib_coefficient_of_friction` is the coefficient of friction used in IB collisions.
+
 ### 5. Fluid Material's {#sec-fluid-materials}
 
 | Parameter | Type   | Description                                    |
@@ -433,6 +445,7 @@ See @ref equations "Equations" for the mathematical models these parameters cont
 | `mp_weno`                  | Logical | Monotonicity preserving WENO |
 | `muscl_order`              | Integer | MUSCL order [1,2] |
 | `muscl_lim`                | Integer | MUSCL Slope Limiter: [1] minmod; [2] monotonized central; [3] Van Albada; [4] Van Leer; [5] SUPERBEE |
+| `muscl_eps`                | Real    | MUSCL limiter slope-product threshold (default: hard-coded thresholds; set to 0 for textbook behavior) |
 | `int_comp`                 | Integer | Interface Compression [0] Off [1] THINC [2] MTHINC (default 0) |
 | `flux_lim`                 | Integer | Flux limiter for post-process: [1] minmod; [2] MUSCL; [3] OSPRE; [4] SUPERBEE |
 | `ic_eps`                   | Real    | Interface compression threshold (default: 1e-4) |
@@ -530,6 +543,10 @@ It is recommended to set `weno_eps` to $10^{-6}$ for WENO-JS, and to $10^{-40}$
 - `muscl_lim` specifies the slope limiter that is used in 2nd order MUSCL Reconstruction by an integer from 1 through 5.
 `muscl_lim = 1`, `2`, `3`, `4`, and `5` correspond to minmod, monotonized central, Van Albada, Van Leer, and SUPERBEE, respectively.
 
+- `muscl_eps` controls the slope-product activation threshold for all MUSCL limiters.
+When not set (default), the threshold is 1e-9 for minmod/MC, and 1e-6 for others.
+Setting `muscl_eps = 0` gives textbook limiter behavior where limiters activate whenever both slopes have the same sign.
+
 - `int_comp` activates interface compression using [1] THINC or [2] MTHINC (default off) used in variable reconstruction, with control parameters (`ic_eps`, and `ic_beta`).
 
 - `riemann_solver` specifies the choice of the Riemann solver that is used in simulation by an integer from 1 through 4.
@@ -1056,8 +1073,20 @@ When ``cyl_coord = 'T'`` is set in 2D the following constraints must be met:
 
 - `cantera_file` specifies the chemical mechanism file. If the file is part of the standard Cantera library, only the filename is required. Otherwise, the file must be located in the same directory as your `case.py` file
 
+### 18. Chemistry-Specific Boundary Conditions
+
+| Parameter          | Type    | Description                                                                 |
+| ---:               | :----:  | :---                                                                        |
+| `bc_[x,y,z]%%isothermal_in`    | Logical | Enable isothermal wall at the domain entrance (minimum coordinate).         |
+| `bc_[x,y,z]%%isothermal_out`   | Logical | Enable isothermal wall at the domain exit (maximum coordinate).             |
+| `bc_[x,y,z]%%Twall_in`         | Real    | Temperature [K] of the entrance isothermal wall.                            |
+| `bc_[x,y,z]%%Twall_out`        | Real    | Temperature [K] of the exit isothermal wall.                                |
+
+This boundary condition can be used for fixed-temperature (isothermal) walls at the domain extremities. It is exclusively available for reacting flows and requires chemistry to be enabled. It properly evaluates heat and species fluxes at the interface when ``chemistry = 'T'``, ``chem_params%%diffusion = 'T'``, and the corresponding domain boundary is set to a slip wall (`bc_[x,y,z]%%[beg,end]` = -15) or a no-slip wall (`bc_[x,y,z]%%[beg,end]` = -16).
+
+
 
-### 18. GPU Performance (NVIDIA UVM)
+### 19. GPU Performance (NVIDIA UVM)
 
 | Parameter                  | Type    | Description                                              |
 | ---:                       | :---:   | :---                                                     |

docs/documentation/visualization.md

@@ -205,7 +205,7 @@ Post-processed data in Silo-HDF5 format (`format=1`) can be opened directly in P
 Paraview 5.11.0 has been confirmed to work with the MFC databases for some parallel environments.
 Nevertheless, the installation and configuration of Paraview can be environment-dependent and are left to the user.
 
-The user can launch Paraview and open the index files under `/silo_hdf5/root`.
+The user can launch Paraview and open the series file `/silo_hdf5/collection.silo.series`.
 Once the database is loaded, flow field variables contained in the database can be added to the render view.
 Further information on Paraview can be found in its [documentation](https://docs.paraview.org/en/latest/).
 The figure below shows the iso-contour of the liquid void fraction (`alpha1`) in the database generated by the example case `3D_sphbubcollapse`.

docs/index.html

@@ -58,6 +58,7 @@
             { name: "Mach 0.3 flow over a corgi (2M STL)", image: "res/simulations/u.png", computer: "HiPerGator", computerUrl: "https://www.rc.ufl.edu/about/hipergator/",                                     accelerators: "2 GPUs",    walltime: "80s",  source: "https://www.youtube.com/watch?v=O8dSRqHLp_o" },
             // Shock-droplet
             { name: "Shedding water droplet",              image: "res/simulations/a.png", computer: "Summit",     computerUrl: "https://www.olcf.ornl.gov/summit/",                                            accelerators: "960 V100s", walltime: "4h",   source: "https://www.youtube.com/watch?v=Gjj-qZkXcrg" },
+            { name: "1K resolved particles",               image: "res/simulations/x.jpg", computer: "Phoenix",    computerUrl: "https://www.pace.gatech.edu/",                                                 accelerators: "4 A100s",   walltime: "~1.5h", source: "https://www.youtube.com/watch?v=ufisHG0KkOU" },
             // Biomedical & acoustics
             { name: "Burstwave lithotripsy",               image: "res/simulations/k.png", computer: "Delta",      computerUrl: "https://www.ncsa.illinois.edu/research/project-highlights/delta/",              accelerators: "128 A100s", walltime: "30m",  source: "https://www.youtube.com/watch?v=XWsUTaJXGF8" },
             { name: "Cavitation fragments kidney stone",   image: "res/simulations/d.png", computer: "Summit",     computerUrl: "https://www.olcf.ornl.gov/summit/",                                            accelerators: "576 V100s", walltime: "30m",  source: "https://doi.org/10.48550/arXiv.2305.09163" },
@@ -74,6 +75,7 @@
             { name: "Breakup of vibrated interface",       image: "res/simulations/f.png", computer: "Summit",     computerUrl: "https://www.olcf.ornl.gov/summit/",                                            accelerators: "128 V100s", walltime: "4h",   source: "https://www.youtube.com/watch?v=XQ3g1oSg8mc" },
             { name: "Viscous Taylor-Green vortex",         image: "res/simulations/h.png", computer: "Delta",      computerUrl: "https://www.ncsa.illinois.edu/research/project-highlights/delta/",              accelerators: "128 A100s", walltime: "17h",  source: "https://www.youtube.com/watch?v=7i2h08dlDQw" },
             { name: "Mach 1.5 shock-helium bubble",        image: "res/simulations/t.png", computer: "Phoenix",    computerUrl: "https://www.pace.gatech.edu/",                                                 accelerators: "1 A100",    walltime: "1h",   source: "https://www.youtube.com/watch?v=zDJoe0NYZsQ" },
+            { name: "Sphere collision and wall rebound",   image: "res/simulations/w.jpg", computer: "Phoenix",    computerUrl: "https://www.pace.gatech.edu/",                                                 accelerators: "4 A100s",   walltime: "~2h",   source: "https://www.youtube.com/watch?v=5ugLuDtV0C0" },
         ];
 
         const scalings = [

docs/module_categories.json

@@ -27,7 +27,8 @@
             "m_chemistry",
             "m_acoustic_src",
             "m_body_forces",
-            "m_pressure_relaxation"
+            "m_pressure_relaxation",
+            "m_collisions"
         ]
     },
     {

docs/res/simulations/w.jpg

@@ -0,0 +1,76 @@
+#!/usr/bin/env python3
+import json
+import math
+
+import cantera as ct
+
+Lx = 0.05
+Ly = 0.05
+
+ctfile = "h2o2.yaml"
+sol_L = ct.Solution(ctfile)
+sol_L.TPY = 1125, ct.one_atm, "O2:0.21,N2:0.79"
+# Configuring case dictionary
+case = {
+    "run_time_info": "T",
+    "x_domain%beg": 0.0,
+    "x_domain%end": Lx,
+    "y_domain%beg": 0.0,
+    "y_domain%end": Ly,
+    "m": 699,
+    "n": 699,
+    "p": 0,
+    "dt": 4.0e-08,
+    "t_step_start": 0,
+    "t_step_stop": 75000,
+    "t_step_save": 4500,
+    "num_patches": 1,
+    "model_eqns": 2,
+    "alt_soundspeed": "F",
+    "num_fluids": 1,
+    "mpp_lim": "F",
+    "mixture_err": "T",
+    "time_stepper": 3,
+    "mp_weno": "F",
+    "weno_order": 5,
+    "weno_eps": 1e-16,
+    "riemann_solver": 2,
+    "wave_speeds": 1,
+    "avg_state": 2,
+    "bc_x%beg": -7,
+    "bc_x%end": -3,
+    "bc_y%beg": -16,
+    "bc_y%end": -3,
+    "bc_y%isothermal_in": "T",
+    "bc_y%Twall_in": 600.0,
+    "format": 1,
+    "precision": 2,
+    "prim_vars_wrt": "T",
+    "parallel_io": "T",
+    "chemistry": "T",
+    "chem_params%diffusion": "T",
+    "chem_params%reactions": "F",
+    "cantera_file": ctfile,
+    "chem_wrt_T": "T",
+    "patch_icpp(1)%geometry": 3,
+    "patch_icpp(1)%hcid": 291,
+    "patch_icpp(1)%x_centroid": Lx / 2,
+    "patch_icpp(1)%y_centroid": Ly / 2,
+    "patch_icpp(1)%length_x": Lx,
+    "patch_icpp(1)%length_y": Ly,
+    "patch_icpp(1)%vel(1)": 0,
+    "patch_icpp(1)%vel(2)": 0,
+    "patch_icpp(1)%pres": 101325,
+    "patch_icpp(1)%alpha_rho(1)": 1.00,
+    "patch_icpp(1)%alpha(1)": 1,
+    "fluid_pp(1)%gamma": 1.0e00 / (1.4e00 - 1.0e00),
+    "fluid_pp(1)%pi_inf": 0.0e00,
+    "viscous": "T",
+    "fluid_pp(1)%Re(1)": 100000,
+}
+for i in range(len(sol_L.Y)):
+    case[f"chem_wrt_Y({i + 1})"] = "T"
+    case[f"patch_icpp(1)%Y({i + 1})"] = sol_L.Y[i]
+
+if __name__ == "__main__":
+    print(json.dumps(case))

docs/res/simulations/x.jpg

@@ -0,0 +1,114 @@
+import json
+import math
+
+# Free stream
+p_inf = 16281
+rho_inf = 0.416
+c_inf = math.sqrt(1.4 * p_inf / rho_inf)
+vel_inf = 2.5 * c_inf
+
+# Jet inlet conditions
+p_J = 1.0 * p_inf
+vel_J = 41.0
+rho_J = 1000.0
+
+# Viscosities
+mu_inf = 1.85e-5
+mu_J = 2.67e-3
+
+d_jet = 1.2e-3
+
+# Simulation parameters
+Ny = 3999
+Nx = 11999
+time_end = 5e-4
+cfl = 0.6
+eps = 1e-8
+
+# Case dictionary
+
+print(
+    json.dumps(
+        {
+            # Logistics
+            "run_time_info": "T",
+            # Computational domain parameters
+            "x_domain%beg": -40 * d_jet,
+            "x_domain%end": 60 * d_jet,
+            "y_domain%beg": 0.0,
+            "y_domain%end": 35 * d_jet,
+            "m": int(Nx),
+            "n": int(Ny),
+            "p": 0,
+            "cfl_adap_dt": "T",
+            "t_stop": time_end,
+            "t_save": time_end / 50,
+            "n_start": 0,
+            "cfl_target": cfl,
+            # Simulation algorithm parameters
+            "num_patches": 2,
+            "model_eqns": 3,
+            "alt_soundspeed": "F",
+            "num_fluids": 2,
+            "mpp_lim": "T",
+            "mixture_err": "F",
+            "time_stepper": 3,
+            "recon_type": 2,
+            "muscl_order": 2,
+            "muscl_lim": 1,
+            "riemann_solver": 2,
+            "wave_speeds": 1,
+            "avg_state": 2,
+            "elliptic_smoothing": "T",
+            "elliptic_smoothing_iters": 10,  # 50,
+            "bc_x%beg": -17,
+            "bc_x%end": -12,
+            "bc_y%beg": -17,
+            "bc_y%end": -3,
+            "num_bc_patches": 0,
+            # Formatted Database File Structures
+            "format": 1,
+            "precision": 2,
+            "prim_vars_wrt": "T",
+            "parallel_io": "T",
+            # Patch 1: Free stream
+            "patch_icpp(1)%geometry": 3,
+            "patch_icpp(1)%x_centroid": 0.0,
+            "patch_icpp(1)%y_centroid": 50 * d_jet,
+            "patch_icpp(1)%length_x": 1000 * d_jet,
+            "patch_icpp(1)%length_y": 1000 * d_jet,
+            "patch_icpp(1)%vel(1)": f"{vel_inf} * tanh(y / {d_jet / 4})",
+            "patch_icpp(1)%vel(2)": 0.0,
+            "patch_icpp(1)%pres": p_inf,
+            "patch_icpp(1)%alpha_rho(1)": rho_inf * (1.0 - eps),
+            "patch_icpp(1)%alpha(1)": 1.0 - eps,
+            "patch_icpp(1)%alpha_rho(2)": eps,
+            "patch_icpp(1)%alpha(2)": eps,
+            # Patch 2: Jet
+            "patch_icpp(2)%geometry": 3,
+            "patch_icpp(2)%alter_patch(1)": "T",
+            "patch_icpp(2)%x_centroid": 0.0,
+            "patch_icpp(2)%y_centroid": 0.0,
+            "patch_icpp(2)%length_x": d_jet,
+            "patch_icpp(2)%length_y": d_jet,
+            "patch_icpp(2)%vel(1)": 0.0,
+            "patch_icpp(2)%vel(2)": vel_J,
+            "patch_icpp(2)%pres": p_J,
+            "patch_icpp(2)%alpha_rho(1)": eps,
+            "patch_icpp(2)%alpha(1)": eps,
+            "patch_icpp(2)%alpha_rho(2)": (1.0 - eps) * rho_J,
+            "patch_icpp(2)%alpha(2)": 1.0 - eps,
+            # Fluid properties
+            "fluid_pp(1)%gamma": 1.00 / (1.4 - 1.0),
+            "fluid_pp(1)%pi_inf": 0.0,
+            "fluid_pp(2)%gamma": 1.0 / (6.12 - 1.0),
+            "fluid_pp(2)%pi_inf": 6.12 * 3.43e8 / (6.12 - 1),
+            "viscous": "T",
+            "fluid_pp(1)%Re(1)": 1 / mu_inf,
+            "fluid_pp(2)%Re(1)": 1 / mu_J,
+            "surface_tension": "T",
+            "cf_wrt": "T",
+            "sigma": 0.072,
+        }
+    )
+)