Changeset 441 for LMDZ.3.3/branches


Ignore:
Timestamp:
Jan 27, 2003, 11:26:36 AM (21 years ago)
Author:
lmdzadmin
Message:

Remplace qsol par qsurf JLD
Rajout zlev, swdown_vrai, albedo_keep JP
IM

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ.3.3/branches/rel-LF/libf/phylmd/interface_surf.F90

    r398 r441  
    5757      & klon, iim, jjm, nisurf, knon, knindex, pctsrf, &
    5858      & rlon, rlat, cufi, cvfi,&
    59       & debut, lafin, ok_veget, soil_model, nsoilmx, tsoil,&
     59      & debut, lafin, ok_veget, soil_model, nsoilmx, tsoil, qsol,&
    6060      & zlev,  u1_lay, v1_lay, temp_air, spechum, epot_air, ccanopy, &
    6161      & tq_cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
    6262      & precip_rain, precip_snow, sollw, sollwdown, swnet, swdown, &
    6363      & fder, taux, tauy, rugos, rugoro, &
    64       & albedo, snow, qsol, &
     64      & albedo, snow, qsurf, &
    6565      & tsurf, p1lay, ps, radsol, &
    6666      & ocean, npas, nexca, zmasq, &
    … …  
    171171  character (len = 6)  :: ocean
    172172  integer              :: npas, nexca ! nombre et pas de temps couplage
    173   real, dimension(klon), intent(INOUT) :: evap, snow, qsol
     173  real, dimension(klon), intent(INOUT) :: evap, snow, qsurf
    174174!! PB ajout pour soil
    175175  logical          :: soil_model
    176176  integer          :: nsoilmx
    177177  REAL, DIMENSION(klon, nsoilmx) :: tsoil
     178  REAL, dimension(klon), intent(INOUT) :: qsol
    178179  REAL, dimension(klon)          :: soilcap
    179180  REAL, dimension(klon)          :: soilflux
    … …  
    193194  integer, save        :: error
    194195  integer              :: ii, index
    195   logical,save              :: check = .false.
     196  logical,save              :: check = .true.
    196197  real, dimension(klon):: cal, beta, dif_grnd, capsol
    197198!!$PB  real, parameter      :: calice=1.0/(5.1444e+06*0.15), tau_gl=86400.*5.
    … …  
    200201  real, dimension(klon):: alb_ice
    201202  real, dimension(klon):: tsurf_temp
    202   real, dimension(klon):: qsol_new
     203  real, dimension(klon):: qsurf_new
    203204!!  real, allocatable, dimension(:), save :: alb_neig_grid
    204205  real, dimension(klon):: alb_neig, alb_eau
    … …  
    309310     & alb_new, z0_new)
    310311! 
    311 ! calcul snow et qsol, hydrol adapté
     312! calcul snow et qsurf, hydrol adapté
    312313!
    313314       CALL calbeta(dtime, nisurf, knon, snow, qsol, beta, capsol, dif_grnd)
    … …  
    323324       CALL calcul_fluxs( klon, knon, nisurf, dtime, &
    324325     &   tsurf, p1lay, cal, beta, tq_cdrag, ps, &
    325      &   precip_rain, precip_snow, snow, qsol,  &
     326     &   precip_rain, precip_snow, snow, qsurf,  &
    326327     &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
    327328     &   petAcoef, peqAcoef, petBcoef, peqBcoef, &
    … …  
    358359     &  evap, fluxsens, fluxlat, &             
    359360     &  tsol_rad, tsurf_new, alb_new, alblw, &
    360      &  emis_new, z0_new, dflux_l, dflux_s, qsol_new)
     361     &  emis_new, z0_new, dflux_l, dflux_s, qsurf_new)
    361362
    362363! 
    … …  
    366367!
    367368! mise a jour de l'humidite saturante calculee par ORCHIDEE
    368       qsol(1:knon) = qsol_new(1:knon)
     369      qsurf(1:knon) = qsurf_new(1:knon)
    369370
    370371    endif   
    … …  
    416417    call calcul_fluxs( klon, knon, nisurf, dtime, &
    417418     &   tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
    418      &   precip_rain, precip_snow, snow, qsol,  &
     419     &   precip_rain, precip_snow, snow, qsurf,  &
    419420     &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
    420421     &   petAcoef, peqAcoef, petBcoef, peqBcoef, &
    … …  
    511512    CALL calcul_fluxs( klon, knon, nisurf, dtime, &
    512513         &   tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
    513          &   precip_rain, precip_snow, snow, qsol,  &
     514         &   precip_rain, precip_snow, snow, qsurf,  &
    514515         &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
    515516         &   petAcoef, peqAcoef, petBcoef, peqBcoef, &
    … …  
    581582    call calcul_fluxs( klon, knon, nisurf, dtime, &
    582583     &   tsurf, p1lay, cal, beta, tq_cdrag, ps, &
    583      &   precip_rain, precip_snow, snow, qsol,  &
     584     &   precip_rain, precip_snow, snow, qsurf,  &
    584585     &   radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
    585586     &   petAcoef, peqAcoef, petBcoef, peqBcoef, &
    … …  
    651652     & knindex, rlon, rlat, cufi, cvfi, iim, jjm, pctsrf, &
    652653     & debut, lafin, ok_veget, &
    653      & zlev,  u1_lay, v1_lay, temp_air, spechum, epot_air, ccanopy, &
     654     & plev,  u1_lay, v1_lay, temp_air, spechum, epot_air, ccanopy, &
    654655     & tq_cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
    655656     & precip_rain, precip_snow, lwdown, swnet, swdown, &
    … …  
    680681!   ok_veget     logical: appel ou non au schema de surface continental
    681682!                     (si false calcul simplifie des fluxs sur les continents)
    682 !   zlev         hauteur de la premiere couche      
     683!   plev         hauteur de la premiere couche (Pa)     
    683684!   u1_lay       vitesse u 1ere couche
    684685!   v1_lay       vitesse v 1ere couche
    … …  
    715716!   emis_new     emissivite
    716717!   z0_new       surface roughness
    717 !   qsurf        saturated air moisture at surface
     718!   qsurf        air moisture at surface
    718719
    719720! Parametres d'entree
    … …  
    730731  real, dimension(klon), intent(IN) :: rlon, rlat
    731732  real, dimension(klon), intent(IN) :: cufi, cvfi
    732   real, dimension(klon), intent(IN) :: zlev
     733  real, dimension(klon), intent(IN) :: plev
    733734  real, dimension(klon), intent(IN) :: u1_lay, v1_lay
    734735  real, dimension(klon), intent(IN) :: temp_air, spechum
    … …  
    739740  real, dimension(klon), intent(IN) :: precip_rain, precip_snow
    740741  real, dimension(klon), intent(IN) :: lwdown, swnet, swdown, ps
     742!IM cf. JP +++
     743  real, dimension(klon) :: swdown_vrai
     744!IM cf. JP ---
    741745  real, dimension(klon), intent(IN) :: tsurf, p1lay
    742746  real, dimension(klon), intent(IN) :: radsol
    … …  
    759763! drapeaux controlant les appels dans SECHIBA
    760764!  type(control_type), save   :: control_in
     765! Preserved albedo
     766!IM cf. JP +++
     767  real, allocatable, dimension(:), save :: albedo_keep, zlev
     768!IM cf. JP ---
    761769! coordonnees geographiques
    762770  real, allocatable, dimension(:,:), save :: lalo
    … …  
    794802
    795803#include "temps.inc"
     804#include "YOMCST.inc"
    796805
    797806  if (check) write(*,*)'Entree ', modname
    … …  
    803812  if (debut) then
    804813
     814  IF ( .NOT. allocated(albedo_keep)) THEN
     815     ALLOCATE(albedo_keep(klon))
     816     ALLOCATE(zlev(klon))
     817  ENDIF
    805818! Pb de correspondances de grilles
    806819   allocate(ig(klon))
    … …  
    9901003  cdrag(1:knon) = tq_cdrag(1:knon)
    9911004
     1005!IM cf. JP +++
     1006! zlev(1:knon) = (100.*plev(1:knon))/((ps(1:knon)/287.05*temp_air(1:knon))*9.80665)
     1007  zlev(1:knon) = (100.*plev(1:knon))/((ps(1:knon)/RD*temp_air(1:knon))*RG)
     1008!IM cf. JP ---
     1009
    9921010   where(cdrag > 0.01)
    9931011     cdrag = 0.01
    … …  
    10081026     & tsol_rad, tsurf_new, qsurf, albedo_out, emis_new, z0_new, &
    10091027     & lon_scat, lat_scat)
     1028
     1029!IM cf. JP +++
     1030    albedo_keep(:) = (albedo_out(:,1)+albedo_out(:,2))/2.
     1031!IM cf. JP ---
     1032
    10101033  endif
     1034
     1035!IM cf. JP +++
     1036  swdown_vrai(:) = swnet(:)/(1. - albedo_keep(:))
     1037!IM cf. JP ---
    10111038
    10121039  call intersurf_main (itime+itau_phy, iim, jjm+1, knon, ktindex, dtime, &
    … …  
    10151042     & zlev,  u1_lay, v1_lay, spechum, temp_air, epot_air, ccanopy, &
    10161043     & cdrag, petA_orc, peqA_orc, petB_orc, peqB_orc, &
    1017      & precip_rain, precip_snow, lwdown, swnet, swdown, ps, &
     1044!IM cf. JP +++
     1045     & precip_rain, precip_snow, lwdown, swnet, swdown_vrai, ps, &
     1046!IM cf. JP ---
    10181047     & evap, fluxsens, fluxlat, coastalflow, riverflow, &
    10191048     & tsol_rad, tsurf_new, qsurf, albedo_out, emis_new, z0_new, &
    10201049     & lon_scat, lat_scat)
     1050
     1051!IM cf. JP +++
     1052    albedo_keep(:) = (albedo_out(:,1)+albedo_out(:,2))/2.
     1053!IM cf. JP ---
    10211054
    10221055    bidule=0.
    … …  
    19271960  character (len = 20),save :: fich ='limit.nc'
    19281961  logical,save     :: newlmt = .false.
    1929   logical,save     :: check = .FALSE.
     1962  logical,save     :: check = .true.
    19301963! Champs lus dans le fichier de CL
    19311964  real, allocatable , save, dimension(:) :: alb_lu, rug_lu
    … …  
    19511984  if (check) write(*,*)modname,':: jour_lu_sur, deja_lu_sur', jour_lu_sur, deja_lu_sur
    19521985  if (check) write(*,*)modname,':: itime, lmt_pas', itime, lmt_pas
    1953   call flush(6)
     1986  if (check) call flush(6)
    19541987
    19551988! Tester d'abord si c'est le moment de lire le fichier
    … …  
    19591992!
    19601993    fich = trim(fich)
     1994    IF (check) WRITE(*,*)modname,' ouverture fichier ',fich
     1995    if (check) CALL flush(6)
    19611996    ierr = NF_OPEN (fich, NF_NOWRITE,nid)
    19621997    if (ierr.NE.NF_NOERR) then
    … …  
    20332068  SUBROUTINE calcul_fluxs( klon, knon, nisurf, dtime, &
    20342069     & tsurf, p1lay, cal, beta, coef1lay, ps, &
    2035      & precip_rain, precip_snow, snow, qsol, &
     2070     & precip_rain, precip_snow, snow, qsurf, &
    20362071     & radsol, dif_grnd, t1lay, q1lay, u1lay, v1lay, &
    20372072     & petAcoef, peqAcoef, petBcoef, peqBcoef, &
    … …  
    20552090!   precip_snow  precipitations solides
    20562091!   snow         champs hauteur de neige
    2057 !   qsol         humidite du sol
    20582092!   runoff       runoff en cas de trop plein
    20592093!   petAcoef     coeff. A de la resolution de la CL pour t
    … …  
    20662100! output:
    20672101!   tsurf_new    temperature au sol
     2102!   qsurf        humidite de l'air au dessus du sol
    20682103!   fluxsens     flux de chaleur sensible
    20692104!   fluxlat      flux de chaleur latente
    … …  
    20862121  real, dimension(klon), intent(IN) :: radsol, dif_grnd
    20872122  real, dimension(klon), intent(IN) :: t1lay, u1lay, v1lay
    2088   real, dimension(klon), intent(INOUT) :: snow, qsol
     2123  real, dimension(klon), intent(INOUT) :: snow, qsurf
    20892124
    20902125! Parametres sorties
    … …  
    21022137  real                  :: bilan_f, fq_fonte
    21032138  REAL                  :: subli, fsno
     2139  REAL                  :: qsat_new, q1_new
    21042140  real, parameter :: t_grnd = 271.35, t_coup = 273.15
    21052141!! PB temporaire en attendant mieux pour le modele de neige
    21062142  REAL, parameter :: chasno = 3.334E+05/(2.3867E+06*0.15)
    21072143!
    2108   logical, save         :: check = .FALSE.
     2144  logical, save         :: check = .true.
    21092145  character (len = 20)  :: modname = 'calcul_fluxs'
    21102146  logical, save         :: fonte_neige = .false.
    … …  
    21142150
    21152151  if (check) write(*,*)'Entree ', modname,' surface = ',nisurf
     2152
     2153  IF (check) THEN
     2154      WRITE(*,*)' radsol (min, max)' &
     2155         &     , MINVAL(radsol(1:knon)), MAXVAL(radsol(1:knon))
     2156      CALL flush(6)
     2157  ENDIF
    21162158
    21172159  if (size(coastalflow) /= knon .AND. nisurf == is_ter) then
    … …  
    21352177!!$      where (precip_rain > 0.) qsol = qsol + (precip_rain - evap) * dtime
    21362178!!$    endif
    2137     IF (nisurf /= is_ter) qsol = max_eau_sol
     2179!!$    IF (nisurf /= is_ter) qsol = max_eau_sol
    21382180
    21392181
    … …  
    22332275    dflux_s(i) = zx_nh(i)
    22342276    dflux_l(i) = (zx_sl(i) * zx_nq(i))
    2235 
     2277! Nouvelle valeure de l'humidite au dessus du sol
     2278    qsat_new=zx_qsat(i) + zx_dq_s_dt(i) * d_ts(i)
     2279    q1_new = peqAcoef(i) - peqBcoef(i)*evap(i)*dtime
     2280    qsurf(i)=q1_new*(1.-beta(i)) + beta(i)*qsat_new
    22362281!
    22372282! en cas de fonte de neige
    … …  
    24122457!   precip_snow  precipitations solides
    24132458!   snow         champs hauteur de neige
    2414 !   qsol         humidite du sol
     2459!   qsol         hauteur d'eau contenu dans le sol
    24152460!   runoff       runoff en cas de trop plein
    24162461!   petAcoef     coeff. A de la resolution de la CL pour t
    … …  
    24592504  real                  :: bilan_f, fq_fonte
    24602505  REAL                  :: subli, fsno
     2506  real, dimension(klon) :: bil_eau_s
    24612507  real, parameter :: t_grnd = 271.35, t_coup = 273.15
    24622508!! PB temporaire en attendant mieux pour le modele de neige
    … …  
    24732519
    24742520! Initialisations
     2521  bil_eau_s(:) = 0.
    24752522  DO i = 1, knon
    24762523    zx_pkh(i) = (ps(i)/ps(i))**RKAPPA
    … …  
    25342581  WHERE (precip_snow > 0.) snow = snow + (precip_snow * dtime)
    25352582  WHERE (evap > 0 ) snow = MAX(0.0, snow - (evap * dtime))
    2536   qsol = qsol + (precip_rain - evap) * dtime
     2583  bil_eau_s = bil_eau_s + (precip_rain - evap) * dtime
    25372584!
    25382585! Y'a-t-il fonte de neige?
    … …  
    25442591      fq_fonte = MIN( MAX((tsurf_new(i)-RTT )/chasno,0.0),snow(i))
    25452592      snow(i) = max(0., snow(i) - fq_fonte)
    2546       qsol(i) = qsol(i) + fq_fonte
     2593      bil_eau_s(i) = bil_eau_s(i) + fq_fonte
    25472594      tsurf_new(i) = tsurf_new(i) - fq_fonte * chasno 
    25482595      IF (nisurf == is_sic .OR. nisurf == is_lic ) tsurf_new(i) = RTT -1.8
    … …  
    25642611!!$      fq_fonte = bilan_f / zx_sl(i)
    25652612    endif
    2566     IF (nisurf == is_ter)  &
    2567        &  run_off(i) = run_off(i) + MAX(qsol(i) - max_eau_sol, 0.0)
    2568     qsol(i) = MIN(qsol(i), max_eau_sol)
     2613
     2614    IF (nisurf == is_ter) then
     2615      qsol(i) = qsol(i) + bil_eau_s(i)
     2616      run_off(i) = run_off(i) + MAX(qsol(i) - max_eau_sol, 0.0)
     2617      qsol(i) = MIN(qsol(i), max_eau_sol)
     2618    else
     2619      run_off(i) = run_off(i) + MAX(bil_eau_s(i), 0.0)
     2620    endif
    25692621  enddo
    25702622
Note: See TracChangeset for help on using the changeset viewer.