************************************************************************ * MASTER MECHANISM - ROUTINE NAME : no3rad * * * * * * PURPOSE : * * perform reaction of NO3 with VOC, i.e. : * * --> abstraction of H from -CHO groups (NO3+aldehyde) (PART I) * * --> addition of NO3 to C=C double bonds (PART II) * * * * INPUT: * * - rdct : name[a6]+formula[a120] of the species for which * * reaction with OH is considered * * - bond(i,j) : carbon-carbon bond matrix of chem * * - group(i) : groups at position (carbon) i * * - nring : min number of rings present * * - temp : reference temperature (more work on this needed) * * - cut_off : ratio below which a pathway is not considered * * - dbrch : NOT USED - MORE WORK ON THIS NEEDED * * - level : number of level (stable + radicals) that were * * necessary to produce the parent of rdct * * - stabl : number of stable level (no radical) that were * * necessary to produce the parent of rdct * * - nfn : total nb. of species having a fixed name * * - namfn(i) : table of the fixed name (6 character) * * - chemfn(i) : formula corresponding the ith species having a * * fixed name * * * * Data for the reactions with NO3 (see no3_rate.dat) * * arrhenius parameter (see no3_rate.dat) * * - nno3dat : total number of known NO3 rate constant * * - no3chem(i) : chemical formula of ith species having a known NO3 * * rate constant * * - no3rate(i,3) : arrhenius parameter of the known NO3 rate constant * * corresponding to the ith species * * known reaction (see no3_prod.dat) * * - nkwno3 : number of VOC having known products for the react * * with NO3 * * - kwno3rct(i) : chemical formula of the VOC having a known * * chemistry for the reaction with NO3 * * - kwno3pd(i,j) : chemical formula of the jth product for the * * reaction of the ith VOC with NO3 * * - kwno3copd(i,j,k): name (6 character) of the k thcoproduct of the * * jth product for the reaction of the ith VOC withNO3 * * - kwno3yld(i,j): yield of the jth product for the reaction of the * * ith VOC with NO3 * * * * * * INPUT/OUTPUT * * - dict(j) : dictionnary line (name + formula + functional * * group info) of species number j * * - namlst(j) : name (lco=6 characters) of the species already * * used at position number j * * - nhldvoc : number of (stable) VOC in the stack * * - holdvoc(i) : list of the VOC in the stack * * - nhldrad : number of radical in the stack * * - holdrad(i) : list of the radicals in the stack * * * ********************************************************************** SUBROUTINE no3rad(rdct,bond,group,nring,brch, & dbrch,dict,namlst, & cut_off, & nhldvoc,holdvoc,nhldrad,holdrad,level,stabl, & nno3dat,no3chem,no3rate,nfn,namfn,chemfn, & nkwno3,kwno3rct,kwno3pd,kwno3copd,kwno3yld) IMPLICIT NONE INCLUDE 'general.h' INCLUDE 'organic.h' INCLUDE 'common.h' * input: CHARACTER(LEN=lcf) rdct CHARACTER(LEN=lgr) group(mca) INTEGER bond(mca,mca) REAL brch REAL cut_off INTEGER nkwno3 CHARACTER(LEN=lfo) no3chem(mrd) REAL no3rate(mrd,3) INTEGER nno3dat CHARACTER(LEN=lfo) kwno3rct(mkr), kwno3pd(mkr,mnr) CHARACTER(LEN=lco) kwno3copd(mkr,mnr,mcp) REAL kwno3yld(mkr,mnr) * input/output CHARACTER(LEN=ldi) dict(mni) REAL dbrch(mni) CHARACTER(LEN=lco) namlst(mni) INTEGER level INTEGER stabl CHARACTER(LEN=lst) holdvoc(mlv) INTEGER nhldvoc CHARACTER(LEN=lst) holdrad(mra) INTEGER nhldrad INTEGER nfn CHARACTER(LEN=lco) namfn(mfn) CHARACTER(LEN=lfo) chemfn(mfn) * internal: allow for up to "mnr" reaction channels CHARACTER(LEN=lfo) chem CHARACTER(LEN=lfo) pchem(mnr,2),tempkc,pchema(mnr,2) CHARACTER(LEN=lgr) tgroup(mca), pold, pnew, tempgr CHARACTER(LEN=lco) coprod(mnr,mca), tprod(mca) INTEGER tbond(mca,mca), flag(mnr), nring INTEGER nr,nch,ich,i,j,k,np,ncd,nc,nca REAL brtio, totr, sum, rate, rno3(mnr) REAL rtot REAL ratio(mnr) REAL garrhc(3),tarrhc(mnr,3),rarrhc(3) REAL rk298 CHARACTER*1 a1,a2,a3,a4 CHARACTER(LEN=lco) r(3), p(mnp) REAL s(mnp), ar1,ar2,ar3,f298,fratio,fac,kstruct INTEGER known_species INTEGER idreac, nlabel REAL xlabel,folow(3),fotroe(4) CHARACTER(LEN=lco) copchem REAL rdtcopchem,A,E,nt INTEGER cdtable(4),tcdtable(4), cdsub(4),conjug, ncdcase INTEGER cdcarbo(4,2),cdeth(4,2) INTEGER cpcdtable(4) INTEGER rjg(mri,2) ! ring-join group pairs CHARACTER(LEN=lfo) rdckprod(mca) CHARACTER(LEN=lco) rdcktprod(mca,mca) INTEGER nip,flag_del(mnr) REAL sc(mca),sc_del(mnr,mca) CHARACTER(LEN=lco) coprod_del(mnr,mca) CHARACTER(LEN=lfo) pchem_del(mnr,2) CHARACTER(LEN=lsb) :: progname='*no3rad* ' CHARACTER(LEN=ler) :: mesg * write info for finding bugs IF (wtflag.NE.0) WRITE(*,*) '*no3_rad*' * initialize: copchem=' ' rdtcopchem=0. DO i=1,mca tgroup(i) = group(i) ENDDO DO i=1,mca DO j=1,mca tbond(i,j) = bond(i,j) ENDDO ENDDO rate = 0. totr = 0. nr = 0 rk298 = 0. DO i=1,mnr pchem(i,1) = ' ' pchem(i,2) = ' ' pchem_del(i,1) = ' ' pchem_del(i,2) = ' ' pchema(i,1) = ' ' pchema(i,2) = ' ' flag(i) = 0 flag_del(i) = 0 rno3(i) = 0. ENDDO DO i=1,mnr DO j=1,mca coprod(i,j) = ' ' coprod_del(i,j) = ' ' ENDDO DO j=1,3 tarrhc(i,j)=0. ENDDO ENDDO DO j=1,3 rarrhc(j)=0. ENDDO DO i=1,4 tcdtable(i)=0 cdtable(i)=0 ENDDO chem=' ' i = lco + INDEX(rdct(lco+1:lcf),' ') chem=rdct(lco+1:i) known_species=0 nca=0 DO i=1,mca IF (tgroup(i).NE.' ') nca=nca+1 ENDDO * ----------------------------------------------------------- * check if the species is known in the rate constant data base * ----------------------------------------------------------- k = lco + index(rdct(lco+1:lcf),' ') DO i=1,nno3dat IF (rdct(lco+1:k).EQ.no3chem(i)) THEN rk298=no3rate(i,1)*298**no3rate(i,2)*exp(-no3rate(i,3)/298.) DO j=1,3 rarrhc(j)=no3rate(i,j) ENDDO ENDIF ENDDO * ----------------------------------------------------------- * check if products are known in no3_prod.dat * WARNING: for those reaction, reaction < 5% are not automatically * eliminated * ----------------------------------------------------------- ctemp IF (rk298.NE.0.) THEN ctemp known_species=0 ctemp DO i=1,nno3prod ctemp IF (rdct(lco+1:o).EQ.reactno3(i)) THEN ctemp nr = nr + 1 ctemp IF (nr.GT.mnr) THEN ctemp WRITE(6,'(a)') '--error--' ctemp WRITE(6,'(a)') 'from MASTER MECHANISM ROUTINE : ho_rad2' ctemp WRITE(6,'(a)') 'too many reactions created for species' ctemp WRITE(6,'(a)') rdct(lco+1:lcf) ctemp WRITE(99,*) 'no3rad',rdct(lco+1:lcf) !STOP ctemp ENDIF ctemp flag(nr) = 1 ctemp known_species = 1 ctemp pchema(nr,1) = prodno3(i,1) ctemp pchema(nr,2) = prodno3(i,2) ctemp ratio(nr) = ratno3(i) * check the products ctemp IF (index(pchema(nr,1),'.').ne.0) THEN ctemp CALL radchk(pchema(nr,1),pchem(nr,1),tprod) ctemp DO j=1,mca ctemp coprod(nr,j) = tprod(j) ctemp ENDDO ctemp ELSE ctemp pchem(nr,1)= pchema(nr,1) ctemp ENDIF ctemp CALL stdchm(pchem(nr,2)) ctemp pchem(nr,2)= pchema(nr,2) * add third product if necessary ctemp IF (coprodno3(i).ne.' ') THEN ctemp DO j=1,mca ctemp IF (coprod(nr,j)(1:2).EQ.' ') THEN ctemp coprod(nr,j) = coprodno3(i) ctemp GOTO 450 ctemp ENDIF ctemp ENDDO ctemp ENDIF ctemp450 CONTINUE ctemp ENDIF ctemp ENDDO ctemp IF (known_species.ne.0) GOTO 800 ctemp ENDIF ************************************************************************ **** H ABSTRACTION *** ************************************************************************ * --------------------------------- * H-ABSTRACTION FROM CH * --------------------------------- * find all C-H bonds - loop over each C DO 110 i=1,nca IF (INDEX(group(i),'CH').EQ.0) GO TO 110 nr = nr + 1 IF (nr.GT.mnr) THEN WRITE(6,'(a)') '--error--' WRITE(6,'(a)') 'from MASTER MECHANISM ROUTINE : ho_rad' WRITE(6,'(a)') 'too many reactions created for species' WRITE(6,'(a)') rdct(lco+1:lcf) WRITE(99,*) 'ho_rad2',chem !STOP ENDIF flag(nr) = 1 * compute associated rate constant: k=A*T**2*exp(-E/T) DO j=1,3 garrhc(j)=0. ENDDO CALL rabsno3(tbond,tgroup,i,garrhc,nring) DO j=1,3 rno3(nr)=garrhc(1) ENDDO * remove one H atom IF (group(i)(1:3).EQ.methyl) THEN pold = methyl pnew = primary ELSE IF (group(i)(1:3).EQ.primary) THEN pold = primary pnew = secondary ELSE IF (group(i)(1:3).EQ.aldehyde) THEN pold = aldehyde pnew = carbonyl ELSE IF (group(i)(1:2).EQ.secondary) THEN pold = secondary pnew = 'C' ENDIF CALL swap(group(i),pold,tgroup(i),pnew) * add radical dot at end of group: nc = INDEX(tgroup(i),' ') tgroup(i)(nc:nc) = '.' * re-build molecule, including ring if present: CALL rebond(tbond,tgroup,tempkc,nring) * check radical and find co-products: CALL radchk(tempkc,rdckprod,rdcktprod,nip,sc) pchem(nr,1) = rdckprod(1) sc_del(nr,1) = sc(1) IF (nip.EQ.2) THEN flag_del(nr) = 1 pchem_del(nr,1) = rdckprod(2) sc_del(nr,2) = sc(2) DO j=1,mca coprod_del(nr,j) = rdcktprod(2,j) ENDDO ENDIF DO j=1,mca coprod(nr,j) = rdcktprod(1,j) ENDDO pchem(nr,2) = 'HNO3' * find standard name: CALL stdchm(pchem(nr,1)) * reset group: tgroup(i) = group(i) * add group rate constant to the sum totr = totr + rno3(nr) 110 CONTINUE *------------------------------------------------ * H ABSTRACTION FROM (OH) - BUT NOT CO(OH): *------------------------------------------------ * find all hydroxy group but not carboxylic acid DO 120 i=1,mca IF (INDEX(group(i),hydroxy).EQ.0) GO TO 120 IF (INDEX(group(i),carboxylic_acid).NE.0) GO TO 120 nr = nr + 1 IF (nr.GT.mnr) THEN mesg = 'too many reactions created for the species' CALL errexit(progname,mesg,chem) ENDIF flag(nr) = 1 * assign rate constant by Kerdouci et al., 2010 rno3(nr)=2.00E-17 * replace hydroxy group by alkoxy group: pold = hydroxy pnew = alkoxy CALL swap(group(i),pold,tgroup(i),pnew) * rebuild, check, and find coproducts: CALL rebond(tbond,tgroup,tempkc,nring) CALL radchk(tempkc,rdckprod,rdcktprod,nip,sc) pchem(nr,1) = rdckprod(1) IF (nip.NE.1) STOP '2 produits no3rad2.f 2' DO j=1,mca coprod(nr,j) = rdcktprod(1,j) ENDDO pchem(nr,2) = 'HNO3' * rename CALL stdchm(pchem(nr,1)) * reset group: tgroup(i) = group(i) 120 CONTINUE * ---------------------------- * abstraction of H from CHO * ---------------------------- * UPDATE RIC : aldehyde H abstraction is now treated before * rate constant for simple aldehyde where measured by Andresen * (2000, Eurotrac symposium) for 14 aldehydes. Except for * CH3CHO (k=~2.5E-15) and C2H5CHO (k=~6E-15), rate constant * where found to be greater than ~1.2E-14 (from 1.17 for butanal * to 4.61 for 2-ethyl butanal). In this preliminary version, * rate constant was set equal to the rate constant for butanal * (since rate constant for CH3CHO and C2H5CHO are known and * estimation are not required). * Based on comparison with the OH reaction rate, it can be expected * that the reactivity is significantly lowered for the -CO-CHO * structures. Based on the SAPRC99 chemical scheme (see discussion * in Carter, 1999), rate constant for glyoxal and methyl glyoxal are * calculated with a correlation from Atkinson 93 * ln kNO3 = 6.498 + 1.611 ln kOH * In SAPRC99, kNO3 are then expressed with temperature dependance of * acetaldehyde (i.e. (kCH3CHO*0.18)*2 per -CHO group, and for methyl * glyoxal (and similar structure) : kCH3CHO*0.89) * In this version temperature dependancy isn't taken into account * so values are 4.86E-16 per -CHO group and 2.40E-15 for methyl glyoxal * Rate constant for -C=C-CHO structure are assumed to be 2.5E-15 * estimated with k acetaldehyde and with the hypothesis of a quasi-total * abstraction for acrolein which rate constant is 2.5E-15. c* find all aldehyde-NO3 reaction channel: c DO i=1,mca c IF ( (group(i)(1:3).EQ.aldehyde) .AND. c & (group(i)(4:).EQ.' ' ) ) THEN c nr = nr + 1 c IF (nr.GT.mnr) THEN c WRITE(6,'(a)') '--error--' c WRITE(6,'(a)') 'from MASTER MECHANISM ROUTINE : no3rad' c WRITE(6,'(a)') 'too many reactions created for species' c WRITE(6,'(a)') rdct(lco+1:lcf) c WRITE(99,*) 'no3rad',rdct(lco+1:lcf) !STOP c ENDIF c flag(nr) = 1 c rno3(nr) = 1.2E-14 * overwrite reaction rate if group in beta is a carbonyl group and * stop if it is a Cd c DO j=1,mca c IF (tbond(i,j).NE.0) THEN c IF (group(j)(1:3).EQ.aldehyde) THEN c rno3(nr) = 4.86E-16 c ENDIF c IF (group(j)(1:2).EQ.carbonyl) THEN c rno3(nr) = 2.40E-15 c ENDIF c IF (INDEX(group(j),'Cd').NE.0) THEN c rno3(nr) = 2.5E-15 c ENDIF c ENDIF c ENDDO * add group rate constant to the sum c totr = totr + rno3(nr) * build product: c tgroup(i) = acyl c CALL rebond(tbond,tgroup,tempkc,nring) * check radical and find co-products: c CALL radchk(tempkc,pchem(nr,1),tprod) c CALL radchk(tempkc,rdckprod,rdcktprod,nip,sc) c IF (nip.EQ.1) pchem(nr,1) = rdckprod(1) c IF (nip.NE.1) STOP 'no3rad.f' c DO j=1,mca c coprod(nr,j) = rdcktprod(1,j) c ENDDO c pchem(nr,2) = 'HNO3' * find standard name: c CALL stdchm(pchem(nr,1)) * reset group: c tgroup(i) = group(i) c ENDIF c ENDDO * ------------------------------------------------------------ * addition of NO3 to Cd=Cd: send to the reaction subroutine * ------------------------------------------------------------ * * check if the molecule contains C=C bond IF (INDEX(rdct(lco+1:lcf),'Cd').EQ.0) GOTO 600 * IF RINGS EXIST remove ring-join characters from groups IF (nring.gt.0) THEN CALL rjgrm(nring,group,rjg) ENDIF * The parameterization used to assign a rate constant for NO3 addition * to C=C bonds depends on whether the C=C bond is conjugated or not * with a C=O bond (i.e. structure of -C=C-C=O type). The subroutine * cdcase return the "case" the species belong. Five cases are * considered : * * CASE 1 : regular Cd molecule : only >C=C< and >C=C-C=C< bonds in * the molecule (i.e, without conjugated C=C-C=O) * CASE 2 : are for structure containing the >C=C-C=O structure * but no C=C-C=C * CASE 3 : are for the -CO-C=C-C=C-C=O structure only (i.e. containing * carbonyl at both side of the conjugated C=C-C=C) * CASE 4 : Two double bonds non conjugated (i.e. C=C-C-C=C) but one * containing at least one C=C-C=O * CASE 5 : are for the -CO-C=C-C=C< structure (i.e. containing * carbonyl at only one side of the conjugated C=C-C=C ! Ludo: ! CASE 7 : reactions of ozone with vinyl ether (-O-C=C) and dihydrofurans CALL cdcase(chem,bond,group, & ncd,ncdcase,conjug, & cdtable,tcdtable,cdsub,cdcarbo,cdeth) * case 1 : C=C bonds conjugated or not, without any C=C-C=O structure * ------- IF (ncdcase.eq.1) THEN CALL no3add_c1(chem,bond,group,ncd,conjug,cdtable,cdsub, & nr,flag,rno3,pchem,coprod,flag_del, & pchem_del,coprod_del,sc_del) * case 2 : C=C-C=O bonds, only 1 C=C (i.e no conjugated Cds) * ------- ELSE IF (ncdcase.eq.2) THEN CALL no3add_c2(chem,bond,group, & ncd,conjug,cdtable,cdsub,cdcarbo, & nr,flag,rno3,pchem,coprod) * case 3 : -CO-C=C-C=C-C=O- structure * ------- ELSE IF (ncdcase.eq.3) THEN CALL no3add_c3(chem,bond,group, & ncd,conjug,cdtable,cdsub,cdcarbo, & nr,flag,rno3,pchem,coprod) * case 4 : -C=C-C-C=C-C=O- * ------- ELSE IF (ncdcase.eq.4) THEN * perform reaction for each C=C, ignoring the other C=C bond DO i=1,3,2 * check everything is ok IF (cdtable(i).eq.0) THEN mesg = 'case 4 encountered but Cd not found for' CALL errexit(progname,mesg,rdct(lco+1:lcf)) ENDIF * make a copy of cdtable and "kill" the Cd bond not considered DO j=1,4 cpcdtable(j)=0 ENDDO cpcdtable(i)=cdtable(i) cpcdtable(i+1)=cdtable(i+1) * send to reaction subroutine IF (cdcarbo(i,1).ne.0) THEN CALL no3add_c2(chem,bond,group, & ncd,conjug,cpcdtable,cdsub,cdcarbo, & nr,flag,rno3,pchem,coprod) ELSE IF (cdcarbo(i+1,1).ne.0) THEN CALL no3add_c2(chem,bond,group, & ncd,conjug,cpcdtable,cdsub,cdcarbo, & nr,flag,rno3,pchem,coprod) ELSE CALL no3add_c1(chem,bond,group,ncd,conjug,cpcdtable,cdsub, & nr,flag,rno3,pchem,coprod,flag_del, & pchem_del,coprod_del,sc_del) ENDIF ENDDO * case 5 : -CO-C=C-C=C< * ------- ELSE IF (ncdcase.eq.5) THEN CALL no3add_c5(chem,bond,group, & ncd,conjug,cdtable,cdsub,cdcarbo, & nr,flag,rno3,pchem,coprod,flag_del,pchem_del, & coprod_del,sc_del) ! case 7 : -C=C-O- ! ---------------- ELSE IF (ncdcase.EQ.7) THEN CALL no3add_c7(chem,bond,group,ncd,conjug,cdtable,cdsub, & cdeth,nr,flag,rno3,pchem,coprod) ENDIF * end of NO3 addition to C=C 600 CONTINUE * ------------------------------------------------------------ * end reactions * ------------------------------------------------------------ * compute total rate constant totr=0. DO i=1,nr totr=totr+rno3(i) ENDDO WRITE(43,'(a60,e9.3)') chem(1:60),totr * if total rate constant less than 1.0E-15, then no reaction: c IF (totr.LT.1.0E-15) THEN c WRITE(99,'(a)') '--warning--' c WRITE(99,'(a)') 'NO3 rate constant estimated less than 1E-15' c WRITE(99,'(a)') 'reaction with NO3 is not considered for:' c WRITE(99,'(a)') rdct(lco+1:lcf) c ENDIF * flag down for duplicate products: IF (nr.gt.1) THEN DO i=1,nr-1 DO 210 j=i+1,nr IF ( (pchem(i,1).EQ.pchem(j,1)) .AND. & (pchem(i,2).EQ.pchem(j,2)) ) THEN DO k=1,mca IF (coprod(i,k).NE.coprod(j,k)) GO TO 210 ENDDO flag(j) = 0 rno3(i) = rno3(i) + rno3(j) ENDIF 210 CONTINUE ENDDO ENDIF * active channel count: 5% cut-off. DO i=1,nr IF (rno3(i).LT.cut_off*totr) flag(i)=0 ENDDO * set branching ratio for each channel nch = 0 sum = 0. DO i=1,nr IF (flag(i).EQ.1) THEN nch=nch+1 sum = sum + rno3(i) ENDIF ENDDO DO i=1,mnr IF (flag(i).EQ.1) THEN ratio(i)=rno3(i)/sum ENDIF ENDDO 800 CONTINUE * set the reaction rate of each channel - prelimary version : if the * rate constant is assigned then multiply the first arrhenius coef. by * the branching ratio else rescale tarrhc(i,1) IF (rk298.NE.0.) THEN DO i=1,nr IF (flag(i).EQ.1) THEN tarrhc(i,1)=ratio(i)*rarrhc(1) tarrhc(i,2)=rarrhc(2) tarrhc(i,3)=rarrhc(3) ENDIF ENDDO ! WRITE(45,'(e9.3,3x,e9.3)') sum,rk298 ! WRITE(46,'(a60,3xe9.3,3x,e9.3)') chem(1:60),sum,rk298 IF (wtflag.NE.0) WRITE(76,*) rk298,sum,' ',chem ELSE DO i=1,mnr IF (flag(i).EQ.1) THEN tarrhc(i,1)=totr*ratio(i) tarrhc(i,2)=0. tarrhc(i,3)=0. ENDIF ENDDO ENDIF * -------------------- * output * -------------------- nt =0. E = 0. ich = 10 A = 0. * write information required for operator IF (wtopeflag.EQ.1) THEN !#Write information required for operator DO 11 i=1,mnr IF (flag(i).NE.0) THEN A = A + tarrhc(i,1) IF ((tarrhc(i,2).NE.nt).OR.(tarrhc(i,3).NE.E)) THEN A = sum nt = 0. E = 0. GOTO 11 ELSE nt=tarrhc(i,2) E=tarrhc(i,3) ENDIF ENDIF 11 CONTINUE WRITE(10,'(A16,A1,A6,A4,1X,E10.3,1X,f4.1,1X,f7.0)') & '****INIT GNO3 + ','G',rdct(1:lco),'****',A,nt,E ENDIF DO 300 i=1,mnr IF (flag(i).EQ.0) GO TO 300 * initialize reaction CALL rxinit3(a1,a2,a3,a4,r,s,p,ar1,ar2,ar3, & f298,fratio,idreac,nlabel,xlabel,folow,fotroe) ich = ich + 1 IF (nch.GT.1) a4 = alfa(ich:ich) a1 = rdct(1:1) ar1=tarrhc(i,1) ar2=tarrhc(i,2) ar3=tarrhc(i,3) IF (rk298.NE.0.) THEN f298=rk298 ELSE f298=totr ENDIF fratio=ratio(i) r(1) = rdct(1:lco) r(2) = 'NO3 ' s(1) = 1. * first main product (organic product) np=1 IF (known_species.eq.0) THEN brtio = brch * rno3(i)/totr ELSE brtio = brch * ratio(nr) ENDIF * CALL bratio(pchem(i,1),brtio,p(1),copchem,rdtcopchem, CALL bratio(pchem(i,1),brtio,p(1), & dbrch,dict,namlst, & nhldvoc,holdvoc,nhldrad,holdrad,level,stabl, & nfn,namfn,chemfn) IF (wtopeflag.EQ.1) & WRITE(10,'(f5.3,2X,A1,A6)') ratio(i), 'G',p(1) * second product (if delocalisation allowed) np = 1 IF (flag_del(i).NE.0) THEN np = np + 1 s(1) = sc_del(i,1) s(np) = sc_del(i,2) CALL bratio(pchem_del(i,1),brtio,p(np), & dbrch,dict,namlst, & nhldvoc,holdvoc,nhldrad,holdrad,level,stabl, & nfn,namfn,chemfn) ENDIF * second main product (if any, must be HNO3 here) IF (pchem(i,2)(1:1).NE.' ') THEN np = np + 1 s(2) = 1. IF (wtopeflag.EQ.1) & WRITE(10,'(f5.3,2X,A1,A6)')s(np)*ratio(i),'G',pchem(i,2) IF (pchem(i,2)(1:5).NE.'HNO3 ') THEN mesg = 'pchem for species 2 is not HNO3' CALL errexit(progname,mesg,rdct(lco+1:lcf)) ELSE p(2)='HNO3 ' ENDIF ENDIF * co-product (linked to the main organic product) DO j=1,mca IF (coprod(i,j)(1:1).NE.' ') THEN np = np + 1 IF (np.GT.mnp) then mesg = 'np > mnp (too many products in the reaction : 1)' CALL errexit(progname,mesg,rdct(lco+1:lcf)) ENDIF s(np) = 1. p(np) = coprod(i,j) IF (wtopeflag.EQ.1) & WRITE(10,'(f5.3,2X,A1,A6)')s(np)*ratio(i),'G',p(np) ENDIF ENDDO DO j=1,mca IF ((coprod_del(i,j)(1:1).NE.' ').AND. & (flag_del(i).NE.0)) THEN np = np + 1 IF (np.GT.mnp) then mesg = 'np > mnp (too many products in the reaction : 2)' CALL errexit(progname,mesg,rdct(lco+1:lcf)) ENDIF s(np) = 1. p(np) = coprod_del(i,j) IF (wtopeflag.EQ.1) write(10,'(20X,A1,A6)') 'G',p(np) ENDIF ENDDO IF (rdtcopchem.GT.0.) THEN np = np + 1 IF (np.GT.mnp) then mesg = 'np > mnp (too many products in the reaction : 3)' CALL errexit(progname,mesg,rdct(lco+1:lcf)) ENDIF s(np) = rdtcopchem p(np) = copchem IF (wtopeflag.EQ.1) & WRITE(10,'(f5.3,2X,A1,A6)')s(np)*ratio(i),'G',p(np) ENDIF * write out: reaction with NO3 are thermal reaction (idreac=0 and does * not require labels) CALL rxwrit3(17,a1,a2,a3,a4,r,s,p,ar1,ar2,ar3, & f298,fratio,idreac,nlabel,xlabel,folow,fotroe) 300 CONTINUE IF (wtopeflag.EQ.1) WRITE(10,*)'end' * end of no3rad RETURN END