!*********************************************************************** module chimiedata_h !*********************************************************************** ! ! subject: ! -------- ! ! data for photochemistry ! ! update 06/03/2021 set to module - add global variable (Yassin Jaziri) ! ! !*********************************************************************** implicit none character*30, save, allocatable :: chemnoms(:) ! name of chemical tracers !-------------------------------------------- ! dimensions of photolysis lookup table !-------------------------------------------- integer, save :: nd ! species integer, save :: nz ! altitude integer, save :: nozo ! ozone integer, save :: nsza ! solar zenith angle integer, save :: ntemp ! temperature integer, save :: ntau ! dust integer, save :: nch4 ! ch4 !-------------------------------------------- real, parameter :: kb = 1.3806e-23 ! photolysis real, save, allocatable :: jphot(:,:,:,:,:,:,:) real, save, allocatable :: ephot(:,:,:,:,:,:,:) real, save, allocatable :: colairtab(:) real, save, allocatable :: szatab(:) real, save, allocatable :: table_ozo(:) real, save, allocatable :: tautab(:) real, save, allocatable :: table_ch4(:) ! deposition real, save, allocatable :: SF_value(:) ! SF_mode=1 (vmr) SF_mode=2 (cm.s-1) real, save, allocatable :: prod_rate(:) ! production flux in molecules/m2/s real, save, allocatable :: surface_flux(:,:) ! Surface flux from deposition molecules/m2/s real, save, allocatable :: surface_flux2(:,:) ! Surface flux mean molecules/m2/s real, save, allocatable :: escape(:,:) ! escape rate in molecules/m2/s integer, save, allocatable :: SF_mode(:) ! SF_mode=1 (fixed mixing ratio) / 2 (fixed sedimentation velocity in cm/s and/or flux in molecules/m^3) !-------------------------------------------- real, save, allocatable :: albedo_snow_chim(:) ! albedo snow on chemistry wavelenght grid real, save, allocatable :: albedo_co2_ice_chim(:) ! albedo co2 ice on chemistry wavelenght grid !-------------------------------------------- ! For hard coding reactions !-------------------------------------------- integer, parameter :: nphot_hard_coding = 0 integer, parameter :: n4_hard_coding = 0 integer, parameter :: n3_hard_coding = 0 contains subroutine indice_HC(nb_phot,nb_reaction_4,nb_reaction_3) use types_asis implicit none integer, intent(inout) :: nb_phot integer, intent(inout) :: nb_reaction_4 integer, intent(inout) :: nb_reaction_3 !!!!!!!!!!! Hard coding reaction !!!!!!!!!!!!!!!!!!!!!!!!!!! !=========================================================== ! d022 : HO2 + NO + M -> HNO3 + M !=========================================================== !nb_reaction_4 = nb_reaction_4 + 1 !indice_4(nb_reaction_4) = z4spec(1.0, indexchim('ho2'), 1.0, indexchim('no'), 1.0, indexchim('hno3'), 0.0, 1) !=========================================================== ! e001 : CO + OH -> CO2 + H !=========================================================== !nb_reaction_4 = nb_reaction_4 + 1 !indice_4(nb_reaction_4) = z4spec(1.0, indexchim('co'), 1.0, indexchim('oh'), 1.0, indexchim('co2'), 1.0, indexchim('h')) !=========================================================== ! f028 : CH + CH4 -> C2H4 + H !=========================================================== !nb_reaction_4 = nb_reaction_4 + 1 !indice_4(nb_reaction_4) = z4spec(1.0, indexchim('ch'), 1.0, indexchim('ch4'), 1.0, indexchim('c2h4'), 1.0, indexchim('h')) !=========================================================== ! r001 : HNO3 + rain -> H2O !=========================================================== !nb_phot = nb_phot + 1 !indice_phot(nb_phot) = z3spec(1.0, indexchim('hno3'), 1.0, indexchim('h2o_vap'), 0.0, 1) !=========================================================== ! e001 : CO + OH -> CO2 + H !=========================================================== !nb_reaction_4 = nb_reaction_4 + 1 !indice_4(nb_reaction_4) = z4spec(1.0, indexchim('co'), 1.0, indexchim('oh'), 1.0, indexchim('co2'), 1.0, indexchim('h')) !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc ! photodissociation of NO : NO + hv -> N + O !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc !nb_phot = nb_phot + 1 !indice_phot(nb_phot) = z3spec(1.0, indexchim('no'), 1.0, indexchim('n'), 1.0, indexchim('o')) !!!!!!!!!!! END Hard coding reaction !!!!!!!!!!!!!!!!!!!!!!! end subroutine indice_HC subroutine reactionrates_HC(nlayer,nesp,dens,t,press,zmmean,sza,c,v_phot,v_4,v_3,nb_phot,nb_reaction_4,nb_reaction_3) use comcstfi_mod, only: g, avocado use types_asis implicit none integer, intent(in) :: nlayer ! number of atmospheric layers integer, intent(in) :: nesp ! number of chemical species integer, intent(inout) :: nb_phot integer, intent(inout) :: nb_reaction_4 integer, intent(inout) :: nb_reaction_3 real, intent(in), dimension(nlayer) :: dens ! total number density (molecule.cm-3) real, intent(in), dimension(nlayer) :: press ! pressure (hPa) real, intent(in), dimension(nlayer) :: t ! temperature (K) real, intent(in), dimension(nlayer) :: zmmean ! mean molar mass (g/mole) real, intent(in) :: sza ! solar zenith angle (deg) real (kind = 8), intent(in), dimension(nlayer,nesp) :: c ! species number density (molecule.cm-3) real (kind = 8), intent(inout), dimension(nlayer, nb_phot_max) :: v_phot real (kind = 8), intent(inout), dimension(nlayer,nb_reaction_4_max) :: v_4 real (kind = 8), intent(inout), dimension(nlayer,nb_reaction_3_max) :: v_3 ! local real, dimension(nlayer) :: colo3 ! ozone columns (molecule.cm-2) integer :: ilev real :: k1a0, k1b0, k1ainf, k1a, k1b, fc, fx, x, y real :: ak0, ak1, rate, xpo, deq, rate1, rate2, xpo1, xpo2 real :: rain_h2o, rain_rate !!!!!!!!!!! Hard coding reaction !!!!!!!!!!!!!!!!!!!!!!!!!!! !---------------------------------------------------------------------- ! nitrogen reactions !---------------------------------------------------------------------- !--- d022: ho2 + no + m -> hno3 + m ! Butkovskaya et al. 2007 !d022(:) = 6.4e-14*exp(1644/T(:)) !nb_reaction_4 = nb_reaction_4 + 1 !v_4(:,nb_reaction_4) = 3.3e-12*exp(270./t(:))*((530/T(:)+6.4*1e-4*press(:)/1.33322-1.73))*1e-2 !---------------------------------------------------------------------- ! carbon reactions !---------------------------------------------------------------------- !--- e001: oh + co -> co2 + h !nb_reaction_4 = nb_reaction_4 + 1 ! jpl 2003 !e001(:) = 1.5e-13*(1 + 0.6*press(:)/1013.) ! mccabe et al., grl, 28, 3135, 2001 !e001(:) = 1.57e-13 + 3.54e-33*dens(:) ! jpl 2006 !ak0 = 1.5e-13*(t(:)/300.)**(0.6) !ak1 = 2.1e-9*(t(:)/300.)**(6.1) !rate1 = ak0/(1. + ak0/(ak1/dens(:))) !xpo1 = 1./(1. + alog10(ak0/(ak1/dens(:)))**2) !ak0 = 5.9e-33*(t(:)/300.)**(-1.4) !ak1 = 1.1e-12*(t(:)/300.)**(1.3) !rate2 = (ak0*dens(:))/(1. + ak0*dens(:)/ak1) !xpo2 = 1./(1. + alog10((ak0*dens(:))/ak1)**2) !e001(:) = rate1*0.6**xpo1 + rate2*0.6**xpo2 ! jpl 2015 !do ilev = 1,nlayer !!oh + co -> h + co2 ! rate1 = 1.5e-13*(t(ilev)/300.)**(0.0) !!oh + co + m -> hoco ! ak0 = 5.9e-33*(t(ilev)/300.)**(-1.0) ! ak1 = 1.1e-12*(t(ilev)/300.)**(1.3) ! rate2 = (ak0*dens(ilev))/(1. + ak0*dens(ilev)/ak1) ! xpo2 = 1./(1. + alog10((ak0*dens(ilev))/ak1)**2) ! ! v_4(ilev,nb_reaction_4) = rate1 + rate2*0.6**xpo2 !end do ! joshi et al., 2006 !do ilev = 1,nlayer ! k1a0 = 1.34*2.5*dens(ilev) & ! *1/(1/(3.62e-26*t(ilev)**(-2.739)*exp(-20./t(ilev))) & ! + 1/(6.48e-33*t(ilev)**(0.14)*exp(-57./t(ilev)))) ! typo in paper corrected ! k1b0 = 1.17e-19*t(ilev)**(2.053)*exp(139./t(ilev)) & ! + 9.56e-12*t(ilev)**(-0.664)*exp(-167./t(ilev)) ! k1ainf = 1.52e-17*t(ilev)**(1.858)*exp(28.8/t(ilev)) & ! + 4.78e-8*t(ilev)**(-1.851)*exp(-318./t(ilev)) ! x = k1a0/(k1ainf - k1b0) ! y = k1b0/(k1ainf - k1b0) ! fc = 0.628*exp(-1223./t(ilev)) + (1. - 0.628)*exp(-39./t(ilev)) & ! + exp(-t(ilev)/255.) ! fx = fc**(1./(1. + (alog(x))**2)) ! typo in paper corrected ! k1a = k1a0*((1. + y)/(1. + x))*fx ! k1b = k1b0*(1./(1.+x))*fx ! ! v_4(ilev,nb_reaction_4) = k1a + k1b !end do !--- f028: ch + ch4 + m -> c2h4 + h + m ! Romani 1993 !nb_reaction_4 = nb_reaction_4 + 1 !v_4(:,nb_reaction_4) = min(2.5e-11*exp(200./t(:)),1.7e-10) !---------------------------------------------------------------------- ! washout r001 : HNO3 + rain -> H2O !---------------------------------------------------------------------- !nb_phot = nb_phot + 1 ! !rain_h2o = 100.e-6 !!rain_rate = 1.e-6 ! 10 days !rain_rate = 1.e-8 ! !do ilev = 1,nlayer ! if (c(ilev,indexchim('h2o_vap'))/dens(ilev) >= rain_h2o) then ! v_phot(ilev,nb_phot) = rain_rate ! else ! v_phot(ilev,nb_phot) = 0. ! end if !end do !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc ! photodissociation of NO !ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc !nb_phot = nb_phot + 1 ! !colo3(nlayer) = 0. !! ozone columns for other levels (molecule.cm-2) !do ilev = nlayer-1,1,-1 ! colo3(ilev) = colo3(ilev+1) + (c(ilev+1,indexchim('o3')) + c(ilev,indexchim('o3')))*0.5*avocado*1e-4*((press(ilev) - press(ilev+1))*100.)/(1.e-3*zmmean(ilev)*g*dens(ilev)) !end do !call jno(nlayer, c(nlayer:1:-1,indexchim('no')), c(nlayer:1:-1,indexchim('o2')), colo3(nlayer:1:-1), dens(nlayer:1:-1), press(nlayer:1:-1), sza, v_phot(nlayer:1:-1,nb_phot)) !!!!!!!!!!! END Hard coding reaction !!!!!!!!!!!!!!!!!!!!!!! end subroutine reactionrates_HC subroutine jno(nivbas, c_no, c_o2, o3t, hnm, pm, sza, tjno) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! ! ! parametrisation de la photodissociation de no ! ! d'apres minschwaner and siskind, a new calculation ! ! of nitric oxide photolysis in the stratosphere, ! ! mesosphere, and lower thermosphere, j. geophys. res., ! ! 98, 20401-20412, 1993 ! ! ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! implicit none ! input integer, intent(in) :: nivbas real (kind = 8), dimension(nivbas), intent(in) :: c_no, c_o2 real, dimension(nivbas), intent(in) :: hnm, o3t, pm real, intent(in) :: sza ! output real, dimension(nivbas), intent(out) :: tjno ! local real, parameter :: factor = 3.141592/180. real, dimension(nivbas) :: colo2, colno, colo3 real, dimension(nivbas) :: r_o2, r_no real, dimension(nivbas) :: p, t_o3_1, t_o3_2, t_o3_3 real, dimension(nivbas) :: bd1, bd2, bd3 real, dimension(6,3) :: sigma_o2, sigma_no_1, sigma_no_2 real, dimension(6,3) :: w_no_1, w_no_2 real, dimension(3) :: delta_lambda, i0, sigma_o3 real :: a, d, kq, n2 real :: sec, zcmt, dp integer :: iniv ! sections efficaces sigma_o2(:,1) = (/1.12e-23, 2.45e-23, 7.19e-23, & 3.04e-22, 1.75e-21, 1.11e-20/) sigma_o2(:,2) = (/1.35e-22, 2.99e-22, 7.33e-22, & 3.07e-21, 1.69e-20, 1.66e-19/) sigma_o2(:,3) = (/2.97e-22, 5.83e-22, 2.05e-21, & 8.19e-21, 4.80e-20, 2.66e-19/) sigma_no_1(:,1) = (/0.00e+00, 1.32e-18, 6.35e-19, & 7.09e-19, 2.18e-19, 4.67e-19/) sigma_no_1(:,2) = (/0.00e+00, 0.00e+00, 3.05e-21, & 5.76e-19, 2.29e-18, 2.21e-18/) sigma_no_1(:,3) = (/1.80e-18, 1.50e-18, 5.01e-19, & 7.20e-20, 6.72e-20, 1.49e-21/) sigma_no_2(:,1) = (/0.00e+00, 4.41e-17, 4.45e-17, & 4.50e-17, 2.94e-17, 4.35e-17/) sigma_no_2(:,2) = (/0.00e+00, 0.00e+00, 3.20e-21, & 5.71e-17, 9.09e-17, 6.00e-17/) sigma_no_2(:,3) = (/1.40e-16, 1.52e-16, 7.00e-17, & 2.83e-17, 2.73e-17, 6.57e-18/) ! facteurs de ponderation w_no_1(:,1) = (/0.00e+00, 5.12e-02, 1.36e-01, & 1.65e-01, 1.41e-01, 4.50e-02/) w_no_1(:,2) = (/0.00e+00, 0.00e+00, 1.93e-03, & 9.73e-02, 9.75e-02, 3.48e-02/) w_no_1(:,3) = (/4.50e-02, 1.80e-01, 2.25e-01, & 2.25e-01, 1.80e-01, 4.50e-02/) w_no_2(:,1) = (/0.00e+00, 5.68e-03, 1.52e-02, & 1.83e-02, 1.57e-02, 5.00e-03/) w_no_2(:,2) = (/0.00e+00, 0.00e+00, 2.14e-04, & 1.08e-02, 1.08e-02, 3.86e-03/) w_no_2(:,3) = (/5.00e-03, 2.00e-02, 2.50e-02, & 2.50e-02, 2.00e-02, 5.00e-03/) ! largeurs spectrales pour les trois bandes considerees (nm) delta_lambda = (/2.3, 1.5, 1.5/) ! flux pour les trois bandes considerees (s-1 cm-2 nm-1) i0 = (/3.98e+11, 2.21e+11, 2.30e+11/) ! sections efficaces de l'ozone pour les trois bandes ! considerees (cm2) d'apres wmo 1985 sigma_o3 = (/4.6e-19, 6.7e-19, 7.1e-19/) ! zcmt: gas constant/boltzmann constant*g zcmt = 287.0/(1.38e-23*9.81) ! d: taux de predissociation spontanee (s-1) d = 1.65e+09 ! a: taux d'emission spontanee (s-1) a = 5.1e+07 ! kq: quenching rate constant (cm3 s-1) kq = 1.5e-09 ! rapports de melange r_no(:) = c_no(:)/hnm(:) r_o2(:) = c_o2(:)/hnm(:) !============================================================= ! calcul des colonnes de o2 et no !============================================================= ! premier niveau du modele (mol/cm2) colo2(1) = 6.8e+05*c_o2(1) colno(1) = 6.8e+05*c_no(1) colo3(1) = o3t(1) ! cas general do iniv = 2,nivbas dp = (pm(iniv) - pm(iniv - 1))*100. dp = max(dp, 0.) colo2(iniv) = colo2(iniv - 1) & + zcmt*(r_o2(iniv-1) + r_o2(iniv))*0.5*dp*1.e-4 colno(iniv) = colno(iniv - 1) & + zcmt*(r_no(iniv-1) + r_no(iniv))*0.5*dp*1.e-4 colo3(iniv) = o3t(iniv) end do !============================================================= ! boucle sur les niveaux !============================================================= do iniv = 1,nivbas if (sza <= 89.0) then sec = 1./cos(sza*factor) colo2(iniv) = colo2(iniv)*sec colno(iniv) = colno(iniv)*sec colo3(iniv) = colo3(iniv)*sec ! facteurs de transmission de l'ozone t_o3_1(iniv) = exp(-sigma_o3(1)*colo3(iniv)) t_o3_2(iniv) = exp(-sigma_o3(2)*colo3(iniv)) t_o3_3(iniv) = exp(-sigma_o3(3)*colo3(iniv)) ! calcul de la probabilite de predissociation n2 = hnm(iniv)*0.78 p(iniv) = d/(a + d + kq*n2) end if end do ! calcul proprement dit, pour les 3 bandes do iniv = 1,nivbas if (sza <= 89.0) then bd1(iniv) = delta_lambda(1)*i0(1)*t_o3_1(iniv)*p(iniv)*( & exp(-sigma_o2(1,1)*colo2(iniv)) & *(w_no_1(1,1)*sigma_no_1(1,1)*exp(-sigma_no_1(1,1)*colno(iniv)) & + w_no_2(1,1)*sigma_no_2(1,1)*exp(-sigma_no_2(1,1)*colno(iniv))) & + exp(-sigma_o2(2,1)*colo2(iniv)) & *(w_no_1(2,1)*sigma_no_1(2,1)*exp(-sigma_no_1(2,1)*colno(iniv)) & + w_no_2(2,1)*sigma_no_2(2,1)*exp(-sigma_no_2(2,1)*colno(iniv))) & + exp(-sigma_o2(3,1)*colo2(iniv)) & *(w_no_1(3,1)*sigma_no_1(3,1)*exp(-sigma_no_1(3,1)*colno(iniv)) & + w_no_2(3,1)*sigma_no_2(3,1)*exp(-sigma_no_2(3,1)*colno(iniv))) & + exp(-sigma_o2(4,1)*colo2(iniv)) & *(w_no_1(4,1)*sigma_no_1(4,1)*exp(-sigma_no_1(4,1)*colno(iniv)) & + w_no_2(4,1)*sigma_no_2(4,1)*exp(-sigma_no_2(4,1)*colno(iniv))) & + exp(-sigma_o2(5,1)*colo2(iniv)) & *(w_no_1(5,1)*sigma_no_1(5,1)*exp(-sigma_no_1(5,1)*colno(iniv)) & + w_no_2(5,1)*sigma_no_2(5,1)*exp(-sigma_no_2(5,1)*colno(iniv))) & + exp(-sigma_o2(6,1)*colo2(iniv)) & *(w_no_1(6,1)*sigma_no_1(6,1)*exp(-sigma_no_1(6,1)*colno(iniv)) & + w_no_2(6,1)*sigma_no_2(6,1)*exp(-sigma_no_2(6,1)*colno(iniv))) & ) bd2(iniv) = delta_lambda(2)*i0(2)*t_o3_2(iniv)*p(iniv)*( & exp(-sigma_o2(1,2)*colo2(iniv)) & *(w_no_1(1,2)*sigma_no_1(1,2)*exp(-sigma_no_1(1,2)*colno(iniv)) & + w_no_2(1,2)*sigma_no_2(1,2)*exp(-sigma_no_2(1,2)*colno(iniv))) & + exp(-sigma_o2(2,2)*colo2(iniv)) & *(w_no_1(2,2)*sigma_no_1(2,2)*exp(-sigma_no_1(2,2)*colno(iniv)) & + w_no_2(2,2)*sigma_no_2(2,2)*exp(-sigma_no_2(2,2)*colno(iniv))) & + exp(-sigma_o2(3,2)*colo2(iniv)) & *(w_no_1(3,2)*sigma_no_1(3,2)*exp(-sigma_no_1(3,2)*colno(iniv)) & + w_no_2(3,2)*sigma_no_2(3,2)*exp(-sigma_no_2(3,2)*colno(iniv))) & + exp(-sigma_o2(4,2)*colo2(iniv)) & *(w_no_1(4,2)*sigma_no_1(4,2)*exp(-sigma_no_1(4,2)*colno(iniv)) & + w_no_2(4,2)*sigma_no_2(4,2)*exp(-sigma_no_2(4,2)*colno(iniv))) & + exp(-sigma_o2(5,2)*colo2(iniv)) & *(w_no_1(5,2)*sigma_no_1(5,2)*exp(-sigma_no_1(5,2)*colno(iniv)) & + w_no_2(5,2)*sigma_no_2(5,2)*exp(-sigma_no_2(5,2)*colno(iniv))) & + exp(-sigma_o2(6,2)*colo2(iniv)) & *(w_no_1(6,2)*sigma_no_1(6,2)*exp(-sigma_no_1(6,2)*colno(iniv)) & + w_no_2(6,2)*sigma_no_2(6,2)*exp(-sigma_no_2(6,2)*colno(iniv))) & ) bd3(iniv) = delta_lambda(3)*i0(3)*t_o3_3(iniv)*p(iniv)*( & exp(-sigma_o2(1,3)*colo2(iniv)) & *(w_no_1(1,3)*sigma_no_1(1,3)*exp(-sigma_no_1(1,3)*colno(iniv)) & + w_no_2(1,3)*sigma_no_2(1,3)*exp(-sigma_no_2(1,3)*colno(iniv))) & + exp(-sigma_o2(2,3)*colo2(iniv)) & *(w_no_1(2,3)*sigma_no_1(2,3)*exp(-sigma_no_1(2,3)*colno(iniv)) & + w_no_2(2,3)*sigma_no_2(2,3)*exp(-sigma_no_2(2,3)*colno(iniv))) & + exp(-sigma_o2(3,3)*colo2(iniv)) & *(w_no_1(3,3)*sigma_no_1(3,3)*exp(-sigma_no_1(3,3)*colno(iniv)) & + w_no_2(3,3)*sigma_no_2(3,3)*exp(-sigma_no_2(3,3)*colno(iniv))) & + exp(-sigma_o2(4,3)*colo2(iniv)) & *(w_no_1(4,3)*sigma_no_1(4,3)*exp(-sigma_no_1(4,3)*colno(iniv)) & + w_no_2(4,3)*sigma_no_2(4,3)*exp(-sigma_no_2(4,3)*colno(iniv))) & + exp(-sigma_o2(5,3)*colo2(iniv)) & *(w_no_1(5,3)*sigma_no_1(5,3)*exp(-sigma_no_1(5,3)*colno(iniv)) & + w_no_2(5,3)*sigma_no_2(5,3)*exp(-sigma_no_2(5,3)*colno(iniv))) & + exp(-sigma_o2(6,3)*colo2(iniv)) & *(w_no_1(6,3)*sigma_no_1(6,3)*exp(-sigma_no_1(6,3)*colno(iniv)) & + w_no_2(6,3)*sigma_no_2(6,3)*exp(-sigma_no_2(6,3)*colno(iniv))) & ) tjno(iniv) = bd1(iniv) + bd2(iniv) + bd3(iniv) else ! night tjno(iniv) = 0. end if end do end subroutine jno function indexchim(traceurname) use tracer_h, only: noms, nqtot, is_chim implicit none !======================================================================= ! Get index of tracer in chemistry with its name using traceur tab ! ! version: April 2019 - Yassin Jaziri (update September 2020) !======================================================================= ! input/output integer :: indexchim ! index of tracer in chemistry character(len=*) :: traceurname ! name of traceur to find ! local variables integer :: iq, iesp iesp = 0 indexchim = 0 do iq = 1,nqtot if (is_chim(iq)==1) then iesp = iesp + 1 if (trim(traceurname)==trim(noms(iq))) then indexchim = iesp exit end if end if end do if (indexchim==0) then print*, 'Error indexchim: wrong traceur name ', trim(traceurname) print*, 'Check if the specie ', trim(traceurname),' is in traceur.def' print*, 'Check if the specie ', trim(traceurname),' is a chemical species' print*, '(option is_chim = 1)' call abort end if return end function indexchim end module chimiedata_h