FLUENTUDF代码_颗粒(煤或生物质)床层多相流燃烧过程数值模拟-【第1部分】

1、Fluent-UDF_Coal-Combustion-Multiphase-Flow-Processes-Codes for various kinds of heterogeneous reactions【第一部分 -Part1】#include "udf.h"#include "stdio.h"#include "time.h"#define SMALL_S 1.e-29#define eps_g_small 1.#define eps_s_small 1.e-6#define spe_small 1.e-6#define spe

2、_small_comb 1.e-6#define TMAX 3000.#define TMIN 280.#define Rgas 1.987 /* cal/mol.K */ /*UNIVERSAL_GAS_CONSTANT;*/#define PR_NUMBER(cp,mu,k) (cp)*(mu)/(k)#define IP_HEAT_COEFF(vof,PvofS,k,nu,d) (vofP)*(vofS)*6.*(k)*(Nu)/(d)/(d)#define Carbon_Init_MF 0.541#define Volatile_Init_MF 0.418#define Moistur

3、e_Init_MF 0.026#define Ash_Init_MF (1.-Carbon_Init_MF-Volatile_Init_MF-Moisture_Init_MF)#define rho_c (1000.)#define rho_ash (2931)#define rho_liq_water (998.2)#define rho_volatile (1000.)#define solid_rho (1./(Carbon_Init_MF/rho_c+Volatile_Init_MF/rho_volatile+Moisture_Init_MF/rho_liq_water+Ash _In

4、it_MF/rho_ash)/* Functions Defination*/void SolidFuel_Reactant(cell_t c, Thread *t, Hetero_Reaction *hr, double* y_carbon, double* mol_weight);void volatile_mass_fractions();void SetSpeciesIndex();void read_c3m_data();double satPressure(double T); double Get_Phase_Index(Hetero_Reaction *hr);double T

5、urbulent_rr(cell_t c, Thread *t, Hetero_Reaction *r, real yiMAX_PHASESMAX_SPE_EQNS);double Mass_Transfer_Coeff(cell_t c, Thread *tp, Thread *ts); double heat_gunn_udf(cell_t c, Thread *ti, Thread *tj);double rr_combustion(cell_t c, Thread *t, Thread *ts, Thread *tp, double yi_O2, double y_ash, doubl

6、e y_carbon);double rr_steam_gasif(cell_t c, Thread *t, Thread *ts, Thread *tp, double p_h2o, double p_co, double p_h2, double y_carbon, double mol_weight, double* direction);double rr_co2_gasif(cell_t c, Thread *t, Thread *ts, Thread *tp, double y_co, double y_co2, double y_carbon, double mol_weight

7、, double* direction);double rr_h2_gasif(cell_t c, Thread *t, Thread *ts, Thread *tp, double y_h2, double y_ch4, double y_carbon, double mol_weight, double* direction);int IP_CH4 = 0, IS_CH4 = 0, IP_CO = 0, IS_CO = 0, IP_CO2 = 0, IS_CO2 = 0, IP_H2 = 0, IS_H2 = 0, IP_H2O = 0, IS_H2O = 0, IP_O2 = 0, IS

8、_O2 = 0, IP_H2S = 0, IS_H2S = 0, IP_CL2 = 0, IS_CL2 = 0, IP_NH3 = 0, IS_NH3 = 0, IP_N2 = 0, IS_N2 = 0, IP_TAR = 0, IS_TAR = 0, IP_C = 0, IS_C = 0, IP_VOL = 0, IS_VOL = 0, IP_MOISTURE = 0, IS_MOISTURE = 0, IP_ASH = 0, IS_ASH = 0, IP_SOOT = 0, IS_SOOT = 0, IP_PAH = 0, IS_PAH = 0, IP_OIL = 0, IS_OIL =

9、0, IP_ASH_R = 0, IS_ASH_R = 0, IP_C_R = 0, IS_C_R = 0, IP_C2H2 = 0, IS_C2H2 = 0, IP_SLURRY = 0, IS_SLURRY = 0, IP_C2H4 = 0, IS_C2H4 = 0, IP_C2H6 = 0, IS_C2H6 = 0, IP_C3H6 = 0, IS_C3H6 = 0, IP_C3H8 = 0, IS_C3H8 = 0, IP_SLURRY_D = 0, IS_SLURRY_D = 0, IP_SAND = 0, IS_SAND = 0, IP_ASH_D = 0, IS_ASH_D =

10、0;int current_c3m_solid_phase = 0;cxboolean init_flag= TRUE;cxboolean PCCL_Devol= FALSE;cxboolean MGAS_Devol= FALSE;cxboolean CPD_Devol= FALSE;cxboolean FGDVC_Devol= FALSE;cxboolean HPTR_Devol= FALSE;cxboolean MGAS_Moisture= FALSE;cxboolean PCCL_Moisture= FALSE;cxboolean MGAS_TarCracking= FALSE;cxbo

11、olean PCCL_2nd_Pyro= FALSE;cxboolean MGAS_Gasif= FALSE;cxboolean PCCL_Gasif= FALSE;cxboolean PCCL_TarCracking= FALSE;cxboolean MGAS_WGS= FALSE;cxboolean PCCL_soot_gasif= FALSE;cxboolean MGAS_char_combustion= FALSE;cxboolean PCCL_char_combustion= FALSE;cxboolean PCCL_soot_oxidation= FALSE;cxboolean T

12、AR_oxidation= FALSE;cxboolean MGAS_gas_phase_oxidation = FALSE;doubleavg_mf_h2o,avg_mf_co,avg_mf_h2,avg_mf_ch4,avg_mf_co2,avg_c,avg_volatile,avg_moisture,av g_ash;/* coal analysis variables */double fc_ar=0.,vm_ar=0.,ash_ar=0.,moist_ar=0.;double f_ep_a = 0.;double mwMAX_PHASESMAX_SPE_EQNS;double A1_

13、devolatilization=0.0, E1_devolatilization=0.0; /* pan : Oct 2012 */ double A2_devolatilization=0.0, E2_devolatilization=0.0; /* pan : Oct 2012 */ double A_tar_cracking=0.0 , E_tar_cracking=0.0;double A_steam_gasification=0.0, E_steam_gasification=0.0;/* pan : Oct 2012 */ /* pan : Oct 2012 */double K

14、_steam_gasification=0.0, N_steam_gasification=0.0;double Annealing_steam_gasification=1.0;double A_co2_gasification=0.0, E_co2_gasification=0.0;/* pan : Oct 2012 */* pan : Oct 2012 */double K_co2_gasification=0.0, N_co2_gasification=0.0;double Annealing_co2_gasification=1.0;double A_h2_gasification=

15、0.0, E_h2_gasification=0.0;double N_h2_gasification=0.0;double Annealing_h2_gasification=1.0;/* pan : Oct 2012 */ /* pan : Oct 2012 */double A_soot_steam_gasification=0.0, E_soot_steam_gasification=0.0; /* pan : oct 2012 */ double K_soot_steam_gasification=0.0, N_soot_steam_gasification=0.0; /* pan

16、: Oct 2012 */ double Annealing_soot_steam_gasification=0.0;/* pan : Oct 2012 */double A_soot_co2_gasification=0.0, E_soot_co2_gasification=0.0;/* pan : Oct 2012 */double K_soot_co2_gasification=0.0, N_soot_co2_gasification=0.0;/* pan : Oct 2012 */double Annealing_soot_co2_gasification=0.0;/* pan : O

17、ct 2012 */double A_soot_h2_gasification=0.0, E_soot_h2_gasification=0.0; /* pan : Oct 2012 */ double N_soot_h2_gasification=0.0; /* pan : Oct 2012 */ double Annealing_soot_h2_gasification=0.0; /* pan : Oct 2012 */double A_Soot_Combustion = 0.0, E_Soot_Combustion = 0.0;double A_c_combustion = 8710. /

18、* g/(atm.cmA2.s) */, E_c_combustion = 27000. /* cal/mole */; double Annealing_c_combustion=0.0 , N_c_combustion=0.0; /* pan : Oct 2012 */double A_moisture_release = 0.0, E_moisture_release = 0.0;double wg3 = 0.014;double Moisture_Flux;DEFINE_ADJUST(gasification,domain)if(init_flag)#if !RP_HOSTif(0)

19、SetSpeciesIndex();#endif #if !RP_NODEvolatile_mass_fractions();#endifhost_to_node_real_5(fc_ar,vm_ar,ash_ar,moist_ar,f_ep_a); host_to_node_int_4(PCCL_Devol,MGAS_Devol,MGAS_Moisture,PCCL_2nd_Pyro); host_to_node_int_3(MGAS_Gasif,PCCL_Gasif,PCCL_Moisture); host_to_node_int_3(CPD_Devol,FGDVC_Devol,HPTR_

20、Devol); /* pan : Oct 2012 */ host_to_node_int_3(PCCL_TarCracking,MGAS_WGS,PCCL_soot_gasif); /* pan : Oct 2012*/host_to_node_int_1(PCCL_char_combustion); /* pan : Oct 2012 */ host_to_node_int_3(MGAS_char_combustion,TAR_oxidation,MGAS_gas_phase_oxidation);/* pan : Oct 2012 */host_to_node_real_2(A_mois

21、ture_release, E_moisture_release);host_to_node_real_2(A1_devolatilization, E1_devolatilization); /* pan : Oct 2012 . added the "1" */host_to_node_real_2(A2_devolatilization, E2_devolatilization); /* pan : Oct 2012 */ host_to_node_real_2(A_tar_cracking, E_tar_cracking);host_to_node_real_2(A

22、_steam_gasification, E_steam_gasification);host_to_node_real_2(K_steam_gasificationN_steam_gasification);/*pan : Oct 2012 */host_to_node_real_1(Annealing_steam_gasification); /* pan : Oct 2012 */host_to_node_real_2(A_co2_gasification, E_co2_gasification);host_to_node_real_2(K_co2_gasification, N_co2

23、_gasification); /* pan : Oct 2012 */host_to_node_real_1(Annealing_co2_gasification);/* pan : Oct 2012 */host_to_node_real_2(A_h2_gasification, E_h2_gasification);host_to_node_real_2(Annealing_h2_gasification, N_h2_gasification); /* pan : Oct 2012 */host_to_node_real_2(A_soot_steam_gasification, E_so

24、ot_steam_gasification);/* pan :Oct 2012 */host_to_node_real_2(K_soot_steam_gasification, N_soot_steam_gasification);/* pan :Oct 2012 */host_to_node_real_1(Annealing_soot_steam_gasification);/* pan : Oct 2012 */host_to_node_real_2(A_soot_co2_gasification, E_soot_co2_gasification);/* pan :Oct 2012 */h

25、ost_to_node_real_2(K_soot_co2_gasification, N_soot_co2_gasification);/* pan :Oct 2012 */host_to_node_real_1(Annealing_soot_co2_gasification);/* pan : Oct 2012 */ host_to_node_real_2(A_soot_h2_gasification, E_soot_h2_gasification); /* pan : Oct 2012 */host_to_node_real_2(Annealing_soot_h2_gasificatio

26、n, N_soot_h2_gasification); /* pan : Oct 2012 */host_to_node_real_2(A_Soot_Combustion, E_Soot_Combustion); host_to_node_real_2(A_c_combustion, E_c_combustion); /* pan : Oct 2012 */host_to_node_real_2(N_c_combustion,Annealing_c_combustion);/* pan : Oct 2012 */host_to_node_real_2(wg3,Moisture_Flux);in

27、it_flag = FALSE; /* end init_flag */DEFINE_ON_DEMAND(Devol_and_Tar_Cracking)#if !RP_NODEvolatile_mass_fractions();#endif void SetSpeciesIndex()Domain *domain = Get_Domain(1);int n, ns, MAX_SPE_EQNS_PRIM = 0;Domain *subdomain;/*int n_phases = DOMAIN_N_DOMAINS(domain);*/* search all the species and sa

28、ved the Molecular Weight */ sub_domain_loop(subdomain, domain, n)Material *m_mat, *s_mat;if (DOMAIN_NSPE(subdomain) > 0) m_matPick_Material(DOMAIN_MATERIAL_NAME(subdomain),NULL); mixture_species_loop(m_mat,s_mat,ns)if (0 = strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"ch4") IP_CH4 = n; IS_CH4 =

29、 ns;if(n = 0) MAX_SPE_EQNS_PRIM +=1; else if (0 strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"co") IP_CO = n; IS_CO = ns;if(n = 0) MAX_SPE_EQNS_PRIM +=1; else if (0 strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"c2h2") IP_C2H2 = n; IS_C2H2 = ns;if(n = 0) MAX_SPE_EQNS_PRIM +=1; else if (0 strcmp(M

30、IXTURE_SPECIE_NAME(m_mat,ns),"c2h4") IP_C2H4 = n;6IS_C2H4 = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;else if strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"c2h6")IP_C2H6 = n;IS_C2H6 = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;else if strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"c3h6")IP_C3H6 = n;IS_C

31、3H6 = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;else if strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"c3h8")IP_C3H8 = n;IS_C3H8 = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;else if strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"co2") IP_CO2 = n;IS_CO2 = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;else if strcmp(MIXTURE_SPE

32、CIE_NAME(m_mat,ns),"h2")IP_H2 = n;IS_H2 = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;(0(0(0(0(08if(0elsestrcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"h2o")IP_H2O = n;IS_H2O = ns;if(n = 0)if(0MAX_SPE_EQNS_PRIM +=1; elsestrcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"o2")IP_O2 = n;IS_O2 = ns;if(n =

33、 0)if(0MAX_SPE_EQNS_PRIM +=1; elsestrcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"h2s")IP_H2S = n;IS_H2S = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;if(0else strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"cl2")IP_CL2 = n;IS_CL2 = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;if(0elsestrcmp(MIXTURE_SPECIE_NAME(m_mat,ns

34、),"nh3")IP_NH3 = n;IS_NH3 = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;if(0else strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"n2")IP_N2 = n;9IS_N2 = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1; else if strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"oil")IP_OIL = n;IS_OIL = ns; if(n = 0)MAX_SPE_EQNS_PRIM +

35、=1; else if strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"pah")IP_PAH = n;IS_PAH = ns; if(n = 0)MAX_SPE_EQNS_PRIM +=1; else if strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"tar")IP_TAR = n;IS_TAR = ns; if(n = 0)MAX_SPE_EQNS_PRIM +=1; else if strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"c")IP_

36、C = n;IS_C = ns;if(n = 0) MAX_SPE_EQNS_PRIM +=1; else if strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"c_recycle")IP_C_R = n;IS_C_R = ns; if(n = 0)MAX_SPE_EQNS_PRIM +=1;(0(0(0(0(011else if (0strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"soot")IP_SOOT = n;IS_SOOT = ns;if(n = 0) MAX_SPE_EQNS_PRIM

37、+=1;else if (0strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"volatile")IP_VOL = n;IS_VOL = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;else if (0strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"h2o<l>")IP_MOISTURE = n;IS_MOISTURE = ns;if(n = 0) MAX_SPE_EQNS_PRIM +=1;else if (0strcmp(MIXTURE_SPECIE_NAME(

38、m_mat,ns),"h2o<l>-slurry")IP_SLURRY = n;IS_SLURRY = ns;mwnns = MATERIAL_PROP(s_mat,PROP_mwi);else if (0 strcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"h2o<l>-dummy")IP_SLURRY_D = n;IS_SLURRY_D = ns;mwnns = MATERIAL_PROP(s_mat,PROP_mwi);else if (0 strcmp(MIXTURE_SPECIE_NAME(m_m

39、at,ns),"ash-coal")IP_ASH = n;12IS_ASH = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;else ifstrcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"ash-dummy") IP_ASH_D = n;IS_ASH_D = ns;else ifstrcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"ash-recycle") IP_ASH_R = n;IS_ASH_R = ns;if(n = 0)MAX_SPE_EQNS_PRIM

40、 +=1;else ifstrcmp(MIXTURE_SPECIE_NAME(m_mat,ns),"si<s>")IP_SAND = n;IS_SAND = ns;if(n = 0)MAX_SPE_EQNS_PRIM +=1;(0(0(0mwnns = MATERIAL_PROP(s_mat,PROP_mwi); elses_mat = Pick_Material(DOMAIN_MATERIAL_NAME(subdomain),NULL); mwn0 = MATERIAL_PROP(s_mat,PROP_mwi);DEFINE_HET_RXN_RATE(mois

41、ture_release,c,t,hr,mw,yi,rr,rr_t) Thread *pt = THREAD_SUB_THREADS(t);int index_phase = Get_Phase_Index(hr);Thread *ts = ptindex_phase; /* solid phase */14double prod = 0.0, Ts = C_T(c,ts);double Pt = MAX(0.1,(op_pres+C_P(c,t);double Tsat = 1./(0.0727/log(Pt/ 611.) - 0.0042) + 273.;*rr = 0;/* Set th

42、e phase and species indices. Ash species index is initialized to zero, with all other indices.Ash species index is used as a flag to execute SetSpeciesIndex only once. This is done by the firstreaction, defined in the heterogeneous reaction panel in FLUENT GUI.*/if(IS_ASH = 0) SetSpeciesIndex();if(C

43、_VOF(c, ts) > 0.0 && Ts >= Tsat) if(MGAS_Moisture)prod =yiIP_MOISTUREIS_MOISTURE*C_R(c,ts)/mwIP_MOISTUREIS_MOISTURE; /* kg-mol/mT*/ A_moisture_release = A1_devolatilization;E_moisture_release = E1_devolatilization;*rr = A_moisture_release*exp(-E_moisture_release/(Rgas*Ts) * prod*C_VOF(

44、c, ts); /* kmol/(m3.s) */if(PCCL_Moisture)*rr = 6.* C_VOF(c, ts) / C_PHASE_DIAMETER(c,ts) * Moisture_Flux/mwIP_MOISTUREIS_MOISTURE ; /* kmol/(m3.s) */DEFINE_HET_RXN_RATE(devolatilization,c,t,hr,mw,yi,rr,rr_t) Thread *pt = THREAD_SUB_THREADS(t);int index_phase = Get_Phase_Index(hr);Thread *ts = ptind

45、ex_phase; /* solid phase */double x0_star = 0., x_star =0.;double Ts = C_T(c,ts);*rr = 0.0;/* Set the phase and species indices. Ash species index is initialized to zero, with all other indices.Ash species index is used as a flag to execute SetSpeciesIndex only once.This is done bythe firstreaction,

46、 defined in the heterogeneous reaction panel in FLUENT GUI.*/if(IS_ASH = 0)SetSpeciesIndex();/*Volatile Matter -> c1 Tar + c2 CO + c3 CO2 + c4 CH4 + c5 H2 + c6 H2O + c7 H2S + c8 NH3 */if(C_YI(c,ts,IS_MOISTURE) < spe_small) if(MGAS_Devol)if(Ts<1223)x0_star = pow(867.2/MAX(Ts-273.),1.),3.914)

47、 / 100.; x_star = solid_rho * (fc_ar/100. + vm_ar/100.)/C_R(c,ts) * x0_star;A2_devolatilization = 0.0;E2_devolatilization = 0.0;if(PCCL_Devol | CPD_Devol | FGDVC_Devol | HPTR_Devol) && yiIP_VOLIS_VOL > x_star)*rr = (A1_devolatilization *exp(-E1_devolatilization/(1.987*Ts) + A2_devolatiliz

48、ation *exp(-E2_devolatilization/(1.987*Ts)*(yiIP_VOLIS_VOL-x_star)*C_VOF(c,ts)*C_R(c,ts)/mwIP_VOLIS_VOL;DEFINE_HET_RXN_RATE(tar_comb,c,t,r,mw,yi,rr,rr_t)Thread *pt = THREAD_SUB_THREADS(t);Thread *tp = pt0; /* gas phase */double T_g = MAX(TMIN, C_T(c,tp);double rho_g = C_R(c,tp)rr = 0.0;1.e-3; /* g/c

49、mA3 */double rr_turb = 1e+20;*rr = 0.0;/* Set the phase and species indices. Ash species index is initialized to zero, with all other indices.Ash species index is used as a flag to execute SetSpeciesIndex only once. This is done by the firstreaction, defined in the heterogeneous reaction panel in FL

50、UENT GUI.*/if(IS_ASH = 0)SetSpeciesIndex();if (rp_ke)rr_turb = Turbulent_rr(c, t, r, yi);if(yiIP_O2IS_O2 > spe_small_comb)double tmp_exp, tmp_o2, tmp_tar;T_g = MIN(T_g, TMAX);tmp_exp = 3.8e11 * exp(-30000./(Rgas*T_g); /* 30000 -> 60000*/tmp_o2 = pow(rho_g*yiIP_O2IS_O2/mwIP_O2IS_O2,1.5); tmp_ta

51、r = pow(rho_g*yiIP_TARIS_TAR/mwIP_TARIS_TAR,0.25);*rr = tmp_exp * tmp_o2 * tmp_tar * C_VOF(c,tp); /* mol/cmA3.s */*rr *= 1000.; /* kmol/(mA3 .s) */*rr = MIN(*rr, rr_turb);DEFINE_HET_RXN_RATE(tar_cracking,c,t,r,mw,yi,rr,rr_t)Thread *pt = THREAD_SUB_THREADS(t);Thread *tp = pt0; /* gas phase */double p

52、rod = 0.0;double T_g = MAX(TMIN, C_T(c,tp);double rr_turb = 1e+20;/* Set the phase and species indices. Ash species index is initialized to zero, with all other indices.Ash species index is used as a flag to execute SetSpeciesIndex only once. This is done by the firstreaction, defined in the heterog

53、eneous reaction panel in FLUENT GUI.*/if(IS_ASH = 0) SetSpeciesIndex();if (rp_ke) rr_turb = Turbulent_rr(c, t, r, yi);if(yiIP_TARIS_TAR > spe_small)prod = yiIP_TARIS_TAR*C_R(c,tp)*C_VOF(c,tp)/mwIP_TARIS_TAR;*rr = A_tar_cracking*exp(-E_tar_cracking/(1.987*T_g)* prod*C_VOF(c, tp); /* kmol/(m3.s) */

54、DEFINE_HET_RXN_RATE(co_comb,c,t,r,mw,yi,rr,rr_t)Thread *pt = THREAD_SUB_THREADS(t);Thread *tp = pt0; /* gas phase */double T_g = MAX(TMIN, C_T(c,tp);double rho_g = C_R(c,tp)*1.e-3; /* g/cmA3 */double rr_turb = 1e+20;*rr = 0;/* Set the phase and species indices. Ash species index is initialized to ze

55、ro, with all other indices.Ash species index is used as a flag to execute SetSpeciesIndex only once. This is done by the firstreaction, defined in the heterogeneous reaction panel in FLUENT GUI.*/if(IS_ASH = 0)SetSpeciesIndex();if(yiIP_O2IS_O2 > spe_small_comb)double tmp_exp, p_o2, p_co, p_h2o;if (rp_ke)rr_turb = Turbulent_rr(c, t, r, yi);T_g = MIN(T_g, TMAX);tmp_exp = 3.98e+14 * exp(-40000./(Rgas*T_g); /*40000 -> 80000 */p_o2 = pow(rho_g*yiIP_O2IS_O2/mwIP_O2IS_O2,0.2