Skip to content

Commit

Permalink
Fixed a bug in the new WWM: use connectivity arrays from SCHISM with
Browse files Browse the repository at this point in the history
extreme caution.
  • Loading branch information
josephzhang8 committed Nov 12, 2020
1 parent a5d90b8 commit f9043e2
Show file tree
Hide file tree
Showing 3 changed files with 26 additions and 23 deletions.
4 changes: 2 additions & 2 deletions cmake/SCHISM.local.build
Original file line number Diff line number Diff line change
Expand Up @@ -27,12 +27,12 @@ set( USE_ICE OFF CACHE BOOLEAN "Use ICE module")
set( USE_GEN OFF CACHE BOOLEAN "Use generic tracer module")
set( USE_AGE OFF CACHE BOOLEAN "Use age module")
set( USE_ECO OFF CACHE BOOLEAN "Use ECO-SIM module")
set( USE_ICM OFF CACHE BOOLEAN "Use ICM module")
set( USE_ICM ON CACHE BOOLEAN "Use ICM module")
##pH model is active only if ICM is on
set( ICM_PH OFF CACHE BOOLEAN "Use pH module inside ICM")
set( USE_COSINE OFF CACHE BOOLEAN "Use CoSiNE module")
set( USE_FIB OFF CACHE BOOLEAN "Use fecal indicating bacteria module")
set( USE_SED OFF CACHE BOOLEAN "Use sediment module")
set( USE_SED ON CACHE BOOLEAN "Use sediment module")
set( USE_FABM OFF CACHE BOOLEAN "FABM BGC model interface")
#If FABM is on, need to set FABM_BASE (after cloning from Joseph's fork: https://github.com/josephzhang8/fabm.git)
set( FABM_BASE /sciclone/home10/yinglong/git/fabm CACHE STRING "Path to FABM base")
Expand Down
44 changes: 23 additions & 21 deletions src/WWMIII/wwm_coupl_selfe.F90
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@ SUBROUTINE RADIATION_STRESS_SCHISM
REAL(rkind) :: WSTMP, DEPLOC

REAL(rkind) :: HTOT
REAL(rkind) :: DSXX3D(2,NVRT,MNS), DSYY3D(2,NVRT,MNS),DSXY3D(2,NVRT,MNS)
REAL(rkind) :: DSXX3D(2,NVRT,nsa), DSYY3D(2,NVRT,nsa),DSXY3D(2,NVRT,nsa)


!GD: imet_dry allows to choose between 2 different methods to compute the
Expand Down Expand Up @@ -164,17 +164,17 @@ SUBROUTINE RADIATION_STRESS_SCHISM

! Computing gradients of the depth-averaged radiation stress
! terms (unit: m^2/s/s)
CALL hgrad_nodes(IMET_DRY,0,NVRT,MNP,MNS,SXX3D,DSXX3D) ! (dSxx/dx , dSxx/dy )
CALL hgrad_nodes(IMET_DRY,0,NVRT,MNP,MNS,SYY3D,DSYY3D) ! (dSyy/dx , dSyy/dy )
CALL hgrad_nodes(IMET_DRY,0,NVRT,MNP,MNS,SXY3D,DSXY3D) ! (dSxy/dx , dSxy/dy )
CALL hgrad_nodes(IMET_DRY,0,NVRT,MNP,nsa,SXX3D,DSXX3D) ! (dSxx/dx , dSxx/dy )
CALL hgrad_nodes(IMET_DRY,0,NVRT,MNP,nsa,SYY3D,DSYY3D) ! (dSyy/dx , dSyy/dy )
CALL hgrad_nodes(IMET_DRY,0,NVRT,MNP,nsa,SXY3D,DSXY3D) ! (dSxy/dx , dSxy/dy )
CALL exchange_s3d_2(DSXX3D)
CALL exchange_s3d_2(DSYY3D)
CALL exchange_s3d_2(DSXY3D)

! Computing the wave forces; these are noted Rsx, Rsy in Rolland
! et al. (2012), see Eq. (9)
! These are stored in wwave_force(:,1:MNS,1:2) (unit: m/s/s)
DO IS = 1, MNS
! These are stored in wwave_force(:,1:nsa,1:2) (unit: m/s/s)
DO IS = 1, nsa
IF(idry_s(IS) == 1) CYCLE

! Total water depth at sides
Expand Down Expand Up @@ -202,7 +202,8 @@ SUBROUTINE STOKES_STRESS_INTEGRAL_SCHISM

USE DATAPOOL
USE schism_glbl, ONLY: errmsg, hmin_radstress, ns, kbs, kbe, nea, idry_e, &
& isdel, indel, elnode, dldxy, zs, area
& isdel, indel, elnode, dldxy, zs, area,nsa,idry_s, &
&isidenode
USE schism_msgp
IMPLICIT NONE

Expand All @@ -213,9 +214,9 @@ SUBROUTINE STOKES_STRESS_INTEGRAL_SCHISM

integer :: IE, isd, j, l, n1, n2, n3, icount
real(rkind) :: tmp0, tmp1, tmp2, ztmp, ubar, vbar, dhdx, dhdy
real(rkind) :: STOKES_WVEL_ELEM(NVRT,MNE), ws_tmp1(NVRT,MNS),ws_tmp2(NVRT,MNS)
real(rkind) :: STOKES_WVEL_ELEM(NVRT,MNE), ws_tmp1(NVRT,nsa),ws_tmp2(NVRT,nsa)

real(rkind) :: dr_dxy_loc(2,NVRT,MNS)
real(rkind) :: dr_dxy_loc(2,NVRT,nsa)


!... Computing Stokes drift horizontal velocities at nodes and
Expand Down Expand Up @@ -307,7 +308,7 @@ SUBROUTINE STOKES_STRESS_INTEGRAL_SCHISM
! The average of the values from vertically adjacent nodes is
! taken
STOKES_HVEL_SIDE = 0.D0; ROLLER_STOKES_HVEL_SIDE = 0.D0
DO IS = 1,MNS
DO IS = 1,nsa
IF(idry_s(IS) == 1) CYCLE

! Indexes of surrounding nodes
Expand All @@ -323,7 +324,7 @@ SUBROUTINE STOKES_STRESS_INTEGRAL_SCHISM
ROLLER_STOKES_HVEL_SIDE(2,k,IS) = (ROLLER_STOKES_HVEL(2,k,n1) + ROLLER_STOKES_HVEL(2,k,n2))/2.D0
END IF
END DO
END DO !MNS
END DO !nsa

!... Compute bottom Stokes drift z-vel. at elements
STOKES_WVEL_ELEM = 0.D0
Expand Down Expand Up @@ -368,12 +369,12 @@ SUBROUTINE STOKES_STRESS_INTEGRAL_SCHISM
!... Compute horizontal gradient of Stokes x and y-vel. (to compute
!Stokes z-vel.)
ws_tmp1 = 0.D0; ws_tmp2 = 0.D0
CALL hgrad_nodes(2,0,NVRT,MNP,MNS,STOKES_HVEL(1,:,:),dr_dxy_loc)
CALL hgrad_nodes(2,0,NVRT,MNP,nsa,STOKES_HVEL(1,:,:),dr_dxy_loc)
ws_tmp1(:,:) = dr_dxy_loc(1,:,:) !valid only in resident
CALL hgrad_nodes(2,0,NVRT,MNP,MNS,STOKES_HVEL(2,:,:),dr_dxy_loc)
CALL hgrad_nodes(2,0,NVRT,MNP,nsa,STOKES_HVEL(2,:,:),dr_dxy_loc)
ws_tmp2(:,:) = dr_dxy_loc(2,:,:)

!... Compute Stokes z-vel. at all levels: STOKES_WVEL_SIDE(NVRT,MNS)
!... Compute Stokes z-vel. at all levels: STOKES_WVEL_SIDE(NVRT,nsa)
STOKES_WVEL_SIDE = 0.D0
DO IS = 1,ns !residents only
IF(idry_s(IS) == 1) CYCLE
Expand Down Expand Up @@ -405,20 +406,20 @@ END SUBROUTINE STOKES_STRESS_INTEGRAL_SCHISM
!**********************************************************************
SUBROUTINE COMPUTE_CONSERVATIVE_VF_TERMS_SCHISM
USE DATAPOOL
USE schism_glbl, ONLY: kbs, ns, idry_e, isdel, elnode, dldxy, cori, zs, su2, sv2
USE schism_glbl, ONLY: kbs, ns, idry_e, isdel, elnode, dldxy, cori, zs, su2, sv2,nsa,idry_s
USE schism_msgp
IMPLICIT NONE

integer :: IS, IE, k, l, icount
real(rkind) :: dJ_dx_loc, dJ_dy_loc, du_loc, dv_loc, dz_loc, Ust_loc, Vst_loc
real(rkind) :: du_dxy(2,NVRT,MNS), dv_dxy(2,NVRT,MNS)
real(rkind) :: du_dxy(2,NVRT,nsa), dv_dxy(2,NVRT,nsa)

!... Initialisation
WWAVE_FORCE = 0.D0

!... Computing the spatial derivative of horizontal velocities
CALL hgrad_nodes(2,0,NVRT,MNP,MNS,uu2,du_dxy)
CALL hgrad_nodes(2,0,NVRT,MNP,MNS,vv2,dv_dxy)
CALL hgrad_nodes(2,0,NVRT,MNP,nsa,uu2,du_dxy)
CALL hgrad_nodes(2,0,NVRT,MNP,nsa,vv2,dv_dxy)

!... Main loop over the sides
DO IS = 1,ns !resident
Expand Down Expand Up @@ -515,7 +516,8 @@ END SUBROUTINE COMPUTE_CONSERVATIVE_VF_TERMS_SCHISM
!**********************************************************************
SUBROUTINE COMPUTE_BREAKING_VF_TERMS_SCHISM
USE DATAPOOL
USE schism_glbl, ONLY: hmin_radstress, kbs, ns, isbs, dps, h0, out_wwm, zs
USE schism_glbl, ONLY: hmin_radstress, kbs, ns, isbs, dps, h0, out_wwm, &
&zs,nsa,idry_s,isidenode
USE schism_msgp
IMPLICIT NONE

Expand Down Expand Up @@ -610,7 +612,7 @@ SUBROUTINE COMPUTE_BREAKING_VF_TERMS_SCHISM
! With this profile, F < 10% of computed F at h < DMIN, and F > 95% of computed F at h > 2.25*DMIN
IF (htot < 3*DMIN) WWAVE_FORCE(:,:,IS) = WWAVE_FORCE(:,:,IS)*tanh((0.5D0*htot/DMIN)**8.D0)

END DO !MNS
END DO !nsa

! Exchange between ghost regions
CALL exchange_s3d_2(WWAVE_FORCE)
Expand Down Expand Up @@ -693,7 +695,7 @@ END SUBROUTINE SHORELINE_WAVE_FORCES
SUBROUTINE APPLY_WAFO_OPBND_RAMP

USE DATAPOOL
USE schism_glbl, only : ns,kbs
USE schism_glbl, only : ns,kbs,idry_s,nsa
USE schism_msgp

IMPLICIT NONE
Expand Down
1 change: 1 addition & 0 deletions src/WWMIII/wwm_datapl.F90
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,7 @@ MODULE DATAPOOL
! CANNOT HANDLE QUADS!
use schism_glbl, only : MNE => nea_wwm, & ! Elements of the augmented domain
& MNP => npa, & ! Nodes in the augmented domain
& NP_RES => np, & ! Local number of resident nodes
& np, &
& npg, & ! number of ghost nodes
& MNEI => mnei_wwm, & ! Max number of neighboring elements surrounding a node, nodes is mnei+1!
Expand Down

0 comments on commit f9043e2

Please sign in to comment.