diff options
| -rw-r--r-- | default.nix | 11 | ||||
| -rw-r--r-- | examples/poiseulle_moving_particle_2d_psm_gpu/sim.cpp | 38 | ||||
| -rw-r--r-- | examples/poiseulle_particles_2d_bgk_gpu/sim.cpp | 14 | ||||
| -rw-r--r-- | examples/poiseulle_particles_2d_psm_gpu/sim.cpp | 37 | ||||
| -rw-r--r-- | lib/core/c++/particle/particle_opa.hpp | 3 | ||||
| -rwxr-xr-x | scripts/python/csv_l2_norm.py | 26 |
6 files changed, 108 insertions, 21 deletions
diff --git a/default.nix b/default.nix index dc75b57..02f9213 100644 --- a/default.nix +++ b/default.nix @@ -151,6 +151,17 @@ in rec { inherit pname version stdenv forstio adaptive-cpp; inherit kel; }; + + poiseulle_particles_2d_all_gpu = pkgs.symlinkJoin { + name = "poiseulle_particles_2d_all_gpu"; + paths = [ + examples.poiseulle_particles_2d_bgk_gpu + examples.poiseulle_particles_2d_psm_gpu + examples.poiseulle_particles_2d_hlbm_gpu + examples.poiseulle_particles_2d_fplbm_gpu + examples.poiseulle_particles_2d_ibm_gpu + ]; + }; poiseulle_moving_particle_2d_psm_gpu = pkgs.callPackage ./examples/poiseulle_moving_particle_2d_psm_gpu/.nix/derivation.nix { inherit pname version stdenv forstio adaptive-cpp; diff --git a/examples/poiseulle_moving_particle_2d_psm_gpu/sim.cpp b/examples/poiseulle_moving_particle_2d_psm_gpu/sim.cpp index 0c10d38..cae8a4d 100644 --- a/examples/poiseulle_moving_particle_2d_psm_gpu/sim.cpp +++ b/examples/poiseulle_moving_particle_2d_psm_gpu/sim.cpp @@ -1,6 +1,7 @@ #include <kel/lbm/sycl/lbm.hpp> #include <kel/lbm/lbm.hpp> #include <kel/lbm/particle.hpp> +#include <kel/lbm/particle_opa.hpp> #include <forstio/io/io.hpp> #include <forstio/remote/filesystem/easy.hpp> @@ -157,25 +158,6 @@ saw::error_or<void> setup_initial_conditions( df_f.get_dims(), {{1u,1u}} ); - - iterator<Desc::D>::apply( - [&](auto& index){ - saw::data<sch::Vector<T,Desc::D>> middle, ind_vec; - middle.at({{0u}}) = dim_x * 0.5; - middle.at({{1u}}) = dim_y * 0.5; - - ind_vec.at({{0u}}) = index.at({{0u}}).template cast_to<T>(); - ind_vec.at({{1u}}) = index.at({{1u}}).template cast_to<T>(); - - auto dist = middle - ind_vec; - auto dist_2 = saw::math::dot(dist,dist); - if(dist_2.at({}).get() < dim_y*dim_y*0.01){ - porous_f.at(index).at({}) = 0.0; - } - }, - {},// 0-index - df_f.get_dims() - ); return saw::make_void(); } @@ -190,6 +172,8 @@ saw::error_or<void> step( auto& q = dev.get_handle(); auto& info_f = fields.template get<"info">(); auto& porous_f = macros.template get<"porosity">(); + sch::data<sch::Scalar<T>> eps; + eps.at({}).set(1.5f); q.submit([&](acpp::sycl::handler& h){ component<T,Desc,cmpt::PSM,encode::Sycl<saw::encode::Native>> collision{0.8}; @@ -199,7 +183,7 @@ saw::error_or<void> step( component<T,Desc,cmpt::ZouHeHorizontal<true>,encode::Sycl<saw::encode::Native>> flow_in{1.0}; component<T,Desc,cmpt::ZouHeHorizontal<false>,encode::Sycl<saw::encode::Native>> flow_out{1.0}; - component<T,Desc,cmpt::OneParticleAt, encode::Sycl<saw::encode::Native>> opa{{},{}}; + component<T,Desc,cmpt::OneParticleAt, encode::Sycl<saw::encode::Native>> opa{{},{},eps}; h.parallel_for(acpp::sycl::range<Desc::D>{dim_x,dim_y}, [=](acpp::sycl::id<Desc::D> idx){ saw::data<sch::FixedArray<sch::UInt64,Desc::D>> index; @@ -371,6 +355,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){ } } */ + /* if(i.get() % 32u == 0u){ { auto eov = dev.copy_to_host(lbm_sycl_macro_data,*lbm_macro_data_ptr); @@ -385,6 +370,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){ } } } + */ // Stream sycl_q.submit([&](acpp::sycl::handler& h){ component<T,Desc,cmpt::Stream,encode::Sycl<saw::encode::Native>> stream; @@ -424,6 +410,18 @@ saw::error_or<void> lbm_main(int argc, char** argv){ return eov; } } + { + auto eov = dev.copy_to_host(lbm_sycl_macro_data,*lbm_macro_data_ptr); + if(eov.is_error()){ + return eov; + } + } + { + auto eov = write_csv_file(out_dir,"m",time_steps.get(), *lbm_macro_data_ptr); + if(eov.is_error()){ + return eov; + } + } sycl_q.wait(); return saw::make_void(); diff --git a/examples/poiseulle_particles_2d_bgk_gpu/sim.cpp b/examples/poiseulle_particles_2d_bgk_gpu/sim.cpp index 9cc77ae..73964b7 100644 --- a/examples/poiseulle_particles_2d_bgk_gpu/sim.cpp +++ b/examples/poiseulle_particles_2d_bgk_gpu/sim.cpp @@ -366,6 +366,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){ } } */ + /* if(i.get() % 32u == 0u){ { auto eov = dev.copy_to_host(lbm_sycl_macro_data,*lbm_macro_data_ptr); @@ -380,6 +381,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){ } } } + */ // Stream sycl_q.submit([&](acpp::sycl::handler& h){ component<T,Desc,cmpt::Stream,encode::Sycl<saw::encode::Native>> stream; @@ -419,6 +421,18 @@ saw::error_or<void> lbm_main(int argc, char** argv){ return eov; } } + { + auto eov = dev.copy_to_host(lbm_sycl_macro_data,*lbm_macro_data_ptr); + if(eov.is_error()){ + return eov; + } + } + { + auto eov = write_csv_file(out_dir,"m",time_steps.get(), *lbm_macro_data_ptr); + if(eov.is_error()){ + return eov; + } + } sycl_q.wait(); return saw::make_void(); diff --git a/examples/poiseulle_particles_2d_psm_gpu/sim.cpp b/examples/poiseulle_particles_2d_psm_gpu/sim.cpp index 0337158..e91686c 100644 --- a/examples/poiseulle_particles_2d_psm_gpu/sim.cpp +++ b/examples/poiseulle_particles_2d_psm_gpu/sim.cpp @@ -1,6 +1,7 @@ #include <kel/lbm/sycl/lbm.hpp> #include <kel/lbm/lbm.hpp> #include <kel/lbm/particle.hpp> +#include <kel/lbm/particle/particle_opa.hpp> #include <forstio/io/io.hpp> #include <forstio/remote/filesystem/easy.hpp> @@ -157,7 +158,27 @@ saw::error_or<void> setup_initial_conditions( df_f.get_dims(), {{1u,1u}} ); + + saw::data<sch::Scalar<T>> eps; + eps.at({}).set(1.5f); + saw::data<sch::Scalar<T>> rad; + rad.at({}).set(dim_y*0.1); + saw::data<sch::Vector<T,Desc::D>> p_pos; + { + p_pos.at({{0u}}) = dim_x * 0.25; + p_pos.at({{1u}}) = dim_y * 0.5; + } + + component<T,Desc,cmpt::OneParticleAt> opa{p_pos,rad,eps}; + iterator<Desc::D>::apply( + [&](auto& index){ + opa.apply(macros,index,{}); + }, + {},// 0-index + df_f.get_dims() + ); + /* iterator<Desc::D>::apply( [&](auto& index){ saw::data<sch::Vector<T,Desc::D>> middle, ind_vec; @@ -176,6 +197,7 @@ saw::error_or<void> setup_initial_conditions( {},// 0-index df_f.get_dims() ); + */ return saw::make_void(); } @@ -203,6 +225,18 @@ saw::error_or<void> step( component<T,Desc,cmpt::Equilibrium,encode::Sycl<saw::encode::Native>> equi{rho_b,vel_b}; + + saw::data<sch::Scalar<T>> eps; + eps.at({}).set(1.5f); + saw::data<sch::Scalar<T>> rad; + rad.at({}).set(dim_y*0.1); + saw::data<sch::Vector<T,Desc::D>> p_pos; + { + p_pos.at({{0u}}) = dim_x * 0.25; + p_pos.at({{1u}}) = dim_y * 0.5; + } + component<T,Desc,cmpt::OneParticleAt,encode::Sycl<saw::encode::Native>> opa{p_pos,rad,eps}; + component<T,Desc,cmpt::ZouHeHorizontal<true>,encode::Sycl<saw::encode::Native>> flow_in{ [&](){ uint64_t target_t_i = 64u; @@ -230,6 +264,7 @@ saw::error_or<void> step( bb.apply(fields,index,t_i); break; case 2u: + opa.apply(macros,index,t_i); collision.apply(fields,macros,index,t_i); break; case 3u: @@ -385,6 +420,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){ } } */ + /* if(i.get() % 32u == 0u){ { auto eov = dev.copy_to_host(lbm_sycl_macro_data,*lbm_macro_data_ptr); @@ -399,6 +435,7 @@ saw::error_or<void> lbm_main(int argc, char** argv){ } } } + */ // Stream sycl_q.submit([&](acpp::sycl::handler& h){ component<T,Desc,cmpt::Stream,encode::Sycl<saw::encode::Native>> stream; diff --git a/lib/core/c++/particle/particle_opa.hpp b/lib/core/c++/particle/particle_opa.hpp index e6396d6..bfd2e3c 100644 --- a/lib/core/c++/particle/particle_opa.hpp +++ b/lib/core/c++/particle/particle_opa.hpp @@ -35,12 +35,13 @@ public: auto& porous = porous_f.at(index); - auto pos_ind = saw::math::vectorize_data(index); auto diff = pos_ind - pos_; // Write out value + auto diff = pos_ind.template cast_to<T>() - pos_; + auto diff_dot = saw::math::dot(diff,diff); porous = particle_porosity<>::calculate(diff, rad_, eps_); } }; diff --git a/scripts/python/csv_l2_norm.py b/scripts/python/csv_l2_norm.py new file mode 100755 index 0000000..88f4817 --- /dev/null +++ b/scripts/python/csv_l2_norm.py @@ -0,0 +1,26 @@ +#!/usr/bin/env python3 + +import argparse +import numpy as np + + +def main(): + parser = argparse.ArgumentParser(description="Compute the L2 Norm"); + parser.add_argument("one"); + parser.add_argument("two"); + args = parser.parse_args(); + + a = np.loadtxt(args.one, delimiter=","); + b = np.loadtxt(args.two, delimiter=","); + + if a.shape != b.shape: + raise ValueError(f"Shape mismatch: {a.shape} vs {b.shape}"); + #endif + + print(np.linalg.norm(a-b)); + +#enddef + +if __name__ == "__main__": + main(); +#enddef |