/[lmdze]/trunk/phylmd/Interface_surf/conf_interface.f
ViewVC logotype

Diff of /trunk/phylmd/Interface_surf/conf_interface.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 14 by guez, Mon Jul 28 14:48:09 2008 UTC revision 38 by guez, Thu Jan 6 17:52:19 2011 UTC
# Line 2  MODULE interface_surf Line 2  MODULE interface_surf
2    
3    ! From phylmd/interface_surf.F90, version 1.8 2005/05/25 13:10:09    ! From phylmd/interface_surf.F90, version 1.8 2005/05/25 13:10:09
4    
5    ! Ce module regroupe toutes les routines gérant l'interface entre le modèle    ! Ce module regroupe toutes les routines gérant l'interface entre le
6    ! atmosphérique et les modèles de surface (sols continentaux,    ! modèle atmosphérique et les modèles de surface (sols continentaux,
7    ! océans, glaces).    ! océans, glaces). Les routines sont les suivantes. "interfsurf_hq" :
8    ! Les routines sont les suivantes:    ! routine d'aiguillage vers les interfaces avec les différents
9    ! interfsurf_hq : routine d'aiguillage vers les interfaces avec les    ! modèles de surface ; "interfoce_*" : routines d'interface proprement
10    !    différents modèles de surface    ! dites.
   ! interfoce_* : routines d'interface proprement dites  
11    
12    ! L. Fairhead, LMD, 02/2000    ! L. Fairhead, LMD, february 2000
13    
14    IMPLICIT none    IMPLICIT none
15    
16    PRIVATE    PRIVATE
17    PUBLIC :: interfsurf_hq    PUBLIC :: interfsurf_hq
18    
19    ! run_off      ruissellement total    ! run_off ruissellement total
20    REAL, ALLOCATABLE, DIMENSION(:), SAVE    :: run_off, run_off_lic    REAL, ALLOCATABLE, DIMENSION(:), SAVE :: run_off, run_off_lic
21    real, allocatable, dimension(:), save    :: coastalflow, riverflow    real, allocatable, dimension(:), save :: coastalflow, riverflow
22    
23    REAL, ALLOCATABLE, DIMENSION(:, :), SAVE :: tmp_rriv, tmp_rcoa, tmp_rlic    REAL, ALLOCATABLE, DIMENSION(:, :), SAVE :: tmp_rriv, tmp_rcoa, tmp_rlic
24    !! pour simuler la fonte des glaciers antarctiques    ! pour simuler la fonte des glaciers antarctiques
25    REAL, ALLOCATABLE, DIMENSION(:, :), SAVE :: coeff_iceberg    REAL, ALLOCATABLE, DIMENSION(:, :), SAVE :: coeff_iceberg
26    real, save                              :: surf_maille    real, save :: surf_maille
27    real, save                              :: cte_flux_iceberg = 6.3e7    real, save :: cte_flux_iceberg = 6.3e7
28    integer, save                           :: num_antarctic = 1    integer, save :: num_antarctic = 1
29    REAL, save                              :: tau_calv    REAL, save :: tau_calv
30    
31  CONTAINS  CONTAINS
32    
33    SUBROUTINE interfsurf_hq(itime, dtime, date0, jour, rmu0, &    SUBROUTINE interfsurf_hq(itime, dtime, date0, jour, rmu0,  &
34         & klon, iim, jjm, nisurf, knon, knindex, pctsrf, &         klon, iim, jjm, nisurf, knon, knindex, pctsrf,  &
35         & rlon, rlat, cufi, cvfi, &         rlon, rlat, cufi, cvfi,  &
36         & debut, lafin, ok_veget, soil_model, nsoilmx, tsoil, qsol, &         debut, lafin, ok_veget, soil_model, nsoilmx, tsoil, qsol,  &
37         & zlev,  u1_lay, v1_lay, temp_air, spechum, epot_air, ccanopy, &         zlev, u1_lay, v1_lay, temp_air, spechum, epot_air, ccanopy, &
38         & tq_cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &         tq_cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef,  &
39         & precip_rain, precip_snow, sollw, sollwdown, swnet, swdown, &         precip_rain, precip_snow, sollw, sollwdown, swnet, swdown,  &
40         & fder, taux, tauy, &         fder, taux, tauy,  &
41         & windsp, &         windsp,  &
42         & rugos, rugoro, &         rugos, rugoro,  &
43         & albedo, snow, qsurf, &         albedo, snow, qsurf,  &
44         & tsurf, p1lay, ps, radsol, &         tsurf, p1lay, ps, radsol,  &
45         & ocean, npas, nexca, zmasq, &         ocean, npas, nexca, zmasq,  &
46         & evap, fluxsens, fluxlat, dflux_l, dflux_s, &                       evap, fluxsens, fluxlat, dflux_l, dflux_s, &
47         & tsol_rad, tsurf_new, alb_new, alblw, emis_new, &         tsol_rad, tsurf_new, alb_new, alblw, emis_new,  &
48         & z0_new, pctsrf_new, agesno, fqcalving, ffonte, run_off_lic_0, &         z0_new, pctsrf_new, agesno, fqcalving, ffonte, run_off_lic_0, &
49         !IM "slab" ocean         flux_o, flux_g, tslab, seaice)
        & flux_o, flux_g, tslab, seaice)  
50    
51      ! Cette routine sert d'aiguillage entre l'atmosphère et la surface      ! Cette routine sert d'aiguillage entre l'atmosphère et la surface
52      ! en général (sols continentaux, océans, glaces) pour les flux de      ! en général (sols continentaux, océans, glaces) pour les flux de
# Line 62  CONTAINS Line 60  CONTAINS
60      use abort_gcm_m, only: abort_gcm      use abort_gcm_m, only: abort_gcm
61      use gath_cpl, only: gath2cpl      use gath_cpl, only: gath2cpl
62      use indicesol      use indicesol
63      use YOMCST      use SUPHEC_M
64      use albsno_m, only: albsno      use albsno_m, only: albsno
65    
66      ! Parametres d'entree      ! Parametres d'entree
67      ! input:      ! input:
68      !   klon         nombre total de points de grille      ! klon nombre total de points de grille
69      !   iim, jjm     nbres de pts de grille      ! iim, jjm nbres de pts de grille
70      !   dtime        pas de temps de la physique (en s)      ! dtime pas de temps de la physique (en s)
71      !   date0        jour initial      ! date0 jour initial
72      !   jour         jour dans l'annee en cours,      ! jour jour dans l'annee en cours,
73      !   rmu0         cosinus de l'angle solaire zenithal      ! rmu0 cosinus de l'angle solaire zenithal
74      !   nexca        pas de temps couplage      ! nexca pas de temps couplage
75      !   nisurf       index de la surface a traiter (1 = sol continental)      ! nisurf index de la surface a traiter (1 = sol continental)
76      !   knon         nombre de points de la surface a traiter      ! knon nombre de points de la surface a traiter
77      !   knindex      index des points de la surface a traiter      ! knindex index des points de la surface a traiter
78      !   pctsrf       tableau des pourcentages de surface de chaque maille      ! pctsrf tableau des pourcentages de surface de chaque maille
79      !   rlon         longitudes      ! rlon longitudes
80      !   rlat         latitudes      ! rlat latitudes
81      !   cufi, cvfi    resolution des mailles en x et y (m)      ! cufi, cvfi resolution des mailles en x et y (m)
82      !   debut        logical: 1er appel a la physique      ! debut logical: 1er appel a la physique
83      !   lafin        logical: dernier appel a la physique      ! lafin logical: dernier appel a la physique
84      !   ok_veget     logical: appel ou non au schema de surface continental      ! ok_veget logical: appel ou non au schema de surface continental
85      !                  (si false calcul simplifie des fluxs sur les continents)      ! (si false calcul simplifie des fluxs sur les continents)
86      !   zlev         hauteur de la premiere couche      ! zlev hauteur de la premiere couche
87      !   u1_lay       vitesse u 1ere couche      ! u1_lay vitesse u 1ere couche
88      !   v1_lay       vitesse v 1ere couche      ! v1_lay vitesse v 1ere couche
89      !   temp_air     temperature de l'air 1ere couche      ! temp_air temperature de l'air 1ere couche
90      !   spechum      humidite specifique 1ere couche      ! spechum humidite specifique 1ere couche
91      !   epot_air     temp potentielle de l'air      ! epot_air temp potentielle de l'air
92      !   ccanopy      concentration CO2 canopee      ! ccanopy concentration CO2 canopee
93      !   tq_cdrag     cdrag      ! tq_cdrag cdrag
94      !   petAcoef     coeff. A de la resolution de la CL pour t      ! petAcoef coeff. A de la resolution de la CL pour t
95      !   peqAcoef     coeff. A de la resolution de la CL pour q      ! peqAcoef coeff. A de la resolution de la CL pour q
96      !   petBcoef     coeff. B de la resolution de la CL pour t      ! petBcoef coeff. B de la resolution de la CL pour t
97      !   peqBcoef     coeff. B de la resolution de la CL pour q      ! peqBcoef coeff. B de la resolution de la CL pour q
98      !   precip_rain  precipitation liquide      ! precip_rain precipitation liquide
99      !   precip_snow  precipitation solide      ! precip_snow precipitation solide
100      !   sollw        flux IR net a la surface      ! sollw flux IR net a la surface
101      !   sollwdown    flux IR descendant a la surface      ! sollwdown flux IR descendant a la surface
102      !   swnet        flux solaire net      ! swnet flux solaire net
103      !   swdown       flux solaire entrant a la surface      ! swdown flux solaire entrant a la surface
104      !   albedo       albedo de la surface      ! albedo albedo de la surface
105      !   tsurf        temperature de surface      ! tsurf temperature de surface
106      !   tslab        temperature slab ocean      ! tslab temperature slab ocean
107      !   pctsrf_slab  pourcentages (0-1) des sous-surfaces dans le slab      ! pctsrf_slab pourcentages (0-1) des sous-surfaces dans le slab
108      !   tmp_pctsrf_slab = pctsrf_slab      ! tmp_pctsrf_slab = pctsrf_slab
109      !   p1lay        pression 1er niveau (milieu de couche)      ! p1lay pression 1er niveau (milieu de couche)
110      !   ps           pression au sol      ! ps pression au sol
111      !   radsol       rayonnement net aus sol (LW + SW)      ! radsol rayonnement net aus sol (LW + SW)
112      !   ocean        type d'ocean utilise ("force" ou "slab" mais pas "couple")      ! ocean type d'ocean utilise ("force" ou "slab" mais pas "couple")
113      !   fder         derivee des flux (pour le couplage)      ! fder derivee des flux (pour le couplage)
114      !   taux, tauy   tension de vents      ! taux, tauy tension de vents
115      !   windsp       module du vent a 10m      ! windsp module du vent a 10m
116      !   rugos        rugosite      ! rugos rugosite
117      !   zmasq        masque terre/ocean      ! zmasq masque terre/ocean
118      !   rugoro       rugosite orographique      ! rugoro rugosite orographique
119      !   run_off_lic_0 runoff glacier du pas de temps precedent      ! run_off_lic_0 runoff glacier du pas de temps precedent
120      integer, intent(IN) :: itime     !  numero du pas de temps      integer, intent(IN) :: itime ! numero du pas de temps
121      integer, intent(IN) :: iim, jjm      integer, intent(IN) :: iim, jjm
122      integer, intent(IN) :: klon      integer, intent(IN) :: klon
123      real, intent(IN) :: dtime      real, intent(IN) :: dtime
124      real, intent(IN) :: date0      real, intent(IN) :: date0
125      integer, intent(IN) :: jour      integer, intent(IN) :: jour
126      real, intent(IN)    :: rmu0(klon)      real, intent(IN) :: rmu0(klon)
127      integer, intent(IN) :: nisurf      integer, intent(IN) :: nisurf
128      integer, intent(IN) :: knon      integer, intent(IN) :: knon
129      integer, dimension(klon), intent(in) :: knindex      integer, dimension(klon), intent(in) :: knindex
# Line 148  CONTAINS Line 146  CONTAINS
146      real, dimension(klon), intent(INOUT) :: tslab      real, dimension(klon), intent(INOUT) :: tslab
147      real, allocatable, dimension(:), save :: tmp_tslab      real, allocatable, dimension(:), save :: tmp_tslab
148      real, dimension(klon), intent(OUT) :: flux_o, flux_g      real, dimension(klon), intent(OUT) :: flux_o, flux_g
149      real, dimension(klon), intent(INOUT)      :: seaice ! glace de mer (kg/m2)      real, dimension(klon), intent(INOUT) :: seaice ! glace de mer (kg/m2)
150      REAL, DIMENSION(klon), INTENT(INOUT) :: radsol, fder      REAL, DIMENSION(klon), INTENT(INOUT) :: radsol, fder
151      real, dimension(klon), intent(IN) :: zmasq      real, dimension(klon), intent(IN) :: zmasq
152      real, dimension(klon), intent(IN) :: taux, tauy, rugos, rugoro      real, dimension(klon), intent(IN) :: taux, tauy, rugos, rugoro
153      real, dimension(klon), intent(IN) :: windsp      real, dimension(klon), intent(IN) :: windsp
154      character(len=*), intent(IN):: ocean      character(len=*), intent(IN):: ocean
155      integer              :: npas, nexca ! nombre et pas de temps couplage      integer :: npas, nexca ! nombre et pas de temps couplage
156      real, dimension(klon), intent(INOUT) :: evap, snow, qsurf      real, dimension(klon), intent(INOUT) :: evap, snow, qsurf
157      !! PB ajout pour soil      !! PB ajout pour soil
158      logical, intent(in):: soil_model      logical, intent(in):: soil_model
159      integer          :: nsoilmx      integer :: nsoilmx
160      REAL, DIMENSION(klon, nsoilmx) :: tsoil      REAL, DIMENSION(klon, nsoilmx) :: tsoil
161      REAL, dimension(klon), intent(INOUT) :: qsol      REAL, dimension(klon), intent(INOUT) :: qsol
162      REAL, dimension(klon)          :: soilcap      REAL, dimension(klon) :: soilcap
163      REAL, dimension(klon)          :: soilflux      REAL, dimension(klon) :: soilflux
164    
165      ! Parametres de sortie      ! Parametres de sortie
166      ! output:      ! output:
167      !   evap         evaporation totale      ! evap evaporation totale
168      !   fluxsens     flux de chaleur sensible      ! fluxsens flux de chaleur sensible
169      !   fluxlat      flux de chaleur latente      ! fluxlat flux de chaleur latente
170      !   tsol_rad          ! tsol_rad
171      !   tsurf_new    temperature au sol      ! tsurf_new temperature au sol
172      !   alb_new      albedo      ! alb_new albedo
173      !   emis_new     emissivite      ! emis_new emissivite
174      !   z0_new       surface roughness      ! z0_new surface roughness
175      !   pctsrf_new   nouvelle repartition des surfaces      ! pctsrf_new nouvelle repartition des surfaces
176      real, dimension(klon), intent(OUT):: fluxsens, fluxlat      real, dimension(klon), intent(OUT):: fluxsens, fluxlat
177      real, dimension(klon), intent(OUT):: tsol_rad, tsurf_new, alb_new      real, dimension(klon), intent(OUT):: tsol_rad, tsurf_new, alb_new
178      real, dimension(klon), intent(OUT):: alblw      real, dimension(klon), intent(OUT):: alblw
# Line 185  CONTAINS Line 183  CONTAINS
183      real, dimension(klon), intent(INOUT):: run_off_lic_0      real, dimension(klon), intent(INOUT):: run_off_lic_0
184    
185      ! Flux thermique utiliser pour fondre la neige      ! Flux thermique utiliser pour fondre la neige
186      !jld a rajouter   real, dimension(klon), intent(INOUT):: ffonte      !jld a rajouter real, dimension(klon), intent(INOUT):: ffonte
187      real, dimension(klon), intent(INOUT):: ffonte      real, dimension(klon), intent(INOUT):: ffonte
188      ! Flux d'eau "perdue" par la surface et nécessaire pour que limiter la      ! Flux d'eau "perdue" par la surface et nécessaire pour que limiter la
189      ! hauteur de neige, en kg/m2/s      ! hauteur de neige, en kg/m2/s
190      !jld a rajouter   real, dimension(klon), intent(INOUT):: fqcalving      !jld a rajouter real, dimension(klon), intent(INOUT):: fqcalving
191      real, dimension(klon), intent(INOUT):: fqcalving      real, dimension(klon), intent(INOUT):: fqcalving
192      !IM: "slab" ocean - Local      !IM: "slab" ocean - Local
193      real, parameter :: t_grnd=271.35      real, parameter :: t_grnd=271.35
# Line 203  CONTAINS Line 201  CONTAINS
201      ! Local      ! Local
202      character (len = 20), save :: modname = 'interfsurf_hq'      character (len = 20), save :: modname = 'interfsurf_hq'
203      character (len = 80) :: abort_message      character (len = 80) :: abort_message
204      logical, save        :: first_call = .true.      logical, save :: first_call = .true.
205      integer, save        :: error      integer, save :: error
206      integer              :: ii      integer :: ii
207      logical, save              :: check = .false.      logical, save :: check = .false.
208      real, dimension(klon):: cal, beta, dif_grnd, capsol      real, dimension(klon):: cal, beta, dif_grnd, capsol
209      real, parameter      :: calice=1.0/(5.1444e+06*0.15), tau_gl=86400.*5.      real, parameter :: calice=1.0/(5.1444e+06*0.15), tau_gl=86400.*5.
210      real, parameter      :: calsno=1./(2.3867e+06*.15)      real, parameter :: calsno=1./(2.3867e+06*.15)
211      real, dimension(klon):: tsurf_temp      real, dimension(klon):: tsurf_temp
212      real, dimension(klon):: alb_neig, alb_eau      real, dimension(klon):: alb_neig, alb_eau
213      real, DIMENSION(klon):: zfra      real, DIMENSION(klon):: zfra
214      logical              :: cumul = .false.      logical :: cumul = .false.
215      INTEGER, dimension(1) :: iloc      INTEGER, dimension(1) :: iloc
216      real, dimension(klon):: fder_prev      real, dimension(klon):: fder_prev
217      REAL, dimension(klon) :: bidule      REAL, dimension(klon) :: bidule
# Line 275  CONTAINS Line 273  CONTAINS
273            call abort_gcm(modname, abort_message, 1)            call abort_gcm(modname, abort_message, 1)
274         endif         endif
275      endif      endif
276      if (.not. allocated(tmp_flux_g)) then        if (.not. allocated(tmp_flux_g)) then
277         allocate(tmp_flux_g(klon), stat = error)         allocate(tmp_flux_g(klon), stat = error)
278         DO i=1, knon         DO i=1, knon
279            tmp_flux_g(knindex(i))=flux_g(i)            tmp_flux_g(knindex(i))=flux_g(i)
# Line 285  CONTAINS Line 283  CONTAINS
283            call abort_gcm(modname, abort_message, 1)            call abort_gcm(modname, abort_message, 1)
284         endif         endif
285      endif      endif
286      if (.not. allocated(tmp_radsol)) then        if (.not. allocated(tmp_radsol)) then
287         allocate(tmp_radsol(klon), stat = error)         allocate(tmp_radsol(klon), stat = error)
288         if (error /= 0) then         if (error /= 0) then
289            abort_message='Pb allocation tmp_radsol'            abort_message='Pb allocation tmp_radsol'
# Line 351  CONTAINS Line 349  CONTAINS
349               abort_message='Pb allocation run_off'               abort_message='Pb allocation run_off'
350               call abort_gcm(modname, abort_message, 1)               call abort_gcm(modname, abort_message, 1)
351            endif            endif
352            !cym                  !cym
353            run_off=0.0            run_off=0.0
354            !cym            !cym
355    
# Line 400  CONTAINS Line 398  CONTAINS
398               CALL soil(dtime, nisurf, knon, snow, tsurf, tsoil, soilcap, &               CALL soil(dtime, nisurf, knon, snow, tsurf, tsoil, soilcap, &
399                    soilflux)                    soilflux)
400               cal(1:knon) = RCPD / soilcap(1:knon)               cal(1:knon) = RCPD / soilcap(1:knon)
401               radsol(1:knon) = radsol(1:knon)  + soilflux(1:knon)               radsol(1:knon) = radsol(1:knon) + soilflux(1:knon)
402            ELSE            ELSE
403               cal = RCPD * capsol               cal = RCPD * capsol
404            ENDIF            ENDIF
405            CALL calcul_fluxs( klon, knon, nisurf, dtime, &            CALL calcul_fluxs( klon, knon, nisurf, dtime, &
406                 &   tsurf, p1lay, cal, beta, tq_cdrag, ps, &                 tsurf, p1lay, cal, beta, tq_cdrag, ps, &
407                 &   precip_rain, precip_snow, snow, qsurf,  &                 precip_rain, precip_snow, snow, qsurf, &
408                 &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &                 radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
409                 &   petAcoef, peqAcoef, petBcoef, peqBcoef, &                 petAcoef, peqAcoef, petBcoef, peqBcoef, &
410                 &   tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)                 tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
411    
412            CALL fonte_neige( klon, knon, nisurf, dtime, &            CALL fonte_neige( klon, knon, nisurf, dtime, &
413                 &   tsurf, p1lay, cal, beta, tq_cdrag, ps, &                 tsurf, p1lay, cal, beta, tq_cdrag, ps, &
414                 &   precip_rain, precip_snow, snow, qsol,  &                 precip_rain, precip_snow, snow, qsol, &
415                 &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &                 radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
416                 &   petAcoef, peqAcoef, petBcoef, peqBcoef, &                 petAcoef, peqAcoef, petBcoef, peqBcoef, &
417                 &   tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &                 tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &
418                 &   fqcalving, ffonte, run_off_lic_0)                 fqcalving, ffonte, run_off_lic_0)
419    
420            call albsno(klon, knon, dtime, agesno, alb_neig, precip_snow)            call albsno(klon, knon, dtime, agesno, alb_neig, precip_snow)
421            where (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.            where (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
422            zfra(1:knon) = max(0.0, min(1.0, snow(1:knon)/(snow(1:knon)+10.0)))            zfra(1:knon) = max(0.0, min(1.0, snow(1:knon)/(snow(1:knon)+10.0)))
423            alb_new(1 : knon)  = alb_neig(1 : knon) *zfra(1:knon) + &            alb_new(1 : knon) = alb_neig(1 : knon) *zfra(1:knon) + &
424                 alb_new(1 : knon)*(1.0-zfra(1:knon))                 alb_new(1 : knon)*(1.0-zfra(1:knon))
425            z0_new = sqrt(z0_new**2+rugoro**2)            z0_new = sqrt(z0_new**2+rugoro**2)
426            alblw(1 : knon) = alb_new(1 : knon)            alblw(1 : knon) = alb_new(1 : knon)
# Line 449  CONTAINS Line 447  CONTAINS
447         agesno = 0.         agesno = 0.
448    
449         call calcul_fluxs( klon, knon, nisurf, dtime, &         call calcul_fluxs( klon, knon, nisurf, dtime, &
450              &   tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &              tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
451              &   precip_rain, precip_snow, snow, qsurf,  &              precip_rain, precip_snow, snow, qsurf, &
452              &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &              radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
453              &   petAcoef, peqAcoef, petBcoef, peqBcoef, &              petAcoef, peqAcoef, petBcoef, peqBcoef, &
454              &   tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)              tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
455    
456         fder_prev = fder             fder_prev = fder
457         fder = fder_prev + dflux_s + dflux_l         fder = fder_prev + dflux_s + dflux_l
458    
459         iloc = maxloc(fder(1:klon))         iloc = maxloc(fder(1:klon))
# Line 463  CONTAINS Line 461  CONTAINS
461            WRITE(*, *)'**** Debug fder****'            WRITE(*, *)'**** Debug fder****'
462            WRITE(*, *)'max fder(', iloc(1), ') = ', fder(iloc(1))            WRITE(*, *)'max fder(', iloc(1), ') = ', fder(iloc(1))
463            WRITE(*, *)'fder_prev, dflux_s, dflux_l', fder_prev(iloc(1)), &            WRITE(*, *)'fder_prev, dflux_s, dflux_l', fder_prev(iloc(1)), &
464                 &                        dflux_s(iloc(1)), dflux_l(iloc(1))                 dflux_s(iloc(1)), dflux_l(iloc(1))
465         endif         endif
466    
467         !IM: flux ocean-atmosphere utile pour le "slab" ocean         !IM: flux ocean-atmosphere utile pour le "slab" ocean
# Line 478  CONTAINS Line 476  CONTAINS
476         ! 2eme appel a interfoce pour le cumul des champs (en particulier         ! 2eme appel a interfoce pour le cumul des champs (en particulier
477         ! fluxsens et fluxlat calcules dans calcul_fluxs)         ! fluxsens et fluxlat calcules dans calcul_fluxs)
478    
479         if (ocean == 'slab  ') then         if (ocean == 'slab ') then
480            seaice=tmp_seaice            seaice=tmp_seaice
481            cumul = .true.            cumul = .true.
482            call interfoce_slab(klon, debut, itime, dtime, jour, &            call interfoce_slab(klon, debut, itime, dtime, jour, &
483                 & tmp_radsol, tmp_flux_o, tmp_flux_g, pctsrf, &                 tmp_radsol, tmp_flux_o, tmp_flux_g, pctsrf, &
484                 & tslab, seaice, pctsrf_new)                 tslab, seaice, pctsrf_new)
485    
486            tmp_pctsrf_slab=pctsrf_new            tmp_pctsrf_slab=pctsrf_new
487            DO i=1, knon            DO i=1, knon
# Line 494  CONTAINS Line 492  CONTAINS
492         ! calcul albedo         ! calcul albedo
493         if ( minval(rmu0) == maxval(rmu0) .and. minval(rmu0) == -999.999 ) then         if ( minval(rmu0) == maxval(rmu0) .and. minval(rmu0) == -999.999 ) then
494            CALL alboc(FLOAT(jour), rlat, alb_eau)            CALL alboc(FLOAT(jour), rlat, alb_eau)
495         else  ! cycle diurne         else ! cycle diurne
496            CALL alboc_cd(rmu0, alb_eau)            CALL alboc_cd(rmu0, alb_eau)
497         endif         endif
498         DO ii =1, knon         DO ii =1, knon
# Line 509  CONTAINS Line 507  CONTAINS
507         ! Surface "glace de mer" appel a l'interface avec l'ocean         ! Surface "glace de mer" appel a l'interface avec l'ocean
508    
509    
510         if (ocean == 'slab  ') then         if (ocean == 'slab ') then
511            pctsrf_new=tmp_pctsrf_slab            pctsrf_new=tmp_pctsrf_slab
512    
513            DO ii = 1, knon            DO ii = 1, knon
# Line 526  CONTAINS Line 524  CONTAINS
524            IF (soil_model) THEN            IF (soil_model) THEN
525               CALL soil(dtime, nisurf, knon, snow, tsurf_new, tsoil, soilcap, soilflux)               CALL soil(dtime, nisurf, knon, snow, tsurf_new, tsoil, soilcap, soilflux)
526               cal(1:knon) = RCPD / soilcap(1:knon)               cal(1:knon) = RCPD / soilcap(1:knon)
527               radsol(1:knon) = radsol(1:knon)  + soilflux(1:knon)               radsol(1:knon) = radsol(1:knon) + soilflux(1:knon)
528            ELSE            ELSE
529               dif_grnd = 1.0 / tau_gl               dif_grnd = 1.0 / tau_gl
530               cal = RCPD * calice               cal = RCPD * calice
# Line 536  CONTAINS Line 534  CONTAINS
534            beta = 1.0            beta = 1.0
535    
536         ELSE         ELSE
537            !                              ! lecture conditions limites            ! ! lecture conditions limites
538            CALL interfoce_lim(itime, dtime, jour, &            CALL interfoce_lim(itime, dtime, jour, &
539                 &  klon, nisurf, knon, knindex, &                 klon, nisurf, knon, knindex, &
540                 &  debut, &                 debut, &
541                 &  tsurf_new, pctsrf_new)                 tsurf_new, pctsrf_new)
542    
543            !IM cf LF            !IM cf LF
544            DO ii = 1, knon            DO ii = 1, knon
# Line 548  CONTAINS Line 546  CONTAINS
546               !IMbad IF (pctsrf_new(ii, nisurf) < EPSFRA) then               !IMbad IF (pctsrf_new(ii, nisurf) < EPSFRA) then
547               IF (pctsrf_new(knindex(ii), nisurf) < EPSFRA) then               IF (pctsrf_new(knindex(ii), nisurf) < EPSFRA) then
548                  snow(ii) = 0.0                  snow(ii) = 0.0
549                  !IM cf LF/JLD         tsurf(ii) = RTT - 1.8                  !IM cf LF/JLD tsurf(ii) = RTT - 1.8
550                  tsurf_new(ii) = RTT - 1.8                  tsurf_new(ii) = RTT - 1.8
551                  IF (soil_model) tsoil(ii, :) = RTT -1.8                  IF (soil_model) tsoil(ii, :) = RTT -1.8
552               endif               endif
# Line 557  CONTAINS Line 555  CONTAINS
555            CALL calbeta(dtime, nisurf, knon, snow, qsol, beta, capsol, dif_grnd)            CALL calbeta(dtime, nisurf, knon, snow, qsol, beta, capsol, dif_grnd)
556    
557            IF (soil_model) THEN            IF (soil_model) THEN
558               !IM cf LF/JLD        CALL soil(dtime, nisurf, knon, snow, tsurf, tsoil, soilcap, soilflux)               !IM cf LF/JLD CALL soil(dtime, nisurf, knon, snow, tsurf, tsoil, soilcap, soilflux)
559               CALL soil(dtime, nisurf, knon, snow, tsurf_new, tsoil, soilcap, soilflux)               CALL soil(dtime, nisurf, knon, snow, tsurf_new, tsoil, soilcap, soilflux)
560               cal(1:knon) = RCPD / soilcap(1:knon)               cal(1:knon) = RCPD / soilcap(1:knon)
561               radsol(1:knon) = radsol(1:knon)  + soilflux(1:knon)               radsol(1:knon) = radsol(1:knon) + soilflux(1:knon)
562               dif_grnd = 0.               dif_grnd = 0.
563            ELSE            ELSE
564               dif_grnd = 1.0 / tau_gl               dif_grnd = 1.0 / tau_gl
# Line 573  CONTAINS Line 571  CONTAINS
571         ENDIF         ENDIF
572    
573         CALL calcul_fluxs( klon, knon, nisurf, dtime, &         CALL calcul_fluxs( klon, knon, nisurf, dtime, &
574              &   tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &              tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
575              &   precip_rain, precip_snow, snow, qsurf,  &              precip_rain, precip_snow, snow, qsurf, &
576              &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &              radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
577              &   petAcoef, peqAcoef, petBcoef, peqBcoef, &              petAcoef, peqAcoef, petBcoef, peqBcoef, &
578              &   tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)              tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
579    
580         !IM: flux entre l'ocean et la glace de mer pour le "slab" ocean         !IM: flux entre l'ocean et la glace de mer pour le "slab" ocean
581         DO i = 1, knon         DO i = 1, knon
582            flux_g(i) = 0.0            flux_g(i) = 0.0
583    
584            !IM: faire dependre le coefficient de conduction de la glace de mer            !IM: faire dependre le coefficient de conduction de la glace de mer
585            !    de l'epaisseur de la glace de mer, dans l'hypothese ou le coeff.            ! de l'epaisseur de la glace de mer, dans l'hypothese ou le coeff.
586            !    actuel correspond a 3m de glace de mer, cf. L.Li            ! actuel correspond a 3m de glace de mer, cf. L.Li
587    
588            !      IF(1.EQ.0) THEN            ! IF(1.EQ.0) THEN
589            !       IF(siceh(i).GT.0.) THEN            ! IF(siceh(i).GT.0.) THEN
590            !        new_dif_grnd(i) = dif_grnd(i)*3./siceh(i)            ! new_dif_grnd(i) = dif_grnd(i)*3./siceh(i)
591            !       ELSE            ! ELSE
592            !        new_dif_grnd(i) = 0.            ! new_dif_grnd(i) = 0.
593            !       ENDIF            ! ENDIF
594            !      ENDIF !(1.EQ.0) THEN            ! ENDIF !(1.EQ.0) THEN
595    
596            IF (cal(i).GT.1.0e-15) flux_g(i)=(tsurf_new(i)-t_grnd) &            IF (cal(i).GT.1.0e-15) flux_g(i)=(tsurf_new(i)-t_grnd) &
597                 &                          * dif_grnd(i) *RCPD/cal(i)                 * dif_grnd(i) *RCPD/cal(i)
598            !   &                          * new_dif_grnd(i) *RCPD/cal(i)            ! & * new_dif_grnd(i) *RCPD/cal(i)
599            tmp_flux_g(knindex(i))=flux_g(i)            tmp_flux_g(knindex(i))=flux_g(i)
600            tmp_radsol(knindex(i))=radsol(i)            tmp_radsol(knindex(i))=radsol(i)
601         ENDDO         ENDDO
602    
603         CALL fonte_neige( klon, knon, nisurf, dtime, &         CALL fonte_neige( klon, knon, nisurf, dtime, &
604              &   tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &              tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
605              &   precip_rain, precip_snow, snow, qsol,  &              precip_rain, precip_snow, snow, qsol, &
606              &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &              radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
607              &   petAcoef, peqAcoef, petBcoef, peqBcoef, &              petAcoef, peqAcoef, petBcoef, peqBcoef, &
608              &   tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &              tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &
609              &   fqcalving, ffonte, run_off_lic_0)              fqcalving, ffonte, run_off_lic_0)
610    
611         !     calcul albedo         ! calcul albedo
612    
613         CALL albsno(klon, knon, dtime, agesno, alb_neig, precip_snow)           CALL albsno(klon, knon, dtime, agesno, alb_neig, precip_snow)
614         WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.         WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
615         zfra(1:knon) = MAX(0.0, MIN(1.0, snow(1:knon)/(snow(1:knon)+10.0)))         zfra(1:knon) = MAX(0.0, MIN(1.0, snow(1:knon)/(snow(1:knon)+10.0)))
616         alb_new(1 : knon) = alb_neig(1 : knon) *zfra(1:knon) + &         alb_new(1 : knon) = alb_neig(1 : knon) *zfra(1:knon) + &
617              0.6 * (1.0-zfra(1:knon))              0.6 * (1.0-zfra(1:knon))
618    
619         fder_prev = fder             fder_prev = fder
620         fder = fder_prev + dflux_s + dflux_l         fder = fder_prev + dflux_s + dflux_l
621    
622         iloc = maxloc(fder(1:klon))         iloc = maxloc(fder(1:klon))
# Line 626  CONTAINS Line 624  CONTAINS
624            WRITE(*, *)'**** Debug fder ****'            WRITE(*, *)'**** Debug fder ****'
625            WRITE(*, *)'max fder(', iloc(1), ') = ', fder(iloc(1))            WRITE(*, *)'max fder(', iloc(1), ') = ', fder(iloc(1))
626            WRITE(*, *)'fder_prev, dflux_s, dflux_l', fder_prev(iloc(1)), &            WRITE(*, *)'fder_prev, dflux_s, dflux_l', fder_prev(iloc(1)), &
627                 &                        dflux_s(iloc(1)), dflux_l(iloc(1))                 dflux_s(iloc(1)), dflux_l(iloc(1))
628         endif         endif
629    
630    
# Line 654  CONTAINS Line 652  CONTAINS
652         IF (soil_model) THEN         IF (soil_model) THEN
653            CALL soil(dtime, nisurf, knon, snow, tsurf, tsoil, soilcap, soilflux)            CALL soil(dtime, nisurf, knon, snow, tsurf, tsoil, soilcap, soilflux)
654            cal(1:knon) = RCPD / soilcap(1:knon)            cal(1:knon) = RCPD / soilcap(1:knon)
655            radsol(1:knon)  = radsol(1:knon) + soilflux(1:knon)            radsol(1:knon) = radsol(1:knon) + soilflux(1:knon)
656         ELSE         ELSE
657            cal = RCPD * calice            cal = RCPD * calice
658            WHERE (snow > 0.0) cal = RCPD * calsno            WHERE (snow > 0.0) cal = RCPD * calsno
# Line 663  CONTAINS Line 661  CONTAINS
661         dif_grnd = 0.0         dif_grnd = 0.0
662    
663         call calcul_fluxs( klon, knon, nisurf, dtime, &         call calcul_fluxs( klon, knon, nisurf, dtime, &
664              &   tsurf, p1lay, cal, beta, tq_cdrag, ps, &              tsurf, p1lay, cal, beta, tq_cdrag, ps, &
665              &   precip_rain, precip_snow, snow, qsurf,  &              precip_rain, precip_snow, snow, qsurf, &
666              &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &              radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
667              &   petAcoef, peqAcoef, petBcoef, peqBcoef, &              petAcoef, peqAcoef, petBcoef, peqBcoef, &
668              &   tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)              tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
669    
670         call fonte_neige( klon, knon, nisurf, dtime, &         call fonte_neige( klon, knon, nisurf, dtime, &
671              &   tsurf, p1lay, cal, beta, tq_cdrag, ps, &              tsurf, p1lay, cal, beta, tq_cdrag, ps, &
672              &   precip_rain, precip_snow, snow, qsol,  &              precip_rain, precip_snow, snow, qsol, &
673              &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &              radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
674              &   petAcoef, peqAcoef, petBcoef, peqBcoef, &              petAcoef, peqAcoef, petBcoef, peqBcoef, &
675              &   tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &              tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &
676              &   fqcalving, ffonte, run_off_lic_0)              fqcalving, ffonte, run_off_lic_0)
677    
678         ! passage du run-off des glaciers calcule dans fonte_neige au coupleur         ! passage du run-off des glaciers calcule dans fonte_neige au coupleur
679         bidule=0.         bidule=0.
680         bidule(1:knon)= run_off_lic(1:knon)             bidule(1:knon)= run_off_lic(1:knon)
681         call gath2cpl(bidule, tmp_rlic, klon, knon, iim, jjm, knindex)         call gath2cpl(bidule, tmp_rlic, klon, knon, iim, jjm, knindex)
682    
683         ! calcul albedo         ! calcul albedo
684    
685         CALL albsno(klon, knon, dtime, agesno, alb_neig, precip_snow)           CALL albsno(klon, knon, dtime, agesno, alb_neig, precip_snow)
686         WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.         WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
687         zfra(1:knon) = MAX(0.0, MIN(1.0, snow(1:knon)/(snow(1:knon)+10.0)))         zfra(1:knon) = MAX(0.0, MIN(1.0, snow(1:knon)/(snow(1:knon)+10.0)))
688         alb_new(1 : knon)  = alb_neig(1 : knon)*zfra(1:knon) + &         alb_new(1 : knon) = alb_neig(1 : knon)*zfra(1:knon) + &
689              &                     0.6 * (1.0-zfra(1:knon))              0.6 * (1.0-zfra(1:knon))
690    
691         !IM: plusieurs choix/tests sur l'albedo des "glaciers continentaux"         !IM: plusieurs choix/tests sur l'albedo des "glaciers continentaux"
692         !       alb_new(1 : knon)  = 0.6 !IM cf FH/GK         ! alb_new(1 : knon) = 0.6 !IM cf FH/GK
693         !       alb_new(1 : knon)  = 0.82         ! alb_new(1 : knon) = 0.82
694         !       alb_new(1 : knon)  = 0.77 !211003 Ksta0.77         ! alb_new(1 : knon) = 0.77 !211003 Ksta0.77
695         !       alb_new(1 : knon)  = 0.8 !KstaTER0.8 & LMD_ARMIP5         ! alb_new(1 : knon) = 0.8 !KstaTER0.8 & LMD_ARMIP5
696         !IM: KstaTER0.77 & LMD_ARMIP6             !IM: KstaTER0.77 & LMD_ARMIP6
697         alb_new(1 : knon)  = 0.77         alb_new(1 : knon) = 0.77
698    
699    
700         ! Rugosite         ! Rugosite
# Line 718  CONTAINS Line 716  CONTAINS
716    
717    !************************    !************************
718    
719    SUBROUTINE interfoce_slab(klon, debut, itap, dtime, ijour, &    SUBROUTINE interfoce_slab(klon, debut, itap, dtime, ijour,  &
720         & radsol, fluxo, fluxg, pctsrf, &         radsol, fluxo, fluxg, pctsrf,  &
721         & tslab, seaice, pctsrf_slab)         tslab, seaice, pctsrf_slab)
722    
723      ! Cette routine calcule la temperature d'un slab ocean, la glace de mer      ! Cette routine calcule la temperature d'un slab ocean, la glace de mer
724      ! et les pourcentages de la maille couverte par l'ocean libre et/ou      ! et les pourcentages de la maille couverte par l'ocean libre et/ou
# Line 729  CONTAINS Line 727  CONTAINS
727      ! I. Musat 04.02.2005      ! I. Musat 04.02.2005
728    
729      ! input:      ! input:
730      !   klon         nombre total de points de grille      ! klon nombre total de points de grille
731      !   debut        logical: 1er appel a la physique      ! debut logical: 1er appel a la physique
732      !   itap         numero du pas de temps      ! itap numero du pas de temps
733      !   dtime        pas de temps de la physique (en s)      ! dtime pas de temps de la physique (en s)
734      !   ijour        jour dans l'annee en cours      ! ijour jour dans l'annee en cours
735      !   radsol       rayonnement net au sol (LW + SW)      ! radsol rayonnement net au sol (LW + SW)
736      !   fluxo        flux turbulent (sensible + latent) sur les mailles oceaniques      ! fluxo flux turbulent (sensible + latent) sur les mailles oceaniques
737      !   fluxg        flux de conduction entre la surface de la glace de mer et l'ocean      ! fluxg flux de conduction entre la surface de la glace de mer et l'ocean
738      !   pctsrf       tableau des pourcentages de surface de chaque maille      ! pctsrf tableau des pourcentages de surface de chaque maille
739      ! output:      ! output:
740      !   tslab        temperature de l'ocean libre      ! tslab temperature de l'ocean libre
741      !   seaice       glace de mer (kg/m2)      ! seaice glace de mer (kg/m2)
742      !   pctsrf_slab  "pourcentages" (valeurs entre 0. et 1.) surfaces issus du slab      ! pctsrf_slab "pourcentages" (valeurs entre 0. et 1.) surfaces issus du slab
743    
744      use indicesol      use indicesol
745      use clesphys      use clesphys
746      use abort_gcm_m, only: abort_gcm      use abort_gcm_m, only: abort_gcm
747      use YOMCST      use SUPHEC_M
748    
749      ! Parametres d'entree      ! Parametres d'entree
750      integer, intent(IN) :: klon      integer, intent(IN) :: klon
# Line 760  CONTAINS Line 758  CONTAINS
758      real, dimension(klon, nbsrf), intent(IN) :: pctsrf      real, dimension(klon, nbsrf), intent(IN) :: pctsrf
759      ! Parametres de sortie      ! Parametres de sortie
760      real, dimension(klon), intent(INOUT) :: tslab      real, dimension(klon), intent(INOUT) :: tslab
761      real, dimension(klon), intent(INOUT)        :: seaice ! glace de mer (kg/m2)      real, dimension(klon), intent(INOUT) :: seaice ! glace de mer (kg/m2)
762      real, dimension(klon, nbsrf), intent(OUT) :: pctsrf_slab      real, dimension(klon, nbsrf), intent(OUT) :: pctsrf_slab
763    
764      ! Variables locales :      ! Variables locales :
# Line 769  CONTAINS Line 767  CONTAINS
767      real, allocatable, dimension(:), save :: tmp_tslab, tmp_seaice      real, allocatable, dimension(:), save :: tmp_tslab, tmp_seaice
768      REAL, allocatable, dimension(:), save :: slab_bils      REAL, allocatable, dimension(:), save :: slab_bils
769      REAL, allocatable, dimension(:), save :: lmt_bils      REAL, allocatable, dimension(:), save :: lmt_bils
770      logical, save              :: check = .false.      logical, save :: check = .false.
771    
772      REAL, parameter :: cyang=50.0 * 4.228e+06 ! capacite calorifique volumetrique de l'eau J/(m2 K)      REAL, parameter :: cyang=50.0 * 4.228e+06 ! capacite calorifique volumetrique de l'eau J/(m2 K)
773      REAL, parameter :: cbing=0.334e+05        ! J/kg      REAL, parameter :: cbing=0.334e+05 ! J/kg
774      real, dimension(klon)                 :: siceh !hauteur de la glace de mer (m)      real, dimension(klon) :: siceh !hauteur de la glace de mer (m)
775      INTEGER :: i      INTEGER :: i
776      integer :: sum_error, error      integer :: sum_error, error
777      REAL :: zz, za, zb      REAL :: zz, za, zb
# Line 798  CONTAINS Line 796  CONTAINS
796    
797         IF (check) THEN         IF (check) THEN
798            PRINT*, 'interfoce_slab klon, debut, itap, dtime, ijour, &            PRINT*, 'interfoce_slab klon, debut, itap, dtime, ijour, &
799                 &          lmt_pas ', klon, debut, itap, dtime, ijour, &                 & lmt_pas ', klon, debut, itap, dtime, ijour, &
800                 &          lmt_pas                 lmt_pas
801         ENDIF !check         ENDIF !check
802    
803         PRINT*, '************************'         PRINT*, '************************'
# Line 824  CONTAINS Line 822  CONTAINS
822    
823      DO i = 1, klon      DO i = 1, klon
824         IF((pctsrf_slab(i, is_oce).GT.epsfra).OR. &         IF((pctsrf_slab(i, is_oce).GT.epsfra).OR. &
825              &  (pctsrf_slab(i, is_sic).GT.epsfra)) THEN              (pctsrf_slab(i, is_sic).GT.epsfra)) THEN
826    
827            ! fabriquer de la glace si congelation atteinte:            ! fabriquer de la glace si congelation atteinte:
828    
829            IF (tmp_tslab(i).LT.(RTT-1.8)) THEN            IF (tmp_tslab(i).LT.(RTT-1.8)) THEN
830               zz =  (RTT-1.8)-tmp_tslab(i)               zz = (RTT-1.8)-tmp_tslab(i)
831               tmp_seaice(i) = tmp_seaice(i) + cyang/cbing * zz               tmp_seaice(i) = tmp_seaice(i) + cyang/cbing * zz
832               seaice(i) = tmp_seaice(i)               seaice(i) = tmp_seaice(i)
833               tmp_tslab(i) = RTT-1.8               tmp_tslab(i) = RTT-1.8
# Line 863  CONTAINS Line 861  CONTAINS
861            ! et pctsrf(i, is_sic) croit lineairement avec seaice de 0. a 20cm d'epaisseur            ! et pctsrf(i, is_sic) croit lineairement avec seaice de 0. a 20cm d'epaisseur
862    
863            pctsrf_slab(i, is_sic)=MIN(siceh(i)/0.20, &            pctsrf_slab(i, is_sic)=MIN(siceh(i)/0.20, &
864                 &                      1.-(pctsrf_slab(i, is_ter)+pctsrf_slab(i, is_lic)))                 1.-(pctsrf_slab(i, is_ter)+pctsrf_slab(i, is_lic)))
865            pctsrf_slab(i, is_oce)=1.0 - &            pctsrf_slab(i, is_oce)=1.0 - &
866                 &      (pctsrf_slab(i, is_ter)+pctsrf_slab(i, is_lic)+pctsrf_slab(i, is_sic))                 (pctsrf_slab(i, is_ter)+pctsrf_slab(i, is_lic)+pctsrf_slab(i, is_sic))
867         ENDIF !pctsrf         ENDIF !pctsrf
868      ENDDO      ENDDO
869    
# Line 875  CONTAINS Line 873  CONTAINS
873         za = radsol(i) + fluxo(i)         za = radsol(i) + fluxo(i)
874         zb = fluxg(i)         zb = fluxg(i)
875         IF((pctsrf_slab(i, is_oce).GT.epsfra).OR. &         IF((pctsrf_slab(i, is_oce).GT.epsfra).OR. &
876              &   (pctsrf_slab(i, is_sic).GT.epsfra)) THEN              (pctsrf_slab(i, is_sic).GT.epsfra)) THEN
877            slab_bils(i)=slab_bils(i)+(za*pctsrf_slab(i, is_oce) &            slab_bils(i)=slab_bils(i)+(za*pctsrf_slab(i, is_oce) &
878                 &             +zb*pctsrf_slab(i, is_sic))/ FLOAT(lmt_pas)                 +zb*pctsrf_slab(i, is_sic))/ FLOAT(lmt_pas)
879         ENDIF         ENDIF
880      ENDDO !klon      ENDDO !klon
881    
# Line 886  CONTAINS Line 884  CONTAINS
884      IF (MOD(itap, lmt_pas).EQ.0) THEN !fin de journee      IF (MOD(itap, lmt_pas).EQ.0) THEN !fin de journee
885         DO i = 1, klon         DO i = 1, klon
886            IF ((pctsrf_slab(i, is_oce).GT.epsfra).OR. &            IF ((pctsrf_slab(i, is_oce).GT.epsfra).OR. &
887                 &    (pctsrf_slab(i, is_sic).GT.epsfra)) THEN                 (pctsrf_slab(i, is_sic).GT.epsfra)) THEN
888               tmp_tslab(i) = tmp_tslab(i) + &               tmp_tslab(i) = tmp_tslab(i) + &
889                    & (slab_bils(i)-lmt_bils(i)) &                    (slab_bils(i)-lmt_bils(i)) &
890                    &                         /cyang*unjour                    /cyang*unjour
891               ! on remet l'accumulation a 0               ! on remet l'accumulation a 0
892               slab_bils(i) = 0.               slab_bils(i) = 0.
893            ENDIF !pctsrf            ENDIF !pctsrf
# Line 902  CONTAINS Line 900  CONTAINS
900    
901    !************************    !************************
902    
903    SUBROUTINE interfoce_lim(itime, dtime, jour, &    SUBROUTINE interfoce_lim(itime, dtime, jour,  &
904         & klon, nisurf, knon, knindex, &         klon, nisurf, knon, knindex,  &
905         & debut,  &         debut,  &
906         & lmt_sst, pctsrf_new)         lmt_sst, pctsrf_new)
907    
908      ! Cette routine sert d'interface entre le modele atmospherique et      ! Cette routine sert d'interface entre le modele atmospherique et
909      ! un fichier de conditions aux limites      ! un fichier de conditions aux limites
# Line 916  CONTAINS Line 914  CONTAINS
914      use indicesol      use indicesol
915    
916      integer, intent(IN) :: itime ! numero du pas de temps courant      integer, intent(IN) :: itime ! numero du pas de temps courant
917      real   , intent(IN) :: dtime ! pas de temps de la physique (en s)      real , intent(IN) :: dtime ! pas de temps de la physique (en s)
918      integer, intent(IN) :: jour ! jour a lire dans l'annee      integer, intent(IN) :: jour ! jour a lire dans l'annee
919      integer, intent(IN) :: nisurf ! index de la surface a traiter (1 = sol continental)      integer, intent(IN) :: nisurf ! index de la surface a traiter (1 = sol continental)
920      integer, intent(IN) :: knon ! nombre de points dans le domaine a traiter      integer, intent(IN) :: knon ! nombre de points dans le domaine a traiter
# Line 926  CONTAINS Line 924  CONTAINS
924    
925      ! Parametres de sortie      ! Parametres de sortie
926      ! output:      ! output:
927      !   lmt_sst      SST lues dans le fichier de CL      ! lmt_sst SST lues dans le fichier de CL
928      !   pctsrf_new   sous-maille fractionnelle      ! pctsrf_new sous-maille fractionnelle
929      real, intent(out), dimension(klon) :: lmt_sst      real, intent(out), dimension(klon) :: lmt_sst
930      real, intent(out), dimension(klon, nbsrf) :: pctsrf_new      real, intent(out), dimension(klon, nbsrf) :: pctsrf_new
931    
932      ! Variables locales      ! Variables locales
933      integer     :: ii      integer :: ii
934      INTEGER, save :: lmt_pas     ! frequence de lecture des conditions limites      INTEGER, save :: lmt_pas ! frequence de lecture des conditions limites
935      ! (en pas de physique)      ! (en pas de physique)
936      logical, save :: deja_lu    ! pour indiquer que le jour a lire a deja      logical, save :: deja_lu ! pour indiquer que le jour a lire a deja
937      ! lu pour une surface precedente      ! lu pour une surface precedente
938      integer, save :: jour_lu      integer, save :: jour_lu
939      integer      :: ierr      integer :: ierr
940      character (len = 20) :: modname = 'interfoce_lim'      character (len = 20) :: modname = 'interfoce_lim'
941      character (len = 80) :: abort_message      character (len = 80) :: abort_message
942      logical, save     :: newlmt = .TRUE.      logical, save :: newlmt = .TRUE.
943      logical, save     :: check = .FALSE.      logical, save :: check = .FALSE.
944      ! Champs lus dans le fichier de CL      ! Champs lus dans le fichier de CL
945      real, allocatable , save, dimension(:) :: sst_lu, rug_lu, nat_lu      real, allocatable , save, dimension(:) :: sst_lu, rug_lu, nat_lu
946      real, allocatable , save, dimension(:, :) :: pct_tmp      real, allocatable , save, dimension(:, :) :: pct_tmp
# Line 950  CONTAINS Line 948  CONTAINS
948      ! quelques variables pour netcdf      ! quelques variables pour netcdf
949    
950      include "netcdf.inc"      include "netcdf.inc"
951      integer              :: nid, nvarid      integer :: nid, nvarid
952      integer, dimension(2) :: start, epais      integer, dimension(2) :: start, epais
953    
954      ! --------------------------------------------------      ! --------------------------------------------------
# Line 1042  CONTAINS Line 1040  CONTAINS
1040               call abort_gcm(modname, abort_message, 1)               call abort_gcm(modname, abort_message, 1)
1041            endif            endif
1042    
1043         else  ! on en est toujours a rnatur         else ! on en est toujours a rnatur
1044    
1045            ierr = NF_INQ_VARID(nid, 'NAT', nvarid)            ierr = NF_INQ_VARID(nid, 'NAT', nvarid)
1046            if (ierr /= NF_NOERR) then            if (ierr /= NF_NOERR) then
# Line 1064  CONTAINS Line 1062  CONTAINS
1062            enddo            enddo
1063    
1064    
1065            !  On se retrouve avec ocean en 1 et terre en 2 alors qu'on veut le contraire            ! On se retrouve avec ocean en 1 et terre en 2 alors qu'on veut le contraire
1066    
1067            pctsrf_new = pct_tmp            pctsrf_new = pct_tmp
1068            pctsrf_new (:, 2)= pct_tmp (:, 1)            pctsrf_new (:, 2)= pct_tmp (:, 1)
# Line 1107  CONTAINS Line 1105  CONTAINS
1105    
1106    !************************    !************************
1107    
1108    SUBROUTINE interfsur_lim(itime, dtime, jour, &    SUBROUTINE interfsur_lim(itime, dtime, jour,  &
1109         & klon, nisurf, knon, knindex, &         klon, nisurf, knon, knindex,  &
1110         & debut,  &         debut,  &
1111         & lmt_alb, lmt_rug)         lmt_alb, lmt_rug)
1112    
1113      ! Cette routine sert d'interface entre le modèle atmosphérique et      ! Cette routine sert d'interface entre le modèle atmosphérique et
1114      ! un fichier de conditions aux limites.      ! un fichier de conditions aux limites.
# Line 1121  CONTAINS Line 1119  CONTAINS
1119    
1120      ! Parametres d'entree      ! Parametres d'entree
1121      ! input:      ! input:
1122      !   itime        numero du pas de temps courant      ! itime numero du pas de temps courant
1123      !   dtime        pas de temps de la physique (en s)      ! dtime pas de temps de la physique (en s)
1124      !   jour         jour a lire dans l'annee      ! jour jour a lire dans l'annee
1125      !   nisurf       index de la surface a traiter (1 = sol continental)      ! nisurf index de la surface a traiter (1 = sol continental)
1126      !   knon         nombre de points dans le domaine a traiter      ! knon nombre de points dans le domaine a traiter
1127      !   knindex      index des points de la surface a traiter      ! knindex index des points de la surface a traiter
1128      !   klon         taille de la grille      ! klon taille de la grille
1129      !   debut        logical: 1er appel a la physique (initialisation)      ! debut logical: 1er appel a la physique (initialisation)
1130      integer, intent(IN) :: itime      integer, intent(IN) :: itime
1131      real   , intent(IN) :: dtime      real , intent(IN) :: dtime
1132      integer, intent(IN) :: jour      integer, intent(IN) :: jour
1133      integer, intent(IN) :: nisurf      integer, intent(IN) :: nisurf
1134      integer, intent(IN) :: knon      integer, intent(IN) :: knon
# Line 1140  CONTAINS Line 1138  CONTAINS
1138    
1139      ! Parametres de sortie      ! Parametres de sortie
1140      ! output:      ! output:
1141      !   lmt_sst      SST lues dans le fichier de CL      ! lmt_sst SST lues dans le fichier de CL
1142      !   lmt_alb      Albedo lu      ! lmt_alb Albedo lu
1143      !   lmt_rug      longueur de rugosité lue      ! lmt_rug longueur de rugosité lue
1144      !   pctsrf_new   sous-maille fractionnelle      ! pctsrf_new sous-maille fractionnelle
1145      real, intent(out), dimension(klon) :: lmt_alb      real, intent(out), dimension(klon) :: lmt_alb
1146      real, intent(out), dimension(klon) :: lmt_rug      real, intent(out), dimension(klon) :: lmt_rug
1147    
1148      ! Variables locales      ! Variables locales
1149      integer     :: ii      integer :: ii
1150      integer, save :: lmt_pas     ! frequence de lecture des conditions limites      integer, save :: lmt_pas ! frequence de lecture des conditions limites
1151      ! (en pas de physique)      ! (en pas de physique)
1152      logical, save :: deja_lu_sur! pour indiquer que le jour a lire a deja      logical, save :: deja_lu_sur! pour indiquer que le jour a lire a deja
1153      ! lu pour une surface precedente      ! lu pour une surface precedente
1154      integer, save :: jour_lu_sur      integer, save :: jour_lu_sur
1155      integer      :: ierr      integer :: ierr
1156      character (len = 20) :: modname = 'interfsur_lim'      character (len = 20) :: modname = 'interfsur_lim'
1157      character (len = 80) :: abort_message      character (len = 80) :: abort_message
1158      logical, save     :: newlmt = .false.      logical, save :: newlmt = .false.
1159      logical, save     :: check = .false.      logical, save :: check = .false.
1160      ! Champs lus dans le fichier de CL      ! Champs lus dans le fichier de CL
1161      real, allocatable , save, dimension(:) :: alb_lu, rug_lu      real, allocatable , save, dimension(:) :: alb_lu, rug_lu
1162    
1163      ! quelques variables pour netcdf      ! quelques variables pour netcdf
1164    
1165      include "netcdf.inc"      include "netcdf.inc"
1166      integer , save             :: nid, nvarid      integer , save :: nid, nvarid
1167      integer, dimension(2), save :: start, epais      integer, dimension(2), save :: start, epais
1168    
1169      !------------------------------------------------------------      !------------------------------------------------------------
# Line 1185  CONTAINS Line 1183  CONTAINS
1183    
1184      ! Tester d'abord si c'est le moment de lire le fichier      ! Tester d'abord si c'est le moment de lire le fichier
1185      if (mod(itime-1, lmt_pas) == 0 .and. .not. deja_lu_sur) then      if (mod(itime-1, lmt_pas) == 0 .and. .not. deja_lu_sur) then
1186      
1187         ! Ouverture du fichier         ! Ouverture du fichier
1188      
1189         ierr = NF_OPEN ('limit.nc', NF_NOWRITE, nid)         ierr = NF_OPEN ('limit.nc', NF_NOWRITE, nid)
1190         if (ierr.NE.NF_NOERR) then         if (ierr.NE.NF_NOERR) then
1191            abort_message &            abort_message &
# Line 1201  CONTAINS Line 1199  CONTAINS
1199         start(2) = jour         start(2) = jour
1200         epais(1) = klon         epais(1) = klon
1201         epais(2) = 1         epais(2) = 1
1202      
1203         ! Lecture Albedo         ! Lecture Albedo
1204      
1205         ierr = NF_INQ_VARID(nid, 'ALB', nvarid)         ierr = NF_INQ_VARID(nid, 'ALB', nvarid)
1206         if (ierr /= NF_NOERR) then         if (ierr /= NF_NOERR) then
1207            abort_message = 'Le champ <ALB> est absent'            abort_message = 'Le champ <ALB> est absent'
# Line 1214  CONTAINS Line 1212  CONTAINS
1212            abort_message = 'Lecture echouee pour <ALB>'            abort_message = 'Lecture echouee pour <ALB>'
1213            call abort_gcm(modname, abort_message, 1)            call abort_gcm(modname, abort_message, 1)
1214         endif         endif
1215      
1216         ! Lecture rugosité         ! Lecture rugosité
1217      
1218         ierr = NF_INQ_VARID(nid, 'RUG', nvarid)         ierr = NF_INQ_VARID(nid, 'RUG', nvarid)
1219         if (ierr /= NF_NOERR) then         if (ierr /= NF_NOERR) then
1220            abort_message = 'Le champ <RUG> est absent'            abort_message = 'Le champ <RUG> est absent'
# Line 1228  CONTAINS Line 1226  CONTAINS
1226            call abort_gcm(modname, abort_message, 1)            call abort_gcm(modname, abort_message, 1)
1227         endif         endif
1228    
1229      
1230         ! Fin de lecture         ! Fin de lecture
1231      
1232         ierr = NF_CLOSE(nid)         ierr = NF_CLOSE(nid)
1233         deja_lu_sur = .true.         deja_lu_sur = .true.
1234         jour_lu_sur = jour         jour_lu_sur = jour
# Line 1238  CONTAINS Line 1236  CONTAINS
1236    
1237      ! Recopie des variables dans les champs de sortie      ! Recopie des variables dans les champs de sortie
1238    
1239  !!$  lmt_alb = 0.0  !!$ lmt_alb = 0.0
1240  !!$  lmt_rug = 0.0  !!$ lmt_rug = 0.0
1241      lmt_alb = 999999.      lmt_alb = 999999.
1242      lmt_rug = 999999.      lmt_rug = 999999.
1243      DO ii = 1, knon      DO ii = 1, knon
# Line 1251  CONTAINS Line 1249  CONTAINS
1249    
1250    !************************    !************************
1251    
1252    SUBROUTINE calcul_fluxs( klon, knon, nisurf, dtime, &    SUBROUTINE calcul_fluxs( klon, knon, nisurf, dtime,  &
1253         & tsurf, p1lay, cal, beta, coef1lay, ps, &         tsurf, p1lay, cal, beta, coef1lay, ps,  &
1254         & precip_rain, precip_snow, snow, qsurf, &         precip_rain, precip_snow, snow, qsurf,  &
1255         & radsol, dif_grnd, t1lay, q1lay, u1lay, v1lay, &         radsol, dif_grnd, t1lay, q1lay, u1lay, v1lay,  &
1256         & petAcoef, peqAcoef, petBcoef, peqBcoef, &         petAcoef, peqAcoef, petBcoef, peqBcoef,  &
1257         & tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
1258    
1259      ! Cette routine calcule les fluxs en h et q a l'interface et eventuellement      ! Cette routine calcule les fluxs en h et q a l'interface et eventuellement
1260      ! une temperature de surface (au cas ou ok_veget = false)      ! une temperature de surface (au cas ou ok_veget = false)
# Line 1264  CONTAINS Line 1262  CONTAINS
1262      ! L. Fairhead 4/2000      ! L. Fairhead 4/2000
1263    
1264      ! input:      ! input:
1265      !   knon         nombre de points a traiter      ! knon nombre de points a traiter
1266      !   nisurf       surface a traiter      ! nisurf surface a traiter
1267      !   tsurf        temperature de surface      ! tsurf temperature de surface
1268      !   p1lay        pression 1er niveau (milieu de couche)      ! p1lay pression 1er niveau (milieu de couche)
1269      !   cal          capacite calorifique du sol      ! cal capacite calorifique du sol
1270      !   beta         evap reelle      ! beta evap reelle
1271      !   coef1lay     coefficient d'echange      ! coef1lay coefficient d'echange
1272      !   ps           pression au sol      ! ps pression au sol
1273      !   precip_rain  precipitations liquides      ! precip_rain precipitations liquides
1274      !   precip_snow  precipitations solides      ! precip_snow precipitations solides
1275      !   snow         champs hauteur de neige      ! snow champs hauteur de neige
1276      !   runoff       runoff en cas de trop plein      ! runoff runoff en cas de trop plein
1277      !   petAcoef     coeff. A de la resolution de la CL pour t      ! petAcoef coeff. A de la resolution de la CL pour t
1278      !   peqAcoef     coeff. A de la resolution de la CL pour q      ! peqAcoef coeff. A de la resolution de la CL pour q
1279      !   petBcoef     coeff. B de la resolution de la CL pour t      ! petBcoef coeff. B de la resolution de la CL pour t
1280      !   peqBcoef     coeff. B de la resolution de la CL pour q      ! peqBcoef coeff. B de la resolution de la CL pour q
1281      !   radsol       rayonnement net aus sol (LW + SW)      ! radsol rayonnement net aus sol (LW + SW)
1282      !   dif_grnd     coeff. diffusion vers le sol profond      ! dif_grnd coeff. diffusion vers le sol profond
1283    
1284      ! output:      ! output:
1285      !   tsurf_new    temperature au sol      ! tsurf_new temperature au sol
1286      !   qsurf        humidite de l'air au dessus du sol      ! qsurf humidite de l'air au dessus du sol
1287      !   fluxsens     flux de chaleur sensible      ! fluxsens flux de chaleur sensible
1288      !   fluxlat      flux de chaleur latente      ! fluxlat flux de chaleur latente
1289      !   dflux_s      derivee du flux de chaleur sensible / Ts      ! dflux_s derivee du flux de chaleur sensible / Ts
1290      !   dflux_l      derivee du flux de chaleur latente  / Ts      ! dflux_l derivee du flux de chaleur latente / Ts
1291    
1292    
1293      use indicesol      use indicesol
1294      use abort_gcm_m, only: abort_gcm      use abort_gcm_m, only: abort_gcm
1295      use yoethf      use yoethf_m
1296      use fcttre, only: thermcep, foeew, qsats, qsatl, foede, dqsats, dqsatl      use fcttre, only: thermcep, foeew, qsats, qsatl, foede, dqsats, dqsatl
1297      use YOMCST      use SUPHEC_M
1298    
1299      ! Parametres d'entree      ! Parametres d'entree
1300      integer, intent(IN) :: knon, nisurf, klon      integer, intent(IN) :: knon, nisurf, klon
1301      real   , intent(IN) :: dtime      real , intent(IN) :: dtime
1302      real, dimension(klon), intent(IN) :: petAcoef, peqAcoef      real, dimension(klon), intent(IN) :: petAcoef, peqAcoef
1303      real, dimension(klon), intent(IN) :: petBcoef, peqBcoef      real, dimension(klon), intent(IN) :: petBcoef, peqBcoef
1304      real, dimension(klon), intent(IN) :: ps, q1lay      real, dimension(klon), intent(IN) :: ps, q1lay
# Line 1321  CONTAINS Line 1319  CONTAINS
1319      real, dimension(klon) :: zx_pkh, zx_dq_s_dt, zx_qsat, zx_coef      real, dimension(klon) :: zx_pkh, zx_dq_s_dt, zx_qsat, zx_coef
1320      real, dimension(klon) :: zx_sl, zx_k1      real, dimension(klon) :: zx_sl, zx_k1
1321      real, dimension(klon) :: zx_q_0 , d_ts      real, dimension(klon) :: zx_q_0 , d_ts
1322      real                  :: zdelta, zcvm5, zx_qs, zcor, zx_dq_s_dh      real :: zdelta, zcvm5, zx_qs, zcor, zx_dq_s_dh
1323      real                  :: bilan_f, fq_fonte      real :: bilan_f, fq_fonte
1324      REAL                  :: subli, fsno      REAL :: subli, fsno
1325      REAL                  :: qsat_new, q1_new      REAL :: qsat_new, q1_new
1326      real, parameter :: t_grnd = 271.35, t_coup = 273.15      real, parameter :: t_grnd = 271.35, t_coup = 273.15
1327      !! PB temporaire en attendant mieux pour le modele de neige      !! PB temporaire en attendant mieux pour le modele de neige
1328      REAL, parameter :: chasno = 3.334E+05/(2.3867E+06*0.15)      REAL, parameter :: chasno = 3.334E+05/(2.3867E+06*0.15)
1329    
1330      logical, save         :: check = .false.      logical, save :: check = .false.
1331      character (len = 20)  :: modname = 'calcul_fluxs'      character (len = 20) :: modname = 'calcul_fluxs'
1332      logical, save         :: fonte_neige = .false.      logical, save :: fonte_neige = .false.
1333      real, save            :: max_eau_sol = 150.0      real, save :: max_eau_sol = 150.0
1334      character (len = 80) :: abort_message      character (len = 80) :: abort_message
1335      logical, save         :: first = .true., second=.false.      logical, save :: first = .true., second=.false.
1336    
1337      if (check) write(*, *)'Entree ', modname, ' surface = ', nisurf      if (check) write(*, *)'Entree ', modname, ' surface = ', nisurf
1338    
1339      IF (check) THEN      IF (check) THEN
1340         WRITE(*, *)' radsol (min, max)' &         WRITE(*, *)' radsol (min, max)' &
1341              &     , MINVAL(radsol(1:knon)), MAXVAL(radsol(1:knon))              , MINVAL(radsol(1:knon)), MAXVAL(radsol(1:knon))
1342         !!CALL flush(6)         !!CALL flush(6)
1343      ENDIF      ENDIF
1344    
# Line 1353  CONTAINS Line 1351  CONTAINS
1351    
1352      ! Traitement neige et humidite du sol      ! Traitement neige et humidite du sol
1353    
1354  !!$  WRITE(*, *)'test calcul_flux, surface ', nisurf  !!$ WRITE(*, *)'test calcul_flux, surface ', nisurf
1355      !!PB test      !!PB test
1356  !!$    if (nisurf == is_oce) then  !!$ if (nisurf == is_oce) then
1357  !!$      snow = 0.  !!$ snow = 0.
1358  !!$      qsol = max_eau_sol  !!$ qsol = max_eau_sol
1359  !!$    else  !!$ else
1360  !!$      where (precip_snow > 0.) snow = snow + (precip_snow * dtime)  !!$ where (precip_snow > 0.) snow = snow + (precip_snow * dtime)
1361  !!$      where (snow > epsilon(snow)) snow = max(0.0, snow - (evap * dtime))  !!$ where (snow > epsilon(snow)) snow = max(0.0, snow - (evap * dtime))
1362  !!$!      snow = max(0.0, snow + (precip_snow - evap) * dtime)  !!$! snow = max(0.0, snow + (precip_snow - evap) * dtime)
1363  !!$      where (precip_rain > 0.) qsol = qsol + (precip_rain - evap) * dtime  !!$ where (precip_rain > 0.) qsol = qsol + (precip_rain - evap) * dtime
1364  !!$    endif  !!$ endif
1365  !!$    IF (nisurf /= is_ter) qsol = max_eau_sol  !!$ IF (nisurf /= is_ter) qsol = max_eau_sol
1366    
1367    
1368      ! Initialisation      ! Initialisation
# Line 1388  CONTAINS Line 1386  CONTAINS
1386            zcor=1./(1.-retv*zx_qs)            zcor=1./(1.-retv*zx_qs)
1387            zx_qs=zx_qs*zcor            zx_qs=zx_qs*zcor
1388            zx_dq_s_dh = FOEDE(tsurf(i), zdelta, zcvm5, zx_qs, zcor) &            zx_dq_s_dh = FOEDE(tsurf(i), zdelta, zcvm5, zx_qs, zcor) &
1389                 &                 /RLVTT / zx_pkh(i)                 /RLVTT / zx_pkh(i)
1390         ELSE         ELSE
1391            IF (tsurf(i).LT.t_coup) THEN            IF (tsurf(i).LT.t_coup) THEN
1392               zx_qs = qsats(tsurf(i)) / ps(i)               zx_qs = qsats(tsurf(i)) / ps(i)
1393               zx_dq_s_dh = dqsats(tsurf(i), zx_qs)/RLVTT &               zx_dq_s_dh = dqsats(tsurf(i), zx_qs)/RLVTT &
1394                    &                    / zx_pkh(i)                    / zx_pkh(i)
1395            ELSE            ELSE
1396               zx_qs = qsatl(tsurf(i)) / ps(i)               zx_qs = qsatl(tsurf(i)) / ps(i)
1397               zx_dq_s_dh = dqsatl(tsurf(i), zx_qs)/RLVTT &               zx_dq_s_dh = dqsatl(tsurf(i), zx_qs)/RLVTT &
1398                    &               / zx_pkh(i)                    / zx_pkh(i)
1399            ENDIF            ENDIF
1400         ENDIF         ENDIF
1401         zx_dq_s_dt(i) = RCPD * zx_pkh(i) * zx_dq_s_dh         zx_dq_s_dt(i) = RCPD * zx_pkh(i) * zx_dq_s_dh
1402         zx_qsat(i) = zx_qs         zx_qsat(i) = zx_qs
1403         zx_coef(i) = coef1lay(i) &         zx_coef(i) = coef1lay(i) &
1404              & * (1.0+SQRT(u1lay(i)**2+v1lay(i)**2)) &              * (1.0+SQRT(u1lay(i)**2+v1lay(i)**2)) &
1405              & * p1lay(i)/(RD*t1lay(i))              * p1lay(i)/(RD*t1lay(i))
1406    
1407      ENDDO      ENDDO
1408    
# Line 1422  CONTAINS Line 1420  CONTAINS
1420         ! Q         ! Q
1421         zx_oq(i) = 1. - (beta(i) * zx_k1(i) * peqBcoef(i) * dtime)         zx_oq(i) = 1. - (beta(i) * zx_k1(i) * peqBcoef(i) * dtime)
1422         zx_mq(i) = beta(i) * zx_k1(i) * &         zx_mq(i) = beta(i) * zx_k1(i) * &
1423              &             (peqAcoef(i) - zx_qsat(i) &              (peqAcoef(i) - zx_qsat(i) &
1424              &                          + zx_dq_s_dt(i) * tsurf(i)) &              + zx_dq_s_dt(i) * tsurf(i)) &
1425              &             / zx_oq(i)              / zx_oq(i)
1426         zx_nq(i) = beta(i) * zx_k1(i) * (-1. * zx_dq_s_dt(i)) &         zx_nq(i) = beta(i) * zx_k1(i) * (-1. * zx_dq_s_dt(i)) &
1427              &                              / zx_oq(i)              / zx_oq(i)
1428    
1429         ! H         ! H
1430         zx_oh(i) = 1. - (zx_k1(i) * petBcoef(i) * dtime)         zx_oh(i) = 1. - (zx_k1(i) * petBcoef(i) * dtime)
# Line 1435  CONTAINS Line 1433  CONTAINS
1433    
1434         ! Tsurface         ! Tsurface
1435         tsurf_new(i) = (tsurf(i) + cal(i)/(RCPD * zx_pkh(i)) * dtime * &         tsurf_new(i) = (tsurf(i) + cal(i)/(RCPD * zx_pkh(i)) * dtime * &
1436              &             (radsol(i) + zx_mh(i) + zx_sl(i) * zx_mq(i)) &              (radsol(i) + zx_mh(i) + zx_sl(i) * zx_mq(i)) &
1437              &                 + dif_grnd(i) * t_grnd * dtime)/ &              + dif_grnd(i) * t_grnd * dtime)/ &
1438              &          ( 1. - dtime * cal(i)/(RCPD * zx_pkh(i)) * ( &              ( 1. - dtime * cal(i)/(RCPD * zx_pkh(i)) * ( &
1439              &                       zx_nh(i) + zx_sl(i) * zx_nq(i)) &                zx_nh(i) + zx_sl(i) * zx_nq(i)) &
1440              &                     + dtime * dif_grnd(i))              + dtime * dif_grnd(i))
1441    
1442    
1443         ! Y'a-t-il fonte de neige?         ! Y'a-t-il fonte de neige?
1444    
1445         !    fonte_neige = (nisurf /= is_oce) .AND. &         ! fonte_neige = (nisurf /= is_oce) .AND. &
1446         !     & (snow(i) > epsfra .OR. nisurf == is_sic .OR. nisurf == is_lic) &         ! & (snow(i) > epsfra .OR. nisurf == is_sic .OR. nisurf == is_lic) &
1447         !     & .AND. (tsurf_new(i) >= RTT)         ! & .AND. (tsurf_new(i) >= RTT)
1448         !    if (fonte_neige) tsurf_new(i) = RTT           ! if (fonte_neige) tsurf_new(i) = RTT
1449         d_ts(i) = tsurf_new(i) - tsurf(i)         d_ts(i) = tsurf_new(i) - tsurf(i)
1450         !    zx_h_ts(i) = tsurf_new(i) * RCPD * zx_pkh(i)         ! zx_h_ts(i) = tsurf_new(i) * RCPD * zx_pkh(i)
1451         !    zx_q_0(i) = zx_qsat(i) + zx_dq_s_dt(i) * d_ts(i)         ! zx_q_0(i) = zx_qsat(i) + zx_dq_s_dt(i) * d_ts(i)
1452         !== flux_q est le flux de vapeur d'eau: kg/(m**2 s)  positive vers bas         !== flux_q est le flux de vapeur d'eau: kg/(m**2 s) positive vers bas
1453         !== flux_t est le flux de cpt (energie sensible): j/(m**2 s)         !== flux_t est le flux de cpt (energie sensible): j/(m**2 s)
1454         evap(i) = - zx_mq(i) - zx_nq(i) * tsurf_new(i)         evap(i) = - zx_mq(i) - zx_nq(i) * tsurf_new(i)
1455         fluxlat(i) = - evap(i) * zx_sl(i)         fluxlat(i) = - evap(i) * zx_sl(i)
# Line 1469  CONTAINS Line 1467  CONTAINS
1467    
1468    !************************    !************************
1469    
1470    SUBROUTINE fonte_neige( klon, knon, nisurf, dtime, &    SUBROUTINE fonte_neige( klon, knon, nisurf, dtime,  &
1471         & tsurf, p1lay, cal, beta, coef1lay, ps, &         tsurf, p1lay, cal, beta, coef1lay, ps,  &
1472         & precip_rain, precip_snow, snow, qsol, &         precip_rain, precip_snow, snow, qsol,  &
1473         & radsol, dif_grnd, t1lay, q1lay, u1lay, v1lay, &         radsol, dif_grnd, t1lay, q1lay, u1lay, v1lay,  &
1474         & petAcoef, peqAcoef, petBcoef, peqBcoef, &         petAcoef, peqAcoef, petBcoef, peqBcoef,  &
1475         & tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l,  &
1476         & fqcalving, ffonte, run_off_lic_0)         fqcalving, ffonte, run_off_lic_0)
1477    
1478      ! Routine de traitement de la fonte de la neige dans le cas du traitement      ! Routine de traitement de la fonte de la neige dans le cas du traitement
1479      ! de sol simplifié      ! de sol simplifié
1480    
1481      ! LF 03/2001      ! LF 03/2001
1482      ! input:      ! input:
1483      !   knon         nombre de points a traiter      ! knon nombre de points a traiter
1484      !   nisurf       surface a traiter      ! nisurf surface a traiter
1485      !   tsurf        temperature de surface      ! tsurf temperature de surface
1486      !   p1lay        pression 1er niveau (milieu de couche)      ! p1lay pression 1er niveau (milieu de couche)
1487      !   cal          capacite calorifique du sol      ! cal capacite calorifique du sol
1488      !   beta         evap reelle      ! beta evap reelle
1489      !   coef1lay     coefficient d'echange      ! coef1lay coefficient d'echange
1490      !   ps           pression au sol      ! ps pression au sol
1491      !   precip_rain  precipitations liquides      ! precip_rain precipitations liquides
1492      !   precip_snow  precipitations solides      ! precip_snow precipitations solides
1493      !   snow         champs hauteur de neige      ! snow champs hauteur de neige
1494      !   qsol         hauteur d'eau contenu dans le sol      ! qsol hauteur d'eau contenu dans le sol
1495      !   runoff       runoff en cas de trop plein      ! runoff runoff en cas de trop plein
1496      !   petAcoef     coeff. A de la resolution de la CL pour t      ! petAcoef coeff. A de la resolution de la CL pour t
1497      !   peqAcoef     coeff. A de la resolution de la CL pour q      ! peqAcoef coeff. A de la resolution de la CL pour q
1498      !   petBcoef     coeff. B de la resolution de la CL pour t      ! petBcoef coeff. B de la resolution de la CL pour t
1499      !   peqBcoef     coeff. B de la resolution de la CL pour q      ! peqBcoef coeff. B de la resolution de la CL pour q
1500      !   radsol       rayonnement net aus sol (LW + SW)      ! radsol rayonnement net aus sol (LW + SW)
1501      !   dif_grnd     coeff. diffusion vers le sol profond      ! dif_grnd coeff. diffusion vers le sol profond
1502    
1503      ! output:      ! output:
1504      !   tsurf_new    temperature au sol      ! tsurf_new temperature au sol
1505      !   fluxsens     flux de chaleur sensible      ! fluxsens flux de chaleur sensible
1506      !   fluxlat      flux de chaleur latente      ! fluxlat flux de chaleur latente
1507      !   dflux_s      derivee du flux de chaleur sensible / Ts      ! dflux_s derivee du flux de chaleur sensible / Ts
1508      !   dflux_l      derivee du flux de chaleur latente  / Ts      ! dflux_l derivee du flux de chaleur latente / Ts
1509      ! in/out:      ! in/out:
1510      !   run_off_lic_0 run off glacier du pas de temps précedent      ! run_off_lic_0 run off glacier du pas de temps précedent
1511    
1512    
1513      use indicesol      use indicesol
1514      use YOMCST      use SUPHEC_M
1515      use yoethf      use yoethf_m
1516      use fcttre      use fcttre
1517      !IM cf JLD      !IM cf JLD
1518    
1519      ! Parametres d'entree      ! Parametres d'entree
1520      integer, intent(IN) :: knon, nisurf, klon      integer, intent(IN) :: knon, nisurf, klon
1521      real   , intent(IN) :: dtime      real , intent(IN) :: dtime
1522      real, dimension(klon), intent(IN) :: petAcoef, peqAcoef      real, dimension(klon), intent(IN) :: petAcoef, peqAcoef
1523      real, dimension(klon), intent(IN) :: petBcoef, peqBcoef      real, dimension(klon), intent(IN) :: petBcoef, peqBcoef
1524      real, dimension(klon), intent(IN) :: ps, q1lay      real, dimension(klon), intent(IN) :: ps, q1lay
# Line 1542  CONTAINS Line 1540  CONTAINS
1540      ! Variables locales      ! Variables locales
1541      ! Masse maximum de neige (kg/m2). Au dessus de ce seuil, la neige      ! Masse maximum de neige (kg/m2). Au dessus de ce seuil, la neige
1542      ! en exces "s'ecoule" (calving)      ! en exces "s'ecoule" (calving)
1543      !  real, parameter :: snow_max=1.      ! real, parameter :: snow_max=1.
1544      !IM cf JLD/GK      !IM cf JLD/GK
1545      real, parameter :: snow_max=3000.      real, parameter :: snow_max=3000.
1546      integer :: i      integer :: i
# Line 1551  CONTAINS Line 1549  CONTAINS
1549      real, dimension(klon) :: zx_pkh, zx_dq_s_dt, zx_qsat, zx_coef      real, dimension(klon) :: zx_pkh, zx_dq_s_dt, zx_qsat, zx_coef
1550      real, dimension(klon) :: zx_sl, zx_k1      real, dimension(klon) :: zx_sl, zx_k1
1551      real, dimension(klon) :: zx_q_0 , d_ts      real, dimension(klon) :: zx_q_0 , d_ts
1552      real                  :: zdelta, zcvm5, zx_qs, zcor, zx_dq_s_dh      real :: zdelta, zcvm5, zx_qs, zcor, zx_dq_s_dh
1553      real                  :: bilan_f, fq_fonte      real :: bilan_f, fq_fonte
1554      REAL                  :: subli, fsno      REAL :: subli, fsno
1555      REAL, DIMENSION(klon) :: bil_eau_s, snow_evap      REAL, DIMENSION(klon) :: bil_eau_s, snow_evap
1556      real, parameter :: t_grnd = 271.35, t_coup = 273.15      real, parameter :: t_grnd = 271.35, t_coup = 273.15
1557      !! PB temporaire en attendant mieux pour le modele de neige      !! PB temporaire en attendant mieux pour le modele de neige
# Line 1563  CONTAINS Line 1561  CONTAINS
1561      REAL, parameter :: chaice = 3.334E+05/(2.3867E+06*0.15)      REAL, parameter :: chaice = 3.334E+05/(2.3867E+06*0.15)
1562      ! fin GKtest      ! fin GKtest
1563    
1564      logical, save         :: check = .FALSE.      logical, save :: check = .FALSE.
1565      character (len = 20)  :: modname = 'fonte_neige'      character (len = 20) :: modname = 'fonte_neige'
1566      logical, save         :: neige_fond = .false.      logical, save :: neige_fond = .false.
1567      real, save            :: max_eau_sol = 150.0      real, save :: max_eau_sol = 150.0
1568      character (len = 80) :: abort_message      character (len = 80) :: abort_message
1569      logical, save         :: first = .true., second=.false.      logical, save :: first = .true., second=.false.
1570      real                 :: coeff_rel      real :: coeff_rel
1571    
1572      if (check) write(*, *)'Entree ', modname, ' surface = ', nisurf      if (check) write(*, *)'Entree ', modname, ' surface = ', nisurf
1573    
# Line 1587  CONTAINS Line 1585  CONTAINS
1585            zcor=1./(1.-retv*zx_qs)            zcor=1./(1.-retv*zx_qs)
1586            zx_qs=zx_qs*zcor            zx_qs=zx_qs*zcor
1587            zx_dq_s_dh = FOEDE(tsurf(i), zdelta, zcvm5, zx_qs, zcor) &            zx_dq_s_dh = FOEDE(tsurf(i), zdelta, zcvm5, zx_qs, zcor) &
1588                 &                 /RLVTT / zx_pkh(i)                 /RLVTT / zx_pkh(i)
1589         ELSE         ELSE
1590            IF (tsurf(i).LT.t_coup) THEN            IF (tsurf(i).LT.t_coup) THEN
1591               zx_qs = qsats(tsurf(i)) / ps(i)               zx_qs = qsats(tsurf(i)) / ps(i)
1592               zx_dq_s_dh = dqsats(tsurf(i), zx_qs)/RLVTT &               zx_dq_s_dh = dqsats(tsurf(i), zx_qs)/RLVTT &
1593                    &                    / zx_pkh(i)                    / zx_pkh(i)
1594            ELSE            ELSE
1595               zx_qs = qsatl(tsurf(i)) / ps(i)               zx_qs = qsatl(tsurf(i)) / ps(i)
1596               zx_dq_s_dh = dqsatl(tsurf(i), zx_qs)/RLVTT &               zx_dq_s_dh = dqsatl(tsurf(i), zx_qs)/RLVTT &
1597                    &               / zx_pkh(i)                    / zx_pkh(i)
1598            ENDIF            ENDIF
1599         ENDIF         ENDIF
1600         zx_dq_s_dt(i) = RCPD * zx_pkh(i) * zx_dq_s_dh         zx_dq_s_dt(i) = RCPD * zx_pkh(i) * zx_dq_s_dh
1601         zx_qsat(i) = zx_qs         zx_qsat(i) = zx_qs
1602         zx_coef(i) = coef1lay(i) &         zx_coef(i) = coef1lay(i) &
1603              & * (1.0+SQRT(u1lay(i)**2+v1lay(i)**2)) &              * (1.0+SQRT(u1lay(i)**2+v1lay(i)**2)) &
1604              & * p1lay(i)/(RD*t1lay(i))              * p1lay(i)/(RD*t1lay(i))
1605      ENDDO      ENDDO
1606    
1607      ! === Calcul de la temperature de surface ===      ! === Calcul de la temperature de surface ===
# Line 1620  CONTAINS Line 1618  CONTAINS
1618         ! Q         ! Q
1619         zx_oq(i) = 1. - (beta(i) * zx_k1(i) * peqBcoef(i) * dtime)         zx_oq(i) = 1. - (beta(i) * zx_k1(i) * peqBcoef(i) * dtime)
1620         zx_mq(i) = beta(i) * zx_k1(i) * &         zx_mq(i) = beta(i) * zx_k1(i) * &
1621              &             (peqAcoef(i) - zx_qsat(i) &              (peqAcoef(i) - zx_qsat(i) &
1622              &                          + zx_dq_s_dt(i) * tsurf(i)) &              + zx_dq_s_dt(i) * tsurf(i)) &
1623              &             / zx_oq(i)              / zx_oq(i)
1624         zx_nq(i) = beta(i) * zx_k1(i) * (-1. * zx_dq_s_dt(i)) &         zx_nq(i) = beta(i) * zx_k1(i) * (-1. * zx_dq_s_dt(i)) &
1625              &                              / zx_oq(i)              / zx_oq(i)
1626    
1627         ! H         ! H
1628         zx_oh(i) = 1. - (zx_k1(i) * petBcoef(i) * dtime)         zx_oh(i) = 1. - (zx_k1(i) * petBcoef(i) * dtime)
# Line 1640  CONTAINS Line 1638  CONTAINS
1638         snow = MAX(0.0, snow)         snow = MAX(0.0, snow)
1639      end where      end where
1640    
1641      !  bil_eau_s = bil_eau_s + (precip_rain * dtime) - (evap - snow_evap) * dtime      ! bil_eau_s = bil_eau_s + (precip_rain * dtime) - (evap - snow_evap) * dtime
1642      bil_eau_s = (precip_rain * dtime) - (evap - snow_evap) * dtime      bil_eau_s = (precip_rain * dtime) - (evap - snow_evap) * dtime
1643    
1644    
# Line 1649  CONTAINS Line 1647  CONTAINS
1647      ffonte=0.      ffonte=0.
1648      do i = 1, knon      do i = 1, knon
1649         neige_fond = ((snow(i) > epsfra .OR. nisurf == is_sic .OR. nisurf == is_lic) &         neige_fond = ((snow(i) > epsfra .OR. nisurf == is_sic .OR. nisurf == is_lic) &
1650              & .AND. tsurf_new(i) >= RTT)              .AND. tsurf_new(i) >= RTT)
1651         if (neige_fond) then         if (neige_fond) then
1652            fq_fonte = MIN( MAX((tsurf_new(i)-RTT )/chasno, 0.0), snow(i))            fq_fonte = MIN( MAX((tsurf_new(i)-RTT )/chasno, 0.0), snow(i))
1653            ffonte(i) = fq_fonte * RLMLT/dtime            ffonte(i) = fq_fonte * RLMLT/dtime
1654            snow(i) = max(0., snow(i) - fq_fonte)            snow(i) = max(0., snow(i) - fq_fonte)
1655            bil_eau_s(i) = bil_eau_s(i) + fq_fonte            bil_eau_s(i) = bil_eau_s(i) + fq_fonte
1656            tsurf_new(i) = tsurf_new(i) - fq_fonte * chasno              tsurf_new(i) = tsurf_new(i) - fq_fonte * chasno
1657            !IM cf JLD OK                !IM cf JLD OK
1658            !IM cf JLD/ GKtest fonte aussi pour la glace            !IM cf JLD/ GKtest fonte aussi pour la glace
1659            IF (nisurf == is_sic .OR. nisurf == is_lic ) THEN            IF (nisurf == is_sic .OR. nisurf == is_lic ) THEN
1660               fq_fonte = MAX((tsurf_new(i)-RTT )/chaice, 0.0)               fq_fonte = MAX((tsurf_new(i)-RTT )/chaice, 0.0)
# Line 1667  CONTAINS Line 1665  CONTAINS
1665            d_ts(i) = tsurf_new(i) - tsurf(i)            d_ts(i) = tsurf_new(i) - tsurf(i)
1666         endif         endif
1667    
1668         !   s'il y a une hauteur trop importante de neige, elle s'coule         ! s'il y a une hauteur trop importante de neige, elle s'coule
1669         fqcalving(i) = max(0., snow(i) - snow_max)/dtime         fqcalving(i) = max(0., snow(i) - snow_max)/dtime
1670         snow(i)=min(snow(i), snow_max)         snow(i)=min(snow(i), snow_max)
1671    
# Line 1676  CONTAINS Line 1674  CONTAINS
1674            run_off(i) = run_off(i) + MAX(qsol(i) - max_eau_sol, 0.0)            run_off(i) = run_off(i) + MAX(qsol(i) - max_eau_sol, 0.0)
1675            qsol(i) = MIN(qsol(i), max_eau_sol)            qsol(i) = MIN(qsol(i), max_eau_sol)
1676         else if (nisurf == is_lic) then         else if (nisurf == is_lic) then
1677            run_off_lic(i) = (coeff_rel *  fqcalving(i)) + &            run_off_lic(i) = (coeff_rel * fqcalving(i)) + &
1678                 &                        (1. - coeff_rel) * run_off_lic_0(i)                 (1. - coeff_rel) * run_off_lic_0(i)
1679            run_off_lic_0(i) = run_off_lic(i)            run_off_lic_0(i) = run_off_lic(i)
1680            run_off_lic(i) = run_off_lic(i) + bil_eau_s(i)/dtime            run_off_lic(i) = run_off_lic(i) + bil_eau_s(i)/dtime
1681         endif         endif
Legend:
Removed from v.14  
changed lines
  Added in v.38
  ViewVC Help
Powered by ViewVC 1.1.21