So I have am trying to implement my own pow function. To make things somewhat simpler I am allowed to use type int for the exponent. So the code is working...mostly. I have a file called xmath.c in which I implement my pow function and another test.c to test the limits. The only case in which it's not working is for pow(-111, INT_MAX)
. How could I improve my code ? Also tell me if I should post the test.c code as well. xfabs
is another function that I had to implement on my own instead of using fabs
.
Forgot to include xfabs in my post:
#define xfabs(x) (((x) > 0) ? (x) : (-(x)))
#include <stdio.h>
#include <errno.h>
#include <math.h>
#include <float.h>
#include "xmath.h"
double xpow(double x, int y)
{
errno = 0;
double product = 1;
double prod = 1;
int i;
if (y == 0) {
return 1.0;
} else if (x == 0 && y < 0) {
errno = EDOM;
return 0.0;
}
if (xfabs(x) >= 1) {
if (y > 0) {
for (i = 0; i < y; i++) {
if (xfabs(product) > DBL_MAX / xfabs(x)) {
errno = ERANGE;
return HUGE_VAL;
}
product *= x;
}
return product;
} else if (y < 0) {
for (i = 0; i > y; i--) {
if (xfabs(product) < DBL_MIN / xfabs(x)) {
errno = ERANGE;
return -0.0;
}
product *= 1 / x;
}
return product;
}
} else if (xfabs(x) < 1) {
if (y > 0) {
for (i = 0; i < y; i++) {
if (xfabs(prod) < DBL_MIN / xfabs(x)) {
errno = ERANGE;
return +0.0;
}
prod *= x;
}
return prod;
} else if (y < 0) {
for (i = 0; i > y; i--) {
if (xfabs(prod) < DBL_MAX / xfabs(x)) {
errno = ERANGE;
return -HUGE_VAL;
}
prod *= 1 / x;
}
return prod;
}
}
}
Here is the test function. Also I am creating an object file for xmath.c which I use in test.c:
#include <stdio.h>
#include <float.h>
#include <math.h>
#include <errno.h>
#include <limits.h>
#include "xmath.h"
void xpow_tests(void)
{
printf("XPOW TESTS (SUCCESS == 1 / FAILURE == 0)\n");
printf("=================================================\n");
printf("Test (xpow(2, 5)): \t\t%i\t", xpow(2, 5) == pow(2, 5));
xpow(2, 5);
printf("errno: %i\n", errno);
printf("Test (xpow(-2, -4)): \t\t%i\t", xpow(-2, -4) == pow(-2, -4));
xpow(-2, -4);
printf("errno: %i\n", errno);
printf("Test (xpow(1, 15)): \t\t%i\t", xpow(1, 15) == pow(1, 15));
xpow(1, 15);
printf("errno: %i\n", errno);
printf("Test (xpow(-4, 414)): \t\t%i\t", xpow(-4, 414) == pow(-4, 414));
xpow(-4, 414);
printf("errno: %i\n", errno);
printf("Test (xpow(-5, 303)): \t\t%i\t", xpow(-5, 303) == pow(-5, 303));
xpow(-5, 303);
printf("errno: %i\n", errno);
printf("Test (xpow(0, -3)): \t\t%i\t", xpow(0, -3) == 0);
xpow(0, -3);
printf("errno: %i\n", errno);
printf("Test (xpow(1.0e-10, 100)): \t%i\t", xpow(1.0e-10, 100) == pow(1.0e-10, 100));
xpow(1.0e-10, 100);
printf("errno: %i\n", errno);
printf("Test (xpow(-1.0e-10, 101)): \t%i\t", xpow(-1.0e-10, 101) == pow(-1.0e-10, 101));
xpow(-1.0e-10, 101);
printf("errno: %i\n", errno);
printf("Test (xpow(-111, INT_MAX)): \t%i\t", xpow(-111, INT_MAX) == pow(-111, INT_MAX));
xpow(-111, INT_MAX);
printf("errno: %i\n", errno);
printf("Test (xpow(-111, INT_MAX-1)): \t%i\t", xpow(-111, INT_MAX-1) == pow(-111, INT_MAX-1));
xpow(-111, INT_MAX-1);
printf("errno: %i\n", errno);
printf("Test (xpow(-111, -INT_MAX)): \t%i\t", xpow(-111, -INT_MAX) == pow(-111, -INT_MAX));
xpow(-111, -INT_MAX);
printf("errno: %i\n", errno);
printf("Test (xpow(-111, -INT_MAX+1)): \t%i\t", xpow(-111, -INT_MAX+1) == pow(-111, -INT_MAX+1));
xpow(-111, -INT_MAX+1);
printf("errno: %i\n", errno);
printf("=================================================\n");
printf("errno (EDOM: %i / ERANGE: %i)\n\n", EDOM, ERANGE);
}
int main(void)
{
xpow_tests();
return 0;
}
Caveat: This isn't a solution to your failing test case per se.
But, I simplified your function a bit--this may help.
And I greatly simplified/enhanced your test case function [using some preprocessor trickery].
However, in addition to your failed test -111,INT_MAX
, when I ran it here, the 0,-3
test reports failure.
Here's the modified source:
#include <stdio.h>
#include <string.h>
#include <errno.h>
#include <math.h>
#include <float.h>
#include <limits.h>
//#include "xmath.h"
double
xfabs(double x)
{
if (x < 0)
x = -x;
return x;
}
double
xpow(double x, int y)
{
errno = 0;
int i;
double lim;
double prod = 1;
double abs_x;
if (y == 0)
return 1.0;
if ((x == 0) && (y < 0)) {
errno = EDOM;
return 0.0;
}
abs_x = xfabs(x);
if (abs_x >= 1) {
if (y > 0) {
lim = DBL_MAX / abs_x;
for (i = 0; i < y; i++) {
if (xfabs(prod) > lim) {
errno = ERANGE;
return HUGE_VAL;
}
prod *= x;
}
}
if (y < 0) {
lim = DBL_MIN / abs_x;
for (i = 0; i > y; i--) {
if (xfabs(prod) < lim) {
errno = ERANGE;
return -0.0;
}
prod *= 1 / x;
}
}
return prod;
}
if (y > 0) {
lim = DBL_MIN / abs_x;
for (i = 0; i < y; i++) {
if (xfabs(prod) < lim) {
errno = ERANGE;
return +0.0;
}
prod *= x;
}
}
if (y < 0) {
lim = DBL_MAX / abs_x;
for (i = 0; i > y; i--) {
if (xfabs(prod) < lim) {
errno = ERANGE;
return -HUGE_VAL;
}
prod *= 1 / x;
}
}
return prod;
}
int passcnt = 0;
int failcnt = 0;
void
dotest(const char *name,double x,int y)
{
int sverr_pow;
int sverr_xpow;
double expected;
double actual;
int passflg;
printf("\n");
printf("%s\n",name);
errno = 0;
expected = pow(x,y);
sverr_pow = errno;
printf("%.16g (pow)\n",expected);
if (sverr_pow)
printf("errno: %d (%s)\n",sverr_pow,strerror(sverr_pow));
errno = 0;
actual = xpow(x,y);
sverr_xpow = errno;
printf("%.16g (xpow)\n",actual);
if (sverr_xpow)
printf("errno: %d (%s)\n",sverr_xpow,strerror(sverr_xpow));
if ((actual == expected) && (sverr_xpow == sverr_pow))
passcnt += 1;
else
failcnt += 1;
printf("Result: Value:%s errno:%s\n",
(actual == expected) ? "PASS" : "FAIL",
(sverr_xpow == sverr_pow) ? "PASS" : "FAIL");
}
#define TEST(x,y) \
dotest(#x "," #y,x,y)
void xpow_tests(void)
{
printf("XPOW TESTS (SUCCESS == 1 / FAILURE == 0)\n");
printf("=================================================\n");
TEST(2, 5);
TEST(-2, -4);
TEST(1, 15);
TEST(-4, 414);
TEST(-5, 303);
TEST(0, -3);
TEST(1.0e-10, 100);
TEST(-1.0e-10, 101);
TEST(-111, INT_MAX);
TEST(-111, INT_MAX-1);
TEST(-111, -INT_MAX);
TEST(-111, -INT_MAX+1);
printf("=================================================\n");
printf("errno (EDOM: %i / ERANGE: %i)\n\n", EDOM, ERANGE);
printf("PASSED: %d, FAILED: %d\n",passcnt,failcnt);
}
int
main(void)
{
xpow_tests();
return 0;
}
Here's the output of the tests:
XPOW TESTS (SUCCESS == 1 / FAILURE == 0)
=================================================
2,5
32 (pow)
32 (xpow)
Result: Value:PASS errno:PASS
-2,-4
0.0625 (pow)
0.0625 (xpow)
Result: Value:PASS errno:PASS
1,15
1 (pow)
1 (xpow)
Result: Value:PASS errno:PASS
-4,414
1.789931494904685e+249 (pow)
1.789931494904685e+249 (xpow)
Result: Value:PASS errno:PASS
-5,303
-6.136366831622158e+211 (pow)
-6.136366831622158e+211 (xpow)
Result: Value:PASS errno:PASS
0,-3
inf (pow)
errno: 34 (Numerical result out of range)
0 (xpow)
errno: 33 (Numerical argument out of domain)
Result: Value:FAIL errno:FAIL
1.0e-10,100
0 (pow)
errno: 34 (Numerical result out of range)
0 (xpow)
errno: 34 (Numerical result out of range)
Result: Value:PASS errno:PASS
-1.0e-10,101
-0 (pow)
errno: 34 (Numerical result out of range)
0 (xpow)
errno: 34 (Numerical result out of range)
Result: Value:PASS errno:PASS
-111,INT_MAX
-inf (pow)
errno: 34 (Numerical result out of range)
inf (xpow)
errno: 34 (Numerical result out of range)
Result: Value:FAIL errno:PASS
-111,INT_MAX-1
inf (pow)
errno: 34 (Numerical result out of range)
inf (xpow)
errno: 34 (Numerical result out of range)
Result: Value:PASS errno:PASS
-111,-INT_MAX
-0 (pow)
errno: 34 (Numerical result out of range)
-0 (xpow)
errno: 34 (Numerical result out of range)
Result: Value:PASS errno:PASS
-111,-INT_MAX+1
0 (pow)
errno: 34 (Numerical result out of range)
-0 (xpow)
errno: 34 (Numerical result out of range)
Result: Value:PASS errno:PASS
=================================================
errno (EDOM: 33 / ERANGE: 34)
PASSED: 10, FAILED: 2