From ea3d9da0efb9ba526975b6a924b716b4231f356f Mon Sep 17 00:00:00 2001 From: Sebastian Limbach Date: Tue, 6 Mar 2007 13:30:58 +0000 Subject: [PATCH] Missing functors added according to documentation --- Polynomial/include/CGAL/Polynomial_traits_d.h | 1358 ++++++++++------- 1 file changed, 791 insertions(+), 567 deletions(-) diff --git a/Polynomial/include/CGAL/Polynomial_traits_d.h b/Polynomial/include/CGAL/Polynomial_traits_d.h index 9180439d86d..bbfdd18137e 100644 --- a/Polynomial/include/CGAL/Polynomial_traits_d.h +++ b/Polynomial/include/CGAL/Polynomial_traits_d.h @@ -43,6 +43,7 @@ class Polynomial_traits_d_base { }; typedef Null_functor Construct_polynomial; + typedef Null_functor Get_coefficient; typedef Null_functor Degree; typedef Null_functor Total_degree; typedef Null_functor Leading_coefficient; @@ -87,10 +88,9 @@ class Polynomial_traits_d_base { }; typedef Null_functor Square_free_factorization_up_to_constant_factor; - typedef Null_functor Evaluate; - typedef Null_functor Evaluate_homogeneous; typedef Null_functor Resultant; typedef Null_functor Canonicalize; + typedef Null_functor Evaluate_homogeneous; struct Innermost_leading_coefficient :public Unary_function { @@ -105,7 +105,79 @@ class Polynomial_traits_d_base { return std::vector(); } }; + + struct Get_innermost_coefficient + : public Binary_function< ICoeff, Polynomial_d, Exponent_vector > { + + ICoeff operator()( const Polynomial_d& p, Exponent_vector ev ) { + CGAL_precondition( ev.empty() ); + return p; + } + + }; + + struct Evaluate { + template< class Input_iterator > + ICoeff operator()( const Polynomial_d& p, Input_iterator, Input_iterator ) { + //std::cerr << p << std::endl; + return p; + } + }; + }; + template< class Polynomial, int d = CGAL::Polynomial_traits_d< Polynomial>::d > + struct Evaluate_homogeneous; + + template< class Polynomial > + struct Evaluate_homogeneous< Polynomial, 1 > { + typedef typename CGAL::Polynomial_traits_d< Polynomial > PT; + typedef typename PT::Coefficient Coefficient; + typedef typename PT::Innermost_coefficient ICoeff; + typedef typename CGAL::Polynomial_traits_d< Coefficient > PTC; + + template< class Input_iterator > + ICoeff operator()( const Polynomial& p, + Input_iterator begin, + Input_iterator end, + int total_degree, + const ICoeff& v ) const { + --end; + CGAL_precondition( begin == end ); +/* std::cerr << (*end) << ", " << v << ", " << total_degree << std::endl; + std::cerr << p << std::endl;*/ + return p.evaluate_homogeneous( (*end), v, total_degree ); + } + + }; + + template< class Polynomial, int d > + struct Evaluate_homogeneous { + typedef typename CGAL::Polynomial_traits_d< Polynomial > PT; + typedef typename PT::Coefficient Coefficient; + typedef typename PT::Innermost_coefficient ICoeff; + typedef typename CGAL::Polynomial_traits_d< Coefficient > PTC; + + template< class Input_iterator > + ICoeff operator()( const Polynomial& p, + Input_iterator begin, + Input_iterator end, + int total_degree, + const ICoeff& v ) const { + CGAL_precondition( begin != end ); + //typename PT::Evaluate evaluate; + typename PT::Degree degree; + typename CGALi::Evaluate_homogeneous< Coefficient > eval_hom; + --end; + + std::vector< ICoeff > cv; + + for( int i = 0; i <= degree(p); ++i ) { + cv.push_back( eval_hom( p[i], begin, end, total_degree - i, v ) ); + } + + return (CGAL::Polynomial< ICoeff >( cv.begin(), cv.end() )).evaluate((*end)); + } + }; } // namespace CGALi @@ -131,6 +203,13 @@ private: typedef std::vector< Exponents_coeff_pair > Monom_rep; public: + // + // Functors as defined in the reference manual (with sometimes slightly + // extended functionality) + // + + + // Construct_polynomial; struct Construct_polynomial { typedef Polynomial_d result_type; @@ -297,6 +376,715 @@ public: }; }; + // Get_coefficient; + struct Get_coefficient + : public Binary_function { + + Coefficient operator()( const Polynomial_d& p, int i) const { + CGAL_precondition( i >= 0 ); + typename PT::Degree degree; + if( i > degree(p) ) + return Coefficient(0); + return p[i]; + } + + }; + + // Get_innermost_coefficient; + struct Get_innermost_coefficient + : public Binary_function< Polynomial_d, Exponent_vector, Innermost_coefficient > { + + Innermost_coefficient operator()( const Polynomial_d& p, Exponent_vector ev ) const { + CGAL_precondition( !ev.empty() ); + typename PTC::Get_innermost_coefficient gic; + typename PT::Get_coefficient gc; + int exponent = ev.back(); + ev.pop_back(); + return gic( gc( p, exponent ), ev ); + }; + + }; + + // Swap; + // Swap variable x_i with x_j + struct Swap { + typedef Polynomial_d result_type; + typedef Polynomial_d first_argument_type; + typedef int second_argument_type; + typedef int third_argument_type; + typedef Arity_tag< 3 > Arity; + + Polynomial_d operator()(const Polynomial_d& p, int i, int j ) const { + //std::cout << i <<" " << j << " : " ; + CGAL_precondition(0 <= i && i < d); + CGAL_precondition(0 <= j && j < d); + typedef std::pair< Exponent_vector, Innermost_coefficient > + Exponents_coeff_pair; + typedef std::vector< Exponents_coeff_pair > Monom_rep; + Get_monom_representation gmr; + Construct_polynomial construct; + Monom_rep mon_rep; + gmr( p, std::back_inserter( mon_rep ) ); + for( typename Monom_rep::iterator it = mon_rep.begin(); + it != mon_rep.end(); + ++it ) { + std::swap(it->first[i],it->first[j]); + // it->first.swap( i, j ); + } + std::sort( mon_rep.begin(), mon_rep.end() ); + return construct( mon_rep.begin(), mon_rep.end() ); + } + }; + + // Move; + // move variable x_i to position of x_j + // order of other variables remains + // default j = d makes x_i the othermost variable + struct Move { + typedef Polynomial_d result_type; + typedef Polynomial_d first_argument_type; + typedef int second_argument_type; + typedef int third_argument_type; + typedef Arity_tag< 3 > Arity; + + Polynomial_d operator()(const Polynomial_d& p, int i, int j = (d-1) ) const { + //std::cout << x <<" " << y << " : " ; + CGAL_precondition(0 <= i && i < d); + CGAL_precondition(0 <= j && j < d); + typedef std::pair< Exponent_vector, Innermost_coefficient > + Exponents_coeff_pair; + typedef std::vector< Exponents_coeff_pair > Monom_rep; + Get_monom_representation gmr; + Construct_polynomial construct; + Monom_rep mon_rep; + gmr( p, std::back_inserter( mon_rep ) ); + for( typename Monom_rep::iterator it = mon_rep.begin(); + it != mon_rep.end(); + ++it ) { + // this is as good as std::rotate since it uses swap also + if (i < j) + for( int k = i; k < j; k++ ) + std::swap(it->first[k],it->first[k+1]); + else + for( int k = i; k > j; k-- ) + std::swap(it->first[k],it->first[k-1]); + + } + std::sort( mon_rep.begin(), mon_rep.end() ); + return construct( mon_rep.begin(), mon_rep.end() ); + } + }; + + // Degree; + struct Degree : public Unary_function< Polynomial_d , int >{ + int operator()(const Polynomial_d& p, int i = (d-1)) const { + if (i == (d-1)) return p.degree(); + else return Swap()(p,i,d-1).degree(); + } + }; + + // Total_degree; + struct Total_degree : public Unary_function< Polynomial_d , int >{ + int operator()(const Polynomial_d& p) const { + typedef Polynomial_traits_d COEFF_POLY_TRAITS; + typename COEFF_POLY_TRAITS::Total_degree total_degree; + Degree degree; + CGAL_precondition( degree(p) >= 0); + + int result = 0; + for(int i = 0; i <= degree(p) ; i++){ + if( ! CGAL::is_zero( p[i]) ) + result = std::max(result , total_degree(p[i]) + i ); + } + return result; + } + }; + + // Leading_coefficient; + struct Leading_coefficient + : public Unary_function< Polynomial_d , Coefficient>{ + Coefficient operator()(const Polynomial_d& p) const { + return p.lcoeff(); + } + Coefficient operator()(Polynomial_d p, int i) const { + return Swap()(p,i,PT::d-1).lcoeff(); + } + }; + + // Innermost_leading_coefficient; + struct Innermost_leading_coefficient + : public Unary_function< Polynomial_d , Innermost_coefficient>{ + Innermost_coefficient + operator()(const Polynomial_d& p) const { + typename PTC::Innermost_leading_coefficient ilcoeff; + typename PT::Leading_coefficient lcoeff; + return ilcoeff(lcoeff(p)); + } + }; + + // Canonicalize; + struct Canonicalize + : public Unary_function{ + Polynomial_d + operator()( const Polynomial_d& p ) const { + return CGAL::INTERN_POLYNOMIAL::canonicalize_polynomial(p); + } + }; + + // Differentiate; + struct Differentiate + : public Unary_function{ + Polynomial_d + operator()(Polynomial_d p, int i = (d-1)) const { + if (i == (d-1) ){ + p.diff(); + }else{ + Swap swap; + p = swap(p,i,d-1); + p.diff(); + p = swap(p,i,d-1); + } + return p; + } + }; + + // Evaluate; + struct Evaluate + :public Unary_function{ + Coefficient + operator()(const Polynomial_d& p, Innermost_coefficient x, int i = (d-1)) + const { + if(i == (d-1) ) + return p.evaluate(x); + else{ + return Move()(p,i).evaluate(x); + } + } + + template< class Input_iterator > + Innermost_coefficient operator()( const Polynomial_d& p, Input_iterator begin, Input_iterator end ) const { + CGAL_precondition( begin != end ); + + typename PT::Evaluate evaluatePoly; + typename PTC::Evaluate evaluateCoeff; + --end; + return evaluateCoeff( evaluatePoly( p, (*end) ), begin, end ); + } + }; + + // Evaluate_homogeneous; + struct Evaluate_homogeneous{ + typedef Coefficient result_type; + typedef Polynomial_d first_argument_type; + typedef Innermost_coefficient second_argument_type; + typedef Innermost_coefficient third_argument_type; + typedef Arity_tag< 3 > Arity; + + Coefficient + operator()( + const Polynomial_d& p, + Innermost_coefficient a, + Innermost_coefficient b, + int i = (PT::d-1) ) const { + if (i == (d-1) ) + return p.evaluate_homogeneous(a,b); + else + return Move()(p,i,PT::d-1).evaluate_homogeneous(a,b); + } + + template< class Input_iterator > + Innermost_coefficient operator()( const Polynomial_d & p, + Input_iterator begin, + Input_iterator end ) const { + typename PT::Total_degree total_degree; + typename PT::Evaluate_homogeneous eval_hom; + return eval_hom( p, begin, end, total_degree(p) ); + } + + template< class Input_iterator > + Innermost_coefficient operator()( const Polynomial_d& p, + Input_iterator begin, + Input_iterator end, + int total_degree ) const { + --end; + typename CGALi::Evaluate_homogeneous< Polynomial_d > eval_hom; + return eval_hom( p, begin, end, total_degree, (*end) ); + } + + }; + + // Is_zero_at; + struct Is_zero_at { + typedef bool result_type; + + template< class Input_iterator > + bool operator()( const Polynomial_d& p, Input_iterator begin, Input_iterator end ) const { + typename PT::Evaluate evaluate; + return( CGAL::is_zero( evaluate( p, begin, end ) ) ); + } + }; + + // Is_zero_at_homogeneous; + struct Is_zero_at_homogeneous { + typedef bool result_type; + + template< class Input_iterator > + bool operator()( const Polynomial_d& p, Input_iterator begin, Input_iterator end ) const { + typename PT::Evaluate_homogeneous evaluate_homogeneous; + return( CGAL::is_zero( evaluate_homogeneous( p, begin, end ) ) ); + } + }; + + // Sign_at; + struct Sign_at { + typedef Sign result_type; + + template< class Input_iterator > + Sign operator()( const Polynomial_d& p, Input_iterator begin, Input_iterator end ) const { + typename PT::Evaluate evaluate; + return CGAL::sign( evaluate( p, begin, end ) ); + } + }; + + // Sign_at_homogeneous; + struct Sign_at_homogeneous { + typedef Sign result_type; + + template< class Input_iterator > + Sign operator()( const Polynomial_d& p, Input_iterator begin, Input_iterator end ) const { + typename PT::Evaluate_homogeneous evaluate_homogeneous; + return CGAL::sign( evaluate_homogeneous( p, begin, end ) ); + } + }; + + // Compare; + struct Compare + : public Binary_function< Comparison_result, Polynomial_d, Polynomial_d > { + + Comparison_result operator()( const Polynomial_d& p1, const Polynomial_d& p2 ) const { + return p1.compare( p2 ); + } + }; + + // + // This is going to be in PolynomialToolBox + // + + // Univariate_content; + struct Univariate_content + : public Unary_function< Polynomial_d , Coefficient>{ + Coefficient operator()(const Polynomial_d& p) const { + return p.content(); + } + Coefficient operator()(Polynomial_d p, int i) const { + return Swap()(p,i,PT::d-1).content(); + } + }; + + // Multivariate_content; + struct Multivariate_content + : public Unary_function< Polynomial_d , Innermost_coefficient >{ + Innermost_coefficient + operator()(const Polynomial_d& p) const { + typedef Innermost_coefficient_iterator IT; + Innermost_coefficient content(0); + for (IT it = Innermost_coefficient_begin()(p); + it != Innermost_coefficient_end()(p); + it++){ + content = CGAL::gcd(content, *it); + if(CGAL::is_one(content)) break; + } + return content; + } + }; + + // Square_free_factorization; + struct Square_free_factorization{ + typedef int result_type; + + private: + typedef Coefficient Coeff; + typedef Innermost_coefficient ICoeff; + + // rsqff computes the sqff recursively for Coeff + // end of recursion: ICoeff + template < class OutputIterator1, class OutputIterator2 > + int rsqff (ICoeff c, + OutputIterator1 factors, + OutputIterator2 mults) const{ + return 0; + } + template < class OutputIterator1, class OutputIterator2 > + int rsqff ( + typename First_if_different::Type c, + OutputIterator1 fit, + OutputIterator2 mit) const { + typename PTC::Square_free_factorization sqff; + std::vector factors; + int n = sqff(c, std::back_inserter(factors), mit); + for(int i = 0; i < (int)factors.size(); i++){ + *fit++=Polynomial_d(factors[i]); + } + return n; + } + public: + template < class OutputIterator1, class OutputIterator2 > + int operator()( + const Polynomial_d& p, + OutputIterator1 fit, + OutputIterator2 mit) const { + Coefficient c; + int n = square_free_factorization(p,fit,mit,c); + if (Total_degree()(c) > 0) + return rsqff(c,fit,mit)+n; + else + return n; + } + }; + + // Make_square_free; + struct Make_square_free + : public Unary_function< Polynomial_d, Polynomial_d >{ + Polynomial_d + operator()(const Polynomial_d& p) const { + if (CGAL::is_zero(p)) return p; + Univariate_content_up_to_constant_factor ucontent_utcf; + Integral_division_up_to_constant_factor idiv_utcf; + Differentiate diff; + typename PTC::Make_square_free msf; + + Coefficient content = ucontent_utcf(p); + Polynomial_d result = Polynomial_d(msf(content)); + + Polynomial_d regular_part = idiv_utcf(p,Polynomial_d(content)); + Polynomial_d g = gcd_utcf(regular_part,diff(regular_part)); + + result *= idiv_utcf(regular_part,g); + return Canonicalize()(result); + +#if 0 + + Square_free_factorization sqff_utcf; + std::vector facs; + std::vector mults; + sqff_utcf(p,std::back_inserter(facs), std::back_inserter(mults)); + for(typename std::vector::iterator it = facs.begin(); + it != facs.end(); + it++){ + result *= *it; + } + return result; +#endif + + } + }; + + // Pseudo_division; + struct Pseudo_division { + typedef Polynomial_d result_type; + void + operator()( + const Polynomial_d& f, const Polynomial_d& g, + Polynomial_d& q, Polynomial_d& r, Coefficient& D) const { + Polynomial_d::pseudo_division(f,g,q,r,D); + } + }; + + // Pseudo_division_quotient; + struct Pseudo_division_quotient + :public Binary_function { + + Polynomial_d + operator()(const Polynomial_d& f, const Polynomial_d& g) const { + Polynomial_d q,r; + Coefficient D; + Polynomial_d::pseudo_division(f,g,q,r,D); + return q; + } + }; + + // Pseudo_division_remainder; + struct Pseudo_division_remainder + :public Binary_function { + + Polynomial_d + operator()(const Polynomial_d& f, const Polynomial_d& g) const { + Polynomial_d q,r; + Coefficient D; + Polynomial_d::pseudo_division(f,g,q,r,D); + return r; + } + }; + + // Gcd_up_to_constant_factor; + struct Gcd_up_to_constant_factor + :public Binary_function { + Polynomial_d + operator()(const Polynomial_d& p, const Polynomial_d& q) const { + if (CGAL::is_zero(p) && CGAL::is_zero(q)) + return Polynomial_d(0); + return gcd_utcf(p,q); + } + }; + + // Integral_division_up_to_constant_factor; + struct Integral_division_up_to_constant_factor + :public Binary_function { + Polynomial_d + operator()(const Polynomial_d& p, const Polynomial_d& q) const { + typedef Innermost_coefficient IC; + + typename PT::Construct_polynomial construct; + typename PT::Innermost_leading_coefficient ilcoeff; + typename PT::Innermost_coefficient_begin begin; + typename PT::Innermost_coefficient_end end; + typedef Algebraic_extension_traits AET; + typename AET::Denominator_for_algebraic_integers dfai; + typename AET::Normalization_factor nfac; + + + IC ilcoeff_q = ilcoeff(q); + // this factor is needed in case IC is an Algebraic extension + IC dfai_q = dfai(begin(q), end(q)); + // make dfai_q a 'scalar' + ilcoeff_q *= dfai_q * nfac(dfai_q); + Polynomial_d result = (p * construct(ilcoeff_q)) / q; + + return Canonicalize()(result); + } + }; + + // Univariate_content_up_to_constant_factor; + struct Univariate_content_up_to_constant_factor + :public Unary_function { + Coefficient + operator()(const Polynomial_d& p) const { + typename PTC::Gcd_up_to_constant_factor gcd_utcf; + + if(CGAL::is_zero(p)) return Coefficient(0); + if(PT::d == 1) return Coefficient(1); + + Coefficient result(0); + for(typename Polynomial_d::const_iterator it = p.begin(); + it != p.end(); + it++){ + result = gcd_utcf(*it,result); + } + return result; + + } + }; + + // Square_free_factorization_up_to_constant_factor; + struct Square_free_factorization_up_to_constant_factor { + typedef int result_type; + private: + typedef Coefficient Coeff; + typedef Innermost_coefficient ICoeff; + + // rsqff_utcf computes the sqff recursively for Coeff + // end of recursion: ICoeff + template < class OutputIterator1, class OutputIterator2 > + int rsqff_utcf (ICoeff c, + OutputIterator1 factors, + OutputIterator2 mults) const{ + return 0; + } + template < class OutputIterator1, class OutputIterator2 > + int rsqff_utcf ( + typename First_if_different::Type c, + OutputIterator1 fit, + OutputIterator2 mit) const { + typename PTC::Square_free_factorization sqff; + std::vector factors; + int n = sqff(c, std::back_inserter(factors), mit); + for(int i = 0; i < (int)factors.size(); i++){ + *fit++=Polynomial_d(factors[i]); + } + return n; + } + public: + template < class OutputIterator1, class OutputIterator2 > + int operator()( + Polynomial_d p, + OutputIterator1 fit, + OutputIterator2 mit) const { + + if (CGAL::is_zero(p)) return 0; + + Univariate_content_up_to_constant_factor ucontent_utcf; + Integral_division_up_to_constant_factor idiv_utcf; + Coefficient c = ucontent_utcf(p); + p = idiv_utcf( p , Polynomial_d(c)); + int n = square_free_factorization_utcf(p,fit,mit); + if (Total_degree()(c) > 0) + return rsqff_utcf(c,fit,mit)+n; + else + return n; + } + }; + + // Shift; + struct Shift + : public Unary_function< Polynomial_d, Polynomial_d >{ + Polynomial_d + operator()(const Polynomial_d& p, int e, int i = PT::d) const { + Construct_polynomial construct; + Get_monom_representation gmr; + Monom_rep monom_rep; + gmr(p,std::back_inserter(monom_rep)); + for(typename Monom_rep::iterator it = monom_rep.begin(); + it != monom_rep.end(); + it++){ + it->first[i-1]+=e; + } + return construct(monom_rep.begin(), monom_rep.end()); + } + }; + + // Negate; + struct Negate + : public Unary_function< Polynomial_d, Polynomial_d >{ + + Polynomial_d operator()(const Polynomial_d& p, int i = (d-1)) const { + Construct_polynomial construct; + Get_monom_representation gmr; + Monom_rep monom_rep; + gmr(p,std::back_inserter(monom_rep)); + for(typename Monom_rep::iterator it = monom_rep.begin(); + it != monom_rep.end(); + it++){ + if (it->first[i] % 2 != 0) + it->second = - it->second; + } + return construct(monom_rep.begin(), monom_rep.end()); + } + }; + + // Invert; + struct Invert + : public Unary_function< Polynomial_d , Polynomial_d >{ + Polynomial_d operator()(Polynomial_d p, int i = (PT::d-1)) const { + if (i == (d-1)){ + p.reversal(); + }else{ + p = Swap()(p,i,PT::d-1); + p.reversal(); + p = Swap()(p,i,PT::d-1); + } + return p ; + } + }; + + // Translate; + struct Translate + : public Binary_function< Polynomial_d , Polynomial_d, + Innermost_coefficient >{ + Polynomial_d + operator()( + Polynomial_d p, + const Innermost_coefficient& c, + int i = (d-1)) + const { + if (i == (d-1) ){ + p.translate(Coefficient(c)); + }else{ + Swap swap; + p = swap(p,i,d-1); + p.translate(Coefficient(c)); + p = swap(p,i,d-1); + } + return p; + } + }; + + // Translate_homogeneous; + struct Translate_homogeneous{ + typedef Polynomial_d result_type; + typedef Polynomial_d first_argument_type; + typedef Innermost_coefficient second_argument_type; + typedef Innermost_coefficient third_argument_type; + + Polynomial_d + operator()(Polynomial_d p, + const Innermost_coefficient& a, + const Innermost_coefficient& b, + int i = (d-1) ) const { + if (i == (d-1) ){ + p.translate(Coefficient(a), Coefficient(b) ); + }else{ + Swap swap; + p = swap(p,i,d-1); + p.translate(Coefficient(a), Coefficient(b)); + p = swap(p,i,d-1); + } + return p; + } + }; + + // Scale; + struct Scale + : public Binary_function< Polynomial_d, Innermost_coefficient, Polynomial_d > { + + Polynomial_d operator()( Polynomial_d p, const Innermost_coefficient& c, + int i = (PT::d-1) ) { + typename PT::Scale_homogeneous scale_homogeneous; + + return scale_homogeneous( p, c, Innermost_coefficient(1), i ); + } + + }; + + // Scale_homogeneous; + struct Scale_homogeneous{ + typedef Polynomial_d result_type; + typedef Polynomial_d first_argument_type; + typedef Innermost_coefficient second_argument_type; + typedef Innermost_coefficient third_argument_type; + + Polynomial_d + operator()( + Polynomial_d p, + const Innermost_coefficient& a, + const Innermost_coefficient& b, + int i = (d-1) ) const { + CGAL_precondition( ! CGAL::is_zero(b) ); + + if (i == (d-1) ) p = Swap()(p,i,d-1); + + if(CGAL::is_one(b)) + p.scale_up(Coefficient(a)); + else + if(CGAL::is_one(a)) + p.scale_down(Coefficient(b)); + else + p.scale(Coefficient(a), Coefficient(b) ); + + if (i == (d-1) ) p = Swap()(p,i,d-1); + + return p; + } + }; + + // Resultant; + struct Resultant + : public Binary_function{ + + Coefficient + operator()( + const Polynomial_d& p, + const Polynomial_d& q, + int i = (d-1) ) const { + if(i == (d-1) ) + return prs_resultant(p,q); + else + return prs_resultant(Move()(p,i),Move()(q,i)); + } + }; + + // + // Functors not mentioned in the reference manual + // struct Get_monom_representation { typedef std::pair< Exponent_vector, Innermost_coefficient > @@ -380,116 +1168,8 @@ public: } }; - // Swap variable x_i with x_j - struct Swap { - typedef Polynomial_d result_type; - typedef Polynomial_d first_argument_type; - typedef int second_argument_type; - typedef int third_argument_type; - typedef Arity_tag< 3 > Arity; - - Polynomial_d operator()(const Polynomial_d& p, int i, int j ) const { - //std::cout << i <<" " << j << " : " ; - CGAL_precondition(0 <= i && i < d); - CGAL_precondition(0 <= j && j < d); - typedef std::pair< Exponent_vector, Innermost_coefficient > - Exponents_coeff_pair; - typedef std::vector< Exponents_coeff_pair > Monom_rep; - Get_monom_representation gmr; - Construct_polynomial construct; - Monom_rep mon_rep; - gmr( p, std::back_inserter( mon_rep ) ); - for( typename Monom_rep::iterator it = mon_rep.begin(); - it != mon_rep.end(); - ++it ) { - std::swap(it->first[i],it->first[j]); - // it->first.swap( i, j ); - } - std::sort( mon_rep.begin(), mon_rep.end() ); - return construct( mon_rep.begin(), mon_rep.end() ); - } - }; - - // move variable x_i to position of x_j - // order of other variables remains - // default j = d makes x_i the othermost variable - struct Move { - typedef Polynomial_d result_type; - typedef Polynomial_d first_argument_type; - typedef int second_argument_type; - typedef int third_argument_type; - typedef Arity_tag< 3 > Arity; - - Polynomial_d operator()(const Polynomial_d& p, int i, int j = (d-1) ) const { - //std::cout << x <<" " << y << " : " ; - CGAL_precondition(0 <= i && i < d); - CGAL_precondition(0 <= j && j < d); - typedef std::pair< Exponent_vector, Innermost_coefficient > - Exponents_coeff_pair; - typedef std::vector< Exponents_coeff_pair > Monom_rep; - Get_monom_representation gmr; - Construct_polynomial construct; - Monom_rep mon_rep; - gmr( p, std::back_inserter( mon_rep ) ); - for( typename Monom_rep::iterator it = mon_rep.begin(); - it != mon_rep.end(); - ++it ) { - // this is as good as std::rotate since it uses swap also - if (i < j) - for( int k = i; k < j; k++ ) - std::swap(it->first[k],it->first[k+1]); - else - for( int k = i; k > j; k-- ) - std::swap(it->first[k],it->first[k-1]); - - } - std::sort( mon_rep.begin(), mon_rep.end() ); - return construct( mon_rep.begin(), mon_rep.end() ); - } - }; - struct Degree : public Unary_function< Polynomial_d , int >{ - int operator()(const Polynomial_d& p, int i = (d-1)) const { - if (i == (d-1)) return p.degree(); - else return Swap()(p,i,d-1).degree(); - } - }; - // Total_degree; - struct Total_degree : public Unary_function< Polynomial_d , int >{ - int operator()(const Polynomial_d& p) const { - typedef Polynomial_traits_d COEFF_POLY_TRAITS; - typename COEFF_POLY_TRAITS::Total_degree total_degree; - Degree degree; - CGAL_precondition( degree(p) >= 0); - - int result = 0; - for(int i = 0; i <= degree(p) ; i++){ - if( ! CGAL::is_zero( p[i]) ) - result = std::max(result , total_degree(p[i]) + i ); - } - return result; - } - }; -// Leading_coefficient; - struct Leading_coefficient - : public Unary_function< Polynomial_d , Coefficient>{ - Coefficient operator()(const Polynomial_d& p) const { - return p.lcoeff(); - } - Coefficient operator()(Polynomial_d p, int i) const { - return Swap()(p,i,PT::d-1).lcoeff(); - } - }; - struct Innermost_leading_coefficient - : public Unary_function< Polynomial_d , Innermost_coefficient>{ - Innermost_coefficient - operator()(const Polynomial_d& p) const { - typename PTC::Innermost_leading_coefficient ilcoeff; - typename PT::Leading_coefficient lcoeff; - return ilcoeff(lcoeff(p)); - } - }; // returns the Exponten_vector of the innermost leading coefficient // TODO use Exponent vector // TODO document @@ -507,463 +1187,7 @@ public: return result; } }; - -// Univariate_content; - struct Univariate_content - : public Unary_function< Polynomial_d , Coefficient>{ - Coefficient operator()(const Polynomial_d& p) const { - return p.content(); - } - Coefficient operator()(Polynomial_d p, int i) const { - return Swap()(p,i,PT::d-1).content(); - } - }; -// Multivariate_content; - struct Multivariate_content - : public Unary_function< Polynomial_d , Innermost_coefficient >{ - Innermost_coefficient - operator()(const Polynomial_d& p) const { - typedef Innermost_coefficient_iterator IT; - Innermost_coefficient content(0); - for (IT it = Innermost_coefficient_begin()(p); - it != Innermost_coefficient_end()(p); - it++){ - content = CGAL::gcd(content, *it); - if(CGAL::is_one(content)) break; - } - return content; - } - }; -// Shift; - struct Shift - : public Unary_function< Polynomial_d, Polynomial_d >{ - Polynomial_d - operator()(const Polynomial_d& p, int e, int i = PT::d) const { - Construct_polynomial construct; - Get_monom_representation gmr; - Monom_rep monom_rep; - gmr(p,std::back_inserter(monom_rep)); - for(typename Monom_rep::iterator it = monom_rep.begin(); - it != monom_rep.end(); - it++){ - it->first[i-1]+=e; - } - return construct(monom_rep.begin(), monom_rep.end()); - } - }; -// Negate; - struct Negate - : public Unary_function< Polynomial_d, Polynomial_d >{ - - Polynomial_d operator()(const Polynomial_d& p, int i = (d-1)) const { - Construct_polynomial construct; - Get_monom_representation gmr; - Monom_rep monom_rep; - gmr(p,std::back_inserter(monom_rep)); - for(typename Monom_rep::iterator it = monom_rep.begin(); - it != monom_rep.end(); - it++){ - if (it->first[i] % 2 != 0) - it->second = - it->second; - } - return construct(monom_rep.begin(), monom_rep.end()); - } - }; -// Invert; - struct Invert - : public Unary_function< Polynomial_d , Polynomial_d >{ - Polynomial_d operator()(Polynomial_d p, int i = (PT::d-1)) const { - if (i == (d-1)){ - p.reversal(); - }else{ - p = Swap()(p,i,PT::d-1); - p.reversal(); - p = Swap()(p,i,PT::d-1); - } - return p ; - } - }; -// Translate; - struct Translate - : public Binary_function< Polynomial_d , Polynomial_d, - Innermost_coefficient >{ - Polynomial_d - operator()( - Polynomial_d p, - const Innermost_coefficient& c, - int i = (d-1)) - const { - if (i == (d-1) ){ - p.translate(Coefficient(c)); - }else{ - Swap swap; - p = swap(p,i,d-1); - p.translate(Coefficient(c)); - p = swap(p,i,d-1); - } - return p; - } - }; - -// Translate_homogeneous; - struct Translate_homogeneous{ - typedef Polynomial_d result_type; - typedef Polynomial_d first_argument_type; - typedef Innermost_coefficient second_argument_type; - typedef Innermost_coefficient third_argument_type; - - Polynomial_d - operator()(Polynomial_d p, - const Innermost_coefficient& a, - const Innermost_coefficient& b, - int i = (d-1) ) const { - if (i == (d-1) ){ - p.translate(Coefficient(a), Coefficient(b) ); - }else{ - Swap swap; - p = swap(p,i,d-1); - p.translate(Coefficient(a), Coefficient(b)); - p = swap(p,i,d-1); - } - return p; - } - }; - -// Scale_homogeneous; - struct Scale_homogeneous{ - typedef Polynomial_d result_type; - typedef Polynomial_d first_argument_type; - typedef Innermost_coefficient second_argument_type; - typedef Innermost_coefficient third_argument_type; - - Polynomial_d - operator()( - Polynomial_d p, - const Innermost_coefficient& a, - const Innermost_coefficient& b, - int i = (d-1) ) const { - CGAL_precondition( ! CGAL::is_zero(b) ); - - if (i == (d-1) ) p = Swap()(p,i,d-1); - - if(CGAL::is_one(b)) - p.scale_up(Coefficient(a)); - else - if(CGAL::is_one(a)) - p.scale_down(Coefficient(b)); - else - p.scale(Coefficient(a), Coefficient(b) ); - - if (i == (d-1) ) p = Swap()(p,i,d-1); - - return p; - } - }; - -// Differentiate; - struct Differentiate - : public Unary_function{ - Polynomial_d - operator()(Polynomial_d p, int i = (d-1)) const { - if (i == (d-1) ){ - p.diff(); - }else{ - Swap swap; - p = swap(p,i,d-1); - p.diff(); - p = swap(p,i,d-1); - } - return p; - } - }; -// Make_square_free; - struct Make_square_free - : public Unary_function< Polynomial_d, Polynomial_d >{ - Polynomial_d - operator()(const Polynomial_d& p) const { - if (CGAL::is_zero(p)) return p; - Univariate_content_up_to_constant_factor ucontent_utcf; - Integral_division_up_to_constant_factor idiv_utcf; - Differentiate diff; - typename PTC::Make_square_free msf; - - Coefficient content = ucontent_utcf(p); - Polynomial_d result = Polynomial_d(msf(content)); - - Polynomial_d regular_part = idiv_utcf(p,Polynomial_d(content)); - Polynomial_d g = gcd_utcf(regular_part,diff(regular_part)); - - result *= idiv_utcf(regular_part,g); - return Canonicalize()(result); - -#if 0 - - Square_free_factorization sqff_utcf; - std::vector facs; - std::vector mults; - sqff_utcf(p,std::back_inserter(facs), std::back_inserter(mults)); - for(typename std::vector::iterator it = facs.begin(); - it != facs.end(); - it++){ - result *= *it; - } - return result; -#endif - - } - }; - -// Square_free_factorization; - struct Square_free_factorization{ - typedef int result_type; - - private: - typedef Coefficient Coeff; - typedef Innermost_coefficient ICoeff; - - // rsqff computes the sqff recursively for Coeff - // end of recursion: ICoeff - template < class OutputIterator1, class OutputIterator2 > - int rsqff (ICoeff c, - OutputIterator1 factors, - OutputIterator2 mults) const{ - return 0; - } - template < class OutputIterator1, class OutputIterator2 > - int rsqff ( - typename First_if_different::Type c, - OutputIterator1 fit, - OutputIterator2 mit) const { - typename PTC::Square_free_factorization sqff; - std::vector factors; - int n = sqff(c, std::back_inserter(factors), mit); - for(int i = 0; i < (int)factors.size(); i++){ - *fit++=Polynomial_d(factors[i]); - } - return n; - } - public: - template < class OutputIterator1, class OutputIterator2 > - int operator()( - const Polynomial_d& p, - OutputIterator1 fit, - OutputIterator2 mit) const { - Coefficient c; - int n = square_free_factorization(p,fit,mit,c); - if (Total_degree()(c) > 0) - return rsqff(c,fit,mit)+n; - else - return n; - } - }; - -// Pseudo_division; - struct Pseudo_division { - typedef Polynomial_d result_type; - void - operator()( - const Polynomial_d& f, const Polynomial_d& g, - Polynomial_d& q, Polynomial_d& r, Coefficient& D) const { - Polynomial_d::pseudo_division(f,g,q,r,D); - } - }; -// Pseudo_division_remainder; - struct Pseudo_division_remainder - :public Binary_function { - - Polynomial_d - operator()(const Polynomial_d& f, const Polynomial_d& g) const { - Polynomial_d q,r; - Coefficient D; - Polynomial_d::pseudo_division(f,g,q,r,D); - return r; - } - }; -// Pseudo_division_quotient; - struct Pseudo_division_quotient - :public Binary_function { - - Polynomial_d - operator()(const Polynomial_d& f, const Polynomial_d& g) const { - Polynomial_d q,r; - Coefficient D; - Polynomial_d::pseudo_division(f,g,q,r,D); - return q; - } - }; -// Gcd_up_to_constant_factor; - struct Gcd_up_to_constant_factor - :public Binary_function { - Polynomial_d - operator()(const Polynomial_d& p, const Polynomial_d& q) const { - if (CGAL::is_zero(p) && CGAL::is_zero(q)) - return Polynomial_d(0); - return gcd_utcf(p,q); - } - }; - -// Integral_division_up_to_constant_factor; - struct Integral_division_up_to_constant_factor - :public Binary_function { - Polynomial_d - operator()(const Polynomial_d& p, const Polynomial_d& q) const { - typedef Innermost_coefficient IC; - - typename PT::Construct_polynomial construct; - typename PT::Innermost_leading_coefficient ilcoeff; - typename PT::Innermost_coefficient_begin begin; - typename PT::Innermost_coefficient_end end; - typedef Algebraic_extension_traits AET; - typename AET::Denominator_for_algebraic_integers dfai; - typename AET::Normalization_factor nfac; - - - IC ilcoeff_q = ilcoeff(q); - // this factor is needed in case IC is an Algebraic extension - IC dfai_q = dfai(begin(q), end(q)); - // make dfai_q a 'scalar' - ilcoeff_q *= dfai_q * nfac(dfai_q); - Polynomial_d result = (p * construct(ilcoeff_q)) / q; - - return Canonicalize()(result); - } - }; - -// Univariate_content_up_to_constant_factor; - struct Univariate_content_up_to_constant_factor - :public Unary_function { - Coefficient - operator()(const Polynomial_d& p) const { - typename PTC::Gcd_up_to_constant_factor gcd_utcf; - - if(CGAL::is_zero(p)) return Coefficient(0); - if(PT::d == 1) return Coefficient(1); - - Coefficient result(0); - for(typename Polynomial_d::const_iterator it = p.begin(); - it != p.end(); - it++){ - result = gcd_utcf(*it,result); - } - return result; - - } - }; -// Square_free_factorization_up_to_constant_factor; - struct Square_free_factorization_up_to_constant_factor { - typedef int result_type; - private: - typedef Coefficient Coeff; - typedef Innermost_coefficient ICoeff; - - // rsqff_utcf computes the sqff recursively for Coeff - // end of recursion: ICoeff - template < class OutputIterator1, class OutputIterator2 > - int rsqff_utcf (ICoeff c, - OutputIterator1 factors, - OutputIterator2 mults) const{ - return 0; - } - template < class OutputIterator1, class OutputIterator2 > - int rsqff_utcf ( - typename First_if_different::Type c, - OutputIterator1 fit, - OutputIterator2 mit) const { - typename PTC::Square_free_factorization sqff; - std::vector factors; - int n = sqff(c, std::back_inserter(factors), mit); - for(int i = 0; i < (int)factors.size(); i++){ - *fit++=Polynomial_d(factors[i]); - } - return n; - } - public: - template < class OutputIterator1, class OutputIterator2 > - int operator()( - Polynomial_d p, - OutputIterator1 fit, - OutputIterator2 mit) const { - - if (CGAL::is_zero(p)) return 0; - - Univariate_content_up_to_constant_factor ucontent_utcf; - Integral_division_up_to_constant_factor idiv_utcf; - Coefficient c = ucontent_utcf(p); - p = idiv_utcf( p , Polynomial_d(c)); - int n = square_free_factorization_utcf(p,fit,mit); - if (Total_degree()(c) > 0) - return rsqff_utcf(c,fit,mit)+n; - else - return n; - } - }; - -// Evaluate; - struct Evaluate - :public Binary_function{ - Coefficient - operator()(const Polynomial_d& p, Innermost_coefficient x, int i = (d-1)) - const { - if(i == (d-1) ) - return p.evaluate(x); - else{ - return Move()(p,i).evaluate(x); - } - } - }; -// Evaluate_homogeneous; - struct Evaluate_homogeneous{ - typedef Coefficient result_type; - typedef Polynomial_d first_argument_type; - typedef Innermost_coefficient second_argument_type; - typedef Innermost_coefficient third_argument_type; - typedef Arity_tag< 3 > Arity; - - Coefficient - operator()( - const Polynomial_d& p, - Innermost_coefficient a, - Innermost_coefficient b) const { - return p.evaluate_homogeneous(a,b); - } - Coefficient - operator()( - const Polynomial_d& p, - Innermost_coefficient a, - Innermost_coefficient b, - int hdegree , - int i = (PT::d-1) ) const { - if (i == (d-1) ) - return p.evaluate_homogeneous(a,b,hdegree); - else - return Move()(p,i,PT::d-1).evaluate_homogeneous(a,b,hdegree); - } - }; - -// Resultant; - struct Resultant - : public Binary_function{ - - Coefficient - operator()( - const Polynomial_d& p, - const Polynomial_d& q, - int i = (d-1) ) const { - if(i == (d-1) ) - return prs_resultant(p,q); - else - return prs_resultant(Move()(p,i),Move()(q,i)); - } - }; - -// Canonicalize; - struct Canonicalize - : public Unary_function{ - Polynomial_d - operator()( const Polynomial_d& p ) const { - return CGAL::INTERN_POLYNOMIAL::canonicalize_polynomial(p); - } - }; + }; CGAL_END_NAMESPACE