From 24bf28a8fb9cc8c3a90b77de9b60728bece7885d Mon Sep 17 00:00:00 2001 From: "Claudius \"keldu\" Holeksa" Date: Sat, 18 Oct 2025 18:01:14 +0200 Subject: Moving project structure for more less compilation --- lib/c++/collision.hpp | 129 ++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 129 insertions(+) create mode 100644 lib/c++/collision.hpp (limited to 'lib/c++/collision.hpp') diff --git a/lib/c++/collision.hpp b/lib/c++/collision.hpp new file mode 100644 index 0000000..9ab542b --- /dev/null +++ b/lib/c++/collision.hpp @@ -0,0 +1,129 @@ +#pragma once + +#include "macroscopic.hpp" +#include "component.hpp" +#include "equilibrium.hpp" + +namespace kel { +namespace lbm { +namespace cmpt { +struct BGK {}; +struct BGKGuo {}; +} + +/** + * Standard BGK collision operator for LBM + */ +template +class component { +private: + typename saw::native_data_type::type relaxation_; + saw::data frequency_; +public: + component(typename saw::native_data_type::type relaxation__): + relaxation_{relaxation__}, + frequency_{typename saw::native_data_type::type(1) / relaxation_} + {} + + using Component = cmpt::BGK; + + /** + * 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 + void apply(saw::data& field, saw::data> index, uint64_t time_step){ + bool is_even = ((time_step % 2) == 0); + 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">(); + + saw::data rho; + saw::data> vel; + compute_rho_u(dfs_old,rho,vel); + auto eq = equilibrium(rho,vel); + + for(uint64_t i = 0u; i < Descriptor::Q; ++i){ + dfs_old({i}) = dfs_old({i}) + frequency_ * (eq.at(i) - dfs_old({i})); + // dfs_old({i}).set(dfs_old({i}).get() + (1.0 / relaxation_) * (eq.at(i).get() - dfs_old({i}).get())); + } + } +}; + +template +class component { +private: + typename saw::native_data_type::type relaxation_; + saw::data frequency_; +public: + component(typename saw::native_data_type::type relaxation__): + relaxation_{relaxation__}, + frequency_{typename saw::native_data_type::type(1) / relaxation_} + {} + + using Component = cmpt::BGKGuo; + + using access = cmpt::access_tuple< + cmpt::access<"dfs", 1, true, true, true>, + cmpt::access<"force", 0, true, false> + >; + + /** + * 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"; + + template + void apply(saw::data& field, saw::data> index, uint64_t time_step){ + bool is_even = ((time_step % 2) == 0); + 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">(); + + saw::data rho; + saw::data> vel; + compute_rho_u(dfs_old,rho,vel); + auto eq = equilibrium(rho,vel); + + using dfi = df_info; + + auto& force = cell.template get<"force">(); + + for(uint64_t i = 0u; i < Descriptor::Q; ++i){ + // saw::data ci_min_u{0}; + saw::data ci_dot_u{0}; + for(uint64_t d = 0u; d < Descriptor::D; ++d){ + ci_dot_u = ci_dot_u + vel.at({d}) * saw::data{static_cast::type>(dfi::directions[i][d])}; + } + saw::data F_i{0};// = saw::data{dfi::weights[i]} * (ci_dot_F * dfi::inv_cs2 + ci_dot_u * ci_dot_F * (dfi::inv_cs2 * dfi::inv_cs2)); + + for(uint64_t d = 0u; d < Descriptor::D; ++d){ + F_i = F_i + + saw::data{dfi::weights[i]} * ((saw::data{static_cast::type>(dfi::directions[i][d])} + - vel.at({d}) ) * dfi::inv_cs2 + ci_dot_u * saw::data{static_cast::type>(dfi::directions[i][d])} * dfi::inv_cs2 * dfi::inv_cs2 ) * force({d}); + } + + + dfs_old({i}) = dfs_old({i}) + frequency_ * (eq.at(i) - dfs_old({i}) ) + F_i * (saw::data{1} - saw::data{0.5f} * frequency_); + } + } +}; +} +} -- cgit v1.2.3