1 #include "burner.hpp"
2 using namespace std;
3
4 /*----------------------------------------------------------*/
5 /* arrondi pour ne conserver que n chiffres significatifs */
6 /*----------------------------------------------------------*/
arrondi(double x,int n)7 double arrondi ( double x , int n ) {
8 if (fabs(x) < EPS)
9 return 0.0;
10 double m = pow ( 10 , ceil(-log10(x)) + n - 1 );
11 return round(m*x)/m;
12 }
13
14 /*---------------------------------------------------------------*/
burner(int nb,chemical ** chem)15 burner::burner ( int nb , chemical ** chem ) {
16
17 rem_nb = nb;
18
19 NO=NO2=N2O=CO=0.0;
20 m = new double[nb];
21 can_burn = new bool[nb];
22
23 // combustion.prop :
24 // -----------------
25 // 64-17-5 3 2 3
26 // 74-82-8 2 1 2
27 // 1333-74-0 0.5 0 1
28 // 100-42-5 10 8 4
29 // 74-85-1 3 2 2
30 // 108-88-3 9 7 4
31 // 100-41-4 10.5 8 5
32 // 71-43-2 7.5 6 3
33
34 for ( i = 0 ; i < nb ; i++ ) {
35
36 can_burn[i] = false;
37
38 if ( chem[i]->CAS == "64-17-5" ||
39 chem[i]->CAS == "74-82-8" ||
40 chem[i]->CAS == "1333-74-0" ||
41 chem[i]->CAS == "100-42-5" ||
42 chem[i]->CAS == "74-85-1" ||
43 chem[i]->CAS == "108-88-3" ||
44 chem[i]->CAS == "100-41-4" ||
45 chem[i]->CAS == "71-43-2" )
46 can_burn[i] = true;
47 }
48
49 O2 = new chemical ("7782-44-7");
50 N2 = new chemical ("7727-37-9");
51 CO2 = new chemical ("124-38-9");
52 H2O = new chemical ("7732-18-5");
53
54 // Construct the rx array;
55 rx = new combrx * [nb];
56 for ( i = 0 ; i < nb ; i++ ) {
57 if ( can_burn[i] )
58 rx[i] = new combrx ( chem[i]->CAS );
59 else
60 rx[i] = NULL;
61 }
62 }
63
64 /*---------------------------------------------------------------*/
~burner(void)65 burner::~burner ( void ) {
66
67 delete [] m;
68 delete [] can_burn;
69
70 for ( i = 0 ; i < rem_nb ; i++ )
71 if (rx[i])
72 delete rx[i];
73 delete [] rx;
74
75 delete O2;
76 delete N2;
77 delete CO2;
78 delete H2O;
79 }
80
81 /*---------------------------------------------------------------*/
solve(double * y)82 bool burner::solve(double * y)
83 {
84 OK=true;
85 //perform mass balance (neglect pollutants flows)
86 out->m = 0.0;
87 for(i=0;i<in->nb;i++)
88 {
89 if (!can_burn[i]) {
90 out->chem[i]->m = in->chem[i]->m;
91 out->m+=out->chem[i]->m;
92 }
93 else {
94 out->chem[i]->m=0.0;
95 O2->m+=rx[i]->O2_flow()*in->chem[i]->n();
96 N2->m+=rx[i]->N2_flow()*in->chem[i]->n();
97 CO2->m+=rx[i]->CO2_flow()*in->chem[i]->n();
98 H2O->m+=rx[i]->H2O_flow()*in->chem[i]->n();
99 }
100 }
101 N2->m*=(1.0+eta);
102 O2->m*=(1.0+eta);
103 //perform energy balance to find Tout
104 T = in->T;
105
106 step=10;
107 Q=1;
108 // find temperature
109 while (fabs(step)>TOL_BURN && fabs(Q)>TOL_BURN && T<MAX_TEMP)
110 {
111 T+=step;
112
113 if(T>MAX_TEMP)
114 T=MAX_TEMP;
115 Q = 0.0;
116 for ( i = 0 ; i < in->nb ; i++ )
117 Q += in->chem[i]->dH ( in->T , T , in->P ) * in->chem[i]->n();
118
119 for ( i = 0 ; i < in->nb ; i++ )
120 if ( can_burn[i] )
121 Q += rx[i]->Hcomb(T) * in->chem[i]->n();
122
123 Q += O2->dH ( 293 , T , in->P ) * O2->n();
124 Q += N2->dH ( 293 , T , in->P ) * N2->n();
125
126
127 if (step/fabs(step)*Q >0)
128 step*= -0.1;
129 else if (fabs(Q)<10)
130 step*=0.25;
131
132 }
133
134
135
136 out->set_thermo(in->thermo);
137 // out->thermo = in->thermo;
138
139 out->set(in->P, T);
140 O2->P=in->P; O2->T=T; O2->state=1; O2->find_v();
141 N2->P=in->P; N2->T=T; N2->state=1; N2->find_v();
142 CO2->P=in->P; CO2->T=T; CO2->state=1; CO2->find_v();
143 H2O->P=in->P; H2O->T=T; H2O->state=1; H2O->find_v();
144 //check if mixture can burn
145 m_can_burn = 0.0;
146 for(i=0;i<in->nb;i++) if(can_burn[i]) m_can_burn+=in->chem[i]->n();
147 LFLmix=0.0;
148 for(i=0;i<in->nb;i++) if(can_burn[i]) LFLmix+=in->chem[i]->n()/m_can_burn*rx[i]->LFL(in->P,T);
149 UFLmix=0.0;
150 for(i=0;i<in->nb;i++) if(can_burn[i]) UFLmix+=in->chem[i]->n()/m_can_burn*rx[i]->UFL(in->P,T);
151 num = 0.0;
152 buff=in->T; in->set(in->P, T);
153 for(i=0;i<in->nb;i++) if(can_burn[i]) num+=in->chem[i]->n()/in->n()*in->v;
154 in->set(in->P, buff);
155 den = O2->v+N2->v+out->v;
156 composition = num/den;
157 if(!(LFLmix<=composition && composition<=UFLmix) || T==MAX_TEMP)
158 {
159 // logf.open(MESSAGES,ios::app);
160 // logf<<" --> Warning <-- Mixture in "<<name<<" can't burn (LFL="<<LFLmix
161 // <<" UFL="<<UFLmix<<" x="<<composition<<").\n";
162 // logf.close();
163 T=in->T;
164 filename = out->name;
165 *out=*in;
166 out->set(filename);
167 // out->write(); // WRITE TOTO
168 OK=false;
169 }
170 else
171 {
172 O2->P=in->P; O2->T=T; O2->state=1; O2->find_v();
173 N2->P=in->P; N2->T=T; N2->state=1; N2->find_v();
174 // out->write(); // WRITE TOTO
175 }
176 if(OK) //compute the pollutants production
177 {
178 fill_K_array();
179 NO = 1e6*sqrt(K[0]*(N2->n()/den)*(O2->n()/den))*den*0.03/(O2->m+N2->m+out->m+H2O->m+CO2->m);
180 N2O = 1e6*K[1]*(N2->n()/den)*sqrt(O2->n()/den)*den*0.044/(O2->m+N2->m+out->m+H2O->m+CO2->m);
181 NO2 = 1e6*K[2]*sqrt(N2->n()/den)*(O2->n()/den)*den*0.046/(O2->m+N2->m+out->m+H2O->m+CO2->m);
182 CO = 1e6*K[3]*(CO2->n()/den)*den/sqrt(O2->n()/den)*0.028/(O2->m+N2->m+out->m+H2O->m+CO2->m);
183 }
184 // logf.open(MESSAGES,ios::app);
185 if (NO>EPS && NO2>EPS && N2O>EPS) {
186 // logf<<" --> Warning <-- Presence of NOx: "<<(NO+NO2+N2O)<<" ppm in "<<name<<".\n";
187 y[7] = NO+NO2+N2O;
188
189 if (ARRONDI)
190 y[7] = arrondi ( y[7] , 6 );
191
192 }
193 if (CO>EPS) {
194 // logf<<" --> Warning <-- Presence of CO: "<<CO<<" ppm in "<<name<<".\n";
195 y[8] = CO;
196 if (ARRONDI)
197 y[8] = arrondi ( y[8] , 6 );
198
199 }
200 // logf.close();
201 return OK;
202 }
203
fill_K_array()204 void burner::fill_K_array()
205 {
206 a[0]=1.0; a[1]=1.0; a[2]=0.5; a[3]=1.0;
207 b[0]=1.0; b[1]=0.5; b[2]=1.0; b[3]=-0.5;
208 c[0]=2.0; c[1]=1.0; c[2]=1.0; c[3]=1.0;
209 K[0] = exp(-120.27*(173.38-T*0.012)/T);
210 K[1] = exp(-120.27*(103.64+T*0.074)/T);
211 K[2] = exp(-120.27*(51.96+T*0.061)/T);
212 K[3] = exp(-120.27*(283.84-T*0.087)/T);
213 for(i=0;i<4;i++)
214 K[i]*=pow(1000, c[i]-a[i]-b[i]);
215 }
216
write()217 void burner::write() {
218
219 cout << setprecision(6);
220
221 cout << "WRITE FILE " << RUNTIME << name << ".unit" << " :\n\tBEGIN\n";
222 cout << "\t>> " << name;
223 cout << endl << "\t>> stream in : " << in->name;
224 cout << endl << "\t>> streams out : " << out->name;
225 cout << endl << "\t>> P = " << in->P << " atm, T(in) = " << in->T
226 << " T(out) = " << T << " K";
227 O2->P = 1;
228 O2->T = 293;
229 O2->state = 1;
230 O2->find_v();
231 N2->P=1;
232 N2->T=293;
233 N2->state=1;
234 N2->find_v();
235 cout << endl << "\t>> Required air flow = "
236 << (O2->m+N2->m) << " kg/s (" << (O2->v+N2->v) << " m3/s)";
237 O2->P=in->P;
238 O2->T=T;
239 O2->state=1;
240 O2->find_v();
241 N2->P=in->P;
242 N2->T=T;
243 N2->state=1;
244 N2->find_v();
245 step=(eta*O2->v/(1+eta)+N2->v+H2O->v+CO2->v+out->v);
246 cout << endl << "\t>> Total flue gases = "
247 << (out->m+CO2->m+H2O->m+N2->m+eta*O2->m/(1+eta))
248 <<" kg/s (" << step << " m3/s)";
249 cout << "\n\tEND\n\n";
250 cost();
251 }
252
253
get_cost(void)254 double burner::get_cost ( void ) {
255
256
257 O2->P = 1;
258 O2->T = 293;
259 O2->state = 1;
260 O2->find_v();
261 N2->P=1;
262 N2->T=293;
263 N2->state=1;
264 N2->find_v();
265 O2->P=in->P;
266 O2->T=T;
267 O2->state=1;
268 O2->find_v();
269 N2->P=in->P;
270 N2->T=T;
271 N2->state=1;
272 N2->find_v();
273 step=(eta*O2->v/(1+eta)+N2->v+H2O->v+CO2->v+out->v);
274
275 buff = 3.1761-0.1373*log10(step) + 0.3414*pow(log10(step),2);
276 buff = 2.7*pow(10, buff);
277 buff = buff*MS_YEAR/MS_2001;
278
279 return buff;
280 }
281
282
cost(void)283 void burner::cost ( void ) {
284 cout << "WRITE FILE " << RUNTIME << name << ".cost" << " :\n\tBEGIN\n";
285 cout << "\t>>" << get_cost();
286 cout << "\n\tEND\n\n";
287 }
288