antioch-0.4.0
cea_thermo_unit.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------bl-
2 //--------------------------------------------------------------------------
3 //
4 // Antioch - A Gas Dynamics Thermochemistry Library
5 //
6 // Copyright (C) 2014-2016 Paul T. Bauman, Benjamin S. Kirk,
7 // Sylvain Plessis, Roy H. Stonger
8 //
9 // Copyright (C) 2013 The PECOS Development Team
10 //
11 // This library is free software; you can redistribute it and/or
12 // modify it under the terms of the Version 2.1 GNU Lesser General
13 // Public License as published by the Free Software Foundation.
14 //
15 // This library is distributed in the hope that it will be useful,
16 // but WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 // Lesser General Public License for more details.
19 //
20 // You should have received a copy of the GNU Lesser General Public
21 // License along with this library; if not, write to the Free Software
22 // Foundation, Inc. 51 Franklin Street, Fifth Floor,
23 // Boston, MA 02110-1301 USA
24 //
25 //-----------------------------------------------------------------------el-
26 //
27 // $Id$
28 //
29 //--------------------------------------------------------------------------
30 //--------------------------------------------------------------------------
31 
32 // C++
33 #include <cmath>
34 #include <limits>
35 
36 // Antioch
39 #include "antioch/cea_thermo.h"
40 
41 template <typename Scalar>
42 int test_cp( const std::string& species_name, unsigned int species, Scalar cp_exact, Scalar T,
44 {
45  using std::abs;
46 
47  int return_flag = 0;
48 
49  const Scalar tol = std::numeric_limits<Scalar>::epsilon() * 5;
50 
51  typedef typename Antioch::CEAThermodynamics<Scalar>::template Cache<Scalar> Cache;
52 
53  const Scalar cp = thermo.cp(Cache(T), species);
54 
55  if( abs( (cp_exact - cp)/cp_exact ) > tol )
56  {
57  std::cerr << "Error: Mismatch in species specific heat."
58  << "\nspecies = " << species_name
59  << "\ncp = " << cp
60  << "\ncp_exact = " << cp_exact
61  << "\ndifference = " << (cp_exact - cp)
62  << "\ntolerance = " << tol
63  << "\nT = " << T << std::endl;
64  return_flag = 1;
65  }
66 
67  return return_flag;
68 }
69 
70 
71 template <typename Scalar>
72 Scalar cp( Scalar T, Scalar a0, Scalar a1, Scalar a2,
73  Scalar a3, Scalar a4, Scalar a5, Scalar a6 )
74 {
75  if( T < 200.1)
76  T = 200.1;
77 
78  return a0/(T*T) + a1/T + a2 + a3*T + a4*(T*T) + a5*(T*T*T) + a6*(T*T*T*T);
79 }
80 
81 
82 template <typename Scalar>
83 int tester()
84 {
85  std::vector<std::string> species_str_list;
86  const unsigned int n_species = 5;
87  species_str_list.reserve(n_species);
88  species_str_list.push_back( "N2" );
89  species_str_list.push_back( "O2" );
90  species_str_list.push_back( "N" );
91  species_str_list.push_back( "O" );
92  species_str_list.push_back( "NO" );
93 
94  const Scalar Mm_N = 14.008e-3L;
95  const Scalar Mm_O = 16.000e-3L;
96  const Scalar Mm_N2 = 2.L * Mm_N;
97  const Scalar Mm_O2 = 2.L * Mm_O;
98  const Scalar Mm_NO = Mm_N + Mm_O;
99 
100  Antioch::ChemicalMixture<Scalar> chem_mixture( species_str_list );
101 
102  Antioch::CEAThermodynamics<Scalar> thermo( chem_mixture );
103 
104  //const Scalar P = 100000.0;
105  const Scalar T1 = 190.0;
106  const Scalar T2 = 1500.0;
107  const Scalar T3 = 10000.0;
108 
109  const Scalar R_N2 = Antioch::Constants::R_universal<Scalar>()/Mm_N2;
110  const Scalar R_O2 = Antioch::Constants::R_universal<Scalar>()/Mm_O2;
111  const Scalar R_N = Antioch::Constants::R_universal<Scalar>()/Mm_N;
112  const Scalar R_O = Antioch::Constants::R_universal<Scalar>()/Mm_O;
113  const Scalar R_NO = Antioch::Constants::R_universal<Scalar>()/Mm_NO;
114 
115  int return_flag = 0;
116 
117  // Test N2 cp
118  {
119  unsigned int index = 0;
120  const Scalar cp_N2_1 = R_N2*cp( T1, Scalar(2.21037122e+04), Scalar(-3.81846145e+02), Scalar(6.08273815e+00),
121  Scalar(-8.53091381e-03), Scalar(1.38464610e-05), Scalar(-9.62579293e-09), Scalar(2.51970560e-12));
122 
123  const Scalar cp_N2_2 = R_N2*cp( T2, Scalar(5.87709908e+05), Scalar(-2.23924255e+03), Scalar(6.06694267e+00),
124  Scalar(-6.13965296e-04), Scalar(1.49179819e-07), Scalar(-1.92309442e-11), Scalar(1.06194871e-15) );
125 
126  const Scalar cp_N2_3 = R_N2*cp( T3, Scalar(8.30971200e+08), Scalar(-6.42048187e+05), Scalar(2.02020507e+02), Scalar(-3.06501961e-02),
127  Scalar(2.48685558e-06), Scalar(-9.70579208e-11), Scalar(1.43751673e-15));
128 
129  const Antioch::Species species = chem_mixture.species_list()[index];
130  const std::string species_name = chem_mixture.species_inverse_name_map().find(species)->second;
131 
132  int return_flag_temp = 0;
133  return_flag_temp = test_cp( species_name, index, cp_N2_1, T1, thermo );
134  if( return_flag_temp != 0 ) return_flag = 1;
135 
136  return_flag = test_cp( species_name, index, cp_N2_2, T2, thermo );
137  if( return_flag_temp != 0 ) return_flag = 1;
138 
139  return_flag = test_cp( species_name, index, cp_N2_3, T3, thermo );
140  if( return_flag_temp != 0 ) return_flag = 1;
141  }
142 
143  // Test O2 cp
144  {
145  unsigned int index = 1;
146  const Scalar cp_1 = R_O2*cp( T1, Scalar(-3.42556269e+04), Scalar(4.84699986e+02), Scalar(1.11901159e+00),
147  Scalar(4.29388743e-03), Scalar(-6.83627313e-07), Scalar(-2.02337478e-09) , Scalar(1.03904064e-12) );
148 
149  const Scalar cp_2 = R_O2*cp( T2, Scalar(-1.03793994e+06), Scalar(2.34483275e+03), Scalar(1.81972949e+00), Scalar(1.26784887e-03),
150  Scalar(-2.18807142e-07), Scalar(2.05372411e-11), Scalar(-8.19349062e-16) );
151 
152  const Scalar cp_3 = R_O2*cp( T3, Scalar(4.97515261e+08), Scalar(-2.86602339e+05), Scalar(6.69015464e+01), Scalar(-6.16971869e-03),
153  Scalar(3.01623757e-07), Scalar(-7.42087888e-12), Scalar(7.27744063e-17));
154 
155  const Antioch::Species species = chem_mixture.species_list()[index];
156  const std::string species_name = chem_mixture.species_inverse_name_map().find(species)->second;
157 
158  int return_flag_temp = 0;
159  return_flag_temp = test_cp( species_name, index, cp_1, T1, thermo );
160  if( return_flag_temp != 0 ) return_flag = 1;
161 
162  return_flag = test_cp( species_name, index, cp_2, T2, thermo );
163  if( return_flag_temp != 0 ) return_flag = 1;
164 
165  return_flag = test_cp( species_name, index, cp_3, T3, thermo );
166  if( return_flag_temp != 0 ) return_flag = 1;
167  }
168 
169  // Test N cp
170  {
171  unsigned int index = 2;
172  const Scalar cp_1 = R_N*cp( T1, Scalar(0.00000000e+00), Scalar(0.00000000e+00), Scalar(2.50000000e+00), Scalar(0.00000000e+00),
173  Scalar(0.00000000e+00), Scalar(0.00000000e+00), Scalar(0.00000000e+00));
174 
175  const Scalar cp_2 = R_N*cp( T2, Scalar(8.87650138e+04), Scalar(-1.07123150e+02), Scalar(2.36218829e+00), Scalar(2.91672008e-04),
176  Scalar(-1.72951510e-07), Scalar(4.01265788e-11), Scalar(-2.67722757e-15) );
177 
178  const Scalar cp_3 = R_N*cp( T3, Scalar(5.47518105e+08), Scalar(-3.10757498e+05), Scalar(6.91678274e+01), Scalar(-6.84798813e-03),
179  Scalar(3.82757240e-07), Scalar(-1.09836771e-11), Scalar(1.27798602e-16));
180 
181  const Antioch::Species species = chem_mixture.species_list()[index];
182  const std::string species_name = chem_mixture.species_inverse_name_map().find(species)->second;
183 
184  int return_flag_temp = 0;
185  return_flag_temp = test_cp( species_name, index, cp_1, T1, thermo );
186  if( return_flag_temp != 0 ) return_flag = 1;
187 
188  return_flag = test_cp( species_name, index, cp_2, T2, thermo );
189  if( return_flag_temp != 0 ) return_flag = 1;
190 
191  return_flag = test_cp( species_name, index, cp_3, T3, thermo );
192  if( return_flag_temp != 0 ) return_flag = 1;
193  }
194 
195 
196  // Test O cp
197  {
198  unsigned int index = 3;
199  const Scalar cp_1 = R_O*cp( T1, Scalar(-7.95361130e+03) , Scalar(1.60717779e+02), Scalar(1.96622644e+00), Scalar(1.01367031e-03),
200  Scalar(-1.11041542e-06), Scalar(6.51750750e-10), Scalar(-1.58477925e-13) );
201 
202  const Scalar cp_2 = R_O*cp( T2, Scalar(2.61902026e+05), Scalar(-7.29872203e+02), Scalar(3.31717727e+00), Scalar(-4.28133436e-04),
203  Scalar(1.03610459e-07), Scalar(-9.43830433e-12), Scalar(2.72503830e-16) );
204 
205  const Scalar cp_3 = R_O*cp( T3, Scalar(1.77900426e+08), Scalar(-1.08232826e+05), Scalar(2.81077837e+01), Scalar(-2.97523226e-03),
206  Scalar(1.85499753e-07), Scalar(-5.79623154e-12), Scalar(7.19172016e-17) );
207 
208  const Antioch::Species species = chem_mixture.species_list()[index];
209  const std::string species_name = chem_mixture.species_inverse_name_map().find(species)->second;
210 
211  int return_flag_temp = 0;
212  return_flag_temp = test_cp( species_name, index, cp_1, T1, thermo );
213  if( return_flag_temp != 0 ) return_flag = 1;
214 
215  return_flag = test_cp( species_name, index, cp_2, T2, thermo );
216  if( return_flag_temp != 0 ) return_flag = 1;
217 
218  return_flag = test_cp( species_name, index, cp_3, T3, thermo );
219  if( return_flag_temp != 0 ) return_flag = 1;
220  }
221 
222 
223  // Test NO cp
224  {
225  unsigned int index = 4;
226  const Scalar cp_1 = R_NO*cp( T1, Scalar(-1.14391658e+04), Scalar(1.53646774e+02), Scalar(3.43146865e+00), Scalar(-2.66859213e-03),
227  Scalar(8.48139877e-06), Scalar(-7.68511079e-09), Scalar(2.38679758e-12) );
228 
229  const Scalar cp_2 = R_NO*cp( T2, Scalar(2.23903708e+05), Scalar(-1.28965624e+03), Scalar(5.43394039e+00), Scalar(-3.65605546e-04),
230  Scalar(9.88101763e-08), Scalar(-1.41608327e-11), Scalar(9.38021642e-16) );
231 
232  const Scalar cp_3 = R_NO*cp( T3, Scalar(-9.57530764e+08), Scalar(5.91243671e+05), Scalar(-1.38456733e+02), Scalar(1.69433998e-02),
233  Scalar(-1.00735146e-06), Scalar(2.91258526e-11), Scalar(-3.29511091e-16) );
234 
235  const Antioch::Species species = chem_mixture.species_list()[index];
236  const std::string species_name = chem_mixture.species_inverse_name_map().find(species)->second;
237 
238  int return_flag_temp = 0;
239  return_flag_temp = test_cp( species_name, index, cp_1, T1, thermo );
240  if( return_flag_temp != 0 ) return_flag = 1;
241 
242  return_flag = test_cp( species_name, index, cp_2, T2, thermo );
243  if( return_flag_temp != 0 ) return_flag = 1;
244 
245  return_flag = test_cp( species_name, index, cp_3, T3, thermo );
246  if( return_flag_temp != 0 ) return_flag = 1;
247  }
248 
249  return return_flag;
250 }
251 
252 
253 int main()
254 {
255 // We're not getting the full long double precision yet?
256  return (tester<double>() ||
257 // tester<long double>() ||
258  tester<float>());
259 }
unsigned int Species
const std::vector< Species > & species_list() const
Scalar cp(Scalar T, Scalar a0, Scalar a1, Scalar a2, Scalar a3, Scalar a4, Scalar a5, Scalar a6)
int main()
int test_cp(const std::string &species_name, unsigned int species, Scalar cp_exact, Scalar T, const Antioch::CEAThermodynamics< Scalar > &thermo)
const std::map< Species, std::string > & species_inverse_name_map() const
Class storing chemical mixture properties.
int tester()
StateType cp(const Cache< StateType > &cache, unsigned int species) const
We currently need different specializations for scalar vs vector inputs here.
Definition: cea_thermo.h:264

Generated on Thu Jul 7 2016 11:09:46 for antioch-0.4.0 by  doxygen 1.8.8