Add test cases for subgrid scale wetting and drying corrections

[sbrus89]

This pull request introduces two test cases for the subgrid scale correction scheme added in this E3SM PR: E3SM-Project/E3SM#6288. These test cases were developed under the ICoM project in collaboration with Dam Wirasaet at U. of Notre Dame.

The subgrid corrections used in these tests account for the nonlinear relationship between ssh and grid-averaged wet volume per unit area for partially wet cells. These relationships are calculated as lookup tables from high-resolution, subgrid scale, bathymetry data. This pre-processing step is implemented in MPAS-Ocean init mode 🙈. Using the subgrid correction scheme increases the accuracy of wetting and drying, and is especially advantageous at coarser resolutions. Currently, this capability is only implemented for single layer configurations.

The two test cases are:

1. A subgrid version of the existing parabolic bowl test case

2. A simulation of Buttermilk Bay, MA with idealized tidal boundary forcing

More information on the subgrid correction methodology can be found in Kennedy et al. 2019. This paper also discusses the Buttermilk Bay test case.

This PR also exercises a new (optional) update to the existing wetting and drying algorithm that takes into account the ssh gradient between cells that share a previously limited edge. This treatment of limited edges improves the wetting and drying convergence for the parabolic bowl. However, this approach is currently limited to a single layer configuration. Thanks to Dam for contributing this improvement!

Checklist

• User's Guide has been updated
• Developer's Guide has been updated
• API documentation in the Developer's Guide (api.rst) has any new or modified class, method and/or functions listed
• Documentation has been built locally and changes look as expected
• The E3SM-Project submodule has been updated with relevant E3SM changes
• Document (in a comment titled Testing in this PR) any testing that was used to verify the changes

[sbrus89]

Testing

These are the results of the parabolic bowl subgrid case:

The standard wetting and drying results for the parabolic bowl:

The convergence comparison between the subgrid and standard cases:

For comparison, the convergence plot without the limited edge ssh treatment is:

[sbrus89]

Testing

The results for the subgrid Buttermilk Bay case are:

The standard wetting and drying Buttermilk Bay results are:

[mark-petersen]

Thanks! Does this file need to be added?

requests.exceptions.HTTPError: 404 Client Error: Not Found for url: https://web.lcrc.anl.gov/public/e3sm/mpas_standalonedata/mpas-ocean/initial_condition_database/buttermilk_bathy.nc

[xylar]

@sbrus89, bathymetry would ideally go in the bathymetry database, not in the initial condition database:
https://web.lcrc.anl.gov/public/e3sm/mpas_standalonedata/mpas-ocean/bathymetry_database/
Picky, I know...

[xylar]

And sorry for not reviewing. I was waiting on the E3SM PR.

[xylar]

Thanks!

[sbrus89]

@mark-petersen, the bathymetry should be available now.

[mark-petersen]

Passes all tests using intel on chrysalis and gnu on perlmutter, both debug and optimized. Tests are as follows:

perlmutter

./conda/configure_compass_env.py --conda ~/miconda3 --compiler gnu -m pm-cpu
source load_dev_compass_1.2.0-alpha.9_pm-cpu_gnu_mpich.sh
compass list|grep -i parabolic_bowl
 370: ocean/parabolic_bowl/standard/ramp_lts
 371: ocean/parabolic_bowl/standard/ramp
 372: ocean/parabolic_bowl/standard/noramp_lts
 373: ocean/parabolic_bowl/standard/noramp
 374: ocean/parabolic_bowl/subgrid/ramp
 375: ocean/parabolic_bowl/subgrid/noramp

compass list|grep -i butter
 107: ocean/buttermilk_bay/standard
 108: ocean/buttermilk_bay/subgrid

compass setup -p ~/repos/E3SM/pr/components/mpas-ocean -w $r/230315_subgrid_sbrus_opt -n 107 108 374 375

salloc --nodes 1 --qos interactive --time 01:00:00 --constraint cpu --account=e3sm
source load_dev_compass_1.2.0-alpha.9_pm-cpu_gnu_mpich.sh
cd  $r/230315_subgrid_sbrus_opt
compass run

chyrsalis

./conda/configure_compass_env.py --conda ~/miconda3 --compiler intel
source load_dev_compass_1.2.0-alpha.9_chrysalis_intel_openmpi.sh
compass setup -p ~/repos/E3SM/pr/components/mpas-ocean -w $r/230315_subgrid_sbrus_opt -n 107 108 374 375

srun -p debug -N 1 -t 1:00:00 --pty bash
source load_dev_compass_1.2.0-alpha.9_pm-cpu_gnu_mpich.sh
cd  $r/230315_subgrid_sbrus_opt
compass run

[mark-petersen]

@sbrus89 this is excellent! Approving based on documentation, output plots, and passing functional testing with two compilers. I did not review the details of the compass code implementation.

[sbrus89]

Thanks for testing, @mark-petersen!

[cbegeman]

This looks fantastic! I ran the standard test cases successfully. Just a few minor suggestions.

[xylar]

#821 will update the E3SM-Project submodule with the changes needed to support this PR. @sbrus89, could you rebase this branch once that goes in?

[sbrus89]

@xylar, thanks for the submodule update. I just rebased.

[xylar]

@sbrus89, I was able to run the Buttermilk Bay test just fine!

I had trouble with the viz step in both parabolic bowl test cases. The ramp test produces:

  in /gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/ocean/tests/parabolic_bowl/viz/__init__.py

      Failed
Exception raised while running the steps of the test case
Traceback (most recent call last):
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/run/serial.py", line 320, in _log_and_run_test
    _run_test(test_case, available_resources)
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/run/serial.py", line 417, in _run_test
    _run_step(test_case, step, test_case.new_step_log_file,
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/run/serial.py", line 468, in _run_step
    step.run()
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/ocean/tests/parabolic_bowl/viz/__init__.py", line 50, in run
    self.rmse_plots()
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/ocean/tests/parabolic_bowl/viz/__init__.py", line 227, in rmse_plots
    rmse[i] = self.compute_rmse(
              ^^^^^^^^^^^^^^^^^^
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/ocean/tests/parabolic_bowl/viz/__init__.py", line 262, in compute_rmse
    ind = time.index(dt.datetime.strptime('0001-01-03_18:00:00',
          ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
ValueError: datetime.datetime(1, 1, 3, 18, 0) is not in list

The noramp test leads to:

Exception raised while running the steps of the test case
Traceback (most recent call last):
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/run/serial.py", line 320, in _log_and_run_test
    _run_test(test_case, available_resources)
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/run/serial.py", line 417, in _run_test
    _run_step(test_case, step, test_case.new_step_log_file,
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/run/serial.py", line 468, in _run_step
    step.run()
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/ocean/tests/parabolic_bowl/viz/__init__.py", line 49, in run
    self.contour_plots(points)
  File "/gpfs/fs1/home/ac.xylar/compass/subgrid_pr/compass/ocean/tests/parabolic_bowl/viz/__init__.py", line 164, in contour_plots
    ds['ssh'][i, :],
    ~~~~~~~~~^^^^^^
  File "/gpfs/fs1/home/ac.xylar/chrysalis/miniforge3/envs/compass_test2/lib/python3.11/site-packages/xarray/core/dataarray.py", line 874, in __getitem__
    return self.isel(indexers=self._item_key_to_dict(key))
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/gpfs/fs1/home/ac.xylar/chrysalis/miniforge3/envs/compass_test2/lib/python3.11/site-packages/xarray/core/dataarray.py", line 1493, in isel
    variable = self._variable.isel(indexers, missing_dims=missing_dims)
               ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/gpfs/fs1/home/ac.xylar/chrysalis/miniforge3/envs/compass_test2/lib/python3.11/site-packages/xarray/core/variable.py", line 1018, in isel
    return self[key]
           ~~~~^^^^^
  File "/gpfs/fs1/home/ac.xylar/chrysalis/miniforge3/envs/compass_test2/lib/python3.11/site-packages/xarray/core/variable.py", line 785, in __getitem__
    data = indexing.apply_indexer(indexable, indexer)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/gpfs/fs1/home/ac.xylar/chrysalis/miniforge3/envs/compass_test2/lib/python3.11/site-packages/xarray/core/indexing.py", line 1002, in apply_indexer
    return indexable[indexer]
           ~~~~~~~~~^^^^^^^^^
  File "/gpfs/fs1/home/ac.xylar/chrysalis/miniforge3/envs/compass_test2/lib/python3.11/site-packages/xarray/core/indexing.py", line 820, in __getitem__
    return type(self)(_wrap_numpy_scalars(self.array[indexer]))
                                          ~~~~~~~~~~^^^^^^^^^
  File "/gpfs/fs1/home/ac.xylar/chrysalis/miniforge3/envs/compass_test2/lib/python3.11/site-packages/xarray/core/indexing.py", line 770, in __getitem__
    return type(self)(_wrap_numpy_scalars(self.array[indexer]))
                                          ~~~~~~~~~~^^^^^^^^^
  File "/gpfs/fs1/home/ac.xylar/chrysalis/miniforge3/envs/compass_test2/lib/python3.11/site-packages/xarray/core/indexing.py", line 645, in __getitem__
    return type(self)(self.array, self._updated_key(indexer))
                                  ^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/gpfs/fs1/home/ac.xylar/chrysalis/miniforge3/envs/compass_test2/lib/python3.11/site-packages/xarray/core/indexing.py", line 600, in _updated_key
    full_key.append(_index_indexer_1d(k, next(iter_new_key), size))
                    ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/gpfs/fs1/home/ac.xylar/chrysalis/miniforge3/envs/compass_test2/lib/python3.11/site-packages/xarray/core/indexing.py", line 287, in _index_indexer_1d
    indexer = _expand_slice(old_indexer, size)[applied_indexer]  # type: ignore[assignment]
              ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~^^^^^^^^^^^^^^^^^
IndexError: index 80 is out of bounds for axis 0 with size 79

[xylar]

Could you see if you are able to reproduce those errors? (This was on Chrysalis with the default Intel and OpenMPI configuration.)

[sbrus89]

@xylar, thanks for testing. I'm looking into it now.

[sbrus89]

@xylar, both the ocean/parabolic_bowl/subgrid/ramp and ocean/parabolic_bowl/subgrid/noramp seemed to pass for me on Chrysalis with the Intel/OpenMPI config. Can you point me to your work directory so I can look more closely at what the issue might be?

[xylar]

@sbrus89, I'll try running again first and point to the errors if they happen again.

[xylar]

@sbrus89, I reran and didn't see the errors from before. I must have messed things up the first time. My job ran out of time and maybe stepped on its own toes somehow.

[sbrus89]

Great, thank you @xylar!

[sbrus89]

@cbegeman, I forgot I made this change in this PR related to #765 .

[cbegeman]

@sbrus89 Great. I just closed the issue.