diff --git a/lib/LeapFrog.py b/lib/LeapFrog.py
index c786abb..f8e937b 100644
--- a/lib/LeapFrog.py
+++ b/lib/LeapFrog.py
@@ -36,7 +36,7 @@ def LP(dyn_syst, dt):
         body.p = body.v * body.m
 
 
-def leapfrog(dyn_syst, duration, dt, recover_param=False, display=False, savename=None):
+def leapfrog(dyn_syst, bin_syst, duration, dt, recover_param=False, display=False, savename=None):
     if display:
         try:
             system("mkdir tmp")
@@ -48,11 +48,15 @@ def leapfrog(dyn_syst, duration, dt, recover_param=False, display=False, savenam
     N = np.ceil(duration / dt).astype(int)
     E = np.zeros(N)
     L = np.zeros((N, 3))
+    sma = np.zeros(N)
+    ecc = np.zeros(N)
     for j in range(N):
         LP(dyn_syst,dt)
 
         E[j] = dyn_syst.E
         L[j] = dyn_syst.L
+        sma[j] = bin_syst.sma
+        ecc[j] = bin_syst.ecc
 
         if display and j % 5 == 0:
             # display progression
@@ -67,4 +71,4 @@ def leapfrog(dyn_syst, duration, dt, recover_param=False, display=False, savenam
             system("convert tmp/temp.gif -fuzz 10% -layers Optimize plots/{0:s}_dynsyst.gif".format(savename))
 
     if recover_param:
-        return E, L
\ No newline at end of file
+        return E, L, sma, ecc
\ No newline at end of file
diff --git a/lib/hermite.py b/lib/hermite.py
index 8a19900..7dc5694 100644
--- a/lib/hermite.py
+++ b/lib/hermite.py
@@ -93,11 +93,16 @@ def hermite(dyn_syst, duration, dt, recover_param=False, display=False, savename
     N = np.ceil(duration / dt).astype(int)
     E = np.zeros(N)
     L = np.zeros((N, 3))
+    sma = np.zeros(N)
+    ecc = np.zeros(N)
+
     for j in range(N):
         HPC(dyn_syst, dt)
 
         E[j] = dyn_syst.E
         L[j] = dyn_syst.L
+        sma[j] = bin_syst.sma
+        ecc[j] = bin_syst.ecc
 
         if display and j % 100 == 0:
             # display progression
@@ -112,4 +117,4 @@ def hermite(dyn_syst, duration, dt, recover_param=False, display=False, savename
             system("convert tmp/temp.gif -fuzz 10% -layers Optimize plots/{0:s}_dynsyst.gif".format(savename))
 
     if recover_param:
-        return E, L
\ No newline at end of file
+        return E, L, sma, ecc
\ No newline at end of file
diff --git a/lib/objects.py b/lib/objects.py
index 4a85b61..3fc48f6 100755
--- a/lib/objects.py
+++ b/lib/objects.py
@@ -23,10 +23,10 @@ class Body:
         self.vp = np.zeros(3)
 
     def __repr__(self): # Called upon "print(body)"
-        return r"Body of mass: {0:.2f} $M_\odot$, position: {1}, velocity: {2}".format(self.m/Ms, self.q, self.v)
+        return r"Body of mass: {0:.1e} $M_\odot$, position: {1}, velocity: {2}".format(self.m/Ms, self.q, self.v)
 
     def __str__(self): # Called upon "str(body)"
-        return r"Body of mass: {0:.2f} $M_\odot$".format(self.m/Ms)
+        return r"Body of mass: {0:.1e} $M_\odot$".format(self.m/Ms)
 
 class System(Body):
 
diff --git a/lib/plots.py b/lib/plots.py
index 320eae8..1622943 100755
--- a/lib/plots.py
+++ b/lib/plots.py
@@ -75,12 +75,12 @@ class DynamicUpdate():
 
     def on_running(self, dyn_syst, step=None, label=None):
         xdata, ydata, zdata = dyn_syst.get_positions
-        values = np.sqrt(np.sum((np.array((xdata,ydata,zdata))**2).T,axis=1))
+        values = np.sqrt(np.sum((np.array((xdata,ydata,zdata))**2).T,axis=1))/au
         self.min_x, self.max_x = -np.max([np.abs(values).max(),self.max_x]), np.max([np.abs(values).max(),self.max_x])
         self.set_lims()
         #Update data (with the new _and_ the old points)
         for i,body in enumerate(dyn_syst.bodylist):
-            x, y, z = body.q
+            x, y, z = body.q/au
             self.lines[i].set_data_3d([x], [y], [z])
         if not label is None:
             if self.blackstyle:
@@ -101,7 +101,7 @@ class DynamicUpdate():
         plt.close()
 
 
-def display_parameters(E,L,parameters,savename=""):
+def display_parameters(E,L,sma,ecc,parameters,savename=""):
     """
     """
     if savename != "":
@@ -110,15 +110,15 @@ def display_parameters(E,L,parameters,savename=""):
     bodies = ""
     for body in dyn_syst.bodylist:
         bodies += str(body)+" ; "
-    title = "Relative difference of the {0:s} "+"for a system composed of {0:s}\n integrated with {1:s} for a duration of {2:.2f} years ".format(bodies, integrator, duration/yr)
+    title1, title2 = "Relative difference of the {0:s} ","for a system composed of {0:s}\n integrated with {1:s} for a duration of {2:.2f} years ".format(bodies, integrator, duration/yr)
 
     fig1 = plt.figure(figsize=(15,7))
     ax1 = fig1.add_subplot(111)
     for i in range(len(E)):
         ax1.plot(np.arange(E[i].shape[0])*step[i]/yr, np.abs((E[i]-E[i][0])/E[i][0]), label="step of {0:.2e}yr".format(step[i]/yr))
-    ax1.set(xlabel=r"$t (yr)$", ylabel=r"$\left|\frac{\delta E_m}{E_m(t=0)}\right|$", yscale='log')
+    ax1.set(xlabel=r"$t [yr]$", ylabel=r"$\left|\frac{\delta E_m}{E_m(t=0)}\right|$", yscale='log')
     ax1.legend()
-    fig1.suptitle(title.format("mechanical energy"))
+    fig1.suptitle(title1.format("mechanical energy")+title2)
     fig1.savefig("plots/{0:s}dEm.png".format(savename),bbox_inches="tight")
 
     fig2 = plt.figure(figsize=(15,7))
@@ -127,9 +127,19 @@ def display_parameters(E,L,parameters,savename=""):
         dL = ((L[i]-L[i][0])/L[i][0])
         dL[np.isnan(dL)] = 0.
         ax2.plot(np.arange(L[i].shape[0])*step[i]/yr, np.abs(np.sum(dL,axis=1)), label="step of {0:.2e}yr".format(step[i]/yr))
-    ax2.set(xlabel=r"$t (yr)$", ylabel=r"$\left|\frac{\delta \vec{L}}{\vec{L}(t=0)}\right|$",yscale='log')
+    ax2.set(xlabel=r"$t [yr]$", ylabel=r"$\left|\frac{\delta \vec{L}}{\vec{L}(t=0)}\right|$",yscale='log')
     ax2.legend()
-    fig2.suptitle(title.format("kinetic moment"))
+    fig2.suptitle(title1.format("kinetic moment")+title2)
     fig2.savefig("plots/{0:s}dL2.png".format(savename),bbox_inches="tight")
+
+    fig3 = plt.figure(figsize=(15,7))
+    ax3 = fig3.add_subplot(111)
+    ax3.plot(np.arange(sma.shape[0])*step[i]/yr, sma, label="a (semi major axis)")
+    ax3.plot(np.arange(ecc.shape[0])*step[i]/yr, ecc, label="e (eccentricity)")
+    ax3.set(xlabel=r"$t [yr]$", ylabel=r"$a [au] or e$")
+    ax3.legend()
+    fig3.suptitle("Semi major axis and eccentricity "+title2)
+    fig3.savefig("plots/{0:s}a_e.png".format(savename),bbox_inches="tight")
+
     plt.show(block=True)
   
\ No newline at end of file
diff --git a/main.py b/main.py
index 133e5e6..b419261 100755
--- a/main.py
+++ b/main.py
@@ -27,30 +27,31 @@ def main():
     v = np.array([v1, v2, v3])
 
     #integration parameters
-    duration, step = 100*yr, np.array([1./(365.25*4.), 1./(365.25*2.), 1./365.25])*yr #integration time and step in years
+    duration, step = 100*yr, np.array([1./(365.25*2.), 1./(365.25*2.), 1./365.25])*yr #integration time and step in years
+
     integrator = "leapfrog"
     n_bodies = 2
-    display = False
+    display = True
     savename = "{0:d}bodies_{1:s}".format(n_bodies, integrator)
 
     #simulation start
     bodylist = []
     for i in range(n_bodies):
         bodylist.append(Body(m[i], q[i], v[i]))
-    dyn_syst = System(bodylist)
-    dyn_syst.COMShift()
+    bin_syst = System(bodylist[0:2])
+    dyn_syst = System(bodylist, main=True)
 
     E, L = [], []
     for step0 in step:
         if integrator.lower() in ['leapfrog', 'frogleap', 'frog']:
-            E0, L0 = leapfrog(dyn_syst,duration, step0, recover_param=True, display=display, savename=savename)
+            E0, L0, sma, ecc = leapfrog(dyn_syst, bin_syst, duration, step0, recover_param=True, display=display, savename=savename)
         elif integrator.lower() in ['hermite','herm']:
-            E0, L0 = hermite(dyn_syst, duration, step0, recover_param=True, display=display, savename=savename)
+            E0, L0, sma, ecc = hermite(dyn_syst, bin_syst, duration, step0, recover_param=True, display=display, savename=savename)
         E.append(E0)
         L.append(L0)
 
     parameters = [duration, step, dyn_syst, integrator]
-    display_parameters(E, L, parameters=parameters, savename=savename)
+    display_parameters(E, L, sma, ecc, parameters=parameters, savename=savename)
 
     return 0