Source code for compass.ocean.tests.overflow.initial_state

import numpy as np
import xarray
from mpas_tools.io import write_netcdf
from mpas_tools.mesh.conversion import convert, cull
from mpas_tools.planar_hex import make_planar_hex_mesh

from compass.ocean.vertical import init_vertical_coord
from compass.step import Step


[docs] class InitialState(Step): """ A step for creating a mesh and initial condition for the nonhydro and the hydro vs nonhydro test cases. """
[docs] def __init__(self, test_case): """ Create the step Parameters ---------- test_case : compass.TestCase """ super().__init__(test_case=test_case, name='initial_state') for file in ['base_mesh.nc', 'culled_mesh.nc', 'culled_graph.info', 'initial_state.nc']: self.add_output_file(file)
[docs] def run(self): """ Run this step of the test case """ config = self.config logger = self.logger section = config['horizontal_grid'] nx = section.getint('nx') ny = section.getint('ny') dc = section.getfloat('dc') dsMesh = make_planar_hex_mesh(nx=nx, ny=ny, dc=dc, nonperiodic_x=True, nonperiodic_y=False) write_netcdf(dsMesh, 'base_mesh.nc') dsMesh = cull(dsMesh, logger=logger) dsMesh = convert(dsMesh, graphInfoFileName='culled_graph.info', logger=logger) write_netcdf(dsMesh, 'culled_mesh.nc') section = config['vertical_grid'] nVertLevels = section.getint('vert_levels') bottom_depth = section.getfloat('bottom_depth') section = config['hydro_vs_nonhydro'] cold_water_range = section.getfloat('cold_water_range') shelf_depth = section.getfloat('shelf_depth') xs = section.getfloat('xs') Ls = section.getfloat('Ls') config_eos_linear_Sref = section.getfloat('eos_linear_Sref') lower_temperature = section.getfloat('lower_temperature') higher_temperature = section.getfloat('higher_temperature') # obtain dimensions and mesh variables # vertical_coordinate = 'uniform' ds = dsMesh.copy() nCells = ds.nCells.size nEdges = ds.nEdges.size nVertices = ds.nVertices.size xCell = ds.xCell # create and initialize variables # bottom depth ds['bottomDepth'] = - (- bottom_depth + 0.5 * (bottom_depth - shelf_depth) * (1.0 + np.tanh((nx * dc - ds.xCell - xs) / Ls))) # ssh ds['ssh'] = xarray.zeros_like(xCell) init_vertical_coord(config, ds) # initial salinity and temperature ds['salinity'] = (config_eos_linear_Sref * xarray.ones_like(ds.zMid)).where(ds.cellMask) # T = Tref - (rho - rhoRef)/alpha ds['temperature'] = xarray.ones_like(ds.zMid).where(ds.cellMask) temperature = ds['temperature'] for iCell in range(0, nCells): temperature[0, iCell, :] = -1.0 for k in range(0, nVertLevels): if (xCell[iCell] < cold_water_range): temperature[0, iCell, k] = lower_temperature else: temperature[0, iCell, k] = higher_temperature # initial velocity on edges ds['normalVelocity'] = (('Time', 'nEdges', 'nVertLevels',), np.zeros([1, nEdges, nVertLevels])) # Coriolis parameter ds['fCell'] = (('nCells', 'nVertLevels',), np.zeros([nCells, nVertLevels])) ds['fEdge'] = (('nEdges', 'nVertLevels',), np.zeros([nEdges, nVertLevels])) ds['fVertex'] = (('nVertices', 'nVertLevels',), np.zeros([nVertices, nVertLevels])) # finalize and write file # If you prefer not to have NaN as the fill value, you should consider # using mpas_tools.io.write_netcdf() instead write_netcdf(ds, 'initial_state.nc')