Skip to content

Commit

Permalink
Removed US argument to find_depth_of_pressure_in_cell
Browse files Browse the repository at this point in the history
  Eliminated the US argument to find_depth_of_pressure_in_cell, which was no
longer being used.  Also stored EOS_domain values in MOM_state_initialization
for reduced overhead.  All answers are bitwise identical.
  • Loading branch information
Hallberg-NOAA committed Apr 19, 2020
1 parent 9c5239e commit 718c3ab
Show file tree
Hide file tree
Showing 2 changed files with 15 additions and 12 deletions.
5 changes: 2 additions & 3 deletions src/equation_of_state/MOM_EOS.F90
Original file line number Diff line number Diff line change
Expand Up @@ -1159,7 +1159,7 @@ end function query_compressible
subroutine EOS_init(param_file, EOS, US)
type(param_file_type), intent(in) :: param_file !< Parameter file structure
type(EOS_type), pointer :: EOS !< Equation of state structure
type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
optional :: US
! Local variables
#include "version_variable.h"
Expand Down Expand Up @@ -1869,7 +1869,7 @@ end subroutine int_density_dz_generic_plm

!> Find the depth at which the reconstructed pressure matches P_tgt
subroutine find_depth_of_pressure_in_cell(T_t, T_b, S_t, S_b, z_t, z_b, P_t, P_tgt, &
rho_ref, G_e, EOS, US, P_b, z_out, z_tol)
rho_ref, G_e, EOS, P_b, z_out, z_tol)
real, intent(in) :: T_t !< Potential temperatue at the cell top [degC]
real, intent(in) :: T_b !< Potential temperatue at the cell bottom [degC]
real, intent(in) :: S_t !< Salinity at the cell top [ppt]
Expand All @@ -1881,7 +1881,6 @@ subroutine find_depth_of_pressure_in_cell(T_t, T_b, S_t, S_b, z_t, z_b, P_t, P_t
real, intent(in) :: rho_ref !< Reference density with which calculation are anomalous to [R ~> kg m-3]
real, intent(in) :: G_e !< Gravitational acceleration [L2 Z-1 T-2 ~> m s-2]
type(EOS_type), pointer :: EOS !< Equation of state structure
type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
real, intent(out) :: P_b !< Pressure at the bottom of the cell [R L2 T-2 ~> Pa]
real, intent(out) :: z_out !< Absolute depth at which anomalous pressure = p_tgt [Z ~> m]
real, optional, intent(in) :: z_tol !< The tolerance in finding z_out [Z ~> m]
Expand Down
22 changes: 13 additions & 9 deletions src/initialization/MOM_state_initialization.F90
Original file line number Diff line number Diff line change
Expand Up @@ -939,6 +939,7 @@ subroutine convert_thickness(h, G, GV, US, tv)
! [H T2 R-1 L-2 ~> s2 m2 kg-1 or s2 m-1]
real :: HR_to_pres ! A conversion factor from the input geometric thicknesses times the layer
! densities into pressure units [L2 T-2 H-1 ~> m s-2 or m4 kg-1 s-2].
integer, dimension(2) :: EOSdom ! The i-computational domain for the equation of state
integer :: i, j, k, is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz
integer :: itt, max_itt

Expand All @@ -956,11 +957,12 @@ subroutine convert_thickness(h, G, GV, US, tv)
do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1
p_bot(i,j) = 0.0 ; p_top(i,j) = 0.0
enddo ; enddo
EOSdom(:) = EOS_domain(G%HI)
do k=1,nz
do j=js,je
do i=is,ie ; p_top(i,j) = p_bot(i,j) ; enddo
call calculate_density(tv%T(:,j,k), tv%S(:,j,k), p_top(:,j), rho, &
tv%eqn_of_state, dom=EOS_domain(G%HI))
tv%eqn_of_state, EOSdom)
do i=is,ie
p_bot(i,j) = p_top(i,j) + HR_to_pres * (h(i,j,k) * rho(i))
enddo
Expand All @@ -971,7 +973,7 @@ subroutine convert_thickness(h, G, GV, US, tv)
tv%eqn_of_state, dz_geo)
if (itt < max_itt) then ; do j=js,je
call calculate_density(tv%T(:,j,k), tv%S(:,j,k), p_bot(:,j), rho, &
tv%eqn_of_state, dom=EOS_domain(G%HI))
tv%eqn_of_state, EOSdom)
! Use Newton's method to correct the bottom value.
! The hydrostatic equation is sufficiently linear that no bounds-checking is needed.
do i=is,ie
Expand Down Expand Up @@ -1215,7 +1217,7 @@ subroutine cut_off_column_top(nk, tv, GV, US, G_earth, depth, min_thickness, T,
e_top = e(1)
do k=1,nk
call find_depth_of_pressure_in_cell(T_t(k), T_b(k), S_t(k), S_b(k), e(K), e(K+1), &
P_t, p_surf, GV%Rho0, G_earth, tv%eqn_of_state, US, &
P_t, p_surf, GV%Rho0, G_earth, tv%eqn_of_state, &
P_b, z_out, z_tol=z_tol)
if (z_out>=e(K)) then
! Imposed pressure was less that pressure at top of cell
Expand Down Expand Up @@ -1601,7 +1603,7 @@ subroutine initialize_temp_salt_fit(T, S, G, GV, US, param_file, eqn_of_state, P
enddo

call calculate_density(T0(1), S0(1), pres(1), rho_guess(1), eqn_of_state)
call calculate_density_derivs(T0, S0, pres, drho_dT, drho_dS, eqn_of_state, dom=(/1,1/))
call calculate_density_derivs(T0, S0, pres, drho_dT, drho_dS, eqn_of_state, (/1,1/) )

if (fit_salin) then
! A first guess of the layers' temperatures.
Expand Down Expand Up @@ -1735,6 +1737,7 @@ subroutine initialize_sponges_file(G, GV, US, use_temperature, tv, param_file, C
real :: Idamp(SZI_(G),SZJ_(G)) ! The inverse damping rate [T-1 ~> s-1].
real :: pres(SZI_(G)) ! An array of the reference pressure [R L2 T-2 ~> Pa].

integer, dimension(2) :: EOSdom ! The i-computational domain for the equation of state
integer :: i, j, k, is, ie, js, je, nz
integer :: isd, ied, jsd, jed
integer, dimension(4) :: siz
Expand Down Expand Up @@ -1864,13 +1867,13 @@ subroutine initialize_sponges_file(G, GV, US, use_temperature, tv, param_file, C
! mixed layer density, which is used in determining which layers can be
! inflated without causing static instabilities.
do i=is-1,ie ; pres(i) = tv%P_Ref ; enddo
EOSdom(:) = EOS_domain(G%HI)

call MOM_read_data(filename, potemp_var, tmp(:,:,:), G%Domain)
call MOM_read_data(filename, salin_var, tmp2(:,:,:), G%Domain)

do j=js,je
call calculate_density(tmp(:,j,1), tmp2(:,j,1), pres, tmp_2d(:,j), tv%eqn_of_state, &
dom=EOS_domain(G%HI))
call calculate_density(tmp(:,j,1), tmp2(:,j,1), pres, tmp_2d(:,j), tv%eqn_of_state, EOSdom)
enddo

call set_up_sponge_ML_density(tmp_2d, G, CSp)
Expand Down Expand Up @@ -1977,6 +1980,7 @@ subroutine MOM_temp_salt_initialize_from_Z(h, tv, G, GV, US, PF, just_read_param
# include "version_variable.h"
character(len=40) :: mdl = "MOM_initialize_layers_from_Z" ! This module's name.

integer, dimension(2) :: EOSdom ! The i-computational domain for the equation of state
integer :: is, ie, js, je, nz ! compute domain indices
integer :: isc,iec,jsc,jec ! global compute domain indices
integer :: isg, ieg, jsg, jeg ! global extent
Expand Down Expand Up @@ -2188,9 +2192,9 @@ subroutine MOM_temp_salt_initialize_from_Z(h, tv, G, GV, US, PF, just_read_param
call convert_temp_salt_for_TEOS10(temp_z, salt_z, G%HI, kd, mask_z, eos)

press(:) = tv%P_Ref
EOSdom(:) = EOS_domain(G%HI)
do k=1,kd ; do j=js,je
call calculate_density(temp_z(:,j,k), salt_z(:,j,k), press, rho_z(:,j,k), eos, &
dom=EOS_domain(G%HI))
call calculate_density(temp_z(:,j,k), salt_z(:,j,k), press, rho_z(:,j,k), eos, EOSdom)
enddo ; enddo

call pass_var(temp_z,G%Domain)
Expand Down Expand Up @@ -2458,7 +2462,7 @@ subroutine MOM_state_init_tests(G, GV, US, tv)
P_t = 0.
do k = 1, nk
call find_depth_of_pressure_in_cell(T_t(k), T_b(k), S_t(k), S_b(k), e(K), e(K+1), P_t, 0.5*P_tot, &
GV%Rho0, GV%g_Earth, tv%eqn_of_state, US, P_b, z_out)
GV%Rho0, GV%g_Earth, tv%eqn_of_state, P_b, z_out)
write(0,*) k, US%RL2_T2_to_Pa*P_t, US%RL2_T2_to_Pa*P_b, 0.5*US%RL2_T2_to_Pa*P_tot, &
US%Z_to_m*e(K), US%Z_to_m*e(K+1), US%Z_to_m*z_out
P_t = P_b
Expand Down

0 comments on commit 718c3ab

Please sign in to comment.