Merge pull request #152 from PhasicFlow/develop

timeInfo is updated

src/Integration/AdamsBashforth2/AdamsBashforth2.cpp

@@ -50,7 +50,7 @@ bool intAllActive(
 		"AdamsBashforth2::correct",
 		rpIntegration (activeRng.start(), activeRng.end()),
 		LAMBDA_HD(uint32 i){
-			d_y[i] = damping*(d_dy[i] + dt*(static_cast<real>(1.5) * d_dy[i] - static_cast<real>(0.5) * d_dy1[i]));
+			d_y[i] = damping*(d_y[i] + dt*(static_cast<real>(1.5) * d_dy[i] - static_cast<real>(0.5) * d_dy1[i]));
 			d_dy1[i] = d_dy[i];
 		});
 	Kokkos::fence();

src/Interaction/Models/contactForce/cGAbsoluteLinearCF.hpp

@@ -25,12 +25,6 @@ Licence:
 #include "symArrays.hpp"
 
 
-
-
-
-
-
-
 namespace pFlow::cfModels
 {
 
@@ -159,7 +153,7 @@ protected:
 		{
 
 			etha_t[i] = -2.0*log( et[i])*sqrt(kt[i]*2/7) /
-					sqrt(log(pow(et[i],2.0))+ pow(Pi,2.0));
+					sqrt(log(pow(et[i],2))+ pow(Pi,2));
 		}
 
 		Vector<linearProperties> prop("prop", nElem);
@@ -285,8 +279,8 @@ public:
 
 		history.overlap_t_ += Vt*dt;
 
-		real mi = 3*Pi/4*pow(Ri,3.0)*rho_[propId_i];
+		real mi = 3*Pi/4*pow(Ri,3)*rho_[propId_i];
-		real mj = 3*Pi/4*pow(Rj,3.0)*rho_[propId_j];
+		real mj = 3*Pi/4*pow(Rj,3)*rho_[propId_j];
 
 		real sqrt_meff = sqrt((mi*mj)/(mi+mj));
 
@@ -299,29 +293,24 @@ public:
 		else if (addDissipationModel_==3)
 		{
 			auto pie =3.14;
-			prop.en_ = exp((pow(f_,1.5)*log(prop.en_)*sqrt( (1-((pow(log(prop.en_),2))/(pow(log(prop.en_),2)+pow(pie,2))))/(1-(pow(f_,3)*(pow(log(prop.en_),2))/(pow(log(prop.en_),2)+pow(pie,2)))) ) ));
+			prop.en_ = exp((pow(f_,static_cast<real>(1.5))*log(prop.en_)*sqrt( (1-((pow(log(prop.en_),2))/(pow(log(prop.en_),2)+pow(pie,2))))/(1-(pow(f_,3)*(pow(log(prop.en_),2))/(pow(log(prop.en_),2)+pow(pie,2)))) ) ));
 		}
 
-
-
-
 		real ethan_  = -2.0*log(prop.en_)*sqrt(prop.kn_)/
-					sqrt(pow(log(prop.en_),2.0)+ pow(Pi,2.0));
+					sqrt(pow(log(prop.en_),2)+ pow(Pi,2));
 
 			//REPORT(0)<<"\n en n is  :  "<<END_REPORT;
 			//REPORT(0)<< prop.en_ <<END_REPORT;
 
-		FCn = (  -pow(f_,3.0)*prop.kn_ * ovrlp_n - sqrt_meff * pow(f_,1.5) * ethan_ * vrn)*Nij;
+		FCn = (  -pow(f_,3)*prop.kn_ * ovrlp_n - sqrt_meff * pow(f_,static_cast<real>(1.5)) * ethan_ * vrn)*Nij;
-		FCt = ( -pow(f_,3.0)*prop.kt_ * history.overlap_t_ - sqrt_meff * pow(f_,1.5) *  prop.ethat_*Vt);
+		FCt = ( -pow(f_,3)*prop.kt_ * history.overlap_t_ - sqrt_meff * pow(f_,static_cast<real>(1.5)) *  prop.ethat_*Vt);
-		
+				
-
-		
 		real ft = length(FCt);
 		real ft_fric = prop.mu_ * length(FCn);
 		
 		if(ft > ft_fric)
 		{
-			if( length(history.overlap_t_) >static_cast<real>(0.0))
+			if( length(history.overlap_t_) >zero)
 			{
 				if constexpr (limited)
 				{

src/Interaction/Models/contactForce/cGNonLinearCF.hpp

@@ -236,7 +236,7 @@ public:
 
 		const auto prop = nonlinearProperties_(propId_i,propId_j);
 
-		const real f = 2.0/( 1.0/cGFi + 1.0/cGFj );
+		const real f = 2/( 1/cGFi + 1/cGFj );
 
 		real vrn = dot(Vr, Nij);	
 		realx3 Vt = Vr - vrn*Nij;
@@ -245,7 +245,7 @@ public:
 
 		real mi = 3*Pi/4*pow(Ri,static_cast<real>(3))*rho_[propId_i];
 		real mj = 3*Pi/4*pow(Rj,static_cast<real>(3))*rho_[propId_j];
-		real Reff = 1.0/(1/Ri + 1/Rj);
+		real Reff = 1/(1/Ri + 1/Rj);
 
 		real K_hertz = 4.0/3.0*prop.Yeff_*sqrt(Reff);
 		real sqrt_meff_K_hertz = sqrt((mi*mj)/(mi+mj) * K_hertz);
@@ -258,13 +258,13 @@ public:
 		else if (addDissipationModel_==3)
 		{
 			
-			en = exp((pow(f,1.5)*log(prop.en_)*sqrt( (1-((pow(log(prop.en_),2))/(pow(log(prop.en_),2)+pow(Pi,2))))/(1-(pow(f,3)*(pow(log(prop.en_),2))/(pow(log(prop.en_),2)+pow(Pi,2)))) ) ));
+			en = exp((pow(f,static_cast<real>(1.5))*log(prop.en_)*sqrt( (1-((pow(log(prop.en_),2))/(pow(log(prop.en_),2)+pow(Pi,2))))/(1-(pow(f,3)*(pow(log(prop.en_),2))/(pow(log(prop.en_),2)+pow(Pi,2)))) ) ));
 		}
 		
 		real Kn = static_cast<real>(4.0/3.0) * prop.Yeff_ * sqrt(Reff*ovrlp_n);
 
 		real ethan_  = -2.0*log(en)*sqrt(Kn)/
-			sqrt(pow(log(en),2.0)+ pow(Pi,2.0));
+			sqrt(pow(log(en),2)+ pow(Pi,2));
 
 		FCn = ( - Kn*ovrlp_n - 
 			   sqrt_meff_K_hertz*ethan_*pow(ovrlp_n,static_cast<real>(0.25))*vrn)*Nij;

src/Interaction/Models/contactForce/cGRelativeLinearCF.hpp

@@ -255,8 +255,8 @@ public:
 
 		history.overlap_t_ += Vt*dt;
 
-		real mi = 3*Pi/4*pow(Ri,3.0)*rho_[propId_i];
+		real mi = 3*Pi/4*pow(Ri,3)*rho_[propId_i];
-		real mj = 3*Pi/4*pow(Rj,3.0)*rho_[propId_j];
+		real mj = 3*Pi/4*pow(Rj,3)*rho_[propId_j];
 
 		real sqrt_meff = sqrt((mi*mj)/(mi+mj));
 		
@@ -268,14 +268,26 @@ public:
 		else if (addDissipationModel_==3)
 		{
 			
-			en = exp((pow(f_,1.5)*log(prop.en_)*sqrt( (1-((pow(log(prop.en_),2))/(pow(log(prop.en_),2)+pow(Pi,2))))/(1-(pow(f_,3)*(pow(log(prop.en_),2))/(pow(log(prop.en_),2)+pow(Pi,2)))) ) ));
+			en = exp(
+				(
+					pow(f_,static_cast<real>(1.5))*
+					log(prop.en_)*
+					sqrt( 
+						(1-((pow(log(prop.en_),static_cast<real>(2.0))
+						)/
+						(
+							pow(log(prop.en_),static_cast<real>(2.0))+
+							pow(Pi,static_cast<real>(2.0)))))/
+							(1-(pow(f_,3)*(pow(log(prop.en_),2))/
+							(pow(log(prop.en_),static_cast<real>(2.0))+
+							pow(Pi,static_cast<real>(2.0))))) ) ));
 		}
 
 		real ethan_  = -2.0*log(en)*sqrt(prop.kn_)/
-			sqrt(pow(log(en),2.0)+ pow(Pi,2.0));
+			sqrt(pow(log(en),2)+ pow(Pi,2));
 		
 
-		FCn = ( -f_*prop.kn_ * ovrlp_n - sqrt_meff * pow(f_,0.5) * ethan_ * vrn)*Nij;
+		FCn = ( -f_*prop.kn_ * ovrlp_n - sqrt_meff * pow(f_,half) * ethan_ * vrn)*Nij;
 		FCt = ( -f_*prop.kt_ * history.overlap_t_);
 		
 

src/Interaction/Models/contactForce/linearCF.hpp

@@ -139,10 +139,10 @@ protected:
 		ForAll(i , kn)
 		{
 			etha_n[i] = -2.0*log(en[i])*sqrt(kn[i])/
-					sqrt(pow(log(en[i]),2.0)+ pow(Pi,2.0));
+					sqrt(pow(log(en[i]),static_cast<real>(2.0))+ pow(Pi,static_cast<real>(2.0)));
 
 			etha_t[i] = -2.0*log( et[i]*sqrt(kt[i]) )/
-					sqrt(pow(log(et[i]),2.0)+ pow(Pi,2.0));
+					sqrt(pow(log(et[i]),static_cast<real>(2.0))+ pow(Pi,static_cast<real>(2.0)));
 		}
 
 		Vector<linearProperties> prop("prop", nElem);
@@ -243,8 +243,8 @@ public:
 
 		history.overlap_t_ += Vt*dt;
 
-		real mi = 3*Pi/4*pow(Ri,3.0)*rho_[propId_i];
+		real mi = 3*Pi/4*pow(Ri,static_cast<real>(3.0))*rho_[propId_i];
-		real mj = 3*Pi/4*pow(Rj,3.0)*rho_[propId_j];
+		real mj = 3*Pi/4*pow(Rj,static_cast<real>(3.0))*rho_[propId_j];
 
 		real sqrt_meff = sqrt((mi*mj)/(mi+mj));
 

src/Interaction/Models/contactForce/nonLinearCF.hpp

@@ -131,7 +131,7 @@ protected:
 			// we take out sqrt(meff*K_hertz) here and then consider this term 
 			// when calculating damping part. 
 			etha_n[i] = -2.2664*log(en[i])/
-					sqrt(pow(log(en[i]),2.0)+ pow(Pi,2.0));
+					sqrt(pow(log(en[i]),static_cast<real>(2.0))+ pow(Pi,static_cast<real>(2.0)));
 
 			// no damping for tangential part 
 
@@ -255,7 +255,7 @@ public:
 		// apply friction 
 		if(ft > ft_fric)
 		{
-			if( length(history.overlap_t_) >0.0)
+			if( length(history.overlap_t_) >zero)
 			{
 				if constexpr (limited)
 				{

src/Interaction/Models/rolling/grainRolling.hpp

@@ -96,10 +96,10 @@ public:
 		realx3 w_hat = wi-wj;
 		real w_hat_mag = length(w_hat);
 
-		if( !equal(w_hat_mag,0.0) )
+		if( !equal(w_hat_mag,zero) )
 			w_hat /= w_hat_mag;
 		else
-			w_hat = 0.0;
+			w_hat = zero;
 
 		auto Reff = (Ri*Rj)/(Ri+Rj);
 

src/Interaction/Models/rolling/normalRolling.hpp

@@ -94,10 +94,10 @@ public:
 		realx3 w_hat = wi-wj;
 		real w_hat_mag = length(w_hat);
 
-		if( !equal(w_hat_mag,0.0) )
+		if( !equal(w_hat_mag,zero) )
 			w_hat /= w_hat_mag;
 		else
-			w_hat = 0.0;
+			w_hat = zero;
 
 		auto Reff = (Ri*Rj)/(Ri+Rj);
 

src/Interaction/contactSearch/boundaries/periodicBoundaryContactSearch/wallBoundaryContactSearch.cpp

@@ -11,7 +11,7 @@ pFlow::wallBoundaryContactSearch::wallBoundaryContactSearch
     const ViewType1D<realx3, memory_space> &normals
 )
 :
-    cellExtent_( max(cellExtent, 0.5 ) ),
+    cellExtent_( max(cellExtent, half ) ),
     numElements_(numElements),
     numPoints_(numPoints),
     vertices_(vertices),

src/Interaction/contactSearch/contactSearch/contactSearch.cpp

@@ -31,6 +31,7 @@ pFlow::contactSearch::contactSearch(
  	Timers& timers)
 :
 	extendedDomainBox_(extDomain),
+	updateInterval_(dict.getValMax<uint32>("updateInterval", 1u)),
 	particles_(prtcl),
 	geometry_(geom),
 	bTimer_("Boundary particles contact search", &timers),
@@ -92,6 +93,7 @@ bool pFlow::contactSearch::boundaryBroadSearch
 		bTimer_.start();
 		for(uint32 i=0u; i<6u; i++)
 		{
+			output<<" boundarySearch "<< i <<" for iter "<< ti.iter()<<endl;
 			if(!BoundaryBroadSearch(
 				i,
 				ti,

src/Interaction/contactSearch/methods/cellBased/NBS/NBS.cpp

@@ -44,8 +44,8 @@ pFlow::NBS::NBS
         position,
         flags,
         diam),
-    sizeRatio_(max(dict.getVal<real>("sizeRatio"), 1.0)),
+    sizeRatio_(max(dict.getVal<real>("sizeRatio"), one)),
-    cellExtent_(max(dict.getVal<real>("cellExtent"), 0.5)),
+    cellExtent_(max(dict.getVal<real>("cellExtent"), half)),
     adjustableBox_(dict.getVal<Logical>("adjustableBox")),
     NBSLevel0_
 	(

src/Interaction/contactSearch/methods/cellBased/NBS/cellsWallLevel0.cpp

@@ -31,7 +31,7 @@ pFlow::cellsWallLevel0::cellsWallLevel0
 	const ViewType1D<realx3, memory_space>& normals
 )
 :
-    cellExtent_( max(cellExtent, 0.5 ) ),
+    cellExtent_( max(cellExtent, half ) ),
     numElements_(numElements),
     numPoints_(numPoints),
     vertices_(vertices),

src/MotionModel/entities/rotatingAxis/rotatingAxis.cpp

@@ -47,7 +47,7 @@ pFlow::rotatingAxis::rotatingAxis
 	timeInterval(),
 	line(p1, p2),
 	omega_(omega),
-	rotating_(!equal(omega,0.0))
+	rotating_(!equal(omega,zero))
 {
 
 }
@@ -58,7 +58,7 @@ pFlow::real pFlow::rotatingAxis::setOmega(real omega)
 {
 	auto tmp = omega_;
 	omega_ = omega;
-	rotating_ = !equal(omega, 0.0);
+	rotating_ = !equal(omega, zero);
 	return tmp;
 }
 
@@ -72,7 +72,7 @@ bool pFlow::rotatingAxis::read
 	if(!timeInterval::read(dict))return false;
 	if(!line::read(dict)) return false;
 		
-	real omega = dict.getValOrSet("omega", 0.0);
+	real omega = dict.getValOrSet("omega", zero);
 	
 	setOmega(omega);
 	return true;

src/Particles/globalDamping/globalDamping.cpp

@@ -23,7 +23,8 @@ Licence:
 
 
 pFlow::globalDamping::globalDamping(const systemControl& control)
-:
+{}
+/*:
     timeControl_(control.settingsDict().subDict("globalDamping"), control.time().dt(), "damping")
 {
     const dictionary& dict = control.settingsDict().subDict("globalDamping");
@@ -40,7 +41,7 @@ pFlow::globalDamping::globalDamping(const systemControl& control)
 	else
 		REPORT(2)<<"Global damping "<<Yellow_Text("is not active")<<"."<<END_REPORT;
 
-}
+}*/
 
 
 /*void pFlow::globalDamping::applyDamping
@@ -68,7 +69,8 @@ pFlow::globalDamping::globalDamping(const systemControl& control)
 
 pFlow::real pFlow::globalDamping::dampingFactor(const timeInfo& ti)const
 {
-	if(!performDamping_) return 1.0;
+  return 1;
+	/*if(!performDamping_) return 1.0;
 	if(!timeControl_.timeEvent(ti.iter(), ti.t(), ti.dt()) )return 1.0;
-	return dampingFactor_;
+	return dampingFactor_;*/
 }

src/Particles/globalDamping/globalDamping.hpp

@@ -37,7 +37,7 @@ private:
 
     real            dampingFactor_;
 
-    baseTimeControl timeControl_;
+  //  baseTimeControl timeControl_;
 
 
 public:

src/phasicFlow/containers/symArrayHD/symArrayHD.hpp

@@ -239,7 +239,7 @@ public:
 		{
 			real nr = 0.5*(sqrt(8.0*nElem+1)-1);
 			n = static_cast<uint32>(nr);
-			if( equal(nr-static_cast<real>(n) , 0.0) ) return true;
+			if( equal(nr-static_cast<real>(n) , zero) ) return true;
 
 			return false;
 		}

src/phasicFlow/ranges/stridedRange.hpp

@@ -111,7 +111,7 @@ public:
 	{
 		
 		if(!isInRange(val)) return false;
-		if(T dist = val-begin_; abs(dist%stride_)<= epsilon) return true;
+		if(const T dist = val-begin_; abs(dist%stride_)<=epsilon) return true;
 		if(equal(val,begin_))return true;
 		if(equal(val,end_))return true;
 		return false;
@@ -151,7 +151,7 @@ public:
 
 template<>
 inline
-bool stridedRange<real>::isMember(real val, real epsilon)const
+bool stridedRange<float>::isMember(float val, float epsilon)const
 	{
 
 	if(!isInRange(val)) return false;
@@ -162,8 +162,31 @@ bool stridedRange<real>::isMember(real val, real epsilon)const
 	if(equal(val,begin_))return true;
 	if(equal(val,end_))return true;
 	return false;
+	/*if(!isInRange(val)) return false;
+	if(const float dist = val-begin_; abs(remainder(dist,stride_)<= epsilon)) return true;	
+	if(equal(val,begin_))return true;
+	if(equal(val,end_))return true;
+	return false;*/
 }
 
+template<>
+inline
+bool stridedRange<double>::isMember(double val, double epsilon)const
+{
+	/*if(!isInRange(val)) return false;
+	if(const double dist = val-begin_; abs(remainder(dist,stride_)<= epsilon)) return true;	
+	if(equal(val,begin_))return true;
+	if(equal(val,end_))return true;
+	return false;*/
+	if(!isInRange(val)) return false;
+	real dist = val-begin_;
+	if(abs(
+		(dist-(static_cast<uint64>((dist+0.01*epsilon)/stride_)*stride_))
+		)<= epsilon) return true;
+	if(equal(val,begin_))return true;
+	if(equal(val,end_))return true;
+	return false;
+}
 
 }
 

src/phasicFlow/repository/Time/timeControl.cpp

@@ -25,7 +25,7 @@ Licence:
 
 bool pFlow::timeControl::screenReport()const
 {
-	return screenReportInterval_.isMember(currentIter_);
+	return screenReportInterval_.isMember(ti_.currentIter());
 }
 
 pFlow::timeControl::timeControl
@@ -33,94 +33,72 @@ pFlow::timeControl::timeControl
 	const dictionary& dict
 )
 :
-	dt_
+	ti_(dict),
-	(
-		dict.getVal<real>("dt")
-	),
 	startTime_
 	(
-		dict.getVal<real>("startTime")
+		dict.getVal<timeValue>("startTime")
 	),
 	endTime_
 	(
-		dict.getVal<real>("endTime")
+		dict.getVal<timeValue>("endTime")
 	),
 	stopAt_(endTime_),
-	currentTime_(startTime_),
 	saveInterval_
 	(
-		dict.getVal<real>("saveInterval")
+		dict.getVal<timeValue>("saveInterval")
 	),
 	lastSaved_(startTime_),
-	currentIter_(0),
-	timePrecision_
-	(
-		dict.getValOrSet("timePrecision", 4)
-	),
-	timersReportInterval_
-	(
-		startTime_,
-		dict.getValOrSet("timersReportInterval", 0.04)
-	),
 	performSorting_
 	(
 		dict.getValOrSet("performSorting", Logical("No"))
-	),
-	sortingInterval_
-	(
-		startTime_,
-		dict.getValOrSet("sortingInterval", static_cast<real>(1.0))
 	)
 {
+	timersReportInterval_ = timeStridedRange(
+		startTime_, 
+		dict.getValOrSet("timersReportInterval", 0.04));
+
+	sortingInterval_ = timeStridedRange(
+		startTime_, 
+		dict.getValOrSet("sortingInterval", 1.0));
 	checkForOutputToFile();
 }
 
 pFlow::timeControl::timeControl(
 		dictionary& dict,
-		real startTime,
+		timeValue startTime,
-		real endTime,
+		timeValue endTime,
-		real saveInterval,
+		timeValue saveInterval,
 		word startTimeName)
 :
-	dt_
+	ti_(startTime, dict),
-	(
-		dict.getVal<real>("dt")
-	),
 	startTime_(startTime),
 	endTime_(endTime),
 	stopAt_(endTime_),
-	currentTime_(startTime_),
 	saveInterval_(saveInterval),
-	lastSaved_(startTime_),
+	lastSaved_(startTime_),	
-	currentIter_(0),
-	timePrecision_
-	(
-		dict.getValOrSet("timePrecision", 4)
-	),
 	managedExternaly_(true),
-	timeName_(startTimeName),
-	timersReportInterval_
-	(
-		startTime_,
-		dict.getValOrSet("timersReportInterval", 0.04)
-	),
 	performSorting_
 	(
 		dict.getValOrSet("performSorting", Logical("No"))
-	),
-	sortingInterval_
-	(
-		startTime_,
-		dict.getValOrSet("sortingInterval", static_cast<real>(1.0))
 	)	
 {
+	timeName_ = startTimeName;
+
+	timersReportInterval_ = timeStridedRange(
+		startTime_, 
+		dict.getValOrSet("timersReportInterval", 0.04));
+
+	sortingInterval_ = timeStridedRange(
+		startTime_, 
+		dict.getValOrSet("sortingInterval", 1.0));
+
 	checkForOutputToFile();
 }
 
-pFlow::real pFlow::timeControl::setTime(real t)
+pFlow::timeValue pFlow::timeControl::setTime(timeValue t)
 {
-	real tmp = currentTime_;
+	timeValue tmp = ti_.currentTime();
-	currentTime_ = t;
+	ti_.currentTime_ = t;
 	lastSaved_ = t;
 	checkForOutputToFile();
 	return tmp;
@@ -136,15 +114,15 @@ pFlow::word pFlow::timeControl::timeName()const
 
 bool pFlow::timeControl::finalTime()const
 {
-	if( currentTime_ >= endTime_ ) return true;
+	if( ti_.currentTime_ >= endTime_ ) return true;
-	if( std::abs(currentTime_-endTime_) < 0.5*dt_ )return true;
+	if( std::abs(ti_.currentTime_-endTime_) < 0.5*ti_.dt_ )return true;
 	return false;	
 }
 
 bool pFlow::timeControl::reachedStopAt()const
 {
-	if( currentTime_ >= stopAt_ ) return true;
+	if( ti_.currentTime_ >= stopAt_ ) return true;
-	if( std::abs(currentTime_-stopAt_) < 0.5*dt_ )return true;
+	if( std::abs(ti_.currentTime_-stopAt_) < 0.5*ti_.dt_ )return true;
 	return false;
 }
 
@@ -154,22 +132,22 @@ void pFlow::timeControl::checkForOutputToFile()
 	bool save = false;
 	if(managedExternaly_)
 	{
-		if( std::abs(currentTime_-writeTime_) < 0.5*dt_)
+		if( std::abs(ti_.currentTime_-writeTime_) < 0.5*ti_.dt_)
 		{
 			save = true;
-			lastSaved_ = currentTime_;
+			lastSaved_ = ti_.currentTime_;
 		}
 	}
 	else
 	{
-		if ( std::abs(currentTime_ - lastSaved_ - saveInterval_) < 0.5 * dt_ )
+		if ( std::abs(ti_.currentTime_ - lastSaved_ - saveInterval_) < 0.5 * ti_.dt_ )
 		{	
-			lastSaved_ = currentTime_;
+			lastSaved_ = ti_.currentTime_;
 			save = true;
 		}
-		else if( std::abs(currentTime_ - lastSaved_) < std::min( pow(10.0,-1.0*timePrecision_), static_cast<real>(0.5 *dt_)) )
+		else if( std::abs(ti_.currentTime_ - lastSaved_) < std::min( std::pow(10.0,-1.0*ti_.precision()), 0.5 *ti_.dt_) )
 		{
-			lastSaved_ = currentTime_;
+			lastSaved_ = ti_.currentTime_;
 			save = true;
 		}
 		
@@ -180,13 +158,13 @@ void pFlow::timeControl::checkForOutputToFile()
 
 bool pFlow::timeControl::timersReportTime()const
 {
-	if(currentIter_<=1)return false;
+	if(ti_.currentIter_<=1)return false;
-	return timersReportInterval_.isMember(currentTime_, 0.55*dt_);
+	return timersReportInterval_.isMember(ti_.currentTime_, 0.55*ti_.dt_);
 }
 
 bool pFlow::timeControl::sortTime()const
 {
-	return performSorting_()&&sortingInterval_.isMember(currentTime_,dt_);
+	return performSorting_()&&sortingInterval_.isMember(ti_.currentTime_,ti_.dt_);
 }
 
 void pFlow::timeControl::setSaveTimeFolder(
@@ -204,10 +182,10 @@ bool pFlow::timeControl::operator ++(int)
 {
 
 	if( reachedStopAt() ) return false;
+	
 	// increament iteration number 
-	currentIter_++;
+	ti_.martchDT();
 
-	currentTime_ += dt_;
 	if(screenReport() && !managedExternaly_)
 	{
 		REPORT(0)<<"Time (s): "<<Cyan_Text( currentTimeWord() )<<END_REPORT;

src/phasicFlow/repository/Time/timeControl.hpp

@@ -36,57 +36,55 @@ class dictionary;
 
 class timeControl
 {
+public:
+	using timeStridedRange = stridedRange<timeValue>;
 private:
 
 
 	//// - Data members
 
-	// - integration time step
+		// - time information 
-	real dt_;
+		timeInfo 		ti_;
-	
+		
-	// - start time of simulation
+		// - start time of simulation
-	real startTime_;
+		timeValue startTime_;
 
-	// - end time of simulation 
+		// - end time of simulation 
-	real endTime_;
+		timeValue endTime_;
 
-	// - stopAt
+		// - stopAt
-	real stopAt_;
+		timeValue stopAt_;	
 
-	// - current time of simulation
+		// - time interval for time folder output 
-	real currentTime_;	
+		timeValue saveInterval_;
 
-	// - time interval for time folder output 
+		// - the last time folder that was saved
-	real saveInterval_;
+		timeValue lastSaved_;
 
-	// - the last time folder that was saved
+		bool  managedExternaly_ = false;
-	real lastSaved_;
 
-	// - current iteration number (for this execution)
+		bool outputToFile_ = false;
-	int32 currentIter_;
 
-	// - time precision for output time folders
+		Logical						performSorting_;
-	int32 timePrecision_;
 
-	bool  managedExternaly_ = false;
+		static 
+		inline int32StridedRagne 	screenReportInterval_ ={0,100};
 
-	word  timeName_ = "wrongSettings"; // for managedExternamly
+		static
+		inline timeStridedRange  	timersReportInterval_{0,100};
 
-	real 	writeTime_ = 0; // for managedExternamly
+		static
+		inline timeStridedRange 	sortingInterval_{0,100};
+		
+		static 
+		inline word  timeName_ = "wrongSettings"; // for managedExternamly
 
-	realStridedRange  	timersReportInterval_;
+		static
+		inline timeValue 	writeTime_ = 0; // for managedExternamly
 
-	Logical				performSorting_;
+		void checkForOutputToFile();
-
+		
-	realStridedRange 	sortingInterval_;
+		bool screenReport()const;
-
-	int32StridedRagne 	screenReportInterval_ ={0,100};
-
-	bool outputToFile_ = false;
-
-	void checkForOutputToFile();
-	
-	bool screenReport()const;
 	
 public:
 
@@ -94,22 +92,22 @@ public:
 
 	timeControl(
 		dictionary& dict,
-		real startTime,
+		timeValue startTime,
-		real endTime,
+		timeValue endTime,
-		real saveInterval,
+		timeValue saveInterval,
 		word startTimeName);
 	
 	virtual ~timeControl() = default;
 	
-	real dt()const
+	timeValue dt()const
 	{
-		return dt_;
+		return ti_.dt();
 	}
 
-	real setTime(real t);
+	timeValue setTime(timeValue t);
 	
 
-	void setStopAt(real sT)
+	void setStopAt(timeValue sT)
 	{
 		if(managedExternaly_)
 		{
@@ -117,37 +115,27 @@ public:
 		}
 	}
 
-	real startTime()const
+	timeValue startTime()const
 	{
 		return startTime_;
 	}
 
 	word timeName()const;	
 
-	real currentTime() const 
+	timeValue currentTime() const 
 	{
-		return currentTime_;
+		return ti_.currentTime();
 	}
 
 	word currentTimeWord(bool forSave = true)const
 	{
-		return real2FixedStripZeros( currentTime(), timePrecision());
+		return ti_.timeName();
-		/*if(forSave)
+		
-		{
-			if(!managedExternaly_)
-				
-			else
-				return timeName_;
-		}
-		else
-		{
-			return real2FixedStripZeros( currentTime(), timePrecision());
-		}*/
 	}
 
 	int32 currentIter()const
 	{
-		return currentIter_;
+		return ti_.currentIter();
 	}
 
 	bool finalTime()const;
@@ -179,17 +167,17 @@ public:
 		bool saveToFile, 
 		const word& timeName = "wrongTimeFolder");
 	
-	int32 timePrecision()const
-	{
-		return timePrecision_;
-	}
-
 	inline
-	timeInfo TimeInfo()const
+	const timeInfo& TimeInfo()const
 	{
-		return {static_cast<uint32>(currentIter_), currentTime_, dt_};
+		return ti_;
 	}
 
+	static
+	int32 timePrecision()
+	{
+		return timeInfo::precision();
+	}
 
 };
 

src/phasicFlow/repository/Time/timeInfo.hpp

@@ -21,6 +21,7 @@ Licence:
 #define __timeInfo__hpp_
 
 #include "types.hpp"
+#include "dictionary.hpp"
 
 namespace pFlow
 {
@@ -28,34 +29,88 @@ namespace pFlow
 class timeInfo
 {
 private:
+	
+	friend class timeControl;
 
-	uint32 currentIter_;
+	// - current iteration number (for this execution)
+	int32 		currentIter_;
 
-	real   currentTime_;
+	// - current time of simulation
+	timeValue   currentTime_;
 
-	real   dt_;
+	// - integration time step
+	timeValue   dt_;
+
+	inline static int32 presicion_ = 5;
 
 public:
 
-	timeInfo(uint32 cIter, real cTime, real dt)
+	timeInfo(int32 cIter, timeValue cTime, timeValue dt)
 	  : currentIter_(cIter),
 	    currentTime_(cTime),
 	    dt_(dt)
 	{
 	}
 
-	inline const real& t() const
+	timeInfo(const dictionary& dict)
+	:
+		currentIter_(0),
+		currentTime_(dict.getVal<timeValue>("startTime")),
+		dt_( dict.getVal<timeValue>("dt"))
+	{
+		presicion_ = dict.getValOrSet<int32>("timePrecision",5);
+	}
+
+	timeInfo(timeValue currentTime, const dictionary& dict)
+	:
+		currentIter_(0),
+		currentTime_(currentTime),
+		dt_( dict.getVal<timeValue>("dt"))
+	{
+		presicion_ = dict.getValOrSet<int32>("timePrecision",5);
+	}
+
+	inline const timeValue& currentTime()const
 	{
 		return currentTime_;
 	}
-	inline const real& dt() const
+
+	inline const timeValue& t() const
+	{
+		return currentTime_;
+	}
+	inline const timeValue& dt() const
 	{
 		return dt_;
 	}
-	inline const uint32& iter() const
+	inline const int32& iter() const
 	{
 		return currentIter_;
 	}
+
+	inline const int32& currentIter() const
+	{
+		return currentIter_;
+	}
+
+	inline
+	void martchDT()
+	{
+		currentIter_++;
+		currentTime_ += dt_;
+	}
+	
+	inline
+	word timeName()const
+	{
+		return real2FixedStripZeros(currentTime_, presicion_);
+	}
+
+	static
+	int32 precision()
+	{
+		return presicion_;
+	}
 };
 
 } // namespace pFlow

src/phasicFlow/structuredData/infinitePlane/infinitePlane.cpp

@@ -29,7 +29,7 @@ pFlow::infinitePlane::infinitePlane
 {
 	auto ln = cross(p2-p1, p3-p1);
 	
-	if( equal(ln.length(),0.0) )
+	if( equal(ln.length(),zero) )
 	{
 		fatalErrorInFunction<<
 		"invalid input to form a infinte wall "<< realx3x3(p1,p2,p3)<<endl;
@@ -60,5 +60,5 @@ bool pFlow::infinitePlane::validPlane3
   const realx3& p3
 )
 {
-  return !equal(cross(p2-p1, p3-p1).length(), 0.0);
+  return !equal(cross(p2-p1, p3-p1).length(), zero);
 }

src/phasicFlow/structuredData/infinitePlane/infinitePlane.hpp

@@ -22,7 +22,6 @@ Licence:
 #define __infinitePlane_hpp__
 
 #include "types.hpp"
-//#include "dictionary.hpp"
 #include "iIstream.hpp"
 #include "iOstream.hpp"
 
@@ -113,20 +112,20 @@ public:
 		bool inPositiveDistance(const realx3& p, real dist)const
 		{	
 			real d = pointFromPlane(p);
-			return d >= 0.0 && d <= dist;
+			return d >= zero && d <= dist;
 		}
 
 		INLINE_FUNCTION_HD
 		bool inNegativeDistance(const realx3& p, real dist)const
 		{	
 			real d = pointFromPlane(p);
-			return d < 0.0 && d >= -dist;
+			return d < zero && d >= -dist;
 		}
 
 		INLINE_FUNCTION_HD
 		bool pointOnPlane(const realx3& p)const
 		{
-			return equal(pointFromPlane(p),0.0);
+			return equal(pointFromPlane(p),zero);
 		}
 
 		INLINE_FUNCTION_HD

src/phasicFlow/triSurface/triangleFunctions.hpp

@@ -38,8 +38,8 @@ INLINE_FUNCTION_HD
 realx3 normal(const realx3& p1, const realx3& p2, const realx3& p3)
 {
 	auto n = cross(p2-p1, p3-p1);
-	if( equal(n.length(), 0.0) )
+	if( equal(n.length(), zero) )
-    	return realx3(0);
+    	return zero3;
 	else
 		return normalize(n);
 }
@@ -52,11 +52,9 @@ bool valid
 	const realx3& p3
 )
 {
-  	return !equal(cross(p2-p1, p3-p1).length(), 0.0);
+  	return !equal(cross(p2-p1, p3-p1).length(), zero);
 }
 
-
-
 }
 
 #endif 

src/phasicFlow/types/basicTypes/bTypesFunctions.cpp

@@ -81,7 +81,7 @@ pFlow::isNo(const word& str)
 }
 
 pFlow::word
-pFlow::real2Fixed(const real& v, int32 numPrecision)
+pFlow::real2Fixed(const float& v, int32 numPrecision)
 {
 	std::stringstream ss;
 
@@ -90,7 +90,32 @@ pFlow::real2Fixed(const real& v, int32 numPrecision)
 }
 
 pFlow::word
-pFlow::real2Word(const real& v, int32 numPrecision)
+pFlow::real2Fixed(const double& v, int32 numPrecision)
+{
+	std::stringstream ss;
+
+	ss << std::fixed << std::setprecision(numPrecision) << v;
+	return ss.str();
+}
+
+pFlow::word
+pFlow::real2Word(const float& v, int32 numPrecision)
+{
+	std::stringstream ss;
+	if (abs(v) < verySmallValue)
+	{
+		ss << "0";
+	}
+	else
+	{
+		ss << std::setprecision(numPrecision) << v;
+	}
+
+	return ss.str();
+}
+
+pFlow::word
+pFlow::real2Word(const double& v, int32 numPrecision)
 {
 	std::stringstream ss;
 	if (abs(v) < verySmallValue)
@@ -145,7 +170,14 @@ pFlow::removeDecimalZeros(const word& str)
 }
 
 pFlow::word
-pFlow::real2FixedStripZeros(const real& v, int32 numPrecision)
+pFlow::real2FixedStripZeros(const float& v, int32 numPrecision)
+{
+	word strVal = real2Fixed(v, numPrecision);
+	return removeDecimalZeros(strVal);
+}
+
+pFlow::word
+pFlow::real2FixedStripZeros(const double& v, int32 numPrecision)
 {
 	word strVal = real2Fixed(v, numPrecision);
 	return removeDecimalZeros(strVal);

src/phasicFlow/types/basicTypes/bTypesFunctions.hpp

@@ -39,6 +39,8 @@ inline const real  zero = 0.0;
 /// one real variable
 inline const real  one = 1.0;
 
+inline const real half = 0.5;
+
 /// zero int32 variable
 inline const int32 zero32 = 0;
 
@@ -76,11 +78,17 @@ isNo(const word& str);
 
 /// Convert floating point variable to string with fixed number of precisions
 word
-real2Fixed(const real& v, int32 numPrecision = 6);
+real2Fixed(const float& v, int32 numPrecision = 6);
+
+word
+real2Fixed(const double& v, int32 numPrecision = 6);
 
 /// Convert floating point variable to string with general format
 word
-real2Word(const real& v, int32 numPrecision = 6);
+real2Word(const float& v, int32 numPrecision = 6);
+
+word
+real2Word(const double& v, int32 numPrecision = 6);
 
 /// Convert int32 to word
 word
@@ -92,7 +100,10 @@ removeDecimalZeros(const word& str);
 
 /// Convert to fixed point variable and remove zeros
 word
-real2FixedStripZeros(const real& v, int32 numPrecision = 6);
+real2FixedStripZeros(const float& v, int32 numPrecision = 6);
+
+word
+real2FixedStripZeros(const double& v, int32 numPrecision = 6);
 
 /// Output <w1,w2>
 word
@@ -225,7 +236,14 @@ readValue(const word& w, bool& val)
 
 INLINE_FUNCTION_HD
 bool
-equal(const real& s1, const real& s2, real tol = smallValue)
+equal(const float& s1, const float& s2, float tol = smallValue)
+{
+	return abs(s1 - s2) <= tol;
+}
+
+INLINE_FUNCTION_HD
+bool
+equal(const double& s1, const double& s2, double tol = smallValue)
 {
 	return abs(s1 - s2) <= tol;
 }

src/phasicFlow/types/basicTypes/builtinTypes.hpp

@@ -63,6 +63,8 @@ using size_t = std::size_t;
 
 using word = std::string;
 
+using timeValue = double;
+
 inline const int numBytesForReal__ = sizeof(real);
 
 inline word

src/phasicFlow/types/basicTypes/math.hpp

@@ -42,8 +42,19 @@ namespace pFlow
 {
 
 INLINE_FUNCTION_HD
-real
+double
-abs(real x)
+abs(double x)
+{
+#ifdef __CUDACC__
+	return ::fabs(x);
+#else
+	return std::fabs(x);
+#endif
+}
+
+INLINE_FUNCTION_HD
+float
+abs(float x)
 {
 #ifdef __CUDACC__
 	return ::fabs(x);
@@ -73,8 +84,18 @@ abs(int32 x)
 #endif
 }
 
-INLINE_FUNCTION_HD real
+INLINE_FUNCTION_HD float
-mod(real x, real y)
+mod(float x, float y)
+{
+#ifdef __CUDACC__
+	return ::fmod(x, y);
+#else
+	return std::fmod(x, y);
+#endif
+}
+
+INLINE_FUNCTION_HD double
+mod(double x, double y)
 {
 #ifdef __CUDACC__
 	return ::fmod(x, y);
@@ -108,8 +129,18 @@ mod(uint32 x, uint32 y)
 	return x % y;
 }
 
-INLINE_FUNCTION_HD real
+INLINE_FUNCTION_HD float
-remainder(real x, real y)
+remainder(float x, float y)
+{
+#ifdef __CUDACC__
+    return ::remainder(x,y);
+#else
+	return std::remainder(x, y);
+#endif
+}
+
+INLINE_FUNCTION_HD double
+remainder(double x, double y)
 {
 #ifdef __CUDACC__
     return ::remainder(x,y);
@@ -120,8 +151,19 @@ remainder(real x, real y)
 
 
 INLINE_FUNCTION_HD
-real
+float
-exp(real x)
+exp(float x)
+{
+#ifdef __CUDACC__
+    return ::exp(x);
+#else
+	return std::exp(x);
+#endif
+}
+
+INLINE_FUNCTION_HD
+double
+exp(double x)
 {
 #ifdef __CUDACC__
     return ::exp(x);
@@ -132,8 +174,17 @@ exp(real x)
 
 
 INLINE_FUNCTION_HD
-real
+float log(float x)
-log(real x)
+{
+#ifdef __CUDACC__
+    return ::log(x);
+#else
+	return std::log(x);
+#endif
+}
+
+INLINE_FUNCTION_HD
+double log(double x)
 {
 #ifdef __CUDACC__
     return ::log(x);
@@ -144,8 +195,7 @@ log(real x)
 
 
 INLINE_FUNCTION_HD
-real
+float log10(float x)
-log10(real x)
 {
 #ifdef __CUDACC__
     return ::log10(x);
@@ -155,8 +205,27 @@ log10(real x)
 }
 
 INLINE_FUNCTION_HD
-real
+double log10(double x)
-pow(real x, real y)
+{
+#ifdef __CUDACC__
+    return ::log10(x);
+#else
+	return std::log10(x);
+#endif
+}
+
+INLINE_FUNCTION_HD
+float pow(float x, float y)
+{
+#ifdef __CUDACC__
+	return ::pow(x,y);
+#else
+	return std::pow(x, y);
+#endif
+}
+
+INLINE_FUNCTION_HD
+double pow(double x, double y)
 {
 #ifdef __CUDACC__
 	return ::pow(x,y);
@@ -167,8 +236,17 @@ pow(real x, real y)
 
 
 INLINE_FUNCTION_HD
-real
+float sqrt(float x)
-sqrt(real x)
+{
+#ifdef __CUDACC__
+    return ::sqrt(x);
+#else
+	return std::sqrt(x);
+#endif
+}
+
+INLINE_FUNCTION_HD
+double sqrt(double x)
 {
 #ifdef __CUDACC__
     return ::sqrt(x);
@@ -180,8 +258,7 @@ sqrt(real x)
 
 
 INLINE_FUNCTION_HD
-real
+float cbrt(float x)
-cbrt(real x)
 {
 #ifdef __CUDACC__
     return ::cbrt(x);
@@ -191,8 +268,17 @@ cbrt(real x)
 }
 
 INLINE_FUNCTION_HD
-real
+double cbrt(double x)
-sin(real x)
+{
+#ifdef __CUDACC__
+    return ::cbrt(x);
+#else
+	return std::cbrt(x);
+#endif 
+}
+
+INLINE_FUNCTION_HD
+float sin(float x)
 {
 #ifdef __CUDACC__
     return ::sin(x);
@@ -202,8 +288,17 @@ sin(real x)
 }
 
 INLINE_FUNCTION_HD
-real
+double sin(double x)
-cos(real x)
+{
+#ifdef __CUDACC__
+    return ::sin(x);
+#else
+	return std::sin(x);
+#endif
+}
+
+INLINE_FUNCTION_HD
+real cos(real x)
 {
 #ifdef __CUDACC__
     return ::cos(x);
@@ -215,8 +310,7 @@ cos(real x)
 
 
 INLINE_FUNCTION_HD
-real
+real tan(real x)
-tan(real x)
 {
 #ifdef __CUDACC__
     return ::tan(x);
@@ -226,8 +320,7 @@ tan(real x)
 }
 
 INLINE_FUNCTION_HD
-real
+real asin(real x)
-asin(real x)
 {
 #ifdef __CUDACC__
     return ::asin(x);
@@ -237,8 +330,7 @@ asin(real x)
 }
 
 INLINE_FUNCTION_HD
-real
+real acos(real x)
-acos(real x)
 {
 #ifdef __CUDACC__
     return ::acos(x);
@@ -248,8 +340,7 @@ acos(real x)
 }
 
 INLINE_FUNCTION_HD
-real
+real atan(real x)
-atan(real x)
 {
 #ifdef __CUDACC__
     return ::atan(x);
@@ -269,8 +360,7 @@ atan2(real y, real x)
 }
 
 INLINE_FUNCTION_HD
-real
+real sinh(real x)
-sinh(real x)
 {
 #ifdef __CUDACC__
     return ::sinh(x);
@@ -280,8 +370,7 @@ sinh(real x)
 }
 
 INLINE_FUNCTION_HD
-real
+real cosh(real x)
-cosh(real x)
 {
 #ifdef __CUDACC__
     return ::cosh(x);
@@ -322,8 +411,7 @@ acosh(real x)
 }
 
 INLINE_FUNCTION_HD
-real
+real atanh(real x)
-atanh(real x)
 {
 #ifdef __CUDACC__
     return ::atanh(x);
@@ -333,8 +421,18 @@ atanh(real x)
 }
 
 
-INLINE_FUNCTION_HD real
+INLINE_FUNCTION_HD 
-min(real x, real y)
+float min(float x, float y)
+{
+#ifdef __CUDACC__
+	return ::fmin(x, y);
+#else
+	return std::fmin(x, y);
+#endif
+}
+
+INLINE_FUNCTION_HD 
+double min(double x, double y)
 {
 #ifdef __CUDACC__
 	return ::fmin(x, y);
@@ -344,8 +442,7 @@ min(real x, real y)
 }
 
 INLINE_FUNCTION_HD
-int64
+int64 min(int32 x, int32 y)
-min(int32 x, int32 y)
 {
 #ifdef __CUDACC__
     return ::min(x,y);
@@ -356,8 +453,7 @@ min(int32 x, int32 y)
 
 
 INLINE_FUNCTION_HD
-int64
+int64 min(int64 x, int64 y)
-min(int64 x, int64 y)
 {
 #ifdef __CUDACC__
     return ::min(x,y);
@@ -369,8 +465,17 @@ min(int64 x, int64 y)
 
 
 INLINE_FUNCTION_HD
-uint64
+uint64 min(uint64 x, uint64 y)
-min(uint64 x, uint64 y)
+{
+#ifdef __CUDACC__
+    return ::min(x,y);
+#else
+	return std::min(x, y);
+#endif
+}
+
+INLINE_FUNCTION_HD 
+uint32 min(uint32 x, uint32 y)
 {
 #ifdef __CUDACC__
     return ::min(x,y);
@@ -380,20 +485,18 @@ min(uint64 x, uint64 y)
 }
 
 
-
+INLINE_FUNCTION_HD 
-INLINE_FUNCTION_HD uint32
+float max(float x, float y)
-min(uint32 x, uint32 y)
 {
 #ifdef __CUDACC__
-    return ::min(x,y);
+	return ::fmax(x, y);
 #else
-	return std::min(x, y);
+	return std::fmax(x, y);
 #endif
 }
 
-
+INLINE_FUNCTION_HD 
-INLINE_FUNCTION_HD real
+double max(double x, double y)
-max(real x, real y)
 {
 #ifdef __CUDACC__
 	return ::fmax(x, y);

utilities/pFlowToVTK/fileSeries.cpp

@@ -1,3 +1,23 @@
+/*------------------------------- phasicFlow ---------------------------------
+      O        C enter of
+     O O       E ngineering and
+    O   O      M ultiscale modeling of
+   OOOOOOO     F luid flow       
+------------------------------------------------------------------------------
+  Copyright (C): www.cemf.ir
+  email: hamid.r.norouzi AT gmail.com
+------------------------------------------------------------------------------  
+Licence:
+  This file is part of phasicFlow code. It is a free software for simulating 
+  granular and multiphase flows. You can redistribute it and/or modify it under
+  the terms of GNU General Public License v3 or any other later versions. 
+ 
+  phasicFlow is distributed to help others in their research in the field of 
+  granular and multiphase flows, but WITHOUT ANY WARRANTY; without even the
+  implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
+
+-----------------------------------------------------------------------------*/
+
 #include "fileSeries.hpp"
 #include "oFstream.hpp"
 #include "streams.hpp"

utilities/pFlowToVTK/fileSeries.hpp

@@ -1,5 +1,22 @@
+/*------------------------------- phasicFlow ---------------------------------
+      O        C enter of
+     O O       E ngineering and
+    O   O      M ultiscale modeling of
+   OOOOOOO     F luid flow       
+------------------------------------------------------------------------------
+  Copyright (C): www.cemf.ir
+  email: hamid.r.norouzi AT gmail.com
+------------------------------------------------------------------------------  
+Licence:
+  This file is part of phasicFlow code. It is a free software for simulating 
+  granular and multiphase flows. You can redistribute it and/or modify it under
+  the terms of GNU General Public License v3 or any other later versions. 
+ 
+  phasicFlow is distributed to help others in their research in the field of 
+  granular and multiphase flows, but WITHOUT ANY WARRANTY; without even the
+  implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
 
-
+-----------------------------------------------------------------------------*/
 
 #ifndef __fileSeries_hpp__
 #define __fileSeries_hpp__

utilities/pFlowToVTK/pointFieldToVTK.cpp

@@ -1,3 +1,23 @@
+/*------------------------------- phasicFlow ---------------------------------
+      O        C enter of
+     O O       E ngineering and
+    O   O      M ultiscale modeling of
+   OOOOOOO     F luid flow       
+------------------------------------------------------------------------------
+  Copyright (C): www.cemf.ir
+  email: hamid.r.norouzi AT gmail.com
+------------------------------------------------------------------------------  
+Licence:
+  This file is part of phasicFlow code. It is a free software for simulating 
+  granular and multiphase flows. You can redistribute it and/or modify it under
+  the terms of GNU General Public License v3 or any other later versions. 
+ 
+  phasicFlow is distributed to help others in their research in the field of 
+  granular and multiphase flows, but WITHOUT ANY WARRANTY; without even the
+  implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
+
+-----------------------------------------------------------------------------*/
+
 #include <regex>
 
 #include "vtkFile.hpp"

utilities/pFlowToVTK/triSurfaceFieldToVTK.cpp

@@ -1,3 +1,23 @@
+/*------------------------------- phasicFlow ---------------------------------
+      O        C enter of
+     O O       E ngineering and
+    O   O      M ultiscale modeling of
+   OOOOOOO     F luid flow       
+------------------------------------------------------------------------------
+  Copyright (C): www.cemf.ir
+  email: hamid.r.norouzi AT gmail.com
+------------------------------------------------------------------------------  
+Licence:
+  This file is part of phasicFlow code. It is a free software for simulating 
+  granular and multiphase flows. You can redistribute it and/or modify it under
+  the terms of GNU General Public License v3 or any other later versions. 
+ 
+  phasicFlow is distributed to help others in their research in the field of 
+  granular and multiphase flows, but WITHOUT ANY WARRANTY; without even the
+  implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
+
+-----------------------------------------------------------------------------*/
+
 #include "vocabs.hpp"
 #include "vtkFile.hpp"
 #include "triSurfaceFieldToVTK.hpp"

utilities/pFlowToVTK/triSurfaceFieldToVTK.hpp

@@ -18,7 +18,6 @@ Licence:
 
 -----------------------------------------------------------------------------*/
 
-
 #ifndef __triSurfaceFieldToVTK_hpp__
 #define __triSurfaceFieldToVTK_hpp__
 

utilities/postprocessPhasicFlow/postprocessPhasicFlow.cpp

@@ -103,7 +103,7 @@ int main(int argc, char** argv )
 			continue;
 		}
 				
-		if( !withZeroFolder && pFlow::equal(folders.time() , 0.0))continue;
+		if( !withZeroFolder && pFlow::equal(folders.time() , pFlow::zero))continue;
 
 		post.processTimeFolder(folders);