diff --git a/analysis/SubPhase.cpp b/analysis/SubPhase.cpp
index c1415ebc..7b466a24 100644
--- a/analysis/SubPhase.cpp
+++ b/analysis/SubPhase.cpp
@@ -28,24 +28,24 @@ SubPhase::SubPhase(std::shared_ptr <Domain> dm):
 
 	//.........................................
 	if (Dm->rank()==0){
-		TIMELOG = fopen("subphase.csv","a+");
-		if (fseek(TIMELOG,0,SEEK_SET) == fseek(TIMELOG,0,SEEK_CUR))
+		SUBPHASE = fopen("subphase.csv","a+");
+		if (ftell(SUBPHASE) == 0)
 		{
 			// If timelog is empty, write a short header to list the averages
-			//fprintf(TIMELOG,"--------------------------------------------------------------------------------------\n");
-			fprintf(TIMELOG,"time rn rw nun nuw Fx Fy Fz iftwn ");				
-			fprintf(TIMELOG,"pwc pwd pnc pnd ");						// pressures 
-			fprintf(TIMELOG,"Mwc Mwd Mwi Mnc Mnd Mni ");				// mass 
-			fprintf(TIMELOG,"Pwc_x Pwd_x Pwi_x Pnc_x Pnd_x Pni_x ");	// momentum 
-			fprintf(TIMELOG,"Pwc_y Pwd_y Pwi_y Pnc_y Pnd_y Pni_y ");			
-			fprintf(TIMELOG,"Pwc_z Pwd_z Pwi_z Pnc_z Pnd_z Pni_z ");			
-			fprintf(TIMELOG,"Kwc Kwd Kwi Knc Knd Kni ");				// kinetic energy
-			fprintf(TIMELOG,"Vwc Awc Hwc Xwc ");					 	// wc region 
-			fprintf(TIMELOG,"Vwd Awd Hwd Xwd Nwd ");					// wd region
-			fprintf(TIMELOG,"Vnc Anc Hnc Xnc ");					 	// nc region
-			fprintf(TIMELOG,"Vnd And Hnd Xnd Nnd ");					// nd region
-			fprintf(TIMELOG,"Vi Ai Hi Xi ");					 		// interface region 
-			fprintf(TIMELOG,"Vic Aic Hic Xic Nic\n");					// interface region 
+			//fprintf(SUBPHASE,"--------------------------------------------------------------------------------------\n");
+			fprintf(SUBPHASE,"time rn rw nun nuw Fx Fy Fz iftwn ");				
+			fprintf(SUBPHASE,"pwc pwd pnc pnd ");						// pressures 
+			fprintf(SUBPHASE,"Mwc Mwd Mwi Mnc Mnd Mni ");				// mass 
+			fprintf(SUBPHASE,"Pwc_x Pwd_x Pwi_x Pnc_x Pnd_x Pni_x ");	// momentum 
+			fprintf(SUBPHASE,"Pwc_y Pwd_y Pwi_y Pnc_y Pnd_y Pni_y ");			
+			fprintf(SUBPHASE,"Pwc_z Pwd_z Pwi_z Pnc_z Pnd_z Pni_z ");			
+			fprintf(SUBPHASE,"Kwc Kwd Kwi Knc Knd Kni ");				// kinetic energy
+			fprintf(SUBPHASE,"Vwc Awc Hwc Xwc ");					 	// wc region 
+			fprintf(SUBPHASE,"Vwd Awd Hwd Xwd Nwd ");					// wd region
+			fprintf(SUBPHASE,"Vnc Anc Hnc Xnc ");					 	// nc region
+			fprintf(SUBPHASE,"Vnd And Hnd Xnd Nnd ");					// nd region
+			fprintf(SUBPHASE,"Vi Ai Hi Xi ");					 		// interface region 
+			fprintf(SUBPHASE,"Vic Aic Hic Xic Nic\n");					// interface region 
 
 			// stress tensor?
 		}
@@ -56,21 +56,21 @@ SubPhase::SubPhase(std::shared_ptr <Domain> dm):
 		sprintf(LocalRankString,"%05d",Dm->rank());
 		char LocalRankFilename[40];
 		sprintf(LocalRankFilename,"%s%s","subphase.csv.",LocalRankString);
-		TIMELOG = fopen(LocalRankFilename,"a+");
-		//fprintf(TIMELOG,"--------------------------------------------------------------------------------------\n");
-		fprintf(TIMELOG,"time rn rw nun nuw Fx Fy Fz iftwn ");				
-		fprintf(TIMELOG,"pwc pwd pnc pnd ");						// pressures 
-		fprintf(TIMELOG,"Mwc Mwd Mwi Mnc Mnd Mni ");				// mass 
-		fprintf(TIMELOG,"Pwc_x Pwd_x Pwi_x Pnc_x Pnd_x Pni_x ");	// momentum 
-		fprintf(TIMELOG,"Pwc_y Pwd_y Pwi_y Pnc_y Pnd_y Pni_y ");			
-		fprintf(TIMELOG,"Pwc_z Pwd_z Pwi_z Pnc_z Pnd_z Pni_z ");			
-		fprintf(TIMELOG,"Kwc Kwd Kwi Knc Knd Kni ");				// kinetic energy
-		fprintf(TIMELOG,"Vwc Awc Hwc Xwc ");					 	// wc region 
-		fprintf(TIMELOG,"Vwd Awd Hwd Xwd Nwd ");					// wd region
-		fprintf(TIMELOG,"Vnc Anc Hnc Xnc ");					 	// nc region
-		fprintf(TIMELOG,"Vnd And Hnd Xnd Nnd ");					// nd region
-		fprintf(TIMELOG,"Vi Ai Hi Xi ");					 		// interface region 
-		fprintf(TIMELOG,"Vic Aic Hic Xic Nic\n");					// interface region 
+		SUBPHASE = fopen(LocalRankFilename,"a+");
+		//fprintf(SUBPHASE,"--------------------------------------------------------------------------------------\n");
+		fprintf(SUBPHASE,"time rn rw nun nuw Fx Fy Fz iftwn ");				
+		fprintf(SUBPHASE,"pwc pwd pnc pnd ");						// pressures 
+		fprintf(SUBPHASE,"Mwc Mwd Mwi Mnc Mnd Mni ");				// mass 
+		fprintf(SUBPHASE,"Pwc_x Pwd_x Pwi_x Pnc_x Pnd_x Pni_x ");	// momentum 
+		fprintf(SUBPHASE,"Pwc_y Pwd_y Pwi_y Pnc_y Pnd_y Pni_y ");			
+		fprintf(SUBPHASE,"Pwc_z Pwd_z Pwi_z Pnc_z Pnd_z Pni_z ");			
+		fprintf(SUBPHASE,"Kwc Kwd Kwi Knc Knd Kni ");				// kinetic energy
+		fprintf(SUBPHASE,"Vwc Awc Hwc Xwc ");					 	// wc region 
+		fprintf(SUBPHASE,"Vwd Awd Hwd Xwd Nwd ");					// wd region
+		fprintf(SUBPHASE,"Vnc Anc Hnc Xnc ");					 	// nc region
+		fprintf(SUBPHASE,"Vnd And Hnd Xnd Nnd ");					// nd region
+		fprintf(SUBPHASE,"Vi Ai Hi Xi ");					 		// interface region 
+		fprintf(SUBPHASE,"Vic Aic Hic Xic Nic\n");					// interface region 
 	}
 }
 
@@ -78,43 +78,43 @@ SubPhase::SubPhase(std::shared_ptr <Domain> dm):
 // Destructor
 SubPhase::~SubPhase()
 {
-    if ( TIMELOG!=NULL ) { fclose(TIMELOG); }
+    if ( SUBPHASE!=NULL ) { fclose(SUBPHASE); }
 
 }
 
 void SubPhase::Write(int timestep)
 {
 	if (Dm->rank()==0){
-		fprintf(TIMELOG,"%i %.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g ",timestep,rho_n,rho_w,nu_n,nu_w,Fx,Fy,Fz,gamma_wn); 
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g ",gwc.p, gwd.p, gnc.p, gnd.p);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g ",gwc.M, gwd.M, giwn.Mw, gnc.M, gnd.M, giwn.Mn);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g ",gwc.Px, gwd.Px, giwn.Px, gnc.Px, gnd.Px, giwn.Px);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g ",gwc.Py, gwd.Py, giwn.Py, gnc.Py, gnd.Py, giwn.Py);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g ",gwc.Pz, gwd.Pz, giwn.Pz, gnc.Pz, gnd.Pz, giwn.Pz);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g ",gwc.K, gwd.K, giwn.Kw, gnc.K, gnd.K, giwn.Kn);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g ",gwc.V, gwc.A, gwc.H, gwc.X);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %i ",gwd.V, gwd.A, gwd.H, gwd.X, gwd.Nc);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g ",gnc.V, gnc.A, gnc.H, gnc.X);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %i ",gnd.V, gnd.A, gnd.H, gnd.X, gnd.Nc);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g ",giwn.V, giwn.A, giwn.H, giwn.X);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %i\n",giwnc.V, giwnc.A, giwnc.H, giwnc.X, giwnc.Nc);
-		fflush(TIMELOG);
+		fprintf(SUBPHASE,"%i %.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g ",timestep,rho_n,rho_w,nu_n,nu_w,Fx,Fy,Fz,gamma_wn); 
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g ",gwc.p, gwd.p, gnc.p, gnd.p);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %.5g %.5g ",gwc.M, gwd.M, giwn.Mw, gnc.M, gnd.M, giwn.Mn);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %.5g %.5g ",gwc.Px, gwd.Px, giwn.Px, gnc.Px, gnd.Px, giwn.Px);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %.5g %.5g ",gwc.Py, gwd.Py, giwn.Py, gnc.Py, gnd.Py, giwn.Py);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %.5g %.5g ",gwc.Pz, gwd.Pz, giwn.Pz, gnc.Pz, gnd.Pz, giwn.Pz);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %.5g %.5g ",gwc.K, gwd.K, giwn.Kw, gnc.K, gnd.K, giwn.Kn);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g ",gwc.V, gwc.A, gwc.H, gwc.X);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %i ",gwd.V, gwd.A, gwd.H, gwd.X, gwd.Nc);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g ",gnc.V, gnc.A, gnc.H, gnc.X);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %i ",gnd.V, gnd.A, gnd.H, gnd.X, gnd.Nc);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g ",giwn.V, giwn.A, giwn.H, giwn.X);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %i\n",giwnc.V, giwnc.A, giwnc.H, giwnc.X, giwnc.Nc);
+		fflush(SUBPHASE);
 	}
 	else{
-		fprintf(TIMELOG,"%i %.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g ",timestep,rho_n,rho_w,nu_n,nu_w,Fx,Fy,Fz,gamma_wn);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g ",wc.p, wd.p, nc.p, nd.p);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g ",wc.M, wd.M, iwn.Mw, nc.M, nd.M, iwn.Mn);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g ",wc.Px, wd.Px, iwn.Px, nc.Px, nd.Px, iwn.Px);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g ",wc.Py, wd.Py, iwn.Py, nc.Py, nd.Py, iwn.Py);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g ",wc.Pz, wd.Pz, iwn.Pz, nc.Pz, nd.Pz, iwn.Pz);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %.5g %.5g ",wc.K, wd.K, iwn.Kw, nc.K, nd.K, iwn.Kn);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g ",wc.V, wc.A, wc.H, wc.X);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %i ",wd.V, wd.A, wd.H, wd.X, wd.Nc);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g ",nc.V, nc.A, nc.H, nc.X);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g %i ",nd.V, nd.A, nd.H, nd.X, nd.Nc);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g ",iwn.V, iwn.A, iwn.H, iwn.X);
-		fprintf(TIMELOG,"%.5g %.5g %.5g %.5g\n",iwnc.V, iwnc.A, iwnc.H, iwnc.X);
-		fflush(TIMELOG);
+		fprintf(SUBPHASE,"%i %.5g %.5g %.5g %.5g %.5g %.5g %.5g %.5g ",timestep,rho_n,rho_w,nu_n,nu_w,Fx,Fy,Fz,gamma_wn);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g ",wc.p, wd.p, nc.p, nd.p);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %.5g %.5g ",wc.M, wd.M, iwn.Mw, nc.M, nd.M, iwn.Mn);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %.5g %.5g ",wc.Px, wd.Px, iwn.Px, nc.Px, nd.Px, iwn.Px);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %.5g %.5g ",wc.Py, wd.Py, iwn.Py, nc.Py, nd.Py, iwn.Py);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %.5g %.5g ",wc.Pz, wd.Pz, iwn.Pz, nc.Pz, nd.Pz, iwn.Pz);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %.5g %.5g ",wc.K, wd.K, iwn.Kw, nc.K, nd.K, iwn.Kn);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g ",wc.V, wc.A, wc.H, wc.X);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %i ",wd.V, wd.A, wd.H, wd.X, wd.Nc);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g ",nc.V, nc.A, nc.H, nc.X);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g %i ",nd.V, nd.A, nd.H, nd.X, nd.Nc);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g ",iwn.V, iwn.A, iwn.H, iwn.X);
+		fprintf(SUBPHASE,"%.5g %.5g %.5g %.5g\n",iwnc.V, iwnc.A, iwnc.H, iwnc.X);
+		fflush(SUBPHASE);
 	}
 
 }
diff --git a/analysis/SubPhase.h b/analysis/SubPhase.h
index b0a68347..738d2188 100644
--- a/analysis/SubPhase.h
+++ b/analysis/SubPhase.h
@@ -104,7 +104,7 @@ public:
 	
 private:
 	FILE *TIMELOG;
-
+	FILE *SUBPHASE;
 };
 
 #endif
diff --git a/analysis/morphology.cpp b/analysis/morphology.cpp
index f86e2de8..07d3a13b 100644
--- a/analysis/morphology.cpp
+++ b/analysis/morphology.cpp
@@ -754,12 +754,20 @@ double MorphGrow(DoubleArray &BoundaryDist, DoubleArray &Dist, Array<char> &id,
 		if (morph_delta / morph_delta_previous > 2.0 ) morph_delta = morph_delta_previous*2.0;
 				
 		//MAX_DISPLACEMENT *= max(TargetGrowth/GrowthEstimate,1.25);
-		if (MAX_DISPLACEMENT > 1.0 ){
-			if (morph_delta > 0.0 ) morph_delta = 1.0;
-			else 					morph_delta = -1.0;
-			
-			//if (COUNT_FOR_LOOP > 2) COUNT_FOR_LOOP = 100;
-			COUNT_FOR_LOOP = 100; // exit loop if displacement is too large
+		
+		if (morph_delta > 0.0 ){
+			// object is growing
+			if (MAX_DISPLACEMENT > 3.0 ){
+				morph_delta = 3.0;
+				COUNT_FOR_LOOP = 100; // exit loop if displacement is too large
+			}
+		}
+		else{
+			// object is shrinking
+			if (MAX_DISPLACEMENT > 1.0 ){
+				morph_delta = -1.0;
+				COUNT_FOR_LOOP = 100; // exit loop if displacement is too large
+			}
 		}
 	}
 	if (rank == 0) printf("Final delta=%f \n",morph_delta);