<        *(Cs + Cl*exp(-5420./T35(no_so)))

     <        *Cf*exp(-5420./T35(no_so))

      exf_BulkqSat = cvapor_fac/exp(cvapor_exp/tkelvin)

              ssq = ssq0*exp( lath*(ssq1-ssq2/Tsf) ) / p0

            swdk(i) =       rfac(jwtype)  * exp(facz/a1(jwtype))

     $                + (1.-rfac(jwtype)) * exp(facz/a2(jwtype))

             flux_u=caexport(i,j)*exp(-(depth_u-zbase)/zca)

                 flux_l=caexport(i,j)*exp(-(depth_l-zbase)/zca)

             lit=lit*exp(-atten)*(1.d0-Fice(i,j,bi,bj))

     &             exp( (24.2-0.085*t)*(pressc-1.0)/(83.143*tk) )

     &             exp( (16.4-0.040*t)*(pressc-1.0)/(83.143*tk) )

           akb(i,j,bi,bj)=exp((-8966.90 - 2890.53*sqrts - 77.942*s +

           akb(i,j,bi,bj) = akb(i,j,bi,bj)*exp(pfactor)

           ak1p(i,j,bi,bj) = exp(-4576.752*invtk + 115.525 - 

           ak2p(i,j,bi,bj) = exp(-8814.715*invtk + 172.0883 - 

           ak3p(i,j,bi,bj) = exp(-3070.75*invtk - 18.141 + 

           aksi(i,j,bi,bj) = exp(-8904.2*invtk + 117.385 - 

           akw(i,j,bi,bj) = exp(-13847.26*invtk + 148.9652 - 

           aks(i,j,bi,bj)=exp(-4276.1*invtk + 141.328 - 

           akf(i,j,bi,bj)=exp(1590.2*invtk - 12.641 + 1.525*sqrtis +

           ff(i,j,bi,bj) = exp(-162.8301 + 218.2968/tk100  +

           ak0(i,j,bi,bj) = exp(93.4517/tk100 - 60.2409 + 

           akb(i,j,bi,bj)=exp((-8966.90 - 2890.53*sqrts - 77.942*s +

           ak1p(i,j,bi,bj) = exp(-4576.752*invtk + 115.525 - 

           ak2p(i,j,bi,bj) = exp(-8814.715*invtk + 172.0883 - 

           ak3p(i,j,bi,bj) = exp(-3070.75*invtk - 18.141 + 

           aksi(i,j,bi,bj) = exp(-8904.2*invtk + 117.385 - 

           akw(i,j,bi,bj) = exp(-13847.26*invtk + 148.9652 - 

           aks(i,j,bi,bj)=exp(-4276.1*invtk + 141.328 - 

           akf(i,j,bi,bj)=exp(1590.2*invtk - 12.641 + 1.525*sqrtis +

           ff(i,j,bi,bj) = exp(-162.8301 + 218.2968/tk100  +

           ak0(i,j,bi,bj) = exp(93.4517/tk100 - 60.2409 + 

                 flux_l = flux_u*exp(dumrate)

           clr(i,j,k)=1.0-(fcld(i,j,k)*(1.0-exp(-1.66*taux)))

         tranal(i,j)=tranal(i,j)*exp(-1.66*taux)

           h2oexp(i,j,k,1) = exp(-xh*xkw(ib))

           conexp(i,j,k,1) = exp(-dcont(i,j,k)*xke(ib))

           co2exp(i,j,k,1,1)=exp(-xc*2.656e-5)

           co2exp(i,j,k,1,2)=exp(-xc*2.656e-3)

           n2oexp(i,j,k,1)=exp(-xc*6.31582e-2)

           n2oexp(i,j,k,1)=exp(-xc*5.35779e-2)

           ch4exp(i,j,k,1)=exp(-xc*5.80708e-3)

           ch4exp(i,j,k,1)=exp(-xc*6.29247e-2)

           comexp(i,j,k,1)=exp(-xc*2.18947e-5)

           comexp(i,j,k,1)=exp(-xc*4.74570e-5)

           cfcexp(i,j,k)=exp(-xf*fk1)

           cfcexp(i,j,k)=exp(-xf*fk2)

           h2oexp(i,j,k,1)=exp(-xx*0.10624)

           co2exp(i,j,k,1,1)=exp(-xx*2.656e-5)

           conexp(i,j,k,1)=exp(-dcont(i,j,k)*1.04995e+2)

           n2oexp(i,j,k,1)=exp(-xx*0.25238)

         dcont(i,j,k) = xx*exp(1800./ta(i,j,k)-6.081)+1.e-10

        n2o(L) = 120.* exp( (26.- hght(L)) / 5.69 ) * 1.e-9

        n2o(L) = 307.* exp( (16.- hght(L)) /10.47 ) * 1.e-9

        methane(L) = 1.7* exp( (14.- hght(L)) /19.16 ) * 1.e-6

      f = 1.0 - min( 1.0 _d 0, max(0.0 _d 0,1.0-exp(-s)) )

            theta(i,j,k,bi,bj)=exp(-0.5*( rD/20.E3 )**2 )

     & * (exp(futureTime/rampTime2)

     &   - exp(-futureTime/rampTime2))

     &   /(exp(futureTime/rampTime2)

     &  + exp(-futureTime/rampTime2))

            hinit=Dinf*(exp(-x)) - Dmax

     &     -  sqrt(gp_inflow*Dinf)*exp(-x)

      rexpt=1.0/exp(futureTime/rampTime2)

     &     +20. _d 0 *exp( -((xC(i,j,bi,bj)-xBump)/dxBump)**2

            smoothAbs_R4 = log(.5*(exp(x*sf)+exp(-x*sf)))*rsf

            smoothAbs_R8 = log(.5*(exp(x*sf)+exp(-x*sf)))*rsf

     &         = exp ( A1_12

     &         = exp ( A1_11

           tmpbulk = cvapor_fac_ice/exp(cvapor_exp_ice/Tsf)

      fswocn = fswpen * exp(-ksolar*hi)

     &       +(albSnow-albOldSnow)*exp(-0.2 _d 0*ageSnow)

     &       +(alb_land(i,j)-albSnow)*exp(-hSnow/hAlbSnow)