Search code examples
cpowmath.h

Pow function in C is outputting zero


I'm writing an integration script, and it uses a non-uniform grid. I'm using the pow function in math.h to calculate by what each value of x should be multiplied to get to the next value of x (this value is xmult in my code). However, for some reason, my pow function outputs zero, even though my calculator says that it should not be. Here are the relevant snippets from my code:

#include "stdio.h"
#include "math.h"
...
int N = 10240;
double recN = 1/N;
double upperlim = 4;
double lowerlim = 0.0001;
double limratio = upperlim/lowerlim;
double xmult;
xmult = pow(limratio,recN);
printf("%e\n",xmult);

The output of that printf function is 1.000000e+00 even though my calculated value does not round to that. I am compiling using gcc on CentOS with the -lm flag, it does not compile without that flag set. If it is any more help, the full code is below (the final output of the program being 0.000000e+00).

#include "stdio.h"
#include "math.h"

int main() {
  int N = 10240;
  double recN = 1/N;
  double pi = M_PI;
  double upperlim = 4;
  double lowerlim = 0.0001;
  double limratio = upperlim/lowerlim;
  double h;
  double f;
  double fnext;
  double ftot = 0;
  double ftot2 = 0;
  double x=lowerlim;
  double xmult;
  double xnext;
  int i;
  double func(double x){
    return sin(x) * exp(-x);
    }
  fnext = func(lowerlim);
  xmult = pow(limratio,recN);
  printf("%e\n",xmult);
  for (i=0;i<N;i++){
    f = fnext;
    fnext = func(x*xmult);
    xnext = x * xmult;
    h = xnext - x;
    ftot = 0.5*h*(f+fnext);
    ftot2 = ftot2  + ftot;
    x = xnext;
  }
  printf("%e\n",ftot2);
 }

Solution

  • More than likely, your problem is a common error when dealing with numbers in C. Consider this:

    double recN = 1/N;

    Here, you define the variable as double, but you then perform integer division. If N>1 then this will result in zero.

    Instead, look for places you've done this type of calculation and convert it to:

    double recN = 1.0/N;