diff --git a/src/Solver/RubarBE.py b/src/Solver/RubarBE.py index 11f571cc..139dc64a 100644 --- a/src/Solver/RubarBE.py +++ b/src/Solver/RubarBE.py @@ -173,6 +173,8 @@ class Rubar3(CommandLineSolver): study.river, hydraulics["edge"], values, trased, hydraulics["warnings"]) saves = self._save_data(values, hydraulics) + laterals = self._lateral_data( + study.river, hydraulics, values, hydraulics["warnings"]) except RubarExportError as error: self._cancel(qlog, str(error)) return False @@ -190,6 +192,8 @@ class Rubar3(CommandLineSolver): self._export_geomac_i(geometry, repertory, qlog, name=name) self._export_mail(study, repertory, qlog, name=name) self._export_devers(study, repertory, qlog, name=name) + if laterals is not None: + self._export_hydevs(laterals, trased, repertory, qlog, name=name) self._export_condin(hydraulics, repertory, qlog, name=name) self._export_stricklers(hydraulics, repertory, qlog, name=name) self._export_boundaries(hydraulics, repertory, qlog, name=name) @@ -529,6 +533,103 @@ class Rubar3(CommandLineSolver): name=bc.name, node=edge.node1.name)) return bc + ############ + # Laterals # + ############ + + # Pamhyr2 lateral contribution types; evaporation is a withdrawal + _lateral_signs = {"LC": 1.0, "RA": 1.0, "EV": -1.0} + # QMU, TMU (variables.f90, LMAX*20): RubarBE copies the series of an + # inflow in each cell of its zone, without checking this size + _max_lateral_points = 3000 * 20 + + def _lateral_data(self, river, hydraulics, values, warnings): + """Lateral inflows of HYDEVS, read when IODEV=O (rubar_init.f90, + init2). Return [(first cell, last cell, [(t, q), ...]), ...] + sorted from upstream, or None when IODEV=N. + + A contribution from section b to section e (1-based) covers the + cells b+1 to e, i.e. from section b to section e (RubarBE divides + the discharge by this length); a contribution at a single section + b covers the cell downstream of it (upstream of it for the last + section). The first and last cells of the reach are never used, + as RubarBE requires. With IODEV=N, the contributions are ignored + with a warning. Otherwise, solid and suspension contributions have + no RubarBE equivalent (their sediment cannot be linked to an + inflow): they are refused""" + trad = RunSolverTranslate() + edge = hydraulics["edge"] + profiles = hydraulics["profiles"] + iodev = values.get("rubarbe_iodev", "n").lower() in ("y", "o") + lst = river.lateral_contribution + + def on_reach(tab): + # A contribution without value carries no data: ignored + kept = [] + for lc in lst.get_tab(tab): + if lc.is_deleted() or lc.reach != edge: + continue + if len(lc) == 0: + warnings.append(trad["lateral_empty_ignored"].format( + name=lc.name, reach=edge.name)) + else: + kept.append(lc) + return kept + + lcs = on_reach("liquid") + sediment = on_reach("solid") + on_reach("suspenssion") + + # Not read by RubarBE: the user is warned, the computation runs + # without them + if not iodev: + for lc in lcs + sediment: + warnings.append(trad["lateral_iodev_disabled"].format( + name=lc.name, reach=edge.name)) + return None + + if sediment: + raise RubarExportError(trad["lateral_solid_unsupported"].format( + name=sediment[0].name, reach=edge.name)) + + index = {id(p): i for i, p in enumerate(profiles, 1)} + count = len(profiles) + blocks = [] + for lc in lcs: + if lc.lctype not in self._lateral_signs: + raise RubarExportError(trad["lateral_type_unsupported"].format( + name=lc.name, reach=edge.name)) + if (id(lc.begin_section) not in index or + id(lc.end_section) not in index): + raise RubarExportError(trad["lateral_sections"].format( + name=lc.name, reach=edge.name)) + b, e = sorted((index[id(lc.begin_section)], + index[id(lc.end_section)])) + first, last = (b + 1, e) if b < e else ((min(b + 1, count),) * 2) + + data = [(float(d[0]), float(d[1])) for d in lc.data] + times = [t for t, _ in data] + if any(t1 < t0 for t0, t1 in zip(times, times[1:])): + raise RubarExportError(trad["lateral_not_sorted"].format( + name=lc.name, reach=edge.name)) + sign = self._lateral_signs[lc.lctype] + blocks.append((first, last, [(t, sign * q) for t, q in data], + lc.name)) + + blocks.sort(key=lambda block: block[0]) + for (_, last, _, name0), (first, _, _, name1) in zip( + blocks, blocks[1:]): + if first <= last: + raise RubarExportError(trad["lateral_overlap"].format( + first=name0, second=name1, reach=edge.name)) + + points = sum((last - first + 1) * len(data) + for first, last, data, _ in blocks) + if points > self._max_lateral_points: + raise RubarExportError(trad["lateral_too_long"].format( + count=points, maximum=self._max_lateral_points)) + + return [block[:3] for block in blocks] + ######### # Saves # ######### @@ -1005,6 +1106,9 @@ class Rubar3(CommandLineSolver): self._write_parameters(repertory, f"tnprof.{name}", lines) def _export_devers(self, study, repertory, qlog, name="0"): + """DEVERS, required by RubarBE when IODEV=O: per cell, weir + level, coefficient and exponent (notice annex 1). Pamhyr2 has no + lateral weir data: a zero coefficient everywhere (no overflow)""" if qlog is not None: qlog.put("Export DEVERS file") @@ -1026,6 +1130,25 @@ class Rubar3(CommandLineSolver): reach_ind += 1 + def _export_hydevs(self, laterals, trased, repertory, qlog, name="0"): + """HYDEVS: per inflow, first and last cells and number of points, + then (Q, t) per point (rubar_init.f90, init2). With sediment + transport, D and S on the first line and (Q, t, Qs, ts) per point: + Pamhyr2 gives no sediment for a liquid inflow, so Qs = 0 (and + D = S = 0, the sediment of the main channel). Written, possibly + empty, so that a previous HYDEVS of the directory is not used""" + if qlog is not None: + qlog.put("Export HYDEVS file") + lines = [] + for first, last, data in laterals: + if trased: + lines.append(f"{first} {last} {len(data)} 0.0 0.0\n") + lines += [f"{q!r} {t!r} 0.0 {t!r}\n" for t, q in data] + else: + lines.append(f"{first} {last} {len(data)}\n") + lines += [f"{q!r} {t!r}\n" for t, q in data] + self._write_parameters(repertory, f"hydevs.{name}", lines) + def _export_stricklers(self, hydraulics, repertory, qlog, name="0"): """FROT (minor bed) and FROT2 (medium bed): LM, then one Strickler per cell (notice, annex 1)""" diff --git a/src/View/RunSolver/Translate.py b/src/View/RunSolver/Translate.py index e0f2902e..600addbd 100644 --- a/src/View/RunSolver/Translate.py +++ b/src/View/RunSolver/Translate.py @@ -129,6 +129,51 @@ class RunSolverTranslate(PamhyrTranslate): "Export cancelled: several {kind} boundary conditions at node " "{node}." ) + self._dict["lateral_empty_ignored"] = _translate( + "Solver", + "Lateral contribution {name} on reach {reach} has no value: it " + "is ignored." + ) + self._dict["lateral_solid_unsupported"] = _translate( + "Solver", + "Export cancelled: solid or suspension lateral contribution " + "{name} on reach {reach}: RubarBE only reads the sediment of a " + "liquid inflow, which Pamhyr2 cannot link to it. Remove it." + ) + self._dict["lateral_iodev_disabled"] = _translate( + "Solver", + "Lateral contribution {name} on reach {reach} is not taken into " + "account: lateral overflows (rubarbe_iodev) are disabled in the " + "solver parameters." + ) + self._dict["lateral_type_unsupported"] = _translate( + "Solver", + "Export cancelled: lateral contribution {name} on reach {reach} " + "has no type. Choose a type (lateral sources, rain or " + "evaporation)." + ) + self._dict["lateral_sections"] = _translate( + "Solver", + "Export cancelled: the sections of lateral contribution {name} " + "are not profiles of reach {reach}." + ) + self._dict["lateral_not_sorted"] = _translate( + "Solver", + "Export cancelled: the times of lateral contribution {name} on " + "reach {reach} must be increasing." + ) + self._dict["lateral_overlap"] = _translate( + "Solver", + "Export cancelled: lateral contributions {first} and {second} " + "on reach {reach} overlap, which RubarBE does not accept. Merge " + "them or move one of them." + ) + self._dict["lateral_too_long"] = _translate( + "Solver", + "Export cancelled: the lateral contributions give {count} values " + "(number of points times number of cells of their zone), " + "RubarBE accepts at most {maximum}." + ) self._dict["boundary_empty_ignored"] = _translate( "Solver", "Boundary condition {name} at node {node} has no value: it is "