antioch-0.4.0
stat_mech_thermo.h
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 
28 #ifndef ANTIOCH_STAT_MECH_THERMO_H
29 #define ANTIOCH_STAT_MECH_THERMO_H
30 
31 // Antioch
32 #include "antioch/chemical_mixture.h"
33 #include "antioch/input_utils.h"
34 #include "antioch/metaprogramming.h"
35 #include "antioch/antioch_exceptions.h"
36 
37 // C++
38 #include <iomanip>
39 #include <vector>
40 #include <cmath>
41 #include <limits>
42 
43 namespace Antioch
44 {
45 
46  template<typename CoeffType=double>
47  class StatMechThermodynamics
48  {
49  public:
50 
51  StatMechThermodynamics( const ChemicalMixture<CoeffType>& chem_mixture );
52 
54 
55  virtual ~StatMechThermodynamics();
56 
65  CoeffType cv_trans( const unsigned int species ) const;
66 
75  CoeffType cv_trans_over_R( const unsigned int species ) const;
76 
88  CoeffType cv_rot( const unsigned int species ) const;
89 
101  CoeffType cv_rot_over_R( const unsigned int species ) const;
102 
107  CoeffType cv_tr (const unsigned int species) const;
108 
113  CoeffType cv_tr_over_R (const unsigned int species) const;
114 
115 
120  template<typename VectorStateType>
121  typename enable_if_c<
122  has_size<VectorStateType>::value,
123  typename Antioch::value_type<VectorStateType>::type
124  >::type
125  cv_tr (const VectorStateType& mass_fractions) const;
126 
131  template<typename StateType>
132  StateType cv_vib (const unsigned int species, const StateType& Tv) const;
133 
138  template<typename StateType>
139  StateType cv_vib_over_R (const unsigned int species, const StateType& Tv) const;
140 
145  template<typename VectorStateType>
146  typename enable_if_c<
147  has_size<VectorStateType>::value,
148  typename Antioch::value_type<VectorStateType>::type
149  >::type
150  cv_vib (const typename Antioch::value_type<VectorStateType>::type& Tv,
151  const VectorStateType& mass_fractions) const;
152 
156  template<typename StateType>
157  StateType cv_el (const unsigned int species, const StateType& Te) const;
158 
162  template<typename VectorStateType>
163  typename enable_if_c<
164  has_size<VectorStateType>::value,
165  typename Antioch::value_type<VectorStateType>::type
166  >::type
167  cv_el (const typename Antioch::value_type<VectorStateType>::type& Te,
168  const VectorStateType& mass_fractions) const;
169 
174  template<typename StateType>
175  StateType cv_ve (const unsigned int species, const StateType& Tv) const;
176 
181  template<typename VectorStateType>
182  typename enable_if_c<
183  has_size<VectorStateType>::value,
184  typename Antioch::value_type<VectorStateType>::type
185  >::type
186  cv_ve (const typename Antioch::value_type<VectorStateType>::type& Tv,
187  const VectorStateType& mass_fractions) const;
188 
195  template<typename StateType>
196  StateType cv (const unsigned int species, const StateType& T, const StateType& Tv) const;
197 
204  template<typename VectorStateType>
205  typename enable_if_c<
206  has_size<VectorStateType>::value,
207  typename Antioch::value_type<VectorStateType>::type
208  >::type
209  cv (const typename Antioch::value_type<VectorStateType>::type& T,
210  const typename Antioch::value_type<VectorStateType>::type& Tv,
211  const VectorStateType& mass_fractions) const;
212 
216  template<typename VectorStateType>
217  typename enable_if_c<
218  has_size<VectorStateType>::value,
219  typename Antioch::value_type<VectorStateType>::type
220  >::type
221  cp (const typename Antioch::value_type<VectorStateType>::type& T,
222  const typename Antioch::value_type<VectorStateType>::type& Tv,
223  const VectorStateType& mass_fractions) const;
224 
228  template<typename StateType>
229  StateType h_tot (const unsigned int species, const StateType& T, const StateType& Tv) const;
230 
235  template<typename StateType>
236  StateType h_tot (const unsigned int species, const StateType& T) const;
237 
241  template<typename VectorStateType>
242  typename enable_if_c<
243  has_size<VectorStateType>::value,
244  typename Antioch::value_type<VectorStateType>::type
245  >::type
246  h_tot (const typename Antioch::value_type<VectorStateType>::type& T,
247  const typename Antioch::value_type<VectorStateType>::type& Tv,
248  const VectorStateType& mass_fractions) const;
249 
254  template<typename VectorStateType>
255  typename enable_if_c<
256  has_size<VectorStateType>::value,
257  typename Antioch::value_type<VectorStateType>::type
258  >::type
259  h_tot (const typename Antioch::value_type<VectorStateType>::type& T,
260  const VectorStateType& mass_fractions) const;
261 
265  template<typename StateType>
266  StateType e_tot (const unsigned int species, const StateType& T, const StateType& Tv) const;
267 
271  template<typename VectorStateType>
272  typename enable_if_c<
273  has_size<VectorStateType>::value,
274  typename Antioch::value_type<VectorStateType>::type
275  >::type
276  e_tot (const typename Antioch::value_type<VectorStateType>::type& T,
277  const typename Antioch::value_type<VectorStateType>::type& Tv,
278  const VectorStateType& mass_fractions) const;
279 
284  template<typename StateType>
285  StateType e_tot (const unsigned int species, const StateType& T) const;
286 
291  template<typename VectorStateType>
292  typename enable_if_c<
293  has_size<VectorStateType>::value,
294  typename Antioch::value_type<VectorStateType>::type
295  >::type
296  e_tot (const typename Antioch::value_type<VectorStateType>::type& T,
297  const VectorStateType& mass_fractions) const;
298 
302  template<typename StateType>
303  StateType e_tr (const unsigned int species, const StateType& T) const;
304 
308  template<typename VectorStateType>
309  typename enable_if_c<
310  has_size<VectorStateType>::value,
311  typename Antioch::value_type<VectorStateType>::type
312  >::type
313  e_tr (const typename Antioch::value_type<VectorStateType>::type& T,
314  const VectorStateType& mass_fractions) const;
315 
319  template<typename StateType>
320  StateType e_vib (const unsigned int species, const StateType& Tv) const;
321 
325  template<typename VectorStateType>
326  typename enable_if_c<
327  has_size<VectorStateType>::value,
328  typename Antioch::value_type<VectorStateType>::type
329  >::type
330  e_vib (const typename Antioch::value_type<VectorStateType>::type& Tv,
331  const VectorStateType& mass_fractions) const;
332 
336  template<typename StateType>
337  StateType e_el (const unsigned int species, const StateType& Te) const;
338 
342  template<typename VectorStateType>
343  typename enable_if_c<
344  has_size<VectorStateType>::value,
345  typename Antioch::value_type<VectorStateType>::type
346  >::type
347  e_el (const typename Antioch::value_type<VectorStateType>::type& Te,
348  const VectorStateType& mass_fractions) const;
349 
353  template<typename StateType>
354  StateType e_ve (const unsigned int species, const StateType& Tv) const;
355 
359  template<typename VectorStateType>
360  typename enable_if_c<
361  has_size<VectorStateType>::value,
362  typename Antioch::value_type<VectorStateType>::type
363  >::type
364  e_ve (const typename Antioch::value_type<VectorStateType>::type& Te,
365  const VectorStateType& mass_fractions) const;
366 
375  CoeffType e_ve_min () const;
376 
382  template<typename VectorStateType>
383  void e_and_cv_ve (const typename Antioch::value_type<VectorStateType>::type& Tv,
384  const VectorStateType& mass_fractions,
385  typename Antioch::value_type<VectorStateType>::type &e_ve,
386  typename Antioch::value_type<VectorStateType>::type &cv_ve) const;
387 
391  CoeffType e_0 (const unsigned int species) const;
392 
396  template<typename VectorStateType>
397  typename enable_if_c<
398  has_size<VectorStateType>::value,
399  typename Antioch::value_type<VectorStateType>::type
400  >::type
401  e_0 (const VectorStateType& mass_fractions) const;
402 
407  template<typename VectorStateType>
408  typename enable_if_c<
409  has_size<VectorStateType>::value,
410  typename Antioch::value_type<VectorStateType>::type
411  >::type
412  Tv_from_e_ve (const typename Antioch::value_type<VectorStateType>::type& e_ve,
413  const VectorStateType& mass_fractions,
414  typename Antioch::value_type<VectorStateType>::type Tv = -1) const;
415 
420  template<typename VectorStateType>
421  typename enable_if_c<
422  has_size<VectorStateType>::value,
423  typename Antioch::value_type<VectorStateType>::type
424  >::type
425  T_from_e_tot (const typename Antioch::value_type<VectorStateType>::type& e_tot,
426  const VectorStateType& mass_fractions,
427  typename Antioch::value_type<VectorStateType>::type T = -1) const;
428 
433  template<typename VectorStateType>
434  typename enable_if_c<
435  has_size<VectorStateType>::value,
436  typename Antioch::value_type<VectorStateType>::type
437  >::type
438  T_from_e_tr (const typename Antioch::value_type<VectorStateType>::type& e_tr,
439  const VectorStateType& mass_fractions,
440  typename Antioch::value_type<VectorStateType>::type T = -1) const;
441 
446  template<typename VectorStateType>
447  typename enable_if_c<
448  has_size<VectorStateType>::value,
449  typename Antioch::value_type<VectorStateType>::type
450  >::type
451  T_from_h_tot (const typename Antioch::value_type<VectorStateType>::type& h_tot,
452  const VectorStateType& mass_fractions,
453  typename Antioch::value_type<VectorStateType>::type T = -1) const;
454 
460  template<typename VectorStateType>
461  typename enable_if_c<
462  has_size<VectorStateType>::value,
463  typename Antioch::value_type<VectorStateType>::type
464  >::type
465  T_from_h_tot_Tv (const typename Antioch::value_type<VectorStateType>::type& h_tot,
466  const typename Antioch::value_type<VectorStateType>::type& Tv,
467  const VectorStateType& mass_fractions,
468  typename Antioch::value_type<VectorStateType>::type T = -1) const;
469 
474  template<typename VectorStateType>
475  typename enable_if_c<
476  has_size<VectorStateType>::value,
477  typename Antioch::value_type<VectorStateType>::type
478  >::type
479  s (const typename Antioch::value_type<VectorStateType>::type& T,
480  const typename Antioch::value_type<VectorStateType>::type& p,
481  const VectorStateType& mass_fractions) const;
482 
483  protected:
484 
485  const ChemicalMixture<CoeffType>& _chem_mixture;
486 
487  private:
488 
490 
491  StatMechThermodynamics();
492  };
493 
494 
495  /* ------------------------- Inline Functions -------------------------*/
496 
497  template<typename CoeffType>
498  inline
499  StatMechThermodynamics<CoeffType>::StatMechThermodynamics( const ChemicalMixture<CoeffType>& chem_mixture )
500  : _chem_mixture(chem_mixture)
501  {
502  // NOP
503  }
504 
505  template<typename CoeffType>
506  inline
507  StatMechThermodynamics<CoeffType>::~StatMechThermodynamics ()
508  {
509  // NOP
510  }
511 
512  template<typename CoeffType>
513  inline
514  CoeffType StatMechThermodynamics<CoeffType>::cv_trans( const unsigned int species ) const
515  {
516  return CoeffType(1.5)*_chem_mixture.R(species);
517  }
518 
519  template<typename CoeffType>
520  inline
521  CoeffType StatMechThermodynamics<CoeffType>::cv_trans_over_R( const unsigned int /*species*/ ) const
522  {
523  return CoeffType(1.5);
524  }
525 
526  template<typename CoeffType>
527  inline
528  CoeffType StatMechThermodynamics<CoeffType>::cv_rot( const unsigned int species ) const
529  {
530  using std::max;
531 
532  return max(this->cv_tr(species) - this->cv_trans(species), CoeffType(0) );
533  }
534 
535  template<typename CoeffType>
536  inline
537  CoeffType StatMechThermodynamics<CoeffType>::cv_rot_over_R( const unsigned int species ) const
538  {
539  using std::max;
540 
541  return max(this->cv_tr_over_R(species) - this->cv_trans_over_R(species), CoeffType(0) );
542  }
543 
544  template<typename CoeffType>
545  inline
546  CoeffType StatMechThermodynamics<CoeffType>::cv_tr (const unsigned int species) const
547  {
548  return _chem_mixture.R(species)*(_chem_mixture.chemical_species()[species])->n_tr_dofs();
549  }
550 
551  template<typename CoeffType>
552  inline
553  CoeffType StatMechThermodynamics<CoeffType>::cv_tr_over_R (const unsigned int species) const
554  {
555  return (_chem_mixture.chemical_species()[species])->n_tr_dofs();
556  }
557 
558  template<typename CoeffType>
559  template<typename VectorStateType>
560  inline
561  typename enable_if_c<
562  has_size<VectorStateType>::value,
563  typename Antioch::value_type<VectorStateType>::type
564  >::type
565  StatMechThermodynamics<CoeffType>::cv_tr (const VectorStateType& mass_fractions) const
566  {
567  typename Antioch::value_type<VectorStateType>::type
568  cv_tr = mass_fractions[0]*this->cv_tr(0);
569 
570  for( unsigned int s = 1; s < _chem_mixture.n_species(); s++ )
571  {
572  cv_tr += mass_fractions[s]*this->cv_tr(s);
573  }
574 
575  return cv_tr;
576  }
577 
578  template<typename CoeffType>
579  template<typename StateType>
580  inline
581  StateType StatMechThermodynamics<CoeffType>::cv_vib (const unsigned int species,
582  const StateType& Tv) const
583  {
584  return this->cv_vib_over_R(species,Tv) * (_chem_mixture.chemical_species()[species])->gas_constant();
585  }
586 
587  template<typename CoeffType>
588  template<typename StateType>
589  inline
590  StateType StatMechThermodynamics<CoeffType>::cv_vib_over_R (const unsigned int species,
591  const StateType& Tv) const
592  {
593  using std::exp;
594 
595  // convenience
596  const ChemicalSpecies<CoeffType>& chem_species = *(_chem_mixture.chemical_species()[species]);
597  const std::vector<CoeffType>& theta_v = chem_species.theta_v();
598  const std::vector<unsigned int>& ndg_v = chem_species.ndg_v();
599 
600  antioch_assert_equal_to(ndg_v.size(), theta_v.size());
601 
602  // Use an input datum to make sure we get the size right
603  StateType cv_vib_over_R = Antioch::zero_clone(Tv);
604 
605  if (theta_v.empty())
606  return cv_vib_over_R;
607 
608  for (unsigned int level=0; level<ndg_v.size(); level++)
609  {
610  typedef typename Antioch::raw_value_type<StateType>::type raw_type;
611  const StateType expval = exp(theta_v[level]/Tv);
612  const StateType expvalminusone = expval - raw_type(1);
613 
614  cv_vib_over_R += (static_cast<CoeffType>(ndg_v[level])*
615  theta_v[level]*theta_v[level]*expval/(expvalminusone*expvalminusone)/(Tv*Tv));
616  }
617 
618  return cv_vib_over_R;
619  }
620 
621  template<typename CoeffType>
622  template<typename VectorStateType>
623  inline
624  typename enable_if_c<
625  has_size<VectorStateType>::value,
626  typename Antioch::value_type<VectorStateType>::type
627  >::type
628  StatMechThermodynamics<CoeffType>::cv_vib (const typename Antioch::value_type<VectorStateType>::type& Tv,
629  const VectorStateType& mass_fractions) const
630  {
631  typename Antioch::value_type<VectorStateType>::type
632  cv_vib = mass_fractions[0]*this->cv_vib(0, Tv);
633 
634  for( unsigned int s = 1; s < _chem_mixture.n_species(); s++ )
635  {
636  cv_vib += mass_fractions[s]*this->cv_vib(s, Tv);
637  }
638 
639  return cv_vib;
640  }
641 
642  template<typename CoeffType>
643  template<typename StateType>
644  inline
645  StateType StatMechThermodynamics<CoeffType>::cv_el (const unsigned int species,
646  const StateType& Te) const
647  {
648  using std::exp;
649 
650  // convenience
651  const ChemicalSpecies<CoeffType>& chem_species = *(_chem_mixture.chemical_species()[species]);
652  const std::vector<CoeffType>& theta_e = chem_species.theta_e();
653  const std::vector<unsigned int>& ndg_e = chem_species.ndg_e();
654 
655  antioch_assert_equal_to(ndg_e.size(), theta_e.size());
656 
657  StateType cv_el = Antioch::zero_clone(Te);
658 
659  // Really < 2? Yes, b/c theta_e[0] = 0.0 always. See
660  // antioch_default_electronic_data.dat
661  if (theta_e.size() < 2) return cv_el;
662 
663  typedef typename Antioch::raw_value_type<StateType>::type raw_type;
664  const raw_type one = static_cast<raw_type>(1);
665 
666  const StateType Teinv = one/Te;
667  const StateType Te2inv = Teinv*Teinv;
668 
669  StateType
670  num = Antioch::zero_clone(Te), dnum = Antioch::zero_clone(Te),
671  den = Antioch::zero_clone(Te), dden = Antioch::zero_clone(Te);
672 
673  for (unsigned int level=0; level<theta_e.size(); level++)
674  {
675  const StateType
676  expval = exp (-theta_e[level] * Teinv),
677  den_l = static_cast<raw_type>(ndg_e[level])*expval,
678  num_l = den_l*theta_e[level],
679  dden_l = num_l*Te2inv,
680  dnum_l = dden_l*theta_e[level];
681 
682  num += num_l;
683  den += den_l;
684 
685  dden += dden_l;
686  dnum += dnum_l;
687  }
688 
689  const StateType invden = one/den;
690 
691  cv_el = chem_species.gas_constant()*(dnum - num*dden*invden) * invden;
692 
693  return cv_el;
694  }
695 
696  template<typename CoeffType>
697  template<typename VectorStateType>
698  inline
699  typename enable_if_c<
700  has_size<VectorStateType>::value,
701  typename Antioch::value_type<VectorStateType>::type
702  >::type
703  StatMechThermodynamics<CoeffType>::cv_el (const typename Antioch::value_type<VectorStateType>::type& Te,
704  const VectorStateType& mass_fractions) const
705  {
706  typename Antioch::value_type<VectorStateType>::type
707  cv_el = mass_fractions[0]*this->cv_el(0, Te);
708 
709  for( unsigned int s = 1; s < _chem_mixture.n_species(); s++ )
710  {
711  cv_el += mass_fractions[s]*this->cv_el(s, Te);
712  }
713 
714  return cv_el;
715  }
716 
717  template<typename CoeffType>
718  template<typename StateType>
719  inline
720  StateType StatMechThermodynamics<CoeffType>::cv_ve (const unsigned int species,
721  const StateType& Tv) const
722  {
723  return (this->cv_vib(species, Tv) + this->cv_el(species, Tv));
724  }
725 
726  template<typename CoeffType>
727  template<typename VectorStateType>
728  inline
729  typename enable_if_c<
730  has_size<VectorStateType>::value,
731  typename Antioch::value_type<VectorStateType>::type
732  >::type
733  StatMechThermodynamics<CoeffType>::cv_ve (const typename Antioch::value_type<VectorStateType>::type& Tv,
734  const VectorStateType& mass_fractions) const
735 
736  {
737  return (this->cv_vib(Tv, mass_fractions) + this->cv_el(Tv, mass_fractions));
738  }
739 
740  template<typename CoeffType>
741  template<typename StateType>
742  inline
743  StateType StatMechThermodynamics<CoeffType>::cv (const unsigned int species,
744  const StateType& /* T */,
745  const StateType& Tv) const
746  {
747  return (this->cv_tr(species) + this->cv_ve(species, Tv));
748  }
749 
750  template<typename CoeffType>
751  template<typename VectorStateType>
752  inline
753  typename enable_if_c<
754  has_size<VectorStateType>::value,
755  typename Antioch::value_type<VectorStateType>::type
756  >::type
757  StatMechThermodynamics<CoeffType>::cv (const typename Antioch::value_type<VectorStateType>::type& /* T */,
758  const typename Antioch::value_type<VectorStateType>::type& Tv,
759  const VectorStateType& mass_fractions) const
760  {
761  return (this->cv_tr(mass_fractions) + this->cv_ve(Tv, mass_fractions));
762  }
763 
764  template<typename CoeffType>
765  template<typename VectorStateType>
766  inline
767  typename enable_if_c<
768  has_size<VectorStateType>::value,
769  typename Antioch::value_type<VectorStateType>::type
770  >::type
771  StatMechThermodynamics<CoeffType>::cp (const typename Antioch::value_type<VectorStateType>::type& T,
772  const typename Antioch::value_type<VectorStateType>::type& Tv,
773  const VectorStateType& mass_fractions) const
774  {
775  return (this->cv(T,Tv,mass_fractions) + this->_chem_mixture.R(mass_fractions));
776  }
777 
778  template<typename CoeffType>
779  template<typename StateType>
780  inline
781  StateType StatMechThermodynamics<CoeffType>::h_tot (const unsigned int species,
782  const StateType& T,
783  const StateType& Tv) const
784  {
785  const ChemicalSpecies<CoeffType>& chem_species = *(_chem_mixture.chemical_species()[species]);
786  return (this->e_tot(species, T, Tv) + chem_species.gas_constant()*T);
787  }
788 
789 
790  template<typename CoeffType>
791  template<typename StateType>
792  inline
793  StateType StatMechThermodynamics<CoeffType>::h_tot (const unsigned int species,
794  const StateType& T) const
795  {
796  return this->h_tot(species, T, T);
797  }
798 
799  template<typename CoeffType>
800  template<typename VectorStateType>
801  inline
802  typename enable_if_c<
803  has_size<VectorStateType>::value,
804  typename Antioch::value_type<VectorStateType>::type
805  >::type
806  StatMechThermodynamics<CoeffType>::h_tot (const typename Antioch::value_type<VectorStateType>::type& T,
807  const typename Antioch::value_type<VectorStateType>::type& Tv,
808  const VectorStateType& mass_fractions) const
809  {
810  typename Antioch::value_type<VectorStateType>::type
811  h_tot = mass_fractions[0]*this->h_tot(0, T, Tv);
812 
813  for( unsigned int s = 1; s < _chem_mixture.n_species(); s++ )
814  {
815  h_tot += mass_fractions[s]*this->h_tot(s, T, Tv);
816  }
817 
818  return h_tot;
819  }
820 
821  template<typename CoeffType>
822  template<typename VectorStateType>
823  inline
824  typename enable_if_c<
825  has_size<VectorStateType>::value,
826  typename Antioch::value_type<VectorStateType>::type
827  >::type
828  StatMechThermodynamics<CoeffType>::h_tot (const typename Antioch::value_type<VectorStateType>::type& T,
829  const VectorStateType& mass_fractions) const
830  {
831  return this->h_tot(T, T, mass_fractions);
832  }
833 
834  template<typename CoeffType>
835  template<typename StateType>
836  inline
837  StateType StatMechThermodynamics<CoeffType>::e_tot (const unsigned int species,
838  const StateType& T,
839  const StateType& Tv) const
840  {
841  return (this->e_tr(species, T) + this->e_ve(species, Tv) + this->e_0(species));
842  }
843 
844 
845  template<typename CoeffType>
846  template<typename VectorStateType>
847  inline
848  typename enable_if_c<
849  has_size<VectorStateType>::value,
850  typename Antioch::value_type<VectorStateType>::type
851  >::type
852  StatMechThermodynamics<CoeffType>::e_tot (const typename Antioch::value_type<VectorStateType>::type& T,
853  const typename Antioch::value_type<VectorStateType>::type& Tv,
854  const VectorStateType& mass_fractions) const
855  {
856  return (this->e_tr(T, mass_fractions) + this->e_ve(Tv, mass_fractions) + this->e_0(mass_fractions));
857  }
858 
859  template<typename CoeffType>
860  template<typename StateType>
861  inline
862  StateType StatMechThermodynamics<CoeffType>::e_tot (const unsigned int species,
863  const StateType& T) const
864  {
865  return this->e_tot(species, T, T);
866  }
867 
868  template<typename CoeffType>
869  template<typename VectorStateType>
870  inline
871  typename enable_if_c<
872  has_size<VectorStateType>::value,
873  typename Antioch::value_type<VectorStateType>::type
874  >::type
875  StatMechThermodynamics<CoeffType>::e_tot (const typename Antioch::value_type<VectorStateType>::type& T,
876  const VectorStateType& mass_fractions) const
877  {
878  return this->e_tot(T, T, mass_fractions);
879  }
880 
881  template<typename CoeffType>
882  template<typename StateType>
883  inline
884  StateType StatMechThermodynamics<CoeffType>::e_tr (const unsigned int species,
885  const StateType& T) const
886  {
887  return this->cv_tr(species)*T;
888  }
889 
890  template<typename CoeffType>
891  template<typename VectorStateType>
892  inline
893  typename enable_if_c<
894  has_size<VectorStateType>::value,
895  typename Antioch::value_type<VectorStateType>::type
896  >::type
897  StatMechThermodynamics<CoeffType>::e_tr (const typename Antioch::value_type<VectorStateType>::type& T,
898  const VectorStateType& mass_fractions) const
899  {
900  typename Antioch::value_type<VectorStateType>::type
901  e_tr = mass_fractions[0]*this->e_tr(0, T);
902 
903  for( unsigned int s = 1; s < _chem_mixture.n_species(); s++ )
904  {
905  e_tr += mass_fractions[s]*this->e_tr(s, T);
906  }
907 
908  return e_tr;
909  }
910 
911  template<typename CoeffType>
912  template<typename StateType>
913  inline
914  StateType StatMechThermodynamics<CoeffType>::e_vib (const unsigned int species,
915  const StateType& Tv) const
916  {
917  using std::exp;
918 
919  // convenience
920  const ChemicalSpecies<CoeffType>& chem_species = *(_chem_mixture.chemical_species()[species]);
921  const std::vector<CoeffType>& theta_v = chem_species.theta_v();
922  const std::vector<unsigned int>& ndg_v = chem_species.ndg_v();
923 
924  antioch_assert_equal_to(ndg_v.size(), theta_v.size());
925 
926  StateType e_vib = 0.0;
927 
928  if (theta_v.empty()) return e_vib;
929 
930  for (unsigned int level=0; level<ndg_v.size(); level++)
931  e_vib += ndg_v[level]*chem_species.gas_constant()*theta_v[level]/(exp(theta_v[level]/Tv) - 1.);
932 
933  return e_vib;
934  }
935 
936  template<typename CoeffType>
937  template<typename VectorStateType>
938  inline
939  typename enable_if_c<
940  has_size<VectorStateType>::value,
941  typename Antioch::value_type<VectorStateType>::type
942  >::type
943  StatMechThermodynamics<CoeffType>::e_vib (const typename Antioch::value_type<VectorStateType>::type& Tv,
944  const VectorStateType& mass_fractions) const
945  {
946  typename Antioch::value_type<VectorStateType>::type
947  e_vib = mass_fractions[0]*this->e_vib(0, Tv);
948 
949  for( unsigned int s = 1; s < _chem_mixture.n_species(); s++ )
950  {
951  e_vib += mass_fractions[s]*this->e_vib(s, Tv);
952  }
953 
954  return e_vib;
955  }
956 
957  template<typename CoeffType>
958  template<typename StateType>
959  inline
960  StateType StatMechThermodynamics<CoeffType>::e_el (const unsigned int species,
961  const StateType& Te) const
962  {
963  using std::exp;
964 
965  // convenience
966  const ChemicalSpecies<CoeffType>& chem_species = *(_chem_mixture.chemical_species()[species]);
967  const std::vector<CoeffType>& theta_e = chem_species.theta_e();
968  const std::vector<unsigned int>& ndg_e = chem_species.ndg_e();
969 
970  antioch_assert_equal_to(ndg_e.size(), theta_e.size());
971 
972  StateType e_el = Antioch::zero_clone(Te);
973 
974  if (theta_e.size() < 2) return e_el;
975 
976  StateType num = Antioch::zero_clone(Te), den = Antioch::zero_clone(Te);
977 
978  for (unsigned int level=0; level<theta_e.size(); level++)
979  {
980  const StateType expval = exp (-theta_e[level] / Te);
981  num += static_cast<StateType>(ndg_e[level])*theta_e[level]*expval;
982  den += static_cast<StateType>(ndg_e[level])*expval;
983  }
984 
985  return chem_species.gas_constant() * num / den;
986  }
987 
988  template<typename CoeffType>
989  template<typename VectorStateType>
990  inline
991  typename enable_if_c<
992  has_size<VectorStateType>::value,
993  typename Antioch::value_type<VectorStateType>::type
994  >::type
995  StatMechThermodynamics<CoeffType>::e_el (const typename Antioch::value_type<VectorStateType>::type& Te,
996  const VectorStateType& mass_fractions) const
997  {
998  typename Antioch::value_type<VectorStateType>::type
999  e_el = mass_fractions[0]*this->e_el(0, Te);
1000 
1001  for( unsigned int s = 1; s < _chem_mixture.n_species(); s++ )
1002  {
1003  e_el += mass_fractions[s]*this->e_el(s, Te);
1004  }
1005 
1006  return e_el;
1007  }
1008 
1009  template<typename CoeffType>
1010  template<typename StateType>
1011  inline
1012  StateType StatMechThermodynamics<CoeffType>::e_ve (const unsigned int species,
1013  const StateType& Tv) const
1014  {
1015  return (this->e_vib(species,Tv) + this->e_el(species,Tv));
1016  }
1017 
1018  template<typename CoeffType>
1019  template<typename VectorStateType>
1020  inline
1021  typename enable_if_c<
1022  has_size<VectorStateType>::value,
1023  typename Antioch::value_type<VectorStateType>::type
1024  >::type
1025  StatMechThermodynamics<CoeffType>::e_ve (const typename Antioch::value_type<VectorStateType>::type& Tv,
1026  const VectorStateType& mass_fractions) const
1027  {
1028  return (this->e_vib(Tv,mass_fractions) + this->e_el(Tv,mass_fractions));
1029  }
1030 
1031  template<typename CoeffType>
1032  inline
1033  CoeffType StatMechThermodynamics<CoeffType>::e_ve_min () const
1034  {
1035  antioch_not_implemented();
1036  }
1037 
1038  template<typename CoeffType>
1039  template<typename VectorStateType>
1040  inline
1041  void StatMechThermodynamics<CoeffType>::e_and_cv_ve (const typename Antioch::value_type<VectorStateType>::type& Tv,
1042  const VectorStateType& mass_fractions,
1043  typename Antioch::value_type<VectorStateType>::type &e_ve,
1044  typename Antioch::value_type<VectorStateType>::type &cv_ve) const
1045  {
1046  e_ve = this->e_ve (Tv, mass_fractions);
1047  cv_ve = this->cv_ve(Tv, mass_fractions);
1048  }
1049 
1050  template<typename CoeffType>
1051  inline
1052  CoeffType StatMechThermodynamics<CoeffType>::e_0 (const unsigned int species) const
1053  {
1054  const ChemicalSpecies<CoeffType>& chem_species = *(_chem_mixture.chemical_species()[species]);
1055  return chem_species.formation_enthalpy();
1056  }
1057 
1058  template<typename CoeffType>
1059  template<typename VectorStateType>
1060  inline
1061  typename enable_if_c<
1062  has_size<VectorStateType>::value,
1063  typename Antioch::value_type<VectorStateType>::type
1064  >::type
1065  StatMechThermodynamics<CoeffType>::e_0 (const VectorStateType& mass_fractions) const
1066  {
1067  typename Antioch::value_type<VectorStateType>::type
1068  e_0 = mass_fractions[0]*this->e_0(0);
1069 
1070  for( unsigned int s = 1; s < _chem_mixture.n_species(); s++ )
1071  {
1072  e_0 += mass_fractions[s]*this->e_0(s);
1073  }
1074 
1075  return e_0;
1076  }
1077 
1078  template<typename CoeffType>
1079  template<typename VectorStateType>
1080  inline
1081  typename enable_if_c<
1082  has_size<VectorStateType>::value,
1083  typename Antioch::value_type<VectorStateType>::type
1084  >::type
1085  StatMechThermodynamics<CoeffType>::Tv_from_e_ve
1086  (const typename Antioch::value_type<VectorStateType>::type& e_ve,
1087  const VectorStateType& mass_fractions,
1088  typename Antioch::value_type<VectorStateType>::type Tv) const
1089  {
1090  antioch_not_implemented();
1091  }
1092 
1093  template<typename CoeffType>
1094  template<typename VectorStateType>
1095  inline
1096  typename enable_if_c<
1097  has_size<VectorStateType>::value,
1098  typename Antioch::value_type<VectorStateType>::type
1099  >::type
1100  StatMechThermodynamics<CoeffType>::T_from_e_tot (const typename Antioch::value_type<VectorStateType>::type& e_tot,
1101  const VectorStateType& mass_fractions,
1102  typename Antioch::value_type<VectorStateType>::type T) const
1103  {
1104  using std::abs;
1105  using std::max;
1106  using std::min;
1107 
1108  typedef typename Antioch::value_type<VectorStateType>::type StateType;
1109 
1110  // Cache the translational/rotational specific heat - this will be used repeatedly
1111  // and involves (2*NS-1) flops to compute, and since this has no functional
1112  // dependence on temperature it will not change throughout the Newton iteration.
1113  const typename Antioch::value_type<VectorStateType>::type
1114  Cv_tr = this->cv_tr(mass_fractions);
1115 
1116  // Similarly for mixture formation energy
1117  const typename Antioch::value_type<VectorStateType>::type
1118  E_0 = this->e_0(mass_fractions);
1119 
1120  // if the user does not provide an initial guess for the temperature
1121  // assume it is all in translation/rotation to compute a starting value.
1122  if (T < 0)
1123  {
1124  T = (e_tot - E_0) / Cv_tr;
1125  T = min(max(T,StateType(10.)),StateType(20000.));
1126 
1127  // FIXME: Use Antioch::Limits or similar? (i.e., don't
1128  // hardcode min and max T)
1129 
1130  // make sure the initial guess is valid
1131  //T = max(T, Limits::CompNSLimits::T_min());
1132  T = max(T, StateType(10.));
1133  T = min(T, StateType(2.e4));
1134  }
1135 
1136  // compute the translational/rotational temperature of the mixture using Newton-Rhapson iteration
1137  CoeffType delta_T = std::numeric_limits<StateType>::max();
1138  const unsigned int max_iterations = 100;
1139 
1140  const CoeffType dT_reltol= std::numeric_limits<CoeffType>::epsilon() * 100;
1141 
1142  // NOTE: FIN-S uses a hardcoded, absolute tolerance on delta_T of
1143  // 1e-8. Using a relative tolerance here of 100*epsilon.
1144  for (unsigned int iter = 0;
1145  abs(delta_T/T) > dT_reltol &&
1146  T >= 0.;
1147  ++iter)
1148  {
1149  if (iter == max_iterations)
1150  throw FailedNewtonTTvInversion ("ERROR: failed to converge T_from_e_tot!");
1151 
1152  // compute the residual, defined as the mismatch between the input e_tot and
1153  // the value corresponding to the current T iterate
1154  typename Antioch::value_type<VectorStateType>::type
1155  Re = 0., dRevdT = 0.;
1156 
1157  this->e_and_cv_ve(T, mass_fractions, Re, dRevdT);
1158  Re += this->e_tr(T, mass_fractions) + E_0 - e_tot;
1159  dRevdT += Cv_tr;
1160 
1161  delta_T = -Re / dRevdT;
1162  T += delta_T;
1163  }
1164 
1165  return T;
1166  }
1167 
1168  template<typename CoeffType>
1169  template<typename VectorStateType>
1170  inline
1171  typename enable_if_c<
1172  has_size<VectorStateType>::value,
1173  typename Antioch::value_type<VectorStateType>::type
1174  >::type
1175  StatMechThermodynamics<CoeffType>::T_from_e_tr
1176  (const typename Antioch::value_type<VectorStateType>::type& e_tr,
1177  const VectorStateType& mass_fractions,
1178  typename Antioch::value_type<VectorStateType>::type T ) const
1179  {
1180  antioch_not_implemented();
1181  }
1182 
1183  template<typename CoeffType>
1184  template<typename VectorStateType>
1185  inline
1186  typename enable_if_c<
1187  has_size<VectorStateType>::value,
1188  typename Antioch::value_type<VectorStateType>::type
1189  >::type
1190  StatMechThermodynamics<CoeffType>::T_from_h_tot
1191  (const typename Antioch::value_type<VectorStateType>::type& h_tot,
1192  const VectorStateType& mass_fractions,
1193  typename Antioch::value_type<VectorStateType>::type T) const
1194  {
1195  antioch_not_implemented();
1196  }
1197 
1198  template<typename CoeffType>
1199  template<typename VectorStateType>
1200  inline
1201  typename enable_if_c<
1202  has_size<VectorStateType>::value,
1203  typename Antioch::value_type<VectorStateType>::type
1204  >::type
1205  StatMechThermodynamics<CoeffType>::T_from_h_tot_Tv
1206  (const typename Antioch::value_type<VectorStateType>::type& h_tot,
1207  const typename Antioch::value_type<VectorStateType>::type& Tv,
1208  const VectorStateType& mass_fractions,
1209  typename Antioch::value_type<VectorStateType>::type T) const
1210  {
1211  antioch_not_implemented();
1212  }
1213 
1214 
1215  template<typename CoeffType>
1216  template<typename VectorStateType>
1217  inline
1218  typename enable_if_c<
1219  has_size<VectorStateType>::value,
1220  typename Antioch::value_type<VectorStateType>::type
1221  >::type
1222  StatMechThermodynamics<CoeffType>::s
1223  (const typename Antioch::value_type<VectorStateType>::type& T,
1224  const typename Antioch::value_type<VectorStateType>::type& p,
1225  const VectorStateType& mass_fractions) const
1226  {
1227  antioch_not_implemented();
1228  }
1229 
1230 } // end namespace Antioch
1231 
1232 #endif // ANTIOCH_STAT_MECH_THERMO_H
Antioch::StatMechThermodynamics::_chem_mixture
const ChemicalMixture< CoeffType > & _chem_mixture
Definition: stat_mech_thermo.h:485
Antioch::max
Antioch::enable_if_c< is_eigen< T >::value, typename value_type< T >::type >::type max(const T &in)
Definition: eigen_utils.h:88
Antioch::StatMechThermodynamics::T_from_e_tr
enable_if_c< has_size< VectorStateType >::value, typename Antioch::value_type< VectorStateType >::type >::type T_from_e_tr(const typename Antioch::value_type< VectorStateType >::type &e_tr, const VectorStateType &mass_fractions, typename Antioch::value_type< VectorStateType >::type T=-1) const
antioch_not_implemented
#define antioch_not_implemented()
Definition: antioch_asserts.h:105
Antioch::ChemicalSpecies::formation_enthalpy
CoeffType formation_enthalpy() const
Returns formation enthalpy in units of [J/kg].
Definition: chemical_species.h:201
Antioch::StatMechThermodynamics::cv_tr
CoeffType cv_tr(const unsigned int species) const
Definition: stat_mech_thermo.h:546
Antioch::StatMechThermodynamics::cp
enable_if_c< has_size< VectorStateType >::value, typename Antioch::value_type< VectorStateType >::type >::type cp(const typename Antioch::value_type< VectorStateType >::type &T, const typename Antioch::value_type< VectorStateType >::type &Tv, const VectorStateType &mass_fractions) const
antioch_assert_equal_to
#define antioch_assert_equal_to(expr1, expr2)
Definition: antioch_asserts.h:74
Antioch::StatMechThermodynamics::s
enable_if_c< has_size< VectorStateType >::value, typename Antioch::value_type< VectorStateType >::type >::type s(const typename Antioch::value_type< VectorStateType >::type &T, const typename Antioch::value_type< VectorStateType >::type &p, const VectorStateType &mass_fractions) const
Antioch::StatMechThermodynamics::cv_tr_over_R
CoeffType cv_tr_over_R(const unsigned int species) const
Definition: stat_mech_thermo.h:553
Antioch::FailedNewtonTTvInversion
A class representing failed Newton convergence in T/Tv inversion.
Definition: antioch_exceptions.h:73
Antioch::StatMechThermodynamics::T_from_e_tot
enable_if_c< has_size< VectorStateType >::value, typename Antioch::value_type< VectorStateType >::type >::type T_from_e_tot(const typename Antioch::value_type< VectorStateType >::type &e_tot, const VectorStateType &mass_fractions, typename Antioch::value_type< VectorStateType >::type T=-1) const
Antioch::StatMechThermodynamics::h_tot
StateType h_tot(const unsigned int species, const StateType &T, const StateType &Tv) const
Definition: stat_mech_thermo.h:781
Antioch::StatMechThermodynamics::e_ve
StateType e_ve(const unsigned int species, const StateType &Tv) const
Definition: stat_mech_thermo.h:1012
Antioch::StatMechThermodynamics::e_vib
StateType e_vib(const unsigned int species, const StateType &Tv) const
Definition: stat_mech_thermo.h:914
Antioch::min
Antioch::enable_if_c< is_eigen< T >::value, typename value_type< T >::type >::type min(const T &in)
Definition: eigen_utils.h:98
Antioch::StatMechThermodynamics::cv_el
StateType cv_el(const unsigned int species, const StateType &Te) const
Definition: stat_mech_thermo.h:645
Antioch::StatMechThermodynamics::cv_rot_over_R
CoeffType cv_rot_over_R(const unsigned int species) const
Definition: stat_mech_thermo.h:537
Antioch::StatMechThermodynamics::StatMechThermodynamics
StatMechThermodynamics()
Default constructor.
Antioch::StatMechThermodynamics::e_el
StateType e_el(const unsigned int species, const StateType &Te) const
Definition: stat_mech_thermo.h:960
Antioch::StatMechThermodynamics::T_from_h_tot_Tv
enable_if_c< has_size< VectorStateType >::value, typename Antioch::value_type< VectorStateType >::type >::type T_from_h_tot_Tv(const typename Antioch::value_type< VectorStateType >::type &h_tot, const typename Antioch::value_type< VectorStateType >::type &Tv, const VectorStateType &mass_fractions, typename Antioch::value_type< VectorStateType >::type T=-1) const
Antioch::ChemicalSpecies::theta_v
const std::vector< CoeffType > & theta_v() const
Characteristic vibrational temperature [K].
Definition: chemical_species.h:229
Antioch::ChemicalSpecies::ndg_v
const std::vector< unsigned int > & ndg_v() const
Degeneracies for each vibrational mode.
Definition: chemical_species.h:234
std::max
max(const _Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > &a, const _Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > &b) ANTIOCH_AUTOFUNC(_Matrix< _Scalar ANTIOCH_COMMA _Rows ANTIOCH_COMMA _Cols ANTIOCH_COMMA _Options ANTIOCH_COMMA _MaxRows ANTIOCH_COMMA _MaxCols >
Antioch::StatMechThermodynamics::cv_rot
CoeffType cv_rot(const unsigned int species) const
Definition: stat_mech_thermo.h:528
Antioch::ChemicalSpecies
Class to encapsulate data for each chemical species.
Definition: chemical_species.h:54
Antioch::StatMechThermodynamics::~StatMechThermodynamics
virtual ~StatMechThermodynamics()
Destructor.
Definition: stat_mech_thermo.h:507
Antioch::StatMechThermodynamics::Tv_from_e_ve
enable_if_c< has_size< VectorStateType >::value, typename Antioch::value_type< VectorStateType >::type >::type Tv_from_e_ve(const typename Antioch::value_type< VectorStateType >::type &e_ve, const VectorStateType &mass_fractions, typename Antioch::value_type< VectorStateType >::type Tv=-1) const
Antioch::StatMechThermodynamics::e_tot
StateType e_tot(const unsigned int species, const StateType &T, const StateType &Tv) const
Definition: stat_mech_thermo.h:837
Antioch::StatMechThermodynamics::e_tr
StateType e_tr(const unsigned int species, const StateType &T) const
Definition: stat_mech_thermo.h:884
Antioch::StatMechThermodynamics::cv
StateType cv(const unsigned int species, const StateType &T, const StateType &Tv) const
Definition: stat_mech_thermo.h:743
Antioch::StatMechThermodynamics::cv_ve
StateType cv_ve(const unsigned int species, const StateType &Tv) const
Definition: stat_mech_thermo.h:720
Antioch::StatMechThermodynamics::cv_trans_over_R
CoeffType cv_trans_over_R(const unsigned int species) const
Definition: stat_mech_thermo.h:521
Antioch::ChemicalSpecies::ndg_e
const std::vector< unsigned int > & ndg_e() const
Degeneracies for each electronic modes.
Definition: chemical_species.h:244
Antioch::StatMechThermodynamics::cv_trans
CoeffType cv_trans(const unsigned int species) const
Definition: stat_mech_thermo.h:514
Antioch::StatMechThermodynamics::cv_vib
StateType cv_vib(const unsigned int species, const StateType &Tv) const
Definition: stat_mech_thermo.h:581
Antioch::ChemicalSpecies::gas_constant
CoeffType gas_constant() const
Returns the species ideal gas constant in [J/kg-K].
Definition: chemical_species.h:196
Antioch::ChemicalMixture
Class storing chemical mixture properties.
Definition: chemical_mixture.h:75
Antioch
The parameters are reduced parameters.
Definition: chemical_mixture.h:47
std::min
typename _Scalar int _Rows int _Cols int _Options int _MaxRows int _MaxCols inline min(const _Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > &a, const _Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > &b) ANTIOCH_AUTOFUNC(_Matrix< _Scalar ANTIOCH_COMMA _Rows ANTIOCH_COMMA _Cols ANTIOCH_COMMA _Options ANTIOCH_COMMA _MaxRows ANTIOCH_COMMA _MaxCols >
Antioch::StatMechThermodynamics::e_0
CoeffType e_0(const unsigned int species) const
Definition: stat_mech_thermo.h:1052
Antioch::zero_clone
_Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > zero_clone(const _Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > &ex)
Definition: eigen_utils.h:145
Antioch::StatMechThermodynamics::T_from_h_tot
enable_if_c< has_size< VectorStateType >::value, typename Antioch::value_type< VectorStateType >::type >::type T_from_h_tot(const typename Antioch::value_type< VectorStateType >::type &h_tot, const VectorStateType &mass_fractions, typename Antioch::value_type< VectorStateType >::type T=-1) const
Antioch::StatMechThermodynamics::e_and_cv_ve
void e_and_cv_ve(const typename Antioch::value_type< VectorStateType >::type &Tv, const VectorStateType &mass_fractions, typename Antioch::value_type< VectorStateType >::type &e_ve, typename Antioch::value_type< VectorStateType >::type &cv_ve) const
Definition: stat_mech_thermo.h:1041
Antioch::ChemicalSpecies::theta_e
const std::vector< CoeffType > & theta_e() const
Characteristic electronic excitation temperatures [K].
Definition: chemical_species.h:239
Antioch::StatMechThermodynamics::cv_vib_over_R
StateType cv_vib_over_R(const unsigned int species, const StateType &Tv) const
Definition: stat_mech_thermo.h:590
Generated on Thu Jul 7 2016 11:09:46 for antioch-0.4.0 by  doxygen 1.8.8