source: LMDZ6/branches/Amaury_dev/libf/phylmd/inlandsis/sisvat_ts2.f90 @ 5153

Last change on this file since 5153 was 5153, checked in by abarral, 7 weeks ago

Revert FCTTRE to INCLUDE to assess impact of inlining
File size: 18.6 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
60  USE VAR_SV
61  USE VARdSV
62
63  USE VARySV
64  USE VARtSV
65  USE VARxSV
66  USE VARphy
67  USE indice_sol_mod
68  USE lmdz_comsoil, ONLY: inertie_sol, inertie_sno, inertie_sic, inertie_lic, iflag_sic, iflag_inertie
69  USE lmdz_yoethf
70
71  USE lmdz_yomcst
72
73  IMPLICIT NONE
74 INCLUDE "FCTTRE.h"
75
76
77  ! +--Global Variables
78  ! +  ================
79
80  ! INCLUDE "indicesol.h"
81  ! include  "LMDZphy.inc"
82
83  ! +--OUTPUT for Stand Alone NetCDF File
84  ! +  ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
85  ! #NC      real*8        SOsoKL(klonv)             ! Absorbed Solar Radiation
86  ! #NC      real*8        IRsoKL(klonv)             ! Absorbed IR    Radiation
87  ! #NC      real*8        HSsoKL(klonv)             ! Absorbed Sensible Heat Flux
88  ! #NC      real*8        HLsoKL(klonv)             ! Absorbed Latent   Heat Flux
89  ! #NC      real*8        HLs_KL(klonv)             ! Evaporation
90  ! #NC      real*8        HLv_KL(klonv)             ! Transpiration
91  ! #NC      common/DumpNC/SOsoKL,IRsoKL
92  ! #NC     .             ,HSsoKL,HLsoKL
93  ! #NC     .             ,HLs_KL,HLv_KL
94
95  ! +--Internal Variables
96  ! +  ==================
97
98  INTEGER :: ig, jk, isl
99  REAL :: mu
100  REAL :: Tsrf(klonv)               ! surface temperature as extrapolated from soil
101  REAL :: mug(klonv)                 !hj coef top layers
102  REAL :: ztherm_i(klonv), zdz2(klonv, -nsol:nsno), z1s
103  REAL :: pfluxgrd(klonv), pcapcal(klonv), cal(klonv)
104  REAL :: beta(klonv), dif_grnd(klonv)
105  REAL :: C_coef(klonv, -nsol:nsno), D_coef(klonv, -nsol:nsno)
106
107  REAL, DIMENSION(klonv) :: zx_mh, zx_nh, zx_oh
108  REAL, DIMENSION(klonv) :: zx_mq, zx_nq, zx_oq
109  REAL, DIMENSION(klonv) :: zx_pkh, zx_dq_s_dt, zx_qsat, zx_coef
110  REAL, DIMENSION(klonv) :: zx_sl, zx_k1
111  REAL, DIMENSION(klonv) :: d_ts
112  REAL :: zdelta, zcvm5, zx_qs, zcor, zx_dq_s_dh
113  REAL :: qsat_new, q1_new
114  ! REAL, PARAMETER          :: t_grnd = 271.35, t_coup = 273.15
115  ! REAL, PARAMETER          :: max_eau_sol = 150.0
116  REAL, DIMENSION(klonv) :: IRs__D, dIRsdT
117
118  REAL :: t_grnd                      ! not used
119  parameter(t_grnd = 271.35)       !
120  REAL :: t_coup                      ! distinguish evap/sublimation
121  parameter(t_coup = 273.15)       !
122  REAL :: max_eau_sol
123  parameter(max_eau_sol = 150.0)
124
125
126  ! WRITE(*,*)'T check'
127
128  ! DO  ig = 1,knonv
129  !     DO  jk = 1,isnoSV(ig) !nsno
130  !       IF (TsisSV(ig,jk) <= 1.) THEN          !hj check
131  !         TsisSV(ig,jk) = TsisSV(ig,isnoSV(ig))
132  !       ENDIF
133
134  !       IF (TsisSV(ig,jk) <= 1.) THEN          !hj check
135  !         TsisSV(ig,jk) = 273.15
136  !       ENDIF
137  !     END DO
138  !        END DO
139
140  !!=======================================================================
141  !! I. First part: corresponds to soil.F90 in LMDZ
142  !!=======================================================================
143
144  DO ig = 1, knonv
145    DO jk = 1, isnoSV(ig)
146      dz2_SV(ig, jk) = dzsnSV(ig, jk)
147      !! use arithmetic center between layers to derive dz1 for snow layers for simplicity:
148      dz1_SV(ig, jk) = 2. / (dzsnSV(ig, jk) + dzsnSV(ig, jk - 1))
149    ENDDO
150  ENDDO
151
152  DO ig = 1, knonv
153    ztherm_i(ig) = inertie_lic
154    IF (isnoSV(ig) > 0) ztherm_i(ig) = inertie_sno
155  ENDDO
156
157  !!-----------------------------------------------------------------------
158  !! 1)
159  !! Calculation of Cgrf and Dgrd coefficients using soil temperature from
160  !! previous time step.
161  !!
162  !! These variables are recalculated on the local compressed grid instead
163  !! of saved in restart file.
164  !!-----------------------------------------------------------------------
165  DO ig = 1, knonv
166    DO jk = -nsol, nsno
167      zdz2(ig, jk) = dz2_SV(ig, jk) / dt__SV                                !ptimestep
168    ENDDO
169  ENDDO
170
171  DO ig = 1, knonv
172    z1s = zdz2(ig, -nsol) + dz1_SV(ig, -nsol + 1)
173    C_coef(ig, -nsol + 1) = zdz2(ig, -nsol) * TsisSV(ig, -nsol) / z1s
174    D_coef(ig, -nsol + 1) = dz1_SV(ig, -nsol + 1) / z1s
175  ENDDO
176
177  DO ig = 1, knonv
178    DO jk = -nsol + 1, isnoSV(ig) - 1, 1
179      z1s = 1. / (zdz2(ig, jk) + dz1_SV(ig, jk + 1) + dz1_SV(ig, jk)           &
180              * (1. - D_coef(ig, jk)))
181      C_coef(ig, jk + 1) = &
182              (TsisSV(ig, jk) * zdz2(ig, jk)                              &
183                      + dz1_SV(ig, jk) * C_coef(ig, jk)) * z1s
184      D_coef(ig, jk + 1) = dz1_SV(ig, jk + 1) * z1s
185    ENDDO
186  ENDDO
187
188  !!-----------------------------------------------------------------------
189  !! 2)
190  !! Computation of the soil temperatures using the Cgrd and Dgrd
191  !! coefficient computed above
192  !!
193  !!-----------------------------------------------------------------------
194  !! Extrapolate surface Temperature                   !hj check
195  mu = 1. / ((2.**1.5 - 1.) / (2.**(0.5) - 1.) - 1.)
196
197  ! IF (knonv>0) THEN
198  !  DO ig=1,8
199  !    WRITE(*,*)ig,'sisvat: Tsis ',TsisSV(ig,isnoSV(ig))
200  !    WRITE(*,*)'max-1            ',TsisSV(ig,isnoSV(ig)-1)
201  !    WRITE(*,*)'max-2            ',TsisSV(ig,isnoSV(ig)-2)
202  !    WRITE(*,*)'0                ',TsisSV(ig,0)
203  !!        WRITE(*,*)min(max(isnoSV(ig),0),1),max(1-isnoSV(ig),0)
204  !  ENDDO
205  ! END IF
206
207  DO ig = 1, knonv
208    IF (isnoSV(ig)>0) THEN
209      IF (isnoSV(ig)>1) THEN
210        mug(ig) = 1. / (1. + dzsnSV(ig, isnoSV(ig) - 1) / dzsnSV(ig, isnoSV(ig))) !mu
211      ELSE
212        mug(ig) = 1. / (1. + dzsnSV(ig, isnoSV(ig) - 1) / dz_dSV(0)) !mu
213      ENDIF
214    ELSE
215      mug(ig) = lambSV
216    ENDIF
217
218    IF (mug(ig)  <= 0.05) THEN
219      WRITE(*, *)'Attention mu low', mug(ig)
220    ENDIF
221    IF (mug(ig)  >= 0.98) THEN
222      WRITE(*, *)'Attention mu high', mug(ig)
223    ENDIF
224
225    Tsrf(ig) = (1.5 * TsisSV(ig, isnoSV(ig)) - 0.5 * TsisSV(ig, isnoSV(ig) - 1))&
226            * min(max(isnoSV(ig), 0), 1) + &
227            ((mug(ig) + 1) * TsisSV(ig, 0) - mug(ig) * TsisSV(ig, -1))       &
228                    * max(1 - isnoSV(ig), 0)
229  ENDDO
230
231
232
233  !! Surface temperature
234  DO ig = 1, knonv
235    TsisSV(ig, isnoSV(ig)) = (mug(ig) * C_coef(ig, isnoSV(ig)) + Tsf_SV(ig)) / &
236            (mug(ig) * (1. - D_coef(ig, isnoSV(ig))) + 1.)
237  ENDDO
238
239  !! Other temperatures
240  DO ig = 1, knonv
241    DO jk = isnoSV(ig), -nsol + 1, -1
242      TsisSV(ig, jk - 1) = C_coef(ig, jk) + D_coef(ig, jk)                   &
243              * TsisSV(ig, jk)
244    ENDDO
245  ENDDO
246  ! WRITE(*,*)ig,'Tsis',TsisSV(ig,0)
247
248  ! IF (indice == is_sic) THEN
249  !   DO ig = 1,knonv
250  !     TsisSV(ig,-nsol) = RTT - 1.8
251  !   END DO
252  ! ENDIF
253
254  !C      !hj new 11 03 2010
255  DO ig = 1, knonv
256    isl = isnoSV(ig)
257    ! dIRsdT(ig) = Eso_sv(ig)* SteBo   * 4.                        & ! - d(IR)/d(T)
258    ! &                             * Tsf_SV(ig)                         & !T TsisSV(ig,isl)           !
259    ! &                             * Tsf_SV(ig)                         & !TsisSV(ig,isl)           !
260    ! &                             * Tsf_SV(ig) !TsisSV(ig,isl)           !
261    !      IRs__D(ig) = dIRsdT(ig)* Tsf_SV(ig) * 0.75 !TsisSV(ig,isl) * 0.75    !:
262    dIRsdT(ig) = Eso_sv(ig) * StefBo * 4.                        & ! - d(IR)/d(T)
263            * TsisSV(ig, isl)                     & !
264            * TsisSV(ig, isl)                     & !
265            * TsisSV(ig, isl)
266    IRs__D(ig) = dIRsdT(ig) * TsisSV(ig, isl) * 0.75                  !:
267  END DO
268  !hj
269  !!-----------------------------------------------------------------------
270  !! 3)
271  !! Calculate the Cgrd and Dgrd coefficient corresponding to actual soil
272  !! temperature
273  !!-----------------------------------------------------------------------
274  DO ig = 1, knonv
275    z1s = zdz2(ig, -nsol) + dz1_SV(ig, -nsol + 1)
276    C_coef(ig, -nsol + 1) = zdz2(ig, -nsol) * TsisSV(ig, -nsol) / z1s
277    D_coef(ig, -nsol + 1) = dz1_SV(ig, -nsol + 1) / z1s
278  ENDDO
279
280  DO ig = 1, knonv
281    DO jk = -nsol + 1, isnoSV(ig) - 1, 1
282      z1s = 1. / (zdz2(ig, jk) + dz1_SV(ig, jk + 1) + dz1_SV(ig, jk)           &
283              * (1. - D_coef(ig, jk)))
284      C_coef(ig, jk + 1) = (TsisSV(ig, jk) * zdz2(ig, jk) + &
285              dz1_SV(ig, jk) * C_coef(ig, jk)) * z1s
286      D_coef(ig, jk + 1) = dz1_SV(ig, jk + 1) * z1s
287    ENDDO
288  ENDDO
289
290  !!-----------------------------------------------------------------------
291  !! 4)
292  !! Computation of the surface diffusive flux from ground and
293  !! calorific capacity of the ground
294  !!-----------------------------------------------------------------------
295  DO ig = 1, knonv
296    !! (pfluxgrd)
297    pfluxgrd(ig) = ztherm_i(ig) * dz1_SV(ig, isnoSV(ig)) * &
298            (C_coef(ig, isnoSV(ig)) + (D_coef(ig, isnoSV(ig)) - 1.)           &
299                    * TsisSV(ig, isnoSV(ig)))
300    !! (pcapcal)
301    pcapcal(ig) = ztherm_i(ig) * &
302            (dz2_SV(ig, isnoSV(ig)) + dt__SV * (1. - D_coef(ig, isnoSV(ig)))    &
303                    * dz1_SV(ig, isnoSV(ig)))
304    z1s = mug(ig) * (1. - D_coef(ig, isnoSV(ig))) + 1.
305    pcapcal(ig) = pcapcal(ig) / z1s
306    pfluxgrd(ig) = (pfluxgrd(ig)                                   &
307            + pcapcal(ig) * (TsisSV(ig, isnoSV(ig)) * z1s               &
308                    - mug(ig) * C_coef(ig, isnoSV(ig))                           &
309                    - Tsf_SV(ig)) / dt__SV)
310  ENDDO
311
312  cal(1:knonv) = RCPD / pcapcal(1:knonv)
313  rsolSV(1:knonv) = rsolSV(1:knonv) + pfluxgrd(1:knonv)
314  !!=======================================================================
315  !! II. Second part: corresponds to calcul_fluxs_mod.F90 in LMDZ
316  !!=======================================================================
317
318  Evp_sv = 0.
319  ! #NC HSsoKL=0.
320  ! #NC HLsoKL=0.
321  dSdTSV = 0.
322  dLdTSV = 0.
323
324  beta(:) = 1.0
325  dif_grnd(:) = 0.0
326
327  !! zx_qs = qsat en kg/kg
328  !!**********************************************************************x***************
329
330  DO ig = 1, knonv
331    IF (ps__SV(ig)<1.) THEN
332      ! WRITE(*,*)'ig',ig,'ps',ps__SV(ig)
333      ps__SV(ig) = max(ps__SV(ig), 1.e-8)
334    ENDIF
335    IF (p1l_SV(ig)<1.) THEN
336      ! WRITE(*,*)'ig',ig,'p1l',p1l_SV(ig)
337      p1l_SV(ig) = max(p1l_SV(ig), 1.e-8)
338    ENDIF
339    IF (TaT_SV(ig)<180.) THEN
340      ! WRITE(*,*)'ig',ig,'TaT',TaT_SV(ig)
341      TaT_SV(ig) = max(TaT_SV(ig), 180.)
342    ENDIF
343    IF (QaT_SV(ig)<1.e-8) THEN
344      ! WRITE(*,*)'ig',ig,'QaT',QaT_SV(ig)
345      QaT_SV(ig) = max(QaT_SV(ig), 1.e-8)
346    ENDIF
347    IF (Tsf_SV(ig)<100.) THEN
348      ! WRITE(*,*)'ig',ig,'Tsf',Tsf_SV(ig)
349      Tsf_SV(ig) = max(Tsf_SV(ig), 180.)
350    ENDIF
351    IF (Tsf_SV(ig)>500.) THEN
352      ! WRITE(*,*)'ig',ig,'Tsf',Tsf_SV(ig)
353      Tsf_SV(ig) = min(Tsf_SV(ig), 400.)
354    ENDIF
355    ! IF (Tsrf(ig).LT.1.) THEN
356    !!          WRITE(*,*)'ig',ig,'Tsrf',Tsrf(ig)
357    !   Tsrf(ig)=max(Tsrf(ig),TaT_SV(ig)-20.)
358    ! ENDIF
359    IF (cdH_SV(ig)<1.e-10) THEN
360      ! IF (ig.le.3)   WRITE(*,*)'ig',ig,'cdH',cdH_SV(ig)
361      cdH_SV(ig) = .5
362    ENDIF
363  ENDDO
364
365  DO ig = 1, knonv
366    zx_pkh(ig) = 1. ! (ps__SV(ig)/ps__SV(ig))**RKAPPA
367    IF (thermcep) THEN
368      zdelta = MAX(0., SIGN(1., rtt - Tsf_SV(ig)))
369      zcvm5 = R5LES * LhvH2O * (1. - zdelta) + R5IES * LhsH2O * zdelta
370      zcvm5 = zcvm5 / RCPD / (1.0 + RVTMP2 * QaT_SV(ig))
371      zx_qs = r2es * FOEEW(Tsf_SV(ig), zdelta) / ps__SV(ig)
372      zx_qs = MIN(0.5, zx_qs)
373      !WRITE(*,*)'zcor',retv*zx_qs
374      zcor = 1. / (1. - retv * zx_qs)
375      zx_qs = zx_qs * zcor
376      zx_dq_s_dh = FOEDE(Tsf_SV(ig), zdelta, zcvm5, zx_qs, zcor)        &
377              / LhvH2O / zx_pkh(ig)
378    ELSE
379      IF (Tsf_SV(ig)<t_coup) THEN
380        zx_qs = qsats(Tsf_SV(ig)) / ps__SV(ig)
381        zx_dq_s_dh = dqsats(Tsf_SV(ig), zx_qs) / LhvH2O               &
382                / zx_pkh(ig)
383      ELSE
384        zx_qs = qsatl(Tsf_SV(ig)) / ps__SV(ig)
385        zx_dq_s_dh = dqsatl(Tsf_SV(ig), zx_qs) / LhvH2O               &
386                / zx_pkh(ig)
387      ENDIF
388    ENDIF
389    zx_dq_s_dt(ig) = RCPD * zx_pkh(ig) * zx_dq_s_dh
390    zx_qsat(ig) = zx_qs
391    ! zx_coef(ig) = cdH_SV(ig) *                                     &
392    ! &       (1.0+SQRT(u1lay(ig)**2+v1lay(ig)**2)) *                   &
393    ! &       p1l_SV(ig)/(RD*t1lay(ig))
394    zx_coef(ig) = cdH_SV(ig) * &
395            (1.0 + VV__SV(ig)) * &
396            p1l_SV(ig) / (RD * TaT_SV(ig))
397
398  ENDDO
399
400
401  !! === Calcul de la temperature de surface ===
402  !! zx_sl = chaleur latente d'evaporation ou de sublimation
403  !!****************************************************************************
404
405  DO ig = 1, knonv
406    zx_sl(ig) = LhvH2O
407    IF (Tsf_SV(ig) < RTT) zx_sl(ig) = LhsH2O
408    zx_k1(ig) = zx_coef(ig)
409  ENDDO
410
411  DO ig = 1, knonv
412    !! Q
413    zx_oq(ig) = 1. - (beta(ig) * zx_k1(ig) * BcoQSV(ig) * dt__SV)
414    zx_mq(ig) = beta(ig) * zx_k1(ig) * &
415            (AcoQSV(ig) - zx_qsat(ig) + &
416                    zx_dq_s_dt(ig) * Tsf_SV(ig))                               &
417            / zx_oq(ig)
418    zx_nq(ig) = beta(ig) * zx_k1(ig) * (-1. * zx_dq_s_dt(ig))       &
419            / zx_oq(ig)
420
421    !! H
422    zx_oh(ig) = 1. - (zx_k1(ig) * BcoHSV(ig) * dt__SV)
423    zx_mh(ig) = zx_k1(ig) * AcoHSV(ig) / zx_oh(ig)
424    zx_nh(ig) = - (zx_k1(ig) * RCPD * zx_pkh(ig)) / zx_oh(ig)
425
426    !! surface temperature
427    TsfnSV(ig) = (Tsf_SV(ig) + cal(ig) / RCPD * zx_pkh(ig) * dt__SV * &
428            (rsolSV(ig) + zx_mh(ig) + zx_sl(ig) * zx_mq(ig))           &
429            + dif_grnd(ig) * t_grnd * dt__SV) / &
430            (1. - dt__SV * cal(ig) / (RCPD * zx_pkh(ig)) * &
431                    (zx_nh(ig) + zx_sl(ig) * zx_nq(ig))                        &
432                    + dt__SV * dif_grnd(ig))
433
434    !hj rajoute 22 11 2010 tuning...
435    TsfnSV(ig) = min(RTT + 0.02, TsfnSV(ig))
436
437    d_ts(ig) = TsfnSV(ig) - Tsf_SV(ig)
438
439
440    !!== flux_q est le flux de vapeur d'eau: kg/(m**2 s)  positive vers bas
441    !!== flux_t est le flux de cpt (energie sensible): j/(m**2 s)
442    Evp_sv(ig) = - zx_mq(ig) - zx_nq(ig) * TsfnSV(ig)
443    HLs_sv(ig) = - Evp_sv(ig) * zx_sl(ig)
444    HSs_sv(ig) = zx_mh(ig) + zx_nh(ig) * TsfnSV(ig)
445
446    !! Derives des flux dF/dTs (W m-2 K-1):
447    dSdTSV(ig) = zx_nh(ig)
448    dLdTSV(ig) = zx_sl(ig) * zx_nq(ig)
449
450
451    !hj  new 11 03 2010
452    isl = isnoSV(ig)
453    ! TsisSV(ig,isl) = TsfnSV(ig)
454    IRs_SV(ig) = IRs__D(ig)                                         & !
455            - dIRsdT(ig) * TsfnSV(ig) !TsisSV(ig,isl)?      !
456
457    ! hj
458    ! #NC   SOsoKL(ig) = sol_SV(ig) * SoSosv(ig)              ! Absorbed Sol.
459    ! #NC   IRsoKL(ig) =               IRs_SV(ig)                           & !Up Surf. IR
460    ! #NC&        +     tau_sv(ig)      *IRd_SV(ig)*Eso_sv(ig)              & !Down Atm IR
461    ! #NC&        -(1.0-tau_sv(ig)) *0.5*IRv_sv(ig)            ! Down Veg IR
462    ! #NC   HLsoKL(ig) = HLs_sv(ig)
463    ! #NC   HSsoKL(ig) = HSs_sv(ig)
464    ! #NC   HLs_KL(ig) = Evp_sv(ig)
465
466    !! Nouvelle valeure de l'humidite au dessus du sol
467    qsat_new = zx_qsat(ig) + zx_dq_s_dt(ig) * d_ts(ig)
468    q1_new = AcoQSV(ig) - BcoQSV(ig) * Evp_sv(ig) * dt__SV
469    QaT_SV(ig) = q1_new * (1. - beta(ig)) + beta(ig) * qsat_new
470
471  ENDDO
472
473END SUBROUTINE sisvat_ts2
Note: See TracBrowser for help on using the repository browser.