Set equilibrium boundary for sim

3 files changed, 26 insertions, 10 deletions

diff --git a/examples/poiseulle_particles_2d_gpu/sim.cpp b/examples/poiseulle_particles_2d_gpu/sim.cpp
index d53a00c..6fdff7b 100644
--- a/examples/poiseulle_particles_2d_gpu/sim.cpp
+++ b/examples/poiseulle_particles_2d_gpu/sim.cpp
@@ -119,6 +119,14 @@ saw::error_or<void> step(
 		// Need nicer things to handle the flow. I see improvement here
 		component<T,Desc,cmpt::BGK,encode::Sycl<saw::encode::Native>> collision{0.6};
 		component<T,Desc,cmpt::BounceBack,encode::Sycl<saw::encode::Native>> bb;
+
+		saw::data<sch::Scalar<T>> rho_b;
+		rho_b.at({}) = 1.0;
+		saw::data<sch::Vector<T,Desc::D>> vel_b;
+		vel_b.at({{0u}}) = 0.01;
+
+		component<T,Desc,cmpt::Equilibrium,encode::Sycl<saw::encode::Native>> equi{rho_b,vel_b};
+
 		component<T,Desc,cmpt::ZouHeHorizontal<true>,encode::Sycl<saw::encode::Native>> flow_in{
 			[&](){
 				uint64_t target_t_i = 256u;
@@ -148,11 +156,13 @@ saw::error_or<void> step(
 					collision.apply(fields,macros,index,t_i);
 				break;
 				case 3u:
-					flow_in.apply(fields,index,t_i);
+					equi.apply(fields,index,t_i);
+					// flow_in.apply(fields,index,t_i);
 					collision.apply(fields,macros,index,t_i);
 				break;
 				case 4u:
-					flow_out.apply(fields,index,t_i);
+					equi.apply(fields,index,t_i);
+					// flow_out.apply(fields,index,t_i);
 					collision.apply(fields,macros,index,t_i);
 				break;
 				default:
diff --git a/lib/core/c++/boundary.hpp b/lib/core/c++/boundary.hpp
index b70ffbf..d5f3022 100644
--- a/lib/core/c++/boundary.hpp
+++ b/lib/core/c++/boundary.hpp
@@ -69,22 +69,28 @@ public:
 template<typename FP, typename Descriptor, typename Encode>
 class component<FP, Descriptor, cmpt::Equilibrium, Encode> final {
 private:
-	saw::data<FP> density_;
-	saw::data<FP> velocity_;
+	saw::data<sch::Scalar<FP>> density_;
+	saw::data<sch::Vector<FP,Descriptor::D>> velocity_;
 public:
 	component(
-		saw::data<FP> density__,
-		saw::Data<FP> velocity__
+		saw::data<sch::Scalar<FP>> density__,
+		saw::data<sch::Vector<FP,Descriptor::D>> velocity__
 	):
 		density_{density__},
 		velocity_{velocity__}
 	{}
 
 	template<typename CellFieldSchema>
-	void apply(const saw::data<CellFieldSchema, Encode>& field, const saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>>& index, saw::data<sch::UInt64> time_step)const{
-		using dfi = df_info<T,Descriptor>;
+	void apply(const saw::data<CellFieldSchema, Encode>& field, const saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>>& index, saw::data<sch::UInt64> time_step) const {
+		
+		bool is_even = ((time_step.get() % 2u) == 0u);
+		using dfi = df_info<FP,Descriptor>;
 		
 		auto& dfs_old_f = (is_even) ? field.template get<"dfs_old">() : field.template get<"dfs">();
+
+		auto eq = equilibrium<FP,Descriptor>(density_,velocity_);
+
+		dfs_old_f.at(index) = eq;
 	}
 };
 
diff --git a/lib/core/c++/collision.hpp b/lib/core/c++/collision.hpp
index a05e263..6c39646 100644
--- a/lib/core/c++/collision.hpp
+++ b/lib/core/c++/collision.hpp
@@ -52,8 +52,8 @@ public:
 		// auto& dfs_f = (is_even) ? field.template get<"dfs">() : field.template get<"dfs_old">();
 		auto& dfs_old_f = (is_even) ? field.template get<"dfs_old">() : field.template get<"dfs">();
 
-		saw::data<T> rho;
-		saw::data<sch::FixedArray<T,Descriptor::D>> vel;
+		saw::data<sch::Scalar<T>> rho;
+		saw::data<sch::Vector<T,Descriptor::D>> vel;
 		compute_rho_u<T,Descriptor>(dfs_old_f.at(index),rho,vel);
 		auto eq = equilibrium<T,Descriptor>(rho,vel);