Changeset 5619 for branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90

Timestamp:

2015-07-20T19:43:15+02:00 (9 years ago)

Author:

mathiot

Message:

ocean/ice sheet coupling: initial commit

File:

• branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90 (modified) (9 diffs)

branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90

@@ -63 +63 @@
 
    REAL(wp), PUBLIC, SAVE ::   rcpi   = 2000.0_wp     ! phycst ?
-   REAL(wp), PUBLIC, SAVE ::   kappa  = 1.54e-6_wp    ! phycst ?
+   REAL(wp), PUBLIC, SAVE ::   kappa  =    0.0_wp    ! phycst ?
    REAL(wp), PUBLIC, SAVE ::   rhoisf = 920.0_wp      ! phycst ?
    REAL(wp), PUBLIC, SAVE ::   tsurf  = -20.0_wp      ! phycst ?
@@ -125 +125 @@
          IF ( lwp ) WRITE(numout,*) '        ln_divisf   = ', ln_divisf 
          IF ( lwp ) WRITE(numout,*) '        nn_gammablk = ', nn_gammablk 
+         IF ( lwp ) WRITE(numout,*) '        rn_gammat0  = ', rn_gammat0  
+         IF ( lwp ) WRITE(numout,*) '        rn_gammas0  = ', rn_gammas0  
          IF ( lwp ) WRITE(numout,*) '        rn_tfri2    = ', rn_tfri2 
          IF (ln_divisf) THEN       ! keep it in the namelist ??? used true anyway as for runoff ? (PM)
@@ -150 +152 @@
             !: read effective lenght (BG03)
             IF (nn_isf == 2) THEN
-               ! Read Data and save some integral values
+               cvarLeff = 'soLeff' 
                CALL iom_open( sn_Leff_isf%clname, inum )
-               cvarLeff  = 'soLeff'               !: variable name for Efficient Length scale
                CALL iom_get( inum, jpdom_data, cvarLeff, risfLeff , 1)
                CALL iom_close(inum)
@@ -237 +238 @@
             CALL sbc_isf_tbl(tsn(:,:,:,jp_tem),ttbl(:,:),'T')
             CALL sbc_isf_tbl(tsn(:,:,:,jp_sal),stbl(:,:),'T')
-            CALL sbc_isf_tbl(un(:,:,:),utbl(:,:),'U')
-            CALL sbc_isf_tbl(vn(:,:,:),vtbl(:,:),'V')
+            CALL sbc_isf_tbl(un(:,:,:)        ,utbl(:,:),'U')
+            CALL sbc_isf_tbl(vn(:,:,:)        ,vtbl(:,:),'V')
             ! iom print
             CALL iom_put('ttbl',ttbl(:,:))
             CALL iom_put('stbl',stbl(:,:))
-            CALL iom_put('utbl',utbl(:,:))
-            CALL iom_put('vtbl',vtbl(:,:))
+            CALL iom_put('utbl',utbl(:,:) * (1._wp - tmask(:,:,1)) * ssmask(:,:))
+            CALL iom_put('vtbl',vtbl(:,:) * (1._wp - tmask(:,:,1)) * ssmask(:,:))
             ! compute fwf and heat flux
             CALL sbc_isf_cav (kt)
@@ -296 +297 @@
          ! 
          ! output
-         CALL iom_put('qisf'  , qisf)
-         IF( iom_use('fwfisf') )   CALL iom_put('fwfisf', fwfisf * stbl(:,:) / soce )
+         IF( iom_use('qisf'  ) )   CALL iom_put('qisf'  , qisf)
+         IF( iom_use('fwfisf') )   CALL iom_put('fwfisf', fwfisf)
       END IF
   
@@ -416 +417 @@
       REAL(wp) ::   zfwflx, zhtflx, zhtflx_b
       REAL(wp) ::   zgammat, zgammas
-      REAL(wp) ::   zeps   =  -1.e-20_wp        !==   Local constant initialization   ==!
+      REAL(wp) ::   zeps = 1.e-20_wp        !==   Local constant initialization   ==!
       INTEGER  ::   ji, jj     ! dummy loop indices
       INTEGER  ::   ii0, ii1, ij0, ij1   ! temporary integers
@@ -506 +507 @@
                      zeps1=rcp*rau0*zgammat
                      zeps2=lfusisf*rau0*zgammas
-                     zeps3=rhoisf*rcpi*kappa/risfdep(ji,jj)
+                     zeps3=rhoisf*rcpi*kappa/MAX(risfdep(ji,jj),zeps)
                      zeps4=zlamb2+zlamb3*risfdep(ji,jj)
                      zeps6=zeps4-zti(ji,jj)
                      zeps7=zeps4-tsurf
                      zaqe=zlamb1 * (zeps1 + zeps3)
-                     zaqer=0.5/zaqe
+                     zaqer=0.5/MIN(zaqe,-zeps)
                      zbqe=zeps1*zeps6+zeps3*zeps7-zeps2
                      zcqe=zeps2*stbl(ji,jj)
@@ -520 +521 @@
                      zfrz(ji,jj)=zeps4+zlamb1*zsfrz
   
-! zfwflx is upward water flux
-                     zfwflx= rau0 * zgammas * ( (zsfrz-stbl(ji,jj)) / zsfrz )
+! zfwflx is upward water flux (MAX usefull if kappa = 0
+                     zfwflx= rau0 * zgammas * ( (zsfrz-stbl(ji,jj)) / MAX(zsfrz,zeps) )
 ! zhtflx is upward heat flux (out of ocean)
 ! If non conservative we have zcfac=0.0 so zhtflx is as ISOMIP but with different zfrz value
@@ -527 +528 @@
 ! zwflx is upward water flux
 ! If non conservative we have zcfac=0.0 so what follows is then zfwflx*sss_m/zsfrz
-                     zfwflx = ( zgammas*rau0 - zcfac*zfwflx ) * (zsfrz - stbl(ji,jj)) / stbl(ji,jj)
+!!!!!!!!                     zfwflx = ( zgammas*rau0 - zcfac*zfwflx ) * (zsfrz - stbl(ji,jj)) / stbl(ji,jj)
 ! test convergence and compute gammat
                      IF (( zhtflx - zhtflx_b) .LE. 0.01 ) lit = .FALSE.