From 2825b172c7bb3d53421d88f3c39b1bac9493ba72 Mon Sep 17 00:00:00 2001
From: Milan <milan@MacBook-Air-von-Milan.local>
Date: Tue, 11 Nov 2025 16:59:34 +0100
Subject: [PATCH 1/2] adding tritium module, migrated from core pathsim

---
 .DS_Store                                 | Bin 0 -> 6148 bytes
 src/.DS_Store                             | Bin 0 -> 6148 bytes
 src/pathsim_chem/.DS_Store                | Bin 0 -> 6148 bytes
 src/pathsim_chem/__init__.py              |   3 +
 src/pathsim_chem/tritium/.DS_Store        | Bin 0 -> 6148 bytes
 src/pathsim_chem/tritium/README.md        |   3 +
 src/pathsim_chem/tritium/__init__.py      |   4 +
 src/pathsim_chem/tritium/bubbler.py       | 218 ++++++++++++++++++++++
 src/pathsim_chem/tritium/residencetime.py | 123 ++++++++++++
 src/pathsim_chem/tritium/splitter.py      |  50 +++++
 src/pathsim_chem/tritium/tcap.py          |  27 +++
 11 files changed, 428 insertions(+)
 create mode 100644 .DS_Store
 create mode 100644 src/.DS_Store
 create mode 100644 src/pathsim_chem/.DS_Store
 create mode 100644 src/pathsim_chem/tritium/.DS_Store
 create mode 100644 src/pathsim_chem/tritium/README.md
 create mode 100644 src/pathsim_chem/tritium/__init__.py
 create mode 100644 src/pathsim_chem/tritium/bubbler.py
 create mode 100644 src/pathsim_chem/tritium/residencetime.py
 create mode 100644 src/pathsim_chem/tritium/splitter.py
 create mode 100644 src/pathsim_chem/tritium/tcap.py

diff --git a/.DS_Store b/.DS_Store
new file mode 100644
index 0000000000000000000000000000000000000000..9a874b5768f336915163bb88cd434575b859f936
GIT binary patch
literal 6148
zcmeH~Jr2S!425ml0g0s}V-^m;4I%_5-~tF3k&vj^b9A16778<}(6eNJu~Vz<8=6`~
zboab&MFtUB!i}=AFfm2m$tVxGT*u4pe81nU<hDp!0Y1oNKeq`gKn17(6`%rCU`7h$
zLB3qg=$ZH^RDcT1Ljn6f6u7Y_ThKop2tERUJ<4ua`z!$#D}XiGf~deWTES>lA49C}
z?O@64YO)2RT{MRe%{!}2F))pG(Sih~)xkgosK7*lF7m<7{{#Hn{6A@7N(HFEpDCdI
z{<z=arSfe3cs;A1vTExF2mNw{x1Rtcb`-DSZrCrj0Bf=ZQGxMCz-3^d0zXyY1+O9z
AbpQYW

literal 0
HcmV?d00001

diff --git a/src/.DS_Store b/src/.DS_Store
new file mode 100644
index 0000000000000000000000000000000000000000..7b69a2719d14a4df45ec0a69cc0f8d3b45ee132f
GIT binary patch
literal 6148
zcmeHK%}&BV5S~TU5;WnUiN{U6k_bc%#7l$t0@mn34Yrg>W9brF<UmMz))(?gd>&_Z
z7o$Or8Z|S??6=I$?9aE&4g&z9Jqhao6##Hh2@54`z7QHGU6O*etRV_LM+y-Pp$k12
zEP1oxH!?tXw*@2cA;UZScz=Bu2GevD#413amtYh}X;!Pfio#lPePh$vbjr@HH&9bA
z^Rr3X_Q%)s>q4nGSoHniGU`vem7P<SWPX(N$2uSi`{;6W6(ymXw$&sFQyuG>0jK1Y
zx|O}z?67rImGx$0UX`;}qgIv2jplq_a(4F*PR_dz@i0+OhE##ypO$Tl3wXxP%JYS3
zFAi|O6ui>u(HkeRO778TvB+X_MrME+U<Upu1NJ#{%71kKylrNH8TdX1=zWl=gs#Qh
zpnf{A(Io(45zSh#PqhTaNQ<t;+#rsi2os8ELWM0cgb7Exw0^F|+@J{uVT%u8R~EKH
z5qfo;Uut&{u0d{@0cPMc162E=R4D(St*-yI6F1BNGw?qd5XFwyY2(>!Ze2Ph<ywh)
rhe|?zxxx1k?C4aCzLbjhP_>|4QU}qsm>Wb33jYz%G;qTV{3ru&Q509O

literal 0
HcmV?d00001

diff --git a/src/pathsim_chem/.DS_Store b/src/pathsim_chem/.DS_Store
new file mode 100644
index 0000000000000000000000000000000000000000..5b72ad7f5cb26099d8b7c6fb06514451f9c2d4b0
GIT binary patch
literal 6148
zcmeHK&1%~~5T141cugockdhn)y&CFZ42JfgI8H9PRi!;N!Lcl3i^^JK$@Re)bWKh<
z^w@J>p^wl92;_M>v%3!d>&=8>W?=Rk&Cjm(liBqEfM|?^Q-BfxNL0dH5r+*zand!(
zSx*5`=s7|d!@!>;LqA%}X2<`?0DZdxq%eXI$A#~we2o0|hf$cMmCBRI?d4zX7i2*m
z$Q!S(CSKZ3$4R3-`btlql#2Yh-S$6+-ASu-bgtsG9md^}E(n7zrd)ps<3LRsY8(fN
zuJuen7G<$jI-X9em3rB!*JfpBTCb5=t<Gje`R47r_a9r|qCu?g43h$XbXs0poWngD
zE1L?EPUPcDQprjv>YaKWH5i)e^TYLwfoeAWNHZdTUUq&xe)w~?DJyu+0?2+`ZOj_>
z7gu|rMOWVEyElp>756YQ-ymhNFayj0Gq61d?CBy8wl^iu#0)S4JHh~+4-P7!YcVsZ
zj}9Dk3xJqMw-&Uim!KNQqH8fTh!GTFN)b(|uq}o#<>+^8o@+5PXv#s@=0n((h3!y;
zULEImbUFyvAdk!dGqB9So_Q?M{eSWA{(rfMpO^t=V5b-m`KH%w;Fj#&y0$sGYaP^E
rR1(U|3_h2jp|4`hrK|V^RSWtZG7w#hnL+fR@E-w90}sr=i!$&RLdbJc

literal 0
HcmV?d00001

diff --git a/src/pathsim_chem/__init__.py b/src/pathsim_chem/__init__.py
index 8a72535..1f4f7c6 100644
--- a/src/pathsim_chem/__init__.py
+++ b/src/pathsim_chem/__init__.py
@@ -11,3 +11,6 @@
     __version__ = "unknown"
 
 __all__ = ["__version__"]
+
+#for direct block import from main package
+from .tritium import *
diff --git a/src/pathsim_chem/tritium/.DS_Store b/src/pathsim_chem/tritium/.DS_Store
new file mode 100644
index 0000000000000000000000000000000000000000..5008ddfcf53c02e82d7eee2e57c38e5672ef89f6
GIT binary patch
literal 6148
zcmeH~Jr2S!425mzP>H1@V-^m;4Wg<&0T*E43hX&L&p$$qDprKhvt+--jT7}7np#A3
zem<@ulZcFPQ@L2!n>{z**<q8>++&mCkOWA81W14cNZ<zv;LbK1Poaz?KmsK2CSc!(
z0ynLxE!0092;Krf2c+FF_Fe*7ECH>lEfg7;MkzE(HCqgga^y>{tEnwC%0;vJ&^%eQ
zLs35+`xjp>T0<F0fCPF1$Cyrb|F7^5{eNG?83~ZUUlGt@xh*qZDeu<Z%US-OSsOPv
j)R!Z4KLME7ReXlK;d!wEw5GODWMKRea10D2@KpjYNUI8I

literal 0
HcmV?d00001

diff --git a/src/pathsim_chem/tritium/README.md b/src/pathsim_chem/tritium/README.md
new file mode 100644
index 0000000..4835a1f
--- /dev/null
+++ b/src/pathsim_chem/tritium/README.md
@@ -0,0 +1,3 @@
+# Tritium Toolbox
+
+This module contains blocks for nuclear fusion modeling. Especially for tritium cycle modeling.
\ No newline at end of file
diff --git a/src/pathsim_chem/tritium/__init__.py b/src/pathsim_chem/tritium/__init__.py
new file mode 100644
index 0000000..4cae89e
--- /dev/null
+++ b/src/pathsim_chem/tritium/__init__.py
@@ -0,0 +1,4 @@
+from .residencetime import *
+from .splitter import *
+from .bubbler import *
+# from .tcap import *
\ No newline at end of file
diff --git a/src/pathsim_chem/tritium/bubbler.py b/src/pathsim_chem/tritium/bubbler.py
new file mode 100644
index 0000000..f272858
--- /dev/null
+++ b/src/pathsim_chem/tritium/bubbler.py
@@ -0,0 +1,218 @@
+#########################################################################################
+##
+##                                      Bubbler Block
+##                                (blocks/fusion/bubbler.py)
+##
+#########################################################################################
+
+# IMPORTS ===============================================================================
+
+import numpy as np
+
+from ..dynsys import DynamicalSystem
+from ...events.schedule import ScheduleList
+
+
+# BLOCK DEFIINITIONS ====================================================================
+
+class Bubbler4(DynamicalSystem):
+    """
+    Tritium bubbling system with sequential vial collection stages.
+
+    This block models a tritium collection system used in fusion reactor blanket 
+    purge gas processing. The system bubbles tritium-containing gas through a series 
+    of liquid-filled vials to capture and concentrate tritium for measurement and 
+    inventory tracking.
+
+
+    Physical Description
+    --------------------
+    The bubbler consists of two parallel processing chains:
+
+    **Soluble Chain (Vials 1-2):** 
+    Tritium already in soluble forms (HTO, HT) flows sequentially through 
+    vials 1 and 2. Each vial has a collection efficiency :math:`\\eta_{vial}`, 
+    representing the fraction of tritium that dissolves into the liquid phase 
+    and is retained.
+
+    **Insoluble Chain (Vials 3-4):** 
+    Tritium in insoluble forms (T₂, organically bound) first undergoes catalytic 
+    conversion to soluble forms with efficiency :math:`\\alpha_{conv}`. The 
+    converted tritium, along with uncaptured soluble tritium from the first chain, 
+    then flows through vials 3 and 4 with the same collection efficiency.
+
+
+    Mathematical Formulation
+    -------------------------
+    The system is governed by the following differential equations for the 
+    vial inventories :math:`x_i`:
+
+    .. math::
+
+        \\frac{dx_1}{dt} &= \\eta_{vial} \\cdot u_{sol}
+
+        \\frac{dx_2}{dt} &= \\eta_{vial} \\cdot (1-\\eta_{vial}) \\cdot u_{sol}
+
+        \\frac{dx_3}{dt} &= \\eta_{vial} \\cdot [\\alpha_{conv} \\cdot u_{insol} + (1-\\eta_{vial})^2 \\cdot u_{sol}]
+
+        \\frac{dx_4}{dt} &= \\eta_{vial} \\cdot (1-\\eta_{vial}) \\cdot [\\alpha_{conv} \\cdot u_{insol} + (1-\\eta_{vial})^2 \\cdot u_{sol}]
+
+
+    The sample output represents uncaptured tritium exiting the system:
+
+    .. math::
+    
+        y_{sample} = (1-\\alpha_{conv}) \\cdot u_{insol} + (1-\\eta_{vial})^2 \\cdot [\\alpha_{conv} \\cdot u_{insol} + (1-\\eta_{vial})^2 \\cdot u_{sol}]
+
+
+    Where:
+        - :math:`u_{sol}` = soluble tritium input flow rate
+        - :math:`u_{insol}` = insoluble tritium input flow rate  
+        - :math:`\\eta_{vial}` = vial collection efficiency
+        - :math:`\\alpha_{conv}` = conversion efficiency from insoluble to soluble
+        - :math:`x_i` = tritium inventory in vial i
+
+    Parameters
+    ----------
+    conversion_efficiency : float 
+       Conversion efficiency from insoluble to soluble forms (:math:`\\alpha_{conv}`), 
+       between 0 and 1.
+    vial_efficiency : float 
+       Collection efficiency of each vial (:math:`\\eta_{vial}`), between 0 and 1.
+    replacement_times : float | list[float] | list[list[float]] 
+        Times at which each vial is replaced with a fresh one. If None, no 
+        replacement events are created. If a single value is provided, it is 
+        used for all vials. If a single list of floats is provided, it will be 
+        used for all vials. If a list of lists is provided, each sublist 
+        corresponds to the replacement times for each vial.
+
+    Notes
+    -----
+    Vial replacement is modeled as instantaneous reset events that set the 
+    corresponding vial inventory to zero, simulating the physical replacement 
+    of a full vial with an empty one.
+    """
+
+    _port_map_out = {
+        "vial1": 0,
+        "vial2": 1,
+        "vial3": 2,
+        "vial4": 3,
+        "sample_out": 4,
+    }
+    _port_map_in = {
+        "sample_in_soluble": 0,
+        "sample_in_insoluble": 1,
+    }
+
+    def __init__(
+        self,
+        conversion_efficiency=0.9,
+        vial_efficiency=0.9,
+        replacement_times=None,
+        ):
+
+        #bubbler parameters
+        self.replacement_times = replacement_times
+        self.vial_efficiency = vial_efficiency
+        self.conversion_efficiency = conversion_efficiency
+
+        #dynamical component, ode rhs
+        def _fn_d(x, u, t):
+
+            #short
+            ve = self.vial_efficiency
+            ce = self.conversion_efficiency
+
+            #unpack inputs
+            sol, ins = u
+
+            #compute vial content change rates
+            dv1 = ve * sol
+            dv2 = dv1 * (1 - ve)
+            dv3 = ve * (ce * ins + (1 - ve)**2 * sol)
+            dv4 = dv3 * (1 - ve)
+
+            return np.array([dv1, dv2, dv3, dv4])
+
+        #algebraic output component
+        def _fn_a(x, u, t):
+
+            #short
+            ve = self.vial_efficiency
+            ce = self.conversion_efficiency
+
+            #unpack inputs
+            sol, ins = u
+
+            sample_out = (1 - ce) * ins + (1 - ve)**2 * (ce * ins + (1 - ve)**2 * sol)
+
+            return np.hstack([x, sample_out])
+
+        #initialization just like `DynamicalSystem` block
+        super().__init__(func_dyn=_fn_d, func_alg=_fn_a, initial_value=np.zeros(4))
+
+        #create internal vial reset events
+        self._create_reset_events()
+
+
+    def _create_reset_event_vial(self, i, reset_times):
+        """Define event action function and return a `ScheduleList` event 
+        per vial `i` that triggers at predefined `reset_times`. 
+        """
+
+        def reset_vial_i(_):
+            #get the full engine state
+            x = self.engine.get()
+            #set index 'i' to zero
+            x[i] = 0.0
+            #set the full engine state 
+            self.engine.set(x)
+
+        return ScheduleList(
+            times_evt=reset_times, 
+            func_act=reset_vial_i
+            )
+
+
+    def _create_reset_events(self):
+        """Create reset events for all vials based on the replacement times.
+
+        Raises
+        ------
+        ValueError : If reset_times is not valid.
+
+        Returns
+        -------
+        events : list[ScheduleList]
+            list of reset events for vials
+        """
+
+        replacement_times = self.replacement_times
+        self.events = []
+
+        # if reset_times is a single list use it for all vials
+        if replacement_times is None:
+            return 
+
+        if isinstance(replacement_times, (int, float)):
+            replacement_times = [replacement_times]
+
+        # if it's a flat list use it for all vials
+        elif isinstance(replacement_times, list) and all(
+            isinstance(t, (int, float)) for t in replacement_times
+            ):
+            replacement_times = [replacement_times] * 4
+
+        elif isinstance(replacement_times, np.ndarray) and replacement_times.ndim == 1:
+            replacement_times = [replacement_times.tolist()] * 4
+
+        elif isinstance(replacement_times, list) and len(replacement_times) != 4:
+            raise ValueError(
+                "replacement_times must be a single value or a list with the same length as the number of vials"
+            )
+
+        #create the internal events
+        self.events = [
+            self._create_reset_event_vial(i, ts) for i, ts in enumerate(replacement_times)
+            ]
diff --git a/src/pathsim_chem/tritium/residencetime.py b/src/pathsim_chem/tritium/residencetime.py
new file mode 100644
index 0000000..814c654
--- /dev/null
+++ b/src/pathsim_chem/tritium/residencetime.py
@@ -0,0 +1,123 @@
+#########################################################################################
+##
+##                           Blocks for residence time modeling
+##                            (blocks/fusion/residencetime.py)
+##
+#########################################################################################
+
+# IMPORTS ===============================================================================
+
+import numpy as np
+
+from ..dynsys import DynamicalSystem
+
+
+# BLOCKS ================================================================================
+
+class ResidenceTime(DynamicalSystem):
+    """Chemical process block with residence time model.
+
+    This block implements an internal 1st order linear ode with 
+    multiple inputs, outputs, an internal constant source term 
+    and no direct passthrough.
+
+    The internal ODE with inputs :math:`u_i` :
+
+    .. math::
+        
+        \\dot{x} = - x / \\tau + \\mathrm{src} + \\sum_i \\beta_i u_i
+
+
+    And the output equation for every output `i` :
+
+    .. math::
+    
+        y_i = \\gamma_i x
+
+
+    Parameters
+    ----------
+    tau : float
+        residence time, inverse natural frequency (eigenvalue)
+    betas: None | list[float] | np.ndarray[float]
+        weights of inputs that are accumulated in state, optional
+    gammas : None | list[float] | np.ndarray[float]
+        weights of states (fractions) for output, optional
+    initial_value : float
+        initial value of state / initial quantity of process
+    source_term : float
+        constant source term / generation term of the process
+    """
+
+    def __init__(self, tau=1, betas=None, gammas=None, initial_value=0, source_term=0):
+
+        #input validation
+        if np.isclose(tau, 0):
+            raise ValueError(f"'tau' must be nonzero but is {tau}")
+    
+        #time constant and input/output weights
+        self.tau = tau
+        self.betas = 1 if betas is None else np.array(betas)
+        self.gammas = 1 if gammas is None else np.array(gammas)
+        self.source_term = source_term
+
+        #rhs of residence time ode
+        def _fn_d(x, u, t):
+            return -x/self.tau + self.source_term + sum(self.betas*u)
+
+        #jacobian of rhs wrt x
+        def _jc_d(x, u, t):
+            return -1/self.tau
+
+        #output function of residence time ode
+        def _fn_a(x, u, t):
+            return self.gammas * x
+
+        #initialization just like `DynamicalSystem` block
+        super().__init__(func_dyn=_fn_d, jac_dyn=_jc_d, func_alg=_fn_a, initial_value=initial_value)
+
+
+class Process(ResidenceTime):
+    """Simplified version of the `ResidenceTime` model block
+    with all inputs being summed equally and only the state 
+    and the flux being returned to the output
+
+    This block implements an internal 1st order linear ode with 
+    multiple inputs, outputs and no direct passthrough.
+
+    The internal ODE with inputs :math:`u_i` :
+
+    .. math::
+        
+        \\dot{x} = - x / \\tau + \\mathrm{src} + \\sum_i u_i
+
+
+    And the output equations for output `i=0` and `i=1`:
+
+    .. math::
+    
+        y_0 = x
+
+    .. math::
+
+        y_1 = x / \\tau
+
+
+    Parameters
+    ----------
+    tau : float
+        residence time, inverse natural frequency (eigenvalue)
+    initial_value : float
+        initial value of state / initial quantity of process
+    source_term : float
+        constant source term / generation term of the process
+    """
+
+    #max number of ports
+    _n_out_max = 2
+
+    #maps for input and output port labels
+    _port_map_out = {"x": 0, "x/tau": 1}
+
+    def __init__(self, tau=1, initial_value=0, source_term=0):
+        super().__init__(tau, 1, [1, 1/tau], initial_value, source_term)
diff --git a/src/pathsim_chem/tritium/splitter.py b/src/pathsim_chem/tritium/splitter.py
new file mode 100644
index 0000000..e034244
--- /dev/null
+++ b/src/pathsim_chem/tritium/splitter.py
@@ -0,0 +1,50 @@
+#########################################################################################
+##
+##                                 Flow Splitter Block 
+##                             (blocks/fusion/splitter.py)
+##
+#########################################################################################
+
+# IMPORTS ===============================================================================
+
+import numpy as np
+
+from ..function import Function
+
+
+# BLOCKS ================================================================================
+
+class Splitter(Function):
+    """Splitter block that splits the input signal into multiple 
+    outputs weighted with the specified fractions.
+
+    Note
+    ----
+    The output fractions must sum to one.
+    
+    Parameters
+    ----------
+    fractions : np.ndarray | list
+        fractions to split the input signal into, 
+        must sum up to one
+    """
+
+    #max number of ports
+    _n_in_max = 1
+
+    #maps for input and output port labels
+    _port_map_in = {"in": 0}
+    
+    def __init__(self, fractions=None):
+
+        self.fractions = np.ones(1) if fractions is None else np.array(fractions)
+
+        #input validation
+        if not np.isclose(sum(self.fractions), 1):
+            raise ValueError(f"'fractions' must sum to one and not {sum(self.fractions)}")
+
+        #initialize like `Function` block
+        super().__init__(func=lambda u: self.fractions*u)
+
+        #dynamically define output port map based on fractions
+        self._port_map_out = {f"out {fr}": i for i, fr in enumerate(self.fractions)}
diff --git a/src/pathsim_chem/tritium/tcap.py b/src/pathsim_chem/tritium/tcap.py
new file mode 100644
index 0000000..8473127
--- /dev/null
+++ b/src/pathsim_chem/tritium/tcap.py
@@ -0,0 +1,27 @@
+#########################################################################################
+##
+##            Blocks for Thermal Cycle Absorption Process (TCAP) modelling
+##                              (blocks/fusion/tcap.py)
+##
+#########################################################################################
+
+# IMPORTS ===============================================================================
+
+from ..ode import ODE
+
+
+# BLOCKS ================================================================================
+
+class TCAP1D(ODE):
+    """This block models the Thermal Cycle Absorption Process (TCAP) in 1d. 
+
+    The model uses a 1d finite difference spatial discretization to construct 
+    a nonlinear ODE internally as proposed in 
+
+        https://doi.org/10.1016/j.ijhydene.2023.03.101
+
+
+    """
+    raise NotImplementedError("TCAP1D block is currently not impolemented!")
+
+

From 1cd5a6fe1207bed998fb6916bbc76f810a504e67 Mon Sep 17 00:00:00 2001
From: Milan <milan@MacBook-Air-von-Milan.local>
Date: Tue, 11 Nov 2025 17:08:03 +0100
Subject: [PATCH 2/2] migrated tests for tritium module

---
 .DS_Store                           | Bin 6148 -> 8196 bytes
 tests/.DS_Store                     | Bin 0 -> 6148 bytes
 tests/tritium/__init__.py           |   0
 tests/tritium/test_bubbler.py       | 174 ++++++++++++++++++++++++++++
 tests/tritium/test_residencetime.py | 115 ++++++++++++++++++
 tests/tritium/test_splitter.py      |  61 ++++++++++
 6 files changed, 350 insertions(+)
 create mode 100644 tests/.DS_Store
 create mode 100644 tests/tritium/__init__.py
 create mode 100644 tests/tritium/test_bubbler.py
 create mode 100644 tests/tritium/test_residencetime.py
 create mode 100644 tests/tritium/test_splitter.py

diff --git a/.DS_Store b/.DS_Store
index 9a874b5768f336915163bb88cd434575b859f936..9d750714b2a0bbd873c8c526d56189c6ca7e55eb 100644
GIT binary patch
literal 8196
zcmeHM&2Jk;6rblL&8$P|q@<8q1#OAcR>TD-q9JJ!xQ<gA6;xE`i#A~PV@OtA?^-)+
zyCRA(s<;56h%0RoTzjDhB-9f(BF+_}{{R;*{0H!6XM#-A04GGnJZa`P@6CI&@BO^q
z@vI3UP-*!Sgp3kGC^m^~2)hS}oR@8f<g~>?1dJzM=*Dq~Lo%?#*%%B1h5^HXVZbn8
z82C>xAUA3hy34+=W=(GxFbw=J8IaEh8=J(k#;(Hb(Scod1VHRZb0#=ua)4-SENkp4
z#2FM&p@=FJ^brG8INELPSJv27sKN>8!w2XqgFc}^c6D6e&fE!<6`I~KU>MlRfIPeR
z5kV5-`J(syo)V7clnA#Ynr;}o==iPM3<HstN<TubuYX|AARVMb^v&9O+^MBb+7cBf
zxg@m(d=21kj~#a@;GH@fnT?~=2_l|k0fEoaWpz1-{J2wzTahoa*qWN47PabZbYtWB
znG;2OVtR73Xm89+mWuY&<n-pIMGuTUcj`j@a@dIC4^>^DpABU1eRY0>`T*9VE83;d
zaBZUI#f@f9a8FPu+M0%3xJ|=P+ca4Dp~Bpm{MkC&-^N<cH`5@F#1c;gVx8?1$jQ7H
zN3SPb&M$j$*lI+g&W1fF^je{l^3!2h3*O~oCFNIAIl3s~B%RG1>x{JozrfKiX!9z1
zi<JPElWrqD>B=Ek_a52%=sq^Q|FMw+PmLWrJ~?&j<rS=}!J||xhK|d_<+jS~l8;Nx
zCLDWzLEvJm((qdB2k!)Exr`spLgX~~ThjchBhvFt9&!Kf<4=qpEaVS8x!ev!;D-DS
zM|eE)gJ^BB6Zq-+A_np+)Eb5&8`R#1pQaY`*??A$97XjU+oRONc#%@}4a;R}zEGO5
zD05h!R>x1E3QEb(oU|+^Sie>ey@V>sltQ60i_2^=O2^92C|<H6tg&)Q(?rQZ@+|IV
zXUQtw4oz~Ed`xbTo8&gRL%t(h<Sw~KekQ+=U&(Jkfk6RmI08rE7>vVln1UBz8cxEC
zFb8Mg6_|%}umFp&1S@b6Jn(_T8U&DlfE3#B9(3V-xC$S@CvXiuh0ow~xDGeq3-}Vg
zf?M!4+*Wzzwr|~eou>x4gM*&Cun||Q`pZbKUC(d*;Cz4jeh>d|Jeuo&b)UCuNfgG>
z8m=mHvNe`9b`@snVrTY$<im(z;9+E-PbEDn&;N6~zyE(2OEWne1`GrLkOAng)~XeB
z&UW<&1j$|A#`ZQgDbla2Fbl!%Pd@}u|Kkrs9NXwBW6B!43ekdMUIb(`n9eZpuQKpE
D>IW0c

delta 146
zcmZp1XfcprU|?W$DortDU=RQ@Ie-{MGpJ516gFW5GME{P8HyN^Q;L&wlJfI&Ktf=^
z1f@X&KrFxjmETx6jeTN+&|-ED4nby!8Ubz~?Fz)I8w<ZPPv)0VWMrIdDx=K<Gy`NI
N$PEme<9X&V0{|pP70&<w

diff --git a/tests/.DS_Store b/tests/.DS_Store
new file mode 100644
index 0000000000000000000000000000000000000000..f53d21c29503cb987e8bd11377318f20bad82408
GIT binary patch
literal 6148
zcmeHK%}T>S5Z-O8-BN@e6g>vK3bquX;3dZT0!H+pQX5m+V45vWY7eE5tG<Nq;KhRv
z;2U`JU7Xq7iv1Bih)9`%*>84sW;gpK>}D8a+!^>4#!8Ga0g9L@K(jz_9Cc1I#*+u+
zYL0*nSl1iIJujS#CdYqd0N-7nC9KZ^Xg|F_>A}d8*9(I<DVJY)CYxJX&5OJ!ic7a6
zhi=kL2C>uZpJVSt%FvtGP46^l4I8EPLm4H_AZqnhLg2R`<l-!dd^vRFAo63C>*@hf
z5QRo*b2O@!YdcnLZ@goTYI|r_tK)G&Y;5oDA2)8oZY1wClN>%eC2IyJa1X}Jru_Ii
z^xz>eWo9UBSKMpa?M<@hBf7F)_d>;pczDY^K07Z5i=zA=MD;q~t#<47qfka|SYfik
zPme|>F+dCu1HZt4IbB5Y7v`cVhyh|?xfsCn!2(6J4W=5^(E*Jr0RVI0RsuHG5}0Ed
zv<;>jVFZNhR6w1|O^LyEI`}OUXB$j4>U73U@xjf?+*BxBtq${BDx7g!Belc;F)+<Q
zR(DHy{~vu_|4%1Tj~E~ZmWlzMtGjgvwq)+sxy|8SYk^*YqF`LA@goHcbrnM_Ud3CW
aO2BWS0caacHG&6(egqT^)DQ!I%D_8olXNox

literal 0
HcmV?d00001

diff --git a/tests/tritium/__init__.py b/tests/tritium/__init__.py
new file mode 100644
index 0000000..e69de29
diff --git a/tests/tritium/test_bubbler.py b/tests/tritium/test_bubbler.py
new file mode 100644
index 0000000..8cb5163
--- /dev/null
+++ b/tests/tritium/test_bubbler.py
@@ -0,0 +1,174 @@
+########################################################################################
+##
+##                                  TESTS FOR 
+##                             'tritium.bubbler.py'
+##
+########################################################################################
+
+# IMPORTS ==============================================================================
+
+import unittest
+import numpy as np
+
+from pathsim_chem.tritium import Bubbler4
+
+from pathsim.solvers import EUF
+from pathsim.events.schedule import ScheduleList
+
+
+# TESTS ================================================================================
+
+class TestFusionBubbler4(unittest.TestCase):
+    """
+    Test the implementation of the 'Bubbler4' block class from the fusion toolbox.
+    The block inherits from `ODE` and models a 4-vial tritium collection system.
+    """
+    
+    def test_init(self):
+        """Test initialization with various parameter combinations"""
+        
+        # Default initialization
+        B = Bubbler4()
+        self.assertEqual(B.conversion_efficiency, 0.9)
+        self.assertEqual(B.vial_efficiency, 0.9)
+        self.assertEqual(B.replacement_times, None)
+        self.assertEqual(len(B.events), 0)  # No events when replacement_times is None
+        
+        # Specific initialization
+        B = Bubbler4(conversion_efficiency=0.8, vial_efficiency=0.7, replacement_times=[100, 200])
+        self.assertEqual(B.conversion_efficiency, 0.8)
+        self.assertEqual(B.vial_efficiency, 0.7)
+        self.assertEqual(B.replacement_times, [100, 200])
+        
+        # Set solver to check internal solver instance
+        B.set_solver(EUF, parent=None)
+        self.assertTrue(B.engine)
+        np.testing.assert_array_equal(B.engine.initial_value, np.zeros(4))
+    
+
+    def test_init_replacement_times(self):
+        """Test different replacement_times configurations"""
+        
+        # Single list - should replicate for all vials
+        B = Bubbler4(replacement_times=[100, 200, 300])
+        B.set_solver(EUF, parent=None)
+        self.assertEqual(len(B.events), 4)
+        
+        # List of lists - one per vial
+        times = [[100, 200], [150, 250], [120, 220], [180, 280]]
+        B = Bubbler4(replacement_times=times)
+        B.set_solver(EUF, parent=None)
+        self.assertEqual(len(B.events), 4)
+        for event in B.events:
+            self.assertIsInstance(event, ScheduleList)
+        
+        # NumPy array
+        B = Bubbler4(replacement_times=np.array([100, 200, 300]))
+        B.set_solver(EUF, parent=None)
+        self.assertEqual(len(B.events), 4)
+        
+        # Invalid case - wrong length
+        with self.assertRaises(ValueError):
+            B = Bubbler4(replacement_times=[[100], [200], [300]])  # Only 3 vials
+    
+
+    def test_update_outputs(self):
+        """Test the update method and output calculations"""
+        
+        # Test with default parameters
+        B = Bubbler4()
+        B.set_solver(EUF, parent=None)
+        
+        # Set some test inputs
+        B.inputs[0] = 10.0  # soluble input
+        B.inputs[1] = 5.0   # insoluble input
+        
+        B.update(0.0)
+        
+        # Check that outputs are set (initial state is zero)
+        self.assertEqual(B.outputs[0], 0.0)  # vial 1
+        self.assertEqual(B.outputs[1], 0.0)  # vial 2
+        self.assertEqual(B.outputs[2], 0.0)  # vial 3
+        self.assertEqual(B.outputs[3], 0.0)  # vial 4
+        
+        # Calculate expected sample_out with default parameters (ve=0.9, ce=0.9)
+        ve, ce = 0.9, 0.9
+        sol, ins = 10.0, 5.0
+        expected_sample_out = (1-ce)*ins + (1-ve)**2 * (ce*ins + (1-ve)**2 * sol)
+        self.assertAlmostEqual(B.outputs[4], expected_sample_out, places=10)
+
+
+    def test_mass_conservation(self):
+        """Test that mass is conserved in the system"""
+        
+        B = Bubbler4(conversion_efficiency=0.8, vial_efficiency=0.7)
+        B.set_solver(EUF, parent=None)
+        
+        # Set test inputs
+        sol_in, ins_in = 10.0, 5.0
+        B.inputs[0] = sol_in
+        B.inputs[1] = ins_in
+        
+        # Get rates and sample_out
+        u = np.array([sol_in, ins_in])
+        x = np.zeros(4)
+        rates = B.func_dyn(x, u, 0.0)
+        
+        B.update(0.0)
+        sample_out = B.outputs[4]
+        
+        # Total input rate
+        total_in = sol_in + ins_in
+        
+        # Total accumulation rate in vials
+        total_vial_rates = np.sum(rates)
+        
+        # Mass conservation: input = vial accumulation + sample output
+        total_out = total_vial_rates + sample_out
+        
+        self.assertAlmostEqual(total_in, total_out, places=10)
+    
+
+    def test_vial_reset_functionality(self):
+        """Test that vial reset events work correctly"""
+        
+        B = Bubbler4(replacement_times=[10.0])
+        B.set_solver(EUF, parent=None)
+        
+        # Manually set some vial inventories
+        x = np.array([1.0, 2.0, 3.0, 4.0])
+        B.engine.set(x)
+        
+        # Test reset function for vial 0
+        reset_event = B.events[0]
+        reset_func = reset_event.func_act
+        
+        # Execute reset
+        reset_func(None)
+        
+        # Check that vial 0 was reset but others remain
+        x_after = B.engine.get()
+        self.assertEqual(x_after[0], 0.0)
+        self.assertEqual(x_after[1], 2.0)
+        self.assertEqual(x_after[2], 3.0)
+        self.assertEqual(x_after[3], 4.0)
+
+        # Test reset function for vial 3
+        reset_event = B.events[3]
+        reset_func = reset_event.func_act
+        
+        # Execute reset
+        reset_func(None)
+        
+        # Check that vial 0 was reset but others remain
+        x_after = B.engine.get()
+        self.assertEqual(x_after[0], 0.0)
+        self.assertEqual(x_after[1], 2.0)
+        self.assertEqual(x_after[2], 3.0)
+        self.assertEqual(x_after[3], 0.0)
+
+
+# RUN TESTS LOCALLY ====================================================================
+
+if __name__ == '__main__':
+    unittest.main(verbosity=2)
\ No newline at end of file
diff --git a/tests/tritium/test_residencetime.py b/tests/tritium/test_residencetime.py
new file mode 100644
index 0000000..447f4ae
--- /dev/null
+++ b/tests/tritium/test_residencetime.py
@@ -0,0 +1,115 @@
+########################################################################################
+##
+##                                  TESTS FOR 
+##                          'tritium.residencetime.py'
+##
+########################################################################################
+
+# IMPORTS ==============================================================================
+
+import unittest
+import numpy as np
+
+from pathsim_chem.tritium import ResidenceTime, Process
+
+from pathsim.solvers import EUF
+
+
+# TESTS ================================================================================
+
+class TestTritiumResidenceTime(unittest.TestCase):
+    """
+    Test the implementation of the 'ResidenceTime' block class from the fusion toolbox
+
+    The block inherits from `ODE`
+    """
+    
+    def test_init(self):
+
+        #default initialization
+        R = ResidenceTime()
+
+        self.assertEqual(R.tau, 1)
+        self.assertEqual(R.betas, 1)
+        self.assertEqual(R.gammas, 1)
+        self.assertEqual(R.source_term, 0)
+
+        #input validation
+        with self.assertRaises(ValueError):
+            R = ResidenceTime(tau=0)
+
+        #specific initialization
+        R = ResidenceTime(tau=0.1, betas=[1, 3, 5], gammas=[2, 4], initial_value=10, source_term=11)
+
+        self.assertEqual(R.tau,0.1)
+        self.assertTrue(np.allclose(R.betas, np.array([1, 3, 5])))
+        self.assertTrue(np.allclose(R.gammas, np.array([2, 4])))
+
+        #set solver to check internal solver instance
+        R.set_solver(EUF, parent=None)
+
+        self.assertTrue(R.engine)
+        self.assertEqual(R.engine.initial_value, 10)
+
+
+    def test_update(self):
+
+        #default
+        R = ResidenceTime()
+        R.set_solver(EUF, parent=None)
+
+        R.update(None)
+
+        self.assertEqual(R.outputs[0], 0)
+
+        #specific
+        R = ResidenceTime(tau=0.1, betas=[1, 3, 5], gammas=[2, 4], initial_value=11, source_term=12)
+        R.set_solver(EUF, parent=None)
+
+        R.update(None)
+
+        self.assertEqual(R.outputs[0], 22) # initial_value * gammas[0]
+        self.assertEqual(R.outputs[1], 44) # initial_value * gammas[1]
+
+
+
+class TestTritiumProcess(unittest.TestCase):
+    """
+    Test the implementation of the 'Process' block class from the fusion toolbox.
+
+    This block inherits from `ResidenceTime`, just testing the initialization.
+    """
+
+    def test_init(self):
+
+        #default initialization
+        P = Process()
+
+        self.assertEqual(P.tau, 1)
+        self.assertEqual(P.betas, 1)
+        self.assertTrue(np.allclose(P.gammas, [1, 1])) # 1, 1/tau
+        self.assertEqual(P.source_term, 0)
+
+        #specific initialization
+        P = Process(tau=0.1, source_term=33, initial_value=44)
+
+        self.assertEqual(P.tau,0.1)
+        self.assertEqual(P.betas, 1)
+        self.assertTrue(np.allclose(P.gammas, np.array([1, 10])))
+        self.assertEqual(P.source_term, 33)
+
+        #set solver to check internal solver instance
+        P.set_solver(EUF, parent=None)
+
+        self.assertTrue(P.engine)
+        self.assertEqual(P.engine.initial_value, 44)
+
+
+
+
+
+
+# RUN TESTS LOCALLY ====================================================================
+
+if __name__ == '__main__':
+    unittest.main(verbosity=2)
diff --git a/tests/tritium/test_splitter.py b/tests/tritium/test_splitter.py
new file mode 100644
index 0000000..71f3e97
--- /dev/null
+++ b/tests/tritium/test_splitter.py
@@ -0,0 +1,61 @@
+########################################################################################
+##
+##                                  TESTS FOR 
+##                            'tritium.splitter.py'
+##
+########################################################################################
+
+# IMPORTS ==============================================================================
+
+import unittest
+import numpy as np
+
+from pathsim_chem.tritium import Splitter
+
+
+# TESTS ================================================================================
+
+class TestTritiumSplitter(unittest.TestCase):
+    """
+    Test the implementation of the 'Splitter' block class from the fusion toolbox
+    """
+
+    def test_init(self):
+        
+        #default initialization
+        S = Splitter()
+        self.assertEqual(S.fractions, np.ones(1))
+
+        #input validation
+        for fracs in [[1, 3], [0.4, 0.6, 0.001], [0.33, 0.33, 0.33]]:
+            with self.assertRaises(ValueError):
+                S = Splitter(fracs)
+
+        #special initialization
+        S = Splitter([0.4, 0.5, 0.1])
+        self.assertEqual(sum(S.fractions - np.array([0.4, 0.5, 0.1])), 0)
+
+        #test the automatic port maps
+        self.assertEqual(S._port_map_out, {"out 0.4":0, "out 0.5":1, "out 0.1":2})
+
+
+    def test_update(self):
+
+        S = Splitter([0.4, 0.5, 0.1])
+
+        #set block inputs
+        S.inputs[0] = 2
+
+        #update block
+        S.update(None)
+
+        #test if update was correct
+        self.assertEqual(S.outputs[0], 0.8)
+        self.assertEqual(S.outputs[1], 1)
+        self.assertEqual(S.outputs[2], 0.2)
+
+
+# RUN TESTS LOCALLY ====================================================================
+
+if __name__ == '__main__':
+    unittest.main(verbosity=2)