From c0549d71b2109f10c1238db8b22362e7826ba61b Mon Sep 17 00:00:00 2001 From: "Claudius \"keldu\" Holeksa" Date: Sun, 5 Jul 2026 15:59:23 +0200 Subject: Just rename from lib to modules --- modules/core/c++/math/math.hpp | 4 ++ modules/core/c++/math/n_closest.hpp | 54 ++++++++++++++ modules/core/c++/math/n_linear.hpp | 139 ++++++++++++++++++++++++++++++++++++ modules/core/c++/math/round.hpp | 26 +++++++ 4 files changed, 223 insertions(+) create mode 100644 modules/core/c++/math/math.hpp create mode 100644 modules/core/c++/math/n_closest.hpp create mode 100644 modules/core/c++/math/n_linear.hpp create mode 100644 modules/core/c++/math/round.hpp (limited to 'modules/core/c++/math') diff --git a/modules/core/c++/math/math.hpp b/modules/core/c++/math/math.hpp new file mode 100644 index 0000000..3920bec --- /dev/null +++ b/modules/core/c++/math/math.hpp @@ -0,0 +1,4 @@ +#pragma once + +#include "n_linear.hpp" +#include "n_closest.hpp" diff --git a/modules/core/c++/math/n_closest.hpp b/modules/core/c++/math/n_closest.hpp new file mode 100644 index 0000000..ac0fe2f --- /dev/null +++ b/modules/core/c++/math/n_closest.hpp @@ -0,0 +1,54 @@ +#pragma once + +#include "../common.hpp" +#include "../iterator.hpp" + +namespace kel { +namespace lbm { + +template +saw::data n_closest_read(const saw::data,Encode>& f, const saw::data>& frac_ind){ + + auto shift_frac_ind = frac_ind; + for(uint64_t i{0u}; i < D; ++i){ + + shift_frac_ind.at({{i}}) = shift_frac_ind.at({{i}}) + saw::data{0.5}; + if(shift_frac_ind.at({{i}}).get() < 0){ + shift_frac_ind.at({{i}}) = {}; + } + } + + saw::data> shift_ind; + for(uint64_t i{0u}; i < D; ++i){ + shift_ind.at({i}) = frac_ind.at({{i}}).template cast_to(); + } + + return f.at(shift_ind); +} + +template +void n_closest_add(const saw::data,Encode>& f, const saw::data>& frac_ind, const saw::data& val){ + auto shift_frac_ind = frac_ind; + for(uint64_t i{0u}; i < D; ++i){ + + shift_frac_ind.at({{i}}) = shift_frac_ind.at({{i}}) + saw::data{0.5}; + if(shift_frac_ind.at({{i}}).get() < 0){ + shift_frac_ind.at({{i}}) = {}; + } + } + + auto f_meta = f.meta(); + saw::data> shift_ind; + for(uint64_t i{0u}; i < D; ++i){ + shift_ind.at({i}) = frac_ind.at({{i}}).template cast_to(); + if(shift_ind.at({i}) < f_meta.at({i})){ + shift_ind.at({i}) = f_meta.at({i}) - 1u; + } + } + auto& f_i = f.at(shift_ind); + + f_i = f_i + val; +} + +} +} diff --git a/modules/core/c++/math/n_linear.hpp b/modules/core/c++/math/n_linear.hpp new file mode 100644 index 0000000..b378440 --- /dev/null +++ b/modules/core/c++/math/n_linear.hpp @@ -0,0 +1,139 @@ +#pragma once + +#include "../common.hpp" +#include "../iterator.hpp" + +namespace kel { +namespace lbm { +namespace impl { +template +struct n_linear_interpolate_helper final { + template + auto apply(const saw::data& field, const saw::data>& pos){ + return pos; + } +}; +} + +template +saw::data,sch::Vector>> position_to_index_and_fraction(const saw::data>& pos){ + saw::data,sch::Vector>> sep; + + auto& ind = sep.template get<0u>(); + auto& frac = sep.template get<1u>(); + + auto pos_cpy = pos; + // Guarantee that the pos is at least 0 + for(uint64_t i = 0u; i < D; ++i){ + pos_cpy.at({{i}}).set(std::max(pos.at({{i}}).get(), static_cast::type>(0))); + } + + // Now we can cast to uint64_t + for(uint64_t i = 0u; i < D; ++i){ + ind.at({{i}}) = pos_cpy.at({{i}}).template cast_to(); + } + + frac = pos_cpy - ind.template cast_to(); + + return sep; +} + +template +auto floor_index_from_position(const saw::data>& pos){ + return position_to_index_and_fraction(pos).template get<0u>(); +} + +template +saw::data,sch::Vector>> position_to_index_and_fraction_bounded( + const saw::data>& pos, + const saw::data>& bound) +{ + auto infr = position_to_index_and_fraction(pos); + auto& ind = infr.template get<0u>(); + auto& fra = infr.template get<1u>(); + for(uint64_t i = 0u; i < D; ++i){ + // If index is higher than bound. Set to bound and reset fraction + if((ind.at({{i}}).get()+1u) >= bound.at({{i}}).get()){ + ind.at({{i}}).set(bound.at({{i}}).get()-1u); + fra.at({{i}}) = {}; + } + } + return infr; +} + + +template +auto n_linear_interpolate( + const saw::data& field, const saw::data>& pos){ + + // Pos + auto pos_bound = pos; + + // Dimensions + auto meta = field.dims(); + + // Lower Index + saw::data> ind; + + for(saw::data i{0u}; i < saw::data{D}; ++i){ + // Native Positive i + auto npos_i = pos.at({i}).get(); + + { + // Ok I want to stay in bounds + npos_i = std::min(npos_i,meta.at(i).get()-1.0); + npos_i = std::max(npos_i,1.0); + } + + // Native Index i + auto nind_i = static_cast(std::floor(npos_i))-1ul; + + // Set index to i + ind.at(i).set(nind_i); + } + saw::data> pos_frac; + for(saw::data i{0u}; i < saw::data{D}; ++i){ + pos_frac.at({i}) = pos_bound.at({i}) - ind.at(i).template cast_to(); + } + + // Base value + saw::data res; + + // constexpr uint64_t d_corners = 1ul << D; + + saw::data> ones_ind; + for(saw::data i{0u}; i < saw::data{D}; ++i){ + ones_ind.at({i}).set(1u); + } + + iterator::apply([&](auto ind){ + // Iterates over (0,0,0) to (1,1,1) + saw::data weight{1.0}; + + for(saw::data d{0u}; d < saw::data{D}; ++d){ + + saw::data t = pos_frac.at({d}); + + if(ind.at(d).get() == 0u){ + weight = weight * (saw::data{1} - t); + }else{ + weight = weight * t; + } + } + + }, {}, ones_ind); + + /// TODO I need to actually calc stuff + return field.at({}); +} + +template +saw::data> bilinear_interpolate(const saw::data& field, const saw::data>& pos){ + saw::data> res; + + { + } + return {}; +} +} +} diff --git a/modules/core/c++/math/round.hpp b/modules/core/c++/math/round.hpp new file mode 100644 index 0000000..d3a2586 --- /dev/null +++ b/modules/core/c++/math/round.hpp @@ -0,0 +1,26 @@ +#pragma once + +#include "../common.hpp" + +namespace kel { +namespace lbm { + +template +saw::data> round_to_unsigned(const saw::data>& inp){ + saw::data> rv; + + auto zero = static_cast::type>(0); + auto half = static_cast::type>(0.5); + + for(uint64_t i{0u}; i < D; ++i){ + auto val = inp.at({{i}}).get()+half; + val = std::max(zero,val); + + rv.at({i}).set(static_cast(val)); + } + + return rv; +} + +} +} -- cgit v1.2.3