Quelle contfrac.gi Sprache: unbekannt

#############################################################################
##
## This file is part of GAP, a system for computational discrete algebra.
## This file's authors include Stefan Kohl.
##
## Copyright of GAP belongs to its developers, whose names are too numerous
## to list here. Please refer to the COPYRIGHT file for details.
##
## SPDX-License-Identifier: GPL-2.0-or-later
##
## This file contains code for computing (with) continued fraction
## expansions of real numbers.
##

#############################################################################
##
#F ContinuedFractionExpansionOfRoot( , <n> )
##
InstallGlobalFunction( ContinuedFractionExpansionOfRoot,

 function ( P, n )

 local a, ai, Pi, Pi_1, pols, d, x0, step, bincoeff, i, j, k;

 if not IsUnivariatePolynomial(P)
 or not (IsPosInt(n) or DegreeOfLaurentPolynomial(P) = 2 and n = 0)
 or not ForAll(CoefficientsOfLaurentPolynomial(P)[1],IsInt)
 or LeadingCoefficient(P) < 0
 then
 Error("usage: ContinuedFractionExpansionOfRoot( , <n> ) ",
 "for a polynomial P with integer coefficients and a ",
 "positive integer <n>");
 fi;
 P := CoefficientsOfLaurentPolynomial(P);
 P := Concatenation(ListWithIdenticalEntries(P[2],0),P[1]);
 d := Length(P) - 1;
 bincoeff := List([0..d],n->List([0..d],k->Binomial(n,k)));
 if ValuePol(P,0) >= 0 then
 Error("the value of at x = 0 has to be negative");
 fi;
 a := []; Pi := ShallowCopy(P); pols := []; i := 1;
 while i <= n or n = 0 do
 if d = 2 and n = 0 then
 Add(pols,Pi);
 fi;
 x0 := 1; step := 1;
 while ValuePol(Pi,x0) < 0 do
 x0 := x0 + step;
 step := 2 * step;
 od;
 step := step/4;
 while step >= 1 do
 if ValuePol(Pi,x0) > 0 then
 x0 := x0 - step;
 else
 x0 := x0 + step;
 fi;
 step := step/2;
 od;
 if ValuePol(Pi,x0) > 0 then
 ai := x0 - 1;
 else
 ai := x0;
 fi;
 a[i] := ai;
 Pi_1 := ShallowCopy(Pi); Pi := ListWithIdenticalEntries(d+1,0);
 for j in [1..d+1] do
 for k in [1..d-j+2] do
 Pi[k] := Pi[k] + Pi_1[d-j+2]*bincoeff[d-j+2][k]*ai^(d+2-j-k);
 od;
 od;
 Pi := -Reversed(Pi);
 if Pi[d+1] = 0 then
 break; # Root is rational.
 fi;
 if d = 2 and n = 0 and Pi in pols then
 break; # One period is done.
 fi;
 i := i + 1;
 od;
 return a;
 end );

#############################################################################
##
#F ContinuedFractionApproximationOfRoot( , <n> )
##
InstallGlobalFunction( ContinuedFractionApproximationOfRoot,

 function ( P, n )

 local a, M;

 if not IsUnivariatePolynomial(P)
 or not (IsPosInt(n) or DegreeOfLaurentPolynomial(P) = 2 and n = 0)
 or not ForAll(CoefficientsOfLaurentPolynomial(P)[1],IsInt)
 or LeadingCoefficient(P) < 0
 then
 Error("usage: ContinuedFractionApproximationOfRoot( , <n> ) ",
 "for a polynomial P with integer coefficients and a ",
 "positive integer <n>");
 fi;
 a := ContinuedFractionExpansionOfRoot(P,n);
 M := Product(a,a_i->[[a_i,1],[1,0]]);
 return M[1][1]/M[2][1];
 end );

Quelle contfrac.gi Sprache: unbekannt

[ Dauer der Verarbeitung: 0.27 Sekunden (vorverarbeitet) ]