Modularize the entropy viscosity regularization

CHANGELOG.md

@@ -1,6 +1,11 @@
 # Changelog
 
 ## Develop
+- *BREAKING* Modularized entropy viscosity in the compressible Navier-Stokes
+  solver. The obsolete `case.numerics.c_avisc_low` and
+  `case.numerics.c_avisc_entropy` options are now rejected. Configure an
+  `artificial_viscosity_model` simulation component with those coefficients
+  and enable its field through the fluid `viscous_regularization` object.
 
 - *BREAKING*, normal_outflow conditions now require specifying `value`, which
   is used to set the value of the tangential components of velocity.
@@ -153,6 +158,9 @@
 - Added `host_array_t` and `device_array_t` temporary array types and support
   for requesting these through `scratch_registry_t`.
 - Added the `cai_sagaut_model_ii` wall model with CPU, CUDA, HIP, and OpenCL.
+- Modularized `entropy_viscosity` into two steps: 1. Calculate the
+  artificial viscosity field in a simcomp; 2. Apply the artificial 
+  viscosity though a viscous regularization object.
 - Added the `create_periodic_zones` contrib utility for converting pairs of
   labeled zones in an existing `.nmsh` mesh into periodic zones.
   backends. This model is based on the work of [Cai and Sagaut (PoF,

doc/pages/user-guide/case-file.md

@@ -373,13 +373,6 @@ The compressible solver accepts the following parameters:
 | `mu`    | Constant physical dynamic viscosity              | Non-negative real | `0.0`         |
 | `kappa` | Constant physical thermal conductivity           | Non-negative real | `0.0`         |
 
-Additional numerics parameters specific to compressible flows:
-
-| Name              | Description                                        | Admissible values | Default value |
-| ----------------- | -------------------------------------------------- | ----------------- | ------------- |
-| `c_avisc_low`     | Coefficient for low-order artificial viscosity     | Positive reals    | `0.5`         |
-| `c_avisc_entropy` | Coefficient for entropy-based artificial viscosity | Positive reals    | `1.0`         |
-
 The compressible solver uses variable time-stepping controlled by the CFL
 number. Set `variable_timestep` to `true` and specify `target_cfl` in the time
 control object.
@@ -403,6 +396,10 @@ Example configuration:
   "fluid": {
     "scheme": "compressible",
     "gamma": 1.4,
+    "viscous_regularization":{
+        "type": "artificial_viscosity",
+        "reg_coeff_name": "entropy_viscosity"
+    },
     "initial_condition": {
       "type": "user"
     },
@@ -426,12 +423,13 @@ Example configuration:
     "output_control": "nsamples",
     "output_value": 20
   },
-  "numerics": {
-    "time_order": 3,
-    "polynomial_order": 5,
-    "c_avisc_low": 0.5,
-    "c_avisc_entropy": 0.5
-  }
+  "simulation_components": [
+    {
+      "type": "artificial_viscosity_model",
+      "model": "entropy_viscosity",
+      "field_name": "entropy_viscosity"
+    }
+  ]
 }
 ~~~~~~~~~~~~~~~
 
@@ -2020,7 +2018,6 @@ concisely directly in the table.
 | `shear_stress.value`                               | The shear stress vector value for `sh` boundaries                                                 | Vector of 3 reals                                           | `[0, 0, 0]`   |
 | `wall_modelling.type`                              | The wall model type for `wm` boundaries. See documentation for additional config parameters.      | `rough_log_law`, `spalding`                                 | -             |
 | `source_terms`                                     | Array of JSON objects, defining additional source terms.                                          | See list of source terms above                              | -             |
-| `gradient_jump_penalty`                            | Array of JSON objects, defining additional gradient jump penalty.                                 | See list of gradient jump penalty above                     | -             |
 | `boundary_types`                                   | Boundary types/conditions labels.                                                                 | Array of strings                                            | -             |
 | `velocity_solver.type`                             | Linear solver for the momentum equation.                                                          | `cg`, `pipecg`, `bicgstab`, `cacg`, `gmres`                 | -             |
 | `velocity_solver.preconditioner.type`              | Linear solver preconditioner for the momentum equation.                                           | `ident`, `hsmg`, `jacobi`                                   | -             |
@@ -2237,3 +2234,43 @@ currently supports 50 regions, with id 1..25 being reserved for internal use.
 | ---------------- | ----------------------------------------------------------- | ----------------- | ------------- |
 | `enabled`        | Whether to enable gathering of runtime statistics           | `true` or `false` | `false`       |
 | `output_profile` | Whether to output all gathered profiling data as a CSV file | `true` or `false` | `false`       |
+
+## Viscous regularization {#case-file_viscous-regularization}
+
+Users can use the viscous regularization object to enhance the smoothness or 
+the numerical stability of the solution via a diffusion term. For example, 
+artificial viscosity can be set up by this object to perform shock capturing.
+Note that some regularization techniques do not have a diffusive mathematical 
+form, for example gradient jump penalty and high-pass filter relaxation terms,
+and they are included in the source terms instead. The viscous regularization
+can be set up with the following options:
+
+* `type`, the viscous regularization type.
+  - `artificial_viscosity`, the standard second-order diffusion term,
+    $\frac{\partial}{\partial x}\left(\mu_\mathrm{artificial}\frac{\partial u}{\partial x}\right)$.
+* `reg_coeff_name`, name of the $\mu_\mathrm{artificial}$ field, usually 
+computed by a simulation component.
+
+Viscous regularization is currently supported by the compressible fluid
+solver. Artificial viscosity requires both the regularization object that
+consumes the coefficient and a simulation component that computes it. The
+`reg_coeff_name` and `field_name` values must match:
+
+~~~~~~~~~~~~~~~{.json}
+{
+  "fluid": {
+    "scheme": "compressible",
+    "viscous_regularization": {
+      "type": "artificial_viscosity",
+      "reg_coeff_name": "entropy_viscosity"
+    }
+  },
+  "simulation_components": [
+    {
+      "type": "artificial_viscosity_model",
+      "model": "entropy_viscosity",
+      "field_name": "entropy_viscosity"
+    }
+  ]
+}
+~~~~~~~~~~~~~~~

doc/pages/user-guide/simcomps.md

@@ -65,6 +65,8 @@ in Neko. The list will be updated as new simcomps are added.
 - Streaming of data for in-situ field manipulation \ref simcomp_data_streamer
 - Sub-sampling of fields by changing polynomial order and masking by point
   zones \ref simcomp_field_subsampler
+- Computation of the viscous coefficient in artificial viscosity method (AVM)
+  via a particular model \ref simcomp_avm
 
 ## Controlling execution and file output
 Each simulation component is, by default, executed once per time step to perform
@@ -764,15 +766,15 @@ keywords:
     the maximum spacing of the GLL nodes within the element. The `les_delta`
   field is added to the registry and written to the .fld files.
 - `nut_field`: The name of the SGS eddy viscosity field added to the registry.
-  Defaults to `nut`. This allows to have two different SGS models active, saved
-  to different fields. For example, one for the scalar and one to the fluid.
-- `extrapolation`: Whether or not extrapolate the velocity to
+  Defaults to `nut`. This allows two different SGS models to be active and saved
+  to different fields, for example one for the scalar and one for the fluid.
+- `extrapolation`: Whether or not to extrapolate the velocity to
   compute the eddy viscosity.
-  - `true`: extrapolate the velocity as the same order as
+  - `true`: extrapolate the velocity to the same order as
   the time scheme.
-  - `false`: the default option, disable the extrapolation.
+  - `false`: the default option, disables extrapolation.
   In this case, the estimation of the eddy viscosity is of first order, while
-  circumvent the risk of unstable extrapolation.
+  circumventing the risk of unstable extrapolation.
 
  ~~~~~~~~~~~~~~~{.json}
  {
@@ -932,4 +934,34 @@ keywords used by the latter can also be specified, with the exception of
    "compute_control": "tsteps",
    "compute_value": 10
  }
- ~~~~~~~~~~~~~~~
+ ~~~~~~~~~~~~~~~
+
+### artificial_viscosity_model {#simcomp_avm}
+Computes an artificial viscosity field using a particular model. **Note:** The simcomp
+*only* computes the artificial viscosity field. You have to select the 
+corresponding `reg_coeff_name` in the compressible fluid solver's
+`viscous_regularization` JSON object to actually enable artificial viscosity;
+see @ref case-file_viscous-regularization. The simcomp is controlled by the
+following keywords:
+
+- `model`: Selects the artificial viscosity model. Currently available models are:
+  - `entropy_viscosity`: An entropy-based viscosity aimed at shock capturing;
+    see Guermond et al. (https://doi.org/10.1016/j.jcp.2010.11.043).
+    - `c_avisc_low`: Coefficient for low-order artificial viscosity, defaults
+    to 0.5
+    - `c_avisc_entropy`: Coefficient for entropy-based artificial viscosity,
+    defaults to 1.0
+
+- `field_name`: The name of the artificial viscosity field added to the registry.
+  Defaults to the model name.
+  It must match the `reg_coeff_name` selected by the consuming regularization
+  object.
+
+ ~~~~~~~~~~~~~~~{.json}
+ {
+   "type": "artificial_viscosity_model",
+   "name": "artificial_viscosity_model",
+   "model": "entropy_viscosity",
+   "field_name": "entropy_viscosity"
+ }
+ ~~~~~~~~~~~~~~~

doc/schemas/numerics.schema.json

@@ -43,12 +43,6 @@
         },
         "oifs_target_cfl": {
           "$ref": "urn:neko:schema:common#/$defs/positiveNumberOrRef"
-        },
-        "c_avisc_low": {
-          "$ref": "urn:neko:schema:common#/$defs/positiveNumberOrRef"
-        },
-        "c_avisc_entropy": {
-          "$ref": "urn:neko:schema:common#/$defs/positiveNumberOrRef"
         }
       },
       "additionalProperties": false

doc/schemas/simcomps/artificial_viscosity_model.schema.json

@@ -0,0 +1,40 @@
+{
+  "$schema": "https://json-schema.org/draft/2020-12/schema",
+  "$id": "urn:neko:schema:simcomps:artificial_viscosity_model",
+  "$defs": {
+    "simulationComponent": {
+      "allOf": [
+        {
+          "$ref": "urn:neko:schema:simulation-components#/$defs/simcompBase"
+        },
+        {
+          "$ref": "urn:neko:schema:simulation-components#/$defs/simcompFieldWriterOptions"
+        },
+        {
+          "type": "object",
+          "required": [
+            "model"
+          ],
+          "properties": {
+            "type": {
+              "const": "artificial_viscosity_model"
+            },
+            "model": {
+              "const": "entropy_viscosity"
+            },
+            "field_name": {
+              "$ref": "urn:neko:schema:simulation-components#/$defs/nonEmptyString"
+            },
+            "c_avisc_low": {
+              "$ref": "urn:neko:schema:simulation-components#/$defs/numberOrRef"
+            },
+            "c_avisc_entropy": {
+              "$ref": "urn:neko:schema:simulation-components#/$defs/numberOrRef"
+            }
+          },
+          "additionalProperties": true
+        }
+      ]
+    }
+  }
+}

doc/schemas/simulation-components.schema.json

@@ -607,6 +607,9 @@
         {
           "$ref": "urn:neko:schema:simcomps:field_subsampler#/$defs/simulationComponent"
         },
+        {
+          "$ref": "urn:neko:schema:simcomps:artificial_viscosity_model#/$defs/simulationComponent"
+        },
         {
           "$ref": "urn:neko:schema:simcomps:boundary_data_writer#/$defs/simulationComponent"
         },
@@ -616,4 +619,4 @@
       ]
     }
   }
-}
+}

examples/cns_2d_shock_tube/2d_shock_tube.case

@@ -13,13 +13,15 @@
         },
         "numerics": {
             "time_order": 3,
-            "polynomial_order": 3,
-            "c_avisc_low": 0.5,
-            "c_avisc_entropy": 0.5
+            "polynomial_order": 3
         },
         "fluid": {
             "scheme": "compressible",
             "gamma": 1.4,
+            "viscous_regularization": {
+                "type": "artificial_viscosity",
+                "reg_coeff_name": "entropy_viscosity"
+            },
             "Re": 1000.0,
             "Pr": 0.73,
             "initial_condition": {
@@ -39,6 +41,12 @@
             "output_value": 100
         },
         "simulation_components": [
+            {
+                "type": "artificial_viscosity_model",
+                "model": "entropy_viscosity",
+                "c_avisc_low": 0.5,
+                "c_avisc_entropy": 0.5
+            },
             {
                 "type": "field_writer",
                 "fields": ["fluid_rho", "S", "artificial_visc"],

examples/euler_1d_sod/sod.case

@@ -13,13 +13,15 @@
         },
         "numerics": {
             "time_order": 4,
-            "polynomial_order": 3,
-            "c_avisc_low": 0.5,
-            "c_avisc_entropy": 1.0
+            "polynomial_order": 3
         },
         "fluid": {
             "scheme": "compressible",
             "gamma": 1.4,
+            "viscous_regularization":{
+                "type": "artificial_viscosity",
+                "reg_coeff_name": "entropy_viscosity"
+            },
             "initial_condition": {
                 "type": "user"
             },
@@ -33,6 +35,10 @@
             "output_value": 100
         },
         "simulation_components": [
+            {
+                "type": "artificial_viscosity_model",
+                "model": "entropy_viscosity"
+            },
             {
                 "type": "field_writer",
                 "fields": ["fluid_rho", "S", "artificial_visc"], // s2, s3, s4 in Paraview

examples/euler_2d_cylinder/euler_2d_cylinder.case

@@ -14,14 +14,16 @@
 
         "numerics": {
             "time_order": 3,
-            "polynomial_order": 5,
-            "c_avisc_low": 0.5,
-            "c_avisc_entropy": 1.0
+            "polynomial_order": 5
         },
 
         "fluid": {
             "scheme": "compressible",
             "gamma": 1.4,
+            "viscous_regularization":{
+                "type": "artificial_viscosity",
+                "reg_coeff_name": "entropy_viscosity"
+            },
 
             "initial_condition": {
                 "type": "user"
@@ -58,6 +60,10 @@
         },
 
         "simulation_components": [
+            {
+                "type": "artificial_viscosity_model",
+                "model": "entropy_viscosity"
+            },
             {
                 "type": "field_writer",
                 "fields": ["fluid_rho", "S", "artificial_visc"], // s2, s3, s4 in Paraview

examples/euler_2d_forward_facing_step/step.case

@@ -14,14 +14,16 @@
 
         "numerics": {
             "time_order": 3,
-            "polynomial_order": 5,
-            "c_avisc_low": 0.5,
-            "c_avisc_entropy": 0.5
+            "polynomial_order": 5
         },
 
         "fluid": {
             "scheme": "compressible",
             "gamma": 1.4,
+            "viscous_regularization":{
+                "type": "artificial_viscosity",
+                "reg_coeff_name": "entropy_viscosity"
+            },
 
             "initial_condition": {
                 "type": "user"
@@ -54,6 +56,10 @@
         },
 
         "simulation_components": [
+            {
+                "type": "artificial_viscosity_model",
+                "model": "entropy_viscosity"
+            },
             {
                 "type": "field_writer",
                 "fields": ["fluid_rho", "S", "artificial_visc", "max_wave_speed"], // s2, s3, s4, s5 in Paraview