xref: /openbsd-src/lib/libm/src/s_cpow.c (revision 2f2c00629eff6a304ebffb255fc56f4fa7a1833b)
1*2f2c0062Sguenther /*	$OpenBSD: s_cpow.c,v 1.7 2016/09/12 19:47:02 guenther Exp $	*/
27b36286aSmartynas /*
37b36286aSmartynas  * Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net>
47b36286aSmartynas  *
57b36286aSmartynas  * Permission to use, copy, modify, and distribute this software for any
67b36286aSmartynas  * purpose with or without fee is hereby granted, provided that the above
77b36286aSmartynas  * copyright notice and this permission notice appear in all copies.
87b36286aSmartynas  *
97b36286aSmartynas  * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
107b36286aSmartynas  * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
117b36286aSmartynas  * MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
127b36286aSmartynas  * ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
137b36286aSmartynas  * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
147b36286aSmartynas  * ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
157b36286aSmartynas  * OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
167b36286aSmartynas  */
177b36286aSmartynas 
187b36286aSmartynas /*							cpow
197b36286aSmartynas  *
207b36286aSmartynas  *	Complex power function
217b36286aSmartynas  *
227b36286aSmartynas  *
237b36286aSmartynas  *
247b36286aSmartynas  * SYNOPSIS:
257b36286aSmartynas  *
267b36286aSmartynas  * double complex cpow();
277b36286aSmartynas  * double complex a, z, w;
287b36286aSmartynas  *
297b36286aSmartynas  * w = cpow (a, z);
307b36286aSmartynas  *
317b36286aSmartynas  *
327b36286aSmartynas  *
337b36286aSmartynas  * DESCRIPTION:
347b36286aSmartynas  *
357b36286aSmartynas  * Raises complex A to the complex Zth power.
367b36286aSmartynas  * Definition is per AMS55 # 4.2.8,
377b36286aSmartynas  * analytically equivalent to cpow(a,z) = cexp(z clog(a)).
387b36286aSmartynas  *
397b36286aSmartynas  * ACCURACY:
407b36286aSmartynas  *
417b36286aSmartynas  *                      Relative error:
427b36286aSmartynas  * arithmetic   domain     # trials      peak         rms
437b36286aSmartynas  *    IEEE      -10,+10     30000       9.4e-15     1.5e-15
447b36286aSmartynas  *
457b36286aSmartynas  */
467b36286aSmartynas 
477b36286aSmartynas #include <complex.h>
48de3697aaSmartynas #include <float.h>
497b36286aSmartynas #include <math.h>
507b36286aSmartynas 
517b36286aSmartynas double complex
cpow(double complex a,double complex z)527b36286aSmartynas cpow(double complex a, double complex z)
537b36286aSmartynas {
547b36286aSmartynas 	double complex w;
557b36286aSmartynas 	double x, y, r, theta, absa, arga;
567b36286aSmartynas 
577b36286aSmartynas 	x = creal (z);
587b36286aSmartynas 	y = cimag (z);
597b36286aSmartynas 	absa = cabs (a);
607b36286aSmartynas 	if (absa == 0.0) {
617b36286aSmartynas 		return (0.0 + 0.0 * I);
627b36286aSmartynas 	}
637b36286aSmartynas 	arga = carg (a);
647b36286aSmartynas 	r = pow (absa, x);
657b36286aSmartynas 	theta = x * arga;
667b36286aSmartynas 	if (y != 0.0) {
677b36286aSmartynas 		r = r * exp (-y * arga);
687b36286aSmartynas 		theta = theta + y * log (absa);
697b36286aSmartynas 	}
707b36286aSmartynas 	w = r * cos (theta) + (r * sin (theta)) * I;
717b36286aSmartynas 	return (w);
727b36286aSmartynas }
73*2f2c0062Sguenther DEF_STD(cpow);
74*2f2c0062Sguenther LDBL_MAYBE_UNUSED_CLONE(cpow);
75