Source code for polaris.ocean.tasks.geostrophic.analysis

import xarray as xr

from polaris.ocean.convergence import ConvergenceAnalysis
from polaris.ocean.tasks.geostrophic.exact_solution import (
    compute_exact_solution,
)


[docs]class Analysis(ConvergenceAnalysis): """ A step for analyzing the output from the geostrophic test case """
[docs] def __init__(self, component, resolutions, subdir, dependencies): """ Create the step Parameters ---------- 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 dependencies : dict of dict of polaris.Steps The dependencies of this step """ convergence_vars = [{'name': 'h', 'title': 'water-column thickness', 'zidx': None}, {'name': 'normalVelocity', 'title': 'normal velocity', 'zidx': 0}] super().__init__(component=component, subdir=subdir, resolutions=resolutions, dependencies=dependencies, convergence_vars=convergence_vars)
[docs] def exact_solution(self, mesh_name, 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 field_name : str The name of the variable of which we evaluate convergence For the default method, we use the same convergence rate for all fields time : float The time at which to evaluate the exact solution in seconds. For the default method, we always use the initial state. zidx : int, optional The z-index for the vertical level at which to evaluate the exact solution Returns ------- solution: xarray.DataArray The exact solution with dimension nCells """ if field_name not in ['h', 'normalVelocity']: print(f'Variable {field_name} not available as an analytic ' 'solution for the geostrophic test case') config = self.config section = config['geostrophic'] alpha = section.getfloat('alpha') vel_period = section.getfloat('vel_period') gh_0 = section.getfloat('gh_0') mesh_filename = f'{mesh_name}_mesh.nc' h, _, _, normalVelocity = compute_exact_solution( alpha, vel_period, gh_0, mesh_filename) if field_name == 'h': return h else: return normalVelocity
[docs] def get_output_field(self, mesh_name, field_name, time, zidx=None): """ Get the model output field at the given time and z index Parameters ---------- mesh_name : str The mesh name which is the prefix for the output file field_name : str The name of the variable of which we evaluate convergence time : float The time at which to evaluate the exact solution in seconds zidx : int, optional The z-index for the vertical level to take the field from Returns ------- field_mpas : xarray.DataArray model output field """ if field_name not in ['h', 'normalVelocity']: print(f'Variable {field_name} not available for analysis in the ' 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) else: ds_init = xr.open_dataset(f'{mesh_name}_init.nc') bottom_depth = ds_init.bottomDepth ssh = super().get_output_field(mesh_name=mesh_name, field_name='ssh', time=time, zidx=None) h = ssh + bottom_depth return h
def convergence_parameters(self, field_name=None): """ Get convergence parameters Parameters ---------- field_name : str The name of the variable of which we evaluate convergence For cosine_bell, we use the same convergence rate for all fields Returns ------- conv_thresh: float The minimum convergence rate conv_thresh: float The maximum convergence rate """ config = self.config conv_thresh = config.getfloat('geostrophic', f'convergence_thresh_{field_name}') error_type = config.get('convergence', 'error_type') return conv_thresh, error_type