"""Shared helpers for heat pump and refrigeration targeting."""
from __future__ import annotations
from typing import Any, Callable
import numpy as np
from CoolProp.CoolProp import PropsSI
from ....classes.stream_collection import StreamCollection
from ....classes.value import Value
from ....lib.config import tol
from ....lib.enums import PT
from ....lib.schemas.hpr import (
HeatPumpTargetInputs,
HPRBackendResult,
HPRParsedState,
HPRPeriodCase,
HPRThermoArtifacts,
MultiPeriodHPRTargetInputs,
SimulatedHPRAnnualizedCostAccounting,
)
from ....utils.blackbox_minimisers import multiminima
from ....utils.costing import (
compute_annual_capital_cost,
compute_annual_energy_cost,
compute_capital_cost,
)
from ...common.miscellaneous import (
g_ineq_penalty,
)
from ...common.problem_table_analysis import (
get_process_heat_cascade,
get_utility_heat_cascade,
)
from ._shared.ambient_preallocation import preallocate_direct_ambient_duties
from ._shared.plotting import plot_multi_hp_profiles_from_results
from ._shared.streams import (
get_ambient_air_stream,
get_carnot_hpr_cycle_streams,
get_Q_vals_at_T_hpr_from_bckgrd_profile,
)
__all__ = [
"PropsSI",
"solve_hpr_placement",
"solve_hpr_multiperiod_placement",
"compute_entropic_mean_temperature",
"calc_carnot_heat_pump_cop",
"calc_carnot_heat_engine_eta",
"get_carnot_hpr_cycle_streams",
"get_Q_vals_at_T_hpr_from_bckgrd_profile",
"validate_vapour_hp_refrigerant_ls",
"calc_hpr_obj",
"calc_simulated_hpr_annualized_costs",
"evaluate_carnot_hpr_result",
"evaluate_vapour_hpr_result",
"plot_multi_hp_profiles_from_results",
"get_process_heat_cascade",
"get_utility_heat_cascade",
"g_ineq_penalty",
"tol",
]
def _verify_x0_ls(x0_ls: list | float | None) -> np.ndarray | None:
if x0_ls is None:
return None
if isinstance(x0_ls, (int, float)):
x0_ls = [x0_ls]
x0_arr = np.asarray(x0_ls, dtype=np.float64)
if x0_arr.size == 0:
return None
if x0_arr.ndim == 0:
return x0_arr.reshape(1, 1)
if x0_arr.ndim == 1:
return x0_arr.reshape(1, -1)
return x0_arr
[docs]
def solve_hpr_placement(
f_obj: Callable,
x0_ls: list | float,
bnds: list,
args: HeatPumpTargetInputs,
) -> HPRBackendResult:
"""Run the configured multistart optimiser and post-process the best result."""
x0_arr = _verify_x0_ls(x0_ls)
try:
local_minima_x, local_minima_f = multiminima(
func=f_obj,
func_kwargs=args,
x0_ls=x0_arr,
bounds=bnds,
optimiser_handle=args.bb_minimiser,
opt_kwargs={"n_runs": max(1, int(args.max_multi_start))},
)
except Exception as exc:
if x0_arr is None or x0_arr.size == 0:
raise ValueError(
"Heat pump and refrigeration targeting "
f"({args.hpr_type}) failed during optimisation and has no "
"initial candidate fallback."
) from exc
local_minima_x, local_minima_f = np.asarray([]), np.asarray([])
candidate_x, candidate_f = _merge_candidate_points(
local_minima_x=local_minima_x,
local_minima_f=local_minima_f,
x0_arr=x0_arr,
f_obj=f_obj,
args=args,
)
if candidate_x.size == 0:
raise ValueError(
"Heat pump and refrigeration targeting "
f"({args.hpr_type}) failed to return any local minima."
)
for candidate_idx in np.argsort(candidate_f):
result = _evaluate_hpr_candidate(f_obj, candidate_x[candidate_idx], args)
if result.success and np.isfinite(float(result.obj)):
return result.with_updates(
success=True,
amb_streams=get_ambient_air_stream(
result.Q_amb_hot, result.Q_amb_cold, args
),
)
raise ValueError(
"Heat pump and refrigeration targeting "
f"({args.hpr_type}) failed to return an optimal result."
)
[docs]
def solve_hpr_multiperiod_placement(
f_obj: Callable,
x0_ls: list | float,
bnds: list,
args: MultiPeriodHPRTargetInputs,
) -> HPRBackendResult:
"""Optimise one shared HPR vector against all prepared period cases."""
x0_arr = _verify_x0_ls(x0_ls)
def multiperiod_objective(
x: np.ndarray,
mp_args: MultiPeriodHPRTargetInputs,
*,
debug: bool = False,
) -> HPRBackendResult:
return _compute_multiperiod_hpr_candidate(f_obj, x, mp_args, debug=debug)
try:
local_minima_x, local_minima_f = multiminima(
func=multiperiod_objective,
func_kwargs=args,
x0_ls=x0_arr,
bounds=bnds,
optimiser_handle=args.bb_minimiser,
opt_kwargs={"n_runs": max(1, int(args.max_multi_start))},
)
except Exception as exc:
if x0_arr is None or x0_arr.size == 0:
raise ValueError(
"Multi-period heat pump and refrigeration targeting "
f"({args.hpr_type}) failed during optimisation and has no "
"initial candidate fallback."
) from exc
local_minima_x, local_minima_f = np.asarray([]), np.asarray([])
candidate_x, candidate_f = _merge_candidate_points(
local_minima_x=local_minima_x,
local_minima_f=local_minima_f,
x0_arr=x0_arr,
f_obj=multiperiod_objective,
args=args,
)
if candidate_x.size == 0:
raise ValueError(
"Multi-period heat pump and refrigeration targeting "
f"({args.hpr_type}) failed to return any local minima."
)
selected_case = _selected_multiperiod_case(args)
for candidate_idx in np.argsort(candidate_f):
result = _evaluate_hpr_candidate(
multiperiod_objective,
candidate_x[candidate_idx],
args,
)
if result.success and np.isfinite(float(result.obj)):
return result.with_updates(
success=True,
amb_streams=get_ambient_air_stream(
result.Q_amb_hot,
result.Q_amb_cold,
selected_case.args,
),
)
raise ValueError(
"Multi-period heat pump and refrigeration targeting "
f"({args.hpr_type}) failed to return an optimal result."
)
def _selected_multiperiod_case(
args: MultiPeriodHPRTargetInputs,
) -> HPRPeriodCase:
for case in args.period_cases:
if str(case.period_id) == str(args.selected_period_id):
return case
for case in args.period_cases:
if int(case.period_idx) == int(args.selected_period_idx):
return case
raise ValueError(
"Selected period is not present in the prepared multi-period HPR cases."
)
def _compute_multiperiod_hpr_candidate(
f_obj: Callable,
x: np.ndarray,
args: MultiPeriodHPRTargetInputs,
*,
debug: bool = False,
) -> HPRBackendResult:
if not args.period_cases:
return HPRBackendResult.failure(reason="No period cases were prepared.")
period_outputs: dict[str, HPRBackendResult] = {}
for case in args.period_cases:
try:
result = f_obj(x, case.args, debug=debug)
except Exception as exc:
return HPRBackendResult.failure(
reason=f"HPR period {case.period_id!r} failed: {exc}"
)
if not isinstance(result, HPRBackendResult):
raise TypeError(
"Heat pump and refrigeration objective functions must return "
"HPRBackendResult."
)
if not result.success or not np.isfinite(float(result.obj)):
reason = result.failure_reason or "candidate failed"
return HPRBackendResult.failure(
reason=f"HPR period {case.period_id!r} failed: {reason}",
Q_amb_hot=result.Q_amb_hot,
Q_amb_cold=result.Q_amb_cold,
)
period_outputs[str(case.period_id)] = result
weights = np.asarray([case.weight for case in args.period_cases], dtype=float)
if not np.isfinite(weights).all() or float(weights.sum()) <= 0.0:
return HPRBackendResult.failure(reason="Period weights must be finite.")
selected = period_outputs[str(_selected_multiperiod_case(args).period_id)]
weighted = _weighted_hpr_backend_result(period_outputs, weights)
try:
shared_objective = _shared_hpr_candidate_objective(
list(period_outputs.values()),
weights,
)
except ValueError as exc:
return HPRBackendResult.failure(reason=str(exc))
weighted = weighted.with_updates(obj=shared_objective)
return selected.with_updates(
obj=shared_objective,
period_outputs=period_outputs,
weighted_output=weighted,
design_vector=np.asarray(x, dtype=float),
period_ids=[str(case.period_id) for case in args.period_cases],
period_weights=[float(case.weight) for case in args.period_cases],
)
def _weighted_hpr_backend_result(
period_outputs: dict[str, HPRBackendResult],
weights: np.ndarray,
) -> HPRBackendResult:
ordered = list(period_outputs.values())
first = ordered[0]
updates: dict[str, Any] = {}
for field in (
"obj",
"utility_tot",
"w_net",
"Q_ext_heat",
"Q_ext_cold",
"feasibility_penalty",
"Q_amb_hot",
"Q_amb_cold",
"w_hpr",
"w_he",
"heat_recovery",
"cop_h",
"eta_he",
"Q_cond",
"Q_evap",
"Q_cond_he",
"Q_evap_he",
"Q_heat",
"Q_cool",
):
weighted = _weighted_hpr_metric(ordered, field, weights)
if weighted is not None:
updates[field] = weighted
for field in ("hpr_operating_cost",):
weighted = _weighted_hpr_metric(ordered, field, weights)
if weighted is not None:
updates[field] = weighted
for field in (
"hpr_capital_cost",
"hpr_annualized_capital_cost",
"hpr_compressor_capital_cost",
"hpr_heat_exchanger_capital_cost",
):
maximum = _max_hpr_metric(ordered, field)
if maximum is not None:
updates[field] = maximum
total_cost = _total_hpr_annualized_cost(updates, ordered, weights)
if total_cost is not None:
updates["hpr_total_annualized_cost"] = total_cost
return first.with_updates(**updates)
def _weighted_hpr_metric(
results: list[HPRBackendResult],
field: str,
weights: np.ndarray,
) -> Any:
values = [getattr(result, field, None) for result in results]
if all(value is None for value in values):
return None
if any(value is None for value in values):
return None
return _aggregate_hpr_values(values, weights=weights, reducer="weighted")
def _max_hpr_metric(results: list[HPRBackendResult], field: str) -> Any:
values = [getattr(result, field, None) for result in results]
if all(value is None for value in values):
return None
if any(value is None for value in values):
return None
return _aggregate_hpr_values(values, weights=None, reducer="max")
def _aggregate_hpr_values(
values: list[Any],
*,
weights: np.ndarray | None,
reducer: str,
) -> Any:
first = values[0]
if isinstance(first, Value):
unit = first.unit
magnitudes = []
for value in values:
if not isinstance(value, Value):
return None
magnitudes.append(float(value.to(unit).value))
aggregate = (
float(np.average(magnitudes, weights=weights))
if reducer == "weighted"
else float(np.max(magnitudes))
)
return Value(aggregate, unit)
try:
arrays = [np.asarray(value, dtype=float) for value in values]
except TypeError, ValueError:
return None
shapes = {array.shape for array in arrays}
if len(shapes) != 1:
return None
stacked = np.stack(arrays, axis=0)
if reducer == "weighted":
aggregate = np.average(stacked, axis=0, weights=weights)
else:
aggregate = np.max(stacked, axis=0)
if aggregate.ndim == 0:
return float(aggregate)
return aggregate
def _total_hpr_annualized_cost(
updates: dict[str, Any],
results: list[HPRBackendResult],
weights: np.ndarray,
) -> Any:
operating = updates.get("hpr_operating_cost")
annualized_capital = updates.get("hpr_annualized_capital_cost")
if operating is not None and annualized_capital is not None:
try:
return operating + annualized_capital
except Exception:
pass
return _weighted_hpr_metric(results, "hpr_total_annualized_cost", weights)
def _shared_hpr_candidate_objective(
results: list[HPRBackendResult],
weights: np.ndarray,
) -> float:
has_cost_breakdown = any(
result.hpr_operating_cost is not None
or result.hpr_annualized_capital_cost is not None
for result in results
)
if not has_cost_breakdown:
fallback = _weighted_hpr_metric(results, "obj", weights)
return float(fallback)
if any(
result.hpr_operating_cost is None or result.hpr_annualized_capital_cost is None
for result in results
):
raise ValueError(
"Shared HPR candidates require a complete operating and annualized "
"capital cost breakdown for every period."
)
operating = _weighted_hpr_metric(results, "hpr_operating_cost", weights)
penalty = _weighted_hpr_metric(results, "feasibility_penalty", weights)
annualized_capital = _max_hpr_metric(
results,
"hpr_annualized_capital_cost",
)
return (
_annual_cost_magnitude(operating)
+ float(penalty)
+ _annual_cost_magnitude(annualized_capital)
)
def _annual_cost_magnitude(value: Any) -> float:
if isinstance(value, Value):
return float(value.to("$/y").value)
return float(value)
def _merge_candidate_points(
local_minima_x: list | np.ndarray,
local_minima_f: list | np.ndarray,
x0_arr: np.ndarray | None,
f_obj: Callable,
args: HeatPumpTargetInputs | MultiPeriodHPRTargetInputs,
) -> tuple[np.ndarray, np.ndarray]:
candidate_blocks = []
objective_blocks = []
local_minima_arr = _normalise_candidate_block(local_minima_x)
if local_minima_arr is not None:
candidate_blocks.append(local_minima_arr)
objective_blocks.append(np.asarray(local_minima_f, dtype=float))
x0_block = _normalise_candidate_block(x0_arr)
if x0_block is not None:
candidate_blocks.append(x0_block)
objective_blocks.append(
np.asarray(
[_score_hpr_candidate_objective(f_obj, x0, args) for x0 in x0_block],
dtype=float,
)
)
if not candidate_blocks:
return np.asarray([]), np.asarray([])
n_cols = x0_block.shape[1] if x0_block is not None else candidate_blocks[0].shape[1]
filtered_blocks = []
filtered_objectives = []
for block, objectives in zip(candidate_blocks, objective_blocks):
if block.shape[1] != n_cols:
continue
filtered_blocks.append(block)
filtered_objectives.append(objectives)
candidate_x = np.vstack(filtered_blocks)
candidate_f = np.concatenate(filtered_objectives)
_, unique_idx = np.unique(candidate_x, axis=0, return_index=True)
unique_idx = np.sort(unique_idx)
return candidate_x[unique_idx], candidate_f[unique_idx]
def _normalise_candidate_block(values: list | np.ndarray | None) -> np.ndarray | None:
if values is None:
return None
block = np.asarray(values, dtype=float)
if block.size == 0:
return None
if block.ndim == 0:
return block.reshape(1, 1)
if block.ndim == 1:
return block.reshape(1, -1)
return block.reshape(block.shape[0], -1)
def _score_hpr_candidate_objective(
f_obj: Callable,
x: np.ndarray,
args: HeatPumpTargetInputs | MultiPeriodHPRTargetInputs,
) -> float:
try:
result = f_obj(x, args, debug=False)
obj = float(result["obj"])
except Exception:
return 1e30
return obj if np.isfinite(obj) else 1e30
def _evaluate_hpr_candidate(
f_obj: Callable,
x: np.ndarray,
args: HeatPumpTargetInputs | MultiPeriodHPRTargetInputs,
) -> HPRBackendResult:
result = f_obj(x, args, debug=args.debug)
if not isinstance(result, HPRBackendResult):
raise TypeError(
"Heat pump and refrigeration objective functions must return "
"HPRBackendResult."
)
return result
def _build_hpr_accounting(
*,
work: float,
Q_ext_heat: float,
Q_ext_cold: float,
args: HeatPumpTargetInputs,
penalty_terms: np.ndarray | None = None,
penalise_external_cold_when_refrigerating: bool = False,
) -> tuple[float, float, float, float]:
"""Standardize external-utility, penalty, and objective semantics."""
positive_penalty_terms = (
np.maximum(penalty_terms, 0.0) if penalty_terms is not None else np.array([])
)
penalty = (
float(
g_ineq_penalty(
positive_penalty_terms,
eta=args.eta_penalty,
rho=args.rho_penalty,
form="square",
)
)
if positive_penalty_terms.size
else 0.0
)
if penalise_external_cold_when_refrigerating and not getattr(
args, "is_heat_pumping", True
):
penalty += float(
g_ineq_penalty(g=Q_ext_cold, rho=args.rho_penalty, form="square")
)
obj = float(
calc_hpr_obj(
work=work,
Q_ext_heat=Q_ext_heat,
Q_ext_cold=Q_ext_cold,
Q_hpr_target=args.Q_hpr_target,
heat_to_power_ratio=args.heat_to_power_ratio,
cold_to_power_ratio=args.cold_to_power_ratio,
penalty=penalty,
)
)
return Q_ext_heat, Q_ext_cold, penalty, obj
def _cycle_penalty(
*,
args: HeatPumpTargetInputs,
cycle_penalty_terms: list[float] | None = None,
) -> float:
if cycle_penalty_terms is None:
cycle_terms = np.array([])
else:
cycle_terms = np.maximum(np.asarray(cycle_penalty_terms, dtype=float), 0.0)
if not cycle_terms.size:
return 0.0
return float(
g_ineq_penalty(
cycle_terms,
eta=float(getattr(args, "eta_penalty", 0.01)),
rho=float(getattr(args, "rho_penalty", 10.0)),
form="square",
)
)
[docs]
def calc_simulated_hpr_annualized_costs(
*,
work: float,
work_arr: np.ndarray | None,
Q_ext_heat: float,
Q_ext_cold: float,
hpr_streams: StreamCollection,
hx_units: int,
penalty_power_equivalent: float,
args: HeatPumpTargetInputs,
) -> SimulatedHPRAnnualizedCostAccounting:
"""Return unit-aware annualized cost accounting for simulated HPR candidates."""
annual_hours = max(float(args.annual_op_time), 0.0)
ele_price = max(float(args.ele_price), 0.0)
heat_price = ele_price * max(float(args.heat_to_power_ratio), 0.0)
cold_price = ele_price * max(float(args.cold_to_power_ratio), 0.0)
operating_cost = (
compute_annual_energy_cost(work, ele_price, annual_hours)
+ compute_annual_energy_cost(Q_ext_heat, heat_price, annual_hours)
+ compute_annual_energy_cost(Q_ext_cold, cold_price, annual_hours)
).to("$/y")
work_values = np.array([]) if work_arr is None else work_arr
compressor_capital = compute_capital_cost(
work,
max(int(np.count_nonzero(work_values > tol)), 1),
args.hpr_comp_fixed_cost,
args.hpr_comp_variable_cost,
args.hpr_comp_cost_exp,
)
hpr_hot_streams = hpr_streams.get_hot_utility_streams()
hpr_cold_streams = hpr_streams.get_cold_utility_streams()
hx_duty = sum(
max(float(streams.sum_stream_attribute("heat_flow", idx=args.period_idx)), 0.0)
for streams in (hpr_hot_streams, hpr_cold_streams)
if len(streams) > 0
)
hx_capital = compute_capital_cost(
hx_duty,
hx_units,
args.hpr_hx_fixed_cost,
args.hpr_hx_variable_cost,
args.hpr_hx_cost_exp,
)
capital_cost = (compressor_capital + hx_capital).to("$")
annualized_capital = compute_annual_capital_cost(
capital_cost,
args.discount_rate,
args.serv_life,
)
total_annualized = (operating_cost + annualized_capital).to("$/y")
feasibility_penalty = compute_annual_energy_cost(
penalty_power_equivalent,
ele_price,
annual_hours,
)
return SimulatedHPRAnnualizedCostAccounting(
hpr_operating_cost=operating_cost,
hpr_capital_cost=capital_cost,
hpr_annualized_capital_cost=annualized_capital,
hpr_total_annualized_cost=total_annualized,
hpr_compressor_capital_cost=compressor_capital,
hpr_heat_exchanger_capital_cost=hx_capital,
feasibility_penalty=feasibility_penalty,
)
[docs]
def evaluate_carnot_hpr_result(
*,
args: HeatPumpTargetInputs,
state: HPRParsedState,
w_net: float,
w_hpr: float | list | np.ndarray,
Q_cond_total: np.ndarray,
Q_evap_total: np.ndarray,
w_he: float | list | np.ndarray | None = None,
heat_recovery: float | list | np.ndarray | None = None,
cop_h: float | list | np.ndarray | None = None,
eta_he: float | list | np.ndarray | None = None,
Q_cond: np.ndarray | None = None,
Q_evap: np.ndarray | None = None,
Q_cond_he: np.ndarray | None = None,
Q_evap_he: np.ndarray | None = None,
penalty_terms: np.ndarray | None = None,
debug: bool = False,
) -> HPRBackendResult:
"""Shared Carnot-family accounting, plotting, and result assembly."""
H_cold_with_amb = args.H_cold + args.z_amb_cold * state.Q_amb_cold
H_hot_with_amb = args.H_hot + args.z_amb_hot * state.Q_amb_hot
Q_ext_heat, Q_ext_cold, penalty, obj = _build_hpr_accounting(
work=float(w_net),
Q_ext_heat=np.abs(H_cold_with_amb[0]) - Q_cond_total.sum(),
Q_ext_cold=np.abs(H_hot_with_amb[-1]) - Q_evap_total.sum(),
args=args,
penalty_terms=penalty_terms,
penalise_external_cold_when_refrigerating=True,
)
hpr_streams = get_carnot_hpr_cycle_streams(
state.T_cond,
Q_cond_total,
state.T_evap,
Q_evap_total,
args,
)
debug_figure = None
if debug:
debug_figure = plot_multi_hp_profiles_from_results(
args.T_hot,
H_hot_with_amb,
args.T_cold,
H_cold_with_amb,
hpr_streams.get_hot_utility_streams(),
hpr_streams.get_cold_utility_streams(),
title=(
f"Obj {obj:.5f} = {(float(w_net) / args.Q_hpr_target):.5f} + "
f"{(Q_ext_heat / args.Q_hpr_target):.5f} + "
f"{(Q_ext_cold / args.Q_hpr_target):.5f} + "
f"{(penalty / args.Q_hpr_target):.5f}"
),
period_idx=args.period_idx,
)
return HPRBackendResult(
obj=obj,
utility_tot=float(w_net + Q_ext_heat + Q_ext_cold),
w_net=float(w_net),
w_hpr=w_hpr,
w_he=w_he,
heat_recovery=heat_recovery,
Q_ext_heat=Q_ext_heat,
Q_ext_cold=Q_ext_cold,
Q_amb_hot=state.Q_amb_hot,
Q_amb_cold=state.Q_amb_cold,
cop_h=cop_h,
eta_he=eta_he,
T_cond=state.T_cond,
T_evap=state.T_evap,
Q_cond=Q_cond_total if Q_cond is None else Q_cond,
Q_evap=Q_evap_total if Q_evap is None else Q_evap,
Q_cond_he=Q_cond_he,
Q_evap_he=Q_evap_he,
artifacts=HPRThermoArtifacts(
hpr_streams=hpr_streams, debug_figure=debug_figure
),
)
[docs]
def evaluate_vapour_hpr_result(
*,
args: HeatPumpTargetInputs,
state: HPRParsedState,
work: float,
work_arr: np.ndarray,
Q_heat: np.ndarray,
Q_cool: np.ndarray,
cop_h: float,
hpr_streams: StreamCollection,
model: Any = None,
penalty_terms: list[float] | None = None,
dT_subcool: np.ndarray | None = None,
dT_superheat: np.ndarray | None = None,
debug: bool = False,
) -> HPRBackendResult:
"""Shared simulated-vapour accounting, plotting, and result assembly."""
ambient_prealloc = preallocate_direct_ambient_duties(
args=args,
Q_amb_hot=state.Q_amb_hot,
Q_amb_cold=state.Q_amb_cold,
)
H_hot_with_amb = ambient_prealloc.H_hot_with_residual_ambient(args)
H_cold_with_amb = ambient_prealloc.H_cold_with_residual_ambient(args)
cond_hot_streams = hpr_streams.get_hot_utility_streams()
cond_cold_streams = ambient_prealloc.bckgrd_cold_streams + get_ambient_air_stream(
Q_amb_cold=ambient_prealloc.Q_amb_cold_residual, args=args
)
if len(cond_hot_streams) or len(cond_cold_streams):
pt_cond = get_process_heat_cascade(
hot_streams=cond_hot_streams,
cold_streams=cond_cold_streams,
is_shifted=True,
period_idx=args.period_idx,
)
Q_ext_heat = float(pt_cond[PT.H_NET][0])
cond_wrong_side = float(pt_cond[PT.H_NET][-1])
else:
Q_ext_heat = 0.0
cond_wrong_side = 0.0
evap_hot_streams = ambient_prealloc.bckgrd_hot_streams + get_ambient_air_stream(
Q_amb_hot=ambient_prealloc.Q_amb_hot_residual, args=args
)
evap_cold_streams = hpr_streams.get_cold_utility_streams()
if len(evap_hot_streams) or len(evap_cold_streams):
pt_evap = get_process_heat_cascade(
hot_streams=evap_hot_streams,
cold_streams=evap_cold_streams,
is_shifted=True,
period_idx=args.period_idx,
)
Q_ext_cold = float(pt_evap[PT.H_NET][-1])
evap_wrong_side = float(pt_evap[PT.H_NET][0])
else:
Q_ext_cold = 0.0
evap_wrong_side = 0.0
hx_units = (
len(cond_hot_streams)
+ len(cond_cold_streams)
+ len(evap_hot_streams)
+ len(evap_cold_streams)
)
all_penalty_terms = [*(penalty_terms or []), cond_wrong_side, evap_wrong_side]
penalty_power_equivalent = _cycle_penalty(
args=args,
cycle_penalty_terms=all_penalty_terms,
)
cost_accounting = calc_simulated_hpr_annualized_costs(
work=float(work),
work_arr=work_arr,
Q_ext_heat=Q_ext_heat,
Q_ext_cold=Q_ext_cold,
hpr_streams=hpr_streams,
hx_units=hx_units,
penalty_power_equivalent=penalty_power_equivalent,
args=args,
)
penalty = float(cost_accounting.feasibility_penalty.to("$/y").value)
obj = float(cost_accounting.hpr_total_annualized_cost.to("$/y").value) + penalty
debug_figure = None
if debug:
debug_figure = plot_multi_hp_profiles_from_results(
ambient_prealloc.T_hot_residual,
H_hot_with_amb,
ambient_prealloc.T_cold_residual,
H_cold_with_amb,
hpr_streams.get_hot_utility_streams(),
hpr_streams.get_cold_utility_streams(),
title=(
f"Obj {obj:.5f} $/y = "
f"{cost_accounting.hpr_total_annualized_cost.value:.5f} $/y + "
f"{penalty:.5f} $/y"
),
period_idx=args.period_idx,
)
return HPRBackendResult(
obj=obj,
utility_tot=float(work + Q_ext_heat + Q_ext_cold),
w_net=float(work),
w_hpr=work_arr,
Q_ext_heat=Q_ext_heat,
Q_ext_cold=Q_ext_cold,
hpr_operating_cost=cost_accounting.hpr_operating_cost,
hpr_capital_cost=cost_accounting.hpr_capital_cost,
hpr_annualized_capital_cost=cost_accounting.hpr_annualized_capital_cost,
hpr_total_annualized_cost=cost_accounting.hpr_total_annualized_cost,
hpr_compressor_capital_cost=cost_accounting.hpr_compressor_capital_cost,
hpr_heat_exchanger_capital_cost=(
cost_accounting.hpr_heat_exchanger_capital_cost
),
feasibility_penalty=penalty,
Q_amb_hot=state.Q_amb_hot,
Q_amb_cold=state.Q_amb_cold,
cop_h=float(cop_h),
T_cond=state.T_cond,
T_evap=state.T_evap,
dT_subcool=dT_subcool,
dT_superheat=dT_superheat,
Q_heat=Q_heat,
Q_cool=Q_cool,
artifacts=HPRThermoArtifacts(
hpr_streams=hpr_streams,
model=model,
debug_figure=debug_figure,
),
)
[docs]
def compute_entropic_mean_temperature(
T_arr: np.ndarray | list,
Q_arr: np.ndarray | list,
*,
input_T_units: str = "C",
) -> float:
"""Return the entropic mean temperature for a distributed heat load."""
T_arr = np.asarray(T_arr, dtype=float)
Q_arr = np.asarray(Q_arr, dtype=float)
unit_offset = 273.15 if input_T_units == "C" else 0
if T_arr.var() < tol:
return T_arr[0] + unit_offset
S_tot = (Q_arr / (T_arr + unit_offset)).sum()
return Q_arr.sum() / S_tot if S_tot > 0 else (T_arr.mean() + unit_offset)
[docs]
def calc_carnot_heat_pump_cop(
T_h: float | np.ndarray,
T_l: float | np.ndarray,
eta_ii: float,
) -> float | np.ndarray:
"""Compute a Carnot-based heating COP with a second-law efficiency factor."""
T_h_arr = np.asarray(T_h, dtype=float)
T_l_arr = np.asarray(T_l, dtype=float)
delta_T = T_h_arr - T_l_arr
with np.errstate(divide="ignore", invalid="ignore"):
cop = np.where(
delta_T > 0.0,
(T_l_arr / delta_T) * eta_ii + 1.0,
np.inf,
)
return cop.item() if cop.ndim == 0 else cop
[docs]
def calc_carnot_heat_engine_eta(
T_h: float | np.ndarray,
T_l: float | np.ndarray,
eta_ii: float,
) -> float | np.ndarray:
"""Compute a Carnot-based heat-engine efficiency with a second-law factor."""
T_h_arr = np.asarray(T_h, dtype=float)
T_l_arr = np.asarray(T_l, dtype=float)
eta = np.where(
T_h_arr > T_l_arr,
(1 - T_l_arr / T_h_arr) * eta_ii,
0.0,
)
return eta.item() if eta.ndim == 0 else eta
[docs]
def validate_vapour_hp_refrigerant_ls(
num_stages: int,
args: HeatPumpTargetInputs,
) -> list:
"""Return one refrigerant name per vapour-compression stage."""
if len(args.refrigerant_ls) > 0:
if args.do_refrigerant_sort:
refrigerants = [
ref
for ref, _ in sorted(
((ref, PropsSI("Tcrit", ref)) for ref in args.refrigerant_ls),
key=lambda x: x[1],
reverse=True,
)
]
else:
refrigerants = args.refrigerant_ls
if num_stages <= len(refrigerants):
return refrigerants[:num_stages]
padding = [refrigerants[-1]] * (num_stages - len(refrigerants))
return refrigerants + padding
return ["water" for _ in range(num_stages)]
[docs]
def calc_hpr_obj(
work: float,
Q_ext_heat: float,
Q_ext_cold: float,
Q_hpr_target: float,
heat_to_power_ratio: float = 1.0,
cold_to_power_ratio: float = 0.0,
penalty: float = 0.0,
) -> float:
"""Return the scalar screening objective used by HPR placement solvers."""
return (
work
+ (Q_ext_heat * heat_to_power_ratio)
+ (Q_ext_cold * cold_to_power_ratio)
+ penalty
) / Q_hpr_target