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

seabed stress - remove if statements #673

Merged
merged 9 commits into from
Feb 23, 2022
Merged

seabed stress - remove if statements #673

Show file tree
Hide file tree
Changes from 1 commit
Commits
Show all changes
9 commits
Select commit Hold shift + click to select a range
5c036cc
refactor seabed_stress. Bit for bit
TillRasmussen Dec 7, 2021
4ab3476
Removed if statement from stepu. Results are binary identical, howeve…
TillRasmussen Dec 12, 2021
0122dd3
changed capping from logical to numeric in order to remove if stateme…
TillRasmussen Dec 20, 2021
34b074c
clean dyn_finish, correct intent(inout) to intent(in) for Cw, resse C…
TillRasmussen Dec 22, 2021
dbbad0d
Merge branch 'main' into neweap
TillRasmussen Feb 19, 2022
22893ca
Reolve conflicts after updating main
TillRasmussen Feb 19, 2022
c833f90
modified environment for Freya to accomodate for additional OMP commands
TillRasmussen Feb 21, 2022
7376807
Requested changes after review. Only changed in seabed stress and not…
TillRasmussen Feb 21, 2022
e312c5b
move deformation to subcycling
TillRasmussen Feb 22, 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
2 changes: 0 additions & 2 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 @@ -541,8 +541,6 @@ subroutine eap (dt)
uvel (:,:,iblk), vvel (:,:,iblk), &
uocn (:,:,iblk), vocn (:,:,iblk), &
aiu (:,:,iblk), fm (:,:,iblk), &
strintx (:,:,iblk), strinty (:,:,iblk), &
strairx (:,:,iblk), strairy (:,:,iblk), &
strocnx (:,:,iblk), strocny (:,:,iblk), &
strocnxT(:,:,iblk), strocnyT(:,:,iblk))

Expand Down
2 changes: 0 additions & 2 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 @@ -616,8 +616,6 @@ subroutine evp (dt)
uvel (:,:,iblk), vvel (:,:,iblk), &
uocn (:,:,iblk), vocn (:,:,iblk), &
aiu (:,:,iblk), fm (:,:,iblk), &
strintx (:,:,iblk), strinty (:,:,iblk), &
strairx (:,:,iblk), strairy (:,:,iblk), &
strocnx (:,:,iblk), strocny (:,:,iblk), &
strocnxT(:,:,iblk), strocnyT(:,:,iblk))

Expand Down
43 changes: 20 additions & 23 deletions cicecore/cicedynB/dynamics/ice_dyn_shared.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -675,7 +675,7 @@ subroutine stepu (nx_block, ny_block, &
taubx , & ! seabed stress, x-direction (N/m^2)
tauby ! seabed stress, y-direction (N/m^2)

real (kind=dbl_kind), dimension (nx_block,ny_block), intent(inout) :: &
real (kind=dbl_kind), dimension (nx_block,ny_block), intent(in) :: &
Cw ! ocean-ice neutral drag coefficient
Comment on lines +678 to 679
Copy link
Member

Choose a reason for hiding this comment

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

nice cleanup


! local variables
Expand Down Expand Up @@ -741,8 +741,8 @@ subroutine stepu (nx_block, ny_block, &

! calculate seabed stress component for outputs
! only needed on last iteration.
taubx(i,j) = -uvel(i,j)*Tbu(i,j) / (sqrt(uold**2 + vold**2) + u0)
tauby(i,j) = -vvel(i,j)*Tbu(i,j) / (sqrt(uold**2 + vold**2) + u0)
taubx(i,j) = -uvel(i,j)*Cb
tauby(i,j) = -vvel(i,j)*Cb
Comment on lines +744 to +745
Copy link
Member

Choose a reason for hiding this comment

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

thanks, I thought we had made this change already somewhere...

enddo ! ij

end subroutine stepu
Expand All @@ -760,8 +760,6 @@ subroutine dyn_finish (nx_block, ny_block, &
uvel, vvel, &
uocn, vocn, &
aiu, fm, &
strintx, strinty, &
strairx, strairy, &
strocnx, strocny, &
strocnxT, strocnyT)

Expand All @@ -779,11 +777,7 @@ subroutine dyn_finish (nx_block, ny_block, &
uocn , & ! ocean current, x-direction (m/s)
vocn , & ! ocean current, y-direction (m/s)
aiu , & ! ice fraction on u-grid
fm , & ! Coriolis param. * mass in U-cell (kg/s)
strintx , & ! divergence of internal ice stress, x (N/m^2)
strinty , & ! divergence of internal ice stress, y (N/m^2)
strairx , & ! stress on ice by air, x-direction
strairy ! stress on ice by air, y-direction
Comment on lines -789 to -792
Copy link
Member

Choose a reason for hiding this comment

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

I'm not sure about those deletions. The computations for Hibler-Bryan stresses are still commented in the subroutine, so if we ever add them back and test them, we would have to add these four arguments back in.

So either we also remove the commented code, or we keep as-is. I do not think that having extra arguments impact performance, as Fortran uses pass-by-reference for arrays, So it should have no impact to have extra arguments passed...

Copy link
Contributor Author

@TillRasmussen TillRasmussen Feb 18, 2022 •
edited
Loading

Choose a reason for hiding this comment

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

I agree that they should be removed or reimplemented, however I dont have an opininion on which is the right to do, however the Hibler-Bryan stresses has been commented out as long as I remember. I would remove the commented code. @eclare108213 any opinion?

Copy link
Contributor

Choose a reason for hiding this comment

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

I would leave them there, but commented out.

fm ! Coriolis param. * mass in U-cell (kg/s)

real (kind=dbl_kind), dimension (nx_block,ny_block), intent(inout) :: &
strocnx , & ! ice-ocean stress, x-direction
Expand All @@ -793,14 +787,15 @@ subroutine dyn_finish (nx_block, ny_block, &
strocnxT, & ! ice-ocean stress, x-direction
strocnyT ! ice-ocean stress, y-direction

real (kind=dbl_kind), dimension (nx_block,ny_block), intent(in) :: &
Cw ! ocean-ice neutral drag coefficient

! local variables

integer (kind=int_kind) :: &
i, j, ij

real (kind=dbl_kind) :: vrel, rhow
real (kind=dbl_kind), dimension (nx_block,ny_block), intent(inout) :: &
Cw ! ocean-ice neutral drag coefficient

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

Expand Down Expand Up @@ -891,9 +886,10 @@ subroutine seabed_stress_factor_LKD (nx_block, ny_block, &
Tbu ! seabed stress factor (N/m^2)

real (kind=dbl_kind) :: &
au, & ! concentration of ice at u location
hu, & ! volume per unit area of ice at u location (mean thickness, m)
hwu, & ! water depth at u location (m)
au , & ! concentration of ice at u location
hu , & ! volume per unit area of ice at u location (mean thickness, m)
hwu , & ! water depth at u location (m)
docalc_tbu, & ! logical as real (C0,C1) decides whether c0 is 0 or
hcu ! critical thickness at u location (m)

integer (kind=int_kind) :: &
Expand All @@ -909,18 +905,19 @@ subroutine seabed_stress_factor_LKD (nx_block, ny_block, &

hwu = min(hwater(i,j),hwater(i+1,j),hwater(i,j+1),hwater(i+1,j+1))

if (hwu < threshold_hw) then
! if (hwu < threshold_hw) then
docalc_tbu = max(sign(c1,threshold_hw-hwu),c0)
Copy link
Member

Choose a reason for hiding this comment

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

maybe this would be easier to read and understand with the merge intrinsic (Intel doc, GFortran doc), which is like a ternary operator:

docalc_tbu = merge(c1, c0, hwu < threshold_hw)

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Done. Will be in a commit a bit later.


au = max(aice(i,j),aice(i+1,j),aice(i,j+1),aice(i+1,j+1))
hu = max(vice(i,j),vice(i+1,j),vice(i,j+1),vice(i+1,j+1))
au = max(aice(i,j),aice(i+1,j),aice(i,j+1),aice(i+1,j+1))
hu = max(vice(i,j),vice(i+1,j),vice(i,j+1),vice(i+1,j+1))

! 1- calculate critical thickness
hcu = au * hwu / k1
! 1- calculate critical thickness
hcu = au * hwu / k1

! 2- calculate seabed stress factor
Tbu(i,j) = k2 * max(c0,(hu - hcu)) * exp(-alphab * (c1 - au))
! 2- calculate seabed stress factor
Tbu(i,j) = docalc_tbu*k2 * max(c0,(hu - hcu)) * exp(-alphab * (c1 - au))

endif
! endif

enddo ! ij

Expand Down
2 changes: 0 additions & 2 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 @@ -642,8 +642,6 @@ subroutine implicit_solver (dt)
uvel (:,:,iblk), vvel (:,:,iblk), &
uocn (:,:,iblk), vocn (:,:,iblk), &
aiu (:,:,iblk), fm (:,:,iblk), &
strintx (:,:,iblk), strinty (:,:,iblk), &
strairx (:,:,iblk), strairy (:,:,iblk), &
strocnx (:,:,iblk), strocny (:,:,iblk), &
strocnxT(:,:,iblk), strocnyT(:,:,iblk))

Expand Down