diff --git a/examples/pstl.py b/examples/pstl.py
index bbf234d..e060681 100644
--- a/examples/pstl.py
+++ b/examples/pstl.py
@@ -47,6 +47,23 @@ def binary_search_settling(ys, save_filename="models/pstl_settling_binary_search
         return np.stack(ϵ_list)
     np.save(save_filename, np.stack(ϵ_list))
 
+def binary_search_settling_vectorize(ys):
+    N = ys.shape[0]
+    y = torch.as_tensor(ys).float().unsqueeze(-1)
+    s = stlcg.Expression('s', torch.abs(y - 1))
+    ϵL = torch.as_tensor(np.zeros([N, 1, 1])).float().requires_grad_(True)
+    ϵU = torch.as_tensor(np.ones([N, 1, 1])).float().requires_grad_(True)
+    ϵ = torch.as_tensor(np.ones([N, 1, 1])).float().requires_grad_(True)
+    ϕL = stlcg.Always(subformula=(s < ϵL), interval=[50, 100])
+    ϕU = stlcg.Always(subformula=(s < ϵU), interval=[50, 100])
+    ϕ =  stlcg.Always(subformula=(s < ϵ), interval=[50, 100])
+    while torch.abs(ϕU.subformula.val - ϕL.subformula.val).max() > 5*1E-3:
+        ϵ = 0.5 * (ϕU.subformula.val + ϕL.subformula.val)
+        ϕ.subformula.val = ϵ
+        r = ϕ.robustness(s)
+        ϕU.subformula.val = torch.where(torch.abs(ϕU.subformula.val - ϕL.subformula.val) > 5*1E-3, torch.where(r > 0, ϵ, ϕU.subformula.val), ϕU.subformula.val)
+        ϕL.subformula.val = torch.where(torch.abs(ϕU.subformula.val - ϕL.subformula.val) > 5*1E-3, torch.where(r <= 0, ϵ, ϕL.subformula.val), ϕL.subformula.val)
+    return ϵ.squeeze().detach().numpy()
 
 def stlcg_settling(ys, save_filename="models/pstl_settling_stlcg.npy"):
     max_epoch = 1000
@@ -97,6 +114,24 @@ def binary_search_overshoot(ys, save_filename="models/pstl_overshoot_binary_sear
     if save_filename is None:
         return np.stack(ϵ_list)
     np.save(save_filename, np.stack(ϵ_list))
+    
+def binary_search_overshoot_vectorize(ys):
+    N = ys.shape[0]
+    y = torch.as_tensor(ys).float().unsqueeze(-1)
+    s = stlcg.Expression('s', y)
+    ϵL = torch.as_tensor(np.zeros([N, 1, 1])).float().requires_grad_(True)
+    ϵU = torch.as_tensor(2*np.ones([N, 1, 1])).float().requires_grad_(True)
+    ϵ = torch.as_tensor(np.ones([N, 1, 1])).float().requires_grad_(True)
+    ϕL = stlcg.Always(subformula=(s < ϵL))
+    ϕU = stlcg.Always(subformula=(s < ϵU))
+    ϕ =  stlcg.Always(subformula=(s < ϵ))
+    while torch.abs(ϕU.subformula.val - ϕL.subformula.val).max() > 5*1E-3:
+        ϵ = 0.5 * (ϕU.subformula.val + ϕL.subformula.val)
+        ϕ.subformula.val = ϵ
+        r = ϕ.robustness(s)
+        ϕU.subformula.val = torch.where(torch.abs(ϕU.subformula.val - ϕL.subformula.val) > 5*1E-3, torch.where(r > 0, ϵ, ϕU.subformula.val), ϕU.subformula.val)
+        ϕL.subformula.val = torch.where(torch.abs(ϕU.subformula.val - ϕL.subformula.val) > 5*1E-3, torch.where(r <= 0, ϵ, ϕL.subformula.val), ϕL.subformula.val)
+    return ϵ.squeeze().detach().numpy()
 
 def stlcg_overshoot(ys, save_filename="models/pstl_overshoot_stlcg.npy"):
     N = ys.shape[0]