AdventOfCode2023/UPolynomialRoots.pas

73 lines
2.0 KiB
Plaintext

{
Solutions to the Advent Of Code.
Copyright (C) 2024 Stefan Müller
This program is free software: you can redistribute it and/or modify it under
the terms of the GNU General Public License as published by the Free Software
Foundation, either version 3 of the License, or (at your option) any later
version.
This program is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with
this program. If not, see <http://www.gnu.org/licenses/>.
}
unit UPolynomialRoots;
{$mode ObjFPC}{$H+}
interface
uses
Classes, SysUtils, UPolynomial, UBigInt;
type
{ TRootIsolation }
TRootIsolation = class
private
function CalcSimpleRootBound(constref APolynomial: TBigIntPolynomial): TBigInt;
public
function Bisect(constref APolynomial: TBigIntPolynomial): Int64;
end;
implementation
{ TRootIsolation }
function TRootIsolation.CalcSimpleRootBound(constref APolynomial: TBigIntPolynomial): TBigInt;
var
i, sign: Integer;
a: TBigInt;
begin
// We need a_n > 0 here, so we use -sign(a_n) instead of actually flipping the polynomial.
// Sign is not 0 because a_n is not 0.
sign := -APolynomial.Coefficient[APolynomial.Degree].Sign;
// This is a simplification of Cauchy's bound to avoid division.
// https://en.wikipedia.org/wiki/Geometrical_properties_of_polynomial_roots#Bounds_of_positive_real_roots
Result := TBigInt.Zero;
for i := 0 to APolynomial.Degree - 1 do begin
a := sign * APolynomial.Coefficient[i];
if Result < a then
Result := a;
end;
Result := Result + 1;
end;
function TRootIsolation.Bisect(constref APolynomial: TBigIntPolynomial): Int64;
var
bound: TBigInt;
p: TBigIntPolynomial;
begin
bound := CalcSimpleRootBound(APolynomial);
p := APolynomial.ScaleVariable(bound);
end;
end.