Source code for dycove.sim.engines.DFM_hydro

###############################################################
#  DFM_hydro.py
###############################################################

import re
import os
import sys
from types import ModuleType
from datetime import datetime
from pathlib import Path
import numpy as np

from dycove.sim.base import HydroSimulationBase, HydroEngineBase
from dycove.utils.simulation_reporting import Reporter
from dycove.constants import H_LIM_VELOCITY


r = Reporter()

def _import_bmi():
    """ Lazy loading of bmi wrapper to avoid import errors when DFM will not be tested/used. """
    try:
         # need to mock this due to incompatibility with updated setuptools versions (python 3.7+)
        sys.modules['pkg_resources'] = ModuleType('pkg_resources')
        from bmi.wrapper import BMIWrapper
        return BMIWrapper
    except ImportError:
        msg = ("The `bmi` package is not installed. "
               "Refer to the documentation for installation instructions.")
        r.report(msg, level="ERROR")
        raise ImportError(msg)
    

[docs] class DFM(HydroSimulationBase): """ Hydrodynamic simulation wrapper for the Delft3D FM model. This class connects the generic :class:`~dycove.sim.base.HydroSimulationBase` to :class:`~dycove.sim.base.engines.DFM_hydro.DFMEngine`, providing a consistent Python interface for running D-Flow FM through its BMI and DIMR interfaces. Notes ----- - All higher-level logic that can be abstracted from the engine classes is handled in :class:`~dycove.sim.base.HydroSimulationBase`; all low-level model interactions are delegated to :class:`~dycove.sim.base.engines.DFM_hydro.DFMEngine`. """
[docs] def __init__(self, dfm_path, config_path, mdu_path, vegetation=None): # build DFM engine engine = DFMEngine(dfm_path, config_path, mdu_path, vegetation) # pass DFM engine to the base class super().__init__(engine)
[docs] class DFMEngine(HydroEngineBase): """ Engine interface for DFM hydro-morphodynamic model. This engine: - Loads and initializes DFM executables (DIMR + D-Flow FM BMI). - Manages exchange of hydrodynamic and vegetation state variables though DFM-specific ``BMI-python`` wrapper. - Ensures that required input files are present and are consistent with simulation settings. Parameters ---------- dfm_path : Path or str Path to the root Delft3D-FM installation directory. Might look like this: 'C:/Program Files (x86)/Deltares/Delft3D Flexible Mesh Suite HM (2021.03)/plugins/DeltaShell.Dimr/kernels/x64' config_path : Path or str Path to DIMR configuration file ``dimr_config.xml`` mdu_path : Path or str Path to DFM MDU file. vegetation : VegetationSpecies or MultipleVegetationSpecies, optional Vegetation object passed down from the base simulation. Notes ----- - Vegetation files (.xyz) required by DFM vegetation module are auto-created if missing. - Parallel mode is not currently implemented. """
[docs] def __init__(self, dfm_path, config_path, mdu_path, vegetation=None): self.dll_dirs = self.add_dll_directories(dfm_path) # do this first to setup PATH before loading self.BMIWrapper = _import_bmi() # lazy bmi loading # Define paths to Delft3D FM .dll files self.dflowfm_path = Path(dfm_path) / ("dflowfm/bin/dflowfm.dll") self.dimr_path = Path(dfm_path) / ("dimr/bin/dimr_dll.dll") self.mdu_path = mdu_path # location of MDU file that contains model directions/inputs self.model_dir = mdu_path.parent # model directory containing MDU and other model files self.config_path = config_path # location of config file used for running DFM using dimr self.veg = vegetation self.open_bmi_wrappers()
[docs] def add_dll_directories(self, dfm_path): """ Add DLL paths to env before calling BMI. Return list to retain handle and avoid accidental garbage collecting. Note that for versions of python 3.7 and earlier, you would need to set the env variables differently: .. code-block:: python os.environ['PATH'] = os.path.join(dfm_path, 'share', 'bin') + ";" + os.path.join(dfm_path, 'dflowfm', 'bin') + ";" + ... ) """ return [os.add_dll_directory(dfm_path / Path("dflowfm/bin")), os.add_dll_directory(dfm_path / Path("dimr/bin")), os.add_dll_directory(dfm_path / Path("share/bin")), os.add_dll_directory(dfm_path / Path("dwaves/bin")), os.add_dll_directory(dfm_path / Path("esmf/scripts")), os.add_dll_directory(dfm_path / Path("swan/scripts")), ]
[docs] def open_bmi_wrappers(self): """ Create BMI wrapper objects for DFM and DIMR """ # BMI wrapper object that interacts with the actual numerical model (e.g., getting and setting variables) self.dflowfm = self.BMIWrapper(engine=str(self.dflowfm_path), configfile=str(self.mdu_path)) # BMI wrapper object that handles the deployment of the model executables self.dimr = self.BMIWrapper(engine=str(self.dimr_path), configfile=str(self.config_path))
[docs] def initialize(self): ### ----- Required for vegetation module ----- ### self.mdu_vars = self.get_model_inputs() self.morphology, self.morph_vars = self.get_morphodynamic_inputs() self.vegetation_file_check() ### ----- Required for numerical model ----- ### self.dimr.initialize()
[docs] def step(self, seconds): self.dimr.update(seconds)
[docs] def cleanup(self): self.dimr.finalize()
[docs] def get_cell_count(self): return int(self.dflowfm.get_var("ndxi")) # number of non-boundary boxes, i.e. within-domain boxes
[docs] def get_refdate(self): # This input has a line in the MDU file, but most other models probably don't care what the date is and we can just hardcode the date refdatestr = self.mdu_vars["RefDate"] return datetime(int(refdatestr[:4]), int(refdatestr[4:6]), int(refdatestr[6:]))
[docs] def get_elevation(self): # DFM returns arrays with boundary values included, slice those out first n_cells = self.get_cell_count() return np.array(self.dflowfm.get_var("bl"))[:n_cells]
[docs] def get_velocity_and_depth(self): n_cells = self.get_cell_count() depth = np.array(self.dflowfm.get_var("hs"))[:n_cells] velocity = np.array(self.dflowfm.get_var("ucmag"))[:n_cells] # Ignore velocities where depth is insufficient velocity = np.where(depth < H_LIM_VELOCITY, 0., velocity) return velocity, depth
[docs] def get_vegetation(self): # Convert to numpy arrays because DFM returns pointers and we don't want to accidentally modify them stemdensity = np.array(self.dflowfm.get_var("rnveg")) stemdiameter = np.array(self.dflowfm.get_var("diaveg")) stemheight = np.array(self.dflowfm.get_var("stemheight")) return stemdensity, stemdiameter, stemheight
[docs] def set_vegetation(self, stemdensity, stemdiameter, stemheight): self.dflowfm.set_var("rnveg", stemdensity) self.dflowfm.set_var("diaveg", stemdiameter) self.dflowfm.set_var("stemheight", stemheight)
[docs] def check_simulation_inputs(self, simstate): """ Compare DYCOVE simulation time to MDU simulation time. All MDU files (DFM models) will have a simulation time specifed, but DYCOVE will run DFM for a period of time based on how many veg years we want to simulate. Basically, the time specified in the MDU needs to be arbitrarily large enough so that we never run into the issue of the model stopping prematurely. """ if simstate.hydro_sim_days*86400 > int(self.mdu_vars["TStop"]): msg = (f"Model simulation time specified in MDU file (TStop = {self.mdu_vars['TStop']}) not long enough based on " "the inputs provided for sim_years, n_ets, and veg_interval, which give simulation length of " f"{simstate.hydro_sim_days*86400}. Please provide an arbitrarily larger number for TStop in MDU file.") r.report(msg, level="ERROR") raise ValueError(msg)
# -------------------------------------------------------- # Some additional required, DFM-specific methods # --------------------------------------------------------
[docs] def get_model_inputs(self): """ Read lines from MDU file into a dictionary """ mdu_lines = self.mdu_path.read_text().splitlines() mdu_vars = {} for line in mdu_lines: if "=" in line: slist = re.split("=|#", line) mdu_vars[slist[0].strip()] = slist[1].strip() return mdu_vars
[docs] def get_morphodynamic_inputs(self): """ Read lines from morphology file into a dictionary, if mprph files exist """ # Same filename as MDU, different extension mor_filepath = self.model_dir / (self.mdu_path.stem + ".mor") morph_vars = {} morphology = True if mor_filepath.exists(): with open(mor_filepath) as f: for line in f: if "=" in line: slist = re.split("=|#", line) morph_vars[slist[0].strip()] = slist[1].strip() else: msg = "Morphology file NOT FOUND; proceeding with morphology off." r.report(msg, level="WARNING") morphology = False # # If morphology is off, ensure that vegetation mor variable is set to 0 (no burial/scour) # # I don't think this is necessary, burial/scour will be zero (nelson) # if self.veg is not None: # self.veg.mor = 0 return morphology, morph_vars
[docs] def vegetation_file_check(self): """ MODIFIES .mdu file if certain vegetation-related lines are not present: - Adds a filename next to 'ExtForceFile' if blank, creates the file too - Adds drag coefficient from VegetationAttributes if [veg] block is present - Adds appropriate Baptist model number (1 if no morph, 2 if morph) - Adds [veg] block if it is not present, including drag coefficient from VegetationAttributes and appropriate Baptist model number Creates empty text files for stem density, stem diameter, and stem height, if they don't already exist. The filenames are those specified in the vegetation .ext file in the model directory (e.g., "FlowFM_veg.ext"). These files can be created beforehand if prior vegetation establishment is desired. Otherwise, blank files are required so that DFM knows to store these variables through time. """ # Read .mdu file lines self.mdu_lines = self.mdu_path.read_text().splitlines() # Track for .mdu file modification self.mdu_modified = False # All only execute if self.veg is not None self.add_extforcefile_to_mdu() self.add_veg_module_to_mdu() self.create_extforcefile() self.create_veg_xyz_files() if self.mdu_modified: self.write_modified_mdu()
[docs] def add_extforcefile_to_mdu(self): """ Add ExtForceFile to .mdu line if it's not there (and if vegetation is active). [external forcing] ExtForceFile = FlowFM.ext # Old format for external forcings file ... """ if self.mdu_vars["ExtForceFile"] == "" and self.veg is not None: self.mdu_modified = True # Get name of model/file based on name of "new" .ext file try: replacement = self.mdu_vars["ExtForceFileNew"].replace("_bnd", "") except: msg = ("Either the 'ExtForceFileNew' file name in the .mdu file does not end in the expected " "'_bnd.ext', or there is no 'ExtForceFileNew' file defined in the .mdu file. If it was " "purposeful that no boundaries were specified for this model, then this check mechanism " "must be updated: please get in touch with us on GitHub.") r.report(msg, level="ERROR") raise NameError(msg) self.mdu_vars["ExtForceFile"] = replacement for i, line in enumerate(self.mdu_lines): # Replace blank space with name of required .ext file if line.startswith("ExtForceFile "): slist = re.split("=|#", line) n_spaces = len(slist[1]) self.mdu_lines[i] = f"{slist[0]}= {replacement}{' '*max(n_spaces - len(replacement) - 1, 1)}#{slist[2]}" # Replace drag coefficient value with the one provided in input .json file (if [veg] block is present) if line.strip().startswith("Cdveg"): slist = re.split("=|#", line) drag = self.veg.get_drag() self.mdu_lines[i] = f"{slist[0]}= {drag:.1f}{' '*13}#{slist[2]}"
[docs] def add_veg_module_to_mdu(self): """ Add [veg] section to .mdu if it's not there (and if vegetation is active). Format: [veg] Vegetationmodelnr = 2 # 1: Baptist, 2: Baptist with morphology correction factor (lambda) Clveg = 0.8 # Stem distance factor, default=0.8 Cdveg = 1.1 # Drag coefficient, pulled from input veg.json file Cbveg = 0.7 # Stem stiffness coefficient, default=0.7 """ veg_block_present = any(line.strip().startswith("[veg]") for line in self.mdu_lines) if not veg_block_present and self.veg is not None: drag = self.veg.get_drag() veg_model_num = 2 if self.veg.mor == 1 else 1 self.mdu_modified = True self.mdu_lines.append("") self.mdu_lines.extend([ "[veg]", f"Vegetationmodelnr = {veg_model_num} # 1: Baptist et al. (2007) equation for calculation of vegetation roughness", "Clveg = 0.8 # Stem distance factor, default=0.8", f"Cdveg = {drag:.1f} # Stem Cd coefficient, default=0.7", "Cbveg = 0.7 # Stem stiffness coefficient, default=0.7", ])
[docs] def create_extforcefile(self): """ Create .ext file in the model directory if it doesn't exist """ ext_force_file = self.model_dir / self.mdu_vars["ExtForceFile"] content = """QUANTITY=stemdensity FILENAME=stemdensity.xyz FILETYPE=7 METHOD=5 OPERAND=O QUANTITY=stemdiameter FILENAME=stemdiameter.xyz FILETYPE=7 METHOD=5 OPERAND=O QUANTITY=stemheight FILENAME=stemheight.xyz FILETYPE=7 METHOD=5 OPERAND=O """ if not ext_force_file.exists() and self.veg is not None: with open(ext_force_file, "w") as f: f.write(content) # It may already exist if other spatially varying parameters are in use; need to append our content elif ext_force_file.exists() and self.veg is not None: with open(ext_force_file, "r") as f: lines = f.read() with open(ext_force_file, "w") as f: f.write(lines) f.write("\n") f.write(content)
[docs] def create_veg_xyz_files(self): """ Create required files for [veg] module to run, even if they are blank """ req_veg_files = ["stemdensity.xyz", "stemdiameter.xyz", "stemheight.xyz"] for filename in req_veg_files: veg_file = self.model_dir / filename if not veg_file.exists() and self.veg is not None: with open(veg_file, "w") as f: f.write("")
[docs] def write_modified_mdu(self): """ Write modified .mdu lines back to file """ self.mdu_path.write_text("\n".join(self.mdu_lines) + "\n") msg = "DFM MDU file updated and rewritten to include required inputs for vegetation module." r.report(msg)
# -------------------------------------------------------- # Parallel methods # --------------------------------------------------------
[docs] def get_rank(self): # TODO: implement parallel processing for DFM return 0
[docs] def is_parallel(self): try: from mpi4py import MPI comm = MPI.COMM_WORLD size = comm.Get_size() return True if size > 1 else False except: return False
[docs] def merge_parallel_veg(self, OutputManager): # DYCOVE-DFM does not currently support parallel processing in this model, but setting up for future use. raise NotImplementedError("Parallel mode not currently implemented for Delft3D FM")