Prep for proper com in circle group

4 files changed, 51 insertions, 3 deletions

diff --git a/examples/poiseulle_particles_2d_hlbm_gpu/sim.cpp b/examples/poiseulle_particles_2d_hlbm_gpu/sim.cpp
index e12b0d8..9375078 100644
--- a/examples/poiseulle_particles_2d_hlbm_gpu/sim.cpp
+++ b/examples/poiseulle_particles_2d_hlbm_gpu/sim.cpp
@@ -176,7 +176,7 @@ saw::error_or<void> step(
 
 	// auto coll_ev = 
 	q.submit([&](acpp::sycl::handler& h){
-		component<T,Desc,cmpt::HLBM,encode::Sycl<saw::encode::Native>> collision{0.65};
+		component<T,Desc,cmpt::Hlbm,encode::Sycl<saw::encode::Native>> collision{0.65};
 		component<T,Desc,cmpt::BounceBack,encode::Sycl<saw::encode::Native>> bb;
 		component<T,Desc,cmpt::AntiBounceBack<0u>,encode::Sycl<saw::encode::Native>> abb;
 
diff --git a/examples/poiseulle_particles_2d_ibm_gpu/sim.cpp b/examples/poiseulle_particles_2d_ibm_gpu/sim.cpp
index ea49d51..e68d7da 100644
--- a/examples/poiseulle_particles_2d_ibm_gpu/sim.cpp
+++ b/examples/poiseulle_particles_2d_ibm_gpu/sim.cpp
@@ -1,6 +1,7 @@
 #include <kel/lbm/sycl/lbm.hpp>
 #include <kel/lbm/lbm.hpp>
 #include <kel/lbm/particle.hpp>
+#include <kel/lbm/math/n_linear.hpp>
 
 #include <forstio/io/io.hpp>
 #include <forstio/remote/filesystem/easy.hpp>
@@ -266,7 +267,39 @@ saw::error_or<void> step(
 	*/
 	q.submit([&](acpp::sycl::handler& h){
 		h.parallel_for(acpp::sycl::range<1u>{1u}, [=](acpp::sycl::id<1u> idx){
-				
+			auto& vel = macros.template get<"velocity">();
+
+			auto& ps = particles;
+			auto& mask = ps.template get<"mask">();
+			auto& dense = ps.template get<"density">().at({});
+			auto& com = ps.template get<"center_of_mass">();
+
+			auto& parts = ps.template get<"particles">();
+
+			auto& p_i = parts.at({{0u}});
+
+			auto& p_i_rb = p_i.template get<"rigid_body">();
+			/// 0. Iterate over mask and calculate position in LBM grid
+			/// In this case it's simple since I'm too lazy to do scaling and rotation
+			/// Technically scale => rotate => translate
+			/// Here it's only translate
+			auto& p_i_rb_pos = p_i_rb.template get<"position">();
+
+			iterator<Desc::D>::apply([&](const auto& index){
+				/// Calculate the shift from the mask
+				saw::data<sch::Vector<T,Desc::D>> index_shift;
+				for(uint64_t i = 0u; i < Desc::D; ++i){
+					index_shift.at({{i}}) = index.at({i}).template cast_to<T>() - com.at({{i}});
+				}
+
+
+				/// TODO 1. Calculate force pickup from neigbouring u_vel cells
+				auto inter_vel = n_linear_interpolate(vel,index_shift);
+
+				/// TODO 3. Distribute force to fluid
+
+
+			}, {}, mask.meta());
 		});
 	}).wait();
 
diff --git a/lib/core/c++/math/n_linear.hpp b/lib/core/c++/math/n_linear.hpp
index 8fb0600..b378440 100644
--- a/lib/core/c++/math/n_linear.hpp
+++ b/lib/core/c++/math/n_linear.hpp
@@ -122,6 +122,9 @@ auto n_linear_interpolate(
 		}
 
 	}, {}, ones_ind);
+
+	/// TODO I need to actually calc stuff
+	return field.at({});
 }
 
 template<typename FieldSchema, typename T>
diff --git a/lib/core/c++/particle/particle.hpp b/lib/core/c++/particle/particle.hpp
index 01429b2..5110893 100644
--- a/lib/core/c++/particle/particle.hpp
+++ b/lib/core/c++/particle/particle.hpp
@@ -55,6 +55,8 @@ template<typename T, uint64_t D, typename CollisionType = ParticleCollisionSpher
 using ParticleGroup = Struct<
 	Member<Array<T,D>, "mask">,
 	Member<FixedArray<Scalar<T>,1u>, "density">,
+	Member<Vector<T,D>, "center_of_mass">,
+	Member<Scalar<T>, "total_mass">,
 	Member<Array<Particle<T,D>>, "particles">
 >;
 }
@@ -76,10 +78,20 @@ saw::data<sch::ParticleGroup<T,D, sch::ParticleCollisionSpheroid<T,radius>>> cre
 	for(uint64_t i = 0u; i < D; ++i){
 		mask_dims.at({i}) = mask_resolution;
 	}
+	saw::data<sch::Scalar<T>> mask_step;
+	saw::data<T> dia_d{radius*2};
+	mask_step.at({}) = dia_d / mask_resolution.template cast_to<T>();
 
 	mask = {mask_dims};
+
+	auto& com = part.template get<"center_of_mass">();
+	for(uint64_t i = 0u; i < D; ++i){
+		com.at({{i}}) = {};
+	}
+	saw::data<sch::UInt64> ele_ctr{0u};
 	iterator<D>::apply([&](const auto& index){
-		
+		++ele_ctr;
+
 	},{},mask_dims);
 
 	return part;