summaryrefslogtreecommitdiff
path: root/algorithms/polynomial.c
blob: ecd90279923aa0f6c98099a54db79dbf11635bf9 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define EPSILON 1E-9


double
eval_poly(double x, double coeffs[], int order)
{
  double res;
  int i;

  res = coeffs[0];

  for ( i = 1; i <= order; i++ ) {
    res = res * x + coeffs[i];
  }

  return res;
}


double
find_zero(double left, double right, double coeffs[], int order)
{
  double c, lp, rp, cp, eps=EPSILON;

  if(eval_poly(left, coeffs, order) * eval_poly(right, coeffs, order) >= 0) {
    exit(1);
  }

  while(1) {
    c = (left+right)/2;

    lp = eval_poly(left, coeffs, order);
    rp = eval_poly(right, coeffs, order);
    cp = eval_poly(c, coeffs, order);

    if(fabs(cp) < eps) {
      break;
    }

    if(cp * lp < 0) {
      right = c;
    }
    else if(cp * rp < 0) {
      left = c;
    }
    else {
      exit(1);
    }
  }

  return c;
}

int
main(void)
{
  double coeffs[] = { 7, 5, 3, 2, 1, 1 };

  printf("eval_poly: %f\n", eval_poly(3, coeffs, 5));
  printf("%f\n", find_zero(-100, 100, coeffs, 5));

  return 0;
}