[Numpy-svn] r5888 - trunk/numpy/core/src
numpy-svn at scipy.org
numpy-svn at scipy.org
Wed Oct 1 14:08:44 EDT 2008
Author: charris
Date: 2008-10-01 13:08:41 -0500 (Wed, 01 Oct 2008)
New Revision: 5888
Modified:
trunk/numpy/core/src/umathmodule.c.src
Log:
Cleanup ufunc loops.
At this point loops are separated into variable kinds, so there is a fair amount
of duplication. I will probably merge loops that look the same in a later
commit. There are no changes to current behavior of loops, this will also be
changed in later work to deal with nans and such.
Modified: trunk/numpy/core/src/umathmodule.c.src
===================================================================
--- trunk/numpy/core/src/umathmodule.c.src 2008-10-01 18:06:04 UTC (rev 5887)
+++ trunk/numpy/core/src/umathmodule.c.src 2008-10-01 18:08:41 UTC (rev 5888)
@@ -1,6 +1,10 @@
/* -*- c -*- */
/*
+ * vim:syntax=c
+ */
+
+/*
*****************************************************************************
** INCLUDES **
*****************************************************************************
@@ -13,24 +17,47 @@
#include "config.h"
#include <math.h>
+#ifndef M_PI
+#define M_PI 3.14159265358979323846264338328
+#endif
+
/*
*****************************************************************************
** BASIC MATH FUNCTIONS **
*****************************************************************************
*/
-/* A whole slew of basic math functions are provided originally
- by Konrad Hinsen. */
+float degreesf(float x) {
+ return x * (float)(180.0/M_PI);
+}
+double degrees(double x) {
+ return x * (180.0/M_PI);
+}
+longdouble degreesl(longdouble x) {
+ return x * (180.0L/M_PI);
+}
+float radiansf(float x) {
+ return x * (float)(M_PI/180.0);
+}
+double radians(double x) {
+ return x * (M_PI/180.0);
+}
+longdouble radiansl(longdouble x) {
+ return x * (M_PI/180.0L);
+}
+
+/*
+ * A whole slew of basic math functions are provided originally
+ * by Konrad Hinsen.
+ */
+
#if !defined(__STDC__) && !defined(_MSC_VER)
extern double fmod (double, double);
extern double frexp (double, int *);
extern double ldexp (double, int);
extern double modf (double, double *);
#endif
-#ifndef M_PI
-#define M_PI 3.14159265358979323846264338328
-#endif
#if defined(DISTUTILS_USE_SDK)
@@ -443,26 +470,6 @@
#define isfinitef(x) (!(isinff((x)) || isnanf((x))))
#define isfinitel(x) (!(isinfl((x)) || isnanl((x))))
-float degreesf(float x) {
- return x * (float)(180.0/M_PI);
-}
-double degrees(double x) {
- return x * (180.0/M_PI);
-}
-longdouble degreesl(longdouble x) {
- return x * (180.0L/M_PI);
-}
-
-float radiansf(float x) {
- return x * (float)(M_PI/180.0);
-}
-double radians(double x) {
- return x * (M_PI/180.0);
-}
-longdouble radiansl(longdouble x) {
- return x * (M_PI/180.0L);
-}
-
/* First, the C functions that do the real work */
/* if C99 extensions not available then define dummy functions that use the
@@ -624,7 +631,62 @@
}
#endif
+/*
+ *****************************************************************************
+ ** PYTHON OBJECT FUNCTIONS **
+ *****************************************************************************
+ */
+static PyObject *
+Py_square(PyObject *o)
+{
+ return PyNumber_Multiply(o, o);
+}
+
+static PyObject *
+Py_get_one(PyObject *o)
+{
+ return PyInt_FromLong(1);
+}
+
+static PyObject *
+Py_reciprocal(PyObject *o)
+{
+ PyObject *one = PyInt_FromLong(1);
+ PyObject *result;
+
+ if (!one) {
+ return NULL;
+ }
+ result = PyNumber_Divide(one, o);
+ Py_DECREF(one);
+ return result;
+}
+
+/**begin repeat
+ * #Kind = Max, Min#
+ * #OP = >=, <=#
+ */
+static PyObject *
+_npy_Object at Kind@(PyObject *i1, PyObject *i2)
+{
+ PyObject *result;
+ int cmp;
+
+ if (PyObject_Cmp(i1, i2, &cmp) < 0) {
+ return NULL;
+ }
+ if (cmp @OP@ 0) {
+ result = i1;
+ }
+ else {
+ result = i2;
+ }
+ Py_INCREF(result);
+ return result;
+}
+/**end repeat**/
+
/*
*****************************************************************************
** COMPLEX FUNCTIONS **
@@ -632,11 +694,12 @@
*/
-/* Don't pass structures between functions (only pointers) because how
- structures are passed is compiler dependent and could cause
- segfaults if ufuncobject.c is compiled with a different compiler
- than an extension that makes use of the UFUNC API
-*/
+/*
+ * Don't pass structures between functions (only pointers) because how
+ * structures are passed is compiler dependent and could cause
+ * segfaults if ufuncobject.c is compiled with a different compiler
+ * than an extension that makes use of the UFUNC API
+ */
/**begin repeat
@@ -650,10 +713,11 @@
static c at typ@ nc_i at c@ = {0., 1.};
static c at typ@ nc_i2 at c@ = {0., 0.5};
/*
- static c at typ@ nc_mi at c@ = {0., -1.};
- static c at typ@ nc_pi2 at c@ = {M_PI/2., 0.};
-*/
+ * static c at typ@ nc_mi at c@ = {0., -1.};
+ * static c at typ@ nc_pi2 at c@ = {M_PI/2., 0.};
+ */
+
static void
nc_sum at c@(c at typ@ *a, c at typ@ *b, c at typ@ *r)
{
@@ -681,7 +745,7 @@
static void
nc_prod at c@(c at typ@ *a, c at typ@ *b, c at typ@ *r)
{
- register @typ@ ar=a->real, br=b->real, ai=a->imag, bi=b->imag;
+ @typ@ ar=a->real, br=b->real, ai=a->imag, bi=b->imag;
r->real = ar*br - ai*bi;
r->imag = ar*bi + ai*br;
return;
@@ -691,8 +755,8 @@
nc_quot at c@(c at typ@ *a, c at typ@ *b, c at typ@ *r)
{
- register @typ@ ar=a->real, br=b->real, ai=a->imag, bi=b->imag;
- register @typ@ d = br*br + bi*bi;
+ @typ@ ar=a->real, br=b->real, ai=a->imag, bi=b->imag;
+ @typ@ d = br*br + bi*bi;
r->real = (ar*br + ai*bi)/d;
r->imag = (ai*br - ar*bi)/d;
return;
@@ -801,8 +865,10 @@
return;
}
}
- /* complexobect.c uses an inline version of this formula
- investigate whether this had better performance or accuracy */
+ /*
+ * complexobect.c uses an inline version of this formula
+ * investigate whether this had better performance or accuracy
+ */
nc_log at c@(a, r);
nc_prod at c@(r, b, r);
nc_exp at c@(r, r);
@@ -823,6 +889,10 @@
static void
nc_acos at c@(c at typ@ *x, c at typ@ *r)
{
+ /*
+ * return nc_neg(nc_prodi(nc_log(nc_sum(x,nc_prod(nc_i,
+ * nc_sqrt(nc_diff(nc_1,nc_prod(x,x))))))));
+ */
nc_prod at c@(x,x,r);
nc_diff at c@(&nc_1 at c@, r, r);
nc_sqrt at c@(r, r);
@@ -832,14 +902,15 @@
nc_prodi at c@(r, r);
nc_neg at c@(r, r);
return;
- /* return nc_neg(nc_prodi(nc_log(nc_sum(x,nc_prod(nc_i,
- nc_sqrt(nc_diff(nc_1,nc_prod(x,x))))))));
- */
}
static void
nc_acosh at c@(c at typ@ *x, c at typ@ *r)
{
+ /*
+ * return nc_log(nc_sum(x,
+ * nc_prod(nc_sqrt(nc_sum(x,nc_1)), nc_sqrt(nc_diff(x,nc_1)))));
+ */
c at typ@ t;
nc_sum at c@(x, &nc_1 at c@, &t);
@@ -850,15 +921,15 @@
nc_sum at c@(x, r, r);
nc_log at c@(r, r);
return;
- /*
- return nc_log(nc_sum(x,
- nc_prod(nc_sqrt(nc_sum(x,nc_1)), nc_sqrt(nc_diff(x,nc_1)))));
- */
}
static void
nc_asin at c@(c at typ@ *x, c at typ@ *r)
{
+ /*
+ * return nc_neg(nc_prodi(nc_log(nc_sum(nc_prod(nc_i,x),
+ * nc_sqrt(nc_diff(nc_1,nc_prod(x,x)))))));
+ */
c at typ@ a, *pa=&a;
nc_prod at c@(x, x, r);
nc_diff at c@(&nc_1 at c@, r, r);
@@ -869,30 +940,29 @@
nc_prodi at c@(r, r);
nc_neg at c@(r, r);
return;
- /*
- return nc_neg(nc_prodi(nc_log(nc_sum(nc_prod(nc_i,x),
- nc_sqrt(nc_diff(nc_1,nc_prod(x,x)))))));
- */
}
static void
nc_asinh at c@(c at typ@ *x, c at typ@ *r)
{
+ /*
+ * return nc_log(nc_sum(nc_sqrt(nc_sum(nc_1,nc_prod(x,x))),x));
+ */
nc_prod at c@(x, x, r);
nc_sum at c@(&nc_1 at c@, r, r);
nc_sqrt at c@(r, r);
nc_sum at c@(r, x, r);
nc_log at c@(r, r);
return;
- /*
- return nc_log(nc_sum(nc_sqrt(nc_sum(nc_1,nc_prod(x,x))),x));
- */
}
static void
nc_atan at c@(c at typ@ *x, c at typ@ *r)
{
+ /*
+ * return nc_prod(nc_i2,nc_log(nc_quot(nc_sum(nc_i,x),nc_diff(nc_i,x))));
+ */
c at typ@ a, *pa=&a;
nc_diff at c@(&nc_i at c@, x, pa);
nc_sum at c@(&nc_i at c@, x, r);
@@ -900,14 +970,14 @@
nc_log at c@(r,r);
nc_prod at c@(&nc_i2 at c@, r, r);
return;
- /*
- return nc_prod(nc_i2,nc_log(nc_quot(nc_sum(nc_i,x),nc_diff(nc_i,x))));
- */
}
static void
nc_atanh at c@(c at typ@ *x, c at typ@ *r)
{
+ /*
+ * return nc_prod(nc_half,nc_log(nc_quot(nc_sum(nc_1,x),nc_diff(nc_1,x))));
+ */
c at typ@ a, *pa=&a;
nc_diff at c@(&nc_1 at c@, x, r);
nc_sum at c@(&nc_1 at c@, x, pa);
@@ -915,9 +985,6 @@
nc_log at c@(r, r);
nc_prod at c@(&nc_half at c@, r, r);
return;
- /*
- return nc_prod(nc_half,nc_log(nc_quot(nc_sum(nc_1,x),nc_diff(nc_1,x))));
- */
}
static void
@@ -1018,1139 +1085,975 @@
*****************************************************************************
*/
+#define OUTPUT_LOOP\
+ char *op = args[1];\
+ intp os = steps[1];\
+ intp n = dimensions[0];\
+ intp i;\
+ for(i = 0; i < n; i++, op += os)
-/**begin repeat
+#define UNARY_LOOP\
+ char *ip1 = args[0], *op = args[1];\
+ intp is1 = steps[0], os = steps[1];\
+ intp n = dimensions[0];\
+ intp i;\
+ for(i = 0; i < n; i++, ip1 += is1, op += os)
- #TYPE=(BOOL, BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG,FLOAT,DOUBLE,LONGDOUBLE)*2#
- #OP=||, +*13, ^, -*13#
- #kind=add*14, subtract*14#
- #typ=(Bool, byte, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble)*2#
-*/
+#define UNARY_LOOP_TWO_OUT\
+ char *ip1 = args[0], *op1 = args[1], *op2 = args[2];\
+ intp is1 = steps[0], os1 = steps[1], os2 = steps[2];\
+ intp n = dimensions[0];\
+ intp i;\
+ for(i = 0; i < n; i++, ip1 += is1, op1 += os1, op2 += os2)
-static void
- at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
-{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((@typ@ *)op)=*((@typ@ *)i1) @OP@ *((@typ@ *)i2);
- }
-}
+#define BINARY_LOOP\
+ char *ip1 = args[0], *ip2 = args[1], *op = args[2];\
+ intp is1 = steps[0], is2 = steps[1], os = steps[2];\
+ intp n = dimensions[0];\
+ intp i;\
+ for(i = 0; i < n; i++, ip1 += is1, ip2 += is2, op += os)
-/**end repeat**/
+#define BINARY_LOOP_TWO_OUT\
+ char *ip1 = args[0], *ip2 = args[1], *op1 = args[2], *op2 = args[3];\
+ intp is1 = steps[0], is2 = steps[1], os1 = steps[2], os2 = steps[3];\
+ intp n = dimensions[0];\
+ intp i;\
+ for(i = 0; i < n; i++, ip1 += is1, ip2 += is2, op1 += os1, op2 += os2)
-/**begin repeat
+/*
+ *****************************************************************************
+ ** BOOLEAN LOOPS **
+ *****************************************************************************
+ */
- #TYPE=(CFLOAT, CDOUBLE, CLONGDOUBLE)*2#
- #OP=+*3,-*3#
- #kind=add*3,subtract*3#
- #typ=(float, double, longdouble)*2#
-*/
+#define BOOL_invert BOOL_logical_not
+#define BOOL_negative BOOL_logical_not
+#define BOOL_add BOOL_logical_or
+#define BOOL_bitwise_and BOOL_logical_and
+#define BOOL_bitwise_or BOOL_logical_or
+#define BOOL_bitwise_xor BOOL_logical_xor
+#define BOOL_multiply BOOL_logical_and
+#define BOOL_subtract BOOL_logical_xor
+/**begin repeat
+ * #kind = equal, not_equal, greater, greater_equal, less, less_equal,
+ * logical_and, logical_or#
+ * #OP = ==, !=, >, >=, <, <=, &&, ||#
+ **/
+
static void
- at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
+BOOL_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- ((@typ@ *)op)[0]=((@typ@ *)i1)[0] @OP@ ((@typ@ *)i2)[0];
- ((@typ@ *)op)[1]=((@typ@ *)i1)[1] @OP@ ((@typ@ *)i2)[1];
+ BINARY_LOOP {
+ Bool in1 = *((Bool *)ip1) != 0;
+ Bool in2 = *((Bool *)ip2) != 0;
+ *((Bool *)op)= in1 @OP@ in2;
}
}
-
/**end repeat**/
-
static void
-BOOL_multiply(char **args, intp *dimensions, intp *steps, void *func) {
- register intp i;
- intp is1=steps[0], is2=steps[1], os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for (i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((Bool *)op) = *((Bool *)i1) && *((Bool *)i2);
+BOOL_logical_xor(char **args, intp *dimensions, intp *steps, void *func)
+{
+ BINARY_LOOP {
+ Bool in1 = *((Bool *)ip1) != 0;
+ Bool in2 = *((Bool *)ip2) != 0;
+ *((Bool *)op)= (in1 && !in2) || (!in1 && in2);
}
}
/**begin repeat
-
- #TYP= BYTE, SHORT, INT, LONG, LONGLONG, UBYTE, USHORT, UINT, ULONG, ULONGLONG, FLOAT, DOUBLE, LONGDOUBLE#
- #typ= byte, short, int, long, longlong, ubyte, ushort, uint, ulong, ulonglong, float, double, longdouble#
-*/
+ * #kind = maximum, minimum#
+ * #OP = >, <#
+ **/
static void
- at TYP@_multiply(char **args, intp *dimensions, intp *steps, void *func)
+BOOL_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0], is2=steps[1], os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for (i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((@typ@ *)op) = (*((@typ@ *)i1)) * (*((@typ@ *)i2));
+ BINARY_LOOP {
+ Bool in1 = *((Bool *)ip1) != 0;
+ Bool in2 = *((Bool *)ip2) != 0;
+ *((Bool *)op) = (in1 @OP@ in2) ? in1 : in2;
}
}
/**end repeat**/
-
/**begin repeat
- #TYP= CFLOAT, CDOUBLE, CLONGDOUBLE#
- #typ= float, double, longdouble#
- #c=f,,l#
-*/
+ * #kind = absolute, logical_not#
+ * #OP = !=, ==#
+ **/
static void
- at TYP@_multiply(char **args, intp *dimensions, intp *steps, void *func)
+BOOL_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0], is2=steps[1], os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for (i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- register @typ@ ar=((c at typ@ *)i1)->real, \
- ai=((c at typ@ *)i1)->imag, \
- br=((c at typ@ *)i2)->real, \
- bi=((c at typ@ *)i2)->imag;
- ((c at typ@ *)op)->real = ar*br - ai*bi;
- ((c at typ@ *)op)->imag = ar*bi + ai*br;
+ UNARY_LOOP {
+ Bool in1 = *(Bool *)ip1;
+ *((Bool *)op) = in1 @OP@ 0;
}
}
+/**end repeat**/
static void
- at TYP@_divide(char **args, intp *dimensions, intp *steps, void *func)
+BOOL_ones_like(char **args, intp *dimensions, intp *steps, void *data)
{
- register intp i;
- intp is1=steps[0], is2=steps[1], os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for (i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- register @typ@ ar=((c at typ@ *)i1)->real, \
- ai=((c at typ@ *)i1)->imag, \
- br=((c at typ@ *)i2)->real, \
- bi=((c at typ@ *)i2)->imag;
- register @typ@ d = br*br + bi*bi;
- ((c at typ@ *)op)->real = (ar*br + ai*bi)/d;
- ((c at typ@ *)op)->imag = (ai*br - ar*bi)/d;
+ OUTPUT_LOOP {
+ *((Bool *)op) = 1;
}
}
+
+/*
+ *****************************************************************************
+ ** INTEGER LOOPS
+ *****************************************************************************
+ */
+
+/**begin repeat
+ * #type = byte, short, int, long, longlong#
+ * #TYPE = BYTE, SHORT, INT, LONG, LONGLONG#
+ * #ftype = float, float, double, double, double#
+ */
+
+/**begin repeat1
+ * both signed and unsigned integer types
+ * #s = , u#
+ * #S = , U#
+ */
+
+#define @S@@TYPE at _floor_divide @S@@TYPE at _divide
+
static void
- at TYP@_floor_divide(char **args, intp *dimensions, intp *steps, void *func)
+ at S@@TYPE at _ones_like(char **args, intp *dimensions, intp *steps, void *data)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2],n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- register @typ@ ar=((c at typ@ *)i1)->real, \
- ai=((c at typ@ *)i1)->imag, \
- br=((c at typ@ *)i2)->real, \
- bi=((c at typ@ *)i2)->imag;
- register @typ@ d = br*br + bi*bi;
- ((c at typ@ *)op)->real = floor at c@((ar*br + ai*bi)/d);
- ((c at typ@ *)op)->imag = 0;
+ OUTPUT_LOOP {
+ *((@s@@type@ *)op) = 1;
}
}
-#define @TYP at _true_divide @TYP at _divide
-/**end repeat**/
-
-
-/**begin repeat
- #TYP=UBYTE,USHORT,UINT,ULONG,ULONGLONG#
- #typ=ubyte, ushort, uint, ulong, ulonglong#
-*/
static void
- at TYP@_divide(char **args, intp *dimensions, intp *steps, void *func)
+ at S@@TYPE at _square(char **args, intp *dimensions, intp *steps, void *data)
{
- register intp i, is1=steps[0],is2=steps[1],os=steps[2],n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- if (*((@typ@ *)i2)==0) {
- generate_divbyzero_error();
- *((@typ@ *)op)=0;
- }
- else {
- *((@typ@ *)op)= *((@typ@ *)i1) / *((@typ@ *)i2);
- }
+ UNARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ *((@s@@type@ *)op) = in1*in1;
}
}
-/**end repeat**/
-
-/**begin repeat
- #TYP=BYTE,SHORT,INT,LONG,LONGLONG#
- #typ=char, short, int, long, longlong#
-*/
static void
- at TYP@_divide(char **args, intp *dimensions, intp *steps, void *func)
+ at S@@TYPE at _reciprocal(char **args, intp *dimensions, intp *steps, void *data)
{
- register intp i, is1=steps[0],is2=steps[1],os=steps[2],n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- @typ@ x, y, tmp;
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- y = *((@typ@ *)i2);
- if (y == 0) {
- generate_divbyzero_error();
- *((@typ@ *)op)=0;
- }
- else {
- x = *((@typ@ *)i1);
- tmp = x / y;
- if (((x > 0) != (y > 0)) && (x % y != 0)) tmp--;
- *((@typ@ *)op)= tmp;
- }
+ UNARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ *((@s@@type@ *)op) = 1.0/in1;
}
}
-/**end repeat**/
-
-/**begin repeat
- #TYP=BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG#
- #typ=char, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong#
- #otyp=float*4, double*6#
-*/
static void
- at TYP@_true_divide(char **args, intp *dimensions, intp *steps, void *func)
+ at S@@TYPE at _conjugate(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i, is1=steps[0],is2=steps[1],os=steps[2],n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- if (*((@typ@ *)i2)==0) {
- generate_divbyzero_error();
- *((@otyp@ *)op)=0;
- }
- else {
- *((@otyp@ *)op)= \
- (@otyp@)((double)*((@typ@ *)i1) / (double)*((@typ@ *)i2));
- }
+ UNARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ *((@s@@type@ *)op) = in1;
}
}
-#define @TYP at _floor_divide @TYP at _divide
-/**end repeat**/
-
-
-/**begin repeat
- #TYP=FLOAT,DOUBLE,LONGDOUBLE#
- #typ=float,double,longdouble#
-*/
static void
- at TYP@_divide(char **args, intp *dimensions, intp *steps, void *func)
+ at S@@TYPE at _negative(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i, is1=steps[0],is2=steps[1],os=steps[2],n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((@typ@ *)op)=*((@typ@ *)i1) / *((@typ@ *)i2);
+ UNARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ *((@s@@type@ *)op) = (@s@@type@)(-(@type@)in1);
}
}
-#define @TYP at _true_divide @TYP at _divide
-/**end repeat**/
-/**begin repeat
-
- #TYP=FLOAT,DOUBLE,LONGDOUBLE#
- #typ=float,double,longdouble#
- #c=f,,l#
-*/
static void
- at TYP@_floor_divide(char **args, intp *dimensions, intp *steps, void *func)
+ at S@@TYPE at _logical_not(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i, is1=steps[0],is2=steps[1],os=steps[2],n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((@typ@ *)op)=floor at c@(*((@typ@ *)i1) / *((@typ@ *)i2));
+ UNARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ *((Bool *)op) = !in1;
}
}
-/**end repeat**/
-/**begin repeat
- #TYP=BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG,FLOAT,DOUBLE,LONGDOUBLE#
- #typ=char, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble#
-*/
static void
- at TYP@_square(char **args, intp *dimensions, intp *steps, void *data)
+ at S@@TYPE at _invert(char **args, intp *dimensions, intp *steps, void *func)
{
- intp i, is1 = steps[0], os = steps[1], n = dimensions[0];
- char *i1 = args[0], *op = args[1];
-
- for (i = 0; i < n; i++, i1 += is1, op += os) {
- @typ@ x = *((@typ@ *)i1);
- *((@typ@ *)op) = x*x;
+ UNARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ *((@s@@type@ *)op) = ~in1;
}
}
-/**end repeat**/
-/**begin repeat
- #TYP=CFLOAT,CDOUBLE,CLONGDOUBLE#
- #typ=float, double, longdouble#
-*/
+/**begin repeat2
+ * Arithmetic
+ * #kind = add, subtract, multiply, bitwise_and, bitwise_or, bitwise_xor,
+ * left_shift, right_shift#
+ * #OP = +, -,*, &, |, ^, <<, >>#
+ */
static void
- at TYP@_square(char **args, intp *dimensions, intp *steps, void *data)
+ at S@@TYPE at _@kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- intp i, is1 = steps[0], os = steps[1], n = dimensions[0];
- char *i1 = args[0], *op = args[1];
- c at typ@ *x, *y;
- @typ@ xr, xi;
-
- for (i = 0; i < n; i++, i1 += is1, op += os) {
- x = (c at typ@ *)i1;
- y = (c at typ@ *)op;
- xr = x->real;
- xi = x->imag;
- y->real = xr*xr - xi*xi;
- y->imag = 2*xr*xi;
+ BINARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ const @s@@type@ in2 = *(@s@@type@ *)ip2;
+ *((@s@@type@ *)op) = in1 @OP@ in2;
}
}
-/**end repeat**/
+/**end repeat2**/
-static PyObject *
-Py_square(PyObject *o)
+/**begin repeat2
+ * #kind = equal, not_equal, greater, greater_equal, less, less_equal,
+ * logical_and, logical_or#
+ * #OP = ==, !=, >, >=, <, <=, &&, ||#
+ */
+static void
+ at S@@TYPE at _@kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- return PyNumber_Multiply(o, o);
+ BINARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ const @s@@type@ in2 = *(@s@@type@ *)ip2;
+ *((Bool *)op) = in1 @OP@ in2;
+ }
}
+/**end repeat2**/
-
-/**begin repeat
- #TYP=BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG,FLOAT,DOUBLE,LONGDOUBLE#
- #typ=char, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble#
-*/
static void
- at TYP@_reciprocal(char **args, intp *dimensions, intp *steps, void *data)
+ at S@@TYPE at _logical_xor(char **args, intp *dimensions, intp *steps, void *func)
{
- intp i, is1 = steps[0], os = steps[1], n = dimensions[0];
- char *i1 = args[0], *op = args[1];
-
- for (i = 0; i < n; i++, i1 += is1, op += os) {
- @typ@ x = *((@typ@ *)i1);
- *((@typ@ *)op) = (@typ@) (1.0 / x);
+ BINARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ const @s@@type@ in2 = *(@s@@type@ *)ip2;
+ *((Bool *)op)= (in1 && !in2) || (!in1 && in2);
}
}
-/**end repeat**/
-/**begin repeat
- #TYP=CFLOAT,CDOUBLE,CLONGDOUBLE#
- #typ=float, double, longdouble#
-*/
+/**begin repeat2
+ * #kind = maximum, minimum#
+ * #OP = >, <#
+ **/
static void
- at TYP@_reciprocal(char **args, intp *dimensions, intp *steps, void *data)
+ at S@@TYPE at _@kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- intp i, is1 = steps[0], os = steps[1], n = dimensions[0];
- char *i1 = args[0], *op = args[1];
- c at typ@ *x, *y;
- @typ@ xr, xi, r, denom;
-
- for (i = 0; i < n; i++, i1 += is1, op += os) {
- x = (c at typ@ *)i1;
- y = (c at typ@ *)op;
- xr = x->real;
- xi = x->imag;
- if (fabs(xi) <= fabs(xr)) {
- r = xi / xr;
- denom = xr + xi * r;
- y->real = 1 / denom;
- y->imag = -r / denom;
- } else {
- r = xr / xi;
- denom = xr * r + xi;
- y->real = r / denom;
- y->imag = -1 / denom;
- }
+ BINARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ const @s@@type@ in2 = *(@s@@type@ *)ip2;
+ *((@s@@type@ *)op) = (in1 @OP@ in2) ? in1 : in2;
}
}
-/**end repeat**/
+/**end repeat2**/
-
-static PyObject *
-Py_reciprocal(PyObject *o)
+static void
+ at S@@TYPE at _true_divide(char **args, intp *dimensions, intp *steps, void *func)
{
- PyObject *one, *result;
- one = PyInt_FromLong(1);
- if (!one) return NULL;
- result = PyNumber_Divide(one, o);
- Py_DECREF(one);
- return result;
+ BINARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ const @s@@type@ in2 = *(@s@@type@ *)ip2;
+ if (in2 == 0) {
+ generate_divbyzero_error();
+ *((@ftype@ *)op) = 0;
+ }
+ else {
+ *((@ftype@ *)op) = (@ftype@)in1 / (@ftype@)in2;
+ }
+ }
}
-static PyObject *
-_npy_ObjectMax(PyObject *i1, PyObject *i2)
+static void
+ at S@@TYPE at _power(char **args, intp *dimensions, intp *steps, void *func)
{
- int cmp;
- PyObject *res;
- if (PyObject_Cmp(i1, i2, &cmp) < 0) return NULL;
-
- if (cmp >= 0) {
- res = i1;
+ BINARY_LOOP {
+ const @ftype@ in1 = (@ftype@)*(@s@@type@ *)ip1;
+ const @ftype@ in2 = (@ftype@)*(@s@@type@ *)ip2;
+ *((@s@@type@ *)op) = (@s@@type@) pow(in1, in2);
}
- else {
- res = i2;
- }
- Py_INCREF(res);
- return res;
}
-static PyObject *
-_npy_ObjectMin(PyObject *i1, PyObject *i2)
+static void
+ at S@@TYPE at _fmod(char **args, intp *dimensions, intp *steps, void *func)
{
- int cmp;
- PyObject *res;
- if (PyObject_Cmp(i1, i2, &cmp) < 0) return NULL;
+ BINARY_LOOP {
+ const @s@@type@ in1 = *(@s@@type@ *)ip1;
+ const @s@@type@ in2 = *(@s@@type@ *)ip2;
+ if (in2 == 0) {
+ generate_divbyzero_error();
+ *((@s@@type@ *)op) = 0;
+ }
+ else {
+ *((@s@@type@ *)op)= in1 % in2;
+ }
- if (cmp <= 0) {
- res = i1;
}
- else {
- res = i2;
- }
- Py_INCREF(res);
- return res;
}
-/* ones_like is defined here because it's used for x**0 */
+/**end repeat1**/
-/**begin repeat
- #TYP=BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG,FLOAT,DOUBLE,LONGDOUBLE#
- #typ=char, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble#
-*/
static void
- at TYP@_ones_like(char **args, intp *dimensions, intp *steps, void *data)
+U at TYPE@_absolute(char **args, intp *dimensions, intp *steps, void *func)
{
- intp i, os = steps[1], n = dimensions[0];
- char *op = args[1];
-
- for (i = 0; i < n; i++, op += os) {
- *((@typ@ *)op) = 1;
+ UNARY_LOOP {
+ const u at type@ in1 = *(u at type@ *)ip1;
+ *((u at type@ *)op) = in1;
}
}
-/**end repeat**/
-/**begin repeat
- #TYP=CFLOAT,CDOUBLE,CLONGDOUBLE#
- #typ=float, double, longdouble#
-*/
static void
- at TYP@_ones_like(char **args, intp *dimensions, intp *steps, void *data)
+ at TYPE@_absolute(char **args, intp *dimensions, intp *steps, void *func)
{
- intp i, is1 = steps[0], os = steps[1], n = dimensions[0];
- char *i1 = args[0], *op = args[1];
- c at typ@ *y;
-
- for (i = 0; i < n; i++, i1 += is1, op += os) {
- y = (c at typ@ *)op;
- y->real = 1.0;
- y->imag = 0.0;
+ UNARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ *((@type@ *)op) = (in1 >= 0) ? in1 : -in1;
}
}
-/**end repeat**/
-static PyObject *
-Py_get_one(PyObject *o)
+static void
+U at TYPE@_sign(char **args, intp *dimensions, intp *steps, void *func)
{
- return PyInt_FromLong(1);
+ UNARY_LOOP {
+ const u at type@ in1 = *(u at type@ *)ip1;
+ *((u at type@ *)op) = in1 > 0 ? 1 : 0;
+ }
}
-
-/**begin repeat
- #TYP=BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG#
- #typ=char, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong#
- #btyp=float*4, double*6#
-*/
static void
- at TYP@_power(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_sign(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i, is1=steps[0],is2=steps[1];
- register intp os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- @btyp@ x, y;
+ UNARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ *((@type@ *)op) = in1 > 0 ? 1 : (in1 < 0 ? -1 : 0);
+ }
+}
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- x = (@btyp@)*((@typ@ *)i1);
- y = (@btyp@)*((@typ@ *)i2);
- *((@typ@ *)op) = (@typ@) pow(x,y);
+static void
+ at TYPE@_divide(char **args, intp *dimensions, intp *steps, void *func)
+{
+ BINARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ const @type@ in2 = *(@type@ *)ip2;
+ if (in2 == 0) {
+ generate_divbyzero_error();
+ *((@type@ *)op) = 0;
+ }
+ else if (((in1 > 0) != (in2 > 0)) && (in1 % in2 != 0)) {
+ *((@type@ *)op) = in1/in2 - 1;
+ }
+ else {
+ *((@type@ *)op) = in1/in2;
+ }
}
}
-/**end repeat**/
-/**begin repeat
- #TYP=UBYTE, BYTE, SHORT, USHORT, INT, UINT, LONG, ULONG, LONGLONG, ULONGLONG, FLOAT, DOUBLE, LONGDOUBLE#
- #typ=ubyte, char, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble#
-*/
static void
- at TYP@_conjugate(char **args, intp *dimensions, intp *steps, void *func)
+U at TYPE@_divide(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i, is1=steps[0], os=steps[1], n=dimensions[0];
- char *i1=args[0], *op=args[1];
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- *((@typ@ *)op)=*((@typ@ *)i1);
+ BINARY_LOOP {
+ const u at type@ in1 = *(u at type@ *)ip1;
+ const u at type@ in2 = *(u at type@ *)ip2;
+ if (in2 == 0) {
+ generate_divbyzero_error();
+ *((u at type@ *)op) = 0;
+ }
+ else {
+ *((u at type@ *)op)= in1/in2;
+ }
}
}
-/**end repeat**/
-/**begin repeat
-
- #TYP=CFLOAT, CDOUBLE, CLONGDOUBLE#
- #typ=float, double, longdouble#
-*/
static void
- at TYP@_conjugate(char **args, intp *dimensions, intp *steps, void *func) {
- register intp i, is1=steps[0], os=steps[1], n=dimensions[0];
- char *i1=args[0], *op=args[1];
+ at TYPE@_remainder(char **args, intp *dimensions, intp *steps, void *func)
+{
+ BINARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ const @type@ in2 = *(@type@ *)ip2;
+ if (in2 == 0) {
+ generate_divbyzero_error();
+ *((@type@ *)op) = 0;
+ }
+ else {
+ /* handle mixed case the way Python does */
+ const @type@ rem = in1 % in2;
+ if ((in1 > 0) == (in2 > 0) || rem == 0) {
+ *((@type@ *)op) = rem;
+ }
+ else {
+ *((@type@ *)op) = rem + in2;
+ }
+ }
+ }
+}
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- ((@typ@ *)op)[0]=((@typ@ *)i1)[0];
- ((@typ@ *)op)[1]=-(((@typ@ *)i1)[1]);
+static void
+U at TYPE@_remainder(char **args, intp *dimensions, intp *steps, void *func)
+{
+ BINARY_LOOP {
+ const u at type@ in1 = *(u at type@ *)ip1;
+ const u at type@ in2 = *(u at type@ *)ip2;
+ if (in2 == 0) {
+ generate_divbyzero_error();
+ *((@type@ *)op) = 0;
+ }
+ else {
+ *((@type@ *)op) = in1 % in2;
+ }
}
}
+
/**end repeat**/
+/*
+ *****************************************************************************
+ ** FLOAT LOOPS **
+ *****************************************************************************
+ */
+
/**begin repeat
+ * Float types
+ * #type = float, double, longdouble#
+ * #TYPE = FLOAT, DOUBLE, LONGDOUBLE#
+ * #c = f, , l#
+ * #C = F, , L#
+ */
- #TYPE=BOOL,UBYTE,USHORT,UINT,ULONG,ULONGLONG#
- #typ=Bool, ubyte, ushort, uint, ulong, ulonglong#
-*/
+/**begin repeat1
+ * Arithmetic
+ * # kind = add, subtract, multiply, divide#
+ * # OP = +, -, *, /#
+ */
static void
- at TYPE@_absolute(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i, n;
- intp is1=steps[0], os=steps[1];
- char *i1=args[0], *op=args[1];
-
- n=dimensions[0];
-
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- *((@typ@ *)op) = *((@typ@*)i1);
+ BINARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ const @type@ in2 = *(@type@ *)ip2;
+ *((@type@ *)op) = in1 @OP@ in2;
}
}
-/**end repeat**/
+/**end repeat1**/
-/**begin repeat
-
- #TYPE=BYTE,SHORT,INT,LONG,LONGLONG,FLOAT,DOUBLE,LONGDOUBLE#
- #typ=byte, short, int, long, longlong, float, double, longdouble#
-*/
+/**begin repeat1
+ * #kind = equal, not_equal, less, less_equal, greater, greater_equal,
+ * logical_and, logical_or#
+ * #OP = ==, !=, <, <=, >, >=, &&, ||#
+ */
static void
- at TYPE@_absolute(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i, n;
- intp is1=steps[0], os=steps[1];
- char *i1=args[0], *op=args[1];
-
- n=dimensions[0];
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- *((@typ@ *)op) = *((@typ@ *)i1) < 0 ? -*((@typ@ *)i1) : *((@typ@ *)i1);
- *((@typ@ *)op) += 0; /* clear sign-bit */
+ BINARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ const @type@ in2 = *(@type@ *)ip2;
+ *((Bool *)op) = in1 @OP@ in2;
}
}
-/**end repeat**/
+/**end repeat1**/
-/**begin repeat
- #TYPE=CFLOAT,CDOUBLE,CLONGDOUBLE#
- #typ= float, double, longdouble#
- #c= f,,l#
-*/
static void
- at TYPE@_absolute(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_logical_xor(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i, n;
- register intp is1=steps[0], os=steps[1];
- char *i1=args[0], *op=args[1];
- n=dimensions[0];
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- *((@typ@ *)op) = (@typ@)sqrt at c@(((@typ@ *)i1)[0]*((@typ@ *)i1)[0] + ((@typ@ *)i1)[1]*((@typ@ *)i1)[1]);
+ BINARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ const @type@ in2 = *(@type@ *)ip2;
+ *((Bool *)op)= (in1 && !in2) || (!in1 && in2);
}
}
-/**end repeat**/
-/**begin repeat
-
- #kind=greater, greater_equal, less, less_equal, equal, not_equal, logical_and, logical_or, bitwise_and, bitwise_or, bitwise_xor#
- #OP=>, >=, <, <=, ==, !=, &&, ||, &, |, ^#
-**/
static void
-BOOL_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_logical_not(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- Bool in1, in2;
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- in1 = (*((Bool *)i1) != 0);
- in2 = (*((Bool *)i2) != 0);
- *((Bool *)op)= in1 @OP@ in2;
+ UNARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ *((Bool *)op) = !in1;
}
}
-/**end repeat**/
-/**begin repeat
-
- #TYPE=(BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG,FLOAT,DOUBLE,LONGDOUBLE)*4#
- #OP= >*13, >=*13, <*13, <=*13#
- #typ=(byte, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble)*4#
- #kind= greater*13, greater_equal*13, less*13, less_equal*13#
-*/
-
+/**begin repeat1
+ * #kind = isnan, isinf, isfinite, signbit#
+ * #func = isnan, isinf, isfinite, signbit#
+ **/
static void
@TYPE at _@kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((Bool *)op)=*((@typ@ *)i1) @OP@ *((@typ@ *)i2);
+ UNARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ *((Bool *)op) = @func@(in1) != 0;
}
}
-/**end repeat**/
+/**end repeat1**/
-
-/**begin repeat
- #TYPE=(CFLOAT,CDOUBLE,CLONGDOUBLE)*4#
- #OP= >*3, >=*3, <*3, <=*3#
- #typ=(cfloat, cdouble, clongdouble)*4#
- #kind= greater*3, greater_equal*3, less*3, less_equal*3#
-*/
-
+/**begin repeat1
+ * #kind = maximum, minimum#
+ * #OP = >, <#
+ **/
static void
@TYPE at _@kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- if (((@typ@ *)i1)->real == ((@typ@ *)i2)->real)
- *((Bool *)op)=((@typ@ *)i1)->imag @OP@ \
- ((@typ@ *)i2)->imag;
- else
- *((Bool *)op)=((@typ@ *)i1)->real @OP@ \
- ((@typ@ *)i2)->real;
+ /* */
+ BINARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ const @type@ in2 = *(@type@ *)ip2;
+ *((@type@ *)op) = in1 @OP@ in2 ? in1 : in2;
}
}
-/**end repeat**/
+/**end repeat1**/
+static void
+ at TYPE@_floor_divide(char **args, intp *dimensions, intp *steps, void *func)
+{
+ BINARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ const @type@ in2 = *(@type@ *)ip2;
+ *((@type@ *)op) = floor at c@(in1/in2);
+ }
+}
-/**begin repeat
- #TYPE=(BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG,FLOAT,DOUBLE,LONGDOUBLE)*4#
- #typ=(byte, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble)*4#
- #OP= ==*13, !=*13, &&*13, ||*13#
- #kind=equal*13, not_equal*13, logical_and*13, logical_or*13#
-*/
static void
- at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_remainder(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((Bool *)op) = *((@typ@ *)i1) @OP@ *((@typ@ *)i2);
+ BINARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ const @type@ in2 = *(@type@ *)ip2;
+ const @type@ res = fmod at c@(in1,in2);
+ if (res && ((in2 < 0) != (res < 0))) {
+ *((@type@ *)op) = res + in2;
+ }
+ else {
+ *((@type@ *)op) = res;
+ }
}
}
-/**end repeat**/
+static void
+ at TYPE@_square(char **args, intp *dimensions, intp *steps, void *data)
+{
+ UNARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ *((@type@ *)op) = in1*in1;
+ }
+}
-/**begin repeat
+static void
+ at TYPE@_reciprocal(char **args, intp *dimensions, intp *steps, void *data)
+{
+ UNARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ *((@type@ *)op) = 1.0/in1;
+ }
+}
- #TYPE=(CFLOAT, CDOUBLE, CLONGDOUBLE)*4#
- #typ=(float, double, longdouble)*4#
- #OP= ==*3, !=*3, &&*3, ||*3#
- #OP2= &&*3, ||*3, &&*3, ||*3#
- #kind=equal*3, not_equal*3, logical_and*3, logical_or*3#
-*/
static void
- at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_ones_like(char **args, intp *dimensions, intp *steps, void *data)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((Bool *)op) = (*((@typ@ *)i1) @OP@ *((@typ@ *)i2)) @OP2@ (*((@typ@ *)i1+1) @OP@ *((@typ@ *)i2+1));
+ OUTPUT_LOOP {
+ *((@type@ *)op) = 1;
}
}
-/**end repeat**/
-
-/** OBJECT comparison for OBJECT arrays **/
-
-/**begin repeat
-
- #kind=greater, greater_equal, less, less_equal, equal, not_equal#
- #op=GT, GE, LT, LE, EQ, NE#
-*/
static void
-OBJECT_ at kind@(char **args, intp *dimensions, intp *steps, void *func) {
- register intp i, is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((Bool *)op)=PyObject_RichCompareBool(*((PyObject **)i1),
- *((PyObject **)i2),
- Py_ at op@);
+ at TYPE@_conjugate(char **args, intp *dimensions, intp *steps, void *func)
+{
+ UNARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ *((@type@ *)op) = in1;
}
}
-/**end repeat**/
-/**begin repeat
-
- #TYPE=BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG,FLOAT,DOUBLE,LONGDOUBLE#
- #typ=byte, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble#
- #styp=byte, byte, short, short, int, int, long, long, longlong, longlong, float, double, longdouble#
-*/
static void
- at TYPE@_negative(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_absolute(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],os=steps[1], n=dimensions[0];
- char *i1=args[0], *op=args[1];
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- *((@typ@ *)op) = (@typ@) (- (@styp@)*((@typ@ *)i1));
+ UNARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ const @type@ tmp = (in1 > 0) ? in1 : -in1;
+ /* add 0 to clear -0.0 */
+ *((@type@ *)op) = tmp + 0;
}
}
-/**end repeat**/
-#define BOOL_negative BOOL_logical_not
-
-#define _SIGN1(x) ((x) > 0 ? 1 : ((x) < 0 ? -1 : 0))
-#define _SIGN2(x) ((x) == 0 ? 0 : 1)
-#define _SIGNC(x) (((x).real > 0) ? 1 : ((x).real < 0 ? -1 : ((x).imag > 0 ? 1 : ((x).imag < 0) ? -1 : 0)))
-/**begin repeat
- #TYPE=BYTE,SHORT,INT,LONG,LONGLONG,FLOAT,DOUBLE,LONGDOUBLE,UBYTE,USHORT,UINT,ULONG,ULONGLONG#
- #typ=byte,short,int,long,longlong,float,double,longdouble,ubyte,ushort,uint,ulong,ulonglong#
- #func=_SIGN1*8,_SIGN2*5#
-*/
static void
- at TYPE@_sign(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_negative(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],os=steps[1], n=dimensions[0];
- char *i1=args[0], *op=args[1];
- @typ@ t1;
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- t1 = *((@typ@ *)i1);
- *((@typ@ *)op) = (@typ@) @func@(t1);
+ UNARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ *((@type@ *)op) = -in1;
}
}
-/**end repeat**/
-/**begin repeat
- #TYPE=CFLOAT,CDOUBLE,CLONGDOUBLE#
- #typ=cfloat,cdouble,clongdouble#
- #rtyp=float,double,longdouble#
-*/
static void
@TYPE at _sign(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],os=steps[1], n=dimensions[0];
- char *i1=args[0], *op=args[1];
- @typ@ t1;
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- t1 = *((@typ@ *)i1);
- (*((@typ@ *)op)).real = (@rtyp@)_SIGNC(t1);
- (*((@typ@ *)op)).imag = (@rtyp@)0;
+ /* */
+ UNARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ if (in1 > 0) {
+ *((@type@ *)op) = 1;
+ }
+ else if (in1 < 0) {
+ *((@type@ *)op) = -1;
+ }
+ else {
+ *((@type@ *)op) = 0;
+ }
}
}
-/**end repeat**/
-#undef _SIGN1
-#undef _SIGN2
-#undef _SIGNC
-
-
static void
-OBJECT_sign(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_modf(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],os=steps[1], n=dimensions[0];
- char *i1=args[0], *op=args[1];
- PyObject *t1, *zero, *res;
- zero = PyInt_FromLong(0);
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- t1 = *((PyObject **)i1);
- res = PyInt_FromLong((long) PyObject_Compare(t1, zero));
- *((PyObject **)op) = res;
+ UNARY_LOOP_TWO_OUT {
+ const @type@ in1 = *(@type@ *)ip1;
+ *(@type@ *)op1 = modf at c@(in1, (@type@ *)op2);
}
- Py_DECREF(zero);
}
-
-/**begin repeat
- #TYPE=BOOL,BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG,FLOAT,DOUBLE,LONGDOUBLE#
- #typ=Bool, byte, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble#
-*/
+#define HAVE_DOUBLE_FUNCS
+#ifdef HAVE_ at TYPE@_FUNCS
static void
- at TYPE@_logical_not(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_frexp(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],os=steps[1], n=dimensions[0];
- char *i1=args[0], *op=args[1];
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- *((Bool *)op) = ! *((@typ@ *)i1);
+ UNARY_LOOP_TWO_OUT {
+ const @type@ in1 = *(@type@ *)ip1;
+ *(@type@ *)op1 = frexp at c@(in1, (int *)op2);
}
}
-/**end repeat**/
-/**begin repeat
- #TYPE=CFLOAT,CDOUBLE,CLONGDOUBLE#
- #typ=cfloat, cdouble, clongdouble#
-*/
static void
- at TYPE@_logical_not(char **args, intp *dimensions, intp *steps, void *func)
+ at TYPE@_ldexp(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],os=steps[1], n=dimensions[0];
- char *i1=args[0], *op=args[1];
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- *((Bool *)op) = ! (((@typ@ *)i1)->real || \
- ((@typ@ *)i1)->imag);
+ BINARY_LOOP {
+ const @type@ in1 = *(@type@ *)ip1;
+ const int in2 = *(int *)ip2;
+ *(@type@ *)op = ldexp at c@(in1, in2);
}
}
-/**end repeat**/
+#endif
+#undef HAVE_DOUBLE_FUNCS
+#define @TYPE at _true_divide @TYPE at _divide
+/**end repeat**/
+/*
+ *****************************************************************************
+ ** COMPLEX LOOPS **
+ *****************************************************************************
+ */
+
/**begin repeat
- #TYPE=BYTE,SHORT,INT,LONG,LONGLONG#
- #typ=byte, short, int, long, longlong#
- #c=f*2,,,l*1#
-*/
+ * complex types
+ * #ctype= cfloat, cdouble, clongdouble#
+ * #CTYPE= CFLOAT, CDOUBLE, CLONGDOUBLE#
+ * #type = float, double, longdouble#
+ * #c = f, , l#
+ */
+
+/**begin repeat1
+ * arithmetic
+ * #kind = add, subtract#
+ * #OP = +, -#
+ */
static void
- at TYPE@_remainder(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- register @typ@ ix,iy, tmp;
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- ix = *((@typ@ *)i1);
- iy = *((@typ@ *)i2);
- if (iy == 0 || ix == 0) {
- if (iy == 0) generate_divbyzero_error();
- *((@typ@ *)op) = 0;
- }
- else if ((ix > 0) == (iy > 0)) {
- *((@typ@ *)op) = ix % iy;
- }
- else { /* handle mixed case the way Python does */
- tmp = ix % iy;
- if (tmp) tmp += iy;
- *((@typ@ *)op)= tmp;
- }
+ BINARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ const @type@ in2r = ((@type@ *)ip2)[0];
+ const @type@ in2i = ((@type@ *)ip2)[1];
+ ((@type@ *)op)[0] = in1r @OP@ in2r;
+ ((@type@ *)op)[1] = in1i @OP@ in2i;
}
}
-/**end repeat**/
+/**end repeat1**/
-/**begin repeat
- #TYPE=UBYTE,USHORT,UINT,ULONG,ULONGLONG#
- #typ=ubyte, ushort, uint, ulong, ulonglong#
-*/
static void
- at TYPE@_remainder(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_multiply(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- register @typ@ ix,iy;
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- ix = *((@typ@ *)i1);
- iy = *((@typ@ *)i2);
- if (iy == 0) {
- generate_divbyzero_error();
- *((@typ@ *)op) = 0;
- }
- *((@typ@ *)op) = ix % iy;
+ BINARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ const @type@ in2r = ((@type@ *)ip2)[0];
+ const @type@ in2i = ((@type@ *)ip2)[1];
+ ((@type@ *)op)[0] = in1r*in2r - in1i*in2i;
+ ((@type@ *)op)[1] = in1r*in2i + in1i*in2r;
}
}
-/**end repeat**/
-/**begin repeat
- #TYPE=FLOAT,DOUBLE,LONGDOUBLE#
- #typ=float,double,longdouble#
- #c=f,,l#
-*/
static void
- at TYPE@_remainder(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_divide(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- @typ@ x, y, res;
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- x = *((@typ@ *)i1);
- y = *((@typ@ *)i2);
- res = fmod at c@(x, y);
- if (res && ((y < 0) != (res < 0))) {
- res += y;
- }
- *((@typ@ *)op)= res;
+ BINARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ const @type@ in2r = ((@type@ *)ip2)[0];
+ const @type@ in2i = ((@type@ *)ip2)[1];
+ @type@ d = in2r*in2r + in2i*in2i;
+ ((@type@ *)op)[0] = (in1r*in2r + in1i*in2i)/d;
+ ((@type@ *)op)[1] = (in1i*in2r - in1r*in2i)/d;
}
}
-/**end repeat**/
+static void
+ at CTYPE@_floor_divide(char **args, intp *dimensions, intp *steps, void *func)
+{
+ BINARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ const @type@ in2r = ((@type@ *)ip2)[0];
+ const @type@ in2i = ((@type@ *)ip2)[1];
+ @type@ d = in2r*in2r + in2i*in2i;
+ ((@type@ *)op)[0] = floor at c@((in1r*in2r + in1i*in2i)/d);
+ ((@type@ *)op)[1] = 0;
+ }
+}
-/**begin repeat
-
- #TYPE=BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG#
- #typ=byte, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong#
+/**begin repeat1
+ #kind = equal, not_equal#
+ #OP1 = ==, !=#
+ #OP2 = &&, ||#
*/
static void
- at TYPE@_fmod(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- @typ@ x, y;
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- x = *((@typ@ *)i1);
- y = *((@typ@ *)i2);
- if (y == 0) {
- generate_divbyzero_error();
- *((@typ@ *)op) = 0;
+ BINARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ const @type@ in2r = ((@type@ *)ip2)[0];
+ const @type@ in2i = ((@type@ *)ip2)[1];
+ *((Bool *)op) = (in1r @OP1@ in2r) @OP2@ (in1i @OP1@ in2i);
+ }
+}
+/**end repeat1**/
+
+/**begin repeat1
+ * #kind= greater, greater_equal, less, less_equal#
+ * #OP = >, >=, <, <=#
+ */
+static void
+ at CTYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
+{
+ BINARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ const @type@ in2r = ((@type@ *)ip2)[0];
+ const @type@ in2i = ((@type@ *)ip2)[1];
+ if (in1r != in2r) {
+ *((Bool *)op) = in1r @OP@ in2r ? 1 : 0;
}
else {
- *((@typ@ *)op)= x % y;
+ *((Bool *)op) = in1i @OP@ in2i ? 1 : 0;
}
-
}
}
-/**end repeat**/
+/**end repeat1**/
-/**begin repeat
-
- #TYPE=(BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG)*5#
- #typ=(byte, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong)*5#
- #OP= &*10, |*10, ^*10, <<*10, >>*10#
- #kind=bitwise_and*10, bitwise_or*10, bitwise_xor*10, left_shift*10, right_shift*10#
-
+/**begin repeat1
+ #kind = logical_and, logical_or#
+ #OP1 = ||, ||#
+ #OP2 = &&, ||#
*/
static void
- at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- register char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((@typ@ *)op)=*((@typ@ *)i1) @OP@ *((@typ@ *)i2);
+ BINARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ const @type@ in2r = ((@type@ *)ip2)[0];
+ const @type@ in2i = ((@type@ *)ip2)[1];
+ *((Bool *)op) = (in1r @OP1@ in1i) @OP2@ (in2r @OP1@ in2i);
}
}
-/**end repeat**/
+/**end repeat1**/
-
-/**begin repeat
- #TYPE=BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG#
- #typ=byte, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong#
-*/
static void
- at TYPE@_invert(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_logical_xor(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0], os=steps[1], n=dimensions[0];
- char *i1=args[0], *op=args[1];
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- *((@typ@ *)op) = ~ *((@typ@*)i1);
+ BINARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ const @type@ in2r = ((@type@ *)ip2)[0];
+ const @type@ in2i = ((@type@ *)ip2)[1];
+ const Bool tmp1 = (in1r || in1i);
+ const Bool tmp2 = (in2r || in2i);
+ *((Bool *)op) = (tmp1 && !tmp2) || (!tmp1 && tmp2);
}
}
-/**end repeat**/
static void
-BOOL_invert(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_logical_not(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0], os=steps[1], n=dimensions[0];
- char *i1=args[0], *op=args[1];
- for(i=0; i<n; i++, i1+=is1, op+=os) {
- *((Bool *)op) = (*((Bool *)i1) ? FALSE : TRUE);
+ UNARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ *((Bool *)op) = !(in1r || in1i);
}
}
-
-/**begin repeat
- #TYPE=BOOL,BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG,FLOAT,DOUBLE,LONGDOUBLE#
- #typ=Bool, byte, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble#
-
-*/
+/**begin repeat1
+ * #kind = isnan, isinf, isfinite#
+ * #func = isnan, isinf, isfinite#
+ * #OP = ||, ||, &&#
+ **/
static void
- at TYPE@_logical_xor(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((Bool *)op)=(*((@typ@ *)i1) || *((@typ@ *)i2)) && !(*((@typ@ *)i1) && *((@typ@ *)i2));
+ UNARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ *((Bool *)op) = @func@(in1r) @OP@ @func@(in1i);
}
}
-/**end repeat**/
+/**end repeat1**/
-
-/**begin repeat
- #TYPE=CFLOAT,CDOUBLE,CLONGDOUBLE#
- #typ=cfloat, cdouble, clongdouble#
-*/
static void
- at TYPE@_logical_xor(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_square(char **args, intp *dimensions, intp *steps, void *data)
{
- Bool p1, p2;
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- p1 = ((@typ@ *)i1)->real || ((@typ@ *)i1)->imag;
- p2 = ((@typ@ *)i2)->real || ((@typ@ *)i2)->imag;
- *((Bool *)op)= (p1 || p2) && !(p1 && p2);
+ UNARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ ((@type@ *)op)[0] = in1r*in1r - in1i*in1i;
+ ((@type@ *)op)[1] = in1r*in1i + in1i*in1r;
}
}
-/**end repeat**/
-
-
-/**begin repeat
-
- #TYPE=(BOOL,BYTE,UBYTE,SHORT,USHORT,INT,UINT,LONG,ULONG,LONGLONG,ULONGLONG,FLOAT,DOUBLE,LONGDOUBLE)*2#
- #OP= >*14, <*14#
- #typ=(Bool, byte, ubyte, short, ushort, int, uint, long, ulong, longlong, ulonglong, float, double, longdouble)*2#
- #kind= maximum*14, minimum*14#
-*/
static void
- at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_reciprocal(char **args, intp *dimensions, intp *steps, void *data)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- *((@typ@ *)op)=*((@typ@ *)i1) @OP@ *((@typ@ *)i2) ? *((@typ@ *)i1) : *((@typ@ *)i2);
+ UNARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ if (fabs at c@(in1i) <= fabs at c@(in1r)) {
+ const @type@ r = in1i/in1r;
+ const @type@ d = in1r + in1i*r;
+ ((@type@ *)op)[0] = 1/d;
+ ((@type@ *)op)[1] = -r/d;
+ } else {
+ const @type@ r = in1r/in1i;
+ const @type@ d = in1r*r + in1i;
+ ((@type@ *)op)[0] = r/d;
+ ((@type@ *)op)[1] = -1/d;
+ }
}
}
-/**end repeat**/
-/**begin repeat
-
- #TYPE=(CFLOAT,CDOUBLE,CLONGDOUBLE)*2#
- #OP= >*3, <*3#
- #typ=(cfloat, cdouble, clongdouble)*2#
- #kind= maximum*3, minimum*3#
-*/
static void
- at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_ones_like(char **args, intp *dimensions, intp *steps, void *data)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2], n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- @typ@ *i1c, *i2c;
- for(i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- i1c = (@typ@ *)i1;
- i2c = (@typ@ *)i2;
- if ((i1c->real @OP@ i2c->real) || \
- ((i1c->real==i2c->real) && (i1c->imag @OP@ i2c->imag)))
- memmove(op, i1, sizeof(@typ@));
- else
- memmove(op, i2, sizeof(@typ@));
+ OUTPUT_LOOP {
+ ((@type@ *)op)[0] = 1;
+ ((@type@ *)op)[1] = 0;
}
}
-/**end repeat**/
+static void
+ at CTYPE@_conjugate(char **args, intp *dimensions, intp *steps, void *func) {
+ UNARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ ((@type@ *)op)[0] = in1r;
+ ((@type@ *)op)[1] = -in1i;
+ }
+}
-
-/*** isinf, isinf, isfinite, signbit ***/
-/**begin repeat
- #kind=isnan*3, isinf*3, isfinite*3, signbit*3#
- #TYPE=(FLOAT, DOUBLE, LONGDOUBLE)*4#
- #typ=(float, double, longdouble)*4#
- #c=(f,,l)*4#
-*/
static void
- at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_absolute(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is=steps[0], os=steps[1], n=dimensions[0];
- char *ip=args[0], *op=args[1];
- for(i=0; i<n; i++, ip+=is, op+=os) {
- *((Bool *)op) = (Bool) (@kind@@c@(*((@typ@ *)ip)) != 0);
+ UNARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ *((@type@ *)op) = sqrt at c@(in1r*in1r + in1i*in1i);
}
}
-/**end repeat**/
-
-/**begin repeat
- #kind=isnan*3, isinf*3, isfinite*3#
- #TYPE=(CFLOAT, CDOUBLE, CLONGDOUBLE)*3#
- #typ=(float, double, longdouble)*3#
- #c=(f,,l)*3#
- #OP=||*6,&&*3#
-*/
static void
- at TYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_sign(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is=steps[0], os=steps[1], n=dimensions[0];
- char *ip=args[0], *op=args[1];
- for(i=0; i<n; i++, ip+=is, op+=os) {
- *((Bool *)op) = @kind@@c@(((@typ@ *)ip)[0]) @OP@ \
- @kind@@c@(((@typ@ *)ip)[1]);
+ UNARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ if (in1r > 0) {
+ ((@type@ *)op)[0] = 1;
+ }
+ else if (in1r < 0) {
+ ((@type@ *)op)[0] = -1;
+ }
+ else {
+ if (in1i > 0) {
+ ((@type@ *)op)[0] = 1;
+ }
+ else if (in1i < 0) {
+ ((@type@ *)op)[0] = -1;
+ }
+ else {
+ ((@type@ *)op)[0] = 0;
+ }
+ }
+ ((@type@ *)op)[1] = 0;
}
}
-/**end repeat**/
-
-
-
-/****** modf ****/
-
-/**begin repeat
- #TYPE=FLOAT, DOUBLE, LONGDOUBLE#
- #typ=float, double, longdouble#
- #c=f,,l#
-*/
+/**begin repeat1
+ * #kind = maximum, minimum#
+ * #OP = >, <#
+ */
static void
- at TYPE@_modf(char **args, intp *dimensions, intp *steps, void *func)
+ at CTYPE@_ at kind@(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],os1=steps[1],os2=steps[2],n=dimensions[0];
- char *i1=args[0], *op1=args[1], *op2=args[2];
- @typ@ x1, y1, y2;
- for (i=0; i<n; i++, i1+=is1, op1+=os1, op2+=os2) {
- x1 = *((@typ@ *)i1);
- y1 = modf at c@(x1, &y2);
- *((@typ@ *)op1) = y1;
- *((@typ@ *)op2) = y2;
+ BINARY_LOOP {
+ const @type@ in1r = ((@type@ *)ip1)[0];
+ const @type@ in1i = ((@type@ *)ip1)[1];
+ const @type@ in2r = ((@type@ *)ip2)[0];
+ const @type@ in2i = ((@type@ *)ip2)[1];
+ if (in1r @OP@ in2r || ((in1r == in2r) && (in1i @OP@ in2i))) {
+ ((@type@ *)op)[0] = in1r;
+ ((@type@ *)op)[1] = in1i;
+ }
+ else {
+ ((@type@ *)op)[0] = in2r;
+ ((@type@ *)op)[1] = in2i;
+ }
}
}
+/**end repeat1**/
+
+#define @CTYPE at _true_divide @CTYPE at _divide
/**end repeat**/
-#define HAVE_DOUBLE_FUNCS
+
+/*
+ *****************************************************************************
+ ** OBJECT LOOPS **
+ *****************************************************************************
+ */
+
/**begin repeat
- #TYPE=FLOAT, DOUBLE, LONGDOUBLE#
- #typ=float, double, longdouble#
- #c=f,,l#
-*/
-#ifdef HAVE_ at TYPE@_FUNCS
+ * #kind = equal, not_equal, greater, greater_equal, less, less_equal#
+ * #OP = EQ, NE, GT, GE, LT, LE#
+ */
static void
- at TYPE@_frexp(char **args, intp *dimensions, intp *steps, void *func)
-{
- register intp i;
- intp is1=steps[0],os1=steps[1],os2=steps[2],n=dimensions[0];
- char *i1=args[0], *op1=args[1], *op2=args[2];
- @typ@ x1, y1;
- int y2;
- for (i=0; i<n; i++, i1+=is1, op1+=os1, op2+=os2) {
- x1 = *((@typ@ *)i1);
- y1 = frexp at c@(x1, &y2);
- *((@typ@ *)op1) = y1;
- *((int *) op2) = y2;
+OBJECT_ at kind@(char **args, intp *dimensions, intp *steps, void *func) {
+ BINARY_LOOP {
+ PyObject *in1 = *(PyObject **)ip1;
+ PyObject *in2 = *(PyObject **)ip2;
+ *(Bool *)op = (Bool) PyObject_RichCompareBool(in1, in2, Py_ at OP@);
}
}
+/**end repeat**/
static void
- at TYPE@_ldexp(char **args, intp *dimensions, intp *steps, void *func)
+OBJECT_sign(char **args, intp *dimensions, intp *steps, void *func)
{
- register intp i;
- intp is1=steps[0],is2=steps[1],os=steps[2],n=dimensions[0];
- char *i1=args[0], *i2=args[1], *op=args[2];
- @typ@ x1, y1;
- int x2;
- for (i=0; i<n; i++, i1+=is1, i2+=is2, op+=os) {
- x1 = *((@typ@ *)i1);
- x2 = *((int *)i2);
- y1 = ldexp at c@(x1, x2);
- *((@typ@ *)op) = y1;
+ PyObject *zero = PyInt_FromLong(0);
+ UNARY_LOOP {
+ PyObject *in1 = *(PyObject **)ip1;
+ *(PyObject **)op = PyInt_FromLong(PyObject_Compare(in1, zero));
}
+ Py_DECREF(zero);
}
-#endif
-/**end repeat**/
-#undef HAVE_DOUBLE_FUNCS
+/*
+ *****************************************************************************
+ ** END LOOPS **
+ *****************************************************************************
+ */
+
static PyUFuncGenericFunction frexp_functions[] = {
#ifdef HAVE_FLOAT_FUNCS
FLOAT_frexp,
@@ -2194,12 +2097,8 @@
};
-
#include "__umath_generated.c"
-
-
#include "ufuncobject.c"
-
#include "__ufunc_api.c"
static double
@@ -2370,6 +2269,7 @@
PyDict_SetItemString(d, "mod", s2);
return;
+
err:
/* Check for errors */
if (!PyErr_Occurred()) {
More information about the Numpy-svn
mailing list