Changeset 2799 for LMDZ5/trunk/libf/phylmd/climb_hq_mod.F90

Timestamp:

Feb 24, 2017, 7:50:33 PM (7 years ago)

Author:

jyg

Message:

New diagnostics to verify energy conservation.
Some corrections to improve energy conservation.
Some bug correction.

New output variables and new flags are

introduced:

(1)New variables: d_xx_col is the budget over

each column of variable "xx".

(2) fl_ebil: 0 -> no diag; 1 -> diags activated
(3) fl_cor_ebil: 0 -> no correction ensuring

energy conservation; 1 -> first set of
corrections.

+ ok_bad_ecmwf_thermo: if true, a bug setting

RVTMP2 to 0 is fixed; for backward compatibility
the default is "false".

File:

• LMDZ5/trunk/libf/phylmd/climb_hq_mod.F90 (modified) (10 diffs)

LMDZ5/trunk/libf/phylmd/climb_hq_mod.F90

@@ -9 +9 @@
   SAVE 
   PRIVATE
-  PUBLIC :: climb_hq_down, climb_hq_up
+  PUBLIC :: climb_hq_down, climb_hq_up, d_h_col_vdf, f_h_bnd
 
   REAL, DIMENSION(:,:), ALLOCATABLE :: gamaq, gamah
@@ -23 +23 @@
   REAL, DIMENSION(:,:), ALLOCATABLE :: Kcoefhq
   !$OMP THREADPRIVATE(Kcoefhq)
+  REAL, SAVE, DIMENSION(:,:), ALLOCATABLE :: h_old ! for diagnostics, h before solving diffusion
+  !$OMP THREADPRIVATE(h_old)
+  REAL, SAVE, DIMENSION(:), ALLOCATABLE :: d_h_col_vdf ! for diagnostics, vertical integral of enthalpy change
+  !$OMP THREADPRIVATE(d_h_col_vdf)
+  REAL, SAVE, DIMENSION(:), ALLOCATABLE :: f_h_bnd ! for diagnostics, enthalpy flux at surface
+  !$OMP THREADPRIVATE(d_h_bnd)
 
 CONTAINS
@@ -71 +77 @@
     LOGICAL, SAVE                            :: first=.TRUE.
     !$OMP THREADPRIVATE(first)
-    REAL, DIMENSION(klon,klev)               :: local_H
+! JLD now renamed h_old and declared in module
+!    REAL, DIMENSION(klon,klev)               :: local_H
     REAL, DIMENSION(klon)                    :: psref 
     REAL                                     :: delz, pkh
     INTEGER                                  :: k, i, ierr
-
 ! Include
 !****************************************************************************************
@@ -113 +119 @@
        ALLOCATE(gamah(1:klon,2:klev), STAT=ierr)
        IF ( ierr /= 0 ) PRINT*,' pb in allloc gamah, ierr=', ierr
+       
+       ALLOCATE(h_old(klon,klev), STAT=ierr)
+       IF ( ierr /= 0 )  PRINT*,' pb in allloc h_old, ierr=', ierr
+       
+       ALLOCATE(d_h_col_vdf(klon), STAT=ierr)
+       IF ( ierr /= 0 )  PRINT*,' pb in allloc d_h_col_vdf, ierr=', ierr
+       
+       ALLOCATE(f_h_bnd(klon), STAT=ierr)
+       IF ( ierr /= 0 )  PRINT*,' pb in allloc f_h_bnd, ierr=', ierr
     END IF
 
@@ -177 +192 @@
 !
 !****************************************************************************************
-    local_H(:,:) = 0.0
+    h_old(:,:) = 0.0
 
     DO k=1,klev
        DO i = 1, knon
           ! convertie la temperature en entalpie potentielle
-          local_H(i,k) = RCPD * temp(i,k) * &
+          h_old(i,k) = RCPD * temp(i,k) * &
                (psref(i)/pplay(i,k))**RKAPPA
        ENDDO
     ENDDO
 
-    CALL calc_coef(knon, Kcoefhq(:,:), gamah(:,:), delp(:,:), local_H(:,:), &
+    CALL calc_coef(knon, Kcoefhq(:,:), gamah(:,:), delp(:,:), h_old(:,:), &
          Ccoef_H(:,:), Dcoef_H(:,:), Acoef_H, Bcoef_H)
  
@@ -349 +364 @@
 ! 1) 
 ! Definition of some variables
+    REAL, DIMENSION(klon,klev)               :: d_h, zairm
 !
 !****************************************************************************************
@@ -355 +371 @@
     d_q(:,:)    = 0.0
     d_t(:,:)    = 0.0
+    d_h(:,:)    = 0.0
+    f_h_bnd(:)= 0.0
 
     psref(1:knon) = paprs(1:knon,1)  
@@ -393 +411 @@
     q_new(1:knon,1) = Acoef_Q(1:knon) + Bcoef_Q(1:knon)*flx_q1(1:knon)*dtime
     h_new(1:knon,1) = Acoef_H(1:knon) + Bcoef_H(1:knon)*flx_h1(1:knon)*dtime
-    
+    f_h_bnd(1:knon) = flx_h1(1:knon)
 !- All the other layers 
     DO k = 2, klev
@@ -427 +445 @@
 !
 !****************************************************************************************
-
+    d_h_col_vdf(:) = 0.0
     DO k = 1, klev
        DO i = 1, knon
           d_t(i,k) = h_new(i,k)/(psref(i)/pplay(i,k))**RKAPPA/RCPD - t_old(i,k)
           d_q(i,k) = q_new(i,k) - q_old(i,k)
-       END DO
+          d_h(i,k) = h_new(i,k) - h_old(i,k)
+!JLD          d_t(i,k) = d_h(i,k)/(psref(i)/pplay(i,k))**RKAPPA/RCPD    !correction a venir
+    ! layer air mass
+          zairm(i, k) = (paprs(i,k)-paprs(i,k+1))/rg
+          d_h_col_vdf(i) = d_h_col_vdf(i) + d_h(i,k)*zairm(i,k)
+        END DO
     END DO
 
@@ -448 +471 @@
        DEALLOCATE(Kcoefhq,stat=ierr)
        IF ( ierr /= 0 )  PRINT*,' pb in dealllocate Kcoefhq, ierr=', ierr
+       DEALLOCATE(h_old, d_h_col_vdf, f_h_bnd, stat=ierr)
+       IF ( ierr /= 0 )  PRINT*,' pb in dealllocate h_old, d_h_col_vdf, f_h_bnd, ierr=', ierr
     END IF
   END SUBROUTINE climb_hq_up