Branch data Line data Source code
1 : : /* $OpenBSD: s_cpowf.c,v 1.2 2010/07/18 18:42:26 guenther Exp $ */
2 : : /*
3 : : * Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net>
4 : : *
5 : : * Permission to use, copy, modify, and distribute this software for any
6 : : * purpose with or without fee is hereby granted, provided that the above
7 : : * copyright notice and this permission notice appear in all copies.
8 : : *
9 : : * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
10 : : * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
11 : : * MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
12 : : * ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
13 : : * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
14 : : * ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
15 : : * OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
16 : : */
17 : :
18 : : /* cpowf
19 : : *
20 : : * Complex power function
21 : : *
22 : : *
23 : : *
24 : : * SYNOPSIS:
25 : : *
26 : : * float complex cpowf();
27 : : * float complex a, z, w;
28 : : *
29 : : * w = cpowf (a, z);
30 : : *
31 : : *
32 : : *
33 : : * DESCRIPTION:
34 : : *
35 : : * Raises complex A to the complex Zth power.
36 : : * Definition is per AMS55 # 4.2.8,
37 : : * analytically equivalent to cpow(a,z) = cexp(z clog(a)).
38 : : *
39 : : * ACCURACY:
40 : : *
41 : : * Relative error:
42 : : * arithmetic domain # trials peak rms
43 : : * IEEE -10,+10 30000 9.4e-15 1.5e-15
44 : : *
45 : : */
46 : :
47 : : #include <openlibm_complex.h>
48 : : #include <openlibm_math.h>
49 : :
50 : : #include "math_private.h"
51 : :
52 : : OLM_DLLEXPORT float complex
53 : 0 : cpowf(float complex a, float complex z)
54 : : {
55 : : float complex w;
56 : : float x, y, r, theta, absa, arga;
57 : :
58 : 0 : x = crealf(z);
59 : 0 : y = cimagf(z);
60 : 0 : absa = cabsf (a);
61 [ # # ]: 0 : if (absa == 0.0f) {
62 : 0 : return (0.0f + 0.0f * I);
63 : : }
64 : 0 : arga = cargf (a);
65 : 0 : r = powf (absa, x);
66 : 0 : theta = x * arga;
67 [ # # ]: 0 : if (y != 0.0f) {
68 : 0 : r = r * expf (-y * arga);
69 : 0 : theta = theta + y * logf (absa);
70 : : }
71 : 0 : w = r * cosf (theta) + (r * sinf (theta)) * I;
72 : 0 : return (w);
73 : : }
|