Changeset 3414 for branches/2012/dev_r3385_NOCS04_HAMF/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90

Timestamp:

2012-06-15T14:53:43+02:00 (12 years ago)

Author:

acc

Message:

Branch: dev_r3385_NOCS04_HAMF; #665. Stage 3 of 2012 development: First working version with (optionally) active embedding. Still some stability questions when running the LIM2-EVP, ORCA2 test-case with full embedding. Works with reduced (halved) timestep. Still investigating.

File:

• branches/2012/dev_r3385_NOCS04_HAMF/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90 (modified) (7 diffs)

branches/2012/dev_r3385_NOCS04_HAMF/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90

@@ -35 +35 @@
    USE dom_ice_2        ! LIM2: ice domain
 #endif
+   USE oce, ONLY :  snwice_mass, snwice_mass_b
 
    IMPLICIT NONE
@@ -124 +125 @@
       REAL(wp) ::   zindb         ! ice (1) or not (0)      
       REAL(wp) ::   zdummy        ! dummy argument
+      REAL(wp) ::   zintb, zintn  ! dummy argument
 
       REAL(wp), POINTER, DIMENSION(:,:) ::   zpresh           ! temporary array for ice strength
@@ -144 +146 @@
       REAL(wp), POINTER, DIMENSION(:,:) ::   zs12             ! Non-diagonal stress tensor component zs12
       REAL(wp), POINTER, DIMENSION(:,:) ::   zu_ice, zv_ice, zresr   ! Local error on velocity
+      REAL(wp), POINTER, DIMENSION(:,:) ::   zpice            ! array used for the calculation of ice surface slope:
+                                                              !   ocean surface (ssh_m) if ice is not embedded
+                                                              !   ice top surface if ice is embedded
       
       !!-------------------------------------------------------------------
@@ -150 +155 @@
       CALL wrk_alloc( jpi,jpj, zc1   , u_oce1, u_oce2, u_ice2, zusw  , v_oce1 , v_oce2, v_ice1                )
       CALL wrk_alloc( jpi,jpj, zf1   , deltat, zu_ice, zf2   , deltac, zv_ice , zdd   , zdt    , zds          )
-      CALL wrk_alloc( jpi,jpj, zdd   , zdt   , zds   , zs1   , zs2   , zs12   , zresr                         )
+      CALL wrk_alloc( jpi,jpj, zdd   , zdt   , zds   , zs1   , zs2   , zs12   , zresr , zpice                 )
 
 #if  defined key_lim2 && ! defined key_lim2_vp
@@ -231 +236 @@
       !  v_oce2: ocean v component on v points                        
 
+      IF( nn_ice_embd == 2 ) THEN             !== embedded sea ice: compute representative ice top surface ==!
+          !                                            
+          ! average interpolation coeff as used in dynspg = (1/nn_fsbc) * {SUM[n/nn_fsbc], n=0,nn_fsbc-1}
+          !                                               = (1/nn_fsbc)^2 * {SUM[n], n=0,nn_fsbc-1}
+         zintn = REAL( nn_fsbc - 1 ) / REAL( nn_fsbc ) * 0.5_wp     
+          !
+          ! average interpolation coeff as used in dynspg = (1/nn_fsbc) * {SUM[1-n/nn_fsbc], n=0,nn_fsbc-1}
+          !                                               = (1/nn_fsbc)^2 * (nn_fsbc^2 - {SUM[n], n=0,nn_fsbc-1})
+         zintb = REAL( nn_fsbc + 1 ) / REAL( nn_fsbc ) * 0.5_wp
+          !
+         zpice(:,:) = ssh_m(:,:) + (  zintn * snwice_mass(:,:) +  zintb * snwice_mass_b(:,:)  ) * r1_rau0
+          !
+      ELSE                                    !== non-embedded sea ice: use ocean surface for slope calculation ==!
+         zpice(:,:) = ssh_m(:,:)
+      ENDIF
+
       DO jj = k_j1+1, k_jpj-1
          DO ji = fs_2, fs_jpim1
@@ -273 +294 @@
             ! include it later
 
-            zdsshx =  ( ssh_m(ji+1,jj) - ssh_m(ji,jj) ) / e1u(ji,jj)
-            zdsshy =  ( ssh_m(ji,jj+1) - ssh_m(ji,jj) ) / e2v(ji,jj)
+            zdsshx =  ( zpice(ji+1,jj) - zpice(ji,jj) ) / e1u(ji,jj)
+            zdsshy =  ( zpice(ji,jj+1) - zpice(ji,jj) ) / e2v(ji,jj)
 
             za1ct(ji,jj) = ztagnx - zmass1(ji,jj) * grav * zdsshx
@@ -746 +767 @@
       CALL wrk_dealloc( jpi,jpj, zc1   , u_oce1, u_oce2, u_ice2, zusw  , v_oce1 , v_oce2, v_ice1                )
       CALL wrk_dealloc( jpi,jpj, zf1   , deltat, zu_ice, zf2   , deltac, zv_ice , zdd   , zdt    , zds          )
-      CALL wrk_dealloc( jpi,jpj, zdd   , zdt   , zds   , zs1   , zs2   , zs12   , zresr                         )
+      CALL wrk_dealloc( jpi,jpj, zdd   , zdt   , zds   , zs1   , zs2   , zs12   , zresr , zpice                 )
 
    END SUBROUTINE lim_rhg