From 52b018bb9adf9c026e8aa3b3bad2e09c75b07472 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Thu, 15 Aug 2019 14:46:33 -0600 Subject: [PATCH 001/225] Add time0 and timenow in MCT --- config_src/mct_driver/ocn_comp_mct.F90 | 14 ++++++++++---- 1 file changed, 10 insertions(+), 4 deletions(-) diff --git a/config_src/mct_driver/ocn_comp_mct.F90 b/config_src/mct_driver/ocn_comp_mct.F90 index 5b581e0427..d7d66b32f6 100644 --- a/config_src/mct_driver/ocn_comp_mct.F90 +++ b/config_src/mct_driver/ocn_comp_mct.F90 @@ -111,6 +111,7 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename ) ! local variables type(time_type) :: time0 !< Model start time + type(time_type) :: timenow !< Model current time type(ESMF_time) :: time_var !< ESMF_time variable to query time type(ESMF_time) :: time_in_ESMF !< Initial time for ocean type(ESMF_timeInterval) :: ocn_cpl_interval !< Ocean coupling interval @@ -202,11 +203,16 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename ) call set_calendar_type(NOLEAP) !TODO: confirm this - ! Get the initial time - call ESMF_ClockGet(EClock, currTime=time_var, rc=rc) + ! Get start time + call ESMF_ClockGet(EClock, StartTime=time_var, rc=rc) call ESMF_TimeGet(time_var, yy=year, mm=month, dd=day, h=hour, m=minute, s=seconds, rc=rc) time0 = set_date(year, month, day, hour, minute, seconds, err_msg=err_msg) + ! Get current time + call ESMF_ClockGet(EClock, currTime=time_var, rc=rc) + call ESMF_TimeGet(time_var, yy=year, mm=month, dd=day, h=hour, m=minute, s=seconds, rc=rc) + timenow = set_date(year, month, day, hour, minute, seconds, err_msg=err_msg) + ! Debugging clocks if (debug .and. is_root_pe()) then write(glb%stdout,*) 'ocn_init_mct, current time: y,m,d-',year,month,day,'h,m,s=',hour,minute,seconds @@ -279,7 +285,7 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename ) runtype = get_runtype() if (runtype == "initial") then ! startup (new run) - 'n' is needed below since we don't specify input_filename in input.nml - call ocean_model_init(glb%ocn_public, glb%ocn_state, time0, time0, input_restart_file = 'n') + call ocean_model_init(glb%ocn_public, glb%ocn_state, time0, timenow, input_restart_file = 'n') else ! hybrid or branch or continuos runs ! get output path root call seq_infodata_GetData( glb%infodata, outPathRoot=restartpath ) @@ -295,7 +301,7 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename ) write(glb%stdout,*) 'Reading restart file: ',trim(restartfile) end if call shr_file_freeUnit(nu) - call ocean_model_init(glb%ocn_public, glb%ocn_state, time0, time0, input_restart_file=trim(restartfile)) + call ocean_model_init(glb%ocn_public, glb%ocn_state, time0, timenow, input_restart_file=trim(restartfile)) endif if (is_root_pe()) then write(glb%stdout,'(/12x,a/)') '======== COMPLETED MOM INITIALIZATION ========' From bd5760550192f146d4cd9e130fb91441c7d85d30 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Thu, 15 Aug 2019 15:18:02 -0600 Subject: [PATCH 002/225] Add time0 and timenow in NUOPC --- config_src/nuopc_driver/mom_cap.F90 | 27 ++++++++++++++++++++++----- 1 file changed, 22 insertions(+), 5 deletions(-) diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index d6a7837c83..b9ac3c0814 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -737,8 +737,9 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) type(ice_ocean_boundary_type), pointer :: Ice_ocean_boundary => NULL() type(ocean_internalstate_wrapper) :: ocean_internalstate type(ocean_grid_type), pointer :: ocean_grid => NULL() - type(time_type) :: Run_len ! length of experiment - type(time_type) :: Time + type(time_type) :: Run_len !< length of experiment + type(time_type) :: time0 !< Model start time + type(time_type) :: timenow !< Model current time type(time_type) :: Time_restart type(time_type) :: DT integer :: DT_OCEAN @@ -841,7 +842,23 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) ! this ocean connector will be driven at set interval DT = set_time (DT_OCEAN, 0) - Time = set_date (YEAR,MONTH,DAY,HOUR,MINUTE,SECOND) + ! get current time + timenow = set_date (YEAR,MONTH,DAY,HOUR,MINUTE,SECOND) + + ! get start/reference time + call ESMF_ClockGet(CLOCK, refTime=MyTime, RC=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & + line=__LINE__, & + file=__FILE__)) & + return ! bail out + + call ESMF_TimeGet (MyTime, YY=YEAR, MM=MONTH, DD=DAY, H=HOUR, M=MINUTE, S=SECOND, RC=rc ) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & + line=__LINE__, & + file=__FILE__)) & + return ! bail out + + time0 = set_date (YEAR,MONTH,DAY,HOUR,MINUTE,SECOND) ! rsd need to figure out how to get this without share code !call shr_nuopc_get_component_instance(gcomp, inst_suffix, inst_index) @@ -971,9 +988,9 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) ocean_public%is_ocean_pe = .true. if (len_trim(restartfile) > 0) then - call ocean_model_init(ocean_public, ocean_state, Time, Time, input_restart_file=trim(restartfile)) + call ocean_model_init(ocean_public, ocean_state, time0, timenow, input_restart_file=trim(restartfile)) else - call ocean_model_init(ocean_public, ocean_state, Time, Time) + call ocean_model_init(ocean_public, ocean_state, time0, timenow) endif call ocean_model_init_sfc(ocean_state, ocean_public) From 4191689d3d213190bef81a1bd8f94a376efc87fc Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Fri, 16 Aug 2019 13:55:06 -0600 Subject: [PATCH 003/225] Deletes unecessary if statement --- config_src/nuopc_driver/mom_cap.F90 | 15 ++++++++++----- 1 file changed, 10 insertions(+), 5 deletions(-) diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index b9ac3c0814..e4afba6524 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -845,6 +845,10 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) ! get current time timenow = set_date (YEAR,MONTH,DAY,HOUR,MINUTE,SECOND) + if (is_root_pe()) then + write(logunit,*) subname//'current time: y,m,d-',year,month,day,'h,m,s=',hour,minute,second + endif + ! get start/reference time call ESMF_ClockGet(CLOCK, refTime=MyTime, RC=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & @@ -860,6 +864,10 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) time0 = set_date (YEAR,MONTH,DAY,HOUR,MINUTE,SECOND) + if (is_root_pe()) then + write(logunit,*) subname//'start time: y,m,d-',year,month,day,'h,m,s=',hour,minute,second + endif + ! rsd need to figure out how to get this without share code !call shr_nuopc_get_component_instance(gcomp, inst_suffix, inst_index) !inst_name = "OCN"//trim(inst_suffix) @@ -987,11 +995,8 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) endif ocean_public%is_ocean_pe = .true. - if (len_trim(restartfile) > 0) then - call ocean_model_init(ocean_public, ocean_state, time0, timenow, input_restart_file=trim(restartfile)) - else - call ocean_model_init(ocean_public, ocean_state, time0, timenow) - endif + + call ocean_model_init(ocean_public, ocean_state, time0, timenow, input_restart_file=trim(restartfile)) call ocean_model_init_sfc(ocean_state, ocean_public) From 52225e34a6d87bd26c7a9f57f79ab3c08d0f07cd Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Fri, 16 Aug 2019 15:35:36 -0600 Subject: [PATCH 004/225] Deletes blank line contains spaces --- config_src/nuopc_driver/mom_cap.F90 | 1 - 1 file changed, 1 deletion(-) diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index e4afba6524..d1681bed25 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -995,7 +995,6 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) endif ocean_public%is_ocean_pe = .true. - call ocean_model_init(ocean_public, ocean_state, time0, timenow, input_restart_file=trim(restartfile)) call ocean_model_init_sfc(ocean_state, ocean_public) From b61053bdf594357b5f6135ff9fde64fd481aa6b4 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Fri, 16 Aug 2019 16:04:23 -0600 Subject: [PATCH 005/225] Replaces timenow with time_start and corrects comments. --- config_src/mct_driver/ocn_comp_mct.F90 | 12 ++++++------ config_src/nuopc_driver/mom_cap.F90 | 8 ++++---- 2 files changed, 10 insertions(+), 10 deletions(-) diff --git a/config_src/mct_driver/ocn_comp_mct.F90 b/config_src/mct_driver/ocn_comp_mct.F90 index d7d66b32f6..a7495d7a7f 100644 --- a/config_src/mct_driver/ocn_comp_mct.F90 +++ b/config_src/mct_driver/ocn_comp_mct.F90 @@ -109,9 +109,9 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename ) type(mct_aVect) , intent(inout) :: o2x_o !< Fluxes from ocean to coupler, computed by ocean character(len=*), optional , intent(in) :: NLFilename !< Namelist filename - ! local variables - type(time_type) :: time0 !< Model start time - type(time_type) :: timenow !< Model current time + ! local variable + type(time_type) :: time0 !< Start time of coupled model's calendar. + type(time_type) :: time_start !< The time at which to initialize the ocean model type(ESMF_time) :: time_var !< ESMF_time variable to query time type(ESMF_time) :: time_in_ESMF !< Initial time for ocean type(ESMF_timeInterval) :: ocn_cpl_interval !< Ocean coupling interval @@ -211,7 +211,7 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename ) ! Get current time call ESMF_ClockGet(EClock, currTime=time_var, rc=rc) call ESMF_TimeGet(time_var, yy=year, mm=month, dd=day, h=hour, m=minute, s=seconds, rc=rc) - timenow = set_date(year, month, day, hour, minute, seconds, err_msg=err_msg) + time_start = set_date(year, month, day, hour, minute, seconds, err_msg=err_msg) ! Debugging clocks if (debug .and. is_root_pe()) then @@ -285,7 +285,7 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename ) runtype = get_runtype() if (runtype == "initial") then ! startup (new run) - 'n' is needed below since we don't specify input_filename in input.nml - call ocean_model_init(glb%ocn_public, glb%ocn_state, time0, timenow, input_restart_file = 'n') + call ocean_model_init(glb%ocn_public, glb%ocn_state, time0, time_start, input_restart_file = 'n') else ! hybrid or branch or continuos runs ! get output path root call seq_infodata_GetData( glb%infodata, outPathRoot=restartpath ) @@ -301,7 +301,7 @@ subroutine ocn_init_mct( EClock, cdata_o, x2o_o, o2x_o, NLFilename ) write(glb%stdout,*) 'Reading restart file: ',trim(restartfile) end if call shr_file_freeUnit(nu) - call ocean_model_init(glb%ocn_public, glb%ocn_state, time0, timenow, input_restart_file=trim(restartfile)) + call ocean_model_init(glb%ocn_public, glb%ocn_state, time0, time_start, input_restart_file=trim(restartfile)) endif if (is_root_pe()) then write(glb%stdout,'(/12x,a/)') '======== COMPLETED MOM INITIALIZATION ========' diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index d1681bed25..aa034e7d36 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -738,8 +738,8 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) type(ocean_internalstate_wrapper) :: ocean_internalstate type(ocean_grid_type), pointer :: ocean_grid => NULL() type(time_type) :: Run_len !< length of experiment - type(time_type) :: time0 !< Model start time - type(time_type) :: timenow !< Model current time + type(time_type) :: time0 !< Start time of coupled model's calendar. + type(time_type) :: time_start !< The time at which to initialize the ocean model type(time_type) :: Time_restart type(time_type) :: DT integer :: DT_OCEAN @@ -843,7 +843,7 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) ! this ocean connector will be driven at set interval DT = set_time (DT_OCEAN, 0) ! get current time - timenow = set_date (YEAR,MONTH,DAY,HOUR,MINUTE,SECOND) + time_start = set_date (YEAR,MONTH,DAY,HOUR,MINUTE,SECOND) if (is_root_pe()) then write(logunit,*) subname//'current time: y,m,d-',year,month,day,'h,m,s=',hour,minute,second @@ -995,7 +995,7 @@ subroutine InitializeAdvertise(gcomp, importState, exportState, clock, rc) endif ocean_public%is_ocean_pe = .true. - call ocean_model_init(ocean_public, ocean_state, time0, timenow, input_restart_file=trim(restartfile)) + call ocean_model_init(ocean_public, ocean_state, time0, time_start, input_restart_file=trim(restartfile)) call ocean_model_init_sfc(ocean_state, ocean_public) From 6e145dfaa49f24c62e7fc1ff6e65aaf48cc0d227 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Tue, 20 Aug 2019 15:45:41 -0600 Subject: [PATCH 006/225] Obsolete parameter ADD_KV_SLOW --- src/core/MOM_variables.F90 | 2 -- src/diagnostics/MOM_obsolete_params.F90 | 3 +++ 2 files changed, 3 insertions(+), 2 deletions(-) diff --git a/src/core/MOM_variables.F90 b/src/core/MOM_variables.F90 index 3748684fd4..a0b31c37de 100644 --- a/src/core/MOM_variables.F90 +++ b/src/core/MOM_variables.F90 @@ -257,8 +257,6 @@ module MOM_variables real, pointer, dimension(:,:,:) :: TKE_turb => NULL() !< The turbulent kinetic energy per unit mass at the interfaces [m2 s-2]. !! This may be at the tracer or corner points - logical :: add_Kv_slow !< If True, add Kv_slow when calculating the 'coupling coefficient' (a_cpl) - !! at the interfaces in find_coupling_coef. end type vertvisc_type !> Container for information about the summed layer transports diff --git a/src/diagnostics/MOM_obsolete_params.F90 b/src/diagnostics/MOM_obsolete_params.F90 index 0e563648f5..4062d04fd9 100644 --- a/src/diagnostics/MOM_obsolete_params.F90 +++ b/src/diagnostics/MOM_obsolete_params.F90 @@ -47,6 +47,9 @@ subroutine find_obsolete_params(param_file) call obsolete_logical(param_file, "BT_CONT_BT_THICK", & hint="Instead use BT_THICK_SCHEME='FROM_BT_CONT'.") + call obsolete_logical(param_file, "ADD_KV_SLOW", & + hint="This option is no longer needed, nor supported.") + call obsolete_logical(param_file, "APPLY_OBC_U", & hint="Instead use OBC_NUMBER_SEGMENTS>0 and use the new segments protocol.") call obsolete_logical(param_file, "APPLY_OBC_V", & From 9f2f709151cc7b0a3023e493198c62fb55658957 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Tue, 20 Aug 2019 15:47:25 -0600 Subject: [PATCH 007/225] Removes visc%Kv_slow from restart files --- .../vertical/MOM_set_viscosity.F90 | 27 +++---------------- 1 file changed, 4 insertions(+), 23 deletions(-) diff --git a/src/parameterizations/vertical/MOM_set_viscosity.F90 b/src/parameterizations/vertical/MOM_set_viscosity.F90 index 1265067ef2..2180a62cc4 100644 --- a/src/parameterizations/vertical/MOM_set_viscosity.F90 +++ b/src/parameterizations/vertical/MOM_set_viscosity.F90 @@ -1735,11 +1735,10 @@ subroutine set_visc_register_restarts(HI, GV, param_file, visc, restart_CS) "Turbulent kinetic energy per unit mass at interfaces", "m2 s-2", z_grid='i') endif - ! MOM_bkgnd_mixing is always used, so always allocate visc%Kv_slow. GMM - call safe_alloc_ptr(visc%Kv_slow, isd, ied, jsd, jed, nz+1) - call register_restart_field(visc%Kv_slow, "Kv_slow", .false., restart_CS, & - "Vertical turbulent viscosity at interfaces due to slow processes", & - "m2 s-1", z_grid='i') + if (useKPP) then + ! MOM_bkgnd_mixing uses Kv_slow when KPP is defined. + call safe_alloc_ptr(visc%Kv_slow, isd, ied, jsd, jed, nz+1) + endif ! visc%MLD is used to communicate the state of the (e)PBL or KPP to the rest of the model call get_param(param_file, mdl, "MLE_USE_PBL_MLD", MLE_use_PBL_MLD, & @@ -1939,24 +1938,11 @@ subroutine set_visc_init(Time, G, GV, US, param_file, diag, visc, CS, restart_CS "The molecular value, ~1e-6 m2 s-1, may be used.", & units="m2 s-1", fail_if_missing=.true.) - call get_param(param_file, mdl, "ADD_KV_SLOW", visc%add_Kv_slow, & - "If true, the background vertical viscosity in the interior "//& - "(i.e., tidal + background + shear + convection) is added "//& - "when computing the coupling coefficient. The purpose of this "//& - "flag is to be able to recover previous answers and it will likely "//& - "be removed in the future since this option should always be true.", & - default=.false.) - call get_param(param_file, mdl, "USE_KPP", use_KPP, & "If true, turns on the [CVMix] KPP scheme of Large et al., 1994, "//& "to calculate diffusivities and non-local transport in the OBL.", & do_not_log=.true., default=.false.) - if (use_KPP .and. visc%add_Kv_slow) call MOM_error(FATAL,"set_visc_init: "//& - "When USE_KPP=True, ADD_KV_SLOW must be false. Otherwise vertical "//& - "viscosity due to slow processes will be double counted. Please set "//& - "ADD_KV_SLOW=False.") - call get_param(param_file, mdl, "KV_BBL_MIN", CS%KV_BBL_min, & "The minimum viscosities in the bottom boundary layer.", & units="m2 s-1", default=Kv_background, scale=US%m_to_Z**2) @@ -2049,11 +2035,6 @@ subroutine set_visc_init(Time, G, GV, US, param_file, diag, visc, CS, restart_CS enddo ; enddo ; enddo endif ; endif - if (associated(visc%Kv_slow)) then ; if (query_initialized(visc%Kv_slow, "Kv_slow", restart_CS)) then - do k=1,nz+1 ; do j=js,je ; do i=is,ie - visc%Kv_slow(i,j,k) = Z_rescale**2 * visc%Kv_slow(i,j,k) - enddo ; enddo ; enddo - endif ; endif endif end subroutine set_visc_init From 93727f8afe45340681c2b3c98ccfae77a8871ded Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Thu, 22 Aug 2019 16:27:47 -0600 Subject: [PATCH 008/225] Fix first four loop indices suggested by Bob NOTE: we need to add a warning that MEKE_EQUILIBRIUM only works if MEKE_KHCOEFF > 0. --- src/parameterizations/lateral/MOM_hor_visc.F90 | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index 5a31405268..e6a8a26631 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -407,7 +407,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, call barotropic_get_tav(BT, ubtav, vbtav, G) call pass_vector(ubtav, vbtav, G%Domain) - do j=js,je ; do i=is,ie + do j=js-1,je+1 ; do i=is-1,ie+1 dudx_bt(i,j) = CS%DY_dxT(i,j)*(G%IdyCu(I,j) * ubtav(I,j) - & G%IdyCu(I-1,j) * ubtav(I-1,j)) dvdy_bt(i,j) = CS%DX_dyT(i,j)*(G%IdxCv(i,J) * vbtav(i,J) - & @@ -417,7 +417,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, call pass_var(dudx_bt, G%Domain, complete=.true.) call pass_var(dvdy_bt, G%Domain, complete=.true.) - do j=Jsq-1,Jeq+2 ; do i=Isq-1,Ieq+2 + do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 sh_xx_bt(i,j) = dudx_bt(i,j) - dvdy_bt(i,j) enddo ; enddo @@ -433,11 +433,11 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, call pass_var(dudy_bt, G%Domain, position=CORNER, complete=.true.) if (CS%no_slip) then - do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 + do J=js-1,Jeq ; do I=is-1,Ieq sh_xy_bt(I,J) = (2.0-G%mask2dBu(I,J)) * ( dvdx_bt(I,J) + dudy_bt(I,J) ) enddo ; enddo else - do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 + do J=js-1,Jeq ; do I=is-1,Ieq sh_xy_bt(I,J) = G%mask2dBu(I,J) * ( dvdx_bt(I,J) + dudy_bt(I,J) ) enddo ; enddo endif From e7781917d054d40973a2848e7ba7215d2334a71e Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Fri, 23 Aug 2019 11:18:39 -0600 Subject: [PATCH 009/225] Deletes max_diss_rate_bt since it is never used --- src/parameterizations/lateral/MOM_hor_visc.F90 | 9 --------- 1 file changed, 9 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index e6a8a26631..a0ff1d96a0 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -251,7 +251,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, grad_vel_mag_h, & ! Magnitude of the velocity gradient tensor squared at h-points [s-2] grad_vel_mag_bt_h, & ! Magnitude of the barotropic velocity gradient tensor squared at h-points [s-2] grad_d2vel_mag_h, & ! Magnitude of the Laplacian of the velocity vector, squared [m-2 s-2] - max_diss_rate_bt, & ! maximum possible energy dissipated by barotropic lateral friction [m2 s-3] boundary_mask ! A mask that zeroes out cells with at least one land edge real, dimension(SZIB_(G),SZJB_(G)) :: & @@ -401,8 +400,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, str_xx_GME(:,:) = 0.0 str_xy_GME(:,:) = 0.0 -! call pass_var(boundary_mask, G%Domain, complete=.true.) - ! Get barotropic velocities and their gradients call barotropic_get_tav(BT, ubtav, vbtav, G) call pass_vector(ubtav, vbtav, G%Domain) @@ -451,12 +448,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, (0.25*(dudy_bt(I,J)+dudy_bt(I-1,J)+dudy_bt(I,J-1)+dudy_bt(I-1,J-1)) )**2) enddo ; enddo - if (associated(MEKE)) then ; if (associated(MEKE%mom_src)) then - do j=Jsq-1,Jeq+2 ; do i=Isq-1,Ieq+2 - max_diss_rate_bt(i,j) = 2.0*MEKE%MEKE(i,j) * grad_vel_mag_bt_h(i,j) - enddo ; enddo - endif ; endif - do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 grad_vel_mag_bt_q(I,J) = boundary_mask(i,j) * (dvdx_bt(i,j)**2 + dudy_bt(i,j)**2 + & (0.25*(dudx_bt(i,j)+dudx_bt(i+1,j)+dudx_bt(i,j+1)+dudx_bt(i+1,j+1)))**2 + & From 3b7cb69a4352f66cadc42599e20821f679658fc2 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Fri, 23 Aug 2019 11:27:04 -0600 Subject: [PATCH 010/225] Fixes the remaining loop indices related to GME --- src/parameterizations/lateral/MOM_hor_visc.F90 | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index a0ff1d96a0..8ad3e7f212 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -713,21 +713,21 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, call pass_var(div_xx, G%Domain, complete=.true.) ! Divergence gradient - do j=Jsq-1,Jeq+2 ; do I=is-2,Ieq+1 + do j=Jsq,Jeq+1 ; do I=Isq-1,Ieq+1 div_xx_dx(I,j) = G%IdxCu(I,j)*(div_xx(i+1,j) - div_xx(i,j)) enddo ; enddo - do J=js-2,Jeq+1 ; do i=Isq-1,Ieq+2 + do J=Jsq-1,Jeq+1 ; do i=Isq,Ieq+1 div_xx_dy(i,J) = G%IdyCv(i,J)*(div_xx(i,j+1) - div_xx(i,j)) enddo ; enddo call pass_vector(div_xx_dx, div_xx_dy, G%Domain) ! Magnitude of divergence gradient - do j=Jsq-1,Jeq+2 ; do i=Isq-1,Ieq+2 + do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 grad_div_mag_h(i,j) =sqrt((0.5*(div_xx_dx(I,j) + div_xx_dx(I-1,j)))**2 + & (0.5 * (div_xx_dy(i,J) + div_xx_dy(i,J-1)))**2) enddo ; enddo - do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 + do J=js-1,Jeq ; do I=is-1,Ieq grad_div_mag_q(I,J) =sqrt((0.5*(div_xx_dx(I,j) + div_xx_dx(I,j+1)))**2 + & (0.5 * (div_xx_dy(i,J) + div_xx_dy(i+1,J)))**2) enddo ; enddo @@ -740,10 +740,10 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, do J=js-2,Jeq+1 ; do i=Isq-1,Ieq+2 div_xx_dy(i,J) = 0.0 enddo ; enddo - do j=Jsq-1,Jeq+2 ; do i=Isq-1,Ieq+2 + do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 grad_div_mag_h(i,j) = 0.0 enddo ; enddo - do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 + do J=js-1,Jeq ; do I=is-1,Ieq grad_div_mag_q(I,J) = 0.0 enddo ; enddo From 32578013e14e5f96fc208672cf7bdb9e014c7a05 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Fri, 23 Aug 2019 13:49:56 -0600 Subject: [PATCH 011/225] Removes unneeded calls to pass_var and add calls to pass_vector --- src/parameterizations/lateral/MOM_hor_visc.F90 | 12 +++--------- 1 file changed, 3 insertions(+), 9 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index 8ad3e7f212..7b1d7e00e2 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -5,7 +5,7 @@ module MOM_hor_visc use MOM_diag_mediator, only : post_data, register_diag_field, safe_alloc_ptr use MOM_diag_mediator, only : diag_ctrl, time_type -use MOM_domains, only : pass_var, CORNER, pass_vector +use MOM_domains, only : pass_var, CORNER, pass_vector, AGRID, BGRID_NE use MOM_error_handler, only : MOM_error, FATAL, WARNING, is_root_pe use MOM_file_parser, only : get_param, log_version, param_file_type use MOM_grid, only : ocean_grid_type @@ -411,8 +411,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, G%IdxCv(i,J-1) * vbtav(i,J-1)) enddo; enddo - call pass_var(dudx_bt, G%Domain, complete=.true.) - call pass_var(dvdy_bt, G%Domain, complete=.true.) + call pass_vector(dudx_bt, dvdy_bt, G%Domain, stagger=BGRID_NE) do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 sh_xx_bt(i,j) = dudx_bt(i,j) - dvdy_bt(i,j) @@ -426,8 +425,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, - ubtav(I,j)*G%IdxCu(I,j)) enddo ; enddo - call pass_var(dvdx_bt, G%Domain, position=CORNER, complete=.true.) - call pass_var(dudy_bt, G%Domain, position=CORNER, complete=.true.) + call pass_vector(dvdx_bt, dudy_bt, G%Domain, stagger=AGRID) if (CS%no_slip) then do J=js-1,Jeq ; do I=is-1,Ieq @@ -710,8 +708,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, (h(i,j,k) + GV%H_subroundoff) enddo ; enddo - call pass_var(div_xx, G%Domain, complete=.true.) - ! Divergence gradient do j=Jsq,Jeq+1 ; do I=Isq-1,Ieq+1 div_xx_dx(I,j) = G%IdxCu(I,j)*(div_xx(i+1,j) - div_xx(i,j)) @@ -720,8 +716,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, div_xx_dy(i,J) = G%IdyCv(i,J)*(div_xx(i,j+1) - div_xx(i,j)) enddo ; enddo - call pass_vector(div_xx_dx, div_xx_dy, G%Domain) - ! Magnitude of divergence gradient do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 grad_div_mag_h(i,j) =sqrt((0.5*(div_xx_dx(I,j) + div_xx_dx(I-1,j)))**2 + & From b339108f879b842a71ca19f6fe32ceda1af66210 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Fri, 23 Aug 2019 13:56:35 -0600 Subject: [PATCH 012/225] Groups terms for rotational invariance and deletes beta_h and beta_q (unused) --- src/parameterizations/lateral/MOM_hor_visc.F90 | 17 ++++------------- 1 file changed, 4 insertions(+), 13 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index 7b1d7e00e2..fe1c069257 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -243,7 +243,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, FrictWorkIntz, & ! depth integrated energy dissipated by lateral friction [W m-2] Leith_Kh_h, & ! Leith Laplacian viscosity at h-points [m2 s-1] Leith_Ah_h, & ! Leith bi-harmonic viscosity at h-points [m4 s-1] - beta_h, & ! Gradient of planetary vorticity at h-points [m-1 s-1] grad_vort_mag_h, & ! Magnitude of vorticity gradient at h-points [m-1 s-1] grad_vort_mag_h_2d, & ! Magnitude of 2d vorticity gradient at h-points [m-1 s-1] grad_div_mag_h, & ! Magnitude of divergence gradient at h-points [m-1 s-1] @@ -264,7 +263,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, vort_xy, & ! Vertical vorticity (dv/dx - du/dy) including metric terms [s-1] Leith_Kh_q, & ! Leith Laplacian viscosity at q-points [m2 s-1] Leith_Ah_q, & ! Leith bi-harmonic viscosity at q-points [m4 s-1] - beta_q, & ! Gradient of planetary vorticity at q-points [m-1 s-1] grad_vort_mag_q, & ! Magnitude of vorticity gradient at q-points [m-1 s-1] grad_vort_mag_q_2d, & ! Magnitude of 2d vorticity gradient at q-points [m-1 s-1] grad_div_mag_q, & ! Magnitude of divergence gradient at q-points [m-1 s-1] @@ -442,14 +440,14 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, do j=Jsq-1,Jeq+2 ; do i=Isq-1,Ieq+2 grad_vel_mag_bt_h(i,j) = boundary_mask(i,j) * (dudx_bt(i,j)**2 + dvdy_bt(i,j)**2 + & - (0.25*(dvdx_bt(I,J)+dvdx_bt(I-1,J)+dvdx_bt(I,J-1)+dvdx_bt(I-1,J-1)) )**2 + & - (0.25*(dudy_bt(I,J)+dudy_bt(I-1,J)+dudy_bt(I,J-1)+dudy_bt(I-1,J-1)) )**2) + (0.25*((dvdx_bt(I,J)+dvdx_bt(I-1,J-1))+(dvdx_bt(I,J-1)+dvdx_bt(I-1,J)) )**2 + & + (0.25*((dudy_bt(I,J)+dudy_bt(I-1,J-1))+(dudy_bt(I,J-1)+dudy_bt(I-1,J)) )**2) enddo ; enddo do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 grad_vel_mag_bt_q(I,J) = boundary_mask(i,j) * (dvdx_bt(i,j)**2 + dudy_bt(i,j)**2 + & - (0.25*(dudx_bt(i,j)+dudx_bt(i+1,j)+dudx_bt(i,j+1)+dudx_bt(i+1,j+1)))**2 + & - (0.25*(dvdy_bt(i,j)+dvdy_bt(i+1,j)+dvdy_bt(i,j+1)+dvdy_bt(i+1,j+1)) )**2) + (0.25*(dudx_bt(i,j)+dudx_bt(i+1,j+1))+(dudx_bt(i,j+1)+dudx_bt(i+1,j)))**2 + & + (0.25*(dvdy_bt(i,j)+dvdy_bt(i+1,j+1))+(dvdy_bt(i,j+1)+dvdy_bt(i+1,j)))**2) enddo ; enddo endif ! use_GME @@ -745,13 +743,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, ! Add in beta for the Leith viscosity if (CS%use_beta_in_Leith) then - do j=Jsq-1,Jeq+2 ; do i=Isq-1,Ieq+2 - beta_h(i,j) = sqrt( G%dF_dx(i,j)**2 + G%dF_dy(i,j)**2 ) - enddo; enddo - do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 - beta_q(I,J) = sqrt( (0.25*(G%dF_dx(i,j)+G%dF_dx(i+1,j)+G%dF_dx(i,j+1)+G%dF_dx(i+1,j+1))**2) + & - (0.25*(G%dF_dy(i,j)+G%dF_dy(i+1,j)+G%dF_dy(i,j+1)+G%dF_dy(i+1,j+1))**2) ) - enddo ; enddo do J=js-2,Jeq+1 ; do i=is-1,Ieq+1 vort_xy_dx(i,J) = vort_xy_dx(i,J) + 0.5 * ( G%dF_dx(i,j) + G%dF_dx(i,j+1)) From 0aeec9de17270c86f19ac2c37bf90e48c7fe182f Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Mon, 26 Aug 2019 15:15:03 -0600 Subject: [PATCH 013/225] Adds correct calculation of GME at q points --- .../lateral/MOM_hor_visc.F90 | 216 +++++++++--------- 1 file changed, 109 insertions(+), 107 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index fe1c069257..44131952a1 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -250,7 +250,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, grad_vel_mag_h, & ! Magnitude of the velocity gradient tensor squared at h-points [s-2] grad_vel_mag_bt_h, & ! Magnitude of the barotropic velocity gradient tensor squared at h-points [s-2] grad_d2vel_mag_h, & ! Magnitude of the Laplacian of the velocity vector, squared [m-2 s-2] - boundary_mask ! A mask that zeroes out cells with at least one land edge + boundary_mask_h ! A mask that zeroes out cells with at least one land edge real, dimension(SZIB_(G),SZJB_(G)) :: & dvdx, dudy, & ! components in the shearing strain [s-1] @@ -269,13 +269,15 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, grad_vel_mag_q, & ! Magnitude of the velocity gradient tensor squared at q-points [s-2] hq, & ! harmonic mean of the harmonic means of the u- & v point thicknesses [H ~> m or kg m-2] ! This form guarantees that hq/hu < 4. - grad_vel_mag_bt_q ! Magnitude of the barotropic velocity gradient tensor squared at q-points [s-2] + grad_vel_mag_bt_q, & ! Magnitude of the barotropic velocity gradient tensor squared at q-points [s-2] + boundary_mask_q ! A mask that zeroes out cells with at least one land edge real, dimension(SZIB_(G),SZJB_(G),SZK_(G)) :: & Ah_q, & ! biharmonic viscosity at corner points [m4 s-1] Kh_q, & ! Laplacian viscosity at corner points [m2 s-1] vort_xy_q, & ! vertical vorticity at corner points [s-1] - GME_coeff_q !< GME coeff. at q-points [m2 s-1] + GME_coeff_q, & !< GME coeff. at q-points [m2 s-1] + max_diss_rate_q ! maximum possible energy dissipated by lateral friction [m2 s-3] real, dimension(SZIB_(G),SZJ_(G),SZK_(G)+1) :: & KH_u_GME !< interface height diffusivities in u-columns [m2 s-1] @@ -285,7 +287,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, Ah_h, & ! biharmonic viscosity at thickness points [m4 s-1] Kh_h, & ! Laplacian viscosity at thickness points [m2 s-1] diss_rate, & ! MKE dissipated by parameterized shear production [m2 s-3] - max_diss_rate, & ! maximum possible energy dissipated by lateral friction [m2 s-3] + max_diss_rate_h, & ! maximum possible energy dissipated by lateral friction [m2 s-3] target_diss_rate_GME, & ! the maximum theoretical dissipation plus the amount spuriously dissipated ! by friction [m2 s-3] FrictWork, & ! work done by MKE dissipation mechanisms [W m-2] @@ -293,7 +295,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, FrictWorkMax, & ! maximum possible work done by MKE dissipation mechanisms [W m-2] FrictWork_GME, & ! work done by GME [W m-2] div_xx_h ! horizontal divergence [s-1] - !real, dimension(SZI_(G),SZJ_(G),SZK_(G)+1) :: & real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: & KH_t_GME, & !< interface height diffusivities in t-columns [m2 s-1] GME_coeff_h !< GME coeff. at h-points [m2 s-1] @@ -312,8 +313,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, real :: hu, hv ! Thicknesses interpolated by arithmetic means to corner ! points; these are first interpolated to u or v velocity ! points where masks are applied [H ~> m or kg m-2]. -! real :: hq ! harmonic mean of the harmonic means of the u- & v- -! ! point thicknesses, in H; This form guarantees that hq/hu < 4. real :: h_neglect ! thickness so small it can be lost in roundoff and so neglected [H ~> m or kg m-2] real :: h_neglect3 ! h_neglect^3 [H3 ~> m3 or kg3 m-6] real :: hrat_min ! minimum thicknesses at the 4 neighboring @@ -382,10 +381,17 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, use_MEKE_Ku = associated(MEKE%Ku) use_MEKE_Au = associated(MEKE%Au) - do j=js,je ; do i=is,ie - boundary_mask(i,j) = (G%mask2dCu(I,j) * G%mask2dCv(i,J) * G%mask2dCu(I-1,j) * G%mask2dCv(i,J-1)) + do j = G%jsc, G%jec ; do i = G%isc, G%iec + boundary_mask_h(i,j) = (G%mask2dCu(I,j) * G%mask2dCv(i,J) * G%mask2dCu(I-1,j) * G%mask2dCv(i,J-1)) enddo ; enddo + do J = G%JscB, G%JecB ; do I = G%IscB, G%IecB + boundary_mask_q(I,J) = (G%mask2dCv(i,J) * G%mask2dCv(i+1,J) * G%mask2dCu(I,j) * G%mask2dCu(I,j-1)) + enddo; enddo + +! call pass_var(boundary_mask_h, G%Domain, complete=.true.) + + if (CS%use_GME) then ! GME tapers off above this depth H0 = 1000.0 @@ -411,6 +417,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, call pass_vector(dudx_bt, dvdy_bt, G%Domain, stagger=BGRID_NE) + do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 sh_xx_bt(i,j) = dudx_bt(i,j) - dvdy_bt(i,j) enddo ; enddo @@ -435,19 +442,16 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, enddo ; enddo endif - ! Get thickness diffusivity for use in GME -! call thickness_diffuse_get_KH(thickness_diffuse, KH_u_GME, KH_v_GME, G) - do j=Jsq-1,Jeq+2 ; do i=Isq-1,Ieq+2 - grad_vel_mag_bt_h(i,j) = boundary_mask(i,j) * (dudx_bt(i,j)**2 + dvdy_bt(i,j)**2 + & - (0.25*((dvdx_bt(I,J)+dvdx_bt(I-1,J-1))+(dvdx_bt(I,J-1)+dvdx_bt(I-1,J)) )**2 + & - (0.25*((dudy_bt(I,J)+dudy_bt(I-1,J-1))+(dudy_bt(I,J-1)+dudy_bt(I-1,J)) )**2) + grad_vel_mag_bt_h(i,j) = boundary_mask_h(i,j) * (dudx_bt(i,j)**2 + dvdy_bt(i,j)**2 + & + (0.25*((dvdx_bt(I,J)+dvdx_bt(I-1,J-1))+(dvdx_bt(I,J-1)+dvdx_bt(I-1,J))))**2 + & + (0.25*((dudy_bt(I,J)+dudy_bt(I-1,J-1))+(dudy_bt(I,J-1)+dudy_bt(I-1,J))))**2) enddo ; enddo do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 - grad_vel_mag_bt_q(I,J) = boundary_mask(i,j) * (dvdx_bt(i,j)**2 + dudy_bt(i,j)**2 + & - (0.25*(dudx_bt(i,j)+dudx_bt(i+1,j+1))+(dudx_bt(i,j+1)+dudx_bt(i+1,j)))**2 + & - (0.25*(dvdy_bt(i,j)+dvdy_bt(i+1,j+1))+(dvdy_bt(i,j+1)+dvdy_bt(i+1,j)))**2) + grad_vel_mag_bt_q(I,J) = boundary_mask_h(I,J) * (dvdx_bt(i,j)**2 + dudy_bt(i,j)**2 + & + (0.25*((dudx_bt(i,j)+dudx_bt(i+1,j+1))+(dudx_bt(i,j+1)+dudx_bt(i+1,j))))**2 + & + (0.25*((dvdy_bt(i,j)+dvdy_bt(i+1,j+1))+(dvdy_bt(i,j+1)+dvdy_bt(i+1,j))))**2) enddo ; enddo endif ! use_GME @@ -555,6 +559,36 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, enddo endif endif + +! if (CS%use_GME) then +! call pass_var(dudx, G%Domain, complete=.true.) +! call pass_var(dvdy, G%Domain, complete=.true.) +! call pass_var(dvdx, G%Domain, position=CORNER, complete=.true.) +! call pass_var(dudy, G%Domain, position=CORNER, complete=.true.) +!! call pass_var(MEKE%MEKE, G%Domain, complete=.true.) +! + +! do j=Jsq-1,Jeq+2 ; do i=Isq-1,Ieq+2 +! do J=Jsq,Jeq+1 ; do i=Isq,Ieq+1 +! do j = js, je ; do i = is, ie +! grad_vel_mag_h(i,j) = boundary_mask_h(i,j) * (dudx(i,j)**2 + dvdy(i,j)**2 + & +! (0.25*((dvdx(I,J)+dvdx(I-1,J-1))+(dvdx(I,J-1)+dvdx(I-1,J))))**2 + & +! (0.25*((dudy(I,J)+dudy(I-1,J-1))+(dudy(I,J-1)+dudy(I-1,J))))**2) +! +! max_diss_rate_h(i,j,k) = 2.0 !*MEKE%MEKE(i,j) * sqrt(grad_vel_mag_h(i,j)) +! enddo ; enddo +! +! do J = G%JscB, G%JecB ; do I = G%IscB, G%IecB +! grad_vel_mag_q(I,J) = boundary_mask_h(I,J) * (dudx(i,j)**2 + dvdy(i,j)**2 + & +! (0.25*((dvdx(I,J)+dvdx(I-1,J-1))+(dvdx(I,J-1)+dvdx(I-1,J))))**2 + & +! (0.25*((dudy(I,J)+dudy(I-1,J-1))+(dudy(I,J-1)+dudy(I-1,J))))**2) +! +! max_diss_rate_q(I,J,k) = 0.5*(MEKE%MEKE(i,j)+MEKE%MEKE(i+1,j)+ & +! MEKE%MEKE(i,j+1)+MEKE%MEKE(i+1,j+1)) * sqrt(grad_vel_mag_q(I,J)) +! enddo ; enddo +! endif + + if (OBC%segment(n)%direction == OBC_DIRECTION_N) then ! There are extra wide halos here to accommodate the cross-corner-point ! OBC projections, but they might not be necessary if the accelerations @@ -659,25 +693,25 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, ! Note this a simple re-calculation of shearing components using the same discretization. ! We will consider using a circulation based calculation of vorticity later. ! Also note this will need OBC boundary conditions re-applied... - do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 - DY_dxBu = G%dyBu(I,J) * G%IdxBu(I,J) - dvdx(I,J) = DY_dxBu * (v(i+1,J,k) * G%IdyCv(i+1,J) - v(i,J,k) * G%IdyCv(i,J)) - DX_dyBu = G%dxBu(I,J) * G%IdyBu(I,J) - dudy(I,J) = DX_dyBu * (u(I,j+1,k) * G%IdxCu(I,j+1) - u(I,j,k) * G%IdxCu(I,j)) - enddo ; enddo +! do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 +! DY_dxBu = G%dyBu(I,J) * G%IdxBu(I,J) +! dvdx(I,J) = DY_dxBu * (v(i+1,J,k) * G%IdyCv(i+1,J) - v(i,J,k) * G%IdyCv(i,J)) +! DX_dyBu = G%dxBu(I,J) * G%IdyBu(I,J) +! dudy(I,J) = DX_dyBu * (u(I,j+1,k) * G%IdxCu(I,j+1) - u(I,j,k) * G%IdxCu(I,j)) +! enddo ; enddo ! Vorticity if (CS%no_slip) then do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 vort_xy(I,J) = (2.0-G%mask2dBu(I,J)) * ( dvdx(I,J) - dudy(I,J) ) - dudy(I,J) = (2.0-G%mask2dBu(I,J)) * dudy(I,J) - dvdx(I,J) = (2.0-G%mask2dBu(I,J)) * dvdx(I,J) +! dudy(I,J) = (2.0-G%mask2dBu(I,J)) * dudy(I,J) +! dvdx(I,J) = (2.0-G%mask2dBu(I,J)) * dvdx(I,J) enddo ; enddo else do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 vort_xy(I,J) = G%mask2dBu(I,J) * ( dvdx(I,J) - dudy(I,J) ) - dudy(I,J) = G%mask2dBu(I,J) * dudy(I,J) - dvdx(I,J) = G%mask2dBu(I,J) * dvdx(I,J) +! dudy(I,J) = G%mask2dBu(I,J) * dudy(I,J) +! dvdx(I,J) = G%mask2dBu(I,J) * dvdx(I,J) enddo ; enddo endif @@ -822,7 +856,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, visc_bound_rem = 0.0 Kh = hrat_min*CS%Kh_Max_xx(i,j) else - !visc_bound_rem = 1.0 - abs(Kh) / (hrat_min*CS%Kh_Max_xx(i,j)) visc_bound_rem = 1.0 - Kh / (hrat_min*CS%Kh_Max_xx(i,j)) endif endif @@ -988,7 +1021,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, visc_bound_rem = 0.0 Kh = hrat_min*CS%Kh_Max_xy(I,J) elseif (CS%Kh_Max_xy(I,J)>0.) then - !visc_bound_rem = 1.0 - abs(Kh) / (hrat_min*CS%Kh_Max_xy(I,J)) visc_bound_rem = 1.0 - Kh / (hrat_min*CS%Kh_Max_xy(I,J)) endif endif @@ -1051,85 +1083,46 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, enddo ; enddo - if (find_FrictWork) then - if (CS%biharmonic) call pass_vector(u0, v0, G%Domain) - call pass_var(dudx, G%Domain, complete=.true.) - call pass_var(dvdy, G%Domain, complete=.true.) - call pass_var(dvdx, G%Domain, position=CORNER, complete=.true.) - call pass_var(dudy, G%Domain, position=CORNER, complete=.true.) + if (CS%use_GME) then +! call pass_var(dudx, G%Domain, complete=.true.) +! call pass_var(dvdy, G%Domain, complete=.true.) +! call pass_var(dvdx, G%Domain, position=CORNER, complete=.true.) +! call pass_var(dudy, G%Domain, position=CORNER, complete=.true.) - if (CS%Laplacian) then - do j=js,je ; do i=is,ie - grad_vel_mag_h(i,j) = boundary_mask(i,j) * (dudx(i,j)**2 + dvdy(i,j)**2 + & - (0.25*(dvdx(I,J)+dvdx(I-1,J)+dvdx(I,J-1)+dvdx(I-1,J-1)) )**2 + & - (0.25*(dudy(I,J)+dudy(I-1,J)+dudy(I,J-1)+dudy(I-1,J-1)) )**2) - enddo ; enddo - else - do j=js,je ; do i=is,ie - grad_vel_mag_h(i,j) = 0.0 - enddo ; enddo - endif - if (CS%biharmonic) then - do j=js,je ; do i=is,ie - grad_d2vel_mag_h(i,j) = boundary_mask(i,j) * ((0.5*(u0(I,j) + u0(I-1,j)))**2 + & - (0.5*(v0(i,J) + v0(i,J-1)))**2) - enddo ; enddo - else - do j=js,je ; do i=is,ie - grad_d2vel_mag_h(i,j) = 0.0 - enddo ; enddo - endif + do j = js, je ; do i = is, ie + grad_vel_mag_h(i,j) = boundary_mask_h(i,j) * (dudx(i,j)**2 + dvdy(i,j)**2 + & + (0.25*((dvdx(I,J)+dvdx(I-1,J-1))+(dvdx(I,J-1)+dvdx(I-1,J))))**2 + & + (0.25*((dudy(I,J)+dudy(I-1,J-1))+(dudy(I,J-1)+dudy(I-1,J))))**2) - do j=js,je ; do i=is,ie - ! Diagnose -Kh * |del u|^2 - Ah * |del^2 u|^2 - diss_rate(i,j,k) = -Kh_h(i,j,k) * grad_vel_mag_h(i,j) - & - Ah_h(i,j,k) * grad_d2vel_mag_h(i,j) - - if (associated(MEKE)) then ; if (associated(MEKE%mom_src)) then - ! This is the maximum possible amount of energy that can be converted - ! per unit time, according to theory (multiplied by h) - max_diss_rate(i,j,k) = 2.0*MEKE%MEKE(i,j) * sqrt(grad_vel_mag_h(i,j)) - FrictWork_diss(i,j,k) = diss_rate(i,j,k) * h(i,j,k) * GV%H_to_kg_m2 - FrictWorkMax(i,j,k) = -max_diss_rate(i,j,k) * h(i,j,k) * GV%H_to_kg_m2 - - ! Determine how much work GME needs to do to reach the "target" ratio between - ! the amount of work actually done and the maximum allowed by theory. Note that - ! we need to add the FrictWork done by the dissipation operators, since this work - ! is done only for numerical stability and is therefore spurious - if (CS%use_GME) then - target_diss_rate_GME(i,j,k) = FWfrac * max_diss_rate(i,j,k) - diss_rate(i,j,k) - endif + max_diss_rate_h(i,j,k) = 2.0 * MEKE%MEKE(i,j) * sqrt(grad_vel_mag_h(i,j)) + enddo ; enddo - else - FrictWork_diss(i,j,k) = diss_rate(i,j,k) * h(i,j,k) * GV%H_to_kg_m2 - endif ; endif + do J = G%JscB, G%JecB ; do I = G%IscB, G%IecB + grad_vel_mag_q(I,J) = boundary_mask_q(I,J) * (dudx(i,j)**2 + dvdy(i,j)**2 + & + (0.25*((dvdx(I,J)+dvdx(I-1,J-1))+(dvdx(I,J-1)+dvdx(I-1,J))))**2 + & + (0.25*((dudy(I,J)+dudy(I-1,J-1))+(dudy(I,J-1)+dudy(I-1,J))))**2) - enddo ; enddo - endif + max_diss_rate_q(I,J,k) = 0.5*(MEKE%MEKE(i,j)+MEKE%MEKE(i+1,j)+ & + MEKE%MEKE(i,j+1)+MEKE%MEKE(i+1,j+1)) * sqrt(grad_vel_mag_q(I,J)) + enddo ; enddo + endif if (CS%use_GME) then - if (.not. (associated(MEKE))) call MOM_error(FATAL, & - "MEKE must be enabled for GME to be used.") - do J=Jsq,Jeq+1 ; do i=Isq,Ieq+1 - - if ((max_diss_rate(i,j,k) > 0) .and. (grad_vel_mag_bt_h(i,j)>0) ) then - GME_coeff = (MIN(G%bathyT(i,j)/H0,1.0)**2) * FWfrac*max_diss_rate(i,j,k) / grad_vel_mag_bt_h(i,j) -! GME_coeff = (MIN(G%bathyT(i,j)/H0,1.0)**2) * FWfrac*target_diss_rate_GME(i,j,k) / grad_vel_mag_bt_h(i,j) + if ((grad_vel_mag_bt_h(i,j)>0) .and. (max_diss_rate_h(i,j,k)>0)) then + GME_coeff = (MIN(G%bathyT(i,j)/H0,1.0)**2) * FWfrac*max_diss_rate_h(i,j,k) / grad_vel_mag_bt_h(i,j) else - GME_coeff = 0.0 + GME_coeff = 1.0 endif ! apply mask - GME_coeff = GME_coeff * boundary_mask(i,j) - + GME_coeff = GME_coeff * boundary_mask_h(i,j) GME_coeff = MIN(GME_coeff, GME_coeff_limiter) if ((CS%id_GME_coeff_h>0) .or. find_FrictWork) GME_coeff_h(i,j,k) = GME_coeff - str_xx_GME(i,j) = GME_coeff * sh_xx_bt(i,j) enddo ; enddo @@ -1137,16 +1130,14 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, do J=js-1,Jeq ; do I=is-1,Ieq - if ((max_diss_rate(i,j,k) > 0) .and. (grad_vel_mag_bt_q(i,j)>0) ) then - GME_coeff = (MIN(G%bathyT(i,j)/H0,1.0)**2) * FWfrac*max_diss_rate(i,j,k) / grad_vel_mag_bt_q(I,J) -! GME_coeff = (MIN(G%bathyT(i,j)/H0,1.0)**2) * FWfrac*target_diss_rate_GME(i,j,k) / grad_vel_mag_bt_q(I,J) + if ((grad_vel_mag_bt_q(I,J)>0) .and. (max_diss_rate_q(I,J,k)>0)) then + GME_coeff = (MIN(G%bathyT(i,j)/H0,1.0)**2) * FWfrac*max_diss_rate_q(I,J,k) / grad_vel_mag_bt_q(I,J) else GME_coeff = 0.0 endif ! apply mask - GME_coeff = GME_coeff * boundary_mask(i,j) - + GME_coeff = GME_coeff * boundary_mask_q(I,J) GME_coeff = MIN(GME_coeff, GME_coeff_limiter) if (CS%id_GME_coeff_q>0) GME_coeff_q(I,J,k) = GME_coeff @@ -1303,16 +1294,16 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, enddo ; enddo endif ! MEKE%backscatter - if (CS%use_GME) then - do j=js,je ; do i=is,ie - ! MEKE%mom_src now is sign definite because it only uses the dissipation - MEKE%mom_src(i,j) = MEKE%mom_src(i,j) + MAX(FrictWork_diss(i,j,k), FrictWorkMax(i,j,k)) - enddo ; enddo - else - do j=js,je ; do i=is,ie - MEKE%mom_src(i,j) = MEKE%mom_src(i,j) + FrictWork(i,j,k) - enddo ; enddo - endif ! CS%use_GME +! if (CS%use_GME) then +! do j=js,je ; do i=is,ie +! ! MEKE%mom_src now is sign definite because it only uses the dissipation +! MEKE%mom_src(i,j) = MEKE%mom_src(i,j) + MAX(FrictWork_diss(i,j,k), FrictWorkMax(i,j,k)) +! enddo ; enddo +! else + do j=js,je ; do i=is,ie + MEKE%mom_src(i,j) = MEKE%mom_src(i,j) + FrictWork(i,j,k) + enddo ; enddo +! endif ! CS%use_GME if (CS%use_GME .and. associated(MEKE)) then if (associated(MEKE%GME_snk)) then @@ -1405,6 +1396,9 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS) ! valid parameters. logical :: split ! If true, use the split time stepping scheme. ! If false and USE_GME = True, issue a FATAL error. + logical :: use_MEKE ! If true, use the MEKE module for calculating eddy kinetic energy. + ! If false and USE_GME = True, issue a FATAL error. + character(len=64) :: inputdir, filename real :: deg2rad ! Converts degrees to radians real :: slat_fn ! sin(lat)**Kh_pwr_of_sine @@ -1659,6 +1653,14 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS) do_not_log=.true.) if (.not. split) call MOM_error(FATAL,"ERROR: Currently, USE_GME = True "// & "cannot be used with SPLIT=False.") + + call get_param(param_file, mdl, "USE_MEKE", use_MEKE, & + "If true, turns on the MEKE scheme which calculates\n"// & + "a sub-grid mesoscale eddy kinetic energy budget.", & + default=.false.) + + if (.not. use_MEKE) call MOM_error(FATAL,"ERROR: Currently, USE_GME = True "// & + "cannot be used with USE_MEKE=False.") endif if (CS%bound_Kh .or. CS%bound_Ah .or. CS%better_bound_Kh .or. CS%better_bound_Ah) & From fa773cdefdcd3eb47b6a9f0e9d085eb91d83b1a3 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Mon, 26 Aug 2019 15:21:35 -0600 Subject: [PATCH 014/225] Created new variables dvdx3/dudy3 for shearing strain on the Laplacian --- src/parameterizations/lateral/MOM_hor_visc.F90 | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index 44131952a1..b83ddd4432 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -254,6 +254,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, real, dimension(SZIB_(G),SZJB_(G)) :: & dvdx, dudy, & ! components in the shearing strain [s-1] + dvdx3, dudy3, & ! components in the shearing strain on the Laplacian [m-2 s-1] dvdx_bt, dudy_bt, & ! components in the barotropic shearing strain [s-1] sh_xy, & ! horizontal shearing strain (du/dy + dv/dx) including metric terms [s-1] sh_xy_bt, & ! barotropic horizontal shearing strain (du/dy + dv/dx) inc. metric terms [s-1] @@ -921,8 +922,8 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if (CS%biharmonic) then ! Gradient of Laplacian, for use in bi-harmonic term do J=js-1,Jeq ; do I=is-1,Ieq - dvdx(I,J) = CS%DY_dxBu(I,J)*(v0(i+1,J)*G%IdyCv(i+1,J) - v0(i,J)*G%IdyCv(i,J)) - dudy(I,J) = CS%DX_dyBu(I,J)*(u0(I,j+1)*G%IdxCu(I,j+1) - u0(I,j)*G%IdxCu(I,j)) + dvdx3(I,J) = CS%DY_dxBu(I,J)*(v0(i+1,J)*G%IdyCv(i+1,J) - v0(i,J)*G%IdyCv(i,J)) + dudy3(I,J) = CS%DX_dyBu(I,J)*(u0(I,j+1)*G%IdxCu(I,j+1) - u0(I,j)*G%IdxCu(I,j)) enddo ; enddo ! Adjust contributions to shearing strain on open boundaries. if (apply_OBC) then ; if (OBC%zero_strain .or. OBC%freeslip_strain) then @@ -931,17 +932,17 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if (OBC%segment(n)%is_N_or_S .and. (J >= js-1) .and. (J <= Jeq)) then do I=OBC%segment(n)%HI%IsdB,OBC%segment(n)%HI%IedB if (OBC%zero_strain) then - dvdx(I,J) = 0. ; dudy(I,J) = 0. + dvdx3(I,J) = 0. ; dudy3(I,J) = 0. elseif (OBC%freeslip_strain) then - dudy(I,J) = 0. + dudy3(I,J) = 0. endif enddo elseif (OBC%segment(n)%is_E_or_W .and. (I >= is-1) .and. (I <= Ieq)) then do J=OBC%segment(n)%HI%JsdB,OBC%segment(n)%HI%JedB if (OBC%zero_strain) then - dvdx(I,J) = 0. ; dudy(I,J) = 0. + dvdx3(I,J) = 0. ; dudy3(I,J) = 0. elseif (OBC%freeslip_strain) then - dvdx(I,J) = 0. + dvdx3(I,J) = 0. endif enddo endif @@ -1072,10 +1073,10 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if (CS%id_Ah_q>0) Ah_q(I,J,k) = Ah - str_xy(I,J) = str_xy(I,J) + Ah * ( dvdx(I,J) + dudy(I,J) ) + str_xy(I,J) = str_xy(I,J) + Ah * ( dvdx3(I,J) + dudy3(I,J) ) ! Keep a copy of the biharmonic contribution for backscatter parameterization - bhstr_xy(I,J) = Ah * ( dvdx(I,J) + dudy(I,J) ) * & + bhstr_xy(I,J) = Ah * ( dvdx3(I,J) + dudy3(I,J) ) * & (hq(I,J) * G%mask2dBu(I,J) * CS%reduction_xy(I,J)) endif ! biharmonic From 2ac12ce6a5f609532eb519fabd18092c06cc2ef3 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Mon, 26 Aug 2019 15:24:33 -0600 Subject: [PATCH 015/225] Deleting comments --- .../lateral/MOM_hor_visc.F90 | 55 ------------------- 1 file changed, 55 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index b83ddd4432..d4fa977cb9 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -561,34 +561,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, endif endif -! if (CS%use_GME) then -! call pass_var(dudx, G%Domain, complete=.true.) -! call pass_var(dvdy, G%Domain, complete=.true.) -! call pass_var(dvdx, G%Domain, position=CORNER, complete=.true.) -! call pass_var(dudy, G%Domain, position=CORNER, complete=.true.) -!! call pass_var(MEKE%MEKE, G%Domain, complete=.true.) -! - -! do j=Jsq-1,Jeq+2 ; do i=Isq-1,Ieq+2 -! do J=Jsq,Jeq+1 ; do i=Isq,Ieq+1 -! do j = js, je ; do i = is, ie -! grad_vel_mag_h(i,j) = boundary_mask_h(i,j) * (dudx(i,j)**2 + dvdy(i,j)**2 + & -! (0.25*((dvdx(I,J)+dvdx(I-1,J-1))+(dvdx(I,J-1)+dvdx(I-1,J))))**2 + & -! (0.25*((dudy(I,J)+dudy(I-1,J-1))+(dudy(I,J-1)+dudy(I-1,J))))**2) -! -! max_diss_rate_h(i,j,k) = 2.0 !*MEKE%MEKE(i,j) * sqrt(grad_vel_mag_h(i,j)) -! enddo ; enddo -! -! do J = G%JscB, G%JecB ; do I = G%IscB, G%IecB -! grad_vel_mag_q(I,J) = boundary_mask_h(I,J) * (dudx(i,j)**2 + dvdy(i,j)**2 + & -! (0.25*((dvdx(I,J)+dvdx(I-1,J-1))+(dvdx(I,J-1)+dvdx(I-1,J))))**2 + & -! (0.25*((dudy(I,J)+dudy(I-1,J-1))+(dudy(I,J-1)+dudy(I-1,J))))**2) -! -! max_diss_rate_q(I,J,k) = 0.5*(MEKE%MEKE(i,j)+MEKE%MEKE(i+1,j)+ & -! MEKE%MEKE(i,j+1)+MEKE%MEKE(i+1,j+1)) * sqrt(grad_vel_mag_q(I,J)) -! enddo ; enddo -! endif - if (OBC%segment(n)%direction == OBC_DIRECTION_N) then ! There are extra wide halos here to accommodate the cross-corner-point @@ -690,29 +662,14 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if ((CS%Leith_Kh) .or. (CS%Leith_Ah)) then - ! Components for the vertical vorticity - ! Note this a simple re-calculation of shearing components using the same discretization. - ! We will consider using a circulation based calculation of vorticity later. - ! Also note this will need OBC boundary conditions re-applied... -! do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 -! DY_dxBu = G%dyBu(I,J) * G%IdxBu(I,J) -! dvdx(I,J) = DY_dxBu * (v(i+1,J,k) * G%IdyCv(i+1,J) - v(i,J,k) * G%IdyCv(i,J)) -! DX_dyBu = G%dxBu(I,J) * G%IdyBu(I,J) -! dudy(I,J) = DX_dyBu * (u(I,j+1,k) * G%IdxCu(I,j+1) - u(I,j,k) * G%IdxCu(I,j)) -! enddo ; enddo - ! Vorticity if (CS%no_slip) then do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 vort_xy(I,J) = (2.0-G%mask2dBu(I,J)) * ( dvdx(I,J) - dudy(I,J) ) -! dudy(I,J) = (2.0-G%mask2dBu(I,J)) * dudy(I,J) -! dvdx(I,J) = (2.0-G%mask2dBu(I,J)) * dvdx(I,J) enddo ; enddo else do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 vort_xy(I,J) = G%mask2dBu(I,J) * ( dvdx(I,J) - dudy(I,J) ) -! dudy(I,J) = G%mask2dBu(I,J) * dudy(I,J) -! dvdx(I,J) = G%mask2dBu(I,J) * dvdx(I,J) enddo ; enddo endif @@ -1085,11 +1042,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if (CS%use_GME) then -! call pass_var(dudx, G%Domain, complete=.true.) -! call pass_var(dvdy, G%Domain, complete=.true.) -! call pass_var(dvdx, G%Domain, position=CORNER, complete=.true.) -! call pass_var(dudy, G%Domain, position=CORNER, complete=.true.) - do j = js, je ; do i = is, ie grad_vel_mag_h(i,j) = boundary_mask_h(i,j) * (dudx(i,j)**2 + dvdy(i,j)**2 + & @@ -1295,16 +1247,9 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, enddo ; enddo endif ! MEKE%backscatter -! if (CS%use_GME) then -! do j=js,je ; do i=is,ie -! ! MEKE%mom_src now is sign definite because it only uses the dissipation -! MEKE%mom_src(i,j) = MEKE%mom_src(i,j) + MAX(FrictWork_diss(i,j,k), FrictWorkMax(i,j,k)) -! enddo ; enddo -! else do j=js,je ; do i=is,ie MEKE%mom_src(i,j) = MEKE%mom_src(i,j) + FrictWork(i,j,k) enddo ; enddo -! endif ! CS%use_GME if (CS%use_GME .and. associated(MEKE)) then if (associated(MEKE%GME_snk)) then From d29a4353fd10bab4e4b1215bb4d3df5b8b414a51 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Tue, 27 Aug 2019 11:05:23 -0600 Subject: [PATCH 016/225] Adds consistency check between mesh and mom6 grid Check for consistency of lat, lon and mask between mesh and mom6 grid --- config_src/nuopc_driver/mom_cap.F90 | 83 ++++++++++++++++++++++++++++- 1 file changed, 81 insertions(+), 2 deletions(-) diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index aa034e7d36..1d0b41b550 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -366,7 +366,7 @@ module MOM_cap_mod use ESMF, only: ESMF_METHOD_INITIALIZE, ESMF_MethodRemove, ESMF_State use ESMF, only: ESMF_LOGMSG_INFO, ESMF_RC_ARG_BAD, ESMF_VM, ESMF_Time use ESMF, only: ESMF_TimeInterval, ESMF_MAXSTR, ESMF_VMGetCurrent -use ESMF, only: ESMF_VMGet, ESMF_TimeGet, ESMF_TimeIntervalGet +use ESMF, only: ESMF_VMGet, ESMF_TimeGet, ESMF_TimeIntervalGet, ESMF_MeshGet use ESMF, only: ESMF_MethodExecute, ESMF_Mesh, ESMF_DeLayout, ESMF_Distgrid use ESMF, only: ESMF_DistGridConnection, ESMF_StateItem_Flag, ESMF_KIND_I4 use ESMF, only: ESMF_KIND_I8, ESMF_FAILURE, ESMF_DistGridCreate, ESMF_MeshCreate @@ -378,7 +378,8 @@ module MOM_cap_mod use ESMF, only: ESMF_FieldCreate, ESMF_LOGMSG_ERROR, ESMF_LOGMSG_WARNING use ESMF, only: ESMF_COORDSYS_SPH_DEG, ESMF_GridCreate, ESMF_INDEX_DELOCAL use ESMF, only: ESMF_MESHLOC_ELEMENT, ESMF_RC_VAL_OUTOFRANGE, ESMF_StateGet -use ESMF, only: ESMF_TimePrint, ESMF_AlarmSet, ESMF_FieldGet +use ESMF, only: ESMF_TimePrint, ESMF_AlarmSet, ESMF_FieldGet, ESMF_Array +use ESMF, only: ESMF_ArrayCreate use ESMF, only: operator(==), operator(/=), operator(+), operator(-) ! TODO ESMF_GridCompGetInternalState does not have an explicit Fortran interface. @@ -1181,6 +1182,14 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) character(len=128) :: fldname character(len=256) :: cvalue character(len=*), parameter :: subname='(MOM_cap:InitializeRealize)' + integer :: spatialDim + integer :: numOwnedElements + type(ESMF_Array) :: elemMaskArray + real(ESMF_KIND_R8) , pointer :: ownedElemCoords(:) + real(ESMF_KIND_R8) , pointer :: lat(:), latMesh(:) + real(ESMF_KIND_R8) , pointer :: lon(:), lonMesh(:) + integer(ESMF_KIND_I4) , pointer :: mask(:), maskMesh(:) + real(ESMF_KIND_R8) :: diff_lon, diff_lat !-------------------------------- rc = ESMF_SUCCESS @@ -1326,6 +1335,76 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) file=__FILE__)) & return + ! Check for consistency of lat, lon and mask between mesh and mom6 grid + call ESMF_MeshGet(Emesh, spatialDim=spatialDim, numOwnedElements=numOwnedElements, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & + line=__LINE__, & + file=__FILE__)) & + return + + allocate(ownedElemCoords(spatialDim*numOwnedElements)) + allocate(lonMesh(numOwnedElements), lon(numOwnedElements)) + allocate(latMesh(numOwnedElements), lat(numOwnedElements)) + allocate(maskMesh(numOwnedElements), mask(numOwnedElements)) + + call ESMF_MeshGet(Emesh, ownedElemCoords=ownedElemCoords, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & + line=__LINE__, & + file=__FILE__)) & + return + do n = 1,numOwnedElements + lonMesh(n) = ownedElemCoords(2*n-1) + latMesh(n) = ownedElemCoords(2*n) + end do + + elemMaskArray = ESMF_ArrayCreate(Distgrid, maskMesh, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & + line=__LINE__, & + file=__FILE__)) & + return + call ESMF_MeshGet(Emesh, elemMaskArray=elemMaskArray, rc=rc) + if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & + line=__LINE__, & + file=__FILE__)) & + return + + call mpp_get_compute_domain(ocean_public%domain, isc, iec, jsc, jec) + n = 0 + do j = jsc, jec + jg = j + ocean_grid%jsc - jsc + do i = isc, iec + ig = i + ocean_grid%isc - isc + n = n+1 + mask(n) = ocean_grid%mask2dT(ig,jg) + lon(n) = ocean_grid%geolonT(ig,jg) + lat(n) = ocean_grid%geolatT(ig,jg) + end do + end do + + do n = 1,numOwnedElements + diff_lon = abs(lonMesh(n) - lon(n)) + if (diff_lon > 1.e-2) then + write(6,100)n,lonMesh(n),lon(n), diff_lon +100 format('ERROR: MOM n, lonMesh(n), lon(n), diff_lon = ',i8,2(f21.13,3x),d21.5) + !call shr_sys_abort() + end if + diff_lat = abs(latMesh(n) - lat(n)) + if (diff_lat > 1.e-2) then + write(6,101)n,latMesh(n),lat(n), diff_lat +101 format('ERROR: CICE n, latMesh(n), lat(n), diff_lat = ',i8,2(f21.13,3x),d21.5) + !call shr_sys_abort() + end if + if (abs(maskMesh(n) - mask(n)) > 0) then + write(6,102)n,maskMesh(n),mask(n) +102 format('ERROR: CICE n, maskMesh(n), mask(n) = ',3(i8,2x)) + !call shr_sys_abort() + end if + end do + + deallocate(ownedElemCoords) + deallocate(lonMesh , lon ) + deallocate(latMesh , lat ) + deallocate(maskMesh, mask) ! realize the import and export fields using the mesh call MOM_RealizeFields(importState, fldsToOcn_num, fldsToOcn, "Ocn import", mesh=Emesh, rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & From fb1133e7d2b3edeadda55941bcdbbab238f84f80 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Tue, 27 Aug 2019 11:14:36 -0600 Subject: [PATCH 017/225] Replaces CICE to MOM in error messages --- config_src/nuopc_driver/mom_cap.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index 1d0b41b550..360e8f833d 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -1391,12 +1391,12 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) diff_lat = abs(latMesh(n) - lat(n)) if (diff_lat > 1.e-2) then write(6,101)n,latMesh(n),lat(n), diff_lat -101 format('ERROR: CICE n, latMesh(n), lat(n), diff_lat = ',i8,2(f21.13,3x),d21.5) +101 format('ERROR: MOM n, latMesh(n), lat(n), diff_lat = ',i8,2(f21.13,3x),d21.5) !call shr_sys_abort() end if if (abs(maskMesh(n) - mask(n)) > 0) then write(6,102)n,maskMesh(n),mask(n) -102 format('ERROR: CICE n, maskMesh(n), mask(n) = ',3(i8,2x)) +102 format('ERROR: MOM n, maskMesh(n), mask(n) = ',3(i8,2x)) !call shr_sys_abort() end if end do From 7ad72961103f4b436579dc34f78ab43d9d774641 Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Wed, 28 Aug 2019 15:30:09 -0600 Subject: [PATCH 018/225] remove labeled format specifiers --- config_src/nuopc_driver/mom_cap.F90 | 13 +++++++------ 1 file changed, 7 insertions(+), 6 deletions(-) diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index 360e8f833d..39ca8c2f66 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -1181,6 +1181,7 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) integer, allocatable :: gindex(:) ! global index space character(len=128) :: fldname character(len=256) :: cvalue + character(len=256) :: frmt ! format specifier for several error msgs character(len=*), parameter :: subname='(MOM_cap:InitializeRealize)' integer :: spatialDim integer :: numOwnedElements @@ -1384,19 +1385,19 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) do n = 1,numOwnedElements diff_lon = abs(lonMesh(n) - lon(n)) if (diff_lon > 1.e-2) then - write(6,100)n,lonMesh(n),lon(n), diff_lon -100 format('ERROR: MOM n, lonMesh(n), lon(n), diff_lon = ',i8,2(f21.13,3x),d21.5) + frmt = "('ERROR: MOM n, lonMesh(n), lon(n), diff_lon = ',i8,2(f21.13,3x),d21.5)" + write(6,frmt)n,lonMesh(n),lon(n), diff_lon !call shr_sys_abort() end if diff_lat = abs(latMesh(n) - lat(n)) if (diff_lat > 1.e-2) then - write(6,101)n,latMesh(n),lat(n), diff_lat -101 format('ERROR: MOM n, latMesh(n), lat(n), diff_lat = ',i8,2(f21.13,3x),d21.5) + frmt = "('ERROR: MOM n, latMesh(n), lat(n), diff_lat = ',i8,2(f21.13,3x),d21.5)" + write(6,frmt)n,latMesh(n),lat(n), diff_lat !call shr_sys_abort() end if if (abs(maskMesh(n) - mask(n)) > 0) then - write(6,102)n,maskMesh(n),mask(n) -102 format('ERROR: MOM n, maskMesh(n), mask(n) = ',3(i8,2x)) + frmt = "('ERROR: MOM n, maskMesh(n), mask(n) = ',3(i8,2x))" + write(6,frmt)n,maskMesh(n),mask(n) !call shr_sys_abort() end if end do From cc9e33d4da0eabacfbeeec826eff9e46c9ae01a6 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Wed, 28 Aug 2019 15:46:31 -0600 Subject: [PATCH 019/225] Adds new parameter, MEKE_EQUILIBRIUM_ALT This will use an alternative formula for computing the equilibrium (initial) value of MEKE. --- src/parameterizations/lateral/MOM_MEKE.F90 | 12 ++++++++++-- 1 file changed, 10 insertions(+), 2 deletions(-) diff --git a/src/parameterizations/lateral/MOM_MEKE.F90 b/src/parameterizations/lateral/MOM_MEKE.F90 index 8f947c4f68..3688c3dfea 100644 --- a/src/parameterizations/lateral/MOM_MEKE.F90 +++ b/src/parameterizations/lateral/MOM_MEKE.F90 @@ -45,6 +45,8 @@ module MOM_MEKE logical :: visc_drag !< If true use the vertvisc_type to calculate bottom drag. logical :: MEKE_GEOMETRIC !< If true, uses the GM coefficient formulation from the GEOMETRIC !! framework (Marshall et al., 2012) + logical :: MEKE_equilibrium_alt !< If true, use an alternative calculation for the + !! equilibrium value of MEKE. logical :: GM_src_alt !< If true, use the GM energy conversion form S^2*N^2*kappa rather !! than the streamfunction for the MEKE GM source term. logical :: Rd_as_max_scale !< If true the length scale can not exceed the @@ -744,8 +746,11 @@ subroutine MEKE_equilibrium(CS, MEKE, G, GV, US, SN_u, SN_v, drag_rate_visc, I_m else EKE = 0. endif - MEKE%MEKE(i,j) = EKE -! MEKE%MEKE(i,j) = (US%Z_to_m*G%bathyT(i,j)*SN / (8*CS%cdrag))**2 + if (CS%MEKE_equilibrium_alt) then + MEKE%MEKE(i,j) = (US%Z_to_m*G%bathyT(i,j)*SN / (8*CS%cdrag))**2 + else + MEKE%MEKE(i,j) = EKE + endif enddo ; enddo end subroutine MEKE_equilibrium @@ -973,6 +978,9 @@ logical function MEKE_init(Time, G, param_file, diag, CS, MEKE, restart_CS) call get_param(param_file, mdl, "MEKE_GEOMETRIC", CS%MEKE_GEOMETRIC, & "If MEKE_GEOMETRIC is true, uses the GM coefficient formulation "//& "from the GEOMETRIC framework (Marshall et al., 2012).", default=.false.) + call get_param(param_file, mdl, "MEKE_EQUILIBRIUM_ALT", CS%MEKE_equilibrium_alt, & + "If true, use an alternative formula for computing the (equilibrium)"//& + "initial value of MEKE.", default=.false.) call get_param(param_file, mdl, "MEKE_FRCOEFF", CS%MEKE_FrCoeff, & "The efficiency of the conversion of mean energy into "//& "MEKE. If MEKE_FRCOEFF is negative, this conversion "//& From d05235f901ca8890c04b69d250a7f3333f143e66 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Wed, 28 Aug 2019 15:48:09 -0600 Subject: [PATCH 020/225] Final GME code fix Changed the loop indices on boundary mask. Added three new runtime parameters for GME (GME_H0, GME_EFFICIENCY, GME_LIMITER) --- .../lateral/MOM_hor_visc.F90 | 40 ++++++++++++------- 1 file changed, 25 insertions(+), 15 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index d4fa977cb9..8074fb2f7d 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -79,6 +79,10 @@ module MOM_hor_visc logical :: res_scale_MEKE !< If true, the viscosity contribution from MEKE is scaled by !! the resolution function. logical :: use_GME !< If true, use GME backscatter scheme. + real :: GME_h0 !< The strength of GME tapers quadratically to zero when the bathymetric + !! depth is shallower than GME_H0 [m] + real :: GME_efficiency !< The nondimensional prefactor multiplying the GME coefficient [nondim] + real :: GME_limiter !< The absolute maximum value the GME coefficient is allowed to take [m2 s-1]. real ALLOCABLE_, dimension(NIMEM_,NJMEM_) :: Kh_bg_xx !< The background Laplacian viscosity at h points [m2 s-1]. @@ -382,23 +386,15 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, use_MEKE_Ku = associated(MEKE%Ku) use_MEKE_Au = associated(MEKE%Au) - do j = G%jsc, G%jec ; do i = G%isc, G%iec + do j=Jsq-1,Jeq+2 ; do i=Isq-1,Ieq+2 boundary_mask_h(i,j) = (G%mask2dCu(I,j) * G%mask2dCv(i,J) * G%mask2dCu(I-1,j) * G%mask2dCv(i,J-1)) enddo ; enddo - do J = G%JscB, G%JecB ; do I = G%IscB, G%IecB + do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 boundary_mask_q(I,J) = (G%mask2dCv(i,J) * G%mask2dCv(i+1,J) * G%mask2dCu(I,j) * G%mask2dCu(I,j-1)) enddo; enddo -! call pass_var(boundary_mask_h, G%Domain, complete=.true.) - - if (CS%use_GME) then - ! GME tapers off above this depth - H0 = 1000.0 - FWfrac = 1.0 - GME_coeff_limiter = 1e7 - ! initialize diag. array with zeros GME_coeff_h(:,:,:) = 0.0 GME_coeff_q(:,:,:) = 0.0 @@ -450,7 +446,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, enddo ; enddo do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 - grad_vel_mag_bt_q(I,J) = boundary_mask_h(I,J) * (dvdx_bt(i,j)**2 + dudy_bt(i,j)**2 + & + grad_vel_mag_bt_q(I,J) = boundary_mask_q(I,J) * (dvdx_bt(i,j)**2 + dudy_bt(i,j)**2 + & (0.25*((dudx_bt(i,j)+dudx_bt(i+1,j+1))+(dudx_bt(i,j+1)+dudx_bt(i+1,j))))**2 + & (0.25*((dvdy_bt(i,j)+dvdy_bt(i+1,j+1))+(dvdy_bt(i,j+1)+dvdy_bt(i+1,j))))**2) enddo ; enddo @@ -1066,14 +1062,14 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, do J=Jsq,Jeq+1 ; do i=Isq,Ieq+1 if ((grad_vel_mag_bt_h(i,j)>0) .and. (max_diss_rate_h(i,j,k)>0)) then - GME_coeff = (MIN(G%bathyT(i,j)/H0,1.0)**2) * FWfrac*max_diss_rate_h(i,j,k) / grad_vel_mag_bt_h(i,j) + GME_coeff = (MIN(G%bathyT(i,j)/CS%GME_h0,1.0)**2) * CS%GME_efficiency*max_diss_rate_h(i,j,k) / grad_vel_mag_bt_h(i,j) else GME_coeff = 1.0 endif ! apply mask GME_coeff = GME_coeff * boundary_mask_h(i,j) - GME_coeff = MIN(GME_coeff, GME_coeff_limiter) + GME_coeff = MIN(GME_coeff, CS%GME_limiter) if ((CS%id_GME_coeff_h>0) .or. find_FrictWork) GME_coeff_h(i,j,k) = GME_coeff str_xx_GME(i,j) = GME_coeff * sh_xx_bt(i,j) @@ -1084,14 +1080,14 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, do J=js-1,Jeq ; do I=is-1,Ieq if ((grad_vel_mag_bt_q(I,J)>0) .and. (max_diss_rate_q(I,J,k)>0)) then - GME_coeff = (MIN(G%bathyT(i,j)/H0,1.0)**2) * FWfrac*max_diss_rate_q(I,J,k) / grad_vel_mag_bt_q(I,J) + GME_coeff = (MIN(G%bathyT(i,j)/CS%GME_h0,1.0)**2) * CS%GME_efficiency*max_diss_rate_q(I,J,k) / grad_vel_mag_bt_q(I,J) else GME_coeff = 0.0 endif ! apply mask GME_coeff = GME_coeff * boundary_mask_q(I,J) - GME_coeff = MIN(GME_coeff, GME_coeff_limiter) + GME_coeff = MIN(GME_coeff, CS%GME_limiter) if (CS%id_GME_coeff_q>0) GME_coeff_q(I,J,k) = GME_coeff str_xy_GME(I,J) = GME_coeff * sh_xy_bt(I,J) @@ -1607,6 +1603,20 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS) if (.not. use_MEKE) call MOM_error(FATAL,"ERROR: Currently, USE_GME = True "// & "cannot be used with USE_MEKE=False.") + + call get_param(param_file, mdl, "GME_H0", CS%GME_h0, & + "The strength of GME tapers quadratically to zero when the bathymetric "//& + "depth is shallower than GME_H0.", units="m", & + default=1000.0) + + call get_param(param_file, mdl, "GME_EFFICIENCY", CS%GME_efficiency, & + "The nondimensional prefactor multiplying the GME coefficient.", & + units="nondim", default=1.0) + + call get_param(param_file, mdl, "GME_LIMITER", CS%GME_limiter, & + "The absolute maximum value the GME coefficient is allowed to take.", & + units="m2 s-1", default=1.0e7) + endif if (CS%bound_Kh .or. CS%bound_Ah .or. CS%better_bound_Kh .or. CS%better_bound_Ah) & From 18c6f1ebb8b3d27c4b666246c4717a9450417196 Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Wed, 28 Aug 2019 16:59:08 -0600 Subject: [PATCH 021/225] replace hard-coded mesh diff. limit with a param: eps_omesh --- config_src/nuopc_driver/mom_ocean_model_nuopc.F90 | 15 ++++++++++++++- 1 file changed, 14 insertions(+), 1 deletion(-) diff --git a/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 b/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 index 95b1bcc6e3..abb6bb2679 100644 --- a/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 +++ b/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 @@ -80,6 +80,7 @@ module MOM_ocean_model_nuopc public ocean_public_type_chksum public ocean_model_data_get public get_ocean_grid +public get_eps_omesh !> This interface extracts a named scalar field or array from the ocean surface or public type interface ocean_model_data_get @@ -182,6 +183,9 @@ module MOM_ocean_model_nuopc logical :: diabatic_first !< If true, apply diabatic and thermodynamic !! processes before time stepping the dynamics. + real :: eps_omesh !< Max allowable difference between ESMF mesh and MOM6 + !! domain coordinates + type(directories) :: dirs !< A structure containing several relevant directory paths. type(mech_forcing) :: forces !< A structure with the driving mechanical surface forces type(forcing) :: fluxes !< A structure containing pointers to @@ -327,6 +331,9 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i else ; call MOM_error(FATAL,"ocean_model_init: OCEAN_SURFACE_STAGGER = "// & trim(stagger)//" is invalid.") ; endif + call get_param(param_file, mdl, "EPS_OMESH",OS%eps_omesh, & + "Maximum allowable difference between ESMF mesh and "//& + "MOM6 domain coordinates.", default=1.e-2) call get_param(param_file, mdl, "RESTORE_SALINITY",OS%restore_salinity, & "If true, the coupled driver will add a globally-balanced "//& "fresh-water flux that drives sea-surface salinity "//& @@ -353,7 +360,7 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i OS%press_to_z = 1.0/(Rho0*G_Earth) - call get_param(param_file, mdl, "HFREEZE", HFrz, & + call get_param(param_file, mdl, "HFREEZE", HFrz, & "If HFREEZE > 0, melt potential will be computed. The actual depth "//& "over which melt potential is computed will be min(HFREEZE, OBLD), "//& "where OBLD is the boundary layer depth. If HFREEZE <= 0 (default), "//& @@ -1174,4 +1181,10 @@ subroutine get_ocean_grid(OS, Gridp) return end subroutine get_ocean_grid +!> Returns eps_omesh read from param file +real function get_eps_omesh(OS) + type(ocean_state_type) :: OS + get_eps_omesh = OS%eps_omesh; return +end function + end module MOM_ocean_model_nuopc From 66f808bf8aec38fef9f5d5c293b8e91d7acd6fa5 Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Wed, 28 Aug 2019 17:20:39 -0600 Subject: [PATCH 022/225] add units to EPS_OMESH --- config_src/nuopc_driver/mom_ocean_model_nuopc.F90 | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 b/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 index abb6bb2679..83dfa0a4f8 100644 --- a/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 +++ b/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 @@ -333,7 +333,8 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i call get_param(param_file, mdl, "EPS_OMESH",OS%eps_omesh, & "Maximum allowable difference between ESMF mesh and "//& - "MOM6 domain coordinates.", default=1.e-2) + "MOM6 domain coordinates in nuopc cap.", & + units="degrees", default=1.e-2) call get_param(param_file, mdl, "RESTORE_SALINITY",OS%restore_salinity, & "If true, the coupled driver will add a globally-balanced "//& "fresh-water flux that drives sea-surface salinity "//& From c698b1a2864b3e2bf6006a5729c876506c5943ae Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Wed, 28 Aug 2019 17:21:34 -0600 Subject: [PATCH 023/225] call get_eps_omesh and remove hardcoded limit --- config_src/nuopc_driver/mom_cap.F90 | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index 39ca8c2f66..642a7eb809 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -344,7 +344,8 @@ module MOM_cap_mod use MOM_grid, only: ocean_grid_type, get_global_grid_size use MOM_ocean_model_nuopc, only: ocean_model_restart, ocean_public_type, ocean_state_type use MOM_ocean_model_nuopc, only: ocean_model_init_sfc -use MOM_ocean_model_nuopc, only: ocean_model_init, update_ocean_model, ocean_model_end, get_ocean_grid +use MOM_ocean_model_nuopc, only: ocean_model_init, update_ocean_model, ocean_model_end +use MOM_ocean_model_nuopc, only: get_ocean_grid, get_eps_omesh use MOM_cap_time, only: AlarmInit use MOM_cap_methods, only: mom_import, mom_export, mom_set_geomtype #ifdef CESMCOUPLED @@ -1191,6 +1192,7 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) real(ESMF_KIND_R8) , pointer :: lon(:), lonMesh(:) integer(ESMF_KIND_I4) , pointer :: mask(:), maskMesh(:) real(ESMF_KIND_R8) :: diff_lon, diff_lat + real :: eps_omesh !-------------------------------- rc = ESMF_SUCCESS @@ -1382,15 +1384,16 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) end do end do + eps_omesh = get_eps_omesh(ocean_state) do n = 1,numOwnedElements diff_lon = abs(lonMesh(n) - lon(n)) - if (diff_lon > 1.e-2) then + if (diff_lon > eps_omesh) then frmt = "('ERROR: MOM n, lonMesh(n), lon(n), diff_lon = ',i8,2(f21.13,3x),d21.5)" write(6,frmt)n,lonMesh(n),lon(n), diff_lon !call shr_sys_abort() end if diff_lat = abs(latMesh(n) - lat(n)) - if (diff_lat > 1.e-2) then + if (diff_lat > eps_omesh) then frmt = "('ERROR: MOM n, latMesh(n), lat(n), diff_lat = ',i8,2(f21.13,3x),d21.5)" write(6,frmt)n,latMesh(n),lat(n), diff_lat !call shr_sys_abort() From d9a00d78175cc034aff47ae4fb0ff4f9f57af897 Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Wed, 28 Aug 2019 17:41:01 -0600 Subject: [PATCH 024/225] call MOM_error if mesh is inconsistent --- config_src/nuopc_driver/mom_cap.F90 | 27 +++++++++++++++------------ 1 file changed, 15 insertions(+), 12 deletions(-) diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index 642a7eb809..25dd9bde76 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -322,7 +322,6 @@ module MOM_cap_mod use mpp_domains_mod, only: mpp_get_ntile_count, mpp_get_pelist, mpp_get_global_domain use mpp_domains_mod, only: mpp_get_domain_npes use mpp_io_mod, only: mpp_open, MPP_RDONLY, MPP_ASCII, MPP_OVERWR, MPP_APPEND, mpp_close, MPP_SINGLE -use mpp_mod, only: input_nml_file, mpp_error, FATAL, NOTE, mpp_pe, mpp_npes, mpp_set_current_pelist use mpp_mod, only: stdlog, stdout, mpp_root_pe, mpp_clock_id use mpp_mod, only: mpp_clock_begin, mpp_clock_end, MPP_CLOCK_SYNC use mpp_mod, only: MPP_CLOCK_DETAILED, CLOCK_COMPONENT, MAXPES @@ -339,7 +338,7 @@ module MOM_cap_mod use MOM_file_parser, only: get_param, log_version, param_file_type, close_param_file use MOM_get_input, only: Get_MOM_Input, directories use MOM_domains, only: pass_var -use MOM_error_handler, only: is_root_pe +use MOM_error_handler, only: MOM_error, FATAL, is_root_pe use MOM_ocean_model_nuopc, only: ice_ocean_boundary_type use MOM_grid, only: ocean_grid_type, get_global_grid_size use MOM_ocean_model_nuopc, only: ocean_model_restart, ocean_public_type, ocean_state_type @@ -1182,7 +1181,8 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) integer, allocatable :: gindex(:) ! global index space character(len=128) :: fldname character(len=256) :: cvalue - character(len=256) :: frmt ! format specifier for several error msgs + character(len=256) :: frmt ! format specifier for several error msgs + character(len=512) :: err_msg ! error messages character(len=*), parameter :: subname='(MOM_cap:InitializeRealize)' integer :: spatialDim integer :: numOwnedElements @@ -1388,20 +1388,23 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) do n = 1,numOwnedElements diff_lon = abs(lonMesh(n) - lon(n)) if (diff_lon > eps_omesh) then - frmt = "('ERROR: MOM n, lonMesh(n), lon(n), diff_lon = ',i8,2(f21.13,3x),d21.5)" - write(6,frmt)n,lonMesh(n),lon(n), diff_lon - !call shr_sys_abort() + frmt = "('ERROR: Inconsistent coords - "//& + "MOM n, lonMesh(n), lon(n), diff_lon = ',i8,2(f21.13,3x),d21.5)" + write(err_msg, frmt)n,lonMesh(n),lon(n), diff_lon + call MOM_error(FATAL, err_msg) end if diff_lat = abs(latMesh(n) - lat(n)) if (diff_lat > eps_omesh) then - frmt = "('ERROR: MOM n, latMesh(n), lat(n), diff_lat = ',i8,2(f21.13,3x),d21.5)" - write(6,frmt)n,latMesh(n),lat(n), diff_lat - !call shr_sys_abort() + frmt = "('ERROR: Inconsistent coords - "//& + "MOM n, latMesh(n), lat(n), diff_lat = ',i8,2(f21.13,3x),d21.5)" + write(err_msg, frmt)n,latMesh(n),lat(n), diff_lat + call MOM_error(FATAL, err_msg) end if if (abs(maskMesh(n) - mask(n)) > 0) then - frmt = "('ERROR: MOM n, maskMesh(n), mask(n) = ',3(i8,2x))" - write(6,frmt)n,maskMesh(n),mask(n) - !call shr_sys_abort() + frmt = "('ERROR: Inconsistent masks - "//& + "MOM n, maskMesh(n), mask(n) = ',3(i8,2x))" + write(err_msg, frmt)n,maskMesh(n),mask(n) + call MOM_error(FATAL, err_msg) end if end do From e02b358ad6cb742de74cb84f7d6033ca9ad01779 Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Thu, 29 Aug 2019 11:46:23 -0600 Subject: [PATCH 025/225] compare mod 360 of diff and eps --- config_src/nuopc_driver/mom_cap.F90 | 2 +- config_src/nuopc_driver/mom_ocean_model_nuopc.F90 | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index 25dd9bde76..d9400e0560 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -1386,7 +1386,7 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) eps_omesh = get_eps_omesh(ocean_state) do n = 1,numOwnedElements - diff_lon = abs(lonMesh(n) - lon(n)) + diff_lon = abs(mod(lonMesh(n) - lon(n),360.0)) if (diff_lon > eps_omesh) then frmt = "('ERROR: Inconsistent coords - "//& "MOM n, lonMesh(n), lon(n), diff_lon = ',i8,2(f21.13,3x),d21.5)" diff --git a/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 b/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 index 83dfa0a4f8..426c7e9922 100644 --- a/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 +++ b/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 @@ -334,7 +334,7 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i call get_param(param_file, mdl, "EPS_OMESH",OS%eps_omesh, & "Maximum allowable difference between ESMF mesh and "//& "MOM6 domain coordinates in nuopc cap.", & - units="degrees", default=1.e-2) + units="degrees", default=1.e-4) call get_param(param_file, mdl, "RESTORE_SALINITY",OS%restore_salinity, & "If true, the coupled driver will add a globally-balanced "//& "fresh-water flux that drives sea-surface salinity "//& From 7f7e8c081b6660031c0c6233c8b82c82607cb44f Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Thu, 29 Aug 2019 15:28:43 -0600 Subject: [PATCH 026/225] Add EPS_OMESH value in error message --- config_src/nuopc_driver/mom_cap.F90 | 16 +++++++++------- 1 file changed, 9 insertions(+), 7 deletions(-) diff --git a/config_src/nuopc_driver/mom_cap.F90 b/config_src/nuopc_driver/mom_cap.F90 index d9400e0560..0caee9510e 100644 --- a/config_src/nuopc_driver/mom_cap.F90 +++ b/config_src/nuopc_driver/mom_cap.F90 @@ -1388,20 +1388,22 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) do n = 1,numOwnedElements diff_lon = abs(mod(lonMesh(n) - lon(n),360.0)) if (diff_lon > eps_omesh) then - frmt = "('ERROR: Inconsistent coords - "//& - "MOM n, lonMesh(n), lon(n), diff_lon = ',i8,2(f21.13,3x),d21.5)" - write(err_msg, frmt)n,lonMesh(n),lon(n), diff_lon + frmt = "('ERROR: Difference between ESMF Mesh and MOM6 domain coords is "//& + "greater than parameter EPS_OMESH. n, lonMesh(n), lon(n), diff_lon, "//& + "EPS_OMESH= ',i8,2(f21.13,3x),2(d21.5))" + write(err_msg, frmt)n,lonMesh(n),lon(n), diff_lon, eps_omesh call MOM_error(FATAL, err_msg) end if diff_lat = abs(latMesh(n) - lat(n)) if (diff_lat > eps_omesh) then - frmt = "('ERROR: Inconsistent coords - "//& - "MOM n, latMesh(n), lat(n), diff_lat = ',i8,2(f21.13,3x),d21.5)" - write(err_msg, frmt)n,latMesh(n),lat(n), diff_lat + frmt = "('ERROR: Difference between ESMF Mesh and MOM6 domain coords is"//& + "greater than parameter EPS_OMESH. n, latMesh(n), lat(n), diff_lat, "//& + "EPS_OMESH= ',i8,2(f21.13,3x),2(d21.5))" + write(err_msg, frmt)n,latMesh(n),lat(n), diff_lat, eps_omesh call MOM_error(FATAL, err_msg) end if if (abs(maskMesh(n) - mask(n)) > 0) then - frmt = "('ERROR: Inconsistent masks - "//& + frmt = "('ERROR: ESMF mesh and MOM6 domain masks are inconsistent! - "//& "MOM n, maskMesh(n), mask(n) = ',3(i8,2x))" write(err_msg, frmt)n,maskMesh(n),mask(n) call MOM_error(FATAL, err_msg) From b810e26be32de73a5c176dd50d91c7e502547576 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Thu, 29 Aug 2019 15:37:15 -0600 Subject: [PATCH 027/225] Fixes line length exceeded --- src/parameterizations/lateral/MOM_hor_visc.F90 | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index 8074fb2f7d..fe9073b7f1 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -1062,7 +1062,8 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, do J=Jsq,Jeq+1 ; do i=Isq,Ieq+1 if ((grad_vel_mag_bt_h(i,j)>0) .and. (max_diss_rate_h(i,j,k)>0)) then - GME_coeff = (MIN(G%bathyT(i,j)/CS%GME_h0,1.0)**2) * CS%GME_efficiency*max_diss_rate_h(i,j,k) / grad_vel_mag_bt_h(i,j) + GME_coeff = (MIN(G%bathyT(i,j)/CS%GME_h0,1.0)**2) * CS%GME_efficiency*max_diss_rate_h(i,j,k) / & + grad_vel_mag_bt_h(i,j) else GME_coeff = 1.0 endif @@ -1080,7 +1081,8 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, do J=js-1,Jeq ; do I=is-1,Ieq if ((grad_vel_mag_bt_q(I,J)>0) .and. (max_diss_rate_q(I,J,k)>0)) then - GME_coeff = (MIN(G%bathyT(i,j)/CS%GME_h0,1.0)**2) * CS%GME_efficiency*max_diss_rate_q(I,J,k) / grad_vel_mag_bt_q(I,J) + GME_coeff = (MIN(G%bathyT(i,j)/CS%GME_h0,1.0)**2) * CS%GME_efficiency*max_diss_rate_q(I,J,k) / & + grad_vel_mag_bt_q(I,J) else GME_coeff = 0.0 endif From 15a03c5aadfc8d2fd968a018dd07171c85451d05 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Tue, 3 Sep 2019 09:25:25 -0600 Subject: [PATCH 028/225] Commented out variables left over from merge These need to be cleaned before submitting a PR to dev/ncar --- .../lateral/MOM_hor_visc.F90 | 23 ++++++++++--------- 1 file changed, 12 insertions(+), 11 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index 0498e34edf..f734bb2788 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -1002,7 +1002,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if (CS%id_Kh_q>0 .or. CS%debug) Kh_q(I,J,k) = Kh if (CS%id_vort_xy_q>0) vort_xy_q(I,J,k) = vort_xy(I,J) - if (CS%debug) sh_xy_3d(I,J,k) = sh_xy(I,J) +! if (CS%debug) sh_xy_3d(I,J,k) = sh_xy(I,J) str_xy(I,J) = -Kh * sh_xy(I,J) else ! not Laplacian @@ -1048,10 +1048,10 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if (CS%id_Ah_q>0 .or. CS%debug) Ah_q(I,J,k) = Ah - str_xy(I,J) = str_xy(I,J) + Ah * ( dvdx3(I,J) + dudy3(I,J) ) + str_xy(I,J) = str_xy(I,J) + Ah * ( dDel2vdx(I,J) + dDel2udy(I,J) ) ! Keep a copy of the biharmonic contribution for backscatter parameterization - bhstr_xy(I,J) = Ah * ( dvdx3(I,J) + dudy3(I,J) ) * & + bhstr_xy(I,J) = Ah * ( dDel2vdx(I,J) + dDel2udy(I,J) ) * & (hq(I,J) * G%mask2dBu(I,J) * CS%reduction_xy(I,J)) endif ! biharmonic @@ -1305,7 +1305,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if (CS%Laplacian) then call hchksum(Kh_h, "Kh_h", G%HI, haloshift=0, scale=US%L_to_m**2*US%s_to_T) call Bchksum(Kh_q, "Kh_q", G%HI, haloshift=0, scale=US%L_to_m**2*US%s_to_T) - call Bchksum(sh_xy_3d, "shear_xy", G%HI, haloshift=0, scale=US%s_to_T) +! call Bchksum(sh_xy_3d, "shear_xy", G%HI, haloshift=0, scale=US%s_to_T) ! call hchksum(sh_xx_3d, "shear_xx", G%HI, haloshift=0, scale=US%s_to_T) endif if (CS%biharmonic) call hchksum(Ah_h, "Ah_h", G%HI, haloshift=0, scale=US%L_to_m**4*US%s_to_T) @@ -1426,13 +1426,14 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS, MEKE) ! parameter spelling checks. call get_param(param_file, mdl, "GET_ALL_PARAMS", get_all, default=.false.) - call get_param(param_file, mdl, "DEFAULT_2018_ANSWERS", default_2018_answers, & - "This sets the default value for the various _2018_ANSWERS parameters.", & - default=.true.) - call get_param(param_file, mdl, "HOR_VISC_2018_ANSWERS", CS%answers_2018, & - "If true, use the order of arithmetic and expressions that recover the "//& - "answers from the end of 2018. Otherwise, use updated and more robust "//& - "forms of the same expressions.", default=default_2018_answers) +! call get_param(param_file, mdl, "DEFAULT_2018_ANSWERS", default_2018_answers, & +! "This sets the default value for the various _2018_ANSWERS parameters.", & +! default=.true.) +! call get_param(param_file, mdl, "HOR_VISC_2018_ANSWERS", CS%answers_2018, & +! "If true, use the order of arithmetic and expressions that recover the "//& +! "answers from the end of 2018. Otherwise, use updated and more robust "//& +! "forms of the same expressions.", default=default_2018_answers) + call get_param(param_file, mdl, "DEBUG", CS%debug, default=.false.) call get_param(param_file, mdl, "LAPLACIAN", CS%Laplacian, & From ab9bc8a3eebfaeefc45478743eeb1664d0c70f9b Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Tue, 3 Sep 2019 09:27:46 -0600 Subject: [PATCH 029/225] Adds OS%US as argument in add_shelf_forces --- config_src/mct_driver/mom_ocean_model_mct.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/config_src/mct_driver/mom_ocean_model_mct.F90 b/config_src/mct_driver/mom_ocean_model_mct.F90 index 4f1c7d963a..94bc15cec8 100644 --- a/config_src/mct_driver/mom_ocean_model_mct.F90 +++ b/config_src/mct_driver/mom_ocean_model_mct.F90 @@ -525,7 +525,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, & if (do_thermo) & call shelf_calc_flux(OS%sfc_state, OS%fluxes, OS%Time, dt_coupling, OS%Ice_shelf_CSp) if (do_dyn) & - call add_shelf_forces(OS%grid, OS%Ice_shelf_CSp, OS%forces) + call add_shelf_forces(OS%grid, OS%US, OS%Ice_shelf_CSp, OS%forces) endif if (OS%icebergs_alter_ocean) then if (do_dyn) & @@ -560,7 +560,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, & if (do_thermo) & call shelf_calc_flux(OS%sfc_state, OS%flux_tmp, OS%Time, dt_coupling, OS%Ice_shelf_CSp) if (do_dyn) & - call add_shelf_forces(OS%grid, OS%Ice_shelf_CSp, OS%forces) + call add_shelf_forces(OS%grid, OS%US, OS%Ice_shelf_CSp, OS%forces) endif if (OS%icebergs_alter_ocean) then if (do_dyn) & From 25b409bdcb3d0b2bc156093f967d5232b8d0e2eb Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Tue, 3 Sep 2019 10:45:15 -0600 Subject: [PATCH 030/225] Adding back 2018_answers and cleanning the code --- .../lateral/MOM_hor_visc.F90 | 22 +++++++------------ 1 file changed, 8 insertions(+), 14 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index f734bb2788..e9ff06d337 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -306,7 +306,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, FrictWork_GME, & ! work done by GME [W m-2] div_xx_h ! horizontal divergence [s-1] real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: & - ! KH_t_GME, & !< interface height diffusivities in t-columns [m2 s-1] GME_coeff_h !< GME coeff. at h-points [L2 T-1 ~> m2 s-1] real :: Ah ! biharmonic viscosity [L4 T-1 ~> m4 s-1] real :: Kh ! Laplacian viscosity [L2 T-1 ~> m2 s-1] @@ -363,11 +362,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, Ah_h(:,:,:) = 0.0 Kh_h(:,:,:) = 0.0 -! if (CS%debug) then -! sh_xx_3d(:,:,:) = 0.0 ; sh_xy_3d(:,:,:) = 0.0 -! Kh_q(:,:,:) = 0.0 ; Ah_q(:,:,:) = 0.0 -! endif - if (present(OBC)) then ; if (associated(OBC)) then ; if (OBC%OBC_pe) then apply_OBC = OBC%Flather_u_BCs_exist_globally .or. OBC%Flather_v_BCs_exist_globally apply_OBC = .true. @@ -1379,6 +1373,7 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS, MEKE) ! If false and USE_GME = True, issue a FATAL error. logical :: use_MEKE ! If true, use the MEKE module for calculating eddy kinetic energy. ! If false and USE_GME = True, issue a FATAL error. + logical :: default_2018_answers character(len=64) :: inputdir, filename real :: deg2rad ! Converts degrees to radians @@ -1426,13 +1421,13 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS, MEKE) ! parameter spelling checks. call get_param(param_file, mdl, "GET_ALL_PARAMS", get_all, default=.false.) -! call get_param(param_file, mdl, "DEFAULT_2018_ANSWERS", default_2018_answers, & -! "This sets the default value for the various _2018_ANSWERS parameters.", & -! default=.true.) -! call get_param(param_file, mdl, "HOR_VISC_2018_ANSWERS", CS%answers_2018, & -! "If true, use the order of arithmetic and expressions that recover the "//& -! "answers from the end of 2018. Otherwise, use updated and more robust "//& -! "forms of the same expressions.", default=default_2018_answers) + call get_param(param_file, mdl, "DEFAULT_2018_ANSWERS", default_2018_answers, & + "This sets the default value for the various _2018_ANSWERS parameters.", & + default=.true.) + call get_param(param_file, mdl, "HOR_VISC_2018_ANSWERS", CS%answers_2018, & + "If true, use the order of arithmetic and expressions that recover the "//& + "answers from the end of 2018. Otherwise, use updated and more robust "//& + "forms of the same expressions.", default=default_2018_answers) call get_param(param_file, mdl, "DEBUG", CS%debug, default=.false.) @@ -2142,7 +2137,6 @@ subroutine smooth_GME(CS,G,GME_flux_h,GME_flux_q) real :: wc, ww, we, wn, ws ! averaging weights for smoothing integer :: i, j, k, s - !do s=1,CS%n_smooth do s=1,1 ! Update halos From 1016be3cba6b07a88764c4aad95ef2253ef4d476 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Wed, 11 Sep 2019 15:04:21 -0400 Subject: [PATCH 031/225] +Pass timestep arguments to vertvisc in [T] Rescaled the timestep arguments to set_viscous_ML, vertvisc, vertvisc_coef, vertvisc_remnant, write_u_accel and write_v_accel to use units of [T]. All answers are bitwise identical, but the units of some public arguments have been rescaled. --- src/core/MOM_dynamics_split_RK2.F90 | 20 +++---- src/core/MOM_dynamics_unsplit.F90 | 17 +++--- src/core/MOM_dynamics_unsplit_RK2.F90 | 17 +++--- src/diagnostics/MOM_PointAccel.F90 | 12 ++-- .../vertical/MOM_set_viscosity.F90 | 10 ++-- .../vertical/MOM_vert_friction.F90 | 58 +++++++++---------- 6 files changed, 66 insertions(+), 68 deletions(-) diff --git a/src/core/MOM_dynamics_split_RK2.F90 b/src/core/MOM_dynamics_split_RK2.F90 index 1f43a699a1..839dcc9f24 100644 --- a/src/core/MOM_dynamics_split_RK2.F90 +++ b/src/core/MOM_dynamics_split_RK2.F90 @@ -56,7 +56,7 @@ module MOM_dynamics_split_RK2 use MOM_tidal_forcing, only : tidal_forcing_init, tidal_forcing_CS use MOM_unit_scaling, only : unit_scale_type use MOM_vert_friction, only : vertvisc, vertvisc_coef, vertvisc_remnant -use MOM_vert_friction, only : vertvisc_limit_vel, vertvisc_init, vertvisc_CS +use MOM_vert_friction, only : vertvisc_init, vertvisc_CS use MOM_vert_friction, only : updateCFLtruncationValue use MOM_verticalGrid, only : verticalGrid_type, get_thickness_units use MOM_verticalGrid, only : get_flux_units, get_tr_flux_units @@ -480,15 +480,15 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, & enddo call enable_averaging(dt, Time_local, CS%diag) - call set_viscous_ML(u, v, h, tv, forces, visc, dt, G, GV, US, & + call set_viscous_ML(u, v, h, tv, forces, visc, dt_in_T, G, GV, US, & CS%set_visc_CSp) call disable_averaging(CS%diag) if (CS%debug) then call uvchksum("before vertvisc: up", up, vp, G%HI, haloshift=0, symmetric=sym) endif - call vertvisc_coef(up, vp, h, forces, visc, dt, G, GV, US, CS%vertvisc_CSp, CS%OBC) - call vertvisc_remnant(visc, CS%visc_rem_u, CS%visc_rem_v, dt, G, GV, US, CS%vertvisc_CSp) + call vertvisc_coef(up, vp, h, forces, visc, dt_in_T, G, GV, US, CS%vertvisc_CSp, CS%OBC) + call vertvisc_remnant(visc, CS%visc_rem_u, CS%visc_rem_v, dt_in_T, G, GV, US, CS%vertvisc_CSp) call cpu_clock_end(id_clock_vertvisc) if (showCallTree) call callTree_wayPoint("done with vertvisc_coef (step_MOM_dyn_split_RK2)") @@ -580,9 +580,9 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, & if (CS%debug) then call uvchksum("0 before vertvisc: [uv]p", up, vp, G%HI,haloshift=0, symmetric=sym, scale=US%L_T_to_m_s) endif - call vertvisc_coef(up, vp, h, forces, visc, US%T_to_s*dt_pred, G, GV, US, CS%vertvisc_CSp, & + call vertvisc_coef(up, vp, h, forces, visc, dt_pred, G, GV, US, CS%vertvisc_CSp, & CS%OBC) - call vertvisc(up, vp, h, forces, visc, US%T_to_s*dt_pred, CS%OBC, CS%ADp, CS%CDp, G, & + call vertvisc(up, vp, h, forces, visc, dt_pred, CS%OBC, CS%ADp, CS%CDp, G, & GV, US, CS%vertvisc_CSp, CS%taux_bot, CS%tauy_bot, waves=waves) if (showCallTree) call callTree_wayPoint("done with vertvisc (step_MOM_dyn_split_RK2)") if (G%nonblocking_updates) then @@ -590,7 +590,7 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, & call start_group_pass(CS%pass_uvp, G%Domain, clock=id_clock_pass) call cpu_clock_begin(id_clock_vertvisc) endif - call vertvisc_remnant(visc, CS%visc_rem_u, CS%visc_rem_v, US%T_to_s*dt_pred, G, GV, US, CS%vertvisc_CSp) + call vertvisc_remnant(visc, CS%visc_rem_u, CS%visc_rem_v, dt_pred, G, GV, US, CS%vertvisc_CSp) call cpu_clock_end(id_clock_vertvisc) call do_group_pass(CS%pass_visc_rem, G%Domain, clock=id_clock_pass) @@ -779,15 +779,15 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, & ! u <- u + dt d/dz visc d/dz u ! u_av <- u_av + dt d/dz visc d/dz u_av call cpu_clock_begin(id_clock_vertvisc) - call vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS%vertvisc_CSp, CS%OBC) - call vertvisc(u, v, h, forces, visc, dt, CS%OBC, CS%ADp, CS%CDp, G, GV, US, & + call vertvisc_coef(u, v, h, forces, visc, dt_in_T, G, GV, US, CS%vertvisc_CSp, CS%OBC) + call vertvisc(u, v, h, forces, visc, dt_in_T, CS%OBC, CS%ADp, CS%CDp, G, GV, US, & CS%vertvisc_CSp, CS%taux_bot, CS%tauy_bot,waves=waves) if (G%nonblocking_updates) then call cpu_clock_end(id_clock_vertvisc) call start_group_pass(CS%pass_uv, G%Domain, clock=id_clock_pass) call cpu_clock_begin(id_clock_vertvisc) endif - call vertvisc_remnant(visc, CS%visc_rem_u, CS%visc_rem_v, dt, G, GV, US, CS%vertvisc_CSp) + call vertvisc_remnant(visc, CS%visc_rem_u, CS%visc_rem_v, dt_in_T, G, GV, US, CS%vertvisc_CSp) call cpu_clock_end(id_clock_vertvisc) if (showCallTree) call callTree_wayPoint("done with vertvisc (step_MOM_dyn_split_RK2)") diff --git a/src/core/MOM_dynamics_unsplit.F90 b/src/core/MOM_dynamics_unsplit.F90 index 9e8be65d7a..58d04cff5a 100644 --- a/src/core/MOM_dynamics_unsplit.F90 +++ b/src/core/MOM_dynamics_unsplit.F90 @@ -94,8 +94,7 @@ module MOM_dynamics_unsplit use MOM_thickness_diffuse, only : thickness_diffuse_CS use MOM_tidal_forcing, only : tidal_forcing_init, tidal_forcing_CS use MOM_unit_scaling, only : unit_scale_type -use MOM_vert_friction, only : vertvisc, vertvisc_coef -use MOM_vert_friction, only : vertvisc_limit_vel, vertvisc_init, vertvisc_CS +use MOM_vert_friction, only : vertvisc, vertvisc_coef, vertvisc_init, vertvisc_CS use MOM_verticalGrid, only : verticalGrid_type, get_thickness_units use MOM_verticalGrid, only : get_flux_units, get_tr_flux_units use MOM_wave_interface, only: wave_parameters_CS @@ -344,13 +343,13 @@ subroutine step_MOM_dyn_unsplit(u, v, h, tv, visc, Time_local, dt, forces, & ! up <- up + dt/2 d/dz visc d/dz up call cpu_clock_begin(id_clock_vertvisc) call enable_averaging(dt, Time_local, CS%diag) - call set_viscous_ML(u, v, h_av, tv, forces, visc, dt*0.5, G, GV, US, & + call set_viscous_ML(u, v, h_av, tv, forces, visc, dt_in_T*0.5, G, GV, US, & CS%set_visc_CSp) call disable_averaging(CS%diag) !### I think that the time steps in the next two calls should be dt_pred. - call vertvisc_coef(up, vp, h_av, forces, visc, dt*0.5, G, GV, US, & + call vertvisc_coef(up, vp, h_av, forces, visc, dt_in_T*0.5, G, GV, US, & CS%vertvisc_CSp, CS%OBC) - call vertvisc(up, vp, h_av, forces, visc, dt*0.5, CS%OBC, CS%ADp, CS%CDp, & + call vertvisc(up, vp, h_av, forces, visc, dt_in_T*0.5, CS%OBC, CS%ADp, CS%CDp, & G, GV, US, CS%vertvisc_CSp, Waves=Waves) call cpu_clock_end(id_clock_vertvisc) call pass_vector(up, vp, G%Domain, clock=id_clock_pass) @@ -412,9 +411,9 @@ subroutine step_MOM_dyn_unsplit(u, v, h, tv, visc, Time_local, dt, forces, & ! upp <- upp + dt/2 d/dz visc d/dz upp call cpu_clock_begin(id_clock_vertvisc) - call vertvisc_coef(upp, vpp, hp, forces, visc, dt*0.5, G, GV, US, & + call vertvisc_coef(upp, vpp, hp, forces, visc, dt_in_T*0.5, G, GV, US, & CS%vertvisc_CSp, CS%OBC) - call vertvisc(upp, vpp, hp, forces, visc, dt*0.5, CS%OBC, CS%ADp, CS%CDp, & + call vertvisc(upp, vpp, hp, forces, visc, dt_in_T*0.5, CS%OBC, CS%ADp, CS%CDp, & G, GV, US, CS%vertvisc_CSp, Waves=Waves) call cpu_clock_end(id_clock_vertvisc) call pass_vector(upp, vpp, G%Domain, clock=id_clock_pass) @@ -483,8 +482,8 @@ subroutine step_MOM_dyn_unsplit(u, v, h, tv, visc, Time_local, dt, forces, & ! u <- u + dt d/dz visc d/dz u call cpu_clock_begin(id_clock_vertvisc) - call vertvisc_coef(u, v, h_av, forces, visc, dt, G, GV, US, CS%vertvisc_CSp, CS%OBC) - call vertvisc(u, v, h_av, forces, visc, dt, CS%OBC, CS%ADp, CS%CDp, & + call vertvisc_coef(u, v, h_av, forces, visc, dt_in_T, G, GV, US, CS%vertvisc_CSp, CS%OBC) + call vertvisc(u, v, h_av, forces, visc, dt_in_T, CS%OBC, CS%ADp, CS%CDp, & G, GV, US, CS%vertvisc_CSp, CS%taux_bot, CS%tauy_bot, Waves=Waves) call cpu_clock_end(id_clock_vertvisc) call pass_vector(u, v, G%Domain, clock=id_clock_pass) diff --git a/src/core/MOM_dynamics_unsplit_RK2.F90 b/src/core/MOM_dynamics_unsplit_RK2.F90 index af33db8011..97ef3ede73 100644 --- a/src/core/MOM_dynamics_unsplit_RK2.F90 +++ b/src/core/MOM_dynamics_unsplit_RK2.F90 @@ -92,8 +92,7 @@ module MOM_dynamics_unsplit_RK2 use MOM_thickness_diffuse, only : thickness_diffuse_CS use MOM_tidal_forcing, only : tidal_forcing_init, tidal_forcing_CS use MOM_unit_scaling, only : unit_scale_type -use MOM_vert_friction, only : vertvisc, vertvisc_coef -use MOM_vert_friction, only : vertvisc_limit_vel, vertvisc_init, vertvisc_CS +use MOM_vert_friction, only : vertvisc, vertvisc_coef, vertvisc_init, vertvisc_CS use MOM_verticalGrid, only : verticalGrid_type, get_thickness_units use MOM_verticalGrid, only : get_flux_units, get_tr_flux_units @@ -342,12 +341,12 @@ subroutine step_MOM_dyn_unsplit_RK2(u_in, v_in, h_in, tv, visc, Time_local, dt, ! up[n-1/2] <- up*[n-1/2] + dt/2 d/dz visc d/dz up[n-1/2] call cpu_clock_begin(id_clock_vertvisc) call enable_averaging(dt, Time_local, CS%diag) - call set_viscous_ML(up, vp, h_av, tv, forces, visc, US%T_to_s*dt_pred, G, GV, US, & + call set_viscous_ML(up, vp, h_av, tv, forces, visc, dt_pred, G, GV, US, & CS%set_visc_CSp) call disable_averaging(CS%diag) - call vertvisc_coef(up, vp, h_av, forces, visc, US%T_to_s*dt_pred, G, GV, US, & + call vertvisc_coef(up, vp, h_av, forces, visc, dt_pred, G, GV, US, & CS%vertvisc_CSp, CS%OBC) - call vertvisc(up, vp, h_av, forces, visc, US%T_to_s*dt_pred, CS%OBC, CS%ADp, CS%CDp, & + call vertvisc(up, vp, h_av, forces, visc, dt_pred, CS%OBC, CS%ADp, CS%CDp, & G, GV, US, CS%vertvisc_CSp) call cpu_clock_end(id_clock_vertvisc) call pass_vector(up, vp, G%Domain, clock=id_clock_pass) @@ -397,13 +396,13 @@ subroutine step_MOM_dyn_unsplit_RK2(u_in, v_in, h_in, tv, visc, Time_local, dt, ! up[n] <- up* + dt d/dz visc d/dz up ! u[n] <- u*[n] + dt d/dz visc d/dz u[n] call cpu_clock_begin(id_clock_vertvisc) - call vertvisc_coef(up, vp, h_av, forces, visc, dt, G, GV, US, & + call vertvisc_coef(up, vp, h_av, forces, visc, dt_in_T, G, GV, US, & CS%vertvisc_CSp, CS%OBC) - call vertvisc(up, vp, h_av, forces, visc, dt, CS%OBC, CS%ADp, CS%CDp, & + call vertvisc(up, vp, h_av, forces, visc, dt_in_T, CS%OBC, CS%ADp, CS%CDp, & G, GV, US, CS%vertvisc_CSp, CS%taux_bot, CS%tauy_bot) - call vertvisc_coef(u_in, v_in, h_av, forces, visc, dt, G, GV, US, & + call vertvisc_coef(u_in, v_in, h_av, forces, visc, dt_in_T, G, GV, US, & CS%vertvisc_CSp, CS%OBC) - call vertvisc(u_in, v_in, h_av, forces, visc, dt, CS%OBC, CS%ADp, CS%CDp,& + call vertvisc(u_in, v_in, h_av, forces, visc, dt_in_T, CS%OBC, CS%ADp, CS%CDp,& G, GV, US, CS%vertvisc_CSp, CS%taux_bot, CS%tauy_bot) call cpu_clock_end(id_clock_vertvisc) call pass_vector(up, vp, G%Domain, clock=id_clock_pass) diff --git a/src/diagnostics/MOM_PointAccel.F90 b/src/diagnostics/MOM_PointAccel.F90 index e0bbd832bb..dd72378671 100644 --- a/src/diagnostics/MOM_PointAccel.F90 +++ b/src/diagnostics/MOM_PointAccel.F90 @@ -66,7 +66,7 @@ module MOM_PointAccel !> This subroutine writes to an output file all of the accelerations !! that have been applied to a column of zonal velocities over the !! previous timestep. This subroutine is called from vertvisc. -subroutine write_u_accel(I, j, um, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, str, a, hv) +subroutine write_u_accel(I, j, um, hin, ADp, CDp, dt_in_T, G, GV, US, CS, vel_rpt, str, a, hv) integer, intent(in) :: I !< The zonal index of the column to be documented. integer, intent(in) :: j !< The meridional index of the column to be documented. type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure. @@ -80,7 +80,7 @@ subroutine write_u_accel(I, j, um, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st !! accelerations in the momentum equations. type(cont_diag_ptrs), intent(in) :: CDp !< A structure with pointers to various terms !! in the continuity equations. - real, intent(in) :: dt !< The ocean dynamics time step [s]. + real, intent(in) :: dt_in_T !< The ocean dynamics time step [T ~> s]. type(PointAccel_CS), pointer :: CS !< The control structure returned by a previous !! call to PointAccel_init. real, intent(in) :: vel_rpt !< The velocity magnitude that triggers a report [L T-1 ~> m s-1]. @@ -95,6 +95,7 @@ subroutine write_u_accel(I, j, um, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st real :: f_eff, CFL real :: Angstrom real :: truncvel, du + real :: dt ! The time step [s] real :: Inorm(SZK_(G)) real :: e(SZK_(G)+1) real :: h_scale, uh_scale @@ -106,6 +107,7 @@ subroutine write_u_accel(I, j, um, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st integer :: file Angstrom = GV%Angstrom_H + GV%H_subroundoff + dt = US%T_to_s*dt_in_T h_scale = GV%H_to_m ; uh_scale = GV%H_to_m ! if (.not.associated(CS)) return @@ -397,7 +399,7 @@ end subroutine write_u_accel !> This subroutine writes to an output file all of the accelerations !! that have been applied to a column of meridional velocities over !! the previous timestep. This subroutine is called from vertvisc. -subroutine write_v_accel(i, J, vm, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, str, a, hv) +subroutine write_v_accel(i, J, vm, hin, ADp, CDp, dt_in_T, G, GV, US, CS, vel_rpt, str, a, hv) integer, intent(in) :: i !< The zonal index of the column to be documented. integer, intent(in) :: J !< The meridional index of the column to be documented. type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure. @@ -411,7 +413,7 @@ subroutine write_v_accel(i, J, vm, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st !! accelerations in the momentum equations. type(cont_diag_ptrs), intent(in) :: CDp !< A structure with pointers to various terms in !! the continuity equations. - real, intent(in) :: dt !< The ocean dynamics time step [s]. + real, intent(in) :: dt_in_T !< The ocean dynamics time step [T ~> s]. type(PointAccel_CS), pointer :: CS !< The control structure returned by a previous !! call to PointAccel_init. real, intent(in) :: vel_rpt !< The velocity magnitude that triggers a report [L T-1 ~> m s-1]. @@ -426,6 +428,7 @@ subroutine write_v_accel(i, J, vm, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st real :: f_eff, CFL real :: Angstrom real :: truncvel, dv + real :: dt ! The time step [s] real :: Inorm(SZK_(G)) real :: e(SZK_(G)+1) real :: h_scale, uh_scale @@ -437,6 +440,7 @@ subroutine write_v_accel(i, J, vm, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st integer :: file Angstrom = GV%Angstrom_H + GV%H_subroundoff + dt = US%T_to_s*dt_in_T h_scale = GV%H_to_m ; uh_scale = GV%H_to_m ! if (.not.associated(CS)) return diff --git a/src/parameterizations/vertical/MOM_set_viscosity.F90 b/src/parameterizations/vertical/MOM_set_viscosity.F90 index 92466266b8..30648c7d61 100644 --- a/src/parameterizations/vertical/MOM_set_viscosity.F90 +++ b/src/parameterizations/vertical/MOM_set_viscosity.F90 @@ -273,7 +273,7 @@ subroutine set_viscous_BBL(u, v, h, tv, visc, G, GV, US, CS, symmetrize) Vol_quit = 0.9*GV%Angstrom_H + h_neglect C2pi_3 = 8.0*atan(1.0)/3.0 - if (.not.associated(CS)) call MOM_error(FATAL,"MOM_vert_friction(BBL): "//& + if (.not.associated(CS)) call MOM_error(FATAL,"MOM_set_viscosity(BBL): "//& "Module must be initialized before it is used.") if (.not.CS%bottomdraglaw) return @@ -1002,7 +1002,7 @@ end function set_u_at_v !! A bulk Richardson criterion or the thickness of the topmost NKML layers (with a bulk mixed layer) !! are currently used. The thicknesses are given in terms of fractional layers, so that this !! thickness will move as the thickness of the topmost layers change. -subroutine set_viscous_ML(u, v, h, tv, forces, visc, dt, G, GV, US, CS, symmetrize) +subroutine set_viscous_ML(u, v, h, tv, forces, visc, dt_in_T, G, GV, US, CS, symmetrize) type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure. type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure. type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type @@ -1018,7 +1018,7 @@ subroutine set_viscous_ML(u, v, h, tv, forces, visc, dt, G, GV, US, CS, symmetri type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces type(vertvisc_type), intent(inout) :: visc !< A structure containing vertical viscosities and !! related fields. - real, intent(in) :: dt !< Time increment [s]. + real, intent(in) :: dt_in_T !< Time increment [T ~> s]. type(set_visc_CS), pointer :: CS !< The control structure returned by a previous !! call to vertvisc_init. logical, optional, intent(in) :: symmetrize !< If present and true, do extra calculations @@ -1125,7 +1125,7 @@ subroutine set_viscous_ML(u, v, h, tv, forces, visc, dt, G, GV, US, CS, symmetri Isq = G%IscB ; Ieq = G%IecB ; Jsq = G%JscB ; Jeq = G%JecB nkmb = GV%nk_rho_varies ; nkml = GV%nkml - if (.not.associated(CS)) call MOM_error(FATAL,"MOM_vert_friction(visc_ML): "//& + if (.not.associated(CS)) call MOM_error(FATAL,"MOM_set_viscosity(visc_ML): "//& "Module must be initialized before it is used.") if (.not.(CS%dynamic_viscous_ML .or. associated(forces%frac_shelf_u) .or. & associated(forces%frac_shelf_v)) ) return @@ -1141,7 +1141,7 @@ subroutine set_viscous_ML(u, v, h, tv, forces, visc, dt, G, GV, US, CS, symmetri OBC => CS%OBC use_EOS = associated(tv%eqn_of_state) - dt_Rho0 = dt/GV%H_to_kg_m2 + dt_Rho0 = US%T_to_s*dt_in_T / GV%H_to_kg_m2 h_neglect = GV%H_subroundoff h_tiny = 2.0*GV%Angstrom_H + h_neglect g_H_Rho0 = (GV%g_Earth*GV%H_to_Z) / GV%Rho0 diff --git a/src/parameterizations/vertical/MOM_vert_friction.F90 b/src/parameterizations/vertical/MOM_vert_friction.F90 index b282995d3f..fe14380617 100644 --- a/src/parameterizations/vertical/MOM_vert_friction.F90 +++ b/src/parameterizations/vertical/MOM_vert_friction.F90 @@ -142,7 +142,7 @@ module MOM_vert_friction !! There is an additional stress term on the right-hand side !! if DIRECT_STRESS is true, applied to the surface layer. -subroutine vertvisc(u, v, h, forces, visc, dt, OBC, ADp, CDp, G, GV, US, CS, & +subroutine vertvisc(u, v, h, forces, visc, dt_in_T, OBC, ADp, CDp, G, GV, US, CS, & taux_bot, tauy_bot, Waves) type(ocean_grid_type), intent(in) :: G !< Ocean grid structure type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure @@ -155,7 +155,7 @@ subroutine vertvisc(u, v, h, forces, visc, dt, OBC, ADp, CDp, G, GV, US, CS, & intent(in) :: h !< Layer thickness [H ~> m or kg m-2] type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces type(vertvisc_type), intent(inout) :: visc !< Viscosities and bottom drag - real, intent(in) :: dt !< Time increment [s] + real, intent(in) :: dt_in_T !< Time increment [T ~> s] type(ocean_OBC_type), pointer :: OBC !< Open boundary condition structure type(accel_diag_ptrs), intent(inout) :: ADp !< Accelerations in the momentum !! equations for diagnostics @@ -185,7 +185,6 @@ subroutine vertvisc(u, v, h, forces, visc, dt, OBC, ADp, CDp, G, GV, US, CS, & real :: Hmix ! The mixed layer thickness over which stress ! is applied with direct_stress [H ~> m or kg m-2]. real :: I_Hmix ! The inverse of Hmix [H-1 ~> m-1 or m2 kg-1]. - real :: dt_in_T ! The timestep [T ~> s] real :: Idt ! The inverse of the time step [T-1 ~> s-1]. real :: dt_Rho0 ! The time step divided by the mean density [L s2 H m T-1 kg-1 ~> s m3 kg-1 or s]. real :: Rho0 ! A density used to convert drag laws into stress in Pa [kg m-3]. @@ -214,7 +213,6 @@ subroutine vertvisc(u, v, h, forces, visc, dt, OBC, ADp, CDp, G, GV, US, CS, & Hmix = CS%Hmix_stress I_Hmix = 1.0 / Hmix endif - dt_in_T = US%s_to_T*dt dt_Rho0 = US%m_s_to_L_T*US%T_to_s * dt_in_T / GV%H_to_kg_m2 dt_Z_to_H = dt_in_T*GV%Z_to_H Rho0 = GV%Rho0 @@ -422,7 +420,7 @@ subroutine vertvisc(u, v, h, forces, visc, dt, OBC, ADp, CDp, G, GV, US, CS, & enddo ! end of v-component J loop - call vertvisc_limit_vel(u, v, h, ADp, CDp, forces, visc, dt, G, GV, US, CS) + call vertvisc_limit_vel(u, v, h, ADp, CDp, forces, visc, dt_in_T, G, GV, US, CS) ! Here the velocities associated with open boundary conditions are applied. if (associated(OBC)) then @@ -459,7 +457,7 @@ end subroutine vertvisc !! after a time-step of viscosity, and the fraction of a time-step's !! worth of barotropic acceleration that a layer experiences after !! viscosity is applied. -subroutine vertvisc_remnant(visc, visc_rem_u, visc_rem_v, dt, G, GV, US, CS) +subroutine vertvisc_remnant(visc, visc_rem_u, visc_rem_v, dt_in_T, G, GV, US, CS) type(ocean_grid_type), intent(in) :: G !< Ocean grid structure type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure type(vertvisc_type), intent(in) :: visc !< Viscosities and bottom drag @@ -471,7 +469,7 @@ subroutine vertvisc_remnant(visc, visc_rem_u, visc_rem_v, dt, G, GV, US, CS) intent(inout) :: visc_rem_v !< Fraction of a time-step's worth of a !! barotopic acceleration that a layer experiences after !! viscosity is applied in the meridional direction [nondim] - real, intent(in) :: dt !< Time increment [s] + real, intent(in) :: dt_in_T !< Time increment [T ~> s] type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(vertvisc_CS), pointer :: CS !< Vertical viscosity control structure @@ -493,7 +491,7 @@ subroutine vertvisc_remnant(visc, visc_rem_u, visc_rem_v, dt, G, GV, US, CS) if (.not.associated(CS)) call MOM_error(FATAL,"MOM_vert_friction(visc): "// & "Module must be initialized before it is used.") - dt_Z_to_H = US%s_to_T*dt*GV%Z_to_H + dt_Z_to_H = dt_in_T*GV%Z_to_H do k=1,nz ; do i=Isq,Ieq ; Ray(i,k) = 0.0 ; enddo ; enddo @@ -567,7 +565,7 @@ end subroutine vertvisc_remnant !> Calculate the coupling coefficients (CS%a_u and CS%a_v) !! and effective layer thicknesses (CS%h_u and CS%h_v) for later use in the !! applying the implicit vertical viscosity via vertvisc(). -subroutine vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS, OBC) +subroutine vertvisc_coef(u, v, h, forces, visc, dt_in_T, G, GV, US, CS, OBC) type(ocean_grid_type), intent(in) :: G !< Ocean grid structure type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type @@ -579,7 +577,7 @@ subroutine vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS, OBC) intent(in) :: h !< Layer thickness [H ~> m or kg m-2] type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces type(vertvisc_type), intent(in) :: visc !< Viscosities and bottom drag - real, intent(in) :: dt !< Time increment [s] + real, intent(in) :: dt_in_T !< Time increment [T ~> s] type(vertvisc_CS), pointer :: CS !< Vertical viscosity control structure type(ocean_OBC_type), pointer :: OBC !< Open boundary condition structure @@ -758,7 +756,7 @@ subroutine vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS, OBC) endif call find_coupling_coef(a_cpl, hvel, do_i, h_harm, bbl_thick, kv_bbl, z_i, h_ml, & - dt, j, G, GV, US, CS, visc, forces, work_on_u=.true., OBC=OBC) + dt_in_T, j, G, GV, US, CS, visc, forces, work_on_u=.true., OBC=OBC) if (allocated(hML_u)) then do i=isq,ieq ; if (do_i(i)) then ; hML_u(I,j) = h_ml(I) ; endif ; enddo endif @@ -773,7 +771,7 @@ subroutine vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS, OBC) if (do_any_shelf) then if (CS%harmonic_visc) then call find_coupling_coef(a_shelf, hvel, do_i_shelf, h_harm, bbl_thick, & - kv_bbl, z_i, h_ml, dt, j, G, GV, US, CS, & + kv_bbl, z_i, h_ml, dt_in_T, j, G, GV, US, CS, & visc, forces, work_on_u=.true., OBC=OBC, shelf=.true.) else ! Find upwind-biased thickness near the surface. ! Perhaps this needs to be done more carefully, via find_eta. @@ -801,7 +799,7 @@ subroutine vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS, OBC) endif ; enddo enddo call find_coupling_coef(a_shelf, hvel_shelf, do_i_shelf, h_harm, & - bbl_thick, kv_bbl, z_i, h_ml, dt, j, G, GV, US, CS, & + bbl_thick, kv_bbl, z_i, h_ml, dt_in_T, j, G, GV, US, CS, & visc, forces, work_on_u=.true., OBC=OBC, shelf=.true.) endif do I=Isq,Ieq ; if (do_i_shelf(I)) CS%a1_shelf_u(I,j) = a_shelf(I,1) ; enddo @@ -927,7 +925,7 @@ subroutine vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS, OBC) endif call find_coupling_coef(a_cpl, hvel, do_i, h_harm, bbl_thick, kv_bbl, z_i, h_ml, & - dt, j, G, GV, US, CS, visc, forces, work_on_u=.false., OBC=OBC) + dt_in_T, j, G, GV, US, CS, visc, forces, work_on_u=.false., OBC=OBC) if ( allocated(hML_v)) then do i=is,ie ; if (do_i(i)) then ; hML_v(i,J) = h_ml(i) ; endif ; enddo endif @@ -941,7 +939,7 @@ subroutine vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS, OBC) if (do_any_shelf) then if (CS%harmonic_visc) then call find_coupling_coef(a_shelf, hvel, do_i_shelf, h_harm, bbl_thick, & - kv_bbl, z_i, h_ml, dt, j, G, GV, US, CS, visc, & + kv_bbl, z_i, h_ml, dt_in_T, j, G, GV, US, CS, visc, & forces, work_on_u=.false., OBC=OBC, shelf=.true.) else ! Find upwind-biased thickness near the surface. ! Perhaps this needs to be done more carefully, via find_eta. @@ -969,7 +967,7 @@ subroutine vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS, OBC) endif ; enddo enddo call find_coupling_coef(a_shelf, hvel_shelf, do_i_shelf, h_harm, & - bbl_thick, kv_bbl, z_i, h_ml, dt, j, G, GV, US, CS, & + bbl_thick, kv_bbl, z_i, h_ml, dt_in_T, j, G, GV, US, CS, & visc, forces, work_on_u=.false., OBC=OBC, shelf=.true.) endif do i=is,ie ; if (do_i_shelf(i)) CS%a1_shelf_v(i,J) = a_shelf(i,1) ; enddo @@ -1034,7 +1032,7 @@ end subroutine vertvisc_coef !! If BOTTOMDRAGLAW is defined, the minimum of Hbbl and half the adjacent !! layer thicknesses are used to calculate a_cpl near the bottom. subroutine find_coupling_coef(a_cpl, hvel, do_i, h_harm, bbl_thick, kv_bbl, z_i, h_ml, & - dt, j, G, GV, US, CS, visc, forces, work_on_u, OBC, shelf) + dt_in_T, j, G, GV, US, CS, visc, forces, work_on_u, OBC, shelf) type(ocean_grid_type), intent(in) :: G !< Ocean grid structure type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type @@ -1054,7 +1052,7 @@ subroutine find_coupling_coef(a_cpl, hvel, do_i, h_harm, bbl_thick, kv_bbl, z_i, !! normalized by the bottom boundary layer thickness real, dimension(SZIB_(G)), intent(out) :: h_ml !< Mixed layer depth [H ~> m or kg m-2] integer, intent(in) :: j !< j-index to find coupling coefficient for - real, intent(in) :: dt !< Time increment [s] + real, intent(in) :: dt_in_T !< Time increment [T ~> s] type(vertvisc_CS), pointer :: CS !< Vertical viscosity control structure type(vertvisc_type), intent(in) :: visc !< Structure containing viscosities and bottom drag type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces @@ -1107,7 +1105,7 @@ subroutine find_coupling_coef(a_cpl, hvel, do_i, h_harm, bbl_thick, kv_bbl, z_i, ! The maximum coupling coefficent was originally introduced to avoid ! truncation error problems in the tridiagonal solver. Effectively, the 1e-10 ! sets the maximum coupling coefficient increment to 1e10 m per timestep. - I_amax = (1.0e-10*US%Z_to_m) * dt*US%s_to_T + I_amax = (1.0e-10*US%Z_to_m) * dt_in_T do_shelf = .false. ; if (present(shelf)) do_shelf = shelf do_OBCs = .false. @@ -1304,10 +1302,10 @@ subroutine find_coupling_coef(a_cpl, hvel, do_i, h_harm, bbl_thick, kv_bbl, z_i, end subroutine find_coupling_coef -!> Velocity components which exceed a threshold for physically -!! reasonable values are truncated. Optionally, any column with excessive -!! velocities may be sent to a diagnostic reporting subroutine. -subroutine vertvisc_limit_vel(u, v, h, ADp, CDp, forces, visc, dt, G, GV, US, CS) +!> Velocity components which exceed a threshold for physically reasonable values +!! are truncated. Optionally, any column with excessive velocities may be sent +!! to a diagnostic reporting subroutine. +subroutine vertvisc_limit_vel(u, v, h, ADp, CDp, forces, visc, dt_in_T, G, GV, US, CS) type(ocean_grid_type), intent(in) :: G !< Ocean grid structure type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type @@ -1321,14 +1319,13 @@ subroutine vertvisc_limit_vel(u, v, h, ADp, CDp, forces, visc, dt, G, GV, US, CS type(cont_diag_ptrs), intent(in) :: CDp !< Continuity diagnostic pointers type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces type(vertvisc_type), intent(in) :: visc !< Viscosities and bottom drag - real, intent(in) :: dt !< Time increment [s] + real, intent(in) :: dt_in_T !< Time increment [T ~> s] type(vertvisc_CS), pointer :: CS !< Vertical viscosity control structure ! Local variables real :: maxvel ! Velocities components greater than maxvel real :: truncvel ! are truncated to truncvel, both [L T-1 ~> m s-1]. - real :: dt_in_T ! The timestep [T ~> s] real :: CFL ! The local CFL number. real :: H_report ! A thickness below which not to report truncations. real :: dt_Rho0 ! The timestep divided by the Boussinesq density [s m3 kg-1]. @@ -1343,8 +1340,7 @@ subroutine vertvisc_limit_vel(u, v, h, ADp, CDp, forces, visc, dt, G, GV, US, CS maxvel = CS%maxvel truncvel = 0.9*maxvel H_report = 6.0 * GV%Angstrom_H - dt_in_T = US%s_to_T*dt - dt_Rho0 = dt / GV%Rho0 + dt_Rho0 = US%T_to_s*dt_in_T / GV%Rho0 if (len_trim(CS%u_trunc_file) > 0) then !$OMP parallel do default(shared) private(trunc_any,CFL) @@ -1415,7 +1411,7 @@ subroutine vertvisc_limit_vel(u, v, h, ADp, CDp, forces, visc, dt, G, GV, US, CS do k=1,nz ; do j=js,je ; do I=Isq,Ieq if (abs(u(I,j,k)) < CS%vel_underflow) then ; u(I,j,k) = 0.0 elseif (abs(u(I,j,k)) > maxvel) then - u(I,j,k) = SIGN(truncvel,u(I,j,k)) + u(I,j,k) = SIGN(truncvel, u(I,j,k)) if (h(i,j,k) + h(i+1,j,k) > H_report) CS%ntrunc = CS%ntrunc + 1 endif enddo ; enddo ; enddo @@ -1426,7 +1422,7 @@ subroutine vertvisc_limit_vel(u, v, h, ADp, CDp, forces, visc, dt, G, GV, US, CS do j=js,je; do I=Isq,Ieq ; if (dowrite(I,j)) then ! Here the diagnostic reporting subroutines are called if ! unphysically large values were found. - call write_u_accel(I, j, u_old, h, ADp, CDp, dt, G, GV, US, CS%PointAccel_CSp, & + call write_u_accel(I, j, u_old, h, ADp, CDp, dt_in_T, G, GV, US, CS%PointAccel_CSp, & vel_report(I,j), forces%taux(I,j)*dt_Rho0, a=CS%a_u, hv=CS%h_u) endif ; enddo ; enddo endif @@ -1500,7 +1496,7 @@ subroutine vertvisc_limit_vel(u, v, h, ADp, CDp, forces, visc, dt, G, GV, US, CS do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie if (abs(v(i,J,k)) < CS%vel_underflow) then ; v(i,J,k) = 0.0 elseif (abs(v(i,J,k)) > maxvel) then - v(i,J,k) = SIGN(truncvel,v(i,J,k)) + v(i,J,k) = SIGN(truncvel, v(i,J,k)) if (h(i,j,k) + h(i,j+1,k) > H_report) CS%ntrunc = CS%ntrunc + 1 endif enddo ; enddo ; enddo @@ -1511,7 +1507,7 @@ subroutine vertvisc_limit_vel(u, v, h, ADp, CDp, forces, visc, dt, G, GV, US, CS do J=Jsq,Jeq; do i=is,ie ; if (dowrite(i,J)) then ! Here the diagnostic reporting subroutines are called if ! unphysically large values were found. - call write_v_accel(i, J, v_old, h, ADp, CDp, dt, G, GV, US, CS%PointAccel_CSp, & + call write_v_accel(i, J, v_old, h, ADp, CDp, dt_in_T, G, GV, US, CS%PointAccel_CSp, & vel_report(i,J), forces%tauy(i,J)*dt_Rho0, a=CS%a_v, hv=CS%h_v) endif ; enddo ; enddo endif From 3704a665a3d9803d3be92c9b21f9b474d3e62844 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Wed, 11 Sep 2019 17:08:08 -0400 Subject: [PATCH 032/225] +Add density rescaling factors in unit_scale_type Added factors for power-of-2 rescaling of density to the unit_scale_type, along with the new run-time parameter R_RESCALE_POWER. All answers are bitwise identical, but there is a new runtime parameter, some new elements in a transparent public type, and a new optional variable in the MOM restart files. This adds a new entry to the MOM_parameter_doc.debugging files. --- src/core/MOM.F90 | 4 +++- src/framework/MOM_unit_scaling.F90 | 23 +++++++++++++++++++---- 2 files changed, 22 insertions(+), 5 deletions(-) diff --git a/src/core/MOM.F90 b/src/core/MOM.F90 index 8a5f9bfe02..0e5ca00823 100644 --- a/src/core/MOM.F90 +++ b/src/core/MOM.F90 @@ -2558,7 +2558,7 @@ subroutine MOM_timing_init(CS) id_clock_diabatic = cpu_clock_id('(Ocean diabatic driver)', grain=CLOCK_MODULE_DRIVER) id_clock_continuity = cpu_clock_id('(Ocean continuity equation *)', grain=CLOCK_MODULE) - id_clock_BBL_visc = cpu_clock_id('(Ocean set BBL viscosity)', grain=CLOCK_MODULE) + id_clock_BBL_visc = cpu_clock_id('(Ocean set BBL viscosity)', grain=CLOCK_MODULE) id_clock_pass = cpu_clock_id('(Ocean message passing *)', grain=CLOCK_MODULE) id_clock_MOM_init = cpu_clock_id('(Ocean MOM_initialize_state)', grain=CLOCK_MODULE) id_clock_pass_init = cpu_clock_id('(Ocean init message passing *)', grain=CLOCK_ROUTINE) @@ -2644,6 +2644,8 @@ subroutine set_restart_fields(GV, US, param_file, CS, restart_CSp) "Length unit conversion factor", "L meter-1") call register_restart_field(US%s_to_T_restart, "s_to_T", .false., restart_CSp, & "Time unit conversion factor", "T second-1") + call register_restart_field(US%kg_m3_to_R_restart, "kg_m3_to_R", .false., restart_CSp, & + "Density unit conversion factor", "R m3 kg-1") end subroutine set_restart_fields diff --git a/src/framework/MOM_unit_scaling.F90 b/src/framework/MOM_unit_scaling.F90 index ca174025bf..fe7f95fc79 100644 --- a/src/framework/MOM_unit_scaling.F90 +++ b/src/framework/MOM_unit_scaling.F90 @@ -16,8 +16,10 @@ module MOM_unit_scaling real :: Z_to_m !< A constant that translates distances in the units of depth to meters. real :: m_to_L !< A constant that translates lengths in meters to the units of horizontal lengths. real :: L_to_m !< A constant that translates lengths in the units of horizontal lengths to meters. - real :: s_to_T !< A constant that time intervals in seconds to the units of time. - real :: T_to_s !< A constant that the units of time to seconds. + real :: s_to_T !< A constant that translates time intervals in seconds to the units of time. + real :: T_to_s !< A constant that translates the units of time to seconds. + real :: R_to_kg_m3 !< A constant that translates the units of density to kilograms per meter cubed. + real :: kg_m3_to_R !< A constant that translates kilograms per meter cubed to the units of density. ! These are useful combinations of the fundamental scale conversion factors above. real :: Z_to_L !< Convert vertical distances to lateral lengths @@ -32,6 +34,7 @@ module MOM_unit_scaling real :: m_to_Z_restart = 0.0 !< A copy of the m_to_Z that is used in restart files. real :: m_to_L_restart = 0.0 !< A copy of the m_to_L that is used in restart files. real :: s_to_T_restart = 0.0 !< A copy of the s_to_T that is used in restart files. + real :: kg_m3_to_R_restart = 0.0 !< A copy of the kg_m3_to_R that is used in restart files. end type unit_scale_type contains @@ -44,8 +47,8 @@ subroutine unit_scaling_init( param_file, US ) ! This routine initializes a unit_scale_type structure (US). ! Local variables - integer :: Z_power, L_power, T_power - real :: Z_rescale_factor, L_rescale_factor, T_rescale_factor + integer :: Z_power, L_power, T_power, R_power + real :: Z_rescale_factor, L_rescale_factor, T_rescale_factor, R_rescale_factor ! This include declares and sets the variable "version". # include "version_variable.h" character(len=16) :: mdl = "MOM_unit_scaling" @@ -69,12 +72,18 @@ subroutine unit_scaling_init( param_file, US ) "An integer power of 2 that is used to rescale the model's "//& "intenal units of time. Valid values range from -300 to 300.", & units="nondim", default=0, debuggingParam=.true.) + call get_param(param_file, mdl, "R_RESCALE_POWER", R_power, & + "An integer power of 2 that is used to rescale the model's "//& + "intenal units of density. Valid values range from -300 to 300.", & + units="nondim", default=0, debuggingParam=.true.) if (abs(Z_power) > 300) call MOM_error(FATAL, "unit_scaling_init: "//& "Z_RESCALE_POWER is outside of the valid range of -300 to 300.") if (abs(L_power) > 300) call MOM_error(FATAL, "unit_scaling_init: "//& "L_RESCALE_POWER is outside of the valid range of -300 to 300.") if (abs(T_power) > 300) call MOM_error(FATAL, "unit_scaling_init: "//& "T_RESCALE_POWER is outside of the valid range of -300 to 300.") + if (abs(R_power) > 300) call MOM_error(FATAL, "unit_scaling_init: "//& + "R_RESCALE_POWER is outside of the valid range of -300 to 300.") Z_rescale_factor = 1.0 if (Z_power /= 0) Z_rescale_factor = 2.0**Z_power @@ -91,6 +100,11 @@ subroutine unit_scaling_init( param_file, US ) US%T_to_s = 1.0 * T_rescale_factor US%s_to_T = 1.0 / T_rescale_factor + R_rescale_factor = 1.0 + if (R_power /= 0) R_rescale_factor = 2.0**R_power + US%R_to_kg_m3 = 1.0 * R_rescale_factor + US%kg_m3_to_R = 1.0 / R_rescale_factor + ! These are useful combinations of the fundamental scale conversion factors set above. US%Z_to_L = US%Z_to_m * US%m_to_L US%L_to_Z = US%L_to_m * US%m_to_Z @@ -111,6 +125,7 @@ subroutine fix_restart_unit_scaling(US) US%m_to_Z_restart = US%m_to_Z US%m_to_L_restart = US%m_to_L US%s_to_T_restart = US%s_to_T + US%kg_m3_to_R_restart = US%kg_m3_to_R end subroutine fix_restart_unit_scaling From 5274403ec806b4187d2e8219ab1dd12048068805 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Thu, 12 Sep 2019 18:00:28 -0400 Subject: [PATCH 033/225] +Add H_to_RZ and RZ_to_H to the verticalGrid_type Added two new dimensional conversion factors, H_to_RZ and RZ_to_H, to the MOM6 vertical grid, in preparation for adding testing of dimensional rescaling of density to the MOM6 code. All answers are bitwise identical, but a transparent type has two new elements. --- src/core/MOM_verticalGrid.F90 | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/src/core/MOM_verticalGrid.F90 b/src/core/MOM_verticalGrid.F90 index c11de0d0dd..66ff737bff 100644 --- a/src/core/MOM_verticalGrid.F90 +++ b/src/core/MOM_verticalGrid.F90 @@ -61,6 +61,8 @@ module MOM_verticalGrid real :: H_to_Pa !< A constant that translates the units of thickness to pressure [Pa]. real :: H_to_Z !< A constant that translates thickness units to the units of depth. real :: Z_to_H !< A constant that translates depth units to thickness units. + real :: H_to_RZ !< A constant that translates thickness units to the units of mass per unit area. + real :: RZ_to_H !< A constant that translates mass per unit area units to thickness units. real :: m_to_H_restart = 0.0 !< A copy of the m_to_H that is used in restart files. end type verticalGrid_type @@ -156,6 +158,9 @@ subroutine verticalGridInit( param_file, GV, US ) GV%Z_to_H = US%Z_to_m * GV%m_to_H GV%Angstrom_Z = US%m_to_Z * GV%Angstrom_m + GV%H_to_RZ = GV%H_to_kg_m2 * US%kg_m3_to_R * US%m_to_Z + GV%RZ_to_H = GV%kg_m2_to_H * US%R_to_kg_m3 * US%Z_to_m + ! Log derivative values. call log_param(param_file, mdl, "M to THICKNESS", GV%m_to_H*H_rescale_factor) call log_param(param_file, mdl, "M to THICKNESS rescaled by 2^-n", GV%m_to_H) From fe4428920a27aa4ca26071a6dfd394bd47f3290c Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Thu, 12 Sep 2019 16:50:56 -0600 Subject: [PATCH 034/225] Fixes a bug when computing PE_release_h PE_release_h was being computed outside of a k-loop and, therefore, it was always zero. --- .../lateral/MOM_thickness_diffuse.F90 | 40 +++++++++---------- 1 file changed, 19 insertions(+), 21 deletions(-) diff --git a/src/parameterizations/lateral/MOM_thickness_diffuse.F90 b/src/parameterizations/lateral/MOM_thickness_diffuse.F90 index 3ebf159e3d..9d170e57c8 100644 --- a/src/parameterizations/lateral/MOM_thickness_diffuse.F90 +++ b/src/parameterizations/lateral/MOM_thickness_diffuse.F90 @@ -1263,26 +1263,27 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV enddo endif + if (find_work) then + do j=js,je ; do i=is,ie + Work_h = 0.5 * G%IareaT(i,j) * & + ((Work_u(I-1,j) + Work_u(I,j)) + (Work_v(i,J-1) + Work_v(i,J))) + if (associated(CS%GMwork)) CS%GMwork(i,j) = Work_h + if (associated(MEKE) .and. associated(MEKE%GM_src) .and. .not. CS%GM_src_alt ) then + MEKE%GM_src(i,j) = MEKE%GM_src(i,j) + Work_h + endif + enddo ; enddo - !if (find_work) then ; do j=js,je ; do i=is,ie ; do k=nz,1,-1 - if (find_work) then ; do j=js,je ; do i=is,ie - ! Note that the units of Work_v and Work_u are W, while Work_h is W m-2. - Work_h = 0.5 * G%IareaT(i,j) * & - ((Work_u(I-1,j) + Work_u(I,j)) + (Work_v(i,J-1) + Work_v(i,J))) - PE_release_h = -0.25*(Kh_u(I,j,k)*(Slope_x_PE(I,j,k)**2) * hN2_x_PE(I,j,k) + & - Kh_u(I-1,j,k)*(Slope_x_PE(I-1,j,k)**2) * hN2_x_PE(I-1,j,k) + & - Kh_v(i,J,k)*(Slope_y_PE(i,J,k)**2) * hN2_y_PE(i,J,k) + & - Kh_v(i,J-1,k)*(Slope_y_PE(i,J-1,k)**2) * hN2_y_PE(i,J-1,k)) - if (associated(CS%GMwork)) CS%GMwork(i,j) = Work_h - if (associated(MEKE)) then ; if (associated(MEKE%GM_src)) then - if (CS%GM_src_alt) then + if (associated(MEKE) .and. associated(MEKE%GM_src) .and. CS%GM_src_alt) then + do j=js,je ; do i=is,ie ; do k=nz,1,-1 + PE_release_h = -0.25*(Kh_u(I,j,k)*(Slope_x_PE(I,j,k)**2) * hN2_x_PE(I,j,k) + & + Kh_u(I-1,j,k)*(Slope_x_PE(I-1,j,k)**2) * hN2_x_PE(I-1,j,k) + & + Kh_v(i,J,k)*(Slope_y_PE(i,J,k)**2) * hN2_y_PE(i,J,k) + & + Kh_v(i,J-1,k)*(Slope_y_PE(i,J-1,k)**2) * hN2_y_PE(i,J-1,k)) + MEKE%GM_src(i,j) = MEKE%GM_src(i,j) + PE_release_h - else - MEKE%GM_src(i,j) = MEKE%GM_src(i,j) + Work_h - endif - endif ; endif - !enddo ; enddo ; enddo ; endif - enddo ; enddo ; endif + enddo ; enddo ; enddo + endif + endif if (CS%id_slope_x > 0) call post_data(CS%id_slope_x, CS%diagSlopeX, CS%diag) if (CS%id_slope_y > 0) call post_data(CS%id_slope_y, CS%diagSlopeY, CS%diag) @@ -2071,9 +2072,6 @@ end subroutine thickness_diffuse_end !! wave-speed or the equivalent barotropic mode wave-speed. !! \f$N_*^2 = \max(N^2,0)\f$ is a non-negative form of the square of the Brunt-Vaisala frequency. !! The parameter \f$\gamma_F\f$ is used to reduce the vertical smoothing length scale. -!! This elliptic form for \f$ \psi \f$ is turned on with the logical KHTH_USE_FGNV_STREAMFUNCTION. -!! -!! Thickness diffusivities are calculated independently at u- and v-points using the following expression !! \f[ !! \kappa_h = \left( \kappa_o + \alpha_{s} L_{s}^2 < S N > + \alpha_{M} \kappa_{M} \right) r(\Delta x,L_d) !! \f] From 584177f56e126218da3c189ad04ee98e40b8421d Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Fri, 13 Sep 2019 09:33:59 -0600 Subject: [PATCH 035/225] Deletes whitespace --- src/parameterizations/lateral/MOM_thickness_diffuse.F90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/parameterizations/lateral/MOM_thickness_diffuse.F90 b/src/parameterizations/lateral/MOM_thickness_diffuse.F90 index 9d170e57c8..728e384f75 100644 --- a/src/parameterizations/lateral/MOM_thickness_diffuse.F90 +++ b/src/parameterizations/lateral/MOM_thickness_diffuse.F90 @@ -1268,9 +1268,9 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV Work_h = 0.5 * G%IareaT(i,j) * & ((Work_u(I-1,j) + Work_u(I,j)) + (Work_v(i,J-1) + Work_v(i,J))) if (associated(CS%GMwork)) CS%GMwork(i,j) = Work_h - if (associated(MEKE) .and. associated(MEKE%GM_src) .and. .not. CS%GM_src_alt ) then + if (associated(MEKE) .and. associated(MEKE%GM_src) .and. .not. CS%GM_src_alt ) then MEKE%GM_src(i,j) = MEKE%GM_src(i,j) + Work_h - endif + endif enddo ; enddo if (associated(MEKE) .and. associated(MEKE%GM_src) .and. CS%GM_src_alt) then @@ -1279,7 +1279,7 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV Kh_u(I-1,j,k)*(Slope_x_PE(I-1,j,k)**2) * hN2_x_PE(I-1,j,k) + & Kh_v(i,J,k)*(Slope_y_PE(i,J,k)**2) * hN2_y_PE(i,J,k) + & Kh_v(i,J-1,k)*(Slope_y_PE(i,J-1,k)**2) * hN2_y_PE(i,J-1,k)) - + MEKE%GM_src(i,j) = MEKE%GM_src(i,j) + PE_release_h enddo ; enddo ; enddo endif From 8eaa01c1764dee1b212ef7f70c374f0947332775 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Fri, 13 Sep 2019 11:52:19 -0400 Subject: [PATCH 036/225] Rescaled density units in MOM_bulk_mixedlayer Rescaled density units in MOM_bulk_mixedlayer for dimensional consistency testing. All answers are bitwise identical. --- .../vertical/MOM_bulk_mixed_layer.F90 | 176 +++++++++--------- 1 file changed, 92 insertions(+), 84 deletions(-) diff --git a/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 b/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 index 2a17bfbd6f..ea7a740df5 100644 --- a/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 +++ b/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 @@ -131,9 +131,9 @@ module MOM_bulk_mixed_layer diag_TKE_mixing, & !< The work done by TKE to deepen the mixed layer. diag_TKE_conv_s2, & !< The convective source of TKE due to to mixing in sigma2. diag_PE_detrain, & !< The spurious source of potential energy due to mixed layer - !! detrainment [kg T-3 Z m-1 ~> W m-2]. + !! detrainment [R Z L2 T-3 ~> W m-2]. diag_PE_detrain2 !< The spurious source of potential energy due to mixed layer only - !! detrainment [kg T-3 Z m-1 ~> W m-2]. + !! detrainment [R Z L2 T-3 ~> W m-2]. logical :: allow_clocks_in_omp_loops !< If true, clocks can be called from inside loops that can !! be threaded. To run with multiple threads, set to False. type(group_pass_type) :: pass_h_sum_hmbl_prev !< For group halo pass @@ -244,8 +244,8 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt_in_T, ea, eb, G, GV, h, & ! The layer thickness [H ~> m or kg m-2]. T, & ! The layer temperatures [degC]. S, & ! The layer salinities [ppt]. - R0, & ! The potential density referenced to the surface [kg m-3]. - Rcv ! The coordinate variable potential density [kg m-3]. + R0, & ! The potential density referenced to the surface [R ~> kg m-3]. + Rcv ! The coordinate variable potential density [R ~> kg m-3]. real, dimension(SZI_(G),SZK_(GV)) :: & u, & ! The zonal velocity [L T-1 ~> m s-1]. v, & ! The meridional velocity [L T-1 ~> m s-1]. @@ -269,9 +269,9 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt_in_T, ea, eb, G, GV, htot, & ! The total depth of the layers being considered for ! entrainment [H ~> m or kg m-2]. R0_tot, & ! The integrated potential density referenced to the surface - ! of the layers which are fully entrained [H kg m-3 ~> kg m-2 or kg2 m-5]. + ! of the layers which are fully entrained [H R ~> kg m-2 or kg2 m-5]. Rcv_tot, & ! The integrated coordinate value potential density of the - ! layers that are fully entrained [H kg m-3 ~> kg m-2 or kg2 m-5]. + ! layers that are fully entrained [H R ~> kg m-2 or kg2 m-5]. Ttot, & ! The integrated temperature of layers which are fully ! entrained [degC H ~> degC m or degC kg m-2]. Stot, & ! The integrated salt of layers which are fully entrained @@ -293,13 +293,13 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt_in_T, ea, eb, G, GV, p_ref_cv, & ! Reference pressure for the potential density which defines ! the coordinate variable, set to P_Ref [Pa]. dR0_dT, & ! Partial derivative of the mixed layer potential density with - ! temperature [kg m-3 degC-1]. + ! temperature [R degC-1 ~> kg m-3 degC-1]. dRcv_dT, & ! Partial derivative of the coordinate variable potential - ! density in the mixed layer with temperature [kg m-3 degC-1]. + ! density in the mixed layer with temperature [R degC-1 ~> kg m-3 degC-1]. dR0_dS, & ! Partial derivative of the mixed layer potential density with - ! salinity [kg m-3 ppt-1]. + ! salinity [R ppt-1 ~> kg m-3 ppt-1]. dRcv_dS, & ! Partial derivative of the coordinate variable potential - ! density in the mixed layer with salinity [kg m-3 ppt-1]. + ! density in the mixed layer with salinity [R ppt-1 ~> kg m-3 ppt-1]. TKE_river ! The source of turbulent kinetic energy available for mixing ! at rivermouths [Z L2 T-3 ~> m3 s-3]. @@ -312,7 +312,7 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt_in_T, ea, eb, G, GV, real :: cMKE(2,SZI_(G)) ! Coefficients of HpE and HpE^2 used in calculating the ! denominator of MKE_rate; the two elements have differing ! units of [H-1 ~> m-1 or m2 kg-1] and [H-2 ~> m-2 or m4 kg-2]. - real :: Irho0 ! 1.0 / rho_0 [m3 kg-1] + real :: Irho0 ! 1.0 / rho_0 [R-1 ~> m3 kg-1] real :: Inkml, Inkmlm1! 1.0 / REAL(nkml) and 1.0 / REAL(nkml-1) real :: Ih ! The inverse of a thickness [H-1 ~> m-1 or m2 kg-1]. real :: Idt_diag ! The inverse of the timestep used for diagnostics [T-1 ~> s-1]. @@ -372,7 +372,7 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt_in_T, ea, eb, G, GV, ! dt_in_T = dt * US%s_to_T - Irho0 = 1.0 / GV%Rho0 + Irho0 = 1.0 / (US%kg_m3_to_R*GV%Rho0) dt__diag = dt_in_T ; if (present(dt_diag)) dt__diag = dt_diag Idt_diag = 1.0 / (dt__diag) write_diags = .true. ; if (present(last_call)) write_diags = last_call @@ -471,11 +471,19 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt_in_T, ea, eb, G, GV, is, ie-is+1, tv%eqn_of_state) call calculate_density_derivs(T(:,1), S(:,1), p_ref_cv, dRcv_dT, dRcv_dS, & is, ie-is+1, tv%eqn_of_state) + if (US%R_to_kg_m3 /= 1.0) then ; do i=is,ie + dR0_dT(i) = US%kg_m3_to_R * dR0_dT(i) ; dR0_dS(i) = US%kg_m3_to_R * dR0_dS(i) + dRcv_dT(i) = US%kg_m3_to_R * dRcv_dT(i) ; dRcv_dS(i) = US%kg_m3_to_R * dRcv_dS(i) + enddo ; endif do k=1,nz call calculate_density(T(:,k), S(:,k), p_ref, R0(:,k), is, ie-is+1, & tv%eqn_of_state) call calculate_density(T(:,k), S(:,k), p_ref_cv, Rcv(:,k), is, & ie-is+1, tv%eqn_of_state) + if (US%kg_m3_to_R /= 1.0) then ; do i=is,ie + R0(i,k) = US%kg_m3_to_R * R0(i,k) + Rcv(i,k) = US%kg_m3_to_R * Rcv(i,k) + enddo ; endif enddo if (id_clock_EOS>0) call cpu_clock_end(id_clock_EOS) @@ -517,7 +525,7 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt_in_T, ea, eb, G, GV, RmixConst = 0.5*CS%rivermix_depth * (GV%g_Earth*US%m_to_Z) * Irho0**2 do i=is,ie TKE_river(i) = max(0.0, RmixConst*dR0_dS(i)* & - US%T_to_s*(fluxes%lrunoff(i,j) + fluxes%frunoff(i,j)) * S(i,1)) + US%T_to_s*US%kg_m3_to_R*(fluxes%lrunoff(i,j) + fluxes%frunoff(i,j)) * S(i,1)) enddo else do i=is,ie ; TKE_river(i) = 0.0 ; enddo @@ -606,7 +614,7 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt_in_T, ea, eb, G, GV, if (CS%ML_resort) then if (id_clock_resort>0) call cpu_clock_begin(id_clock_resort) - call resort_ML(h(:,0:), T(:,0:), S(:,0:), R0(:,0:), Rcv(:,0:), GV%Rlay, eps, & + call resort_ML(h(:,0:), T(:,0:), S(:,0:), R0(:,0:), Rcv(:,0:), US%kg_m3_to_R*GV%Rlay(:), eps, & d_ea, d_eb, ksort, G, GV, CS, dR0_dT, dR0_dS, dRcv_dT, dRcv_dS) if (id_clock_resort>0) call cpu_clock_end(id_clock_resort) endif @@ -642,11 +650,11 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt_in_T, ea, eb, G, GV, if (id_clock_detrain>0) call cpu_clock_begin(id_clock_detrain) if (CS%nkbl == 1) then call mixedlayer_detrain_1(h(:,0:), T(:,0:), S(:,0:), R0(:,0:), Rcv(:,0:), & - GV%Rlay, dt_in_T, dt__diag, d_ea, d_eb, j, G, GV, US, CS, & + US%kg_m3_to_R*GV%Rlay(:), dt_in_T, dt__diag, d_ea, d_eb, j, G, GV, US, CS, & dRcv_dT, dRcv_dS, max_BL_det) elseif (CS%nkbl == 2) then call mixedlayer_detrain_2(h(:,0:), T(:,0:), S(:,0:), R0(:,0:), Rcv(:,0:), & - GV%Rlay, dt_in_T, dt__diag, d_ea, j, G, GV, US, CS, & + US%kg_m3_to_R*GV%Rlay(:), dt_in_T, dt__diag, d_ea, j, G, GV, US, CS, & dR0_dT, dR0_dS, dRcv_dT, dRcv_dS, max_BL_det) else ! CS%nkbl not = 1 or 2 ! This code only works with 1 or 2 buffer layers. @@ -814,9 +822,9 @@ subroutine convective_adjustment(h, u, v, R0, Rcv, T, S, eps, d_eb, & real, dimension(SZI_(G),SZK_(GV)), intent(inout) :: T !< Layer temperatures [degC]. real, dimension(SZI_(G),SZK_(GV)), intent(inout) :: S !< Layer salinities [ppt]. real, dimension(SZI_(G),SZK_(GV)), intent(inout) :: R0 !< Potential density referenced to - !! surface pressure [kg m-3]. + !! surface pressure [R ~> kg m-3]. real, dimension(SZI_(G),SZK_(GV)), intent(inout) :: Rcv !< The coordinate defining potential - !! density [kg m-3]. + !! density [R ~> kg m-3]. real, dimension(SZI_(G),SZK_(GV)), intent(inout) :: d_eb !< The downward increase across a layer !! in the entrainment from below [H ~> m or kg m-2]. !! Positive values go with mass gain by @@ -845,9 +853,9 @@ subroutine convective_adjustment(h, u, v, R0, Rcv, T, S, eps, d_eb, & htot, & ! The total depth of the layers being considered for ! entrainment [H ~> m or kg m-2]. R0_tot, & ! The integrated potential density referenced to the surface - ! of the layers which are fully entrained [H kg m-3 ~> kg m-2 or kg2 m-5]. + ! of the layers which are fully entrained [H R ~> kg m-2 or kg2 m-5]. Rcv_tot, & ! The integrated coordinate value potential density of the - ! layers that are fully entrained [H kg m-3 ~> kg m-2 or kg2 m-5]. + ! layers that are fully entrained [H R ~> kg m-2 or kg2 m-5]. Ttot, & ! The integrated temperature of layers which are fully ! entrained [degC H ~> degC m or degC kg m-2]. Stot, & ! The integrated salt of layers which are fully entrained @@ -861,11 +869,11 @@ subroutine convective_adjustment(h, u, v, R0, Rcv, T, S, eps, d_eb, & real :: Ih ! The inverse of a thickness [H-1 ~> m-1 or m2 kg-1]. real :: g_H2_2Rho0 ! Half the gravitational acceleration times the square of ! the conversion from H to Z divided by the mean density, - ! in [L2 Z m3 T-3 H-2 kg-1 ~> m4 s-3 kg-1 or m10 s-3 kg-3]. + ! in [L2 Z T-3 H-2 R-1 ~> m4 s-3 kg-1 or m10 s-3 kg-3]. integer :: is, ie, nz, i, k, k1, nzc, nkmb is = G%isc ; ie = G%iec ; nz = GV%ke - g_H2_2Rho0 = (GV%g_Earth * GV%H_to_Z**2) / (2.0 * GV%Rho0) + g_H2_2Rho0 = (GV%g_Earth * GV%H_to_Z**2) / (2.0 * US%kg_m3_to_R*GV%Rho0) nzc = nz ; if (present(nz_conv)) nzc = nz_conv nkmb = CS%nkml+CS%nkbl @@ -959,9 +967,9 @@ subroutine mixedlayer_convection(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & real, dimension(SZI_(G)), intent(out) :: vhtot !< The integrated mixed layer meridional !! velocity [H L T-1 ~> m2 s-1 or kg m-1 s-1]. real, dimension(SZI_(G)), intent(out) :: R0_tot !< The integrated mixed layer potential density referenced - !! to 0 pressure [H kg m-2 ~> kg m-1 or kg2 m-4]. + !! to 0 pressure [H R ~> kg m-2 or kg2 m-5]. real, dimension(SZI_(G)), intent(out) :: Rcv_tot !< The integrated mixed layer coordinate - !! variable potential density [H kg m-2 ~> kg m-1 or kg2 m-4]. + !! variable potential density [H R ~> kg m-2 or kg2 m-5]. real, dimension(SZI_(G),SZK_(GV)), & intent(in) :: u !< Zonal velocities interpolated to h points [L T-1 ~> m s-1]. real, dimension(SZI_(G),SZK_(GV)), & @@ -972,21 +980,21 @@ subroutine mixedlayer_convection(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & intent(in) :: S !< Layer salinities [ppt]. real, dimension(SZI_(G),SZK_(GV)), & intent(in) :: R0 !< Potential density referenced to - !! surface pressure [kg m-3]. + !! surface pressure [R ~> kg m-3]. real, dimension(SZI_(G),SZK_(GV)), & intent(in) :: Rcv !< The coordinate defining potential - !! density [kg m-3]. + !! density [R ~> kg m-3]. real, dimension(SZI_(G),SZK_(GV)), & intent(in) :: eps !< The negligibly small amount of water !! that will be left in each layer [H ~> m or kg m-2]. real, dimension(SZI_(G)), intent(in) :: dR0_dT !< The partial derivative of R0 with respect to - !! temperature [kg m-3 degC-1]. + !! temperature [R degC-1 ~> kg m-3 degC-1]. real, dimension(SZI_(G)), intent(in) :: dRcv_dT !< The partial derivative of Rcv with respect to - !! temperature [kg m-3 degC-1]. + !! temperature [R degC-1 ~> kg m-3 degC-1]. real, dimension(SZI_(G)), intent(in) :: dR0_dS !< The partial derivative of R0 with respect to - !! salinity [kg m-3 ppt-1]. + !! salinity [R ppt-1 ~> kg m-3 ppt-1]. real, dimension(SZI_(G)), intent(in) :: dRcv_dS !< The partial derivative of Rcv with respect to - !! salinity [kg m-3 ppt-1]. + !! salinity [R ppt-1 ~> kg m-3 ppt-1]. real, dimension(SZI_(G)), intent(in) :: netMassInOut !< The net mass flux (if non-Boussinesq) !! or volume flux (if Boussinesq) into the ocean !! within a time step [H ~> m or kg m-2]. (I.e. P+R-E.) @@ -1043,9 +1051,9 @@ subroutine mixedlayer_convection(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & real :: T_precip ! The temperature of the precipitation [degC]. real :: C1_3, C1_6 ! 1/3 and 1/6. real :: En_fn, Frac, x1 ! Nondimensional temporary variables. - real :: dr, dr0 ! Temporary variables [kg m-3 H ~> kg m-2 or kg2 m-5]. - real :: dr_ent, dr_comp ! Temporary variables [kg m-3 H ~> kg m-2 or kg2 m-5]. - real :: dr_dh ! The partial derivative of dr_ent with h_ent [kg m-3]. + real :: dr, dr0 ! Temporary variables [R H ~> kg m-2 or kg2 m-5]. + real :: dr_ent, dr_comp ! Temporary variables [R H ~> kg m-2 or kg2 m-5]. + real :: dr_dh ! The partial derivative of dr_ent with h_ent [R ~> kg m-3]. real :: h_min, h_max ! The minimum, maximum, and previous estimates for real :: h_prev ! h_ent [H ~> m or kg m-2]. real :: h_evap ! The thickness that is evaporated [H ~> m or kg m-2]. @@ -1053,22 +1061,22 @@ subroutine mixedlayer_convection(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & ! h_ent between iterations [H ~> m or kg m-2]. real :: g_H2_2Rho0 ! Half the gravitational acceleration times the square of ! the conversion from H to Z divided by the mean density, - ! [L2 Z m3 T-3 H-2 kg-1 ~> m4 s-3 kg-1 or m10 s-3 kg-3]. + ! [L2 Z T-3 H-2 R-1 ~> m4 s-3 kg-1 or m10 s-3 kg-3]. real :: Angstrom ! The minimum layer thickness [H ~> m or kg m-2]. real :: opacity ! The opacity converted to inverse thickness units [H-1 ~> m-1 or m2 kg-1] real :: sum_Pen_En ! The potential energy change due to penetrating ! shortwave radiation, integrated over a layer - ! [H kg m-3 ~> kg m-2 or kg2 m-5]. + ! [H R ~> kg m-2 or kg2 m-5]. real :: Idt ! 1.0/dt [T-1 ~> s-1] real :: netHeatOut ! accumulated heat content of mass leaving ocean integer :: is, ie, nz, i, k, ks, itt, n real, dimension(max(nsw,1)) :: & - C2, & ! Temporary variable [kg m-3 H-1 ~> kg m-4 or m-1]. - r_SW_top ! Temporary variables [H kg m-3 ~> kg m-2 or kg2 m-5]. + C2, & ! Temporary variable R H-1 ~> kg m-4 or m-1]. + r_SW_top ! Temporary variables [H R ~> kg m-2 or kg2 m-5]. Angstrom = GV%Angstrom_H C1_3 = 1.0/3.0 ; C1_6 = 1.0/6.0 - g_H2_2Rho0 = (GV%g_Earth * GV%H_to_Z**2) / (2.0 * GV%Rho0) + g_H2_2Rho0 = (GV%g_Earth * GV%H_to_Z**2) / (2.0 * US%kg_m3_to_R*GV%Rho0) Idt = 1.0 / dt_in_T is = G%isc ; ie = G%iec ; nz = GV%ke @@ -1514,9 +1522,9 @@ subroutine mechanical_entrainment(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & real, dimension(SZI_(G)), intent(inout) :: vhtot !< The integrated mixed layer meridional !! velocity [H L T-1 ~> m2 s-1 or kg m-1 s-1]. real, dimension(SZI_(G)), intent(inout) :: R0_tot !< The integrated mixed layer potential density - !! referenced to 0 pressure [H kg m-3 ~> kg m-2 or kg2 m-5]. + !! referenced to 0 pressure [H R ~> kg m-2 or kg2 m-5]. real, dimension(SZI_(G)), intent(inout) :: Rcv_tot !< The integrated mixed layer coordinate variable - !! potential density [H kg m-3 ~> kg m-2 or kg2 m-5]. + !! potential density [H R ~> kg m-2 or kg2 m-5]. real, dimension(SZI_(G),SZK_(GV)), & intent(in) :: u !< Zonal velocities interpolated to h points [L T-1 ~> m s-1]. real, dimension(SZI_(G),SZK_(GV)), & @@ -1527,17 +1535,17 @@ subroutine mechanical_entrainment(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & intent(in) :: S !< Layer salinities [ppt]. real, dimension(SZI_(G),SZK_(GV)), & intent(in) :: R0 !< Potential density referenced to - !! surface pressure [kg m-3]. + !! surface pressure [R ~> kg m-3]. real, dimension(SZI_(G),SZK_(GV)), & intent(in) :: Rcv !< The coordinate defining potential - !! density [kg m-3]. + !! density [R ~> kg m-3]. real, dimension(SZI_(G),SZK_(GV)), & intent(in) :: eps !< The negligibly small amount of water !! that will be left in each layer [H ~> m or kg m-2]. real, dimension(SZI_(G)), intent(in) :: dR0_dT !< The partial derivative of R0 with respect to - !! temperature [kg m-3 degC-1]. + !! temperature [R degC-1 ~> kg m-3 degC-1]. real, dimension(SZI_(G)), intent(in) :: dRcv_dT !< The partial derivative of Rcv with respect to - !! temperature [kg m-3 degC-1]. + !! temperature [R degC-1 ~> kg m-3 degC-1]. real, dimension(2,SZI_(G)), intent(in) :: cMKE !< Coefficients of HpE and HpE^2 used in calculating the !! denominator of MKE_rate; the two elements have differing !! units of [H-1 ~> m-1 or m2 kg-1] and [H-2 ~> m-2 or m4 kg-2]. @@ -1577,7 +1585,7 @@ subroutine mechanical_entrainment(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & real :: HpE ! The current thickness plus entrainment [H ~> m or kg m-2]. real :: g_H_2Rho0 ! Half the gravitational acceleration times the ! conversion from H to m divided by the mean density, - ! in [L2 m3 T-2 H-1 kg-1 ~> m4 s-2 kg-1 or m7 s-2 kg-2]. + ! in [L2 T-2 H-1 R-1 ~> m4 s-2 kg-1 or m7 s-2 kg-2]. real :: TKE_full_ent ! The TKE remaining if a layer is fully entrained ! [Z L2 T-2 ~> m3 s-2]. real :: dRL ! Work required to mix water from the next layer @@ -1611,7 +1619,7 @@ subroutine mechanical_entrainment(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & integer :: is, ie, nz, i, k, ks, itt, n C1_3 = 1.0/3.0 ; C1_6 = 1.0/6.0 ; C1_24 = 1.0/24.0 - g_H_2Rho0 = (GV%g_Earth * GV%H_to_Z) / (2.0 * GV%Rho0) + g_H_2Rho0 = (GV%g_Earth * GV%H_to_Z) / (2.0 * US%kg_m3_to_R*GV%Rho0) Hmix_min = CS%Hmix_min h_neglect = GV%H_subroundoff is = G%isc ; ie = G%iec ; nz = GV%ke @@ -1837,7 +1845,7 @@ subroutine sort_ML(h, R0, eps, G, GV, CS, ksort) type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure. real, dimension(SZI_(G),SZK_(GV)), intent(in) :: h !< Layer thickness [H ~> m or kg m-2]. real, dimension(SZI_(G),SZK_(GV)), intent(in) :: R0 !< The potential density used to sort - !! the layers [kg m-3]. + !! the layers [R ~> kg m-3]. real, dimension(SZI_(G),SZK_(GV)), intent(in) :: eps !< The (small) thickness that must !! remain in each layer [H ~> m or kg m-2]. type(bulkmixedlayer_CS), pointer :: CS !< The control structure returned by a @@ -1893,11 +1901,11 @@ subroutine resort_ML(h, T, S, R0, Rcv, RcvTgt, eps, d_ea, d_eb, ksort, G, GV, CS real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: T !< Layer temperatures [degC]. real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: S !< Layer salinities [ppt]. real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: R0 !< Potential density referenced to - !! surface pressure [kg m-3]. + !! surface pressure [R ~> kg m-3]. real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: Rcv !< The coordinate defining - !! potential density [kg m-3]. + !! potential density [R ~> kg m-3]. real, dimension(SZK_(GV)), intent(in) :: RcvTgt !< The target value of Rcv for each - !! layer [kg m-3]. + !! layer [R ~> kg m-3]. real, dimension(SZI_(G),SZK_(GV)), intent(inout) :: eps !< The (small) thickness that must !! remain in each layer [H ~> m or kg m-2]. real, dimension(SZI_(G),SZK_(GV)), intent(inout) :: d_ea !< The upward increase across a @@ -1915,19 +1923,19 @@ subroutine resort_ML(h, T, S, R0, Rcv, RcvTgt, eps, d_ea, d_eb, ksort, G, GV, CS real, dimension(SZI_(G)), intent(in) :: dR0_dT !< The partial derivative of !! potential density referenced !! to the surface with potential - !! temperature [kg m-3 degC-1]. + !! temperature [R degC-1 ~> kg m-3 degC-1]. real, dimension(SZI_(G)), intent(in) :: dR0_dS !< The partial derivative of !! cpotential density referenced !! to the surface with salinity, - !! [kg m-3 ppt-1]. + !! [R ppt-1 ~> kg m-3 ppt-1]. real, dimension(SZI_(G)), intent(in) :: dRcv_dT !< The partial derivative of !! coordinate defining potential !! density with potential - !! temperature [kg m-3 degC-1]. + !! temperature [R degC-1 ~> kg m-3 degC-1]. real, dimension(SZI_(G)), intent(in) :: dRcv_dS !< The partial derivative of !! coordinate defining potential !! density with salinity, - !! [kg m-3 ppt-1]. + !! [R ppt-1 ~> kg m-3 ppt-1]. ! If there are no massive light layers above the deepest of the mixed- and ! buffer layers, do nothing (except perhaps to reshuffle these layers). @@ -2213,11 +2221,11 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: T !< Potential temperature [degC]. real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: S !< Salinity [ppt]. real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: R0 !< Potential density referenced to - !! surface pressure [kg m-3]. + !! surface pressure [R ~> kg m-3]. real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: Rcv !< The coordinate defining potential - !! density [kg m-3]. + !! density [R ~> kg m-3]. real, dimension(SZK_(GV)), intent(in) :: RcvTgt !< The target value of Rcv for each - !! layer [kg m-3]. + !! layer [R ~> kg m-3]. real, intent(in) :: dt_in_T !< Time increment [T ~> s]. real, intent(in) :: dt_diag !< The diagnostic time step [T ~> s]. real, dimension(SZI_(G),SZK_(GV)), intent(inout) :: d_ea !< The upward increase across a layer in @@ -2231,18 +2239,18 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea real, dimension(SZI_(G)), intent(in) :: dR0_dT !< The partial derivative of !! potential density referenced to the !! surface with potential temperature, - !! [kg m-3 degC-1]. + !! [R degC-1 ~> kg m-3 degC-1]. real, dimension(SZI_(G)), intent(in) :: dR0_dS !< The partial derivative of !! cpotential density referenced to the !! surface with salinity - !! [kg m-3 ppt-1]. + !! [R ppt-1 ~> kg m-3 ppt-1]. real, dimension(SZI_(G)), intent(in) :: dRcv_dT !< The partial derivative of !! coordinate defining potential density !! with potential temperature, - !! [kg m-3 degC-1]. + !! [R degC-1 ~> kg m-3 degC-1]. real, dimension(SZI_(G)), intent(in) :: dRcv_dS !< The partial derivative of !! coordinate defining potential density - !! with salinity [kg m-3 ppt-1]. + !! with salinity [R ppt-1 ~> kg m-3 ppt-1]. real, dimension(SZI_(G)), intent(in) :: max_BL_det !< If non-negative, the maximum !! detrainment permitted from the buffer !! layers [H ~> m or kg m-2]. @@ -2255,9 +2263,9 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea real :: h_to_bl ! The total thickness detrained to the buffer ! layers [H ~> m or kg m-2]. real :: R0_to_bl ! The depth integrated amount of R0 that is detrained to the - ! buffer layer [H kg m-3 ~> kg m-2 or kg2 m-5] + ! buffer layer [H R ~> kg m-2 or kg2 m-5] real :: Rcv_to_bl ! The depth integrated amount of Rcv that is detrained to the - ! buffer layer [H kg m-3 ~> kg m-2 or kg2 m-5] + ! buffer layer [H R ~> kg m-2 or kg2 m-5] real :: T_to_bl ! The depth integrated amount of T that is detrained to the ! buffer layer [degC H ~> degC m or degC kg m-2] real :: S_to_bl ! The depth integrated amount of S that is detrained to the @@ -2282,7 +2290,7 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea ! layer that remains [H ~> m or kg m-2]. real :: stays_min_merge ! The minimum allowed value of stays_merge [H ~> m or kg m-2]. - real :: dR0_2dz, dRcv_2dz ! Half the vertical gradients of R0 and Rcv [kg m-3 H-1 ~> kg m-4 or m-1] + real :: dR0_2dz, dRcv_2dz ! Half the vertical gradients of R0 and Rcv [R H-1 ~> kg m-4 or m-1] ! real :: dT_2dz, dS_2dz ! Half the vertical gradients of T and S, in degC H-1, and ppt H-1. real :: scale_slope ! A nondimensional number < 1 used to scale down ! the slope within the upper buffer layer when @@ -2293,7 +2301,7 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea ! rho_0*g [H2 ~> m2 or kg2 m-4]. real :: dPE_det, dPE_merge ! The energy required to mix the detrained water ! into the buffer layer or the merge the two - ! buffer layers [kg H2 Z T-2 L-2 m-1 ~> J m-2 or J kg2 m-8]. + ! buffer layers [R H2 L2 Z-1 T-2 ~> J m-2 or J kg2 m-8]. real :: h_from_ml ! The amount of additional water that must be ! drawn from the mixed layer [H ~> m or kg m-2]. @@ -2308,18 +2316,18 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea ! interior layers that are just lighter and ! just denser than the lower buffer layer. - real :: R0_det, T_det, S_det ! Detrained values of R0 [kg m-3], T [degC], and S [ppt]. + real :: R0_det, T_det, S_det ! Detrained values of R0 [R ~> kg m-3], T [degC], and S [ppt]. real :: Rcv_stays, R0_stays ! Values of Rcv and R0 that stay in a layer. real :: T_stays, S_stays ! Values of T and S that stay in a layer. real :: dSpice_det, dSpice_stays! The spiciness difference between an original ! buffer layer and the water that moves into ! an interior layer or that stays in that - ! layer [kg m-3]. + ! layer [R ~> kg m-3]. real :: dSpice_lim, dSpice_lim2 ! Limits to the spiciness difference between ! the lower buffer layer and the water that - ! moves into an interior layer [kg m-3]. + ! moves into an interior layer [R ~> kg m-3]. real :: dSpice_2dz ! The vertical gradient of spiciness used for - ! advection [kg m-3 H-1 ~> kg m-4 or m-1]. + ! advection [R H-1 ~> kg m-4 or m-1]. real :: dPE_ratio ! Multiplier of dPE_det at which merging is ! permitted - here (detrainment_per_day/dt)*30 @@ -2333,8 +2341,8 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea real :: I_denom ! A work variable with units of [ppt2 m6 kg-2]. real :: g_2 ! 1/2 g_Earth [L2 Z-1 T-2 ~> m s-2]. - real :: Rho0xG ! Rho0 times G_Earth [kg L2 m-3 Z-1 T-2 ~> kg m-2 s-2]. - real :: I2Rho0 ! 1 / (2 Rho0) [m3 kg-1]. + real :: Rho0xG ! Rho0 times G_Earth [R L2 Z-1 T-2 ~> kg m-2 s-2]. + real :: I2Rho0 ! 1 / (2 Rho0) [R-1 ~> m3 kg-1]. real :: Idt_H2 ! The square of the conversion from thickness to Z ! divided by the time step [Z2 H-2 T-1 ~> s-1 or m6 kg-2 s-1]. logical :: stable_Rcv ! If true, the buffer layers are stable with @@ -2349,7 +2357,7 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea real :: Ih, Ihdet, Ih1f, Ih2f ! Assorted inverse thickness work variables, real :: Ihk0, Ihk1, Ih12 ! all in [H-1 ~> m-1 or m2 kg-1]. real :: dR1, dR2, dR2b, dRk1 ! Assorted density difference work variables, - real :: dR0, dR21, dRcv ! all in [kg m-3]. + real :: dR0, dR21, dRcv ! all in [R ~> kg m-3]. real :: dRcv_stays, dRcv_det, dRcv_lim real :: Angstrom ! The minumum layer thickness [H ~> m or kg m-2]. @@ -2362,9 +2370,9 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea nkmb = CS%nkml+CS%nkbl h_neglect = GV%H_subroundoff g_2 = 0.5 * GV%g_Earth - Rho0xG = GV%Rho0 * GV%g_Earth + Rho0xG = US%kg_m3_to_R*GV%Rho0 * GV%g_Earth Idt_H2 = GV%H_to_Z**2 / dt_diag - I2Rho0 = 0.5 / GV%Rho0 + I2Rho0 = 0.5 / (US%kg_m3_to_R*GV%Rho0) Angstrom = GV%Angstrom_H ! This is hard coding of arbitrary and dimensional numbers. @@ -2802,7 +2810,7 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea h_det_to_h2*( (R0(i,kb1)-R0_det)*h1 + (R0(i,kb2)-R0_det)*h2 ) + & h_ml_to_h2*( (R0(i,kb2)-R0(i,0))*h2 + (R0(i,kb1)-R0(i,0))*h1 + & (R0_det-R0(i,0))*h_det_to_h2 ) + & - h_det_to_h1*h_ml_to_h1*(R0_det-R0(i,0))) - 2.0*GV%Rho0*dPE_extrap ) + h_det_to_h1*h_ml_to_h1*(R0_det-R0(i,0))) - 2.0*US%kg_m3_to_R*GV%Rho0*dPE_extrap ) if (allocated(CS%diag_PE_detrain)) & CS%diag_PE_detrain(i,j) = CS%diag_PE_detrain(i,j) + s1en @@ -3104,11 +3112,11 @@ subroutine mixedlayer_detrain_1(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: T !< Potential temperature [degC]. real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: S !< Salinity [ppt]. real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: R0 !< Potential density referenced to - !! surface pressure [kg m-3]. + !! surface pressure [R ~> kg m-3]. real, dimension(SZI_(G),SZK0_(GV)), intent(inout) :: Rcv !< The coordinate defining potential - !! density [kg m-3]. + !! density [R ~> kg m-3]. real, dimension(SZK_(GV)), intent(in) :: RcvTgt !< The target value of Rcv for each - !! layer [kg m-3]. + !! layer [R ~> kg m-3]. real, intent(in) :: dt_in_T !< Time increment [T ~> s]. real, intent(in) :: dt_diag !< The accumulated time interval for !! diagnostics [T ~> s]. @@ -3127,10 +3135,10 @@ subroutine mixedlayer_detrain_1(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea real, dimension(SZI_(G)), intent(in) :: dRcv_dT !< The partial derivative of !! coordinate defining potential density !! with potential temperature - !! [kg m-3 degC-1]. + !! [R degC-1 ~> kg m-3 degC-1]. real, dimension(SZI_(G)), intent(in) :: dRcv_dS !< The partial derivative of !! coordinate defining potential density - !! with salinity [kg m-3 ppt-1]. + !! with salinity [R ppt-1 ~> kg m-3 ppt-1]. real, dimension(SZI_(G)), intent(in) :: max_BL_det !< If non-negative, the maximum !! detrainment permitted from the buffer !! layers [H ~> m or kg m-2]. @@ -3148,7 +3156,7 @@ subroutine mixedlayer_detrain_1(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea real :: dt_Time ! The timestep divided by the detrainment timescale [nondim]. real :: g_H2_2Rho0dt ! Half the gravitational acceleration times the square of the ! conversion from H to m divided by the mean density times the time - ! step [L2 Z m3 T-3 H-2 kg-1 ~> m4 s-3 kg-1 or m10 s-3 kg-3]. + ! step [L2 Z T-3 H-2 R-1 ~> m4 s-3 kg-1 or m10 s-3 kg-3]. real :: g_H2_2dt ! Half the gravitational acceleration times the square of the ! conversion from H to Z divided by the diagnostic time step ! [L2 Z H-2 T-3 ~> m s-3 or m7 kg-2 s-3]. @@ -3163,7 +3171,7 @@ subroutine mixedlayer_detrain_1(h, T, S, R0, Rcv, RcvTgt, dt_in_T, dt_diag, d_ea "CS%nkbl must be 1 in mixedlayer_detrain_1.") dt_Time = dt_in_T / CS%BL_detrain_time - g_H2_2Rho0dt = (GV%g_Earth * GV%H_to_Z**2) / (2.0 * GV%Rho0 * dt_diag) + g_H2_2Rho0dt = (GV%g_Earth * GV%H_to_Z**2) / (2.0 * US%kg_m3_to_R*GV%Rho0 * dt_diag) g_H2_2dt = (GV%g_Earth * GV%H_to_Z**2) / (2.0 * dt_diag) ! Move detrained water into the buffer layer. @@ -3606,10 +3614,10 @@ subroutine bulkmixedlayer_init(Time, G, GV, US, param_file, diag, CS) 'm3 s-3', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_PE_detrain = register_diag_field('ocean_model', 'PE_detrain', diag%axesT1, & Time, 'Spurious source of potential energy from mixed layer detrainment', & - 'W m-2', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) + 'W m-2', conversion=US%R_to_kg_m3*US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_PE_detrain2 = register_diag_field('ocean_model', 'PE_detrain2', diag%axesT1, & Time, 'Spurious source of potential energy from mixed layer only detrainment', & - 'W m-2', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) + 'W m-2', conversion=US%R_to_kg_m3*US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_h_mismatch = register_diag_field('ocean_model', 'h_miss_ML', diag%axesT1, & Time, 'Summed absolute mismatch in entrainment terms', 'm', conversion=US%Z_to_m) CS%id_Hsfc_used = register_diag_field('ocean_model', 'Hs_used', diag%axesT1, & From 466903bd8a85b63ce7c849f9c5da8e700c867970 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Fri, 13 Sep 2019 10:58:44 -0600 Subject: [PATCH 037/225] Deletes more whitespace --- src/parameterizations/lateral/MOM_thickness_diffuse.F90 | 1 - 1 file changed, 1 deletion(-) diff --git a/src/parameterizations/lateral/MOM_thickness_diffuse.F90 b/src/parameterizations/lateral/MOM_thickness_diffuse.F90 index 728e384f75..9a5cfede6f 100644 --- a/src/parameterizations/lateral/MOM_thickness_diffuse.F90 +++ b/src/parameterizations/lateral/MOM_thickness_diffuse.F90 @@ -1279,7 +1279,6 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV Kh_u(I-1,j,k)*(Slope_x_PE(I-1,j,k)**2) * hN2_x_PE(I-1,j,k) + & Kh_v(i,J,k)*(Slope_y_PE(i,J,k)**2) * hN2_y_PE(i,J,k) + & Kh_v(i,J-1,k)*(Slope_y_PE(i,J-1,k)**2) * hN2_y_PE(i,J-1,k)) - MEKE%GM_src(i,j) = MEKE%GM_src(i,j) + PE_release_h enddo ; enddo ; enddo endif From 185235666d2e8799021cc88396849e91ec264bc9 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Mon, 16 Sep 2019 11:13:23 -0600 Subject: [PATCH 038/225] Removes white space --- src/parameterizations/lateral/MOM_hor_visc.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index e9ff06d337..213e9a485f 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -263,7 +263,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, real, dimension(SZIB_(G),SZJB_(G)) :: & dvdx, dudy, & ! components in the shearing strain [s-1] - dDel2vdx, dDel2udy, & ! Components in the biharmonic equivalent of the shearing strain [L-2 T-1 ~> m-2 s-1] + dDel2vdx, dDel2udy, & ! Components in the biharmonic equivalent of the shearing strain [L-2 T-1 ~> m-2 s-1] dvdx_bt, dudy_bt, & ! components in the barotropic shearing strain [s-1] sh_xy, & ! horizontal shearing strain (du/dy + dv/dx) including metric terms [s-1] sh_xy_bt, & ! barotropic horizontal shearing strain (du/dy + dv/dx) inc. metric terms [s-1] From ba2f186aedb18320ca4110729ad65574c6c53e06 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Tue, 17 Sep 2019 09:43:33 -0600 Subject: [PATCH 039/225] Removed diagnostics related to FrictWork_diss and FrictWork_Max * We deleted these terms previously but forgot to remove posting the diagnostics. Thanks to Travis and .testing/tc0 we were able to identidy this mistake. --- .../lateral/MOM_hor_visc.F90 | 19 +------------------ 1 file changed, 1 insertion(+), 18 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index 213e9a485f..b798fb4d86 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -180,8 +180,7 @@ module MOM_hor_visc integer :: id_GME_coeff_h = -1, id_GME_coeff_q = -1 integer :: id_vort_xy_q = -1, id_div_xx_h = -1 integer :: id_FrictWork = -1, id_FrictWorkIntz = -1 - integer :: id_FrictWorkMax = -1 - integer :: id_FrictWork_diss = -1, id_FrictWork_GME = -1 + integer :: id_FrictWork_GME = -1 !!@} @@ -301,8 +300,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, target_diss_rate_GME, & ! the maximum theoretical dissipation plus the amount spuriously dissipated ! by friction [m2 s-3] FrictWork, & ! work done by MKE dissipation mechanisms [W m-2] - FrictWork_diss, & ! negative definite work done by MKE dissipation mechanisms [W m-2] - FrictWorkMax, & ! maximum possible work done by MKE dissipation mechanisms [W m-2] FrictWork_GME, & ! work done by GME [W m-2] div_xx_h ! horizontal divergence [s-1] real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: & @@ -1283,8 +1280,6 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if (CS%id_diffu>0) call post_data(CS%id_diffu, diffu, CS%diag) if (CS%id_diffv>0) call post_data(CS%id_diffv, diffv, CS%diag) if (CS%id_FrictWork>0) call post_data(CS%id_FrictWork, FrictWork, CS%diag) - if (CS%id_FrictWorkMax>0) call post_data(CS%id_FrictWorkMax, FrictWorkMax, CS%diag) - if (CS%id_FrictWork_diss>0) call post_data(CS%id_FrictWork_diss, FrictWork_diss, CS%diag) if (CS%id_FrictWork_GME>0) call post_data(CS%id_FrictWork_GME, FrictWork_GME, CS%diag) if (CS%id_Ah_h>0) call post_data(CS%id_Ah_h, Ah_h, CS%diag) if (CS%id_div_xx_h>0) call post_data(CS%id_div_xx_h, div_xx_h, CS%diag) @@ -2077,18 +2072,6 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS, MEKE) CS%id_FrictWork = register_diag_field('ocean_model','FrictWork',diag%axesTL,Time,& 'Integral work done by lateral friction terms', 'W m-2', conversion=US%s_to_T**3*US%L_to_m**2) - CS%id_FrictWork_diss = register_diag_field('ocean_model','FrictWork_diss',diag%axesTL,Time,& - 'Integral work done by lateral friction terms (excluding diffusion of energy)', & - 'W m-2', conversion=US%s_to_T**3*US%L_to_m**2) - - if (associated(MEKE)) then - if (associated(MEKE%mom_src)) then - CS%id_FrictWorkMax = register_diag_field('ocean_model', 'FrictWorkMax', diag%axesTL, Time,& - 'Maximum possible integral work done by lateral friction terms', & - 'W m-2', conversion=US%s_to_T**3*US%L_to_m**2) - endif - endif - CS%id_FrictWorkIntz = register_diag_field('ocean_model','FrictWorkIntz',diag%axesT1,Time, & 'Depth integrated work done by lateral friction', 'W m-2', conversion=US%s_to_T**3*US%L_to_m**2, & cmor_field_name='dispkexyfo', & From 992e82652f73190b330eeec83f2375291e7f7c68 Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Tue, 17 Sep 2019 10:23:21 -0600 Subject: [PATCH 040/225] Avoids posting Kv_slow when it is not initialized --- src/parameterizations/vertical/MOM_vert_friction.F90 | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/parameterizations/vertical/MOM_vert_friction.F90 b/src/parameterizations/vertical/MOM_vert_friction.F90 index 1bed36e75e..d286bc815a 100644 --- a/src/parameterizations/vertical/MOM_vert_friction.F90 +++ b/src/parameterizations/vertical/MOM_vert_friction.F90 @@ -1015,7 +1015,8 @@ subroutine vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS, OBC) endif ! Offer diagnostic fields for averaging. - if (CS%id_Kv_slow > 0) call post_data(CS%id_Kv_slow, visc%Kv_slow, CS%diag) + if (associated(visc%Kv_slow) .and. (CS%id_Kv_slow > 0)) & + call post_data(CS%id_Kv_slow, visc%Kv_slow, CS%diag) if (CS%id_Kv_u > 0) call post_data(CS%id_Kv_u, Kv_u, CS%diag) if (CS%id_Kv_v > 0) call post_data(CS%id_Kv_v, Kv_v, CS%diag) if (CS%id_au_vv > 0) call post_data(CS%id_au_vv, CS%a_u, CS%diag) From 353cea1d24689686480b3356f725666d291bc3f3 Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Tue, 17 Sep 2019 13:22:08 -0600 Subject: [PATCH 041/225] remove unused modules --- config_src/mct_driver/mom_ocean_model_mct.F90 | 1 - config_src/nuopc_driver/mom_ocean_model_nuopc.F90 | 1 - 2 files changed, 2 deletions(-) diff --git a/config_src/mct_driver/mom_ocean_model_mct.F90 b/config_src/mct_driver/mom_ocean_model_mct.F90 index ec894f1ebb..3dd302cf2d 100644 --- a/config_src/mct_driver/mom_ocean_model_mct.F90 +++ b/config_src/mct_driver/mom_ocean_model_mct.F90 @@ -56,7 +56,6 @@ module MOM_ocean_model_mct use coupler_types_mod, only : coupler_type_set_diags, coupler_type_send_data use mpp_domains_mod, only : domain2d, mpp_get_layout, mpp_get_global_domain use mpp_domains_mod, only : mpp_define_domains, mpp_get_compute_domain, mpp_get_data_domain -use atmos_ocean_fluxes_mod, only : aof_set_coupler_flux use fms_mod, only : stdout use mpp_mod, only : mpp_chksum use MOM_EOS, only : gsw_sp_from_sr, gsw_pt_from_ct diff --git a/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 b/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 index 426c7e9922..a99c1b60eb 100644 --- a/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 +++ b/config_src/nuopc_driver/mom_ocean_model_nuopc.F90 @@ -52,7 +52,6 @@ module MOM_ocean_model_nuopc use coupler_types_mod, only : coupler_type_set_diags, coupler_type_send_data use mpp_domains_mod, only : domain2d, mpp_get_layout, mpp_get_global_domain use mpp_domains_mod, only : mpp_define_domains, mpp_get_compute_domain, mpp_get_data_domain -use atmos_ocean_fluxes_mod, only : aof_set_coupler_flux use fms_mod, only : stdout use mpp_mod, only : mpp_chksum use MOM_EOS, only : gsw_sp_from_sr, gsw_pt_from_ct From c295ce1bdb5760deaea84667eeb4ee6df98ebb2d Mon Sep 17 00:00:00 2001 From: alperaltuntas Date: Tue, 17 Sep 2019 13:24:27 -0600 Subject: [PATCH 042/225] add ncar version of Doxyfile --- docs/Doxyfile_ncar_rtd | 2437 ++++++++++++++++++++++++++++++++++++++++ 1 file changed, 2437 insertions(+) create mode 100644 docs/Doxyfile_ncar_rtd diff --git a/docs/Doxyfile_ncar_rtd b/docs/Doxyfile_ncar_rtd new file mode 100644 index 0000000000..46e06ef109 --- /dev/null +++ b/docs/Doxyfile_ncar_rtd @@ -0,0 +1,2437 @@ +# Doxyfile 1.8.12 + +# This file describes the settings to be used by the documentation system +# doxygen (www.doxygen.org) for a project. +# +# All text after a double hash (##) is considered a comment and is placed in +# front of the TAG it is preceding. +# +# All text after a single hash (#) is considered a comment and will be ignored. +# The format is: +# TAG = value [value, ...] +# For lists, items can also be appended using: +# TAG += value [value, ...] +# Values that contain spaces should be placed between quotes (\" \"). + +#--------------------------------------------------------------------------- +# Project related configuration options +#--------------------------------------------------------------------------- + +# This tag specifies the encoding used for all characters in the config file +# that follow. The default is UTF-8 which is also the encoding used for all text +# before the first occurrence of this tag. Doxygen uses libiconv (or the iconv +# built into libc) for the transcoding. See http://www.gnu.org/software/libiconv +# for the list of possible encodings. +# The default value is: UTF-8. + +DOXYFILE_ENCODING = UTF-8 + +# The PROJECT_NAME tag is a single word (or a sequence of words surrounded by +# double-quotes, unless you are using Doxywizard) that should identify the +# project for which the documentation is generated. This name is used in the +# title of most generated pages and in a few other places. +# The default value is: My Project. + +PROJECT_NAME = "MOM6" + +# The PROJECT_NUMBER tag can be used to enter a project or revision number. This +# could be handy for archiving the generated documentation or if some version +# control system is used. + +PROJECT_NUMBER = + +# Using the PROJECT_BRIEF tag one can provide an optional one line description +# for a project that appears at the top of each page and should give viewer a +# quick idea about the purpose of the project. Keep the description short. + +PROJECT_BRIEF = + +# With the PROJECT_LOGO tag one can specify a logo or an icon that is included +# in the documentation. The maximum height of the logo should not exceed 55 +# pixels and the maximum width should not exceed 200 pixels. Doxygen will copy +# the logo to the output directory. + +PROJECT_LOGO = + +# The OUTPUT_DIRECTORY tag is used to specify the (relative or absolute) path +# into which the generated documentation will be written. If a relative path is +# entered, it will be relative to the location where doxygen was started. If +# left blank the current directory will be used. + +#OUTPUT_DIRECTORY = + +# If the CREATE_SUBDIRS tag is set to YES then doxygen will create 4096 sub- +# directories (in 2 levels) under the output directory of each output format and +# will distribute the generated files over these directories. Enabling this +# option can be useful when feeding doxygen a huge amount of source files, where +# putting all generated files in the same directory would otherwise causes +# performance problems for the file system. +# The default value is: NO. + +CREATE_SUBDIRS = NO + +# If the ALLOW_UNICODE_NAMES tag is set to YES, doxygen will allow non-ASCII +# characters to appear in the names of generated files. If set to NO, non-ASCII +# characters will be escaped, for example _xE3_x81_x84 will be used for Unicode +# U+3044. +# The default value is: NO. + +ALLOW_UNICODE_NAMES = NO + +# The OUTPUT_LANGUAGE tag is used to specify the language in which all +# documentation generated by doxygen is written. Doxygen will use this +# information to generate all constant output in the proper language. +# Possible values are: Afrikaans, Arabic, Armenian, Brazilian, Catalan, Chinese, +# Chinese-Traditional, Croatian, Czech, Danish, Dutch, English (United States), +# Esperanto, Farsi (Persian), Finnish, French, German, Greek, Hungarian, +# Indonesian, Italian, Japanese, Japanese-en (Japanese with English messages), +# Korean, Korean-en (Korean with English messages), Latvian, Lithuanian, +# Macedonian, Norwegian, Persian (Farsi), Polish, Portuguese, Romanian, Russian, +# Serbian, Serbian-Cyrillic, Slovak, Slovene, Spanish, Swedish, Turkish, +# Ukrainian and Vietnamese. +# The default value is: English. + +OUTPUT_LANGUAGE = English + +# If the BRIEF_MEMBER_DESC tag is set to YES, doxygen will include brief member +# descriptions after the members that are listed in the file and class +# documentation (similar to Javadoc). Set to NO to disable this. +# The default value is: YES. + +BRIEF_MEMBER_DESC = YES + +# If the REPEAT_BRIEF tag is set to YES, doxygen will prepend the brief +# description of a member or function before the detailed description +# +# Note: If both HIDE_UNDOC_MEMBERS and BRIEF_MEMBER_DESC are set to NO, the +# brief descriptions will be completely suppressed. +# The default value is: YES. + +REPEAT_BRIEF = YES + +# This tag implements a quasi-intelligent brief description abbreviator that is +# used to form the text in various listings. Each string in this list, if found +# as the leading text of the brief description, will be stripped from the text +# and the result, after processing the whole list, is used as the annotated +# text. Otherwise, the brief description is used as-is. If left blank, the +# following values are used ($name is automatically replaced with the name of +# the entity):The $name class, The $name widget, The $name file, is, provides, +# specifies, contains, represents, a, an and the. + +ABBREVIATE_BRIEF = + +# If the ALWAYS_DETAILED_SEC and REPEAT_BRIEF tags are both set to YES then +# doxygen will generate a detailed section even if there is only a brief +# description. +# The default value is: NO. + +ALWAYS_DETAILED_SEC = NO + +# If the INLINE_INHERITED_MEMB tag is set to YES, doxygen will show all +# inherited members of a class in the documentation of that class as if those +# members were ordinary class members. Constructors, destructors and assignment +# operators of the base classes will not be shown. +# The default value is: NO. + +INLINE_INHERITED_MEMB = NO + +# If the FULL_PATH_NAMES tag is set to YES, doxygen will prepend the full path +# before files name in the file list and in the header files. If set to NO the +# shortest path that makes the file name unique will be used +# The default value is: YES. + +FULL_PATH_NAMES = YES + +# The STRIP_FROM_PATH tag can be used to strip a user-defined part of the path. +# Stripping is only done if one of the specified strings matches the left-hand +# part of the path. The tag can be used to show relative paths in the file list. +# If left blank the directory from which doxygen is run is used as the path to +# strip. +# +# Note that you can specify absolute paths here, but also relative paths, which +# will be relative from the directory where doxygen is started. +# This tag requires that the tag FULL_PATH_NAMES is set to YES. + +STRIP_FROM_PATH = + +# The STRIP_FROM_INC_PATH tag can be used to strip a user-defined part of the +# path mentioned in the documentation of a class, which tells the reader which +# header file to include in order to use a class. If left blank only the name of +# the header file containing the class definition is used. Otherwise one should +# specify the list of include paths that are normally passed to the compiler +# using the -I flag. + +STRIP_FROM_INC_PATH = + +# If the SHORT_NAMES tag is set to YES, doxygen will generate much shorter (but +# less readable) file names. This can be useful is your file systems doesn't +# support long names like on DOS, Mac, or CD-ROM. +# The default value is: NO. + +SHORT_NAMES = NO + +# If the JAVADOC_AUTOBRIEF tag is set to YES then doxygen will interpret the +# first line (until the first dot) of a Javadoc-style comment as the brief +# description. If set to NO, the Javadoc-style will behave just like regular Qt- +# style comments (thus requiring an explicit @brief command for a brief +# description.) +# The default value is: NO. + +JAVADOC_AUTOBRIEF = NO + +# If the QT_AUTOBRIEF tag is set to YES then doxygen will interpret the first +# line (until the first dot) of a Qt-style comment as the brief description. If +# set to NO, the Qt-style will behave just like regular Qt-style comments (thus +# requiring an explicit \brief command for a brief description.) +# The default value is: NO. + +QT_AUTOBRIEF = NO + +# The MULTILINE_CPP_IS_BRIEF tag can be set to YES to make doxygen treat a +# multi-line C++ special comment block (i.e. a block of //! or /// comments) as +# a brief description. This used to be the default behavior. The new default is +# to treat a multi-line C++ comment block as a detailed description. Set this +# tag to YES if you prefer the old behavior instead. +# +# Note that setting this tag to YES also means that rational rose comments are +# not recognized any more. +# The default value is: NO. + +MULTILINE_CPP_IS_BRIEF = NO + +# If the INHERIT_DOCS tag is set to YES then an undocumented member inherits the +# documentation from any documented member that it re-implements. +# The default value is: YES. + +INHERIT_DOCS = YES + +# If the SEPARATE_MEMBER_PAGES tag is set to YES then doxygen will produce a new +# page for each member. If set to NO, the documentation of a member will be part +# of the file/class/namespace that contains it. +# The default value is: NO. + +SEPARATE_MEMBER_PAGES = NO + +# The TAB_SIZE tag can be used to set the number of spaces in a tab. Doxygen +# uses this value to replace tabs by spaces in code fragments. +# Minimum value: 1, maximum value: 16, default value: 4. + +TAB_SIZE = 2 + +# This tag can be used to specify a number of aliases that act as commands in +# the documentation. An alias has the form: +# name=value +# For example adding +# "sideeffect=@par Side Effects:\n" +# will allow you to put the command \sideeffect (or @sideeffect) in the +# documentation, which will result in a user-defined paragraph with heading +# "Side Effects:". You can put \n's in the value part of an alias to insert +# newlines. + +ALIASES = + +# This tag can be used to specify a number of word-keyword mappings (TCL only). +# A mapping has the form "name=value". For example adding "class=itcl::class" +# will allow you to use the command class in the itcl::class meaning. + +TCL_SUBST = + +# Set the OPTIMIZE_OUTPUT_FOR_C tag to YES if your project consists of C sources +# only. Doxygen will then generate output that is more tailored for C. For +# instance, some of the names that are used will be different. The list of all +# members will be omitted, etc. +# The default value is: NO. + +OPTIMIZE_OUTPUT_FOR_C = NO + +# Set the OPTIMIZE_OUTPUT_JAVA tag to YES if your project consists of Java or +# Python sources only. Doxygen will then generate output that is more tailored +# for that language. For instance, namespaces will be presented as packages, +# qualified scopes will look different, etc. +# The default value is: NO. + +OPTIMIZE_OUTPUT_JAVA = NO + +# Set the OPTIMIZE_FOR_FORTRAN tag to YES if your project consists of Fortran +# sources. Doxygen will then generate output that is tailored for Fortran. +# The default value is: NO. + +OPTIMIZE_FOR_FORTRAN = YES + +# Set the OPTIMIZE_OUTPUT_VHDL tag to YES if your project consists of VHDL +# sources. Doxygen will then generate output that is tailored for VHDL. +# The default value is: NO. + +OPTIMIZE_OUTPUT_VHDL = NO + +# Doxygen selects the parser to use depending on the extension of the files it +# parses. With this tag you can assign which parser to use for a given +# extension. Doxygen has a built-in mapping, but you can override or extend it +# using this tag. The format is ext=language, where ext is a file extension, and +# language is one of the parsers supported by doxygen: IDL, Java, Javascript, +# C#, C, C++, D, PHP, Objective-C, Python, Fortran (fixed format Fortran: +# FortranFixed, free formatted Fortran: FortranFree, unknown formatted Fortran: +# Fortran. In the later case the parser tries to guess whether the code is fixed +# or free formatted code, this is the default for Fortran type files), VHDL. For +# instance to make doxygen treat .inc files as Fortran files (default is PHP), +# and .f files as C (default is Fortran), use: inc=Fortran f=C. +# +# Note: For files without extension you can use no_extension as a placeholder. +# +# Note that for custom extensions you also need to set FILE_PATTERNS otherwise +# the files are not read by doxygen. + +EXTENSION_MAPPING = + +# If the MARKDOWN_SUPPORT tag is enabled then doxygen pre-processes all comments +# according to the Markdown format, which allows for more readable +# documentation. See http://daringfireball.net/projects/markdown/ for details. +# The output of markdown processing is further processed by doxygen, so you can +# mix doxygen, HTML, and XML commands with Markdown formatting. Disable only in +# case of backward compatibilities issues. +# The default value is: YES. + +MARKDOWN_SUPPORT = YES + +# When the TOC_INCLUDE_HEADINGS tag is set to a non-zero value, all headings up +# to that level are automatically included in the table of contents, even if +# they do not have an id attribute. +# Note: This feature currently applies only to Markdown headings. +# Minimum value: 0, maximum value: 99, default value: 0. +# This tag requires that the tag MARKDOWN_SUPPORT is set to YES. + +TOC_INCLUDE_HEADINGS = 0 + +# When enabled doxygen tries to link words that correspond to documented +# classes, or namespaces to their corresponding documentation. Such a link can +# be prevented in individual cases by putting a % sign in front of the word or +# globally by setting AUTOLINK_SUPPORT to NO. +# The default value is: YES. + +AUTOLINK_SUPPORT = YES + +# If you use STL classes (i.e. std::string, std::vector, etc.) but do not want +# to include (a tag file for) the STL sources as input, then you should set this +# tag to YES in order to let doxygen match functions declarations and +# definitions whose arguments contain STL classes (e.g. func(std::string); +# versus func(std::string) {}). This also make the inheritance and collaboration +# diagrams that involve STL classes more complete and accurate. +# The default value is: NO. + +BUILTIN_STL_SUPPORT = NO + +# If you use Microsoft's C++/CLI language, you should set this option to YES to +# enable parsing support. +# The default value is: NO. + +CPP_CLI_SUPPORT = NO + +# Set the SIP_SUPPORT tag to YES if your project consists of sip (see: +# http://www.riverbankcomputing.co.uk/software/sip/intro) sources only. Doxygen +# will parse them like normal C++ but will assume all classes use public instead +# of private inheritance when no explicit protection keyword is present. +# The default value is: NO. + +SIP_SUPPORT = NO + +# For Microsoft's IDL there are propget and propput attributes to indicate +# getter and setter methods for a property. Setting this option to YES will make +# doxygen to replace the get and set methods by a property in the documentation. +# This will only work if the methods are indeed getting or setting a simple +# type. If this is not the case, or you want to show the methods anyway, you +# should set this option to NO. +# The default value is: YES. + +IDL_PROPERTY_SUPPORT = YES + +# If member grouping is used in the documentation and the DISTRIBUTE_GROUP_DOC +# tag is set to YES then doxygen will reuse the documentation of the first +# member in the group (if any) for the other members of the group. By default +# all members of a group must be documented explicitly. +# The default value is: NO. + +DISTRIBUTE_GROUP_DOC = YES + +# If one adds a struct or class to a group and this option is enabled, then also +# any nested class or struct is added to the same group. By default this option +# is disabled and one has to add nested compounds explicitly via \ingroup. +# The default value is: NO. + +GROUP_NESTED_COMPOUNDS = NO + +# Set the SUBGROUPING tag to YES to allow class member groups of the same type +# (for instance a group of public functions) to be put as a subgroup of that +# type (e.g. under the Public Functions section). Set it to NO to prevent +# subgrouping. Alternatively, this can be done per class using the +# \nosubgrouping command. +# The default value is: YES. + +SUBGROUPING = YES + +# When the INLINE_GROUPED_CLASSES tag is set to YES, classes, structs and unions +# are shown inside the group in which they are included (e.g. using \ingroup) +# instead of on a separate page (for HTML and Man pages) or section (for LaTeX +# and RTF). +# +# Note that this feature does not work in combination with +# SEPARATE_MEMBER_PAGES. +# The default value is: NO. + +INLINE_GROUPED_CLASSES = NO + +# When the INLINE_SIMPLE_STRUCTS tag is set to YES, structs, classes, and unions +# with only public data fields or simple typedef fields will be shown inline in +# the documentation of the scope in which they are defined (i.e. file, +# namespace, or group documentation), provided this scope is documented. If set +# to NO, structs, classes, and unions are shown on a separate page (for HTML and +# Man pages) or section (for LaTeX and RTF). +# The default value is: NO. + +INLINE_SIMPLE_STRUCTS = NO + +# When TYPEDEF_HIDES_STRUCT tag is enabled, a typedef of a struct, union, or +# enum is documented as struct, union, or enum with the name of the typedef. So +# typedef struct TypeS {} TypeT, will appear in the documentation as a struct +# with name TypeT. When disabled the typedef will appear as a member of a file, +# namespace, or class. And the struct will be named TypeS. This can typically be +# useful for C code in case the coding convention dictates that all compound +# types are typedef'ed and only the typedef is referenced, never the tag name. +# The default value is: NO. + +TYPEDEF_HIDES_STRUCT = NO + +# The size of the symbol lookup cache can be set using LOOKUP_CACHE_SIZE. This +# cache is used to resolve symbols given their name and scope. Since this can be +# an expensive process and often the same symbol appears multiple times in the +# code, doxygen keeps a cache of pre-resolved symbols. If the cache is too small +# doxygen will become slower. If the cache is too large, memory is wasted. The +# cache size is given by this formula: 2^(16+LOOKUP_CACHE_SIZE). The valid range +# is 0..9, the default is 0, corresponding to a cache size of 2^16=65536 +# symbols. At the end of a run doxygen will report the cache usage and suggest +# the optimal cache size from a speed point of view. +# Minimum value: 0, maximum value: 9, default value: 0. + +LOOKUP_CACHE_SIZE = 0 + +#--------------------------------------------------------------------------- +# Build related configuration options +#--------------------------------------------------------------------------- + +# If the EXTRACT_ALL tag is set to YES, doxygen will assume all entities in +# documentation are documented, even if no documentation was available. Private +# class members and static file members will be hidden unless the +# EXTRACT_PRIVATE respectively EXTRACT_STATIC tags are set to YES. +# Note: This will also disable the warnings about undocumented members that are +# normally produced when WARNINGS is set to YES. +# The default value is: NO. + +EXTRACT_ALL = YES + +# If the EXTRACT_PRIVATE tag is set to YES, all private members of a class will +# be included in the documentation. +# The default value is: NO. + +EXTRACT_PRIVATE = YES + +# If the EXTRACT_PACKAGE tag is set to YES, all members with package or internal +# scope will be included in the documentation. +# The default value is: NO. + +EXTRACT_PACKAGE = YES + +# If the EXTRACT_STATIC tag is set to YES, all static members of a file will be +# included in the documentation. +# The default value is: NO. + +EXTRACT_STATIC = YES + +# If the EXTRACT_LOCAL_CLASSES tag is set to YES, classes (and structs) defined +# locally in source files will be included in the documentation. If set to NO, +# only classes defined in header files are included. Does not have any effect +# for Java sources. +# The default value is: YES. + +EXTRACT_LOCAL_CLASSES = YES + +# This flag is only useful for Objective-C code. If set to YES, local methods, +# which are defined in the implementation section but not in the interface are +# included in the documentation. If set to NO, only methods in the interface are +# included. +# The default value is: NO. + +EXTRACT_LOCAL_METHODS = YES + +# If this flag is set to YES, the members of anonymous namespaces will be +# extracted and appear in the documentation as a namespace called +# 'anonymous_namespace{file}', where file will be replaced with the base name of +# the file that contains the anonymous namespace. By default anonymous namespace +# are hidden. +# The default value is: NO. + +EXTRACT_ANON_NSPACES = YES + +# If the HIDE_UNDOC_MEMBERS tag is set to YES, doxygen will hide all +# undocumented members inside documented classes or files. If set to NO these +# members will be included in the various overviews, but no documentation +# section is generated. This option has no effect if EXTRACT_ALL is enabled. +# The default value is: NO. + +HIDE_UNDOC_MEMBERS = NO + +# If the HIDE_UNDOC_CLASSES tag is set to YES, doxygen will hide all +# undocumented classes that are normally visible in the class hierarchy. If set +# to NO, these classes will be included in the various overviews. This option +# has no effect if EXTRACT_ALL is enabled. +# The default value is: NO. + +HIDE_UNDOC_CLASSES = NO + +# If the HIDE_FRIEND_COMPOUNDS tag is set to YES, doxygen will hide all friend +# (class|struct|union) declarations. If set to NO, these declarations will be +# included in the documentation. +# The default value is: NO. + +HIDE_FRIEND_COMPOUNDS = NO + +# If the HIDE_IN_BODY_DOCS tag is set to YES, doxygen will hide any +# documentation blocks found inside the body of a function. If set to NO, these +# blocks will be appended to the function's detailed documentation block. +# The default value is: NO. + +HIDE_IN_BODY_DOCS = NO + +# The INTERNAL_DOCS tag determines if documentation that is typed after a +# \internal command is included. If the tag is set to NO then the documentation +# will be excluded. Set it to YES to include the internal documentation. +# The default value is: NO. + +INTERNAL_DOCS = YES + +# If the CASE_SENSE_NAMES tag is set to NO then doxygen will only generate file +# names in lower-case letters. If set to YES, upper-case letters are also +# allowed. This is useful if you have classes or files whose names only differ +# in case and if your file system supports case sensitive file names. Windows +# and Mac users are advised to set this option to NO. +# The default value is: system dependent. + +CASE_SENSE_NAMES = YES + +# If the HIDE_SCOPE_NAMES tag is set to NO then doxygen will show members with +# their full class and namespace scopes in the documentation. If set to YES, the +# scope will be hidden. +# The default value is: NO. + +HIDE_SCOPE_NAMES = NO + +# If the HIDE_COMPOUND_REFERENCE tag is set to NO (default) then doxygen will +# append additional text to a page's title, such as Class Reference. If set to +# YES the compound reference will be hidden. +# The default value is: NO. + +HIDE_COMPOUND_REFERENCE= NO + +# If the SHOW_INCLUDE_FILES tag is set to YES then doxygen will put a list of +# the files that are included by a file in the documentation of that file. +# The default value is: YES. + +SHOW_INCLUDE_FILES = YES + +# If the SHOW_GROUPED_MEMB_INC tag is set to YES then Doxygen will add for each +# grouped member an include statement to the documentation, telling the reader +# which file to include in order to use the member. +# The default value is: NO. + +SHOW_GROUPED_MEMB_INC = NO + +# If the FORCE_LOCAL_INCLUDES tag is set to YES then doxygen will list include +# files with double quotes in the documentation rather than with sharp brackets. +# The default value is: NO. + +FORCE_LOCAL_INCLUDES = NO + +# If the INLINE_INFO tag is set to YES then a tag [inline] is inserted in the +# documentation for inline members. +# The default value is: YES. + +INLINE_INFO = YES + +# If the SORT_MEMBER_DOCS tag is set to YES then doxygen will sort the +# (detailed) documentation of file and class members alphabetically by member +# name. If set to NO, the members will appear in declaration order. +# The default value is: YES. + +SORT_MEMBER_DOCS = YES + +# If the SORT_BRIEF_DOCS tag is set to YES then doxygen will sort the brief +# descriptions of file, namespace and class members alphabetically by member +# name. If set to NO, the members will appear in declaration order. Note that +# this will also influence the order of the classes in the class list. +# The default value is: NO. + +SORT_BRIEF_DOCS = NO + +# If the SORT_MEMBERS_CTORS_1ST tag is set to YES then doxygen will sort the +# (brief and detailed) documentation of class members so that constructors and +# destructors are listed first. If set to NO the constructors will appear in the +# respective orders defined by SORT_BRIEF_DOCS and SORT_MEMBER_DOCS. +# Note: If SORT_BRIEF_DOCS is set to NO this option is ignored for sorting brief +# member documentation. +# Note: If SORT_MEMBER_DOCS is set to NO this option is ignored for sorting +# detailed member documentation. +# The default value is: NO. + +SORT_MEMBERS_CTORS_1ST = NO + +# If the SORT_GROUP_NAMES tag is set to YES then doxygen will sort the hierarchy +# of group names into alphabetical order. If set to NO the group names will +# appear in their defined order. +# The default value is: NO. + +SORT_GROUP_NAMES = NO + +# If the SORT_BY_SCOPE_NAME tag is set to YES, the class list will be sorted by +# fully-qualified names, including namespaces. If set to NO, the class list will +# be sorted only by class name, not including the namespace part. +# Note: This option is not very useful if HIDE_SCOPE_NAMES is set to YES. +# Note: This option applies only to the class list, not to the alphabetical +# list. +# The default value is: NO. + +SORT_BY_SCOPE_NAME = NO + +# If the STRICT_PROTO_MATCHING option is enabled and doxygen fails to do proper +# type resolution of all parameters of a function it will reject a match between +# the prototype and the implementation of a member function even if there is +# only one candidate or it is obvious which candidate to choose by doing a +# simple string match. By disabling STRICT_PROTO_MATCHING doxygen will still +# accept a match between prototype and implementation in such cases. +# The default value is: NO. + +STRICT_PROTO_MATCHING = NO + +# The GENERATE_TODOLIST tag can be used to enable (YES) or disable (NO) the todo +# list. This list is created by putting \todo commands in the documentation. +# The default value is: YES. + +GENERATE_TODOLIST = NO + +# The GENERATE_TESTLIST tag can be used to enable (YES) or disable (NO) the test +# list. This list is created by putting \test commands in the documentation. +# The default value is: YES. + +GENERATE_TESTLIST = YES + +# The GENERATE_BUGLIST tag can be used to enable (YES) or disable (NO) the bug +# list. This list is created by putting \bug commands in the documentation. +# The default value is: YES. + +GENERATE_BUGLIST = YES + +# The GENERATE_DEPRECATEDLIST tag can be used to enable (YES) or disable (NO) +# the deprecated list. This list is created by putting \deprecated commands in +# the documentation. +# The default value is: YES. + +GENERATE_DEPRECATEDLIST= YES + +# The ENABLED_SECTIONS tag can be used to enable conditional documentation +# sections, marked by \if ... \endif and \cond +# ... \endcond blocks. + +ENABLED_SECTIONS = + +# The MAX_INITIALIZER_LINES tag determines the maximum number of lines that the +# initial value of a variable or macro / define can have for it to appear in the +# documentation. If the initializer consists of more lines than specified here +# it will be hidden. Use a value of 0 to hide initializers completely. The +# appearance of the value of individual variables and macros / defines can be +# controlled using \showinitializer or \hideinitializer command in the +# documentation regardless of this setting. +# Minimum value: 0, maximum value: 10000, default value: 30. + +MAX_INITIALIZER_LINES = 30 + +# Set the SHOW_USED_FILES tag to NO to disable the list of files generated at +# the bottom of the documentation of classes and structs. If set to YES, the +# list will mention the files that were used to generate the documentation. +# The default value is: YES. + +SHOW_USED_FILES = YES + +# Set the SHOW_FILES tag to NO to disable the generation of the Files page. This +# will remove the Files entry from the Quick Index and from the Folder Tree View +# (if specified). +# The default value is: YES. + +SHOW_FILES = YES + +# Set the SHOW_NAMESPACES tag to NO to disable the generation of the Namespaces +# page. This will remove the Namespaces entry from the Quick Index and from the +# Folder Tree View (if specified). +# The default value is: YES. + +SHOW_NAMESPACES = YES + +# The FILE_VERSION_FILTER tag can be used to specify a program or script that +# doxygen should invoke to get the current version for each file (typically from +# the version control system). Doxygen will invoke the program by executing (via +# popen()) the command command input-file, where command is the value of the +# FILE_VERSION_FILTER tag, and input-file is the name of an input file provided +# by doxygen. Whatever the program writes to standard output is used as the file +# version. For an example see the documentation. + +FILE_VERSION_FILTER = + +# The LAYOUT_FILE tag can be used to specify a layout file which will be parsed +# by doxygen. The layout file controls the global structure of the generated +# output files in an output format independent way. To create the layout file +# that represents doxygen's defaults, run doxygen with the -l option. You can +# optionally specify a file name after the option, if omitted DoxygenLayout.xml +# will be used as the name of the layout file. +# +# Note that if you run doxygen from a directory containing a file called +# DoxygenLayout.xml, doxygen will parse it automatically even if the LAYOUT_FILE +# tag is left empty. + +LAYOUT_FILE = + +# The CITE_BIB_FILES tag can be used to specify one or more bib files containing +# the reference definitions. This must be a list of .bib files. The .bib +# extension is automatically appended if omitted. This requires the bibtex tool +# to be installed. See also http://en.wikipedia.org/wiki/BibTeX for more info. +# For LaTeX the style of the bibliography can be controlled using +# LATEX_BIB_STYLE. To use this feature you need bibtex and perl available in the +# search path. See also \cite for info how to create references. + +CITE_BIB_FILES = + +#--------------------------------------------------------------------------- +# Configuration options related to warning and progress messages +#--------------------------------------------------------------------------- + +# The QUIET tag can be used to turn on/off the messages that are generated to +# standard output by doxygen. If QUIET is set to YES this implies that the +# messages are off. +# The default value is: NO. + +QUIET = NO + +# The WARNINGS tag can be used to turn on/off the warning messages that are +# generated to standard error (stderr) by doxygen. If WARNINGS is set to YES +# this implies that the warnings are on. +# +# Tip: Turn warnings on while writing the documentation. +# The default value is: YES. + +WARNINGS = YES + +# If the WARN_IF_UNDOCUMENTED tag is set to YES then doxygen will generate +# warnings for undocumented members. If EXTRACT_ALL is set to YES then this flag +# will automatically be disabled. +# The default value is: YES. + +WARN_IF_UNDOCUMENTED = YES + +# If the WARN_IF_DOC_ERROR tag is set to YES, doxygen will generate warnings for +# potential errors in the documentation, such as not documenting some parameters +# in a documented function, or documenting parameters that don't exist or using +# markup commands wrongly. +# The default value is: YES. + +WARN_IF_DOC_ERROR = YES + +# This WARN_NO_PARAMDOC option can be enabled to get warnings for functions that +# are documented, but have no documentation for their parameters or return +# value. If set to NO, doxygen will only warn about wrong or incomplete +# parameter documentation, but not about the absence of documentation. +# The default value is: NO. + +WARN_NO_PARAMDOC = NO + +# If the WARN_AS_ERROR tag is set to YES then doxygen will immediately stop when +# a warning is encountered. +# The default value is: NO. + +WARN_AS_ERROR = NO + +# The WARN_FORMAT tag determines the format of the warning messages that doxygen +# can produce. The string should contain the $file, $line, and $text tags, which +# will be replaced by the file and line number from which the warning originated +# and the warning text. Optionally the format may contain $version, which will +# be replaced by the version of the file (if it could be obtained via +# FILE_VERSION_FILTER) +# The default value is: $file:$line: $text. + +WARN_FORMAT = "$file:$line: $text" + +# The WARN_LOGFILE tag can be used to specify a file to which warning and error +# messages should be written. If left blank the output is written to standard +# error (stderr). + +WARN_LOGFILE = doxygen.log + +#--------------------------------------------------------------------------- +# Configuration options related to the input files +#--------------------------------------------------------------------------- + +# The INPUT tag is used to specify the files and/or directories that contain +# documented source files. You may enter file names like myfile.cpp or +# directories like /usr/src/myproject. Separate the files or directories with +# spaces. +# Note: If this tag is empty the current directory is searched. + +INPUT = ../src \ + ../config_src/solo_driver \ + ../config_src/mct_driver \ + ../config_src/nuopc_driver \ + +# This tag can be used to specify the character encoding of the source files +# that doxygen parses. Internally doxygen uses the UTF-8 encoding. Doxygen uses +# libiconv (or the iconv built into libc) for the transcoding. See the libiconv +# documentation (see: http://www.gnu.org/software/libiconv) for the list of +# possible encodings. +# The default value is: UTF-8. + +INPUT_ENCODING = UTF-8 + +# If the value of the INPUT tag contains directories, you can use the +# FILE_PATTERNS tag to specify one or more wildcard patterns (like *.cpp and +# *.h) to filter out the source-files in the directories. +# +# Note that for custom extensions or not directly supported extensions you also +# need to set EXTENSION_MAPPING for the extension otherwise the files are not +# read by doxygen. +# +# If left blank the following patterns are tested:*.c, *.cc, *.cxx, *.cpp, +# *.c++, *.java, *.ii, *.ixx, *.ipp, *.i++, *.inl, *.idl, *.ddl, *.odl, *.h, +# *.hh, *.hxx, *.hpp, *.h++, *.cs, *.d, *.php, *.php4, *.php5, *.phtml, *.inc, +# *.m, *.markdown, *.md, *.mm, *.dox, *.py, *.pyw, *.f90, *.f, *.for, *.tcl, +# *.vhd, *.vhdl, *.ucf and *.qsf. + +FILE_PATTERNS = *.c \ + *.cc \ + *.cxx \ + *.cpp \ + *.c++ \ + *.h \ + *.hh \ + *.hxx \ + *.hpp \ + *.h++ \ + *.inc \ + *.m \ + *.markdown \ + *.md \ + *.mm \ + *.dox \ + *.f90 \ + *.f \ + *.for \ + *.F90 + +# The RECURSIVE tag can be used to specify whether or not subdirectories should +# be searched for input files as well. +# The default value is: NO. + +RECURSIVE = YES + +# The EXCLUDE tag can be used to specify files and/or directories that should be +# excluded from the INPUT source files. This way you can easily exclude a +# subdirectory from a directory tree whose root is specified with the INPUT tag. +# +# Note that relative paths are relative to the directory from which doxygen is +# run. + +EXCLUDE = ../src/equation_of_state/TEOS10 + +# The EXCLUDE_SYMLINKS tag can be used to select whether or not files or +# directories that are symbolic links (a Unix file system feature) are excluded +# from the input. +# The default value is: NO. + +EXCLUDE_SYMLINKS = NO + +# If the value of the INPUT tag contains directories, you can use the +# EXCLUDE_PATTERNS tag to specify one or more wildcard patterns to exclude +# certain files from those directories. +# +# Note that the wildcards are matched against the file with absolute path, so to +# exclude all test directories for example use the pattern */test/* + +EXCLUDE_PATTERNS = makedep.py Makefile INSTALL + +# The EXCLUDE_SYMBOLS tag can be used to specify one or more symbol names +# (namespaces, classes, functions, etc.) that should be excluded from the +# output. The symbol name can be a fully qualified name, a word, or if the +# wildcard * is used, a substring. Examples: ANamespace, AClass, +# AClass::ANamespace, ANamespace::*Test +# +# Note that the wildcards are matched against the file with absolute path, so to +# exclude all test directories use the pattern */test/* + +EXCLUDE_SYMBOLS = + +# The EXAMPLE_PATH tag can be used to specify one or more files or directories +# that contain example code fragments that are included (see the \include +# command). + +EXAMPLE_PATH = ../src + +# If the value of the EXAMPLE_PATH tag contains directories, you can use the +# EXAMPLE_PATTERNS tag to specify one or more wildcard pattern (like *.cpp and +# *.h) to filter out the source-files in the directories. If left blank all +# files are included. + +EXAMPLE_PATTERNS = * + +# If the EXAMPLE_RECURSIVE tag is set to YES then subdirectories will be +# searched for input files to be used with the \include or \dontinclude commands +# irrespective of the value of the RECURSIVE tag. +# The default value is: NO. + +EXAMPLE_RECURSIVE = NO + +# The IMAGE_PATH tag can be used to specify one or more files or directories +# that contain images that are to be included in the documentation (see the +# \image command). + +IMAGE_PATH = images ../src + +# The INPUT_FILTER tag can be used to specify a program that doxygen should +# invoke to filter for each input file. Doxygen will invoke the filter program +# by executing (via popen()) the command: +# +# +# +# where is the value of the INPUT_FILTER tag, and is the +# name of an input file. Doxygen will then use the output that the filter +# program writes to standard output. If FILTER_PATTERNS is specified, this tag +# will be ignored. +# +# Note that the filter must not add or remove lines; it is applied before the +# code is scanned, but not when the output code is generated. If lines are added +# or removed, the anchors will not be placed correctly. +# +# Note that for custom extensions or not directly supported extensions you also +# need to set EXTENSION_MAPPING for the extension otherwise the files are not +# properly processed by doxygen. + +INPUT_FILTER = + +# The FILTER_PATTERNS tag can be used to specify filters on a per file pattern +# basis. Doxygen will compare the file name with each pattern and apply the +# filter if there is a match. The filters are a list of the form: pattern=filter +# (like *.cpp=my_cpp_filter). See INPUT_FILTER for further information on how +# filters are used. If the FILTER_PATTERNS tag is empty or if none of the +# patterns match the file name, INPUT_FILTER is applied. +# +# Note that for custom extensions or not directly supported extensions you also +# need to set EXTENSION_MAPPING for the extension otherwise the files are not +# properly processed by doxygen. + +FILTER_PATTERNS = + +# If the FILTER_SOURCE_FILES tag is set to YES, the input filter (if set using +# INPUT_FILTER) will also be used to filter the input files that are used for +# producing the source files to browse (i.e. when SOURCE_BROWSER is set to YES). +# The default value is: NO. + +FILTER_SOURCE_FILES = NO + +# The FILTER_SOURCE_PATTERNS tag can be used to specify source filters per file +# pattern. A pattern will override the setting for FILTER_PATTERN (if any) and +# it is also possible to disable source filtering for a specific pattern using +# *.ext= (so without naming a filter). +# This tag requires that the tag FILTER_SOURCE_FILES is set to YES. + +FILTER_SOURCE_PATTERNS = + +# If the USE_MDFILE_AS_MAINPAGE tag refers to the name of a markdown file that +# is part of the input, its contents will be placed on the main page +# (index.html). This can be useful if you have a project on for instance GitHub +# and want to reuse the introduction page also for the doxygen output. + +USE_MDFILE_AS_MAINPAGE = ../README.md + +#--------------------------------------------------------------------------- +# Configuration options related to source browsing +#--------------------------------------------------------------------------- + +# If the SOURCE_BROWSER tag is set to YES then a list of source files will be +# generated. Documented entities will be cross-referenced with these sources. +# +# Note: To get rid of all source code in the generated output, make sure that +# also VERBATIM_HEADERS is set to NO. +# The default value is: NO. + +SOURCE_BROWSER = YES + +# Setting the INLINE_SOURCES tag to YES will include the body of functions, +# classes and enums directly into the documentation. +# The default value is: NO. + +INLINE_SOURCES = YES + +# Setting the STRIP_CODE_COMMENTS tag to YES will instruct doxygen to hide any +# special comment blocks from generated source code fragments. Normal C, C++ and +# Fortran comments will always remain visible. +# The default value is: YES. + +STRIP_CODE_COMMENTS = NO + +# If the REFERENCED_BY_RELATION tag is set to YES then for each documented +# function all documented functions referencing it will be listed. +# The default value is: NO. + +REFERENCED_BY_RELATION = YES + +# If the REFERENCES_RELATION tag is set to YES then for each documented function +# all documented entities called/used by that function will be listed. +# The default value is: NO. + +REFERENCES_RELATION = YES + +# If the REFERENCES_LINK_SOURCE tag is set to YES and SOURCE_BROWSER tag is set +# to YES then the hyperlinks from functions in REFERENCES_RELATION and +# REFERENCED_BY_RELATION lists will link to the source code. Otherwise they will +# link to the documentation. +# The default value is: YES. + +REFERENCES_LINK_SOURCE = NO + +# If SOURCE_TOOLTIPS is enabled (the default) then hovering a hyperlink in the +# source code will show a tooltip with additional information such as prototype, +# brief description and links to the definition and documentation. Since this +# will make the HTML file larger and loading of large files a bit slower, you +# can opt to disable this feature. +# The default value is: YES. +# This tag requires that the tag SOURCE_BROWSER is set to YES. + +SOURCE_TOOLTIPS = YES + +# If the USE_HTAGS tag is set to YES then the references to source code will +# point to the HTML generated by the htags(1) tool instead of doxygen built-in +# source browser. The htags tool is part of GNU's global source tagging system +# (see http://www.gnu.org/software/global/global.html). You will need version +# 4.8.6 or higher. +# +# To use it do the following: +# - Install the latest version of global +# - Enable SOURCE_BROWSER and USE_HTAGS in the config file +# - Make sure the INPUT points to the root of the source tree +# - Run doxygen as normal +# +# Doxygen will invoke htags (and that will in turn invoke gtags), so these +# tools must be available from the command line (i.e. in the search path). +# +# The result: instead of the source browser generated by doxygen, the links to +# source code will now point to the output of htags. +# The default value is: NO. +# This tag requires that the tag SOURCE_BROWSER is set to YES. + +USE_HTAGS = NO + +# If the VERBATIM_HEADERS tag is set the YES then doxygen will generate a +# verbatim copy of the header file for each class for which an include is +# specified. Set to NO to disable this. +# See also: Section \class. +# The default value is: YES. + +VERBATIM_HEADERS = YES + +#--------------------------------------------------------------------------- +# Configuration options related to the alphabetical class index +#--------------------------------------------------------------------------- + +# If the ALPHABETICAL_INDEX tag is set to YES, an alphabetical index of all +# compounds will be generated. Enable this if the project contains a lot of +# classes, structs, unions or interfaces. +# The default value is: YES. + +ALPHABETICAL_INDEX = YES + +# The COLS_IN_ALPHA_INDEX tag can be used to specify the number of columns in +# which the alphabetical index list will be split. +# Minimum value: 1, maximum value: 20, default value: 5. +# This tag requires that the tag ALPHABETICAL_INDEX is set to YES. + +COLS_IN_ALPHA_INDEX = 1 + +# In case all classes in a project start with a common prefix, all classes will +# be put under the same header in the alphabetical index. The IGNORE_PREFIX tag +# can be used to specify a prefix (or a list of prefixes) that should be ignored +# while generating the index headers. +# This tag requires that the tag ALPHABETICAL_INDEX is set to YES. + +IGNORE_PREFIX = + +#--------------------------------------------------------------------------- +# Configuration options related to the HTML output +#--------------------------------------------------------------------------- + +# If the GENERATE_HTML tag is set to YES, doxygen will generate HTML output +# The default value is: YES. + +GENERATE_HTML = NO + +# The HTML_OUTPUT tag is used to specify where the HTML docs will be put. If a +# relative path is entered the value of OUTPUT_DIRECTORY will be put in front of +# it. +# The default directory is: html. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_OUTPUT = APIs + +# The HTML_FILE_EXTENSION tag can be used to specify the file extension for each +# generated HTML page (for example: .htm, .php, .asp). +# The default value is: .html. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_FILE_EXTENSION = .html + +# The HTML_HEADER tag can be used to specify a user-defined HTML header file for +# each generated HTML page. If the tag is left blank doxygen will generate a +# standard header. +# +# To get valid HTML the header file that includes any scripts and style sheets +# that doxygen needs, which is dependent on the configuration options used (e.g. +# the setting GENERATE_TREEVIEW). It is highly recommended to start with a +# default header using +# doxygen -w html new_header.html new_footer.html new_stylesheet.css +# YourConfigFile +# and then modify the file new_header.html. See also section "Doxygen usage" +# for information on how to generate the default header that doxygen normally +# uses. +# Note: The header is subject to change so you typically have to regenerate the +# default header when upgrading to a newer version of doxygen. For a description +# of the possible markers and block names see the documentation. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_HEADER = dox_rtd_header.html + +# The HTML_FOOTER tag can be used to specify a user-defined HTML footer for each +# generated HTML page. If the tag is left blank doxygen will generate a standard +# footer. See HTML_HEADER for more information on how to generate a default +# footer and what special commands can be used inside the footer. See also +# section "Doxygen usage" for information on how to generate the default footer +# that doxygen normally uses. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_FOOTER = dox_rtd_footer.html + +# The HTML_STYLESHEET tag can be used to specify a user-defined cascading style +# sheet that is used by each HTML page. It can be used to fine-tune the look of +# the HTML output. If left blank doxygen will generate a default style sheet. +# See also section "Doxygen usage" for information on how to generate the style +# sheet that doxygen normally uses. +# Note: It is recommended to use HTML_EXTRA_STYLESHEET instead of this tag, as +# it is more robust and this tag (HTML_STYLESHEET) will in the future become +# obsolete. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_STYLESHEET = + +# The HTML_EXTRA_STYLESHEET tag can be used to specify additional user-defined +# cascading style sheets that are included after the standard style sheets +# created by doxygen. Using this option one can overrule certain style aspects. +# This is preferred over using HTML_STYLESHEET since it does not replace the +# standard style sheet and is therefore more robust against future updates. +# Doxygen will copy the style sheet files to the output directory. +# Note: The order of the extra style sheet files is of importance (e.g. the last +# style sheet in the list overrules the setting of the previous ones in the +# list). For an example see the documentation. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_EXTRA_STYLESHEET = + +# The HTML_EXTRA_FILES tag can be used to specify one or more extra images or +# other source files which should be copied to the HTML output directory. Note +# that these files will be copied to the base HTML output directory. Use the +# $relpath^ marker in the HTML_HEADER and/or HTML_FOOTER files to load these +# files. In the HTML_STYLESHEET file, use the file name only. Also note that the +# files will be copied as-is; there are no commands or markers available. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_EXTRA_FILES = + +# The HTML_COLORSTYLE_HUE tag controls the color of the HTML output. Doxygen +# will adjust the colors in the style sheet and background images according to +# this color. Hue is specified as an angle on a colorwheel, see +# http://en.wikipedia.org/wiki/Hue for more information. For instance the value +# 0 represents red, 60 is yellow, 120 is green, 180 is cyan, 240 is blue, 300 +# purple, and 360 is red again. +# Minimum value: 0, maximum value: 359, default value: 220. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_COLORSTYLE_HUE = 220 + +# The HTML_COLORSTYLE_SAT tag controls the purity (or saturation) of the colors +# in the HTML output. For a value of 0 the output will use grayscales only. A +# value of 255 will produce the most vivid colors. +# Minimum value: 0, maximum value: 255, default value: 100. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_COLORSTYLE_SAT = 100 + +# The HTML_COLORSTYLE_GAMMA tag controls the gamma correction applied to the +# luminance component of the colors in the HTML output. Values below 100 +# gradually make the output lighter, whereas values above 100 make the output +# darker. The value divided by 100 is the actual gamma applied, so 80 represents +# a gamma of 0.8, The value 220 represents a gamma of 2.2, and 100 does not +# change the gamma. +# Minimum value: 40, maximum value: 240, default value: 80. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_COLORSTYLE_GAMMA = 80 + +# If the HTML_TIMESTAMP tag is set to YES then the footer of each generated HTML +# page will contain the date and time when the page was generated. Setting this +# to YES can help to show when doxygen was last run and thus if the +# to NO can help when comparing the output of multiple runs. +# The default value is: YES. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_TIMESTAMP = NO + +# If the HTML_DYNAMIC_SECTIONS tag is set to YES then the generated HTML +# documentation will contain sections that can be hidden and shown after the +# page has loaded. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_DYNAMIC_SECTIONS = NO + +# With HTML_INDEX_NUM_ENTRIES one can control the preferred number of entries +# shown in the various tree structured indices initially; the user can expand +# and collapse entries dynamically later on. Doxygen will expand the tree to +# such a level that at most the specified number of entries are visible (unless +# a fully collapsed tree already exceeds this amount). So setting the number of +# entries 1 will produce a full collapsed tree by default. 0 is a special value +# representing an infinite number of entries and will result in a full expanded +# tree by default. +# Minimum value: 0, maximum value: 9999, default value: 100. +# This tag requires that the tag GENERATE_HTML is set to YES. + +HTML_INDEX_NUM_ENTRIES = 900 + +# If the GENERATE_DOCSET tag is set to YES, additional index files will be +# generated that can be used as input for Apple's Xcode 3 integrated development +# environment (see: http://developer.apple.com/tools/xcode/), introduced with +# OSX 10.5 (Leopard). To create a documentation set, doxygen will generate a +# Makefile in the HTML output directory. Running make will produce the docset in +# that directory and running make install will install the docset in +# ~/Library/Developer/Shared/Documentation/DocSets so that Xcode will find it at +# startup. See http://developer.apple.com/tools/creatingdocsetswithdoxygen.html +# for more information. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTML is set to YES. + +GENERATE_DOCSET = NO + +# This tag determines the name of the docset feed. A documentation feed provides +# an umbrella under which multiple documentation sets from a single provider +# (such as a company or product suite) can be grouped. +# The default value is: Doxygen generated docs. +# This tag requires that the tag GENERATE_DOCSET is set to YES. + +DOCSET_FEEDNAME = "Doxygen generated docs" + +# This tag specifies a string that should uniquely identify the documentation +# set bundle. This should be a reverse domain-name style string, e.g. +# com.mycompany.MyDocSet. Doxygen will append .docset to the name. +# The default value is: org.doxygen.Project. +# This tag requires that the tag GENERATE_DOCSET is set to YES. + +DOCSET_BUNDLE_ID = org.doxygen.Project + +# The DOCSET_PUBLISHER_ID tag specifies a string that should uniquely identify +# the documentation publisher. This should be a reverse domain-name style +# string, e.g. com.mycompany.MyDocSet.documentation. +# The default value is: org.doxygen.Publisher. +# This tag requires that the tag GENERATE_DOCSET is set to YES. + +DOCSET_PUBLISHER_ID = org.doxygen.Publisher + +# The DOCSET_PUBLISHER_NAME tag identifies the documentation publisher. +# The default value is: Publisher. +# This tag requires that the tag GENERATE_DOCSET is set to YES. + +DOCSET_PUBLISHER_NAME = Publisher + +# If the GENERATE_HTMLHELP tag is set to YES then doxygen generates three +# additional HTML index files: index.hhp, index.hhc, and index.hhk. The +# index.hhp is a project file that can be read by Microsoft's HTML Help Workshop +# (see: http://www.microsoft.com/en-us/download/details.aspx?id=21138) on +# Windows. +# +# The HTML Help Workshop contains a compiler that can convert all HTML output +# generated by doxygen into a single compiled HTML file (.chm). Compiled HTML +# files are now used as the Windows 98 help format, and will replace the old +# Windows help format (.hlp) on all Windows platforms in the future. Compressed +# HTML files also contain an index, a table of contents, and you can search for +# words in the documentation. The HTML workshop also contains a viewer for +# compressed HTML files. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTML is set to YES. + +GENERATE_HTMLHELP = NO + +# The CHM_FILE tag can be used to specify the file name of the resulting .chm +# file. You can add a path in front of the file if the result should not be +# written to the html output directory. +# This tag requires that the tag GENERATE_HTMLHELP is set to YES. + +CHM_FILE = + +# The HHC_LOCATION tag can be used to specify the location (absolute path +# including file name) of the HTML help compiler (hhc.exe). If non-empty, +# doxygen will try to run the HTML help compiler on the generated index.hhp. +# The file has to be specified with full path. +# This tag requires that the tag GENERATE_HTMLHELP is set to YES. + +HHC_LOCATION = + +# The GENERATE_CHI flag controls if a separate .chi index file is generated +# (YES) or that it should be included in the master .chm file (NO). +# The default value is: NO. +# This tag requires that the tag GENERATE_HTMLHELP is set to YES. + +GENERATE_CHI = NO + +# The CHM_INDEX_ENCODING is used to encode HtmlHelp index (hhk), content (hhc) +# and project file content. +# This tag requires that the tag GENERATE_HTMLHELP is set to YES. + +CHM_INDEX_ENCODING = + +# The BINARY_TOC flag controls whether a binary table of contents is generated +# (YES) or a normal table of contents (NO) in the .chm file. Furthermore it +# enables the Previous and Next buttons. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTMLHELP is set to YES. + +BINARY_TOC = NO + +# The TOC_EXPAND flag can be set to YES to add extra items for group members to +# the table of contents of the HTML help documentation and to the tree view. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTMLHELP is set to YES. + +TOC_EXPAND = NO + +# If the GENERATE_QHP tag is set to YES and both QHP_NAMESPACE and +# QHP_VIRTUAL_FOLDER are set, an additional index file will be generated that +# can be used as input for Qt's qhelpgenerator to generate a Qt Compressed Help +# (.qch) of the generated HTML documentation. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTML is set to YES. + +GENERATE_QHP = NO + +# If the QHG_LOCATION tag is specified, the QCH_FILE tag can be used to specify +# the file name of the resulting .qch file. The path specified is relative to +# the HTML output folder. +# This tag requires that the tag GENERATE_QHP is set to YES. + +QCH_FILE = + +# The QHP_NAMESPACE tag specifies the namespace to use when generating Qt Help +# Project output. For more information please see Qt Help Project / Namespace +# (see: http://qt-project.org/doc/qt-4.8/qthelpproject.html#namespace). +# The default value is: org.doxygen.Project. +# This tag requires that the tag GENERATE_QHP is set to YES. + +QHP_NAMESPACE = org.doxygen.Project + +# The QHP_VIRTUAL_FOLDER tag specifies the namespace to use when generating Qt +# Help Project output. For more information please see Qt Help Project / Virtual +# Folders (see: http://qt-project.org/doc/qt-4.8/qthelpproject.html#virtual- +# folders). +# The default value is: doc. +# This tag requires that the tag GENERATE_QHP is set to YES. + +QHP_VIRTUAL_FOLDER = doc + +# If the QHP_CUST_FILTER_NAME tag is set, it specifies the name of a custom +# filter to add. For more information please see Qt Help Project / Custom +# Filters (see: http://qt-project.org/doc/qt-4.8/qthelpproject.html#custom- +# filters). +# This tag requires that the tag GENERATE_QHP is set to YES. + +QHP_CUST_FILTER_NAME = + +# The QHP_CUST_FILTER_ATTRS tag specifies the list of the attributes of the +# custom filter to add. For more information please see Qt Help Project / Custom +# Filters (see: http://qt-project.org/doc/qt-4.8/qthelpproject.html#custom- +# filters). +# This tag requires that the tag GENERATE_QHP is set to YES. + +QHP_CUST_FILTER_ATTRS = + +# The QHP_SECT_FILTER_ATTRS tag specifies the list of the attributes this +# project's filter section matches. Qt Help Project / Filter Attributes (see: +# http://qt-project.org/doc/qt-4.8/qthelpproject.html#filter-attributes). +# This tag requires that the tag GENERATE_QHP is set to YES. + +QHP_SECT_FILTER_ATTRS = + +# The QHG_LOCATION tag can be used to specify the location of Qt's +# qhelpgenerator. If non-empty doxygen will try to run qhelpgenerator on the +# generated .qhp file. +# This tag requires that the tag GENERATE_QHP is set to YES. + +QHG_LOCATION = + +# If the GENERATE_ECLIPSEHELP tag is set to YES, additional index files will be +# generated, together with the HTML files, they form an Eclipse help plugin. To +# install this plugin and make it available under the help contents menu in +# Eclipse, the contents of the directory containing the HTML and XML files needs +# to be copied into the plugins directory of eclipse. The name of the directory +# within the plugins directory should be the same as the ECLIPSE_DOC_ID value. +# After copying Eclipse needs to be restarted before the help appears. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTML is set to YES. + +GENERATE_ECLIPSEHELP = NO + +# A unique identifier for the Eclipse help plugin. When installing the plugin +# the directory name containing the HTML and XML files should also have this +# name. Each documentation set should have its own identifier. +# The default value is: org.doxygen.Project. +# This tag requires that the tag GENERATE_ECLIPSEHELP is set to YES. + +ECLIPSE_DOC_ID = org.doxygen.Project + +# If you want full control over the layout of the generated HTML pages it might +# be necessary to disable the index and replace it with your own. The +# DISABLE_INDEX tag can be used to turn on/off the condensed index (tabs) at top +# of each HTML page. A value of NO enables the index and the value YES disables +# it. Since the tabs in the index contain the same information as the navigation +# tree, you can set this option to YES if you also set GENERATE_TREEVIEW to YES. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTML is set to YES. + +DISABLE_INDEX = YES + +# The GENERATE_TREEVIEW tag is used to specify whether a tree-like index +# structure should be generated to display hierarchical information. If the tag +# value is set to YES, a side panel will be generated containing a tree-like +# index structure (just like the one that is generated for HTML Help). For this +# to work a browser that supports JavaScript, DHTML, CSS and frames is required +# (i.e. any modern browser). Windows users are probably better off using the +# HTML help feature. Via custom style sheets (see HTML_EXTRA_STYLESHEET) one can +# further fine-tune the look of the index. As an example, the default style +# sheet generated by doxygen has an example that shows how to put an image at +# the root of the tree instead of the PROJECT_NAME. Since the tree basically has +# the same information as the tab index, you could consider setting +# DISABLE_INDEX to YES when enabling this option. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTML is set to YES. + +GENERATE_TREEVIEW = NO + +# The ENUM_VALUES_PER_LINE tag can be used to set the number of enum values that +# doxygen will group on one line in the generated HTML documentation. +# +# Note that a value of 0 will completely suppress the enum values from appearing +# in the overview section. +# Minimum value: 0, maximum value: 20, default value: 4. +# This tag requires that the tag GENERATE_HTML is set to YES. + +ENUM_VALUES_PER_LINE = 4 + +# If the treeview is enabled (see GENERATE_TREEVIEW) then this tag can be used +# to set the initial width (in pixels) of the frame in which the tree is shown. +# Minimum value: 0, maximum value: 1500, default value: 250. +# This tag requires that the tag GENERATE_HTML is set to YES. + +TREEVIEW_WIDTH = 250 + +# If the EXT_LINKS_IN_WINDOW option is set to YES, doxygen will open links to +# external symbols imported via tag files in a separate window. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTML is set to YES. + +EXT_LINKS_IN_WINDOW = YES + +# Use this tag to change the font size of LaTeX formulas included as images in +# the HTML documentation. When you change the font size after a successful +# doxygen run you need to manually remove any form_*.png images from the HTML +# output directory to force them to be regenerated. +# Minimum value: 8, maximum value: 50, default value: 10. +# This tag requires that the tag GENERATE_HTML is set to YES. + +FORMULA_FONTSIZE = 10 + +# Use the FORMULA_TRANPARENT tag to determine whether or not the images +# generated for formulas are transparent PNGs. Transparent PNGs are not +# supported properly for IE 6.0, but are supported on all modern browsers. +# +# Note that when changing this option you need to delete any form_*.png files in +# the HTML output directory before the changes have effect. +# The default value is: YES. +# This tag requires that the tag GENERATE_HTML is set to YES. + +FORMULA_TRANSPARENT = YES + +# Enable the USE_MATHJAX option to render LaTeX formulas using MathJax (see +# http://www.mathjax.org) which uses client side Javascript for the rendering +# instead of using pre-rendered bitmaps. Use this if you do not have LaTeX +# installed or if you want to formulas look prettier in the HTML output. When +# enabled you may also need to install MathJax separately and configure the path +# to it using the MATHJAX_RELPATH option. +# The default value is: NO. +# This tag requires that the tag GENERATE_HTML is set to YES. + +USE_MATHJAX = YES + +# When MathJax is enabled you can set the default output format to be used for +# the MathJax output. See the MathJax site (see: +# http://docs.mathjax.org/en/latest/output.html) for more details. +# Possible values are: HTML-CSS (which is slower, but has the best +# compatibility), NativeMML (i.e. MathML) and SVG. +# The default value is: HTML-CSS. +# This tag requires that the tag USE_MATHJAX is set to YES. + +MATHJAX_FORMAT = HTML-CSS + +# When MathJax is enabled you need to specify the location relative to the HTML +# output directory using the MATHJAX_RELPATH option. The destination directory +# should contain the MathJax.js script. For instance, if the mathjax directory +# is located at the same level as the HTML output directory, then +# MATHJAX_RELPATH should be ../mathjax. The default value points to the MathJax +# Content Delivery Network so you can quickly see the result without installing +# MathJax. However, it is strongly recommended to install a local copy of +# MathJax from http://www.mathjax.org before deployment. +# The default value is: http://cdn.mathjax.org/mathjax/latest. +# This tag requires that the tag USE_MATHJAX is set to YES. + +MATHJAX_RELPATH = http://cdn.mathjax.org/mathjax/latest + +# The MATHJAX_EXTENSIONS tag can be used to specify one or more MathJax +# extension names that should be enabled during MathJax rendering. For example +# MATHJAX_EXTENSIONS = TeX/AMSmath TeX/AMSsymbols +# This tag requires that the tag USE_MATHJAX is set to YES. + +MATHJAX_EXTENSIONS = + +# The MATHJAX_CODEFILE tag can be used to specify a file with javascript pieces +# of code that will be used on startup of the MathJax code. See the MathJax site +# (see: http://docs.mathjax.org/en/latest/output.html) for more details. For an +# example see the documentation. +# This tag requires that the tag USE_MATHJAX is set to YES. + +MATHJAX_CODEFILE = + +# When the SEARCHENGINE tag is enabled doxygen will generate a search box for +# the HTML output. The underlying search engine uses javascript and DHTML and +# should work on any modern browser. Note that when using HTML help +# (GENERATE_HTMLHELP), Qt help (GENERATE_QHP), or docsets (GENERATE_DOCSET) +# there is already a search function so this one should typically be disabled. +# For large projects the javascript based search engine can be slow, then +# enabling SERVER_BASED_SEARCH may provide a better solution. It is possible to +# search using the keyboard; to jump to the search box use + S +# (what the is depends on the OS and browser, but it is typically +# , /