Skip to content

Commit

Permalink
ww3_gse2: fix severall bugs and improve scheme for filtered points
Browse files Browse the repository at this point in the history
  • Loading branch information
aronroland committed Aug 3, 2024
1 parent 6a972bb commit afd8ae3
Show file tree
Hide file tree
Showing 4 changed files with 54 additions and 25 deletions.
13 changes: 10 additions & 3 deletions model/src/w3gridmd.F90
Original file line number Diff line number Diff line change
Expand Up @@ -2467,8 +2467,9 @@ SUBROUTINE W3GRID()
JGS_NLEVEL = 0
JGS_SOURCE_NONLINEAR = .FALSE.
! read data from the unstructured devoted namelist
CALL READNL ( NDSS, 'UNST', STATUS )

CALL READNL ( NDSS, 'UNST', STATUS )

B_JGS_USE_JACOBI = JGS_USE_JACOBI
B_JGS_TERMINATE_MAXITER = JGS_TERMINATE_MAXITER
B_JGS_TERMINATE_DIFFERENCE = JGS_TERMINATE_DIFFERENCE
Expand Down Expand Up @@ -3346,7 +3347,10 @@ SUBROUTINE W3GRID()
JGS_DIFF_THR, &
JGS_NORM_THR, &
JGS_NLEVEL, &
JGS_SOURCE_NONLINEAR
JGS_SOURCE_NONLINEAR, &
JGS_LGSE, &
JGS_GSE_METHOD, &
JGS_GSE_TS
!
WRITE (NDSO,2976) P2SF, I1P2SF, I2P2SF, &
US3D, I1US3D, I2US3D, &
Expand Down Expand Up @@ -6672,7 +6676,10 @@ SUBROUTINE W3GRID()
', JGS_DIFF_THR=', F8.3, &
', JGS_NORM_THR=', F8.3, &
', JGS_NLEVEL=', I3, &
', JGS_SOURCE_NONLINEAR=', L3 / )
', JGS_SOURCE_NONLINEAR=', L3, &
', JGS_LGSE=', L3, &
', JGS_GSE_METHOD=', I3, &
', JGS_GSE_TS=', F15.3/)
!
960 FORMAT (/' Miscellaneous ',A/ &
' --------------------------------------------------')
Expand Down
9 changes: 6 additions & 3 deletions model/src/w3iogrmd.F90
Original file line number Diff line number Diff line change
Expand Up @@ -796,7 +796,8 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT &
B_JGS_DIFF_THR, &
B_JGS_NORM_THR, &
B_JGS_NLEVEL, &
B_JGS_SOURCE_NONLINEAR
B_JGS_SOURCE_NONLINEAR, &
B_JGS_LGSE, B_JGS_GSE_TS, B_JGS_GSE_METHOD
#ifdef W3_ASCII
WRITE (NDSA,*) &
'FSN, FSPSI,FSFCT,FSNIMP,FSTOTALIMP,FSTOTALEXP, &
Expand Down Expand Up @@ -830,7 +831,8 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT &
B_JGS_DIFF_THR, &
B_JGS_NORM_THR, &
B_JGS_NLEVEL, &
B_JGS_SOURCE_NONLINEAR
B_JGS_SOURCE_NONLINEAR, &
B_JGS_LGSE, B_JGS_GSE_TS, B_JGS_GSE_METHOD
#endif
!Init COUNTCON and IOBDP to zero, it needs to be set somewhere or
!removed
Expand Down Expand Up @@ -995,7 +997,8 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT &
B_JGS_DIFF_THR, &
B_JGS_NORM_THR, &
B_JGS_NLEVEL, &
B_JGS_SOURCE_NONLINEAR
B_JGS_SOURCE_NONLINEAR, &
B_JGS_LGSE, B_JGS_GSE_TS, B_JGS_GSE_METHOD
IF (.NOT. GUGINIT) THEN
CALL W3DIMUG ( IGRD, NTRI, NX, COUNTOT, NNZ, NDSE, NDST )
END IF
Expand Down
53 changes: 34 additions & 19 deletions model/src/w3profsmd_pdlib.F90
Original file line number Diff line number Diff line change
Expand Up @@ -7722,7 +7722,7 @@ SUBROUTINE BLOCK_SOLVER_DIFFUSION(DTG)
INTEGER NB_ITER, iIter, ip_global
REAL*8 DeltaTmax, eDeltaT, CLATSMN, DFAC, RFAC, eDiffNorm
REAL*8 eNorm, DTquot, diffc, dcell, XWIND, DIFFTOT
REAL*8 DVDX(NPA), DVDY(NPA), DV2DXX(NPA), DV2DXY(NPA), KH, SWFAC
REAL*8 DVDX(NPA), DVDY(NPA), DV2DXX(NPA), DV2DXY(NPA), KH, SWFAC(NPA)
REAL eRealA(1), eRealB(1)


Expand All @@ -7747,23 +7747,26 @@ SUBROUTINE BLOCK_SOLVER_DIFFUSION(DTG)
DFAC = 1.
END IF
CLATMN = COS ( 70 * DERA )

DO ITH = 1, NTH
DO IK = 1, NK
ISP = ITH + (IK-1)*NTH
DO JSEA = 1, NPA
CALL INIT_GET_ISEA(ISEA, JSEA)
TFAC = MIN ( 1. , (CLATS(ISEA)/CLATMN)**2 )
CGD = 0.5 * GRAV / SIG(IK) * IOBDP_LOC(JSEA)
DSS = ( CGD * (XFR-1.))**2 * DTME / 12.! * TFAC
DNN = ( CGD * DTH )**2 * DTME / 12.! * TFAC
DCELL = CGD / 10.0 ! -> CELLP needs probably redifinition ...
KH = WN(IK,ISEA)*DW(JSEA)
IF (DW(JSEA) .gt. 1000.) THEN
SWFAC = 1.
IF (DW(ISEA) .gt. 1000.) THEN
SWFAC(JSEA) = 1.
ELSE
SWFAC = 0.!(1.-0.5*(1+(2*KH/SINH(MIN(20.d0,2*KH)))))**8
SWFAC(JSEA) = 0.!(1.-0.5*(1+(2*KH/SINH(MIN(20.d0,2*KH)))))**8
ENDIF
IF (SWFAC(JSEA) .LT. TINY(1.)) THEN
DIFFVEC(:,JSEA) = 0.
CYCLE
ENDIF
CGD = 0.5 * GRAV / SIG(IK) * IOBDP_LOC(JSEA)
DSS = ( CGD * (XFR-1.))**2 * DTME / 12.
DNN = ( CGD * DTH )**2 * DTME / 12.
DCELL = CGD / 10.0 ! -> CELLP needs probably redifinition ...
KH = WN(IK,ISEA)*DW(ISEA)
XWIND = 3.3 * U10(ISEA)*WN(IK,ISEA)/SIG(IK) - 2.3
XWIND = MAX ( 0. , MIN ( 1. , XWIND ) )
#ifdef W3_XW0
Expand All @@ -7773,16 +7776,16 @@ SUBROUTINE BLOCK_SOLVER_DIFFUSION(DTG)
XWIND = 1.
#endif
TFAC = MIN ( 1. , (CLATS(ISEA)/CLATMN)**2 )
DSS = SWFAC * (XWIND * DCELL + (1.-XWIND) * DSS * TFAC)
DSS = SWFAC(JSEA) * (XWIND * DCELL + (1.-XWIND) * DSS * TFAC)
#ifdef W3_DSS0
DSS = 0.
#endif
DNN = SWFAC * (XWIND * DCELL + (1.-XWIND) * DNN * TFAC)
DNN = SWFAC(JSEA) * (XWIND * DCELL + (1.-XWIND) * DNN * TFAC)

DIFFVEC(1,JSEA) = (DSS*ECOS(ITH)**2+DNN*ESIN(ITH)**2)
DIFFVEC(2,JSEA) = (DSS*ESIN(ITH)**2+DNN*ECOS(ITH)**2) / CLATS(ISEA)**2
DIFFVEC(3,JSEA) = ((DSS-DNN) * ESIN(ITH)*ECOS(ITH)) / CLATS(ISEA)
!write(3000+myrank,'(I10,20F20.10)') IK, SIG(IK), CGD, CLATS(ISEA), DSS, DNN, DIFFVEC(1,JSEA), DIFFVEC(2,JSEA)
!IF (DW(ISEA) .gt. 1000.) write(3000+myrank,'(I10,20F20.10)') IK, SIG(IK), CGD, CLATS(ISEA), DSS, DNN, DIFFVEC(1,JSEA), DIFFVEC(2,JSEA), DCELL, XWIND
END DO

CALL DIFFERENTIATE_XYDIR(DBLE(VA(ISP,:)),DVDX,DVDY) ! AR: this two lines can be optimized, however seems fast enough by now.
Expand All @@ -7791,7 +7794,10 @@ SUBROUTINE BLOCK_SOLVER_DIFFUSION(DTG)
DeltaTmax = 1./TINY(1.)
DO IE=1,NE
NI = INE(:,IE)
eDet = 2. * PDLIB_TRIA(IE)
IF (SUM(SWFAC(NI)) .EQ. 0) THEN
CYCLE
ENDIF
eDet = 2. * PDLIB_TRIA(IE)
DO IDX=1,3
V(1) = PDLIB_IEN(2*IDX-1,IE)
V(2) = PDLIB_IEN(2*IDX ,IE)
Expand Down Expand Up @@ -7819,12 +7825,15 @@ SUBROUTINE BLOCK_SOLVER_DIFFUSION(DTG)
DT_DIFF = DTG/NB_ITER
PHI_V = 0.

!WRITE(5000+myrank,*) 'NUMBER OF SUB ITERATIONS', ITH, IK, NB_ITER, DT_DIFF, DeltaTmax
WRITE(5000+myrank,*) 'NUMBER OF SUB ITERATIONS', ITH, IK, NB_ITER, DT_DIFF, DeltaTmax
!CALL FLUSH(5000+myrank)

DO IT = 1, NB_ITER
DO IE = 1, NE
NI = INE(:,IE)
IF (SUM(SWFAC(NI)) .EQ. 0) THEN
CYCLE
ENDIF
eDet = 2. * PDLIB_TRIA(IE)
DEDX(1) = PDLIB_IEN(1,IE)
DEDX(2) = PDLIB_IEN(3,IE)
Expand All @@ -7847,7 +7856,11 @@ SUBROUTINE BLOCK_SOLVER_DIFFUSION(DTG)
END DO
CALL PDLIB_exchange1DREAL(PHI_V)
DO JSEA =1, NSEAL
VA(ISP,JSEA) = MAX(0.,VA(ISP,JSEA) - (DT_DIFF * PHI_V(JSEA) / PDLIB_SI(JSEA) + 2 * DV2DXY(JSEA) * DIFFVEC(3,JSEA)) * DFAC * IOBDP_LOC(JSEA))
DIFFTOT = - (DT_DIFF * PHI_V(JSEA) / PDLIB_SI(JSEA) + DT_DIFF * 2 * DV2DXY(JSEA) * DIFFVEC(3,JSEA)) * DFAC * IOBDP_LOC(JSEA)
VA(ISP,JSEA) = MAX(0.,VA(ISP,JSEA) + DIFFTOT )
!IF (ABS(DIFFTOT) .GT. 0.) THEN
! WRITE(50000+myrank,'(2I10,10F20.10)') JSEA, ISP, VA(ISP,JSEA), DT_DIFF, PHI_V(JSEA), PDLIB_SI(JSEA), DV2DXY(JSEA), DIFFVEC(3,JSEA), DFAC, IOBDP_LOC(JSEA)
!ENDIF
END DO
END DO
END DO
Expand All @@ -7865,6 +7878,7 @@ SUBROUTINE DIFFERENTIATE_XYDIR(VAR, DVDX, DVDY)
USE yowExchangeModule, only : PDLIB_exchange1DREAL
USE yowNodepool, only : PDLIB_IEN, PDLIB_TRIA, NP, NPA
USE yowElementpool, only : NE, INE
USE W3PARALL, only: INIT_GET_ISEA
IMPLICIT NONE
REAL*8, INTENT(IN) :: VAR(NPA)
REAL*8, INTENT(INOUT) :: DVDX(NPA), DVDY(NPA)
Expand All @@ -7873,7 +7887,7 @@ SUBROUTINE DIFFERENTIATE_XYDIR(VAR, DVDX, DVDY)
REAL*8 :: DVDXIE, DVDYIE
REAL*8 :: WEI(NPA)
INTEGER :: NI(3)
INTEGER :: IE, JSEA
INTEGER :: IE, JSEA, ISEA

WEI(:) = 0.d0
DVDX(:) = 0.d0
Expand Down Expand Up @@ -7905,7 +7919,8 @@ SUBROUTINE DIFFERENTIATE_XYDIR(VAR, DVDX, DVDY)
DVDY = DVDY/WEI

DO JSEA = 1, NP
IF (DW(JSEA) .LT. DMIN) THEN
CALL INIT_GET_ISEA(ISEA, JSEA)
IF (DW(ISEA) .LT. DMIN) THEN
DVDX(JSEA) = 0.
DVDY(JSEA) = 0.
END IF
Expand Down
4 changes: 4 additions & 0 deletions model/src/w3wavemd.F90
Original file line number Diff line number Diff line change
Expand Up @@ -1812,6 +1812,10 @@ SUBROUTINE W3WAVE ( IMOD, ODAT, TEND, STAMP, NO_OUT &
END IF

IF (LPDLIB) THEN

!WRITE(*,*) B_JGS_LGSE, FSTOTALIMP, FSTOTALEXP
!STOP 'TEST TEST TEST'

!
#ifdef W3_PDLIB
IF (FLCX .or. FLCY) THEN
Expand Down

0 comments on commit afd8ae3

Please sign in to comment.