diff options
| -rw-r--r-- | examples/schaefer_turek_durst_krause_rannbacher/sim.cpp | 2 | ||||
| -rw-r--r-- | examples/schaefer_turek_durst_krause_rannbacher/step.hpp | 23 | ||||
| -rw-r--r-- | modules/core/c++/fplbm.hpp | 2 | ||||
| -rw-r--r-- | modules/core/c++/macroscopic.hpp | 25 |
4 files changed, 44 insertions, 8 deletions
diff --git a/examples/schaefer_turek_durst_krause_rannbacher/sim.cpp b/examples/schaefer_turek_durst_krause_rannbacher/sim.cpp index 4d015c9..ff6a43a 100644 --- a/examples/schaefer_turek_durst_krause_rannbacher/sim.cpp +++ b/examples/schaefer_turek_durst_krause_rannbacher/sim.cpp @@ -48,7 +48,7 @@ saw::error_or<void> lbm_main(const saw::data<args::LbmArgs>& args){ // delta_x {{0.41 / dim_y}}, // delta_t - {{1.067e-4}} + {{1.067e-4 * 0.5}} }; print_lbm_meta<T,Desc>(conv,{1e-3},{1e-4},{0.5 * dim_y}); diff --git a/examples/schaefer_turek_durst_krause_rannbacher/step.hpp b/examples/schaefer_turek_durst_krause_rannbacher/step.hpp index cbdfd10..9748ec8 100644 --- a/examples/schaefer_turek_durst_krause_rannbacher/step.hpp +++ b/examples/schaefer_turek_durst_krause_rannbacher/step.hpp @@ -8,6 +8,7 @@ namespace lbm { template<typename T, typename Desc, typename Coll> saw::error_or<void> step( const converter<T>& conv, + const saw::data<sch::Scalar<T>>& tau, 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::Ptr<sch::ParticleSpheroidGroup<T,Desc>>,encode::Sycl<saw::encode::Native>>& particles, @@ -48,7 +49,7 @@ saw::error_or<void> step( }).wait(); q.submit([&](acpp::sycl::handler& h){ - component<T,Desc,cmpt::Hlbm,encode::Sycl<saw::encode::Native>> collision{1.0}; + component<T,Desc,cmpt::Hlbm,encode::Sycl<saw::encode::Native>> collision{0.75}; component<T,Desc,cmpt::BounceBack,encode::Sycl<saw::encode::Native>> bb; component<T,Desc,cmpt::ZouHeHorizontal<false>,encode::Sycl<saw::encode::Native>> flow_out{1.0}; @@ -74,9 +75,19 @@ saw::error_or<void> step( component<T,Desc,cmpt::ZouHeVelocityX<true>,encode::Sycl<saw::encode::Native>> flow_in{ [&]() -> saw::data<sch::Vector<T,Desc::D>> { saw::data<sch::Vector<T,Desc::D>> vel; + saw::data<T> scale; + { + scale.set([&](){ + if(t_i.get() > 2048u){ + return 1.0f; + } + + return t_i.get() / 2048.0f; + }()); + } { auto y_si = conv.meter_lbm_to_si(index.at({{1u}}).template cast_to<T>()).handle(); - vel.at({{0u}}) = conv.velocity_si_to_lbm({{static_cast<typename saw::native_data_type<T>::type>((1.2 / 0.41) * y_si.get() - (1.2 / (0.41*0.41)) * y_si.get() * y_si.get())}}).handle(); + vel.at({{0u}}) = conv.velocity_si_to_lbm({{static_cast<typename saw::native_data_type<T>::type>((1.2 / 0.41) * y_si.get() - (1.2 / (0.41*0.41)) * y_si.get() * y_si.get())}}).handle() * scale; vel.at({{1u}}) = 0.0f; } return vel; @@ -123,8 +134,8 @@ saw::error_or<void> step( // auto coll_ev = q.submit([&](acpp::sycl::handler& h){ - component<T,Desc,cmpt::BGK,encode::Sycl<saw::encode::Native>> bgk{1.0}; - component<T,Desc,cmpt::FpLbm,encode::Sycl<saw::encode::Native>> collision{1.0}; + component<T,Desc,cmpt::BGK,encode::Sycl<saw::encode::Native>> bgk{0.75}; + component<T,Desc,cmpt::FpLbm,encode::Sycl<saw::encode::Native>> collision{0.75}; component<T,Desc,cmpt::BounceBack,encode::Sycl<saw::encode::Native>> bb; component<T,Desc,cmpt::ZouHeHorizontal<false>,encode::Sycl<saw::encode::Native>> flow_out{1.0}; @@ -199,8 +210,8 @@ saw::error_or<void> step( // auto coll_ev = q.submit([&](acpp::sycl::handler& h){ - component<T,Desc,cmpt::BGK,encode::Sycl<saw::encode::Native>> bgk{1.0}; - component<T,Desc,cmpt::Psm,encode::Sycl<saw::encode::Native>> collision{1.0}; + component<T,Desc,cmpt::BGK,encode::Sycl<saw::encode::Native>> bgk{0.75}; + component<T,Desc,cmpt::Psm,encode::Sycl<saw::encode::Native>> collision{0.75}; component<T,Desc,cmpt::BounceBack,encode::Sycl<saw::encode::Native>> bb; component<T,Desc,cmpt::ZouHeHorizontal<false>,encode::Sycl<saw::encode::Native>> flow_out{1.0}; diff --git a/modules/core/c++/fplbm.hpp b/modules/core/c++/fplbm.hpp index c9ed6e9..d8c075b 100644 --- a/modules/core/c++/fplbm.hpp +++ b/modules/core/c++/fplbm.hpp @@ -92,7 +92,7 @@ public: auto& force_f = macros.template get<"force">(); auto& force = force_f.at(index); - vel = vel + force * half; + vel = vel + force * half / rho; // compute_rho_u<T,Descriptor>(dfs,rho,vel); auto eq = equilibrium<T,Descriptor>(rho,vel); diff --git a/modules/core/c++/macroscopic.hpp b/modules/core/c++/macroscopic.hpp index 768fef0..9db1538 100644 --- a/modules/core/c++/macroscopic.hpp +++ b/modules/core/c++/macroscopic.hpp @@ -33,5 +33,30 @@ void compute_rho_u ( vel().at({{i}}) = vel().at({{i}}) / rho().at({}); } } + +/** + * Calculate the macroscopic variables rho and u in Lattice Units. + */ +template<typename T, typename Desc> +void compute_rho_momentum ( + const saw::data<sch::FixedArray<T,Desc::Q>>& dfs, + saw::ref<saw::data<sch::Scalar<T>>> rho, + saw::ref<saw::data<sch::Vector<T,Desc::D>>> mom + ) +{ + using dfi = df_info<T, Desc>; + + rho().at({}).set(0); + for(uint64_t i = 0; i < Desc::D; ++i){ + mom().at({{i}}).set(0); + } + + for(size_t j = 0; j < Desc::Q; ++j){ + rho().at({}) = rho().at({}) + dfs.at({j}); + for(size_t i = 0; i < Desc::D; ++i){ + mom().at({{i}}) = mom().at({{i}}) + saw::data<T>{static_cast<typename saw::native_data_type<T>::type>(dfi::directions[j][i])} * dfs.at({j}); + } + } +} } } |