A.3 Example

The following compiled UDF, named arrh, defines an Arrhenius-type reaction rate. The rate exponents are assumed to be same as the stoichiometric coefficients.

#include "udf.h"

static const real Arrhenius = 1.e15;
static const real E_Activation = 1.e6;
#define SMALL_S 1.e-29

DEFINE_HET_RXN_RATE(arrh,c,t,hr,mw,yi,rr,rr_t)
{
     Domain **domain_reactant = hr->domain_reactant;
     real *stoich_reactant = hr->stoich_reactant;
     int *reactant = hr->reactant;
     int i;
     int sp_id;
     int dindex;
     Thread *t_reactant;
     real ci;
     real T = 1200.; /* should obtain from cell */

     /* instead of compute rr directly, compute log(rr) and then
        take exp */

     *rr = 0;
     for (i=0; i < hr->n_reactants; i++)
       {
         sp_id = reactant[i]; /* species ID to access mw and yi */

         if (sp_id == -1) sp_id = 0; /* if phase does not have species,
                                  mw, etc. will be stored at index 0 */

         dindex = DOMAIN_INDEX(domain_reactant[i]);
                              /* domain index to access mw & yi */

         t_reactant = THREAD_SUB_THREAD(t,dindex);

         /* get conc. */
         ci = yi[dindex][sp_id]*C_R(c,t_reactant)/mw[dindex][sp_id];

         ci = MAX(ci,SMALL_S);

         *rr += stoich_reactant[i]*log(ci);
       }

     *rr += log(Arrhenius + SMALL_S) -
                  E_Activation/(UNIVERSAL_GAS_CONSTANT*T);

     /* 1.e-40 < rr < 1.e40 */
     *rr = MAX(*rr,-40);
     *rr = MIN(*rr,40);

     *rr = exp(*rr);
}