Changeset 486 for LMDZ.3.3/branches/rel-LF/libf/phylmd/clmain.F

Timestamp:

Dec 15, 2003, 6:50:41 PM (21 years ago)

Author:

lmdzadmin

Message:

Phasage avec la version de Ionela
IM/LF

File:

• LMDZ.3.3/branches/rel-LF/libf/phylmd/clmain.F (modified) (10 diffs)

LMDZ.3.3/branches/rel-LF/libf/phylmd/clmain.F

@@ -7 +7 @@
      .                  jour, rmu0, co2_ppm,
      .                  ok_veget, ocean, npas, nexca, ts,
-     .                  soil_model,cdmmax, cdhmax, ftsoil,qsol,
+     .                  soil_model,cdmmax, cdhmax,
+     .                  ksta, ksta_ter, ok_kzmin, ftsoil,qsol,
      .                  paprs,pplay,radsol,snow,qsurf,evap,albe,alblw,
      .                  fluxlat,
@@ -143 +144 @@
 cIM ajout seuils cdrm, cdrh
       REAL cdmmax, cdhmax
+cIM: 261103
+      REAL ksta, ksta_ter
+      LOGICAL ok_kzmin
+cIM: 261103
       REAL ftsoil(klon,nsoilmx,nbsrf)
       REAL ytsoil(klon,nsoilmx)
@@ -472 +477 @@
 c calculer Cdrag et les coefficients d'echange
       CALL coefkz(nsrf, knon, ypaprs, ypplay,
+cIM 261103
+     .     ksta, ksta_ter,
+cIM 261103
      .            yts, yrugos, yu, yv, yt, yq,
      .            yqsurf, 
@@ -494 +502 @@
       endif
 
+c
+cIM: 261103
+      if (ok_kzmin) THEN
+cIM cf FH: 201103 BEG
+c   Calcul d'une diffusion minimale pour les conditions tres stables.
+      call coefkzmin(knon,ypaprs,ypplay,yu,yv,yt,yq,ycoefm
+     .   ,ycoefm0,ycoefh0)
+c      call dump2d(iim,jjm-1,ycoefm(2:klon-1,2), 'KZ         ')
+c      call dump2d(iim,jjm-1,ycoefm0(2:klon-1,2),'KZMIN      ')
+ 
+       if ( 1.eq.1 ) then
+       DO k = 1, klev
+       DO i = 1, knon
+          ycoefm(i,k) = MAX(ycoefm(i,k),ycoefm0(i,k))
+          ycoefh(i,k) = MAX(ycoefh(i,k),ycoefh0(i,k))
+       ENDDO
+       ENDDO
+       endif
+cIM cf FH: 201103 END
+      endif !ok_kzmin
+cIM: 261103
+
+c
 c calculer la diffusion des vitesses "u" et "v"
       CALL clvent(knon,dtime,yu1,yv1,ycoefm,yt,yu,ypaprs,ypplay,ydelp,
@@ -503 +534 @@
       ytaux = y_flux_u(:,1)
       ytauy = y_flux_v(:,1)
-
-      if (nsrf.eq.is_oce) then
-         do j=1,knon
-            ycoefm(j,1)=1.e-3
-c            ycoefh(j,1)=0.8*ycoefm(j,1)
-         enddo
-      endif
 
 c FH modif sur le cdrag temperature
@@ -737 +761 @@
 #else
        DO j=1, knon
+         i = ni(j) 
          t2m(i,nsrf)=0.
          q2m(i,nsrf)=0.
@@ -1219 +1244 @@
       END
       SUBROUTINE coefkz(nsrf, knon, paprs, pplay,
+cIM 261103
+     .                  ksta, ksta_ter,
+cIM 261103
      .                  ts, rugos,
      .                  u,v,t,q,
@@ -1279 +1307 @@
       PARAMETER (prandtl=0.4)
       REAL kstable ! diffusion minimale (situation stable)
-      PARAMETER (kstable=1.0e-10)
+      ! GKtest
+      ! PARAMETER (kstable=1.0e-10)
+      REAL ksta, ksta_ter
+cIM: 261103     REAL kstable_ter, kstable_sinon
+cIM: 211003 cf GK   PARAMETER (kstable_ter = 1.0e-6)
+cIM: 261103     PARAMETER (kstable_ter = 1.0e-8)
+cIM: 261103   PARAMETER (kstable_ter = 1.0e-10)
+cIM: 261103   PARAMETER (kstable_sinon = 1.0e-10)
+      ! fin GKtest
       REAL mixlen ! constante controlant longueur de melange
       PARAMETER (mixlen=35.0)
@@ -1371 +1407 @@
          gamt(2) = -2.5E-03
       ENDIF
+cIM cf JLD/ GKtest
+      IF ( nsrf .EQ. is_ter ) THEN
+cIM 261103     kstable = kstable_ter
+        kstable = ksta_ter
+      ELSE
+cIM 261103     kstable = kstable_sinon
+        kstable = ksta
+      ENDIF
+cIM cf JLD/ GKtest fin
 c
 c Calculer les geopotentiels de chaque couche
@@ -1614 +1659 @@
       DO k = 2, klev
       DO i = 1, knon
-      IF ( (nsrf.NE.is_oce) .OR.  ! si ce n'est pas sur l'ocean
-     .     (invb(i).EQ.klev) .OR. ! s'il n'y a pas d'inversion
-     .     (zdthmin(i).GT.seuil) )THEN ! si l'inversion est trop faible
+cIM cf FH/GK   IF ( (nsrf.NE.is_oce) .OR.  ! si ce n'est pas sur l'ocean
+cIM cf FH/GK  .     (invb(i).EQ.klev) .OR. ! s'il n'y a pas d'inversion
+      !IM cf JLD/ GKtest TERkz2
+      ! IF ( (nsrf.EQ.is_ter) .OR.  ! si on est sur la terre
+      ! fin GKtest
+      IF ( (nsrf.EQ.is_oce) .AND.  ! si on est sur ocean et si 
+     .     ( (invb(i).EQ.klev) .OR.      ! s'il n'y a pas d'inversion
+     .     (zdthmin(i).GT.seuil) ) )THEN ! si l'inversion est trop faible
          zl2(i)=(mixlen*MAX(0.0,(paprs(i,k)-paprs(i,klev+1))
      .                       /(paprs(i,2)-paprs(i,klev+1)) ))**2