xref: /freebsd/lib/msun/src/s_cpowf.c (revision 5a4c3b831b2f27edc11dd6ee7d0f8036dad76a87)
16813d08fSMatt Macy /*-
26813d08fSMatt Macy  * Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net>
36813d08fSMatt Macy  *
46813d08fSMatt Macy  * Permission to use, copy, modify, and distribute this software for any
56813d08fSMatt Macy  * purpose with or without fee is hereby granted, provided that the above
66813d08fSMatt Macy  * copyright notice and this permission notice appear in all copies.
76813d08fSMatt Macy  *
86813d08fSMatt Macy  * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
96813d08fSMatt Macy  * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
106813d08fSMatt Macy  * MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
116813d08fSMatt Macy  * ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
126813d08fSMatt Macy  * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
136813d08fSMatt Macy  * ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
146813d08fSMatt Macy  * OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
156813d08fSMatt Macy  */
166813d08fSMatt Macy 
176813d08fSMatt Macy /*							cpowf
186813d08fSMatt Macy  *
196813d08fSMatt Macy  *	Complex power function
206813d08fSMatt Macy  *
216813d08fSMatt Macy  *
226813d08fSMatt Macy  *
236813d08fSMatt Macy  * SYNOPSIS:
246813d08fSMatt Macy  *
256813d08fSMatt Macy  * float complex cpowf();
266813d08fSMatt Macy  * float complex a, z, w;
276813d08fSMatt Macy  *
286813d08fSMatt Macy  * w = cpowf (a, z);
296813d08fSMatt Macy  *
306813d08fSMatt Macy  *
316813d08fSMatt Macy  *
326813d08fSMatt Macy  * DESCRIPTION:
336813d08fSMatt Macy  *
346813d08fSMatt Macy  * Raises complex A to the complex Zth power.
356813d08fSMatt Macy  * Definition is per AMS55 # 4.2.8,
366813d08fSMatt Macy  * analytically equivalent to cpow(a,z) = cexp(z clog(a)).
376813d08fSMatt Macy  *
386813d08fSMatt Macy  * ACCURACY:
396813d08fSMatt Macy  *
406813d08fSMatt Macy  *                      Relative error:
416813d08fSMatt Macy  * arithmetic   domain     # trials      peak         rms
426813d08fSMatt Macy  *    IEEE      -10,+10     30000       9.4e-15     1.5e-15
436813d08fSMatt Macy  *
446813d08fSMatt Macy  */
456813d08fSMatt Macy 
466813d08fSMatt Macy #include <sys/cdefs.h>
476813d08fSMatt Macy __FBSDID("$FreeBSD$");
486813d08fSMatt Macy 
496813d08fSMatt Macy #include <complex.h>
506813d08fSMatt Macy #include <math.h>
51*5a4c3b83SDimitry Andric #include "math_private.h"
526813d08fSMatt Macy 
536813d08fSMatt Macy float complex
546813d08fSMatt Macy cpowf(float complex a, float complex z)
556813d08fSMatt Macy {
566813d08fSMatt Macy 	float complex w;
576813d08fSMatt Macy 	float x, y, r, theta, absa, arga;
586813d08fSMatt Macy 
596813d08fSMatt Macy 	x = crealf(z);
606813d08fSMatt Macy 	y = cimagf(z);
616813d08fSMatt Macy 	absa = cabsf (a);
626813d08fSMatt Macy 	if (absa == 0.0f) {
63*5a4c3b83SDimitry Andric 		return (CMPLXF(0.0f, 0.0f));
646813d08fSMatt Macy 	}
656813d08fSMatt Macy 	arga = cargf (a);
666813d08fSMatt Macy 	r = powf (absa, x);
676813d08fSMatt Macy 	theta = x * arga;
686813d08fSMatt Macy 	if (y != 0.0f) {
696813d08fSMatt Macy 		r = r * expf (-y * arga);
706813d08fSMatt Macy 		theta = theta + y * logf (absa);
716813d08fSMatt Macy 	}
72*5a4c3b83SDimitry Andric 	w = CMPLXF(r * cosf (theta), r * sinf (theta));
736813d08fSMatt Macy 	return (w);
746813d08fSMatt Macy }
75