1 /* $OpenBSD: polevll.c,v 1.1 2011/07/06 00:02:42 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 /* 66 * Polynomial evaluator: 67 * P[0] x^n + P[1] x^(n-1) + ... + P[n] 68 */ 69 long double 70 __polevll(long double x, void *PP, int n) 71 { 72 long double y; 73 long double *P; 74 75 P = (long double *)PP; 76 y = *P++; 77 do { 78 y = y * x + *P++; 79 } while (--n); 80 81 return (y); 82 } 83 84 /* 85 * Polynomial evaluator: 86 * x^n + P[0] x^(n-1) + P[1] x^(n-2) + ... + P[n] 87 */ 88 long double 89 __p1evll(long double x, void *PP, int n) 90 { 91 long double y; 92 long double *P; 93 94 P = (long double *)PP; 95 n -= 1; 96 y = x + *P++; 97 do { 98 y = y * x + *P++; 99 } while (--n); 100 101 return (y); 102 } 103