diff --git a/analysis/python/astar/data/read.py b/analysis/python/astar/data/read.py
index 270263f..2dd0d45 100644
--- a/analysis/python/astar/data/read.py
+++ b/analysis/python/astar/data/read.py
@@ -116,6 +116,7 @@ def parse_ts(fdir, fname):
         line[i] = line[i].replace('VTXBUF_', "")
         line[i] = line[i].replace('UU', "uu")
         line[i] = line[i].replace('_total', "tot")
+        line[i] = line[i].replace('ACCRETION', "acc")
         line[i] = line[i].replace('A', "aa")
         line[i] = line[i].replace('LNRHO', "lnrho")
         line[i] = line[i].replace('ENTROPY', "ss")
diff --git a/analysis/python/astar/visual/lineplot.py b/analysis/python/astar/visual/lineplot.py
index f33d75d..fe750ed 100644
--- a/analysis/python/astar/visual/lineplot.py
+++ b/analysis/python/astar/visual/lineplot.py
@@ -32,7 +32,7 @@ def plot_min_man_rms(ts, xaxis, yaxis1, yaxis2, yaxis3):
     plt.xlabel(xaxis)
     plt.legend()
 
-def plot_ts(ts, show_all=False, lnrho=False, uutot=False, uux=False, uuy=False, uuz=False, ss=False):
+def plot_ts(ts, show_all=False, lnrho=False, uutot=False, uux=False, uuy=False, uuz=False, ss=False, acc=False):
 
     if show_all:
         lnrho=True
@@ -41,6 +41,7 @@ def plot_ts(ts, show_all=False, lnrho=False, uutot=False, uux=False, uuy=False,
         uuy=True
         uuz=True
         ss=True
+        acc=True
 
     if lnrho:
         plt.figure()
@@ -89,5 +90,13 @@ def plot_ts(ts, show_all=False, lnrho=False, uutot=False, uux=False, uuy=False,
         yaxis2 = 'ss_min'
         yaxis3 = 'ss_max'
         plot_min_man_rms(ts, xaxis, yaxis1, yaxis2, yaxis3)
+
+    if acc:   
+        plt.figure()
+        xaxis = 't_step'
+        yaxis1 = 'acc_rms'
+        yaxis2 = 'acc_min'
+        yaxis3 = 'acc_max'
+        plot_min_man_rms(ts, xaxis, yaxis1, yaxis2, yaxis3)
   
     plt.show()
diff --git a/src/standalone/simulation.cc b/src/standalone/simulation.cc
index 911dfb8..c16bdee 100644
--- a/src/standalone/simulation.cc
+++ b/src/standalone/simulation.cc
@@ -252,17 +252,9 @@ run_simulation(void)
 #if LSINK
 
         const AcReal sum_mass = acReduceScal(RTYPE_MAX, VTXBUF_ACCRETION);
-//	if (i > 1000) {
-            accreted_mass = accreted_mass + sum_mass;
-//        } else {
-//            accreted_mass = 0.0;
-//        }
+        accreted_mass = accreted_mass + sum_mass;
         AcReal sink_mass = 0.0;
-        //if (i > 1000 ) {
-	    sink_mass = mesh_info.real_params[AC_M_sink_init] + accreted_mass;
-	//}
-        printf("sink mass is: %e \n", sink_mass); 
-        printf("accreted mass is: %e \n", accreted_mass); 
+	sink_mass = mesh_info.real_params[AC_M_sink_init] + accreted_mass;
         acLoadDeviceConstant(AC_M_sink, sink_mass);
         vertex_buffer_set(VTXBUF_ACCRETION, 0.0, mesh);
         
@@ -315,6 +307,8 @@ run_simulation(void)
             */
 
             print_diagnostics(i, dt, t_step, diag_file);
+            printf("sink mass is: %.15e \n", sink_mass); 
+            printf("accreted mass is: %.15e \n", accreted_mass); 
 
             /*
                 We would also might want an XY-average calculating funtion,