From 38152bf4443fae3f161f77003fbb3c24136f08d0 Mon Sep 17 00:00:00 2001 From: sarangbhagwat Date: Sun, 23 Aug 2026 00:00:12 -0700 Subject: [PATCH 1/3] add synthetic HXN regression tests guarding energy balance, MER targets, and documented heat recovery Ten synthetic systems of increasing complexity (2 to 10 process streams, built from simulated HXutility units; every case from the third on has phase-changing streams: total and partial condensers, boiling and partially boiling cold streams, mixed pressures, a threshold problem and a near-degenerate cascade). For each, the synthesized network must (i) close its energy balance to |error| < 0.1 % with no numerical RuntimeWarning escaping simulate(); (ii) never beat the problem-table MER targets computed on the same streams (a network below MER is a bookkeeping error), and satisfy actual_heat - actual_cool == net process duty; (iii) use no more hot or cold utility than documented in CASES (recorded at 1ab689ff by running the file directly), so a future change to hxn cannot silently make the synthesizer recover less heat. The recorded baselines also document where the greedy heuristic currently falls short of MER: cases 1-4 and 7 hit the target exactly; cases 8 and 9 are within 1 %; case 5 (boiling cold stream) uses 5.8e6 kJ/hr of hot utility where MER is 0, and case 6 uses 7.1e6 against an MER of 2.1e6. Case 10 is a strict xfail: synthesis emits a zero-duty exchanger (HX_8_2_cs, hot 355 K against cold 365 K, Q = 0) whose cold inlet is hotter than its hot inlet, and compute_LMTD then raises FloatingPointError on the log of a negative ratio. Fixing that flips the xfail into a failure that prompts recording the case's baseline. Validation: tests/test_hxn_regression.py -> 9 passed, 1 xfailed (5.7 s, not marked slow). Full CI-style suite -> 74 failed, 488 passed, 62 skipped, 1 xfailed; the 4 failures beyond the documented 70 (test_oilcane_O6/O8/O9: AttributeError inside the sibling biorefineries clone; test_tire: value drift, fails identically on unmodified master) are pre-existing stack drift since the baseline was measured against the pip thermosteam wheel. Co-Authored-By: Claude Fable 5 --- tests/test_hxn_regression.py | 239 +++++++++++++++++++++++++++++++++++ 1 file changed, 239 insertions(+) create mode 100644 tests/test_hxn_regression.py diff --git a/tests/test_hxn_regression.py b/tests/test_hxn_regression.py new file mode 100644 index 00000000..de0fb93f --- /dev/null +++ b/tests/test_hxn_regression.py @@ -0,0 +1,239 @@ +# -*- coding: utf-8 -*- +# BioSTEAM: The Biorefinery Simulation and Techno-Economic Analysis Modules +# Copyright (C) 2020-, Yoel Cortes-Pena +# Copyright (C) 2026-, Sarang Bhagwat +# +# This module is under the UIUC open-source license. See +# github.com/BioSTEAMDevelopmentGroup/biosteam/blob/master/LICENSE.txt +# for license details. +""" +Regression tests for heat exchanger network synthesis on synthetic systems. + +Ten synthetic systems of increasing complexity (all with phase-changing +streams from case 3 on). For each, the synthesized network must + +(i) close its energy balance (|error| < 0.1 %) without RuntimeWarnings, +(ii) never beat the minimum-energy-requirement (MER) targets of the problem + table computed on the same streams, and +(iii) recover at least as much heat as documented in ``CASES`` below, so that + no future change to ``hxn`` silently makes the synthesizer perform worse. + +The documented utility loads were recorded by running this file directly +(``python tests/test_hxn_regression.py`` prints them) at commit +``1ab689ff`` (branch ``hxn-pinch-diagram``). Improvements leave slack; a +maintainer lowers the numbers deliberately when a better network is +intended. Never raise them to make a failing test pass. +""" +import warnings +import pytest +import biosteam as bst +from numpy.testing import assert_allclose +from biosteam.facilities.hxn.hxn_synthesis import problem_table + +EB_TOLERANCE = 0.1 # percent +MER_RTOL = 1e-3 # network may not beat the MER target by more than this +DOC_RTOL = 1e-3 # network may not be worse than documented by more than this + +def utility_hx(ID, T, P, phase, T_out, rigorous=None, **flow): + """A simulated HXutility acting as one process stream (kmol/hr flows).""" + s = bst.Stream(ID + '_in', T=T, P=P, phase=phase, units='kmol/hr', **flow) + if rigorous is None: rigorous = phase == 'g' + hx = bst.HXutility(ID, ins=s, T=T_out, rigorous=rigorous) + hx.simulate() + return hx + +def boiling_hx(ID, T, P, T_out, **flow): + """A cold liquid stream heated past its bubble point (rigorous VLE).""" + return utility_hx(ID, T, P, 'l', T_out, rigorous=True, **flow) + +def setup(name, chemicals=('Water', 'Ethanol')): + bst.settings.set_thermo(list(chemicals), cache=True) + bst.main_flowsheet.set_flowsheet('test_hxn_regression_' + name) + +# --------------------------------------------------------------------------- +# Cases +# --------------------------------------------------------------------------- + +def case_01_two_liquids(): + """Hot liquid, cold liquid; trivial counter-current match.""" + setup('01') + return [utility_hx('H1', 400., 5e5, 'l', 320., Water=1000.), + utility_hx('C1', 300., 101325., 'l', 360., Water=1000.)], 5. + +def case_02_pinch_limited(): + """Cold target above the hot inlet: part of the heating must be utility.""" + setup('02') + return [utility_hx('H1', 360., 5e5, 'l', 320., Water=1000.), + utility_hx('C1', 300., 5e5, 'l', 380., Water=800.)], 5. + +def case_03_condenser_two_colds(): + """Condensing ethanol vapor against two cold liquids.""" + setup('03') + return [utility_hx('H1', 355., 101325., 'g', 340., Ethanol=400.), + utility_hx('C1', 300., 101325., 'l', 345., Water=1500.), + utility_hx('C2', 310., 101325., 'l', 340., Water=300., Ethanol=300.)], 5. + +def case_04_report_case(): + """The 4-stream report case: 2 colds, condensing hot, hot liquid.""" + setup('04') + return [utility_hx('C1', 300., 101325., 'l', 390., Water=2000.), + utility_hx('C2', 310., 101325., 'l', 345., Water=500., Ethanol=500.), + utility_hx('H1', 352., 101325., 'g', 340., Ethanol=300.), + utility_hx('H2', 420., 5e5, 'l', 320., Water=800.)], 5. + +def case_05_boiling_cold(): + """A cold stream that boils (water -> steam) and a condensing hot stream.""" + setup('05') + return [boiling_hx('C1', 330., 101325., 380., Water=300.), + utility_hx('C2', 300., 101325., 'l', 350., Water=1000.), + utility_hx('H1', 420., 5e5, 'g', 330., Water=250.), + utility_hx('H2', 400., 5e5, 'l', 310., Water=1500.)], 5. + +def case_06_mixed_pressures(): + """5-bar condensing hot against 1-atm boiling cold; T_min_app = 10.""" + setup('06') + return [utility_hx('H1', 430., 5e5, 'g', 400., Water=200.), + utility_hx('H2', 380., 5e5, 'l', 320., Water=1200.), + boiling_hx('C1', 340., 101325., 375., Water=150.), + utility_hx('C2', 300., 101325., 'l', 360., Ethanol=800.), + utility_hx('C3', 320., 101325., 'l', 390., Water=600.)], 10. + +def case_07_threshold(): + """Threshold problem: heating dominates so the cold target is ~0; + one hot stream condenses and subcools.""" + setup('07') + return [utility_hx('H1', 375., 101325., 'g', 310., Water=80.), + utility_hx('H2', 360., 101325., 'l', 330., Water=200.), + utility_hx('C1', 300., 101325., 'l', 370., Water=1500.), + utility_hx('C2', 305., 101325., 'l', 350., Ethanol=800.), + boiling_hx('C3', 340., 101325., 355., Ethanol=200.), + utility_hx('C4', 320., 101325., 'l', 360., Water=700.)], 5. + +def case_08_two_condensers(): + """Two condensers at different temperatures (ethanol 1 atm, water 2 bar) + against three colds, one of which boils.""" + setup('08') + return [utility_hx('H1', 355., 101325., 'g', 335., Ethanol=300.), + utility_hx('H2', 400., 2e5, 'g', 360., Water=150.), + utility_hx('H3', 390., 5e5, 'l', 330., Water=900.), + utility_hx('C1', 300., 101325., 'l', 345., Water=1200.), + boiling_hx('C2', 330., 101325., 370., Ethanol=250.), + utility_hx('C3', 310., 101325., 'l', 380., Water=700.)], 5. + +def case_09_near_degenerate(): + """Eight streams with two near-equal pinch candidates, a partially + condensing hot stream (wet outlet) and a partially boiling cold stream.""" + setup('09') + return [utility_hx('H1', 380., 101325., 'g', 372., Water=120.), # partial condensation + utility_hx('H2', 365., 101325., 'g', 330., Ethanol=250.), + utility_hx('H3', 410., 5e5, 'l', 340., Water=700.), + utility_hx('H4', 345., 101325., 'l', 305., Ethanol=900.), + boiling_hx('C1', 350., 101325., 373.5, Water=200.), # partial boiling + utility_hx('C2', 300., 101325., 'l', 340., Water=1500.), + boiling_hx('C3', 320., 101325., 352., Ethanol=300.), + utility_hx('C4', 335., 101325., 'l', 395., Water=500.)], 5. + +def case_10_ten_streams(): + """Ten streams mixing liquids, condensers, boilers, and pressures.""" + setup('10') + return [utility_hx('H1', 355., 101325., 'g', 320., Ethanol=300.), + utility_hx('H2', 420., 5e5, 'g', 340., Water=150.), + utility_hx('H3', 395., 5e5, 'l', 330., Water=1000.), + utility_hx('H4', 370., 101325., 'l', 310., Ethanol=700.), + utility_hx('H5', 380., 2e5, 'g', 375., Water=100.), # partial condensation + utility_hx('C1', 300., 101325., 'l', 360., Water=2000.), + boiling_hx('C2', 330., 101325., 380., Water=200.), + boiling_hx('C3', 320., 101325., 352., Ethanol=400.), + utility_hx('C4', 310., 101325., 'l', 345., Water=400., Ethanol=400.), + utility_hx('C5', 340., 101325., 'l', 390., Water=600.)], 5. + +# name -> (builder, documented hot utility load [kJ/hr], documented cold utility load [kJ/hr]) +# Documented values recorded at 1ab689ff (see module docstring). +CASES = { + 'case_01_two_liquids': (case_01_two_liquids, 0, 1.53912e+06), + 'case_02_pinch_limited': (case_02_pinch_limited, 1.81522e+06, 0), + 'case_03_condenser_two_colds': (case_03_condenser_two_colds, 0, 9.49905e+06), + 'case_04_report_case': (case_04_report_case, 2.37319e+06, 3.56871e+06), + 'case_05_boiling_cold': (case_05_boiling_cold, 5.8433e+06, 1.04356e+07), + 'case_06_mixed_pressures': (case_06_mixed_pressures, 7.05541e+06, 4.93079e+06), + 'case_07_threshold': (case_07_threshold, 1.85977e+07, 0), + 'case_08_two_condensers': (case_08_two_condensers, 3.02237e+06, 7.36431e+06), + 'case_09_near_degenerate': (case_09_near_degenerate, 1.40965e+07, 9.66427e+06), + 'case_10_ten_streams': (case_10_ten_streams, None, None), +} + +# --------------------------------------------------------------------------- +# Harness +# --------------------------------------------------------------------------- + +def synthesize(builder): + units, T_min_app = builder() + HXN = bst.HeatExchangerNetwork('HXN', T_min_app=T_min_app) + sys = bst.System.from_units('sys', units=[*units, HXN]) + with warnings.catch_warnings(): + warnings.simplefilter('error', RuntimeWarning) + # thermosteam registry bookkeeping on temporary stream copies; not numerical + warnings.filterwarnings('ignore', message='.*has been replaced in registry', + category=RuntimeWarning) + sys.simulate() + return units, HXN, T_min_app + +def mer_targets(units, T_min_app): + hus = [hx.heat_utilities[0] for hx in units] + hus.sort(key=lambda hu: hu.duty) + streams_inlet = [hu.unit.ins[0].copy() for hu in hus] + streams_quenched = [hu.unit.outs[0].copy() for hu in hus] + for s in streams_quenched: s.vle(H=s.H, P=s.P) + is_hot = [hu.duty < 0 for hu in hus] + table = problem_table(streams_inlet, streams_quenched, is_hot, T_min_app) + return table.hot_util_load, table.cold_util_load + +def actual_loads(HXN): + hus = [hu for hx in HXN.new_HX_utils for hu in hx.heat_utilities] + heat = sum(hu.unit_duty for hu in hus if hu.unit_duty > 0) + cool = -sum(hu.unit_duty for hu in hus if hu.unit_duty < 0) + return heat, cool + +# Cases that currently crash inside synthesis. Strict xfail: once the bug is +# fixed the test fails, prompting the maintainer to record the case's baseline. +KNOWN_FAILURES = { + # The synthesizer emits a zero-duty exchanger (HX_8_2_cs: H1 at 355 K + # against C3 at 365 K, Q = 0) whose cold inlet is hotter than its hot + # inlet; compute_LMTD then takes log of a negative ratio. + 'case_10_ten_streams': FloatingPointError, +} + +@pytest.mark.parametrize('name', [ + pytest.param(name, marks=pytest.mark.xfail(raises=KNOWN_FAILURES[name], strict=True, + reason='synthesizer emits a zero-duty exchanger with crossed temperatures')) + if name in KNOWN_FAILURES else name + for name in CASES +]) +def test_hxn_regression(name): + builder, doc_heat, doc_cool = CASES[name] + units, HXN, T_min_app = synthesize(builder) + # (i) energy balance + assert abs(HXN.energy_balance_percent_error) < EB_TOLERANCE, name + # (ii) MER targets are a lower bound; energy identity holds + heat, cool = actual_loads(HXN) + hot_target, cold_target = mer_targets(units, T_min_app) + net_duty = sum(hx.heat_utilities[0].unit_duty for hx in units) + assert heat >= hot_target * (1 - MER_RTOL), (name, heat, hot_target) + assert cool >= cold_target * (1 - MER_RTOL), (name, cool, cold_target) + assert_allclose(heat - cool, net_duty, rtol=1e-3, err_msg=name) + # (iii) never worse than documented + assert doc_heat is not None and doc_cool is not None, f'{name}: baseline not recorded' + assert heat <= doc_heat * (1 + DOC_RTOL) + 1e-9, (name, heat, doc_heat) + assert cool <= doc_cool * (1 + DOC_RTOL) + 1e-9, (name, cool, doc_cool) + +if __name__ == '__main__': + for name, (builder, *_) in CASES.items(): + if name in KNOWN_FAILURES: + print(f"{name}: known failure ({KNOWN_FAILURES[name].__name__})") + continue + units, HXN, T_min_app = synthesize(builder) + heat, cool = actual_loads(HXN) + hot_target, cold_target = mer_targets(units, T_min_app) + print(f"{name}: heat={heat:.6g} cool={cool:.6g} " + f"(MER hot={hot_target:.6g} cold={cold_target:.6g}; " + f"EB error={HXN.energy_balance_percent_error:.4f}%)") From 75f1402ac280f444fd60dea7661c1b971c6a729f Mon Sep 17 00:00:00 2001 From: sarangbhagwat Date: Sun, 23 Aug 2026 00:37:19 -0700 Subject: [PATCH 2/3] fix HXN synthesis for non-equilibrium inlets, stale network ordering, and H_lim at the bubble point Regression case 10 (tests/test_hxn_regression.py) crashed in compute_LMTD with a zero-duty exchanger whose cold inlet was hotter than its hot inlet. Tracing it back exposed three independent defects, each fixed at its cause: 1. Phantom pinch loads for non-equilibrium inlets (hxn_synthesis.py). load_duties() and the transient pinch-state setup in synthesize_network() flashed each stream at its pinch temperature and took that equilibrium enthalpy at face value. A non-rigorous HXutility can carry a superheated liquid (ethanol at 370 K, 1 atm, H = 6.4e6 kJ/hr); flashing it at the 355 K pinch yields vapor with H = 32.0e6, and abs(H_pinch - H_in) turned that into a 25.6e6 kJ/hr "cooling load" above the pinch for a stream whose whole duty is 5.45e6. Exchangers were then synthesized against a transient stream carrying latent heat the real stream does not have; rewired with the real inlet they behaved differently, overheating a cold stream and leaving a downstream match with crossed, zero-duty inlets. New pinch_state() flashes at the pinch T and, when the result lies outside [min(H_in, H_out), max(H_in, H_out)], re-flashes at the clipped enthalpy, so the hot-side and cold-side loads split |H_in - H_out| exactly by construction (the abs() is gone). This is the synthesizer's counterpart of the clipping the problem table got in 08e873c5. 2. HXN_sys path in synthesis order (_heat_exchanger_network.py). _cost built bst.System(path=new_HXs + new_HX_utils). Hot-side exchangers are synthesized before the cold-side exchangers that feed them after the life-cycle rewiring, and with no recycle declared the system ran the path once, leaving those exchangers with stale inlets (an energy-balance residual equal to the upstream duty change). The path is now built with tmo.Network.from_units(all_units, interaction=False) and System._from_network, which orders units by their connections and detects the genuine loops in larger networks for the fixed-point convergence that was already configured. interaction=False is required: HXprocess is an interaction unit and the default strips those out and disconnects them. 3. heat_exchange_to_condition ignored H_lim at the bubble point (heat_transfer.py). When the temperature limit lands within 1e-3 K of the stream's bubble point (e.g. T_hot_in - dT = 373.124 K for water at 1 atm) the outlet was set to the saturated phase and the enthalpy limit was never applied, so a cold stream's apparent capacity became its full vaporization and it overshot its enthalpy target; its "heating" utility then had to cool it back. The limit is now applied in that branch too, by VLE since the clipped state is two-phase. Also in counter_current_heat_exchange, both no-exchange early exits did s1_in.copy_like(s1_out), overwriting an inlet with the already modified outlet (observed rewriting a superheated-liquid inlet to vapor); corrected to s1_out.copy_like(s1_in). Tests: tests/test_heat_exchange.py (new) covers H_lim at the bubble point for heat_exchange_to_condition and HXprocess, and inlet immutability; tests/test_hxn.py gains a load_duties/pinch_state energy-conservation test for a non-equilibrium inlet. Case 10 loses its xfail and records its baseline (hot 1.40742e7, cold 8.06488e6 kJ/hr, both above MER, energy balance error 0.0000 %); the H5 input is corrected to 1 atm (it was a subcooled vapor at 2 bar). Cases 1-9 are bit-identical to their baselines. Validation: targeted HX/HXN tests and doctests -> 37 passed. Full CI-style suite -> 74 failed, 493 passed, 62 skipped: exactly the same 74 pre-existing failures as the previous run, +5 passed. Co-Authored-By: Claude Fable 5 --- .../facilities/hxn/_heat_exchanger_network.py | 11 ++- biosteam/facilities/hxn/hxn_synthesis.py | 42 ++++++++-- biosteam/units/design_tools/heat_transfer.py | 13 ++- tests/test_heat_exchange.py | 82 +++++++++++++++++++ tests/test_hxn.py | 30 ++++++- tests/test_hxn_regression.py | 32 ++------ 6 files changed, 175 insertions(+), 35 deletions(-) create mode 100644 tests/test_heat_exchange.py diff --git a/biosteam/facilities/hxn/_heat_exchanger_network.py b/biosteam/facilities/hxn/_heat_exchanger_network.py index 8e732a45..55be67ff 100644 --- a/biosteam/facilities/hxn/_heat_exchanger_network.py +++ b/biosteam/facilities/hxn/_heat_exchanger_network.py @@ -10,6 +10,7 @@ @author: sarangbhagwat and yoelcp """ import biosteam as bst +import thermosteam as tmo import numpy as np from .hxn_synthesis import synthesize_network, StreamLifeCycle, plot_pinch_diagram from warnings import warn @@ -225,7 +226,15 @@ def _cost(self): unit = i.unit if s_out: unit.ins[i.index] = s_out s_out = unit.outs[i.index] - self.HXN_sys = sys = bst.System(ID=None, path=all_units) + # Order the path by the rewired stream connections (and + # detect recycle loops) rather than using synthesis order: a + # hot-side exchanger is synthesized before the cold-side + # exchangers that feed it, and a single pass in synthesis + # order would leave it with stale inlets. HXprocess units are + # interaction units, which Network.from_units strips out (and + # disconnects) by default; keep them with interaction=False. + network = tmo.Network.from_units(all_units, interaction=False) + self.HXN_sys = sys = bst.System._from_network(None, network) sys.set_tolerance(method='fixedpoint', subsystems=True) original_purchase_costs = [hx.purchase_cost for hx in hxs] diff --git a/biosteam/facilities/hxn/hxn_synthesis.py b/biosteam/facilities/hxn/hxn_synthesis.py index 1041c205..8e633449 100644 --- a/biosteam/facilities/hxn/hxn_synthesis.py +++ b/biosteam/facilities/hxn/hxn_synthesis.py @@ -368,16 +368,39 @@ def temperature_interval_pinch_analysis(hus, streams_quenched +def pinch_state(stream_in, H_in, H_out, T_pinch): + """ + Return a copy of `stream_in` in the state it has when it crosses the + pinch, with enthalpy guaranteed to lie within [min(H_in, H_out), + max(H_in, H_out)]. + + The copy is first flashed at `T_pinch`. If the resulting equilibrium + enthalpy lies outside the stream's own enthalpy range the stream never + actually passes through that state: a non-equilibrium inlet (e.g. a + superheated liquid from a non-rigorous HXutility) or outlet has less + (or more) enthalpy than the equilibrium fluid at the pinch. The copy is + then flashed at the nearer end of the range instead, so the hot-side + and cold-side loads split `|H_in - H_out|` exactly and the transient + stream used for matching never carries heat the real stream does not + have. This is the synthesizer's counterpart of the clipping done by + `_stream_H_at_boundaries` for the problem table. + """ + stream = stream_in.copy() + stream.vle(T=T_pinch, P=stream.P) + H_lo, H_hi = sorted((H_in, H_out)) + H = stream.H + if H < H_lo or H > H_hi: + stream.vle(H=min(max(H, H_lo), H_hi), P=stream.P) + return stream + def load_duties(streams, streams_quenched, pinch_T_arr, T_out_arr, indices, is_cold, Q_hot_side, Q_cold_side): for index in indices: - stream = streams[index].copy() - H_in = stream.H - stream.vle(T = pinch_T_arr[index], P = stream.P) - H_pinch = stream.H + H_in = streams[index].H H_out = streams_quenched[index].H + H_pinch = pinch_state(streams[index], H_in, H_out, pinch_T_arr[index]).H if not is_cold(index): - dH1 = abs(H_pinch - H_in) - dH2 = abs(H_out - H_pinch) + dH1 = H_in - H_pinch + dH2 = H_pinch - H_out if abs(dH1)<0.01: dH1 = 0 if abs(dH2)<0.01: dH2 = 0 Q_hot_side[index] = ['cool', dH1] @@ -418,14 +441,17 @@ def synthesize_network(hus, T_min_app=5., Qmin=1e-3, force_ideal_thermo=False, HXs_cold_side = [] streams_transient_cold_side = streams_inlet streams_transient_hot_side = [i.copy() for i in streams_inlet] + # Hot streams enter the cold-side design at their pinch state and cold + # streams enter the hot-side design at theirs; the enthalpy of that + # state is clipped to the stream's real range (see `pinch_state`). for i in hot_indices: s = streams_transient_cold_side[i] if s.T != pinch_T_arr[i]: - s.vle(T=pinch_T_arr[i], P=s.P) + streams_transient_cold_side[i] = pinch_state(s, s.H, H_out_arr[i], pinch_T_arr[i]) for i in cold_indices: s = streams_transient_hot_side[i] if s.T != pinch_T_arr[i]: - s.vle(T=pinch_T_arr[i], P=s.P) + streams_transient_hot_side[i] = pinch_state(s, s.H, H_out_arr[i], pinch_T_arr[i]) def get_stream_at_H_max(cold): s_cs = streams_transient_cold_side[cold] diff --git a/biosteam/units/design_tools/heat_transfer.py b/biosteam/units/design_tools/heat_transfer.py index 8ade3361..af53ebe9 100644 --- a/biosteam/units/design_tools/heat_transfer.py +++ b/biosteam/units/design_tools/heat_transfer.py @@ -77,7 +77,16 @@ def heat_exchange_to_condition(s_in, s_out, T=None, phase=None, else: if s_out.H < H_lim: s_out.H = H_lim else: + # At the bubble point: the most the stream can absorb (release) + # without leaving T is full vaporization (condensation). The + # enthalpy limit still applies; a limit short of the full phase + # change lands in the two-phase region, so solve it by VLE. s_out.phase = 'g' if heating else 'l' + if H_lim_given: + if heating: + if s_out.H > H_lim: s_out.vle(H=H_lim, P=s_out.P) + else: + if s_out.H < H_lim: s_out.vle(H=H_lim, P=s_out.P) else: s_out.vle(T=T, P=s_out.P) if H_lim_given: @@ -168,14 +177,14 @@ def counter_current_heat_exchange(s0_in, s1_in, s0_out, s1_out, if Q_hot_stream == Q_cold_stream == 0.: s0_out.copy_like(s0_in) - s1_in.copy_like(s1_out) + s1_out.copy_like(s1_in) return 0. if Q_hot_stream > 0 or Q_cold_stream < 0: # Sanity check if Q_hot_stream / s_hot_in.C < 0.1 or Q_cold_stream / s_cold_in.C > -0.1: s0_out.copy_like(s0_in) - s1_in.copy_like(s1_out) + s1_out.copy_like(s1_in) return 0. raise RuntimeError('inlet stream not in vapor-liquid equilibrium') diff --git a/tests/test_heat_exchange.py b/tests/test_heat_exchange.py new file mode 100644 index 00000000..5843498c --- /dev/null +++ b/tests/test_heat_exchange.py @@ -0,0 +1,82 @@ +# -*- coding: utf-8 -*- +# BioSTEAM: The Biorefinery Simulation and Techno-Economic Analysis Modules +# Copyright (C) 2020-, Yoel Cortes-Pena +# Copyright (C) 2026-, Sarang Bhagwat +# +# This module is under the UIUC open-source license. See +# github.com/BioSTEAMDevelopmentGroup/biosteam/blob/master/LICENSE.txt +# for license details. +""" +Tests for heat exchanger units and the counter-current heat exchange solver. +""" +import biosteam as bst +from numpy.testing import assert_allclose +from biosteam.units.design_tools.heat_transfer import heat_exchange_to_condition + +def test_heat_exchange_to_condition_respects_H_lim_at_bubble_point(): + # When the temperature limit coincides with the stream's bubble point + # (within the solver's 1e-3 K tolerance), the outlet is set to the + # saturated phase; the enthalpy limit must still be honored. + bst.settings.set_thermo(['Water', 'Ethanol'], cache=True) + s_in = bst.Stream('w_in', Water=2000., T=350., P=101325., phase='l', + units='kmol/hr') + T_bp = s_in.bubble_point_at_P().T + s_lim = s_in.copy(); s_lim.T = 360. + H_lim = s_lim.H + # heating: limit below full vaporization at the bubble point + s_out = s_in.copy() + Q = heat_exchange_to_condition(s_in, s_out, T=T_bp, H_lim=H_lim, heating=True) + assert_allclose(s_out.H, H_lim, rtol=1e-9) + assert_allclose(Q, H_lim - s_in.H, rtol=1e-9) + assert s_out.phase == 'l' and s_out.T < T_bp + # cooling: a saturated vapor whose limit is above full condensation + v_in = s_in.copy('w_vap'); v_in.phase = 'g'; v_in.T = T_bp + s_out = v_in.copy() + H_lim = v_in.H - 0.5 * (v_in.H - s_in.H) # half-way to liquid at 350 K + Q = heat_exchange_to_condition(v_in, s_out, T=T_bp, H_lim=H_lim, heating=False) + assert_allclose(s_out.H, H_lim, rtol=1e-9) + assert_allclose(Q, H_lim - v_in.H, rtol=1e-9) + +def test_HXprocess_H_lim_when_pinch_is_at_bubble_point(): + # Hot liquid water at exactly T_bp + dT against cold liquid water with an + # enthalpy limit below vaporization: the cold outlet may not exceed its + # enthalpy limit just because its temperature limit lands on the + # bubble point. + bst.settings.set_thermo(['Water', 'Ethanol'], cache=True) + cold = bst.Stream('cold', Water=2000., T=350., P=101325., phase='l', + units='kmol/hr') + T_bp = cold.bubble_point_at_P().T + hot = bst.Stream('hot', Water=1000., T=T_bp + 5., P=5e5, phase='l', + units='kmol/hr') + s_lim = cold.copy(); s_lim.T = 360. + H_lim = s_lim.H + hx = bst.HXprocess('hx', ins=(cold, hot), H_lim0=H_lim, T_lim1=355., dT=5.) + hx.simulate() + cold_out, hot_out = hx.outs + assert_allclose(cold_out.H, H_lim, rtol=1e-9) + assert_allclose(hx.Q, H_lim - cold.H, rtol=1e-9) + assert_allclose(hot_out.H, hot.H - hx.Q, rtol=1e-9) + assert hot_out.T > 355. # the hot stream was not the limiting side + +def test_HXprocess_never_modifies_inlets(): + # A superheated-liquid hot inlet (ethanol at 370 K, 1 atm; bp 351.4 K) + # against a two-phase cold inlet of the same fluid: the solver finds no + # feasible exchange. It must leave both inlet streams untouched rather + # than overwrite one with its (already modified) outlet. + bst.settings.set_thermo(['Water', 'Ethanol'], cache=True) + hot = bst.Stream('hot', Ethanol=700., T=370., P=101325., phase='l', + units='kmol/hr') + cold = bst.Stream('cold', Ethanol=400., T=348.6, P=101325., phase='l', + units='kmol/hr') + cold.vle(H=cold.H + 4e6, P=101325.) # two-phase at the bubble point + H_hot, H_cold = hot.H, cold.H + hx = bst.HXprocess('hx', ins=(cold, hot), dT=5.) + hx._run() + assert hx.Q == 0. + assert_allclose([hot.H, cold.H], [H_hot, H_cold], rtol=1e-12) + assert hot.phase == 'l' and hot.T == 370. + +if __name__ == '__main__': + test_heat_exchange_to_condition_respects_H_lim_at_bubble_point() + test_HXprocess_H_lim_when_pinch_is_at_bubble_point() + test_HXprocess_never_modifies_inlets() diff --git a/tests/test_hxn.py b/tests/test_hxn.py index 214372d5..151db4b3 100644 --- a/tests/test_hxn.py +++ b/tests/test_hxn.py @@ -96,7 +96,7 @@ def test_energy_balance_error_contributions_ignored_none(): # --------------------------------------------------------------------------- from biosteam.facilities.hxn.hxn_synthesis import ( - temperature_interval_pinch_analysis, problem_table, + temperature_interval_pinch_analysis, problem_table, load_duties, pinch_state, ) def utility_hx(ID, T, P, phase, T_out, **flow): @@ -236,6 +236,34 @@ def test_problem_table_point_load_cannot_heat_above_itself(): assert table.pinch_T == 395. assert_allclose(table.cold_util_load, P_hot - (H_395 - H_392), rtol=1e-9) +def test_load_duties_non_equilibrium_inlet_conserves_energy(): + # A non-rigorous HXutility can carry a superheated liquid (ethanol at + # 370 K, 1 atm; bp 351.4 K). Flashing it at the 355 K pinch gives vapor + # with far more enthalpy than the inlet has; the pinch split must clip + # to the stream's real enthalpy range so that the hot-side and cold-side + # loads sum to the stream's duty (and are not phantom latent heat). + bst.settings.set_thermo(['Water', 'Ethanol'], cache=True) + s_in = bst.Stream('ne_in', Ethanol=700., T=370., P=101325., phase='l', + units='kmol/hr') + s_out = bst.Stream('ne_out', Ethanol=700., T=310., P=101325., phase='l', + units='kmol/hr') + duty = s_in.H - s_out.H + flashed = s_in.copy(); flashed.vle(T=355., P=101325.) + assert flashed.H > s_in.H # the trap: equilibrium enthalpy above H_in + Q_hot_side, Q_cold_side = {}, {} + load_duties([s_in], [s_out], np.array([355.]), np.array([310.]), [0], + lambda i: False, Q_hot_side, Q_cold_side) + assert Q_hot_side[0][0] == Q_cold_side[0][0] == 'cool' + assert_allclose(Q_hot_side[0][1] + Q_cold_side[0][1], duty, rtol=1e-9) + # the stream cannot give any heat above the pinch: its whole enthalpy + # range lies below the equilibrium state at 355 K + assert_allclose(Q_hot_side[0][1], 0., atol=1e-6) + assert_allclose(Q_cold_side[0][1], duty, rtol=1e-9) + # and the transient state used for matching carries the clipped enthalpy + s = pinch_state(s_in, s_in.H, s_out.H, 355.) + assert_allclose(s.H, s_in.H, rtol=1e-9) + assert s.T < 355. + def test_synthetic_network_reaches_MER(): units = synthetic_units() HXN = bst.HeatExchangerNetwork('HXN', T_min_app=5.) diff --git a/tests/test_hxn_regression.py b/tests/test_hxn_regression.py index de0fb93f..173321e3 100644 --- a/tests/test_hxn_regression.py +++ b/tests/test_hxn_regression.py @@ -19,8 +19,11 @@ no future change to ``hxn`` silently makes the synthesizer perform worse. The documented utility loads were recorded by running this file directly -(``python tests/test_hxn_regression.py`` prints them) at commit -``1ab689ff`` (branch ``hxn-pinch-diagram``). Improvements leave slack; a +(``python tests/test_hxn_regression.py`` prints them): cases 1-9 at commit +``1ab689ff`` (branch ``hxn-pinch-diagram``), case 10 after the synthesizer +fixes on ``hxn-regression-tests`` (non-equilibrium inlets clipped to the +stream's enthalpy range; network path ordered by its connections; H_lim +honored at the bubble point). Improvements leave slack; a maintainer lowers the numbers deliberately when a better network is intended. Never raise them to make a failing test pass. """ @@ -140,7 +143,7 @@ def case_10_ten_streams(): utility_hx('H2', 420., 5e5, 'g', 340., Water=150.), utility_hx('H3', 395., 5e5, 'l', 330., Water=1000.), utility_hx('H4', 370., 101325., 'l', 310., Ethanol=700.), - utility_hx('H5', 380., 2e5, 'g', 375., Water=100.), # partial condensation + utility_hx('H5', 380., 101325., 'g', 372.5, Water=100.), # partial condensation utility_hx('C1', 300., 101325., 'l', 360., Water=2000.), boiling_hx('C2', 330., 101325., 380., Water=200.), boiling_hx('C3', 320., 101325., 352., Ethanol=400.), @@ -148,7 +151,7 @@ def case_10_ten_streams(): utility_hx('C5', 340., 101325., 'l', 390., Water=600.)], 5. # name -> (builder, documented hot utility load [kJ/hr], documented cold utility load [kJ/hr]) -# Documented values recorded at 1ab689ff (see module docstring). +# Documented values: see module docstring for provenance. CASES = { 'case_01_two_liquids': (case_01_two_liquids, 0, 1.53912e+06), 'case_02_pinch_limited': (case_02_pinch_limited, 1.81522e+06, 0), @@ -159,7 +162,7 @@ def case_10_ten_streams(): 'case_07_threshold': (case_07_threshold, 1.85977e+07, 0), 'case_08_two_condensers': (case_08_two_condensers, 3.02237e+06, 7.36431e+06), 'case_09_near_degenerate': (case_09_near_degenerate, 1.40965e+07, 9.66427e+06), - 'case_10_ten_streams': (case_10_ten_streams, None, None), + 'case_10_ten_streams': (case_10_ten_streams, 1.40742e+07, 8.06488e+06), } # --------------------------------------------------------------------------- @@ -194,21 +197,7 @@ def actual_loads(HXN): cool = -sum(hu.unit_duty for hu in hus if hu.unit_duty < 0) return heat, cool -# Cases that currently crash inside synthesis. Strict xfail: once the bug is -# fixed the test fails, prompting the maintainer to record the case's baseline. -KNOWN_FAILURES = { - # The synthesizer emits a zero-duty exchanger (HX_8_2_cs: H1 at 355 K - # against C3 at 365 K, Q = 0) whose cold inlet is hotter than its hot - # inlet; compute_LMTD then takes log of a negative ratio. - 'case_10_ten_streams': FloatingPointError, -} - -@pytest.mark.parametrize('name', [ - pytest.param(name, marks=pytest.mark.xfail(raises=KNOWN_FAILURES[name], strict=True, - reason='synthesizer emits a zero-duty exchanger with crossed temperatures')) - if name in KNOWN_FAILURES else name - for name in CASES -]) +@pytest.mark.parametrize('name', list(CASES)) def test_hxn_regression(name): builder, doc_heat, doc_cool = CASES[name] units, HXN, T_min_app = synthesize(builder) @@ -228,9 +217,6 @@ def test_hxn_regression(name): if __name__ == '__main__': for name, (builder, *_) in CASES.items(): - if name in KNOWN_FAILURES: - print(f"{name}: known failure ({KNOWN_FAILURES[name].__name__})") - continue units, HXN, T_min_app = synthesize(builder) heat, cool = actual_loads(HXN) hot_target, cold_target = mer_targets(units, T_min_app) From 1d882c4b3040bbbc9049273ca9318abe968a27e3 Mon Sep 17 00:00:00 2001 From: sarangbhagwat Date: Sun, 23 Aug 2026 03:27:24 -0700 Subject: [PATCH 3/3] address review: exact pinch states at end temperatures, HXN-specific path warning, tighter regression tolerances Follow-up to 75f1402a from a code review of the branch. pinch_state (hxn_synthesis.py) now takes the stream's real end states and, when the pinch temperature coincides with an end temperature or the flash at the pinch leaves the stream's enthalpy range, returns the *equilibrium state at that end enthalpy* (new helper _end_state). Two defects closed: - A non-equilibrium inlet with pinch_T == T_in still had its duty split across the pinch: flashing a superheated liquid at its own temperature does not reproduce H_in, and the result can lie inside the range so the clip never fired (superheated liquid water 380 -> 400 K put 98 % of the duty on the side the stream is "unavailable" on). Now H_pinch == H_in exactly, so the whole duty lands on one side as the pinch analysis promises. - A phase-mislabelled non-condensable (N2 as 'l' at 400 K) re-flashed at the clipped enthalpy came back as gas at ~2500 K and drove matching. _end_state falls back to the stream as given when the equilibrium state leaves the stream's own [T_out, T_in] range. The review suggested returning the stream as given at end points; that was tried and rejected: it lets the synthesizer offer a non-equilibrium inlet's heat at its fictitious temperature (liquid ethanol "at 370 K") and case 10 then beat the problem-table MER by 0.3 %. The equilibrium state at the end enthalpy keeps synthesis consistent with the table, which is the invariant the regression suite enforces. _heat_exchanger_network.py: thermosteam's Network.sort can warn 'network path could not be determined' when its ordering heuristic does not settle; since 75f1402a that could surface from HXN with no context. It is caught and re-issued as an HXN-specific RuntimeWarning (tested by monkeypatching Network.from_units). tests/test_hxn_regression.py: converged networks close their energy balance to ~1e-10 %, so the 0.1 % tolerance could not catch the stale-ordering defect (ablating it gives 0.049 %). EB_TOLERANCE is now 1e-6 % and the heat - cool == net duty identity uses rtol 1e-8. Case 5's documented baseline is lowered (hot 5.8433e6 -> 3.2224e6, cold 1.04356e7 -> 7.81466e6): its 420 K, 5 bar vapor feed is below water's boiling point there, a non-equilibrium inlet now handled consistently. Cases 1-4 and 6-10 are unchanged. Validation: HX/HXN tests and doctests -> 37 passed. Full CI-style suite -> 74 failed, 495 passed, 62 skipped; the same 74 pre-existing failures as the previous two runs, +2 passed. Co-Authored-By: Claude Fable 5 --- .../facilities/hxn/_heat_exchanger_network.py | 13 +++- biosteam/facilities/hxn/hxn_synthesis.py | 68 ++++++++++++++----- tests/test_hxn.py | 59 +++++++++++++++- tests/test_hxn_regression.py | 21 +++--- 4 files changed, 130 insertions(+), 31 deletions(-) diff --git a/biosteam/facilities/hxn/_heat_exchanger_network.py b/biosteam/facilities/hxn/_heat_exchanger_network.py index 55be67ff..fbb3c447 100644 --- a/biosteam/facilities/hxn/_heat_exchanger_network.py +++ b/biosteam/facilities/hxn/_heat_exchanger_network.py @@ -14,6 +14,7 @@ import numpy as np from .hxn_synthesis import synthesize_network, StreamLifeCycle, plot_pinch_diagram from warnings import warn +import warnings __all__ = ('HeatExchangerNetwork',) @@ -233,7 +234,17 @@ def _cost(self): # order would leave it with stale inlets. HXprocess units are # interaction units, which Network.from_units strips out (and # disconnects) by default; keep them with interaction=False. - network = tmo.Network.from_units(all_units, interaction=False) + with warnings.catch_warnings(record=True) as caught: + warnings.simplefilter('always', RuntimeWarning) + network = tmo.Network.from_units(all_units, interaction=False) + for w in caught: + if 'network path could not be determined' in str(w.message): + warn('heat exchanger network path could not be fully ' + 'ordered from its stream connections; exchangers ' + 'fed by later ones in the path may be simulated ' + 'with stale inlets until convergence', RuntimeWarning) + else: + warn(w.message, w.category) self.HXN_sys = sys = bst.System._from_network(None, network) sys.set_tolerance(method='fixedpoint', subsystems=True) diff --git a/biosteam/facilities/hxn/hxn_synthesis.py b/biosteam/facilities/hxn/hxn_synthesis.py index 8e633449..c47dde15 100644 --- a/biosteam/facilities/hxn/hxn_synthesis.py +++ b/biosteam/facilities/hxn/hxn_synthesis.py @@ -368,36 +368,68 @@ def temperature_interval_pinch_analysis(hus, streams_quenched -def pinch_state(stream_in, H_in, H_out, T_pinch): +def _end_state(stream_end, T_lo, T_hi): """ - Return a copy of `stream_in` in the state it has when it crosses the + Return a copy of `stream_end` at equilibrium at its own enthalpy, or the + stream as given if that equilibrium state lies outside the stream's own + temperature range [T_lo, T_hi] (e.g. a non-condensable mislabelled as a + liquid, whose equilibrium state at the same enthalpy is a gas at an + absurd temperature). Either way the enthalpy is exactly `stream_end.H`. + """ + stream = stream_end.copy() + try: + stream.vle(H=stream_end.H, P=stream.P) + except Exception: + return stream_end.copy() + if T_lo <= stream.T <= T_hi: return stream + return stream_end.copy() + +def pinch_state(stream_in, stream_out, T_pinch): + """ + Return a copy of the stream in the state it has when it crosses the pinch, with enthalpy guaranteed to lie within [min(H_in, H_out), max(H_in, H_out)]. - The copy is first flashed at `T_pinch`. If the resulting equilibrium - enthalpy lies outside the stream's own enthalpy range the stream never - actually passes through that state: a non-equilibrium inlet (e.g. a - superheated liquid from a non-rigorous HXutility) or outlet has less - (or more) enthalpy than the equilibrium fluid at the pinch. The copy is - then flashed at the nearer end of the range instead, so the hot-side - and cold-side loads split `|H_in - H_out|` exactly and the transient - stream used for matching never carries heat the real stream does not - have. This is the synthesizer's counterpart of the clipping done by - `_stream_H_at_boundaries` for the problem table. + `stream_in` and `stream_out` are the stream's real end states (the + outlet quenched to equilibrium at its own enthalpy). For an interior + pinch the inlet copy is flashed at `T_pinch`; the result is used as is + when its enthalpy lies within the stream's own range. Otherwise the + stream never passes through that equilibrium state: a non-equilibrium + inlet (e.g. a superheated liquid from a non-rigorous HXutility) has + less enthalpy than the equilibrium fluid at the pinch, and the state + returned is instead the equilibrium state at the nearer end enthalpy. + The same end state is returned when `T_pinch` coincides with an end + temperature, because flashing a non-equilibrium inlet at its own + temperature does not reproduce `H_in` (and the result may even lie + inside the range). Using the *equilibrium* state at the end enthalpy, + rather than the stream as given, keeps the synthesizer consistent with + the problem table: the heat is offered at the temperature the + equilibrium model says it is available, not at a fictitious one; see + `_end_state` for the fallback when that state is unphysical. + + Either way the hot-side and cold-side loads split `|H_in - H_out|` + exactly and the transient stream used for matching never carries heat + the real stream does not have. This is the synthesizer's counterpart of + the clipping done by `_stream_H_at_boundaries` for the problem table. """ + T_lo, T_hi = sorted((stream_in.T, stream_out.T)) + if T_pinch == stream_in.T: return _end_state(stream_in, T_lo, T_hi) + if T_pinch == stream_out.T: return _end_state(stream_out, T_lo, T_hi) stream = stream_in.copy() stream.vle(T=T_pinch, P=stream.P) + H_in = stream_in.H + H_out = stream_out.H H_lo, H_hi = sorted((H_in, H_out)) H = stream.H - if H < H_lo or H > H_hi: - stream.vle(H=min(max(H, H_lo), H_hi), P=stream.P) - return stream + if H_lo <= H <= H_hi: return stream + H_clipped = H_lo if H < H_lo else H_hi + return _end_state(stream_in if H_clipped == H_in else stream_out, T_lo, T_hi) def load_duties(streams, streams_quenched, pinch_T_arr, T_out_arr, indices, is_cold, Q_hot_side, Q_cold_side): for index in indices: H_in = streams[index].H H_out = streams_quenched[index].H - H_pinch = pinch_state(streams[index], H_in, H_out, pinch_T_arr[index]).H + H_pinch = pinch_state(streams[index], streams_quenched[index], pinch_T_arr[index]).H if not is_cold(index): dH1 = H_in - H_pinch dH2 = H_pinch - H_out @@ -447,11 +479,11 @@ def synthesize_network(hus, T_min_app=5., Qmin=1e-3, force_ideal_thermo=False, for i in hot_indices: s = streams_transient_cold_side[i] if s.T != pinch_T_arr[i]: - streams_transient_cold_side[i] = pinch_state(s, s.H, H_out_arr[i], pinch_T_arr[i]) + streams_transient_cold_side[i] = pinch_state(s, streams_quenched[i], pinch_T_arr[i]) for i in cold_indices: s = streams_transient_hot_side[i] if s.T != pinch_T_arr[i]: - streams_transient_hot_side[i] = pinch_state(s, s.H, H_out_arr[i], pinch_T_arr[i]) + streams_transient_hot_side[i] = pinch_state(s, streams_quenched[i], pinch_T_arr[i]) def get_stream_at_H_max(cold): s_cs = streams_transient_cold_side[cold] diff --git a/tests/test_hxn.py b/tests/test_hxn.py index 151db4b3..4e7a147b 100644 --- a/tests/test_hxn.py +++ b/tests/test_hxn.py @@ -259,10 +259,63 @@ def test_load_duties_non_equilibrium_inlet_conserves_energy(): # range lies below the equilibrium state at 355 K assert_allclose(Q_hot_side[0][1], 0., atol=1e-6) assert_allclose(Q_cold_side[0][1], duty, rtol=1e-9) - # and the transient state used for matching carries the clipped enthalpy - s = pinch_state(s_in, s_in.H, s_out.H, 355.) + # and the transient state used for matching carries exactly H_in, at + # equilibrium (two-phase at the bubble point), not at the fictitious + # 370 K of the superheated liquid, so matching is consistent with the + # problem table + s = pinch_state(s_in, s_out, 355.) assert_allclose(s.H, s_in.H, rtol=1e-9) - assert s.T < 355. + assert s.T < 355. and 'g' in s.phase and 'l' in s.phase + +def test_pinch_state_at_endpoints_uses_real_states(): + # pinch_T == T_in must put the whole duty on one side even when the + # inlet is not at equilibrium: flashing a superheated liquid at its own + # T does not reproduce H_in, and the resulting enthalpy can lie inside + # the stream's range so that clipping alone would not catch it. + bst.settings.set_thermo(['Water', 'Ethanol'], cache=True) + s_in = bst.Stream('sh_in', Water=100., T=380., P=101325., phase='l', + units='kmol/hr') + s_out = bst.Stream('sh_out', Water=100., T=400., P=101325., phase='g', + units='kmol/hr') + duty = s_out.H - s_in.H + flashed = s_in.copy(); flashed.vle(T=380., P=101325.) + assert s_in.H < flashed.H < s_out.H # the trap: in range, but not H_in + Q_hot_side, Q_cold_side = {}, {} + load_duties([s_in], [s_out], np.array([380.]), np.array([400.]), [0], + lambda i: True, Q_hot_side, Q_cold_side) + assert_allclose(Q_hot_side[0][1], duty, rtol=1e-9) + assert_allclose(Q_cold_side[0][1], 0., atol=1e-6) + assert_allclose(pinch_state(s_in, s_out, 380.).H, s_in.H, rtol=1e-9) + assert_allclose(pinch_state(s_in, s_out, 400.).H, s_out.H, rtol=1e-9) + # a phase-mislabelled non-condensable: clipping must return the real + # end state, never an equilibrium re-flash (N2 "liquid" at 400 K would + # otherwise come back as gas at thousands of K) + bst.settings.set_thermo(['Water', 'N2'], cache=True) + h_in = bst.Stream('n2_in', N2=100., T=400., P=101325., phase='l', + units='kmol/hr') + h_out = bst.Stream('n2_out', N2=100., T=320., P=101325., phase='l', + units='kmol/hr') + s = pinch_state(h_in, h_out, 355.) + assert h_out.H <= s.H <= h_in.H + assert 320. <= s.T <= 400. + +def test_unordered_network_path_warns_with_context(monkeypatch): + # thermosteam's Network.sort warns 'network path could not be determined' + # when its ordering heuristic does not settle; HXN re-raises that as its + # own warning so the user knows which facility it concerns. + import thermosteam as tmo + original = tmo.Network.from_units + def from_units_with_warning(*args, **kwargs): + network = original(*args, **kwargs) + warnings.warn('network path could not be determined', RuntimeWarning) + return network + monkeypatch.setattr(tmo.Network, 'from_units', from_units_with_warning) + units = synthetic_units() + HXN = bst.HeatExchangerNetwork('HXN', T_min_app=5.) + sys = bst.System.from_units('sys_unordered', units=[*units, HXN]) + with pytest.warns(RuntimeWarning, match='heat exchanger network path could not be fully ordered'): + sys.simulate() + assert abs(HXN.energy_balance_percent_error) < 1e-6 def test_synthetic_network_reaches_MER(): units = synthetic_units() diff --git a/tests/test_hxn_regression.py b/tests/test_hxn_regression.py index 173321e3..7c2d1864 100644 --- a/tests/test_hxn_regression.py +++ b/tests/test_hxn_regression.py @@ -12,18 +12,21 @@ Ten synthetic systems of increasing complexity (all with phase-changing streams from case 3 on). For each, the synthesized network must -(i) close its energy balance (|error| < 0.1 %) without RuntimeWarnings, +(i) close its energy balance (|error| < 1e-6 %) without RuntimeWarnings, (ii) never beat the minimum-energy-requirement (MER) targets of the problem table computed on the same streams, and (iii) recover at least as much heat as documented in ``CASES`` below, so that no future change to ``hxn`` silently makes the synthesizer perform worse. The documented utility loads were recorded by running this file directly -(``python tests/test_hxn_regression.py`` prints them): cases 1-9 at commit -``1ab689ff`` (branch ``hxn-pinch-diagram``), case 10 after the synthesizer -fixes on ``hxn-regression-tests`` (non-equilibrium inlets clipped to the -stream's enthalpy range; network path ordered by its connections; H_lim -honored at the bubble point). Improvements leave slack; a +(``python tests/test_hxn_regression.py`` prints them): cases 1-4 and 6-9 +at commit ``1ab689ff`` (branch ``hxn-pinch-diagram``); case 10 after the +synthesizer fixes on ``hxn-regression-tests`` (non-equilibrium inlets +clipped to the stream's enthalpy range; network path ordered by its +connections; H_lim honored at the bubble point); case 5 lowered after the +pinch state at an end temperature became the equilibrium state at that +end enthalpy (its 420 K, 5 bar vapor feed is below water's boiling point +there, a non-equilibrium inlet). Improvements leave slack; a maintainer lowers the numbers deliberately when a better network is intended. Never raise them to make a failing test pass. """ @@ -33,7 +36,7 @@ from numpy.testing import assert_allclose from biosteam.facilities.hxn.hxn_synthesis import problem_table -EB_TOLERANCE = 0.1 # percent +EB_TOLERANCE = 1e-6 # percent; converged networks close to ~1e-10 % MER_RTOL = 1e-3 # network may not beat the MER target by more than this DOC_RTOL = 1e-3 # network may not be worse than documented by more than this @@ -157,7 +160,7 @@ def case_10_ten_streams(): 'case_02_pinch_limited': (case_02_pinch_limited, 1.81522e+06, 0), 'case_03_condenser_two_colds': (case_03_condenser_two_colds, 0, 9.49905e+06), 'case_04_report_case': (case_04_report_case, 2.37319e+06, 3.56871e+06), - 'case_05_boiling_cold': (case_05_boiling_cold, 5.8433e+06, 1.04356e+07), + 'case_05_boiling_cold': (case_05_boiling_cold, 3.2224e+06, 7.81466e+06), 'case_06_mixed_pressures': (case_06_mixed_pressures, 7.05541e+06, 4.93079e+06), 'case_07_threshold': (case_07_threshold, 1.85977e+07, 0), 'case_08_two_condensers': (case_08_two_condensers, 3.02237e+06, 7.36431e+06), @@ -209,7 +212,7 @@ def test_hxn_regression(name): net_duty = sum(hx.heat_utilities[0].unit_duty for hx in units) assert heat >= hot_target * (1 - MER_RTOL), (name, heat, hot_target) assert cool >= cold_target * (1 - MER_RTOL), (name, cool, cold_target) - assert_allclose(heat - cool, net_duty, rtol=1e-3, err_msg=name) + assert_allclose(heat - cool, net_duty, rtol=1e-8, err_msg=name) # (iii) never worse than documented assert doc_heat is not None and doc_cool is not None, f'{name}: baseline not recorded' assert heat <= doc_heat * (1 + DOC_RTOL) + 1e-9, (name, heat, doc_heat)