Ignore:
Timestamp:
Apr 7, 2021, 3:16:53 PM (4 years ago)
Author:
cmathe
Message:

Mars GCM:
co2_ice as scatterer in radiative transfert. Need co2clouds and

activeco2ice .eqv. True. Files involved:

  • aeropacity_mod.F
  • callradite_mod.F
  • physiq_mod.F
  • updatereffrad_mod.F
  • suaer.F90
  • determine co2_ice density from temperature. Used in riceco2 computation.

Files involved:

  • co2cloud.F90
  • improvedco2clouds_mod.F90
  • updaterad.F90
  • updatereffrad_mod.F
  • co2condens_mod4micro.F: variable initialization
  • initracer.F: add nuiceco2_ref = 0.2
  • phyredem.F: remove co2_ice from qsurf since co2_ice => co2ice
  • watercloud_mod.F: tiny typo

CM

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

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.MARS/libf/phymars/aeropacity_mod.F

    r2449 r2494  
    1010
    1111      SUBROUTINE aeropacity(ngrid,nlayer,nq,zday,pplay,pplev,ls,
    12      &    pq,tauscaling,dust_rad_adjust,tau_pref_scenario,tau_pref_gcm,
    13      &    tau,taucloudtes,aerosol,dsodust,reffrad,
     12     &    pq,pt,tauscaling,dust_rad_adjust,tau_pref_scenario,
     13     &    tau_pref_gcm,tau,taucloudtes,aerosol,dsodust,reffrad,
    1414     &    QREFvis3d,QREFir3d,omegaREFir3d,
    1515     &    totstormfract,clearatm,dsords,dsotop,
     
    3737      use dust_param_mod, only: odpref, freedust
    3838      use dust_scaling_mod, only: compute_dustscaling
     39      use density_co2_ice_mod, only: density_co2_ice
    3940       IMPLICIT NONE
    4041c=======================================================================
     
    7475                         ! of the atmospheric layers
    7576      REAL,INTENT(IN) ::  pq(ngrid,nlayer,nq) ! tracers
     77      REAL,INTENT(IN) ::  pt(ngrid,nlayer) !temperature
    7678      REAL,INTENT(OUT) :: tau_pref_scenario(ngrid) ! prescribed dust column
    7779                          ! visible opacity at odpref from scenario
     
    132134      REAL topdust0(ngrid)
    133135
     136!      -- CO2 clouds
     137       real CLFtotco2
     138       real taucloudco2vis(ngrid)
     139       real taucloudco2tes(ngrid)
     140       real totcloudco2frac(ngrid) ! a mettre en (in) [CM]
     141       double precision :: rho_ice_co2
    134142#ifdef DUSTSTORM
    135143!! Local dust storms
     
    174182        ! identify scatterers that are dust
    175183        naerdust=0
     184        iaerdust(1:naerkind) = 0
     185        nqdust(1:nq) = 0
    176186        DO iaer=1,naerkind
    177187          txt=name_iaer(iaer)
     
    192202          ENDIF
    193203        ENDDO
    194 
    195204        IF (water.AND.activice) THEN
    196205          i_ice=igcm_h2o_ice
     
    236245                                       !avoid unrealistic values due to constant lifting
    237246        ENDIF
    238 
    239247
    240248#ifndef DUSTSTORM
     
    488496c       1. Initialization
    489497        aerosol(1:ngrid,1:nlayer,iaer) = 0.
    490         taucloudvis(1:ngrid) = 0.
    491         taucloudtes(1:ngrid) = 0.
     498        taucloudco2vis(1:ngrid) = 0.
     499        taucloudco2tes(1:ngrid) = 0.
    492500c       2. Opacity calculation
    493501        ! NO CLOUDS
    494502        IF (clearsky) THEN
    495             aerosol(1:ngrid,1:nlayer,iaer) =1.e-9
     503            aerosol(1:ngrid,1:nlayer,iaer) = 1.e-9
    496504        ! CLOUDSs
    497505        ELSE ! else (clearsky)
    498           DO ig=1, ngrid
    499             DO l=1,nlayer
     506          DO ig = 1, ngrid
     507            DO l = 1, nlayer
     508              call density_co2_ice(dble(pt(ig,l)), rho_ice_co2)
     509
    500510              aerosol(ig,l,iaer) = max(1E-20,
    501511     &          (  0.75 * QREFvis3d(ig,l,iaer) /
    502      &          ( rho_ice * reffrad(ig,l,iaer) )  ) *
     512     &          ( rho_ice_co2 * reffrad(ig,l,iaer) )  ) *
    503513     &          pq(ig,l,i_co2ice) *
    504514     &          ( pplev(ig,l) - pplev(ig,l+1) ) / g
    505515     &                              )
    506               taucloudvis(ig) = taucloudvis(ig) + aerosol(ig,l,iaer)
    507               taucloudtes(ig) = taucloudtes(ig) + aerosol(ig,l,iaer)*
     516              taucloudco2vis(ig) = taucloudco2vis(ig)
     517     &                             + aerosol(ig,l,iaer)
     518              taucloudco2tes(ig) = taucloudco2tes(ig)
     519     &                             + aerosol(ig,l,iaer) *
    508520     &          QREFir3d(ig,l,iaer) / QREFvis3d(ig,l,iaer) *
    509521     &          ( 1.E0 - omegaREFir3d(ig,l,iaer) )
     
    511523          ENDDO
    512524          ! SUB-GRID SCALE CLOUDS
    513           IF (CLFvarying) THEN
     525          IF (CLFvaryingCO2) THEN
    514526             DO ig=1, ngrid
    515                 DO l=1,nlayer-1
    516                    CLFtot  = max(totcloudfrac(ig),0.01)
     527                DO l= 1, nlayer-1
     528                   CLFtotco2  = max(totcloudco2frac(ig),0.01)
    517529                   aerosol(ig,l,iaer)=
    518      &                    aerosol(ig,l,iaer)/CLFtot
     530     &                    aerosol(ig,l,iaer)/CLFtotco2
    519531                   aerosol(ig,l,iaer) =
    520532     &                    max(aerosol(ig,l,iaer),1.e-9)
    521533                ENDDO
    522534             ENDDO
    523           ENDIF ! end (CLFvarying)             
     535          ENDIF ! end (CLFvaryingCO2)
    524536        ENDIF ! end (clearsky)
    525537
  • trunk/LMDZ.MARS/libf/phymars/callradite_mod.F

    r2459 r2494  
    427427     &                rdust,rstormdust,rtopdust,rice,nuice,
    428428     &                reffrad,nueffrad, riceco2, nuiceco2,
    429      &                pq,tauscaling,tau,pplay)
     429     &                pq,tauscaling,tau,pplay, pt)
    430430c     Computing 3D scattering parameters:
     431      gVIS3d(:,:,:,:) = 0.
    431432      CALL aeroptproperties(ngrid,nlayer,reffrad,nueffrad,
    432433     &                      QVISsQREF3d,omegaVIS3d,gVIS3d,
     
    436437c     Computing aerosol optical depth in each layer:
    437438      CALL aeropacity(ngrid,nlayer,nq,zday,pplay,pplev,ls,
    438      &    pq,tauscaling,dust_rad_adjust,tau_pref_scenario,tau_pref_gcm,
    439      &    tau,taucloudtes,aerosol,dsodust,reffrad,
     439     &    pq,pt,tauscaling,dust_rad_adjust,tau_pref_scenario,
     440     &    tau_pref_gcm,tau,taucloudtes,aerosol,dsodust,reffrad,
    440441     &    QREFvis3d,QREFir3d,omegaREFir3d,
    441442     &    totstormfract,clearatm,dsords,dsotop,
    442443     &    alpha_hmons,nohmons,
    443444     &    clearsky,totcloudfrac)
    444 
    445445c     Starting loop on sub-domain
    446446c     ----------------------------
    447 
     447      zgVIS3d(:,:,:,:) = 0.
     448      zfluxd_sw(:,:,:) = 0.
     449      zfluxu_sw(:,:,:) = 0.
     450      zQVISsQREF3d(:,:,:,:) = 0.
     451      zomegaVIS3d(:,:,:,:) = 0.
    448452      DO jd=1,ndomain
    449453        ig0=(jd-1)*ndomainsz
     
    479483         enddo
    480484        enddo
    481 
     485        zplev(:,:) = 0.
    482486        do l=1,nlaylte+1
    483487         do ig = 1,nd
     
    485489         enddo
    486490        enddo
    487 
     491        zdp(:,:) = 0.
     492       
    488493        do l=1,nlaylte
    489494         do ig = 1,nd
     
    494499         enddo
    495500        enddo
    496 
     501        zaerosol(:,:,:) = 0.
    497502        do n=1,naerkind
    498503          do l=1,nlaylte
     
    502507          enddo
    503508        enddo
    504 
     509        zalbedo(:,:) = 0.
    505510        do j=1,2
    506511          do ig = 1,nd
     
    546551c          1370 W.m-2 is the solar constant at 1 AU.
    547552           cste_mars=1370./(dist_sol*dist_sol)
    548 
     553           zzdtsw(:,:) = 0.
    549554           call swmain ( nd, nflev,
    550555     S     cste_mars, zalbedo,
     
    552557     S     zzdtsw, zfluxd_sw, zfluxu_sw,
    553558     &     zQVISsQREF3d,zomegaVIS3d,zgVIS3d)
    554 
    555559c       ------------------------------------------------------------
    556560c       Un-spliting output variable from sub-domain input variables
     
    564568        enddo
    565569
     570        ptlev(:, :) = 0.
    566571        do l=1,nlaylte+1
    567572         do ig = 1,nd
     
    588593         enddo
    589594      endif
    590 
    591595c     Output for debugging if lwrite=T
    592596c     --------------------------------
  • trunk/LMDZ.MARS/libf/phymars/co2cloud.F90

    r2460 r2494  
    108108  use tracer_mod,      only: igcm_co2, igcm_co2_ice, igcm_dust_mass, igcm_dust_number, igcm_h2o_ice, &
    109109                             igcm_ccn_mass, igcm_ccn_number, igcm_ccnco2_mass, igcm_ccnco2_number, rho_dust, &
    110                              nuiceco2_sed, nuiceco2_ref, rho_ice_co2, r3n_q, rho_ice, nuice_sed
     110                             nuiceco2_sed, nuiceco2_ref, r3n_q, rho_ice, nuice_sed
    111111
    112112  use newsedim_mod,    only: newsedim
     
    114114  use datafile_mod,    only: datadir
    115115
    116   use improvedCO2clouds_mod, only: improvedCO2clouds, density_co2_ice
     116  use improvedCO2clouds_mod, only: improvedCO2clouds
    117117
    118118#ifndef MESOSCALE
     
    209209     dev2,                   &! 1. / ( sqrt(2.) * sigma_iceco2 )
    210210     Qext1bins(nbinco2_cld), &! Extinction coefficients for rb_cldco2 radius of CO2 ice particles
     211     Qextv1mic(var_dim_qext), &
     212     radv(var_dim_qext),      &    ! radius of CO2 ice at 1 µm (read from file_qext)
    211213     rb_cldco2(nbinco2_cld+1) ! boundary values of each rad_cldco2 bin (m)
    212 
    213214  logical, save :: &
    214215     firstcall = .true. ! Used to compute saved variables
     
    273274     tau1mic(ngrid),                &! CO2 ice column opacity at 1µm
    274275     Qext1bins2(ngrid,nlay),        &! CO2 ice opacities
    275      radv(var_dim_qext),            &! radius of CO2 ice at 1 µm (read from file_qext)
    276      Qextv1mic(var_dim_qext),       &! extinction coefficient of CO2 ice at 1 µm (read from file_qext)
    277276! ---For Saturation Index computation
    278277     rho,                           &! background density
     
    285284! ---Misc
    286285     myT,                           &! Temperature scalar for co2 density computation
    287      tcond(ngrid,nlay),             &! CO2 condensation temperature
    288      rho_ice_co2T(ngrid,nlay)        ! T-dependant CO2 ice density
     286     tcond(ngrid,nlay)               ! CO2 condensation temperature
    289287
    290288  logical :: &
     
    419417  epaisseur(1:ngrid,1:nlay) = 0.
    420418  masse(1:ngrid,1:nlay) = 0.
    421 
     419  riceco2(1:ngrid, 1:nlay) = 0.
    422420  zqsed0(1:ngrid,1:nlay,1:nq) = 0.
    423421  sum_subpdqs_sedco2(1:ngrid) = 0.
     
    611609                          subpdqcloudco2, subpdtcloudco2, nq, tauscaling, mem_Mccn_co2, mem_Mh2o_co2, mem_Nccn_co2, &
    612610                          rb_cldco2, sigma_iceco2, dev2)
     611
    613612   do l = 1, nlay
    614613     do ig = 1, ngrid
     
    682681          zqsed(ig,l,:) = pqeff(ig,l,:) + sum_subpdq(ig,l,:) * microtimestep
    683682
    684           ! density of co2 ice
    685           call density_co2_ice(dble(ztsed(ig,l)), rho_ice_co2T(ig,l))
    686 
    687683          ! assure positive value of co2_ice mmr, ccnco2 number, ccnco2 mass
    688684          Niceco2 = max(zqsed(ig,l,igcm_co2_ice), threshold)
     
    692688          ! Get density cloud and co2 ice particle radius
    693689          if (Niceco2.ne.0d0) then
    694             call updaterice_microco2(dble(Niceco2), dble(Qccnco2), dble(Nccnco2), tauscaling(ig), riceco2(ig,l), rhocloudco2t(ig,l))
     690            call updaterice_microco2(dble(Niceco2), dble(Qccnco2), dble(Nccnco2), ztsed(ig,l), tauscaling(ig), &
     691                                     riceco2(ig,l), rhocloudco2t(ig,l))
    695692          else
    696693            riceco2(ig,l) = 0.
     
    776773    pdqs_sedco2(ig) = sum_subpdqs_sedco2(ig) / real(imicroco2)
    777774  end do
    778 
    779775  ! temperature tendency (T.s-1)
    780776  do l = 1, nlay
     
    835831      myT = pteff(ig,l) + (pdt(ig,l)+pdtcloudco2(ig,l))*ptimestep
    836832
    837       ! compute density of co2 ice
    838       call density_co2_ice(myT, rho_ice_co2T(ig,l))
    839 
    840       rho_ice_co2 = rho_ice_co2T(ig,l) ! rho_ice_co2 is shared by tracer_mod and used in updaterice
    841 
    842833      ! Compute particle size
    843       call updaterice_microco2(dble(Niceco2), dble(Qccnco2), dble(Nccnco2), tauscaling(ig), riceco2(ig,l), rhocloudco2(ig,l))
     834      call updaterice_microco2(dble(Niceco2), dble(Qccnco2), dble(Nccnco2), myT, tauscaling(ig), riceco2(ig,l), &
     835                               rhocloudco2(ig,l))
    844836
    845837     ! Compute opacities
  • trunk/LMDZ.MARS/libf/phymars/co2condens_mod4micro.F90

    r2456 r2494  
    302302
    303303! Compute without microphysics
     304 diff_zcondicea(1:ngrid, 1:nlayer) = 0.
     305 diff_zfallice(1:ngrid) = 0.
    304306 do l =1, nlayer
    305307   do ig = 1, ngrid
     
    319321
    320322 ztcond_2(:,nlayer+1)=ztcond_2(:,nlayer)
    321  zfallice_2(:,:) = 0
    322 
     323 zfallice_2(:,:) = 0.
     324 zcondicea_2(:,:) = 0.
    323325 do l = nlayer , 1, -1
    324326   do ig = 1, ngrid
  • trunk/LMDZ.MARS/libf/phymars/improvedco2clouds_mod.F90

    r2456 r2494  
    77! Contains subroutines:
    88!     - improvedco2clouds_mod: nucleation
    9 !
    10 !     - density_co2_ice: compute density of co2 ice particle
    119!======================================================================================================================!
    1210module improvedco2clouds_mod
     
    7270  use tracer_mod,   only: igcm_dust_mass, igcm_dust_number, rho_dust, igcm_h2o_ice, igcm_ccn_mass, igcm_ccn_number, &
    7371                          nuice_sed, igcm_co2, igcm_co2_ice, igcm_ccnco2_mass, igcm_ccnco2_number, nuiceco2_sed, &
    74                           rho_ice_co2, nuiceco2_ref
     72                          nuiceco2_ref
    7573
    7674  use conc_mod,     only: mmean
    7775
    7876  use datafile_mod, only: datadir
     77
     78  use density_co2_ice_mod, only: density_co2_ice
    7979
    8080  implicit none
     
    211211    rate(nbinco2_cld),                 &! nucleation rate
    212212    rateh2o(nbinco2_cld),              &! nucleation rate for h2o
    213     rho_ice_co2T(ngrid,nlay),          &! density of co2 ice
     213    rho_ice_co2T,                      &! density of co2 ice Temperature-dependent
    214214    riceco2(ngrid,nlay),               &! CO2 ice mean radius (m)
    215215    vrat_cld,                          &! Volume ratio
     
    232232    firstcall = .false.
    233233
    234 !   Volume of a co2 molecule (m3)
    235     vo1co2 = m0co2 / dble(rho_ice_co2)
    236 
    237234!   Variance of the ice and CCN distributions
    238235    sigma_ice = sqrt(log(1.+nuice_sed))
     
    287284
    288285  ! pteff temperature layer; sum_subpdt dT.s-1
     286  zt(1:ngrid,1:nlay) = 0.
    289287  zt(:,:) = pteff(:,:) + sum_subpdt(:,:) * microtimestep
    290288
     
    294292
    295293  zq0(:,:,:) = zq(:,:,:)
     294  zqsat(:,:) = 0.
    296295!----------------------------------------------------------------------------------------------------------------------!
    297296! 2. Compute saturation
     
    338337
    339338      !T-dependant CO2 ice density
    340       call density_co2_ice(dble(zt(ig,l)), rho_ice_co2T(ig,l))
    341 
    342       vo2co2 = m0co2 / dble(rho_ice_co2T(ig,l))
    343       rho_ice_co2 = rho_ice_co2T(ig,l)
     339      call density_co2_ice(dble(zt(ig,l)), rho_ice_co2T)
     340
     341      vo2co2 = m0co2 / rho_ice_co2T
    344342!----------------------------------------------------------------------------------------------------------------------!
    345343! 4.1 Nucleation
     
    421419          ! Call to nucleation routine
    422420          call nucleaCO2(dble(pco2), zt(ig,l), dble(satu), n_aer, rate, n_aer_h2oice, rad_h2oice, rateh2o, vo2co2)
     421
    423422          dN = 0.
    424423          dM = 0.
     
    492491      if (zq(ig,l,igcm_ccnco2_number) * tauscaling(ig) + threshold >= 1) then
    493492
    494         call updaterice_microco2(dble(zq(ig,l,igcm_co2_ice)), dble(zq(ig,l,igcm_ccnco2_mass)), dble(zq(ig,l,igcm_ccnco2_number)), &
    495                                  tauscaling(ig), riceco2(ig,l), rhocloudco2(ig,l))
    496 
     493        call updaterice_microco2(dble(zq(ig,l,igcm_co2_ice)), dble(zq(ig,l,igcm_ccnco2_mass)), &
     494                                 dble(zq(ig,l,igcm_ccnco2_number)), zt(ig,l), tauscaling(ig), riceco2(ig,l), &
     495                                 rhocloudco2(ig,l))
    497496        Ic_rice = 0.
    498497
     
    594593  end subroutine improvedCO2clouds
    595594
    596 
    597 !**********************************************************************************************************************!
    598 !**********************************************************************************************************************!
    599 
    600 
    601 !======================================================================================================================!
    602 ! SUBROUTINE: density_co2_ice =========================================================================================!
    603 !======================================================================================================================!
    604 ! Subject:
    605 !---------
    606 !   Compute co2 ice particles density
    607 !======================================================================================================================!
    608   subroutine density_co2_ice(temperature, density)
    609 
    610   implicit none
    611 
    612   double precision, intent(in) :: &
    613      temperature
    614 
    615   double precision, intent(out) :: &
    616      density
    617 
    618   density = 1000. * (1.72391 - 2.53e-4*temperature - 2.87e-6*temperature*temperature)
    619 
    620   end subroutine density_co2_ice
    621 
    622595end module improvedCO2clouds_mod
    623596
  • trunk/LMDZ.MARS/libf/phymars/initracer.F

    r2463 r2494  
    497497                      ! of the water-ice size distribution
    498498
     499      nuiceco2_ref=0.2    !C.M. Effective variance nueff of the
     500                          ! co2-ice size distribution
    499501      if (doubleq) then
    500502c       "doubleq" technique
  • trunk/LMDZ.MARS/libf/phymars/phyredem.F90

    r2417 r2494  
    302302      endif ! of if (txt.eq."hdo_vap")
    303303
     304      ! co2_ice has been added to co2ice in co2condens4micro
     305      if (txt.eq."co2_ice") then
     306        write(*,*)"physdem1: skipping co2_ice tracer"
     307        cycle
     308      end if     
     309
    304310      call put_field(trim(txt),"tracer on surface",qsurf(:,iq),time)
    305311    enddo
     
    317323     call put_field("dv_nonoro_gwd","Meridional wind tendency due to GW",dv_nonoro_gwd,time)
    318324  endif
     325
    319326  ! Close file
    320327  call close_restartphy
  • trunk/LMDZ.MARS/libf/phymars/physiq_mod.F

    r2491 r2494  
    13201320
    13211321c        Calling vdif (Martian version WITH CO2 condensation)
     1322         dwatercap_dif(:) = 0.
     1323         zcdh(:) = 0.
     1324         zcdv(:) = 0.
    13221325         CALL vdifc(ngrid,nlayer,nq,co2ice,zpopsk,
    13231326     $        ptimestep,capcal,lwrite,
     
    16401643
    16411644c Temperature variation due to latent heat release
    1642 c            if (activice) then  !Maybe create activice_co2 ?
    16431645               pdt(1:ngrid,1:nlayer) =
    16441646     &              pdt(1:ngrid,1:nlayer) +
    16451647     &              zdtcloudco2(1:ngrid,1:nlayer)! --> in co2condens
    1646 c            endif
    16471648           
    16481649
     
    19421943     $                  zdqssed_co2,zcondicea_co2microp,
    19431944     &                  zdtcloudco2)
     1945         ! no scavenging yet
     1946         zdqsc(:,:) = 0.
    19441947        else
    19451948          CALL co2condens(ngrid,nlayer,nq,ptimestep,
     
    28102813c        WRITEDIAGFI can ALSO be called from any other subroutines
    28112814c        for any variables !!
    2812 c         call WRITEDIAGFI(ngrid,"emis","Surface emissivity","w.m-1",2,
    2813 c     &                  emis)
     2815         call WRITEDIAGFI(ngrid,"emis","Surface emissivity","w.m-1",2,
     2816     &                  emis)
    28142817c        call WRITEDIAGFI(ngrid,"pplay","Pressure","Pa",3,zplay)
    28152818c        call WRITEDIAGFI(ngrid,"pplev","Pressure","Pa",3,zplev)
     
    33463349            call WRITEDIAGFI(ngrid,"q15um","15 um cooling","K/s",
    33473350     $           3,zdtnlte)
    3348             call WRITEDIAGFI(ngrid,"quv","UV heating","K/s",
    3349      $           3,zdteuv)
    3350             call WRITEDIAGFI(ngrid,"cond","Thermal conduction","K/s",
    3351      $           3,zdtconduc)
    33523351            call WRITEDIAGFI(ngrid,"qnir","NIR heating","K/s",
    33533352     $           3,zdtnirco2)
     
    34863485          endif
    34873486
    3488 c         CALL writeg1d(ngrid,nlayer,zt,'temp','K')
    3489 c         CALL writeg1d(ngrid,nlayer,riceco2,'riceco2','m')
    3490 c         CALL writeg1d(ngrid,nlayer,satuco2,'satuco2','satu')
    3491          
    3492          
    3493 c         call WRITEDIAGFI(ngrid,"satuco2","vap in satu","kg/kg",1,
    3494 c     &        satuco2)
    3495 c         call WRITEDIAGFI(ngrid,"riceco2","ice radius","m"
    3496 c     &        ,1,riceco2)
    3497 ! or output in diagfi.nc (for testphys1d)
    34983487         call WRITEDIAGFI(ngrid,'ps','Surface pressure','Pa',0,ps)
    34993488         call WRITEDIAGFI(ngrid,'temp','Temperature ',
  • trunk/LMDZ.MARS/libf/phymars/suaer.F90

    r2447 r2494  
    164164        longrefvis(iaer) = 0.67E-6  ! 1.5um OMEGA/MEx
    165165!       Reference wavelength in the infrared:
    166         longrefir(iaer) = 1E-6   ! 1 micron a vérifier si on prends cette valeur
     166        longrefir(iaer) = 4.26E-6  ! 2347 cm-1 OMEGA/MEx
    167167!==================================================================
    168168        CASE("stormdust_doubleq") aerkind   ! Two-moment scheme for stormdust - radiative properties
  • trunk/LMDZ.MARS/libf/phymars/updaterad.F90

    r1996 r2494  
    100100!============================================================================
    101101!============================================================================
    102 
    103 subroutine updaterice_microco2(qice,qccn,nccn,coeff,rice,rhocloudco2)
    104 use tracer_mod, only: rho_dust, rho_ice_co2
    105 USE comcstfi_h, only:  pi
     102subroutine updaterice_microco2(qice, qccn, nccn, temperature, coeff, rice, rhocloudco2)
     103use tracer_mod, only: rho_dust
     104use comcstfi_h, only: pi
     105use density_co2_ice_mod, only: density_co2_ice
     106
    106107implicit none
    107108!CO2 clouds parameter update by CL and JA 09/16
     
    109110DOUBLE PRECISION, intent(in)  :: qice,qccn,nccn
    110111real, intent(in)  :: coeff         ! this coeff is tauscaling if microphy = T (possibly ccn_factor^-1 otherwise)
     112real, intent(in)  :: temperature ! temperature effective for density co2_ice computation
    111113real, intent(out) :: rhocloudco2 ! rhocloud is needed for sedimentation and is also a good diagnostic variable
    112114double precision, intent(out) :: rice
    113115real nccn_true,qccn_true ! nombre et masse de CCN
    114    
     116double precision :: rho_ice_co2T ! density co2_ice Temperature-dependent
     117
    115118nccn_true = max(nccn * coeff, 1.e-30)
    116119qccn_true = max(qccn * coeff, 1.e-30)
    117120
    118 
    119   rhocloudco2 = (qice *rho_ice_co2 + qccn_true*rho_dust) / (qice + qccn_true)
    120 
    121   rhocloudco2 = min(max(rhocloudco2,rho_ice_co2),rho_dust)
    122 
    123   rice = (qice + qccn_true) * 0.75 / pi / rhocloudco2 / nccn_true
    124 
    125   if (rice .le. r3iceco2min) then !r3icemin radius power 3 ?
    126     rice = riceco2min
    127   else if (rice .ge. r3iceco2max) then !r3icemin radius power 3 ?
    128     rice = riceco2max
    129   else
    130     rice = rice**(1./3.) ! here rice is always positive
    131   endif
     121call density_co2_ice(dble(temperature), rho_ice_co2T)
     122
     123rhocloudco2 = (qice * rho_ice_co2T + qccn_true*rho_dust) / (qice + qccn_true)
     124
     125rhocloudco2 = min(max(rhocloudco2,rho_ice_co2T), rho_dust)
     126
     127rice = (qice + qccn_true) * 0.75 / pi / rhocloudco2 / nccn_true
     128
     129rhocloudco2 = (qice * rho_ice_co2T + qccn_true*rho_dust) / (qice + qccn_true)
     130
     131if (rice .le. r3iceco2min) then !r3icemin radius power 3 ?
     132  rice = riceco2min
     133else if (rice .ge. r3iceco2max) then !r3icemin radius power 3 ?
     134  rice = riceco2max
     135else
     136  rice = rice**(1./3.) ! here rice is always positive
     137endif
    132138
    133139
  • trunk/LMDZ.MARS/libf/phymars/updatereffrad_mod.F

    r2459 r2494  
    88     &                rdust,rstormdust,rtopdust,rice,nuice,
    99     &                reffrad,nueffrad, riceco2, nuiceco2,
    10      &                pq,tauscaling,tau,pplay)
     10     &                pq,tauscaling,tau,pplay, pt)
    1111       USE updaterad, ONLY: updaterdust, updaterice_micro,
    1212     &                      updaterice_microco2, updaterice_typ
     
    1818     &                       ref_r0, igcm_dust_submicron,
    1919     &                       igcm_stormdust_mass,igcm_stormdust_number,
    20      &                       igcm_topdust_mass,igcm_topdust_number
     20     &                       igcm_topdust_mass,igcm_topdust_number,
     21     &                       rho_ice
    2122       USE dimradmars_mod, only: nueffdust,naerkind,
    2223     &            name_iaer,
     
    7374      double precision, INTENT(out) :: riceco2(ngrid,nlayer) ! co2 ice radius
    7475      REAL, INTENT(out) :: nuiceco2(ngrid,nlayer)
    75 
     76      REAL, INTENT(in) :: pt(ngrid,nlayer) ! temperature
    7677     
    7778c     Local variables:
     
    185186
    186187
    187 c       1.3 CO2-ice particles --- in progress.....
    188 c       ----------------------- j'ai juste copié-collé la partie 1.2
     188c       1.3 CO2-ice particles
    189189
    190190        IF (co2clouds.AND.activeco2ice) THEN
    191 
    192 c    At firstcall, the true number and true mass of cloud condensation nuclei are not known.
    193 c    Indeed it is scaled on the prescribed dust opacity via a 'tauscaling' coefficient
    194 c    computed after radiative transfer. If tauscaling is not in startfi, we make an assumption for its value.
    195 
    196           IF (firstcall) THEN
    197             !IF (minval(tauscaling).lt.0) tauscaling(:) = 1.e-3 ! default value when non-read in startfi is -1
    198             !IF (freedust)                tauscaling(:) = 1.    ! if freedust, enforce no rescaling at all
    199             firstcall = .false.
    200           ENDIF
    201191          DO l=1,nlayer
    202192            DO ig=1,ngrid
     
    204194     &                              dble(pq(ig,l,igcm_ccnco2_mass)),
    205195     &                              dble(pq(ig,l,igcm_ccnco2_number)),
     196     &                              pt(ig,l),
    206197     &                              tauscaling(ig),riceco2(ig,l),
    207198     &                              rhocloudco2(ig,l))
     
    266257          DO l=1,nlayer
    267258            DO ig=1,ngrid
    268               reffrad(ig,l,iaer)=riceco2(ig,l)*(1.+nuiceco2_ref)
     259              reffrad(ig,l,iaer)=real(riceco2(ig,l))*(1.+nuiceco2_ref)
    269260              nueffrad(ig,l,iaer)=nuiceco2_ref
    270261            ENDDO
  • trunk/LMDZ.MARS/libf/phymars/watercloud_mod.F

    r2437 r2494  
    681681     &      ,"kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
    682682     &      igcm_ccn_mass))
    683       call WRITEDIAGFI(ngrid,"pdqccnN2","pdqcloudccnN apres
    684      &      microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
     683      call WRITEDIAGFI(ngrid,"pdqccnN2","pdqcloudccnN apres "//
     684     &      "microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
    685685     &      igcm_ccn_number))
    686       call WRITEDIAGFI(ngrid,"pdqdust2", "pdqclouddust apres
    687      &      microphysique","kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
     686      call WRITEDIAGFI(ngrid,"pdqdust2", "pdqclouddust apres "//
     687     &      "microphysique","kg/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
    688688     &      igcm_dust_mass))
    689       call WRITEDIAGFI(ngrid,"pdqdustN2", "pdqclouddustN apres
    690      &      microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
     689      call WRITEDIAGFI(ngrid,"pdqdustN2", "pdqclouddustN apres "//
     690     &      "microphysique","nb/kg.s-1",3,pdqcloud(1:ngrid,1:nlay,
    691691     &      igcm_dust_number))
    692692c=======================================================================
Note: See TracChangeset for help on using the changeset viewer.