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

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

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

revision 175 by guez, Fri Feb 5 16:02:34 2016 UTC revision 228 by guez, Fri Nov 3 12:38:47 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         co2_ppm, 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
     USE dimens_m, ONLY: iim, jjm  
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):: co2_ppm ! taux CO2 atmosphere      REAL, INTENT(IN):: ftsol(:, :) ! (klon, nbsrf) temp\'erature du sol (en K)
     REAL, INTENT(IN):: ts(klon, nbsrf) ! temperature du sol (en Kelvin)  
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 65  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 101  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 à 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      !IM cf. AM : pbl, hbtm (Comme les autres diagnostics on cumule ds  
115      ! physiq ce qui permet de sortir les grdeurs par sous surface)      REAL, INTENT(inout):: u10m_srf(:, :), v10m_srf(:, :) ! (klon, nbsrf)
116      REAL pblh(klon, nbsrf)      ! composantes du vent \`a 10m sans spirale d'Ekman
117      ! pblh------- HCL  
118        ! 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 146  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)  
     ! -- LOOP  
     REAL yu10mx(klon)  
     REAL yu10my(klon)  
     REAL ywindsp(klon)  
     ! -- LOOP  
185    
186      REAL yt10m(klon), yq10m(klon)      REAL yt10m(klon), yq10m(klon)
187      REAL ypblh(klon)      REAL ypblh(klon)
# Line 229  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 262  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 282  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.
     yu10mx = 0.  
     yu10my = 0.  
     ywindsp = 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 308  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 331  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)
              yu10mx(j) = u10m(i, nsrf)  
              yu10my(j) = v10m(i, nsrf)  
              ywindsp(j) = sqrt(yu10mx(j)*yu10mx(j)+yu10my(j)*yu10my(j))  
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 379  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 401  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 430  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))
              IF (prt_level > 9) PRINT *, 'USTAR = ', yustar  
372    
373               ! iflag_pbl peut \^etre utilis\'e comme longueur de m\'elange               ! iflag_pbl peut \^etre utilis\'e comme longueur de m\'elange
374    
375               IF (iflag_pbl >= 11) THEN               CALL yamada4(dtime, rg, zlev(:knon, :), yzlay(:knon, :), &
376                  CALL vdif_kcay(knon, dtime, rg, ypaprs, yzlev, yzlay, yu, yv, &                    yu(:knon, :), yv(:knon, :), yteta(:knon, :), &
377                       yteta, coefm(:knon, 1), yq2, q2diag, ykmm, ykmn, yustar, &                    coefm(:knon, 1), yq2(:knon, :), ykmm(:knon, :), &
378                       iflag_pbl)                    ykmn(:knon, :), ykmq(:knon, :), ustar(:knon), 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  
379    
380               coefm(:knon, 2:) = ykmm(:knon, 2:klev)               coefm(:knon, 2:) = ykmm(:knon, 2:klev)
381               coefh(:knon, 2:) = ykmn(:knon, 2:klev)               coefh(:knon, 2:) = ykmn(:knon, 2:klev)
382            END IF            END IF
383    
384            ! calculer la diffusion des vitesses "u" et "v"            ! calculer la diffusion des vitesses "u" et "v"
385            CALL clvent(knon, dtime, yu1, yv1, coefm(:knon, :), yt, yu, ypaprs, &            CALL clvent(knon, dtime, yu(:knon, 1), yv(:knon, 1), &
386                 ypplay, ydelp, y_d_u, y_flux_u)                 coefm(:knon, :), yt, yu, ypaprs, ypplay, ydelp, y_d_u, &
387            CALL clvent(knon, dtime, yu1, yv1, coefm(:knon, :), yt, yv, ypaprs, &                 y_flux_u(:knon))
388                 ypplay, ydelp, y_d_v, y_flux_v)            CALL clvent(knon, dtime, yu(:knon, 1), yv(:knon, 1), &
389                   coefm(:knon, :), yt, yv, ypaprs, ypplay, ydelp, y_d_v, &
390                   y_flux_v(:knon))
391    
392            ! calculer la diffusion de "q" et de "h"            ! calculer la diffusion de "q" et de "h"
393            CALL clqh(dtime, itap, jour, debut, rlat, knon, nsrf, ni(:knon), &            CALL clqh(dtime, julien, firstcal, nsrf, ni(:knon), &
394                 pctsrf, ytsoil, yqsol, rmu0, co2_ppm, yrugos, yrugoro, yu1, &                 ytsoil(:knon, :), yqsol(:knon), mu0, yrugos, yrugoro, &
395                 yv1, coefh(:knon, :), yt, yq, yts, ypaprs, ypplay, ydelp, &                 yu(:knon, 1), yv(:knon, 1), coefh(:knon, :), yt, yq, &
396                 yrads, yalb(:knon), ysnow, yqsurf, yrain_f, ysnow_f, yfder, &                 yts(:knon), ypaprs, ypplay, ydelp, yrads(:knon), yalb(:knon), &
397                 yfluxlat, pctsrf_new, yagesno(:knon), y_d_t, y_d_q, &                 snow(:knon), yqsurf, yrain_f, ysnow_f, yfluxlat(:knon), &
398                 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), &
399                 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), &
400                   y_dflux_q(:knon), y_fqcalving, y_ffonte, y_run_off_lic_0)
401    
402            ! calculer la longueur de rugosite sur ocean            ! calculer la longueur de rugosite sur ocean
403            yrugm = 0.            yrugm = 0.
404            IF (nsrf == is_oce) THEN            IF (nsrf == is_oce) THEN
405               DO j = 1, knon               DO j = 1, knon
406                  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) &
407                       0.11*14E-6/sqrt(coefm(j, 1)*(yu1(j)**2+yv1(j)**2))                       / rg + 0.11 * 14E-6 &
408                         / sqrt(coefm(j, 1) * (yu(j, 1)**2 + yv(j, 1)**2))
409                  yrugm(j) = max(1.5E-05, yrugm(j))                  yrugm(j) = max(1.5E-05, yrugm(j))
410               END DO               END DO
411            END IF            END IF
412            DO j = 1, knon            DO j = 1, knon
413               y_dflux_t(j) = y_dflux_t(j)*ypct(j)               y_dflux_t(j) = y_dflux_t(j) * ypct(j)
414               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)  
415            END DO            END DO
416    
417            DO k = 1, klev            DO k = 1, klev
418               DO j = 1, knon               DO j = 1, knon
419                  i = ni(j)                  i = ni(j)
420                  coefh(j, k) = coefh(j, k)*ypct(j)                  coefh(j, k) = coefh(j, k) * ypct(j)
421                  coefm(j, k) = coefm(j, k)*ypct(j)                  coefm(j, k) = coefm(j, k) * ypct(j)
422                  y_d_t(j, k) = y_d_t(j, k)*ypct(j)                  y_d_t(j, k) = y_d_t(j, k) * ypct(j)
423                  y_d_q(j, k) = y_d_q(j, k)*ypct(j)                  y_d_q(j, k) = y_d_q(j, k) * ypct(j)
424                  flux_t(i, k, nsrf) = y_flux_t(j, k)                  y_d_u(j, k) = y_d_u(j, k) * ypct(j)
425                  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)  
426               END DO               END DO
427            END DO            END DO
428    
429            evap(:, nsrf) = -flux_q(:, 1, nsrf)            flux_t(ni(:knon), nsrf) = y_flux_t(:knon)
430              flux_q(ni(:knon), nsrf) = y_flux_q(:knon)
431              flux_u(ni(:knon), nsrf) = y_flux_u(:knon)
432              flux_v(ni(:knon), nsrf) = y_flux_v(:knon)
433    
434              evap(:, nsrf) = -flux_q(:, nsrf)
435    
436            falbe(:, nsrf) = 0.            falbe(:, nsrf) = 0.
437            snow(:, nsrf) = 0.            fsnow(:, nsrf) = 0.
438            qsurf(:, nsrf) = 0.            qsurf(:, nsrf) = 0.
439            rugos(:, nsrf) = 0.            frugs(:, nsrf) = 0.
           fluxlat(:, nsrf) = 0.  
440            DO j = 1, knon            DO j = 1, knon
441               i = ni(j)               i = ni(j)
442               d_ts(i, nsrf) = y_d_ts(j)               d_ts(i, nsrf) = y_d_ts(j)
443               falbe(i, nsrf) = yalb(j)               falbe(i, nsrf) = yalb(j)
444               snow(i, nsrf) = ysnow(j)               fsnow(i, nsrf) = snow(j)
445               qsurf(i, nsrf) = yqsurf(j)               qsurf(i, nsrf) = yqsurf(j)
446               rugos(i, nsrf) = yz0_new(j)               frugs(i, nsrf) = yz0_new(j)
447               fluxlat(i, nsrf) = yfluxlat(j)               fluxlat(i, nsrf) = yfluxlat(j)
448               IF (nsrf == is_oce) THEN               IF (nsrf == is_oce) THEN
449                  rugmer(i) = yrugm(j)                  rugmer(i) = yrugm(j)
450                  rugos(i, nsrf) = yrugm(j)                  frugs(i, nsrf) = yrugm(j)
451               END IF               END IF
452               agesno(i, nsrf) = yagesno(j)               agesno(i, nsrf) = yagesno(j)
453               fqcalving(i, nsrf) = y_fqcalving(j)               fqcalving(i, nsrf) = y_fqcalving(j)
# Line 521  contains Line 456  contains
456               cdragm(i) = cdragm(i) + coefm(j, 1)               cdragm(i) = cdragm(i) + coefm(j, 1)
457               dflux_t(i) = dflux_t(i) + y_dflux_t(j)               dflux_t(i) = dflux_t(i) + y_dflux_t(j)
458               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)  
459            END DO            END DO
460            IF (nsrf == is_ter) THEN            IF (nsrf == is_ter) THEN
461               qsol(ni(:knon)) = yqsol(:knon)               qsol(ni(:knon)) = yqsol(:knon)
# Line 534  contains Line 467  contains
467            END IF            END IF
468    
469            ftsoil(:, :, nsrf) = 0.            ftsoil(:, :, nsrf) = 0.
470            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  
471    
472            DO j = 1, knon            DO j = 1, knon
473               i = ni(j)               i = ni(j)
# Line 556  contains Line 484  contains
484    
485            DO j = 1, knon            DO j = 1, knon
486               i = ni(j)               i = ni(j)
487               uzon(j) = yu(j, 1) + y_d_u(j, 1)               u1(j) = yu(j, 1) + y_d_u(j, 1)
488               vmer(j) = yv(j, 1) + y_d_v(j, 1)               v1(j) = yv(j, 1) + y_d_v(j, 1)
489               tair1(j) = yt(j, 1) + y_d_t(j, 1)               tair1(j) = yt(j, 1) + y_d_t(j, 1)
490               qair1(j) = yq(j, 1) + y_d_q(j, 1)               qair1(j) = yq(j, 1) + y_d_q(j, 1)
491               zgeo1(j) = rd*tair1(j)/(0.5*(ypaprs(j, 1)+ypplay(j, &               zgeo1(j) = rd * tair1(j) / (0.5 * (ypaprs(j, 1) + ypplay(j, &
492                    1)))*(ypaprs(j, 1)-ypplay(j, 1))                    1))) * (ypaprs(j, 1)-ypplay(j, 1))
493               tairsol(j) = yts(j) + y_d_ts(j)               tairsol(j) = yts(j) + y_d_ts(j)
494               rugo1(j) = yrugos(j)               rugo1(j) = yrugos(j)
495               IF (nsrf == is_oce) THEN               IF (nsrf == is_oce) THEN
496                  rugo1(j) = rugos(i, nsrf)                  rugo1(j) = frugs(i, nsrf)
497               END IF               END IF
498               psfce(j) = ypaprs(j, 1)               psfce(j) = ypaprs(j, 1)
499               patm(j) = ypplay(j, 1)               patm(j) = ypplay(j, 1)
# Line 573  contains Line 501  contains
501               qairsol(j) = yqsurf(j)               qairsol(j) = yqsurf(j)
502            END DO            END DO
503    
504            CALL stdlevvar(klon, knon, nsrf, zxli, uzon, vmer, tair1, qair1, &            CALL stdlevvar(klon, knon, nsrf, u1(:knon), v1(:knon), tair1(:knon), &
505                 zgeo1, tairsol, qairsol, rugo1, psfce, patm, yt2m, yq2m, &                 qair1, zgeo1, tairsol, qairsol, rugo1, psfce, patm, yt2m, &
506                 yt10m, yq10m, yu10m, yustar)                 yq2m, yt10m, yq10m, wind10m(:knon), ustar)
507    
508            DO j = 1, knon            DO j = 1, knon
509               i = ni(j)               i = ni(j)
510               t2m(i, nsrf) = yt2m(j)               t2m(i, nsrf) = yt2m(j)
511               q2m(i, nsrf) = yq2m(j)               q2m(i, nsrf) = yq2m(j)
512    
513               ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman               u10m_srf(i, nsrf) = (wind10m(j) * u1(j)) &
514               u10m(i, nsrf) = (yu10m(j)*uzon(j))/sqrt(uzon(j)**2+vmer(j)**2)                    / sqrt(u1(j)**2 + v1(j)**2)
515               v10m(i, nsrf) = (yu10m(j)*vmer(j))/sqrt(uzon(j)**2+vmer(j)**2)               v10m_srf(i, nsrf) = (wind10m(j) * v1(j)) &
516                      / sqrt(u1(j)**2 + v1(j)**2)
517            END DO            END DO
518    
519            CALL hbtm(knon, ypaprs, ypplay, yt2m, yq2m, yustar, &            CALL hbtm(ypaprs, ypplay, yt2m, yq2m, ustar(:knon), y_flux_t(:knon), &
520                 y_flux_t, y_flux_q, yu, yv, yt, yq, ypblh, ycapcl, yoliqcl, &                 y_flux_q(:knon), yu, yv, yt, yq, ypblh(:knon), ycapcl, &
521                 ycteicl, ypblt, ytherm, ytrmb1, ytrmb2, ytrmb3, ylcl)                 yoliqcl, ycteicl, ypblt, ytherm, ytrmb1, ytrmb2, ytrmb3, ylcl)
522    
523            DO j = 1, knon            DO j = 1, knon
524               i = ni(j)               i = ni(j)
# Line 612  contains Line 540  contains
540                  q2(i, k, nsrf) = yq2(j, k)                  q2(i, k, nsrf) = yq2(j, k)
541               END DO               END DO
542            END DO            END DO
543           else
544              fsnow(:, nsrf) = 0.
545         end IF if_knon         end IF if_knon
546      END DO loop_surface      END DO loop_surface
547    
548      ! On utilise les nouvelles surfaces      ! On utilise les nouvelles surfaces
549        frugs(:, is_oce) = rugmer
550        pctsrf(:, is_oce) = pctsrf_new_oce
551        pctsrf(:, is_sic) = pctsrf_new_sic
552    
553      rugos(:, is_oce) = rugmer      firstcal = .false.
     pctsrf = pctsrf_new  
554    
555    END SUBROUTINE clmain    END SUBROUTINE clmain
556    
Legend:
Removed from v.175  
changed lines
  Added in v.228
  ViewVC Help
Powered by ViewVC 1.1.21