Skip to content

Commit

Permalink
Speed up redsea gulf bay fix (mom-ocean#247)
Browse files Browse the repository at this point in the history 
* Clean up red_sea_fix routine

* remove global fields from red sea fix

* Pass Salt directly into redsea_gulfbay_hmix_s

* Fixed typos in interface to gulfbay fix. Also fixed references to Salt%field with incorrect rank, and
call to mpp_sum.
  • Loading branch information
@aidanheerdegen
aidanheerdegen authored Sep 19, 2018
1 parent afe80bf commit 4af4189
Show file tree
Hide file tree
Showing 2 changed files with 69 additions and 94 deletions.
2 changes: 1 addition & 1 deletion bin/mkmf.template.travis
View file
Original file line number Diff line number Diff line change
Expand Up @@ -24,7 +24,7 @@ INCLUDE := -I$(shell nc-config --includedir) \

FPPFLAGS := $(INCLUDE)

FFLAGS := -O2 -fcray-pointer -fdefault-real-8 -ffree-line-length-none -fno-range-check -W
FFLAGS := -O2 -fcray-pointer -fdefault-real-8 -ffree-line-length-none -fno-range-check -w

CFLAGS := -D__IFC -O2 $(shell nc-config --cflags) -I/usr/include/openmpi $(INCLUDE)

Expand Down
161 changes: 68 additions & 93 deletions src/mom5/ocean_core/ocean_model.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2052,8 +2052,8 @@ subroutine update_ocean_model(Ice_ocean_boundary, Ocean_state, Ocean_sfc, &
if (mpp_pe() == mpp_root_pe()) then
call print_time(time_start_update, 'Calling redsea_gulfbay_hmix_s at runtime = ')
endif
call redsea_gulfbay_hmix_s(Time, Grid, Thickness, &
T_prog(1:num_prog_tracers), Ocean_sfc)
call redsea_gulfbay_hmix_s(Time, Grid, Domain, Thickness, &
T_prog(index_salt))
call mpp_clock_end(id_sfix)
endif
#endif
Expand Down Expand Up @@ -2083,140 +2083,115 @@ end subroutine update_ocean_model
! </SUBROUTINE> NAME="update_ocean_model"

#if defined(ACCESS)
subroutine redsea_gulfbay_hmix_s(Time, Grid, Thickness, T_prog, Ocean_sfc)
subroutine redsea_gulfbay_hmix_s(Time, Grid, Domain, Thickness, Salt)

use mpp_domains_mod, only : mpp_global_field, mpp_get_data_domain
use mpp_mod, only : mpp_broadcast
use mpp_mod, only : mpp_sum

use auscom_ice_parameters_mod, only : irs1, ire1, jrs1, jre1, irs2, ire2,jrs2, jre2, &
igs, ige, jgs, jge, ksmax

implicit none

type(ocean_time_type), intent(in) :: Time
type(ocean_grid_type), target :: Grid ! domain and grid information for ocean model
type(ocean_grid_type), intent(in) :: Grid ! grid information for ocean model
type(ocean_domain_type), intent(in) :: Domain ! domain grid information for ocean model
type(ocean_thickness_type), intent(in) :: Thickness
type(ocean_prog_tracer_type), intent(inout) :: T_prog(:)
type(ocean_public_type), intent(in) :: Ocean_sfc
type(ocean_prog_tracer_type), intent(inout) :: Salt

real, dimension(:,:,:), allocatable :: global_tmask ! for global mask
real, dimension(:,:,:), allocatable :: global_dzt ! for global dzt
real, dimension(:,:,:), allocatable :: global_sp ! for global salinity
real, dimension(:,:) , allocatable :: global_dat ! for global area

real :: volume = 0.0
real :: wetvolume = 0.0
real :: tot_sp = 0.0
real :: ave_sp = 0.0

real, dimension(:,:,:), allocatable :: salt_vol_sums ! total salt and volume sums. Red sea salt(:,1,1) vol(:,2,1)
! Gulf salt(:,1,2) vol(:,2,2)

integer :: tau, taup1
integer :: i, j, k

integer :: nx, ny, nz
integer :: iisd, iied, jjsd, jjed

nx = Grid%ni
ny = Grid%nj
nz = Grid%nk
integer :: isc, iec, jsc, jec
integer :: isd, ied, jsd, jed

tau = Time%tau
taup1 = Time%taup1

allocate (global_tmask(nx,ny,nz)) ; global_tmask=0.0
call mpp_global_field(Domain%domain2d, Grid%tmask, global_tmask)
allocate (global_dat(nx,ny)) ; global_dat=0.0
call mpp_global_field(Domain%domain2d, Grid%dat, global_dat)

allocate (global_dzt(nx,ny,nz)) ; global_dzt=0.0
call mpp_global_field(Domain%domain2d, Thickness%dzt(:,:,:), global_dzt)
allocate (global_sp(nx,ny,nz)) ; global_sp=0.0
call mpp_global_field(Domain%domain2d, T_prog(index_salt)%field(:,:,:,taup1),global_sp)

call mpp_get_data_domain(Ocean_sfc%domain, iisd, iied, jjsd, jjed)
call get_local_indices(Domain, isd, ied, jsd, jed, isc, iec, jsc, jec)

allocate(salt_vol_sums(ksmax,2,2))
salt_vol_sums = 0.0
do k = 1, ksmax
!
!for Red Sea
!
wetvolume = 0.0
tot_sp = 0.0
do j=jrs1,jre1
do i=irs1,ire1
if(global_tmask(i,j,k) == 1.0) then
volume = global_dat(i,j) * global_dzt(i,j,k)
wetvolume = wetvolume + volume
tot_sp = tot_sp + global_sp(i,j,k) * volume
do j=max(jrs1,jsc),min(jre1,jec)
do i=max(irs1,isc),min(ire1,iec)
if(Grid%tmask(i,j,k) == 1.0) then
volume = Grid%dat(i,j) * Thickness%dzt(i,j,k)
salt_vol_sums(k,1,1) = salt_vol_sums(k,1,1) + Salt%field(i,j,k,taup1) * volume
salt_vol_sums(k,2,1) = salt_vol_sums(k,2,1) + volume
endif
enddo
enddo
do j=jrs2,jre2
do i=irs2,ire2
if(global_tmask(i,j,k) == 1.0) then
volume = global_dat(i,j) * global_dzt(i,j,k)
wetvolume = wetvolume + volume
tot_sp = tot_sp + global_sp(i,j,k) * volume
do j=max(jrs2,jsc),min(jre2,jec)
do i=max(irs2,isc),min(ire2,iec)
if(Grid%tmask(i,j,k) == 1.0) then
volume = Grid%dat(i,j) * Thickness%dzt(i,j,k)
salt_vol_sums(k,1,1) = salt_vol_sums(k,1,1) + Salt%field(i,j,k,taup1) * volume
salt_vol_sums(k,2,1) = salt_vol_sums(k,2,1) + volume
endif
enddo
enddo
if (wetvolume /= 0.0) then
ave_sp = tot_sp/wetvolume
do j=jrs1,jre1
do i=irs1,ire1
if(global_tmask(i,j,k) == 1.0) then
global_sp(i,j,k) = ave_sp
endif
enddo
enddo
do j=jrs2,jre2
do i=irs2,ire2
if(global_tmask(i,j,k) == 1.0) then
global_sp(i,j,k) = ave_sp
endif
enddo
enddo
endif
!
!for Gulf Bay
!
wetvolume = 0.0
tot_sp = 0.0
do j=jgs,jge
do i=igs,ige
if(global_tmask(i,j,k) == 1.0) then
volume = global_dat(i,j) * global_dzt(i,j,k)
wetvolume = wetvolume + volume
tot_sp = tot_sp + global_sp(i,j,k) * volume
do j=max(jgs,jsc),min(jge,jec)
do i=max(igs,isc),min(ige,iec)
if(Grid%tmask(i,j,k) == 1.0) then
volume = Grid%dat(i,j) * Thickness%dzt(i,j,k)
salt_vol_sums(k,1,2) = salt_vol_sums(k,1,2) + Salt%field(i,j,k,taup1) * volume
salt_vol_sums(k,2,2) = salt_vol_sums(k,2,2) + volume
endif
enddo
enddo
if (wetvolume /= 0.0) then
ave_sp = tot_sp/wetvolume
do j=jgs,jge
do i=igs,ige
if(global_tmask(i,j,k) == 1.0) then
global_sp(i,j,k) = ave_sp
endif
enddo
enddo
endif

call mpp_broadcast(global_sp(:,:,k),nx*ny,mpp_root_pe())
T_prog(index_salt)%field(iisd:iied,jjsd:jjed,k,taup1) = global_sp(iisd:iied,jjsd:jjed,k)
enddo !k=1,ksmax

call mpp_sum(salt_vol_sums, size(salt_vol_sums))

! Replace by sums
! Note that we fill to domain edge. No need to update halo

enddo !k=1,kdmax
! Red Sea
do k=1,ksmax
if ( salt_vol_sums(k,2,1) /= 0.0 ) then
ave_sp = salt_vol_sums(k,1,1) / salt_vol_sums(k,2,1)
do j=max(jrs1,jsd),min(jre1,jed)
do i=max(irs1,isd),min(ire1,ied)
Salt%field(i,j,k,taup1) = ave_sp * Grid%tmask(i,j,k)
enddo
enddo
do j=max(jrs2,jsd),min(jre2,jed)
do i=max(irs2,isd),min(ire2,ied)
Salt%field(i,j,k,taup1) = ave_sp * Grid%tmask(i,j,k)
enddo
enddo
endif

! if(mpp_pe() == mpp_root_pe()) then
! write(115,'(10e12.5)') global_sp
! endif
!Persian gulf

! global_sp = 0.0
! call mpp_global_field(Domain%domain2d,
! T_prog(index_salt)%field(:,:,:,taup1),global_sp)
! if(mpp_pe() == mpp_root_pe()) then
! write(116,'(10e12.5)') global_sp
! endif
if ( salt_vol_sums(k,2,2) /= 0.0 ) then
ave_sp = salt_vol_sums(k,1,2) / salt_vol_sums(k,2,2)
do j=max(jgs,jsd),min(jge,jed)
do i=max(igs,isd),min(ige,ied)
Salt%field(i,j,k,taup1) = ave_sp * Grid%tmask(i,j,k)
enddo
enddo
endif
enddo




deallocate (global_tmask, global_dzt, global_sp, global_dat)
deallocate (salt_vol_sums)

end subroutine redsea_gulfbay_hmix_s
#endif
Expand Down

0 comments on commit 4af4189

Please sign in to comment.