Add configurable inequality tolerance to the VMCON solver

process/data_structure/numerics.py

@@ -427,15 +427,25 @@
 
 force_vmcon_inequality_satisfication: int = None
 """If 1, adds an additional convergence criteria to the VMCON solver
-that enforces the inequality constraints with zero margin/tolerance.
-I.e. VMCON cannot converge until all inequality constraints are
-exactly satisfied.
+that enforces a margin on the inequality constraints.
+I.e. VMCON cannot converge until all inequality constraints are satisfied
+to within a tolerance of `force_vmcon_inequality_tolerance`.
 
-Default is 0.
+Default is 1 (enabled).
 
 NOTE: this only affects the VMCON solver.
 """
 
+force_vmcon_inequality_tolerance: float = None
+"""The relative tolerance for the additional VMCON convergence criteria
+that forces inequality constraints to be satisfied within a tolerance.
+
+Default is 1e-8.
+
+NOTE: has no effect if `force_vmcon_inequality_satisfication` is 0
+NOTE: this only affects the VMCON solver.
+"""
+
 
 def init_numerics():
     global ipnvars
@@ -480,6 +490,7 @@ def init_numerics():
     global xcs
     global vlam
     global force_vmcon_inequality_satisfication
+    global force_vmcon_inequality_tolerance
 
     """Initialise module variables"""
     ioptimz = 1
@@ -639,4 +650,5 @@ def init_numerics():
     xcs = np.array([0.0] * ipnvars)
     vlam = np.array([0.0] * ipnvars)
     name_xc = [""] * ipnvars
-    force_vmcon_inequality_satisfication = 0
+    force_vmcon_inequality_satisfication = 1
+    force_vmcon_inequality_tolerance = 1e-8

process/input.py

@@ -2339,6 +2339,9 @@ def __post_init__(self):
         int,
         choices=(0, 1),
     ),
+    "force_vmcon_inequality_tolerance": InputVariable(
+        data_structure.numerics, float, range=(0.0, 1e10)
+    ),
 }
 
 

process/solver.py

@@ -204,8 +204,10 @@ def _ineq_cons_satisfied(
             :return: True if inequality constraints satisfied
             :rtype: bool
             """
-            # Check all ineqs positive, i.e. satisfied
-            return bool(np.all(result.ie >= 0.0))
+            # negative constraint value = violated
+            # Check all ineqs are satisfied to within the tolerance
+            # E.g. the relative violations are no more than v=0-tolerance
+            return bool(np.all(result.ie >= -numerics.force_vmcon_inequality_tolerance))
 
         try:
             x, _, _, res = solve(