#include #include #include #include #include #include namespace kel { namespace lbm { constexpr uint64_t dim_x = 1024u; constexpr uint64_t dim_y = 512u; namespace sch { using namespace saw::schema; using InfoChunk = Chunk; template using DfChunk = Chunk, 1u, dim_x, dim_y>; template using ChunkStruct = Struct< Member, Member, "dfs">, Member, "dfs_old"> >; template using VelChunk = Chunk, 0u, dim_x, dim_y>; template using RhoChunk = Chunk, 0u, dim_x, dim_y>; template using MacroStruct = Struct< Member, "velocity">, Member, "density"> >; // template // using Particles = Particle< } template saw::error_or setup_initial_conditions( saw::data>& fields, saw::data>& macros ){ auto& info_f = fields.template get<"info">(); // Set everything as walls iterator::apply( [&](auto& index){ info_f.at(index).set(1u); }, {}, info_f.get_dims(), {} ); // Fluid iterator::apply( [&](auto& index){ info_f.at(index).set(2u); }, {}, info_f.get_dims(), {{1u,1u}} ); // Inflow iterator::apply( [&](auto& index){ info_f.at(index).set(3u); }, {{0u,0u}}, {{1u, dim_y}}, {{0u,1u}} ); // Outflow iterator::apply( [&](auto& index){ info_f.at(index).set(4u); }, {{dim_x-1u,0u}}, {{dim_x, dim_y}}, {{0u,1u}} ); // auto& df_f = fields.template get<"dfs_old">(); auto& rho_f = macros.template get<"density">(); auto& vel_f = macros.template get<"velocity">(); iterator::apply( [&](auto& index){ auto& df = df_f.at(index); auto& rho = rho_f.at(index); rho.at({}) = {1}; auto& vel = vel_f.at(index); auto eq = equilibrium(rho,vel); df = eq; }, {},// 0-index df_f.get_dims() ); iterator::apply( [&](auto& index){ auto& df = df_f.at(index); auto& rho = rho_f.at(index); rho.at({}) = {1}; auto& vel = vel_f.at(index); vel.at({{0u}}) = 0.1; auto eq = equilibrium(rho,vel); df = eq; }, {},// 0-index df_f.get_dims(), {{1u,1u}} ); return saw::make_void(); } template saw::error_or step( saw::data>,encode::Sycl>& fields, saw::data>,encode::Sycl>& macros, saw::data t_i, device& dev ){ auto& q = dev.get_handle(); auto& info_f = fields.template get<"info">(); // auto coll_ev = q.submit([&](acpp::sycl::handler& h){ // Need nicer things to handle the flow. I see improvement here component> collision{0.6}; component> bb; saw::data> rho_b; rho_b.at({}) = 1.0; saw::data> vel_b; vel_b.at({{0u}}) = 0.1; component> equi{rho_b,vel_b}; component,encode::Sycl> flow_in{ [&](){ uint64_t target_t_i = 256u; if(t_i.get() < target_t_i){ return 1.0 + (0.001 / target_t_i) * t_i.get(); } return 1.001; }() }; component,encode::Sycl> flow_out{1.0}; h.parallel_for(acpp::sycl::range{dim_x,dim_y}, [=](acpp::sycl::id idx){ saw::data> index; for(uint64_t i = 0u; i < Desc::D; ++i){ index.at({{i}}).set(idx[i]); } auto info = info_f.at(index); switch(info.get()){ case 0u: break; case 1u: bb.apply(fields,index,t_i); break; case 2u: collision.apply(fields,macros,index,t_i); break; case 3u: equi.apply(fields,index,t_i); // flow_in.apply(fields,index,t_i); collision.apply(fields,macros,index,t_i); break; case 4u: equi.apply(fields,index,t_i); // flow_out.apply(fields,index,t_i); collision.apply(fields,macros,index,t_i); break; default: break; } }); }).wait(); q.submit([&](acpp::sycl::handler& h){ component> stream; h.parallel_for(acpp::sycl::range{dim_x,dim_y}, [=](acpp::sycl::id idx){ saw::data> index; for(uint64_t i = 0u; i < Desc::D; ++i){ index.at({{i}}).set(idx[i]); } auto info = info_f.at(index); if(info.get() > 0u){ stream.apply(fields,index,t_i); } }); }).wait(); // Step /* q.submit([&](acpp::sycl::handler& h){ // h.depends_on(collision_ev); }).wait(); */ return saw::make_void(); } } } template saw::error_or lbm_main(int argc, char** argv){ using namespace kel::lbm; using dfi = df_info; auto eo_lbm_dir = output_directory(); if(eo_lbm_dir.is_error()){ return std::move(eo_lbm_dir.get_error()); } auto& lbm_dir = eo_lbm_dir.get_value(); auto out_dir = lbm_dir / "poiseulle_particles_2d_gpu"; converter conv { // delta_x {{1.0}}, // delta_t {{1.0}} }; // saw::data> meta{{dim_x,dim_y}}; auto lbm_data_ptr = saw::heap>>(); auto lbm_macro_data_ptr = saw::heap>>(); auto eo_aio = saw::setup_async_io(); if(eo_aio.is_error()){ return std::move(eo_aio.get_error()); } auto& aio = eo_aio.get_value(); saw::wait_scope wait{aio.event_loop}; bool krun = true; bool print_status = false; aio.event_port.on_signal(saw::Signal::Terminate).then([&](){ krun = false; }).detach(); aio.event_port.on_signal(saw::Signal::User1).then([&](){ print_status = true; }).detach(); device dev; auto& sycl_q = dev.get_handle(); sycl_q.wait(); { auto eov = setup_initial_conditions(*lbm_data_ptr,*lbm_macro_data_ptr); if(eov.is_error()){ return eov; } } saw::data, encode::Sycl> lbm_sycl_data{sycl_q}; saw::data, encode::Sycl> lbm_sycl_macro_data{sycl_q}; sycl_q.wait(); auto lsd_view = make_chunk_struct_view(lbm_sycl_data); auto lsdm_view = make_chunk_struct_view(lbm_sycl_macro_data); { auto eov = dev.copy_to_device(*lbm_data_ptr,lbm_sycl_data); if(eov.is_error()){ return eov; } } { auto eov = dev.copy_to_device(*lbm_macro_data_ptr,lbm_sycl_macro_data); if(eov.is_error()){ return eov; } } sycl_q.wait(); saw::data time_steps{4096ul}; for(saw::data i{0u}; i < time_steps and krun; ++i){ { std::string file_name = "tmp/t_"; file_name += std::to_string(i.get()); file_name += ".vtk"; { auto eov = write_vtk_file(file_name, *lbm_macro_data_ptr); if(eov.is_error()){ return eov; } } { auto eov = saw::easy::encode_and_write_file,saw::encode::KelSimple>(file_name, *lbm_macro_data_ptr); if(eov.is_error()){ return eov; } } } { std::string file_name = "tmp/df_"; file_name += std::to_string(i.get()); file_name += ".vtk"; auto eov = write_vtk_file(file_name, *lbm_data_ptr); if(eov.is_error()){ return eov; } } { auto eov = step(lsd_view,lsdm_view,i,dev); if(eov.is_error()){ return eov; } } { auto eov = dev.copy_to_host(lbm_sycl_macro_data,*lbm_macro_data_ptr); if(eov.is_error()){ return eov; } } { auto eov = dev.copy_to_host(lbm_sycl_data,*lbm_data_ptr); if(eov.is_error()){ return eov; } } wait.poll(); if(print_status){ std::cout<<"Status: "<().get() * 100 / time_steps.get())<<"%"<::apply( [&](auto& index){ std::cout<().at(index).template cast_to().get()<<"\n"; }, {{0u,0u}}, {{dim_x, dim_y}} ); */ sycl_q.wait(); return saw::make_void(); } using FloatT = kel::lbm::sch::Float32; int main(int argc, char** argv){ auto eov = lbm_main(argc, argv); if(eov.is_error()){ auto& err = eov.get_error(); std::cerr<<"[Error] "< 0u){ std::cerr<<" - "<