source: LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/surf_land_bucket_mod.F90 @ 5441

Last change on this file since 5441 was 4143, checked in by dcugnet, 3 years ago
  • Some variables are renamed or replaced by direct equivalents:
    • iso_indnum -> tracers(:)%iso_iName
    • niso_possibles -> niso
    • iqiso -> iqIsoPha ; index_trac -> itZonIso
    • ok_iso_verif -> isoCheck
    • ntraceurs_zone -> nzone ; ntraciso -> ntiso
    • qperemin -> min_qparent ; masseqmin -> min_qmass ; ratiomin -> min_ratio
  • Some renamed variables are only aliased with the older name (using USE <module>, ONLY: <oldName> => <newName>) in routines where they are repeated many times.
  • Few hard-coded indexes are now computed (examples: ilic, iso, ivap, irneb, iq_vap, iq_liq, iso_H2O, iso_HDO, iso_HTO, iso_O17, iso_O18).
  • The IF(isoCheck) test is now embedded in the check_isotopes_seq and check_isotopes_loc routines (lighter calling).
File size: 12.1 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
111#ifdef ISO
112    integer ixt
113    REAL, DIMENSION(niso,klon) :: xtsnow_prec,xtsol_prec
114    REAL, DIMENSION(klon) :: snow_prec,qsol_prec
115    real, parameter :: t_coup = 273.15
116    real, dimension(klon) :: fq_fonte_diag
117    real, dimension(klon) :: fqfonte_diag
118    real, dimension(klon) ::  snow_evap_diag
119    real, dimension(klon) ::  fqcalving_diag
120    real max_eau_sol_diag 
121    real, dimension(klon) ::  run_off_lic_diag
122    real :: coeff_rel_diag
123!    real, dimension(klon), intent(out) ::  runoff_diag   
124#endif       
125!
126!****************************************************************************************
127
128#ifdef ISO
129#ifdef ISOVERIF
130        !write(*,*) 'surf_land_bucket 152'
131        do i=1,knon
132          if (iso_eau.gt.0) then
133            call iso_verif_egalite_choix(precip_snow(i), &
134     &          xtprecip_snow(iso_eau,i),'surf_land_bucket 131', &
135     &          errmax,errmaxrel)
136            call iso_verif_egalite_choix(qsol(i), &
137     &          xtsol(iso_eau,i),'surf_land_bucket 134', &
138     &          errmax,errmaxrel)
139          endif 
140        enddo
141#endif
142#ifdef ISOVERIF
143        do i=1,knon
144         do ixt=1,niso
145          call iso_verif_noNaN(xtsol(ixt,i),'surf_land_mod_bucket 142')
146         enddo !do ixt=1,niso
147        enddo !do i=1,knon
148        !write(*,*) 'surf_land_bucket 152'
149#endif
150#endif
151
152!
153!* Read from limit.nc : albedo(alb_lim), length of rugosity(z0_new)
154!
155    CALL limit_read_rug_alb(itime, dtime, jour,&
156         knon, knindex, &
157         z0_new, alb_lim)
158!        write(*,*) 'surf_land_bucket 166'
159!
160!* Calcultaion of fluxes
161!
162
163! calculate total absorbed radiance at surface
164       radsol(:) = 0.0
165       radsol(1:knon) = swnet(1:knon) + lwnet(1:knon)
166
167! calculate constants
168!        write(*,*) 'surf_land_bucket 176'
169    CALL calbeta(dtime, is_ter, knon, snow, qsol, beta, capsol, dif_grnd)
170    if (type_veget=='betaclim') then
171       CALL calbeta_clim(knon,jour,latitude(knindex(1:knon)),beta)
172    endif
173       
174! calculate temperature, heat capacity and conduction flux in soil
175!        write(*,*) 'surf_land_bucket 183: soil_model=',soil_model
176    IF (soil_model) THEN
177!       write(*,*) 'surf_land_bucket 185'
178       CALL soil(dtime, is_ter, knon, snow, tsurf, qsol,  &
179      & longitude(knindex(1:knon)), latitude(knindex(1:knon)), tsoil, soilcap, soilflux)
180
181!       write(*,*) 'surf_land_bucket 187'
182       DO i=1, knon
183          cal(i) = RCPD / soilcap(i)
184          radsol(i) = radsol(i)  + soilflux(i)
185       END DO
186    ELSE
187       cal(:) = RCPD * capsol(:)
188       IF (klon_glo .EQ. 1) THEN
189         cal(:) = 0.
190       ENDIF
191    ENDIF
192   
193! Suppose zero surface speed
194!        write(*,*) 'surf_land_bucket 198'
195    u0(:)=0.0
196    v0(:)=0.0
197    u1_lay(:) = u1(:) - u0(:)
198    v1_lay(:) = v1(:) - v0(:)
199
200!        write(*,*) 'surf_land_bucket 201'
201    CALL calcul_fluxs(knon, is_ter, dtime, &
202         tsurf, p1lay, cal, beta, tq_cdrag, tq_cdrag, pref, &
203         precip_rain, precip_snow, snow, qsurf,  &
204         radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, gustiness, &
205         1.,petAcoef, peqAcoef, petBcoef, peqBcoef, &
206         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
207   
208
209
210#ifdef ISO
211   ! verif
212#ifdef ISOVERIF
213    !write(*,*) 'surf_land_bucket 211'
214    do i=1,knon
215      if (iso_eau.gt.0) then
216        call iso_verif_egalite_choix(xtsnow(iso_eau,i), &
217     &           snow(i),'surf_land_bucket 522', &
218     &           errmax,errmaxrel)
219       endif !if (iso_eau.gt.0) then
220    enddo !do i=1,knon
221#endif
222   ! end verif
223#endif         
224#ifdef ISO
225    do i=1,knon
226      snow_prec(i)=snow(i)
227      qsol_prec(i)=qsol(i)
228      do ixt=1,niso
229        xtsnow_prec(ixt,i)=xtsnow(ixt,i)
230        xtsol_prec(ixt,i)=xtsol(ixt,i)
231      enddo !do ixt=1,niso
232      ! initialisation:
233      fqfonte_diag(i)=0.0
234      fq_fonte_diag(i)=0.0
235      snow_evap_diag(i)=0.0
236   enddo !do i=1,knon
237#ifdef ISOVERIF
238       ! write(*,*) 'surf_land_bucket 235'
239       do i=1,knon 
240        if (iso_eau.gt.0) then
241            call iso_verif_egalite(qsol_prec(i),xtsol_prec(iso_eau,i), &
242    &            'surf_land_bucket 141')
243        endif
244      enddo !do i=1,knon
245#endif   
246#endif   
247!
248!* Calculate snow height, run_off, age of snow
249!     
250    CALL fonte_neige( knon, is_ter, knindex, dtime, &
251         tsurf, precip_rain, precip_snow, &
252         snow, qsol, tsurf_new, evap &
253#ifdef ISO   
254     & ,fq_fonte_diag,fqfonte_diag,snow_evap_diag,fqcalving_diag   &
255     & ,max_eau_sol_diag,runoff_diag,run_off_lic_diag,coeff_rel_diag   &
256#endif
257     &   )
258
259#ifdef ISO
260#ifdef ISOVERIF
261        do i=1,knon
262         do ixt=1,niso
263          call iso_verif_noNaN(xtsol_prec(ixt,i),'surf_land_burcket 237')
264         enddo
265        enddo
266#endif
267#ifdef ISOVERIF
268        !write(*,*) 'surf_land_bucket 235'
269        do i=1,knon
270          if (iso_eau.gt.0) then
271            call iso_verif_egalite_choix(qsol_prec(i), &
272     &          xtsol_prec(iso_eau,i),'surf_land_bucket 628', &
273     &          errmax,errmaxrel)
274            call iso_verif_egalite_choix(precip_snow(i), &
275     &          xtprecip_snow(iso_eau,i),'surf_land_bucket 227', &
276     &          errmax,errmaxrel)
277             ! attention, dans fonte_neige, on modifie snow sans modifier
278             ! xtsnow
279             ! c'est fait plus tard dans gestion_neige
280!            write(*,*) 'surf_land_bucket 287: i=',i
281!            write(*,*) 'snow(i)=',snow(i)
282            call iso_verif_egalite_choix(xtsnow(iso_eau,i), &
283     &           snow_prec(i),'surf_land_bucket 245', &
284     &           errmax,errmaxrel)
285          endif 
286          if ((iso_O17.gt.0).and.(iso_O18.gt.0)) then
287              if (qsol_prec(i).gt.ridicule_qsol) then
288                call iso_verif_aberrant_o17(xtsol_prec(iso_O17,i) &
289     &           /qsol_prec(i),xtsol_prec(iso_O18,i) &
290     &           /qsol_prec(i),'surf_land_bucket 642')
291              endif !if ((qsol_prec(i).gt.ridicule_qsol) &
292          endif !if ((iso_O17.gt.0).and.(iso_O18.gt.0)) then
293        enddo  !do i=1,knon
294        !write(*,*) 'surf_land_mod 291'
295        !write(*,*) 'snow_evap_diag(1)=',snow_evap_diag(1)
296#endif         
297        call calcul_iso_surf_ter_vectall(klon,knon, &
298     &           evap,snow_evap_diag,snow, &
299     &           fq_fonte_diag,fqfonte_diag,dtime,precip_rain,xtprecip_rain, &
300     &           precip_snow,xtprecip_snow, snow_prec,xtsnow_prec, &
301     &           tsurf_new,xtspechum,pref,spechum,t_coup,u1_lay,v1_lay,p1lay, &
302     &           qsol,xtsol,qsol_prec,xtsol_prec, &
303     &           max_eau_sol_diag, &
304     &           xtevap,xtsnow,h1,runoff_diag,xtrunoff_diag,fqcalving_diag, &
305     &           knindex,is_ter,run_off_lic_diag,coeff_rel_diag,Rland_ice &
306     &   )
307!#ifdef ISOVERIF
308!        write(*,*) 'surf_land_bucket 303'
309!#endif
310#endif
311
312!
313!* Calculate the age of snow
314!
315    CALL albsno(klon,knon,dtime,agesno(:),alb_neig(:), precip_snow(:)) 
316   
317    WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
318   
319    DO i=1, knon
320       zfra(i) = MAX(0.0,MIN(1.0, snow(i)/(snow(i)+10.0)))
321       alb_lim(i)  = alb_neig(i) *zfra(i) + alb_lim(i)*(1.0-zfra(i))
322    END DO
323
324!
325!* Return albedo :
326!    alb1_new and alb2_new are here given the same values
327!
328    alb1_new(:) = 0.0
329    alb2_new(:) = 0.0
330    alb1_new(1:knon) = alb_lim(1:knon)
331    alb2_new(1:knon) = alb_lim(1:knon)
332       
333!
334!* Calculate the rugosity
335!
336    DO i = 1, knon
337       z0_new(i) = MAX(1.5e-05,SQRT(z0_new(i)**2 + rugoro(i)**2))
338    END DO
339
340!* Send to coupler
341!  The run-off from river and coast are not calculated in the bucket modele.
342!  For testing purpose of the coupled modele we put the run-off to zero.
343    IF (type_ocean=='couple') THEN
344       dummy_riverflow(:)   = 0.0
345       dummy_coastalflow(:) = 0.0
346       CALL cpl_send_land_fields(itime, knon, knindex, &
347            dummy_riverflow, dummy_coastalflow)
348    ENDIF
349
350!
351!* End
352!
353  END SUBROUTINE surf_land_bucket
354!
355!****************************************************************************************
356!
357END MODULE surf_land_bucket_mod
Note: See TracBrowser for help on using the repository browser.