pythongh-61103: Support double complex (_Complex) type in ctypes (python#120894)

Example:

```pycon
>>> import ctypes
>>> ctypes.__STDC_IEC_559_COMPLEX__
1
>>> libm = ctypes.CDLL('libm.so.6')
>>> libm.clog.argtypes = [ctypes.c_double_complex]
>>> libm.clog.restype = ctypes.c_double_complex
>>> libm.clog(1+1j)
(0.34657359027997264+0.7853981633974483j)
```

Co-authored-by: Nice Zombies <[email protected]>
Co-authored-by: Bénédikt Tran <[email protected]>
Co-authored-by: Victor Stinner <[email protected]>

Author: 4 people

Doc/library/ctypes.rst

@@ -266,6 +266,16 @@ Fundamental data types
 (1)
    The constructor accepts any object with a truth value.
 
+Additionally, if IEC 60559 compatible complex arithmetic (Annex G) is supported, the following
+complex types are available:
+
++----------------------------------+---------------------------------+-----------------+
+| ctypes type                      | C type                          | Python type     |
++==================================+=================================+=================+
+| :class:`c_double_complex`        | :c:expr:`double complex`        | complex         |
++----------------------------------+---------------------------------+-----------------+
+
+
 All these types can be created by calling them with an optional initializer of
 the correct type and value::
 
@@ -2284,6 +2294,14 @@ These are the fundamental ctypes data types:
    optional float initializer.
 
 
+.. class:: c_double_complex
+
+   Represents the C :c:expr:`double complex` datatype, if available.  The
+   constructor accepts an optional :class:`complex` initializer.
+
+   .. versionadded:: 3.14
+
+
 .. class:: c_int
 
    Represents the C :c:expr:`signed int` datatype.  The constructor accepts an

Lib/ctypes/__init__.py

@@ -205,6 +205,12 @@ class c_longdouble(_SimpleCData):
 if sizeof(c_longdouble) == sizeof(c_double):
     c_longdouble = c_double
 
+try:
+    class c_double_complex(_SimpleCData):
+        _type_ = "C"
+except AttributeError:
+    pass
+
 if _calcsize("l") == _calcsize("q"):
     # if long and long long have the same size, make c_longlong an alias for c_long
     c_longlong = c_long

Lib/test/test_ctypes/test_libc.py

@@ -1,3 +1,4 @@
+import ctypes
 import math
 import unittest
 from ctypes import (CDLL, CFUNCTYPE, POINTER, create_string_buffer, sizeof,
@@ -21,6 +22,17 @@ def test_sqrt(self):
         self.assertEqual(lib.my_sqrt(4.0), 2.0)
         self.assertEqual(lib.my_sqrt(2.0), math.sqrt(2.0))
 
+    @unittest.skipUnless(hasattr(ctypes, "c_double_complex"),
+                         "requires C11 complex type")
+    def test_csqrt(self):
+        lib.my_csqrt.argtypes = ctypes.c_double_complex,
+        lib.my_csqrt.restype = ctypes.c_double_complex
+        self.assertEqual(lib.my_csqrt(4), 2+0j)
+        self.assertAlmostEqual(lib.my_csqrt(-1+0.01j),
+                               0.004999937502734214+1.0000124996093955j)
+        self.assertAlmostEqual(lib.my_csqrt(-1-0.01j),
+                               0.004999937502734214-1.0000124996093955j)
+
     def test_qsort(self):
         comparefunc = CFUNCTYPE(c_int, POINTER(c_char), POINTER(c_char))
         lib.my_qsort.argtypes = c_void_p, c_size_t, c_size_t, comparefunc

Lib/test/test_ctypes/test_numbers.py

@@ -1,7 +1,10 @@
 import array
+import ctypes
 import struct
 import sys
 import unittest
+from itertools import combinations
+from math import copysign, isnan
 from operator import truth
 from ctypes import (byref, sizeof, alignment,
                     c_char, c_byte, c_ubyte, c_short, c_ushort, c_int, c_uint,
@@ -38,8 +41,55 @@ def valid_ranges(*types):
 signed_ranges = valid_ranges(*signed_types)
 bool_values = [True, False, 0, 1, -1, 5000, 'test', [], [1]]
 
+class IntLike:
+    def __int__(self):
+        return 2
+
+class IndexLike:
+    def __index__(self):
+        return 2
+
+class FloatLike:
+    def __float__(self):
+        return 2.0
+
+class ComplexLike:
+    def __complex__(self):
+        return 1+1j
+
+
+INF = float("inf")
+NAN = float("nan")
+
 
 class NumberTestCase(unittest.TestCase):
+    # from Lib/test/test_complex.py
+    def assertFloatsAreIdentical(self, x, y):
+        """assert that floats x and y are identical, in the sense that:
+        (1) both x and y are nans, or
+        (2) both x and y are infinities, with the same sign, or
+        (3) both x and y are zeros, with the same sign, or
+        (4) x and y are both finite and nonzero, and x == y
+
+        """
+        msg = 'floats {!r} and {!r} are not identical'
+
+        if isnan(x) or isnan(y):
+            if isnan(x) and isnan(y):
+                return
+        elif x == y:
+            if x != 0.0:
+                return
+            # both zero; check that signs match
+            elif copysign(1.0, x) == copysign(1.0, y):
+                return
+            else:
+                msg += ': zeros have different signs'
+        self.fail(msg.format(x, y))
+
+    def assertComplexesAreIdentical(self, x, y):
+        self.assertFloatsAreIdentical(x.real, y.real)
+        self.assertFloatsAreIdentical(x.imag, y.imag)
 
     def test_default_init(self):
         # default values are set to zero
@@ -86,28 +136,39 @@ def test_byref(self):
     def test_floats(self):
         # c_float and c_double can be created from
         # Python int and float
-        class FloatLike:
-            def __float__(self):
-                return 2.0
         f = FloatLike()
         for t in float_types:
             self.assertEqual(t(2.0).value, 2.0)
             self.assertEqual(t(2).value, 2.0)
             self.assertEqual(t(2).value, 2.0)
             self.assertEqual(t(f).value, 2.0)
 
+    @unittest.skipUnless(hasattr(ctypes, "c_double_complex"),
+                         "requires C11 complex type")
+    def test_complex(self):
+        for t in [ctypes.c_double_complex]:
+            self.assertEqual(t(1).value, 1+0j)
+            self.assertEqual(t(1.0).value, 1+0j)
+            self.assertEqual(t(1+0.125j).value, 1+0.125j)
+            self.assertEqual(t(IndexLike()).value, 2+0j)
+            self.assertEqual(t(FloatLike()).value, 2+0j)
+            self.assertEqual(t(ComplexLike()).value, 1+1j)
+
+    @unittest.skipUnless(hasattr(ctypes, "c_double_complex"),
+                         "requires C11 complex type")
+    def test_complex_round_trip(self):
+        # Ensure complexes transformed exactly.  The CMPLX macro should
+        # preserve special components (like inf/nan or signed zero).
+        values = [complex(*_) for _ in combinations([1, -1, 0.0, -0.0, 2,
+                                                     -3, INF, -INF, NAN], 2)]
+        for z in values:
+            with self.subTest(z=z):
+                z2 = ctypes.c_double_complex(z).value
+                self.assertComplexesAreIdentical(z, z2)
+
     def test_integers(self):
-        class FloatLike:
-            def __float__(self):
-                return 2.0
         f = FloatLike()
-        class IntLike:
-            def __int__(self):
-                return 2
         d = IntLike()
-        class IndexLike:
-            def __index__(self):
-                return 2
         i = IndexLike()
         # integers cannot be constructed from floats,
         # but from integer-like objects

Makefile.pre.in

@@ -3125,7 +3125,7 @@ MODULE_MATH_DEPS=$(srcdir)/Modules/_math.h
 MODULE_PYEXPAT_DEPS=@LIBEXPAT_INTERNAL@
 MODULE_UNICODEDATA_DEPS=$(srcdir)/Modules/unicodedata_db.h $(srcdir)/Modules/unicodename_db.h
 MODULE__BLAKE2_DEPS=$(srcdir)/Modules/_blake2/impl/blake2-config.h $(srcdir)/Modules/_blake2/impl/blake2-impl.h $(srcdir)/Modules/_blake2/impl/blake2.h $(srcdir)/Modules/_blake2/impl/blake2b-load-sse2.h $(srcdir)/Modules/_blake2/impl/blake2b-load-sse41.h $(srcdir)/Modules/_blake2/impl/blake2b-ref.c $(srcdir)/Modules/_blake2/impl/blake2b-round.h $(srcdir)/Modules/_blake2/impl/blake2b.c $(srcdir)/Modules/_blake2/impl/blake2s-load-sse2.h $(srcdir)/Modules/_blake2/impl/blake2s-load-sse41.h $(srcdir)/Modules/_blake2/impl/blake2s-load-xop.h $(srcdir)/Modules/_blake2/impl/blake2s-ref.c $(srcdir)/Modules/_blake2/impl/blake2s-round.h $(srcdir)/Modules/_blake2/impl/blake2s.c $(srcdir)/Modules/_blake2/blake2module.h $(srcdir)/Modules/hashlib.h
-MODULE__CTYPES_DEPS=$(srcdir)/Modules/_ctypes/ctypes.h
+MODULE__CTYPES_DEPS=$(srcdir)/Modules/_ctypes/ctypes.h $(srcdir)/Modules/_complex.h
 MODULE__CTYPES_TEST_DEPS=$(srcdir)/Modules/_ctypes/_ctypes_test_generated.c.h
 MODULE__CTYPES_MALLOC_CLOSURE=@MODULE__CTYPES_MALLOC_CLOSURE@
 MODULE__DECIMAL_DEPS=$(srcdir)/Modules/_decimal/docstrings.h @LIBMPDEC_INTERNAL@

Misc/NEWS.d/next/Library/2024-06-23-07-23-08.gh-issue-61103.ca_U_l.rst

@@ -0,0 +1,3 @@
+Support :c:expr:`double complex` C type in :mod:`ctypes` via
+:class:`~ctypes.c_double_complex` if compiler has C11 complex
+arithmetic.  Patch by Sergey B Kirpichev.

Modules/_complex.h

@@ -0,0 +1,34 @@
+/* Workarounds for buggy complex number arithmetic implementations. */
+
+#ifndef Py_HAVE_C_COMPLEX
+#  error "this header file should only be included if Py_HAVE_C_COMPLEX is defined"
+#endif
+
+#include <complex.h>
+
+/* Other compilers (than clang), that claims to
+   implement C11 *and* define __STDC_IEC_559_COMPLEX__ - don't have
+   issue with CMPLX().  This is specific to glibc & clang combination:
+   https://sourceware.org/bugzilla/show_bug.cgi?id=26287
+
+   Here we fallback to using __builtin_complex(), available in clang
+   v12+.  Else CMPLX implemented following C11 6.2.5p13: "Each complex type
+   has the same representation and alignment requirements as an array
+   type containing exactly two elements of the corresponding real type;
+   the first element is equal to the real part, and the second element
+   to the imaginary part, of the complex number.
+ */
+#if !defined(CMPLX)
+#  if defined(__clang__) && __has_builtin(__builtin_complex)
+#    define CMPLX(x, y) __builtin_complex ((double) (x), (double) (y))
+#  else
+static inline double complex
+CMPLX(double real, double imag)
+{
+    double complex z;
+    ((double *)(&z))[0] = real;
+    ((double *)(&z))[1] = imag;
+    return z;
+}
+#  endif
+#endif

Modules/_ctypes/_ctypes.c

@@ -1750,7 +1750,11 @@ class _ctypes.c_void_p "PyObject *" "clinic_state_sub()->PyCSimpleType_Type"
 [clinic start generated code]*/
 /*[clinic end generated code: output=da39a3ee5e6b4b0d input=dd4d9646c56f43a9]*/
 
+#if defined(Py_HAVE_C_COMPLEX) && defined(FFI_TARGET_HAS_COMPLEX_TYPE)
+static const char SIMPLE_TYPE_CHARS[] = "cbBhHiIlLdCfuzZqQPXOv?g";
+#else
 static const char SIMPLE_TYPE_CHARS[] = "cbBhHiIlLdfuzZqQPXOv?g";
+#endif
 
 /*[clinic input]
 _ctypes.c_wchar_p.from_param as c_wchar_p_from_param
@@ -2226,7 +2230,17 @@ PyCSimpleType_init(PyObject *self, PyObject *args, PyObject *kwds)
         goto error;
     }
 
-    stginfo->ffi_type_pointer = *fmt->pffi_type;
+    if (!fmt->pffi_type->elements) {
+        stginfo->ffi_type_pointer = *fmt->pffi_type;
+    }
+    else {
+        stginfo->ffi_type_pointer.size = fmt->pffi_type->size;
+        stginfo->ffi_type_pointer.alignment = fmt->pffi_type->alignment;
+        stginfo->ffi_type_pointer.type = fmt->pffi_type->type;
+        stginfo->ffi_type_pointer.elements = PyMem_Malloc(2 * sizeof(ffi_type));
+        memcpy(stginfo->ffi_type_pointer.elements,
+               fmt->pffi_type->elements, 2 * sizeof(ffi_type));
+    }
     stginfo->align = fmt->pffi_type->alignment;
     stginfo->length = 0;
     stginfo->size = fmt->pffi_type->size;

Modules/_ctypes/_ctypes_test.c

@@ -13,6 +13,12 @@
 
 #include <Python.h>
 
+#include <ffi.h>                  // FFI_TARGET_HAS_COMPLEX_TYPE
+
+#if defined(Py_HAVE_C_COMPLEX) && defined(FFI_TARGET_HAS_COMPLEX_TYPE)
+#  include "../_complex.h"        // csqrt()
+#  undef I                        // for _ctypes_test_generated.c.h
+#endif
 #include <stdio.h>                // printf()
 #include <stdlib.h>               // qsort()
 #include <string.h>               // memset()
@@ -443,6 +449,13 @@ EXPORT(double) my_sqrt(double a)
     return sqrt(a);
 }
 
+#if defined(Py_HAVE_C_COMPLEX) && defined(FFI_TARGET_HAS_COMPLEX_TYPE)
+EXPORT(double complex) my_csqrt(double complex a)
+{
+    return csqrt(a);
+}
+#endif
+
 EXPORT(void) my_qsort(void *base, size_t num, size_t width, int(*compare)(const void*, const void*))
 {
     qsort(base, num, width, compare);

Modules/_ctypes/callproc.c

@@ -105,6 +105,10 @@ module _ctypes
 #include "pycore_global_objects.h"// _Py_ID()
 #include "pycore_traceback.h"     // _PyTraceback_Add()
 
+#if defined(Py_HAVE_C_COMPLEX) && defined(FFI_TARGET_HAS_COMPLEX_TYPE)
+#include "../_complex.h"          // complex
+#endif
+
 #include "clinic/callproc.c.h"
 
 #define CTYPES_CAPSULE_NAME_PYMEM "_ctypes pymem"
@@ -651,6 +655,9 @@ union result {
     double d;
     float f;
     void *p;
+#if defined(Py_HAVE_C_COMPLEX) && defined(FFI_TARGET_HAS_COMPLEX_TYPE)
+    double complex C;
+#endif
 };
 
 struct argument {

Modules/_ctypes/cfield.c

@@ -14,6 +14,9 @@
 #include <ffi.h>
 #include "ctypes.h"
 
+#if defined(Py_HAVE_C_COMPLEX) && defined(FFI_TARGET_HAS_COMPLEX_TYPE)
+#  include "../_complex.h"        // complex
+#endif
 
 #define CTYPES_CFIELD_CAPSULE_NAME_PYMEM "_ctypes/cfield.c pymem"
 
@@ -1087,6 +1090,30 @@ d_get(void *ptr, Py_ssize_t size)
     return PyFloat_FromDouble(val);
 }
 
+#if defined(Py_HAVE_C_COMPLEX) && defined(FFI_TARGET_HAS_COMPLEX_TYPE)
+static PyObject *
+C_set(void *ptr, PyObject *value, Py_ssize_t size)
+{
+    Py_complex c = PyComplex_AsCComplex(value);
+
+    if (c.real == -1 && PyErr_Occurred()) {
+        return NULL;
+    }
+    double complex x = CMPLX(c.real, c.imag);
+    memcpy(ptr, &x, sizeof(x));
+    _RET(value);
+}
+
+static PyObject *
+C_get(void *ptr, Py_ssize_t size)
+{
+    double complex x;
+
+    memcpy(&x, ptr, sizeof(x));
+    return PyComplex_FromDoubles(creal(x), cimag(x));
+}
+#endif
+
 static PyObject *
 d_set_sw(void *ptr, PyObject *value, Py_ssize_t size)
 {
@@ -1592,6 +1619,9 @@ static struct fielddesc formattable[] = {
     { 'B', B_set, B_get, NULL},
     { 'c', c_set, c_get, NULL},
     { 'd', d_set, d_get, NULL, d_set_sw, d_get_sw},
+#if defined(Py_HAVE_C_COMPLEX) && defined(FFI_TARGET_HAS_COMPLEX_TYPE)
+    { 'C', C_set, C_get, NULL},
+#endif
     { 'g', g_set, g_get, NULL},
     { 'f', f_set, f_get, NULL, f_set_sw, f_get_sw},
     { 'h', h_set, h_get, NULL, h_set_sw, h_get_sw},
@@ -1642,6 +1672,9 @@ _ctypes_init_fielddesc(void)
         case 'B': fd->pffi_type = &ffi_type_uchar; break;
         case 'c': fd->pffi_type = &ffi_type_schar; break;
         case 'd': fd->pffi_type = &ffi_type_double; break;
+#if defined(Py_HAVE_C_COMPLEX) && defined(FFI_TARGET_HAS_COMPLEX_TYPE)
+        case 'C': fd->pffi_type = &ffi_type_complex_double; break;
+#endif
         case 'g': fd->pffi_type = &ffi_type_longdouble; break;
         case 'f': fd->pffi_type = &ffi_type_float; break;
         case 'h': fd->pffi_type = &ffi_type_sshort; break;