RadLib
Loading...
Searching...
No Matches
rad_rcslw.cc
Go to the documentation of this file.
1
6#include <algorithm>
7#include <cmath>
8#include <fstream>
9#include <iostream>
10#include <cstdlib> // exit
11#include "multilinear_interpolation.h"
12#include "rad_rcslw.h"
13
14using namespace std;
15
16#define GETSTRINGEDVALUE(themacro) STRINGIFY(themacro)
17#define STRINGIFY(thestring) #thestring
18
31
32rad_rcslw::rad_rcslw(const int p_nGG,
33 const double TbTref,
34 const double p_P,
35 const double fvsoot,
36 const double xH2O,
37 const double xCO2,
38 const double xCO) :
39 rad(p_nGG, p_nGG + 1){
40
41 P = p_P/101325; // atm
42 Tref = TbTref;
43 Tb = TbTref;
44
45 Cmin = 0.0001;
46 Cmax = 1000.0;
47
48 P_table = vector<double>{0.1, 0.25, 0.5, 1,2,4,8,15,30,50};
49 C_table = vector<double>(71);
50 for(int i=0; i<71; i++)
51 C_table[i] = Cmin*pow(Cmax/Cmin, i/70.0);
52
53 Tg_table = vector<double>(28);
54 for(int i=0; i<28; i++)
55 Tg_table[i] = 300 + i*(3000.0-300.0)/27;
56
57 Tb_table = vector<double>(28);
58 for(int i=0; i<28; i++)
59 Tb_table[i] = 300 + i*(3000.0-300.0)/27;
60
61 xH2O_table = vector<double>{0, 0.05, 0.1, 0.2, 0.3, 0.4, 0.6, 0.8, 1.0};
62
63 nP = P_table.size(); // 10
64 nC = C_table.size(); // 71
65 nTg = Tg_table.size(); // 28
66 nTb = Tb_table.size(); // 28
67 ny_H2O = xH2O_table.size(); // 9
68
69 set_Falbdf_CO2_CO_H2O_at_P();
70
71 Fmin = get_F_albdf(Cmin, Tref, Tref, xCO2, xCO, xH2O, fvsoot);
72 Fmax = get_F_albdf(Cmax, Tref, Tref, xCO2, xCO, xH2O, fvsoot);
73
74 set_Fpts();
75
76}
77
98
99void rad_rcslw::get_k_a_oneband(double &kabs,
100 double &awts,
101 const int iband,
102 const double T,
103 const double P_not_used,
104 const double fvsoot,
105 const double xH2O,
106 const double xCO2,
107 const double xCO,
108 const double xCH4_not_used){
109
110 if(iband < 0 || iband >= nGGa) {
111 cerr << "\n\n***** rad_rcslw::get_k_a_oneband: iband out of range *****\n" << endl;
112 exit(0);
113 }
114
115 //--------------- kabs
116
117 if(iband==0)
118 kabs = 0.0;
119 else{
120 double Nconc = P*101325/8.31446/T; // mol/m3
121 double C = get_FI_albdf(F_pts[iband-1], T, Tb, xCO2, xCO, xH2O, fvsoot);
122 kabs = Nconc * C;
123 }
124
125 //--------------- awts
126
127 if(iband==0){
128 double Ct = get_FI_albdf(Ft_pts[0], T, Tb, xCO2, xCO, xH2O, fvsoot);
129 awts = get_F_albdf(Ct, T, T, xCO2, xCO, xH2O, fvsoot);
130 }
131 else{
132 double Ct_i = get_FI_albdf(Ft_pts[iband], T, Tb, xCO2, xCO, xH2O, fvsoot);
133 double Ct_im1 = get_FI_albdf(Ft_pts[iband-1], T, Tb, xCO2, xCO, xH2O, fvsoot);
134 awts = get_F_albdf(Ct_i, T, T, xCO2, xCO, xH2O, fvsoot) -
135 get_F_albdf(Ct_im1, T, T, xCO2, xCO, xH2O, fvsoot);
136 }
137
138}
139
159
160void rad_rcslw::get_k_a(vector<double> &kabs,
161 vector<double> &awts,
162 const double T,
163 const double P_not_used,
164 const double fvsoot,
165 const double xH2O,
166 const double xCO2,
167 const double xCO,
168 const double xCH4_not_used){
169
170 vector<double> C(nGG);
171 vector<double> Ct(nGGa);
172 vector<double> FCt(nGGa);
173
174 for(int j=0; j < nGG; j++)
175 C[j] = get_FI_albdf(F_pts[j], T, Tb, xCO2, xCO, xH2O, fvsoot);
176
177 for(int j=0; j < nGGa; j++)
178 Ct[j] = get_FI_albdf(Ft_pts[j], T, Tb, xCO2, xCO, xH2O, fvsoot);
179
180 double Nconc = P*101325/8.31446/T; // mol/m3
181
182 kabs.resize(nGGa);
183 kabs[0] = 0.0;
184 for(int j=1; j < nGGa; j++)
185 kabs[j] = Nconc * C[j-1];
186
187 for(int j=0; j < nGGa; j++)
188 FCt[j] = get_F_albdf(Ct[j], T, T, xCO2, xCO, xH2O, fvsoot);
189
190 awts.resize(nGGa);
191 awts[0] = FCt[0];
192 for(int j=1; j < nGGa; j++)
193 awts[j] = FCt[j] - FCt[j-1];
194
195 //------------------ contents below could be used instead for this whole function, but slower
196
197 //kabs.resize(nGGa);
198 //awts.resize(nGGa);
199
200 //double k, a;
201 //for(int i=0; i<nGGa; i++){
202 // get_k_a_oneband(k, a, i, T, P, fvsoot, xH2O, xCO2, xCO, xCH4_not_used);
203 // kabs[i] = k;
204 // awts[i] = a;
205 //}
206}
207
221
222double rad_rcslw::get_F_albdf(const double C, double Tg, double Tb,
223 double xCO2, double xCO, double xH2O,
224 double fvsoot){
225
226
227 if (xCO2 <= 1E-20) xCO2 = 1E-20;
228 if (xCO <= 1e-20) xCO = 1e-20;
229 if (xH2O <= 1e-20) xH2O = 1e-20;
230 if (fvsoot < 0.0) fvsoot = 0.0;
231
232 if (Tg < Tg_table[0]) Tg = Tg_table[0];
233 if (Tg > Tg_table.back()) Tg = Tg_table.back();
234
235 if (Tb < Tb_table[0]) Tb = Tb_table[0];
236 if (Tb > Tb_table.back()) Tb = Tb_table.back();
237
238 if (xH2O < xH2O_table[0]) xH2O = xH2O_table[0];
239 if (xH2O > xH2O_table.back()) xH2O = xH2O_table.back();
240
241 double CxCO2 = C/xCO2;
242 double CxCO = C/xCO;
243 double CxH2O = C/xH2O;
244
245 if (CxCO2 < C_table[0]) CxCO2 = C_table[0];
246 if (CxCO2 > C_table.back()) CxCO2 = C_table.back();
247 if (CxCO < C_table[0]) CxCO = C_table[0];
248 if (CxCO > C_table.back()) CxCO = C_table.back();
249 if (CxH2O < C_table[0]) CxH2O = C_table[0];
250 if (CxH2O > C_table.back()) CxH2O = C_table.back();
251
252 double F_CO2, F_CO, F_H2O;
253
254 F_CO2 = LI_3D(nTg, nTb, nC, &Tg_table[0], &Tb_table[0], &C_table[0], &Falbdf_CO2[0], Tg, Tb, CxCO2);
255 F_CO = LI_3D(nTg, nTb, nC, &Tg_table[0], &Tb_table[0], &C_table[0], &Falbdf_CO[0], Tg, Tb, CxCO);
256 F_H2O = LI_4D(ny_H2O, nTg, nTb, nC, &xH2O_table[0], &Tg_table[0], &Tb_table[0], &C_table[0], &Falbdf_H2O[0], xH2O, Tg, Tb, CxH2O);
257
258 return F_CO2 * F_CO * F_H2O * F_albdf_soot(C, Tg, Tb, fvsoot);
259}
260
261
275
276double rad_rcslw::get_FI_albdf(const double F, const double Tg, const double Tb,
277 const double xCO2, const double xCO, const double xH2O,
278 const double fvsoot){
279
280 if (F <= get_F_albdf(Cmin, Tg, Tb, xCO2, xCO, xH2O, fvsoot))
281 return(Cmin);
282 if (F >= get_F_albdf(Cmax, Tg, Tb, xCO2, xCO, xH2O, fvsoot))
283 return(Cmax);
284
285 // Find table location
286
287 int iLo = 0;
288 int iHi = nC - 1;
289
290 int maxit = 100;
291 int it;
292 for (it=0; it < maxit; it++){
293 if (iHi-iLo == 1)
294 break;
295 int iMi = (iLo+iHi)/2;
296 if (F > get_F_albdf(C_table[iMi], Tg, Tb, xCO2, xCO, xH2O, fvsoot))
297 iLo = iMi;
298 else
299 iHi = iMi;
300 }
301
302 int iMi = (iLo+iHi)/2;
303 if(it == maxit){
304 cout << "WARNING, NO CONVERGENCE IN" << maxit << "iterations" << "\n";
305 return C_table[iMi] ;
306 }
307
308 // Now interpolate to the solution
309
310 double Flo = get_F_albdf(C_table[iLo], Tg, Tb, xCO2, xCO, xH2O, fvsoot);
311 double Fhi = get_F_albdf(C_table[iHi], Tg, Tb, xCO2, xCO, xH2O, fvsoot);
312 double cc = C_table[iLo] + (F-Flo) * (C_table[iHi]-C_table[iLo])/(Fhi-Flo);
313
314 // cout << "relative error:" << (F - get_F_albdf(cc, Tg, Tb, xCO2, xCO, xH2O, fvsoot))/F << "\n";
315 // return cc; //doldb
316
317 //----------- but F(cc) is not equal to F! So give it another interpolation:
318 // rel error ~ 0.1% --> 1E-6
319
320 double FF = get_F_albdf(cc, Tg, Tb, xCO2, xCO, xH2O, fvsoot);
321 double ccc = C_table[iLo] + (F-Flo)*(cc-C_table[iLo])/(FF-Flo);
322
323 return ccc;
324
327
328 //double FFF = get_F_albdf(ccc, Tg, Tb, xCO2, xCO, xH2O, fvsoot);
329 //double cccc = cc + (F-FF)*(ccc-cc)/(FFF-FF);
330
332
333 //return cccc;
334
335}
336
342
343void rad_rcslw::set_Fpts(){
344
345 if(nGG > 25){
346 cout << endl << "ERROR: nGG too high for quadrature: nGG max is currently 25" << endl;
347 exit(0);
348 }
349
350 vector<vector<double> > W(26);
351
352 W[0] = vector<double> {0.0,0.0};
353 W[1] = vector<double> {1.0000000000000000E+00, 1.0000000000000000E+00};
354 W[2] = vector<double> {3.4785484513745368E-01, 6.5214515486254621E-01, 6.5214515486254621E-01, 3.4785484513745368E-01};
355 W[3] = vector<double> {1.7132449237916975E-01, 3.6076157304813894E-01, 4.6791393457269137E-01, 4.6791393457269137E-01, 3.6076157304813894E-01, 1.7132449237916975E-01};
356 W[4] = vector<double> {1.0122853629037669E-01, 2.2238103445337434E-01, 3.1370664587788705E-01, 3.6268378337836177E-01, 3.6268378337836177E-01, 3.1370664587788705E-01, 2.2238103445337434E-01, 1.0122853629037669E-01};
357 W[5] = vector<double> {6.6671344308688069E-02, 1.4945134915058036E-01, 2.1908636251598201E-01, 2.6926671930999652E-01, 2.9552422471475298E-01, 2.9552422471475298E-01, 2.6926671930999652E-01, 2.1908636251598201E-01, 1.4945134915058036E-01, 6.6671344308688069E-02};
358 W[6] = vector<double> {4.7175336386512022E-02, 1.0693932599531888E-01, 1.6007832854334611E-01, 2.0316742672306565E-01, 2.3349253653835464E-01, 2.4914704581340269E-01, 2.4914704581340269E-01, 2.3349253653835464E-01, 2.0316742672306565E-01, 1.6007832854334611E-01, 1.0693932599531888E-01, 4.7175336386512022E-02};
359 W[7] = vector<double> {3.5119460331752374E-02, 8.0158087159760305E-02, 1.2151857068790296E-01, 1.5720316715819341E-01, 1.8553839747793763E-01, 2.0519846372129555E-01, 2.1526385346315766E-01, 2.1526385346315766E-01, 2.0519846372129555E-01, 1.8553839747793763E-01, 1.5720316715819341E-01, 1.2151857068790296E-01, 8.0158087159760305E-02, 3.5119460331752374E-02};
360 W[8] = vector<double> {2.7152459411754037E-02, 6.2253523938647706E-02, 9.5158511682492591E-02, 1.2462897125553403E-01, 1.4959598881657676E-01, 1.6915651939500262E-01, 1.8260341504492361E-01, 1.8945061045506859E-01, 1.8945061045506859E-01, 1.8260341504492361E-01, 1.6915651939500262E-01, 1.4959598881657676E-01, 1.2462897125553403E-01, 9.5158511682492591E-02, 6.2253523938647706E-02, 2.7152459411754037E-02};
361 W[9] = vector<double> {2.1616013526484130E-02, 4.9714548894969221E-02, 7.6425730254889246E-02, 1.0094204410628699E-01, 1.2255520671147836E-01, 1.4064291467065063E-01, 1.5468467512626521E-01, 1.6427648374583273E-01, 1.6914238296314363E-01, 1.6914238296314363E-01, 1.6427648374583273E-01, 1.5468467512626521E-01, 1.4064291467065063E-01, 1.2255520671147836E-01, 1.0094204410628699E-01, 7.6425730254889246E-02, 4.9714548894969221E-02, 2.1616013526484130E-02};
362 W[10] = vector<double> {1.7614007139153273E-02, 4.0601429800386217E-02, 6.2672048334109443E-02, 8.3276741576704671E-02, 1.0193011981724026E-01, 1.1819453196151825E-01, 1.3168863844917653E-01, 1.4209610931838187E-01, 1.4917298647260366E-01, 1.5275338713072578E-01, 1.5275338713072578E-01, 1.4917298647260366E-01, 1.4209610931838187E-01, 1.3168863844917653E-01, 1.1819453196151825E-01, 1.0193011981724026E-01, 8.3276741576704671E-02, 6.2672048334109443E-02, 4.0601429800386217E-02, 1.7614007139153273E-02};
363 W[11] = vector<double> {1.4627995298274705E-02, 3.3774901584815178E-02, 5.2293335152682870E-02, 6.9796468424520197E-02, 8.5941606217067396E-02, 1.0041414444288072E-01, 1.1293229608053883E-01, 1.2325237681051199E-01, 1.3117350478706188E-01, 1.3654149834601478E-01, 1.3925187285563156E-01, 1.3925187285563156E-01, 1.3654149834601478E-01, 1.3117350478706188E-01, 1.2325237681051199E-01, 1.1293229608053883E-01, 1.0041414444288072E-01, 8.5941606217067396E-02, 6.9796468424520197E-02, 5.2293335152682870E-02, 3.3774901584815178E-02, 1.4627995298274705E-02};
364 W[12] = vector<double> {1.2341229799987091E-02, 2.8531388628933743E-02, 4.4277438817419551E-02, 5.9298584915436742E-02, 7.3346481411080411E-02, 8.6190161531953288E-02, 9.7618652104114065E-02, 1.0744427011596561E-01, 1.1550566805372561E-01, 1.2167047292780342E-01, 1.2583745634682830E-01, 1.2793819534675221E-01, 1.2793819534675221E-01, 1.2583745634682830E-01, 1.2167047292780342E-01, 1.1550566805372561E-01, 1.0744427011596561E-01, 9.7618652104114065E-02, 8.6190161531953288E-02, 7.3346481411080411E-02, 5.9298584915436742E-02, 4.4277438817419551E-02, 2.8531388628933743E-02, 1.2341229799987091E-02};
365 W[13] = vector<double> {1.0551372617343395E-02, 2.4417851092631938E-02, 3.7962383294363120E-02, 5.0975825297148079E-02, 6.3274046329574674E-02, 7.4684149765659763E-02, 8.5045894313485068E-02, 9.4213800355914160E-02, 1.0205916109442532E-01, 1.0847184052857647E-01, 1.1336181654631956E-01, 1.1666044348529646E-01, 1.1832141527926213E-01, 1.1832141527926213E-01, 1.1666044348529646E-01, 1.1336181654631956E-01, 1.0847184052857647E-01, 1.0205916109442532E-01, 9.4213800355914160E-02, 8.5045894313485068E-02, 7.4684149765659763E-02, 6.3274046329574674E-02, 5.0975825297148079E-02, 3.7962383294363120E-02, 2.4417851092631938E-02, 1.0551372617343395E-02};
366 W[14] = vector<double> {9.1242825930943974E-03, 2.1132112592771271E-02, 3.2901427782304517E-02, 4.4272934759003985E-02, 5.5107345675716936E-02, 6.5272923966999755E-02, 7.4646214234568811E-02, 8.3113417228900935E-02, 9.0571744393032852E-02, 9.6930657997929923E-02, 1.0211296757806078E-01, 1.0605576592284637E-01, 1.0871119225829413E-01, 1.1004701301647524E-01, 1.1004701301647524E-01, 1.0871119225829413E-01, 1.0605576592284637E-01, 1.0211296757806078E-01, 9.6930657997929923E-02, 9.0571744393032852E-02, 8.3113417228900935E-02, 7.4646214234568811E-02, 6.5272923966999755E-02, 5.5107345675716936E-02, 4.4272934759003985E-02, 3.2901427782304517E-02, 2.1132112592771271E-02, 9.1242825930943974E-03};
367 W[15] = vector<double> {7.9681924961695228E-03, 1.8466468311091087E-02, 2.8784707883322873E-02, 3.8799192569626793E-02, 4.8402672830594434E-02, 5.7493156217619093E-02, 6.5974229882180324E-02, 7.3755974737704802E-02, 8.0755895229419811E-02, 8.6899787201082698E-02, 9.2122522237785789E-02, 9.6368737174643990E-02, 9.9593420586794934E-02, 1.0176238974840521E-01, 1.0285265289355848E-01, 1.0285265289355848E-01, 1.0176238974840521E-01, 9.9593420586794934E-02, 9.6368737174643990E-02, 9.2122522237785789E-02, 8.6899787201082698E-02, 8.0755895229419811E-02, 7.3755974737704802E-02, 6.5974229882180324E-02, 5.7493156217619093E-02, 4.8402672830594434E-02, 3.8799192569626793E-02, 2.8784707883322873E-02, 1.8466468311091087E-02, 7.9681924961695228E-03};
368 W[16] = vector<double> {7.0186100094692984E-03, 1.6274394730905965E-02, 2.5392065309262427E-02, 3.4273862913021626E-02, 4.2835898022226426E-02, 5.0998059262376244E-02, 5.8684093478535704E-02, 6.5822222776361752E-02, 7.2345794108848449E-02, 7.8193895787070311E-02, 8.3311924226946846E-02, 8.7652093004403908E-02, 9.1173878695763863E-02, 9.3844399080804566E-02, 9.5638720079274833E-02, 9.6540088514727812E-02, 9.6540088514727812E-02, 9.5638720079274833E-02, 9.3844399080804566E-02, 9.1173878695763863E-02, 8.7652093004403908E-02, 8.3311924226946846E-02, 7.8193895787070311E-02, 7.2345794108848449E-02, 6.5822222776361752E-02, 5.8684093478535704E-02, 5.0998059262376244E-02, 4.2835898022226426E-02, 3.4273862913021626E-02, 2.5392065309262427E-02, 1.6274394730905965E-02, 7.0186100094692984E-03};
369 W[17] = vector<double> {6.2291405559100326E-03, 1.4450162748594548E-02, 2.2563721985495038E-02, 3.0491380638445777E-02, 3.8166593796386906E-02, 4.5525611523353778E-02, 5.2507414572678060E-02, 5.9054135827524654E-02, 6.5111521554076457E-02, 7.0629375814255727E-02, 7.5561974660031797E-02, 7.9868444339771819E-02, 8.3513099699845522E-02, 8.6465739747035669E-02, 8.8701897835693738E-02, 9.0203044370640612E-02, 9.0956740330259786E-02, 9.0956740330259786E-02, 9.0203044370640612E-02, 8.8701897835693738E-02, 8.6465739747035669E-02, 8.3513099699845522E-02, 7.9868444339771819E-02, 7.5561974660031797E-02, 7.0629375814255727E-02, 6.5111521554076457E-02, 5.9054135827524654E-02, 5.2507414572678060E-02, 4.5525611523353778E-02, 3.8166593796386906E-02, 3.0491380638445777E-02, 2.2563721985495038E-02, 1.4450162748594548E-02, 6.2291405559100326E-03};
370 W[18] = vector<double> {5.5657196642477837E-03, 1.2915947284064104E-02, 2.0181515297735174E-02, 2.7298621498568355E-02, 3.4213810770307482E-02, 4.0875750923645232E-02, 4.7235083490266047E-02, 5.3244713977759678E-02, 5.8860144245324549E-02, 6.4039797355015429E-02, 6.8745323835736297E-02, 7.2941885005653004E-02, 7.6598410645870627E-02, 7.9687828912071559E-02, 8.2187266704339651E-02, 8.4078218979661792E-02, 8.5346685739338499E-02, 8.5983275670394627E-02, 8.5983275670394627E-02, 8.5346685739338499E-02, 8.4078218979661792E-02, 8.2187266704339651E-02, 7.9687828912071559E-02, 7.6598410645870627E-02, 7.2941885005653004E-02, 6.8745323835736297E-02, 6.4039797355015429E-02, 5.8860144245324549E-02, 5.3244713977759678E-02, 4.7235083490266047E-02, 4.0875750923645232E-02, 3.4213810770307482E-02, 2.7298621498568355E-02, 2.0181515297735174E-02, 1.2915947284064104E-02, 5.5657196642477837E-03};
371 W[19] = vector<double> {5.0028807496389329E-03, 1.1613444716468455E-02, 1.8156577709613410E-02, 2.4579739738232007E-02, 3.0839500545175719E-02, 3.6894081594024859E-02, 4.2703158504674703E-02, 4.8228061860758530E-02, 5.3432019910332196E-02, 5.8280399146997154E-02, 6.2740933392133172E-02, 6.6783937979140381E-02, 7.0382507066898942E-02, 7.3512692584743411E-02, 7.6153663548446479E-02, 7.8287844658210981E-02, 7.9901033243527819E-02, 8.0982493770597061E-02, 8.1525029280385797E-02, 8.1525029280385797E-02, 8.0982493770597061E-02, 7.9901033243527819E-02, 7.8287844658210981E-02, 7.6153663548446479E-02, 7.3512692584743411E-02, 7.0382507066898942E-02, 6.6783937979140381E-02, 6.2740933392133172E-02, 5.8280399146997154E-02, 5.3432019910332196E-02, 4.8228061860758530E-02, 4.2703158504674703E-02, 3.6894081594024859E-02, 3.0839500545175719E-02, 2.4579739738232007E-02, 1.8156577709613410E-02, 1.1613444716468455E-02, 5.0028807496389329E-03};
372 W[20] = vector<double> {4.5212770985300181E-03, 1.0498284531151609E-02, 1.6421058381907345E-02, 2.2245849194166653E-02, 2.7937006980023528E-02, 3.3460195282547678E-02, 3.8782167974472377E-02, 4.3870908185673324E-02, 4.8695807635072405E-02, 5.3227846983937115E-02, 5.7439769099391892E-02, 6.1306242492929319E-02, 6.4804013456601486E-02, 6.7912045815234398E-02, 7.0611647391287169E-02, 7.2886582395804478E-02, 7.4723169057968677E-02, 7.6110361900626741E-02, 7.7039818164248389E-02, 7.7505947978425332E-02, 7.7505947978425332E-02, 7.7039818164248389E-02, 7.6110361900626741E-02, 7.4723169057968677E-02, 7.2886582395804478E-02, 7.0611647391287169E-02, 6.7912045815234398E-02, 6.4804013456601486E-02, 6.1306242492929319E-02, 5.7439769099391892E-02, 5.3227846983937115E-02, 4.8695807635072405E-02, 4.3870908185673324E-02, 3.8782167974472377E-02, 3.3460195282547678E-02, 2.7937006980023528E-02, 2.2245849194166653E-02, 1.6421058381907345E-02, 1.0498284531151609E-02, 4.5212770985300181E-03};
373
374 W[21] = vector<double> {4.105998604646913e-03, 9.536220301748407e-03, 1.4922443697357294e-02, 2.0227869569052214e-02, 2.5422959526113627e-02, 3.0479240699603335e-02, 3.5369071097592e-02, 4.0065735180692286e-02, 4.454357777196598e-02, 4.877814079280347e-02, 5.274629569917415e-02, 5.642636935801852e-02, 5.979826222758681e-02, 6.284355804500275e-02, 6.554562436490913e-02, 6.788970337652217e-02, 6.986299249259438e-02, 7.14547142651712e-02, 7.265617524380438e-02, 7.346081345346782e-02, 7.386423423217317e-02, 7.386423423217317e-02, 7.346081345346782e-02, 7.265617524380438e-02, 7.14547142651712e-02, 6.986299249259438e-02, 6.788970337652217e-02, 6.554562436490913e-02, 6.284355804500275e-02, 5.979826222758681e-02, 5.642636935801852e-02, 5.274629569917415e-02, 4.877814079280347e-02, 4.454357777196598e-02, 4.0065735180692286e-02, 3.5369071097592e-02, 3.0479240699603335e-02, 2.5422959526113627e-02, 2.0227869569052214e-02, 1.4922443697357294e-02, 9.536220301748407e-03, 4.105998604646913e-03};
375 W[22] = vector<double> {3.7454048031164908e-03, 8.700481367523676e-03, 1.3619586755579383e-02, 1.8471481736814985e-02, 2.323148190201923e-02, 2.787578282128097e-02, 3.238122281206993e-02, 3.672534781380876e-02, 4.0886512310346054e-02, 4.484398408197005e-02, 4.857804644835181e-02, 5.207009609170433e-02, 5.53027355637279e-02, 5.825985987759541e-02, 6.092673670156192e-02, 6.329007973320361e-02, 6.533811487918127e-02, 6.706063890629348e-02, 6.844907026936646e-02, 6.949649186157239e-02, 7.019768547355804e-02, 7.054915778935386e-02, 7.054915778935386e-02, 7.019768547355804e-02, 6.949649186157239e-02, 6.844907026936646e-02, 6.706063890629348e-02, 6.533811487918127e-02, 6.329007973320361e-02, 6.092673670156192e-02, 5.825985987759541e-02, 5.53027355637279e-02, 5.207009609170433e-02, 4.857804644835181e-02, 4.484398408197005e-02, 4.0886512310346054e-02, 3.672534781380876e-02, 3.238122281206993e-02, 2.787578282128097e-02, 2.323148190201923e-02, 1.8471481736814985e-02, 1.3619586755579383e-02, 8.700481367523676e-03, 3.7454048031164908e-03};
376 W[23] = vector<double> {3.4303008681091456e-03, 7.969898229724492e-03, 1.2479883770988784e-02, 1.6933514007836395e-02, 2.1309998754136483e-02, 2.5589286397129797e-02, 2.975182955220257e-02, 3.37786279991068e-02, 3.765130535738592e-02, 4.135219010967875e-02, 4.486439527731805e-02, 4.817189510171219e-02, 5.1259598007142824e-02, 5.411341538585653e-02, 5.6720325843991094e-02, 5.906843459554615e-02, 6.1147027724650325e-02, 6.294662106439443e-02, 6.445900346713893e-02, 6.567727426778112e-02, 6.659587476845476e-02, 6.721061360067808e-02, 6.751868584903632e-02, 6.751868584903632e-02, 6.721061360067808e-02, 6.659587476845476e-02, 6.567727426778112e-02, 6.445900346713893e-02, 6.294662106439443e-02, 6.1147027724650325e-02, 5.906843459554615e-02, 5.6720325843991094e-02, 5.411341538585653e-02, 5.1259598007142824e-02, 4.817189510171219e-02, 4.486439527731805e-02, 4.135219010967875e-02, 3.765130535738592e-02, 3.37786279991068e-02, 2.975182955220257e-02, 2.5589286397129797e-02, 2.1309998754136483e-02, 1.6933514007836395e-02, 1.2479883770988784e-02, 7.969898229724492e-03, 3.4303008681091456e-03};
377 W[24] = vector<double> {3.1533460523091796e-03, 7.327553901276492e-03, 1.1477234579234974e-02, 1.5579315722942928e-02, 1.9616160457355297e-02, 2.3570760839324092e-02, 2.7426509708356882e-02, 3.116722783279834e-02, 3.477722256477066e-02, 3.8241351065830674e-02, 4.1545082943464554e-02, 4.46745608566941e-02, 4.7616658492490284e-02, 5.035903555385428e-02, 5.289018948519349e-02, 5.5199503699984054e-02, 5.727729210040293e-02, 5.911483969839548e-02, 6.070443916589358e-02, 6.2039423159892464e-02, 6.311419228625378e-02, 6.392423858464795e-02, 6.446616443594984e-02, 6.473769681268368e-02, 6.473769681268368e-02, 6.446616443594984e-02, 6.392423858464795e-02, 6.311419228625378e-02, 6.2039423159892464e-02, 6.070443916589358e-02, 5.911483969839548e-02, 5.727729210040293e-02, 5.5199503699984054e-02, 5.289018948519349e-02, 5.035903555385428e-02, 4.7616658492490284e-02, 4.46745608566941e-02, 4.1545082943464554e-02, 3.8241351065830674e-02, 3.477722256477066e-02, 3.116722783279834e-02, 2.7426509708356882e-02, 2.3570760839324092e-02, 1.9616160457355297e-02, 1.5579315722942928e-02, 1.1477234579234974e-02, 7.327553901276492e-03, 3.1533460523091796e-03};
378 W[25] = vector<double> {2.9086225531578685e-03, 6.759799195745505e-03, 1.059054838365167e-02, 1.438082276148571e-02, 1.811556071348957e-02, 2.1780243170124582e-02, 2.536067357001259e-02, 2.8842993580534923e-02, 3.2213728223577896e-02, 3.5459835615146026e-02, 3.8568756612587844e-02, 4.15284630901474e-02, 4.43275043388031e-02, 4.695505130394827e-02, 4.9400938449466136e-02, 5.165570306958092e-02, 5.371062188899592e-02, 5.555774480621228e-02, 5.718992564772815e-02, 5.8600849813222215e-02, 5.9785058704265225e-02, 6.073797084176991e-02, 6.145589959031641e-02, 6.193606742068298e-02, 6.2176616655347e-02, 6.2176616655347e-02, 6.193606742068298e-02, 6.145589959031641e-02, 6.073797084176991e-02, 5.9785058704265225e-02, 5.8600849813222215e-02, 5.718992564772815e-02, 5.555774480621228e-02, 5.371062188899592e-02, 5.165570306958092e-02, 4.9400938449466136e-02, 4.695505130394827e-02, 4.43275043388031e-02, 4.15284630901474e-02, 3.8568756612587844e-02, 3.5459835615146026e-02, 3.2213728223577896e-02, 2.8842993580534923e-02, 2.536067357001259e-02, 2.1780243170124582e-02, 1.811556071348957e-02, 1.438082276148571e-02, 1.059054838365167e-02, 6.759799195745505e-03, 2.9086225531578685e-03};
379
380 vector<vector<double> > X(26);
381
382 X[0] = vector<double> {0.0, 0.0};
383 X[1] = vector<double> {-5.7735026918962573E-01, 5.7735026918962573E-01};
384 X[2] = vector<double> {-8.6113631159405257E-01, -3.3998104358485626E-01, 3.3998104358485626E-01, 8.6113631159405257E-01};
385 X[3] = vector<double> {-9.3246951420315205E-01, -6.6120938646626448E-01, -2.3861918608319693E-01, 2.3861918608319693E-01, 6.6120938646626448E-01, 9.3246951420315205E-01};
386 X[4] = vector<double> {-9.6028985649753618E-01, -7.9666647741362673E-01, -5.2553240991632899E-01, -1.8343464249564978E-01, 1.8343464249564978E-01, 5.2553240991632899E-01, 7.9666647741362673E-01, 9.6028985649753618E-01};
387 X[5] = vector<double> {-9.7390652851717174E-01, -8.6506336668898454E-01, -6.7940956829902444E-01, -4.3339539412924721E-01, -1.4887433898163122E-01, 1.4887433898163122E-01, 4.3339539412924721E-01, 6.7940956829902444E-01, 8.6506336668898454E-01, 9.7390652851717174E-01};
388 X[6] = vector<double> {-9.8156063424671924E-01, -9.0411725637047480E-01, -7.6990267419430469E-01, -5.8731795428661748E-01, -3.6783149899818018E-01, -1.2523340851146891E-01, 1.2523340851146891E-01, 3.6783149899818018E-01, 5.8731795428661748E-01, 7.6990267419430469E-01, 9.0411725637047480E-01, 9.8156063424671924E-01};
389 X[7] = vector<double> {-9.8628380869681231E-01, -9.2843488366357352E-01, -8.2720131506976502E-01, -6.8729290481168548E-01, -5.1524863635815410E-01, -3.1911236892788974E-01, -1.0805494870734367E-01, 1.0805494870734367E-01, 3.1911236892788974E-01, 5.1524863635815410E-01, 6.8729290481168548E-01, 8.2720131506976502E-01, 9.2843488366357352E-01, 9.8628380869681231E-01};
390 X[8] = vector<double> {-9.8940093499164994E-01, -9.4457502307323260E-01, -8.6563120238783176E-01, -7.5540440835500300E-01, -6.1787624440264377E-01, -4.5801677765722737E-01, -2.8160355077925892E-01, -9.5012509837637454E-02, 9.5012509837637454E-02, 2.8160355077925892E-01, 4.5801677765722737E-01, 6.1787624440264377E-01, 7.5540440835500300E-01, 8.6563120238783176E-01, 9.4457502307323260E-01, 9.8940093499164994E-01};
391 X[9] = vector<double> {-9.9156516842093090E-01, -9.5582394957139782E-01, -8.9260246649755570E-01, -8.0370495897252314E-01, -6.9168704306035322E-01, -5.5977083107394754E-01, -4.1175116146284263E-01, -2.5188622569150548E-01, -8.4775013041735292E-02, 8.4775013041735292E-02, 2.5188622569150548E-01, 4.1175116146284263E-01, 5.5977083107394754E-01, 6.9168704306035322E-01, 8.0370495897252314E-01, 8.9260246649755570E-01, 9.5582394957139782E-01, 9.9156516842093090E-01};
392 X[10] = vector<double> {-9.9312859918509488E-01, -9.6397192727791381E-01, -9.1223442825132584E-01, -8.3911697182221878E-01, -7.4633190646015080E-01, -6.3605368072651502E-01, -5.1086700195082713E-01, -3.7370608871541955E-01, -2.2778585114164510E-01, -7.6526521133497338E-02, 7.6526521133497338E-02, 2.2778585114164510E-01, 3.7370608871541955E-01, 5.1086700195082713E-01, 6.3605368072651502E-01, 7.4633190646015080E-01, 8.3911697182221878E-01, 9.1223442825132584E-01, 9.6397192727791381E-01, 9.9312859918509488E-01};
393 X[11] = vector<double> {-9.9429458548239924E-01, -9.7006049783542869E-01, -9.2695677218717398E-01, -8.6581257772030018E-01, -7.8781680597920811E-01, -6.9448726318668275E-01, -5.8764040350691160E-01, -4.6935583798675706E-01, -3.4193582089208424E-01, -2.0786042668822130E-01, -6.9739273319722211E-02, 6.9739273319722211E-02, 2.0786042668822130E-01, 3.4193582089208424E-01, 4.6935583798675706E-01, 5.8764040350691160E-01, 6.9448726318668275E-01, 7.8781680597920811E-01, 8.6581257772030018E-01, 9.2695677218717398E-01, 9.7006049783542869E-01, 9.9429458548239924E-01};
394 X[12] = vector<double> {-9.9518721999702131E-01, -9.7472855597130947E-01, -9.3827455200273280E-01, -8.8641552700440096E-01, -8.2000198597390295E-01, -7.4012419157855436E-01, -6.4809365193697555E-01, -5.4542147138883956E-01, -4.3379350762604513E-01, -3.1504267969616340E-01, -1.9111886747361631E-01, -6.4056892862605630E-02, 6.4056892862605630E-02, 1.9111886747361631E-01, 3.1504267969616340E-01, 4.3379350762604513E-01, 5.4542147138883956E-01, 6.4809365193697555E-01, 7.4012419157855436E-01, 8.2000198597390295E-01, 8.8641552700440096E-01, 9.3827455200273280E-01, 9.7472855597130947E-01, 9.9518721999702131E-01};
395 X[13] = vector<double> {-9.9588570114561692E-01, -9.7838544595647092E-01, -9.4715906666171423E-01, -9.0263786198430707E-01, -8.4544594278849805E-01, -7.7638594882067880E-01, -6.9642726041995728E-01, -6.0669229301761807E-01, -5.0844071482450570E-01, -4.0305175512348629E-01, -2.9200483948595690E-01, -1.7685882035689018E-01, -5.9230093429313208E-02, 5.9230093429313208E-02, 1.7685882035689018E-01, 2.9200483948595690E-01, 4.0305175512348629E-01, 5.0844071482450570E-01, 6.0669229301761807E-01, 6.9642726041995728E-01, 7.7638594882067880E-01, 8.4544594278849805E-01, 9.0263786198430707E-01, 9.4715906666171423E-01, 9.7838544595647092E-01, 9.9588570114561692E-01};
396 X[14] = vector<double> {-9.9644249757395442E-01, -9.8130316537087281E-01, -9.5425928062893817E-01, -9.1563302639213207E-01, -8.6589252257439497E-01, -8.0564137091717913E-01, -7.3561087801363179E-01, -6.5665109403886501E-01, -5.6972047181140173E-01, -4.7587422495511827E-01, -3.7625151608907870E-01, -2.7206162763517810E-01, -1.6456928213338079E-01, -5.5079289884034266E-02, 5.5079289884034266E-02, 1.6456928213338079E-01, 2.7206162763517810E-01, 3.7625151608907870E-01, 4.7587422495511827E-01, 5.6972047181140173E-01, 6.5665109403886501E-01, 7.3561087801363179E-01, 8.0564137091717913E-01, 8.6589252257439497E-01, 9.1563302639213207E-01, 9.5425928062893817E-01, 9.8130316537087281E-01, 9.9644249757395442E-01};
397 X[15] = vector<double> {-9.9689348407464951E-01, -9.8366812327974729E-01, -9.6002186496830755E-01, -9.2620004742927431E-01, -8.8256053579205263E-01, -8.2956576238276836E-01, -7.6777743210482619E-01, -6.9785049479331585E-01, -6.2052618298924289E-01, -5.3662414814201986E-01, -4.4703376953808915E-01, -3.5270472553087812E-01, -2.5463692616788985E-01, -1.5386991360858354E-01, -5.1471842555317698E-02, 5.1471842555317698E-02, 1.5386991360858354E-01, 2.5463692616788985E-01, 3.5270472553087812E-01, 4.4703376953808915E-01, 5.3662414814201986E-01, 6.2052618298924289E-01, 6.9785049479331585E-01, 7.6777743210482619E-01, 8.2956576238276836E-01, 8.8256053579205263E-01, 9.2620004742927431E-01, 9.6002186496830755E-01, 9.8366812327974729E-01, 9.9689348407464951E-01};
398 X[16] = vector<double> {-9.9726386184948157E-01, -9.8561151154526838E-01, -9.6476225558750639E-01, -9.3490607593773967E-01, -8.9632115576605220E-01, -8.4936761373256997E-01, -7.9448379596794239E-01, -7.3218211874028971E-01, -6.6304426693021523E-01, -5.8771575724076230E-01, -5.0689990893222936E-01, -4.2135127613063533E-01, -3.3186860228212767E-01, -2.3928736225213706E-01, -1.4447196158279649E-01, -4.8307665687738310E-02, 4.8307665687738310E-02, 1.4447196158279649E-01, 2.3928736225213706E-01, 3.3186860228212767E-01, 4.2135127613063533E-01, 5.0689990893222936E-01, 5.8771575724076230E-01, 6.6304426693021523E-01, 7.3218211874028971E-01, 7.9448379596794239E-01, 8.4936761373256997E-01, 8.9632115576605220E-01, 9.3490607593773967E-01, 9.6476225558750639E-01, 9.8561151154526838E-01, 9.9726386184948157E-01};
399 X[17] = vector<double> {-9.9757175379084195E-01, -9.8722781640630952E-01, -9.6870826253334430E-01, -9.4216239740510710E-01, -9.0780967771832455E-01, -8.6593463833456441E-01, -8.1688422790093362E-01, -7.6106487662987299E-01, -6.9893911321626290E-01, -6.3102172708052851E-01, -5.5787550066974667E-01, -4.8010654519032703E-01, -3.9835927775864594E-01, -3.1331108133946323E-01, -2.2566669161644948E-01, -1.3615235725918298E-01, -4.5509821953102533E-02, 4.5509821953102533E-02, 1.3615235725918298E-01, 2.2566669161644948E-01, 3.1331108133946323E-01, 3.9835927775864594E-01, 4.8010654519032703E-01, 5.5787550066974667E-01, 6.3102172708052851E-01, 6.9893911321626290E-01, 7.6106487662987299E-01, 8.1688422790093362E-01, 8.6593463833456441E-01, 9.0780967771832455E-01, 9.4216239740510710E-01, 9.6870826253334430E-01, 9.8722781640630952E-01, 9.9757175379084195E-01};
400 X[18] = vector<double> {-9.9783046248408580E-01, -9.8858647890221230E-01, -9.7202769104969799E-01, -9.4827298439950758E-01, -9.1749777451565906E-01, -8.7992980089039707E-01, -8.3584716699247530E-01, -7.8557623013220657E-01, -7.2948917159355664E-01, -6.6800123658552102E-01, -6.0156765813598057E-01, -5.3068028592624517E-01, -4.5586394443342027E-01, -3.7767254711968923E-01, -2.9668499534402826E-01, -2.1350089231686559E-01, -1.2873610380938480E-01, -4.3018198473708608E-02, 4.3018198473708608E-02, 1.2873610380938480E-01, 2.1350089231686559E-01, 2.9668499534402826E-01, 3.7767254711968923E-01, 4.5586394443342027E-01, 5.3068028592624517E-01, 6.0156765813598057E-01, 6.6800123658552102E-01, 7.2948917159355664E-01, 7.8557623013220657E-01, 8.3584716699247530E-01, 8.7992980089039707E-01, 9.1749777451565906E-01, 9.4827298439950758E-01, 9.7202769104969799E-01, 9.8858647890221230E-01, 9.9783046248408580E-01};
401 X[19] = vector<double> {-9.9804993053568758E-01, -9.8973945426638554E-01, -9.7484632859015352E-01, -9.5346633093352962E-01, -9.2574133204858433E-01, -8.9185573900463222E-01, -8.5203502193236214E-01, -8.0654416760531689E-01, -7.5568590375397071E-01, -6.9979868037918436E-01, -6.3925441582968168E-01, -5.7445602104780713E-01, -5.0583471792793111E-01, -4.3384716943237650E-01, -3.5897244047943500E-01, -2.8170880979016527E-01, -2.0257045389211670E-01, -1.2208402533786741E-01, -4.0785147904578239E-02, 4.0785147904578239E-02, 1.2208402533786741E-01, 2.0257045389211670E-01, 2.8170880979016527E-01, 3.5897244047943500E-01, 4.3384716943237650E-01, 5.0583471792793111E-01, 5.7445602104780713E-01, 6.3925441582968168E-01, 6.9979868037918436E-01, 7.5568590375397071E-01, 8.0654416760531689E-01, 8.5203502193236214E-01, 8.9185573900463222E-01, 9.2574133204858433E-01, 9.5346633093352962E-01, 9.7484632859015352E-01, 9.8973945426638554E-01, 9.9804993053568758E-01};
402 X[20] = vector<double> {-9.9823770971055925E-01, -9.9072623869945708E-01, -9.7725994998377430E-01, -9.5791681921379168E-01, -9.3281280827867652E-01, -9.0209880696887434E-01, -8.6595950321225956E-01, -8.2461223083331170E-01, -7.7830565142651942E-01, -7.2731825518992710E-01, -6.7195668461417957E-01, -6.1255388966798030E-01, -5.4946712509512818E-01, -4.8307580168617870E-01, -4.1377920437160498E-01, -3.4199409082575849E-01, -2.6815218500725369E-01, -1.9269758070137111E-01, -1.1608407067525521E-01, -3.8772417506050816E-02, 3.8772417506050816E-02, 1.1608407067525521E-01, 1.9269758070137111E-01, 2.6815218500725369E-01, 3.4199409082575849E-01, 4.1377920437160498E-01, 4.8307580168617870E-01, 5.4946712509512818E-01, 6.1255388966798030E-01, 6.7195668461417957E-01, 7.2731825518992710E-01, 7.7830565142651942E-01, 8.2461223083331170E-01, 8.6595950321225956E-01, 9.0209880696887434E-01, 9.3281280827867652E-01, 9.5791681921379168E-01, 9.7725994998377430E-01, 9.9072623869945708E-01, 9.9823770971055925E-01};
403
404 X[21] = vector<double> {-9.983996189900625e-01, -9.915772883408609e-01, -9.793425080637482e-01, -9.617593653382045e-01, -9.389235573549881e-01, -9.109597249041275e-01, -8.780205698121728e-01, -8.402859832618169e-01, -7.979620532554874e-01, -7.512799356894805e-01, -7.004945905561712e-01, -6.458833888692479e-01, -5.877445974851093e-01, -5.263957499311923e-01, -4.6217191207042196e-01, -3.9542385204297503e-01, -3.265161244654115e-01, -2.5582507934287907e-01, -1.8373680656485455e-01, -1.1064502720851986e-01, -3.694894316535177e-02, 3.694894316535177e-02, 1.1064502720851986e-01, 1.8373680656485455e-01, 2.5582507934287907e-01, 3.265161244654115e-01, 3.9542385204297503e-01, 4.6217191207042196e-01, 5.263957499311923e-01, 5.877445974851093e-01, 6.458833888692479e-01, 7.004945905561712e-01, 7.512799356894805e-01, 7.979620532554874e-01, 8.402859832618169e-01, 8.780205698121728e-01, 9.109597249041275e-01, 9.389235573549881e-01, 9.617593653382045e-01, 9.793425080637482e-01, 9.915772883408609e-01, 9.983996189900625e-01};
405 X[22] = vector<double> {-9.985402006367742e-01, -9.923163921385159e-01, -9.81151833077914e-01, -9.650996504224931e-01, -9.442395091181941e-01, -9.186752599841758e-01, -8.885342382860432e-01, -8.539665950047104e-01, -8.15144539645135e-01, -7.722614792487559e-01, -7.25531053660717e-01, -6.751860706661224e-01, -6.214773459035758e-01, -5.646724531854708e-01, -5.050543913882023e-01, -4.429201745254115e-01, -3.7857935201470716e-01, -3.123524665027858e-01, -2.4456945692820126e-01, -1.7556801477551678e-01, -1.0569190170865325e-01, -3.5289236964135356e-02, 3.5289236964135356e-02, 1.0569190170865325e-01, 1.7556801477551678e-01, 2.4456945692820126e-01, 3.123524665027858e-01, 3.7857935201470716e-01, 4.429201745254115e-01, 5.050543913882023e-01, 5.646724531854708e-01, 6.214773459035758e-01, 6.751860706661224e-01, 7.25531053660717e-01, 7.722614792487559e-01, 8.15144539645135e-01, 8.539665950047104e-01, 8.885342382860432e-01, 9.186752599841758e-01, 9.442395091181941e-01, 9.650996504224931e-01, 9.81151833077914e-01, 9.923163921385159e-01, 9.985402006367742e-01};
406 X[23] = vector<double> {-9.98663042133818e-01, -9.929623489061743e-01, -9.827336698041669e-01, -9.680213918539919e-01, -9.488923634460897e-01, -9.25433798806754e-01, -8.97752711533942e-01, -8.65975394866858e-01, -8.302468370660661e-01, -7.907300570752742e-01, -7.47605359615666e-01, -7.010695120204057e-01, -6.513348462019977e-01, -5.986282897127152e-01, -5.431903302618026e-01, -4.852739183881647e-01, -4.251433132828284e-01, -3.630728770209957e-01, -2.9934582270187005e-01, -2.3425292220626975e-01, -1.6809117946710353e-01, -1.0116247530558424e-01, -3.377219001605204e-02, 3.377219001605204e-02, 1.0116247530558424e-01, 1.6809117946710353e-01, 2.3425292220626975e-01, 2.9934582270187005e-01, 3.630728770209957e-01, 4.251433132828284e-01, 4.852739183881647e-01, 5.431903302618026e-01, 5.986282897127152e-01, 6.513348462019977e-01, 7.010695120204057e-01, 7.47605359615666e-01, 7.907300570752742e-01, 8.302468370660661e-01, 8.65975394866858e-01, 8.97752711533942e-01, 9.25433798806754e-01, 9.488923634460897e-01, 9.680213918539919e-01, 9.827336698041669e-01, 9.929623489061743e-01, 9.98663042133818e-01};
407 X[24] = vector<double> {-9.987710072524261e-01, -9.935301722663508e-01, -9.841245837228269e-01, -9.705915925462473e-01, -9.529877031604309e-01, -9.313866907065543e-01, -9.058791367155696e-01, -8.765720202742479e-01, -8.435882616243935e-01, -8.070662040294426e-01, -7.671590325157404e-01, -7.240341309238146e-01, -6.778723796326639e-01, -6.288673967765136e-01, -5.772247260839727e-01, -5.23160974722233e-01, -4.669029047509584e-01, -4.086864819907167e-01, -3.4875588629216075e-01, -2.8736248735545555e-01, -2.2476379039468905e-01, -1.612223560688917e-01, -9.70046992094627e-02, -3.238017096286937e-02, 3.238017096286937e-02, 9.70046992094627e-02, 1.612223560688917e-01, 2.2476379039468905e-01, 2.8736248735545555e-01, 3.4875588629216075e-01, 4.086864819907167e-01, 4.669029047509584e-01, 5.23160974722233e-01, 5.772247260839727e-01, 6.288673967765136e-01, 6.778723796326639e-01, 7.240341309238146e-01, 7.671590325157404e-01, 8.070662040294426e-01, 8.435882616243935e-01, 8.765720202742479e-01, 9.058791367155696e-01, 9.313866907065543e-01, 9.529877031604309e-01, 9.705915925462473e-01, 9.841245837228269e-01, 9.935301722663508e-01, 9.987710072524261e-01};
408 X[25] = vector<double> {-9.98866404420071e-01, -9.940319694320907e-01, -9.853540840480058e-01, -9.72864385106692e-01, -9.566109552428079e-01, -9.36656618944878e-01, -9.130785566557919e-01, -8.859679795236131e-01, -8.554297694299461e-01, -8.21582070859336e-01, -7.845558329003992e-01, -7.444943022260686e-01, -7.015524687068222e-01, -6.558964656854394e-01, -6.077029271849502e-01, -5.5715830451465e-01, -5.044581449074642e-01, -4.4980633497403877e-01, -3.9341431189756515e-01, -3.3550024541943735e-01, -2.7628819377953195e-01, -2.1600723687604176e-01, -1.548905899981459e-01, -9.317470156008613e-02, -3.1098338327188876e-02, 3.1098338327188876e-02, 9.317470156008613e-02, 1.548905899981459e-01, 2.1600723687604176e-01, 2.7628819377953195e-01, 3.3550024541943735e-01, 3.9341431189756515e-01, 4.4980633497403877e-01, 5.044581449074642e-01, 5.5715830451465e-01, 6.077029271849502e-01, 6.558964656854394e-01, 7.015524687068222e-01, 7.444943022260686e-01, 7.845558329003992e-01, 8.21582070859336e-01, 8.554297694299461e-01, 8.859679795236131e-01, 9.130785566557919e-01, 9.36656618944878e-01, 9.566109552428079e-01, 9.72864385106692e-01, 9.853540840480058e-01, 9.940319694320907e-01, 9.98866404420071e-01};
409
410
411 vector<double> cumsum_w(nGG, W[nGG][nGG]);
412 for(int i=1; i<nGG; i++)
413 cumsum_w[i] = cumsum_w[i-1] + W[nGG][i+nGG];
414
415 Ft_pts.resize(nGGa); // \tilde{F} grid
416 Ft_pts[0] = Fmin;
417 for(int i=1; i<nGGa; i++)
418 Ft_pts[i] = Fmin + (Fmax-Fmin)*cumsum_w[i-1];
419
420 F_pts.resize(nGG);
421 for(auto i=nGG; i < X[nGG].size(); i++)
422 F_pts[i-nGG] = Fmin + X[nGG][i]*(Fmax - Fmin); // F grid (vals bet. \tild{F} pnts)
423
424}
425
434
435void rad_rcslw::set_Falbdf_CO2_CO_H2O_at_P(){
436
437 double P1;
438 double P2;
439 double f;
440
441 if(P < P_table[0] || P > P_table.back() ) {
442 cout << "Pressure = " << P << " atm is out of range\n";
443 exit(0);
444 }
445 else if(P==P_table[0]){
446 P1 = P_table[0];
447 P2 = P_table[1];
448 }
449 else if(P==P_table.back()){
450 P1 = P_table[nP-2];
451 P2 = P_table[nP-1];
452 }
453 else{
454 int i;
455 for(i=0; P_table[i] < P; i++)
456 P1 = P_table[i];
457 P2 = P_table[i];
458 }
459
460 f = (P-P1)/(P2-P1); // Interpolation factor
461
462 string Pres_1 = to_string(P1);
463 string Pres_2 = to_string(P2);
464 Pres_1.erase(Pres_1.find_last_not_of('0')+1, string::npos);
465 Pres_2.erase(Pres_2.find_last_not_of('0')+2, string::npos);
466 if(Pres_1.back() == '.') Pres_1.push_back('0');
467 if(Pres_2.back() == '.') Pres_2.push_back('0');
468 replace(Pres_1.begin(), Pres_1.end(), '.', '_');
469 replace(Pres_2.begin(), Pres_2.end(), '.', '_');
470
471 //--------------------- CO2
472
473 string CO2_file1 = string(GETSTRINGEDVALUE(RCSLW_DATA_DIR)) + "/co2_p" + Pres_1 + ".bin";
474 string CO2_file2 = string(GETSTRINGEDVALUE(RCSLW_DATA_DIR)) + "/co2_p" + Pres_2 + ".bin";
475 //string CO2_file1 = "ALBDF_Tables/co2_p" + Pres_1 + ".bin";
476 //string CO2_file2 = "ALBDF_Tables/co2_p" + Pres_2 + ".bin";
477
478 vector<double> CO2_F1(nTg*nTb*nC);
479 vector<double> CO2_F2(nTg*nTb*nC);
480 Falbdf_CO2.resize(nTg*nTb*nC);
481
482 get_FI_albdf_tables(CO2_file1, nTg, nTb, nC, CO2_F1);
483 get_FI_albdf_tables(CO2_file2, nTg, nTb, nC, CO2_F2);
484
485 for(auto ind=0; ind < Falbdf_CO2.size(); ++ind)
486 Falbdf_CO2[ind] = CO2_F1[ind]*(1-f) + CO2_F2[ind]*f;
487
488 //---------------------CO
489
490 string CO_file1 = string(GETSTRINGEDVALUE(RCSLW_DATA_DIR)) + "/co_p" + Pres_1 + ".bin";
491 string CO_file2 = string(GETSTRINGEDVALUE(RCSLW_DATA_DIR)) + "/co_p" + Pres_2 + ".bin";
492 //string CO_file1 = "ALBDF_Tables/co_p" + Pres_1 + ".bin";
493 //string CO_file2 = "ALBDF_Tables/co_p" + Pres_2 + ".bin";
494
495 vector<double> CO_F1(nTg*nTb*nC);
496 vector<double> CO_F2(nTg*nTb*nC);
497 Falbdf_CO.resize(nTg*nTb*nC);
498
499 get_FI_albdf_tables(CO_file1, nTg, nTb, nC, CO_F1);
500 get_FI_albdf_tables(CO_file2, nTg, nTb, nC, CO_F2);
501
502 for(auto ind=0; ind < Falbdf_CO.size(); ++ind)
503 Falbdf_CO[ind] = CO_F1[ind]*(1-f) + CO_F2[ind]*f;
504
505 //---------------------H2O
506
507 string H2O_file1 = string(GETSTRINGEDVALUE(RCSLW_DATA_DIR)) + "/h2o_p" + Pres_1 + ".bin";
508 string H2O_file2 = string(GETSTRINGEDVALUE(RCSLW_DATA_DIR)) + "/h2o_p" + Pres_2 + ".bin";
509 //string H2O_file1 = "ALBDF_Tables/h2o_p" + Pres_1 + ".bin";
510 //string H2O_file2 = "ALBDF_Tables/h2o_p" + Pres_2 + ".bin";
511
512 vector<double> H2O_F1(ny_H2O*nTg*nTb*nC);
513 vector<double> H2O_F2(ny_H2O*nTg*nTb*nC);
514 Falbdf_H2O.resize(ny_H2O*nTg*nTb*nC);
515
516 get_FI_albdf_tables(H2O_file1, ny_H2O, nTg, nTb, nC, H2O_F1);
517 get_FI_albdf_tables(H2O_file2, ny_H2O, nTg, nTb, nC, H2O_F2);
518
519 for(auto ind=0; ind < Falbdf_H2O.size(); ++ind)
520 Falbdf_H2O[ind] = H2O_F1[ind]*(1-f) + H2O_F2[ind]*f;
521
522}
523
534
535void rad_rcslw::get_FI_albdf_tables(const string Ptable_file_name,
536 const int nx, const int ny, const int nz,
537 vector<double> &myarray){
538
539 ifstream ifile;
540 ifile.open(Ptable_file_name, ios::in | ios::binary);
541 if(!ifile){
542 cout << endl << "error opening file: " << Ptable_file_name << endl;
543 exit(0);
544 }
545 float data;
546 int ntot = nx*ny*nz;
547 for(int ind=0; ind<ntot; ++ind){
548 ifile.read((char *) &data, sizeof(data));
549 myarray[ind] = static_cast<double>(data);
550 }
551 ifile.close();
552
553}
554
566
567void rad_rcslw::get_FI_albdf_tables(const string Ptable_file_name,
568 const int nx, const int ny, const int nz, const int nw,
569 vector<double> &myarray){
570
571 ifstream ifile;
572 ifile.open(Ptable_file_name, ios::in | ios::binary);
573 if(!ifile){
574 cout << endl << "error opening file: " << Ptable_file_name << endl;
575 exit(0);
576 }
577 float data;
578 int ntot = nx*ny*nz*nw;
579 for(int ind=0; ind<ntot; ++ind){
580 ifile.read((char *) &data, sizeof(data));
581 myarray[ind] = static_cast<double>(data);
582 }
583 ifile.close();
584
585}
587
618
619double rad_rcslw::F_albdf_soot(const double C, const double Tg, const double Tb, const double fvsoot){
620
621 if (fvsoot < 1E-12)
622 return 1.0;
623
624 double ksoot = 1.03; // real part of complex refractive index
625 double nsoot = 1.75; // imag part of complex refractive index
626 double csoot = 36*M_PI*nsoot*ksoot/(pow((nsoot*nsoot - ksoot*ksoot + 2),2) + 4*pow(nsoot*ksoot,2));
627
628 double hCokb = 0.01438777354; // h*Co/kb = C2 = PlancK*lightSpeed/Boltzmann = 6.62607004E-34 m2*kg/s * 299792458 m/s / 1.38064852E-23 m2*kg/s2*K
629
630 double Nconc = P*101325/8.31446/Tg; // mol/m3
631
632 double x = hCokb*C*Nconc / (csoot*fvsoot*Tb);
633 double sum = 0.0;
634 double nx;
635 for(int n=1; n<=3; n++){ // sum from n=1 to oo, but n=1 to 3 is enough
636 nx = n*x;
637 sum += exp(-nx)/pow(n,4)*(nx*(nx*(nx+3)+6)+6);
638 }
639 return 1.0 - 15.0/pow(M_PI,4)*sum;
640}
641
rad_rcslw::F_albdf_soot
double F_albdf_soot(const double C, const double Tg, const double Tb, const double fvsoot)
Definition rad_rcslw.cc:619
rad_rcslw::Tg_table
std::vector< double > Tg_table
gas temperatures (K, abscissas) in table
Definition rad_rcslw.h:33
rad_rcslw::rad_rcslw
rad_rcslw(const int p_nGG, const double TbTref, const double p_P, const double fvsoot, const double xH2O, const double xCO2, const double xCO)
Definition rad_rcslw.cc:32
rad_rcslw::C_table
std::vector< double > C_table
cross sections (m2/mol abscissas) in table
Definition rad_rcslw.h:32
rad_rcslw::P
double P
system pressure (atm); assumed constant
Definition rad_rcslw.h:27
rad_rcslw::Fmax
double Fmax
maximum albdf, corresponding to Cmax
Definition rad_rcslw.h:53
rad_rcslw::Tb_table
std::vector< double > Tb_table
blackbody temperature (K, abscissas) in table
Definition rad_rcslw.h:34
rad_rcslw::P_table
std::vector< double > P_table
pressures (atm, abscissas) in table
Definition rad_rcslw.h:31
rad_rcslw::xH2O_table
std::vector< double > xH2O_table
mole fractions H2O (abscissas) in table
Definition rad_rcslw.h:35
rad_rcslw::Falbdf_H2O
std::vector< double > Falbdf_H2O
H2O albdf table values.
Definition rad_rcslw.h:42
rad_rcslw::Cmax
double Cmax
maximum absorption cross section
Definition rad_rcslw.h:51
rad_rcslw::get_k_a_oneband
void get_k_a_oneband(double &kabs, double &awts, const int iband, const double T, const double P_not_used, const double fvsoot, const double xH2O, const double xCO2, const double xCO, const double xCH4_not_used)
Definition rad_rcslw.cc:99
rad_rcslw::set_Fpts
void set_Fpts(void)
Definition rad_rcslw.cc:343
rad_rcslw::get_F_albdf
double get_F_albdf(const double C, double Tg, double Tb, double xCO2, double xCO, double xH2O, double fvsoot)
Definition rad_rcslw.cc:222
rad_rcslw::ny_H2O
int ny_H2O
number of h2o values (abscissas) in table
Definition rad_rcslw.h:48
rad_rcslw::Ft_pts
std::vector< double > Ft_pts
albdf grid (t for tilde)
Definition rad_rcslw.h:38
rad_rcslw::get_k_a
void get_k_a(std::vector< double > &kabs, std::vector< double > &awts, const double T, const double P_not_used, const double fvsoot, const double xH2O, const double xCO2, const double xCO, const double xCH4_not_used)
Definition rad_rcslw.cc:160
rad_rcslw::nTg
int nTg
number of Tg values (abscissas) in table
Definition rad_rcslw.h:46
rad_rcslw::F_pts
std::vector< double > F_pts
albdf grid
Definition rad_rcslw.h:37
rad_rcslw::nTb
int nTb
number of Tb values (abscissas) in table
Definition rad_rcslw.h:47
rad_rcslw::set_Falbdf_CO2_CO_H2O_at_P
void set_Falbdf_CO2_CO_H2O_at_P()
Definition rad_rcslw.cc:435
rad_rcslw::Tref
double Tref
reference temperature (K) for setting F grid (normally same as Tb)
Definition rad_rcslw.h:28
rad_rcslw::get_FI_albdf_tables
void get_FI_albdf_tables(const std::string Ptable_file_name, const int nx, const int ny, const int nz, std::vector< double > &myarray)
rad_rcslw::Falbdf_CO2
std::vector< double > Falbdf_CO2
CO2 albdf table values (with abscissas above)
Definition rad_rcslw.h:40
rad_rcslw::Fmin
double Fmin
minimum albdf, corresponding to Cmin
Definition rad_rcslw.h:52
rad_rcslw::Tb
double Tb
black temperature (K)
Definition rad_rcslw.h:29
rad_rcslw::get_FI_albdf
double get_FI_albdf(const double F, const double Tg, const double Tb, const double xCO2, const double xCO, const double xH2O, const double fvsoot)
Definition rad_rcslw.cc:276
rad_rcslw::nC
int nC
number of C values (abscissas) in table
Definition rad_rcslw.h:45
rad_rcslw::Falbdf_CO
std::vector< double > Falbdf_CO
CO albdf table values.
Definition rad_rcslw.h:41
rad_rcslw::nP
int nP
number of P values (abscissas) in table
Definition rad_rcslw.h:44
rad_rcslw::Cmin
double Cmin
minimum absorption cross section
Definition rad_rcslw.h:50
rad
Definition rad.h:20
rad::nGGa
int nGGa
number of grey gases including the clear gas
Definition rad.h:32
rad::nGG
int nGG
number of gray gases, not including the clear gas
Definition rad.h:31
multilinear_interpolation.h
Multilinear interpolation (header only).
LI_4D
double LI_4D(const int nx, const int ny, const int nz, const int nw, const double *x, const double *y, const double *z, const double *w, const double *f, const double xP, const double yP, const double zP, const double wP)
Definition multilinear_interpolation.h:121
LI_3D
double LI_3D(const int nx, const int ny, const int nz, const double *x, const double *y, const double *z, const double *f, const double xP, const double yP, const double zP)
Definition multilinear_interpolation.h:80
GETSTRINGEDVALUE
#define GETSTRINGEDVALUE(themacro)
Definition rad_rcslw.cc:16
rad_rcslw.h
Header file for child class rad_rcslw.