c*******************************************************************************

	subroutine getptrates(tc,rc,ysl,temp1tab,odtab,temptab,ntmx,
	1                     cgsptab,wcheml)

	IMPLICIT NONE
        include 'parameter.par' 
	include 'chem.com'


	integer ntmx

	real*8 temptab(ntmx)
	real*8 cgsptab(nsmx,ntmx)

	double precision temp1tab,odtab


	real*8 cgspl(nsmx),wspl(nsmx),cspl(nsmx),ropl(idim)
	real*8 qfac(idim),rr_f(idim)
	real*8 hspl(nsmx),ysl(nsmx),wcheml(nsmx)

	integer ksp,k,isp,ip3,nyp,nnn,kthb,ispthb,kfl,kfc,ifp
	integer ir,ngs,jtab

	double precision one,otc,tc,vlntemp,rc,rr_r,xik,crprod
	double precision cvp,cpprod,ctot,vt,tempref,vot,cthb,rklow
	double precision pr,prlog,fc,xp,fcent,fclog,xn,cprlog
	double precision flog,clog,pcor,rtab


c-------------------------------------------------------------------------------
	one = 1.0
c-------------------------------------------------------------------------------

	otc= 1./tc
	vlntemp = log(tc)

c ... get gibbs free energies

	call getcgsp(tc,temp1tab,odtab,temptab,ntmx,cgsptab,cgspl)

c-------------------------------------------------------------------------------
c-------------------------------------------------------------------------------
c-------------------------------------------------------------------------------
c ... sum wcheml(ksp=1,2,...,ns) over all reactions k=1,2,...,numreac

c ... initialize

	do ksp=1,ns
	   wspl(ksp)  = 0.0
	   cspl(ksp)  = rc*ysl(ksp)*ospwt(ksp)	! concentrations
	end do

c-------------------------------------------------------------------------------
c-------------------------------------------------------------------------------
c ... loop over all reactions to evaluate rate-of-progress: ropl(k)

	do k = 1, numreac
	      
c-------------------------------------------------------------------------------
c ... initialize qfac array, product of concentrations to be used in 
c ... third body and fall off reactions
	
	   qfac(k) = 1.0

c-------------------------------------------------------------------------------
c ... evaluate forward reaction rate, rr_f

	   if(itfdep(k).ne.0)then
	      rr_f(k) = capb(k)*exp(par(2,k)*vlntemp - cape(k)*otc)
	   else
	      rr_f(k) = capb(k)
	   endif

c-------------------------------------------------------------------------------
c ... nspec(k) > 0 iff reaction is reversible
c ... xik = SUM_i=1,ns nu(i,k)*G_ik
c ... evaluate (1/equilibrium_reaction_rate) = orr_e
c ... evaluate reverse reaction rate, rr_r
c ... kflagrev(k) = 1 -> reverse Arrhenius rate parameters are defined
c ...                    use these rather than orr_e to find rr_r

	   if(nspec(k).gt.0)then
	      if(kflagrev(k).eq.1)then
		 if(itrdep(k).ne.0)then
		   rr_r = rcapb(k)*exp(rexpo(k)*vlntemp-rcape(k)*otc)
		 else
		   rr_r = rcapb(k) 
		 endif
	      else
	         xik = 0.0
	         do isp=1,nreac(k)
	            xik = xik + anu(isp,k)*cgspl(nunk(isp,k))
	         end do
	         do ip3=4,nprod(k)+3
	            xik = xik + anu(ip3,k)*cgspl(nunk(ip3,k))
	         end do
	         rr_r  = rr_f(k)*exp(xik*otc+anusum(k)*vlntemp) 
	      endif
	   else
	      rr_r  = 0.0			! reaction irreversible-> rr_r=0
	   endif

c-------------------------------------------------------------------------------
c ... find reactants concentrations product PROD(c_j**nu_jk^prime)_j=1,ns
c ... it is assumed below that nyp is a positive integer (non-zero)
c ... which should be true in general

	   crprod = 1.0
	   do isp = 1,nreac(k)
	      nyp = - nu(isp,k)
	      cvp = cspl(nunk(isp,k))
	      do nnn=1,nyp
	         crprod = crprod * cvp
	      end do
	   end do

c-------------------------------------------------------------------------------
c ... for reverse reaction,
c ... find products concentrations product PROD(c_j**nu_jk^double-prime)_j=1,ns
c ... nspec(k) > 0 if reaction is reversible
c ... it is assumed below that nu(ip3,k) is a positive integer (non-zero)
c ... which should be true in general

	   if(nspec(k).gt.0)then
	      cpprod = 1.0
	      do ip3 = 4,nprod(k)+3
		 cvp = cspl(nunk(ip3,k))
		 do nnn=1,nu(ip3,k)
	            cpprod = cpprod * cvp
		 end do
	      end do
	   else
	      cpprod = 0.0
	   endif

c-------------------------------------------------------------------------------
c ... evaluate reaction rate-of-progress

	   ropl(k) = rr_f(k)*crprod - rr_r*cpprod

c-------------------------------------------------------------------------------

	end do

c ... end of loop over all reactions for rate-of-progress
c-------------------------------------------------------------------------------
c-------------------------------------------------------------------------------
c ... evaluate third-body-scaled reaction rates-of-progress
c ... nb. qfac is the product of concentrations, with possible 3rd
c ... body efficiencies.
c ... if the reaction is also a fall-off reaction with kfal=0, then
c ... qfac is just used to find pr below, and then reset to 1.0
c ... otherwise, it is actually used to scale the rate of progress 
c ... later on below

	if(nthb.gt.0)then

	   ctot = 0.0
	   do ksp=1,ns
	      ctot = ctot + cspl(ksp)
	   end do

	   do kthb = 1,nthb
	      k = ithb(kthb)
	      qfac(k) = ctot
	      do ispthb = 1,ntbs(kthb)
		 ksp  = nktb(ispthb,kthb)
		 qfac(k) = qfac(k)+(aik(ispthb,kthb)-one)*cspl(ksp)
	      end do
	   end do
	endif

c-------------------------------------------------------------------------------
c ... evaluate fall-off pressure corrections
c ... vt is dimensional temperature (K), vot = 1/vt 
c ... these are used in the computation of the dimensionless F

	if(nfal.gt.0)then

	   vt = tc*tempref
	   vot= 1.0/vt

	   do kfl = 1, nfal

	      k   = ifal(kfl)
              kfc = kfal(kfl)
	      ifp = ifop(kfl)

	      if(kfc.eq.0)then
		 cthb    = qfac(k)
		 qfac(k) = 1.0
	      else
		 cthb    = cspl(kfc)
	      endif

	      rklow = capbl(k)*exp(fpar(2,kfl)*vlntemp-capel(k)*otc)

	      pr    = rklow * cthb / rr_f(k)
	      prlog = log10(pr)

	      if(ifp.eq.1)then	                    ! Lindemann form

	         fc = 1.0	

	      else			            ! ifop = 2, 3, or 4

		 if(ifp.eq.2)then                   ! SRI form

		    xp = 1.0/(1.0+prlog*prlog)
		    fc = (
     +		          ( fpar(4,kfl) * exp(-fpar(5,kfl)*vot) +
     +                      exp(-vt/fpar(6,kfl)) 
     +                    )**xp 
     +                   ) * fpar(7,kfl) * vt**fpar(8,kfl)

	         else 			            ! 6 or 7-param TROE form 

		    fcent = (1.0-fpar(4,kfl))*exp(-vt/fpar(5,kfl)) +
     +                         fpar(4,kfl)*exp(-vt/fpar(6,kfl))

		    if(ifp.eq.4)then                ! 7-parameter TROE form

		       fcent = fcent + exp(-fpar(7,kfl)*vot)

		    endif

		    fclog  = log10(fcent)
		    xn     = 0.75 - 1.27*fclog
		    cprlog = prlog - (0.4+0.67*fclog)
		    clog   = cprlog/(xn-0.14*cprlog)
		    flog   = fclog/(1.0+clog*clog)
		    fc     = 10.0**flog

		 endif

	      endif

	      pcor = fc*pr/(1.0+pr)

	      ropl(k) = ropl(k) * pcor

	   end do 

	endif

c-------------------------------------------------------------------------------
c ... modify ropl due to 3rd bodies, given possible mod to qfac bec falloff

	if(nthb.gt.0)then
	   do kthb = 1,nthb
	      k = ithb(kthb)
	      ropl(k) = qfac(k) * ropl(k) 
	   end do
	endif

c-------------------------------------------------------------------------------
c ... evaluate contributions to reactants and products reaction rates

	do k=1,numreac   
	   do ir = 1, nreac(k)
	      ksp    = nunk(ir,k)
	      wspl(ksp) = wspl(ksp) + anu(ir,k)*ropl(k)	! reactants
	   end do
	   do ip3 = 4, nprod(k)+3
	      ksp    = nunk(ip3,k)
	      wspl(ksp) = wspl(ksp) + anu(ip3,k)*ropl(k) 	! products
	   end do
	end do

c-------------------------------------------------------------------------------
c ... final scaling

	do ksp=1,ns
	   wcheml(ksp) = wspl(ksp)*spwt(ksp)
	end do

c-------------------------------------------------------------------------------
c-------------------------------------------------------------------------------

        return
	end

c*******************************************************************************


c ... get G for each of the pure species at given (i,j)

	subroutine getcgsp(tc,temp1tab,odtab,temptab,ntmx,cgsptab,cgspl)

	IMPLICIT NONE 

        include 'parameter.par'
	include 'chem.com'


	integer jtab,jtp,k,ntmx

	real*8 cgspl(nsmx)

	real*8 temptab(ntmx)
	real*8 cgsptab(nsmx,ntmx)

	double precision tc,rtab,cgsv,temp1tab,odtab

c-------------------------------------------------------------------------------
c ... linear interpolation from table of g_k(T)
	    
	    jtab = int((tc-temp1tab)*odtab) + 1
	    rtab = (tc-temptab(jtab))*odtab
	    jtp  = jtab + 1

	    do k=1,ns
	       cgsv    = cgsptab(k,jtab)
	       cgspl(k) = cgsv + rtab*(cgsptab(k,jtp)-cgsv)     
	    end do
c-------------------------------------------------------------------------------

	return
	end