Skip to content

Remove US from most MOM_EOS function calls #2

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
Hallberg-NOAA merged 5 commits into from
Apr 19, 2020
Merged

Remove US from most MOM_EOS function calls #2

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
415a6bc
Adds unit conversions to EOS type
adcroft Apr 17, 2020
18d520e
Remove optional US from most MOM_EOS functions
adcroft Apr 17, 2020
abecf02
Merge branch 'Hallberg-NOAA-rescale_pressure' into rescale_pressure
adcroft Apr 17, 2020
0caf9cc
Missed a conflict resolution
adcroft Apr 17, 2020
371ed85
Initialize scaling params without EOS_init()
adcroft Apr 17, 2020
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 3 additions & 3 deletions src/ALE/MOM_regridding.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1386,7 +1386,7 @@ subroutine build_rho_grid( G, GV, US, h, tv, dzInterface, remapCS, CS )
! Local depth (G%bathyT is positive)
nominalDepth = G%bathyT(i,j)*GV%Z_to_H

call build_rho_column(CS%rho_CS, US, nz, nominalDepth, h(i, j, :), &
call build_rho_column(CS%rho_CS, nz, nominalDepth, h(i, j, :), &
tv%T(i, j, :), tv%S(i, j, :), tv%eqn_of_state, zNew, &
h_neglect=h_neglect, h_neglect_edge=h_neglect_edge)

Expand Down Expand Up @@ -1501,7 +1501,7 @@ subroutine build_grid_HyCOM1( G, GV, US, h, tv, h_new, dzInterface, CS )
( 0.5 * ( z_col(K) + z_col(K+1) ) * (GV%H_to_RZ*GV%g_Earth) - tv%P_Ref )
enddo

call build_hycom1_column(CS%hycom_CS, US, tv%eqn_of_state, GV%ke, depth, &
call build_hycom1_column(CS%hycom_CS, tv%eqn_of_state, GV%ke, depth, &
h(i,j,:), tv%T(i,j,:), tv%S(i,j,:), p_col, &
z_col, z_col_new, zScale=GV%Z_to_H, &
h_neglect=h_neglect, h_neglect_edge=h_neglect_edge)
Expand Down Expand Up @@ -1635,7 +1635,7 @@ subroutine build_grid_SLight(G, GV, US, h, tv, dzInterface, CS)
( 0.5 * ( z_col(K) + z_col(K+1) ) * (GV%H_to_RZ*GV%g_Earth) - tv%P_Ref )
enddo

call build_slight_column(CS%slight_CS, US, tv%eqn_of_state, GV%H_to_RZ*GV%g_Earth, &
call build_slight_column(CS%slight_CS, tv%eqn_of_state, GV%H_to_RZ*GV%g_Earth, &
GV%H_subroundoff, nz, depth, h(i, j, :), &
tv%T(i, j, :), tv%S(i, j, :), p_col, z_col, z_col_new, &
h_neglect=h_neglect, h_neglect_edge=h_neglect_edge)
Expand Down
6 changes: 2 additions & 4 deletions src/ALE/coord_hycom.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,6 @@ module coord_hycom
! This file is part of MOM6. See LICENSE.md for the license.

use MOM_error_handler, only : MOM_error, FATAL
use MOM_unit_scaling, only : unit_scale_type
use MOM_EOS, only : EOS_type, calculate_density
use regrid_interp, only : interp_CS_type, build_and_interpolate_grid

Expand Down Expand Up @@ -96,10 +95,9 @@ subroutine set_hycom_params(CS, max_interface_depths, max_layer_thickness, inter
end subroutine set_hycom_params

!> Build a HyCOM coordinate column
subroutine build_hycom1_column(CS, US, eqn_of_state, nz, depth, h, T, S, p_col, &
subroutine build_hycom1_column(CS, eqn_of_state, nz, depth, h, T, S, p_col, &
z_col, z_col_new, zScale, h_neglect, h_neglect_edge)
type(hycom_CS), intent(in) :: CS !< Coordinate control structure
type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
type(EOS_type), pointer :: eqn_of_state !< Equation of state structure
integer, intent(in) :: nz !< Number of levels
real, intent(in) :: depth !< Depth of ocean bottom (positive [H ~> m or kg m-2])
Expand Down Expand Up @@ -133,7 +131,7 @@ subroutine build_hycom1_column(CS, US, eqn_of_state, nz, depth, h, T, S, p_col,
z_scale = 1.0 ; if (present(zScale)) z_scale = zScale

! Work bottom recording potential density
call calculate_density(T, S, p_col, rho_col, eqn_of_state, US=US)
call calculate_density(T, S, p_col, rho_col, eqn_of_state)
! This ensures the potential density profile is monotonic
! although not necessarily single valued.
do k = nz-1, 1, -1
Expand Down
11 changes: 4 additions & 7 deletions src/ALE/coord_rho.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,6 @@ module coord_rho

use MOM_error_handler, only : MOM_error, FATAL
use MOM_remapping, only : remapping_CS, remapping_core_h
use MOM_unit_scaling, only : unit_scale_type
use MOM_EOS, only : EOS_type, calculate_density
use regrid_interp, only : interp_CS_type, build_and_interpolate_grid, DEGREE_MAX

Expand Down Expand Up @@ -88,10 +87,9 @@ end subroutine set_rho_params
!!
!! 1. Density profiles are calculated on the source grid.
!! 2. Positions of target densities (for interfaces) are found by interpolation.
subroutine build_rho_column(CS, US, nz, depth, h, T, S, eqn_of_state, z_interface, &
subroutine build_rho_column(CS, nz, depth, h, T, S, eqn_of_state, z_interface, &
h_neglect, h_neglect_edge)
type(rho_CS), intent(in) :: CS !< coord_rho control structure
type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
integer, intent(in) :: nz !< Number of levels on source grid (i.e. length of h, T, S)
real, intent(in) :: depth !< Depth of ocean bottom (positive in m)
real, dimension(nz), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
Expand Down Expand Up @@ -126,7 +124,7 @@ subroutine build_rho_column(CS, US, nz, depth, h, T, S, eqn_of_state, z_interfac

! Compute densities on source column
pres(:) = CS%ref_pressure
call calculate_density(T, S, pres, densities, eqn_of_state, US=US)
call calculate_density(T, S, pres, densities, eqn_of_state)
do k = 1,count_nonzero_layers
densities(k) = densities(mapping(k))
enddo
Expand Down Expand Up @@ -185,10 +183,9 @@ end subroutine build_rho_column
!! 4. T & S are remapped onto the new grid.
!! 5. Return to step 1 until convergence or until the maximum number of
!! iterations is reached, whichever comes first.
subroutine build_rho_column_iteratively(CS, US, remapCS, nz, depth, h, T, S, eqn_of_state, &
subroutine build_rho_column_iteratively(CS, remapCS, nz, depth, h, T, S, eqn_of_state, &
zInterface, h_neglect, h_neglect_edge, dev_tol)
type(rho_CS), intent(in) :: CS !< Regridding control structure
type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
type(remapping_CS), intent(in) :: remapCS !< Remapping parameters and options
integer, intent(in) :: nz !< Number of levels
real, intent(in) :: depth !< Depth of ocean bottom [Z ~> m]
Expand Down Expand Up @@ -250,7 +247,7 @@ subroutine build_rho_column_iteratively(CS, US, remapCS, nz, depth, h, T, S, eqn
enddo

! Compute densities within current water column
call calculate_density( T_tmp, S_tmp, pres, densities, 1, nz, eqn_of_state, US=US)
call calculate_density( T_tmp, S_tmp, pres, densities, 1, nz, eqn_of_state)

do k = 1,count_nonzero_layers
densities(k) = densities(mapping(k))
Expand Down
12 changes: 5 additions & 7 deletions src/ALE/coord_slight.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,6 @@ module coord_slight
! This file is part of MOM6. See LICENSE.md for the license.

use MOM_error_handler, only : MOM_error, FATAL
use MOM_unit_scaling, only : unit_scale_type
use MOM_EOS, only : EOS_type, calculate_compress
use MOM_EOS, only : calculate_density, calculate_density_derivs
use regrid_interp, only : interp_CS_type, regridding_set_ppolys
Expand Down Expand Up @@ -178,11 +177,10 @@ subroutine set_slight_params(CS, max_interface_depths, max_layer_thickness, &
end subroutine set_slight_params

!> Build a SLight coordinate column
subroutine build_slight_column(CS, US, eqn_of_state, H_to_pres, H_subroundoff, &
subroutine build_slight_column(CS, eqn_of_state, H_to_pres, H_subroundoff, &
nz, depth, h_col, T_col, S_col, p_col, z_col, z_col_new, &
h_neglect, h_neglect_edge)
type(slight_CS), intent(in) :: CS !< Coordinate control structure
type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
type(EOS_type), pointer :: eqn_of_state !< Equation of state structure
real, intent(in) :: H_to_pres !< A conversion factor from thicknesses to
!! scaled pressure [R L2 T-2 H-1 ~> Pa m-1 or Pa m2 kg-1]
Expand Down Expand Up @@ -253,7 +251,7 @@ subroutine build_slight_column(CS, US, eqn_of_state, H_to_pres, H_subroundoff, &
dz = (z_col(nz+1) - z_col(1)) / real(nz)
do K=2,nz ; z_col_new(K) = z_col(1) + dz*real(K-1) ; enddo
else
call calculate_density(T_col, S_col, p_col, rho_col, eqn_of_state, US=US)
call calculate_density(T_col, S_col, p_col, rho_col, eqn_of_state)

! Find the locations of the target potential densities, flagging
! locations in apparently unstable regions as not reliable.
Expand Down Expand Up @@ -375,11 +373,11 @@ subroutine build_slight_column(CS, US, eqn_of_state, H_to_pres, H_subroundoff, &
T_int(nz+1) = T_f(nz) ; S_int(nz+1) = S_f(nz)
p_IS(nz+1) = z_col(nz+1) * H_to_pres
call calculate_density_derivs(T_int, S_int, p_IS, drhoIS_dT, drhoIS_dS, &
eqn_of_state, US, dom=(/2,nz/))
eqn_of_state, dom=(/2,nz/))
call calculate_density_derivs(T_int, S_int, p_R, drhoR_dT, drhoR_dS, &
eqn_of_state, US, dom=(/2,nz/))
eqn_of_state, dom=(/2,nz/))
if (CS%compressibility_fraction > 0.0) then
call calculate_compress(T_int, S_int, p_R(:), rho_tmp, drho_dp, 2, nz-1, eqn_of_state, US)
call calculate_compress(T_int, S_int, p_R(:), rho_tmp, drho_dp, 2, nz-1, eqn_of_state)
else
do K=2,nz ; drho_dp(K) = 0.0 ; enddo
endif
Expand Down
4 changes: 2 additions & 2 deletions src/core/MOM.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2242,7 +2242,7 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
! Use the Wright equation of state by default, unless otherwise specified
! Note: this line and the following block ought to be in a separate
! initialization routine for tv.
if (use_EOS) call EOS_init(param_file, CS%tv%eqn_of_state)
if (use_EOS) call EOS_init(param_file, CS%tv%eqn_of_state, US)
if (use_temperature) then
allocate(CS%tv%TempxPmE(isd:ied,jsd:jed)) ; CS%tv%TempxPmE(:,:) = 0.0
if (use_geothermal) then
Expand Down Expand Up @@ -2905,7 +2905,7 @@ subroutine adjust_ssh_for_p_atm(tv, G, GV, US, ssh, p_atm, use_EOS)
do j=js,je
if (calc_rho) then
call calculate_density(tv%T(:,j,1), tv%S(:,j,1), 0.5*p_atm(:,j), Rho_conv, &
tv%eqn_of_state, US, dom=EOS_domain(G%HI))
tv%eqn_of_state, dom=EOS_domain(G%HI))
do i=is,ie
IgR0 = US%Z_to_m / (Rho_conv(i) * GV%g_Earth)
ssh(i,j) = ssh(i,j) + p_atm(i,j) * IgR0
Expand Down
20 changes: 10 additions & 10 deletions src/core/MOM_PressureForce_Montgomery.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -188,7 +188,7 @@ subroutine PressureForce_Mont_nonBouss(h, tv, PFu, PFv, G, GV, US, CS, p_atm, pb
!$OMP parallel do default(shared)
do k=1,nz
call int_specific_vol_dp(tv%T(:,:,k), tv%S(:,:,k), p(:,:,k), p(:,:,k+1), &
0.0, G%HI, tv%eqn_of_state, dz_geo(:,:,k), halo_size=1, US=US)
0.0, G%HI, tv%eqn_of_state, dz_geo(:,:,k), halo_size=1)
enddo
!$OMP parallel do default(shared)
do j=Jsq,Jeq+1 ; do k=1,nz ; do i=Isq,Ieq+1
Expand Down Expand Up @@ -230,7 +230,7 @@ subroutine PressureForce_Mont_nonBouss(h, tv, PFu, PFv, G, GV, US, CS, p_atm, pb
tv_tmp%T(i,j,k) = tv%T(i,j,k) ; tv_tmp%S(i,j,k) = tv%S(i,j,k)
enddo ; enddo
call calculate_density(tv%T(:,j,nkmb), tv%S(:,j,nkmb), p_ref, Rho_cv_BL(:), &
tv%eqn_of_state, US=US, dom=EOSdom)
tv%eqn_of_state, dom=EOSdom)
do k=nkmb+1,nz ; do i=Isq,Ieq+1
if (GV%Rlay(k) < Rho_cv_BL(i)) then
tv_tmp%T(i,j,k) = tv%T(i,j,nkmb) ; tv_tmp%S(i,j,k) = tv%S(i,j,nkmb)
Expand All @@ -247,7 +247,7 @@ subroutine PressureForce_Mont_nonBouss(h, tv, PFu, PFv, G, GV, US, CS, p_atm, pb
!$OMP parallel do default(shared) private(rho_in_situ)
do k=1,nz ; do j=Jsq,Jeq+1
call calculate_density(tv_tmp%T(:,j,k), tv_tmp%S(:,j,k), p_ref, rho_in_situ, &
tv%eqn_of_state, US=US, dom=EOSdom)
tv%eqn_of_state, dom=EOSdom)
do i=Isq,Ieq+1 ; alpha_star(i,j,k) = 1.0 / rho_in_situ(i) ; enddo
enddo ; enddo
endif ! use_EOS
Expand Down Expand Up @@ -487,7 +487,7 @@ subroutine PressureForce_Mont_Bouss(h, tv, PFu, PFv, G, GV, US, CS, p_atm, pbce,
tv_tmp%T(i,j,k) = tv%T(i,j,k) ; tv_tmp%S(i,j,k) = tv%S(i,j,k)
enddo ; enddo
call calculate_density(tv%T(:,j,nkmb), tv%S(:,j,nkmb), p_ref, Rho_cv_BL(:), &
tv%eqn_of_state, US=US, dom=EOSdom)
tv%eqn_of_state, dom=EOSdom)

do k=nkmb+1,nz ; do i=Isq,Ieq+1
if (GV%Rlay(k) < Rho_cv_BL(i)) then
Expand All @@ -508,7 +508,7 @@ subroutine PressureForce_Mont_Bouss(h, tv, PFu, PFv, G, GV, US, CS, p_atm, pbce,
!$OMP parallel do default(shared)
do k=1,nz+1 ; do j=Jsq,Jeq+1
call calculate_density(tv_tmp%T(:,j,k), tv_tmp%S(:,j,k), p_ref, rho_star(:,j,k), &
tv%eqn_of_state, US=US, dom=EOSdom)
tv%eqn_of_state, dom=EOSdom)
do i=Isq,Ieq+1 ; rho_star(i,j,k) = G_Rho0*rho_star(i,j,k) ; enddo
enddo ; enddo
endif ! use_EOS
Expand Down Expand Up @@ -666,7 +666,7 @@ subroutine Set_pbce_Bouss(e, tv, G, GV, US, Rho0, GFS_scale, pbce, rho_star)
press(i) = -Rho0xG*e(i,j,1)
enddo
call calculate_density(tv%T(:,j,1), tv%S(:,j,1), press, rho_in_situ, &
tv%eqn_of_state, US=US, dom=EOSdom)
tv%eqn_of_state, dom=EOSdom)
do i=Isq,Ieq+1
pbce(i,j,1) = G_Rho0*(GFS_scale * rho_in_situ(i)) * GV%H_to_Z
enddo
Expand All @@ -677,7 +677,7 @@ subroutine Set_pbce_Bouss(e, tv, G, GV, US, Rho0, GFS_scale, pbce, rho_star)
S_int(i) = 0.5*(tv%S(i,j,k-1)+tv%S(i,j,k))
enddo
call calculate_density_derivs(T_int, S_int, press, dR_dT, dR_dS, &
tv%eqn_of_state, US=US, dom=EOSdom)
tv%eqn_of_state, dom=EOSdom)
do i=Isq,Ieq+1
pbce(i,j,k) = pbce(i,j,k-1) + G_Rho0 * &
((e(i,j,K) - e(i,j,nz+1)) * Ihtot(i)) * &
Expand Down Expand Up @@ -764,7 +764,7 @@ subroutine Set_pbce_nonBouss(p, tv, G, GV, US, GFS_scale, pbce, alpha_star)
!$OMP parallel do default(shared) private(T_int,S_int,dR_dT,dR_dS,rho_in_situ)
do j=Jsq,Jeq+1
call calculate_density(tv%T(:,j,nz), tv%S(:,j,nz), p(:,j,nz+1), rho_in_situ, &
tv%eqn_of_state, US=US, dom=EOSdom)
tv%eqn_of_state, dom=EOSdom)
do i=Isq,Ieq+1
C_htot(i,j) = dP_dH / ((p(i,j,nz+1)-p(i,j,1)) + dp_neglect)
pbce(i,j,nz) = dP_dH / (rho_in_situ(i))
Expand All @@ -774,9 +774,9 @@ subroutine Set_pbce_nonBouss(p, tv, G, GV, US, GFS_scale, pbce, alpha_star)
T_int(i) = 0.5*(tv%T(i,j,k)+tv%T(i,j,k+1))
S_int(i) = 0.5*(tv%S(i,j,k)+tv%S(i,j,k+1))
enddo
call calculate_density(T_int, S_int, p(:,j,k+1), rho_in_situ, tv%eqn_of_state, US=US, dom=EOSdom)
call calculate_density(T_int, S_int, p(:,j,k+1), rho_in_situ, tv%eqn_of_state, dom=EOSdom)
call calculate_density_derivs(T_int, S_int, p(:,j,k+1), dR_dT, dR_dS, &
tv%eqn_of_state, US=US, dom=EOSdom)
tv%eqn_of_state, dom=EOSdom)
do i=Isq,Ieq+1
pbce(i,j,k) = pbce(i,j,k+1) + ((p(i,j,K+1)-p(i,j,1))*C_htot(i,j)) * &
((dR_dT(i)*(tv%T(i,j,k+1)-tv%T(i,j,k)) + &
Expand Down
Loading