Ignore:
Timestamp:
Jan 12, 2022, 10:54:09 PM (2 years ago)
Author:
dcugnet
Message:

Most of the changes are intended to help to eventually remove the constraints about the tracers assumptions, in particular water tracers.

  • Remove index tables itr_indice and niadv, replaced by tracers(:)%isAdvected and tracers(:)%isH2OFamily. Most of the loops are now from 1 to nqtot:
    • DO iq=nqo+1,nqtot loops are replaced with: DO iq=1,nqtot

IF(tracers(iq)%isH2Ofamily) CYCLE

  • DO it=1,nbtr; iq=niadv(it+nqo)

and DO it=1,nqtottr; iq=itr_indice(it) loops are replaced with:

it = 0
DO iq = 1, nqtot

IF(.NOT.tracers(iq)%isAdvected .OR. tracers(iq)%isH2Ofamily) CYCLE
it = it+1

  • Move some StratAer? related code from infotrac to infotrac_phy
  • Remove "nqperes" variable:

DO iq=1,nqpere loops are replaced with:
DO iq=1,nqtot

IF(tracers(iq)%parent/='air') CYCLE

  • Cosmetic changes (justification, SELECT CASE instead of multiple IF...) mostly in advtrac* routines.
File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/phylmdiso/physiq_mod.F90

    r4050 r4056  
    3939    USE ioipsl_getin_p_mod, ONLY : getin_p
    4040    USE indice_sol_mod
    41     USE infotrac_phy, ONLY: nqtot, nbtr, nqo, type_trac,ok_isotopes, &
    42         nqtottr,itr_indice ! C Risi
    43 
     41    USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, type_trac, nqCO2, ok_isotopes
     42    USE readTracFiles_mod, ONLY: phases_sep
     43    USE strings_mod,  ONLY: strIdx
    4444    USE iophy
    4545    USE limit_read_mod, ONLY : init_limit_read
    … …  
    6161    USE phys_output_mod
    6262    USE phys_output_ctrlout_mod
    63     USE print_control_mod, ONLY: mydebug=>debug , lunout, prt_level
     63    USE print_control_mod, ONLY: mydebug=>debug , lunout, prt_level, &
     64         alert_first_call, call_alert, prt_alerte
    6465    USE readaerosol_mod, ONLY : init_aero_fromfile
    6566    USE readaerosolstrato_m, ONLY : init_readaerosolstrato
    … …  
    123124#ifdef ISO
    124125    USE infotrac_phy, ONLY:  &
    125         iqiso,iso_indnum,tracers,ok_isotrac, &
    126         niso,ntraciso,nqtottr,itr_indice ! ajout C Risi pour isos
     126        iqiso,iso_indnum,ok_isotrac,niso, ntraciso
    127127     USE isotopes_mod, ONLY: iso_eau,iso_HDO,iso_O18,iso_O17,iso_HTO, &
    128128        & bidouille_anti_divergence,ok_bidouille_wake, &
    … …  
    188188       !
    189189       d_t_vdf,d_q_vdf,d_u_vdf,d_v_vdf,d_t_diss, &
    190        d_t_vdf_w,d_q_vdf_w, &
    191        d_t_vdf_x,d_q_vdf_x, &
     190       d_t_vdf_x, d_t_vdf_w, &
     191       d_q_vdf_x, d_q_vdf_w, &
    192192       d_ts, &
    193193       !
    … …  
    262262       zxfluxlat_x, zxfluxlat_w, &
    263263       !
    264        d_t_vdf_x, d_t_vdf_w, &
    265        d_q_vdf_x, d_q_vdf_w, &
    266264       pbl_tke_input, tke_dissip, l_mix, wprime,&
    267265       t_therm, q_therm, u_therm, v_therm, &
    … …  
    939937    real zqsat(klon,klev)
    940938    !
    941     INTEGER i, k, iq, j, nsrf, ll, l
    942     INTEGER itr ! C Risi
     939    INTEGER i, k, iq, j, nsrf, ll, l, itr
    943940#ifdef ISO
    944941    real zxt_apres(ntraciso,klon)
    … …  
    11331130!JLD    REAL zstophy, zout
    11341131
    1135     CHARACTER*20 modname
     1132    CHARACTER (LEN=20) :: modname='physiq_mod'
    11361133    CHARACTER*80 abort_message
    11371134    LOGICAL, SAVE ::  ok_sync, ok_sync_omp
    … …  
    13061303    pi = 4. * ATAN(1.)
    13071304
     1305    ! set-up call to alerte function
     1306    call_alert = (alert_first_call .AND. is_master)
     1307   
    13081308    ! Ehouarn: set value of jjmp1 since it is no longer a "fixed parameter"
    13091309    jjmp1=nbp_lat
    … …  
    14241424    forall (k=1: nbp_lev) zmasse(:, k) = (paprs(:, k)-paprs(:, k+1)) / rg
    14251425
    1426     modname = 'physiq'
    14271426
    14281427    IF (debut) THEN
    … …  
    18531852
    18541853       CALL iniradia(klon,klev,paprs(1,1:klev+1))
    1855 
    1856        ! Initialisation des champs dans phytrac* qui sont utilisés par phys_output_write*
     1854       !
     1855!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1856       ! Initialisation des champs dans phytrac* qui sont utilises par phys_output_write*
     1857       !
     1858!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1859
    18571860#ifdef CPP_Dust
    18581861       ! Quand on utilise SPLA, on force iflag_phytrac=1
    … …  
    18991902            ENDDO
    19001903          ENDDO
    1901         ELSE
     1904       ELSE
    19021905          pbl_tke(:,:,is_ave) = 0. !ym missing init : maybe must be initialized in the same way that for klon_glo==1 ??
    19031906!>jyg
    … …  
    19491952          CALL abort_physic(modname,abort_message,1)
    19501953       ENDIF
     1954
     1955!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1956       ! Initialisation pour la convection de K.E. et pour les poches froides
     1957       !
     1958!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1959
    19511960       WRITE(lunout,*)"Clef pour la convection, iflag_con=", iflag_con
    1952        WRITE(lunout,*)"Clef pour le driver de la convection, ok_cvl=", &
    1953             ok_cvl
     1961       WRITE(lunout,*)"Clef pour le driver de la convection, ok_cvl=", ok_cvl
    19541962       !
    19551963       !KE43
    … …  
    20042012             d_deltaxt_ajs_cv(:,:,:) = 0.
    20052013#endif
    2006           ENDIF
     2014          ENDIF  !  (iflag_wake>=1)
    20072015
    20082016          !        do i = 1,klon
    … …  
    20152023       !   ALLOCATE(lonGCM(0), latGCM(0))
    20162024       !   ALLOCATE(iGCM(0), jGCM(0))
    2017        ENDIF
    2018 
     2025       ENDIF  !  (iflag_con.GE.3)
     2026       !
    20192027       DO i=1,klon
    20202028          rugoro(i) = f_rugoro * MAX(1.0e-05, zstd(i)*zsig(i)/2.0)
    … …  
    20852093       !$OMP BARRIER
    20862094       missing_val=missing_val_omp
     2095       !
     2096       ! Now we activate some double radiation call flags only if some
     2097       ! diagnostics are requested, otherwise there is no point in doing this
     2098       IF (is_master) THEN
     2099         !--setting up swaero_diag to TRUE in XIOS case
     2100         IF (xios_field_is_active("topswad").OR.xios_field_is_active("topswad0").OR. &
     2101            xios_field_is_active("solswad").OR.xios_field_is_active("solswad0").OR. &
     2102            xios_field_is_active("topswai").OR.xios_field_is_active("solswai").OR.  &
     2103              (iflag_rrtm==1.AND.(xios_field_is_active("toplwad").OR.xios_field_is_active("toplwad0").OR. &
     2104                                  xios_field_is_active("sollwad").OR.xios_field_is_active("sollwad0"))))  &
     2105            !!!--for now these fields are not in the XML files so they are omitted
     2106            !!!  xios_field_is_active("toplwai").OR.xios_field_is_active("sollwai") !))) &
     2107            swaero_diag=.TRUE.
     2108 
     2109         !--setting up swaerofree_diag to TRUE in XIOS case
     2110         IF (xios_field_is_active("SWdnSFCcleanclr").OR.xios_field_is_active("SWupSFCcleanclr").OR. &
     2111            xios_field_is_active("SWupTOAcleanclr").OR.xios_field_is_active("rsucsaf").OR.   &
     2112            xios_field_is_active("rsdcsaf") .OR. xios_field_is_active("LWdnSFCcleanclr").OR. &
     2113            xios_field_is_active("LWupTOAcleanclr")) &
     2114            swaerofree_diag=.TRUE.
     2115 
     2116         !--setting up dryaod_diag to TRUE in XIOS case
     2117         DO naero = 1, naero_tot-1
     2118          IF (xios_field_is_active("dryod550_"//name_aero_tau(naero))) dryaod_diag=.TRUE.
     2119         ENDDO
     2120         !
     2121         !--setting up ok_4xCO2atm to TRUE in XIOS case
     2122         IF (xios_field_is_active("rsut4co2").OR.xios_field_is_active("rlut4co2").OR. &
     2123            xios_field_is_active("rsutcs4co2").OR.xios_field_is_active("rlutcs4co2").OR. &
     2124            xios_field_is_active("rsu4co2").OR.xios_field_is_active("rsucs4co2").OR. &
     2125            xios_field_is_active("rsd4co2").OR.xios_field_is_active("rsdcs4co2").OR. &
     2126            xios_field_is_active("rlu4co2").OR.xios_field_is_active("rlucs4co2").OR. &
     2127            xios_field_is_active("rld4co2").OR.xios_field_is_active("rldcs4co2")) &
     2128            ok_4xCO2atm=.TRUE.
     2129       ENDIF
     2130       !$OMP BARRIER
     2131       CALL bcast(swaero_diag)
     2132       CALL bcast(swaerofree_diag)
     2133       CALL bcast(dryaod_diag)
     2134       CALL bcast(ok_4xCO2atm)
    20872135#endif
    2088 
    2089 
     2136       !
    20902137       CALL printflag( tabcntr0,radpas,ok_journe, &
    20912138            ok_instan, ok_region )
    20922139       !
    20932140       !
    2094        !
    20952141       ! Prescrire l'ozone dans l'atmosphere
    2096        !
    20972142       !
    20982143       !c         DO i = 1, klon
    … …  
    21502195#endif
    21512196       ENDIF
     2197       !
    21522198       IF (type_trac == 'repr') THEN
    21532199#ifdef REPROBUS
    … …  
    21982244          SFRWL(6)=3.02191470E-02
    21992245       END SELECT
    2200 
    2201 
    22022246       !albedo SB <<<
    22032247
    … …  
    23312375      ! RomP <<<
    23322376    ENDIF
    2333 
    23342377    !
    23352378    ! Ne pas affecter les valeurs entrees de u, v, h, et q
    … …  
    24082451
    24092452    tke0(:,:)=pbl_tke(:,:,is_ave)
    2410     !C Risi:Nombre de traceurs de l'eau: nqo
    2411     !  IF (nqtot.GE.3) THEN
    2412     IF (nqtot.GE.(nqo+1)) THEN
    2413        !     DO iq = 3, nqtot       
    2414 !       DO iq = nqo+1, nqtot 
    2415        ! CR: on modifie directement le code ici.
    2416        ! les isotopes ne sont pas dispatchés dans tr_seri, il faut les enlever.
    2417        ! on a prévu pour ça un tableau d'indice dans infotrac
    2418 #ifdef ISOVERIF
    2419        write(*,*) 'physiq 1971: nqtottr=',nqtottr
    2420 #endif
    2421        do itr=1,nqtottr
    2422          iq=itr_indice(itr)
     2453    IF (nqtot > nqo) THEN
     2454       ! water isotopes are not included in tr_seri
     2455       itr = 0
     2456       DO iq = 1, nqtot
     2457         IF(tracers(iq)%isH2Ofamily) CYCLE
     2458         itr = itr+1
    24232459!#ifdef ISOVERIF
    24242460!         write(*,*) 'physiq 1973: itr,iq=',itr,iq
    … …  
    24272463         DO  k = 1, klev
    24282464             DO  i = 1, klon
    2429                 tr_seri(i,k,itr) = qx(i,k,iq) ! modif C Risi
     2465                tr_seri(i,k,itr) = qx(i,k,iq)
    24302466             ENDDO
    2431           ENDDO !DO  k = 1, klev
    2432           !write(*,*) 'physiq 1980'
    2433        enddo !do itr=1,nqtottr
    2434 
    2435     ELSE !IF (nqtot.GE.(nqo+1)) THEN
    2436        DO k = 1, klev
    2437           DO i = 1, klon
    2438              tr_seri(i,k,1) = 0.0
    24392467          ENDDO
    24402468       ENDDO
    2441     ENDIF !IF (nqtot.GE.(nqo+1)) THEN
     2469    ELSE
     2470! DC: make sure the final "1" index was meant for 1st H2O phase (vapor) !!!
     2471!       tr_seri(:,:,strIdx(tracers(:)%name,'H2O'//phases_sep//'g')) = 0.0
     2472       tr_seri(:,:,strIdx(tracers(:)%name,'H2Ov')) = 0.0
     2473    ENDIF
    24422474!
    24432475! Temporary solutions adressing ticket #104 and the non initialisation of tr_ancien
    … …  
    24452477    IF (debut) THEN
    24462478      WRITE(lunout,*)' WARNING: tr_ancien initialised to tr_seri'
    2447 !      DO iq = nqo+1, nqtot
    2448 !           tr_ancien(:,:,iq-nqo)=tr_seri(:,:,iq-nqo)
    2449 !      ENDDO
    2450        ! modif CRisi:
    2451        do itr=1,nqtottr
     2479       itr = 0
     2480       do iq = 1, nqtot
     2481         IF(tracers(iq)%isH2Ofamily) CYCLE
     2482         itr = itr+1
    24522483         tr_ancien(:,:,itr)=tr_seri(:,:,itr)       
    24532484       enddo
    … …  
    25192550       d_qs_dyn2d(:)=(zx_tmp_fi2d(:)-prsw_ancien(:))/phys_tstep
    25202551       ! !! RomP >>>   td dyn traceur
    2521        IF (nqtot.GT.nqo) THEN     ! jyg
    2522 !          DO iq = nqo+1, nqtot      ! jyg
    2523           DO itr=1,nqtottr      ! C Risi modif directe
    2524               d_tr_dyn(:,:,itr)=(tr_seri(:,:,itr)-tr_ancien(:,:,itr))/phys_tstep ! jyg
    2525           ENDDO
    2526        ENDIF
     2552       IF (nqtot > nqo) d_tr_dyn(:,:,:)=(tr_seri(:,:,:)-tr_ancien(:,:,:))/phys_tstep
    25272553       ! !! RomP <<<
    25282554
    … …  
    26272653
    26282654       ! !! RomP >>>   td dyn traceur
    2629        IF (nqtot.GT.nqo) THEN                                       ! jyg
    2630 !          DO iq = nqo+1, nqtot                                      ! jyg
    2631 !              d_tr_dyn(:,:,iq-nqo)= 0.0                             ! jyg
    2632 ! Modif C Risi:
    2633           DO itr=1,nqtottr
    2634                 d_tr_dyn(:,:,itr)= 0.0
    2635           ENDDO
    2636        ENDIF
     2655       IF (nqtot > nqo) d_tr_dyn(:,:,:)= 0.0
    26372656       ! !! RomP <<<
    26382657       ancien_ok = .TRUE.
    … …  
    33633382            ENDDO
    33643383         ENDDO
    3365        ELSE !IF (iflag_wake>=1) THEN
     3384       ELSE
    33663385                t_w(:,:) = t_seri(:,:)
    33673386                q_w(:,:) = q_seri(:,:)
    … …  
    37523771          ! où i n'est pas un point convectif et donc ibas_con(i)=0
    37533772          ! c'est un pb indépendant des isotopes
    3754           if (ibas_con(i).gt.0) then
    3755             ema_pcb(i)  = paprs(i,ibas_con(i))
    3756           else ! if (ibas_con(i).gt.0) then
    3757               ema_pcb(i)  = 0.0
    3758           endif ! if (ibas_con(i).gt.0) then
    3759 
     3773          if (ibas_con(i) > 0) then
     3774             ema_pcb(i)  = paprs(i,ibas_con(i))
     3775          else
     3776             ema_pcb(i)  = 0.0
     3777          endif
    37603778       ENDDO
    37613779       DO i = 1, klon
    … …  
    44844502    ENDDO
    44854503
    4486    CALL  calcratqs(klon,klev,prt_level,lunout,        &
     4504    CALL  calcratqs(klon,klev,prt_level,lunout,        &
    44874505         iflag_ratqs,iflag_con,iflag_cld_th,pdtphys,  &
    44884506         ratqsbas,ratqshaut,ratqsp0, ratqsdp, &
    … …  
    44924510         pbl_tke(:,:,is_ave),tke_dissip_ave,l_mix_ave,wprime_ave,t2m,q2m,fm_therm, &
    44934511         ratqs,ratqsc,ratqs_inter)
    4494 
    44954512
    44964513    !
    … …  
    56085625                     ZLWFT0_i, ZFLDN0, ZFLUP0, &
    56095626                     ZSWFT0_i, ZFSDN0, ZFSUP0)
    5610           endif !ok_4xCO2atm
     5627          ENDIF !ok_4xCO2atm
    56115628       ENDIF ! aerosol_couple
    56125629       itaprad = 0
    … …  
    64856502#endif
    64866503! #ifdef ISO
    6487     !
    6488     !CR: nb de traceurs eau: nqo
    6489     !  IF (nqtot.GE.3) THEN
    6490     IF (nqtot.GE.(nqo+1)) THEN
    6491        !     DO iq = 3, nqtot
    6492 !       DO iq = nqo+1, nqtot ! modif C Risi
    6493         do itr=1,nqtottr
    6494          iq=itr_indice(itr)
    6495           DO  k = 1, klev
    6496              DO  i = 1, klon
    6497                 ! d_qx(i,k,iq) = ( tr_seri(i,k,iq-2) - qx(i,k,iq) ) / phys_tstep
    6498                 d_qx(i,k,iq) = ( tr_seri(i,k,iq-nqo) - qx(i,k,iq) ) / phys_tstep
    6499              ENDDO
     6504    ! DC: All iterations are cycled if nqtot==nqo, so no nqtot>nqo condition required
     6505    itr = 0
     6506    DO iq = 1, nqtot
     6507       IF(tracers(iq)%isH2Ofamily) CYCLE
     6508       itr = itr+1
     6509       DO  k = 1, klev
     6510          DO  i = 1, klon
     6511             d_qx(i,k,iq) = ( tr_seri(i,k,itr) - qx(i,k,iq) ) / phys_tstep
    65006512          ENDDO
    6501        ENDDO ! !do itr=1,nqtottr
    6502     ENDIF
     6513       ENDDO
     6514    ENDDO
    65036515    !
    65046516    !IM rajout diagnostiques bilan KP pour analyse MJO par Jun-Ichi Yano
    … …  
    65586570    CALL water_int(klon,klev,qs_ancien,zmasse,prsw_ancien)
    65596571    ! !! RomP >>>
    6560     !CR: nb de traceurs eau: nqo
    6561     IF (nqtot.GT.nqo) THEN
    6562        ! DO iq = nqo+1, nqtot ! modif C Risi
    6563        do itr=1,nqtottr
    6564           tr_ancien(:,:,itr) = tr_seri(:,:,itr)
    6565        ENDDO
    6566     ENDIF
     6572    IF (nqtot > nqo) tr_ancien(:,:,:) = tr_seri(:,:,:)
    65676573    ! !! RomP <<<
    65686574    !==========================================================================
    … …  
    67956801! ISO
    67966802
     6803    ! Disabling calls to the prt_alerte function
     6804    alert_first_call = .FALSE.
     6805   
    67976806    IF (lafin) THEN
    67986807       itau_phy = itau_phy + itap
Note: See TracChangeset for help on using the changeset viewer.