Ignore:
Timestamp:
Mar 29, 2023, 3:14:27 PM (15 months ago)
Author:
lguez
Message:

Sync latest trunk changes to branch LMDZ_ECRad

Location:
LMDZ6/branches/LMDZ_ECRad
Files:
1 deleted
15 edited
9 copied

Legend:

Unmodified
Added
Removed
  • LMDZ6/branches/LMDZ_ECRad

  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/fisrtilp.F90

    r4143 r4482  
    143143  ! Coeffients de fraction lessivee : pour OFF-LINE
    144144  !
    145   REAL, DIMENSION(klon,klev),      INTENT(OUT)  :: pfrac_nucl
    146   REAL, DIMENSION(klon,klev),      INTENT(OUT)  :: pfrac_1nucl
    147   REAL, DIMENSION(klon,klev),      INTENT(OUT)  :: pfrac_impa
     145  REAL, DIMENSION(klon,klev),      INTENT(INOUT)  :: pfrac_nucl
     146  REAL, DIMENSION(klon,klev),      INTENT(INOUT)  :: pfrac_1nucl
     147  REAL, DIMENSION(klon,klev),      INTENT(INOUT)  :: pfrac_impa
    148148  !
    149149  ! Fraction d'aerosols lessivee par impaction et par nucleation
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/isotopes_mod.F90

    r4158 r4482  
    133133
    134134SUBROUTINE iso_init()
    135    USE ioipsl_getin_p_mod, ONLY: getin_p
    136135   USE infotrac_phy,       ONLY: ntiso, niso, getKey
    137136    USE strings_mod,       ONLY: maxlen
     
    175174
    176175   !--- Type of water isotopes:
    177    iso_eau = strIdx(isoName, 'H2[16]O'); CALL msg('iso_eau='//int2str(iso_eau), modname)
    178    iso_HDO = strIdx(isoName, 'H[2]HO');  CALL msg('iso_HDO='//int2str(iso_HDO), modname)
    179    iso_O18 = strIdx(isoName, 'H2[18]O'); CALL msg('iso_O18='//int2str(iso_O18), modname)
    180    iso_O17 = strIdx(isoName, 'H2[17]O'); CALL msg('iso_O17='//int2str(iso_O17), modname)
    181    iso_HTO = strIdx(isoName, 'H[3]HO');  CALL msg('iso_HTO='//int2str(iso_HTO), modname)
    182 
    183    ! initialisation
    184    ! lecture des parametres isotopiques:
    185    ! pour que ca marche en openMP, il faut utiliser getin_p. Car le getin ne peut
    186    ! etre appele que par un thread a la fois, et ca pose tout un tas de problemes,
    187    ! d'ou tout un tas de magouilles comme dans conf_phys_m. A terme, tout le monde
    188    ! lira par getin_p.
     176   iso_eau = strIdx(isoName, 'H216O'); CALL msg('iso_eau='//int2str(iso_eau), modname)
     177   iso_HDO = strIdx(isoName, 'HDO');   CALL msg('iso_HDO='//int2str(iso_HDO), modname)
     178   iso_O18 = strIdx(isoName, 'H218O'); CALL msg('iso_O18='//int2str(iso_O18), modname)
     179   iso_O17 = strIdx(isoName, 'H217O'); CALL msg('iso_O17='//int2str(iso_O17), modname)
     180   iso_HTO = strIdx(isoName, 'HTO');   CALL msg('iso_HTO='//int2str(iso_HTO), modname)
     181
     182   !--- Initialiaation: reading the isotopic parameters.
    189183   CALL get_in('lambda',     lambda_sursat, 0.004)
    190184   CALL get_in('thumxt1',    thumxt1,       0.75*1.2)
     
    252246   ! bugs quand temperature dans ascendances convs est mal calculee
    253247   CALL get_in('cond_temp_env',        cond_temp_env,        .FALSE.)
    254    IF(ANY(isoName == 'H[3]HO')) &
     248   IF(ANY(isoName == 'HTO')) &
    255249   CALL get_in('ok_prod_nucl_tritium', ok_prod_nucl_tritium, .FALSE., .FALSE.)
    256250
     
    339333   USE mod_phys_lmdz_omp_data, ONLY :  is_omp_root
    340334   USE mod_phys_lmdz_transfert_para, ONLY : bcast
    341    CHARACTER(LEN=*),  INTENT(IN)    :: nam
    342    CHARACTER(LEN=*),  INTENT(INOUT) :: val
    343    CHARACTER(LEN=*), INTENT(IN)    :: def
    344    LOGICAL, OPTIONAL, INTENT(IN)    :: lDisp
     335   CHARACTER(LEN=*),           INTENT(IN)    :: nam
     336   CHARACTER(LEN=*),           INTENT(INOUT) :: val
     337   CHARACTER(LEN=*), OPTIONAL, INTENT(IN)    :: def
     338   LOGICAL,          OPTIONAL, INTENT(IN)    :: lDisp
    345339   LOGICAL :: lD
    346340!$OMP BARRIER
    347    IF(.NOT.(is_mpi_root.AND.is_omp_root)) RETURN
    348    val=def; CALL getin(nam,val); CALL bcast(val)
    349    lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp
    350    IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(val))
     341   IF(is_mpi_root.AND.is_omp_root) THEN
     342      IF(PRESENT(def)) val=def; CALL getin(nam,val)
     343      lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp
     344      IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(val))
     345  END IF
     346  CALL bcast(val)
    351347END SUBROUTINE getinp_s
    352348
     
    358354   CHARACTER(LEN=*),  INTENT(IN)    :: nam
    359355   INTEGER,           INTENT(INOUT) :: val
    360    INTEGER,           INTENT(IN)    :: def
     356   INTEGER, OPTIONAL, INTENT(IN)    :: def
    361357   LOGICAL, OPTIONAL, INTENT(IN)    :: lDisp
    362358   LOGICAL :: lD
    363359!$OMP BARRIER
    364    IF(.NOT.(is_mpi_root.AND.is_omp_root)) RETURN
    365    val=def; CALL getin(nam,val); CALL bcast(val)
    366    lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp
    367    IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(int2str(val)))
     360   IF(is_mpi_root.AND.is_omp_root) THEN
     361      IF(PRESENT(def)) val=def; CALL getin(nam,val)
     362      lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp
     363      IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(int2str(val)))
     364  END IF
     365  CALL bcast(val)
    368366END SUBROUTINE getinp_i
    369367
     
    375373   CHARACTER(LEN=*),  INTENT(IN)    :: nam
    376374   REAL,              INTENT(INOUT) :: val
    377    REAL,              INTENT(IN)    :: def
     375   REAL,    OPTIONAL, INTENT(IN)    :: def
    378376   LOGICAL, OPTIONAL, INTENT(IN)    :: lDisp
    379377   LOGICAL :: lD
    380378!$OMP BARRIER
    381    IF(.NOT.(is_mpi_root.AND.is_omp_root)) RETURN
    382    val=def; CALL getin(nam,val); CALL bcast(val)
    383    lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp
    384    IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(real2str(val)))
     379   IF(is_mpi_root.AND.is_omp_root) THEN
     380      IF(PRESENT(def)) val=def; CALL getin(nam,val)
     381      lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp
     382      IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(real2str(val)))
     383  END IF
     384  CALL bcast(val)
    385385END SUBROUTINE getinp_r
    386386
     
    392392   CHARACTER(LEN=*),  INTENT(IN)    :: nam
    393393   LOGICAL,           INTENT(INOUT) :: val
    394    LOGICAL,           INTENT(IN)    :: def
     394   LOGICAL, OPTIONAL, INTENT(IN)    :: def
    395395   LOGICAL, OPTIONAL, INTENT(IN)    :: lDisp
    396396   LOGICAL :: lD
    397397!$OMP BARRIER
    398    IF(.NOT.(is_mpi_root.AND.is_omp_root)) RETURN
    399    val=def; CALL getin(nam,val); CALL bcast(val)
    400    lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp
    401    IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(bool2str(val)))
     398   IF(is_mpi_root.AND.is_omp_root) THEN
     399      IF(PRESENT(def)) val=def; CALL getin(nam,val)
     400      lD=.TRUE.; IF(PRESENT(lDisp)) lD=lDisp
     401      IF(lD) CALL msg(TRIM(nam)//' = '//TRIM(bool2str(val)))
     402  END IF
     403  CALL bcast(val)
    402404END SUBROUTINE getinp_l
    403405
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/isotopes_routines_mod.F90

    r4149 r4482  
    17311731       else
    17321732           write(*,*) 'iso_fractcal 1734: non prévu: ixt=',ixt
    1733 #ifdef ISOVERIF
    1734            stop
    1735 #endif
     1733!#ifdef ISOVERIF
     1734           CALL abort_physic('isotopes_routines_mod', 'iso_fractcal 1734', 1)
     1735!#endif
    17361736       endif
    17371737
     
    19491949       else
    19501950           write(*,*) 'iso_fractcal 1954: non prévu: ixt=',ixt
    1951 #ifdef ISOVERIF
    1952            stop
    1953 #endif
     1951!#ifdef ISOVERIF
     1952          CALL abort_physic('isotopes_routines_mod', 'iso_fractcal 1954', 1)
     1953!#endif
    19541954       endif
    19551955
     
    21852185        write(*,*) 'qsol(i)=',qsol(i),' mais evap(i)=',evap(i)
    21862186        write(*,*) 'q1lay(i)=',q1lay(i)
    2187         stop
     2187        CALL abort_physic('isotopes_routines_mod', 'iso_surf 2187', 1)
    21882188 endif !if (q1lay(i).gt.ridicule) then
    21892189endif !if (qsol(i).gt.ridicule_rain) then
     
    24142414    write(*,*) 'iso_surf>iso_rosee_givre 3189: evap=',evap(i)
    24152415    write(*,*) 'q1lay(i)=',q1lay(i)
    2416     stop
     2416    CALL abort_physic('isotopes_routines_mod', 'iso_surf 2416', 1)
    24172417endif !if (q1lay(i).gt.0) then
    24182418
     
    25522552  write(*,*) 'il,Pqisup(il),Eqi(il),Pqiinf(il)=', &
    25532553&          il,Pqisup(il),Eqi(il),Pqiinf(il)
    2554   stop
     2554  CALL abort_physic('isotopes_routines_mod', 'stewart 2554', 1)
    25552555endif !if (iso_verif_egalite
    25562556enddo !do il=1,ncas
     
    1650816508      end subroutine phyiso_etat0_dur
    1650916509
    16510       subroutine phyiso_etat0_fichier( &
    16511      &           snow,run_off_lic_0, &
    16512      &           xtsnow,xtrun_off_lic_0, &
    16513      &           Rland_ice)
    16514       USE dimphy, only: klon,klev
    16515       !USE mod_grid_phy_lmdz
    16516       !USE mod_phys_lmdz_para
    16517       USE iophy
    16518       USE phys_state_var_mod, ONLY: q_ancien,xt_ancien,wake_deltaq,wake_deltaxt, &
    16519 #ifdef ISOVERIF
    16520         rain_fall,snow_fall,fevap,qsol, &
    16521 #endif
    16522         xtrain_fall,xtsnow_fall,ql_ancien,xtl_ancien,qs_ancien,xts_ancien, &
    16523         fxtevap,xtsol
    16524       !USE iostart
    16525       !USE write_field_phy
    16526       USE indice_sol_mod, only: nbsrf 
    16527   USE isotopes_mod, ONLY: isoName,iso_HDO,iso_eau
    16528 #ifdef ISOVERIF
    16529   USE isotopes_verif_mod
     16510SUBROUTINE phyiso_etat0_fichier(snow, run_off_lic_0, xtsnow, xtrun_off_lic_0, Rland_ice)
     16511   USE dimphy,             ONLY: klon,klev
     16512   USE iophy
     16513   USE phys_state_var_mod, ONLY: q_ancien, xt_ancien, wake_deltaq, wake_deltaxt, &
     16514#ifdef ISOVERIF
     16515     rain_fall, snow_fall, fevap,qsol, &
     16516#endif
     16517     xtrain_fall, xtsnow_fall, ql_ancien, xtl_ancien, qs_ancien, xts_ancien, fxtevap, xtsol
     16518   USE indice_sol_mod,    ONLY: nbsrf
     16519   USE isotopes_mod,      ONLY: isoName,iso_HDO,iso_eau
     16520   USE phyetat0_get_mod,  ONLY: phyetat0_get, phyetat0_srf
     16521   USE readTracFiles_mod, ONLY: new2oldH2O
     16522   USE strings_mod,       ONLY: strIdx, strTail, maxlen, msg, int2str
     16523#ifdef ISOVERIF
     16524   USE isotopes_verif_mod
    1653016525#endif
    1653116526#ifdef ISOTRAC
    16532  USE isotrac_mod, ONLY: strtrac,initialisation_isotrac,index_iso, &
    16533 &       index_zone,izone_init
    16534 #endif
    16535         implicit none
     16527   USE isotrac_mod, ONLY: strtrac, initialisation_isotrac, index_iso, index_zone, izone_init
     16528#endif
     16529   IMPLICIT NONE
    1653616530
    1653716531#include "netcdf.inc"
    1653816532#include "dimsoil.h"
    1653916533#include "clesphys.h"
    16540 ! #include "thermcell.h"
    1654116534#include "compbl.h"   
    1654216535
    16543         ! inputs
    16544         !REAL qsol(klon)
    16545         REAL snow(klon,nbsrf)
    16546         !REAL evap(klon,nbsrf)
    16547         REAL run_off_lic_0(klon)
    16548         ! outputs   
    16549         !REAL xtsol(niso,klon)
    16550         REAL xtsnow(niso,klon,nbsrf)
    16551         !REAL xtevap(ntraciso,klon,nbsrf)     
    16552         REAL xtrun_off_lic_0(niso,klon)
    16553         REAL Rland_ice(niso,klon)
    16554 
    16555         ! locals
    16556         real iso_tmp(klon)
    16557         real iso_tmp_lonlev(klon,klev)
    16558         real iso_tmp_lonsrf(klon,nbsrf)
    16559         INTEGER ierr
    16560         integer i,ixt,k,nsrf
    16561         INTEGER nid, nvarid
    16562         CHARACTER*2 str2
    16563         CHARACTER*5 str5
    16564         real xmin,xmax   
    16565         CHARACTER*50 outiso 
    16566         integer lnblnk
    16567         LOGICAL :: found,phyetat0_get,phyetat0_srf
    16568 
    16569 !#ifdef ISOVERIF
    16570 !      integer iso_verif_egalite_nostop
    16571 !#endif
    16572 !#ifdef ISOVERIF
    16573 !        real deltaD
    16574 !        integer iso_verif_noNaN_nostop
    16575 !#endif
     16536   REAL, INTENT(IN) ::             snow     (klon,nbsrf)
     16537   REAL, INTENT(IN) ::    run_off_lic_0     (klon)
     16538   REAL, INTENT(OUT) ::          xtsnow(niso,klon,nbsrf)
     16539   REAL, INTENT(OUT) :: xtrun_off_lic_0(niso,klon)
     16540   REAL, INTENT(OUT) ::       Rland_ice(niso,klon)
     16541
     16542   INTEGER :: ierr, i, ixt, k, nsrf, nid, nvarid, lnblnk
     16543   CHARACTER(LEN=2) :: str2
     16544   CHARACTER(LEN=5) :: str5
     16545   CHARACTER(LEN=maxlen) :: outiso, oldIso, modname, nam(2)
     16546   REAL :: xmin, xmax
     16547   LOGICAL :: found
    1657616548#ifdef ISOTRAC
    16577         integer iiso,izone
    16578 #endif
    16579 
    16580 
    16581    write(*,*) 'phyiso_etat0_fichier 3'
    16582    write(*,*) 'niso=',niso
    16583    write(*,*) 'isoName(1)='//TRIM(isoName(1))
    16584 
    16585    do ixt=1,ntraciso
    16586 
    16587       outiso=TRIM(isoName(ixt))
    16588       i = INDEX(outiso, '_', .TRUE.)
    16589       outiso = outiso(1:i-1)//outiso(i+1:LEN_TRIM(outiso))
    16590       write(*,*) 'phyiso_etat0_fichier 16621: ixt,outiso=',ixt,TRIM(outiso)
    16591 
    16592            
    16593       ! on lit seulement si ixt<=niso ou si on initialise les traceurs d'après
    16594       ! fichier:
     16549   INTEGER :: iiso, izone
     16550#endif
     16551
     16552   modname = 'phyiso_etat0_fichier'
     16553   CALL msg('3', modname)
     16554   CALL msg('niso = '//TRIM(int2str(niso)), modname)
     16555   CALL msg('isoName(1) = '//TRIM(isoName(1)), modname)
     16556
     16557   DO ixt = 1, ntraciso
     16558
     16559      outiso = isoName(ixt)
     16560      oldIso = strTail(new2oldH2O(outiso), '_')            !--- Remove "H2O_" from "H2O_<iso>[_<tag>]"
     16561      ! on lit seulement si ixt<=niso ou si on initialise les traceurs d'après fichier:
    1659516562#ifdef ISOTRAC
    16596       if ((ixt.le.niso).or.(initialisation_isotrac.eq.0)) then
    16597 #endif
    16598 
    16599       found=phyetat0_srf(1,iso_tmp_lonsrf,"XTSNOW"//TRIM(outiso),"Surface snow",0.)
    16600       if (.NOT.found) CALL abort_physic('isotopes_routines_mod', &
    16601                             'phyiso_etat0_fichier 16581: variable isotopique not found',1)
    16602       xtsnow(ixt,:,:)=iso_tmp_lonsrf(:,:)
    16603      
    16604       found=phyetat0_srf(1,iso_tmp_lonsrf,"XTEVAP"//TRIM(outiso),"evaporation",0.)
    16605       fxtevap(ixt,:,:)=iso_tmp_lonsrf(:,:)
    16606 
    16607       found=phyetat0_get(1,iso_tmp,"xtrain_f"//TRIM(outiso),"xrain fall",0.)
    16608       xtrain_fall(ixt,:)=iso_tmp(:)
    16609 
    16610       found=phyetat0_get(1,iso_tmp,"xtsnow_f"//TRIM(outiso),"snow fall",0.)
    16611       xtsnow_fall(ixt,:)=iso_tmp(:)
    16612 
    16613       found=phyetat0_get(klev,iso_tmp_lonlev,"XTANCIEN"//TRIM(outiso),"QANCIEN",0.)
    16614       xt_ancien(ixt,:,:)=iso_tmp_lonlev(:,:)
    16615 
    16616       found=phyetat0_get(klev,iso_tmp_lonlev,"XTLANCIEN"//TRIM(outiso),"QLANCIEN",0.)
    16617       xtl_ancien(ixt,:,:)=iso_tmp_lonlev(:,:)
    16618 
    16619       found=phyetat0_get(klev,iso_tmp_lonlev,"XTSANCIEN"//TRIM(outiso),"QSANCIEN",0.)
    16620       xts_ancien(ixt,:,:)=iso_tmp_lonlev(:,:)
    16621 
    16622       found=phyetat0_get(1,iso_tmp,"XTRUNOFFLIC0"//TRIM(outiso),"RUNOFFLIC0",0.) 
    16623       xtrun_off_lic_0(ixt,:)=iso_tmp(:)
    16624 
    16625       found=phyetat0_get(klev,iso_tmp_lonlev,"WAKE_DELTAXT"//TRIM(outiso),"Delta hum. wake/env",0.) 
    16626       wake_deltaxt(ixt,:,:)=iso_tmp_lonlev(:,:)
    16627 
    16628 #ifdef ISOVERIF           
    16629       if ((ixt.eq.iso_eau).and.(iso_eau.gt.0)) then
    16630         do i=1,klon
    16631          call iso_verif_egalite(xtrain_fall(iso_eau,i),rain_fall(i), &
    16632      &           'phyisoetat0_fichier 231a')
    16633          call iso_verif_egalite(xtsnow_fall(iso_eau,i),snow_fall(i), &
    16634      &           'phyisoetat0_fichier 231b')
    16635          DO nsrf = 1, nbsrf
    16636          call iso_verif_egalite(fxtevap(iso_eau,i,nsrf),fevap(i,nsrf), &
    16637      &           'phyisoetat0_fichier 231c')
    16638          call iso_verif_egalite(xtsnow(iso_eau,i,nsrf),snow(i,nsrf), &
    16639      &           'phyisoetat0_fichier 231d')
    16640          enddo !DO nsrf = 1, nbsrf
    16641         enddo !do i=1,klon       
    16642       endif !if (iso_eau.gt.0) then 
    16643         if ((iso_HDO.gt.0).and.(ixt.eq.iso_HDO)) then
    16644               do k=1,klev
    16645                do i=1,klon
    16646                 if (q_ancien(i,k).gt.2e-3) then
    16647                 call iso_verif_aberrant(xt_ancien(iso_hdo,i,k) &
    16648      &           /q_ancien(i,k),'phyisoetat0_fichier 312')
    16649                 endif !if (q_ancien(i,k).gt.2e-3) then
    16650                enddo !do i=1,klon
    16651               enddo !do k=1,klev
    16652       endif !if ((iso_HDO.gt.0).and.(ixt.eq.iso_HDO)) then
    16653       if (iso_eau.gt.0) then   
    16654         do i=1,klon
    16655           if (iso_verif_egalite_nostop(run_off_lic_0(i), &
    16656      &           xtrun_off_lic_0(iso_eau,i), &
    16657      &          'phyiso_etat0_fichier 326').eq.1) then
    16658             write(*,*) 'i=',i
    16659             stop
    16660           endif !if (iso_verif_egalite_nostop(run_off_lic_0(i),
    16661         enddo !do i=1,klon
    16662       endif !if (iso_eau.gt.0) then
    16663 #endif
    16664 
    16665        ! ces variables n'ont pas de traceurs:
    16666        if (ixt.le.niso) then
    16667         found=phyetat0_get(1,iso_tmp,"XTSOL"//TRIM(outiso),"Surface hmidity / bucket",0.) 
    16668         xtsol(ixt,:)=iso_tmp(:)
    16669 
    16670         found=phyetat0_get(1,iso_tmp,"Rland_ice"//TRIM(outiso),"R land ice",0.)
    16671         Rland_ice(ixt,:)=iso_tmp(:)
    16672 
    16673 #ifdef ISOVERIF
    16674       do i=1,klon
    16675           if (iso_verif_noNaN_nostop(xtsol(ixt,i), &
    16676      &          'phyiso_etat0_fichier 95').eq.1) then
    16677             write(*,*) 'ixt,i=',ixt,i
    16678             stop
    16679           endif       
    16680       enddo !do i=1,klon
    16681 #endif
    16682 
    16683        endif
    16684 
     16563      IF(ixt <= niso .OR. initialisation_isotrac == 0) THEN
     16564#endif
     16565      found = phyetat0iso_srf3(xtsnow,      "XTSNOW", "Surface snow", 0.)
     16566      if (.NOT.found) CALL abort_physic('isotopes_routines_mod', 'phyiso_etat0_fichier 16581: unfound isotopic variable',1)
     16567      found = phyetat0iso_srf3(fxtevap,     "XTEVAP", "evaporation",  0.)
     16568      found = phyetat0iso_get2(xtrain_fall, "xtrain_f", "xrain fall", 0.)
     16569      found = phyetat0iso_get2(xtsnow_fall, "xtsnow_f", "xsnow fall", 0.)
     16570      found = phyetat0iso_get3(xt_ancien,   "XTANCIEN",  "QANCIEN",   0.)
     16571      found = phyetat0iso_get3(xtl_ancien,  "XTLANCIEN", "QLANCIEN",  0.)
     16572      found = phyetat0iso_get3(xts_ancien,  "XTSANCIEN", "QSANCIEN",  0.)
     16573      found = phyetat0iso_get2(xtrun_off_lic_0, "XTRUNOFFLIC0", "RUNOFFLIC0", 0.)
     16574      found = phyetat0iso_get3(wake_deltaxt,  "WAKE_DELTAXT", "Delta hum. wake/env",  0.)
     16575#ifdef ISOVERIF
     16576      IF(ixt == iso_eau .AND. iso_eau > 0) THEN
     16577         DO i=1,klon
     16578            CALL iso_verif_egalite(xtrain_fall(iso_eau,i),rain_fall(i),TRIM(modname)//' 231a')
     16579            CALL iso_verif_egalite(xtsnow_fall(iso_eau,i),snow_fall(i),TRIM(modname)//' 231b')
     16580            DO nsrf = 1, nbsrf
     16581               CALL iso_verif_egalite(fxtevap(iso_eau,i,nsrf),fevap(i,nsrf),TRIM(modname)//' 231c')
     16582               CALL iso_verif_egalite( xtsnow(iso_eau,i,nsrf), snow(i,nsrf),TRIM(modname)//' 231d')
     16583            END DO
     16584         END DO
     16585      END IF
     16586      IF(ixt == iso_HDO .AND. iso_HDO > 0) THEN
     16587         DO k=1,klev
     16588            DO i=1,klon
     16589               IF(q_ancien(i,k) > 2e-3) &
     16590                  CALL iso_verif_aberrant(xt_ancien(iso_hdo,i,k)/q_ancien(i,k),TRIM(modname)//' 312')
     16591            END DO
     16592         END DO
     16593      END IF
     16594      IF(iso_eau > 0 .AND. ixt == iso_eau) THEN
     16595         DO i=1,klon
     16596            IF(iso_verif_egalite_nostop(run_off_lic_0(i),xtrun_off_lic_0(iso_eau,i),TRIM(modname)//' 326') == 1) THEN
     16597               WRITE(*,*) 'i=',i
     16598               STOP
     16599            END IF
     16600         END DO
     16601      END IF
     16602#endif
     16603      ! ces variables n'ont pas de traceurs:
     16604      IF(ixt <= niso) THEN
     16605         found = phyetat0iso_get2(xtsol, "XTSOL", "Surface humidity / bucket", 0.)
     16606         found = phyetat0iso_get2(Rland_ice, "Rland_ice", "SR land ice", 0.)
     16607#ifdef ISOVERIF
     16608
     16609         DO i=1,klon
     16610            IF(iso_verif_noNaN_nostop(xtsol(ixt,i),TRIM(modname)//' 95') == 1) THEN
     16611               WRITE(*,*) 'ixt,i=',ixt,i
     16612               STOP
     16613            END IF
     16614         END DO
     16615#endif
     16616      END IF
    1668516617#ifdef ISOTRAC
    16686      endif !if ((ixt.le.niso).or.(initialisation_isotrac.eq.0)) then
    16687 #endif
    16688 
    16689   enddo !do ixt=1,ntraciso
     16618      END IF ! IF(ixt > niso .OR. initialisation_isotrac == 0))
     16619#endif
     16620
     16621   END DO
    1669016622
    1669116623#ifdef ISOTRAC
    16692         if (initialisation_isotrac.ne.0) then
    16693         ! on n'initialise pas d'après le fichier
    16694         ! l'eau normale est mise dans la zone izone_init
    16695 
    16696         do ixt=niso+1,ntraciso
    16697 
    16698              iiso=index_iso(ixt)
    16699 
    16700              if (index_zone(ixt).eq.izone_init) then
    16701                 do i=1,klon
    16702                  do nsrf = 1, nbsrf
    16703                   fxtevap(ixt,i,nsrf)=fxtevap(iiso,i,nsrf)
    16704                  enddo !do nsrf = 1, nbsrf
    16705                  xtsnow_fall(ixt,i)=xtsnow_fall(iiso,i)
    16706                  xtrain_fall(ixt,i)=xtrain_fall(iiso,i)
    16707                  do k=1,klev
    16708                    xt_ancien(ixt,i,k)=xt_ancien(iiso,i,k)
    16709                    xtl_ancien(ixt,i,k)=xtl_ancien(iiso,i,k)
    16710                    xts_ancien(ixt,i,k)=xts_ancien(iiso,i,k)
    16711                    wake_deltaxt(ixt,i,k)= wake_deltaxt(iiso,i,k)   
    16712                  enddo
    16713                 enddo !do i=1,klon
    16714              else !if (index_zone(ixt).eq.izone_init) then
    16715                 do i=1,klon
    16716                  do nsrf = 1, nbsrf
    16717                   fxtevap(ixt,i,nsrf)=0.0
    16718                  enddo !do nsrf = 1, nbsrf
    16719                  xtsnow_fall(ixt,i)=0.0
    16720                  xtrain_fall(ixt,i)=0.0
    16721                  do k=1,klev
    16722                    xt_ancien(ixt,i,k)=0.0
    16723                    xtl_ancien(ixt,i,k)=0.0
    16724                    xts_ancien(ixt,i,k)=0.0
    16725                  enddo
    16726                 enddo !do i=1,klon
    16727              endif !if (index_zone(ixt).eq.izone_init) then
    16728 
    16729          enddo  !do ixt=1,niso
    16730       endif !if (initialisation_isotrac.eq.0) then
    16731 
    16732 
    16733 #ifdef ISOVERIF
    16734         DO nsrf = 1, nbsrf
    16735          do i=1,klon
    16736                call iso_verif_traceur(fxtevap(1,i,nsrf), &
    16737      &                   'phyiso_etat0_fichier 426')
    16738          enddo !do i=1,klon
    16739         enddo !DO nsrf = 1, nbsrf
    16740         do i=1,klon
    16741            call iso_verif_traceur(xtrain_fall(1,i), &
    16742      &                   'phyiso_etat0_fichier 466')
    16743            call iso_verif_traceur(xtsnow_fall(1,i), &
    16744      &                   'phyiso_etat0_fichier 468')
    16745         enddo !do i=1,klon
    16746         do k=1,klev
    16747           do i=1,klon
    16748                call iso_verif_traceur(xt_ancien(1,i,k), &
    16749      &                   'phyiso_etat0_fichier 591')
    16750           enddo !do i=1,klon
    16751         enddo !do k=1,klev             
     16624   IF(initialisation_isotrac /= 0) THEN
     16625      ! On n'initialise pas d'apres le fichier. L'eau normale est mise dans la zone izone_init
     16626      DO ixt=niso+1,ntraciso
     16627         iiso=index_iso(ixt)
     16628         IF(index_zone(ixt) == izone_init) THEN
     16629            DO i = 1, klon
     16630               fxtevap(ixt,i,1:nbsrf) = fxtevap(iiso,i,1:nbsrf)
     16631               xtsnow_fall(ixt,i) = xtsnow_fall(iiso,i)
     16632               xtrain_fall(ixt,i) = xtrain_fall(iiso,i)
     16633               DO k = 1, klev
     16634                  xt_ancien   (ixt,i,k) = xt_ancien   (iiso,i,k)
     16635                  xtl_ancien  (ixt,i,k) = xtl_ancien  (iiso,i,k)
     16636                  xts_ancien  (ixt,i,k) = xts_ancien  (iiso,i,k)
     16637                  wake_deltaxt(ixt,i,k) = wake_deltaxt(iiso,i,k)   
     16638               END DO
     16639            END DO
     16640         ELSE
     16641            DO i = 1, klon
     16642               fxtevap(ixt,i,1:nbsrf)=0.0
     16643               xtsnow_fall(ixt,i)=0.0
     16644               xtrain_fall(ixt,i)=0.0
     16645               xt_ancien (ixt,i,1:klev) = 0.0
     16646               xtl_ancien(ixt,i,1:klev) = 0.0
     16647               xts_ancien(ixt,i,1:klev) = 0.0
     16648            END DO
     16649         END IF
     16650      END DO
     16651   END IF
     16652
     16653#ifdef ISOVERIF
     16654   DO nsrf = 1, nbsrf
     16655      DO i = 1, klon
     16656         CALL iso_verif_traceur(fxtevap(1,i,nsrf), 'phyiso_etat0_fichier 426')
     16657      END DO
     16658   END DO
     16659   DO i=1,klon
     16660      CALL iso_verif_traceur(xtrain_fall(1,i), 'phyiso_etat0_fichier 466')
     16661      CALL iso_verif_traceur(xtsnow_fall(1,i), 'phyiso_etat0_fichier 468')
     16662   END DO
     16663   DO k = 1, klev
     16664      DO i = 1, klon
     16665         CALL iso_verif_traceur(xt_ancien(1,i,k), 'phyiso_etat0_fichier 591')
     16666      END DO
     16667   END DO
    1675216668#endif
    1675316669        ! endif ISOVERIF       
     
    1675516671        ! endif ISOTRAC     
    1675616672
    16757 ! on ferme le fichier
    16758 !      CALL close_startphy
    16759 ! déjà fermé dans phyetat0
     16673CONTAINS
     16674
     16675LOGICAL FUNCTION phyetat0iso_get2(field, pref, descr, default) RESULT(lFound)
     16676  REAL,             INTENT(INOUT) :: field(:,:)
     16677  CHARACTER(LEN=*), INTENT(IN)    :: pref, descr
     16678  REAL,             INTENT(IN)    :: default
     16679  REAL :: iso_tmp(klon)
     16680  nam(1) = TRIM(pref)//TRIM(outiso)
     16681  nam(2) = TRIM(pref)//TRIM(oldIso)
     16682  lFound = phyetat0_get(iso_tmp, nam, descr, default)
     16683  field(ixt,:) = iso_tmp
     16684END FUNCTION phyetat0iso_get2
     16685
     16686
     16687LOGICAL FUNCTION phyetat0iso_get3(field, pref, descr, default) RESULT(lFound)
     16688  REAL,             INTENT(INOUT) :: field(:,:,:)
     16689  CHARACTER(LEN=*), INTENT(IN)    :: pref, descr
     16690  REAL,             INTENT(IN)    :: default
     16691  REAL :: iso_tmp_lonlev(klon,klev)
     16692  nam(1) = TRIM(pref)//TRIM(outiso)
     16693  nam(2) = TRIM(pref)//TRIM(oldIso)
     16694  lFound = phyetat0_get(iso_tmp_lonlev, nam, descr, default)
     16695  field(ixt,:,:) = iso_tmp_lonlev(:,:)
     16696END FUNCTION phyetat0iso_get3
     16697
     16698LOGICAL FUNCTION phyetat0iso_srf3(field, pref, descr, default) RESULT(lFound)
     16699  REAL,             INTENT(INOUT) :: field(:,:,:)
     16700  CHARACTER(LEN=*), INTENT(IN)    :: pref, descr
     16701  REAL,             INTENT(IN)    :: default
     16702  REAL :: iso_tmp_lonsrf(klon,nbsrf)
     16703  nam(1) = TRIM(pref)//TRIM(outiso)
     16704  nam(2) = TRIM(pref)//TRIM(oldIso)
     16705  lFound = phyetat0_srf(iso_tmp_lonsrf, nam, descr, default)
     16706  field(ixt,:,:) = iso_tmp_lonsrf
     16707END FUNCTION phyetat0iso_srf3
    1676016708
    1676116709        end subroutine phyiso_etat0_fichier
     16710
     16711
    1676216712
    1676316713
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/isotopes_verif_mod.F90

    r4149 r4482  
    55MODULE isotopes_verif_mod
    66!use isotopes_mod, ONLY:
    7 #ifdef ISOTRAC
    8    USE isotrac_mod, ONLY: nzone
    9 #endif
    10 USE infotrac_phy, ONLY: ntraciso=>ntiso, niso, itZonIso
     7!#ifdef ISOTRAC
     8!   USE isotrac_mod, ONLY: nzone
     9!#endif
     10USE infotrac_phy, ONLY: ntraciso=>ntiso, niso, itZonIso, nzone
    1111implicit none
    1212save
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/isotrac_mod.F90

    r4143 r4482  
    11#ifdef ISO
    22#ifdef ISOTRAC
    3 ! $Id: $
    43
    54MODULE isotrac_mod
    6 use infotrac_phy, ONLY: niso,ntiso,ntraceurs_zone=>nzone
    7 use isotopes_mod, only: ridicule
    8 
    9 IMPLICIT NONE
    10 SAVE
    11 
    12 ! contient toutes les variables traceurs isotopiques + les routines specifiquement
    13 ! traceurs isotopiques
    14 
    15       real ridicule_trac
    16       parameter (ridicule_trac=ridicule*1e4)
    17 
    18 integer, save ::  option_traceurs
    19 integer, save ::  ntraceurs_zone_opt ! ntraceurs_zone propre à l'option
    20 ! on vérifie que ça correspond bien à ntraceurs_zone d'infotrac
    21 integer, save ::  ntraceurs_zoneOR
    22 !$OMP THREADPRIVATE(option_traceurs,ntraceurs_zone_opt,ntraceurs_zoneOR)
    23 integer, save ::  initialisation_isotrac
    24                 ! 1 pour idéalisé
    25                 ! 0 pour lecture dans fichier
    26 !$OMP THREADPRIVATE(initialisation_isotrac)
    27 
    28 ! variables spécifiques aux différentes options, mais necessaires au
    29 ! calcul du nombre de zones de traceurs
    30 ! si option=3
    31 integer, save :: use_bassin_atlantic
    32 !$OMP THREADPRIVATE(use_bassin_atlantic)
    33 integer, save :: use_bassin_medit
    34 !$OMP THREADPRIVATE(use_bassin_medit)
    35 integer, save :: use_bassin_indian
    36 !$OMP THREADPRIVATE(use_bassin_indian)
    37 integer, save :: use_bassin_austral
    38 !$OMP THREADPRIVATE(use_bassin_austral)
    39 integer, save :: use_bassin_pacific
    40 !$OMP THREADPRIVATE(use_bassin_pacific)
    41 integer, save :: use_bassin_merarabie
    42 !$OMP THREADPRIVATE(use_bassin_merarabie)
    43 integer, save :: use_bassin_golfebengale
    44 !$OMP THREADPRIVATE(use_bassin_golfebengale)
    45 integer, save :: use_bassin_indiansud
    46 !$OMP THREADPRIVATE(use_bassin_indiansud)
    47 integer, save :: use_bassin_tropics
    48 !$OMP THREADPRIVATE(use_bassin_tropics)
    49 integer, save :: use_bassin_midlats
    50 !$OMP THREADPRIVATE(use_bassin_midlats)
    51 integer, save :: use_bassin_hauteslats
    52 !$OMP THREADPRIVATE(use_bassin_hauteslats)
    53 integer, save :: bassin_atlantic
    54 !$OMP THREADPRIVATE(bassin_atlantic)
    55 integer, save :: bassin_medit
    56 !$OMP THREADPRIVATE(bassin_medit)
    57 integer, save :: bassin_indian
    58 !$OMP THREADPRIVATE(bassin_indian)
    59 integer, save :: bassin_austral
    60 !$OMP THREADPRIVATE(bassin_austral)
    61 integer, save :: bassin_pacific
    62 !$OMP THREADPRIVATE(bassin_pacific)
    63 integer, save :: bassin_merarabie
    64 !$OMP THREADPRIVATE(bassin_merarabie)
    65 integer, save :: bassin_golfebengale
    66 !$OMP THREADPRIVATE(bassin_golfebengale)
    67 integer, save :: bassin_indiansud
    68 !$OMP THREADPRIVATE(bassin_indiansud)
    69 integer, save :: bassin_tropics
    70 !$OMP THREADPRIVATE(bassin_tropics)
    71 integer, save :: bassin_midlats
    72 !$OMP THREADPRIVATE(bassin_midlats)
    73 integer, save :: bassin_hauteslats
    74 !$OMP THREADPRIVATE(bassin_hauteslats)
    75 ! si option=4
    76 integer nzone_temp
    77 parameter (nzone_temp=1)
    78 real, save :: zone_temp1,zone_tempf,zone_tempa 
    79 !$OMP THREADPRIVATE(zone_temp1,zone_tempf,zone_tempa)
    80 ! si option 14
    81 integer nzone_lat
    82 parameter (nzone_lat=4)
    83 integer nzone_pres
    84 parameter (nzone_pres=3)
    85 real, save :: zone_pres1,zone_presf,zone_presa
    86 !$OMP THREADPRIVATE(zone_pres1,zone_presf,zone_presa)
    87 real, save :: dlattag,lattag_min
    88 !$OMP THREADPRIVATE(dlattag,lattag_min)
    89 
    90 
    91 ! option 1: on trace evap ocean et continent séparement 
     5  USE infotrac_phy,      ONLY: niso, ntiso, nzone
     6  USE readTracFiles_mod, ONLY: delPhase
     7  USE isotopes_mod,      ONLY: ridicule, get_in
     8
     9  IMPLICIT NONE
     10  SAVE
     11
     12!=== CONTENT: ALL THE ISOTOPIC TRACERS RELATED VARIABLES ===
     13!
     14! option 1: on trace evap ocean et continent separement 
    9215! option 2: on trace evap ocean, continent et evap precip
    93 ! option 3: on trace evap différents bassins océaniques
    94 !       + continents + résidu
    95 !       attention, choisir dans ce cas les bassins océaniques
     16! option 3: on trace evap differents bassins oceaniques
     17!       + continents + residu
     18!       attention, choisir dans ce cas les bassins oceaniques
    9619!       dans iso_traceurs_opt3F90.h
    97 ! option 4: tracage par température minimale
    98 !       dans ce cas, on définit des bins dans iso_traceurs_opt4.h
    99 ! option 5: pour AMMA: on taggue résidu/AEJ/flux mousson/Harmattan
     20! option 4: tracage par temperature minimale
     21!       dans ce cas, on definit des bins dans iso_traceurs_opt4.h
     22! option 5: pour AMMA: on taggue residu/AEJ/flux mousson/Harmattan
    10023! option 6: taggage des ddfts
    101 ! option 7: pour Sandrine: taggage de la vapeur à 700hPa pour omega500<-20 TODO
    102 ! option 8: pour Sandrine: taggage de la vapeur entre 950 et 800hPa, omega de 0 à 25 hPa et de l'évaoration en omega<-20. TODO
     24! option 7: pour Sandrine: taggage de la vapeur a 700hPa pour omega500<-20 TODO
     25! option 8: pour Sandrine: taggage de la vapeur entre 950 et 800hPa, omega de 0 a 25 hPa et de l'evaoration en omega<-20. TODO
    10326! option 9: taggage du condensat et de la revap precip
    10427! option 10: taggage evap oce, transpiration et evaporation
     
    10730! option 12: taggage evap oce, sol nu, canop et reste evap cont.
    10831! A utiliser quand on couple avec ORCHIDEE
    109 ! option 13: taggage température minimale + revap precip
    110 ! option 14: taggage lat et altitude de dernière saturation (niveaux de pression) + evap surf
     32! option 13: taggage temperature minimale + revap precip
     33! option 14: taggage lat et altitude de derniere saturation (niveaux de pression) + evap surf
    11134! otion 15: taggage irrigation
    11235! option 16: taggage precip selon saisons et fonte neige: seulement pour ORCHIDEE
    113 ! option 17: taggage température minimum de condensation directement dans la convection et la cond LS, + evap sfc, condensat et precipitation
     36! option 17: taggage temperature minimum de condensation directement dans la convection et la cond LS, + evap sfc, condensat et precipitation
    11437! option 18: idem 17 mais on tague qsmin au lieu de Tmin
    11538! option 19: on tag vap residuelle, vap residuelle dans ddfts, sfc, cond, rev
    11639! option 20: on taggue vapeur tropicale vs vapeur extratropicale
    11740! option 21: taggage de 2 boites 3D: extratropiques (>35°) et UT tropicale (15-15°, > 500hPa)
    118 ! option 22: tagage de la vapeur proccessée dans les zones très convectives
     41! option 22: tagage de la vapeur proccessee dans les zones tres convectives
    11942               
    120         ! ces variables sont initialisées dans traceurs_init
     43   !--- nzone_opt (value of nzone for the selected option) must be equal to nzone as defined in onfotrac
     44   REAL, PARAMETER :: ridicule_trac = ridicule * 1e4
     45   INTEGER, SAVE :: option_traceurs, nzone_opt, nzoneOR
     46!$OMP THREADPRIVATE(option_traceurs,nzone_opt,nzoneOR)
     47   INTEGER, SAVE :: initialisation_isotrac
     48!$OMP THREADPRIVATE(initialisation_isotrac)     
     49                ! 1 pour idealise
     50                ! 0 pour lecture dans fichier
     51
     52   !=== VARIABLES SPECIFIC TO THE SELECTED OPTION, BUT NEEDED FOR THE COMPUTATION OF THE NUMBER OF ZONES ; TO BE INITIALIZED IN traceurs_init
     53
     54   !--- option 3
     55   LOGICAL, SAVE :: use_bassin_Austral, use_bassin_Atlantic, use_bassin_MidLats, use_bassin_SouthIndian, use_bassin_MerArabie
     56!$OMP THREADPRIVATE(use_bassin_Austral, use_bassin_Atlantic, use_bassin_MidLats, use_bassin_SouthIndian, use_bassin_MerArabie)
     57   INTEGER, SAVE ::     bassin_Austral,     bassin_Atlantic,     bassin_MidLats,     bassin_SouthIndian,     bassin_MerArabie
     58!$OMP THREADPRIVATE(    bassin_Austral,     bassin_Atlantic,     bassin_MidLats,     bassin_SouthIndian,     bassin_MerArabie)
     59   LOGICAL, SAVE :: use_bassin_Pacific, use_bassin_Indian,   use_bassin_Tropics, use_bassin_BengalGolf,  use_bassin_HighLats, use_bassin_Medit
     60!$OMP THREADPRIVATE(use_bassin_Pacific, use_bassin_Indian,   use_bassin_Tropics, use_bassin_BengalGolf,  use_bassin_HighLats, use_bassin_Medit)
     61   INTEGER, SAVE ::     bassin_Pacific,     bassin_Indian,       bassin_Tropics,     bassin_BengalGolf,      bassin_HighLats,     bassin_Medit
     62!$OMP THREADPRIVATE(    bassin_Pacific,     bassin_Indian,       bassin_Tropics,     bassin_BengalGolf,      bassin_HighLats,     bassin_Medit)
     63
     64   !--- option 4
     65   INTEGER, PARAMETER :: nzone_temp = 1
     66   REAL,   SAVE ::  zone_temp1, zone_tempf, zone_tempa 
     67!$OMP THREADPRIVATE(zone_temp1, zone_tempf, zone_tempa)
     68   REAL,   SAVE ::  zone_temp(nzone_temp-1)
     69!$OMP THREADPRIVATE(zone_temp)
     70
     71   !--- option 5
     72   INTEGER, SAVE :: izone_aej, izone_harmattan, izone_mousson
     73!$OMP THREADPRIVATE(izone_aej, izone_harmattan, izone_mousson)
     74
     75   !--- option 6
     76   INTEGER, SAVE :: izone_ddft
     77!$OMP THREADPRIVATE(izone_ddft)
     78
     79   !--- option 10
     80   INTEGER, SAVE :: izone_contfrac
     81!$OMP THREADPRIVATE(izone_contfrac)
     82
     83   !--- option 12       
     84   INTEGER, SAVE :: izone_contcanop
     85!$OMP THREADPRIVATE(izone_contcanop)
     86
     87   !--- option 13
     88   INTEGER, PARAMETER :: nzone_pres = 3
     89   REAL, SAVE ::  zone_pres(nzone_pres-1)
     90!$OMP THREADPRIVATE(zone_pres)
     91
     92   !--- option 14
     93   INTEGER, PARAMETER :: nzone_lat = 4
     94   REAL,    SAVE :: zone_pres1, zone_presf, zone_presa
     95!$OMP THREADPRIVATE(zone_pres1, zone_presf, zone_presa)
     96   REAL,    SAVE :: dlattag, lattag_min, zone_lat(nzone_lat-1)
     97!$OMP THREADPRIVATE(dlattag, lattag_min, zone_lat)
     98
     99   !--- option 15
     100   INTEGER, SAVE :: izone_irrig
     101!$OMP THREADPRIVATE(izone_irrig)
     102
     103   !--- option 17
     104   REAL,    SAVE :: seuil_tag_tmin, seuil_tag_tmin_ls
     105!$OMP THREADPRIVATE(seuil_tag_tmin, seuil_tag_tmin_ls)
     106  INTEGER,  SAVE :: option_seuil_tag_tmin
     107!$OMP THREADPRIVATE(option_seuil_tag_tmin)
     108
     109   !--- option 20
     110   INTEGER, SAVE :: izone_trop, izone_extra
     111!$OMP THREADPRIVATE(izone_trop, izone_extra)
     112   REAL,    SAVE :: lim_tag20
     113!$OMP THREADPRIVATE(lim_tag20)
     114
     115   !--- option 21: on garde izone_trop, izone_extra 
     116
     117   !--- option 22
     118   INTEGER, SAVE :: izone_conv_BT, izone_conv_UT
     119!$OMP THREADPRIVATE(izone_conv_BT, izone_conv_UT)
     120   REAL,    SAVE :: lim_precip_tag22
     121!$OMP THREADPRIVATE(lim_precip_tag22)
     122
    121123       
    122 integer, ALLOCATABLE, DIMENSION(:), save :: index_iso
    123 !$OMP THREADPRIVATE(index_iso)
    124 integer, ALLOCATABLE, DIMENSION(:), save ::  index_zone
    125 !$OMP THREADPRIVATE(index_zone)
    126 integer, ALLOCATABLE, DIMENSION(:,:), save ::  itZonIso_loc ! il y a déjà un itZonIso dans infotrac: vérifier que c'est le même
    127 !$OMP THREADPRIVATE(itZonIso_loc)
    128 character*3, ALLOCATABLE, DIMENSION(:), save :: strtrac
    129 !$OMP THREADPRIVATE(strtrac)
    130 ! -> tout ça passe maintenant par infotrac
    131 
    132 integer, ALLOCATABLE, DIMENSION(:), save :: bassin_map
    133 integer, ALLOCATABLE, DIMENSION(:,:), save :: boite_map
    134 !$OMP THREADPRIVATE(bassin_map,boite_map)
    135 
    136 
    137         ! traitement recyclage et evap
    138 integer, save :: izone_cont ! pour le recyclage continental
    139 !$OMP THREADPRIVATE(izone_cont)
    140 integer, save :: izone_oce ! pour l'océan
    141 !$OMP THREADPRIVATE(izone_oce)
    142 integer, save :: izone_poubelle ! pour les petits résidus numériques
     124  INTEGER, ALLOCATABLE, SAVE :: index_iso(:), index_zone(:), itZonIso_loc(:,:)
     125!$OMP THREADPRIVATE(            index_iso,    index_zone,    itZonIso_loc)
     126  CHARACTER(LEN=3), ALLOCATABLE :: strtrac(:)
     127!$OMP THREADPRIVATE(               strtrac)
     128  INTEGER, ALLOCATABLE, SAVE :: bassin_map(:), boite_map(:,:)
     129!$OMP THREADPRIVATE(            bassin_map,    boite_map)
     130
     131   !=== RECYCLING AND EVAPORATION TREATMENT
     132   INTEGER, SAVE :: izone_cont, izone_oce        !--- For land and ocean recycling
     133!$OMP THREADPRIVATE(izone_cont, izone_oce)
     134   INTEGER, SAVE :: izone_poubelle               !--- For small numerical residues
    143135!$OMP THREADPRIVATE(izone_poubelle)
    144 integer, save :: izone_init ! pour l'initialisation par défaut
     136   INTEGER, SAVE :: izone_init                   !--- For default initialization
    145137!$OMP THREADPRIVATE(izone_init)
    146 integer, save :: izone_revap ! pour l'évap des gouttes
     138   INTEGER, SAVE :: izone_revap                  !--- For droplets evaporation
    147139!$OMP THREADPRIVATE(izone_revap)
    148 integer, save :: option_revap
    149 !$OMP THREADPRIVATE(option_revap)
    150 integer, save :: option_tmin
    151 !$OMP THREADPRIVATE(option_tmin)
    152 integer, save :: option_cond
    153 !$OMP THREADPRIVATE(option_cond)
    154 integer, save :: izone_cond
    155 !$OMP THREADPRIVATE(izone_cond)
    156 real evap_franche
    157 parameter (evap_franche=1e-6) ! en kg/m2/s
    158 
    159 ! specifique à option 4:
    160 real, save ::  zone_temp(nzone_temp-1)
    161 !$OMP THREADPRIVATE(zone_temp)
    162 ! si option 5
    163 integer, save :: izone_aej
    164 !$OMP THREADPRIVATE(izone_aej)
    165 integer, save :: izone_harmattan
    166 !$OMP THREADPRIVATE(izone_harmattan)
    167 integer, save :: izone_mousson
    168 !$OMP THREADPRIVATE(izone_mousson)
    169 ! si option 6
    170 integer, save :: izone_ddft
    171 !$OMP THREADPRIVATE(izone_ddft)
    172 ! si option 10
    173 integer, save :: izone_contfrac
    174 !$OMP THREADPRIVATE(izone_contfrac)
    175 ! si option 12 
    176 integer, save :: izone_contcanop
    177 !$OMP THREADPRIVATE(izone_contcanop)
    178 ! specifique à option 13:
    179 real, save ::  zone_pres(nzone_pres-1)
    180 !$OMP THREADPRIVATE(zone_pres)
    181 ! si option 14
    182 real, save ::  zone_lat(nzone_lat-1)
    183 !$OMP THREADPRIVATE(zone_lat)
    184 ! si option 15
    185 integer, save :: izone_irrig
    186 !$OMP THREADPRIVATE(izone_irrig)
    187 ! si option 17
    188 real, save ::  seuil_tag_tmin
    189 !$OMP THREADPRIVATE(seuil_tag_tmin)
    190 real, save ::  seuil_tag_tmin_ls
    191 !$OMP THREADPRIVATE(seuil_tag_tmin_ls)
    192 integer, save :: option_seuil_tag_tmin
    193 !$OMP THREADPRIVATE(option_seuil_tag_tmin)
    194 ! si option 20
    195 integer, save :: izone_trop,izone_extra
    196 real, save ::  lim_tag20
    197 !$OMP THREADPRIVATE(izone_trop,izone_extra,lim_tag20)
    198 ! si option 21: on garde izone_trop,izone_extra 
    199 ! si opt 22
    200 integer, save :: izone_conv_BT,izone_conv_UT
    201 real, save ::  lim_precip_tag22
    202 !$OMP THREADPRIVATE(izone_conv_BT,izone_conv_UT,lim_precip_tag22)
    203 
     140   INTEGER, SAVE :: option_revap, option_tmin, option_cond, izone_cond
     141!$OMP THREADPRIVATE(option_revap, option_tmin, option_cond, izone_cond)
     142   REAL, PARAMETER :: evap_franche = 1e-6        !--- In kg/m2/s
    204143
    205144CONTAINS
    206145
    207       subroutine iso_traceurs_init()
    208 
    209       use IOIPSL ! getin
    210       USE infotrac_phy, ONLY: itZonIso
    211       USE isotopes_mod, ONLY: iso_eau,ntracisoOR,initialisation_iso
    212       USE dimphy, only: klon,klev
    213 
    214         implicit none
    215 
    216 
    217         ! définition de quelles zones et quelles isotopes représentent
    218         ! les traceurs
    219 
    220         ! inputs, outputs
    221         ! ! c'est les variables dans traceurs.h qui sont modifiées
    222 
    223         ! locals
    224         integer itrac,izone,ixt,k
    225         integer izone_pres,izone_lat
    226         character*2 strz,strz_preslat
    227         character*1 strz_pres,strz_lat
    228         integer ntraceurs_zone_opt
    229 
    230         ! vérifier que on a bien l'eau comme traceurs
    231         if (iso_eau.eq.0) then
    232             write(*,*) 'traceurs_init 18: isotrac ne marche que si ', &
    233      &            'on met l''eau comme isotope'
    234             stop
    235         endif
    236 
    237         ! initialiser
    238         option_traceurs=0
    239         initialisation_isotrac=0
    240 
    241         ! allouer
    242         allocate (index_iso(ntiso))
    243         allocate (index_zone(ntiso))
    244         allocate (itZonIso_loc(ntraceurs_zone,niso))
    245         allocate (strtrac(ntraceurs_zone))
    246         allocate (bassin_map(klon))
    247         allocate (boite_map(klon,klev))
    248 
    249         if (initialisation_iso.eq.0) then
    250           call getin('initialisation_isotrac',initialisation_isotrac)
    251           write(*,*) 'initialisation_isotrac=',initialisation_isotrac
    252         endif !if (initialisation_iso.eq.0) then
    253 
    254         ! lire l'option de traçage
    255         call getin('option_traceurs',option_traceurs)
    256         write(*,*) 'option_traceurs=',option_traceurs
    257 
    258         ! cas général: pas de traceurs dans ORCHIDEE
    259         ntracisoOR=niso
    260 
    261         ! partie à éditer ! pour définir les différentes zones
    262         if (option_traceurs.eq.1) then
    263           ! on trace continents/ocean 
    264 
    265           ntraceurs_zone_opt=2
    266           izone_cont=1
    267           izone_oce=2         
    268           izone_poubelle=2 ! zone où on met les flux non physiques, de
    269                 ! réajustement
    270           izone_init=2 ! zone d'initialisation par défaut         
    271           option_revap=0
    272           option_tmin=0
    273           izone_revap=0
    274           option_cond=0
    275 
    276           strtrac(izone_cont)='con'
    277           strtrac(izone_oce)='oce'
    278 
    279         elseif (option_traceurs.eq.2) then
    280           ! on trace continent/ ocean/reevap des gouttes
    281 
    282           ntraceurs_zone_opt=3
    283           izone_cont=1
    284           izone_oce=2
    285           izone_poubelle=2 ! zone où on met les flux non physiques, de
    286                 ! réajustement
    287           izone_init=2 ! zone d'initialisation par défaut
    288           option_revap=1
    289           option_tmin=0
    290           izone_revap=3
    291           option_cond=0
    292 
    293           strtrac(izone_cont)='con'
    294           strtrac(izone_oce)='oce'
    295           strtrac(izone_revap)='rev'
    296          
    297 
    298         else if (option_traceurs.eq.3) then
    299             ! on trace des bassins océaniques + un résidu. On ne trace
    300             ! pas l'évap des gouttes à part
    301             ! le résidu est la dernère dimension
    302            
    303           ! lire les use_bassin
    304           call getin('use_bassin_atlantic',use_bassin_atlantic)     
    305           call getin('use_bassin_medit',use_bassin_medit)     
    306           call getin('use_bassin_indian',use_bassin_indian)     
    307           call getin('use_bassin_austral',use_bassin_austral)     
    308           call getin('use_bassin_pacific',use_bassin_pacific)     
    309           call getin('use_bassin_merarabie',use_bassin_merarabie)     
    310           call getin('use_bassin_golfebengale',use_bassin_golfebengale)     
    311           call getin('use_bassin_indiansud',use_bassin_indiansud)     
    312           call getin('use_bassin_tropics',use_bassin_tropics)     
    313           call getin('use_bassin_midlats',use_bassin_midlats)     
    314           call getin('use_bassin_hauteslats',use_bassin_hauteslats)
    315 
    316           write(*,*) 'use_bassin_atlantic=' ,use_bassin_atlantic 
    317           write(*,*) 'use_bassin_medit=' ,use_bassin_medit
    318           write(*,*) 'use_bassin_indian=' ,use_bassin_indian
    319           write(*,*) 'use_bassin_austral=' ,use_bassin_austral
    320           write(*,*) 'use_bassin_merarabie=' ,use_bassin_merarabie
    321           write(*,*) 'use_bassin_golfebengale=' ,use_bassin_golfebengale
    322           write(*,*) 'use_bassin_indiansud=' ,use_bassin_indiansud
    323           write(*,*) 'use_bassin_tropics=' ,use_bassin_tropics
    324           write(*,*) 'use_bassin_midlats=' ,use_bassin_midlats
    325           write(*,*) 'use_bassin_hauteslats=' ,use_bassin_hauteslats
    326 
    327        
    328           ntraceurs_zone_opt=2 &
    329      &                   +use_bassin_atlantic &
    330      &                   +use_bassin_medit &
    331      &                   +use_bassin_indian &
    332      &                   +use_bassin_austral &
    333      &                   +use_bassin_pacific &
    334      &                   +use_bassin_merarabie &
    335      &                   +use_bassin_golfebengale &
    336      &                   +use_bassin_indiansud &
    337      &                   +use_bassin_tropics &
    338      &                   +use_bassin_midlats &
    339      &                   +use_bassin_hauteslats
    340 
    341           izone_cont=ntraceurs_zone
    342           izone_oce=0 ! pas de sens car séparée en bassins         
    343           izone_poubelle=ntraceurs_zone-1 ! zone où on met les flux non physiques, de
    344                 ! réajustement
    345           izone_init=ntraceurs_zone-1 ! zone d'initialisation par défaut
    346           option_revap=0 ! on ne trace pas les gouttes
    347           option_tmin=0
    348           izone_revap=0 ! pas de sens car on taggue pas les gouttes séparemment 
    349           option_cond=0
    350 
    351           ! si on a use_bassin_indian, on n'a pas le découpage détaillé
    352           ! de l'indian:
     146   SUBROUTINE iso_traceurs_init()
     147
     148   USE infotrac_phy, ONLY: itZonIso, isoName, isoZone
     149   USE isotopes_mod, ONLY: iso_eau, ntracisoOR, initialisation_iso
     150   USE dimphy,       ONLY: klon, klev
     151   USE  strings_mod, ONLY: int2str, strStack, strTail, strHead, strIdx, fmsg
     152
     153   IMPLICIT NONE
     154   ! Define which zones and isotopes correspond to isotopic tagging tracers
     155   ! Modify traceurs.h variables
     156   INTEGER :: izone, ixt, k
     157   INTEGER :: izone_pres, izone_lat
     158   INTEGER :: nzone_opt
     159
     160   IF(fmsg("traceurs_init 18: isotrac ne marche que si on met l'eau comme isotope", 'iso_traceurs_init', iso_eau==0)) STOP
     161
     162   !--- Initialize
     163   option_traceurs = 0
     164   initialisation_isotrac = 0
     165
     166   !--- Allocate
     167   ALLOCATE(index_iso (ntiso))
     168   ALLOCATE(index_zone(ntiso))
     169   ALLOCATE(itZonIso_loc(nzone,niso))
     170   ALLOCATE(strtrac(nzone))
     171   ALLOCATE(bassin_map(klon))
     172   ALLOCATE( boite_map(klon,klev))
     173
     174   IF(initialisation_iso == 0) CALL get_in('initialisation_isotrac', initialisation_isotrac)
     175
     176   !--- Read tracing option
     177   CALL get_in('option_traceurs', option_traceurs)
     178
     179   !--- Genral case: no traceurs in ORCHIDEE
     180   ntracisoOR=niso
     181
     182   ! partie a editer ! pour definir les differentes zones
     183   SELECT CASE(option_traceurs)
     184      !========================================================================================================================
     185      CASE(1)      !=== TRACING LAND/OCEAN
     186      !========================================================================================================================
     187         nzone_opt=2
     188         izone_cont=1
     189         izone_oce=2
     190         izone_poubelle=2    ! zone ou on met les flux non physiques, de reajustement
     191         izone_init=2        ! zone d'initialisation par defaut         
     192         option_revap=0
     193         option_tmin=0
     194         izone_revap=0
     195         option_cond=0
     196         strtrac(izone_cont) = 'con'
     197         strtrac(izone_oce)  = 'oce'
     198      !========================================================================================================================
     199      CASE(2)      !=== TRACING LAND/OCEAN/DROPLETS REEVAPORATION
     200      !========================================================================================================================
     201         nzone_opt=3
     202         izone_cont=1
     203         izone_oce=2
     204         izone_poubelle=2    ! zone ou on met les flux non physiques, de reajustement
     205         izone_init=2        ! zone d'initialisation par defaut         
     206         option_revap=1
     207         option_tmin=0
     208         izone_revap=3
     209         option_cond=0
     210         strtrac(izone_cont) = 'con'
     211         strtrac(izone_oce)  = 'oce'
     212         strtrac(izone_revap)= 'rev'
     213      !========================================================================================================================
     214      CASE(3)      !=== TRACING OCEANS BASINS + RESIDUE (LAST DIMENSION). NO DROPLETS EVAPORATION TRACING.
     215      !========================================================================================================================
     216         ! lire les use_bassin
     217         CALL get_in('use_bassin_Atlantic',   use_bassin_Atlantic)
     218         CALL get_in('use_bassin_Medit',      use_bassin_Medit)
     219         CALL get_in('use_bassin_Indian',     use_bassin_Indian)
     220         CALL get_in('use_bassin_Austral',    use_bassin_Austral)
     221         CALL get_in('use_bassin_Pacific',    use_bassin_Pacific)
     222         CALL get_in('use_bassin_MerArabie',  use_bassin_MerArabie)
     223         CALL get_in('use_bassin_BengalGolf', use_bassin_BengalGolf)
     224         CALL get_in('use_bassin_SouthIndian',use_bassin_SouthIndian)
     225         CALL get_in('use_bassin_Tropics',    use_bassin_Tropics)
     226         CALL get_in('use_bassin_Midlats',    use_bassin_Midlats)
     227         CALL get_in('use_bassin_HighLats',   use_bassin_HighLats)
     228         nzone_opt  =  2  +  COUNT([use_bassin_Atlantic, use_bassin_Medit,     use_bassin_Indian,     &
     229            use_bassin_Austral,     use_bassin_Pacific,  use_bassin_MerArabie, use_bassin_BengalGolf, &
     230            use_bassin_SouthIndian, use_bassin_Tropics,  use_bassin_Midlats,   use_bassin_HighLats])
     231         izone_cont=nzone
     232         izone_oce=0             ! pas de sens car separee en bassins         
     233         izone_poubelle=nzone-1  ! zone ou on met les flux non physiques, de reajustement
     234         izone_init=nzone-1      ! zone d'initialisation par defaut
     235         option_revap=0          ! on ne trace pas les gouttes
     236         option_tmin=0
     237         izone_revap=0           ! pas de sens car on taggue pas les gouttes separemment 
     238         option_cond=0
    353239#ifdef ISOVERIF
    354           if (use_bassin_indian.eq.1) then
    355 !              call iso_verif_egalite(float(use_bassin_merarabie), &
    356 !     &            0.0,'iso_traceurs_init 73: revoir def des bassins')
    357                if ((use_bassin_merarabie.ne.0).or. &
    358       &            (use_bassin_indiansud.ne.0).or. &
    359       &            (use_bassin_golfebengale.ne.0)) then
    360                 write(*,*) 'traceurs_init 73'
    361                 stop
    362                endif
    363 !              call iso_verif_egalite(float(use_bassin_golfebengale), &
    364 !     &            0.0,'iso_traceurs_init 73: revoir def des bassins')
    365 !              call iso_verif_egalite(float(use_bassin_indiansud), &
    366 !     &            0.0,'iso_traceurs_init 73: revoir def des bassins')
    367           endif
     240         IF(use_bassin_Indian) THEN   !=== NON COMPATIBLE WITH A DETAILED INDIAN CUTTING
     241            IF(use_bassin_MerArabie .OR. use_bassin_SouthIndian .OR. use_bassin_BengalGolf) THEN
     242               WRITE(*,*) 'traceurs_init 73'; STOP
     243            END IF
     244!           CALL iso_verif_egalite(float(use_bassin_MerArabie),   0.0, 'iso_traceurs_init 73: revoir def des bassins')
     245!           CALL iso_verif_egalite(float(use_bassin_BengalGolf),  0.0, 'iso_traceurs_init 73: revoir def des bassins')
     246!           CALL iso_verif_egalite(float(use_bassin_SouthIndian), 0.0, 'iso_traceurs_init 73: revoir def des bassins')
     247         END IF
    368248#endif   
    369          
    370           bassin_atlantic= max(use_bassin_atlantic,1)
    371           bassin_medit=max(use_bassin_atlantic &
    372      &           +use_bassin_medit,1)
    373           bassin_indian=max(use_bassin_atlantic &
    374      &           +use_bassin_medit &
    375      &           +use_bassin_indian,1)
    376           bassin_austral=max(use_bassin_atlantic &
    377      &           +use_bassin_medit &
    378      &           +use_bassin_indian &
    379      &           +use_bassin_austral,1)
    380           bassin_pacific=max(use_bassin_atlantic &
    381      &           +use_bassin_medit &
    382      &           +use_bassin_indian &
    383      &           +use_bassin_austral &
    384      &           +use_bassin_pacific,1)
    385           bassin_merarabie=max(use_bassin_atlantic &
    386      &           +use_bassin_medit &
    387      &           +use_bassin_indian &
    388      &           +use_bassin_austral &
    389      &           +use_bassin_pacific &
    390      &           +use_bassin_merarabie,1)
    391           bassin_golfebengale=max(use_bassin_atlantic&
    392      &           +use_bassin_medit &
    393      &           +use_bassin_indian &
    394      &           +use_bassin_austral &
    395      &           +use_bassin_pacific &
    396      &           +use_bassin_merarabie &
    397      &           +use_bassin_golfebengale,1)
    398           bassin_indiansud=max(use_bassin_atlantic &
    399      &           +use_bassin_medit &
    400      &           +use_bassin_indian &
    401      &           +use_bassin_austral &
    402      &           +use_bassin_pacific &
    403      &           +use_bassin_merarabie &
    404      &           +use_bassin_golfebengale &
    405      &           +use_bassin_indiansud,1)
    406           bassin_tropics=max(use_bassin_atlantic &
    407      &                       +use_bassin_medit &
    408      &                       +use_bassin_indian &
    409      &                       +use_bassin_austral &
    410      &                       +use_bassin_pacific &
    411      &                       +use_bassin_merarabie &
    412      &                       +use_bassin_golfebengale &
    413      &                       +use_bassin_indiansud, &
    414      &                       +use_bassin_tropics,1)
    415           bassin_midlats=max(use_bassin_atlantic &
    416      &                       +use_bassin_medit &
    417      &                       +use_bassin_indian &
    418      &                       +use_bassin_austral &
    419      &                       +use_bassin_pacific &
    420      &                       +use_bassin_merarabie &
    421      &                       +use_bassin_golfebengale &
    422      &                       +use_bassin_indiansud &
    423      &                       +use_bassin_tropics &
    424      &                       +use_bassin_midlats,1)
    425           bassin_hauteslats=max(use_bassin_atlantic &
    426      &                       +use_bassin_medit &
    427      &                       +use_bassin_indian &
    428      &                       +use_bassin_austral &
    429      &                       +use_bassin_pacific &
    430      &                       +use_bassin_merarabie &
    431      &                       +use_bassin_golfebengale &
    432      &                       +use_bassin_indiansud &
    433      &                       +use_bassin_tropics &
    434      &                       +use_bassin_midlats &
    435      &                       +use_bassin_hauteslats,1)
    436 
    437           write(*,*) 'bassin_atlantic=' ,bassin_atlantic 
    438           write(*,*) 'bassin_medit=' ,bassin_medit
    439           write(*,*) 'bassin_indian=' ,bassin_indian
    440           write(*,*) 'bassin_austral=' ,bassin_austral
    441           write(*,*) 'bassin_merarabie=' ,bassin_merarabie
    442           write(*,*) 'bassin_golfebengale=' ,bassin_golfebengale
    443           write(*,*) 'bassin_indiansud=' ,bassin_indiansud
    444           write(*,*) 'bassin_tropics=' ,bassin_tropics
    445           write(*,*) 'bassin_midlats=' ,bassin_midlats
    446           write(*,*) 'bassin_hauteslats=' ,bassin_hauteslats
    447 
    448           if (use_bassin_atlantic.eq.1) then
    449             strtrac(bassin_atlantic)='atl'
    450           endif
    451           if (use_bassin_medit.eq.1) then
    452             strtrac(bassin_medit)='med'
    453           endif
    454           if (use_bassin_indian.eq.1) then
    455             strtrac(bassin_indian)='ind'
    456           endif
    457           if (use_bassin_austral.eq.1) then
    458             strtrac(bassin_austral)='aus'
    459           endif
    460           if (use_bassin_pacific.eq.1) then
    461             strtrac(bassin_pacific)='pac'
    462           endif
    463           if (use_bassin_merarabie.eq.1) then
    464             strtrac(bassin_merarabie)='ara'
    465           endif
    466           if (use_bassin_golfebengale.eq.1) then
    467             strtrac(bassin_golfebengale)='ben'
    468           endif
    469           if (use_bassin_indiansud.eq.1) then
    470             strtrac(bassin_indiansud)='ins'
    471           endif
    472           if (use_bassin_tropics.eq.1) then
    473             strtrac(bassin_tropics)='tro'
    474           endif
    475           if (use_bassin_midlats.eq.1) then
    476             strtrac(bassin_midlats)='mid'
    477           endif
    478           if (use_bassin_hauteslats.eq.1) then
    479             strtrac(bassin_hauteslats)='hau'
    480           endif
    481           strtrac(ntraceurs_zone-1)='res'
    482           strtrac(ntraceurs_zone)='con'
    483 
    484         else if (option_traceurs.eq.4) then
    485           ! on trace les température minimales vécues
    486           ! comme dans article sur LdG sauf pas de revap
    487            
    488           zone_temp1=293.0 ! en K
    489 !          zone_tempf=223.0 ! en K
    490           zone_tempf=243.0 ! en K
    491  ! courbure de la relation entre l'indice et la température: 0 pour linéaire, <0 pour plus de détal en bas
    492 
     249         bassin_Atlantic   = 1
     250         bassin_Medit      = bassin_Atlantic    + COUNT([use_bassin_Medit]);       WRITE(*,*) 'bassin_Atlantic    =' ,bassin_Atlantic
     251         bassin_Indian     = bassin_Medit       + COUNT([use_bassin_Indian]);      WRITE(*,*) 'bassin_Medit       =' ,bassin_Medit
     252         bassin_Austral    = bassin_Indian      + COUNT([use_bassin_Austral]);     WRITE(*,*) 'bassin_Indian      =' ,bassin_Indian
     253         bassin_Pacific    = bassin_Austral     + COUNT([use_bassin_Pacific]);     WRITE(*,*) 'bassin_Austral     =' ,bassin_Austral
     254         bassin_MerArabie  = bassin_Pacific     + COUNT([use_bassin_MerArabie]);   WRITE(*,*) 'bassin_MerArabie   =' ,bassin_MerArabie
     255         bassin_BengalGolf = bassin_MerArabie   + COUNT([use_bassin_BengalGolf]);  WRITE(*,*) 'bassin_BengalGolf  =' ,bassin_BengalGolf
     256         bassin_SouthIndian= bassin_BengalGolf  + COUNT([use_bassin_SouthIndian]); WRITE(*,*) 'bassin_SouthIndian =' ,bassin_SouthIndian
     257         bassin_Tropics    = bassin_SouthIndian + COUNT([use_bassin_Tropics]);     WRITE(*,*) 'bassin_Tropics     =' ,bassin_Tropics
     258         bassin_MidLats    = bassin_Tropics     + COUNT([use_bassin_MidLats]);     WRITE(*,*) 'bassin_MidLats     =' ,bassin_MidLats
     259         bassin_HighLats   = bassin_MidLats     + COUNT([use_bassin_HighLats]);    WRITE(*,*) 'bassin_HighLats    =' ,bassin_HighLats
     260         IF(use_bassin_atlantic   ) strtrac(bassin_atlantic)   = 'atl'
     261         IF(use_bassin_medit      ) strtrac(bassin_medit)      = 'med'
     262         IF(use_bassin_indian     ) strtrac(bassin_indian)     = 'ind'
     263         IF(use_bassin_austral    ) strtrac(bassin_austral)    = 'aus'
     264         IF(use_bassin_pacific    ) strtrac(bassin_pacific)    = 'pac'
     265         IF(use_bassin_merarabie  ) strtrac(bassin_merarabie)  = 'ara'
     266         IF(use_bassin_BengalGolf ) strtrac(bassin_BengalGolf) = 'ben'
     267         IF(use_bassin_SouthIndian) strtrac(bassin_SouthIndian)= 'ins'
     268         IF(use_bassin_tropics    ) strtrac(bassin_tropics)    = 'tro'
     269         IF(use_bassin_midlats    ) strtrac(bassin_midlats)    = 'mid'
     270         IF(use_bassin_HighLats   ) strtrac(bassin_HighLats)   = 'hau'
     271         strtrac(nzone-1)='res'
     272         strtrac(nzone)='con'
     273      !========================================================================================================================
     274      CASE(4)      !=== TRACING MINIMAL EXPERIENCED TEMPERATURE AS IN THE ARTICLE ON LfG, EXCEPT NO REVAPORATION
     275      !========================================================================================================================
     276         zone_temp1 = 293.0  ! en K
     277!        zone_tempf = 223.0  ! en K
     278         zone_tempf = 243.0  ! en K
     279        ! courbure de la relation entre l'indice et la temperature: 0 pour lineaire, <0 pour plus de detal en bas
    493280        ! zone 1: >= zone_temp1
    494         ! zone 2 à 4: intermédiaire,
     281        ! zone 2 a 4: intermediaire,
    495282        ! zone 5: <zone_tempf
    496        
    497           ntraceurs_zone_opt=nzone_temp+1
    498 
    499           zone_tempa=-4.0 ! en K
    500           izone_cont=ntraceurs_zone
    501           izone_oce=ntraceurs_zone 
    502           izone_poubelle=ntraceurs_zone
    503           izone_init=ntraceurs_zone ! zone d'initialisation par défaut
    504           option_revap=0
    505           option_tmin=0 
    506           izone_revap=0
    507           option_cond=0
    508           do izone=1,nzone_temp
    509             write(strz,'(i2.2)') izone
    510             strtrac(izone)='t'//strz
    511             write(*,*) 'izone,strz,strtrac=',izone,strz,strtrac(izone)
    512           enddo
    513           strtrac(izone_poubelle)='pou'
    514 
    515           ! initialisation des zones de tempéarture
    516           do izone=1,nzone_temp-1
    517             zone_temp(izone)=zone_temp1+float(izone-1) &
    518      &                      *(zone_tempa*float(izone-nzone_temp+1) &
    519      &                      +(zone_tempf-zone_temp1)/float(nzone_temp-2))
    520           enddo
    521           write(*,*) 'iso_trac_init 183: zone_temp=',zone_temp         
    522 
    523         elseif (option_traceurs.eq.5) then
    524           ! on trace AEJ/flux de mousson/Harmattan
    525 !          write(*,*) 'iso_traceurs_init 129'
    526 
    527           ntraceurs_zone_opt=4
    528           izone_cont=1
    529           izone_oce=1
    530           izone_poubelle=1 ! zone où on met les flux non physiques, de
    531                 ! réajustement
    532           izone_init=1 ! zone d'initialisation par défaut
    533           option_revap=0
    534           option_tmin=0
    535           izone_revap=0
    536           izone_aej=2
    537           izone_mousson=3
    538           izone_harmattan=4
    539           option_cond=0
    540 
    541           strtrac(izone_poubelle)='res'
    542           strtrac(izone_aej)='aej'
    543           strtrac(izone_mousson)='mou'
    544           strtrac(izone_harmattan)='sah'
    545 
    546         elseif (option_traceurs.eq.6) then
    547           ! on trace les ddfts
    548 
    549           ntraceurs_zone_opt=2
    550           izone_cont=1
    551           izone_oce=1
    552           izone_poubelle=1 ! zone où on met les flux non physiques, de
    553                 ! réajustement
    554           izone_init=1 ! zone d'initialisation par défaut
    555           option_revap=0
    556           option_tmin=0
    557           izone_revap=0
    558           izone_ddft=2
    559           option_cond=0
    560 
    561           strtrac(izone_poubelle)='res'
    562           strtrac(izone_ddft)='dft'
    563 
    564         elseif (option_traceurs.eq.9) then
    565           ! on trace le condensat
    566 
    567           ntraceurs_zone_opt=3
    568           izone_cont=1
    569           izone_oce=1
    570           izone_poubelle=1 ! zone où on met les flux non physiques, de
    571                 ! réajustement
    572           izone_init=1 ! zone d'initialisation par défaut
    573           option_revap=1
    574           option_tmin=0
    575           izone_revap=2
    576           izone_cond=3
    577           option_cond=1
    578 
    579           ! 1 par défaut pour colorier à la fois condensat LS et
    580           ! condensat convectif. Mais on peut mettre 2 si on ne veut que
    581           ! collorier que le condensat convectif.
    582           call getin('option_cond',option_cond)
    583           write(*,*) 'option_cond=',option_cond
    584 
    585           strtrac(izone_poubelle)='res'
    586           strtrac(izone_cond)='con'
    587           strtrac(izone_revap)='rev'
    588 
    589         elseif (option_traceurs.eq.10) then
    590           ! on trace l'évap venant de ocean/continent no frac/continent frac
    591           !  utilse seulement si couplé avec ORCHIDEE
    592 #ifdef CPP_VEGET
    593 #else
    594           write(*,*) 'iso_traceurs_init 219: option_traceurs=10 ', &
    595      &                      'inutile si on ne couple pas avec ORCHIDEE'
    596           stop
     283         nzone_opt=nzone_temp+1
     284         zone_tempa=-4.0     ! en K
     285         izone_cont=nzone
     286         izone_oce=nzone 
     287         izone_poubelle=nzone
     288         izone_init=nzone    ! zone d'initialisation par defaut
     289         option_revap=0
     290         option_tmin=0 
     291         izone_revap=0
     292         option_cond=0
     293         DO izone=1,nzone_temp
     294            strtrac(izone) = 't'//TRIM(int2str(izone))
     295            WRITE(*,*) 'izone, strtrac=', izone, strtrac(izone)
     296         END DO
     297         strtrac(izone_poubelle)='pou'
     298         ! Initialization of temperatures zones
     299         DO izone=1,nzone_temp-1
     300            zone_temp(izone) = zone_temp1+float(izone-1)            &
     301                            * (zone_tempa*float(izone-nzone_temp+1) &
     302                            + (zone_tempf-zone_temp1)/float(nzone_temp-2))
     303         END DO
     304         WRITE(*,*) 'iso_trac_init 183: zone_temp=', zone_temp
     305      !========================================================================================================================
     306      CASE(5)      !=== TRACING AEJ/MOONSOON FLUX/Harmattan
     307      !========================================================================================================================
     308!        WRITE*,*) 'iso_traceurs_init 129'
     309         nzone_opt=4
     310         izone_cont=1
     311         izone_oce=1
     312         izone_poubelle=1    ! zone ou on met les flux non physiques, de reajustement
     313         izone_init=1        ! zone d'initialisation par defaut         
     314         option_revap=0
     315         option_tmin=0
     316         izone_revap=0
     317         izone_aej=2
     318         izone_mousson=3
     319         izone_harmattan=4
     320         option_cond=0
     321         strtrac(izone_poubelle) = 'res'
     322         strtrac(izone_aej)      = 'aej'
     323         strtrac(izone_mousson)  = 'mou'
     324         strtrac(izone_harmattan)= 'sah'
     325      !========================================================================================================================
     326      CASE(6)      !=== TRACING DDFTS
     327      !========================================================================================================================
     328         nzone_opt=2
     329         izone_cont=1
     330         izone_oce=1
     331         izone_poubelle=1    ! zone ou on met les flux non physiques, de reajustement
     332         izone_init=1        ! zone d'initialisation par defaut         
     333         option_revap=0
     334         option_tmin=0
     335         izone_revap=0
     336         izone_ddft=2
     337         option_cond=0
     338         strtrac(izone_poubelle)='res'
     339         strtrac(izone_ddft)='dft'
     340      !========================================================================================================================
     341      CASE(9)      !=== TRACING CONDENSATION
     342      !========================================================================================================================
     343         nzone_opt=3
     344         izone_cont=1
     345         izone_oce=1
     346         izone_poubelle=1    ! zone ou on met les flux non physiques, de reajustement
     347         izone_init=1        ! zone d'initialisation par defaut         
     348         option_revap=1
     349         option_tmin=0
     350         izone_revap=2
     351         izone_cond=3
     352         option_cond=1
     353         ! 1 par defaut pour colorier a la fois condensat LS et condensat convectif.
     354         ! Mais on peut mettre 2 si on ne veut que colorier que le condensat convectif.
     355         CALL get_in('option_cond',option_cond)
     356         strtrac(izone_poubelle)='res'
     357         strtrac(izone_cond)='con'
     358         strtrac(izone_revap)='rev'
     359      !========================================================================================================================
     360      CASE(10)     !=== TRACING EVAPORATION FROM OCEAN/LAND, NON FRAC/LAND FRAC ; ONLY WHEN COUPLED WITH ORCHIDEE
     361      !========================================================================================================================
     362#ifndef CPP_VEGET
     363         WRITE(*,*) 'iso_traceurs_init 219: option_traceurs=10 inutile si on ne couple pas avec ORCHIDEE'; STOP
    597364#endif         
    598 
    599           ntraceurs_zone_opt=3
    600           izone_cont=1 ! sous-entendu non fractionnant
    601           izone_oce=2
    602           izone_poubelle=2 ! zone où on met les flux non physiques, de
    603                 ! réajustement
    604           izone_init=2 ! zone d'initialisation par défaut
    605           option_revap=0
    606           option_tmin=0
    607           izone_revap=0
    608           izone_contfrac=3
    609           izone_contcanop=3
    610           izone_irrig=0
    611           option_cond=0
    612 
    613           strtrac(izone_oce)='oce'
    614           strtrac(izone_cont)='con' 
    615           strtrac(izone_contfrac)='enu'  ! evap sol nu
    616 
    617         elseif (option_traceurs.eq.11) then
    618           ! on trace reevap des gouttes et le reste
    619 
    620           ntraceurs_zone_opt=2
    621           izone_cont=1
    622           izone_oce=1
    623           izone_poubelle=1 ! zone où on met les flux non physiques, de
    624                 ! réajustement
    625           izone_init=1 ! zone d'initialisation par défaut
    626           option_revap=1
    627           option_tmin=0
    628           izone_revap=2
    629           izone_irrig=0
    630           option_cond=0
    631 
    632           strtrac(izone_poubelle)='res'
    633           strtrac(izone_revap)='rev'
    634 
    635         elseif (option_traceurs.eq.12) then
    636           ! on trace evap du sol nu, evap de la canopée, reste de l'evap cont et
    637           ! evap oce
    638 #ifdef CPP_VEGET
    639 #else
    640           write(*,*) 'iso_traceurs_init 257: option_traceurs=10 ', &
    641      &                      'inutile si on ne couple pas avec ORCHIDEE'
    642           stop
     365         nzone_opt=3
     366         izone_cont=1        ! sous-entendu non fractionnant
     367         izone_oce=2
     368         izone_poubelle=2    ! zone ou on met les flux non physiques, de reajustement
     369         izone_init=2        ! zone d'initialisation par defaut
     370         option_revap=0
     371         option_tmin=0
     372         izone_revap=0
     373         izone_contfrac=3
     374         izone_contcanop=3
     375         izone_irrig=0
     376         option_cond=0
     377         strtrac(izone_oce)='oce'
     378         strtrac(izone_cont)='con' 
     379         strtrac(izone_contfrac)='enu'  ! evap sol nu
     380      !========================================================================================================================
     381      CASE(11)     !=== TRACING DROPLETS REEVAPORATION + REST
     382      !========================================================================================================================
     383         nzone_opt=2
     384         izone_cont=1
     385         izone_oce=1
     386         izone_poubelle=1    ! zone ou on met les flux non physiques, de reajustement
     387         izone_init=1        ! zone d'initialisation par defaut
     388         option_revap=1
     389         option_tmin=0
     390         izone_revap=2
     391         izone_irrig=0
     392         option_cond=0
     393         strtrac(izone_poubelle)='res'
     394         strtrac(izone_revap)='rev'
     395      !========================================================================================================================
     396      CASE(12)     !=== TRACING NAKED GROUND EVAPORATION, CANOPY EVAPORATION, REST OF LAND EVAPORATION AND OCEAN EVAPORATION
     397      !========================================================================================================================
     398#ifndef CPP_VEGET
     399         WRITE(*,*) 'iso_traceurs_init 257: option_traceurs=10 inutile si on ne couple pas avec ORCHIDEE'; STOP
    643400#endif           
    644 
    645           ntraceurs_zone_opt=2
    646           izone_cont=1
    647           izone_oce=2
    648           izone_poubelle=2 ! zone où on met les flux non physiques, de
    649                 ! réajustement
    650           izone_init=2 ! zone d'initialisation par défaut
    651           option_revap=0
    652           option_tmin=0
    653           izone_revap=0
    654           izone_contfrac=3
    655           izone_contcanop=4
    656           izone_irrig=0   
    657           option_cond=0
    658 
    659           strtrac(izone_oce)='oce'
    660           strtrac(izone_cont)='con'
    661           strtrac(izone_contfrac)='enu'  ! evap sol nu
    662           strtrac(izone_contcanop)='eca'  ! evap canop
    663 
    664        else if (option_traceurs.eq.13) then
    665           ! on trace les température minimales vécues + la revap
    666           ! comme dans article sur LdG
    667            
    668         zone_temp1=293.0         ! en K       
    669 !        parameter (zone_tempf=223.0) ! en K
    670         zone_tempf=243.0 ! en K
    671         zone_tempa=-4.0 ! courbure de la relation entre l'indice et la température: 0 pour linéaire, <0 pour plus de détal en bas
    672 
    673         ! zone 1: >= zone_temp1
    674         ! zone 2 à 4: intermédiaire,
    675         ! zone 5: <zone_tempf
    676        
    677           ntraceurs_zone_opt=nzone_temp+1
    678          
    679           izone_cont=1
    680           izone_oce=1 
    681           izone_poubelle=1
    682           izone_init=1 ! zone d'initialisation par défaut
    683           option_revap=1   
    684           option_tmin=0
    685           izone_revap=ntraceurs_zone
    686           izone_irrig=0
    687           option_cond=0
    688           do izone=1,nzone_temp
    689             write(strz,'(i2.2)') izone
    690             strtrac(izone)='t'//strz
    691             write(*,*) 'izone,strz,strtrac=',izone,strz,strtrac(izone)
    692           enddo
    693           strtrac(izone_revap)='rev'
    694 
    695           ! initialisation des zones de tempéarture
    696           do izone=1,nzone_temp-1
    697             zone_temp(izone)=zone_temp1+float(izone-1) &
    698      &                      *(zone_tempa*float(izone-nzone_temp+1) &
    699      &                      +(zone_tempf-zone_temp1)/float(nzone_temp-2))
    700           enddo
    701           write(*,*) 'zone_temp=',zone_temp
    702 
    703        else if (option_traceurs.eq.14) then
    704           ! on trace les pres et lat de dernière saturation définies
    705           ! comme rh>90%
    706            
    707         zone_pres1=600.0*100.0 ! en Pa       
    708         zone_presf=300.0*100.0 ! en Pa
    709         zone_presa=0.0 ! courbure de la relation entre l'indice et la température: 0 pour linéaire, <0 pour plus de détal en bas
    710 
    711         lattag_min=10.0 ! en degrès
    712         dlattag=15.0
    713 
    714         ! zone 1: >= zone_pres1
    715         ! zone 2 à 4: intermédiaire,
    716         ! zone 5: <zone_presf
    717        
    718          ntraceurs_zone_opt=nzone_pres*nzone_lat+1         
    719           izone_cont=ntraceurs_zone
    720           izone_oce=ntraceurs_zone
    721           izone_poubelle=ntraceurs_zone
    722           izone_init=ntraceurs_zone ! zone d'initialisation par défaut
    723           option_revap=0 
    724           option_tmin=0
    725           izone_revap=0
    726           izone_irrig=0
    727           option_cond=0
    728           do izone_pres=1,nzone_pres
    729            do izone_lat=1,nzone_lat
    730             write(strz_pres,'(i1.1)') izone_pres
    731             write(strz_lat,'(i1.1)') izone_lat
    732             strz_preslat=strz_pres//strz_lat
    733             izone=izone_lat+(izone_pres-1)*nzone_lat
    734             strtrac(izone)='t'//strz_preslat
    735             write(*,*) 'izone_pres,izone_lat,strtrac=', &
    736      &                        izone_pres,izone_lat,izone,strtrac(izone)
    737            enddo !do izone_lat=1,nzone_lat
    738           enddo !do izone_pres=1,nzone_pres
    739           strtrac(ntraceurs_zone)='sfc'
    740 
    741           ! initialisation des zones de tempéarture
    742           do izone=1,nzone_pres-1
    743             zone_pres(izone)=zone_pres1+float(izone-1) &
    744      &                      *(zone_presa*float(izone-nzone_pres+1) &
    745      &                      +(zone_presf-zone_pres1)/float(nzone_pres-2))
    746           enddo !do izone=1,nzone_pres-1
    747           write(*,*) 'traceurs_init 332: zone_pres=',zone_pres
    748 !          stop
    749 !
    750        elseif (option_traceurs.eq.15) then
    751           ! on trace l'irrigation dans ORCHIDEE
    752 #ifdef CPP_VEGET
    753 #else
    754           write(*,*) 'iso_traceurs_init 257: option_traceurs=15 ', &
    755      &                      'inutile si on ne couple pas avec ORCHIDEE'
    756           stop
     401         nzone_opt=2
     402         izone_cont=1
     403         izone_oce=2
     404         izone_poubelle=2    ! zone ou on met les flux non physiques, de reajustement
     405         izone_init=2        ! zone d'initialisation par defaut
     406         option_revap=0
     407         option_tmin=0
     408         izone_revap=0
     409         izone_contfrac=3
     410         izone_contcanop=4
     411         izone_irrig=0   
     412         option_cond=0
     413         strtrac(izone_oce)='oce'
     414         strtrac(izone_cont)='con'
     415         strtrac(izone_contfrac)='enu' ! evap sol nu
     416         strtrac(izone_contcanop)='eca'! evap canop
     417      !========================================================================================================================
     418      CASE(13)     !=== TRACING MINIMUM EXPERIENCED TEMPERATIRES + REEVAPORATION AS IN THE ARTICLE ON LdG
     419      !========================================================================================================================
     420         zone_temp1=293.0    ! en K       
     421!        zone_tempf=223.0    ! en K
     422         zone_tempf=243.0    ! en K
     423         zone_tempa=-4.0     ! courbure de la relation entre l'indice et la temperature: 0 pour lineaire, <0 pour plus de detal en bas
     424         ! zone 1: >= zone_temp1
     425         ! zone 2 a 4: intermediaire,
     426         ! zone 5: <zone_tempf
     427         nzone_opt=nzone_temp+1
     428         izone_cont=1
     429         izone_oce=1 
     430         izone_poubelle=1
     431         izone_init=1        ! zone d'initialisation par defaut
     432         option_revap=1   
     433         option_tmin=0
     434         izone_revap=nzone
     435         izone_irrig=0
     436         option_cond=0
     437         DO izone=1,nzone_temp
     438            strtrac(izone) = 't'//TRIM(int2str(izone))
     439            WRITE(*,*) 'izone, strtrac = ', izone, strtrac(izone)
     440         END DO
     441         strtrac(izone_revap)='rev'
     442         ! initialisation des zones de tempearture
     443         DO izone=1,nzone_temp-1
     444            zone_temp(izone) = zone_temp1+float(izone-1) &
     445                             *(zone_tempa*float(izone-nzone_temp+1) &
     446                             +(zone_tempf-zone_temp1)/float(nzone_temp-2))
     447         END DO
     448         WRITE(*,*) 'zone_temp=',zone_temp
     449      !========================================================================================================================
     450      CASE(14)     !=== TRACING PRES AND LAT OF LAST SATURATION DEFINED AS rh>90%
     451      !========================================================================================================================
     452         zone_pres1=600.0*100.0   ! en Pa       
     453         zone_presf=300.0*100.0   ! en Pa
     454         zone_presa=0.0           ! courbure de la relation entre l'indice et la temperature: 0 pour lineaire
     455         lattag_min=10.0          ! en degres
     456         dlattag=15.0
     457         ! zone 1: >= zone_pres1
     458         ! zone 2 a 4: intermediaire,
     459         ! zone 5: <zone_presf
     460         nzone_opt=nzone_pres*nzone_lat+1         
     461         izone_cont=nzone
     462         izone_oce=nzone
     463         izone_poubelle=nzone
     464         izone_init=nzone         ! zone d'initialisation par defaut
     465         option_revap=0 
     466         option_tmin=0
     467         izone_revap=0
     468         izone_irrig=0
     469         option_cond=0
     470         DO izone_pres=1,nzone_pres
     471            DO izone_lat=1,nzone_lat
     472               izone=izone_lat+(izone_pres-1)*nzone_lat
     473               strtrac(izone) = 't'//TRIM(int2str(izone_pres))//TRIM(int2str(izone_lat))
     474               write(*,*) 'izone_pres, izone_lat, izone, strtrac = ',izone_pres, izone_lat, izone, strtrac(izone)
     475            END DO
     476         END DO
     477         strtrac(nzone)='sfc'
     478         ! initialisation des zones de temperature
     479         DO izone=1,nzone_pres-1
     480            zone_pres(izone) = zone_pres1+float(izone-1) &
     481                             *(zone_presa*float(izone-nzone_pres+1) &
     482                             +(zone_presf-zone_pres1)/float(nzone_pres-2))
     483         END DO
     484         WRITE(*,*) 'traceurs_init 332: zone_pres=',zone_pres
     485      !========================================================================================================================
     486      CASE(15)     !=== TRACING IRRIGATION IN ORCHIDEE
     487      !========================================================================================================================
     488#ifndef CPP_VEGET
     489         WRITE(*,*) 'iso_traceurs_init 257: option_traceurs=15 inutile si on ne couple pas avec ORCHIDEE'; STOP
    757490#endif
    758 
    759           ntraceurs_zone_opt=1
    760           izone_cont=1
    761           izone_oce=1
    762           izone_poubelle=1 ! zone où on met les flux non physiques, de
    763                 ! réajustement
    764           izone_init=1 ! zone d'initialisation par défaut
    765           option_revap=0
    766           option_tmin=0
    767           izone_revap=0
    768           izone_contfrac=0
    769           izone_contcanop=0
    770           izone_irrig=2
    771           option_cond=0
    772          
    773           strtrac(izone_poubelle)='res'
    774           strtrac(izone_irrig)='irrig'
    775 
    776           ! dans ce cas particulier, il y a des traceurs dans ORCHIDEE
    777           ntracisoOR=ntiso
    778 
    779         else if ((option_traceurs.eq.17).or. &
    780      &           (option_traceurs.eq.18)) then
    781           ! on trace les température minimales vécues
    782           ! comme dans article sur LdG sauf pas de revap
    783            
    784         zone_temp1=12.0e-3 ! en kg/kg       
    785         zone_tempf=0.2e-3 ! en kg/kg
    786         zone_tempa=1.2e-3 ! courbure de la relation entre l'indice et la température: 0 pour linéaire, <0 pour plus de détail en bas
    787 
    788 !       parameter (zone_temp1=14.0e-3) ! en kg/kg       
    789 !       parameter (zone_tempf=0.2e-3) ! en kg/kg
    790 !       parameter (zone_tempa=0.5e-3)       
    791 
    792 !        parameter (zone_temp1=10.0e-3) ! en kg/kg
    793 !       parameter (zone_tempf=0.5e-3) ! en kg/kg
    794 !       parameter (zone_tempa=0.5e-3)
    795 
    796         ! zone 1: >= zone_temp1
    797         ! zone 2 à 4: intermédiaire,
    798         ! zone 5: <zone_tempf
    799        
    800         ntraceurs_zone_opt=nzone_temp+3
    801        
    802           izone_cont=nzone_temp+1
    803           izone_oce=nzone_temp+1
    804           izone_poubelle=nzone_temp+1
    805           izone_init=nzone_temp+1 ! zone d'initialisation par défaut
    806           option_revap=1 
    807           option_tmin=1
    808           option_cond=1
    809 
    810           izone_revap=nzone_temp+3
    811           izone_cond=nzone_temp+2
    812           do izone=1,nzone_temp
    813             write(strz,'(i2.2)') izone
    814             strtrac(izone)='t'//strz
    815             write(*,*) 'izone,strz,strtrac=',izone,strz,strtrac(izone)
    816           enddo !do izone=1,nzone_temp
    817           strtrac(izone_poubelle)='sfc'
    818           strtrac(izone_cond)='con'
    819           strtrac(izone_revap)='rev'
    820 
    821           ! initialisation des zones de tempéarture
    822           do izone=1,nzone_temp-1
    823             zone_temp(izone)=zone_temp1+float(izone-1) &
    824      &                      *(zone_tempa*float(izone-nzone_temp+1) &
    825      &             +(zone_tempf-zone_temp1)/float(nzone_temp-2))
    826           enddo
    827          write(*,*) 'zone_temp1,zone_tempf,zone_tempa=', &
    828      &              zone_temp1,zone_tempf,zone_tempa
    829           write(*,*) 'zone_temp=',zone_temp
    830 !          stop         
    831 
    832         else if (option_traceurs.eq.19) then
    833 
    834         zone_temp1=12.0e-3 ! en kg/kg       
    835         zone_tempf=0.2e-3 ! en kg/kg
    836         zone_tempa=1.2e-3 ! courbure de la relation entre l'indice et la température: 0 pour linéaire, <0 pour plus de détail en bas
    837 
    838 !       parameter (zone_temp1=14.0e-3) ! en kg/kg       
    839 !       parameter (zone_tempf=0.2e-3) ! en kg/kg
    840 !       parameter (zone_tempa=0.5e-3)       
    841 
    842 !        parameter (zone_temp1=10.0e-3) ! en kg/kg
    843 !       parameter (zone_tempf=0.5e-3) ! en kg/kg
    844 !       parameter (zone_tempa=0.5e-3)
    845 
    846         ! zone 1: >= zone_temp1
    847         ! zone 2 à 4: intermédiaire,
    848         ! zone 5: <zone_tempf
    849        
    850         ntraceurs_zone_opt=nzone_temp+4
    851        
    852           izone_cont=nzone_temp+1
    853           izone_oce=nzone_temp+1
    854           izone_poubelle=nzone_temp+1
    855           if (option_seuil_tag_tmin.eq.1) then
    856             izone_init=nzone_temp+1 ! zone d'initialisation par défaut
    857           else
     491         nzone_opt=1
     492         izone_cont=1
     493         izone_oce=1
     494         izone_poubelle=1    ! zone ou on met les flux non physiques, de reajustement
     495         izone_init=1        ! zone d'initialisation par defaut
     496         option_revap=0
     497         option_tmin=0
     498         izone_revap=0
     499         izone_contfrac=0
     500         izone_contcanop=0
     501         izone_irrig=2
     502         option_cond=0
     503         strtrac(izone_poubelle)='res'
     504         strtrac(izone_irrig)='irrig'
     505         ! dans ce cas particulier, il y a des traceurs dans ORCHIDEE
     506         ntracisoOR=ntiso
     507      !========================================================================================================================
     508      CASE(17,18)  !=== TRACING MINIMAL EXPERIENCES TEMPERATURES AS IN THE ARTICLE ABOUT LdG, BUT NO EVAPORATION
     509      !========================================================================================================================
     510         zone_temp1=12.0e-3  ! en kg/kg       
     511         zone_tempf=0.2e-3   ! en kg/kg
     512         zone_tempa=1.2e-3   ! courbure de la relation entre l'indice et la temperature: 0 pour lineaire
     513!        zone_temp1=14.0e-3  ! en kg/kg       
     514!        zone_tempf=0.2e-3   ! en kg/kg
     515!        zone_tempa=0.5e-3       
     516!        zone_temp1=10.0e-3  ! en kg/kg
     517!        zone_tempf=0.5e-3   ! en kg/kg
     518!        zone_tempa=0.5e-3
     519         ! zone 1: >= zone_temp1
     520         ! zone 2 a 4: intermediaire,
     521         ! zone 5: <zone_tempf
     522         nzone_opt=nzone_temp+3
     523         izone_cont=nzone_temp+1
     524         izone_oce=nzone_temp+1
     525         izone_poubelle=nzone_temp+1
     526         izone_init=nzone_temp+1 ! zone d'initialisation par defaut
     527         option_revap=1 
     528         option_tmin=1
     529         option_cond=1
     530         izone_revap=nzone_temp+3
     531         izone_cond=nzone_temp+2
     532         DO izone=1,nzone_temp
     533            strtrac(izone) = 't'//TRIM(int2str(izone))
     534            WRITE(*,*) 'izone, strtrac = ', izone, strtrac(izone)
     535         END DO !do izone=1,nzone_temp
     536         strtrac(izone_poubelle)='sfc'
     537         strtrac(izone_cond)='con'
     538         strtrac(izone_revap)='rev'
     539         ! initialisation des zones de tempearture
     540         DO izone=1,nzone_temp-1
     541            zone_temp(izone) = zone_temp1+float(izone-1) &
     542                             *(zone_tempa*float(izone-nzone_temp+1) &
     543                             +(zone_tempf-zone_temp1)/float(nzone_temp-2))
     544         END DO
     545         WRITE(*,*) 'zone_temp1,zone_tempf,zone_tempa=',zone_temp1,zone_tempf,zone_tempa
     546         WRITE(*,*) 'zone_temp=',zone_temp
     547!        STOP         
     548      !========================================================================================================================
     549      CASE(19)     !=== TRACING TROPICAL AND EXTRATROPICAL VAPOUR
     550      !========================================================================================================================
     551         zone_temp1=12.0e-3  ! en kg/kg       
     552         zone_tempf=0.2e-3   ! en kg/kg
     553         zone_tempa=1.2e-3   ! courbure de la relation entre l'indice et la temperature: 0 pour lineaire, <0 pour plus de detail en bas
     554!        zone_temp1=14.0e-3  ! en kg/kg       
     555!        zone_tempf=0.2e-3   ! en kg/kg
     556!        zone_tempa=0.5e-3
     557!        zone_temp1=10.0e-3  ! en kg/kg       
     558!        zone_tempf=0.5e-3   ! en kg/kg
     559!        zone_tempa=0.5e-3
     560         ! zone 1: >= zone_temp1
     561         ! zone 2 a 4: intermediaire,
     562         ! zone 5: <zone_tempf
     563         nzone_opt=nzone_temp+4
     564         izone_cont=nzone_temp+1
     565         izone_oce=nzone_temp+1
     566         izone_poubelle=nzone_temp+1
     567         IF(option_seuil_tag_tmin == 1) THEN
     568            izone_init=nzone_temp+1 ! zone d'initialisation par defaut
     569         ELSE
    858570            izone_init=nzone_temp
    859           endif
    860           option_revap=1   
    861           izone_revap=nzone_temp+3
    862           izone_cond=nzone_temp+2
    863           izone_ddft=nzone_temp+4
    864           option_tmin=1         
    865           option_cond=1
    866           do izone=1,nzone_temp
    867             write(strz,'(i2.2)') izone
    868             strtrac(izone)='t'//strz
    869             write(*,*) 'izone,strz,strtrac=',izone,strz,strtrac(izone)
    870           enddo !do izone=1,nzone_temp
    871           strtrac(izone_poubelle)='sfc'
    872           strtrac(izone_cond)='con'
    873           strtrac(izone_revap)='rev'
    874           strtrac(izone_ddft)='dft'
    875 
    876         elseif (option_traceurs.eq.20) then
    877           ! on vapeur tropical/extractropicale/recyclage extractropical
    878           ! pour comprendre controles humidité et isotopes subtropicaux.       
    879          
    880           lim_tag20=35.0
    881           call getin('lim_tag20',lim_tag20)
    882           write(*,*) 'lim_tag20=',lim_tag20
    883 
    884           ntraceurs_zone_opt=3
    885           izone_cont=1
    886           izone_oce=1
    887           izone_poubelle=2 ! zone où on met les flux non physiques, de
    888                 ! réajustement
    889           izone_init=2 ! zone d'initialisation par défaut
    890           option_revap=0
    891           option_tmin=0
    892           izone_revap=0
    893           izone_trop=2
    894           izone_extra=3
    895 
    896           strtrac(izone_trop)='tro' ! vapeur tropicale
    897           strtrac(izone_extra)='ext' ! vapeur extractropicale evaporée
    898                 ! dans les tropiques
    899           strtrac(izone_cont)='rec' ! recyclage
    900 
    901         elseif (option_traceurs.eq.21) then
    902           ! on trace 2 boites 3D: UT tropicale et extratropiques
    903           ! fonctionnement similaire à option 5 pour taggage des zones
    904           ! AMMA
    905 !          write(*,*) 'iso_traceurs_init 129'
    906 
    907           ntraceurs_zone_opt=3
    908           izone_cont=1
    909           izone_oce=1
    910           izone_poubelle=1 ! zone où on met les flux non physiques, de
    911                 ! réajustement
    912           izone_init=1 ! zone d'initialisation par défaut
    913           option_revap=0
    914           option_tmin=0
    915           izone_revap=0
    916           izone_trop=2
    917           izone_extra=3
    918           option_cond=0
    919 
    920           strtrac(izone_poubelle)='res'
    921           strtrac(izone_trop)='tro'
    922           strtrac(izone_extra)='ext'
    923 
    924         elseif (option_traceurs.eq.22) then
    925           ! on trace la vapeur qui a été processée dans zones de
    926           ! convections à 3 niveaux: BT, MT et UT
    927 
    928           lim_precip_tag22=20.0
    929           call getin('lim_precip_tag22',lim_precip_tag22)
    930           write(*,*) 'lim_precip_tag22=',lim_precip_tag22
    931 
    932           ntraceurs_zone_opt=3
    933           izone_cont=1
    934           izone_oce=1
    935           izone_poubelle=1 ! zone où on met les flux non physiques, de
    936                 ! réajustement
    937           izone_init=1 ! zone d'initialisation par défaut
    938           option_revap=0
    939           option_tmin=0
    940           izone_revap=0
    941           izone_conv_BT=2
    942           izone_conv_UT=3
    943           option_cond=0
    944 
    945           strtrac(izone_poubelle)='res'
    946           strtrac(izone_conv_BT)='cbt'
    947           strtrac(izone_conv_UT)='cut'
    948 
    949         else
    950             write(*,*) 'traceurs_init 36: option pas encore prévue'
    951             stop
    952         endif
    953 
    954        
    955           if (ntraceurs_zone_opt.ne.ntraceurs_zone) then
    956                 write(*,*) 'ntraceurs_zone_opt,ntraceurs_zone=', &
    957                         & ntraceurs_zone_opt,ntraceurs_zone
    958                 call abort_physic ('isotrac_mod','ntraceurs_zone incoherent',1)
    959           endif
    960 
    961        
    962         ! seuil sur le taux de condensation
    963         if (option_tmin.eq.1) then
    964           seuil_tag_tmin=0.01
    965           call getin('seuil_tag_tmin',seuil_tag_tmin)
    966           write(*,*) 'seuil_tag_tmin=',seuil_tag_tmin
    967 
    968           seuil_tag_tmin_ls=seuil_tag_tmin
    969           call getin('seuil_tag_tmin_ls',seuil_tag_tmin_ls)
    970           write(*,*) 'seuil_tag_tmin_ls=',seuil_tag_tmin_ls
    971 
    972           option_seuil_tag_tmin=1
    973           call getin('option_seuil_tag_tmin',option_seuil_tag_tmin)
    974           write(*,*) 'option_seuil_tag_tmin=',option_seuil_tag_tmin
    975         endif
    976 
    977 
    978         do ixt=1,niso
    979            index_zone(ixt)=0
    980            index_iso(ixt)=ixt
    981         enddo
    982         itrac=niso       
    983         do izone=1,ntraceurs_zone
    984           do ixt=1,niso
    985             itrac=itrac+1
    986             index_zone(itrac)=izone
    987             index_iso(itrac)=ixt
    988             itZonIso_loc(izone,ixt)=itrac
    989             if (itZonIso(izone,ixt).ne.itZonIso_loc(izone,ixt)) then
    990                 write(*,*) 'isotrac 989: izone,ixt,itrac=',izone,ixt,itrac
    991                 CALL abort_physic ('isotrac','isotrac 989',1)
    992             endif
    993           enddo
    994         enddo
     571         END IF
     572         option_revap=1   
     573         izone_revap=nzone_temp+3
     574         izone_cond=nzone_temp+2
     575         izone_ddft=nzone_temp+4
     576         option_tmin=1         
     577         option_cond=1
     578         DO izone=1,nzone_temp
     579            strtrac(izone) = 't'//TRIM(int2str(izone))
     580            WRITE(*,*) 'izone, strtrac = ', izone, strtrac(izone)
     581         END DO
     582         strtrac(izone_poubelle)='sfc'
     583         strtrac(izone_cond)='con'
     584         strtrac(izone_revap)='rev'
     585         strtrac(izone_ddft)='dft'
     586      !========================================================================================================================
     587      CASE(20)     !=== TRACING TROPICAL/EXTRATROPICAL/EXTRATROPICAL RECYCLING TO STUDY HUMIDITY AND SUBTROPICAL ISOTOPES CONTROL
     588      !========================================================================================================================
     589         CALL get_in('lim_tag20', lim_tag20, 35.0)
     590         nzone_opt=3
     591         izone_cont=1
     592         izone_oce=1
     593         izone_poubelle=2    ! zone ou on met les flux non physiques, de reajustement
     594         izone_init=2        ! zone d'initialisation par defaut
     595         option_revap=0
     596         option_tmin=0
     597         izone_revap=0
     598         izone_trop=2
     599         izone_extra=3
     600         strtrac(izone_trop)='tro'     ! tropical vapour
     601         strtrac(izone_extra)='ext'    ! extratropical vapour evaporated in the tropics
     602         strtrac(izone_cont)='rec'     ! recycling
     603      !========================================================================================================================
     604      CASE(21)     !=== TRACING TWO 3D BOXES: TROPICAL UT AND EXTRATROPICS ; SIMILAR TO 5 FOR AMMA ZONES TAGGING
     605      !========================================================================================================================
     606!        WRITE(*,*) 'iso_traceurs_init 129'
     607         nzone_opt=3
     608         izone_cont=1
     609         izone_oce=1
     610         izone_poubelle=1    ! zone ou on met les flux non physiques, de reajustement
     611         izone_init=1        ! zone d'initialisation par defaut
     612         option_revap=0
     613         option_tmin=0
     614         izone_revap=0
     615         izone_trop=2
     616         izone_extra=3
     617         option_cond=0
     618         strtrac(izone_poubelle)='res'
     619         strtrac(izone_trop)='tro'
     620         strtrac(izone_extra)='ext'
     621      !========================================================================================================================
     622      CASE(22)     !=== TRACING WATER VAPOUR PROCESSED IN THE 3-LEVELS SCONVECTION ZONES BT, MT AND UT
     623      !========================================================================================================================
     624         CALL get_in('lim_precip_tag22', lim_precip_tag22, 20.0)
     625         nzone_opt=3
     626         izone_cont=1
     627         izone_oce=1
     628         izone_poubelle=1    ! zone ou on met les flux non physiques, de reajustement
     629         izone_init=1        ! zone d'initialisation par defaut
     630         option_revap=0
     631         option_tmin=0
     632         izone_revap=0
     633         izone_conv_BT=2
     634         izone_conv_UT=3
     635         option_cond=0
     636         strtrac(izone_poubelle)='res'
     637         strtrac(izone_conv_BT)='cbt'
     638         strtrac(izone_conv_UT)='cut'
     639      CASE DEFAULT
     640         WRITE(*,*) 'traceurs_init 36: option pas encore prevue' ; STOP
     641   END SELECT
     642
     643   IF(nzone_opt /= nzone) THEN
     644      WRITE(*,*) 'nzone_opt, nzone=', nzone_opt, nzone
     645      CALL abort_physic ('isotrac_mod','nzone incoherent',1)
     646   END IF
     647
     648   !--- Condensation rate threshold
     649   IF(option_tmin == 1) THEN
     650      seuil_tag_tmin = 0.01
     651      CALL get_in('seuil_tag_tmin',        seuil_tag_tmin,        0.01)
     652      CALL get_in('seuil_tag_tmin_ls',     seuil_tag_tmin_ls,     seuil_tag_tmin)
     653      CALL get_in('option_seuil_tag_tmin', option_seuil_tag_tmin, 1)
     654   END IF
     655
     656   index_zone = [(strIdx(isoZone, strTail(isoName(ixt) ,'_',.TRUE.)), ixt=1, ntiso)]
     657   index_iso  = [(strIdx(isoName, strHead(isoName(ixt) ,'_',.TRUE.)), ixt=1, ntiso)]
     658   itZonIso_loc = itZonIso(:,:)
    995659#ifdef ISOVERIF
    996 !        call iso_verif_egalite(float(itrac),float(ntiso), &
    997 !     &           'traceurs_init 50')
    998         if (itrac.ne.ntiso) then
    999           write(*,*) 'traceurs_init 50'
    1000           stop
    1001         endif
    1002      
    1003         write(*,*) 'traceurs_init 65: bilan de l''init:'
    1004         write(*,*) 'index_zone=',index_zone(1:ntiso)
    1005         write(*,*) 'index_iso=',index_iso(1:ntiso)
    1006         write(*,*) 'itZonIso=',itZonIso(1:ntraceurs_zone,1:niso)
    1007         do izone=1,ntraceurs_zone
    1008           write(*,*) 'strtrac(',izone,')=',strtrac(izone)
    1009         enddo !do izone=1,ntraceurs_zone
    1010         write(*,*) 'ntracisoOR=',ntracisoOR
     660   WRITE(*,*) 'traceurs_init 65: bilan de l''init:'
     661   WRITE(*,*) 'index_zone = '//TRIM(strStack(int2str(index_zone(1:ntiso))))
     662   WRITE(*,*) 'index_iso  = '//TRIM(strStack(int2str(index_iso (1:ntiso))))
     663   DO izone=1,nzone
     664      WRITE(*,*)'itZonIso('//TRIM(int2str(izone))//',:) = '//strStack(int2str(itZonIso(izone,:)))
     665   END DO
     666   DO izone=1,nzone
     667      WRITE(*,*)'strtrac('//TRIM(int2str(izone))//',:) = '//TRIM(strtrac(izone))
     668   END DO
     669   WRITE(*,*) 'ntracisoOR=',ntracisoOR
    1011670#endif 
    1012671
    1013         end subroutine iso_traceurs_init
    1014 
     672END SUBROUTINE iso_traceurs_init
    1015673
    1016674END MODULE isotrac_mod
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/isotrac_routines_mod.F90

    r4143 r4482  
    11131113      USE isotrac_mod, only: use_bassin_atlantic,use_bassin_medit, &
    11141114&       use_bassin_indian,use_bassin_austral,use_bassin_pacific, &
    1115 &       use_bassin_merarabie,use_bassin_golfebengale,use_bassin_indiansud, &
    1116 &       use_bassin_tropics,use_bassin_midlats,use_bassin_hauteslats, &
     1115&       use_bassin_MerArabie,use_bassin_BengalGolf,use_bassin_SouthIndian, &
     1116&       use_bassin_tropics,use_bassin_midlats,use_bassin_HighLats, &
    11171117&       bassin_atlantic,bassin_medit, &
    11181118&       bassin_indian,bassin_austral,bassin_pacific, &
    1119 &       bassin_merarabie,bassin_golfebengale,bassin_indiansud, &
    1120 &       bassin_tropics,bassin_midlats,bassin_hauteslats
     1119&       bassin_MerArabie,bassin_BengalGolf,bassin_SouthIndian, &
     1120&       bassin_tropics,bassin_midlats,bassin_HighLats
    11211121      implicit none
    11221122      ! répond true si lat,lon se trouve dans le bassin numéroté bassin
     
    11371137      write(*,*) 'is_in_basin 84: entree,bassin=',bassin
    11381138#endif
    1139       if ((use_bassin_atlantic.eq.1).and. &
    1140      &           (bassin.eq.bassin_atlantic)) then
     1139      if (use_bassin_atlantic .and. bassin==bassin_atlantic) then
    11411140#ifdef ISOVERIF           
    11421141          write(*,*) 'bassin Atlantique?'
     
    11691168          endif
    11701169
    1171       else if ((use_bassin_medit.eq.1).and. &
    1172      &           (bassin.eq.bassin_medit)) then
     1170      else if (use_bassin_medit .and. bassin==bassin_medit) then
    11731171#ifdef ISOVERIF           
    11741172          write(*,*) 'bassin Medit?'
     
    11831181          endif
    11841182
    1185       else if ((use_bassin_indian.eq.1).and. &
    1186      &           (bassin.eq.bassin_indian)) then
     1183      else if (use_bassin_indian .and. bassin==bassin_indian) then
    11871184#ifdef ISOVERIF           
    11881185          write(*,*) 'bassin indian?'
     
    11991196          endif   
    12001197
    1201       else if ((use_bassin_indiansud.eq.1).and. &
    1202      &           (bassin.eq.bassin_indiansud)) then
     1198      else if (use_bassin_SouthIndian .and. bassin==bassin_SouthIndian) then
    12031199#ifdef ISOVERIF           
    12041200          write(*,*) 'bassin indian hemisphere Sud?'
     
    12091205          endif
    12101206         
    1211       else if ((use_bassin_merarabie.eq.1).and. &
    1212      &           (bassin.eq.bassin_merarabie)) then
     1207      else if (use_bassin_MerArabie .and. bassin==bassin_MerArabie) then
    12131208#ifdef ISOVERIF           
    12141209          write(*,*) 'bassin Mer d''Arabie?'
     
    12191214          endif
    12201215
    1221       else if ((use_bassin_golfebengale.eq.1).and. &
    1222      &           (bassin.eq.bassin_golfebengale)) then
     1216      else if (use_bassin_BengalGolf .and. bassin==bassin_BengalGolf) then
    12231217#ifdef ISOVERIF           
    12241218          write(*,*) 'bassin Golfe du Bengale?'
     
    12291223          endif         
    12301224
    1231       else if ((use_bassin_pacific.eq.1).and. &
    1232      &           (bassin.eq.bassin_pacific)) then
     1225      else if (use_bassin_pacific .and. bassin==bassin_pacific) then
    12331226#ifdef ISOVERIF           
    12341227          write(*,*) 'bassin Pacific?'
     
    12671260          endif
    12681261
    1269       else if ((use_bassin_austral.eq.1).and. &
    1270      &           (bassin.eq.bassin_austral)) then 
     1262      else if (use_bassin_austral .and. bassin==bassin_austral) then 
    12711263#ifdef ISOVERIF           
    12721264          write(*,*) 'bassin austral?'
     
    12771269          endif 
    12781270
    1279       else if ((use_bassin_hauteslats.eq.1).and. &
    1280      &           (bassin.eq.bassin_hauteslats)) then 
     1271      else if (use_bassin_HighLats .and. bassin==bassin_HighLats) then 
    12811272#ifdef ISOVERIF           
    12821273          write(*,*) 'bassin hautes lats?'
     
    12871278          endif
    12881279
    1289       else if ((use_bassin_tropics.eq.1).and. &
    1290      &           (bassin.eq.bassin_tropics)) then 
     1280      else if (use_bassin_tropics .and. bassin==bassin_tropics) then 
    12911281#ifdef ISOVERIF           
    12921282          write(*,*) 'bassin tropics?'
     
    12971287          endif
    12981288
    1299        else if ((use_bassin_midlats.eq.1).and. &
    1300      &           (bassin.eq.bassin_midlats)) then 
     1289       else if (use_bassin_midlats .and. bassin==bassin_midlats) then 
    13011290#ifdef ISOVERIF           
    13021291          write(*,*) 'bassin mid lats?'
     
    13141303          write(*,*) 'bassin_indian=' ,bassin_indian
    13151304          write(*,*) 'bassin_austral=' ,bassin_austral
    1316           write(*,*) 'bassin_merarabie=' ,bassin_merarabie
    1317           write(*,*) 'bassin_golfebengale=' ,bassin_golfebengale
    1318           write(*,*) 'bassin_indiansud=' ,bassin_indiansud
     1305          write(*,*) 'bassin_MerArabie=' ,bassin_MerArabie
     1306          write(*,*) 'bassin_BengalGolf=' ,bassin_BengalGolf
     1307          write(*,*) 'bassin_SouthIndian=' ,bassin_SouthIndian
    13191308          write(*,*) 'use_bassin_atlantic=' ,use_bassin_atlantic 
    13201309          write(*,*) 'use_bassin_medit=' ,use_bassin_medit
    13211310          write(*,*) 'use_bassin_indian=' ,use_bassin_indian
    13221311          write(*,*) 'use_bassin_austral=' ,use_bassin_austral
    1323           write(*,*) 'use_bassin_merarabie=' ,use_bassin_merarabie
    1324           write(*,*) 'use_bassin_golfebengale=' ,use_bassin_golfebengale
    1325           write(*,*) 'use_bassin_indiansud=' ,use_bassin_indiansud
     1312          write(*,*) 'use_bassin_MerArabie=' ,use_bassin_MerArabie
     1313          write(*,*) 'use_bassin_BengalGolf=' ,use_bassin_BengalGolf
     1314          write(*,*) 'use_bassin_SouthIndian=' ,use_bassin_SouthIndian
    13261315          stop
    13271316      endif
     
    23332322        subroutine iso_verif_traceur_jbid_vect(x,n,m)
    23342323        USE isotopes_mod, ONLY: bidouille_anti_divergence,iso_eau,ridicule
    2335         use isotrac_mod, only: ntraceurs_zone=>nzone
     2324        !use isotrac_mod, only: ntraceurs_zone=>nzone
     2325        USE infotrac_phy, ONLY: ntraceurs_zone=>nzone
    23362326        implicit none
    23372327       
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/pbl_surface_mod.F90

    r4143 r4482  
    1414  USE mod_grid_phy_lmdz,   ONLY : klon_glo
    1515  USE ioipsl
    16   USE surface_data,        ONLY : type_ocean, ok_veget
     16  USE surface_data,        ONLY : type_ocean, ok_veget, landice_opt
    1717  USE surf_land_mod,       ONLY : surf_land
    1818  USE surf_landice_mod,    ONLY : surf_landice
     
    2323  USE climb_wind_mod,      ONLY : climb_wind_down, climb_wind_up
    2424  USE coef_diff_turb_mod,  ONLY : coef_diff_turb
     25  USE atke_exchange_coeff_mod, ONLY :  atke_compute_km_kh
    2526  USE ioipsl_getin_p_mod,  ONLY : getin_p
    2627  USE cdrag_mod
     
    406407#endif
    407408    USE ioipsl_getin_p_mod, ONLY : getin_p
    408     use phys_state_var_mod, only: ds_ns, dt_ns, delta_sst, delta_sal, zsig, zmea
    409     use phys_output_var_mod, only: dter, dser, tkt, tks, taur, sss
     409    use phys_state_var_mod, only: ds_ns, dt_ns, delta_sst, delta_sal, dter, &
     410         dser, dt_ds, zsig, zmea
     411    use phys_output_var_mod, only: tkt, tks, taur, sss
    410412#ifdef CPP_XIOS
    411413    USE wxios, ONLY: missing_val
     
    10281030    REAL, DIMENSION(klon)              :: yrmu0
    10291031    ! Martin
    1030 
    1031     REAL, DIMENSION(klon):: ydelta_sst, ydelta_sal, yds_ns, ydt_ns, ydter, ydser, &
    1032          ytkt, ytks, ytaur, ysss
    1033     ! compression of delta_sst, delta_sal, ds_ns, dt_ns, dter, dser, tkt, tks,
    1034     ! taur, sss on ocean points
     1032    REAL, DIMENSION(klon)              :: yri0
     1033    REAL, DIMENSION(klon):: ydelta_sst, ydelta_sal, yds_ns, ydt_ns, ydter, &
     1034         ydser, ydt_ds, ytkt, ytks, ytaur, ysss
     1035    ! compression of delta_sst, delta_sal, ds_ns, dt_ns, dter, dser,
     1036    ! dt_ds, tkt, tks, taur, sss on ocean points
    10351037
    10361038#ifdef ISO
     
    12771279
    12781280    ytke=0.
     1281    yri0(:)=0.
    12791282!FC
    12801283    y_treedrg=0.
     
    18131816             ydelta_sal(:knon) = delta_sal(ni(:knon))
    18141817             ydelta_sst(:knon) = delta_sst(ni(:knon))
     1818             ydter(:knon) = dter(ni(:knon))
     1819             ydser(:knon) = dser(ni(:knon))
     1820             ydt_ds(:knon) = dt_ds(ni(:knon))
    18151821          end if
    18161822         
     
    18441850        CALL cdrag(knon, nsrf, &
    18451851            speed, yt(:,1), yq(:,1), zgeo1, ypaprs(:,1),&
    1846             yts, yqsurf, yz0m, yz0h, &
     1852            yts, yqsurf, yz0m, yz0h, yri0, 0, &
    18471853            ycdragm, ycdragh, zri1, pref )
    18481854
     
    18781884            CALL cdrag(knon, nsrf, &
    18791885            speed_x, yt_x(:,1), yq_x(:,1), zgeo1_x, ypaprs(:,1),&
    1880             yts_x, yqsurf_x, yz0m, yz0h, &
     1886            yts_x, yqsurf_x, yz0m, yz0h, yri0, 0, &
    18811887            ycdragm_x, ycdragh_x, zri1_x, pref_x )
    18821888
     
    19051911        CALL cdrag(knon, nsrf, &
    19061912            speed_w, yt_w(:,1), yq_w(:,1), zgeo1_w, ypaprs(:,1),&
    1907             yts_w, yqsurf_w, yz0m, yz0h, &
     1913            yts_w, yqsurf_w, yz0m, yz0h, yri0, 0, &
    19081914            ycdragm_w, ycdragh_w, zri1_w, pref_w )
    19091915!
     
    19471953
    19481954       ENDIF
     1955
     1956        IF (iflag_pbl>=50) THEN
     1957
     1958        CALL atke_compute_km_kh(knon,klev,yu,yv,yt, &
     1959             ypplay,ypaprs,ytke,ycoefm, ycoefh)
     1960
     1961        ELSE
     1962
     1963
    19491964        CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
    19501965            ypaprs, ypplay, yu, yv, yq, yt, yts, yqsurf, ycdragm, &
     
    19621977       ENDDO
    19631978       ENDIF
     1979
     1980       ENDIF ! iflag_pbl >= 50
     1981
    19641982        IF (prt_level >=10) print *,'coef_diff_turb -> ycoefh ',ycoefh
    19651983!
     
    19771995
    19781996       ENDIF
     1997
     1998        IF (iflag_pbl>=50) THEN
     1999
     2000        CALL atke_compute_km_kh(knon,klev,yu_x,yv_x,yt_x, &
     2001             ypplay,ypaprs,ytke_x,ycoefm_x, ycoefh_x)
     2002
     2003        ELSE
     2004
     2005
    19792006        CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
    19802007            ypaprs, ypplay, yu_x, yv_x, yq_x, yt_x, yts_x, yqsurf_x, ycdragm_x, &
     
    19922019       ENDDO
    19932020       ENDIF
     2021
     2022       ENDIF ! iflag_pbl >= 50
     2023
    19942024        IF (prt_level >=10) print *,'coef_diff_turb -> ycoefh_x ',ycoefh_x
    19952025!
     
    20052035      print *,' args coef_diff_turb: ytke_w ', ytke_w
    20062036       ENDIF
     2037
     2038        IF (iflag_pbl>=50) THEN
     2039
     2040        CALL atke_compute_km_kh(knon,klev,yu_w,yv_w,yt_w, &
     2041             ypplay,ypaprs,ytke_w,ycoefm_w, ycoefh_w)
     2042
     2043        ELSE
     2044
     2045
    20072046        CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
    20082047            ypaprs, ypplay, yu_w, yv_w, yq_w, yt_w, yts_w, yqsurf_w, ycdragm_w, &
     
    20192058       ENDDO
    20202059       ENDIF
     2060
     2061       ENDIF ! iflag_pbl >= 50
     2062
    20212063        IF (prt_level >=10) print *,'coef_diff_turb -> ycoefh_w ',ycoefh_w
    20222064!
     
    23822424       CASE(is_lic)
    23832425          ! Martin
    2384           CALL surf_landice(itap, dtime, knon, ni, &
    2385                rlon, rlat, debut, lafin, &
    2386                yrmu0, ylwdown, yalb, zgeo1, &
    2387                ysolsw, ysollw, yts, ypplay(:,1), &
    2388 !!jyg               ycdragh, ycdragm, yrain_f, ysnow_f, yt(:,1), yq(:,1),&
    2389                ycdragh, ycdragm, yrain_f, ysnow_f, yt1, yq1,&
    2390                AcoefH, AcoefQ, BcoefH, BcoefQ, &
    2391                AcoefU, AcoefV, BcoefU, BcoefV, &
    2392                ypsref, yu1, yv1, ygustiness, yrugoro, pctsrf, &
    2393                ysnow, yqsurf, yqsol, yagesno, &
    2394                ytsoil, yz0m, yz0h, SFRWL, yalb_dir_new, yalb_dif_new, yevap,yfluxsens,yfluxlat, &
    2395                ytsurf_new, y_dflux_t, y_dflux_q, &
    2396                yzmea, yzsig, ycldt, &
    2397                ysnowhgt, yqsnow, ytoice, ysissnow, &
    2398                yalb3_new, yrunoff, &
    2399                y_flux_u1, y_flux_v1 &
    2400 #ifdef ISO
    2401            &    ,yxtrain_f, yxtsnow_f,yxt1,yRland_ice &
    2402            &    ,yxtsnow,yxtsol,yxtevap &
     2426          IF (landice_opt .LT. 2) THEN
     2427             ! Land ice is treated by LMDZ and not by ORCHIDEE
     2428             
     2429             CALL surf_landice(itap, dtime, knon, ni, &
     2430                  rlon, rlat, debut, lafin, &
     2431                  yrmu0, ylwdown, yalb, zgeo1, &
     2432                  ysolsw, ysollw, yts, ypplay(:,1), &
     2433                  !!jyg               ycdragh, ycdragm, yrain_f, ysnow_f, yt(:,1), yq(:,1),&
     2434                  ycdragh, ycdragm, yrain_f, ysnow_f, yt1, yq1,&
     2435                  AcoefH, AcoefQ, BcoefH, BcoefQ, &
     2436                  AcoefU, AcoefV, BcoefU, BcoefV, &
     2437                  ypsref, yu1, yv1, ygustiness, yrugoro, pctsrf, &
     2438                  ysnow, yqsurf, yqsol, yagesno, &
     2439                  ytsoil, yz0m, yz0h, SFRWL, yalb_dir_new, yalb_dif_new, yevap,yfluxsens,yfluxlat, &
     2440                  ytsurf_new, y_dflux_t, y_dflux_q, &
     2441                  yzmea, yzsig, ycldt, &
     2442                  ysnowhgt, yqsnow, ytoice, ysissnow, &
     2443                  yalb3_new, yrunoff, &
     2444                  y_flux_u1, y_flux_v1 &
     2445#ifdef ISO
     2446                  &    ,yxtrain_f, yxtsnow_f,yxt1,yRland_ice &
     2447                  &    ,yxtsnow,yxtsol,yxtevap &
    24032448#endif             
    2404            &    )
    2405 
    2406 !jyg<
    2407 !!          alb3_lic(:)=0.
    2408 !>jyg
    2409           DO j = 1, knon
    2410              i = ni(j)
    2411              alb3_lic(i) = yalb3_new(j)
    2412              snowhgt(i)   = ysnowhgt(j)
    2413              qsnow(i)     = yqsnow(j)
    2414              to_ice(i)    = ytoice(j)
    2415              sissnow(i)   = ysissnow(j)
    2416              runoff(i)    = yrunoff(j)
    2417           ENDDO
    2418           ! Martin
    2419 ! Special DICE MPL 05082013 puis BOMEX MPL 20150410
    2420        IF (ok_prescr_ust) THEN
    2421           DO j=1,knon
    2422           y_flux_u1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yu(j,1)*ypplay(j,1)/RD/yt(j,1)
    2423           y_flux_v1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yv(j,1)*ypplay(j,1)/RD/yt(j,1)
    2424           ENDDO
    2425       ENDIF
    2426 
     2449                  &    )
     2450             
     2451             !jyg<
     2452             !!          alb3_lic(:)=0.
     2453             !>jyg
     2454             DO j = 1, knon
     2455                i = ni(j)
     2456                alb3_lic(i) = yalb3_new(j)
     2457                snowhgt(i)   = ysnowhgt(j)
     2458                qsnow(i)     = yqsnow(j)
     2459                to_ice(i)    = ytoice(j)
     2460                sissnow(i)   = ysissnow(j)
     2461                runoff(i)    = yrunoff(j)
     2462             ENDDO
     2463             ! Martin
     2464             ! Special DICE MPL 05082013 puis BOMEX MPL 20150410
     2465             IF (ok_prescr_ust) THEN
     2466                DO j=1,knon
     2467                   y_flux_u1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yu(j,1)*ypplay(j,1)/RD/yt(j,1)
     2468                   y_flux_v1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yv(j,1)*ypplay(j,1)/RD/yt(j,1)
     2469                ENDDO
     2470             ENDIF
     2471             
    24272472#ifdef ISOVERIF
    2428         do j=1,knon
    2429           do ixt=1,ntraciso
    2430             call iso_verif_noNaN(yxtevap(ixt,j), &
    2431          &      'pbl_surface 1095a: apres surf_landice')
    2432             call iso_verif_noNaN(yxtsol(ixt,j), &
    2433          &      'pbl_surface 1095b: apres surf_landice')
    2434           enddo
    2435         enddo
     2473             do j=1,knon
     2474                do ixt=1,ntraciso
     2475                   call iso_verif_noNaN(yxtevap(ixt,j), &
     2476                        &      'pbl_surface 1095a: apres surf_landice')
     2477                   call iso_verif_noNaN(yxtsol(ixt,j), &
     2478                        &      'pbl_surface 1095b: apres surf_landice')
     2479                enddo
     2480             enddo
    24362481#endif
    24372482#ifdef ISOVERIF
    2438         !write(*,*) 'pbl_surface_mod 1060: sortie surf_landice'
    2439         do j=1,knon
    2440           if (iso_eau.gt.0) then     
    2441                  call iso_verif_egalite(yxtsnow(iso_eau,j), &
    2442      &                  ysnow(j),'pbl_surf_mod 1064')
    2443            endif !if (iso_eau.gt.0) then
    2444         enddo !do i=1,klon
    2445 #endif
    2446          
     2483             !write(*,*) 'pbl_surface_mod 1060: sortie surf_landice'
     2484             do j=1,knon
     2485                if (iso_eau.gt.0) then     
     2486                   call iso_verif_egalite(yxtsnow(iso_eau,j), &
     2487                        &                  ysnow(j),'pbl_surf_mod 1064')
     2488                endif !if (iso_eau.gt.0) then
     2489             enddo !do i=1,klon
     2490#endif
     2491          END IF
    24472492       CASE(is_oce)
    24482493           CALL surf_ocean(rlon, rlat, ysolsw, ysollw, yalb_vis, &
     
    24592504               y_flux_u1, y_flux_v1, ydelta_sst(:knon), ydelta_sal(:knon), &
    24602505               yds_ns(:knon), ydt_ns(:knon), ydter(:knon), ydser(:knon), &
    2461                ytkt(:knon), ytks(:knon), ytaur(:knon), ysss &
     2506               ydt_ds(:knon), ytkt(:knon), ytks(:knon), ytaur(:knon), ysss &
    24622507#ifdef ISO
    24632508         &      ,yxtrain_f, yxtsnow_f,yxt1,Roce, &
     
    33913436          taur(ni(:knon)) = ytaur(:knon)
    33923437          sss(ni(:knon)) = ysss(:knon)
     3438
     3439          if (activate_ocean_skin == 2 .and. type_ocean == "couple") then
     3440             dt_ds = missing_val
     3441             dt_ds(ni(:knon)) = ydt_ds(:knon)
     3442          end if
    33933443       end if
    33943444
     
    41284178
    41294179    USE indice_sol_mod
    4130     use phys_state_var_mod, only: delta_sal, ds_ns, dt_ns, delta_sst
     4180    use phys_state_var_mod, only: delta_sal, ds_ns, dt_ns, delta_sst, dter, &
     4181         dser, dt_ds
    41314182    use config_ocean_skin_m, only: activate_ocean_skin
    41324183
     
    42344285                         delta_sal(i) = 0.
    42354286                         delta_sst(i) = 0.
     4287                         dter(i) = 0.
     4288                         dser(i) = 0.
     4289                         dt_ds(i) = 0.
    42364290                      end if
    42374291                     
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/phyredem.F90

    r4170 r4482  
    3030                                du_gwd_rando, du_gwd_front, u10m, v10m, &
    3131                                treedrg, solswfdiff, delta_sal, ds_ns, dt_ns, &
    32                                 delta_sst, ratqs_inter
     32                                delta_sst, ratqs_inter, dter, dser, dt_ds
    3333#ifdef ISO
    3434  USE phys_state_var_mod, ONLY: xtsol, fxtevap,xtrain_fall, xtsnow_fall,     &
     
    3939  USE iostart, ONLY: open_restartphy, close_restartphy, enddef_restartphy, put_field, put_var
    4040  USE traclmdz_mod, ONLY : traclmdz_to_restart
    41   USE infotrac_phy, ONLY: types_trac, nqtot, tracers, nbtr, niso
     41  USE infotrac_phy, ONLY: type_trac, nqtot, tracers, nbtr, niso
    4242#ifdef ISO
    4343#ifdef ISOVERIF
     
    4545#endif
    4646#endif
    47   USE carbon_cycle_mod, ONLY : carbon_cycle_cpl, co2_send
     47  USE carbon_cycle_mod, ONLY : carbon_cycle_cpl, co2_send, carbon_cycle_rad, RCO2_glo
    4848  USE indice_sol_mod, ONLY: nbsrf, is_oce, is_sic, is_ter, is_lic, epsfra
    4949  USE surface_data, ONLY: type_ocean, version_ocean
     
    130130
    131131  ! co2_ppm0 : initial value of atmospheric CO2
    132   tab_cntrl(16) = co2_ppm0
     132  ! tab_cntrl(16) = co2_ppm0
     133
     134  !  PC -- initial value of RCO2 for the radiation scheme
     135  !  tab_cntrl(17) = co2_ppm * 1.0e-06 * RMCO2 / RMD
     136  IF (carbon_cycle_rad) tab_cntrl(17) = RCO2_glo
     137  !PRINT*, "PC : phyredem RCO2_glo =",RCO2_glo
    133138
    134139  DO pass=1,2   ! pass=1 netcdf definition ; pass=2 netcdf write
     
    171176    CALL put_field(pass,"FSIC", "fraction glace mer", pctsrf(:, is_sic))
    172177
    173     IF(nbsrf>99) THEN
    174       PRINT*, "Trop de sous-mailles";  CALL abort_physic("phyredem", "", 1)
    175     END IF
    176     IF(nsoilmx>99) THEN
    177       PRINT*, "Trop de sous-surfaces"; CALL abort_physic("phyredem", "", 1)
    178     END IF
    179     IF(nsw>99) THEN
    180       PRINT*, "Trop de bandes"; CALL abort_physic("phyredem", "", 1)
    181     END IF
     178    IF(nbsrf  >99) CALL abort_physic("phyredem", "Trop de sous-mailles", 1)
     179    IF(nsoilmx>99) CALL abort_physic("phyredem", "Trop de sous-mailles", 1)
     180    IF(nsw    >99) CALL abort_physic("phyredem", "Trop de bandes", 1)
    182181
    183182!    Surface variables
     
    345344
    346345
    347     IF (ANY(types_trac == 'co2i') .OR. ANY(types_trac == 'inco')) THEN
     346    IF (ANY(type_trac == ['co2i','inco'])) THEN
    348347       IF (carbon_cycle_cpl) THEN
    349348          IF (.NOT. ALLOCATED(co2_send)) THEN
     
    356355
    357356    ! trs from traclmdz_mod
    358     ELSE IF (ANY(types_trac == 'lmdz')) THEN
     357    ELSE IF (type_trac == 'lmdz') THEN
    359358       CALL traclmdz_to_restart(trs)
    360359       it = 0
     
    394393          CALL put_field(pass, "delta_SST", &
    395394               "ocean-air interface temperature minus bulk SST", delta_sst)
     395          CALL put_field(pass, "dter", &
     396               "ocean-air interface temperature minus subskin temperature", &
     397               dter)
     398          CALL put_field(pass, "dser", &
     399               "ocean-air interface salinity minus subskin salinity", dser)
     400          CALL put_field(pass, "dt_ds", &
     401               "(tks / tkt) * dTer", dt_ds)
    396402       end if
    397403       
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/phys_output_var_mod.F90

    r3940 r4482  
    135135  ! Ocean-atmosphere interface, subskin ocean and near-surface ocean:
    136136 
    137   REAL, ALLOCATABLE, SAVE:: dter(:)
    138   ! Temperature variation in the diffusive microlayer, that is
    139   ! ocean-air interface temperature minus subskin temperature. In K.
    140      
    141   REAL, SAVE, ALLOCATABLE:: dser(:)
    142   ! Temperature variation in the diffusive microlayer, that is
    143   ! subskin temperature minus ocean-air interface temperature. In K.
    144 
    145   REAL, SAVE, ALLOCATABLE:: tkt(:)
     137  REAL, SAVE, ALLOCATABLE:: tkt(:) ! (klon)
    146138  ! épaisseur (m) de la couche de diffusion thermique (microlayer)
    147139  ! cool skin thickness
    148140
    149   REAL, SAVE, ALLOCATABLE:: tks(:)
     141  REAL, SAVE, ALLOCATABLE:: tks(:) ! (klon)
    150142  ! épaisseur (m) de la couche de diffusion de masse (microlayer)
    151143 
    152   REAL, SAVE, ALLOCATABLE:: taur(:) ! momentum flux due to rain, in Pa
    153 
    154   REAL, SAVE, ALLOCATABLE:: sss(:)
     144  REAL, SAVE, ALLOCATABLE:: taur(:) ! (klon) momentum flux due to rain, in Pa
     145
     146  REAL, SAVE, ALLOCATABLE:: sss(:) ! (klon)
    155147  ! bulk salinity of the surface layer of the ocean, in ppt
    156148 
    157   !$OMP THREADPRIVATE(dter, dser, tkt, tks, taur, sss)
     149  !$OMP THREADPRIVATE(tkt, tks, taur, sss)
    158150
    159151CONTAINS
     
    216208    IF (ok_gwd_rando) allocate(zustr_gwd_rando(klon), zvstr_gwd_rando(klon))
    217209
    218     if (activate_ocean_skin >= 1) allocate(dter(klon), dser(klon), tkt(klon), &
    219          tks(klon), taur(klon), sss(klon))
     210    if (activate_ocean_skin >= 1) allocate(tkt(klon), tks(klon), taur(klon), &
     211         sss(klon))
    220212
    221213  END SUBROUTINE phys_output_var_init
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/physiq_mod.F90

    r4170 r4482  
    3939    USE ioipsl_getin_p_mod, ONLY : getin_p
    4040    USE indice_sol_mod
    41     USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, types_trac, nqCO2
     41    USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, type_trac
    4242    USE readTracFiles_mod, ONLY: addPhase
    43     USE strings_mod,  ONLY: strIdx, strStack, int2str
     43    USE strings_mod,  ONLY: strIdx
    4444    USE iophy
    4545    USE limit_read_mod, ONLY : init_limit_read
     
    5353    USE pbl_surface_mod, ONLY : pbl_surface
    5454    USE phyaqua_mod, only: zenang_an
     55    USE phyetat0_mod, only: phyetat0
    5556    USE phystokenc_mod, ONLY: offline, phystokenc
    5657    USE phys_cal_mod, only: year_len, mth_len, days_elapsed, jh_1jan, &
    57          year_cur, mth_cur,jD_cur, jH_cur, jD_ref, day_cur, hour
     58         year_cur, mth_cur,jD_cur, jH_cur, jD_ref, day_cur, hour, calend
    5859!!  USE phys_local_var_mod, ONLY : a long list of variables
    5960!!              ==> see below, after "CPP Keys" section
     
    6869    USE regr_horiz_time_climoz_m, ONLY: regr_horiz_time_climoz
    6970    USE regr_pr_time_av_m, only: regr_pr_time_av
    70     USE surface_data,     ONLY : type_ocean, ok_veget, landice_opt
    71     USE time_phylmdz_mod, only: annee_ref, current_time, day_ini, day_ref, &
    72           day_step_phy, itau_phy, pdtphys, raz_date, start_time, update_time
     71    USE surface_data,     ONLY : type_ocean, ok_veget
     72    USE time_phylmdz_mod, only: current_time, itau_phy, pdtphys, raz_date, update_time
    7373    USE tracinca_mod, ONLY: config_inca
    7474    USE tropopause_m,     ONLY: dyn_tropopause
    7575    USE ice_sursat_mod,  ONLY: flight_init, airplane
    7676    USE vampir
    77     USE VERTICAL_LAYERS_MOD, ONLY: aps,bps, ap, bp
    7877    USE write_field_phy
     78#ifdef CPP_XIOS
     79    USE wxios, ONLY: g_ctx, wxios_set_context
     80#endif
     81    USE atke_turbulence_ini_mod, ONLY : atke_ini
     82    USE lscp_ini_mod, ONLY : lscp_ini
    7983    USE lscp_mod, ONLY : lscp
     84    USE wake_ini_mod, ONLY : wake_ini
     85    USE yamada_ini_mod, ONLY : yamada_ini
    8086    USE thermcell_ini_mod, ONLY : thermcell_ini
    8187
     
    96102
    97103
     104#ifdef INCA
     105    USE geometry_mod,      ONLY: longitude, latitude, boundslon, boundslat, ind_cell_glo
     106    USE time_phylmdz_mod,  ONLY: ndays
     107    USE infotrac_phy,      ONLY: nqCO2
     108#endif
    98109#ifdef REPROBUS
    99     USE CHEM_REP, ONLY : Init_chem_rep_xjour, &
    100          d_q_rep,d_ql_rep,d_qi_rep,ptrop,ttrop, &
    101          ztrop, gravit,itroprep, Z1,Z2,fac,B
     110    USE chem_rep, ONLY: Init_chem_rep_xjour, d_q_rep, d_ql_rep, d_qi_rep, &
     111                        ptrop, ttrop, ztrop, gravit, itroprep, Z1, Z2, fac, B
     112#endif
     113#if defined INCA || defined REPROBUS
     114    USE time_phylmdz_mod,    ONLY: annee_ref, day_ini, day_ref, start_time
     115    USE vertical_layers_mod, ONLY: aps, bps, ap, bp
    102116#endif
    103117
     
    105119#ifdef CPP_RRTM
    106120    USE YOERAD, ONLY : NRADLP
    107     USE YOESW, ONLY : RSUN
     121!    USE YOESW, ONLY : RSUN
    108122#endif
    109123
     
    115129
    116130#ifdef CPP_XIOS
    117     USE xios, ONLY: xios_update_calendar, xios_context_finalize, &
    118             xios_get_field_attr, xios_field_is_active
     131    USE xios, ONLY: xios_update_calendar, xios_context_finalize
     132    USE xios, ONLY: xios_get_field_attr, xios_field_is_active, xios_context
     133    USE xios, ONLY: xios_set_current_context
    119134    USE wxios, ONLY: missing_val, missing_val_omp
    120135#endif
     
    179194       d_t_ajsb,d_q_ajsb, &
    180195       d_t_ajs,d_q_ajs,d_u_ajs,d_v_ajs, &
    181        d_t_ajs_w,d_q_ajs_w, &
    182        d_t_ajs_x,d_q_ajs_x, &
     196!       d_t_ajs_w,d_q_ajs_w, &
     197!       d_t_ajs_x,d_q_ajs_x, &
    183198       !
    184199       d_t_eva,d_q_eva,d_ql_eva,d_qi_eva, &
     
    193208       d_ts, &
    194209       !
    195        d_t_oli,d_u_oli,d_v_oli, &
     210!       d_t_oli,d_u_oli,d_v_oli, &
    196211       d_t_oro,d_u_oro,d_v_oro, &
    197212       d_t_oro_gw,d_u_oro_gw,d_v_oro_gw, &
     
    312327       sij, &
    313328       !
     329       rneblsvol, &
     330       zqsatl, zqsats, &
     331       qclr, qcld, qss, qvc, rnebclr, rnebss, gamma_ss, &
     332       Tcontr, qcontr, qcontr2, fcontrN, fcontrP, &
    314333       cldemi,  &
    315334       cldfra, cldtau, fiwc,  &
     
    365384#endif
    366385       !
    367 
    368386
    369387    IMPLICIT NONE
     
    578596    !
    579597    !
    580     INTEGER debug
    581598    INTEGER n
    582599    !ym      INTEGER npoints
     
    635652    ! Upmost level reached by deep convection and related variable (jyg)
    636653    !
    637     INTEGER izero
     654!    INTEGER izero
    638655    INTEGER k_upper_cv
    639656    !------------------------------------------------------------------
     
    813830    REAL beta_prec_fisrt(klon,klev) ! taux de conv de l'eau cond (fisrt)
    814831    ! RomP <<<
    815     REAL          :: calday
    816832
    817833    !IM cf FH pour Tiedtke 080604
     
    869885    !C      EXTERNAL o3cm      ! initialiser l'ozone
    870886    EXTERNAL orbite    ! calculer l'orbite terrestre
    871     EXTERNAL phyetat0  ! lire l'etat initial de la physique
    872887    EXTERNAL phyredem  ! ecrire l'etat de redemarrage de la physique
    873888    EXTERNAL suphel    ! initialiser certaines constantes
     
    934949    real zqsat(klon,klev)
    935950    !
    936     INTEGER i, k, iq, j, nsrf, ll, l, itr
     951    INTEGER i, k, iq, nsrf, l, itr
    937952#ifdef ISO
    938953    real zxt_apres(ntraciso,klon)
     
    10561071
    10571072    REAL picefra(klon,klev)
     1073    REAL zrel_oro(klon)
    10581074    !IM cf. AM 081204 END
    10591075    !
     
    12651281
    12661282#ifdef INCA
     1283    REAL :: calday, zxsnow_dummy(klon)
    12671284    ! set de variables utilisees pour l'initialisation des valeurs provenant de INCA
    12681285    REAL, DIMENSION(klon,klev,naero_grp,nbands) :: init_tauinca
     
    13131330    phys_tstep=NINT(pdtphys)
    13141331#ifdef CPP_XIOS
    1315     IF (.NOT. debut .AND. is_omp_master) CALL xios_update_calendar(itap+1)
     1332! switch to XIOS LMDZ physics context
     1333    IF (.NOT. debut .AND. is_omp_master) THEN
     1334       CALL wxios_set_context()
     1335       CALL xios_update_calendar(itap+1)
     1336    ENDIF
    13161337#endif
    13171338
     
    14031424          WRITE (lunout, *) ' ok_ice_sursat=y requires 4 H2O tracers ', &
    14041425               '(H2O_g, H2O_l, H2O_s, H2O_r) but nqo=', nqo, '. Might as well stop here.'
     1426          abort_message='see above'
     1427          CALL abort_physic(modname,abort_message,1)
     1428       ENDIF
     1429
     1430       IF (ok_plane_h2o.AND..NOT.ok_ice_sursat) THEN
     1431          WRITE (lunout, *) ' ok_plane_h2o=y requires ok_ice_sursat=y '
     1432          abort_message='see above'
     1433          CALL abort_physic(modname,abort_message,1)
     1434       ENDIF
     1435
     1436       IF (ok_plane_contrail.AND..NOT.ok_ice_sursat) THEN
     1437          WRITE (lunout, *) ' ok_plane_contrail=y requires ok_ice_sursat=y '
    14051438          abort_message='see above'
    14061439          CALL abort_physic(modname,abort_message,1)
     
    15391572       tau_overturning_th(:)=0.
    15401573
    1541        IF (ANY(types_trac == 'inca') .OR. ANY(types_trac == 'inco')) THEN
     1574       IF (ANY(type_trac == ['inca','inco'])) THEN
    15421575          ! jg : initialisation jusqu'au ces variables sont dans restart
    15431576          ccm(:,:,:) = 0.
     
    18471880!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    18481881       ! Nouvelle initialisation pour le rayonnement RRTM
    1849        !
    18501882!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    18511883
    18521884       CALL iniradia(klon,klev,paprs(1,1:klev+1))
    18531885
    1854 
    18551886!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1887       CALL wake_ini(rg,rd,rv,prt_level)
     1888       CALL yamada_ini(klon,lunout,prt_level)
     1889       CALL atke_ini(prt_level, lunout, RG, RD, RPI, RCPD)
    18561890       CALL thermcell_ini(iflag_thermals,prt_level,tau_thermals,lunout, &
    18571891   &    RG,RD,RCPD,RKAPPA,RLVTT,RETV)
     1892       IF (ok_new_lscp) then
     1893           CALL lscp_ini(pdtphys,ok_ice_sursat)
     1894       endif
     1895
     1896!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1897
    18581898       !
    18591899!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     
    21512191       !c         ENDDO
    21522192       !
    2153        IF (ANY(types_trac == 'inca') .OR. ANY(types_trac == 'inco')) THEN ! ModThL
     2193       IF (ANY(type_trac == ['inca','inco'])) THEN ! ModThL
    21542194#ifdef INCA
    21552195          CALL VTe(VTphysiq)
     
    21582198          WRITE(lunout,*) 'initial time chemini', days_elapsed, calday
    21592199
    2160           CALL chemini(  &
    2161                rg, &
    2162                ra, &
    2163                cell_area, &
    2164                latitude_deg, &
    2165                longitude_deg, &
    2166                presnivs, &
    2167                calday, &
    2168                klon, &
    2169                nqtot, &
    2170                nqo+nqCO2, &
    2171                pdtphys, &
    2172                annee_ref, &
    2173                year_cur, &
    2174                day_ref,  &
    2175                day_ini, &
    2176                start_time, &
    2177                itau_phy, &
    2178                date0, &
    2179                io_lon, &
    2180                io_lat, &
    2181                chemistry_couple, &
    2182                init_source, &
    2183                init_tauinca, &
    2184                init_pizinca, &
    2185                init_cginca, &
    2186                init_ccminca)
     2200          call init_const_lmdz( &
     2201          ndays, nbsrf, is_oce,is_sic, is_ter,is_lic, calend, &
     2202          config_inca)
     2203
     2204          CALL init_inca_geometry( &
     2205               longitude, latitude, &
     2206               boundslon, boundslat, &
     2207               cell_area, ind_cell_glo)
     2208
     2209          if (grid_type==unstructured) THEN
     2210             CALL chemini(  pplay, &
     2211                  nbp_lon, nbp_lat, &
     2212                  latitude_deg, &
     2213                  longitude_deg, &
     2214                  presnivs, &
     2215                  calday, &
     2216                  klon, &
     2217                  nqtot, &
     2218                  nqo+nqCO2, &
     2219                  pdtphys, &
     2220                  annee_ref, &
     2221                  year_cur, &
     2222                  day_ref,  &
     2223                  day_ini, &
     2224                  start_time, &
     2225                  itau_phy, &
     2226                  date0, &
     2227                  chemistry_couple, &
     2228                  init_source, &
     2229                  init_tauinca, &
     2230                  init_pizinca, &
     2231                  init_cginca, &
     2232                  init_ccminca)
     2233          ELSE
     2234             CALL chemini(  pplay, &
     2235                  nbp_lon, nbp_lat, &
     2236                  latitude_deg, &
     2237                  longitude_deg, &
     2238                  presnivs, &
     2239                  calday, &
     2240                  klon, &
     2241                  nqtot, &
     2242                  nqo+nqCO2, &
     2243                  pdtphys, &
     2244                  annee_ref, &
     2245                  year_cur, &
     2246                  day_ref,  &
     2247                  day_ini, &
     2248                  start_time, &
     2249                  itau_phy, &
     2250                  date0, &
     2251                  chemistry_couple, &
     2252                  init_source, &
     2253                  init_tauinca, &
     2254                  init_pizinca, &
     2255                  init_cginca, &
     2256                  init_ccminca, &
     2257                  io_lon, &
     2258                  io_lat)
     2259          ENDIF
    21872260
    21882261
     
    22002273       ENDIF
    22012274       !
    2202        IF (ANY(types_trac == 'repr')) THEN
     2275       IF (type_trac == 'repr') THEN
    22032276#ifdef REPROBUS
    22042277          CALL chemini_rep(  &
     
    22862359
    22872360
     2361
    22882362    ENDIF
    22892363    !
     
    23162390
    23172391    ! Update time and other variables in Reprobus
    2318     IF (ANY(types_trac == 'repr')) THEN
     2392    IF (type_trac == 'repr') THEN
    23192393#ifdef REPROBUS
    23202394       CALL Init_chem_rep_xjour(jD_cur-jD_ref+day_ref)
     
    23902464          ql_seri(i,k) = qx(i,k,iliq)
    23912465          !CR: ATTENTION, on rajoute la variable glace
    2392           IF (nqo.eq.2) THEN
     2466          IF (nqo.EQ.2) THEN             !--vapour and liquid only
    23932467             qs_seri(i,k) = 0.
    2394           ELSE IF (nqo.eq.3) THEN
     2468             rneb_seri(i,k) = 0.
     2469          ELSE IF (nqo.EQ.3) THEN        !--vapour, liquid and ice
    23952470             qs_seri(i,k) = qx(i,k,isol)
    2396           ELSE IF (nqo.eq.4) THEN
     2471             rneb_seri(i,k) = 0.
     2472          ELSE IF (nqo.EQ.4) THEN        !--vapour, liquid, ice and rneb
    23972473             qs_seri(i,k) = qx(i,k,isol)
    23982474             rneb_seri(i,k) = qx(i,k,irneb)
     
    26062682             DO i = 1, klon
    26072683              if (q_seri(i,k).gt.3e-3) then
    2608               call iso_verif_positif(deltaD(xt_seri(iso_eau,i,k) &
     2684              call iso_verif_positif(deltaD(xt_seri(iso_HDO,i,k) &
    26092685     &           /q_seri(i,k))+400.0,'physiq 2045a')
    2610               call iso_verif_positif(deltaD(xt_ancien(iso_eau,i,k) &
     2686              call iso_verif_positif(deltaD(xt_ancien(iso_HDO,i,k) &
    26112687     &            /q_ancien(i,k))+400.0,'physiq 2045b')
    26122688!              call iso_verif_egalite_choix(d_xt_dyn(iso_hdo,i,k) &
     
    26582734       ! !! RomP >>>   td dyn traceur
    26592735       IF (nqtot > nqo) d_tr_dyn(:,:,:)= 0.0
     2736       ! !! RomP <<<
    26602737       d_rneb_dyn(:,:)=0.0
    2661        ! !! RomP <<<
    26622738       ancien_ok = .TRUE.
    26632739    ENDIF
     
    35263602           endif ! if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then
    35273603           if ((iso_O17.gt.0).and.(iso_O18.gt.0)) then
    3528              if ((q_seri(i,k).gt.ridicule).and.(l.lt.nlevmaxO17)) then
     3604             if ((q_seri(i,k).gt.ridicule).and.(k.lt.nlevmaxO17)) then
    35293605              call iso_verif_aberrant_o17(xt_seri(iso_o17,i,k) &
    35303606     &           /q_seri(i,k),xt_seri(iso_o18,i,k) &
     
    35823658          !
    35833659          !>jyg
    3584           IF (ANY(types_trac == 'repr')) THEN
     3660          IF (type_trac == 'repr') THEN
    35853661             nbtr_tmp=ntra
    35863662          ELSE
     
    46084684    IF (ok_new_lscp) THEN
    46094685
    4610     CALL lscp(phys_tstep,missing_val,paprs,pplay, &
     4686    !--mise à jour de flight_m et flight_h2o dans leur module
     4687    IF (ok_plane_h2o .OR. ok_plane_contrail) THEN
     4688      CALL airplane(debut,pphis,pplay,paprs,t_seri)
     4689    ENDIF
     4690
     4691    CALL lscp(klon,klev,phys_tstep,missing_val,paprs,pplay, &
    46114692         t_seri, q_seri,ptconv,ratqs, &
    4612          d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, rneb_seri, &
     4693         d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, rneblsvol, rneb_seri, &
    46134694         cldliq, picefra, rain_lsc, snow_lsc, &
    4614          pfrac_impa, pfrac_nucl, pfrac_1nucl, &
    46154695         frac_impa, frac_nucl, beta_prec_fisrt, &
    46164696         prfl, psfl, rhcl,  &
    46174697         zqasc, fraca,ztv,zpspsk,ztla,zthl,iflag_cld_th, &
    4618          iflag_ice_thermo, ok_ice_sursat)
     4698         iflag_ice_thermo, ok_ice_sursat, zqsatl, zqsats, &
     4699         qclr, qcld, qss, qvc, rnebclr, rnebss, gamma_ss, &
     4700         Tcontr, qcontr, qcontr2, fcontrN, fcontrP )
    46194701
    46204702    ELSE
     4703
    46214704    CALL fisrtilp(phys_tstep,paprs,pplay, &
    46224705         t_seri, q_seri,ptconv,ratqs, &
     
    46904773       ENDDO
    46914774    ENDDO
    4692     IF (nqo==3) THEN
     4775    IF (nqo >= 3) THEN
    46934776    DO k = 1, klev
    46944777       DO i = 1, klon
     
    51125195    ENDDO
    51135196
    5114     IF (ANY(types_trac == 'inca') .OR. ANY(types_trac == 'inco')) THEN ! ModThL
     5197    IF (ANY(type_trac == ['inca','inco'])) THEN ! ModThL
    51155198#ifdef INCA
    51165199       CALL VTe(VTphysiq)
     
    51685251#endif
    51695252    ENDIF !type_trac = inca or inco
    5170     IF (ANY(types_trac == 'repr')) THEN
     5253    IF (type_trac == 'repr') THEN
    51715254#ifdef REPROBUS
    51725255    !CALL chemtime_rep(itap+itau_phy-1, date0, dtime, itap)
     
    54975580          !
    54985581          !--interactive CO2 in ppm from carbon cycle
    5499           IF (carbon_cycle_rad.AND..NOT.debut) THEN
    5500             RCO2=RCO2_glo
    5501           ENDIF
     5582          IF (carbon_cycle_rad) RCO2=RCO2_glo
    55025583          !
    55035584          IF (prt_level .GE.10) THEN
     
    55555636 
    55565637#ifndef CPP_XIOS
    5557           !--OB 30/05/2016 modified 21/10/2016
    5558           !--here we return swaero_diag and dryaod_diag to FALSE
    5559           !--and histdef will switch it back to TRUE if necessary
    5560           !--this is necessary to get the right swaero at first step
    5561           !--but only in the case of no XIOS as XIOS is covered elsewhere
    5562           IF (debut) swaerofree_diag = .FALSE.
    5563           IF (debut) swaero_diag = .FALSE.
    5564           IF (debut) dryaod_diag = .FALSE.
    5565           !--IM 15/09/2017 here we return ok_4xCO2atm to FALSE
    5566           !--as for swaero_diag, see above
    5567           IF (debut) ok_4xCO2atm = .FALSE.
    5568 
    55695638          !
    55705639          !IM 2eme calcul radiatif pour le cas perturbe ou au moins un
     
    57345803       DO i=1,klon
    57355804          itest(i)=0
    5736           !        IF ((zstd(i).gt.10.0)) THEN
    5737           IF (((zpic(i)-zmea(i)).GT.100.).AND.(zstd(i).GT.10.0)) THEN
     5805          zrel_oro(i)=zstd(i)/(max(zsig(i),1.E-8)*sqrt(cell_area(i)))
     5806          !zrel_oro: relative mountain height wrt relief explained by mean slope
     5807          ! -> condition on zrel_oro can deactivate the drag on tilted planar terrains
     5808          !    such as ice sheets (work by V. Wiener)
     5809          ! zpmm_orodr_t and zstd_orodr_t are activation thresholds set by F. Lott to
     5810          ! earn computation time but they are not physical.
     5811          IF (((zpic(i)-zmea(i)).GT.zpmm_orodr_t).AND.(zstd(i).GT.zstd_orodr_t).AND.(zrel_oro(i).LE.zrel_oro_t)) THEN
    57385812             itest(i)=1
    57395813             igwd=igwd+1
     
    57885862       DO i=1,klon
    57895863          itest(i)=0
    5790           IF ((zpic(i)-zmea(i)).GT.100.) THEN
     5864          !zrel_oro: relative mountain height wrt relief explained by mean slope
     5865          ! -> condition on zrel_oro can deactivate the lifting on tilted planar terrains
     5866          !    such as ice sheets (work by V. Wiener)
     5867          zrel_oro(i)=zstd(i)/(max(zsig(i),1.E-8)*sqrt(cell_area(i)))
     5868          IF (((zpic(i)-zmea(i)).GT.zpmm_orolf_t).AND.(zrel_oro(i).LE.zrel_oro_t)) THEN
    57915869             itest(i)=1
    57925870             igwd=igwd+1
     
    60046082! car on peut s'attendre a ce que les petites echelles produisent aussi de la TKE
    60056083! Mais attention, cela ne va pas dans le sens de la conservation de l'energie!
    6006           IF (zstd(i).GT.1.0) THEN
     6084          IF ((zstd(i).GT.1.0).AND.(zrel_oro(i).LE.zrel_oro_t)) THEN
    60076085             itest(i)=1
    60086086             igwd=igwd+1
     
    60166094       DO i=1,klon
    60176095          itest(i)=0
    6018         IF (((zpic(i)-zmea(i)).GT.100.).AND.(zstd(i).GT.10.0)) THEN
     6096        IF (((zpic(i)-zmea(i)).GT.zpmm_orodr_t).AND.(zstd(i).GT.zstd_orodr_t).AND.(zrel_oro(i).LE.zrel_oro_t)) THEN
    60196097             itest(i)=1
    60206098             igwd=igwd+1
     
    62646342    !
    62656343
    6266     IF (ANY(types_trac=='repr')) THEN
     6344    IF (type_trac=='repr') THEN
    62676345!MM pas d'impact, car on recupere q_seri,tr_seri,t_seri via phys_local_var_mod
    62686346!MM                               dans Reprobus
     
    62756353    ELSE
    62766354       sh_in(:,:) = qx(:,:,ivap)
    6277        IF (nqo .EQ. 3) THEN
     6355       IF (nqo >= 3) THEN
    62786356          ch_in(:,:) = qx(:,:,iliq) + qx(:,:,isol)
    62796357       ELSE
     
    63626440    ! Calculer le transport de l'eau et de l'energie (diagnostique)
    63636441    !
    6364     CALL transp (paprs,zxtsol, &
    6365          t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, zphi, &
    6366          ve, vq, ue, uq, vwat, uwat)
     6442    CALL transp (paprs,zxtsol, t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, zphi, &
     6443                 ue, ve, uq, vq, uwat, vwat)
    63676444    !
    63686445    !IM global posePB BEG
    63696446    IF(1.EQ.0) THEN
    63706447       !
    6371        CALL transp_lay (paprs,zxtsol, &
    6372             t_seri, q_seri, u_seri, v_seri, zphi, &
     6448       CALL transp_lay (paprs,zxtsol, t_seri, q_seri, u_seri, v_seri, zphi, &
    63736449            ve_lay, vq_lay, ue_lay, uq_lay)
    63746450       !
    63756451    ENDIF !(1.EQ.0) THEN
    63766452    !IM global posePB END
     6453    !
    63776454    ! Accumuler les variables a stocker dans les fichiers histoire:
    63786455    !
     
    63856462    d_t_ec(:,:)=0.
    63866463    forall (k=1: nbp_lev) exner(:, k) = (pplay(:, k)/paprs(:,1))**RKAPPA
    6387     CALL ener_conserv(klon,klev,pdtphys,u,v,t,qx(:,:,ivap),qx(:,:,iliq),qx(:,:,isol), &
     6464    CALL ener_conserv(klon,klev,pdtphys,u,v,t,qx,ivap,iliq,isol, &
    63886465         u_seri,v_seri,t_seri,q_seri,ql_seri,qs_seri,pbl_tke(:,:,is_ave)-tke0(:,:), &
    63896466         zmasse,exner,d_t_ec)
     
    64306507#endif
    64316508    !
    6432     IF (ANY(types_trac == 'inca') .OR. ANY(types_trac == 'inco')) THEN
     6509    IF (ANY(type_trac == ['inca','inco'])) THEN
    64336510#ifdef INCA
    64346511       CALL VTe(VTphysiq)
     
    64476524            pphis, &
    64486525            zx_rh, &
    6449             aps, bps, ap, bp)
     6526            aps, bps, ap, bp, lafin)
    64506527
    64516528       CALL VTe(VTinca)
     
    64546531    ENDIF
    64556532
     6533    IF (type_trac == 'repr') THEN
     6534#ifdef REPROBUS
     6535        CALL coord_hyb_rep(paprs, pplay, aps, bps, ap, bp, cell_area)
     6536#endif
     6537    ENDIF
    64566538
    64576539    !
     
    64776559          d_qx(i,k,iliq) = ( ql_seri(i,k) - qx(i,k,iliq) ) / phys_tstep
    64786560          !CR: on ajoute le contenu en glace
    6479           IF (nqo.ge.3) THEN
     6561          IF (nqo >= 3) THEN
    64806562             d_qx(i,k,isol) = ( qs_seri(i,k) - qx(i,k,isol) ) / phys_tstep
    64816563          ENDIF
    64826564          !--ice_sursat: nqo=4, on ajoute rneb
    6483           IF (nqo.eq.4) THEN
     6565          IF (nqo == 4) THEN
    64846566             d_qx(i,k,irneb) = ( rneb_seri(i,k) - qx(i,k,irneb) ) / phys_tstep
    64856567          ENDIF
     
    65036585      enddo ! DO k = 1, klev
    65046586    enddo !do ixt=1,ntraciso
    6505 !#ifdef ISOVERIF
    6506 !        write(*,*) 'physiq 6120: d_qx(1,1,:)=',d_qx(1,1,:)
    6507 !        write(*,*) 'qx(1,1,:)=',qx(1,1,:)
    6508 !        write(*,*) 'xt_seri(:,1,1)=',xt_seri(:,1,1)
    6509 !#endif
    6510 #endif
    6511 ! #ifdef ISO
     6587#endif
    65126588    ! DC: All iterations are cycled if nqtot==nqo, so no nqtot>nqo condition required
    65136589    itr = 0
     
    65696645    ql_ancien(:,:) = ql_seri(:,:)
    65706646    qs_ancien(:,:) = qs_seri(:,:)
     6647    rneb_ancien(:,:) = rneb_seri(:,:)
    65716648#ifdef ISO
    65726649    xt_ancien(:,:,:)=xt_seri(:,:,:)
     
    68286905#ifdef CPP_XIOS
    68296906       IF (is_omp_master) CALL xios_context_finalize
     6907
     6908#ifdef INCA
     6909       if (type_trac == 'inca') then
     6910          IF (is_omp_master .and. grid_type==unstructured) THEN
     6911             CALL finalize_inca
     6912          ENDIF
     6913       endif
     6914#endif
     6915
    68306916#endif
    68316917       WRITE(lunout,*) ' physiq fin, nombre de steps ou cvpas = 1 : ', Ncvpaseq1
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/surf_land_mod.F90

    r4143 r4482  
    3030    USE dimphy
    3131    USE surface_data, ONLY    : ok_veget
    32 ! >> PC
    3332    USE carbon_cycle_mod
    34 ! << PC
    3533
    3634    ! See comments in each module surf_land_orchidee_xxx for compatiblity with ORCHIDEE
     
    5149    USE surf_land_orchidee_nounstruct_mod
    5250#else
     51#if ORCHIDEE_NOLIC
     52    ! Compilation with cpp key ORCHIDEE_NOLIC
     53    USE surf_land_orchidee_nolic_mod
     54#else
     55    ! Default version#else
    5356    USE surf_land_orchidee_mod
     57#endif
    5458#endif
    5559#endif
     
    6771#endif
    6872#endif
    69 
    70 ! >> PC
     73   
    7174    USE print_control_mod, ONLY: lunout
    72 ! << PC
    7375
    7476    INCLUDE "dimsoil.h"
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/surf_landice_mod.F90

    r4143 r4482  
    144144    INTEGER                  :: i,j,nt
    145145    REAL, DIMENSION(klon)    :: fqfonte,ffonte
     146    REAL, DIMENSION(klon)    :: run_off_lic_frac
    146147#ifdef ISO       
    147148      real, parameter :: t_coup = 273.15
     
    484485! Send run-off on land-ice to coupler if coupled ocean.
    485486! run_off_lic has been calculated in fonte_neige or surf_inlandsis
    486 !
    487 !****************************************************************************************
    488     IF (type_ocean=='couple') THEN
    489        CALL cpl_send_landice_fields(itime, knon, knindex, run_off_lic)
     487! If landice_opt>=2, corresponding call is done from surf_land_orchidee
     488!****************************************************************************************
     489    IF (type_ocean=='couple' .AND. landice_opt .LT. 2) THEN
     490       ! Compress fraction where run_off_lic is active (here all pctsrf(is_lic))
     491       run_off_lic_frac(:)=0.0
     492       DO j = 1, knon
     493          i = knindex(j)
     494          run_off_lic_frac(j) = pctsrf(i,is_lic)
     495       ENDDO
     496       
     497       CALL cpl_send_landice_fields(itime, knon, knindex, run_off_lic, run_off_lic_frac)
    490498    ENDIF
    491499
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/surf_ocean_mod.F90

    r4143 r4482  
    2020       z0m, z0h, SFRWL, alb_dir_new, alb_dif_new, evap, fluxsens, fluxlat, &
    2121       tsurf_new, dflux_s, dflux_l, lmt_bils, &
    22        flux_u1, flux_v1, delta_sst, delta_sal, ds_ns, dt_ns, dter, dser, tkt, tks, &
    23        taur, sss &
     22       flux_u1, flux_v1, delta_sst, delta_sal, ds_ns, dt_ns, dter, dser, &
     23       dt_ds, tkt, tks, taur, sss &
    2424#ifdef ISO
    2525        &       ,xtprecip_rain, xtprecip_snow,xtspechum,Roce, &
     
    114114    ! minus foundation temperature. (Can be negative.) In K.
    115115
     116    REAL, intent(inout):: dter(:) ! (knon)
     117    ! Temperature variation in the diffusive microlayer, that is
     118    ! ocean-air interface temperature minus subskin temperature. In
     119    ! K.
     120
     121    REAL, intent(inout):: dser(:) ! (knon)
     122    ! Salinity variation in the diffusive microlayer, that is
     123    ! ocean-air interface salinity minus subskin salinity. In ppt.
     124
     125    real, intent(inout):: dt_ds(:) ! (knon)
     126    ! (tks / tkt) * dTer, in K
     127
    116128    ! Output variables
    117129    !**************************************************************************
     
    129141    REAL, DIMENSION(klon), INTENT(OUT)       :: flux_u1, flux_v1
    130142
    131     REAL, intent(out):: dter(:) ! (knon)
    132     ! Temperature variation in the diffusive microlayer, that is
    133     ! ocean-air interface temperature minus subskin temperature. In
    134     ! K.
    135 
    136     REAL, intent(out):: dser(:) ! (knon)
    137     ! Salinity variation in the diffusive microlayer, that is
    138     ! ocean-air interface salinity minus subskin salinity. In ppt.
    139 
    140143    REAL, intent(out):: tkt(:) ! (knon)
    141144    ! épaisseur (m) de la couche de diffusion thermique (microlayer)
     
    152155    ! defined for subscripts 1:knon, but we have to declare it with
    153156    ! size klon because of the coupling machinery.)
     157
    154158#ifdef ISO
    155159    REAL, DIMENSION(ntraciso,klon), INTENT(out)        :: xtevap ! isotopes in surface evaporation flux
     
    172176    real s_int(knon) ! ocean-air interface salinity, in ppt
    173177
    174     !******************************************************************************
     178    !**************************************************************************
    175179
    176180#ifdef ISO
     
    213217    ENDIF
    214218
    215 
    216219    rhoa = PS(:KNON) / (Rd * temp_air(:knon) * (1. + retv * spechum(:knon)))
    217 
    218220    !******************************************************************************
    219221    ! Switch according to type of ocean (couple, slab or forced)
     
    232234            qsurf, evap, fluxsens, fluxlat, flux_u1, flux_v1, &
    233235            tsurf_new, dflux_s, dflux_l, sens_prec_liq, sss, delta_sal, rhoa, &
    234             delta_sst)
     236            delta_sst, dTer, dSer, dt_ds)
    235237
    236238    CASE('slab')
     
    376378       delta_sst = t_int - tsurf_new(:knon)
    377379       delta_sal = s_int - sss(:knon)
    378        if (activate_ocean_skin >= 2) tsurf_new(:knon) = t_int
     380
     381       if (activate_ocean_skin == 2) then
     382          tsurf_new(:knon) = t_int
     383          if (type_ocean == 'couple') dt_ds = (tks / tkt) * dter
     384       end if
    379385    end if
    380 
     386   
    381387  END SUBROUTINE surf_ocean
    382388  !****************************************************************************
  • LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/wake.F90

    r4143 r4482  
    25502550    d_deltaqw, sigmaw, d_sigmaw, alpha)
    25512551  ! ------------------------------------------------------
    2552   ! Dtermination du coefficient alpha tel que les tendances
     2552  ! D\'etermination du coefficient alpha tel que les tendances
    25532553  ! corriges alpha*d_G, pour toutes les grandeurs G, correspondent
    25542554  ! a une humidite positive dans la zone (x) et dans la zone (w).
Note: See TracChangeset for help on using the changeset viewer.