Changeset 2558 for LMDZ5/trunk/libf

Timestamp:

Jun 10, 2016, 1:03:29 PM (8 years ago)

Author:

jyg

Message:

Introduction of the Effective Lifting Power (ELP),
which takes into account the area (cv_feed_area)
covered by thermals contributing to each
cumulonimbus.

The use of ELP prevents singularities when the

trigger probability tends to zero. It is
activated by iflag_clos_bl = 3.

The ELP values are stored in the ALP_bl variable.

File:

• LMDZ5/trunk/libf/phylmd/alpale_th.F90 (modified) (4 diffs)

LMDZ5/trunk/libf/phylmd/alpale_th.F90

@@ -1 @@
-SUBROUTINE alpale_th ( dtime, lmax_th, t_seri,  &
+SUBROUTINE alpale_th ( dtime, lmax_th, t_seri, cell_area,  &
                        cin, s2, n2,  &
                        ale_bl_trig, ale_bl_stat, ale_bl,  &
@@ -28 +28 @@
 !----------------
   REAL, INTENT(IN)                                           :: dtime
+  REAL, DIMENSION(klon), INTENT(IN)                          :: cell_area
   INTEGER, DIMENSION(klon), INTENT(IN)                       :: lmax_th
   REAL, DIMENSION(klon,klev), INTENT(IN)                     :: t_seri
@@ -47 +48 @@
   LOGICAL, SAVE                                              :: first = .TRUE.
   REAL, SAVE                                                 :: random_notrig_max=1.
+  REAL, SAVE                                                 :: cv_feed_area
+  REAL                                                       :: birth_number
+  REAL, DIMENSION(klon)                                      :: ale_bl_ref
   REAL, DIMENSION(klon)                                      :: tau_trig
+  REAL, DIMENSION(klon)                                      :: birth_rate
 !
     !$OMP THREADPRIVATE(random_notrig_max) 
+    !$OMP THREADPRIVATE(cv_feed_area) 
     !$OMP THREADPRIVATE(first)
 !
+ REAL umexp  ! expression of (1.-exp(-x))/x valid for all x, especially when x->0
+ REAL x
+ umexp(x) = max(sign(1.,x-1.e-3),0.)*(1.-exp(-x))/max(x,1.e-3) + &
+            (1.-max(sign(1.,x-1.e-3),0.))*(1.-0.5*x*(1.-x/3.*(1.-0.25*x)))
+!
+!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+!  JYG, 20160513 : Introduction of the Effective Lifting Power (ELP), which 
+! takes into account the area (cv_feed_area) covered by thermals contributing 
+! to each cumulonimbus.
+!   The use of ELP prevents singularities when the trigger probability tends to
+! zero. It is activated by iflag_clos_bl = 3.
+!   The ELP values are stored in the ALP_bl variable.
+!   
+!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+!
+!---------------------------------------
+  IF (iflag_clos_bl .LT. 3) THEN
+!---------------------------------------
+!
+!      Original code (Nicolas Rochetin)
+!     --------------------------------
+
     IF (first) THEN
+       random_notrig_max=1.
        CALL getin_p('random_notrig_max',random_notrig_max)
        first=.FALSE.
@@ -158 +187 @@
           endif !(iflag_clos_bl)
 
+!
+!---------------------------------------
+  ELSEIF (iflag_clos_bl .EQ. 3) THEN  ! (iflag_clos_bl .LT. 3)
+!---------------------------------------
+!
+!      New code with Effective Lifting Power
+!     -------------------------------------
+    IF (first) THEN
+       cv_feed_area = 1.e10   ! m2
+       CALL getin_p('cv_feed_area', cv_feed_area)
+       first=.FALSE.
+    ENDIF
+
+          !-----------Stochastic triggering-----------
+     if (iflag_trig_bl.ge.1) then
+        !
+        IF (prt_level .GE. 10) THEN
+           print *,'cin, ale_bl_stat, alp_bl_stat ', &
+                cin, ale_bl_stat, alp_bl_stat
+        ENDIF
+
+        ! Use ale_bl_stat (Rochetin's code) or ale_bl (old code) according to 
+        ! iflag_trig_bl value.
+        IF (iflag_trig_bl.eq.1) then         ! use ale_bl_stat (Rochetin computation)
+         do i=1,klon
+              ale_bl_ref(i)=ale_bl_stat(i)
+         enddo
+        ELSE IF (iflag_trig_bl.ge.2) then    ! use ale_bl (old computation)
+         do i=1,klon
+              ale_bl_ref(i)=Ale_bl(i)
+         enddo
+        ENDIF ! (iflag_trig_bl.eq.1)
+
+
+        !----Initializations and random number generation
+        do i=1,klon
+           proba_notrig(i)=1.
+           random_notrig(i)=1e6*ale_bl_stat(i)-int(1e6*ale_bl_stat(i))
+           if ( ale_bl_trig(i) .lt. abs(cin(i))+1.e-10 ) then 
+              tau_trig(i)=tau_trig_shallow
+           else
+              tau_trig(i)=tau_trig_deep
+           endif
+        enddo
+        !
+        IF (prt_level .GE. 10) THEN
+           print *,'random_notrig, tau_trig ', &
+                random_notrig, tau_trig
+           print *,'s_trig,s2,n2 ', &
+                s_trig,s2,n2
+        ENDIF
+
+        !----alp_bl computation
+        do i=1,klon
+           if ( (ale_bl_ref(i) .gt. abs(cin(i))+1.e-10) )  then
+              birth_number = n2(i)*exp(-s_trig/s2(i))
+              birth_rate(i) = birth_number/(tau_trig(i)*cell_area(i))
+              proba_notrig(i)=exp(-birth_number*dtime/tau_trig(i))
+              Alp_bl(i) = Alp_bl(i)* &
+                          umexp(-birth_number*cv_feed_area/cell_area(i))/ &
+                          umexp(-birth_number*dtime/tau_trig(i))*  &
+                          tau_trig(i)*cv_feed_area/(dtime*cell_area(i))
+          else 
+              proba_notrig(i)=1.
+              random_notrig(i)=0.
+              alp_bl(i)=0.
+           endif
+        enddo
+
+        !----ale_bl_trig computation
+         do i=1,klon
+           if (random_notrig(i) .ge. proba_notrig(i)) then 
+              ale_bl_trig(i)=ale_bl_ref(i)
+           else
+              ale_bl_trig(i)=0.
+           endif
+         enddo
+
+        !
+        IF (prt_level .GE. 10) THEN
+           print *,'proba_notrig, ale_bl_trig ', &
+                proba_notrig, ale_bl_trig
+        ENDIF
+
+     endif !(iflag_trig_bl .ge. 1)
+
+!---------------------------------------
+  ENDIF ! (iflag_clos_bl .LT. 3)
+!---------------------------------------
+
           IF (prt_level .GE. 10) THEN
              print *,'ale_bl_trig, alp_bl_stat ',ale_bl_trig, alp_bl_stat