source: LMDZ6/trunk/libf/phylmd/surf_land_bucket_mod.F90 @ 5051

Last change on this file since 5051 was 5022, checked in by Sebastien Nguyen, 4 months ago

include ISO keys in pbl_surface and associated routines in phylmd

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 12.0 KB
Line 
1!
2MODULE surf_land_bucket_mod
3!
4! Surface land bucket module
5!
6! This module is used when no external land model is choosen.
7!
8  IMPLICIT NONE
9
10CONTAINS
11
12  SUBROUTINE surf_land_bucket(itime, jour, knon, knindex, debut, dtime,&
13       tsurf, p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, &
14       spechum, petAcoef, peqAcoef, petBcoef, peqBcoef, pref, &
15       u1, v1, gustiness, rugoro, swnet, lwnet, &
16       snow, qsol, agesno, tsoil, &
17       qsurf, z0_new, alb1_new, alb2_new, evap, &
18       fluxsens, fluxlat, tsurf_new, dflux_s, dflux_l &
19#ifdef ISO
20       ,xtprecip_rain, xtprecip_snow,xtspechum, &
21       xtsnow, xtsol,xtevap,h1, &
22       runoff_diag,xtrunoff_diag,Rland_ice &
23#endif           
24            )
25
26    USE limit_read_mod
27    USE surface_data
28    USE fonte_neige_mod
29    USE calcul_fluxs_mod
30    USE cpl_mod
31    USE dimphy
32    USE geometry_mod, ONLY: longitude,latitude
33    USE mod_grid_phy_lmdz
34    USE mod_phys_lmdz_para
35    USE indice_sol_mod
36#ifdef ISO
37    use infotrac_phy, ONLY: ntiso,niso
38    USE isotopes_mod, ONLY: iso_eau, iso_HDO, iso_O18, iso_O17, &
39        ridicule_qsol
40    USE isotopes_routines_mod, ONLY: calcul_iso_surf_ter_vectall
41#ifdef ISOVERIF
42    USE isotopes_verif_mod, ONLY: errmax,errmaxrel,iso_verif_noNaN, &
43        iso_verif_aberrant_o17,iso_verif_egalite_choix,iso_verif_egalite
44#endif
45#endif
46!****************************************************************************************
47! Bucket calculations for surface.
48!
49    INCLUDE "clesphys.h"
50    INCLUDE "dimsoil.h"
51    INCLUDE "YOMCST.h"
52
53! Input variables 
54!****************************************************************************************
55    INTEGER, INTENT(IN)                     :: itime, jour, knon
56    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex
57    LOGICAL, INTENT(IN)                     :: debut
58    REAL, INTENT(IN)                        :: dtime
59    REAL, DIMENSION(klon), INTENT(IN)       :: tsurf
60    REAL, DIMENSION(klon), INTENT(IN)       :: p1lay
61    REAL, DIMENSION(klon), INTENT(IN)       :: tq_cdrag
62    REAL, DIMENSION(klon), INTENT(IN)       :: precip_rain, precip_snow
63    REAL, DIMENSION(klon), INTENT(IN)       :: temp_air, spechum
64    REAL, DIMENSION(klon), INTENT(IN)       :: petAcoef, peqAcoef
65    REAL, DIMENSION(klon), INTENT(IN)       :: petBcoef, peqBcoef
66    REAL, DIMENSION(klon), INTENT(IN)       :: pref
67    REAL, DIMENSION(klon), INTENT(IN)       :: u1, v1, gustiness
68    REAL, DIMENSION(klon), INTENT(IN)       :: rugoro
69    REAL, DIMENSION(klon), INTENT(IN)       :: swnet, lwnet
70#ifdef ISO
71    REAL, DIMENSION(ntiso,klon), INTENT(IN) :: xtprecip_rain, xtprecip_snow
72    REAL, DIMENSION(ntiso,klon), INTENT(IN) :: xtspechum   
73#endif
74
75! In/Output variables
76!****************************************************************************************
77    REAL, DIMENSION(klon), INTENT(INOUT)          :: snow, qsol
78    REAL, DIMENSION(klon), INTENT(INOUT)          :: agesno
79    REAL, DIMENSION(klon, nsoilmx), INTENT(INOUT) :: tsoil
80#ifdef ISO
81    REAL, DIMENSION(niso,klon), INTENT(INOUT)       :: xtsnow,xtsol
82#endif
83
84! Output variables
85!****************************************************************************************
86    REAL, DIMENSION(klon), INTENT(OUT)       :: qsurf
87    REAL, DIMENSION(klon), INTENT(OUT)       :: z0_new
88    REAL, DIMENSION(klon), INTENT(OUT)       :: alb1_new, alb2_new
89    REAL, DIMENSION(klon), INTENT(OUT)       :: evap, fluxsens, fluxlat
90    REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
91    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l     
92#ifdef ISO
93    REAL, DIMENSION(ntiso,klon), INTENT(OUT) :: xtevap
94    REAL, DIMENSION(klon),       INTENT(OUT) :: h1
95    REAL, DIMENSION(niso,klon),  INTENT(OUT) :: xtrunoff_diag
96    REAL, DIMENSION(klon),       INTENT(OUT) :: runoff_diag
97    REAL, DIMENSION(niso,klon),  INTENT(IN)  :: Rland_ice
98#endif
99
100! Local variables
101!****************************************************************************************
102    REAL, DIMENSION(klon) :: soilcap, soilflux
103    REAL, DIMENSION(klon) :: cal, beta, dif_grnd, capsol
104    REAL, DIMENSION(klon) :: alb_neig, alb_lim
105    REAL, DIMENSION(klon) :: zfra
106    REAL, DIMENSION(klon) :: radsol       ! total net radiance at surface
107    REAL, DIMENSION(klon) :: u0, v0, u1_lay, v1_lay
108    REAL, DIMENSION(klon) :: dummy_riverflow,dummy_coastalflow
109    INTEGER               :: i
110#ifdef ISO
111    INTEGER               :: ixt
112    REAL, DIMENSION(niso,klon) :: xtsnow_prec,xtsol_prec
113    REAL, DIMENSION(klon) :: snow_prec,qsol_prec
114    REAL, PARAMETER       :: t_coup = 273.15
115    REAL, DIMENSION(klon) :: fq_fonte_diag
116    REAL, DIMENSION(klon) :: fqfonte_diag
117    REAL, DIMENSION(klon) :: snow_evap_diag
118    REAL, DIMENSION(klon) :: fqcalving_diag
119    REAL                  :: max_eau_sol_diag 
120    REAL, DIMENSION(klon) :: run_off_lic_diag
121    REAL :: coeff_rel_diag
122#endif
123!
124!****************************************************************************************
125
126#ifdef ISO
127#ifdef ISOVERIF
128        !write(*,*) 'surf_land_bucket 152'
129        DO i=1,knon
130          IF (iso_eau > 0) THEN
131            CALL iso_verif_egalite_choix(precip_snow(i), &
132     &                                   xtprecip_snow(iso_eau,i),'surf_land_bucket 131', &
133     &                                   errmax,errmaxrel)
134            CALL iso_verif_egalite_choix(qsol(i), &
135     &                                   xtsol(iso_eau,i),'surf_land_bucket 134', &
136     &                                   errmax,errmaxrel)
137          ENDIF
138        ENDDO
139#endif
140#ifdef ISOVERIF
141        DO i=1,knon
142          DO ixt=1,niso
143            CALL iso_verif_noNaN(xtsol(ixt,i),'surf_land_mod_bucket 142')
144          ENDDO !do ixt=1,niso
145        ENDDO !do i=1,knon
146        !write(*,*) 'surf_land_bucket 152'
147#endif
148#endif
149
150!
151!* Read from limit.nc : albedo(alb_lim), length of rugosity(z0_new)
152!
153    CALL limit_read_rug_alb(itime, dtime, jour,&
154         knon, knindex, &
155         z0_new, alb_lim)
156!
157!* Calcultaion of fluxes
158!
159
160! calculate total absorbed radiance at surface
161       radsol(:) = 0.0
162       radsol(1:knon) = swnet(1:knon) + lwnet(1:knon)
163
164! calculate constants
165    CALL calbeta(dtime, is_ter, knon, snow, qsol, beta, capsol, dif_grnd)
166    if (type_veget=='betaclim') then
167       CALL calbeta_clim(knon,jour,latitude(knindex(1:knon)),beta)
168    endif
169       
170! calculate temperature, heat capacity and conduction flux in soil
171    IF (soil_model) THEN
172       CALL soil(dtime, is_ter, knon, snow, tsurf, qsol,  &
173        & longitude(knindex(1:knon)), latitude(knindex(1:knon)), tsoil, soilcap, soilflux)
174
175       DO i=1, knon
176          cal(i) = RCPD / soilcap(i)
177          radsol(i) = radsol(i)  + soilflux(i)
178       END DO
179    ELSE
180       cal(:) = RCPD * capsol(:)
181       IF (klon_glo .EQ. 1) THEN
182         cal(:) = 0.
183       ENDIF
184    ENDIF
185   
186! Suppose zero surface speed
187    u0(:)=0.0
188    v0(:)=0.0
189    u1_lay(:) = u1(:) - u0(:)
190    v1_lay(:) = v1(:) - v0(:)
191
192    CALL calcul_fluxs(knon, is_ter, dtime, &
193         tsurf, p1lay, cal, beta, tq_cdrag, tq_cdrag, pref, &
194         precip_rain, precip_snow, snow, qsurf,  &
195         radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, gustiness, &
196         1.,petAcoef, peqAcoef, petBcoef, peqBcoef, &
197         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
198   
199#ifdef ISO
200   ! verif
201#ifdef ISOVERIF
202    !write(*,*) 'surf_land_bucket 211'
203    DO i=1,knon
204      IF (iso_eau > 0) THEN
205        CALL iso_verif_egalite_choix(xtsnow(iso_eau,i), &
206     &           snow(i),'surf_land_bucket 522', &
207     &           errmax,errmaxrel)
208      ENDIF !IF (iso_eau > 0) then
209    ENDDO !DO i=1,knon
210#endif
211   ! end verif
212#endif         
213#ifdef ISO
214    DO i=1,knon
215      snow_prec(i)=snow(i)
216      qsol_prec(i)=qsol(i)
217      DO ixt=1,niso
218        xtsnow_prec(ixt,i)=xtsnow(ixt,i)
219        xtsol_prec(ixt,i) =xtsol(ixt,i)
220      ENDDO !DO ixt=1,niso
221      ! initialisation:
222      fqfonte_diag(i)  =0.0
223      fq_fonte_diag(i) =0.0
224      snow_evap_diag(i)=0.0
225    ENDDO !DO i=1,knon
226#ifdef ISOVERIF
227    ! write(*,*) 'surf_land_bucket 235'
228    DO i=1,knon 
229      IF (iso_eau > 0) THEN
230        CALL iso_verif_egalite(qsol_prec(i),xtsol_prec(iso_eau,i), &
231    &                              'surf_land_bucket 141')
232      ENDIF
233    ENDDO !DO i=1,knon
234#endif   
235#endif   
236!
237!* Calculate snow height, run_off, age of snow
238!     
239    CALL fonte_neige( knon, is_ter, knindex, dtime, &
240         tsurf, precip_rain, precip_snow, &
241         snow, qsol, tsurf_new, evap &
242#ifdef ISO   
243     & ,fq_fonte_diag,fqfonte_diag,snow_evap_diag,fqcalving_diag   &
244     & ,max_eau_sol_diag,runoff_diag,run_off_lic_diag,coeff_rel_diag   &
245#endif
246     &   )
247
248#ifdef ISO
249#ifdef ISOVERIF
250        DO i=1,knon
251          DO ixt=1,niso
252            CALL iso_verif_noNaN(xtsol_prec(ixt,i),'surf_land_burcket 237')
253          ENDDO
254        ENDDO
255#endif
256#ifdef ISOVERIF
257        !write(*,*) 'surf_land_bucket 235'
258        DO i=1,knon
259          IF (iso_eau > 0) THEN
260            CALL iso_verif_egalite_choix(qsol_prec(i), &
261     &                                   xtsol_prec(iso_eau,i),'surf_land_bucket 628', &
262     &                                   errmax,errmaxrel)
263            CALL iso_verif_egalite_choix(precip_snow(i), &
264     &                                   xtprecip_snow(iso_eau,i),'surf_land_bucket 227', &
265     &                                   errmax,errmaxrel)
266             ! attention, dans fonte_neige, on modifie snow sans modifier
267             ! xtsnow
268             ! c'est fait plus tard dans gestion_neige
269!            write(*,*) 'surf_land_bucket 287: i=',i
270!            write(*,*) 'snow(i)=',snow(i)
271            CALL iso_verif_egalite_choix(xtsnow(iso_eau,i), &
272     &                                   snow_prec(i),'surf_land_bucket 245', &
273     &                                   errmax,errmaxrel)
274          ENDIF 
275          IF ((iso_O17 > 0).AND.(iso_O18 > 0)) THEN
276              IF (qsol_prec(i) > ridicule_qsol) THEN
277                CALL iso_verif_aberrant_o17(xtsol_prec(iso_O17,i)/qsol_prec(i) &
278     &                                     ,xtsol_prec(iso_O18,i)/qsol_prec(i) &
279     &                                     ,'surf_land_bucket 642')
280              ENDIF !IF ((qsol_prec(i) > ridicule_qsol) &
281          ENDIF !IF ((iso_O17 > 0).AND.(iso_O18 > 0)) THEN
282        ENDDO  !DO i=1,knon
283        !write(*,*) 'surf_land_mod 291'
284        !write(*,*) 'snow_evap_diag(1)=',snow_evap_diag(1)
285#endif         
286        CALL calcul_iso_surf_ter_vectall(klon,knon, &
287     &           evap,snow_evap_diag,snow, &
288     &           fq_fonte_diag,fqfonte_diag,dtime,precip_rain,xtprecip_rain, &
289     &           precip_snow,xtprecip_snow, snow_prec,xtsnow_prec, &
290     &           tsurf_new,xtspechum,pref,spechum,t_coup,u1_lay,v1_lay,p1lay, &
291     &           qsol,xtsol,qsol_prec,xtsol_prec, &
292     &           max_eau_sol_diag, &
293     &           xtevap,xtsnow,h1,runoff_diag,xtrunoff_diag,fqcalving_diag, &
294     &           knindex,is_ter,run_off_lic_diag,coeff_rel_diag,Rland_ice &
295     &   )
296!#ifdef ISOVERIF
297!        write(*,*) 'surf_land_bucket 303'
298!#endif
299#endif
300
301!
302!* Calculate the age of snow
303!
304    CALL albsno(klon,knon,dtime,agesno(:),alb_neig(:), precip_snow(:)) 
305   
306    WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
307   
308    DO i=1, knon
309       zfra(i) = MAX(0.0,MIN(1.0, snow(i)/(snow(i)+10.0)))
310       alb_lim(i)  = alb_neig(i) *zfra(i) + alb_lim(i)*(1.0-zfra(i))
311    END DO
312
313!
314!* Return albedo :
315!    alb1_new and alb2_new are here given the same values
316!
317    alb1_new(:) = 0.0
318    alb2_new(:) = 0.0
319    alb1_new(1:knon) = alb_lim(1:knon)
320    alb2_new(1:knon) = alb_lim(1:knon)
321       
322!
323!* Calculate the rugosity
324!
325    DO i = 1, knon
326       z0_new(i) = MAX(1.5e-05,SQRT(z0_new(i)**2 + rugoro(i)**2))
327    END DO
328
329!* Send to coupler
330!  The run-off from river and coast are not calculated in the bucket modele.
331!  For testing purpose of the coupled modele we put the run-off to zero.
332    IF (type_ocean=='couple') THEN
333       dummy_riverflow(:)   = 0.0
334       dummy_coastalflow(:) = 0.0
335       CALL cpl_send_land_fields(itime, knon, knindex, &
336            dummy_riverflow, dummy_coastalflow)
337    ENDIF
338
339!
340!* End
341!
342  END SUBROUTINE surf_land_bucket
343!
344!****************************************************************************************
345!
346END MODULE surf_land_bucket_mod
Note: See TracBrowser for help on using the repository browser.