From bf2c83f3d6f9c45bdaf7132b74fbddaeadfafccb Mon Sep 17 00:00:00 2001 From: WenhaoChen89 <96131003+WenhaoChen89@users.noreply.github.com> Date: Thu, 7 Jul 2022 12:23:09 -0400 Subject: [PATCH 1/9] Change dumbbell initialization --- .DS_Store | Bin 0 -> 6148 bytes src/.DS_Store | Bin 0 -> 8196 bytes .../MOM_state_initialization.F90 | 2 +- src/user/dumbbell_initialization.F90 | 74 ++++++++++++++++-- 4 files changed, 69 insertions(+), 7 deletions(-) create mode 100644 .DS_Store create mode 100644 src/.DS_Store diff --git a/.DS_Store b/.DS_Store new file mode 100644 index 0000000000000000000000000000000000000000..967887c2aff4da305da78dd97e26ae40dcba8b4e GIT binary patch literal 6148 zcmeHKL5tHs6n>LNH|Zkupr8j~z-w6>#fo^zcDsiK^{^2=sMMsXHEcFhns$3Ah1~VO zc=m_rfAOU6&Aft5>q!tP^We?*X6Ah_ov&mj0|25oPXmA!0C?zxZWEgwM)Txb=J1G* zD5yq~s{~^zcFc1YtO8bn|5pKecO8gfc$TYS?fpIasZIa*CJaGA0UscPF$|!D=@#b; zX8z7`Aw|v*%FWv1ZY{3~&)^JB;1pg!gzXVrLR#Z|s&N`HfeRSdvVPREd^pE?$Iyr8 zyX!jR8l$mmWzM(VS;?ZDux4MANuHInPUpLD8t%302cG90dheo-Y7x!i*}Uw<(>K(6 zsZ^e<<2ZShjTVE}%@bA3;;a}=H6cq!nDX{@R-|gtQ}ZG%H8(O7p6~mE)}7_@LHBW6 z22X-jTP}mn!$)n|4Z5q9@7=n0|LOVQBEKvY35OLG4Ckz=YgL>=0MqdFCT^!biF@63MM>2+tl|h_A z6NZXtsDiH;!q9PCa$d+-88mbdeEAUkWWhHSp`VWRrA`MC8nmrdz$&n*K*MxfbpD_I z`ToC2vOTMSRba0a5Ng2H+b2D}E%C`QFgHtrq-Pe$~h5@#~m-C-xwaV9Zf2;{EM z;466cLG)cb=~vyoF*At^B8swI&{h3YUw6+}#dLZ?M4~+oJ)#{Va!}aTs!(f$$GNP? zhUvKoD}X223t}0dRa^139hL#hfMvikU>UFs{2L75%w}1l+NxIvl)3`ImeH*Qb@k&P*5wY^ma%hf(E|}03ag>QTrq^XI|vQOc+2?DwbgJE zYGur0Ru<-lBFxo8xXPVGOIxcg1D1hR2Keq?qFs7K`w!#X{e2fB*S{%Iurl>XQcC@T z-^YSql|qWhr$b8MpV0`d9%XryxLQDMQ=eYb7&!Kqz;>Yzpm_%FEvzfv6=tZwJB8-} zGfmV^j|*#lT|w~9`f80M*q9FTh>gPn;sW-Tpzt1j_($~i$Bqlw&%c5VCfU$@JNXn0 z<0u<78ehalseEp8%W<5l^Rja!C!LW$8fR@kea@{XQpUka-4C8d{YkfW>48i}ew6f6 z1&G2vTwXkjl2A_Ca-4*j!gln4Q>j$CwJX!<&E}oD>)rEab$9ADZr!fCO|LndRh-M$ zZrp#|J&p&7WWpiByn*SwZZGtE)J`SnO(#uanJ}q52fi06`h6AbP!Z`P7(POtf!aZy zP$VhQhoywvrXhS!LXHGMJr|H-(h;ZX@(h~@B;z>=()?;FL6?2YNX~4HAG5O7@6?cP z{lSnrCZryR^lsyKL%MkI+H4h%<5EM^Xz=eyC| ppt] real :: S_light, S_dense ! The lightest and densest salinities in the sponges [S ~> ppt]. real :: eta_IC_quanta ! The granularity of quantization of intial interface heights [Z-1 ~> m-1]. + logical :: dbrotate ! If true, rotate the domain. + logical :: use_ALE ! True if ALE is being used, False if in layered mode + ! This include declares and sets the variable "version". # include "version_variable.h" character(len=20) :: verticalCoordinate integer :: i, j, k, is, ie, js, je, nz + real :: x, y is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke @@ -128,6 +132,8 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, units='m', default=1.0e-3, scale=US%m_to_Z, do_not_log=just_read) call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE", verticalCoordinate, & default=DEFAULT_COORDINATE_MODE, do_not_log=just_read) + call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) + if(.not. use_ALE) verticalCoordinate = "LAYER" ! WARNING: this routine specifies the interface heights so that the last layer ! is vanished, even at maximum depth. In order to have a uniform @@ -141,8 +147,38 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, !enddo select case ( coordinateMode(verticalCoordinate) ) - - case ( REGRIDDING_LAYER, REGRIDDING_RHO ) ! Initial thicknesses for isopycnal coordinates + case ( REGRIDDING_LAYER) ! Initial thicknesses for isopycnal coordinates + call get_param(param_file, mdl, "DUMBBELL_ROTATION", dbrotate, & + 'Logical for rotation of dumbbell domain.',& + units='nondim', default=.false., do_not_log=just_read) + do j=js,je + do i=is,ie + ! Compute normalized zonal coordinates (x,y=0 at center of domain) + if (dbrotate) then + ! This is really y in the rotated case + x = G%geoLatT(i,j) + else + x = G%geoLonT(i,j) + endif + + eta1D(1) = 0.0 + eta1D(nz+1) = -depth_tot(i,j) + if (x<0.0) then + do k=nz,2, -1 + eta1D(k) = eta1D(k+1) + min_thickness + enddo + else + do k=2,nz + eta1D(k) = eta1D(k-1) - min_thickness + enddo + endif + + do k=1,nz + h(i,j,k) = GV%Z_to_H * (eta1D(k) - eta1D(k+1)) + enddo + enddo; enddo + + case ( REGRIDDING_RHO, REGRIDDING_HYCOM1) ! Initial thicknesses for isopycnal coordinates call get_param(param_file, mdl, "INITIAL_SSS", S_surf, & units='1e-3', default=34., scale=US%ppt_to_S, do_not_log=.true.) call get_param(param_file, mdl, "INITIAL_S_RANGE", S_range, & @@ -231,12 +267,18 @@ subroutine dumbbell_initialize_temperature_salinity ( T, S, h, G, GV, US, param_ real :: x ! The fractional position in the domain [nondim] real :: dblen ! The size of the dumbbell test case [axis_units] logical :: dbrotate ! If true, rotate the domain. + logical :: use_ALE ! If false, use layer mode. character(len=20) :: verticalCoordinate, density_profile is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke T_surf = 20.0*US%degC_to_C + ! layer mode + call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) + if(.not. use_ALE) call MOM_error(FATAL, "dumbbell_initialize_temperature_salinity: "//& + "Please use 'fit' for 'TS_CONFIG' in the LAYER mode.") + call get_param(param_file, mdl, "REGRIDDING_COORDINATE_MODE", verticalCoordinate, & default=DEFAULT_COORDINATE_MODE, do_not_log=just_read) call get_param(param_file, mdl, "INITIAL_DENSITY_PROFILE", density_profile, & @@ -288,11 +330,12 @@ subroutine dumbbell_initialize_temperature_salinity ( T, S, h, G, GV, US, param_ end subroutine dumbbell_initialize_temperature_salinity !> Initialize the restoring sponges for the dumbbell test case -subroutine dumbbell_initialize_sponges(G, GV, US, tv, depth_tot, param_file, use_ALE, CSp, ACSp) +subroutine dumbbell_initialize_sponges(G, GV, US, tv, h_in, depth_tot, param_file, use_ALE, CSp, ACSp) type(ocean_grid_type), intent(in) :: G !< Horizontal grid control structure type(verticalGrid_type), intent(in) :: GV !< Vertical grid control structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamic variables + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h_in !< Layer thickness [H ~> m or kg m-2] real, dimension(SZI_(G),SZJ_(G)), & intent(in) :: depth_tot !< The nominal total depth of the ocean [Z ~> m] type(param_file_type), intent(in) :: param_file !< Parameter file structure @@ -306,6 +349,7 @@ subroutine dumbbell_initialize_sponges(G, GV, US, tv, depth_tot, param_file, use real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h ! sponge thicknesses [H ~> m or kg m-2] real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: S ! sponge salinities [S ~> ppt] real, dimension(SZK_(GV)+1) :: eta1D ! interface positions for ALE sponge + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)+1) :: eta ! A temporary array for interface heights [Z ~> m]. integer :: i, j, k, nz real :: x ! The fractional position in the domain [nondim] @@ -404,9 +448,27 @@ subroutine dumbbell_initialize_sponges(G, GV, US, tv, depth_tot, param_file, use enddo endif enddo ; enddo - endif - if (associated(tv%S)) call set_up_ALE_sponge_field(S, G, GV, tv%S, ACSp) + + else + do j=G%jsc,G%jec ; do i=G%isc,G%iec + eta(i,j,1) = 0.0 + do k=2,nz + eta(i,j,k) = eta(i,j,k-1)-h_in(i,j,k-1) + enddo + eta(i,j,nz+1) = -depth_tot(i,j) + do k=1,nz + S(i,j,k)= tv%S(i,j,k) + enddo + enddo ; enddo + + ! This call sets up the damping rates and interface heights. + ! This sets the inverse damping timescale fields in the sponges. ! + call initialize_sponge(Idamp, eta, G, param_file, CSp, GV) + + ! The remaining calls to set_up_sponge_field can be in any order. ! + if ( associated(tv%S) ) call set_up_sponge_field(S, tv%S, G, GV, nz, CSp) + endif end subroutine dumbbell_initialize_sponges From 63d7624e48892be36e770dce0df0c0b9b70663ee Mon Sep 17 00:00:00 2001 From: WenhaoChen89 <96131003+WenhaoChen89@users.noreply.github.com> Date: Thu, 7 Jul 2022 12:34:50 -0400 Subject: [PATCH 2/9] Change in Dumbbell Layer Mode --- .DS_Store | Bin 6148 -> 0 bytes src/.DS_Store | Bin 8196 -> 0 bytes 2 files changed, 0 insertions(+), 0 deletions(-) delete mode 100644 .DS_Store delete mode 100644 src/.DS_Store diff --git a/.DS_Store b/.DS_Store deleted file mode 100644 index 967887c2aff4da305da78dd97e26ae40dcba8b4e..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 6148 zcmeHKL5tHs6n>LNH|Zkupr8j~z-w6>#fo^zcDsiK^{^2=sMMsXHEcFhns$3Ah1~VO zc=m_rfAOU6&Aft5>q!tP^We?*X6Ah_ov&mj0|25oPXmA!0C?zxZWEgwM)Txb=J1G* zD5yq~s{~^zcFc1YtO8bn|5pKecO8gfc$TYS?fpIasZIa*CJaGA0UscPF$|!D=@#b; zX8z7`Aw|v*%FWv1ZY{3~&)^JB;1pg!gzXVrLR#Z|s&N`HfeRSdvVPREd^pE?$Iyr8 zyX!jR8l$mmWzM(VS;?ZDux4MANuHInPUpLD8t%302cG90dheo-Y7x!i*}Uw<(>K(6 zsZ^e<<2ZShjTVE}%@bA3;;a}=H6cq!nDX{@R-|gtQ}ZG%H8(O7p6~mE)}7_@LHBW6 z22X-jTP}mn!$)n|4Z5q9@7=n0|LOVQBEKvY35OLG4Ckz=YgL>=0MqdFCT^!biF@63MM>2+tl|h_A z6NZXtsDiH;!q9PCa$d+-88mbdeEAUkWWhHSp`VWRrA`MC8nmrdz$&n*K*MxfbpD_I z`ToC2vOTMSRba0a5Ng2H+b2D}E%C`QFgHtrq-Pe$~h5@#~m-C-xwaV9Zf2;{EM z;466cLG)cb=~vyoF*At^B8swI&{h3YUw6+}#dLZ?M4~+oJ)#{Va!}aTs!(f$$GNP? zhUvKoD}X223t}0dRa^139hL#hfMvikU>UFs{2L75%w}1l+NxIvl)3`ImeH*Qb@k&P*5wY^ma%hf(E|}03ag>QTrq^XI|vQOc+2?DwbgJE zYGur0Ru<-lBFxo8xXPVGOIxcg1D1hR2Keq?qFs7K`w!#X{e2fB*S{%Iurl>XQcC@T z-^YSql|qWhr$b8MpV0`d9%XryxLQDMQ=eYb7&!Kqz;>Yzpm_%FEvzfv6=tZwJB8-} zGfmV^j|*#lT|w~9`f80M*q9FTh>gPn;sW-Tpzt1j_($~i$Bqlw&%c5VCfU$@JNXn0 z<0u<78ehalseEp8%W<5l^Rja!C!LW$8fR@kea@{XQpUka-4C8d{YkfW>48i}ew6f6 z1&G2vTwXkjl2A_Ca-4*j!gln4Q>j$CwJX!<&E}oD>)rEab$9ADZr!fCO|LndRh-M$ zZrp#|J&p&7WWpiByn*SwZZGtE)J`SnO(#uanJ}q52fi06`h6AbP!Z`P7(POtf!aZy zP$VhQhoywvrXhS!LXHGMJr|H-(h;ZX@(h~@B;z>=()?;FL6?2YNX~4HAG5O7@6?cP z{lSnrCZryR^lsyKL%MkI+H4h%<5EM^Xz=eyC| Date: Thu, 7 Jul 2022 13:04:26 -0400 Subject: [PATCH 3/9] Revert "Change dumbbell initialization" This reverts commit bf2c83f3d6f9c45bdaf7132b74fbddaeadfafccb. --- .../MOM_state_initialization.F90 | 2 +- src/user/dumbbell_initialization.F90 | 74 ++----------------- 2 files changed, 7 insertions(+), 69 deletions(-) diff --git a/src/initialization/MOM_state_initialization.F90 b/src/initialization/MOM_state_initialization.F90 index cd6ddbeca7..0757fd887f 100644 --- a/src/initialization/MOM_state_initialization.F90 +++ b/src/initialization/MOM_state_initialization.F90 @@ -585,7 +585,7 @@ subroutine MOM_initialize_state(u, v, h, tv, Time, G, GV, US, PF, dirs, & case ("USER"); call user_initialize_sponges(G, GV, use_temperature, tv, PF, sponge_CSp, h) case ("BFB"); call BFB_initialize_sponges_southonly(G, GV, US, use_temperature, tv, depth_tot, PF, & sponge_CSp, h) - case ("DUMBBELL"); call dumbbell_initialize_sponges(G, GV, US, tv, h, depth_tot, PF, useALE, & + case ("DUMBBELL"); call dumbbell_initialize_sponges(G, GV, US, tv, depth_tot, PF, useALE, & sponge_CSp, ALE_sponge_CSp) case ("phillips"); call Phillips_initialize_sponges(G, GV, US, tv, PF, sponge_CSp, h) case ("dense"); call dense_water_initialize_sponges(G, GV, US, tv, depth_tot, PF, useALE, & diff --git a/src/user/dumbbell_initialization.F90 b/src/user/dumbbell_initialization.F90 index 820ce69049..c197274067 100644 --- a/src/user/dumbbell_initialization.F90 +++ b/src/user/dumbbell_initialization.F90 @@ -16,7 +16,7 @@ module dumbbell_initialization use MOM_verticalGrid, only : verticalGrid_type use regrid_consts, only : coordinateMode, DEFAULT_COORDINATE_MODE use regrid_consts, only : REGRIDDING_LAYER, REGRIDDING_ZSTAR -use regrid_consts, only : REGRIDDING_RHO, REGRIDDING_SIGMA, REGRIDDING_HYCOM1 +use regrid_consts, only : REGRIDDING_RHO, REGRIDDING_SIGMA use MOM_ALE_sponge, only : ALE_sponge_CS, set_up_ALE_sponge_field, initialize_ALE_sponge implicit none ; private @@ -112,14 +112,10 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, real :: S_range ! The range of salinities in this test case [S ~> ppt] real :: S_light, S_dense ! The lightest and densest salinities in the sponges [S ~> ppt]. real :: eta_IC_quanta ! The granularity of quantization of intial interface heights [Z-1 ~> m-1]. - logical :: dbrotate ! If true, rotate the domain. - logical :: use_ALE ! True if ALE is being used, False if in layered mode - ! This include declares and sets the variable "version". # include "version_variable.h" character(len=20) :: verticalCoordinate integer :: i, j, k, is, ie, js, je, nz - real :: x, y is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke @@ -132,8 +128,6 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, units='m', default=1.0e-3, scale=US%m_to_Z, do_not_log=just_read) call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE", verticalCoordinate, & default=DEFAULT_COORDINATE_MODE, do_not_log=just_read) - call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) - if(.not. use_ALE) verticalCoordinate = "LAYER" ! WARNING: this routine specifies the interface heights so that the last layer ! is vanished, even at maximum depth. In order to have a uniform @@ -147,38 +141,8 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, !enddo select case ( coordinateMode(verticalCoordinate) ) - case ( REGRIDDING_LAYER) ! Initial thicknesses for isopycnal coordinates - call get_param(param_file, mdl, "DUMBBELL_ROTATION", dbrotate, & - 'Logical for rotation of dumbbell domain.',& - units='nondim', default=.false., do_not_log=just_read) - do j=js,je - do i=is,ie - ! Compute normalized zonal coordinates (x,y=0 at center of domain) - if (dbrotate) then - ! This is really y in the rotated case - x = G%geoLatT(i,j) - else - x = G%geoLonT(i,j) - endif - - eta1D(1) = 0.0 - eta1D(nz+1) = -depth_tot(i,j) - if (x<0.0) then - do k=nz,2, -1 - eta1D(k) = eta1D(k+1) + min_thickness - enddo - else - do k=2,nz - eta1D(k) = eta1D(k-1) - min_thickness - enddo - endif - - do k=1,nz - h(i,j,k) = GV%Z_to_H * (eta1D(k) - eta1D(k+1)) - enddo - enddo; enddo - - case ( REGRIDDING_RHO, REGRIDDING_HYCOM1) ! Initial thicknesses for isopycnal coordinates + + case ( REGRIDDING_LAYER, REGRIDDING_RHO ) ! Initial thicknesses for isopycnal coordinates call get_param(param_file, mdl, "INITIAL_SSS", S_surf, & units='1e-3', default=34., scale=US%ppt_to_S, do_not_log=.true.) call get_param(param_file, mdl, "INITIAL_S_RANGE", S_range, & @@ -267,18 +231,12 @@ subroutine dumbbell_initialize_temperature_salinity ( T, S, h, G, GV, US, param_ real :: x ! The fractional position in the domain [nondim] real :: dblen ! The size of the dumbbell test case [axis_units] logical :: dbrotate ! If true, rotate the domain. - logical :: use_ALE ! If false, use layer mode. character(len=20) :: verticalCoordinate, density_profile is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke T_surf = 20.0*US%degC_to_C - ! layer mode - call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) - if(.not. use_ALE) call MOM_error(FATAL, "dumbbell_initialize_temperature_salinity: "//& - "Please use 'fit' for 'TS_CONFIG' in the LAYER mode.") - call get_param(param_file, mdl, "REGRIDDING_COORDINATE_MODE", verticalCoordinate, & default=DEFAULT_COORDINATE_MODE, do_not_log=just_read) call get_param(param_file, mdl, "INITIAL_DENSITY_PROFILE", density_profile, & @@ -330,12 +288,11 @@ subroutine dumbbell_initialize_temperature_salinity ( T, S, h, G, GV, US, param_ end subroutine dumbbell_initialize_temperature_salinity !> Initialize the restoring sponges for the dumbbell test case -subroutine dumbbell_initialize_sponges(G, GV, US, tv, h_in, depth_tot, param_file, use_ALE, CSp, ACSp) +subroutine dumbbell_initialize_sponges(G, GV, US, tv, depth_tot, param_file, use_ALE, CSp, ACSp) type(ocean_grid_type), intent(in) :: G !< Horizontal grid control structure type(verticalGrid_type), intent(in) :: GV !< Vertical grid control structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamic variables - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h_in !< Layer thickness [H ~> m or kg m-2] real, dimension(SZI_(G),SZJ_(G)), & intent(in) :: depth_tot !< The nominal total depth of the ocean [Z ~> m] type(param_file_type), intent(in) :: param_file !< Parameter file structure @@ -349,7 +306,6 @@ subroutine dumbbell_initialize_sponges(G, GV, US, tv, h_in, depth_tot, param_fil real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h ! sponge thicknesses [H ~> m or kg m-2] real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: S ! sponge salinities [S ~> ppt] real, dimension(SZK_(GV)+1) :: eta1D ! interface positions for ALE sponge - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)+1) :: eta ! A temporary array for interface heights [Z ~> m]. integer :: i, j, k, nz real :: x ! The fractional position in the domain [nondim] @@ -448,27 +404,9 @@ subroutine dumbbell_initialize_sponges(G, GV, US, tv, h_in, depth_tot, param_fil enddo endif enddo ; enddo - if (associated(tv%S)) call set_up_ALE_sponge_field(S, G, GV, tv%S, ACSp) - - else - do j=G%jsc,G%jec ; do i=G%isc,G%iec - eta(i,j,1) = 0.0 - do k=2,nz - eta(i,j,k) = eta(i,j,k-1)-h_in(i,j,k-1) - enddo - eta(i,j,nz+1) = -depth_tot(i,j) - do k=1,nz - S(i,j,k)= tv%S(i,j,k) - enddo - enddo ; enddo - - ! This call sets up the damping rates and interface heights. - ! This sets the inverse damping timescale fields in the sponges. ! - call initialize_sponge(Idamp, eta, G, param_file, CSp, GV) + endif - ! The remaining calls to set_up_sponge_field can be in any order. ! - if ( associated(tv%S) ) call set_up_sponge_field(S, tv%S, G, GV, nz, CSp) - endif + if (associated(tv%S)) call set_up_ALE_sponge_field(S, G, GV, tv%S, ACSp) end subroutine dumbbell_initialize_sponges From 83f5777e0ba29a9c4c5eb64bc912f674908714ed Mon Sep 17 00:00:00 2001 From: WenhaoChen89 <96131003+WenhaoChen89@users.noreply.github.com> Date: Thu, 7 Jul 2022 13:09:41 -0400 Subject: [PATCH 4/9] Revert "Revert "Change dumbbell initialization"" This reverts commit aaf446d116382a86bceb6f28c5fb4b4fc66a70fc. --- .../MOM_state_initialization.F90 | 2 +- src/user/dumbbell_initialization.F90 | 74 +++++++++++++++++-- 2 files changed, 69 insertions(+), 7 deletions(-) diff --git a/src/initialization/MOM_state_initialization.F90 b/src/initialization/MOM_state_initialization.F90 index 0757fd887f..cd6ddbeca7 100644 --- a/src/initialization/MOM_state_initialization.F90 +++ b/src/initialization/MOM_state_initialization.F90 @@ -585,7 +585,7 @@ subroutine MOM_initialize_state(u, v, h, tv, Time, G, GV, US, PF, dirs, & case ("USER"); call user_initialize_sponges(G, GV, use_temperature, tv, PF, sponge_CSp, h) case ("BFB"); call BFB_initialize_sponges_southonly(G, GV, US, use_temperature, tv, depth_tot, PF, & sponge_CSp, h) - case ("DUMBBELL"); call dumbbell_initialize_sponges(G, GV, US, tv, depth_tot, PF, useALE, & + case ("DUMBBELL"); call dumbbell_initialize_sponges(G, GV, US, tv, h, depth_tot, PF, useALE, & sponge_CSp, ALE_sponge_CSp) case ("phillips"); call Phillips_initialize_sponges(G, GV, US, tv, PF, sponge_CSp, h) case ("dense"); call dense_water_initialize_sponges(G, GV, US, tv, depth_tot, PF, useALE, & diff --git a/src/user/dumbbell_initialization.F90 b/src/user/dumbbell_initialization.F90 index c197274067..820ce69049 100644 --- a/src/user/dumbbell_initialization.F90 +++ b/src/user/dumbbell_initialization.F90 @@ -16,7 +16,7 @@ module dumbbell_initialization use MOM_verticalGrid, only : verticalGrid_type use regrid_consts, only : coordinateMode, DEFAULT_COORDINATE_MODE use regrid_consts, only : REGRIDDING_LAYER, REGRIDDING_ZSTAR -use regrid_consts, only : REGRIDDING_RHO, REGRIDDING_SIGMA +use regrid_consts, only : REGRIDDING_RHO, REGRIDDING_SIGMA, REGRIDDING_HYCOM1 use MOM_ALE_sponge, only : ALE_sponge_CS, set_up_ALE_sponge_field, initialize_ALE_sponge implicit none ; private @@ -112,10 +112,14 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, real :: S_range ! The range of salinities in this test case [S ~> ppt] real :: S_light, S_dense ! The lightest and densest salinities in the sponges [S ~> ppt]. real :: eta_IC_quanta ! The granularity of quantization of intial interface heights [Z-1 ~> m-1]. + logical :: dbrotate ! If true, rotate the domain. + logical :: use_ALE ! True if ALE is being used, False if in layered mode + ! This include declares and sets the variable "version". # include "version_variable.h" character(len=20) :: verticalCoordinate integer :: i, j, k, is, ie, js, je, nz + real :: x, y is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke @@ -128,6 +132,8 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, units='m', default=1.0e-3, scale=US%m_to_Z, do_not_log=just_read) call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE", verticalCoordinate, & default=DEFAULT_COORDINATE_MODE, do_not_log=just_read) + call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) + if(.not. use_ALE) verticalCoordinate = "LAYER" ! WARNING: this routine specifies the interface heights so that the last layer ! is vanished, even at maximum depth. In order to have a uniform @@ -141,8 +147,38 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, !enddo select case ( coordinateMode(verticalCoordinate) ) - - case ( REGRIDDING_LAYER, REGRIDDING_RHO ) ! Initial thicknesses for isopycnal coordinates + case ( REGRIDDING_LAYER) ! Initial thicknesses for isopycnal coordinates + call get_param(param_file, mdl, "DUMBBELL_ROTATION", dbrotate, & + 'Logical for rotation of dumbbell domain.',& + units='nondim', default=.false., do_not_log=just_read) + do j=js,je + do i=is,ie + ! Compute normalized zonal coordinates (x,y=0 at center of domain) + if (dbrotate) then + ! This is really y in the rotated case + x = G%geoLatT(i,j) + else + x = G%geoLonT(i,j) + endif + + eta1D(1) = 0.0 + eta1D(nz+1) = -depth_tot(i,j) + if (x<0.0) then + do k=nz,2, -1 + eta1D(k) = eta1D(k+1) + min_thickness + enddo + else + do k=2,nz + eta1D(k) = eta1D(k-1) - min_thickness + enddo + endif + + do k=1,nz + h(i,j,k) = GV%Z_to_H * (eta1D(k) - eta1D(k+1)) + enddo + enddo; enddo + + case ( REGRIDDING_RHO, REGRIDDING_HYCOM1) ! Initial thicknesses for isopycnal coordinates call get_param(param_file, mdl, "INITIAL_SSS", S_surf, & units='1e-3', default=34., scale=US%ppt_to_S, do_not_log=.true.) call get_param(param_file, mdl, "INITIAL_S_RANGE", S_range, & @@ -231,12 +267,18 @@ subroutine dumbbell_initialize_temperature_salinity ( T, S, h, G, GV, US, param_ real :: x ! The fractional position in the domain [nondim] real :: dblen ! The size of the dumbbell test case [axis_units] logical :: dbrotate ! If true, rotate the domain. + logical :: use_ALE ! If false, use layer mode. character(len=20) :: verticalCoordinate, density_profile is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke T_surf = 20.0*US%degC_to_C + ! layer mode + call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) + if(.not. use_ALE) call MOM_error(FATAL, "dumbbell_initialize_temperature_salinity: "//& + "Please use 'fit' for 'TS_CONFIG' in the LAYER mode.") + call get_param(param_file, mdl, "REGRIDDING_COORDINATE_MODE", verticalCoordinate, & default=DEFAULT_COORDINATE_MODE, do_not_log=just_read) call get_param(param_file, mdl, "INITIAL_DENSITY_PROFILE", density_profile, & @@ -288,11 +330,12 @@ subroutine dumbbell_initialize_temperature_salinity ( T, S, h, G, GV, US, param_ end subroutine dumbbell_initialize_temperature_salinity !> Initialize the restoring sponges for the dumbbell test case -subroutine dumbbell_initialize_sponges(G, GV, US, tv, depth_tot, param_file, use_ALE, CSp, ACSp) +subroutine dumbbell_initialize_sponges(G, GV, US, tv, h_in, depth_tot, param_file, use_ALE, CSp, ACSp) type(ocean_grid_type), intent(in) :: G !< Horizontal grid control structure type(verticalGrid_type), intent(in) :: GV !< Vertical grid control structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamic variables + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h_in !< Layer thickness [H ~> m or kg m-2] real, dimension(SZI_(G),SZJ_(G)), & intent(in) :: depth_tot !< The nominal total depth of the ocean [Z ~> m] type(param_file_type), intent(in) :: param_file !< Parameter file structure @@ -306,6 +349,7 @@ subroutine dumbbell_initialize_sponges(G, GV, US, tv, depth_tot, param_file, use real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h ! sponge thicknesses [H ~> m or kg m-2] real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: S ! sponge salinities [S ~> ppt] real, dimension(SZK_(GV)+1) :: eta1D ! interface positions for ALE sponge + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)+1) :: eta ! A temporary array for interface heights [Z ~> m]. integer :: i, j, k, nz real :: x ! The fractional position in the domain [nondim] @@ -404,9 +448,27 @@ subroutine dumbbell_initialize_sponges(G, GV, US, tv, depth_tot, param_file, use enddo endif enddo ; enddo - endif - if (associated(tv%S)) call set_up_ALE_sponge_field(S, G, GV, tv%S, ACSp) + + else + do j=G%jsc,G%jec ; do i=G%isc,G%iec + eta(i,j,1) = 0.0 + do k=2,nz + eta(i,j,k) = eta(i,j,k-1)-h_in(i,j,k-1) + enddo + eta(i,j,nz+1) = -depth_tot(i,j) + do k=1,nz + S(i,j,k)= tv%S(i,j,k) + enddo + enddo ; enddo + + ! This call sets up the damping rates and interface heights. + ! This sets the inverse damping timescale fields in the sponges. ! + call initialize_sponge(Idamp, eta, G, param_file, CSp, GV) + + ! The remaining calls to set_up_sponge_field can be in any order. ! + if ( associated(tv%S) ) call set_up_sponge_field(S, tv%S, G, GV, nz, CSp) + endif end subroutine dumbbell_initialize_sponges From dc82d3ad952befc9bd555259ae184479a3711f2e Mon Sep 17 00:00:00 2001 From: WenhaoChen89 <96131003+WenhaoChen89@users.noreply.github.com> Date: Thu, 7 Jul 2022 17:32:46 -0400 Subject: [PATCH 5/9] Update dumbbell_initialization.F90 --- src/user/dumbbell_initialization.F90 | 57 +++++++++++++--------------- 1 file changed, 27 insertions(+), 30 deletions(-) diff --git a/src/user/dumbbell_initialization.F90 b/src/user/dumbbell_initialization.F90 index 820ce69049..9390f26dca 100644 --- a/src/user/dumbbell_initialization.F90 +++ b/src/user/dumbbell_initialization.F90 @@ -113,8 +113,8 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, real :: S_light, S_dense ! The lightest and densest salinities in the sponges [S ~> ppt]. real :: eta_IC_quanta ! The granularity of quantization of intial interface heights [Z-1 ~> m-1]. logical :: dbrotate ! If true, rotate the domain. - logical :: use_ALE ! True if ALE is being used, False if in layered mode - + logical :: use_ALE ! True if ALE is being used, False if in layered mode + ! This include declares and sets the variable "version". # include "version_variable.h" character(len=20) :: verticalCoordinate @@ -132,7 +132,7 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, units='m', default=1.0e-3, scale=US%m_to_Z, do_not_log=just_read) call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE", verticalCoordinate, & default=DEFAULT_COORDINATE_MODE, do_not_log=just_read) - call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) + call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) if(.not. use_ALE) verticalCoordinate = "LAYER" ! WARNING: this routine specifies the interface heights so that the last layer @@ -148,36 +148,34 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, select case ( coordinateMode(verticalCoordinate) ) case ( REGRIDDING_LAYER) ! Initial thicknesses for isopycnal coordinates - call get_param(param_file, mdl, "DUMBBELL_ROTATION", dbrotate, & + call get_param(param_file, mdl, "DUMBBELL_ROTATION", dbrotate, & 'Logical for rotation of dumbbell domain.',& units='nondim', default=.false., do_not_log=just_read) do j=js,je do i=is,ie ! Compute normalized zonal coordinates (x,y=0 at center of domain) if (dbrotate) then - ! This is really y in the rotated case - x = G%geoLatT(i,j) + ! This is really y in the rotated case + x = G%geoLatT(i,j) else - x = G%geoLonT(i,j) + x = G%geoLonT(i,j) endif - - eta1D(1) = 0.0 - eta1D(nz+1) = -depth_tot(i,j) + eta1D(1) = 0.0 + eta1D(nz+1) = -depth_tot(i,j) if (x<0.0) then do k=nz,2, -1 eta1D(k) = eta1D(k+1) + min_thickness - enddo + enddo else do k=2,nz eta1D(k) = eta1D(k-1) - min_thickness enddo endif - do k=1,nz - h(i,j,k) = GV%Z_to_H * (eta1D(k) - eta1D(k+1)) + h(i,j,k) = GV%Z_to_H * (eta1D(k) - eta1D(k+1)) enddo - enddo; enddo - + enddo; enddo + case ( REGRIDDING_RHO, REGRIDDING_HYCOM1) ! Initial thicknesses for isopycnal coordinates call get_param(param_file, mdl, "INITIAL_SSS", S_surf, & units='1e-3', default=34., scale=US%ppt_to_S, do_not_log=.true.) @@ -278,7 +276,7 @@ subroutine dumbbell_initialize_temperature_salinity ( T, S, h, G, GV, US, param_ call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) if(.not. use_ALE) call MOM_error(FATAL, "dumbbell_initialize_temperature_salinity: "//& "Please use 'fit' for 'TS_CONFIG' in the LAYER mode.") - + call get_param(param_file, mdl, "REGRIDDING_COORDINATE_MODE", verticalCoordinate, & default=DEFAULT_COORDINATE_MODE, do_not_log=just_read) call get_param(param_file, mdl, "INITIAL_DENSITY_PROFILE", density_profile, & @@ -449,26 +447,25 @@ subroutine dumbbell_initialize_sponges(G, GV, US, tv, h_in, depth_tot, param_fil endif enddo ; enddo if (associated(tv%S)) call set_up_ALE_sponge_field(S, G, GV, tv%S, ACSp) - else do j=G%jsc,G%jec ; do i=G%isc,G%iec eta(i,j,1) = 0.0 do k=2,nz - eta(i,j,k) = eta(i,j,k-1)-h_in(i,j,k-1) - enddo - eta(i,j,nz+1) = -depth_tot(i,j) + eta(i,j,k) = eta(i,j,k-1)-h_in(i,j,k-1) + enddo + eta(i,j,nz+1) = -depth_tot(i,j) do k=1,nz - S(i,j,k)= tv%S(i,j,k) + S(i,j,k)= tv%S(i,j,k) enddo - enddo ; enddo - - ! This call sets up the damping rates and interface heights. - ! This sets the inverse damping timescale fields in the sponges. ! - call initialize_sponge(Idamp, eta, G, param_file, CSp, GV) - - ! The remaining calls to set_up_sponge_field can be in any order. ! - if ( associated(tv%S) ) call set_up_sponge_field(S, tv%S, G, GV, nz, CSp) - endif + enddo ; enddo + + ! This call sets up the damping rates and interface heights. + ! This sets the inverse damping timescale fields in the sponges. ! + call initialize_sponge(Idamp, eta, G, param_file, CSp, GV) + + ! The remaining calls to set_up_sponge_field can be in any order. ! + if ( associated(tv%S) ) call set_up_sponge_field(S, tv%S, G, GV, nz, CSp) + endif end subroutine dumbbell_initialize_sponges From 16375ad48bdebe192b3436ea1919e818f73a8455 Mon Sep 17 00:00:00 2001 From: WenhaoChen89 <96131003+WenhaoChen89@users.noreply.github.com> Date: Thu, 7 Jul 2022 17:57:45 -0400 Subject: [PATCH 6/9] Update dumbbell_initialization.F90 --- src/user/dumbbell_initialization.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/user/dumbbell_initialization.F90 b/src/user/dumbbell_initialization.F90 index 9390f26dca..7162e4bc81 100644 --- a/src/user/dumbbell_initialization.F90 +++ b/src/user/dumbbell_initialization.F90 @@ -133,7 +133,7 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE", verticalCoordinate, & default=DEFAULT_COORDINATE_MODE, do_not_log=just_read) call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) - if(.not. use_ALE) verticalCoordinate = "LAYER" + if(.not. use_ALE) verticalCoordinate = "LAYER" ! WARNING: this routine specifies the interface heights so that the last layer ! is vanished, even at maximum depth. In order to have a uniform @@ -162,7 +162,7 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, endif eta1D(1) = 0.0 eta1D(nz+1) = -depth_tot(i,j) - if (x<0.0) then + if (x<0.0) then do k=nz,2, -1 eta1D(k) = eta1D(k+1) + min_thickness enddo From 3e56e422a2c394f6f1f6f56d8a85fa7b006197f3 Mon Sep 17 00:00:00 2001 From: WenhaoChen89 <96131003+WenhaoChen89@users.noreply.github.com> Date: Thu, 7 Jul 2022 18:12:47 -0400 Subject: [PATCH 7/9] Update dumbbell_initialization.F90 --- src/user/dumbbell_initialization.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/user/dumbbell_initialization.F90 b/src/user/dumbbell_initialization.F90 index 7162e4bc81..54e5862ea7 100644 --- a/src/user/dumbbell_initialization.F90 +++ b/src/user/dumbbell_initialization.F90 @@ -458,11 +458,11 @@ subroutine dumbbell_initialize_sponges(G, GV, US, tv, h_in, depth_tot, param_fil S(i,j,k)= tv%S(i,j,k) enddo enddo ; enddo - + ! This call sets up the damping rates and interface heights. ! This sets the inverse damping timescale fields in the sponges. ! call initialize_sponge(Idamp, eta, G, param_file, CSp, GV) - + ! The remaining calls to set_up_sponge_field can be in any order. ! if ( associated(tv%S) ) call set_up_sponge_field(S, tv%S, G, GV, nz, CSp) endif From 9de382d75ddc85d4b5fa7f0ff989fd5b782d9398 Mon Sep 17 00:00:00 2001 From: WenhaoChen89 <96131003+WenhaoChen89@users.noreply.github.com> Date: Sun, 10 Jul 2022 14:59:48 -0400 Subject: [PATCH 8/9] Update dumbbell_initialization.F90 --- src/user/dumbbell_initialization.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/user/dumbbell_initialization.F90 b/src/user/dumbbell_initialization.F90 index 54e5862ea7..e3f0c6dad3 100644 --- a/src/user/dumbbell_initialization.F90 +++ b/src/user/dumbbell_initialization.F90 @@ -133,7 +133,7 @@ subroutine dumbbell_initialize_thickness ( h, depth_tot, G, GV, US, param_file, call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE", verticalCoordinate, & default=DEFAULT_COORDINATE_MODE, do_not_log=just_read) call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) - if(.not. use_ALE) verticalCoordinate = "LAYER" + if (.not. use_ALE) verticalCoordinate = "LAYER" ! WARNING: this routine specifies the interface heights so that the last layer ! is vanished, even at maximum depth. In order to have a uniform @@ -274,7 +274,7 @@ subroutine dumbbell_initialize_temperature_salinity ( T, S, h, G, GV, US, param_ ! layer mode call get_param(param_file, mdl, "USE_REGRIDDING", use_ALE, do_not_log = .true.) - if(.not. use_ALE) call MOM_error(FATAL, "dumbbell_initialize_temperature_salinity: "//& + if (.not. use_ALE) call MOM_error(FATAL, "dumbbell_initialize_temperature_salinity: "//& "Please use 'fit' for 'TS_CONFIG' in the LAYER mode.") call get_param(param_file, mdl, "REGRIDDING_COORDINATE_MODE", verticalCoordinate, & From 433ee56da89c59b8ad5af1f6e12c29af13acc5fd Mon Sep 17 00:00:00 2001 From: WenhaoChen89 <96131003+WenhaoChen89@users.noreply.github.com> Date: Fri, 15 Jul 2022 12:03:26 -0400 Subject: [PATCH 9/9] Update dumbbell_initialization.F90 --- src/user/dumbbell_initialization.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/user/dumbbell_initialization.F90 b/src/user/dumbbell_initialization.F90 index e3f0c6dad3..570e638465 100644 --- a/src/user/dumbbell_initialization.F90 +++ b/src/user/dumbbell_initialization.F90 @@ -451,7 +451,7 @@ subroutine dumbbell_initialize_sponges(G, GV, US, tv, h_in, depth_tot, param_fil do j=G%jsc,G%jec ; do i=G%isc,G%iec eta(i,j,1) = 0.0 do k=2,nz - eta(i,j,k) = eta(i,j,k-1)-h_in(i,j,k-1) + eta(i,j,k) = eta(i,j,k-1)- GV%H_to_Z * h_in(i,j,k-1) enddo eta(i,j,nz+1) = -depth_tot(i,j) do k=1,nz