research/cantera
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

*** empty log message ***

Browse source
This commit is contained in:
Dave Goodwin 2007-12-17 07:51:37 +00:00
parent fbf5953af7
commit dc94ba0dec
4 changed files with 104 additions and 21 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

74
Cantera/src/spectra/rotor.cpp
View file

@ -1,4 +1,7 @@
/**
* @file rotor.cpp
*
*/
#ifdef WIN32
#pragma warning(disable:4786)
#pragma warning(disable:4503)
@ -9,10 +12,13 @@
namespace Cantera {
/*
* @param mu reduced mass in kg
* @param re bond length in meters
* @dipoleMoment permanent dipole moment in ...
/**
* @param Bv Rotational constant, wavenumbers
* @param Dv Coefficient describing centrifugal
* effects on the bond length. For a rigid rotor, Bv = 0.
* @param Hv Coefficient describing higher-order vibration-rotation
* interactions. For a rigid rotor, Hv = 0.
* @dipoleMoment permanent dipole moment.
*/
Rotor::Rotor(doublereal Bv, doublereal dipoleMoment,
doublereal Dv, doublereal Hv ) : m_Bv(Bv),
@ -20,30 +26,72 @@ namespace Cantera {
m_Hv(Hv),
m_dipole(dipoleMoment) {}
// energy in wavenumbers
/**
* The energy of the level with rotational quantum number J,
* in wavenumber units.
* \f[
* E(J) = J(J+1)B - [J(J+1)]^2 D + [J(J+1)]^3H
* \f]
* For a rigid rotor, only B is non-zero. The parameters B, D, and H
* are set in the constructor.
*/
doublereal Rotor::energy_w(int J) {
int jjp1 = J*(J + 1);
return jjp1*(m_Bv + jjp1*(m_Hv*jjp1 - m_Dv));
}
doublereal Rotor::degeneracy(int J) {
/**
* The number of quantum states with the same J.
*/
int Rotor::degeneracy(int J) {
return 2*J + 1;
}
/**
* The rotational partition function.
*
* If T/Trot > 100, then the classical value (T/Trot) is
* is returned. Otherwise, it is computed as a sum
* \f[
* z = \sum_{J=0}^{J_{max} (2J + 1) \exp(-E(J)/kT)
* \]
*/
doublereal Rotor::partitionFunction(doublereal T, int cutoff) {
int j = 0;
if (cutoff < 0) cutoff = 100;
doublereal dsum = 0.0, sum = 0.0;
for (j = 0; j < cutoff; j++) {
dsum = degeneracy(j)*exp(-Planck*energy_w(j)/(Boltzmann * T));
sum += dsum;
doublereal T_Trot = wnum_to_J(m_Bv)/(Boltzmann*T);
if (T_Trot > 100.0)
return T_Trot;
else {
if (cutoff < 0) cutoff = 3.0*sqrt(T/m_Bv);
doublereal dsum = 0.0, sum = 0.0;
for (j = 0; j < cutoff; j++) {
dsum = degeneracy(j)*exp(-wnum_to_J(energy_w(j))/(Boltzmann * T));
sum += dsum;
}
return sum;
}
return sum;
}
/**
* Ratio of the population of all states with rotational quantum
* number J to the ground state population.
*/
doublereal Rotor::relPopulation(int J, doublereal T) {
return degeneracy(J)*exp(-wnum_to_J(energy_w(J))/(Boltzmann*T));
}
/** The difference in the energies of an upper and a lower state.
*
*/
doublereal Rotor::frequency(int J_lower, int J_upper) {
return (energy_w(J_upper) - energy_w(J_lower));
}
doublereal Rotor::intensity(int J_lower, int J_upper, doublereal T) {
int dJ = J_upper - J_lower;
if (dJ > 1 || dJ < -1) return 0;
return relPopulation(J_lower, T);
}
}

21
Cantera/src/spectra/rotor.h
View file

@ -1,3 +1,8 @@
/**
* @file rotor.h
* Header file for class Rotor.
*/
#ifndef CT_ROTOR
#define CT_ROTOR
@ -5,6 +10,10 @@
namespace Cantera {
/**
* Class Rotor represents a non-rigid quantum-mechanical rotor.
* @ingroup Spectroscopy
*/
class Rotor {
public:
Rotor() {}
@ -16,15 +25,23 @@ namespace Cantera {
doublereal Dv = 0.0, doublereal Hv = 0.0);
// energy in wavenumbers
/// energy in wavenumbers
doublereal energy_w(int J);
doublereal degeneracy(int J);
/// degeneracy
int degeneracy(int J);
doublereal partitionFunction(doublereal T, int cutoff=-1);
doublereal frequency(int J_lower, int J_upper);
doublereal relPopulation(int J, doublereal T);
doublereal population(int J, doublereal T) {
return relPopulation(J,T)/partitionFunction(T);
}
doublereal intensity(int J_lower, int J_upper, doublereal T);
protected:
doublereal m_Bv;

9
test_problems/spectroscopy/Makefile.in
View file

@ -53,6 +53,13 @@ CANTERA_LIBDIR=@buildlib@
# required Cantera libraries
CANTERA_LIBS = @LOCAL_LIBS@ -lctcxx
# Cantera library Dependencies
CANTERA_LIB_DEPS = $(CANTERA_LIBDIR)/libtransport.a \
$(CANTERA_LIBDIR)/libthermo.a \
$(CANTERA_LIBDIR)/libctspectra.a \
$(CANTERA_LIBDIR)/libctnumerics.a \
$(CANTERA_LIBDIR)/libctbase.a
# the directory where Cantera include files may be found.
ifeq ($(src_dir_tree), 1)
CANTERA_INCDIR=../../Cantera/src
@ -77,7 +84,7 @@ PROGRAM = $(PROG_NAME)$(EXE_EXT)
all: $(PROGRAM)
$(PROGRAM): $(OBJS) $(CANTERA_LIBDIR)/libctbase.a
$(PROGRAM): $(OBJS) $(CANTERA_LIB_DEPS)
$(CXX) -o $(PROGRAM) $(OBJS) $(LCXX_FLAGS) $(LINK_OPTIONS) \
$(CANTERA_LIBS) @LIBS@ $(FORT_LIBS) \
$(LCXX_END_LIBS)

21
test_problems/spectroscopy/spectratest.cpp
View file

@ -6,9 +6,20 @@ using namespace std;
using namespace Cantera;
int main() {
Rotor* r = new Rotor(1.0);
double w = 8065.0;
cout << "eV: " << wnum_to_eV(w) << endl;
delete r;
double B;
double T;
cout << "enter B, T: ";
cin >> B >> T;
Rotor* r = new Rotor(B);
double theta = wnum_to_J(B)/Boltzmann;
cout << "theta = " << theta << endl;
double pop, cpop = 0.0;
for (int j = 0; j < 50; j++) {
pop = r->population(j, T);
cpop += pop;
if (cpop > 0.999) break;
cout << j << ", " << r->energy_w(j) << ", "
<< r->frequency(j, j+1) << ", "
<< pop << ", " << cpop << ", " << r->partitionFunction(T) << ", " << T/theta << endl;
}
}