/[lmdze]/trunk/Sources/phylmd/clmain.f

Diff of /trunk/Sources/phylmd/clmain.f

Parent Directory | Revision Log | View Patch Patch

-revision 221 by guez,
Thu Apr 20 14:44:47 2017 UTC
+revision 248 by guez,
Fri Jan  5 16:40:13 2018 UTC
 Line 6 
 contains
    SUBROUTINE clmain(dtime, pctsrf, t, q, u, v, julien, mu0, ftsol, cdmmax, &
         cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, qsol, paprs, pplay, fsnow, &
-        qsurf, evap, falbe, fluxlat, rain_fall, snow_f, solsw, sollw, fder, &
+        qsurf, evap, falbe, fluxlat, rain_fall, snow_f, fsolsw, fsollw, frugs, &
-        rugos, agesno, rugoro, d_t, d_q, d_u, d_v, d_ts, flux_t, flux_q, &
+        agesno, rugoro, d_t, d_q, d_u, d_v, d_ts, flux_t, flux_q, flux_u, &
-        flux_u, flux_v, cdragh, cdragm, q2, dflux_t, dflux_q, ycoefh, zu1, &
+        flux_v, cdragh, cdragm, q2, dflux_t, dflux_q, coefh, t2m, q2m, &
-        zv1, t2m, q2m, u10m, v10m, pblh, capcl, oliqcl, cteicl, pblt, therm, &
+        u10m_srf, v10m_srf, pblh, capcl, oliqcl, cteicl, pblt, therm, trmb1, &
-        trmb1, trmb2, trmb3, plcl, fqcalving, ffonte, run_off_lic_0)
+        trmb2, trmb3, plcl, fqcalving, ffonte, run_off_lic_0)
      ! From phylmd/clmain.F, version 1.6, 2005/11/16 14:47:19
      ! Author: Z. X. Li (LMD/CNRS), date: 1993/08/18
 Line 21 
 contains
      ! ne tient pas compte de la diff\'erentiation des sous-fractions
      ! de sol.
-     ! Pour pouvoir extraire les coefficients d'\'echanges et le vent
+     use clcdrag_m, only: clcdrag
-     ! dans la premi\`ere couche, trois champs ont \'et\'e cr\'e\'es : "ycoefh",
-     ! "zu1" et "zv1". Nous avons moyenn\'e les valeurs de ces trois
-     ! champs sur les quatre sous-surfaces du mod\`ele.
      use clqh_m, only: clqh
      use clvent_m, only: clvent
      use coefkz_m, only: coefkz
      use coefkzmin_m, only: coefkzmin
-     USE conf_gcm_m, ONLY: prt_level, lmt_pas
+     use coefkz2_m, only: coefkz2
+     USE conf_gcm_m, ONLY: lmt_pas
      USE conf_phys_m, ONLY: iflag_pbl
      USE dimphy, ONLY: klev, klon, zmasq
      USE dimsoil, ONLY: nsoilmx
-Line 41 
 contains
+Line 38 
 contains
      USE suphec_m, ONLY: rd, rg, rkappa
      use time_phylmdz, only: itap
      use ustarhb_m, only: ustarhb
-     use vdif_kcay_m, only: vdif_kcay
      use yamada4_m, only: yamada4
      REAL, INTENT(IN):: dtime ! interval du temps (secondes)
-Line 50 
 contains
+Line 46 
 contains
      ! tableau des pourcentages de surface de chaque maille
      REAL, INTENT(IN):: t(klon, klev) ! temperature (K)
-     REAL, INTENT(IN):: q(klon, klev) ! vapeur d'eau (kg/kg)
+     REAL, INTENT(IN):: q(klon, klev) ! vapeur d'eau (kg / kg)
      REAL, INTENT(IN):: u(klon, klev), v(klon, klev) ! vitesse
      INTEGER, INTENT(IN):: julien ! jour de l'annee en cours
      REAL, intent(in):: mu0(klon) ! cosinus de l'angle solaire zenithal
-     REAL, INTENT(IN):: ftsol(klon, nbsrf) ! temp\'erature du sol (en K)
+     REAL, INTENT(IN):: ftsol(:, :) ! (klon, nbsrf) temp\'erature du sol (en K)
      REAL, INTENT(IN):: cdmmax, cdhmax ! seuils cdrm, cdrh
      REAL, INTENT(IN):: ksta, ksta_ter
      LOGICAL, INTENT(IN):: ok_kzmin
-Line 62 
 contains
+Line 58 
 contains
      REAL, INTENT(inout):: ftsoil(klon, nsoilmx, nbsrf)
      ! soil temperature of surface fraction
-     REAL, INTENT(inout):: qsol(klon)
+     REAL, INTENT(inout):: qsol(:) ! (klon)
      ! column-density of water in soil, in kg m-2
-     REAL, INTENT(IN):: paprs(klon, klev+1) ! pression a intercouche (Pa)
+     REAL, INTENT(IN):: paprs(klon, klev + 1) ! pression a intercouche (Pa)
      REAL, INTENT(IN):: pplay(klon, klev) ! pression au milieu de couche (Pa)
      REAL, INTENT(inout):: fsnow(:, :) ! (klon, nbsrf) \'epaisseur neigeuse
      REAL qsurf(klon, nbsrf)
-Line 74 
 contains
+Line 70 
 contains
      REAL, intent(out):: fluxlat(:, :) ! (klon, nbsrf)
      REAL, intent(in):: rain_fall(klon)
-     ! liquid water mass flux (kg/m2/s), positive down
+     ! liquid water mass flux (kg / m2 / s), positive down
      REAL, intent(in):: snow_f(klon)
-     ! solid water mass flux (kg/m2/s), positive down
+     ! solid water mass flux (kg / m2 / s), positive down
-     REAL, INTENT(IN):: solsw(klon, nbsrf), sollw(klon, nbsrf)
+     REAL, INTENT(IN):: fsolsw(klon, nbsrf), fsollw(klon, nbsrf)
-     REAL, intent(in):: fder(klon)
+     REAL, intent(inout):: frugs(klon, nbsrf) ! longueur de rugosit\'e (en m)
-     REAL, intent(inout):: rugos(klon, nbsrf) ! longueur de rugosit\'e (en m)
      real agesno(klon, nbsrf)
      REAL, INTENT(IN):: rugoro(klon)
-Line 95 
 contains
+Line 90 
 contains
      REAL, intent(out):: d_ts(:, :) ! (klon, nbsrf) variation of ftsol
      REAL, intent(out):: flux_t(klon, nbsrf)
-     ! flux de chaleur sensible (Cp T) (W/m2) (orientation positive vers
+     ! flux de chaleur sensible (Cp T) (W / m2) (orientation positive vers
      ! le bas) à la surface
      REAL, intent(out):: flux_q(klon, nbsrf)
-     ! flux de vapeur d'eau (kg/m2/s) à la surface
+     ! flux de vapeur d'eau (kg / m2 / s) à la surface
      REAL, intent(out):: flux_u(klon, nbsrf), flux_v(klon, nbsrf)
-     ! tension du vent à la surface, en Pa
+     ! tension du vent (flux turbulent de vent) à la surface, en Pa
      REAL, INTENT(out):: cdragh(klon), cdragm(klon)
-     real q2(klon, klev+1, nbsrf)
+     real q2(klon, klev + 1, nbsrf)
      REAL, INTENT(out):: dflux_t(klon), dflux_q(klon)
      ! dflux_t derive du flux sensible
      ! dflux_q derive du flux latent
      ! IM "slab" ocean
-     REAL, intent(out):: ycoefh(klon, klev)
+     REAL, intent(out):: coefh(:, 2:) ! (klon, 2:klev)
-     REAL, intent(out):: zu1(klon)
+     ! Pour pouvoir extraire les coefficients d'\'echange, le champ
-     REAL zv1(klon)
+     ! "coefh" a \'et\'e cr\'e\'e. Nous avons moyenn\'e les valeurs de
+     ! ce champ sur les quatre sous-surfaces du mod\`ele.
      REAL, INTENT(inout):: t2m(klon, nbsrf), q2m(klon, nbsrf)
-     REAL u10m(klon, nbsrf), v10m(klon, nbsrf)
-     ! Ionela Musat cf. Anne Mathieu : planetary boundary layer, hbtm
+     REAL, INTENT(inout):: u10m_srf(:, :), v10m_srf(:, :) ! (klon, nbsrf)
-     ! (Comme les autres diagnostics on cumule dans physiq ce qui
+     ! composantes du vent \`a 10m sans spirale d'Ekman
-     ! permet de sortir les grandeurs par sous-surface)
+     ! Ionela Musat. Cf. Anne Mathieu : planetary boundary layer, hbtm.
+     ! Comme les autres diagnostics on cumule dans physiq ce qui permet
+     ! de sortir les grandeurs par sous-surface.
      REAL pblh(klon, nbsrf) ! height of planetary boundary layer
      REAL capcl(klon, nbsrf)
      REAL oliqcl(klon, nbsrf)
-Line 137 
 contains
+Line 136 
 contains
      REAL fqcalving(klon, nbsrf), ffonte(klon, nbsrf)
      ! ffonte----Flux thermique utilise pour fondre la neige
      ! fqcalving-Flux d'eau "perdue" par la surface et necessaire pour limiter la
-     !           hauteur de neige, en kg/m2/s
+     !           hauteur de neige, en kg / m2 / s
      REAL run_off_lic_0(klon)
      ! Local:
-Line 152 
 contains
+Line 151 
 contains
      real y_run_off_lic_0(klon)
      REAL rugmer(klon)
      REAL ytsoil(klon, nsoilmx)
-     REAL yts(klon), yrugos(klon), ypct(klon), yz0_new(klon)
+     REAL yts(klon), ypct(klon), yz0_new(klon)
+     real yrugos(klon) ! longeur de rugosite (en m)
      REAL yalb(klon)
-     REAL yu1(klon), yv1(klon)
-     ! On ajoute en output yu1 et yv1 qui sont les vents dans
-     ! la premi\`ere couche.
      REAL snow(klon), yqsurf(klon), yagesno(klon)
+     real yqsol(klon) ! column-density of water in soil, in kg m-2
-     real yqsol(klon)
+     REAL yrain_f(klon) ! liquid water mass flux (kg / m2 / s), positive down
-     ! column-density of water in soil, in kg m-2
+     REAL ysnow_f(klon) ! solid water mass flux (kg / m2 / s), positive down
-     REAL yrain_f(klon)
-     ! liquid water mass flux (kg/m2/s), positive down
-     REAL ysnow_f(klon)
-     ! solid water mass flux (kg/m2/s), positive down
-     REAL yfder(klon)
      REAL yrugm(klon), yrads(klon), yrugoro(klon)
      REAL yfluxlat(klon)
      REAL y_d_ts(klon)
-Line 179 
 contains
+Line 166 
 contains
      REAL y_flux_t(klon), y_flux_q(klon)
      REAL y_flux_u(klon), y_flux_v(klon)
      REAL y_dflux_t(klon), y_dflux_q(klon)
-     REAL coefh(klon, klev), coefm(klon, klev)
+     REAL ycoefh(klon, 2:klev), ycoefm(klon, 2:klev)
+     real ycdragh(klon), ycdragm(klon)
      REAL yu(klon, klev), yv(klon, klev)
      REAL yt(klon, klev), yq(klon, klev)
-     REAL ypaprs(klon, klev+1), ypplay(klon, klev), ydelp(klon, klev)
+     REAL ypaprs(klon, klev + 1), ypplay(klon, klev), ydelp(klon, klev)
+     REAL ycoefm0(klon, 2:klev), ycoefh0(klon, 2:klev)
-     REAL ycoefm0(klon, klev), ycoefh0(klon, klev)
+     REAL yzlay(klon, klev), zlev(klon, klev + 1), yteta(klon, klev)
+     REAL yq2(klon, klev + 1)
-     REAL yzlay(klon, klev), yzlev(klon, klev+1), yteta(klon, klev)
-     REAL ykmm(klon, klev+1), ykmn(klon, klev+1)
-     REAL ykmq(klon, klev+1)
-     REAL yq2(klon, klev+1)
-     REAL q2diag(klon, klev+1)
-     REAL u1lay(klon), v1lay(klon)
      REAL delp(klon, klev)
      INTEGER i, k, nsrf
      INTEGER ni(klon), knon, j
      REAL pctsrf_pot(klon, nbsrf)
      ! "pourcentage potentiel" pour tenir compte des \'eventuelles
      ! apparitions ou disparitions de la glace de mer
-     REAL zx_alf1, zx_alf2 ! valeur ambiante par extrapolation
+     REAL yt2m(klon), yq2m(klon), wind10m(klon)
+     REAL ustar(klon)
-     REAL yt2m(klon), yq2m(klon), yu10m(klon)
-     REAL yustar(klon)
      REAL yt10m(klon), yq10m(klon)
      REAL ypblh(klon)
-Line 218 
 contains
+Line 196 
 contains
      REAL ytrmb1(klon)
      REAL ytrmb2(klon)
      REAL ytrmb3(klon)
-     REAL uzon(klon), vmer(klon)
+     REAL u1(klon), v1(klon)
      REAL tair1(klon), qair1(klon), tairsol(klon)
      REAL psfce(klon), patm(klon)
      REAL qairsol(klon), zgeo1(klon)
      REAL rugo1(klon)
+     REAL zgeop(klon, klev)
-     ! utiliser un jeu de fonctions simples
-     LOGICAL zxli
-     PARAMETER (zxli=.FALSE.)
      !------------------------------------------------------------
-Line 235 
 contains
+Line 210 
 contains
      DO k = 1, klev ! epaisseur de couche
         DO i = 1, klon
-           delp(i, k) = paprs(i, k) - paprs(i, k+1)
+           delp(i, k) = paprs(i, k) - paprs(i, k + 1)
         END DO
      END DO
-     DO i = 1, klon ! vent de la premiere couche
-        zx_alf1 = 1.0
-        zx_alf2 = 1.0 - zx_alf1
-        u1lay(i) = u(i, 1)*zx_alf1 + u(i, 2)*zx_alf2
-        v1lay(i) = v(i, 1)*zx_alf1 + v(i, 2)*zx_alf2
-     END DO
      ! Initialization:
      rugmer = 0.
-Line 251 
 contains
+Line 220 
 contains
      cdragm = 0.
      dflux_t = 0.
      dflux_q = 0.
-     zu1 = 0.
-     zv1 = 0.
      ypct = 0.
      yqsurf = 0.
      yrain_f = 0.
      ysnow_f = 0.
-     yfder = 0.
      yrugos = 0.
-     yu1 = 0.
-     yv1 = 0.
-     yrads = 0.
      ypaprs = 0.
      ypplay = 0.
      ydelp = 0.
-Line 282 
 contains
+Line 245 
 contains
      d_q = 0.
      d_u = 0.
      d_v = 0.
-     ycoefh = 0.
+     coefh = 0.
      ! Initialisation des "pourcentages potentiels". On consid\`ere ici qu'on
      ! peut avoir potentiellement de la glace sur tout le domaine oc\'eanique
-Line 324 
 contains
+Line 287 
 contains
               yrain_f(j) = rain_fall(i)
               ysnow_f(j) = snow_f(i)
               yagesno(j) = agesno(i, nsrf)
-              yfder(j) = fder(i)
+              yrugos(j) = frugs(i, nsrf)
-              yrugos(j) = rugos(i, nsrf)
               yrugoro(j) = rugoro(i)
-              yu1(j) = u1lay(i)
+              yrads(j) = fsolsw(i, nsrf) + fsollw(i, nsrf)
-              yv1(j) = v1lay(i)
+              ypaprs(j, klev + 1) = paprs(i, klev + 1)
-              yrads(j) = solsw(i, nsrf) + sollw(i, nsrf)
-              ypaprs(j, klev+1) = paprs(i, klev+1)
               y_run_off_lic_0(j) = run_off_lic_0(i)
            END DO
            ! For continent, copy soil water content
-           IF (nsrf == is_ter) THEN
+           IF (nsrf == is_ter) yqsol(:knon) = qsol(ni(:knon))
-              yqsol(:knon) = qsol(ni(:knon))
-           ELSE
-              yqsol = 0.
-           END IF
            ytsoil(:knon, :) = ftsoil(ni(:knon), :, nsrf)
-Line 356 
 contains
+Line 312 
 contains
               END DO
            END DO
-           ! calculer Cdrag et les coefficients d'echange
+           ! Calculer les géopotentiels de chaque couche:
-           CALL coefkz(nsrf, ypaprs, ypplay, ksta, ksta_ter, yts(:knon), &
-                yrugos, yu, yv, yt, yq, yqsurf(:knon), coefm(:knon, :), &
+           zgeop(:knon, 1) = RD * yt(:knon, 1) / (0.5 * (ypaprs(:knon, 1) &
-                coefh(:knon, :))
+                + ypplay(:knon, 1))) * (ypaprs(:knon, 1) - ypplay(:knon, 1))
+           DO k = 2, klev
+              zgeop(:knon, k) = zgeop(:knon, k - 1) + RD * 0.5 &
+                   * (yt(:knon, k - 1) + yt(:knon, k)) / ypaprs(:knon, k) &
+                   * (ypplay(:knon, k - 1) - ypplay(:knon, k))
+           ENDDO
+           CALL clcdrag(nsrf, yu(:knon, 1), yv(:knon, 1), yt(:knon, 1), &
+                yq(:knon, 1), zgeop(:knon, 1), yts(:knon), yqsurf(:knon), &
+                yrugos(:knon), ycdragm(:knon), ycdragh(:knon))
+           CALL coefkz(nsrf, ypaprs(:knon, :), ypplay(:knon, :), ksta, &
+                ksta_ter, yts(:knon), yu(:knon, :), yv(:knon, :), yt(:knon, :), &
+                yq(:knon, :), zgeop(:knon, :), ycoefm(:knon, :), &
+                ycoefh(:knon, :))
            IF (iflag_pbl == 1) THEN
-              CALL coefkz2(nsrf, knon, ypaprs, ypplay, yt, ycoefm0, ycoefh0)
+              CALL coefkz2(nsrf, knon, ypaprs, ypplay, yt, ycoefm0(:knon, :), &
-              coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))
+                   ycoefh0(:knon, :))
-              coefh(:knon, :) = max(coefh(:knon, :), ycoefh0(:knon, :))
+              ycoefm(:knon, :) = max(ycoefm(:knon, :), ycoefm0(:knon, :))
+              ycoefh(:knon, :) = max(ycoefh(:knon, :), ycoefh0(:knon, :))
+              ycdragm(:knon) = max(ycdragm(:knon), 0.)
+              ycdragh(:knon) = max(ycdragh(:knon), 0.)
            END IF
-           ! on met un seuil pour coefm et coefh
+           ! on met un seuil pour ycdragm et ycdragh
            IF (nsrf == is_oce) THEN
-              coefm(:knon, 1) = min(coefm(:knon, 1), cdmmax)
+              ycdragm(:knon) = min(ycdragm(:knon), cdmmax)
-              coefh(:knon, 1) = min(coefh(:knon, 1), cdhmax)
+              ycdragh(:knon) = min(ycdragh(:knon), cdhmax)
            END IF
            IF (ok_kzmin) THEN
               ! Calcul d'une diffusion minimale pour les conditions tres stables
               CALL coefkzmin(knon, ypaprs, ypplay, yu, yv, yt, yq, &
-                   coefm(:knon, 1), ycoefm0, ycoefh0)
+                   ycdragm(:knon), ycoefh0(:knon, :))
-              coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))
+              ycoefm0(:knon, :) = ycoefh0(:knon, :)
-              coefh(:knon, :) = max(coefh(:knon, :), ycoefh0(:knon, :))
+              ycoefm(:knon, :) = max(ycoefm(:knon, :), ycoefm0(:knon, :))
+              ycoefh(:knon, :) = max(ycoefh(:knon, :), ycoefh0(:knon, :))
            END IF
-           IF (iflag_pbl >= 3) THEN
+           IF (iflag_pbl >= 6) THEN
               ! Mellor et Yamada adapt\'e \`a Mars, Richard Fournier et
               ! Fr\'ed\'eric Hourdin
               yzlay(:knon, 1) = rd * yt(:knon, 1) / (0.5 * (ypaprs(:knon, 1) &
                    + ypplay(:knon, 1))) &
                    * (ypaprs(:knon, 1) - ypplay(:knon, 1)) / rg
               DO k = 2, klev
-                 yzlay(1:knon, k) = yzlay(1:knon, k-1) &
+                 yzlay(:knon, k) = yzlay(:knon, k-1) &
                       + rd * 0.5 * (yt(1:knon, k-1) + yt(1:knon, k)) &
                       / ypaprs(1:knon, k) &
                       * (ypplay(1:knon, k-1) - ypplay(1:knon, k)) / rg
               END DO
               DO k = 1, klev
-                 yteta(1:knon, k) = yt(1:knon, k)*(ypaprs(1:knon, 1) &
+                 yteta(1:knon, k) = yt(1:knon, k) * (ypaprs(1:knon, 1) &
-                      / ypplay(1:knon, k))**rkappa * (1.+0.61*yq(1:knon, k))
+                      / ypplay(1:knon, k))**rkappa * (1. + 0.61 * yq(1:knon, k))
               END DO
-              yzlev(1:knon, 1) = 0.
-              yzlev(:knon, klev+1) = 2. * yzlay(:knon, klev) &
+              zlev(:knon, 1) = 0.
+              zlev(:knon, klev + 1) = 2. * yzlay(:knon, klev) &
                    - yzlay(:knon, klev - 1)
               DO k = 2, klev
-                 yzlev(1:knon, k) = 0.5*(yzlay(1:knon, k)+yzlay(1:knon, k-1))
+                 zlev(:knon, k) = 0.5 * (yzlay(:knon, k) + yzlay(:knon, k-1))
               END DO
               DO k = 1, klev + 1
                  DO j = 1, knon
                     i = ni(j)
-Line 409 
 contains
+Line 390 
 contains
                  END DO
               END DO
-              CALL ustarhb(knon, yu, yv, coefm(:knon, 1), yustar)
+              ustar(:knon) = ustarhb(yu(:knon, 1), yv(:knon, 1), ycdragm(:knon))
-              IF (prt_level > 9) PRINT *, 'USTAR = ', yustar
+              CALL yamada4(dtime, rg, zlev(:knon, :), yzlay(:knon, :), &
+                   yu(:knon, :), yv(:knon, :), yteta(:knon, :), yq2(:knon, :), &
-              ! iflag_pbl peut \^etre utilis\'e comme longueur de m\'elange
+                   ycoefm(:knon, :), ycoefh(:knon, :), ustar(:knon))
-              IF (iflag_pbl >= 11) THEN
-                 CALL vdif_kcay(knon, dtime, rg, ypaprs, yzlev, yzlay, yu, yv, &
-                      yteta, coefm(:knon, 1), yq2, q2diag, ykmm, ykmn, yustar, &
-                      iflag_pbl)
-              ELSE
-                 CALL yamada4(knon, dtime, rg, yzlev, yzlay, yu, yv, yteta, &
-                      coefm(:knon, 1), yq2, ykmm, ykmn, ykmq, yustar, iflag_pbl)
-              END IF
-              coefm(:knon, 2:) = ykmm(:knon, 2:klev)
-              coefh(:knon, 2:) = ykmn(:knon, 2:klev)
            END IF
-           ! calculer la diffusion des vitesses "u" et "v"
+           CALL clvent(dtime, yu(:knon, 1), yv(:knon, 1), ycoefm(:knon, :), &
-           CALL clvent(knon, dtime, yu1, yv1, coefm(:knon, :), yt, yu, ypaprs, &
+                ycdragm(:knon), yt(:knon, :), yu(:knon, :), ypaprs(:knon, :), &
-                ypplay, ydelp, y_d_u, y_flux_u(:knon))
+                ypplay(:knon, :), ydelp(:knon, :), y_d_u(:knon, :), &
-           CALL clvent(knon, dtime, yu1, yv1, coefm(:knon, :), yt, yv, ypaprs, &
+                y_flux_u(:knon))
-                ypplay, ydelp, y_d_v, y_flux_v(:knon))
+           CALL clvent(dtime, yu(:knon, 1), yv(:knon, 1), ycoefm(:knon, :), &
+                ycdragm(:knon), yt(:knon, :), yv(:knon, :), ypaprs(:knon, :), &
+                ypplay(:knon, :), ydelp(:knon, :), y_d_v(:knon, :), &
+                y_flux_v(:knon))
            ! calculer la diffusion de "q" et de "h"
            CALL clqh(dtime, julien, firstcal, nsrf, ni(:knon), &
-                ytsoil(:knon, :), yqsol, mu0, yrugos, yrugoro, yu1, yv1, &
+                ytsoil(:knon, :), yqsol(:knon), mu0, yrugos, yrugoro, &
-                coefh(:knon, :), yt, yq, yts(:knon), ypaprs, ypplay, ydelp, &
+                yu(:knon, 1), yv(:knon, 1), ycoefh(:knon, :), ycdragh(:knon), &
-                yrads, yalb(:knon), snow(:knon), yqsurf, yrain_f, ysnow_f, &
+                yt, yq, yts(:knon), ypaprs, ypplay, ydelp, yrads(:knon), &
-                yfder, yfluxlat(:knon), pctsrf_new_sic, yagesno(:knon), y_d_t, &
+                yalb(:knon), snow(:knon), yqsurf, yrain_f, ysnow_f, &
-                y_d_q, y_d_ts(:knon), yz0_new, y_flux_t(:knon), &
+                yfluxlat(:knon), pctsrf_new_sic, yagesno(:knon), y_d_t, y_d_q, &
-                y_flux_q(:knon), y_dflux_t, y_dflux_q, y_fqcalving, y_ffonte, &
+                y_d_ts(:knon), yz0_new, y_flux_t(:knon), y_flux_q(:knon), &
+                y_dflux_t(:knon), y_dflux_q(:knon), y_fqcalving, y_ffonte, &
                 y_run_off_lic_0)
            ! calculer la longueur de rugosite sur ocean
            yrugm = 0.
            IF (nsrf == is_oce) THEN
               DO j = 1, knon
-                 yrugm(j) = 0.018*coefm(j, 1)*(yu1(j)**2+yv1(j)**2)/rg + &
+                 yrugm(j) = 0.018 * ycdragm(j) * (yu(j, 1)**2 + yv(j, 1)**2) &
-.11*14E-6/sqrt(coefm(j, 1)*(yu1(j)**2+yv1(j)**2))
+                      / rg + 0.11 * 14E-6 &
+                      / sqrt(ycdragm(j) * (yu(j, 1)**2 + yv(j, 1)**2))
                  yrugm(j) = max(1.5E-05, yrugm(j))
               END DO
            END IF
            DO j = 1, knon
-              y_dflux_t(j) = y_dflux_t(j)*ypct(j)
+              y_dflux_t(j) = y_dflux_t(j) * ypct(j)
-              y_dflux_q(j) = y_dflux_q(j)*ypct(j)
+              y_dflux_q(j) = y_dflux_q(j) * ypct(j)
-              yu1(j) = yu1(j)*ypct(j)
-              yv1(j) = yv1(j)*ypct(j)
            END DO
            DO k = 1, klev
               DO j = 1, knon
                  i = ni(j)
-                 coefh(j, k) = coefh(j, k)*ypct(j)
+                 y_d_t(j, k) = y_d_t(j, k) * ypct(j)
-                 coefm(j, k) = coefm(j, k)*ypct(j)
+                 y_d_q(j, k) = y_d_q(j, k) * ypct(j)
-                 y_d_t(j, k) = y_d_t(j, k)*ypct(j)
+                 y_d_u(j, k) = y_d_u(j, k) * ypct(j)
-                 y_d_q(j, k) = y_d_q(j, k)*ypct(j)
+                 y_d_v(j, k) = y_d_v(j, k) * ypct(j)
-                 y_d_u(j, k) = y_d_u(j, k)*ypct(j)
-                 y_d_v(j, k) = y_d_v(j, k)*ypct(j)
               END DO
            END DO
-Line 481 
 contains
+Line 451 
 contains
            falbe(:, nsrf) = 0.
            fsnow(:, nsrf) = 0.
            qsurf(:, nsrf) = 0.
-           rugos(:, nsrf) = 0.
+           frugs(:, nsrf) = 0.
            DO j = 1, knon
               i = ni(j)
               d_ts(i, nsrf) = y_d_ts(j)
               falbe(i, nsrf) = yalb(j)
               fsnow(i, nsrf) = snow(j)
               qsurf(i, nsrf) = yqsurf(j)
-              rugos(i, nsrf) = yz0_new(j)
+              frugs(i, nsrf) = yz0_new(j)
               fluxlat(i, nsrf) = yfluxlat(j)
               IF (nsrf == is_oce) THEN
                  rugmer(i) = yrugm(j)
-                 rugos(i, nsrf) = yrugm(j)
+                 frugs(i, nsrf) = yrugm(j)
               END IF
               agesno(i, nsrf) = yagesno(j)
               fqcalving(i, nsrf) = y_fqcalving(j)
               ffonte(i, nsrf) = y_ffonte(j)
-              cdragh(i) = cdragh(i) + coefh(j, 1)
+              cdragh(i) = cdragh(i) + ycdragh(j) * ypct(j)
-              cdragm(i) = cdragm(i) + coefm(j, 1)
+              cdragm(i) = cdragm(i) + ycdragm(j) * ypct(j)
               dflux_t(i) = dflux_t(i) + y_dflux_t(j)
               dflux_q(i) = dflux_q(i) + y_dflux_q(j)
-              zu1(i) = zu1(i) + yu1(j)
-              zv1(i) = zv1(i) + yv1(j)
            END DO
            IF (nsrf == is_ter) THEN
               qsol(ni(:knon)) = yqsol(:knon)
-Line 523 
 contains
+Line 491 
 contains
                  d_q(i, k) = d_q(i, k) + y_d_q(j, k)
                  d_u(i, k) = d_u(i, k) + y_d_u(j, k)
                  d_v(i, k) = d_v(i, k) + y_d_v(j, k)
-                 ycoefh(i, k) = ycoefh(i, k) + coefh(j, k)
               END DO
            END DO
+           forall (k = 2:klev) coefh(ni(:knon), k) &
+                = coefh(ni(:knon), k) + ycoefh(:knon, k) * ypct(:knon)
            ! diagnostic t, q a 2m et u, v a 10m
            DO j = 1, knon
               i = ni(j)
-              uzon(j) = yu(j, 1) + y_d_u(j, 1)
+              u1(j) = yu(j, 1) + y_d_u(j, 1)
-              vmer(j) = yv(j, 1) + y_d_v(j, 1)
+              v1(j) = yv(j, 1) + y_d_v(j, 1)
               tair1(j) = yt(j, 1) + y_d_t(j, 1)
               qair1(j) = yq(j, 1) + y_d_q(j, 1)
-              zgeo1(j) = rd*tair1(j)/(0.5*(ypaprs(j, 1)+ypplay(j, &
+              zgeo1(j) = rd * tair1(j) / (0.5 * (ypaprs(j, 1) + ypplay(j, &
-)))*(ypaprs(j, 1)-ypplay(j, 1))
+))) * (ypaprs(j, 1)-ypplay(j, 1))
               tairsol(j) = yts(j) + y_d_ts(j)
               rugo1(j) = yrugos(j)
               IF (nsrf == is_oce) THEN
-                 rugo1(j) = rugos(i, nsrf)
+                 rugo1(j) = frugs(i, nsrf)
               END IF
               psfce(j) = ypaprs(j, 1)
               patm(j) = ypplay(j, 1)
-Line 548 
 contains
+Line 518 
 contains
               qairsol(j) = yqsurf(j)
            END DO
-           CALL stdlevvar(klon, knon, nsrf, zxli, uzon, vmer, tair1, qair1, &
+           CALL stdlevvar(klon, knon, nsrf, u1(:knon), v1(:knon), tair1(:knon), &
-                zgeo1, tairsol, qairsol, rugo1, psfce, patm, yt2m, yq2m, &
+                qair1, zgeo1, tairsol, qairsol, rugo1, psfce, patm, yt2m, &
-                yt10m, yq10m, yu10m, yustar)
+                yq2m, yt10m, yq10m, wind10m(:knon), ustar(:knon))
            DO j = 1, knon
               i = ni(j)
               t2m(i, nsrf) = yt2m(j)
               q2m(i, nsrf) = yq2m(j)
-              ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman
+              u10m_srf(i, nsrf) = (wind10m(j) * u1(j)) &
-              u10m(i, nsrf) = (yu10m(j)*uzon(j))/sqrt(uzon(j)**2+vmer(j)**2)
+                   / sqrt(u1(j)**2 + v1(j)**2)
-              v10m(i, nsrf) = (yu10m(j)*vmer(j))/sqrt(uzon(j)**2+vmer(j)**2)
+              v10m_srf(i, nsrf) = (wind10m(j) * v1(j)) &
+                   / sqrt(u1(j)**2 + v1(j)**2)
            END DO
-           CALL hbtm(ypaprs, ypplay, yt2m, yq2m, yustar, y_flux_t(:knon), &
+           CALL hbtm(ypaprs, ypplay, yt2m, yq2m, ustar(:knon), y_flux_t(:knon), &
                 y_flux_q(:knon), yu, yv, yt, yq, ypblh(:knon), ycapcl, &
                 yoliqcl, ycteicl, ypblt, ytherm, ytrmb1, ytrmb2, ytrmb3, ylcl)
-Line 592 
 contains
+Line 563 
 contains
      END DO loop_surface
      ! On utilise les nouvelles surfaces
-     rugos(:, is_oce) = rugmer
+     frugs(:, is_oce) = rugmer
      pctsrf(:, is_oce) = pctsrf_new_oce
      pctsrf(:, is_sic) = pctsrf_new_sic

 Legend:



Removed from v.221
 


changed lines


 
Added in v.248
 Legend:



Removed from v.221
 


changed lines


 
Added in v.248
-Removed from v.221
+Added in v.248

	ViewVC Help
Powered by ViewVC 1.1.21