diff --git a/.github/workflows/unit_tests.yaml b/.github/workflows/unit_tests.yaml
index 0e7aa421..46d3ca53 100644
--- a/.github/workflows/unit_tests.yaml
+++ b/.github/workflows/unit_tests.yaml
@@ -19,7 +19,18 @@ jobs:
     strategy:
       matrix:
         python-version: ['3.11', '3.12', '3.13']
-    name: Python ${{ matrix.python-version }}
+        extra: ['test,zarr']
+        marker-serial: ['not parallel and not gpu and not pyfms']
+        marker-parallel: ['parallel and not gpu and not pyfms']
+        include:
+          # add pyfms tests for 3.12
+          - extra: 'test,pyfms'
+            marker-serial: 'pyfms and not parallel'
+            marker-parallel: 'pyfms and parallel'
+            python-version: '3.12'
+      # don't cancel other jobs if one fails
+      fail-fast: false
+    name: Python ${{ matrix.python-version }}${{ contains(matrix.extra, 'pyfms') && ' (pyFMS)' || '' }}
     steps:
       - name: Checkout repository
         uses: actions/checkout@v6
@@ -31,17 +42,24 @@ jobs:
         with:
           python-version: ${{ matrix.python-version }}
 
-      - name: Install mpi (MPICH flavor)
-        run: pip3 install mpich
+      - name: Install MPI
+        if: ${{!contains(matrix.extra, 'pyfms')}}
+        run: pip3 install openmpi
+
+      - name: Install pyFMS dependencies (includes system MPI)
+        if: contains(matrix.extra, 'pyfms')
+        run: |
+          sudo apt-get update
+          sudo apt-get install libopenmpi-dev netcdf-bin libnetcdf-dev libnetcdff-dev nco libyaml-dev diffutils
 
       - name: Install Python packages
-        run: pip3 install .[test,zarr]
+        run: pip3 install .[${{matrix.extra}}]
 
       - name: Run serial cpu tests
-        run: coverage run --rcfile=pyproject.toml -m pytest -m "not parallel and not gpu" tests
+        run: coverage run --rcfile=pyproject.toml -m pytest -m "${{matrix.marker-serial}}" tests
 
       - name: Run parallel cpu tests
-        run: mpiexec -np 6 coverage run --rcfile=pyproject.toml -m mpi4py -m pytest -m "parallel and not gpu" tests
+        run: mpiexec -np 6 --oversubscribe coverage run --rcfile=pyproject.toml -m mpi4py -m pytest -m "${{matrix.marker-parallel}}" tests
 
       - name: Output code coverage
         run: |
diff --git a/docs/docstrings/monitor/diag_manager_monitor.md b/docs/docstrings/monitor/diag_manager_monitor.md
new file mode 100644
index 00000000..36b92e81
--- /dev/null
+++ b/docs/docstrings/monitor/diag_manager_monitor.md
@@ -0,0 +1,3 @@
+# diag_manager_monitor
+
+::: monitor.diag_manager_monitor
diff --git a/mkdocs.yml b/mkdocs.yml
index 7af6e5ac..d9230943 100644
--- a/mkdocs.yml
+++ b/mkdocs.yml
@@ -75,6 +75,7 @@ nav:
           - "subtile_grid_sizer": docstrings/initialization/subtile_grid_sizer.md
       - monitor:
           - "convert": docstrings/monitor/convert.md
+          - "diag_manager_monitor": docstrings/monitor/diag_manager_monitor.md
           - "netcdf_monitor": docstrings/monitor/netcdf_monitor.md
           - "protocol": docstrings/monitor/protocol.md
           - "zarr_monitor": docstrings/monitor/zarr_monitor.md
diff --git a/ndsl/optional_imports.py b/ndsl/optional_imports.py
index 4cc81ea6..af2ec8d0 100644
--- a/ndsl/optional_imports.py
+++ b/ndsl/optional_imports.py
@@ -17,6 +17,11 @@ def __call__(self, *args: Any, **kwargs: dict) -> None:
 except ModuleNotFoundError as err:
     zarr = RaiseWhenAccessed(err)
 
+try:
+    import pyfms
+except ModuleNotFoundError as err:
+    pyfms = RaiseWhenAccessed(err)
+
 try:
     import cupy
 except ImportError:
diff --git a/pyproject.toml b/pyproject.toml
index 8bce607f..1a7ff3dc 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -124,6 +124,8 @@ module = [
 markers = [
   # tests running on a cpu (e.g. with cupy)
   "gpu",
+  # tests relying on the optional pyFMS dependency
+  "pyfms",
   # tests relying on at least two MPI processes
   "parallel",
   # tests relying on the optional zarr dependency
diff --git a/tests/monitor/__init__.py b/tests/monitor/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/tests/test_dm_monitor_cubed.py b/tests/monitor/test_dm_monitor_cubed.py
similarity index 81%
rename from tests/test_dm_monitor_cubed.py
rename to tests/monitor/test_dm_monitor_cubed.py
index a1ec8c59..e1bdb557 100644
--- a/tests/test_dm_monitor_cubed.py
+++ b/tests/monitor/test_dm_monitor_cubed.py
@@ -8,8 +8,8 @@
 import cftime
 import numpy as np
 import pytest
-import xarray as xr
 import yaml
+from netCDF4 import Dataset, num2date
 
 from ndsl import (
     CubedSphereCommunicator,
@@ -21,9 +21,7 @@
 )
 from ndsl.config import Backend
 from ndsl.initialization import SubtileGridSizer
-
-
-pyfms = pytest.importorskip("pyfms")
+from ndsl.optional_imports import pyfms
 
 
 # init fms mpi and set up a simple domain
@@ -52,7 +50,7 @@ def fms_mpp_init():
     return domain_id
 
 
-def _create_input(reduction: str = "none"):
+def _create_input(reduction: str = "none") -> None:
     diag_config = {
         "title": "ndsl_diag_manager_test",
         "base_date": "1 1 1 0 0 0",
@@ -89,9 +87,9 @@ def _create_input(reduction: str = "none"):
 
 
 # Simple test, uses a lat/lon grid and (1, npes) layout
+@pytest.mark.pyfms
 @pytest.mark.parallel
-def test_dm_monitor():
-
+def test_dm_monitor() -> None:
     npes = MPIComm()._comm.Get_size()
     if npes % 6 != 0:
         raise RuntimeError("this test requires npes to be a multiple of 6 to run")
@@ -203,24 +201,28 @@ def test_dm_monitor():
     pe = MPIComm()._comm.Get_rank() + 1
     filename = "diag_manager_cubed_sphere.tile" + str(pe) + ".nc"
     assert Path(filename).exists()
-    ds = xr.open_mfdataset(filename, decode_times=True)
-    assert "var1" in ds
-    np.testing.assert_array_equal(ds["var1"].shape, (ntimesteps, ny, nx))
-    assert "var2" in ds
-    np.testing.assert_array_equal(ds["var2"].shape, (ntimesteps, nz, ny, nx))
-    assert ds["var1"].dims == ("time", "y", "x")
-    assert ds["var2"].dims == ("time", "z", "y", "x")
-    assert ds["time"].shape == (ntimesteps,)
-    assert ds["time"].dims == ("time",)
-    assert ds["time"].values[0] == cftime.DatetimeNoLeap(1, 1, 1, 0, 0, 15)
-    assert ds["time"].values[1] == cftime.DatetimeNoLeap(1, 1, 1, 0, 0, 30)
-    assert ds["time"].values[2] == cftime.DatetimeNoLeap(1, 1, 1, 0, 0, 45)
+    ds = Dataset(filename)
+    assert "var1" in ds.variables
+    var1_ds = ds.variables["var1"]
+    np.testing.assert_array_equal(var1_ds.shape, (ntimesteps, ny, nx))
+    assert "var2" in ds.variables
+    var2_ds = ds.variables["var2"]
+    np.testing.assert_array_equal(var2_ds.shape, (ntimesteps, nz, ny, nx))
+    assert var1_ds.dimensions == ("time", "y", "x")
+    assert var2_ds.dimensions == ("time", "z", "y", "x")
+    time_var = ds["time"]
+    dates = num2date(time_var[:], units=time_var.units, calendar=time_var.calendar)
+    assert time_var.shape == (ntimesteps,)
+    assert time_var.dimensions == ("time",)
+    assert dates[0] == cftime.DatetimeNoLeap(1, 1, 1, 0, 0, 15)
+    assert dates[1] == cftime.DatetimeNoLeap(1, 1, 1, 0, 0, 30)
+    assert dates[2] == cftime.DatetimeNoLeap(1, 1, 1, 0, 0, 45)
     # data is just the timestep number
-    np.testing.assert_array_equal(ds["var1"].values[0, :, :], 1)
-    np.testing.assert_array_equal(ds["var1"].values[1, :, :], 2)
-    np.testing.assert_array_equal(ds["var1"].values[2, :, :], 3)
-    np.testing.assert_array_equal(ds["var2"].values[0, :, :, :], 2)
-    np.testing.assert_array_equal(ds["var2"].values[1, :, :, :], 4)
-    np.testing.assert_array_equal(ds["var2"].values[2, :, :, :], 6)
+    np.testing.assert_array_equal(var1_ds[0, :, :], 1)
+    np.testing.assert_array_equal(var1_ds[1, :, :], 2)
+    np.testing.assert_array_equal(var1_ds[2, :, :], 3)
+    np.testing.assert_array_equal(var2_ds[0, :, :, :], 2)
+    np.testing.assert_array_equal(var2_ds[1, :, :, :], 4)
+    np.testing.assert_array_equal(var2_ds[2, :, :, :], 6)
 
     pyfms.fms.end()
diff --git a/tests/test_dm_monitor_single.py b/tests/monitor/test_dm_monitor_single.py
similarity index 81%
rename from tests/test_dm_monitor_single.py
rename to tests/monitor/test_dm_monitor_single.py
index f8c86328..bc9e4a55 100644
--- a/tests/test_dm_monitor_single.py
+++ b/tests/monitor/test_dm_monitor_single.py
@@ -9,8 +9,8 @@
 import cftime
 import numpy as np
 import pytest
-import xarray as xr
 import yaml
+from netCDF4 import Dataset, num2date
 
 from ndsl import (
     DiagManagerMonitor,
@@ -22,12 +22,10 @@
 )
 from ndsl.config import Backend
 from ndsl.initialization import SubtileGridSizer
+from ndsl.optional_imports import pyfms
 
 
-pyfms = pytest.importorskip("pyfms")
-
-
-def _create_input(reduction: str = "none"):
+def _create_input() -> None:
     diag_config = {
         "title": "ndsl_diag_manager_test",
         "base_date": "2 1 1 1 1 1",
@@ -63,10 +61,10 @@ def _create_input(reduction: str = "none"):
         f.write(text_content)
 
 
-def test_dm_monitor_single_tile():
+@pytest.mark.pyfms
+def test_dm_monitor_single_tile() -> None:
     # mpi info
     npes = MPIComm()._comm.Get_size()
-    pe = MPIComm()._comm.Get_rank()
     # tile parameters for quantities/domains
     nx = 8
     ny = 8
@@ -220,30 +218,30 @@ def test_dm_monitor_single_tile():
 
     # check output!
     assert Path("diag_manager_single_tile.nc").exists()
-    ds = xr.open_mfdataset("diag_manager_single_tile.nc", decode_times=True)
-    assert "var_2d" in ds
-    np.testing.assert_array_equal(ds["var_2d"].shape, (ntimesteps, nx, ny))
-    assert ds["var_2d"].dims == ("time", "y", "x")
-    assert ds["var_2d"].attrs["units"] == "muntin"
-    assert ds["var_3d"].dims == ("time", "z", "y", "x")
-    assert ds["var_3d"].attrs["units"] == "muntin"
-    assert ds["time"].shape == (ntimesteps,)
-    assert ds["time"].dims == ("time",)
-    assert ds["time"].values[0] == cftime.DatetimeNoLeap(2, 1, 1, 2, 1, 1)
-    assert ds["time"].values[1] == cftime.DatetimeNoLeap(2, 1, 1, 3, 1, 1)
-    assert ds["time"].values[2] == cftime.DatetimeNoLeap(2, 1, 1, 4, 1, 1)
-    np.testing.assert_array_equal(ds["var_2d"].values[0, :, :], var2_global.transpose())
-    np.testing.assert_array_equal(ds["var_2d"].values[1, :, :], var2_global.transpose())
-    np.testing.assert_array_equal(ds["var_2d"].values[2, :, :], var2_global.transpose())
+    ds = Dataset("diag_manager_single_tile.nc")
+    assert "var_2d" in ds.variables
+    assert "time" in ds.variables
+    assert "var_3d" in ds.variables
+    var2_ds = ds.variables["var_2d"]
+    time_var = ds.variables["time"]
+    var3_ds = ds.variables["var_3d"]
+    np.testing.assert_array_equal(var2_ds.shape, (ntimesteps, nx, ny))
+    assert var2_ds.dimensions == ("time", "y", "x")
+    assert var2_ds.units == "muntin"
+    assert var3_ds.dimensions == ("time", "z", "y", "x")
+    assert var3_ds.units == "muntin"
+    assert time_var.shape == (ntimesteps,)
+    assert time_var.dimensions == ("time",)
+    dates = num2date(time_var[:], units=time_var.units, calendar=time_var.calendar)
+    assert dates[0] == cftime.DatetimeNoLeap(2, 1, 1, 2, 1, 1)
+    assert dates[1] == cftime.DatetimeNoLeap(2, 1, 1, 3, 1, 1)
+    assert dates[2] == cftime.DatetimeNoLeap(2, 1, 1, 4, 1, 1)
+    np.testing.assert_array_equal(var2_ds[0, :, :], var2_global.transpose())
+    np.testing.assert_array_equal(var2_ds[1, :, :], var2_global.transpose())
+    np.testing.assert_array_equal(var2_ds[2, :, :], var2_global.transpose())
     # data is transposed when passed into fortran
-    np.testing.assert_array_equal(
-        ds["var_3d"].values[0, :, :, :], var3_global.transpose()
-    )
-    np.testing.assert_array_equal(
-        ds["var_3d"].values[1, :, :, :], var3_global.transpose()
-    )
-    np.testing.assert_array_equal(
-        ds["var_3d"].values[2, :, :, :], var3_global.transpose()
-    )
+    np.testing.assert_array_equal(var3_ds[0, :, :, :], var3_global.transpose())
+    np.testing.assert_array_equal(var3_ds[1, :, :, :], var3_global.transpose())
+    np.testing.assert_array_equal(var3_ds[2, :, :, :], var3_global.transpose())
 
     pyfms.fms.end()
diff --git a/tests/test_netcdf_monitor.py b/tests/monitor/test_netcdf_monitor.py
similarity index 100%
rename from tests/test_netcdf_monitor.py
rename to tests/monitor/test_netcdf_monitor.py
diff --git a/tests/test_zarr_monitor.py b/tests/monitor/test_zarr_monitor.py
similarity index 100%
rename from tests/test_zarr_monitor.py
rename to tests/monitor/test_zarr_monitor.py