Streamlining and thoughts about changes in particle dynamicsdev

4 files changed, 87 insertions, 9 deletions

diff --git a/c++/descriptor.hpp b/c++/descriptor.hpp
index 0c6707a..d04d217 100644
--- a/c++/descriptor.hpp
+++ b/c++/descriptor.hpp
@@ -1,6 +1,7 @@
 #pragma once
 
 #include <forstio/codec/data.hpp>
+#include <forstio/codec/data_math.hpp>
 #include <forstio/codec/schema_factory.hpp>
 
 namespace kel {
diff --git a/c++/macroscopic.hpp b/c++/macroscopic.hpp
index 8a25248..51368e9 100644
--- a/c++/macroscopic.hpp
+++ b/c++/macroscopic.hpp
@@ -59,5 +59,34 @@ void compute_rho_u (
 		vel().at({i}) = vel().at({i}) / rho();
 	}
 }
+
+/**
+ * Calculate the macroscopic variables rho and u in Lattice Units.
+ */
+template<typename T, typename Desc>
+void compute_rho_u (
+		const saw::data<sch::Cell<T, Desc, 0, 0, 1>>& dfs,
+		saw::ref<saw::data<T>> rho,
+		saw::ref<saw::data<sch::Vector<T,Desc::D>>> vel
+	)
+{
+	using dfi = df_info<T, Desc>;
+
+	rho().set(0);
+	for(uint64_t i = 0; i < Desc::D; ++i){
+		vel().at({{i}}).set(0);
+	}
+
+	for(size_t j = 0; j < Desc::Q; ++j){
+		rho() = rho() + dfs(j);
+		for(size_t i = 0; i < Desc::D; ++i){
+			vel().at({{i}}) = vel().at({{i}}) + saw::data<T>{static_cast<typename saw::native_data_type<T>::type>(dfi::directions[j][i])} * dfs(j);
+		}
+	}
+
+	for(size_t i = 0; i < Desc::D; ++i){
+		vel().at({{i}}) = vel().at({{i}}) / rho();
+	}
+}
 }
 }
diff --git a/c++/util.hpp b/c++/util.hpp
index 441ffc8..c5b6a05 100644
--- a/c++/util.hpp
+++ b/c++/util.hpp
@@ -4,6 +4,8 @@
 #include <forstio/codec/data.hpp>
 
 #include <iostream>
+#include <iomanip>
+#include <sys/ioctl.h>
 
 namespace kel {
 namespace lbm {
@@ -46,18 +48,37 @@ public:
 };
 */
 
-void print_progress_bar(uint32_t progress, uint32_t progress_target){
-	std::cout
-		<<"\r"
-		<<"Progress: "
-		<<((100 * progress) / progress_target)
-		<<"% [";
+void print_progress_bar(const uint32_t progress, const uint32_t progress_target){
+	std::cout<<"\r";
+	//	<<"Progress: "
+	//	<<((100 * progress) / progress_target)
+	//	<<"% [";
 
-	uint64_t i{0u};
+	uint64_t perc_prog = (100ul*progress) / progress_target;
+	constexpr uint64_t max_perc_progress = 100u;
+	perc_prog = std::min(perc_prog, max_perc_progress);
 
-	uint64_t max_display = 64u;
+	std::string_view prog_str{"Progress: "};
+	// Progress string + (100 width and perc char) + ([] brackets) + space + pointer
+	uint64_t i{prog_str.size() + 4u + 2u + 1u + 1u};
+
+	std::cout<<prog_str;
+	std::cout<<std::setw(3u)<<perc_prog<<"%";
+	std::cout<<" ";
+	std::cout<<"[";
+
+	uint64_t max_display = []() -> uint64_t{
+		struct winsize w;
+		if(ioctl(STDOUT_FILENO, TIOCGWINSZ,&w) == -1){
+			// Standardized Terminal size
+			return 80u;
+		}
+		return w.ws_col;
+	}();
+	max_display = std::max(max_display,i) - i;
 	uint64_t progress_display = (max_display * progress) / progress_target;
-	for(; i < progress_display; ++i){
+
+	for(i = 0u; i < progress_display; ++i){
 		std::cout<<"#";
 	}
 	for(; i < max_display; ++i){
@@ -65,6 +86,7 @@ void print_progress_bar(uint32_t progress, uint32_t progress_target){
 	}
 
 	std::cout<<"]";
+
 	std::cout<<std::flush;
 }
 }
diff --git a/c++/write_vtk.hpp b/c++/write_vtk.hpp
index c9ecc99..0647db5 100644
--- a/c++/write_vtk.hpp
+++ b/c++/write_vtk.hpp
@@ -3,6 +3,7 @@
 #include <forstio/error.hpp>
 
 #include <forstio/codec/data.hpp>
+#include <forstio/codec/data_math.hpp>
 
 #include "descriptor.hpp"
 
@@ -55,6 +56,31 @@ struct lbm_vtk_writer<sch::FixedArray<T,D>> {
 	}
 };
 
+template<typename T, uint64_t D>
+struct lbm_vtk_writer<sch::Vector<T,D>> {
+	static saw::error_or<void> apply_header(std::ostream& vtk_file, std::string_view name){
+		vtk_file<<"VECTORS "<<name<<" float\n";
+		return saw::make_void();
+	}
+	static saw::error_or<void> apply(std::ostream& vtk_file, const saw::data<sch::Vector<T,D>>& field){
+		static_assert(D > 0, "Non-dimensionality is bad for velocity.");
+		static_assert(D <= 3, "4th dimension as well. Mostly due to vtk.");
+
+		// vtk_file<<"VECTORS "<<name<<" float\n";
+		for(uint64_t i = 0u; i < D; ++i){
+			if(i > 0){
+				vtk_file<<" ";
+			}
+			vtk_file<<field.at({{i}}).get();
+		}
+		for(uint64_t i = D; i < 3; ++i){
+			vtk_file<<" 0";
+		}
+		vtk_file<<"\n";
+		return saw::make_void();
+	}
+};
+
 template<uint64_t Dim, typename... StructT, saw::string_literal... StructN>
 struct lbm_vtk_writer<sch::Array<sch::Struct<sch::Member<StructT,StructN>...> , Dim>> {