Add singly degenerate Cahn-Hilliard mobility

.github/workflows/linux.yml

@@ -234,8 +234,12 @@ jobs:
         echo "## 📄 Report: test-2d-24-04-clang" >> $GITHUB_STEP_SUMMARY
         echo "" >> $GITHUB_STEP_SUMMARY
         echo "[⬇️ Click here to download](${{ steps.upload_report.outputs.artifact-url }})" >> $GITHUB_STEP_SUMMARY
-  test-2d-24-04-clang-fft:
+  test-24-04-clang-fft:
+    name: test-${{ matrix.dim }}d-24-04-clang-fft
     runs-on: ubuntu-24.04
+    strategy:
+      matrix:
+        dim: [2, 3]
     steps:
     - uses: actions/checkout@v4
     - name: installing preliminaries
@@ -244,19 +248,19 @@ jobs:
       run: sudo apt install clang
     - name: installing fftw3
       run: sudo apt install libfftw3-dev
-    - name: configure (2d)
-      run: ./configure --dim=2 --comp=clang++ --fft
+    - name: configure (${{ matrix.dim }}d)
+      run: ./configure --dim=${{ matrix.dim }} --comp=clang++ --fft
     - name: make
       run: make
     - name: regression tests
-      run: scripts/runtests.py --serial --no-coverage --comp=clang++  --timeout=200 --no-backspace --fft-only
+      run: scripts/runtests.py --serial --no-coverage --dim=${{ matrix.dim }} --comp=clang++  --timeout=200 --no-backspace --fft-only
     - name: create report
       if: always()
-      run: python3 scripts/package_report.py --name test-2d-24-04-clang-fft
+      run: python3 scripts/package_report.py --name test-${{ matrix.dim }}d-24-04-clang-fft
     - name: report artifact summary
       if: always()
       run: |
-        name="test-2d-24-04-clang-fft"
+        name="test-${{ matrix.dim }}d-24-04-clang-fft"
         echo "### Report bundle: $name" >> $GITHUB_STEP_SUMMARY
         if [ -d "$name" ]; then
           echo "Files: $(find "$name" -type f | wc -l)" >> $GITHUB_STEP_SUMMARY
@@ -270,11 +274,11 @@ jobs:
       id: upload_report
       uses: actions/upload-artifact@v4
       with:
-        name: test-2d-24-04-clang-fft
-        path: test-2d-24-04-clang-fft
+        name: test-${{ matrix.dim }}d-24-04-clang-fft
+        path: test-${{ matrix.dim }}d-24-04-clang-fft
     - name: Add report link to summary
       if: always()
       run: |
-        echo "## 📄 Report: test-2d-24-04-clang-fft" >> $GITHUB_STEP_SUMMARY
+        echo "## 📄 Report: test-${{ matrix.dim }}d-24-04-clang-fft" >> $GITHUB_STEP_SUMMARY
         echo "" >> $GITHUB_STEP_SUMMARY
         echo "[⬇️ Click here to download](${{ steps.upload_report.outputs.artifact-url }})" >> $GITHUB_STEP_SUMMARY

scripts/runtests.py

@@ -145,6 +145,7 @@ def readMetadata(metadatafile):
 parser.add_argument('--mpirun-flags',dest="mpirun_flags",default="",help="Extra arguments to pass to mpirun (like --oversubscribe). All arguments must be in a string.")
 parser.add_argument('--fft',dest="fft",default=False,action='store_true',help="Enable fft-based tests")
 parser.add_argument('--fft-only',dest="fft_only",default=False,action='store_true',help="Run fft tests only")
+parser.add_argument('--hdf5', default=False, action='store_true', help='Use HDF5-enabled executables')
 parser.add_argument('--cuda', const='auto', default=None, nargs='?', type=str,
                     help="Run CUDA tests, optionally specifying the executable architecture suffix")
 parser.add_argument('--post-timeout', dest="post_timeout", default=10000, help='How long to wait before skipping results posting')
@@ -445,6 +446,7 @@ def test(testdir):
             if args.profile: exestr += "-profile"
             if args.perf: exestr += "-perf"
             if coverage: exestr += "-coverage"
+            if args.hdf5: exestr += "-hdf5"
             if args.cuda:
                 cuda_arch = args.cuda
                 if cuda_arch == "auto":

src/Integrator/CahnHilliard.H

@@ -14,9 +14,36 @@
 //
 //    \frac{\partial\eta}{\partial t} = L\nabla^2\Big(\eta^3 - \eta - \gamma\nabla^2\eta\Big)
 //
-// which is integrated using a forward Euler scheme.
+// Constant mobility is integrated using forward Euler in realspace or a
+// semi-implicit spectral scheme. With :code:`mobility = singly_degenerate`,
+// the equation is instead
 //
-// This is tested in :ref:`CahnHilliard`
+// .. math::
+//
+//    \partial_t\eta = \nabla\cdot(M(\eta)\nabla\mu),\qquad
+//    M(\eta) = L\alpha^2(1-\alpha)^2 + M_0,\qquad
+//    \alpha = \mathrm{clip}((\eta+1)/2,0,1),
+//
+// where :math:`\mu=\eta^3-\eta-\gamma\nabla^2\eta` and
+// :math:`M_0` is :code:`mobility_floor` (default zero).
+// Both SDCH methods use conservative face fluxes and the same discrete Laplacian.
+// The spectral method applies the Fourier multiplier :math:`1/(1+\Delta t S D^2)`
+// to the increment, where :math:`D` is the negative discrete Laplacian symbol
+// and :math:`S` is :code:`spectral_stabilization`, default
+// :math:`\gamma(L/16+M_0)`. SDCH requires periodic boundaries. Its spectral
+// method assembles a uniform finest-level grid and restricts the update to
+// the AMR hierarchy; realspace AMR has no flux refluxing at coarse/fine interfaces.
+// By default the solution is not clipped, so timesteps must resolve the
+// interface dynamics. Clipping can change the conserved mean.
+// Optional input and field names are :code:`input_name`, :code:`field_name`,
+// :code:`old_field_name`, and :code:`intermediate_name`. SDCH can optionally
+// clip its updated field with :code:`clipping` (default false). Derived integrators can override
+// :code:`FillMobilityMask` to apply a spatial multiplier to the face mobility.
+// The hook receives each level's native box layout and valid cells; the caller
+// fills periodic ghosts. Set :code:`diagnostics = true` for SDCH mobility, mask, flux, RHS, raw update,
+// and update delta fields alongside the solution.
+//
+// This is tested in :ref:`CahnHilliard` and :ref:`SinglyDegenerateCahnHilliard`
 
 #ifndef INTEGRATOR_CAHNHILLIARD_H
 #define INTEGRATOR_CAHNHILLIARD_H
@@ -32,7 +59,7 @@
 namespace Integrator
 {
 
-class CahnHilliard : public Integrator
+class CahnHilliard : virtual public Integrator
 {
 public:
     static constexpr const char* name = "cahnhilliard";
@@ -60,23 +87,49 @@ protected:
     void AdvanceSpectral (int lev, Set::Scalar time, Set::Scalar dt);
     /// realspace version of advance
     void AdvanceReal (int lev, Set::Scalar time, Set::Scalar dt);
+    /// Conservative realspace update with variable mobility
+    void AdvanceDegenerateReal (int lev, Set::Scalar time, Set::Scalar dt);
+    /// Stabilized spectral update with variable mobility
+    void AdvanceDegenerateSpectral (int lev, Set::Scalar time, Set::Scalar dt);
     /// Mark any cells that need to be refined
     void TagCellsForRefinement (int lev, amrex::TagBoxArray& tags, amrex::Real time, int ngrow) override;
 
-private:
+    /// Fill a cell-centered multiplier for the singly-degenerate mobility.
+    /// The default is one everywhere; derived integrators may override it.
+    virtual void FillMobilityMask(int lev, amrex::MultiFab& mask);
+    void InitializeDiagnostics(int lev);
+
+protected:
 
     Set::Field<Set::Scalar> etanew_mf;     /// The new value for eta this timestep
     Set::Field<Set::Scalar> etaold_mf;     /// Last timestep's value for eta
     Set::Field<Set::Scalar> intermediate;  /// Intermediate field used for CH kinetics
+    Set::Field<Set::Scalar> diagnostic_mobility;
+    Set::Field<Set::Scalar> diagnostic_mobility_mask;
+    Set::Field<Set::Scalar> diagnostic_flux;
+    Set::Field<Set::Scalar> diagnostic_rhs;
+    Set::Field<Set::Scalar> diagnostic_raw;
+    Set::Field<Set::Scalar> diagnostic_delta;
 
-    BC::BC<Set::Scalar> *bc; /// eta's bc object
-    IC::IC<Set::Scalar> *ic; /// eta's ic object
+    BC::BC<Set::Scalar> *bc = nullptr; /// eta's bc object
+    IC::IC<Set::Scalar> *ic = nullptr; /// eta's ic object
 
     Set::Scalar gamma = NAN; 
     Set::Scalar L     = NAN; 
     Set::Scalar refinement_threshold = NAN; 
 
-    std::string method; 
+    std::string method;
+    std::string mobility = "constant";
+    Set::Scalar mobility_floor = 0.0;
+    Set::Scalar spectral_stabilization = NAN;
+    Set::Scalar tstart = 1E100;
+    Set::Scalar tstart_timestep = -1.0;
+    bool clipping = false;
+    bool output_diagnostics = false;
+    std::string input_name = "eta";
+    std::string field_name = "eta";
+    std::string old_field_name = "eta_old";
+    std::string intermediate_name = "int";
 };
 }
 #endif