diff --git a/src/Solver/RubarBE.py b/src/Solver/RubarBE.py index 2f1ce90d..e34b9fe3 100644 --- a/src/Solver/RubarBE.py +++ b/src/Solver/RubarBE.py @@ -189,6 +189,9 @@ class Rubar3(CommandLineSolver): self._export_boundaries(hydraulics, repertory, qlog, name=name) if solid is not None: self._export_hydros(solid, repertory, qlog, name=name) + if trased: + self._export_condins(geometry, hydraulics["tinit"], repertory, + qlog, name=name) # self._export_abshyd(study, repertory, qlog, name=name) return True @@ -302,17 +305,11 @@ class Rubar3(CommandLineSolver): name=names[marker], profile=self._profile_label(profile), reach=edge.name)) - # Defaults of RubarBE: first and last points for G and D, G and D - # for g and d. Erosion only applies strictly between g and d - def index(marker, default): - return markers.index(marker) if marker in markers else default - - left = index("G", 0) - right = index("D", len(points) - 1) + left, right, active_left, active_right = self._bed_limits(markers) limits = [("main_channel", "rg", "rd", left, right)] if trased: limits.append(("active_bed", "g", "d", - index("g", left), index("d", right))) + active_left, active_right)) for bed, lname, rname, i, j in limits: if j - i < 2: raise RubarExportError(trad["bed_limits_invalid"].format( @@ -320,6 +317,18 @@ class Rubar3(CommandLineSolver): profile=self._profile_label(profile), reach=edge.name)) return markers + @staticmethod + def _bed_limits(markers): + """Indices of the G, D, g and d points of a section. Defaults of + RubarBE: first and last points for G and D, G and D for g and d. + Erosion only applies strictly between g and d""" + def index(marker, default): + return markers.index(marker) if marker in markers else default + + left = index("G", 0) + right = index("D", len(markers) - 1) + return left, right, index("g", left), index("d", right) + def _compartments(self, edge, profile, number, point): """Sediment layers of a point as RubarBE compartments, from the top: (top elevation, D, S, critical stress). A non erodible @@ -965,6 +974,27 @@ class Rubar3(CommandLineSolver): ] self._write_parameters(repertory, f"hydros.{name}", lines) + def _export_condins(self, geometry, tinit, repertory, qlog, name="0"): + """CONDINS: TINIT, then for each section (intermaille) its number, + Qs, D and S of the active layer (notice annex 1, rubar_init.f90, + inits). Pamhyr2 holds no initial solid condition: Qs = 0 (empty + active layer) and D, S are the mean of the top compartments of the + points strictly between g and d, as advised by the notice ยง3.5""" + if qlog is not None: + qlog.put("Export CONDINS file") + lines = [f"{tinit!r}\n"] + for i, (_, rows) in enumerate(geometry, 1): + markers = [marker for marker, _, _ in rows] + _, _, left, right = self._bed_limits(markers) + tops = [compartments[0] for _, _, compartments + in rows[left + 1:right]] + # Rounded to 10 significant digits, so that a mean of equal + # values is written as these values + d = float(f"{sum(c[1] for c in tops) / len(tops):.10g}") + s = float(f"{sum(c[2] for c in tops) / len(tops):.10g}") + lines.append(f"{i} 0.0 {d!r} {s!r}\n") + self._write_parameters(repertory, f"condins.{name}", lines) + def read_profil(self, study, fname, results, qlog=None, name="0"): logger.info(f"read: profil.{name}") with open(fname, "r") as f: