diff --git a/UB4/UB4.py b/UB4/UB4.py
index 4c6ed980b2da31579a3accbc4aa6983489baffed..1d81dfccd2463fdb96a4d70c2fdc2629f144456c 100644
--- a/UB4/UB4.py
+++ b/UB4/UB4.py
@@ -1,9 +1,10 @@
 import numpy as np
 import matplotlib.pyplot  as plt
+from tqdm import trange
 
 # Ex1: Numerical integrators 
 
-def implicit_euler(f, y0: np.ndarray, t: float, dt: float, **kwargs) -> np.ndarray:
+def implicit_euler(f, y0: np.ndarray, t: float, dt: float, eps = 1e-5, **kwargs) -> np.ndarray:
     """
     A-stable
     :param f: function to be integrated
@@ -14,10 +15,15 @@ def implicit_euler(f, y0: np.ndarray, t: float, dt: float, **kwargs) -> np.ndarr
     no_of_steps = int(t // dt)
     y = np.zeros((no_of_steps, *y0.shape))
     y[0] = y0
-    
-    for i in range(1, no_of_steps):
-        y[i] = y[i-1] + f(y[i-1], **kwargs) * dt
-        y[i] = y[i-1] + f(y[i], **kwargs) * dt
+
+    for i in trange(1, no_of_steps):
+        y_0 = y[i-1] + dt * f(y[i-1]) # predicted state obtained through EE
+        y_1 = y[i-1] + dt * f(y_0)
+
+        while abs(y_1 - y_0) > eps:
+            y_0 = y_1
+            y_1 = y[i-1] + dt * f(y_0)
+        y[i] = y_1
 
     return y
 
@@ -32,14 +38,13 @@ def explicit_euler(f, y0: np.ndarray, t: float, dt: float, **kwargs) -> np.ndarr
     y = np.zeros((no_of_steps, *y0.shape))
     y[0] = y0
 
-    for i in range(1, no_of_steps):
+    for i in trange(1, no_of_steps):
         y[i] = y[i-1] + f(y[i-1], **kwargs) * dt
     
     return y
 
-def implicit_midpoint(f, y0: np.ndarray, t: float, dt: float, **kwargs) -> np.ndarray:
+def implicit_midpoint(f, y0: np.ndarray, t: float, dt: float, eps = 1e-5, **kwargs) -> np.ndarray:
     """
-    Not L-stable, doesn't decay properly - oscillation
     :param f: function to be integrated
     :param y0: initial value
     :param t: time interval
@@ -48,18 +53,24 @@ def implicit_midpoint(f, y0: np.ndarray, t: float, dt: float, **kwargs) -> np.nd
     no_of_steps = int(t // dt)
     y = np.zeros((no_of_steps, *y0.shape))
     y[0] = y0
-    
-    for i in range(1, no_of_steps):
-        y[i] = y[i-1] + dt * f(y[i-1] + dt/2 * f(y[i-1], **kwargs), **kwargs)
-        y[i] = y[i-1] + dt * f(1/2 *  (y[i-1] + y[i]), **kwargs)
-    
+
+    for i in trange(1, no_of_steps):
+        y_0 = y[i-1] + dt * f(y[i-1] + dt/2 * f(y[i-1], **kwargs), **kwargs) # predicted state obtained through IM
+        y_1 = y[i-1] + dt * f( 1/2 * (y[i-1] + y_0) ,**kwargs)
+        
+        while abs(y_1 - y_0) > eps:
+            y_0 = y_1
+            y_1 = y[i-1] + dt * f( 1/2 * (y[i-1] + y_0) ,**kwargs)
+        y[i] = y_1
+
     return y
 
+
 # Ex2: dy/dt = lambda * y, y(0) = 1, lambda = -1
-y0 = np.array([1])
+y0 = np.array(1)
 lamda = -1
 t = [10, 100]
-dt = [1.5, 1, 1e-1, 1e-2]
+dt = [0.9, 1e-1, 1e-2]
 f = lambda y: lamda * y
 exact_solution = lambda t : np.exp(lamda * t) * y0
 
@@ -84,14 +95,8 @@ exact_solution = lambda t : np.exp(lamda * t) * y0
 
 #### For A-stability, dt < 1
 #### Observation
-### t = 10, dt = 1.5 <-- Strange behaviour !
-# IE grows and EE decays!
-# IM is robust
-
-### t = 10, dt = 1 <-- Strange behaviour !
-# IE gives a constant function (but it doesn't blow up as time progresses!).
-# However, this means it's not A-stable because it does not appropriately damp the solution. 
-# EE & IM decay faster than ES -> L-stable. EE gives a worse approximation in comparison to IM
+### t = 10, dt = 0.9  <-- Strange behaviour ! 
+# IM decas as rapidly as ES, IE not as fast, EE too fast
 ### t = 10, dt = 1e-1  <-- Strange behaviour ! 
 # They all behave similarly. They do not decay as rapidly as ES.
 ### t = 10, dt = 1e-2 <-- Strange behaviour !
@@ -110,17 +115,14 @@ for i in range(len(dt)):
     plt.show()
 
 #### Observation
-### t = 100, dt = 1.5 <-- Strange behaviour !
-# IE blows up at the end.
-### t = 100, dt = 1 <-- Strange behaviour !
-# Same as in the case t = 10, dt = 1
-### t = 100, dt = 1e-1  <-- Strange behaviour ! 
-# Same as in the case t = 100, dt = 1e-1
-### t = 100, dt = 1e-2 <-- Strange behaviour !
+### t = 100, dt = 0.9  <-- Strange behaviour ! 
+# Same as in the case t = 10, dt = 0.9
+### t = 100, dt = 0.1 <-- Strange behaviour !
+# Same as in the case t = 100, dt = 0.1
+### t = 100, dt = 1e-2  <-- Strange behaviour ! 
 # Same as in the case t = 100, dt = 1e-2
 
 
-
 # Ex3: 
 # def f(y: np.ndarray, k: float) -> np.ndarray:
 #     assert y.shape[0] == 3, 'y has the wrong dimension. It should be 3'