diff --git a/process/core/solver/constraints.py b/process/core/solver/constraints.py
index c5d786d0f8..f6e3d057cd 100644
--- a/process/core/solver/constraints.py
+++ b/process/core/solver/constraints.py
@@ -18,7 +18,6 @@
 from process.models.physics.exhaust import PlasmaExhaust
 from process.models.physics.physics import (
     BetaComponentLimits,
-    PlasmaBeta,
     PlasmaIgnitionModel,
 )
 from process.models.tfcoil.base import TFConductorModel
@@ -216,43 +215,6 @@ def eq(value: float, bound: float, registration: ConstraintRegistration):
     )
 
 
-@ConstraintManager.register_constraint(1, "", "=")
-def constraint_equation_1(constraint_registration, data):
-    """Relationship between beta, temperature (keV) and density
-
-    beta_total_vol_avg: total plasma beta
-    beta_{ft}: fast alpha beta component
-    beta_{NBI}: neutral beam beta component
-    n_e: electron density [/m3]
-    n_i: total ion density [/m3]
-    T_e: density weighted average electron temperature [keV]
-    T_i: density weighted average ion temperature [keV]
-    B_{tot}: total toroidal + poloidal field [T]
-    """
-    # Density weighted temperature is used here as 〈nT〉 != 〈n〉_V * 〈T〉_V
-    beta_thermal_total_vol_avg = PlasmaBeta.calculate_plasma_beta(
-        pres_plasma=(
-            constants.KILOELECTRON_VOLT
-            * (
-                data.physics.nd_plasma_electrons_vol_avg
-                * data.physics.temp_plasma_electron_density_weighted_kev
-                + data.physics.nd_plasma_ions_total_vol_avg
-                * data.physics.temp_plasma_ion_density_weighted_kev
-            )
-        ),
-        b_field=data.physics.b_plasma_total,
-    )
-    return eq(
-        (
-            data.physics.beta_fast_alpha
-            + data.physics.beta_beam
-            + beta_thermal_total_vol_avg
-        ),
-        data.physics.beta_total_vol_avg,
-        constraint_registration,
-    )
-
-
 @ConstraintManager.register_constraint(2, "MW/m³", "=")
 def constraint_equation_2(constraint_registration, data):
     """
diff --git a/process/core/solver/iteration_variables.py b/process/core/solver/iteration_variables.py
index fd468c5537..1e4665bb6b 100644
--- a/process/core/solver/iteration_variables.py
+++ b/process/core/solver/iteration_variables.py
@@ -46,7 +46,6 @@ class IterationVariable:
     2: IterationVariable("b_plasma_toroidal_on_axis", "physics", 0.010, 30.00),
     3: IterationVariable("rmajor", "physics", 0.1, 50.00),
     4: IterationVariable("temp_plasma_electron_vol_avg_kev", "physics", 5.0, 150.0),
-    5: IterationVariable("beta_total_vol_avg", "physics", 0.001, 1.0),
     6: IterationVariable("nd_plasma_electrons_vol_avg", "physics", 2.0e19, 1.0e21),
     7: IterationVariable("f_nd_beam_electron", "physics", 1.0e-6, 1.0),
     10: IterationVariable("hfact", "physics", 0.1, 3.0),
diff --git a/process/data_structure/numerics.py b/process/data_structure/numerics.py
index 85785576a4..d2d318666e 100644
--- a/process/data_structure/numerics.py
+++ b/process/data_structure/numerics.py
@@ -200,7 +200,7 @@ class NumericsData:
 
     lablcc: list[str] = field(
         default_factory=lambda: [
-            "⟨β⟩ consistency                   ",
+            "NOT USED",
             "Global power balance consistency ",
             "Ion power balance                ",
             "Electron power balance           ",
@@ -295,11 +295,11 @@ class NumericsData:
         ]
     )
     """Labels describing constraint equations (corresponding itvs)<UL>
-    * ( 1) Beta (consistency equation) (itv 5)
+    * ( 1) NOT USED
     * ( 2) Global power balance (consistency equation) (itv 10,1,2,3,4,6,11)
     * ( 3) Ion power balance DEPRECATED (itv 10,1,2,3,4,6,11)
     * ( 4) Electron power balance DEPRECATED (itv 10,1,2,3,4,6,11)
-    * ( 5) Density upper limit (itv 9,1,2,3,4,5,6)
+    * ( 5) Density upper limit (itv 9,1,2,3,4,6)
     * ( 6) (Epsilon x beta poloidal) upper limit (itv 8,1,2,3,4,6)
     * ( 7) Beam ion density (NBI) (consistency equation) (itv 7)
     * ( 8) Neutron wall load upper limit (itv 14,1,2,3,4,6)
@@ -409,7 +409,6 @@ class NumericsData:
     * ( 2) b_plasma_toroidal_on_axis
     * ( 3) rmajor
     * ( 4) temp_plasma_electron_vol_avg_kev
-    * ( 5) beta_total_vol_avg
     * ( 6) nd_plasma_electrons_vol_avg
     * ( 7) f_nd_beam_electron
     * ( 8) NOT USED
diff --git a/process/models/physics/physics.py b/process/models/physics/physics.py
index 43975e3f0e..9189f30d20 100644
--- a/process/models/physics/physics.py
+++ b/process/models/physics/physics.py
@@ -430,12 +430,6 @@ def run(self):
 
         # ============================================
 
-        # -----------------------------------------------------
-        # Beta Components
-        # -----------------------------------------------------
-
-        self.beta.run()
-
         # Stored thermal energies in the plasma
 
         (
@@ -694,6 +688,9 @@ def run(self):
             f_plasma_fuel_deuterium=self.data.physics.f_plasma_fuel_deuterium,
         )
 
+        # Calculate beta after the fast-particle contributions have been updated.
+        self.beta.run()
+
         # Calculate ion/electron equilibration power
 
         self.data.physics.pden_ion_electron_equilibration_mw = rether(
@@ -3506,6 +3503,27 @@ def run(self):
              If an illegal value of i_beta_norm_max is provided.
 
         """
+        # Thermal beta is derived directly from the current plasma state.
+        self.data.physics.beta_thermal_vol_avg = self.calculate_plasma_beta(
+            pres_plasma=(
+                constants.KILOELECTRON_VOLT
+                * (
+                    self.data.physics.nd_plasma_electrons_vol_avg
+                    * self.data.physics.temp_plasma_electron_density_weighted_kev
+                    + self.data.physics.nd_plasma_ions_total_vol_avg
+                    * self.data.physics.temp_plasma_ion_density_weighted_kev
+                )
+            ),
+            b_field=self.data.physics.b_plasma_total,
+        )
+
+        # Total beta is derived from the thermal and fast-particle components.
+        self.data.physics.beta_total_vol_avg = (
+            self.data.physics.beta_thermal_vol_avg
+            + self.data.physics.beta_fast_alpha
+            + self.data.physics.beta_beam
+        )
+
         # -----------------------------------------------------
         # Normalised Beta Limit
         # -----------------------------------------------------
@@ -3583,12 +3601,6 @@ def run(self):
             beta=self.data.physics.beta_total_vol_avg,
         )
 
-        self.data.physics.beta_thermal_vol_avg = (
-            self.data.physics.beta_total_vol_avg
-            - self.data.physics.beta_fast_alpha
-            - self.data.physics.beta_beam
-        )
-
         self.data.physics.beta_poloidal_eps = (
             self.data.physics.beta_poloidal_vol_avg * self.data.physics.eps
         )
diff --git a/tests/regression/input_files/large_tokamak_eval.IN.DAT b/tests/regression/input_files/large_tokamak_eval.IN.DAT
index 66a44b7c0a..ea66f6dd51 100644
--- a/tests/regression/input_files/large_tokamak_eval.IN.DAT
+++ b/tests/regression/input_files/large_tokamak_eval.IN.DAT
@@ -3,10 +3,9 @@ i_process_run_mode  = -2
 
 *---------------Constraint Equations---------------*
 * Define number of equality constraints
-n_equality_constraints = 2
+n_equality_constraints = 1
 
 * Equalities
-icc = 1 * Beta
 icc = 2 * Global power balance
 
 * Inequalities
diff --git a/tests/regression/input_files/large_tokamak_nof.IN.DAT b/tests/regression/input_files/large_tokamak_nof.IN.DAT
index 6bdc4128be..8acbca19e2 100644
--- a/tests/regression/input_files/large_tokamak_nof.IN.DAT
+++ b/tests/regression/input_files/large_tokamak_nof.IN.DAT
@@ -19,7 +19,7 @@ epsvmc = 1e-7
 
 * Number of equality constraints
 * (the first n_equality_constraints icc=... statements are equality constraints)
-n_equality_constraints = 3
+n_equality_constraints = 2
 
 * Maximum number of solver iterations
 maxcal = 100
@@ -29,8 +29,6 @@ maxcal = 100
 
 * Beta consistency *
 *------------------*
-icc = 1
-ixc = 5 * beta_total_vol_avg
 beta_total_vol_avg = 0.03
 
 * Global power balance *
diff --git a/tests/regression/input_files/low_aspect_ratio_DEMO.IN.DAT b/tests/regression/input_files/low_aspect_ratio_DEMO.IN.DAT
index fc7557a2ca..f103b7046f 100644
--- a/tests/regression/input_files/low_aspect_ratio_DEMO.IN.DAT
+++ b/tests/regression/input_files/low_aspect_ratio_DEMO.IN.DAT
@@ -33,15 +33,12 @@ i_figure_merit   = -14 * Switch for figure-of-merit (see lablmm for descriptions
 * DESCRIPTION:   Switch for Figure-of-Merit (1: Minimise major radius)
 * JUSTIFICATION: Aim to minimise major radius
 
-n_equality_constraints = 4
+n_equality_constraints = 3
 * DESCRIPTION:   No of equality constraints
 * JUSTIFICATION: 
 
 *--------------------------Constraint Equations-------------------------*
 
-icc = 1 
-* DESCRIPTION:   Beta Consistency
-* JUSTIFICATION: Consistency equations should always be on
 
 icc = 2 
 * DESCRIPTION:   Global Power Balance Consistency
@@ -173,7 +170,6 @@ boundu(4) = 150.0
 * DESCRIPTION:   Volume averaged electron temperature (keV)
 * JUSTIFICATION: Iterate temperature to get solution
 
-ixc = 5 
 beta_total_vol_avg = 3.71965626913446160e-02
 * DESCRIPTION:   Total Plasma Beta
 * JUSTIFICATION: Beta is limited by dnbeta
diff --git a/tests/regression/input_files/spherical_tokamak_eval.IN.DAT b/tests/regression/input_files/spherical_tokamak_eval.IN.DAT
index f79f323f4e..f204a891a0 100644
--- a/tests/regression/input_files/spherical_tokamak_eval.IN.DAT
+++ b/tests/regression/input_files/spherical_tokamak_eval.IN.DAT
@@ -3,10 +3,9 @@ i_process_run_mode = -2
 
 *---------------Constraint Equations---------------*
 * Define number of equality constraints
-n_equality_constraints = 3
+n_equality_constraints = 2
 
 * Equalities
-icc = 1 * Beta
 icc = 2 * Global power balance
 icc = 11 * Radial build
 
diff --git a/tests/regression/input_files/st_regression.IN.DAT b/tests/regression/input_files/st_regression.IN.DAT
index 5b2f404a4c..a112c0d424 100644
--- a/tests/regression/input_files/st_regression.IN.DAT
+++ b/tests/regression/input_files/st_regression.IN.DAT
@@ -55,7 +55,7 @@ i_figure_merit = -5
 * DESCRIPTION:   Switch for Figure-of-Merit 
 * JUSTIFICATION: Selected optimisation
 
-n_equality_constraints = 3
+n_equality_constraints = 2
 * DESCRIPTION:   Number of equality constraints (the first n_equality_constraints icc = ... statements are considered equality)
 * JUSTIFICATION: 3 equality (consistency) equations
 
@@ -70,10 +70,6 @@ itart = 1
 *_______________________________________________________________________________________*
 *------------------------------- *Consistency Equations* -------------------------------*
 *‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾‾*
-icc = 1
-* DESCRIPTION:   Beta Consistency
-* JUSTIFICATION: Consistency equations should always be on
-* VARIABLES:     Calculated in-situ, no iteration variables
 
 icc = 2
 * DESCRIPTION:   Global Power Balance Consistency
@@ -295,13 +291,6 @@ icc = 24
 * DESCRIPTION:   Plasma beta upper limit
 * JUSTIFICATION: Turned off, not in use. icc=24. Allow to be calculated in-situ
 
-ixc = 5
-beta_total_vol_avg = 0.2
-boundl(5) = 0.01
-boundu(5) = 1.0
-* DESCRIPTION:   Total Plasma Beta, 
-* JUSTIFICATION: Beta is limited by beta_norm_max
-
 
 *~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~