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