/************************************************************************** * Test program for Bauhuber's method * * ----------------------------------------------------------------------- * * This program uses Bauhuber's Method to find all real or complex * * roots of a polynomial of degree n: * * n-1 n * * P(x) = a[0] + a[1] * x + ... + a[n-1] * x + a[n] * x , * * * * where a[i], i=0, ..., n, can be complex. * * ----------------------------------------------------------------------- * * SAMPLE RUN: * * (Find all roots of polynomial of degree 4: * * P(x) = 0.000089248 - 0.04368 X + 2.9948 X2 - 6.798 X3 + X4) * * * * -------------------------------------------------------------------- * * Complex and Real Roots of a Polynomial (Bauhuber's method) * * -------------------------------------------------------------------- * * * * Polynomial coefficients: * * * * a[ 0] = 0.0000892480 0.0000000000 * * a[ 1] = -0.0436800000 0.0000000000 * * a[ 2] = 2.9948000000 0.0000000000 * * a[ 3] = -6.7980000000 0.0000000000 * * a[ 4] = 1.0000000000 0.0000000000 * * * * Roots (without scaling) * * No. Real part Imaginary part function value * * * * 0 2.4537212366448121e-003 0.0000000000000000e+000 1.355253e-020 * * 1 1.2573281756029101e-002 0.0000000000000000e+000 4.065758e-020 * * 2 4.5731893744203350e-001 0.0000000000000000e+000 4.540097e-018 * * 3 6.3256540595652933e+000 0.0000000000000000e+000 1.618457e-013 * * * * * * Roots (with scaling) * * No. Real part Imaginary part function value * * * * 0 2.4537212366448130e-003 0.0000000000000000e+000 0.000000e+000 * * 1 1.2573281756029099e-002 0.0000000000000000e+000 2.710505e-020 * * 2 4.5731893744203356e-001 0.0000000000000000e+000 1.382358e-018 * * 3 6.3256540595652924e+000 0.0000000000000000e+000 5.143054e-014 * * * * -------------------------------------------------------------------- * * Ref.: "Numerical algorithms with C, By Gisela Engeln-Muellges and * * Frank Uhlig, Springer-Verlag, 1996" [BIBLI 11]. * * * * C++ Release By J-P Moreau, Paris. * * (www.jpmoreau.fr) * *************************************************************************** Note: to link with basis_r.cpp, fbauhube.cpp, vmblock.cpp --------------------------------------------------------- */ #include // for general declarations #include // for memory dynamic allocations int WriteHead1( char *); int WriteEnd1(void); void SetVec (int, REAL *, REAL); int bauhub /* Bauhuber's method for complex polynomials .*/ ( int real, /* Are the coefficients real ? .....*/ int scale, /* Scaling ? .......................*/ int n, /* degree of polynomial ............*/ REAL ar[], /* Real parts of coefficients ......*/ REAL ai[], /* Imaginary parts, coefficients ...*/ REAL rootr[], /* Real parts of roots .............*/ REAL rooti[], /* Imaginary parts of roots ........*/ REAL absf[] /* Absolute value of function values*/ ); void main () { REAL *ar, *ai, *rootr, *rooti, *val; int i, j, n, rc, skala; void *vmblock; WriteHead1(" Complex and Real Roots of a Polynomial (Bauhuber's method)"); n = 4; //order of polynomial vmblock = vminit(); ar = (REAL *) vmalloc(vmblock, VEKTOR, n + 1, 0); ai = (REAL *) vmalloc(vmblock, VEKTOR, n + 1, 0); rootr = (REAL *) vmalloc(vmblock, VEKTOR, n + 1, 0); rooti = (REAL *) vmalloc(vmblock, VEKTOR, n + 1, 0); val = (REAL *) vmalloc(vmblock, VEKTOR, n + 1, 0); if (! vmcomplete(vmblock)) LogError ("No Memory", 0, __FILE__, __LINE__); //define ar vector ar[0] = 0.000089248; ar[1] = -0.04368; ar[2] = 2.9948; ar[3] = -6.798; ar[4] = ONE; //ai vector is null SetVec (n + 1, ai, ZERO); printf ("Polynomial coefficients:\n"); for (i = 0; i <= n; i++) { printf("\n a[%2d] = ", i); printf (FORMAT_2010LF, ar[i]); printf (FORMAT_2010LF, ai[i]); } for (j = 0; j < 2; j++) { skala = j; rc = bauhub (0, skala, n, ar, ai, rootr, rooti, val); if (rc == 0) { if (skala == 0) printf ("\n\nRoots (without scaling)\n"); else printf ("\n\nRoots (with scaling)\n"); printf (" No. Real part Imaginary part " "function value\n"); for (i = 0; i < n; i++) { printf (" \n %4d ", i); printf (FORMAT_2016LE, rootr[i]); printf (FORMAT_2016LE, rooti[i]); printf (FORMAT_LE, val[i]); } printf ("\n"); } else LogError ("bauhube", rc, __FILE__, __LINE__); } WriteEnd1(); } //end of file tbauhube.cpp