diff --git a/compass/ocean/suites/wetdry.txt b/compass/ocean/suites/wetdry.txt index 42169c7f29..0a1307c10c 100644 --- a/compass/ocean/suites/wetdry.txt +++ b/compass/ocean/suites/wetdry.txt @@ -1,13 +1,21 @@ -ocean/drying_slope/250m/sigma/default -ocean/drying_slope/250m/sigma/ramp -ocean/drying_slope/250m/single_layer/default -ocean/drying_slope/250m/single_layer/ramp -ocean/drying_slope/1km/sigma/default -ocean/drying_slope/1km/sigma/ramp -ocean/drying_slope/1km/sigma/decomp -ocean/drying_slope/1km/single_layer/default -ocean/drying_slope/1km/single_layer/ramp -ocean/drying_slope/1km/single_layer/decomp +ocean/drying_slope/sigma/standard/250m/default +ocean/drying_slope/sigma/standard/1km/default +ocean/drying_slope/sigma/standard/250m/decomp +ocean/drying_slope/sigma/standard/1km/decomp +ocean/drying_slope/single_layer/standard/250m/default +ocean/drying_slope/single_layer/standard/1km/default +ocean/drying_slope/sigma/ramp/250m/default +ocean/drying_slope/sigma/ramp/1km/default +ocean/drying_slope/sigma/ramp/250m/decomp +ocean/drying_slope/sigma/ramp/1km/decomp +ocean/drying_slope/single_layer/ramp/250m/default +ocean/drying_slope/single_layer/ramp/1km/default +ocean/drying_slope/sigma/standard/250m/loglaw +ocean/drying_slope/sigma/standard/1km/loglaw +ocean/drying_slope/single_layer/standard/250m/loglaw +ocean/drying_slope/single_layer/standard/1km/loglaw +ocean/drying_slope/sigma/standard/convergence +ocean/drying_slope/sigma/ramp/convergence ocean/dam_break/40cm/default ocean/dam_break/40cm/ramp ocean/dam_break/120cm/default diff --git a/compass/ocean/tests/drying_slope/__init__.py b/compass/ocean/tests/drying_slope/__init__.py index 75b0cd8c3e..5aa6d0d8fb 100644 --- a/compass/ocean/tests/drying_slope/__init__.py +++ b/compass/ocean/tests/drying_slope/__init__.py @@ -1,7 +1,7 @@ +from compass.ocean.tests.drying_slope.convergence import Convergence from compass.ocean.tests.drying_slope.decomp import Decomp from compass.ocean.tests.drying_slope.default import Default from compass.ocean.tests.drying_slope.loglaw import LogLaw -from compass.ocean.tests.drying_slope.ramp import Ramp from compass.testgroup import TestGroup @@ -17,17 +17,22 @@ def __init__(self, mpas_core): """ super().__init__(mpas_core=mpas_core, name='drying_slope') - for resolution in [0.25, 1.]: + for method in ['standard', 'ramp']: for coord_type in ['sigma', 'single_layer']: + for resolution in [0.25, 1.]: + self.add_test_case( + Default(test_group=self, resolution=resolution, + coord_type=coord_type, method=method)) + self.add_test_case( + Decomp(test_group=self, resolution=resolution, + coord_type=coord_type, method=method)) + for coord_type in ['sigma']: self.add_test_case( - Default(test_group=self, resolution=resolution, - coord_type=coord_type)) - self.add_test_case( - Decomp(test_group=self, resolution=resolution, - coord_type=coord_type)) - self.add_test_case( - Ramp(test_group=self, resolution=resolution, - coord_type=coord_type)) + Convergence(test_group=self, + coord_type=coord_type, + method=method)) + for coord_type in ['sigma', 'single_layer']: + for resolution in [0.25, 1.]: self.add_test_case( LogLaw(test_group=self, resolution=resolution, - coord_type=coord_type)) + coord_type=coord_type, method='standard')) diff --git a/compass/ocean/tests/drying_slope/analysis.py b/compass/ocean/tests/drying_slope/analysis.py new file mode 100644 index 0000000000..e852c11f30 --- /dev/null +++ b/compass/ocean/tests/drying_slope/analysis.py @@ -0,0 +1,168 @@ +import os + +import matplotlib.pyplot as plt +import numpy as np +import pandas as pd +import xarray as xr +from scipy.interpolate import interp1d + +from compass.step import Step + + +class Analysis(Step): + """ + A step for analyzing the convergence of drying slope results and producing + a convergence plot. + + Attributes + ---------- + damping_coeff : float + The Rayleigh damping coefficient used for the forward runs + + resolutions : float + The resolution of the test case + + times : list of float + The times at which to compare to the analytical solution + """ + def __init__(self, test_case, resolutions, damping_coeff): + """ + Create a forward step + + Parameters + ---------- + test_case : compass.TestCase + The test case this step belongs to + + resolutions : list of floats + The resolution of the test case + + damping_coeff: float + the value of the rayleigh damping coefficient + """ + + super().__init__(test_case=test_case, name='analysis') + + self.damping_coeff = damping_coeff + self.resolutions = resolutions + self.times = ['0.05', '0.15', '0.25', '0.30', '0.40', '0.50'] + + for resolution in resolutions: + if resolution < 1.: + res_name = f'{int(resolution*1e3)}m' + else: + res_name = f'{int(resolution)}km' + self.add_input_file(filename=f'output_{res_name}.nc', + target=f'../forward_{res_name}/output.nc') + for time in self.times: + filename = f'r{damping_coeff}d{time}-analytical'\ + '.csv' + self.add_input_file(filename=filename, target=filename, + database='drying_slope') + + self.add_output_file(filename='convergence.png') + + def run(self): + """ + Run this step of the test case + """ + self._plot_convergence() + + def _compute_rmse(self, ds, t): + """ + Get the time step + + Parameters + ---------- + ds : xarray.Dataset + the MPAS dataset containing output from a forward step + t : float + the time to evaluate the RMSE at + + Returns + ------- + rmse : float + the root-mean-squared-error of the MPAS output relative to the + analytical solution at time t + """ + + x_exact, ssh_exact = self._exact_solution(t) + tidx = int((float(t) / 0.2 + 1e-16) * 24.0) + ds = ds.drop_vars(np.setdiff1d([j for j in ds.variables], + ['yCell', 'ssh'])) + ds = ds.isel(Time=tidx) + x_mpas = ds.yCell.values / 1000.0 + ssh_mpas = ds.ssh.values + # Interpolate mpas solution to the points at which we have an exact + # solution + idx_min = np.argwhere(x_exact - x_mpas[0] >= 0.).item(0) + idx_max = np.argwhere(x_exact - x_mpas[-1] <= 0.).item(-1) + f = interp1d(x_mpas, ssh_mpas) + ssh_mpas_interp = f(x_exact[idx_min:idx_max]) + rmse = np.sqrt(np.mean(np.square(ssh_mpas_interp - + ssh_exact[idx_min:idx_max]))) + return rmse + + def _plot_convergence(self): + """ + Plot convergence curves + """ + comparisons = [] + cases = {'standard': '../../../standard/convergence/analysis', + 'ramp': '../../../ramp/convergence/analysis'} + for case in cases: + include = True + for resolution in self.resolutions: + if resolution < 1.: + res_name = f'{int(resolution*1e3)}m' + else: + res_name = f'{int(resolution)}km' + if not os.path.exists(f'{cases[case]}/output_{res_name}.nc'): + include = False + if include: + comparisons.append(case) + + fig, ax = plt.subplots(nrows=1, ncols=1) + + max_rmse = 0 + for k, comp in enumerate(comparisons): + rmse = np.zeros((len(self.resolutions), len(self.times))) + for i, resolution in enumerate(self.resolutions): + if resolution < 1.: + res_name = f'{int(resolution*1e3)}m' + else: + res_name = f'{int(resolution)}km' + ds = xr.open_dataset( + f'{cases[comp]}/output_{res_name}.nc') + for j, time in enumerate(self.times): + rmse[i, j] = self._compute_rmse(ds, time) + rmse_tav = np.mean(rmse, axis=1) + if np.max(rmse_tav) > max_rmse: + max_rmse = np.max(rmse_tav) + ax.loglog(self.resolutions, rmse_tav, + linestyle='-', marker='o', label=comp) + + rmse_1st_order = np.zeros(len(self.resolutions)) + rmse_1st_order[0] = max_rmse + for i in range(len(self.resolutions) - 1): + rmse_1st_order[i + 1] = rmse_1st_order[i] / 2.0 + + ax.loglog(self.resolutions, np.flip(rmse_1st_order), + linestyle='-', color='k', alpha=.25, label='1st order') + + ax.set_xlabel('Cell size (km)') + ax.set_ylabel('RMSE (m)') + + ax.legend(loc='lower right') + ax.set_title('SSH convergence') + fig.tight_layout() + fig.savefig('convergence.png') + + def _exact_solution(self, time): + """ + Returns distance, ssh of the analytic solution + """ + datafile = f'./r{self.damping_coeff}d{time}-'\ + f'analytical.csv' + data = pd.read_csv(datafile, header=None) + return data[0], -data[1] diff --git a/compass/ocean/tests/drying_slope/convergence/__init__.py b/compass/ocean/tests/drying_slope/convergence/__init__.py new file mode 100644 index 0000000000..5e09e7c8f0 --- /dev/null +++ b/compass/ocean/tests/drying_slope/convergence/__init__.py @@ -0,0 +1,126 @@ +from numpy import ceil + +from compass.config import CompassConfigParser +from compass.ocean.tests.drying_slope.analysis import Analysis +from compass.ocean.tests.drying_slope.forward import Forward +from compass.ocean.tests.drying_slope.initial_state import InitialState +from compass.ocean.tests.drying_slope.viz import Viz +from compass.testcase import TestCase +from compass.validate import compare_variables + + +class Convergence(TestCase): + """ + The convergence drying_slope test case + + Attributes + ---------- + resolution : float + The resolution of the test case in km + + coord_type : str + The type of vertical coordinate (``sigma``, ``single_layer``, etc.) + + damping_coeffs: list of float + The damping coefficients at which to evaluate convergence. Must be of + length 1. + """ + + def __init__(self, test_group, coord_type, method): + """ + Create the test case + + Parameters + ---------- + test_group : compass.ocean.tests.drying_slope.DryingSlope + The test group that this test case belongs to + + coord_type : str + The type of vertical coordinate (``sigma``, ``single_layer``) + + method: str + The wetting-and-drying method (``standard``, ``ramp``) + """ + name = 'convergence' + + self.coord_type = coord_type + damping_coeffs = [0.01] + self.damping_coeffs = damping_coeffs + subdir = f'{coord_type}/{method}/{name}' + super().__init__(test_group=test_group, name=name, + subdir=subdir) + self.resolutions = None + # add the steps with default resolutions so they can be listed + config = CompassConfigParser() + config.add_from_package('compass.ocean.tests.drying_slope', + 'drying_slope.cfg') + self._setup_steps(config, subdir=subdir, method=method) + + def _setup_steps(self, config, subdir, method): + """ setup steps given resolutions """ + + default_resolutions = '0.25, 0.5, 1, 2' + + # set the default values that a user may change before setup + config.set('drying_slope_convergence', 'resolutions', + default_resolutions, + comment='a list of resolutions (km) to test') + + # get the resolutions back, perhaps with values set in the user's + # config file + resolutions = config.getlist('drying_slope_convergence', + 'resolutions', dtype=float) + + if self.resolutions is not None and self.resolutions == resolutions: + return + + # start fresh with no steps + self.steps = dict() + self.steps_to_run = list() + + self.resolutions = resolutions + section = config['drying_slope'] + ntasks_baseline = section.getint('ntasks_baseline') + min_tasks = section.getint('min_tasks') + + for resolution in self.resolutions: + + if resolution < 1.: + res_name = f'{int(resolution*1e3)}m' + else: + res_name = f'{int(resolution)}km' + init_name = f'initial_state_{res_name}' + self.add_step(InitialState(test_case=self, + name=init_name, + resolution=resolution, + coord_type=self.coord_type)) + ntasks = max(min_tasks, + int(ceil(ntasks_baseline / resolution**2.))) + forward_step = Forward(test_case=self, resolution=resolution, + name=f'forward_{res_name}', + input_path=f'../{init_name}', + ntasks=ntasks, min_tasks=min_tasks, + openmp_threads=1, + damping_coeff=self.damping_coeffs[0], + coord_type=self.coord_type) + if method == 'ramp': + forward_step.add_namelist_options( + {'config_zero_drying_velocity_ramp': ".true."}) + self.add_step(forward_step) + self.add_step(Analysis(test_case=self, + resolutions=resolutions, + damping_coeff=self.damping_coeffs[0])) + + def validate(self): + """ + Validate variables against a baseline + """ + super().validate() + variables = ['layerThickness', 'normalVelocity'] + for resolution in self.resolutions: + if resolution < 1.: + res_name = f'{int(resolution*1e3)}m' + else: + res_name = f'{int(resolution)}km' + compare_variables(test_case=self, variables=variables, + filename1=f'forward_{res_name}/output.nc') diff --git a/compass/ocean/tests/drying_slope/decomp/__init__.py b/compass/ocean/tests/drying_slope/decomp/__init__.py index aa969f4cd5..5b99b7c170 100644 --- a/compass/ocean/tests/drying_slope/decomp/__init__.py +++ b/compass/ocean/tests/drying_slope/decomp/__init__.py @@ -1,4 +1,3 @@ -from compass.ocean.tests import drying_slope from compass.ocean.tests.drying_slope.forward import Forward from compass.ocean.tests.drying_slope.initial_state import InitialState from compass.testcase import TestCase @@ -14,19 +13,28 @@ class Decomp(TestCase): ---------- resolution : str The resolution of the test case + + coord_type : str + The type of vertical coordinate (``sigma``, ``single_layer``, etc.) """ - def __init__(self, test_group, resolution, coord_type): + def __init__(self, test_group, resolution, coord_type, method): """ Create the test case Parameters ---------- - test_group : compass.ocean.tests.baroclinic_channel.BaroclinicChannel + test_group : compass.ocean.tests.drying_slope.DryingSlope The test group that this test case belongs to - resolution : str - The resolution of the test case + resolution : float + The resolution of the test case in km + + coord_type : str + The type of vertical coordinate (``sigma``, ``single_layer``) + + method : str + The type of wetting-and-drying algorithm """ name = 'decomp' self.resolution = resolution @@ -35,11 +43,12 @@ def __init__(self, test_group, resolution, coord_type): res_name = f'{int(resolution*1e3)}m' else: res_name = f'{int(resolution)}km' - subdir = f'{res_name}/{coord_type}/{name}' + subdir = f'{coord_type}/{method}/{res_name}/{name}' super().__init__(test_group=test_group, name=name, subdir=subdir) + self.add_step(InitialState(test_case=self, resolution=resolution, + coord_type=coord_type)) - self.add_step(InitialState(test_case=self, coord_type=coord_type)) if coord_type == 'single_layer': damping_coeff = None else: @@ -51,27 +60,11 @@ def __init__(self, test_group, resolution, coord_type): ntasks=procs, openmp_threads=1, damping_coeff=damping_coeff, coord_type=coord_type) + if method == 'ramp': + forward_step.add_namelist_options( + {'config_zero_drying_velocity_ramp': ".true."}) self.add_step(forward_step) - def configure(self): - """ - Modify the configuration options for this test case. - """ - - resolution = self.resolution - config = self.config - ny = round(28 / resolution) - if resolution < 1.: - ny += 2 - dc = 1e3 * resolution - - config.set('drying_slope', 'ny', f'{ny}', comment='the number of ' - 'mesh cells in the y direction') - config.set('drying_slope', 'dc', f'{dc}', comment='the distance ' - 'between adjacent cell centers') - - # no run() method is needed - def validate(self): """ Test cases can override this method to perform validation of variables diff --git a/compass/ocean/tests/drying_slope/default/__init__.py b/compass/ocean/tests/drying_slope/default/__init__.py index 977af6920b..3af74a6bba 100644 --- a/compass/ocean/tests/drying_slope/default/__init__.py +++ b/compass/ocean/tests/drying_slope/default/__init__.py @@ -1,7 +1,10 @@ -from compass.testcase import TestCase -from compass.ocean.tests.drying_slope.initial_state import InitialState +from numpy import ceil + +from compass.config import CompassConfigParser from compass.ocean.tests.drying_slope.forward import Forward +from compass.ocean.tests.drying_slope.initial_state import InitialState from compass.ocean.tests.drying_slope.viz import Viz +from compass.testcase import TestCase from compass.validate import compare_variables @@ -18,7 +21,7 @@ class Default(TestCase): The type of vertical coordinate (``sigma``, ``single_layer``, etc.) """ - def __init__(self, test_group, resolution, coord_type): + def __init__(self, test_group, resolution, coord_type, method): """ Create the test case @@ -32,6 +35,9 @@ def __init__(self, test_group, resolution, coord_type): coord_type : str The type of vertical coordinate (``sigma``, ``single_layer``) + + method : str + The type of wetting-and-drying algorithm """ name = 'default' @@ -41,42 +47,44 @@ def __init__(self, test_group, resolution, coord_type): res_name = f'{int(resolution*1e3)}m' else: res_name = f'{int(resolution)}km' - subdir = f'{res_name}/{coord_type}/{name}' + subdir = f'{coord_type}/{method}/{res_name}/{name}' super().__init__(test_group=test_group, name=name, subdir=subdir) - self.add_step(InitialState(test_case=self, coord_type=coord_type)) + self.add_step(InitialState(test_case=self, resolution=resolution, + coord_type=coord_type)) + damping_coeffs = None + config = CompassConfigParser() + config.add_from_package('compass.ocean.tests.drying_slope', + 'drying_slope.cfg') + section = config['drying_slope'] + ntasks_baseline = section.getint('ntasks_baseline') + min_tasks = section.getint('min_tasks') + ntasks = max(min_tasks, int(ceil(ntasks_baseline / resolution**2.))) if coord_type == 'single_layer': - self.add_step(Forward(test_case=self, resolution=resolution, - ntasks=4, openmp_threads=1, - coord_type=coord_type)) - damping_coeffs = None + forward_step = Forward(test_case=self, resolution=resolution, + ntasks=ntasks, min_tasks=min_tasks, + openmp_threads=1, + coord_type=coord_type) + if method == 'ramp': + forward_step.add_namelist_options( + {'config_zero_drying_velocity_ramp': ".true."}) + self.add_step(forward_step) else: damping_coeffs = [0.0025, 0.01] for damping_coeff in damping_coeffs: - self.add_step(Forward(test_case=self, resolution=resolution, - ntasks=4, openmp_threads=1, - damping_coeff=damping_coeff, - coord_type=coord_type)) + forward_step = Forward(test_case=self, resolution=resolution, + name=f'forward_{damping_coeff}', + ntasks=ntasks, min_tasks=min_tasks, + openmp_threads=1, + damping_coeff=damping_coeff, + coord_type=coord_type) + if method == 'ramp': + forward_step.add_namelist_options( + {'config_zero_drying_velocity_ramp': ".true."}) + self.add_step(forward_step) self.damping_coeffs = damping_coeffs self.add_step(Viz(test_case=self, damping_coeffs=damping_coeffs)) - def configure(self): - """ - Modify the configuration options for this test case. - """ - - resolution = self.resolution - config = self.config - ny = round(28 / resolution) - if resolution < 1.: - ny += 2 - dc = 1e3 * resolution - - config.set('drying_slope', 'ny', f'{ny}', comment='the number of ' - 'mesh cells in the y direction') - config.set('drying_slope', 'dc', f'{dc}', comment='the distance ' - 'between adjacent cell centers') - def validate(self): """ Validate variables against a baseline diff --git a/compass/ocean/tests/drying_slope/drying_slope.cfg b/compass/ocean/tests/drying_slope/drying_slope.cfg index 31c514c747..a9d4c58c9d 100644 --- a/compass/ocean/tests/drying_slope/drying_slope.cfg +++ b/compass/ocean/tests/drying_slope/drying_slope.cfg @@ -4,12 +4,29 @@ # Number of vertical levels vert_levels = 10 +# Thickness of each layer in the thin film region +thin_film_thickness = 1.0e-3 + # config options for drying slope test cases [drying_slope] # the number of grid cells in x nx = 6 +# time step in s per km of horizontal resolution +dt_per_km = 30 + +# Number of tasks at 1km resolution +ntasks_baseline = 4 + +# Minimum number of tasks +min_tasks = 1 + +# config options for visualizing drying slope ouptut +[drying_slope_convergence] + +resolutions = 0.25, 0.5, 1, 2 + # config options for visualizing drying slope ouptut [drying_slope_viz] diff --git a/compass/ocean/tests/drying_slope/forward.py b/compass/ocean/tests/drying_slope/forward.py index 64bd2e20ec..36c9c23527 100644 --- a/compass/ocean/tests/drying_slope/forward.py +++ b/compass/ocean/tests/drying_slope/forward.py @@ -1,3 +1,5 @@ +import time + from compass.model import run_model from compass.step import Step @@ -8,10 +10,11 @@ class Forward(Step): test cases. """ def __init__(self, test_case, resolution, name='forward', subdir=None, + input_path='../initial_state', ntasks=1, min_tasks=None, openmp_threads=1, damping_coeff=None, coord_type='sigma'): """ - Create a new test case + Create a forward step Parameters ---------- @@ -48,23 +51,17 @@ def __init__(self, test_case, resolution, name='forward', subdir=None, """ if min_tasks is None: min_tasks = ntasks - if damping_coeff is not None: - name = f'{name}_{damping_coeff}' super().__init__(test_case=test_case, name=name, subdir=subdir, ntasks=ntasks, min_tasks=min_tasks, openmp_threads=openmp_threads) + self.resolution = resolution self.add_namelist_file('compass.ocean.tests.drying_slope', 'namelist.forward') - if resolution < 1.: - res_name = f'{int(resolution*1e3)}m' - else: - res_name = f'{int(resolution)}km' - self.add_namelist_file('compass.ocean.tests.drying_slope', - f'namelist.{res_name}.forward') - self.add_namelist_file('compass.ocean.tests.drying_slope', - f'namelist.{coord_type}.forward') + if coord_type == 'single_layer' or coord_type == 'sigma': + self.add_namelist_file('compass.ocean.tests.drying_slope', + f'namelist.{coord_type}.forward') if damping_coeff is not None: # update the Rayleigh damping coeff to the requested value options = {'config_Rayleigh_damping_coeff': f'{damping_coeff}'} @@ -73,7 +70,6 @@ def __init__(self, test_case, resolution, name='forward', subdir=None, self.add_streams_file('compass.ocean.tests.drying_slope', 'streams.forward') - input_path = '../initial_state' self.add_input_file(filename='mesh.nc', target=f'{input_path}/culled_mesh.nc') @@ -96,5 +92,31 @@ def run(self): """ Run this step of the test case """ + dt = self.get_dt() + section = self.config['vertical_grid'] + thin_film_thickness = section.getfloat('thin_film_thickness') + options = {'config_dt': f"'{dt}'", + 'config_drying_min_cell_height': f'{thin_film_thickness}'} + self.update_namelist_at_runtime(options=options, + out_name='namelist.ocean') run_model(self) + + def get_dt(self): + """ + Get the time step + + Returns + ------- + dt : str + the time step in HH:MM:SS + """ + config = self.config + # dt is proportional to resolution + dt_per_km = config.getfloat('drying_slope', 'dt_per_km') + + dt = dt_per_km * self.resolution + # https://stackoverflow.com/a/1384565/7728169 + dt = time.strftime('%H:%M:%S', time.gmtime(dt)) + + return dt diff --git a/compass/ocean/tests/drying_slope/initial_state.py b/compass/ocean/tests/drying_slope/initial_state.py index ea2be441fe..97d312af42 100644 --- a/compass/ocean/tests/drying_slope/initial_state.py +++ b/compass/ocean/tests/drying_slope/initial_state.py @@ -1,12 +1,9 @@ -import xarray -import numpy - -from mpas_tools.planar_hex import make_planar_hex_mesh from mpas_tools.io import write_netcdf from mpas_tools.mesh.conversion import convert, cull +from mpas_tools.planar_hex import make_planar_hex_mesh -from compass.step import Step from compass.model import run_model +from compass.step import Step class InitialState(Step): @@ -14,7 +11,8 @@ class InitialState(Step): A step for creating a mesh and initial condition for drying slope test cases """ - def __init__(self, test_case, coord_type='sigma'): + def __init__(self, test_case, resolution, name='initial_state', + coord_type='sigma'): """ Create the step @@ -23,11 +21,13 @@ def __init__(self, test_case, coord_type='sigma'): test_case : compass.ocean.tests.drying_slope.default.Default The test case this step belongs to """ - super().__init__(test_case=test_case, name='initial_state', ntasks=1, + super().__init__(test_case=test_case, name=name, ntasks=1, min_tasks=1, openmp_threads=1) self.coord_type = coord_type + self.resolution = resolution + self.add_namelist_file('compass.ocean.tests.drying_slope', 'namelist.init', mode='init') @@ -47,22 +47,28 @@ def run(self): config = self.config logger = self.logger - config = self.config section = config['vertical_grid'] coord_type = self.coord_type + thin_film_thickness = section.getfloat('thin_film_thickness') + 1.0e-9 if coord_type == 'single_layer': - options = {'config_tidal_boundary_vert_levels': '1'} + options = {'config_tidal_boundary_vert_levels': '1', + 'config_drying_min_cell_height': + f'{thin_film_thickness}'} self.update_namelist_at_runtime(options) else: - vert_levels = section.get('vert_levels') + vert_levels = section.getint('vert_levels') options = {'config_tidal_boundary_vert_levels': f'{vert_levels}', - 'config_tidal_boundary_layer_type': f"'{coord_type}'"} + 'config_tidal_boundary_layer_type': f"'{coord_type}'", + 'config_drying_min_cell_height': + f'{thin_film_thickness}'} self.update_namelist_at_runtime(options) - section = config['drying_slope'] - nx = section.getint('nx') - ny = section.getint('ny') - dc = section.getfloat('dc') + # Determine mesh parameters + nx = config.getint('drying_slope', 'nx') + ny = round(28 / self.resolution) + if self.resolution < 1.: + ny += 2 + dc = 1e3 * self.resolution logger.info(' * Make planar hex mesh') dsMesh = make_planar_hex_mesh(nx=nx, ny=ny, dc=dc, nonperiodic_x=False, diff --git a/compass/ocean/tests/drying_slope/loglaw/__init__.py b/compass/ocean/tests/drying_slope/loglaw/__init__.py index 1e941d41d9..0cbd50b721 100644 --- a/compass/ocean/tests/drying_slope/loglaw/__init__.py +++ b/compass/ocean/tests/drying_slope/loglaw/__init__.py @@ -1,3 +1,6 @@ +from numpy import ceil + +from compass.config import CompassConfigParser from compass.ocean.tests.drying_slope.forward import Forward from compass.ocean.tests.drying_slope.initial_state import InitialState from compass.ocean.tests.drying_slope.viz import Viz @@ -18,7 +21,7 @@ class LogLaw(TestCase): The type of vertical coordinate (``sigma``, ``single_layer``, etc.) """ - def __init__(self, test_group, resolution, coord_type): + def __init__(self, test_group, resolution, coord_type, method): """ Create the test case @@ -41,35 +44,27 @@ def __init__(self, test_group, resolution, coord_type): res_name = f'{int(resolution*1e3)}m' else: res_name = f'{int(resolution)}km' - subdir = f'{res_name}/{coord_type}/{name}' + subdir = f'{coord_type}/{method}/{res_name}/{name}' super().__init__(test_group=test_group, name=name, subdir=subdir) - self.add_step(InitialState(test_case=self, coord_type=coord_type)) + self.add_step(InitialState(test_case=self, coord_type=coord_type, + resolution=resolution)) + config = CompassConfigParser() + config.add_from_package('compass.ocean.tests.drying_slope', + 'drying_slope.cfg') + section = config['drying_slope'] + ntasks_baseline = section.getint('ntasks_baseline') + min_tasks = section.getint('min_tasks') + ntasks = max(min_tasks, int(ceil(ntasks_baseline / resolution**2.))) forward_step = Forward(test_case=self, resolution=resolution, - ntasks=4, openmp_threads=1, + ntasks=ntasks, min_tasks=min_tasks, + openmp_threads=1, coord_type=coord_type) forward_step.add_namelist_options( {'config_implicit_bottom_drag_type': "'loglaw'"}) self.add_step(forward_step) self.add_step(Viz(test_case=self, damping_coeffs=None)) - def configure(self): - """ - Modify the configuration options for this test case. - """ - - resolution = self.resolution - config = self.config - ny = round(28 / resolution) - if resolution < 1.: - ny += 2 - dc = 1e3 * resolution - - config.set('drying_slope', 'ny', f'{ny}', comment='the number of ' - 'mesh cells in the y direction') - config.set('drying_slope', 'dc', f'{dc}', comment='the distance ' - 'between adjacent cell centers') - def validate(self): """ Validate variables against a baseline diff --git a/compass/ocean/tests/drying_slope/namelist.1km.forward b/compass/ocean/tests/drying_slope/namelist.1km.forward deleted file mode 100644 index ef4481cc0a..0000000000 --- a/compass/ocean/tests/drying_slope/namelist.1km.forward +++ /dev/null @@ -1 +0,0 @@ -config_dt='0000_00:00:30' diff --git a/compass/ocean/tests/drying_slope/namelist.250m.forward b/compass/ocean/tests/drying_slope/namelist.250m.forward deleted file mode 100644 index 6198aa4969..0000000000 --- a/compass/ocean/tests/drying_slope/namelist.250m.forward +++ /dev/null @@ -1 +0,0 @@ -config_dt='0000_00:00:08' diff --git a/compass/ocean/tests/drying_slope/namelist.forward b/compass/ocean/tests/drying_slope/namelist.forward index ecde58b0d0..3fcac6b181 100644 --- a/compass/ocean/tests/drying_slope/namelist.forward +++ b/compass/ocean/tests/drying_slope/namelist.forward @@ -13,7 +13,7 @@ config_vert_coord_movement='impermeable_interfaces' config_ALE_thickness_proportionality='weights_only' config_use_wetting_drying=.true. config_prevent_drying=.true. -config_drying_min_cell_height=1.0e-3 +config_zero_drying_velocity_ramp_hmax = 1e-1 config_zero_drying_velocity=.true. config_verify_not_dry=.true. config_thickness_flux_type='upwind' diff --git a/compass/ocean/tests/drying_slope/namelist.init b/compass/ocean/tests/drying_slope/namelist.init index f1ff98bfc2..a15464bd7c 100644 --- a/compass/ocean/tests/drying_slope/namelist.init +++ b/compass/ocean/tests/drying_slope/namelist.init @@ -6,6 +6,5 @@ config_tidal_boundary_left_bottom_depth=0.0 config_tidal_boundary_use_distances=.true. config_tidal_forcing_monochromatic_baseline=10.0 config_use_wetting_drying=.true. -config_drying_min_cell_height=1.000001e-3 config_use_tidal_forcing=.true. config_write_cull_cell_mask=.false. diff --git a/compass/ocean/tests/drying_slope/namelist.single_layer.forward b/compass/ocean/tests/drying_slope/namelist.single_layer.forward index 7226997da9..96506c1159 100644 --- a/compass/ocean/tests/drying_slope/namelist.single_layer.forward +++ b/compass/ocean/tests/drying_slope/namelist.single_layer.forward @@ -1,4 +1,4 @@ -config_bottom_drag_mode = 'implicit' +config_implicit_bottom_drag_type = 'constant' config_implicit_constant_bottom_drag_coeff = 3.0e-3 config_disable_thick_vadv = .true. config_disable_thick_sflux = .true. diff --git a/compass/ocean/tests/drying_slope/ramp/__init__.py b/compass/ocean/tests/drying_slope/ramp/__init__.py deleted file mode 100644 index 39f1380311..0000000000 --- a/compass/ocean/tests/drying_slope/ramp/__init__.py +++ /dev/null @@ -1,99 +0,0 @@ -from compass.ocean.tests.drying_slope.forward import Forward -from compass.ocean.tests.drying_slope.initial_state import InitialState -from compass.ocean.tests.drying_slope.viz import Viz -from compass.testcase import TestCase -from compass.validate import compare_variables - - -class Ramp(TestCase): - """ - The drying_slope test case with ramped wetting-velocity factor - - Attributes - ---------- - resolution : float - The resolution of the test case in km - - coord_type : str - The type of vertical coordinate (``sigma``, ``single_layer``, etc.) - """ - - def __init__(self, test_group, resolution, coord_type): - """ - Create the test case - - Parameters - ---------- - test_group : compass.ocean.tests.drying_slope.DryingSlope - The test group that this test case belongs to - - resolution : float - The resolution of the test case in km - - coord_type : str - The type of vertical coordinate (``sigma``, ``single_layer``) - """ - name = 'ramp' - - self.resolution = resolution - self.coord_type = coord_type - if resolution < 1.: - res_name = f'{int(resolution*1e3)}m' - else: - res_name = f'{int(resolution)}km' - subdir = f'{res_name}/{coord_type}/{name}' - super().__init__(test_group=test_group, name=name, - subdir=subdir) - self.add_step(InitialState(test_case=self, coord_type=coord_type)) - if coord_type == 'single_layer': - forward_step = Forward(test_case=self, resolution=resolution, - ntasks=4, openmp_threads=1, - coord_type=coord_type) - damping_coeffs = None - forward_step.add_namelist_options( - {'config_zero_drying_velocity_ramp': ".true."}) - self.add_step(forward_step) - else: - damping_coeffs = [0.0025, 0.01] - for damping_coeff in damping_coeffs: - forward_step = Forward(test_case=self, resolution=resolution, - ntasks=4, openmp_threads=1, - damping_coeff=damping_coeff, - coord_type=coord_type) - forward_step.add_namelist_options( - {'config_zero_drying_velocity_ramp': ".true."}) - self.add_step(forward_step) - self.damping_coeffs = damping_coeffs - self.add_step(Viz(test_case=self, damping_coeffs=damping_coeffs)) - - def configure(self): - """ - Modify the configuration options for this test case. - """ - - resolution = self.resolution - config = self.config - ny = round(28 / resolution) - if resolution < 1.: - ny += 2 - dc = 1e3 * resolution - - config.set('drying_slope', 'ny', f'{ny}', comment='the number of ' - 'mesh cells in the y direction') - config.set('drying_slope', 'dc', f'{dc}', comment='the distance ' - 'between adjacent cell centers') - - def validate(self): - """ - Validate variables against a baseline - """ - damping_coeffs = self.damping_coeffs - variables = ['layerThickness', 'normalVelocity'] - if damping_coeffs is not None: - for damping_coeff in damping_coeffs: - compare_variables(test_case=self, variables=variables, - filename1=f'forward_{damping_coeff}/' - 'output.nc') - else: - compare_variables(test_case=self, variables=variables, - filename1='forward/output.nc') diff --git a/compass/ocean/tests/drying_slope/viz/__init__.py b/compass/ocean/tests/drying_slope/viz/__init__.py index 46a96ba7a8..875859d94f 100644 --- a/compass/ocean/tests/drying_slope/viz/__init__.py +++ b/compass/ocean/tests/drying_slope/viz/__init__.py @@ -1,9 +1,10 @@ import os -import xarray -import numpy +import subprocess + import matplotlib.pyplot as plt +import numpy as np import pandas as pd -import subprocess +import xarray as xr from compass.step import Step @@ -43,6 +44,8 @@ def __init__(self, test_case, damping_coeffs=None): self.times = times self.datatypes = datatypes + self.add_input_file(filename='init.nc', + target='../initial_state/ocean.nc') if damping_coeffs is None: self.add_input_file(filename='output.nc', target='../forward/output.nc') @@ -62,13 +65,11 @@ def run(self): Run this step of the test case """ section = self.config['paths'] - datapath = section.get('ocean_database_root') section = self.config['vertical_grid'] - vert_levels = section.get('vert_levels') section = self.config['drying_slope_viz'] generate_movie = section.getboolean('generate_movie') - self._plot_ssh_validation() + self._plot_ssh_validation(times=self.times) self._plot_ssh_time_series() if generate_movie: frames_per_second = section.getint('frames_per_second') @@ -79,10 +80,16 @@ def run(self): os.makedirs(os.path.join(os.getcwd(), outFolder)) except OSError: pass - self._plot_ssh_validation_for_movie(outFolder=outFolder) + for tidx, itime in enumerate(np.linspace(0, 0.5, 5 * 12 + 1)): + self._plot_ssh_validation(times=[itime], tidx=tidx, + outFolder=outFolder) self._images_to_movies(framesPerSecond=frames_per_second, outFolder=outFolder, extension=movie_format) + def _forcing(self, t): + ssh = 10. * np.sin(t * np.pi / 12.) - 10. + return ssh + def _plot_ssh_time_series(self, outFolder='.'): """ Plot ssh forcing on the right x boundary as a function of time against @@ -91,11 +98,10 @@ def _plot_ssh_time_series(self, outFolder='.'): (2013) test case. """ colors = {'MPAS-O': 'k', 'analytical': 'b', 'ROMS': 'g'} - xSsh = numpy.linspace(0, 12.0, 100) - ySsh = 10.0*numpy.sin(xSsh*numpy.pi/12.0) - 10.0 + xSsh = np.linspace(0, 12.0, 100) + ySsh = 10.0 * np.sin(xSsh * np.pi / 12.0) - 10.0 figsize = [6.4, 4.8] - markersize = 20 damping_coeffs = self.damping_coeffs if damping_coeffs is None: @@ -108,11 +114,10 @@ def _plot_ssh_time_series(self, outFolder='.'): fig, _ = plt.subplots(nrows=naxes, ncols=1, figsize=figsize, dpi=100) for i in range(naxes): - ax = plt.subplot(naxes, 1, i+1) - ds = xarray.open_dataset(ncFilename[i]) - ssh = ds.ssh + ax = plt.subplot(naxes, 1, i + 1) + ds = xr.open_dataset(ncFilename[i]) ympas = ds.ssh.where(ds.tidalInputMask).mean('nCells').values - xmpas = numpy.linspace(0, 1.0, len(ds.xtime))*12.0 + xmpas = np.linspace(0, 1.0, len(ds.xtime)) * 12.0 ax.plot(xmpas, ympas, marker='o', label='MPAS-O forward', color=colors['MPAS-O']) ax.plot(xSsh, ySsh, lw=3, label='analytical', @@ -128,7 +133,7 @@ def _plot_ssh_time_series(self, outFolder='.'): plt.close(fig) - def _plot_ssh_validation(self, outFolder='.'): + def _plot_ssh_validation(self, times, tidx=None, outFolder='.'): """ Plot ssh as a function of along-channel distance for all times for which there is validation data @@ -136,10 +141,9 @@ def _plot_ssh_validation(self, outFolder='.'): colors = {'MPAS-O': 'k', 'analytical': 'b', 'ROMS': 'g'} locs = [7.2, 2.2, 0.2, 1.2, 4.2, 9.3] - locs = 0.92 - numpy.divide(locs, 11.) + locs = 0.92 - np.divide(locs, 11.) damping_coeffs = self.damping_coeffs - times = self.times datatypes = self.datatypes if damping_coeffs is None: @@ -152,16 +156,22 @@ def _plot_ssh_validation(self, outFolder='.'): ncFilename = [f'output_{damping_coeff}.nc' for damping_coeff in damping_coeffs] - xBed = numpy.linspace(0, 25, 100) - yBed = 10.0/25.0*xBed + ds_mesh = xr.open_dataset('init.nc') + mesh_ymean = ds_mesh.isel(Time=0).groupby('yCell').mean( + dim=xr.ALL_DIMS) + bottom_depth = mesh_ymean.bottomDepth.values + x = mesh_ymean.yCell.values / 1000.0 + + xBed = np.linspace(0, 25, 100) + yBed = 10.0 / 25.0 * xBed fig, _ = plt.subplots(nrows=naxes, ncols=1, sharex=True) for i in range(naxes): - ax = plt.subplot(naxes, 1, i+1) - ds = xarray.open_dataset(ncFilename[i]) - ds = ds.drop_vars(numpy.setdiff1d([j for j in ds.variables], - ['yCell', 'ssh'])) + ax = plt.subplot(naxes, 1, i + 1) + ds = xr.open_dataset(ncFilename[i]) + ds = ds.drop_vars(np.setdiff1d([j for j in ds.variables], + ['yCell', 'ssh'])) ax.plot(xBed, yBed, '-k', lw=3) ax.set_xlim(0, 25) @@ -176,24 +186,42 @@ def _plot_ssh_validation(self, outFolder='.'): # Plot MPAS-O data # factor of 1e- needed to account for annoying round-off issue # to get right time slices - plottime = int((float(atime)/0.2 + 1e-16)*24.0) + plottime = int((float(atime) / 0.2 + 1e-16) * 24.0) ymean = ds.isel(Time=plottime).groupby('yCell').mean( - dim=xarray.ALL_DIMS) - x = ymean.yCell.values/1000.0 + dim=xr.ALL_DIMS) y = ymean.ssh.values - mpas = ax.plot(x, -y, label='MPAS-O', color=colors['MPAS-O']) - ax.text(1, ay, atime + ' days', size=8, - transform=ax.transAxes) + ax.plot(x, -y, label='MPAS-O', color=colors['MPAS-O']) if damping_coeffs is not None: ax.text(0.5, 5, 'r = ' + str(damping_coeffs[i])) # Plot comparison data - for datatype in datatypes: - datafile = f'./r{damping_coeffs[i]}d{atime}-'\ - f'{datatype.lower()}.csv' - data = pd.read_csv(datafile, header=None) - ax.scatter(data[0], data[1], marker='.', - color=colors[datatype], label=datatype) + if tidx is not None: + plt.title(f'{atime:03f} days') + for atime, ay in zip(self.times, locs): + ax.text(1, ay, f'{atime} days', size=8, + transform=ax.transAxes) + for datatype in datatypes: + datafile = f'./r{damping_coeffs[i]}d{atime}-'\ + f'{datatype.lower()}.csv' + if os.path.exists(datafile): + data = pd.read_csv(datafile, header=None) + ax.scatter(data[0], data[1], marker='.', + color=colors[datatype], + label=datatype) + else: + ax.text(1, ay, f'{atime} days', size=8, + transform=ax.transAxes) + for datatype in datatypes: + datafile = f'./r{damping_coeffs[i]}d{atime}-'\ + f'{datatype.lower()}.csv' + if os.path.exists(datafile): + data = pd.read_csv(datafile, header=None) + ax.scatter(data[0], data[1], marker='.', + color=colors[datatype], + label=datatype) + # Plot bottom depth, but line will not be visible unless bottom + # depth is incorrect + ax.plot(x, bottom_depth, ':k') ax.legend(frameon=False, loc='lower left') ds.close() @@ -206,101 +234,12 @@ def _plot_ssh_validation(self, outFolder='.'): rotation='vertical') fig.text(0.5, 0.02, 'Along channel distance (km)', ha='center') - fig.savefig(f'{outFolder}/ssh_depth_section.png', dpi=200) + filename = f'{outFolder}/ssh_depth_section' + if tidx is not None: + filename = f'{filename}_t{tidx:03d}' + fig.savefig(f'{filename}.png', dpi=200, format='png') plt.close(fig) - def _plot_ssh_validation_for_movie(self, outFolder='.'): - """ - Compare ssh along the channel at different time slices with the - analytical solution and ROMS results. - - Parameters - ---------- - - """ - colors = {'MPAS-O': 'k', 'analytical': 'b', 'ROMS': 'g'} - - locs = [7.2, 2.2, 0.2, 1.2, 4.2, 9.3] - locs = 0.92 - numpy.divide(locs, 11.) - - damping_coeffs = self.damping_coeffs - if damping_coeffs is None: - naxes = 1 - nhandles = 1 - ncFilename = ['output.nc'] - else: - naxes = len(damping_coeffs) - nhandles = naxes + 2 - ncFilename = [f'output_{damping_coeff}.nc' - for damping_coeff in damping_coeffs] - - times = self.times - datatypes = self.datatypes - - xBed = numpy.linspace(0, 25, 100) - yBed = 10.0/25.0*xBed - - ii = 0 - # Plot profiles over the 12h simulation duration - for itime in numpy.linspace(0, 0.5, 5*12+1): - - plottime = int((float(itime)/0.2 + 1e-16)*24.0) - - fig, _ = plt.subplots(nrows=naxes, ncols=1, sharex=True) - - for i in range(naxes): - ax = plt.subplot(naxes, 1, i+1) - ds = xarray.open_dataset(ncFilename[i]) - ds = ds.drop_vars(numpy.setdiff1d([j for j in ds.variables], - ['yCell', 'ssh'])) - - # Plot MPAS-O snapshots - # factor of 1e- needed to account for annoying round-off issue - # to get right time slices - ymean = ds.isel(Time=plottime).groupby('yCell').mean( - dim=xarray.ALL_DIMS) - x = ymean.yCell.values/1000.0 - y = ymean.ssh.values - ax.plot(xBed, yBed, '-k', lw=3) - mpas = ax.plot(x, -y, label='MPAS-O', color=colors['MPAS-O']) - - ax.set_ylim(-1, 11) - ax.set_xlim(0, 25) - ax.invert_yaxis() - ax.spines['top'].set_visible(False) - ax.spines['right'].set_visible(False) - ax.legend(frameon=False, loc='lower left') - ax.set_title(f't = {itime:.3f} days') - if damping_coeffs is not None: - ax.text(0.5, 5, 'r = ' + str(damping_coeffs[i])) - # Plot comparison data - for atime, ay in zip(times, locs): - ax.text(1, ay, f'{atime} days', size=8, - transform=ax.transAxes) - - for datatype in datatypes: - datafile = f'./r{damping_coeffs[i]}d{atime}-'\ - f'{datatype.lower()}.csv' - data = pd.read_csv(datafile, header=None) - ax.scatter(data[0], data[1], marker='.', - color=colors[datatype], label=datatype) - - h, l0 = ax.get_legend_handles_labels() - ax.legend(h[0:nhandles], l0[0:nhandles], frameon=False, - loc='lower left') - ax.set_title(f't = {itime:.3f} days') - - ds.close() - - fig.text(0.04, 0.5, 'Channel depth (m)', va='center', - rotation='vertical') - fig.text(0.5, 0.02, 'Along channel distance (km)', ha='center') - - fig.savefig(f'{outFolder}/ssh_depth_section_{ii:03d}.png', dpi=200) - - plt.close(fig) - ii += 1 - def _images_to_movies(self, outFolder='.', framesPerSecond=30, extension='mp4', overwrite=True): """ @@ -316,7 +255,7 @@ def _images_to_movies(self, outFolder='.', framesPerSecond=30, outFileName = f'{outFolder}/{prefix}.{extension}' if overwrite or not os.path.exists(outFileName): - imageFileTemplate = f'{outFolder}/{prefix}_%03d.png' + imageFileTemplate = f'{outFolder}/{prefix}_t%03d.png' logFileName = f'{outFolder}/logs/{prefix}.log' with open(logFileName, 'w') as logFile: args = ['ffmpeg', '-y', '-r', framesPerSecond, diff --git a/docs/developers_guide/ocean/api.rst b/docs/developers_guide/ocean/api.rst index 06ad6950b3..b6926e0e8b 100644 --- a/docs/developers_guide/ocean/api.rst +++ b/docs/developers_guide/ocean/api.rst @@ -87,12 +87,16 @@ drying_slope DryingSlope + analysis.Analysis + analysis.Analysis.run + + convergence.Convergence + convergence.Convergence.validate + decomp.Decomp - decomp.Decomp.configure decomp.Decomp.validate default.Default - default.Default.configure default.Default.validate forward.Forward @@ -105,14 +109,9 @@ drying_slope loglaw.LogLaw.configure loglaw.LogLaw.validate - ramp.Ramp - ramp.Ramp.configure - ramp.Ramp.validate - viz.Viz viz.Viz.run - global_convergence ~~~~~~~~~~~~~~~~~~ diff --git a/docs/developers_guide/ocean/test_groups/drying_slope.rst b/docs/developers_guide/ocean/test_groups/drying_slope.rst index 5b58fbdbba..aef6c4cd42 100644 --- a/docs/developers_guide/ocean/test_groups/drying_slope.rst +++ b/docs/developers_guide/ocean/test_groups/drying_slope.rst @@ -22,8 +22,8 @@ duration, bottom drag, and tidal forcing options, as well as shared ``streams.init`` and ``streams.forward`` files that defines ``mesh``, ``input``, ``restart``, ``forcing`` and ``output`` streams. -Namelist options specific to resolutions and vertical coordinates are given in -``namelist.${RES}.*`` and ``namelist.${COORD}*`` files. +Namelist options specific vertical coordinates are given in +``namelist.${COORD}*`` files. initial_state ~~~~~~~~~~~~~ @@ -66,6 +66,13 @@ between the analytical solution, MPAS-Ocean and ROMS. Similar plots are used to create a movie showing the solution from MPAS-Ocean at more fine-grained time intervals. +analysis +~~~~~~~~ + +The class :py:class:`compass.ocean.tests.drying_slope.analysis.Analysis` +produces a convergence plot for a series of forward steps at different +resolutions. It uses the analytical solution available at 5 discrete times t +compute the RMSE. .. _dev_ocean_drying_slope_default: @@ -81,6 +88,14 @@ two cases at different values of ``config_Rayleigh_damping_coeff``, 0.0025 and 0.01, for which there is comparison data. The ``single_layer`` case runs at one value of the implicit bottom drag coefficient. +.. _dev_ocean_drying_slope_convergence: + +convergence +----------- + +The :py:class:`compass.ocean.tests.drying_slope.convergence.Convergence` expands +on the default class to include initial and forward steps for multiple +resolutions and an analysis step to generate a convergence plot. .. _dev_ocean_drying_slope_decomp: @@ -105,8 +120,6 @@ ramp The :py:class:`compass.ocean.tests.drying_slope.ramp.Ramp` is identical to the default class except it sets ``ramp`` to ``True`` for the forward step to enable the ramp feature for wetting and drying. - - .. _dev_ocean_drying_slope_log_law: loglaw diff --git a/docs/users_guide/ocean/test_groups/drying_slope.rst b/docs/users_guide/ocean/test_groups/drying_slope.rst index 37c40b12f5..b3723053f2 100644 --- a/docs/users_guide/ocean/test_groups/drying_slope.rst +++ b/docs/users_guide/ocean/test_groups/drying_slope.rst @@ -64,6 +64,14 @@ The config options for this test case are: # the number of grid cells in x and y nx = 6 + # time step in s per km of horizontal resolution + dt_per_km = 30 + + # config options for visualizing drying slope ouptut + [drying_slope_convergence] + + resolutions = 0.25, 0.5, 1, 2 + # config options for visualizing drying slope ouptut [drying_slope_viz] @@ -81,12 +89,24 @@ All units are mks. default ------- -``ocean/drying_slope/${RES}/${COORD}/default`` is the default version of the -drying slope test case for two short (12h) test runs with two different drag -coefficients and validation of sea surface height through visual inspection +``ocean/drying_slope/${COORD}/${METHOD}/${RES}/default`` is the default version +of the drying slope test case for two short (12h) test runs with two different +drag coefficients and validation of sea surface height through visual inspection against analytic and ROMS solutions. ``RES`` is either 250m or 1km. ``COORD`` -is either ``single_layer`` or ``sigma``. Rayleigh drag is not compatible with -``single_layer`` so implicit drag with a constant coefficient is used. +is either ``single_layer`` or ``sigma``. The wetting and drying ``METHOD`` is +either ``standard`` or ``ramp``. The ``ramp`` method ramps velocities and +velocity tendencies over a given layer thickness range rather than a binary +switch at the minimum thickness. + +convergence +----------- + +``ocean/drying_slope/${COORD}/${METHOD}/convergence`` is a convergence test in +horizontal resolution and time. It produces a convergence plot for the +resolutions specified in the config file. ``COORD`` is either ``single_layer`` +or ``sigma``. The wetting and drying ``METHOD`` is either ``standard`` or +``ramp``. If the other of the two methods has already been run, its convergence +curve will be included in the plot. decomp ------