Ignore:
Timestamp:
Oct 19, 2023, 4:02:57 PM (8 months ago)
Author:
idelkadi
Message:

Merged trunk changes -r4488:4726 LMDZ_ECRad branch

Location:
LMDZ6/branches/LMDZ_ECRad
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/branches/LMDZ_ECRad

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

    r4482 r4727  
    154154&                   'revap_ilp 131')
    155155       enddo
     156#ifdef ISOTRAC
     157     call iso_verif_traceur(zxtrfl_ancien(1,i), &
     158&           'iso_revap_fisrtilp 158: debut')
     159#endif
    156160#endif
    157161    endif !if (zrfln(i).gt.ridicule*1e-2) then 
    … …  
    193197        endif !if (zq(i).gt.ridicule) then
    194198     endif !if ((iso_HDO.gt.0.and.(iso_O18.gt.0) then
     199#ifdef ISOTRAC
     200     call iso_verif_traceur(zxtrfl_ancien(1,i), &
     201&           'iso_revap_fisrtilp 201: debut quand pas de precip')
     202#endif
    195203!     write(*,*) 'iso_routines tmp 184'
    196204#endif               
    … …  
    504512!             stop
    505513  endif
    506   enddo
     514  enddo !do iiso=1,niso
    507515enddo !do i=1,ncas_evap_liq
    508516#endif         
    … …  
    514522&       izone,zqevfl(1),Exi(1,1),fac_fluxtomixratio(1), &
    515523&       xtrevap_tag(1,1),1,hdiag(1))
     524        ! dans cette routine, zxtrfl reçoit zxtrfln_cas
    516525
    517526  enddo !do izone=1,ntraceurs_zone
     527
    518528#ifdef ISOVERIF
    519529do i=1,ncas_evap_liq
    … …  
    522532&              0.0,'revap_ilp 414')
    523533enddo
     534     call iso_verif_traceur(zxtrfl(1,cas_evap_liq(i)), &
     535&           'iso_revap_fisrtilp 470a: apres stewart_explicite_vectall')
    524536enddo !do i=1,ncas
    525537#endif           
    … …  
    536548     call iso_verif_traceur(zxt(1,cas_evap_liq(i)), &
    537549&           'iso_revap_fisrtilp 282')
     550     call iso_verif_traceur(zxtrfl(1,cas_evap_liq(i)), &
     551&           'iso_revap_fisrtilp 804a')
    538552     call iso_verif_traceur(zxtrfln(1,cas_evap_liq(i)), &
    539 &           'iso_revap_fisrtilp 804')
     553&           'iso_revap_fisrtilp 804b')
    540554     do ixt=1,ntraciso
    541555        call iso_verif_positif_choix(zxt(ixt,cas_evap_liq(i)), &
    … …  
    676690   if (iso_HDO.gt.0) then
    677691      do i=1,ncas_evap_glace
    678          if (zrfln_cas(i).gt.ridicule_rain) then
    679                 call iso_verif_aberrant( &
    680 &                  (zxtrfln_cas(iso_HDO,i) &
    681 &                 /zrfln_cas(i)), 'iso_revap_fisrtilp 4563')
    682          endif
     692                call iso_verif_aberrant_choix(zxtrfln_cas(iso_HDO,i), zrfln_cas(i), &
     693                  ridicule_rain,deltalim_snow, 'iso_revap_fisrtilp 4563')
    683694       enddo !do i=1,ncas_evap_glace
    684695   endif !if ((iso_HDO.gt.0).and.(ixt.eq.iso_HDO)) then
    … …  
    789800     call iso_verif_traceur(zxt(1,cas_evap_glace(i)), &
    790801&           'iso_revap_fisrtilp 1033')
     802     call iso_verif_traceur(zxtrfl(1,cas_evap_glace(i)), &
     803&           'iso_revap_fisrtilp 1035a')
    791804     call iso_verif_traceur(zxtrfln(1,cas_evap_glace(i)), &
    792 &           'iso_revap_fisrtilp 1035')
     805&           'iso_revap_fisrtilp 1035b')
    793806  enddo
    794807#endif         
    … …  
    830843    call iso_verif_traceur(zxt(1,i),'iso_revap_fisrtilp 532')
    831844    call iso_verif_traceur(zxtrfln(1,i), &
    832 &           'iso_revap_fisrtilp 533')
     845&           'iso_revap_fisrtilp 533a')
     846    call iso_verif_traceur(zxtrfl(1,i), &
     847&           'iso_revap_fisrtilp 533b')
    833848    do ixt=1,ntraciso
    834849        call iso_verif_positif_choix(zxt(ixt,i),0.0, &
    … …  
    836851    enddo
    837852enddo !do i=1,klon
    838         write(*,*) 'revap_ilp 814: sortie'
     853        !write(*,*) 'revap_ilp 814: sortie'
    839854#endif       
    840855#endif         
    … …  
    23842399    endif !if (iso_eau.gt.0) then
    23852400    if (abs(evap(i)).gt.ridicule_rain) then
    2386       if (iso_HDO.gt.0) then                 
    2387          if (iso_verif_aberrant_nostop(xtevap(iso_HDO,i)/evap(i), &
    2388 &           'iso_surf>iso_rosee_givre 3193').eq.1) then   
     2401      if (iso_HDO.gt.0) then                     
     2402         if (iso_verif_aberrant_choix_nostop(xtevap(iso_HDO,i),evap(i), &
     2403&           ridicule_rain,deltalim_snow,'iso_surf>iso_rosee_givre 3193').eq.1) then   
    23892404                write(*,*) 'zxtalphai(iso_HDO)=',zxtalphai(iso_HDO)
    23902405                write(*,*) 'deltaD1eff=',deltaD(xt1lay(iso_HDO,i)/q1lay(i))
    23912406                write(*,*) 'tsurf(i)=',tsurf(i)
    23922407                write(*,*) 'q1lay(i)=',q1lay(i)   
    2393                 stop 
     2408                !stop 
    23942409           endif !if (iso_verif_aberrant_nostop
    23952410      endif  !if (iso_HDO.gt.0) then 
    … …  
    25342549
    25352550! quelques verifs de bilan d'eau
    2536 #ifdef ISOVERIF      
     2551#ifdef ISOVERIF   
    25372552do il=1,ncas
    25382553  do ixt=1,niso
    … …  
    26472662&             'stewart_explicite_vectall 220')
    26482663    enddo
    2649 #endif
    2650 #ifdef ISOVERIF
    2651       if (iso_eau.gt.0) then
     2664    if (iso_eau.gt.0) then
    26522665        call iso_verif_egalite_choix( &
    26532666&              (Exi(iso_eau,il)*fac_ftmr(il)), &
    … …  
    26692682&                   *fac_ftmr(il))),'stewart_explicite 214')
    26702683       endif !if ((iso_HDO.gt.0).and.
    2671        if ((debug.eq.1).and.(il.eq.il_debug)) then
    2672            write(*,*) 'stewart_explicit 224: cas Pqisup<=0'
    2673            write(*,*) 'Eqi(il),deltaD=',Eqi(il), &
    2674 &                   deltaD((Exi(iso_HDO,il)/Eqi(il)))
    2675        endif
     2684
    26762685#endif             
    26772686else !if (Pqisup.eq.0) then
    … …  
    26792688h(il)=qeff(il)/qs(il)       
    26802689h(il)= MAX(MIN(h(il),1.0),0.0)
    2681 #ifdef ISOVERIF       
     2690#ifdef ISOVERIF
    26822691call iso_verif_positif(h(il)-thumxt1,'stewart_explicit 209')
    26832692#endif       
    … …  
    27392748      endif !if ((iso_HDO.gt.0).and.
    27402749    endif !if (iso_HDO.gt.0) then
    2741     if ((debug.eq.1).and.(il.eq.il_debug)) then
     2750    if ((debug.eq.1).and.(il.eq.il_debug).and.(Eqi(il).gt.0.)) then
    27422751           write(*,*) 'stewart_explicit 302: cas evap~0'
    27432752           write(*,*) 'deltaDv est inchangé:',deltaD( &
    … …  
    27792788f(il)=m(il)/m0(il)   
    27802789! verifs
    2781 #ifdef ISOVERIF        
     2790#ifdef ISOVERIF       
    27822791call iso_verif_positif((m(il)), &
    27832792&           'stewart_explicite 173')
    … …  
    28362845    ! rajout verif 4 sept 2009
    28372846    if (iso_HDO.gt.0) then
    2838      if (Pqisup(il).gt.ridicule) then
    2839         call iso_verif_aberrant((Rl0(iso_HDO,il)), &
    2840 &                   'stewart_explicite 368')
    2841      endif
    2842     endif
     2847        call iso_verif_aberrant_choix(Rl0(iso_HDO,il)*Pqisup(il),Pqisup(il), &
     2848&                   ridicule_rain,deltalim_snow,'stewart_explicite 368')
     2849    endif !if (iso_HDO.gt.0) then
    28432850  endif !(iso_eau.gt.0)
    28442851#endif
    … …  
    28682875&                   'stewart_explicite 271')
    28692876   enddo !do ixt=1,niso   
    2870 #endif
    2871 #ifdef ISOVERIF 
    28722877   if (iso_eau.gt.0) then
    28732878      call iso_verif_egalite_choix( &
    … …  
    29232928       endif !if (qp0(il)+Eqi(il)*fac_ftmr(il).gt.ridicule) then
    29242929     endif !if (iso_HDO.gt.0)
    2925      if ((debug.eq.1).and.(il.eq.il_debug)) then
     2930     if ((debug.eq.1).and.(il.eq.il_debug).and.(Eqi(il).gt.0.)) then
    29262931         write(*,*) 'stewart_explicit 442: tout se réévapore'
    29272932         write(*,*) 'Eqi(il),deltaD=',Eqi(il), &
    … …  
    29432948if ((h(il).gt.0.99).or. &
    29442949&           (h(il).gt.0.98).and.(f(il).lt.1e-3)) then
    2945 #ifdef ISOVERIF         
    2946 !            write(*,*) 'stewart_explicit 191: cas h=1: il=',il
    2947 #endif           
    29482950    do ixt=1,niso
    29492951    interm(ixt,il)=alphap(ixt,il)*(1.0-h(il)) &
    … …  
    29812983&           'stewart_explicite 261')
    29822984    enddo !do ixt=1,niso
    2983 #endif           
    2984 #ifdef ISOVERIF
    29852985    if (iso_eau.gt.0) then
    29862986        call iso_verif_egalite_choix( &
    … …  
    30073007     
    30083008     if (iso_HDO.gt.0) then
    3009        if (Pqiinf(il).gt.ridicule_rain) then
    3010            if (iso_verif_aberrant_nostop( &
    3011 &                  (Pxtiinf(iso_HDO,il)/Pqiinf(il)), &
    3012 &                   'stewart_explicie 248').eq.1) then
     3009       if (iso_verif_aberrant_choix_nostop(Pxtiinf(iso_HDO,il),Pqiinf(il), &
     3010&                ridicule_rain,deltalim_snow,'stewart_explicite 248').eq.1) then
    30133011             write(*,*) 'cas reeq totale, il=',il
    30143012             write(*,*) 'deltaDl0=',deltaD( &
    … …  
    30193017&                   (Rb(iso_hdo,il)))
    30203018             stop
    3021            endif
    3022        endif  !if ((iso_HDO.gt.0).and.(ixt.eq.iso_HDO).and.
     3019       endif !if (iso_verif_aberrant_choix_nostop
    30233020       if (qp0(il)+Eqi(il)*fac_ftmr(il).gt.ridicule) then
    30243021         call iso_verif_aberrant(( &
    … …  
    30273024       endif !if ((iso_HDO.gt.0).and.
    30283025     endif  !if (iso_HDO.gt.0) then
    3029     if ((debug.eq.1).and.(il.eq.il_debug)) then
     3026     if ((debug.eq.1).and.(il.eq.il_debug).and.(Eqi(il).gt.0.)) then
    30303027         write(*,*) 'stewart_explicit 526: cas h~1: rééq'
    30313028         write(*,*) 'Eqi(il),deltaD=',Eqi(il), &
    … …  
    30343031&                   (Rb0(iso_hdo,il))),deltaD( &
    30353032&                   (Rl0(iso_hdo,il)))
    3036     endif
     3033     endif !if ((debug.eq.1).and.(il.eq.il_debug)) then
    30373034#endif
    30383035    ! end verifs
    30393036
    30403037else if ((f(il).gt.0.998).and. &
    3041 &           (Eqi(il)*fac_ftmr(il).lt.1e-2*qp0(il))) then
     3038&           (Eqi(il)*fac_ftmr(il).lt.1e-2*qp0(il))) then ! if ((h(il).gt.0.99).or.
    30423039
    30433040!*** cas particulier pour éviter imprécisions numériques:
    … …  
    30833080&           'stewart_explicite 397')
    30843081   enddo !do ixt=1,niso
    3085 #endif
    3086 #ifdef ISOVERIF
    30873082   if (iso_eau.gt.0) then
    30883083        call iso_verif_egalite_choix( &
    … …  
    31673162           stop
    31683163         endif !if (iso_verif_aberrant_nostop((
    3169        endif !if ((iso_HDO.gt.0).and.
     3164       endif !if (qp0(il)+Eqi(il)*fac_ftmr(il).gt.ridicule) then
    31703165     endif !!if ((iso_HDO.gt.0)
    3171      if ((debug.eq.1).and.(il.eq.il_debug)) then
     3166     if ((debug.eq.1).and.(il.eq.il_debug).and.(Eqi(il).gt.0.)) then
    31723167         write(*,*) 'stewart_explicit 663: cas où réévap faible'
    31733168         write(*,*) 'ordre 1 pour la vapeur et le liquide'
    … …  
    32713266             endif !if (qp0(il)+Eqi(il)*fac_ftmr(il).gt.ridicule) then
    32723267          endif !if ((iso_HDO.gt.0)
    3273      if ((debug.eq.1).and.(il.eq.il_debug)) then
     3268     if ((debug.eq.1).and.(il.eq.il_debug).and.(Eqi(il).gt.0.)) then
    32743269         write(*,*) 'stewart_explicit 767: cas de réévap sèche'
    32753270         write(*,*) 'distill de Rayleigh'
    … …  
    33173312&                   'stewart_explicite 467')
    33183313    enddo !do ixt=1,niso   
    3319 #endif         
    3320 #ifdef ISOVERIF 
    33213314          if (iso_eau.gt.0) then
    33223315            call iso_verif_egalite_choix( &
    … …  
    33603353            endif !if (qp0(il)+Eqi(il)*fac_ftmr(il).gt.ridicule) then
    33613354          endif !if (iso_HDO.gt.0)
    3362       if ((debug.eq.1).and.(il.eq.il_debug)) then
     3355      if ((debug.eq.1).and.(il.eq.il_debug).and.(Eqi(il).gt.0.)) then
    33633356         write(*,*) 'stewart_explicit 831: flux de masse vap~0'
    33643357         write(*,*) 'Eqi(il),deltaD=',Eqi(il), &
    … …  
    34293422        enddo !do ixt=1,niso   
    34303423#endif
    3431 #ifdef ISOVERIF   
     3424#ifdef ISOVERIF
    34323425          if (iso_eau.gt.0) then
    34333426            call iso_verif_egalite_choix( &
    … …  
    35463539        enddo !do ixt=1,niso
    35473540#endif
    3548 #ifdef ISOVERIF      
     3541#ifdef ISOVERIF   
    35493542          if (iso_eau.gt.0) then
    35503543            call iso_verif_egalite_choix( &
    … …  
    37503743       ! pour meilleure convergence numérique:
    37513744       !xtnew=qp0+Eqi*fac_ftmr
    3752    endif   ! if (iso_eau.gt.0).and.(ixt.eq.iso_eau)   
    3753 
     3745   endif   ! if (iso_eau.gt.0).and.(ixt.eq.iso_eau) 
    37543746   if (iso_HDO.gt.0) then
    37553747     if (Pqiinf(il).gt.ridicule_rain) then
    … …  
    37783770&                   xtnew(iso_HDO,il)/(qp0(il)+Eqi(il) &
    37793771&                   *fac_ftmr(il))),'stewart_explicite 912')
    3780         endif !if (qp0(il)+Eqi(il)*fac_ftmr(il).gt.ridicule) then
    3781       endif !if (iso_HDO.gt.0)
    3782       if ((debug.eq.1).and.(il.eq.il_debug)) then
     3772       endif !if (qp0(il)+Eqi(il)*fac_ftmr(il).gt.ridicule) then
     3773       if ((debug.eq.1).and.(il.eq.il_debug).and.(Eqi(il).gt.0.)) then
    37833774         write(*,*) 'stewart_explicit 991: fcas général'
    37843775         write(*,*) 'mais avec formule simplifiée'
    3785          write(*,*) 'il,Eqi(il),deltaD=',il,Eqi(il), &
    3786 &                   deltaD((Exi(iso_HDO,il)/Eqi(il)))
    3787       endif
    3788 #endif     
     3776         write(*,*) 'il,Eqi(il)=',il,Eqi(il)
     3777         write(*,*) 'deltaD=',deltaD((Exi(iso_HDO,il)/Eqi(il)))
     3778       endif
     3779      endif !if (iso_HDO.gt.0) 
     3780#endif
    37893781  ! end verifs 
    37903782   
    … …  
    38383830! compression
    38393831if (ncas_Jsimple+ncas_rieman.gt.0) then
    3840 !#ifdef ISOVERIF          
     3832!#ifdef ISOVERIF         
    38413833!        write(*,*) 'stewart_explicite_vectall 873:compression_calculJ'
    38423834!#endif       
    … …  
    39273919    ! ******* traitement vectoriel du cas Rieman (=2533)
    39283920if (ncas_rieman.gt.0) then
    3929 !#ifdef ISOVERIF
    3930 !!          write(*,*) 'traitement vectoriel rieman: x',ncas_rieman
    3931 !          do icas_rieman=1+ncas_Jsimple,ncas_rieman+ncas_Jsimple
    3932 !!            write(*,*) 'ntot_cas(icas_rieman)=',ntot_cas(icas_rieman)
    3933 !          call iso_verif_positif(float(ntot_cas(icas_rieman))-1.0,
    3934 !     :           'stewart_expl 984: ntot faux')
    3935 !          enddo   !do icas_rieman=1,ncas_rieman       
    3936 !#endif           
     3921           
    39373922   icas_rieman=1+ncas_Jsimple
    3938 !           write(*,*) 'stewart_expl 988 tmp: icas_rieman=',icas_rieman
    3939 !           write(*,*) 'qp0_cas(1)=',qp0_cas(1)
    3940 !           write(*,*) 'A_cas(1)=',A_cas(1)
    3941 !           write(*,*) 'm0_cas(1)=',m0_cas(1)
    3942 !           write(*,*) 'm_cas(1)=',m_cas(1)
    3943 !           if (iso_eau.gt.0) then
    3944 !           write(*,*) 'beta_cas(iso_eau,1)=',beta_cas(iso_eau,1)
    3945 !           write(*,*) 'gama_cas(iso_eau,1)=',gama_cas(iso_eau,1)
    3946 !           endif
    3947 !           write(*,*) 'f_cas(1)=',f_cas(1)
    3948 !           write(*,*) 'g_cas(1)=',g_cas(1)
    3949 !#ifdef rieman           
    3950 !           call integrale_rieman_vectall
    3951 !     :       (ncas_rieman,m_cas(icas_rieman),
    3952 !     :       J(1,icas_rieman),e(1,icas_rieman),
    3953 !     :       qp0_cas(icas_rieman),A_cas(icas_rieman),
    3954 !     :       m0_cas(icas_rieman),beta_cas(1,icas_rieman),
    3955 !     :       gama_cas(1,icas_rieman),f_cas(icas_rieman),
    3956 !     :       g_cas(icas_rieman),ntot_cas(icas_rieman))
    3957 !#else
     3923
    39583924   call integrale_gauss_vectall &
    39593925&       (ncas_rieman,m_cas(icas_rieman), &
    … …  
    39643930!     :       g_cas(icas_rieman),ntot_cas(icas_rieman))
    39653931&       g_cas(icas_rieman))     
    3966 !#endif
    3967 !#ifdef ISOVERIF
    3968 !           do il=1+ncas_Jsimple,ncas_rieman+ncas_Jsimple
    3969 !             do ixt=1,niso
    3970 !               call integrale_Rieman_precision(m_cas(il),m0_cas(il),
    3971 !     :          Jtmp,etmp,ntot_cas(il)*1e2,
    3972 !     :           qp0_cas(il),A_cas(il),m0_cas(il),
    3973 !     :          beta_cas(ixt,il),gama_cas(ixt,il),f_cas(il),g_cas(il))
    3974 !               call iso_verif_egalite_choix((Jtmp),
    3975 !     :           (J(ixt,il)),
    3976 !     :           'stewart_exp 999: test rieman',errmax,errmaxrel)
    3977 !               if ((iso_eau.gt.0).and.(ixt.eq.iso_eau)) then
    3978 !                 write(*,*) 'stew exp tmp 1005: il,J(iso_eau,il),Jtmp=',
    3979 !     :           il,J(iso_eau,il),Jtmp
    3980 !               endif
    3981 !             enddo             
    3982 !           enddo
    3983 !#endif           
    39843932
    39853933  endif !if (ncas_rieman.gt.0) then
    … …  
    40083956&          -beta_cas(ixt,il)*gama_cas(ixt,il)*r_jqp0(ixt,il) &
    40093957&                   /g_cas(il)/g_cas(il))
    4010 !#ifdef ISOVERIF           
    4011 !            if ((iso_eau.gt.0).and.(iso_eau.eq.ixt)) then
    4012 !              write(*,*) 'stewart_explicite tmp 1071: il=',il
    4013 !              if (il.le.ncas_Jsimple) then
    4014 !                write(*,*) 'cas_Jsimple(il)=',cas_Jsimple(il)
    4015 !              else !if (il.le.ncas_Jsimple) then
    4016 !                write(*,*) 'cas_rieman(il)=',cas_rieman(il)
    4017 !              endif !if (il.le.ncas_Jsimple) then
    4018 !              write(*,*) 'f_cas(il),beta_cas(ixt,il),gama_cas(ixt,il)=',
    4019 !     :                   f_cas(il),beta_cas(ixt,il),gama_cas(ixt,il)
    4020 !              write(*,*) 'g_cas(il),r_jqp0(ixt,il),r_jl0(ixt,il)=',
    4021 !     ;          g_cas(il),r_jqp0(ixt,il),r_jl0(ixt,il)
    4022 !              write(*,*) 'A_cas(il)=',A_cas(il)
    4023 !              write(*,*) 'pond Rl0=',(f_cas(il)**beta_cas(ixt,il))
    4024 !     :          *(g_cas(il)**(-beta_cas(ixt,il)*gama_cas(ixt,il)))
    4025 !     :          +beta(ixt,il)*gama_cas(ixt,il)*r_jqp0(ixt,il)
    4026 !     :                   /f_cas(il)/g_cas(il)
    4027 !              write(*,*) 'pond Rb0=',gama_cas(ixt,il)*beta_cas(ixt,il)
    4028 !     :                   *r_jl0(ixt,il)/f_cas(il)/g_cas(il)
    4029 !              write(*,*) 'pondRl0=fac1*fac2+t3'
    4030 !              write(*,*) 'fac1=',
    4031 !     :           f_cas(il)**beta_cas(ixt,il)
    4032 !              write(*,*) 'fac2=',
    4033 !     :             g_cas(il)**(-beta_cas(ixt,il)*gama_cas(ixt,il))
    4034 !              write(*,*) 't3=',
    4035 !     :           beta_cas(ixt,il)*gama_cas(ixt,il)*r_jqp0(ixt,il)
    4036 !     :                   /f_cas(il)/g_cas(il)     
    4037 !            endif
    4038 !#endif
    40393958
    40403959    Pxtiinf_cas(ixt,il)=Pqiinf_cas(il)*Rl(ixt,il)
    … …  
    40974016    enddo
    40984017#endif
    4099 #ifdef ISOVERIF            
     4018#ifdef ISOVERIF         
    41004019    if (iso_eau.gt.0) then
    41014020      if (iso_verif_egalite_choix_nostop( &
    … …  
    44114330   if (iso_HDO.gt.0) then
    44124331    if (Pqiinf(il).gt.ridicule_rain) then
    4413       if (iso_verif_aberrant_nostop( &
    4414 &                  (Pxtiinf(iso_HDO,il)/Pqiinf(il)), &
    4415 &                   'stewart_explicie 871').eq.1)  then
     4332      if (iso_verif_aberrant_choix_nostop(Pxtiinf(iso_HDO,il),Pqiinf(il),ridicule_rain,deltalim_snow, &
     4333&                   'stewart_explicite 871').eq.1)  then
    44164334      write(*,*) 'deltaDl0=',deltaD( &
    44174335&           (Rl0(iso_HDO,il)))
    … …  
    44284346    endif !if (iso_HDO.gt.0)
    44294347
    4430     if ((debug.eq.1).and.(il.eq.il_debug)) then
     4348    if ((debug.eq.1).and.(il.eq.il_debug).and.(Eqi(il).gt.0.)) then
    44314349          write(*,*) 'stewart_explicit 1558: cas avec calcul J'
    44324350          write(*,*) 'Eqi(il),deltaD=',Eqi(il), &
    … …  
    44784396
    44794397#ifdef ISOVERIF
    4480 !          write(*,*) 'stewart_explicite vectall 1179: fin'
     4398          write(*,*) 'stewart_explicite vectall 1179: fin'
    44814399#endif         
    44824400
    … …  
    47024620&           'stewart_sublim_nofrac 39')             
    47034621endif !if ((iso_eau.gt.0).and.(ixt.eq.iso_eau)) then
    4704 if ((iso_HDO.gt.0).and. &
    4705 &           (Pqisup(il).gt.ridicule_rain)) then
    4706     call iso_verif_aberrant( &
    4707 &           (Pxtisup(iso_HDO,il)/Pqisup(il)), &
    4708 &           'stewart_sublim_nofrac 40')     
     4622if (iso_HDO.gt.0) then ! Camille 9 mars 2023: moins stricte pour condensat
     4623    call iso_verif_aberrant_choix(Pxtisup(iso_HDO,il),Pqisup(il), &
     4624&           ridicule_rain,deltalim_snow, 'stewart_sublim_nofrac 40')
    47094625endif !if ((iso_HDO.gt.0).and.(ixt.eq.iso_HDO).and.
    47104626enddo !do il=1,ncas
    … …  
    48344750        enddo ! do ixt=1,niso
    48354751        ! verif que deltaD(Pqiinf) raisonable
    4836         if ((iso_HDO.gt.0).and. &
    4837      &           (Pqiinf(il).gt.ridicule_rain)) then
    4838             call iso_verif_aberrant( &
    4839      &           (Pxtiinf(iso_HDO,il)/Pqiinf(il)), &
    4840      &           'stewart_sublim 175')
     4752        if (iso_HDO.gt.0) then
     4753            call iso_verif_aberrant_choix(Pxtiinf(iso_HDO,il),Pqiinf(il), &
     4754     &           ridicule_rain,deltalim_snow, 'stewart_sublim 175')
    48414755        endif !if ((iso_HDO.gt.0).and.(ixt.eq.iso_HDO).and.
    48424756        if (iso_eau.gt.0) then
    … …  
    49234837!         integer ntot_cas(ncas)
    49244838         integer il,ixt       
    4925          integer coeff_precision
    4926          parameter (coeff_precision=0.2)
    49274839
    49284840
    … …  
    1082710739     &          Exi_cas(1,1),Exi(1,1), &
    1082810740#endif       
    10829      &          xtp_avantevaptrac_cas,0,hdiag(1)) ! hdiag ne sera pas utilisé
     10741     &          xtp_avantevaptrac_cas,0,hdiag(1)) ! hdiag ne sera pas utilise
    1083010742
    1083110743       enddo ! do izone=1,ntraceurs_zone
    … …  
    1261512527                  if (iso_verif_aberrant_choix_nostop( &
    1261612528     &               zxtliq(iso_HDO,i),cond(i), &
    12617      &               ridicule,deltalim, &
     12529     &               ridicule,deltalim_snow, &
    1261812530     &               'condiso_liq_ice_vectall 32b').eq.1) then
    1261912531                    write(*,*) 'deltaDvap=',deltaD(xt(iso_hdo,i)/qt(i))
    … …  
    1263612548                        write(*,*) 'deltaO18(zxtliq/cond)=',deltaO(zxtliq(iso_O18,i)/cond(i))
    1263712549                        write(*,*) 'tcond(i)=',tcond(i)-t_coup,'°C'
    12638                         stop
     12550                        !stop ! Camille 9 mars 2023: trop strict
    1263912551                    endif  !if (iso_verif_O18_aberrant_nostop(     
    1264012552                    endif ! if (iso_O18.gt.0) then     
    … …  
    1271412626                        write(*,*) 'zxtalphai(iso_O18,i)=',zxtalphai(iso_O18,i)
    1271512627                        write(*,*) 'xt(1:niso,i)=',xt(1:niso,i)
    12716                         stop
     12628                        !stop ! Camille 9 mars 2023: trop strict
    1271712629                    endif  !if (iso_verif_O18_aberrant_nostop(     
    1271812630                    endif ! if (iso_O18.gt.0) then
    … …  
    1276512677          if (zfice(i).gt.0.9) then
    1276612678              if (iso_verif_aberrant_choix_nostop( &
    12767      &          zxtice(iso_HDO,i),cond(i),ridicule,deltalim_snow, &
     12679     &          zxtice(iso_HDO,i),cond(i),ridicule,deltalim_snow, &
     12680                ! Camille 9 mars 2023: pour le condensat, on laisse plus de
     12681                ! marge
    1276812682     &          'condiso_liq_ice_vect 412').eq.1) then
    1276912683                write(*,*) 'debug condiso_liq_ice_vect 449: i,zfice=', &
    … …  
    1347013384          endif !if (iso_eau.gt.0) then
    1347113385#ifdef ISOTRAC
    13472           call iso_verif_traceur(xtsnow_evap(1,i), &
     13386!          call iso_verif_traceur(xtsnow_evap(1,i), &
     13387!     &           'gestion neige 2146') ! attention car snow_evap parfois
     13388!     négatif -> il ne faut pas passer dans les verifs de positivité.
     13389          call iso_verif_traceur_justmass(xtsnow_evap(1,i), &
    1347313390     &           'gestion neige 2146')
    1347413391#endif         
    … …  
    1414114058        if (iso_HDO.gt.0) then
    1414214059            call iso_verif_aberrant_choix(-xtsol_evap(iso_HDO,i), &
    14143      &           sol_evap(i),ridicule_evap,deltalim, &
    14144      &           'calcul_iso_surf_sic 257')
     14060     &           sol_evap(i),ridicule_evap,deltalim_snow, &
     14061     &           'calcul_iso_surf_sic 257_sol_evap')
    1414514062        endif
    1414614063#endif           
    … …  
    1416214079              endif !if (iso_eau.gt.0) then
    1416314080              if (iso_HDO.gt.0) then
    14164                 if (evap(i).gt.ridicule_evap) then 
    14165                    call iso_verif_aberrant(xtevap(iso_HDO,i)/evap(i), &
    14166      &                  'calcul_iso_surf_sic 257')
    14167                endif !if (evap(i).gt.ridicule_evap) then 
     14081                   call iso_verif_aberrant_choix(xtevap(iso_HDO,i),evap(i), &
     14082     &                  ridicule_evap,deltalim_snow,'calcul_iso_surf_sic 257_evap')
    1416814083              endif !if (iso_eau.gt.0) then
    1416914084#ifdef ISOTRAC
    … …  
    1432314238        enddo !do i=1,knon
    1432414239        endif !if (iso_HDO.gt.0) then
    14325        
    14326         if (iso_eau.gt.0) then 
     14240           
    1432714241        do i=1,knon
     14242        if (iso_eau.gt.0) then
    1432814243          call iso_verif_egalite_choix(xtsnow(iso_eau,i),snow(i),  &
    14329      &         'calcul_iso_surf_lic_vectall 587',errmax,errmaxrel)
    14330         enddo
     14244     &         'calcul_iso_surf_lic_vectall 587a',errmax,errmaxrel)
    1433114245        endif
     14246        if (iso_HDO.gt.0) then
     14247          call iso_verif_aberrant_choix(xtsnow(iso_HDO,i), &
     14248     &               snow(i),ridicule_snow,deltalim_snow, &
     14249     &               'calcul_iso_surf_lic 587b')
     14250        endif
     14251        enddo !do i=1,knon
    1433214252#endif             
    1433314253       
    … …  
    1447614396     &                  'calcul_iso_surf_lic 363',errmax,errmaxrel) 
    1447714397              endif !if (iso_eau.gt.0) then
    14478               if (iso_HDO.gt.0) then 
    14479                 if (snow(i).gt.ridicule_snow) then 
    14480                  call iso_verif_aberrant(xtsnow(iso_HDO,i)/snow(i), &
    14481      &                  'calcul_iso_surf_lic 367')
    14482                 endif !if (snow(i).gt.ridicule_snow) then
     14398              if (iso_HDO.gt.0) then   
    1448314399                call iso_verif_aberrant_choix(xtsnow(iso_HDO,i), &
    14484      &               snow(i),ridicule,deltalim_snow, &
     14400     &               snow(i),ridicule_snow,deltalim_snow, &
    1448514401     &               'calcul_iso_surf_lic 797')
    14486                 if (abs(evap(i)).gt.ridicule_evap) then 
    14487                  call iso_verif_aberrant(xtevap(iso_HDO,i)/evap(i), &
    14488      &                  'calcul_iso_surf_lic 369')
    14489                 endif ! if (evap(i).gt.ridicule_evap) then
     14402                call iso_verif_aberrant_choix(xtevap(iso_HDO,i),evap(i), &
     14403     &                  ridicule_evap,deltalim_snow, 'calcul_iso_surf_lic 369')
    1449014404              endif !if (iso_eau.gt.0) then
    1449114405#ifdef ISOTRAC
    … …  
    1579515709              endif !if (iso_eau.gt.0) then
    1579615710              if (iso_HDO.gt.0) then 
    15797                  if (snow(i).gt.ridicule_snow) then
    15798                     call iso_verif_aberrant(xtsnow(iso_HDO,i)/snow(i), &
    15799      &                  'calcul_iso_surf_ter 749')
    15800                  endif !if (snow(i).gt.ridicule_snow) then
     15711                 call iso_verif_aberrant_choix(xtsnow(iso_HDO,i),snow(i), &
     15712     &                  ridicule_snow,deltalim_snow, 'calcul_iso_surf_ter 749')
    1580115713                 call iso_verif_aberrant_choix(xtsnow(iso_HDO,i), &
    1580215714     &               snow(i),ridicule,deltalim_snow, &
    … …  
    1580915721                     write(*,*) 'sol_evap,snow_evap=', &
    1581015722     &                   sol_evap(i),snow_evap(i)
    15811                      write(*,*) 'deltaDsol_evap=', &
     15723                     if (sol_evap(i).gt.ridicule_evap)write(*,*) 'deltaDsol_evap=', &
    1581215724     &                   deltaD(xtsol_evap(iso_hdo,i)/sol_evap(i))
    15813                      write(*,*) 'deltaDsnow_evap=', &
     15725                     if (snow_evap(i).gt.ridicule_evap)write(*,*) 'deltaDsnow_evap=', &
    1581415726     &                   deltaD(xtsnow_evap(iso_hdo,i)/snow_evap(i))
    1581515727                     write(*,*) 'deltaD1new=',deltaD( &
    … …  
    1609416006          endif !if (iso_eau.gt.0) then
    1609516007          if (iso_HDO.gt.0) then
    16096             if (snow(i,nsrf).gt.ridicule_snow) then
    16097             call iso_verif_aberrant(xtsnow(iso_hdo,i,nsrf)/snow(i,nsrf), &
    16098      &                'phyisoetat0 117')
    16099             endif
     16008            call iso_verif_aberrant_choix(xtsnow(iso_hdo,i,nsrf),snow(i,nsrf), &
     16009     &                ridicule_snow, deltalim_snow, 'phyisoetat0 117')
    1610016010          endif !if (iso_eau.gt.0) then
    1610116011        enddo !do nsrf=1,nbsrf
    … …  
    1618016090        real deltaD_snow_fall_O18,deltaD_rain_fall_O18
    1618116091        real alpha(niso),kcin(niso)
    16182         real q0,h0
    16183         parameter (q0=20e-3,h0=0.7)
    1618416092!      character*50 text
    1618516093       
    … …  
    1631416222        do k=1,klev
    1631516223          do ixt=1,niso
    16316            RMerlivat(ixt)=toce(ixt)/alpha(ixt) &
    16317      &                  *(1-kcin(ixt))/(1.0-kcin(ixt)*h0)
    16318            xt_ancien(ixt,i,k)=q_ancien(i,k)*RMerlivat(ixt) &
    16319      &                  *(min(q0,q_ancien(i,k))/q0)**(alpha(ixt)-1.0)
     16224           call iso_init_ideal(q_ancien(i,k),xt_ancien(ixt,i,k),ixt, &
     16225                alpha(ixt),kcin(ixt),toce(ixt))
     16226
    1632016227           if (q_ancien(i,k).gt.ridicule) then
    1632116228           xtl_ancien(ixt,i,k)=ql_ancien(i,k)*alpha(ixt) &
    … …  
    1633816245     &           'phyisoetat0 16067') 
    1633916246        enddo !do ixt=1,niso
    16340         if ((k.eq.1).and.(iso_HDO.gt.0).and.(iso_O18.gt.0) &
    16341                 .and.(abs(q_ancien(i,k)-q0).lt.1e-3)) then
    16342                 ! vérifier qu'on est proche de la fermeture de Merlivat
    16343             write(*,*) 'i,k=',i,k
    16344             write(*,*) 'q_ancien(i,k)=',q_ancien(i,k)
    16345             write(*,*) 'deltaD=',deltaD(xt_ancien(iso_HDO,i,k)/q_ancien(i,k))
    16346             write(*,*) 'deltaDM=',deltaD(RMerlivat(iso_HDO))
    16347             write(*,*) 'deltaO=',deltaO(xt_ancien(iso_O18,i,k)/q_ancien(i,k))
    16348             write(*,*) 'deltaOM=',deltaO(RMerlivat(iso_O18))
    16349             write(*,*) 'dexcess=',dexcess(xt_ancien(iso_HDO,i,k)/q_ancien(i,k), &
    16350                     xt_ancien(iso_O18,i,k)/q_ancien(i,k))
    16351             write(*,*) 'dexcessM=',dexcess(RMerlivat(iso_HDO),RMerlivat(iso_O18))
    16352             write(*,*) 'kcin=',kcin
    16353             write(*,*) 'toce=',toce
    16354             write(*,*) 'alpha=',alpha
    16355             call iso_verif_positif(20.0-abs(-80.0-deltaD(xt_ancien(iso_HDO,i,k)/q_ancien(i,k))), &
    16356                     'phyisoetat0 16398a')   
    16357             call iso_verif_positif(5.0-abs(10.0-dexcess(xt_ancien(iso_HDO,i,k)/q_ancien(i,k), &
    16358                     xt_ancien(iso_O18,i,k)/q_ancien(i,k))),'phyisoetat0 16398b')
    16359         endif
     16247
     16248        ! Camille 7 mars 2023: ajout d'un check
     16249        if ((i.eq.1).and.(k.eq.1).and.(iso_HDO.gt.0)) then
     16250        write(*,*) 'phyisoetat0 16362: q_ancien(1,1)=',q_ancien(1,1)
     16251        write(*,*) 'deltaD_ancien=',deltaD(xt_ancien(iso_HDO,i,k)/q_ancien(i,k))
     16252        write(*,*) 'xt_ancien(:,i,k)=',xt_ancien(:,i,k)
     16253        endif !if ((i.eq.1).and.(k.eq.1)) then
     16254
     16255        if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then
     16256            if (q_ancien(i,k).gt.ridicule) then 
     16257               if (iso_verif_o18_aberrant_nostop( &
     16258     &              xt_ancien(iso_HDO,i,k)/q_ancien(i,k), &
     16259     &              xt_ancien(iso_O18,i,k)/q_ancien(i,k), &
     16260     &              'phyisoetat0 16366 q_ancien').eq.1) then
     16261                  write(*,*) 'phyisoetat0 16367: i,k,q_ancien(i,k)=',i,k,q_ancien(i,k)
     16262                  write(*,*) 'xt_ancien(:,i,k)=',xt_ancien(:,i,k)
     16263                  stop
     16264              endif !  if (iso_verif_o18_aberrant_nostop
     16265            endif !if (q_seri(i,k).gt.errmax) then
     16266        endif !if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then
    1636016267#endif
    1636116268
    … …  
    1654316450   CHARACTER(LEN=2) :: str2
    1654416451   CHARACTER(LEN=5) :: str5
    16545    CHARACTER(LEN=maxlen) :: outiso, oldIso, modname, nam(2)
     16452   CHARACTER(LEN=maxlen) :: outiso, oldIso, modname, nam(3), oldIso2
    1654616453   REAL :: xmin, xmax
    1654716454   LOGICAL :: found
    … …  
    1655916466      outiso = isoName(ixt)
    1656016467      oldIso = strTail(new2oldH2O(outiso), '_')            !--- Remove "H2O_" from "H2O_<iso>[_<tag>]"
     16468      i = INDEX(outiso, '_', .TRUE.)
     16469      oldIso2 = outiso(1:i-1)//outiso(i+1:LEN_TRIM(outiso)) ! CR 2023: on ajoute cette possibilité aussi, elle correspond au cas le plus récent.
     16470!      write(*,*) 'tmp 16541:'
     16471!      write(*,*) 'outiso=',outiso
     16472!      write(*,*) 'oldIso=',oldIso
     16473!      write(*,*) 'oldIso2=',oldIso2
     16474
    1656116475      ! on lit seulement si ixt<=niso ou si on initialise les traceurs d'après fichier:
    1656216476#ifdef ISOTRAC
    1656316477      IF(ixt <= niso .OR. initialisation_isotrac == 0) THEN
    1656416478#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)
    1656716479      found = phyetat0iso_srf3(fxtevap,     "XTEVAP", "evaporation",  0.)
     16480      if (.NOT.found) CALL abort_physic('isotopes_routines_mod', 'phyiso_etat0_fichier 16581a: unfound isotopic variable',1)
    1656816481      found = phyetat0iso_get2(xtrain_fall, "xtrain_f", "xrain fall", 0.)
    1656916482      found = phyetat0iso_get2(xtsnow_fall, "xtsnow_f", "xsnow fall", 0.)
    … …  
    1657116484      found = phyetat0iso_get3(xtl_ancien,  "XTLANCIEN", "QLANCIEN",  0.)
    1657216485      found = phyetat0iso_get3(xts_ancien,  "XTSANCIEN", "QSANCIEN",  0.)
    16573       found = phyetat0iso_get2(xtrun_off_lic_0, "XTRUNOFFLIC0", "RUNOFFLIC0", 0.)
    1657416486      found = phyetat0iso_get3(wake_deltaxt,  "WAKE_DELTAXT", "Delta hum. wake/env",  0.)
    1657516487#ifdef ISOVERIF
    … …  
    1658016492            DO nsrf = 1, nbsrf
    1658116493               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')
    1658316494            END DO
    1658416495         END DO
    … …  
    1659216503         END DO
    1659316504      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
    1660216505#endif
    1660316506      ! ces variables n'ont pas de traceurs:
    1660416507      IF(ixt <= niso) THEN
    1660516508         found = phyetat0iso_get2(xtsol, "XTSOL", "Surface humidity / bucket", 0.)
     16509         if (.NOT.found) CALL abort_physic('isotopes_routines_mod', 'phyiso_etat0_fichier 16581b: unfound isotopic variable',1)
    1660616510         found = phyetat0iso_get2(Rland_ice, "Rland_ice", "SR land ice", 0.)
    16607 #ifdef ISOVERIF
    16608 
     16511         found = phyetat0iso_srf3(xtsnow,      "XTSNOW", "Surface snow", 0.) ! CR avril 2023: deplacer ici
     16512         found = phyetat0iso_get2(xtrun_off_lic_0, "XTRUNOFFLIC0", "RUNOFFLIC0", 0.)
     16513#ifdef ISOVERIF
    1660916514         DO i=1,klon
    1661016515            IF(iso_verif_noNaN_nostop(xtsol(ixt,i),TRIM(modname)//' 95') == 1) THEN
    … …  
    1661216517               STOP
    1661316518            END IF
    16614          END DO
     16519            IF(ixt == iso_eau .AND. iso_eau > 0) THEN
     16520             DO nsrf = 1, nbsrf
     16521               CALL iso_verif_egalite(fxtevap(iso_eau,i,nsrf),fevap(i,nsrf),TRIM(modname)//' 231c')
     16522               CALL iso_verif_egalite( xtsnow(iso_eau,i,nsrf), snow(i,nsrf),TRIM(modname)//' 231d')
     16523             END DO
     16524             CALL iso_verif_egalite( xtrun_off_lic_0(iso_eau,i), run_off_lic_0(i),TRIM(modname)//' 231e')
     16525            ENDIF !IF(ixt == iso_eau .AND. iso_eau > 0) THEN
     16526         END DO !DO i=1,klon
    1661516527#endif
    1661616528      END IF
    … …  
    1668016592  nam(1) = TRIM(pref)//TRIM(outiso)
    1668116593  nam(2) = TRIM(pref)//TRIM(oldIso)
     16594  nam(3) = TRIM(pref)//TRIM(oldIso2)
    1668216595  lFound = phyetat0_get(iso_tmp, nam, descr, default)
    1668316596  field(ixt,:) = iso_tmp
    … …  
    1669216605  nam(1) = TRIM(pref)//TRIM(outiso)
    1669316606  nam(2) = TRIM(pref)//TRIM(oldIso)
     16607  nam(3) = TRIM(pref)//TRIM(oldIso2)
    1669416608  lFound = phyetat0_get(iso_tmp_lonlev, nam, descr, default)
    1669516609  field(ixt,:,:) = iso_tmp_lonlev(:,:)
    … …  
    1670316617  nam(1) = TRIM(pref)//TRIM(outiso)
    1670416618  nam(2) = TRIM(pref)//TRIM(oldIso)
     16619  nam(3) = TRIM(pref)//TRIM(oldIso2)
    1670516620  lFound = phyetat0_srf(iso_tmp_lonsrf, nam, descr, default)
    1670616621  field(ixt,:,:) = iso_tmp_lonsrf
    … …  
    1812918044
    1813018045      do iessai=1,nessai
    18131          day_nucl(iessai)   = 0.
    18132          month_nucl(iessai) = 0.
    18133          year_nucl(iessai)  = 0.
     18046         day_nucl(iessai)   = 0
     18047         month_nucl(iessai) = 0
     18048         year_nucl(iessai)  = 0
    1813418049         lat_nucl(iessai)   = 0.
    1813518050         lon_nucl(iessai)   = 0.
    … …  
    1881618731          do i=1,n
    1881718732            qttrac(i)=xt(ieau,i)
    18818 !            if (qt(i).gt.0.0) then ! modif C Risi juillet 2020
    18819             if ((qt(i).gt.0.0).and.(xt(ieau,i).gt.0.0)) then
    18820                zcondtrac(i)=(zcond(i)/qt(i))*xt(ieau,i)
     18733            if (qt(i).gt.0.0) then ! modif C Risi juillet 2020 ! remodif Camille 9 mars 2023
     18734!            if ((qt(i).gt.0.0).and.(xt(ieau,i).gt.0.0)) then
     18735               zcondtrac(i)=(zcond(i)/qt(i))*qttrac(i)
    1882118736            else !if (qt(i).eq.0) then
    1882218737#ifdef ISOVERIF             
    … …  
    1883618751            endif
    1883718752            if (iso_HDO.gt.0) then
    18838 !              if (qttrac(i).gt.ridicule_trac) then
    1883918753                call iso_verif_aberrant_choix(xttrac(iso_HDO,i), &
    1884018754     &           qttrac(i),ridicule_trac,deltalimtrac, &
    1884118755     &           'condisotrac 205')
    18842 !              endif
    1884318756            endif
    1884418757            call iso_verif_positif(qt(i)-cond(i), &
    … …  
    1885318766          call condiso_liq_ice_vectall(xttrac,qttrac,zcondtrac, &
    1885418767     &           tcond,zfice,zxticetrac,zxtliqtrac,n)
    18855 #ifdef ISOVERIF         
    18856           write(*,*) 'condisotrac 167: après condiso'
    18857 #endif           
     18768
    1885818769          do i=1,n
    1885918770          do iiso=1,niso         
    … …  
    1888118792#endif
    1888218793
     18794subroutine iso_init_ideal(q,xt,ixt,alpha,kcin,toce)
     18795
     18796        USE isotopes_mod, ONLY: iso_eau,iso_HDO,ridicule
     18797#ifdef ISOVERIF
     18798        USE isotopes_verif_mod
     18799#endif
     18800        implicit none
     18801
     18802        ! inputs
     18803        real q ! humidité spec
     18804        integer ixt ! indice isotopique
     18805        real alpha ! coef frac à l'eq
     18806        real kcin ! coef frac cinétique
     18807        real toce ! rapport iso ds ocean surface
     18808
     18809        ! outputs
     18810        real xt ! equivalent iso de l'humidité spec, même unité.
     18811
     18812        ! locals
     18813        real RMerlivat
     18814        real q0,h0 ! conditions initiales de la distill de Rayleigh
     18815        parameter (q0=20e-3,h0=0.7)
     18816
     18817        ! verifier que ixt est un isotope et pas un tagging
     18818        if (ixt.gt.niso) then
     18819          CALL abort_physic('isotopes_routines_mod', 'iso_init_ideal, ixt>niso', 1)
     18820        endif
     18821
     18822        ! R selon Merlivat:
     18823        RMerlivat=toce/alpha *(1.0-kcin)/(1.0-kcin*h0)
     18824
     18825        ! R d'après Rayleigh
     18826        xt=q*RMerlivat*(min(q0,q)/q0)**(alpha-1.0)
     18827
     18828#ifdef ISOVERIF
     18829        call iso_verif_noNaN(xt, 'isotopes_routines_mod 18930a: iso_init_ideal')
     18830        if ((iso_HDO.gt.0).and.(ixt.eq.iso_HDO)) then
     18831            if (q.gt.ridicule) then
     18832                call iso_verif_aberrant_encadre( xt/q, 'isotopes_routines_mod 18930b: iso_init_ideal')
     18833            endif
     18834        endif
     18835        if ((iso_eau.gt.0).and.(ixt.eq.iso_eau)) then
     18836             call iso_verif_egalite(xt,q, 'isotopes_routines_mod 18930c: iso_init_ideal')
     18837        endif
     18838#endif
     18839       
     18840
     18841end subroutine iso_init_ideal
     18842
     18843
     18844subroutine appel_stewart_debug(lwork,nloc,inb,na,i, &
     18845                evap,water,rpprec,rr,wdtrain, &
     18846                xtevap,xtwater,xtp,xt,xtwdtrain)
     18847USE isotopes_mod, ONLY: iso_eau, iso_HDO,thumxt1, &
     18848&       bidouille_anti_divergence,ridicule,Rdefault
     18849use infotrac_phy, ONLY: ntraciso=>ntiso, niso
     18850#ifdef ISOTRAC
     18851    use isotrac_mod, only: option_cond,izone_cond,index_iso,index_zone,izone_poubelle
     18852#endif
     18853#ifdef ISOVERIF
     18854USE isotopes_verif_mod
     18855#endif
     18856implicit none
     18857
     18858
     18859! inputs
     18860integer nloc,na,i ! dimension horiz effective
     18861logical lwork(nloc)
     18862real wdtrain(nloc),xtwdtrain(ntraciso,nloc)
     18863real xt(ntraciso,nloc,na)
     18864real evap(nloc,na),water(nloc,na),rpprec(nloc,na),rr(nloc,na)
     18865integer inb(nloc)
     18866
     18867! outputs
     18868real xtevap(ntraciso,nloc,na),xtwater(ntraciso,nloc,na),xtp(ntraciso,nloc,na)
     18869
     18870! locals
     18871integer il,ixt
     18872
     18873     do il=1,nloc
     18874       if (i.le.inb(il) .and. lwork(il)) then
     18875          if (wdtrain(il).gt.0.) then
     18876            do ixt=1,ntraciso
     18877             xtwater(ixt,il,i)= xtwdtrain(ixt,il)/wdtrain(il)*water(il,i)
     18878             xtevap(ixt,il,i)= xtwdtrain(ixt,il)/wdtrain(il)*evap(il,i)
     18879            enddo
     18880          else !if (wdtrain(il).gt.0.) then
     18881            do ixt=1,niso
     18882             xtwater(ixt,il,i)= Rdefault(ixt)*water(il,i)
     18883             xtevap(ixt,il,i)= Rdefault(ixt)*evap(il,i)
     18884            enddo
     18885#ifdef ISOTRAC
     18886            do ixt=1+niso,ntraciso
     18887             if (index_zone(ixt).eq.izone_poubelle) then
     18888               xtwater(ixt,il,i)= Rdefault(index_iso(ixt))*water(il,i)
     18889               xtevap(ixt,il,i)= Rdefault(index_iso(ixt))*evap(il,i)
     18890             else
     18891               xtwater(ixt,il,i)= 0.
     18892               xtevap(ixt,il,i)=0.
     18893             endif
     18894            enddo ! do ixt=1+niso,ntraciso
     18895#endif
     18896         endif !if (wdtrain(il).gt.0.) then
     18897         do ixt=1,ntraciso
     18898             xtp(ixt,il,i)= xt(ixt,il,i)/rr(il,i)*rpprec(il,i)
     18899         enddo !do ixt=1,ntraciso
     18900      endif
     18901    enddo ! do il=1,ncum
     18902end subroutine appel_stewart_debug
     18903
    1888318904END MODULE isotopes_routines_mod
    1888418905#endif
Note: See TracChangeset for help on using the changeset viewer.