src/math/__polevll.c

   1 /* origin: OpenBSD /usr/src/lib/libm/src/polevll.c */
   2 /*
   3  * Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net>
   4  *
   5  * Permission to use, copy, modify, and distribute this software for any
   6  * purpose with or without fee is hereby granted, provided that the above
   7  * copyright notice and this permission notice appear in all copies.
   8  *
   9  * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
  10  * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
  11  * MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
  12  * ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
  13  * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
  14  * ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
  15  * OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
  16  */
  17 /*
  18  *      Evaluate polynomial
  19  *
  20  *
  21  * SYNOPSIS:
  22  *
  23  * int N;
  24  * long double x, y, coef[N+1], polevl[];
  25  *
  26  * y = polevll( x, coef, N );
  27  *
  28  *
  29  * DESCRIPTION:
  30  *
  31  * Evaluates polynomial of degree N:
  32  *
  33  *                     2          N
  34  * y  =  C  + C x + C x  +...+ C x
  35  *        0    1     2          N
  36  *
  37  * Coefficients are stored in reverse order:
  38  *
  39  * coef[0] = C  , ..., coef[N] = C  .
  40  *            N                   0
  41  *
  42  *  The function p1evll() assumes that coef[N] = 1.0 and is
  43  * omitted from the array.  Its calling arguments are
  44  * otherwise the same as polevll().
  45  *
  46  *
  47  * SPEED:
  48  *
  49  * In the interest of speed, there are no checks for out
  50  * of bounds arithmetic.  This routine is used by most of
  51  * the functions in the library.  Depending on available
  52  * equipment features, the user may wish to rewrite the
  53  * program in microcode or assembly language.
  54  *
  55  */
  56
  57 #include "libm.h"
  58
  59 /*
  60  * Polynomial evaluator:
  61  *  P[0] x^n  +  P[1] x^(n-1)  +  ...  +  P[n]
  62  */
  63 long double __polevll(long double x, const long double *P, int n)
  64 {
  65         long double y;
  66
  67         y = *P++;
  68         do {
  69                 y = y * x + *P++;
  70         } while (--n);
  71
  72         return y;
  73 }
  74
  75 /*
  76  * Polynomial evaluator:
  77  *  x^n  +  P[0] x^(n-1)  +  P[1] x^(n-2)  +  ...  +  P[n]
  78  */
  79 long double __p1evll(long double x, const long double *P, int n)
  80 {
  81         long double y;
  82
  83         n -= 1;
  84         y = x + *P++;
  85         do {
  86                 y = y * x + *P++;
  87         } while (--n);
  88
  89         return y;
  90 }