15 January 2017 09:35:18 AM DIVDIF_PRB C version Test the DIVDIF library. TEST01 DATA_TO_DIF_DISPLAY sets up a difference table and displays intermediate calculations; DIF_APPEND appends a new data point; DIF_ANTIDERIV computes the antiderivative; DIF_DERIV_TABLE computes the derivative; DIF_SHIFT_ZERO shifts all the abscissas to 0; DIF_VAL evaluates at a point; The divided difference table: 1.000000 2.000000 3.000000 4.000000 0 1.000000 4.000000 9.000000 16.000000 1 3.000000 5.000000 7.000000 2 1.000000 1.000000 3 0.000000 The divided difference polynomial: p(x) = 1.000000 + ( x - 1.000000 ) * ( 3.000000 + ( x - 2.000000 ) * ( 1.000000 + ( x - 3.000000 ) * ( 0.000000 ))) DIF_APPEND can add the data (5,25) to the table. The updated divided difference polynomial: p(x) = 25.000000 + ( x - 5.000000 ) * ( 6.000000 + ( x - 1.000000 ) * ( 1.000000 + ( x - 2.000000 ) * ( -0.000000 + ( x - 3.000000 ) * ( -0.000000 )))) DIF_VAL can evaluate the table at a point. DIF_VAL reports P(2.500000) = 6.250000 The divided difference table after DIF_SHIFT_ZERO: p(x) = 0.000000 + ( x - 0.000000 ) * ( 0.000000 + ( x - 0.000000 ) * ( 1.000000 + ( x - 0.000000 ) * ( 0.000000 + ( x - 0.000000 ) * ( -0.000000 )))) The divided difference table for the derivative: p(x) = 0.000000 + ( x - 0.000000 ) * ( 2.000000 + ( x - 0.000000 ) * ( 0.000000 + ( x - 0.000000 ) * ( -0.000000 ))) DIF_VAL reports P'(2.500000) = 5.000000 The divided difference table for the antiderivative: p(x) = 0.000000 + ( x - 0.000000 ) * ( 0.000000 + ( x - 0.000000 ) * ( 0.000000 + ( x - 0.000000 ) * ( 0.333333 + ( x - 0.000000 ) * ( 0.000000 + ( x - 0.000000 ) * ( -0.000000 ))))) DIF_VAL reports (Anti)P(2.500000) = 5.208333 TEST02 Approximate Y = EXP(X) using orders 1 to 8. Original data: X Y 0.000000 1.000000 1.000000 2.718282 2.000000 7.389056 3.000000 20.085537 4.000000 54.598150 5.000000 148.413159 6.000000 403.428793 7.000000 1096.633158 Evaluate at X = 2.500000 where EXP(X) = 12.182494 Order Approximate Y Error 1 1.000000 -11.182494 2 5.295705 -6.886789 3 10.831628 -1.350866 4 12.417007 0.234513 5 12.076491 -0.106003 6 12.252022 0.069528 7 12.126351 -0.056143 8 12.234320 0.051826 TEST03 DIF_BASIS computes Lagrange basis polynomials in difference form. The base points: 1 1.000000 2 2.000000 3 3.000000 4 4.000000 5 5.000000 The table of difference vectors defining the basis polynomials. Each ROW represents a polynomial. 1.000000 -1.000000 0.500000 -0.166667 0.041667 0.000000 1.000000 -1.000000 0.500000 -0.166667 0.000000 0.000000 0.500000 -0.500000 0.250000 0.000000 0.000000 0.000000 0.166667 -0.166667 0.000000 0.000000 0.000000 0.000000 0.041667 Evaluate basis polynomial #3 at a set of points. X Y 1.000000 0.000000 1.500000 -0.546875 2.000000 0.000000 2.500000 0.703125 3.000000 1.000000 3.500000 0.703125 4.000000 0.000000 4.500000 -0.546875 5.000000 0.000000 TEST05 DIF_TO_R8POLY converts a difference table to a polynomial; DIF_SHIFT_ZERO shifts a divided difference table to 0; These are equivalent operations The divided difference table: 1.000000 2.000000 3.000000 4.000000 0 -2.000000 2.000000 14.000000 40.000000 1 4.000000 12.000000 26.000000 2 4.000000 7.000000 3 1.000000 The divided difference table: 1.000000 2.000000 3.000000 4.000000 0 -2.000000 2.000000 14.000000 40.000000 1 4.000000 12.000000 26.000000 2 4.000000 7.000000 3 1.000000 The divided difference table: p(x) = -2.000000 + ( x - 1.000000 ) * ( 4.000000 + ( x - 2.000000 ) * ( 4.000000 + ( x - 3.000000 ) * ( 1.000000 ))) The polynomial using DIF_SHIFT_ZERO: p(x) = 1.000000 * x ^ 3 - 2.000000 * x ^ 2 + 3.000000 * x - 4.000000 The polynomial using DIF_TO_R8POLY: p(x) = 1.000000 * x ^ 3 - 2.000000 * x ^ 2 + 3.000000 * x - 4.000000 TEST06 R8POLY_ANT_COF computes the coefficients of the antiderivative of a polynomial; R8POLY_ANT_VAL evaluates the antiderivative of a polynomial; R8POLY_DER_COF computes the coefficients of the derivative of a polynomial; R8POLY_DER_VAL evaluates the derivative of a polynomial; R8POLY_PRINT prints a polynomial; R8POLY_VAL evaluates a polynomial. The initial polynomial: p(x) = 5.000000 * x ^ 4 + 4.000000 * x ^ 3 + 3.000000 * x ^ 2 + 2.000000 * x + 1.000000 The antiderivative polynomial: p(x) = 1.000000 * x ^ 5 + 1.000000 * x ^ 4 + 1.000000 * x ^ 3 + 1.000000 * x ^ 2 + 1.000000 * x The derivative polynomial: p(x) = 20.000000 * x ^ 3 + 12.000000 * x ^ 2 + 6.000000 * x + 2.000000 Evaluate the polynomial, antiderivative and derivative, using only the original polynomial coefficients: X P(X) Anti_P(X) P'(X) 0.000000 30.000000 0.000000 2.000000 1.000000 44.000000 34.000000 40.000000 2.000000 158.000000 120.000000 222.000000 TEST07 R8POLY_BASIS computes Lagrange basis polynomials in standard form. The table of difference vectors defining the basis polynomials. Each ROW represents a polynomial. 5.000000 -6.416667 2.958333 -0.583333 0.041667 -10.000000 17.833333 -9.833333 2.166667 -0.166667 10.000000 -19.500000 12.250000 -3.000000 0.250000 -5.000000 10.166667 -6.833333 1.833333 -0.166667 1.000000 -2.083333 1.458333 -0.416667 0.041667 One basis polynomial in standard form: p(x) = 0.250000 * x ^ 4 - 3.000000 * x ^ 3 + 12.250000 * x ^ 2 - 19.500000 * x + 10.000000 Evaluate basis polynomial #3 at a set of points. X Y 1.000000 0.000000 1.500000 -0.546875 2.000000 0.000000 2.500000 0.703125 3.000000 1.000000 3.500000 0.703125 4.000000 0.000000 4.500000 -0.546875 5.000000 0.000000 TEST08 R8POLY_SHIFT shifts polynomial coefficients. Polynomial coefficients for argument X 0 6.000000 1 -1.000000 2 2.000000 SCALE = 2.000000 SHIFT = 3.000000 Polynomial coefficients for argument Z = SCALE * X + SHIFT: 0 12.000000 1 -3.500000 2 0.500000 TEST16 NCC_RULE computes closed Newton Cotes formulas; Newton-Cotes Closed Quadrature Rule: Abscissa Weight 1 -1.000000 0.086921 2 -0.714286 0.414005 3 -0.428571 0.153125 4 -0.142857 0.345949 5 0.142857 0.345949 6 0.428571 0.153125 7 0.714286 0.414005 8 1.000000 0.086921 TEST17 NCO_RULE computes open Newton Cotes formulas. Newton-Cotes Open Quadrature Rule: Abscissa Weight 1 -0.777778 0.797768 2 -0.555556 -1.251339 3 -0.333333 2.217411 4 -0.111111 -0.763839 5 0.111111 -0.763839 6 0.333333 2.217411 7 0.555556 -1.251339 8 0.777778 0.797768 TEST18 ROOTS_TO_DIF computes a divided difference polynomial with the given roots; DIF_TO_R8POLY converts it to a standard form polynomial. The roots: 1 3.000000 The polynomial: p(x) = 1.000000 * x - 3.000000 The roots: 1 3.000000 2 1.000000 The polynomial: p(x) = 1.000000 * x ^ 2 - 4.000000 * x + 3.000000 The roots: 1 3.000000 2 1.000000 3 2.000000 The polynomial: p(x) = 1.000000 * x ^ 3 - 6.000000 * x ^ 2 + 11.000000 * x - 6.000000 The roots: 1 3.000000 2 1.000000 3 2.000000 4 4.000000 The polynomial: p(x) = 1.000000 * x ^ 4 - 10.000000 * x ^ 3 + 35.000000 * x ^ 2 - 50.000000 * x + 24.000000 TEST19 ROOTS_TO_R8POLY computes polynomial coefficients from roots. The roots: 1 3.000000 The polynomial: p(x) = 1.000000 * x - 3.000000 The roots: 1 3.000000 2 1.000000 The polynomial: p(x) = 1.000000 * x ^ 2 - 4.000000 * x + 3.000000 The roots: 1 3.000000 2 1.000000 3 2.000000 The polynomial: p(x) = 1.000000 * x ^ 3 - 6.000000 * x ^ 2 + 11.000000 * x - 6.000000 The roots: 1 3.000000 2 1.000000 3 2.000000 4 4.000000 The polynomial: p(x) = 1.000000 * x ^ 4 - 10.000000 * x ^ 3 + 35.000000 * x ^ 2 - 50.000000 * x + 24.000000 TEST20 For a divided difference polynomial: DIF_DERIVK_TABLE computes the K-th derivative; The divided difference polynomial P0: p(x) = 0.333333 + ( x - -2.000000 ) * ( 0.041667 + ( x - -1.000000 ) * ( 0.291667 + ( x - 0.000000 ) * ( 0.083333 + ( x - 1.000000 ) * ( 0.041667 )))) Using DIF_TO_R8POLY p(x) = 0.041667 * x ^ 4 + 0.166667 * x ^ 3 + 0.500000 * x ^ 2 + 1.000000 * x + 1.000000 Evaluate difference tables for the function P0 and its first four derivatives, P1...P4. X P0 P1 P2 P3 P4 0.0000 1.0000 1.0000 1.0000 1.0000 1.0000 0.2000 1.2214 1.2213 1.2200 1.2000 1.0000 0.4000 1.4917 1.4907 1.4800 1.4000 1.0000 0.6000 1.8214 1.8160 1.7800 1.6000 1.0000 0.8000 2.2224 2.2053 2.1200 1.8000 1.0000 1.0000 2.7083 2.6667 2.5000 2.0000 1.0000 1.2000 3.2944 3.2080 2.9200 2.2000 1.0000 1.4000 3.9974 3.8373 3.3800 2.4000 1.0000 1.6000 4.8357 4.5627 3.8800 2.6000 1.0000 1.8000 5.8294 5.3920 4.4200 2.8000 1.0000 2.0000 7.0000 6.3333 5.0000 3.0000 1.0000 TEST21 DIF_BASIS_DERIV computes difference tables for the first derivative of each Lagrange basis. Lagrange basis derivative polynomial coefficients: Row: 0 1 Col 0: -0.285714 0.095238 1: 0.250000 -0.166667 2: 0.035714 0.071429 P1'=-(2x-6)/21 p(x) = 0.095238 * x - 0.285714 P2'=-(2x-3)/12 p(x) = - 0.166667 * x + 0.250000 P3'=(2x+1)/28 p(x) = 0.071429 * x + 0.035714 TEST22 DIF_BASIS_DERIVK computes difference tables for the K-th derivative of each Lagrange basis. Lagrange basis K-th derivative polynomial coefficients: Row: 0 1 2 Col 0: 5.916667 -3.500000 0.500000 1: -19.666667 13.000000 -2.000000 2: 24.500000 -18.000000 3.000000 3: -13.666667 11.000000 -2.000000 4: 2.916667 -2.500000 0.500000 P1''=(12x^2-84x+142)/24 p(x) = 0.500000 * x ^ 2 - 3.500000 * x + 5.916667 P2''=-2x^2+13x-59/3 p(x) = - 2.000000 * x ^ 2 + 13.000000 * x - 19.666667 P3''=3x^2-18x+49/2 p(x) = 3.000000 * x ^ 2 - 18.000000 * x + 24.500000 P4''=-2x^2+11x-41/3 p(x) = - 2.000000 * x ^ 2 + 11.000000 * x - 13.666667 P5''=(6x^2-30x+35)/12 p(x) = 0.500000 * x ^ 2 - 2.500000 * x + 2.916667 DIVDIF_PRB Normal end of execution. 15 January 2017 09:35:18 AM