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

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

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

trunk/Sources/phylmd/clmain.f revision 207 by guez, Thu Sep 1 10:30:53 2016 UTC trunk/phylmd/Interface_surf/pbl_surface.f90 revision 330 by guez, Wed Jul 31 14:55:23 2019 UTC
# Line 1  Line 1 
1  module clmain_m  module pbl_surface_m
2    
3    IMPLICIT NONE    IMPLICIT NONE
4    
5  contains  contains
6    
7    SUBROUTINE clmain(dtime, pctsrf, t, q, u, v, jour, rmu0, ftsol, cdmmax, &    SUBROUTINE pbl_surface(pctsrf, t, q, u, v, julien, mu0, ftsol, cdmmax, &
8         cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, qsol, paprs, pplay, snow, &         cdhmax, ftsoil, qsol, paprs, play, fsnow, fqsurf, falbe, fluxlat, &
9         qsurf, evap, falbe, fluxlat, rain_fall, snow_f, solsw, sollw, fder, &         rain_fall, snow_fall, frugs, agesno, rugoro, d_t, d_q, d_u, d_v, &
10         rlat, rugos, agesno, rugoro, d_t, d_q, d_u, d_v, d_ts, flux_t, flux_q, &         flux_t, flux_q, flux_u, flux_v, cdragh, cdragm, q2, dflux_t, dflux_q, &
11         flux_u, flux_v, cdragh, cdragm, q2, dflux_t, dflux_q, ycoefh, zu1, &         coefh, t2m, q2m, u10m_srf, v10m_srf, pblh, capcl, oliqcl, cteicl, pblt, &
12         zv1, t2m, q2m, u10m, v10m, pblh, capcl, oliqcl, cteicl, pblt, therm, &         therm, plcl, fqcalving, ffonte, run_off_lic_0, albsol, sollw, solsw, &
13         trmb1, trmb2, trmb3, plcl, fqcalving, ffonte, run_off_lic_0)         tsol)
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)
17        ! Date: Aug. 18th, 1993
18      ! Objet : interface de couche limite (diffusion verticale)      ! Objet : interface de couche limite (diffusion verticale)
19    
20      ! Tout ce qui a trait aux traceurs est dans "phytrac". Le calcul      ! Tout ce qui a trait aux traceurs est dans "phytrac". Le calcul
# Line 21  contains Line 22  contains
22      ! ne tient pas compte de la diff\'erentiation des sous-fractions      ! ne tient pas compte de la diff\'erentiation des sous-fractions
23      ! de sol.      ! de sol.
24    
25      ! Pour pouvoir extraire les coefficients d'\'echanges et le vent      use cdrag_m, only: cdrag
     ! 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.  
   
26      use clqh_m, only: clqh      use clqh_m, only: clqh
27      use clvent_m, only: clvent      use clvent_m, only: clvent
28      use coefkz_m, only: coefkz      use coef_diff_turb_m, only: coef_diff_turb
29      use coefkzmin_m, only: coefkzmin      USE conf_gcm_m, ONLY: lmt_pas
     USE conf_gcm_m, ONLY: prt_level, lmt_pas  
30      USE conf_phys_m, ONLY: iflag_pbl      USE conf_phys_m, ONLY: iflag_pbl
31      USE dimphy, ONLY: klev, klon, zmasq      USE dimphy, ONLY: klev, klon
32      USE dimsoil, ONLY: nsoilmx      USE dimsoil, ONLY: nsoilmx
33      use hbtm_m, only: hbtm      use hbtm_m, only: hbtm
34        USE histwrite_phy_m, ONLY: histwrite_phy
35      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
36      USE interfoce_lim_m, ONLY: interfoce_lim      USE interfoce_lim_m, ONLY: interfoce_lim
37        use phyetat0_m, only: masque
38      use stdlevvar_m, only: stdlevvar      use stdlevvar_m, only: stdlevvar
39      USE suphec_m, ONLY: rd, rg, rkappa      USE suphec_m, ONLY: rd, rg, rsigma
40      use time_phylmdz, only: itap      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)  
41    
42      REAL, INTENT(inout):: pctsrf(klon, nbsrf)      REAL, INTENT(inout):: pctsrf(:, :) ! (klon, nbsrf)
43      ! tableau des pourcentages de surface de chaque maille      ! pourcentages de surface de chaque maille
44    
45      REAL, INTENT(IN):: t(klon, klev) ! temperature (K)      REAL, INTENT(IN):: t(klon, klev) ! temperature (K)
46      REAL, INTENT(IN):: q(klon, klev) ! vapeur d'eau (kg/kg)      REAL, INTENT(IN):: q(klon, klev) ! vapeur d'eau (kg / kg)
47      REAL, INTENT(IN):: u(klon, klev), v(klon, klev) ! vitesse      REAL, INTENT(IN):: u(klon, klev), v(klon, klev) ! vitesse
48      INTEGER, INTENT(IN):: jour ! jour de l'annee en cours      INTEGER, INTENT(IN):: julien ! jour de l'annee en cours
49      REAL, intent(in):: rmu0(klon) ! cosinus de l'angle solaire zenithal          REAL, intent(in):: mu0(klon) ! cosinus de l'angle solaire zenithal    
50      REAL, INTENT(IN):: ftsol(klon, nbsrf) ! temperature du sol (en Kelvin)  
51        REAL, INTENT(INout):: ftsol(:, :) ! (klon, nbsrf)
52        ! skin temperature of surface fraction, in K
53    
54      REAL, INTENT(IN):: cdmmax, cdhmax ! seuils cdrm, cdrh      REAL, INTENT(IN):: cdmmax, cdhmax ! seuils cdrm, cdrh
     REAL, INTENT(IN):: ksta, ksta_ter  
     LOGICAL, INTENT(IN):: ok_kzmin  
55    
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):: play(klon, klev) ! pression au milieu de couche (Pa)
     REAL, INTENT(inout):: snow(klon, nbsrf)  
     REAL qsurf(klon, nbsrf)  
     REAL evap(klon, nbsrf)  
     REAL, intent(inout):: falbe(klon, nbsrf)  
64    
65      REAL fluxlat(klon, nbsrf)      REAL, INTENT(inout):: fsnow(:, :) ! (klon, nbsrf)
66        ! column-density of mass of snow, in kg m-2
67    
68      REAL, intent(in):: rain_fall(klon)      REAL, INTENT(inout):: fqsurf(klon, nbsrf)
69      ! liquid water mass flux (kg/m2/s), positive down      REAL, intent(inout):: falbe(klon, nbsrf)
70    
71      REAL, intent(in):: snow_f(klon)      REAL, intent(out):: fluxlat(:, :) ! (klon, nbsrf)
72      ! solid water mass flux (kg/m2/s), positive down      ! flux de chaleur latente, en W m-2
73    
74      REAL, INTENT(IN):: solsw(klon, nbsrf), sollw(klon, nbsrf)      REAL, intent(in):: rain_fall(klon)
75      REAL, intent(in):: fder(klon)      ! liquid water mass flux (kg / m2 / s), positive down
     REAL, INTENT(IN):: rlat(klon) ! latitude en degr\'es  
76    
77      REAL, intent(inout):: rugos(klon, nbsrf) ! longueur de rugosit\'e (en m)      REAL, intent(in):: snow_fall(klon)
78        ! solid water mass flux (kg / m2 / s), positive down
79    
80      real agesno(klon, nbsrf)      REAL, intent(inout):: frugs(klon, nbsrf) ! longueur de rugosit\'e (en m)
81        real, intent(inout):: agesno(:, :) ! (klon, nbsrf)
82      REAL, INTENT(IN):: rugoro(klon)      REAL, INTENT(IN):: rugoro(klon)
83    
84      REAL d_t(klon, klev), d_q(klon, klev)      REAL, intent(out):: d_t(:, :), d_q(:, :) ! (klon, klev)
85      ! d_t------output-R- le changement pour "t"      ! changement pour t et q
     ! d_q------output-R- le changement pour "q"  
86    
87      REAL, intent(out):: d_u(klon, klev), d_v(klon, klev)      REAL, intent(out):: d_u(klon, klev), d_v(klon, klev)
88      ! changement pour "u" et "v"      ! changement pour "u" et "v"
89    
     REAL, intent(out):: d_ts(klon, nbsrf) ! le changement pour "ftsol"  
   
90      REAL, intent(out):: flux_t(klon, nbsrf)      REAL, intent(out):: flux_t(klon, nbsrf)
91      ! flux de chaleur sensible (Cp T) (W/m2) (orientation positive vers      ! flux de chaleur sensible (c_p T) (W / m2) (orientation positive
92      ! le bas) à la surface      ! vers le bas) à la surface
93    
94      REAL, intent(out):: flux_q(klon, nbsrf)      REAL, intent(out):: flux_q(klon, nbsrf)
95      ! flux de vapeur d'eau (kg/m2/s) à la surface      ! flux de vapeur d'eau (kg / m2 / s) à la surface
96    
97      REAL, intent(out):: flux_u(klon, nbsrf), flux_v(klon, nbsrf)      REAL, intent(out):: flux_u(:, :), flux_v(:, :) ! (klon, nbsrf)
98      ! tension du vent à la surface, en Pa      ! 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        ! Ocean slab:
104        REAL, INTENT(out):: dflux_t(klon) ! derive du flux sensible
105        REAL, INTENT(out):: dflux_q(klon) ! derive du flux latent
106    
107        REAL, intent(out):: coefh(:, 2:) ! (klon, 2:klev)
108        ! Pour pouvoir extraire les coefficients d'\'echange, le champ
109        ! "coefh" a \'et\'e cr\'e\'e. Nous avons moyenn\'e les valeurs de
110        ! ce champ sur les quatre sous-surfaces du mod\`ele.
111    
112      REAL, INTENT(out):: dflux_t(klon), dflux_q(klon)      REAL, INTENT(inout):: t2m(klon, nbsrf), q2m(klon, nbsrf)
113      ! dflux_t derive du flux sensible  
114      ! dflux_q derive du flux latent      REAL, INTENT(inout):: u10m_srf(:, :), v10m_srf(:, :) ! (klon, nbsrf)
115      ! IM "slab" ocean      ! composantes du vent \`a 10m sans spirale d'Ekman
116    
117      REAL, intent(out):: ycoefh(klon, klev)      ! Ionela Musat. Cf. Anne Mathieu : planetary boundary layer, hbtm.
118      REAL, intent(out):: zu1(klon)      ! Comme les autres diagnostics on cumule dans physiq ce qui permet
119      REAL zv1(klon)      ! de sortir les grandeurs par sous-surface.
     REAL t2m(klon, nbsrf), q2m(klon, nbsrf)  
     REAL u10m(klon, nbsrf), v10m(klon, nbsrf)  
   
     ! 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)  
120      REAL pblh(klon, nbsrf) ! height of planetary boundary layer      REAL pblh(klon, nbsrf) ! height of planetary boundary layer
121      REAL capcl(klon, nbsrf)      REAL capcl(klon, nbsrf)
122      REAL oliqcl(klon, nbsrf)      REAL oliqcl(klon, nbsrf)
123      REAL cteicl(klon, nbsrf)      REAL cteicl(klon, nbsrf)
124      REAL pblt(klon, nbsrf)      REAL, INTENT(inout):: pblt(klon, nbsrf) ! T au nveau HCL
     ! pblT------- T au nveau HCL  
125      REAL therm(klon, nbsrf)      REAL therm(klon, nbsrf)
     REAL trmb1(klon, nbsrf)  
     ! trmb1-------deep_cape  
     REAL trmb2(klon, nbsrf)  
     ! trmb2--------inhibition  
     REAL trmb3(klon, nbsrf)  
     ! trmb3-------Point Omega  
126      REAL plcl(klon, nbsrf)      REAL plcl(klon, nbsrf)
127      REAL fqcalving(klon, nbsrf), ffonte(klon, nbsrf)  
128      ! ffonte----Flux thermique utilise pour fondre la neige      REAL, intent(out):: fqcalving(klon, nbsrf)
129      ! fqcalving-Flux d'eau "perdue" par la surface et necessaire pour limiter la      ! flux d'eau "perdue" par la surface et necessaire pour limiter la
130      !           hauteur de neige, en kg/m2/s      ! hauteur de neige, en kg / m2 / s
131      REAL run_off_lic_0(klon)  
132        real ffonte(klon, nbsrf) ! flux thermique utilise pour fondre la neige
133        REAL, intent(inout):: run_off_lic_0(:) ! (klon)
134    
135        REAL, intent(out):: albsol(:) ! (klon)
136        ! albedo du sol total, visible, moyen par maille
137    
138        REAL, intent(in):: sollw(:) ! (klon)
139        ! surface net downward longwave flux, in W m-2
140    
141        REAL, intent(in):: solsw(:) ! (klon)
142        ! surface net downward shortwave flux, in W m-2
143    
144        REAL, intent(in):: tsol(:) ! (klon)
145    
146      ! Local:      ! Local:
147    
148      LOGICAL:: firstcal = .true.      REAL d_ts(klon, nbsrf) ! variation of ftsol
149        REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous-surface
150        REAL fsolsw(klon, nbsrf) ! flux solaire absorb\'e pour chaque sous-surface
151    
152      ! la nouvelle repartition des surfaces sortie de l'interface      ! la nouvelle repartition des surfaces sortie de l'interface
153      REAL, save:: pctsrf_new_oce(klon)      REAL, save:: pctsrf_new_oce(klon)
154      REAL, save:: pctsrf_new_sic(klon)      REAL, save:: pctsrf_new_sic(klon)
155    
156      REAL y_fqcalving(klon), y_ffonte(klon)      REAL y_fqcalving(klon), y_ffonte(klon)
157      real y_run_off_lic_0(klon)      real y_run_off_lic_0(klon), y_run_off_lic(klon)
158        REAL run_off_lic(klon) ! ruissellement total
159      REAL rugmer(klon)      REAL rugmer(klon)
160      REAL ytsoil(klon, nsoilmx)      REAL ytsoil(klon, nsoilmx)
161      REAL yts(klon), yrugos(klon), ypct(klon), yz0_new(klon)      REAL yts(klon), ypctsrf(klon), yz0_new(klon)
162        real yrugos(klon) ! longueur de rugosite (en m)
163      REAL yalb(klon)      REAL yalb(klon)
164      REAL yu1(klon), yv1(klon)      REAL snow(klon) ! column-density of mass of snow, in kg m-2
165      ! on rajoute en output yu1 et yv1 qui sont les vents dans      real yqsurf(klon), yagesno(klon)
166      ! la premiere couche      real yqsol(klon) ! column-density of water in soil, in kg m-2
167      REAL ysnow(klon), yqsurf(klon), yagesno(klon)      REAL yrain_fall(klon) ! liquid water mass flux (kg / m2 / s), positive down
168        REAL ysnow_fall(klon) ! solid water mass flux (kg / m2 / s), positive down
169      real yqsol(klon)      REAL yrugm(klon), radsol(klon), yrugoro(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)  
     REAL yrugm(klon), yrads(klon), yrugoro(klon)  
   
170      REAL yfluxlat(klon)      REAL yfluxlat(klon)
   
171      REAL y_d_ts(klon)      REAL y_d_ts(klon)
172      REAL y_d_t(klon, klev), y_d_q(klon, klev)      REAL y_d_t(klon, klev), y_d_q(klon, klev)
173      REAL y_d_u(klon, klev), y_d_v(klon, klev)      REAL y_d_u(klon, klev), y_d_v(klon, klev)
174      REAL y_flux_t(klon), y_flux_q(klon)      REAL y_flux_t(klon), y_flux_q(klon)
175      REAL y_flux_u(klon), y_flux_v(klon)      REAL y_flux_u(klon), y_flux_v(klon)
176      REAL y_dflux_t(klon), y_dflux_q(klon)      REAL y_dflux_t(klon), y_dflux_q(klon)
177      REAL coefh(klon, klev), coefm(klon, klev)      REAL ycoefh(klon, 2:klev), ycoefm(klon, 2:klev)
178        real ycdragh(klon), ycdragm(klon)
179      REAL yu(klon, klev), yv(klon, klev)      REAL yu(klon, klev), yv(klon, klev)
180      REAL yt(klon, klev), yq(klon, klev)      REAL yt(klon, klev), yq(klon, klev)
181      REAL ypaprs(klon, klev+1), ypplay(klon, klev), ydelp(klon, klev)      REAL ypaprs(klon, klev + 1), ypplay(klon, klev), ydelp(klon, klev)
182        REAL yq2(klon, klev + 1)
     REAL ycoefm0(klon, klev), ycoefh0(klon, klev)  
   
     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)  
183      REAL delp(klon, klev)      REAL delp(klon, klev)
184      INTEGER i, k, nsrf      INTEGER i, k, nsrf
   
185      INTEGER ni(klon), knon, j      INTEGER ni(klon), knon, j
186    
187      REAL pctsrf_pot(klon, nbsrf)      REAL pctsrf_pot(klon, nbsrf)
188      ! "pourcentage potentiel" pour tenir compte des \'eventuelles      ! "pourcentage potentiel" pour tenir compte des \'eventuelles
189      ! apparitions ou disparitions de la glace de mer      ! apparitions ou disparitions de la glace de mer
190    
191      REAL zx_alf1, zx_alf2 ! valeur ambiante par extrapolation      REAL yt2m(klon), yq2m(klon), wind10m(klon)
192        REAL ustar(klon)
     REAL yt2m(klon), yq2m(klon), yu10m(klon)  
     REAL yustar(klon)  
193    
194      REAL yt10m(klon), yq10m(klon)      REAL yt10m(klon), yq10m(klon)
195      REAL ypblh(klon)      REAL ypblh(klon)
# Line 220  contains Line 199  contains
199      REAL ycteicl(klon)      REAL ycteicl(klon)
200      REAL ypblt(klon)      REAL ypblt(klon)
201      REAL ytherm(klon)      REAL ytherm(klon)
202      REAL ytrmb1(klon)      REAL u1(klon), v1(klon)
     REAL ytrmb2(klon)  
     REAL ytrmb3(klon)  
     REAL uzon(klon), vmer(klon)  
203      REAL tair1(klon), qair1(klon), tairsol(klon)      REAL tair1(klon), qair1(klon), tairsol(klon)
204      REAL psfce(klon), patm(klon)      REAL psfce(klon), patm(klon)
205        REAL zgeo1(klon)
     REAL qairsol(klon), zgeo1(klon)  
206      REAL rugo1(klon)      REAL rugo1(klon)
207        REAL zgeop(klon, klev)
     ! utiliser un jeu de fonctions simples                
     LOGICAL zxli  
     PARAMETER (zxli=.FALSE.)  
208    
209      !------------------------------------------------------------      !------------------------------------------------------------
210    
211        albsol = sum(falbe * pctsrf, dim = 2)
212    
213        ! R\'epartition sous maille des flux longwave et shortwave
214        ! R\'epartition du longwave par sous-surface lin\'earis\'ee
215    
216        forall (nsrf = 1:nbsrf)
217           fsollw(:, nsrf) = sollw + 4. * RSIGMA * tsol**3 &
218                * (tsol - ftsol(:, nsrf))
219           fsolsw(:, nsrf) = solsw * (1. - falbe(:, nsrf)) / (1. - albsol)
220        END forall
221    
222      ytherm = 0.      ytherm = 0.
223    
224      DO k = 1, klev ! epaisseur de couche      DO k = 1, klev ! epaisseur de couche
225         DO i = 1, klon         DO i = 1, klon
226            delp(i, k) = paprs(i, k) - paprs(i, k+1)            delp(i, k) = paprs(i, k) - paprs(i, k + 1)
227         END DO         END DO
228      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  
229    
230      ! Initialization:      ! Initialization:
231      rugmer = 0.      rugmer = 0.
# Line 256  contains Line 233  contains
233      cdragm = 0.      cdragm = 0.
234      dflux_t = 0.      dflux_t = 0.
235      dflux_q = 0.      dflux_q = 0.
     zu1 = 0.  
     zv1 = 0.  
     ypct = 0.  
     yts = 0.  
     ysnow = 0.  
     yqsurf = 0.  
     yrain_f = 0.  
     ysnow_f = 0.  
     yfder = 0.  
236      yrugos = 0.      yrugos = 0.
     yu1 = 0.  
     yv1 = 0.  
     yrads = 0.  
237      ypaprs = 0.      ypaprs = 0.
238      ypplay = 0.      ypplay = 0.
239      ydelp = 0.      ydelp = 0.
     yu = 0.  
     yv = 0.  
     yt = 0.  
     yq = 0.  
     y_dflux_t = 0.  
     y_dflux_q = 0.  
     ytsoil = 999999.  
240      yrugoro = 0.      yrugoro = 0.
241      d_ts = 0.      d_ts = 0.
     yfluxlat = 0.  
242      flux_t = 0.      flux_t = 0.
243      flux_q = 0.      flux_q = 0.
244      flux_u = 0.      flux_u = 0.
245      flux_v = 0.      flux_v = 0.
246        fluxlat = 0.
247      d_t = 0.      d_t = 0.
248      d_q = 0.      d_q = 0.
249      d_u = 0.      d_u = 0.
250      d_v = 0.      d_v = 0.
251      ycoefh = 0.      coefh = 0.
252        fqcalving = 0.
253        run_off_lic = 0.
254    
255      ! Initialisation des "pourcentages potentiels". On consid\`ere ici qu'on      ! Initialisation des "pourcentages potentiels". On consid\`ere ici qu'on
256      ! peut avoir potentiellement de la glace sur tout le domaine oc\'eanique      ! peut avoir potentiellement de la glace sur tout le domaine oc\'eanique
257      ! (\`a affiner)      ! (\`a affiner).
258    
259      pctsrf_pot(:, is_ter) = pctsrf(:, is_ter)      pctsrf_pot(:, is_ter) = pctsrf(:, is_ter)
260      pctsrf_pot(:, is_lic) = pctsrf(:, is_lic)      pctsrf_pot(:, is_lic) = pctsrf(:, is_lic)
261      pctsrf_pot(:, is_oce) = 1. - zmasq      pctsrf_pot(:, is_oce) = 1. - masque
262      pctsrf_pot(:, is_sic) = 1. - zmasq      pctsrf_pot(:, is_sic) = 1. - masque
263    
264      ! Tester si c'est le moment de lire le fichier:      ! Tester si c'est le moment de lire le fichier:
265      if (mod(itap - 1, lmt_pas) == 0) then      if (mod(itap - 1, lmt_pas) == 0) then
266         CALL interfoce_lim(jour, pctsrf_new_oce, pctsrf_new_sic)         CALL interfoce_lim(julien, pctsrf_new_oce, pctsrf_new_sic)
267      endif      endif
268    
269      ! Boucler sur toutes les sous-fractions du sol:      ! Boucler sur toutes les sous-fractions du sol:
270    
271      loop_surface: DO nsrf = 1, nbsrf      loop_surface: DO nsrf = 1, nbsrf
272         ! Chercher les indices :         ! Define ni and knon:
273    
274         ni = 0         ni = 0
275         knon = 0         knon = 0
276    
277         DO i = 1, klon         DO i = 1, klon
278            ! Pour d\'eterminer le domaine \`a traiter, on utilise les surfaces            ! Pour d\'eterminer le domaine \`a traiter, on utilise les surfaces
279            ! "potentielles"            ! "potentielles"
# Line 322  contains Line 284  contains
284         END DO         END DO
285    
286         if_knon: IF (knon /= 0) then         if_knon: IF (knon /= 0) then
287            DO j = 1, knon            ypctsrf(:knon) = pctsrf(ni(:knon), nsrf)
288               i = ni(j)            yts(:knon) = ftsol(ni(:knon), nsrf)
289               ypct(j) = pctsrf(i, nsrf)            snow(:knon) = fsnow(ni(:knon), nsrf)
290               yts(j) = ftsol(i, nsrf)            yqsurf(:knon) = fqsurf(ni(:knon), nsrf)
291               ysnow(j) = snow(i, nsrf)            yalb(:knon) = falbe(ni(:knon), nsrf)
292               yqsurf(j) = qsurf(i, nsrf)            yrain_fall(:knon) = rain_fall(ni(:knon))
293               yalb(j) = falbe(i, nsrf)            ysnow_fall(:knon) = snow_fall(ni(:knon))
294               yrain_f(j) = rain_fall(i)            yagesno(:knon) = agesno(ni(:knon), nsrf)
295               ysnow_f(j) = snow_f(i)            yrugos(:knon) = frugs(ni(:knon), nsrf)
296               yagesno(j) = agesno(i, nsrf)            yrugoro(:knon) = rugoro(ni(:knon))
297               yfder(j) = fder(i)            radsol(:knon) = fsolsw(ni(:knon), nsrf) + fsollw(ni(:knon), nsrf)
298               yrugos(j) = rugos(i, nsrf)            ypaprs(:knon, klev + 1) = paprs(ni(:knon), klev + 1)
299               yrugoro(j) = rugoro(i)            y_run_off_lic_0(:knon) = run_off_lic_0(ni(:knon))
              yu1(j) = u1lay(i)  
              yv1(j) = v1lay(i)  
              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  
300    
301            ! For continent, copy soil water content            ! For continent, copy soil water content
302            IF (nsrf == is_ter) THEN            IF (nsrf == is_ter) yqsol(:knon) = qsol(ni(:knon))
              yqsol(:knon) = qsol(ni(:knon))  
           ELSE  
              yqsol = 0.  
           END IF  
303    
304            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  
305    
306            DO k = 1, klev            DO k = 1, klev
307               DO j = 1, knon               DO j = 1, knon
308                  i = ni(j)                  i = ni(j)
309                  ypaprs(j, k) = paprs(i, k)                  ypaprs(j, k) = paprs(i, k)
310                  ypplay(j, k) = pplay(i, k)                  ypplay(j, k) = play(i, k)
311                  ydelp(j, k) = delp(i, k)                  ydelp(j, k) = delp(i, k)
312                  yu(j, k) = u(i, k)                  yu(j, k) = u(i, k)
313                  yv(j, k) = v(i, k)                  yv(j, k) = v(i, k)
# Line 369  contains Line 316  contains
316               END DO               END DO
317            END DO            END DO
318    
319            ! calculer Cdrag et les coefficients d'echange            ! Calculer les géopotentiels de chaque couche:
           CALL coefkz(nsrf, knon, ypaprs, ypplay, ksta, ksta_ter, yts, yrugos, &  
                yu, yv, yt, yq, yqsurf, coefm(:knon, :), coefh(:knon, :))  
           IF (iflag_pbl == 1) THEN  
              CALL coefkz2(nsrf, knon, ypaprs, ypplay, yt, ycoefm0, ycoefh0)  
              coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))  
              coefh(:knon, :) = max(coefh(:knon, :), ycoefh0(:knon, :))  
           END IF  
320    
321            ! on met un seuil pour coefm et coefh            zgeop(:knon, 1) = RD * yt(:knon, 1) / (0.5 * (ypaprs(:knon, 1) &
322            IF (nsrf == is_oce) THEN                 + ypplay(:knon, 1))) * (ypaprs(:knon, 1) - ypplay(:knon, 1))
              coefm(:knon, 1) = min(coefm(:knon, 1), cdmmax)  
              coefh(:knon, 1) = min(coefh(:knon, 1), cdhmax)  
           END IF  
323    
324            IF (ok_kzmin) THEN            DO k = 2, klev
325               ! Calcul d'une diffusion minimale pour les conditions tres stables               zgeop(:knon, k) = zgeop(:knon, k - 1) + RD * 0.5 &
326               CALL coefkzmin(knon, ypaprs, ypplay, yu, yv, yt, yq, &                    * (yt(:knon, k - 1) + yt(:knon, k)) / ypaprs(:knon, k) &
327                    coefm(:knon, 1), ycoefm0, ycoefh0)                    * (ypplay(:knon, k - 1) - ypplay(:knon, k))
328               coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))            ENDDO
329               coefh(:knon, :) = max(coefh(:knon, :), ycoefh0(:knon, :))  
330            END IF            CALL cdrag(nsrf, sqrt(yu(:knon, 1)**2 + yv(:knon, 1)**2), &
331                   yt(:knon, 1), yq(:knon, 1), zgeop(:knon, 1), ypaprs(:knon, 1), &
332                   yts(:knon), yqsurf(:knon), yrugos(:knon), ycdragm(:knon), &
333                   ycdragh(:knon))
334    
335            IF (iflag_pbl >= 3) THEN            IF (iflag_pbl == 1) THEN
336               ! Mellor et Yamada adapt\'e \`a Mars, Richard Fournier et               ycdragm(:knon) = max(ycdragm(:knon), 0.)
337               ! Fr\'ed\'eric Hourdin               ycdragh(:knon) = max(ycdragh(:knon), 0.)
338               yzlay(:knon, 1) = rd * yt(:knon, 1) / (0.5 * (ypaprs(:knon, 1) &            end IF
                   + ypplay(:knon, 1))) &  
                   * (ypaprs(:knon, 1) - ypplay(:knon, 1)) / rg  
              DO k = 2, klev  
                 yzlay(1:knon, k) = yzlay(1: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) &  
                      / 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) &  
                   - yzlay(:knon, klev - 1)  
              DO k = 2, klev  
                 yzlev(1:knon, k) = 0.5*(yzlay(1:knon, k)+yzlay(1:knon, k-1))  
              END DO  
              DO k = 1, klev + 1  
                 DO j = 1, knon  
                    i = ni(j)  
                    yq2(j, k) = q2(i, k, nsrf)  
                 END DO  
              END DO  
   
              CALL ustarhb(knon, yu, yv, coefm(:knon, 1), yustar)  
              IF (prt_level > 9) PRINT *, 'USTAR = ', yustar  
   
              ! iflag_pbl peut \^etre utilis\'e comme longueur de m\'elange  
   
              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  
339    
340               coefm(:knon, 2:) = ykmm(:knon, 2:klev)            ! on met un seuil pour ycdragm et ycdragh
341               coefh(:knon, 2:) = ykmn(:knon, 2:klev)            IF (nsrf == is_oce) THEN
342                 ycdragm(:knon) = min(ycdragm(:knon), cdmmax)
343                 ycdragh(:knon) = min(ycdragh(:knon), cdhmax)
344            END IF            END IF
345    
346            ! calculer la diffusion des vitesses "u" et "v"            IF (iflag_pbl >= 6) yq2(:knon, :) = q2(ni(:knon), :, nsrf)
347            CALL clvent(knon, dtime, yu1, yv1, coefm(:knon, :), yt, yu, ypaprs, &            call coef_diff_turb(nsrf, ni(:knon), ypaprs(:knon, :), &
348                 ypplay, ydelp, y_d_u, y_flux_u(:knon))                 ypplay(:knon, :), yu(:knon, :), yv(:knon, :), yq(:knon, :), &
349            CALL clvent(knon, dtime, yu1, yv1, coefm(:knon, :), yt, yv, ypaprs, &                 yt(:knon, :), yts(:knon), ycdragm(:knon), zgeop(:knon, :), &
350                 ypplay, ydelp, y_d_v, y_flux_v(:knon))                 ycoefm(:knon, :), ycoefh(:knon, :), yq2(:knon, :))
351    
352            ! calculer la diffusion de "q" et de "h"            CALL clvent(yu(:knon, 1), yv(:knon, 1), ycoefm(:knon, :), &
353            CALL clqh(dtime, jour, firstcal, rlat, nsrf, ni(:knon), ytsoil, &                 ycdragm(:knon), yt(:knon, :), yu(:knon, :), ypaprs(:knon, :), &
354                 yqsol, rmu0, yrugos, yrugoro, yu1, yv1, coefh(:knon, :), yt, &                 ypplay(:knon, :), ydelp(:knon, :), y_d_u(:knon, :), &
355                 yq, yts(:knon), ypaprs, ypplay, ydelp, yrads, yalb(:knon), &                 y_flux_u(:knon))
356                 ysnow, yqsurf, yrain_f, ysnow_f, yfder, yfluxlat, &            CALL clvent(yu(:knon, 1), yv(:knon, 1), ycoefm(:knon, :), &
357                 pctsrf_new_sic, yagesno(:knon), y_d_t, y_d_q, y_d_ts(:knon), &                 ycdragm(:knon), yt(:knon, :), yv(:knon, :), ypaprs(:knon, :), &
358                 yz0_new, y_flux_t(:knon), y_flux_q(:knon), y_dflux_t, &                 ypplay(:knon, :), ydelp(:knon, :), y_d_v(:knon, :), &
359                 y_dflux_q, y_fqcalving, y_ffonte, y_run_off_lic_0)                 y_flux_v(:knon))
360    
361              CALL clqh(julien, nsrf, ni(:knon), ytsoil(:knon, :), yqsol(:knon), &
362                   mu0(ni(:knon)), yrugos(:knon), yrugoro(:knon), yu(:knon, 1), &
363                   yv(:knon, 1), ycoefh(:knon, :), ycdragh(:knon), yt(:knon, :), &
364                   yq(:knon, :), yts(:knon), ypaprs(:knon, :), ypplay(:knon, :), &
365                   ydelp(:knon, :), radsol(:knon), yalb(:knon), snow(:knon), &
366                   yqsurf(:knon), yrain_fall(:knon), ysnow_fall(:knon), &
367                   yfluxlat(:knon), pctsrf_new_sic(ni(:knon)), yagesno(:knon), &
368                   y_d_t(:knon, :), y_d_q(:knon, :), y_d_ts(:knon), &
369                   yz0_new(:knon), y_flux_t(:knon), y_flux_q(:knon), &
370                   y_dflux_t(:knon), y_dflux_q(:knon), y_fqcalving(:knon), &
371                   y_ffonte(:knon), y_run_off_lic_0(:knon), y_run_off_lic(:knon))
372    
373            ! calculer la longueur de rugosite sur ocean            ! calculer la longueur de rugosite sur ocean
374    
375            yrugm = 0.            yrugm = 0.
376    
377            IF (nsrf == is_oce) THEN            IF (nsrf == is_oce) THEN
378               DO j = 1, knon               DO j = 1, knon
379                  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) &
380                       0.11*14E-6/sqrt(coefm(j, 1)*(yu1(j)**2+yv1(j)**2))                       / rg + 0.11 * 14E-6 &
381                         / sqrt(ycdragm(j) * (yu(j, 1)**2 + yv(j, 1)**2))
382                  yrugm(j) = max(1.5E-05, yrugm(j))                  yrugm(j) = max(1.5E-05, yrugm(j))
383               END DO               END DO
384            END IF            END IF
           DO j = 1, knon  
              y_dflux_t(j) = y_dflux_t(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  
385    
386            DO k = 1, klev            DO k = 1, klev
387               DO j = 1, knon               DO j = 1, knon
388                  i = ni(j)                  i = ni(j)
389                  coefh(j, k) = coefh(j, k)*ypct(j)                  y_d_t(j, k) = y_d_t(j, k) * ypctsrf(j)
390                  coefm(j, k) = coefm(j, k)*ypct(j)                  y_d_q(j, k) = y_d_q(j, k) * ypctsrf(j)
391                  y_d_t(j, k) = y_d_t(j, k)*ypct(j)                  y_d_u(j, k) = y_d_u(j, k) * ypctsrf(j)
392                  y_d_q(j, k) = y_d_q(j, k)*ypct(j)                  y_d_v(j, k) = y_d_v(j, k) * ypctsrf(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)  
393               END DO               END DO
394            END DO            END DO
395    
396            DO j = 1, knon            flux_t(ni(:knon), nsrf) = y_flux_t(:knon)
397               i = ni(j)            flux_q(ni(:knon), nsrf) = y_flux_q(:knon)
398               flux_t(i, nsrf) = y_flux_t(j)            flux_u(ni(:knon), nsrf) = y_flux_u(:knon)
399               flux_q(i, nsrf) = y_flux_q(j)            flux_v(ni(:knon), nsrf) = y_flux_v(:knon)
              flux_u(i, nsrf) = y_flux_u(j)  
              flux_v(i, nsrf) = y_flux_v(j)  
           END DO  
   
           evap(:, nsrf) = -flux_q(:, nsrf)  
400    
401            falbe(:, nsrf) = 0.            falbe(:, nsrf) = 0.
402            snow(:, nsrf) = 0.            fsnow(:, nsrf) = 0.
403            qsurf(:, nsrf) = 0.            fqsurf(:, nsrf) = 0.
404            rugos(:, nsrf) = 0.            frugs(:, nsrf) = 0.
           fluxlat(:, nsrf) = 0.  
405            DO j = 1, knon            DO j = 1, knon
406               i = ni(j)               i = ni(j)
407               d_ts(i, nsrf) = y_d_ts(j)               d_ts(i, nsrf) = y_d_ts(j)
408               falbe(i, nsrf) = yalb(j)               falbe(i, nsrf) = yalb(j)
409               snow(i, nsrf) = ysnow(j)               fsnow(i, nsrf) = snow(j)
410               qsurf(i, nsrf) = yqsurf(j)               fqsurf(i, nsrf) = yqsurf(j)
411               rugos(i, nsrf) = yz0_new(j)               frugs(i, nsrf) = yz0_new(j)
412               fluxlat(i, nsrf) = yfluxlat(j)               fluxlat(i, nsrf) = yfluxlat(j)
413               IF (nsrf == is_oce) THEN               IF (nsrf == is_oce) THEN
414                  rugmer(i) = yrugm(j)                  rugmer(i) = yrugm(j)
415                  rugos(i, nsrf) = yrugm(j)                  frugs(i, nsrf) = yrugm(j)
416               END IF               END IF
417               agesno(i, nsrf) = yagesno(j)               agesno(i, nsrf) = yagesno(j)
418               fqcalving(i, nsrf) = y_fqcalving(j)               fqcalving(i, nsrf) = y_fqcalving(j)
419               ffonte(i, nsrf) = y_ffonte(j)               ffonte(i, nsrf) = y_ffonte(j)
420               cdragh(i) = cdragh(i) + coefh(j, 1)               cdragh(i) = cdragh(i) + ycdragh(j) * ypctsrf(j)
421               cdragm(i) = cdragm(i) + coefm(j, 1)               cdragm(i) = cdragm(i) + ycdragm(j) * ypctsrf(j)
422               dflux_t(i) = dflux_t(i) + y_dflux_t(j)               dflux_t(i) = dflux_t(i) + y_dflux_t(j) * ypctsrf(j)
423               dflux_q(i) = dflux_q(i) + y_dflux_q(j)               dflux_q(i) = dflux_q(i) + y_dflux_q(j) * ypctsrf(j)
              zu1(i) = zu1(i) + yu1(j)  
              zv1(i) = zv1(i) + yv1(j)  
424            END DO            END DO
425            IF (nsrf == is_ter) THEN            IF (nsrf == is_ter) THEN
426               qsol(ni(:knon)) = yqsol(:knon)               qsol(ni(:knon)) = yqsol(:knon)
# Line 525  contains Line 428  contains
428               DO j = 1, knon               DO j = 1, knon
429                  i = ni(j)                  i = ni(j)
430                  run_off_lic_0(i) = y_run_off_lic_0(j)                  run_off_lic_0(i) = y_run_off_lic_0(j)
431                    run_off_lic(i) = y_run_off_lic(j)
432               END DO               END DO
433            END IF            END IF
434    
435            ftsoil(:, :, nsrf) = 0.            ftsoil(:, :, nsrf) = 0.
436            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  
437    
438            DO j = 1, knon            DO j = 1, knon
439               i = ni(j)               i = ni(j)
# Line 543  contains Line 442  contains
442                  d_q(i, k) = d_q(i, k) + y_d_q(j, k)                  d_q(i, k) = d_q(i, k) + y_d_q(j, k)
443                  d_u(i, k) = d_u(i, k) + y_d_u(j, k)                  d_u(i, k) = d_u(i, k) + y_d_u(j, k)
444                  d_v(i, k) = d_v(i, k) + y_d_v(j, k)                  d_v(i, k) = d_v(i, k) + y_d_v(j, k)
                 ycoefh(i, k) = ycoefh(i, k) + coefh(j, k)  
445               END DO               END DO
446            END DO            END DO
447    
448              forall (k = 2:klev) coefh(ni(:knon), k) &
449                   = coefh(ni(:knon), k) + ycoefh(:knon, k) * ypctsrf(:knon)
450    
451            ! diagnostic t, q a 2m et u, v a 10m            ! diagnostic t, q a 2m et u, v a 10m
452    
453            DO j = 1, knon            DO j = 1, knon
454               i = ni(j)               i = ni(j)
455               uzon(j) = yu(j, 1) + y_d_u(j, 1)               u1(j) = yu(j, 1) + y_d_u(j, 1)
456               vmer(j) = yv(j, 1) + y_d_v(j, 1)               v1(j) = yv(j, 1) + y_d_v(j, 1)
457               tair1(j) = yt(j, 1) + y_d_t(j, 1)               tair1(j) = yt(j, 1) + y_d_t(j, 1)
458               qair1(j) = yq(j, 1) + y_d_q(j, 1)               qair1(j) = yq(j, 1) + y_d_q(j, 1)
459               zgeo1(j) = rd*tair1(j)/(0.5*(ypaprs(j, 1)+ypplay(j, &               zgeo1(j) = rd * tair1(j) / (0.5 * (ypaprs(j, 1) + ypplay(j, &
460                    1)))*(ypaprs(j, 1)-ypplay(j, 1))                    1))) * (ypaprs(j, 1)-ypplay(j, 1))
461               tairsol(j) = yts(j) + y_d_ts(j)               tairsol(j) = yts(j) + y_d_ts(j)
462               rugo1(j) = yrugos(j)               rugo1(j) = yrugos(j)
463               IF (nsrf == is_oce) THEN               IF (nsrf == is_oce) THEN
464                  rugo1(j) = rugos(i, nsrf)                  rugo1(j) = frugs(i, nsrf)
465               END IF               END IF
466               psfce(j) = ypaprs(j, 1)               psfce(j) = ypaprs(j, 1)
467               patm(j) = ypplay(j, 1)               patm(j) = ypplay(j, 1)
   
              qairsol(j) = yqsurf(j)  
468            END DO            END DO
469    
470            CALL stdlevvar(klon, knon, nsrf, zxli, uzon, vmer, tair1, qair1, &            CALL stdlevvar(nsrf, u1(:knon), v1(:knon), tair1(:knon), qair1, &
471                 zgeo1, tairsol, qairsol, rugo1, psfce, patm, yt2m, yq2m, &                 zgeo1, tairsol, yqsurf(:knon), rugo1, psfce, patm, yt2m, yq2m, &
472                 yt10m, yq10m, yu10m, yustar)                 yt10m, yq10m, wind10m(:knon), ustar(:knon))
473    
474            DO j = 1, knon            DO j = 1, knon
475               i = ni(j)               i = ni(j)
476               t2m(i, nsrf) = yt2m(j)               t2m(i, nsrf) = yt2m(j)
477               q2m(i, nsrf) = yq2m(j)               q2m(i, nsrf) = yq2m(j)
478    
479               ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman               u10m_srf(i, nsrf) = (wind10m(j) * u1(j)) &
480               u10m(i, nsrf) = (yu10m(j)*uzon(j))/sqrt(uzon(j)**2+vmer(j)**2)                    / sqrt(u1(j)**2 + v1(j)**2)
481               v10m(i, nsrf) = (yu10m(j)*vmer(j))/sqrt(uzon(j)**2+vmer(j)**2)               v10m_srf(i, nsrf) = (wind10m(j) * v1(j)) &
482                      / sqrt(u1(j)**2 + v1(j)**2)
483            END DO            END DO
484    
485            CALL hbtm(ypaprs, ypplay, yt2m, yq2m, yustar, y_flux_t(:knon), &            CALL hbtm(ypaprs, ypplay, yt2m, yq2m, ustar(:knon), y_flux_t(:knon), &
486                 y_flux_q(:knon), yu, yv, yt, yq, ypblh(:knon), ycapcl, &                 y_flux_q(:knon), yu(:knon, :), yv(:knon, :), yt(:knon, :), &
487                 yoliqcl, ycteicl, ypblt, ytherm, ytrmb1, ytrmb2, ytrmb3, ylcl)                 yq(:knon, :), ypblh(:knon), ycapcl, yoliqcl, ycteicl, ypblt, &
488                   ytherm, ylcl)
489    
490            DO j = 1, knon            DO j = 1, knon
491               i = ni(j)               i = ni(j)
# Line 596  contains Line 496  contains
496               cteicl(i, nsrf) = ycteicl(j)               cteicl(i, nsrf) = ycteicl(j)
497               pblt(i, nsrf) = ypblt(j)               pblt(i, nsrf) = ypblt(j)
498               therm(i, nsrf) = ytherm(j)               therm(i, nsrf) = ytherm(j)
              trmb1(i, nsrf) = ytrmb1(j)  
              trmb2(i, nsrf) = ytrmb2(j)  
              trmb3(i, nsrf) = ytrmb3(j)  
499            END DO            END DO
500    
501            DO j = 1, knon            IF (iflag_pbl >= 6) q2(ni(:knon), :, nsrf) = yq2(:knon, :)
502               DO k = 1, klev + 1         else
503                  i = ni(j)            fsnow(:, nsrf) = 0.
                 q2(i, k, nsrf) = yq2(j, k)  
              END DO  
           END DO  
504         end IF if_knon         end IF if_knon
505      END DO loop_surface      END DO loop_surface
506    
507      ! On utilise les nouvelles surfaces      ! On utilise les nouvelles surfaces
508      rugos(:, is_oce) = rugmer      frugs(:, is_oce) = rugmer
509      pctsrf(:, is_oce) = pctsrf_new_oce      pctsrf(:, is_oce) = pctsrf_new_oce
510      pctsrf(:, is_sic) = pctsrf_new_sic      pctsrf(:, is_sic) = pctsrf_new_sic
511    
512      firstcal = .false.      CALL histwrite_phy("run_off_lic", run_off_lic)
513        ftsol = ftsol + d_ts ! update surface temperature
514        CALL histwrite_phy("dtsvdfo", d_ts(:, is_oce))
515        CALL histwrite_phy("dtsvdft", d_ts(:, is_ter))
516        CALL histwrite_phy("dtsvdfg", d_ts(:, is_lic))
517        CALL histwrite_phy("dtsvdfi", d_ts(:, is_sic))
518    
519    END SUBROUTINE clmain    END SUBROUTINE pbl_surface
520    
521  end module clmain_m  end module pbl_surface_m

Legend:
Removed from v.207  
changed lines
  Added in v.330

  ViewVC Help
Powered by ViewVC 1.1.21