diff --git a/constants.py b/constants.py index 45340ec7..6a45a99b 100644 --- a/constants.py +++ b/constants.py @@ -50,14 +50,23 @@ ' REMESHING OPTIONS' -remesh_method = 'log_profile' # Remeshing (log_profile or adaptive_profile) +remesh_method = 'log_profile' # possibilities: 'log_profile', 'adaptive_profile', 'uniform_profile' + +' (a) Logarithmic ' first_layer_height = 0.01 # The first layer will always have the defined height (m) layer_stretching = 1.20 # Stretching factor used by the log_profile method (e.g. 1.1 mean the subsequent layer is 10% greater than the previous +' (b) Adaptive ' merge_max = 1 # How many mergings are allowed per time step density_threshold_merging = 5 # If merging is true threshold for layer densities difference two layer try: 5-10 (kg m^-3) temperature_threshold_merging = 0.01 # If mering is true threshold for layer temperatures to merge try: 0.05-0.1 (K) +' (c) Uniform ' +layer_height_threshold = 0.1 # Height threshold for sub-surface layers (m) - layer height will stay approx. uniform except for compaction differences try 0.1 (m) + +dual_layer_height_profile = False # Enable to create a coarser sub-surface mesh below a user-defined threshold (decreases computation speed of uniform remeshing technique) +coarse_layer_depth_threshold = 20 # Simulation depth threshold at which layers are merged to form a coarser mesh. +coarse_layer_height_threshold = 0.5 # Height threshold for coarser sub-surface layers beneath the depth threshold try: 0.5 (m) ' PHYSICAL CONSTANTS ' constant_density = 300. # constant density of freshly fallen snow [kg m-3], if densification_method is set to 'constant' diff --git a/cosipy/cpkernel/grid.py b/cosipy/cpkernel/grid.py index 9ec5e375..67b002d5 100644 --- a/cosipy/cpkernel/grid.py +++ b/cosipy/cpkernel/grid.py @@ -20,13 +20,14 @@ spec['new_snow_height'] = float64 spec['new_snow_timestamp'] = float64 spec['old_snow_timestamp'] = float64 +spec['fresh_snow_timestamp'] = float64 spec['grid'] = types.ListType(node_type) @jitclass(spec) class Grid: def __init__(self, layer_heights, layer_densities, layer_temperatures, layer_liquid_water_content, layer_ice_fraction=None, - new_snow_height=None, new_snow_timestamp=None, old_snow_timestamp=None): + new_snow_height=None, new_snow_timestamp=None, old_snow_timestamp=None, fresh_snow_timestamp=None): """ The Grid-class controls the numerical mesh. The grid consists of a list of nodes (layers) that store the information @@ -63,14 +64,16 @@ def __init__(self, layer_heights, layer_densities, layer_temperatures, layer_liq # snow layer age (old_snow_timestamp) if (new_snow_height is not None) and (new_snow_timestamp is not None) and \ (old_snow_timestamp is not None): - self.new_snow_height = new_snow_height # meter snow accumulation - self.new_snow_timestamp = new_snow_timestamp # seconds since snowfall - self.old_snow_timestamp = old_snow_timestamp # snow age below fresh snow layer + self.new_snow_height = new_snow_height # meter snow accumulation + self.new_snow_timestamp = new_snow_timestamp # seconds since snowfall (reset if uppermost layer melts) + self.old_snow_timestamp = old_snow_timestamp # snow age below fresh snow layer + self.fresh_snow_timestamp = fresh_snow_timestamp # seconds since snowfall (not reset) else: #TO DO: pick better initialization values self.new_snow_height = 0.0 self.new_snow_timestamp = 0.0 self.old_snow_timestamp = 0.0 + self.fresh_snow_timestamp = 0.0 # Do the grid initialization self.grid = typed.List.empty_list(node_type) @@ -364,7 +367,33 @@ def adaptive_profile(self): # Correct first layer self.correct_layer(0 ,first_layer_height) + + def uniform_profile(self): + """ Remesh in order to try to preserve approximately uniform layer heights. + The user has the option to remesh deeper layers into coarser layer heights in order to speed up computation time. + + The layer_height_threshold, coarse_layer_height_threshold and coarse_layer_depth_threshold variables in the + configuration file (constants.py) define the corresponding thresholds. + """ + # Merge fresh snow with the uppermost layer unless it exceeds the maximum layer height + if ((self.get_node_height(0) + self.get_node_height(1) <= layer_height_threshold) & (self.fresh_snow_timestamp == 0)): + self.merge_nodes(0) + + # Merge internal snow layers if they subsceed the minimum snow layer height + idx = 1 + while ((idx < self.get_number_snow_layers()-1)): + if ((self.get_node_height(idx) <= minimum_snow_layer_height)): + self.merge_nodes(idx-1) + idx += 1 + + if dual_layer_height_profile == True: + # Merge the lower layers if they go beyond the user-defined depth threshold to make a coarser mesh: + idx = 0 + while ((idx < self.get_number_snow_layers()-1)): + if ((self.get_node_depth(idx) > coarse_layer_depth_threshold) & ((self.get_node_height(idx) + self.get_node_height(idx + 1) <= coarse_layer_height_threshold))): + self.merge_nodes(idx) + idx += 1 def split_node(self, pos): """ Split node at position pos @@ -438,21 +467,28 @@ def check(self, name): def update_grid(self): """ Re-meshes the layers (numerical grid). - Two algorithms are currently implemented to re-mesh the layers: + Three algorithms are currently implemented to re-mesh the layers: (i) log_profile (ii) adaptive_profile - - (i) The log-profile algorithm arranges the mesh logarithmically. - The user provides a stretching factor (layer_stretching in the configuration file) - that determines the increase in layer heights. - - (ii) The adjustment of the profile is done on the basis of the similarity of layers. - Layers with very similar states (temperature and density) are joined together. The - similarity is determined by user-specified threshold values - (temperature_threshold_merging, density_threshold_merging). In - addition, the maximum number of merging steps per time step - can be specified (merge_max). + (iii) uniform_profile + + (i) The log-profile algorithm arranges the mesh logarithmically. + The user provides a stretching factor (layer_stretching in the configuration file) + that determines the increase in layer heights. + + (ii) The adjustment of the profile is done on the basis of the similarity of layers. + Layers with very similar states (temperature and density) are joined together. The + similarity is determined by user-specified threshold values + (temperature_threshold_merging, density_threshold_merging). In + addition, the maximum number of merging steps per time step + can be specified (merge_max). + + (iii) The uniform algorithm aims to maintain layers of approximately uniform height. + Fresh fresh snowfall events are merged into a common layer until a user-specified height + threshold is reached and a new layer is created. Note that layers are not strictly forced + to be exactly uniform - compaction and differences in snowfall mean they will always be + slighly different in height. """ #------------------------------------------------------------------------- @@ -462,12 +498,13 @@ def update_grid(self): self.log_profile() elif (remesh_method=='adaptive_profile'): self.adaptive_profile() + elif (remesh_method=='uniform_profile'): + self.uniform_profile() # if first layer becomes very small, remove it if (self.get_node_height(0)