before checking parallelIO for dictionary

src/phasicFlow/CMakeLists.txt

@@ -29,6 +29,7 @@ fileSystem/fileSystem.cpp
 processors/processors.cpp
 
 dictionary/dictionary.cpp
+dictionary/fileDictionary.cpp
 dictionary/entry/iEntry.cpp
 dictionary/entry/dataEntry.cpp 
 dictionary/twoPartEntry/twoPartEntry.cpp
@@ -36,6 +37,7 @@ dictionary/twoPartEntry/twoPartEntry.cpp
 containers/Vector/Vectors.cpp
 containers/Field/Fields.cpp
 containers/List/anyList/anyList.cpp
+containers/pointField/pointFields.cpp
 
 Timer/Timer.cpp
 Timer/Timers.cpp
@@ -53,6 +55,8 @@ repository/systemControl/dynamicLinkLibs.cpp
 eventManagement/subscriber.cpp
 eventManagement/observer.cpp
 
+demComponent/demComponent.cpp
+
 structuredData/box/box.cpp
 structuredData/line/line.cpp
 structuredData/infinitePlane/infinitePlane.cpp

src/phasicFlow/containers/VectorHD/VectorSingle.cpp

@@ -0,0 +1,596 @@
+/*------------------------------- 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.
+
+-----------------------------------------------------------------------------*/
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+pFlow::uint32 pFlow::VectorSingle<T,MemorySpace>::changeSize
+(
+    uint32 n, 
+    bool withInit
+)
+{
+    if(n > this->capacity() )
+    {
+        uint32 newCap = evalCapacity(n);
+        changeCapacity(newCap, withInit);
+    }
+    return setSize(n);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+pFlow::uint32 pFlow::VectorSingle<T,MemorySpace>::changeCapacitySize
+(
+    uint32 actualCap, 
+    uint32 n, 
+    bool withInit
+)
+{
+    changeCapacity(actualCap, withInit);
+    return setSize(n);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+void pFlow::VectorSingle<T,MemorySpace>::changeCapacity
+(
+    uint32 actualCap, 
+    bool withInit
+)
+{
+    if(withInit)
+    {
+        resizeInit(view_, actualCap);
+    }
+    else
+    {
+        resizeNoInit(view_, actualCap);
+    }
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+pFlow::uint32 pFlow::VectorSingle<T,MemorySpace>::reallocateCapacitySize
+(
+    uint32 cap, 
+    uint32 s
+)
+{
+    reallocNoInit(view_, cap);
+    return setSize(s);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+pFlow::uint32 pFlow::VectorSingle<T,MemorySpace>::setSize(uint32 n)
+{
+    auto os = size_;
+    size_ = n;
+    return os;
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *//
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>::VectorSingle()
+:
+    VectorSingle("VectorSingle")
+{}
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>::VectorSingle(const word& name)
+:
+    size_(0),
+    view_(name,2)
+{}
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>::VectorSingle
+(
+    const word& name, 
+    uint32 n
+)
+:
+    size_(n),
+    view_(name, evalCapacity(n))
+{}
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>::VectorSingle
+(
+    const word& name, 
+    uint32 n, 
+    const T& val
+)
+:
+    VectorSingle(name, n)
+{
+    assign(n, val); 
+}
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>::VectorSingle
+(
+    const word& name, 
+    uint32 cap, 
+    uint32 n, 
+    RESERVE
+)
+:
+    size_(n),
+    view_(name, cap)
+{}
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>::VectorSingle
+(
+    const word& name, 
+    const std::vector<T> & src
+) 
+:
+    VectorSingle(name)
+{
+    assign(src);
+}
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>::VectorSingle
+(
+    const word& name, 
+    const std::vector<T> & src, 
+    uint32 cap
+) 
+:
+    VectorSingle(name)
+{
+    assign(src, cap);
+}
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>::VectorSingle
+(
+    const VectorSingle& src
+)
+:
+    VectorSingle(src.name(), src)
+{}
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>::VectorSingle
+(
+    const word& name, 
+    const VectorSingle& src
+)
+:
+    VectorSingle(name, src.capacity(), src.size(), RESERVE())
+{
+    copy(deviceVector(), src.deviceVector());	
+}
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>& 
+    pFlow::VectorSingle<T,MemorySpace>::operator = (const VectorSingle& rhs) 
+{
+    if(&rhs == this) return *this;
+    VectorSingle temp(rhs);
+
+    size_ 	= temp.size();	
+    view_   = temp.view_;
+
+    return *this;
+}
+
+template<typename T, typename MemorySpace>
+pFlow::VectorSingle<T,MemorySpace>::~VectorSingle()
+{
+    view_ = viewType();
+    size_ = 0;
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+pFlow::uniquePtr<pFlow::VectorSingle<T,MemorySpace>>
+ pFlow::VectorSingle<T,MemorySpace>::clone() const  
+{
+    return makeUnique<VectorSingle>( this->name(), *this);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+pFlow::VectorSingle<T,MemorySpace>*
+    pFlow::VectorSingle<T,MemorySpace>::clonePtr()const   
+{
+    return new VectorSingle(this->name(), *this);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+pFlow::VectorSingle<T,MemorySpace>& 
+    pFlow::VectorSingle<T,MemorySpace>::VectorField()
+{
+    return *this;
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+const pFlow::VectorSingle<T,MemorySpace>& 
+pFlow::VectorSingle<T,MemorySpace>::VectorField()const
+{
+    return *this;
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+auto pFlow::VectorSingle<T,MemorySpace>::deviceVectorAll() const
+{
+    return view_;
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+auto pFlow::VectorSingle<T,MemorySpace>::deviceVector()const
+{
+    return Kokkos::subview(view_, Kokkos::make_pair(0,int(size_)));
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+auto pFlow::VectorSingle<T,MemorySpace>::hostVectorAll()const
+{
+    auto hView = Kokkos::create_mirror_view(view_);
+    copy(hView, view_);
+    return hView;
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+auto pFlow::VectorSingle<T,MemorySpace>::hostVector()const
+{
+    auto hView = Kokkos::create_mirror_view(deviceVector());
+    copy(hView, deviceVector());
+    return hView;
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+pFlow::word pFlow::VectorSingle<T,MemorySpace>::name()const
+{
+    return view_.label();
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+pFlow::uint32 pFlow::VectorSingle<T,MemorySpace>::size()const
+{
+    return size_;
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+pFlow::uint32 pFlow::VectorSingle<T,MemorySpace>::capacity()const
+{
+    return view_.extent(0);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+bool pFlow::VectorSingle<T,MemorySpace>::empty()const
+{
+    return size_==0;
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+void pFlow::VectorSingle<T,MemorySpace>::reserve(uint32 cap) 
+{
+    changeCapacity(cap);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+void pFlow::VectorSingle<T,MemorySpace>::reallocate(uint32 cap) 
+{
+    reallocateCapacitySize(cap, 0);			
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+void pFlow::VectorSingle<T,MemorySpace>::reallocate
+(
+    uint32 cap, 
+    uint32 newSize
+) 
+{
+    reallocateCapacitySize(cap, newSize);			
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H   
+void pFlow::VectorSingle<T,MemorySpace>::resize(uint32 n)
+{
+    changeSize(n);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H  
+void pFlow::VectorSingle<T,MemorySpace>::resize(uint32 n, const T& val) 
+{
+    assign(n, val);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+void pFlow::VectorSingle<T,MemorySpace>::clear() 
+{
+    size_ = 0;
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+void pFlow::VectorSingle<T,MemorySpace>::fill(const T& val)
+{
+    if(empty())return;
+    pFlow::fill(view_, rangeU32(0 ,size_) ,val);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+void pFlow::VectorSingle<T,MemorySpace>::assign
+(
+    size_t n, 
+    const T& val
+)
+{
+    if( n > capacity() )
+    {
+        reallocateCapacitySize(evalCapacity(n), n);
+    }
+    else
+    {
+        setSize(n);	
+    }
+    this->fill(val);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+void pFlow::VectorSingle<T,MemorySpace>::assign
+(
+    const std::vector<T>& src, 
+    uint32 cap
+) 
+{
+    uint32 srcSize = src.size();
+    
+    if(cap != capacity())
+    {
+        reallocateCapacitySize(cap, srcSize);
+    }
+    else
+    {
+        setSize(srcSize);
+    }
+    
+    // - unmanaged view in the host
+    hostViewType1D<const T> temp(src.data(), srcSize );
+    copy(deviceVector(), temp);
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+void pFlow::VectorSingle<T,MemorySpace>::assign
+(
+    const std::vector<T>& src
+) 
+{
+    assign(src, src.capacity());
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+auto pFlow::VectorSingle<T,MemorySpace>::getSpan()
+{
+    return span<T>(view_.data(), size());
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H 
+auto pFlow::VectorSingle<T,MemorySpace>::getSpan()const
+{
+    return span<const T>(view_.data(), this->size());
+}
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+bool pFlow::VectorSingle<T,MemorySpace>::insertSetElement(uint32IndexContainer indices, const T& val)
+{
+	if(indices.empty()) return true;
+	
+	auto maxInd = indices.max();
+
+	if(this->empty() || maxInd > size()-1 )
+	{
+		resize(maxInd+1);
+	}
+
+	using policy = Kokkos::RangePolicy<execution_space,Kokkos::IndexType<uint32>>;
+
+	if constexpr( isDeviceAccessible_ )
+	{
+		auto v = view_;
+		auto ind = indices.deviceView();
+
+		Kokkos::parallel_for(
+			"VectorSingle::insertSetElement",
+			policy(0,indices.size()),
+			LAMBDA_HD(uint32 i){
+				v[ind[i]]= val;
+			});
+
+		Kokkos::fence();
+	}
+	else
+	{
+
+		auto v = view_;
+		auto ind = indices.hostView();
+		
+		Kokkos::parallel_for(
+			"VectorSingle::insertSetElement",
+			policy(0,indices.size()),
+			LAMBDA_HD(uint32 i){
+				v[ind[i]]= val;
+			});
+
+		Kokkos::fence();
+
+	}
+
+	return true;
+}
+
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+bool pFlow::VectorSingle<T,MemorySpace>::insertSetElement
+(
+	const uint32IndexContainer indices, 
+	const std::vector<T>& vals
+)
+{
+
+	if(indices.size() == 0)return true;
+	if(indices.size() != vals.size())return false;
+
+	auto maxInd = indices.max(); 
+	
+	if(this->empty() || maxInd > size()-1 )
+	{
+		resize(maxInd+1);
+	} 
+	
+	using policy = Kokkos::RangePolicy<execution_space,Kokkos::IndexType<uint32>>;
+
+	hostViewType1D<const T> hVals( vals.data(), vals.size());
+
+	if constexpr( isDeviceAccessible_ )
+	{
+		deviceViewType1D<T> dVals("dVals", indices.size());
+		copy(dVals, hVals);
+		auto dVec = view_;
+		auto ind = indices.deviceView();
+
+		Kokkos::parallel_for(
+			"VectorSingle::insertSetElement",
+			policy(0,indices.size()), LAMBDA_HD(int32 i){	
+				dVec(ind(i))= dVals(i);}
+			);
+
+		Kokkos::fence();
+
+	}
+	else
+	{
+		auto dVec = view_;
+		auto ind = indices.hostView();
+
+		Kokkos::parallel_for(
+			"VectorSingle::insertSetElement",
+			policy(0,indices.size()), LAMBDA_HD(int32 i){	
+				dVec(ind(i))= hVals(i);}
+			);
+
+		Kokkos::fence();
+
+	}
+		
+	return true;
+}
+
+
+template<typename T, typename MemorySpace>
+INLINE_FUNCTION_H
+bool pFlow::VectorSingle<T,MemorySpace>::reorderItems(uint32IndexContainer indices)
+{
+	if(indices.size() == 0)
+	{
+		setSize(0);
+		return true;
+	}
+	
+	auto maxInd = indices.max();
+
+	if(maxInd >= this->size())
+	{
+		fatalErrorInFunction<<"In reordering the VectorSingle ("
+		<< this->name()<< ") maximum index ("<< maxInd << 
+		") exceeds the size of the vector ("	<< this->size()<<")"<<endl;
+		return false;
+	}
+
+	uint32 		newSize = indices.size();
+	
+	setSize(newSize);
+
+	viewType 	sortedView(this->name(), newSize);
+
+	using policy = Kokkos::RangePolicy< execution_space,Kokkos::IndexType<uint32>>;
+
+	if constexpr( isDeviceAccessible_)
+	{
+		auto d_indices = indices.deviceView();
+		auto d_view = view_;
+		
+		Kokkos::parallel_for
+		(
+			"VectorSingle::sortItems",
+			policy(0,newSize), 
+			LAMBDA_HD(uint32 i)
+			{	
+				sortedView(i) = d_view(d_indices(i));
+			}
+		);
+
+		Kokkos::fence();
+
+	}
+	else
+	{
+		auto h_indices = indices.hostView();
+		auto d_view = view_;
+		
+		Kokkos::parallel_for
+		(
+			"VectorSingle::sortItems",
+			policy(0,newSize), 
+			LAMBDA_HD(uint32 i)
+			{	
+				sortedView(i) = d_view(h_indices(i));
+			}
+		);
+		
+		Kokkos::fence();
+	}
+
+	copy(deviceVector(), sortedView);
+	
+	return true;
+}

src/phasicFlow/containers/VectorHD/VectorSingle.hpp

@@ -17,23 +17,18 @@ Licence:
   implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
 
 -----------------------------------------------------------------------------*/
-
-
 #ifndef __VectorSingle_hpp__
 #define __VectorSingle_hpp__
 
 #include <vector>
 
-
 #include "globalSettings.hpp"
+#include "phasicFlowKokkos.hpp"
 #include "types.hpp"
 #include "error.hpp"
 #include "indexContainer.hpp"
-#include "iOstream.hpp"
+#include "streams.hpp"
-#include "iIstream.hpp"
 #include "span.hpp"
-#include "Vector.hpp"
-#include "phasicFlowKokkos.hpp"
 #include "dataIO.hpp"
 
 
@@ -49,8 +44,6 @@ namespace pFlow
 template<typename T, typename MemorySpace>
 class VectorSingle;
 
-
-
 template<typename T, typename MemorySpace=void>
 class VectorSingle
 {
@@ -103,7 +96,6 @@ protected:
 		/// Is the memory of this vector accessible from Host
 		static constexpr 
 		bool  isHostAccessible_ = isHostAccessible<execution_space>();
-			//Kokkos::SpaceAccessibility<execution_space,HostSpace>::accessible;
 
 		/// Is the memory of this vector accessiple from Divce
 		static constexpr
@@ -123,144 +115,60 @@ protected:
 			return static_cast<uint32>(n*growthFactor_+1);
 		}
 
-	/// Change size to n and preserve the conetent if realloc 
+	/// @brief Change size to n and preserve the conetent if realloc 
 	/// occurs 
 	INLINE_FUNCTION_H 
-	uint32 changeSize(uint32 n, bool withInit= false)
+	uint32 changeSize(uint32 n, bool withInit= false);
-	{
+
-		if(n > this->capacity() )
+    /// @brief Change the size and capacity of Vector
-		{
+	INLINE_FUNCTION_H
-			uint32 newCap = evalCapacity(n);
+	uint32 changeCapacitySize(uint32 actualCap, uint32 n, bool withInit= false);
-			changeCapacity(newCap, withInit);
-		}
-		return setSize(n);
-	}
 	
 	INLINE_FUNCTION_H
-	uint32 changeCapacitySize(uint32 actualCap, uint32 n, bool withInit= false)
+	void changeCapacity(uint32 actualCap, bool withInit= false);
-	{
-		changeCapacity(actualCap, withInit);
-		return setSize(n);
-	}
 	
 	INLINE_FUNCTION_H
-	void changeCapacity(uint32 actualCap, bool withInit= false)
+	uint32 reallocateCapacitySize(uint32 cap, uint32 s);
-	{
-		if(withInit)
-		{
-			resizeInit(view_, actualCap);
-		}
-		else
-		{
-			resizeNoInit(view_, actualCap);
-		}
-	}
 	
 	INLINE_FUNCTION_H 
-	uint32 reallocateCapacitySize(uint32 cap, uint32 s)
+	uint32 setSize(uint32 n);
-	{
-		reallocNoInit(view_, cap);
-		return setSize(s);
-	}
-
-	INLINE_FUNCTION_H 
-	uint32 setSize(uint32 n)
-	{
-		auto os = size_;
-		size_ = n;
-		return os;
-	}
 	
 public:
 
-	// - type info
+	/// Type info
 	TypeInfoTemplateNV2("VectorSingle", T, memoerySpaceName());
 
 	//// - Constructors
 
 		/// Empty vector
-		VectorSingle()
+		VectorSingle();
-		:
-			VectorSingle("VectorSingle")
-		{}
-
-		/// Empty vector with a name 
-		VectorSingle(const word& name)
-		:
-			size_(0),
-			view_(name,2)
-		{}
 		
+		/// Empty vector with a name (capacity = 2) 
+		VectorSingle(const word& name);
 		
 		/// Vector with name and size n
-		VectorSingle(const word& name, uint32 n)
+		VectorSingle(const word& name, uint32 n);
-		:
-			size_(n),
-			view_(name, evalCapacity(n))
-		{}
-
-		
 		
 		/// Vector with name, size and value
-		VectorSingle(const word& name, uint32 n, const T& val)
+		VectorSingle(const word& name, uint32 n, const T& val);
-		:
-			VectorSingle(name, n)
-		{
-			assign(n, val); 
-		}
-
 
 		/// Vector with name, size (n) and reserved capacity 
-		VectorSingle(const word& name, uint32 cap, uint32 n, RESERVE )
+		VectorSingle(const word& name, uint32 cap, uint32 n, RESERVE );
-		:
-			size_(n),
-			view_(name, cap)
-		{}
 		
 		/// Construct with a name and form std::vector (host memory)
-		VectorSingle(const word& name, const std::vector<T> & src) 
+		VectorSingle(const word& name, const std::vector<T> & src); 
-		:
-			VectorSingle(name)
-		{
-			assign(src);
-		}
 
-		/// Construct with a name and form std::vector (host memory)
+		/// Construct with a name and form std::vector (host memory) and with a desired capacity. 
-		/// and with a desired capacity. 
+		VectorSingle(const word& name, const std::vector<T> & src, uint32 cap);
-		VectorSingle(const word& name, const std::vector<T> & src, uint32 cap) 
-		:
-			VectorSingle(name)
-		{
-			assign(src, cap);
-		}
 		
 		/// Copy construct (performs deep copy)
-		VectorSingle(const VectorSingle& src)
+		VectorSingle(const VectorSingle& src);
-		:
-			VectorSingle(src.name(), src)
-		{	
-			//copy(deviceVector(), src.deviceVector());
-		}
 		
     	/// Copy construct with a new name (perform deep copy)
-		VectorSingle(const word& name, const VectorSingle& src)
+		VectorSingle(const word& name, const VectorSingle& src);
-		:
-			VectorSingle(name, src.capacity(), src.size(), RESERVE())
-		{
-			copy(deviceVector(), src.deviceVector());	
-		}
 		
 		/// Copy assignment (perform deep copy from rhs to *this)
-		VectorSingle& operator = (const VectorSingle& rhs) 
+		VectorSingle& operator = (const VectorSingle& rhs) ;
-		{
-			if(&rhs == this) return *this;
-			VectorSingle temp(rhs);
-
-			size_ 	= temp.size();	
-			view_   = temp.view_;
-
-			return *this;
-		}
 		
 		/// Move construct
 		VectorSingle(VectorSingle&&) = default;   
@@ -268,210 +176,115 @@ public:
 		/// Move assignment
 		VectorSingle& operator= (VectorSingle&&) = default; 
 
+		/// @brief Descructor
         /// This may not destroy the underlying memory, sice view is 
 		/// shared_ptr and maybe referenced by another object too 
-		~VectorSingle()
+		~VectorSingle();
-		{
-			view_ = viewType();
-			size_ = 0;
-		}
 
 		/// Clone as a uniquePtr (perform deep copy)
 		INLINE_FUNCTION_H 
-		uniquePtr<VectorSingle> clone() const  
+		uniquePtr<VectorSingle> clone() const;
-		{
-			return makeUnique<VectorSingle>( this->name(), *this);
-		}
 
 		/// Clone as a pointer (perform deep copy)
 		INLINE_FUNCTION_H
-		VectorSingle* clonePtr()const   
+		VectorSingle* clonePtr()const;
-		{
-			return new VectorSingle(this->name(), *this);
-		}
 	
 	//// - Methods
 
 		/// Return *this
 		INLINE_FUNCTION_H
-		VectorType& VectorField()
+		VectorType& VectorField();
-		{
-			return *this;
-		}
 
 		/// Return *this
 		INLINE_FUNCTION_H
-		const VectorType& VectorField()const
+		const VectorType& VectorField()const;
-		{
-			return *this;
-		}
 
 		/// Device vector range [0,capcity)
 		INLINE_FUNCTION_H 
-		viewType deviceVectorAll() const{
+		auto deviceVectorAll() const;
-			return view_;
-		}
 
 		///  Device vector range [0, size)
 		INLINE_FUNCTION_H 
-		viewType deviceVector()const
+		auto deviceVector()const;
-		{
-			return Kokkos::subview(view_, Kokkos::make_pair(0,int(size_)));
-		}
 
 		/// Return a vector accessible on Host in range [0,capacity)
 		INLINE_FUNCTION_H 
-		auto hostVectorAll()const
+		auto hostVectorAll()const;
-		{
+		
-			auto hView = Kokkos::create_mirror_view(view_);
-			copy(hView, view_);
-			return hView;
-		}
 
 		/// Return a vector accessible on Host in range [0,size)
 		INLINE_FUNCTION_H 
-		auto hostVector()const
+		auto hostVector()const;
-		{
-			auto hView = Kokkos::create_mirror_view(deviceVector());
-			copy(hView, deviceVector());
-			return hView;
-		}
 
 		/// Name of the vector 
 		INLINE_FUNCTION_H 
-		const word name()const
+		word name()const;
-		{
+		
-			return view_.label();
-		}
 		
 		/// Size of the vector 
 		INLINE_FUNCTION_H 
-		uint32 size()const
+		uint32 size()const;
-		{
+		
-			return size_;
-		}
 
 		// Capcity of the vector 
 		INLINE_FUNCTION_H 
-		uint32 capacity()const
+		uint32 capacity()const;
-		{
-			return view_.extent(0);
-		}
 		
 		/// If vector is empty 
 		INLINE_FUNCTION_H 
-		bool empty()const
+		bool empty()const;
-		{
-			return size_==0;
-		}
 		
 		/// Reserve capacity for vector 
 		/// Preserve the content.
 		INLINE_FUNCTION_H 
-		void reserve(uint32 cap) 
+		void reserve(uint32 cap); 
-		{
-			changeCapacity(cap);
-		}
 
 		/// Reallocate memory to new cap and set size to 0. 
 		INLINE_FUNCTION_H 
-		void reallocate(uint32 cap) 
+		void reallocate(uint32 cap); 
-		{
-			reallocateCapacitySize(cap, 0);			
-		}
 		
 		/// Reallocate memory to new cap and set size to newSize. 
 		/// Do not preserve the content
 		INLINE_FUNCTION_H 
-		void reallocate(uint32 cap, uint32 newSize) 
+		void reallocate(uint32 cap, uint32 newSize); 
-		{
-			reallocateCapacitySize(cap, newSize);			
-		}
 		
 		/// Resize the vector and preserve the content
 		INLINE_FUNCTION_H   
-		void resize(uint32 n){
+		void resize(uint32 n);
-			changeSize(n);
-		}
 
 		/// Resize the vector and assign the value to it. 
 		INLINE_FUNCTION_H  
-		void resize(uint32 n, const T& val) {
+		void resize(uint32 n, const T& val);
-			assign(n, val);
-		}
 
 		/// Clear the vector, but keep the allocated memory unchanged 
 		INLINE_FUNCTION_H 
-		void clear() 
+		void clear();
-		{
-			 size_ = 0;
-		}
 		
 		/// Fill the range [0,size) with val
 		INLINE_FUNCTION_H 
-		void fill(const T& val)
+		void fill(const T& val);
-		{
-			if(empty())return;
-			pFlow::fill(view_, rangeU32(0 ,size_) ,val);
-		}
 		
 		/// Change size of the vector and assign val to vector and 
 		INLINE_FUNCTION_H
-		void assign(size_t n, const T& val)
+		void assign(size_t n, const T& val);
-		{
-			if( n > capacity() )
-			{
-				reallocateCapacitySize(evalCapacity(n), n);
-			}
-			else
-			{
-				setSize(n);	
-			}
-			this->fill(val);
-		}
-
 		
 		/// Assign source vector with specified capacity.
 		/// The size of *this becomes the size of src. 
 		INLINE_FUNCTION_H 
-		void assign(const std::vector<T>& src, uint32 cap) 
+		void assign(const std::vector<T>& src, uint32 cap); 
-		{
-			uint32 srcSize = src.size();
 		
-			if(cap != capacity())
-			{
-				reallocateCapacitySize(cap, srcSize);
-			}
-			else
-			{
-				setSize(srcSize);
-			}
-			
-			// - unmanaged view in the host
-			hostViewType1D<const T> temp(src.data(), srcSize );
-			copy(deviceVector(), temp);
-		}
 
 		/// Assign source vector.
 		/// The size of *this becomes the size of src. 
 		/// The capacity of *this becomes the capacity of src.
 		INLINE_FUNCTION_H 
-		void assign(const std::vector<T>& src) 
+		void assign(const std::vector<T>& src);
-		{
-			assign(src, src.capacity());
-		}
 		
 		INLINE_FUNCTION_H
-		auto getSpan()
+		auto getSpan();
-		{
-			return span<T>(view_.data(), size());
-		}
 				
 		INLINE_FUNCTION_H 
-		auto getSpan()const
+		auto getSpan()const;
-		{
-			return span<const T>(view_.data(), this->size());
-		}
 		
 		INLINE_FUNCTION_H
 		bool insertSetElement(uint32IndexContainer indices, const T& val);
@@ -482,40 +295,8 @@ public:
 		INLINE_FUNCTION_H
 		bool reorderItems(uint32IndexContainer indices);
 
-
+		/// @brief push a new element at the end (host call only)
-		/*INLINE_FUNCTION_H
+		///  resize if necessary and works on host accessible vector.
-		bool append(const deviceViewType1D<T>& dVec, size_t numElems)
-		{
-
-			if(numElems == 0 )return true;
-			auto oldSize = size_;
-			auto newSize = size_ + numElems;
-
-			if(this->empty() || newSize > capacity() )
-			{
-				resize(newSize);
-			}
-			else
-			{
-				size_ = size_+numElems;
-			}
-			
-			auto dSubView = Kokkos::subview(view_, Kokkos::make_pair(oldSize, newSize)); 
-			copy(dSubView, dVec);
-
-			return true;
-		}
-
-		INLINE_FUNCTION_H
-		bool append(const VectorSingle& Vec)
-		{
-			return append(Vec.deviceVector(), Vec.size());
-		}*/
-
-		// - host calls only
-		//   push a new element at the end
-		//   resize if necessary
-		//   works on host accessible vector 
 		template<bool Enable = true>
 		INLINE_FUNCTION_H
 		typename std::enable_if<isHostAccessible_ && Enable, void>::type
@@ -535,7 +316,7 @@ public:
 			return view_.data();
 		}
 
-		/// Return begin iterator. it works when host is accessible.
+		/// Return begin iterator. It works when devices is host accessible.
 		template<bool Enable = true>
 		INLINE_FUNCTION_H 
 		typename std::enable_if_t<isHostAccessible_ && Enable,iterator>
@@ -666,8 +447,39 @@ inline iOstream& operator << (iOstream& os, const VectorSingle<T, MemorySpace>&
 } // - pFlow
 
 #include "VectorSingleAlgorithms.hpp"
-#include "VectorSingleI.hpp"
+#include "VectorSingle.cpp"
+
+
 
 
 #endif //__VectorSingle_hpp__
 
+
+/*INLINE_FUNCTION_H
+		bool append(const deviceViewType1D<T>& dVec, size_t numElems)
+		{
+
+			if(numElems == 0 )return true;
+			auto oldSize = size_;
+			auto newSize = size_ + numElems;
+
+			if(this->empty() || newSize > capacity() )
+			{
+				resize(newSize);
+			}
+			else
+			{
+				size_ = size_+numElems;
+			}
+			
+			auto dSubView = Kokkos::subview(view_, Kokkos::make_pair(oldSize, newSize)); 
+			copy(dSubView, dVec);
+
+			return true;
+		}
+
+		INLINE_FUNCTION_H
+		bool append(const VectorSingle& Vec)
+		{
+			return append(Vec.deviceVector(), Vec.size());
+		}*/

src/phasicFlow/containers/VectorHD/VectorSingleAlgorithms.hpp

@@ -17,8 +17,6 @@ Licence:
   implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
 
 -----------------------------------------------------------------------------*/
-
-
 #ifndef __VectorSingleMath_hpp__
 #define __VectorSingleMath_hpp__
 

src/phasicFlow/containers/VectorHD/VectorSingleI.hpp

@@ -1,197 +0,0 @@
-/*------------------------------- 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.
-
------------------------------------------------------------------------------*/
-
-
-template<typename T, typename MemorySpace>
-INLINE_FUNCTION_H
-bool pFlow::VectorSingle<T,MemorySpace>::insertSetElement(uint32IndexContainer indices, const T& val)
-{
-	if(indices.empty()) return true;
-	
-	auto maxInd = indices.max();
-
-	if(this->empty() || maxInd > size()-1 )
-	{
-		resize(maxInd+1);
-	}
-
-	using policy = Kokkos::RangePolicy<execution_space,Kokkos::IndexType<uint32>>;
-
-	if constexpr( isDeviceAccessible_ )
-	{
-		auto v = view_;
-		auto ind = indices.deviceView();
-
-		Kokkos::parallel_for(
-			"VectorSingle::insertSetElement",
-			policy(0,indices.size()),
-			LAMBDA_HD(uint32 i){
-				v[ind[i]]= val;
-			});
-
-		Kokkos::fence();
-	}
-	else
-	{
-
-		auto v = view_;
-		auto ind = indices.hostView();
-		
-		Kokkos::parallel_for(
-			"VectorSingle::insertSetElement",
-			policy(0,indices.size()),
-			LAMBDA_HD(uint32 i){
-				v[ind[i]]= val;
-			});
-
-		Kokkos::fence();
-
-	}
-
-	return true;
-}
-
-
-template<typename T, typename MemorySpace>
-INLINE_FUNCTION_H
-bool pFlow::VectorSingle<T,MemorySpace>::insertSetElement
-(
-	const uint32IndexContainer indices, 
-	const std::vector<T>& vals
-)
-{
-
-	if(indices.size() == 0)return true;
-	if(indices.size() != vals.size())return false;
-
-	auto maxInd = indices.max(); 
-	
-	if(this->empty() || maxInd > size()-1 )
-	{
-		resize(maxInd+1);
-	} 
-	
-	using policy = Kokkos::RangePolicy<execution_space,Kokkos::IndexType<uint32>>;
-
-	hostViewType1D<const T> hVals( vals.data(), vals.size());
-
-	if constexpr( isDeviceAccessible_ )
-	{
-		deviceViewType1D<T> dVals("dVals", indices.size());
-		copy(dVals, hVals);
-		auto dVec = view_;
-		auto ind = indices.deviceView();
-
-		Kokkos::parallel_for(
-			"VectorSingle::insertSetElement",
-			policy(0,indices.size()), LAMBDA_HD(int32 i){	
-				dVec(ind(i))= dVals(i);}
-			);
-
-		Kokkos::fence();
-
-	}
-	else
-	{
-		auto dVec = view_;
-		auto ind = indices.hostView();
-
-		Kokkos::parallel_for(
-			"VectorSingle::insertSetElement",
-			policy(0,indices.size()), LAMBDA_HD(int32 i){	
-				dVec(ind(i))= hVals(i);}
-			);
-
-		Kokkos::fence();
-
-	}
-		
-	return true;
-}
-
-
-template<typename T, typename MemorySpace>
-INLINE_FUNCTION_H
-bool pFlow::VectorSingle<T,MemorySpace>::reorderItems(uint32IndexContainer indices)
-{
-	if(indices.size() == 0)
-	{
-		setSize(0);
-		return true;
-	}
-	
-	auto maxInd = indices.max();
-
-	if(maxInd >= this->size())
-	{
-		fatalErrorInFunction<<"In reordering the VectorSingle ("
-		<< this->name()<< ") maximum index ("<< maxInd << 
-		") exceeds the size of the vector ("	<< this->size()<<")"<<endl;
-		return false;
-	}
-
-	uint32 		newSize = indices.size();
-	
-	setSize(newSize);
-
-	viewType 	sortedView(this->name(), newSize);
-
-	using policy = Kokkos::RangePolicy< execution_space,Kokkos::IndexType<uint32>>;
-
-	if constexpr( isDeviceAccessible_)
-	{
-		auto d_indices = indices.deviceView();
-		auto d_view = view_;
-		
-		Kokkos::parallel_for
-		(
-			"VectorSingle::sortItems",
-			policy(0,newSize), 
-			LAMBDA_HD(uint32 i)
-			{	
-				sortedView(i) = d_view(d_indices(i));
-			}
-		);
-
-		Kokkos::fence();
-
-	}
-	else
-	{
-		auto h_indices = indices.hostView();
-		auto d_view = view_;
-		
-		Kokkos::parallel_for
-		(
-			"VectorSingle::sortItems",
-			policy(0,newSize), 
-			LAMBDA_HD(uint32 i)
-			{	
-				sortedView(i) = d_view(h_indices(i));
-			}
-		);
-		
-		Kokkos::fence();
-	}
-
-	copy(deviceVector(), sortedView);
-	
-	return true;
-}

src/phasicFlow/containers/pointField/boundaryField.cpp

@@ -0,0 +1,32 @@
+/*------------------------------- 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.
+
+-----------------------------------------------------------------------------*/
+
+template<template<class, class> class VectorField, class T, class MemorySpace>
+pFlow::boundaryField<VectorField, T, MemorySpace>::boundaryField
+(
+	const boundaryBase& boundary,
+	FieldType &internal
+)
+:
+	observer(),
+	boundary_(boundary),
+	internal_(internal)
+{}
+

src/phasicFlow/containers/pointField/boundaryField.hpp

@@ -0,0 +1,65 @@
+/*------------------------------- 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 __boundaryField_hpp__
+#define __boundaryField_hpp__
+
+#include "observer.hpp"
+#include "boundaryBase.hpp"
+#include "Field.hpp"
+
+namespace pFlow
+{
+
+template< template<class, class> class VectorField, class T, class MemorySpace = void>
+class boundaryField
+:
+    public observer 
+{
+public:
+	using FieldType = Field<VectorField, T, MemorySpace>;
+protected:
+
+    const boundaryBase&     boundary_;
+
+    /// @brief a ref to the internal field 
+    FieldType& internal_;
+
+public:
+
+   boundaryField(const boundaryBase& boundary, FieldType& internal);
+
+	bool update
+	(
+		const message& msg, 
+    	const anyList& varList
+	) override
+    {
+		notImplementedFunction;
+		return false;
+	}
+
+};
+
+}
+
+#include "boundaryField.cpp"
+
+#endif
+

src/phasicFlow/containers/pointField/boundaryFieldList.hpp

@@ -0,0 +1,72 @@
+/*------------------------------- 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 __boundaryFieldList_hpp__
+#define __boundaryFieldList_hpp__
+
+#include "boundaryField.hpp"
+#include "boundaryList.hpp"
+#include "ListPtr.hpp"
+
+namespace pFlow
+{
+
+template< template<class, class> class VectorField, class T, class MemorySpace = void >
+class boundaryFieldList
+:
+    ListPtr<boundaryField<VectorField, T, MemorySpace>>
+{
+public:
+
+    using boundaryFieldType = boundaryField<VectorField, T, MemorySpace>;
+
+    using FieldType = typename boundaryFieldType::FieldType;
+
+protected:
+
+    const boundaryList&         boundaries_;
+
+public:
+
+    boundaryFieldList(const boundaryList& boundaries, FieldType& internal)
+    :
+        ListPtr<boundaryFieldType>(boundaries.size()),
+        boundaries_(boundaries)
+    {
+        for(auto i=0; i<boundaries.size(); i++)
+        {
+            this->set
+            (
+                i, 
+                makeUnique<boundaryFieldType>
+                (
+                    boundaries_[i], 
+                    internal
+                )
+            );
+        }
+    }
+
+
+
+};
+
+}
+
+#endif

src/phasicFlow/containers/pointField/pointField.cpp

@@ -18,8 +18,20 @@ Licence:
 
 -----------------------------------------------------------------------------*/
 
-
 template<template<class, class> class VectorField, class T, class MemorySpace>
+pFlow::pointField<VectorField, T, MemorySpace>::pointField
+(
+    const pointStructure& pStruct, 
+    const T& defVal, 
+    message msg
+)
+:
+	boundaryFieldList_(pStruct.boundaries(), *this),
+	pStruct_(pStruct),
+	defaultValue_(defVal)
+{}
+
+/*template<template<class, class> class VectorField, class T, class MemorySpace>
 bool pFlow::pointField<VectorField, T, MemorySpace>::readPointField
 (
 	iIstream& is
@@ -113,9 +125,7 @@ bool pFlow::pointField<VectorField, T, MemorySpace>::update(const eventMessage&
  
 	if( msg.isDeleted() )
 	{
-		/*const auto& dp = pStruct_.markedDeletePoints();
+		
-		return this->deleteElement(dp);
-		notImplementedFunction;*/
 	}
 	else if( msg.isInsert())
 	{
@@ -132,5 +142,5 @@ bool pFlow::pointField<VectorField, T, MemorySpace>::update(const eventMessage&
 	}
 	
 	return true;
-}
+}*/
 

src/phasicFlow/containers/pointField/pointField.hpp

@@ -21,20 +21,18 @@ Licence:
 #ifndef __pointField_hpp__
 #define __pointField_hpp__
 
-
+#include "boundaryFieldList.hpp"
-#include "Field.hpp"
 #include "pointStructure.hpp"
-#include "error.hpp"
+#include "Field.hpp"
+#include "observer.hpp"
 
 namespace pFlow
 {
 
-
-
 template<template<class, class> class VectorField, class T, class MemorySpace=void>
 class pointField
 :
-	public eventObserver,
+	public observer,
 	public Field<VectorField, T, MemorySpace>
 {
 public:
@@ -43,28 +41,35 @@ public:
 	
 	using FieldType         = Field<VectorField, T, MemorySpace>;
 
+    using boundaryFieldListType = boundaryFieldList<VectorField, T, MemorySpace>;
+  	
   	using VectorType  		= typename FieldType::VectorType;
 
 	using iterator        	= typename FieldType::iterator;
 
-  	using constIterator   	= typename FieldType::constIterator;
+  	using const_iterator   	= typename FieldType::const_iterator;
 
 	using reference       	= typename FieldType::reference;
   	
-  	using constReference  	= typename FieldType::constReference;
+  	using const_reference  	= typename FieldType::const_reference;
 
-	using valueType       	= typename FieldType::valueType;
+	using value_type       	= typename FieldType::value_type;
   	
   	using pointer         	= typename FieldType::pointer;
   	
-  	using constPointer    	= typename FieldType::constPointer;
+  	using const_pointer    	= typename FieldType::const_pointer;
 
 protected:
 
-	////- data members
+    //// - data members 
+
+        /// @brief list of boundaries
+        boundaryFieldListType       boundaryFieldList_;
+        
+		/// @brief  refrence to point structure 
 		const pointStructure&       pStruct_;
 
-		// - value when a new item is added to field
+		/// @brief value when a new item is added to field
 		T 	defaultValue_;
 
 
@@ -77,10 +82,13 @@ public:
 	//// - Constructors
 
 		// - construct a field from pointStructure and set defaultValue_ and field value to defVal
-		pointField( const pointStructure& pStruct, const T& defVal, bool subscribe = true);
+		pointField( 
+            const pointStructure& pStruct, 
+            const T& defVal, 
+            message msg);
 
 		// - construct from iIOEntity, pointStructure and a value
-		pointField( const pointStructure& pStruct, const T& val, const T& defVal, bool subscribe = true);
+		/*pointField( const pointStructure& pStruct, const T& val, const T& defVal, bool subscribe = true);
 
 		// - construct from another pointField
 		//   subscribe to events if true
@@ -159,10 +167,10 @@ public:
 		bool write(iOstream& os)const
 		{
 			return writePointField(os);
-		}
+		}*/
 };
 
-template<template<class, class> class VectorField, class T, class MemorySpace>
+/*template<template<class, class> class VectorField, class T, class MemorySpace>
 iIstream& operator >> (iIstream & is, pointField<VectorField, T, MemorySpace> & pF )
 {
 	if( !pF.read(is))
@@ -186,11 +194,11 @@ iOstream& operator << (iOstream& os, const pointField<VectorField, T, MemorySpac
 	}
 
 	return os;
-}
+}*/
 
 }
 
 #include "pointField.cpp"
-#include "pointFieldAlgorithms.hpp"
+//#include "pointFieldAlgorithms.hpp"
 
 #endif //  __pointField_hpp__

src/phasicFlow/containers/pointField/pointField_old.hpp

@@ -0,0 +1,196 @@
+/*------------------------------- 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 __pointField_hpp__
+#define __pointField_hpp__
+
+
+#include "Field.hpp"
+#include "pointStructure.hpp"
+#include "error.hpp"
+
+namespace pFlow
+{
+
+
+
+template<template<class, class> class VectorField, class T, class MemorySpace=void>
+class pointField
+:
+	public eventObserver,
+	public Field<VectorField, T, MemorySpace>
+{
+public:
+	
+	using pointFieldType 	= pointField<VectorField, T, MemorySpace>;
+	
+	using FieldType         = Field<VectorField, T, MemorySpace>;
+  	
+  	using VectorType  		= typename FieldType::VectorType;
+
+	using iterator        	= typename FieldType::iterator;
+
+  	using constIterator   	= typename FieldType::constIterator;
+
+	using reference       	= typename FieldType::reference;
+  	
+  	using constReference  	= typename FieldType::constReference;
+
+	using valueType       	= typename FieldType::valueType;
+  	
+  	using pointer         	= typename FieldType::pointer;
+  	
+  	using constPointer    	= typename FieldType::constPointer;
+
+protected:
+
+	////- data members
+		const pointStructure& pStruct_;
+
+		// - value when a new item is added to field
+		T 	defaultValue_;
+
+
+public:
+ 
+	// - type info
+	TypeInfoTemplateNV2("pointField", T, VectorType::memoerySpaceName());
+
+	
+	//// - Constructors
+
+		// - construct a field from pointStructure and set defaultValue_ and field value to defVal
+		pointField( const pointStructure& pStruct, const T& defVal, bool subscribe = true);
+
+		// - construct from iIOEntity, pointStructure and a value
+		pointField( const pointStructure& pStruct, const T& val, const T& defVal, bool subscribe = true);
+
+		// - construct from another pointField
+		//   subscribe to events if true
+		pointField( const pointField& src, bool subscribe);
+
+
+		// - copy construct 
+		pointField(const pointField& src);
+
+		// - no move construct
+		pointField(pointField&& src) = delete;
+
+
+		// assignment, only assign the VectorField and preserve other parts of this 
+		pointField& operator = (const pointField& rhs);
+
+		// no move assignment 
+		pointField& operator = (pointField&&) = delete;		
+
+
+		inline uniquePtr<pointFieldType> clone() const
+		{
+			return makeUnique<pointFieldType>(*this);
+		}
+
+		inline pointFieldType* clonePtr()const
+		{
+			return new pointFieldType(*this);
+		}
+
+	//// - Methods
+
+		// - reference to pointStructure
+		inline const pointStructure& pStruct()const {
+			return pStruct_;
+		}
+
+		// if all points are active 
+		INLINE_FUNCTION_H
+		bool allActive()const {	
+			return pStruct_.allActive();
+		}
+		
+
+		INLINE_FUNCTION_H
+		bool isActive(label i)const {
+			return pStruct_.isActive(i);
+		}
+
+		const auto& pointFlag()const
+		{
+			return pStruct_.pointFlag();
+		}
+
+		range activeRange()const
+		{
+			return pStruct_.activeRange();
+		}
+
+		// - update the field if any changes occure in pStruct 
+		//   for now it checks for deleted points
+		bool update(const eventMessage& msg);
+
+
+	//// -  IO operations
+		bool readPointField(iIstream& is);
+
+		bool writePointField(iOstream& os)const;
+
+
+		bool read(iIstream& is)
+		{
+			return readPointField(is);
+		}
+
+		bool write(iOstream& os)const
+		{
+			return writePointField(os);
+		}
+};
+
+template<template<class, class> class VectorField, class T, class MemorySpace>
+iIstream& operator >> (iIstream & is, pointField<VectorField, T, MemorySpace> & pF )
+{
+	if( !pF.read(is))
+	{
+		ioErrorInFile( is.name(), is.lineNumber() ) <<
+		"error in reading pointField from file. \n";
+		fatalExit;
+	}
+
+	return is;
+}
+
+template<template<class, class> class VectorField, class T, class MemorySpace>
+iOstream& operator << (iOstream& os, const pointField<VectorField, T, MemorySpace>& pF )
+{
+	if(! pF.write(os) )
+	{
+		ioErrorInFile( os.name(), os.lineNumber() )<<
+		"error in writing pointField into file. \n";
+		fatalExit;
+	}
+
+	return os;
+}
+
+}
+
+#include "pointField.cpp"
+#include "pointFieldAlgorithms.hpp"
+
+#endif //  __pointField_hpp__

src/phasicFlow/containers/pointField/pointFields.cpp

@@ -24,7 +24,7 @@ Licence:
 
 template class pFlow::pointField<pFlow::VectorSingle, pFlow::int8>; 
 
-template class pFlow::pointField<pFlow::VectorSingle, pFlow::int8, pFlow::HostSpace>; 
+/*template class pFlow::pointField<pFlow::VectorSingle, pFlow::int8, pFlow::HostSpace>; 
 
 template class pFlow::pointField<pFlow::VectorSingle, pFlow::int16>; 
 
@@ -52,7 +52,7 @@ template class pFlow::pointField<pFlow::VectorSingle, pFlow::real, pFlow::HostSp
 
 template class pFlow::pointField<pFlow::VectorSingle, pFlow::realx3>; 
 
-template class pFlow::pointField<pFlow::VectorSingle, pFlow::realx3, pFlow::HostSpace>;
+template class pFlow::pointField<pFlow::VectorSingle, pFlow::realx3, pFlow::HostSpace>;*/
 
 
 

src/phasicFlow/containers/pointField/pointFields.hpp

@@ -35,21 +35,17 @@ template<typename T>
 template<typename T>
   using pointField_D = pointField<VectorSingle, T>;
 
-template<typename T>
+/*template<typename T>
-  using pointField_HD = pointField<VectorDual, T>;
+  using pointField_HD = pointField<VectorDual, T>;*/
 
 
 using int8PointField_D 		= pointField<VectorSingle, int8>;
 
 using int8PointField_H 		= pointField<VectorSingle, int8, HostSpace>;
 
-using int16PointField_D 	= pointField<VectorSingle, int16>;
-
-using int16PointField_H 	= pointField<VectorSingle, int16, HostSpace>;
-
 using int32PointField_D 	= pointField<VectorSingle, int32>;
 
-using int32PointField_H 	= pointField<VectorSingle, int32, HostSpace>;
+/*using int32PointField_H 	= pointField<VectorSingle, int32, HostSpace>;
 
 using int64PointField_D 	= pointField<VectorSingle, int64>;
 
@@ -89,7 +85,7 @@ using realPointField_HD 	= pointField<VectorDual, real>;
 
 using realx3PointField_HD 	= pointField<VectorDual, realx3>;
 
-using wordPointField 		= pointField<Vector, word, vecAllocator<word>>;
+using wordPointField 		= pointField<Vector, word, vecAllocator<word>>;*/
 
 }
 

src/phasicFlow/demComponent/demComponent.cpp

@@ -0,0 +1,31 @@
+/*------------------------------- 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 "demComponent.hpp"
+
+
+
+pFlow::demComponent::demComponent(const word& name, systemControl& control)
+:
+	componentName_(name),
+	control_(control),
+	time_(control.time()),
+	timers_(name, &control.timers())
+{}

src/phasicFlow/demComponent/demComponent.hpp

@@ -22,7 +22,7 @@ Licence:
 #define __demComponent_hpp__
 
 #include "systemControl.hpp"
-
+#include "Time.hpp"
 
 
 namespace pFlow
@@ -48,6 +48,8 @@ protected:
 		/// Reference to systemControl 
 		systemControl& 	control_;
 
+		Time& 			time_;
+
 		/// All timers (if any) of this component
 		Timers 	timers_;
 
@@ -59,12 +61,7 @@ public:
 	// - Constructors
 	
 		/// construct from components
-		demComponent(const word& name, systemControl& control)
+		demComponent(const word& name, systemControl& control);
-		:
-			componentName_(name),
-			control_(control),
-			timers_(name, &control.timers())
-		{}
 		
 		/// No copy constructor
 		demComponent(const demComponent&) = delete;
@@ -102,16 +99,29 @@ public:
 		inline
 		real dt()const
 		{
-			return control_.time().dt();
+			return time_.dt();
 		}
 
 		/// Current simulation time 
 		inline
 		real currentTime()const
 		{
-			return control_.time().currentTime();
+			return time_.currentTime();
 		}
 
+		inline
+		const auto& time()const
+		{
+			return time_;
+		}
+
+		inline 
+		auto& time()
+		{
+			return time_;
+		}
+
+
 		/// Const ref to timers
 		inline
 		const auto& timers() const

src/phasicFlow/dictionary/dictionary.hpp

@@ -84,7 +84,7 @@ protected:
 		/// ref to parrent dictionary
 		const dictionary& parDict_;
 
-		bool isGlobal_;
+		bool isGlobal_ = false;
 
 
 	///// protected methods 
@@ -107,6 +107,12 @@ protected:
 		/// write dictionary to stream - with keyword
 		bool writeDictionary(iOstream& os, bool withBlock = true)const;
 
+        /// construct an empty dictionary with keyword and make it global/fileDictionary (if true)
+		dictionary(const word& keyword, bool global);
+
+		/// construct a dictionary with name and read it from file 
+		dictionary(const word& keyword, const fileSystem& file);
+	
 
 public:
 
@@ -123,12 +129,6 @@ public:
 		/// construct an empty dictionary with keyword
 		dictionary(const word& keyword);
 
-		/// construct an empty dictionary with keyword and make it global/fileDictionary (if true)
-		dictionary(const word& keyword, bool global);
-
-		/// construct a dictionary with name and read it from file 
-		dictionary(const word& keyword, const fileSystem& file);
-
 		/// cunstruct an empty dictionary with keyword and parDict
 		dictionary(const word& keyword, const dictionary& parDict);
 

src/phasicFlow/dictionary/fileDictionary.cpp

@@ -0,0 +1,26 @@
+#include "fileDictionary.hpp"
+#include "IOPattern.hpp"
+
+
+pFlow::fileDictionary::fileDictionary(const word &keyword)
+:
+    dictionary(keyword, true)
+{
+}
+
+pFlow::fileDictionary::fileDictionary(const word &keyword, const fileSystem &file)
+:
+    dictionary(keyword, file)
+{
+}
+
+bool pFlow::fileDictionary::read(iIstream &is, const IOPattern &iop)
+{
+    
+    return dictionary::read(is);
+}
+
+bool pFlow::fileDictionary::write(iOstream &os, const IOPattern &iop) const
+{
+    return dictionary::write(os);
+}

src/phasicFlow/dictionary/fileDictionary.hpp

@@ -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 "dictionary.hpp"
+
+
+namespace pFlow
+{
+
+class IOPattern;
+
+class fileDictionary
+:
+    public dictionary
+{
+public:
+
+    TypeInfo("fileDictionary");
+
+    /// construct an empty dictionary with keyword and make it global/fileDictionary
+	fileDictionary(const word& keyword);
+
+		/// construct a dictionary with name and read it from file 
+	fileDictionary(const word& keyword, const fileSystem& file);
+
+
+
+    /// read from stream
+	virtual bool read(iIstream& is, const IOPattern& iop);
+
+	/// write to stream
+	virtual bool write(iOstream& os, const IOPattern& iop) const;
+
+};
+
+}

src/phasicFlow/repository/IOobject/IOobject.cpp

@@ -123,7 +123,7 @@ bool pFlow::IOobject::read(iIstream& is, bool rdHdr)
 	{
 		if(!readHeader(is))return false;
 	}
-	return object_->read_object_t(is);
+	return object_->read_object_t(is, ioPattern_);
 }
 
 
@@ -132,5 +132,5 @@ bool pFlow::IOobject::write(iOstream& os) const
 	if(this->readWriteHeader())
 		writeHeader(os, typeName());
 	
-	return (object_->write_object_t(os) && writeSeparator(os));
+	return (object_->write_object_t(os, ioPattern_) && writeSeparator(os));
 }

src/phasicFlow/repository/IOobject/IOobject.hpp

@@ -50,9 +50,9 @@ private:
 
 		virtual word typeName()const = 0;
 
-		virtual bool read_object_t(iIstream& is) = 0;
+		virtual bool read_object_t(iIstream& is, IOPattern iop) = 0;
 
-		virtual bool write_object_t(iOstream& os)const = 0;
+		virtual bool write_object_t(iOstream& os, IOPattern iop)const = 0;
 
 	};
 
@@ -92,14 +92,14 @@ private:
 			return data_.typeName();
 		}
 
-		virtual bool read_object_t(iIstream& is)
+		bool read_object_t(iIstream& is, IOPattern iop) override
 		{
-			return data_.read(is);
+			return data_.read(is, iop);
 		}
 
-		virtual bool write_object_t(iOstream& os)const
+		bool write_object_t(iOstream& os, IOPattern iop)const override
 		{
-			return data_.write(os);
+			return data_.write(os, iop);
 		}
 
 		auto& data()

src/phasicFlow/repository/systemControl/systemControl.cpp

@@ -122,17 +122,17 @@ pFlow::systemControl::systemControl
 	),
 	settingsDict_
 	(
-		settings().emplaceObject<dictionary>
+		settings().emplaceObject<fileDictionary>
 		(
 			objectFile
 			(
 				settingsFile__,
 				"",
 				objectFile::READ_ALWAYS,
-				objectFile::WRITE_NEVER
+				objectFile::WRITE_NEVER,
+                IOPattern::AllProcessorsSimilar
 			),
-			settingsFile__,
+			settingsFile__
-			true
 		)
 	),
 	libs_(settingsDict_),
@@ -193,7 +193,7 @@ pFlow::systemControl::systemControl(
 	),
 	settingsDict_
 	(
-		settings().emplaceObject<dictionary>
+		settings().emplaceObject<fileDictionary>
 		(
 			objectFile
 			(
@@ -202,8 +202,7 @@ pFlow::systemControl::systemControl(
 				objectFile::READ_ALWAYS,
 				objectFile::WRITE_NEVER
 			),
-			settingsFile__,
+			settingsFile__
-			true
 		)
 	),
 	libs_(settingsDict_),
@@ -231,10 +230,10 @@ pFlow::systemControl::systemControl(
 	writeToFileTimer_("Write to file", &timers_)	
 {}
 
-pFlow::dictionary& pFlow::systemControl::domainDict()
+pFlow::fileDictionary& pFlow::systemControl::domainDict()
 {
 	return 
-	settings().emplaceObjectOrGet<dictionary>
+	settings().emplaceObjectOrGet<fileDictionary>
 	(
 		objectFile
 		(
@@ -243,8 +242,7 @@ pFlow::dictionary& pFlow::systemControl::domainDict()
 			objectFile::READ_ALWAYS,
 			objectFile::WRITE_NEVER
 		),
-		domainFile__,
+		domainFile__
-		true
 	);
 }
 

src/phasicFlow/repository/systemControl/systemControl.hpp

@@ -30,7 +30,7 @@ Licence:
 
 #include "types.hpp"
 #include "Time.hpp"
-#include "dictionary.hpp"
+#include "fileDictionary.hpp"
 #include "box.hpp"
 #include "Timers.hpp"
 #include "dynamicLinkLibs.hpp"
@@ -59,7 +59,7 @@ protected:
 	repository 			caseSetup_;
 
 	/// settingsDict fileDictionary
-	dictionary& 		settingsDict_;
+	fileDictionary& 		settingsDict_;
 
 	/// extra libs to be loaded
 	dynamicLinkLibs 	libs_;
@@ -155,15 +155,15 @@ public:
 		return timersReport_();
 	}
 
-	const dictionary& settingsDict()const{
+	const fileDictionary& settingsDict()const{
 		return settingsDict_;
 	}
 
-	dictionary& settingsDict(){
+	fileDictionary& settingsDict(){
 		return settingsDict_;
 	}
 
-	dictionary& domainDict();
+	fileDictionary& domainDict();
 	
 
 	virtual word runName() const

src/phasicFlow/setFieldList/setFieldEntryTemplates.cpp

@@ -50,8 +50,8 @@ void* pFlow::setFieldEntry::setPointFieldDefaultValueNew
 	Type defValue = entry_.secondPartVal<Type>();
 	
 	if(verbose)
-		REPORT(2)<<"Creating pointField " << greenText(fieldName())<< " with default value " << cyanText(defValue)<<
+		REPORT(2)<<"Creating pointField " << Green_Text(fieldName())<< " with default value " << cyanText(defValue)<<
-	             " in repository "<< owner.name() <<endREPORT;
+	             " in repository "<< owner.name() <<END_REPORT;
 
 	
 	auto& field = 
@@ -85,8 +85,8 @@ void* pFlow::setFieldEntry::setPointFieldDefaultValueStdNew
 	Type defValue = entry_.secondPartVal<Type>();
 	
 	if(verbose)
-		REPORT(2)<<"Creating pointField " << greenText(fieldName())<< " with default value " << cyanText(defValue)<<
+		REPORT(2)<<"Creating pointField " << Green_Text(fieldName())<< " with default value " << cyanText(defValue)<<
-	             " in repository "<< owner.name() <<endREPORT;
+	             " in repository "<< owner.name() <<END_REPORT;
 
 	// by default we perform operations on host
 	auto& field = 
@@ -129,8 +129,8 @@ void* pFlow::setFieldEntry::setPointFieldSelected
 	Type value = entry_.secondPartVal<Type>();
 
 	if(verbose)
-		REPORT(2)<< "Setting selected points of " << greenText(fName)
+		REPORT(2)<< "Setting selected points of " << Green_Text(fName)
-		 		 << " to value " << cyanText(value) <<endREPORT;
+		 		 << " to value " << cyanText(value) <<END_REPORT;
 	
 	
 	auto fieldTypeName = owner.lookupObjectTypeName(fName);
@@ -155,7 +155,7 @@ void* pFlow::setFieldEntry::setPointFieldSelected
 			return nullptr;
 	}
 
-	if( pointField<VectorDual,Type>::TYPENAME() == fieldTypeName )
+	/*if( pointField<VectorDual,Type>::TYPENAME() == fieldTypeName )
 	{
 		
 		
@@ -164,7 +164,7 @@ void* pFlow::setFieldEntry::setPointFieldSelected
 			return &field;
 		else
 			return nullptr;
-	}
+	}*/
 
 	fatalErrorInFunction<<
 	fieldTypeName<< " is not a supported field type for setFieldEntry.\n";
@@ -197,8 +197,8 @@ void* pFlow::setFieldEntry::setPointFieldSelectedStd
 	Type value = entry_.secondPartVal<Type>();
 
 	if(verbose)
-		REPORT(2)<< "Setting selected points of " << greenText(fName)
+		REPORT(2)<< "Setting selected points of " << Green_Text(fName)
-		 		 << " to value " << cyanText(value) <<endREPORT;
+		 		 << " to value " << cyanText(value) <<END_REPORT;
 	
 	
 	auto fieldTypeName = owner.lookupObjectTypeName(fName);

src/phasicFlow/smartPointers/uniquePtr.hpp

@@ -105,7 +105,4 @@ inline uniquePtr<T> makeUnique(Args&&... args)
 
 }
 
-
-
-
 #endif

src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp

@@ -20,6 +20,7 @@ Licence:
 
 #include "boundaryBase.hpp"
 #include "dictionary.hpp"
+#include "internalPoints.hpp"
 
 /*pFlow::boundaryBase::boundaryBase
 (
@@ -59,13 +60,30 @@ pFlow::boundaryBase::boundaryBase
 
 }
 
+typename pFlow::boundaryBase::pointFieldAccessType 
+    pFlow::boundaryBase::thisPoints()
+{
     
-pFlow::uniquePtr<pFlow::boundaryBase> pFlow::boundaryBase::create
+    return pointFieldAccessType
-(
+    (
-	const dictionary& 	dict,
+        indexList_.size(),
-	const plane& 		bplane,
+        indexList_.deviceVectorAll(),
-	internalPoints& 	internal
+        internal_.pointPosition().deviceVectorAll()
-)
+    );
+   
+}
+
+typename pFlow::boundaryBase::pointFieldAccessType
+    pFlow::boundaryBase::neighborPoints()
+{
+    notImplementedFunction;
+    return pointFieldAccessType();
+}
+
+pFlow::uniquePtr<pFlow::boundaryBase> pFlow::boundaryBase::create(
+    const dictionary &dict,
+    const plane &bplane,
+    internalPoints &internal)
 {
 	word type = dict.getVal<word>("type");
 	word bType = angleBracketsNames("boundary", type);

src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp

@@ -25,6 +25,7 @@ Licence:
 #include "subscriber.hpp"
 #include "VectorSingles.hpp"
 #include "plane.hpp"
+#include "scatterFieldAccess.hpp"
 
 #include "streams.hpp"
 
@@ -42,6 +43,9 @@ class boundaryBase
 {
 public:
 
+    using pointFieldAccessType = 
+        scatterFieldAccess<realx3,DefaultExecutionSpace>;
+
 	enum DIRECTION: int8
 	{
 		
@@ -106,11 +110,20 @@ public:
 		(dict, bplane, internal)
 	);
 
-	virtual bool beforeIteratoin(uint32 iterNum, real t) = 0 ;
+	virtual 
+    bool beforeIteratoin(uint32 iterNum, real t) = 0 ;
 
-	virtual bool iterate(uint32 iterNum, real t) = 0;
+	virtual 
+    bool iterate(uint32 iterNum, real t) = 0;
+
+	virtual 
+    bool afterIteration(uint32 iterNum, real t) = 0;
+
+    pointFieldAccessType thisPoints();
+
+    virtual
+    pointFieldAccessType neighborPoints();
     
-	virtual bool afterIteration(uint32 iterNum, real t) = 0;
 
 	static
 	uniquePtr<boundaryBase> create

src/phasicFlow/structuredData/boundaries/boundaryBase/scatterFieldAccess.hpp

@@ -9,14 +9,18 @@ namespace pFlow
 {
 
 
-template<typename T, typename ExecutionSpace>
+template<typename T, typename MemorySpace>
 class scatterFieldAccess
 {
 public:
 	
-	using execution_space 	= ExecutionSpace;
+	using viewType 			= ViewType1D<T, MemorySpace>;
 
-	using memory_space 		= typename execution_space::memory_space;
+    using device_type 		= typename viewType::device_type;
+
+    using memory_space 		= typename viewType::memory_space;
+
+    using execution_space 	= typename viewType::execution_space;
 
 protected:
 
@@ -24,12 +28,17 @@ protected:
 
 	ViewType1D<uint32, memory_space> 	indices_;
 
-	ViewType1D<T, memory_space> 		fieldVals_;
+	viewType 		fieldVals_;
 
 public:
 
+    scatterFieldAccess():
+        indices_(),
+        fieldVals_()    
+    {}
+
 	scatterFieldAccess(
-		uint32 	sz
+		uint32 	sz,
 		ViewType1D<uint32, memory_space> ind,
 		ViewType1D<T, memory_space> fVals)
 	:
@@ -50,31 +59,37 @@ public:
 
 	// - Methods
 
+    INLINE_FUNCTION_HD
 	T& operator()(uint32 i)
 	{
 		return fieldVals_(indices_(i));
 	}
 
+    INLINE_FUNCTION_HD
 	const T& operator()(uint32 i)const
 	{
 		return fieldVals_(indices_(i));
 	}
 
+    INLINE_FUNCTION_HD
 	T& operator[](uint32 i)
 	{
 		return fieldVals_(indices_(i));
 	}
 
+    INLINE_FUNCTION_HD
 	const T& operator[](uint32 i)const
 	{
 		return fieldVals_(indices_(i));
 	}
 
+    INLINE_FUNCTION_HD
 	uint32 	size()const
 	{
 		return size_;
 	}
 
+    INLINE_FUNCTION_HD
 	bool empty()const
 	{
 		return size_ == 0;

src/phasicFlow/structuredData/pointStructure/boundaryList.cpp

@@ -27,9 +27,9 @@ pFlow::boundaryList::boundaryList
 	internalPoints& internal
 )
 :
+    ListPtr<boundaryBase>(simD.sizeOfBoundaries()),
 	simDomain_(simD),
-	internal_(internal),
+	internal_(internal)
-	boundaryList_(simD.sizeOfBoundaries())
 {}
 
 bool pFlow::boundaryList::updateLists()
@@ -37,7 +37,7 @@ bool pFlow::boundaryList::updateLists()
 
 	for(auto i=0; i<simDomain_.sizeOfBoundaries();i++)
 	{
-		boundaryList_.set
+		this->set
 		(
 			i,
 			boundaryBase::create

src/phasicFlow/structuredData/pointStructure/boundaryList.hpp

@@ -33,6 +33,8 @@ namespace pFlow
 class simulationDomain;
 
 class boundaryList
+:
+    public ListPtr<boundaryBase>
 {
 
 protected:
@@ -42,8 +44,6 @@ protected:
 		
         const simulationDomain& simDomain_;
 
-        ListPtr<boundaryBase> 	boundaryList_;
-
         bool                    listSet_ = false;
 
     
@@ -65,6 +65,16 @@ public:
   ~boundaryList() = default;
 
   bool updateLists();
+
+  auto& boundary(size_t i)
+  {
+    return ListPtr<boundaryBase>::operator[](i);
+  }
+
+  const auto& boundary(size_t i)const
+  {
+    return ListPtr<boundaryBase>::operator[](i);
+  }
 };
 
 } // pFlow

src/phasicFlow/structuredData/pointStructure/internalPoints.cpp

@@ -187,6 +187,12 @@ const pFlow::internalPoints::pFlagTypeHost&
 	return pFlagsH_;
 }
 
+FUNCTION_H
+pFlow::realx3Field_D &pFlow::internalPoints::pointPosition()
+{
+    return pointPosition_;
+}
+
 FUNCTION_H
 void pFlow::internalPoints::updateFlag
 (

src/phasicFlow/structuredData/pointStructure/internalPoints.hpp

@@ -117,7 +117,8 @@ public:
 		FUNCTION_H
 		const realx3Field_D& pointPosition()const;
 		
-		
+		FUNCTION_H
+        realx3Field_D& pointPosition();
 
 		/*INLINE_FUNCTION_H
 		auto pointPositionHostAll()

src/phasicFlow/structuredData/pointStructure/pointStructure.cpp

@@ -23,23 +23,6 @@ Licence:
 
 #include "streams.hpp"
 
-pFlow::pointStructure::pointStructure
-(
-  const dictionary& domainDict
-)
-:
-  internalPoints(),
-  simulationDomain_
-  (
-    simulationDomain::create(domainDict)
-  ),
-  boundaries_
-  (
-    simulationDomain_(),
-    *this
-  )
-{}
-
 bool pFlow::pointStructure::distributePoints
 (
 	const realx3Vector &points
@@ -85,12 +68,33 @@ bool pFlow::pointStructure::distributePoints
 
 }
 
+pFlow::pointStructure::pointStructure
+(
+    systemControl& control
+)
+:
+	demComponent("particlesCM", control),
+  	internalPoints(),
+	simulationDomain_
+	(
+		simulationDomain::create(control.domainDict())
+	),
+	boundaries_
+	(
+		simulationDomain_(),
+		*this
+	)
+{}
+
+
 pFlow::pointStructure::pointStructure(
-    const dictionary &domainDict,
+    systemControl& control,
     const realx3Vector &posVec)
-    : internalPoints(),
+:
+	demComponent("particlesCM", control), 
+	internalPoints(),
 	simulationDomain_(
-          simulationDomain::create(domainDict)),
+		simulationDomain::create(control.domainDict())),
 	boundaries_(
 		simulationDomain_(),
 		*this)
@@ -113,12 +117,24 @@ pFlow::pointStructure::pointStructure(
     
 }
 
-FUNCTION_H
+bool pFlow::pointStructure::beforeIteration()
-bool pFlow::pointStructure::read
+{
-(
+    return true;
-    iIstream& is, 
+}
-    IOPattern::IOType iotype
+
-)
+bool pFlow::pointStructure::iterate()
+{
+    return true;
+}
+
+bool pFlow::pointStructure::afterIteration()
+{
+    return true;
+}
+
+bool pFlow::pointStructure::read(
+    iIstream &is,
+    IOPattern::IOType iotype)
 {
     Field<Vector, realx3 , vecAllocator<realx3>> fRead("file_internalPoints", "internalPoints");
     

src/phasicFlow/structuredData/pointStructure/pointStructure.hpp

@@ -17,25 +17,23 @@ Licence:
   implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
 
 -----------------------------------------------------------------------------*/
-
-
 #ifndef __pointStructure_hpp__
 #define __pointStructure_hpp__
 
+#include "demComponent.hpp"
 #include "internalPoints.hpp"
 #include "simulationDomain.hpp"
 #include "boundaryList.hpp"
-#include "uniquePtr.hpp"
+
 
 
 namespace pFlow
 {
 
-class simulationDomain;
-
 
 class pointStructure
 :
+    public demComponent,
 	public internalPoints
 {
 protected:
@@ -57,15 +55,12 @@ public:
 
 	//// - Constructors 
 
-		// - an empty pointStructure, good for reading from file 
+		/// an empty pointStructure, good for reading from file 
-		pointStructure(const dictionary& domainDict);
+		pointStructure(systemControl& control);
-
-		// - construct from components
-		//pointStructure(const int8Vector& flgVec, const realx3Vector& posVec);
 
 		/// construct from point positions, assume all points are active
 		pointStructure(
-			const dictionary& domainDict, 
+            systemControl& control,
 			const realx3Vector& posVec);
 
 		// - copy construct
@@ -89,8 +84,46 @@ public:
 		// - destructor 
 		virtual ~pointStructure() = default;
 
+	// - Iteration methods
 
-	/// Read 
+		/// In the time loop before iterate
+		bool beforeIteration() override;
+
+		/// @brief This is called in time loop. Perform the main calculations 
+		/// when the component should evolve along time.
+		bool iterate() override;
+
+		/// This is called in time loop, after iterate.
+		bool afterIteration() override;
+
+    // - Access methods 
+
+        auto& boundaries()
+        {
+            return boundaries_;
+        }
+
+        const auto& boundaries()const
+        {
+            return boundaries_;
+        }
+
+        auto& boundary(size_t i )
+        {
+            return boundaries_[i];
+        }
+
+        auto& boundary(size_t i)const
+        {
+            return boundaries_[i];
+        }
+
+	// - IO methods 
+
+		/// @brief read the point structure from the input stream.
+		/// @param is: input stream
+		/// @param iotype: type of input pattern
+		/// @return true on success
 		FUNCTION_H
 		bool read(iIstream& is, IOPattern::IOType iotype);
 		

utilities/particlesPhasicFlow/CMakeLists.txt

@@ -3,7 +3,7 @@ set(source_files
 particlesPhasicFlow.cpp
 positionParticles/positionParticles.cpp 
 positionOrdered/positionOrdered.cpp 
-positionRandom/positionRandom.cpp
+#positionRandom/positionRandom.cpp
 empty/empty.cpp
 )
 set(link_lib phasicFlow Kokkos::kokkos Interaction Utilities)

utilities/particlesPhasicFlow/empty/empty.cpp

@@ -30,7 +30,7 @@ pFlow::empty::empty(
 	positionParticles(control, dict),
 	position_
 	(
-		maxNumberOfParticles_, 0, RESERVE()
+		"empty",maxNumberOfParticles_, 0, RESERVE()
 	)
 {
 }

utilities/particlesPhasicFlow/particlesPhasicFlow.cpp

@@ -24,6 +24,8 @@ Licence:
 #include "setFields.hpp"
 #include "systemControl.hpp"
 #include "commandLine.hpp"
+#include "vocabs.hpp"
+
 //#include "readControlDict.hpp"
 
 using pFlow::output;
@@ -74,7 +76,7 @@ int main( int argc, char* argv[] )
 	if(setOnly && positionOnly)
 	{
 		ERR<<
-		"Options --positionParticles-only and --setFields-only cannot be used simeltanuously. \n"<<endERR;
+		"Options --positionParticles-only and --setFields-only cannot be used simeltanuously. \n"<<END_ERR;
 		return 1;
 	}	
 
@@ -103,8 +105,8 @@ int main( int argc, char* argv[] )
 	{
 
 		// position particles based on the dict content 
-		REPORT(0)<< "Positioning points . . . \n"<<endREPORT;
+		REPORT(0)<< "Positioning points . . . \n"<<END_REPORT;
-		auto pointPosition = positionParticles::create(cpDict.subDict("positionParticles"));
+		auto pointPosition = positionParticles::create(Control, cpDict.subDict("positionParticles"));
 
 		fileSystem pStructPath = Control.time().path()+pointStructureFile__;
 
@@ -120,6 +122,7 @@ int main( int argc, char* argv[] )
 				objectFile::READ_NEVER,
 				objectFile::WRITE_ALWAYS			
 			),
+			Control,
             finalPos
 		);
 
@@ -127,9 +130,9 @@ int main( int argc, char* argv[] )
 
 
 		REPORT(1)<< "Created pStruct with "<< pStruct.size() << " points and capacity "<<
-		pStruct.capacity()<<" . . ."<< endREPORT;
+		pStruct.capacity()<<" . . ."<< END_REPORT;
 
-		REPORT(1)<< "Writing pStruct to " << Control.time().path()+ pointStructureFile__<< endREPORT<<endl<<endl;
+		REPORT(1)<< "Writing pStruct to " << Control.time().path()+ pointStructureFile__<< END_REPORT<<endl<<endl;
 
 		if( !Control.time().write())
 		{
@@ -148,7 +151,8 @@ int main( int argc, char* argv[] )
 				Control.time().path(),
 				objectFile::READ_NEVER,
 				objectFile::WRITE_ALWAYS			
-			)
+			),
+            Control
 		);
 
 		pStructPtr = &pStruct;
@@ -170,7 +174,7 @@ int main( int argc, char* argv[] )
 		{
 			if( !sfEntry.setPointFieldDefaultValueNewAll(Control.time(), pStruct, true))
 			{
-				ERR<< "\n error occured in setting default value fields.\n"<<endERR;
+				ERR<< "\n error occured in setting default value fields.\n"<<END_ERR;
 				return 1;
 			}
 		}
@@ -183,13 +187,13 @@ int main( int argc, char* argv[] )
 
 		for(auto name: selNames)
 		{
-			REPORT(1)<< "Applying selector " << greenText(name) <<endREPORT;
+			REPORT(1)<< "Applying selector " << Green_Text(name) <<END_REPORT;
 			
 			if(
 				!applySelector(Control, pStruct, selectorsDict.subDict(name))
 			 )
 			{
-				ERR<<"\n error occured in setting selector. \n"<<endERR;
+				ERR<<"\n error occured in setting selector. \n"<<END_ERR;
 				return 1;
 			}
 			output<<endl;
@@ -197,7 +201,7 @@ int main( int argc, char* argv[] )
 	}
 
 	Control.time().write(true);
-	REPORT(0)<< greenText("\nFinished successfully.\n")<<endREPORT;
+	REPORT(0)<< Green_Text("\nFinished successfully.\n")<<END_REPORT;
 
 
 // this should be palced in each main 

utilities/particlesPhasicFlow/positionOrdered/positionOrdered.cpp

@@ -127,7 +127,7 @@ pFlow::positionOrdered::positionOrdered
 	const dictionary& dict
 )
 :
-	positionParticles(dict),
+	positionParticles(control, dict),
 	poDict_
 	(  
 		dict.subDict("positionOrderedInfo")
@@ -138,7 +138,7 @@ pFlow::positionOrdered::positionOrdered
 	),
 	numPoints_
 	(
-		poDict_.getVal<size_t>("numPoints")
+		poDict_.getVal<uint64>("numPoints")
 	),
 	axisOrder_
 	(
@@ -146,7 +146,7 @@ pFlow::positionOrdered::positionOrdered
 	),
 	position_
 	(
-		maxNumberOfParticles_, RESERVE()
+		"positionOrdered", maxNumberOfParticles_, numPoints_ ,RESERVE()
 	)
 {
 	

utilities/particlesPhasicFlow/positionOrdered/positionOrdered.hpp

@@ -76,12 +76,12 @@ public:
 
 	//// - Methods 
 
-	virtual label numPoints()const
+	virtual uint64 numPoints()const
 	{
 		return position_.size();
 	}
 
-	virtual label size()const
+	virtual uint64 size()const
 	{
 		return position_.size();
 	}

utilities/particlesPhasicFlow/positionParticles/positionParticles.cpp

@@ -54,7 +54,7 @@ pFlow::realx3Vector pFlow::positionParticles::sortByMortonCode(realx3Vector& pos
 			ind++});
 	}
 
-	INFORMATION<<"Performing morton sorting."<<endINFO;
+	INFORMATION<<"Performing morton sorting."<<END_INFO;
 	std::sort(
 		indMor.begin(),
 		indMor.end(),
@@ -120,7 +120,11 @@ pFlow::realx3Vector pFlow::positionParticles::getFinalPosition()
 }
 
 pFlow::uniquePtr<pFlow::positionParticles> 
-	pFlow::positionParticles::create(const dictionary & dict)
+	pFlow::positionParticles::create
+(
+    systemControl& control,
+    const dictionary & dict
+)
 {
 
 	word method = dict.getVal<word>("method");

utilities/particlesPhasicFlow/positionParticles/positionParticles.hpp

@@ -127,8 +127,8 @@ public:
 	(
 		positionParticles,
 		dictionary,
-		(const dictionary& dict),
+		(systemControl& control, const dictionary& dict),
-		(dict)
+		(control, dict)
 	);
 
 	virtual ~positionParticles() = default;
@@ -151,7 +151,7 @@ public:
 	virtual realx3Vector getFinalPosition();
 
 	static
-	uniquePtr<positionParticles> create(const dictionary & dict);
+	uniquePtr<positionParticles> create(systemControl& control, const dictionary & dict);
 	
 };