/[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

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

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

  ViewVC Help
Powered by ViewVC 1.1.21