NewtonPolynomialForwardTest.php 3.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149
  1. <?php
  2. namespace MathPHP\Tests\NumericalAnalysis\Interpolation;
  3. use MathPHP\NumericalAnalysis\Interpolation\NewtonPolynomialForward;
  4. class NewtonPolynomialForwardTest extends \PHPUnit\Framework\TestCase
  5. {
  6. /**
  7. * @test Interpolated piecewise function computes expected values: p(x) = expected
  8. * @dataProvider dataProviderForPolynomialAgrees
  9. * @param int $x
  10. * @param int $expected
  11. * @throws \Exception
  12. */
  13. public function testPolynomialAgrees(int $x, int $expected)
  14. {
  15. // Given
  16. $points = [[0, 0], [1, 5], [3, 2], [7, 10], [10, -4]];
  17. // And
  18. $p = NewtonPolynomialForward::interpolate($points);
  19. // When
  20. $evaluated = $p($x);
  21. // Then
  22. $this->assertEqualsWithDelta($expected, $evaluated, 0.00001);
  23. }
  24. /**
  25. * @return array (x, expected)
  26. */
  27. public function dataProviderForPolynomialAgrees(): array
  28. {
  29. return [
  30. [0, 0], // p(0) = 0
  31. [1, 5], // p(1) = 5
  32. [3, 2], // p(3) = 2
  33. [7, 10], // p(7) = 10
  34. [10, -4], // p(10) = -4
  35. ];
  36. }
  37. /**
  38. * @test Solve
  39. * @dataProvider dataProviderForSolve
  40. * @param float $x
  41. * @throws \Exception
  42. *
  43. * f(x) = x⁴ + 8x³ -13x² -92x + 96
  44. *
  45. * Given n points, the error in the Newton Polynomials is proportional
  46. * to the max value of the nth derivative. Thus, if we if interpolate n at
  47. * 6 points, the 5th derivative of our original function f(x) = 0, and so
  48. * our resulting polynomial will have no error.
  49. *
  50. * p(x) agrees with f(x) at x = $_
  51. */
  52. public function testSolve($x)
  53. {
  54. // f(x) = x⁴ + 8x³ -13x² -92x + 96
  55. $f = function ($x) {
  56. return $x ** 4 + 8 * $x ** 3 - 13 * $x ** 2 - 92 * $x + 96;
  57. };
  58. // And
  59. $a = 0;
  60. $b = 10;
  61. $n = 5;
  62. // And
  63. $p = NewtonPolynomialForward::interpolate($f, $a, $b, $n);
  64. $expected = $f($x);
  65. // When
  66. $actual = $p($x);
  67. // Then
  68. $this->assertEquals($expected, $actual);
  69. }
  70. /**
  71. * @test Solve
  72. * @dataProvider dataProviderForSolve
  73. * @param float $x
  74. * @throws \Exception
  75. *
  76. * f(x) = x⁴ + 8x³ -13x² -92x + 96
  77. *
  78. * The error is bounded by:
  79. * |f(x)-p(x)| = tol <= (max f⁽ⁿ⁺¹⁾(x))*(x-x₀)*...*(x-xn)/(n+1)!
  80. *
  81. * f'(x) = 4x³ +24x² -26x - 92
  82. * f''(x) = 12x² - 48x - 26
  83. * f'''(x) = 24x - 48
  84. * f⁽⁴⁾(x) = 24
  85. *
  86. * So, tol <= 24*(x-x₀)*...*(x-xn)/(4!) = (x-x₀)*...*(x-xn) where
  87. *
  88. * p(x) agrees with f(x) at x = $_
  89. */
  90. public function testSolveNonZeroError($x)
  91. {
  92. // f(x) = x⁴ + 8x³ -13x² -92x + 96
  93. $f = function ($x) {
  94. return $x ** 4 + 8 * $x ** 3 - 13 * $x ** 2 - 92 * $x + 96;
  95. };
  96. // And
  97. $a = 0;
  98. $b = 10;
  99. $n = 5;
  100. // and
  101. $x₀ = 0;
  102. $x₁ = 3;
  103. $x₂ = 6;
  104. $x₃ = 9;
  105. $tol = \abs(($x - $x₀) * ($x - $x₁) * ($x - $x₂) * ($x - $x₃));
  106. // And
  107. $p = NewtonPolynomialForward::interpolate($f, $a, $b, $n);
  108. $expected = $f($x);
  109. // When
  110. $actual = $p($x);
  111. // Then
  112. $this->assertEqualsWithDelta($expected, $actual, $tol);
  113. }
  114. /**
  115. * @return array p(x) agrees with f(x) at x = $_
  116. */
  117. public function dataProviderForSolve(): array
  118. {
  119. return [
  120. [0],
  121. [2],
  122. [4],
  123. [6],
  124. [8],
  125. [10],
  126. [-90],
  127. [-99],
  128. ];
  129. }
  130. }