Skip to content

Commit

Permalink
Changes in ice_dyn* are interpolation of uvelE/vvelN to B grid. ice_t… (
Browse files Browse the repository at this point in the history 
#17)

* Changes in ice_dyn* are interpolation of uvelE/vvelN to B grid. ice_transport files are changed in a way so that velocities are interpolated to b grid for depature point function and kept at E or N grid possible.

* changed according to comments. changed average from F to S and. commented out in vp and eap

* comment out grid_system, uvelE and vvelN
  • Loading branch information
@TillRasmussen
TillRasmussen authored Nov 17, 2021
1 parent 0aeee2c commit 838bf30
Show file tree
Hide file tree
Showing 5 changed files with 96 additions and 43 deletions.
13 changes: 10 additions & 3 deletions cicecore/cicedynB/dynamics/ice_dyn_eap.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -134,9 +134,9 @@ subroutine eap (dt)
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
tarear, uarear, grid_average_X2Y!, grid_system commented out until implementation of cd-grid
use ice_state, only: aice, vice, vsno, uvel, vvel, divu, shear, &
aice_init, aice0, aicen, vicen, strength
aice_init, aice0, aicen, vicen, strength !, uvelE, vvelN grid_system commented out until implementation of cd-grid
! use ice_timers, only: timer_dynamics, timer_bound, &
! ice_timer_start, ice_timer_stop, &
! timer_tmp1, timer_tmp2, timer_tmp3
Expand Down Expand Up @@ -558,7 +558,14 @@ subroutine eap (dt)
field_loc_NEcorner, field_type_vector)
call grid_average_X2Y('U2TF',strocnxT) ! shift
call grid_average_X2Y('U2TF',strocnyT)

! 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
! if (grid_system == 'CD') then
! call grid_average_X2Y('E2US',uvelE,uvel)
! call grid_average_X2Y('N2US',vvelN,vvel)
! endif
!end comment out
call ice_timer_stop(timer_dynamics) ! dynamics

end subroutine eap
Expand Down
6 changes: 5 additions & 1 deletion cicecore/cicedynB/dynamics/ice_dyn_evp.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -729,7 +729,11 @@ subroutine evp (dt)
field_loc_NEcorner, field_type_vector)
call grid_average_X2Y('U2TF',strocnxT) ! shift
call grid_average_X2Y('U2TF',strocnyT)

! shift velocity components from CD grid locations (N, E) to B grid location (U) for transport
if (grid_system == 'CD') then
call grid_average_X2Y('E2US',uvelE,uvel)
call grid_average_X2Y('N2US',vvelN,vvel)
endif
call ice_timer_stop(timer_dynamics) ! dynamics

end subroutine evp
Expand Down
13 changes: 10 additions & 3 deletions cicecore/cicedynB/dynamics/ice_dyn_vp.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -200,9 +200,9 @@ subroutine implicit_solver (dt)
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
tarear, grid_type, grid_average_X2Y !, grid_system commented out until implementation of c grid
use ice_state, only: aice, vice, vsno, uvel, vvel, divu, shear, &
aice_init, aice0, aicen, vicen, strength
aice_init, aice0, aicen, vicen, strength!, uvelE, vvelN ommented out until implementation of c grid
use ice_timers, only: timer_dynamics, timer_bound, &
ice_timer_start, ice_timer_stop

Expand Down Expand Up @@ -657,7 +657,14 @@ subroutine implicit_solver (dt)
field_loc_NEcorner, field_type_vector)
call grid_average_X2Y('U2TF',strocnxT) ! shift
call grid_average_X2Y('U2TF',strocnyT)

! 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
! if (grid_system == 'CD') then
! call grid_average_X2Y('E2US',uvelE,uvel)
! call grid_average_X2Y('N2US',vvelN,vvel)
! endif
!end comment out
call ice_timer_stop(timer_dynamics) ! dynamics

end subroutine implicit_solver
Expand Down
77 changes: 45 additions & 32 deletions cicecore/cicedynB/dynamics/ice_transport_driver.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -260,8 +260,8 @@ subroutine transport_remap (dt)
use ice_domain_size, only: ncat, max_blocks
use ice_blocks, only: nx_block, ny_block, block, get_block, nghost
use ice_state, only: aice0, aicen, vicen, vsnon, trcrn, &
uvel, vvel, bound_state
use ice_grid, only: tarea
uvel, vvel, bound_state, uvelE, vvelN
use ice_grid, only: tarea, grid_system
use ice_calendar, only: istep1
use ice_timers, only: ice_timer_start, ice_timer_stop, &
timer_advect, timer_bound
Expand Down Expand Up @@ -538,14 +538,24 @@ subroutine transport_remap (dt)
!-------------------------------------------------------------------
! Main remapping routine: Step ice area and tracers forward in time.
!-------------------------------------------------------------------

call horizontal_remap (dt, ntrace, &
if (grid_system == 'CD') then
call horizontal_remap (dt, ntrace, &
uvel (:,:,:), vvel (:,:,:), &
aim (:,:,:,:), trm (:,:,:,:,:), &
l_fixed_area, &
tracer_type, depend, &
has_dependents, integral_order, &
l_dp_midpt)
l_dp_midpt, grid_system, &
uvelE(:,:,:),vvelN(:,:,:))
else
call horizontal_remap (dt, ntrace, &
uvel (:,:,:), vvel (:,:,:), &
aim (:,:,:,:), trm (:,:,:,:,:), &
l_fixed_area, &
tracer_type, depend, &
has_dependents, integral_order, &
l_dp_midpt, grid_system)
endif

!-------------------------------------------------------------------
! Given new fields, recompute state variables.
Expand Down Expand Up @@ -709,8 +719,8 @@ subroutine transport_upwind (dt)
use ice_domain_size, only: ncat, max_blocks
use ice_state, only: aice0, aicen, vicen, vsnon, trcrn, &
uvel, vvel, trcr_depend, bound_state, trcr_base, &
n_trcr_strata, nt_strata
use ice_grid, only: HTE, HTN, tarea
n_trcr_strata, nt_strata, uvelE, vvelN
use ice_grid, only: HTE, HTN, tarea, grid_system
use ice_timers, only: ice_timer_start, ice_timer_stop, &
timer_bound, timer_advect

Expand Down Expand Up @@ -761,31 +771,34 @@ subroutine transport_upwind (dt)
!-------------------------------------------------------------------
! Average corner velocities to edges.
!-------------------------------------------------------------------

!$OMP PARALLEL DO PRIVATE(iblk,i,j,ilo,ihi,jlo,jhi,this_block)
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi

do j = jlo, jhi
do i = ilo, ihi
uee(i,j,iblk) = p5*(uvel(i,j,iblk) + uvel(i,j-1,iblk))
vnn(i,j,iblk) = p5*(vvel(i,j,iblk) + vvel(i-1,j,iblk))
enddo
enddo
enddo
!$OMP END PARALLEL DO

call ice_timer_start(timer_bound)
call ice_HaloUpdate (uee, halo_info, &
field_loc_Eface, field_type_vector)
call ice_HaloUpdate (vnn, halo_info, &
field_loc_Nface, field_type_vector)
call ice_timer_stop(timer_bound)

if (grid_system == 'CD') then
uee(:,:,:)=uvelE(:,:,:)
vnn(:,:,:)=vvelN(:,:,:)
else
!$OMP PARALLEL DO PRIVATE(iblk,i,j,ilo,ihi,jlo,jhi,this_block)
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi

do j = jlo, jhi
do i = ilo, ihi
uee(i,j,iblk) = p5*(uvel(i,j,iblk) + uvel(i,j-1,iblk))
vnn(i,j,iblk) = p5*(vvel(i,j,iblk) + vvel(i-1,j,iblk))
enddo
enddo
enddo
!$OMP END PARALLEL DO

call ice_timer_start(timer_bound)
call ice_HaloUpdate (uee, halo_info, &
field_loc_Eface, field_type_vector)
call ice_HaloUpdate (vnn, halo_info, &
field_loc_Nface, field_type_vector)
call ice_timer_stop(timer_bound)
endif
!$OMP PARALLEL DO PRIVATE(iblk,ilo,ihi,jlo,jhi,this_block)
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
Expand Down
30 changes: 26 additions & 4 deletions cicecore/cicedynB/dynamics/ice_transport_remap.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -319,7 +319,8 @@ subroutine horizontal_remap (dt, ntrace, &
tracer_type, depend, &
has_dependents, &
integral_order, &
l_dp_midpt)
l_dp_midpt, grid_system, &
uvelE, vvelN)

use ice_boundary, only: ice_halo, ice_HaloMask, ice_HaloUpdate, &
ice_HaloDestroy
Expand All @@ -339,15 +340,21 @@ subroutine horizontal_remap (dt, ntrace, &
ntrace ! number of tracers in use

real (kind=dbl_kind), intent(in), dimension(nx_block,ny_block,max_blocks) :: &
uvel ,&! x-component of velocity (m/s)
vvel ! y-component of velocity (m/s)
uvel ,&! x-component of velocity (m/s) ugrid
vvel ! y-component of velocity (m/s) ugrid

real (kind=dbl_kind), intent(in), optional, dimension(nx_block,ny_block,max_blocks) :: &
uvelE ,&! x-component of velocity (m/s) egrid
vvelN ! y-component of velocity (m/s) ngrid

real (kind=dbl_kind), intent(inout), dimension (nx_block,ny_block,0:ncat,max_blocks) :: &
mm ! mean mass values in each grid cell

real (kind=dbl_kind), intent(inout), dimension (nx_block,ny_block,ntrace,ncat,max_blocks) :: &
tm ! mean tracer values in each grid cell

character (len=char_len_long), intent(in) :: grid_system

!-------------------------------------------------------------------
! If l_fixed_area is true, the area of each departure region is
! computed in advance (e.g., by taking the divergence of the
Expand Down Expand Up @@ -663,6 +670,20 @@ subroutine horizontal_remap (dt, ntrace, &
enddo

if (l_fixed_area) then
if (grid_system == 'CD') then ! velocities are already on the center
do j = jlo, jhi
do i = ilo-1, ihi
edgearea_e(i,j) = uvelE(i,j,iblk) * HTE(i,j,iblk) * dt
enddo
enddo

do j = jlo-1, jhi
do i = ilo, ihi
edgearea_n(i,j) = vvelN(i,j,iblk)*HTN(i,j,iblk) * dt
enddo
enddo

else
do j = jlo, jhi
do i = ilo-1, ihi
edgearea_e(i,j) = (uvel(i,j,iblk) + uvel(i,j-1,iblk)) &
Expand All @@ -676,7 +697,8 @@ subroutine horizontal_remap (dt, ntrace, &
* p5 * HTN(i,j,iblk) * dt
enddo
enddo
endif
endif
endif

!-------------------------------------------------------------------
! Transports for east cell edges.
Expand Down

0 comments on commit 838bf30

Please sign in to comment.