Utility postProcess is modified to be used in Version 1.0

- counting dictionary is added to postProcess
- Code need to be modified to cleaned (fields in the repository should be removed onces postProcess is passed the next time step)

CMakeLists.txt

@@ -76,7 +76,7 @@ add_subdirectory(solvers)
 add_subdirectory(utilities)
 
 #add_subdirectory(DEMSystems)
-add_subdirectory(testIO)
+#add_subdirectory(testIO)
 
 install(FILES "${PROJECT_BINARY_DIR}/phasicFlowConfig.H"
   DESTINATION include)

src/phasicFlow/repository/repository/repository.cpp

@@ -145,6 +145,7 @@ bool pFlow::repository::addToRepository(IOobject* io)
 	{
 		warningInFunction<<
     	"Failed to add object " << io->name() << 
+		" with type "<< lookupObjectTypeName(io->name())<<
     	" to repository " << this->name() << 
     	". It is already in this repository. \n";
     	return false;

utilities/CMakeLists.txt

@@ -9,5 +9,5 @@ add_subdirectory(pFlowToVTK)
 
 add_subdirectory(Utilities)
 
-#add_subdirectory(postprocessPhasicFlow)
+add_subdirectory(postprocessPhasicFlow)
 

utilities/Utilities/readControlDict.hpp

@@ -46,8 +46,8 @@ protected:
 
 	int32 precision_;
 
-	inline static fileSystem defaultRootPath = CWD();
+	inline static fileSystem defaultRootPath = pFlow::CWD();
-	inline static fileSystem defaultCDPath = CWD()/"system"+"controlDict";
+	inline static fileSystem defaultCDPath = pFlow::CWD()/word("system")+word("controlDict");
 	
 	word convertTimeToName(const real t, const int32 precision)const;
 

utilities/Utilities/readFromTimeFolder.hpp

@@ -113,7 +113,7 @@ public:
 	}
 
 	template<typename T>
-	uniquePtr<pointField_H<T>>  createUniformPointField_H(word const& name, T value)
+	uniquePtr<pointField_H<T>>  createUniformPointField_H(word const& name, T value, bool regRep = true)
 	{
 		if( !checkForPointStructure() )
 		{
@@ -124,7 +124,8 @@ public:
 
 		auto& pStruct = repository_.lookupObject<pointStructure>(pointStructureFile__);
 
-
+		if(regRep)
+		{
 			auto Ptr = makeUnique<pointField_H<T>>
 			(
 				objectFile
@@ -140,6 +141,27 @@ public:
 			);
 			return Ptr;
 		}
+		else
+		{
+			auto Ptr = makeUnique<pointField_H<T>>
+			(
+				objectFile
+				(
+					name,
+					"",
+					objectFile::READ_NEVER,
+					objectFile::WRITE_NEVER
+				),
+				nullptr,
+				pStruct,
+				value
+			);
+			Ptr().fill(value);
+			return Ptr;
+		}
+		
+		return nullptr;
+	}
 
 	template<typename T>
 	uniquePtr<pointField_H<T>> readPointField_H(word const& name) 

utilities/postprocessPhasicFlow/CMakeLists.txt

@@ -1,12 +1,15 @@
 
 set(source_files 
-postprocessPhasicFlow.cpp 
+cellMapper.cpp
+rectangleMesh.cpp
+countField.cpp
+countFields.cpp
 postprocess.cpp
 processField.cpp
 ProcessFields.cpp
 includeMask.cpp
 IncludeMasks.cpp
-
+postprocessPhasicFlow.cpp
 )
 set(link_lib phasicFlow Interaction Kokkos::kokkos Utilities)
 

utilities/postprocessPhasicFlow/IncludeMask.hpp

@@ -186,11 +186,12 @@ protected:
 		
 	Operator 				operator_;
 	
-	pointField_H<T> 		field_;
+	uniquePtr<pointField_H<T>> 		fieldPtr_;
 
+	hostViewType1D<T> 			field_;
 public:
 
-	TypeInfoTemplate2("IncludeMask", T, Operator);
+	TypeInfoTemplate12("IncludeMask", T, Operator);
 
 	IncludeMask(
 		const dictionary& dict,
@@ -199,8 +200,10 @@ public:
 	:
 		includeMask(dict, opType, timeFolder),
 		operator_(dict),
-		field_(timeFolder.readPointField_H<T>(this->fieldName()))
+		fieldPtr_(timeFolder.readPointField_H<T>(this->fieldName())),
-	{}
+		field_(fieldPtr_().hostView())
+	{
+	}
 
 	add_vCtor(
 		includeMask,
@@ -212,6 +215,11 @@ public:
 		return operator_(field_[n]);
 	}
 
+	uint32 size()const override
+	{
+		return field_.size();
+	}
+
 };
 
 
@@ -221,7 +229,7 @@ class IncludeMask<T,allOp<T>>
 	public includeMask
 {
 public:
-	TypeInfoTemplate2("IncludeMask", T, allOp<int8>);
+	TypeInfoTemplate12("IncludeMask", T, allOp<int8>);
 
 	IncludeMask(
 		const dictionary& dict,
@@ -240,6 +248,11 @@ public:
 	{
 		return true;
 	}
+
+	uint32 size()const override
+	{
+		return 0;
+	}
 };
 
 

utilities/postprocessPhasicFlow/ProcessField.hpp

@@ -40,14 +40,14 @@ class ProcessField
 
 protected:
 
-	pointField_H<T>& 	field_;
+	uniquePtr<pointField_H<T>> 	field_;
 
 
-	rectMeshField_H<T>& 	processedField_;
+	rectMeshField_H<T> 	processedField_;
 
 public:
 
-	TypeInfoTemplate("ProcessField", T);
+	TypeInfoTemplate11("ProcessField", T);
 
 
 	ProcessField(
@@ -58,25 +58,15 @@ public:
 		processField(dict, pToCell, rep),
 		field_(
 			this->isUniform()?
-			timeFolder().createUniformPointField_H(this->fieldName(), getUniformValue() ):
+			timeFolder().createUniformPointField_H(this->fieldName(), getUniformValue(), false ):
 			timeFolder().readPointField_H<T>(this->fieldName()) 
 			),
 		processedField_
 		(
-			processedRepository().emplaceObject<rectMeshField_H<T>>
-			(
-				objectFile
-				(
-					processedFieldName(),
-					"",
-					objectFile::READ_NEVER,
-					objectFile::WRITE_ALWAYS
-				),
 			mesh(),
 			processedFieldName(),
 			T{}
 		)
-		)
 	{
 		
 	}
@@ -99,32 +89,58 @@ public:
 		
 		const includeMask& incMask = includeMask_();
 		
-		auto numerator = sumMaksOp( field_ , this->pointToCell(), incMask);
+		auto numeratorPtr =  sumMaksOp( field_() , this->pointToCell(), incMask);
-		
+		uniquePtr<rectMeshField_H<real>> denomeratorPtr;
-		rectMeshField_H<real> denomerator( this->mesh(), real{} );
 		
 		if(operation() == "sum")
 		{
-			denomerator = rectMeshField_H<real>(this->mesh(), static_cast<real>(1.0));
+			denomeratorPtr = makeUnique<rectMeshField_H<real>>(this->mesh(), static_cast<real>(1.0));
 
 		}else if(operation() == "average")
 		{
 
-			pointField_H<real> oneFld(field_.pStruct(), static_cast<real>(1.0), static_cast<real>(1.0));
+			pointField_H<real> oneFld(
+				objectFile
+				(
+					"oneField",
+					"", 
+					objectFile::READ_NEVER, 
+					objectFile::WRITE_NEVER
+				),
+				const_cast<pointStructure&>(field_().pStruct()), 
+				static_cast<real>(1.0), 
+				static_cast<real>(1.0)
+			);
 
-			denomerator = sumOp(oneFld, this->pointToCell());
+			denomeratorPtr = sumOp(oneFld, this->pointToCell());
 
 		}else if(operation() == "averageMask")
 		{
-			pointField_H<real> oneFld(field_.pStruct(), static_cast<real>(1.0), static_cast<real>(1.0));
+			//pointField_H<real> oneFld(field_().pStruct(), static_cast<real>(1.0), static_cast<real>(1.0));
+			pointField_H<real> oneFld(
+				objectFile
+				(
+					"oneField",
+					"", 
+					objectFile::READ_NEVER, 
+					objectFile::WRITE_NEVER
+				),
+				const_cast<pointStructure&>(field_().pStruct()), 
+				static_cast<real>(1.0), 
+				static_cast<real>(1.0)
+			);
 			
-			denomerator = sumMaksOp(oneFld, this->pointToCell(), incMask);
+			denomeratorPtr = sumMaksOp(oneFld, this->pointToCell(), incMask);
-		}else
+
+		}
+		else
 		{
 			fatalErrorInFunction<<"operation is not known: "<< operation()<<endl;
 			fatalExit;
 		}
 
+		auto& numerator = numeratorPtr();
+		auto& denomerator = denomeratorPtr();
 
 		for(int32 i=0; i<this->mesh().nx(); i++ )
 		{

utilities/postprocessPhasicFlow/cellMapper.cpp

@@ -0,0 +1,91 @@
+/*------------------------------- 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 "cellMapper.hpp"
+//#include "cellMapperKernels.hpp"
+#include "streams.hpp"
+
+
+void pFlow::cellMapper::allocateArrays(rangeU32 nextRng)
+{
+    checkAllocateNext(nextRng);
+    nullifyNext(nextRng);
+    reallocFill(head_, domainCells_.nx(), domainCells_.ny(), domainCells_.nz(), NoPos);
+}
+
+void pFlow::cellMapper::checkAllocateNext(rangeU32 nextRng)
+{
+
+    auto newCap = nextRng.end();
+    
+    if( nextCapacity_ < newCap)
+    {
+        nextCapacity_ = newCap;	
+        reallocNoInit(next_, nextCapacity_);
+    }
+}
+
+void pFlow::cellMapper::nullifyHead()
+{
+    fill(head_, NoPos);	
+}
+
+void pFlow::cellMapper::nullifyNext(rangeU32 nextRng)
+{
+	fill(next_, nextRng, NoPos);
+}
+
+pFlow::cellMapper::cellMapper(
+    const rectangleMesh& rectMesh,
+	const hostViewType1D<realx3>& pointPos,
+	const pFlagTypeHost& flags
+)
+: 
+   	domainCells_(rectMesh)
+{
+
+    allocateArrays(flags.activeRange());
+}
+
+bool pFlow::cellMapper::build
+(
+    const hostViewType1D<realx3>& pointPos, 
+    const pFlagTypeHost & flags
+)
+{
+    auto aRange = flags.activeRange();
+    
+    checkAllocateNext(aRange);
+    nullifyHead();
+    nullifyNext(aRange);
+	int32x3 ind;
+	for(uint32 i=aRange.start(); i<aRange.end(); i++)
+	{
+		if( domainCells_.isInsideIndex(pointPos[i], ind) )
+		{
+			next_[i] = head_(ind.x(), ind.y(), ind.z());
+			head_(ind.x(), ind.y(), ind.z()) = i;
+		}
+	}     
+    
+    return true;
+
+}
+

utilities/postprocessPhasicFlow/cellMapper.hpp

@@ -0,0 +1,137 @@
+/*------------------------------- 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 __cellMapper_hpp__
+#define __cellMapper_hpp__
+
+#include "phasicFlowKokkos.hpp"
+#include "pointFlag.hpp"
+#include "rectangleMesh.hpp"
+
+
+namespace pFlow
+{
+
+class cellMapper
+{
+public:
+
+	using HeadType 			= hostViewType3D<uint32>;
+
+	using NextType 			= hostViewType1D<uint32>;	
+
+	
+	static constexpr uint32 NoPos 	= 0xFFFFFFFF;
+
+	class CellIterator
+	{
+	private:
+		HeadType 		head_;
+
+		NextType 		next_;
+
+	public:
+
+		CellIterator(const HeadType& head, const NextType& next)
+		:
+			head_(head),
+			next_(next)
+		{}
+		
+		static constexpr uint32 NoPos = 0xFFFFFFFF;
+
+		INLINE_FUNCTION_H
+		int32x3 numCells()const {
+			return int32x3(head_.extent(0), head_.extent(1), head_.extent(2));}
+
+		INLINE_FUNCTION_H
+		uint32 start(int32 i, int32 j, int32 k)const {
+			return head_(i,j,k); }
+
+		INLINE_FUNCTION_H
+		uint32 getNext(uint32 n)const {
+			if(n == NoPos ) return NoPos;
+			return next_(n); }
+
+		INLINE_FUNCTION_H
+		uint32 next(uint32 n)const{
+			return next_(n);}
+	};
+
+private:
+	
+	const rectangleMesh&  	domainCells_;
+
+	HeadType 			head_{"NBS::head",1,1,1};
+
+	NextType 			next_{"NBS::next", 1};
+
+	uint32 				nextCapacity_ = 0;
+
+	void allocateArrays(rangeU32 nextRng);
+
+	void checkAllocateNext(rangeU32 nextRng);
+	
+	void nullifyHead();
+	
+	void nullifyNext(rangeU32 nextRng);
+	
+public:
+
+	TypeInfoNV("cellMapper");
+
+		
+		cellMapper(
+			const rectangleMesh& rectMesh,
+			const hostViewType1D<realx3>& pointPos,
+			const pFlagTypeHost& flags);
+					
+		INLINE_FUNCTION_HD
+		cellMapper(const cellMapper&) = default;
+
+		INLINE_FUNCTION_HD
+		cellMapper(cellMapper&&) = default;
+		
+		INLINE_FUNCTION_HD
+		cellMapper& operator = (const cellMapper&) = default;
+
+		INLINE_FUNCTION_HD
+		cellMapper& operator = (cellMapper&&) = default;
+
+		INLINE_FUNCTION_HD
+		~cellMapper()=default;
+
+	//// - Methods
+
+		auto getCellIterator()const
+		{
+			return CellIterator(head_, next_);
+		}
+
+		
+		bool build(
+			const hostViewType1D<realx3>& pointPos,
+			const pFlagTypeHost& flags);
+};
+
+} // pFlow
+
+#endif // __cellMapper_hpp__

utilities/postprocessPhasicFlow/countField.cpp

@@ -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.
+
+-----------------------------------------------------------------------------*/
+
+#include "countField.hpp"
+#include "repository.hpp"
+#include "twoPartEntry.hpp"
+
+bool pFlow::countField::getFieldType(
+	const dictionary& dict,
+	readFromTimeFolder& timeFolder,
+	word& fieldName,
+	word& fieldType)
+{
+	if(dict.containsDataEntry("field"))
+	{
+		const dataEntry& entry = dict.dataEntryRef("field");
+
+		if( isTwoPartEntry(entry))
+		{
+			twoPartEntry tpEntry(entry);
+			fieldName = "uniformField";
+			fieldType = tpEntry.firstPart();
+		}
+		else
+		{
+			fieldName = dict.getVal<word>("field");
+			if( !timeFolder.pointFieldFileGetType(fieldName, fieldType) )
+			{
+				fatalErrorInFunction<<"error in reading field type from file "<< fieldName<<
+				"in folder "<< timeFolder.path()<<endl;
+				return false;
+			}
+		}
+	}
+	else
+	{
+		fatalErrorInFunction<< "dictionary "<< dict.globalName()<<
+		"does not contain field keyword"<<endl;
+		return false;
+	}
+
+	return true;
+}
+
+
+pFlow::countField::countField(const dictionary& dict, repository& rep)
+:
+	dict_(dict),
+	timeFolder_(rep)
+{
+
+	word includeMaskType = dict_.getVal<word>("includeMask");
+
+	auto& incDict = dict_.subDictOrCreate(includeMaskType+"Info");
+	
+	includeMask_ = includeMask::create(incDict, includeMaskType, timeFolder_);
+
+}
+
+
+bool pFlow::countField::process(uint32& countedValue)
+{
+	auto& incMask = includeMask_();
+
+	countedValue = 0;
+	uint32 n = incMask.size();
+
+	for(uint32 i=0; i<n; i++)
+	{
+		if( incMask(i) )
+		{
+			countedValue++;
+		}
+	}
+
+	return true;
+}

utilities/postprocessPhasicFlow/countField.hpp

@@ -0,0 +1,102 @@
+/*------------------------------- 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 __countField_hpp__
+#define __countField_hpp__
+
+
+#include "virtualConstructor.hpp"
+#include "dictionary.hpp"
+#include "readFromTimeFolder.hpp"
+#include "includeMask.hpp"
+
+
+namespace pFlow
+{
+
+
+class repository;
+
+class countField
+{
+protected:
+
+	dictionary 	dict_;
+	
+	mutable readFromTimeFolder 	timeFolder_;
+
+	uniquePtr<includeMask>  	includeMask_ = nullptr;
+	
+	bool static getFieldType(
+			const dictionary& dict,
+			readFromTimeFolder& timeFolder,
+			word& fieldName,
+			word& fieldType);
+
+public:
+
+	TypeInfo("countField");
+
+	countField(const dictionary& dict, repository& rep);
+
+	
+	auto& dict()
+	{
+		return dict_;
+	}
+
+	const auto& dict()const
+	{
+		return dict_;
+	}
+
+	auto& timeFolderRepository()
+	{
+		return timeFolder_.db();
+	}
+	
+	
+	auto& timeFolder()
+	{
+		return timeFolder_;
+	}
+
+	word variableName()const
+	{
+		return dict_.name();
+	}
+
+	// requires a class to read pointField from timefolder 
+	bool process(uint32& countedValue);
+
+
+
+	//virtual bool writeToFile(iOstream& is) const = 0;
+
+	/*static 
+	uniquePtr<countField> create(
+		const dictionary& dict,
+		repository& rep);*/
+};
+
+
+}
+
+
+#endif //__countField_hpp__

utilities/postprocessPhasicFlow/countFields.cpp

@@ -0,0 +1,53 @@
+/*------------------------------- 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 "countFields.hpp"
+#include "countField.hpp"
+#include "repository.hpp"
+#include "includeMask.hpp"
+
+
+pFlow::countFields::countFields(const dictionary& dict)
+:
+	dict_(dict),
+	variableNames_(dict.dictionaryKeywords()),
+	countedValues_(variableNames_.size(), 0u)
+{
+
+}
+
+
+bool pFlow::countFields::process(repository& rep)
+{
+	uint32List counted;
+
+	for(const auto& name: variableNames_)
+	{
+		const dictionary& varDict = dict_.subDict(name);
+		uint32 count;
+		countField cField(varDict, rep);
+		cField.process(count);
+		counted.push_back(count);
+	}
+
+	countedValues_ = counted;
+
+	return true;
+}

utilities/postprocessPhasicFlow/countFields.hpp

@@ -0,0 +1,77 @@
+/*------------------------------- 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 __countFields_hpp__
+#define __countFields_hpp__
+
+
+#include "dictionary.hpp"
+
+
+namespace pFlow
+{
+
+
+class repository;
+
+class countFields
+{
+protected:
+
+	dictionary 	dict_;
+	
+	wordList 	variableNames_;
+
+	uint32List 	countedValues_;
+
+public:
+
+	TypeInfo("countFields");
+
+	countFields(const dictionary& dict);
+
+	auto& dict()
+	{
+		return dict_;
+	}
+
+	const auto& dict()const
+	{
+		return dict_;
+	}
+
+	const wordList& variableNames()const
+	{
+		return variableNames_;
+	}
+	const uint32List& countedValues()const
+	{
+		return countedValues_;
+	}
+	
+	// requires a class to read pointField from timefolder 
+	bool process(repository& rep);
+
+};
+
+
+}
+
+
+#endif //__countFields_hpp__

utilities/postprocessPhasicFlow/fieldOperations.hpp

@@ -31,14 +31,15 @@ namespace pFlow
 
 
 template<typename T>
-rectMeshField_H<T> sumOp( const pointField_H<T> field, const pointRectCell& pointToCell)
+uniquePtr<rectMeshField_H<T>> sumOp( pointField_H<T>& field, pointRectCell& pointToCell)
 {
 	// create field
-	const auto& mesh = pointToCell.mesh();
+	auto& mesh = pointToCell.mesh();
 	auto iterator = pointToCell.getCellIterator();
+	auto f = field.deviceView();
 
-	rectMeshField_H<T> results(mesh, T(0));
+	auto resultsPtr = makeUnique<rectMeshField_H<T>>(mesh, T(0));
-
+	auto& results = resultsPtr();
 	for(int32 i=0; i<mesh.nx(); i++)
 	{
 		for(int32 j=0; j<mesh.ny(); j++)
@@ -49,7 +50,7 @@ rectMeshField_H<T> sumOp( const pointField_H<T> field, const pointRectCell& poin
 				T res (0);
 				while(n>-1)
 				{
-					res += field[n];
+					res += f[n];
 					n = iterator.getNext(n);
 				}
 
@@ -58,17 +59,19 @@ rectMeshField_H<T> sumOp( const pointField_H<T> field, const pointRectCell& poin
 		}
 	}
 
-	return results;
+	return resultsPtr;
 }
 
 template<typename T, typename incMask>
-rectMeshField_H<T> sumMaksOp( const pointField_H<T> field, const pointRectCell& pointToCell, const incMask& mask)
+uniquePtr<rectMeshField_H<T>> sumMaksOp( pointField_H<T>& field, pointRectCell& pointToCell, const incMask& mask)
 {
 	// create field
-	const auto& mesh = pointToCell.mesh();
+	auto& mesh = pointToCell.mesh();
 	auto iterator = pointToCell.getCellIterator();
+	auto f = field.deviceView();
 
-	rectMeshField_H<T> results(mesh, T(0));
+	auto resultsPtr = makeUnique<rectMeshField_H<T>>(mesh, T(0));
+	auto& results = resultsPtr();
 
 	for(int32 i=0; i<mesh.nx(); i++)
 	{
@@ -85,7 +88,7 @@ rectMeshField_H<T> sumMaksOp( const pointField_H<T> field, const pointRectCell&
 					
 					if(mask(n))
 					{
-						res += field[n];
+						res += f[n];
 					}
 
 					n = iterator.getNext(n);
@@ -96,7 +99,7 @@ rectMeshField_H<T> sumMaksOp( const pointField_H<T> field, const pointRectCell&
 		}
 	}
 
-	return results;
+	return resultsPtr;
 }
 
 

utilities/postprocessPhasicFlow/includeMask.cpp

@@ -85,7 +85,7 @@ pFlow::uniquePtr<pFlow::includeMask> pFlow::includeMask::create(
 		auto objPtr = 
 			dictionaryvCtorSelector_[method]
 			(dict, opType, timeFolder);
-		REPORT(2)<< dict.name()<< " with model "<<greenText(method)<<" is processed."<<endREPORT;
+		REPORT(2)<< dict.name()<< " with model "<<Green_Text(method)<<" is processed."<<END_REPORT;
 		return objPtr;
 	}
 	else

utilities/postprocessPhasicFlow/includeMask.hpp

@@ -85,6 +85,8 @@ public:
 
 	virtual bool isIncluded(int32 n) const = 0;
 
+	virtual uint32 size()const = 0;
+
 	bool operator()(int32 n) const 
 	{
 		return isIncluded(n);

utilities/postprocessPhasicFlow/pointRectCell.hpp

@@ -21,7 +21,7 @@ Licence:
 #ifndef __pointRectCell_hpp__
 #define __pointRectCell_hpp__
 
-#include "mapperNBS.hpp"
+#include "cellMapper.hpp"
 #include "rectMeshFields.hpp"
 #include "pointStructure.hpp"
 
@@ -33,20 +33,19 @@ class pointRectCell
 {
 public:
 
-	using mapType = mapperNBS<DefaultHostExecutionSpace>;
+	using mapType = cellMapper;
 
-	using memory_space = typename mapType::memory_space;
+	//using memory_space = typename mapType::memory_space;
 
 protected:
 
 	repository&				processedRepository_;
 
-	rectangleMesh& 			mesh_;
+	rectangleMesh 			mesh_;
-
 
 	const pointStructure& 	pStruct_;
 
-	ViewType1D<realx3, memory_space> 		pointPosition_;
+	hostViewType1D<realx3> 	pointPosition_;
 
 	mapType 				map_;
 
@@ -62,32 +61,28 @@ public:
 		processedRepository_(rep),
 		mesh_
 		(
-			processedRepository_.emplaceObject<rectangleMesh>
+			dictMesh,
-			(
+			&rep
-				objectFile
-				(
-					"rectMesh",
-					"",
-					objectFile::READ_NEVER,
-					objectFile::WRITE_NEVER
-				),
-				dictMesh
-			)
 		),
 		pStruct_(pStruct),
-		pointPosition_(pStruct_.pointPosition().hostVectorAll()),
+		pointPosition_(pStruct_.pointPosition().hostViewAll()),
-		map_( 
+		map_
-			mesh_.domain(), 
+		( 
-			mesh_.nx(),
+			mesh_,
-			mesh_.ny(),
+			pointPosition_,
-			mesh_.nz(),
+			pStruct_.activePointsMaskHost()
-			pointPosition_),
+		),
-		nPointInCell_(mesh_, 0)
+		nPointInCell_(mesh_,"nPointInCell", 0)
 	{
 
 		mapPOints();
 	}
 
+	auto& mesh()
+	{
+		return mesh_;
+	}
+
 	const auto& mesh()const
 	{
 		return mesh_;
@@ -100,20 +95,19 @@ public:
 
 	void mapPOints()
 	{
-		range activeRange = pStruct_.activeRange();
+		auto points = pStruct_.pointPositionHost();
-		auto activeMask  = pStruct_.activePointsMaskH();
+		auto activeMask  = pStruct_.activePointsMaskHost();
 
-		
+		map_.build(points, activeMask);
-		map_.buildCheckInDomain(activeRange, activeMask);
 		
 		auto iterator = map_.getCellIterator();
 		
 	
-		for(int32 i=0; i<map_.nx(); i++)
+		for(int32 i=0; i<mesh_.nx(); i++)
 		{
-			for(int32 j=0; j<map_.ny(); j++)
+			for(int32 j=0; j<mesh_.ny(); j++)
 			{
-				for(int32 k=0; k<map_.nz(); k++)
+				for(int32 k=0; k<mesh_.nz(); k++)
 				{
 
 					int32 res = 0;

utilities/postprocessPhasicFlow/postprocess.cpp

@@ -21,6 +21,7 @@ Licence:
 #include "postprocess.hpp" 
 #include "timeFolder.hpp"
 #include "pointStructure.hpp"
+#include "countFields.hpp"
 #include "vocabs.hpp"
 #include "vtkFile.hpp"
 
@@ -29,22 +30,26 @@ pFlow::postprocess::postprocess(systemControl& control)
 	control_(control),
 	dict_(postprocessFile__, control_.settings().path()+postprocessFile__)
 {
-	REPORT(1)<<"Reading numberBased dictionary ..."<<endREPORT;
+	REPORT(1)<<"Reading numberBased dictionary ..."<<END_REPORT;
 	auto nbDict = dict_.subDictOrCreate("numberBased");
 
 	numberBasedDictNames_ = dict_.subDictOrCreate("numberBased").dictionaryKeywords();
 	if(!numberBasedDictNames_.empty())
 	{
-		REPORT(1)<< "numberBased dictionary contains " << yellowText(numberBasedDictNames_)<<endREPORT<<endl;	
+		REPORT(1)<< "numberBased dictionary contains " << Yellow_Text(numberBasedDictNames_)<<END_REPORT<<endl;	
 	}
 	
-
 	weightBasedDictNames_ = dict_.subDictOrCreate("weightBased").dictionaryKeywords();
 	if(!weightBasedDictNames_.empty())
 	{
-		REPORT(1)<< "numberBased dictionary contains " << yellowText(weightBasedDictNames_)<<endREPORT<<endl;	
+		REPORT(1)<< "numberBased dictionary contains " << Yellow_Text(weightBasedDictNames_)<<END_REPORT<<endl;	
 	}
 
+	countDictNames_ = dict_.subDictOrCreate("counting").dictionaryKeywords();
+	if(!countDictNames_.empty())
+	{
+		REPORT(1)<< "counting dictionary contains " << Yellow_Text(countDictNames_)<<END_REPORT<<endl;	
+	}
 
 }
 
@@ -55,43 +60,32 @@ bool pFlow::postprocess::processTimeFolder(real time, const word& tName, const f
 	
 	time_ = time;
 
-	REPORT(0)<<"Working on time folder "<< cyanText(time)<<endREPORT; 
+	REPORT(0)<<"Working on time folder "<< Cyan_Text(time)<<END_REPORT; 
-	timeFolderReposiory_ = 
-		makeUnique<repository>
-		(
-			"timeFolder-"+tName,
-			localFolder,
-			&control_
-		);
 	
-	REPORT(1)<<"Reading pointStructure"<<endREPORT;
+	control_.time().setTime(time);
-	timeFolderReposiory().emplaceObject<pointStructure>(
-		objectFile
-		(
-			pFlow::pointStructureFile__,
-			"",
-			objectFile::READ_ALWAYS,
-			objectFile::WRITE_NEVER
-		));	
 
+	REPORT(1)<<"Reading pointStructure"<<END_REPORT;
+	pointStructure pStruct(control_);
 
-	REPORT(1)<<"Creating mesh and point to cell mapper"<<endREPORT;
+	// first delete the object to remove fields from repository 
+	pointToCell_.reset(nullptr);
+	
+	REPORT(1)<<"Creating mesh and point to cell mapper"<<END_REPORT;
 	pointToCell_ = makeUnique<pointRectCell>(
 		dict_.subDict("rectMesh"),
-		timeFolderReposiory().lookupObject<pointStructure>(pointStructureFile__),
+		pStruct,
-		timeFolderReposiory());
+		control_.time());
 
 	// first numberbased dict 
+	
 	processedFields_.clear();
 	for(word& dictName:numberBasedDictNames_)
 	{
-		
-		
 		auto fieldDict  = dict_.subDict("numberBased").subDict(dictName);
 		auto ppFieldPtr = processField::create(
 			fieldDict,
 			pointToCell_(),
-			timeFolderReposiory());
+			control_.time());
 
 		if(!ppFieldPtr->process())
 		{
@@ -99,11 +93,28 @@ bool pFlow::postprocess::processTimeFolder(real time, const word& tName, const f
 		}
 
 		processedFields_.push_back( ppFieldPtr.release() );
-
 		output<<endl;
-
 	}
 
+	countedVariableNamesList_.clear();
+	countedVairablesLists_.clear();
+
+	for(const auto& countDictName:countDictNames_)
+	{
+		REPORT(1)<<"Processing "<< Yellow_Text("counting."<<countDictName)<<END_REPORT;
+		const dictionary& countDict = dict_.subDict("counting").subDict(countDictName);
+
+		countFields cFields(countDict);
+
+		cFields.process(timeFolderReposiory());
+
+		countedVariableNamesList_.push_back(
+			cFields.variableNames());
+
+		countedVairablesLists_.push_back(cFields.countedValues());
+
+	}
+	output<<"\n";
 
 
 	return true;
@@ -125,7 +136,7 @@ bool pFlow::postprocess::writeToVTK(fileSystem destPath, word bName)const
 
 	if(!vtk) return false;
 
-	REPORT(1)<<"Writing processed fields to vtk file..."<<endREPORT;
+	REPORT(1)<<"Writing processed fields to vtk file..."<<END_REPORT;
 	// mesh
 	pointToCell_->mesh().writeToVtk(vtk());
 

utilities/postprocessPhasicFlow/postprocess.hpp

@@ -46,6 +46,12 @@ protected:
 
 	wordList 				weightBasedDictNames_;
 
+	wordList 				countDictNames_;
+
+	List<wordList>  		countedVariableNamesList_;
+
+	List<uint32List> 		countedVairablesLists_;
+
 	uniquePtr<repository> 	timeFolderReposiory_ {nullptr};
 
 	uniquePtr<pointRectCell>	pointToCell_ {nullptr};

utilities/postprocessPhasicFlow/postprocessPhasicFlow.cpp

@@ -74,7 +74,7 @@ int main(int argc, char** argv )
 	#include "initialize_Control.hpp"
 
 	
-	pFlow::postprocess post(Control);
+	
 		
 	// time folders in case 
 	timeFolder folders(Control);
@@ -93,11 +93,15 @@ int main(int argc, char** argv )
 
 	pFlow::fileSystem destFolder = pFlow::fileSystem(outFolder);
 
+	pFlow::postprocess post(Control);
+
 	do 
 	{
 
-
+		if( !validRange.isMember( folders.time() ) )
-		if( !validRange.isMember( folders.time() ) )continue;
+		{
+			continue;
+		}
 				
 		if( !withZeroFolder && pFlow::equal(folders.time() , 0.0))continue;
 
@@ -109,6 +113,7 @@ int main(int argc, char** argv )
 		}
 
 		
+			
 	}while (folders++);
 
 	#include "finalize.hpp"

utilities/postprocessPhasicFlow/processField.cpp

@@ -107,7 +107,6 @@ pFlow::processField::create(
 		return nullptr;
 	}
 	
-
 	auto method = angleBracketsNames("ProcessField", fType);
 
 	if( dictionaryvCtorSelector_.search(method) )
@@ -115,7 +114,7 @@ pFlow::processField::create(
 		auto objPtr = 
 			dictionaryvCtorSelector_[method]
 			(dict, pToCell, rep);
-		REPORT(2)<<"Processing/creating " << yellowText(dict.name())<< " with model "<<greenText(method)<<"."<<endREPORT;
+		REPORT(2)<<"Processing/creating " << Yellow_Text(dict.name())<< " with model "<<Green_Text(method)<<"."<<END_REPORT;
 		return objPtr;
 	}
 	else

utilities/postprocessPhasicFlow/processField.hpp

@@ -70,6 +70,7 @@ public:
 
 	processField(const dictionary& dict, pointRectCell& pToCell, repository& rep);
 
+	virtual ~processField() = default;
 
 	create_vCtor(
 			processField,
@@ -85,11 +86,21 @@ public:
 		return pointToCell_.mesh();
 	}
 
+	auto& mesh()
+	{
+		return pointToCell_.mesh();
+	}
+
 	const auto& pointToCell()const
 	{
 		return pointToCell_;
 	}
 
+	auto& pointToCell()
+	{
+		return pointToCell_;
+	}
+
 	auto& dict()
 	{
 		return dict_;

utilities/postprocessPhasicFlow/rectMeshField.hpp

@@ -22,17 +22,19 @@ Licence:
 #define __rectMeshField_hpp__
 
 #include "rectangleMesh.hpp"
-#include "baseAlgorithms.hpp"
+#include "ViewAlgorithms.hpp"
 
 namespace pFlow
 {
 
-template<typename T, typename MemorySpace=void>
+template<typename T>
 class rectMeshField
+:
+	public IOobject
 {
 public:
 
-	using viewType = ViewType3D<T,MemorySpace>;
+	using viewType = ViewType3D<T,HostSpace>;
 
 	using memory_space = typename viewType::memory_space;
 
@@ -54,22 +56,49 @@ protected:
 public:
 
 
-	TypeInfoTemplateNV2("rectMeshField", T, memoerySpaceName());
+	TypeInfoTemplateNV111("rectMeshField", T, memoerySpaceName());
 	
-	rectMeshField(const rectangleMesh& mesh, const word& name ,const T& defVal)
+	rectMeshField(rectangleMesh& mesh, const word& name ,const T& defVal)
 	:
+		IOobject(
+			objectFile
+			(
+				name,
+				"",
+				objectFile::READ_NEVER,
+				objectFile::WRITE_NEVER
+			),
+			IOPattern::MasterProcessorOnly,
+			nullptr
+		),
 		mesh_(&mesh),
 		name_(name),
-		field_("pFlow::reactMeshField", mesh_->nx(), mesh_->ny(), mesh_->nz()),
+		field_("reactMeshField."+name, mesh_->nx(), mesh_->ny(), mesh_->nz()),
 		defaultValue_(defVal)
 	{
 		this->fill(defaultValue_);
 	}
 
-	rectMeshField(const rectangleMesh& mesh, const T& defVal)
+	rectMeshField(rectangleMesh& mesh, const T& defVal)
 	:
-		rectMeshField(mesh, "noName", defVal)
+		IOobject(
-	{}
+			objectFile
+			(
+				"noName",
+				"",
+				objectFile::READ_NEVER,
+				objectFile::WRITE_NEVER
+			),
+			IOPattern::MasterProcessorOnly,
+			nullptr
+		),
+		mesh_(&mesh),
+		name_("noName"),
+		field_("reactMeshField.noName", mesh_->nx(), mesh_->ny(), mesh_->nz()),
+		defaultValue_(defVal)
+	{
+		this->fill(defaultValue_);
+	}
 
 	rectMeshField(const rectMeshField&) = default;
 
@@ -144,13 +173,23 @@ public:
 	{
 		pFlow::fill(
 			field_,
-			range(0,mesh_->nx()),
-			range(0,mesh_->ny()),
-			range(0,mesh_->nz()),
 			val
 			);
 	}
 
+	 
+	bool write(iOstream& is, const IOPattern& iop)const override
+	{
+		notImplementedFunction;
+		return true;
+	}
+
+	bool read(iIstream& is, const IOPattern& iop) override
+	{
+		notImplementedFunction;
+		return true;
+	}
+
 	bool read(iIstream& is)
 	{
 		notImplementedFunction;

utilities/postprocessPhasicFlow/rectMeshFieldToVTK.hpp

@@ -26,7 +26,7 @@ namespace pFlow
 {
 
 template<typename T>
-bool convertRectMeshField(iOstream& os, rectMeshField_H<T>& field)
+bool convertRectMeshField(iOstream& os, const rectMeshField_H<T>& field)
 {
 	fatalErrorInFunction<< "this type is not supported "<<
 	field.typeName()<<endl;
@@ -36,7 +36,7 @@ bool convertRectMeshField(iOstream& os, rectMeshField_H<T>& field)
 
 
 template<>
-bool convertRectMeshField(iOstream& os, rectMeshField_H<real>& field)
+bool convertRectMeshField(iOstream& os, const rectMeshField_H<real>& field)
 {
 	
 	os<<"FIELD FieldData 1 " << field.name() << " 1 "<< field.size() << " float\n";
@@ -55,7 +55,7 @@ bool convertRectMeshField(iOstream& os, rectMeshField_H<real>& field)
 }
 
 template<>
-bool convertRectMeshField(iOstream& os, rectMeshField_H<realx3>& field)
+bool convertRectMeshField(iOstream& os, const rectMeshField_H<realx3>& field)
 {
 	
 	os<<"FIELD FieldData 1 " << field.name() << " 3 "<< field.size() << " float\n";
@@ -74,7 +74,7 @@ bool convertRectMeshField(iOstream& os, rectMeshField_H<realx3>& field)
 }
 
 template<>
-bool convertRectMeshField(iOstream& os, rectMeshField_H<int32>& field)
+bool convertRectMeshField(iOstream& os, const rectMeshField_H<int32>& field)
 {
 	
 	os<<"FIELD FieldData 1 " << field.name() << " 1 "<< field.size() << " int\n";

utilities/postprocessPhasicFlow/rectMeshFields.hpp

@@ -27,17 +27,17 @@ namespace pFlow
 {
 
 template<typename T>
-using rectMeshField_H 		= rectMeshField<T,HostSpace>;
+using rectMeshField_H 		= rectMeshField<T>;
 
-using int8RectMeshField_H  = rectMeshField<int8, HostSpace>;
+using int8RectMeshField_H  = rectMeshField<int8>;
 
-using int32RectMeshField_H 	= rectMeshField<int32, HostSpace>;
+using int32RectMeshField_H 	= rectMeshField<int32>;
 
-using int64RectMeshField_H  = rectMeshField<int64, HostSpace>;
+using int64RectMeshField_H  = rectMeshField<int64>;
 
-using realRectMeshField_H	= rectMeshField<real, HostSpace>;
+using realRectMeshField_H	= rectMeshField<real>;
 
-using realx3RectMeshField_H	= rectMeshField<realx3, HostSpace>;
+using realx3RectMeshField_H	= rectMeshField<realx3>;
 
 }
 

utilities/postprocessPhasicFlow/rectangleMesh.cpp

@@ -0,0 +1,52 @@
+#include "rectangleMesh.hpp"
+
+
+pFlow::rectangleMesh::rectangleMesh
+(
+    const box& mshBox, 
+    int32 nx, 
+    int32 ny, 
+    int32 nz,
+    repository* rep
+)
+:
+    IOobject(
+        objectFile
+        (
+            "rectMesh",
+            "",
+            objectFile::READ_NEVER,
+            objectFile::WRITE_NEVER
+        ),
+        IOPattern::MasterProcessorOnly,
+        rep
+    ),
+    meshBox_(mshBox),
+    numCells_(nx, ny, nz)
+{
+    if(mshBox.minPoint()>= mshBox.maxPoint())
+    {
+        fatalErrorInFunction<<"Lower corner point of mesh "<<mshBox.minPoint()<<
+        " confilicts with upper corner point of mesh "<<mshBox.maxPoint()<<endl;
+        fatalExit;
+    }
+
+    numCells_ = max( numCells_ , int32x3(1) );
+
+    dx_ = (mshBox.maxPoint() - mshBox.minPoint())/
+        realx3(numCells_.x_, numCells_.y_, numCells_.z_);
+
+}
+
+pFlow::rectangleMesh::rectangleMesh(const dictionary &dict, repository* rep)
+:
+    rectangleMesh(
+        box(dict),
+        dict.getVal<int32>("nx"),
+		dict.getVal<int32>("ny"),
+		dict.getVal<int32>("nz"),
+        rep
+    )
+{
+
+}

utilities/postprocessPhasicFlow/rectangleMesh.hpp

@@ -21,100 +21,112 @@ Licence:
 #ifndef __rectangleMesh_hpp__
 #define __rectangleMesh_hpp__
 
-#include "cells.hpp"
+#include "types.hpp"
+#include "error.hpp"
+#include "box.hpp"
+#include "IOobject.hpp"
+
+class repository;
 
 namespace pFlow
 {
 
-
-
 class rectangleMesh
 :
-	public cells<int32>
+	public IOobject
 {
+private:
+
+	box 		meshBox_;
+
+	int32x3 	numCells_;
+
+	realx3 		dx_;
 
 public:
 
-	TypeInfoNV("rectangleMesh");
+	TypeInfo("rectangleMesh");
 
-
-	INLINE_FUNCTION_HD
-	rectangleMesh(){};
-
-	INLINE_FUNCTION_HD
 	rectangleMesh(
-		const realx3& minP,
+		const box& mshBox,
-		const realx3& maxP,
 		int32 nx,
 		int32 ny,
-		int32 nz)
+		int32 nz,
-	:
+		repository* rep);
-	cells(
+	
-		box(minP, maxP),
+	rectangleMesh(const dictionary & dict, repository* rep);
-		nx, ny, nz)
-	{}
 
-	INLINE_FUNCTION_H
-	rectangleMesh(const dictionary & dict)
-	:
-	cells(
-		box(
-			dict.getVal<realx3>("min"),
-			dict.getVal<realx3>("max")),
-		dict.getVal<int32>("nx"),
-		dict.getVal<int32>("ny"),
-		dict.getVal<int32>("nz")
-		)
-	{}
 	
-	INLINE_FUNCTION_HD
 	rectangleMesh(const rectangleMesh&) = default;
 
-	INLINE_FUNCTION_HD
+	
 	rectangleMesh& operator = (const rectangleMesh&) = default;
 
-	INLINE_FUNCTION_HD
 	rectangleMesh(rectangleMesh&&) = default;
 
-	INLINE_FUNCTION_HD
+	
 	rectangleMesh& operator = (rectangleMesh&&) = default;
 
-	INLINE_FUNCTION_HD
+	~rectangleMesh() override = default;
-	~rectangleMesh() = default;
 
-
-	INLINE_FUNCTION_HD
 	int64 size()const
 	{
-		return this->totalCells();
+		return numCells_.x_*numCells_.y_*numCells_.z_;
+	}
+
+	int32 nx()const
+	{
+		return numCells_.x();
+	}
+
+	int32 ny()const
+	{
+		return numCells_.y();
+	}
+
+	int32 nz()const
+	{
+		return numCells_.z();
 	}
 
-	INLINE_FUNCTION_HD
 	real cellVol()const
 	{
-		auto [dx,dy,dz] = this->cellSize();
+		return dx_.x_*dx_.y_*dx_.z_;
-		return dx*dy*dz;
 	}
 	
-	INLINE_FUNCTION_H
+	realx3 minPoint()const
-	auto minPoint()const
 	{
-		return domain().minPoint();
+		return meshBox_.minPoint();
 	}
 	
-	INLINE_FUNCTION_H
+	realx3 maxPoint()const
-	auto maxPoint()const
 	{
-		return domain().maxPoint();
+		return meshBox_.maxPoint();
 	}
 
-	bool read(iIstream& is)
+	inline
+	bool isInsideIndex(const realx3 p, int32x3 & ind )const
 	{
+		if(meshBox_.isInside(p))
+		{
+			ind = (p - meshBox_.minPoint())/dx_ ;
+			return true;
+		}
+		else
+		{
+			return false;
+		}
+	}
+
+	bool write(iOstream& is, const IOPattern& iop)const override
+	{
+		notImplementedFunction;
 		return true;
 	}
 
-	bool write(iOstream& os)const
+	bool read(iIstream& is, const IOPattern& iop) override
 	{
+		notImplementedFunction;
 		return true;
 	}
 
@@ -124,7 +136,7 @@ public:
 		os<<"DIMENSIONS "<<nx()+1<<" "<< ny()+1 << " "<< nz()+1 <<endl;
 		
 		auto [x, y , z] = this->minPoint();
-		auto [dx, dy, dz] = this->cellSize();
+		auto [dx, dy, dz] = dx_;
 
 		os<<"X_COORDINATES "<< nx()+1 <<" float\n";
 		for(int32 i=0; i<nx()+1; i++)