extelpf_edge.f90

Go to the documentation of this file.

!/===========================================================================/
! Copyright (c) 2007, The University of Massachusetts Dartmouth 
! Produced at the School of Marine Science & Technology 
! Marine Ecosystem Dynamics Modeling group
! All rights reserved.
!
! FVCOM has been developed by the joint UMASSD-WHOI research team. For 
! details of authorship and attribution of credit please see the FVCOM
! technical manual or contact the MEDM group.
!
! 
! This file is part of FVCOM. For details, see http://fvcom.smast.umassd.edu 
! The full copyright notice is contained in the file COPYRIGHT located in the 
! root directory of the FVCOM code. This original header must be maintained
! in all distributed versions.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" 
! AND ANY EXPRESS OR  IMPLIED WARRANTIES, INCLUDING,  BUT NOT  LIMITED TO,
! THE IMPLIED WARRANTIES OF MERCHANTABILITY AND  FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED.  
!
!/---------------------------------------------------------------------------/
! CVS VERSION INFORMATION
! $Id$
! $Name$
! $Revision$
!/===========================================================================/

!==============================================================================|
!  CALCULATE n+1 time step perturbation water surface elevation                |
!==============================================================================|
   SUBROUTINE extelpf_edge(KKT)       
!==============================================================================|
   USE all_vars
   USE bcs

   USE mod_obcs2
   USE mod_obcs3

   IMPLICIT NONE
   INTEGER, INTENT(IN) :: KKT
   REAL(SP) :: XFLUX(0:MT)
   REAL(SP) :: DIJ,UIJ,VIJ,DTK,EXFLUX,FXFLUX,FFLUX(1:IOBCN),TP,CC,CP
   INTEGER  :: I,J,I1,IA,IB,JJ,J1,J2,J3,II,JN

!---------ACCUMULATE FLUX BY LOOPING OVER CONTROL VOLUME HALF EDGES------------!

   xflux = 0.0_sp
   DO ii=1,nobcv
     i = nobedge_lst(ii)
     i1  = ntrg(i)
     ia  = niec(i,1)
     ib  = niec(i,2)
     j1  = i_obc_node(nv(i1,1))
     j2  = i_obc_node(nv(i1,2))
     j3  = i_obc_node(nv(i1,3))
     dij = h1(i1) + one_third*(elt(j1)+elt(j2)+elt(j3))
     uij = uat(i_obc_cell(i1))
     vij = vat(i_obc_cell(i1))
     exflux = dij*(-uij*dltye(i) + vij*dltxe(i))  


     xflux(ia) = xflux(ia)-exflux
     xflux(ib) = xflux(ib)+exflux
   END DO

   DO i = 1, iobcn
      j = i_obc_n(i)
      j1= i_obc_node(j)
      fxflux = -(eltf(j1)-elrkt(j1))*art1(j)/(alpha_rk(kkt)*dte) - xflux(j)
      fflux(i)= fxflux / (h(j)+elt(j1)) * elp(j1)
   END DO
       
   xflux = 0.0_sp
   DO ii=1,nobcv
     i = nobedge_lst(ii)
     i1  = ntrg(i)
     ia  = niec(i,1)
     ib  = niec(i,2)
     j1  = i_obc_node(nv(i1,1))
     j2  = i_obc_node(nv(i1,2))
     j3  = i_obc_node(nv(i1,3))
     dij = one_third * (elp(j1)+elp(j2)+elp(j3))
     uij = uat(i_obc_cell(i1))
     vij = vat(i_obc_cell(i1))
     exflux = dij*(-uij*dltye(i) + vij*dltxe(i))

  
     xflux(ia) = xflux(ia)-exflux
     xflux(ib) = xflux(ib)+exflux
   END DO

   DO i = 1, iobcn
      j = i_obc_n(i)
      fflux(i) = fflux(i) + xflux(j)
   END DO

   xflux = 0.0_sp
   DO ii=1,nobcv
     i = nobedge_lst(ii)
     i1  = ntrg(i)
     ia  = niec(i,1)
     ib  = niec(i,2)
     j1  = i_obc_node(nv(i1,1))
     j2  = i_obc_node(nv(i1,2))
     j3  = i_obc_node(nv(i1,3))
     dij = h1(i1) + one_third*(elt(j1)+elt(j2)+elt(j3)+ &
           elp(j1)+ elp(j2)+elp(j3))
     uij = uap(i_obc_cell(i1))
     vij = vap(i_obc_cell(i1))
     exflux = dij*(-uij*dltye(i) + vij*dltxe(i))  


     xflux(ia) = xflux(ia)-exflux
     xflux(ib) = xflux(ib)+exflux
   END DO

   DO i = 1, iobcn
      j = i_obc_n(i)
      fflux(i) = fflux(i) + xflux(j) + fluxobn2(i)
   END DO
   fluxobn2 = 0.0_sp

   IF (ibcn(1) > 0) THEN
   DO i  = 1, ibcn(1)
      jn = obc_lst(1,i)
      j  = i_obc_n(jn)
      j1 = i_obc_node(j)
      elpf(j1) = elrkp(j1) - alpha_rk(kkt)*dte*fflux(jn)/art1(j)
   END DO
   END IF

   tp = 10800.0_sp       ! 3 hours
   IF (ibcn(4) > 0) THEN
   DO i  = 1, ibcn(4)
      jn = obc_lst(4,i)
      j  = i_obc_n(jn)
      j1 = i_obc_node(j)

      i1 = i_obc_node(next_obc(jn))
      cc = sqrt(grav_n(j)*d(j))*alpha_rk(kkt)*dte/dltn_obc(jn)
      cp = cc + 1.0_sp

      elpf(j1) = (cc*elpf(i1) + elrkp(j1)*(1.0_sp-alpha_rk(kkt)*dte/tp))/cp
   END DO
   END IF
 
   IF((ibcn(2)>0).or.(ibcn(3)>0).or.(ibcn(5)>0)) print *,"error in EXTELPF_EDGE.F"

   RETURN
   END SUBROUTINE extelpf_edge
!==============================================================================|