| 8.1 | Introduction | ||||
| 8.2 | Univariate Algebraic Kernel | ||||
| 8.2.1 Models | |||||
| 8.2.2 Examples | |||||
| 8.3 | Design and Implementation History | ||||
So far the package only provides models for the univariate kernel. Nevertheless, we also provide the concepts for the bivariate kernel, since this settles the interface for upcoming implementations.
The kernel provides two different function objects to construct an AlgebraicKernel_d_1::Algebraic_real_1. The most general way is to use AlgebraicKernel_d_1::Isolate_1; the function object takes a univariate polynomial and writes all real roots into a given output iterator. It is also possible to retrieve the multiplicity of each root. The second option is to construct one particular algebraic real using AlgebraicKernel_d_1::Construct_algebraic_real_1. This function object provides construction from the native int type, the coefficient type as well as the bound type. Moreover, it is possible to construct an algebraic real by giving a polynomial and either an isolating interval or the index of the root. A related function object is AlgebraicKernel_d_1::Number_of_solutions_1 computing the number of real roots of a polynomial.
An AlgebraicKernel_d_1::Algebraic_real_1 is model of RealEmbeddable, for instance, it is possible to compare two algebraic reals, to determine the sign of an algebraic real or to ask for its double approximation, see also section 4.3. Moreover, there is AlgebraicKernel_d_1::Compare_1 which provides comparison with int, the coefficient type and the bound type.
There are several ways to approximate an AlgebraicKernel_d_1::Algebraic_real_1:
AlgebraicKernel_d_1::Approximate_absolute_1 - provides an approximation that is
better than the passed absolute error bound.
AlgebraicKernel_d_1::Approximate_relative_1 - provides an approximation that is
better than the passed relative error bound.
AlgebraicKernel_d_1::Isolate_1 - returns an isolating interval with respect to
a given univariate polynomial.
A related function object is AlgebraicKernel_d_1::Bound_between_1, which computes
a number that isolates two algebraic real numbers.
Though the AlgebraicKernel_d_1 does not provide arithmetic operations on AlgebraicKernel_d_1::Algebraic_real_1, it is possible to compute the sign of a polynomial at a given algebraic real using AlgebraicKernel_d_1::Sign_at_1. Or alternatively, just compute whether the polynomial is zero at an algebraic real number using AlgebraicKernel_d_1::Is_zero_at_1. Note that this operation can be significantly less expensive, in particular, if the polynomial is not zero at the given algebraic real.
First of all the type AlgebraicKernel_d_1::Polynomial_1 is required
to be a model of the concept Polynomial_1, which is defined
in the Polynomial package (see chapter 7).
This implies that all essential functionality is provided via
CGAL::Polynomial_traits_d. However, the algebraic kernel also
provides several function objects to handle polynomials:
AlgebraicKernel_d_1::Is_square_free_1
- determines whether a polynomial is square free
AlgebraicKernel_d_1::Make_square_free_1
- computes the square free part of a polynomial
AlgebraicKernel_d_1::Square_free_factorize_1
- computes a square free factorization of a polynomial
AlgebraicKernel_d_1::Is_coprime_1
- Computes whether a pair of polynomials is square free
AlgebraicKernel_d_1::Make_coprime_1
- decompose two polynomials into the coprime factors and their common factor.
Though the polynomial package provides similar functionality we suggest to use the function objects provided by the kernel, since the design of the algebraic kernel allows for instance internal caching by the kernel.
Also note that AlgebraicKernel_d_1::Square_free_factorize_1 only computes the square free
factorization up to a constant factor. This is a slight modification with respect to its
counter part in CGAL::Polynomial_traits_d. In this way it was possible that the concepts just require
the coefficient type to be a model of IntegralDomain, instead of Field or UniqueFactorizationDomain.
For more details see also:
PolynomialTraits_d::SquareFreeFactorize
PolynomialTraits_d::SquareFreeFactorizeUpToConstantFactor
Most implementations of an AlgebraicKernel_d_1 will represent an algebraic real number by the root of a square free polynomial and an isolating interval, that is, the number is defined as the only root of the polynomial within the interval. Usually, one will refrain from computing the minimal polynomial since the computation of the minimal polynomial is much more expensive and does not pay of. However, besides the representation by a polynomial and an isolating interval one can also imagine the representation by a polynomial and the index of the root, e.g., as the ith real root when enumerated from minus to plus infinity. Moreover, it may very well be that the kernel just computes an approximation of the number, whereas the representing polynomial is not computed yet. This is in particular relevant in relation to the AlgebraicKernel_d_2, where AlgebraicKernel_d_1::Algebraic_real_1 is used to represent coordinates of solutions of bivariate systems. Hence, the design does not allow a direct access to any, seemingly obvious, members of an AlgebraicKernel_d_1::Algebraic_real_1. Instead there is, e.g., AlgebraicKernel_d_1::Compute_polynomial_1 which emphasizes that the requested polynomial may not be computed yet. Similarly, there is no way to directly ask for the refinement of the current isolating interval since this would impose a state to every object of an AlgebraicKernel_d_1::Algebraic_real_1.
The package offers two univariate algebraic kernels that are based on the library Rs [RS], namely CGAL::Algebraic_kernel_d_1_RS_Gmpz and CGAL::Algebraic_kernel_d_1_RS_Gmpq. As the names indicate, the kernels are based on the library Rs [RS] and support univariate polynomials over CGAL::Gmpz or CGAL::Gmpq, respectively.
In general we encourage to use CGAL::Algebraic_kernel_d_1_RS_Gmpz instead of CGAL::Algebraic_kernel_d_1_RS_Gmpq. This is caused by the fact that the most efficient way to compute operations (such as gcd) on polynomials with rational coefficients is to use the corresponding implementation for polynomials with integer coefficients. That is, the CGAL::Algebraic_kernel_d_1_RS_Gmpq is slightly slower due to overhead caused by the necessary conversions. However, since this may not always be a major issue the CGAL::Algebraic_kernel_d_1_RS_Gmpq is provided for convenience.
The core of both kernels is the implementation of the interval Descartes algorithm [RZ04] of the library RS [RS], which is used to isolate the roots of the polynomial. The RS library restricts its attention to univariate integer polynomials and some substantial gain of efficiency can be made by using a kernel that does not follow the generic programming paradigm, by avoiding interfaces between layers. Specifically, the fact of working with only a number type allows to optimize some polynomial operations as well as memory handling. The implementation of these kernels make heavy use of the MPFR [MPFb] and MPFI [MPFa] libraries, and of their CGAL interfaces, Gmpfr and Gmpfi. The algebraic numbers (roots of the polynomials) are represented in the two RS kernels by a Gmpfi interval and a pointer to the polynomial of which they are roots. See [LPT09] for more details on the implementation, tests of these kernels, comparisons with other algebraic kernels and discussions about the efficiency.
File: examples/Algebraic_kernel_d/Construct_algebraic_real_1.cpp
#include <CGAL/Algebraic_kernel_d_1_RS_Gmpz.h>
#include <vector>
typedef CGAL::Algebraic_kernel_d_1_RS_Gmpz AK;
typedef AK::Polynomial_1 Polynomial_1;
typedef AK::Algebraic_real_1 Algebraic_real_1;
typedef AK::Coefficient Coefficient;
typedef AK::Bound Bound;
typedef AK::Multiplicity_type Multiplicity_type;
int main(){
AK ak; // an object of Algebraic_kernel_d_1_RS_Gmpz
AK::Construct_algebraic_real_1 construct_algreal_1 = ak.construct_algebraic_real_1_object();
std::cout << "Construct from int : " << construct_algreal_1(int(2)) << "\n";
std::cout << "Construct from Coefficient : " << construct_algreal_1(Coefficient(2)) << "\n";
std::cout << "Construct from Bound : " << construct_algreal_1(Bound(2)) << "\n\n";
Polynomial_1 x = CGAL::shift(AK::Polynomial_1(1),1); // the monomial x
std::cout << "Construct by index : "
<< construct_algreal_1(x*x-2,1) << "\n"
<< to_double(construct_algreal_1(x*x-2,1)) << "\n";
std::cout << "Construct by isolating interval : "
<< construct_algreal_1(x*x-2,Bound(0),Bound(2)) << "\n"
<< to_double(construct_algreal_1(x*x-2,Bound(0),Bound(2))) << "\n\n";
return 0;
}
File: examples/Algebraic_kernel_d/Solve_1.cpp
#include <CGAL/Algebraic_kernel_d_1_RS_Gmpz.h>
#include <vector>
typedef CGAL::Algebraic_kernel_d_1_RS_Gmpz AK;
typedef AK::Polynomial_1 Polynomial_1;
typedef AK::Algebraic_real_1 Algebraic_real_1;
typedef AK::Bound Bound;
typedef AK::Multiplicity_type Multiplicity_type;
int main(){
AK ak; // an object of Algebraic_kernel_d_1_RS_Gmpz
AK::Solve_1 solve_1 = ak.solve_1_object();
Polynomial_1 x = CGAL::shift(AK::Polynomial_1(1),1); // the monomial x
// variant using a bool indicating a square free polynomial
// multiplicities are not computed
std::vector<Algebraic_real_1> roots;
solve_1(x*x-2,true, std::back_inserter(roots));
std::cout << "Number of roots is : " << roots.size() << "\n";
std::cout << "First root should be -sqrt(2): " << CGAL::to_double(roots[0]) << "\n";
std::cout << "Second root should be sqrt(2): " << CGAL::to_double(roots[1]) << "\n\n";
roots.clear();
// variant for roots in a given range of a square free polynomial
solve_1((x*x-2)*(x*x-3),true, Bound(0),Bound(10),std::back_inserter(roots));
std::cout << "Number of roots is : " << roots.size() << "\n";
std::cout << "First root should be sqrt(2): " << CGAL::to_double(roots[0]) << "\n";
std::cout << "Second root should be sqrt(3): " << CGAL::to_double(roots[1]) << "\n\n";
roots.clear();
// variant computing all roots with multiplicities
std::vector<std::pair<Algebraic_real_1,Multiplicity_type> > mroots;
solve_1((x*x-2), std::back_inserter(mroots));
std::cout << "Number of roots is : " << mroots.size() << "\n";
std::cout << "First root should be -sqrt(2): " << CGAL::to_double(mroots[0].first) << ""
<< " with multiplicity " << mroots[0].second << "\n";
std::cout << "Second root should be sqrt(2): " << CGAL::to_double(mroots[1].first) << ""
<< " with multiplicity " << mroots[1].second << "\n\n";
mroots.clear();
// variant computing roots with multiplicities for a range
solve_1((x*x-2)*(x*x-3),Bound(0),Bound(10),std::back_inserter(mroots));
std::cout << "Number of roots is : " << mroots.size() << "\n";
std::cout << "First root should be sqrt(2): " << CGAL::to_double(mroots[0].first) << ""
<< " with multiplicity " << mroots[0].second << "\n";
std::cout << "Second root should be sqrt(3): " << CGAL::to_double(mroots[1].first) << ""
<< " with multiplicity " << mroots[1].second << "\n\n";
return 0;
}
The following example illustrates the comparison of AlgebraicKernel_d_1::Algebraic_real_1 numbers:
File: examples/Algebraic_kernel_d/Compare_1.cpp
#include <CGAL/Algebraic_kernel_d_1_RS_Gmpz.h>
#include <vector>
typedef CGAL::Algebraic_kernel_d_1_RS_Gmpz AK;
typedef AK::Coefficient Coefficient;
typedef AK::Polynomial_1 Polynomial_1;
typedef AK::Algebraic_real_1 Algebraic_real_1;
typedef AK::Bound Bound;
typedef std::pair<Bound,Bound> Interval;
int main(){
AK ak;
AK::Construct_algebraic_real_1 construct_algebraic_real_1 = ak.construct_algebraic_real_1_object();
Polynomial_1 x = CGAL::shift(AK::Polynomial_1(1),1); // the monomial x
Algebraic_real_1 a = construct_algebraic_real_1(x*x-2,1); // sqrt(2)
Algebraic_real_1 b = construct_algebraic_real_1(x*x-3,1); // sqrt(3)
// Algebraic_real_1 is RealEmbeddable (just some functions:)
std::cout << "sign of a is : " << CGAL::sign(a) << "\n";
std::cout << "double approximation of a is : " << CGAL::to_double(a) << "\n";
std::cout << "double approximation of b is : " << CGAL::to_double(b) << "\n";
std::cout << "double lower bound of a : " << CGAL::to_interval(a).first << "\n";
std::cout << "double upper bound of a : " << CGAL::to_interval(a).second << "\n";
std::cout << "LessThanComparable (a<b) : " << (a<b) << "\n\n";
// use compare_1 with int, Bound, Coefficient, Algebraic_real_1
AK::Compare_1 compare_1 = ak.compare_1_object();
std::cout << " compare with an int : " << compare_1(a ,int(2)) << "\n";
std::cout << " compare with an Coefficient : " << compare_1(a ,Coefficient(2)) << "\n";
std::cout << " compare with an Bound : " << compare_1(a ,Bound(2)) << "\n";
std::cout << " compare with another Algebraic_real_1: " << compare_1(a ,b) << "\n\n";
// get a value between two roots
AK::Bound_between_1 bound_between_1 = ak.bound_between_1_object();
std::cout << " value between sqrt(2) and sqrt(3) " << bound_between_1(a,b) << "\n";
std::cout << " is larger than sqrt(2) " << compare_1(bound_between_1(a,b),a) << "\n";
std::cout << " is less than sqrt(3) " << compare_1(bound_between_1(a,b),b) << "\n\n";
// approximate with relative precision
AK::Approximate_relative_1 approx_r = ak.approximate_relative_1_object();
std::cout << " lower bound of a with at least 100 bits: "<< approx_r(a,100).first << "\n";
std::cout << " upper bound of a with at least 100 bits: "<< approx_r(a,100).second << "\n\n";
// approximate with absolute error
AK::Approximate_absolute_1 approx_a = ak.approximate_absolute_1_object();
std::cout << " lower bound of b with error less than 2^-100: "<< approx_a(b,100).first << "\n";
std::cout << " upper bound of b with error less than 2^-100: "<< approx_a(b,100).second << "\n\n";
return 0;
}
The following example illustrates the isolation of AlgebraicKernel_d_1::Algebraic_real_1 numbers:
File: examples/Algebraic_kernel_d/Isolate_1.cpp
#include <CGAL/Algebraic_kernel_d_1_RS_Gmpz.h>
#include <vector>
typedef CGAL::Algebraic_kernel_d_1_RS_Gmpz AK;
typedef AK::Polynomial_1 Polynomial_1;
typedef AK::Algebraic_real_1 Algebraic_real_1;
typedef AK::Coefficient Coefficient;
typedef AK::Bound Bound;
typedef AK::Multiplicity_type Multiplicity_type;
int main(){
AK ak; // an object of Algebraic_kernel_d_1_RS_Gmpz
AK::Construct_algebraic_real_1 construct_algreal_1 = ak.construct_algebraic_real_1_object();
AK::Isolate_1 isolate_1 = ak.isolate_1_object();
AK::Compute_polynomial_1 compute_polynomial_1 = ak.compute_polynomial_1_object();
// construct an algebraic number from an integer
Algebraic_real_1 frominteger=construct_algreal_1(int(2));
std::cout << "Construct from int: " << frominteger << "\n";
// the constructed algebraic number is root of a polynomial
Polynomial_1 pol=compute_polynomial_1(frominteger);
std::cout << "The constructed number is root of: " << pol << "\n";
// construct an algebraic number from a polynomial and an isolating interval
Polynomial_1 x = CGAL::shift(AK::Polynomial_1(1),1); // the monomial x
Algebraic_real_1 frominterval=construct_algreal_1(x*x-2,Bound(0),Bound(2));
std::cout << "Construct from isolating interval: " << frominterval << "\n";
// isolate the second algebraic number from the first: this is to say,
// isolating the second algebraic number with respect to the polynomial
// of which the first constructed number is root
std::pair<Bound,Bound> isolation1 = isolate_1(frominterval,pol);
std::cout << "Isolating the second algebraic number gives: ["
<< isolation1.first << "," << isolation1.second << "]\n";
// isolate again the same algebraic number, this time with respect to
// the polynomial 10*x-14 (which has root 1.4, close to this algebraic
// number)
std::pair<Bound,Bound> isolation2 = isolate_1(frominterval,10*x-14);
std::cout << "Isolating again the second algebraic number gives: ["
<< isolation2.first << "," << isolation2.second << "]\n";
return 0;
}
The following example illustrates the sign evaluation of AlgebraicKernel_d_1::Algebraic_real_1 numbers in polynomials:
File: examples/Algebraic_kernel_d/Sign_at_1.cpp
#include <CGAL/Algebraic_kernel_d_1_RS_Gmpz.h>
#include <vector>
typedef CGAL::Algebraic_kernel_d_1_RS_Gmpz AK;
typedef AK::Polynomial_1 Polynomial_1;
typedef AK::Algebraic_real_1 Algebraic_real_1;
typedef AK::Coefficient Coefficient;
typedef AK::Bound Bound;
typedef AK::Multiplicity_type Multiplicity_type;
int main(){
AK ak;
AK::Construct_algebraic_real_1 construct_algreal_1 = ak.construct_algebraic_real_1_object();
AK::Solve_1 solve_1 = ak.solve_1_object();
AK::Sign_at_1 sign_at_1 = ak.sign_at_1_object();
AK::Is_zero_at_1 is_zero_at_1 = ak.is_zero_at_1_object();
// construct the polynomials p=x^2-5 and q=x-2
Polynomial_1 x = CGAL::shift(AK::Polynomial_1(1),1); // the monomial x
Polynomial_1 p = x*x-5;
std::cout << "Polynomial p: " << p << "\n";
Polynomial_1 q = x-2;
std::cout << "Polynomial q: " << q << "\n";
// find the roots of p (it has two roots) and q (one root)
std::vector<Algebraic_real_1> roots_p,roots_q;
solve_1(p,true, std::back_inserter(roots_p));
solve_1(q,true, std::back_inserter(roots_q));
// evaluate the second root of p in q
std::cout << "Sign of the evaluation of root 2 of p in q: "
<< sign_at_1(q,roots_p[1]) << "\n";
// evaluate the root of q in p
std::cout << "Sign of the evaluation of root 1 of q in p: "
<< sign_at_1(p,roots_q[0]) << "\n";
// check whether the evaluation of the first root of p in p is zero
std::cout << "Is zero the evaluation of root 1 of p in p? "
<< is_zero_at_1(p,roots_p[0]) << "\n";
return 0;
}
This package is clearly split into a univariate and bivariate kernel. However, with respect to its history the package splits into a design part and an implementation part. The concepts, which make up the design part, where written by Eric Berberich, Michael Hemmer, and Monique Teillaud. The implementation part is currently comprised of two univariate kernels that interface the library RS [RS]. These models were written by Luis Peñaranda and Sylvain Lazard.
The design history of the package is fairly old and several ideas that influenced this package can already be found in [BHKT07]. Since then, the initial design underwent considerable changes. For instance, it was decided that the algebraic numbers should be under the control of the algebraic kernel. On the other hand the initial support for polynomials was extended to a separate and independent package that is not restricted to a certain number of variables. Thus, the authors want to thank for all the useful feedback and ideas that was brought to them throughout the last years. In particular, they want to thank Menelaos Karavelas and Elias Tsigaridas for their initial contributions.
The implementation history shall be considered as unclosed, since there are some more models in the pipeline. So far, the package provides two models of a univariate algebraic kernel. Both models interface the library RS [RS] by Fabrice Rouillier. The authors want to thank Fabrice Rouillier and Elias Tsigaridas for strong support and many useful discussions that lead to the integration of RS.