Skip to content

Commit

Permalink
more updates to get nuopc and mct closer
Browse files Browse the repository at this point in the history
  • Loading branch information
Mariana Vertenstein committed May 28, 2018
1 parent f482863 commit 9722934
Show file tree
Hide file tree
Showing 4 changed files with 142 additions and 159 deletions.
121 changes: 4 additions & 117 deletions config_src/mct_driver/MOM_ocean_model.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -94,11 +94,9 @@ module MOM_ocean_model
public ocean_model_save_restart, Ocean_stock_pe
public ocean_model_init_sfc, ocean_model_flux_init
public ocean_model_restart
public ice_ocean_boundary_type
public ice_ocn_bnd_type_chksum
public ocean_public_type_chksum
public ocean_model_data_get
public ocn_export
public ice_ocn_bnd_type_chksum

interface ocean_model_data_get
module procedure ocean_model_data1D_get
Expand Down Expand Up @@ -394,8 +392,7 @@ end subroutine ocean_model_init
!! of Ocean_coupling_time_step, returning the publicly visible ocean surface properties in
!! Ocean_sfc and storing the new ocean properties in Ocean_state.
subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, &
time_start_update, Ocean_coupling_time_step, &
x2o_o, ind)
time_start_update, Ocean_coupling_time_step)

type(ice_ocean_boundary_type), &
intent(in) :: Ice_ocean_boundary !< A structure containing the
Expand All @@ -415,9 +412,6 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, &
type(time_type), intent(in) :: Ocean_coupling_time_step !< The amount of time over
!! which to advance the ocean.

real(kind=8), intent(in) :: x2o_o(:,:) !< Fluxes from coupler to ocean, computed by ocean
type(cpl_indices_type), intent(inout) :: ind !< Structure with MCT attribute vectors and indices

! local variables
type(time_type) :: Master_time !< This allows step_MOM to temporarily change
!! the time that is seen by internal modules.
Expand Down Expand Up @@ -459,8 +453,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, &

! Import fluxes from coupler to ocean. Also, perform do SST and SSS restoring, if needed.
call convert_IOB_to_fluxes_and_forces(Ice_ocean_boundary, OS%fluxes, OS%Time, OS%grid, OS%forcing_CSp, &
OS%sfc_state, OS%restore_salinity, OS%restore_temp, &
forces=OS%forces, x2o=x2o_o, ind=ind)
OS%sfc_state, OS%restore_salinity, OS%restore_temp, forces=OS%forces)

! Fields that exist in both the forcing and mech_forcing types must be copied.
call copy_common_forcing_fields(OS%forces, OS%fluxes, OS%grid)
Expand Down Expand Up @@ -491,8 +484,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, &

! Import fluxes from coupler to ocean. Also, perform do SST and SSS restoring, if needed.
call convert_IOB_to_fluxes_and_forces(Ice_ocean_boundary, OS%fluxes, OS%Time, OS%grid, OS%forcing_CSp, &
OS%sfc_state, OS%restore_salinity, OS%restore_temp, &
forces=OS%forces, x2o=x2o_o, ind=ind)
OS%sfc_state, OS%restore_salinity, OS%restore_temp, forces=OS%forces)

if (OS%use_ice_shelf) then
call shelf_calc_flux(OS%sfc_state, OS%forces, OS%flux_tmp, OS%Time, time_step, OS%Ice_shelf_CSp)
Expand Down Expand Up @@ -1037,109 +1029,4 @@ subroutine get_ocean_grid(OS, Gridp)
end subroutine get_ocean_grid
! </SUBROUTINE> NAME="get_ocean_grid"

!=======================================================================
! Routines that are specific to MCT driver
!=======================================================================


!=======================================================================

!> Maps outgoing ocean data to MCT buffer.
!! See \ref section_ocn_export for a summary of the data
!! that is transferred from MOM6 to MCT.
subroutine ocn_export(ind, ocn_public, grid, o2x)
type(cpl_indices_type), intent(inout) :: ind !< Structure with coupler indices and vectors
type(ocean_public_type), intent(in) :: ocn_public !< Ocean surface state
type(ocean_grid_type), intent(in) :: grid !< Ocean model grid
real(kind=8), intent(inout) :: o2x(:,:) !< MCT outgoing bugger

! Local variables
real, dimension(grid%isd:grid%ied,grid%jsd:grid%jed) :: ssh !< Local copy of sea_lev with updated halo
integer :: i, j, n, ig, jg !< Grid indices
real :: slp_L, slp_R, slp_C, slope, u_min, u_max

! Copy from ocn_public to o2x. ocn_public uses global indexing with no halos.
! The mask comes from "grid" that uses the usual MOM domain that has halos
! and does not use global indexing.
n = 0
do j=grid%jsc, grid%jec
jg = j + grid%jdg_offset
do i=grid%isc,grid%iec
n = n+1
ig = i + grid%idg_offset
! surface temperature in Kelvin
o2x(ind%o2x_So_t, n) = ocn_public%t_surf(ig,jg) * grid%mask2dT(i,j)
o2x(ind%o2x_So_s, n) = ocn_public%s_surf(ig,jg) * grid%mask2dT(i,j)
o2x(ind%o2x_So_u, n) = ocn_public%u_surf(ig,jg) * grid%mask2dT(i,j)
o2x(ind%o2x_So_v, n) = ocn_public%v_surf(ig,jg) * grid%mask2dT(i,j)
! Make a copy of ssh in order to do a halo update. We use the usual MOM domain
! in order to update halos. i.e. does not use global indexing.
ssh(i,j) = ocn_public%sea_lev(ig,jg)
end do
end do

! Update halo of ssh so we can calculate gradients
call pass_var(ssh, grid%domain)

! d/dx ssh
n = 0
do j=grid%jsc, grid%jec ; do i=grid%isc,grid%iec
n = n+1
! This is a simple second-order difference
!o2x(ind%o2x_So_dhdx, n) = 0.5 * (ssh(i+1,j) - ssh(i-1,j)) * grid%IdxT(i,j) * grid%mask2dT(i,j)
! This is a PLM slope which might be less prone to the A-grid null mode
slp_L = (ssh(I,j) - ssh(I-1,j)) * grid%mask2dCu(I-1,j)
if (grid%mask2dCu(I-1,j)==0.) slp_L = 0.
slp_R = (ssh(I+1,j) - ssh(I,j)) * grid%mask2dCu(I,j)
if (grid%mask2dCu(I+1,j)==0.) slp_R = 0.
slp_C = 0.5 * (slp_L + slp_R)
if ( (slp_L * slp_R) > 0.0 ) then
! This limits the slope so that the edge values are bounded by the
! two cell averages spanning the edge.
u_min = min( ssh(i-1,j), ssh(i,j), ssh(i+1,j) )
u_max = max( ssh(i-1,j), ssh(i,j), ssh(i+1,j) )
slope = sign( min( abs(slp_C), 2.*min( ssh(i,j) - u_min, u_max - ssh(i,j) ) ), slp_C )
else
! Extrema in the mean values require a PCM reconstruction avoid generating
! larger extreme values.
slope = 0.0
end if
o2x(ind%o2x_So_dhdx, n) = slope * grid%IdxT(i,j) * grid%mask2dT(i,j)
if (grid%mask2dT(i,j)==0.) o2x(ind%o2x_So_dhdx, n) = 0.0
end do; end do

! d/dy ssh
n = 0
do j=grid%jsc, grid%jec ; do i=grid%isc,grid%iec
n = n+1
! This is a simple second-order difference
! o2x(ind%o2x_So_dhdy, n) = 0.5 * (ssh(i,j+1) - ssh(i,j-1)) * grid%IdyT(i,j) * grid%mask2dT(i,j)
! This is a PLM slope which might be less prone to the A-grid null mode
slp_L = ssh(i,J) - ssh(i,J-1) * grid%mask2dCv(i,J-1)
if (grid%mask2dCv(i,J-1)==0.) slp_L = 0.

slp_R = ssh(i,J+1) - ssh(i,J) * grid%mask2dCv(i,J)
if (grid%mask2dCv(i,J+1)==0.) slp_R = 0.

slp_C = 0.5 * (slp_L + slp_R)
!write(6,*)'slp_L, slp_R,i,j,slp_L*slp_R', slp_L, slp_R,i,j,slp_L*slp_R
if ((slp_L * slp_R) > 0.0) then
! This limits the slope so that the edge values are bounded by the
! two cell averages spanning the edge.
u_min = min( ssh(i,j-1), ssh(i,j), ssh(i,j+1) )
u_max = max( ssh(i,j-1), ssh(i,j), ssh(i,j+1) )
slope = sign( min( abs(slp_C), 2.*min( ssh(i,j) - u_min, u_max - ssh(i,j) ) ), slp_C )
else
! Extrema in the mean values require a PCM reconstruction avoid generating
! larger extreme values.
slope = 0.0
end if
o2x(ind%o2x_So_dhdy, n) = slope * grid%IdyT(i,j) * grid%mask2dT(i,j)
if (grid%mask2dT(i,j)==0.) o2x(ind%o2x_So_dhdy, n) = 0.0
end do; end do

end subroutine ocn_export

!=======================================================================

end module MOM_ocean_model
51 changes: 24 additions & 27 deletions config_src/mct_driver/MOM_surface_forcing.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -150,10 +150,12 @@ module MOM_surface_forcing
! the elements, units, and conventions that exactly conform to the use for
! MOM-based coupled models.
type, public :: ice_ocean_boundary_type
real, pointer, dimension(:,:) :: latent_flux =>NULL() !< latent flux (W/m2)
real, pointer, dimension(:,:) :: rofl_flux =>NULL() !< liquid runoff (W/m2)
real, pointer, dimension(:,:) :: rofi_flux =>NULL() !< ice runoff (W/m2)
real, pointer, dimension(:,:) :: u_flux =>NULL() !< i-direction wind stress (Pa)
real, pointer, dimension(:,:) :: v_flux =>NULL() !< j-direction wind stress (Pa)
real, pointer, dimension(:,:) :: t_flux =>NULL() !< sensible heat flux (W/m2)
real, pointer, dimension(:,:) :: latent_flux =>NULL() !< latent flux (W/m2)
real, pointer, dimension(:,:) :: q_flux =>NULL() !< specific humidity flux (kg/m2/s)
real, pointer, dimension(:,:) :: salt_flux =>NULL() !< salt flux (kg/m2/s)
real, pointer, dimension(:,:) :: lw_flux =>NULL() !< long wave radiation (W/m2)
Expand Down Expand Up @@ -200,8 +202,7 @@ module MOM_surface_forcing
!! passed from MCT to MOM6, including fluxes that need to be included in
!! the future.
subroutine convert_IOB_to_fluxes_and_forces(IOB, fluxes, Time, G, CS, &
sfc_state, restore_salt, restore_temp, &
forces, x2o, ind)
sfc_state, restore_salt, restore_temp, forces)

type(ice_ocean_boundary_type), &
target, intent(in) :: IOB !< An ice-ocean boundary type with fluxes to drive
Expand All @@ -222,8 +223,6 @@ subroutine convert_IOB_to_fluxes_and_forces(IOB, fluxes, Time, G, CS, &
logical, optional, intent(in) :: restore_temp !< If true, temperature is restored to a target value.

type(mech_forcing), intent(inout) :: forces !< Driving mechanical forces
real(kind=8), intent(in) :: x2o(:,:) !< Fluxes from coupler to ocean, computed by ocean
type(cpl_indices_type), intent(inout) :: ind !< Structure with MCT attribute vectors and indices

! local variables
real, dimension(SZIB_(G),SZJB_(G)) :: &
Expand Down Expand Up @@ -269,8 +268,6 @@ subroutine convert_IOB_to_fluxes_and_forces(IOB, fluxes, Time, G, CS, &

real :: C_p ! heat capacity of seawater ( J/(K kg) )

real :: sw_vis_dir, sw_vis_dif, sw_nir_dir, sw_nir_dif

call cpu_clock_begin(id_clock_forcing)

is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec
Expand Down Expand Up @@ -447,21 +444,17 @@ subroutine convert_IOB_to_fluxes_and_forces(IOB, fluxes, Time, G, CS, &
k = k + 1 ! Increment position within gindex

if (wind_stagger == BGRID_NE) then
taux_at_q(I,J) = x2o(ind%x2o_Foxx_taux,k) * CS%wind_stress_multiplier
tauy_at_q(I,J) = x2o(ind%x2o_Foxx_tauy,k) * CS%wind_stress_multiplier
taux_at_q(I,J) = IOB%u_flux(i,j) * CS%wind_stress_multiplier
tauy_at_q(I,J) = IOB%v_flux(i,j) * CS%wind_stress_multiplier
! GMM, cime uses AGRID
elseif (wind_stagger == AGRID) then
taux_at_h(i,j) = x2o(ind%x2o_Foxx_taux,k) * CS%wind_stress_multiplier
tauy_at_h(i,j) = x2o(ind%x2o_Foxx_tauy,k) * CS%wind_stress_multiplier
taux_at_h(i,j) = IOB%u_flux(i,j) * CS%wind_stress_multiplier
tauy_at_h(i,j) = IOB%v_flux(i,j) * CS%wind_stress_multiplier
else ! C-grid wind stresses.
forces%taux(I,j) = x2o(ind%x2o_Foxx_taux,k) * CS%wind_stress_multiplier
forces%tauy(i,J) = x2o(ind%x2o_Foxx_tauy,k) * CS%wind_stress_multiplier
forces%taux(I,j) = IOB%u_flux(i,j) * CS%wind_stress_multiplier
forces%tauy(i,J) = IOB%v_flux(i,j) * CS%wind_stress_multiplier
endif

! NOTE: in convert_IOB_to_fluxes x2o below is replace by
! IOB%flux_quantity where flux_quantity is what we use to
! compute fluxes%flux_quantity

! liquid precipitation (rain)
if (associated(fluxes%lprec)) &
fluxes%lprec(i,j) = G%mask2dT(i,j) * IOB%lprec(i-i0,j-j0)
Expand All @@ -472,15 +465,15 @@ subroutine convert_IOB_to_fluxes_and_forces(IOB, fluxes, Time, G, CS, &

! evaporation
if (associated(fluxes%evap)) &
fluxes%evap(i,j) = x2o(ind%x2o_Foxx_evap,k) * G%mask2dT(i,j)
fluxes%evap(i,j) = G%mask2dT(i,j) * IOB%q_flux(i-i0,j-j0)

! river runoff flux
if (associated(fluxes%lrunoff)) &
fluxes%lrunoff(i,j) = x2o(ind%x2o_Foxx_rofl,k) * G%mask2dT(i,j)
fluxes%lrunoff(i,j) = G%mask2dT(i,j) * IOB%rofl_flux(i-i0,j-j0)

! ice runoff flux
if (associated(fluxes%frunoff)) &
fluxes%frunoff(i,j) = x2o(ind%x2o_Foxx_rofi,k) * G%mask2dT(i,j)
fluxes%frunoff(i,j) = G%mask2dT(i,j) * IOB%rofi_flux(i-i0,j-j0)

! GMM, we don't have an icebergs yet so the following is not needed
!if (((associated(IOB%ustar_berg) .and. (.not. associated(fluxes%ustar_berg))) &
Expand Down Expand Up @@ -529,7 +522,8 @@ subroutine convert_IOB_to_fluxes_and_forces(IOB, fluxes, Time, G, CS, &
! applied surface pressure from atmosphere and cryosphere
! sea-level pressure (Pa)
if (associated(forces%p_surf_full) .and. associated(forces%p_surf)) then
forces%p_surf_full(i,j) = G%mask2dT(i,j) * x2o(ind%x2o_Sa_pslv,k)
forces%p_surf_full(i,j) = G%mask2dT(i,j) * IOB%p(i-i0,j-j0)

if (CS%max_p_surf >= 0.0) then
forces%p_surf(i,j) = MIN(forces%p_surf_full(i,j),CS%max_p_surf)
else
Expand All @@ -541,16 +535,15 @@ subroutine convert_IOB_to_fluxes_and_forces(IOB, fluxes, Time, G, CS, &
else
forces%p_surf_SSH => forces%p_surf_full
endif

endif

! salt flux
! more salt restoring logic
if (associated(fluxes%salt_flux)) &
fluxes%salt_flux(i,j) = G%mask2dT(i,j)*(x2o(ind%x2o_Fioi_salt,k) + fluxes%salt_flux(i,j))
fluxes%salt_flux(i,j) = G%mask2dT(i,j)*(IOB%salt_flux(i-i0,j-j0) + fluxes%salt_flux(i,j))

if (associated(fluxes%salt_flux_in)) &
fluxes%salt_flux_in(i,j) = G%mask2dT(i,j)*x2o(ind%x2o_Fioi_salt,k)
fluxes%salt_flux_in(i,j) = G%mask2dT(i,j)*IOB%salt_flux(i-i0,j-j0)

enddo; enddo
! ############################ END OF MCT to MOM ##############################
Expand Down Expand Up @@ -1006,10 +999,12 @@ subroutine IOB_allocate(IOB, isc, iec, jsc, jec)
type(ice_ocean_boundary_type), intent(inout) :: IOB !< An ice-ocean boundary type with fluxes to drive
integer, intent(in) :: isc, iec, jsc, jec !< The ocean's local grid size

allocate ( IOB% u_flux (isc:iec,jsc:jec), &
allocate ( IOB% latent_flux (isc:iec,jsc:jec), &
IOB% rofl_flux (isc:iec,jsc:jec), &
IOB% rofi_flux (isc:iec,jsc:jec), &
IOB% u_flux (isc:iec,jsc:jec), &
IOB% v_flux (isc:iec,jsc:jec), &
IOB% t_flux (isc:iec,jsc:jec), &
IOB% latent_flux (isc:iec,jsc:jec), &
IOB% q_flux (isc:iec,jsc:jec), &
IOB% salt_flux (isc:iec,jsc:jec), &
IOB% lw_flux (isc:iec,jsc:jec), &
Expand All @@ -1029,10 +1024,12 @@ subroutine IOB_allocate(IOB, isc, iec, jsc, jec)
IOB% mi (isc:iec,jsc:jec), &
IOB% p (isc:iec,jsc:jec))

IOB%latent_flux = 0.0
IOB%rofl_flux = 0.0
IOB%rofi_flux = 0.0
IOB%u_flux = 0.0
IOB%v_flux = 0.0
IOB%t_flux = 0.0
IOB%latent_flux = 0.0
IOB%q_flux = 0.0
IOB%salt_flux = 0.0
IOB%lw_flux = 0.0
Expand Down
Loading

0 comments on commit 9722934

Please sign in to comment.