Commit 4be82470 authored by eric pellegrini's avatar eric pellegrini
Browse files

fixe conflicts

parents 0d899064 d3439bc8
/* mconf.h
*
* Common include file for math routines
*
*
*
* SYNOPSIS:
*
* #include "mconf.h"
*
*
*
* DESCRIPTION:
*
* This file contains definitions for error codes that are
* passed to the common error handling routine mtherr()
* (which see).
*
* The file also includes a conditional assembly definition
* for the type of computer arithmetic (IEEE, DEC, Motorola
* IEEE, or UNKnown).
*
* For Digital Equipment PDP-11 and VAX computers, certain
* IBM systems, and others that use numbers with a 56-bit
* significand, the symbol DEC should be defined. In this
* mode, most floating point constants are given as arrays
* of octal integers to eliminate decimal to binary conversion
* errors that might be introduced by the compiler.
*
* For little-endian computers, such as IBM PC, that follow the
* IEEE Standard for Binary Floating Point Arithmetic (ANSI/IEEE
* Std 754-1985), the symbol IBMPC should be defined. These
* numbers have 53-bit significands. In this mode, constants
* are provided as arrays of hexadecimal 16 bit integers.
*
* Big-endian IEEE format is denoted MIEEE. On some RISC
* systems such as Sun SPARC, double precision constants
* must be stored on 8-byte address boundaries. Since integer
* arrays may be aligned differently, the MIEEE configuration
* may fail on such machines.
*
* To accommodate other types of computer arithmetic, all
* constants are also provided in a normal decimal radix
* which one can hope are correctly converted to a suitable
* format by the available C language compiler. To invoke
* this mode, define the symbol UNK.
*
* An important difference among these modes is a predefined
* set of machine arithmetic constants for each. The numbers
* MACHEP (the machine roundoff error), MAXNUM (largest number
* represented), and several other parameters are preset by
* the configuration symbol. Check the file const.c to
* ensure that these values are correct for your computer.
*
* Configurations NANS, INFINITIES, MINUSZERO, and DENORMAL
* may fail on many systems. Verify that they are supposed
* to work on your computer.
*/
/*
Cephes Math Library Release 2.3: June, 1995
Copyright 1984, 1987, 1989, 1995 by Stephen L. Moshier
*/
#ifndef _MCONF_H
#define _MCONF_H
/* Constant definitions for math error conditions
*/
#define DOMAIN 1 /* argument domain error */
#define SING 2 /* argument singularity */
#define OVERFLOW 3 /* overflow range error */
#define UNDERFLOW 4 /* underflow range error */
#define TLOSS 5 /* total loss of precision */
#define PLOSS 6 /* partial loss of precision */
#define EDOM 33
#define ERANGE 34
/* Complex numeral. */
typedef struct
{
double r;
double i;
} cmplx;
/* Long double complex numeral. */
/*
typedef struct
{
long double r;
long double i;
} cmplxl;
*/
/* Type of computer arithmetic */
/* PDP-11, Pro350, VAX:
*/
#undef DEC
/* Intel IEEE, low order words come first:
*/
#undef IBMPC
/* Motorola IEEE, high order words come first
* (Sun 680x0 workstation):
*/
#undef MIEEE
/* UNKnown arithmetic, invokes coefficients given in
* normal decimal format. Beware of range boundary
* problems (MACHEP, MAXLOG, etc. in const.c) and
* roundoff problems in pow.c:
* (Sun SPARCstation)
*/
#define UNK 1
/* If you define UNK, then be sure to set BIGENDIAN properly. */
#define BIGENDIAN 0
/* Define this `volatile' if your compiler thinks
* that floating point arithmetic obeys the associative
* and distributive laws. It will defeat some optimizations
* (but probably not enough of them).
*
* #define VOLATILE volatile
*/
#define VOLATILE
/* For 12-byte long doubles on an i386, pad a 16-bit short 0
* to the end of real constants initialized by integer arrays.
*
* #define XPD 0,
*
* Otherwise, the type is 10 bytes long and XPD should be
* defined blank (e.g., Microsoft C).
*
* #define XPD
*/
#define XPD 0,
/* Define to support tiny denormal numbers, else undefine. */
#define DENORMAL 1
/* Define to ask for infinity support, else undefine. */
#define INFINITIES 1
/* Define to ask for support of numbers that are Not-a-Number,
else undefine. This may automatically define INFINITIES in some files. */
#define NANS 1
/* Define to distinguish between -0.0 and +0.0. */
#define MINUSZERO 1
/* Define 1 for ANSI C atan2() function
See atan.c and clog.c. */
#define ANSIC 1
/* ANSI function prototypes */
#define ANSIPROT
//extern double polevl ( double x, double *P, int N );
//extern double p1evl ( double x, double *P, int N );
//extern double ndtr ( double a );
//extern double erfc ( double a );
//extern double erf ( double x );
/* No error reporting */
#define mtherr(function, type)
/* Constants */
#define MAXLOG 8.8029691931113054295988E1
#define SQRTH 7.07106781186547524401E-1
#endif
/* ndtr.c
*
* Normal distribution function
*
*
*
* SYNOPSIS:
*
* double x, y, ndtr();
*
* y = ndtr( x );
*
*
*
* DESCRIPTION:
*
* Returns the area under the Gaussian probability density
* function, integrated from minus infinity to x:
*
* x
* -
* 1 | | 2
* ndtr(x) = --------- | exp( - t /2 ) dt
* sqrt(2pi) | |
* -
* -inf.
*
* = ( 1 + erf(z) ) / 2
* = erfc(z) / 2
*
* where z = x/sqrt(2). Computation is via the functions
* erf and erfc.
*
*
* ACCURACY:
*
* Relative error:
* arithmetic domain # trials peak rms
* DEC -13,0 8000 2.1e-15 4.8e-16
* IEEE -13,0 30000 3.4e-14 6.7e-15
*
*
* ERROR MESSAGES:
*
* message condition value returned
* erfc underflow x > 37.519379347 0.0
*
*/
/* erf.c
*
* Error function
*
*
*
* SYNOPSIS:
*
* double x, y, erf();
*
* y = erf( x );
*
*
*
* DESCRIPTION:
*
* The integral is
*
* x
* -
* 2 | | 2
* erf(x) = -------- | exp( - t ) dt.
* sqrt(pi) | |
* -
* 0
*
* The magnitude of x is limited to 9.231948545 for DEC
* arithmetic; 1 or -1 is returned outside this range.
*
* For 0 <= |x| < 1, erf(x) = x * P4(x**2)/Q5(x**2); otherwise
* erf(x) = 1 - erfc(x).
*
*
*
* ACCURACY:
*
* Relative error:
* arithmetic domain # trials peak rms
* DEC 0,1 14000 4.7e-17 1.5e-17
* IEEE 0,1 30000 3.7e-16 1.0e-16
*
*/
/* erfc.c
*
* Complementary error function
*
*
*
* SYNOPSIS:
*
* double x, y, erfc();
*
* y = erfc( x );
*
*
*
* DESCRIPTION:
*
*
* 1 - erf(x) =
*
* inf.
* -
* 2 | | 2
* erfc(x) = -------- | exp( - t ) dt
* sqrt(pi) | |
* -
* x
*
*
* For small x, erfc(x) = 1 - erf(x); otherwise rational
* approximations are computed.
*
*
*
* ACCURACY:
*
* Relative error:
* arithmetic domain # trials peak rms
* DEC 0, 9.2319 12000 5.1e-16 1.2e-16
* IEEE 0,26.6417 30000 5.7e-14 1.5e-14
*
*
* ERROR MESSAGES:
*
* message condition value returned
* erfc underflow x > 9.231948545 (DEC) 0.0
*
*
*/
/*
Cephes Math Library Release 2.2: June, 1992
Copyright 1984, 1987, 1988, 1992 by Stephen L. Moshier
Direct inquiries to 30 Frost Street, Cambridge, MA 02140
*/
extern double ndtr(double a);
extern double erfc(double a);
extern double erf(double x);
/* polevl.c
* p1evl.c
*
* Evaluate polynomial
*
*
*
* SYNOPSIS:
*
* int N;
* double x, y, coef[N+1], polevl[];
*
* y = polevl( x, coef, N );
*
*
*
* DESCRIPTION:
*
* Evaluates polynomial of degree N:
*
* 2 N
* y = C + C x + C x +...+ C x
* 0 1 2 N
*
* Coefficients are stored in reverse order:
*
* coef[0] = C , ..., coef[N] = C .
* N 0
*
* The function p1evl() assumes that coef[N] = 1.0 and is
* omitted from the array. Its calling arguments are
* otherwise the same as polevl().
*
*
* SPEED:
*
* In the interest of speed, there are no checks for out
* of bounds arithmetic. This routine is used by most of
* the functions in the library. Depending on available
* equipment features, the user may wish to rewrite the
* program in microcode or assembly language.
*
*/
/*
Cephes Math Library Release 2.1: December, 1988
Copyright 1984, 1987, 1988 by Stephen L. Moshier
Direct inquiries to 30 Frost Street, Cambridge, MA 02140
*/
extern double polevl( double x, double coef[], int N );
extern double p1evl( double x, double coef[], int N );
......@@ -131,7 +131,7 @@ readDCD(PyObject *dummy, PyObject *args)
goto error;
}
if( atoms != dcdAtoms ){
snprintf(buffer, sizeof(buffer),
PyOS_snprintf(buffer, sizeof(buffer),
"number of atoms in DCD file (%d) doesn't "
"match universe (%d)", dcdAtoms, atoms);
PyErr_SetString(PyExc_ValueError, buffer);
......
......@@ -29,7 +29,7 @@ add_point(const TPoint p, TPoint *r, int r_index, PyObject *pt_dict)
{
PyObject *p_orig;
char key[512];
snprintf(key, sizeof(key), "%.6f,%.6f,%.6f", p[0], p[1], p[2]);
PyOS_snprintf(key, sizeof(key), "%.6f,%.6f,%.6f", p[0], p[1], p[2]);
p_orig = PyDict_GetItemString(pt_dict, key);
if(p_orig)
Py_DECREF(p_orig);
......@@ -178,7 +178,7 @@ nbor_data_1_atom(PyObject *nbors, int i, PyObject *atom_data,
{
char key[128];
PyObject *alist;
snprintf(key, sizeof(key), "%d %d %d", key0+nbor_boxes[boxn][0],
PyOS_snprintf(key, sizeof(key), "%d %d %d", key0+nbor_boxes[boxn][0],
key1+nbor_boxes[boxn][1], key2+nbor_boxes[boxn][2]);
alist = PyDict_GetItemString(boxes, key);
if(!alist && i == -1) printf("none in list at %s, atom %d\n", key, i);
......@@ -421,7 +421,7 @@ FindNeighbors(PyObject *dummy, PyObject *args)
PyObject *v, *tmp;
char key[128];
PyObject *pos = PyList_GetItem(atom_data, (Py_ssize_t)i);
snprintf(key, sizeof(key), "%d %d %d",
PyOS_snprintf(key, sizeof(key), "%d %d %d",
(int)floor(PyFloat_AsDouble(PyTuple_GetItem(pos,
(Py_ssize_t)0))/box_size),
(int)floor(PyFloat_AsDouble(PyTuple_GetItem(pos,
......@@ -484,7 +484,7 @@ FindNeighbors(PyObject *dummy, PyObject *args)
{
char key[128];
PyObject *alist;
snprintf(key, sizeof(key), "%d %d %d", key0+nbor_boxes[boxn][0],
PyOS_snprintf(key, sizeof(key), "%d %d %d", key0+nbor_boxes[boxn][0],
key1+nbor_boxes[boxn][1], key2+nbor_boxes[boxn][2]);
alist = PyDict_GetItemString(boxes, key);
if(!alist && i == -1) printf("none in list at %s\n", key);
......@@ -590,7 +590,7 @@ FindNeighborsOfAtom(PyObject *dummy, PyObject *args)
{
char key[128];
PyObject *alist;
snprintf(key, sizeof(key),
PyOS_snprintf(key, sizeof(key),
"%d %d %d", key0+nbor_boxes[boxn][0],
key1+nbor_boxes[boxn][1], key2+nbor_boxes[boxn][2]);
alist = PyDict_GetItemString(boxes, key);
......
......@@ -769,7 +769,7 @@ PyTrajectory_TimeStamp(PyTrajectoryObject *self, char *text)
{
time_t now = time(NULL);
static char time_stamp[200];
snprintf(time_stamp, sizeof(time_stamp), text, ctime(&now));
PyOS_snprintf(time_stamp, sizeof(time_stamp), text, ctime(&now));
time_stamp[strlen(time_stamp)-1] = '\0';
return PyNetCDFFile_AddHistoryLine(self->file, time_stamp);
}
......@@ -1397,14 +1397,14 @@ PyTrajectory_Output(PyTrajectoryOutputSpec *spec, int step,
PyEval_RestoreThread(*thread);
if (step > 0)
PyFile_WriteString("\n", spec->destination);
snprintf(buffer, sizeof(buffer), "Step %d\n", step);
PyOS_snprintf(buffer, sizeof(buffer), "Step %d\n", step);
PyFile_WriteString(buffer, spec->destination);
for (var = data; var->name != NULL; var++)
if (spec->what & var->class) {
switch (var->type) {
case PyTrajectory_Scalar:
snprintf(buffer, sizeof(buffer), var->text, *var->value.dp);
PyOS_snprintf(buffer, sizeof(buffer), var->text, *var->value.dp);
PyFile_WriteString(buffer, spec->destination);
break;
case PyTrajectory_ParticleScalar:
......@@ -1412,7 +1412,7 @@ PyTrajectory_Output(PyTrajectoryOutputSpec *spec, int step,
PyFile_WriteString(buffer, spec->destination);
array_data = (double *)var->value.array->data;
for (i = 0; i < var->value.array->dimensions[0]; i++) {
snprintf(buffer, sizeof(buffer),
PyOS_snprintf(buffer, sizeof(buffer),
" %5d: %g\n", i, *array_data++);
PyFile_WriteString(buffer, spec->destination);
}
......@@ -1422,7 +1422,7 @@ PyTrajectory_Output(PyTrajectoryOutputSpec *spec, int step,
PyFile_WriteString(buffer, spec->destination);
array_data = (double *)var->value.array->data;
for (i = 0; i < var->value.array->dimensions[0]; i++) {
snprintf(buffer, sizeof(buffer)," %5d: %g, %g, %g\n", i,
PyOS_snprintf(buffer, sizeof(buffer)," %5d: %g, %g, %g\n", i,
array_data[0], array_data[1], array_data[2]);
PyFile_WriteString(buffer, spec->destination);
array_data += 3;
......@@ -1433,12 +1433,12 @@ PyTrajectory_Output(PyTrajectoryOutputSpec *spec, int step,
PyFile_WriteString(buffer, spec->destination);
array_data = var->value.dp;
for (i = 0; i < var->length; i++) {
snprintf(buffer, sizeof(buffer), " %g", *array_data++);
PyOS_snprintf(buffer, sizeof(buffer), " %g", *array_data++);
PyFile_WriteString(buffer, spec->destination);
if (i == var->length-1)
snprintf(buffer, sizeof(buffer), "\n");
PyOS_snprintf(buffer, sizeof(buffer), "\n");
else
snprintf(buffer, sizeof(buffer), ",");
PyOS_snprintf(buffer, sizeof(buffer), ",");
PyFile_WriteString(buffer, spec->destination);
}
break;
......
......@@ -607,7 +607,7 @@ int write_dcdheader(FILE *fd, char *filename, int N, int NSET, int ISTART,
out_integer = 2;
fwrite((char *) & out_integer, sizeof(int), 1, fd);
snprintf(title_string, sizeof(title_string),
PyOS_snprintf(title_string, sizeof(title_string),
"REMARKS FILENAME=%s CREATED BY VMD", filename);
pad(title_string, 80);
fwrite(title_string, sizeof(char), 80, fd);
......@@ -616,7 +616,7 @@ int write_dcdheader(FILE *fd, char *filename, int N, int NSET, int ISTART,
tmbuf=localtime(&cur_time);
strftime(time_str, 10, "%m/%d/%y", tmbuf);
snprintf(title_string, sizeof(title_string),
PyOS_snprintf(title_string, sizeof(title_string),
"REMARKS DATE: %s CREATED BY MMTK.", time_str);
pad(title_string, 80);
fwrite(title_string, sizeof(char), 80, fd);
......
......@@ -16,6 +16,12 @@
/*
* Include erfc code, unless user wants to use an erfc from libm.
*/
<<<<<<< HEAD
=======
#ifndef LIBM_HAS_ERFC
#include "MMTK/ndtr.h"
#endif
>>>>>>> feature-remove-old-numeric-2
#ifndef M_PI
#define M_PI 3.14159265358979323846264338327950288
#endif
......
/* ndtr.c
*
* Normal distribution function
*
*
*
* SYNOPSIS:
*
* double x, y, ndtr();
*
* y = ndtr( x );
*
*
*
* DESCRIPTION:
*
* Returns the area under the Gaussian probability density
* function, integrated from minus infinity to x:
*
* x
* -
* 1 | | 2
* ndtr(x) = --------- | exp( - t /2 ) dt
* sqrt(2pi) | |
* -
* -inf.
*
* = ( 1 + erf(z) ) / 2
* = erfc(z) / 2
*
* where z = x/sqrt(2). Computation is via the functions
* erf and erfc.
*
*
* ACCURACY:
*
* Relative error:
* arithmetic domain # trials peak rms
* DEC -13,0 8000 2.1e-15 4.8e-16
* IEEE -13,0 30000 3.4e-14 6.7e-15
*
*
* ERROR MESSAGES:
*
* message condition value returned
* erfc underflow x > 37.519379347 0.0
*
*/
/* erf.c
*
* Error function
*
*
*
* SYNOPSIS:
*
* double x, y, erf();
*
* y = erf( x );
*
*
*
* DESCRIPTION:
*
* The integral is
*
* x
* -
* 2 | | 2
* erf(x) = -------- | exp( - t ) dt.
* sqrt(pi) | |
* -
* 0
*
* The magnitude of x is limited to 9.231948545 for DEC
* arithmetic; 1 or -1 is returned outside this range.
*
* For 0 <= |x| < 1, erf(x) = x * P4(x**2)/Q5(x**2); otherwise
* erf(x) = 1 - erfc(x).
*
*
*
* ACCURACY:
*
* Relative error:
* arithmetic domain # trials peak rms
* DEC 0,1 14000 4.7e-17 1.5e-17
* IEEE 0,1 30000 3.7e-16 1.0e-16
*
*/
/* erfc.c
*
* Complementary error function
*
*
*
* SYNOPSIS:
*
* double x, y, erfc();
*
* y = erfc( x );
*
*
*
* DESCRIPTION:
*
*
* 1 - erf(x) =
*
* inf.
* -
* 2 | | 2
* erfc(x) = -------- | exp( - t ) dt
* sqrt(pi) | |