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

Moved the "checkFinite" function into the Cantera namespace

Browse source
This commit is contained in:
Ray Speth 2013-04-29 18:05:22 +00:00
parent 1d40085abc
commit 2c58b7237a
5 changed files with 26 additions and 157 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

40
include/cantera/base/mdp_allo.h
View file

@ -687,46 +687,6 @@ extern void mdp_init_dbl_2(double** const v, const double value,
*/
extern void mdp_init_int_1(int* const v, const int value, const int len);
/*
* Utility routines to check that a number is finite
*/
//! Utility routine to check to see that a number is neither zero
//! nor indefinite.
/*!
* This check can be used before using the number in a denominator.
*
* @param tmp number to be checked
*/
extern void checkZeroFinite(const double tmp);
//! Utility routine to check to see that a number is finite.
/*!
* @param tmp number to be checked
*/
extern void checkFinite(const double tmp);
//! Utility routine to link checkFinte() to fortran program
/*!
* This routine is accessible from fortran, usually
*
* @param tmp Pointer to the number to check
*
* @todo link it into the usual way Cantera handles Fortran calls
*/
extern "C" void checkfinite_(double* tmp);
//! utility routine to check that a double stays bounded
/*!
* This routine checks to see if a number stays bounded. The absolute
* value of the number is required to stay below the trigger.
*
* @param tmp Number to be checked
* @param trigger bounds on the number. Defaults to 1.0E20
*/
extern void checkMagnitude(const double tmp, const double trigger = 1.0E20);
} /* end of mdp namespace */
/****************************************************************************/
#endif

3
include/cantera/base/utilities.h
View file

@ -678,6 +678,9 @@ void deepStdVectorPointerCopy(const std::vector<D*> &fromVec, std::vector<D*> &t
}
//@}
//! Check to see that a number is finite (not NaN, +Inf or -Inf)
void checkFinite(const double tmp);
}
#endif

110
src/base/checkFinite.cpp
View file

@ -1,7 +1,6 @@
/**
* @file checkFinite.cpp
* Declarations for Multi Dimensional Pointer (mdp) routines that
* check for the presence of NaNs in the code.
* @file checkFinite.cpp Declarations for routines that check for the
* presence of NaNs in the code.
*/
/*
* Copyright 2004 Sandia Corporation. Under the terms of Contract
@ -13,7 +12,6 @@
#include "cantera/base/ct_defs.h"
#include <stdexcept>
#include <string>
#include <cmath>
#include <cstdlib>
@ -33,16 +31,11 @@
using namespace std;
namespace mdp
{
namespace Cantera {
// Utility routine to check to see that a number is finite.
/*
* @param tmp number to be checked
*/
#ifdef _WIN32
void checkFinite(const double tmp)
{
#ifdef _WIN32
if (_finite(tmp)) {
if (_isnan(tmp)) {
printf("checkFinite() ERROR: we have encountered a nan!\n");
@ -51,14 +44,10 @@ void checkFinite(const double tmp)
} else {
printf("checkFinite() ERROR: we have encountered a neg inf!\n");
}
const std::string s = "checkFinite()";
throw std::range_error(s);
throw std::range_error("checkFinite()");
}
}
#else
void checkFinite(const double tmp)
{
if (! finite(tmp)) {
if (!finite(tmp)) {
if (isnan(tmp)) {
printf("checkFinite() ERROR: we have encountered a nan!\n");
} else if (isinf(tmp) == 1) {
@ -66,90 +55,9 @@ void checkFinite(const double tmp)
} else {
printf("checkFinite() ERROR: we have encountered a neg inf!\n");
}
const std::string s = "checkFinite()";
throw std::range_error(s);
throw std::range_error("checkFinite()");
}
}
#endif
// Utility routine to link checkFinte() to fortran program
/*
* This routine is accessible from fortran, usually
*
* @param tmp Pointer to the number to check
*
* @todo link it into the usual way Cantera handles Fortran calls
*/
extern "C" void checkfinite_(double* tmp)
{
checkFinite(*tmp);
}
// Utility routine to check that a double stays bounded
/*
* This routine checks to see if a number stays bounded. The absolute
* value of the number is required to stay below the trigger.
*
* @param tmp Number to be checked
* @param trigger bounds on the number. Defaults to 1.0E20
*/
void checkMagnitude(const double tmp, const double trigger)
{
checkFinite(tmp);
if (fabs(tmp) >= trigger) {
char sbuf[64];
sprintf(sbuf, "checkMagnitude() ERROR: Trigger %g exceeded: %g\n", trigger,
tmp);
throw std::range_error(sbuf);
}
}
// Utility routine to check to see that a number is neither zero
// nor indefinite.
/*
* This check can be used before using the number in a denominator.
*
* @param tmp number to be checked
*/
#ifdef _WIN32
void checkZeroFinite(const double tmp)
{
if ((tmp == 0.0) || (! _finite(tmp))) {
if (tmp == 0.0) {
printf("checkZeroFinite() ERROR: we have encountered a zero!\n");
} else if (_isnan(tmp)) {
printf("checkZeroFinite() ERROR: we have encountered a nan!\n");
} else if (_fpclass(tmp) == _FPCLASS_PINF) {
printf("checkZeroFinite() ERROR: we have encountered a pos inf!\n");
} else {
printf("checkZeroFinite() ERROR: we have encountered a neg inf!\n");
}
char sbuf[64];
sprintf(sbuf, "checkZeroFinite() ERROR: zero or indef exceeded: %g\n",
tmp);
throw std::range_error(sbuf);
}
}
#else
void checkZeroFinite(const double tmp)
{
if ((tmp == 0.0) || (! finite(tmp))) {
if (tmp == 0.0) {
printf("checkZeroFinite() ERROR: we have encountered a zero!\n");
} else if (isnan(tmp)) {
printf("checkZeroFinite() ERROR: we have encountered a nan!\n");
} else if (isinf(tmp) == 1) {
printf("checkZeroFinite() ERROR: we have encountered a pos inf!\n");
} else {
printf("checkZeroFinite() ERROR: we have encountered a neg inf!\n");
}
char sbuf[64];
sprintf(sbuf, "checkZeroFinite() ERROR: zero or indef exceeded: %g\n",
tmp);
throw std::range_error(sbuf);
}
}
#endif
}
}

22
src/numerics/NonlinearSolver.cpp
View file

@ -574,17 +574,17 @@ doublereal NonlinearSolver::residErrorNorm(const doublereal* const resid, const
for (size_t i = 0; i < neq_; i++) {
#ifdef DEBUG_MODE
mdp::checkFinite(resid[i]);
checkFinite(resid[i]);
#endif
error = resid[i] / m_residWts[i];
#ifdef DEBUG_MODE
mdp::checkFinite(error);
checkFinite(error);
#endif
sum_norm += (error * error);
}
sum_norm = sqrt(sum_norm / neq_);
#ifdef DEBUG_MODE
mdp::checkFinite(sum_norm);
checkFinite(sum_norm);
#endif
if (printLargest) {
const int num_entries = printLargest;
@ -820,7 +820,7 @@ void NonlinearSolver::scaleMatrix(GeneralMatrix& jac, doublereal* const y_comm,
m_rowWtScales[irow] += fabs(*jptr) * m_ewt[jcol];
}
#ifdef DEBUG_MODE
mdp::checkFinite(m_rowWtScales[irow]);
checkFinite(m_rowWtScales[irow]);
#endif
jptr++;
}
@ -844,7 +844,7 @@ void NonlinearSolver::scaleMatrix(GeneralMatrix& jac, doublereal* const y_comm,
m_rowWtScales[irow] += vv * m_ewt[jcol];
}
#ifdef DEBUG_MODE
mdp::checkFinite(m_rowWtScales[irow]);
checkFinite(m_rowWtScales[irow]);
#endif
}
}
@ -1441,7 +1441,7 @@ doublereal NonlinearSolver::doCauchyPointSolve(GeneralMatrix& jac)
deltaX_CP_[j] -= m_resid[i] * jac(i,j) * colFac * rowFac * m_ewt[j] * m_ewt[j]
/ (m_residWts[i] * m_residWts[i]);
#ifdef DEBUG_MODE
mdp::checkFinite(deltaX_CP_[j]);
checkFinite(deltaX_CP_[j]);
#endif
}
}
@ -3773,7 +3773,7 @@ int NonlinearSolver::beuler_jac(GeneralMatrix& J, doublereal* const f,
m_nJacEval++;
#ifdef DEBUG_MODE
for (int ii = 0; ii < neq_; ii++) {
mdp::checkFinite(f[ii]);
checkFinite(f[ii]);
}
#endif
@ -3844,7 +3844,7 @@ int NonlinearSolver::beuler_jac(GeneralMatrix& J, doublereal* const f,
throw CanteraError("NonlinearSolver::beuler_jac", "dy is equal to zero");
}
for (int ii = 0; ii < neq_; ii++) {
mdp::checkFinite(m_wksp[ii]);
checkFinite(m_wksp[ii]);
}
#endif
@ -3929,7 +3929,7 @@ int NonlinearSolver::beuler_jac(GeneralMatrix& J, doublereal* const f,
throw CanteraError("NonlinearSolver::beuler_jac", "dy is equal to zero");
}
for (int ii = 0; ii < neq_; ii++) {
mdp::checkFinite(m_wksp[ii]);
checkFinite(m_wksp[ii]);
}
#endif
if (info != 1) {
@ -4094,7 +4094,7 @@ NonlinearSolver::computeResidWts()
for (size_t i = 0; i < neq_; i++) {
m_residWts[i] = userResidAtol_[i] + userResidRtol_ * m_rowWtScales[i] / rtol_;
#ifdef DEBUG_MODE
mdp::checkFinite(m_residWts[i]);
checkFinite(m_residWts[i]);
#endif
}
} else {
@ -4102,7 +4102,7 @@ NonlinearSolver::computeResidWts()
for (size_t i = 0; i < neq_; i++) {
m_residWts[i] = m_rowWtScales[i];
#ifdef DEBUG_MODE
mdp::checkFinite(m_residWts[i]);
checkFinite(m_residWts[i]);
#endif
sum += m_residWts[i];
}

8
src/numerics/RootFind.cpp
View file

@ -18,9 +18,7 @@
#endif
#include "cantera/base/global.h"
#ifdef DEBUG_MODE
#include "cantera/base/mdp_allo.h"
#endif
#include "cantera/base/utilities.h"
#include "cantera/base/stringUtils.h"
/* Standard include files */
@ -1183,11 +1181,11 @@ doublereal RootFind::func(doublereal x)
{
doublereal r;
#ifdef DEBUG_MODE
mdp::checkFinite(x);
checkFinite(x);
#endif
m_residFunc->evalSS(0.0, &x, &r);
#ifdef DEBUG_MODE
mdp::checkFinite(r);
checkFinite(r);
#endif
doublereal ff = r - m_funcTargetValue;
if (x >= x_maxTried_) {