From 0c790f2dcb975f0aeae024bc8651cfa4096716cc Mon Sep 17 00:00:00 2001 From: Marshall Ward Date: Fri, 10 May 2024 15:08:21 -0400 Subject: [PATCH] SAL: Explicitly allocate sht variable in SAL_CS For reasons not entirely clear, deallocation of the RK2 control structure caused a segmentation fault when compiled by Nvidia. A deeper investigation suggested that the compiler was attempting to deallocate a nullified pointer during automatic cleanup of RK2->SAL->sht. Apparently deallocation of a NULL pointer is an error (even though free() is not supposed to care). By redefining `sht` as an allocatable component rather than a local component of the SAL_CS instance, it seemed to satisfy the compiler that nothing was needed and the error went away. There are still some lingering questions behind the cause of the segfault, but for now I am going to put this under "compiler bug". This patch also initializes some of the logical flags in the SAL_CS type. This is possibly unnecessary, but it is consistent with the general rules of safety followed in other MOM6 derived types. --- .../lateral/MOM_self_attr_load.F90 | 16 +++++++++++----- 1 file changed, 11 insertions(+), 5 deletions(-) diff --git a/src/parameterizations/lateral/MOM_self_attr_load.F90 b/src/parameterizations/lateral/MOM_self_attr_load.F90 index 7f7215c9d8..e7f8a73ab2 100644 --- a/src/parameterizations/lateral/MOM_self_attr_load.F90 +++ b/src/parameterizations/lateral/MOM_self_attr_load.F90 @@ -20,14 +20,17 @@ module MOM_self_attr_load !> The control structure for the MOM_self_attr_load module type, public :: SAL_CS ; private - logical :: use_sal_scalar !< If true, use the scalar approximation to calculate SAL. - logical :: use_sal_sht !< If true, use online spherical harmonics to calculate SAL - logical :: use_tidal_sal_prev !< If true, read the tidal SAL from the previous iteration of - !! the tides to facilitate convergence. + logical :: use_sal_scalar = .false. + !< If true, use the scalar approximation to calculate SAL. + logical :: use_sal_sht = .false. + !< If true, use online spherical harmonics to calculate SAL + logical :: use_tidal_sal_prev = .false. + !< If true, read the tidal SAL from the previous iteration of the tides to + !! facilitate convergence. real :: sal_scalar_value !< The constant of proportionality between sea surface height !! (really it should be bottom pressure) anomalies and bottom !! geopotential anomalies [nondim]. - type(sht_CS) :: sht !< Spherical harmonic transforms (SHT) control structure + type(sht_CS), allocatable :: sht !< Spherical harmonic transforms (SHT) control structure integer :: sal_sht_Nd !< Maximum degree for SHT [nodim] real, allocatable :: Love_Scaling(:) !< Love number for each SHT mode [nodim] real, allocatable :: Snm_Re(:), & !< Real SHT coefficient for SHT SAL [Z ~> m] @@ -218,6 +221,8 @@ subroutine SAL_init(G, US, param_file, CS) allocate(CS%Love_Scaling(lmax)); CS%Love_Scaling(:) = 0.0 call calc_love_scaling(CS%sal_sht_Nd, rhoW, rhoE, CS%Love_Scaling) + + allocate(CS%sht) call spherical_harmonics_init(G, param_file, CS%sht) endif @@ -234,6 +239,7 @@ subroutine SAL_end(CS) if (allocated(CS%Snm_Re)) deallocate(CS%Snm_Re) if (allocated(CS%Snm_Im)) deallocate(CS%Snm_Im) call spherical_harmonics_end(CS%sht) + deallocate(CS%sht) endif end subroutine SAL_end