diff options
Diffstat (limited to 'c++')
| -rw-r--r-- | c++/boundary.hpp | 2 | ||||
| -rw-r--r-- | c++/collision.hpp | 32 | ||||
| -rw-r--r-- | c++/particle/particle.hpp | 3 | ||||
| -rw-r--r-- | c++/statistics.hpp | 11 | ||||
| -rw-r--r-- | c++/write_vtk.hpp | 2 |
5 files changed, 31 insertions, 19 deletions
diff --git a/c++/boundary.hpp b/c++/boundary.hpp index 3cf6f3f..bc4859c 100644 --- a/c++/boundary.hpp +++ b/c++/boundary.hpp @@ -51,7 +51,7 @@ public: 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 = (not is_even) ? cell.template get<"dfs_old">() : cell.template get<"dfs">(); using dfi = df_info<T,Descriptor>; diff --git a/c++/collision.hpp b/c++/collision.hpp index df9a734..cd6cf14 100644 --- a/c++/collision.hpp +++ b/c++/collision.hpp @@ -17,11 +17,11 @@ template<typename T, typename Descriptor> class component<T, Descriptor, cmpt::BGK> { private: typename saw::native_data_type<T>::type relaxation_; - typename saw::native_data_type<T>::type frequency_; + saw::data<T> frequency_; public: component(typename saw::native_data_type<T>::type relaxation__): relaxation_{relaxation__}, - frequency_{1 / relaxation_} + frequency_{typename saw::native_data_type<T>::type(1) / relaxation_} {} using Component = cmpt::BGK; @@ -57,7 +57,8 @@ public: auto eq = equilibrium<T,Descriptor>(rho,vel); for(uint64_t i = 0u; i < Descriptor::Q; ++i){ - dfs_old({i}).set(dfs_old({i}).get() + (1.0 / relaxation_) * (eq.at(i).get() - dfs_old({i}).get())); + 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())); } } }; @@ -66,11 +67,11 @@ template<typename T, typename Descriptor> class component<T, Descriptor, cmpt::BGKGuo> { private: typename saw::native_data_type<T>::type relaxation_; - typename saw::native_data_type<T>::type frequency_; + saw::data<T> frequency_; public: component(typename saw::native_data_type<T>::type relaxation__): relaxation_{relaxation__}, - frequency_{1 / relaxation_} + frequency_{typename saw::native_data_type<T>::type(1) / relaxation_} {} using Component = cmpt::BGKGuo; @@ -89,7 +90,7 @@ public: template<typename CellFieldSchema> void apply(saw::data<CellFieldSchema>& field, saw::data<sch::FixedArray<sch::UInt64, Descriptor::D>> index, uint64_t time_step){ - bool is_even = ((time_step & 2) == 0); + 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">(); @@ -105,22 +106,21 @@ public: auto& force = cell.template get<"force">(); for(uint64_t i = 0u; i < Descriptor::Q; ++i){ + saw::data<T> ci_min_u{0}; saw::data<T> ci_dot_u{0}; - saw::data<T> ci_dot_F{0}; - saw::data<T> u_dot_F{0}; for(uint64_t d = 0u; d < Descriptor::D; ++d){ - ci_dot_u = ci_dot_u + vel.at({d}) * dfi::directions[i][d]; - ci_dot_F = ci_dot_F + force({d}) * dfi::directions[i][d]; - u_dot_F = u_dot_F + vel.at({d}) * force({d}); + ci_dot_u = ci_dot_u + vel.at({d}) * saw::data<T>{static_cast<typename saw::native_data_type<T>::type>(dfi::directions[i][d])}; } + saw::data<T> F_i{0};// = saw::data<T>{dfi::weights[i]} * (ci_dot_F * dfi::inv_cs2 + ci_dot_u * ci_dot_F * (dfi::inv_cs2 * dfi::inv_cs2)); - // auto df_old = dfs_old({i}); - saw::data<T> F_i = saw::data<T>{dfi::weights[i]} * (ci_dot_F / dfi::cs2 + ci_dot_u * ci_dot_F / (dfi::cs2 * dfi::cs2)); - if(F_i.get() != 0){ - std::cerr<<"AAAH"; + for(uint64_t d = 0u; d < Descriptor::D; ++d){ + F_i = F_i + + saw::data<T>{dfi::weights[i]} * ((saw::data<T>{static_cast<typename saw::native_data_type<T>::type>(dfi::directions[i][d])} + - vel.at({d}) ) * dfi::inv_cs2 + ci_dot_u * saw::data<T>{static_cast<typename saw::native_data_type<T>::type>(dfi::directions[i][d])} * dfi::inv_cs2 * dfi::inv_cs2 ) * force({d}); } - dfs_old({i}) = dfs_old({i}) + saw::data<T>{frequency_} * (eq.at(i) - dfs_old({i}) ) - F_i * saw::data<T>{1.0f - 0.5f * frequency_}; + + dfs_old({i}) = dfs_old({i}) + frequency_ * (eq.at(i) - dfs_old({i}) ) + F_i * (saw::data<T>{1} - saw::data<T>{0.5f} * frequency_); } } }; diff --git a/c++/particle/particle.hpp b/c++/particle/particle.hpp index a55b06e..ae7769a 100644 --- a/c++/particle/particle.hpp +++ b/c++/particle/particle.hpp @@ -27,7 +27,7 @@ using ParticleMask = Struct< template<typename T, uint64_t D> using Particle = Struct< Member<ParticleRigidBody<T,D>, "rigid_body">, - Member<ParticleMask<Float32,D>, "mask">, + Member<ParticleMask<T,D>, "mask">, Member<T, "size"> >; } @@ -84,6 +84,7 @@ public: template<typename LbmLattice> void update_particle_border(saw::data<LbmLattice>& latt){ + (void) latt; for(auto& iter : particles_){ auto& par = iter; diff --git a/c++/statistics.hpp b/c++/statistics.hpp new file mode 100644 index 0000000..c07ccb7 --- /dev/null +++ b/c++/statistics.hpp @@ -0,0 +1,11 @@ +#pragma once + +namespace kel { +namespace lbm { +template<typename T> +class statistics { +private: +public: +}; +} +} diff --git a/c++/write_vtk.hpp b/c++/write_vtk.hpp index eb2545a..c9ecc99 100644 --- a/c++/write_vtk.hpp +++ b/c++/write_vtk.hpp @@ -81,7 +81,7 @@ struct lbm_vtk_writer<sch::Array<sch::Struct<sch::Member<StructT,StructN>...> , template<uint64_t i> static saw::error_or<void> write_i(std::ostream& vtk_file, const saw::data<sch::Array<sch::Struct<sch::Member<StructT,StructN>...>,Dim>>& field){ - auto meta = field.get_dims(); + // auto meta = field.get_dims(); saw::data<sch::FixedArray<sch::UInt64,Dim>> index; for(saw::data<sch::UInt64> it{0}; it.get() < Dim; ++it){ index.at({0u}).set(0u); |