Add model_base class

CMakeLists.txt

@@ -92,6 +92,7 @@ if(USE_CXX)
     )
 
     set(SOURCES_CXX 
+            srcCXX/model_base.cpp
             srcCXX/sfx_esm.cpp
             srcCXX/sfx_sheba.cpp
             srcCXX/sfx_compute_sheba.cpp
@@ -101,6 +102,8 @@ if(USE_CXX)
             includeCU/sfx_surface.cuh
             includeCU/sfx_math.cuh
             includeCU/sfx_esm_compute_subfunc.cuh
+
+            includeCXX/model_base.h
             includeCXX/sfx_esm.h
             includeCXX/sfx_sheba.h
             includeCXX/sfx_compute_sheba.h

includeC/sfx_data.h

@@ -54,28 +54,6 @@ extern "C" {
         float *Pr_t_inv;
     };
 
-    // use sfx_esm_param
-    struct sfx_esm_param
-    {
-        float kappa;
-        float Pr_t_0_inv;
-        float Pr_t_inf_inv;
-
-        float alpha_m;
-        float alpha_h;
-        float alpha_h_fix;
-        float beta_m;
-        float beta_h;
-        float Rib_max;
-    };
-
-    struct sfx_esm_numericsTypeC
-    {
-        int maxiters_convection;
-        int maxiters_charnock;
-    };
-    
-
     struct sfx_surface_param
     {
         int surface_ocean;
@@ -98,19 +76,6 @@ extern "C" {
         float B_max_land;
     };
 
-
-    // TODO: add later after esm implementation
-    // struct sfx_sheba_param
-    // {
-    //     float alpha_m;
-    //     float alpha_h;
-    //     float a_m;
-    //     float b_m;
-    //     float a_h;
-    //     float b_h;
-    //     float c_h;
-    // };
-
     struct sfx_phys_constants
     {
         float Pr_m;
@@ -118,6 +83,47 @@ extern "C" {
         float nu_air;
     };
 
+    // use sfx_esm_param
+    struct sfx_esm_param
+    {
+        float kappa;
+        float Pr_t_0_inv;
+        float Pr_t_inf_inv;
+
+        float alpha_m;
+        float alpha_h;
+        float alpha_h_fix;
+        float beta_m;
+        float beta_h;
+        float Rib_max;
+    };
+
+    struct sfx_esm_numericsTypeC
+    {
+        int maxiters_convection;
+        int maxiters_charnock;
+    };
+
+    // use sfx_sheba_param
+    struct sfx_sheba_param
+    {
+        float kappa;
+        float Pr_t_0_inv;
+
+        float alpha_m;
+        float alpha_h;
+        float a_m;
+        float b_m;
+        float a_h;
+        float b_h;
+        float c_h;
+    };
+
+    struct sfx_sheba_numericsTypeC
+    {
+        int maxiters_charnock;
+    };
+    
 #ifdef __cplusplus
 }
 #endif
\ No newline at end of file

includeCXX/model_base.h

+#pragma once
+#include <cstddef>
+#include "sfx_template_parameters.h"
+#include "../includeC/sfx_data.h"
+
+template<typename T, MemType memIn, MemType memOut, MemType RunMem >
+class ModelBase
+{
+public:
+    sfxDataVecTypeC* res_sfx;
+    sfxDataVecTypeC sfx;
+    meteoDataVecTypeC meteo;
+    sfx_surface_param surface_param;
+    sfx_phys_constants phys_constants;
+
+    int grid_size;
+    bool ifAllocated;
+    size_t allocated_size;
+    
+    ModelBase(sfxDataVecTypeC* sfx,
+                meteoDataVecTypeC* meteo,
+                const sfx_surface_param surface_param,
+                const sfx_phys_constants phys_constants,
+                const int grid_size);
+    ~ModelBase();
+};
\ No newline at end of file

includeCXX/sfx_call_class_func.h

@@ -13,21 +13,13 @@ extern "C" {
                            const struct sfx_phys_constants* constants,
                            const int grid_size);
 
-    void surf_flux_sheba_CXX (float *zeta_, float *Rib_, float *Re_, float *B_, float *z0_m_, float *z0_t_, float *Rib_conv_lim_, float *Cm_, float *Ct_, float *Km_, float *Pr_t_inv_,
-    float *U_, float *dT_, float *Tsemi_, float *dQ_, float *h_, float *in_z0_m_,
-    const float kappa, const float Pr_t_0_inv, 
-    const float alpha_m, const float alpha_h, 
-    const float a_m, const float a_h, 
-    const float b_m, const float b_h,
-    const float c_h,
-    const float Re_rough_min, 
-    const float B1_rough, const float B2_rough,
-    const float B_max_land, const float B_max_ocean, const float B_max_lake,
-    const float gamma_c, const float Re_visc_min,
-    const float Pr_m, const float nu_air, const float g, 
-    const int maxiters_charnock, 
+    void surf_flux_sheba_CXX(struct sfxDataVecTypeC* sfx,
+                           struct meteoDataVecTypeC* meteo,
+                           const struct sfx_sheba_param* model_param, 
+                           const struct sfx_surface_param* surface_param,
+                           const struct sfx_sheba_numericsTypeC* numerics,
+                           const struct sfx_phys_constants* constants,
                            const int grid_size);
-
 #ifdef __cplusplus
 }
 #endif
\ No newline at end of file

includeCXX/sfx_esm.h

 #pragma once
 #include "sfx_template_parameters.h"
+#include "../includeCXX/model_base.h"
 #include "../includeC/sfx_data.h"
-// #include <cstddef>
 
 template<typename T, MemType memIn, MemType memOut, MemType RunMem >
-class FluxEsmBase
+class FluxEsmBase : public ModelBase<T, memIn, memOut, RunMem>
 {
-protected:
-    struct sfxDataVecTypeC sfx;
-    struct meteoDataVecTypeC meteo;
-    struct sfx_esm_param model_param;
-    struct sfx_surface_param surface_param;
-    struct sfx_esm_numericsTypeC numerics;
-    struct sfx_phys_constants phys_constants;
+public:
+    using ModelBase<T, memIn, memOut, RunMem>::res_sfx;
+    using ModelBase<T, memIn, memOut, RunMem>::sfx;
+    using ModelBase<T, memIn, memOut, RunMem>::meteo;
+    using ModelBase<T, memIn, memOut, RunMem>::surface_param;
+    using ModelBase<T, memIn, memOut, RunMem>::phys_constants;
+    using ModelBase<T, memIn, memOut, RunMem>::grid_size;
+    using ModelBase<T, memIn, memOut, RunMem>::ifAllocated;
+    using ModelBase<T, memIn, memOut, RunMem>::allocated_size;
 
-    struct sfxDataVecTypeC* res_sfx;
+    sfx_esm_param model_param;
+    sfx_esm_numericsTypeC numerics;
 
-    int grid_size;
-    bool ifAllocated;
-    size_t allocated_size;
-public:
-    FluxEsmBase() = default;
-    FluxEsmBase(struct sfxDataVecTypeC* sfx,
-                struct meteoDataVecTypeC* meteo,
-                const struct sfx_esm_param model_param, 
-                const struct sfx_surface_param surface_param,
-                const struct sfx_esm_numericsTypeC numerics,
-                const struct sfx_phys_constants phys_constants,
+    FluxEsmBase(sfxDataVecTypeC* sfx,
+                meteoDataVecTypeC* meteo,
+                const sfx_esm_param model_param, 
+                const sfx_surface_param surface_param,
+                const sfx_esm_numericsTypeC numerics,
+                const sfx_phys_constants phys_constants,
                 const int grid_size);
     ~FluxEsmBase();
 };
@@ -41,23 +39,22 @@ class FluxEsm<T, memIn, memOut, MemType::CPU> : public FluxEsmBase<T, memIn, mem
     using FluxEsmBase<T, memIn, memOut, MemType::CPU>::res_sfx;
     using FluxEsmBase<T, memIn, memOut, MemType::CPU>::sfx;
     using FluxEsmBase<T, memIn, memOut, MemType::CPU>::meteo;
-    using FluxEsmBase<T, memIn, memOut, MemType::CPU>::model_param;
     using FluxEsmBase<T, memIn, memOut, MemType::CPU>::surface_param;
-    using FluxEsmBase<T, memIn, memOut, MemType::CPU>::numerics;
     using FluxEsmBase<T, memIn, memOut, MemType::CPU>::phys_constants;
     using FluxEsmBase<T, memIn, memOut, MemType::CPU>::grid_size;
     using FluxEsmBase<T, memIn, memOut, MemType::CPU>::ifAllocated;
     using FluxEsmBase<T, memIn, memOut, MemType::CPU>::allocated_size;
+    using FluxEsmBase<T, memIn, memOut, MemType::CPU>::model_param;
+    using FluxEsmBase<T, memIn, memOut, MemType::CPU>::numerics;
 public:
-    FluxEsm() = default;
-    FluxEsm(struct sfxDataVecTypeC* sfx,
-                struct meteoDataVecTypeC* meteo,
-                const struct sfx_esm_param model_param, 
-                const struct sfx_surface_param surface_param,
-                const struct sfx_esm_numericsTypeC numerics,
-                const struct sfx_phys_constants phys_constants,
-                const int grid_size) : FluxEsmBase<T, memIn, memOut, MemType::CPU>(sfx, meteo, 
-                model_param, surface_param, numerics, phys_constants, grid_size) {};
+    FluxEsm(sfxDataVecTypeC* sfx,
+                meteoDataVecTypeC* meteo,
+                const sfx_esm_param model_param, 
+                const sfx_surface_param surface_param,
+                const sfx_esm_numericsTypeC numerics,
+                const sfx_phys_constants phys_constants,
+                const int grid_size) : FluxEsmBase<T, memIn, memOut, MemType::CPU>(sfx, meteo, model_param, 
+                                       surface_param, numerics, phys_constants, grid_size) {}
     ~FluxEsm() = default;
     void compute_flux();
 };
@@ -69,23 +66,22 @@ class FluxEsm<T, memIn, memOut, MemType::GPU> : public FluxEsmBase<T, memIn, mem
     using FluxEsmBase<T, memIn, memOut, MemType::GPU>::res_sfx;
     using FluxEsmBase<T, memIn, memOut, MemType::GPU>::sfx;
     using FluxEsmBase<T, memIn, memOut, MemType::GPU>::meteo;
-    using FluxEsmBase<T, memIn, memOut, MemType::GPU>::model_param;
     using FluxEsmBase<T, memIn, memOut, MemType::GPU>::surface_param;
-    using FluxEsmBase<T, memIn, memOut, MemType::GPU>::numerics;
     using FluxEsmBase<T, memIn, memOut, MemType::GPU>::phys_constants;
     using FluxEsmBase<T, memIn, memOut, MemType::GPU>::grid_size;
     using FluxEsmBase<T, memIn, memOut, MemType::GPU>::ifAllocated;
     using FluxEsmBase<T, memIn, memOut, MemType::GPU>::allocated_size;
+    using FluxEsmBase<T, memIn, memOut, MemType::GPU>::model_param;
+    using FluxEsmBase<T, memIn, memOut, MemType::GPU>::numerics;
 public:
-    FluxEsm() = default;
-    FluxEsm(struct sfxDataVecTypeC* sfx,
-                struct meteoDataVecTypeC* meteo,
-                const struct sfx_esm_param model_param, 
-                const struct sfx_surface_param surface_param,
-                const struct sfx_esm_numericsTypeC numerics,
-                const struct sfx_phys_constants phys_constants,
-                const int grid_size) : FluxEsmBase<T, memIn, memOut, MemType::GPU>(sfx, meteo, 
-                model_param, surface_param, numerics, phys_constants, grid_size) {};
+    FluxEsm(sfxDataVecTypeC* sfx,
+                meteoDataVecTypeC* meteo,
+                const sfx_esm_param model_param, 
+                const sfx_surface_param surface_param,
+                const sfx_esm_numericsTypeC numerics,
+                const sfx_phys_constants phys_constants,
+                const int grid_size) : FluxEsmBase<T, memIn, memOut, MemType::GPU>(sfx, meteo, model_param, 
+                                       surface_param, numerics, phys_constants, grid_size) {}
     ~FluxEsm() = default;
     void compute_flux();
 };

srcC/sfx_call_cxx.c

@@ -15,33 +15,13 @@ void get_surface_fluxes_esm (struct sfxDataVecTypeC* sfx,
   surf_flux_esm_CXX(sfx, meteo, model_param, surface_param, numerics, constants, *grid_size); 
 }
 
-void get_surface_fluxes_sheba (float *zeta_, float *Rib_, float *Re_, float *B_, float *z0_m_, float *z0_t_, float *Rib_conv_lim_, float *Cm_, float *Ct_, float *Km_, float *Pr_t_inv_,
-    float *U_, float *dT_, float *Tsemi_, float *dQ_, float *h_, float *in_z0_m_,
-    const float *kappa, const float *Pr_t_0_inv, 
-    const float *alpha_m, const float *alpha_h, 
-    const float *a_m, const float *a_h, 
-    const float *b_m, const float *b_h,
-    const float *c_h,
-    const float *Re_rough_min, 
-    const float *B1_rough, const float *B2_rough,
-    const float *B_max_land, const float *B_max_ocean, const float *B_max_lake,
-    const float *gamma_c, const float *Re_visc_min,
-    const float *Pr_m, const float *nu_air, const float *g, 
-    const int *maxiters_charnock,
+void get_surface_fluxes_sheba (struct sfxDataVecTypeC* sfx,
+                                 struct meteoDataVecTypeC* meteo,
+                                 const struct sfx_sheba_param* model_param, 
+                                 const struct sfx_surface_param* surface_param,
+                                 const struct sfx_sheba_numericsTypeC* numerics,
+                                 const struct sfx_phys_constants* constants,
                                  const int *grid_size)
 {
-	surf_flux_sheba_CXX ( zeta_, Rib_, Re_, B_, z0_m_, z0_t_, Rib_conv_lim_, Cm_, Ct_, Km_, Pr_t_inv_,
-  U_, dT_, Tsemi_, dQ_, h_, in_z0_m_,
-  *kappa, *Pr_t_0_inv, 
-  *alpha_m, *alpha_h,
-  *a_m, *a_h, 
-  *b_m, *b_h,
-  *c_h,
-  *Re_rough_min, 
-  *B1_rough, *B2_rough,
-  *B_max_land, *B_max_ocean, *B_max_lake,
-  *gamma_c, *Re_visc_min,
-  *Pr_m, *nu_air, *g, 
-  *maxiters_charnock,
-  *grid_size);
+  surf_flux_sheba_CXX(sfx, meteo, model_param, surface_param, numerics, constants, *grid_size); 
 }
\ No newline at end of file

srcCXX/model_base.cpp

+#include "../includeCXX/model_base.h"
+
+#ifdef INCLUDE_CUDA
+    #include "../includeCU/sfx_memory_processing.cuh"
+#endif
+#include "../includeCXX/sfx_memory_processing.h"
+
+template<typename T, MemType memIn, MemType memOut, MemType RunMem >
+ModelBase<T, memIn, memOut, RunMem>::ModelBase(sfxDataVecTypeC* sfx_in,
+                meteoDataVecTypeC* meteo_in,
+                const sfx_surface_param surface_param_in,
+                const sfx_phys_constants phys_constants_in,
+                const int grid_size_in)
+{
+    ifAllocated = false;
+    grid_size = grid_size_in;
+    surface_param = surface_param_in;
+    phys_constants = phys_constants_in;
+
+    if(RunMem != memOut)
+        res_sfx = sfx_in;
+    else
+        res_sfx = nullptr;
+    
+    if(RunMem != memIn)
+    {
+        const size_t new_size = grid_size * sizeof(T);
+
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(meteo.U), allocated_size, new_size);
+        memproc::memcopy<RunMem, memIn>(meteo.U, meteo_in->U, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(meteo.dT), allocated_size, new_size);
+        memproc::memcopy<RunMem, memIn>(meteo.dT, meteo_in->dT, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(meteo.Tsemi), allocated_size, new_size);
+        memproc::memcopy<RunMem, memIn>(meteo.Tsemi, meteo_in->Tsemi, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(meteo.dQ), allocated_size, new_size);
+        memproc::memcopy<RunMem, memIn>(meteo.dQ, meteo_in->dQ, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(meteo.h), allocated_size, new_size);
+        memproc::memcopy<RunMem, memIn>(meteo.h, meteo_in->h, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(meteo.z0_m), allocated_size, new_size);
+        memproc::memcopy<RunMem, memIn>(meteo.z0_m, meteo_in->z0_m, new_size);
+
+        ifAllocated = true;
+    }
+    else
+    {
+        meteo = *meteo_in;
+    }
+
+    if(RunMem != memOut)
+    {
+        const size_t new_size = grid_size * sizeof(T);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.zeta), allocated_size, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.Rib), allocated_size, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.Re), allocated_size, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.Rib_conv_lim), allocated_size, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.z0_m), allocated_size, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.z0_t), allocated_size, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.B), allocated_size, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.Cm), allocated_size, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.Ct), allocated_size, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.Km), allocated_size, new_size);
+        allocated_size = 0;
+        memproc::realloc<RunMem>((void *&)(sfx.Pr_t_inv), allocated_size, new_size);
+
+        ifAllocated = true;
+    }
+    else
+    {
+        sfx = *sfx_in;
+    }
+}
+
+template<typename T, MemType memIn, MemType memOut, MemType RunMem >
+ModelBase<T, memIn, memOut, RunMem>::~ModelBase()
+{
+    if(ifAllocated == true)
+    {
+        if(RunMem != memIn)
+        {    
+            memproc::dealloc<RunMem>((void*&)meteo.U);
+            memproc::dealloc<RunMem>((void*&)meteo.dT);
+            memproc::dealloc<RunMem>((void*&)meteo.Tsemi);
+            memproc::dealloc<RunMem>((void*&)meteo.dQ);
+            memproc::dealloc<RunMem>((void*&)meteo.h);
+            memproc::dealloc<RunMem>((void*&)meteo.z0_m);
+        }
+        if(RunMem != memOut)
+        {
+            memproc::dealloc<RunMem>((void*&)sfx.zeta);
+            memproc::dealloc<RunMem>((void*&)sfx.Rib);
+            memproc::dealloc<RunMem>((void*&)sfx.Re);
+            memproc::dealloc<RunMem>((void*&)sfx.Rib_conv_lim);
+            memproc::dealloc<RunMem>((void*&)sfx.z0_m);
+            memproc::dealloc<RunMem>((void*&)sfx.z0_t);
+            memproc::dealloc<RunMem>((void*&)sfx.B);
+            memproc::dealloc<RunMem>((void*&)sfx.Cm);
+            memproc::dealloc<RunMem>((void*&)sfx.Ct);
+            memproc::dealloc<RunMem>((void*&)sfx.Km);
+            memproc::dealloc<RunMem>((void*&)sfx.Pr_t_inv);
+        }
+
+        ifAllocated = false;
+        allocated_size = 0;
+    }
+}
+
+template class ModelBase<float, MemType::CPU, MemType::CPU, MemType::CPU>;
+#ifdef INCLUDE_CUDA
+    template class ModelBase<float, MemType::GPU, MemType::GPU, MemType::GPU>;
+    template class ModelBase<float, MemType::GPU, MemType::GPU, MemType::CPU>;
+    template class ModelBase<float, MemType::GPU, MemType::CPU, MemType::GPU>;
+    template class ModelBase<float, MemType::CPU, MemType::GPU, MemType::GPU>;
+    template class ModelBase<float, MemType::CPU, MemType::CPU, MemType::GPU>;
+    template class ModelBase<float, MemType::CPU, MemType::GPU, MemType::CPU>;
+    template class ModelBase<float, MemType::GPU, MemType::CPU, MemType::CPU>;
+#endif 
\ No newline at end of file

srcCXX/sfx_call_class_func.cpp

@@ -6,12 +6,12 @@
 #include <vector>
 
 // -------------------------------------------------------------------------- //
-void surf_flux_esm_CXX (struct sfxDataVecTypeC* sfx,
-                           struct meteoDataVecTypeC* meteo,
-                           const struct sfx_esm_param* model_param, 
-                           const struct sfx_surface_param* surface_param,
-                           const struct sfx_esm_numericsTypeC* numerics,
-                           const struct sfx_phys_constants* constants,
+void surf_flux_esm_CXX (sfxDataVecTypeC* sfx,
+                           meteoDataVecTypeC* meteo,
+                           const sfx_esm_param* model_param, 
+                           const sfx_surface_param* surface_param,
+                           const sfx_esm_numericsTypeC* numerics,
+                           const sfx_phys_constants* constants,
                            const int grid_size)
 {
 #ifdef INCLUDE_CUDA
@@ -23,38 +23,19 @@ void surf_flux_esm_CXX (struct sfxDataVecTypeC* sfx,
 #endif
 }
 
-void surf_flux_sheba_CXX (float *zeta_, float *Rib_, float *Re_, float *B_, float *z0_m_, float *z0_t_, float *Rib_conv_lim_, float *Cm_, float *Ct_, float *Km_, float *Pr_t_inv_,
-    float *U_, float *dT_, float *Tsemi_, float *dQ_, float *h_, float *in_z0_m_,
-    const float kappa, const float Pr_t_0_inv, 
-    const float alpha_m, const float alpha_h, 
-    const float a_m, const float a_h, 
-    const float b_m, const float b_h,
-    const float c_h,
-    const float Re_rough_min, 
-    const float B1_rough, const float B2_rough,
-    const float B_max_land, const float B_max_ocean, const float B_max_lake,
-    const float gamma_c, const float Re_visc_min,
-    const float Pr_m, const float nu_air, const float g, 
-    const int maxiters_charnock, 
+void surf_flux_sheba_CXX (sfxDataVecTypeC* sfx,
+                           meteoDataVecTypeC* meteo,
+                           const sfx_sheba_param* model_param, 
+                           const sfx_surface_param* surface_param,
+                           const sfx_sheba_numericsTypeC* numerics,
+                           const sfx_phys_constants* constants,
                            const int grid_size)
 {
-#ifdef INCLUDE_CUDA
-    static FluxSheba<float, MemType::CPU, MemType::CPU, MemType::GPU> F;
-#else
-    static FluxSheba<float, MemType::CPU, MemType::CPU, MemType::CPU> F;
-#endif
-    F.set_params(grid_size, kappa, Pr_t_0_inv, 
-    alpha_m, alpha_h, 
-    a_m, a_h, 
-    b_m, b_h,
-    c_h,
-    Re_rough_min, 
-    B1_rough, B2_rough,
-    B_max_land, B_max_ocean, B_max_lake,
-    gamma_c, Re_visc_min,
-    Pr_m, nu_air, g);
-
-    F.compute_flux_sheba(zeta_, Rib_, Re_, B_, z0_m_, z0_t_, Rib_conv_lim_, Cm_, Ct_, Km_, Pr_t_inv_,
-    U_, dT_, Tsemi_, dQ_, h_, in_z0_m_,
-    maxiters_charnock);
+// #ifdef INCLUDE_CUDA
+//     static FluxSheba<float, MemType::CPU, MemType::CPU, MemType::GPU> F(sfx, meteo, *model_param, *surface_param, *numerics, *constants, grid_size);
+//     F.compute_flux();
+// #else
+//     static FluxSheba<float, MemType::CPU, MemType::CPU, MemType::CPU> F(sfx, meteo, *model_param, *surface_param, *numerics, *constants, grid_size);
+//     F.compute_flux();
+// #endif
 }
\ No newline at end of file

srcCXX/sfx_esm.cpp

@@ -11,124 +11,21 @@
 #include "../includeCU/sfx_esm_compute_subfunc.cuh"
 
 template<typename T, MemType memIn, MemType memOut, MemType RunMem >
-FluxEsmBase<T, memIn, memOut, RunMem>::FluxEsmBase(struct sfxDataVecTypeC* sfx_in,
-                struct meteoDataVecTypeC* meteo_in,
-                const struct sfx_esm_param model_param_in, 
-                const struct sfx_surface_param surface_param_in,
-                const struct sfx_esm_numericsTypeC numerics_in,
-                const struct sfx_phys_constants phys_constants_in,
-                const int grid_size_in)
+FluxEsmBase<T, memIn, memOut, RunMem>::FluxEsmBase(sfxDataVecTypeC* sfx_in,
+                meteoDataVecTypeC* meteo_in,
+                const sfx_esm_param model_param_in, 
+                const sfx_surface_param surface_param_in,
+                const sfx_esm_numericsTypeC numerics_in,
+                const sfx_phys_constants phys_constants_in,
+                const int grid_size_in) : ModelBase<T, memIn, memOut, RunMem>(sfx_in, meteo_in, 
+                    surface_param_in, phys_constants_in, grid_size_in)
 {
-    ifAllocated = false;
-    grid_size = grid_size_in;
-
     model_param = model_param_in;
-    surface_param = surface_param_in;
     numerics = numerics_in;
-    phys_constants = phys_constants_in;
-
-    if(RunMem != memOut)
-        res_sfx = sfx_in;
-    else
-        res_sfx = nullptr;
-    
-    if(RunMem != memIn)
-    {
-        const size_t new_size = grid_size * sizeof(T);
-
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(meteo.U), allocated_size, new_size);
-        memproc::memcopy<RunMem, memIn>(meteo.U, meteo_in->U, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(meteo.dT), allocated_size, new_size);
-        memproc::memcopy<RunMem, memIn>(meteo.dT, meteo_in->dT, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(meteo.Tsemi), allocated_size, new_size);
-        memproc::memcopy<RunMem, memIn>(meteo.Tsemi, meteo_in->Tsemi, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(meteo.dQ), allocated_size, new_size);
-        memproc::memcopy<RunMem, memIn>(meteo.dQ, meteo_in->dQ, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(meteo.h), allocated_size, new_size);
-        memproc::memcopy<RunMem, memIn>(meteo.h, meteo_in->h, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(meteo.z0_m), allocated_size, new_size);
-        memproc::memcopy<RunMem, memIn>(meteo.z0_m, meteo_in->z0_m, new_size);
-
-        ifAllocated = true;
-    }
-    else
-    {
-        meteo = *meteo_in;
-    }
-
-    if(RunMem != memOut)
-    {
-        const size_t new_size = grid_size * sizeof(T);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.zeta), allocated_size, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.Rib), allocated_size, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.Re), allocated_size, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.Rib_conv_lim), allocated_size, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.z0_m), allocated_size, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.z0_t), allocated_size, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.B), allocated_size, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.Cm), allocated_size, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.Ct), allocated_size, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.Km), allocated_size, new_size);
-        allocated_size = 0;
-        memproc::realloc<RunMem>((void *&)(sfx.Pr_t_inv), allocated_size, new_size);
-
-        ifAllocated = true;
-    }
-    else
-    {
-        sfx = *sfx_in;
-    }
 }
 
 template<typename T, MemType memIn, MemType memOut, MemType RunMem >
-FluxEsmBase<T, memIn, memOut, RunMem>::~FluxEsmBase()
-{
-    if(ifAllocated == true)
-    {
-        if(RunMem != memIn)
-        {    
-            memproc::dealloc<RunMem>((void*&)meteo.U);
-            memproc::dealloc<RunMem>((void*&)meteo.dT);
-            memproc::dealloc<RunMem>((void*&)meteo.Tsemi);
-            memproc::dealloc<RunMem>((void*&)meteo.dQ);
-            memproc::dealloc<RunMem>((void*&)meteo.h);
-            memproc::dealloc<RunMem>((void*&)meteo.z0_m);
-        }
-        if(RunMem != memOut)
-        {
-            memproc::dealloc<RunMem>((void*&)sfx.zeta);
-            memproc::dealloc<RunMem>((void*&)sfx.Rib);
-            memproc::dealloc<RunMem>((void*&)sfx.Re);
-            memproc::dealloc<RunMem>((void*&)sfx.Rib_conv_lim);
-            memproc::dealloc<RunMem>((void*&)sfx.z0_m);
-            memproc::dealloc<RunMem>((void*&)sfx.z0_t);
-            memproc::dealloc<RunMem>((void*&)sfx.B);
-            memproc::dealloc<RunMem>((void*&)sfx.Cm);
-            memproc::dealloc<RunMem>((void*&)sfx.Ct);
-            memproc::dealloc<RunMem>((void*&)sfx.Km);
-            memproc::dealloc<RunMem>((void*&)sfx.Pr_t_inv);
-        }
-
-        ifAllocated = false;
-        allocated_size = 0;
-    }
-}
+FluxEsmBase<T, memIn, memOut, RunMem>::~FluxEsmBase() {}
 
 template<typename T, MemType memIn, MemType memOut >
 void FluxEsm<T, memIn, memOut, MemType::CPU>::compute_flux()
@@ -258,12 +155,4 @@ template class FluxEsmBase<float, MemType::CPU, MemType::CPU, MemType::CPU>;
     template class FluxEsm<float, MemType::CPU, MemType::CPU, MemType::GPU>;
     template class FluxEsm<float, MemType::CPU, MemType::GPU, MemType::CPU>;
     template class FluxEsm<float, MemType::GPU, MemType::CPU, MemType::CPU>;
-
-    // template class FluxEsm<double, MemType::GPU, MemType::GPU, MemType::GPU>;
-    // template class FluxEsm<double, MemType::GPU, MemType::GPU, MemType::CPU>;
-    // template class FluxEsm<double, MemType::GPU, MemType::CPU, MemType::GPU>;
-    // template class FluxEsm<double, MemType::CPU, MemType::GPU, MemType::GPU>;
-    // template class FluxEsm<double, MemType::CPU, MemType::CPU, MemType::GPU>;
-    // template class FluxEsm<double, MemType::CPU, MemType::GPU, MemType::CPU>;
-    // template class FluxEsm<double, MemType::GPU, MemType::CPU, MemType::CPU>;
 #endif 
\ No newline at end of file

srcF/sfx_data.f90

@@ -132,19 +132,6 @@ module sfx_data
         type(C_PTR) :: Pr_t_inv        !< inverse turbulent Prandtl number at h [n/d]
     end type
 
-    type, BIND(C), public :: sfx_esm_param 
-        real(C_FLOAT) :: kappa           
-        real(C_FLOAT) :: Pr_t_0_inv
-        real(C_FLOAT) :: Pr_t_inf_inv
-
-        real(C_FLOAT) :: alpha_m           
-        real(C_FLOAT) :: alpha_h
-        real(C_FLOAT) :: alpha_h_fix
-        real(C_FLOAT) :: beta_m
-        real(C_FLOAT) :: beta_h
-        real(C_FLOAT) :: Rib_max
-    end type
-
     type, BIND(C), public :: sfx_surface_param 
         integer(C_INT) :: surface_ocean           
         integer(C_INT) :: surface_land
@@ -165,15 +152,47 @@ module sfx_data
         real(C_FLOAT)  :: B_max_land;
     end type
 
+    type, BIND(C), public :: sfx_phys_constants 
+        real(C_FLOAT)  :: Pr_m;
+        real(C_FLOAT)  :: g;
+        real(C_FLOAT)  :: nu_air;
+    end type
+
+    ! ESM strcuctures
+    type, BIND(C), public :: sfx_esm_param 
+        real(C_FLOAT) :: kappa           
+        real(C_FLOAT) :: Pr_t_0_inv
+        real(C_FLOAT) :: Pr_t_inf_inv
+
+        real(C_FLOAT) :: alpha_m           
+        real(C_FLOAT) :: alpha_h
+        real(C_FLOAT) :: alpha_h_fix
+        real(C_FLOAT) :: beta_m
+        real(C_FLOAT) :: beta_h
+        real(C_FLOAT) :: Rib_max
+    end type
+
     type, BIND(C), public :: sfx_esm_numericsTypeC 
         integer(C_INT) :: maxiters_convection           
         integer(C_INT) :: maxiters_charnock 
     end type
 
-    type, BIND(C), public :: sfx_phys_constants 
-        real(C_FLOAT)  :: Pr_m;
-        real(C_FLOAT)  :: g;
-        real(C_FLOAT)  :: nu_air;
+    ! SHEBA strcuctures
+    type, BIND(C), public :: sfx_sheba_param 
+        real(C_FLOAT) :: kappa           
+        real(C_FLOAT) :: Pr_t_0_inv
+
+        real(C_FLOAT) :: alpha_m           
+        real(C_FLOAT) :: alpha_h
+        real(C_FLOAT) :: a_m
+        real(C_FLOAT) :: b_m
+        real(C_FLOAT) :: a_h
+        real(C_FLOAT) :: b_h
+        real(C_FLOAT) :: c_h
+    end type
+
+    type, BIND(C), public :: sfx_sheba_numericsTypeC 
+        integer(C_INT) :: maxiters_charnock 
     end type
 #endif
     ! --------------------------------------------------------------------------------

srcF/sfx_fc_wrapper.F90

@@ -16,30 +16,17 @@ module C_FUNC
             type(sfx_phys_constants) :: constants
         END SUBROUTINE get_surface_fluxes_esm
 
-        SUBROUTINE get_surface_fluxes_sheba(zeta, Rib, Re, B, z0_m, z0_t, Rib_conv_lim, &
-            Cm, Ct, Km, Prt_inv, & 
-            U, dT, Tsemi, dQ, h, in_z0_m, & 
-            kappa, Pr_t_0_inv, & 
-            alpha_m, alpha_h, & 
-            a_m, a_h,         &
-            b_m, b_h,         &
-            c_h,              &
-            Re_rough_min, & 
-            B1_rough, B2_rough, & 
-            B_max_land, B_max_ocean, B_max_lake, & 
-            gamma_c, Re_visc_min, & 
-            Pr_m, nu_air, g, & 
-            maxiters_charnock, & 
-            grid_size) BIND(C)
-            USE, INTRINSIC :: ISO_C_BINDING, ONLY: C_INT
-            USE, INTRINSIC :: ISO_C_BINDING, ONLY: C_FLOAT
+        SUBROUTINE get_surface_fluxes_sheba(sfx, meteo, model_param, surface_param, numerics, constants, grid_size) BIND(C)
+            use sfx_data
+            USE, INTRINSIC :: ISO_C_BINDING, ONLY: C_INT, C_PTR
             IMPLICIT NONE
-            INTEGER(C_INT) :: grid_size, maxiters_charnock
-            REAL(C_FLOAT)  :: kappa, Pr_t_0_inv, Pr_t_inf_inv, alpha_m, alpha_h, a_m, a_h, b_m, b_h, &
-            c_h, Rib_max, Re_rough_min, B1_rough, B2_rough, B_max_land, B_max_ocean, & 
-            B_max_lake, gamma_c, Re_visc_min, Pr_m, nu_air, g
-            REAL(C_FLOAT), dimension(grid_size) :: U, dT, Tsemi, dQ, h, in_z0_m, zeta, Rib, Re, &
-            Rib_conv_lim, z0_m, z0_t, B, Cm, Ct, Km, Prt_inv
+            INTEGER(C_INT) :: grid_size
+            type(C_PTR), value :: sfx
+            type(C_PTR), value :: meteo
+            type(sfx_sheba_param) :: model_param
+            type(sfx_surface_param) :: surface_param
+            type(sfx_sheba_numericsTypeC) :: numerics
+            type(sfx_phys_constants) :: constants
         END SUBROUTINE get_surface_fluxes_sheba
     END INTERFACE
 

srcF/sfx_sheba.f90

@@ -13,6 +13,7 @@ module sfx_sheba
     use sfx_sheba_param
 
 #if defined(INCLUDE_CUDA) || defined(INCLUDE_CXX)
+    use iso_c_binding, only: c_loc
     use C_FUNC
 #endif
     ! --------------------------------------------------------------------------------
@@ -37,6 +38,22 @@ module sfx_sheba
 
 contains
 
+#if defined(INCLUDE_CXX)
+    subroutine set_c_struct_sfx_sheba_param_values(sfx_model_param)
+        type (sfx_sheba_param), intent(inout) :: sfx_model_param
+        sfx_model_param%kappa = kappa
+        sfx_model_param%Pr_t_0_inv = Pr_t_0_inv
+
+        sfx_model_param%alpha_m = alpha_m
+        sfx_model_param%alpha_h = alpha_h
+        sfx_model_param%a_m = a_m
+        sfx_model_param%b_m = b_m
+        sfx_model_param%a_h = a_h
+        sfx_model_param%b_h = b_h
+        sfx_model_param%c_h = c_h
+    end subroutine set_c_struct_sfx_sheba_param_values
+#endif
+
     ! --------------------------------------------------------------------------------
     subroutine get_surface_fluxes_vec(sfx, meteo, numerics, n)
         !< @brief surface flux calculation for array data
@@ -54,22 +71,29 @@ contains
         type (sfxDataType) sfx_cell
         integer i
         ! ----------------------------------------------------------------------------
-#if defined(INCLUDE_CUDA) || defined(INCLUDE_CXX)
-        call get_surface_fluxes_sheba(sfx%zeta, sfx%Rib, sfx%Re, sfx%B, sfx%z0_m, sfx%z0_t, &
-        sfx%Rib_conv_lim, sfx%Cm, sfx%Ct, sfx%Km, sfx%Pr_t_inv, &
-        meteo%U, meteo%dT, meteo%Tsemi, meteo%dQ, meteo%h, meteo%z0_m, &
-        kappa, Pr_t_0_inv, & 
-        alpha_m, alpha_h, & 
-        a_m, a_h,         &
-        b_m, b_h,         &
-        c_h,              &
-        Re_rough_min, & 
-        B1_rough, B2_rough, & 
-        B_max_land, B_max_ocean, B_max_lake, & 
-        gamma_c, Re_visc_min, & 
-        Pr_m, nu_air, g, & 
-        numerics%maxiters_charnock, & 
-        n)
+#if defined(INCLUDE_CXX)
+        type (meteoDataVecTypeC), pointer :: meteo_c         !< meteorological data (input)
+        type (sfxDataVecTypeC), pointer :: sfx_c             !< surface flux data (output)
+        type (sfx_sheba_param) :: model_param
+        type (sfx_surface_param) :: surface_param
+        type (sfx_sheba_numericsTypeC) :: numerics_c
+        type (sfx_phys_constants) :: phys_constants
+
+        numerics_c%maxiters_charnock = numerics%maxiters_charnock
+
+        phys_constants%Pr_m = Pr_m;
+        phys_constants%nu_air = nu_air;
+        phys_constants%g = g;
+
+        call set_c_struct_sfx_sheba_param_values(model_param)
+        call set_c_struct_sfx_surface_param_values(surface_param)
+        call set_meteo_vec_c(meteo, meteo_c)
+        call set_sfx_vec_c(sfx, sfx_c)
+
+        ! call get_surface_fluxes_sheba(c_loc(sfx_c), c_loc(meteo_c), model_param, surface_param, numerics_c, phys_constants, n)
+
+        deallocate(meteo_c)
+        deallocate(sfx_c)
 #else
         do i = 1, n
             meteo_cell = meteoDataType(&