FSCHAPMAN and SSH_TIDES, ADD_FSOBC

General scientific issues regarding ROMS

Moderators: arango, robertson

Post Reply
Message
Author
longmtm
Posts: 55
Joined: Tue Dec 13, 2005 7:23 pm
Location: Univ of Maryland Center for Environmental Science

FSCHAPMAN and SSH_TIDES, ADD_FSOBC

#1 Unread post by longmtm »

Dear all,

Please help, I'm a bit confused here. I wonder how open boundary tidal water level is picked up by the model when EAST_FSCHAPMAN condition is used.

Code: https://www.myroms.org/svn/src/trunk/RO ... r/zetabc.F

In the following section of Nonlinear/zetabc.F

Code: Select all


# elif defined EAST_FSCHAPMAN
!
!  Eastern edge, Chapman boundary condition.
!
        DO j=Jstr,Jend
          cff=dt2d*GRID(ng)%pm(Iend,j)
          cff1=SQRT(g*(GRID(ng)%h(Iend,j)+                              &
     &                 zeta(Iend,j,know)))
          Cx=cff*cff1
          cff2=1.0_r8/(1.0_r8+Cx)
          zeta(Iend+1,j,kout)=cff2*(zeta(Iend+1,j,know)+                &
     &                              Cx*zeta(Iend,j,kout))
#  ifdef MASKING
          zeta(Iend+1,j,kout)=zeta(Iend+1,j,kout)*                      &
     &                        GRID(ng)%rmask(Iend+1,j)
#  endif
        END DO

In the upper section, I can NOT see how "BOUNDARY % zeta_east" is used although in Nonlinear/set_tides.F https://www.myroms.org/svn/src/trunk/RO ... et_tides.F , BOUNDARY % zeta_east is calculated using tidal forcing:

Code: Select all

#  if defined EAST_FSOBC || defined EAST_M2OBC
      update=.FALSE.
      IF (EASTERN_EDGE) THEN
        DO j=JstrR,JendR
#   ifdef ADD_FSOBC
          zeta_east(j)=zeta_east(j)+                                    &
     &                 0.5_r8*(Etide(Iend,j)+Etide(Iend+1,j))
#   else
          zeta_east(j)=0.5_r8*(Etide(Iend,j)+Etide(Iend+1,j))
#   endif

Anything special or weird here?

Thanks!

Wen

User avatar
arango
Site Admin
Posts: 1361
Joined: Wed Feb 26, 2003 4:41 pm
Location: DMCS, Rutgers University
Contact:

Re: FSCHAPMAN and SSH_TIDES, ADD_FSOBC

#2 Unread post by arango »

Well, this have been discussed many times in this forum. Please use the search capabilities of the forum to answer your question(s). I will advise you also to read the Chapman and Flather papers.

If you want to force ROMS with tides, you need to force the prognostic, vertically integrated equations for ubar and vbar and not the vertically integrated continuity equation for zeta. The zeta equation in ROMS is forced by the 2D-momentum divergence, div(ubar, vbar). That's how the tidal forcing enter into ROMS and not as a boundary conditions in zeta :!:

Therefore, you need to look at the 1D-radiation Flather boundary conditions in u2dbc_im.F and v3dbc_im.F. Check those routines instead and review the Flather condition.

There is nothing weird here but misunderstanding.

longmtm
Posts: 55
Joined: Tue Dec 13, 2005 7:23 pm
Location: Univ of Maryland Center for Environmental Science

Re: FSCHAPMAN and SSH_TIDES, ADD_FSOBC

#3 Unread post by longmtm »

Thanks.

The u2dbc_im.F file https://www.myroms.org/svn/src/trunk/RO ... u2dbc_im.F treated EAST_M2FLATHER so that BOUNDARY(ng)%zeta_east is used to specify ubar.

Code: Select all


# elif defined EAST_M2FLATHER
!
!  Eastern edge, Flather boundary condition.
!
        DO j=Jstr,Jend
#  if defined SSH_TIDES && !defined UV_TIDES
#   ifdef FSOBC_REDUCED
          bry_pgr=-g*(BOUNDARY(ng)%zeta_east(j)-                        &
     &                zeta(Iend,j,know))*                               &
     &            0.5_r8*GRID(ng)%pm(Iend,j)
#   else
          bry_pgr=-g*(zeta(Iend+1,j,know)-                              &
     &                zeta(Iend  ,j,know))*                             &
     &            0.5_r8*(GRID(ng)%pm(Iend  ,j)+                        &
     &                    GRID(ng)%pm(Iend+1,j))
#   endif
#   ifdef UV_COR
          bry_cor=0.125_r8*(vbar(Iend  ,j  ,know)+                      &
     &                      vbar(Iend  ,j+1,know)+                      &
     &                      vbar(Iend+1,j  ,know)+                      &
     &                      vbar(Iend+1,j+1,know))*                     &
     &                     (GRID(ng)%f(Iend  ,j)+                       &
     &                      GRID(ng)%f(Iend+1,j))
#   else
          bry_cor=0.0_r8
#   endif
          cff1=1.0_r8/(0.5_r8*(GRID(ng)%h(Iend  ,j)+                    &
     &                         zeta(Iend  ,j,know)+                     &
     &                         GRID(ng)%h(Iend+1,j)+                    &
     &                         zeta(Iend+1,j,know)))
          bry_str=cff1*(FORCES(ng)%sustr(Iend+1,j)-                     &
     &                  FORCES(ng)%bustr(Iend+1,j))
          Cx=1.0_r8/SQRT(g*0.5_r8*(GRID(ng)%h(Iend+1,j)+                &
     &                             zeta(Iend+1,j,know)+                 &
     &                             GRID(ng)%h(Iend  ,j)+                &
     &                             zeta(Iend  ,j,know)))
          cff2=GRID(ng)%om_u(Iend+1,j)*Cx
!!        cff2=dt2d
          bry_val=ubar(Iend,j,know)+                                    &
     &            cff2*(bry_pgr+                                        &
     &                  bry_cor+                                        &
     &                  bry_str)
#  else
          bry_val=BOUNDARY(ng)%ubar_east(j)
#  endif
          cff=1.0_r8/(0.5_r8*(GRID(ng)%h(Iend  ,j)+                     &
     &                        zeta(Iend  ,j,know)+                      &
     &                        GRID(ng)%h(Iend+1,j)+                     &
     &                        zeta(Iend+1,j,know)))
          Cx=SQRT(g*cff)
          ubar(Iend+1,j,kout)=bry_val+                                  &
     &                        Cx*(0.5_r8*(zeta(Iend  ,j,know)+          &
     &                                    zeta(Iend+1,j,know))-         &
     &                            BOUNDARY(ng)%zeta_east(j))
#  ifdef MASKING
          ubar(Iend+1,j,kout)=ubar(Iend+1,j,kout)*                      &
     &                        GRID(ng)%umask(Iend+1,j)
#  endif
        END DO

Now I see how boundary forcing values of both ubar_east and zeta_east are used:

Code: Select all


   #  else
          bry_val=BOUNDARY(ng)%ubar_east(j)
#  endif
          cff=1.0_r8/(0.5_r8*(GRID(ng)%h(Iend  ,j)+                     &
     &                        zeta(Iend  ,j,know)+                      &
     &                        GRID(ng)%h(Iend+1,j)+                     &
     &                        zeta(Iend+1,j,know)))
          Cx=SQRT(g*cff)
          ubar(Iend+1,j,kout)=bry_val+                                  &
     &                        Cx*(0.5_r8*(zeta(Iend  ,j,know)+          &
     &                                    zeta(Iend+1,j,know))-         &
     &                            BOUNDARY(ng)%zeta_east(j))

The caveat is that if I use EAST_M2RADIATION+EAST_M2NUDGING for (ubar,vbar) and EAST_FSCHAPMAN for zeta
then the boundary forcing value for zeta, BOUNDARY(ng)%zeta_east, will NOT be used because the EAST_M2RADIATION+EAST_M2NUDGING section (below) in u2dbc_im.F does not involve BOUNDARY(ng)%zeta_east, neither does the EAST_FSCHAPMAN section of zetabc.F https://www.myroms.org/svn/src/trunk/RO ... r/zetabc.F

Code: Select all

# if defined EAST_M2RADIATION
!
!  Eastern edge, implicit upstream radiation condition.
!
        DO j=Jstr,Jend+1
          grad(Iend  ,j)=ubar(Iend  ,j  ,know)-                         &
     &                   ubar(Iend  ,j-1,know)
          grad(Iend+1,j)=ubar(Iend+1,j  ,know)-                         &
     &                   ubar(Iend+1,j-1,know)
        END DO
        DO j=Jstr,Jend
          dUdt=ubar(Iend,j,know)-ubar(Iend  ,j,kout)
          dUdx=ubar(Iend,j,kout)-ubar(Iend-1,j,kout)
#  ifdef EAST_M2NUDGING
          IF ((dUdt*dUdx).lt.0.0_r8) THEN
            tau=M2obc_in(ng,ieast)
          ELSE
            tau=M2obc_out(ng,ieast)
          END IF
          tau=tau*dt2d
#  endif
          IF ((dUdt*dUdx).lt.0.0_r8) dUdt=0.0_r8
          IF ((dUdt*(grad(Iend,j)+grad(Iend,j+1))).gt.0.0_r8) THEN
            dUde=grad(Iend,j)
          ELSE
            dUde=grad(Iend,j+1)
          END IF
          cff=MAX(dUdx*dUdx+dUde*dUde,eps)
          Cx=dUdt*dUdx
#  ifdef RADIATION_2D
          Ce=MIN(cff,MAX(dUdt*dUde,-cff))
#  else
          Ce=0.0_r8
#  endif
#  if defined CELERITY_WRITE && defined FORWARD_WRITE
          BOUNDARY(ng)%ubar_east_Cx(j)=Cx
          BOUNDARY(ng)%ubar_east_Ce(j)=Ce
          BOUNDARY(ng)%ubar_east_C2(j)=cff
#  endif
          ubar(Iend+1,j,kout)=(cff*ubar(Iend+1,j,know)+                 &
     &                         Cx *ubar(Iend  ,j,kout)-                 &
     &                         MAX(Ce,0.0_r8)*grad(Iend+1,j  )-         &
     &                         MIN(Ce,0.0_r8)*grad(Iend+1,j+1))/        &
     &                        (cff+Cx)
#  ifdef EAST_M2NUDGING
          ubar(Iend+1,j,kout)=ubar(Iend+1,j,kout)+                      &
     &                        tau*(BOUNDARY(ng)%ubar_east(j)-           &
     &                             ubar(Iend+1,j,know))
#  endif
#  ifdef MASKING
          ubar(Iend+1,j,kout)=ubar(Iend+1,j,kout)*                      &
     &                        GRID(ng)%umask(Iend+1,j)
#  endif
        END DO

Post Reply