Coverage for backend/ahuora-builder/src/ahuora_builder/custom/thermal_utility_systems/header_OLD.py: 17%
231 statements
« prev ^ index » next coverage.py v7.10.7, created at 2026-07-22 05:22 +0000
1# Pyomo core
2##Previous header implementation, new header implementation in header.py, kept for reference in case of any issues with the new implementation
3import pyomo.environ as pyo
4from pyomo.environ import (
5 Constraint,
6 Expression,
7 Param,
8 PositiveReals,
9 RangeSet,
10 Suffix,
11 Var,
12 value,
13 units as UNIT,
14)
15from pyomo.core.base.reference import Reference
16from pyomo.common.config import ConfigBlock, ConfigValue, In, Bool
18# IDAES core
19from idaes.core import (
20 declare_process_block_class,
21 UnitModelBlockData,
22 useDefault,
23 StateBlock,
24)
25from idaes.core.util import scaling
26from idaes.core.util.config import is_physical_parameter_block
27from idaes.core.util.math import smooth_min, smooth_max
28from idaes.core.util.tables import create_stream_table_dataframe
29from idaes.core.solvers import get_solver
30from idaes.core.initialization import ModularInitializerBase
31from idaes.core.util.model_statistics import degrees_of_freedom
33# Logger
34import idaes.logger as idaeslog
36# Typing
37from typing import List
40__author__ = "Ahuora Centre for Smart Energy Systems, University of Waikato, New Zealand"
42# Set up logger
43_log = idaeslog.getLogger(__name__)
45class SimpleHeaderInitializer(ModularInitializerBase):
46 """Initialize a Header unit block with staged seeding and solves.
48 This routine performs a two-stage initialization:
49 1) Seed inlet and internal state variables, relax selected constraints, and
50 perform a first solve.
51 2) Reactivate/tighten constraints and perform a second solve.
53 Args:
54 blk: The Header unit model block to initialize.
55 **kwargs: Optional keyword arguments:
56 solver: A Pyomo/IDAES solver object. If not provided, uses ``get_solver()``.
57 solver_options (dict): Options to set on the solver, e.g. tolerances.
58 outlvl: IDAES log level (e.g., ``idaeslog.WARNING``).
60 Returns:
61 pyomo.opt.results.results_.SolverResults: The result object from the final solve.
63 Notes:
64 - Inlet state blocks are initialized via their own ``initialize`` if available.
65 - Mixed state is seeded from inlet totals/minimums (pressure) and average
66 enthalpy; works with temperature- or enthalpy-based property packages.
67 - Temporary seeds/relaxations are undone, leaving original DOF intact.
68 """
70 def initialize(self, blk, **kwargs):
71 # --- Solver setup
72 solver = kwargs.get("solver", None) or get_solver()
73 solver_options = kwargs.get("solver_options", {})
74 for k, v in solver_options.items():
75 solver.options[k] = v
77 outlvl = kwargs.get("outlvl", idaeslog.WARNING)
78 log = idaeslog.getLogger(__name__)
80 # --- Time index
81 t0 = blk.flowsheet().time.first()
83 # --- 1) Initialize inlet state blocks
84 inlet_blocks = list(blk.inlet_blocks)
85 if len(inlet_blocks) < 1:
86 raise ValueError("No inlet added to header.")
88 for sb in inlet_blocks:
89 if hasattr(sb, "initialize"):
90 sb.initialize(outlvl=outlvl)
92 # --- 2) Aggregate inlet info for seeding mixed state block
93 F_mixed = sum(
94 value(sb[t0].flow_mol)
95 for sb in inlet_blocks
96 )
97 P_mixed = min(
98 value(sb[t0].pressure)
99 for sb in inlet_blocks
100 )
101 E_mixed = sum(
102 value(sb[t0].flow_mol * sb[t0].enth_mol, 0.0)
103 for sb in inlet_blocks
104 )
105 if F_mixed > 0:
106 h_mixed = E_mixed / F_mixed
107 else:
108 # Seed from the first inlet’s enthalpy (no double subscripting)
109 first_inlet = inlet_blocks[0]
110 h_mixed = value(first_inlet[t0].enth_mol)
112 # --- 3) Seed mixed_state: flow, pressure, enthalpy
113 ms = blk.mixed_state
114 ms[t0].flow_mol.set_value(
115 F_mixed
116 )
117 ms[t0].pressure.set_value(
118 P_mixed
119 )
120 ms[t0].enth_mol.set_value(
121 h_mixed
122 )
123 ms.initialize(outlvl=outlvl)
125 # --- 4) Seed outlet_states with pressure, enthalpy
126 flow_undefined = []
127 defined_flow = 0
128 for sb in blk.outlet_blocks:
129 sb[t0].pressure.set_value(
130 value(ms[t0].pressure)
131 )
132 sb[t0].enth_mol.set_value(
133 value(ms[t0].enth_mol)
134 )
135 if sb in [blk.outlet_condensate_state, blk.outlet_vent_state]:
136 sb[t0].flow_mol.set_value(
137 0.0
138 )
139 else:
140 if value(sb[t0].flow_mol, exception=False) is None:
141 flow_undefined.append(sb)
142 else:
143 defined_flow += value(sb[t0].flow_mol)
145 tot_undefined_flow = max(sum(value(sb[t0].flow_mol) for sb in blk.inlet_blocks) - defined_flow, 0)
146 for sb in flow_undefined:
147 sb[t0].flow_mol.set_value(
148 tot_undefined_flow / len(flow_undefined)
149 )
151 for sb in blk.outlet_blocks:
152 sb.initialize(outlvl=outlvl)
154 res2 = solver.solve(blk, tee=False)
155 log.info(f"Header init status: {res2.solver.termination_condition}")
157 return res2
159def _make_config_block(config):
160 """Declare configuration options for the Header unit.
162 Declares property package references and integer counts for inlets and outlets.
164 Args:
165 config (ConfigBlock): The mutable configuration block to populate.
167 Raises:
168 ValueError: If invalid option values are provided by the caller (via IDAES).
169 """
171 config.declare(
172 "property_package",
173 ConfigValue(
174 default=useDefault,
175 domain=is_physical_parameter_block,
176 description="Property package to use for control volume",
177 ),
178 )
179 config.declare(
180 "property_package_args",
181 ConfigBlock(
182 implicit=True,
183 description="Arguments to use for constructing property packages",
184 ),
185 )
186 config.declare(
187 "num_inlets",
188 ConfigValue(
189 default=1,
190 domain=In(list(range(1, 100))),
191 description="Number of utility providers at inlets.",
192 ),
193 )
194 config.declare(
195 "num_outlets",
196 ConfigValue(
197 default=1,
198 domain=In(list(range(0, 100))),
199 description="Number of utility users at outlets." \
200 "Excludes outlets associated with condensate and vent flows.",
201 ),
202 )
203 config.declare(
204 "is_liquid_header",
205 ConfigValue(
206 default=False,
207 domain=Bool,
208 description="Flag for selecting liquid or vapour (including steam and other gases).",
209 ),
210 )
211@declare_process_block_class("simple_header")
212class SimpleHeaderData(UnitModelBlockData):
213 """Thermal utility header unit operation.
215 The Header aggregates multiple inlet providers and distributes utility to
216 multiple users, with optional venting, condensate removal (or liquid overflow), heat loss, and
217 pressure loss. A mixed (intermediate) state is used for balances and
218 pressure/enthalpy coupling across outlets.
220 Key features:
221 - Material, energy, and momentum balances with smooth min/max functions.
222 - Vapour/liquid equilibrium calculation for mixed state.
223 - Shared mixed enthalpy across outlets of the same phase.
224 - Computed excess flow from an overall flow balance.
225 - Optional heat and pressure losses.
227 Attributes:
228 inlet_list (list[str]): Names for inlet ports.
229 outlet_list (list[str]): Names for outlet ports (incl. condensate/ and vent).
230 inlet_blocks (list): StateBlocks for all inlets.
231 outlet_blocks (list): StateBlocks for all outlets.
232 mixed_state: Intermediate mixture StateBlock.
233 heat_loss (Var): Heat loss from the header (W).
234 pressure_loss (Var): Pressure drop from inlet minimum to mixed state (Pa).
235 makeup_flow_mol (Var): Required inlet makeup molar flow (mol/s).
236 """
238 default_initializer=SimpleHeaderInitializer
239 CONFIG = UnitModelBlockData.CONFIG()
240 _make_config_block(CONFIG)
242 def build(self) -> None:
243 # 1. Inherit standard UnitModelBlockData properties and functions
244 super().build()
246 # 2. Validate input parameters are valid
247 self._validate_model_config()
249 # 3. Create lists of ports with state blocks to add
250 self.inlet_list = self._create_inlet_port_name_list()
251 self.outlet_list = self._create_outlet_port_name_list()
253 # 4. Declare ports, state blocks and state property bounds
254 self.inlet_blocks = self._add_ports_with_state_blocks(
255 stream_list=self.inlet_list,
256 is_inlet=True,
257 has_phase_equilibrium=False,
258 is_defined_state=True,
259 )
260 self.outlet_blocks = self._add_ports_with_state_blocks(
261 stream_list=self.outlet_list,
262 is_inlet=False,
263 has_phase_equilibrium=False,
264 is_defined_state=False
265 )
266 self._internal_blocks = self._add_internal_state_blocks()
267 self._add_bounds_to_state_properties()
268 self._outlet_supply_blocks = self._create_custom_state_lists()
270 # 4. Declare references, variables and expressions for external and internal use
271 self._create_references()
272 self._create_variables()
273 self._create_expressions()
275 # 5. Set balance equations
276 self._add_material_balances()
277 self._add_energy_balances()
278 self._add_momentum_balances()
279 self._add_additional_constraints()
281 # 6. Other
282 self.scaling_factor = Suffix(direction=Suffix.EXPORT)
283 self.split_flow = self._create_flow_map_references()
285 def _validate_model_config(self) -> bool:
286 """Validate configuration for inlet and outlet counts.
288 Raises:
289 ValueError: If ``num_inlets < 1`` or ``num_outlets < 1``.
290 """
291 if self.config.num_inlets < 1:
292 raise ValueError("Header requires at least one provider (num_inlets >= 1).")
293 if self.config.num_outlets < 1:
294 raise ValueError("Header requires at least one user (num_outlets >= 1).")
295 return True
297 def _create_inlet_port_name_list(self) -> List[str]:
298 """Build ordered inlet port names.
300 Returns:
301 list[str]: Names ``["inlet_1", ..., "inlet_N"]`` based on ``num_inlets``.
302 """
303 return [
304 f"inlet_{i+1}" for i in range(self.config.num_inlets)
305 ]
307 def _create_outlet_port_name_list(self) -> List[str]:
308 """Build ordered outlet port names.
310 Returns:
311 list[str]: Names ``["outlet_1", ..., "outlet_n", "outlet_condensate", "outlet_vent"]``.
312 """
313 return [
314 f"outlet_{i+1}" for i in range(self.config.num_outlets)
315 ] + ["outlet_condensate"] + ["outlet_vent"]
317 def _add_ports_with_state_blocks(self,
318 stream_list: List[str],
319 is_inlet: List[str],
320 has_phase_equilibrium: bool=False,
321 is_defined_state: bool=None,
322 ) -> List[StateBlock]:
323 """Construct StateBlocks and expose them as ports.
325 Creates a StateBlock per named stream and attaches a corresponding inlet or
326 outlet Port. Inlet blocks are defined states; outlet blocks are calculated states.
328 Args:
329 stream_list (list[str]): Port/StateBlock base names to create.
330 is_inlet (bool): If True, create inlet ports with ``defined_state=True``;
331 otherwise create outlet ports with ``defined_state=False``.
332 has_phase_equilibrium (bool)
334 Returns:
335 list: The created StateBlocks, in the same order as ``stream_list``.
336 """
337 # Create empty list to hold StateBlocks for return
338 state_block_ls = []
340 # Setup StateBlock argument dict
341 tmp_dict = dict(**self.config.property_package_args)
342 tmp_dict["has_phase_equilibrium"] = has_phase_equilibrium
343 if is_defined_state == None:
344 tmp_dict["defined_state"] = True if is_inlet else False
345 else:
346 tmp_dict["defined_state"] = is_defined_state
348 # Create an instance of StateBlock for all streams
349 for s in stream_list:
350 sb = self.config.property_package.build_state_block(
351 self.flowsheet().time, doc=f"Thermophysical properties at {s}", **tmp_dict
352 )
353 setattr(
354 self, s + "_state",
355 sb
356 )
357 state_block_ls.append(sb)
358 add_fn = self.add_inlet_port if is_inlet else self.add_outlet_port
359 add_fn(
360 name=s,
361 block=sb,
362 )
364 return state_block_ls
366 def _add_internal_state_blocks(self) -> List[StateBlock]:
367 """Create the intermediate (mixed) StateBlock.
369 The mixed state:
370 - Has phase equilibrium enabled.
371 - Is not a defined state (solved from balances).
372 """
373 tmp_dict = dict(**self.config.property_package_args)
374 tmp_dict["has_phase_equilibrium"] = True
375 tmp_dict["defined_state"] = False
377 self.mixed_state = self.config.property_package.build_state_block(
378 self.flowsheet().time,
379 doc=f"Thermophysical properties at intermediate mixed state.",
380 **tmp_dict
381 )
382 return [
383 self.mixed_state
384 ]
386 def _add_bounds_to_state_properties(self) -> None:
387 """Add lower and/or upper bounds to state properties.
389 - Set nonnegativity lower bounds on all inlet/outlet molar flows.
390 """
391 for sb in (self.inlet_blocks + self.outlet_blocks):
392 for t in sb:
393 sb[t].flow_mol.setlb(0.0)
395 def _create_custom_state_lists(self) -> List[StateBlock]:
396 """Partition outlet names into vapour outlets and capture their StateBlocks.
398 Populates:
399 - ``_outlet_supply_list``: Outlet names excluding condensate and vent.
400 - ``_outlet_supply_blocks``: Corresponding StateBlocks.
401 """
402 self._outlet_supply_list = [
403 v for v in self.outlet_list
404 if not v in ["outlet_condensate", "outlet_vent"]
405 ]
406 return [
407 getattr(self, n + "_state")
408 for n in self._outlet_supply_list
409 ]
411 def _create_references(self) -> None:
412 """Create convenient References.
414 Creates references to mixed_state properties:
415 - ``total_flow_mol``
416 - ``total_flow_mass``
417 - ``pressure``
418 - ``temperature``
419 - ``enth_mol``
420 - ``enth_mass``
421 - ``vapor_frac``
422 """
423 self.total_flow_mol = Reference(
424 self.mixed_state[:].flow_mol
425 )
426 self.total_flow_mass = Reference(
427 self.mixed_state[:].flow_mass
428 )
429 self.pressure = Reference(
430 self.mixed_state[:].pressure
431 )
432 self.temperature = Reference(
433 self.mixed_state[:].temperature
434 )
435 self.enth_mol = Reference(
436 self.mixed_state[:].enth_mol
437 )
438 self.enth_mass = Reference(
439 self.mixed_state[:].enth_mass
440 )
441 self.vapor_frac = Reference(
442 self.mixed_state[:].vapor_frac
443 )
445 def _create_variables(self) -> None:
446 """Create required variables.
448 Creates:
449 - ``heat_loss`` (W)
450 - ``pressure_loss`` (Pa)
451 """
452 self.heat_loss = Var(
453 self.flowsheet().time,
454 initialize=0.0,
455 doc="Heat loss",
456 units=UNIT.W
457 )
458 self.pressure_loss = Var(
459 self.flowsheet().time,
460 initialize=0.0,
461 doc="Pressure loss",
462 units=UNIT.Pa
463 )
465 def _create_expressions(self) -> None:
466 """Create convenient Expressions.
468 Creates:
469 - ``balance_flow_mol`` (mol/s)
470 - ``degree_of_superheat`` (K)
471 - ``makeup_flow_mol`` (mol/s)
472 - ``_partial_total_flow_mol`` (mol/s): used for scaling purposes in a material balance
473 """
474 self.degree_of_superheat = Expression(
475 self.flowsheet().time,
476 rule=lambda b, t: b.temperature[t] - b.outlet_condensate_state[t].temperature
477 )
478 self._partial_total_flow_mol = Expression(
479 self.flowsheet().time,
480 rule=lambda b, t: (
481 sum(
482 o[t].flow_mol
483 for o in (b.inlet_blocks + b._outlet_supply_blocks)
484 )
485 )
486 )
487 self.balance_flow_mol = Expression(
488 self.flowsheet().time,
489 rule=lambda b, t: (
490 sum(
491 i[t].flow_mol
492 for i in b.inlet_blocks
493 )
494 -
495 sum(
496 o[t].flow_mol
497 for o in (
498 b._outlet_supply_blocks +
499 [
500 b.outlet_vent_state
501 if self.config.is_liquid_header
502 else b.outlet_condensate_state
503 ]
504 )
505 )
506 )
507 )
508 self.makeup_flow_mol = Expression(
509 self.flowsheet().time,
510 rule=lambda b, t: (
511 (
512 b.outlet_condensate_state[t].flow_mol
513 if self.config.is_liquid_header
514 else b.outlet_vent_state[t].flow_mol
515 )
516 -
517 b.balance_flow_mol[t]
518 )
519 )
521 def _add_material_balances(self) -> None:
522 """Material balance equations summary.
524 Introduces:
525 - ``_partial_total_flow_mol``: Sum of known inlet and vapour outlet flows,
526 used for scaling a smooth vent calculation.
528 Constraints:
529 - ``mixed_state_material_balance``: Mixed flow equals total inlet flow.
530 - ``vent_flow_balance``: Depends on the header's primary phase: liquid vs gas
531 If gas header, smoothly enforces nonnegative vent flow.
532 If liquid header, determines flow from mixed-state vapour fraction
533 - ``condensate_flow_balance``: Depends on the header's primary phase: liquid vs gas
534 If gas header, determines flow from mixed-state vapour fraction
535 If liquid header, smoothly enforces nonnegative condensate flow.
536 """
538 @self.Constraint(
539 self.flowsheet().time,
540 doc="Mixed state material balance",
541 )
542 def mixed_state_material_balance(b, t):
543 return (
544 b.mixed_state[t].flow_mol
545 ==
546 sum(
547 i[t].flow_mol
548 for i in b.inlet_blocks
549 )
550 )
552 eps = 1e-5 # smoothing parameter; smaller = closer to exact max, larger = smoother
553 if self.config.is_liquid_header:
554 # Assigns excess liquid flow to outlet_condensate
555 @self.Constraint(
556 self.flowsheet().time,
557 doc="Condensate flow balance." \
558 "Determines the positive amount of excess flow that exits through outlet_condensate"
559 )
560 def condensate_flow_balance(b, t):
561 return (
562 b.outlet_condensate_state[t].flow_mol
563 ==
564 smooth_max(
565 b.balance_flow_mol[t] / (b._partial_total_flow_mol[t] + 1e-6),
566 0.0,
567 eps,
568 ) * (b._partial_total_flow_mol[t] + 1e-6)
569 )
571 # Removes any gas/vapour from a liquid header
572 @self.Constraint(
573 self.flowsheet().time,
574 doc="Vent balance."
575 )
576 def vent_flow_balance(b, t):
577 return b.outlet_vent_state[t].flow_mol == (
578 b.mixed_state[t].flow_mol * b.mixed_state[t].vapor_frac
579 )
580 else:
581 # Assigns excess steam/vapour flow to outlet_vent
582 @self.Constraint(
583 self.flowsheet().time,
584 doc="Vent flow balance." \
585 "Determines the positive amount of excess flow that exits through the vent"
586 )
587 def vent_flow_balance(b, t):
588 return (
589 b.outlet_vent_state[t].flow_mol
590 ==
591 smooth_max(
592 b.balance_flow_mol[t] / (b._partial_total_flow_mol[t] + 1e-6),
593 0.0,
594 eps,
595 ) * (b._partial_total_flow_mol[t] + 1e-6)
596 )
598 # Removes any condensate/liquid from a steam/gas header
599 @self.Constraint(
600 self.flowsheet().time,
601 doc="Condensate balance."
602 )
603 def condensate_flow_balance(b, t):
604 return (
605 b.outlet_condensate_state[t].flow_mol
606 ==
607 b.mixed_state[t].flow_mol * (1 - b.mixed_state[t].vapor_frac)
608 )
610 def _add_energy_balances(self) -> None:
611 """Energy balance equations summary.
613 Introduces:
614 - ``_liq_out_enth_mol``: Shared molar enthalpy for all liquid outlets,
615 including the condensate.
616 - ``_vap_out_enth_mol``: Shared molar enthalpy for all vapour outlets,
617 including the vent.
619 Constraints:
620 - ``inlets_to_mixed_state_energy_balance``: Inlet energy to mixed state (+ heat loss).
621 - ``mixed_state_to_outlets_energy_balance``: Mixed state to all outlets.
622 - ``molar_enthalpy_equality_eqn``: Common vapour enthalpy across vapour outlets and vent.
623 """
624 @self.Constraint(self.flowsheet().time, doc="Inlets to mixed state energy balance including heat loss")
625 def inlets_to_mixed_state_energy_balance(b, t):
626 return (
627 b.mixed_state[t].flow_mol * b.mixed_state[t].enth_mol
628 + b.heat_loss[t]
629 ==
630 sum(
631 i[t].flow_mol * i[t].enth_mol
632 for i in b.inlet_blocks
633 )
634 )
635 @self.Constraint(
636 self.flowsheet().time,
637 doc="Mixed state to outlets energy balance"
638 )
639 def mixed_state_to_outlets_energy_balance(b, t):
640 return (
641 b.mixed_state[t].enth_mol
642 *
643 sum(
644 o[t].flow_mol
645 for o in b.outlet_blocks
646 )
647 ==
648 sum(
649 o[t].flow_mol * o[t].enth_mol
650 for o in b.outlet_blocks
651 )
652 )
653 if self.config.is_liquid_header:
654 self._liq_out_enth_mol = Var(
655 self.flowsheet().time,
656 initialize=42.0 * 18,
657 doc="Molar enthalpy of the liquid outlets",
658 units=UNIT.J / UNIT.mol
659 )
660 @self.Constraint(
661 self.flowsheet().time,
662 self._outlet_supply_blocks + [self.outlet_condensate_state], # exclude vent outlet
663 doc="All liquid outlets (incl. condensate) share a common liquid enthalpy",
664 )
665 def molar_enthalpy_equality_eqn(b, t, o):
666 return (
667 o[t].enth_mol
668 ==
669 b._liq_out_enth_mol[t]
670 )
671 else:
672 self._vap_out_enth_mol = Var(
673 self.flowsheet().time,
674 initialize=2700.0 * 18,
675 doc="Molar enthalpy of the vapour outlets",
676 units=UNIT.J / UNIT.mol
677 )
678 @self.Constraint(
679 self.flowsheet().time,
680 self._outlet_supply_blocks + [self.outlet_vent_state], # exclude condensate outlet
681 doc="All vapour outlets (incl. vent) share a common vapour enthalpy",
682 )
683 def molar_enthalpy_equality_eqn(b, t, o):
684 return (
685 o[t].enth_mol
686 ==
687 b._vap_out_enth_mol[t]
688 )
690 def _add_momentum_balances(self) -> None:
691 """Momentum balance equations summary.
693 Computes the minimum inlet pressure via a sequential smooth minimum and
694 sets the mixed-state pressure to that minimum minus ``pressure_loss``,
695 then enforces equality to every outlet pressure.
697 Notes:
698 - Uses IDAES ``smooth_min`` for differentiable minimum pressure.
699 - ``_eps_pressure`` is a smoothing parameter (units of pressure).
700 """
701 inlet_idx = RangeSet(len(self.inlet_blocks))
702 # Get units metadata
703 units = self.mixed_state.params.get_metadata()
704 # Add variables
705 self._minimum_pressure = Var(
706 self.flowsheet().time,
707 inlet_idx,
708 doc="Variable for calculating minimum inlet pressure",
709 units=units.get_derived_units("pressure"),
710 )
711 self._eps_pressure = Param(
712 mutable=True,
713 initialize=1e-3,
714 domain=PositiveReals,
715 doc="Smoothing term for minimum inlet pressure",
716 units=units.get_derived_units("pressure"),
717 )
718 # Calculate minimum inlet pressure
719 @self.Constraint(
720 self.flowsheet().time,
721 inlet_idx,
722 doc="Calculation for minimum inlet pressure",
723 )
724 def minimum_pressure_constraint(b, t, i):
725 if i == inlet_idx.first():
726 return (
727 b._minimum_pressure[t, i]
728 ==
729 (b.inlet_blocks[i - 1][t].pressure)
730 )
731 else:
732 return (
733 b._minimum_pressure[t, i]
734 ==
735 smooth_min(
736 b._minimum_pressure[t, i - 1],
737 b.inlet_blocks[i - 1][t].pressure,
738 b._eps_pressure,
739 )
740 )
741 # Set mixed pressure to minimum inlet pressure minus any pressure loss
742 @self.Constraint(
743 self.flowsheet().time,
744 doc="Pressure equality constraint from minimum inlet to mixed state",
745 )
746 def mixture_pressure(b, t):
747 return (
748 b.mixed_state[t].pressure
749 ==
750 b._minimum_pressure[t, inlet_idx.last()] - b.pressure_loss[t]
751 )
752 # Set outlet pressures to mixed pressure
753 @self.Constraint(
754 self.flowsheet().time,
755 self.outlet_blocks,
756 doc="Pressure equality constraint from mixed state to outlets",
757 )
758 def pressure_equality_eqn(b, t, o):
759 return (
760 b.mixed_state[t].pressure
761 ==
762 o[t].pressure
763 )
765 def _add_additional_constraints(self) -> None:
766 """Add auxiliary constraints and bounds.
768 - Fix vent vapour fraction to near one (near 100% vapour).
769 OR
770 - Fix condensate vapour fraction to a small value (near 100% liquid).
771 """
772 if self.config.is_liquid_header:
773 @self.Constraint(self.flowsheet().time, doc="Vent vapour fraction.")
774 def vent_vapour_fraction(b, t):
775 return (
776 b.outlet_vent_state[t].vapor_frac
777 ==
778 1 #1 - 1e-6
779 )
780 else:
781 @self.Constraint(self.flowsheet().time, doc="Condensate vapour fraction.")
782 def condensate_vapour_fraction(b, t):
783 return (
784 b.outlet_condensate_state[t].vapor_frac
785 ==
786 0 #1e-6
787 )
789 def _create_flow_map_references(self):
790 """Create a two-key Reference for outlet flows over time and outlet name.
792 Builds a mapping ``(t, outlet_name) -> outlet_state[t].flow_mol`` and exposes it
793 as a Reference for compact access to outlet flow splits.
795 Returns:
796 pyomo.core.base.reference.Reference: A Reference indexed by ``(time, outlet)``.
797 """
798 self.outlet_idx = pyo.Set(initialize=self.outlet_list)
799 # Map each (t, o) to the outlet state's flow var
800 ref_map = {}
801 for o in self.outlet_list:
802 if o != "vent":
803 outlet_state_block = getattr(self, f"{o}_state")
804 for t in self.flowsheet().time:
805 ref_map[(t, o)] = outlet_state_block[t].flow_mol
807 return Reference(ref_map)
809 def calculate_scaling_factors(self):
810 """Assign scaling factors to improve numerical conditioning.
812 Sets scaling factors for performance and auxiliary variables. If present,
813 also scales the shared vapour enthalpy variable ``_vap_out_enth_mol``.
814 """
815 super().calculate_scaling_factors()
816 scaling.set_scaling_factor(self.heat_loss, 1e-6)
817 scaling.set_scaling_factor(self.pressure_loss, 1e-6)
818 scaling.set_scaling_factor(self.balance_flow_mol, 1e-3)
819 scaling.set_scaling_factor(self._partial_total_flow_mol, 1e-3)
820 if hasattr(self, "_vap_out_enth_mol"):
821 scaling.set_scaling_factor(self._vap_out_enth_mol, 1e-6)
823 def _get_stream_table_contents(self, time_point=0):
824 """Create a stream table for all inlets and outlets.
826 Args:
827 time_point (int | float): Time index at which to extract stream data.
829 Returns:
830 pandas.DataFrame: A tabular view suitable for reporting via
831 ``create_stream_table_dataframe``.
832 """
833 io_dict = {}
835 for inlet_name in self.inlet_list:
836 io_dict[inlet_name] = getattr(self, inlet_name)
838 for outlet_name in self.outlet_list:
839 io_dict[outlet_name] = getattr(self, outlet_name)
841 return create_stream_table_dataframe(io_dict, time_point=time_point)
843 def _get_performance_contents(self, time_point=0, is_full_report=True):
844 """Collect performance results for reporting.
846 Args:
847 time_point (int | float): Time index at which to report values.
848 is_full_report (bool): Flag for full or partial performance report.
850 Returns:
851 dict: A report of internal unit model results.
852 """
853 return (
854 {
855 "vars": {
856 "Heat Loss": self.heat_loss[time_point],
857 "Pressure Drop": self.pressure_loss[time_point],
858 "Mass Flow": self.mixed_state[time_point].flow_mass,
859 "Molar Flow": self.mixed_state[time_point].flow_mol,
860 "Balance Flow": self.balance_flow_mol[time_point],
861 "Pressure": self.mixed_state[time_point].pressure,
862 "Temperature": self.mixed_state[time_point].temperature,
863 "Degree of Superheat": self.degree_of_superheat[time_point],
864 "Vapour Fraction": self.mixed_state[time_point].vapor_frac,
865 "Mass Specific Enthalpy": self.mixed_state[time_point].enth_mass,
866 "Molar Specific Enthalpy": self.mixed_state[time_point].enth_mol,
867 }
868 } if is_full_report else {
869 "vars": {
870 "Balance Flow": self.balance_flow_mol[time_point],
871 "Pressure": self.mixed_state[time_point].pressure,
872 "Temperature": self.mixed_state[time_point].temperature,
873 "Degree of Superheat": self.degree_of_superheat[time_point],
874 }
875 }
877 )
879 def initialize(self, *args, **kwargs):
880 """Initialize the Header unit using :class:`SimpleHeaderInitializer`.
882 Args:
883 *args: Forwarded to ``SimpleHeaderInitializer.initialize``.
884 **kwargs: Forwarded to ``SimpleHeaderInitializer.initialize`` (e.g., solver, options).
886 Returns:
887 pyomo.opt.results.results_.SolverResults: Results from the initializer's solve.
888 """
889 init = SimpleHeaderInitializer()
890 return init.initialize(self, *args, **kwargs)