diff --git a/src/psimpy/simulator/__init__.py b/src/psimpy/simulator/__init__.py index 75652d1cbde76dd875545a431733b67abf9dbc42..78b86034ef6baba0aa10c877f89de4f7e2a44904 100644 --- a/src/psimpy/simulator/__init__.py +++ b/src/psimpy/simulator/__init__.py @@ -1,3 +1,4 @@ from psimpy.simulator.mass_point_model import MassPointModel from psimpy.simulator.ravaflow24 import Ravaflow24Mixture -from psimpy.simulator.run_simulator import RunSimulator \ No newline at end of file +from psimpy.simulator.ravaflow3G import Ravaflow3GMixture +from psimpy.simulator.run_simulator import RunSimulator diff --git a/src/psimpy/simulator/ravaflow3G.py b/src/psimpy/simulator/ravaflow3G.py index 74b75820aaddc95ec445fc049351ce94d62dcdf0..0c0910b8980dcd5dcb12851e43dccdd7a80cfbb0 100644 --- a/src/psimpy/simulator/ravaflow3G.py +++ b/src/psimpy/simulator/ravaflow3G.py @@ -19,12 +19,11 @@ class Ravaflow3GMixture: curvature_control: str = "1", surface_control: str = "0", entrainment_control: str = "0", - shear_velocity_coef: float = 0.05, - basal_friction_diff: float = 0.0, stopping_control: str = "0", stopping_threshold: float = 0.0, friction_control: str = "0", non_hydro_control: str = "0", + phase_sep_control: str = "0", hydrograph_control: str = "2", deceleration_control: str = "0", ) -> None: @@ -56,29 +55,22 @@ class Ravaflow3GMixture: `'1'`: apply balancing of forces. entrainment_control : str `'0'`: no entrainment (default). - `'1'`: the entrainment coefficient is multiplied with flow momentum. + `'1'`: the entrainment coefficient (see `preprocess`) is multiplied with flow momentum. `'2'`: simplified entrainment and deposition model is used. `'3'`: combination of `'1'` for entrainmnet and `'2'` for deposition. `'4'`: acceleration-deceleration entrainment and deposition model. - shear_velocity_coef : float - Only used with ``entrainment_control`` being `'2'` or `'3'`, range - :math:`[0,1]`. - basal_friction_diff : float - Difference between the basal friction angles of the basal surface - and the flow, in degrees. Only used with ``entrainment_control`` - being `'2'` or `'3'`. stopping_control : str `'0'`: no stopping control is used (default). - `'1'`: stop if the flow kinetic energy is equal or smaller than a - given threshold. - `'2'`: stop if the flow momentum is equal or smaller than a given - threshold. + `'1'`: stop if the flow kinetic energy ratio is equal or smaller than + given threshold ``stopping_threshold``. + `'2'`: stop if the flow momentum ratio is equal or smaller than + given threshold ``stopping_threshold``. `'3'`: stop if the dynamic flow pressure of all raster cells is euqal or smaller than a given threshold. stopping_threshold : float - Threshold value for ``stopping_control``. - If ``stopping_control`` is `'1'` or `'2'`, ``stopping_threshold`` - has to be given as the maximum value reached during the flow. + Threshold value for `stopping_control`. + If `stopping_control` is `'1'` or `'2'`, a ``stopping_threshold`` + has to be given with respect to maximal flow kinetic energy/momentum. If ``stopping_control`` is `'3'`, the pressure threshold has to be specified. friction_control : str, optional @@ -90,6 +82,9 @@ class Ravaflow3GMixture: `'1'`: Only dispersion is considered. `'2'`: Only enhanced gravity is considered. `'3'`: Both dispersion and enhanced gravity are considered. + phase_sep_control: str + `'0'`: Phase separation is neglected (default). + `'1'`: Phase separation is considered. hydrograph_control: str `'0'`: Impose hydrograph on the flow. `'1'`: Reset flow at the hydrograph profile before imposing the hydrograph. @@ -124,9 +119,6 @@ class Ravaflow3GMixture: raise ValueError("entrainment_control must be '0', '1', '2', '3', or '4'") self.entrainment_control = entrainment_control - self.shear_velocity_coef = shear_velocity_coef - self.basal_friction_diff = basal_friction_diff - if stopping_control not in ["0", "1", "2", "3"]: raise ValueError("stopping_control must be '0', '1', '2', or '3'") self.stopping_control = stopping_control @@ -141,6 +133,10 @@ class Ravaflow3GMixture: raise ValueError("non_hydro_control must be '0', '1', '2', or '3'") self.non_hydro_control = non_hydro_control + if phase_sep_control not in ["0", "1"]: + raise ValueError("phase_sep_control must be '0' or '1'.") + self.phase_sep_control = phase_sep_control + if hydrograph_control not in ["0", "1", "2"]: raise ValueError("hydrograph_control must be '0', '1', or '2'") self.hydrograph_control = hydrograph_control @@ -164,7 +160,7 @@ class Ravaflow3GMixture: EPSG: Union[str, None] = None, ) -> tuple[str, str]: """ - Preprocess simulation input data simulation. + Preprocess simulation input data. Parameters ---------- @@ -195,9 +191,9 @@ class Ravaflow3GMixture: except for :math:`0` meaning no entrainment. EPSG : str, optional `EPSG` (European Petroleum Survey Group) code to create - `GRASS Location <https://grass.osgeo.org/grass82/manuals/grass_database.html>`_. + GRASS Location <https://grass.osgeo.org/grass82/manuals/grass_database.html>. If `None`, ``elevation`` must be a georeferenced file which has - metadata to create the `GRASS Location`. + metadata to create the GRASS Location. Returns ------- @@ -205,7 +201,7 @@ class Ravaflow3GMixture: Name of the `GRASS Location` (including path). sh_file : str Name of the shell file (including path), which will be called by - `GRASS <https://grass.osgeo.org/>`_ to run the simulation. + `GRASS <https://grass.osgeo.org/>` to run the simulation. """ sh_file = os.path.join(self.dir_sim, f"{prefix}_shell.sh") @@ -230,18 +226,19 @@ class Ravaflow3GMixture: # create grass location if EPSG is None: - os.system(f"grass -e -c {elevation} {grass_location}") + os.system(f"grass -c XY -e {grass_location}") else: os.system(f"grass -e -c EPSG:{EPSG} {grass_location}") # create shell file with open(sh_file, "w") as sh: sh.write("# import elevation raster \n") + sh.write("g.region -d \n") sh.write(f"r.in.gdal -o input={elevation} output=elev --overwrite") sh.write("\n\n") sh.write("# set region based on elev \n") - sh.write("g.region raster=elev \n\n") + sh.write("g.region -s raster=elev \n\n") sh.write("# import hrelease raster \n") sh.write(f"r.in.gdal -o input={hrelease} output=raw_hrel --overwrite") @@ -255,12 +252,11 @@ class Ravaflow3GMixture: f"elevation=elev hrelease=hrel " f"friction={internal_friction},{basal_friction}," f"{turbulent_friction} " - f"basal={entrainment_coef},{self.shear_velocity_coef}," - f"{self.basal_friction_diff},{self.stopping_threshold} " + f"basal={entrainment_coef},{self.stopping_threshold} " f"controls=0,{self.diffusion_control},{self.curvature_control}," f"{self.surface_control},{self.entrainment_control}," f"{self.stopping_control},{self.friction_control},{self.non_hydro_control}," - f"{self.hydrograph_control},{self.deceleration_control} " + f"{self.phase_sep_control},{self.hydrograph_control},{self.deceleration_control} " f"time={self.time_step},{self.time_end} " f"cfl={self.cfl},{self.time_step_length}" ) @@ -311,7 +307,7 @@ class Ravaflow3GMixture: `'t'`: maximum flow kinetic energy, in :math:`J`. threshold : float or int Threshold of ``qoi`` to determine impact area. Areas where ``qoi`` - is larger than ``threshold`` are regarded as impar area. + is larger than ``threshold`` are regarded as impact area. Returns ------- @@ -356,9 +352,9 @@ class Ravaflow3GMixture: `'p'`: maximum flow pressure, in :math:`Pa`. `'t'`: maximum flow kinetic energy, in :math:`J`. aggregate : bool - If `True`, returns the overall maximum value over all spatio-tempo + If `True`, returns the overall maximum value over all spatio-temporal grids. If `False`, returns the maximum values over all time steps at - each spatio location, namely a raster of maximum values. + each spatial location, namely a raster of maximum values. Returns -------