diff options
| -rw-r--r-- | c++/macroscopic.hpp | 31 | ||||
| -rw-r--r-- | c++/write_vtk.hpp | 90 | ||||
| -rw-r--r-- | tests/particles.cpp | 8 | ||||
| -rw-r--r-- | tests/vtk_write.cpp | 16 |
4 files changed, 118 insertions, 27 deletions
diff --git a/c++/macroscopic.hpp b/c++/macroscopic.hpp index 43c727b..e126bbe 100644 --- a/c++/macroscopic.hpp +++ b/c++/macroscopic.hpp @@ -9,7 +9,7 @@ namespace lbm { */ template<typename T, typename Desc> void compute_rho_u ( - saw::data<sch::Cell<T, Desc, 0, 0, 1>>& dfs, + const saw::data<sch::Cell<T, Desc, 0, 0, 1>>& dfs, typename saw::native_data_type<T>::type& rho, std::array<typename saw::native_data_type<T>::type, 2>& vel ) @@ -30,5 +30,34 @@ void compute_rho_u ( vel[i] /= rho; } } + +/** + * Calculate the macroscopic variables rho and u in Lattice Units. + */ +template<typename T, typename Desc> +void compute_rho_u ( + const saw::data<sch::Cell<T, Desc, 0, 0, 1>>& dfs, + saw::ref<saw::data<T>> rho, + saw::ref<saw::data<sch::FixedArray<T,Desc::D>>> vel + ) +{ + using dfi = df_info<T, Desc>; + + rho().set(0); + for(uint64_t i = 0; i < Desc::D; ++i){ + vel().at({i}).set(0); + } + + for(size_t j = 0; j < Desc::Q; ++j){ + rho = rho + dfs(j); + for(size_t i = 0; i < Desc::D; ++i){ + vel[i] = vel[i] + saw::data<T>{dfi::directions[j][i]} * dfs(j); + } + } + + for(size_t i = 0; i < Desc::D; ++i){ + vel[i] = vel[i] / rho; + } +} } } diff --git a/c++/write_vtk.hpp b/c++/write_vtk.hpp index fec4ea5..5cbc6c0 100644 --- a/c++/write_vtk.hpp +++ b/c++/write_vtk.hpp @@ -16,50 +16,91 @@ template<typename CellFieldSchema> struct lbm_vtk_writer { }; -template<typename... MemberTypes, saw::string_literal... MemberNames> -struct lbm_vtk_writer<sch::Struct<sch::Member<MemberTypes,MemberNames>...>> { - - static saw::error_or<void> apply(){ +template<typename T, uint64_t D> +struct lbm_vtk_writer<sch::Primitive<T,D>> { + static saw::error_or<void> apply_header(std::ostream& vtk_file, std::string_view name){ + vtk_file<<"SCALARS "<<name<<" float\n"; + vtk_file<<"LOOKUP_TABLE default\n"; + return saw::make_void(); + } + static saw::error_or<void> apply(std::ostream& vtk_file, const saw::data<sch::Primitive<T,D>>& field){ + vtk_file<<field.get()<<"\n"; return saw::make_void(); } }; template<typename T, uint64_t D> struct lbm_vtk_writer<sch::FixedArray<T,D>> { - static saw::error_or<void> apply(std::ostream& vtk_file, std::string_view name){ + static saw::error_or<void> apply_header(std::ostream& vtk_file, std::string_view name){ vtk_file<<"VECTORS "<<name<<" float\n"; + return saw::make_void(); + } + static saw::error_or<void> apply(std::ostream& vtk_file, const saw::data<sch::FixedArray<T,D>>& field){ + static_assert(D > 0, "Non-dimensionality is bad for velocity."); + static_assert(D <= 3, "4th dimension as well. Mostly due to vtk."); + + // vtk_file<<"VECTORS "<<name<<" float\n"; + for(uint64_t i = 0u; i < D; ++i){ + if(i > 0){ + vtk_file<<" "; + } + vtk_file<<field.at({i}).get(); + } + for(uint64_t i = D; i < 3; ++i){ + vtk_file<<" 0"; + } + vtk_file<<"\n"; + return saw::make_void(); } }; -template<typename Desc, typename... StructT, saw::string_literal... StructN> -struct lbm_vtk_writer<sch::CellField<Desc,sch::Struct<sch::Member<StructT,StructN>...>>> { - template<uint64_t i, uint64_t Dep> - static saw::error_or<void> write_i_iterate_d(std::ostream& vtk_file, const saw::data<sch::CellField<Desc,sch::Struct<sch::Member<StructT,StructN>...>>>& field, saw::data<sch::FixedArray<sch::UInt64,Desc::D>>& index){ - if constexpr (Dep == Desc::D){ +template<uint64_t Dim, typename... StructT, saw::string_literal... StructN> +struct lbm_vtk_writer<sch::Array<sch::Struct<sch::Member<StructT,StructN>...> , Dim>> { + template<uint64_t i, uint64_t Dep> + static saw::error_or<void> write_i_iterate_d(std::ostream& vtk_file, const saw::data<sch::Array<sch::Struct<sch::Member<StructT,StructN>...>,Dim>>& field, saw::data<sch::FixedArray<sch::UInt64,Dim>>& index){ + constexpr auto Lit = saw::parameter_key_pack_type<i, StructN...>::literal; + using Type = typename saw::parameter_pack_type<i,StructT...>::type; + + if constexpr (Dep == Dim){ + return lbm_vtk_writer<Type>::apply(vtk_file, field.at(index).template get<Lit>()); + } else { + // Dep < Dim // I hope + static_assert(Dep < Dim, "Don't fall into this case"); + for(index.at({Dep}) = 0; index.at({Dep}) < field.get_dims().at({Dep}); ++index.at({Dep})){ + auto eov = write_i_iterate_d<i,Dep+1>(vtk_file, field, index); + if(eov.is_error()){ + return eov; + } + } } + return saw::make_void(); } template<uint64_t i> - static saw::error_or<void> write_i(std::ostream& vtk_file, const saw::data<sch::CellField<Desc,sch::Struct<sch::Member<StructT,StructN>...>>>& field){ + 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.meta(); - saw::data<sch::FixedArray<sch::UInt64,Desc::D>> index; - for(saw::data<sch::UInt64> it{0}; it.get() < Desc::D; ++it){ + 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); } // vtk_file write? // Data header { - auto eov = lbm_vtk_writer<typename saw::parameter_pack_type<i,StructT...>::type>::apply(vtk_file, saw::parameter_key_pack_type<i, StructN...>::literal.view()); + + auto eov = lbm_vtk_writer<typename saw::parameter_pack_type<i,StructT...>::type>::apply_header(vtk_file, saw::parameter_key_pack_type<i, StructN...>::literal.view()); + if(eov.is_error()){ + return eov; + } + } return write_i_iterate_d<i,0u>(vtk_file, field, index); - - return saw::make_void(); } template<uint64_t i> - static saw::error_or<void> iterate_i(std::ostream& vtk_file, const saw::data<sch::CellField<Desc,sch::Struct<sch::Member<StructT,StructN>...>>>& field){ + static saw::error_or<void> iterate_i(std::ostream& vtk_file, const saw::data<sch::Array<sch::Struct<sch::Member<StructT,StructN>...>, Dim>>& field){ + constexpr auto Lit = saw::parameter_key_pack_type<i, StructN...>::literal; { auto eov = write_i<i>(vtk_file, field); if(eov.is_error()){ @@ -73,27 +114,28 @@ struct lbm_vtk_writer<sch::CellField<Desc,sch::Struct<sch::Member<StructT,Struct } static saw::error_or<void> apply(std::ostream& vtk_file, - const saw::data<sch::CellField<Desc,sch::Struct<sch::Member<StructT,StructN>...>>>& field){ + const saw::data<sch::Array<sch::Struct<sch::Member<StructT,StructN>...>, Dim>>& field){ vtk_file <<"# vtk DataFile Version 3.0\n" <<"LBM File\n" <<"ASCII\n" - <<"DATASET STRUCTURED POINTS\n" + <<"DATASET STRUCTURED_POINTS\n" <<"SPACING 1.0 1.0 1.0\n" <<"ORIGIN 0.0 0.0 0.0\n" ; - auto meta = field.meta(); + auto meta = field.get_dims(); saw::data<sch::UInt64> pd_size{1u}; // DIMENSIONS + { - vtk_file << "DIMENSIONS "; - for(saw::data<sch::UInt64> i{0u}; i.get() < Desc::D; ++i){ + vtk_file << "DIMENSIONS"; + for(saw::data<sch::UInt64> i{0u}; i.get() < Dim; ++i){ pd_size = pd_size * meta.at(i); vtk_file << " " << meta.at(i).get(); } - for(saw::data<sch::UInt64> i{Desc::D}; i.get() < 3u; ++i){ + for(saw::data<sch::UInt64> i{Dim}; i.get() < 3u; ++i){ vtk_file << " 1"; } diff --git a/tests/particles.cpp b/tests/particles.cpp index 6332195..873b8ad 100644 --- a/tests/particles.cpp +++ b/tests/particles.cpp @@ -31,4 +31,12 @@ SAW_TEST("Minor Test for mask"){ //saw::data<sch::Array<sch::T,2>> reference_mask{{{4+2,4+2}}}; //reference_mask.at({{0,0}}); } + +SAW_TEST("Verlet integration test"){ + using namespace kel; + lbm::particle_system<sch::T,2,sch::Particle<sch::T,2>> system; + + // system.step(); + +} } diff --git a/tests/vtk_write.cpp b/tests/vtk_write.cpp index 0fcc169..49c9c52 100644 --- a/tests/vtk_write.cpp +++ b/tests/vtk_write.cpp @@ -29,6 +29,12 @@ using CellStruct = Struct< Member<CellInfo<Desc>, "info"> >; +template<typename T, uint64_t D> +using MacroStruct = Struct< + Member<FixedArray<T,D>, "velocity">, + Member<T, "pressure"> +>; + } SAW_TEST("VTK Write test example"){ @@ -38,11 +44,17 @@ SAW_TEST("VTK Write test example"){ std::stringstream sstream; - saw::data<sch::CellField<sch::D2Q5,sch::CellStruct<sch::D2Q5>>> cells; + saw::data<sch::Array<sch::MacroStruct<sch::T,2>, 2>> cells{{{2u,2u}}}; - auto eov = lbm::impl::lbm_vtk_writer<sch::CellField<sch::D2Q5,sch::CellStruct<sch::D2Q5>>>::apply(sstream, cells); + auto& cell_0 = cells.at({{{0,0}}}); + cell_0.template get<"velocity">()= {{0.5,-0.1}}; + cell_0.template get<"pressure">().set({1.1}); + + auto eov = lbm::impl::lbm_vtk_writer<sch::Array<sch::MacroStruct<sch::T,2>, 2>>::apply(sstream, cells); SAW_EXPECT(eov.is_value(), "vtk writer failed to write"); + std::cout<<sstream.str()<<std::endl; + // using Type = typename parameter_pack_type<i,T...>::type; } } |