| /* origin: OpenBSD /usr/src/lib/libm/src/polevll.c */ |
| /* |
| * Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net> |
| * |
| * Permission to use, copy, modify, and distribute this software for any |
| * purpose with or without fee is hereby granted, provided that the above |
| * copyright notice and this permission notice appear in all copies. |
| * |
| * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES |
| * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF |
| * MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR |
| * ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES |
| * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN |
| * ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF |
| * OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE. |
| */ |
| /* |
| * Evaluate polynomial |
| * |
| * |
| * SYNOPSIS: |
| * |
| * int N; |
| * long double x, y, coef[N+1], polevl[]; |
| * |
| * y = polevll( x, coef, N ); |
| * |
| * |
| * DESCRIPTION: |
| * |
| * Evaluates polynomial of degree N: |
| * |
| * 2 N |
| * y = C + C x + C x +...+ C x |
| * 0 1 2 N |
| * |
| * Coefficients are stored in reverse order: |
| * |
| * coef[0] = C , ..., coef[N] = C . |
| * N 0 |
| * |
| * The function p1evll() assumes that coef[N] = 1.0 and is |
| * omitted from the array. Its calling arguments are |
| * otherwise the same as polevll(). |
| * |
| * |
| * SPEED: |
| * |
| * In the interest of speed, there are no checks for out |
| * of bounds arithmetic. This routine is used by most of |
| * the functions in the library. Depending on available |
| * equipment features, the user may wish to rewrite the |
| * program in microcode or assembly language. |
| * |
| */ |
| |
| #include "libm.h" |
| |
| #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024 |
| #else |
| /* |
| * Polynomial evaluator: |
| * P[0] x^n + P[1] x^(n-1) + ... + P[n] |
| */ |
| long double __polevll(long double x, const long double *P, int n) |
| { |
| long double y; |
| |
| y = *P++; |
| do { |
| y = y * x + *P++; |
| } while (--n); |
| |
| return y; |
| } |
| |
| /* |
| * Polynomial evaluator: |
| * x^n + P[0] x^(n-1) + P[1] x^(n-2) + ... + P[n] |
| */ |
| long double __p1evll(long double x, const long double *P, int n) |
| { |
| long double y; |
| |
| n -= 1; |
| y = x + *P++; |
| do { |
| y = y * x + *P++; |
| } while (--n); |
| |
| return y; |
| } |
| #endif |