Skip to content

Wave-Ice Coupling for WW3v6.07 and E3SM v3 #7

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

Open
wants to merge 4 commits into
base: e3sm
Choose a base branch
from
Open

Wave-Ice Coupling for WW3v6.07 and E3SM v3 #7

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
1cc49e7
add IC4M10
erinethomas May 29, 2025
5aadc04
IC NUMERIC
erinethomas May 29, 2025
8575c2b
add wavice switch file
erinethomas May 29, 2025
9793971
method 8
erinethomas Jun 11, 2025
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
1 change: 1 addition & 0 deletions model/bin/switch_E3SM_wavice
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
F90 NOGRB NC4 DIST MPI PR3 UQ FLX0 LN1 ST4 STAB0 NL1 BT1 DB1 MLIM TR0 BS0 IC4 IS0 REF0 XX0 WNT1 WNX1 CRT1 CRX1 O0 O1 O2 O3 O4 O5 O6 O7 SCRIPNC SCRIP RTD RWND UOST
6 changes: 5 additions & 1 deletion model/src/w3fld1md.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1134,7 +1134,11 @@ SUBROUTINE APPENDTAIL(INSPC, WN2, NKT, KA1, KA2, KA3, WNDDIR,SAT)
ENDDO
AVG=SUM(NORMSPC)/MAX(REAL(NTH),1.)
DO T=1, NTH
INSPC(K,T) = SAT * NORMSPC(T)/TPI/(WN2(K)**3.0)/AVG
IF (AVG /= 0.0) THEN
INSPC(K,T)=BT(K)*INSPC(K,T)/TPI/(WN2(K)**3.0)/AVG
ELSE
INSPC(K,T)=0.0
ENDIF
ENDDO
ENDDO
DO T=1, NTH
Expand Down
5 changes: 4 additions & 1 deletion model/src/w3gdatmd.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -680,6 +680,7 @@ MODULE W3GDATMD

LOGICAL :: GINIT, FLDRY, FLCX, FLCY, FLCTH, FLCK, FLSOU, IICEDISP,&
IICESMOOTH
LOGICAL :: IC_NUMERICS
LOGICAL :: FLAGLL
LOGICAL :: CMPRTRCK
LOGICAL, POINTER :: FLAGST(:)
Expand Down Expand Up @@ -1182,6 +1183,7 @@ MODULE W3GDATMD

LOGICAL, POINTER :: GINIT, FLDRY, FLCX, FLCY, FLCTH, FLCK, FLSOU, IICEDISP,&
IICESMOOTH
LOGICAL, POINTER :: IC_NUMERICS
LOGICAL, POINTER :: FLAGLL
LOGICAL, POINTER :: CMPRTRCK
LOGICAL, POINTER :: FLAGST(:)
Expand Down Expand Up @@ -2359,7 +2361,8 @@ SUBROUTINE W3SETG ( IMOD, NDSE, NDST )
FLSOU => GRIDS(IMOD)%FLSOU
IICEDISP => GRIDS(IMOD)%IICEDISP
IICESMOOTH => GRIDS(IMOD)%IICESMOOTH
!
IC_NUMERICS => GRIDS(IMOD)%IC_NUMERICS
!
GNAME => GRIDS(IMOD)%GNAME
FILEXT => GRIDS(IMOD)%FILEXT
TRIGP => GRIDS(IMOD)%TRIGP
Expand Down
9 changes: 6 additions & 3 deletions model/src/w3gridmd.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -808,6 +808,7 @@ MODULE W3GRIDMD
!
REAL(8) :: GSHIFT ! see notes in WMGHGH
LOGICAL :: FLC, ICEDISP, TRCKCMPR
LOGICAL :: ICNUMERICS
INTEGER :: PTM ! Partitioning method
REAL :: PTFC ! Part. cut off freq (for method 5)
REAL :: AIRCMIN, AIRGB
Expand Down Expand Up @@ -1101,7 +1102,7 @@ MODULE W3GRIDMD
STDX, STDY, STDT, ICEHMIN, ICEHINIT, ICEDISP, &
ICESLN, ICEWIND, ICESNL, ICESDS, ICEHFAC, &
ICEHDISP, ICEDDISP, ICEFDISP, CALTYPE, &
TRCKCMPR, PTM, PTFC, BTBET
TRCKCMPR, PTM, PTFC, BTBET, ICNUMERICS
NAMELIST /OUTS/ P2SF, I1P2SF, I2P2SF, &
US3D, I1US3D, I2US3D, &
USSP, IUSSP, STK_WN, &
Expand Down Expand Up @@ -2744,6 +2745,7 @@ SUBROUTINE W3GRID(MDS)
STDY = -1.
STDT = -1.
ICEDISP = .FALSE.
ICNUMERICS=.FALSE.
CALTYPE = 'standard'
! Variables for 3D array output
E3D=0
Expand Down Expand Up @@ -3017,6 +3019,7 @@ SUBROUTINE W3GRID(MDS)
IICEHDISP = ICEHDISP
IICEDDISP = ICEDDISP
IICEFDISP = ICEFDISP
IC_NUMERICS=ICNUMERICS
PMOVE = MAX ( 0. , PMOVE )
PFMOVE = PMOVE
!
Expand Down Expand Up @@ -3404,7 +3407,7 @@ SUBROUTINE W3GRID(MDS)
ICEHINIT, ICEDISP, ICEHDISP, &
ICESLN, ICEWIND, ICESNL, ICESDS, &
ICEDDISP,ICEFDISP, CALTYPE, TRCKCMPR, &
BTBETA
BTBETA,ICNUMERICS
ELSE
WRITE (NDSO,2966) CICE0, CICEN, LICE, PMOVE, XSEED, FLAGTR, &
XP, XR, XFILT, IHMAX, HSPMIN, WSMULT, &
Expand All @@ -3414,7 +3417,7 @@ SUBROUTINE W3GRID(MDS)
ICEHINIT, ICEDISP, ICEHDISP, &
ICESLN, ICEWIND, ICESNL, ICESDS, &
ICEDDISP, ICEFDISP, CALTYPE, TRCKCMPR,&
BTBETA
BTBETA,ICNUMERICS
END IF
!
#ifdef W3_FLD1
Expand Down
4 changes: 2 additions & 2 deletions model/src/w3iogrmd.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -787,7 +787,7 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT )
FLCK, FLSOU, FLBPI, FLBPO, CLATS, CLATIS, CTHG0S, &
STEXU, STEYU, STEDU, IICEHMIN, IICEHINIT, IICEDISP, &
ICESCALES(1:4), CALTYPE, CMPRTRCK, IICEHFAC, IICEHDISP,&
IICEDDISP, IICEFDISP, BTBETA, &
IICEDDISP, IICEFDISP, BTBETA, IC_NUMERICS, &
AAIRCMIN, AAIRGB

WRITE(NDSM)GRIDSHIFT
Expand Down Expand Up @@ -980,7 +980,7 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT )
FLCTH, FLCK, FLSOU, FLBPI, FLBPO, CLATS, CLATIS, &
CTHG0S, STEXU, STEYU, STEDU, IICEHMIN, IICEHINIT, &
IICEDISP, ICESCALES(1:4), CALTYPE, CMPRTRCK, IICEHFAC, &
IICEDDISP, IICEHDISP, IICEFDISP, BTBETA, &
IICEDDISP, IICEHDISP, IICEFDISP, BTBETA, IC_NUMERICS, &
AAIRCMIN, AAIRGB
#ifdef W3_DEBUGIOGR
WRITE(740+IAPROC,*) 'W3IOGR, step 7.14'
Expand Down
40 changes: 39 additions & 1 deletion model/src/w3sic4md.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -322,7 +322,8 @@ SUBROUTINE W3SIC4 (A, DEPTH, CG, IX, IY, S, D)
REAL, ALLOCATABLE :: ALPHA(:) ! exponential decay rate for energy
REAL, ALLOCATABLE :: MARG1(:), MARG2(:) ! Arguments for M2
REAL, ALLOCATABLE :: KARG1(:), KARG2(:), KARG3(:) !Arguments for M3

REAL :: x1,x2,x3,x1sqr,x2sqr,x3sqr !Arguments for M10
REAL :: perfour,amhb,bmhb !Arguments for M10
!/
!/ ------------------------------------------------------------------- /
!/
Expand Down Expand Up @@ -513,6 +514,43 @@ SUBROUTINE W3SIC4 (A, DEPTH, CG, IX, IY, S, D)
ALPHA(IK) = 0.2*(FREQ**2.13)*HICE
END DO
WN_I= 0.5 * ALPHA

CASE (8)
! Cubic fit to Meylan, Horvat & Bitz 2021
! ICECOEF1 is thickness
! ICECOEF5 is floe size
! TPI/SIG is period
x3=min(ICECOEF1,3.5) ! limit thickness to 3.5 m
x3=max(x3,0.1) ! limit thickness >0.1 m since I make fit below
x2=min(ICECOEF5*0.5,100.0) ! convert dia to radius, limit to 100m
x2=max(2.5,x2)
x2sqr=x2*x2
x3sqr=x3*x3
amhb = 2.12e-3
bmhb = 4.59e-2

DO IK=1, NK
x1=TPI/SIG(IK) ! period
x1sqr=x1*x1
KARG1(ik)=-0.26982 + 1.5043*x3 - 0.70112*x3sqr + 0.011037*x2 + &
(-0.0073178)*x2*x3 + 0.00036604*x2*x3sqr + &
(-0.00045789)*x2sqr + 1.8034e-05*x2sqr*x3 + &
(-0.7246)*x1 + 0.12068*x1*x3 + &
(-0.0051311)*x1*x3sqr + 0.0059241*x1*x2 + &
0.00010771*x1*x2*x3 - 1.0171e-05*x1*x2sqr + &
0.0035412*x1sqr - 0.0031893*x1sqr*x3 + &
(-0.00010791)*x1sqr*x2 + &
0.00031073*x1**3 + 1.5996e-06*x2**3 + 0.090994*x3**3
KARG1(IK)=min(KARG1(IK),0.0)
ALPHA(IK) = 10.0**KARG1(IK)
perfour=x1sqr*x1sqr
if ((x1.gt.5.0) .and. (x1.lt.20.0)) then
ALPHA(IK) = ALPHA(IK) + amhb/x1sqr+bmhb/perfour
else if (x1.gt.20.0) then
ALPHA(IK) = amhb/x1sqr+bmhb/perfour
endif
WN_I(IK) = ALPHA(IK) * 0.5
end do

CASE DEFAULT
WN_I = ICECOEF1 !Default to IC1: Uniform in k
Expand Down
66 changes: 65 additions & 1 deletion model/src/w3srcemd.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -366,6 +366,7 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, &
XFC, XFLT, XREL, XFT, FXFM, FXPM, DDEN, &
FTE, FTF, FHMAX, ECOS, ESIN, IICEDISP, &
ICESCALES, IICESMOOTH
USE W3GDATMD, ONLY: IC_NUMERICS
USE W3GDATMD, ONLY: FSSOURCE, optionCall
USE W3GDATMD, ONLY: B_JGS_NLEVEL, B_JGS_SOURCE_NONLINEAR
#ifdef W3_REF1
Expand Down Expand Up @@ -1227,6 +1228,32 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, &
CALL UOST_SRCTRMCOMPUTE(IX, IY, SPEC, CG1, DT, &
U10ABS, U10DIR, VSUO, VDUO)
#endif
! Sea Ice Source Terms if IC_NUMERICS namelist flag = True
IF (IC_NUMERICS) THEN
#ifdef W3_IC1
IF (ICE .GT. 0) CALL W3SIC1 ( SPEC,DEPTH, CG1, IX, IY, VSIC, VDIC )
#endif
#ifdef W3_IS2
IF (ICE .GT. 0) CALL W3SIS2 ( SPEC, DEPTH, ICE, ICEH, ICEF, ICEDMAX, IX, IY, &
VSIR, VDIR, VDIR2, WN1, CG1, WN_R, CG_ICE, R )
#endif
#ifdef W3_IC2
IF (ICE .GT. 0) CALL W3SIC2 ( SPEC, DEPTH, ICEH, ICEF, CG1, WN1,&
IX, IY, VSIC, VDIC, WN_R, CG_ICE, ALPHA_LIU, R)
#endif
#ifdef W3_IC3
IF (ICE .GT. 0) CALL W3SIC3 ( SPEC,DEPTH, CG1, WN1, IX, IY, VSIC, VDIC )
#endif
#ifdef W3_IC4
IF (ICE .GT. 0) CALL W3SIC4 ( SPEC,DEPTH, CG1, IX, IY, VSIC, VDIC )
#endif
#ifdef W3_IC5
IF (ICE .GT. 0) CALL W3SIC5 ( SPEC,DEPTH, CG1, WN1, IX, IY, VSIC, VDIC )
#endif
#ifdef W3_IS1
IF (ICE .GT. 0) CALL W3SIS1 ( SPEC, ICE, VSIR )
#endif
ENDIF
!
! 2.g Dump training data if necessary
!
Expand Down Expand Up @@ -1288,6 +1315,12 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, &
VDIN(1:NSPECH) = ICESCALEIN * VDIN(1:NSPECH)
VSDS(1:NSPECH) = ICESCALEDS * VSDS(1:NSPECH)
VDDS(1:NSPECH) = ICESCALEDS * VDDS(1:NSPECH)
IF(IC_NUMERICS) THEN
#if defined(W3_IC1) || defined(W3_IC2) || defined(W3_IC3) || defined(W3_IC4) || defined(W3_IC5)
VSIC(1:NSPECH) = ICE * VSIC(1:NSPECH) ! (see Rogers et al 2016)
VDIC(1:NSPECH) = ICE * VDIC(1:NSPECH)
#endif
ENDIF
END IF
!
VS = 0
Expand All @@ -1307,6 +1340,11 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, &
#ifdef W3_UOST
VS(IS) = VS(IS) + VSUO(IS)
#endif
IF ( IC_NUMERICS .AND. ICE.GT.0. ) THEN
#if defined(W3_IC1) || defined(W3_IC2) || defined(W3_IC3) || defined(W3_IC4) || defined(W3_IC5)
VS(IS) = VS(IS) + VSIC(IS)
#endif
ENDIF
VD(IS) = VDIN(IS) + VDNL(IS) &
+ VDDS(IS) + VDBT(IS)
#ifdef W3_ST6
Expand All @@ -1321,6 +1359,11 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, &
#ifdef W3_UOST
VD(IS) = VD(IS) + VDUO(IS)
#endif
IF ( IC_NUMERICS .AND. ICE.GT.0. ) THEN
#if defined(W3_IC1) || defined(W3_IC2) || defined(W3_IC3) || defined(W3_IC4) || defined(W3_IC5)
VD(IS) = VD(IS) + VDIC(IS)
#endif
ENDIF
DAMAX = MIN ( DAM(IS) , MAX ( XREL*SPECINIT(IS) , AFILT ) )
AFAC = 1. / MAX( 1.E-10 , ABS(VS(IS)/DAMAX) )
#ifdef W3_NL5
Expand Down Expand Up @@ -1603,6 +1646,14 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, &
/ MAX ( 1. , (1.-HDT*VDBT(IS))) ! semi-implict integration scheme
PHINL = PHINL + VSNL(IS)* DT * FACTOR &
/ MAX ( 1. , (1.-HDT*VDNL(IS))) ! semi-implict integration scheme
IF ( IC_NUMERICS .AND. ICE.GT.0 ) THEN
#if defined(W3_IC1) || defined(W3_IC2) || defined(W3_IC3) || defined(W3_IC4) || defined(W3_IC5)
PHICE = PHICE + VSIC(IS) * DT * FACTOR &
/ MAX ( 1. , (1.-HDT*VDIC(IS))) ! semi-implicit integration
TAUICE(:) = TAUICE(:) - FACTOR2*COSI(:)*VSIC(IS) * DT &
/ MAX ( 1. , (1.-HDT*VDIC(IS)))
#endif
ENDIF
IF (VSIN(IS).GT.0.) WHITECAP(3) = WHITECAP(3) + SPEC(IS) * FACTOR
HSTOT = HSTOT + SPEC(IS) * FACTOR
END DO
Expand Down Expand Up @@ -1867,6 +1918,13 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, &
!
TAUOX=(GRAV*MWXFINISH+TAUWIX-TAUBBL(1))/DTG
TAUOY=(GRAV*MWYFINISH+TAUWIY-TAUBBL(2))/DTG
IF (IC_NUMERICS) THEN
#if defined(W3_IC1) || defined(W3_IC2) || defined(W3_IC3) || defined(W3_IC4) || defined(W3_IC5)
TAUICE(:)=TAUICE(:)/DTG
TAUOX = TAUOX - TAUICE(1)
TAUOY = TAUOY - TAUICE(2)
#endif
ENDIF
TAUWIX=TAUWIX/DTG
TAUWIY=TAUWIY/DTG
TAUWNX=TAUWNX/DTG
Expand All @@ -1881,6 +1939,11 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, &
PHIAW =DWAT*GRAV*PHIAW /DTG
PHINL =DWAT*GRAV*PHINL /DTG
PHIBBL=DWAT*GRAV*PHIBBL/DTG
IF (IC_NUMERICS) THEN
#if defined(W3_IC1) || defined(W3_IC2) || defined(W3_IC3) || defined(W3_IC4) || defined(W3_IC5)
PHICE =-1.*DWAT*GRAV*PHICE/DTG
#endif
ENDIF
!
! 10.1 Adds ice scattering and dissipation: implicit integration---------------- *
! INFLAGS2(4) is true if ice concentration was ever read during
Expand All @@ -1893,7 +1956,7 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, &
#endif

IF ( INFLAGS2(4).AND.ICE.GT.0 ) THEN

IF (.NOT. IC_NUMERICS ) THEN
IF (IICEDISP) THEN
ICECOEF2 = 1E-6
CALL LIU_FORWARD_DISPERSION (ICEH,ICECOEF2,DEPTH, &
Expand Down Expand Up @@ -2021,6 +2084,7 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, &
END DO
PHICE =-1.*DWAT*GRAV*PHICE /DTG
TAUICE(:)=TAUICE(:)/DTG
ENDIF ! end if IC_NUMERICS
ELSE
#ifdef W3_IS2
IF (IS2PARS(10).LT.0.5) THEN
Expand Down
Loading