Skip to content

Commit

Permalink
Added dxgrow, dygrow to facilitate variable spaced grid. Modified rec… (
Browse files Browse the repository at this point in the history
#746)

* Added dxgrow, dygrow to facilitate variable spaced grid. Modified rectgrid to generate grid from center outward using growth factor

* Adding vargrid namelist options.

* Refactored rectgrid to compute dx/dy first. Then ULON/ULAT. Added scale_dxdy flag to check of want grid spacing scaled. Renamed dxgrow/dygrow to dxscale/dyscale.

* Added method to check for odd nx_global/ny_global when applying grid spacing scale factors

* Update comments before computing dx/dy

* Update comments when checking for even/odd

* made grid_lonref, grid_latref namelist varaibles. Removed vargrid_suite.ts. Updated the box nml to specify the default grid_lonref and grid_latref for future reference.

* Change grid_lonref/grid_latref to lonrefrect,latrefrect. Reduce default vargrid tests to 3 per B,C,CD grid.

* Make new subroutine rectgrid_scale_dxdy to implemnet grid scaling. Remove explicit latrefrec/lonrefrect from set_nml. Make dxscale,dyscale,latrefrec,lonrefrec double precition in ice_in

* Add set_nml.scale1 to test vargrid with dxscale/dyscale = 1.d0

* Remove lonrefrec/latrefrec from boxnodyn

* Add lonrefrect, latrafrect to documentation

* Inserted new scaled grid varibles in alphabetical order in documentation
  • Loading branch information
daveh150 authored Aug 22, 2022
1 parent fea412a commit c87dcd3
Show file tree
Hide file tree
Showing 9 changed files with 310 additions and 69 deletions.
16 changes: 14 additions & 2 deletions cicecore/cicedynB/general/ice_init.F90
Original file line number Diff line number Diff line change
Expand Up @@ -110,8 +110,8 @@ subroutine input_data
grid_ice, grid_ice_thrm, grid_ice_dynu, grid_ice_dynv, &
grid_ocn, grid_ocn_thrm, grid_ocn_dynu, grid_ocn_dynv, &
grid_atm, grid_atm_thrm, grid_atm_dynu, grid_atm_dynv, &
dxrect, dyrect, &
pgl_global_ext
dxrect, dyrect, dxscale, dyscale, scale_dxdy, &
lonrefrect, latrefrect, pgl_global_ext
use ice_dyn_shared, only: ndte, kdyn, revised_evp, yield_curve, &
evp_algorithm, visc_method, &
seabed_stress, seabed_stress_method, &
Expand Down Expand Up @@ -212,6 +212,8 @@ subroutine input_data
bathymetry_file, use_bathymetry, nfsd, bathymetry_format, &
ncat, nilyr, nslyr, nblyr, &
kcatbound, gridcpl_file, dxrect, dyrect, &
dxscale, dyscale, lonrefrect, latrefrect, &
scale_dxdy, &
close_boundaries, orca_halogrid, grid_ice, kmt_type, &
grid_atm, grid_ocn

Expand Down Expand Up @@ -398,6 +400,11 @@ subroutine input_data
ksno = 0.3_dbl_kind ! snow thermal conductivity
dxrect = 0.0_dbl_kind ! user defined grid spacing in cm in x direction
dyrect = 0.0_dbl_kind ! user defined grid spacing in cm in y direction
lonrefrect = -156.50_dbl_kind ! lower left corner lon for rectgrid
latrefrect = 71.35_dbl_kind ! lower left corner lat for rectgrid
scale_dxdy = .false. ! apply dxscale, dyscale to rectgrid
dxscale = 1.0_dbl_kind ! user defined rectgrid x-grid scale factor (e.g., 1.02)
dyscale = 1.0_dbl_kind ! user defined rectgrid y-grid scale factor (e.g., 1.02)
close_boundaries = .false. ! true = set land on edges of grid
seabed_stress= .false. ! if true, seabed stress for landfast is on
seabed_stress_method = 'LKD'! LKD = Lemieux et al 2015, probabilistic = Dupont et al. in prep
Expand Down Expand Up @@ -853,6 +860,11 @@ subroutine input_data
call broadcast_scalar(grid_format, master_task)
call broadcast_scalar(dxrect, master_task)
call broadcast_scalar(dyrect, master_task)
call broadcast_scalar(scale_dxdy, master_task)
call broadcast_scalar(dxscale, master_task)
call broadcast_scalar(dyscale, master_task)
call broadcast_scalar(lonrefrect, master_task)
call broadcast_scalar(latrefrect, master_task)
call broadcast_scalar(close_boundaries, master_task)
call broadcast_scalar(grid_type, master_task)
call broadcast_scalar(grid_ice, master_task)
Expand Down
298 changes: 236 additions & 62 deletions cicecore/cicedynB/infrastructure/ice_grid.F90
Original file line number Diff line number Diff line change
Expand Up @@ -129,6 +129,15 @@ module ice_grid
dxrect, & ! user_specified spacing (cm) in x-direction (uniform HTN)
dyrect ! user_specified spacing (cm) in y-direction (uniform HTE)

! growth factor for variable spaced grid
real (kind=dbl_kind), public :: &
dxscale, & ! scale factor for grid spacing in x direction (e.g., 1.02)
dyscale ! scale factor for gird spacing in y direction (e.g., 1.02)

real (kind=dbl_kind), public :: &
lonrefrect, & ! lower left lon for rectgrid
latrefrect ! lower left lat for rectgrid

! Corners of grid boxes for history output
real (kind=dbl_kind), dimension (:,:,:,:), allocatable, public :: &
lont_bounds, & ! longitude of gridbox corners for T point
Expand Down Expand Up @@ -171,7 +180,8 @@ module ice_grid

logical (kind=log_kind), public :: &
use_bathymetry, & ! flag for reading in bathymetry_file
pgl_global_ext ! flag for init primary grid lengths (global ext.)
pgl_global_ext, & ! flag for init primary grid lengths (global ext.)
scale_dxdy ! flag to apply scale factor to vary dx/dy in rectgrid

logical (kind=log_kind), dimension (:,:,:), allocatable, public :: &
tmask , & ! land/boundary mask, thickness (T-cell)
Expand Down Expand Up @@ -1360,7 +1370,7 @@ subroutine rectgrid
imid, jmid

real (kind=dbl_kind) :: &
length, &
length, &
rad_to_deg

real (kind=dbl_kind), dimension(:,:), allocatable :: &
Expand All @@ -1383,69 +1393,71 @@ subroutine rectgrid

allocate(work_g1(nx_global,ny_global))

! Weddell Sea
! lower left corner of grid is 55W, 75S

! Barrow AK
! lower left corner of grid is 156.5W, 71.35N

if (my_task == master_task) then
work_g1 = c0
length = dxrect*cm_to_m/radius*rad_to_deg

! work_g1(1,:) = -55._dbl_kind ! Weddell Sea
work_g1(1,:) = -156.5_dbl_kind ! Barrow AK

do j = 1, ny_global
do i = 2, nx_global
work_g1(i,j) = work_g1(i-1,j) + length ! ULON
enddo
enddo
work_g1(:,:) = work_g1(:,:) / rad_to_deg
endif
call gridbox_verts(work_g1,lont_bounds)
call scatter_global(ULON, work_g1, master_task, distrb_info, &
field_loc_NEcorner, field_type_scalar)
call ice_HaloExtrapolate(ULON, distrb_info, &
ew_boundary_type, ns_boundary_type)

if (my_task == master_task) then
work_g1 = c0
length = dyrect*cm_to_m/radius*rad_to_deg

! work_g1(:,1) = -75._dbl_kind ! Weddell Sea
work_g1(:,1) = 71.35_dbl_kind ! Barrow AK
if (scale_dxdy) then
! scale grid spacing from center outward.
! this different than original method in it
! needs to define grid spacing before lat/lon.
! original rectgrid defines latlon first
call rectgrid_scale_dxdy
else
! rectgrid no grid spacing.
! original method with addition to use namelist lat/lon reference

if (my_task == master_task) then
work_g1 = c0
length = dxrect*cm_to_m/radius*rad_to_deg

work_g1(1,:) = lonrefrect ! reference lon from namelist

do j = 1, ny_global
do i = 2, nx_global
work_g1(i,j) = work_g1(i-1,j) + length ! ULON
enddo
enddo
work_g1(:,:) = work_g1(:,:) / rad_to_deg
endif
call scatter_global(ULON, work_g1, master_task, distrb_info, &
field_loc_NEcorner, field_type_scalar)
call ice_HaloExtrapolate(ULON, distrb_info, &
ew_boundary_type, ns_boundary_type)

if (my_task == master_task) then
work_g1 = c0
length = dyrect*cm_to_m/radius*rad_to_deg

work_g1(:,1) = latrefrect ! reference latitude from namelist

do i = 1, nx_global
do j = 2, ny_global
work_g1(i,j) = work_g1(i,j-1) + length ! ULAT
enddo
enddo
work_g1(:,:) = work_g1(:,:) / rad_to_deg
endif
call scatter_global(ULAT, work_g1, master_task, distrb_info, &
field_loc_NEcorner, field_type_scalar)
call ice_HaloExtrapolate(ULAT, distrb_info, &
ew_boundary_type, ns_boundary_type)

do i = 1, nx_global
do j = 2, ny_global
work_g1(i,j) = work_g1(i,j-1) + length ! ULAT
enddo
enddo
work_g1(:,:) = work_g1(:,:) / rad_to_deg
endif
call gridbox_verts(work_g1,latt_bounds)
call scatter_global(ULAT, work_g1, master_task, distrb_info, &
field_loc_NEcorner, field_type_scalar)
call ice_HaloExtrapolate(ULAT, distrb_info, &
ew_boundary_type, ns_boundary_type)
if (my_task == master_task) then
do j = 1, ny_global
do i = 1, nx_global
work_g1(i,j) = dxrect ! HTN
enddo
enddo
endif
call primary_grid_lengths_HTN(work_g1) ! dxU, dxT, dxN, dxE

if (my_task == master_task) then
do j = 1, ny_global
do i = 1, nx_global
work_g1(i,j) = dxrect ! HTN
enddo
enddo
endif
call primary_grid_lengths_HTN(work_g1) ! dxU, dxT, dxN, dxE
if (my_task == master_task) then
do j = 1, ny_global
do i = 1, nx_global
work_g1(i,j) = dyrect ! HTE
enddo
enddo
endif
call primary_grid_lengths_HTE(work_g1) ! dyU, dyT, dyN, dyE

if (my_task == master_task) then
do j = 1, ny_global
do i = 1, nx_global
work_g1(i,j) = dyrect ! HTE
enddo
enddo
endif
call primary_grid_lengths_HTE(work_g1) ! dyU, dyT, dyN, dyE
endif ! scale_dxdy

!-----------------------------------------------------------------
! Construct T-cell land mask
Expand Down Expand Up @@ -1526,6 +1538,168 @@ subroutine rectgrid

end subroutine rectgrid

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

subroutine rectgrid_scale_dxdy

! generate a variable spaced rectangluar grid.
! extend spacing from center of grid outward.
use ice_constants, only: c0, c1, c2, radius, cm_to_m, &
field_loc_center, field_loc_NEcorner, field_type_scalar

integer (kind=int_kind) :: &
i, j, iblk, &
imid, jmid, &
center1, center2 ! array centers for expanding dx, dy

real (kind=dbl_kind) :: &
length, &
rad_to_deg

real (kind=dbl_kind), dimension(:,:), allocatable :: &
work_g1

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

call icepack_query_parameters(rad_to_deg_out=rad_to_deg)

allocate(work_g1(nx_global,ny_global))

! determine dx spacing
! strategy: initialize with dxrect.
! if want to scale the grid, work from center outwards,
! multplying neighbor cell by scale factor.
! this assumes dx varies in x direction only.
! (i.e, dx is the same across same y location)
if (my_task == master_task) then

! initialize with initial dxrect
work_g1(:,:) = dxrect

! check if nx is even or odd
! if even, middle 2 columns are center
! of odd, middle 1 column is center
if (mod(nx_global,2) == 0) then ! nx_global is even

! with even number of x locatons,
! the center two y columns are center
center1 = nx_global/2 ! integer math
center2 = center1 + 1 ! integer math

else ! nx_global = odd
! only one center index. set center2=center1
center1 = ceiling(real(nx_global/2),int_kind)
center2 = center1
endif

! note loop over only half the x grid points (center1)-1
! working from the center outward.
do j = 1, ny_global
do i = 1, center1-1
! work from center1 to left
work_g1(center1-i,j) = dxscale*work_g1(center1-i+1,j)

! work from center2 to right
work_g1(center2+i,j) = dxscale*work_g1(center2+i-1,j)
enddo ! i
enddo ! j

endif ! my_task == master_task


! note work_g1 is converted to meters in primary_grid_lengths_HTN
call primary_grid_lengths_HTN(work_g1) ! dxU, dxT, dxN, dxE

! make ULON array
if (my_task == master_task) then

! make first column reference lon in radians.
! the remaining work_g1 is still dx in meters
work_g1(1,:) = lonrefrect/rad_to_deg ! radians

! loop over remaining points and add spacing to successive
! x locations
do j = 1, ny_global
do i = 2, nx_global ! start from i=2. i=1 is lonrefrect
length = work_g1(i,j)/radius ! grid spacing in radians
work_g1(i,j) = work_g1(i-1,j) + length ! ULON
enddo ! i
enddo ! j
endif ! mytask == master_task
call scatter_global(ULON, work_g1, master_task, distrb_info, &
field_loc_NEcorner, field_type_scalar)
call ice_HaloExtrapolate(ULON, distrb_info, &
ew_boundary_type, ns_boundary_type)

! determine dy spacing
! strategy: initialize with dyrect.
! if want to scale the grid, work from center outwards,
! multplying neighbor cell by scale factor.
! this assumes dy varies in y direction only.
! (i.e, dy is the same across same x location)
if (my_task == master_task) then

! initialize with initial dxrect
work_g1(:,:) = dyrect

! check if ny is even or odd
! if even, middle 2 rows are center
! of odd, middle 1 row is center
if (mod(ny_global,2) == 0) then ! ny_global is even

! with even number of x locatons,
! the center two y columns are center
center1 = ny_global/2 ! integer math
center2 = center1 + 1 ! integer math

else ! ny_global = odd
! only one center index. set center2=center1
center1 = ceiling(real(ny_global/2),int_kind)
center2 = center1
endif

! note loop over only half the y grid points (center1)-1
! working from the center outward.
do i = 1, nx_global
do j = 1, center1-1
! work from center1 to bottom
work_g1(i,center1-j) = dyscale*work_g1(i,center1-j+1)

! work from center2 to top
work_g1(i,center2+j) = dyscale*work_g1(i,center2+j-1)
enddo ! i
enddo ! j
endif ! mytask == master_task
! note work_g1 is converted to meters primary_grid_lengths_HTE
call primary_grid_lengths_HTE(work_g1) ! dyU, dyT, dyN, dyE

! make ULAT array
if (my_task == master_task) then

! make first row reference lat in radians.
! the remaining work_g1 is still dy in meters
work_g1(:,1) = latrefrect/rad_to_deg ! radians


! loop over remaining points and add spacing to successive
! x locations
do j = 2, ny_global ! start from j=2. j=1 is latrefrect
do i = 1, nx_global
length = work_g1(i,j)/radius ! grid spacing in radians
work_g1(i,j) = work_g1(i,j-1) + length ! ULAT
enddo ! i
enddo ! j
endif ! mytask == master_task
call scatter_global(ULAT, work_g1, master_task, distrb_info, &
field_loc_NEcorner, field_type_scalar)
call ice_HaloExtrapolate(ULAT, distrb_info, &
ew_boundary_type, ns_boundary_type)


deallocate(work_g1)

end subroutine rectgrid_scale_dxdy

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

! Complex land mask for testing box cases
Expand Down
4 changes: 2 additions & 2 deletions configuration/scripts/cice.batch.csh
Original file line number Diff line number Diff line change
Expand Up @@ -90,8 +90,8 @@ else if (${ICE_MACHINE} =~ nrlssc*) then
# nrlssc queue system has nodes with different task per node
if (${taskpernode} <= 12) set tpnstr = 'twelve'
if (${taskpernode} == 20) set tpnstr = 'twenty'
if (${taskpernode} == 24) set tpnstr = 'twentyfour'
if (${taskpernode} == 28) set tpnstr = 'twentyeight'
if (${taskpernode} >= 24) set tpnstr = 'twentyfour'
#if (${taskpernode} == 28) set tpnstr = 'twentyeight'

cat >> ${jobfile} <<EOFB
#PBS -N ${shortcase}
Expand Down
Loading

0 comments on commit c87dcd3

Please sign in to comment.