Enable Kokkos OpenACC build

Author: mhoemmen

cpp-mdspan/CMakeLists.txt

@@ -106,6 +106,9 @@ endif()
 # Users can set FETCHCONTENT_SOURCE_DIR_KOKKOS
 # to a local path to Kokkos, to override automatic downloading.
 
+set(Kokkos_ENABLE_SERIAL ON CACHE INTERNAL "")
+set(Kokkos_ENABLE_OPENACC ON CACHE INTERNAL "")
+
 FetchContent_Declare(
   Kokkos
   GIT_REPOSITORY https://github.com/kokkos/kokkos.git
@@ -181,13 +184,36 @@ target_compile_options(miniWeather_serial PRIVATE
 )
 
 if (kokkos_POPULATED)
-  add_executable(miniWeather_kokkos miniWeather_mdspan.cpp miniWeather_common.cpp)
-  target_include_directories(miniWeather_kokkos PRIVATE "${PROJECT_SOURCE_DIR}")
+  add_executable(miniWeather_kokkos_serial miniWeather_mdspan.cpp miniWeather_common.cpp)
+  target_include_directories(miniWeather_kokkos_serial PRIVATE "${PROJECT_SOURCE_DIR}")
+
+  target_link_libraries(miniWeather_kokkos_serial PRIVATE Kokkos::kokkos)
+  target_compile_definitions(miniWeather_kokkos_serial PRIVATE MINIWEATHER_KOKKOS)
+  target_compile_definitions(miniWeather_kokkos_serial PRIVATE MINIWEATHER_KOKKOS_SERIAL)
+
+  # We got mdspan from Kokkos.
+  target_compile_options(miniWeather_kokkos_serial PRIVATE
+    $<$<OR:$<CXX_COMPILER_ID:Clang>,$<CXX_COMPILER_ID:AppleClang>,$<CXX_COMPILER_ID:GNU>>:
+      -Wall>
+    $<$<CXX_COMPILER_ID:MSVC>:
+      /W4>
+  )
+endif()
+
+# FIXME (mfh 2025/04/04) Figure out how to pass OpenACC flags to Kokkos.
+if (kokkos_POPULATED AND OpenACC_FOUND)
+  set(Kokkos_ENABLE_OPENACC ON)
+
+  add_executable(miniWeather_kokkos_openacc miniWeather_mdspan.cpp miniWeather_common.cpp)
+  target_include_directories(miniWeather_kokkos_openacc PRIVATE "${PROJECT_SOURCE_DIR}")
+
+  target_link_libraries(miniWeather_kokkos_openacc PRIVATE Kokkos::kokkos)
+  target_compile_definitions(miniWeather_kokkos_openacc PRIVATE MINIWEATHER_KOKKOS)
+  target_compile_definitions(miniWeather_kokkos_openacc PRIVATE MINIWEATHER_KOKKOS_OPENACC)
 
-  target_link_libraries(miniWeather_kokkos PRIVATE Kokkos::kokkos)
-  target_compile_definitions(miniWeather_kokkos PRIVATE MINIWEATHER_KOKKOS)
   # We got mdspan from Kokkos.
-  target_compile_options(miniWeather_kokkos PRIVATE
+  target_link_libraries(miniWeather_kokkos_openacc PRIVATE OpenACC::OpenACC_CXX)
+  target_compile_options(miniWeather_kokkos_openacc PRIVATE
     $<$<OR:$<CXX_COMPILER_ID:Clang>,$<CXX_COMPILER_ID:AppleClang>,$<CXX_COMPILER_ID:GNU>>:
       -Wall>
     $<$<CXX_COMPILER_ID:MSVC>:

cpp-mdspan/miniWeather_common.cpp

@@ -10,102 +10,9 @@ make_unique_array_1d(host_memory_space, int X) {
   return std::make_unique<double[]>(X);
 }
 
-test_case injection(double x , double z) {
-  auto [hr, ht] = hydro_const_theta(z);
-  double r = 0.0;
-  double t = 0.0;
-  double u = 0.0;
-  double w = 0.0;
-  return {r, u, w, t, hr, ht};
-}
-
-test_case density_current(double x , double z) {
-  auto [hr, ht] = hydro_const_theta(z);
-  double r = 0.0;
-  double t = sample_ellipse_cosine(x, z, -20.0, xlen/2, 5000.0, 4000.0, 2000.0);
-  double u = 0.0;
-  double w = 0.0;
-  return {r, u, w, t, hr, ht};
-}
-
-test_case gravity_waves(double x, double z) {
-  auto [hr, ht] = hydro_const_bvfreq(z, 0.02);
-  double r = 0.0;
-  double t = 0.0;
-  double u = 15.0;
-  double w = 0.0;
-  return {r, u, w, t, hr, ht};
-}
-
-test_case thermal(double x, double z) {
-  auto [hr, ht] = hydro_const_theta(z);
-  double r = 0.0;
-  double t = sample_ellipse_cosine(x, z, 3.0, xlen/2,2000.0, 2000.0, 2000.0);
-  double u = 0.0;
-  double w = 0.0;
-  return {r, u, w, t, hr, ht};
-}
-
-test_case collision(double x , double z) {
-  auto [hr, ht] = hydro_const_theta(z);
-  double r = 0.0;
-  double t = 0.0;
-  double u = 0.0;
-  double w = 0.0;
-  t = t + sample_ellipse_cosine(x, z,  20.0, xlen/2,2000.0, 2000.0, 2000.0);
-  t = t + sample_ellipse_cosine(x, z, -20.0, xlen/2,8000.0, 2000.0, 2000.0);
-  return {r, u, w, t, hr, ht};
-}
-
-test_case get_test_case(int data_spec, double x_, double z_) {
-  if (data_spec == DATA_SPEC_COLLISION      ) { return collision(x_, z_); }
-  if (data_spec == DATA_SPEC_THERMAL        ) { return thermal(x_, z_); }
-  if (data_spec == DATA_SPEC_GRAVITY_WAVES  ) { return gravity_waves(x_, z_); }
-  if (data_spec == DATA_SPEC_DENSITY_CURRENT) { return density_current(x_, z_); }
-  if (data_spec == DATA_SPEC_INJECTION      ) { return injection(x_, z_); }
-  assert(false);
-  return test_case{};
-}
-
-r_t_pair hydro_const_theta(double z) {
-  const double theta0 = 300.;  //Background potential temperature
-  const double exner0 = 1.;    //Surface-level Exner pressure
-  double       p,exner,rt;
-  //Establish hydrostatic balance first using Exner pressure
-  double t = theta0;                         //Potential Temperature at z
-  exner = exner0 - grav * z / (cp * theta0); //Exner pressure at z
-  p = p0 * pow(exner,(cp/rd));               //Pressure at z
-  rt = pow((p / C0),(1. / gamm));            //rho*theta at z
-  double r = rt / t;                         //Density at z
-
-  return {r, t};
-}
-
-r_t_pair hydro_const_bvfreq(double z, double bv_freq0) {
-  const double theta0 = 300.;  //Background potential temperature
-  const double exner0 = 1.;    //Surface-level Exner pressure
-  double       p, exner, rt;
-  double t = theta0 * exp( bv_freq0*bv_freq0 / grav * z );                                    //Pot temp at z
-  exner = exner0 - grav*grav / (cp * bv_freq0*bv_freq0) * (t - theta0) / (t * theta0); //Exner pressure at z
-  p = p0 * pow(exner,(cp/rd));                                                         //Pressure at z
-  rt = pow((p / C0), (1. / gamm));                                                  //rho*theta at z
-  double r = rt / t;                                                                          //Density at z
-
-  return {r, t};
-}
-
-double sample_ellipse_cosine( double x , double z , double amp , double x0 , double z0 , double xrad , double zrad ) {
-  double dist;
-  //Compute distance from bubble center
-  dist = sqrt( ((x-x0)/xrad)*((x-x0)/xrad) + ((z-z0)/zrad)*((z-z0)/zrad) ) * pi / 2.0;
-  //If the distance from bubble center is less than the radius, create a cos**2 profile
-  if (dist <= pi / 2.0) {
-    return amp * pow(cos(dist), 2.0);
-  } else {
-    return 0.;
-  }
-}
-
 void finalize() {
+#if defined(MINIWEATHER_KOKKOS)
+  Kokkos::finalize();
+#endif
   (void) MPI_Finalize();
 }

cpp-mdspan/miniWeather_common.hpp

@@ -20,6 +20,13 @@
 #include "mdspan/mdspan.hpp"
 #include "unique_mdarray.hpp"
 
+#if defined(MINIWEATHER_KOKKOS)
+#  include "Kokkos_Core.hpp"
+#  define MINIWEATHER_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
+#else
+#  define MINIWEATHER_INLINE_FUNCTION inline
+#endif
+
 constexpr double pi        = 3.14159265358979323846264338327;   //Pi
 constexpr double grav      = 9.8;                               //Gravitational acceleration (m / s^2)
 constexpr double cp        = 1004.;                             //Specific heat of dry air at constant pressure
@@ -251,6 +258,64 @@ struct init_result {
   global_arrays<MemorySpace> arrays;
 };
 
+struct r_t_pair {
+  double r;
+  double t;
+};
+
+// Establish hydrostatic balance using constant potential temperature
+// (thermally neutral atmosphere)
+// z is the input coordinate
+// r and t are the output background hydrostatic density and potential temperature
+MINIWEATHER_INLINE_FUNCTION
+r_t_pair hydro_const_theta(double z) {
+  const double theta0 = 300.;  //Background potential temperature
+  const double exner0 = 1.;    //Surface-level Exner pressure
+  double       p,exner,rt;
+  //Establish hydrostatic balance first using Exner pressure
+  double t = theta0;                         //Potential Temperature at z
+  exner = exner0 - grav * z / (cp * theta0); //Exner pressure at z
+  p = p0 * pow(exner,(cp/rd));               //Pressure at z
+  rt = pow((p / C0),(1. / gamm));            //rho*theta at z
+  double r = rt / t;                         //Density at z
+
+  return {r, t};
+}
+
+//Establish hydrostatic balance using constant Brunt-Vaisala frequency
+//z is the input coordinate
+//bv_freq0 is the constant Brunt-Vaisala frequency
+//r and t are the output background hydrostatic density and potential temperature
+MINIWEATHER_INLINE_FUNCTION
+r_t_pair hydro_const_bvfreq(double z, double bv_freq0) {
+  const double theta0 = 300.;  //Background potential temperature
+  const double exner0 = 1.;    //Surface-level Exner pressure
+  double       p, exner, rt;
+  double t = theta0 * exp( bv_freq0*bv_freq0 / grav * z );                                    //Pot temp at z
+  exner = exner0 - grav*grav / (cp * bv_freq0*bv_freq0) * (t - theta0) / (t * theta0); //Exner pressure at z
+  p = p0 * pow(exner,(cp/rd));                                                         //Pressure at z
+  rt = pow((p / C0), (1. / gamm));                                                  //rho*theta at z
+  double r = rt / t;                                                                          //Density at z
+
+  return {r, t};
+}
+
+//Sample from an ellipse of a specified center, radius, and amplitude at a specified location
+//x and z are input coordinates
+//amp,x0,z0,xrad,zrad are input amplitude, center, and radius of the ellipse
+MINIWEATHER_INLINE_FUNCTION
+double sample_ellipse_cosine( double x , double z , double amp , double x0 , double z0 , double xrad , double zrad ) {
+  double dist;
+  //Compute distance from bubble center
+  dist = sqrt( ((x-x0)/xrad)*((x-x0)/xrad) + ((z-z0)/zrad)*((z-z0)/zrad) ) * pi / 2.0;
+  //If the distance from bubble center is less than the radius, create a cos**2 profile
+  if (dist <= pi / 2.0) {
+    return amp * pow(cos(dist), 2.0);
+  } else {
+    return 0.;
+  }
+}
+
 struct test_case {
   double r;
   double u;
@@ -265,43 +330,71 @@ struct test_case {
 // hr and ht are output background hydrostatic density and potential temperature at that location.
 
 //This test case is initially balanced but injects fast, cold air from the left boundary near the model top
-test_case injection(double x, double z);
+MINIWEATHER_INLINE_FUNCTION
+test_case injection(double x , double z) {
+  auto [hr, ht] = hydro_const_theta(z);
+  double r = 0.0;
+  double t = 0.0;
+  double u = 0.0;
+  double w = 0.0;
+  return {r, u, w, t, hr, ht};
+}
 
 //Initialize a density current (falling cold thermal that propagates along the model bottom)
-test_case density_current(double x, double z);
-
-test_case gravity_waves(double x, double z);
+MINIWEATHER_INLINE_FUNCTION
+test_case density_current(double x , double z) {
+  auto [hr, ht] = hydro_const_theta(z);
+  double r = 0.0;
+  double t = sample_ellipse_cosine(x, z, -20.0, xlen/2, 5000.0, 4000.0, 2000.0);
+  double u = 0.0;
+  double w = 0.0;
+  return {r, u, w, t, hr, ht};
+}
+
+MINIWEATHER_INLINE_FUNCTION
+test_case gravity_waves(double x, double z) {
+  auto [hr, ht] = hydro_const_bvfreq(z, 0.02);
+  double r = 0.0;
+  double t = 0.0;
+  double u = 15.0;
+  double w = 0.0;
+  return {r, u, w, t, hr, ht};
+}
 
 //Rising thermal
-test_case thermal(double x, double z);
+MINIWEATHER_INLINE_FUNCTION
+test_case thermal(double x, double z) {
+  auto [hr, ht] = hydro_const_theta(z);
+  double r = 0.0;
+  double t = sample_ellipse_cosine(x, z, 3.0, xlen/2,2000.0, 2000.0, 2000.0);
+  double u = 0.0;
+  double w = 0.0;
+  return {r, u, w, t, hr, ht};
+}
 
 //Colliding thermals
-test_case collision(double x, double z);
-
-test_case get_test_case(int data_spec, double x_, double z_);
-
-struct r_t_pair {
-  double r;
-  double t;
-};
-
-// Establish hydrostatic balance using constant potential temperature
-// (thermally neutral atmosphere)
-// z is the input coordinate
-// r and t are the output background hydrostatic density and potential temperature
-r_t_pair hydro_const_theta(double z);
-
-//Establish hydrostatic balance using constant Brunt-Vaisala frequency
-//z is the input coordinate
-//bv_freq0 is the constant Brunt-Vaisala frequency
-//r and t are the output background hydrostatic density and potential temperature
-r_t_pair hydro_const_bvfreq(double z, double bv_freq0);
-
-//Sample from an ellipse of a specified center, radius, and amplitude at a specified location
-//x and z are input coordinates
-//amp,x0,z0,xrad,zrad are input amplitude, center, and radius of the ellipse
-double sample_ellipse_cosine(double x, double z, double amp, double x0, double z0,
-                             double xrad, double zrad);
+MINIWEATHER_INLINE_FUNCTION
+test_case collision(double x , double z) {
+  auto [hr, ht] = hydro_const_theta(z);
+  double r = 0.0;
+  double t = 0.0;
+  double u = 0.0;
+  double w = 0.0;
+  t = t + sample_ellipse_cosine(x, z,  20.0, xlen/2,2000.0, 2000.0, 2000.0);
+  t = t + sample_ellipse_cosine(x, z, -20.0, xlen/2,8000.0, 2000.0, 2000.0);
+  return {r, u, w, t, hr, ht};
+}
+
+MINIWEATHER_INLINE_FUNCTION
+test_case get_test_case(int data_spec, double x, double z) {
+  if (data_spec == DATA_SPEC_COLLISION      ) { return collision(x, z); }
+  if (data_spec == DATA_SPEC_THERMAL        ) { return thermal(x, z); }
+  if (data_spec == DATA_SPEC_GRAVITY_WAVES  ) { return gravity_waves(x, z); }
+  if (data_spec == DATA_SPEC_DENSITY_CURRENT) { return density_current(x, z); }
+  if (data_spec == DATA_SPEC_INJECTION      ) { return injection(x, z); }
+  assert(false);
+  return test_case{};
+}
 
 struct reduction_result {
   double mass;

cpp-mdspan/miniWeather_generic_algs.hpp

@@ -1,6 +1,9 @@
 #pragma once
 
 #include "miniWeather_common.hpp"
+#if defined(MINIWEATHER_KOKKOS)
+#  include "Kokkos_Core.hpp"
+#endif
 
 // Perform a single time step.
 // Time steps are dimensionally split and
@@ -87,7 +90,9 @@ init_result<MemorySpace> init(
   int *argc , char ***argv)
 {
   (void) MPI_Init(argc,argv);
-
+#if defined(MINIWEATHER_KOKKOS)
+  Kokkos::initialize(*argc, *argv);
+#endif
   /////////////////////////////////////////////////////////////
   // BEGIN MPI DUMMY SECTION
   // TODO: (1) GET NUMBER OF MPI RANKS
@@ -205,6 +210,3 @@ reduction_result reductions(
     .te = glob[1]
   };
 }
-
-
-