summaryrefslogtreecommitdiff
path: root/lib/core/c++
diff options
context:
space:
mode:
Diffstat (limited to 'lib/core/c++')
-rw-r--r--lib/core/c++/abstract/templates.hpp3
-rw-r--r--lib/core/c++/collision.hpp54
-rw-r--r--lib/core/c++/environment.hpp54
-rw-r--r--lib/core/c++/particle/particle.hpp5
4 files changed, 101 insertions, 15 deletions
diff --git a/lib/core/c++/abstract/templates.hpp b/lib/core/c++/abstract/templates.hpp
index d4d4ace..f675a99 100644
--- a/lib/core/c++/abstract/templates.hpp
+++ b/lib/core/c++/abstract/templates.hpp
@@ -1,4 +1,7 @@
#pragma once
namespace kel {
+namespace lbm {
+
+}
}
diff --git a/lib/core/c++/collision.hpp b/lib/core/c++/collision.hpp
index ef6bfe2..df45bbe 100644
--- a/lib/core/c++/collision.hpp
+++ b/lib/core/c++/collision.hpp
@@ -95,10 +95,21 @@ class component<T, Descriptor, cmpt::BGKGuo, Encode> {
private:
typename saw::native_data_type<T>::type relaxation_;
saw::data<T> frequency_;
+ saw::data<sch::Vector<T,Descriptor::D>> global_force_;
public:
component(typename saw::native_data_type<T>::type relaxation__):
relaxation_{relaxation__},
- frequency_{typename saw::native_data_type<T>::type(1) / relaxation_}
+ frequency_{typename saw::native_data_type<T>::type(1) / relaxation_},
+ global_force_{}
+ {}
+
+ component(
+ typename saw::native_data_type<T>::type relaxation__,
+ saw::data<sch::Vector<T,Descriptor::D>>& global_force__
+ ):
+ relaxation_{relaxation__},
+ frequency_{typename saw::native_data_type<T>::type(1) / relaxation_},
+ global_force_{global_force__}
{}
using Component = cmpt::BGKGuo;
@@ -115,39 +126,56 @@ public:
static constexpr saw::string_literal after = "";
static constexpr saw::string_literal before = "streaming";
- template<typename CellFieldSchema>
- void apply(saw::data<CellFieldSchema, Encode>& field, saw::data<sch::FixedArray<sch::UInt64, Descriptor::D>> index, saw::data<sch::UInt64> time_step){
+ template<typename CellFieldSchema, typename MacroFieldSchema>
+ void apply(const saw::data<CellFieldSchema, Encode>& field, const saw::data<MacroFieldSchema,Encode>& macros, saw::data<sch::FixedArray<sch::UInt64,Descriptor::D>> index, saw::data<sch::UInt64> time_step) const {
+ // void apply(saw::data<CellFieldSchema, Encode>& field, saw::data<sch::FixedArray<sch::UInt64, Descriptor::D>> index, saw::data<sch::UInt64> time_step){
bool is_even = ((time_step.get() % 2) == 0);
auto& dfs_old_f = (is_even) ? field.template get<"dfs_old">() : field.template get<"dfs">();
auto& dfs = dfs_old_f.at(index);
- saw::data<sch::Scalar<T>> rho;
- saw::data<sch::Vector<T,Descriptor::D>> vel;
+ auto& rho_f = macros.template get<"density">();
+ auto& vel_f = macros.template get<"velocity">();
+
+ saw::data<sch::Scalar<T>>& rho = rho_f.at(index);
+ saw::data<sch::Vector<T,Descriptor::D>>& vel = vel_f.at(index);
+
compute_rho_u<T,Descriptor>(dfs_old_f.at(index),rho,vel);
auto eq = equilibrium<T,Descriptor>(rho,vel);
using dfi = df_info<T,Descriptor>;
- auto& force_f = field.template get<"force">();
+ auto& force_f = macros.template get<"force">();
auto& force = force_f.at(index);
+
+ saw::data<sch::Scalar<T>> dfi_inv_cs2;
+ dfi_inv_cs2.at({}).set(dfi::inv_cs2);
for(uint64_t i = 0u; i < Descriptor::Q; ++i){
// saw::data<T> ci_min_u{0};
- saw::data<T> ci_dot_u{0};
+ saw::data<sch::Vector<T,Descriptor::D>> ci;
for(uint64_t d = 0u; d < Descriptor::D; ++d){
- ci_dot_u = ci_dot_u + vel.at({{d}}) * saw::data<T>{static_cast<typename saw::native_data_type<T>::type>(dfi::directions[i][d])};
+ ci.at({{d}}).set(static_cast<typename saw::native_data_type<T>::type>(dfi::directions[i][d]));
}
- saw::data<T> F_i{0};// = saw::data<T>{dfi::weights[i]} * (ci_dot_F * dfi::inv_cs2 + ci_dot_u * ci_dot_F * (dfi::inv_cs2 * dfi::inv_cs2));
+ auto ci_dot_u = saw::math::dot(ci,vel);
+
+ // saw::data<sch::Vector<T,Descriptor::D>> F_i;
+ // F_i = f * (c_i - u * ics2 + <c_i,u> * c_i * ics2 * ics2) * w_i;
+ saw::data<sch::Scalar<T>> w;
+ w.at({}).set(dfi::weights[i]);
+ auto F_i_d = saw::math::dot((force+global_force_) * w, (ci - vel * dfi_inv_cs2 + ci * ci_dot_u * dfi_inv_cs2 * dfi_inv_cs2 ));
+ /*
+ saw::data<sch::Scalar<T>> F_i_sum;
for(uint64_t d = 0u; d < Descriptor::D; ++d){
- F_i = F_i +
- saw::data<T>{dfi::weights[i]} * ((saw::data<T>{static_cast<typename saw::native_data_type<T>::type>(dfi::directions[i][d])}
- - vel.at({d}) ) * dfi::inv_cs2 + ci_dot_u * saw::data<T>{static_cast<typename saw::native_data_type<T>::type>(dfi::directions[i][d])} * dfi::inv_cs2 * dfi::inv_cs2 ) * force.at({d});
+ saw::data<sch::Scalar<T>> F_i_d;
+ F_i_d.at({}) = F_i.at({{d}});
+ F_i_sum = F_i_sum + F_i_d;
}
+ */
- dfs.at({i}) = dfs.at({i}) + frequency_ * (eq.at(i) - dfs.at({i}) ) + F_i * (saw::data<T>{1} - saw::data<T>{0.5f} * frequency_);
+ dfs.at({i}) = dfs.at({i}) + frequency_ * (eq.at(i) - dfs.at({i}) ) + F_i_d.at({}) * (saw::data<T>{1} - saw::data<T>{0.5f} * frequency_);
}
}
};
diff --git a/lib/core/c++/environment.hpp b/lib/core/c++/environment.hpp
index 9e587d0..27d8f3a 100644
--- a/lib/core/c++/environment.hpp
+++ b/lib/core/c++/environment.hpp
@@ -7,18 +7,70 @@
namespace kel {
namespace lbm {
+struct environment {
+ std::filesystem::path lbm_dir;
+ std::filesystem::path data_dir;
+};
+
+saw::error_or<environment> lbm_directory(){
+ namespace fs = std::filesystem;
+
+ const char* home_dir = std::getenv("HOME");
+ if(not home_dir){
+ return saw::make_error<saw::err::not_found>("Couldn't find home dir");
+ }
+
+ environment env;
+
+ env.lbm_dir = std::filesystem::path{home_dir} / ".lbm";
+
+ {
+ env.data_dir = env.lbm_dir / "data";
+ // LBM Data Location
+ auto lbm_dir_config = env.lbm_dir / "data_dir_location.txt";
+ if( fs::exists(lbm_dir_config) && fs::is_regular_file(lbm_dir_config) ){
+ std::ifstream file{lbm_dir_config};
+
+ if(file.is_open()){
+ std::stringstream buffer;
+ buffer <<file.rdbuf();
+ env.data_dir = fs::path{buffer.str()};
+ }
+ }
+ }
+
+ return env;
+}
+
/**
* Returns the default output directory.
* Located outside the project dir because dispatching build jobs with output data in the git directory
* also copies simulated data which takes a long time.
*/
saw::error_or<std::filesystem::path> output_directory(){
+ namespace fs = std::filesystem;
+
const char* home_dir = std::getenv("HOME");
if(not home_dir){
return saw::make_error<saw::err::not_found>("Couldn't find home dir");
}
- return std::filesystem::path{home_dir} / ".lbm";
+ auto lbm_dir = std::filesystem::path{home_dir} / ".lbm";
+
+ {
+ auto lbm_dir_config = lbm_dir / "data_dir_location.txt";
+ if( fs::exists(lbm_dir_config) && fs::is_regular_file(lbm_dir_config) ){
+ std::ifstream file{lbm_dir_config};
+
+ if(file.is_open()){
+ std::stringstream buffer;
+ buffer <<file.rdbuf();
+ return fs::path{buffer.str()};
+ }
+ }
+ }
+
+ return lbm_dir / "data";
}
}
diff --git a/lib/core/c++/particle/particle.hpp b/lib/core/c++/particle/particle.hpp
index 8ef590d..f8104fd 100644
--- a/lib/core/c++/particle/particle.hpp
+++ b/lib/core/c++/particle/particle.hpp
@@ -43,12 +43,15 @@ using ParticleCollisionSpheroid = Struct<
Member<Scalar<T>, "radius">
>;
+tem
+
template<typename T, uint64_t D, typename CollisionType = ParticleCollisionSpheroid<T>>
using Particle = Struct<
Member<ParticleRigidBody<T,D>, "rigid_body">,
Member<CollisionType, "collision">,
// Problem is that dynamic data would two layered
- // Member<FixedArray<Float64,D,D,D>, "mask">,
+ // Member<Array<Float64,D>, "mask">,
+ Member<UInt64, "mask_id">,
Member<Scalar<T>, "mass">
>;
}