1/* complex/gsl_complex.h
2 *
3 * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Gerard Jungman, Brian Gough
4 * Copyright (C) 2020, 2021 Patrick Alken
5 *
6 * This program is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 3 of the License, or (at
9 * your option) any later version.
10 *
11 * This program is distributed in the hope that it will be useful, but
12 * WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 * General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program; if not, write to the Free Software
18 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
19 */
20
21#ifndef __GSL_COMPLEX_H__
22#define __GSL_COMPLEX_H__
23
24#undef __BEGIN_DECLS
25#undef __END_DECLS
26#ifdef __cplusplus
27# define __BEGIN_DECLS extern "C" {
28# define __END_DECLS }
29#else
30# define __BEGIN_DECLS /* empty */
31# define __END_DECLS /* empty */
32#endif
33
34__BEGIN_DECLS
35
36
37/* two consecutive built-in types as a complex number */
38typedef double * gsl_complex_packed ;
39typedef float * gsl_complex_packed_float ;
40typedef long double * gsl_complex_packed_long_double ;
41
42typedef const double * gsl_const_complex_packed ;
43typedef const float * gsl_const_complex_packed_float ;
44typedef const long double * gsl_const_complex_packed_long_double ;
45
46
47/* 2N consecutive built-in types as N complex numbers */
48typedef double * gsl_complex_packed_array ;
49typedef float * gsl_complex_packed_array_float ;
50typedef long double * gsl_complex_packed_array_long_double ;
51
52typedef const double * gsl_const_complex_packed_array ;
53typedef const float * gsl_const_complex_packed_array_float ;
54typedef const long double * gsl_const_complex_packed_array_long_double ;
55
56
57/* Yes... this seems weird. Trust us. The point is just that
58 sometimes you want to make it obvious that something is
59 an output value. The fact that it lacks a 'const' may not
60 be enough of a clue for people in some contexts.
61 */
62typedef double * gsl_complex_packed_ptr ;
63typedef float * gsl_complex_packed_float_ptr ;
64typedef long double * gsl_complex_packed_long_double_ptr ;
65
66typedef const double * gsl_const_complex_packed_ptr ;
67typedef const float * gsl_const_complex_packed_float_ptr ;
68typedef const long double * gsl_const_complex_packed_long_double_ptr ;
69
70/*
71 * If <complex.h> is included, use the C99 complex type. Otherwise
72 * define a type bit-compatible with C99 complex. The GSL_REAL and GSL_IMAG
73 * macros require C11 functionality also (_Generic)
74 */
75
76/* older gcc compilers claim to be C11 compliant but do not support _Generic */
77#if defined(__GNUC__) && (__GNUC__ < 7)
78# define GSL_COMPLEX_LEGACY 1
79#endif
80
81#if !defined(GSL_COMPLEX_LEGACY) && \
82 defined(_Complex_I) && \
83 defined(complex) && \
84 defined(I) && \
85 defined(__STDC__) && (__STDC__ == 1) && \
86 defined(__STDC_VERSION__) && (__STDC_VERSION__ >= 201112L) /* C11 */
87
88# define GSL_COMPLEX_DEFINE(R, C) typedef R _Complex C ;
89
90# define GSL_COMPLEX_P(zp) (&(zp))
91# define GSL_COMPLEX_EQ(z1,z2) ((z1) == (z2))
92# define GSL_SET_COMPLEX(zp,x,y) (*(zp) = (x) + I * (y))
93
94# define GSL_REAL(z) (_Generic((z), \
95 complex float : ((float *) &(z)), \
96 complex double : ((double *) &(z)), \
97 complex long double : ((long double *) &(z)))[0])
98
99# define GSL_IMAG(z) (_Generic((z), \
100 complex float : ((float *) &(z)), \
101 complex double : ((double *) &(z)), \
102 complex long double : ((long double *) &(z)))[1])
103
104# define GSL_COMPLEX_P_REAL(zp) GSL_REAL(*(zp))
105# define GSL_COMPLEX_P_IMAG(zp) GSL_IMAG(*(zp))
106# define GSL_SET_REAL(zp,x) do { GSL_REAL(*(zp)) = (x); } while(0)
107# define GSL_SET_IMAG(zp,x) do { GSL_IMAG(*(zp)) = (x); } while(0)
108
109#else /* legacy complex definitions */
110
111/*
112 * According to the C17 standard, 6.2.5 paragraph 13:
113 *
114 * "Each complex type has the same representation and alignment requirements
115 * as an array type containing exactly two elements of the corresponding real
116 * type; the first element is equal to the real part, and the second element to
117 * the imaginary part, of the complex number."
118 */
119
120/*# define GSL_COMPLEX_DEFINE(R, C) typedef R C[2]*/
121# define GSL_COMPLEX_DEFINE(R, C) typedef struct { R dat[2]; } C ;
122
123# define GSL_REAL(z) ((z).dat[0])
124# define GSL_IMAG(z) ((z).dat[1])
125# define GSL_COMPLEX_P(zp) ((zp)->dat)
126# define GSL_COMPLEX_P_REAL(zp) ((zp)->dat[0])
127# define GSL_COMPLEX_P_IMAG(zp) ((zp)->dat[1])
128# define GSL_COMPLEX_EQ(z1,z2) (((z1).dat[0] == (z2).dat[0]) && ((z1).dat[1] == (z2).dat[1]))
129
130# define GSL_SET_COMPLEX(zp,x,y) do {(zp)->dat[0]=(x); (zp)->dat[1]=(y);} while(0)
131# define GSL_SET_REAL(zp,x) do {(zp)->dat[0]=(x);} while(0)
132# define GSL_SET_IMAG(zp,y) do {(zp)->dat[1]=(y);} while(0)
133
134#endif
135
136GSL_COMPLEX_DEFINE(double, gsl_complex)
137GSL_COMPLEX_DEFINE(long double, gsl_complex_long_double)
138GSL_COMPLEX_DEFINE(float, gsl_complex_float)
139
140#define GSL_SET_COMPLEX_PACKED(zp,n,x,y) do {*((zp)+2*(n))=(x); *((zp)+(2*(n)+1))=(y);} while(0)
141
142__END_DECLS
143
144#endif /* __GSL_COMPLEX_H__ */
145