Working on making cavity code available to everyone else :)

4 files changed, 80 insertions, 23 deletions

diff --git a/c++/SConscript b/c++/SConscript
index 57091af..3b79d1a 100644
--- a/c++/SConscript
+++ b/c++/SConscript
@@ -13,7 +13,7 @@ dir_path = Dir('.').abspath
 core_env = env.Clone();
 
 core_env.sources = sorted(glob.glob(dir_path + "/*.cpp"));
-core_env.headers = sorted(glob.glob(dir_path + "/*.h"));
+core_env.headers = sorted(glob.glob(dir_path + "/*.hpp"));
 
 env.sources += core_env.sources;
 env.headers += core_env.headers;
diff --git a/c++/boundary.hpp b/c++/boundary.hpp
new file mode 100644
index 0000000..77542c5
--- /dev/null
+++ b/c++/boundary.hpp
@@ -0,0 +1,59 @@
+#pragma once
+
+#include "component.hpp"
+#include "equilibrium.hpp"
+
+namespace kel {
+namespace lbm {
+namespace cmpt {
+struct BounceBack {};
+}
+
+template<typename T, typename Descriptor>
+class component<T, Descriptor, cmpt::BounceBack> {
+private:
+	typename saw::native_data_type<T>::type relaxation_;
+public:
+	component(typename saw::native_data_type<T>::type relaxation__):
+		relaxation_{relaxation__}
+	{}
+
+	using Component = cmpt::BounceBack;
+
+	/**
+	 * Thoughts regarding automating order setup
+	 */
+	using access = cmpt::access_tuple<
+		cmpt::access<"dfs", 1, true, true, true>
+	>;
+
+	/**
+	 * Thoughts regarding automating order setup
+	 */
+	static constexpr saw::string_literal name = "collision";
+	static constexpr saw::string_literal after = "";
+	static constexpr saw::string_literal before = "streaming";
+
+	/**
+	 * Raw setup
+	 */
+	template<typename CellFieldSchema>
+	void apply(saw::data<CellFieldSchema>& field, saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>> index, uint64_t time_step){
+		bool is_even = ((time_step % 2) == 0);
+		auto& cell = field.at(index);
+
+		auto& dfs_old = is_even ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
+		auto& dfs = (!is_even) ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
+
+		using dfi = df_info<sch::T,Desc>;
+
+		// Technically use .copy()
+		auto df_cpy = dfs;
+
+		for(uint64_t i = 1u; i < Desc::Q; ++i){
+			dfs({i}) = df_cpy({dfi::opposite_index.at(i)});
+		}
+	}
+};
+}
+}
diff --git a/c++/collision.hpp b/c++/collision.hpp
index 9143862..7ca7c9f 100644
--- a/c++/collision.hpp
+++ b/c++/collision.hpp
@@ -11,7 +11,13 @@ struct BGK {};
 
 template<typename T, typename Descriptor>
 class component<T, Descriptor, cmpt::BGK> {
+private:
+	typename saw::native_data_type<T>::type relaxation_;
 public:
+	component(typename saw::native_data_type<T>::type relaxation__):
+		relaxation_{relaxation__}
+	{}
+
 	using Component = cmpt::BGK;
 
 	/**
@@ -31,14 +37,22 @@ public:
 	/**
 	 * Raw setup
 	 */
-	void apply(saw::data<sch::CellField>& field, saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>> index, uint64_t time_step){
+	template<typename CellFieldSchema>
+	void apply(saw::data<CellFieldSchema>& field, saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>> index, uint64_t time_step){
 		bool is_even = ((time_step % 2) == 0);
 		auto& cell = field.at(index);
 
 		auto& dfs_old = is_even ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
 		auto& dfs = (!is_even) ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
 
+		typename saw::native_data_type<T>::type rho;
+		std::array<typename saw::native_data_type<T>::type, Descriptor::D> vel;
+		compute_rho_u<T,Descriptor>(dfs,rho,vel);
+		auto eq = equilibrium<Descriptor>(rho,vel);
 
+		for(uint64_t i = 0u; i < Descriptor::Q; ++i){
+			dfs({i}).set(dfs({i}).get() + (1.0 / relaxation_) * (eq[i] - dfs({i}).get()));
+		}
 	}
 
 
diff --git a/examples/cavity_2d.cpp b/examples/cavity_2d.cpp
index 654a9ac..5683d85 100644
--- a/examples/cavity_2d.cpp
+++ b/examples/cavity_2d.cpp
@@ -2,6 +2,7 @@
 #include "../c++/macroscopic.hpp"
 #include "../c++/lbm.hpp"
 #include "../c++/component.hpp"
+#include "../c++/collision.hpp"
 
 /**
  */
@@ -106,8 +107,6 @@ public:
 namespace cmpt {
 struct BounceBack{};
 struct MovingWall {};
-struct BGK {};
-struct ConstRhoBGK {};
 }
 
 
@@ -155,23 +154,6 @@ public:
 		*/
 	}
 };
-
-template<typename Desc>
-class component<sch::T, Desc, cmpt::BGK> {
-public:
-	typename saw::native_data_type<sch::T>::type relaxation_;
-public:
-	void apply(saw::data<sch::DfCell<Desc>>& dfs){
-		typename saw::native_data_type<sch::T>::type rho;
-		std::array<typename saw::native_data_type<sch::T>::type, Desc::D> vel;
-		compute_rho_u<sch::T,Desc>(dfs,rho,vel);
-		auto eq = equilibrium<Desc>(rho,vel);
-		
-		for(uint64_t i = 0u; i < Desc::Q; ++i){
-			dfs({i}).set(dfs({i}).get() + (1.0 / relaxation_) * (eq[i] - dfs({i}).get()));
-		}
-	}
-};
 }
 }
 
@@ -260,8 +242,10 @@ void lbm_step(
 	using namespace kel::lbm;
   using dfi = df_info<sch::T,sch::D2Q9>;
 
-	component<sch::T, sch::D2Q9, cmpt::BGK> coll;
-	coll.relaxation_ = 0.5384;
+  /**
+	* Relaxation parameter \tau
+	*/
+	component<sch::T, sch::D2Q9, cmpt::BGK> coll{0.5384};
 
 	component<sch::T, sch::D2Q9, cmpt::BounceBack> bb;
 	component<sch::T, sch::D2Q9, cmpt::MovingWall> bb_lid;