source: LMDZ6/branches/Ocean_skin/libf/phylmdiso/surf_land_mod.F90 @ 5453

Last change on this file since 5453 was 4368, checked in by lguez, 2 years ago

Sync latest trunk changes to Ocean_skin

File size: 15.0 KB
Line 
1!
2MODULE surf_land_mod
3
4  IMPLICIT NONE
5
6CONTAINS
7!
8!****************************************************************************************
9
10  SUBROUTINE surf_land(itime, dtime, date0, jour, knon, knindex, &
11       rlon, rlat, yrmu0, &
12       debut, lafin, zlev, ccanopy, swnet, lwnet, albedo, &
13       tsurf, p1lay, cdragh, cdragm, precip_rain, precip_snow, temp_air, spechum, &
14       AcoefH, AcoefQ, BcoefH, BcoefQ, &
15       AcoefU, AcoefV, BcoefU, BcoefV, &
16       pref, u1, v1, gustiness, rugoro, pctsrf, &
17       lwdown_m, q2m, t2m, &
18       snow, qsol, agesno, tsoil, &
19       z0m, z0h, SFRWL, alb_dir_new, alb_dif_new, evap, fluxsens, fluxlat, &   
20       qsurf, tsurf_new, dflux_s, dflux_l, &
21       flux_u1, flux_v1 , &
22       veget,lai,height &
23#ifdef ISO
24       ,xtprecip_rain, xtprecip_snow,xtspechum, &
25       xtsnow, xtsol,xtevap,h1, &
26       runoff_diag,xtrunoff_diag,Rland_ice &
27#endif               
28               )
29
30    USE dimphy
31    USE surface_data, ONLY    : ok_veget
32    USE carbon_cycle_mod
33
34    ! See comments in each module surf_land_orchidee_xxx for compatiblity with ORCHIDEE
35#ifdef ORCHIDEE_NOOPENMP
36    ! Compilation with cpp key ORCHIDEE NOOPENMP
37    USE surf_land_orchidee_noopenmp_mod
38#else
39#if ORCHIDEE_NOZ0H
40    ! Compilation with cpp key ORCHIDEE NOZ0H
41    USE surf_land_orchidee_noz0h_mod
42#else
43#if ORCHIDEE_NOFREIN
44    ! Compilation with cpp key ORCHIDEE_NOFREIN
45    USE surf_land_orchidee_nofrein_mod
46#else
47#if ORCHIDEE_NOUNSTRUCT
48    ! Compilation with cpp key ORCHIDEE_NOUNSTRUCT
49    USE surf_land_orchidee_nounstruct_mod
50#else
51#if ORCHIDEE_NOLIC
52    ! Compilation with cpp key ORCHIDEE_NOLIC
53    USE surf_land_orchidee_nolic_mod
54#else
55    ! Default version#else
56    USE surf_land_orchidee_mod
57#endif
58#endif
59#endif
60#endif
61#endif
62
63    USE surf_land_bucket_mod
64    USE calcul_fluxs_mod
65    USE indice_sol_mod
66#ifdef ISO
67    use infotrac_phy, ONLY: ntiso,niso
68    use isotopes_mod, ONLY: nudge_qsol, iso_eau
69#ifdef ISOVERIF
70    use isotopes_verif_mod
71#endif
72#endif
73   
74    USE print_control_mod, ONLY: lunout
75
76    INCLUDE "dimsoil.h"
77    INCLUDE "YOMCST.h"
78    INCLUDE "clesphys.h"
79    INCLUDE "dimpft.h"
80
81! Input variables 
82!****************************************************************************************
83    INTEGER, INTENT(IN)                     :: itime, jour, knon
84    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex
85    REAL, INTENT(IN)                        :: date0
86    REAL, DIMENSION(klon), INTENT(IN)       :: rlon, rlat
87    REAL, DIMENSION(klon), INTENT(IN)       :: yrmu0  ! cosine of solar zenith angle
88    LOGICAL, INTENT(IN)                     :: debut, lafin
89    REAL, INTENT(IN)                        :: dtime
90    REAL, DIMENSION(klon), INTENT(IN)       :: zlev, ccanopy
91    REAL, DIMENSION(klon), INTENT(IN)       :: swnet, lwnet
92    REAL, DIMENSION(klon), INTENT(IN)       :: albedo  ! albedo for whole short-wave interval
93    REAL, DIMENSION(klon), INTENT(IN)       :: tsurf
94    REAL, DIMENSION(klon), INTENT(IN)       :: p1lay
95    REAL, DIMENSION(klon), INTENT(IN)       :: cdragh, cdragm
96    REAL, DIMENSION(klon), INTENT(IN)       :: precip_rain, precip_snow
97    REAL, DIMENSION(klon), INTENT(IN)       :: temp_air, spechum
98    REAL, DIMENSION(klon), INTENT(IN)       :: AcoefH, AcoefQ, BcoefH, BcoefQ
99    REAL, DIMENSION(klon), INTENT(IN)       :: AcoefU, AcoefV, BcoefU, BcoefV
100    REAL, DIMENSION(klon), INTENT(IN)       :: pref   ! pressure reference
101    REAL, DIMENSION(klon), INTENT(IN)       :: u1, v1, gustiness
102    REAL, DIMENSION(klon), INTENT(IN)       :: rugoro
103    REAL, DIMENSION(klon,nbsrf), INTENT(IN) :: pctsrf
104    REAL, DIMENSION(klon), INTENT(IN)       :: lwdown_m  ! downwelling longwave radiation at mean surface
105                                                         ! corresponds to previous sollwdown
106    REAL, DIMENSION(klon), INTENT(IN)       :: q2m, t2m
107#ifdef ISO
108    REAL, DIMENSION(ntiso,klon), INTENT(IN)       :: xtprecip_rain, xtprecip_snow
109    REAL, DIMENSION(ntiso,klon), INTENT(IN)       :: xtspechum
110#endif
111
112! In/Output variables
113!****************************************************************************************
114    REAL, DIMENSION(klon), INTENT(INOUT)          :: snow, qsol
115    REAL, DIMENSION(klon), INTENT(INOUT)          :: agesno
116    REAL, DIMENSION(klon, nsoilmx), INTENT(INOUT) :: tsoil
117#ifdef ISO
118    REAL, DIMENSION(niso,klon), INTENT(INOUT)    :: xtsnow, xtsol
119#endif
120
121! Output variables
122!****************************************************************************************
123    REAL, DIMENSION(klon), INTENT(OUT)       :: z0m, z0h
124!albedo SB >>>
125!    REAL, DIMENSION(klon), INTENT(OUT)       :: alb1_new ! albdeo for shortwave interval 1(visible)
126!    REAL, DIMENSION(klon), INTENT(OUT)       :: alb2_new ! albedo for shortwave interval 2(near infrared)
127    REAL, DIMENSION(6), INTENT(IN) :: SFRWL
128    REAL, DIMENSION(klon,nsw), INTENT(OUT)       :: alb_dir_new,alb_dif_new
129!albedo SB <<<
130    REAL, DIMENSION(klon), INTENT(OUT)       :: evap
131    REAL, DIMENSION(klon), INTENT(OUT)       :: fluxsens, fluxlat
132    REAL, DIMENSION(klon), INTENT(OUT)       :: qsurf
133    REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
134    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l     
135    REAL, DIMENSION(klon), INTENT(OUT)       :: flux_u1, flux_v1  ! flux for U and V at first model level
136    REAL, DIMENSION(klon,nvm_lmdz), INTENT(OUT) :: veget,lai
137    REAL, DIMENSION(klon,nvm_lmdz), INTENT(OUT) :: height
138#ifdef ISO
139    REAL, DIMENSION(ntiso,klon), INTENT(OUT)      :: xtevap
140    REAL, DIMENSION(klon), INTENT(OUT)      :: h1
141    REAL, DIMENSION(niso,klon), INTENT(OUT)      :: xtrunoff_diag
142    REAL, DIMENSION(klon), INTENT(OUT)      :: runoff_diag
143    REAL, DIMENSION(niso,klon), INTENT(IN)        :: Rland_ice
144#endif
145
146! Local variables
147!****************************************************************************************
148    REAL, DIMENSION(klon) :: p1lay_tmp
149    REAL, DIMENSION(klon) :: pref_tmp
150    REAL, DIMENSION(klon) :: swdown     ! downwelling shortwave radiation at land surface
151    REAL, DIMENSION(klon) :: epot_air           ! potential air temperature
152    REAL, DIMENSION(klon) :: tsol_rad, emis_new ! output from interfsol not used
153    REAL, DIMENSION(klon) :: u0, v0     ! surface speed
154    INTEGER               :: i
155
156!albedo SB >>>
157    REAL, DIMENSION(klon)      :: alb1_new,alb2_new
158!albedo SB <<<
159
160#ifdef ISO       
161      real, parameter :: t_coup = 273.15
162      real, dimension(klon) :: fqfonte_diag
163      real, dimension(klon) ::  snow_evap_diag
164      real, dimension(klon) ::  fqcalving_diag
165      integer :: ixt
166#endif 
167
168#ifdef ISO
169#ifdef ISOVERIF
170!        write(*,*) 'surf_land_mod 162'
171        do i=1,knon
172          if (iso_eau.gt.0) then
173            call iso_verif_egalite_choix(precip_snow(i), &
174     &          xtprecip_snow(iso_eau,i),'surf_land_mod 129', &
175     &          errmax,errmaxrel)
176            call iso_verif_egalite_choix(qsol(i), &
177     &          xtsol(iso_eau,i),'surf_land_mod 139', &
178     &          errmax,errmaxrel)
179          endif 
180        enddo
181#endif
182#ifdef ISOVERIF
183!       write(*,*) 'surf_land 169: ok_veget=',ok_veget
184        do i=1,knon
185         do ixt=1,ntiso
186           call iso_verif_noNaN(xtprecip_snow(ixt,i),'surf_land 146')
187         enddo
188        enddo
189#endif
190#endif
191
192
193!****************************************************************************************
194! Choice between call to vegetation model (ok_veget=true) or simple calculation below
195!
196!****************************************************************************************
197   IF (ok_veget) THEN
198!****************************************************************************************
199!  Call model sechiba in model ORCHIDEE
200!
201!****************************************************************************************
202       p1lay_tmp(:)      = 0.0
203       pref_tmp(:)       = 0.0
204       p1lay_tmp(1:knon) = p1lay(1:knon)/100.
205       pref_tmp(1:knon)  = pref(1:knon)/100.
206!
207!* Calculate incoming flux for SW and LW interval: swdown
208!
209       swdown(:) = 0.0
210       DO i = 1, knon
211          swdown(i) = swnet(i)/(1-albedo(i))
212       END DO
213!
214!* Calculate potential air temperature
215!
216       epot_air(:) = 0.0
217       DO i = 1, knon
218          epot_air(i) = RCPD*temp_air(i)*(pref(i)/p1lay(i))**RKAPPA
219       END DO
220
221#ifdef ISO
222      CALL abort_gcm('surf_land_mod 220','isos pas prevus dans orchidee',1)
223#endif
224       ! temporary for keeping same results using lwdown_m instead of lwdown
225       CALL surf_land_orchidee(itime, dtime, date0, knon, &
226            knindex, rlon, rlat, yrmu0, pctsrf, &
227            debut, lafin, &
228            zlev,  u1, v1, gustiness, temp_air, spechum, epot_air, ccanopy, &
229            cdragh, AcoefH, AcoefQ, BcoefH, BcoefQ, &
230            precip_rain, precip_snow, lwdown_m, swnet, swdown, &
231            pref_tmp, q2m, t2m, &
232            evap, fluxsens, fluxlat, &             
233            tsol_rad, tsurf_new, alb1_new, alb2_new, &
234            emis_new, z0m, z0h, qsurf, &
235            veget, lai, height &
236!#ifdef ISO
237!            , xtprecip_rain, xtprecip_snow, xtspechum, xtevap &
238!#endif
239            )                 
240
241#ifdef ISO
242#ifdef ISOVERIF
243     write(*,*) 'surf_land 193: apres surf_land_orchidee'   
244     do i=1,knon
245        if (iso_eau.gt.0) then
246             call iso_verif_egalite_choix(xtevap(iso_eau,i),evap(i), &
247    &            'surf_land 197',errmax,errmaxrel)
248        endif !if (iso_eau.gt.0) then     
249      enddo !do i=1,knon 
250#endif
251#endif
252
253!* Add contribution of relief to surface roughness
254
255       DO i=1,knon
256          z0m(i) = MAX(1.5e-05,SQRT(z0m(i)**2 + rugoro(i)**2))
257       ENDDO
258
259    ELSE  ! not ok_veget
260!****************************************************************************************
261! No extern vegetation model choosen, call simple bucket calculations instead.
262!
263!****************************************************************************************
264#ifdef ISO
265#ifdef ISOVERIF
266!       write(*,*) 'surf_land 247'
267        call iso_verif_egalite_vect1D( &
268     &           xtsnow,snow,'surf_land_mod 207',niso,klon)
269#endif
270#endif
271
272#ifdef ISO
273        if (nudge_qsol.eq.1) then
274          call surf_land_nudge_qsol(knon,rlat,rlon,qsol,xtsol,knindex)
275        endif
276        !write(*,*) 'surf_land 258'
277#endif
278       CALL surf_land_bucket(itime, jour, knon, knindex, debut, dtime,&
279            tsurf, p1lay, cdragh, precip_rain, precip_snow, temp_air, &
280            spechum, AcoefH, AcoefQ, BcoefH, BcoefQ, pref, &
281            u1, v1, gustiness, rugoro, swnet, lwnet, &
282            snow, qsol, agesno, tsoil, &
283            qsurf, z0m, alb1_new, alb2_new, evap, &
284            fluxsens, fluxlat, tsurf_new, dflux_s, dflux_l &
285#ifdef ISO
286            ,xtprecip_rain, xtprecip_snow,xtspechum, &
287            xtsnow, xtsol,xtevap,h1, &
288     &      runoff_diag, xtrunoff_diag,Rland_ice &
289#endif           
290     &       )
291        z0h(1:knon)=z0m(1:knon) ! En attendant mieux
292
293    ENDIF ! ok_veget
294
295!****************************************************************************************
296! Calculation for all land models
297! - Flux calculation at first modele level for U and V
298!****************************************************************************************
299! Suppose zero surface speed
300    u0(:)=0.0
301    v0(:)=0.0
302    CALL calcul_flux_wind(knon, dtime, &
303         u0, v0, u1, v1, gustiness, cdragm, &
304         AcoefU, AcoefV, BcoefU, BcoefV, &
305         p1lay, temp_air, &
306         flux_u1, flux_v1)
307
308#ifdef ISO
309#ifdef ISOVERIF
310!     write(*,*) 'surf_land 237: sortie'   
311     do i=1,knon
312        if (iso_eau.gt.0) then
313             call iso_verif_egalite_choix(xtsnow(iso_eau,i),snow(i), &
314    &            'surf_land 241',errmax,errmaxrel)
315        endif !if (iso_eau.gt.0) then     
316      enddo !do i=1,knon 
317#endif
318#endif
319
320!albedo SB >>>
321     SELECT CASE(NSW)
322     CASE(2)
323       alb_dir_new(1:knon,1)=alb1_new(1:knon)
324       alb_dir_new(1:knon,2)=alb2_new(1:knon)
325     CASE(4)
326       alb_dir_new(1:knon,1)=alb1_new(1:knon)
327       alb_dir_new(1:knon,2)=alb2_new(1:knon)
328       alb_dir_new(1:knon,3)=alb2_new(1:knon)
329       alb_dir_new(1:knon,4)=alb2_new(1:knon)
330     CASE(6)
331       alb_dir_new(1:knon,1)=alb1_new(1:knon)
332       alb_dir_new(1:knon,2)=alb1_new(1:knon)
333       alb_dir_new(1:knon,3)=alb1_new(1:knon)
334       alb_dir_new(1:knon,4)=alb2_new(1:knon)
335       alb_dir_new(1:knon,5)=alb2_new(1:knon)
336       alb_dir_new(1:knon,6)=alb2_new(1:knon)
337     END SELECT
338
339     alb_dif_new=alb_dir_new
340!albedo SB <<<
341   
342  END SUBROUTINE surf_land
343
344
345#ifdef ISO
346    SUBROUTINE surf_land_nudge_qsol(knon,rlat,rlon,qsol,xtsol,knindex)
347
348    USE dimphy   
349    USE infotrac_phy, ONLY: niso
350    use isotopes_mod, ONLY: region_nudge_qsol   
351    INTEGER, INTENT(IN)                     ::  knon         
352    REAL, DIMENSION(klon), INTENT(IN)       :: rlon, rlat
353    REAL, DIMENSION(klon), INTENT(INOUT)    :: qsol
354    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex   
355    REAL, DIMENSION(niso,klon), INTENT(INOUT)    :: xtsol
356    real lat_min_nudge_qsol,lat_max_nudge_qsol
357    real lon_min_nudge_qsol,lon_max_nudge_qsol
358    integer i,ixt
359    real qsol_new
360
361  if (region_nudge_qsol.eq.1) then
362        ! Aamzonie du Sud
363        lat_min_nudge_qsol=-15.0
364        lat_max_nudge_qsol=-5.0
365        lon_min_nudge_qsol=-70.0
366        lon_max_nudge_qsol=-50.0
367  else if (region_nudge_qsol.eq.2) then
368        ! Aamzonie du Nord
369        lat_min_nudge_qsol=-5.0
370        lat_max_nudge_qsol=5.0
371        lon_min_nudge_qsol=-70.0
372        lon_max_nudge_qsol=-50.0
373  else
374        write(*,*) 'surf_land 298: cas pas prevu'
375        write(*,*) 'region_nudge_qsol=',region_nudge_qsol
376        stop
377  endif
378
379!  write(*,*) 'surf_land 314: knon=',knon
380!  write(*,*) 'rlat=',rlat
381!  write(*,*) 'rlon=',rlon
382!  write(*,*) 'region_nudge_qsol=',region_nudge_qsol
383
384    do i=1,knon
385     if ((rlat(knindex(i)).ge.lat_min_nudge_qsol).and. &
386  &    (rlat(knindex(i)).le.lat_max_nudge_qsol).and. &
387  &    (rlon(knindex(i)).ge.lon_min_nudge_qsol).and. &
388  &    (rlon(knindex(i)).le.lon_max_nudge_qsol)) then
389!        write(*,*) 'surf_land 324: bon domaine: rlat,rlon,qsol=', &
390!  &             rlat(knindex(i)),rlon(knindex(i)),qsol(knindex(i))
391        qsol_new=qsol(i)
392        if (region_nudge_qsol.eq.1) then       
393           qsol_new=max(qsol(i),50.0)   
394        else if (region_nudge_qsol.eq.2) then     
395           qsol_new=max(qsol(i),120.0)
396        else !if (region_nudge_qsol.eq.1) then
397          write(*,*) 'surf_land 317: cas pas prevu'
398          write(*,*) 'region_nudge_qsol=',region_nudge_qsol
399          stop
400        endif !if (region_nudge_qsol.eq.1) then
401        if (qsol(i).gt.0.0) then
402          do ixt=1,niso
403           xtsol(ixt,i)=xtsol(ixt,i)*qsol_new/qsol(i)
404          enddo
405        else !if (qsol(i).gt.0.0) then
406          do ixt=1,niso
407           xtsol(ixt,i)=0.0
408          enddo
409        endif !if (qsol(i).gt.0.0) then
410        qsol(i)=qsol_new
411        write(*,*) 'surf_land 346: qsol_new=',qsol(i)     
412     endif ! if ((rlat(i).ge.lat_min_nudge_qsol).and.
413  enddo !do i=1,knon
414
415  END SUBROUTINE surf_land_nudge_qsol
416#endif
417
418!
419!****************************************************************************************
420
421END MODULE surf_land_mod
422!
423!****************************************************************************************
424
Note: See TracBrowser for help on using the repository browser.