summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
Diffstat
-rw-r--r--examples/poiseulle_particles_2d_gpu/sim.cpp52
-rw-r--r--lib/core/c++/boundary.hpp90
-rw-r--r--lib/core/c++/chunk.hpp13
-rw-r--r--lib/core/c++/collision.hpp1
-rw-r--r--lib/core/c++/lbm.hpp1
-rw-r--r--lib/core/c++/stream.hpp25
6 files changed, 109 insertions, 73 deletions
diff --git a/examples/poiseulle_particles_2d_gpu/sim.cpp b/examples/poiseulle_particles_2d_gpu/sim.cpp
index fc9fc61..e655315 100644
--- a/examples/poiseulle_particles_2d_gpu/sim.cpp
+++ b/examples/poiseulle_particles_2d_gpu/sim.cpp
@@ -32,7 +32,7 @@ saw::error_or<void> setup_initial_conditions(saw::data<sch::ChunkStruct<T,Desc>>
// Set everything as walls
iterator<Desc::D>::apply(
[&](auto& index){
- info_f.at(index).set(2u);
+ info_f.at(index).set(1u);
},
{},
info_f.get_dims(),
@@ -41,7 +41,7 @@ saw::error_or<void> setup_initial_conditions(saw::data<sch::ChunkStruct<T,Desc>>
// Fluid
iterator<Desc::D>::apply(
[&](auto& index){
- info_f.at(index).set(1u);
+ info_f.at(index).set(2u);
},
{},
info_f.get_dims(),
@@ -91,15 +91,59 @@ saw::error_or<void> step(saw::data<sch::Ptr<sch::ChunkStruct<T,Desc>>,encode::Sy
auto& q = dev.get_handle();
auto& info_f = fields.template get<"info">();
+ // auto coll_ev =
q.submit([&](acpp::sycl::handler& h){
+ // 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;
+ component<T,Desc,cmpt::ZouHeHorizontal<true>,encode::Sycl<saw::encode::Native>> flow_in{1.01};
+ component<T,Desc,cmpt::ZouHeHorizontal<true>,encode::Sycl<saw::encode::Native>> flow_out{1.0};
+
h.parallel_for(acpp::sycl::range<Desc::D>{dim_x,dim_y}, [=](acpp::sycl::id<Desc::D> idx){
saw::data<sch::FixedArray<sch::UInt64,Desc::D>> index;
for(uint64_t i = 0u; i < Desc::D; ++i){
index.at({{i}}).set(idx[i]);
}
- collision.apply(fields,index,t_i);
+ auto info = info_f.at(index);
+
+ switch(info.get()){
+ case 0u:
+ break;
+ case 1u:
+ bb.apply(fields,index,t_i);
+ break;
+ case 2u:
+ collision.apply(fields,index,t_i);
+ break;
+ case 3u:
+ flow_in.apply(fields,index,t_i);
+ collision.apply(fields,index,t_i);
+ break;
+ case 4u:
+ flow_out.apply(fields,index,t_i);
+ collision.apply(fields,index,t_i);
+ break;
+ default:
+ break;
+ }
+ });
+ }).wait();
+
+ q.submit([&](acpp::sycl::handler& h){
+ component<T,Desc,cmpt::Stream,encode::Sycl<saw::encode::Native>> stream;
+
+ h.parallel_for(acpp::sycl::range<Desc::D>{dim_x,dim_y}, [=](acpp::sycl::id<Desc::D> idx){
+ saw::data<sch::FixedArray<sch::UInt64,Desc::D>> index;
+ for(uint64_t i = 0u; i < Desc::D; ++i){
+ index.at({{i}}).set(idx[i]);
+ }
+
+ auto info = info_f.at(index);
+
+ if(info.get() > 0u){
+ stream.apply(fields,index,t_i);
+ }
});
}).wait();
@@ -155,7 +199,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){
}
sycl_q.wait();
- for(saw::data<sch::UInt64> i{0u}; i < saw::data<sch::UInt64>{32ul}; ++i){
+ for(saw::data<sch::UInt64> i{0u}; i < saw::data<sch::UInt64>{1024ul}; ++i){
auto eov = step<T,Desc>(lsd_view,i,dev);
if(eov.is_error()){
return eov;
diff --git a/lib/core/c++/boundary.hpp b/lib/core/c++/boundary.hpp
index 8ab2457..ec16edf 100644
--- a/lib/core/c++/boundary.hpp
+++ b/lib/core/c++/boundary.hpp
@@ -44,47 +44,22 @@ public:
* Raw setup
*/
template<typename CellFieldSchema>
- void apply(saw::data<CellFieldSchema, Encode>& field, const saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>>& index, saw::data<sch::UInt64> time_step){
+ 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);
// This is a ref
- auto& cell = field(index);
- auto& dfs_old = (is_even) ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
- // auto& dfs = (not is_even) ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
+ // 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">();
using dfi = df_info<T,Descriptor>;
// Technically use .copy()
- auto df_cpy = dfs_old.copy();
+ auto& dfs_old = dfs_old_f.at(index);
+ auto df_cpy = dfs_old;
- for(uint64_t i = 1u; i < Descriptor::Q; ++i){
- dfs_old({i}) = df_cpy({dfi::opposite_index.at(i)});
- }
- }
-
- template<typename CellStructSchema>
- void apply(
- saw::data<CellStructSchema, Encode>* field,
- const saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>>& meta,
- const saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>>& index,
- saw::data<sch::UInt64> time_step
- ){
- bool is_even = ((time_step.get() % 2u) == 0u);
-
- // This is a ref
- auto& cell = field[0u];
-
- auto& dfs_old = (is_even) ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
- // auto& dfs = (not is_even) ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
-
- using dfi = df_info<T,Descriptor>;
-
- // Technically use .copy()
- auto df_cpy = dfs_old.copy();
-
- for(uint64_t i = 1u; i < Descriptor::Q; ++i){
- dfs_old({i}) = df_cpy({dfi::opposite_index.at(i)});
+ for(uint64_t i = 0u; i < Descriptor::Q; ++i){
+ dfs_old.at({i}) = df_cpy.at({dfi::opposite_index.at(i)});
}
}
};
@@ -104,59 +79,46 @@ public:
template<typename FP, typename Descriptor, bool Dir, typename Encode>
class component<FP, Descriptor, cmpt::ZouHeHorizontal<Dir>, Encode> final {
private:
- saw::data<FP> pressure_setting_;
saw::data<FP> rho_setting_;
public:
- component(const saw::data<FP>& pressure_setting__):
- pressure_setting_{pressure_setting__},
- rho_setting_{pressure_setting__ * df_info<FP,Descriptor>::inv_cs2}
+ component(const saw::data<FP>& density_setting__):
+ rho_setting_{density_setting__}
{}
template<typename CellFieldSchema>
- void apply(saw::data<CellFieldSchema, Encode>& field, saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>> index, saw::data<sch::UInt64> time_step){
+ void apply(const saw::data<CellFieldSchema, Encode>& field, saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>> index, saw::data<sch::UInt64> time_step) const {
using dfi = df_info<FP,Descriptor>;
+ constexpr int known_dir = Dir ? 1 : -1;
bool is_even = ((time_step.get() % 2) == 0);
- auto& cell = field(index);
- auto& info = cell.template get<"info">();
- if(info({0u}).get() == 0u){
- return;
- }
- auto& dfs_old = (is_even) ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
- // auto& dfs = (not is_even) ? cell.template get<"dfs_old">() : cell.template get<"dfs">();
+ auto& info_f = field.template get<"info">();
+
+ // 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">();
+ auto& dfs_old = dfs_old_f.at(index);
/**
* Sum all known DFs
*/
saw::data<FP> sum_df{0};
for(saw::data<sch::UInt64> k{0u}; k < saw::data<sch::UInt64>{Descriptor::Q}; ++k){
auto c_k = dfi::directions[k.get()];
- auto& cell_n = field({{index.at({0u})+c_k[0u], index.at({1u})+c_k[1u]}});
- auto& info_n = cell_n.template get<"info">();
- auto info_n_val = info_n({0u});
- auto k_opp = dfi::opposite_index[k.get()];
- if(info_n_val.get() > 0u){
- sum_df += dfs_old({k_opp});
+ if(c_k[0u]*known_dir <= 0){
+ sum_df += dfs_old.at({k});
}
}
/**
* Get the sum of the unknown dfs and precalculate the direction
*/
- constexpr int known_dir = Dir ? 1 : -1;
auto sum_unknown_dfs = (rho_setting_ - sum_df) * saw::data<FP>{known_dir};
for(saw::data<sch::UInt64> k{0u}; k < saw::data<sch::UInt64>{Descriptor::Q}; ++k){
auto c_k = dfi::directions[k.get()];
- auto& cell_n = field({{index.at({0u})+c_k[0u], index.at({1u})+c_k[1u]}});
- auto& info_n = cell_n.template get<"info">();
- auto info_n_val = info_n({0u});
- // auto k_opp = dfi::opposite_index[k.get()];
-
- if(info_n_val.get() > 0u){
- sum_unknown_dfs += dfs_old({k}) * c_k[0u];
+ if(c_k[0u]*known_dir > 0){
+ sum_unknown_dfs += dfs_old.at({k}) * c_k[0u];
}
}
@@ -165,13 +127,13 @@ public:
static_assert(Descriptor::D == 2u and Descriptor::Q == 9u, "Some parts are hard coded sadly");
if constexpr (Dir) {
- dfs_old({2u}) = dfs_old({1u}) + saw::data<FP>{2.0 / 3.0} * rho_setting_ * vel_x;
- dfs_old({6u}) = dfs_old({5u}) + saw::data<FP>{1.0 / 6.0} * rho_setting_ * vel_x + saw::data<FP>{0.5} * (dfs_old({4u}) - dfs_old({3u}));
- dfs_old({8u}) = dfs_old({7u}) + saw::data<FP>{1.0 / 6.0} * rho_setting_ * vel_x + saw::data<FP>{0.5} * (dfs_old({3u}) - dfs_old({4u}));
+ dfs_old.at({2u}) = dfs_old.at({1u}) + saw::data<FP>{2.0 / 3.0} * rho_setting_ * vel_x;
+ dfs_old.at({6u}) = dfs_old.at({5u}) + saw::data<FP>{1.0 / 6.0} * rho_setting_ * 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_setting_ * vel_x + saw::data<FP>{0.5} * (dfs_old.at({3u}) - dfs_old.at({4u}));
}else if constexpr (not Dir){
- dfs_old({1u}) = dfs_old({2u}) - saw::data<FP>{2.0 / 3.0} * rho_setting_ * vel_x;
- dfs_old({5u}) = dfs_old({6u}) - saw::data<FP>{1.0 / 6.0} * rho_setting_ * vel_x + saw::data<FP>{0.5} * (dfs_old({3u}) - dfs_old({4u}));
- dfs_old({7u}) = dfs_old({8u}) - saw::data<FP>{1.0 / 6.0} * rho_setting_ * vel_x + saw::data<FP>{0.5} * (dfs_old({4u}) - dfs_old({3u}));
+ dfs_old.at({1u}) = dfs_old.at({2u}) - saw::data<FP>{2.0 / 3.0} * rho_setting_ * vel_x;
+ dfs_old.at({5u}) = dfs_old.at({6u}) - saw::data<FP>{1.0 / 6.0} * rho_setting_ * 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_setting_ * vel_x + saw::data<FP>{0.5} * (dfs_old.at({4u}) - dfs_old.at({3u}));
}
}
};
diff --git a/lib/core/c++/chunk.hpp b/lib/core/c++/chunk.hpp
index 5d20faa..a1f2451 100644
--- a/lib/core/c++/chunk.hpp
+++ b/lib/core/c++/chunk.hpp
@@ -35,6 +35,7 @@ using SuperChunk = Array<ChunkSchema,Dim>;
}
namespace saw {
+
template<typename Sch, uint64_t Ghost, uint64_t... Sides, typename Encode>
class data<kel::lbm::sch::Chunk<Sch,Ghost,Sides...>,Encode> final {
public:
@@ -57,19 +58,21 @@ public:
return data<InnerSchema,Encode>::get_dims();
}
- data<ValueSchema, Encode>& at(const data<schema::FixedArray<schema::UInt64,sizeof...(Sides)>>& index){
+ static constexpr auto to_ghost_index(const data<schema::FixedArray<schema::UInt64,sizeof...(Sides)>>& index){
std::decay_t<decltype(index)> ind;
for(uint64_t i = 0u; i < sizeof...(Sides); ++i){
ind.at({i}) = index.at({i}) + Ghost;
}
+ return ind;
+ }
+
+ data<ValueSchema, Encode>& at(const data<schema::FixedArray<schema::UInt64,sizeof...(Sides)>>& index){
+ std::decay_t<decltype(index)> ind = to_ghost_index(index);
return values_.at(ind);
}
const data<ValueSchema, Encode>& at(const data<schema::FixedArray<schema::UInt64,sizeof...(Sides)>>& index) const {
- std::decay_t<decltype(index)> ind;
- for(uint64_t i = 0u; i < sizeof...(Sides); ++i){
- ind.at({i}) = index.at({i}) + Ghost;
- }
+ std::decay_t<decltype(index)> ind = to_ghost_index(index);
return values_.at(ind);
}
diff --git a/lib/core/c++/collision.hpp b/lib/core/c++/collision.hpp
index 85c3af9..ed26a08 100644
--- a/lib/core/c++/collision.hpp
+++ b/lib/core/c++/collision.hpp
@@ -3,6 +3,7 @@
#include "macroscopic.hpp"
#include "component.hpp"
#include "equilibrium.hpp"
+#include "chunk.hpp"
namespace kel {
namespace lbm {
diff --git a/lib/core/c++/lbm.hpp b/lib/core/c++/lbm.hpp
index 00f153a..5334e4e 100644
--- a/lib/core/c++/lbm.hpp
+++ b/lib/core/c++/lbm.hpp
@@ -13,6 +13,7 @@
#include "equilibrium.hpp"
#include "iterator.hpp"
#include "macroscopic.hpp"
+#include "stream.hpp"
#include "write_vtk.hpp"
#include "util.hpp"
diff --git a/lib/core/c++/stream.hpp b/lib/core/c++/stream.hpp
new file mode 100644
index 0000000..d217373
--- /dev/null
+++ b/lib/core/c++/stream.hpp
@@ -0,0 +1,25 @@
+#pragma once
+
+#include "component.hpp"
+
+namespace kel {
+namespace lbm {
+namespace cmpt {
+struct Stream {};
+}
+
+template<typename T, typename Descriptor, typename Encode>
+class component<T,Descriptor, cmpt::Stream, Encode> final {
+private:
+public:
+ static constexpr saw::string_literal name = "streaming";
+ static constexpr saw::string_literal after = "collide";
+ static constexpr saw::string_literal before = "";
+
+ 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 {
+
+ }
+};
+}
+}
generated by cgit v1.2.3 (git 2.25.1) at 2026-02-17 02:28:51 +0000