1570af302Sopenharmony_ci/* origin: OpenBSD /usr/src/lib/libm/src/polevll.c */ 2570af302Sopenharmony_ci/* 3570af302Sopenharmony_ci * Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net> 4570af302Sopenharmony_ci * 5570af302Sopenharmony_ci * Permission to use, copy, modify, and distribute this software for any 6570af302Sopenharmony_ci * purpose with or without fee is hereby granted, provided that the above 7570af302Sopenharmony_ci * copyright notice and this permission notice appear in all copies. 8570af302Sopenharmony_ci * 9570af302Sopenharmony_ci * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES 10570af302Sopenharmony_ci * WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF 11570af302Sopenharmony_ci * MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR 12570af302Sopenharmony_ci * ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES 13570af302Sopenharmony_ci * WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN 14570af302Sopenharmony_ci * ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF 15570af302Sopenharmony_ci * OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE. 16570af302Sopenharmony_ci */ 17570af302Sopenharmony_ci/* 18570af302Sopenharmony_ci * Evaluate polynomial 19570af302Sopenharmony_ci * 20570af302Sopenharmony_ci * 21570af302Sopenharmony_ci * SYNOPSIS: 22570af302Sopenharmony_ci * 23570af302Sopenharmony_ci * int N; 24570af302Sopenharmony_ci * long double x, y, coef[N+1], polevl[]; 25570af302Sopenharmony_ci * 26570af302Sopenharmony_ci * y = polevll( x, coef, N ); 27570af302Sopenharmony_ci * 28570af302Sopenharmony_ci * 29570af302Sopenharmony_ci * DESCRIPTION: 30570af302Sopenharmony_ci * 31570af302Sopenharmony_ci * Evaluates polynomial of degree N: 32570af302Sopenharmony_ci * 33570af302Sopenharmony_ci * 2 N 34570af302Sopenharmony_ci * y = C + C x + C x +...+ C x 35570af302Sopenharmony_ci * 0 1 2 N 36570af302Sopenharmony_ci * 37570af302Sopenharmony_ci * Coefficients are stored in reverse order: 38570af302Sopenharmony_ci * 39570af302Sopenharmony_ci * coef[0] = C , ..., coef[N] = C . 40570af302Sopenharmony_ci * N 0 41570af302Sopenharmony_ci * 42570af302Sopenharmony_ci * The function p1evll() assumes that coef[N] = 1.0 and is 43570af302Sopenharmony_ci * omitted from the array. Its calling arguments are 44570af302Sopenharmony_ci * otherwise the same as polevll(). 45570af302Sopenharmony_ci * 46570af302Sopenharmony_ci * 47570af302Sopenharmony_ci * SPEED: 48570af302Sopenharmony_ci * 49570af302Sopenharmony_ci * In the interest of speed, there are no checks for out 50570af302Sopenharmony_ci * of bounds arithmetic. This routine is used by most of 51570af302Sopenharmony_ci * the functions in the library. Depending on available 52570af302Sopenharmony_ci * equipment features, the user may wish to rewrite the 53570af302Sopenharmony_ci * program in microcode or assembly language. 54570af302Sopenharmony_ci * 55570af302Sopenharmony_ci */ 56570af302Sopenharmony_ci 57570af302Sopenharmony_ci#include "libm.h" 58570af302Sopenharmony_ci 59570af302Sopenharmony_ci#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024 60570af302Sopenharmony_ci#else 61570af302Sopenharmony_ci/* 62570af302Sopenharmony_ci * Polynomial evaluator: 63570af302Sopenharmony_ci * P[0] x^n + P[1] x^(n-1) + ... + P[n] 64570af302Sopenharmony_ci */ 65570af302Sopenharmony_cilong double __polevll(long double x, const long double *P, int n) 66570af302Sopenharmony_ci{ 67570af302Sopenharmony_ci long double y; 68570af302Sopenharmony_ci 69570af302Sopenharmony_ci y = *P++; 70570af302Sopenharmony_ci do { 71570af302Sopenharmony_ci y = y * x + *P++; 72570af302Sopenharmony_ci } while (--n); 73570af302Sopenharmony_ci 74570af302Sopenharmony_ci return y; 75570af302Sopenharmony_ci} 76570af302Sopenharmony_ci 77570af302Sopenharmony_ci/* 78570af302Sopenharmony_ci * Polynomial evaluator: 79570af302Sopenharmony_ci * x^n + P[0] x^(n-1) + P[1] x^(n-2) + ... + P[n] 80570af302Sopenharmony_ci */ 81570af302Sopenharmony_cilong double __p1evll(long double x, const long double *P, int n) 82570af302Sopenharmony_ci{ 83570af302Sopenharmony_ci long double y; 84570af302Sopenharmony_ci 85570af302Sopenharmony_ci n -= 1; 86570af302Sopenharmony_ci y = x + *P++; 87570af302Sopenharmony_ci do { 88570af302Sopenharmony_ci y = y * x + *P++; 89570af302Sopenharmony_ci } while (--n); 90570af302Sopenharmony_ci 91570af302Sopenharmony_ci return y; 92570af302Sopenharmony_ci} 93570af302Sopenharmony_ci#endif 94