Isolating the gpu setup

4 files changed, 473 insertions, 0 deletions

diff --git a/examples/cavity_2d_gpu/.nix/derivation.nix b/examples/cavity_2d_gpu/.nix/derivation.nix
new file mode 100644
index 0000000..a98795b
--- /dev/null
+++ b/examples/cavity_2d_gpu/.nix/derivation.nix
@@ -0,0 +1,36 @@
+{ lib
+, stdenv
+, scons
+, clang-tools
+, forstio
+, pname
+, version
+}:
+
+stdenv.mkDerivation {
+	inherit pname version;
+  src = ./..;
+
+  nativeBuildInputs = [
+    scons
+    clang-tools
+  ];
+
+  buildInputs = [
+    forstio.core
+    forstio.async
+    forstio.codec
+    forstio.codec-unit
+    forstio.codec-json
+  ];
+
+  doCheck = true;
+  checkPhase = ''
+    scons test
+    ./bin/tests
+	'';
+
+	preferLocalBuild = true;
+
+  outputs = [ "out" "dev" ];
+}
diff --git a/examples/cavity_2d_gpu/SConscript b/examples/cavity_2d_gpu/SConscript
new file mode 100644
index 0000000..077fb99
--- /dev/null
+++ b/examples/cavity_2d_gpu/SConscript
@@ -0,0 +1,41 @@
+#!/bin/false
+
+import os
+import os.path
+import glob
+
+
+Import('env')
+
+dir_path = Dir('.').abspath
+
+# Environment for base library
+cavity2d_env = examples_env.Clone();
+
+examples_env.sources = sorted(glob.glob(dir_path + "/*.cpp"))
+examples_env.headers = sorted(glob.glob(dir_path + "/*.hpp"))
+
+env.sources += examples_env.sources;
+env.headers += examples_env.headers;
+
+# Cavity2D
+examples_objects = [];
+examples_env.add_source_files(examples_objects, ['cavity_2d.cpp'], shared=False);
+examples_env.cavity_2d = examples_env.Program('#bin/cavity_2d', [env.library_static, examples_objects]);
+
+# Set Alias
+env.examples = [
+	examples_env.meta_2d,
+#	examples_env.cavity_2d,
+#	examples_env.cavity_3d,
+#	examples_env.particle_ibm_2d
+	examples_env.poiseulle_2d,
+	examples_env.poiseulle_channel_2d,
+	examples_env.poiseulle_particles_channel_2d
+];
+env.Alias('examples', env.examples);
+
+if env["build_examples"]:
+	env.targets += ['examples'];
+	env.Install('$prefix/bin/', env.examples);
+#endif
diff --git a/examples/cavity_2d_gpu/SConstruct b/examples/cavity_2d_gpu/SConstruct
new file mode 100644
index 0000000..fc60882
--- /dev/null
+++ b/examples/cavity_2d_gpu/SConstruct
@@ -0,0 +1,80 @@
+#!/usr/bin/env python3
+
+import sys
+import os
+import os.path
+import glob
+import re
+
+
+if sys.version_info < (3,):
+    def isbasestring(s):
+        return isinstance(s,basestring)
+else:
+    def isbasestring(s):
+        return isinstance(s, (str,bytes))
+
+def add_kel_source_files(self, sources, filetype, lib_env=None, shared=False, target_post=""):
+
+    if isbasestring(filetype):
+        dir_path = self.Dir('.').abspath
+        filetype = sorted(glob.glob(dir_path+"/"+filetype))
+
+    for path in filetype:
+        target_name = re.sub( r'(.*?)(\.cpp|\.c\+\+)', r'\1' + target_post, path )
+        if shared:
+            target_name+='.os'
+            sources.append( self.SharedObject( target=target_name, source=path ) )
+        else:
+            target_name+='.o'
+            sources.append( self.StaticObject( target=target_name, source=path ) )
+    pass
+
+def isAbsolutePath(key, dirname, env):
+	assert os.path.isabs(dirname), "%r must have absolute path syntax" % (key,)
+
+env_vars = Variables(
+	args=ARGUMENTS
+)
+
+env_vars.Add('prefix',
+	help='Installation target location of build results and headers',
+	default='/usr/local/',
+	validator=isAbsolutePath
+)
+
+env_vars.Add('build_examples',
+	help='If examples should be built',
+	default="true"
+)
+
+env=Environment(ENV=os.environ, variables=env_vars, CPPPATH=[],
+    CPPDEFINES=['SAW_UNIX'],
+    CXXFLAGS=[
+        '-std=c++20',
+        '-g',
+        '-Wall',
+        '-Wextra'
+    ],
+    LIBS=[
+        'forstio-core',
+        'forstio-codec'
+    ]
+);
+env.__class__.add_source_files = add_kel_source_files
+env.Tool('compilation_db');
+env.cdb = env.CompilationDatabase('compile_commands.json');
+
+env.objects = [];
+env.sources = [];
+env.headers = [];
+env.targets = [];
+
+Export('env')
+SConscript('SConscript')
+
+env.Alias('cdb', env.cdb);
+env.Alias('all', [env.targets]);
+env.Default('all');
+
+env.Alias('install', '$prefix')
diff --git a/examples/cavity_2d_gpu/cavity_2d_gpu.cpp b/examples/cavity_2d_gpu/cavity_2d_gpu.cpp
new file mode 100644
index 0000000..a72cc6a
--- /dev/null
+++ b/examples/cavity_2d_gpu/cavity_2d_gpu.cpp
@@ -0,0 +1,316 @@
+#include "../c++/descriptor.hpp"
+#include "../c++/macroscopic.hpp"
+#include "../c++/lbm.hpp"
+#include "../c++/component.hpp"
+#include "../c++/collision.hpp"
+#include "../c++/boundary.hpp"
+
+#include <forstio/codec/data.hpp>
+// #include <forstio/remote/
+
+#include <iostream>
+#include <fstream>
+#include <cmath>
+
+namespace kel {
+namespace lbm {
+namespace sch {
+using namespace saw::schema;
+
+/**
+ * Basic distribution function
+ * Base type
+ * D
+ * Q
+ * Scalar factor
+ * D factor
+ * Q factor
+ */
+using T = Float32;
+using D2Q5 = Descriptor<2u,5u>;
+using D2Q9 = Descriptor<2u,9u>;
+
+template<typename Desc>
+using DfCell = Cell<T, Desc, 0u, 0u, 1u>;
+
+template<typename Desc>
+using CellInfo = Cell<UInt8, D2Q9, 1u, 0u, 0u>;
+
+/**
+ * Basic type for simulation
+ */
+template<typename Desc>
+using CellStruct = Struct<
+	Member<DfCell<Desc>, "dfs">,
+	Member<DfCell<Desc>, "dfs_old">,
+	Member<CellInfo<Desc>, "info">
+>;
+
+template<typename T, uint64_t D>
+using MacroStruct = Struct<
+	Member<FixedArray<T,D>, "velocity">,
+	Member<T, "pressure">
+>;
+
+using CavityFieldD2Q9 = CellField<D2Q9, CellStruct<D2Q9>>;
+}
+
+
+/*
+template<typename T, typename Encode>
+class df_cell_view;
+*/
+/**
+ * Minor helper for the AA-Pull Pattern, so I can use only one lattice
+ *
+ * Am I sure I want to use AA this way?
+ * Esoteric Twist technically reduces the needed memory access footprint
+ */
+/*
+template<typename Desc, size_t SN, size_t DN, size_t QN, typename Encode>
+class df_cell_view<sch::Cell<sch::T, Desc, SN, DN, QN>, Encode> {
+public:
+	using Schema = sch::Cell<sch::T,Desc,SN,DN,QN>;
+private:
+		std::array<std::decay_t<typename saw::native_data_type<sch::T>::type>*, QN> view_;
+public:
+		df_cell_view(const std::array<std::decay_t<typename saw::native_data_type<sch::T>::type>*, QN>& view):
+				view_{view}
+		{}
+};
+*/
+namespace cmpt {
+struct MovingWall {};
+}
+
+/**
+ * Full-Way moving wall Bounce back, something is not right here.
+ * Technically it should reflect properly.
+ */
+template<typename Desc>
+class component<sch::T, Desc, cmpt::MovingWall> {
+public:
+	std::array<typename saw::native_data_type<sch::T>::type, Desc::D> lid_vel;
+
+public:
+	void apply(
+		saw::data<sch::DfCell<Desc>>& dfs
+	){
+		using dfi = df_info<sch::T,Desc>;
+
+		// Technically use .copy()
+		/*
+		auto dfs_cpy = dfs;
+
+		for(uint64_t i = 0u; i < Desc::Q; ++i){
+			dfs({dfi::opposite_index.at(i)}) = dfs_cpy({i}) - 2.0 * dfi::weights[i] * 1.0 * ( lid_vel[0] * dfi::directions[i][0] + lid_vel[1] * dfi::directions[i][1]) * dfi::inv_cs2;
+		}
+		*/
+	}
+};
+}
+}
+
+constexpr size_t dim_size = 2;
+constexpr size_t dim_x = 128;
+constexpr size_t dim_y = 128;
+
+void set_geometry(saw::data<kel::lbm::sch::CavityFieldD2Q9>& latt){
+	using namespace kel::lbm;
+	/**
+	 * Set ghost
+	 */
+	iterate_over([&](const saw::data<sch::FixedArray<sch::UInt64,2u>>& index){
+		auto& cell = latt(index);
+		auto& info = cell.template get<"info">();
+
+		info({0u}).set(0u);
+
+	}, {{0u,0u}}, meta);
+
+	/**
+	 * Set wall
+	 */
+	iterate_over([&](const saw::data<sch::FixedArray<sch::UInt64,2u>>& index){
+		auto& cell = latt(index);
+		auto& info = cell.template get<"info">();
+
+		info({0u}).set(2u);
+
+	}, {{0u,0u}}, meta, {{1u,1u}});
+
+	/**
+	 * Set fluid
+	 */
+	iterate_over([&](const saw::data<sch::FixedArray<sch::UInt64,2u>>& index){
+		auto& cell = latt(index);
+		auto& info = cell.template get<"info">();
+
+		info({0u}).set(1u);
+
+	}, {{0u,0u}}, meta, {{2u,2u}});
+
+	/**
+	 * Set top lid
+	 */
+	iterate_over([&](const saw::data<sch::FixedArray<sch::UInt64,2u>>& index){
+		auto& cell = latt(index);
+		auto& info = cell.template get<"info">();
+
+		info({0u}).set(3u);
+
+	}, {{0u,1u}}, {{meta.at({0u}), 2u}}, {{2u,0u}});
+}
+
+void set_initial_conditions(saw::data<kel::lbm::sch::CavityFieldD2Q9>& latt){
+	using namespace kel::lbm;
+
+	saw::data<sch::T> rho{1.0};
+	saw::data<sch::FixedArray<sch::T,sch::D2Q9::D>> vel{{0.0,0.0}};
+	auto eq = equilibrium<sch::T,sch::D2Q9>(rho, vel);
+
+	auto meta = latt.meta();
+
+	/**
+	 * Set distribution
+	 */
+	iterate_over([&](const saw::data<sch::FixedArray<sch::UInt64,2u>>& index){
+		auto& cell = latt(index);
+		auto& dfs = cell.template get<"dfs">();
+		auto& dfs_old = cell.template get<"dfs_old">();
+
+		for(saw::data<sch::UInt64> k = 0; k < saw::data<sch::UInt64>{sch::D2Q9::Q}; ++k){
+			dfs(k) = eq.at(k);
+			dfs_old(k) = eq.at(k);
+		}
+
+	}, {{0u,0u}}, meta);
+}
+
+void lbm_step(
+	saw::data<kel::lbm::sch::CavityFieldD2Q9>& latt,
+	uint64_t time_step,
+	sycl::queue& sycl_q
+){
+	using namespace kel::lbm;
+	using dfi = df_info<sch::T,sch::D2Q9>;
+
+  /**
+	* 1. Relaxation parameter \tau
+	*/
+	component<sch::T, sch::D2Q9, cmpt::BGK> coll{0.59};
+	component<sch::T, sch::D2Q9, cmpt::BounceBack> bb;
+	component<sch::T, sch::D2Q9, cmpt::MovingWall> bb_lid;
+	bb_lid.lid_vel = {0.1,0.0};
+
+	auto meta = latt.meta();
+
+	/**
+	 * Collision
+	 */
+	iterate_over([&](const saw::data<sch::FixedArray<sch::UInt64,2u>>& index){
+		auto& cell = latt(index);
+		auto& info = cell.template get<"info">();
+
+		auto& dfs = cell.template get<"dfs">();
+		auto& dfs_old = cell.template get<"dfs_old">();
+
+		switch(info({0u}).get()){
+			case 1u: {
+				coll.apply(latt, index, time_step);
+				break;
+			}
+			case 2u: {
+				bb.apply(latt, index, time_step);
+				break;
+			}
+			default:
+				break;
+		}
+
+	}, {{0u,0u}}, meta);
+
+	// Stream
+	for(uint64_t i = 1u; (i+1u) < latt.template get_dim_size<0>().get(); ++i){
+		for(uint64_t j = 1u; (j+1u) < latt.template get_dim_size<1>().get(); ++j){
+			auto& cell = latt({{i,j}});
+			auto& df_new = even_step ? cell.template get<"dfs">() : cell.template get<"dfs_old">();
+			auto& info_new = cell.template get<"info">();
+
+			if(info_new({0u}).get() > 0u && info_new({0u}).get() != 3u){
+				for(uint64_t k = 0u; k < sch::D2Q9::Q; ++k){
+					auto dir = dfi::directions[dfi::opposite_index[k]];
+					auto& cell_dir_old = latt({{i+dir[0],j+dir[1]}});
+
+					auto& df_old = even_step ? cell_dir_old.template get<"dfs_old">() : cell_dir_old.template get<"dfs">();
+					auto& info_old = cell_dir_old.template get<"info">();
+
+					if( info_old({0}).get() == 3u ){
+						auto& df_old_loc = even_step ? latt({{i,j}}).template get<"dfs_old">() : latt({{i,j}}).template get<"dfs">();
+						df_new({k}) = df_old_loc({dfi::opposite_index.at(k)}) - 2.0 * dfi::inv_cs2 * dfi::weights.at(k) * 1.0 * ( bb_lid.lid_vel[0] * dir[0] + bb_lid.lid_vel[1] * dir[1]);
+						// dfs({dfi::opposite_index.at(i)}) = dfs_cpy({i}) - 2.0 * dfi::weights[i] * 1.0 * ( lid_vel[0] * dfi::directions[i][0] + lid_vel[1] * dfi::directions[i][1]) * dfi::inv_cs2;
+					} else {
+						df_new({k}) = df_old({k});
+					}
+				}
+			}
+		}
+	}
+}
+
+int main(){
+	using namespace kel::lbm;
+
+	saw::data<sch::FixedArray<sch::UInt64,sch::D2Q9::D>> dim{{dim_x, dim_y}};
+
+	saw::data<sch::CavityFieldD2Q9, saw::encode::Native> lattice{dim};
+
+	converter<sch::T> conv{
+		{0.1},
+		{0.1}
+	};
+
+	print_lbm_meta<sch::T, sch::D2Q9>(conv, {1e-3});
+
+	auto eo_lbm_dir = output_directory();
+	if(eo_lbm_dir.is_error()){
+		return -1;
+	}
+	auto& lbm_dir = eo_lbm_dir.get_value();
+	auto out_dir = lbm_dir / "cavity_gpu_2d";
+
+	/**
+	 * Set meta information describing what this cell is
+	 */
+	set_geometry(lattice);
+
+	/**
+	 *
+	 */
+	set_initial_conditions(lattice);
+
+	sycl::queue sycl_q{sycl::default_selector_v, sycl::property::queue::in_order{}};
+
+	/**
+	 * Timeloop
+	 */
+
+	uint64_t lattice_steps = 512000u;
+	bool even_step = true;
+
+	uint64_t print_every = 256u;
+	uint64_t file_no = 0u;
+
+	saw::data<sch::Array<sch::MacroStruct<sch::T,sch::D2Q9::D>,sch::D2Q9::D>> macros{dim};
+
+	for(uint64_t i = 0u; i < 256u; ++i){
+		{
+			std::string vtk_f_name{"tmp/poiseulle_2d_"};
+			vtk_f_name += std::to_string(i) + ".vtk";
+			write_vtk_file(vtk_f_name, macros);
+		}
+
+		lbm_step(lattice, i, sycl_q);
+	}
+	return 0;
+}