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.

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

Refactor strocnxT, strocnyT implementation #764

Merged
merged 3 commits into from
Oct 10, 2022
Merged

Refactor strocnxT, strocnyT implementation #764

Show file tree
Hide file tree
Changes from 2 commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
84ce3b2
Refactor strocnxT, strocnyT implementation
apcraig Oct 5, 2022
98fc557
Update computation of strocnxT, strocnyT passed into icepack_step_therm1
apcraig Oct 6, 2022
de0f7cd
Rename strocn[x,y]T_sf to strocn[x,y]T_iavg
apcraig Oct 8, 2022
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
View file
Open in desktop
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
View file
Open in desktop
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
View file
Open in desktop
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
12 changes: 8 additions & 4 deletions cicecore/cicedynB/general/ice_flux.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -56,12 +56,14 @@ 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_sf, & ! ice-ocean stress, x-direction at T points, per ice fraction (scaled flux)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The "scaled flux" terminology is confusing. Scaled how? "per ice fraction" means that the quantity is the mean (of whatever) over just the ice fraction of the grid cell. That should be enough information, but if it's not clear, then say "averaged over the ice fraction". This comment implies that the stress terms have been divided by the ice fractional area, so that they are not averaged over the full grid cell area. Instead of _sf, maybe use _ice (although that's confusing too, with the history output standard of using _ai to mean the quantity includes a multiplicative factor of aice)

strocnyT_sf ! ice-ocean stress, y-direction at T points, per ice fraction (scaled flux)

! diagnostic

real (kind=dbl_kind), dimension (:,:,:), allocatable, public :: &
strocnxT, & ! ice-ocean stress, x-direction at T points, computed at end of dynamics
strocnyT, & ! ice-ocean stress, y-direction at T points, computed at end of dynamics
sig1 , & ! normalized principal stress component
sig2 , & ! normalized principal stress component
sigP , & ! internal ice pressure (N/m)
Expand Down Expand Up @@ -391,6 +393,8 @@ subroutine alloc_flux
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_sf(nx_block,ny_block,max_blocks), & ! ice-ocean stress, x-direction, per ice area
strocnyT_sf(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 +769,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_sf (:,:,:) = c0 ! ice-ocean stress, x-direction (T-cell)
strocnyT_sf (:,:,:) = 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
View file
Open in desktop
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