Changeset 648 for codes/icosagcm/trunk/src

Timestamp:

12/21/17 22:38:57 (6 years ago)

Author:

aslmd

Message:

added an option for simple Rayleigh friction at bottom (giant planets, but could be used for any case really)

File:

• codes/icosagcm/trunk/src/dissip/dissip_gcm.f90 (modified) (7 diffs)

codes/icosagcm/trunk/src/dissip/dissip_gcm.f90

@@ -35 +35 @@
 !$OMP THREADPRIVATE(rayleigh_friction_type)
 !$OMP THREADPRIVATE(rayleigh_shear)
-  REAL, SAVE    :: rayleigh_tau
+  REAL, SAVE    :: rayleigh_tau, rayleigh_limlat
 !$OMP THREADPRIVATE(rayleigh_tau)
+!$OMP THREADPRIVATE(rayleigh_limlat)
   
   REAL,SAVE    :: dtdissip
@@ -103 +104 @@
       rayleigh_friction_type=2
       IF (is_master) PRINT *, 'Rayleigh friction : Schaer-like mountain with shear DCMIP2.2'
+   CASE('giant_liu_schneider') 
+      rayleigh_friction_type=99
+      IF (is_master) PRINT *, 'Rayleigh friction : giant planets Liu Schneider 2010'
    CASE DEFAULT
       IF (is_master) PRINT *, 'Bad selector : rayleigh_friction_type =', TRIM(rayleigh_friction_key), &
@@ -117 +121 @@
          STOP
       END IF
+      IF(rayleigh_friction_type == 99) THEN
+         rayleigh_limlat=0.
+         CALL getin("rayleigh_limlat",rayleigh_limlat)
+         rayleigh_limlat = rayleigh_limlat*3.14159/180.
+      ENDIF
    END IF
 
@@ -479 +488 @@
 
     INTEGER :: ind, shear
-    INTEGER :: l,ij
+    INTEGER :: l,ij,nn
 
 !$OMP BARRIER
@@ -515 +524 @@
 
       IF(rayleigh_friction_type>0) THEN
+       IF(rayleigh_friction_type<99) THEN
          phi=f_geopot(ind)
          ue=f_ue(ind)
@@ -524 +534 @@
             ENDDO
          END DO
+       ELSE
+         ue=f_ue(ind)
+            DO ij=ij_begin,ij_end
+              nn = ij+u_right
+              IF (ABS(lat_e(nn)) .gt. rayleigh_limlat) THEN
+                !print*, "latitude", lat_e(nn)*180./3.14159
+                due(nn,ll_begin:ll_begin+1) = due(nn,ll_begin:ll_begin+1) - (ue(nn,ll_begin:ll_begin+1)/rayleigh_tau)
+              ENDIF
+              nn = ij+u_lup
+              IF (ABS(lat_e(nn)) .gt. rayleigh_limlat) THEN
+                due(nn,ll_begin:ll_begin+1) = due(nn,ll_begin:ll_begin+1) - (ue(nn,ll_begin:ll_begin+1)/rayleigh_tau)
+              ENDIF
+              nn = ij+u_ldown
+              IF (ABS(lat_e(nn)) .gt. rayleigh_limlat) THEN
+                due(nn,ll_begin:ll_begin+1) = due(nn,ll_begin:ll_begin+1) - (ue(nn,ll_begin:ll_begin+1)/rayleigh_tau)
+              ENDIF
+            ENDDO
+       ENDIF
       END IF
    END DO
@@ -532 +560 @@
 
     CONTAINS
+
       SUBROUTINE relax(shift_t, shift_u)
         USE dcmip_initial_conditions_test_1_2_3