From e2c81e91690b7b9511ee439ebc05164b2c39b72f Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Wed, 8 Dec 2021 15:13:51 -0500 Subject: [PATCH] +Rescale forcing arguments and revise ctrl_forcing (#25) * +Rescale forcing arguments and revise ctrl_forcing This commit revisits the units of the input arguments to the various ocean- only surfaces forcing routines, including: - Rescaled the units of the time intervals passed to the various forcing routines to [T ~> s] - Applied dimensional scaling to MOM_controlled_forcing.F90. This code is not yet in active use, so these changes can not change answers, but it is now much closer to compliance with modern MOM6 standards, including improved documentation, and could be ready to try without too much more effort. - Documented the remaining real variables in benchmark_initialization.F90, along with their units. All answers are bitwise identical, but there are changes to the units of some arguments in public interfaces. --- .../FMS_cap/MOM_surface_forcing_gfdl.F90 | 8 +- .../solo_driver/MESO_surface_forcing.F90 | 6 +- .../solo_driver/MOM_surface_forcing.F90 | 20 +- .../solo_driver/user_surface_forcing.F90 | 6 +- src/user/BFB_surface_forcing.F90 | 7 +- src/user/MOM_controlled_forcing.F90 | 347 +++++++++++------- src/user/benchmark_initialization.F90 | 53 +-- src/user/dumbbell_surface_forcing.F90 | 17 +- src/user/user_revise_forcing.F90 | 8 +- 9 files changed, 281 insertions(+), 191 deletions(-) diff --git a/config_src/drivers/FMS_cap/MOM_surface_forcing_gfdl.F90 b/config_src/drivers/FMS_cap/MOM_surface_forcing_gfdl.F90 index 84beb3fcf4..d7c483ce49 100644 --- a/config_src/drivers/FMS_cap/MOM_surface_forcing_gfdl.F90 +++ b/config_src/drivers/FMS_cap/MOM_surface_forcing_gfdl.F90 @@ -404,7 +404,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G, delta_sst = data_restore(i,j)- sfc_state%SST(i,j) delta_sst = sign(1.0,delta_sst)*min(abs(delta_sst),CS%max_delta_trestore) fluxes%heat_added(i,j) = G%mask2dT(i,j) * CS%trestore_mask(i,j) * & - rhoXcp * delta_sst * CS%Flux_const_temp ! W m-2 + rhoXcp * delta_sst * CS%Flux_const_temp ! [Q R Z T-1 ~> W m-2] enddo ; enddo endif @@ -568,8 +568,8 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G, !#CTRL# SSS_anom(i,j) = sfc_state%SSS(i,j) - CS%S_Restore(i,j) !#CTRL# SSS_mean(i,j) = 0.5*(sfc_state%SSS(i,j) + CS%S_Restore(i,j)) !#CTRL# enddo ; enddo -!#CTRL# call apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, fluxes%heat_restore, & -!#CTRL# fluxes%vprec, day, dt, G, CS%ctrl_forcing_CSp) +!#CTRL# call apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, fluxes%heat_added, & +!#CTRL# fluxes%vprec, day, US%s_to_T*valid_time, G, US, CS%ctrl_forcing_CSp) !#CTRL# endif ! adjust the NET fresh-water flux to zero, if flagged @@ -1601,7 +1601,7 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, wind_stagger) endif endif -!#CTRL# call controlled_forcing_init(Time, G, param_file, diag, CS%ctrl_forcing_CSp) +!#CTRL# call controlled_forcing_init(Time, G, US, param_file, diag, CS%ctrl_forcing_CSp) call user_revise_forcing_init(param_file, CS%urf_CS) diff --git a/config_src/drivers/solo_driver/MESO_surface_forcing.F90 b/config_src/drivers/solo_driver/MESO_surface_forcing.F90 index f39dff2a0b..300c736802 100644 --- a/config_src/drivers/solo_driver/MESO_surface_forcing.F90 +++ b/config_src/drivers/solo_driver/MESO_surface_forcing.F90 @@ -61,7 +61,7 @@ subroutine MESO_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) type(forcing), intent(inout) :: fluxes !< A structure containing thermodynamic forcing fields type(time_type), intent(in) :: day !< The time of the fluxes real, intent(in) :: dt !< The amount of time over which - !! the fluxes apply [s] + !! the fluxes apply [T ~> s] type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(MESO_surface_forcing_CS), pointer :: CS !< A pointer to the control structure returned by @@ -215,8 +215,8 @@ subroutine MESO_surface_forcing_init(Time, G, US, param_file, diag, CS) type(MESO_surface_forcing_CS), pointer :: CS !< A pointer that is set to point to the !! control structure for this module -! This include declares and sets the variable "version". -#include "version_variable.h" + ! This include declares and sets the variable "version". +# include "version_variable.h" character(len=40) :: mdl = "MESO_surface_forcing" ! This module's name. if (associated(CS)) then diff --git a/config_src/drivers/solo_driver/MOM_surface_forcing.F90 b/config_src/drivers/solo_driver/MOM_surface_forcing.F90 index bf3d517b3d..7c26c2f194 100644 --- a/config_src/drivers/solo_driver/MOM_surface_forcing.F90 +++ b/config_src/drivers/solo_driver/MOM_surface_forcing.F90 @@ -233,7 +233,7 @@ subroutine set_forcing(sfc_state, forces, fluxes, day_start, day_interval, G, US type(surface_forcing_CS), pointer :: CS !< pointer to control structure returned by !! a previous surface_forcing_init call ! Local variables - real :: dt ! length of time over which fluxes applied [s] + real :: dt ! length of time over which fluxes applied [T ~> s] type(time_type) :: day_center ! central time of the fluxes. integer :: isd, ied, jsd, jed isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed @@ -242,7 +242,7 @@ subroutine set_forcing(sfc_state, forces, fluxes, day_start, day_interval, G, US call callTree_enter("set_forcing, MOM_surface_forcing.F90") day_center = day_start + day_interval/2 - dt = time_type_to_real(day_interval) + dt = US%s_to_T * time_type_to_real(day_interval) if (CS%first_call_set_forcing) then ! Allocate memory for the mechanical and thermodynamic forcing fields. @@ -899,7 +899,7 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS) type(forcing), intent(inout) :: fluxes !< A structure containing thermodynamic forcing fields type(time_type), intent(in) :: day !< The time of the fluxes real, intent(in) :: dt !< The amount of time over which - !! the fluxes apply [s] + !! the fluxes apply [T ~> s] type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(surface_forcing_CS), pointer :: CS !< pointer to control structure returned by @@ -1162,7 +1162,7 @@ subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS) !#CTRL# SSS_mean(i,j) = 0.5*(sfc_state%SSS(i,j) + CS%S_Restore(i,j)) !#CTRL# enddo ; enddo !#CTRL# call apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, fluxes%heat_added, & -!#CTRL# fluxes%vprec, day, dt, G, CS%ctrl_forcing_CSp) +!#CTRL# fluxes%vprec, day, dt, G, US, CS%ctrl_forcing_CSp) !#CTRL# endif call callTree_leave("buoyancy_forcing_from_files") @@ -1175,7 +1175,7 @@ subroutine buoyancy_forcing_from_data_override(sfc_state, fluxes, day, dt, G, US type(forcing), intent(inout) :: fluxes !< A structure containing thermodynamic forcing fields type(time_type), intent(in) :: day !< The time of the fluxes real, intent(in) :: dt !< The amount of time over which - !! the fluxes apply [s] + !! the fluxes apply [T ~> s] type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(surface_forcing_CS), pointer :: CS !< pointer to control structure returned by @@ -1289,7 +1289,7 @@ subroutine buoyancy_forcing_from_data_override(sfc_state, fluxes, day, dt, G, US !#CTRL# SSS_mean(i,j) = 0.5*(sfc_state%SSS(i,j) + CS%S_Restore(i,j)) !#CTRL# enddo ; enddo !#CTRL# call apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, fluxes%heat_added, & -!#CTRL# fluxes%vprec, day, dt, G, CS%ctrl_forcing_CSp) +!#CTRL# fluxes%vprec, day, US%T_to_s*dt, G, US, CS%ctrl_forcing_CSp) !#CTRL# endif call callTree_leave("buoyancy_forcing_from_data_override") @@ -1302,7 +1302,7 @@ subroutine buoyancy_forcing_zero(sfc_state, fluxes, day, dt, G, CS) type(forcing), intent(inout) :: fluxes !< A structure containing thermodynamic forcing fields type(time_type), intent(in) :: day !< The time of the fluxes real, intent(in) :: dt !< The amount of time over which - !! the fluxes apply [s] + !! the fluxes apply [T ~> s] type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure type(surface_forcing_CS), pointer :: CS !< pointer to control structure returned by !! a previous surface_forcing_init call @@ -1345,7 +1345,7 @@ subroutine buoyancy_forcing_const(sfc_state, fluxes, day, dt, G, US, CS) type(forcing), intent(inout) :: fluxes !< A structure containing thermodynamic forcing fields type(time_type), intent(in) :: day !< The time of the fluxes real, intent(in) :: dt !< The amount of time over which - !! the fluxes apply [s] + !! the fluxes apply [T ~> s] type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(surface_forcing_CS), pointer :: CS !< pointer to control structure returned by @@ -1388,7 +1388,7 @@ subroutine buoyancy_forcing_linear(sfc_state, fluxes, day, dt, G, US, CS) type(forcing), intent(inout) :: fluxes !< A structure containing thermodynamic forcing fields type(time_type), intent(in) :: day !< The time of the fluxes real, intent(in) :: dt !< The amount of time over which - !! the fluxes apply [s] + !! the fluxes apply [T ~> s] type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(surface_forcing_CS), pointer :: CS !< pointer to control structure returned by @@ -1908,7 +1908,7 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, tracer_flow_C if (trim(CS%wind_config) == "file") & CS%wind_nlev = num_timelevels(CS%wind_file, CS%stress_x_var, min_dims=3) -!#CTRL# call controlled_forcing_init(Time, G, param_file, diag, CS%ctrl_forcing_CSp) +!#CTRL# call controlled_forcing_init(Time, G, US, param_file, diag, CS%ctrl_forcing_CSp) call user_revise_forcing_init(param_file, CS%urf_CS) diff --git a/config_src/drivers/solo_driver/user_surface_forcing.F90 b/config_src/drivers/solo_driver/user_surface_forcing.F90 index 960cddf2ac..0af6b126e1 100644 --- a/config_src/drivers/solo_driver/user_surface_forcing.F90 +++ b/config_src/drivers/solo_driver/user_surface_forcing.F90 @@ -104,7 +104,7 @@ subroutine USER_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) type(forcing), intent(inout) :: fluxes !< A structure containing thermodynamic forcing fields type(time_type), intent(in) :: day !< The time of the fluxes real, intent(in) :: dt !< The amount of time over which - !! the fluxes apply [s] + !! the fluxes apply [T ~> s] type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(user_surface_forcing_CS), pointer :: CS !< A pointer to the control structure returned @@ -242,8 +242,8 @@ subroutine USER_surface_forcing_init(Time, G, US, param_file, diag, CS) type(user_surface_forcing_CS), pointer :: CS !< A pointer that is set to point to !! the control structure for this module -! This include declares and sets the variable "version". -#include "version_variable.h" + ! This include declares and sets the variable "version". +# include "version_variable.h" character(len=40) :: mdl = "user_surface_forcing" ! This module's name. if (associated(CS)) then diff --git a/src/user/BFB_surface_forcing.F90 b/src/user/BFB_surface_forcing.F90 index 6b17d64697..6214f2d095 100644 --- a/src/user/BFB_surface_forcing.F90 +++ b/src/user/BFB_surface_forcing.F90 @@ -52,7 +52,7 @@ subroutine BFB_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) !! have NULL ptrs. type(time_type), intent(in) :: day !< Time of the fluxes. real, intent(in) :: dt !< The amount of time over which - !! the fluxes apply [s] + !! the fluxes apply [T ~> s] type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(BFB_surface_forcing_CS), pointer :: CS !< A pointer to the control structure @@ -177,8 +177,9 @@ subroutine BFB_surface_forcing_init(Time, G, US, param_file, diag, CS) type(diag_ctrl), target, intent(in) :: diag !< A structure that is used to !! regulate diagnostic output. type(BFB_surface_forcing_CS), pointer :: CS !< A pointer to the control structure for this module -! This include declares and sets the variable "version". -#include "version_variable.h" + + ! This include declares and sets the variable "version". +# include "version_variable.h" character(len=40) :: mdl = "BFB_surface_forcing" ! This module's name. if (associated(CS)) then diff --git a/src/user/MOM_controlled_forcing.F90 b/src/user/MOM_controlled_forcing.F90 index 4d44e580e0..f783a271a6 100644 --- a/src/user/MOM_controlled_forcing.F90 +++ b/src/user/MOM_controlled_forcing.F90 @@ -8,20 +8,19 @@ module MOM_controlled_forcing ! This file is part of MOM6. See LICENSE.md for the license. -use MOM_diag_mediator, only : post_data, query_averaging_enabled +use MOM_diag_mediator, only : post_data, query_averaging_enabled, enable_averages, disable_averaging use MOM_diag_mediator, only : register_diag_field, diag_ctrl, safe_alloc_ptr -use MOM_domains, only : pass_var, pass_vector, AGRID, To_South, To_West, To_All +use MOM_domains, only : pass_var, pass_vector, AGRID, To_South, To_West, To_All use MOM_error_handler, only : MOM_error, FATAL, WARNING, MOM_mesg, is_root_pe -use MOM_file_parser, only : read_param, get_param, log_param, log_version, param_file_type -use MOM_forcing_type, only : forcing -use MOM_grid, only : ocean_grid_type -use MOM_io, only : vardesc, var_desc -use MOM_restart, only : register_restart_field, MOM_restart_CS -use MOM_time_manager, only : time_type, operator(+), operator(/), operator(-) -use MOM_time_manager, only : get_date, set_date -use MOM_time_manager, only : time_type_to_real, real_to_time -use MOM_unit_scaling, only : unit_scale_type -use MOM_variables, only : surface +use MOM_file_parser, only : read_param, get_param, log_param, log_version, param_file_type +use MOM_forcing_type, only : forcing +use MOM_grid, only : ocean_grid_type +use MOM_restart, only : register_restart_field, MOM_restart_CS +use MOM_time_manager, only : time_type, operator(+), operator(/), operator(-) +use MOM_time_manager, only : get_date, set_date +use MOM_time_manager, only : time_type_to_real, real_to_time +use MOM_unit_scaling, only : unit_scale_type +use MOM_variables, only : surface implicit none ; private @@ -32,89 +31,103 @@ module MOM_controlled_forcing !> Control structure for MOM_controlled_forcing type, public :: ctrl_forcing_CS ; private - logical :: use_temperature !< If true, temperature and salinity are used as - !! state variables. - logical :: do_integrated !< If true, use time-integrated anomalies to control - !! the surface state. - integer :: num_cycle !< The number of elements in the forcing cycle. - real :: heat_int_rate !< The rate at which heating anomalies accumulate [s-1]. - real :: prec_int_rate !< The rate at which precipitation anomalies accumulate [s-1]. - real :: heat_cyc_rate !< The rate at which cyclical heating anomaliess - !! accumulate [s-1]. - real :: prec_cyc_rate !< The rate at which cyclical precipitation anomaliess - !! accumulate [s-1]. + logical :: use_temperature !< If true, temperature and salinity are used as state variables. + logical :: do_integrated !< If true, use time-integrated anomalies to control the surface state. + integer :: num_cycle !< The number of elements in the forcing cycle. + real :: heat_int_rate !< The rate at which heating anomalies accumulate [T-1 ~> s-1] + real :: prec_int_rate !< The rate at which precipitation anomalies accumulate [T-1 ~> s-1] + real :: heat_cyc_rate !< The rate at which cyclical heating anomalies accumulate [T-1 ~> s-1] + real :: prec_cyc_rate !< The rate at which cyclical precipitation anomalies + !! accumulate [T-1 ~> s-1] real :: Len2 !< The square of the length scale over which the anomalies - !! are smoothed via a Laplacian filter [m2]. + !! are smoothed via a Laplacian filter [L2 ~> m2] real :: lam_heat !< A constant of proportionality between SST anomalies - !! and heat fluxes [W m-2 degC-1]. + !! and heat fluxes [Q R Z T-1 degC-1 ~> W m-2 degC-1] real :: lam_prec !< A constant of proportionality between SSS anomalies - !! (normalised by mean SSS) and precipitation [kg m-2]. + !! (normalised by mean SSS) and precipitation [R Z T-1 ~> kg m-2 s-1] real :: lam_cyc_heat !< A constant of proportionality between cyclical SST - !! anomalies and corrective heat fluxes [W m-2 degC-1]. + !! anomalies and corrective heat fluxes [W m-2 degC-1] real :: lam_cyc_prec !< A constant of proportionality between cyclical SSS !! anomalies (normalised by mean SSS) and corrective - !! precipitation [kg m-2]. + !! precipitation [R Z T-1 ~> kg m-2 s-1] - !>@{ Pointers for data. - !! \todo Needs more complete documentation. - real, pointer, dimension(:) :: & - avg_time => NULL() real, pointer, dimension(:,:) :: & - heat_0 => NULL(), & - precip_0 => NULL() + heat_0 => NULL(), & !< The non-periodic integrative corrective heat flux that has been + !! evolved to control mean SST anomalies [Q R Z T-1 ~> W m-2] + precip_0 => NULL() !< The non-periodic integrative corrective precipitation that has been + !! evolved to control mean SSS anomalies [R Z T-1 ~> kg m-2 s-1] + + ! The final dimension of each of the six variables that follow is for the periodic bins. + real, pointer, dimension(:,:,:) :: & + heat_cyc => NULL(), & !< The periodic integrative corrective heat flux that has been evolved + !! to control periodic (seasonal) SST anomalies [Q R Z T-1 ~> W m-2]. + !! The third dimension is the periodic bins. + precip_cyc => NULL() !< The non-periodic integrative corrective precipitation that has been + !! evolved to control periodic (seasonal) SSS anomalies [R Z T-1 ~> kg m-2 s-1]. + !! The third dimension is the periodic bins. + real, pointer, dimension(:) :: & + avg_time => NULL() !< The accumulated averaging time in each part of the cycle [T ~> s] or + !! a negative value to indicate that the variables like avg_SST_anom are + !! the actual averages, and not time integrals. + !! The dimension is the periodic bins. real, pointer, dimension(:,:,:) :: & - heat_cyc => NULL(), & - precip_cyc => NULL(), & - avg_SST_anom => NULL(), & - avg_SSS_anom => NULL(), & - avg_SSS => NULL() - !>@} + avg_SST_anom => NULL(), & !< The time-averaged periodic sea surface temperature anomalies [degC], + !! or (at some points in the code), the time-integrated periodic + !! temperature anomalies [T degC ~> s degC]. + !! The third dimension is the periodic bins. + avg_SSS_anom => NULL(), & !< The time-averaged periodic sea surface salinity anomalies [ppt], + !! or (at some points in the code), the time-integrated periodic + !! salinity anomalies [T ppt ~> s ppt]. + !! The third dimension is the periodic bins. + avg_SSS => NULL() !< The time-averaged periodic sea surface salinities [ppt], or (at + !! some points in the code), the time-integrated periodic + !! salinities [T ppt ~> s ppt]. + !! The third dimension is the periodic bins. + type(diag_ctrl), pointer :: diag => NULL() !< A structure that is used to !! regulate the timing of diagnostic output. - integer :: id_heat_0 = -1 !< Diagnostic handle + integer :: id_heat_0 = -1 !< Diagnostic handle for the steady heat flux + integer :: id_prec_0 = -1 !< Diagnostic handle for the steady precipitation end type ctrl_forcing_CS contains -!> This subroutine calls any of the other subroutines in this file -!! that are needed to specify the current surface forcing fields. +!> This subroutine determines corrective surface forcing fields using simple control theory. subroutine apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, virt_heat, virt_precip, & day_start, dt, G, US, CS) - type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure. - real, dimension(SZI_(G),SZJ_(G)), intent(in) :: SST_anom !< The sea surface temperature - !! anomalies [degC]. - real, dimension(SZI_(G),SZJ_(G)), intent(in) :: SSS_anom !< The sea surface salinity - !! anomlies [ppt]. - real, dimension(SZI_(G),SZJ_(G)), intent(in) :: SSS_mean !< The mean sea surface - !! salinity [ppt]. + type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure + real, dimension(SZI_(G),SZJ_(G)), intent(in) :: SST_anom !< The sea surface temperature anomalies [degC] + real, dimension(SZI_(G),SZJ_(G)), intent(in) :: SSS_anom !< The sea surface salinity anomlies [ppt] + real, dimension(SZI_(G),SZJ_(G)), intent(in) :: SSS_mean !< The mean sea surface salinity [ppt] real, dimension(SZI_(G),SZJ_(G)), intent(inout) :: virt_heat !< Virtual (corrective) heat - !! fluxes that are augmented - !! in this subroutine [W m-2]. + !! fluxes that are augmented in this + !! subroutine [Q R Z T-1 ~> W m-2] real, dimension(SZI_(G),SZJ_(G)), intent(inout) :: virt_precip !< Virtual (corrective) - !! precipitation fluxes that - !! are augmented in this - !! subroutine [kg m-2 s-1]. - type(time_type), intent(in) :: day_start !< Start time of the fluxes. - real, intent(in) :: dt !< Length of time over which these - !! fluxes will be applied [s]. - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type - type(ctrl_forcing_CS), pointer :: CS !< A pointer to the control structure - !! returned by a previous call to - !! ctrl_forcing_init. -! + !! precipitation fluxes that are augmented + !! in this subroutine [R Z T-1 ~> kg m-2 s-1] + type(time_type), intent(in) :: day_start !< Start time of the fluxes. + real, intent(in) :: dt !< Length of time over which these fluxes + !! will be applied [T ~> s] + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + type(ctrl_forcing_CS), pointer :: CS !< A pointer to the control structure returned + !! by a previous call to ctrl_forcing_init. + + ! Local variables real, dimension(SZIB_(G),SZJ_(G)) :: & - flux_heat_x, & - flux_prec_x + flux_heat_x, & ! Zonal smoothing flux of the virtual heat fluxes [L2 Q R Z T-1 ~> W] + flux_prec_x ! Zonal smoothing flux of the virtual precipitation [L2 R Z T-1 ~> kg s-1] real, dimension(SZI_(G),SZJB_(G)) :: & - flux_heat_y, & - flux_prec_y + flux_heat_y, & ! Meridional smoothing flux of the virtual heat fluxes [L2 Q R Z T-1 ~> W] + flux_prec_y ! Meridional smoothing flux of the virtual precipitation [L2 R Z T-1 ~> kg s-1] type(time_type) :: day_end - real :: coef ! A heat-flux coefficient [m2]. - real :: mr_st, mr_end, mr_mid, mr_prev, mr_next - real :: dt_wt, dt_heat_rate, dt_prec_rate - real :: dt1_heat_rate, dt1_prec_rate, dt2_heat_rate, dt2_prec_rate - real :: wt_per1, wt_st, wt_end, wt_mid - integer :: m_st, m_end, m_mid, m_u1, m_u2, m_u3 + real :: coef ! A heat-flux coefficient [L2 ~> m2] + real :: mr_st, mr_end, mr_mid ! Position of various times in the periodic cycle [nondim] + real :: mr_prev, mr_next ! Position of various times in the periodic cycle [nondim] + real :: dt_wt ! The timestep times a fractional weight used to accumulate averages [T ~> s] + real :: dt_heat_rate, dt_prec_rate ! Timestep times the flux accumulation rate [nondim] + real :: dt1_heat_rate, dt1_prec_rate, dt2_heat_rate, dt2_prec_rate ! [nondim] + real :: wt_per1, wt_st, wt_end, wt_mid ! Averaging weights [nondim] + integer :: m_st, m_end, m_mid, m_u1, m_u2, m_u3 ! Indices (nominally months) in the periodic cycle integer :: yr, mon, day, hr, min, sec integer :: i, j, is, ie, js, je @@ -123,7 +136,7 @@ subroutine apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, virt_heat, virt_prec if (.not.associated(CS)) return if ((CS%num_cycle <= 0) .and. (.not.CS%do_integrated)) return - day_end = day_start + real_to_time(dt) + day_end = day_start + real_to_time(US%T_to_s*dt) do j=js,je ; do i=is,ie virt_heat(i,j) = 0.0 ; virt_precip(i,j) = 0.0 @@ -148,12 +161,12 @@ subroutine apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, virt_heat, virt_prec do j=js,je ; do i=is,ie CS%heat_0(i,j) = CS%heat_0(i,j) + dt_heat_rate * ( & -CS%lam_heat*G%mask2dT(i,j)*SST_anom(i,j) + & - (US%m_to_L**2*G%IareaT(i,j) * ((flux_heat_x(I-1,j) - flux_heat_x(I,j)) + & + (G%IareaT(i,j) * ((flux_heat_x(I-1,j) - flux_heat_x(I,j)) + & (flux_heat_y(i,J-1) - flux_heat_y(i,J))) ) ) CS%precip_0(i,j) = CS%precip_0(i,j) + dt_prec_rate * ( & CS%lam_prec * G%mask2dT(i,j)*(SSS_anom(i,j) / SSS_mean(i,j)) + & - (US%m_to_L**2*G%IareaT(i,j) * ((flux_prec_x(I-1,j) - flux_prec_x(I,j)) + & + (G%IareaT(i,j) * ((flux_prec_x(I-1,j) - flux_prec_x(I,j)) + & (flux_prec_y(i,J-1) - flux_prec_y(i,J))) ) ) virt_heat(i,j) = virt_heat(i,j) + CS%heat_0(i,j) @@ -257,6 +270,7 @@ subroutine apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, virt_heat, virt_prec ! Accumulate the average anomalies for this period. dt_wt = wt_per1 * dt CS%avg_time(m_mid) = CS%avg_time(m_mid) + dt_wt + ! These loops temporarily change the units of the CS%avg_ variables to [degC s] or [ppt s]. do j=js,je ; do i=is,ie CS%avg_SST_anom(i,j,m_mid) = CS%avg_SST_anom(i,j,m_mid) + & dt_wt * G%mask2dT(i,j) * SST_anom(i,j) @@ -281,6 +295,7 @@ subroutine apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, virt_heat, virt_prec m_u2 = periodic_int(m_st - 3.0, CS%num_cycle) m_u3 = periodic_int(m_st - 2.0, CS%num_cycle) + ! These loops restore the units of the CS%avg variables to [degC] or [ppt] if (CS%avg_time(m_u1) > 0.0) then do j=js,je ; do i=is,ie CS%avg_SST_anom(i,j,m_u1) = CS%avg_SST_anom(i,j,m_u1) / CS%avg_time(m_u1) @@ -332,13 +347,13 @@ subroutine apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, virt_heat, virt_prec do j=js,je ; do i=is,ie CS%heat_cyc(i,j,m_u1) = CS%heat_cyc(i,j,m_u1) + dt1_heat_rate * ( & -CS%lam_cyc_heat*(CS%avg_SST_anom(i,j,m_u2) - CS%avg_SST_anom(i,j,m_u1)) + & - (US%m_to_L**2*G%IareaT(i,j) * ((flux_heat_x(I-1,j) - flux_heat_x(I,j)) + & + (G%IareaT(i,j) * ((flux_heat_x(I-1,j) - flux_heat_x(I,j)) + & (flux_heat_y(i,J-1) - flux_heat_y(i,J))) ) ) CS%precip_cyc(i,j,m_u1) = CS%precip_cyc(i,j,m_u1) + dt1_prec_rate * ( & - CS%lam_cyc_prec * (CS%avg_SSS_anom(i,j,m_u2) - CS%avg_SSS_anom(i,j,m_u1)) / & + CS%lam_prec * (CS%avg_SSS_anom(i,j,m_u2) - CS%avg_SSS_anom(i,j,m_u1)) / & (0.5*(CS%avg_SSS(i,j,m_u2) + CS%avg_SSS(i,j,m_u1))) + & - (US%m_to_L**2*G%IareaT(i,j) * ((flux_prec_x(I-1,j) - flux_prec_x(I,j)) + & + (G%IareaT(i,j) * ((flux_prec_x(I-1,j) - flux_prec_x(I,j)) + & (flux_prec_y(i,J-1) - flux_prec_y(i,J))) ) ) enddo ; enddo endif @@ -357,19 +372,26 @@ subroutine apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, virt_heat, virt_prec do j=js,je ; do i=is,ie CS%heat_cyc(i,j,m_u2) = CS%heat_cyc(i,j,m_u2) + dt1_heat_rate * ( & -CS%lam_cyc_heat*(CS%avg_SST_anom(i,j,m_u3) - CS%avg_SST_anom(i,j,m_u2)) + & - (US%m_to_L**2*G%IareaT(i,j) * ((flux_heat_x(I-1,j) - flux_heat_x(I,j)) + & + (G%IareaT(i,j) * ((flux_heat_x(I-1,j) - flux_heat_x(I,j)) + & (flux_heat_y(i,J-1) - flux_heat_y(i,J))) ) ) CS%precip_cyc(i,j,m_u2) = CS%precip_cyc(i,j,m_u2) + dt1_prec_rate * ( & - CS%lam_cyc_prec * (CS%avg_SSS_anom(i,j,m_u3) - CS%avg_SSS_anom(i,j,m_u2)) / & + CS%lam_prec * (CS%avg_SSS_anom(i,j,m_u3) - CS%avg_SSS_anom(i,j,m_u2)) / & (0.5*(CS%avg_SSS(i,j,m_u3) + CS%avg_SSS(i,j,m_u2))) + & - (US%m_to_L**2*G%IareaT(i,j) * ((flux_prec_x(I-1,j) - flux_prec_x(I,j)) + & + (G%IareaT(i,j) * ((flux_prec_x(I-1,j) - flux_prec_x(I,j)) + & (flux_prec_y(i,J-1) - flux_prec_y(i,J))) ) ) enddo ; enddo endif endif ! (CS%num_cycle > 0) + if (CS%do_integrated .and. ((CS%id_heat_0 > 0) .or. (CS%id_prec_0 > 0))) then + call enable_averages(dt, day_start + real_to_time(US%T_to_s*dt), CS%diag) + if (CS%id_heat_0 > 0) call post_data(CS%id_heat_0, CS%heat_0, CS%diag) + if (CS%id_prec_0 > 0) call post_data(CS%id_prec_0, CS%precip_0, CS%diag) + call disable_averaging(CS%diag) + endif + end subroutine apply_ctrl_forcing !> This function maps rval into an integer in the range from 1 to num_period. @@ -415,7 +437,6 @@ subroutine register_ctrl_forcing_restarts(G, param_file, CS, restart_CS) logical :: controlled, use_temperature character (len=8) :: period_str - type(vardesc) :: vd integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB @@ -443,47 +464,44 @@ subroutine register_ctrl_forcing_restarts(G, param_file, CS, restart_CS) call read_param(param_file, "CTRL_FORCE_NUM_CYCLE", CS%num_cycle) if (CS%do_integrated) then - call safe_alloc_ptr(CS%heat_0,isd,ied,jsd,jed) ; CS%heat_0(:,:) = 0.0 - call safe_alloc_ptr(CS%precip_0,isd,ied,jsd,jed) ; CS%precip_0(:,:) = 0.0 - vd = var_desc("Ctrl_heat","W m-2","Control Integrative Heating",z_grid='1') - call register_restart_field(CS%heat_0, vd, .false., restart_CS) - vd = var_desc("Ctrl_precip","kg m-2 s-1","Control Integrative Precipitation",z_grid='1') - call register_restart_field(CS%precip_0, vd, .false., restart_CS) + allocate(CS%heat_0(isd:ied,jsd:jed), source=0.0) + allocate(CS%precip_0(isd:ied,jsd:jed), source=0.0) + + call register_restart_field(CS%heat_0, "Ctrl_heat", .false., restart_CS, & + longname="Control Integrative Heating", units="W m-2", z_grid='1') + call register_restart_field(CS%precip_0, "Ctrl_precip", .false., restart_CS, & + longname="Control Integrative Precipitation", units="kg m-2 s-1", z_grid='1') endif if (CS%num_cycle > 0) then + allocate(CS%heat_cyc(isd:ied,jsd:jed,CS%num_cycle), source=0.0) + allocate(CS%precip_cyc(isd:ied,jsd:jed,CS%num_cycle), source=0.0) + allocate(CS%avg_time(CS%num_cycle), source=0.0) + allocate(CS%avg_SST_anom(isd:ied,jsd:jed,CS%num_cycle), source=0.0) + allocate(CS%avg_SSS_anom(isd:ied,jsd:jed,CS%num_cycle), source=0.0) + write (period_str, '(i8)') CS%num_cycle period_str = trim('p ')//trim(adjustl(period_str)) - call safe_alloc_ptr(CS%heat_cyc,isd,ied,jsd,jed,CS%num_cycle) ; CS%heat_cyc(:,:,:) = 0.0 - call safe_alloc_ptr(CS%precip_cyc,isd,ied,jsd,jed,CS%num_cycle) ; CS%precip_cyc(:,:,:) = 0.0 - vd = var_desc("Ctrl_heat_cycle", "W m-2","Cyclical Control Heating",& - z_grid='1', t_grid=period_str) - call register_restart_field(CS%heat_cyc, vd, .false., restart_CS) - vd = var_desc("Ctrl_precip_cycle","kg m-2 s-1","Cyclical Control Precipitation", & - z_grid='1', t_grid=period_str) - call register_restart_field(CS%precip_cyc, vd, .false., restart_CS) - - call safe_alloc_ptr(CS%avg_time,CS%num_cycle) ; CS%avg_time(:) = 0.0 - vd = var_desc("avg_time","sec","Cyclical accumulated averaging time", & - '1',z_grid='1',t_grid=period_str) - call register_restart_field(CS%avg_time, vd, .false., restart_CS) - - call safe_alloc_ptr(CS%avg_SST_anom,isd,ied,jsd,jed,CS%num_cycle) ; CS%avg_SST_anom(:,:,:) = 0.0 - call safe_alloc_ptr(CS%avg_SSS_anom,isd,ied,jsd,jed,CS%num_cycle) ; CS%avg_SSS_anom(:,:,:) = 0.0 - vd = var_desc("avg_SST_anom","deg C","Cyclical average SST Anomaly", & - z_grid='1',t_grid=period_str) - call register_restart_field(CS%avg_SST_anom, vd, .false., restart_CS) - vd = var_desc("avg_SSS_anom","g kg-1","Cyclical average SSS Anomaly", & - z_grid='1',t_grid=period_str) - call register_restart_field(CS%avg_SSS_anom, vd, .false., restart_CS) + + call register_restart_field(CS%heat_cyc, "Ctrl_heat_cycle", .false., restart_CS, & + longname="Cyclical Control Heating", units="W m-2", z_grid='1', t_grid=period_str) + call register_restart_field(CS%precip_cyc, "Ctrl_precip_cycle", .false., restart_CS, & + longname="Cyclical Control Precipitation", units="kg m-2 s-1", z_grid='1', t_grid=period_str) + call register_restart_field(CS%avg_time, "avg_time", .false., restart_CS, & + longname="Cyclical accumulated averaging time", units="sec", z_grid='1', t_grid=period_str) + call register_restart_field(CS%avg_SST_anom, "avg_SST_anom", .false., restart_CS, & + longname="Cyclical average SST Anomaly", units="deg C", z_grid='1', t_grid=period_str) + call register_restart_field(CS%avg_SSS_anom, "avg_SSS_anom", .false., restart_CS, & + longname="Cyclical average SSS Anomaly", units="g kg-1", z_grid='1', t_grid=period_str) endif end subroutine register_ctrl_forcing_restarts !> Set up this modules control structure. -subroutine controlled_forcing_init(Time, G, param_file, diag, CS) +subroutine controlled_forcing_init(Time, G, US, param_file, diag, CS) type(time_type), intent(in) :: Time !< The current model time. type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure. + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(param_file_type), intent(in) :: param_file !< A structure indicating the !! open file to parse for model !! parameter values. @@ -491,13 +509,20 @@ subroutine controlled_forcing_init(Time, G, param_file, diag, CS) !! diagnostic output. type(ctrl_forcing_CS), pointer :: CS !< A pointer that is set to point to the !! control structure for this module. - real :: smooth_len + + ! Local variables + real :: smooth_len ! A smoothing lengthscale [L ~> m] + real :: RZ_T_rescale ! Unit conversion factor for precipiation [T kg m-2 s-1 R-1 Z-1 ~> 1] + real :: QRZ_T_rescale ! Unit conversion factor for head fluxes [T W m-2 Q-1 R-1 Z-1 ~> 1] logical :: do_integrated integer :: num_cycle -! This include declares and sets the variable "version". -#include "version_variable.h" + integer :: i, j, isc, iec, jsc, jec, m + ! This include declares and sets the variable "version". +# include "version_variable.h" character(len=40) :: mdl = "MOM_controlled_forcing" ! This module's name. + isc = G%isc ; iec = G%iec ; jsc = G%jsc ; jec = G%jec + ! These should have already been called. ! call read_param(param_file, "CTRL_FORCE_INTEGRATED", CS%do_integrated) ! call read_param(param_file, "CTRL_FORCE_NUM_CYCLE", CS%num_cycle) @@ -523,40 +548,96 @@ subroutine controlled_forcing_init(Time, G, param_file, diag, CS) CS%diag => diag call get_param(param_file, mdl, "CTRL_FORCE_HEAT_INT_RATE", CS%heat_int_rate, & - "The integrated rate at which heat flux anomalies are "//& - "accumulated.", units="s-1", default=0.0) + "The integrated rate at which heat flux anomalies are accumulated.", & + units="s-1", default=0.0, scale=US%T_to_s) call get_param(param_file, mdl, "CTRL_FORCE_PREC_INT_RATE", CS%prec_int_rate, & - "The integrated rate at which precipitation anomalies "//& - "are accumulated.", units="s-1", default=0.0) + "The integrated rate at which precipitation anomalies are accumulated.", & + units="s-1", default=0.0, scale=US%T_to_s) call get_param(param_file, mdl, "CTRL_FORCE_HEAT_CYC_RATE", CS%heat_cyc_rate, & - "The integrated rate at which cyclical heat flux "//& - "anomalies are accumulated.", units="s-1", default=0.0) + "The integrated rate at which cyclical heat flux anomalies are accumulated.", & + units="s-1", default=0.0, scale=US%T_to_s) call get_param(param_file, mdl, "CTRL_FORCE_PREC_CYC_RATE", CS%prec_cyc_rate, & - "The integrated rate at which cyclical precipitation "//& - "anomalies are accumulated.", units="s-1", default=0.0) + "The integrated rate at which cyclical precipitation anomalies are accumulated.", & + units="s-1", default=0.0, scale=US%T_to_s) call get_param(param_file, mdl, "CTRL_FORCE_SMOOTH_LENGTH", smooth_len, & - "The length scales over which controlled forcing "//& - "anomalies are smoothed.", units="m", default=0.0) + "The length scales over which controlled forcing anomalies are smoothed.", & + units="m", default=0.0, scale=US%m_to_L) call get_param(param_file, mdl, "CTRL_FORCE_LAMDA_HEAT", CS%lam_heat, & "A constant of proportionality between SST anomalies "//& - "and controlling heat fluxes", "W m-2 K-1", default=0.0) + "and controlling heat fluxes", & + units="W m-2 K-1", default=0.0, scale=US%W_m2_to_QRZ_T) call get_param(param_file, mdl, "CTRL_FORCE_LAMDA_PREC", CS%lam_prec, & "A constant of proportionality between SSS anomalies "//& "(normalised by mean SSS) and controlling precipitation.", & - "kg m-2", default=0.0) + units="kg m-2 s-1", default=0.0, scale=US%kg_m2s_to_RZ_T) call get_param(param_file, mdl, "CTRL_FORCE_LAMDA_CYC_HEAT", CS%lam_cyc_heat, & "A constant of proportionality between SST anomalies "//& - "and cyclical controlling heat fluxes", "W m-2 K-1", default=0.0) + "and cyclical controlling heat fluxes", & + units="W m-2 K-1", default=0.0, scale=US%W_m2_to_QRZ_T) call get_param(param_file, mdl, "CTRL_FORCE_LAMDA_CYC_PREC", CS%lam_cyc_prec, & "A constant of proportionality between SSS anomalies "//& - "(normalised by mean SSS) and cyclical controlling "//& - "precipitation.", "kg m-2", default=0.0) + "(normalised by mean SSS) and cyclical controlling precipitation.", & + units="kg m-2 s-1", default=0.0, scale=US%kg_m2s_to_RZ_T) CS%Len2 = smooth_len**2 -! ### REPLACE THIS WITH ANY DIAGNOSTICS FROM THIS MODULE. -! CS%id_taux = register_diag_field('ocean_model', 'taux', diag%axesu1, Time, & -! 'Zonal Wind Stress', 'Pascal') + if (CS%do_integrated) then + CS%id_heat_0 = register_diag_field('ocean_model', 'Ctrl_heat', diag%axesT1, Time, & + 'Control Corrective Heating', 'W m-2', conversion=US%QRZ_T_to_W_m2) + CS%id_prec_0 = register_diag_field('ocean_model', 'Ctrl_prec', diag%axesT1, Time, & + 'Control Corrective Precipitation', 'kg m-2 s-1', conversion=US%RZ_T_to_kg_m2s) + endif + + ! Rescale if there are differences between the dimensional scaling of variables in + ! restart files from those in use for this run. + if ((US%J_kg_to_Q_restart*US%kg_m3_to_R_restart*US%m_to_Z_restart*US%s_to_T_restart /= 0.0) .and. & + ((US%J_kg_to_Q * US%kg_m3_to_R * US%m_to_Z * US%s_to_T_restart) /= & + (US%J_kg_to_Q_restart * US%kg_m3_to_R_restart * US%m_to_Z_restart * US%s_to_T)) ) then + ! Redo the scaling of the corrective heat fluxes to [Q R Z T-1 ~> W m-2] + QRZ_T_rescale = (US%J_kg_to_Q * US%kg_m3_to_R * US%m_to_Z * US%s_to_T_restart) / & + (US%J_kg_to_Q_restart * US%kg_m3_to_R_restart * US%m_to_Z_restart * US%s_to_T) + + if (associated(CS%heat_0)) then + do j=jsc,jec ; do i=isc,iec + CS%heat_0(i,j) = QRZ_T_rescale * CS%heat_0(i,j) + enddo ; enddo + endif + + if ((CS%num_cycle > 0) .and. associated(CS%heat_cyc)) then + do m=1,CS%num_cycle ; do j=jsc,jec ; do i=isc,iec + CS%heat_cyc(i,j,m) = QRZ_T_rescale * CS%heat_cyc(i,j,m) + enddo ; enddo ; enddo + endif + endif + + if ((US%kg_m3_to_R_restart * US%m_to_Z_restart * US%s_to_T_restart /= 0.0) .and. & + ((US%kg_m3_to_R * US%m_to_Z * US%s_to_T_restart) /= & + (US%kg_m3_to_R_restart * US%m_to_Z_restart * US%s_to_T)) ) then + ! Redo the scaling of the corrective precipitation to [R Z T-1 ~> kg m-2 s-1] + RZ_T_rescale = (US%kg_m3_to_R * US%m_to_Z * US%s_to_T_restart) / & + (US%kg_m3_to_R_restart * US%m_to_Z_restart * US%s_to_T) + + if (associated(CS%precip_0)) then + do j=jsc,jec ; do i=isc,iec + CS%precip_0(i,j) = RZ_T_rescale * CS%precip_0(i,j) + enddo ; enddo + endif + + if ((CS%num_cycle > 0) .and. associated(CS%precip_cyc)) then + do m=1,CS%num_cycle ; do j=jsc,jec ; do i=isc,iec + CS%precip_cyc(i,j,m) = RZ_T_rescale * CS%precip_cyc(i,j,m) + enddo ; enddo ; enddo + endif + endif + + if ((CS%num_cycle > 0) .and. associated(CS%avg_time) .and. & + ((US%s_to_T_restart /= 0.0) .and. ((US%s_to_T_restart) /= US%s_to_T)) ) then + ! Redo the scaling of the accumulated times to [T ~> s] + do m=1,CS%num_cycle + CS%avg_time(m) = (US%s_to_T / US%s_to_T_restart) * CS%avg_time(m) + enddo + endif + end subroutine controlled_forcing_init diff --git a/src/user/benchmark_initialization.F90 b/src/user/benchmark_initialization.F90 index b955f75a32..30c30f264a 100644 --- a/src/user/benchmark_initialization.F90 +++ b/src/user/benchmark_initialization.F90 @@ -34,18 +34,18 @@ module benchmark_initialization subroutine benchmark_initialize_topography(D, G, param_file, max_depth, US) type(dyn_horgrid_type), intent(in) :: G !< The dynamic horizontal grid type real, dimension(G%isd:G%ied,G%jsd:G%jed), & - intent(out) :: D !< Ocean bottom depth in m or [Z ~> m] if US is present + intent(out) :: D !< Ocean bottom depth in [m] or [Z ~> m] if US is present type(param_file_type), intent(in) :: param_file !< Parameter file structure - real, intent(in) :: max_depth !< Maximum model depth in the units of D + real, intent(in) :: max_depth !< Maximum model depth in the units of D, [m] or [Z ~> m] type(unit_scale_type), optional, intent(in) :: US !< A dimensional unit scaling type ! Local variables - real :: min_depth ! The minimum and maximum depths [Z ~> m]. - real :: PI ! 3.1415926... calculated as 4*atan(1) - real :: D0 ! A constant to make the maximum ! - ! basin depth MAXIMUM_DEPTH. ! - real :: m_to_Z ! A dimensional rescaling factor. - real :: x, y + real :: min_depth ! The minimum basin depth [m] or [Z ~> m] + real :: PI ! 3.1415926... calculated as 4*atan(1) + real :: D0 ! A constant to make the maximum basin depth MAXIMUM_DEPTH [m] or [Z ~> m] + real :: m_to_Z ! A dimensional rescaling factor [Z m-1 ~> 1] + real :: x ! Longitude relative to the domain edge, normalized by its extent [nondim] + real :: y ! Latitude relative to the domain edge, normalized by its extent [nondim] ! This include declares and sets the variable "version". # include "version_variable.h" character(len=40) :: mdl = "benchmark_initialize_topography" ! This subroutine's name. @@ -118,10 +118,13 @@ subroutine benchmark_initialize_thickness(h, depth_tot, G, GV, US, param_file, e real :: a_exp ! The fraction of the overall stratification that is exponential. real :: I_ts, I_md ! Inverse lengthscales [Z-1 ~> m-1]. real :: T_frac ! A ratio of the interface temperature to the range - ! between SST and the bottom temperature. - real :: err, derr_dz ! The error between the profile's temperature and the - ! interface temperature for a given z and its derivative. - real :: pi, z + ! between SST and the bottom temperature [nondim]. + real :: err ! The normalized error between the profile's temperature and the + ! interface temperature for a given z [nondim] + real :: derr_dz ! The derivative of the normalized error between the profile's + ! temperature and the interface temperature with z [Z-1 ~> m-1] + real :: pi ! 3.1415926... calculated as 4*atan(1) + real :: z ! A work variable for the interface position [Z ~> m] ! This include declares and sets the variable "version". # include "version_variable.h" character(len=40) :: mdl = "benchmark_initialize_thickness" ! This subroutine's name. @@ -178,9 +181,10 @@ subroutine benchmark_initialize_thickness(h, depth_tot, G, GV, US, param_file, e do k=1,nz ; e_pert(K) = 0.0 ; enddo ! This sets the initial thickness (in [H ~> m or kg m-2]) of the layers. The thicknesses - ! are set to insure that: 1. each layer is at least Gv%Angstrom_m thick, and - ! 2. the interfaces are where they should be based on the resting depths and interface - ! height perturbations, as long at this doesn't interfere with 1. + ! are set to insure that: + ! 1. each layer is at least GV%Angstrom_H thick, and + ! 2. the interfaces are where they should be based on the resting depths and + ! interface height perturbations, as long at this doesn't interfere with 1. eta1D(nz+1) = -depth_tot(i,j) do k=nz,2,-1 @@ -214,8 +218,8 @@ end subroutine benchmark_initialize_thickness !> Initializes layer temperatures and salinities for benchmark subroutine benchmark_init_temperature_salinity(T, S, G, GV, US, param_file, & eqn_of_state, P_Ref, just_read) - type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure. - type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure. + type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure + type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(out) :: T !< The potential temperature !! that is being initialized [degC] real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(out) :: S !< The salinity that is being @@ -226,19 +230,18 @@ subroutine benchmark_init_temperature_salinity(T, S, G, GV, US, param_file, & !! model parameter values. type(EOS_type), intent(in) :: eqn_of_state !< Equation of state structure real, intent(in) :: P_Ref !< The coordinate-density - !! reference pressure [R L2 T-2 ~> Pa]. + !! reference pressure [R L2 T-2 ~> Pa] logical, intent(in) :: just_read !< If true, this call will only read !! parameters without changing T & S. ! Local variables real :: T0(SZK_(GV)) ! A profile of temperatures [degC] real :: S0(SZK_(GV)) ! A profile of salinities [ppt] - real :: pres(SZK_(GV)) ! Reference pressure [R L2 T-2 ~> Pa]. - real :: drho_dT(SZK_(GV)) ! Derivative of density with temperature [R degC-1 ~> kg m-3 degC-1]. - real :: drho_dS(SZK_(GV)) ! Derivative of density with salinity [R ppt-1 ~> kg m-3 ppt-1]. - real :: rho_guess(SZK_(GV)) ! Potential density at T0 & S0 [R ~> kg m-3]. - real :: PI ! 3.1415926... calculated as 4*atan(1) - real :: SST ! The initial sea surface temperature [degC]. - real :: lat + real :: pres(SZK_(GV)) ! Reference pressure [R L2 T-2 ~> Pa] + real :: drho_dT(SZK_(GV)) ! Derivative of density with temperature [R degC-1 ~> kg m-3 degC-1] + real :: drho_dS(SZK_(GV)) ! Derivative of density with salinity [R ppt-1 ~> kg m-3 ppt-1] + real :: rho_guess(SZK_(GV)) ! Potential density at T0 & S0 [R ~> kg m-3] + real :: PI ! 3.1415926... calculated as 4*atan(1) + real :: SST ! The initial sea surface temperature [degC] character(len=40) :: mdl = "benchmark_init_temperature_salinity" ! This subroutine's name. integer :: i, j, k, k1, is, ie, js, je, nz, itt diff --git a/src/user/dumbbell_surface_forcing.F90 b/src/user/dumbbell_surface_forcing.F90 index 693d2b5ceb..a1d8bf4b52 100644 --- a/src/user/dumbbell_surface_forcing.F90 +++ b/src/user/dumbbell_surface_forcing.F90 @@ -52,7 +52,7 @@ subroutine dumbbell_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) !! have NULL ptrs. type(time_type), intent(in) :: day !< Time of the fluxes. real, intent(in) :: dt !< The amount of time over which - !! the fluxes apply [s] + !! the fluxes apply [T ~> s] type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(dumbbell_surface_forcing_CS), pointer :: CS !< A control structure returned by a previous @@ -126,7 +126,7 @@ subroutine dumbbell_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS) end subroutine dumbbell_buoyancy_forcing !> Dynamic forcing for the dumbbell test case -subroutine dumbbell_dynamic_forcing(sfc_state, fluxes, day, dt, G, CS) +subroutine dumbbell_dynamic_forcing(sfc_state, fluxes, day, dt, G, US, CS) type(surface), intent(inout) :: sfc_state !< A structure containing fields that !! describe the surface state of the ocean. type(forcing), intent(inout) :: fluxes !< A structure containing pointers to any @@ -134,15 +134,17 @@ subroutine dumbbell_dynamic_forcing(sfc_state, fluxes, day, dt, G, CS) !! have NULL ptrs. type(time_type), intent(in) :: day !< Time of the fluxes. real, intent(in) :: dt !< The amount of time over which - !! the fluxes apply [s] + !! the fluxes apply [T ~> s] type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(dumbbell_surface_forcing_CS), pointer :: CS !< A control structure returned by a previous !! call to dumbbell_surface_forcing_init ! Local variables integer :: i, j, is, ie, js, je integer :: isd, ied, jsd, jed integer :: idays, isecs - real :: deg_rad, rdays + real :: deg_rad ! A conversion factor from degrees to radians [nondim] + real :: rdays ! The elapsed time [days] is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec @@ -178,11 +180,12 @@ subroutine dumbbell_surface_forcing_init(Time, G, US, param_file, diag, CS) type(dumbbell_surface_forcing_CS), & pointer :: CS !< A pointer to the control structure for this module ! Local variables - real :: S_surf, S_range - real :: x, y + real :: S_surf ! Initial surface salinity [ppt] + real :: S_range ! Range of the initial vertical distribution of salinity [ppt] + real :: x, y ! Latitude and longitude normalized by the domain size [nondim] integer :: i, j logical :: dbrotate ! If true, rotate the domain. -#include "version_variable.h" +# include "version_variable.h" character(len=40) :: mdl = "dumbbell_surface_forcing" ! This module's name. if (associated(CS)) then diff --git a/src/user/user_revise_forcing.F90 b/src/user/user_revise_forcing.F90 index bf31ca02f8..eb9694a091 100644 --- a/src/user/user_revise_forcing.F90 +++ b/src/user/user_revise_forcing.F90 @@ -24,9 +24,6 @@ module user_revise_forcing real :: cdrag !< The quadratic bottom drag coefficient. end type user_revise_forcing_CS -! This include declares and sets the variable "version". -#include "version_variable.h" - character(len=40) :: mdl = "user_revise_forcing" !< This module's name. contains !> This subroutine sets the surface wind stresses. @@ -41,6 +38,7 @@ subroutine user_alter_forcing(sfc_state, fluxes, day, G, CS) type(user_revise_forcing_CS), pointer :: CS !< A pointer to the control structure !! returned by a previous call to !! surface_forcing_init. + return end subroutine user_alter_forcing @@ -52,6 +50,10 @@ subroutine user_revise_forcing_init(param_file,CS) !! returned by a previous call to !! surface_forcing_init. + ! This include declares and sets the variable "version". +# include "version_variable.h" + character(len=40) :: mdl = "user_revise_forcing" !< This module's name. + call log_version(param_file, mdl, version) end subroutine user_revise_forcing_init