Skip to content

Commit

Permalink
UKMO Staging jan2021 (#327)
Browse files Browse the repository at this point in the history 
* First set of changes intended to fix the bug (#19)

Fixes: #314
* Interpolation weights now correctly calculated on points next to land and BC locations.
* Changes to improve the code: the possibility of reading zero values from
the input is considered, and points that should not be taken into
account in the interpolation are identified by the netcdf fill value; a
subroutine is created to avoid code duplication

* Bug fix and small simplification/optimization change (#18)

* Fixes #290 (ww3_multi hanging when generating restart with IOSTYP >= 2)
* Also fixes out-of-bounds array access error.
* Includes some MPI optimizations

* Correction to the bug fix in branch bf_multi_hang
to take into account the coupled
    configurations, that are also affected

* Small correction to the multi_hang branch: revert changes to JSEA index
in w3iorsmd

Co-authored-by: Juan Manuel Castillo Sanchez <[email protected]>
Co-authored-by: ukmo-juan.castillo <[email protected]>
  • Loading branch information
@ukmo-ccbunney @ukmo-juan-castillo
3 people authored Mar 19, 2021
1 parent c26cd2a commit 6d95c33
Show file tree
Hide file tree
Showing 5 changed files with 165 additions and 61 deletions.
6 changes: 5 additions & 1 deletion model/ftn/w3gsrumd.ftn
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3726,10 +3726,14 @@
IF ( K .EQ. 1 ) THEN
SIGN1 = SIGN(ONE,CROSS)
ELSE
IF ( SIGN(ONE,CROSS) .NE. SIGN1 ) THEN
! If point lies along a border, the cross product
! is zero and its sign is not well defined
IF ( ABS(CROSS) .GT. LEPS ) THEN
IF ( SIGN(ONE,CROSS) .NE. SIGN1 ) THEN
LSBC = .FALSE.
CYCLE SUBCELL_LOOP
END IF
END IF
END IF
END DO !K
IF ( LSBC ) RETURN
Expand Down
4 changes: 2 additions & 2 deletions model/ftn/w3odatmd.ftn
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1591,9 +1591,9 @@
!/MPI RSBLKS => OUTPTS(IMOD)%OUT4%RSBLKS
!/MPI IF ( NRQRS .NE. 0 ) THEN
!/MPI IRQRS => OUTPTS(IMOD)%OUT4%IRQRS
!/MPI IRQRSS => OUTPTS(IMOD)%OUT4%IRQRSS
!/MPI VAAUX => OUTPTS(IMOD)%OUT4%VAAUX
!/MPI END IF
!/MPI IRQRSS => OUTPTS(IMOD)%OUT4%IRQRSS
!/MPI VAAUX => OUTPTS(IMOD)%OUT4%VAAUX
!
NBI => OUTPTS(IMOD)%OUT5%NBI
NBI2 => OUTPTS(IMOD)%OUT5%NBI2
Expand Down
209 changes: 154 additions & 55 deletions model/ftn/ww3_prnc.ftn
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1649,7 +1649,6 @@
! forces to 0 values that are undefined
WHERE(XC.NE.XC) XC = FILLVALUE
WHERE (XC.NE.FILLVALUE) XC=XC*XCFAC+XCOFF
WHERE (XC.EQ.FILLVALUE) XC=0.

!
!/T2 WRITE (NDST,9060) 1
Expand Down Expand Up @@ -1683,7 +1682,6 @@
CALL CHECK_ERR(IRET)
WHERE(YC.NE.YC) YC = FILLVALUE
WHERE (YC.NE.FILLVALUE) YC=YC*YCFAC+YCOFF
WHERE (YC.EQ.FILLVALUE) YC=0.
!
!/T2 WRITE (NDST,9060) 2
!/T2 IXP0 = 1
Expand Down Expand Up @@ -1802,72 +1800,44 @@
!/O3 IF ( IAPROC .EQ. NAPOUT ) WRITE (NDSO,976) ' '
IF (( IFLD.LE.2 ).AND.( .NOT. FLBERG )) THEN
!
DO IY=1,NY
DO IX=1,NX
FA(IX,IY) &
= RD11(IX,IY) * XC(IX21(IX,IY),IY21(IX,IY)) &
+ RD21(IX,IY) * XC(IX22(IX,IY),IY21(IX,IY)) &
+ RD12(IX,IY) * XC(IX21(IX,IY),IY22(IX,IY)) &
+ RD22(IX,IY) * XC(IX22(IX,IY),IY22(IX,IY))
END DO
END DO
CALL INTERP(MXM, MYM, XC, IX21, IX22, IY21, IY22, &
RD11, RD12, RD21, RD22, FILLVALUE, FA)
!
IF (NFCOMP.EQ.2) THEN
!/O3 IF ( IAPROC .EQ. NAPOUT ) WRITE (NDSO,976) ' (2) '
DO IY=1,NY
DO IX=1,NX
FA(IX,IY) = FA(IX,IY) &
+ XD11(IX,IY) * YC(JX21(IX,IY),JY21(IX,IY)) &
+ XD21(IX,IY) * YC(JX22(IX,IY),JY21(IX,IY)) &
+ XD12(IX,IY) * YC(JX21(IX,IY),JY22(IX,IY)) &
+ XD22(IX,IY) * YC(JX22(IX,IY),JY22(IX,IY))
END DO
END DO
END IF
CALL INTERP(YC, JX21, JX22, JY21, JY22, XD11, XD12,&
XD21, XD22, FILLVALUE, FA)
END IF
!
! ... Two-component fields
!
ELSE !so if IFLD.GT.2
!
CALL INTERP(MXM, MYM, XC, IX21, IX22, IY21, IY22, &
RD11, RD12, RD21, RD22, FILLVALUE, FX)

CALL INTERP(MXM, MYM, YC, IX21, IX22, IY21, IY22, &
RD11, RD12, RD21, RD22, FILLVALUE, FY)

CALL INTERP(MXM, MYM, AC, IX21, IX22, IY21, IY22, &
RD11, RD12, RD21, RD22, FILLVALUE, FA)

CALL INTERP(MXM, MYM, SQRT(XC**2+YC**2), IX21, IX22, &
IY21, IY22, RD11, RD12, RD21, RD22, &
SQRT(2.0)*FILLVALUE, A2)

CALL INTERP(MXM, MYM, XC**2+YC**2, IX21, IX22, &
IY21, IY22, RD11, RD12, RD21, RD22, &
2.0*FILLVALUE*FILLVALUE, A3)

DO IY=1,NY
DO IX=1,NX
FX(IX,IY) &
= RD11(IX,IY) * XC(IX21(IX,IY),IY21(IX,IY)) &
+ RD21(IX,IY) * XC(IX22(IX,IY),IY21(IX,IY)) &
+ RD12(IX,IY) * XC(IX21(IX,IY),IY22(IX,IY)) &
+ RD22(IX,IY) * XC(IX22(IX,IY),IY22(IX,IY))
FY(IX,IY) &
= RD11(IX,IY) * YC(IX21(IX,IY),IY21(IX,IY)) &
+ RD21(IX,IY) * YC(IX22(IX,IY),IY21(IX,IY)) &
+ RD12(IX,IY) * YC(IX21(IX,IY),IY22(IX,IY)) &
+ RD22(IX,IY) * YC(IX22(IX,IY),IY22(IX,IY))
FA(IX,IY) &
= RD11(IX,IY) * AC(IX21(IX,IY),IY21(IX,IY)) &
+ RD21(IX,IY) * AC(IX22(IX,IY),IY21(IX,IY)) &
+ RD12(IX,IY) * AC(IX21(IX,IY),IY22(IX,IY)) &
+ RD22(IX,IY) * AC(IX22(IX,IY),IY22(IX,IY))
A1(IX,IY) = MAX ( 1.E-10 , &
SQRT( FX(IX,IY)**2 + FY(IX,IY)**2 ) )
A2(IX,IY) &
= RD11(IX,IY) * SQRT(XC(IX21(IX,IY),IY21(IX,IY))**2 &
+YC(IX21(IX,IY),IY21(IX,IY))**2) &
+ RD21(IX,IY) * SQRT(XC(IX22(IX,IY),IY21(IX,IY))**2 &
+YC(IX22(IX,IY),IY21(IX,IY))**2) &
+ RD12(IX,IY) * SQRT(XC(IX21(IX,IY),IY22(IX,IY))**2 &
+YC(IX21(IX,IY),IY22(IX,IY))**2) &
+ RD22(IX,IY) * SQRT(XC(IX22(IX,IY),IY22(IX,IY))**2 &
+YC(IX22(IX,IY),IY22(IX,IY))**2)
A3(IX,IY) = SQRT ( &
RD11(IX,IY) * ( XC(IX21(IX,IY),IY21(IX,IY))**2 &
+ YC(IX21(IX,IY),IY21(IX,IY))**2 ) &
+ RD21(IX,IY) * ( XC(IX22(IX,IY),IY21(IX,IY))**2 &
+ YC(IX22(IX,IY),IY21(IX,IY))**2 ) &
+ RD12(IX,IY) * ( XC(IX21(IX,IY),IY22(IX,IY))**2 &
+ YC(IX21(IX,IY),IY22(IX,IY))**2 ) &
+ RD22(IX,IY) * ( XC(IX22(IX,IY),IY22(IX,IY))**2 &
+ YC(IX22(IX,IY),IY22(IX,IY))**2 ) )
END DO

A3(IX,IY) = SQRT( A3(IX,IY) )
END DO
END DO
!
! ... Winds, correct for velocity or energy conservation
!
Expand Down Expand Up @@ -2238,6 +2208,135 @@

END PROGRAM W3PRNC

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

SUBROUTINE INTERP(MXM, MYM, XC, IX21, IX22, IY21, IY22, &
RD11, RD12, RD21, RD22, FILLVALUE, FA)
!/
!/ +-----------------------------------+
!/ | WAVEWATCH III NOAA/NCEP |
!/ | J. M. Castillo |
!/ | FORTRAN 90 |
!/ | Last update : 23-Feb-2021 |
!/ +-----------------------------------+
!/
!/ 23-Feb-2021 : First version ( version 7.xx )
!/
! 1. Purpose :
!
! Interpolate from a field read from file to the wave grid
!
! 2. Method :
!
! Invalid points are identified by the fill value read from the
! netcdf input, and interpolation does not take into account
! these points. The valid interpolation coefficients are scaled
! so that the sum is one, otherwise unphysical values can be
! generated.
!
! When one point is on the boundary but is not an ocean grid point,
! the interpolation coefficients are zero, and in this case we
! provide a sensible value - the value as read, not interpolated
!
! 3. Parameters :
!
! Parameter list
! ----------------------------------------------------------------
! MxM I I Dimensions of the XC variable
! XC R.A. I Field to be interpolated, as read from the
! input netcdf
! IXxx I.A. I List of x-index to convert from the original
! field to the model grid
! IYxx I.A. I List of y-index to convert from the original
! field to the model grid
! RDxx R.A. I Interpolation factors
! FILLVALUE R I Fill value identifying non valid input
! FA F O Result of the interpolation
! ----------------------------------------------------------------
!
! 4. Subroutines used :
!
! None
!
! 5. Called by :
!
! Name Type Module Description
! ----------------------------------------------------------------
! WW3_PRNC Prog. N/A Input data preprocessor.
! ----------------------------------------------------------------
!
! 6. Error messages :
!
! None
!
! 7. Remarks :
!
! 8. Structure :
!
! See source code.
!
! 9. Switches :
!
! 10. Source code :
!
!/ ------------------------------------------------------------------- /
USE W3GDATMD, ONLY: NX, NY

IMPLICIT NONE
!/
!/ ------------------------------------------------------------------- /
!/ Parameter list
!/
INTEGER, INTENT(IN) :: MXM, MYM
REAL, DIMENSION(MXM,MYM), INTENT(IN) :: XC
INTEGER, DIMENSION(NX,NY), INTENT(IN) :: IX21, IX22, IY21, IY22
REAL, DIMENSION(NX,NY), INTENT(IN) :: RD11, RD12, RD21, RD22
REAL, INTENT(IN) :: FILLVALUE
REAL, DIMENSION(NX,NY), INTENT(OUT) :: FA
!/
!/ ------------------------------------------------------------------- /
!/ Local variables
!/
INTEGER :: IX, IY
REAL :: FACTOR
!/ ------------------------------------------------------------------- /

DO IY=1,NY
DO IX=1,NX
FACTOR = 0.0
FA(IX,IY) = 0.0

IF(XC(IX21(IX,IY),IY21(IX,IY)).NE.FILLVALUE) THEN
FACTOR = FACTOR + RD11(IX,IY)
FA(IX,IY) = RD11(IX,IY) * XC(IX21(IX,IY),IY21(IX,IY))
ENDIF
IF(XC(IX22(IX,IY),IY21(IX,IY)).NE.FILLVALUE) THEN
FACTOR = FACTOR + RD21(IX,IY)
FA(IX,IY) = FA(IX,IY) + RD21(IX,IY) * XC(IX22(IX,IY),IY21(IX,IY))
ENDIF
IF(XC(IX21(IX,IY),IY22(IX,IY)).NE.FILLVALUE) THEN
FACTOR = FACTOR + RD12(IX,IY)
FA(IX,IY) = FA(IX,IY) + RD12(IX,IY) * XC(IX21(IX,IY),IY22(IX,IY))
ENDIF
IF(XC(IX22(IX,IY),IY22(IX,IY)).NE.FILLVALUE) THEN
FACTOR = FACTOR + RD22(IX,IY)
FA(IX,IY) = FA(IX,IY) + RD22(IX,IY) * XC(IX22(IX,IY),IY22(IX,IY))
ENDIF

IF(FACTOR.GT.0.0) THEN
FA(IX,IY) = FA(IX,IY) / FACTOR
ELSE
IF(XC(IX,IY).EQ.FILLVALUE) THEN
FA(IX,IY) = 0.0
ELSE
FA(IX,IY) = XC(IX,IY)
ENDIF
END IF
END DO
END DO

END SUBROUTINE INTERP

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

SUBROUTINE CHECK_ERR(IRET)
Expand Down
4 changes: 2 additions & 2 deletions regtests/mww3_test_04/input/ww3_multi_grdset_d.inp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,7 @@ $
19680606 000000 1200 19680608 000000
$
N
DPT CUR HS FP DP DIR SPR
DPT CUR HS FP DP DIR SPR WND UST LM
$
19680606 000000 3600 19680608 000000
-100.E3 50.E3 'A1'
Expand Down Expand Up @@ -43,7 +43,7 @@ $
0.E3 0.E3 'STOPSTRING'
$
19680606 000000 0 19680608 000000
19680606 011200 0 19680606 011200
19680606 020000 3600 19680606 020000
19680606 000000 0 19680608 000000
19680606 000000 0 19680608 000000
$
Expand Down
3 changes: 2 additions & 1 deletion regtests/mww3_test_04/input/ww3_multi_grdset_d.nml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -164,7 +164,7 @@
! ITYPE(3)%TRACK%FORMAT = F
! -------------------------------------------------------------------- !
&OUTPUT_TYPE_NML
ALLTYPE%FIELD%LIST = 'DPT CUR HS FP DP DIR SPR'
ALLTYPE%FIELD%LIST = 'DPT CUR HS FP DP DIR SPR WND UST LM'
ALLTYPE%POINT%NAME = 'points'
ALLTYPE%POINT%FILE = '../input/points.list'
/
Expand Down Expand Up @@ -208,6 +208,7 @@
&OUTPUT_DATE_NML
ALLDATE%FIELD = '19680606 000000' '1200' '19680608 000000'
ALLDATE%POINT = '19680606 000000' '3600' '19680608 000000'
ALLDATE%RESTART = '19680606 020000' '3600' '19680606 020000'
/

! -------------------------------------------------------------------- !
Expand Down

0 comments on commit 6d95c33

Please sign in to comment.