Skip to content
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
@@ -0,0 +1,112 @@
from libecalc.common.chart_type import ChartType
from libecalc.domain.component_validation_error import (
ProcessChartTypeValidationException,
)
from libecalc.domain.process.compressor.core.train.utils.enthalpy_calculations import (
calculate_enthalpy_change_head_iteration,
)
from libecalc.domain.process.process_solver.boundary import Boundary
from libecalc.domain.process.process_system.pressure_ratio_unit import PressureRatioUnit
from libecalc.domain.process.process_system.process_unit import ProcessUnitId
from libecalc.domain.process.value_objects.chart.chart import ChartData
from libecalc.domain.process.value_objects.chart.compressor import CompressorChart
from libecalc.domain.process.value_objects.fluid_stream import FluidService, FluidStream
from libecalc.domain.process.value_objects.fluid_stream.fluid import Fluid

_ALLOWED_CHART_TYPES = (ChartType.GENERIC_FROM_INPUT, ChartType.GENERIC_FROM_DESIGN_POINT)


class PressureRatioCompressor(PressureRatioUnit):
"""Compressor stage that operates at a given pressure ratio.

Provides the two methods OutletPressureSolverRatio needs:
- propagate_stream_at_pressure_ratio()
- get_recirculation_range_at_pressure_ratio()

propagate_stream() raises NotImplementedError — a pressure ratio is always
required to evaluate this unit.
"""

def __init__(
self,
process_unit_id: ProcessUnitId,
compressor_chart: ChartData,
fluid_service: FluidService,
) -> None:
if compressor_chart.origin_of_chart_data not in _ALLOWED_CHART_TYPES:
raise ProcessChartTypeValidationException(
f"PressureRatioCompressor requires a generic chart, "
f"got {compressor_chart.origin_of_chart_data.value}."
)
self._id = process_unit_id
self._compressor_chart = CompressorChart(compressor_chart)
self._fluid_service = fluid_service

def get_id(self) -> ProcessUnitId:
return self._id

@property
def compressor_chart(self) -> CompressorChart:
return self._compressor_chart

def get_recirculation_range_at_pressure_ratio(self, inlet_stream: FluidStream, pressure_ratio: float) -> Boundary:
"""Surge and stonewall boundaries at the given pressure ratio.

Uses head-based surge/stonewall limits — no shaft speed required.

Returns:
Boundary where:
min = additional rate (sm³/day) needed to reach the surge line
max = additional rate (sm³/day) available before stonewall
"""
outlet_pressure = inlet_stream.pressure_bara * pressure_ratio
enthalpy_change, polytropic_efficiency = calculate_enthalpy_change_head_iteration(
outlet_pressure=outlet_pressure,
polytropic_efficiency_vs_rate_and_head_function=self._compressor_chart.efficiency_as_function_of_rate_and_head,
inlet_streams=inlet_stream,
fluid_service=self._fluid_service,
)
head_joule_per_kg = float(enthalpy_change * polytropic_efficiency)

min_actual_rate = float(self._compressor_chart.minimum_rate_as_function_of_head(head_joule_per_kg))
max_actual_rate = float(self._compressor_chart.maximum_rate_as_function_of_head(head_joule_per_kg))

density = inlet_stream.density
min_sm3_per_day = self._fluid_service.mass_rate_to_standard_rate(
fluid_model=inlet_stream.fluid_model,
mass_rate_kg_per_h=min_actual_rate * density,
)
max_sm3_per_day = self._fluid_service.mass_rate_to_standard_rate(
fluid_model=inlet_stream.fluid_model,
mass_rate_kg_per_h=max_actual_rate * density,
)
inlet_sm3_per_day = inlet_stream.standard_rate_sm3_per_day
return Boundary(
min=max(0.0, min_sm3_per_day - inlet_sm3_per_day),
max=max(0.0, max_sm3_per_day - inlet_sm3_per_day),
)

def propagate_stream_at_pressure_ratio(self, inlet_stream: FluidStream, pressure_ratio: float) -> FluidStream:
Copy link
Copy Markdown
Contributor

@frodehk frodehk Apr 16, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggestion: Have a setter for pressure_ratio, and use propagate_stream() to fulfill the ProcessUnit contract. Other process units are already doing this, e.g. Choke (set_pressure_change).

Alternatively, split the ProcessUnit interface? e.g. RatioBasedUnit

"""Propagate stream to an outlet pressure defined by a pressure ratio.

The operating point is determined purely from inlet conditions and the
target pressure ratio — no shaft speed required.

Args:
inlet_stream: The inlet fluid stream.
pressure_ratio: Target outlet/inlet pressure ratio (e.g. 2.0 means double the pressure).

Returns:
The outlet FluidStream at target_pressure = inlet.pressure * pressure_ratio.
"""
outlet_pressure = inlet_stream.pressure_bara * pressure_ratio
enthalpy_change_joule_per_kg, polytropic_efficiency = calculate_enthalpy_change_head_iteration(
outlet_pressure=outlet_pressure,
polytropic_efficiency_vs_rate_and_head_function=self._compressor_chart.efficiency_as_function_of_rate_and_head,
inlet_streams=inlet_stream,
fluid_service=self._fluid_service,
)
target_enthalpy = float(inlet_stream.enthalpy_joule_per_kg + enthalpy_change_joule_per_kg)
props = self._fluid_service.flash_ph(inlet_stream.fluid_model, float(outlet_pressure), target_enthalpy)
outlet_fluid = Fluid(fluid_model=inlet_stream.fluid_model, properties=props)
return inlet_stream.with_new_fluid(outlet_fluid)
87 changes: 54 additions & 33 deletions src/libecalc/domain/process/process_solver/feasibility_solver.py
Original file line number Diff line number Diff line change
@@ -1,64 +1,62 @@
import abc

from libecalc.domain.process.entities.process_units.compressor import Compressor
from libecalc.domain.process.process_solver.float_constraint import FloatConstraint
from libecalc.domain.process.process_solver.outlet_pressure_solver import OutletPressureSolver
from libecalc.domain.process.process_solver.outlet_pressure_solver_ratio import OutletPressureSolverRatio
from libecalc.domain.process.process_solver.process_runner import ProcessRunner
from libecalc.domain.process.process_solver.process_system_solver import ProcessSystemSolver
from libecalc.domain.process.value_objects.fluid_stream import FluidStream


class FeasibilitySolver:
"""
Calculates how much of a given inlet rate exceeds what a compressor train
class FeasibilitySolver(abc.ABC):
"""Calculates how much of a given inlet rate exceeds what a compressor train
can handle for a target pressure.

Orchestrates OutletPressureSolver and queries compressor charts to find
the feasible rate — the excess is redirected via stream distribution.
"""

def __init__(
self,
outlet_pressure_solver: OutletPressureSolver,
compressors: list[Compressor],
runner: ProcessRunner,
):
self._solver = outlet_pressure_solver
self._compressors = compressors
self._runner = runner

@abc.abstractmethod
def get_excess_rate(
self,
inlet_stream: FluidStream,
target_pressure: FloatConstraint,
) -> float:
"""
Rate [sm³/day] that exceeds what this train can handle.
"""Rate [sm³/day] that exceeds what this train can handle.

This is the amount that must be redirected (e.g. via overflow)
to another train in the stream distribution.
"""
feasible = self._find_feasible_rate(inlet_stream, target_pressure)
excess_rate = max(0.0, inlet_stream.standard_rate_sm3_per_day - feasible)
return excess_rate
...


class RunnerBasedFeasibility(FeasibilitySolver):
"""Uses runner + stonewall limits to find the maximum feasible rate.

When the solver fails to meet the target pressure, applies the partial
solution and queries each compressor's stonewall limit to determine the
bottleneck rate.
"""

def __init__(
self,
solver: ProcessSystemSolver,
compressors: list[Compressor],
runner: ProcessRunner,
):
self._solver = solver
self._compressors = compressors
self._runner = runner

def _find_feasible_rate(
def get_excess_rate(
self,
inlet_stream: FluidStream,
target_pressure: FloatConstraint,
) -> float:
"""Highest standard rate [sm³/day] for which the train can meet target_pressure.

Returns the full inlet rate if the solver succeeds, otherwise finds the
bottleneck compressor's stone wall limit at the current operating point.
"""
solution = self._solver.find_solution(target_pressure, inlet_stream)

if solution.success:
# The train can handle the full rate — no need to search for a bottleneck.
return inlet_stream.standard_rate_sm3_per_day
return 0.0

# Apply the configuration before querying compressor charts.
self._runner.apply_configurations(solution.configuration)

# Search for compressor with the lowest max rate
min_max_rate = float("inf")
for compressor in self._compressors:
compressor_inlet = self._runner.run(
Expand All @@ -69,5 +67,28 @@ def _find_feasible_rate(
min_max_rate = min(min_max_rate, max_rate)

feasible_rate = max(0.0, min_max_rate)
return max(0.0, inlet_stream.standard_rate_sm3_per_day - feasible_rate)


class RatioBasedFeasibility(FeasibilitySolver):
"""Computes excess rate from stonewall limits at each stage.

Uses OutletPressureSolverRatio.get_max_standard_rate() to find the
bottleneck stage's stonewall limit, then returns the excess above that.
No runner needed — chart data lives on PressureRatioCompressor.
"""

def __init__(self, solver: OutletPressureSolverRatio):
self._solver = solver

def get_excess_rate(
self,
inlet_stream: FluidStream,
target_pressure: FloatConstraint,
) -> float:
solution = self._solver.find_solution(target_pressure, inlet_stream)
if solution.success:
return 0.0

return feasible_rate
max_rate = self._solver.get_max_standard_rate(target_pressure, inlet_stream)
return max(0.0, inlet_stream.standard_rate_sm3_per_day - max_rate)
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,9 @@
from libecalc.domain.component_validation_error import DomainValidationException
from libecalc.domain.process.process_solver.configuration import Configuration, OperatingConfiguration
from libecalc.domain.process.process_solver.float_constraint import FloatConstraint
from libecalc.domain.process.process_solver.outlet_pressure_solver import OutletPressureSolver
from libecalc.domain.process.process_solver.outlet_pressure_solver_speed import (
OutletPressureSolverSpeed,
)
from libecalc.domain.process.process_solver.pressure_control.downstream_choke import (
DownstreamChokePressureControlStrategy,
)
Expand All @@ -22,7 +24,7 @@ class MultiPressureSolver:
"""Tries to find the shaft speed satisfying N ordered pressure targets, one per segment. Segments share a
single physical shaft. Each segment's runner covers a disjoint sub-sequence of the stream propagation chain.

Independent OutletPressureSolver per segment (sequential) against its own pressure target. This gives the
Independent OutletPressureSolverSpeed per segment (sequential) against its own pressure target. This gives the
speed each segment would require if unconstrained. The segment requiring the highest speed is the binding
constraint. At binding speed, non-binding segments produce more pressure than needed and must have it reduced
by their pressure-control strategy.
Expand All @@ -35,7 +37,7 @@ class MultiPressureSolver:

def __init__(
self,
segments: list[OutletPressureSolver],
segments: list[OutletPressureSolverSpeed],
) -> None:
if len(segments) < 2:
raise DomainValidationException("MultiPressureSolver requires at least 2 segments.")
Expand All @@ -47,7 +49,7 @@ def __init__(
self._segments: Final = segments

@staticmethod
def _validate_pressure_control_placement(segments: list[OutletPressureSolver]) -> None:
def _validate_pressure_control_placement(segments: list[OutletPressureSolverSpeed]) -> None:
"""Upstream choke is only valid on the first segment; downstream choke only on the last."""
for i, segment in enumerate(segments):
strategy = segment.pressure_control_strategy
Expand Down
Loading
Loading