Skip to content
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.

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

Refactor strocnxT, strocnyT implementation #764

Merged
merged 3 commits into from
Oct 10, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
30 changes: 2 additions & 28 deletions cicecore/cicedynB/dynamics/ice_dyn_eap.F90
Original file line number Diff line number Diff line change
Expand Up @@ -138,15 +138,15 @@ subroutine eap (dt)
use ice_flux, only: rdg_conv, strairxT, strairyT, &
strairxU, strairyU, uocn, vocn, ss_tltx, ss_tlty, fmU, &
strtltxU, strtltyU, strocnxU, strocnyU, strintxU, strintyU, taubxU, taubyU, &
strocnxT, strocnyT, strax, stray, &
strax, stray, &
TbU, hwater, &
stressp_1, stressp_2, stressp_3, stressp_4, &
stressm_1, stressm_2, stressm_3, stressm_4, &
stress12_1, stress12_2, stress12_3, stress12_4
use ice_grid, only: tmask, umask, dxT, dyT, dxhy, dyhx, cxp, cyp, cxm, cym, &
tarear, uarear, grid_average_X2Y, iceumask, &
grid_atm_dynu, grid_atm_dynv, grid_ocn_dynu, grid_ocn_dynv
use ice_state, only: aice, vice, vsno, uvel, vvel, divu, shear, &
use ice_state, only: aice, aiU, vice, vsno, uvel, vvel, divu, shear, &
aice_init, aice0, aicen, vicen, strength
use ice_timers, only: timer_dynamics, timer_bound, &
ice_timer_start, ice_timer_stop
Expand Down Expand Up @@ -182,7 +182,6 @@ subroutine eap (dt)
wateryU , & ! for ocean stress calculation, y (m/s)
forcexU , & ! work array: combined atm stress and ocn tilt, x
forceyU , & ! work array: combined atm stress and ocn tilt, y
aiU , & ! ice fraction on u-grid
umass , & ! total mass of ice and snow (u grid)
umassdti ! mass of U-cell/dte (kg/m^2 s)

Expand All @@ -205,10 +204,6 @@ subroutine eap (dt)
type (block) :: &
this_block ! block information for current block

real (kind=dbl_kind), dimension (nx_block,ny_block,max_blocks) :: &
work1 , & ! temporary
work2 ! temporary

character(len=*), parameter :: subname = '(eap)'

call ice_timer_start(timer_dynamics) ! dynamics
Expand Down Expand Up @@ -567,27 +562,6 @@ subroutine eap (dt)
enddo
!$OMP END PARALLEL DO

! strocn computed on U, N, E as needed. Map strocn U divided by aiU to T
! TODO: This should be done elsewhere as part of generalization?
! conservation requires aiU be divided before averaging
work1 = c0
work2 = c0
!$OMP PARALLEL DO PRIVATE(iblk,ij,i,j)
do iblk = 1, nblocks
do ij = 1, icellu(iblk)
i = indxui(ij,iblk)
j = indxuj(ij,iblk)
work1(i,j,iblk) = strocnxU(i,j,iblk)/aiU(i,j,iblk)
work2(i,j,iblk) = strocnyU(i,j,iblk)/aiU(i,j,iblk)
enddo
enddo
call ice_HaloUpdate (work1, halo_info, &
field_loc_NEcorner, field_type_vector)
call ice_HaloUpdate (work2, halo_info, &
field_loc_NEcorner, field_type_vector)
call grid_average_X2Y('F', work1, 'U', strocnxT, 'T') ! shift
call grid_average_X2Y('F', work2, 'U', strocnyT, 'T')

call ice_timer_stop(timer_dynamics) ! dynamics

end subroutine eap
Expand Down
32 changes: 2 additions & 30 deletions cicecore/cicedynB/dynamics/ice_dyn_evp.F90
Original file line number Diff line number Diff line change
Expand Up @@ -85,7 +85,7 @@ subroutine evp (dt)
use ice_flux, only: rdg_conv, rdg_shear, strairxT, strairyT, &
strairxU, strairyU, uocn, vocn, ss_tltx, ss_tlty, fmU, &
strtltxU, strtltyU, strocnxU, strocnyU, strintxU, strintyU, taubxU, taubyU, &
strocnxT, strocnyT, strax, stray, &
strax, stray, &
TbU, hwater, &
strairxN, strairyN, fmN, &
strtltxN, strtltyN, strocnxN, strocnyN, strintxN, strintyN, taubxN, taubyN, &
Expand All @@ -106,7 +106,7 @@ subroutine evp (dt)
tarear, uarear, earear, narear, grid_average_X2Y, uarea, &
grid_type, grid_ice, &
grid_atm_dynu, grid_atm_dynv, grid_ocn_dynu, grid_ocn_dynv
use ice_state, only: aice, vice, vsno, uvel, vvel, uvelN, vvelN, &
use ice_state, only: aice, aiU, vice, vsno, uvel, vvel, uvelN, vvelN, &
uvelE, vvelE, divu, shear, &
aice_init, aice0, aicen, vicen, strength
use ice_timers, only: timer_dynamics, timer_bound, &
Expand Down Expand Up @@ -155,7 +155,6 @@ subroutine evp (dt)
wateryU , & ! for ocean stress calculation, y (m/s)
forcexU , & ! work array: combined atm stress and ocn tilt, x
forceyU , & ! work array: combined atm stress and ocn tilt, y
aiU , & ! ice fraction on u-grid
umass , & ! total mass of ice and snow (u grid)
umassdti ! mass of U-cell/dte (kg/m^2 s)

Expand Down Expand Up @@ -217,10 +216,6 @@ subroutine evp (dt)
type (block) :: &
this_block ! block information for current block

real (kind=dbl_kind), dimension (nx_block,ny_block,max_blocks) :: &
work1, & ! temporary
work2 ! temporary

logical (kind=log_kind), save :: &
first_time = .true. ! first time logical

Expand Down Expand Up @@ -1326,29 +1321,6 @@ subroutine evp (dt)

endif

! strocn computed on U, N, E as needed. Map strocn U divided by aiU to T
! TODO: This should be done elsewhere as part of generalization?
! TODO: Rename strocn[x,y]T since it's different than strocn[x,y][U,N,E]
! conservation requires aiU be divided before averaging
work1 = c0
work2 = c0
!$OMP PARALLEL DO PRIVATE(iblk,ij,i,j) SCHEDULE(runtime)
do iblk = 1, nblocks
do ij = 1, icellu(iblk)
i = indxui(ij,iblk)
j = indxuj(ij,iblk)
work1(i,j,iblk) = strocnxU(i,j,iblk)/aiU(i,j,iblk)
work2(i,j,iblk) = strocnyU(i,j,iblk)/aiU(i,j,iblk)
enddo
enddo
!$OMP END PARALLEL DO
call ice_HaloUpdate (work1, halo_info, &
field_loc_NEcorner, field_type_vector)
call ice_HaloUpdate (work2, halo_info, &
field_loc_NEcorner, field_type_vector)
call grid_average_X2Y('F', work1, 'U', strocnxT, 'T') ! shift
call grid_average_X2Y('F', work2, 'U', strocnyT, 'T')

if (grid_ice == 'CD' .or. grid_ice == 'C') then
call grid_average_X2Y('S', strintxE, 'E', strintxU, 'U') ! diagnostic
call grid_average_X2Y('S', strintyN, 'N', strintyU, 'U') ! diagnostic
Expand Down
30 changes: 2 additions & 28 deletions cicecore/cicedynB/dynamics/ice_dyn_vp.F90
Original file line number Diff line number Diff line change
Expand Up @@ -171,15 +171,15 @@ subroutine implicit_solver (dt)
use ice_flux, only: rdg_conv, rdg_shear, strairxT, strairyT, &
strairxU, strairyU, uocn, vocn, ss_tltx, ss_tlty, fmU, &
strtltxU, strtltyU, strocnxU, strocnyU, strintxU, strintyU, taubxU, taubyU, &
strocnxT, strocnyT, strax, stray, &
strax, stray, &
TbU, hwater, &
stressp_1, stressp_2, stressp_3, stressp_4, &
stressm_1, stressm_2, stressm_3, stressm_4, &
stress12_1, stress12_2, stress12_3, stress12_4
use ice_grid, only: tmask, umask, dxT, dyT, cxp, cyp, cxm, cym, &
tarear, grid_type, grid_average_X2Y, iceumask, &
grid_atm_dynu, grid_atm_dynv, grid_ocn_dynu, grid_ocn_dynv
use ice_state, only: aice, vice, vsno, uvel, vvel, divu, shear, &
use ice_state, only: aice, aiU, vice, vsno, uvel, vvel, divu, shear, &
aice_init, aice0, aicen, vicen, strength
use ice_timers, only: timer_dynamics, timer_bound, &
ice_timer_start, ice_timer_stop
Expand Down Expand Up @@ -210,7 +210,6 @@ subroutine implicit_solver (dt)
Cb , & ! seabed stress coefficient
fpresx , & ! fixed point residual vector, x components: fx = uvel - uprev_k
fpresy , & ! fixed point residual vector, y components: fy = vvel - vprev_k
aiU , & ! ice fraction on u-grid
umass , & ! total mass of ice and snow (u grid)
umassdti ! mass of U-cell/dte (kg/m^2 s)

Expand All @@ -234,10 +233,6 @@ subroutine implicit_solver (dt)
real (kind=dbl_kind), allocatable :: &
sol(:) ! solution vector

real (kind=dbl_kind), dimension (nx_block,ny_block,max_blocks) :: &
work1, & ! temporary
work2 ! temporary

character(len=*), parameter :: subname = '(implicit_solver)'

call ice_timer_start(timer_dynamics) ! dynamics
Expand Down Expand Up @@ -640,27 +635,6 @@ subroutine implicit_solver (dt)
enddo
!$OMP END PARALLEL DO

! strocn computed on U, N, E as needed. Map strocn U divided by aiU to T
! TODO: This should be done elsewhere as part of generalization?
! conservation requires aiU be divided before averaging
work1 = c0
work2 = c0
!$OMP PARALLEL DO PRIVATE(iblk,ij,i,j)
do iblk = 1, nblocks
do ij = 1, icellu(iblk)
i = indxui(ij,iblk)
j = indxuj(ij,iblk)
work1(i,j,iblk) = strocnxU(i,j,iblk)/aiU(i,j,iblk)
work2(i,j,iblk) = strocnyU(i,j,iblk)/aiU(i,j,iblk)
enddo
enddo
call ice_HaloUpdate (work1, halo_info, &
field_loc_NEcorner, field_type_vector)
call ice_HaloUpdate (work2, halo_info, &
field_loc_NEcorner, field_type_vector)
call grid_average_X2Y('F',work1,'U',strocnxT,'T') ! shift
call grid_average_X2Y('F',work2,'U',strocnyT,'T')

! shift velocity components from CD grid locations (N, E) to B grid location (U) for transport
! commented out in order to focus on EVP for now within the cdgrid
! should be used when routine is ready
Expand Down
22 changes: 16 additions & 6 deletions cicecore/cicedynB/general/ice_flux.F90
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,16 @@ module ice_flux
! Dynamics component
! All variables are assumed to be on the atm or ocn thermodynamic
! grid except as noted
!
! scale_fluxes divides several of these by aice "in place", so
! the state of some of these variables is not well defined. In the
! future, we need to refactor and add "_iavg" versions of the
! fields to clearly differentiate fields that have been divided
! by aice and others that are not. The challenge is that we need
! to go thru each field carefully to see which version is used.
! For instance, in diagnostics, there are places where these
! fields are multiplied by aice to compute things properly.
! strocn[x,y]T_iavg is the first field defined using _iavg.
!-----------------------------------------------------------------

real (kind=dbl_kind), dimension (:,:,:), allocatable, public :: &
Expand All @@ -56,8 +66,8 @@ module ice_flux

! out to ocean T-cell (kg/m s^2)
! Note, CICE_IN_NEMO uses strocnx and strocny for coupling
strocnxT, & ! ice-ocean stress, x-direction at T points, per ice fraction
strocnyT ! ice-ocean stress, y-direction at T points, per ice fraction
strocnxT_iavg, & ! ice-ocean stress, x-direction at T points, per ice fraction (scaled flux)
strocnyT_iavg ! ice-ocean stress, y-direction at T points, per ice fraction (scaled flux)

! diagnostic

Expand Down Expand Up @@ -389,8 +399,8 @@ subroutine alloc_flux
hwater (nx_block,ny_block,max_blocks), & ! water depth for seabed stress calc (landfast ice)
strairxT (nx_block,ny_block,max_blocks), & ! stress on ice by air, x-direction
strairyT (nx_block,ny_block,max_blocks), & ! stress on ice by air, y-direction
strocnxT (nx_block,ny_block,max_blocks), & ! ice-ocean stress, x-direction
strocnyT (nx_block,ny_block,max_blocks), & ! ice-ocean stress, y-direction
strocnxT_iavg(nx_block,ny_block,max_blocks), & ! ice-ocean stress, x-direction, per ice area
strocnyT_iavg(nx_block,ny_block,max_blocks), & ! ice-ocean stress, y-direction, per ice area
sig1 (nx_block,ny_block,max_blocks), & ! normalized principal stress component
sig2 (nx_block,ny_block,max_blocks), & ! normalized principal stress component
sigP (nx_block,ny_block,max_blocks), & ! internal ice pressure (N/m)
Expand Down Expand Up @@ -765,8 +775,8 @@ subroutine init_coupler_flux
! fluxes sent to ocean
!-----------------------------------------------------------------

strocnxT(:,:,:) = c0 ! ice-ocean stress, x-direction (T-cell)
strocnyT(:,:,:) = c0 ! ice-ocean stress, y-direction (T-cell)
strocnxT_iavg (:,:,:) = c0 ! ice-ocean stress, x-direction (T-cell)
strocnyT_iavg (:,:,:) = c0 ! ice-ocean stress, y-direction (T-cell)
fresh (:,:,:) = c0
fsalt (:,:,:) = c0
fpond (:,:,:) = c0
Expand Down
6 changes: 4 additions & 2 deletions cicecore/cicedynB/general/ice_state.F90
Original file line number Diff line number Diff line change
Expand Up @@ -54,7 +54,8 @@ module ice_state

real (kind=dbl_kind), dimension(:,:,:), allocatable, &
public :: &
aice , & ! concentration of ice
aice , & ! concentration of ice on T grid
aiU , & ! concentration of ice on U grid
vice , & ! volume per unit area of ice (m)
vsno ! volume per unit area of snow (m)

Expand Down Expand Up @@ -151,7 +152,8 @@ subroutine alloc_state
file=__FILE__, line=__LINE__)

allocate ( &
aice (nx_block,ny_block,max_blocks) , & ! concentration of ice
aice (nx_block,ny_block,max_blocks) , & ! concentration of ice T grid
aiU (nx_block,ny_block,max_blocks) , & ! concentration of ice U grid
vice (nx_block,ny_block,max_blocks) , & ! volume per unit area of ice (m)
vsno (nx_block,ny_block,max_blocks) , & ! volume per unit area of snow (m)
aice0 (nx_block,ny_block,max_blocks) , & ! concentration of open water
Expand Down
Loading