CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   FLUENT (http://www.cfd-online.com/Forums/fluent/)
-   -   NOx UDF (http://www.cfd-online.com/Forums/fluent/44944-nox-udf.html)

Abhishek Asthana June 6, 2007 16:54

NOx UDF
 
Hi everyone,

I want to use my own rate for NOx formation in a UDF.

my NO rate= abs(NO Thermal)+abs(NO Prompt)+abs(NO Fuel)-abs(NO Reburn)

I do not understand how to get the NO thermal rate (and the other rates too)in Fluent.. I can see in result s but do not know how to have that value. Can anyonehelp me with this macro? Any pointers given will be helpful too. thanks in advance.

Allan Walsh June 7, 2007 17:33

Re: NOx UDF
 
As you are seeing, NOx predictions are quite complicated. I'm guessing you are looking at NOx from coal. If you are just looking for something representative, why not use the built in NOx stuff.

We have done our own NOx macros using the species transport equations but it has taken several years - not continuous of course. But, you need to start with reaction rates from experimental data, work out all of the kinetics (and likely turbulence reaction parmaters), write your code, model it, check against measured data, and then keep repeating as necessary.

Good luck.

cg June 9, 2007 07:25

Re: NOx UDF
 
Here is a UDF that will reproduce Fluent's default NOx model:

/************************************************** **********

UDF example of User-Defined NOx Rate

If used with the "replace with udf " radio buttons activated,

this udf will exactly reproduce the default fluent NOx

rates for thermal, prompt, fuel and n2o-interm pathways.

The flag "NOx->nox_io_pdf == IN_PDF" should allways be used

for rates other than that from char N, so that if requested,

the contributions will be pdf integrated. Any contribution

from char must be included within an switch statement of the

form "NOx->nox_io_pdf == OUT_PDF".

Prepared by: Anura Perera, FLUENT(UK) ************************************************** ***********/

#include "udf.h"

#define SMALL_S 1.e-29 #define SMALL_N 1.e-6

/* const real R_Gas_C = 1.9872; Universal gas constant in cal/mole/K */ const real P_Atm = 101325.; /* unit atm in Pa */

DEFINE_NOX_RATE(user_nox, c, t, NOx) {

NOx->fwdrate = 0.0;

NOx->revrate = 0.0;

switch (NOx->nox_io_pdf) {

case IN_PDF:

/* Source terms other than those from char must be included here */

if (NOX_EQN(NOx) == EQ_NO) {

/* Thermal NOx */

if (NOx->thermal_nox && NOx->thermal_udf_replace) {

real kf1, kr1, kf2, kr2, kf3;

real o_eq=0.0, oh_eq=0.0;

real s1, s2, s3, rf, rr;

Rate_Const K_Z_F[3] = {{1.80e8, 0.0, 38370.0},

{1.80e4, 1.0, 4680.0},

{7.10e7, 0.0, 450.0}};

Rate_Const K_Z_R[3] = {{3.80e7, 0.0, 425.0},

{3.81e3, 1.0, 20820.0},

{1.70e8, 0.0, 24560.0}};

Rate_Const K_O[2] = {{36.64, 0.5, 27123.0},

{3.97e5, -0.5, 31090.0}};

Rate_Const K_OH = {2.129e2, -0.57, 4595.0};

kf1 = ARRH(NOx, K_Z_F[0]);

kr1 = ARRH(NOx, K_Z_R[0]);

kf2 = ARRH(NOx, K_Z_F[1]);

kr2 = ARRH(NOx, K_Z_R[1]);

kf3 = ARRH(NOx, K_Z_F[2]);

/* determine O concentration */

if (NOx->o_mode == PARTIAL_EQUILIBRIUM) {

/* use partial equilibrium O model (Warnatz data) */

o_eq = ARRH(NOx, K_O[0]);

}

else if (NOx->o_mode == FULL_EQUILIBRIUM) {

/* use full equilibrium O model (Westenburg) */

o_eq = ARRH(NOx, K_O[1]);

}

else if (NOx->o_mode == LOCAL_MASS_FRACTION) {

/* use local [O] mass fraction */

o_eq = MOLECON(NOx, O);

}

if ((NOx->o_mode == PARTIAL_EQUILIBRIUM || NOx->o_mode == FULL_EQUILIBRIUM))

o_eq *= sqrt(MOLECON(NOx, O2));

/* determine OH concentration */

/* if local OH concentration is zero then calculate it from [O] and [H2O] */

if (NOx->oh_mode == PARTIAL_EQUILIBRIUM) {

/* (Westbrook&Baulch -- Andre Peters) */

oh_eq = ARRH(NOx, K_OH);

oh_eq *= sqrt(o_eq*MOLECON(NOx, H2O));

}

else if (NOx->oh_mode == LOCAL_MASS_FRACTION) {

/* use local OH mass fraction */

oh_eq = MOLECON(NOx, OH);

}

else if (NOx->oh_mode == 0) {

/* do not include OH in NO calculation */

oh_eq = 0.0;

}

rf = rr = 0.0;

s2 = s3 = 0.0;

s1 = kf2*MOLECON(NOx, O2);

if (s1 < SMALL_N) {

rf = 2.*o_eq*kf1*MOLECON(NOx, N2);

rr = 0.0;

}

else {

s2 = 1.+ kf3*oh_eq/s1;

s3 = s1 + kr1*MOLECON(NOx, NO)/s2;

rf = 2.*o_eq*kf1*MOLECON(NOx, N2)*s1/s3;

rr = 2.*o_eq*kr1*MOLECON(NOx, NO)*kr2*MOLECON(NOx, NO)/s3;

}

NOx->fwdrate += rf;

NOx->revrate += -rr;

}

/* Prompt NOx */

if (NOx->prompt_nox && NOx->prompt_udf_replace) {

real f,rf;

Rate_Const K_PM = {6.4e6, 0.0, 36483.49436};

f = 4.75 + 0.0819*NOx->c_number -

23.2*NOx->equiv_ratio + 32.0*pow(NOx->equiv_ratio,2.) - 12.2*pow(NOx->equiv_ratio,3.);

rf = ARRH(NOx, K_PM);

rf *= pow((NOx->uni_R*NOx->temp_m/NOx->press), (1+NOx->oxy_order));

rf *= pow(MOLECON(NOx, O2), NOx->oxy_order);

rf *= MOLECON(NOx, N2)*MOLECON(NOx, FUEL);

NOx->fwdrate += f*rf;

}

/* Fuel NOx */

if (NOx->fuel_nox && NOx->fuel_udf_replace) {

real rf=0.0, rr=0.0;

real k1, k2, k3, k4;

real s1=0.0, s2=0.0;

real rf_sncr=0.0, rr_sncr=0.0;

Rate_Const K_FL[5] = {{1.0e10, 0.0, 33715.781},

{3.0e12, 0.0, 30193.23671},

{4.0e6, 0.0, 16103.05958},

{1.8e8, 0.0, 13586.95652},

{2.3e2, 0.0, 17159.82287}};

k1 = ARRH(NOx, K_FL[0]);

k2 = ARRH(NOx, K_FL[1]);

k3 = ARRH(NOx, K_FL[2]);

k4 = ARRH(NOx, K_FL[3]);

if(NOx->solve_hcn) {

/* NO source from HCN */

s1 = k1*MOLECON(NOx, HCN)*pow(MOLECON(NOx, O2), NOx->oxy_order)*

pow((NOx->uni_R*NOx->temp_m/NOx->press), (1.+NOx->oxy_order));

s2 = k2*MOLECON(NOx, HCN)*MOLECON(NOx, NO)*

pow((NOx->uni_R*NOx->temp_m/NOx->press), 2.);

}

if(NOx->solve_nh3) {

/* NO source from NH3 */

if ( NOx->nh3_switch == 1 ) {

if (!NOx->sncr_nox) nox_nh3_rate(NOx, &rf_sncr, &rr_sncr);

}

else {

s1 += k3*MOLECON(NOx, NH3)*pow(MOLECON(NOx, O2), NOx->oxy_order)*

pow((NOx->uni_R*NOx->temp_m/NOx->press), (1.+NOx->oxy_order));

s2 += k4*MOLECON(NOx, NH3)*MOLECON(NOx, NO)*

pow((NOx->uni_R*NOx->temp_m/NOx->press), 2.);

}

}

rf = s1*(NOx->press/(NOx->uni_R*NOx->temp_m)) + rf_sncr;

rr = -s2*(NOx->press/(NOx->uni_R*NOx->temp_m)) + rr_sncr;

if (NOx->fuel_n_conv == 2) {

switch (NOx->nox_type) {

case FUEL_G:

rf += NOx->r_fuel*NOx->Yn_fuel*NOx->fueln_nox_frac*1000./NOx->sp[N].mw;

break;

case FUEL_L:

rf += NOx->r_gas*NOx->Yn_fuel*NOx->fueln_nox_frac*1000./(NOx->sp[N].mw*NOx->cell_V);

break;

case FUEL_S:

rf += NOx->r_volatile*NOx->Yn_volatile*NOx->fueln_nox_frac*1000./(NOx->sp[N].mw*NOx->cell_V);

break;

}

}

NOx->fwdrate += rf;

NOx->revrate += rr;

}

/* N2O Path NOx */

if (NOx->n2opath_nox && NOx->n2opath_udf_replace) {

real rf, rr;

real k3, k4;

real s1, s2;

Rate_Const K_N2O_F[2] = {{4.44e32, -8.358, 28234.0},

{2.90e7, 0.000, 11651.0}};

Rate_Const K_N2O_R[2] = {{4.00e8, 0.000, 28234.0},

{1.45e-29, 9.259, 11651.0}};

k3 = ARRH(NOx, K_N2O_F[1]);

k4 = ARRH(NOx, K_N2O_R[1]);

s1 = 2.0*k3*MOLECON(NOx, N2O)*MOLECON(NOx, O2);

s2 = 2.0*k4*pow(MOLECON(NOx, NO), 2.);

rf = s1;

rr = -s2;

NOx->fwdrate += rf;

NOx->revrate += rr;

}

}

else if (NOX_EQN(NOx) == EQ_HCN) {

/* Fuel NOx */

if (NOx->fuel_nox && NOx->fuel_udf_replace) {

real rf=0.0, rr=0.0;

real k1, k2;

real s1=0.0, s2=0.0;

Rate_Const K_FL[5] = {{1.0e10, 0.0, 33715.781},

{3.0e12, 0.0, 30193.23671},

{4.0e6, 0.0, 16103.05958},

{1.8e8, 0.0, 13586.95652},

{2.3e2, 0.0, 17159.82287}};

k1 = ARRH(NOx, K_FL[0]);

k2 = ARRH(NOx, K_FL[1]);

/* HCN production */

switch (NOx->nox_type) {

case FUEL_G:

rf = NOx->r_fuel*NOx->Yn_fuel*NOx->fueln_hcn_frac*1000./NOx->sp[N].mw;

break;

case FUEL_L:

rf = NOx->r_gas*NOx->Yn_fuel*NOx->fueln_hcn_frac*1000./(NOx->sp[N].mw*NOx->cell_V);

break;

case FUEL_S:

if (NOx->fuel_n_conv == 0 || NOx->fuel_n_conv == 2) {

rf = NOx->r_volatile*NOx->Yn_volatile*NOx->fueln_hcn_frac*1000./(NOx->sp[N].mw*NOx->cell_V);

}

break;

default:

rf = 0.0;

}

/* HCN depletion */

s1 = k1*MOLECON(NOx, HCN)*pow(MOLECON(NOx, O2), NOx->oxy_order)*

pow((NOx->uni_R*NOx->temp_m/NOx->press), (1.+NOx->oxy_order));

s2 = k2*MOLECON(NOx, HCN)*MOLECON(NOx, NO)*

pow((NOx->uni_R*NOx->temp_m/NOx->press), 2.);

rr = -(s1+s2)*(NOx->press/(NOx->uni_R*NOx->temp_m));

NOx->fwdrate += rf;

NOx->revrate += rr;

}

}

else if (NOX_EQN(NOx) == EQ_NH3) {

/* Fuel NOx */

if (NOx->fuel_nox && NOx->fuel_udf_replace) {

real rf=0.0, rr=0.0;

real k3, k4;

real s3=0.0, s4=0.0;

real rf_sncr=0.0, rr_sncr=0.0;

Rate_Const K_FL[5] = {{1.0e10, 0.0, 33715.781},

{3.0e12, 0.0, 30193.23671},

{4.0e6, 0.0, 16103.05958},

{1.8e8, 0.0, 13586.95652},

{2.3e2, 0.0, 17159.82287}};

k3 = ARRH(NOx, K_FL[2]);

k4 = ARRH(NOx, K_FL[3]);

/* NH3 production */

switch (NOx->nox_type) {

case FUEL_G:

rf = NOx->r_fuel*NOx->Yn_fuel*NOx->fueln_nh3_frac*1000./NOx->sp[N].mw;

break;

case FUEL_L:

rf = NOx->r_gas*NOx->Yn_fuel*NOx->fueln_nh3_frac*1000./(NOx->sp[N].mw*NOx->cell_V);

break;

case FUEL_S:

if (NOx->fuel_n_conv == 1 || NOx->fuel_n_conv == 2) {

rf = NOx->r_volatile*NOx->Yn_volatile*NOx->fueln_nh3_frac*1000./(NOx->sp[N].mw*NOx->cell_V);

}

break;

default:

rf = 0.0;

}

/* NH3 depletion */

if ( NOx->nh3_switch == 1 ) {

if (!NOx->sncr_nox) nox_nh3_rate(NOx, &rf_sncr, &rr_sncr);

rr = rr_sncr;

}

else {

s3 = k3*MOLECON(NOx, NH3)*pow(MOLECON(NOx, O2), NOx->oxy_order)*

pow((NOx->uni_R*NOx->temp_m/NOx->press), (1+NOx->oxy_order));

s4 = k4*MOLECON(NOx, NH3)*MOLECON(NOx, NO)*

pow((NOx->uni_R*NOx->temp_m/NOx->press), 2.);

rr = -(s3+s4)*(NOx->press/(NOx->uni_R*NOx->temp_m));

}

NOx->fwdrate += rf;

NOx->revrate += rr;

}

}

else if (NOX_EQN(NOx) == EQ_N2O) {

if (NOx->n2opath_udf_replace) {

real rf, rr;

real k1, k2, k3, k4;

real o_eq=0.0, oh_eq=0.0, thirdb, n2o_qs=0.0;

real s1, s2, s3, s4;

Rate_Const K_N2O_F[2] = {{4.44e32, -8.358, 28234.0},

{2.90e7, 0.000, 11651.0}};

Rate_Const K_N2O_R[2] = {{4.00e8, 0.000, 28234.0},

{1.45e-29, 9.259, 11651.0}};

Rate_Const K_O[2] = {{36.64, 0.5, 27123.0},

{3.97e5, -0.5, 31090.0}};

Rate_Const K_OH = {2.129e2, -0.57, 4595.0};

k1 = ARRH(NOx, K_N2O_F[0]);

k2 = ARRH(NOx, K_N2O_R[0]);

k3 = ARRH(NOx, K_N2O_F[1]);

k4 = ARRH(NOx, K_N2O_R[1]);

/* determine O concentration */

if (NOx->n2opath_o_mode == PARTIAL_EQUILIBRIUM) {

/* use partial equilibrium O model (Warnatz data) */

o_eq = ARRH(NOx, K_O[0]);

}

else if (NOx->n2opath_o_mode == FULL_EQUILIBRIUM) {

/* use full equilibrium O model (Westenburg) */

o_eq = ARRH(NOx, K_O[1]);

}

else if (NOx->n2opath_o_mode == LOCAL_MASS_FRACTION) {

/* use local [O] mass fraction */

o_eq = MOLECON(NOx, O);

}

if ((NOx->n2opath_o_mode == PARTIAL_EQUILIBRIUM || NOx->n2opath_o_mode == FULL_EQUILIBRIUM))

o_eq *= sqrt(MOLECON(NOx, O2));

/* determine OH concentration */

/* if local OH concentration is zero then calculate it from [O] and [H2O] */

if (NOx->n2opath_oh_mode == PARTIAL_EQUILIBRIUM) {

/* (Westbrook&Baulch -- Andre Peters) */

oh_eq = ARRH(NOx, K_OH);

oh_eq *= sqrt(o_eq*MOLECON(NOx, H2O));

}

else if (NOx->n2opath_oh_mode == LOCAL_MASS_FRACTION) {

/* use local OH mass fraction */

oh_eq = MOLECON(NOx, OH);

}

else if (NOx->n2opath_oh_mode == 0) {

/* do not include OH in NO calculation */

oh_eq = 0.0;

}

thirdb = (NOx->press/(NOx->uni_R*NOx->temp_m)); /* general third-body concentration mol/m^3 */

if (NOx->n2opath_n2o_model == 0) {

s1 = k1*MOLECON(NOx, N2)*o_eq*thirdb;

s2 = k4*pow(MOLECON(NOx, NO), 2.);

s3 = k2*thirdb;

s4 = k3*o_eq;

if ( (s3+s4) != 0.0 ) {

/* N2O mole concentration (moles/m^3) */

n2o_qs = (s1+s2)/(s3+s4);

}

/* Update N2O mass fraction when Quasi-steady assumption is made */

C_POLLUT(c,t,3) = n2o_qs*NOx->sp[N2O].mw/NOx->den/1000.0;

}

else {

s1 = k1*MOLECON(NOx, N2)*o_eq*thirdb;

s2 = k4*pow(MOLECON(NOx, NO), 2.);

s3 = k2*MOLECON(NOx, N2O)*thirdb;

s4 = k3*MOLECON(NOx, N2O)*o_eq;

rf = s1+s2;

rr = -(s3+s4);

NOx->fwdrate += rf;

NOx->revrate += rr;

}

}

}

break;

case OUT_PDF:

/* Char Contributions, that do not go into pdf loop must be included here */

if (NOX_EQN(NOx) == EQ_NO) {

/* Fuel NOx */

if (NOx->fuel_udf_replace) {

real rf=0.0, rr=0.0;

real k5;

real s5;

Rate_Const K_FL[5] = {{1.0e10, 0.0, 33715.781},

{3.0e12, 0.0, 30193.23671},

{4.0e6, 0.0, 16103.05958},

{1.8e8, 0.0, 13586.95652},

{2.3e2, 0.0, 17159.82287}};

if (NOx->char_n_conv == 0 || NOx->char_n_conv == 3) {

rf += NOx->r_char*NOx->Yn_char*NOx->charn_nox_frac*1000./(NOx->sp[N].mw*NOx->cell_V);

}

k5 = (K_FL[4].A/P_Atm)*pow(NOx->temp_m, K_FL[4].B)*exp(-K_FL[4].C/NOx->temp_m);

s5 = k5*NOx->den*NOx->uni_R*NOx->temp_m*NOx->sp[NO].yi;

rr += -s5*NOx->bet_area*NOx->part_con;

NOx->fwdrate += rf;

NOx->revrate += rr;

}

}

else if (NOX_EQN(NOx) == EQ_HCN) {

/* Fuel NOx */

if (NOx->fuel_udf_replace) {

real rf=0.0, rr=0.0;

if (NOx->char_n_conv == 1 || NOx->char_n_conv == 3) {

rf = NOx->r_char*NOx->Yn_char*NOx->charn_hcn_frac*1000./(NOx->sp[N].mw*NOx->cell_V);

}

NOx->fwdrate += rf;

NOx->revrate += rr;

}

}

else if (NOX_EQN(NOx) == EQ_NH3) {

/* Fuel NOx */

if (NOx->fuel_udf_replace) {

real rf=0.0, rr=0.0;

rf = rr = 0.0; /* char contributions are not included in PDF integral */

if (NOx->char_n_conv == 2 || NOx->char_n_conv == 3) {

rf = NOx->r_char*NOx->Yn_char*NOx->charn_nh3_frac*1000./(NOx->sp[N].mw*NOx->cell_V);

}

NOx->fwdrate += rf;

NOx->revrate += rr;

}

}

break;

default:

;

}

}


Abhishek Asthana June 11, 2007 06:08

Re: NOx UDF
 
Thanks CG,

ur macro has been a great help. i'll let u know how its workin n the results.

Abhishek


All times are GMT -4. The time now is 00:04.