From 457c77c96c4234e39dff5cc70accf50b454b4cba Mon Sep 17 00:00:00 2001 From: Carolyn Begeman Date: Mon, 7 Oct 2024 15:03:16 -0500 Subject: [PATCH 1/9] Extend convergence framework to space and time --- polaris/ocean/convergence/__init__.py | 65 ++++++- polaris/ocean/convergence/analysis.py | 119 +++++++----- polaris/ocean/convergence/convergence.cfg | 11 ++ polaris/ocean/convergence/forward.py | 37 ++-- .../ocean/convergence/spherical/forward.py | 8 +- .../ocean/convergence/spherical/spherical.cfg | 12 +- polaris/ocean/tasks/cosine_bell/__init__.py | 169 +++++++++++------- polaris/ocean/tasks/cosine_bell/analysis.py | 17 +- polaris/ocean/tasks/cosine_bell/forward.py | 13 +- polaris/ocean/tasks/geostrophic/__init__.py | 163 ++++++++++------- polaris/ocean/tasks/geostrophic/analysis.py | 27 ++- polaris/ocean/tasks/geostrophic/forward.py | 13 +- .../tasks/inertial_gravity_wave/__init__.py | 106 ++++++++--- .../tasks/inertial_gravity_wave/analysis.py | 12 +- .../tasks/inertial_gravity_wave/forward.py | 14 +- .../inertial_gravity_wave.cfg | 2 + .../ocean/tasks/inertial_gravity_wave/init.py | 6 +- .../tasks/manufactured_solution/__init__.py | 107 ++++++++--- .../tasks/manufactured_solution/analysis.py | 12 +- .../tasks/manufactured_solution/forward.py | 13 +- .../ocean/tasks/manufactured_solution/init.py | 4 +- .../manufactured_solution.cfg | 2 + .../ocean/tasks/sphere_transport/__init__.py | 120 ++++++++----- .../ocean/tasks/sphere_transport/analysis.py | 13 +- .../ocean/tasks/sphere_transport/forward.py | 13 +- 25 files changed, 693 insertions(+), 385 deletions(-) diff --git a/polaris/ocean/convergence/__init__.py b/polaris/ocean/convergence/__init__.py index 7fb71d53f..77971171d 100644 --- a/polaris/ocean/convergence/__init__.py +++ b/polaris/ocean/convergence/__init__.py @@ -1,2 +1,63 @@ -from polaris.ocean.convergence.analysis import ConvergenceAnalysis -from polaris.ocean.convergence.forward import ConvergenceForward +import numpy as np + + +def get_resolution_for_task(config, refinement_factor, + refinement='both'): + + base_resolution = config.getfloat('convergence', 'base_resolution') + refinement_factors = config.getlist('convergence', + 'refinement_factors', + dtype=float) + + if refinement_factor not in refinement_factors: + raise ValueError( + 'refinement_factor not found in config option refinement_factors') + + if refinement == 'time': + resolution = base_resolution + else: + resolution = refinement_factor * base_resolution + + return resolution + + +def get_timestep_for_task(config, refinement_factor, + refinement='both'): + + base_resolution = config.getfloat('convergence', 'base_resolution') + refinement_factors = config.getlist('convergence', + 'refinement_factors', + dtype=float) + + if refinement_factor not in refinement_factors: + raise ValueError( + 'refinement_factor not found in config option refinement_factors') + + resolution = get_resolution_for_task( + config, refinement_factor, refinement=refinement) + + 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') + btr_timestep = 0. + else: + dt_per_km = section.getfloat('split_dt_per_km') + btr_dt_per_km = section.getfloat('btr_dt_per_km') + for idx, refinement_factor in enumerate(refinement_factors): + if refinement == 'time': + btr_timestep = btr_dt_per_km * base_resolution * \ + refinement_factor + elif refinement == 'space': + btr_timestep = btr_dt_per_km * base_resolution + else: + btr_timestep = btr_dt_per_km * resolution + if refinement == 'time': + timestep = dt_per_km * refinement_factor * base_resolution + elif refinement == 'space': + timestep = dt_per_km * base_resolution + else: + timestep = dt_per_km * resolution + + return timestep, btr_timestep diff --git a/polaris/ocean/convergence/analysis.py b/polaris/ocean/convergence/analysis.py index f2bbf85a1..0586fceec 100644 --- a/polaris/ocean/convergence/analysis.py +++ b/polaris/ocean/convergence/analysis.py @@ -6,7 +6,10 @@ from polaris.mpas import area_for_field, time_index_from_xtime from polaris.ocean.model import OceanIOStep -from polaris.ocean.resolution import resolution_to_subdir +from polaris.ocean.convergence import ( + get_resolution_for_task, + get_timestep_for_task, +) from polaris.viz import use_mplstyle @@ -53,8 +56,8 @@ class ConvergenceAnalysis(OceanIOStep): 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, subdir, dependencies, - convergence_vars): + def __init__(self, component, subdir, dependencies, convergence_vars, + refinement='both'): """ Create the step @@ -63,9 +66,6 @@ def __init__(self, component, resolutions, subdir, dependencies, component : polaris.Component The component the step belongs to - resolutions : list of float - The resolutions of the meshes that have been run - subdir : str The subdirectory that the step resides in @@ -102,11 +102,15 @@ def __init__(self, component, resolutions, subdir, dependencies, zidx : int The z-index to use for variables that have an nVertLevels dimension, which should be None for variables that don't + + refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ super().__init__(component=component, name='analysis', subdir=subdir) - self.resolutions = resolutions self.dependencies_dict = dependencies self.convergence_vars = convergence_vars + self.refinement = refinement for var in convergence_vars: self.add_output_file(f'convergence_{var["name"]}.png') @@ -116,21 +120,23 @@ def setup(self): Add input files based on resolutions, which may have been changed by user config options """ + config = self.config dependencies = self.dependencies_dict - - for resolution in self.resolutions: - mesh_name = resolution_to_subdir(resolution) - base_mesh = dependencies['mesh'][resolution] - init = dependencies['init'][resolution] - forward = dependencies['forward'][resolution] + refinement_factors = config.getlist('convergence', + 'refinement_factors', + dtype=float) + for refinement_factor in refinement_factors: + base_mesh = dependencies['mesh'][refinement_factor] + init = dependencies['init'][refinement_factor] + forward = dependencies['forward'][refinement_factor] self.add_input_file( - filename=f'{mesh_name}_mesh.nc', + filename=f'mesh_r{refinement_factor:02g}.nc', work_dir_target=f'{base_mesh.path}/base_mesh.nc') self.add_input_file( - filename=f'{mesh_name}_init.nc', + filename=f'init_r{refinement_factor:02g}.nc', work_dir_target=f'{init.path}/initial_state.nc') self.add_input_file( - filename=f'{mesh_name}_output.nc', + filename=f'output_r{refinement_factor:02g}.nc', work_dir_target=f'{forward.path}/output.nc') def run(self): @@ -162,37 +168,54 @@ def plot_convergence(self, variable_name, title, zidx): The z-index to use for variables that have an nVertLevels dimension, which should be None for variables that don't """ - resolutions = self.resolutions + config = self.config logger = self.logger conv_thresh, error_type = self.convergence_parameters( field_name=variable_name) error = [] - for resolution in resolutions: - mesh_name = resolution_to_subdir(resolution) + refinement_factors = config.getlist('convergence', + 'refinement_factors', + dtype=float) + resolutions = list() + timesteps = list() + for refinement_factor in refinement_factors: + timestep, _ = get_timestep_for_task( + config, refinement_factor, refinement=self.refinement) + timesteps.append(timestep) + resolution = get_resolution_for_task( + config, refinement_factor, refinement=self.refinement) + resolutions.append(resolution) error_res = self.compute_error( - mesh_name=mesh_name, + refinement_factor=refinement_factor, variable_name=variable_name, zidx=zidx, error_type=error_type) error.append(error_res) - res_array = np.array(resolutions) + if self.refinement == 'time': + refinement_array = np.array(timesteps) + x_label = 'Time (s)' + else: + refinement_array = np.array(resolutions) + x_label = 'Horizontal resolution (km)' error_array = np.array(error) filename = f'convergence_{variable_name}.csv' - data = np.stack((res_array, error_array), axis=1) + data = np.stack((refinement_array, error_array), axis=1) df = pd.DataFrame(data, columns=['resolution', error_type]) df.to_csv(f'convergence_{variable_name}.csv', index=False) convergence_failed = False - poly = np.polyfit(np.log10(res_array), np.log10(error_array), 1) + poly = np.polyfit(np.log10(refinement_array), np.log10(error_array), 1) convergence = poly[0] conv_round = np.round(convergence, 3) - fit = res_array ** poly[0] * 10 ** poly[1] + fit = refinement_array ** poly[0] * 10 ** poly[1] - order1 = 0.5 * error_array[-1] * (res_array / res_array[-1]) - order2 = 0.5 * error_array[-1] * (res_array / res_array[-1])**2 + order1 = 0.5 * error_array[-1] * \ + (refinement_array / refinement_array[-1]) + order2 = 0.5 * error_array[-1] * \ + (refinement_array / refinement_array[-1])**2 use_mplstyle() fig = plt.figure() @@ -205,9 +228,9 @@ def plot_convergence(self, variable_name, title, zidx): alpha=0.3) ax.loglog(resolutions, order2, 'k', label='second order', alpha=0.3) - ax.loglog(res_array, fit, 'k', + ax.loglog(refinement_array, fit, 'k', label=f'linear fit (order={conv_round})') - ax.loglog(res_array, error_array, 'o', label='numerical') + ax.loglog(refinement_array, error_array, 'o', label='numerical') if self.baseline_dir is not None: baseline_filename = os.path.join(self.baseline_dir, filename) @@ -217,20 +240,20 @@ def plot_convergence(self, variable_name, title, zidx): raise ValueError( f'{error_type} not available for baseline') else: - res_array = data.resolution.to_numpy() + refinement_array = data.resolution.to_numpy() error_array = data[error_type].to_numpy() poly = np.polyfit( - np.log10(res_array), np.log10(error_array), 1) + np.log10(refinement_array), np.log10(error_array), 1) base_convergence = poly[0] conv_round = np.round(base_convergence, 3) - fit = res_array ** poly[0] * 10 ** poly[1] - ax.loglog(res_array, error_array, 'o', color='#ff7f0e', - label='baseline') - ax.loglog(res_array, fit, color='#ff7f0e', + fit = refinement_array ** poly[0] * 10 ** poly[1] + ax.loglog(refinement_array, error_array, 'o', + color='#ff7f0e', label='baseline') + ax.loglog(refinement_array, fit, color='#ff7f0e', label=f'linear fit, baseline ' f'(order={conv_round})') - ax.set_xlabel('resolution (km)') + ax.set_xlabel(x_label) ax.set_ylabel(f'{error_title}') ax.set_title(f'Error Convergence of {title}') ax.legend(loc='lower left') @@ -251,15 +274,15 @@ def plot_convergence(self, variable_name, title, zidx): if convergence_failed: raise ValueError('Convergence rate below minimum tolerance.') - def compute_error(self, mesh_name, variable_name, zidx=None, + def compute_error(self, refinement_factor, variable_name, zidx=None, error_type='l2'): """ Compute the error for a given resolution Parameters ---------- - mesh_name : str - The name of the mesh + refinement_factor : float + The factor by which step is refined in space, time or both variable_name : str The variable name of which to evaluate error with respect to the @@ -277,16 +300,16 @@ def compute_error(self, mesh_name, variable_name, zidx=None, The error of the variable given by variable_name """ norm_type = {'l2': None, 'inf': np.inf} - ds_mesh = self.open_model_dataset(f'{mesh_name}_mesh.nc') + ds_mesh = self.open_model_dataset(f'mesh_r{refinement_factor:02g}.nc') config = self.config section = config['convergence'] eval_time = section.getfloat('convergence_eval_time') s_per_hour = 3600.0 - field_exact = self.exact_solution(mesh_name, variable_name, + field_exact = self.exact_solution(refinement_factor, variable_name, time=eval_time * s_per_hour, zidx=zidx) - field_mpas = self.get_output_field(mesh_name, variable_name, + field_mpas = self.get_output_field(refinement_factor, variable_name, time=eval_time * s_per_hour, zidx=zidx) diff = field_exact - field_mpas @@ -304,14 +327,14 @@ def compute_error(self, mesh_name, variable_name, zidx=None, return error - def exact_solution(self, mesh_name, field_name, time, zidx=None): + def exact_solution(self, refinement_factor, field_name, time, zidx=None): """ Get the exact solution Parameters ---------- - mesh_name : str - The mesh name which is the prefix for the initial condition file + refinement_factor : float + The factor by which step is refined in space, time or both field_name : str The name of the variable of which we evaluate convergence @@ -332,20 +355,20 @@ def exact_solution(self, mesh_name, field_name, time, zidx=None): The exact solution as derived from the initial condition """ - ds_init = self.open_model_dataset(f'{mesh_name}_init.nc') + ds_init = self.open_model_dataset(f'init_r{refinement_factor:02g}.nc') ds_init = ds_init.isel(Time=0) if zidx is not None: ds_init = ds_init.isel(nVertLevels=zidx) return ds_init[field_name] - def get_output_field(self, mesh_name, field_name, time, zidx=None): + def get_output_field(self, refinement_factor, field_name, time, zidx=None): """ Get the model output field at the given time and z index Parameters ---------- - mesh_name : str + e str The mesh name which is the prefix for the output file field_name : str @@ -362,7 +385,7 @@ def get_output_field(self, mesh_name, field_name, time, zidx=None): field_mpas : xarray.DataArray model output field """ - ds_out = self.open_model_dataset(f'{mesh_name}_output.nc') + ds_out = self.open_model_dataset(f'output_r{refinement_factor:02g}.nc') tidx = time_index_from_xtime(ds_out.xtime.values, time) ds_out = ds_out.isel(Time=tidx) diff --git a/polaris/ocean/convergence/convergence.cfg b/polaris/ocean/convergence/convergence.cfg index 5f8b0752d..f93d368c4 100644 --- a/polaris/ocean/convergence/convergence.cfg +++ b/polaris/ocean/convergence/convergence.cfg @@ -9,6 +9,15 @@ convergence_thresh = 1.0 # Type of error to compute error_type = l2 +# the base mesh resolution (km) to which refinement_factors +# are applied if refinement is 'space' or 'both' on a planar mesh +# base resolutions for spherical meshes are given in section spherical_convergence +base_resolution = 120 + +# refinement factors for a planar mesh applied to either space or time +# refinement factors for a spherical mesh given in section spherical_convergence +refinement_factors = 4., 2., 1., 0.5 + # config options for convergence forward steps [convergence_forward] @@ -18,9 +27,11 @@ error_type = l2 time_integrator = RK4 # RK4 time step per resolution (s/km), since dt is proportional to resolution +# if using convergence in time only, this is used for the largest resolution rk4_dt_per_km = 3.0 # split time step per resolution (s/km), since dt is proportional to resolution +# if using convergence in time only, this is used for the largest resolution split_dt_per_km = 30.0 # the barotropic time step (s/km) for simulations using split time stepping, diff --git a/polaris/ocean/convergence/forward.py b/polaris/ocean/convergence/forward.py index 2866853e2..710d1f3d8 100644 --- a/polaris/ocean/convergence/forward.py +++ b/polaris/ocean/convergence/forward.py @@ -1,3 +1,4 @@ +from polaris.ocean.convergence import get_timestep_for_task from polaris.ocean.model import OceanModelStep, get_time_interval_string @@ -8,9 +9,6 @@ class ConvergenceForward(OceanModelStep): Attributes ---------- - resolution : float - The resolution of the (uniform) mesh in km - package : Package The package name or module object that contains the YAML file @@ -19,10 +17,10 @@ class ConvergenceForward(OceanModelStep): """ - def __init__(self, component, name, subdir, resolution, mesh, init, - package, yaml_filename='forward.yaml', options=None, - graph_target=None, output_filename='output.nc', - validate_vars=None): + def __init__(self, component, name, subdir, refinement_factor, mesh, init, + package, yaml_filename='forward.yaml', + options=None, graph_target=None, output_filename='output.nc', + validate_vars=None, refinement='both'): """ Create a new step @@ -37,9 +35,6 @@ def __init__(self, component, name, subdir, resolution, mesh, init, 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 @@ -64,7 +59,8 @@ def __init__(self, component, name, subdir, resolution, mesh, init, super().__init__(component=component, name=name, subdir=subdir, openmp_threads=1, graph_target=graph_target) - self.resolution = resolution + self.refinement = refinement + self.refinement_factor = refinement_factor self.package = package self.yaml_filename = yaml_filename @@ -109,25 +105,16 @@ def dynamic_model_config(self, at_setup): super().dynamic_model_config(at_setup=at_setup) config = self.config + refinement_factor = self.refinement_factor 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) + timestep, btr_timestep = get_timestep_for_task( + config, refinement_factor, refinement=self.refinement) + dt_str = get_time_interval_string(seconds=timestep) + btr_dt_str = get_time_interval_string(seconds=btr_timestep) s_per_hour = 3600. run_duration = section.getfloat('run_duration') diff --git a/polaris/ocean/convergence/spherical/forward.py b/polaris/ocean/convergence/spherical/forward.py index 88be45217..26f61ac47 100644 --- a/polaris/ocean/convergence/spherical/forward.py +++ b/polaris/ocean/convergence/spherical/forward.py @@ -1,4 +1,5 @@ -from polaris.ocean.convergence import ConvergenceForward +from polaris.ocean.convergence import get_resolution_for_task +from polaris.ocean.convergence.forward import ConvergenceForward class SphericalConvergenceForward(ConvergenceForward): @@ -30,5 +31,8 @@ def compute_cell_count(self): The approximate number of cells in the mesh """ # use a heuristic based on QU30 (65275 cells) and QU240 (10383 cells) - cell_count = 6e8 / self.resolution**2 + resolution = get_resolution_for_task( + self.config, self.refinement_factor, + refinement=self.refinement) + cell_count = 6e8 / resolution**2 return cell_count diff --git a/polaris/ocean/convergence/spherical/spherical.cfg b/polaris/ocean/convergence/spherical/spherical.cfg index 69ea2dd22..0dc5fbbb3 100644 --- a/polaris/ocean/convergence/spherical/spherical.cfg +++ b/polaris/ocean/convergence/spherical/spherical.cfg @@ -1,8 +1,16 @@ # config options for spherical convergence tests [spherical_convergence] +# The base resolution for the icosahedral mesh to which the refinement +# factors are applied +icos_base_resolution = 60. + # a list of icosahedral mesh resolutions (km) to test -icos_resolutions = 60, 120, 240, 480 +icos_refinement_factors = 8., 4., 2., 1. + +# The base resolution for the quasi-uniform mesh to which the refinement +# factors are applied +qu_base_resolution = 120. # a list of quasi-uniform mesh resolutions (km) to test -qu_resolutions = 60, 90, 120, 150, 180, 210, 240 +qu_refinement_factors = 0.5, 0.75, 1., 1.25, 1.5, 1.75, 2. diff --git a/polaris/ocean/tasks/cosine_bell/__init__.py b/polaris/ocean/tasks/cosine_bell/__init__.py index 89124c644..55ca91d82 100644 --- a/polaris/ocean/tasks/cosine_bell/__init__.py +++ b/polaris/ocean/tasks/cosine_bell/__init__.py @@ -1,8 +1,14 @@ +from math import ceil from typing import Dict from polaris import Step, Task from polaris.config import PolarisConfigParser +from polaris.ocean.convergence import ( + get_resolution_for_task, + get_timestep_for_task, +) from polaris.ocean.mesh.spherical import add_spherical_base_mesh_step +from polaris.ocean.resolution import resolution_to_subdir from polaris.ocean.tasks.cosine_bell.analysis import Analysis from polaris.ocean.tasks.cosine_bell.forward import Forward from polaris.ocean.tasks.cosine_bell.init import Init @@ -29,10 +35,12 @@ def add_cosine_bell_tasks(component): 'cosine_bell.cfg') for include_viz in [False, True]: - component.add_task(CosineBell(component=component, - config=config, - icosahedral=icosahedral, - include_viz=include_viz)) + for refinement in ['space', 'time', 'both']: + component.add_task(CosineBell(component=component, + config=config, + icosahedral=icosahedral, + include_viz=include_viz, + refinement=refinement)) class CosineBell(Task): @@ -50,7 +58,8 @@ class CosineBell(Task): include_viz : bool Include VizMap and Viz steps for each resolution """ - def __init__(self, component, config, icosahedral, include_viz): + def __init__(self, component, config, icosahedral, include_viz, + refinement='both'): """ Create the convergence test @@ -68,29 +77,29 @@ def __init__(self, component, config, icosahedral, include_viz): include_viz : bool Include VizMap and Viz steps for each resolution + + refinement : str, optional """ if icosahedral: prefix = 'icos' else: prefix = 'qu' - subdir = f'spherical/{prefix}/cosine_bell' - name = f'{prefix}_cosine_bell' + subdir = f'spherical/{prefix}/cosine_bell/convergence_{refinement}' + name = f'{prefix}_cosine_bell_convergence_{refinement}' if include_viz: subdir = f'{subdir}/with_viz' name = f'{name}_with_viz' - link = 'cosine_bell.cfg' - else: - # config options live in the task already so no need for a symlink - link = None + link = 'cosine_bell.cfg' super().__init__(component=component, name=name, subdir=subdir) self.resolutions = list() + self.refinement = refinement self.icosahedral = icosahedral self.include_viz = include_viz self.set_shared_config(config, link=link) - self._setup_steps() + self._setup_steps(refinement=refinement) def configure(self): """ @@ -99,10 +108,17 @@ def configure(self): super().configure() # set up the steps again in case a user has provided new resolutions - self._setup_steps() + self._setup_steps(self.refinement) + + def _setup_steps(self, refinement): + """ + setup steps given resolutions - def _setup_steps(self): - """ setup steps given resolutions """ + Parameters + ---------- + refinement : str + refinement type. One of 'space', 'time' or 'both' + """ icosahedral = self.icosahedral config = self.config config_filename = self.config_filename @@ -112,85 +128,106 @@ def _setup_steps(self): else: prefix = 'qu' - resolutions = config.getlist('spherical_convergence', - f'{prefix}_resolutions', dtype=float) - - if self.resolutions == resolutions: - return + refinement_factors = config.getlist('spherical_convergence', + f'{prefix}_refinement_factors', + dtype=str) + refinement_factors = ', '.join(refinement_factors) + config.set('convergence', 'refinement_factors', + value=refinement_factors) + refinement_factors = config.getlist('convergence', + 'refinement_factors', + dtype=float) + + base_resolution = config.getfloat('spherical_convergence', + f'{prefix}_base_resolution') + config.set('convergence', 'base_resolution', + value=f'{base_resolution:03g}') # start fresh with no steps for step in list(self.steps.values()): self.remove_step(step) - self.resolutions = resolutions - component = self.component analysis_dependencies: Dict[str, Dict[str, Step]] = ( dict(mesh=dict(), init=dict(), forward=dict())) - for resolution in resolutions: + + resolutions = list() + timesteps = list() + + case_dir = f'spherical/{prefix}/cosine_bell' + + for refinement_factor in refinement_factors: + resolution = get_resolution_for_task( + config, refinement_factor, refinement=refinement) + base_mesh_step, mesh_name = add_spherical_base_mesh_step( component, resolution, icosahedral) - self.add_step(base_mesh_step, symlink=f'base_mesh/{mesh_name}') - analysis_dependencies['mesh'][resolution] = base_mesh_step - - cos_bell_dir = f'spherical/{prefix}/cosine_bell' + analysis_dependencies['mesh'][refinement_factor] = base_mesh_step name = f'{prefix}_init_{mesh_name}' - subdir = f'{cos_bell_dir}/init/{mesh_name}' - if self.include_viz: - symlink = f'init/{mesh_name}' - else: - symlink = None + subdir = f'{case_dir}/init/{mesh_name}' if subdir in component.steps: init_step = component.steps[subdir] else: - init_step = Init(component=component, name=name, subdir=subdir, - base_mesh=base_mesh_step) + init_step = Init(component=component, name=name, + subdir=subdir, base_mesh=base_mesh_step) init_step.set_shared_config(config, link=config_filename) - self.add_step(init_step, symlink=symlink) - analysis_dependencies['init'][resolution] = init_step + analysis_dependencies['init'][refinement_factor] = init_step - name = f'{prefix}_forward_{mesh_name}' - subdir = f'{cos_bell_dir}/forward/{mesh_name}' - if self.include_viz: - symlink = f'forward/{mesh_name}' - else: - symlink = None + if resolution not in resolutions: + self.add_step(base_mesh_step, symlink=f'base_mesh/{mesh_name}') + self.add_step(init_step, symlink=f'init/{mesh_name}') + resolutions.append(resolution) + + timestep, _ = get_timestep_for_task( + config, refinement_factor, refinement=refinement) + timestep = ceil(timestep) + timesteps.append(timestep) + + subdir = f'{case_dir}/forward/{mesh_name}_{timestep}s' + symlink = f'forward/{mesh_name}_{timestep}s' if subdir in component.steps: forward_step = component.steps[subdir] else: - forward_step = Forward(component=component, name=name, - subdir=subdir, resolution=resolution, - mesh=base_mesh_step, - init=init_step) - forward_step.set_shared_config(config, link=config_filename) + name = f'{prefix}_forward_{mesh_name}_{timestep}s' + forward_step = Forward( + component=component, name=name, + subdir=subdir, + refinement_factor=refinement_factor, + mesh=base_mesh_step, + init=init_step, refinement=refinement) + forward_step.set_shared_config( + config, link=config_filename) self.add_step(forward_step, symlink=symlink) - analysis_dependencies['forward'][resolution] = forward_step + analysis_dependencies['forward'][refinement_factor] = forward_step if self.include_viz: - with_viz_dir = f'spherical/{prefix}/cosine_bell/with_viz' - name = f'{prefix}_viz_{mesh_name}' - subdir = f'{with_viz_dir}/viz/{mesh_name}' - step = Viz(component=component, name=name, - subdir=subdir, base_mesh=base_mesh_step, - init=init_step, forward=forward_step, - mesh_name=mesh_name) - step.set_shared_config(config, link=config_filename) - self.add_step(step) - - subdir = f'spherical/{prefix}/cosine_bell/analysis' - if self.include_viz: - symlink = 'analysis' - else: - symlink = None + name = f'{prefix}_viz_{mesh_name}_{timestep}s' + subdir = f'{case_dir}/viz/{mesh_name}_{timestep}s' + if subdir in component.steps: + viz_step = component.steps[subdir] + else: + viz_step = Viz(component=component, name=name, + subdir=subdir, base_mesh=base_mesh_step, + init=init_step, forward=forward_step, + mesh_name=mesh_name) + viz_step.set_shared_config(config, link=config_filename) + self.add_step(viz_step) + + subdir = f'{case_dir}/convergence_{refinement}/analysis' if subdir in component.steps: step = component.steps[subdir] step.resolutions = resolutions step.dependencies_dict = analysis_dependencies else: - step = Analysis(component=component, resolutions=resolutions, + step = Analysis(component=component, subdir=subdir, - dependencies=analysis_dependencies) + dependencies=analysis_dependencies, + refinement=refinement) step.set_shared_config(config, link=config_filename) - self.add_step(step, symlink=symlink) + if self.include_viz: + symlink = f'analysis_{refinement}' + self.add_step(step, symlink=symlink) + else: + self.add_step(step) diff --git a/polaris/ocean/tasks/cosine_bell/analysis.py b/polaris/ocean/tasks/cosine_bell/analysis.py index 7ba98f8ae..79106a655 100644 --- a/polaris/ocean/tasks/cosine_bell/analysis.py +++ b/polaris/ocean/tasks/cosine_bell/analysis.py @@ -3,7 +3,7 @@ from mpas_tools.transects import lon_lat_to_cartesian from mpas_tools.vector import Vector -from polaris.ocean.convergence import ConvergenceAnalysis +from polaris.ocean.convergence.analysis import ConvergenceAnalysis from polaris.ocean.tasks.cosine_bell.init import cosine_bell @@ -11,7 +11,7 @@ class Analysis(ConvergenceAnalysis): """ A step for analyzing the output from the cosine bell test case """ - def __init__(self, component, resolutions, subdir, dependencies): + def __init__(self, component, subdir, dependencies, refinement='both'): """ Create the step @@ -20,9 +20,6 @@ def __init__(self, component, resolutions, subdir, dependencies): component : polaris.Component The component the step belongs to - resolutions : list of float - The resolutions of the meshes that have been run - subdir : str The subdirectory that the step resides in @@ -33,11 +30,11 @@ def __init__(self, component, resolutions, subdir, dependencies): 'title': 'tracer1', 'zidx': 0}] super().__init__(component=component, subdir=subdir, - resolutions=resolutions, dependencies=dependencies, - convergence_vars=convergence_vars) + convergence_vars=convergence_vars, + refinement=refinement) - def exact_solution(self, mesh_name, field_name, time, zidx=None): + def exact_solution(self, refinement_factor, field_name, time, zidx=None): """ Get the exact solution @@ -76,10 +73,10 @@ def exact_solution(self, mesh_name, field_name, time, zidx=None): psi0 = config.getfloat('cosine_bell', 'psi0') vel_pd = config.getfloat('cosine_bell', 'vel_pd') - ds_mesh = xr.open_dataset(f'{mesh_name}_mesh.nc') + ds_mesh = xr.open_dataset(f'mesh_r{refinement_factor:02g}.nc') sphere_radius = ds_mesh.sphere_radius - ds_init = xr.open_dataset(f'{mesh_name}_init.nc') + ds_init = xr.open_dataset(f'init_r{refinement_factor:02g}.nc') latCell = ds_init.latCell.values lonCell = ds_init.lonCell.values diff --git a/polaris/ocean/tasks/cosine_bell/forward.py b/polaris/ocean/tasks/cosine_bell/forward.py index 19ae04ee0..275b0b8a8 100644 --- a/polaris/ocean/tasks/cosine_bell/forward.py +++ b/polaris/ocean/tasks/cosine_bell/forward.py @@ -7,7 +7,8 @@ class Forward(SphericalConvergenceForward): bell test case """ - def __init__(self, component, name, subdir, resolution, mesh, init): + def __init__(self, component, name, subdir, mesh, init, + refinement_factor, refinement='both'): """ Create a new step @@ -30,13 +31,17 @@ def __init__(self, component, name, subdir, resolution, mesh, init): init : polaris.Step The init step + + time_refinement_factor : float + The amount by which to scale dt_per_km and btr_dt_per_km """ package = 'polaris.ocean.tasks.cosine_bell' validate_vars = ['normalVelocity', 'tracer1'] super().__init__(component=component, name=name, subdir=subdir, - resolution=resolution, mesh=mesh, - init=init, package=package, + mesh=mesh, init=init, package=package, yaml_filename='forward.yaml', output_filename='output.nc', validate_vars=validate_vars, - graph_target=f'{mesh.path}/graph.info') + graph_target=f'{mesh.path}/graph.info', + refinement_factor=refinement_factor, + refinement=refinement) diff --git a/polaris/ocean/tasks/geostrophic/__init__.py b/polaris/ocean/tasks/geostrophic/__init__.py index d29ad36e7..2ee90f021 100644 --- a/polaris/ocean/tasks/geostrophic/__init__.py +++ b/polaris/ocean/tasks/geostrophic/__init__.py @@ -1,7 +1,12 @@ +from math import ceil from typing import Dict from polaris import Step, Task from polaris.config import PolarisConfigParser +from polaris.ocean.convergence import ( + get_resolution_for_task, + get_timestep_for_task, +) from polaris.ocean.mesh.spherical import add_spherical_base_mesh_step from polaris.ocean.tasks.geostrophic.analysis import Analysis from polaris.ocean.tasks.geostrophic.forward import Forward @@ -29,10 +34,12 @@ def add_geostrophic_tasks(component): 'geostrophic.cfg') for include_viz in [False, True]: - component.add_task(Geostrophic(component=component, - config=config, - icosahedral=icosahedral, - include_viz=include_viz)) + for refinement in ['space', 'time', 'both']: + component.add_task(Geostrophic(component=component, + config=config, + icosahedral=icosahedral, + include_viz=include_viz, + refinement=refinement)) class Geostrophic(Task): @@ -50,7 +57,8 @@ class Geostrophic(Task): include_viz : bool Include VizMap and Viz steps for each resolution """ - def __init__(self, component, config, icosahedral, include_viz): + def __init__(self, component, config, icosahedral, include_viz, + refinement='both'): """ Create the convergence test @@ -74,23 +82,22 @@ def __init__(self, component, config, icosahedral, include_viz): else: prefix = 'qu' - subdir = f'spherical/{prefix}/geostrophic' - name = f'{prefix}_geostrophic' + subdir = f'spherical/{prefix}/geostrophic/convergence_{refinement}' + name = f'{prefix}_geostrophic_convergence_{refinement}' if include_viz: subdir = f'{subdir}/with_viz' name = f'{name}_with_viz' - link = 'geostrophic.cfg' - else: - # config options live in the task already so no need for a symlink - link = None + link = 'geostrophic.cfg' + super().__init__(component=component, name=name, subdir=subdir) self.resolutions = list() + self.refinement = refinement self.icosahedral = icosahedral self.include_viz = include_viz self.set_shared_config(config, link=link) - self._setup_steps() + self._setup_steps(refinement=refinement) def configure(self): """ @@ -99,96 +106,126 @@ def configure(self): super().configure() # set up the steps again in case a user has provided new resolutions - self._setup_steps() + self._setup_steps(self.refinement) + + def _setup_steps(self, refinement): + """ + setup steps given resolutions - def _setup_steps(self): - """ setup steps given resolutions """ + Parameters + ---------- + refinement : str + refinement type. One of 'space', 'time' or 'both' + """ icosahedral = self.icosahedral config = self.config config_filename = self.config_filename + if icosahedral: prefix = 'icos' else: prefix = 'qu' - resolutions = config.getlist('spherical_convergence', - f'{prefix}_resolutions', dtype=float) - - if self.resolutions == resolutions: - return + refinement_factors = config.getlist('spherical_convergence', + f'{prefix}_refinement_factors', + dtype=str) + refinement_factors = ', '.join(refinement_factors) + config.set('convergence', 'refinement_factors', + value=refinement_factors) + refinement_factors = config.getlist('convergence', + 'refinement_factors', + dtype=float) + + base_resolution = config.getfloat('spherical_convergence', + f'{prefix}_base_resolution') + config.set('convergence', 'base_resolution', + value=f'{base_resolution:03g}') # start fresh with no steps for step in list(self.steps.values()): self.remove_step(step) - self.resolutions = resolutions - component = self.component + resolutions = list() + timesteps = list() + case_dir = f'spherical/{prefix}/geostrophic' analysis_dependencies: Dict[str, Dict[float, Step]] = ( dict(mesh=dict(), init=dict(), forward=dict())) - for resolution in resolutions: + + for refinement_factor in refinement_factors: + resolution = get_resolution_for_task( + config, refinement_factor, refinement=refinement) + base_mesh_step, mesh_name = add_spherical_base_mesh_step( component, resolution, icosahedral) - self.add_step(base_mesh_step, symlink=f'base_mesh/{mesh_name}') - analysis_dependencies['mesh'][resolution] = base_mesh_step + analysis_dependencies['mesh'][refinement_factor] = base_mesh_step name = f'{prefix}_init_{mesh_name}' subdir = f'{case_dir}/init/{mesh_name}' - if self.include_viz: - symlink = f'init/{mesh_name}' - else: - symlink = None if subdir in component.steps: init_step = component.steps[subdir] else: - init_step = Init(component=component, name=name, subdir=subdir, - base_mesh=base_mesh_step) + init_step = Init(component=component, name=name, + subdir=subdir, base_mesh=base_mesh_step) init_step.set_shared_config(config, link=config_filename) - self.add_step(init_step, symlink=symlink) - analysis_dependencies['init'][resolution] = init_step + analysis_dependencies['init'][refinement_factor] = init_step - name = f'{prefix}_forward_{mesh_name}' - subdir = f'{case_dir}/forward/{mesh_name}' - if self.include_viz: - symlink = f'forward/{mesh_name}' - else: - symlink = None + if resolution not in resolutions: + self.add_step(base_mesh_step, symlink=f'base_mesh/{mesh_name}') + self.add_step(init_step, symlink=f'init/{mesh_name}') + resolutions.append(resolution) + + timestep, _ = get_timestep_for_task( + config, refinement_factor, refinement=refinement) + timestep = ceil(timestep) + timesteps.append(timestep) + + subdir = f'{case_dir}/forward/{mesh_name}_{timestep}s' + symlink = f'forward/{mesh_name}_{timestep}s' if subdir in component.steps: forward_step = component.steps[subdir] else: - forward_step = Forward(component=component, name=name, - subdir=subdir, resolution=resolution, - mesh=base_mesh_step, - init=init_step) - forward_step.set_shared_config(config, link=config_filename) + name = f'{prefix}_forward_{mesh_name}_{timestep}s' + forward_step = Forward( + component=component, name=name, + subdir=subdir, + refinement_factor=refinement_factor, + mesh=base_mesh_step, + init=init_step, refinement=refinement) + forward_step.set_shared_config( + config, link=config_filename) self.add_step(forward_step, symlink=symlink) - analysis_dependencies['forward'][resolution] = forward_step + analysis_dependencies['forward'][refinement_factor] = forward_step if self.include_viz: - with_viz_dir = f'{case_dir}/with_viz' - name = f'{prefix}_viz_{mesh_name}' - subdir = f'{with_viz_dir}/viz/{mesh_name}' - step = Viz(component=component, name=name, - subdir=subdir, base_mesh=base_mesh_step, - init=init_step, forward=forward_step, - mesh_name=mesh_name) - step.set_shared_config(config, link=config_filename) - self.add_step(step) - - subdir = f'{case_dir}/analysis' - if self.include_viz: - symlink = 'analysis' - else: - symlink = None + name = f'{prefix}_viz_{mesh_name}_{timestep}s' + subdir = f'{case_dir}/viz/{mesh_name}_{timestep}s' + if subdir in component.steps: + viz_step = component.steps[subdir] + else: + viz_step = Viz(component=component, name=name, + subdir=subdir, base_mesh=base_mesh_step, + init=init_step, forward=forward_step, + mesh_name=mesh_name) + viz_step.set_shared_config(config, link=config_filename) + self.add_step(viz_step) + + subdir = f'{case_dir}/convergence_{refinement}/analysis' if subdir in component.steps: step = component.steps[subdir] step.resolutions = resolutions step.dependencies_dict = analysis_dependencies else: - step = Analysis(component=component, resolutions=resolutions, - subdir=subdir, dependencies=analysis_dependencies) + step = Analysis(component=component, + subdir=subdir, + dependencies=analysis_dependencies, + refinement=refinement) step.set_shared_config(config, link=config_filename) - self.add_step(step, symlink=symlink) + if self.include_viz: + symlink = f'analysis_{refinement}' + self.add_step(step, symlink=symlink) + else: + self.add_step(step) diff --git a/polaris/ocean/tasks/geostrophic/analysis.py b/polaris/ocean/tasks/geostrophic/analysis.py index f56847cb8..e4fab12d4 100644 --- a/polaris/ocean/tasks/geostrophic/analysis.py +++ b/polaris/ocean/tasks/geostrophic/analysis.py @@ -1,6 +1,6 @@ import xarray as xr -from polaris.ocean.convergence import ConvergenceAnalysis +from polaris.ocean.convergence.analysis import ConvergenceAnalysis from polaris.ocean.tasks.geostrophic.exact_solution import ( compute_exact_solution, ) @@ -10,7 +10,7 @@ class Analysis(ConvergenceAnalysis): """ A step for analyzing the output from the geostrophic test case """ - def __init__(self, component, resolutions, subdir, dependencies): + def __init__(self, component, subdir, dependencies, refinement='both'): """ Create the step @@ -19,9 +19,6 @@ def __init__(self, component, resolutions, subdir, dependencies): component : polaris.Component The component the step belongs to - resolutions : list of float - The resolutions of the meshes that have been run - subdir : str The subdirectory that the step resides in @@ -35,9 +32,9 @@ def __init__(self, component, resolutions, subdir, dependencies): 'title': 'normal velocity', 'zidx': 0}] super().__init__(component=component, subdir=subdir, - resolutions=resolutions, dependencies=dependencies, - convergence_vars=convergence_vars) + convergence_vars=convergence_vars, + refinement=refinement) def exact_solution(self, mesh_name, field_name, time, zidx=None): """ @@ -88,7 +85,7 @@ def exact_solution(self, mesh_name, field_name, time, zidx=None): else: return normalVelocity - def get_output_field(self, mesh_name, field_name, time, zidx=None): + def get_output_field(self, refinement_factor, field_name, time, zidx=None): """ Get the model output field at the given time and z index @@ -117,15 +114,15 @@ def get_output_field(self, mesh_name, field_name, time, zidx=None): f'geostrophic test case') if field_name == 'normalVelocity': - return super().get_output_field(mesh_name=mesh_name, - field_name=field_name, - time=time, zidx=zidx) + return super().get_output_field( + refinement_factor=refinement_factor, + field_name=field_name, time=time, zidx=zidx) else: - ds_init = xr.open_dataset(f'{mesh_name}_init.nc') + ds_init = xr.open_dataset(f'init_r{refinement_factor:02g}.nc') bottom_depth = ds_init.bottomDepth - ssh = super().get_output_field(mesh_name=mesh_name, - field_name='ssh', - time=time, zidx=None) + ssh = super().get_output_field( + refinement_factor=refinement_factor, + field_name='ssh', time=time, zidx=None) h = ssh + bottom_depth return h diff --git a/polaris/ocean/tasks/geostrophic/forward.py b/polaris/ocean/tasks/geostrophic/forward.py index 1e788b468..13249a1c7 100644 --- a/polaris/ocean/tasks/geostrophic/forward.py +++ b/polaris/ocean/tasks/geostrophic/forward.py @@ -7,7 +7,8 @@ class Forward(SphericalConvergenceForward): bell test case """ - def __init__(self, component, name, subdir, resolution, mesh, init): + def __init__(self, component, name, subdir, mesh, init, + refinement_factor, refinement='both'): """ Create a new step @@ -22,9 +23,6 @@ def __init__(self, component, name, subdir, resolution, mesh, init): subdir : str The subdirectory for the step - resolution : float - The resolution of the (uniform) mesh in km - mesh : polaris.Step The base mesh step @@ -37,9 +35,10 @@ def __init__(self, component, name, subdir, resolution, mesh, init): 'layerThickness', 'normalVelocity'] super().__init__(component=component, name=name, subdir=subdir, - resolution=resolution, mesh=mesh, - init=init, package=package, + mesh=mesh, init=init, package=package, yaml_filename='forward.yaml', output_filename='output.nc', validate_vars=validate_vars, - graph_target=f'{mesh.path}/graph.info') + graph_target=f'{mesh.path}/graph.info', + refinement_factor=refinement_factor, + refinement=refinement) diff --git a/polaris/ocean/tasks/inertial_gravity_wave/__init__.py b/polaris/ocean/tasks/inertial_gravity_wave/__init__.py index a3023a74e..69f9c6a54 100644 --- a/polaris/ocean/tasks/inertial_gravity_wave/__init__.py +++ b/polaris/ocean/tasks/inertial_gravity_wave/__init__.py @@ -1,6 +1,12 @@ +from math import ceil from typing import Dict from polaris import Step, Task +from polaris.config import PolarisConfigParser +from polaris.ocean.convergence import ( + get_resolution_for_task, + get_timestep_for_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 @@ -15,7 +21,19 @@ def add_inertial_gravity_wave_tasks(component): component : polaris.ocean.Ocean the ocean component that the task will be added to """ - component.add_task(InertialGravityWave(component=component)) + basedir = 'planar/inertial_gravity_wave' + config_filename = 'inertial_gravity_wave.cfg' + filepath = f'{basedir}/{config_filename}' + config = PolarisConfigParser(filepath=filepath) + config.add_from_package('polaris.ocean.convergence', + 'convergence.cfg') + config.add_from_package( + 'polaris.ocean.tasks.inertial_gravity_wave', + config_filename) + for refinement in ['space', 'time', 'both']: + component.add_task(InertialGravityWave(component=component, + config=config, + refinement=refinement)) class InertialGravityWave(Task): @@ -23,7 +41,7 @@ class InertialGravityWave(Task): The convergence test case for inertial gravity waves """ - def __init__(self, component): + def __init__(self, component, config, refinement='both'): """ Create the test case @@ -31,39 +49,73 @@ def __init__(self, component): ---------- component : polaris.ocean.Ocean The ocean component that this task belongs to + + config : polaris.config.PolarisConfigParser + A shared config parser + + refinement : str, optional + Whether to refine in space, time or both space and time """ - name = 'inertial_gravity_wave' - subdir = f'planar/{name}' + name = f'inertial_gravity_wave_convergence_{refinement}' + basedir = 'planar/inertial_gravity_wave' + subdir = f'{basedir}/convergence_{refinement}' + config_filename = 'inertial_gravity_wave.cfg' super().__init__(component=component, name=name, subdir=subdir) + self.set_shared_config(config, link=config_filename) - self.resolutions = [200., 100., 50., 25.] analysis_dependencies: Dict[str, Dict[float, Step]] = ( dict(mesh=dict(), init=dict(), forward=dict())) - for resolution in self.resolutions: + + refinement_factors = self.config.getlist( + 'convergence', 'refinement_factors', dtype=float) + timesteps = list() + resolutions = list() + for refinement_factor in refinement_factors: + resolution = get_resolution_for_task( + self.config, refinement_factor, refinement=refinement) 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 + + subdir = f'{basedir}/init/{mesh_name}' + symlink = f'init/{mesh_name}' + if subdir in component.steps: + init_step = component.steps[subdir] + else: + init_step = Init(component=component, resolution=resolution, + subdir=subdir) + init_step.set_shared_config(self.config, link=config_filename) + if resolution not in resolutions: + self.add_step(init_step, symlink=symlink) + resolutions.append(resolution) + + timestep, _ = get_timestep_for_task( + config, refinement_factor, refinement=refinement) + timestep = ceil(timestep) + timesteps.append(timestep) + + subdir = f'{basedir}/forward/{mesh_name}_{timestep}s' + symlink = f'forward/{mesh_name}_{timestep}s' + if subdir in component.steps: + forward_step = component.steps[subdir] + else: + forward_step = Forward( + component=component, + refinement=refinement, + refinement_factor=refinement_factor, + name=f'forward_{mesh_name}_{timestep}s', + subdir=subdir, + init=init_step) + forward_step.set_shared_config( + config, link=config_filename) + self.add_step(forward_step, symlink=symlink) + + analysis_dependencies['mesh'][refinement_factor] = init_step + analysis_dependencies['init'][refinement_factor] = init_step + analysis_dependencies['forward'][refinement_factor] = forward_step self.add_step(Analysis(component=component, - resolutions=self.resolutions, - subdir=f'{subdir}/analysis', + subdir=f'{self.subdir}/analysis', + refinement=refinement, dependencies=analysis_dependencies)) - self.add_step(Viz(component=component, resolutions=self.resolutions, + self.add_step(Viz(component=component, resolutions=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 c9383358e..852af6723 100644 --- a/polaris/ocean/tasks/inertial_gravity_wave/analysis.py +++ b/polaris/ocean/tasks/inertial_gravity_wave/analysis.py @@ -1,6 +1,6 @@ import xarray as xr -from polaris.ocean.convergence import ConvergenceAnalysis +from polaris.ocean.convergence.analysis import ConvergenceAnalysis from polaris.ocean.tasks.inertial_gravity_wave.exact_solution import ( ExactSolution, ) @@ -16,7 +16,7 @@ class Analysis(ConvergenceAnalysis): resolutions : list of float The resolutions of the meshes that have been run """ - def __init__(self, component, resolutions, subdir, dependencies): + def __init__(self, component, subdir, dependencies, refinement='both'): """ Create the step @@ -38,11 +38,11 @@ def __init__(self, component, resolutions, subdir, dependencies): 'title': 'SSH', 'zidx': None}] super().__init__(component=component, subdir=subdir, - resolutions=resolutions, dependencies=dependencies, - convergence_vars=convergence_vars) + convergence_vars=convergence_vars, + refinement=refinement) - def exact_solution(self, mesh_name, field_name, time, zidx=None): + def exact_solution(self, refinement_factor, field_name, time, zidx=None): """ Get the exact solution @@ -69,7 +69,7 @@ def exact_solution(self, mesh_name, field_name, time, zidx=None): solution : xarray.DataArray The exact solution as derived from the initial condition """ - init = xr.open_dataset(f'{mesh_name}_init.nc') + init = xr.open_dataset(f'init_r{refinement_factor:02g}.nc') exact = ExactSolution(init, self.config) if field_name != 'ssh': raise ValueError(f'{field_name} is not currently supported') diff --git a/polaris/ocean/tasks/inertial_gravity_wave/forward.py b/polaris/ocean/tasks/inertial_gravity_wave/forward.py index abca61c8f..fecbb37d2 100644 --- a/polaris/ocean/tasks/inertial_gravity_wave/forward.py +++ b/polaris/ocean/tasks/inertial_gravity_wave/forward.py @@ -1,7 +1,8 @@ import numpy as np from polaris.mesh.planar import compute_planar_hex_nx_ny -from polaris.ocean.convergence import ConvergenceForward +from polaris.ocean.convergence import get_resolution_for_task +from polaris.ocean.convergence.forward import ConvergenceForward class Forward(ConvergenceForward): @@ -14,7 +15,8 @@ class Forward(ConvergenceForward): resolution : float The resolution of the test case in km """ - def __init__(self, component, name, resolution, subdir, init): + def __init__(self, component, name, refinement_factor, subdir, + init, refinement='both'): """ Create a new test case @@ -34,7 +36,9 @@ def __init__(self, component, name, resolution, subdir, init): """ super().__init__(component=component, name=name, subdir=subdir, - resolution=resolution, mesh=init, init=init, + refinement=refinement, + refinement_factor=refinement_factor, + mesh=init, init=init, package='polaris.ocean.tasks.inertial_gravity_wave', yaml_filename='forward.yaml', graph_target=f'{init.path}/culled_graph.info', @@ -52,8 +56,10 @@ def compute_cell_count(self): The approximate number of cells in the mesh """ section = self.config['inertial_gravity_wave'] + resolution = get_resolution_for_task( + self.config, self.refinement_factor, refinement=self.refinement) lx = section.getfloat('lx') ly = np.sqrt(3.0) / 2.0 * lx - nx, ny = compute_planar_hex_nx_ny(lx, ly, self.resolution) + nx, ny = compute_planar_hex_nx_ny(lx, ly, resolution) cell_count = nx * ny return cell_count 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 c33e82230..3a319ac7c 100644 --- a/polaris/ocean/tasks/inertial_gravity_wave/inertial_gravity_wave.cfg +++ b/polaris/ocean/tasks/inertial_gravity_wave/inertial_gravity_wave.cfg @@ -57,6 +57,8 @@ convergence_thresh = ${inertial_gravity_wave:conv_thresh} # Type of error to compute error_type = l2 +base_resolution = 100. +refinement_factors = 2., 1., 0.5, 0.25 # config options for spherical convergence tests [convergence_forward] diff --git a/polaris/ocean/tasks/inertial_gravity_wave/init.py b/polaris/ocean/tasks/inertial_gravity_wave/init.py index 4935e912a..f19d89563 100644 --- a/polaris/ocean/tasks/inertial_gravity_wave/init.py +++ b/polaris/ocean/tasks/inertial_gravity_wave/init.py @@ -23,7 +23,7 @@ class Init(Step): resolution : float The resolution of the test case in km """ - def __init__(self, component, resolution, taskdir): + def __init__(self, component, resolution, subdir): """ Create the step @@ -35,13 +35,13 @@ def __init__(self, component, resolution, taskdir): resolution : float The resolution of the test case in km - taskdir : str + subdir : str The subdirectory that the task belongs to """ mesh_name = resolution_to_subdir(resolution) super().__init__(component=component, name=f'init_{mesh_name}', - subdir=f'{taskdir}/init/{mesh_name}') + subdir=subdir) self.resolution = resolution for filename in ['culled_mesh.nc', 'initial_state.nc', 'culled_graph.info']: diff --git a/polaris/ocean/tasks/manufactured_solution/__init__.py b/polaris/ocean/tasks/manufactured_solution/__init__.py index 2c95b634e..7471602c2 100644 --- a/polaris/ocean/tasks/manufactured_solution/__init__.py +++ b/polaris/ocean/tasks/manufactured_solution/__init__.py @@ -1,6 +1,12 @@ +from math import ceil from typing import Dict from polaris import Step, Task +from polaris.config import PolarisConfigParser +from polaris.ocean.convergence import ( + get_resolution_for_task, + get_timestep_for_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 @@ -16,19 +22,22 @@ def add_manufactured_solution_tasks(component): component : polaris.ocean.Ocean the ocean component that the task will be added to """ - component.add_task(ManufacturedSolution(component=component)) + basedir = 'planar/manufactured_solution' + config_filename = 'manufactured_solution.cfg' + filepath = f'{basedir}/{config_filename}' + config = PolarisConfigParser(filepath=filepath) + for refinement in ['space', 'time', 'both']: + component.add_task(ManufacturedSolution(component=component, + config=config, + refinement=refinement)) class ManufacturedSolution(Task): """ The convergence test case using the method of manufactured solutions - - Attributes - ---------- - resolutions : list of floats - The resolutions of the test case in km """ - def __init__(self, component): + + def __init__(self, component, config, refinement='both'): """ Create the test case @@ -36,38 +45,78 @@ def __init__(self, component): ---------- component : polaris.ocean.Ocean The ocean component that this task belongs to + + config : polaris.config.PolarisConfigParser + A shared config parser + + refinement : str, optional + Whether to refine in space, time or both space and time """ - name = 'manufactured_solution' - subdir = f'planar/{name}' + name = f'manufactured_solution_convergence_{refinement}' + basedir = 'planar/manufactured_solution' + subdir = f'{basedir}/convergence_{refinement}' + config_filename = 'manufactured_solution.cfg' super().__init__(component=component, name=name, subdir=subdir) + self.set_shared_config(config, link=config_filename) - self.resolutions = [200., 100., 50., 25.] analysis_dependencies: Dict[str, Dict[float, Step]] = ( dict(mesh=dict(), init=dict(), forward=dict())) - for resolution in self.resolutions: + + refinement_factors = self.config.getlist( + 'convergence', 'refinement_factors', dtype=float) + timesteps = list() + resolutions = list() + for refinement_factor in refinement_factors: + resolution = get_resolution_for_task( + self.config, refinement_factor, refinement=refinement) 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 + + subdir = f'{basedir}/init/{mesh_name}' + symlink = f'init/{mesh_name}' + if subdir in component.steps: + init_step = component.steps[subdir] + else: + init_step = Init(component=component, resolution=resolution, + subdir=subdir) + init_step.set_shared_config(self.config, link=config_filename) + if resolution not in resolutions: + self.add_step(init_step, symlink=symlink) + resolutions.append(resolution) + + timestep, _ = get_timestep_for_task( + config, refinement_factor, refinement=refinement) + timestep = ceil(timestep) + timesteps.append(timestep) + + subdir = f'{basedir}/forward/{mesh_name}_{timestep}s' + symlink = f'forward/{mesh_name}_{timestep}s' + if subdir in component.steps: + forward_step = component.steps[subdir] + else: + forward_step = Forward( + component=component, + refinement=refinement, + refinement_factor=refinement_factor, + name=f'forward_{mesh_name}_{timestep}s', + subdir=subdir, + init=init_step) + forward_step.set_shared_config( + config, link=config_filename) + self.add_step(forward_step, symlink=symlink) + + analysis_dependencies['mesh'][refinement_factor] = init_step + analysis_dependencies['init'][refinement_factor] = init_step + analysis_dependencies['forward'][refinement_factor] = forward_step self.add_step(Analysis(component=component, - resolutions=self.resolutions, subdir=f'{self.subdir}/analysis', + refinement=refinement, dependencies=analysis_dependencies)) - self.add_step(Viz(component=component, resolutions=self.resolutions, + self.add_step(Viz(component=component, resolutions=resolutions, taskdir=self.subdir), run_by_default=False) - - self.config.add_from_package('polaris.ocean.convergence', - 'convergence.cfg') - self.config.add_from_package( + config.add_from_package('polaris.ocean.convergence', + 'convergence.cfg') + config.add_from_package( 'polaris.ocean.tasks.manufactured_solution', - 'manufactured_solution.cfg') + config_filename) diff --git a/polaris/ocean/tasks/manufactured_solution/analysis.py b/polaris/ocean/tasks/manufactured_solution/analysis.py index e7947192e..19fca84aa 100644 --- a/polaris/ocean/tasks/manufactured_solution/analysis.py +++ b/polaris/ocean/tasks/manufactured_solution/analysis.py @@ -1,4 +1,4 @@ -from polaris.ocean.convergence import ConvergenceAnalysis +from polaris.ocean.convergence.analysis import ConvergenceAnalysis from polaris.ocean.tasks.manufactured_solution.exact_solution import ( ExactSolution, ) @@ -14,7 +14,7 @@ class Analysis(ConvergenceAnalysis): resolutions : list of float The resolutions of the meshes that have been run """ - def __init__(self, component, resolutions, subdir, dependencies): + def __init__(self, component, subdir, dependencies, refinement='both'): """ Create the step @@ -36,11 +36,11 @@ def __init__(self, component, resolutions, subdir, dependencies): 'title': 'SSH', 'zidx': None}] super().__init__(component=component, subdir=subdir, - resolutions=resolutions, dependencies=dependencies, - convergence_vars=convergence_vars) + convergence_vars=convergence_vars, + refinement=refinement) - def exact_solution(self, mesh_name, field_name, time, zidx=None): + def exact_solution(self, refinement_factor, field_name, time, zidx=None): """ Get the exact solution @@ -67,7 +67,7 @@ def exact_solution(self, mesh_name, field_name, time, zidx=None): solution : xarray.DataArray The exact solution as derived from the initial condition """ - init = self.open_model_dataset(f'{mesh_name}_init.nc') + init = self.open_model_dataset(f'init_r{refinement_factor:02g}.nc') exact = ExactSolution(self.config, init) if field_name != 'ssh': raise ValueError(f'{field_name} is not currently supported') diff --git a/polaris/ocean/tasks/manufactured_solution/forward.py b/polaris/ocean/tasks/manufactured_solution/forward.py index 62bd76a4d..9643c48f6 100644 --- a/polaris/ocean/tasks/manufactured_solution/forward.py +++ b/polaris/ocean/tasks/manufactured_solution/forward.py @@ -1,7 +1,8 @@ import numpy as np from polaris.mesh.planar import compute_planar_hex_nx_ny -from polaris.ocean.convergence import ConvergenceForward +from polaris.ocean.convergence import get_resolution_for_task +from polaris.ocean.convergence.forward import ConvergenceForward from polaris.ocean.tasks.manufactured_solution.exact_solution import ( ExactSolution, ) @@ -17,7 +18,8 @@ class Forward(ConvergenceForward): resolution : float The resolution of the test case in km """ - def __init__(self, component, name, resolution, subdir, init): + def __init__(self, component, name, refinement_factor, subdir, + init, refinement='both'): """ Create a new test case @@ -37,7 +39,8 @@ def __init__(self, component, name, resolution, subdir, init): """ super().__init__(component=component, name=name, subdir=subdir, - resolution=resolution, mesh=init, init=init, + refinement_factor=refinement_factor, + mesh=init, init=init, refinement=refinement, package='polaris.ocean.tasks.manufactured_solution', yaml_filename='forward.yaml', graph_target=f'{init.path}/culled_graph.info', @@ -66,10 +69,12 @@ def compute_cell_count(self): The approximate number of cells in the mesh """ # no file to read from, so we'll compute it based on config options + resolution = get_resolution_for_task( + self.config, self.refinement_factor, refinement=self.refinement) section = self.config['manufactured_solution'] lx = section.getfloat('lx') ly = np.sqrt(3.0) / 2.0 * lx - nx, ny = compute_planar_hex_nx_ny(lx, ly, self.resolution) + nx, ny = compute_planar_hex_nx_ny(lx, ly, resolution) cell_count = nx * ny return cell_count diff --git a/polaris/ocean/tasks/manufactured_solution/init.py b/polaris/ocean/tasks/manufactured_solution/init.py index 63cf166f2..d08cf15fc 100644 --- a/polaris/ocean/tasks/manufactured_solution/init.py +++ b/polaris/ocean/tasks/manufactured_solution/init.py @@ -22,7 +22,7 @@ class Init(OceanIOStep): resolution : float The resolution of the test case in km """ - def __init__(self, component, resolution, taskdir): + def __init__(self, component, resolution, subdir): """ Create the step @@ -40,7 +40,7 @@ def __init__(self, component, resolution, taskdir): mesh_name = resolution_to_subdir(resolution) super().__init__(component=component, name=f'init_{mesh_name}', - subdir=f'{taskdir}/init/{mesh_name}') + subdir=subdir) self.resolution = resolution def setup(self): diff --git a/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg b/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg index 079efdd92..9771a16e8 100644 --- a/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg +++ b/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg @@ -56,6 +56,8 @@ min_pc_fraction = 0.1 # config options for spherical convergence tests [convergence] +base_resolution = 50. +refinement_factors = 4., 2., 1., 0.5 # Evaluation time for convergence analysis (in hours) convergence_eval_time = ${manufactured_solution:run_duration} diff --git a/polaris/ocean/tasks/sphere_transport/__init__.py b/polaris/ocean/tasks/sphere_transport/__init__.py index 972199fad..1c3154848 100644 --- a/polaris/ocean/tasks/sphere_transport/__init__.py +++ b/polaris/ocean/tasks/sphere_transport/__init__.py @@ -1,7 +1,12 @@ +from math import ceil from typing import Dict from polaris import Step, Task from polaris.config import PolarisConfigParser +from polaris.ocean.convergence import ( + get_resolution_for_task, + get_timestep_for_task, +) from polaris.ocean.mesh.spherical import add_spherical_base_mesh_step from polaris.ocean.tasks.sphere_transport.analysis import Analysis from polaris.ocean.tasks.sphere_transport.filament_analysis import ( @@ -56,7 +61,8 @@ class SphereTransport(Task): include_viz : bool Include VizMap and Viz steps for each resolution """ - def __init__(self, component, config, case_name, icosahedral, include_viz): + def __init__(self, component, config, case_name, icosahedral, include_viz, + refinement='both'): """ Create test case for creating a global MPAS-Ocean mesh @@ -95,6 +101,7 @@ def __init__(self, component, config, case_name, icosahedral, include_viz): link = None super().__init__(component=component, name=name, subdir=subdir) self.resolutions = list() + self.refinement = refinement self.icosahedral = icosahedral self.case_name = case_name self.include_viz = include_viz @@ -102,7 +109,7 @@ def __init__(self, component, config, case_name, icosahedral, include_viz): self.set_shared_config(config, link=link) # add the steps with default resolutions so they can be listed - self._setup_steps() + self._setup_steps(refinement=refinement) def configure(self): """ @@ -111,9 +118,9 @@ def configure(self): super().configure() # set up the steps again in case a user has provided new resolutions - self._setup_steps() + self._setup_steps(self.refinement) - def _setup_steps(self): + def _setup_steps(self, refinement): """ setup steps given resolutions """ case_name = self.case_name icosahedral = self.icosahedral @@ -125,17 +132,26 @@ def _setup_steps(self): else: prefix = 'qu' - resolutions = config.getlist('spherical_convergence', - f'{prefix}_resolutions', dtype=float) - - if self.resolutions == resolutions: - return + refinement_factors = config.getlist('spherical_convergence', + f'{prefix}_refinement_factors', + dtype=str) + refinement_factors = ', '.join(refinement_factors) + config.set('convergence', 'refinement_factors', + value=refinement_factors) + refinement_factors = config.getlist('convergence', + 'refinement_factors', + dtype=float) + + base_resolution = config.getfloat('spherical_convergence', + f'{prefix}_base_resolution') + config.set('convergence', 'base_resolution', + value=f'{base_resolution:03g}') # start fresh with no steps for step in list(self.steps.values()): self.remove_step(step) - self.resolutions = resolutions + resolutions = self.resolutions component = self.component @@ -143,49 +159,62 @@ def _setup_steps(self): analysis_dependencies: Dict[str, Dict[str, Step]] = ( dict(mesh=dict(), init=dict(), forward=dict())) - for resolution in resolutions: - base_mesh_step, mesh_name = add_spherical_base_mesh_step( - component, resolution, icosahedral) - self.add_step(base_mesh_step, symlink=f'base_mesh/{mesh_name}') - analysis_dependencies['mesh'][resolution] = base_mesh_step - - name = f'{prefix}_init_{mesh_name}' - subdir = f'{sph_trans_dir}/init/{mesh_name}' - if self.include_viz: - symlink = f'init/{mesh_name}' - else: - symlink = None - if subdir in component.steps: - init_step = component.steps[subdir] - else: - init_step = Init(component=component, name=name, subdir=subdir, - base_mesh=base_mesh_step, case_name=case_name) - init_step.set_shared_config(config, link=config_filename) - self.add_step(init_step, symlink=symlink) - analysis_dependencies['init'][resolution] = init_step - - name = f'{prefix}_forward_{mesh_name}' - subdir = f'{sph_trans_dir}/forward/{mesh_name}' + timesteps = list() + for idx, refinement_factor in enumerate(refinement_factors): + resolution = get_resolution_for_task( + config, refinement_factor, refinement=refinement) + if resolution not in self.resolutions: + resolutions.append(resolution) + base_mesh_step, mesh_name = add_spherical_base_mesh_step( + component, resolution, icosahedral) + self.add_step(base_mesh_step, symlink=f'base_mesh/{mesh_name}') + analysis_dependencies['mesh'][resolution] = base_mesh_step + + name = f'{prefix}_init_{mesh_name}' + subdir = f'{sph_trans_dir}/init/{mesh_name}' + if self.include_viz: + symlink = f'init/{mesh_name}' + else: + symlink = None + if subdir in component.steps: + init_step = component.steps[subdir] + else: + init_step = Init(component=component, name=name, + subdir=subdir, base_mesh=base_mesh_step, + case_name=case_name) + init_step.set_shared_config(config, link=config_filename) + self.add_step(init_step, symlink=symlink) + analysis_dependencies['init'][resolution] = init_step + + timestep, _ = get_timestep_for_task( + config, refinement_factor, refinement=refinement) + timestep = ceil(timestep) + timesteps.append(timestep) + name = f'{prefix}_forward_{mesh_name}_{timestep}s' + subdir = f'{sph_trans_dir}/forward/{mesh_name}_{timestep}s' if self.include_viz: - symlink = f'forward/{mesh_name}' + symlink = f'forward/{mesh_name}_{timestep}s' else: symlink = None if subdir in component.steps: forward_step = component.steps[subdir] else: - forward_step = Forward(component=component, name=name, - subdir=subdir, resolution=resolution, - base_mesh=base_mesh_step, - init=init_step, - case_name=case_name) + forward_step = Forward( + component=component, name=name, + subdir=subdir, + base_mesh=base_mesh_step, + init=init_step, + case_name=case_name, + refinement_factor=refinement_factor, + refinement=refinement) forward_step.set_shared_config(config, link=config_filename) self.add_step(forward_step, symlink=symlink) analysis_dependencies['forward'][resolution] = forward_step if self.include_viz: with_viz_dir = f'{sph_trans_dir}/with_viz' - name = f'{prefix}_viz_{mesh_name}' - subdir = f'{with_viz_dir}/viz/{mesh_name}' + name = f'{prefix}_viz_{mesh_name}_{timestep}s' + subdir = f'{with_viz_dir}/viz/{mesh_name}_{timestep}s' step = Viz(component=component, name=name, subdir=subdir, base_mesh=base_mesh_step, init=init_step, forward=forward_step, @@ -193,9 +222,9 @@ def _setup_steps(self): step.set_shared_config(config, link=config_filename) self.add_step(step) - subdir = f'{sph_trans_dir}/analysis' + subdir = f'{sph_trans_dir}/analysis/{refinement}' if self.include_viz: - symlink = 'analysis' + symlink = f'analysis_{refinement}' else: symlink = None if subdir in component.steps: @@ -203,9 +232,10 @@ def _setup_steps(self): step.resolutions = resolutions step.dependencies_dict = analysis_dependencies else: - step = Analysis(component=component, resolutions=resolutions, + step = Analysis(component=component, subdir=subdir, case_name=case_name, - dependencies=analysis_dependencies) + dependencies=analysis_dependencies, + refinement=refinement) step.set_shared_config(config, link=config_filename) self.add_step(step, symlink=symlink) diff --git a/polaris/ocean/tasks/sphere_transport/analysis.py b/polaris/ocean/tasks/sphere_transport/analysis.py index 620a80651..00659bef2 100644 --- a/polaris/ocean/tasks/sphere_transport/analysis.py +++ b/polaris/ocean/tasks/sphere_transport/analysis.py @@ -1,4 +1,4 @@ -from polaris.ocean.convergence import ConvergenceAnalysis +from polaris.ocean.convergence.analysis import ConvergenceAnalysis class Analysis(ConvergenceAnalysis): @@ -13,8 +13,8 @@ class Analysis(ConvergenceAnalysis): case_name : str The name of the test case """ - def __init__(self, component, resolutions, subdir, case_name, - dependencies): + def __init__(self, component, subdir, case_name, + dependencies, refinement='both'): """ Create the step @@ -23,9 +23,6 @@ def __init__(self, component, resolutions, subdir, case_name, component : polaris.Component The component the step belongs to - resolutions : list of float - The resolutions of the meshes that have been run - subdir : str The subdirectory that the step resides in @@ -46,9 +43,9 @@ def __init__(self, component, resolutions, subdir, case_name, 'title': 'tracer3', 'zidx': 1}] super().__init__(component=component, subdir=subdir, - resolutions=resolutions, dependencies=dependencies, - convergence_vars=convergence_vars) + convergence_vars=convergence_vars, + refinement=refinement) # Note: there is no need to overwrite the default method exact_solution # which uses the initial condition diff --git a/polaris/ocean/tasks/sphere_transport/forward.py b/polaris/ocean/tasks/sphere_transport/forward.py index c795993af..0b2250881 100644 --- a/polaris/ocean/tasks/sphere_transport/forward.py +++ b/polaris/ocean/tasks/sphere_transport/forward.py @@ -7,8 +7,8 @@ class Forward(SphericalConvergenceForward): transport test case """ - def __init__(self, component, name, subdir, resolution, base_mesh, init, - case_name): + def __init__(self, component, name, subdir, base_mesh, init, + case_name, refinement_factor, refinement): """ Create a new step @@ -23,9 +23,6 @@ def __init__(self, component, name, subdir, resolution, base_mesh, init, subdir : str The subdirectory for the step - resolution : float - The resolution of the (uniform) mesh in km - base_mesh : polaris.Step The base mesh step @@ -47,10 +44,12 @@ def __init__(self, component, name, subdir, resolution, base_mesh, init, } validate_vars = ['normalVelocity', 'tracer1', 'tracer2', 'tracer3'] super().__init__(component=component, name=name, subdir=subdir, - resolution=resolution, mesh=base_mesh, + mesh=base_mesh, init=init, package=package, yaml_filename='forward.yaml', output_filename='output.nc', validate_vars=validate_vars, options=namelist_options, - graph_target=f'{base_mesh.path}/graph.info') + graph_target=f'{base_mesh.path}/graph.info', + refinement_factor=refinement_factor, + refinement=refinement) From ad8e4d681f67b154dcfaec8628c15b9f12e897b7 Mon Sep 17 00:00:00 2001 From: Carolyn Begeman Date: Wed, 9 Oct 2024 15:32:45 -0500 Subject: [PATCH 2/9] Update docs for new convergence tests --- docs/developers_guide/ocean/framework.md | 51 ++++++++++++++++--- docs/users_guide/ocean/tasks/cosine_bell.md | 42 +++++++++++---- .../ocean/tasks/inertial_gravity_wave.md | 8 +++ .../ocean/tasks/manufactured_solution.md | 10 +++- 4 files changed, 93 insertions(+), 18 deletions(-) diff --git a/docs/developers_guide/ocean/framework.md b/docs/developers_guide/ocean/framework.md index cca5f7986..d1c324f0a 100644 --- a/docs/developers_guide/ocean/framework.md +++ b/docs/developers_guide/ocean/framework.md @@ -219,16 +219,35 @@ 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 key config options that control the convergence test are `base_resolution` +and `refinement_factors`. The `base_resolution` is multipled by the +`refinement_factors` to determine which resolutions to test when the +convergence is being tested in space (or space and time together). The +`base_resolution` is applied to all steps when convergence in time is tested. +`base_resolution` times `dt_per_km` determines the base timestep in that case +and is then multiplied by the `refinement_factors` to determine which time steps +to test. When spherical meshes are being tested, the values in the +`convergence` section are overridden by their values in the +`spherical_convergence` section with a prefix indicating the mesh type. + The shared config options are: ```cfg # config options for spherical convergence tests [spherical_convergence] +# The base resolution for the icosahedral mesh to which the refinement +# factors are applied +icos_base_resolution = 60. + # a list of icosahedral mesh resolutions (km) to test -icos_resolutions = 60, 120, 240, 480 +icos_refinement_factors = 8., 4., 2., 1. + +# The base resolution for the quasi-uniform mesh to which the refinement +# factors are applied +qu_base_resolution = 120. # a list of quasi-uniform mesh resolutions (km) to test -qu_resolutions = 60, 90, 120, 150, 180, 210, 240 +qu_refinement_factors = 0.5, 0.75, 1., 1.25, 1.5, 1.75, 2. [convergence] @@ -241,6 +260,15 @@ convergence_thresh = 1.0 # Type of error to compute error_type = l2 +# the base mesh resolution (km) to which refinement_factors +# are applied if refinement is 'space' or 'both' on a planar mesh +# base resolutions for spherical meshes are given in section spherical_convergence +base_resolution = 120 + +# refinement factors for a planar mesh applied to either space or time +# refinement factors for a spherical mesh given in section spherical_convergence +refinement_factors = 4., 2., 1., 0.5 + # config options for convergence forward steps [convergence_forward] @@ -325,7 +353,8 @@ class Forward(SphericalConvergenceForward): bell test case """ - def __init__(self, component, name, subdir, resolution, mesh, init): + def __init__(self, component, name, subdir, mesh, init, + refinement_factor, refinement='both'): """ Create a new step @@ -348,6 +377,9 @@ class Forward(SphericalConvergenceForward): init : polaris.Step The init step + + refinement_factor : float + The factor by which to scale space, time or both """ package = 'polaris.ocean.tasks.cosine_bell' validate_vars = ['normalVelocity', 'tracer1'] @@ -356,7 +388,10 @@ class Forward(SphericalConvergenceForward): init=init, package=package, yaml_filename='forward.yaml', output_filename='output.nc', - validate_vars=validate_vars) + validate_vars=validate_vars, + graph_target=f'{init.path}/graph.info', + refinement_factor=refinement_factor, + refinement=refinement) ``` Each convergence test must define a YAML file with model config options, called `forward.yaml` by default. The `package` parameter is the location of this @@ -419,7 +454,7 @@ class Analysis(ConvergenceAnalysis): """ A step for analyzing the output from the cosine bell test case """ - def __init__(self, component, resolutions, subdir, dependencies): + def __init__(self, component, subdir, dependencies, refinement='both'): """ Create the step @@ -436,6 +471,9 @@ class Analysis(ConvergenceAnalysis): dependencies : dict of dict of polaris.Steps The dependencies of this step + + refinement : str, optional + Whether to refine in space, time or both space and time """ convergence_vars = [{'name': 'tracer1', 'title': 'tracer1', @@ -443,7 +481,8 @@ class Analysis(ConvergenceAnalysis): super().__init__(component=component, subdir=subdir, resolutions=resolutions, dependencies=dependencies, - convergence_vars=convergence_vars) + convergence_vars=convergence_vars, + refinement=refinement) ``` Many tasks will also need to override the diff --git a/docs/users_guide/ocean/tasks/cosine_bell.md b/docs/users_guide/ocean/tasks/cosine_bell.md index f49afa97a..5230fa36b 100644 --- a/docs/users_guide/ocean/tasks/cosine_bell.md +++ b/docs/users_guide/ocean/tasks/cosine_bell.md @@ -33,34 +33,56 @@ These tasks support only MPAS-Ocean. (ocean-cosine-bell-mesh)= ## mesh -Two global mesh variants are tested, quasi-uniform (QU) and icosohydral. Thus, -there are 4 variants of the task: +Two global mesh variants are tested, quasi-uniform (QU) and icosohydral. There +are also variants to test convergence in space, time, or both space and time. +In addition, the tests can be set up with or without the viz step. Thus, there +are 12 variants of the task: ``` -ocean/spherical/icos/cosine_bell -ocean/spherical/icos/cosine_bell/with_viz -ocean/spherical/qu/cosine_bell -ocean/spherical/qu/cosine_bell/with_viz +ocean/spherical/icos/cosine_bell/convergence_both/ +ocean/spherical/icos/cosine_bell/convergence_both/with_viz +ocean/spherical/qu/cosine_bell/convergence_both/ +ocean/spherical/qu/cosine_bell/convergence_both/with_viz +ocean/spherical/icos/cosine_bell/convergence_space/ +ocean/spherical/icos/cosine_bell/convergence_space/with_viz +ocean/spherical/qu/cosine_bell/convergence_space/ +ocean/spherical/qu/cosine_bell/convergence_space/with_viz +ocean/spherical/icos/cosine_bell/convergence_time/ +ocean/spherical/icos/cosine_bell/convergence_time/with_viz +ocean/spherical/qu/cosine_bell/convergence_time/ +ocean/spherical/qu/cosine_bell/convergence_time/with_viz ``` The default resolutions used in the task depends on the mesh type. -For the `icos` mesh type, the defaults are: +For the `icos` mesh type, the defaults are 60, 120, 240, 480 km, as determined +by the following config options. See {ref}`dev-ocean-convergence` for more +details. ```cfg # config options for spherical convergence tests [spherical_convergence] +# The base resolution for the icosahedral mesh to which the refinement +# factors are applied +icos_base_resolution = 60. + +# a list of icosahedral mesh resolutions (km) to test +icos_refinement_factors = 8., 4., 2., 1. # a list of icosahedral mesh resolutions (km) to test -icos_resolutions = 60, 120, 240, 480 ``` -for the `qu` mesh type, they are: +For the `qu` mesh type, they are 60, 90, 120, 150, 180, 210, 240 km as +determined by the following config options: ```cfg # config options for spherical convergence tests [spherical_convergence] +# The base resolution for the quasi-uniform mesh to which the refinement +# factors are applied +qu_base_resolution = 120. + # a list of quasi-uniform mesh resolutions (km) to test -qu_resolutions = 60, 90, 120, 150, 180, 210, 240 +qu_refinement_factors = 0.5, 0.75, 1., 1.25, 1.5, 1.75, 2. ``` To alter the resolutions used in this task, you will need to create your own diff --git a/docs/users_guide/ocean/tasks/inertial_gravity_wave.md b/docs/users_guide/ocean/tasks/inertial_gravity_wave.md index f858f232e..709166e4c 100644 --- a/docs/users_guide/ocean/tasks/inertial_gravity_wave.md +++ b/docs/users_guide/ocean/tasks/inertial_gravity_wave.md @@ -11,6 +11,14 @@ The implementation is from ## description +There are 3 versions of the convergence test case, `convergence_space`, +`convergence_time`, and `convergence_both` corresponding to space, time, and +space and time convergence tests. Tests involving spatial convergence run the +inertial gravity wave simulation for 4 different resolutions: 200, 100, 50, +and 25 km. Tests involving temporal convergence use the parameter `dt_per_km` +at the `base_resolution` multiplied by `refinement_factors` (see +{ref}`dev-ocean-convergence` for more details on how to change resolutions or +time steps tested). The `inertial_gravity_wave` task runs the inertial gravity wave simulation for 4 different resolutions: 200, 100, 50, and 25 km. diff --git a/docs/users_guide/ocean/tasks/manufactured_solution.md b/docs/users_guide/ocean/tasks/manufactured_solution.md index f9e840e0f..866f2d85d 100644 --- a/docs/users_guide/ocean/tasks/manufactured_solution.md +++ b/docs/users_guide/ocean/tasks/manufactured_solution.md @@ -22,8 +22,14 @@ These tasks support both MPAS-Ocean and Omega. ### description -The `convergence` test case runs the manufactured solution simulation for 4 -different resolutions: 200, 100, 50, and 25 km. +There are 3 versions of the convergence test case, `convergence_space`, +`convergence_time`, and `convergence_both` corresponding to space, time, and +space and time convergence tests. Tests involving spatial convergence run the +manufactured solution simulation for 4 different resolutions: 200, 100, 50, +and 25 km. Tests involving temporal convergence use the parameter `dt_per_km` +at the `base_resolution` multiplied by `refinement_factors` (see +{ref}`dev-ocean-convergence` for more details on how to change resolutions or +time steps tested). The forward step for each resolution runs the simulation for 10 hours. The model is configured without vertical advection and mixing. No tracers are enabled From c661b9397c00eb4faf931877fd7323b3bbb4ed7d Mon Sep 17 00:00:00 2001 From: Carolyn Begeman Date: Wed, 9 Oct 2024 15:34:38 -0500 Subject: [PATCH 3/9] Update doc strings --- polaris/ocean/convergence/__init__.py | 39 +++++++++++++++++++ polaris/ocean/convergence/analysis.py | 13 ++++--- polaris/ocean/convergence/forward.py | 15 +++++++ .../ocean/convergence/spherical/spherical.cfg | 4 +- polaris/ocean/tasks/cosine_bell/__init__.py | 8 ++-- polaris/ocean/tasks/cosine_bell/analysis.py | 9 +++-- polaris/ocean/tasks/cosine_bell/forward.py | 11 +++--- polaris/ocean/tasks/geostrophic/__init__.py | 10 +++-- polaris/ocean/tasks/geostrophic/analysis.py | 11 ++++-- polaris/ocean/tasks/geostrophic/forward.py | 7 ++++ .../tasks/inertial_gravity_wave/analysis.py | 7 +++- .../tasks/inertial_gravity_wave/forward.py | 20 +++++++--- .../inertial_gravity_wave.cfg | 4 ++ .../tasks/manufactured_solution/analysis.py | 7 +++- .../tasks/manufactured_solution/forward.py | 19 +++++++-- .../manufactured_solution.cfg | 8 +++- .../ocean/tasks/sphere_transport/__init__.py | 8 ++++ .../ocean/tasks/sphere_transport/analysis.py | 4 ++ .../ocean/tasks/sphere_transport/forward.py | 7 ++++ .../sphere_transport/resources/flow_types.py | 4 +- 20 files changed, 172 insertions(+), 43 deletions(-) diff --git a/polaris/ocean/convergence/__init__.py b/polaris/ocean/convergence/__init__.py index 77971171d..a66abb2e4 100644 --- a/polaris/ocean/convergence/__init__.py +++ b/polaris/ocean/convergence/__init__.py @@ -3,7 +3,26 @@ def get_resolution_for_task(config, refinement_factor, refinement='both'): + """ + Get the resolution for a step in a convergence task + Parameters + ---------- + config : polaris.Config + The config options for this task + + refinement_factor : float + The factor by which either resolution or time is refined for this step + + refinement : str, optional + Whether to refine in space, time or both + + Returns + ------- + resolution : float + The resolution corresponding to the refinement_factor and convergence + test type + """ base_resolution = config.getfloat('convergence', 'base_resolution') refinement_factors = config.getlist('convergence', 'refinement_factors', @@ -23,6 +42,26 @@ def get_resolution_for_task(config, refinement_factor, def get_timestep_for_task(config, refinement_factor, refinement='both'): + """ + Get the time step for a forward step in a convergence task + + Parameters + ---------- + config : polaris.Config + The config options for this task + + refinement_factor : float + The factor by which either resolution or time is refined for this step + + refinement : str, optional + Whether to refine in space, time or both + + Returns + ------- + resolution : float + The resolution corresponding to the refinement_factor and convergence + test type + """ base_resolution = config.getfloat('convergence', 'base_resolution') refinement_factors = config.getlist('convergence', diff --git a/polaris/ocean/convergence/analysis.py b/polaris/ocean/convergence/analysis.py index 0586fceec..f44f719e1 100644 --- a/polaris/ocean/convergence/analysis.py +++ b/polaris/ocean/convergence/analysis.py @@ -19,9 +19,6 @@ class ConvergenceAnalysis(OceanIOStep): Attributes ---------- - resolutions : list of float - The resolutions of the meshes that have been run - dependencies_dict : dict of dict of polaris.Steps The dependencies of this step must be given as separate keys in the dict: @@ -55,6 +52,10 @@ class ConvergenceAnalysis(OceanIOStep): zidx : int The z-index to use for variables that have an nVertLevels dimension, which should be None for variables that don't + + refinement : str + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ def __init__(self, component, subdir, dependencies, convergence_vars, refinement='both'): @@ -74,17 +75,17 @@ def __init__(self, component, subdir, dependencies, convergence_vars, dict: mesh : dict of polaris.Steps - Keys of the dict correspond to `resolutions` + Keys of the dict correspond to `refinement_factors` Values of the dict are polaris.Steps, which must have the attribute `path`, the path to `base_mesh.nc` of that resolution init : dict of polaris.Steps - Keys of the dict correspond to `resolutions` + Keys of the dict correspond to `refinement_factors` Values of the dict are polaris.Steps, which must have the attribute `path`, the path to `initial_state.nc` of that resolution forward : dict of polaris.Steps - Keys of the dict correspond to `resolutions` + Keys of the dict correspond to `refinement_factors` Values of the dict are polaris.Steps, which must have the attribute `path`, the path to `forward.nc` of that resolution diff --git a/polaris/ocean/convergence/forward.py b/polaris/ocean/convergence/forward.py index 710d1f3d8..1ea1c47e2 100644 --- a/polaris/ocean/convergence/forward.py +++ b/polaris/ocean/convergence/forward.py @@ -15,6 +15,13 @@ class ConvergenceForward(OceanModelStep): yaml_filename : str A YAML file that is a Jinja2 template with model config options + refinement : str + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time + + refinement_factor : float + Refinement factor use to scale space, time or both depending on + refinement option """ def __init__(self, component, name, subdir, refinement_factor, mesh, init, @@ -35,6 +42,10 @@ def __init__(self, component, name, subdir, refinement_factor, mesh, init, subdir : str The subdirectory for the step + refinement_factor : float + Refinement factor use to scale space, time or both depending on + refinement option + package : Package The package name or module object that contains the YAML file @@ -55,6 +66,10 @@ def __init__(self, component, name, subdir, refinement_factor, mesh, init, validate_vars : list of str, optional A list of variables to validate against a baseline if requested + + refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ super().__init__(component=component, name=name, subdir=subdir, openmp_threads=1, graph_target=graph_target) diff --git a/polaris/ocean/convergence/spherical/spherical.cfg b/polaris/ocean/convergence/spherical/spherical.cfg index 0dc5fbbb3..fd769f8ba 100644 --- a/polaris/ocean/convergence/spherical/spherical.cfg +++ b/polaris/ocean/convergence/spherical/spherical.cfg @@ -5,12 +5,12 @@ # factors are applied icos_base_resolution = 60. -# a list of icosahedral mesh resolutions (km) to test +# The factors by which to scale space or time based on icos_base_resolution icos_refinement_factors = 8., 4., 2., 1. # The base resolution for the quasi-uniform mesh to which the refinement # factors are applied qu_base_resolution = 120. -# a list of quasi-uniform mesh resolutions (km) to test +# The factors by which to scale space or time based on qu_base_resolution qu_refinement_factors = 0.5, 0.75, 1., 1.25, 1.5, 1.75, 2. diff --git a/polaris/ocean/tasks/cosine_bell/__init__.py b/polaris/ocean/tasks/cosine_bell/__init__.py index 55ca91d82..34e527211 100644 --- a/polaris/ocean/tasks/cosine_bell/__init__.py +++ b/polaris/ocean/tasks/cosine_bell/__init__.py @@ -49,8 +49,9 @@ class CosineBell(Task): Attributes ---------- - resolutions : list of float - A list of mesh resolutions + refinement : str + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time icosahedral : bool Whether to use icosahedral, as opposed to less regular, JIGSAW meshes @@ -79,6 +80,8 @@ def __init__(self, component, config, icosahedral, include_viz, Include VizMap and Viz steps for each resolution refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ if icosahedral: prefix = 'icos' @@ -92,7 +95,6 @@ def __init__(self, component, config, icosahedral, include_viz, name = f'{name}_with_viz' link = 'cosine_bell.cfg' super().__init__(component=component, name=name, subdir=subdir) - self.resolutions = list() self.refinement = refinement self.icosahedral = icosahedral self.include_viz = include_viz diff --git a/polaris/ocean/tasks/cosine_bell/analysis.py b/polaris/ocean/tasks/cosine_bell/analysis.py index 79106a655..28fd7576e 100644 --- a/polaris/ocean/tasks/cosine_bell/analysis.py +++ b/polaris/ocean/tasks/cosine_bell/analysis.py @@ -25,6 +25,9 @@ def __init__(self, component, subdir, dependencies, refinement='both'): dependencies : dict of dict of polaris.Steps The dependencies of this step + + refinement : str, optional + Whether to refine in space, time or both space and time """ convergence_vars = [{'name': 'tracer1', 'title': 'tracer1', @@ -40,8 +43,8 @@ def exact_solution(self, refinement_factor, field_name, time, zidx=None): Parameters ---------- - mesh_name : str - The mesh name which is the prefix for the initial condition file + refinement_factor : float + The factor by which to scale space, time or both field_name : str The name of the variable of which we evaluate convergence @@ -58,7 +61,7 @@ def exact_solution(self, refinement_factor, field_name, time, zidx=None): Returns ------- - solution: xarray.DataArray + solution : xarray.DataArray The exact solution with dimension nCells """ diff --git a/polaris/ocean/tasks/cosine_bell/forward.py b/polaris/ocean/tasks/cosine_bell/forward.py index 275b0b8a8..8c477d9e7 100644 --- a/polaris/ocean/tasks/cosine_bell/forward.py +++ b/polaris/ocean/tasks/cosine_bell/forward.py @@ -23,17 +23,18 @@ def __init__(self, component, name, subdir, mesh, init, subdir : str The subdirectory for the step - resolution : float - The resolution of the (uniform) mesh in km - mesh : polaris.Step The base mesh step init : polaris.Step The init step - time_refinement_factor : float - The amount by which to scale dt_per_km and btr_dt_per_km + refinement_factor : float + The factor by which to scale space, time or both + + refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ package = 'polaris.ocean.tasks.cosine_bell' validate_vars = ['normalVelocity', 'tracer1'] diff --git a/polaris/ocean/tasks/geostrophic/__init__.py b/polaris/ocean/tasks/geostrophic/__init__.py index 2ee90f021..a3b2e941c 100644 --- a/polaris/ocean/tasks/geostrophic/__init__.py +++ b/polaris/ocean/tasks/geostrophic/__init__.py @@ -48,8 +48,9 @@ class Geostrophic(Task): Attributes ---------- - resolutions : list of float - A list of mesh resolutions + refinement : str + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time icosahedral : bool Whether to use icosahedral, as opposed to less regular, JIGSAW meshes @@ -76,6 +77,10 @@ def __init__(self, component, config, icosahedral, include_viz, include_viz : bool Include VizMap and Viz steps for each resolution + + refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ if icosahedral: prefix = 'icos' @@ -90,7 +95,6 @@ def __init__(self, component, config, icosahedral, include_viz, link = 'geostrophic.cfg' super().__init__(component=component, name=name, subdir=subdir) - self.resolutions = list() self.refinement = refinement self.icosahedral = icosahedral self.include_viz = include_viz diff --git a/polaris/ocean/tasks/geostrophic/analysis.py b/polaris/ocean/tasks/geostrophic/analysis.py index e4fab12d4..a26f01fdb 100644 --- a/polaris/ocean/tasks/geostrophic/analysis.py +++ b/polaris/ocean/tasks/geostrophic/analysis.py @@ -24,6 +24,9 @@ def __init__(self, component, subdir, dependencies, refinement='both'): dependencies : dict of dict of polaris.Steps The dependencies of this step + + refinement : str, optional + Whether to refine in space, time or both space and time """ convergence_vars = [{'name': 'h', 'title': 'water-column thickness', @@ -36,14 +39,14 @@ def __init__(self, component, subdir, dependencies, refinement='both'): convergence_vars=convergence_vars, refinement=refinement) - def exact_solution(self, mesh_name, field_name, time, zidx=None): + def exact_solution(self, refinement_factor, field_name, time, zidx=None): """ Get the exact solution Parameters ---------- - mesh_name : str - The mesh name which is the prefix for the initial condition file + refinement_factor : float + The factor by which to scale space, time or both field_name : str The name of the variable of which we evaluate convergence @@ -75,7 +78,7 @@ def exact_solution(self, mesh_name, field_name, time, zidx=None): vel_period = section.getfloat('vel_period') gh_0 = section.getfloat('gh_0') - mesh_filename = f'{mesh_name}_mesh.nc' + mesh_filename = f'mesh_r{refinement_factor:02g}.nc' h, _, _, normalVelocity = compute_exact_solution( alpha, vel_period, gh_0, mesh_filename) diff --git a/polaris/ocean/tasks/geostrophic/forward.py b/polaris/ocean/tasks/geostrophic/forward.py index 13249a1c7..8e8188ea8 100644 --- a/polaris/ocean/tasks/geostrophic/forward.py +++ b/polaris/ocean/tasks/geostrophic/forward.py @@ -28,6 +28,13 @@ def __init__(self, component, name, subdir, mesh, init, init : polaris.Step The init step + + refinement_factor : float + The factor by which to scale space, time or both + + refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ package = 'polaris.ocean.tasks.geostrophic' validate_vars = ['temperature', diff --git a/polaris/ocean/tasks/inertial_gravity_wave/analysis.py b/polaris/ocean/tasks/inertial_gravity_wave/analysis.py index 852af6723..4cdb4e59a 100644 --- a/polaris/ocean/tasks/inertial_gravity_wave/analysis.py +++ b/polaris/ocean/tasks/inertial_gravity_wave/analysis.py @@ -33,6 +33,9 @@ def __init__(self, component, subdir, dependencies, refinement='both'): dependencies : dict of dict of polaris.Steps The dependencies of this step + + refinement : str, optional + Whether to refine in space, time or both space and time """ convergence_vars = [{'name': 'ssh', 'title': 'SSH', @@ -48,8 +51,8 @@ def exact_solution(self, refinement_factor, field_name, time, zidx=None): Parameters ---------- - mesh_name : str - The mesh name which is the prefix for the initial condition file + refinement_factor : float + The factor by which to scale space, time or both field_name : str The name of the variable of which we evaluate convergence diff --git a/polaris/ocean/tasks/inertial_gravity_wave/forward.py b/polaris/ocean/tasks/inertial_gravity_wave/forward.py index fecbb37d2..68930dde6 100644 --- a/polaris/ocean/tasks/inertial_gravity_wave/forward.py +++ b/polaris/ocean/tasks/inertial_gravity_wave/forward.py @@ -12,8 +12,12 @@ class Forward(ConvergenceForward): Attributes ---------- - resolution : float - The resolution of the test case in km + refinement_factor : float + The factor by which to scale space, time or both + + refinement : str + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ def __init__(self, component, name, refinement_factor, subdir, init, refinement='both'): @@ -25,14 +29,20 @@ def __init__(self, component, name, refinement_factor, subdir, component : polaris.Component The component the step belongs to - resolution : km - The resolution of the test case in km + name : str + The name of the step + refinement_factor : float + The factor by which to scale space, time or both subdir : str - The subdirectory that the step belongs to + The subdirectory that the task belongs to init : polaris.Step The step which generates the mesh and initial condition + + refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ super().__init__(component=component, name=name, subdir=subdir, 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 3a319ac7c..b690a3de1 100644 --- a/polaris/ocean/tasks/inertial_gravity_wave/inertial_gravity_wave.cfg +++ b/polaris/ocean/tasks/inertial_gravity_wave/inertial_gravity_wave.cfg @@ -57,7 +57,11 @@ convergence_thresh = ${inertial_gravity_wave:conv_thresh} # Type of error to compute error_type = l2 +# the base mesh resolution (km) to which refinement_factors +# are applied if refinement is 'space' or 'both' on a planar mesh base_resolution = 100. + +# refinement factors for a planar mesh applied to either space or time refinement_factors = 2., 1., 0.5, 0.25 # config options for spherical convergence tests diff --git a/polaris/ocean/tasks/manufactured_solution/analysis.py b/polaris/ocean/tasks/manufactured_solution/analysis.py index 19fca84aa..06ba74046 100644 --- a/polaris/ocean/tasks/manufactured_solution/analysis.py +++ b/polaris/ocean/tasks/manufactured_solution/analysis.py @@ -31,6 +31,9 @@ def __init__(self, component, subdir, dependencies, refinement='both'): dependencies : dict of dict of polaris.Steps The dependencies of this step + + refinement : str, optional + Whether to refine in space, time or both space and time """ convergence_vars = [{'name': 'ssh', 'title': 'SSH', @@ -46,8 +49,8 @@ def exact_solution(self, refinement_factor, field_name, time, zidx=None): Parameters ---------- - mesh_name : str - The mesh name which is the prefix for the initial condition file + refinement_factor : float + The factor by which to scale space, time or both field_name : str The name of the variable of which we evaluate convergence diff --git a/polaris/ocean/tasks/manufactured_solution/forward.py b/polaris/ocean/tasks/manufactured_solution/forward.py index 9643c48f6..012bea98f 100644 --- a/polaris/ocean/tasks/manufactured_solution/forward.py +++ b/polaris/ocean/tasks/manufactured_solution/forward.py @@ -15,8 +15,12 @@ class Forward(ConvergenceForward): Attributes ---------- - resolution : float - The resolution of the test case in km + refinement_factor : float + The factor by which to scale space, time or both + + refinement : str + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ def __init__(self, component, name, refinement_factor, subdir, init, refinement='both'): @@ -28,14 +32,21 @@ def __init__(self, component, name, refinement_factor, subdir, component : polaris.Component The component the step belongs to - resolution : km - The resolution of the test case in km + name : str + The name of the step + + refinement_factor : float + The factor by which to scale space, time or both subdir : str The subdirectory that the task belongs to init : polaris.Step The step which generates the mesh and initial condition + + refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ super().__init__(component=component, name=name, subdir=subdir, diff --git a/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg b/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg index 9771a16e8..be89113c0 100644 --- a/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg +++ b/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg @@ -56,8 +56,6 @@ min_pc_fraction = 0.1 # config options for spherical convergence tests [convergence] -base_resolution = 50. -refinement_factors = 4., 2., 1., 0.5 # Evaluation time for convergence analysis (in hours) convergence_eval_time = ${manufactured_solution:run_duration} @@ -67,6 +65,12 @@ convergence_thresh = ${manufactured_solution:conv_thresh} # Type of error to compute error_type = l2 +# the base mesh resolution (km) to which refinement_factors +# are applied if refinement is 'space' or 'both' on a planar mesh +base_resolution = 50. + +# refinement factors for a planar mesh applied to either space or time +refinement_factors = 4., 2., 1., 0.5 # config options for spherical convergence tests [convergence_forward] diff --git a/polaris/ocean/tasks/sphere_transport/__init__.py b/polaris/ocean/tasks/sphere_transport/__init__.py index 1c3154848..3596c5e95 100644 --- a/polaris/ocean/tasks/sphere_transport/__init__.py +++ b/polaris/ocean/tasks/sphere_transport/__init__.py @@ -55,6 +55,10 @@ class SphereTransport(Task): resolutions : list of float A list of mesh resolutions + refinement : str + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time + icosahedral : bool Whether to use icosahedral, as opposed to less regular, JIGSAW meshes @@ -84,6 +88,10 @@ def __init__(self, component, config, case_name, icosahedral, include_viz, include_viz : bool Include VizMap and Viz steps for each resolution + + refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ if icosahedral: prefix = 'icos' diff --git a/polaris/ocean/tasks/sphere_transport/analysis.py b/polaris/ocean/tasks/sphere_transport/analysis.py index 00659bef2..17983cce3 100644 --- a/polaris/ocean/tasks/sphere_transport/analysis.py +++ b/polaris/ocean/tasks/sphere_transport/analysis.py @@ -31,6 +31,10 @@ def __init__(self, component, subdir, case_name, dependencies : dict of dict of polaris.Steps The dependencies of this step + + refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ self.case_name = case_name convergence_vars = [{'name': 'tracer1', diff --git a/polaris/ocean/tasks/sphere_transport/forward.py b/polaris/ocean/tasks/sphere_transport/forward.py index 0b2250881..6ce2195b8 100644 --- a/polaris/ocean/tasks/sphere_transport/forward.py +++ b/polaris/ocean/tasks/sphere_transport/forward.py @@ -31,6 +31,13 @@ def __init__(self, component, name, subdir, base_mesh, init, case_name: str The name of the test case + + refinement_factor : float + The factor by which to scale space, time or both + + refinement : str, optional + Refinement type. One of 'space', 'time' or 'both' indicating both + space and time """ package = 'polaris.ocean.tasks.sphere_transport' flow_id = {'rotation_2d': 1, diff --git a/polaris/ocean/tasks/sphere_transport/resources/flow_types.py b/polaris/ocean/tasks/sphere_transport/resources/flow_types.py index 78a9adf95..3645d2d27 100644 --- a/polaris/ocean/tasks/sphere_transport/resources/flow_types.py +++ b/polaris/ocean/tasks/sphere_transport/resources/flow_types.py @@ -43,7 +43,7 @@ def flow_nondivergent(t, lon, lat, u_0, tau): def flow_divergent(t, lon, lat, u_0, tau): """ - Compute a nondivergent velocity field + Compute a divergent velocity field Parameters ---------- @@ -81,7 +81,7 @@ def flow_divergent(t, lon, lat, u_0, tau): def flow_rotation(lon, lat, omega, tau, sphere_radius): """ - Compute a nondivergent velocity field + Compute a rotational velocity field Parameters ---------- From d29a1990746ac995b09321161a323d9868b83cd2 Mon Sep 17 00:00:00 2001 From: Carolyn Begeman Date: Thu, 10 Oct 2024 10:34:01 -0500 Subject: [PATCH 4/9] Update suites --- polaris/ocean/suites/convergence.txt | 12 ++++++------ polaris/ocean/suites/geostrophic.txt | 4 ++-- polaris/ocean/suites/nightly.txt | 2 +- polaris/ocean/suites/pr.txt | 2 +- 4 files changed, 10 insertions(+), 10 deletions(-) diff --git a/polaris/ocean/suites/convergence.txt b/polaris/ocean/suites/convergence.txt index 56a8e63c7..9a92ef47e 100644 --- a/polaris/ocean/suites/convergence.txt +++ b/polaris/ocean/suites/convergence.txt @@ -1,22 +1,22 @@ -ocean/planar/inertial_gravity_wave -ocean/planar/manufactured_solution +ocean/planar/inertial_gravity_wave/convergence_both +ocean/planar/manufactured_solution/convergence_both ocean/spherical/icos/correlated_tracers_2d cached: icos_base_mesh_60km icos_base_mesh_120km icos_base_mesh_240km icos_base_mesh_480km ocean/spherical/qu/correlated_tracers_2d cached: qu_base_mesh_60km qu_base_mesh_90km qu_base_mesh_120km qu_base_mesh_150km cached: qu_base_mesh_180km qu_base_mesh_210km qu_base_mesh_240km -ocean/spherical/icos/cosine_bell +ocean/spherical/icos/cosine_bell/convergence_both 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 +ocean/spherical/qu/cosine_bell/convergence_both 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 -ocean/spherical/icos/geostrophic +ocean/spherical/icos/geostrophic/convergence_both 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/geostrophic +ocean/spherical/qu/geostrophic/convergence_both 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 diff --git a/polaris/ocean/suites/geostrophic.txt b/polaris/ocean/suites/geostrophic.txt index a8a573e24..7a88f6da9 100644 --- a/polaris/ocean/suites/geostrophic.txt +++ b/polaris/ocean/suites/geostrophic.txt @@ -1,2 +1,2 @@ -ocean/spherical/icos/geostrophic -ocean/spherical/qu/geostrophic +ocean/spherical/icos/geostrophic/convergence_both +ocean/spherical/qu/geostrophic/convergence_both diff --git a/polaris/ocean/suites/nightly.txt b/polaris/ocean/suites/nightly.txt index 725b504e1..e2b2ef4aa 100644 --- a/polaris/ocean/suites/nightly.txt +++ b/polaris/ocean/suites/nightly.txt @@ -3,5 +3,5 @@ ocean/planar/baroclinic_channel/10km/decomp ocean/planar/baroclinic_channel/10km/restart ocean/planar/ice_shelf_2d/5km/z-star/default/with_restart ocean/planar/ice_shelf_2d/5km/z-level/default/with_restart -ocean/planar/inertial_gravity_wave +ocean/planar/inertial_gravity_wave/convergence_both # ocean/planar/manufactured_solution diff --git a/polaris/ocean/suites/pr.txt b/polaris/ocean/suites/pr.txt index 3d63ada67..4bb3cf52a 100644 --- a/polaris/ocean/suites/pr.txt +++ b/polaris/ocean/suites/pr.txt @@ -3,7 +3,7 @@ ocean/planar/baroclinic_channel/10km/decomp ocean/planar/baroclinic_channel/10km/restart ocean/planar/ice_shelf_2d/5km/z-star/default/with_restart ocean/planar/ice_shelf_2d/5km/z-level/default/with_restart -ocean/planar/inertial_gravity_wave +ocean/planar/inertial_gravity_wave/convergence_both ocean/planar/internal_wave/standard/default ocean/planar/internal_wave/vlr/default # ocean/planar/manufactured_solution From e46197a5f69f4482e2a06d5219025f41feef7ead Mon Sep 17 00:00:00 2001 From: Carolyn Begeman Date: Thu, 10 Oct 2024 11:11:47 -0500 Subject: [PATCH 5/9] Update docs for new convergence tests --- docs/developers_guide/ocean/api.md | 27 ++++++++++++--------- docs/users_guide/ocean/tasks/geostrophic.md | 15 ++++++++---- 2 files changed, 25 insertions(+), 17 deletions(-) diff --git a/docs/developers_guide/ocean/api.md b/docs/developers_guide/ocean/api.md index 99da871be..c1f2fa9fa 100644 --- a/docs/developers_guide/ocean/api.md +++ b/docs/developers_guide/ocean/api.md @@ -313,18 +313,21 @@ .. 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 + get_resolution_for_tasks + get_timestep_for_tasks + + forward.ConvergenceForward + forward.ConvergenceForward.compute_cell_count + forward.ConvergenceForward.dynamic_model_config + + analysis.ConvergenceAnalysis + analysis.ConvergenceAnalysis.compute_error + analysis.ConvergenceAnalysis.convergence_parameters + analysis.ConvergenceAnalysis.exact_solution + analysis.ConvergenceAnalysis.get_output_field + analysis.ConvergenceAnalysis.plot_convergence + analysis.ConvergenceAnalysis.run + analysis.ConvergenceAnalysis.setup ``` ### Spherical Convergence Tests diff --git a/docs/users_guide/ocean/tasks/geostrophic.md b/docs/users_guide/ocean/tasks/geostrophic.md index 82301a154..614eeacce 100644 --- a/docs/users_guide/ocean/tasks/geostrophic.md +++ b/docs/users_guide/ocean/tasks/geostrophic.md @@ -4,20 +4,25 @@ ## description -The `geostrophic` and `geostrophic/with_viz` tasks implement the "Global Steady +The `geostrophic` tasks implement the "Global Steady State Nonlinear Zonal Geostrophic Flow" test case described in [Williamson et al. 1992]() -The task is a convergence test with time step varying proportionately to grid -size. The result of the `analysis` step of the task are plots like the -following showing convergence of water-column thickness and normal velocity as -functions of the mesh resolution: +The task `geostrophic/convergence_both` is a space-and-time convergence test +with time step varying proportionately to grid size. Convergence tests in +either space or time are also available as `convergence_space` or +`convergence_time`. The result of the `analysis` step of the task are plots +like the following showing convergence of water-column thickness and normal +velocity as functions of the mesh resolution: ```{image} images/geostrophic_convergence.png :align: center :width: 500 px ``` +Visualizations of the fields themselves can be added in viz steps by appending +`with_viz` to the task name at setup. + ## suppported models These tasks support only MPAS-Ocean. From f39c33db7b36fefeb0f816a4163aa1d6d8e8ebf8 Mon Sep 17 00:00:00 2001 From: Carolyn Begeman Date: Thu, 10 Oct 2024 21:18:28 -0500 Subject: [PATCH 6/9] Update sphere_transport for convergence --- .../ocean/tasks/sphere_transport/__init__.py | 54 ++++++++++--------- .../sphere_transport/filament_analysis.py | 35 ++++++------ .../ocean/tasks/sphere_transport/forward.yaml | 2 + .../tasks/sphere_transport/mixing_analysis.py | 35 ++++++------ 4 files changed, 70 insertions(+), 56 deletions(-) diff --git a/polaris/ocean/tasks/sphere_transport/__init__.py b/polaris/ocean/tasks/sphere_transport/__init__.py index 3596c5e95..e1aabb1e1 100644 --- a/polaris/ocean/tasks/sphere_transport/__init__.py +++ b/polaris/ocean/tasks/sphere_transport/__init__.py @@ -171,28 +171,28 @@ def _setup_steps(self, refinement): for idx, refinement_factor in enumerate(refinement_factors): resolution = get_resolution_for_task( config, refinement_factor, refinement=refinement) - if resolution not in self.resolutions: + base_mesh_step, mesh_name = add_spherical_base_mesh_step( + component, resolution, icosahedral) + analysis_dependencies['mesh'][refinement_factor] = base_mesh_step + + name = f'{prefix}_init_{mesh_name}' + subdir = f'{sph_trans_dir}/init/{mesh_name}' + if self.include_viz: + symlink = f'init/{mesh_name}' + else: + symlink = None + if subdir in component.steps: + init_step = component.steps[subdir] + else: + init_step = Init(component=component, name=name, + subdir=subdir, base_mesh=base_mesh_step, + case_name=case_name) + init_step.set_shared_config(config, link=config_filename) + analysis_dependencies['init'][refinement_factor] = init_step + self.add_step(base_mesh_step, symlink=f'base_mesh/{mesh_name}') + self.add_step(init_step, symlink=symlink) + if resolution not in resolutions: resolutions.append(resolution) - base_mesh_step, mesh_name = add_spherical_base_mesh_step( - component, resolution, icosahedral) - self.add_step(base_mesh_step, symlink=f'base_mesh/{mesh_name}') - analysis_dependencies['mesh'][resolution] = base_mesh_step - - name = f'{prefix}_init_{mesh_name}' - subdir = f'{sph_trans_dir}/init/{mesh_name}' - if self.include_viz: - symlink = f'init/{mesh_name}' - else: - symlink = None - if subdir in component.steps: - init_step = component.steps[subdir] - else: - init_step = Init(component=component, name=name, - subdir=subdir, base_mesh=base_mesh_step, - case_name=case_name) - init_step.set_shared_config(config, link=config_filename) - self.add_step(init_step, symlink=symlink) - analysis_dependencies['init'][resolution] = init_step timestep, _ = get_timestep_for_task( config, refinement_factor, refinement=refinement) @@ -217,7 +217,7 @@ def _setup_steps(self, refinement): refinement=refinement) forward_step.set_shared_config(config, link=config_filename) self.add_step(forward_step, symlink=symlink) - analysis_dependencies['forward'][resolution] = forward_step + analysis_dependencies['forward'][refinement_factor] = forward_step if self.include_viz: with_viz_dir = f'{sph_trans_dir}/with_viz' @@ -257,10 +257,11 @@ def _setup_steps(self, refinement): step = component.steps[subdir] else: step = MixingAnalysis(component=component, - resolutions=resolutions, icosahedral=icosahedral, subdir=subdir, + refinement_factors=refinement_factors, case_name=case_name, - dependencies=analysis_dependencies) + dependencies=analysis_dependencies, + refinement=refinement) step.set_shared_config(config, link=config_filename) self.add_step(step, symlink=symlink) @@ -274,9 +275,10 @@ def _setup_steps(self, refinement): step = component.steps[subdir] else: step = FilamentAnalysis(component=component, - resolutions=resolutions, + refinement_factors=refinement_factors, icosahedral=icosahedral, subdir=subdir, case_name=case_name, - dependencies=analysis_dependencies) + dependencies=analysis_dependencies, + refinement=refinement) step.set_shared_config(config, link=config_filename) self.add_step(step, symlink=symlink) diff --git a/polaris/ocean/tasks/sphere_transport/filament_analysis.py b/polaris/ocean/tasks/sphere_transport/filament_analysis.py index a3f196439..f5a3d6895 100644 --- a/polaris/ocean/tasks/sphere_transport/filament_analysis.py +++ b/polaris/ocean/tasks/sphere_transport/filament_analysis.py @@ -5,6 +5,7 @@ from polaris import Step from polaris.mpas import time_index_from_xtime +from polaris.ocean.convergence import get_resolution_for_task from polaris.ocean.resolution import resolution_to_subdir from polaris.viz import use_mplstyle @@ -25,8 +26,8 @@ class FilamentAnalysis(Step): case_name : str The name of the test case """ - def __init__(self, component, resolutions, icosahedral, subdir, - case_name, dependencies): + def __init__(self, component, refinement_factors, icosahedral, subdir, + case_name, dependencies, refinement='both'): """ Create the step @@ -53,22 +54,22 @@ def __init__(self, component, resolutions, icosahedral, subdir, """ super().__init__(component=component, name='filament_analysis', subdir=subdir) - self.resolutions = resolutions + self.refinement_factors = refinement_factors + self.refinement = refinement self.case_name = case_name - for resolution in resolutions: - mesh_name = resolution_to_subdir(resolution) - base_mesh = dependencies['mesh'][resolution] - init = dependencies['init'][resolution] - forward = dependencies['forward'][resolution] + for refinement_factor in refinement_factors: + base_mesh = dependencies['mesh'][refinement_factor] + init = dependencies['init'][refinement_factor] + forward = dependencies['forward'][refinement_factor] self.add_input_file( - filename=f'{mesh_name}_mesh.nc', + filename=f'mesh_r{refinement_factor:02g}.nc', work_dir_target=f'{base_mesh.path}/base_mesh.nc') self.add_input_file( - filename=f'{mesh_name}_init.nc', + filename=f'init_r{refinement_factor:02g}.nc', work_dir_target=f'{init.path}/initial_state.nc') self.add_input_file( - filename=f'{mesh_name}_output.nc', + filename=f'output_r{refinement_factor:02g}.nc', work_dir_target=f'{forward.path}/output.nc') self.add_output_file('filament.png') @@ -77,7 +78,11 @@ def run(self): Run this step of the test case """ plt.switch_backend('Agg') - resolutions = self.resolutions + resolutions = list() + for refinement_factor in self.refinement_factors: + resolution = get_resolution_for_task( + self.config, refinement_factor, self.refinement) + resolutions.append(resolution) config = self.config section = config[self.case_name] eval_time = section.getfloat('filament_evaluation_time') @@ -89,9 +94,9 @@ def run(self): filament_norm = np.zeros((len(resolutions), num_tau)) use_mplstyle() fig, ax = plt.subplots() - for i, resolution in enumerate(resolutions): - mesh_name = resolution_to_subdir(resolution) - ds = xr.open_dataset(f'{mesh_name}_output.nc') + for i, refinement_factor in enumerate(self.refinement_factors): + mesh_name = resolution_to_subdir(resolutions[i]) + ds = xr.open_dataset(f'output_r{refinement_factor:02g}.nc') tidx = time_index_from_xtime(ds.xtime.values, eval_time * s_per_day) tracer = ds[variable_name] diff --git a/polaris/ocean/tasks/sphere_transport/forward.yaml b/polaris/ocean/tasks/sphere_transport/forward.yaml index e8811b654..be8a53ec0 100644 --- a/polaris/ocean/tasks/sphere_transport/forward.yaml +++ b/polaris/ocean/tasks/sphere_transport/forward.yaml @@ -49,6 +49,8 @@ mpas-ocean: - tracers - mesh - xtime + - velocityZonal + - velocityMeridional - normalVelocity - layerThickness - refZMid diff --git a/polaris/ocean/tasks/sphere_transport/mixing_analysis.py b/polaris/ocean/tasks/sphere_transport/mixing_analysis.py index 750518a7f..88ee06702 100644 --- a/polaris/ocean/tasks/sphere_transport/mixing_analysis.py +++ b/polaris/ocean/tasks/sphere_transport/mixing_analysis.py @@ -7,6 +7,7 @@ from polaris import Step from polaris.mpas import time_index_from_xtime +from polaris.ocean.convergence import get_resolution_for_task from polaris.ocean.resolution import resolution_to_subdir from polaris.viz import use_mplstyle @@ -27,8 +28,8 @@ class MixingAnalysis(Step): case_name : str The name of the test case """ - def __init__(self, component, resolutions, icosahedral, subdir, - case_name, dependencies): + def __init__(self, component, refinement_factors, icosahedral, subdir, + case_name, dependencies, refinement='both'): """ Create the step @@ -55,22 +56,22 @@ def __init__(self, component, resolutions, icosahedral, subdir, """ super().__init__(component=component, name='mixing_analysis', subdir=subdir) - self.resolutions = resolutions + self.refinement_factors = refinement_factors + self.refinement = refinement self.case_name = case_name - for resolution in resolutions: - mesh_name = resolution_to_subdir(resolution) - base_mesh = dependencies['mesh'][resolution] - init = dependencies['init'][resolution] - forward = dependencies['forward'][resolution] + for refinement_factor in refinement_factors: + base_mesh = dependencies['mesh'][refinement_factor] + init = dependencies['init'][refinement_factor] + forward = dependencies['forward'][refinement_factor] self.add_input_file( - filename=f'{mesh_name}_mesh.nc', + filename=f'mesh_r{refinement_factor:02g}.nc', work_dir_target=f'{base_mesh.path}/base_mesh.nc') self.add_input_file( - filename=f'{mesh_name}_init.nc', + filename=f'init_r{refinement_factor:02g}.nc', work_dir_target=f'{init.path}/initial_state.nc') self.add_input_file( - filename=f'{mesh_name}_output.nc', + filename=f'output_r{refinement_factor:02g}.nc', work_dir_target=f'{forward.path}/output.nc') self.add_output_file('triplots.png') @@ -79,7 +80,11 @@ def run(self): Run this step of the test case """ plt.switch_backend('Agg') - resolutions = self.resolutions + resolutions = list() + for refinement_factor in self.refinement_factors: + resolution = get_resolution_for_task( + self.config, refinement_factor, self.refinement) + resolutions.append(resolution) config = self.config section = config[self.case_name] eval_time = section.getfloat('mixing_evaluation_time') @@ -89,12 +94,12 @@ def run(self): use_mplstyle() fig, axes = plt.subplots(nrows=nrows, ncols=2, sharex=True, sharey=True, figsize=(5.5, 7)) - for i, resolution in enumerate(resolutions): + for i, refinement_factor in enumerate(self.refinement_factors): ax = axes[int(i / 2), i % 2] _init_triplot_axes(ax) - mesh_name = resolution_to_subdir(resolution) + mesh_name = resolution_to_subdir(resolutions[i]) ax.set(title=mesh_name) - ds = xr.open_dataset(f'{mesh_name}_output.nc') + ds = xr.open_dataset(f'output_r{refinement_factor:02g}.nc') if i % 2 == 0: ax.set_ylabel("tracer3") if int(i / 2) == nrows - 1: From fa5d6c148e414be358ae9a3369de91c350013dcc Mon Sep 17 00:00:00 2001 From: Carolyn Begeman Date: Tue, 15 Oct 2024 14:56:25 -0500 Subject: [PATCH 7/9] fixup imports --- polaris/ocean/tasks/cosine_bell/__init__.py | 1 - 1 file changed, 1 deletion(-) diff --git a/polaris/ocean/tasks/cosine_bell/__init__.py b/polaris/ocean/tasks/cosine_bell/__init__.py index 34e527211..6e26b2498 100644 --- a/polaris/ocean/tasks/cosine_bell/__init__.py +++ b/polaris/ocean/tasks/cosine_bell/__init__.py @@ -8,7 +8,6 @@ get_timestep_for_task, ) from polaris.ocean.mesh.spherical import add_spherical_base_mesh_step -from polaris.ocean.resolution import resolution_to_subdir from polaris.ocean.tasks.cosine_bell.analysis import Analysis from polaris.ocean.tasks.cosine_bell.forward import Forward from polaris.ocean.tasks.cosine_bell.init import Init From 97ff0dd0302a4a28a3c0e42e1e666ab4c0378c45 Mon Sep 17 00:00:00 2001 From: Carolyn Begeman Date: Tue, 12 Nov 2024 22:30:13 -0600 Subject: [PATCH 8/9] Separate refinement factors config options by refinement type --- polaris/ocean/convergence/__init__.py | 18 ++++++++++------- polaris/ocean/convergence/analysis.py | 20 +++++++++++++------ polaris/ocean/convergence/convergence.cfg | 14 ++++++++----- polaris/ocean/convergence/forward.py | 2 ++ .../ocean/convergence/spherical/spherical.cfg | 6 ++++-- polaris/ocean/tasks/cosine_bell/__init__.py | 13 ++++++------ polaris/ocean/tasks/geostrophic/__init__.py | 13 ++++++------ .../tasks/inertial_gravity_wave/__init__.py | 11 +++++++++- .../inertial_gravity_wave.cfg | 3 ++- .../tasks/manufactured_solution/__init__.py | 11 +++++++++- .../manufactured_solution.cfg | 3 ++- .../ocean/tasks/sphere_transport/__init__.py | 15 +++++++------- 12 files changed, 86 insertions(+), 43 deletions(-) diff --git a/polaris/ocean/convergence/__init__.py b/polaris/ocean/convergence/__init__.py index a66abb2e4..0a04e50e7 100644 --- a/polaris/ocean/convergence/__init__.py +++ b/polaris/ocean/convergence/__init__.py @@ -23,10 +23,12 @@ def get_resolution_for_task(config, refinement_factor, The resolution corresponding to the refinement_factor and convergence test type """ + if refinement == 'both': + option = 'refinement_factors_space' + else: + option = f'refinement_factors_{refinement}' base_resolution = config.getfloat('convergence', 'base_resolution') - refinement_factors = config.getlist('convergence', - 'refinement_factors', - dtype=float) + refinement_factors = config.getlist('convergence', option, dtype=float) if refinement_factor not in refinement_factors: raise ValueError( @@ -63,10 +65,12 @@ def get_timestep_for_task(config, refinement_factor, test type """ + if refinement == 'both': + option = 'refinement_factors_space' + else: + option = f'refinement_factors_{refinement}' base_resolution = config.getfloat('convergence', 'base_resolution') - refinement_factors = config.getlist('convergence', - 'refinement_factors', - dtype=float) + refinement_factors = config.getlist('convergence', option, dtype=float) if refinement_factor not in refinement_factors: raise ValueError( @@ -93,7 +97,7 @@ def get_timestep_for_task(config, refinement_factor, else: btr_timestep = btr_dt_per_km * resolution if refinement == 'time': - timestep = dt_per_km * refinement_factor * base_resolution + timestep = dt_per_km * refinement_factor * resolution elif refinement == 'space': timestep = dt_per_km * base_resolution else: diff --git a/polaris/ocean/convergence/analysis.py b/polaris/ocean/convergence/analysis.py index f44f719e1..54d3b5ec4 100644 --- a/polaris/ocean/convergence/analysis.py +++ b/polaris/ocean/convergence/analysis.py @@ -5,11 +5,11 @@ import pandas as pd from polaris.mpas import area_for_field, time_index_from_xtime -from polaris.ocean.model import OceanIOStep from polaris.ocean.convergence import ( get_resolution_for_task, get_timestep_for_task, ) +from polaris.ocean.model import OceanIOStep from polaris.viz import use_mplstyle @@ -123,8 +123,11 @@ def setup(self): """ config = self.config dependencies = self.dependencies_dict - refinement_factors = config.getlist('convergence', - 'refinement_factors', + if self.refinement == 'time': + option = 'refinement_factors_time' + else: + option = 'refinement_factors_space' + refinement_factors = config.getlist('convergence', option, dtype=float) for refinement_factor in refinement_factors: base_mesh = dependencies['mesh'][refinement_factor] @@ -175,8 +178,13 @@ def plot_convergence(self, variable_name, title, zidx): field_name=variable_name) error = [] - refinement_factors = config.getlist('convergence', - 'refinement_factors', + if self.refinement == 'time': + option = 'refinement_factors_time' + header = 'time step' + else: + option = 'refinement_factors_space' + header = 'resolution' + refinement_factors = config.getlist('convergence', option, dtype=float) resolutions = list() timesteps = list() @@ -203,7 +211,7 @@ def plot_convergence(self, variable_name, title, zidx): error_array = np.array(error) filename = f'convergence_{variable_name}.csv' data = np.stack((refinement_array, error_array), axis=1) - df = pd.DataFrame(data, columns=['resolution', error_type]) + df = pd.DataFrame(data, columns=[header, error_type]) df.to_csv(f'convergence_{variable_name}.csv', index=False) convergence_failed = False diff --git a/polaris/ocean/convergence/convergence.cfg b/polaris/ocean/convergence/convergence.cfg index f93d368c4..9957bc0d2 100644 --- a/polaris/ocean/convergence/convergence.cfg +++ b/polaris/ocean/convergence/convergence.cfg @@ -4,19 +4,23 @@ convergence_eval_time = 24.0 # Convergence threshold below which a test fails -convergence_thresh = 1.0 +convergence_thresh_space = 1.0 +convergence_thresh_time = 1.0 # Type of error to compute error_type = l2 -# the base mesh resolution (km) to which refinement_factors -# are applied if refinement is 'space' or 'both' on a planar mesh +# the base mesh resolution (km) to which refinement_factors are applied # base resolutions for spherical meshes are given in section spherical_convergence base_resolution = 120 -# refinement factors for a planar mesh applied to either space or time +# refinement factors for a planar mesh applied to either space or both space and time # refinement factors for a spherical mesh given in section spherical_convergence -refinement_factors = 4., 2., 1., 0.5 +refinement_factors_space = 4., 2., 1., 0.5 + +# refinement factors for a planar mesh applied to time with the base timestep +# determined by base_resolution * dt_per_km +refinement_factors_time = 1., 0.5, 0.25 # config options for convergence forward steps [convergence_forward] diff --git a/polaris/ocean/convergence/forward.py b/polaris/ocean/convergence/forward.py index 1ea1c47e2..7a7ecbc2e 100644 --- a/polaris/ocean/convergence/forward.py +++ b/polaris/ocean/convergence/forward.py @@ -132,6 +132,7 @@ def dynamic_model_config(self, at_setup): btr_dt_str = get_time_interval_string(seconds=btr_timestep) s_per_hour = 3600. + section = config['convergence_forward'] run_duration = section.getfloat('run_duration') run_duration_str = get_time_interval_string( seconds=run_duration * s_per_hour) @@ -140,6 +141,7 @@ def dynamic_model_config(self, at_setup): output_interval_str = get_time_interval_string( seconds=output_interval * s_per_hour) + time_integrator = section.get('time_integrator') time_integrator_map = dict([('RK4', 'RungeKutta4')]) model = config.get('ocean', 'model') if model == 'omega': diff --git a/polaris/ocean/convergence/spherical/spherical.cfg b/polaris/ocean/convergence/spherical/spherical.cfg index fd769f8ba..cb4b8d4fd 100644 --- a/polaris/ocean/convergence/spherical/spherical.cfg +++ b/polaris/ocean/convergence/spherical/spherical.cfg @@ -6,11 +6,13 @@ icos_base_resolution = 60. # The factors by which to scale space or time based on icos_base_resolution -icos_refinement_factors = 8., 4., 2., 1. +icos_refinement_factors_space = 8., 4., 2., 1. +icos_refinement_factors_time = 1., 0.5, 0.25 # The base resolution for the quasi-uniform mesh to which the refinement # factors are applied qu_base_resolution = 120. # The factors by which to scale space or time based on qu_base_resolution -qu_refinement_factors = 0.5, 0.75, 1., 1.25, 1.5, 1.75, 2. +qu_refinement_factors_space = 2., 1.75, 1.5, 1.25, 1., 0.75, 0.5 +qu_refinement_factors_time = 1., 0.5, 0.25 diff --git a/polaris/ocean/tasks/cosine_bell/__init__.py b/polaris/ocean/tasks/cosine_bell/__init__.py index 6e26b2498..a26af8d35 100644 --- a/polaris/ocean/tasks/cosine_bell/__init__.py +++ b/polaris/ocean/tasks/cosine_bell/__init__.py @@ -129,15 +129,16 @@ def _setup_steps(self, refinement): else: prefix = 'qu' + if refinement == 'time': + option = 'refinement_factors_time' + else: + option = 'refinement_factors_space' refinement_factors = config.getlist('spherical_convergence', - f'{prefix}_refinement_factors', - dtype=str) + f'{prefix}_{option}', dtype=str) refinement_factors = ', '.join(refinement_factors) - config.set('convergence', 'refinement_factors', - value=refinement_factors) + config.set('convergence', option, value=refinement_factors) refinement_factors = config.getlist('convergence', - 'refinement_factors', - dtype=float) + option, dtype=float) base_resolution = config.getfloat('spherical_convergence', f'{prefix}_base_resolution') diff --git a/polaris/ocean/tasks/geostrophic/__init__.py b/polaris/ocean/tasks/geostrophic/__init__.py index a3b2e941c..f5cb1ca72 100644 --- a/polaris/ocean/tasks/geostrophic/__init__.py +++ b/polaris/ocean/tasks/geostrophic/__init__.py @@ -130,15 +130,16 @@ def _setup_steps(self, refinement): else: prefix = 'qu' + if refinement == 'time': + option = 'refinement_factors_time' + else: + option = 'refinement_factors_space' refinement_factors = config.getlist('spherical_convergence', - f'{prefix}_refinement_factors', - dtype=str) + f'{prefix}_{option}', dtype=str) refinement_factors = ', '.join(refinement_factors) - config.set('convergence', 'refinement_factors', - value=refinement_factors) + config.set('convergence', option, value=refinement_factors) refinement_factors = config.getlist('convergence', - 'refinement_factors', - dtype=float) + option, dtype=float) base_resolution = config.getfloat('spherical_convergence', f'{prefix}_base_resolution') diff --git a/polaris/ocean/tasks/inertial_gravity_wave/__init__.py b/polaris/ocean/tasks/inertial_gravity_wave/__init__.py index 69f9c6a54..f9b2de3a0 100644 --- a/polaris/ocean/tasks/inertial_gravity_wave/__init__.py +++ b/polaris/ocean/tasks/inertial_gravity_wave/__init__.py @@ -66,8 +66,12 @@ def __init__(self, component, config, refinement='both'): analysis_dependencies: Dict[str, Dict[float, Step]] = ( dict(mesh=dict(), init=dict(), forward=dict())) + if refinement == 'time': + option = 'refinement_factors_time' + else: + option = 'refinement_factors_space' refinement_factors = self.config.getlist( - 'convergence', 'refinement_factors', dtype=float) + 'convergence', option, dtype=float) timesteps = list() resolutions = list() for refinement_factor in refinement_factors: @@ -119,3 +123,8 @@ def __init__(self, component, config, refinement='both'): self.add_step(Viz(component=component, resolutions=resolutions, taskdir=self.subdir), run_by_default=False) + config.add_from_package('polaris.ocean.convergence', + 'convergence.cfg') + config.add_from_package( + 'polaris.ocean.tasks.inertial_gravity_wave', + config_filename) 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 b690a3de1..aebb28655 100644 --- a/polaris/ocean/tasks/inertial_gravity_wave/inertial_gravity_wave.cfg +++ b/polaris/ocean/tasks/inertial_gravity_wave/inertial_gravity_wave.cfg @@ -62,7 +62,8 @@ error_type = l2 base_resolution = 100. # refinement factors for a planar mesh applied to either space or time -refinement_factors = 2., 1., 0.5, 0.25 +refinement_factors_space = 2., 1., 0.5, 0.25 +refinement_factors_time = 1., 0.5, 0.25 # config options for spherical convergence tests [convergence_forward] diff --git a/polaris/ocean/tasks/manufactured_solution/__init__.py b/polaris/ocean/tasks/manufactured_solution/__init__.py index 7471602c2..f85917698 100644 --- a/polaris/ocean/tasks/manufactured_solution/__init__.py +++ b/polaris/ocean/tasks/manufactured_solution/__init__.py @@ -26,6 +26,11 @@ def add_manufactured_solution_tasks(component): config_filename = 'manufactured_solution.cfg' filepath = f'{basedir}/{config_filename}' config = PolarisConfigParser(filepath=filepath) + config.add_from_package('polaris.ocean.convergence', + 'convergence.cfg') + config.add_from_package( + 'polaris.ocean.tasks.manufactured_solution', + config_filename) for refinement in ['space', 'time', 'both']: component.add_task(ManufacturedSolution(component=component, config=config, @@ -62,8 +67,12 @@ def __init__(self, component, config, refinement='both'): analysis_dependencies: Dict[str, Dict[float, Step]] = ( dict(mesh=dict(), init=dict(), forward=dict())) + if refinement == 'time': + option = 'refinement_factors_time' + else: + option = 'refinement_factors_space' refinement_factors = self.config.getlist( - 'convergence', 'refinement_factors', dtype=float) + 'convergence', option, dtype=float) timesteps = list() resolutions = list() for refinement_factor in refinement_factors: diff --git a/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg b/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg index be89113c0..a4edc7e49 100644 --- a/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg +++ b/polaris/ocean/tasks/manufactured_solution/manufactured_solution.cfg @@ -70,7 +70,8 @@ error_type = l2 base_resolution = 50. # refinement factors for a planar mesh applied to either space or time -refinement_factors = 4., 2., 1., 0.5 +refinement_factors_space = 4., 2., 1., 0.5 +refinement_factors_time = 1., 0.5, 0.25 # config options for spherical convergence tests [convergence_forward] diff --git a/polaris/ocean/tasks/sphere_transport/__init__.py b/polaris/ocean/tasks/sphere_transport/__init__.py index e1aabb1e1..599bedf7b 100644 --- a/polaris/ocean/tasks/sphere_transport/__init__.py +++ b/polaris/ocean/tasks/sphere_transport/__init__.py @@ -128,7 +128,7 @@ def configure(self): # set up the steps again in case a user has provided new resolutions self._setup_steps(self.refinement) - def _setup_steps(self, refinement): + def _setup_steps(self, refinement): # noqa: C901 """ setup steps given resolutions """ case_name = self.case_name icosahedral = self.icosahedral @@ -140,15 +140,16 @@ def _setup_steps(self, refinement): else: prefix = 'qu' + if refinement == 'time': + option = 'refinement_factors_time' + else: + option = 'refinement_factors_space' refinement_factors = config.getlist('spherical_convergence', - f'{prefix}_refinement_factors', - dtype=str) + f'{prefix}_{option}', dtype=str) refinement_factors = ', '.join(refinement_factors) - config.set('convergence', 'refinement_factors', - value=refinement_factors) + config.set('convergence', option, value=refinement_factors) refinement_factors = config.getlist('convergence', - 'refinement_factors', - dtype=float) + option, dtype=float) base_resolution = config.getfloat('spherical_convergence', f'{prefix}_base_resolution') From 9fe597a27bee2f470afa13b8b266cee44af7724e Mon Sep 17 00:00:00 2001 From: Carolyn Begeman Date: Mon, 25 Nov 2024 10:33:57 -0600 Subject: [PATCH 9/9] Change task list order --- polaris/ocean/tasks/cosine_bell/__init__.py | 4 ++-- polaris/ocean/tasks/geostrophic/__init__.py | 4 ++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/polaris/ocean/tasks/cosine_bell/__init__.py b/polaris/ocean/tasks/cosine_bell/__init__.py index a26af8d35..07b4e20bc 100644 --- a/polaris/ocean/tasks/cosine_bell/__init__.py +++ b/polaris/ocean/tasks/cosine_bell/__init__.py @@ -33,8 +33,8 @@ def add_cosine_bell_tasks(component): config.add_from_package('polaris.ocean.tasks.cosine_bell', 'cosine_bell.cfg') - for include_viz in [False, True]: - for refinement in ['space', 'time', 'both']: + for refinement in ['space', 'time', 'both']: + for include_viz in [False, True]: component.add_task(CosineBell(component=component, config=config, icosahedral=icosahedral, diff --git a/polaris/ocean/tasks/geostrophic/__init__.py b/polaris/ocean/tasks/geostrophic/__init__.py index f5cb1ca72..9d67179ce 100644 --- a/polaris/ocean/tasks/geostrophic/__init__.py +++ b/polaris/ocean/tasks/geostrophic/__init__.py @@ -33,8 +33,8 @@ def add_geostrophic_tasks(component): config.add_from_package('polaris.ocean.tasks.geostrophic', 'geostrophic.cfg') - for include_viz in [False, True]: - for refinement in ['space', 'time', 'both']: + for refinement in ['space', 'time', 'both']: + for include_viz in [False, True]: component.add_task(Geostrophic(component=component, config=config, icosahedral=icosahedral,