Homme(xx)/SL: Finish C++/Kokkos for ETM; modify vertical discretization.

The atmosphere semi-Lagrangian tracer transport method's enhanced trajectory
method (ETM) is a stealth feature intended to make tracer transport time steps
even longer and, in any case, more flexible.

This PR completes the C++/Kokkos version of the ETM, in particular support for
semi_lagrange_trajectory_nvelocity > 2.

It also modifies the vertical discretization relative to the original ETM
implementation, in light of issues at and around the tropopause when using the
ETM in test v3 simulations. This new approach resolves the tropopause issue.

Finally, this PR strengthens the homme_integration test suite as follows:
* add a new planar tracer transport test with analytical solution at the end
  time, with space- and time-varying surface pressure;
* add instances of this test to the homme_integration test suite;
* add a cmake tool to check error output for the recently added sphere and new
  plane tracer transport tests when running this suite;
* activate the Kokkos builds in the non-BFB part of homme_integration, turning
  on the parts of the unit tests that are meaningful in the non-BFB build and
  the recently added and new tracer transport tests.

Author: ambrad

components/eamxx/cime_config/namelist_defaults_eamxx.xml

@@ -897,6 +897,9 @@ be lost if SCREAM_HACK_XML is not enabled.
     <hypervis_subcycle_q hgrid=".*pg2">6</hypervis_subcycle_q>
     <transport_alg hgrid=".*pg2">12</transport_alg>
     <semi_lagrange_trajectory_nsubstep>0</semi_lagrange_trajectory_nsubstep>
+    <semi_lagrange_trajectory_nvelocity>-1</semi_lagrange_trajectory_nvelocity>
+    <semi_lagrange_halo>-1</semi_lagrange_halo>
+    <semi_lagrange_diagnostics>0</semi_lagrange_diagnostics>
     <!-- Other settings that we'll trigger based on pg2 for convenience -->
     <se_ftype valid_values="0,2" hgrid=".*pg2">2</se_ftype>
     <mesh_file type="file">none</mesh_file>

components/eamxx/cime_config/testdefs/testmods_dirs/eamxx/sl_nsubstep2/shell_commands

@@ -1,3 +1,4 @@
 ATMCHANGE=$CIMEROOT/../components/eamxx/scripts/atmchange
 
 $ATMCHANGE semi_lagrange_trajectory_nsubstep=2 -b
+$ATMCHANGE semi_lagrange_trajectory_nvelocity=3 -b

components/homme/cmake/HommeMacros.cmake

@@ -829,6 +829,24 @@ MACRO(CREATE_CXX_VS_F90_TESTS_WITH_PROFILE TESTS_LIST testProfile)
   ENDFOREACH ()
 ENDMACRO(CREATE_CXX_VS_F90_TESTS_WITH_PROFILE)
 
+function(check_transport_error_norms
+    TEST_NAME TCEN_ERROR_ANCHOR TCEN_FILENAME TCEN_UPPER_BOUNDS)
+  # Check prescribed-wind tracer transport test error-norm output, which is
+  # available with some of these tests. This function encapsulates the setup
+  # steps to add a check after a tracer transport test runs. See
+  # TransportCheckErrorNorms.cmake.in for details.
+  configure_file(
+    ${HOMME_SOURCE_DIR}/cmake/TransportCheckErrorNorms.cmake.in
+    ${HOMME_BINARY_DIR}/tests/${TEST_NAME}/check.cmake
+    @ONLY)
+  add_test(
+    NAME "${TEST_NAME}_l2err"
+    COMMAND ${CMAKE_COMMAND} -P check.cmake
+    WORKING_DIRECTORY ${HOMME_BINARY_DIR}/tests/${TEST_NAME})
+  set_tests_properties(
+    "${TEST_NAME}_l2err" PROPERTIES DEPENDS "${TEST_NAME}")
+endfunction()
+
 macro(testQuadPrec HOMME_QUAD_PREC)
 
   TRY_COMPILE(COMPILE_RESULT_VAR

components/homme/cmake/TransportCheckErrorNorms.cmake.in

@@ -0,0 +1,39 @@
+# This utility checks output of the form
+#   planar_conv>                          l2                    linf
+#   planar_conv> Q 1  3.7612654153016604E-03  9.8620728873592828E-03
+#   planar_conv> Q 2  3.9526251584022527E-03  8.2639849437505641E-03
+# against caller-provided upper bounds on the l2 norms. Output of this sort
+# appears on stdout at the end of some tracer transport tests. 'planar_conv' can
+# be something else; it is caller-configurable. Example usage:
+#   set(TCEN_ERROR_ANCHOR "planar_conv")
+#   set(TCEN_FILENAME "out.txt")
+#   set(TCEN_UPPER_BOUNDS "3.9e-3;4.1e-3")
+#   configure_file(TransportCheckErrorNorms.cmake.in check.cmake @ONLY)
+#   add_test(
+#     NAME "${TEST_NAME}_l2err"
+#     COMMAND ${CMAKE_COMMAND} -P check.cmake
+#     WORKING_DIRECTORY ${HOMME_BINARY_DIR}/tests/${TEST_NAME})
+#   set_tests_properties(
+#     "${TEST_NAME}_l2err" PROPERTIES DEPENDS "${TEST_NAME}")
+# Function check_transport_error_norms in HommeMacros.cmake encapsulates the
+# above commands.
+
+function(check_error_norms test_name filename upper_bounds)
+  execute_process(COMMAND grep "${test_name}> Q" ${filename} OUTPUT_VARIABLE grep_stdout)
+  separate_arguments(tokens NATIVE_COMMAND "${grep_stdout}")
+  list(LENGTH upper_bounds nchk)
+  math(EXPR nchk "${nchk} - 1")
+  foreach(idx RANGE 0 ${nchk})
+    math(EXPR tidx "3 + 5 * ${idx}")
+    list(GET tokens ${tidx} e)
+    list(GET upper_bounds ${idx} ub)
+    if(e GREATER ub)
+      message(FATAL_ERROR "${e} should be <= ${ub} but is not")
+    endif()
+  endforeach()
+endfunction()
+
+check_error_norms(
+  "@TCEN_ERROR_ANCHOR@"
+  "@TCEN_FILENAME@"
+  "@TCEN_UPPER_BOUNDS@")

components/homme/cmake/machineFiles/aurora-aot.cmake

@@ -44,7 +44,8 @@ SET(CMAKE_Fortran_COMPILER "mpifort" CACHE STRING "")
 SET(CMAKE_CXX_COMPILER "mpicxx" CACHE STRING "")
 
 #AOT flags
-SET(SYCL_COMPILE_FLAGS "-std=c++17 -fsycl -fsycl-device-code-split=per_kernel -fno-sycl-id-queries-fit-in-int -fsycl-unnamed-lambda -Xclang -fsycl-allow-virtual-functions")
+#SET(SYCL_COMPILE_FLAGS "-std=c++17 -fsycl -fsycl-device-code-split=per_kernel -fno-sycl-id-queries-fit-in-int -fsycl-unnamed-lambda -Xclang -fsycl-allow-virtual-functions")
+SET(SYCL_COMPILE_FLAGS "-std=c++17 -fsycl -fsycl-device-code-split=per_kernel -fno-sycl-id-queries-fit-in-int -fsycl-unnamed-lambda")
 SET(SYCL_LINK_FLAGS "-Wl,--no-relax -flink-huge-device-code -fsycl-max-parallel-link-jobs=32 -fsycl -fsycl-device-code-split=per_kernel -fsycl-targets=intel_gpu_pvc")
 
 SET(ADD_Fortran_FLAGS "-fc=ifx -fpscomp logicals -O3 -DNDEBUG -DCPRINTEL -g" CACHE STRING "")

components/homme/cmake/machineFiles/chrysalis.cmake

@@ -54,6 +54,7 @@ IF (${IFORT_RESULT} EQUAL 0)
   SET (ADD_Fortran_FLAGS "-traceback" CACHE STRING "")
   SET (ADD_C_FLAGS       "-traceback" CACHE STRING "")
   SET (ADD_CXX_FLAGS     "-traceback" CACHE STRING "")
+  SET (BUILD_HOMME_THETA_KOKKOS TRUE CACHE BOOL "")
 ELSE()
   SET (MKLROOT $ENV{MKLROOT} CACHE FILEPATH "")
   SET (HOMME_FIND_BLASLAPACK TRUE CACHE BOOL "")

components/homme/src/preqx_kokkos/CMakeLists.txt

@@ -97,6 +97,7 @@ MACRO(PREQX_KOKKOS_SETUP)
     ${SRC_SHARE_DIR}/cxx/prim_cxx_driver_base.F90
     ${SRC_SHARE_DIR}/cxx/utilities/bfb_mod.F90
     ${SRC_SHARE_DIR}/planar_mod.F90
+    ${TEST_SRC_DIR}/planar_transport.F90
     ${SRC_SHARE_DIR}/geometry_mod.F90
     ${SRC_SHARE_DIR}/planar_mesh_mod.F90
   )

components/homme/src/share/compose/compose_hommexx.cpp

@@ -24,8 +24,8 @@ using View = typename TracerArrays<ko::MachineTraits>::View<DataType>;
 void set_views (const SetView<HommexxReal***>& spheremp,
                 const SetView<HommexxReal****>& dp, const SetView5& dp3d,
                 const SetView<HommexxReal******>& qdp, const SetView5& q,
-                const SetView5& dep_points, const SetView5& vnode,
-                const SetView5& vdep, const Int ndim) {
+                const SetView5& dep_points, const SetView5& vnode, const Int ndim,
+                const SetView5& vdep, const Int vdep_ndim) {
   static_assert(std::is_same<Real, HommexxReal>::value,
                 "Hommexx and Compose real types must be the same.");
 #ifdef COMPOSE_PORT
@@ -42,23 +42,22 @@ void set_views (const SetView<HommexxReal***>& spheremp,
   if (vnode.data())
     ta.vnode = View<Real****>(vnode.data(), nel, vnode.extent_int(1), np2, ndim);
   if (vdep.data())
-    ta.vdep  = View<Real****>(vdep.data(),  nel, vdep .extent_int(1), np2, ndim);
+    ta.vdep  = View<Real****>(vdep.data(),  nel, vdep .extent_int(1), np2, vdep_ndim);
 #else
   slmm_throw_if(true, "Running a Hommexx code path with the non-Hommexx build"
                 " is not supported.\n");
 #endif
 }
 
-void set_hvcoord (const HommexxReal etai_beg, const HommexxReal etai_end,
-                  const HommexxReal* etam) {
+void set_hvcoord (const HommexxReal* etai, const HommexxReal* etam) {
   auto& cm = *get_isl_mpi_singleton();
-  islmpi::set_hvcoord(cm, etai_beg, etai_end, etam);
+  islmpi::set_hvcoord(cm, etai, etam);
 }
 
-void calc_v_departure (const int step, const HommexxReal dtsub) {
+void interp_v_update (const int step, const HommexxReal dtsub) {
   auto& cm = *get_isl_mpi_singleton();
-  islmpi::calc_v_departure<>(cm, 0, cm.nelemd - 1, step, dtsub,
-                             nullptr, nullptr, nullptr);
+  islmpi::interp_v_update<>(cm, 0, cm.nelemd - 1, step, dtsub,
+                            nullptr, nullptr, nullptr);
 }
 
 void advect (const int np1, const int n0_qdp, const int np1_qdp) {

components/homme/src/share/compose/compose_hommexx.hpp

@@ -16,12 +16,11 @@ typedef SetView<HommexxReal*****> SetView5;
 void set_views(const SetView<HommexxReal***>& spheremp,
                const SetView<HommexxReal****>& dp, const SetView5& dp3d,
                const SetView<HommexxReal******>& qdp, const SetView5& q,
-               const SetView5& dep_points, const SetView5& vnode,
-               const SetView5& vdep, const int trajectory_ndim);
+               const SetView5& dep_points, const SetView5& vnode, const int ndim,
+               const SetView5& vdep, const int vdep_ndim);
 
-void set_hvcoord(const HommexxReal etai_beg, const HommexxReal etai_end,
-                 const HommexxReal* etam);
-void calc_v_departure(const int step, const HommexxReal dtsub);
+void set_hvcoord(const HommexxReal* etai, const HommexxReal* etam);
+void interp_v_update(const int step, const HommexxReal dtsub);
 
 void advect(const int np1, const int n0_qdp, const int np1_qdp);
 

components/homme/src/share/compose/compose_slmm.cpp

@@ -124,41 +124,47 @@ void finalize_init_phase (IslMpi<MT>& cm, typename IslMpi<MT>::Advecter& advecte
 }
 
 template <typename MT>
-void set_hvcoord (IslMpi<MT>& cm, const Real etai_beg, const Real etai_end,
-                  const Real* etam) {
+void set_hvcoord (IslMpi<MT>& cm, const Real* etai, const Real* etam) {
   if (cm.etam.size() > 0) return;
 #if defined COMPOSE_HORIZ_OPENMP
 # pragma omp barrier
 # pragma omp master
 #endif
   {
     slmm_assert(cm.nlev > 0);
-    cm.etai_beg = etai_beg;
-    cm.etai_end = etai_end;
+    cm.etai_beg = etai[0];
+    cm.etai_end = etai[cm.nlev];
+    cm.etai = typename IslMpi<MT>::template ArrayD<Real*>("etai", cm.nlev+1);
     cm.etam = typename IslMpi<MT>::template ArrayD<Real*>("etam", cm.nlev);
-    const auto h = ko::create_mirror_view(cm.etam);
+    const auto hi = ko::create_mirror_view(cm.etai);
+    const auto hm = ko::create_mirror_view(cm.etam);
+    for (int k = 0; k <= cm.nlev; ++k) {
+      hi(k) = etai[k];
+      slmm_assert(k == 0 or hi(k) > hi(k-1));
+      slmm_assert(hi(k) >= 0 && hi(k) <= 1);
+    }
     for (int k = 0; k < cm.nlev; ++k) {
-      h(k) = etam[k];
-      slmm_assert(k == 0 || h(k) > h(k-1));
-      slmm_assert(h(k) > 0 && h(k) < 1);
+      hm(k) = etam[k];
+      slmm_assert(k == 0 or hm(k) > hm(k-1));
+      slmm_assert(hm(k) > 0 && hm(k) < 1);
     }
-    ko::deep_copy(cm.etam, h);
+    ko::deep_copy(cm.etai, hi);
+    ko::deep_copy(cm.etam, hm);
   }
 #if defined COMPOSE_HORIZ_OPENMP
 # pragma omp barrier
 #endif
 }
 
 template void set_hvcoord(
-  IslMpi<ko::MachineTraits>& cm, const Real etai_beg, const Real etai_end,
-  const Real* etam);
+  IslMpi<ko::MachineTraits>& cm, const Real* etai, const Real* etam);
 
 // Set pointers to HOMME data arrays.
 template <typename MT>
 void set_elem_data (IslMpi<MT>& cm, const Int ie, Real* qdp, const Int n0_qdp,
                     const Real* dp, Real* q, const Int nelem_in_patch) {
   slmm_assert(ie < cm.ed_h.size());
-  slmm_assert(cm.halo > 1 || cm.ed_h(ie).nbrs.size() == nelem_in_patch);
+  slmm_assert(cm.halo > 1 or cm.ed_h(ie).nbrs.size() == nelem_in_patch);
   auto& e = cm.ed_h(ie);
 #if defined COMPOSE_PORT
   cm.tracer_arrays->pqdp.set_ie_ptr(ie, qdp);
@@ -395,15 +401,14 @@ void slmm_check_ref2sphere (homme::Int ie, homme::Cartesian3D* p) {
   amb::dev_fin_threads();
 }
 
-void slmm_set_hvcoord (const homme::Real etai_beg, const homme::Real etai_end,
-                       const homme::Real* etam) {
+void slmm_set_hvcoord (const homme::Real* etai, const homme::Real* etam) {
   amb::dev_init_threads();
   slmm_assert(homme::g_csl_mpi);
-  homme::islmpi::set_hvcoord(*homme::g_csl_mpi, etai_beg, etai_end, etam);
+  homme::islmpi::set_hvcoord(*homme::g_csl_mpi, etai, etam);
   amb::dev_fin_threads();
 }
 
-void slmm_calc_v_departure (
+void slmm_interp_v_update (
   homme::Int nets, homme::Int nete, homme::Int step, homme::Real dtsub,
   homme::Real* dep_points, homme::Int dep_points_ndim, homme::Real* vnode,
   homme::Real* vdep, homme::Int* info)
@@ -419,8 +424,8 @@ void slmm_calc_v_departure (
     homme::sl_traj_h2d(*cm.tracer_arrays, dep_points, vnode, vdep,
                        cm.dep_points_ndim);
   }
-  homme::islmpi::calc_v_departure(cm, nets - 1, nete - 1, step - 1,
-                                  dtsub, dep_points, vnode, vdep);
+  homme::islmpi::interp_v_update(cm, nets - 1, nete - 1, step - 1,
+                                 dtsub, dep_points, vnode, vdep);
   *info = 0;
   {
     slmm::Timer timer("d2h");