Skip to content

Commit efc51e9

Browse files
committed
fixed quad2quad casting
1 parent 1c339f9 commit efc51e9

2 files changed

Lines changed: 84 additions & 83 deletions

File tree

quaddtype/numpy_quaddtype/src/casts.cpp

Lines changed: 54 additions & 83 deletions
Original file line numberDiff line numberDiff line change
@@ -36,31 +36,39 @@ quad_to_quad_resolve_descriptors(PyObject *NPY_UNUSED(self),
3636
QuadPrecDTypeObject *given_descrs[2],
3737
QuadPrecDTypeObject *loop_descrs[2], npy_intp *view_offset)
3838
{
39-
NPY_CASTING casting = NPY_NO_CASTING;
40-
4139
Py_INCREF(given_descrs[0]);
4240
loop_descrs[0] = given_descrs[0];
4341

4442
if (given_descrs[1] == NULL) {
4543
Py_INCREF(given_descrs[0]);
4644
loop_descrs[1] = given_descrs[0];
45+
*view_offset = 0;
46+
return NPY_NO_CASTING;
4747
}
48-
else {
49-
Py_INCREF(given_descrs[1]);
50-
loop_descrs[1] = given_descrs[1];
51-
if (given_descrs[0]->backend != given_descrs[1]->backend) {
52-
casting = NPY_UNSAFE_CASTING;
48+
49+
Py_INCREF(given_descrs[1]);
50+
loop_descrs[1] = given_descrs[1];
51+
52+
if (given_descrs[0]->backend != given_descrs[1]->backend) {
53+
// Different backends require actual conversion, no view possible
54+
*view_offset = NPY_MIN_INTP;
55+
if (given_descrs[0]->backend == BACKEND_SLEEF) {
56+
// SLEEF -> long double may lose precision
57+
return NPY_SAME_KIND_CASTING;
5358
}
59+
// long double -> SLEEF preserves value exactly
60+
return NPY_SAFE_CASTING;
5461
}
5562

5663
*view_offset = 0;
57-
return casting;
64+
return NPY_NO_CASTING;
5865
}
5966

67+
template <bool Aligned>
6068
static int
61-
quad_to_quad_strided_loop_unaligned(PyArrayMethod_Context *context, char *const data[],
62-
npy_intp const dimensions[], npy_intp const strides[],
63-
void *NPY_UNUSED(auxdata))
69+
quad_to_quad_strided_loop(PyArrayMethod_Context *context, char *const data[],
70+
npy_intp const dimensions[], npy_intp const strides[],
71+
void *NPY_UNUSED(auxdata))
6472
{
6573
npy_intp N = dimensions[0];
6674
char *in_ptr = data[0];
@@ -70,93 +78,56 @@ quad_to_quad_strided_loop_unaligned(PyArrayMethod_Context *context, char *const
7078

7179
QuadPrecDTypeObject *descr_in = (QuadPrecDTypeObject *)context->descriptors[0];
7280
QuadPrecDTypeObject *descr_out = (QuadPrecDTypeObject *)context->descriptors[1];
81+
QuadBackendType backend_in = descr_in->backend;
82+
QuadBackendType backend_out = descr_out->backend;
7383

7484
// inter-backend casting
75-
if (descr_in->backend != descr_out->backend) {
85+
if (backend_in != backend_out) {
7686
while (N--) {
77-
quad_value in_val, out_val;
78-
if (descr_in->backend == BACKEND_SLEEF) {
79-
memcpy(&in_val.sleef_value, in_ptr, sizeof(Sleef_quad));
80-
out_val.longdouble_value = cast_sleef_to_double(in_val.sleef_value);
87+
quad_value in_val;
88+
load_quad<Aligned>(in_ptr, backend_in, &in_val);
89+
quad_value out_val;
90+
if (backend_in == BACKEND_SLEEF)
91+
{
92+
out_val.longdouble_value = static_cast<long double>(cast_sleef_to_double(in_val.sleef_value));
8193
}
82-
else {
83-
memcpy(&in_val.longdouble_value, in_ptr, sizeof(long double));
84-
out_val.sleef_value = Sleef_cast_from_doubleq1(in_val.longdouble_value);
94+
else
95+
{
96+
long double ld = in_val.longdouble_value;
97+
if (std::isnan(ld)) {
98+
out_val.sleef_value = (!ld_signbit(&ld)) ? QUAD_PRECISION_NAN : QUAD_PRECISION_NEG_NAN;
99+
}
100+
else if (std::isinf(ld)) {
101+
out_val.sleef_value = (ld > 0) ? QUAD_PRECISION_INF : QUAD_PRECISION_NINF;
102+
}
103+
else
104+
{
105+
// to prevent compiler optimizations, ABI handling issues with __float128 on x86-64 machines
106+
// won't be expensive as for fixed size compiler can optimize memcpy with movq
107+
Sleef_quad temp = Sleef_cast_from_doubleq1(static_cast<double>(ld));
108+
std::memcpy(&out_val.sleef_value, &temp, sizeof(Sleef_quad));
109+
}
85110
}
86-
memcpy(out_ptr, &out_val,
87-
(descr_out->backend == BACKEND_SLEEF) ? sizeof(Sleef_quad)
88-
: sizeof(long double));
111+
112+
store_quad<Aligned>(out_ptr, &out_val, backend_out);
89113
in_ptr += in_stride;
90114
out_ptr += out_stride;
91115
}
92-
93116
return 0;
94117
}
95118

96-
size_t elem_size =
97-
(descr_in->backend == BACKEND_SLEEF) ? sizeof(Sleef_quad) : sizeof(long double);
98-
99-
while (N--) {
100-
memcpy(out_ptr, in_ptr, elem_size);
119+
// same backend: direct copy
120+
// same_value casting not needed here as values are identical
121+
while(N--) {
122+
quad_value val;
123+
load_quad<Aligned>(in_ptr, backend_in, &val);
124+
store_quad<Aligned>(out_ptr, &val, backend_out);
101125
in_ptr += in_stride;
102126
out_ptr += out_stride;
103127
}
104128
return 0;
105129
}
106130

107-
static int
108-
quad_to_quad_strided_loop_aligned(PyArrayMethod_Context *context, char *const data[],
109-
npy_intp const dimensions[], npy_intp const strides[],
110-
void *NPY_UNUSED(auxdata))
111-
{
112-
npy_intp N = dimensions[0];
113-
char *in_ptr = data[0];
114-
char *out_ptr = data[1];
115-
npy_intp in_stride = strides[0];
116-
npy_intp out_stride = strides[1];
117-
118-
QuadPrecDTypeObject *descr_in = (QuadPrecDTypeObject *)context->descriptors[0];
119-
QuadPrecDTypeObject *descr_out = (QuadPrecDTypeObject *)context->descriptors[1];
120-
121-
// inter-backend casting
122-
if (descr_in->backend != descr_out->backend) {
123-
if (descr_in->backend == BACKEND_SLEEF) {
124-
while (N--) {
125-
Sleef_quad in_val = *(Sleef_quad *)in_ptr;
126-
*(long double *)out_ptr = cast_sleef_to_double(in_val);
127-
in_ptr += in_stride;
128-
out_ptr += out_stride;
129-
}
130-
}
131-
else {
132-
while (N--) {
133-
long double in_val = *(long double *)in_ptr;
134-
*(Sleef_quad *)out_ptr = Sleef_cast_from_doubleq1(in_val);
135-
in_ptr += in_stride;
136-
out_ptr += out_stride;
137-
}
138-
}
139-
140-
return 0;
141-
}
142-
143-
if (descr_in->backend == BACKEND_SLEEF) {
144-
while (N--) {
145-
*(Sleef_quad *)out_ptr = *(Sleef_quad *)in_ptr;
146-
in_ptr += in_stride;
147-
out_ptr += out_stride;
148-
}
149-
}
150-
else {
151-
while (N--) {
152-
*(long double *)out_ptr = *(long double *)in_ptr;
153-
in_ptr += in_stride;
154-
out_ptr += out_stride;
155-
}
156-
}
157-
158-
return 0;
159-
}
160131

161132
static NPY_CASTING
162133
void_to_quad_resolve_descriptors(PyObject *NPY_UNUSED(self), PyArray_DTypeMeta *dtypes[2],
@@ -1443,8 +1414,8 @@ init_casts_internal(void)
14431414
PyArray_DTypeMeta **quad2quad_dtypes = new PyArray_DTypeMeta *[2]{nullptr, nullptr};
14441415
PyType_Slot *quad2quad_slots = new PyType_Slot[4]{
14451416
{NPY_METH_resolve_descriptors, (void *)&quad_to_quad_resolve_descriptors},
1446-
{NPY_METH_strided_loop, (void *)&quad_to_quad_strided_loop_aligned},
1447-
{NPY_METH_unaligned_strided_loop, (void *)&quad_to_quad_strided_loop_unaligned},
1417+
{NPY_METH_strided_loop, (void *)&quad_to_quad_strided_loop<true>},
1418+
{NPY_METH_unaligned_strided_loop, (void *)&quad_to_quad_strided_loop<false>},
14481419
{0, nullptr}};
14491420

14501421
PyArrayMethod_Spec *quad2quad_spec = new PyArrayMethod_Spec{

quaddtype/tests/test_quaddtype.py

Lines changed: 30 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -5368,3 +5368,33 @@ def test_hash_backends(self, backend):
53685368
"""Test hash works for both backends."""
53695369
quad_val = QuadPrecision(1.5, backend=backend)
53705370
assert hash(quad_val) == hash(1.5)
5371+
5372+
@pytest.mark.parametrize("src_backend,dst_backend", [
5373+
("sleef", "longdouble"),
5374+
("longdouble", "sleef"),
5375+
("sleef", "sleef"),
5376+
("longdouble", "longdouble"),
5377+
])
5378+
@pytest.mark.parametrize("value", [
5379+
"0.0", "-0.0", "1.0", "-1.0", "3.14159265358979323846",
5380+
"inf", "-inf", "nan", "1e100", "1e-100", "-nan"
5381+
])
5382+
def test_quad_to_quad_backend_casting(src_backend, dst_backend, value):
5383+
"""Test casting between QuadPrecDType with different backends."""
5384+
5385+
src_arr = np.array([value], dtype=QuadPrecDType(backend=src_backend))
5386+
dst_arr = src_arr.astype(QuadPrecDType(backend=dst_backend))
5387+
res_arr = np.array([value], dtype=QuadPrecDType(backend=dst_backend))
5388+
5389+
expected_backend = 0 if dst_backend == 'sleef' else 1
5390+
assert dst_arr.dtype.backend == expected_backend
5391+
5392+
assert np.signbit(src_arr[0]) == np.signbit(dst_arr[0])
5393+
if np.isnan(src_arr[0]):
5394+
assert np.isnan(dst_arr[0])
5395+
elif np.isinf(src_arr[0]):
5396+
assert np.isinf(dst_arr[0])
5397+
elif src_backend != dst_backend:
5398+
np.testing.assert_allclose(dst_arr, res_arr, rtol=1e-15)
5399+
else:
5400+
np.testing.assert_array_equal(dst_arr, res_arr)

0 commit comments

Comments
 (0)