codekingpro/portable-devtools
114k
1/*
2 * Copyright 1993-2012 NVIDIA Corporation. All rights reserved.
3 *
4 * NOTICE TO LICENSEE:
5 *
6 * This source code and/or documentation ("Licensed Deliverables") are
7 * subject to NVIDIA intellectual property rights under U.S. and
8 * international Copyright laws.
9 *
10 * These Licensed Deliverables contained herein is PROPRIETARY and
11 * CONFIDENTIAL to NVIDIA and is being provided under the terms and
12 * conditions of a form of NVIDIA software license agreement by and
13 * between NVIDIA and Licensee ("License Agreement") or electronically
14 * accepted by Licensee. Notwithstanding any terms or conditions to
15 * the contrary in the License Agreement, reproduction or disclosure
16 * of the Licensed Deliverables to any third party without the express
17 * written consent of NVIDIA is prohibited.
18 *
19 * NOTWITHSTANDING ANY TERMS OR CONDITIONS TO THE CONTRARY IN THE
20 * LICENSE AGREEMENT, NVIDIA MAKES NO REPRESENTATION ABOUT THE
21 * SUITABILITY OF THESE LICENSED DELIVERABLES FOR ANY PURPOSE. IT IS
22 * PROVIDED "AS IS" WITHOUT EXPRESS OR IMPLIED WARRANTY OF ANY KIND.
23 * NVIDIA DISCLAIMS ALL WARRANTIES WITH REGARD TO THESE LICENSED
24 * DELIVERABLES, INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY,
25 * NONINFRINGEMENT, AND FITNESS FOR A PARTICULAR PURPOSE.
26 * NOTWITHSTANDING ANY TERMS OR CONDITIONS TO THE CONTRARY IN THE
27 * LICENSE AGREEMENT, IN NO EVENT SHALL NVIDIA BE LIABLE FOR ANY
28 * SPECIAL, INDIRECT, INCIDENTAL, OR CONSEQUENTIAL DAMAGES, OR ANY
29 * DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS,
30 * WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS
31 * ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE
32 * OF THESE LICENSED DELIVERABLES.
33 *
34 * U.S. Government End Users. These Licensed Deliverables are a
35 * "commercial item" as that term is defined at 48 C.F.R. 2.101 (OCT
36 * 1995), consisting of "commercial computer software" and "commercial
37 * computer software documentation" as such terms are used in 48
38 * C.F.R. 12.212 (SEPT 1995) and is provided to the U.S. Government
39 * only as a commercial end item. Consistent with 48 C.F.R.12.212 and
40 * 48 C.F.R. 227.7202-1 through 227.7202-4 (JUNE 1995), all
41 * U.S. Government End Users acquire the Licensed Deliverables with
42 * only those rights set forth herein.
43 *
44 * Any use of the Licensed Deliverables in individual and commercial
45 * software must include, in the user documentation and internal
46 * comments to the code, the above Disclaimer and U.S. Government End
47 * Users Notice.
48 */
49
50#if !defined(CU_COMPLEX_H_)
51#define CU_COMPLEX_H_
52
53#if !defined(__CUDACC_RTC__)
54#if defined(__GNUC__)
55#if defined(__clang__) || (!defined(__PGIC__) && (__GNUC__ > 4 || (__GNUC__ == 4 && __GNUC_MINOR__ >= 2)))
56#pragma GCC diagnostic ignored "-Wunused-function"
57#endif
58#endif
59#endif
60
61/* When trying to include C header file in C++ Code extern "C" is required
62 * But the Standard QNX headers already have ifdef extern in them when compiling C++ Code
63 * extern "C" cannot be nested
64 * Hence keep the header out of extern "C" block
65 */
66
67#if !defined(__CUDACC__)
68#include <math.h> /* import fabsf, sqrt */
69#endif /* !defined(__CUDACC__) */
70
71#if defined(__cplusplus)
72extern "C" {
73#endif /* __cplusplus */
74
75#include "vector_types.h"
76
77typedef float2 cuFloatComplex;
78
79__host__ __device__ static __inline__ float cuCrealf (cuFloatComplex x)
80{
81 return x.x;
82}
83
84__host__ __device__ static __inline__ float cuCimagf (cuFloatComplex x)
85{
86 return x.y;
87}
88
89__host__ __device__ static __inline__ cuFloatComplex make_cuFloatComplex
90 (float r, float i)
91{
92 cuFloatComplex res;
93 res.x = r;
94 res.y = i;
95 return res;
96}
97
98__host__ __device__ static __inline__ cuFloatComplex cuConjf (cuFloatComplex x)
99{
100 return make_cuFloatComplex (cuCrealf(x), -cuCimagf(x));
101}
102__host__ __device__ static __inline__ cuFloatComplex cuCaddf (cuFloatComplex x,
103 cuFloatComplex y)
104{
105 return make_cuFloatComplex (cuCrealf(x) + cuCrealf(y),
106 cuCimagf(x) + cuCimagf(y));
107}
108
109__host__ __device__ static __inline__ cuFloatComplex cuCsubf (cuFloatComplex x,
110 cuFloatComplex y)
111{
112 return make_cuFloatComplex (cuCrealf(x) - cuCrealf(y),
113 cuCimagf(x) - cuCimagf(y));
114}
115
116/* This implementation could suffer from intermediate overflow even though
117 * the final result would be in range. However, various implementations do
118 * not guard against this (presumably to avoid losing performance), so we
119 * don't do it either to stay competitive.
120 */
121__host__ __device__ static __inline__ cuFloatComplex cuCmulf (cuFloatComplex x,
122 cuFloatComplex y)
123{
124 cuFloatComplex prod;
125 prod = make_cuFloatComplex ((cuCrealf(x) * cuCrealf(y)) -
126 (cuCimagf(x) * cuCimagf(y)),
127 (cuCrealf(x) * cuCimagf(y)) +
128 (cuCimagf(x) * cuCrealf(y)));
129 return prod;
130}
131
132/* This implementation guards against intermediate underflow and overflow
133 * by scaling. Such guarded implementations are usually the default for
134 * complex library implementations, with some also offering an unguarded,
135 * faster version.
136 */
137__host__ __device__ static __inline__ cuFloatComplex cuCdivf (cuFloatComplex x,
138 cuFloatComplex y)
139{
140 cuFloatComplex quot;
141 float s = fabsf(cuCrealf(y)) + fabsf(cuCimagf(y));
142 float oos = 1.0f / s;
143 float ars = cuCrealf(x) * oos;
144 float ais = cuCimagf(x) * oos;
145 float brs = cuCrealf(y) * oos;
146 float bis = cuCimagf(y) * oos;
147 s = (brs * brs) + (bis * bis);
148 oos = 1.0f / s;
149 quot = make_cuFloatComplex (((ars * brs) + (ais * bis)) * oos,
150 ((ais * brs) - (ars * bis)) * oos);
151 return quot;
152}
153
154/*
155 * We would like to call hypotf(), but it's not available on all platforms.
156 * This discrete implementation guards against intermediate underflow and
157 * overflow by scaling. Otherwise we would lose half the exponent range.
158 * There are various ways of doing guarded computation. For now chose the
159 * simplest and fastest solution, however this may suffer from inaccuracies
160 * if sqrt and division are not IEEE compliant.
161 */
162__host__ __device__ static __inline__ float cuCabsf (cuFloatComplex x)
163{
164 float a = cuCrealf(x);
165 float b = cuCimagf(x);
166 float v, w, t;
167 a = fabsf(a);
168 b = fabsf(b);
169 if (a > b) {
170 v = a;
171 w = b;
172 } else {
173 v = b;
174 w = a;
175 }
176 t = w / v;
177 t = 1.0f + t * t;
178 t = v * sqrtf(t);
179 if ((v == 0.0f) || (v > 3.402823466e38f) || (w > 3.402823466e38f)) {
180 t = v + w;
181 }
182 return t;
183}
184
185/* Double precision */
186typedef double2 cuDoubleComplex;
187
188__host__ __device__ static __inline__ double cuCreal (cuDoubleComplex x)
189{
190 return x.x;
191}
192
193__host__ __device__ static __inline__ double cuCimag (cuDoubleComplex x)
194{
195 return x.y;
196}
197
198__host__ __device__ static __inline__ cuDoubleComplex make_cuDoubleComplex
199 (double r, double i)
200{
201 cuDoubleComplex res;
202 res.x = r;
203 res.y = i;
204 return res;
205}
206
207__host__ __device__ static __inline__ cuDoubleComplex cuConj(cuDoubleComplex x)
208{
209 return make_cuDoubleComplex (cuCreal(x), -cuCimag(x));
210}
211
212__host__ __device__ static __inline__ cuDoubleComplex cuCadd(cuDoubleComplex x,
213 cuDoubleComplex y)
214{
215 return make_cuDoubleComplex (cuCreal(x) + cuCreal(y),
216 cuCimag(x) + cuCimag(y));
217}
218
219__host__ __device__ static __inline__ cuDoubleComplex cuCsub(cuDoubleComplex x,
220 cuDoubleComplex y)
221{
222 return make_cuDoubleComplex (cuCreal(x) - cuCreal(y),
223 cuCimag(x) - cuCimag(y));
224}
225
226/* This implementation could suffer from intermediate overflow even though
227 * the final result would be in range. However, various implementations do
228 * not guard against this (presumably to avoid losing performance), so we
229 * don't do it either to stay competitive.
230 */
231__host__ __device__ static __inline__ cuDoubleComplex cuCmul(cuDoubleComplex x,
232 cuDoubleComplex y)
233{
234 cuDoubleComplex prod;
235 prod = make_cuDoubleComplex ((cuCreal(x) * cuCreal(y)) -
236 (cuCimag(x) * cuCimag(y)),
237 (cuCreal(x) * cuCimag(y)) +
238 (cuCimag(x) * cuCreal(y)));
239 return prod;
240}
241
242/* This implementation guards against intermediate underflow and overflow
243 * by scaling. Such guarded implementations are usually the default for
244 * complex library implementations, with some also offering an unguarded,
245 * faster version.
246 */
247__host__ __device__ static __inline__ cuDoubleComplex cuCdiv(cuDoubleComplex x,
248 cuDoubleComplex y)
249{
250 cuDoubleComplex quot;
251 double s = (fabs(cuCreal(y))) + (fabs(cuCimag(y)));
252 double oos = 1.0 / s;
253 double ars = cuCreal(x) * oos;
254 double ais = cuCimag(x) * oos;
255 double brs = cuCreal(y) * oos;
256 double bis = cuCimag(y) * oos;
257 s = (brs * brs) + (bis * bis);
258 oos = 1.0 / s;
259 quot = make_cuDoubleComplex (((ars * brs) + (ais * bis)) * oos,
260 ((ais * brs) - (ars * bis)) * oos);
261 return quot;
262}
263
264/* This implementation guards against intermediate underflow and overflow
265 * by scaling. Otherwise we would lose half the exponent range. There are
266 * various ways of doing guarded computation. For now chose the simplest
267 * and fastest solution, however this may suffer from inaccuracies if sqrt
268 * and division are not IEEE compliant.
269 */
270__host__ __device__ static __inline__ double cuCabs (cuDoubleComplex x)
271{
272 double a = cuCreal(x);
273 double b = cuCimag(x);
274 double v, w, t;
275 a = fabs(a);
276 b = fabs(b);
277 if (a > b) {
278 v = a;
279 w = b;
280 } else {
281 v = b;
282 w = a;
283 }
284 t = w / v;
285 t = 1.0 + t * t;
286 t = v * sqrt(t);
287 if ((v == 0.0) ||
288 (v > 1.79769313486231570e+308) || (w > 1.79769313486231570e+308)) {
289 t = v + w;
290 }
291 return t;
292}
293
294#if defined(__cplusplus)
295}
296#endif /* __cplusplus */
297
298/* aliases */
299typedef cuFloatComplex cuComplex;
300__host__ __device__ static __inline__ cuComplex make_cuComplex (float x,
301 float y)
302{
303 return make_cuFloatComplex (x, y);
304}
305
306/* float-to-double promotion */
307__host__ __device__ static __inline__ cuDoubleComplex cuComplexFloatToDouble
308 (cuFloatComplex c)
309{
310 return make_cuDoubleComplex ((double)cuCrealf(c), (double)cuCimagf(c));
311}
312
313__host__ __device__ static __inline__ cuFloatComplex cuComplexDoubleToFloat
314(cuDoubleComplex c)
315{
316 return make_cuFloatComplex ((float)cuCreal(c), (float)cuCimag(c));
317}
318
319
320__host__ __device__ static __inline__ cuComplex cuCfmaf( cuComplex x, cuComplex y, cuComplex d)
321{
322 float real_res;
323 float imag_res;
324
325 real_res = (cuCrealf(x) * cuCrealf(y)) + cuCrealf(d);
326 imag_res = (cuCrealf(x) * cuCimagf(y)) + cuCimagf(d);
327
328 real_res = -(cuCimagf(x) * cuCimagf(y)) + real_res;
329 imag_res = (cuCimagf(x) * cuCrealf(y)) + imag_res;
330
331 return make_cuComplex(real_res, imag_res);
332}
333
334__host__ __device__ static __inline__ cuDoubleComplex cuCfma( cuDoubleComplex x, cuDoubleComplex y, cuDoubleComplex d)
335{
336 double real_res;
337 double imag_res;
338
339 real_res = (cuCreal(x) * cuCreal(y)) + cuCreal(d);
340 imag_res = (cuCreal(x) * cuCimag(y)) + cuCimag(d);
341
342 real_res = -(cuCimag(x) * cuCimag(y)) + real_res;
343 imag_res = (cuCimag(x) * cuCreal(y)) + imag_res;
344
345 return make_cuDoubleComplex(real_res, imag_res);
346}
347
348#endif /* !defined(CU_COMPLEX_H_) */
349 