Icb into v2.5

CMakeLists.txt

@@ -2,8 +2,8 @@ cmake_minimum_required(VERSION 3.9)
 
 # set default build type cache entry (do so before project(...) is called, which would create this cache entry on its own)
 if(NOT CMAKE_BUILD_TYPE)
-	message(STATUS "setting default build type: Release")
-	set(CMAKE_BUILD_TYPE "Release" CACHE STRING "Choose the type of build, options are: None(CMAKE_CXX_FLAGS or CMAKE_C_FLAGS used) Debug Release RelWithDebInfo MinSizeRel.")
+  message(STATUS "setting default build type: Release")
+  set(CMAKE_BUILD_TYPE "Release" CACHE STRING "Choose the type of build, options are: None(CMAKE_CXX_FLAGS or CMAKE_C_FLAGS used) Debug Release RelWithDebInfo MinSizeRel.")
 endif()
 
 project(FESOM2.0)
@@ -14,7 +14,8 @@ set(OIFS_COUPLED OFF CACHE BOOL "compile fesom coupled to OpenIFS. (Also needs F
 set(CRAY OFF CACHE BOOL "compile with cray ftn")
 set(USE_ICEPACK OFF CACHE BOOL "compile fesom with the Iceapck modules for sea ice column physics.")
 set(OPENMP_REPRODUCIBLE OFF CACHE BOOL "serialize OpenMP loops that are critical for reproducible results")
+set(ASYNC_ICEBERGS ON CACHE BOOL "compile fesom with or without support for asynchronous iceberg computations")
 
-#set(VERBOSE OFF CACHE BOOL "toggle debug output")
+set(VERBOSE OFF CACHE BOOL "toggle debug output")
 #add_subdirectory(oasis3-mct/lib/psmile)
 add_subdirectory(src)

config/namelist.config

@@ -65,3 +65,11 @@ flag_debug=.false.
 n_levels=2
 n_part= 12, 36          ! 432 number of partitions on each hierarchy level
 /
+
+&icebergs
+use_icesheet_coupling=.false.
+ib_num=1
+use_icebergs=.false.
+steps_per_ib_step=8
+ib_async_mode=0
+/

config/namelist.io

@@ -12,7 +12,7 @@ ldiag_extflds    =.false.
 /
 
 &nml_general
-io_listsize    =100 !number of streams to allocate. shallbe large or equal to the number of streams in &nml_list
+io_listsize    =120 !number of streams to allocate. shallbe large or equal to the number of streams in &nml_list
 vec_autorotate =.false.
 /
 
@@ -48,4 +48,7 @@ io_list =  'sst       ',1, 'm', 4,
            'bolus_u   ',1, 'y', 4,
            'bolus_v   ',1, 'y', 4,
            'bolus_w   ',1, 'y', 4,
+           'icb       ',1, 'm', 4,
+           'fw        ',1, 'm', 4,
+           'fh        ',1, 'm', 4,
 /

env/levante.dkrz.de/shell

@@ -0,0 +1 @@
+shell.intel

env/levante.dkrz.de/shell.gnu

@@ -0,0 +1,45 @@
+# make the contents as shell agnostic as possible so we can include them with bash, zsh and others
+export LC_ALL=en_US.UTF-8
+export CPU_MODEL=AMD_EPYC_ZEN3
+
+module --force purge
+
+module load git
+
+module load gcc/11.2.0-gcc-11.2.0
+
+module load openmpi/4.1.2-gcc-11.2.0
+
+module load netcdf-c/4.8.1-gcc-11.2.0
+module load netcdf-fortran/4.5.3-gcc-11.2.0
+
+export FC=mpif90 CC=mpicc CXX=mpicxx 
+
+# following is only needed for libblas which is needed by params lib and often provided by lapack
+#module load intel-oneapi-mkl/2022.0.1-gcc-11.2.0
+# so use the LD_LIBRARY_PATH or other paths like prefix paths etc for cmake
+#export LD_LIBRARY_PATH=/sw/spack-levante/intel-oneapi-mkl-2022.0.1-ttdktf/mkl/2022.0.1/lib/intel64:$LD_LIBRARY_PATH
+spack load intel-oneapi-mkl@2022.0.1%gcc@11.2.0
+
+#other alternative blas
+#spack load netlib-lapack@3.9.1%gcc@11.2.0
+
+ulimit -s unlimited # without setting the stack size we get a segfault from the levante netcdf library at runtime
+ulimit -c 0 # do not create a coredump after a crash
+
+#OPENMPI specific runtime settings some deviation from ref:https://docs.dkrz.de/doc/levante/running-jobs/runtime-settings.html
+export OMPI_MCA_pml="ucx"
+export OMPI_MCA_btl=self
+export OMPI_MCA_osc="pt2pt"
+export UCX_IB_ADDR_TYPE=ib_global
+# for most runs one may or may not want to disable HCOLL
+export OMPI_MCA_coll="^ml,hcoll"
+export OMPI_MCA_coll_hcoll_enable="0"
+export HCOLL_ENABLE_MCAST_ALL="0"
+export HCOLL_MAIN_IB=mlx5_0:1
+export UCX_NET_DEVICES=mlx5_0:1
+export UCX_TLS=mm,knem,cma,dc_mlx5,dc_x,self
+export UCX_UNIFIED_MODE=y
+export HDF5_USE_FILE_LOCKING=FALSE
+export OMPI_MCA_io="romio321"
+export UCX_HANDLE_ERRORS=bt

env/levante.dkrz.de/shell.intel

@@ -0,0 +1,36 @@
+# make the contents as shell agnostic as possible so we can include them with bash, zsh and others
+export LC_ALL=en_US.UTF-8
+export CPU_MODEL=AMD_EPYC_ZEN3
+
+module --force purge 
+
+module load intel-oneapi-compilers/2022.0.1-gcc-11.2.0
+module load openmpi/4.1.2-intel-2021.5.0
+export FC=mpif90 CC=mpicc CXX=mpicxx ;
+spack load intel-oneapi-mkl@2022.0.1%gcc@11.2.0 # this handles adding to path elegantly then using hardcoded path below
+#module load intel-oneapi-mkl/2022.0.1-gcc-11.2.0
+#export LD_LIBRARY_PATH=/sw/spack-levante/intel-oneapi-mkl-2022.0.1-ttdktf/mkl/2022.0.1/lib/intel64:$LD_LIBRARY_PATH
+
+module load netcdf-c/4.8.1-openmpi-4.1.2-intel-2021.5.0
+module load netcdf-fortran/4.5.3-openmpi-4.1.2-intel-2021.5.0
+module load git # to be able to determine the fesom git SHA when compiling
+
+ulimit -s unlimited # without setting the stack size we get a segfault from the levante netcdf library at runtime
+ulimit -c 0 # do not create a coredump after a crash
+
+# environment for Open MPI 4.0.0 and later from https://docs.dkrz.de/doc/levante/running-jobs/runtime-settings.html
+export OMPI_MCA_pml="ucx"
+export OMPI_MCA_btl=self
+export OMPI_MCA_osc="pt2pt"
+export UCX_IB_ADDR_TYPE=ib_global
+# for most runs one may or may not want to disable HCOLL
+export OMPI_MCA_coll="^ml,hcoll"
+export OMPI_MCA_coll_hcoll_enable="0"
+export HCOLL_ENABLE_MCAST_ALL="0"
+export HCOLL_MAIN_IB=mlx5_0:1
+export UCX_NET_DEVICES=mlx5_0:1
+export UCX_TLS=mm,knem,cma,dc_mlx5,dc_x,self
+export UCX_UNIFIED_MODE=y
+export HDF5_USE_FILE_LOCKING=FALSE
+export OMPI_MCA_io="romio321"
+export UCX_HANDLE_ERRORS=bt

env/levante.dkrz.de/shell.nvhpc

@@ -0,0 +1,29 @@
+# make the contents as shell agnostic as possible so we can include them with bash, zsh and others
+export LC_ALL=en_US.UTF-8
+export CPU_MODEL=AMD_EPYC_ZEN3
+
+# module load intel-oneapi-compilers/2022.0.1-gcc-11.2.0
+# module load openmpi/4.1.2-intel-2021.5.0
+spack load openmpi@4.1.2%nvhpc@22.5
+export FC=mpif90 CC=mpicc CXX=mpicxx;
+
+module load netcdf-c/4.8.1-openmpi-4.1.2-intel-2021.5.0
+module load netcdf-fortran/4.5.3-openmpi-4.1.2-intel-2021.5.0
+module load git # to be able to determine the fesom git SHA when compiling
+
+# environment for Open MPI 4.0.0 and later from https://docs.dkrz.de/doc/levante/running-jobs/runtime-settings.html
+export OMPI_MCA_pml="ucx"
+export OMPI_MCA_btl=self
+export OMPI_MCA_osc="pt2pt"
+export UCX_IB_ADDR_TYPE=ib_global
+# for most runs one may or may not want to disable HCOLL
+export OMPI_MCA_coll="^ml,hcoll"
+export OMPI_MCA_coll_hcoll_enable="0"
+export HCOLL_ENABLE_MCAST_ALL="0"
+export HCOLL_MAIN_IB=mlx5_0:1
+export UCX_NET_DEVICES=mlx5_0:1
+export UCX_TLS=mm,knem,cma,dc_mlx5,dc_x,self
+export UCX_UNIFIED_MODE=y
+export HDF5_USE_FILE_LOCKING=FALSE
+export OMPI_MCA_io="romio321"
+export UCX_HANDLE_ERRORS=bt

env/mistral.dkrz.de/shell-intel+openmpi

@@ -42,9 +42,5 @@ export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:${NETCDFFROOT}/lib:${HDF5ROOT}/lib:${NET
 
 export FESOM_USE_CPLNG='active'
 
-export DR_HOOK=1
-export DR_HOOK_IGNORE_SIGNALS='-1'
-export DR_HOOK_OPT=prof
-export DR_HOOK_PROFILE_LIMIT=0.5
 export OIFS_DUMMY_ACTION=ABORT
 export HDF5_DISABLE_VERSION_CHECK=1

src/CMakeLists.txt

@@ -126,6 +126,9 @@ if(${FESOM_COUPLED})
    include(${CMAKE_CURRENT_LIST_DIR}/../cmake/FindOASIS.cmake)
    target_compile_definitions(${PROJECT_NAME} PRIVATE __oasis)
 endif()
+if(${ASYNC_ICEBERGS})
+   target_compile_definitions(${PROJECT_NAME} PRIVATE __async_icebergs)
+endif()
 if(${OIFS_COUPLED})
    target_compile_definitions(${PROJECT_NAME} PRIVATE __oifs)
 endif()

src/MOD_DYN.F90

@@ -72,6 +72,10 @@ MODULE MOD_DYN
     ! sea surface height arrays
     real(kind=WP), allocatable, dimension(:)    :: eta_n, d_eta, ssh_rhs, ssh_rhs_old
 
+    ! LA: 2023-05-17 iceberg arrays
+    real(kind=WP), allocatable, dimension(:)    :: eta_n_ib ! kh 18.03.21 additional array for asynchronous iceberg computations
+    real(kind=WP), allocatable, dimension(:,:,:):: uv_ib    ! kh 18.03.21 additional array for asynchronous iceberg computations
+
     !___________________________________________________________________________
     ! summarizes solver input parameter
     type(t_solverinfo)                          :: solverinfo

src/MOD_ICE.F90

@@ -131,8 +131,8 @@ MODULE MOD_ICE
 
     !___________________________________________________________________________
     ! zonal & merdional ice velocity
-    real(kind=WP), allocatable, dimension(:)    :: uice, uice_rhs, uice_old, uice_aux
-    real(kind=WP), allocatable, dimension(:)    :: vice, vice_rhs, vice_old, vice_aux
+    real(kind=WP), allocatable, dimension(:)    :: uice, uice_rhs, uice_old, uice_aux, uice_ib
+    real(kind=WP), allocatable, dimension(:)    :: vice, vice_rhs, vice_old, vice_aux, vice_ib
 
     ! surface stess atm<-->ice, oce<-->ice
     real(kind=WP), allocatable, dimension(:)    :: stress_atmice_x, stress_iceoce_x
@@ -153,9 +153,18 @@ MODULE MOD_ICE
     ! total number of ice tracers (default=3, 1=area, 2=mice, 3=msnow, (4=ice_temp)
 #if defined (__oifs) || defined (__ifsinterface)
     integer                                     :: num_itracers=4
+#else
+!    integer                                     :: num_itracers=3
+    !------------------------------
+    ! LA 2023-01-31 add icebergs
+#if defined(__async_icebergs)
+    integer                                     :: num_itracers=5
 #else
     integer                                     :: num_itracers=3
 #endif 
+    !------------------------------
+#endif 
+
 
     ! put ice tracers data arrays
     type(t_ice_data), allocatable, dimension(:) :: data
@@ -769,15 +778,17 @@ subroutine ice_init(ice, partit, mesh)
     ! to here since namelist.ice is now read in ice_init where whichEVP is not available
     ! when  mesh_auxiliary_arrays is called
     !array of 2D boundary conditions is used in ice_maEVP
-    if (ice%whichEVP > 0) then
+
+    ! LA 2023-05-24 initiate bc_index_nod2D also for whichEVP==0
+    !if (ice%whichEVP > 0) then
         allocate(mesh%bc_index_nod2D(myDim_nod2D+eDim_nod2D))
         mesh%bc_index_nod2D=1._WP
         do n=1, myDim_edge2D
             ed=mesh%edges(:, n)
             if (myList_edge2D(n) <= mesh%edge2D_in) cycle
             mesh%bc_index_nod2D(ed)=0._WP
         end do
-    end if
+    !end if
 
 end subroutine ice_init  
 !

src/MOD_MESH.F90

@@ -107,7 +107,13 @@ MODULE MOD_MESH
 ! Arrays added for ALE implementation:
 ! --> layer thinkness at node and depthlayer for t=n and t=n+1
 real(kind=WP), allocatable,dimension(:,:)   :: hnode, hnode_new, zbar_3d_n, Z_3d_n
+!------------------------------
+! LA 2023-01-31 add icebergs
+!#if defined(__async_icebergs)
+        real(kind=WP), allocatable,dimension(:,:)   :: Z_3d_n_ib
+!#endif 
 
+!------------------------------
 ! --> layer thinkness at elements, interpolated from hnode
 real(kind=WP), allocatable,dimension(:,:)   :: helem
 
@@ -135,6 +141,7 @@ MODULE MOD_MESH
 
 character(:), allocatable :: representative_checksum
 
+
 contains
 #if defined(__PGI)
   private
@@ -225,6 +232,10 @@ subroutine write_t_mesh(mesh, unit, iostat, iomsg)
     call write_bin_array(mesh%hnode_new,               unit, iostat, iomsg)
     call write_bin_array(mesh%zbar_3d_n,               unit, iostat, iomsg)
     call write_bin_array(mesh%Z_3d_n,                  unit, iostat, iomsg)
+! LA 2023-01-31 add icebergs
+!#if defined(__async_icebergs)
+    call write_bin_array(mesh%Z_3d_n_ib,               unit, iostat, iomsg)
+!#endif 
     call write_bin_array(mesh%helem,                   unit, iostat, iomsg)
     call write_bin_array(mesh%bottom_elem_thickness,   unit, iostat, iomsg)
     call write_bin_array(mesh%bottom_node_thickness,   unit, iostat, iomsg)
@@ -322,6 +333,10 @@ subroutine read_t_mesh(mesh, unit, iostat, iomsg)
     call read_bin_array(mesh%hnode_new,               unit, iostat, iomsg)
     call read_bin_array(mesh%zbar_3d_n,               unit, iostat, iomsg)
     call read_bin_array(mesh%Z_3d_n,                  unit, iostat, iomsg)
+! LA 2023-01-31 add icebergs
+!#if defined(__async_icebergs)
+    call read_bin_array(mesh%Z_3d_n_ib,               unit, iostat, iomsg)
+!#endif 
     call read_bin_array(mesh%helem,                   unit, iostat, iomsg)
     call read_bin_array(mesh%bottom_elem_thickness,   unit, iostat, iomsg)
     call read_bin_array(mesh%bottom_node_thickness,   unit, iostat, iomsg)

src/MOD_PARTIT.F90

@@ -39,6 +39,12 @@ module MOD_PARTIT
 
 TYPE T_PARTIT
   integer              :: MPI_COMM_FESOM ! FESOM communicator (for ocean only runs if often a copy of MPI_COMM_WORLD)
+
+  !---------------------------------------------------
+  !LA 2023-01-31 add asynchronous icebergs
+  ! kh 10.02.21 communicator for async iceberg computations based on OpenMP
+  integer              :: MPI_COMM_FESOM_IB
+  !---------------------------------------------------
 
   type(com_struct) :: com_nod2D
   type(com_struct) :: com_elem2D
@@ -61,6 +67,11 @@ module MOD_PARTIT
 
   ! general MPI part
   integer            :: MPIERR
+  !---------------------------------------------------
+  !LA 2023-01-31 add asynchronous icebergs
+  ! kh 11.02.21
+  integer            :: MPIERR_IB
+  !---------------------------------------------------
   integer            :: npes
   integer            :: mype
   integer            :: maxPEnum=100

src/associate_mesh_ass.h

@@ -54,6 +54,9 @@ hnode(1:mesh%nl-1, 1:myDim_nod2D+eDim_nod2D)               => mesh%hnode(:,:)
 hnode_new(1:mesh%nl-1, 1:myDim_nod2D+eDim_nod2D)           => mesh%hnode_new(:,:)
 zbar_3d_n(1:mesh%nl, 1:myDim_nod2D+eDim_nod2D)             => mesh%zbar_3d_n(:,:)
 Z_3d_n(1:mesh%nl-1, 1:myDim_nod2D+eDim_nod2D)              => mesh%Z_3d_n(:,:)
+#if defined(__async_icebergs)
+Z_3d_n_ib(1:mesh%nl-1, 1:myDim_nod2D+eDim_nod2D)           => mesh%Z_3d_n_ib(:,:)
+#endif
 helem(1:mesh%nl-1, 1:myDim_elem2D)                         => mesh%helem(:,:)
 bottom_elem_thickness(1:myDim_elem2D)                      => mesh%bottom_elem_thickness(:)
 bottom_node_thickness(1:myDim_nod2D+eDim_nod2D)            => mesh%bottom_node_thickness(:)

src/associate_mesh_def.h

@@ -38,6 +38,7 @@ real(kind=WP), dimension(:,:), pointer :: hnode
 real(kind=WP), dimension(:,:), pointer :: hnode_new
 real(kind=WP), dimension(:,:), pointer :: zbar_3d_n
 real(kind=WP), dimension(:,:), pointer :: Z_3d_n
+real(kind=WP), dimension(:,:), pointer :: Z_3d_n_ib
 real(kind=WP), dimension(:,:), pointer :: helem
 real(kind=WP), dimension(:)  , pointer :: bottom_elem_thickness
 real(kind=WP), dimension(:)  , pointer :: bottom_node_thickness

src/associate_part_ass.h

@@ -1,4 +1,5 @@
 MPI_COMM_FESOM  => partit%MPI_COMM_FESOM
+MPI_COMM_FESOM_IB       => partit%MPI_COMM_FESOM_IB
 com_nod2D       => partit%com_nod2D
 com_elem2D      => partit%com_elem2D
 com_elem2D_full => partit%com_elem2D_full
@@ -12,6 +13,7 @@ eDim_edge2D     => partit%eDim_edge2D
 pe_status       => partit%pe_status
 elem_full_flag  => partit%elem_full_flag
 MPIERR          => partit%MPIERR
+MPIERR_IB               => partit%MPIERR_IB
 npes            => partit%npes
 mype            => partit%mype
 maxPEnum        => partit%maxPEnum

src/associate_part_def.h

@@ -1,5 +1,6 @@
 
-  integer,          pointer     :: MPI_COMM_FESOM ! FESOM communicator (for ocean only runs if often a copy of MPI_COMM_WORLD)
+  integer,          pointer     :: MPI_COMM_FESOM       ! FESOM communicator (for ocean only runs if often a copy of MPI_COMM_WORLD)
+  integer,          pointer     :: MPI_COMM_FESOM_IB    ! FESOM communicator copy for icebergs LA: 2023-05-22
   type(com_struct), pointer     :: com_nod2D
   type(com_struct), pointer     :: com_elem2D
   type(com_struct), pointer     :: com_elem2D_full
@@ -17,6 +18,7 @@
   integer, dimension(:,:,:), pointer :: s_mpitype_nod3D, r_mpitype_nod3D
 
   integer, pointer   :: MPIERR
+  integer, pointer   :: MPIERR_IB       ! copy for icebergs LA: 2023-05-22
   integer, pointer   :: npes
   integer, pointer   :: mype
   integer, pointer   :: maxPEnum