source: LMDZ6/trunk/libf/phylmd/inlandsis/sisvat_ts2.f90 @ 5284

Last change on this file since 5284 was 5284, checked in by abarral, 6 hours ago

Turn alpale.h alpale.f90 YOETHF.h into modules
File size: 19.2 KB
Line 
1subroutine SISVAT_TS2
2  ! #ES.                     (ETSo_0,ETSo_1,ETSo_d)
3
4  ! +------------------------------------------------------------------------+
5  ! | MAR          SISVAT_TS2                            Mon 16-08-2009  MAR |
6  ! |   SubRoutine SISVAT_TS2 computes the Soil/Snow temperature and fluxes  |
7  ! |   using the same method as in LMDZ for consistency.                    |
8  ! |   The corresponding LMDZ routines are soil (soil.F90) and calcul_fluxs |
9  ! |   (calcul_fluxs_mod.F90).                                              |
10  ! +------------------------------------------------------------------------+
11  ! |                                                                        |
12  ! |                                                                        |
13  ! |   PARAMETERS:  klonv: Total Number of columns =                        |
14  ! |   ^^^^^^^^^^        = Total Number of grid boxes of surface type       |
15  ! |                       (land ice for now)                               |
16  ! |                                                                        |
17  ! |   INPUT:   isnoSV   = total Nb of Ice/Snow Layers                      |
18  ! |   ^^^^^    sol_SV   : Downward Solar Radiation                  [W/m2] |
19  ! |            IRd_SV   : Surface Downward  Longwave   Radiation    [W/m2] |
20  ! |            VV__SV   : SBL Top    Wind Speed                      [m/s] |
21  ! |            TaT_SV   : SBL Top    Temperature                       [K] |
22  ! |            QaT_SV   : SBL Top    Specific  Humidity            [kg/kg] |
23  ! |            dzsnSV   : Snow Layer Thickness                         [m] |
24  ! |            dt__SV   : Time Step                                    [s] |
25  ! |                                                                        |
26  ! |            SoSosv   : Absorbed Solar Radiation by Surfac.(Normaliz)[-] |
27  ! |            Eso_sv   : Soil+Snow       Emissivity                   [-] |
28  ! |   ?        rah_sv   : Aerodynamic Resistance for Heat            [s/m] |
29  ! |                                                                        |
30  ! |            dz1_SV    : "inverse" layer thickness (centered)      [1/m] |
31  ! |            dz2_SV    : layer thickness (layer above (?))           [m] |
32  ! |            AcoHSV    : coefficient for Enthalpy evolution, from atm.   |
33  ! |            AcoHSV    : coefficient for Enthalpy evolution, from atm.   |
34  ! |            AcoQSV    : coefficient for Humidity evolution, from atm.   |
35  ! |            BcoQSV    : coefficient for Humidity evolution, from atm.   |
36  ! |            ps__SV    : surface pressure                           [Pa] |
37  ! |            p1l_SV    : 1st atmospheric layer pressure             [Pa] |
38  ! |            cdH_SV    : drag coeff Energy (?)                           |
39  ! |            rsolSV    : Radiation balance surface                [W/m2] |
40  ! |            lambSV    : Coefficient for soil layer geometry         [-] |
41  ! |                                                                        |
42  ! |   INPUT /  TsisSV   : Soil/Ice Temperatures (layers -nsol,-nsol+1,..,0)|
43  ! |   OUTPUT:           & Snow     Temperatures (layers  1,2,...,nsno) [K] |
44  ! |   ^^^^^^   rsolSV   : Radiation balance surface                 [W/m2] |
45  ! |                                                                        |
46  ! |   OUTPUT:  IRs_SV   : Soil      IR Radiation                    [W/m2] |
47  ! |   ^^^^^^   HSs_sv   : Sensible  Heat Flux                       [W/m2] |
48  ! |            HLs_sv   : Latent    Heat Flux                       [W/m2] |
49  ! |            TsfnSV   : new surface temperature                      [K] |
50  ! |            Evp_sv   : Evaporation                              [kg/m2] |
51  ! |            dSdTSV   : Sensible Heat Flux temp. derivative     [W/m2/K] |
52  ! |            dLdTSV   : Latent Heat Flux temp. derivative       [W/m2/K] |
53  ! |                                                                        |
54  ! |   ?        ETSo_0   : Snow/Soil Energy Power, before Forcing    [W/m2] |
55  ! |   ?        ETSo_1   : Snow/Soil Energy Power, after  Forcing    [W/m2] |
56  ! |   ?        ETSo_d   : Snow/Soil Energy Power         Forcing    [W/m2] |
57  ! |                                                                        |
58  ! |________________________________________________________________________|
59
60USE yoethf_mod_h
61    USE VAR_SV
62  USE VARdSV
63
64  USE VARySV
65  USE VARtSV
66  USE VARxSV
67  USE VARphy
68  USE indice_sol_mod
69
70
71  USE yomcst_mod_h, ONLY: RPI, RCLUM, RHPLA, RKBOL, RNAVO                   &
72          , RDAY, REA, REPSM, RSIYEA, RSIDAY, ROMEGA                  &
73          , R_ecc, R_peri, R_incl                                      &
74          , RA, RG, R1SA                                         &
75          , RSIGMA                                                     &
76          , R, RMD, RMV, RD, RV, RCPD                    &
77          , RMO3, RMCO2, RMC, RMCH4, RMN2O, RMCFC11, RMCFC12        &
78          , RCPV, RCVD, RCVV, RKAPPA, RETV, eps_w                    &
79          , RCW, RCS                                                 &
80          , RLVTT, RLSTT, RLMLT, RTT, RATM                           &
81          , RESTT, RALPW, RBETW, RGAMW, RALPS, RBETS, RGAMS            &
82          , RALPD, RBETD, RGAMD
83IMPLICIT NONE
84
85
86  ! +--Global Variables
87  ! +  ================
88
89
90  INCLUDE "FCTTRE.h"
91   ! INCLUDE "indicesol.h"
92  INCLUDE "comsoil.h"
93   ! include  "LMDZphy.inc"
94
95  ! +--OUTPUT for Stand Alone NetCDF File
96  ! +  ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
97  ! #NC      real*8        SOsoKL(klonv)             ! Absorbed Solar Radiation
98  ! #NC      real*8        IRsoKL(klonv)             ! Absorbed IR    Radiation
99  ! #NC      real*8        HSsoKL(klonv)             ! Absorbed Sensible Heat Flux
100  ! #NC      real*8        HLsoKL(klonv)             ! Absorbed Latent   Heat Flux
101  ! #NC      real*8        HLs_KL(klonv)             ! Evaporation
102  ! #NC      real*8        HLv_KL(klonv)             ! Transpiration
103  ! #NC      common/DumpNC/SOsoKL,IRsoKL
104  ! #NC     .             ,HSsoKL,HLsoKL
105  ! #NC     .             ,HLs_KL,HLv_KL
106
107  ! +--Internal Variables
108  ! +  ==================
109
110  integer :: ig,jk,isl
111  real :: mu
112  real :: Tsrf(klonv)               ! surface temperature as extrapolated from soil
113  real :: mug(klonv)                 !hj coef top layers
114  real :: ztherm_i(klonv),zdz2(klonv,-nsol:nsno),z1s
115  real :: pfluxgrd(klonv), pcapcal(klonv), cal(klonv)
116  real :: beta(klonv), dif_grnd(klonv)
117  real :: C_coef(klonv,-nsol:nsno),D_coef(klonv,-nsol:nsno)
118
119  REAL, DIMENSION(klonv)   :: zx_mh, zx_nh, zx_oh
120  REAL, DIMENSION(klonv)   :: zx_mq, zx_nq, zx_oq
121  REAL, DIMENSION(klonv)   :: zx_pkh, zx_dq_s_dt, zx_qsat, zx_coef
122  REAL, DIMENSION(klonv)   :: zx_sl, zx_k1
123  REAL, DIMENSION(klonv)   :: d_ts
124  REAL                     :: zdelta, zcvm5, zx_qs, zcor, zx_dq_s_dh
125  REAL                     :: qsat_new, q1_new
126   ! REAL, PARAMETER          :: t_grnd = 271.35, t_coup = 273.15
127   ! REAL, PARAMETER          :: max_eau_sol = 150.0
128  REAL, DIMENSION(klonv)   :: IRs__D, dIRsdT
129
130
131  REAL :: t_grnd                      ! not used
132  parameter(t_grnd = 271.35)       !
133  REAL :: t_coup                      ! distinguish evap/sublimation
134  parameter(t_coup = 273.15)       !
135  REAL :: max_eau_sol
136  parameter(max_eau_sol = 150.0)
137
138
139     ! write(*,*)'T check'
140  !
141  !    DO  ig = 1,knonv
142  !        DO  jk = 1,isnoSV(ig) !nsno
143  !          IF (TsisSV(ig,jk) <= 1.) THEN          !hj check
144  !            TsisSV(ig,jk) = TsisSV(ig,isnoSV(ig))
145  !          ENDIF
146  !
147  !          IF (TsisSV(ig,jk) <= 1.) THEN          !hj check
148  !            TsisSV(ig,jk) = 273.15
149  !          ENDIF
150  !        END DO
151  !        END DO
152
153  !!=======================================================================
154  !! I. First part: corresponds to soil.F90 in LMDZ
155  !!=======================================================================
156
157  DO ig = 1,knonv
158    DO jk =1,isnoSV(ig)
159      dz2_SV(ig,jk)=dzsnSV(ig,jk)
160  !! use arithmetic center between layers to derive dz1 for snow layers for simplicity:
161      dz1_SV(ig,jk)=2./(dzsnSV(ig,jk)+dzsnSV(ig,jk-1))
162    ENDDO
163  ENDDO
164
165  DO ig = 1,knonv
166    ztherm_i(ig)   = inertie_lic
167    IF (isnoSV(ig) > 0) ztherm_i(ig)   = inertie_sno
168  ENDDO
169
170  !!-----------------------------------------------------------------------
171  !! 1)
172  !! Calculation of Cgrf and Dgrd coefficients using soil temperature from
173  !! previous time step.
174  !!
175  !! These variables are recalculated on the local compressed grid instead
176  !! of saved in restart file.
177  !!-----------------------------------------------------------------------
178  DO ig=1,knonv
179    DO jk=-nsol,nsno
180      zdz2(ig,jk)=dz2_SV(ig,jk)/dt__SV                                !ptimestep
181    ENDDO
182  ENDDO
183
184  DO ig=1,knonv
185    z1s = zdz2(ig,-nsol)+dz1_SV(ig,-nsol+1)
186    C_coef(ig,-nsol+1)=zdz2(ig,-nsol)*TsisSV(ig,-nsol)/z1s
187    D_coef(ig,-nsol+1)=dz1_SV(ig,-nsol+1)/z1s
188  ENDDO
189
190  DO ig=1,knonv
191    DO jk=-nsol+1,isnoSV(ig)-1,1
192      z1s = 1./(zdz2(ig,jk)+dz1_SV(ig,jk+1)+dz1_SV(ig,jk)           &
193            *(1.-D_coef(ig,jk)))
194      C_coef(ig,jk+1)=                                              &
195            (TsisSV(ig,jk)*zdz2(ig,jk)                              &
196            +dz1_SV(ig,jk)*C_coef(ig,jk)) * z1s
197      D_coef(ig,jk+1)=dz1_SV(ig,jk+1)*z1s
198    ENDDO
199  ENDDO
200
201  !!-----------------------------------------------------------------------
202  !! 2)
203  !! Computation of the soil temperatures using the Cgrd and Dgrd
204  !! coefficient computed above
205  !!
206  !!-----------------------------------------------------------------------
207  !! Extrapolate surface Temperature                   !hj check
208  mu=1./((2.**1.5-1.)/(2.**(0.5)-1.)-1.)
209
210  ! IF (knonv>0) THEN
211  !  DO ig=1,8
212  !    write(*,*)ig,'sisvat: Tsis ',TsisSV(ig,isnoSV(ig))
213  !    write(*,*)'max-1            ',TsisSV(ig,isnoSV(ig)-1)
214  !    write(*,*)'max-2            ',TsisSV(ig,isnoSV(ig)-2)
215  !    write(*,*)'0                ',TsisSV(ig,0)
216  !!        write(*,*)min(max(isnoSV(ig),0),1),max(1-isnoSV(ig),0)
217  !  ENDDO
218  ! END IF
219
220  DO ig=1,knonv
221    IF (isnoSV(ig).GT.0) THEN
222      IF (isnoSV(ig).GT.1) THEN
223       mug(ig)=1./(1.+dzsnSV(ig,isnoSV(ig)-1)/dzsnSV(ig,isnoSV(ig))) !mu
224      ELSE
225       mug(ig) = 1./(1.+dzsnSV(ig,isnoSV(ig)-1)/dz_dSV(0)) !mu
226      ENDIF
227    ELSE
228      mug(ig) = lambSV
229    ENDIF
230
231    IF (mug(ig)  .LE. 0.05) THEN
232     write(*,*)'Attention mu low', mug(ig)
233    ENDIF
234    IF (mug(ig)  .GE. 0.98) THEN
235     write(*,*)'Attention mu high', mug(ig)
236    ENDIF
237
238    Tsrf(ig)=(1.5*TsisSV(ig,isnoSV(ig))-0.5*TsisSV(ig,isnoSV(ig)-1))&
239          *min(max(isnoSV(ig),0),1)+                             &
240          ((mug(ig)+1)*TsisSV(ig,0)-mug(ig)*TsisSV(ig,-1))       &
241          *max(1-isnoSV(ig),0)
242  ENDDO
243
244
245
246  !! Surface temperature
247  DO ig=1,knonv
248  TsisSV(ig,isnoSV(ig))=(mug(ig)*C_coef(ig,isnoSV(ig))+Tsf_SV(ig))/ &
249        (mug(ig)*(1.-D_coef(ig,isnoSV(ig)))+1.)
250  ENDDO
251
252  !! Other temperatures
253  DO ig=1,knonv
254    DO jk=isnoSV(ig),-nsol+1,-1
255      TsisSV(ig,jk-1)=C_coef(ig,jk)+D_coef(ig,jk)                   &
256            *TsisSV(ig,jk)
257    ENDDO
258  ENDDO
259   ! write(*,*)ig,'Tsis',TsisSV(ig,0)
260
261   ! IF (indice == is_sic) THEN
262   !   DO ig = 1,knonv
263   !     TsisSV(ig,-nsol) = RTT - 1.8
264   !   END DO
265   ! ENDIF
266
267  !C      !hj new 11 03 2010
268    DO ig=1,knonv
269      isl         = isnoSV(ig)
270       ! dIRsdT(ig) = Eso_sv(ig)* SteBo   * 4.                        & ! - d(IR)/d(T)
271  ! &                             * Tsf_SV(ig)                         & !T TsisSV(ig,isl)           !
272  ! &                             * Tsf_SV(ig)                         & !TsisSV(ig,isl)           !
273  ! &                             * Tsf_SV(ig) !TsisSV(ig,isl)           !
274  !      IRs__D(ig) = dIRsdT(ig)* Tsf_SV(ig) * 0.75 !TsisSV(ig,isl) * 0.75    !:
275      dIRsdT(ig) = Eso_sv(ig)* StefBo   * 4.                        & ! - d(IR)/d(T)
276            * TsisSV(ig,isl)                     &
277            * TsisSV(ig,isl)                     &
278            * TsisSV(ig,isl)
279      IRs__D(ig) = dIRsdT(ig)* TsisSV(ig,isl) * 0.75                  !:
280     END DO
281   ! !hj
282  !!-----------------------------------------------------------------------
283  !! 3)
284  !! Calculate the Cgrd and Dgrd coefficient corresponding to actual soil
285  !! temperature
286  !!-----------------------------------------------------------------------
287  DO ig=1,knonv
288    z1s = zdz2(ig,-nsol)+dz1_SV(ig,-nsol+1)
289    C_coef(ig,-nsol+1) = zdz2(ig,-nsol)*TsisSV(ig,-nsol)/z1s
290    D_coef(ig,-nsol+1) = dz1_SV(ig,-nsol+1)/z1s
291  ENDDO
292
293  DO ig=1,knonv
294    DO jk=-nsol+1,isnoSV(ig)-1,1
295      z1s = 1./(zdz2(ig,jk)+dz1_SV(ig,jk+1)+dz1_SV(ig,jk)           &
296            *(1.-D_coef(ig,jk)))
297      C_coef(ig,jk+1) = (TsisSV(ig,jk)*zdz2(ig,jk)+                 &
298            dz1_SV(ig,jk)*C_coef(ig,jk)) * z1s
299      D_coef(ig,jk+1) = dz1_SV(ig,jk+1)*z1s
300    ENDDO
301  ENDDO
302
303  !!-----------------------------------------------------------------------
304  !! 4)
305  !! Computation of the surface diffusive flux from ground and
306  !! calorific capacity of the ground
307  !!-----------------------------------------------------------------------
308  DO ig=1,knonv
309  !! (pfluxgrd)
310    pfluxgrd(ig) = ztherm_i(ig)*dz1_SV(ig,isnoSV(ig))*              &
311          (C_coef(ig,isnoSV(ig))+(D_coef(ig,isnoSV(ig))-1.)           &
312          *TsisSV(ig,isnoSV(ig)))
313  !! (pcapcal)
314    pcapcal(ig)  = ztherm_i(ig)*                                    &
315          (dz2_SV(ig,isnoSV(ig))+dt__SV*(1.-D_coef(ig,isnoSV(ig)))    &
316          *dz1_SV(ig,isnoSV(ig)))
317    z1s = mug(ig)*(1.-D_coef(ig,isnoSV(ig)))+1.
318    pcapcal(ig)  = pcapcal(ig)/z1s
319    pfluxgrd(ig) = ( pfluxgrd(ig)                                   &
320          + pcapcal(ig) * (TsisSV(ig,isnoSV(ig)) * z1s               &
321          - mug(ig)* C_coef(ig,isnoSV(ig))                           &
322          - Tsf_SV(ig))       /dt__SV )
323  ENDDO
324
325
326  cal(1:knonv) = RCPD / pcapcal(1:knonv)
327  rsolSV(1:knonv)  = rsolSV(1:knonv) + pfluxgrd(1:knonv)
328  !!=======================================================================
329  !! II. Second part: corresponds to calcul_fluxs_mod.F90 in LMDZ
330  !!=======================================================================
331
332  Evp_sv = 0.
333  ! #NC HSsoKL=0.
334  ! #NC HLsoKL=0.
335  dSdTSV = 0.
336  dLdTSV = 0.
337
338  beta(:) = 1.0
339  dif_grnd(:) = 0.0
340
341  !! zx_qs = qsat en kg/kg
342  !!**********************************************************************x***************
343
344  DO ig = 1,knonv
345    IF (ps__SV(ig).LT.1.) THEN
346       ! write(*,*)'ig',ig,'ps',ps__SV(ig)
347      ps__SV(ig)=max(ps__SV(ig),1.e-8)
348    ENDIF
349    IF (p1l_SV(ig).LT.1.) THEN
350       ! write(*,*)'ig',ig,'p1l',p1l_SV(ig)
351      p1l_SV(ig)=max(p1l_SV(ig),1.e-8)
352    ENDIF
353    IF (TaT_SV(ig).LT.180.) THEN
354       ! write(*,*)'ig',ig,'TaT',TaT_SV(ig)
355      TaT_SV(ig)=max(TaT_SV(ig),180.)
356    ENDIF
357    IF (QaT_SV(ig).LT.1.e-8) THEN
358       ! write(*,*)'ig',ig,'QaT',QaT_SV(ig)
359      QaT_SV(ig)=max(QaT_SV(ig),1.e-8)
360    ENDIF
361    IF (Tsf_SV(ig).LT.100.) THEN
362       ! write(*,*)'ig',ig,'Tsf',Tsf_SV(ig)
363      Tsf_SV(ig)=max(Tsf_SV(ig),180.)
364    ENDIF
365    IF (Tsf_SV(ig).GT.500.) THEN
366       ! write(*,*)'ig',ig,'Tsf',Tsf_SV(ig)
367      Tsf_SV(ig)=min(Tsf_SV(ig),400.)
368    ENDIF
369      ! IF (Tsrf(ig).LT.1.) THEN
370  !!          write(*,*)'ig',ig,'Tsrf',Tsrf(ig)
371      !   Tsrf(ig)=max(Tsrf(ig),TaT_SV(ig)-20.)
372      ! ENDIF
373     IF (cdH_SV(ig).LT.1.e-10) THEN
374        ! IF (ig.le.3)   write(*,*)'ig',ig,'cdH',cdH_SV(ig)
375       cdH_SV(ig)=.5
376     ENDIF
377  ENDDO
378
379
380  DO ig = 1,knonv
381    zx_pkh(ig) = 1. ! (ps__SV(ig)/ps__SV(ig))**RKAPPA
382    IF (thermcep) THEN
383      zdelta=MAX(0.,SIGN(1.,rtt-Tsf_SV(ig)))
384      zcvm5 = R5LES*LhvH2O*(1.-zdelta) + R5IES*LhsH2O*zdelta
385      zcvm5 = zcvm5 / RCPD / (1.0+RVTMP2*QaT_SV(ig))
386      zx_qs= r2es * FOEEW(Tsf_SV(ig),zdelta)/ps__SV(ig)
387      zx_qs=MIN(0.5,zx_qs)
388      ! !write(*,*)'zcor',retv*zx_qs
389      zcor=1./(1.-retv*zx_qs)
390      zx_qs=zx_qs*zcor
391      zx_dq_s_dh = FOEDE(Tsf_SV(ig),zdelta,zcvm5,zx_qs,zcor)        &
392            /LhvH2O / zx_pkh(ig)
393    ELSE
394      IF (Tsf_SV(ig).LT.t_coup) THEN
395         zx_qs = qsats(Tsf_SV(ig)) / ps__SV(ig)
396         zx_dq_s_dh = dqsats(Tsf_SV(ig),zx_qs)/LhvH2O               &
397               / zx_pkh(ig)
398      ELSE
399         zx_qs = qsatl(Tsf_SV(ig)) / ps__SV(ig)
400         zx_dq_s_dh = dqsatl(Tsf_SV(ig),zx_qs)/LhvH2O               &
401               / zx_pkh(ig)
402      ENDIF
403    ENDIF
404    zx_dq_s_dt(ig) = RCPD * zx_pkh(ig) * zx_dq_s_dh
405    zx_qsat(ig) = zx_qs
406     ! zx_coef(ig) = cdH_SV(ig) *                                     &
407  ! &       (1.0+SQRT(u1lay(ig)**2+v1lay(ig)**2)) *                   &
408  ! &       p1l_SV(ig)/(RD*t1lay(ig))
409    zx_coef(ig) = cdH_SV(ig) *                                      &
410          (1.0+VV__SV(ig)) *                                         &
411          p1l_SV(ig)/(RD*TaT_SV(ig))
412
413  ENDDO
414
415
416  !! === Calcul de la temperature de surface ===
417  !! zx_sl = chaleur latente d'evaporation ou de sublimation
418  !!****************************************************************************
419
420  DO ig = 1,knonv
421    zx_sl(ig) = LhvH2O
422    IF (Tsf_SV(ig) .LT. RTT) zx_sl(ig) = LhsH2O
423    zx_k1(ig) = zx_coef(ig)
424  ENDDO
425
426
427  DO ig = 1,knonv
428  !! Q
429    zx_oq(ig) = 1. - (beta(ig) * zx_k1(ig) * BcoQSV(ig) * dt__SV)
430    zx_mq(ig) = beta(ig) * zx_k1(ig) *                              &
431          (AcoQSV(ig) - zx_qsat(ig) +                                &
432          zx_dq_s_dt(ig) * Tsf_SV(ig))                               &
433          / zx_oq(ig)
434    zx_nq(ig) = beta(ig) * zx_k1(ig) * (-1. * zx_dq_s_dt(ig))       &
435          / zx_oq(ig)
436
437  !! H
438    zx_oh(ig) = 1. - (zx_k1(ig) * BcoHSV(ig) * dt__SV)
439    zx_mh(ig) = zx_k1(ig) * AcoHSV(ig) / zx_oh(ig)
440    zx_nh(ig) = - (zx_k1(ig) * RCPD * zx_pkh(ig))/ zx_oh(ig)
441
442  !! surface temperature
443    TsfnSV(ig) = (Tsf_SV(ig) + cal(ig)/RCPD * zx_pkh(ig) * dt__SV * &
444          (rsolSV(ig) + zx_mh(ig) + zx_sl(ig) * zx_mq(ig))           &
445          + dif_grnd(ig) * t_grnd * dt__SV)/                         &
446          ( 1. - dt__SV * cal(ig)/(RCPD * zx_pkh(ig)) *              &
447          (zx_nh(ig) + zx_sl(ig) * zx_nq(ig))                        &
448          + dt__SV * dif_grnd(ig))
449
450  !hj rajoute 22 11 2010 tuning...
451    TsfnSV(ig) = min(RTT+0.02,TsfnSV(ig))
452
453    d_ts(ig) = TsfnSV(ig) - Tsf_SV(ig)
454
455
456  !!== flux_q est le flux de vapeur d'eau: kg/(m**2 s)  positive vers bas
457  !!== flux_t est le flux de cpt (energie sensible): j/(m**2 s)
458    Evp_sv(ig) = - zx_mq(ig) - zx_nq(ig) * TsfnSV(ig)
459    HLs_sv(ig) = - Evp_sv(ig) * zx_sl(ig)
460    HSs_sv(ig) = zx_mh(ig) + zx_nh(ig) * TsfnSV(ig)
461
462  !! Derives des flux dF/dTs (W m-2 K-1):
463    dSdTSV(ig) = zx_nh(ig)
464    dLdTSV(ig) = zx_sl(ig) * zx_nq(ig)
465
466
467  !hj  new 11 03 2010
468    isl         = isnoSV(ig)
469     ! TsisSV(ig,isl) = TsfnSV(ig)
470    IRs_SV(ig) = IRs__D(ig)                                         & !
471          - dIRsdT(ig) * TsfnSV(ig) !TsisSV(ig,isl)?      !
472
473  ! hj
474  ! #NC   SOsoKL(ig) = sol_SV(ig) * SoSosv(ig)              ! Absorbed Sol.
475  ! #NC   IRsoKL(ig) =               IRs_SV(ig)                           & !Up Surf. IR
476  ! #NC&        +     tau_sv(ig)      *IRd_SV(ig)*Eso_sv(ig)              & !Down Atm IR
477  ! #NC&        -(1.0-tau_sv(ig)) *0.5*IRv_sv(ig)            ! Down Veg IR
478  ! #NC   HLsoKL(ig) = HLs_sv(ig)
479  ! #NC   HSsoKL(ig) = HSs_sv(ig)
480  ! #NC   HLs_KL(ig) = Evp_sv(ig)
481
482  !! Nouvelle valeure de l'humidite au dessus du sol
483    qsat_new=zx_qsat(ig) + zx_dq_s_dt(ig) * d_ts(ig)
484    q1_new = AcoQSV(ig) - BcoQSV(ig)* Evp_sv(ig)*dt__SV
485    QaT_SV(ig)=q1_new*(1.-beta(ig)) + beta(ig)*qsat_new
486
487  ENDDO
488
489end subroutine sisvat_ts2
Note: See TracBrowser for help on using the repository browser.