Ignore:
Timestamp:
Apr 13, 2020, 4:48:01 PM (5 years ago)
Author:
emillour
Message:

Mars GCM:

  • Adding HCO+ ion chemistry.

FGG

Location:
trunk/LMDZ.MARS/libf/aeronomars
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.MARS/libf/aeronomars/calchim_mod.F90

    r2273 r2284  
    2222                            igcm_coplus, igcm_cplus, igcm_nplus,          &
    2323                            igcm_noplus, igcm_n2plus, igcm_hplus,         &
    24                             igcm_hco2plus, igcm_elec, mmol
     24                            igcm_hco2plus, igcm_hcoplus, igcm_elec, mmol
    2525
    2626      use conc_mod, only: mmean ! mean molecular mass of the atmosphere
    … …  
    141141      integer,save :: i_hplus=0
    142142      integer,save :: i_hco2plus=0
     143      integer,save :: i_hcoplus=0
    143144      integer,save :: i_elec=0
    144145
    … …  
    542543            if (i_hco2plus /= 0) then
    543544               write(*,*) "calchim: Error: HCO2+ is present, but O2+ is not!!!"
     545               write(*,*) "Both must be in traceur.def if ion chemistry wanted"
     546            endif
     547         endif
     548         i_hcoplus=igcm_hcoplus
     549         if(ionchem) then
     550            if (i_hcoplus == 0) then
     551               write(*,*) "calchim: Error, no HCO+ tracer !!!"
     552               write(*,*) "HCO+ is needed if O2+ is in traceur.def"
     553               stop
     554            else
     555               nbq = nbq + 1
     556               niq(nbq) = i_hcoplus
     557            end if
     558         else
     559            if (i_hcoplus /= 0) then
     560               write(*,*) "calchim: Error: HCO+ is present, but O2+ is not!!!"
    544561               write(*,*) "Both must be in traceur.def if ion chemistry wanted"
    545562            endif
    … …  
    630647            ! set number of reactions, depending on ion chemistry or not
    631648            if (ionchem) then
    632                nb_reaction_4_max = 61   ! set number of bimolecular reactions
     649               nb_reaction_4_max = 67   ! set number of bimolecular reactions
    633650               nb_reaction_3_max = 6    ! set number of quadratic reactions
    634651               nquench           = 9    ! set number of quenching + heterogeneous reactions

  • trunk/LMDZ.MARS/libf/aeronomars/inichim_newstart.F90

    r1918 r2284  
    117117      igcm_hplus=0
    118118      igcm_hco2plus=0
     119      igcm_hcoplus=0
    119120      igcm_elec=0
    120121
    … …  
    330331           igcm_hco2plus = iq
    331332           mmol(igcm_hco2plus) = 45.
     333           count = count + 1
     334        end if
     335        if (noms(iq) == "hcoplus") then
     336           igcm_hcoplus = iq
     337           mmol(igcm_hcoplus) = 29.
    332338           count = count + 1
    333339        end if
    … …  
    593599              .or. igcm_cplus == 0 .or. igcm_nplus == 0 .or. igcm_noplus == 0    &
    594600              .or. igcm_n2plus == 0 .or. igcm_hplus == 0 .or. igcm_hco2plus == 0 &
    595               .or. igcm_elec == 0) then
     601              .or. igcm_hcoplus == 0 .or. igcm_elec == 0) then
    596602            write(*,*)'inichim_newstart error:'
    597603            write(*,*)'if co2plus is in traceur.def, all other ions must also be'
    598604            write(*,*)'o2plus, oplus, coplus, cplus, nplus, noplus, n2plus'
    599             write(*,*)'hplus, hco2plus and elec'
     605            write(*,*)'hplus, hco2plus, hcoplus, and elec'
    600606            write(*,*)'stop'
    601607            stop
    … …  
    616622                  pq(i,j,l,igcm_hplus)    = 0.
    617623                  pq(i,j,l,igcm_hco2plus) = 0.
     624                  pq(i,j,l,igcm_hcoplus)  = 0.
    618625                  pq(i,j,l,igcm_elec)     = 0.
    619626               end do
    … …  
    633640         qsurf(1:ngrid,igcm_hplus)    = 0.
    634641         qsurf(1:ngrid,igcm_hco2plus) = 0.
     642         qsurf(1:ngrid,igcm_hcoplus)  = 0.
    635643         qsurf(1:ngrid,igcm_elec)     = 0.
    636644
    … …  
    640648              .or. igcm_cplus /= 0 .or. igcm_nplus /= 0 .or. igcm_noplus /= 0    &
    641649              .or. igcm_n2plus /= 0 .or. igcm_hplus /= 0 .or. igcm_hco2plus /= 0 &
    642               .or. igcm_elec /= 0) then
     650              .or. igcm_hcoplus /= 0 .or. igcm_elec /= 0) then
    643651            write(*,*)'inichim_newstart error:'
    644652            write(*,*)'some ions are in traceur.def, but not co2plus'

  • trunk/LMDZ.MARS/libf/aeronomars/photochemistry.F90

    r2273 r2284  
    112112integer,parameter :: i_hplus   = 26
    113113integer,parameter :: i_hco2plus= 27
    114 integer,parameter :: i_elec    = 28
     114integer,parameter :: i_hcoplus = 28
     115integer,parameter :: i_elec    = 29
    115116
    116117integer :: ilay
    … …  
    154155               i_n, i_n2d, i_no, i_no2, i_n2, i_co2plus,            &
    155156               i_oplus, i_o2plus, i_noplus, i_coplus, i_cplus,      &
    156                i_n2plus, i_nplus, i_hplus, i_hco2plus, i_elec)
     157               i_n2plus, i_nplus, i_hplus, i_hco2plus, i_hcoplus,   &
     158               i_elec)
    157159   firstcall = .false.
    158160end if
    … …  
    168170               i_oplus, i_o2plus, i_noplus, i_coplus, i_cplus, &
    169171               i_n2plus, i_nplus, i_hplus, i_hco2plus,         &
    170                i_elec, dens, rm, c)
     172               i_hcoplus, i_elec, dens, rm, c)
    171173
    172174!===================================================================
    … …  
    363365      if(c(ilev,i_elec).ne.c(ilev,i_co2plus)+c(ilev,i_oplus)+c(ilev,i_o2plus)+&
    364366           c(ilev,i_noplus)+c(ilev,i_coplus)+c(ilev,i_cplus)+c(ilev,i_n2plus)+&
    365            c(ilev,i_nplus)+c(ilev,i_hplus)+c(ilev,i_hco2plus)) then
     367           c(ilev,i_nplus)+c(ilev,i_hplus)+c(ilev,i_hco2plus)+                &
     368           c(ilev,i_hcoplus)) then
    366369         c(ilev,i_elec) = c(ilev,i_co2plus)+c(ilev,i_oplus)+c(ilev,i_o2plus)+ &
    367370              c(ilev,i_noplus)+c(ilev,i_coplus)+c(ilev,i_cplus)+              &
    368371              c(ilev,i_n2plus)+c(ilev,i_nplus)+c(ilev,i_hplus)+               &
    369               c(ilev,i_hco2plus)
     372              c(ilev,i_hco2plus)+c(ilev,i_hcoplus)
    370373         !      write(*,*)'photochemistry/359'
    371374         !      write(*,*)'Forcing charge neutrality at ilev,',ilev,' ig=',ig
    … …  
    393396               i_oplus, i_o2plus, i_noplus, i_coplus, i_cplus, &
    394397               i_n2plus, i_nplus, i_hplus, i_hco2plus,         &
    395                i_elec, dens, c)
     398               i_hcoplus, i_elec, dens, c)
    396399contains
    397400
    … …  
    592595                           i013, i014, i015, i016, i017, i018, i019,   &
    593596                           i020, i021, i022, i023, i024, i025, i026,   &
    594                            i027, i028, i029, i030,                     &
     597                           i027, i028, i029, i030, i031, i032, i033,   &
     598                           i034, i035, i036,                           &
    595599                           h001, h002, h003, h004, h005
    596600
    … …  
    13601364         nb_reaction_4 = nb_reaction_4 + 1
    13611365         v_4(:,nb_reaction_4) = i030(:)
     1366
     1367!---     i031: HCO2+ + O -> HCO+ + O2
     1368
     1369!        UMIST
     1370
     1371         i031(:) = 1.e-9
     1372         nb_reaction_4 = nb_reaction_4 + 1
     1373         v_4(:,nb_reaction_4) = i031(:)
     1374
     1375!---     i032: HCO2+ + CO -> HCO+ + CO2
     1376
     1377!        UMIST, from Prassad & Huntress 1980
     1378
     1379         i032(:) = 7.8e-10
     1380         nb_reaction_4 = nb_reaction_4 + 1
     1381         v_4(:,nb_reaction_4) = i032(:)
     1382
     1383!---     i033: H+ + CO2 -> HCO+ + O
     1384
     1385!        UMIST, from Smith et al., Int. J. Mass Spectrom. Ion Proc., 117, 457-473(1992)
     1386
     1387         i033(:) = 3.5e-9
     1388         nb_reaction_4 = nb_reaction_4 + 1
     1389         v_4(:,nb_reaction_4) = i033(:)
     1390
     1391
     1392!---     i034: CO2+ + H -> HCO+ + O
     1393
     1394!        Seen in Fox 2015, from Borodi et al., Int. J. Mass Spectrom. 280, 218-225, 2009
     1395
     1396         i034(:) = 4.5e-10
     1397         nb_reaction_4 = nb_reaction_4 + 1
     1398         v_4(:,nb_reaction_4) = i034(:)
     1399
     1400!---     i035: CO+ + H2 -> HCO+ + H
     1401
     1402         !UMIST, from Scott et al., J. Chem. Phys., 106, 3982-3987(1997)
     1403
     1404         i035(:) = 7.5e-10
     1405         nb_reaction_4 = nb_reaction_4 + 1
     1406         v_4(:,nb_reaction_4) = i035(:)
     1407
     1408!---     i036: HCO+ + e- -> CO + H
     1409
     1410         !UMIST, from Mitchell, Phys. Rep., 186, 215 (1990)
     1411
     1412         i036(:) = 2.4e-7 *((300./t_elect(:))**0.69)
     1413         nb_reaction_4 = nb_reaction_4 + 1
     1414         v_4(:,nb_reaction_4) = i036(:)
    13621415
    13631416      end if   !ionchem
    … …  
    15721625                   i_n, i_n2d, i_no, i_no2, i_n2, i_co2plus,            &
    15731626                   i_oplus, i_o2plus, i_noplus, i_coplus, i_cplus,      &
    1574                    i_n2plus, i_nplus, i_hplus, i_hco2plus, i_elec)
     1627                   i_n2plus, i_nplus, i_hplus, i_hco2plus, i_hcoplus,   &
     1628                   i_elec)
    15751629
    15761630!================================================================
    … …  
    15961650           i_n, i_n2d, i_no, i_no2, i_n2,                   &
    15971651           i_co2plus, i_oplus, i_o2plus, i_noplus, i_coplus, &
    1598            i_cplus, i_n2plus, i_nplus, i_hplus, i_hco2plus, i_elec
     1652           i_cplus, i_n2plus, i_nplus, i_hplus, i_hco2plus, &
     1653           i_hcoplus, i_elec
    15991654integer, intent(in) :: nb_reaction_3_max
    16001655                       ! number of quadratic reactions
    … …  
    24412496   indice_4(nb_reaction_4) = z4spec(1.0, i_hco2plus, 1.0, i_elec, 1.0, i_h, 1.0, i_co2)
    24422497
     2498
     2499!===========================================================
     2500!      i031 : HCO2+ + O -> HCO+ + O2
     2501!===========================================================
     2502
     2503   nb_reaction_4 = nb_reaction_4 + 1
     2504
     2505   indice_4(nb_reaction_4) = z4spec(1.0, i_hco2plus, 1.0, i_o, 1.0, i_hcoplus, 1.0, i_o2)
     2506
     2507
     2508!===========================================================
     2509!      i032 : HCO2+ + CO -> HCO+ + CO2
     2510!===========================================================
     2511
     2512   nb_reaction_4 = nb_reaction_4 + 1
     2513   indice_4(nb_reaction_4) = z4spec(1.0, i_hco2plus, 1.0, i_co, 1.0, i_hcoplus, 1.0, i_co2)
     2514
     2515
     2516!===========================================================
     2517!      i033 : H+ + CO2 -> HCO+ + O
     2518!===========================================================
     2519
     2520   nb_reaction_4 = nb_reaction_4 + 1
     2521   indice_4(nb_reaction_4) = z4spec(1.0, i_hplus, 1.0, i_co2, 1.0, i_hcoplus, 1.0, i_o)
     2522
     2523
     2524!===========================================================
     2525!      i034 : CO2+ + H -> HCO+ + O
     2526!===========================================================
     2527
     2528   nb_reaction_4 = nb_reaction_4 + 1
     2529   indice_4(nb_reaction_4) = z4spec(1.0, i_co2plus, 1.0, i_h, 1.0, i_hcoplus, 1.0, i_o)
     2530
     2531
     2532!===========================================================
     2533!      i035 : CO+ + H2 -> HCO+ + H
     2534!===========================================================
     2535
     2536   nb_reaction_4 = nb_reaction_4 + 1
     2537   indice_4(nb_reaction_4) = z4spec(1.0, i_coplus, 1.0, i_h2, 1.0, i_hcoplus, 1.0, i_h)
     2538
     2539
     2540!===========================================================
     2541!      i036 : HCO+ + e- -> CO + H
     2542!===========================================================
     2543
     2544   nb_reaction_4 = nb_reaction_4 + 1
     2545   indice_4(nb_reaction_4) = z4spec(1.0, i_hcoplus, 1.0, i_elec, 1.0, i_co, 1.0, i_h)
     2546
     2547
    24432548end if    !ionchem
    24442549
    … …  
    25182623                           i_co2plus, i_oplus, i_o2plus, i_noplus,   &
    25192624                           i_coplus, i_cplus, i_n2plus, i_nplus,     &
    2520                            i_hplus, i_hco2plus, i_elec, dens, rm, c)
     2625                           i_hplus, i_hco2plus, i_hcoplus, i_elec,   &
     2626                           dens, rm, c)
    25212627       
    25222628!*****************************************************************
    … …  
    25292635     &                      igcm_noplus, igcm_coplus, igcm_cplus,        &
    25302636     &                      igcm_n2plus, igcm_nplus, igcm_hplus,         &
    2531      &                      igcm_hco2plus, igcm_elec
     2637     &                      igcm_hco2plus, igcm_hcoplus, igcm_elec
    25322638
    25332639      implicit none
    … …  
    25492655                 i_n, i_n2d, i_no, i_no2, i_n2,                      &
    25502656                 i_co2plus, i_oplus, i_o2plus, i_noplus, i_coplus,   &
    2551                  i_cplus, i_n2plus, i_nplus, i_hplus, i_hco2plus, i_elec
     2657                 i_cplus, i_n2plus, i_nplus, i_hplus, i_hco2plus,    &
     2658                 i_hcoplus, i_elec
    25522659
    25532660      real :: zycol(nlayer,nq)      ! volume mixing ratios in the gcm
    … …  
    26832790               stop
    26842791            endif
     2792            if (igcm_hcoplus == 0) then
     2793               write(*,*) "gcmtochim: Error; no HCO+ tracer !!!"
     2794               stop
     2795            endif
    26852796            if (igcm_elec == 0) then
    26862797               write(*,*) "gcmtochim: Error; no e- tracer !!!"
    … …  
    27272838            rm(l,i_hplus)    = zycol(l, igcm_hplus)
    27282839            rm(l,i_hco2plus) = zycol(l, igcm_hco2plus)
     2840            rm(l,i_hcoplus)  = zycol(l, igcm_hcoplus)
    27292841            rm(l,i_elec)     = zycol(l, igcm_elec)
    27302842         end do
    … …  
    27552867                           i_co2plus, i_oplus, i_o2plus, i_noplus,    &
    27562868                           i_coplus, i_cplus, i_n2plus, i_nplus,      &
    2757                            i_hplus, i_hco2plus, i_elec, dens, c)
     2869                           i_hplus, i_hco2plus, i_hcoplus, i_elec,    &
     2870                           dens, c)
    27582871 
    27592872!*****************************************************************
    … …  
    27662879                            igcm_noplus, igcm_coplus, igcm_cplus,         &
    27672880                            igcm_n2plus, igcm_nplus, igcm_hplus,          &
    2768                             igcm_hco2plus, igcm_elec
     2881                            igcm_hco2plus, igcm_hcoplus, igcm_elec
    27692882
    27702883      implicit none
    … …  
    27862899                 i_co2plus, i_oplus, i_o2plus, i_noplus,         &
    27872900                 i_coplus, i_cplus, i_n2plus, i_nplus,           &
    2788                  i_hplus, i_hco2plus, i_elec
     2901                 i_hplus, i_hco2plus, i_hcoplus, i_elec
    27892902
    27902903      real :: dens(nlayer)     ! total number density (molecule.cm-3)
    … …  
    28392952            zycol(l, igcm_hplus)   = c(l,i_hplus)/dens(l)
    28402953            zycol(l, igcm_hco2plus)= c(l,i_hco2plus)/dens(l)
     2954            zycol(l, igcm_hcoplus) = c(l,i_hcoplus)/dens(l)
    28412955            zycol(l, igcm_elec)    = c(l,i_elec)/dens(l)
    28422956         end do
Note: See TracChangeset for help on using the changeset viewer.