utilities added

utilities/particlesPhasicFlow/CMakeLists.txt

@@ -0,0 +1,11 @@
+
+set(source_files
+particlesPhasicFlow.C
+positionParticles/positionParticles.C 
+positionOrdered/positionOrdered.C 
+positionRandom/positionRandom.C
+empty/empty.C
+)
+set(link_lib phasicFlow Kokkos::kokkos Interaction)
+
+pFlow_make_executable_install(particlesPhasicFlow source_files link_lib)

utilities/particlesPhasicFlow/empty/empty.C

@@ -0,0 +1,36 @@
+/*------------------------------- 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 "empty.H"
+
+
+
+pFlow::empty::empty(
+	const dictionary& dict
+)
+:
+	positionParticles(dict),
+	position_
+	(
+		maxNumberOfParticles_, RESERVE()
+	)
+{
+	position_.clear();
+}

utilities/particlesPhasicFlow/empty/empty.H

@@ -0,0 +1,94 @@
+/*------------------------------- 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 __empty_H__
+#define __empty_H__
+
+#include "positionParticles.H"
+#include "box.H"
+
+namespace pFlow
+{
+
+
+class empty 
+:
+	public positionParticles
+{
+protected:
+
+	dictionary 	emptyDict_;
+
+
+	realx3Vector 	position_;
+
+public:
+
+	// - type Info
+	TypeName("empty");
+
+	empty(const dictionary& dict);
+
+	// - add this class to vCtor selection table 
+	add_vCtor(
+		positionParticles,
+		empty,
+		dictionary);
+
+	virtual ~empty() = default;
+
+	//// - Methods 
+
+	virtual label numPoints()const
+	{
+		return 0;
+	}
+
+	virtual label size()const
+	{
+		return 0;
+	}
+
+	real maxDiameter() const override
+	{
+		return 1.0;
+	}
+
+	// - const access to position
+	virtual const realx3Vector& position()const 
+	{
+		return position_;
+	}
+
+	// - access to position 
+	virtual realx3Vector& position()
+	{
+		return position_;
+	}
+
+	
+};	
+
+
+}
+
+
+
+#endif // __empety_H__

utilities/particlesPhasicFlow/particlesPhasicFlow.C

@@ -0,0 +1,195 @@
+/*------------------------------- 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 "positionParticles.H"
+#include "pointStructure.H"
+#include "setFields.H"
+#include "systemControl.H"
+#include "commandLine.H"
+
+using pFlow::output;
+using pFlow::endl;
+using pFlow::dictionary;
+using pFlow::uniquePtr;
+using pFlow::IOobject;
+using pFlow::objectFile;
+using pFlow::fileSystem;
+using pFlow::pointStructure;
+using pFlow::pointStructureFile__;
+using pFlow::positionParticles;
+using pFlow::commandLine;
+using pFlow::setFieldList;
+
+int main( int argc, char* argv[] )
+{
+
+	commandLine cmds(
+		"createParticles",
+		"Read the dictionary createParticles and create particles"
+		" based on the two sub-dictionaries positionParticles and setFields.\n"
+		"First executespositionParticles and then setFields, except "
+		"otherwise selected in the command line.");
+
+			
+	bool positionOnly = false;
+	cmds.add_flag(
+		"--positionParticles-only",
+		positionOnly,
+		"Exectue the positionParticles part only and store the created "
+		"pointStructure in the time folder.");
+	
+	bool setOnly = false;
+	cmds.add_flag("--setFields-only",
+		setOnly,
+		"Exectue the setFields part only. Read the pointStructure from "
+		"time folder and setFields and save the result in the same time folder.");
+
+
+	if(!cmds.parse(argc, argv)) return 0;
+	
+	if(setOnly && positionOnly)
+	{
+		Err<<
+		"Options --positionParticles-only and --setFields-only cannot be used simeltanuously. \n"<<endErr;
+		return 1;
+	}	
+
+// this should be palced in each main 
+#include "initialize_Control.H"
+
+	auto objCPDict = IOobject::make<dictionary>
+	(
+		objectFile
+		(
+			"createParticles",
+			Control.settings().path(),
+			objectFile::READ_ALWAYS,
+			objectFile::WRITE_ALWAYS
+		),
+		"createParticles",
+		true
+	);
+
+	auto& cpDict = objCPDict().getObject<dictionary>();
+
+	uniquePtr<IOobject> pStructObj{nullptr};
+
+	if(!setOnly)
+	{
+
+		// position particles based on the dict content 
+		Report(0)<< "Positioning points . . . \n"<<endReport;
+		auto pointPosition = positionParticles::create(cpDict.subDict("positionParticles"));
+
+		fileSystem pStructPath = Control.time().path()+pointStructureFile__;
+
+		auto finalPos = pointPosition().getFinalPosition();
+
+		
+		pStructObj = IOobject::make<pointStructure>
+		(
+			objectFile
+			(
+				pointStructureFile__,
+				Control.time().path(),
+				objectFile::READ_NEVER,
+				objectFile::WRITE_ALWAYS			
+			),
+			finalPos
+		);
+
+		auto& pSruct = pStructObj().getObject<pointStructure>();
+
+		Report(1)<< "Created pStruct with "<< pSruct.size() << " points and capacity "<<
+		pSruct.capacity()<<" . . ."<< endReport;
+
+		Report(1)<< "Writing pStruct to " << pStructObj().path() << endReport<<endl<<endl;
+
+		if( !pStructObj().write())
+		{
+			fatalErrorInFunction<<
+			"Error in writing to file. \n ";
+			return 1; 
+		}
+	}else
+	{
+		pStructObj = IOobject::make<pointStructure>
+		(
+			objectFile
+			(
+				pointStructureFile__,
+				Control.time().path(),
+				objectFile::READ_ALWAYS,
+				objectFile::WRITE_NEVER			
+			)
+		);
+	}
+
+	
+
+	if(!positionOnly)
+	{
+
+		auto& pStruct = pStructObj().getObject<pointStructure>();
+
+		auto& sfDict = cpDict.subDict("setFields");
+
+		setFieldList defValueList(sfDict.subDict("defaultValue"));
+
+		for(auto& sfEntry: defValueList)
+		{
+			if( !sfEntry.setPointFieldDefaultValueNewAll(Control.time(), pStruct, true))
+			{
+				Err<< "\n error occured in setting default value fields.\n"<<endErr;
+				return 1;
+			}
+		}
+		
+		output<<endl;
+
+		auto& selectorsDict = sfDict.subDict("selectors");
+
+		auto selNames = selectorsDict.dictionaryKeywords();
+
+		for(auto name: selNames)
+		{
+			Report(1)<< "Applying selector " << greenText(name) <<endReport;
+			
+			if(
+				!applySelector(Control, pStruct, selectorsDict.subDict(name))
+			 )
+			{
+				Err<<"\n error occured in setting selector. \n"<<endErr;
+				return 1;
+			}
+			output<<endl;
+		}
+	}
+
+	Control.time().write(true);
+	Report(0)<< greenText("\nFinished successfully.\n")<<endReport;
+
+
+// this should be palced in each main 
+#include "finalize.H"
+
+	return 0;
+} 

utilities/particlesPhasicFlow/positionOrdered/positionOrdered.C

@@ -0,0 +1,161 @@
+/*------------------------------- 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 "positionOrdered.H"
+#include "error.H"
+
+
+#include "streams.H"
+
+
+bool pFlow::positionOrdered::findAxisIndex()
+{
+	if( axisOrder_.size() != 3 )
+	{
+		fatalErrorInFunction <<
+		"axisOrder should have 3 components, but " << axisOrder_.size() << 
+		" has been provided. \n";
+		return false;
+	}
+
+
+	if( axisOrder_[0] == axisOrder_[1] || 
+		axisOrder_[0] == axisOrder_[2] ||
+		axisOrder_[1] == axisOrder_[2] )
+	{
+		fatalErrorInFunction<<
+		"invalid/repeated axis names in axisOrder. This is provided: " << axisOrder_ <<endl;
+		return false;
+	}
+
+	realx3 uV[3];
+	size_t i=0;
+	for(auto& ca: axisOrder_)
+	{
+		if(ca == "x")
+		{
+			uV[i] = realx3(1.0, 0.0, 0.0);
+		}
+		else if(ca == "y")
+		{
+			uV[i] = realx3(0.0, 1.0, 0.0);
+		}
+		else if(ca == "z")
+		{
+			uV[i] = realx3(0.0, 0.0, 1.0);
+		}
+		else
+		{
+			fatalErrorInFunction << 
+			"unknown name for axis in axisOrder: " << ca <<endl; 
+			return false;
+		}
+		i++;
+	}
+
+	uVector1_ = uV[0];
+	uVector2_ = uV[1];
+	uVector3_ = uV[2];
+
+	return true;
+}
+
+bool pFlow::positionOrdered::positionPointsOrdered()
+{
+	position_.clear();
+
+	realx3 dl(diameter_);
+	auto minP = static_cast<real>(0.5)*dl + box_.minPoint();
+	auto maxP = static_cast<real>(0.5)*dl + box_.maxPoint();
+
+	auto cntr = minP;
+
+	size_t n = 0;
+	while( n < numPoints_ )
+	{
+		position_.push_back(cntr);
+
+		cntr += dl*uVector1_;
+		
+		if( dot(uVector1_, cntr) > dot(uVector1_, maxP) )
+		{
+			cntr = (minP*uVector1_) + ( (cntr+dl) * uVector2_) + (cntr*uVector3_);
+
+			if( dot(uVector2_, cntr) > dot(uVector2_, maxP) )
+			{
+				cntr = (cntr*uVector1_) + (minP*uVector2_) + ((cntr+dl)*uVector3_);
+
+				if( dot(uVector3_,cntr) > dot(uVector3_, maxP) )
+				{
+					fatalErrorInFunction<<
+					"positioned " << n << " points in the domain and it is full. \n" <<
+					"request to position "<< numPoints_<< " points has failed.\n";
+					return false;
+				}
+			}
+		}
+		n++;
+	}
+
+	return true;
+}
+
+pFlow::positionOrdered::positionOrdered
+(
+	const dictionary& dict
+)
+:
+	positionParticles(dict),
+	poDict_
+	(  
+		dict.subDict("positionOrderedInfo")
+	),
+	diameter_
+	(
+		poDict_.getVal<real>("diameter")
+	),
+	numPoints_
+	(
+		poDict_.getVal<size_t>("numPoints")
+	),
+	axisOrder_
+	(
+		poDict_.getValOrSet("axisOrder", wordList{"x", "y", "z"})
+	),
+	box_
+	(
+		poDict_.subDict("box")
+	),
+	position_
+	(
+		maxNumberOfParticles_, RESERVE()
+	)
+{
+	
+	if( !findAxisIndex() )
+	{
+		fatalExit;
+	}
+
+	if(!positionPointsOrdered())
+	{
+		fatalExit;
+	}
+}

utilities/particlesPhasicFlow/positionOrdered/positionOrdered.H

@@ -0,0 +1,115 @@
+/*------------------------------- 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 __positionOrdered_H__
+#define __positionOrdered_H__
+
+#include "positionParticles.H"
+#include "box.H"
+
+namespace pFlow
+{
+
+
+class positionOrdered 
+:
+	public positionParticles
+{
+protected:
+
+	dictionary 	poDict_;
+
+	real 		diameter_;
+
+	size_t 		numPoints_;
+	
+	wordList 	axisOrder_;
+
+	box 		box_;
+
+	// - unit vector of the first axis
+	realx3 		uVector1_;
+
+	// - unit vector of the second axis
+	realx3 		uVector2_;
+
+	// - unit vector of the third axis
+	realx3 		uVector3_;
+
+
+	realx3Vector 	position_;
+
+	bool findAxisIndex();
+
+	bool positionPointsOrdered();
+
+public:
+
+	// - type Info
+	TypeName("positionOrdered");
+
+	positionOrdered(const dictionary& dict);
+
+	// - add this class to vCtor selection table 
+	add_vCtor(
+		positionParticles,
+		positionOrdered,
+		dictionary);
+
+	virtual ~positionOrdered() = default;
+
+	//// - Methods 
+
+	virtual label numPoints()const
+	{
+		return position_.size();
+	}
+
+	virtual label size()const
+	{
+		return position_.size();
+	}
+
+	real maxDiameter() const override
+	{
+		return diameter_;
+	}
+	
+	// - const access to position
+	virtual const realx3Vector& position()const 
+	{
+		return position_;
+	}
+
+	// - access to position 
+	virtual realx3Vector& position()
+	{
+		return position_;
+	}
+
+	
+};	
+
+
+}
+
+
+
+#endif // __positionOrdered_H__

utilities/particlesPhasicFlow/positionParticles/positionParticles.C

@@ -0,0 +1,110 @@
+/*------------------------------- 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 "positionParticles.H"
+#include "cells.H"
+#include "contactSearchFunctions.H"
+
+#include "streams.H"
+
+
+pFlow::realx3Vector pFlow::positionParticles::sortByMortonCode(realx3Vector& position)const
+{
+	struct indexMorton
+	{
+		size_t morton;
+		size_t index;
+	};
+
+	realx3 minP = min(position);
+	realx3 maxP = max(position);
+	real cellsize = maxDiameter();
+	cells<size_t> allCells( box(minP, maxP), cellsize);
+
+	Vector<indexMorton> indMor(position.size(),RESERVE());
+
+	indMor.clear();
+
+	size_t ind=0;
+	for(const auto& p:position)
+	{
+		auto cellInd = allCells.pointIndex(p);
+		indMor.push_back( 
+			{ xyzToMortonCode64(cellInd.x(), cellInd.y(), cellInd.z()),
+			ind++});
+	}
+
+	Info<<"Performing morton sorting."<<endInfo;
+	std::sort(
+		indMor.begin(),
+		indMor.end(),
+		[]( const indexMorton &lhs, const indexMorton &rhs){
+			return lhs.morton < rhs.morton; } );
+
+	realx3Vector sortedPos(position.capacity(), RESERVE());
+	sortedPos.clear();
+
+
+	for(auto& ind:indMor)
+	{
+		sortedPos.push_back( position[ind.index] );
+	}
+
+	return sortedPos;
+}
+
+
+pFlow::positionParticles::positionParticles
+(
+	const dictionary& dict
+)
+{
+	maxNumberOfParticles_ = dict.getValOrSet("maxNumberOfParticles", static_cast<size_t>(10000));
+	
+	mortonSorting_ = dict.getValOrSet("mortonSorting", Logical("Yes"));	
+}
+
+
+
+
+pFlow::uniquePtr<pFlow::positionParticles> pFlow::positionParticles::create(const dictionary & dict)
+{
+
+	word method = dict.getVal<word>("method");
+
+
+	if( dictionaryvCtorSelector_.search(method) )
+	{
+		return dictionaryvCtorSelector_[method] (dict);
+	}
+	else
+	{
+		printKeys
+		( 
+			fatalError << "Ctor Selector "<< method << " dose not exist. \n"
+			<<"Avaiable ones are: \n\n"
+			,
+			dictionaryvCtorSelector_
+		);
+		fatalExit;
+	}
+
+	return nullptr;
+}

utilities/particlesPhasicFlow/positionParticles/positionParticles.H

@@ -0,0 +1,98 @@
+/*------------------------------- 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 __positionParticles_H__
+#define __positionParticles_H__
+
+#include "virtualConstructor.H"
+#include "Vectors.H"
+#include "dictionary.H"
+
+namespace pFlow
+{
+
+
+class positionParticles
+{
+protected:
+
+	size_t  		maxNumberOfParticles_ = 10000;
+
+	Logical 		mortonSorting_;
+
+	static const size_t numReports_ = 40;
+
+	realx3Vector sortByMortonCode(realx3Vector& position)const;
+	
+public:
+
+	// - type Info
+	TypeName("positionParticles");
+
+	positionParticles(const dictionary& dict);
+
+	create_vCtor
+	(
+		positionParticles,
+		dictionary,
+		(const dictionary& dict),
+		(dict)
+	);
+
+	virtual ~positionParticles() = default;
+
+	//// - Methods 
+
+	virtual label numPoints()const = 0;
+
+	virtual label size()const = 0;
+
+	virtual real maxDiameter() const = 0;
+
+	// - const access to position
+	virtual const realx3Vector& position()const = 0;
+	
+
+	// - access to position 
+	virtual realx3Vector& position() = 0;
+	
+	virtual realx3Vector getFinalPosition()
+	{
+		if(mortonSorting_)
+		{
+			return sortByMortonCode(position());
+		}
+		else
+		{
+			return position();
+		}
+	}
+
+
+	static
+	uniquePtr<positionParticles> create(const dictionary & dict);
+};	
+
+
+}
+
+
+
+#endif // __positionParticles_H__

utilities/particlesPhasicFlow/positionRandom/positionRandom.C

@@ -0,0 +1,260 @@
+/*------------------------------- 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 "positionRandom.H"
+#include "uniformRandomReal.H"
+#include "VectorSingles.H"
+#include "VectorDuals.H"
+#include "NBS.H"
+#include "unsortedPairs.H"
+
+
+
+namespace pFlow
+{
+
+using SearchType = NBS<DefaultExecutionSpace, int32> ;
+using ContainerType = unsortedPairs<DefaultExecutionSpace, int32>;
+
+
+void fillPoints(
+		uint numPoints,
+		realx3 minP,
+		realx3 maxP,
+		realx3Vector_HD& points,
+		int32Vector_HD& flags );
+
+int32 findCollisions(
+		ContainerType& pairs,
+		int32Vector_HD& flags);
+
+int32 findCollisions(int32 num, realx3Vector_HD& points, real diam)
+{
+	int32 res =0;
+	for(auto i=0; i<num;i++)
+	{
+		for(auto j=i+1; j<num; j++)
+		{
+			if(sphereSphereCheck(points[i],points[j],diam,diam))res++;
+		}
+	}
+
+	return res;
+}
+
+}
+
+
+bool pFlow::positionRandom::positionOnePass(int32 pass, int32 startNum)
+{
+
+	realVector_D diameter(startNum , diameter_);
+	int32Vector_HD  flagHD(startNum, 0);
+	realx3Vector_HD positionHD(startNum);
+	
+	auto minP = box_.minPoint() + static_cast<real>(0.5)*realx3(diameter_);
+	auto maxP = box_.maxPoint() - static_cast<real>(0.5)*realx3(diameter_);
+
+	SearchType search(
+					box_,
+					diameter_,
+					positionHD.deviceVectorAll(),
+					diameter.deviceVectorAll());
+
+	ContainerType pairs(3*startNum);
+		
+	Report(1)<< "Positioning "<<
+		greenText("(Pass #"<< pass+1<<")")<< 
+		": started with "<< startNum <<" points."<<endReport;
+	
+	fillPoints(startNum, minP, maxP, positionHD, flagHD);
+		
+	search.broadSearch(pairs, range(0, startNum), true);
+
+	
+	int32 numCollisions = findCollisions(pairs, flagHD);
+	
+	
+	Report(2)<< "Positioned " << cyanText(startNum - numCollisions) <<
+	" without collision \n"<<endReport; 
+
+	if(startNum-numCollisions >= numPoints_ )
+	{
+		
+		Report(1)<<"Selected "<< cyanText(numPoints_)<< " for the final field.\n"<<endReport;
+
+		positionHD.syncViews();
+		position_.clear();
+		int32 n=0;
+		for(int32 i=0; i<startNum; i++)
+		{
+			if(flagHD[i] == 0 )
+			{
+				position_.push_back( positionHD[i]);
+				n++;
+				if(n==numPoints_)break;	
+			}
+			
+		}
+		
+		return true;
+	}
+	
+	
+	return false;
+	
+	
+}
+
+bool pFlow::positionRandom::positionPointsRandom()
+{
+	
+	position_.clear();
+
+	if(numPoints_ == 0)return true;
+	
+	size_t pass = 0;
+	int32 startNum = numPoints_;
+
+	while ( pass <maxIterations_)
+	{
+		if( positionOnePass(pass, startNum) )return true;
+		startNum = 1.1*startNum+1;
+		pass++;
+	}	
+	
+	
+	fatalErrorInFunction<< 
+	"  cannot position "<< numPoints_ << " in the domain in " << maxIterations_ << " iterations.\n" <<
+	"  you may increase maxIterations for positioning points.\n";
+
+	return false;
+}
+
+bool pFlow::positionRandom::inCollision
+(
+	const realx3 &cntr,
+	real diam
+)
+{
+	for(const auto& cp: position_)
+	{
+		if( length(cp-cntr) <= diam ) return true;
+	}
+
+	return false;
+}
+
+pFlow::positionRandom::positionRandom
+(
+	const dictionary& dict
+)
+:
+	positionParticles(dict),
+	prDict_
+	(  
+		dict.subDict("positionRandomInfo")
+	),
+	diameter_
+	(
+		prDict_.getVal<real>("diameter")
+	),
+	numPoints_
+	(
+		prDict_.getVal<size_t>("numPoints")
+	),
+	box_
+	(
+		prDict_.subDict("box")
+	),
+	maxIterations_
+	(
+		prDict_.getValOrSet("maxIterations", 10)
+	),
+	position_
+	(
+		maxNumberOfParticles_, RESERVE()
+	)
+{
+
+	reportInterval_ = max(numPoints_/numReports_, static_cast<size_t>(2));
+
+	if( !positionPointsRandom() )
+	{
+		fatalExit;
+	}
+}
+
+void pFlow::fillPoints(
+	uint numPoints,
+	realx3 minP,
+	realx3 maxP,
+	realx3Vector_HD& points,
+	int32Vector_HD& flags )
+{
+
+	uniformRandomReal rand;
+	
+
+	for(size_t i=0; i<numPoints; i++)
+	{
+		if(flags[i] == 0)
+		{
+			points[i] =rand(minP, maxP);	
+		}
+	}
+	points.modifyOnHost();
+	points.syncViews();
+}
+
+pFlow::int32 pFlow::findCollisions(
+		ContainerType& pairs,
+		int32Vector_HD& flags)
+{
+	auto allPairs = pairs.getPairs();
+	auto num = pairs.capacity();
+	auto dFlags = flags.deviceVector();
+	
+
+	int32 numCollisions = 0;
+
+	Kokkos::parallel_reduce(
+		"positionRandom::findCollisions",
+		num,
+		LAMBDA_HD(int32 i, int32& valueToUpdate){
+			if(allPairs.isValid(i))
+			{
+				auto pair = allPairs.getPair(i);
+				if( dFlags[pair.first] ==0 )
+				{
+					dFlags[pair.first] = 1;
+					valueToUpdate++;	
+				}
+			}
+		}, numCollisions);
+
+	flags.modifyOnDevice();
+	flags.syncViews();	
+
+	return numCollisions;
+}
+
+

utilities/particlesPhasicFlow/positionRandom/positionRandom.H

@@ -0,0 +1,113 @@
+/*------------------------------- 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 __positionOrdered_H__
+#define __positionOrdered_H__
+
+#include "positionParticles.H"
+#include "VectorSingles.H"
+
+#include "box.H"
+
+namespace pFlow
+{
+
+class positionRandom
+:
+	public positionParticles
+{
+protected:
+
+	dictionary 	prDict_;
+
+	real 		diameter_;
+
+	size_t 		numPoints_;
+	
+	box 		box_;
+
+	size_t 		maxIterations_;
+	
+	realx3Vector 	position_;
+
+	size_t 		reportInterval_;
+
+	bool positionOnePass(int32 pass, int32 startNum);
+	
+	bool positionPointsRandom();
+
+	bool inCollision(const realx3 &cntr, real diam);
+
+public:
+
+	// - type Info
+	TypeName("positionRandom");
+
+	positionRandom(const dictionary& dict);
+
+	// - add this class to vCtor selection table 
+	add_vCtor(
+		positionParticles,
+		positionRandom,
+		dictionary);
+
+	virtual ~positionRandom() = default;
+
+	//// - Methods 
+
+	virtual label numPoints()const
+	{
+		return position_.size();
+	}
+
+	virtual label size()const
+	{
+		return position_.size();
+	}
+
+	real maxDiameter() const override
+	{
+		return diameter_;
+	}
+
+	// - const access to position
+	virtual const realx3Vector& position()const 
+	{
+		return position_;
+	}
+
+	// - access to position 
+	virtual realx3Vector& position()
+	{
+		return position_;
+	}
+
+
+
+	
+};	
+
+
+}
+
+
+
+#endif // __positionOrdered_H__

utilities/particlesPhasicFlow/setFields.H

@@ -0,0 +1,56 @@
+/*------------------------------- 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 __setFields_H__
+#define __setFields_H__
+
+#include "pStructSelector.H"
+#include "pointFields.H"
+#include "setFieldList.H"
+#include "systemControl.H"
+
+namespace pFlow
+{
+
+
+bool applySelector(systemControl& control, const pointStructure& pStruct, const dictionary& selDict)
+{
+
+	
+	auto selector = pStructSelector::create(pStruct, selDict);
+
+	auto& selected = selector().selectedPoinsts();
+
+	int32IndexContainer selIndex(selected.data(), selected.size());
+
+	setFieldList sfList(selDict.subDict("fieldValue"));
+
+	for(auto& sfEntry:sfList)
+	{
+		if(!sfEntry.setPointFieldSelectedAll(control.time(), selIndex, true ))return false;
+	}
+	
+	return true;
+}
+
+} // pFlow
+
+#endif //__setFields_H__