Fix SVT implementation and add tests

include/iblgf/domain/decomposition/server.hpp

@@ -105,7 +105,7 @@ class Server : public ServerBase<ServerClientTraits<Domain>>
     };
 
   public: //Ctors
-    using super_type::ServerBase;
+    Server() = default;
 
     Server(const Server& other) = default;
     Server(Server&& other) = default;

include/iblgf/domain/mpi/server_base.hpp

@@ -40,14 +40,13 @@ class ServerBase : public ServerClientBase<Traits>
     using super_type::task_manager_;
 
   public: // ctors
-    using super_type::ServerClientBase;
+    ServerBase() = default;
 
     ServerBase(const ServerBase&) = default;
     ServerBase(ServerBase&&) = default;
     ServerBase& operator=(const ServerBase&) & = default;
     ServerBase& operator=(ServerBase&&) & = default;
     ~ServerBase() = default;
-    ServerBase() = default;
 
   protected: //helper struct
     class ClientInfo

include/iblgf/domain/mpi/task.hpp

@@ -285,7 +285,7 @@ class Task : public Task_base<TaskBuffer<Tag, T, ID>, BufferPolicy, ID>
 {
   public:
     using super_type = Task_base<TaskBuffer<Tag, T, ID>, BufferPolicy, ID>;
-    using super_type::Task_base;
+    using typename super_type::Task_base;
     using id_type = ID;
     using buffer_type = TaskBuffer<Tag, T, ID>;
     using buffer_container_type = typename buffer_type::container_t;

include/iblgf/domain/mpi/task_communicator.hpp

@@ -452,7 +452,7 @@ struct SendTaskCommunicator
         SendTaskCommunicator<typename TaskType::inplace_task_type>;
 
   public: //Ctors
-    using super_type::TaskCommunicator;
+    SendTaskCommunicator() = default;
 
   public: //Memebers:
     int tag_rank_impl(int rank_other) const noexcept
@@ -514,7 +514,7 @@ class RecvTaskCommunicator
         RecvTaskCommunicator<typename TaskType::inplace_task_type>;
 
   public: //Ctors
-    using  super_type::TaskCommunicator;
+    RecvTaskCommunicator() = default;
 
   public: //members
     int  tag_rank_impl(int _rank_other) const noexcept { return _rank_other; }

include/iblgf/domain/octree/iterators.hpp

@@ -196,8 +196,6 @@ struct IteratorBfs : public IteratorBase<Node, IteratorBfs<Node, Dim>>
     using iterator_category = typename iterator_base_type::iterator_category;
 
   public: //ctors
-    using iterator_base_type::IteratorBase;
-
     IteratorBfs(const IteratorBfs&) = default;
     IteratorBfs(IteratorBfs&&) = default;
     IteratorBfs& operator=(const IteratorBfs&) & = default;
@@ -258,8 +256,6 @@ struct IteratorDfs : public IteratorBase<Node, IteratorDfs<Node, Dim>>
     using iterator_category = typename iterator_base_type::iterator_category;
 
   public: //ctors
-    using iterator_base_type::IteratorBase;
-
     IteratorDfs(const IteratorDfs&) = default;
     IteratorDfs& operator=(const IteratorDfs&) & = default;
     ~IteratorDfs() = default;

include/iblgf/domain/svt.hpp

@@ -13,12 +13,17 @@
 #ifndef IBLGF_DOMAIN_STARTING_VORTEX_THEORY_HPP
 #define IBLGF_DOMAIN_STARTING_VORTEX_THEORY_HPP
 #include <iblgf/dictionary/dictionary.hpp>
+#include <iblgf/types.hpp>
+#include <cmath>
+#include <stdexcept>
 
 namespace iblgf
 {
 namespace ib
 {
 using namespace dictionary;
+using namespace types;
+
 template<int Dim>
 class SVT
 {
@@ -43,76 +48,102 @@ class SVT
         auto d2 = dict->get_dictionary("svt");
         p = d2->template get_or<float_type>("p", 2);
         m = d2->template get_or<float_type>("m", 1);
-        beta_hat = d2->template get_or<float_type>("beta_hat", 0.5);
-        rotV = d2->template get_or<bool>("rotV", false);
+        omega0 = d2->template get_or<float_type>("Omega0", 0.0);
+        U0 = d2->template get_or<float_type>("U0", 1.0);
+        alpha0 = d2->template get_or<float_type>("alpha0", 0.0);
         x0 = d2->template get_or<float_type>("x0", 0.0);
+        // beta_hat = d2->template get_or<float_type>("beta_hat", 0.5);
+        // beta = d2->template get_or<float_type>("beta", 1.0);
+        // x0 = d2->template get_or<float_type>("x0", 0.0);
     }
 
     float_type operator()(std::size_t idx, float_type t, real_coordinate_type coord = real_coordinate_type{})
     {
-        if (!rotV) // seperate U and V power laws
+        float_type alpha = alpha0;
+        if (p >= 0) alpha = alpha0 + 1 / (p + 1) * omega0 * std::pow(t, p + 1); // for pure translation in y set m for power and alpha=pi/2,p=-1
+        if constexpr (Dim == 2)
         {
-            if (Dim == 3 && idx == 1) return 0.0;
-            if ((Dim == 2 && idx == 0) || (Dim == 3 && idx == 2)) //tangent direction, for just tangent set p<0 
+            // float_type f_alpha;
+            // float_type f_u;
+            // if (idx == 0) f_alpha=-coord[1]*omega0; //tangent direction
+            // else f_alpha = (coord[0]-x0)*omega0; //normal direction
+            if (idx == 0) // tangent direction
             {
-                if (m < 0) return 0.0;
-                float_type h1;
-                if (p < 0) h1=1.0; //beta hat not used because no normal component
-                else h1 = 1.0 / std::sqrt(4 * std::pow(beta_hat, 2) + 1);
-                return -h1 * std::pow(t, m);
+                float_type u_tan;
+                if (m < 0) u_tan = 0.0;
+                else
+                    u_tan = U0 * std::cos(alpha) * std::pow(t, m);
+
+                float_type u_rot;
+                if (p < 0) u_rot = 0.0;
+                else
+                    u_rot = -coord[1] * omega0 * std::pow(t, p);
+
+                return -(u_tan + u_rot);
+            }
+            else if (idx == 1) //normal direction
+            {
+                float_type u_tan;
+                if (m < 0) u_tan = 0.0;
+                else
+                    u_tan = -U0 * std::sin(alpha) * std::pow(t, m);
+
+                float_type u_rot;
+                if (p < 0) u_rot = 0.0;
+                else
+                    u_rot = (coord[0] - x0) * omega0 * std::pow(t, p);
+
+                return -(u_tan + u_rot);
             }
-            else if ((Dim == 2 && idx == 1) ||
-                     (Dim == 3 && idx == 0)) //normal direction, for just normal set m<0
+            else
             {
-                if (p < 0) return 0.0;
-                float_type h1;
-                if (m < 0 ) { h1 = 1.0; }
-                else { h1 = 2.0 * beta_hat / std::sqrt(4 * std::pow(beta_hat, 2) + 1); }
-                return -h1 * std::pow(t, p);
+                throw std::runtime_error("idx should be 0 or 1 for 2D problems");
             }
         }
-        else // combined rotation
+        else if constexpr (Dim == 3)
         {
-            float_type f_alpha = 0.0;
-            float_type f_u = 0.0;
-            // if (Dim != 2) throw std::runtime_error("rotV option is only implemented for 2-D problems");
-            if ((Dim==2&&idx == 0) || (Dim==3&&idx==2)) //tangent direction
+            if (idx == 1) return 0.0; //no velocity in y direction
+            if (idx == 2)             //tangent direction
             {
-                if (m < 0) return 0.0;
-                float_type h1;
-                if (p < 0) h1=1.0; //beta hat not used because no normal component
-                else h1 = 1.0 / std::sqrt(4 * std::pow(beta_hat, 2) + 1);
-                return -h1 * std::pow(t, m);
+                float_type u_tan;
+                if (m < 0) u_tan = 0.0;
+                else
+                    u_tan = U0 * std::cos(alpha) * std::pow(t, m);
+
+                float_type u_rot;
+                if (p < 0) u_rot = 0.0;
+                else
+                    u_rot = -coord[0] * omega0 * std::pow(t, p);
+
+                return -(u_tan + u_rot);
             }
-            else if ((Dim==2&&idx == 1) || (Dim==3&&idx==0)) //normal directions
+            else if (idx == 0) //normal direction
             {
-                if (p < 0) return 0.0;
-                float_type h1;
-                if (m < 0 ) { h1 = 2.0; }
-                else { h1 = (4 * beta_hat) / std::sqrt(4 * std::pow(beta_hat, 2) + 1); }
-                float_type coord_n;
-                if(Dim==2)
-                    coord_n=coord[0]-x0;
-                else //3D
-                    coord_n=coord[2]-x0; //
-                return -h1 * coord_n * std::pow(t, p);
-                // f_alpha = ((coord[0] - x0) * beta_hat * 4)/ std::sqrt(4 * std::pow(beta_hat, 2) + 1);
-                // f_u = 0.0;
+                float_type u_tan;
+                if (m < 0) u_tan = 0.0;
+                else
+                    u_tan = -U0 * std::sin(alpha) * std::pow(t, m);
+
+                float_type u_rot;
+                if (p < 0) u_rot = 0.0;
+                else
+                    u_rot = (coord[2] - x0) * omega0 * std::pow(t, p);
+
+                return -(u_tan + u_rot);
+            }
+            else
+            {
+                throw std::runtime_error("idx should be 0,1, 2 for 3D problems");
             }
-			// float_type val_u=std::pow(t,m);
-			// float_type val_a=std::pow(t,p);
-			// return -(f_u * val_u + f_alpha * val_a);
         }
-        return 0.0;
     }
-
-  private:
-    float_type        p;
-    float_type        m;
-    float_type beta_hat;
-    float_type d;
+    float_type p;
+    float_type m;
     bool       rotV;
     float_type x0;
+    float_type omega0;
+    float_type U0;
+    float_type alpha0;
 };
 } // namespace ib
 

tests/CMakeLists.txt

@@ -32,6 +32,7 @@ add_subdirectory (triple_decomposition)
 #add_subdirectory (helmholtz)
 add_subdirectory (Poisson2D)
 add_subdirectory(conv_test)
+add_subdirectory(svt)
 #add_subdirectory (operators_jac)
 #add_subdirectory (helmholtz_FFT)
 #add_subdirectory (operators_helmholtz)

tests/server_client/client.hpp

@@ -46,7 +46,11 @@ class Client : public ClientBase<ClientTraits>
     Client& operator=(Client&&) & = default;
     ~Client() = default;
 
-    using super_type::ClientBase;
+    Client() = default;
+    Client(int _server_rank = 0)
+    : super_type(_server_rank)
+    {
+    }
 
   public:
     auto test()

tests/svt/CMakeLists.txt

@@ -0,0 +1,27 @@
+#      ▄▄▄▄▄▄▄▄▄▄▄  ▄▄▄▄▄▄▄▄▄▄   ▄            ▄▄▄▄▄▄▄▄▄▄▄  ▄▄▄▄▄▄▄▄▄▄▄
+#     ▐░░░░░░░░░░░▌▐░░░░░░░░░░▌ ▐░▌          ▐░░░░░░░░░░░▌▐░░░░░░░░░░░▌
+#      ▀▀▀▀█░█▀▀▀▀ ▐░█▀▀▀▀▀▀▀█░▌▐░▌          ▐░█▀▀▀▀▀▀▀▀▀ ▐░█▀▀▀▀▀▀▀▀▀
+#          ▐░▌     ▐░▌       ▐░▌▐░▌          ▐░▌          ▐░▌
+#          ▐░▌     ▐░█▄▄▄▄▄▄▄█░▌▐░▌          ▐░▌ ▄▄▄▄▄▄▄▄ ▐░█▄▄▄▄▄▄▄▄▄
+#          ▐░▌     ▐░░░░░░░░░░▌ ▐░▌          ▐░▌▐░░░░░░░░▌▐░░░░░░░░░░░▌
+#          ▐░▌     ▐░█▀▀▀▀▀▀▀█░▌▐░▌          ▐░▌ ▀▀▀▀▀▀█░▌▐░█▀▀▀▀▀▀▀▀▀
+#          ▐░▌     ▐░▌       ▐░▌▐░▌          ▐░▌       ▐░▌▐░▌
+#      ▄▄▄▄█░█▄▄▄▄ ▐░█▄▄▄▄▄▄▄█░▌▐░█▄▄▄▄▄▄▄▄▄ ▐░█▄▄▄▄▄▄▄█░▌▐░▌
+#     ▐░░░░░░░░░░░▌▐░░░░░░░░░░▌ ▐░░░░░░░░░░░▌▐░░░░░░░░░░░▌▐░▌
+#      ▀▀▀▀▀▀▀▀▀▀▀  ▀▀▀▀▀▀▀▀▀▀   ▀▀▀▀▀▀▀▀▀▀▀  ▀▀▀▀▀▀▀▀▀▀▀  ▀
+
+set(svt_test svt)
+
+foreach (item ${svt_test})
+    add_executable(${item} ${item}.cpp ${SOURCES})
+    target_link_libraries(${item} iblgf_test_main)
+    add_test(
+        NAME ${item}
+        COMMAND ${MPIEXEC_EXECUTABLE} ${MPIEXEC_PREFLAGS} ${MPIEXEC_NUMPROC_FLAG} 1 $<TARGET_FILE:${item}>
+        # COMMAND ${CMAKE_COMMAND} -E env
+        #     ASAN_OPTIONS=detect_leaks=0
+        #     LSAN_OPTIONS=detect_leaks=0
+        #     ${MPIEXEC_EXECUTABLE} ${MPIEXEC_NUMPROC_FLAG} 1 ${MPIEXEC_PREFLAGS}
+        #     $<TARGET_FILE:datafields>
+    )
+endforeach()