Merge pull request #1630 from knutfrode/dev

Removed gls_tke functionality from OceanDrift and reader_ROMS_native,…

Author: knutfrode

opendrift/legacy.py

@@ -0,0 +1,43 @@
+# This is a collection of methods that have been removed from various places in the OpenDrift codebase.
+
+
+# See earlier version of reader_ROMS_native how the GLS parameters were obtained
+def gls_tke(windstress, depth, sea_water_density,
+            tke, generic_length_scale, gls_parameters=None):
+    '''From LADIM model, based on ROMS files.'''
+
+    g = 9.81
+    f0 = 0.1  # mean wave frequency
+    c_w = 4.0  # wave mixing parameter
+    c_i = 0.2  # coefficient for the interior
+    if gls_parameters is None:
+        # GLS parameters from ROMS, k-omega closure (see ocean.in)
+        p = 0.0
+        m = 1.0
+        n = 1.0
+        cmu0 = 0.5477  # for KANTHA_CLAYSON stability function
+    else:
+        p = gls_parameters['gls_p']
+        m = gls_parameters['gls_m']
+        n = gls_parameters['gls_n']
+        cmu0 = gls_parameters['gls_cmu0']
+
+    phi = 100. * (windstress/sea_water_density)**(3./2.)
+
+    # dissipation and turbulent length scale for interiour of mixed layer
+    eps = cmu0**(3.+p/n)*tke**(3./2.+m/n)*generic_length_scale**(-1./n)
+    l_i = c_i * tke**(3./2.) * eps**(-1.)
+
+    # diffusivity for interior of mixed layer
+    # c_i = sqrt(2.) * cmu0**3
+    ki = c_i * (2.*tke)**0.5 * l_i
+
+    # length scale and diffusivity of wave-enhanced layer
+    l_w = np.sqrt(phi / (g*f0))
+    kwave = c_w * (2*tke)**0.5 * l_w
+    kmix = ki + kwave
+
+    K, N = np.meshgrid(kmix, depths)
+
+    return K
+

opendrift/models/oceandrift.py

@@ -21,7 +21,7 @@
 import logging; logger = logging.getLogger(__name__)
 from opendrift.models.basemodel import OpenDriftSimulation
 from opendrift.elements import LagrangianArray
-from opendrift.models.physics_methods import verticaldiffusivity_Large1994, verticaldiffusivity_Sundby1983, gls_tke
+from opendrift.models.physics_methods import verticaldiffusivity_Large1994, verticaldiffusivity_Sundby1983
 from opendrift.config import CONFIG_LEVEL_ESSENTIAL, CONFIG_LEVEL_BASIC, CONFIG_LEVEL_ADVANCED
 
 # Defining the oil element properties
@@ -92,10 +92,6 @@ class OceanDrift(OpenDriftSimulation):
         'sea_surface_wind_wave_to_direction': {'fallback': 0, 'important': False},
         'sea_surface_wind_wave_mean_period': {'fallback': 0, 'important': False},
         'sea_surface_wind_wave_significant_height': {'fallback': 0, 'important': False},
-        'surface_downward_x_stress': {'fallback': 0},
-        'surface_downward_y_stress': {'fallback': 0},
-        'turbulent_kinetic_energy': {'fallback': 0},
-        'turbulent_generic_length_scale': {'fallback': 0},
         'ocean_mixed_layer_thickness': {'fallback': 50},
         'sea_floor_depth_below_sea_level': {'fallback': 10000},
         'land_binary_mask': {'fallback': None},
@@ -142,7 +138,7 @@ def __init__(self, *args, **kwargs):
                 'Time step used for inner loop of vertical mixing.'},
             'vertical_mixing:diffusivitymodel': {'type': 'enum', 'default': 'environment',
                 'enum': ['environment', 'stepfunction', 'windspeed_Sundby1983',
-                 'windspeed_Large1994', 'gls_tke','constant'], 'level': CONFIG_LEVEL_ADVANCED,
+                 'windspeed_Large1994', 'constant'], 'level': CONFIG_LEVEL_ADVANCED,
                  'units': 'seconds', 'description': 'Algorithm/source used for profile of vertical diffusivity. Environment means that diffusivity is aquired from readers or environment constants/fallback.'},
             'vertical_mixing:background_diffusivity': {'type': 'float', 'min': 0, 'max': 1, 'default': 1.2e-5,
                 'level': CONFIG_LEVEL_ADVANCED, 'units': 'm2s-1', 'description':
@@ -363,24 +359,6 @@ def get_diffusivity_profile(self, model, depths):
             return verticaldiffusivity_Large1994(wind, depth, MLD, background_diffusivity)
         elif model == 'windspeed_Sundby1983':
             return verticaldiffusivity_Sundby1983(wind, depth, MLD, background_diffusivity)
-        elif model == 'gls_tke':
-            if not hasattr(self, 'gls_parameters'):
-                logger.info('Searching readers for GLS parameters...')
-                for reader_name, reader in self.readers.items():
-                    if hasattr(reader, 'gls_parameters'):
-                        self.gls_parameters = reader.gls_parameters
-                        logger.info('Found gls-parameters in ' + reader_name)
-                        break  # Success
-                if not hasattr(self, 'gls_parameters'):
-                    logger.info('Did not find gls-parameters in any readers.')
-                    self.gls_parameters = None
-            windstress = np.sqrt(self.environment.surface_downward_x_stress**2 +
-                                 self.environment.surface_downward_y_stress**2)
-            return gls_tke(windstress, depth, self.sea_water_density(),
-                           self.environment.turbulent_kinetic_energy,
-                           self.environment.turbulent_generic_length_scale,
-                           gls_parameters)
-
         else:
             raise ValueError('Unknown diffusivity model: ' + model)
 

opendrift/models/physics_methods.py

@@ -215,45 +215,6 @@ def verticaldiffusivity_stepfunction(depth, MLD=20,
     K[depth>MLD] = k_below
     return K
 
-def gls_tke(windstress, depth, sea_water_density,
-            tke, generic_length_scale, gls_parameters=None):
-    '''From LADIM model.'''
-
-    g = 9.81
-    f0 = 0.1  # mean wave frequency
-    c_w = 4.0  # wave mixing parameter
-    c_i = 0.2  # coefficient for the interior
-    if gls_parameters is None:
-        # GLS parameters from ROMS, k-omega closure (see ocean.in)
-        p = 0.0
-        m = 1.0
-        n = 1.0
-        cmu0 = 0.5477  # for KANTHA_CLAYSON stability function
-    else:
-        p = gls_parameters['gls_p']
-        m = gls_parameters['gls_m']
-        n = gls_parameters['gls_n']
-        cmu0 = gls_parameters['gls_cmu0']
-
-    phi = 100. * (windstress/sea_water_density)**(3./2.)
-
-    # dissipation and turbulent length scale for interiour of mixed layer
-    eps = cmu0**(3.+p/n)*tke**(3./2.+m/n)*generic_length_scale**(-1./n)
-    l_i = c_i * tke**(3./2.) * eps**(-1.)
-
-    # diffusivity for interior of mixed layer
-    # c_i = sqrt(2.) * cmu0**3
-    ki = c_i * (2.*tke)**0.5 * l_i
-
-    # length scale and diffusivity of wave-enhanced layer
-    l_w = np.sqrt(phi / (g*f0))
-    kwave = c_w * (2*tke)**0.5 * l_w
-    kmix = ki + kwave
-
-    K, N = np.meshgrid(kmix, depths)
-
-    return K
-
 def plot_stokes_profile(profiles, view=['vertical', 'birdseye', 'u', 'v'], filename=None):
     '''Plot vertical profile of Stokes drift
 

opendrift/readers/reader_ROMS_native.py

@@ -135,13 +135,11 @@ def __init__(self, filename=None, name=None, gridfile=None, standard_name_mappin
             -2000, -2500, -3000, -3500, -4000, -4500, -5000, -5500, -6000,
             -6500, -7000, -7500, -8000])
 
-        gls_param = ['gls_cmu0', 'gls_p', 'gls_m', 'gls_n']
-        
         self.name = name or 'roms native'
 
         def drop_non_essential_vars_pop(ds):
             dropvars = [v for v in ds.variables if v not in
-                        list(self.ROMS_variable_mapping.keys()) + gls_param +
+                        list(self.ROMS_variable_mapping.keys()) +
                         ['ocean_time', 'time', 'bulk_time', 's_rho',
                          'Cs_r', 'Cs_rho', 'hc', 'angle', 'Vtransform']
                         and v[0:3] not in ['lon', 'lat', 'mas']]
@@ -232,16 +230,6 @@ def drop_non_essential_vars_pop(ds):
             if var not in self.Dataset.variables:
                 del self.ROMS_variable_mapping[var]
 
-        try:  # Check for GLS parameters (diffusivity)
-            self.gls_parameters = {}
-            for gls_par in gls_param:
-                self.gls_parameters[gls_par] = \
-                    self.Dataset.variables[gls_par][()]
-            logger.info('Read GLS parameters from file.')
-        except Exception as e:
-            logger.info(e)
-            logger.info('Did not find complete set of GLS parameters')
-
         # Get time coverage
         ocean_time = None
         for tv in ['ocean_time', 'time', 'bulk_time']:

tests/models/test_basemodel.py

@@ -10,7 +10,7 @@
 def test_logging(tmpdir, capsys):
     # Accepting small variations in log output,
     # depending on machine, and from which folder test is run
-    accepted = (278, 281)
+    accepted = (262, 265)
 
     # Logging to console
     logfile = None