array-api icon indicating copy to clipboard operation
array-api copied to clipboard

RFC: Adding complex number support to the specification

Open kgryte opened this issue 3 years ago • 3 comments

Complex Number Support

Plan for complex number support in the array API specification.

What follows is a plan for adding complex number support to the 2022 array API specification. This RFC is comprised of the following sections:

  • Prior discusssions
  • High-level summary concerns and questions
  • Individual API changes

Prior Discussions

General Concerns

New APIs

abs2 (x*x^H)
real
conj/conjugate (complex conjugate)
imag
arg/angle/phase (phase angle)
linalg.eig
linalg.eigvals

Decisions

  • initial data types: complex64 and complex128?

  • rounding complex numbers (component-wise)

    • NumPy supports round, but not ceil, floor, and trunc
  • ordering (lexicographic?)

  • branch cut policy?

  • casting from complex to real

    • discard imaginary component (as in C, NumPy)
    • other choice would to require user to explicitly invoke real before casting
  • complex numbers with components which are infinity and/or NaN

    • in Python, a complex number can be both infinite and NaN (according to cmath.isinf, cmath.isnan)

      In [1]: x = complex(float('inf'), float('nan'))
      
      In [2]: x
      Out[2]: (inf+nanj)
      
      In [3]: cmath.isinf(x)
      Out[3]: True
      
      In [4]: cmath.isnan(x)
      Out[4]: True
      
    • in NumPy, if any component is NaN, the complex number is NaN

      In [1]: np.isnan([np.inf+np.nan*1j])
      Out[1]: array([ True], dtype=bool)
      
      In [2]: np.isinf([np.inf+np.nan*1j])
      Out[2]: array([False], dtype=bool)
      
    • in C99 (see "complex floating types" section), one infinity model (e.g., inf + nan*j === inf)

      • "In order to support the one-infinity model of complex number arithmetic, C regards any complex value with at least one infinite part as an infinity even if its other part is a NaN..."
      • See C99 rationale 7.3.9.4 and SO

Creation Functions

arange

asarray

empty

  • No changes. Output array data type has no restrictions.

empty_like

  • No changes. Output array data type has no restrictions.

eye

from_dlpack

  • No changes.

full

full_like

linspace

  • start

    • add support for complex
  • stop

    • add support for complex
  • when either start or stop is complex, the result must be complex.

  • status: todo

    • PR: TODO

meshgrid

ones

ones_like

tril

  • No changes necessary. No restrictions on input array data types.

triu

  • No changes necessary. No restrictions on input array data types.

zeros

  • No changes necessary. No restrictions on output array data types.

zeros_like

  • No changes necessary. No restrictions on output array data types.

Data Type Functions

astype

  • Add note stating that, when casting from a complex data type to a real data type stating, the imaginary component is discarded and casting to integral data types is unspecified and thus implementation-dependent.

    • Discarding the imaginary component follows C and is done in NumPy.
    • Could also disallow/omit/not explicitly support casting from complex to real and require users to explicitly call real if they want to discard the imaginary component.
  • status: wip

broadcast_arrays

  • No changes necessary.

broadcast_to

  • No changes necessary.

can_cast

  • No changes necessary. Follows promotion rules for complex numbers.

finfo

  • No changes necessary.

iinfo

  • No changes necessary.

result_type

  • No changes necessary. Follows promotion rules for complex numbers.

Data Types

  • Add complex64 and complex128 following precedent where the numeric suffix specifies the number of bits.

    • The real and imaginary components should be IEEE 754 floating-point numbers.
  • status: implemented

Default Data Types

  • The default complex number data type is dependent on the default floating-point data type. If the latter is float32, the default complex data type must be complex64.

    • PyTorch already does this.
  • status: implemented

Data Type Categories

  • Add complex data types to list of numeric data types.

  • Add "Real Data Types" category which should be equal to the current list of numeric data types.

  • Add "Complex Data Types" category which only includes the complex number data types.

  • status: implemented


Element-wise Functions

abs

  • Add support by returning the complex magnitude.
  • Add complex number special cases: same as hypot(creal(z), cimag(z)).

acos

acosh

add

  • Add complex number support. Complex number addition is well-defined.
  • Some care should be given to infinities and nans.

asin

asinh

atan

atan2

  • No complex number support.
  • No changes necessary as already limited to real number data types.

atanh

bitwise_and

  • No changes. Complex dtypes are not allowed.

bitwise_left_shift

  • No changes. Complex dtypes are not allowed.

bitwise_invert

  • No changes. Complex dtypes are not allowed.

bitwise_or

  • No changes. Complex dtypes are not allowed.

bitwise_right_shift

  • No changes. Complex dtypes are not allowed.

bitwise_xor

  • No changes. Complex dtypes are not allowed.

ceil

  • ~~Add complex number support. Independently round components (e.g., as in MATLAB).~~

    • NumPy does not support in np.ceil.
    • No complex number support due to implementation ambiguity.

cos

cosh

divide

  • Add complex number support.
  • Define special cases.

equal

  • Add complex number support.
  • Care should be taken with infinities and nans.

exp

expm1

floor

  • ~~Add complex number support. Independently round components (e.g., as in MATLAB).~~

    • NumPy does not support in np.floor.
    • No complex number support due to implementation ambiguity.

floor_divide

  • No complex number support.
  • No changes necessary as already limited to real number data types.

greater

  • Add support for complex numbers by ordering lexicographic order. Similar to sorting functions.

greater_equal

  • Same as greater.

isfinite

  • Add support for complex numbers.
  • Care should be taken for infinities and nans.

isinf

  • Add support for complex numbers.
  • Care should be taken for infinities and nans.
    • how to classify a complex number as infinite?

isnan

  • Add support for complex numbers.
  • Care should be taken for infinities and nans.
    • how to classify a complex number as NaN?

less

  • Same as greater.

less_equal

  • Same as less.

log

log1p

log2

log10

logaddexp

  • No complex number support. This function is mainly useful for stats.
  • No changes necessary as already limited to real-valued data types.

logical_and

  • No changes necessary. Input arrays expected to have boolean data type.

logical_not

  • No changes necessary. Input arrays expected to have boolean data type.

logical_or

  • No changes necessary. Input arrays expected to have boolean data type.

logical_xor

  • No changes necessary. Input arrays expected to have boolean data type.

multiply

  • Add complex number support. Complex number multiplication is well-defined.
  • Define special cases. Care should be taken with regard to infinities and nans.

negative

not_equal

  • Add complex number support.
  • Care should be taken with infinities and nans.

positive

pow

remainder

  • No complex number support.
  • No changes necessary as already limited to real-valued data types.

round

sign

sin

sinh

square

sqrt

subtract

  • Add complex number support. Complex number addition is well-defined.

tan

tanh

trunc

  • ~~Add complex number support. Independently truncate components (e.g., as in MATLAB).~~

    • NumPy does not support in np.trunc.
    • No complex number support due to implementation ambiguity.

Linear Algebra Functions

matmul

  • No changes necessary. x1 and x2 are allowed to be any numeric data type. Complex number multiplication and addition is well-defined.

matrix_transpose

  • No changes necessary.

tensordot

  • No changes necessary. x1 and x2 are allowed to be any numeric data type. Complex number multiplication and addition is well-defined.

vecdot


Manipulation Functions

concat

  • No changes necessary.

expand_dims

  • No changes necessary.

flip

  • No changes necessary.

permute_dims

  • No changes necessary.

reshape

  • No changes necessary.

roll

  • No changes necessary.

squeeze

  • No changes necessary.

stack

  • No changes necessary.

Searching Functions

argmax

  • Require that complex numbers be ordered in lexicographic order?

    • While not compatible with multiplication, lexicographic order is a common total order. E.g., NumPy uses lexicographic order.
    • Another alternative is by magnitude and then phase angle. E.g., MATLAB supports this order.
    • The specification could default to lexicographic, and revisit the ordering relation in the future with possible support for optionally specifying an alternative comparison method.
    • Some care needs to be given to sorting when one or more components is NaN.

argmin

  • Same as argmax.

nonzero

where

  • No changes necessary.

Set Functions

unique_all

  • Specify complex number value equality.

  • In C, "In order to support the one-infinity model of complex number arithmetic, C regards any complex value with at least one infinite part as an infinity even if its other part is a NaN..."

  • NumPy (and Python) does not follow C. If one component is NaN, even if the other component is infinite, np.equal is false.

    • For Python's cmath, a complex number can be both infinite (cmath.isinf) and NaN (cmath.isnan).

unique_counts

  • Same as unique_all.

unique_inverse

  • Same as unique_all.

unique_values

  • Same as unique_all.

Sorting Functions

argsort

  • Require that complex numbers be ordered in lexicographic order?

    • While not compatible with multiplication, lexicographic order is a common total order. E.g., NumPy uses lexicographic order.
    • Another alternative is by magnitude and then phase angle. E.g., MATLAB supports this order.
    • The specification could default to lexicographic, and revisit the ordering relation in the future with possible support for optionally specifying an alternative comparison method.
    • Some care needs to be given to sorting when one or more components is NaN.

sort

  • Same as argsort

Statistical Functions

max

  • Similar to elsewhere, require that complex numbers be ordered in lexicographic order?

mean

  • No complex number support.
  • Change to only real number data types.

min

  • Same as max.

prod

  • Add support for complex numbers. Complex number multiplication is well-defined.
  • The dtype option needs to be updated to accommodate complex numbers.

std

  • No complex number support.
  • Change to only real number data types.

sum

  • Add support for complex numbers. Complex number addition is well-defined.
  • The dtype option needs to be updated to accommodate complex numbers.

var

  • No complex number support.
  • Change to only real number data types.

Type Promotion

Rules

  • Add complex number type promotion table.
  • Related complex number type promotion to floating-point number type promotion (e.g., f8 and c16)

Mixing arrays with Python Scalars

  • Add complex Python scalars.
  • Update note concerning mixed "kind" operations.

Utility Functions

all

any


Linear Algebra Extension

cholesky

cross

  • Add complex number support.

determinant

  • Add complex number support.

diagonal

  • No changes necessary.

eigh

  • Add complex number support.
  • Require that each matrix must be Hermitian and the returned Q unitary.
  • May want to make the dtype of the eigevalues unconditionally real.

eigvalsh

  • Add complex number support.
  • Require that each matrix must be Hermitian.

inv

  • Add complex number support.

matmul

  • Support added in main namespace.

matrix_norm

  • Add complex number support.

matrix_power

  • Add complex number support.

matrix_rank

  • Add complex number support.

matrix_transpose

  • No changes necessary.

outer

  • No changes necessary.

pinv

  • Add complex number support.

qr

  • Add complex number support.

slogdet

  • Add complex number support.

solve

  • Add complex number support.

svd

  • Add coplex number support.
  • Require that each U, Vh must be uniary.
  • May want to require that the returned dtype of S is unconditionally real.

svdvals

  • Add complex number support.

trace

  • No changes necessary.

vecdot

  • Changes made in main namespace.

vector_norm

  • Add complex number support.

kgryte avatar Jan 20 '22 11:01 kgryte

Would we want a data type function for complex dtypes i.e. cinfo()? Currently NumPy and PyTorch don't have a dedicate method, but complex dtypes can fallback on finfo().

Interestingly there was a little discussion of a torch.cinfo() at https://github.com/pytorch/pytorch/issues/35954#issuecomment-620223025

honno avatar Jan 21 '22 10:01 honno

  • branch cut policy?

I don't know why I missed it... CuPy uses Thrust (and it's likely being used by many other projects that depends on CUDA; IIUC the implementation was ported to libcudacxx) so it's worth mentioning.

Thrust's complex math funcs were based on FreeBSD's. In each FreeBSD man page the branch cut is clearly documented so it's a good starting point. However, as Thrust didn't document it clearly (I asked internally and was told to inspect the source directly, which is challenging given my lack of bandwidth) and it's unclear to me if the Thrust impl was a one-to-one translation of FreeBSD's since I don't have access to the FreeBSD source code, I'd proceed with caution.

leofang avatar Apr 18 '22 19:04 leofang

Tracking changes to existing APIs within the array API specification:

New APIs

kgryte avatar May 23 '22 00:05 kgryte

At this point, PRs have been submitted for all existing APIs for which we currently expect complex number support.

kgryte avatar Dec 13 '22 10:12 kgryte

All PRs have been merged - adding complex number support to the API is complete. Thanks a lot to everyone who contributed - and in particular @kgryte for all the spec writing.

rgommers avatar Dec 14 '22 21:12 rgommers