first commit

include/dmt/fdmt_base.hpp

@@ -46,15 +46,18 @@ class FDMT {
     virtual ~FDMT()              = default;
     float get_df() const;
     float get_correction() const;
+    int get_m_nsamps() const;
     SizeType get_niters() const;
     const FDMTPlan& get_plan() const;
     const DtGridType& get_dt_grid_final() const;
     std::vector<float> get_dm_grid_final() const;
-    static void set_log_level(int level);
-    virtual void execute(const float* waterfall, size_t waterfall_size,
-                         float* dmt, size_t dmt_size)             = 0;
-    virtual void initialise(const float* waterfall, float* state) = 0;
+    //static void set_log_level(int level);
+    virtual void execute(const float* __restrict waterfall, size_t waterfall_size,
+                         float* __restrict dmt, size_t dmt_size)             = 0;
+    virtual void initialise(const float* __restrict waterfall, float* __restrict state) = 0;
 
+
+
 protected:
     void check_inputs(size_t waterfall_size, size_t dmt_size) const;
 

include/dmt/fdmt_cpu.hpp

@@ -1,20 +1,20 @@
 #pragma once
 
-#include <dmt/fdmt_base.hpp>
+#include <fdmt_base.hpp>
 
 class FDMTCPU : public FDMT {
 public:
     FDMTCPU(float f_min, float f_max, size_t nchans, size_t nsamps, float tsamp,
             size_t dt_max, size_t dt_step = 1, size_t dt_min = 0);
     static void set_num_threads(int nthreads);
-    void execute(const float* waterfall, size_t waterfall_size, float* dmt,
+    void execute(const float* __restrict waterfall, size_t waterfall_size, float* __restrict dmt,
                  size_t dmt_size) override;
-    void initialise(const float* waterfall, float* state) override;
+    void initialise(const float* __restrict waterfall, float* __restrict state) override;
 
 private:
     // Buffers
     std::vector<float> m_state_in;
     std::vector<float> m_state_out;
 
-    void execute_iter(const float* state_in, float* state_out, SizeType i_iter);
+    void execute_iter(const float* __restrict state_in, float* __restrict state_out, SizeType i_iter);
 };

include/dmt/fdmt_gpu.cuh

@@ -0,0 +1,126 @@
+#pragma once
+
+#include <thrust/device_vector.h>
+#include <fdmt_base.hpp>
+
+using StShapeTypeD   = int4;
+using FDMTCoordTypeD = int2;
+
+struct FDMTCoordMappingD {
+    FDMTCoordTypeD head;
+    FDMTCoordTypeD tail;
+    int offset;
+};
+
+struct FDMTPlanD
+{  
+    // is analogue of std::vector<std::vector<FDMTCoordType>> coordinates;
+    // "coordinates_d" is flattened vector "coordinates" 
+    // in vector "len_inner_vects_coordinates_cumsum" we store cummulative sums of elements  inner vectors of "coordinates" 
+    // So, 
+    // len_inner_vects_coordinates_cumsum[0] = 0
+    // len_inner_vects_coordinates_cumsum[1] = len_inner_vects_coordinates_cumsum[0] + coordinates[0].size() *2
+    // ...
+    // len_inner_vects_coordinates_cumsum[n] = len_inner_vects_coordinates_cumsum[n-1] + coordinates[n-1].size() *2
+    // ...
+    // Remember that always:  len_inner_vects_coordinates_cumsum.size() = NUmIter +1
+    thrust::device_vector<int> coordinates_d;
+    std::vector<int> lenof_innerVects_coords_cumsum_h;
+
+    // Is an analogues as previous
+    thrust::device_vector<int> coordinates_to_copy_d;
+    std::vector<int> lenof_innerVects_coords_to_copy_cumsum_h;
+
+
+    // It is analogue of:   std::vector<std::vector<FDMTCoordMapping>> mappings;
+    // each FDMTCoordMapping consists of 5 elements
+    // "mappings_d" is flattened vector "mappings" 
+    // in vector "len_mappings_cumsum" we store cummulative sums of elements  inner vectors of "mappings" 
+    // So, 
+    // len_mappings_cumsum[0] = 0
+    // len_mappings_cumsum[1] = len_mappings_cumsum[0] + mappings[0].size() *5
+    // ...
+    // len_mappings_cumsum[n] = len_mappings_cumsum[n-1] + mappings[n-1].size() *5
+    // ...
+    // Remember that always:  len_mappings_cumsum.size() = NUmIter +1
+    thrust::device_vector<int> mappings_d;
+    std::vector<int> len_mappings_cumsum_h;
+
+
+    // Is an analogues as previous
+    thrust::device_vector<int> mappings_to_copy_d;
+    std::vector<int> len_mappings_to_copy_cumsum_h;
+
+
+    // It is analogue of state_sub_idx
+    // Has size: state_sub_idx_d.size() = m_niters +1
+    thrust::device_vector<int>state_sub_idx_d;
+    std::vector<int> len_state_sub_idx_cumsum_h;
+
+
+    // It is analogue of dt_grid
+    // Has size: dt_grid_d.size() = m_niters +1
+    thrust::device_vector<int>dt_grid_d;
+    thrust::device_vector<int> pos_gridInnerVects_d;
+    std::vector<int> pos_gridSubVects_h;
+
+};
+class FDMTGPU : public FDMT {
+public:
+    FDMTGPU(float f_min, float f_max, size_t nchans, size_t nsamps, float tsamp,
+            size_t dt_max, size_t dt_step = 1, size_t dt_min = 0);
+    void execute(const float* __restrict  waterfall, size_t waterfall_size, float* __restrict  dmt,
+                 size_t dmt_size) override;
+
+
+
+
+private:
+    thrust::device_vector<float> m_state_in_d;
+    thrust::device_vector<float> m_state_out_d;
+
+
+    FDMTPlanD m_fdmt_plan_d;
+
+    static void transfer_plan_to_device(const FDMTPlan& plan,
+                                        FDMTPlanD& plan_d);
+    void initialise(const float* __restrict waterfall, float* __restrict  state) override;
+
+};
+
+std::vector<SizeType>  flatten_mappings(const std::vector<std::vector<FDMTCoordMapping>>& mappings);
+
+__global__
+void  kernel_init_fdmt(const float* __restrict waterfall, int* __restrict p_state_sub_idx
+    , int* __restrict p_dt_grid, int* __restrict p_pos_gridInnerVects, float* __restrict state, const int nsamps);
+
+__global__
+void  kernel_init_fdmt_v1(const float* __restrict waterfall, int* __restrict p_state_sub_idx
+    , int* __restrict p_dt_grid, int* __restrict p_pos_gridInnerVects, float* __restrict state, const int nsamps);
+
+__global__
+void kernel_execute_iter(const float* __restrict state_in, float* __restrict state_out
+    , int* __restrict pcoords_cur//+   
+    , int* __restrict pmappings_cur //+
+    , int* __restrict pcoords_copy_cur //+
+    , int* __restrict pmappings_copy_cur //+
+    , int* __restrict pstate_sub_idx_cur //+
+    , int* __restrict pstate_sub_idx_prev//+
+    , int nsamps //+
+    , int  coords_cur_size//+
+    , int  coords_copy_cur_size//+
+);
+
+
+__global__
+void kernel_execute_iter_v1(const float* __restrict state_in, float* __restrict state_out
+    , int* __restrict pcoords_cur//+   
+    , int* __restrict pmappings_cur //+
+    , int* __restrict pcoords_copy_cur //+
+    , int* __restrict pmappings_copy_cur //+
+    , int* __restrict pstate_sub_idx_cur //+
+    , int* __restrict pstate_sub_idx_prev//+
+    , int nsamps //+
+    , int  coords_cur_size//+
+    , int  coords_copy_cur_size//+
+);

lib/fdmt_base.cpp

@@ -3,10 +3,10 @@
 #include <cstddef>
 #include <cstdint>
 #include <numeric>
-#include <spdlog/spdlog.h>
-
-#include "dmt/fdmt_utils.hpp"
-#include <dmt/fdmt_base.hpp>
+//#include <spdlog/spdlog.h>
+#include <stdexcept>  // For std::invalid_argument
+#include "fdmt_utils.hpp"
+#include <fdmt_base.hpp>
 
 size_t FDMTPlan::calculate_memory_usage() const {
     size_t mem_use = 0;
@@ -48,13 +48,15 @@ FDMT::FDMT(float f_min, float f_max, SizeType nchans, SizeType nsamps,
       m_correction(m_df / 2),
       m_niters(calculate_niters(m_nchans)) {
     configure_fdmt_plan();
-    spdlog::debug("FDMT: df={}, dt_max={}, dt_min={}, dt_step={}, niters={}",
-                  m_df, m_dt_max, m_dt_min, m_dt_step, m_niters);
+   /* spdlog::debug("FDMT: df={}, dt_max={}, dt_min={}, dt_step={}, niters={}",
+                  m_df, m_dt_max, m_dt_min, m_dt_step, m_niters);*/
 }
 
 // Getters
 float FDMT::get_df() const { return m_df; }
 float FDMT::get_correction() const { return m_correction; }
+int FDMT::get_m_nsamps() const { return static_cast<int>(m_nsamps); };
+
 SizeType FDMT::get_niters() const { return m_niters; }
 const FDMTPlan& FDMT::get_plan() const { return m_fdmt_plan; }
 const DtGridType& FDMT::get_dt_grid_final() const {
@@ -73,13 +75,13 @@ std::vector<float> FDMT::get_dm_grid_final() const {
 }
 
 // Setters
-void FDMT::set_log_level(int level) {
-    if (level < static_cast<int>(spdlog::level::trace) ||
-        level > static_cast<int>(spdlog::level::off)) {
-        spdlog::set_level(spdlog::level::info);
-    }
-    spdlog::set_level(static_cast<spdlog::level::level_enum>(level));
-}
+//void FDMT::set_log_level(int level) {
+//    if (level < static_cast<int>(spdlog::level::trace) ||
+//        level > static_cast<int>(spdlog::level::off)) {
+//        spdlog::set_level(spdlog::level::info);
+//    }
+//    spdlog::set_level(static_cast<spdlog::level::level_enum>(level));
+//}
 
 // Private methods
 SizeType FDMT::calculate_niters(SizeType nchans) {
@@ -111,7 +113,7 @@ void FDMT::check_inputs(size_t waterfall_size, size_t dmt_size) const {
         plan.state_shape[m_niters][3] * plan.state_shape[m_niters][4]) {
         throw std::invalid_argument("Invalid size of dmt");
     }
-    spdlog::debug("FDMT: Input dimensions: {}x{}", m_nchans, m_nsamps);
+   // spdlog::debug("FDMT: Input dimensions: {}x{}", m_nchans, m_nsamps);
 }
 
 void FDMT::configure_fdmt_plan() {
@@ -132,7 +134,7 @@ void FDMT::configure_fdmt_plan() {
     for (SizeType i_iter = 1; i_iter < m_niters + 1; ++i_iter) {
         make_fdmt_plan(i_iter);
     }
-    spdlog::debug("FDMT: configured fdmt plan");
+   // spdlog::debug("FDMT: configured fdmt plan");
 }
 
 void FDMT::make_fdmt_plan_iter0() {