129 lines
3.9 KiB
PHP
129 lines
3.9 KiB
PHP
<?php
|
||
|
||
namespace App\Support;
|
||
|
||
/**
|
||
* Regresja wielomianowa - wierny port legacy polynomial/PolynomialRegression.php
|
||
* (Andrew Que, GPL). Ten sam algorytm (równania normalne + eliminacja
|
||
* wierszowa), arytmetyka natywna (float) zamiast bcmath - wynik identyczny dla
|
||
* zakresów danych pomp, bez zależności od ext-bcmath.
|
||
*
|
||
* Użycie jak w legacy: new PolynomialRegression(4) → 4 współczynniki (wielomian
|
||
* stopnia 3); addData(x,y) ×N; getCoefficients(); interpolate(coeffs, x).
|
||
*/
|
||
class PolynomialRegression
|
||
{
|
||
/** @var float[] */
|
||
private array $xPowers = [];
|
||
/** @var float[] */
|
||
private array $xyPowers = [];
|
||
|
||
public function __construct(private int $degree = 3)
|
||
{
|
||
$this->reset();
|
||
}
|
||
|
||
public function reset(): void
|
||
{
|
||
$squares = ($this->degree - 1) * 2;
|
||
for ($i = 0; $i <= $squares; $i++) {
|
||
$this->xPowers[$i] = 0.0;
|
||
$this->xyPowers[$i] = 0.0;
|
||
}
|
||
}
|
||
|
||
public function addData(float $x, float $y): void
|
||
{
|
||
$squares = ($this->degree - 1) * 2;
|
||
for ($i = 0; $i <= $squares; $i++) {
|
||
$this->xPowers[$i] += $x ** $i;
|
||
$this->xyPowers[$i] += $y * ($x ** $i);
|
||
}
|
||
}
|
||
|
||
/**
|
||
* Współczynniki [c0, c1, ...] (c0 = wyraz wolny).
|
||
*
|
||
* @return float[]
|
||
*/
|
||
public function getCoefficients(?int $degree = null): array
|
||
{
|
||
$degree ??= $this->degree;
|
||
$eps = 1e-12;
|
||
|
||
// Macierz równań normalnych + kolumna augmentacji.
|
||
$matrix = [];
|
||
for ($row = 0; $row < $degree; $row++) {
|
||
for ($col = 0; $col < $degree; $col++) {
|
||
$matrix[$row][$col] = $this->xPowers[$row + $col];
|
||
}
|
||
$matrix[$row][$degree] = $this->xyPowers[$row];
|
||
}
|
||
|
||
$rows = count($matrix);
|
||
$isDone = array_fill(0, $rows, false);
|
||
$order = [];
|
||
|
||
// Eliminacja wierszowa do macierzy górnotrójkątnej z 1 na przekątnej.
|
||
for ($column = 0; $column < $rows; $column++) {
|
||
$activeRow = 0;
|
||
while ($activeRow < $rows
|
||
&& (abs($matrix[$activeRow][$column]) < $eps || $isDone[$activeRow])) {
|
||
$activeRow++;
|
||
}
|
||
|
||
if ($activeRow < $rows) {
|
||
$order[$column] = $activeRow;
|
||
|
||
$firstTerm = $matrix[$activeRow][$column];
|
||
for ($sub = $column; $sub <= $rows; $sub++) {
|
||
$matrix[$activeRow][$sub] /= $firstTerm;
|
||
}
|
||
$isDone[$activeRow] = true;
|
||
|
||
for ($row = 0; $row < $rows; $row++) {
|
||
if (! $isDone[$row] && abs($matrix[$row][$column]) >= $eps) {
|
||
$factor = $matrix[$row][$column];
|
||
for ($sub = $column; $sub <= $rows; $sub++) {
|
||
$matrix[$row][$sub] -= $factor * $matrix[$activeRow][$sub];
|
||
}
|
||
}
|
||
}
|
||
}
|
||
}
|
||
|
||
// Podstawianie wsteczne.
|
||
$isDone = array_fill(0, $rows, false);
|
||
$coefficients = [];
|
||
for ($column = $rows - 1; $column >= 0; $column--) {
|
||
$activeRow = $order[$column];
|
||
$isDone[$activeRow] = true;
|
||
|
||
for ($row = 0; $row < $rows; $row++) {
|
||
if (! $isDone[$row]) {
|
||
$factor = $matrix[$row][$column];
|
||
for ($sub = $column; $sub <= $rows; $sub++) {
|
||
$matrix[$row][$sub] -= $factor * $matrix[$activeRow][$sub];
|
||
}
|
||
}
|
||
}
|
||
$coefficients[$column] = $matrix[$activeRow][$rows];
|
||
}
|
||
|
||
ksort($coefficients);
|
||
|
||
return $coefficients;
|
||
}
|
||
|
||
/** @param float[] $coefficients */
|
||
public function interpolate(array $coefficients, float $x): float
|
||
{
|
||
$y = 0.0;
|
||
foreach ($coefficients as $i => $c) {
|
||
$y += $c * ($x ** $i);
|
||
}
|
||
|
||
return $y;
|
||
}
|
||
}
|