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

Jump to bottom

First round of housekeeping on ice_grid #921

Merged
merged 8 commits into from
Jan 11, 2024
Merged

First round of housekeeping on ice_grid #921

Show file tree
Hide file tree
Changes from 5 commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
0ef1ff6
removal of unused variables.
TillRasmussen Dec 15, 2023
32e7e26
moved xav to transport. Could remove commented code. Could remove xav…
TillRasmussen Dec 15, 2023
6eac1a6
Move derived parameters and only allocate if needed
TillRasmussen Dec 17, 2023
18e7161
bugfixes for cxp, cyp...
TillRasmussen Dec 19, 2023
461713f
Merge branch 'CICE-Consortium:main' into newgrid
TillRasmussen Dec 21, 2023
1cc2b24
fix index and remove commented code in ice_grid
TillRasmussen Jan 1, 2024
a618aef
new version of transport_remap. xav, yav array where needed. xxav, yy…
TillRasmussen Jan 1, 2024
f848367
Removed comments rom ice_transport_remap and arrays for nonuniform grids
TillRasmussen Jan 6, 2024
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
4 changes: 2 additions & 2 deletions cicecore/cicedyn/dynamics/ice_dyn_eap.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -109,7 +109,7 @@ subroutine eap (dt)
seabed_stress_factor_LKD, seabed_stress_factor_prob, &
seabed_stress_method, seabed_stress, &
stack_fields, unstack_fields, iceTmask, iceUmask, &
fld2, fld3, fld4
fld2, fld3, fld4, dxhy, dyhx, cxp, cyp, cxm, cym
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, &
Expand All @@ -118,7 +118,7 @@ subroutine eap (dt)
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, &
use ice_grid, only: tmask, umask, dxT, dyT, &
tarear, uarear, grid_average_X2Y, &
grid_atm_dynu, grid_atm_dynv, grid_ocn_dynu, grid_ocn_dynv
use ice_state, only: aice, aiU, vice, vsno, uvel, vvel, divu, shear, &
Expand Down
54 changes: 47 additions & 7 deletions cicecore/cicedyn/dynamics/ice_dyn_evp.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -84,6 +84,12 @@ module ice_dyn_evp
emass (:,:,:) , & ! total mass of ice and snow (E grid)
emassdti (:,:,:) ! mass of E-cell/dte (kg/m^2 s)

real (kind=dbl_kind), dimension (:,:,:), allocatable, public :: &
ratiodxN , & ! - dxN(i+1,j) / dxN(i,j)
ratiodyE , & ! - dyE(i ,j+1) / dyE(i,j)
ratiodxNr , & ! 1 / ratiodxN
ratiodyEr ! 1 / ratiodyE

real (kind=dbl_kind), allocatable :: &
strengthU(:,:,:) , & ! strength averaged to U points
divergU (:,:,:) , & ! div array on U points, differentiate from divu
Expand Down Expand Up @@ -119,9 +125,10 @@ module ice_dyn_evp
! Elastic-viscous-plastic dynamics driver
!
subroutine init_evp
use ice_blocks, only: nx_block, ny_block, nghost
use ice_domain_size, only: max_blocks, nx_global, ny_global
use ice_grid, only: grid_ice, dyT, dxT, uarear, tmask, G_HTE, G_HTN
use ice_blocks, only: get_block, nx_block, ny_block, nghost, block
use ice_domain_size, only: max_blocks
use ice_domain, only: nblocks, blocks_ice
use ice_grid, only: grid_ice, dyT, dxT, uarear, tmask, G_HTE, G_HTN, dxN, dyE
use ice_calendar, only: dt_dyn
use ice_dyn_shared, only: init_dyn_shared, evp_algorithm
use ice_dyn_evp1d, only: dyn_evp1d_init
Expand All @@ -131,6 +138,14 @@ subroutine init_evp

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

type (block) :: &
this_block ! block information for current block

integer (kind=int_kind) :: &
i, j, iblk , & ! block index
ilo,ihi,jlo,jhi ! beginning and end of physical domain


call init_dyn_shared(dt_dyn)

if (evp_algorithm == "shared_mem_1d" ) then
Expand Down Expand Up @@ -197,6 +212,33 @@ subroutine init_evp
stat=ierr)
if (ierr/=0) call abort_ice(subname//' ERROR: Out of memory E evp')

allocate( &
ratiodxN (nx_block,ny_block,max_blocks), &
ratiodyE (nx_block,ny_block,max_blocks), &
ratiodxNr(nx_block,ny_block,max_blocks), &
ratiodyEr(nx_block,ny_block,max_blocks), &
stat=ierr)
if (ierr/=0) call abort_ice(subname//'ERROR: Out of memory ratio')

!$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
ratiodxN (i,j,iblk) = - dxN(i+1,j ,iblk) / dxN(i,j,iblk)
ratiodyE (i,j,iblk) = - dyE(i ,j+1,iblk) / dyE(i,j,iblk)
ratiodxNr(i,j,iblk) = c1 / ratiodxN(i,j,iblk)
ratiodyEr(i,j,iblk) = c1 / ratiodyE(i,j,iblk)
enddo
enddo
enddo ! iblk
!$OMP END PARALLEL DO

endif

end subroutine init_evp
Expand All @@ -219,7 +261,7 @@ subroutine evp (dt)
ice_HaloDestroy, ice_HaloUpdate_stress
use ice_blocks, only: block, get_block, nx_block, ny_block, nghost
use ice_domain, only: nblocks, blocks_ice, halo_info, maskhalo_dyn
use ice_domain_size, only: max_blocks, ncat, nx_global, ny_global
use ice_domain_size, only: max_blocks, ncat
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, &
Expand All @@ -238,8 +280,6 @@ subroutine evp (dt)
stresspU, stressmU, stress12U
use ice_grid, only: tmask, umask, umaskCD, nmask, emask, uvm, epm, npm, &
dxE, dxN, dxT, dxU, dyE, dyN, dyT, dyU, &
ratiodxN, ratiodxNr, ratiodyE, ratiodyEr, &
dxhy, dyhx, cxp, cyp, cxm, cym, &
tarear, uarear, earear, narear, grid_average_X2Y, uarea, &
grid_type, grid_ice, &
grid_atm_dynu, grid_atm_dynv, grid_ocn_dynu, grid_ocn_dynv
Expand All @@ -250,7 +290,7 @@ subroutine evp (dt)
ice_timer_start, ice_timer_stop, timer_evp
use ice_dyn_shared, only: evp_algorithm, stack_fields, unstack_fields, &
DminTarea, visc_method, deformations, deformationsC_T, deformationsCD_T, &
strain_rates_U, &
strain_rates_U, dxhy, dyhx, cxp, cyp, cxm, cym, &
iceTmask, iceUmask, iceEmask, iceNmask, &
dyn_haloUpdate, fld2, fld3, fld4
use ice_dyn_evp1d, only: dyn_evp1d_run
Expand Down
86 changes: 79 additions & 7 deletions cicecore/cicedyn/dynamics/ice_dyn_shared.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@ module ice_dyn_shared

use ice_kinds_mod
use ice_communicate, only: my_task, master_task, get_num_procs
use ice_constants, only: c0, c1, c2, c3, c4, c6
use ice_constants, only: c0, c1, c2, c3, c4, c6, c1p5
use ice_constants, only: omega, spval_dbl, p01, p001, p5
use ice_blocks, only: nx_block, ny_block
use ice_domain_size, only: max_blocks
Expand Down Expand Up @@ -119,6 +119,14 @@ module ice_dyn_shared
real (kind=dbl_kind), allocatable, public :: &
DminTarea(:,:,:) ! deltamin * tarea (m^2/s)

real (kind=dbl_kind), dimension (:,:,:), allocatable, public :: &
cyp , & ! 1.5*HTE(i,j)-0.5*HTW(i,j) = 1.5*HTE(i,j)-0.5*HTE(i-1,j)
cxp , & ! 1.5*HTN(i,j)-0.5*HTS(i,j) = 1.5*HTN(i,j)-0.5*HTN(i,j-1)
cym , & ! 0.5*HTE(i,j)-1.5*HTW(i,j) = 0.5*HTE(i,j)-1.5*HTE(i-1,j)
cxm , & ! 0.5*HTN(i,j)-1.5*HTS(i,j) = 0.5*HTN(i,j)-1.5*HTN(i,j-1)
dxhy , & ! 0.5*(HTE(i,j) - HTW(i,j)) = 0.5*(HTE(i,j) - HTE(i-1,j))
dyhx ! 0.5*(HTN(i,j) - HTS(i,j)) = 0.5*(HTN(i,j) - HTN(i,j-1))

! ice isotropic tensile strength parameter
real (kind=dbl_kind), public :: &
Ktens ! T=Ktens*P (tensile strength: see Konig and Holland, 2010)
Expand Down Expand Up @@ -192,6 +200,23 @@ subroutine alloc_dyn_shared
stat=ierr)
if (ierr/=0) call abort_ice(subname//': Out of memory')

! if (evp_algorithm == "standard_2d") then
allocate( &
cyp (nx_block,ny_block,max_blocks), & ! 1.5*HTE - 0.5*HTW
cxp (nx_block,ny_block,max_blocks), & ! 1.5*HTN - 0.5*HTS
cym (nx_block,ny_block,max_blocks), & ! 0.5*HTE - 1.5*HTW
cxm (nx_block,ny_block,max_blocks), & ! 0.5*HTN - 1.5*HTS
stat=ierr)
if (ierr/=0) call abort_ice(subname//': Out of memory')

if (grid_ice == 'B' .and. evp_algorithm == "standard_2d") then
allocate( &
dxhy (nx_block,ny_block,max_blocks), & ! 0.5*(HTE - HTW)
dyhx (nx_block,ny_block,max_blocks), & ! 0.5*(HTN - HTS)
stat=ierr)
if (ierr/=0) call abort_ice(subname//': Out of memory')
endif

if (grid_ice == 'CD' .or. grid_ice == 'C') then
allocate( &
uvelE_init (nx_block,ny_block,max_blocks), & ! x-component of velocity (m/s), beginning of timestep
Expand All @@ -214,8 +239,9 @@ end subroutine alloc_dyn_shared

subroutine init_dyn_shared (dt)

use ice_blocks, only: nx_block, ny_block
use ice_domain, only: nblocks, halo_dynbundle
use ice_blocks, only: block, get_block
use ice_boundary, only: ice_halo, ice_haloUpdate
use ice_domain, only: nblocks, halo_dynbundle, blocks_ice, halo_info
use ice_domain_size, only: max_blocks
use ice_flux, only: &
stressp_1, stressp_2, stressp_3, stressp_4, &
Expand All @@ -224,17 +250,22 @@ subroutine init_dyn_shared (dt)
stresspT, stressmT, stress12T, &
stresspU, stressmU, stress12U
use ice_state, only: uvel, vvel, uvelE, vvelE, uvelN, vvelN
use ice_grid, only: ULAT, NLAT, ELAT, tarea
use ice_grid, only: ULAT, NLAT, ELAT, tarea, HTE, HTN
use ice_constants, only: field_loc_center, field_type_vector

real (kind=dbl_kind), intent(in) :: &
dt ! time step

! local variables

integer (kind=int_kind) :: &
i, j , & ! indices
nprocs, & ! number of processors
iblk ! block index
i, j , & ! indices
ilo, ihi, jlo, jhi, & !min and max index for interior of blocks
nprocs, & ! number of processors
iblk ! block index

type (block) :: &
this_block ! block information for current block

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

Expand Down Expand Up @@ -333,6 +364,47 @@ subroutine init_dyn_shared (dt)
enddo ! iblk
!$OMP END PARALLEL DO

if (grid_ice == 'B' .and. evp_algorithm == "standard_2d") then
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
dxhy(i,j,iblk) = p5*(HTE(i,j,iblk) - HTE(i-1,j,iblk))
dyhx(i,j,iblk) = p5*(HTN(i,j,iblk) - HTN(i,j-1,iblk))
enddo
enddo
enddo
call ice_HaloUpdate (dxhy, halo_info, &
field_loc_center, field_type_vector, &
fillValue=c1)
call ice_HaloUpdate (dyhx, halo_info, &
field_loc_center, field_type_vector, &
fillValue=c1)

endif

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+1
do i = ilo, ihi+1
cyp(i,j,iblk) = (c1p5*HTE(i,j,iblk) - p5*HTE(i-1,j,iblk))
cxp(i,j,iblk) = (c1p5*HTN(i,j,iblk) - p5*HTN(i,j-1,iblk))
! match order of operations in cyp, cxp for tripole grids
cym(i,j,iblk) = -(c1p5*HTE(i-1,j,iblk) - p5*HTE(i,j,iblk))
cxm(i,j,iblk) = -(c1p5*HTN(i,j-1,iblk) - p5*HTN(i,j,iblk))
enddo
enddo
enddo

end subroutine init_dyn_shared

!=======================================================================
Expand Down
26 changes: 9 additions & 17 deletions cicecore/cicedyn/dynamics/ice_dyn_vp.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -52,7 +52,7 @@ module ice_dyn_vp
use ice_fileunits, only: nu_diag
use ice_flux, only: fmU
use ice_global_reductions, only: global_sum
use ice_grid, only: dxT, dyT, dxhy, dyhx, cxp, cyp, cxm, cym, uarear
use ice_grid, only: dxT, dyT, uarear
use ice_exit, only: abort_ice
use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted
use icepack_intfc, only: icepack_ice_strength, icepack_query_parameters
Expand Down Expand Up @@ -120,19 +120,9 @@ subroutine init_vp
use ice_boundary, only: ice_HaloUpdate
use ice_constants, only: c1, &
field_loc_center, field_type_scalar
use ice_domain, only: blocks_ice, halo_info
use ice_domain, only: blocks_ice
use ice_calendar, only: dt_dyn
use ice_dyn_shared, only: init_dyn_shared
! use ice_grid, only: tarea

! local variables

integer (kind=int_kind) :: &
i, j, iblk, &
ilo,ihi,jlo,jhi ! beginning and end of physical domain

type (block) :: &
this_block ! block information for current block

call init_dyn_shared(dt_dyn)

Expand Down Expand Up @@ -167,7 +157,8 @@ subroutine implicit_solver (dt)
use ice_blocks, only: block, get_block, nx_block, ny_block
use ice_domain, only: blocks_ice, halo_info, maskhalo_dyn
use ice_domain_size, only: max_blocks, ncat
use ice_dyn_shared, only: deformations, iceTmask, iceUmask
use ice_dyn_shared, only: deformations, iceTmask, iceUmask, &
cxp, cyp, cxm, cym
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, &
Expand All @@ -176,7 +167,7 @@ subroutine implicit_solver (dt)
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, &
use ice_grid, only: tmask, umask, dxT, dyT, &
tarear, grid_type, grid_average_X2Y, &
grid_atm_dynu, grid_atm_dynv, grid_ocn_dynu, grid_ocn_dynv
use ice_state, only: aice, aiU, vice, vsno, uvel, vvel, divu, shear, &
Expand Down Expand Up @@ -686,9 +677,8 @@ subroutine anderson_solver (icellT , icellU , &
use ice_domain, only: maskhalo_dyn, halo_info
use ice_domain_size, only: max_blocks
use ice_flux, only: fmU, TbU
use ice_grid, only: dxT, dyT, dxhy, dyhx, cxp, cyp, cxm, cym, &
uarear
use ice_dyn_shared, only: DminTarea
use ice_grid, only: dxT, dyT, uarear
use ice_dyn_shared, only: DminTarea, dxhy, dyhx, cxp, cyp, cxm, cym
use ice_state, only: uvel, vvel, strength
use ice_timers, only: ice_timer_start, ice_timer_stop, timer_bound

Expand Down Expand Up @@ -2753,6 +2743,7 @@ subroutine fgmres (zetax2 , etax2 , &
use ice_boundary, only: ice_HaloUpdate
use ice_domain, only: maskhalo_dyn, halo_info
use ice_timers, only: ice_timer_start, ice_timer_stop, timer_bound
use ice_dyn_shared, only: dxhy, dyhx, cxp, cyp, cxm, cym

real (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks,4), intent(in) :: &
zetax2 , & ! zetax2 = 2*zeta (bulk viscosity)
Expand Down Expand Up @@ -3154,6 +3145,7 @@ subroutine pgmres (zetax2 , etax2 , &
use ice_boundary, only: ice_HaloUpdate
use ice_domain, only: maskhalo_dyn, halo_info
use ice_timers, only: ice_timer_start, ice_timer_stop, timer_bound
use ice_dyn_shared, only: dyhx, dxhy, cxp, cyp, cxm, cym

real (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks,4), intent(in) :: &
zetax2 , & ! zetax2 = 2*zeta (bulk viscosity)
Expand Down
Loading
Loading