From 34bcca5a8d9234febbed9c11c2c801c3e5d63581 Mon Sep 17 00:00:00 2001 From: "Claudius \"keldu\" Holeksa" Date: Thu, 23 Jul 2026 23:20:05 +0200 Subject: Fixed HLBM galilean issue --- examples/moving_poiseulle_particles_2d_hlbm_gpu/sim.cpp | 5 ++--- modules/core/c++/hlbm.hpp | 10 ++++++++-- 2 files changed, 10 insertions(+), 5 deletions(-) diff --git a/examples/moving_poiseulle_particles_2d_hlbm_gpu/sim.cpp b/examples/moving_poiseulle_particles_2d_hlbm_gpu/sim.cpp index 97ecfde..f89f5da 100644 --- a/examples/moving_poiseulle_particles_2d_hlbm_gpu/sim.cpp +++ b/examples/moving_poiseulle_particles_2d_hlbm_gpu/sim.cpp @@ -201,7 +201,7 @@ saw::error_or step( index.at({{i}}).set(idx[i]); } - hlbm_one_part.apply(fields,macros,particles,index,t_i,{{128u}}); + hlbm_one_part.apply(fields,macros,particles,index,t_i,{16u}); }); }).wait(); @@ -380,14 +380,13 @@ saw::error_or lbm_main(int argc, char** argv){ } sycl_q.wait(); - if( i.get() % 32u == 33u){ + if( i.get() % 32u == 0u){ { auto eov = dev.copy_to_host(lbm_sycl_macro_data,*lbm_macro_data_ptr); if(eov.is_error()){ return eov; } } - { auto eov = write_vtk_file(out_dir,"m",i.get(), *lbm_macro_data_ptr); if(eov.is_error()){ diff --git a/modules/core/c++/hlbm.hpp b/modules/core/c++/hlbm.hpp index 4b1f2b4..bb21ae5 100644 --- a/modules/core/c++/hlbm.hpp +++ b/modules/core/c++/hlbm.hpp @@ -123,6 +123,7 @@ public: auto& part_spheroid_group = part_group; auto& mvel_f = macros.template get<"velocity">(); + auto& mrho_f = macros.template get<"density">(); auto& mpor_f = macros.template get<"porosity">(); auto& particle_N_f = field.template get<"particle_N">(); @@ -160,6 +161,7 @@ public: for(uint64_t i{0u}; i < Desc::D; ++i){ force_p.at({{i}}) = 0.0; } + auto vel_p_old = (pirb_pos-pirb_pos_old) / ts; iterator::apply([&](const auto& index_f) -> void{ // ask for the d_k value here. @@ -193,19 +195,23 @@ public: uint64_t i_opp = dfi::opposite_index[i]; - saw::data dfs_added = dfs.at({i}) - dfs_old_f.at(n_ind_i).at({i_opp}); + auto u_p_e = saw::math::dot(e_i,vel_p_old); + + saw::data dfs_added = dfs.at({i})*(saw::data{1}-u_p_e.at({})) + dfs_old_f.at(n_ind_i).at({i_opp})*(saw::data{1}+u_p_e.at({})); saw::data> dfs_added_v; dfs_added_v.at({}) = dfs_added; auto ei_dfs = e_i * dfs_added_v; momentum = momentum + ei_dfs; } + // technically needs to adjust for rotation as well + auto& rho = mrho_f.at(index_f); force_p = force_p + momentum; },start,stop); auto& pirb_acc = pirb.template get<"acceleration">(); - pirb_acc = force_p / part_spheroid_group.template get<"total_mass">().at({}); + pirb_acc = - force_p / part_spheroid_group.template get<"total_mass">().at({}); for(saw::data i{0u}; i < sub_steps; ++i){ verlet_step_lambda(pi,ts); -- cgit v1.2.3