Reflective boundary condition is added and tested.

It requires messaging integration when changing velocity.

src/Integration/AdamsBashforth2/AdamsBashforth2.hpp

@@ -41,11 +41,11 @@ class AdamsBashforth2
 {
 private:
 
-
 	auto& dy1()
 	{
 		return static_cast<realx3PointField_D&>(*this);
 	}
+	
 public:
 
 	/// Type info

src/phasicFlow/CMakeLists.txt

@@ -78,6 +78,7 @@ structuredData/boundaries/boundaryBase/boundaryBaseKernels.cpp
 structuredData/boundaries/boundaryExit/boundaryExit.cpp
 structuredData/boundaries/boundaryNone/boundaryNone.cpp
 structuredData/boundaries/boundaryPeriodic/boundaryPeriodic.cpp
+structuredData/boundaries/boundaryReflective/boundaryReflective.cpp
 structuredData/boundaries/boundaryList.cpp
 structuredData/pointStructure/pointStructure/pointStructure.cpp
 structuredData/pointStructure/selectors/pStructSelector/pStructSelector.cpp

src/phasicFlow/containers/pointField/boundaryField/boundaryField/boundaryField.cpp

@@ -68,7 +68,7 @@ pFlow::uniquePtr<pFlow::boundaryField<T, MemorySpace>>
 		printKeys
 		( 
 			fatalError << "Ctor Selector "<< bType << "for type "<< 
-			Yellow_Text(getTypeName<T>())	<< " dose not exist.\n"
+			Yellow_Text(getTypeName<T>())	<< " does not exist.\n"
 			<<"Avaiable ones are: \n\n"
 			,
 			boundaryBasevCtorSelector_

src/phasicFlow/containers/pointField/boundaryField/createBoundaryFields.hpp

@@ -21,15 +21,23 @@ Licence:
 #ifndef __createBoundaryFields_hpp__
 #define __createBoundaryFields_hpp__
 
-#include "boundaryField.hpp"
+#include "boundaryField/boundaryField.hpp"
-#include "exitBoundaryField.hpp"
+#include "exitBoundaryField/exitBoundaryField.hpp"
+#include "periodicBoundaryField/periodicBoundaryField.hpp"
+#include "reflectiveBoundaryField/reflectiveBoundaryField.hpp"
+
+#define createDerivedBoundary(DataType, MemorySpaceType)          \
+	template class pFlow::exitBoundaryField<DataType, MemorySpaceType>; \
+	template class pFlow::periodicBoundaryField<DataType, MemorySpaceType>; \
+	template class pFlow::reflectiveBoundaryField<DataType, MemorySpaceType>;
 
 #define createBaseBoundary(DataType, MemorySpaceType)          \
-    template class pFlow::boundaryField<DataType, MemorySpaceType>;
+	template class pFlow::boundaryField<DataType, MemorySpaceType>;
     
 
-#define createBoundary(DataType, MemorySpaceType, BoundaryType)    \
+#define createBoundaryFields(DataType, MemorySpaceType)    \
-    template class pFlow::BoundaryType##BoundaryField<DataType, MemorySpaceType>;
+    createBaseBoundary(DataType, MemorySpaceType); \
+	createDerivedBoundary(DataType, MemorySpaceType);
 
 
 #endif //__createBoundaryFields_hpp__

src/phasicFlow/containers/pointField/boundaryField/reflectiveBoundaryField/reflectiveBoundaryField.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.
+-----------------------------------------------------------------------------*/
+
+template<class T, class MemorySpace>
+    pFlow::reflectiveBoundaryField<T, MemorySpace>::reflectiveBoundaryField
+(
+	const boundaryBase& boundary, 
+	const pointStructure& pStruct,
+	InternalFieldType& internal
+)
+:
+	BoundaryFieldType(boundary, pStruct, internal)
+{
+	//this->addEvent(message::BNDR_DELETE);
+}

src/phasicFlow/containers/pointField/boundaryField/reflectiveBoundaryField/reflectiveBoundaryField.hpp

@@ -0,0 +1,90 @@
+/*------------------------------- 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 __reflectiveBoundaryField_hpp__
+#define __reflectiveBoundaryField_hpp__
+
+#include "boundaryField.hpp"
+
+namespace pFlow
+{
+
+template< class T, class MemorySpace = void>
+class reflectiveBoundaryField
+:
+    public boundaryField<T, MemorySpace> 
+{
+public:
+	
+    using ReflectiveBoundaryFieldType = reflectiveBoundaryField<T, MemorySpace>;
+
+	using BoundaryFieldType = boundaryField<T, MemorySpace>;
+
+	using InternalFieldType = typename BoundaryFieldType::InternalFieldType;
+
+	using memory_space 		= typename BoundaryFieldType::memory_space;
+
+	using execution_space 	= typename BoundaryFieldType::execution_space;
+
+   
+
+public:
+
+	TypeInfoTemplate211("boundaryField","reflective", T, memory_space::name());
+
+	reflectiveBoundaryField(
+		const boundaryBase& boundary, 
+		const pointStructure& pStruct,
+		InternalFieldType& internal);
+		
+
+	add_vCtor
+	(
+		BoundaryFieldType,
+		ReflectiveBoundaryFieldType,
+		boundaryBase
+	);
+
+
+	bool hearChanges
+	(
+		real t,
+		real dt,
+		uint32 iter,
+		const message& msg, 
+    	const anyList& varList
+	) override
+    {
+		BoundaryFieldType::hearChanges(t,dt,iter, msg,varList);
+		
+		if(msg.equivalentTo(message::BNDR_DELETE))
+		{
+			// do nothing;
+		}
+		
+		return true;
+	}
+
+};
+
+}
+
+#include "reflectiveBoundaryField.cpp"
+
+#endif

src/phasicFlow/containers/pointField/pointField/pointFields.cpp

@@ -18,76 +18,37 @@ Licence:
 
 -----------------------------------------------------------------------------*/
 
-
 #include "pointFields.hpp"
 #include "createBoundaryFields.hpp"
-#include "periodicBoundaryField.hpp"
-
-#define createAllBoundary(DataType, MemorySpaceType)          \
-    template class pFlow::exitBoundaryField<DataType, MemorySpaceType>; \
-    template class pFlow::periodicBoundaryField<DataType, MemorySpaceType>;
 
 
+#define createPointFields(DataType) 							\
+  template class pFlow::pointField<DataType, pFlow::HostSpace>; \
+  createBoundaryFields(DataType, pFlow::HostSpace); 			\
+              													\
+  template class pFlow::pointField<DataType>;   				\
+  createBoundaryFields(DataType, void);
 
 // uint8
-template class pFlow::pointField<pFlow::uint8, pFlow::HostSpace>; 
+createPointFields(pFlow::uint8);
-createBaseBoundary(pFlow::uint8, pFlow::HostSpace);
-createAllBoundary(pFlow::uint8, pFlow::HostSpace);
-
-template class pFlow::pointField<pFlow::uint8>; 
-createBaseBoundary(pFlow::uint8, void);
-createAllBoundary(pFlow::uint8, void);
 
 /// uint32
-template class pFlow::pointField<pFlow::uint32, pFlow::HostSpace>; 
+createPointFields(pFlow::uint32);
-createBaseBoundary(pFlow::uint32, pFlow::HostSpace);
-createAllBoundary(pFlow::uint32, pFlow::HostSpace);
-
-template class pFlow::pointField<pFlow::uint32>; 
-createBaseBoundary(pFlow::uint32, void);
-createAllBoundary(pFlow::uint32, void);
 
 /// uint64
-template class pFlow::pointField<pFlow::uint64, pFlow::HostSpace>; 
+createPointFields(pFlow::uint64);
-createBaseBoundary(pFlow::uint64, pFlow::HostSpace);
-createAllBoundary(pFlow::uint64, pFlow::HostSpace);
-
-template class pFlow::pointField<pFlow::uint64>; 
-createBaseBoundary(pFlow::uint64, void);
-createAllBoundary(pFlow::uint64, void);
 
 /// real 
-template class pFlow::pointField<pFlow::real, pFlow::HostSpace>; 
+createPointFields(pFlow::real);
-createBaseBoundary(pFlow::real, pFlow::HostSpace);
-createAllBoundary(pFlow::real, pFlow::HostSpace);
-
-template class pFlow::pointField<pFlow::real>; 
-createBaseBoundary(pFlow::real, void);
-createAllBoundary(pFlow::real, void);
 
 /// realx3
-template class pFlow::pointField<pFlow::realx3, pFlow::HostSpace>; 
+createPointFields(pFlow::realx3);
-createBaseBoundary(pFlow::realx3, pFlow::HostSpace);
-createAllBoundary(pFlow::realx3, pFlow::HostSpace);
-
-
-template class pFlow::pointField<pFlow::realx3>; 
-createBaseBoundary(pFlow::realx3, void);
-createAllBoundary(pFlow::realx3, void);
 
 
 /// realx4
-template class pFlow::pointField<pFlow::realx4, pFlow::HostSpace>; 
+createPointFields(pFlow::realx4);
-createBaseBoundary(pFlow::realx4, pFlow::HostSpace);
-createAllBoundary(pFlow::realx4, pFlow::HostSpace);
 
 
-template class pFlow::pointField<pFlow::realx4>; 
+/// word, only on host
-createBaseBoundary(pFlow::realx4, void);
-createAllBoundary(pFlow::realx4, void);
-
-
-/// word 
 template class pFlow::pointField<pFlow::word, pFlow::HostSpace>; 
-createBaseBoundary(pFlow::word, pFlow::HostSpace);
+createBoundaryFields(pFlow::word, pFlow::HostSpace);
-createAllBoundary(pFlow::word, pFlow::HostSpace);

src/phasicFlow/repository/repository/repository.hpp

@@ -33,7 +33,7 @@ namespace pFlow
 
 class repository
 {
-protected:
+private:
 	
 	// - repository name
 	word 		name_;
@@ -53,10 +53,11 @@ protected:
 
 
 	template <typename Type1>
-		word reportTypeError (IOobject& object);
+	word reportTypeError (IOobject& object)const;
 
 	template <typename Type>
-		bool checkForObjectType(IOobject& object);
+	static
+	bool checkForObjectType(IOobject& object);
 
 public:
 
@@ -167,10 +168,14 @@ public:
 			}
 		}
 		
-		// - return a ref to the underlaying data in the object 
+		/// return a ref to the underlaying data in the object 
 		template<typename T>
 		T& lookupObject(const word& name);
 
+		/// return a const ref to the underlaying data in the object 
+		template<typename T>
+		const T& lookupObject(const word& name)const;
+		
 		// - search the name and return a ref to repository
 		repository& lookupRepository(const word& name);
 

src/phasicFlow/repository/repository/repositoryTemplates.cpp

@@ -20,7 +20,7 @@ Licence:
 
 
 template <typename Type1>
-pFlow::word pFlow::repository::reportTypeError(IOobject& object)
+pFlow::word pFlow::repository::reportTypeError(IOobject& object)const
 {
 	word err;
 	err = "Object " + object.name() + " with type " + Type1::TYPENAME() + 
@@ -40,28 +40,57 @@ template<typename T>
 T& pFlow::repository::lookupObject(const word& name)
 {
 	if( auto [iter, success] = objects_.findIf(name); success )
-    {
+	{
 
-    	if( checkType<T>(iter->second) )
+		if( checkType<T>(iter->second) )
-    	{
+		{
-    		return static_cast<T&>(*iter->second);
+			return static_cast<T&>(*iter->second);
 			
 
-    	}else
+		}else
-    	{
+		{
-    		fatalErrorInFunction << 
+			fatalErrorInFunction << 
-    		reportTypeError<T>(*iter->second)<<endl;
+			reportTypeError<T>(*iter->second)<<endl;
-    		fatalExit;
+			fatalExit;
-    		return static_cast<T&>(*iter->second);
+			return static_cast<T&>(*iter->second);
-    	}
+		}
 		
-    }
+	}
-    else
+	else
-    {
+	{
-    	fatalErrorInFunction << 
+		fatalErrorInFunction << 
 		"Object with name " << name << " is not found in repository " << this->name()<<endl <<
 		"list of avaiable objest is \n" << objectNames();
-        fatalExit;
+		fatalExit;
-        return static_cast<T&>(*iter->second);
+		return static_cast<T&>(*iter->second);
-    }
+	}
+}
+
+template<typename T>
+const T& pFlow::repository::lookupObject(const word& name)const
+{
+	if( auto [iter, success] = objects_.findIf(name); success )
+	{
+
+		if( checkType<T>(iter->second) )
+		{
+			return static_cast<const T&>(*iter->second);
+		}
+		else
+		{
+			fatalErrorInFunction << 
+			reportTypeError<T>(*iter->second)<<endl;
+			fatalExit;
+			return static_cast<T&>(*iter->second);
+		}
+		
+	}
+	else
+	{
+		fatalErrorInFunction << 
+		"Object with name " << name << " is not found in repository " << this->name()<<endl <<
+		"list of avaiable objest is \n" << objectNames();
+		fatalExit;
+		return static_cast<T&>(*iter->second);
+	}
 }

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

@@ -24,6 +24,7 @@ Licence:
 #include "internalPoints.hpp"
 #include "anyList.hpp"
 #include "Time.hpp"
+#include "pointStructure.hpp"
 #include "boundaryBaseKernels.hpp"
 
 void pFlow::boundaryBase::setSize(uint32 newSize)
@@ -114,9 +115,9 @@ bool pFlow::boundaryBase::removeIndices
 	
 	message msgBndry 	= message::BNDR_RESET;
 
-	uint32 iter = internal_.time().currentIter();
+	uint32 iter = time().currentIter();
-	real t = internal_.time().currentTime();
+	real t = time().currentTime();
-	real dt = internal_.time().dt();
+	real dt = time().dt();
 
 	if( !this->notify(iter, t, dt, msgBndry, aList) )
 	{
@@ -197,6 +198,16 @@ pFlow::boundaryBase::boundaryBase
 	unSyncLists();
 }
 
+const pFlow::pointStructure &pFlow::boundaryBase::pStruct() const
+{
+    return boundaries_.pStruct();
+}
+
+const pFlow::Time &pFlow::boundaryBase::time() const
+{
+    return boundaries_.pStruct().time();
+}
+
 pFlow::boundaryBase &pFlow::boundaryBase::mirrorBoundary()
 {
     return boundaries_[mirrorBoundaryIndex()];
@@ -256,7 +267,7 @@ pFlow::uniquePtr<pFlow::boundaryBase> pFlow::boundaryBase::create
 	{
 		printKeys
 		( 
-			fatalError << "Ctor Selector "<< bType << " dose not exist. \n"
+			fatalError << "Ctor Selector "<< bType << " does not exist. \n"
 			<<"Avaiable ones are: \n\n"
 			,
 			dictionaryvCtorSelector_

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

@@ -37,6 +37,8 @@ namespace pFlow
 class internalPoints;
 class dictionary;
 class boundaryList;
+class pointStructure;
+class Time;
 
 class boundaryBase
 :
@@ -217,6 +219,10 @@ public:
 		return internal_;
 	}
 
+	const pointStructure& pStruct()const;
+
+	const Time& time()const;
+
 	inline
 	const auto& indexList()const
 	{

src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.cpp

@@ -0,0 +1,149 @@
+/*------------------------------- 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 "boundaryReflective.hpp"
+#include "pointFields.hpp"
+#include "dictionary.hpp"
+#include "Time.hpp"
+
+pFlow::boundaryReflective::boundaryReflective
+(
+    const dictionary &dict, 
+    const plane &bplane, 
+    internalPoints &internal, 
+    boundaryList &bndrs, 
+    uint32 thisIndex
+)
+:
+    boundaryBase
+    (
+        dict,
+        bplane,
+        internal,
+        bndrs,
+        thisIndex
+    )
+{
+    restitution_ = dict.getValOrSet("restitution", restitution_);
+    velocityName_ = dict.getValOrSet("velocityName",velocityName_);
+    diameterName_ = dict.getValOrSet("diameterName", diameterName_);
+}
+
+bool pFlow::boundaryReflective::beforeIteration(
+    uint32 iterNum,
+    real t,
+    real dt)
+{
+    return true;
+}
+
+bool pFlow::boundaryReflective::iterate
+(
+    uint32 iterNum, 
+    real t, 
+    real dt
+)
+{
+    return true;
+}
+
+bool pFlow::boundaryReflective::afterIteration
+(
+    uint32 iterNum, 
+    real t, 
+    real dt
+)
+{
+    if(empty())return true;
+    
+    uint32 s = size();
+	uint32Vector_D inContactFlags("inContactFlags",s+1, s+1, RESERVE()); 
+	inContactFlags.fill(0u);
+	auto inContactFlagsD = inContactFlags.deviceViewAll();
+
+	auto points = thisPoints();
+    
+	auto p = boundaryPlane().infPlane();
+    const auto &diam = time().lookupObject<realPointField_D>(diameterName_);
+	auto diams = diam.BoundaryField(thisBoundaryIndex()).thisField();
+
+	uint32 numInContact = 0;
+
+	Kokkos::parallel_reduce
+	(
+		"pFlow::boundaryReflective::afterIteration",
+		deviceRPolicyStatic(0u,s),
+		LAMBDA_HD(uint32 i, uint32& nContactToUpdate)
+		{
+			if(p.inPositiveDistance(points(i), 0.5*diams(i)))
+			{
+				inContactFlagsD(i)=1;
+				nContactToUpdate++;
+			}
+		}, 
+		numInContact
+	);
+
+    // no particle in contact 
+	if(numInContact == 0 )
+	{
+		return true;
+	}
+
+    uint32Vector_D inContactList("inContactList", numInContact);
+    const auto& inContactListD = inContactList.deviceViewAll();
+
+    exclusiveScan(inContactFlagsD, 0u, s+1, inContactFlagsD, 0u);
+    
+    Kokkos::parallel_for
+    (
+        "pFlow::boundaryReflective::afterIteration",
+        deviceRPolicyStatic(0, s),
+        LAMBDA_HD(uint32 i)
+		{
+			if(inContactFlagsD(i)!= inContactFlagsD(i+1))
+				inContactListD(inContactFlagsD(i)) = points.index(i);
+		}
+    );
+    Kokkos::fence(); 
+
+    const auto& velocity = time().lookupObject<realx3PointField_D>(velocityName_);
+    const auto& velocityD = velocity.deviceViewAll();
+
+    Kokkos::parallel_for(
+        "pFlow::boundaryReflective::velocityChange",
+        deviceRPolicyStatic(0,numInContact),
+        LAMBDA_HD(uint32 i)
+        {
+            auto& vel = velocityD(inContactListD(i));
+            real vn = dot(p.normal(), vel);
+            
+            if(vn < 0)
+            {
+                realx3 vt = vel - vn*p.normal();
+                vel = restitution_*(vt - vn*p.normal());
+            }
+        }
+    );
+    Kokkos::fence();
+    
+    // TODO: notify integration for changes in the velocity 
+
+    return true;
+}

src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.hpp

@@ -0,0 +1,73 @@
+/*------------------------------- 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 __boundaryReflective_hpp__
+#define __boundaryReflective_hpp__
+
+
+#include "boundaryBase.hpp"
+
+namespace pFlow
+{
+
+class boundaryReflective
+:
+ 	public boundaryBase
+{
+private:
+
+	real        restitution_ = 0.95;
+
+    word        velocityName_{"velocity"};
+	
+    word        diameterName_{"diameter"};
+
+public:
+
+	TypeInfo("boundary<reflective>");
+
+	boundaryReflective(
+		const dictionary &dict,
+		const plane 	&bplane,
+		internalPoints 	&internal,
+		boundaryList	&bndrs,
+		uint32 			thisIndex);
+
+	 
+	~boundaryReflective()override = default;
+	
+	add_vCtor
+	(
+		boundaryBase,
+		boundaryReflective,
+		dictionary
+	);
+
+	bool beforeIteration(uint32 iterNum, real t, real dt) override;
+
+	bool iterate(uint32 iterNum, real t, real dt) override;
+
+	bool afterIteration(uint32 iterNum, real t, real dt) override;
+
+
+};
+
+}
+
+#endif