Skip to content

Commit

Permalink
Add probabilistic grounding scheme for landfast ice (#565)
Browse files Browse the repository at this point in the history 
* In the process of adding seabed2 parameterization

* Added calls foe icepack_query_parameters

* Modified modules used in subroutine to match latest code (it now compiles)

* Changed expression basal stress to seabed stress

* Replaced basal by better expression seabed at many places

* seabed everywhere...basalstress logical replaced by seabedstress

* Added kseabed in namelist for choosing seabed stress method

* Changed gacc to gravit (from icepack) and use of alphab from ice_in

* Cleaned up the seabed2 code, added comments

* Changing the doc for new grounding scheme

* Almost done with the doc

* Minor changes to the doc

* Modifs to the doc

* Value of sigmab was not the right one with respect to Dupont et el in prep

* For some rebase did not work completely...replced basalstress logical by seabedstress

* Other issue with the rebase for ice_init...fixed

* Minor mistake corrected...it compiles

* changed basalstress to seabedstress in options for tests

* Addressed some of the changes required for the PR

* Replace kseabed by character string for choice of method

* Updated the doc

* Cosmetic change ice_init

* Update doc/source/user_guide/ug_testing.rst

Co-authored-by: Philippe Blain <[email protected]>

* Cosmetic changes to ice_init to align prints

* Some modifs requested for PR: same taub after 4 months (gx1)

* Done with requested changes to new grounding method: same taub after 4 months (gx1)

* Modifs to doc and new test of prob method

Co-authored-by: Philippe Blain <[email protected]>
  • Loading branch information
@JFLemieux73 @phil-blain
JFLemieux73 and phil-blain authored Feb 22, 2021
1 parent d9a9f6d commit 3c516c8
Show file tree
Hide file tree
Showing 18 changed files with 658 additions and 335 deletions.
4 changes: 2 additions & 2 deletions cicecore/cicedynB/analysis/ice_history.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -925,12 +925,12 @@ subroutine init_hist (dt)
ns1, f_strinty)

call define_hist_field(n_taubx,"taubx","N/m^2",ustr2D, ucstr, &
"basal (seabed) stress (x)", &
"seabed (basal) stress (x)", &
"positive is x direction on U grid", c1, c0, &
ns1, f_taubx)

call define_hist_field(n_tauby,"tauby","N/m^2",ustr2D, ucstr, &
"basal (seabed) stress (y)", &
"seabed (basal) stress (y)", &
"positive is y direction on U grid", c1, c0, &
ns1, f_tauby)

Expand Down
31 changes: 23 additions & 8 deletions cicecore/cicedynB/dynamics/ice_dyn_eap.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -122,7 +122,8 @@ subroutine eap (dt)
use ice_dyn_shared, only: fcor_blk, ndte, dtei, &
denom1, uvel_init, vvel_init, arlx1i, &
dyn_prep1, dyn_prep2, stepu, dyn_finish, &
basal_stress_coeff, basalstress, &
seabed_stress_factor_LKD, seabed_stress_factor_prob, &
seabed_stress_method, seabed_stress, &
stack_velocity_field, unstack_velocity_field
use ice_flux, only: rdg_conv, strairxT, strairyT, &
strairx, strairy, uocn, vocn, ss_tltx, ss_tlty, iceumask, fm, &
Expand Down Expand Up @@ -383,17 +384,31 @@ subroutine eap (dt)
endif

!-----------------------------------------------------------------
! basal stress coefficients (landfast ice)
! seabed stress factor Tbu (Tbu is part of Cb coefficient)
!-----------------------------------------------------------------

if (basalstress) then
if (seabed_stress) then

!$OMP PARALLEL DO PRIVATE(iblk)
do iblk = 1, nblocks
call basal_stress_coeff (nx_block, ny_block, &
icellu (iblk), &
indxui(:,iblk), indxuj(:,iblk), &
vice(:,:,iblk), aice(:,:,iblk), &
hwater(:,:,iblk), Tbu(:,:,iblk))

if ( seabed_stress_method == 'LKD' ) then

call seabed_stress_factor_LKD (nx_block, ny_block, &
icellu (iblk), &
indxui(:,iblk), indxuj(:,iblk), &
vice(:,:,iblk), aice(:,:,iblk), &
hwater(:,:,iblk), Tbu(:,:,iblk))

elseif ( seabed_stress_method == 'probabilistic' ) then

call seabed_stress_factor_prob (nx_block, ny_block, &
icellt(iblk), indxti(:,iblk), indxtj(:,iblk), &
icellu(iblk), indxui(:,iblk), indxuj(:,iblk), &
aicen(:,:,:,iblk), vicen(:,:,:,iblk), &
hwater(:,:,iblk), Tbu(:,:,iblk))
endif

enddo
!$OMP END PARALLEL DO
endif
Expand Down
35 changes: 25 additions & 10 deletions cicecore/cicedynB/dynamics/ice_dyn_evp.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -40,8 +40,9 @@ module ice_dyn_evp
use ice_constants, only: c0, p027, p055, p111, p166, &
p222, p25, p333, p5, c1
use ice_dyn_shared, only: stepu, dyn_prep1, dyn_prep2, dyn_finish, &
ndte, yield_curve, ecci, denom1, arlx1i, fcor_blk, uvel_init, &
vvel_init, basal_stress_coeff, basalstress, Ktens, revp
ndte, yield_curve, ecci, denom1, arlx1i, fcor_blk, uvel_init, vvel_init, &
seabed_stress_factor_LKD, seabed_stress_factor_prob, seabed_stress_method, &
seabed_stress, Ktens, revp
use ice_fileunits, only: nu_diag
use ice_exit, only: abort_ice
use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted
Expand Down Expand Up @@ -325,17 +326,31 @@ subroutine evp (dt)
endif

!-----------------------------------------------------------------
! basal stress coefficients (landfast ice)
! seabed stress factor Tbu (Tbu is part of Cb coefficient)
!-----------------------------------------------------------------

if (basalstress) then
if (seabed_stress) then

!$OMP PARALLEL DO PRIVATE(iblk)
do iblk = 1, nblocks
call basal_stress_coeff (nx_block, ny_block, &
icellu (iblk), &
indxui(:,iblk), indxuj(:,iblk), &
vice(:,:,iblk), aice(:,:,iblk), &
hwater(:,:,iblk), Tbu(:,:,iblk))
do iblk = 1, nblocks

if ( seabed_stress_method == 'LKD' ) then

call seabed_stress_factor_LKD (nx_block, ny_block, &
icellu (iblk), &
indxui(:,iblk), indxuj(:,iblk), &
vice(:,:,iblk), aice(:,:,iblk), &
hwater(:,:,iblk), Tbu(:,:,iblk))

elseif ( seabed_stress_method == 'probabilistic' ) then

call seabed_stress_factor_prob (nx_block, ny_block, &
icellt(iblk), indxti(:,iblk), indxtj(:,iblk), &
icellu(iblk), indxui(:,iblk), indxuj(:,iblk), &
aicen(:,:,:,iblk), vicen(:,:,:,iblk), &
hwater(:,:,iblk), Tbu(:,:,iblk))
endif

enddo
!$OMP END PARALLEL DO
endif
Expand Down
10 changes: 5 additions & 5 deletions cicecore/cicedynB/dynamics/ice_dyn_evp_1d.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -818,7 +818,7 @@ subroutine stepu_iter(NA_len,rhow, &
real (kind=dbl_kind) :: tmp_str2_nw,tmp_str3_se,tmp_str4_sw, tmp_strintx
real (kind=dbl_kind) :: tmp_str6_se,tmp_str7_nw,tmp_str8_sw, tmp_strinty
real (kind=dbl_kind) :: waterx,watery
real (kind=dbl_kind) :: u0 = 5.e-5_dbl_kind ! residual velocity for basal stress (m/s)
real (kind=dbl_kind) :: u0 = 5.e-5_dbl_kind ! residual velocity for seabed stress (m/s)

character(len=*), parameter :: subname = '(stepu_iter)'
!---------------------------------------
Expand Down Expand Up @@ -881,7 +881,7 @@ subroutine stepu_last(NA_len, rhow, &

use ice_kinds_mod
use ice_constants, only: c0, c1
use ice_dyn_shared, only: brlx, revp, basalstress
use ice_dyn_shared, only: brlx, revp, seabed_stress

implicit none

Expand All @@ -908,7 +908,7 @@ subroutine stepu_last(NA_len, rhow, &
real (kind=dbl_kind) :: tmp_str2_nw,tmp_str3_se,tmp_str4_sw
real (kind=dbl_kind) :: tmp_str6_se,tmp_str7_nw,tmp_str8_sw
real (kind=dbl_kind) :: waterx,watery
real (kind=dbl_kind) :: u0 = 5.e-5_dbl_kind ! residual velocity for basal stress (m/s)
real (kind=dbl_kind) :: u0 = 5.e-5_dbl_kind ! residual velocity for seabed stress (m/s)

character(len=*), parameter :: subname = '(stepu_last)'
!---------------------------------------
Expand Down Expand Up @@ -954,8 +954,8 @@ subroutine stepu_last(NA_len, rhow, &
+ umassdti(iw)*(brlx*vold + revp*vvel_init(iw))
uvel(iw) = (cca*cc1 + ccb*cc2) / ab2
vvel(iw) = (cca*cc2 - ccb*cc1) / ab2
! calculate basal stress component for outputs
if ( basalstress ) then
! calculate seabed stress component for outputs
if ( seabed_stress ) then
taubx(iw) = -uvel(iw)*Tbu(iw) / (sqrt(uold**2 + vold**2) + u0)
tauby(iw) = -vvel(iw)*Tbu(iw) / (sqrt(uold**2 + vold**2) + u0)
endif
Expand Down
Loading

0 comments on commit 3c516c8

Please sign in to comment.