From fb588d873b60a0eb455e328a043534c00315f56d Mon Sep 17 00:00:00 2001 From: Ray Speth Date: Fri, 30 May 2014 19:23:37 +0000 Subject: [PATCH] Remove deprecated / unfinished 'spectra' classes --- doc/doxygen/Doxyfile | 1 - include/cantera/radiation.h | 7 - include/cantera/spectra.h | 7 - include/cantera/spectra/LineBroadener.h | 151 ------------------ include/cantera/spectra/rotor.h | 99 ------------ src/SConscript | 1 - src/spectra/DiatomicMolecule.h | 18 --- src/spectra/LineBroadener.cpp | 174 --------------------- src/spectra/Nuclei.h | 131 ---------------- src/spectra/rotor.cpp | 109 ------------- src/spectra/spectralUtilities.cpp | 12 -- src/spectra/spectralUtilities.h | 6 - test_problems/SConscript | 1 - test_problems/spectroscopy/spectratest.cpp | 72 --------- 14 files changed, 789 deletions(-) delete mode 100644 include/cantera/radiation.h delete mode 100644 include/cantera/spectra.h delete mode 100644 include/cantera/spectra/LineBroadener.h delete mode 100644 include/cantera/spectra/rotor.h delete mode 100644 src/spectra/DiatomicMolecule.h delete mode 100644 src/spectra/LineBroadener.cpp delete mode 100644 src/spectra/Nuclei.h delete mode 100644 src/spectra/rotor.cpp delete mode 100644 src/spectra/spectralUtilities.cpp delete mode 100644 src/spectra/spectralUtilities.h delete mode 100644 test_problems/spectroscopy/spectratest.cpp diff --git a/doc/doxygen/Doxyfile b/doc/doxygen/Doxyfile index 493f3377e..abfd5170f 100644 --- a/doc/doxygen/Doxyfile +++ b/doc/doxygen/Doxyfile @@ -577,7 +577,6 @@ INPUT = src/apps \ src/kinetics \ src/numerics \ src/oneD \ - src/spectra \ src/thermo \ src/tpx \ src/transport \ diff --git a/include/cantera/radiation.h b/include/cantera/radiation.h deleted file mode 100644 index 3c4d20d13..000000000 --- a/include/cantera/radiation.h +++ /dev/null @@ -1,7 +0,0 @@ -#ifndef CT_RAD_H_INCL -#define CT_RAD_H_INCL - -#include "spectra/rotor.h" -#include "spectra/LineBroadener.h" - -#endif diff --git a/include/cantera/spectra.h b/include/cantera/spectra.h deleted file mode 100644 index 9be836294..000000000 --- a/include/cantera/spectra.h +++ /dev/null @@ -1,7 +0,0 @@ -#ifndef CT_SPECTRA_H_INCL -#define CT_SECTRA_H_INCL - -#include "spectra/LineBroadener.h" -#include "spectra/rotor.h" - -#endif diff --git a/include/cantera/spectra/LineBroadener.h b/include/cantera/spectra/LineBroadener.h deleted file mode 100644 index 0ffdecde3..000000000 --- a/include/cantera/spectra/LineBroadener.h +++ /dev/null @@ -1,151 +0,0 @@ -/** - * @file LineBroadener.h - * Header file for class LineBroadener - * @ingroup spectroscopy - */ - -#include "cantera/base/ct_defs.h" -#include "cantera/base/ctexceptions.h" -#include "cantera/base/global.h" - -namespace Cantera -{ - -/** - * Base class for classes implementing line shapes of - * various types. - * @ingroup spectroscopy - * @deprecated incomplete / abandoned - */ -class LineBroadener -{ - -public: - - /// Default constructor - LineBroadener() { - warn_deprecated("class LineBroadener"); - } - - /// Destructor - virtual ~LineBroadener() {} - - /** - * The line shape profile - * \f[ - * P(\Delta\nu) - *\f] - * as a function of distance from line - * center \f$ \Delta\nu \f$. - * This function must have total area = 1.0. - * Note that this method must be overloaded in each - * derived class. If the base class method is called, - * an exception will be thrown. - */ - virtual doublereal profile(doublereal deltaFreq) { - throw CanteraError("LineBroadener::profile", - "base class method called!"); - } - - doublereal operator()(doublereal deltaFreq) { - return profile(deltaFreq); - } - - /** - * The cumulative profile, defined as - * \f[ - * C(\Delta \nu) = \int_{-\infty}^{\Delta \nu} P(x) dx - * \f] - */ - virtual doublereal cumulative(doublereal deltaFreq) { - throw CanteraError("LineBroadener::cumulative", - "base class method called!"); - } - - virtual doublereal width() { - return 0.0; - } -}; - -/** - * The line shape for pure collisional broadening. The Lorentzian line - * shape is - * \f[ - * L(\Delta\nu) = \frac{1}{\pi}\frac{\gamma}{\Delta\nu^2 + \gamma^2} - * \f] - * where \f$ \gamma = {\mbox{FWHM}/2} \f$. - */ -class LorentzianProfile : public LineBroadener -{ -public: - LorentzianProfile(doublereal FWHM); - virtual doublereal profile(doublereal deltaFreq); - virtual doublereal cumulative(doublereal deltaFreq); - virtual double width(); - -protected: - doublereal m_hwhm; - doublereal m_hwhm2; -}; - -/** - * A Gaussian line profile. This profile results when Doppler - * broadening is dominant. - */ -class GaussianProfile : public LineBroadener -{ -public: - - /** - * Constructor. - */ - GaussianProfile(doublereal sigma); - virtual doublereal profile(doublereal deltaFreq); - virtual doublereal cumulative(doublereal deltaFreq); - virtual doublereal width(); - - doublereal standardDev() { - return m_sigma; - } - -protected: - doublereal m_sigma; - doublereal m_sigma2; - doublereal m_width; -}; - - -/** - * A Voigt profile is the convolution of a Lorentzian and a - * Gaussian profile. This profile results when Doppler - * broadening and collisional broadening both are important. - */ -class Voigt : public LineBroadener -{ -public: - - /** - * Constructor. - */ - Voigt(doublereal sigma, doublereal gamma); - virtual doublereal profile(doublereal deltaFreq); - //virtual doublereal cumulative(doublereal deltaFreq) - //virtual doublereal width() - - void testv(); - -protected: - - doublereal F(doublereal x); - - doublereal m_sigma; - doublereal m_gamma_lor; - doublereal m_sigma2; - doublereal m_width; - doublereal m_gamma; - doublereal m_sigsqrt2; - doublereal m_a; - doublereal m_eps; -}; - -} diff --git a/include/cantera/spectra/rotor.h b/include/cantera/spectra/rotor.h deleted file mode 100644 index fb0efa110..000000000 --- a/include/cantera/spectra/rotor.h +++ /dev/null @@ -1,99 +0,0 @@ -#ifndef CT_ROTOR -#define CT_ROTOR - -/** - * @file rotor.h - * Header file for class Rotor. - */ - -/** - * @defgroup spectroscopy Spectroscopic Models - * - * These classes are used to simulate the absorption and emission spectra of - * molecules. - */ - -#include "cantera/base/ct_defs.h" -#include "cantera/base/global.h" - -/** - * Namespace for spectroscopic functions and classes. - */ -namespace Cantera -{ - -/** - * Class Rotor represents a non-rigid quantum-mechanical rotor. - * @ingroup spectroscopy - * @deprecated incomplete / abandoned - */ -class Rotor -{ -public: - - /// Default Constructor. - Rotor() { - warn_deprecated("class Rotor"); - } - - /// Destructor. - virtual ~Rotor() {} - - /// Full Constructor. - Rotor(doublereal Bv, doublereal dipoleMoment = 0.0, - doublereal Dv = 0.0, doublereal Hv = 0.0); - - doublereal energy_w(int J); - - 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; - doublereal m_Dv; - doublereal m_Hv; - doublereal m_dipole; -}; - -/** convert from Hz to wavenmbers */ -inline doublereal hz_to_wnum(doublereal freq) -{ - return freq/(100.0*Cantera::lightSpeed); -} - -/** Convert from wavenumbers to Joules. */ -inline doublereal wnum_to_J(doublereal w) -{ - return Cantera::Planck * w * 100.0 * Cantera::lightSpeed; -} - -inline doublereal J_to_wnum(doublereal e) -{ - return e /(Cantera::Planck * 100.0 * Cantera::lightSpeed); -} - -inline doublereal wnum_to_eV(doublereal w) -{ - return Cantera::Planck * w * 100.0 * Cantera::lightSpeed / Cantera::ElectronCharge; -} - -inline doublereal eV_to_wnum(doublereal e) -{ - return e * Cantera::ElectronCharge / (Cantera::Planck * 100.0 * Cantera::lightSpeed); -} -} - -#endif - - diff --git a/src/SConscript b/src/SConscript index 190b9c8d3..0d1f3b766 100644 --- a/src/SConscript +++ b/src/SConscript @@ -30,7 +30,6 @@ libs = [('base', ['cpp'], baseSetup), ('numerics', ['cpp'], numericsSetup), ('kinetics', ['cpp'], defaultSetup), ('transport', ['cpp'], defaultSetup), - ('spectra', ['cpp'], defaultSetup), ('oneD', ['cpp'], defaultSetup), ('zeroD', ['cpp'], defaultSetup), ('clib', ['cpp'], defaultSetup), diff --git a/src/spectra/DiatomicMolecule.h b/src/spectra/DiatomicMolecule.h deleted file mode 100644 index 36b252962..000000000 --- a/src/spectra/DiatomicMolecule.h +++ /dev/null @@ -1,18 +0,0 @@ -#ifndef CT_DIATOMIC_H -#define CT_DIATOMIC_H - -#include "Nuclei.h" - -namespace Cantera -{ - -class DiatomicMolecule -{ -public: - DiatomicMolecule(Nucleus n1, Nucleus n2, doublereal - doublereal absorptionCrossSection(doublereal nu); -}; -} - -#endif - diff --git a/src/spectra/LineBroadener.cpp b/src/spectra/LineBroadener.cpp deleted file mode 100644 index 81258d8e0..000000000 --- a/src/spectra/LineBroadener.cpp +++ /dev/null @@ -1,174 +0,0 @@ -#include "cantera/base/ct_defs.h" -#include -#include - -#ifdef USE_BOOST_MATH -#include -using boost::math::erf; -#endif - -#include "cantera/spectra/LineBroadener.h" - -namespace Cantera -{ - -LorentzianProfile::LorentzianProfile(doublereal gamma) -{ - m_hwhm = gamma; - m_hwhm2 = m_hwhm*m_hwhm; -} - -/** - * The Lorentzian profile for collision-broadened lines. - * - *\f[ - * \frac{1}{\pi} \frac{\gamma}{ (\Delta\nu)^2 + \gamma^2} - *\f] - * Units: 1/wavenumber (or cm). - */ -doublereal LorentzianProfile::profile(doublereal deltaFreq) -{ - return (1.0/Cantera::Pi) *m_hwhm/(deltaFreq*deltaFreq + m_hwhm2); -} - -/** - * - * The cumulative profile, given by - * \f[ - * \frac{1}{\pi} \tan^{-1}\left(\frac{\Delta\nu}{gamma}\right) + 0.5 - * \f] - */ -doublereal LorentzianProfile::cumulative(doublereal deltaFreq) -{ - return (1.0/Pi) * atan(deltaFreq/m_hwhm) + 0.5; -} - -doublereal LorentzianProfile::width() -{ - return 2.0*m_hwhm; -} - -GaussianProfile::GaussianProfile(doublereal sigma) -{ - m_sigma = sigma; - m_sigma2 = m_sigma*m_sigma; -} - -doublereal GaussianProfile::profile(doublereal deltaFreq) -{ - //cout << "entered Gaussian::profile" << endl; - //cout << "deltaFreq = " << deltaFreq << endl; - //cout << "m_sigma = " << m_sigma << endl; - return 1.0/(m_sigma * Cantera::SqrtTwo *Cantera::SqrtPi) * - exp(-deltaFreq*deltaFreq/(2.0*m_sigma2)); -} - -doublereal GaussianProfile::cumulative(doublereal deltaFreq) -{ - return 0.5*(1.0 + erf(deltaFreq/(m_sigma*SqrtTwo))); -} - -doublereal GaussianProfile::width() -{ - return 2.0*m_sigma*sqrt(log(4.0)); -} - - - -/** - * @param sigma The standard deviation of the Gaussian - * @param gamma The half-width of the Lorentzian. - */ -Voigt::Voigt(doublereal sigma, doublereal gamma) -{ - m_sigma = sigma; - m_sigma2 = m_sigma*m_sigma; - m_gamma_lor = gamma; - m_sigsqrt2 = SqrtTwo*m_sigma; - m_gamma = gamma/m_sigsqrt2; - m_eps = 1.0e-20; -} - -void Voigt::testv() -{ - m_gamma = 1.0e1; - std::cout << F(1.0) << std::endl; - m_gamma = 0.5; - std::cout << F(1.0) << std::endl; - m_gamma = 0.0001; - std::cout << F(10.0) << std::endl; -} -/** - * This method evaluates the function - * \f[ - * F(x, y) = \frac{y}{\pi}\int_{-\infty}^{+\infty} \frac{e^{-z^2}} - * {(x - z)^2 + y^2} dz - * \f] - * The algorithm used to cmpute this function is described in the - * reference below. @see F. G. Lether and P. R. Wenston, "The - * numerical computation of the %Voigt function by a corrected - * midpoint quadrature rule for \f$ (-\infty, \infty) \f$. Journal - * of Computational and Applied Mathematics}, 34 (1):75--92, 1991. - */ -doublereal Voigt::F(doublereal x) -{ - - if (x < 0.0) { - x = -x; - } - double y = m_gamma; - - double c3 = log(Pi*m_eps/2.0); - double tau = sqrt(-log(y) - c3); - double b = (tau + x)/y; - double t = b*y; - double f1, f2, f3; - const double c0 = 2.0/(Pi*exp(0.0)); - const double c1 = 1.0/SqrtTwo; - const double c2 = 2.0/SqrtPi; - - if (y > c0/m_eps) { - return 0.0; - } - double f0, ef0; - while (1 > 0) { - f0 = Pi*Pi/(t*t); - ef0 = exp(-f0); - f1 = c2*y*ef0; - f2 = fabs(y*y - Pi*Pi/(t*t)); - f3 = 1.0 - ef0*ef0; - t *= c1; - if (f1/(f2*f3) < 0.5*m_eps) { - break; - } - } - double h = t/y; - int N = int(0.5 + b/h); - double S = 0.0; - double u = h/2; - for (int i = 0; i < N; i++) { - S += (1.0 + exp(-4.0*x*y*u))*exp(-pow(y*u-x,2))/(u*u+1.0); - u += h; - } - double Q = h*S/Pi; - double C = 0.0; - if (y*y < Pi/h) { - C = 2.0*exp(y*y - x*x)*cos(2*x*y)/(1.0 + exp(2*Pi/h)); - } else { - return 0.0; - } - return Q + C; -} - -/** - * Voigt profile. - * - * Not sure that constant is right. - */ -doublereal Voigt::profile(doublereal deltaFreq) -{ - const double ff = 1.0/(m_sigsqrt2*SqrtPi); - return ff*F(deltaFreq/m_sigsqrt2); -} - -} diff --git a/src/spectra/Nuclei.h b/src/spectra/Nuclei.h deleted file mode 100644 index 130e42fc3..000000000 --- a/src/spectra/Nuclei.h +++ /dev/null @@ -1,131 +0,0 @@ -/** - * @file Nuclei.h Provides class Nucleus - */ - -#ifndef CT_NUCL_H -#define CT_NUCL_H - -#include "cantera/base/ct_defs.h" -#include "cantera/base/global.h" - -namespace Cantera -{ - -/** - * Represents atomic nuclei. These classes only provide minimal - * information, and are designed only to handle nuclear statistics - * effects on spectra. - * @ingroup spectroscopy - * @deprecated incomplete / abandoned - */ -class Nucleus -{ -public: - Nucleus(const std::string& symbol, - int nP, int nN, doublereal spin) : m_np(nP), - m_nn(nN), m_spin(spin), - m_sym(symbol) { - warn_deprecated("class Nucleus"); - } - virtual ~Nucleus() {} - int nProtons() { - return m_np; - } - int mNeutrons() { - return m_nn; - } - doublereal spin() { - return m_spin; - } - int multiplicity() { - return (int)(2*m_spin) + 1; - } - int atomicNumber() { - return m_np; - } - std::string symbol() { - return m_sym; - } - - bool operator==(Nucleus& b) { - if (m_np == b.m_np && m_nn == b.m_nn) { - return true; - } else { - return false; - } - } - - bool operator!=(Nucleus& b) { - return !(*this == b); - } - - bool isBoson() { - return (m_spin - std::floor(m_spin) < 0.001); - } - - -protected: - - int m_np; //< Number of protons - int m_nn; //< Number of electrons - doublereal m_spin; //< Spin. - std::string m_sym; //< Symbol. -}; - -inline Nucleus* HydrogenNucleus() -{ - return new Nucleus("H", 1, 0, 0.5); -} -inline Nucleus* DeuteriumNucleus() -{ - return new Nucleus("D", 1, 1, 1.0); -} -inline Nucleus* TritiumNucleus() -{ - return new Nucleus("T", 1, 2, 0.5); -} -inline Nucleus* He3Nucleus() -{ - return new Nucleus("He3", 2, 1, 0.5); -} -inline Nucleus* He4Nucleus() -{ - return new Nucleus("He3", 2, 2, 0.0); -} -inline Nucleus* C12nucleus() -{ - return new Nucleus("C12", 6, 6, 0.0); -} -inline Nucleus* C13nucleus() -{ - return new Nucleus("C13", 6, 7, 0.5); -} -inline Nucleus* N14nucleus() -{ - return new Nucleus("N14", 7, 7, 1.0); -} -inline Nucleus* N15nucleus() -{ - return new Nucleus("N15", 7, 8, 0.5); -} -inline Nucleus* O16nucleus() -{ - return new Nucleus("O16", 8, 8, 0.0); -} -inline Nucleus* O17nucleus() -{ - return new Nucleus("O17", 8, 9, 2.5); -} -inline Nucleus* O18nucleus() -{ - return new Nucleus("O18", 8, 10, 0.0); -} -inline Nucleus* F19nucleus() -{ - return new Nucleus("F19", 9, 10, 0.5); -} - - -} // Cantera - -#endif diff --git a/src/spectra/rotor.cpp b/src/spectra/rotor.cpp deleted file mode 100644 index 9ae5180dc..000000000 --- a/src/spectra/rotor.cpp +++ /dev/null @@ -1,109 +0,0 @@ -/** - * @file rotor.cpp - * -*/ -#include "cantera/base/ct_defs.h" -#include "cantera/spectra/rotor.h" - -namespace Cantera -{ - -/** - * Constructor. - * - * @param Bv Rotational constant, wavenumbers. - * @dipoleMoment permanent dipole moment. - * @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. - */ -Rotor::Rotor(doublereal Bv, doublereal dipoleMoment, - doublereal Dv, doublereal Hv) : m_Bv(Bv), - m_Dv(Dv), - m_Hv(Hv), - m_dipole(dipoleMoment) {} - -/** - * 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)); -} - -/** - * The number of quantum states with the same J. For a - * quantum-mechanical rotor, this is simply 2J+1. - */ -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) - * \f] - */ -doublereal Rotor::partitionFunction(doublereal T, int cutoff) -{ - int j = 0; - doublereal T_Trot = wnum_to_J(m_Bv)/(Boltzmann*T); - if (T_Trot > 100.0) { - return T_Trot; - } else { - if (cutoff < 0) { - cutoff = (int)(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; - } -} - -/** - * 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 frequency at which radiation is absorbed by a transition - * from the lower to the upper state in wavenumber units. - */ -doublereal Rotor::frequency(int J_lower, int J_upper) -{ - return energy_w(J_upper) - energy_w(J_lower); -} - -/** - * The spectral intensity of a rotational transition. - */ -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); -} - -} diff --git a/src/spectra/spectralUtilities.cpp b/src/spectra/spectralUtilities.cpp deleted file mode 100644 index 9f5e553b3..000000000 --- a/src/spectra/spectralUtilities.cpp +++ /dev/null @@ -1,12 +0,0 @@ -#include "spectralUtilities.h" -#include "Nuclei.h" - -namespace Cantera -{ - -//Nucleus_syms = {"H":1, "D":2, "T":3, "He3":4, "He4":5, -// "C12":6, "C13":7, "N14":8, "N15":9, -// "O16":10, "O17":11, "O18":12, "F19":13}; -} - - diff --git a/src/spectra/spectralUtilities.h b/src/spectra/spectralUtilities.h deleted file mode 100644 index 3b4bb4130..000000000 --- a/src/spectra/spectralUtilities.h +++ /dev/null @@ -1,6 +0,0 @@ -#ifndef CT_SPEC_UTILS_H -#define CT_SPEC_UTILS_H - -#include "Nuclei.h" - -#endif diff --git a/test_problems/SConscript b/test_problems/SConscript index e560bf784..8561ee708 100644 --- a/test_problems/SConscript +++ b/test_problems/SConscript @@ -281,7 +281,6 @@ CompileAndTest('pureFluid', 'pureFluidTest', 'testPureWater', 'output_blessed.tx if haveConverters: CompileAndTest('rankine_democxx', 'rankine_democxx', 'rankine', 'output_blessed.txt') CompileAndTest('silane_equil', 'silane_equil', 'silane_equi', 'output_blessed.txt') -# spectroscopy is incomplete CompileAndTest('simpleTransport', 'simpleTransport', 'simpleTransport', 'output_blessed.txt') CompileAndTest('stoichSolidKinetics', 'stoichSolidKinetics', diff --git a/test_problems/spectroscopy/spectratest.cpp b/test_problems/spectroscopy/spectratest.cpp deleted file mode 100644 index ef4a45a04..000000000 --- a/test_problems/spectroscopy/spectratest.cpp +++ /dev/null @@ -1,72 +0,0 @@ -#include "cantera/spectra.h" -#include "spectra/Nuclei.h" -#include -using namespace std; -using namespace Cantera; - -int main() -{ - - Nucleus* a = HydrogenNucleus(); - Nucleus* b = HydrogenNucleus(); - if (*a == *b) { - cout << "a and b and indistinguishable" << endl; - } else { - cout << "\nwhy are a and b not indistinguishable?\n"; - return 1; - } - - // test line broading classes - double gam = 2.0e0; - double sigma = 5.0; - - LineBroadener* lor = new LorentzianProfile(gam); - LineBroadener* gaus = new GaussianProfile(sigma); - Voigt* voig = new Voigt(sigma, gam); - //voig->testv(); - - double dnu = 0.1; - double nu; - double sum = 0.0, sumg = 0.0, sumlor = 0.0; - for (int n = -2000; n < 2000; n++) { - //cout << n << endl; - nu = n*dnu; - sumg += gaus->profile(nu)*dnu; - sum += voig->profile(nu)*dnu; - sumlor += lor->profile(nu)*dnu; - //cout << nu << ", " << (*lor)(nu) << ", " << (*gaus)(nu) - // << ", " << (*voig)(nu) << endl; - } - - /* old output - cout << "Voigt area = " << sum << endl; - cout << "Gaussian area = " << sumg << endl; - cout << "Lorentzian area = " << sumlor << endl; - */ - - // 'blessed' output: - // Voigt area = 0.99363 - // Gaussian area = 1 - // Lorentzian area = 0.993634 - - // guessing a sane tolerance - double TOL = .0001; - - if (abs(sum-.99363) > TOL) { - cout << "\nVOIGT AREA REGRESSION TEST FAILURE\n"; - return 1; - } - - if (abs(sumg-1.0) > TOL) { - cout << "\nGAUSSIAN AREA REGRESSION TEST FAILURE\n"; - return 1; - } - if (abs(sumlor-.993634) > TOL) { - cout << "\nLORENTZIAN AREA REGRESSION TEST FAILURE\n"; - return 1; - } - - // steady as she goes - return 0; - -}