diff --git a/quaddtype/.gitignore b/.gitignore similarity index 100% rename from quaddtype/.gitignore rename to .gitignore diff --git a/quaddtype/pyproject.toml b/quaddtype/pyproject.toml index 99b5d9c3..687b2637 100644 --- a/quaddtype/pyproject.toml +++ b/quaddtype/pyproject.toml @@ -4,7 +4,7 @@ requires = [ "meson-python", "patchelf", "wheel", - "numpy @ file:///home/pdmurray/Desktop/numpy-1.24.0.dev0+1120.gf30af6acd-cp311-cp311-manylinux_2_17_x86_64.manylinux2014_x86_64.whl" + "numpy @ file:///home/pdmurray/Desktop/numpy-1.25.0.dev0+167.g4ca8204c5-cp311-cp311d-linux_x86_64.whl" ] build-backend = "mesonpy" @@ -16,7 +16,7 @@ readme = 'README.md' author = "Peyton Murray" requires-python = ">=3.9.0" dependencies = [ - "numpy @ file:///home/pdmurray/Desktop/numpy-1.24.0.dev0+1120.gf30af6acd-cp311-cp311-manylinux_2_17_x86_64.manylinux2014_x86_64.whl" + "numpy @ file:///home/pdmurray/Desktop/numpy-1.25.0.dev0+167.g4ca8204c5-cp311-cp311d-linux_x86_64.whl" ] [project.optional-dependencies] diff --git a/quaddtype/quaddtype/src/casts.c b/quaddtype/quaddtype/src/casts.c index 3df6529f..bc0ee4c2 100644 --- a/quaddtype/quaddtype/src/casts.c +++ b/quaddtype/quaddtype/src/casts.c @@ -29,24 +29,69 @@ quad_to_quad_resolve_descriptors(PyObject *NPY_UNUSED(self), return NPY_SAME_KIND_CASTING; } +// Each element is a __float128 element; no casting needed static int quad_to_quad_contiguous(PyArrayMethod_Context *NPY_UNUSED(context), char *const data[], - npy_intp const dimensions[], npy_intp const strides[], void *auxdata) + npy_intp const dimensions[], npy_intp const strides[], + NpyAuxData *NPY_UNUSED(auxdata)) { + npy_intp N = dimensions[0]; + __float128 *in = (__float128 *)data[0]; + __float128 *out = (__float128 *)data[1]; + + while (N--) { + *out = *in; + out++; + in++; + } + return 0; } +// Elements are strided, e.g. +// +// x = np.linspace(40) +// x[::3] +// +// Therefore the stride needs to be used to increment the pointers inside the loop. static int quad_to_quad_strided(PyArrayMethod_Context *NPY_UNUSED(context), char *const data[], - npy_intp const dimensions[], npy_intp const strides[], void *auxdata) + npy_intp const dimensions[], npy_intp const strides[], + NpyAuxData *NPY_UNUSED(auxdata)) { + npy_intp N = dimensions[0]; + char *in = data[0]; + char *out = data[1]; + npy_intp in_stride = strides[0]; + npy_intp out_stride = strides[1]; + + while (N--) { + *(__float128 *)out = *(__float128 *)in; + in += in_stride; + out += out_stride; + } + return 0; } +// Arrays are unaligned. static int quad_to_quad_unaligned(PyArrayMethod_Context *NPY_UNUSED(context), char *const data[], - npy_intp const dimensions[], npy_intp const strides[], void *auxdata) + npy_intp const dimensions[], npy_intp const strides[], + NpyAuxData *NPY_UNUSED(auxdata)) { + npy_intp N = dimensions[0]; + char *in = data[0]; + char *out = data[1]; + npy_intp in_stride = strides[0]; + npy_intp out_stride = strides[1]; + + while (N--) { + memcpy(out, in, sizeof(__float128)); // NOLINT + in += in_stride; + out += out_stride; + } + return 0; } diff --git a/quaddtype/quaddtype/src/dtype.c b/quaddtype/quaddtype/src/dtype.c index b70da987..470245be 100644 --- a/quaddtype/quaddtype/src/dtype.c +++ b/quaddtype/quaddtype/src/dtype.c @@ -1,14 +1,5 @@ -#include - -#define PY_ARRAY_UNIQUE_SYMBOL quaddtype_ARRAY_API -#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION -#define NO_IMPORT_ARRAY -#include "numpy/arrayobject.h" -#include "numpy/experimental_dtype_api.h" -#include "numpy/ndarraytypes.h" - -#include "casts.h" #include "dtype.h" +#include "casts.h" PyTypeObject *QuadScalar_Type = NULL; @@ -17,8 +8,9 @@ new_quaddtype_instance(void) { QuadDTypeObject *new = (QuadDTypeObject *)PyArrayDescr_Type.tp_new((PyTypeObject *)&QuadDType, NULL, NULL); - if (new == NULL) + if (new == NULL) { return NULL; + } new->base.elsize = sizeof(__float128); new->base.alignment = _Alignof(__float128); @@ -82,14 +74,34 @@ common_dtype(PyArray_DTypeMeta *self, PyArray_DTypeMeta *other) return (PyArray_DTypeMeta *)Py_NotImplemented; } -static PyType_Slot QuadDType_Slots[] = {{NPY_DT_common_instance, &common_instance}, - {NPY_DT_common_dtype, &common_dtype}, - // {NPY_DT_discover_descr_from_pyobject, - // &unit_discover_descriptor_from_pyobject}, - /* The header is wrong on main :(, so we add 1 */ - {NPY_DT_setitem, &quad_setitem}, - {NPY_DT_getitem, &quad_getitem}, - {0, NULL}}; +static QuadDTypeObject * +quaddtype_ensure_canonical(QuadDTypeObject *self) +{ + Py_INCREF(self); + return self; +} + +static PyArray_Descr * +quad_discover_descriptor_from_pyobject(PyArray_DTypeMeta *NPY_UNUSED(cls), PyObject *obj) +{ + if (Py_TYPE(obj) != QuadScalar_Type) { + PyErr_SetString(PyExc_TypeError, "Can only store QuadScalars in a QuadDType array."); + return NULL; + } + + // Get the dtype attribute from the object. + return (PyArray_Descr *)PyObject_GetAttrString(obj, "dtype"); +} + +static PyType_Slot QuadDType_Slots[] = { + {NPY_DT_common_instance, &common_instance}, + {NPY_DT_common_dtype, &common_dtype}, + {NPY_DT_discover_descr_from_pyobject, &quad_discover_descriptor_from_pyobject}, + /* The header is wrong on main :(, so we add 1 */ + {NPY_DT_setitem, &quad_setitem}, + {NPY_DT_getitem, &quad_getitem}, + {NPY_DT_ensure_canonical, &quaddtype_ensure_canonical}, + {0, NULL}}; /* * The following defines everything type object related (i.e. not NumPy @@ -114,7 +126,8 @@ quaddtype_dealloc(QuadDTypeObject *self) static PyObject * quaddtype_repr(QuadDTypeObject *self) { - return PyUnicode_FromString("This is a quad (128-bit float) dtype."); + PyObject *res = PyUnicode_FromString("This is a quad (128-bit float) dtype."); + return res; } // These are the basic things that you need to create a Python Type/Class in C. @@ -140,7 +153,6 @@ init_quad_dtype(void) // do it. You first have to create a static type, but see the note there! PyArrayMethod_Spec *casts[] = { &QuadToQuadCastSpec, - &QuadToFloat128CastSpec, NULL, }; diff --git a/quaddtype/quaddtype/src/dtype.h b/quaddtype/quaddtype/src/dtype.h index 2427beb5..fdce9398 100644 --- a/quaddtype/quaddtype/src/dtype.h +++ b/quaddtype/quaddtype/src/dtype.h @@ -3,6 +3,10 @@ #include +#define PY_ARRAY_UNIQUE_SYMBOL quaddtype_ARRAY_API +#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION +#define NO_IMPORT_ARRAY +#include "numpy/ndarraytypes.h" #include "numpy/arrayobject.h" #include "numpy/experimental_dtype_api.h" @@ -15,6 +19,7 @@ extern PyTypeObject *QuadScalar_Type; QuadDTypeObject * new_quaddtype_instance(void); + int init_quad_dtype(void); diff --git a/quaddtype/quaddtype/src/umath.c b/quaddtype/quaddtype/src/umath.c index f8fed77c..86c65171 100644 --- a/quaddtype/quaddtype/src/umath.c +++ b/quaddtype/quaddtype/src/umath.c @@ -25,7 +25,7 @@ quad_multiply_strided_loop(PyArrayMethod_Context *context, char *const data[], char *in1 = data[0], *in2 = data[1]; char *out = data[2]; npy_intp in1_stride = strides[0]; - npy_intp in2_stride = strides[0]; + npy_intp in2_stride = strides[1]; npy_intp out_stride = strides[2]; while (N--) { @@ -56,9 +56,18 @@ quad_multiply_resolve_descriptors(PyObject *self, PyArray_DTypeMeta *dtypes[], // The operand units can be used as-is; no casting required for quad types. Py_INCREF(given_descrs[0]); loop_descrs[0] = given_descrs[0]; + + if (given_descrs[1] == NULL) { + Py_INCREF(given_descrs[0]); + loop_descrs[1] = given_descrs[0]; + return NPY_NO_CASTING; + } + Py_INCREF(given_descrs[1]); + loop_descrs[1] = given_descrs[1]; + Py_INCREF(given_descrs[1]); loop_descrs[1] = given_descrs[1]; - return NPY_NO_CASTING; + return NPY_SAFE_CASTING; } // Function that adds our multiply loop to NumPy's multiply ufunc.