Dangling changes. Restructuring particle setup

6 files changed, 51 insertions, 90 deletions

diff --git a/examples/poiseulle_particles_2d_bgk_gpu/.nix/derivation.nix b/examples/poiseulle_particles_2d_bgk_gpu/.nix/derivation.nix
index 517b6b8..8aca15e 100644
--- a/examples/poiseulle_particles_2d_bgk_gpu/.nix/derivation.nix
+++ b/examples/poiseulle_particles_2d_bgk_gpu/.nix/derivation.nix
@@ -11,7 +11,7 @@
 }:
 
 stdenv.mkDerivation {
-	pname = pname + "-examples-" + "poiseulle_2d_gpu";
+	pname = pname + "-examples-" + "poiseulle_particles_2d_bgk_gpu";
 	inherit version;
 	src = ./..;
 
diff --git a/examples/poiseulle_particles_2d_bgk_gpu/sim.cpp b/examples/poiseulle_particles_2d_bgk_gpu/sim.cpp
index 6bf4f84..7d1a624 100644
--- a/examples/poiseulle_particles_2d_bgk_gpu/sim.cpp
+++ b/examples/poiseulle_particles_2d_bgk_gpu/sim.cpp
@@ -57,8 +57,7 @@ using MacroStruct = Struct<
 template<typename T, typename Desc>
 saw::error_or<void> setup_initial_conditions(
 		saw::data<sch::ChunkStruct<T,Desc>>& fields,
-		saw::data<sch::MacroStruct<T,Desc>>& macros,
-		saw::data<sch::FixedArray<sch::Particle<T,Desc::D>, particle_amount>>& particles
+		saw::data<sch::MacroStruct<T,Desc>>& macros
 ){
 	auto& info_f = fields.template get<"info">();
 	// Set everything as walls
@@ -162,21 +161,12 @@ template<typename T, typename Desc>
 saw::error_or<void> step(
 		saw::data<sch::Ptr<sch::ChunkStruct<T,Desc>>,encode::Sycl<saw::encode::Native>>& fields,
 		saw::data<sch::Ptr<sch::MacroStruct<T,Desc>>,encode::Sycl<saw::encode::Native>>& macros,
-		saw::data<sch::FixedArray<sch::Particle<T,Desc::D>, particle_amount>>& particles,
 		saw::data<sch::UInt64> t_i,
 		device& dev
 ){
 	auto& q = dev.get_handle();
 	auto& info_f = fields.template get<"info">();
 
-	{
-		auto& p = particles.at({{0u}});
-
-		auto& p_coll = p.template get<"collision">();
-		auto& p_rad = p_coll.template get<"radius">();
-	}
-
-
 	// auto coll_ev = 
 	q.submit([&](acpp::sycl::handler& h){
 		// Need nicer things to handle the flow. I see improvement here
@@ -283,7 +273,6 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 	// saw::data<sch::FixedArray<sch::UInt64,Desc::D>> meta{{dim_x,dim_y}};
 	auto lbm_data_ptr = saw::heap<saw::data<sch::ChunkStruct<T,Desc>>>();
 	auto lbm_macro_data_ptr = saw::heap<saw::data<sch::MacroStruct<T,Desc>>>();
-	auto lbm_particle_data_ptr = saw::heap<saw::data<sch::FixedArray<sch::Particle<T,Desc::D>, particle_amount>>>();
 	
 	std::cout<<"Estimated Bytes: "<<memory_estimate<sch::ChunkStruct<T,Desc>,sch::MacroStruct<T,Desc>>().get()<<std::endl;
 
@@ -309,7 +298,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 
 	sycl_q.wait();
 	{
-		auto eov = setup_initial_conditions<T,Desc>(*lbm_data_ptr,*lbm_macro_data_ptr,*lbm_particle_data_ptr);
+		auto eov = setup_initial_conditions<T,Desc>(*lbm_data_ptr,*lbm_macro_data_ptr);
 		if(eov.is_error()){
 			return eov;
 		}
@@ -323,7 +312,6 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 
 	saw::data<sch::ChunkStruct<T,Desc>, encode::Sycl<saw::encode::Native>> lbm_sycl_data{sycl_q};
 	saw::data<sch::MacroStruct<T,Desc>, encode::Sycl<saw::encode::Native>> lbm_sycl_macro_data{sycl_q};
-	saw::data<sch::FixedArray<sch::Particle<T,Desc::D>, particle_amount>, encode::Sycl<saw::encode::Native>> lbm_sycl_particle_data{sycl_q};
 	sycl_q.wait();
 
 	auto lsd_view = make_chunk_struct_view(lbm_sycl_data);
@@ -340,12 +328,6 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 			return eov;
 		}
 	}
-	{
-		auto eov = dev.copy_to_device(*lbm_particle_data_ptr,lbm_sycl_particle_data);
-		if(eov.is_error()){
-			return eov;
-		}
-	}
 	sycl_q.wait();
 	saw::data<sch::UInt64> time_steps{16u*4096ul};
 
@@ -354,7 +336,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 	for(saw::data<sch::UInt64> i{0u}; i < time_steps and krun; ++i){
 		// BC + Collision
 		{
-			auto eov = step<T,Desc>(lsd_view,lsdm_view,*lbm_particle_data_ptr,i,dev);
+			auto eov = step<T,Desc>(lsd_view,lsdm_view,i,dev);
 			if(eov.is_error()){
 				return eov;
 			}
diff --git a/examples/poiseulle_particles_2d_psm_gpu/sim.cpp b/examples/poiseulle_particles_2d_psm_gpu/sim.cpp
index 7702808..7e76da6 100644
--- a/examples/poiseulle_particles_2d_psm_gpu/sim.cpp
+++ b/examples/poiseulle_particles_2d_psm_gpu/sim.cpp
@@ -60,8 +60,7 @@ using MacroStruct = Struct<
 template<typename T, typename Desc>
 saw::error_or<void> setup_initial_conditions(
 		saw::data<sch::ChunkStruct<T,Desc>>& fields,
-		saw::data<sch::MacroStruct<T,Desc>>& macros,
-		saw::data<sch::FixedArray<sch::Particle<T,Desc::D>, particle_amount>>& particles
+		saw::data<sch::MacroStruct<T,Desc>>& macros
 ){
 	auto& info_f = fields.template get<"info">();
 	auto& porous_f = macros.template get<"porosity">();
@@ -178,23 +177,6 @@ saw::error_or<void> setup_initial_conditions(
 		df_f.get_dims()
 	);
 
-	for(saw::data<sch::UInt64> i{0u}; i < saw::data<sch::UInt64>{particle_amount}; ++i){
-		auto& part = particles.at(i);
-
-		saw::data<sch::Vector<T,Desc::D>> pos;
-		pos.at({{0u}}) = dim_x * 0.25;
-		pos.at({{1u}}) = dim_y * 0.5;
-		saw::data<sch::Scalar<T>> rad, dense, dt;
-		rad.at({}) = dim_y * 0.1;
-		dense.at({}) = 1.0;
-		dt.at({}) = 1.0;
-		part = create_spheroid_particle(
-			pos,{},{},
-			{},{},{},
-			rad, dense,dt
-		);
-	}
-
 	return saw::make_void();
 }
 
@@ -202,7 +184,6 @@ template<typename T, typename Desc>
 saw::error_or<void> step(
 		saw::data<sch::Ptr<sch::ChunkStruct<T,Desc>>,encode::Sycl<saw::encode::Native>>& fields,
 		saw::data<sch::Ptr<sch::MacroStruct<T,Desc>>,encode::Sycl<saw::encode::Native>>& macros,
-		saw::data<sch::FixedArray<sch::Particle<T,Desc::D>, particle_amount>>& particles,
 		saw::data<sch::UInt64> t_i,
 		device& dev
 ){
@@ -210,14 +191,6 @@ saw::error_or<void> step(
 	auto& info_f = fields.template get<"info">();
 	auto& porous_f = macros.template get<"porosity">();
 
-	{
-		auto& p = particles.at({{0u}});
-
-		auto& p_coll = p.template get<"collision">();
-		auto& p_rad = p_coll.template get<"radius">();
-	}
-
-
 	q.submit([&](acpp::sycl::handler& h){
 		component<T,Desc,cmpt::PSM,encode::Sycl<saw::encode::Native>> collision{0.8};
 		component<T,Desc,cmpt::BounceBack,encode::Sycl<saw::encode::Native>> bb;
@@ -326,7 +299,6 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 	// saw::data<sch::FixedArray<sch::UInt64,Desc::D>> meta{{dim_x,dim_y}};
 	auto lbm_data_ptr = saw::heap<saw::data<sch::ChunkStruct<T,Desc>>>();
 	auto lbm_macro_data_ptr = saw::heap<saw::data<sch::MacroStruct<T,Desc>>>();
-	auto lbm_particle_data_ptr = saw::heap<saw::data<sch::FixedArray<sch::Particle<T,Desc::D>, particle_amount>>>();
 	
 	std::cout<<"Estimated Bytes: "<<memory_estimate<sch::ChunkStruct<T,Desc>,sch::MacroStruct<T,Desc>>().get()<<std::endl;
 
@@ -352,7 +324,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 
 	sycl_q.wait();
 	{
-		auto eov = setup_initial_conditions<T,Desc>(*lbm_data_ptr,*lbm_macro_data_ptr,*lbm_particle_data_ptr);
+		auto eov = setup_initial_conditions<T,Desc>(*lbm_data_ptr,*lbm_macro_data_ptr);
 		if(eov.is_error()){
 			return eov;
 		}
@@ -366,7 +338,6 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 
 	saw::data<sch::ChunkStruct<T,Desc>, encode::Sycl<saw::encode::Native>> lbm_sycl_data{sycl_q};
 	saw::data<sch::MacroStruct<T,Desc>, encode::Sycl<saw::encode::Native>> lbm_sycl_macro_data{sycl_q};
-	saw::data<sch::FixedArray<sch::Particle<T,Desc::D>, particle_amount>, encode::Sycl<saw::encode::Native>> lbm_sycl_particle_data{sycl_q};
 	sycl_q.wait();
 
 	auto lsd_view = make_chunk_struct_view(lbm_sycl_data);
@@ -383,12 +354,6 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 			return eov;
 		}
 	}
-	{
-		auto eov = dev.copy_to_device(*lbm_particle_data_ptr,lbm_sycl_particle_data);
-		if(eov.is_error()){
-			return eov;
-		}
-	}
 	sycl_q.wait();
 	saw::data<sch::UInt64> time_steps{16u*4096ul};
 
@@ -397,13 +362,13 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 	for(saw::data<sch::UInt64> i{0u}; i < time_steps and krun; ++i){
 		// BC + Collision
 		{
-			auto eov = step<T,Desc>(lsd_view,lsdm_view,*lbm_particle_data_ptr,i,dev);
+			auto eov = step<T,Desc>(lsd_view,lsdm_view,i,dev);
 			if(eov.is_error()){
 				return eov;
 			}
 		}
 		sycl_q.wait();
-		if(i.get() %  1 == 0u){
+		if(i.get() %  32u == 0u){
 			{
 				auto eov = dev.copy_to_host(lbm_sycl_macro_data,*lbm_macro_data_ptr);
 				if(eov.is_error()){
diff --git a/examples/settling_cubes_2d_ibm_gpu/sim.cpp b/examples/settling_cubes_2d_ibm_gpu/sim.cpp
index bb6ab59..aec8bb9 100644
--- a/examples/settling_cubes_2d_ibm_gpu/sim.cpp
+++ b/examples/settling_cubes_2d_ibm_gpu/sim.cpp
@@ -49,15 +49,10 @@ using MacroStruct = Struct<
 	Member<VectorChunk<T,Desc>, "force">
 >;
 
-//template<typename T, typename Desc>
-//using ParticleArray = Array<
-//	Particle<T,Desc::D>
-//>;
-
 template<typename T, typename Desc>
 using ParticleGroups = Tuple<
 	ParticleGroup<
-		T,Desc::D
+		T,Desc::D,sch::ParticleCollisionSpheroid<T>
 	>
 >;
 }
@@ -65,7 +60,8 @@ using ParticleGroups = Tuple<
 template<typename T, typename Desc>
 saw::error_or<void> setup_initial_conditions(
 		saw::data<sch::ChunkStruct<T,Desc>>& fields,
-		saw::data<sch::MacroStruct<T,Desc>>& macros
+		saw::data<sch::MacroStruct<T,Desc>>& macros,
+		saw::data<sch::ParticleGroups<T,Desc>>& particles
 ){
 	auto& info_f = fields.template get<"info">();
 	// Set everything as walls
@@ -105,12 +101,20 @@ saw::error_or<void> setup_initial_conditions(
 		df_f.get_dims()
 	);
 
-	return saw::make_void();
-}
+	// Particles
+	{
+		saw::data<sch::Scalar<T>> rad_p, dense_p;
+		rad_p.at({}).set(2.0);
+		dense_p.at({}).set(1);
+		auto& spheroid_group = particles.template get<0u>();
+		spheroid_group = create_spheroid_particle_group<T,Desc::D>(
+			rad_p,
+			dense_p,
+			{64u}
+		);
+	}
 
-template<typename T, typename Desc>
-void add_particles(saw::data<sch::Array<sch::Particle<T,Desc::D>>>& particles){
-	////// TODO
+	return saw::make_void();
 }
 
 template<typename T, typename Desc>
@@ -129,6 +133,8 @@ saw::error_or<void> step(
 
 	// auto coll_ev = 
 	q.submit([&](acpp::sycl::handler& h){
+		saw::data<sch::Vector<T,Desc::D>> force;
+		force.at({{1}}).set(-1.0);
 		// Need nicer things to handle the flow. I see improvement here
 		component<T,Desc,cmpt::BGKGuo, encode::Sycl<saw::encode::Native>> collision{0.8};
 		component<T,Desc,cmpt::BounceBack,encode::Sycl<saw::encode::Native>> bb;
@@ -205,7 +211,8 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 	// saw::data<sch::FixedArray<sch::UInt64,Desc::D>> meta{{dim_x,dim_y}};
 	auto lbm_data_ptr = saw::heap<saw::data<sch::ChunkStruct<T,Desc>>>();
 	auto lbm_macro_data_ptr = saw::heap<saw::data<sch::MacroStruct<T,Desc>>>();
-	
+	auto lbm_particles_data = saw::data<sch::ParticleGroups<T,Desc>>();
+
 	std::cout<<"Estimated Bytes: "<<memory_estimate<sch::ChunkStruct<T,Desc>,sch::MacroStruct<T,Desc>>().get()<<std::endl;
 
 	auto eo_aio = saw::setup_async_io();
@@ -230,7 +237,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){
 
 	sycl_q.wait();
 	{
-		auto eov = setup_initial_conditions<T,Desc>(*lbm_data_ptr,*lbm_macro_data_ptr);
+		auto eov = setup_initial_conditions<T,Desc>(*lbm_data_ptr,*lbm_macro_data_ptr,lbm_particles_data);
 		if(eov.is_error()){
 			return eov;
 		}
diff --git a/lib/core/c++/boundary.hpp b/lib/core/c++/boundary.hpp
index 4dbbdf8..01ae7b5 100644
--- a/lib/core/c++/boundary.hpp
+++ b/lib/core/c++/boundary.hpp
@@ -211,12 +211,12 @@ public:
 
 		if constexpr (Dir) {
 			dfs_old.at({2u}) = dfs_old.at({1u}) + saw::data<FP>{2.0 / 3.0} * rho_vel_x;
-			dfs_old.at({6u}) = dfs_old.at({5u}) + saw::data<FP>{1.0 / 6.0} * rho_vel_x + saw::data<FP>{0.5} * (dfs_old.at({4u}) - dfs_old.at({3u}));
-			dfs_old.at({8u}) = dfs_old.at({7u}) + saw::data<FP>{1.0 / 6.0} * rho_vel_x + saw::data<FP>{0.5} * (dfs_old.at({3u}) - dfs_old.at({4u}));
+			dfs_old.at({6u}) = dfs_old.at({5u}) + saw::data<FP>{1.0 / 6.0} * rho_vel_x + saw::data<FP>{0.5} * (dfs_old.at({3u}) - dfs_old.at({4u}));
+			dfs_old.at({8u}) = dfs_old.at({7u}) + saw::data<FP>{1.0 / 6.0} * rho_vel_x + saw::data<FP>{0.5} * (dfs_old.at({4u}) - dfs_old.at({3u}));
 		}else if constexpr (not Dir){
 			dfs_old.at({1u}) = dfs_old.at({2u}) - saw::data<FP>{2.0 / 3.0} * rho_vel_x;
-			dfs_old.at({5u}) = dfs_old.at({6u}) - saw::data<FP>{1.0 / 6.0} * rho_vel_x + saw::data<FP>{0.5} * (dfs_old.at({3u}) - dfs_old.at({4u}));
-			dfs_old.at({7u}) = dfs_old.at({8u}) - saw::data<FP>{1.0 / 6.0} * rho_vel_x + saw::data<FP>{0.5} * (dfs_old.at({4u}) - dfs_old.at({3u}));
+			dfs_old.at({5u}) = dfs_old.at({6u}) - saw::data<FP>{1.0 / 6.0} * rho_vel_x + saw::data<FP>{0.5} * (dfs_old.at({4u}) - dfs_old.at({3u}));
+			dfs_old.at({7u}) = dfs_old.at({8u}) - saw::data<FP>{1.0 / 6.0} * rho_vel_x + saw::data<FP>{0.5} * (dfs_old.at({3u}) - dfs_old.at({4u}));
 		}
 	}
 };
diff --git a/lib/core/c++/particle/particle.hpp b/lib/core/c++/particle/particle.hpp
index a95a173..c0d115f 100644
--- a/lib/core/c++/particle/particle.hpp
+++ b/lib/core/c++/particle/particle.hpp
@@ -1,8 +1,10 @@
 #pragma once
 
-#include <forstio/codec/data.hpp>
-#include <forstio/codec/data_math.hpp>
 #include <forstio/codec/math.hpp>
+#include <forstio/codec/data_math.hpp>
+#include <forstio/codec/data.hpp>
+
+#include "../iterator.hpp"
 
 namespace kel {
 namespace lbm {
@@ -52,25 +54,20 @@ using Particle = Struct<
 
 template<typename T, uint64_t D, typename CollisionType = ParticleCollisionSpheroid<T>>
 using ParticleGroup = Struct<
-	Member<Array<Float64,D>, "mask">,
+	Member<Array<T,D>, "mask">,
 	Member<CollisionType, "collision">,
 	Member<Scalar<T>, "mass">,
 	Member<Array<Particle<T,D>>, "particles">
 >;
-
-template<typename... PartGroups>
-using ParticleGroupTuple = Tuple<
-	PartGroups...
->;
-
 }
 
 template<typename T, uint64_t D>
 saw::data<sch::ParticleGroup<T,D, sch::ParticleCollisionSpheroid<T>>> create_spheroid_particle_group(
 		saw::data<sch::Scalar<T>> rad_p,
-		saw::data<sch::Scalar<T>> density_p
+		saw::data<sch::Scalar<T>> density_p,
+		const saw::data<sch::UInt64>& mask_resolution
 ){
-	saw::data<sch::ParticleGroup<T,D>,sch::ParticleCollisionSpheroid<T>> part;
+	saw::data<sch::ParticleGroup<T,D,sch::ParticleCollisionSpheroid<T>>> part;
 
 	auto& mask = part.template get<"mask">();
 	auto& collision = part.template get<"collision">();
@@ -92,6 +89,16 @@ saw::data<sch::ParticleGroup<T,D, sch::ParticleCollisionSpheroid<T>>> create_sph
 		static_assert(D == 0u or D > 3u, "Dimensions only supported for Dim 1,2 & 3.");
 	}
 
+	saw::data<sch::FixedArray<sch::UInt64,D>> mask_dims;
+	for(uint64_t i = 0u; i < D; ++i){
+		mask_dims.at({i}) = mask_resolution;
+	}
+
+	mask = {mask_dims};
+	iterator<D>::apply([&](const auto& index){
+		
+	},{},mask_dims);
+
 	return part;
 }