Ignore:
Timestamp:
Nov 30, 2021, 10:19:18 AM (4 years ago)
Author:
cmathe
Message:

GCM Mars: add meteoritic particles as CN for CO2 microphysics

Location:
trunk/LMDZ.MARS/libf/phymars
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.MARS/libf/phymars/callsedim_mod.F

    r2562 r2589  
    1919     &                      igcm_h2o_ice, igcm_hdo_ice,
    2020     &                      nuice_sed, nuice_ref,
    21      &                      igcm_ccnco2_mass,igcm_ccnco2_number,
    22      &                      igcm_co2_ice, igcm_stormdust_mass,
     21     &                      igcm_co2_ice, igcm_stormdust_mass,
    2322     &                      igcm_stormdust_number,igcm_topdust_mass,
    2423     &                      igcm_topdust_number,
    … …  
    162161      INTEGER,SAVE :: iccnco2_h2o_mass_ice ! index of tracer containing CCN number
    163162      INTEGER,SAVE :: iccnco2_h2o_mass_ccn ! index of tracer containing CCN number
     163      INTEGER,SAVE :: iccnco2_meteor_number ! index of tracer containing CCN number
     164      INTEGER,SAVE :: iccnco2_meteor_mass ! index of tracer containing CCN number
    164165      INTEGER,SAVE :: ico2_ice ! index of tracer containing CCN number
    165166
    … …  
    247248          iccnco2_h2o_mass_ccn=0
    248249          iccnco2_h2o_number=0
     250          iccnco2_meteor_mass=0
     251          iccnco2_meteor_number=0
    249252          ico2_ice=0
    250253          do iq=1,nq
    … …  
    275278              write(*,*)"callsedim: iccnco2_h2o_mass_ccn=",
    276279     &                 iccnco2_h2o_mass_ccn
     280            endif
     281            if (noms(iq).eq."ccnco2_meteor_number") then
     282              iccnco2_meteor_number=iq
     283              write(*,*)"callsedim: iccnco2_meteor_number=",
     284     &                   iccnco2_meteor_number
     285            endif
     286            if (noms(iq).eq."ccnco2_meteor_mass") then
     287              iccnco2_meteor_mass=iq
     288              write(*,*)"callsedim: iccnco2_meteor_mass=",
     289     &                 iccnco2_meteor_mass
    277290            endif
    278291          enddo
    … …  
    294307           end if
    295308          end if
     309          if (meteo_flux) then
     310            if((iccnco2_meteor_mass.eq.0).or.
     311     &         (iccnco2_meteor_number.eq.0)) then
     312            write(*,*) 'callsedim: error! could not identify'
     313            write(*,*) ' tracers for ccn co2 mass and number mixing'
     314            write(*,*) ' ratio and co2clouds are activated!'
     315            call abort_physic(modname,"missing co2 ccn tracers",1)
     316           end if
     317          end if
     318
    296319       ENDIF                    !of if (co2clouds)
    297320
    … …  
    432455     &        (iq.ne.iccnco2_h2o_number).and.
    433456     &        (iq.ne.iccnco2_h2o_mass_ice).and.
    434      &        (iq.ne.iccnco2_h2o_mass_ccn).and.  ! no sedim for gaz or CO2 clouds  (done in microtimestep)
     457     &        (iq.ne.iccnco2_h2o_mass_ccn).and.
     458     &        (iq.ne.iccnco2_meteor_mass).and.
     459     &        (iq.ne.iccnco2_meteor_number).and.  ! no sedim for gaz or CO2 clouds  (done in microtimestep)
    435460     &        iq .ne. igcm_hdo_ice) then !MVals: hdo is transported by h2o
    436461c -----------------------------------------------------------------

  • trunk/LMDZ.MARS/libf/phymars/co2cloud.F90

    r2562 r2589  
    100100                             igcm_ccn_mass, igcm_ccn_number, igcm_ccnco2_mass, igcm_ccnco2_number, &
    101101                             igcm_ccnco2_h2o_number, igcm_ccnco2_h2o_mass_ice, igcm_ccnco2_h2o_mass_ccn, rho_dust, &
    102                              nuiceco2_sed, nuiceco2_ref, r3n_q, rho_ice, nuice_sed
     102                             nuiceco2_sed, nuiceco2_ref, r3n_q, rho_ice, nuice_sed, igcm_ccnco2_meteor_mass, &
     103                             igcm_ccnco2_meteor_number
    103104
    104105  use newsedim_mod,    only: newsedim
    … …  
    591592       sum_subpdq(ig,l,igcm_co2) = sum_subpdq(ig,l,igcm_co2) + pdq(ig,l,igcm_co2)
    592593
     594       if (meteo_flux) then
     595        sum_subpdq(ig,l,igcm_ccnco2_meteor_number) = sum_subpdq(ig,l,igcm_ccnco2_meteor_number) + &
     596                                                  pdq(ig,l,igcm_ccnco2_meteor_number)
     597
     598        sum_subpdq(ig,l,igcm_ccnco2_meteor_mass) = sum_subpdq(ig,l,igcm_ccnco2_meteor_mass) + &
     599                                                    pdq(ig,l,igcm_ccnco2_meteor_mass)
     600       end if
    593601       if (co2useh2o) then
    594602         sum_subpdq(ig,l,igcm_h2o_ice) = sum_subpdq(ig,l,igcm_h2o_ice) + pdq(ig,l,igcm_h2o_ice)
    … …  
    598606         sum_subpdq(ig,l,igcm_ccn_number) = sum_subpdq(ig,l,igcm_ccn_number) + pdq(ig,l,igcm_ccn_number)
    599607
    600         sum_subpdq(ig,l,igcm_ccnco2_h2o_number) = sum_subpdq(ig,l,igcm_ccnco2_h2o_number) + &
    601                                                   pdq(ig,l,igcm_ccnco2_h2o_number)
    602 
    603         sum_subpdq(ig,l,igcm_ccnco2_h2o_mass_ice) = sum_subpdq(ig,l,igcm_ccnco2_h2o_mass_ice) + &
    604                                                     pdq(ig,l,igcm_ccnco2_h2o_mass_ice)
    605 
    606         sum_subpdq(ig,l,igcm_ccnco2_h2o_mass_ccn) = sum_subpdq(ig,l,igcm_ccnco2_h2o_mass_ccn) + &
    607                                                     pdq(ig,l,igcm_ccnco2_h2o_mass_ccn)
     608         sum_subpdq(ig,l,igcm_ccnco2_h2o_number) = sum_subpdq(ig,l,igcm_ccnco2_h2o_number) + &
     609                                                   pdq(ig,l,igcm_ccnco2_h2o_number)
     610
     611         sum_subpdq(ig,l,igcm_ccnco2_h2o_mass_ice) = sum_subpdq(ig,l,igcm_ccnco2_h2o_mass_ice) + &
     612                                                     pdq(ig,l,igcm_ccnco2_h2o_mass_ice)
     613
     614         sum_subpdq(ig,l,igcm_ccnco2_h2o_mass_ccn) = sum_subpdq(ig,l,igcm_ccnco2_h2o_mass_ccn) + &
     615                                                     pdq(ig,l,igcm_ccnco2_h2o_mass_ccn)
    608616       end if
    609617     end do ! ngrid
    … …  
    642650       end if
    643651
    644    ! ccnco2_h2o_number and masses
     652       ! ccnco2_meteor_number and ccnco2_meteor_mass
     653       if (meteo_flux) then
     654         if (((pq(ig,l,igcm_ccnco2_meteor_number)+(sum_subpdq(ig,l,igcm_ccnco2_meteor_number)+ &
     655               subpdqcloudco2(ig,l,igcm_ccnco2_meteor_number))*microtimestep)<=1.) .or. &
     656             (pq(ig,l,igcm_ccnco2_meteor_mass)+(sum_subpdq(ig,l,igcm_ccnco2_meteor_mass)+ &
     657               subpdqcloudco2(ig,l,igcm_ccnco2_meteor_mass))*microtimestep<=1e-20)) then
     658           subpdqcloudco2(ig,l,igcm_ccnco2_meteor_number) = - pq(ig,l,igcm_ccnco2_meteor_number)/microtimestep + 1. &
     659                                                    - sum_subpdq(ig,l,igcm_ccnco2_meteor_number)
     660             subpdqcloudco2(ig,l,igcm_ccnco2_meteor_mass) = - pq(ig,l,igcm_ccnco2_meteor_mass)/microtimestep + 1e-20 &
     661                                                  - sum_subpdq(ig,l,igcm_ccnco2_meteor_mass)
     662         end if
     663       end if
     664
     665       ! ccnco2_h2o_number and masses
    645666       if (co2useh2o) then
    646667         if (((pq(ig,l,igcm_ccnco2_h2o_number) + (sum_subpdq(ig,l,igcm_ccnco2_h2o_number) + &
    … …  
    688709        sum_subpdq(ig,l,igcm_co2) = sum_subpdq(ig,l,igcm_co2) + subpdqcloudco2(ig,l,igcm_co2)
    689710
     711        if (meteo_flux) then
     712          sum_subpdq(ig,l,igcm_ccnco2_meteor_mass) = sum_subpdq(ig,l,igcm_ccnco2_meteor_mass) + &
     713                                                     subpdqcloudco2(ig,l,igcm_ccnco2_meteor_mass)
     714
     715          sum_subpdq(ig,l,igcm_ccnco2_meteor_number) = sum_subpdq(ig,l,igcm_ccnco2_meteor_number) + &
     716                                                       subpdqcloudco2(ig,l,igcm_ccnco2_meteor_number)
     717        end if
    690718        if (co2useh2o) then
    691719          sum_subpdq(ig,l,igcm_h2o_ice) = sum_subpdq(ig,l,igcm_h2o_ice) + subpdqcloudco2(ig,l,igcm_h2o_ice)
    … …  
    723751
    724752          ! assure positive value of co2_ice mmr, ccnco2 number, ccnco2 mass
     753          ! meteoritic particle are considered like dust, rho_dust
    725754          Niceco2 = max(zqsed(ig,l,igcm_co2_ice), threshold)
    726755          Nccnco2 = max(zqsed(ig,l,igcm_ccnco2_number), threshold)
    … …  
    753782      zqsed0(:,:,igcm_ccnco2_number) = zqsed(:,:,igcm_ccnco2_number)
    754783
     784      if (meteo_flux) then
     785        zqsed0(:,:,igcm_ccnco2_meteor_mass) = zqsed(:,:,igcm_ccnco2_meteor_mass)
     786        zqsed0(:,:,igcm_ccnco2_meteor_number) = zqsed(:,:,igcm_ccnco2_meteor_number)
     787      end if
     788
    755789      if (co2useh2o) then
    756790        zqsed0(:,:,igcm_ccnco2_h2o_number) = zqsed(:,:,igcm_ccnco2_h2o_number)
    … …  
    794828      end do
    795829
     830      if (meteo_flux) then
     831        wq(:,:) = 0.
     832        ! for ccnco2_meteor_mass
     833        call newsedim(ngrid, nlay, ngrid*nlay, ngrid*nlay, microtimestep, pplev, masse, epaisseur, ztsed, &
     834                     rsedcloudco2, rhocloudco2t, zqsed(:,:,igcm_ccnco2_meteor_mass), wq, beta)
     835        do ig = 1, ngrid
     836          sum_subpdqs_sedccn(ig,igcm_ccnco2_meteor_mass) = sum_subpdqs_sedccn(ig,igcm_ccnco2_meteor_mass) + &
     837                                                           wq(ig,1)/ microtimestep
     838        end do
     839
     840        wq(:,:) = 0.
     841        ! for ccnco2_meteor_number
     842        call newsedim(ngrid, nlay, ngrid*nlay, ngrid*nlay,microtimestep, pplev, masse, epaisseur, ztsed, &
     843                     rsedcloudco2, rhocloudco2t, zqsed(:,:,igcm_ccnco2_meteor_number), wq, beta)
     844        do ig = 1, ngrid
     845          sum_subpdqs_sedccn(ig,igcm_ccnco2_meteor_number) = sum_subpdqs_sedccn(ig,igcm_ccnco2_meteor_number) + &
     846                                                             wq(ig,1) / microtimestep
     847        end do
     848      end if
    796849      ! for ccnco2_h2o_mass_ice
    797850      if (co2useh2o) then
    … …  
    835888          subpdqsed(ig,l,igcm_co2_ice) = ( zqsed(ig,l,igcm_co2_ice) - zqsed0(ig,l,igcm_co2_ice) ) / microtimestep
    836889
     890          if (meteo_flux) then
     891            subpdqsed(ig,l,igcm_ccnco2_meteor_mass) = ( zqsed(ig,l,igcm_ccnco2_meteor_mass) - &
     892                                                        zqsed0(ig,l,igcm_ccnco2_meteor_mass) ) / microtimestep
     893
     894            subpdqsed(ig,l,igcm_ccnco2_meteor_number) = ( zqsed(ig,l,igcm_ccnco2_meteor_number) - &
     895                                                          zqsed0(ig,l,igcm_ccnco2_meteor_number) ) / microtimestep
     896          end if
    837897          if (co2useh2o) then
    838898          subpdqsed(ig,l,igcm_ccnco2_h2o_number) = ( zqsed(ig,l,igcm_ccnco2_h2o_number) - &
    … …  
    855915
    856916          sum_subpdq(ig,l,igcm_co2_ice) = sum_subpdq(ig,l,igcm_co2_ice) + subpdqsed(ig,l,igcm_co2_ice)
     917          if (meteo_flux) then
     918            sum_subpdq(ig,l,igcm_ccnco2_meteor_mass) = sum_subpdq(ig,l,igcm_ccnco2_meteor_mass) + &
     919                                                       subpdqsed(ig,l,igcm_ccnco2_meteor_mass)
     920
     921            sum_subpdq(ig,l,igcm_ccnco2_meteor_number) = sum_subpdq(ig,l,igcm_ccnco2_meteor_number) + &
     922                                                         subpdqsed(ig,l,igcm_ccnco2_meteor_number)
     923          end if
    857924          if (co2useh2o) then
    858925            sum_subpdq(ig,l,igcm_ccnco2_h2o_mass_ice) = sum_subpdq(ig,l,igcm_ccnco2_h2o_mass_ice) + &
    … …  
    911978                                           pdq(ig,l,igcm_dust_number)
    912979
     980      if (meteo_flux) then
     981        pdqcloudco2(ig,l,igcm_ccnco2_meteor_mass) = ( sum_subpdq(ig,l,igcm_ccnco2_meteor_mass)/real(imicroco2) ) &
     982                                                    - pdq(ig,l,igcm_ccnco2_meteor_mass)
     983
     984        pdqcloudco2(ig,l,igcm_ccnco2_meteor_number) = ( sum_subpdq(ig,l,igcm_ccnco2_meteor_number) / real(imicroco2) ) &
     985                                                       - pdq(ig,l,igcm_ccnco2_meteor_number)
     986      end if
    913987      if (co2useh2o) then
    914988        pdqcloudco2(ig,l,igcm_h2o_ice) = ( sum_subpdq(ig,l,igcm_h2o_ice) / real(imicroco2) ) - &
    … …  
    9431017      Niceco2 = max (Niceco2, threshold)
    9441018
     1019      ! meteoritic particles are considered like dust, => rho_dust
    9451020      Nccnco2 = max( (pqeff(ig,l,igcm_ccnco2_number) + (pdq(ig,l,igcm_ccnco2_number) + &
    9461021                       pdqcloudco2(ig,l, igcm_ccnco2_number)) * ptimestep) &
    … …  
    10591134        Qext1bins2(ig,l) = Qext1bins2(ig,l) * co2cloudfrac(ig,l)
    10601135
     1136        if (meteo_flux) then
     1137          pdqcloudco2(ig,l,igcm_ccnco2_meteor_mass) = pdqcloudco2(ig,l,igcm_ccnco2_meteor_mass) * co2cloudfrac(ig,l)
     1138
     1139          pdqcloudco2(ig,l,igcm_ccnco2_meteor_number) = pdqcloudco2(ig,l,igcm_ccnco2_meteor_number) * co2cloudfrac(ig,l)
     1140        end if
    10611141        if (co2useh2o) then
    10621142          pdqcloudco2(ig,l,igcm_h2o_ice) = pdqcloudco2(ig,l,igcm_h2o_ice) * co2cloudfrac(ig,l)
    … …  
    10881168!----------------------------------------------------------------------------------------------------------------------!
    10891169  call WRITEDIAGFI(ngrid, "satuco2", "vap in satu", " ", 3, satuco2)
    1090 
    1091   call WRITEDIAGFI(ngrid, "precip_co2_ice", "surface deposition of co2 ice", "kg.m-2", 2,  pdqs_sedco2(:)*ptimestep)
    10921170
    10931171  call WRITEDIAGFI(ngrid, "precip_co2_ice_rate", "surface deposition rate of co2 ice", "kg.m-2.s-1", 2, pdqs_sedco2(:))

  • trunk/LMDZ.MARS/libf/phymars/improvedco2clouds_mod.F90

    r2562 r2589  
    7070                          nuice_sed, igcm_co2, igcm_co2_ice, igcm_ccnco2_mass, igcm_ccnco2_number, &
    7171                          igcm_ccnco2_h2o_mass_ice, igcm_ccnco2_h2o_mass_ccn, igcm_ccnco2_h2o_number, nuiceco2_sed, &
    72                           nuiceco2_ref
     72                          nuiceco2_ref, igcm_ccnco2_meteor_mass, igcm_ccnco2_meteor_number
    7373
    7474  use conc_mod,     only: mmean
    75 
     75  use time_phylmdz_mod, only: daysec
    7676  use nucleaco2_mod, only: nucleaco2
    7777  use datafile_mod, only: datadir
    … …  
    135135     threshold = 1e-30 ! limit value
    136136
    137   character(len=20), parameter:: &
    138      file_meteoritic_flux = 'Meteo_flux_Plane.dat'
     137  character(len=23), parameter:: &
     138     file_meteoritic_flux = 'Meteo_flux_Plane_v2.dat'
    139139!----------------------------------------------------------------------------------------------------------------------!
    140140!----2) Saved:
    … …  
    144144
    145145  double precision, save :: &
     146     pression_meteor(nlev_meteor),    &! pressure from meteoritic flux input file
    146147     meteor(nlev_meteor,nbin_meteor), &! Meteoritic flux read from file uMeteor
    147148     dev3                              ! 1. / ( sqrt(2.) * sigma_ice )
    … …  
    158159! ---Variables for meteoritic flux
    159160     ibin,   &! loop on nbin_meteor
    160      nelem,  &! nb of points during interpolation of meteoritic flux
    161      lebon1, &! index where P_meteor is the nearest of pplev(ig,l)
    162      lebon2, &! index where P_meteor is the nearest of pplev(ig,l+1)
     161     idx_min,&! index of min(diff_pressure)
    163162     read_ok  ! file uMeteor iostat
    164163
    … …  
    196195    n_aer(nbinco2_cld),                &! Radius used by the microphysical scheme (m)
    197196    m_aer(nbinco2_cld),                &! number concentration V-1 of particle/each size bin
     197    n_aer_meteor(nbinco2_cld),                &! Radius used by the microphysical scheme (m)
     198    m_aer_meteor(nbinco2_cld),                &! number concentration V-1 of particle/each size bin
    198199    n_aer_h2oice(nbinco2_cld),         &! mass mixing ratio of particle/each size bin
    199200    m_aer_h2oice(nbinco2_cld),         &! Same - for CO2 nucleation
    … …  
    203204    dN,                                &! number of particle of dust used as ccn
    204205    dM,                                &! mass of dN
     206    dN_meteor,                                &! number of particle of dust used as ccn
     207    dM_meteor,                                &! mass of dN
    205208    dNh2o,                             &! number of particle of h2o ice used as ccn
    206209    dMh2o,                             &! mass of dNh2o
    207210    dNN,                               &! min(dN,zq(ig,l,igcm_dust_number))
    208211    dMM,                               &! min(dM,zq(ig,l,igcm_dust_mass))
     212    dNN_meteor,                                &! number of particle of dust used as ccn
     213    dMM_meteor,                                &! mass of dN
    209214    dNNh2o,                            &! min(dNNh2o,zq(ig,l,igcm_ccn_number))
    210215    dMh2o_ice,                         &! min(dMh2o*zq(ig,l,igcm_h2o_ice)*ratioh2o_ccn, zq(ig,l,igcm_h2o_ice))
    211216    dMh2o_ccn,                         &! min(dMh2o_ccn,zq(ig,l,igcm_ccn_mass))
    212217    rate(nbinco2_cld),                 &! nucleation rate
     218    rate_meteor(nbinco2_cld),                 &! nucleation rate
    213219    rateh2o(nbinco2_cld),              &! nucleation rate for h2o
    214220    rho_ice_co2T,                      &! density of co2 ice Temperature-dependent
    … …  
    217223    Proba,                             &! 1.0 - exp(-1.*microtimestep*rate(i))
    218224    Probah2o,                          &! 1.0 - exp(-1.*microtimestep*rateh2o(i))
    219     mtemp(nbinco2_cld),                &! sum(meteor(lebon1:lebon2,ibin))
    220     pression_meteor(nlev_meteor),      &! pressure from meteoritic flux input file
    221     ltemp1(nlev_meteor),               &! abs(pression_meteor(:)-pplev(ig,l))
    222     ltemp2(nlev_meteor),               &! abs(pression_meteor(:)-pplev(ig,l+1))
    223     meteor_ccn(ngrid,nlay,nbinco2_cld)  ! nbinco2_cld = 100
     225    Proba_meteor,                          &! 1.0 - exp(-1.*microtimestep*rateh2o(i))
     226    diff_pressure(nlev_meteor),        &! abs(pression_meteor(:)-pplev(ig,l))
     227    meteor_ccn(ngrid,nlay,nbinco2_cld)  !
    224228
    225229  logical :: &
    … …  
    264268        ! les mêmes 100 bins size que la distri nuclea : on touche pas
    265269        do ibin = 1, nbin_meteor
    266           read(uMeteor,'(F12.6)') meteor(i,ibin)
     270          read(uMeteor,'(F12.6)')meteor(i,ibin)
    267271        end do
    268272      end do
    … …  
    276280!----------------------------------------------------------------------------------------------------------------------!
    277281  rdust(:,:) = 0.
    278   meteor_ccn(:,:,:) = 0.
     282  meteor_ccn(:,:,:) = 0d0
    279283  rice(:,:) = 1.e-8
    280284  riceco2(:,:) = 1.e-11
    … …  
    302306! 3. Bonus: additional meteoritic particles for nucleation
    303307!----------------------------------------------------------------------------------------------------------------------!
    304   ! TODO: instead of intepolation, used only the nearest pplev(ig,l)
    305308  if (meteo_flux) then
    306309    do l = 1, nlay
    … …  
    308311        masse(ig,l) = (pplev(ig,l) - pplev(ig,l+1)) / g
    309312
    310         ltemp1 = abs(pression_meteor(:)-pplev(ig,l))
    311         ltemp2 = abs(pression_meteor(:)-pplev(ig,l+1))
    312 
    313         lebon1 = minloc(ltemp1,DIM=1)
    314         lebon2 = minloc(ltemp2,DIM=1)
    315 
    316         nelem = lebon2-lebon1+1.
    317 
    318         mtemp(:) = 0d0
    319 
    320         do ibin = 1, nbin_meteor
    321           mtemp(ibin) = sum(meteor(lebon1:lebon2,ibin))
    322         end do
    323 
     313        diff_pressure(:) = abs(pression_meteor(:)-pplev(ig,l))
     314
     315        idx_min = minloc(diff_pressure, DIM=1)
    324316        ! Par kg air csi par m carre, x epaisseur/masse pour par kg/air. Check original unit with J. Plane
    325         meteor_ccn(ig,l,:)=mtemp(:)/nelem/masse(ig,l)
     317        meteor_ccn(ig,l,:) = meteor(idx_min,:)/masse(ig,l)
    326318      end do
    327319    end do
    … …  
    387379
    388380          ! Now increment CCN tracers and update dust tracers
    389           dNN = min(dN,zq(ig,l,igcm_dust_number)) ! dNN est devenu DN
    390           dMM = min(dM,zq(ig,l,igcm_dust_mass))  ! idem dans le min
     381          dNN = min(dN, zq(ig,l,igcm_dust_number)) ! dNN est devenu DN
     382          dMM = min(dM, zq(ig,l,igcm_dust_mass))  ! idem dans le min
    391383
    392384          zq(ig,l,igcm_ccnco2_mass) = zq(ig,l,igcm_ccnco2_mass) + dMM /tauscaling(ig)
    393 
    394385          zq(ig,l,igcm_ccnco2_number) = zq(ig,l,igcm_ccnco2_number) + dNN /tauscaling(ig)
    395386
    396387          zq(ig,l,igcm_dust_mass) = zq(ig,l,igcm_dust_mass) - dMM /tauscaling(ig)
    397 
    398388          zq(ig,l,igcm_dust_number) = zq(ig,l,igcm_dust_number) - dNN /tauscaling(ig)
    399 
    400           ! Ajout meteor_ccn particles aux particules de poussière background
    401           if (meteo_flux) then
     389        end if ! of if No > 1e-30
     390
     391        ! Ajout meteor_ccn particles aux particules de poussière background
     392        if (meteo_flux) then
     393            n_aer_meteor(1:nbinco2_cld) = 0d0
     394            m_aer_meteor(1:nbinco2_cld) = 0d0
    402395            do i = 1, nbinco2_cld
    403               n_aer(i) = n_aer(i) + meteor_ccn(ig,l,i)
    404 
    405               m_aer(i) = m_aer(i) + (4./3.) * pi * rho_dust * meteor_ccn(ig,l,i) * rad_cldco2(i)**3
    406              end do
     396              n_aer_meteor(i) = meteor_ccn(ig,l,i)
     397              m_aer_meteor(i) = (4./3.) * pi * rho_dust * meteor_ccn(ig,l,i) * rad_cldco2(i)**3
     398            end do
     399            ! Call to nucleation routine
     400            rate_meteor(1:nbinco2_cld) = 0d0
     401            call nucleaco2(dble(pco2), zt(ig,l), dble(satu), n_aer_meteor, rate_meteor, vo2co2, mtetaco2)
     402
     403            dN_meteor = 0.
     404            dM_meteor = 0.
     405            do i = 1, nbinco2_cld
     406              Proba_meteor = 1.0 - exp(-1.*microtimestep*rate_meteor(i))
     407              dN_meteor = dN_meteor + n_aer_meteor(i) * Proba_meteor
     408              dM_meteor = dM_meteor + m_aer_meteor(i) * Proba_meteor
     409            end do
     410            ! Now increment CCN tracers and update dust tracers
     411            zq(ig,l,igcm_ccnco2_meteor_mass)   = zq(ig,l,igcm_ccnco2_meteor_mass) + dM_meteor
     412            zq(ig,l,igcm_ccnco2_meteor_number) = zq(ig,l,igcm_ccnco2_meteor_number) + dN_meteor
     413
     414            zq(ig,l,igcm_ccnco2_mass)   = zq(ig,l,igcm_ccnco2_mass) + dM_meteor
     415            zq(ig,l,igcm_ccnco2_number) = zq(ig,l,igcm_ccnco2_number) + dN_meteor
    407416          end if
    408         end if ! of if No > 1e-30
    409417
    410418        ! Same but with h2o particles as CCN only if co2useh2o = .true.
    … …  
    539547          zq(ig,l,igcm_ccnco2_mass) = 0.
    540548          zq(ig,l,igcm_ccnco2_number) = 0.
     549          if (meteo_flux) then
     550            zq(ig,l,igcm_dust_mass) = zq(ig,l,igcm_dust_mass) - zq(ig,l,igcm_ccnco2_meteor_mass)
     551            zq(ig,l,igcm_dust_number) = zq(ig,l,igcm_dust_number) - zq(ig,l,igcm_ccnco2_meteor_number)
     552            zq(ig,l,igcm_ccnco2_meteor_mass) = 0.
     553            zq(ig,l,igcm_ccnco2_meteor_number) = 0.
     554          end if
    541555          if (co2useh2o) then
    542556            zq(ig,l,igcm_dust_mass) = zq(ig,l,igcm_dust_mass) - zq(ig,l,igcm_ccnco2_h2o_mass_ccn) - &
    … …  
    572586  subpdqcloudco2(:,:,igcm_dust_number) = ( zq(:,:,igcm_dust_number) - zq0(:,:,igcm_dust_number) )   / microtimestep
    573587
     588  if (meteo_flux) then
     589    subpdqcloudco2(:,:,igcm_ccnco2_meteor_mass) = (zq(:,:,igcm_ccnco2_meteor_mass)-zq0(:,:,igcm_ccnco2_meteor_mass) &
     590                                                   )/microtimestep
     591
     592    subpdqcloudco2(:,:,igcm_ccnco2_meteor_number) = ( zq(:,:,igcm_ccnco2_meteor_number) - &
     593                                                      zq0(:,:,igcm_ccnco2_meteor_number) )/microtimestep
     594  end if
     595
    574596  if (co2useh2o) then
    575597    subpdqcloudco2(:,:,igcm_h2o_ice) = ( zq(:,:,igcm_h2o_ice) - zq0(:,:,igcm_h2o_ice) )  / microtimestep

  • trunk/LMDZ.MARS/libf/phymars/initracer.F

    r2562 r2589  
    6464      igcm_ccnco2_mass=0
    6565      igcm_ccnco2_number=0
     66      igcm_ccnco2_meteor_mass=0
     67      igcm_ccnco2_meteor_number=0
    6668      igcm_ccnco2_h2o_mass_ice=0
    6769      igcm_ccnco2_h2o_mass_ccn=0
    … …  
    428430              count=count+1
    429431           endif
     432           if (meteo_flux) then
     433             if (noms(iq).eq."ccnco2_meteor_mass") then
     434                igcm_ccnco2_meteor_mass=iq
     435                count=count+1
     436             endif
     437             if (noms(iq).eq."ccnco2_meteor_number") then
     438                igcm_ccnco2_meteor_number=iq
     439                count=count+1
     440             endif
     441           end if
    430442           if (co2useh2o) then
    431443           if (noms(iq).eq."ccnco2_h2o_number") then
    … …  
    821833     &            "a co2 tracer !"
    822834             call abort_physic("initracer","co2 clouds issue",1)
    823           endif
     835          end if
    824836          if (igcm_co2_ice .eq. 0) then
    825837             write(*,*) "initracer: error !!"
    … …  
    827839     &            "a co2_ice tracer !"
    828840             call abort_physic("initracer","co2 clouds issue",1)
    829           endif
     841          end if
     842          ! TODO: ajouter les tests avec co2useh2o, et meteo_flux, ainsi que les traceurs ccnco2_number
    830843       endif
    831844 

  • trunk/LMDZ.MARS/libf/phymars/physiq_mod.F

    r2578 r2589  
    3333     &                      igcm_ccnco2_h2o_mass_ccn,
    3434     &                      igcm_ccnco2_h2o_number,
     35     &                      igcm_ccnco2_meteor_mass,
     36     &                      igcm_ccnco2_meteor_number,
    3537     &                      rho_ice_co2,nuiceco2_sed,nuiceco2_ref,
    3638     &                      igcm_dust_mass, igcm_dust_number, igcm_h2o2,
    … …  
    16761678     &              zdtcloudco2(1:ngrid,1:nlayer)
    16771679
    1678 
    16791680! increment dust and ccn masses and numbers
    16801681! We need to check that we have Nccn & Ndust > 0
    … …  
    16991700     &                             + zdqcloudco2(:,:,igcm_ccnco2_number)
    17001701
     1702             if (meteo_flux) then
     1703               pdq(:,:,igcm_ccnco2_meteor_mass) =
     1704     &                 pdq(:,:,igcm_ccnco2_meteor_mass) +
     1705     &                 zdqcloudco2(:,:,igcm_ccnco2_meteor_mass)
     1706
     1707               pdq(:,:,igcm_ccnco2_meteor_number) =
     1708     &                 pdq(:,:,igcm_ccnco2_meteor_number)
     1709     &                 + zdqcloudco2(:,:,igcm_ccnco2_meteor_number)
     1710             end if
    17011711!Update water ice clouds values as well
    17021712             if (co2useh2o) then
    … …  
    18001810     &           - pq(:,:,igcm_dust_number)/ptimestep + 1.e-30
    18011811             end where
     1812             if (meteo_flux) then
     1813             where (pq(:,:,igcm_ccnco2_meteor_mass) +
     1814     &              ptimestep*pdq(:,:,igcm_ccnco2_meteor_mass) < 0.)
     1815               pdq(:,:,igcm_ccnco2_meteor_mass) =
     1816     &            - pq(:,:,igcm_ccnco2_meteor_mass)/ptimestep + 1.e-30
     1817               pdq(:,:,igcm_ccnco2_meteor_number) =
     1818     &            - pq(:,:,igcm_ccnco2_meteor_number)/ptimestep + 1.e-30
     1819             end where
     1820             where (pq(:,:,igcm_ccnco2_meteor_number) +
     1821     &              ptimestep*pdq(:,:,igcm_ccnco2_meteor_number) < 0.)
     1822               pdq(:,:,igcm_ccnco2_meteor_mass) =
     1823     &            - pq(:,:,igcm_ccnco2_meteor_mass)/ptimestep + 1.e-30
     1824               pdq(:,:,igcm_ccnco2_meteor_number) =
     1825     &            - pq(:,:,igcm_ccnco2_meteor_number)/ptimestep + 1.e-30
     1826             end where
     1827             end if
    18021828      END IF ! of IF (co2clouds)
    18031829
    … …  
    30233049         call WRITEDIAGFI(ngrid,"emis","Surface emissivity","w.m-1",2,
    30243050     &                  emis)
     3051         if (co2useh2o) then
     3052           call WRITEDIAGFI(ngrid,'ccnqco2_h2o_m_ice',
     3053     &                      'CCNco2_h2o_mass_ice mmr',
     3054     &                    'kg.kg-1',3,zq(:,:,igcm_ccnco2_h2o_mass_ice))
     3055
     3056           call WRITEDIAGFI(ngrid,'ccnqco2_h2o_m_ccn',
     3057     &                      'CCNco2_h2o_mass_ccn mmr',
     3058     &                   'kg.kg-1',3,zq(:,:,igcm_ccnco2_h2o_mass_ccn))
     3059
     3060           call WRITEDIAGFI(ngrid,'ccnNco2_h2o','CCNco2_h2o number',
     3061     &                   'part.kg-1',3,zq(:,:,igcm_ccnco2_h2o_number))
     3062         end if
     3063
     3064           if (meteo_flux) then
     3065          call WRITEDIAGFI(ngrid,'ccnqco2_meteor','CCNco2_meteor mmr',
     3066     &                  'kg.kg-1',3,zq(:,:,igcm_ccnco2_meteor_mass))
     3067
     3068        call WRITEDIAGFI(ngrid,'ccnNco2_meteor','CCNco2_meteor number',
     3069     &                'part.kg-1',3,zq(:,:,igcm_ccnco2_meteor_number))
     3070           end if
     3071
    30253072        end if ! of if (co2clouds)
    30263073      end if ! of if (tracer.and.(igcm_co2.ne.0))

  • trunk/LMDZ.MARS/libf/phymars/tracer_mod.F90

    r2578 r2589  
    6565      integer,save :: igcm_ccnco2_mass   ! CCN (dust and/or water ice) for CO2 mass mixing ratio
    6666      integer,save :: igcm_ccnco2_number ! CCN (dust and/or water ice) for CO2 number mixing ratio
     67      integer,save :: igcm_ccnco2_meteor_mass   ! CCN (dust and/or water ice) for CO2 mass mixing ratio
     68      integer,save :: igcm_ccnco2_meteor_number ! CCN (dust and/or water ice) for CO2 number mixing ratio
    6769      integer,save :: igcm_ccnco2_h2o_mass_ice   ! CCN (dust and/or water ice) for CO2 mass mixing ratio
    6870      integer,save :: igcm_ccnco2_h2o_mass_ccn   ! CCN (dust and/or water ice) for CO2 mass mixing ratio
Note: See TracChangeset for help on using the changeset viewer.