/[lmdze]/trunk/phylmd/Interface_surf/pbl_surface.f
ViewVC logotype

Diff of /trunk/phylmd/Interface_surf/pbl_surface.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 186 by guez, Mon Mar 21 15:36:26 2016 UTC revision 229 by guez, Mon Nov 6 17:20:45 2017 UTC
# Line 4  module clmain_m Line 4  module clmain_m
4    
5  contains  contains
6    
7    SUBROUTINE clmain(dtime, itap, pctsrf, pctsrf_new, t, q, u, v, jour, rmu0, &    SUBROUTINE clmain(dtime, pctsrf, t, q, u, v, julien, mu0, ftsol, cdmmax, &
8         ts, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, qsol, &         cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, qsol, paprs, pplay, fsnow, &
9         paprs, pplay, snow, qsurf, evap, falbe, fluxlat, rain_fall, snow_f, &         qsurf, evap, falbe, fluxlat, rain_fall, snow_f, fsolsw, fsollw, frugs, &
10         solsw, sollw, fder, rlat, rugos, debut, agesno, rugoro, d_t, d_q, d_u, &         agesno, rugoro, d_t, d_q, d_u, d_v, d_ts, flux_t, flux_q, flux_u, &
11         d_v, d_ts, flux_t, flux_q, flux_u, flux_v, cdragh, cdragm, q2, &         flux_v, cdragh, cdragm, q2, dflux_t, dflux_q, ycoefh, t2m, q2m, &
12         dflux_t, dflux_q, ycoefh, zu1, zv1, t2m, q2m, u10m, v10m, pblh, capcl, &         u10m_srf, v10m_srf, pblh, capcl, oliqcl, cteicl, pblt, therm, trmb1, &
13         oliqcl, cteicl, pblt, therm, trmb1, trmb2, trmb3, plcl, fqcalving, &         trmb2, trmb3, plcl, fqcalving, ffonte, run_off_lic_0)
        ffonte, run_off_lic_0)  
14    
15      ! From phylmd/clmain.F, version 1.6, 2005/11/16 14:47:19      ! From phylmd/clmain.F, version 1.6, 2005/11/16 14:47:19
16      ! Author: Z. X. Li (LMD/CNRS), date: 1993/08/18      ! Author: Z. X. Li (LMD/CNRS), date: 1993/08/18
# Line 22  contains Line 21  contains
21      ! ne tient pas compte de la diff\'erentiation des sous-fractions      ! ne tient pas compte de la diff\'erentiation des sous-fractions
22      ! de sol.      ! de sol.
23    
     ! Pour pouvoir extraire les coefficients d'\'echanges et le vent  
     ! 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.  
   
24      use clqh_m, only: clqh      use clqh_m, only: clqh
25      use clvent_m, only: clvent      use clvent_m, only: clvent
26      use coefkz_m, only: coefkz      use coefkz_m, only: coefkz
27      use coefkzmin_m, only: coefkzmin      use coefkzmin_m, only: coefkzmin
28      USE conf_gcm_m, ONLY: prt_level      USE conf_gcm_m, ONLY: lmt_pas
29      USE conf_phys_m, ONLY: iflag_pbl      USE conf_phys_m, ONLY: iflag_pbl
30      USE dimphy, ONLY: klev, klon, zmasq      USE dimphy, ONLY: klev, klon, zmasq
31      USE dimsoil, ONLY: nsoilmx      USE dimsoil, ONLY: nsoilmx
32      use hbtm_m, only: hbtm      use hbtm_m, only: hbtm
33      USE indicesol, ONLY: epsfra, is_lic, is_oce, is_sic, is_ter, nbsrf      USE indicesol, ONLY: epsfra, is_lic, is_oce, is_sic, is_ter, nbsrf
34        USE interfoce_lim_m, ONLY: interfoce_lim
35      use stdlevvar_m, only: stdlevvar      use stdlevvar_m, only: stdlevvar
36      USE suphec_m, ONLY: rd, rg, rkappa      USE suphec_m, ONLY: rd, rg, rkappa
37        use time_phylmdz, only: itap
38      use ustarhb_m, only: ustarhb      use ustarhb_m, only: ustarhb
     use vdif_kcay_m, only: vdif_kcay  
39      use yamada4_m, only: yamada4      use yamada4_m, only: yamada4
40    
41      REAL, INTENT(IN):: dtime ! interval du temps (secondes)      REAL, INTENT(IN):: dtime ! interval du temps (secondes)
     INTEGER, INTENT(IN):: itap ! numero du pas de temps  
     REAL, INTENT(inout):: pctsrf(klon, nbsrf)  
42    
43      ! la nouvelle repartition des surfaces sortie de l'interface      REAL, INTENT(inout):: pctsrf(klon, nbsrf)
44      REAL, INTENT(out):: pctsrf_new(klon, nbsrf)      ! tableau des pourcentages de surface de chaque maille
45    
46      REAL, INTENT(IN):: t(klon, klev) ! temperature (K)      REAL, INTENT(IN):: t(klon, klev) ! temperature (K)
47      REAL, INTENT(IN):: q(klon, klev) ! vapeur d'eau (kg/kg)      REAL, INTENT(IN):: q(klon, klev) ! vapeur d'eau (kg / kg)
48      REAL, INTENT(IN):: u(klon, klev), v(klon, klev) ! vitesse      REAL, INTENT(IN):: u(klon, klev), v(klon, klev) ! vitesse
49      INTEGER, INTENT(IN):: jour ! jour de l'annee en cours      INTEGER, INTENT(IN):: julien ! jour de l'annee en cours
50      REAL, intent(in):: rmu0(klon) ! cosinus de l'angle solaire zenithal          REAL, intent(in):: mu0(klon) ! cosinus de l'angle solaire zenithal    
51      REAL, INTENT(IN):: ts(klon, nbsrf) ! temperature du sol (en Kelvin)      REAL, INTENT(IN):: ftsol(:, :) ! (klon, nbsrf) temp\'erature du sol (en K)
52      REAL, INTENT(IN):: cdmmax, cdhmax ! seuils cdrm, cdrh      REAL, INTENT(IN):: cdmmax, cdhmax ! seuils cdrm, cdrh
53      REAL, INTENT(IN):: ksta, ksta_ter      REAL, INTENT(IN):: ksta, ksta_ter
54      LOGICAL, INTENT(IN):: ok_kzmin      LOGICAL, INTENT(IN):: ok_kzmin
# Line 63  contains Line 56  contains
56      REAL, INTENT(inout):: ftsoil(klon, nsoilmx, nbsrf)      REAL, INTENT(inout):: ftsoil(klon, nsoilmx, nbsrf)
57      ! soil temperature of surface fraction      ! soil temperature of surface fraction
58    
59      REAL, INTENT(inout):: qsol(klon)      REAL, INTENT(inout):: qsol(:) ! (klon)
60      ! column-density of water in soil, in kg m-2      ! column-density of water in soil, in kg m-2
61    
62      REAL, INTENT(IN):: paprs(klon, klev+1) ! pression a intercouche (Pa)      REAL, INTENT(IN):: paprs(klon, klev + 1) ! pression a intercouche (Pa)
63      REAL, INTENT(IN):: pplay(klon, klev) ! pression au milieu de couche (Pa)      REAL, INTENT(IN):: pplay(klon, klev) ! pression au milieu de couche (Pa)
64      REAL snow(klon, nbsrf)      REAL, INTENT(inout):: fsnow(:, :) ! (klon, nbsrf) \'epaisseur neigeuse
65      REAL qsurf(klon, nbsrf)      REAL qsurf(klon, nbsrf)
66      REAL evap(klon, nbsrf)      REAL evap(klon, nbsrf)
67      REAL, intent(inout):: falbe(klon, nbsrf)      REAL, intent(inout):: falbe(klon, nbsrf)
68        REAL, intent(out):: fluxlat(:, :) ! (klon, nbsrf)
     REAL fluxlat(klon, nbsrf)  
69    
70      REAL, intent(in):: rain_fall(klon)      REAL, intent(in):: rain_fall(klon)
71      ! liquid water mass flux (kg/m2/s), positive down      ! liquid water mass flux (kg / m2 / s), positive down
72    
73      REAL, intent(in):: snow_f(klon)      REAL, intent(in):: snow_f(klon)
74      ! 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):: fder(klon)  
     REAL, INTENT(IN):: rlat(klon) ! latitude en degr\'es  
   
     REAL rugos(klon, nbsrf)  
     ! rugos----input-R- longeur de rugosite (en m)  
75    
76      LOGICAL, INTENT(IN):: debut      REAL, INTENT(IN):: fsolsw(klon, nbsrf), fsollw(klon, nbsrf)
77        REAL, intent(inout):: frugs(klon, nbsrf) ! longueur de rugosit\'e (en m)
78      real agesno(klon, nbsrf)      real agesno(klon, nbsrf)
79      REAL, INTENT(IN):: rugoro(klon)      REAL, INTENT(IN):: rugoro(klon)
80    
# Line 99  contains Line 85  contains
85      REAL, intent(out):: d_u(klon, klev), d_v(klon, klev)      REAL, intent(out):: d_u(klon, klev), d_v(klon, klev)
86      ! changement pour "u" et "v"      ! changement pour "u" et "v"
87    
88      REAL, intent(out):: d_ts(klon, nbsrf) ! le changement pour "ts"      REAL, intent(out):: d_ts(:, :) ! (klon, nbsrf) variation of ftsol
89    
90      REAL flux_t(klon, klev, nbsrf), flux_q(klon, klev, nbsrf)      REAL, intent(out):: flux_t(klon, nbsrf)
91      ! flux_t---output-R- flux de chaleur sensible (CpT) J/m**2/s (W/m**2)      ! flux de chaleur sensible (Cp T) (W / m2) (orientation positive vers
92      !                    (orientation positive vers le bas)      ! le bas) à la surface
93      ! flux_q---output-R- flux de vapeur d'eau (kg/m**2/s)  
94        REAL, intent(out):: flux_q(klon, nbsrf)
95      REAL flux_u(klon, klev, nbsrf), flux_v(klon, klev, nbsrf)      ! flux de vapeur d'eau (kg / m2 / s) à la surface
96      ! flux_u---output-R- tension du vent X: (kg m/s)/(m**2 s) ou Pascal  
97      ! flux_v---output-R- tension du vent Y: (kg m/s)/(m**2 s) ou Pascal      REAL, intent(out):: flux_u(klon, nbsrf), flux_v(klon, nbsrf)
98        ! tension du vent (flux turbulent de vent) à la surface, en Pa
99    
100      REAL, INTENT(out):: cdragh(klon), cdragm(klon)      REAL, INTENT(out):: cdragh(klon), cdragm(klon)
101      real q2(klon, klev+1, nbsrf)      real q2(klon, klev + 1, nbsrf)
102    
103      REAL, INTENT(out):: dflux_t(klon), dflux_q(klon)      REAL, INTENT(out):: dflux_t(klon), dflux_q(klon)
104      ! dflux_t derive du flux sensible      ! dflux_t derive du flux sensible
105      ! dflux_q derive du flux latent      ! dflux_q derive du flux latent
106      !IM "slab" ocean      ! IM "slab" ocean
107    
108      REAL, intent(out):: ycoefh(klon, klev)      REAL, intent(out):: ycoefh(klon, klev)
109      REAL, intent(out):: zu1(klon)      ! Pour pouvoir extraire les coefficients d'\'echange, le champ
110      REAL zv1(klon)      ! "ycoefh" a \'et\'e cr\'e\'e. Nous avons moyenn\'e les valeurs de
111      REAL t2m(klon, nbsrf), q2m(klon, nbsrf)      ! ce champ sur les quatre sous-surfaces du mod\`ele.
112      REAL u10m(klon, nbsrf), v10m(klon, nbsrf)  
113        REAL, INTENT(inout):: t2m(klon, nbsrf), q2m(klon, nbsrf)
114      ! Ionela Musat cf. Anne Mathieu : pbl, hbtm (Comme les autres  
115      ! diagnostics on cumule dans physiq ce qui permet de sortir les      REAL, INTENT(inout):: u10m_srf(:, :), v10m_srf(:, :) ! (klon, nbsrf)
116      ! grandeurs par sous-surface)      ! composantes du vent \`a 10m sans spirale d'Ekman
117      REAL pblh(klon, nbsrf)  
118      ! pblh------- HCL      ! Ionela Musat. Cf. Anne Mathieu : planetary boundary layer, hbtm.
119        ! Comme les autres diagnostics on cumule dans physiq ce qui permet
120        ! de sortir les grandeurs par sous-surface.
121        REAL pblh(klon, nbsrf) ! height of planetary boundary layer
122      REAL capcl(klon, nbsrf)      REAL capcl(klon, nbsrf)
123      REAL oliqcl(klon, nbsrf)      REAL oliqcl(klon, nbsrf)
124      REAL cteicl(klon, nbsrf)      REAL cteicl(klon, nbsrf)
125      REAL pblt(klon, nbsrf)      REAL, INTENT(inout):: pblt(klon, nbsrf) ! T au nveau HCL
     ! pblT------- T au nveau HCL  
126      REAL therm(klon, nbsrf)      REAL therm(klon, nbsrf)
127      REAL trmb1(klon, nbsrf)      REAL trmb1(klon, nbsrf)
128      ! trmb1-------deep_cape      ! trmb1-------deep_cape
# Line 145  contains Line 134  contains
134      REAL fqcalving(klon, nbsrf), ffonte(klon, nbsrf)      REAL fqcalving(klon, nbsrf), ffonte(klon, nbsrf)
135      ! ffonte----Flux thermique utilise pour fondre la neige      ! ffonte----Flux thermique utilise pour fondre la neige
136      ! fqcalving-Flux d'eau "perdue" par la surface et necessaire pour limiter la      ! fqcalving-Flux d'eau "perdue" par la surface et necessaire pour limiter la
137      !           hauteur de neige, en kg/m2/s      !           hauteur de neige, en kg / m2 / s
138      REAL run_off_lic_0(klon)      REAL run_off_lic_0(klon)
139    
140      ! Local:      ! Local:
141    
142        LOGICAL:: firstcal = .true.
143    
144        ! la nouvelle repartition des surfaces sortie de l'interface
145        REAL, save:: pctsrf_new_oce(klon)
146        REAL, save:: pctsrf_new_sic(klon)
147    
148      REAL y_fqcalving(klon), y_ffonte(klon)      REAL y_fqcalving(klon), y_ffonte(klon)
149      real y_run_off_lic_0(klon)      real y_run_off_lic_0(klon)
   
150      REAL rugmer(klon)      REAL rugmer(klon)
   
151      REAL ytsoil(klon, nsoilmx)      REAL ytsoil(klon, nsoilmx)
   
152      REAL yts(klon), yrugos(klon), ypct(klon), yz0_new(klon)      REAL yts(klon), yrugos(klon), ypct(klon), yz0_new(klon)
153      REAL yalb(klon)      REAL yalb(klon)
154      REAL yu1(klon), yv1(klon)      REAL snow(klon), yqsurf(klon), yagesno(klon)
155      ! on rajoute en output yu1 et yv1 qui sont les vents dans      real yqsol(klon) ! column-density of water in soil, in kg m-2
156      ! la premiere couche      REAL yrain_f(klon) ! liquid water mass flux (kg / m2 / s), positive down
157      REAL ysnow(klon), yqsurf(klon), yagesno(klon)      REAL ysnow_f(klon) ! solid water mass flux (kg / m2 / s), positive down
   
     real yqsol(klon)  
     ! column-density of water in soil, in kg m-2  
   
     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)  
158      REAL yrugm(klon), yrads(klon), yrugoro(klon)      REAL yrugm(klon), yrads(klon), yrugoro(klon)
   
159      REAL yfluxlat(klon)      REAL yfluxlat(klon)
   
160      REAL y_d_ts(klon)      REAL y_d_ts(klon)
161      REAL y_d_t(klon, klev), y_d_q(klon, klev)      REAL y_d_t(klon, klev), y_d_q(klon, klev)
162      REAL y_d_u(klon, klev), y_d_v(klon, klev)      REAL y_d_u(klon, klev), y_d_v(klon, klev)
163      REAL y_flux_t(klon, klev), y_flux_q(klon, klev)      REAL y_flux_t(klon), y_flux_q(klon)
164      REAL y_flux_u(klon, klev), y_flux_v(klon, klev)      REAL y_flux_u(klon), y_flux_v(klon)
165      REAL y_dflux_t(klon), y_dflux_q(klon)      REAL y_dflux_t(klon), y_dflux_q(klon)
166      REAL coefh(klon, klev), coefm(klon, klev)      REAL coefh(klon, klev), coefm(klon, klev)
167      REAL yu(klon, klev), yv(klon, klev)      REAL yu(klon, klev), yv(klon, klev)
168      REAL yt(klon, klev), yq(klon, klev)      REAL yt(klon, klev), yq(klon, klev)
169      REAL ypaprs(klon, klev+1), ypplay(klon, klev), ydelp(klon, klev)      REAL ypaprs(klon, klev + 1), ypplay(klon, klev), ydelp(klon, klev)
   
170      REAL ycoefm0(klon, klev), ycoefh0(klon, klev)      REAL ycoefm0(klon, klev), ycoefh0(klon, klev)
171        REAL yzlay(klon, klev), zlev(klon, klev + 1), yteta(klon, klev)
172      REAL yzlay(klon, klev), yzlev(klon, klev+1), yteta(klon, klev)      REAL ykmm(klon, klev + 1), ykmn(klon, klev + 1)
173      REAL ykmm(klon, klev+1), ykmn(klon, klev+1)      REAL ykmq(klon, klev + 1)
174      REAL ykmq(klon, klev+1)      REAL yq2(klon, klev + 1)
     REAL yq2(klon, klev+1)  
     REAL q2diag(klon, klev+1)  
   
     REAL u1lay(klon), v1lay(klon)  
175      REAL delp(klon, klev)      REAL delp(klon, klev)
176      INTEGER i, k, nsrf      INTEGER i, k, nsrf
   
177      INTEGER ni(klon), knon, j      INTEGER ni(klon), knon, j
178    
179      REAL pctsrf_pot(klon, nbsrf)      REAL pctsrf_pot(klon, nbsrf)
180      ! "pourcentage potentiel" pour tenir compte des \'eventuelles      ! "pourcentage potentiel" pour tenir compte des \'eventuelles
181      ! apparitions ou disparitions de la glace de mer      ! apparitions ou disparitions de la glace de mer
182    
183      REAL zx_alf1, zx_alf2 !valeur ambiante par extrapola.      REAL yt2m(klon), yq2m(klon), wind10m(klon)
184        REAL ustar(klon)
     REAL yt2m(klon), yq2m(klon), yu10m(klon)  
     REAL yustar(klon)  
185    
186      REAL yt10m(klon), yq10m(klon)      REAL yt10m(klon), yq10m(klon)
187      REAL ypblh(klon)      REAL ypblh(klon)
# Line 223  contains Line 194  contains
194      REAL ytrmb1(klon)      REAL ytrmb1(klon)
195      REAL ytrmb2(klon)      REAL ytrmb2(klon)
196      REAL ytrmb3(klon)      REAL ytrmb3(klon)
197      REAL uzon(klon), vmer(klon)      REAL u1(klon), v1(klon)
198      REAL tair1(klon), qair1(klon), tairsol(klon)      REAL tair1(klon), qair1(klon), tairsol(klon)
199      REAL psfce(klon), patm(klon)      REAL psfce(klon), patm(klon)
200    
201      REAL qairsol(klon), zgeo1(klon)      REAL qairsol(klon), zgeo1(klon)
202      REAL rugo1(klon)      REAL rugo1(klon)
203    
     ! utiliser un jeu de fonctions simples                
     LOGICAL zxli  
     PARAMETER (zxli=.FALSE.)  
   
204      !------------------------------------------------------------      !------------------------------------------------------------
205    
206      ytherm = 0.      ytherm = 0.
207    
208      DO k = 1, klev ! epaisseur de couche      DO k = 1, klev ! epaisseur de couche
209         DO i = 1, klon         DO i = 1, klon
210            delp(i, k) = paprs(i, k) - paprs(i, k+1)            delp(i, k) = paprs(i, k) - paprs(i, k + 1)
211         END DO         END DO
212      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  
213    
214      ! Initialization:      ! Initialization:
215      rugmer = 0.      rugmer = 0.
# Line 256  contains Line 217  contains
217      cdragm = 0.      cdragm = 0.
218      dflux_t = 0.      dflux_t = 0.
219      dflux_q = 0.      dflux_q = 0.
     zu1 = 0.  
     zv1 = 0.  
220      ypct = 0.      ypct = 0.
     yts = 0.  
     ysnow = 0.  
221      yqsurf = 0.      yqsurf = 0.
222      yrain_f = 0.      yrain_f = 0.
223      ysnow_f = 0.      ysnow_f = 0.
     yfder = 0.  
224      yrugos = 0.      yrugos = 0.
     yu1 = 0.  
     yv1 = 0.  
     yrads = 0.  
225      ypaprs = 0.      ypaprs = 0.
226      ypplay = 0.      ypplay = 0.
227      ydelp = 0.      ydelp = 0.
# Line 276  contains Line 229  contains
229      yv = 0.      yv = 0.
230      yt = 0.      yt = 0.
231      yq = 0.      yq = 0.
     pctsrf_new = 0.  
     y_flux_u = 0.  
     y_flux_v = 0.  
232      y_dflux_t = 0.      y_dflux_t = 0.
233      y_dflux_q = 0.      y_dflux_q = 0.
     ytsoil = 999999.  
234      yrugoro = 0.      yrugoro = 0.
235      d_ts = 0.      d_ts = 0.
     yfluxlat = 0.  
236      flux_t = 0.      flux_t = 0.
237      flux_q = 0.      flux_q = 0.
238      flux_u = 0.      flux_u = 0.
239      flux_v = 0.      flux_v = 0.
240        fluxlat = 0.
241      d_t = 0.      d_t = 0.
242      d_q = 0.      d_q = 0.
243      d_u = 0.      d_u = 0.
# Line 299  contains Line 248  contains
248      ! peut avoir potentiellement de la glace sur tout le domaine oc\'eanique      ! peut avoir potentiellement de la glace sur tout le domaine oc\'eanique
249      ! (\`a affiner)      ! (\`a affiner)
250    
251      pctsrf_pot = pctsrf      pctsrf_pot(:, is_ter) = pctsrf(:, is_ter)
252        pctsrf_pot(:, is_lic) = pctsrf(:, is_lic)
253      pctsrf_pot(:, is_oce) = 1. - zmasq      pctsrf_pot(:, is_oce) = 1. - zmasq
254      pctsrf_pot(:, is_sic) = 1. - zmasq      pctsrf_pot(:, is_sic) = 1. - zmasq
255    
256        ! Tester si c'est le moment de lire le fichier:
257        if (mod(itap - 1, lmt_pas) == 0) then
258           CALL interfoce_lim(julien, pctsrf_new_oce, pctsrf_new_sic)
259        endif
260    
261      ! Boucler sur toutes les sous-fractions du sol:      ! Boucler sur toutes les sous-fractions du sol:
262    
263      loop_surface: DO nsrf = 1, nbsrf      loop_surface: DO nsrf = 1, nbsrf
# Line 322  contains Line 277  contains
277            DO j = 1, knon            DO j = 1, knon
278               i = ni(j)               i = ni(j)
279               ypct(j) = pctsrf(i, nsrf)               ypct(j) = pctsrf(i, nsrf)
280               yts(j) = ts(i, nsrf)               yts(j) = ftsol(i, nsrf)
281               ysnow(j) = snow(i, nsrf)               snow(j) = fsnow(i, nsrf)
282               yqsurf(j) = qsurf(i, nsrf)               yqsurf(j) = qsurf(i, nsrf)
283               yalb(j) = falbe(i, nsrf)               yalb(j) = falbe(i, nsrf)
284               yrain_f(j) = rain_fall(i)               yrain_f(j) = rain_fall(i)
285               ysnow_f(j) = snow_f(i)               ysnow_f(j) = snow_f(i)
286               yagesno(j) = agesno(i, nsrf)               yagesno(j) = agesno(i, nsrf)
287               yfder(j) = fder(i)               yrugos(j) = frugs(i, nsrf)
              yrugos(j) = rugos(i, nsrf)  
288               yrugoro(j) = rugoro(i)               yrugoro(j) = rugoro(i)
289               yu1(j) = u1lay(i)               yrads(j) = fsolsw(i, nsrf) + fsollw(i, nsrf)
290               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)  
291               y_run_off_lic_0(j) = run_off_lic_0(i)               y_run_off_lic_0(j) = run_off_lic_0(i)
292            END DO            END DO
293    
294            ! For continent, copy soil water content            ! For continent, copy soil water content
295            IF (nsrf == is_ter) THEN            IF (nsrf == is_ter) yqsol(:knon) = qsol(ni(:knon))
              yqsol(:knon) = qsol(ni(:knon))  
           ELSE  
              yqsol = 0.  
           END IF  
296    
297            DO k = 1, nsoilmx            ytsoil(:knon, :) = ftsoil(ni(:knon), :, nsrf)
              DO j = 1, knon  
                 i = ni(j)  
                 ytsoil(j, k) = ftsoil(i, k, nsrf)  
              END DO  
           END DO  
298    
299            DO k = 1, klev            DO k = 1, klev
300               DO j = 1, knon               DO j = 1, knon
# Line 367  contains Line 310  contains
310            END DO            END DO
311    
312            ! calculer Cdrag et les coefficients d'echange            ! calculer Cdrag et les coefficients d'echange
313            CALL coefkz(nsrf, knon, ypaprs, ypplay, ksta, ksta_ter, yts, yrugos, &            CALL coefkz(nsrf, ypaprs, ypplay, ksta, ksta_ter, yts(:knon), &
314                 yu, yv, yt, yq, yqsurf, coefm(:knon, :), coefh(:knon, :))                 yrugos, yu, yv, yt, yq, yqsurf(:knon), coefm(:knon, :), &
315                   coefh(:knon, :))
316    
317            IF (iflag_pbl == 1) THEN            IF (iflag_pbl == 1) THEN
318               CALL coefkz2(nsrf, knon, ypaprs, ypplay, yt, ycoefm0, ycoefh0)               CALL coefkz2(nsrf, knon, ypaprs, ypplay, yt, ycoefm0, ycoefh0)
319               coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))               coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))
# Line 389  contains Line 334  contains
334               coefh(:knon, :) = max(coefh(:knon, :), ycoefh0(:knon, :))               coefh(:knon, :) = max(coefh(:knon, :), ycoefh0(:knon, :))
335            END IF            END IF
336    
337            IF (iflag_pbl >= 3) THEN            IF (iflag_pbl >= 6) THEN
338               ! Mellor et Yamada adapt\'e \`a Mars, Richard Fournier et               ! Mellor et Yamada adapt\'e \`a Mars, Richard Fournier et
339               ! Fr\'ed\'eric Hourdin               ! Fr\'ed\'eric Hourdin
340               yzlay(:knon, 1) = rd * yt(:knon, 1) / (0.5 * (ypaprs(:knon, 1) &               yzlay(:knon, 1) = rd * yt(:knon, 1) / (0.5 * (ypaprs(:knon, 1) &
341                    + ypplay(:knon, 1))) &                    + ypplay(:knon, 1))) &
342                    * (ypaprs(:knon, 1) - ypplay(:knon, 1)) / rg                    * (ypaprs(:knon, 1) - ypplay(:knon, 1)) / rg
343    
344               DO k = 2, klev               DO k = 2, klev
345                  yzlay(1:knon, k) = yzlay(1:knon, k-1) &                  yzlay(:knon, k) = yzlay(:knon, k-1) &
346                       + rd * 0.5 * (yt(1:knon, k-1) + yt(1:knon, k)) &                       + rd * 0.5 * (yt(1:knon, k-1) + yt(1:knon, k)) &
347                       / ypaprs(1:knon, k) &                       / ypaprs(1:knon, k) &
348                       * (ypplay(1:knon, k-1) - ypplay(1:knon, k)) / rg                       * (ypplay(1:knon, k-1) - ypplay(1:knon, k)) / rg
349               END DO               END DO
350    
351               DO k = 1, klev               DO k = 1, klev
352                  yteta(1:knon, k) = yt(1:knon, k)*(ypaprs(1:knon, 1) &                  yteta(1:knon, k) = yt(1:knon, k) * (ypaprs(1:knon, 1) &
353                       / ypplay(1:knon, k))**rkappa * (1.+0.61*yq(1:knon, k))                       / ypplay(1:knon, k))**rkappa * (1. + 0.61 * yq(1:knon, k))
354               END DO               END DO
355               yzlev(1:knon, 1) = 0.  
356               yzlev(:knon, klev+1) = 2. * yzlay(:knon, klev) &               zlev(:knon, 1) = 0.
357                 zlev(:knon, klev + 1) = 2. * yzlay(:knon, klev) &
358                    - yzlay(:knon, klev - 1)                    - yzlay(:knon, klev - 1)
359    
360               DO k = 2, klev               DO k = 2, klev
361                  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))
362               END DO               END DO
363    
364               DO k = 1, klev + 1               DO k = 1, klev + 1
365                  DO j = 1, knon                  DO j = 1, knon
366                     i = ni(j)                     i = ni(j)
# Line 418  contains Line 368  contains
368                  END DO                  END DO
369               END DO               END DO
370    
371               CALL ustarhb(knon, yu, yv, coefm(:knon, 1), yustar)               ustar(:knon) = ustarhb(yu(:knon, 1), yv(:knon, 1), coefm(:knon, 1))
372               IF (prt_level > 9) PRINT *, 'USTAR = ', yustar               CALL yamada4(dtime, rg, zlev(:knon, :), yzlay(:knon, :), &
373                      yu(:knon, :), yv(:knon, :), yteta(:knon, :), &
374               ! iflag_pbl peut \^etre utilis\'e comme longueur de m\'elange                    coefm(:knon, 1), yq2(:knon, :), ykmm(:knon, :), &
375                      ykmn(:knon, :), ykmq(: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  
   
376               coefm(:knon, 2:) = ykmm(:knon, 2:klev)               coefm(:knon, 2:) = ykmm(:knon, 2:klev)
377               coefh(:knon, 2:) = ykmn(:knon, 2:klev)               coefh(:knon, 2:) = ykmn(:knon, 2:klev)
378            END IF            END IF
379    
380            ! calculer la diffusion des vitesses "u" et "v"            CALL clvent(dtime, yu(:knon, 1), yv(:knon, 1), coefm(:knon, :), &
381            CALL clvent(knon, dtime, yu1, yv1, coefm(:knon, :), yt, yu, ypaprs, &                 yt(:knon, :), yu(:knon, :), ypaprs(:knon, :), &
382                 ypplay, ydelp, y_d_u, y_flux_u)                 ypplay(:knon, :), ydelp(:knon, :), y_d_u(:knon, :), &
383            CALL clvent(knon, dtime, yu1, yv1, coefm(:knon, :), yt, yv, ypaprs, &                 y_flux_u(:knon))
384                 ypplay, ydelp, y_d_v, y_flux_v)            CALL clvent(dtime, yu(:knon, 1), yv(:knon, 1), coefm(:knon, :), &
385                   yt(:knon, :), yv(:knon, :), ypaprs(:knon, :), &
386                   ypplay(:knon, :), ydelp(:knon, :), y_d_v(:knon, :), &
387                   y_flux_v(:knon))
388    
389            ! calculer la diffusion de "q" et de "h"            ! calculer la diffusion de "q" et de "h"
390            CALL clqh(dtime, itap, jour, debut, rlat, knon, nsrf, ni(:knon), &            CALL clqh(dtime, julien, firstcal, nsrf, ni(:knon), &
391                 pctsrf, ytsoil, yqsol, rmu0, yrugos, yrugoro, yu1, &                 ytsoil(:knon, :), yqsol(:knon), mu0, yrugos, yrugoro, &
392                 yv1, coefh(:knon, :), yt, yq, yts, ypaprs, ypplay, ydelp, &                 yu(:knon, 1), yv(:knon, 1), coefh(:knon, :), yt, yq, &
393                 yrads, yalb(:knon), ysnow, yqsurf, yrain_f, ysnow_f, yfder, &                 yts(:knon), ypaprs, ypplay, ydelp, yrads(:knon), yalb(:knon), &
394                 yfluxlat, pctsrf_new, yagesno(:knon), y_d_t, y_d_q, &                 snow(:knon), yqsurf, yrain_f, ysnow_f, yfluxlat(:knon), &
395                 y_d_ts(:knon), yz0_new, y_flux_t, y_flux_q, y_dflux_t, &                 pctsrf_new_sic, yagesno(:knon), y_d_t, y_d_q, y_d_ts(:knon), &
396                 y_dflux_q, y_fqcalving, y_ffonte, y_run_off_lic_0)                 yz0_new, y_flux_t(:knon), y_flux_q(:knon), y_dflux_t(:knon), &
397                   y_dflux_q(:knon), y_fqcalving, y_ffonte, y_run_off_lic_0)
398    
399            ! calculer la longueur de rugosite sur ocean            ! calculer la longueur de rugosite sur ocean
400            yrugm = 0.            yrugm = 0.
401            IF (nsrf == is_oce) THEN            IF (nsrf == is_oce) THEN
402               DO j = 1, knon               DO j = 1, knon
403                  yrugm(j) = 0.018*coefm(j, 1)*(yu1(j)**2+yv1(j)**2)/rg + &                  yrugm(j) = 0.018 * coefm(j, 1) * (yu(j, 1)**2 + yv(j, 1)**2) &
404                       0.11*14E-6/sqrt(coefm(j, 1)*(yu1(j)**2+yv1(j)**2))                       / rg + 0.11 * 14E-6 &
405                         / sqrt(coefm(j, 1) * (yu(j, 1)**2 + yv(j, 1)**2))
406                  yrugm(j) = max(1.5E-05, yrugm(j))                  yrugm(j) = max(1.5E-05, yrugm(j))
407               END DO               END DO
408            END IF            END IF
409            DO j = 1, knon            DO j = 1, knon
410               y_dflux_t(j) = y_dflux_t(j)*ypct(j)               y_dflux_t(j) = y_dflux_t(j) * ypct(j)
411               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)  
412            END DO            END DO
413    
414            DO k = 1, klev            DO k = 1, klev
415               DO j = 1, knon               DO j = 1, knon
416                  i = ni(j)                  i = ni(j)
417                  coefh(j, k) = coefh(j, k)*ypct(j)                  coefh(j, k) = coefh(j, k) * ypct(j)
418                  coefm(j, k) = coefm(j, k)*ypct(j)                  coefm(j, k) = coefm(j, k) * ypct(j)
419                  y_d_t(j, k) = y_d_t(j, k)*ypct(j)                  y_d_t(j, k) = y_d_t(j, k) * ypct(j)
420                  y_d_q(j, k) = y_d_q(j, k)*ypct(j)                  y_d_q(j, k) = y_d_q(j, k) * ypct(j)
421                  flux_t(i, k, nsrf) = y_flux_t(j, k)                  y_d_u(j, k) = y_d_u(j, k) * ypct(j)
422                  flux_q(i, k, nsrf) = y_flux_q(j, k)                  y_d_v(j, k) = y_d_v(j, k) * ypct(j)
                 flux_u(i, k, nsrf) = y_flux_u(j, k)  
                 flux_v(i, k, nsrf) = y_flux_v(j, k)  
                 y_d_u(j, k) = y_d_u(j, k)*ypct(j)  
                 y_d_v(j, k) = y_d_v(j, k)*ypct(j)  
423               END DO               END DO
424            END DO            END DO
425    
426            evap(:, nsrf) = -flux_q(:, 1, nsrf)            flux_t(ni(:knon), nsrf) = y_flux_t(:knon)
427              flux_q(ni(:knon), nsrf) = y_flux_q(:knon)
428              flux_u(ni(:knon), nsrf) = y_flux_u(:knon)
429              flux_v(ni(:knon), nsrf) = y_flux_v(:knon)
430    
431              evap(:, nsrf) = -flux_q(:, nsrf)
432    
433            falbe(:, nsrf) = 0.            falbe(:, nsrf) = 0.
434            snow(:, nsrf) = 0.            fsnow(:, nsrf) = 0.
435            qsurf(:, nsrf) = 0.            qsurf(:, nsrf) = 0.
436            rugos(:, nsrf) = 0.            frugs(:, nsrf) = 0.
           fluxlat(:, nsrf) = 0.  
437            DO j = 1, knon            DO j = 1, knon
438               i = ni(j)               i = ni(j)
439               d_ts(i, nsrf) = y_d_ts(j)               d_ts(i, nsrf) = y_d_ts(j)
440               falbe(i, nsrf) = yalb(j)               falbe(i, nsrf) = yalb(j)
441               snow(i, nsrf) = ysnow(j)               fsnow(i, nsrf) = snow(j)
442               qsurf(i, nsrf) = yqsurf(j)               qsurf(i, nsrf) = yqsurf(j)
443               rugos(i, nsrf) = yz0_new(j)               frugs(i, nsrf) = yz0_new(j)
444               fluxlat(i, nsrf) = yfluxlat(j)               fluxlat(i, nsrf) = yfluxlat(j)
445               IF (nsrf == is_oce) THEN               IF (nsrf == is_oce) THEN
446                  rugmer(i) = yrugm(j)                  rugmer(i) = yrugm(j)
447                  rugos(i, nsrf) = yrugm(j)                  frugs(i, nsrf) = yrugm(j)
448               END IF               END IF
449               agesno(i, nsrf) = yagesno(j)               agesno(i, nsrf) = yagesno(j)
450               fqcalving(i, nsrf) = y_fqcalving(j)               fqcalving(i, nsrf) = y_fqcalving(j)
# Line 509  contains Line 453  contains
453               cdragm(i) = cdragm(i) + coefm(j, 1)               cdragm(i) = cdragm(i) + coefm(j, 1)
454               dflux_t(i) = dflux_t(i) + y_dflux_t(j)               dflux_t(i) = dflux_t(i) + y_dflux_t(j)
455               dflux_q(i) = dflux_q(i) + y_dflux_q(j)               dflux_q(i) = dflux_q(i) + y_dflux_q(j)
              zu1(i) = zu1(i) + yu1(j)  
              zv1(i) = zv1(i) + yv1(j)  
456            END DO            END DO
457            IF (nsrf == is_ter) THEN            IF (nsrf == is_ter) THEN
458               qsol(ni(:knon)) = yqsol(:knon)               qsol(ni(:knon)) = yqsol(:knon)
# Line 522  contains Line 464  contains
464            END IF            END IF
465    
466            ftsoil(:, :, nsrf) = 0.            ftsoil(:, :, nsrf) = 0.
467            DO k = 1, nsoilmx            ftsoil(ni(:knon), :, nsrf) = ytsoil(:knon, :)
              DO j = 1, knon  
                 i = ni(j)  
                 ftsoil(i, k, nsrf) = ytsoil(j, k)  
              END DO  
           END DO  
468    
469            DO j = 1, knon            DO j = 1, knon
470               i = ni(j)               i = ni(j)
# Line 544  contains Line 481  contains
481    
482            DO j = 1, knon            DO j = 1, knon
483               i = ni(j)               i = ni(j)
484               uzon(j) = yu(j, 1) + y_d_u(j, 1)               u1(j) = yu(j, 1) + y_d_u(j, 1)
485               vmer(j) = yv(j, 1) + y_d_v(j, 1)               v1(j) = yv(j, 1) + y_d_v(j, 1)
486               tair1(j) = yt(j, 1) + y_d_t(j, 1)               tair1(j) = yt(j, 1) + y_d_t(j, 1)
487               qair1(j) = yq(j, 1) + y_d_q(j, 1)               qair1(j) = yq(j, 1) + y_d_q(j, 1)
488               zgeo1(j) = rd*tair1(j)/(0.5*(ypaprs(j, 1)+ypplay(j, &               zgeo1(j) = rd * tair1(j) / (0.5 * (ypaprs(j, 1) + ypplay(j, &
489                    1)))*(ypaprs(j, 1)-ypplay(j, 1))                    1))) * (ypaprs(j, 1)-ypplay(j, 1))
490               tairsol(j) = yts(j) + y_d_ts(j)               tairsol(j) = yts(j) + y_d_ts(j)
491               rugo1(j) = yrugos(j)               rugo1(j) = yrugos(j)
492               IF (nsrf == is_oce) THEN               IF (nsrf == is_oce) THEN
493                  rugo1(j) = rugos(i, nsrf)                  rugo1(j) = frugs(i, nsrf)
494               END IF               END IF
495               psfce(j) = ypaprs(j, 1)               psfce(j) = ypaprs(j, 1)
496               patm(j) = ypplay(j, 1)               patm(j) = ypplay(j, 1)
# Line 561  contains Line 498  contains
498               qairsol(j) = yqsurf(j)               qairsol(j) = yqsurf(j)
499            END DO            END DO
500    
501            CALL stdlevvar(klon, knon, nsrf, zxli, uzon, vmer, tair1, qair1, &            CALL stdlevvar(klon, knon, nsrf, u1(:knon), v1(:knon), tair1(:knon), &
502                 zgeo1, tairsol, qairsol, rugo1, psfce, patm, yt2m, yq2m, &                 qair1, zgeo1, tairsol, qairsol, rugo1, psfce, patm, yt2m, &
503                 yt10m, yq10m, yu10m, yustar)                 yq2m, yt10m, yq10m, wind10m(:knon), ustar)
504    
505            DO j = 1, knon            DO j = 1, knon
506               i = ni(j)               i = ni(j)
507               t2m(i, nsrf) = yt2m(j)               t2m(i, nsrf) = yt2m(j)
508               q2m(i, nsrf) = yq2m(j)               q2m(i, nsrf) = yq2m(j)
509    
510               ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman               u10m_srf(i, nsrf) = (wind10m(j) * u1(j)) &
511               u10m(i, nsrf) = (yu10m(j)*uzon(j))/sqrt(uzon(j)**2+vmer(j)**2)                    / sqrt(u1(j)**2 + v1(j)**2)
512               v10m(i, nsrf) = (yu10m(j)*vmer(j))/sqrt(uzon(j)**2+vmer(j)**2)               v10m_srf(i, nsrf) = (wind10m(j) * v1(j)) &
513                      / sqrt(u1(j)**2 + v1(j)**2)
514            END DO            END DO
515    
516            CALL hbtm(knon, ypaprs, ypplay, yt2m, yq2m, yustar, y_flux_t, &            CALL hbtm(ypaprs, ypplay, yt2m, yq2m, ustar(:knon), y_flux_t(:knon), &
517                 y_flux_q, yu, yv, yt, yq, ypblh(:knon), ycapcl, yoliqcl, &                 y_flux_q(:knon), yu, yv, yt, yq, ypblh(:knon), ycapcl, &
518                 ycteicl, ypblt, ytherm, ytrmb1, ytrmb2, ytrmb3, ylcl)                 yoliqcl, ycteicl, ypblt, ytherm, ytrmb1, ytrmb2, ytrmb3, ylcl)
519    
520            DO j = 1, knon            DO j = 1, knon
521               i = ni(j)               i = ni(j)
# Line 600  contains Line 537  contains
537                  q2(i, k, nsrf) = yq2(j, k)                  q2(i, k, nsrf) = yq2(j, k)
538               END DO               END DO
539            END DO            END DO
540           else
541              fsnow(:, nsrf) = 0.
542         end IF if_knon         end IF if_knon
543      END DO loop_surface      END DO loop_surface
544    
545      ! On utilise les nouvelles surfaces      ! On utilise les nouvelles surfaces
546        frugs(:, is_oce) = rugmer
547        pctsrf(:, is_oce) = pctsrf_new_oce
548        pctsrf(:, is_sic) = pctsrf_new_sic
549    
550      rugos(:, is_oce) = rugmer      firstcal = .false.
     pctsrf = pctsrf_new  
551    
552    END SUBROUTINE clmain    END SUBROUTINE clmain
553    
Legend:
Removed from v.186  
changed lines
  Added in v.229
  ViewVC Help
Powered by ViewVC 1.1.21