Changeset 4292 for branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn.F90

Timestamp:

2013-11-20T17:28:04+01:00 (11 years ago)

Author:

cetlod

Message:

dev_MERGE_2013 : 1st step of the merge, see ticket #1185

File:

• branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn.F90 (modified) (6 diffs)

branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn.F90

@@ -30 +30 @@
    USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
    USE in_out_manager  !
-   USE domvvl          ! variable volume
+   USE domvvl
 
    IMPLICIT NONE
@@ -57 +57 @@
       LOGICAL, INTENT( in ), OPTIONAL :: dyn3d_only       ! T => only update baroclinic velocities
       !!
-      INTEGER               :: jk,ii,ij,ib,igrd     ! Loop counter
-      LOGICAL               :: ll_dyn2d, ll_dyn3d  
-      !!
+      INTEGER               :: jk,ii,ij,ib_bdy,ib,igrd     ! Loop counter
+      LOGICAL               :: ll_dyn2d, ll_dyn3d, ll_orlanski
+      !!
+      REAL(wp), POINTER, DIMENSION(:,:) :: phur1, phvr1     ! inverse depth at u and v points
 
       IF( nn_timing == 1 ) CALL timing_start('bdy_dyn')
@@ -70 +71 @@
       ENDIF
 
+      ll_orlanski = .false.
+      DO ib_bdy = 1, nb_bdy
+         IF ( cn_dyn2d(ib_bdy) == 'orlanski' .or. cn_dyn2d(ib_bdy) == 'orlanski_npo' &
+     &   .or. cn_dyn3d(ib_bdy) == 'orlanski' .or. cn_dyn3d(ib_bdy) == 'orlanski_npo') ll_orlanski = .true.
+      ENDDO
+
       !-------------------------------------------------------
       ! Set pointers
@@ -77 +84 @@
       phur => hur
       phvr => hvr
-      CALL wrk_alloc(jpi,jpj,pu2d,pv2d) 
+      CALL wrk_alloc(jpi,jpj,pua2d,pva2d) 
+      IF ( ll_orlanski ) CALL wrk_alloc(jpi,jpj,pub2d,pvb2d,phur1,phvr1) 
 
       !-------------------------------------------------------
@@ -83 +91 @@
       !-------------------------------------------------------
 
-      pu2d(:,:) = 0.e0
-      pv2d(:,:) = 0.e0
+      ! "After" velocities: 
+
+      pua2d(:,:) = 0.e0
+      pva2d(:,:) = 0.e0
+      
       IF (lk_vvl) THEN
-         DO jk = 1, jpkm1   !! Vertically integrated momentum trends
-            pu2d(:,:) = pu2d(:,:) + fse3u_a(:,:,jk) * umask(:,:,jk) * ua(:,:,jk)
-            pv2d(:,:) = pv2d(:,:) + fse3v_a(:,:,jk) * vmask(:,:,jk) * va(:,:,jk)
-         END DO
-         pu2d(:,:) = pu2d(:,:) / ( hu_0(:,:) + sshu_a(:,:) + 1._wp - umask(:,:,1) )
-         pv2d(:,:) = pv2d(:,:) / ( hv_0(:,:) + sshv_a(:,:) + 1._wp - vmask(:,:,1) )
+         phur1(:,:) = 0.
+         phvr1(:,:) = 0.
+         DO jk = 1, jpkm1
+            phur1(:,:) = phur1(:,:) + fse3u_a(:,:,jk) * umask(:,:,jk)
+            phvr1(:,:) = phvr1(:,:) + fse3v_a(:,:,jk) * vmask(:,:,jk)
+            pua2d(:,:) = pua2d(:,:) + fse3u_a(:,:,jk) * umask(:,:,jk) * ua(:,:,jk)
+            pva2d(:,:) = pva2d(:,:) + fse3v_a(:,:,jk) * vmask(:,:,jk) * va(:,:,jk)
+         END DO
+         phur1(:,:) = umask(:,:,1) / ( phur1(:,:) + 1. - umask(:,:,1) )
+         phvr1(:,:) = vmask(:,:,1) / ( phvr1(:,:) + 1. - vmask(:,:,1) )
+         pua2d(:,:) = pua2d(:,:) * phur1(:,:)
+         pva2d(:,:) = pva2d(:,:) * phvr1(:,:)
       ELSE
-         DO jk = 1, jpkm1   !! Vertically integrated momentum trends
-            pu2d(:,:) = pu2d(:,:) + fse3u(:,:,jk) * umask(:,:,jk) * ua(:,:,jk)
-            pv2d(:,:) = pv2d(:,:) + fse3v(:,:,jk) * vmask(:,:,jk) * va(:,:,jk)
-         END DO
-         pu2d(:,:) = pu2d(:,:) * phur(:,:)
-         pv2d(:,:) = pv2d(:,:) * phvr(:,:)
+         DO jk = 1, jpkm1
+            pua2d(:,:) = pua2d(:,:) + fse3u(:,:,jk) * umask(:,:,jk) * ua(:,:,jk)
+            pva2d(:,:) = pva2d(:,:) + fse3v(:,:,jk) * vmask(:,:,jk) * va(:,:,jk)
+         END DO
+         pua2d(:,:) = pua2d(:,:) * phur(:,:)
+         pva2d(:,:) = pva2d(:,:) * phvr(:,:)
       ENDIF
+
       DO jk = 1 , jpkm1
-         ua(:,:,jk) = ua(:,:,jk) - pu2d(:,:) * umask(:,:,jk)
-         va(:,:,jk) = va(:,:,jk) - pv2d(:,:) * vmask(:,:,jk)
+         ua(:,:,jk) = ua(:,:,jk) - pua2d(:,:)
+         va(:,:,jk) = va(:,:,jk) - pva2d(:,:)
       END DO
+
+      ! "Before" velocities (required for Orlanski condition): 
+
+      IF ( ll_orlanski ) THEN          
+         pub2d(:,:) = 0.e0
+         pvb2d(:,:) = 0.e0
+
+         IF (lk_vvl) THEN
+            phur1(:,:) = 0.
+            phvr1(:,:) = 0.
+            DO jk = 1, jpkm1   !! Vertically integrated momentum trends
+               phur1(:,:) = phur1(:,:) + fse3u_b(:,:,jk) * umask(:,:,jk)
+               phvr1(:,:) = phvr1(:,:) + fse3v_b(:,:,jk) * vmask(:,:,jk)
+               pub2d(:,:) = pub2d(:,:) + fse3u_b(:,:,jk) * umask(:,:,jk) * ub(:,:,jk)
+               pvb2d(:,:) = pvb2d(:,:) + fse3v_b(:,:,jk) * vmask(:,:,jk) * vb(:,:,jk)
+            END DO
+            phur1(:,:) = umask(:,:,1) / ( phur1(:,:) + 1. - umask(:,:,1) )
+            phvr1(:,:) = vmask(:,:,1) / ( phvr1(:,:) + 1. - vmask(:,:,1) )
+            pub2d(:,:) = pub2d(:,:) * phur1(:,:)
+            pvb2d(:,:) = pvb2d(:,:) * phvr1(:,:)
+         ELSE
+            DO jk = 1, jpkm1   !! Vertically integrated momentum trends
+               pub2d(:,:) = pub2d(:,:) + fse3u(:,:,jk) * umask(:,:,jk) * ub(:,:,jk)
+               pvb2d(:,:) = pvb2d(:,:) + fse3v(:,:,jk) * vmask(:,:,jk) * vb(:,:,jk)
+            END DO
+            pub2d(:,:) = pub2d(:,:) * phur(:,:)
+            pvb2d(:,:) = pvb2d(:,:) * phvr(:,:)
+         ENDIF
+
+         DO jk = 1 , jpkm1
+            ub(:,:,jk) = ub(:,:,jk) - pub2d(:,:)
+            vb(:,:,jk) = vb(:,:,jk) - pvb2d(:,:)
+         END DO
+      END IF
 
       !-------------------------------------------------------
@@ -119 +171 @@
 
       DO jk = 1 , jpkm1
-         ua(:,:,jk) = ( ua(:,:,jk) + pu2d(:,:) ) * umask(:,:,jk)
-         va(:,:,jk) = ( va(:,:,jk) + pv2d(:,:) ) * vmask(:,:,jk)
+         ua(:,:,jk) = ( ua(:,:,jk) + pua2d(:,:) ) * umask(:,:,jk)
+         va(:,:,jk) = ( va(:,:,jk) + pva2d(:,:) ) * vmask(:,:,jk)
       END DO
 
-      CALL wrk_dealloc(jpi,jpj,pu2d,pv2d) 
+      IF ( ll_orlanski ) THEN
+         DO jk = 1 , jpkm1
+            ub(:,:,jk) = ( ub(:,:,jk) + pub2d(:,:) ) * umask(:,:,jk)
+            vb(:,:,jk) = ( vb(:,:,jk) + pvb2d(:,:) ) * vmask(:,:,jk)
+         END DO
+      END IF
+
+      CALL wrk_dealloc(jpi,jpj,pua2d,pva2d) 
+      IF ( ll_orlanski ) CALL wrk_dealloc(jpi,jpj,pub2d,pvb2d) 
 
       IF( nn_timing == 1 ) CALL timing_stop('bdy_dyn')