source: LMDZ6/trunk/libf/phymar/PHY_SISVAT_RUN.f90 @ 3472

Last change on this file since 3472 was 2089, checked in by Laurent Fairhead, 10 years ago

Inclusion de la physique de MAR


Integration of MAR physics

File size: 54.4 KB
Line 
1      subroutine PHY_SISVAT_RUN(                                        &
2!------------------------------------------------------------------------------+
3!                                                                              |
4!     SubRoutine PHY_SISVAT_RUN interfaces 3d-grid        Wed 26-Jun-2013 MAR  |
5!               with SISVAT                2d-grid                             |
6!                    SISVAT as  Soil                                           |
7!                               Ice                                            |
8!                               Snow                                           |
9!                               Vegetation                                     |
10!                               Atmosphere                                     |
11!                               Transfer   Scheme                              |
12!                                                                              |
13!     version 3.p.4.1 created by H. Gallee,               Tue 26-Feb-2013      |
14!           Last Modification by H. Gallee,               Wed 26-Jun-2013      |
15!                                                                              |
16!------------------------------------------------------------------------------+
17! #TC&                        qxTC,uqTC  ,qsTC  ,                       &      ! Aerosols: Atm.Conc., Surf.Flux
18     &                     )
19
20!------------------------------------------------------------------------------+
21!                                                                              |
22!     INPUT:    it_RUN: Run Iterations Counter                                 |
23!     ^^^^^                                                                    |
24!                                                                              |
25!     INPUT    (via module Mod_SISVAT_ctr)                                     |
26!     ^^^^^     VegMod: SISVAT    is set up when .T.                           |
27!               SnoMod: Snow Pack is set up when .T.                           |
28!               inpSBC: Update Bound.Condit.when .T.                           |
29!                                                                              |
30!     INPUT    (via common block)                                              |
31!     ^^^^^     xxxxTV: SISVAT/MAR interfacing variables                       |
32!                                                                              |
33!     Preprocessing  Option: SISVAT PHYSICS                                    |
34!     ^^^^^^^^^^^^^^^^^^^^^  ^^^^^^^^^^^^^^                                    |
35!   #                        #SN: Snow         Model                           |
36!   #                        #GP  LAI and GLF  Variations not specified        |
37!   #                        #OP  SST       is interactive                     |
38!                                                                              |
39!   #                        #DS: diffuse radiation differing from direct      |
40!                                (variable RADsod must still be included)      |
41!                                                                              |
42!     Preprocessing  Option: SISVAT PHYSICS: Col de Porte                      |
43!     ^^^^^^^^^^^^^^^^^^^^^  ^^^^^^^^^^^^^^^^^^^^^^^^^^^^                      |
44!   #                        #CP: SBL,                       Col de Porte      |
45!   #                        #cp  Solar Radiation,           Col de Porte      |
46!   #                        #AG: Snow Ageing,               Col de Porte      |
47!                                                                              |
48!                                                                              |
49!     Preprocessing  Option: STANDARD Possibility                              |
50!     ^^^^^^^^^^^^^^^^^^^^^  ^^^^^^^^^^^^^^^^^^^^                              |
51!     #HY: Explicit Cloud MICROPHYSICS              may be turned ON           |
52!     #BS: Explicit Cloud MICROPHYSICS: Blow. *(Snow)         Model)           |
53!     #SI: SISVAT: Sea-Ice Fraction calculated from prescribed SST             |
54!     #IP: SISVAT: Sea-Ice Fraction prescribed from SMMR and SSM/I             |
55!     #SN: SNOW Model                               may be turned ON           |
56!     #AO: COUPLING with  NEMO Ocean-Sea-Ice Model using OASIS                 |
57!     #AE: TURBULENCE: Aerosols Erosion / Turbulent Diffusion Coeff.           |
58!     #BD: TraCer   Aeolian Erosion  Submodel           is turned ON           |
59!     #OR: SBL: Orography Roughness included from z0__SV (MARdom)              |
60!     #SZ: SBL: Mom.: Roughn.Length= F(u*) Andreas &al.(1995)  Snow            |
61!     #ZM: SBL: M/H   Roughn.Length: Box Moving Average (in Time)              |
62!     #IB: OUTPUT: Ice-Sheet Surface Mass Balance  (on NetCDF File )           |
63!     #vL: PORTABILITY: Vectorization enhanced                                 |
64!     #vS: PORTABILITY: Vectorization enhanced: Snow    Model *                |
65!     #vD: PORTABILITY: Vectorization enhanced: Blowing  Dust .                |
66!     #vZ: PORTABILITY: Vectorization enhanced: Av.Roughness Length            |
67!     #AA  TURBULENCE: SBL  Time Mean (BOX Moving Average)                     |
68!     #AW  TURBULENCE: Wind Time Mean (BOX Moving Average)                     |
69!     #AH  TURBULENCE: Ta-T Time Mean (BOX Moving Average)                     |
70!                                                                              |
71!                                                                              |
72!     Preprocessing  Option: OBSOLETE                                          |
73!     ^^^^^^^^^^^^^^^^^^^^^  ^^^^^^^^                                          |
74!     #AM  TURBULENCE: u*   Time Mean (BOX Moving Average)                     |
75!     #AT  TURBULENCE: u*T* Time Mean (BOX Moving Average)                     |
76!     #AS  TURBULENCE: u*s* Time Mean (BOX Moving Average)                     |
77!                                                                              |
78!                                                                              |
79!     Preprocessing  Option:                                                   |
80!     ^^^^^^^^^^^^^^^^^^^^^                                                    |
81!     #TC: TraCer   Advection-Diffusion Equation        is turned ON           |
82!     #CM: Update Green Leaf Fraction may be performed (in routine INIglf)     |
83!     #VM: TURBULENCE: Ua   minimum is prescribed                              |
84!     #GP: Soil /Vegetation Model: LAI, GLF Variations NOT prescrib.           |
85!     #OP: SISVAT: Interactive Sea Surface Temperature     turned ON           |
86!     #op: SISVAT: SST Nudging -->   prescribed values     turned ON           |
87!     #RE: SISVAT: SST Condition     modified   by         epsilon             |
88!     #PO: POLYNYA Model                            may be turned ON           |
89!                                                                              |
90!                                                                              |
91!     Preprocessing  Option: SISVAT IO (not always a standard preprocess.)     |
92!     ^^^^^^^^^^^^^^^^^^^^^  ^^^^^^^^^                                         |
93!     FILE                 |      CONTENT                                      |
94!     ~~~~~~~~~~~~~~~~~~~~~+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~     |
95!   # SISVAT_iii_jjj_n     | #e0: OUTPUT on ASCII  File (SISVAT Variables)     |
96!   #                      |(#e0  MUST BE PREPROCESSED BEFORE #e1 & #e2 !)     |
97!   # SISVAT_iii_jjj_n     | #e1: OUTPUT/Verification: Energy Conservation     |
98!   # SISVAT_iii_jjj_n     | #e2: OUTPUT/Verification: Energy Consrv.2e pt.    |
99!                          |                           (no premature stop)     |
100!   # SISVAT_iii_jjj_n     | #e3: OUTPUT/Verification: Energy Consrv. (HL)     |
101!   # SISVAT_iii_jjj_n     | #e4: OUTPUT/Verification: Energy Consrv. (HLS)    |
102!   # SISVAT_iii_jjj_n     | #el: OUTPUT Format      : Addition.Frame Line     |
103!                          |                                                   |
104!   # SISVAT_iii_jjj_n     | #m0: OUTPUT/Verification: H2O    Conservation     |
105!   # SISVAT_iii_jjj_n     | #m1: OUTPUT/Verification: * Mass Conservation     |
106!   # SISVAT_iii_jjj_n     | #m2: OUTPUT/Verification: SeaIce Conservation     |
107!   # stdout               | #mw: OUTPUT/Verification: H2O    Conservation     |
108!                          |      unit  6, SubRoutine  SISVAT_qSo **ONLY**     |
109!                          |                                                   |
110!   # SISVAT_iii_jjj_n     | #mb: OUTPUT/Verification: Blowing Snow            |
111!                          |               SubRoutine  PHY_SISVAT **ONLY**     |
112!                          |                                                   |
113!   # SISVAT_zSn.vz        | #vz: OUTPUT/Verification: Snow Layers Agrega.     |
114!                          |      unit 41, SubRoutine  SISVAT_zSn **ONLY**     |
115!   # SISVAT_qSo.vw        | #vw: OUTPUT/Verif+Detail: H2O    Conservation     |
116!                          |      unit 42, SubRoutine  SISVAT_qSo **ONLY**     |
117!   # SISVAT_qSn.vm        | #vm: OUTPUT/Verification: Energy/Water Budget     |
118!                          |      unit 43, SubRoutine  SISVAT_qSn **ONLY**     |
119!   # SISVAT_qSn.vu        | #vu: OUTPUT/Verification: Slush  Parameteriz.     |
120!                          |      unit 44, SubRoutine  SISVAT_qSn **ONLY**     |
121!   # SISVAT_wEq.ve        | #ve: OUTPUT/Verification: Snow/Ice Water Eqv.     |
122!                          |      unit 45, SubRoutine  SISVAT_wEq **ONLY**     |
123!   # SnOptP____.va        | #va: OUTPUT/Verification: Albedo Parameteriz.     |
124!                          |      unit 46, SubRoutine  SnOptP     **ONLY**     |
125!   # SISVAT_GSn.vp        | #vp: OUTPUT/Verification: Snow   Properties       |
126!                          |      unit 47, SubRoutines SISVAT_zSn, _GSn        |
127!   # PHY_SISVAT.v0        | #v0: OUTPUT/Verification: DUMP                    |
128!                          |      unit 50, SubRoutine  PHY_SISVAT **ONLY**     |
129!                          |                                                   |
130!   # stdout               | #b0: OUTPUT of Snow Erosion Variables             |
131!                          |      unit  6, SubRoutine  SISVAT_BSn **ONLY**     |
132!   # stdout               | #b1: OUTPUT of Snow Erosion Turbulence            |
133!                          |      unit  6, SubRoutine  SISVATeSBL **ONLY**     |
134!   # stdout               | #b2: OUTPUT of Snow Erosion Turbulence (2)        |
135!                          |      unit  6, SubRoutine  SISVATeSBL **ONLY**     |
136!                          |                                                   |
137!   # stdout               | #s0: OUTPUT of Snow Buffer Layer                  |
138!                          |      unit  6, SubRoutine  SISVAT     **ONLY**     |
139!   # stdout               | #s1: OUTPUT of Snow Layers Agregation             |
140!                          |      unit  6, SubRoutine  SISVAT_zSn, _zAg        |
141!   # stdout               | #s2: OUTPUT of SnowFall, Snow Buffer              |
142!                          |      unit  6, SubRoutine  SISVAT_BSn, _qSn        |
143!   # stdout               | #sb: OUTPUT/Verification: SISVAT_SBL              |
144!                          |      unit  6, SubRoutine  SISVAT_SBL **ONLY**     |
145!   # stdout               | #sd: OUTPUT/Verification: Snow Thinest Layer      |
146!                          |      unit  6, SubRoutine  SISVAT_zSn **ONLY**     |
147!   # stdout               | #sf: OUTPUT of SnowFall, Z0 and Drag Coeff.       |
148!                          |      unit  6, SubRoutines PHY_SISVAT, SISVAT      |
149!   # stdout               | #sg: OUTPUT/Verification: Gravitational Front     |
150!                          |      unit  6, SubRoutine  SISVAT_qSo **ONLY**     |
151!   # stdout               | #sh: OUTPUT/Verification: SBL Turbulent Fluxes    |
152!                          |      unit  6, SubRoutine  PHY_SISVAT **ONLY**     |
153!   # stdout               | #si: OUTPUT of Sea Ice Fraction, Temperature      |
154!                          |      unit  6, SubRoutine  PHY_SISVAT **ONLY**     |
155!   # stdout               | #sm: OUTPUT/Verification: Mosaic Fractions Sum    |
156!                          |      unit  6, SubRoutine  PHY_SISVAT **ONLY**     |
157!   # stdout               | #so: OUTPUT/Verification: Soil Vertic.Discret.    |
158!                          |      unit  6, SubRoutine  SISVAT_ini **ONLY**     |
159!   # stdout               | #ss: OUTPUT of Blowing Snow Variables             |
160!                          |      unit  6, SubRoutine  SISVATeSBL **ONLY**     |
161!   # stdout               | #sz: OUTPUT of Roughness Length & Drag Coeff.     |
162!                          |      unit  6, SubRoutine  SISVAT     **ONLY**     |
163!                          |                                                   |
164!------------------------------------------------------------------------------+
165
166
167! Global Variables
168! ================
169
170      use Mod_Real
171      use Mod_PHY____dat
172      use Mod_PHY____grd
173      use Mod_SISVAT_dim
174      use Mod_SISVAT_grd
175      use Mod_SISVAT_ctr
176      use Mod_SISVAT_gpt
177
178      use Mod_PHY_S0_dat
179      use Mod_PHY_S0_kkl
180      use Mod_PHY_RT_kkl
181      use Mod_PHY_DY_kkl
182      use Mod_PHY_AT_kkl
183      use Mod_PHY_CM_kkl
184
185
186! INTERFACE Variables
187! ===================
188
189      use Mod_SISVAT_dat
190      use Mod_SISVAT_dzS
191      use Mod_SISVAT_kkl
192      use Mod_SISVAT_loc
193! #AW use Mod_SISVAT_xAW
194! #AH use Mod_SISVAT_xAH
195! #ZM use Mod_SISVAT_xZM
196! #WL use Mod_SISVAT_xWL
197
198
199
200      IMPLICIT NONE
201
202
203
204      INTEGER    ntrac
205      PARAMETER (ntrac=1)
206
207
208! Tracers
209! ^^^^^^^
210! #TC real       uqTC  ( mx, my,    ntrac)  ! q*     (Tracer)
211! #TC real       qsTC  ( mx, my,    ntrac)  ! qSurf  (Tracer)
212! #TC real       qxTC  ( mx, my, mz,ntrac)  ! q      (Tracer)
213
214! Auxiliary Variables
215! ^^^^^^^^^^^^^^^^^^^
216! #CM integer                                          ::  newglfSIS    !
217! #IP integer                                          ::  newsicSI     !
218
219
220
221
222! Internal  Variables
223! ===================
224
225      logical            glfFIX
226
227      character(len=1)   cha  ,chb
228      integer            jjtime
229      integer            iwr
230! #SI integer            l
231      integer            isl  ,isn
232      integer            i    ,j    ,k     ,ikp
233      integer            n    ,nm
234
235      real(kind=real8)   rhAir                ! Air     Densitity
236      real(kind=real8)   Ua_min               ! Minimum Air Velocity
237      real(kind=real8)   d_snow,SnowOK        ! Snow Precip.: Total
238! #BS real(kind=real8)   dbsnow               ! Snow Precip.: from Drift
239! #WL real(kind=real8)   uqstar               ! u*q*
240
241! #SI real(kind=real8)   SrfSIC,SIc0OK        ! Oceanic Fraction: previous
242! #SI real(kind=real8)   FraOcn,SIceOK        ! Oceanic Fraction
243
244! #Za real(kind=real8)   S_Eros               ! Snow Erosion (status) = (1,0)
245! #Za real(kind=real8)   SnEros               ! Snow Erosion (status) = (1,0)
246
247! OUTPUT/Verification:   Blowing Snow
248! #mb integer            noUNIT
249! #mb real(kind=real8)   BlowST,SnowSB
250
251
252! DATA
253! ====
254
255      data     glfFIX   /  .false./
256      data     cha      /'-'/
257      data     chb      /':'/
258
259
260
261! SISVAT Time
262! ===========
263
264          jjtime = HourTU*3600+minuTU*60+Sec_TU
265
266
267
268
269! SISVAT Variables:      (x,y,z) Domain
270! =====================================
271
272      DO ikp=1,kcolp
273      DO nm =1,mwp
274
275
276! SISVAT Variables: Interpolation of V at 10 m a.g.l.
277! ---------------------------------------------------
278
279        IF (k_SL.EQ.mzp)                                            THEN
280          VV10SV(ikp,nm) = r_SL10 * WindSV(ikp,nm,k_SL)                 !
281        ELSE
282          VV10SV(ikp,nm) =          WindSV(ikp,nm,k_SL)                &!
283     &                   + r_SL10 *(WindSV(ikp,nm,k_SL-1)              &!
284     &                             -WindSV(ikp,nm,k_SL))                !
285        END IF
286
287
288! Sastrugi Height decreased by Precipitation if V < 6 m/s (Kotlyakov, 1961)
289! -------------------------------------------------------------------------
290
291! #ZS   dz0_SV(ikp,nm)  =      max(zer0,(SnowCM(ikp)-Sno0CM(ikp))      &!
292! #ZS&    /max(0.05,0.104*sqrt(max(zer0, WindSV(ikp,nm,mzp)-6.00))))    !
293! #ZS   dz0_SV(ikp,nm)  =                dz0_SV(ikp,nm) *0.01*max(2-n,0)! dz0(Sastrugi dh)
294
295
296! Influence of the Angle(Wind,Sastrugi) (Andreas, 1995, CCREL report 95-16)
297! -------------------------------------------------------------------------
298
299! #Za   S_Eros     = max(zer0,sign(un_1,-uss_CM(ikp) -eps9))
300! #Za   SnEros     = max(zer0,sign(un_1, uss_CM(ikp) +eps9))
301! #Za   VVs_SV(ikp,nm)    = VVs_SV(ikp,nm) * RRsxSV(ikp,nm)
302! #Za   VVs_SV(ikp,nm)    =                                            &
303! #Za&           SnEros* VVs_SV(ikp,nm)                                &
304! #Za&         + S_Eros*(VVs_SV(ikp,nm) * Adz0dt     + VV10SV(ikp,nm))
305! #Za   RRsxSV(ikp,nm)    =                                            &
306! #Za&           SnEros* RRsxSV(ikp,nm)                                &
307! #Za&         + S_Eros*(RRsxSV(ikp,nm) * Adz0dt     + 1.0 )
308! #Za   DDsxSV(ikp,nm)    =                                            &
309! #Za&           SnEros* DDsxSV(ikp,nm)                                &
310! #Za&         + S_Eros* DDsxSV(ikp,nm) * Adz0dt                       &
311! #Za&         +        (Va__SV(ikp,nm) *(Ua__SV(ikp,nm) -Ua0_SV(ikp,nm))  &
312! #Za&                  -Ua__SV(ikp,nm) *(Va__SV(ikp,nm) -Va0_SV(ikp,nm))) &
313! #Za&   /(degrad*max(.3,WindSV(ikp,nm,mzp) * WindSV(ikp,nm,mzp)))
314! #Za    IF (DDsxSV(ikp,nm)   .GT.360.) DDsxSV(ikp,nm)  = DDsxSV(ikp,nm) - 360.
315! #Za    IF (DDsxSV(ikp,nm)   .LT.  0.) DDsxSV(ikp,nm)  = DDsxSV(ikp,nm) + 360.
316      END DO
317      END DO
318
319
320! Water Vapor Flux Limitor
321! ------------------------
322
323! #WL DO ikp=1,kcolp
324! #WL DO nm =1,mwp
325
326! #WL  DO n=1,nLimi-1
327! #WL    WL_mem(ikp,nm,n) = WL_mem(ikp,nm,n+1)
328! #WL  END DO
329
330! #WL    uqstar  =  max(abs(uqs_SV_gpt(ikp),eps6)*sign(1.,uqs_SV_gpt(ikp))
331! #WL    WL_mem(ikp,nm,n)                                              &
332! #WL& = Kzh_AT(ikp,nm,k1m(mzp))
333! #WL&    *(qv__SV(ikp,nm,mzpp-k2m(mzp))-qv__SV(ikp,nm,mzpp-k1m(mzp))) &
334! #WL&    /(uqstar    *(hsigma(k2m(mzp))            -hsigma(k1m(mzp))))
335
336! #WL     WLmmem(ikp,nm) = 0.
337! #WL  DO n=1,nLimi
338! #WL     WLmmem(ikp,nm) = WLmmem(ikp,nm)+WL_mem(ikp,nm,n)
339! #WL  ENDDO
340! #WL     WLmmem(ikp,nm) = WLmmem(ikp,nm)/nLimit
341
342! #WL END DO
343! #WL END DO
344
345
346
347! SWD / SWA INPUT
348! ---------------
349
350      DO ikp = 1,kcolp
351
352      IF (.NOT.InpSWD)                                              THEN
353        SWDsRT(ikp) = SWAsRT(ikp)/(1.d0-Alb_SV_gpt(ikp))                ! => downward Solar
354      END IF
355
356
357! Grid Averages
358! -------------
359
360        Alb_SV_gpt(ikp)  = 0.
361        EmisSV_gpt(ikp)  = 0.
362        LWUsRT    (ikp)  = 0.
363        LMO_SV_gpt(ikp)  = 0.
364        us__SV_gpt(ikp)  = 0.
365        uts_SV_gpt(ikp)  = 0.
366        uqs_SV_gpt(ikp)  = 0.
367        uss_CM(ikp)      = 0.
368        qs__CM(ikp,mzpp) = 0.
369        hFraSV_gpt(ikp)  = 0.
370        Tas_SV_gpt(ikp)  = 0.
371        qv__DY(ikp,mzpp) = 0.
372! #TC     uqTC(i,j,1)    = 0.
373! #TC     qsTC(i,j,1)    = 0.
374      END DO
375
376
377
378! SISVAT Variables: from (x,y,z) Domain to (vector,mosaic) Domain
379! ===============================================================
380
381! SISVAT Variables Update
382! ^^^^^^^^^^^^^^^^^^^^^^^
383      DO      ikp         = 1,kcolp
384      DO      nm          = 1,mwp
385              i           = ii__AP(ikp)
386              j           = jj__AP(ikp)
387
388! OUTPUT of SnowFall, Roughness Length and Drag Coefficients
389! #sf         IF (ikp.EQ.1 .AND. Sec_TU.EQ.0)
390! #sf&        write(6,6659)
391! #sf 6659    format(20x,'   dsn_SV   us__SV   Z0SaSi   Z0Sa_N'        &
392! #sf&                  ,'   Z0SaSV   Z0m_Sn   Z0m_SV')
393
394! Atmospheric Forcing                                       (INPUT)
395! ^^^^^^^^^^^^^^^^^^^                                        ^^^^^
396            DO isn=1,mzp
397              isl =  mzp + 1  - isn                                     !
398              Kz__SV(ikp,nm,isn) = Kzh_AT(ikp,isl)                      !
399              pktaSV(ikp,nm,isn) = pkt0SV(ikp,nm,isn)                   !
400            ENDDO
401              za__SV(ikp,nm) =    zza_SV(ikp,nm,1)                      ! [m]
402                                                                        !
403              VV__SV(ikp,nm) =    WindSV(ikp,nm,mzp)                    !
404
405              Ua_min         = eps6                                     !
406! #VM         Ua_min         = 0.2 * sqrt(za__SV(ikp,nm) )              !
407              VV__SV(ikp,nm) = max(Ua_min,VV__SV(ikp,nm) )              !
408                                                                        !
409! #Za         VVs_SV(ikp,nm) = VVs_SV(ikp,nm) / RRsxSV(ikp,nm)          !
410! #Za         DDs_SV(ikp,nm) =         max(zer0,DDsxSV(ikp,nm)-180.)   &
411! #Za&        +180.*min(un_1,zer0-min(zer0,DDsxSV(ikp,nm)     -180.))  &
412! #Za&                           +min(zer0,DDsxSV(ikp,nm)     -180.)
413              rhT_SV(ikp,nm) = pkPaSV(ikp,nm)   *1.e3                  &! [kg/m3]
414     &                       /(TaT_SV(ikp,nm)   *R_DAir)                !
415              QaT_SV(ikp,nm) = qv__SV(ikp,nm,1)                         !
416! #WL         dQa_SV(ikp,nm) =  max(0.,1.-WLmmem(ikp,nm))              &! Water  Vapor
417! #WL&                        * dt__SV/hsigma(mzp)                      ! Flux Limitor
418              qsnoSV(ikp,nm) =      0.                                 &!
419     &                       + qs__CM(ikp,mzp)                         &!
420! #TC&                       +   qxTC(i,j,mzp,1)                       &!
421     &                       +   0.
422
423! Energy Fluxes                                             (INPUT)
424! ^^^^^^^^^^^^^                                              ^^^^^
425              coszSV(ikp,nm) = max(cszEPS,csz_S0(ikp))                  ! cos(zenith.Dist.)
426              sol_SV(ikp,nm) =            SWDsRT(ikp)                   !    downward Solar
427              IRd_SV(ikp,nm) =            LWDsRT(ikp)                   !    downward IR
428
429! Water  Fluxes                                             (INPUT)
430! ^^^^^^^^^^^^^                                              ^^^^^
431              drr_SV(ikp,nm) =(RainCM(ikp)-Rai0CM(ikp)) *1.e3 /dt__SV   ! [m/s] -> [mm/s] .OR. [kg/m2/s]
432              d_snow         = SnowCM(ikp)-SnobCM(ikp)                  ! Only SnowFall
433              dsn_SV(ikp,nm) = d_snow                   *1.e3 /dt__SV   ! Erosion NOT incl.
434              SnowOK         =                                         &! Correction
435     &         max(zer0,sign(un_1,       qs__CM(ikp,mzp)-eps6))        &!
436     &        *max(zer0,sign(un_1,Tf_Sno-TaT_SV(ikp,nm) -eps6))         !
437              dsn_SV(ikp,nm) = dsn_SV(ikp,nm)+drr_SV(ikp,nm)*SnowOK     !
438              drr_SV(ikp,nm) = drr_SV(ikp,nm)           *(1.-SnowOK)    !
439! #BS         dbsnow         =-ussxSV(ikp,nm)                          &! Erosion
440! #BS&                        *dt__SV     *rhT_SV(ikp,nm)               !
441! #bS         dsnbSV(ikp,nm) =min(max(zer0,dbsnow)                     &!
442! #bS&                    /       max(eps6,d_snow),un_1)                !
443! #BS         dsnbSV(ikp,nm) =1.0-min(qs__CM(ikp,k_zb)                 &! k_zb level ~ 25 magl
444! #BS&                           /max(qs__CM(ikp,mzp ),epsn),un_1)      !(BS negligib.at k_zb)
445!             dsnbSV is used and modified in SISVAT_BSn,
446!                  then used for Buffer Layer Update
447! #BS         dbs_SV(ikp,nm) =        ussbSV(ikp,nm)                    !          [kg/m2]
448
449! Soil/Canopy                                               (INPUT)
450! ^^^^^^^^^^^                                                ^^^^^
451              slorSV(ikp,nm) = atan(slopSV(ikp,nm))                     ! Fall Line Slope       [radian]
452
453! Energy Fluxes                                             (INPUT/OUTPUT)
454! ^^^^^^^^^^^^^                                              ^^^^^^^^^^^^
455              cld_SV(ikp,nm) =      ClouRT(ikp)                         ! Cloudiness
456
457! Water  Fluxes                                             (INPUT/OUTPUT)
458! ^^^^^^^^^^^^^                                              ^^^^^^^^^^^^
459! #BS         uss_SV(ikp,nm) =      ussxSV(ikp,nm)                      ! u*s*
460
461! #vL ENDDO
462! #vL ENDDO
463
464! Snow Roughness                                            (INPUT/OUTPUT)
465! ^^^^^^^^^^^^^^                                             ^^^^^^^^^^^^
466! #vL DO      ikp = 1,kcolp
467! #vL DO      nm  = 1,mwp
468
469! #ZM         Z0mmSV(ikp,nm) =     0.                                   !
470! #ZM         Z0emSV(ikp,nm) =     0.                                   !
471! #ZM         Z0hmSV(ikp,nm) =     0.                                   !
472! #ZM       DO k  =   1,ntavz                                           !
473! #ZM         Z0mmSV(ikp,nm) =     Z0mmSV(ikp,nm)                      &!
474! #ZM&                       +     z0_mem(ikp,nm,k)                     !
475! #ZM         Z0emSV(ikp,nm) =     Z0emSV(ikp,nm)                      &!
476! #ZM&                       +     b0_mem(ikp,nm,k)                     !
477! #ZM         Z0hmSV(ikp,nm) =     Z0hmSV(ikp,nm)                       !
478! #ZM&                       +     r0_mem(ikp,nm,k)                     !
479! #ZM       ENDDO                                                       !
480! #ZM         Z0mmSV(ikp,nm) = min(Z0mmSV(ikp,nm)  /ntavz              &!  z0(Mom., Box Av.)
481! #ZM&                            ,hsigma(mzp)     /    3.)             !
482! #ZM         Z0emSV(ikp,nm) =     Z0emSV(ikp,nm)  /ntavz               !  b0(Eros, Box Av.)
483! #ZM         Z0hmSV(ikp,nm) =     Z0hmSV(ikp,nm)  /ntavz               !  r0(Heat, Box Av.)
484
485! V,  dT(a-s)    Time Moving Averages
486! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
487! #AA       DO k =  1,ntaver-1
488! #AW         V__mem(ikp,nm,k    ) = V__mem(ikp,nm,k+1)
489! #AH         T__mem(ikp,nm,k    ) = T__mem(ikp,nm,k+1)
490! #AA       ENDDO
491! #AW         V__mem(ikp,nm,ntave) = VV__SV(ikp,nm)
492! #AH         T__mem(ikp,nm,ntave) = TaT_SV(ikp,nm)-Ts__SV(ikp,nm)
493
494! #AW         VVmmem(ikp,nm)       = 0.0
495! #AH         dTmmem(ikp,nm)       = 0.0
496! #AA       DO k=1,ntaver
497! #AW         VVmmem(ikp,nm)       = VVmmem(ikp,nm)+V__mem(ikp,nm,k)
498! #AH         dTmmem(ikp,nm)       = dTmmem(ikp,nm)+T__mem(ikp,nm,k)
499! #AA       ENDDO
500! #AW         VVmmem(ikp,nm)        = VVmmem(ikp,nm)/ntaver
501! #AH         dTmmem(ikp,nm)        = dTmmem(ikp,nm)/ntaver
502
503! #vL ENDDO
504! #vL ENDDO
505
506! Snow Pack                                                 (INPUT/OUTPUT)
507! ^^^^^^^^^                                                  ^^^^^^^^^^^^
508! #vS DO      ikp = 1,kcolp
509! #vS DO      nm  = 1,mwp
510              dsn_SV(ikp,nm) =     dsn_SV(ikp,nm)                      &!
511     &                        +max(BufsSV(ikp,nm)-SMndSV,0.)           &!
512     &                        /    dt__SV                               !
513              BufsSV(ikp,nm) = min(BufsSV(ikp,nm),SMndSV   )            !
514! #vS ENDDO
515! #vS ENDDO
516
517      DO      isn = 1,nsnow
518! #vS DO      ikp = 1,kcolp
519! #vS DO      nm  = 1,mwp
520              G1snSV(ikp,nm,isn) = max(-G1_dSV,min(G1_dSV,             &!
521     &                                  G1snSV(ikp,nm,isn)))            ! [-]        [-]
522              G2snSV(ikp,nm,isn) = max(-G1_dSV,min(G1_dSV,             &!
523     &                                  G2snSV(ikp,nm,isn)))            ! [-] [0.0001 m]
524! #vS END DO
525! #vS END DO
526      END DO
527
528! #vL DO      ikp = 1,kcolp
529! #vL DO      nm  = 1,mwp
530! #vL         i   = ii__AP(ikp)   
531! #vL         j   = jj__AP(ikp)   
532
533! #SI         HFraSV(ikp,nm)     = 0.                         ! Frazil Thickness
534
535! RunOFF Intensity                                          (INPUT/OUTPUT)
536! ^^^^^^^^^^^^^^^^                                           ^^^^^^^^^^^^
537              RnofSV(ikp,nm)     = 0.                         ! RunOFF Intensity
538
539! Blown Snow OUTPUT for a specific Grid Point
540! ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
541! #mb    IF  (lwriSV(ikp,nm).ne.0.AND.it_RUN.gt.0)                     THEN
542! #mb         noUNIT =   no__SV(lwriSV(ikp,nm))
543! #mb         write(noUNIT,5012)
544! #mb 5012    format(/,1x)
545! #mb         write(noUNIT,5013)
546! #mb 5013    format(' -----+--------+--------+--------+--------+',    &
547! #mb&             '--------+--------+--------+--------+--------+',    &
548! #mb&             '--------+')
549! #mb         write(noUNIT,5014)
550! #mb 5014    format('    n |     z  |     qs |      V |        |',    &
551! #mb&             '     T  | TKE^0.5|        |        |        |',    &
552! #mb&             '        |',                                        &
553! #mb&             /,'      |    [m] | [g/kg] |  [m/s] |        |',    &
554! #mb&             '    [K] |  [m/s] |        |        |        |',    &
555! #mb&             '        |')
556! #mb           BlowST=0.
557! #mb           k=0
558! #mb 5011    CONTINUE
559! #mb           k=k+1
560! #mb         IF (               k             .gt.mzp )      GO TO 5010
561! #mb         IF (Z___DY    (ikp,k)-sh__AP(ikp).lt.100.)            THEN
562! #mb           BlowST=BlowST+WindSV(ikp,nm,k)*qs__CM    (ikp,k)       &
563! #mb&                 *1.e3*(pkPaSV(ikp,nm)  -pt__DY)*d1Sigm(k) *Grav_I
564! #mb           write(noUNIT,5015) mzpp-k                              &
565! #mb&           ,     Z___DY(ikp,k)   -sh__AP(ikp)                    &
566! #mb&           ,     qs__CM(ikp,k)   *1.e3                           &
567! #mb&           ,     WindSV(ikp,nm,k),Ta__DY(ikp,k)                  &
568! #mb&           ,sqrt(TKE_AT(ikp,k))         
569! #mb 5015      format(i5,' |',f7.2,' |',f7.3,' |',f7.2,' |',          &
570! #mb&                  8x,'|',f7.2,' |',f7.3,' |', 4(8x,'|'))
571! #mb         END IF
572! #mb         GO TO 5011
573! #mb 5010    CONTINUE
574
575! #mb           SnowSB =   BufsSV(ikp,nm)     
576! #mb         IF        (  isnoSV(ikp,nm)   .GT.     0      )       THEN
577! #mb         DO isn=max(0,isnoSV(ikp,nm)  ) ,isnoSV(ikp,nm) 
578! #mb           SnowSB =   SnowSB                                      &
579! #mb&                +    dzsnSV(ikp,nm,isn)*ro__SV(ikp,nm,isn)
580! #mb         END DO
581! #mb         END IF
582
583! #mb         write(noUNIT,5016) BlowST,SnowSB
584! #mb 5016    format(' * TRANSPORT = ',e12.3,' kg/m/s'  ,  8x,'|',     &
585! #mb&               ' * BUDGET    = ',f12.6,' mm w.e.|',2(8x,'|'))
586! #mb         write(noUNIT,5013)
587! #mb    END IF
588
589! OUTPUT, for Stand-Alone VERIFICATION
590! ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
591! #v0    IF  (ikp .le. kcolp .and. nm .le. mxpp)                    THEN
592! #v0         write(50,5001) it_RUN,ikp,nm                             &
593! #v0&                   za__SV(ikp,nm),VV__SV(ikp,nm),TaT_SV(ikp,nm), &
594! #v0&                   rhT_SV(ikp,nm),QaT_SV(ikp,nm),qsnoSV(ikp,nm), &
595! #v0&                   coszSV(ikp,nm),sol_SV(ikp,nm),IRd_SV(ikp,nm), &
596! #v0&                   drr_SV(ikp,nm),dsn_SV(ikp,nm),dbs_SV(ikp,nm), &
597! #v0&                   LSmask(ikp,nm),isotSV(ikp,nm),alb0SV(ikp,nm), &
598! #v0&                   IRs_SV(ikp,nm),                               &
599! #v0&                   ivgtSV(ikp,nm),LAI0SV(ikp,nm),glf0SV(ikp,nm), &
600! #v0&                   TvegSV(ikp,nm),LMO_SV(ikp,nm),us__SV(ikp,nm), &
601! #v0&                   uqs_SV(ikp,nm),uts_SV(ikp,nm),uss_SV(ikp,nm), &
602! #v0&                   snCaSV(ikp,nm),rrCaSV(ikp,nm),psivSV(ikp,nm)
603! #v0 5001    format(/,'c #INFO   it_RUN             = ',i15  ,        &
604! #v0&               /,'c #INFO   ikp                = ',i15  ,        &
605! #v0&               /,'c #INFO   nm                 = ',i15  ,        &
606! #v0&               /,'          za__SV(ikp,nm)     = ',e15.6,        &
607! #v0&               /,'          VV__SV(ikp,nm)     = ',e15.6,        &
608! #v0&               /,'          TaT_SV(ikp,nm)     = ',e15.6,        &
609! #v0&               /,'          rhT_SV(ikp,nm)     = ',e15.6,        &
610! #v0&               /,'          QaT_SV(ikp,nm)     = ',e15.6,        &
611! #v0&               /,'          qsnoSV(ikp,nm)     = ',e15.6,        &
612! #v0&               /,'          coszSV(ikp,nm)     = ',e15.6,        &
613! #v0&               /,'          sol_SV(ikp,nm)     = ',e15.6,        &
614! #v0&               /,'          IRd_SV(ikp,nm)     = ',e15.6,        &
615! #v0&               /,'          drr_SV(ikp,nm)     = ',e15.6,        &
616! #v0&               /,'          dsn_SV(ikp,nm)     = ',e15.6,        &
617! #v0&               /,'          dbs_SV(ikp,nm)     = ',e15.6,        &
618! #v0&               /,'          LSmask(ikp,nm)     = ',i15  ,        &
619! #v0&               /,'          isotSV(ikp,nm)     = ',i15  ,        &
620! #v0&               /,'          alb0SV(ikp,nm)     = ',e15.6,        &
621! #v0&               /,'          IRs_SV(ikp,nm)     = ',e15.6,        &
622! #v0&               /,'          ivgtSV(ikp,nm)     = ',i15  ,        &
623! #v0&               /,'          LAI0SV(ikp,nm)     = ',e15.6,        &
624! #v0&               /,'          glf0SV(ikp,nm)     = ',e15.6,        &
625! #v0&               /,'          TvegSV(ikp,nm)     = ',e15.6,        &
626! #v0&               /,'          LMO_SV(ikp,nm)     = ',e15.6,        &
627! #v0&               /,'          us__SV(ikp,nm)     = ',e15.6,        &
628! #v0&               /,'          uqs_SV(ikp,nm)     = ',e15.6,        &
629! #v0&               /,'          uts_SV(ikp,nm)     = ',e15.6,        &
630! #v0&               /,'          uss_SV(ikp,nm)     = ',e15.6,        &
631! #v0&               /,'          snCaSV(ikp,nm)     = ',e15.6,        &
632! #v0&               /,'          rrCaSV(ikp,nm)     = ',e15.6,        &
633! #v0&               /,'          psivSV(ikp,nm)     = ',e15.6)
634! #v0         DO isl = -nsoil,0
635! #v0           write(50,5002)    isl,TsisSV(ikp,nm,isl),              &
636! #v0&                            isl,eta_SV(ikp,nm,isl)
637! #v0 5002      format('          TsisSV(ikp,nm,',i2,')  = ',e15.6,    &
638! #v0&                 '          eta_SV(ikp,nm,',i2,')  = ',e15.6)
639! #v0         END DO
640! #v0         DO isl =       1,nsnow
641! #v0           write(50,5003)    isl,TsisSV(ikp,nm,isl),              &
642! #v0&                            isl,dzsnSV(ikp,nm,isl)
643! #v0 5003      format('          TsisSV(ikp,nm,',i2,')  = ',e15.6,    &
644! #v0&                 '          dzsnSV(ikp,nm,',i2,')  = ',e15.6)
645! #v0         END DO
646! #v0    END IF
647      END DO
648      END DO
649
650! SISVAT Execution
651! ^^^^^^^^^^^^^^^^
652          write(daHost,'(i2,a3,i4,i3,2(a1,i2))')                        &
653     &          Day_TU,LabMon(Mon_TU),YearTU,                           &
654     &          HourTU,chb,minuTU,chb,Sec_TU
655
656! OUTPUT of SnowFall, Roughness Length and Drag Coefficients
657! #sf     write(6,666) Day_TU,Mon_TU,YearTU,HourTU,minuTU,Sec_TU
658! #sf 666 format(2(i2,'-'),2i4,2(':',i2),3x,$)
659
660
661
662!           ******
663      call  SISVAT(jjtime,kcolv)
664!           ******
665
666
667
668! MAR    Variables Update
669! ^^^^^^^^^^^^^^^^^^^^^^^
670      DO      ikp = 1,kcolp
671      DO      k   =    1,mzp
672              dpktSV_gpt(ikp,k) = 0.
673      ENDDO
674      ENDDO
675
676      DO      ikp = 1,kcolp
677      DO      nm  = 1,mwp
678
679! Atmospheric Forcing                                       (OUTPUT)
680! ^^^^^^^^^^^^^^^^^^^                                        ^^^^^^
681            DO isn=1,mzp
682              k   =  mzp + 1  - isn                                     !
683              dpktSV_gpt(ikp,k)     =                                  &!
684      &       dpktSV_gpt(ikp,k)     +                                  &!
685      &       FracSV(ikp,nm)        *                                  &!
686      &      (pktaSV(ikp,nm,isn)-pkt0SV(ikp,nm,isn)) / dt__SV           !
687            ENDDO
688
689! Water  Fluxes                                             (INPUT/OUTPUT)
690! ^^^^^^^^^^^^^                                              ^^^^^^^^^^^^
691! #BS     ussxSV(ikp,nm)         = 0.                                   !
692! #BS     ussxSV(ikp,nm)         =                                     &! MAX        Erosion
693! #BS&   (ussbSV(ikp,nm)  -        dbs_SV(ikp,nm))                     &!-unconsumed Erosion
694! #BS&                    /(dt__SV*rhT_SV(ikp,nm))                      ! ==> actual u*s*
695! #vL ENDDO
696! #vL ENDDO
697
698! #vL DO      ikp = 1,kcolp
699! #vL DO      nm  = 1,mwp 
700! #BS     ussbSV(ikp,nm)  = dt__SV*uss_SV(ikp,nm)                      &! NEW  MAX   Erosion
701! #BS&                            *rhT_SV(ikp,nm)                       ! rho u*s* dt[kg/m2]
702
703! Soil/Canopy                                               (INPUT/OUTPUT)
704! ^^^^^^^^^^^                                                ^^^^^^^^^^^^
705          Z0m_SV(ikp,nm)  =    min(Z0m_SV(ikp,nm)                      &! Moment.Roughn.L.
706      &                           ,hsigma(mzp)/3.)                      !
707! #vL END DO
708! #vL END DO
709
710! Snow Roughness                                            (INPUT/OUTPUT)
711! ^^^^^^^^^^^^^^                                             ^^^^^^^^^^^^
712! #vZ DO      ikp = 1,kcolp
713! #vZ DO      nm  = 1,mwp 
714! #ZM       DO k=1,ntavz-1                                              !
715! #ZM         z0_mem(ikp,nm,k)       = z0_mem(ikp,nm,k+1)               !
716! #ZM         b0_mem(ikp,nm,k)       = b0_mem(ikp,nm,k+1)               !
717! #ZM         r0_mem(ikp,nm,k)       = r0_mem(ikp,nm,k+1)               !
718! #ZM       ENDDO                                                       !
719! #vZ ENDDO
720! #vZ ENDDO
721
722! #vL DO      ikp = 1,kcolp
723! #vL DO      nm  = 1,mwp 
724! #ZM         z0_mem(ikp,nm,ntavz)   = Z0mnSV(ikp,nm)                   ! z0(Momentum)
725! #ZM         b0_mem(ikp,nm,ntavz)   = Z0enSV(ikp,nm)                   ! z0(Mom., Erosion)
726! #ZM         r0_mem(ikp,nm,ntavz)   = Z0hnSV(ikp,nm)                   ! z0(Heat)
727! #vL END DO
728! #vL END DO
729
730! Dust   Fluxes                                             (INPUT/OUTPUT)
731! ^^^^^^^^^^^^^
732! #vD DO      ikp = 1,kcolp
733! #vD DO      nm  = 1,mwp 
734! #BD         uds_SV(ikp,nm)         =(1-min(1,isnoSV(ikp,nm)))        &! Snow Free  Surface
735! #BD&                                        *uss_SV(ikp,nm)          &! DUST       Erosion
736! #BD&                                  *max(1,(2-my)*3)                ! Tuning Factor (2D)
737! #vD END DO
738! #vD END DO
739
740! Snow Pack                                                 (INPUT/OUTPUT)
741! ^^^^^^^^^                                                  ^^^^^^^^^^^^
742! #vL DO      ikp = 1,kcolp
743! #vL DO      nm  = 1,mwp 
744! #IB         wes0SV(ikp,nm)         = dwesSV(ikp,nm)                   ! Depo. / Subli.
745! #IB         wem0SV(ikp,nm)         = dwemSV(ikp,nm)                   ! Melting
746! #IB         wer0SV(ikp,nm)         = dwerSV(ikp,nm)                   ! Refreezing
747
748! Radiative Properties                                            (OUTPUT)
749! ^^^^^^^^^^^^^^^^^^^^                                             ^^^^^^
750          albcSV(ikp,nm)             = alb_SV(ikp,nm)                  &! Mosaic Albedo
751! #AO&                             *   LSmask(ikp,nm)                  &!
752! #AO&                             +   Alb_AO(ikp,nm)                  &! Mosaic AlbedoNEMO
753! #AO&                             *(1-LSmask(ikp,nm))                 &!
754     &                             + 0.
755          ROF_SV(ikp,nm)             = RnofSV(ikp,nm)   *dt__SV        &!
756     &                             +   ROF_SV(ikp,nm)                   ! Integrated Run OFF
757
758      END DO
759      END DO
760
761
762
763! Surface Temperature: Prescription of relevant Medium (Snow, precribed SST)
764! ==========================================================================
765
766! Control Martin
767!PRINT*, 'nsoil=',nsoil
768!PRINT*, 'FixSST=',FixSST
769!PRINT*, 'LSMask=',LSMask
770! Control Martin
771
772      DO ikp=1,kcolp
773            DO isl =   -nsoil,0                                          !
774
775
776! Open Ocean
777! ----------
778
779              eta_SV(ikp,1,isl) =                                      &!
780     &        eta_SV(ikp,1,isl) *   LSMask(ikp,1)                      &!
781     &                          +(1-LSMask(ikp,1)    )                  ! Sea: Humidity:=1
782
783              TsisSV(ikp,1,isl) =                                      &!
784     &       (TsisSV(ikp,1,isl) *   LSMask(ikp,1)                      &! Soil Temperature
785     &       +sst_SB(ikp)       *(1-LSMask(ikp,1)    )) * FixSST       &! Prescribed   SST
786! #OP&                                                  * FixSST       &!
787! #OP&      +(TsisSV(ikp,1,isl)                                        &!
788! #op&      +(sst_SB(ikp)                                              &!~Prescribed   SST
789! #op&       -TsisSV(ikp,1,isl))*(1-LSMask(ikp,1)    )  * SSTnud       &! (Nudging)
790! #OP&                                               )  * VarSST       &! Interactive  SST
791     &      + 0.
792
793
794! Sea Ice
795! -------
796
797! #AO         eta_SV(ikp,2,isl) =                                      &!
798! #AO&        eta_SV(ikp,2,isl) *   LSMask(ikp,2)                       ! Sea: Humidity:=0
799
800! #AO         TsisSV(ikp,2,isl) =                                      &!
801! #AO&       (TsisSV(ikp,2,isl) =   LSMask(ikp,2)                      &! Soil Temperature
802! #AO&       +  271.2           *(1-LSMask(ikp,2)    ))                 ! Prescribed    ST
803            END DO
804
805! #AO       DO isl =         0,nsnow                                    !
806! #AO         TsisSV(ikp,2,isl) =                                      &!
807! #AO&        s_T_AO(ikp,2)     *(1-LSMask(ikp,2)    )                  ! Prescribed   SIT
808! #AO       END DO
809      ENDDO
810
811
812      DO      ikp = 1,kcolp
813      DO      nm  = 1,mwp 
814          Ts__SV(ikp,nm)         =                                     &! Surf.Temperature
815     &    TsisSV(ikp,nm,0)                                             &!
816     &                    *(1-min(1,isnoSV    (ikp,nm)  ))             &!
817     &   +TsisSV(ikp,nm,      max(1,isnoSV    (ikp,nm)  ))             &!
818     &                    *   min(1,isnoSV    (ikp,nm)  )              &!
819     &       +0.
820      ENDDO
821      ENDDO
822
823
824
825! Mosaic Cleaning
826! ===============
827
828      DO      ikp = 1,kcolp
829        IF (LSMask(ikp,1)    .EQ.   0)                              THEN
830            isnoSV    (ikp,1)     = 0
831            ispiSV    (ikp,1)     = 0
832            iiceSV    (ikp,1)     = 0
833          DO ISL=1,nsnow
834            TsisSV    (ikp,1,isl) = 0.
835            dzsnSV    (ikp,1,isl) = 0.
836            ro__SV    (ikp,1,isl) = 0.
837            eta_SV    (ikp,1,isl) = 0.
838            G1snSV    (ikp,1,isl) = 0.
839            G2snSV    (ikp,1,isl) = 0.
840            agsnSV    (ikp,1,isl) = 0.
841            istoSV    (ikp,1,isl) = 0.
842          ENDDO
843        ENDIF
844        IF (FracSV(ikp,n2) .LT. epsn .AND. n2 .GT. 1)               THEN
845            Ts__SV(ikp,n2)         = Ts__SV(ikp,1) 
846          DO isl=-nsoil,0
847            TsisSV    (ikp,n2,isl) = TsisSV(ikp,1,isl)
848          ENDDO    ! #n2 
849            isnoSV    (ikp,n2    ) = isnoSV(ikp,1    ) * LSMask(ikp,n2)
850            ispiSV    (ikp,n2    ) = ispiSV(ikp,1    ) * LSMask(ikp,n2)
851            iiceSV    (ikp,n2    ) = iiceSV(ikp,1    ) * LSMask(ikp,n2)
852          DO isl=1,nsnow      !   
853            TsisSV    (ikp,n2,isl) = TsisSV(ikp,1,isl) * LSMask(ikp,n2)
854            dzsnSV    (ikp,n2,isl) = dzsnSV(ikp,1,isl) * LSMask(ikp,n2)
855            ro__SV    (ikp,n2,isl) = ro__SV(ikp,1,isl) * LSMask(ikp,n2)
856            eta_SV    (ikp,n2,isl) = eta_SV(ikp,1,isl) * LSMask(ikp,n2)
857            G1snSV    (ikp,n2,isl) = G1snSV(ikp,1,isl) * LSMask(ikp,n2)
858            G2snSV    (ikp,n2,isl) = G2snSV(ikp,1,isl) * LSMask(ikp,n2)
859            agsnSV    (ikp,n2,isl) = agsnSV(ikp,1,isl) * LSMask(ikp,n2)
860            istoSV    (ikp,n2,isl) = istoSV(ikp,1,isl) * LSMask(ikp,n2)
861          ENDDO    ! #n2 
862        ENDIF      ! #n2 
863      ENDDO
864
865
866
867! Grid Averages / Diagnostics
868! ===========================
869
870! Grid Averages                                                   (OUTPUT)
871! ^^^^^^^^^^^^^                                                    ^^^^^^
872      DO ikp = 1,kcolp
873      DO nm  = 1,mwp
874! #IB         wes_SV    (ikp,nm)= - wes0SV(ikp,nm) + wes_SV(ikp,nm)     ! Depo. / Subli.
875! #IB         wem_SV    (ikp,nm)=   wem0SV(ikp,nm) + wem_SV(ikp,nm)     ! Melting
876! #IB         wer_SV    (ikp,nm)=   wer0SV(ikp,nm) + wer_SV(ikp,nm)     ! Refreezing
877! #IB         wee_SV    (ikp,nm)= - uqs_SV(ikp,nm) + wee_SV(ikp,nm)    &! Evapotranspiration
878! #IB&                            * dt__SV         * roa_SV(ikp,nm,1)   !
879
880              Alb_SV_gpt(ikp)   =   Alb_SV_gpt(ikp)                    &! Grid   Albedo
881     &      + FracSV    (ikp,nm)  * albcSV    (ikp,nm)                  ! Mosaic Albedo
882              EmisSV_gpt(ikp)   =   EmisSV_gpt(ikp)                    &! Grid   Emissivity
883     &      + FracSV    (ikp,nm)  * emi_SV    (ikp,nm)                  ! Mosaic Emissivity
884              LWUsRT    (ikp)   =   LWUsRT    (ikp)                    &!
885     &      + FracSV    (ikp,nm)  * IRu_SV    (ikp,nm)                  ! Mosaic Upw.IR
886              LMO_SV_gpt(ikp)   =   LMO_SV_gpt(ikp)                    &!
887     &      + FracSV    (ikp,nm)  * LMO_SV    (ikp,nm)                  ! Mosaic Mon.Ob.
888              us__SV_gpt(ikp)   =   us__SV_gpt(ikp)                    &! Grid   u*
889     &      + FracSV    (ikp,nm)  * us__SV    (ikp,nm)                  ! Mosaic u*
890              uts_SV_gpt(ikp)   =   uts_SV_gpt(ikp)                    &! Grid   u*T*
891     &      + FracSV    (ikp,nm)  * uts_SV    (ikp,nm)                  ! Mosaic u*T*
892              uqs_SV_gpt(ikp)   =   uqs_SV_gpt(ikp)                    &! Grid   u*q*
893     &      + FracSV    (ikp,nm)  * uqs_SV    (ikp,nm)                  ! Mosaic u*q*
894! #BS         uss_CM    (ikp)   =   uss_CM    (ikp)                    &! Grid   u*s*
895! #BS&      + FracSV    (ikp,nm)  * ussxSV    (ikp,nm)                  ! Mosaic u*s*
896!     NO !    ussxSV    (ikp,nm)  = uss_SV    (ikp,nm)                  !        u*s*
897!        !    Upper Update = wrong Source of Atmospher.Snow             !
898! #BS         qs__CM    (ikp,mzpp) =                                   &! Salt.Part.Concent.
899! #BS&        qs__CM    (ikp,mzpp)                                     &!
900! #BS&      + FracSV    (ikp,nm)  * qSalSV    (ikp,nm)                 &!
901! #BS&                       *min(1,isnoSV    (ikp,nm)  )               !
902! #PO         hFraSV_gpt(ikp)     = hFraSV_gpt(ikp)                    &! Frazil  Thickness
903! #PO&      + FracSV    (ikp,nm)  * HFraSV    (ikp,nm)                  !
904              Tas_SV_gpt(ikp)     = Tas_SV_gpt(ikp)                    &! Surface Air
905     &      + FracSV    (ikp,nm)  * Ts__SV    (ikp,nm)                  !         Temperatur
906              qv__DY    (ikp,mzpp)= qv__DY    (ikp,mzpp)               &!
907     &      + FracSV    (ikp,nm)  * SHumSV    (ikp,nm)                  ! Surf.Specif.Humid.
908                                                                        ! to adapt over soil
909! #TC         uqTC      (i,j,1)   = uqTC      (i,j,1)                  &! Grid   u*d*
910! #TC&      + FracSV    (ikp,nm)  * uds_SV    (ikp,nm)                  ! Mosaic u*d*
911! #TC         qsTC      (i,j,1)   = qsTC      (i,j,1)                  &! Salt.Part.Concent.
912! #TC&      + FracSV    (ikp,nm)  * qSalSV    (ikp,nm)                 &!
913! #TC&                    *(1-min(1,isnoSV    (ikp,nm)  ))              !
914      ENDDO
915      ENDDO
916
917
918
919! De la grille "SISVAT" vers la grille "AtmPHY"
920! =============================================
921
922      DO ikp = 1,kcolp
923        Sno0CM    (ikp) = SnowCM    (ikp)                               !
924         rhAir          = roa_SV    (ikp,1,1)                           ! Air    Densitity
925        Ta__DY    (ikp,mzpp) =                                         &!
926     &                    Tas_SV_gpt(ikp)                               !
927        pkt_DY    (ikp,mzpp) =                                         &!
928     &                    Tas_SV_gpt(ikp) /pkPaSV(ikp,1) ** RCp         !
929        HsenSV_gpt(ikp) =-uts_SV_gpt(ikp) *rhAir    *CpdAir             ! Sensible Heat Flux
930        HLatSV_gpt(ikp) =-uqs_SV_gpt(ikp) *rhAir    *LhvH2O             ! Latent   Heat Flux
931        WE2aSV_gpt(ikp) = WE2aSV_gpt(ikp)                              &! Total    Evaporat.
932     &                   -uqs_SV_gpt(ikp) *rhAir    *dt__SV             ! [mm w.e.]
933      END DO
934
935
936
937! Blown Snow/Dust Accumulation
938! ============================
939
940! #AE DO ikp=1,kcolp
941! #AE DO nm =1,mwp
942
943! #BS    SnowCM(ikp) = SnowCM(ikp)                                     &!
944! #BS&               + dt__SV *1.e-3 *roa_SV(ikp,1) *uss_CM(ikp)
945
946! #AE END DO
947! #AE END DO
948
949
950
951! Formation of Lakes                                                                 \VER
952! ==================
953
954
955
956! Sea-Ice Ice Floe Size
957! =====================
958
959! Prescription from SST
960! ---------------------
961
962! #SI IF (VarSST.le.eps6)                                           THEN
963! #SI   DO ikp=1,kcolp
964! #SI        FraOcn          =(TsisSV    (ikp,1,0)-TOF_SV)/TSIdSV       ! Prescribed
965!                                                                       ! from SST
966! #IP        FraOcn          =  1.-sif_SB(ikp)                          ! Prescribed
967!                                                                       ! from SSM/I
968! #SI        FraOcn          = min(  un_1,FraOcn)                       ! UpperLimit
969! #SI        FraOcn          = max(OcnMin,FraOcn)                       ! LowerLimit
970! #SI        FracSV(ikp,1)   =     LSMask(ikp,1)     * FracSV(ikp,1)   &! New Ocean
971! #SI&                        + (1-LSMask(ikp,1)    )* FraOcn           !
972! #SI        SrfSIC          =                         FracSV(ikp,2)    ! Old Sea Ice
973! #SI        SIc0OK          = max(zer0, sign(un_1,SrfSIC-eps6))        !
974! #SI        FracSV(ikp,2)     *   LSMask(ikp,2)     * FracSV(ikp,2)   &! New Sea Ice
975! #SI&                        + (1-LSMask(ikp,2)    )*(1.-FraOcn)       !
976! #SI        SIceOK          = max(zer0, sign(un_1,FracSV(ikp,2)       &!
977! #SI&                                                   -eps6))        !
978
979! Sea-Ice Vertical Discretization
980! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
981! #SI        isnoSV(ikp,2)     =                                       &
982! #SI&       isnoSV(ikp,2)                     *   LSMask(ikp,2)       &
983! #SI&     +(max(1                                                     &
984! #SI&      ,isnoSV(ikp,2)  )  *SIc0OK                                 &
985! #SI&      +     3        *(1.-SIc0OK)*SIceOK)*(1-LSMask(ikp,2))
986
987! #SI        iiceSV(ikp,2)     =                                       &
988! #SI&       iiceSV(ikp,2)                     *   LSMask(ikp,2)       &
989! #SI&     +(max(1                                                     &
990! #SI&      ,iiceSV(ikp,2)  )  *SIc0OK                                 &
991! #SI&      +     3        *(1.-SIc0OK)*SIceOK)*(1-LSMask(ikp,2))
992! #SI        ispiSV(ikp,2)     =                   iiceSV(ikp,2) 
993
994! #SI     DO l=1,nsnow
995! #SI        dzsnSV(ikp,2,l) =                                         &
996! #SI&       dzsnSV(ikp,2,l)                   *   LSMask(ikp,2)       &
997! #SI&     +(max                                                       &
998! #SI&      (SIc_OK(min(2,l))*    SIcMIN                               &
999! #SI&      ,dzsnSV(ikp,2,l))*SIc0OK                                   &
1000! #SI&      +dzSIce(min(4,l))*(1.-SIc0OK)*SIceOK)*(1-LSMask(ikp,2))
1001
1002! #SI        TsisSV(ikp,2,l) =                                         &
1003! #SI&       TsisSV(ikp,2,l)                 *   LSMask(ikp,2)         &
1004! #SI&     +(TsisSV(ikp,2,l) *    SIc0OK                               &
1005! #SI&      +TsisSV(ikp,1,0) *(1.-SIc0OK)   )*(1-LSMask(ikp,2))
1006
1007! #SI        ro__SV(ikp,2,l) =                                         &
1008! #SI&       ro__SV(ikp,2,l)                 *   LSMask(ikp,2)         &
1009! #SI&     +(max                                                       &
1010! #SI&      (SIc_OK(min(2,l))*rhoIce                                   &
1011! #SI&      ,ro__SV(ikp,2,l))*SIc0OK                                   &
1012! #SI&      +rhoIce      *(1.-SIc0OK)*SIceOK)*(1-LSMask(ikp,2))
1013
1014! #SI        G1snSV(ikp,2,l) =                                         &
1015! #SI&       G1snSV(ikp,2,l)                 *   LSMask(ikp,2)         &
1016! #SI&     +(G1snSV(ikp,2,l) *SIc0OK                                   &
1017! #SI&      +G1_dSV      *(1.-SIc0OK)*SIceOK)*(1-LSMask(ikp,2))
1018
1019! #SI        G2snSV(ikp,2,l) =                                         &
1020! #SI&       G2snSV(ikp,2,l)                 *   LSMask(ikp,2)         &
1021! #SI&     +(G2snSV(ikp,2,l) *SIc0OK                                   &
1022! #SI&      +30.         *(1.-SIc0OK)*SIceOK)*(1-LSMask(ikp,2))
1023
1024! #SI        istoSV(ikp,2,l) =                                         &
1025! #SI&       istoSV(ikp,2,l)                 *   LSMask(ikp,2)         &
1026!!#SI&new  +(istoSV(ikp,2,l) *SIc0OK                                   &
1027!!#SI&new   +istdSV(2)   *(1.-SIc0OK)*SIceOK)*(1-LSMask(ikp,2))        &
1028! #SI&     + istdSV(2)   *                    (1-LSMask(ikp,2))
1029! #SI     END DO
1030! #SI     DO l=-nsoil,0
1031! #SI        TsisSV(ikp,2,l) =                                         &
1032! #SI&       TsisSV(ikp,2,l)                 *   LSMask(ikp,2)         &
1033! #SI&     +(TsisSV(ikp,2,l) *    SIc0OK                               &
1034! #SI&      +TsisSV(ikp,1,l) *(1.-SIc0OK)   )*(1-LSMask(ikp,2))
1035
1036! #SI        eta_SV(ikp,2,l) =                                         &
1037! #SI&       eta_SV(ikp,2,l)                 *   LSMask(ikp,2)         &
1038! #SI&     + eta_SV(ikp,2,l) *    SIc0OK     *(1-LSMask(ikp,2))
1039!                                 No Pore in Ice => No Water
1040! #SI     END DO
1041
1042! OUTPUT of Sea Ice Fraction, Temperature, Discretization
1043! #si        write(6,6001) Day_TU,LabMon(Mon_TU),YearTU                &
1044! #si&                    ,HourTU,minuTU,Sec_TU ,TsisSV(ikp,1,0)       &
1045! #si&                    ,FraOcn,FracSV(ikp,1) ,TsisSV(ikp,2,0)       &
1046! #si&                    ,       iiceSV(ikp,2) ,isnoSV(ikp,2) 
1047
1048! #SI   END DO
1049! #SI END IF
1050
1051! Otherwise SST and FrLead have been computed in the Sea-Ice Polynya Model
1052! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1053
1054
1055! Rainfall, Snowfall Time Integral at previous Time Step
1056! ------------------------------------------------------
1057
1058        DO ikp = 1,kcolp
1059           Rai0CM(ikp) = RainCM(ikp)                                    ! Rainfall Time Integral
1060           SnobCM(ikp) = SnowCM(ikp)                                    ! Snowfall Time Integral
1061                                                                        ! Erosion  skipped
1062
1063
1064!  Wind Horizontal Components         at previous Time Step
1065!  --------------------------------------------------------
1066
1067        DO nm=1,mwp
1068! #Za      Ua0_SV(ikp,nm) = Ua__SV(ikp,nm)                              !
1069! #Za      Va0_SV(ikp,nm) = Va__SV(ikp,nm)                              !
1070        END DO
1071        END DO
1072
1073
1074! Additional OUTPUT for VERIFICATION
1075! ----------------------------------
1076
1077! OUTPUT/Verification: SBL Turbulent Fluxes
1078! #sh i   = i_x0 + 10.*111.111e3/dxHOST
1079! #sh j   = j_y0
1080! #sh nm  = 1
1081! #sh ikp = ikl_AP(i,j)
1082! #sh write(6,6060) it_EXP,Day_TU,LabMon(Mon_TU),YearTU                    &
1083! #sh&                    ,HourTU,minuTU        ,lat__r(ikp)/Dg2Rad        &
1084! #sh&         ,TaT_SV    (ikp,nm)              ,roa_SV(ikp,nm,1)          &
1085! #sh&         ,HsenSV_gpt(ikp),HLatSV_gpt(ikp),-86400.*uqs_SV_gpt(ikp)    &
1086! #sh&        ,1.e3*RainCM(ikp),WE2aSV_gpt(ikp),        ROF_SV    (ikp,nm) 
1087! #sh 6060 format(i6,i3,'-',a3,'-',i4,':',i2,':',i2,f6.2,' N',             &
1088! #sh&            f9.3,' K'     ,f6.3,' kg/m3',2(f6.1,' W/m2'),            &
1089! #sh&            f6.3,' mm/day',3(f9.3,' mm'))
1090
1091
1092
1093! SISVAT OUTPUTs
1094! ==============
1095
1096!     IF (WRIsbc)                                                   THEN
1097
1098!            **************
1099!       CALL PHY_SISVAT_IOs
1100!            **************
1101
1102!     END IF
1103
1104
1105
1106      return
1107      end subroutine PHY_SISVAT_RUN
Note: See TracBrowser for help on using the repository browser.