diff --git a/docs/developers_guide/ocean/api.md b/docs/developers_guide/ocean/api.md
index c4ee9fd5c..7043f7642 100644
--- a/docs/developers_guide/ocean/api.md
+++ b/docs/developers_guide/ocean/api.md
@@ -66,7 +66,7 @@
    InertialGravityWave
 
    analysis.Analysis
-   analysis.Analysis.run
+   analysis.Analysis.exact_solution
 
    exact_solution.ExactSolution
    exact_solution.ExactSolution.ssh
@@ -74,7 +74,6 @@
 
    forward.Forward
    forward.Forward.compute_cell_count
-   forward.Forward.dynamic_model_config
 
    init.Init
    init.Init.run
@@ -127,7 +126,7 @@
    ManufacturedSolution
 
    analysis.Analysis
-   analysis.Analysis.run
+   analysis.Analysis.exact_solution
 
    exact_solution.ExactSolution
    exact_solution.ExactSolution.ssh
@@ -135,7 +134,6 @@
 
    forward.Forward
    forward.Forward.compute_cell_count
-   forward.Forward.dynamic_model_config
 
    init.Init
    init.Init.run
@@ -169,6 +167,27 @@
 
 ## Ocean Framework
 
+### Convergence Tests
+
+```{eval-rst}
+.. currentmodule:: polaris.ocean.convergence
+
+.. autosummary::
+   :toctree: generated/
+
+   ConvergenceForward
+   ConvergenceForward.compute_cell_count
+   ConvergenceForward.dynamic_model_config
+
+   ConvergenceAnalysis
+   ConvergenceAnalysis.compute_error
+   ConvergenceAnalysis.convergence_parameters
+   ConvergenceAnalysis.exact_solution
+   ConvergenceAnalysis.get_output_field
+   ConvergenceAnalysis.plot_convergence
+   ConvergenceAnalysis.run
+   ConvergenceAnalysis.setup
+```
 ### Spherical Convergence Tests
 
 ```{eval-rst}
@@ -178,13 +197,7 @@
    :toctree: generated/
 
    SphericalConvergenceForward
-   SphericalConvergenceAnalysis
-   SphericalConvergenceAnalysis.compute_error
-   SphericalConvergenceAnalysis.convergence_parameters
-   SphericalConvergenceAnalysis.exact_solution
-   SphericalConvergenceAnalysis.get_output_field
-   SphericalConvergenceAnalysis.run
-   SphericalConvergenceAnalysis.plot_convergence
+   SphericalConvergenceForward.compute_cell_count
 ```
 
 ### Ocean Model
diff --git a/docs/developers_guide/ocean/framework.md b/docs/developers_guide/ocean/framework.md
index ab0b524d3..dfcea88ca 100644
--- a/docs/developers_guide/ocean/framework.md
+++ b/docs/developers_guide/ocean/framework.md
@@ -136,14 +136,14 @@ The function {py:func}`polaris.ocean.mesh.spherical.add_spherical_base_mesh_step
 returns a step for for a spherical `qu` or `icos` mesh of a given resolution 
 (in km).  The step can be shared between tasks.
 
-(dev-ocean-spherical-convergence)=
+(dev-ocean-convergence)=
 
-## Spherical Convergence Tests
+## Convergence Tests
 
 Several tests that are in Polaris or which we plan to add are convergence
-tests on {ref}`dev-ocean-spherical-meshes`. The ocean framework includes
-shared config options and a base class for forward steps that are expected
-to be useful across these tests.
+tests on {ref}`dev-ocean-spherical-meshes` and planar meshes.
+The ocean framework includes shared config options and base classes for
+forward and analysis steps that are expected to be useful across these tests.
 
 The shared config options are:
 ```cfg
@@ -156,8 +156,10 @@ icos_resolutions = 60, 120, 240, 480
 # a list of quasi-uniform mesh resolutions (km) to test
 qu_resolutions = 60, 90, 120, 150, 180, 210, 240
 
-# Evaluation time for convergence analysis (in days)
-convergence_eval_time = 1.0
+[convergence]
+
+# Evaluation time for convergence analysis (in hours)
+convergence_eval_time = 24.0
 
 # Convergence threshold below which a test fails
 convergence_thresh = 1.0
@@ -165,9 +167,8 @@ convergence_thresh = 1.0
 # Type of error to compute
 error_type = l2
 
-
-# config options for spherical convergence forward steps
-[spherical_convergence_forward]
+# config options for convergence forward steps
+[convergence_forward]
 
 # time integrator: {'split_explicit', 'RK4'}
 time_integrator = RK4
@@ -182,17 +183,19 @@ split_dt_per_km = 30.0
 # since btr_dt is proportional to resolution
 btr_dt_per_km = 1.5
 
-# Run duration in days
-run_duration = ${spherical_convergence:convergence_eval_time}
+# Run duration in hours
+run_duration = ${convergence:convergence_eval_time}
 
-# Output interval in days
+# Output interval in hours
 output_interval = ${run_duration}
 ```
 The first 2 are the default resolutions for icosahedral and quasi-uniform
-base meshes, respectively.  The `convergence_eval_time` will generally be
-modified by each test case. The `convergence_thresh` will also be modified by
-each test case, and will depend on the numerical methods being tested. The
-`error_type` is the L2 norm by default.
+base meshes, respectively.
+
+The `convergence_eval_time` will generally be modified by each test case. The
+`convergence_thresh` will also be modified by each test case, and will depend
+on the numerical methods being tested. The `error_type` is the L2 norm by
+default. The L-infty norm, `inf`, is also supported.
 
 `time_integrator` will typically be overridden by the specific convergence
 task's config options, and indicates which time integrator to use for the
@@ -201,7 +204,7 @@ forward run.  Depending on the time integrator, either `rk4_dt_per_km` or
 mesh resolution (proportional to the cell size). For split time integrators,
 `btr_dt_per_km` will be used to compute the barotropic time step in a similar
 way.  The `run_duration` and `output_interval` are typically the same, and
-they are given in days.
+they are given in hours.
 
 Each convergence test can override these defaults with its own defaults by 
 defining them in its own config file.  Convergence tests should bring in this
@@ -218,6 +221,8 @@ def add_cosine_bell_tasks(component):
 
         filepath = f'spherical/{prefix}/cosine_bell/cosine_bell.cfg'
         config = PolarisConfigParser(filepath=filepath)
+        config.add_from_package('polaris.ocean.convergence',
+                                'convergence.cfg')
         config.add_from_package('polaris.ocean.convergence.spherical',
                                 'spherical.cfg')
         config.add_from_package('polaris.ocean.tasks.cosine_bell',
@@ -228,9 +233,13 @@ In addition, the {py:class}`polaris.ocean.convergence.spherical.SphericalConverg
 step can serve as a parent class for forward steps in convergence tests.  This
 parent class takes care of setting the time step based on the `dt_per_km`
 config option and computes the approximate number of cells in the mesh, used
-for determining the computational resources required, using a heuristic 
-appropriate for approximately uniform spherical meshes.  A convergence test's
-`Forward` step should descend from this class like in this example:
+for determining the computational resources required. When convergence tests
+are run on spherical meshes,
+the {py:class}`polaris.ocean.convergence.spherical.SphericalConvergenceForward`
+should be invoked and overrides the `compute_cell_count` method with a
+heuristic appropriate for approximately uniform spherical meshes.  A
+convergence test's `Forward` step should descend from this class like in this
+example:
 
 ```python
 from polaris.ocean.convergence.spherical import SphericalConvergenceForward
@@ -242,7 +251,7 @@ class Forward(SphericalConvergenceForward):
     bell test case
     """
 
-    def __init__(self, component, name, subdir, resolution, base_mesh, init):
+    def __init__(self, component, name, subdir, resolution, mesh, init):
         """
         Create a new step
 
@@ -260,7 +269,7 @@ class Forward(SphericalConvergenceForward):
         resolution : float
             The resolution of the (uniform) mesh in km
 
-        base_mesh : polaris.Step
+        mesh : polaris.Step
             The base mesh step
 
         init : polaris.Step
@@ -269,7 +278,7 @@ class Forward(SphericalConvergenceForward):
         package = 'polaris.ocean.tasks.cosine_bell'
         validate_vars = ['normalVelocity', 'tracer1']
         super().__init__(component=component, name=name, subdir=subdir,
-                         resolution=resolution, base_mesh=base_mesh,
+                         resolution=resolution, mesh=mesh,
                          init=init, package=package,
                          yaml_filename='forward.yaml',
                          output_filename='output.nc',
@@ -286,9 +295,9 @@ analyze.  The `validate_vars` are a list of variables to compare against a
 baseline (if one is provided), and can be `None` if baseline validation should
 not be performed.
 
-The `base_mesh` step should be created with the function described in
+The `mesh` step should be created with the function described in
 {ref}`dev-ocean-spherical-meshes`, and the `init` step should produce a file
-`initial_state.nc` that will be the initial condition for the forward run.
+`init.nc` that will be the initial condition for the forward run.
 
 The `forward.yaml` file should be a YAML file with Jinja templating for the 
 time integrator, time step, run duration and output interval, e.g.:
@@ -318,11 +327,11 @@ omega:
       - normalVelocity
       - layerThickness
 ```
-`SphericalConvergenceForward` takes care of filling in the template based
+`ConvergenceForward` takes care of filling in the template based
 on the associated config options (first at setup and again at runtime in case
 the config options have changed).
 
-In addition, the {py:class}`polaris.ocean.convergence.spherical.SphericalConvergenceAnalysis`
+In addition, the {py:class}`polaris.ocean.convergence.ConvergenceAnalysis`
 step can serve as a parent class for analysis steps in convergence tests.  This
 parent class computes the error norm for the output from each resolution's
 forward step. It also produces the convergence plot.
@@ -331,12 +340,11 @@ This is an example of how a task's analysis step can descend from the parent
 class:
 
 ```python
-class Analysis(SphericalConvergenceAnalysis):
+class Analysis(ConvergenceAnalysis):
     """
     A step for analyzing the output from the cosine bell test case
     """
-    def __init__(self, component, resolutions, icosahedral, subdir,
-                 dependencies):
+    def __init__(self, component, resolutions, subdir, dependencies):
         """
         Create the step
 
@@ -348,10 +356,6 @@ class Analysis(SphericalConvergenceAnalysis):
         resolutions : list of float
             The resolutions of the meshes that have been run
 
-        icosahedral : bool
-            Whether to use icosahedral, as opposed to less regular, JIGSAW
-            meshes
-
         subdir : str
             The subdirectory that the step resides in
 
@@ -363,18 +367,17 @@ class Analysis(SphericalConvergenceAnalysis):
                              'zidx': 0}]
         super().__init__(component=component, subdir=subdir,
                          resolutions=resolutions,
-                         icosahedral=icosahedral,
                          dependencies=dependencies,
                          convergence_vars=convergence_vars)
 ```
 
 Many tasks will also need to override the 
-{py:meth}`polaris.ocean.convergence.spherical.SphericalConvergenceAnalysis.exact_solution()` 
+{py:meth}`polaris.ocean.convergence.ConvergenceAnalysis.exact_solution()` 
 method. If not overridden, the analysis step will compute the difference of the 
 output from the initial state.
 
 In some cases, the child class will also need to override the 
-{py:meth}`polaris.ocean.convergence.spherical.SphericalConvergenceAnalysis.get_output_field()`
+{py:meth}`polaris.ocean.convergence.ConvergenceAnalysis.get_output_field()`
 method if the requested field is not available directly from the output put
 rather needs to be computed.  The default behavior is to read the requested
 variable (the value associate the `'name'` key) at the time index closest to
diff --git a/docs/developers_guide/ocean/tasks/cosine_bell.md b/docs/developers_guide/ocean/tasks/cosine_bell.md
index e825e39f4..4d4cdef2f 100644
--- a/docs/developers_guide/ocean/tasks/cosine_bell.md
+++ b/docs/developers_guide/ocean/tasks/cosine_bell.md
@@ -18,7 +18,7 @@ model config options related to drag and default horizontal and
 vertical momentum and tracer diffusion, as well as defining `mesh`, `input`,
 `restart`, and `output` streams.  This file has Jinja templating that is
 used to update model config options based on Polaris config options, see
-{ref}`dev-ocean-spherical-convergence`.
+{ref}`dev-ocean-convergence`.
 
 ### base_mesh
 
@@ -36,16 +36,16 @@ tracer distributed in a cosine-bell shape.
 The class {py:class}`polaris.ocean.tasks.cosine_bell.forward.Forward`
 descends from {py:class}`polaris.ocean.convergence.spherical.SphericalConvergenceForward`,
 and defines a step for running MPAS-Ocean from an initial condition produced in
-an `init` step. See {ref}`dev-ocean-spherical-convergence` for some relevant
+an `init` step. See {ref}`dev-ocean-convergence` for some relevant
 discussion of the parent class. The time step is determined from the resolution
-based on the `dt_per_km` config option in the `[spherical_convergences]` 
+based on the `dt_per_km` config option in the `[convergence_forward]` 
 section.  Other model config options are taken from `forward.yaml`.
 
 ### analysis
 
 The class {py:class}`polaris.ocean.tasks.cosine_bell.analysis.Analysis`
 descends from
-{py:class}`polaris.ocean.convergence.spherical.SphericalConvergenceAnalysis`,
+{py:class}`polaris.ocean.convergence.ConvergenceAnalysis`,
 and defines a step for computing the error norm (L2) for the results
 at each resolution, saving them in `convergence_tracer1.csv` and plotting them
 in `convergence_tracer1.png`.
diff --git a/docs/developers_guide/ocean/tasks/inertial_gravity_wave.md b/docs/developers_guide/ocean/tasks/inertial_gravity_wave.md
index 100576758..751be935f 100644
--- a/docs/developers_guide/ocean/tasks/inertial_gravity_wave.md
+++ b/docs/developers_guide/ocean/tasks/inertial_gravity_wave.md
@@ -48,10 +48,11 @@ the exact solution. The tracer and coriolis fields are uniform in space.
 ### forward
 
 The class {py:class}`polaris.ocean.tasks.inertial_gravity_wave.forward.Forward`
-defines a step for running MPAS-Ocean from the initial condition produced in the
-`init` step.  Namelist and streams files are updated in
-{py:meth}`polaris.ocean.tasks.inertial_gravity_wave.forward.Forward.dynamic_model_config()`
-with time steps determined algorithmically based on config options.  The number
+descends from
+ {py:class}`polaris.ocean.convergence.ConvergenceForward`
+and runs MPAS-Ocean from the initial condition produced in the `init` step.
+Namelist and streams files are updated by the parent class with time steps
+determined algorithmically based on config options. The number
 of cells is approximated from config options in
 {py:meth}`polaris.ocean.tasks.inertial_gravity_wave.forward.Forward.compute_cell_count()`
 so that this can be used to constrain the number of MPI tasks that Polaris 
@@ -65,14 +66,14 @@ performed against a baseline if one is provided when calling
 ### analysis
 
 The class
-{py:class}`polaris.ocean.tasks.inertial_gravity_wave.analysis.Analysis` defines
-a step for computing the root mean-square-error from the final simulated field
+{py:class}`polaris.ocean.tasks.inertial_gravity_wave.analysis.Analysis`
+descends from {py:class}`polaris.ocean.convergence.ConvergenceAnalysis`
+a step for computing the error from the final simulated field
 and the exact solution. It uses the config options to determine whether the
 convergence rate falls within acceptable bounds.
 
 ### viz
 
 The class {py:class}`polaris.ocean.tasks.inertial_gravity_wave.viz.Viz` defines
-a step for visualization. It produces two plots: the convergence of the RMSE
-with resolution and a plan-view of the simulated, exact, and (simulated - exact)
-SSH fields.
+a step for visualization. It produces a plan-view figure of the simulated,
+exact, and (simulated - exact) SSH fields.
diff --git a/docs/developers_guide/ocean/tasks/manufactured_solution.md b/docs/developers_guide/ocean/tasks/manufactured_solution.md
index 14f00e49a..0e5bee95d 100644
--- a/docs/developers_guide/ocean/tasks/manufactured_solution.md
+++ b/docs/developers_guide/ocean/tasks/manufactured_solution.md
@@ -45,31 +45,32 @@ the exact solution. The tracer and coriolis fields are uniform in space.
 ### forward
 
 The class {py:class}`polaris.ocean.tasks.manufactured_solution.forward.Forward`
-defines a step for running MPAS-Ocean from the initial condition produced in
-the `init` step.  Namelist and streams files are updated in
-{py:meth}`polaris.ocean.tasks.manufactured_solution.forward.Forward.dynamic_model_config()`
-with time steps determined algorithmically based on config options.  The
-number of cells is approximated from config options in
-{py:meth}`polaris.ocean.tasks.manufactured_solution.forward.Forward.compute_cell_count()`
-so that this can be used to constrain the number of MPI tasks that Polaris tasks
-have as their target and minimum (if the resources are not explicitly
-prescribed).  For MPAS-Ocean, PIO namelist options are modified and a
-graph partition is generated as part of `runtime_setup()`.  Next, the ocean 
-model is run.  Finally, validation of `temperature`, `layerThickness` and 
+descends from
+ {py:class}`polaris.ocean.convergence.ConvergenceForward`
+and runs MPAS-Ocean from the initial condition produced in the `init` step.
+Namelist and streams files are updated by the parent class with time steps
+determined algorithmically based on config options. The number
+of cells is approximated from config options in
+{py:meth}`polaris.ocean.tasks.inertial_gravity_wave.forward.Forward.compute_cell_count()`
+so that this can be used to constrain the number of MPI tasks that Polaris 
+tasks have as  their target and minimum (if the resources are not explicitly 
+prescribed).  For MPAS-Ocean, PIO namelist options are modified and a graph 
+partition is generated as part of `runtime_setup()`.  Then, the ocean model 
+is run.  Finally, validation of `temperature`, `layerThickness` and 
 `normalVelocity` are performed against a baseline if one is provided when 
 calling {ref}`dev-polaris-setup`.
 
 ### analysis
 
 The class {py:class}`polaris.ocean.tasks.manufactured_solution.analysis.Analysis`
-defines a step for computing the root mean-square-error from the final
-simulated field and the exact solution. It uses the config options to determine
-whether the convergence rate falls within acceptable bounds.
+descends from {py:class}`polaris.ocean.convergence.ConvergenceAnalysis`
+a step for computing the error from the final simulated field
+and the exact solution. It uses the config options to determine whether the
+convergence rate falls within acceptable bounds.
 
 ### viz
 
 The class {py:class}`polaris.ocean.tasks.manufactured_solution.viz.Viz`
-defines a step for visualization. It produces two plots: the convergence of the
-RMSE with resolution and a plan-view of the simulated, exact, and (simulated -
-exact) SSH fields.
+defines a step for visualization. It produces a plan-view figure of the
+simulated, exact, and (simulated - exact) SSH fields.
 
diff --git a/docs/users_guide/ocean/tasks/cosine_bell.md b/docs/users_guide/ocean/tasks/cosine_bell.md
index 18b1898d7..7955e139f 100644
--- a/docs/users_guide/ocean/tasks/cosine_bell.md
+++ b/docs/users_guide/ocean/tasks/cosine_bell.md
@@ -167,8 +167,8 @@ alter this before setup (in a user config file) or before running the task (in
 the config file in the work directory). 
 
 ```cfg
-# config options for spherical convergence tests
-[spherical_convergence_forward]
+# config options for convergence tests
+[convergence_forward]
 
 # time integrator: {'split_explicit', 'RK4'}
 time_integrator = RK4
@@ -182,12 +182,12 @@ period for advection to make a full rotation around the globe, 24 days:
 
 ```cfg
 # config options for spherical convergence tests
-[spherical_convergence_forward]
+[convergence_forward]
 
-# Run duration in days
+# Run duration in hours
 run_duration = ${cosine_bell:vel_pd}
 
-# Output interval in days
+# Output interval in hours
 output_interval = ${cosine_bell:vel_pd}
 ```
 
@@ -220,8 +220,8 @@ radius = 2123666.6667
 # hill max of tracer
 psi0 = 1.0
 
-# time (days) for bell to transit equator once
-vel_pd = 24.0
+# time (hours) for bell to transit equator once
+vel_pd = 576.0
 
 # convergence threshold below which the test fails
 convergence_thresh = 1.8
@@ -254,8 +254,8 @@ norm_args = {'vmin': 0., 'vmax': 1.}
 The 7 options from `temperature` to `vel_pd` are used to control properties of
 the cosine bell and the rest of the sphere, as well as the advection.
 
-The options `qu_conv_thresh` to `icos_conv_max` are thresholds for determining
-when the convergence rates are not within the expected range.
+The option `convergence_thresh` is a threshold for determining
+when the convergence rates are not above a minimum convergence rate.
 
 The options in the `cosine_bell_viz` section are used in visualizing the
 initial and final states on a lon-lat grid for `cosine_bell/with_viz` tasks.
@@ -266,11 +266,14 @@ of these config options can be changed here:
 
 ```cfg
 # config options for spherical convergence tests
-[spherical_convergence]
+[convergence]
 
 # Evaluation time for convergence analysis (in days)
 convergence_eval_time = ${cosine_bell:vel_pd}
 
+# Convergence threshold below which a test fails
+convergence_thresh = ${cosine_bell:convergence_thresh}
+
 # Type of error to compute
 error_type = l2
 ```
diff --git a/docs/users_guide/ocean/tasks/images/inertial_gravity_wave_convergence.png b/docs/users_guide/ocean/tasks/images/inertial_gravity_wave_convergence.png
index 320ad561c..a07e76624 100644
Binary files a/docs/users_guide/ocean/tasks/images/inertial_gravity_wave_convergence.png and b/docs/users_guide/ocean/tasks/images/inertial_gravity_wave_convergence.png differ
diff --git a/docs/users_guide/ocean/tasks/inertial_gravity_wave.md b/docs/users_guide/ocean/tasks/inertial_gravity_wave.md
index 2783f6426..7c27b5992 100644
--- a/docs/users_guide/ocean/tasks/inertial_gravity_wave.md
+++ b/docs/users_guide/ocean/tasks/inertial_gravity_wave.md
@@ -20,11 +20,11 @@ and the pressure gradient used is the gradient of the sea surface height.
 Horizontal mixing and bottom friction are also neglected. The nonlinear momentum
 terms are not included and the layer thickness equation is linearized.
 
-The analysis step computes the root mean-square-error of the difference between
+The analysis step computes the L2-norm of the difference between
 the simulated SSH field and the exact solution at the end of the simulation. It
 also computes the convergence rate with resolution.
 
-The visualization step produces two plots: the convergence of the RMSE with
+The visualization step produces two plots: the convergence of the L2-norm with
 resolution and a plan-view of the simulated, exact, and (simulated-analytical)
 SSH fields.
 
@@ -114,9 +114,6 @@ n_wavelengths_y = 2
 # Convergence threshold below which the test fails
 conv_thresh = 1.8
 
-# Convergence rate above which a warning is issued
-conv_max = 2.2
-
 # time step per resolution (s/km), since dt is proportional to resolution
 dt_per_km = 3.0
 
diff --git a/docs/users_guide/ocean/tasks/manufactured_solution.md b/docs/users_guide/ocean/tasks/manufactured_solution.md
index 0016e894f..ae97febe0 100644
--- a/docs/users_guide/ocean/tasks/manufactured_solution.md
+++ b/docs/users_guide/ocean/tasks/manufactured_solution.md
@@ -26,11 +26,11 @@ model is configured without vertical advection and mixing. No tracers are enable
 and the pressure gradient used is the gradient of the sea surface height.
 Horizontal mixing and bottom friction are also neglected.
 
-The analysis step computes the root mean-square-error of the difference between
+The analysis step computes the L2-norm of the difference between
 the simulated SSH field and the exact solution at the end of the simulation. It
-also computes the convergence rate with resolution
+also computes the convergence rate with resolution.
 
-The visualization step produces two plots: the convergence of the RMSE with
+The visualization step produces two plots: the convergence of the L2-norm with
 resolution and a plan-view of the simulated, exact, and (simulated-exact)
 SSH fields.
 
@@ -111,9 +111,6 @@ dt_per_km = 3.0
 
 # Convergence threshold below which the test fails
 conv_thresh = 1.8
-
-# Convergence rate above which a warning is issued
-conv_max = 2.2
 ```
 
 ### cores
diff --git a/polaris/ocean/convergence/__init__.py b/polaris/ocean/convergence/__init__.py
index e69de29bb..7fb71d53f 100644
--- a/polaris/ocean/convergence/__init__.py
+++ b/polaris/ocean/convergence/__init__.py
@@ -0,0 +1,2 @@
+from polaris.ocean.convergence.analysis import ConvergenceAnalysis
+from polaris.ocean.convergence.forward import ConvergenceForward
diff --git a/polaris/ocean/convergence/spherical/analysis.py b/polaris/ocean/convergence/analysis.py
similarity index 92%
rename from polaris/ocean/convergence/spherical/analysis.py
rename to polaris/ocean/convergence/analysis.py
index f2bdebc39..2e39fcc45 100644
--- a/polaris/ocean/convergence/spherical/analysis.py
+++ b/polaris/ocean/convergence/analysis.py
@@ -11,19 +11,15 @@
 from polaris.viz import use_mplstyle
 
 
-class SphericalConvergenceAnalysis(Step):
+class ConvergenceAnalysis(Step):
     """
-    A step for analyzing the output from the geostrophic convergence test
+    A step for analyzing the output from convergence tests
 
     Attributes
     ----------
     resolutions : list of float
         The resolutions of the meshes that have been run
 
-    icosahedral : bool
-        Whether to use icosahedral, as opposed to less regular, JIGSAW
-        meshes
-
     dependencies_dict : dict of dict of polaris.Steps
         The dependencies of this step must be given as separate keys in the
         dict:
@@ -58,8 +54,8 @@ class SphericalConvergenceAnalysis(Step):
                 The z-index to use for variables that have an nVertLevels
                 dimension, which should be None for variables that don't
     """
-    def __init__(self, component, resolutions, icosahedral, subdir,
-                 dependencies, convergence_vars):
+    def __init__(self, component, resolutions, subdir, dependencies,
+                 convergence_vars):
         """
         Create the step
 
@@ -71,10 +67,6 @@ def __init__(self, component, resolutions, icosahedral, subdir,
         resolutions : list of float
             The resolutions of the meshes that have been run
 
-        icosahedral : bool
-            Whether to use icosahedral, as opposed to less regular, JIGSAW
-            meshes
-
         subdir : str
             The subdirectory that the step resides in
 
@@ -114,7 +106,6 @@ def __init__(self, component, resolutions, icosahedral, subdir,
         """
         super().__init__(component=component, name='analysis', subdir=subdir)
         self.resolutions = resolutions
-        self.icosahedral = icosahedral
         self.dependencies_dict = dependencies
         self.convergence_vars = convergence_vars
 
@@ -252,7 +243,7 @@ def plot_convergence(self, variable_name, title, zidx):
         logger.info(f'Order of convergence for {title}: '
                     f'{conv_round}')
 
-        if convergence < conv_thresh:
+        if conv_round < conv_thresh:
             logger.error(f'Error: order of convergence for {title}\n'
                          f'  {conv_round} < min tolerance '
                          f'{conv_thresh}')
@@ -286,30 +277,31 @@ def compute_error(self, mesh_name, variable_name, zidx=None,
         error : float
             The error of the variable given by variable_name
         """
+        norm_type = {'l2': None, 'inf': np.inf}
         ds_mesh = xr.open_dataset(f'{mesh_name}_mesh.nc')
         config = self.config
-        section = config['spherical_convergence']
+        section = config['convergence']
         eval_time = section.getfloat('convergence_eval_time')
-        s_per_day = 86400.0
+        s_per_hour = 3600.0
 
         field_exact = self.exact_solution(mesh_name, variable_name,
-                                          time=eval_time * s_per_day,
+                                          time=eval_time * s_per_hour,
                                           zidx=zidx)
         field_mpas = self.get_output_field(mesh_name, variable_name,
-                                           time=eval_time * s_per_day,
+                                           time=eval_time * s_per_hour,
                                            zidx=zidx)
         diff = field_exact - field_mpas
 
+        # Only the L2 norm is area-weighted
         if error_type == 'l2':
             area = area_for_field(ds_mesh, diff)
-            total_area = np.sum(area)
-            den_l2 = np.sum(field_exact**2 * area) / total_area
-            num_l2 = np.sum(diff**2 * area) / total_area
-            error = np.sqrt(num_l2) / np.sqrt(den_l2)
-        elif error_type == 'inf':
-            error = np.amax(diff) / np.amax(np.abs(field_exact))
-        else:
-            raise ValueError(f'Unsupported error type {error_type}')
+            field_exact = field_exact * area
+            diff = diff * area
+        error = np.linalg.norm(diff, ord=norm_type[error_type])
+
+        # Normalize the error norm by the vector norm of the exact solution
+        den = np.linalg.norm(field_exact, ord=norm_type[error_type])
+        error = np.divide(error, den)
 
         return error
 
@@ -401,7 +393,7 @@ def convergence_parameters(self, field_name=None):
             The error norm to compute
         """
         config = self.config
-        section = config['spherical_convergence']
+        section = config['convergence']
         conv_thresh = section.getfloat('convergence_thresh')
         error_type = section.get('error_type')
         return conv_thresh, error_type
diff --git a/polaris/ocean/convergence/convergence.cfg b/polaris/ocean/convergence/convergence.cfg
new file mode 100644
index 000000000..e4e202244
--- /dev/null
+++ b/polaris/ocean/convergence/convergence.cfg
@@ -0,0 +1,32 @@
+[convergence]
+
+# Evaluation time for convergence analysis (in hours)
+convergence_eval_time = 24.0
+
+# Convergence threshold below which a test fails
+convergence_thresh = 1.0
+
+# Type of error to compute
+error_type = l2
+
+# config options for convergence forward steps
+[convergence_forward]
+
+# time integrator: {'split_explicit', 'RK4'}
+time_integrator = RK4
+
+# RK4 time step per resolution (s/km), since dt is proportional to resolution
+rk4_dt_per_km = 3.0
+
+# split time step per resolution (s/km), since dt is proportional to resolution
+split_dt_per_km = 30.0
+
+# the barotropic time step (s/km) for simulations using split time stepping,
+# since btr_dt is proportional to resolution
+btr_dt_per_km = 1.5
+
+# Run duration in hours
+run_duration = ${convergence:convergence_eval_time}
+
+# Output interval in hours
+output_interval = ${run_duration}
diff --git a/polaris/ocean/convergence/forward.py b/polaris/ocean/convergence/forward.py
new file mode 100644
index 000000000..d1754cc85
--- /dev/null
+++ b/polaris/ocean/convergence/forward.py
@@ -0,0 +1,148 @@
+from polaris.ocean.model import OceanModelStep, get_time_interval_string
+
+
+class ConvergenceForward(OceanModelStep):
+    """
+    A step for performing forward ocean component runs as part of a spherical
+    convergence test
+
+    Attributes
+    ----------
+    resolution : float
+        The resolution of the (uniform) mesh in km
+
+    package : Package
+        The package name or module object that contains the YAML file
+
+    yaml_filename : str
+        A YAML file that is a Jinja2 template with model config options
+
+    """
+
+    def __init__(self, component, name, subdir, resolution, mesh, init,
+                 package, yaml_filename='forward.yaml', options=None,
+                 graph_filename='graph.info', output_filename='output.nc',
+                 validate_vars=None):
+        """
+        Create a new step
+
+        Parameters
+        ----------
+        component : polaris.Component
+            The component the step belongs to
+
+        name : str
+            The name of the step
+
+        subdir : str
+            The subdirectory for the step
+
+        resolution : float
+            The resolution of the (uniform) mesh in km
+
+        package : Package
+            The package name or module object that contains the YAML file
+
+        yaml_filename : str, optional
+            A YAML file that is a Jinja2 template with model config options
+
+        options : dict, optional
+            A dictionary of options and value to replace model config options
+            with new values
+
+        output_filename : str, optional
+            The output file that will be written out at the end of the forward
+            run
+
+        validate_vars : list of str, optional
+            A list of variables to validate against a baseline if requested
+        """
+        super().__init__(component=component, name=name, subdir=subdir,
+                         openmp_threads=1)
+
+        self.resolution = resolution
+        self.package = package
+        self.yaml_filename = yaml_filename
+
+        # make sure output is double precision
+        self.add_yaml_file('polaris.ocean.config', 'output.yaml')
+
+        if options is not None:
+            self.add_model_config_options(options=options)
+
+        self.add_input_file(
+            filename='init.nc',
+            work_dir_target=f'{init.path}/initial_state.nc')
+        self.add_input_file(
+            filename='graph.info',
+            work_dir_target=f'{mesh.path}/{graph_filename}')
+
+        self.add_output_file(filename=output_filename,
+                             validate_vars=validate_vars)
+
+    def compute_cell_count(self):
+        """
+        Compute the approximate number of cells in the mesh, used to constrain
+        resources
+
+        Returns
+        -------
+        cell_count : int or None
+            The approximate number of cells in the mesh
+        """
+        raise ValueError('compute_cell_count method must be overridden by '
+                         'spherical or planar method')
+
+    def dynamic_model_config(self, at_setup):
+        """
+        Set the model time step from config options at setup and runtime
+
+        Parameters
+        ----------
+        at_setup : bool
+            Whether this method is being run during setup of the step, as
+            opposed to at runtime
+        """
+        super().dynamic_model_config(at_setup=at_setup)
+
+        config = self.config
+
+        vert_levels = config.getfloat('vertical_grid', 'vert_levels')
+        if not at_setup and vert_levels == 1:
+            self.add_yaml_file('polaris.ocean.config', 'single_layer.yaml')
+
+        section = config['convergence_forward']
+
+        time_integrator = section.get('time_integrator')
+
+        # dt is proportional to resolution: default 30 seconds per km
+        if time_integrator == 'RK4':
+            dt_per_km = section.getfloat('rk4_dt_per_km')
+        else:
+            dt_per_km = section.getfloat('split_dt_per_km')
+        dt_str = get_time_interval_string(seconds=dt_per_km * self.resolution)
+
+        # btr_dt is also proportional to resolution: default 1.5 seconds per km
+        btr_dt_per_km = section.getfloat('btr_dt_per_km')
+        btr_dt_str = get_time_interval_string(
+            seconds=btr_dt_per_km * self.resolution)
+
+        s_per_hour = 3600.
+        run_duration = section.getfloat('run_duration')
+        run_duration_str = get_time_interval_string(
+            seconds=run_duration * s_per_hour)
+
+        output_interval = section.getfloat('output_interval')
+        output_interval_str = get_time_interval_string(
+            seconds=output_interval * s_per_hour)
+
+        replacements = dict(
+            time_integrator=time_integrator,
+            dt=dt_str,
+            btr_dt=btr_dt_str,
+            run_duration=run_duration_str,
+            output_interval=output_interval_str,
+        )
+
+        self.add_yaml_file(self.package, self.yaml_filename,
+                           template_replacements=replacements)
diff --git a/polaris/ocean/convergence/spherical/__init__.py b/polaris/ocean/convergence/spherical/__init__.py
index 418252ded..702a1352f 100644
--- a/polaris/ocean/convergence/spherical/__init__.py
+++ b/polaris/ocean/convergence/spherical/__init__.py
@@ -1,6 +1,3 @@
-from polaris.ocean.convergence.spherical.analysis import (
-    SphericalConvergenceAnalysis,
-)
 from polaris.ocean.convergence.spherical.forward import (
     SphericalConvergenceForward,
 )
diff --git a/polaris/ocean/convergence/spherical/forward.py b/polaris/ocean/convergence/spherical/forward.py
index aa14f26d8..88be45217 100644
--- a/polaris/ocean/convergence/spherical/forward.py
+++ b/polaris/ocean/convergence/spherical/forward.py
@@ -1,7 +1,7 @@
-from polaris.ocean.model import OceanModelStep, get_time_interval_string
+from polaris.ocean.convergence import ConvergenceForward
 
 
-class SphericalConvergenceForward(OceanModelStep):
+class SphericalConvergenceForward(ConvergenceForward):
     """
     A step for performing forward ocean component runs as part of a spherical
     convergence test
@@ -19,66 +19,6 @@ class SphericalConvergenceForward(OceanModelStep):
 
     """
 
-    def __init__(self, component, name, subdir, resolution, base_mesh, init,
-                 package, yaml_filename='forward.yaml', options=None,
-                 output_filename='output.nc', validate_vars=None):
-        """
-        Create a new step
-
-        Parameters
-        ----------
-        component : polaris.Component
-            The component the step belongs to
-
-        name : str
-            The name of the step
-
-        subdir : str
-            The subdirectory for the step
-
-        resolution : float
-            The resolution of the (uniform) mesh in km
-
-        package : Package
-            The package name or module object that contains the YAML file
-
-        yaml_filename : str, optional
-            A YAML file that is a Jinja2 template with model config options
-
-        options : dict, optional
-            A dictionary of options and value to replace model config options
-            with new values
-
-        output_filename : str, optional
-            The output file that will be written out at the end of the forward
-            run
-
-        validate_vars : list of str, optional
-            A list of variables to validate against a baseline if requested
-        """
-        super().__init__(component=component, name=name, subdir=subdir,
-                         openmp_threads=1)
-
-        self.resolution = resolution
-        self.package = package
-        self.yaml_filename = yaml_filename
-
-        # make sure output is double precision
-        self.add_yaml_file('polaris.ocean.config', 'output.yaml')
-
-        if options is not None:
-            self.add_model_config_options(options=options)
-
-        self.add_input_file(
-            filename='init.nc',
-            work_dir_target=f'{init.path}/initial_state.nc')
-        self.add_input_file(
-            filename='graph.info',
-            work_dir_target=f'{base_mesh.path}/graph.info')
-
-        self.add_output_file(filename=output_filename,
-                             validate_vars=validate_vars)
-
     def compute_cell_count(self):
         """
         Compute the approximate number of cells in the mesh, used to constrain
@@ -92,54 +32,3 @@ def compute_cell_count(self):
         # use a heuristic based on QU30 (65275 cells) and QU240 (10383 cells)
         cell_count = 6e8 / self.resolution**2
         return cell_count
-
-    def dynamic_model_config(self, at_setup):
-        """
-        Set the model time step from config options at setup and runtime
-
-        Parameters
-        ----------
-        at_setup : bool
-            Whether this method is being run during setup of the step, as
-            opposed to at runtime
-        """
-        super().dynamic_model_config(at_setup=at_setup)
-
-        config = self.config
-
-        vert_levels = config.getfloat('vertical_grid', 'vert_levels')
-        if not at_setup and vert_levels == 1:
-            self.add_yaml_file('polaris.ocean.config', 'single_layer.yaml')
-
-        section = config['spherical_convergence_forward']
-
-        time_integrator = section.get('time_integrator')
-
-        # dt is proportional to resolution: default 30 seconds per km
-        if time_integrator == 'RK4':
-            dt_per_km = section.getfloat('rk4_dt_per_km')
-        else:
-            dt_per_km = section.getfloat('split_dt_per_km')
-        dt_str = get_time_interval_string(seconds=dt_per_km * self.resolution)
-
-        # btr_dt is also proportional to resolution: default 1.5 seconds per km
-        btr_dt_per_km = section.getfloat('btr_dt_per_km')
-        btr_dt_str = get_time_interval_string(
-            seconds=btr_dt_per_km * self.resolution)
-
-        run_duration = section.getfloat('run_duration')
-        run_duration_str = get_time_interval_string(days=run_duration)
-
-        output_interval = section.getfloat('output_interval')
-        output_interval_str = get_time_interval_string(days=output_interval)
-
-        replacements = dict(
-            time_integrator=time_integrator,
-            dt=dt_str,
-            btr_dt=btr_dt_str,
-            run_duration=run_duration_str,
-            output_interval=output_interval_str,
-        )
-
-        self.add_yaml_file(self.package, self.yaml_filename,
-                           template_replacements=replacements)
diff --git a/polaris/ocean/convergence/spherical/spherical.cfg b/polaris/ocean/convergence/spherical/spherical.cfg
index 630250da4..69ea2dd22 100644
--- a/polaris/ocean/convergence/spherical/spherical.cfg
+++ b/polaris/ocean/convergence/spherical/spherical.cfg
@@ -6,35 +6,3 @@ icos_resolutions = 60, 120, 240, 480
 
 # a list of quasi-uniform mesh resolutions (km) to test
 qu_resolutions = 60, 90, 120, 150, 180, 210, 240
-
-# Evaluation time for convergence analysis (in days)
-convergence_eval_time = 1.0
-
-# Convergence threshold below which a test fails
-convergence_thresh = 1.0
-
-# Type of error to compute
-error_type = l2
-
-
-# config options for spherical convergence forward steps
-[spherical_convergence_forward]
-
-# time integrator: {'split_explicit', 'RK4'}
-time_integrator = RK4
-
-# RK4 time step per resolution (s/km), since dt is proportional to resolution
-rk4_dt_per_km = 3.0
-
-# split time step per resolution (s/km), since dt is proportional to resolution
-split_dt_per_km = 30.0
-
-# the barotropic time step (s/km) for simulations using split time stepping,
-# since btr_dt is proportional to resolution
-btr_dt_per_km = 1.5
-
-# Run duration in days
-run_duration = ${spherical_convergence:convergence_eval_time}
-
-# Output interval in days
-output_interval = ${run_duration}
diff --git a/polaris/ocean/suites/convergence.txt b/polaris/ocean/suites/convergence.txt
new file mode 100644
index 000000000..6c01e8ccf
--- /dev/null
+++ b/polaris/ocean/suites/convergence.txt
@@ -0,0 +1,10 @@
+ocean/planar/inertial_gravity_wave
+ocean/planar/manufactured_solution
+ocean/spherical/icos/cosine_bell
+  cached: icos_base_mesh_60km icos_init_60km icos_base_mesh_120km icos_init_120km
+  cached: icos_base_mesh_240km icos_init_240km icos_base_mesh_480km icos_init_480km
+ocean/spherical/qu/cosine_bell
+  cached: qu_base_mesh_60km qu_init_60km qu_base_mesh_90km qu_init_90km
+  cached: qu_base_mesh_120km qu_init_120km qu_base_mesh_150km qu_init_150km
+  cached: qu_base_mesh_180km qu_init_180km qu_base_mesh_210km qu_init_210km
+  cached: qu_base_mesh_240km qu_init_240km
diff --git a/polaris/ocean/tasks/cosine_bell/__init__.py b/polaris/ocean/tasks/cosine_bell/__init__.py
index c3144c534..64433d1c1 100644
--- a/polaris/ocean/tasks/cosine_bell/__init__.py
+++ b/polaris/ocean/tasks/cosine_bell/__init__.py
@@ -21,6 +21,8 @@ def add_cosine_bell_tasks(component):
 
         filepath = f'spherical/{prefix}/cosine_bell/cosine_bell.cfg'
         config = PolarisConfigParser(filepath=filepath)
+        config.add_from_package('polaris.ocean.convergence',
+                                'convergence.cfg')
         config.add_from_package('polaris.ocean.convergence.spherical',
                                 'spherical.cfg')
         config.add_from_package('polaris.ocean.tasks.cosine_bell',
@@ -160,7 +162,7 @@ def _setup_steps(self):
             else:
                 forward_step = Forward(component=component, name=name,
                                        subdir=subdir, resolution=resolution,
-                                       base_mesh=base_mesh_step,
+                                       mesh=base_mesh_step,
                                        init=init_step)
                 forward_step.set_shared_config(config, link=config_filename)
             self.add_step(forward_step, symlink=symlink)
@@ -197,7 +199,7 @@ def _setup_steps(self):
             step.dependencies_dict = analysis_dependencies
         else:
             step = Analysis(component=component, resolutions=resolutions,
-                            icosahedral=icosahedral, subdir=subdir,
+                            subdir=subdir,
                             dependencies=analysis_dependencies)
             step.set_shared_config(config, link=config_filename)
         self.add_step(step, symlink=symlink)
diff --git a/polaris/ocean/tasks/cosine_bell/analysis.py b/polaris/ocean/tasks/cosine_bell/analysis.py
index 2d77ded8e..7ba98f8ae 100644
--- a/polaris/ocean/tasks/cosine_bell/analysis.py
+++ b/polaris/ocean/tasks/cosine_bell/analysis.py
@@ -3,16 +3,15 @@
 from mpas_tools.transects import lon_lat_to_cartesian
 from mpas_tools.vector import Vector
 
-from polaris.ocean.convergence.spherical import SphericalConvergenceAnalysis
+from polaris.ocean.convergence import ConvergenceAnalysis
 from polaris.ocean.tasks.cosine_bell.init import cosine_bell
 
 
-class Analysis(SphericalConvergenceAnalysis):
+class Analysis(ConvergenceAnalysis):
     """
     A step for analyzing the output from the cosine bell test case
     """
-    def __init__(self, component, resolutions, icosahedral, subdir,
-                 dependencies):
+    def __init__(self, component, resolutions, subdir, dependencies):
         """
         Create the step
 
@@ -24,10 +23,6 @@ def __init__(self, component, resolutions, icosahedral, subdir,
         resolutions : list of float
             The resolutions of the meshes that have been run
 
-        icosahedral : bool
-            Whether to use icosahedral, as opposed to less regular, JIGSAW
-            meshes
-
         subdir : str
             The subdirectory that the step resides in
 
@@ -39,7 +34,6 @@ def __init__(self, component, resolutions, icosahedral, subdir,
                              'zidx': 0}]
         super().__init__(component=component, subdir=subdir,
                          resolutions=resolutions,
-                         icosahedral=icosahedral,
                          dependencies=dependencies,
                          convergence_vars=convergence_vars)
 
@@ -91,7 +85,7 @@ def exact_solution(self, mesh_name, field_name, time, zidx=None):
 
         # distance that the cosine bell center traveled in radians
         # based on equatorial velocity
-        distance = 2.0 * np.pi * time / (86400.0 * vel_pd)
+        distance = 2.0 * np.pi * time / (3600.0 * vel_pd)
 
         # new location of blob center
         lon_new = lon_center + distance
diff --git a/polaris/ocean/tasks/cosine_bell/cosine_bell.cfg b/polaris/ocean/tasks/cosine_bell/cosine_bell.cfg
index c2b8b5ff7..fb17734cb 100644
--- a/polaris/ocean/tasks/cosine_bell/cosine_bell.cfg
+++ b/polaris/ocean/tasks/cosine_bell/cosine_bell.cfg
@@ -20,10 +20,10 @@ partial_cell_type = None
 min_pc_fraction = 0.1
 
 
-# config options for spherical convergence tests
-[spherical_convergence]
+# config options for convergence tests
+[convergence]
 
-# Evaluation time for convergence analysis (in days)
+# Evaluation time for convergence analysis (in hours)
 convergence_eval_time = ${cosine_bell:vel_pd}
 
 # Convergence threshold below which a test fails
@@ -33,8 +33,8 @@ convergence_thresh = ${cosine_bell:convergence_thresh}
 error_type = l2
 
 
-# config options for spherical convergence tests
-[spherical_convergence_forward]
+# config options for convergence forward steps
+[convergence_forward]
 
 # time integrator: {'split_explicit', 'RK4'}
 time_integrator = RK4
@@ -42,10 +42,10 @@ time_integrator = RK4
 # RK4 time step per resolution (s/km), since dt is proportional to resolution
 rk4_dt_per_km = 3.0
 
-# Run duration in days
+# Run duration in hours
 run_duration = ${cosine_bell:vel_pd}
 
-# Output interval in days
+# Output interval in hours
 output_interval = ${cosine_bell:vel_pd}
 
 
@@ -70,8 +70,8 @@ radius = 2123666.6667
 # hill max of tracer
 psi0 = 1.0
 
-# time (days) for bell to transit equator once
-vel_pd = 24.0
+# time (hours) for bell to transit equator once
+vel_pd = 576.
 
 # convergence threshold below which the test fails
 convergence_thresh = 1.8
diff --git a/polaris/ocean/tasks/cosine_bell/forward.py b/polaris/ocean/tasks/cosine_bell/forward.py
index 5212332d0..3c2804a4d 100644
--- a/polaris/ocean/tasks/cosine_bell/forward.py
+++ b/polaris/ocean/tasks/cosine_bell/forward.py
@@ -7,7 +7,7 @@ class Forward(SphericalConvergenceForward):
     bell test case
     """
 
-    def __init__(self, component, name, subdir, resolution, base_mesh, init):
+    def __init__(self, component, name, subdir, resolution, mesh, init):
         """
         Create a new step
 
@@ -25,7 +25,7 @@ def __init__(self, component, name, subdir, resolution, base_mesh, init):
         resolution : float
             The resolution of the (uniform) mesh in km
 
-        base_mesh : polaris.Step
+        mesh : polaris.Step
             The base mesh step
 
         init : polaris.Step
@@ -34,7 +34,7 @@ def __init__(self, component, name, subdir, resolution, base_mesh, init):
         package = 'polaris.ocean.tasks.cosine_bell'
         validate_vars = ['normalVelocity', 'tracer1']
         super().__init__(component=component, name=name, subdir=subdir,
-                         resolution=resolution, base_mesh=base_mesh,
+                         resolution=resolution, mesh=mesh,
                          init=init, package=package,
                          yaml_filename='forward.yaml',
                          output_filename='output.nc',
diff --git a/polaris/ocean/tasks/cosine_bell/init.py b/polaris/ocean/tasks/cosine_bell/init.py
index 251918a04..4c3b26ef9 100644
--- a/polaris/ocean/tasks/cosine_bell/init.py
+++ b/polaris/ocean/tasks/cosine_bell/init.py
@@ -97,9 +97,9 @@ def run(self):
         ds['tracer3'] = ds.tracer1
 
         # Initialize velocity
-        seconds_per_day = 86400.0
+        s_per_hour = 3600.
         velocity = (2.0 * np.pi * np.cos(angleEdge) * sphere_radius *
-                    np.cos(latEdge) / (seconds_per_day * vel_pd))
+                    np.cos(latEdge) / (s_per_hour * vel_pd))
         velocity_array, _ = xr.broadcast(velocity, ds.refZMid)
         ds['normalVelocity'] = velocity_array.expand_dims(dim='Time', axis=0)
 
diff --git a/polaris/ocean/tasks/cosine_bell/viz.py b/polaris/ocean/tasks/cosine_bell/viz.py
index 65b313692..c115629ba 100644
--- a/polaris/ocean/tasks/cosine_bell/viz.py
+++ b/polaris/ocean/tasks/cosine_bell/viz.py
@@ -148,5 +148,5 @@ def run(self):
             ds_init.lon.values, ds_init.lat.values, ds_out.tracer1.values,
             out_filename='final.png', config=config,
             colormap_section='cosine_bell_viz',
-            title=f'{mesh_name} tracer after {run_duration:g} days',
+            title=f'{mesh_name} tracer after {run_duration/24.:g} days',
             plot_land=False)
diff --git a/polaris/ocean/tasks/inertial_gravity_wave/__init__.py b/polaris/ocean/tasks/inertial_gravity_wave/__init__.py
index ab61bae4c..a3023a74e 100644
--- a/polaris/ocean/tasks/inertial_gravity_wave/__init__.py
+++ b/polaris/ocean/tasks/inertial_gravity_wave/__init__.py
@@ -1,4 +1,7 @@
-from polaris import Task
+from typing import Dict
+
+from polaris import Step, Task
+from polaris.ocean.resolution import resolution_to_subdir
 from polaris.ocean.tasks.inertial_gravity_wave.analysis import Analysis
 from polaris.ocean.tasks.inertial_gravity_wave.forward import Forward
 from polaris.ocean.tasks.inertial_gravity_wave.init import Init
@@ -34,19 +37,33 @@ def __init__(self, component):
         super().__init__(component=component, name=name, subdir=subdir)
 
         self.resolutions = [200., 100., 50., 25.]
-        for res in self.resolutions:
-            self.add_step(Init(component=component, resolution=res,
-                               taskdir=self.subdir))
-            self.add_step(Forward(component=component, resolution=res,
-                                  taskdir=self.subdir))
+        analysis_dependencies: Dict[str, Dict[float, Step]] = (
+            dict(mesh=dict(), init=dict(), forward=dict()))
+        for resolution in self.resolutions:
+            mesh_name = resolution_to_subdir(resolution)
+            init_step = Init(component=component, resolution=resolution,
+                             taskdir=self.subdir)
+            self.add_step(init_step)
+            forward_step = Forward(component=component,
+                                   name=f'forward_{mesh_name}',
+                                   resolution=resolution,
+                                   subdir=f'{self.subdir}/forward/{mesh_name}',
+                                   init=init_step)
+            self.add_step(forward_step)
+            analysis_dependencies['mesh'][resolution] = init_step
+            analysis_dependencies['init'][resolution] = init_step
+            analysis_dependencies['forward'][resolution] = forward_step
 
         self.add_step(Analysis(component=component,
                                resolutions=self.resolutions,
-                               taskdir=self.subdir))
+                               subdir=f'{subdir}/analysis',
+                               dependencies=analysis_dependencies))
         self.add_step(Viz(component=component, resolutions=self.resolutions,
                           taskdir=self.subdir),
                       run_by_default=False)
 
+        self.config.add_from_package('polaris.ocean.convergence',
+                                     'convergence.cfg')
         self.config.add_from_package(
             'polaris.ocean.tasks.inertial_gravity_wave',
             'inertial_gravity_wave.cfg')
diff --git a/polaris/ocean/tasks/inertial_gravity_wave/analysis.py b/polaris/ocean/tasks/inertial_gravity_wave/analysis.py
index b5b344d11..c9383358e 100644
--- a/polaris/ocean/tasks/inertial_gravity_wave/analysis.py
+++ b/polaris/ocean/tasks/inertial_gravity_wave/analysis.py
@@ -1,18 +1,12 @@
-import datetime
-import warnings
-
-import cmocean  # noqa: F401
-import numpy as np
 import xarray as xr
 
-from polaris import Step
-from polaris.ocean.resolution import resolution_to_subdir
+from polaris.ocean.convergence import ConvergenceAnalysis
 from polaris.ocean.tasks.inertial_gravity_wave.exact_solution import (
     ExactSolution,
 )
 
 
-class Analysis(Step):
+class Analysis(ConvergenceAnalysis):
     """
     A step for analysing the output from the inertial gravity wave
     test case
@@ -22,7 +16,7 @@ class Analysis(Step):
     resolutions : list of float
         The resolutions of the meshes that have been run
     """
-    def __init__(self, component, resolutions, taskdir):
+    def __init__(self, component, resolutions, subdir, dependencies):
         """
         Create the step
 
@@ -34,54 +28,49 @@ def __init__(self, component, resolutions, taskdir):
         resolutions : list of float
             The resolutions of the meshes that have been run
 
-        taskdir : str
-            The subdirectory that the task belongs to
-        """
-        super().__init__(component=component, name='analysis', indir=taskdir)
-        self.resolutions = resolutions
-
-        for resolution in resolutions:
-            mesh_name = resolution_to_subdir(resolution)
-            self.add_input_file(
-                filename=f'init_{mesh_name}.nc',
-                target=f'../init/{mesh_name}/initial_state.nc')
-            self.add_input_file(
-                filename=f'output_{mesh_name}.nc',
-                target=f'../forward/{mesh_name}/output.nc')
+        subdir : str
+            The subdirectory that the step resides in
 
-    def run(self):
+        dependencies : dict of dict of polaris.Steps
+            The dependencies of this step
         """
-        Run this step of the test case
+        convergence_vars = [{'name': 'ssh',
+                             'title': 'SSH',
+                             'zidx': None}]
+        super().__init__(component=component, subdir=subdir,
+                         resolutions=resolutions,
+                         dependencies=dependencies,
+                         convergence_vars=convergence_vars)
+
+    def exact_solution(self, mesh_name, field_name, time, zidx=None):
         """
-        config = self.config
-        resolutions = self.resolutions
-
-        section = config['inertial_gravity_wave']
-        conv_thresh = section.getfloat('conv_thresh')
-        conv_max = section.getfloat('conv_max')
-
-        rmse = []
-        for i, res in enumerate(resolutions):
-            mesh_name = f'{res:g}km'
-            init = xr.open_dataset(f'init_{mesh_name}.nc')
-            ds = xr.open_dataset(f'output_{mesh_name}.nc')
-            exact = ExactSolution(init, config)
-
-            t0 = datetime.datetime.strptime(ds.xtime.values[0].decode(),
-                                            '%Y-%m-%d_%H:%M:%S')
-            tf = datetime.datetime.strptime(ds.xtime.values[-1].decode(),
-                                            '%Y-%m-%d_%H:%M:%S')
-            t = (tf - t0).total_seconds()
-            ssh_model = ds.ssh.values[-1, :]
-            rmse.append(np.sqrt(np.mean((ssh_model - exact.ssh(t).values)**2)))
+        Get the exact solution
 
-        p = np.polyfit(np.log10(resolutions), np.log10(rmse), 1)
-        conv = p[0]
-
-        if conv < conv_thresh:
-            raise ValueError(f'order of convergence '
-                             f' {conv} < min tolerence {conv_thresh}')
-
-        if conv > conv_max:
-            warnings.warn(f'order of convergence '
-                          f'{conv} > max tolerence {conv_max}')
+        Parameters
+        ----------
+        mesh_name : str
+            The mesh name which is the prefix for the initial condition file
+
+        field_name : str
+            The name of the variable of which we evaluate convergence
+            For the default method, we use the same convergence rate for all
+            fields
+
+        time : float
+            The time at which to evaluate the exact solution in seconds.
+            For the default method, we always use the initial state.
+
+        zidx : int, optional
+            The z-index for the vertical level at which to evaluate the exact
+            solution
+
+        Returns
+        -------
+        solution : xarray.DataArray
+            The exact solution as derived from the initial condition
+        """
+        init = xr.open_dataset(f'{mesh_name}_init.nc')
+        exact = ExactSolution(init, self.config)
+        if field_name != 'ssh':
+            raise ValueError(f'{field_name} is not currently supported')
+        return exact.ssh(time)
diff --git a/polaris/ocean/tasks/inertial_gravity_wave/exact_solution.py b/polaris/ocean/tasks/inertial_gravity_wave/exact_solution.py
index 5fe57c214..d117f1d05 100644
--- a/polaris/ocean/tasks/inertial_gravity_wave/exact_solution.py
+++ b/polaris/ocean/tasks/inertial_gravity_wave/exact_solution.py
@@ -1,4 +1,5 @@
 import numpy as np
+from mpas_tools.cime.constants import constants
 
 
 class ExactSolution():
@@ -64,7 +65,7 @@ def __init__(self, ds, config):
         npx = section.getfloat('n_wavelengths_x')
         npy = section.getfloat('n_wavelengths_y')
 
-        self.g = 9.80616
+        self.g = constants['SHR_CONST_G']
         ly = np.sqrt(3.0) / 2.0 * lx
         self.kx = npx * 2.0 * np.pi / (lx * 1e3)
         self.ky = npy * 2.0 * np.pi / (ly * 1e3)
diff --git a/polaris/ocean/tasks/inertial_gravity_wave/forward.py b/polaris/ocean/tasks/inertial_gravity_wave/forward.py
index 9727b2f78..9d50aba11 100644
--- a/polaris/ocean/tasks/inertial_gravity_wave/forward.py
+++ b/polaris/ocean/tasks/inertial_gravity_wave/forward.py
@@ -1,13 +1,10 @@
-import time
-
 import numpy as np
 
 from polaris.mesh.planar import compute_planar_hex_nx_ny
-from polaris.ocean.model import OceanModelStep
-from polaris.ocean.resolution import resolution_to_subdir
+from polaris.ocean.convergence import ConvergenceForward
 
 
-class Forward(OceanModelStep):
+class Forward(ConvergenceForward):
     """
     A step for performing forward ocean component runs as part of inertial
     gravity wave test cases.
@@ -17,8 +14,7 @@ class Forward(OceanModelStep):
     resolution : float
         The resolution of the test case in km
     """
-    def __init__(self, component, resolution, taskdir,
-                 ntasks=None, min_tasks=None, openmp_threads=1):
+    def __init__(self, component, name, resolution, subdir, init):
         """
         Create a new test case
 
@@ -30,42 +26,20 @@ def __init__(self, component, resolution, taskdir,
         resolution : km
             The resolution of the test case in km
 
-        taskdir : str
-            The subdirectory that the task belongs to
-
-        ntasks : int, optional
-            the number of tasks the step would ideally use.  If fewer tasks
-            are available on the system, the step will run on all available
-            tasks as long as this is not below ``min_tasks``
+        subdir : str
+            The subdirectory that the step belongs to
 
-        min_tasks : int, optional
-            the number of tasks the step requires.  If the system has fewer
-            than this number of tasks, the step will fail
-
-        openmp_threads : int, optional
-            the number of OpenMP threads the step will use
+        init : polaris.Step
+            The step which generates the mesh and initial condition
         """
-        self.resolution = resolution
-        mesh_name = resolution_to_subdir(resolution)
         super().__init__(component=component,
-                         name=f'forward_{mesh_name}',
-                         subdir=f'{taskdir}/forward/{mesh_name}',
-                         ntasks=ntasks, min_tasks=min_tasks,
-                         openmp_threads=openmp_threads)
-
-        self.add_input_file(filename='initial_state.nc',
-                            target=f'../../init/{mesh_name}/initial_state.nc')
-        self.add_input_file(filename='graph.info',
-                            target=f'../../init/{mesh_name}/culled_graph.info')
-
-        self.add_output_file(
-            filename='output.nc',
-            validate_vars=['layerThickness', 'normalVelocity'])
-
-        self.add_yaml_file('polaris.ocean.config',
-                           'single_layer.yaml')
-        self.add_yaml_file('polaris.ocean.tasks.inertial_gravity_wave',
-                           'forward.yaml')
+                         name=name, subdir=subdir,
+                         resolution=resolution, mesh=init, init=init,
+                         package='polaris.ocean.tasks.inertial_gravity_wave',
+                         yaml_filename='forward.yaml',
+                         graph_filename='culled_graph.info',
+                         output_filename='output.nc',
+                         validate_vars=['layerThickness', 'normalVelocity'])
 
     def compute_cell_count(self):
         """
@@ -83,25 +57,3 @@ def compute_cell_count(self):
         nx, ny = compute_planar_hex_nx_ny(lx, ly, self.resolution)
         cell_count = nx * ny
         return cell_count
-
-    def dynamic_model_config(self, at_setup):
-        """
-        Set the model time step from config options at setup and runtime
-
-        Parameters
-        ----------
-        at_setup : bool
-            Whether this method is being run during setup of the step, as
-            opposed to at runtime
-        """
-        super().dynamic_model_config(at_setup=at_setup)
-
-        # dt is proportional to resolution
-        config = self.config
-        section = config['inertial_gravity_wave']
-        dt_per_km = section.getfloat('dt_per_km')
-        dt = dt_per_km * self.resolution
-        # https://stackoverflow.com/a/1384565/7728169
-        dt_str = time.strftime('%H:%M:%S', time.gmtime(dt))
-        options = {'config_dt': dt_str}
-        self.add_model_config_options(options)
diff --git a/polaris/ocean/tasks/inertial_gravity_wave/forward.yaml b/polaris/ocean/tasks/inertial_gravity_wave/forward.yaml
index be242a139..13dae39fa 100644
--- a/polaris/ocean/tasks/inertial_gravity_wave/forward.yaml
+++ b/polaris/ocean/tasks/inertial_gravity_wave/forward.yaml
@@ -1,8 +1,9 @@
 omega:
   time_management:
-    config_run_duration: 10:00:00
+    config_run_duration: {{ run_duration }}
   time_integration:
-    config_time_integrator: RK4
+    config_dt: {{ dt }}
+    config_time_integrator: {{ time_integrator }}
   bottom_drag:
     config_bottom_drag_mode: implicit
     config_implicit_bottom_drag_type: constant
@@ -14,14 +15,14 @@ omega:
     config_disable_vel_hmix: true
   streams:
     mesh:
-      filename_template: initial_state.nc
+      filename_template: init.nc
     input:
-      filename_template: initial_state.nc
+      filename_template: init.nc
     restart: {}
     output:
       type: output
       filename_template: output.nc
-      output_interval: 0000_00:20:00
+      output_interval: {{ output_interval }}
       clobber_mode: truncate
       contents:
       - xtime
diff --git a/polaris/ocean/tasks/inertial_gravity_wave/inertial_gravity_wave.cfg b/polaris/ocean/tasks/inertial_gravity_wave/inertial_gravity_wave.cfg
index f08ad30a0..c33e82230 100644
--- a/polaris/ocean/tasks/inertial_gravity_wave/inertial_gravity_wave.cfg
+++ b/polaris/ocean/tasks/inertial_gravity_wave/inertial_gravity_wave.cfg
@@ -19,12 +19,12 @@ n_wavelengths_y = 2
 # Convergence threshold below which the test fails
 conv_thresh = 1.8
 
-# Convergence rate above which a warning is issued
-conv_max = 2.2
-
 # time step per resolution (s/km), since dt is proportional to resolution
 dt_per_km = 3.0
 
+# Run duration in hours
+run_duration = 10.0
+
 [vertical_grid]
 
 # The type of vertical grid
@@ -45,6 +45,35 @@ partial_cell_type = None
 # The minimum fraction of a layer for partial cells
 min_pc_fraction = 0.1
 
+# config options for spherical convergence tests
+[convergence]
+
+# Evaluation time for convergence analysis (in hours)
+convergence_eval_time = ${inertial_gravity_wave:run_duration}
+
+# Convergence threshold below which a test fails
+convergence_thresh = ${inertial_gravity_wave:conv_thresh}
+
+# Type of error to compute
+error_type = l2
+
+
+# config options for spherical convergence tests
+[convergence_forward]
+
+# time integrator: {'split_explicit', 'RK4'}
+time_integrator = RK4
+
+# RK4 time step per resolution (s/km), since dt is proportional to resolution
+rk4_dt_per_km = ${inertial_gravity_wave:dt_per_km}
+
+# Run duration in hours
+run_duration = ${inertial_gravity_wave:run_duration}
+
+# Output interval in hours
+output_interval = ${inertial_gravity_wave:run_duration}
+
+
 [ocean]
 # the number of cells per core to aim for
 goal_cells_per_core = 200
diff --git a/polaris/ocean/tasks/inertial_gravity_wave/viz.py b/polaris/ocean/tasks/inertial_gravity_wave/viz.py
index 73dbf27a0..703b2bdc2 100644
--- a/polaris/ocean/tasks/inertial_gravity_wave/viz.py
+++ b/polaris/ocean/tasks/inertial_gravity_wave/viz.py
@@ -53,7 +53,6 @@ def __init__(self, component, resolutions, taskdir):
                 filename=f'output_{mesh_name}.nc',
                 target=f'../forward/{mesh_name}/output.nc')
 
-        self.add_output_file('convergence.png')
         self.add_output_file('comparison.png')
 
     def run(self):
@@ -73,9 +72,10 @@ def run(self):
         rmse = []
         error_range = None
         for i, res in enumerate(resolutions):
-            ds_mesh = xr.open_dataset(f'mesh_{res}km.nc')
-            ds_init = xr.open_dataset(f'init_{res}km.nc')
-            ds = xr.open_dataset(f'output_{res}km.nc')
+            mesh_name = resolution_to_subdir(res)
+            ds_mesh = xr.open_dataset(f'mesh_{mesh_name}.nc')
+            ds_init = xr.open_dataset(f'init_{mesh_name}.nc')
+            ds = xr.open_dataset(f'output_{mesh_name}.nc')
             ds['maxLevelCell'] = ds_init.maxLevelCell
             exact = ExactSolution(ds_init, config)
 
@@ -116,24 +116,3 @@ def run(self):
                         size='large', ha='right', va='center')
 
         fig.savefig('comparison.png', bbox_inches='tight', pad_inches=0.1)
-
-        # Convergence plots
-        fig = plt.figure()
-        ax = fig.add_subplot(111)
-        p = np.polyfit(np.log10(resolutions), np.log10(rmse), 1)
-        conv = np.round(p[0], 3)
-        ax.loglog(resolutions, rmse, '-ok', label=f'numerical (order={conv})')
-
-        c = rmse[0] * 1.5 / resolutions[0]
-        order1 = c * np.power(resolutions, 1)
-        c = rmse[0] * 1.5 / resolutions[0]**2
-        order2 = c * np.power(resolutions, 2)
-
-        ax.loglog(resolutions, order1, '--k', label='first order', alpha=0.3)
-        ax.loglog(resolutions, order2, 'k', label='second order', alpha=0.3)
-        ax.set_xlabel('resolution (km)')
-        ax.set_ylabel('RMS error (m)')
-        ax.invert_xaxis()
-        ax.set_title('Error Convergence')
-        ax.legend(loc='lower left')
-        fig.savefig('convergence.png', bbox_inches='tight', pad_inches=0.1)
diff --git a/polaris/ocean/tasks/manufactured_solution/__init__.py b/polaris/ocean/tasks/manufactured_solution/__init__.py
index 056e89047..2c95b634e 100644
--- a/polaris/ocean/tasks/manufactured_solution/__init__.py
+++ b/polaris/ocean/tasks/manufactured_solution/__init__.py
@@ -1,4 +1,7 @@
-from polaris import Task
+from typing import Dict
+
+from polaris import Step, Task
+from polaris.ocean.resolution import resolution_to_subdir
 from polaris.ocean.tasks.manufactured_solution.analysis import Analysis
 from polaris.ocean.tasks.manufactured_solution.forward import Forward
 from polaris.ocean.tasks.manufactured_solution.init import Init
@@ -39,19 +42,32 @@ def __init__(self, component):
         super().__init__(component=component, name=name, subdir=subdir)
 
         self.resolutions = [200., 100., 50., 25.]
-        for res in self.resolutions:
-            self.add_step(Init(component=component, resolution=res,
-                               taskdir=self.subdir))
-            self.add_step(Forward(component=component, resolution=res,
-                                  taskdir=self.subdir))
+        analysis_dependencies: Dict[str, Dict[float, Step]] = (
+            dict(mesh=dict(), init=dict(), forward=dict()))
+        for resolution in self.resolutions:
+            mesh_name = resolution_to_subdir(resolution)
+            init_step = Init(component=component, resolution=resolution,
+                             taskdir=self.subdir)
+            self.add_step(init_step)
+            forward_step = Forward(component=component, resolution=resolution,
+                                   name=f'forward_{mesh_name}',
+                                   subdir=f'{self.subdir}/forward/{mesh_name}',
+                                   init=init_step)
+            self.add_step(forward_step)
+            analysis_dependencies['mesh'][resolution] = init_step
+            analysis_dependencies['init'][resolution] = init_step
+            analysis_dependencies['forward'][resolution] = forward_step
 
         self.add_step(Analysis(component=component,
                                resolutions=self.resolutions,
-                               taskdir=self.subdir))
+                               subdir=f'{self.subdir}/analysis',
+                               dependencies=analysis_dependencies))
         self.add_step(Viz(component=component, resolutions=self.resolutions,
                           taskdir=self.subdir),
                       run_by_default=False)
 
+        self.config.add_from_package('polaris.ocean.convergence',
+                                     'convergence.cfg')
         self.config.add_from_package(
             'polaris.ocean.tasks.manufactured_solution',
             'manufactured_solution.cfg')
diff --git a/polaris/ocean/tasks/manufactured_solution/analysis.py b/polaris/ocean/tasks/manufactured_solution/analysis.py
index a1748049e..25a17e717 100644
--- a/polaris/ocean/tasks/manufactured_solution/analysis.py
+++ b/polaris/ocean/tasks/manufactured_solution/analysis.py
@@ -1,18 +1,12 @@
-import datetime
-import warnings
-
-import cmocean  # noqa: F401
-import numpy as np
 import xarray as xr
 
-from polaris import Step
-from polaris.ocean.resolution import resolution_to_subdir
+from polaris.ocean.convergence import ConvergenceAnalysis
 from polaris.ocean.tasks.manufactured_solution.exact_solution import (
     ExactSolution,
 )
 
 
-class Analysis(Step):
+class Analysis(ConvergenceAnalysis):
     """
     A step for analysing the output from the manufactured solution
     test case
@@ -22,7 +16,7 @@ class Analysis(Step):
     resolutions : list of float
         The resolutions of the meshes that have been run
     """
-    def __init__(self, component, resolutions, taskdir):
+    def __init__(self, component, resolutions, subdir, dependencies):
         """
         Create the step
 
@@ -34,54 +28,49 @@ def __init__(self, component, resolutions, taskdir):
         resolutions : list of float
             The resolutions of the meshes that have been run
 
-        taskdir : str
-            The subdirectory that the task belongs to
-        """
-        super().__init__(component=component, name='analysis', indir=taskdir)
-        self.resolutions = resolutions
-
-        for resolution in resolutions:
-            mesh_name = resolution_to_subdir(resolution)
-            self.add_input_file(
-                filename=f'init_{mesh_name}.nc',
-                target=f'../init/{mesh_name}/initial_state.nc')
-            self.add_input_file(
-                filename=f'output_{mesh_name}.nc',
-                target=f'../forward/{mesh_name}/output.nc')
+        subdir : str
+            The subdirectory that the step resides in
 
-    def run(self):
+        dependencies : dict of dict of polaris.Steps
+            The dependencies of this step
         """
-        Run this step of the test case
+        convergence_vars = [{'name': 'ssh',
+                             'title': 'SSH',
+                             'zidx': None}]
+        super().__init__(component=component, subdir=subdir,
+                         resolutions=resolutions,
+                         dependencies=dependencies,
+                         convergence_vars=convergence_vars)
+
+    def exact_solution(self, mesh_name, field_name, time, zidx=None):
         """
-        config = self.config
-        resolutions = self.resolutions
-
-        section = config['manufactured_solution']
-        conv_thresh = section.getfloat('conv_thresh')
-        conv_max = section.getfloat('conv_max')
-
-        rmse = []
-        for i, res in enumerate(resolutions):
-            mesh_name = f'{res:g}km'
-            init = xr.open_dataset(f'init_{mesh_name}.nc')
-            ds = xr.open_dataset(f'output_{mesh_name}.nc')
-            exact = ExactSolution(config, init)
-
-            t0 = datetime.datetime.strptime(ds.xtime.values[0].decode(),
-                                            '%Y-%m-%d_%H:%M:%S')
-            tf = datetime.datetime.strptime(ds.xtime.values[-1].decode(),
-                                            '%Y-%m-%d_%H:%M:%S')
-            t = (tf - t0).total_seconds()
-            ssh_model = ds.ssh.values[-1, :]
-            rmse.append(np.sqrt(np.mean((ssh_model - exact.ssh(t).values)**2)))
+        Get the exact solution
 
-        p = np.polyfit(np.log10(resolutions), np.log10(rmse), 1)
-        conv = p[0]
-
-        if conv < conv_thresh:
-            raise ValueError(f'order of convergence '
-                             f' {conv} < min tolerence {conv_thresh}')
-
-        if conv > conv_max:
-            warnings.warn(f'order of convergence '
-                          f'{conv} > max tolerence {conv_max}')
+        Parameters
+        ----------
+        mesh_name : str
+            The mesh name which is the prefix for the initial condition file
+
+        field_name : str
+            The name of the variable of which we evaluate convergence
+            For the default method, we use the same convergence rate for all
+            fields
+
+        time : float
+            The time at which to evaluate the exact solution in seconds.
+            For the default method, we always use the initial state.
+
+        zidx : int, optional
+            The z-index for the vertical level at which to evaluate the exact
+            solution
+
+        Returns
+        -------
+        solution : xarray.DataArray
+            The exact solution as derived from the initial condition
+        """
+        init = xr.open_dataset(f'{mesh_name}_init.nc')
+        exact = ExactSolution(self.config, init)
+        if field_name != 'ssh':
+            raise ValueError(f'{field_name} is not currently supported')
+        return exact.ssh(time)
diff --git a/polaris/ocean/tasks/manufactured_solution/forward.py b/polaris/ocean/tasks/manufactured_solution/forward.py
index d8f56c92b..cd79f3e52 100644
--- a/polaris/ocean/tasks/manufactured_solution/forward.py
+++ b/polaris/ocean/tasks/manufactured_solution/forward.py
@@ -1,16 +1,13 @@
-import time
-
 import numpy as np
 
 from polaris.mesh.planar import compute_planar_hex_nx_ny
-from polaris.ocean.model import OceanModelStep
-from polaris.ocean.resolution import resolution_to_subdir
+from polaris.ocean.convergence import ConvergenceForward
 from polaris.ocean.tasks.manufactured_solution.exact_solution import (
     ExactSolution,
 )
 
 
-class Forward(OceanModelStep):
+class Forward(ConvergenceForward):
     """
     A step for performing forward ocean component runs as part of manufactured
     solution test cases.
@@ -20,8 +17,7 @@ class Forward(OceanModelStep):
     resolution : float
         The resolution of the test case in km
     """
-    def __init__(self, component, resolution, taskdir,
-                 ntasks=None, min_tasks=None, openmp_threads=1):
+    def __init__(self, component, name, resolution, subdir, init):
         """
         Create a new test case
 
@@ -33,42 +29,20 @@ def __init__(self, component, resolution, taskdir,
         resolution : km
             The resolution of the test case in km
 
-        taskdir : str
+        subdir : str
             The subdirectory that the task belongs to
 
-        ntasks : int, optional
-            the number of tasks the step would ideally use.  If fewer tasks
-            are available on the system, the step will run on all available
-            tasks as long as this is not below ``min_tasks``
-
-        min_tasks : int, optional
-            the number of tasks the step requires.  If the system has fewer
-            than this number of tasks, the step will fail
-
-        openmp_threads : int, optional
-            the number of OpenMP threads the step will use
+        init : polaris.Step
+            The step which generates the mesh and initial condition
         """
-        self.resolution = resolution
-        mesh_name = resolution_to_subdir(resolution)
         super().__init__(component=component,
-                         name=f'forward_{mesh_name}',
-                         subdir=f'{taskdir}/forward/{mesh_name}',
-                         ntasks=ntasks, min_tasks=min_tasks,
-                         openmp_threads=openmp_threads)
-
-        self.add_input_file(filename='initial_state.nc',
-                            target=f'../../init/{mesh_name}/initial_state.nc')
-        self.add_input_file(filename='graph.info',
-                            target=f'../../init/{mesh_name}/culled_graph.info')
-
-        self.add_output_file(
-            filename='output.nc',
-            validate_vars=['layerThickness', 'normalVelocity'])
-
-        self.add_yaml_file('polaris.ocean.config',
-                           'single_layer.yaml')
-        self.add_yaml_file('polaris.ocean.tasks.manufactured_solution',
-                           'forward.yaml')
+                         name=name, subdir=subdir,
+                         resolution=resolution, mesh=init, init=init,
+                         package='polaris.ocean.tasks.manufactured_solution',
+                         yaml_filename='forward.yaml',
+                         graph_filename='culled_graph.info',
+                         output_filename='output.nc',
+                         validate_vars=['layerThickness', 'normalVelocity'])
 
     def compute_cell_count(self):
         """
@@ -100,16 +74,8 @@ def dynamic_model_config(self, at_setup):
         """
         super().dynamic_model_config(at_setup=at_setup)
 
-        # dt is proportional to resolution
-        config = self.config
-        section = config['manufactured_solution']
-        dt_per_km = section.getfloat('dt_per_km')
-        dt = dt_per_km * self.resolution
-        # https://stackoverflow.com/a/1384565/7728169
-        dt_str = time.strftime('%H:%M:%S', time.gmtime(dt))
-        exact_solution = ExactSolution(config)
-        options = {'config_dt': dt_str,
-                   'config_manufactured_solution_amplitude':
+        exact_solution = ExactSolution(self.config)
+        options = {'config_manufactured_solution_amplitude':
                    exact_solution.eta0,
                    'config_manufactured_solution_wavelength_x':
                    exact_solution.lambda_x,
diff --git a/polaris/ocean/tasks/manufactured_solution/forward.yaml b/polaris/ocean/tasks/manufactured_solution/forward.yaml
index 6b3502fce..aa40ed4b1 100644
--- a/polaris/ocean/tasks/manufactured_solution/forward.yaml
+++ b/polaris/ocean/tasks/manufactured_solution/forward.yaml
@@ -1,8 +1,9 @@
 omega:
   time_management:
-    config_run_duration: 10:00:00
+    config_run_duration: {{ run_duration }}
   time_integration:
-    config_time_integrator: RK4
+    config_dt: {{ dt }}
+    config_time_integrator: {{ time_integrator }}
   bottom_drag:
     config_bottom_drag_mode: implicit
     config_implicit_bottom_drag_type: constant
@@ -13,14 +14,14 @@ omega:
     config_disable_vel_hmix: true
   streams:
     mesh:
-      filename_template: initial_state.nc
+      filename_template: init.nc
     input:
-      filename_template: initial_state.nc
+      filename_template: init.nc
     restart: {}
     output:
       type: output
       filename_template: output.nc
-      output_interval: 0000_00:20:00
+      output_interval: {{ output_interval }}
       clobber_mode: truncate
       contents:
       - xtime
diff --git a/polaris/ocean/tasks/manufactured_solution/init.py b/polaris/ocean/tasks/manufactured_solution/init.py
index 73c997d48..af6353a1d 100644
--- a/polaris/ocean/tasks/manufactured_solution/init.py
+++ b/polaris/ocean/tasks/manufactured_solution/init.py
@@ -43,6 +43,9 @@ def __init__(self, component, resolution, taskdir):
                          name=f'init_{mesh_name}',
                          subdir=f'{taskdir}/init/{mesh_name}')
         self.resolution = resolution
+        for filename in ['culled_mesh.nc', 'initial_state.nc',
+                         'culled_graph.info']:
+            self.add_output_file(filename=filename)
 
     def run(self):
         """
diff --git a/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg b/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg
index 14d42ca80..4dc4195b5 100644
--- a/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg
+++ b/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg
@@ -30,8 +30,8 @@ dt_per_km = 3.0
 # Convergence threshold below which the test fails
 conv_thresh = 1.8
 
-# Convergence rate above which a warning is issued
-conv_max = 2.2
+# Run duration in hours
+run_duration = 10.0
 
 [vertical_grid]
 
@@ -52,3 +52,31 @@ partial_cell_type = None
 
 # The minimum fraction of a layer for partial cells
 min_pc_fraction = 0.1
+
+# config options for spherical convergence tests
+[convergence]
+
+# Evaluation time for convergence analysis (in hours)
+convergence_eval_time = ${manufactured_solution:run_duration}
+
+# Convergence threshold below which a test fails
+convergence_thresh = ${manufactured_solution:conv_thresh}
+
+# Type of error to compute
+error_type = l2
+
+
+# config options for spherical convergence tests
+[convergence_forward]
+
+# time integrator: {'split_explicit', 'RK4'}
+time_integrator = RK4
+
+# RK4 time step per resolution (s/km), since dt is proportional to resolution
+rk4_dt_per_km = ${manufactured_solution:dt_per_km}
+
+# Run duration in hours
+run_duration = ${manufactured_solution:run_duration}
+
+# Output interval in hours
+output_interval = ${manufactured_solution:run_duration}
diff --git a/polaris/ocean/tasks/manufactured_solution/viz.py b/polaris/ocean/tasks/manufactured_solution/viz.py
index 5ebbbf4e0..b485334e7 100644
--- a/polaris/ocean/tasks/manufactured_solution/viz.py
+++ b/polaris/ocean/tasks/manufactured_solution/viz.py
@@ -53,7 +53,7 @@ def __init__(self, component, resolutions, taskdir):
                 filename=f'output_{mesh_name}.nc',
                 target=f'../forward/{mesh_name}/output.nc')
 
-        self.add_output_file('convergence.png')
+        self.add_output_file('comparison.png')
 
     def run(self):
         """
@@ -114,24 +114,3 @@ def run(self):
                         size='large', ha='right', va='center')
 
         fig.savefig('comparison.png', bbox_inches='tight', pad_inches=0.1)
-
-        # Convergence polts
-        fig = plt.figure()
-        ax = fig.add_subplot(111)
-        p = np.polyfit(np.log10(resolutions), np.log10(rmse), 1)
-        conv = np.round(p[0], 3)
-        ax.loglog(resolutions, rmse, '-ok', label=f'numerical (order={conv})')
-
-        c = rmse[0] * 1.5 / resolutions[0]
-        order1 = c * np.power(resolutions, 1)
-        c = rmse[0] * 1.5 / resolutions[0]**2
-        order2 = c * np.power(resolutions, 2)
-
-        ax.loglog(resolutions, order1, '--k', label='first order', alpha=0.3)
-        ax.loglog(resolutions, order2, 'k', label='second order', alpha=0.3)
-        ax.set_xlabel('resolution (km)')
-        ax.set_ylabel('RMS error (m)')
-        ax.set_title('Error Convergence')
-        ax.invert_xaxis()
-        ax.legend(loc='lower left')
-        fig.savefig('convergence.png', bbox_inches='tight', pad_inches=0.1)