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 
64 /*
65  * Polynomial evaluator:
66  *  P[0] x^n  +  P[1] x^(n-1)  +  ...  +  P[n]
67  */
68 long double
__polevll(long double x,const long double * P,int n)69 __polevll(long double x, const long double *P, int n)
70 {
71 	long double y;
72 
73 	y = *P++;
74 	do {
75 		y = y * x + *P++;
76 	} while (--n);
77 
78 	return (y);
79 }
80 
81 /*
82  * Polynomial evaluator:
83  *  x^n  +  P[0] x^(n-1)  +  P[1] x^(n-2)  +  ...  +  P[n]
84  */
85 long double
__p1evll(long double x,const long double * P,int n)86 __p1evll(long double x, const long double *P, int n)
87 {
88 	long double y;
89 
90 	n -= 1;
91 	y = x + *P++;
92 	do {
93 		y = y * x + *P++;
94 	} while (--n);
95 
96 	return (y);
97 }
98