ThreePointFormula.php 5.4 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135
  1. <?php
  2. namespace MathPHP\NumericalAnalysis\NumericalDifferentiation;
  3. use MathPHP\Exception;
  4. /**
  5. * Three Point Formula
  6. *
  7. * In numerical analysis, the three point formula is used for approximating
  8. * the derivative of a function at a point in its domain.
  9. *
  10. * We can either directly supply a set of inputs and their corresponding outputs
  11. * for said function, or if we explicitly know the function, we can define it as a
  12. * callback function and then generate a set of points by evaluating that function
  13. * at 3 points between a start and end point.
  14. */
  15. class ThreePointFormula extends NumericalDifferentiation
  16. {
  17. /**
  18. * Use the Three Point Formula to approximate the derivative of a function at
  19. * our $target. Our input can support either a set of arrays, or a callback
  20. * function with arguments (to produce a set of arrays). Each array in our
  21. * input contains two numbers which correspond to coordinates (x, y) or
  22. * equivalently, (x, f(x)), of the function f(x) whose derivative we are
  23. * approximating.
  24. *
  25. * The Three Point Formula requires we supply 3 points that are evenly spaced
  26. * apart, and that our target equals the x-components of one of our 3 points.
  27. *
  28. * Example: differentiation(2, function($x) {return $x**2;}, 0, 4 ,3) will produce
  29. * a set of arrays by evaluating the callback at 3 evenly spaced points
  30. * between 0 and 4. Then, this array will be used in our approximation.
  31. *
  32. * Three Point Formula:
  33. *
  34. * - If the 2nd point is our $target, use the Midpoint Formula:
  35. *
  36. * 1 h²
  37. * f′(x₀) = - [f(x₀+h)-f(x₀-h)] - - f⁽³⁾(ζ₁)
  38. * 2h 6
  39. *
  40. * where ζ₁ lies between x₀ - h and x₀ + h
  41. *
  42. * - If the 1st or 3rd point is our $target, use the Endpoint Formula:
  43. * - Note that when the 3rd point is our $target, we use a negative h.
  44. *
  45. * 1 h²
  46. * f′(x₀) = - [-3f(x₀)+4f(x₀+h)-f(x₀+2h)] + - f⁽³⁾(ζ₀)
  47. * 2h 3
  48. *
  49. * where ζ₀ lies between x₀ and x₀ + 2h
  50. *
  51. * @param float $target
  52. * The value at which we are approximating the derivative
  53. * @param callable|array<array{int|float, int|float}> $source
  54. * The source of our approximation. Should be either
  55. * a callback function or a set of arrays. Each array
  56. * (point) contains precisely two numbers, an x and y.
  57. * Example array: [[1,2], [2,3], [3,4]].
  58. * Example callback: function($x) {return $x**2;}
  59. * @param int|float ...$args
  60. * The arguments of our callback function: start,
  61. * end, and n. Example: differentiate($target, $source, 0, 8, 3).
  62. * If $source is a set of points, do not input any
  63. * $args. Example: approximate($source).
  64. *
  65. * @return float
  66. * The approximation of f'($target), i.e. the derivative
  67. * of our input at our target point
  68. *
  69. * @throws Exception\BadDataException
  70. */
  71. public static function differentiate(float $target, $source, ...$args): float
  72. {
  73. // Get an array of points from our $source argument
  74. $points = self::getPoints($source, $args);
  75. // Validate input, sort points, make sure spacing is constant, and make
  76. // sure our target is contained in an interval supplied by our $source
  77. self::validate($points, $degree = 3);
  78. $sorted = self::sort($points);
  79. self::assertSpacingConstant($sorted);
  80. self::assertTargetInPoints($target, $sorted);
  81. // Descriptive constants
  82. $x = self::X;
  83. $y = self::Y;
  84. // Initialize
  85. $h = ($sorted[2][$x] - $sorted[0][$x]) / 2;
  86. /*
  87. * If the 2nd point is our $target, use the Midpoint Formula:
  88. *
  89. * 1 h²
  90. * f′(x₀) = - [f(x₀+h)-f(x₀-h)] - - f⁽³⁾(ζ₁)
  91. * 2h 6
  92. *
  93. * where ζ₁ lies between x₀ - h and x₀ + h
  94. */
  95. if ($sorted[1][$x] == $target) {
  96. $f⟮x₀⧿h⟯ = $sorted[0][$y];
  97. $f⟮x₀⧾h⟯ = $sorted[2][$y];
  98. $derivative = ($f⟮x₀⧾h⟯ - $f⟮x₀⧿h⟯) / (2 * $h);
  99. return $derivative;
  100. }
  101. /*
  102. * If the 1st or 3rd point is our $target, use the Endpoint Formula:
  103. * Note that when the 3rd point is our $target, we use a negative h.
  104. *
  105. * 1 h²
  106. * f′(x₀) = - [-3f(x₀)+4f(x₀+h)-f(x₀+2h)] + - f⁽³⁾(ζ₀)
  107. * 2h 3
  108. *
  109. * where ζ₀ lies between x₀ and x₀ + 2h
  110. */
  111. if ($sorted[0][$x] == $target) { // The 1st point is our $target
  112. $f⟮x₀⟯ = $sorted[0][$y];
  113. $f⟮x₀⧾h⟯ = $sorted[1][$y];
  114. $f⟮x₀⧾2h⟯ = $sorted[2][$y];
  115. } else { // The 3rd point is our $target, use negative h
  116. $h = -$h;
  117. $f⟮x₀⟯ = $sorted[2][$y];
  118. $f⟮x₀⧾h⟯ = $sorted[1][$y];
  119. $f⟮x₀⧾2h⟯ = $sorted[0][$y];
  120. }
  121. $derivative = (-3 * $f⟮x₀⟯ + 4 * $f⟮x₀⧾h⟯ - $f⟮x₀⧾2h⟯) / (2 * $h);
  122. return $derivative;
  123. }
  124. }