Changeset 1750

Timestamp:

Jul 24, 2017, 5:57:02 PM (7 years ago)

Author:

mlefevre

Message:

Update of module_initialize_real : implementation of Venus pressure computation and potential temperature using cp(T)

File:

• trunk/MESOSCALE/LMD_MM_MARS/SRC/WRFV2/dyn_em/module_initialize_real.F (modified) (11 diffs)

trunk/MESOSCALE/LMD_MM_MARS/SRC/WRFV2/dyn_em/module_initialize_real.F

@@ -118 +118 @@
       REAL :: max_dz
 
-!      INTEGER , PARAMETER :: nl_max = 1000
+!     INTEGER , PARAMETER :: nl_max = 1000
 !      REAL , DIMENSION(nl_max) :: grid%em_dn
 
@@ -149 +149 @@
       !LOGICAL :: interp_theta = .false. ! Wee et al. 2012 correction
       REAL :: pfu, pfd, phm
-
+      REAL :: tpot
 
 #ifdef DM_PARALLEL
@@ -1080 +1080 @@
 !!END DO
 
-
+IF ( planet == "mars" ) then
 !--get vertical size of the GCM input array and allocate new stuff
-sizegcm=SIZE(grid%em_rh_gc(1,:,1))
-ALLOCATE(sig(MIN(ite,ide-1)-its+1,sizegcm, MIN(jte,jde-1)-jts+1))
-!ALLOCATE(ap(MIN(ite,ide-1)-its+1,sizegcm, MIN(jte,jde-1)-jts+1))
-ALLOCATE(bp(MIN(ite,ide-1)-its+1,sizegcm, MIN(jte,jde-1)-jts+1))
-
-DO j = jts , MIN ( jde-1 , jte )
-DO i = its , MIN (ide-1 , ite )
+  sizegcm=SIZE(grid%em_rh_gc(1,:,1))
+  ALLOCATE(sig(MIN(ite,ide-1)-its+1,sizegcm, MIN(jte,jde-1)-jts+1))
+  !ALLOCATE(ap(MIN(ite,ide-1)-its+1,sizegcm, MIN(jte,jde-1)-jts+1))
+  ALLOCATE(bp(MIN(ite,ide-1)-its+1,sizegcm, MIN(jte,jde-1)-jts+1))
+
+  DO j = jts , MIN ( jde-1 , jte )
+  DO i = its , MIN (ide-1 , ite )
 
 !!! Define old sigma levels for each column
-        sig(i,:,j)=grid%em_p_gc(i,:,j)/grid%em_psfc_gc(i,j)
+          sig(i,:,j)=grid%em_p_gc(i,:,j)/grid%em_psfc_gc(i,j)
 
 !!! Compute new sigma levels from old sigma levels with GCM (low-res) and WRF (hi-res) surface pressure 
 !!!                        (dimlevs,sigma_gcm, ps_gcm,                ps_hr,         sigma_hr)
-        CALL build_sigma_hr(sizegcm,sig(i,:,j),grid%em_psfc_gc(i,j),grid%psfc(i,j),bp(i,:,j))
+          CALL build_sigma_hr(sizegcm,sig(i,:,j),grid%em_psfc_gc(i,j),grid%psfc(i,j),bp(i,:,j))
 
 !!! Calculate new pressure levels
-        grid%em_pd_gc(i,:,j)=bp(i,:,j)*grid%psfc(i,j)
-
-END DO
-END DO
-
-DEALLOCATE(sig)
-DEALLOCATE(bp)
+          grid%em_pd_gc(i,:,j)=bp(i,:,j)*grid%psfc(i,j)
+
+  END DO
+  END DO
+
+  DEALLOCATE(sig)
+  DEALLOCATE(bp)
 
 !!****MARS who knows...
-grid%em_rh_gc(:,:,:)=0.
+  grid%em_rh_gc(:,:,:)=0.
 
 
@@ -1113 +1113 @@
 !grid%em_pd_gc=grid%em_p_gc
 !!****MARS
-
+ELSE ! VENUS
 
 
@@ -1120 +1120 @@
          !!  dry top pressure (constant).
          !
-         !CALL p_dts ( grid%em_mu0 , grid%em_intq_gc , grid%psfc , grid%p_top , &
-         !             ids , ide , jds , jde , 1   , num_metgrid_levels , &
-         !             ims , ime , jms , jme , 1   , num_metgrid_levels , &
-         !             its , ite , jts , jte , 1   , num_metgrid_levels )
-
+           CALL p_dts ( grid%em_mu0 , grid%em_intq_gc , grid%psfc , grid%p_top , &
+                        ids , ide , jds , jde , 1   , num_metgrid_levels , &
+                        ims , ime , jms , jme , 1   , num_metgrid_levels , &
+                        its , ite , jts , jte , 1   , num_metgrid_levels )
+ENDIF
+IF ( planet == "mars" ) then
 !!****MARS
-DO j = jts , MIN ( jde-1 , jte )
-DO i = its , MIN (ide-1 , ite )
-
-   grid%em_mu0(i,j) = grid%psfc(i,j) - grid%p_top
-
-END DO
-END DO
+  DO j = jts , MIN ( jde-1 , jte )
+  DO i = its , MIN (ide-1 , ite )
+
+     grid%em_mu0(i,j) = grid%psfc(i,j) - grid%p_top
+
+  END DO
+  END DO
 !!****MARS
-
+ELSE ! VENUS
    
          !!  Compute the dry, hydrostatic surface pressure.
          !
-         !CALL p_dhs ( grid%em_pdhs , grid%ht , p00 , t00 , a , &
-         !             ids , ide , jds , jde , kds , kde , &
-         !             ims , ime , jms , jme , kms , kme , &
-         !             its , ite , jts , jte , kts , kte )
+           CALL p_dhs ( grid%em_pdhs , grid%ht , p00 , t00 , a , &
+                        ids , ide , jds , jde , kds , kde , &
+                        ims , ime , jms , jme , kms , kme , &
+                        its , ite , jts , jte , kts , kte )
+ENDIF
 !!****MARS: voir remarques dans la routine
 !!****MARS: dry hydrostatic pressure comes from the GCM ...
@@ -2294 +2296 @@
 IF (( i .EQ. its ) .AND. ( j .EQ. jts )) print *, temp, k
 !!! MODIF WRFV3.1 - parameter tiso
-               grid%em_t_init(i,k,j) = temp*(p00/grid%em_pb(i,k,j))**(r_d/cp) - t0
+               IF (planet .eq. "mars" ) THEN
+                 grid%em_t_init(i,k,j) = temp*(p00/grid%em_pb(i,k,j))**(r_d/cp) - t0
+               ELSE
+                 grid%em_t_init(i,k,j) = (temp**nu + nu*(TT00**nu)*log((p00/grid%em_pb(i,k,j))**rcp))**(1/nu) -t0
+               ENDIF
                grid%em_alb(i,k,j) = (r_d/p1000mb)*(grid%em_t_init(i,k,j)+t0)*(grid%em_pb(i,k,j)/p1000mb)**cvpm
             END DO
@@ -2780 +2786 @@
 temp1 = MAX(tiso,t00+A*LOG(grid%em_pb(i,k,j)/p00))
 temp2 = MAX(tiso,t00+A*LOG(           pb_int/p00))
-grid%em_t_init(i,k,j) = temp1*(p00/grid%em_pb(i,k,j))**(r_d/cp) - t0
-t_init_int(i,k,j)     = temp2*(p00/pb_int           )**(r_d/cp) - t0
+IF (planet .eq. "mars" ) THEN
+  grid%em_t_init(i,k,j) = temp1*(p00/grid%em_pb(i,k,j))**(r_d/cp) - t0
+  t_init_int(i,k,j)     = temp2*(p00/pb_int           )**(r_d/cp) - t0
+ELSE
+  grid%em_t_init(i,k,j) = (temp1**nu + nu*(TT00**nu)*log((p00/grid%em_pb(i,k,j))**(rcp)))**(1/nu) - t0
+  t_init_int(i,k,j)     = (temp2**nu + nu*(TT00**nu)*log((p00/pb_int)**(rcp)))**(1/nu) - t0
+ENDIF
                grid%em_alb(i,k,j) = (r_d/p1000mb)*(grid%em_t_init(i,k,j)+t0)*(grid%em_pb(i,k,j)/p1000mb)**cvpm
             END DO
@@ -4969 +4980 @@
 !            temp =             t00 + A*LOG(pb/p00)
 temp = MAX ( tiso, t00 + A*LOG(pb/p00) )
-            t_init = temp*(p00/pb)**(r_d/cp) - t0
+           IF (planet .eq. "mars" ) THEN
+              t_init = temp*(p00/pb)**(r_d/cp) - t0
+            ELSE
+              t_init = (temp**nu + nu*(TT00**nu)*log((p00/pb)**(rcp)))**(1/nu) - t0
+            ENDIF
             alb(k) = (r_d/p1000mb)*(t_init+t0)*(pb/p1000mb)**cvpm
          END DO
@@ -5005 +5020 @@
 !            temp =             t00 + A*LOG(pb/p00)
 temp = MAX ( tiso, t00 + A*LOG(pb/p00) )
-            t_init = temp*(p00/pb)**(r_d/cp) - t0
+           IF (planet .eq. "mars" ) THEN
+              t_init = temp*(p00/pb)**(r_d/cp) - t0
+            ELSE
+              t_init = (temp**nu + nu*(TT00**nu)*log((p00/pb)**(rcp)))**(1/nu) -t0
+            ENDIF
             alb(k) = (r_d/p1000mb)*(t_init+t0)*(pb/p1000mb)**cvpm
             znw(k+1) = znw(k) - dz*g / ( mub*alb(k) )
@@ -5025 +5044 @@
 !                  temp =             t00 + A*LOG(pb/p00)
 temp = MAX ( tiso, t00 + A*LOG(pb/p00) )
-                  t_init = temp*(p00/pb)**(r_d/cp) - t0
+           IF (planet .eq. "mars" ) THEN
+              t_init = temp*(p00/pb)**(r_d/cp) - t0
+            ELSE
+              t_init = (temp**nu + nu*(TT00**nu)*log((p00/pb)**(rcp)))**(1/nu) -t0
+            ENDIF
                   alb(k) = (r_d/p1000mb)*(t_init+t0)*(pb/p1000mb)**cvpm
                   znw(k+1) = znw(k) - dz*g / ( mub*alb(k) )
@@ -5052 +5075 @@
 !               temp =             t00 + A*LOG(pb/p00)
 temp = MAX ( tiso, t00 + A*LOG(pb/p00) )
-               t_init = temp*(p00/pb)**(r_d/cp) - t0
+           IF (planet .eq. "mars" ) THEN
+              t_init = temp*(p00/pb)**(r_d/cp) - t0
+            ELSE
+              t_init = (temp**nu + nu*(TT00**nu)*log((p00/pb)**(rcp)))**(1/nu) -t0
+            ENDIF
                alb(k) = (r_d/p1000mb)*(t_init+t0)*(pb/p1000mb)**cvpm
             END DO