Skip to content

Commit

Permalink
Update of all the kinetics model tests.
Browse files Browse the repository at this point in the history 
Now all the tests are standardized, testing both KineticsConditions and
StateType for the temperature. All the vectorizations are also tested.
  • Loading branch information
@SylvainPlessis
SylvainPlessis committed Apr 15, 2015
1 parent 90d791c commit c3059ed
Show file tree
Hide file tree
Showing 13 changed files with 1,356 additions and 661 deletions.
94 changes: 50 additions & 44 deletions test/arrhenius_rate_unit.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -37,67 +37,73 @@
#include "antioch/units.h"

template <typename Scalar>
int test_values(const Scalar & Cf, const Scalar & Ea, const Scalar & R, const Antioch::ArrheniusRate<Scalar> & arrhenius_rate)
int check_rate_and_derivative(const Scalar & rate_exact, const Scalar & derive_exact,
const Scalar & rate, const Scalar & derive, const Scalar & T)
{
using std::abs;
using std::exp;
int return_flag = 0;

const Scalar tol = std::numeric_limits<Scalar>::epsilon() * 100;

for(Scalar T = 300.1; T <= 2500.1; T += 10.)
{
const Scalar rate_exact = Cf*exp(-Ea/(R*T));
const Scalar derive_exact = Ea/(R*T*T) * Cf * exp(-Ea/(R*T));

Scalar rate1 = arrhenius_rate(T);
Scalar deriveRate1 = arrhenius_rate.derivative(T);
Scalar rate;
Scalar deriveRate;

arrhenius_rate.rate_and_derivative(T,rate,deriveRate);

if( abs( (rate1 - rate_exact)/rate_exact ) > tol )
{
std::cout << std::scientific << std::setprecision(16)
<< "Error: Mismatch in rate values." << std::endl
<< "T = " << T << " K" << std::endl
<< "rate(T) = " << rate1 << std::endl
<< "rate_exact = " << rate_exact << std::endl;

return_flag = 1;
}
const Scalar tol = std::numeric_limits<Scalar>::epsilon() * 2;
int return_flag(0);
if( abs( (rate - rate_exact)/rate_exact ) > tol )
{
std::cout << std::scientific << std::setprecision(16)
<< "Error: Mismatch in rate values." << std::endl
<< "T = " << T << " K" << std::endl
<< "rate(T) = " << rate << std::endl
<< "rate_exact = " << rate_exact << std::endl;

<< "rate_exact = " << rate_exact << std::endl
<< "relative difference = " << abs( (rate - rate_exact)/rate_exact ) << std::endl
<< "tolerance = " << tol << std::endl;

return_flag = 1;
}
if( abs( (deriveRate1 - derive_exact)/derive_exact ) > tol )
if( abs( (derive - derive_exact)/derive_exact ) > tol )
{
std::cout << std::scientific << std::setprecision(16)
<< "Error: Mismatch in rate derivative values." << std::endl
<< "T = " << T << " K" << std::endl
<< "drate_dT(T) = " << deriveRate1 << std::endl
<< "derive_exact = " << derive_exact << std::endl;
<< "drate_dT(T) = " << derive << std::endl
<< "derive_exact = " << derive_exact << std::endl
<< "relative difference = " << abs( (derive - derive_exact)/derive_exact ) << std::endl
<< "tolerance = " << tol << std::endl;

return_flag = 1;
}
if( abs( (deriveRate - derive_exact)/derive_exact ) > tol )
{
std::cout << std::scientific << std::setprecision(16)
<< "Error: Mismatch in rate derivative values." << std::endl
<< "T = " << T << " K" << std::endl
<< "drate_dT(T) = " << deriveRate << std::endl
<< "derive_exact = " << derive_exact << std::endl;

return_flag = 1;
}
if(return_flag)break;
return return_flag;
}

template <typename Scalar>
int test_values(const Scalar & Cf, const Scalar & Ea, const Scalar & R, const Antioch::ArrheniusRate<Scalar> & arrhenius_rate)
{
using std::abs;
using std::exp;
int return_flag = 0;

for(Scalar T = 300.1; T <= 2500.1; T += 10.)
{
const Scalar rate_exact = Cf*exp(-Ea/(R*T));
const Scalar derive_exact = Ea/(R*T*T) * Cf * exp(-Ea/(R*T));
Antioch::KineticsConditions<Scalar> cond(T);

// KineticsConditions

Scalar rate = arrhenius_rate(cond);
Scalar deriveRate = arrhenius_rate.derivative(cond);

return_flag = check_rate_and_derivative(rate_exact,derive_exact,rate,deriveRate,T) || return_flag;

arrhenius_rate.rate_and_derivative(cond,rate,deriveRate);

return_flag = check_rate_and_derivative(rate_exact,derive_exact,rate,deriveRate,T) || return_flag;

// T

rate = arrhenius_rate(T);
deriveRate = arrhenius_rate.derivative(T);

return_flag = check_rate_and_derivative(rate_exact,derive_exact,rate,deriveRate,T) || return_flag;

arrhenius_rate.rate_and_derivative(T,rate,deriveRate);

return_flag = check_rate_and_derivative(rate_exact,derive_exact,rate,deriveRate,T) || return_flag;
}

return return_flag;
Expand Down
161 changes: 110 additions & 51 deletions test/arrhenius_rate_vec_unit.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -68,6 +68,69 @@ GRVY::GRVY_Timer_Class gt;
#include <cmath>
#include <limits>

template <typename PairScalars>
int check_rate_and_derivative(const PairScalars & rate_exact, const PairScalars & derive_exact,
const PairScalars & rate, const PairScalars & derive, const PairScalars & T)
{
typedef typename Antioch::value_type<PairScalars>::type Scalar;
const Scalar tol = std::numeric_limits<Scalar>::epsilon() * 2;

int return_flag(0);
for (unsigned int tuple=0; tuple != ANTIOCH_N_TUPLES; ++tuple)
{
if( abs( (rate[2*tuple] - rate_exact[2*tuple])/rate_exact[2*tuple] ) > tol )
{
std::cout << std::scientific << std::setprecision(16)
<< "Error: Mismatch in rate values." << std::endl
<< "T = " << T << " K" << std::endl
<< "rate(T) = " << rate[2*tuple] << std::endl
<< "rate_exact = " << rate_exact[2*tuple] << std::endl
<< "relative difference = " << abs( (rate[2*tuple] - rate_exact[2*tuple])/rate_exact[2*tuple] ) << std::endl
<< "tolerance = " << tol << std::endl;

return_flag = 1;
}
if( abs( (rate[2*tuple+1] - rate_exact[2*tuple+1])/rate_exact[2*tuple+1] ) > tol )
{
std::cout << std::scientific << std::setprecision(16)
<< "Error: Mismatch in rate values." << std::endl
<< "T = " << T << " K" << std::endl
<< "rate(T) = " << rate[2*tuple] << std::endl
<< "rate_exact = " << rate_exact[2*tuple+1] << std::endl
<< "relative difference = " << abs( (rate[2*tuple] - rate_exact[2*tuple+1])/rate_exact[2*tuple+1] ) << std::endl
<< "tolerance = " << tol << std::endl;

return_flag = 1;
}
if( abs( (derive[2*tuple] - derive_exact[2*tuple])/derive_exact[2*tuple] ) > tol )
{
std::cout << std::scientific << std::setprecision(16)
<< "Error: Mismatch in rate derivative values." << std::endl
<< "T = " << T << " K" << std::endl
<< "drate_dT(T) = " << derive[2*tuple] << std::endl
<< "derive_exact = " << derive_exact[2*tuple] << std::endl
<< "relative difference = " << abs( (derive[2*tuple] - derive_exact[2*tuple])/derive_exact[2*tuple] ) << std::endl
<< "tolerance = " << tol << std::endl;

return_flag = 1;
}
if( abs( (derive[2*tuple+1] - derive_exact[2*tuple+1])/derive_exact[2*tuple+1] ) > tol )
{
std::cout << std::scientific << std::setprecision(16)
<< "Error: Mismatch in rate derivative values." << std::endl
<< "T = " << T << " K" << std::endl
<< "drate_dT(T) = " << derive[2*tuple+1] << std::endl
<< "derive_exact = " << derive_exact[2*tuple+1] << std::endl
<< "relative difference = " << abs( (derive[2*tuple+1] - derive_exact[2*tuple+1])/derive_exact[2*tuple+1] ) << std::endl
<< "tolerance = " << tol << std::endl;

return_flag = 1;
}

}
return return_flag;
}

template <typename PairScalars>
int vectester(const PairScalars& example, const std::string& testname)
{
Expand All @@ -83,74 +146,72 @@ int vectester(const PairScalars& example, const std::string& testname)

// Construct from example to avoid resizing issues
PairScalars T = example;
PairScalars rate_exact = example;
PairScalars derive_exact = example;
for (unsigned int tuple=0; tuple != ANTIOCH_N_TUPLES; ++tuple)
{
T[2*tuple] = 1500.1;
T[2*tuple+1] = 1600.1;
rate_exact[2*tuple] = Cf*exp(-Ea/1500.1);
rate_exact[2*tuple+1] = Cf*exp(-Ea/1600.1);
derive_exact[2*tuple] = Ea/(Scalar(1500.1)*Scalar(1500.1)) * Cf * exp(-Ea/Scalar(1500.1));
derive_exact[2*tuple+1] = Ea/(Scalar(1600.1)*Scalar(1600.1)) * Cf * exp(-Ea/Scalar(1600.1));
}

const Scalar rate_exact0 = Cf*exp(-Ea/1500.1);
const Scalar rate_exact1 = Cf*exp(-Ea/1600.1);
const Scalar derive_exact0 = Ea/(Scalar(1500.1)*Scalar(1500.1)) * Cf * exp(-Ea/Scalar(1500.1));
const Scalar derive_exact1 = Ea/(Scalar(1600.1)*Scalar(1600.1)) * Cf * exp(-Ea/Scalar(1600.1));
Antioch::KineticsConditions<PairScalars> cond(T);

int return_flag = 0;

// KineticsConditions
#ifdef ANTIOCH_HAVE_GRVY
gt.BeginTimer(testname);
#endif

const PairScalars rate = arrhenius_rate(T);
const PairScalars deriveRate = arrhenius_rate.derivative(T);
PairScalars rate = arrhenius_rate(cond);
PairScalars derive = arrhenius_rate.derivative(cond);

#ifdef ANTIOCH_HAVE_GRVY
gt.EndTimer(testname);
#endif

const Scalar tol = std::numeric_limits<Scalar>::epsilon()*10;
return_flag = check_rate_and_derivative(rate_exact, derive_exact, rate, derive,T) || return_flag;

for (unsigned int tuple=0; tuple != ANTIOCH_N_TUPLES; ++tuple)
{
if( abs( (rate[2*tuple] - rate_exact0)/rate_exact0 ) > tol )
{
std::cout << "Error: Mismatch in rate values." << std::endl
<< "rate(T0) = " << rate[2*tuple] << std::endl
<< "rate_exact = " << rate_exact0 << std::endl
<< "difference = " << rate[2*tuple] - rate_exact0 << std::endl;

return_flag = 1;
break;
}

if( abs( (rate[2*tuple+1] - rate_exact1)/rate_exact1 ) > tol )
{
std::cout << "Error: Mismatch in rate values." << std::endl
<< "rate(T1) = " << rate[2*tuple+1] << std::endl
<< "rate_exact = " << rate_exact1 << std::endl
<< "difference = " << rate[2*tuple+1] - rate_exact1 << std::endl;

return_flag = 1;
break;
}
}
if( abs( (deriveRate[0] - derive_exact0)/derive_exact0 ) > tol )
{
std::cout << std::scientific << std::setprecision(16)
<< "Error: Mismatch in rate derivative values." << std::endl
<< "drate_dT(T0) = " << deriveRate[0] << std::endl
<< "derive_exact = " << derive_exact0 << std::endl;
#ifdef ANTIOCH_HAVE_GRVY
gt.BeginTimer(testname);
#endif

return_flag = 1;
}
if( abs( (deriveRate[1] - derive_exact1)/derive_exact1 ) > tol )
{
std::cout << std::scientific << std::setprecision(16)
<< "Error: Mismatch in rate derivative values." << std::endl
<< "drate_dT(T1) = " << deriveRate[1] << std::endl
<< "derive_exact = " << derive_exact1 << std::endl;
arrhenius_rate.rate_and_derivative(cond,rate,derive);

return_flag = 1;
}
#ifdef ANTIOCH_HAVE_GRVY
gt.EndTimer(testname);
#endif

return_flag = check_rate_and_derivative(rate_exact, derive_exact, rate, derive,T) || return_flag;

// T
#ifdef ANTIOCH_HAVE_GRVY
gt.BeginTimer(testname);
#endif

rate = arrhenius_rate(T);
derive = arrhenius_rate.derivative(T);

#ifdef ANTIOCH_HAVE_GRVY
gt.EndTimer(testname);
#endif

return_flag = check_rate_and_derivative(rate_exact, derive_exact, rate, derive,T) || return_flag;

#ifdef ANTIOCH_HAVE_GRVY
gt.BeginTimer(testname);
#endif

arrhenius_rate.rate_and_derivative(T,rate,derive);

#ifdef ANTIOCH_HAVE_GRVY
gt.EndTimer(testname);
#endif

return_flag = check_rate_and_derivative(rate_exact, derive_exact, rate, derive,T) || return_flag;

std::cout << "Arrhenius rate: " << arrhenius_rate << std::endl;

Expand Down Expand Up @@ -185,13 +246,11 @@ int main()
vectester (MetaPhysicL::NumberArray<2*ANTIOCH_N_TUPLES, long double> (0), "NumberArray<ld>");
#endif
#ifdef ANTIOCH_HAVE_VEXCL
std::cout << "vexcl start" << std::endl;
vex::Context ctx_f (vex::Filter::All);
std::cout << "vexcl start" << std::endl;
if (!ctx_f.empty())
returnval = returnval ||
vectester (vex::vector<float> (ctx_f, 2*ANTIOCH_N_TUPLES), "vex::vector<float>");
std::cout << "vexcl float ok" << std::endl;

vex::Context ctx_d (vex::Filter::DoublePrecision);
if (!ctx_d.empty())
returnval = returnval ||
Expand Down
Loading

0 comments on commit c3059ed

Please sign in to comment.