Horner Method for polynomial Computation

Source: Internet
Author: User
The Horner algorithm is a fast algorithm named after William George Horner, a British mathematician. In fact, it has been used by Paolo ruffini before. Chinese mathematician Qin Jiuyi proposed this algorithm more than William George Horner 600 years earlier.

P (x) is a polynomial:

We want to calculate the polynomial value p (x0) when X is a special value x0 ).

Construct a sequence:

Then the B0 value of this sequence is the polynomial value.

The program implementation is as follows:

double horner(double p[], int n, double x){    double sum;    sum = p[--n];    while ( n > 0 )    {        sum = p[--n] + sum * x;    }    return sum;}

I often need to design FIR and IIR filters. After designing the filters, I need to verify whether the frequency characteristics are correct. In this case, we need to calculate the result of the real system number Polynomial in the case of X to obtain the complex value. Therefore, the following code is available:

 double _Complex horner_C(double p[], int n, double _Complex x) {    double _Complex sum;    sum = p[--n];    while ( n > 0 )    {        sum = p[--n] + sum * x;    }    return sum; }

This Code uses the new features supported by the plural in c99 and requires the compiler to support the c99 standard.

Finally, a test code is provided. The polynomial P constitutes a FIR low-pass filter, which is designed by using the ffilt function in SCILAB.

The code for constructing a filter is as follows:

ffilt('lp',20 , 0.2 , 0.5);
The following code calculates the frequency response of the filter.

#include <stdio.h>#include <stdlib.h>#include <math.h>#include <complex.h>double horner(double p[], int n, double x);double _Complex horner_C(double p[], int n, double _Complex x);int main(){    int i;    double p[] = {    -0.0196945,    -0.0356154,    1.559E-17,    0.0465740,    0.0340178,    -0.0415773,    -0.0864945,    1.559E-17,    0.2018205,    0.3741957,    0.3741957,    0.2018205,    1.559E-17,    -0.0864945,    -0.0415773,    0.0340178,    0.0465740,    1.559E-17,    -0.0356154,    -0.0196945    };    double x[201], mag[201], pha[201];    double _Complex xx;    for(i = 0; i <= 200; i++)    {        x[i] = 0.5 * i / 200.0;        xx = cexp(-2 * M_PI * x[i] * I);        mag[i] = cabs(horner_C(p, 20, xx));        pha[i] = carg(horner_C(p, 20, xx));        printf("%f, %f, %f\n", x[i] , mag[i], pha[i]);    }    return 0;}

The output results show the following amplitude-frequency characteristics:

The amplitude-frequency characteristics calculated by SCILAB are as follows:

We can see that the two images are exactly the same, indicating that our code calculation results are correct.

Contact Us

The content source of this page is from Internet, which doesn't represent Alibaba Cloud's opinion; products and services mentioned on that page don't have any relationship with Alibaba Cloud. If the content of the page makes you feel confusing, please write us an email, we will handle the problem within 5 days after receiving your email.

If you find any instances of plagiarism from the community, please send an email to: info-contact@alibabacloud.com and provide relevant evidence. A staff member will contact you within 5 working days.

A Free Trial That Lets You Build Big!

Start building with 50+ products and up to 12 months usage for Elastic Compute Service

  • Sales Support

    1 on 1 presale consultation

  • After-Sales Support

    24/7 Technical Support 6 Free Tickets per Quarter Faster Response

  • Alibaba Cloud offers highly flexible support services tailored to meet your exact needs.