From e2c6323d890aff0fc57aec56ecab5028a2aa42d3 Mon Sep 17 00:00:00 2001 From: Claude Date: Wed, 17 Jun 2026 13:24:58 +0000 Subject: [PATCH 1/3] Fix surpyval imports for surpyval >= 0.11 compatibility non_repairable.py imported NonParametric, Parametric and ExactEventTime from submodule paths that were reorganised in surpyval 0.11. Import them from the top-level surpyval namespace instead, which works on both the old and new layouts, so RePyability can be used alongside the current surpyval that the Reliafy backend depends on. Co-Authored-By: Claude Opus 4.8 Claude-Session: https://claude.ai/code/session_013LT4iFs4HRQraP52HwVmgR --- repyability/non_repairable.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/repyability/non_repairable.py b/repyability/non_repairable.py index 358ca4d..9469d35 100644 --- a/repyability/non_repairable.py +++ b/repyability/non_repairable.py @@ -2,9 +2,7 @@ from scipy.integrate import quad from scipy.interpolate import interp1d from scipy.optimize import minimize -from surpyval.nonparametric import NonParametric -from surpyval.parametric import Parametric -from surpyval.parametric.exact_event_time import ExactEventTime +from surpyval import ExactEventTime, NonParametric, Parametric from repyability.rbd.standby_node import StandbyModel From c850c1c46c9a3d07edda4523e7d776d32988b5ec Mon Sep 17 00:00:00 2001 From: Claude Date: Wed, 17 Jun 2026 21:16:25 +0000 Subject: [PATCH 2/3] Upgrade to surpyval 0.11.1 (and NumPy 2) Bump the pinned dependencies to the latest surpyval (0.11.1), which requires numpy>=2.0,<3 and scipy>=1.13, and fix the incompatibilities that the newer surpyval/NumPy stack exposed: - NonRepairable._cost_rate: surpyval 0.11 returns 1-element arrays from sf/mean/qf, and scipy.optimize.minimize passes the time as a 1-element array. Coerce it to a float so the scipy.integrate.quad bound in avg_replacement_time() is a scalar (quad rejects array bounds). Fixes find_optimal_replacement() for parametric models. - NonRepairableRBD.random(): take the scalar from .random(1) before ordering the event queue and assigning into the output array (NumPy 2 rejects assigning a 1-element array to a scalar slot). Fixes the Monte-Carlo mean()/mean_time_to_failure(). - NonRepairableRBD.ff_by_node(): add the missing @check_x so it returns arrays like sf_by_node() rather than a scalar. Tests: relax the exact ff == 1 - sf assertion to approx (surpyval 0.11 computes ff(t) directly, so it can differ from 1 - sf(t) by a ULP), and add a system-level MTTF regression test. Full suite: 201 passed. Co-Authored-By: Claude Opus 4.8 Claude-Session: https://claude.ai/code/session_013LT4iFs4HRQraP52HwVmgR --- repyability/non_repairable.py | 5 +++++ repyability/rbd/non_repairable_rbd.py | 9 ++++++--- repyability/tests/test_rbd_misc.py | 18 +++++++++++++++++- repyability/tests/test_rbd_sf.py | 4 +++- requirements.txt | 8 ++++---- 5 files changed, 35 insertions(+), 9 deletions(-) diff --git a/repyability/non_repairable.py b/repyability/non_repairable.py index 9469d35..16b5588 100644 --- a/repyability/non_repairable.py +++ b/repyability/non_repairable.py @@ -61,6 +61,11 @@ def avg_replacement_time(self, t): return out def _cost_rate(self, t): + # surpyval >= 0.11 returns 1-element arrays from sf/mean/qf, and t can + # arrive as a 1-element array from scipy.optimize.minimize. Coerce to a + # plain float so the quadrature bound in avg_replacement_time() is a + # scalar (scipy.integrate.quad rejects array bounds). + t = float(np.asarray(t).item()) planned_costs = self.reliability_function(t) * self.cp unplanned_costs = (1 - self.reliability_function(t)) * self.cu avg_repl_time = self.avg_replacement_time(t) diff --git a/repyability/rbd/non_repairable_rbd.py b/repyability/rbd/non_repairable_rbd.py index 439d2f7..9883d7e 100644 --- a/repyability/rbd/non_repairable_rbd.py +++ b/repyability/rbd/non_repairable_rbd.py @@ -340,11 +340,10 @@ def sf_by_node( node_sf[node_name] = node.sf(x) return node_sf + @check_x def ff_by_node( self, x: Optional[ArrayLike] = None, *args, **kwargs ) -> Dict[Any, np.ndarray]: - - # The return dict node_ff: Dict[Any, np.ndarray] = {} # Cache the component reliabilities for efficiency @@ -437,7 +436,11 @@ def random(self, size): for i in range(size): event_queue: PriorityQueue = PriorityQueue() for node in self.G.nodes: - time = self.reliabilities[node].random(1) + # .random(1) returns a 1-element array; take the scalar so the + # event time orders the PriorityQueue and assigns into ``out`` + # (NumPy >= 2 rejects assigning a 1-element array to a scalar). + draw = np.asarray(self.reliabilities[node].random(1)) + time = float(draw.reshape(-1)[0]) event_queue.put(NodeFailure(time, node)) working_nodes = {k: True for k in self.G.nodes} diff --git a/repyability/tests/test_rbd_misc.py b/repyability/tests/test_rbd_misc.py index 72d92b7..7f456e2 100644 --- a/repyability/tests/test_rbd_misc.py +++ b/repyability/tests/test_rbd_misc.py @@ -3,12 +3,28 @@ Uses pytest fixtures located in conftest.py in the tests/ directory. """ +import numpy as np import pytest -from surpyval import FixedEventProbability +from surpyval import Exponential, FixedEventProbability from repyability.rbd.non_repairable_rbd import NonRepairableRBD +def test_rbd_mean_time_to_failure_series(): + # Two exponential components in series fail at the first failure, so the + # system MTTF is 1 / (lambda_1 + lambda_2). Exercises the Monte-Carlo + # random()/mean() path. + rbd = NonRepairableRBD( + [("input", "a"), ("a", "b"), ("b", "output")], + { + "a": Exponential.from_params([0.01]), + "b": Exponential.from_params([0.02]), + }, + ) + np.random.seed(0) + assert rbd.mean(mc_samples=20_000) == pytest.approx(1 / 0.03, rel=5e-2) + + # Check components are correct lengths def test_rbd_components(rbd1: NonRepairableRBD, rbd2: NonRepairableRBD): assert len(rbd1.nodes) == 3 diff --git a/repyability/tests/test_rbd_sf.py b/repyability/tests/test_rbd_sf.py index 76b35db..c898ec9 100644 --- a/repyability/tests/test_rbd_sf.py +++ b/repyability/tests/test_rbd_sf.py @@ -15,7 +15,9 @@ def test_rbd_sf_ff_by_node(rbd_series: NonRepairableRBD): node_sf = rbd_series.sf_by_node(t) node_ff = rbd_series.ff_by_node(t) for k in node_sf.keys(): - assert node_ff[k] == 1 - node_sf[k] + # surpyval >= 0.11 computes ff(t) directly (more precisely than + # 1 - sf(t)), so the two can differ by a floating-point ULP. + assert node_ff[k] == pytest.approx(1 - node_sf[k]) def test_rel_unrel(rbd_series: NonRepairableRBD): diff --git a/requirements.txt b/requirements.txt index 6ef5fd1..6209e49 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,6 +1,6 @@ # Production requirements # For development, use requirements_dev.txt instead -networkx==3.0 -numpy==1.23.5 -surpyval==0.10.10 -tqdm==4.64.1 +networkx>=3.0 +numpy>=2.0,<3 +surpyval==0.11.1 +tqdm>=4.64 From 3f5d2cac910f03a488e9b2213c400b6509608815 Mon Sep 17 00:00:00 2001 From: Claude Date: Wed, 17 Jun 2026 23:28:58 +0000 Subject: [PATCH 3/3] Add conditional survival cs() across the RBD package MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Implement cs(x, X) = sf(X + x) / sf(X) — the probability of surviving a further x given survival to X — as a first-class method: - a shared conditional_survival(model, x, X) helper (utils.wrappers) that works for anything exposing sf, clipping to [0, 1] and returning 0 where the item has all but surely failed by X; - NonRepairableRBD.cs (passes through sf() args like working/broken nodes and method); - StandbyModel, RepeatedNode, RepeatedStandbyNode; - PerfectReliability (always 1) and PerfectUnreliability (always 0). This matches SurPyval's own model.cs definition, so conditional survival can be computed uniformly across systems and every node type. Full suite: 201 passing plus a new cs test. Co-Authored-By: Claude Opus 4.8 Claude-Session: https://claude.ai/code/session_013LT4iFs4HRQraP52HwVmgR --- repyability/rbd/helper_classes.py | 10 +++++++ repyability/rbd/non_repairable_rbd.py | 24 +++++++++++++++ repyability/rbd/repeated_node.py | 6 ++++ repyability/rbd/repeated_standby_node.py | 6 ++++ repyability/rbd/standby_node.py | 6 ++++ repyability/tests/test_rbd_sf.py | 20 +++++++++++++ repyability/utils/wrappers.py | 38 ++++++++++++++++++++++++ 7 files changed, 110 insertions(+) diff --git a/repyability/rbd/helper_classes.py b/repyability/rbd/helper_classes.py index 5741c1b..d5cd23b 100644 --- a/repyability/rbd/helper_classes.py +++ b/repyability/rbd/helper_classes.py @@ -10,6 +10,11 @@ def sf(cls, x): def ff(cls, x): return np.zeros_like(x).astype(float) + @classmethod + def cs(cls, x, X): + # Always working -> conditional survival is always 1. + return np.ones_like(np.atleast_1d(x)).astype(float) + @classmethod def random(cls, size): return np.ones(size) * np.inf @@ -24,6 +29,11 @@ def sf(cls, x): def ff(cls, x): return np.ones_like(x).astype(float) + @classmethod + def cs(cls, x, X): + # Never working -> conditional survival is always 0. + return np.zeros_like(np.atleast_1d(x)).astype(float) + @classmethod def random(self, size): return np.zeros(size) diff --git a/repyability/rbd/non_repairable_rbd.py b/repyability/rbd/non_repairable_rbd.py index 9883d7e..bb4491e 100644 --- a/repyability/rbd/non_repairable_rbd.py +++ b/repyability/rbd/non_repairable_rbd.py @@ -10,6 +10,8 @@ from numpy.typing import ArrayLike from surpyval import NonParametric +from repyability.utils.wrappers import conditional_survival + from .helper_classes import PerfectReliability, PerfectUnreliability from .rbd import RBD from .repeated_node import RepeatedNode @@ -327,6 +329,28 @@ def reliability(self, x: Optional[ArrayLike] = None, *args, **kwargs): """ return self.sf(x, *args, **kwargs) + def cs(self, x: ArrayLike, X: ArrayLike, *args, **kwargs) -> np.ndarray: + """Returns the conditional survival of the system. + + That is, the probability the system survives a *further* ``x`` given it + has already survived to ``X``: ``R(x | X) = sf(X + x) / sf(X)``. + + Parameters + ---------- + x : ArrayLike + The further duration/s at which conditional survival is evaluated. + X : ArrayLike + The age/s the system is known to have survived to. + *args, **kwargs : + Any sf() arguments (e.g. working_nodes, broken_nodes, method). + + Returns + ------- + np.ndarray + The conditional survival probability/ies. + """ + return conditional_survival(self, x, X, *args, **kwargs) + @check_x def sf_by_node( self, x: Optional[ArrayLike] = None, *args, **kwargs diff --git a/repyability/rbd/repeated_node.py b/repyability/rbd/repeated_node.py index 14ff116..74a8c8c 100644 --- a/repyability/rbd/repeated_node.py +++ b/repyability/rbd/repeated_node.py @@ -44,3 +44,9 @@ def ff(self, x): else: ff = self.model.ff(x) ** self.repeats return ff + + def cs(self, x, X): + """Conditional survival ``R(x | X) = sf(X + x) / sf(X)``.""" + from repyability.utils.wrappers import conditional_survival + + return conditional_survival(self, x, X) diff --git a/repyability/rbd/repeated_standby_node.py b/repyability/rbd/repeated_standby_node.py index 4b38825..ce833ac 100644 --- a/repyability/rbd/repeated_standby_node.py +++ b/repyability/rbd/repeated_standby_node.py @@ -60,3 +60,9 @@ def sf(self, *args, **kwargs): def ff(self, *args, **kwargs): return self._sf_model.ff(*args, **kwargs) + + def cs(self, x, X): + """Conditional survival ``R(x | X) = sf(X + x) / sf(X)``.""" + from repyability.utils.wrappers import conditional_survival + + return conditional_survival(self, x, X) diff --git a/repyability/rbd/standby_node.py b/repyability/rbd/standby_node.py index f96eb18..613ba1f 100644 --- a/repyability/rbd/standby_node.py +++ b/repyability/rbd/standby_node.py @@ -179,3 +179,9 @@ def ff(self, *args, **kwargs): if self._sf_model is not None: return self._sf_model.ff(*args, **kwargs) return self.model.ff(*args, **kwargs) + + def cs(self, x, X): + """Conditional survival ``R(x | X) = sf(X + x) / sf(X)``.""" + from repyability.utils.wrappers import conditional_survival + + return conditional_survival(self, x, X) diff --git a/repyability/tests/test_rbd_sf.py b/repyability/tests/test_rbd_sf.py index c898ec9..bcfb77a 100644 --- a/repyability/tests/test_rbd_sf.py +++ b/repyability/tests/test_rbd_sf.py @@ -146,3 +146,23 @@ def test_rbd_sf_NonRepairableRBD_as_node(rbd1: NonRepairableRBD, method: str): * rbd.reliabilities["NonRepairableRBD"].ff(t) * rbd.reliabilities[4].ff(t) ) == rbd.ff(t, method=method) + + +def test_rbd_conditional_survival(): + """cs(x, X) returns sf(X + x) / sf(X), clipped to [0, 1].""" + import numpy as np + from surpyval import Exponential, Weibull + + rbd = NonRepairableRBD( + [("input", "c1"), ("c1", "c2"), ("c2", "output")], + { + "c1": Weibull.from_params([100, 2]), + "c2": Exponential.from_params([0.01]), + }, + ) + x = np.array([10.0, 20.0, 50.0]) + X = 30.0 + expected = rbd.sf(x + X) / rbd.sf(X) + assert np.allclose(rbd.cs(x, X), expected) + # Surviving zero further time is certain. + assert rbd.cs(0.0, X)[0] == 1.0 diff --git a/repyability/utils/wrappers.py b/repyability/utils/wrappers.py index cdf4c54..21024c2 100644 --- a/repyability/utils/wrappers.py +++ b/repyability/utils/wrappers.py @@ -1,3 +1,6 @@ +import numpy as np + + def check_probability(func): """checks probability is between 0 and 1""" @@ -10,3 +13,38 @@ def wrap(obj, target: float, *args, **kwargs): return func(obj, target, *args, **kwargs) return wrap + + +def conditional_survival(model, x, X, *args, **kwargs): + """Conditional survival from any model exposing ``sf``. + + Returns the probability of surviving a *further* ``x`` given the item has + already survived to ``X``: + + .. math:: + R(x \\mid X) = \\frac{R(X + x)}{R(X)} + + Parameters + ---------- + model : object + Anything with an ``sf(x, ...)`` method (a distribution, a standby + arrangement, an RBD, ...). + x : array_like or scalar + The further duration(s) at which conditional survival is evaluated. + X : array_like or scalar + The age(s) the item is known to have survived to. + + Returns + ------- + numpy.ndarray + The conditional survival probability, clipped to ``[0, 1]``; where the + item has all but surely failed by ``X`` (``R(X) ≈ 0``) it is ``0``. + """ + x = np.atleast_1d(np.asarray(x, dtype=float)) + X = np.atleast_1d(np.asarray(X, dtype=float)) + denom = np.asarray(model.sf(X, *args, **kwargs), dtype=float) + numer = np.asarray(model.sf(x + X, *args, **kwargs), dtype=float) + with np.errstate(divide="ignore", invalid="ignore"): + out = numer / denom + out = np.where(np.isfinite(out), out, 0.0) + return np.clip(out, 0.0, 1.0)