Home | History | Annotate | Download | only in src
      1 /*-
      2  * Copyright (c) 2008 Stephen L. Moshier <steve (at) moshier.net>
      3  *
      4  * Permission to use, copy, modify, and distribute this software for any
      5  * purpose with or without fee is hereby granted, provided that the above
      6  * copyright notice and this permission notice appear in all copies.
      7  *
      8  * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
      9  * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
     10  * MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
     11  * ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
     12  * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
     13  * ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
     14  * OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
     15  */
     16 
     17 /*							cpowl
     18  *
     19  *	Complex power function
     20  *
     21  *
     22  *
     23  * SYNOPSIS:
     24  *
     25  * long double complex cpowl();
     26  * long double complex a, z, w;
     27  *
     28  * w = cpowl (a, z);
     29  *
     30  *
     31  *
     32  * DESCRIPTION:
     33  *
     34  * Raises complex A to the complex Zth power.
     35  * Definition is per AMS55 # 4.2.8,
     36  * analytically equivalent to cpow(a,z) = cexp(z clog(a)).
     37  *
     38  * ACCURACY:
     39  *
     40  *                      Relative error:
     41  * arithmetic   domain     # trials      peak         rms
     42  *    IEEE      -10,+10     30000       9.4e-15     1.5e-15
     43  *
     44  */
     45 
     46 #include <sys/cdefs.h>
     47 __FBSDID("$FreeBSD: head/lib/msun/src/s_cpowl.c 336563 2018-07-20 18:27:30Z dim $");
     48 
     49 #include <complex.h>
     50 #include <math.h>
     51 #include "math_private.h"
     52 
     53 long double complex
     54 cpowl(long double complex a, long double complex z)
     55 {
     56 	long double complex w;
     57 	long double x, y, r, theta, absa, arga;
     58 
     59 	x = creall(z);
     60 	y = cimagl(z);
     61 	absa = cabsl(a);
     62 	if (absa == 0.0L) {
     63 		return (CMPLXL(0.0L, 0.0L));
     64 	}
     65 	arga = cargl(a);
     66 	r = powl(absa, x);
     67 	theta = x * arga;
     68 	if (y != 0.0L) {
     69 		r = r * expl(-y * arga);
     70 		theta = theta + y * logl(absa);
     71 	}
     72 	w = CMPLXL(r * cosl(theta), r * sinl(theta));
     73 	return (w);
     74 }
     75