diff --git a/Polynomial/include/CGAL/Polynomial/fwd.h b/Polynomial/include/CGAL/Polynomial/fwd.h
index 272ffdf6d4d..8a6317da537 100644
--- a/Polynomial/include/CGAL/Polynomial/fwd.h
+++ b/Polynomial/include/CGAL/Polynomial/fwd.h
@@ -117,6 +117,15 @@ template<typename OutputIterator, typename NT> OutputIterator
                                           OutputIterator2 out_f,
                                           OutputIterator3 out_fx);
 
+// eliminates outermost variable
+template <class Coeff> 
+inline Coeff new_resultant( 
+    const CGAL::Polynomial<Coeff>&, const CGAL::Polynomial<Coeff>&);
+// eliminates innermost variable 
+template <class Coeff> 
+inline Coeff new_resultant_( 
+    const CGAL::Polynomial<Coeff>&, const CGAL::Polynomial<Coeff>&);
+  
 
 } // namespace CGALi
 
diff --git a/Polynomial/include/CGAL/Polynomial/resultant.h b/Polynomial/include/CGAL/Polynomial/resultant.h
index af1ae07349f..5461a9e4423 100644
--- a/Polynomial/include/CGAL/Polynomial/resultant.h
+++ b/Polynomial/include/CGAL/Polynomial/resultant.h
@@ -1,109 +1,396 @@
-// TODO: Add licence
+// Copyright (c) 1997-2000  Max-Planck-Institute Saarbruecken (Germany).
+// All rights reserved.
+//
+// This file is part of CGAL (www.cgal.org); you may redistribute it under
+// the terms of the Q Public License version 1.0.
+// See the file LICENSE.QPL distributed with CGAL.
+//
+// Licensees holding a valid commercial license may use this file in
+// accordance with the commercial license agreement provided with the software.
 //
 // This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
 // WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
 //
 // $URL:$
-// $Id: $
-// 
+// $Id:$
 //
-// Author(s)     :  
 //
-// ============================================================================
+// Author(s)     : Michael Hemmer <hemmer@mpi-inf.mpg.de>
+
+
+#ifndef CGAL_NEW_RESULTANT_H
+#define CGAL_NEW_RESULTANT_H
+
+// Modular arithmetic is slower, hence the default is 0
+#ifndef CGAL_RESULTANT_USE_MODULAR_ARITHMETIC
+#define CGAL_RESULTANT_USE_MODULAR_ARITHMETIC 0
+#endif
+
+#ifndef CGAL_RESULTANT_USE_DECOMPOSE
+#define CGAL_RESULTANT_USE_DECOMPOSE 1
+#endif 
+
+#ifndef CGAL_RESULTANT_USE_INTERPOLATION
+#define CGAL_RESULTANT_USE_INTERPOLATION 1
+#endif 
 
-// TODO: The comments are all original EXACUS comments and aren't adapted. So
-//         they may be wrong now.
 
-#ifndef CGAL_POLYNOMIAL_RESULTANT_H
-#define CGAL_POLYNOMIAL_RESULTANT_H
 
 #include <CGAL/basic.h>
 #include <CGAL/Polynomial.h>
+
+#include <CGAL/Polynomial_traits_d.h>
+#include <CGAL/Polynomial/Interpolator.h>
 #include <CGAL/Polynomial/prs_resultant.h>
 #include <CGAL/Polynomial/bezout_matrix.h>
+#include <CGAL/Polynomial/subresultants.h>
 
-#include <boost/type_traits.hpp>
-#include <boost/mpl/if.hpp>
+#include <CGAL/Modular.h>
+#include <CGAL/Modular_traits.h>
+#include <CGAL/Chinese_remainder_traits.h>
+#include <CGAL/primes.h>
+#include <CGAL/Polynomial/Cached_extended_euclidean_algorithm.h>
 
 CGAL_BEGIN_NAMESPACE
 
-namespace CGALi {
 
-    template <class NT> inline
-    NT resultant_(Polynomial<NT> A, Polynomial<NT> B) {
-        typedef typename Algebraic_structure_traits<NT>::Algebraic_category Algebraic_category;
-        return resultant_(A,B,Algebraic_category());     
-    }
-    
-    // the general function for CGAL::Integral_domain_without_div_tag
-    template < class NT> inline 
-    NT resultant_(Polynomial<NT> A, 
-                  Polynomial<NT> B,
-                  Integral_domain_without_division_tag){
-        return hybrid_bezout_subresultant(A,B,0);
-    }
-    
-    // the specialization for CGAL::Unique_factorization_domain_tag
-    template < class NT> inline
-    NT resultant_(Polynomial<NT> A, 
-                  Polynomial<NT> B,
-                  Unique_factorization_domain_tag ){
-        return prs_resultant_ufd(A,B);
-    }
-        
-    template <class NT> inline
-    NT resultant_field(Polynomial<NT> A, Polynomial<NT> B, ::CGAL::Tag_false) {       
-        return prs_resultant_field(A,B);       
-    } 
-    template <class NT> inline
-    NT resultant_field(Polynomial<NT> A, Polynomial<NT> B, ::CGAL::Tag_true) {       
-        return prs_resultant_decompose(A,B);       
-    } 
-    
-    // the specialization for CGAL::Field_tag
-    template < class NT> inline
-    NT resultant_(Polynomial<NT> A, Polynomial<NT> B, Field_tag){
-        // prs_resultant_decompose can only be used if 
-        // NT is decomposable && the resulting type is at least a UFDomain
-        typedef typename Fraction_traits<NT>::Is_fraction Is_decomposable;
-        const bool is_decomposable 
-            = ::boost::is_same<Is_decomposable, ::CGAL::Tag_true>::value;
-        
-        
-        typedef typename Fraction_traits<NT>::Numerator_type NUM;
-        typedef typename Algebraic_structure_traits<NUM>::Algebraic_category ALG_TYPE;
-        const bool is_inheritance 
-            = ::boost::is_base_and_derived<Unique_factorization_domain_tag,ALG_TYPE>::value ||
-              ::boost::is_same<Unique_factorization_domain_tag,ALG_TYPE>::value;
+// The main function provided within this file is CGAL::CGALi::resultant(F,G),
+// all other functions are used for dispatching. 
+// The implementation uses interpolatation for multivariate polynomials
+// Due to the recursive structuture of CGAL::Polynomial<Coeff> it is better 
+// to write the function such that the inner most variabel is eliminated. 
+// However,  CGAL::CGALi::resultant(F,G) eliminates the outer most variabel.
+// This is due to backward compatibility issues with code base on EXACUS. 
+// In turn CGAL::CGALi::resultant_(F,G) eliminates the innermost variable. 
+ 
+// Dispatching
+// CGAL::CGALi::new_resultant_decompose applies if Coeff is a Fraction
+// CGAL::CGALi::new_resultant_modularize applies if Coeff is Modularizable
+// CGAL::CGALi::new_resultant_interpolate applies for multivairate polynomials
+// CGAL::CGALi::new_resultant_univariate selects the proper algorithm for IC 
 
-        return resultant_field(
-          A, B,
-          typename ::boost::mpl::if_c< is_decomposable && is_inheritance,
-                                 ::CGAL::Tag_true,
-                                 ::CGAL::Tag_false >::type() );
-    }    
+namespace CGALi{
 
-    /*! \ingroup CGAL_Polynomial
-     *  \relates CGAL::Polynomial
-     *  \brief compute the resultant of the polynomials \c A and \c B
-     *
-     *  The way the resultant is computed depends on the Algebraic_category. 
-     *  In general the resultant will be computed by the function
-     *  CGAL::hybrid_bezout_subresultant, but if possible the function
-     *  CGAL::prs_resultant_ufd or CGAL::prs_resultant_field are used. 
-     *  
-     *  Up to now it is not clear, that the functions based on the polynomial
-     *  remainder sequence are faster than the one based on the bezoutian. 
-     *  Thus you can use CGAL::hybrid_bezout_subresultant instead, which will 
-     *  work for any Algebraic_category
-     */
-    template <class NT> inline
-    NT resultant(Polynomial<NT> A, Polynomial<NT> B) {
-        return CGALi::resultant_(A, B);
-    }   
+template <class Coeff> 
+inline Coeff new_resultant_interpolate( 
+    const CGAL::Polynomial<Coeff>&, const CGAL::Polynomial<Coeff>& );
+template <class Coeff> 
+inline Coeff new_resultant_modularize( 
+    const CGAL::Polynomial<Coeff>&, 
+    const CGAL::Polynomial<Coeff>&, CGAL::Tag_true);
+template <class Coeff> 
+inline Coeff new_resultant_modularize( 
+    const CGAL::Polynomial<Coeff>&, 
+    const CGAL::Polynomial<Coeff>&, CGAL::Tag_false);
+template <class Coeff> 
+inline Coeff new_resultant_decompose( 
+    const CGAL::Polynomial<Coeff>&, 
+    const CGAL::Polynomial<Coeff>&, CGAL::Tag_true);
+template <class Coeff> 
+inline Coeff new_resultant_decompose( 
+    const CGAL::Polynomial<Coeff>&, 
+    const CGAL::Polynomial<Coeff>&, CGAL::Tag_false);
+template <class Coeff> 
+inline Coeff new_resultant_( 
+    const CGAL::Polynomial<Coeff>&, const CGAL::Polynomial<Coeff>&);
+
+template <class Coeff> 
+inline Coeff new_resultant_univariate( 
+    const CGAL::Polynomial<Coeff>& A, 
+    const CGAL::Polynomial<Coeff>& B, CGAL::Integral_domain_without_division_tag){ 
+  return hybrid_bezout_subresultant(A,B,0);
+}
+template <class Coeff> 
+inline Coeff new_resultant_univariate( 
+    const CGAL::Polynomial<Coeff>& A, 
+    const CGAL::Polynomial<Coeff>& B, CGAL::Unique_factorization_domain_tag){
+  return prs_resultant_ufd(A,B);
+}
+
+template <class Coeff> 
+inline Coeff new_resultant_univariate( 
+    const CGAL::Polynomial<Coeff>& A, 
+    const CGAL::Polynomial<Coeff>& B, CGAL::Field_tag){
+  return prs_resultant_field(A,B);  
+}
 
 } // namespace CGALi
 
+namespace CGALi{
+
+
+template <class IC> 
+inline IC 
+new_resultant_interpolate( 
+    const CGAL::Polynomial<IC>& F, 
+    const CGAL::Polynomial<IC>& G){
+  CGAL_precondition(CGAL::Polynomial_traits_d<CGAL::Polynomial<IC> >::d == 1);
+    typedef CGAL::Algebraic_structure_traits<IC> AST_IC;
+    typedef typename AST_IC::Algebraic_category Algebraic_category;
+    return CGALi::new_resultant_univariate(F,G,Algebraic_category()); 
+}
+
+#if CGAL_RESULTANT_USE_INTERPOLATION
+template <class Coeff_2> 
+inline
+CGAL::Polynomial<Coeff_2>  new_resultant_interpolate(
+        const CGAL::Polynomial<CGAL::Polynomial<Coeff_2> >& F, 
+        const CGAL::Polynomial<CGAL::Polynomial<Coeff_2> >& G){
+    
+    typedef CGAL::Polynomial<Coeff_2> Coeff_1;
+    typedef CGAL::Polynomial<Coeff_1> POLY;
+    typedef CGAL::Polynomial_traits_d<POLY> PT;
+    typedef typename PT::Innermost_coefficient IC; 
+
+    CGAL_precondition(PT::d >= 2);
+    
+    typename PT::Degree degree; 
+    typename CGAL::Polynomial_traits_d<Coeff_1>::Degree_vector degree_vector; 
+
+    int maxdegree = degree(F,0)*degree(G,PT::d-1) + degree(F,PT::d-1)*degree(G,0); 
+
+    typedef std::pair<IC,Coeff_2> Point; 
+    std::vector<Point> points; // interpolation points  
+    
+   
+    int i(0);
+    CGAL::Exponent_vector ev_f(degree_vector(Coeff_1()));
+    CGAL::Exponent_vector ev_g(degree_vector(Coeff_1()));
+    
+   
+    while((int) points.size() <= maxdegree + 1){
+        i++;
+        // timer1.start();
+        Coeff_1 Fat_i = F.evaluate(Coeff_1(i));
+        Coeff_1 Gat_i = G.evaluate(Coeff_1(i));
+        // timer1.stop();
+        
+        if(degree_vector(Fat_i) >  ev_f || degree_vector(Gat_i) >  ev_g){
+            points.clear();
+            ev_f  = degree_vector(Fat_i);
+            ev_g  = degree_vector(Gat_i);
+            CGAL_postcondition(points.size() == 0);
+        }
+        if(degree_vector(Fat_i) ==  ev_f && degree_vector(Gat_i) ==  ev_g){
+            // timer2.start();
+            Coeff_2 res_at_i = new_resultant_interpolate(Fat_i, Gat_i);
+            // timer2.stop();
+            points.push_back(Point(IC(i),res_at_i));
+        }      
+    }
+   
+    // timer3.start();
+    CGAL::CGALi::Interpolator<Coeff_1> interpolator(points.begin(),points.end());
+    Coeff_1 result = interpolator.get_interpolant();
+    // timer3.stop();
+
+#ifndef CGAL_NDEBUG
+    while((int) points.size() <= maxdegree + 3){
+        i++;
+        Coeff_1 Fat_i = typename PT::Evaluate()(F,IC(i));
+        Coeff_1 Gat_i = typename PT::Evaluate()(G,IC(i));
+        
+        assert(degree_vector(Fat_i) <= ev_f);
+        assert(degree_vector(Gat_i) <= ev_g);
+        
+        if(degree_vector(Fat_i) ==  ev_f && degree_vector(Gat_i) ==  ev_g){
+            Coeff_2 res_at_i = new_resultant_interpolate(Fat_i, Gat_i);
+            points.push_back(Point(IC(i), res_at_i));
+        }
+    }
+    CGAL::CGALi::Interpolator<Coeff_1> 
+      interpolator_(points.begin(),points.end());
+    Coeff_1 result_= interpolator_.get_interpolant();
+    
+     // the interpolate polynomial has to be stable !
+    assert(result_ == result); 
+#endif 
+    return result; 
+}
+#endif
+
+
+
+template <class Coeff> 
+inline
+Coeff new_resultant_modularize( 
+        const CGAL::Polynomial<Coeff>& F, 
+        const CGAL::Polynomial<Coeff>& G, 
+        CGAL::Tag_false){
+    return new_resultant_interpolate(F,G);
+};
+
+template <class Coeff> 
+inline
+Coeff new_resultant_modularize( 
+        const CGAL::Polynomial<Coeff>& F, 
+        const CGAL::Polynomial<Coeff>& G, 
+        CGAL::Tag_true){
+    
+    typedef Polynomial_traits_d<CGAL::Polynomial<Coeff> > PT;
+    typedef typename PT::Polynomial_d Polynomial;
+    typedef typename PT::Innermost_coefficient IC;
+    
+    
+    typedef Chinese_remainder_traits<Coeff> CRT;
+    typedef typename CRT::Scalar_type Scalar;
+
+
+    typedef typename CGAL::Modular_traits<Polynomial>::Modular_NT MPolynomial; 
+    typedef typename CGAL::Modular_traits<Coeff>::Modular_NT      MCoeff; 
+    typedef typename CGAL::Modular_traits<Scalar>::Modular_NT     MScalar;  
+        
+    typename CRT::Chinese_remainder chinese_remainder; 
+    typename CGAL::Modular_traits<Coeff>::Modular_image_inv inv_map;
+
+
+    typename PT::Degree_vector                                     degree_vector; 
+    typename CGAL::Polynomial_traits_d<MPolynomial>::Degree_vector mdegree_vector;
+
+    bool solved = false; 
+    int prime_index = 0; 
+    int n = 0;
+    Scalar p,q,pq,s,t; 
+    Coeff R, R_old; 
+    
+    // CGAL::Timer timer_evaluate, timer_resultant, timer_cr; 
+    
+    do{
+        MPolynomial mF, mG;
+        MCoeff mR;
+        //timer_evaluate.start();
+        do{
+            // select a prime number
+            int current_prime = -1;
+            prime_index++;
+            if(prime_index >= 2000){
+                std::cerr<<"primes in the array exhausted"<<std::endl;
+                assert(false);
+                current_prime = CGALi::get_next_lower_prime(current_prime);
+            } else{
+                current_prime = CGALi::primes[prime_index];
+            }
+            CGAL::Modular::set_current_prime(current_prime);
+            
+            mF = CGAL::modular_image(F);
+            mG = CGAL::modular_image(G);
+            
+        }while( degree_vector(F) != mdegree_vector(mF) || 
+                degree_vector(G) != mdegree_vector(mG));
+        //timer_evaluate.stop();
+        
+        //timer_resultant.start();
+        n++;
+        mR = new_resultant_interpolate(mF,mG);
+        //timer_resultant.stop();
+        //timer_cr.start();
+        if(n == 1){ 
+            // init chinese remainder
+            q =  CGAL::Modular::get_current_prime(); // implicit ! 
+            R = inv_map(mR);
+        }else{
+            // continue chinese remainder
+            p = CGAL::Modular::get_current_prime(); // implicit!  
+            R_old  = R ;
+//            chinese_remainder(q,Gs ,p,inv_map(mG_),pq,Gs);             
+//            cached_extended_euclidean_algorithm(q,p,s,t);
+            CGALi::Cached_extended_euclidean_algorithm
+                <typename CRT::Scalar_type> ceea;
+            ceea(q,p,s,t);
+            pq =p*q;
+            chinese_remainder(q,p,pq,s,t,R_old,inv_map(mR),R);
+            q=pq;
+        }
+        solved = (R==R_old);
+        //timer_cr.stop();       
+    } while(!solved);
+        
+    //std::cout << "Time Evaluate   : " << timer_evaluate.time() << std::endl; 
+    //std::cout << "Time Resultant  : " << timer_resultant.time() << std::endl; 
+    //std::cout << "Time Chinese R  : " << timer_cr.time() << std::endl; 
+    // CGAL_postcondition(R == new_resultant_interpolate(F,G));
+    return R;
+    // return new_resultant_interpolate(F,G);
+}
+
+
+template <class Coeff> 
+inline
+Coeff new_resultant_decompose( 
+    const CGAL::Polynomial<Coeff>& F,
+    const CGAL::Polynomial<Coeff>& G, 
+    CGAL::Tag_false){
+#if CGAL_RESULTANT_USE_MODULAR_ARITHMETIC
+  typedef CGAL::Polynomial<Coeff> Polynomial; 
+  typedef typename Modular_traits<Polynomial>::Is_modularizable Is_modularizable; 
+  return new_resultant_modularize(F,G,Is_modularizable());
+#else
+  return  new_resultant_modularize(F,G,CGAL::Tag_false());
+#endif
+}
+
+template <class Coeff> 
+inline
+Coeff new_resultant_decompose( 
+        const CGAL::Polynomial<Coeff>& F, 
+        const CGAL::Polynomial<Coeff>& G, 
+        CGAL::Tag_true){  
+    
+    typedef Polynomial<Coeff> POLY;
+    typedef typename Fraction_traits<POLY>::Numerator_type Numerator;
+    typedef typename Fraction_traits<POLY>::Denominator_type Denominator;
+    typename Fraction_traits<POLY>::Decompose decompose;
+    typedef typename Numerator::NT RES;
+    
+    Denominator a, b;
+    // F.simplify_coefficients(); not const 
+    // G.simplify_coefficients(); not const 
+    Numerator F0; decompose(F,F0,a);
+    Numerator G0; decompose(G,G0,b);
+    Denominator c = CGAL::ipower(a, G.degree()) * CGAL::ipower(b, F.degree());
+    typedef Algebraic_structure_traits<RES> AST_RES;
+    typedef typename AST_RES::Algebraic_category Algebraic_category;
+    RES res0 =  CGAL::CGALi::new_resultant_(F0, G0);
+    typename Fraction_traits<Coeff>::Compose comp_frac;
+    Coeff res = comp_frac(res0, c);
+    typename Algebraic_structure_traits<Coeff>::Simplify simplify;
+    simplify(res);
+    return res;
+}
+
+
+template <class Coeff> 
+inline
+Coeff new_resultant_( 
+        const CGAL::Polynomial<Coeff>& F, 
+        const CGAL::Polynomial<Coeff>& G){
+#if CGAL_RESULTANT_USE_DECOMPOSE
+    typedef CGAL::Fraction_traits<Polynomial<Coeff > > FT;
+    typedef typename FT::Is_fraction Is_fraction; 
+    return new_resultant_decompose(F,G,Is_fraction());
+#else
+    return new_resultant_decompose(F,G,CGAL::Tag_false());
+#endif
+}
+
+
+
+template <class Coeff> 
+inline
+Coeff  new_resultant( 
+        const CGAL::Polynomial<Coeff>& F_, 
+        const CGAL::Polynomial<Coeff>& G_){
+  // make the variable to be elimnated the innermost one.
+    typedef CGAL::Polynomial_traits_d<CGAL::Polynomial<Coeff> > PT;
+    CGAL::Polynomial<Coeff> F = typename PT::Move()(F_, PT::d-1, 0);
+    CGAL::Polynomial<Coeff> G = typename PT::Move()(G_, PT::d-1, 0);
+    return CGALi::new_resultant_(F,G);
+}
+
+} // namespace CGALi    
 CGAL_END_NAMESPACE
 
-#endif// CGAL_POLYNOMIAL_RESULTANT_H
+
+
+#endif // CGAL_NEW_RESULTANT_H
+
diff --git a/Polynomial/test/Polynomial/Polynomial_traits_d.cpp b/Polynomial/test/Polynomial/Polynomial_traits_d.cpp
index ac484d9d027..07aacb5134f 100644
--- a/Polynomial/test/Polynomial/Polynomial_traits_d.cpp
+++ b/Polynomial/test/Polynomial/Polynomial_traits_d.cpp
@@ -1145,7 +1145,7 @@ void test_compare(const Polynomial_traits_d&) {
     
   CGAL_SNAP_CGALi_TRAITS_D(Polynomial_traits_d);
     
-  typename PT::Compare compare;
+  typename PT::Compare compare; (void) compare;
     
   Polynomial_d p0(Coeff(0));
   Polynomial_d pp2(Coeff(2));
diff --git a/Polynomial/test/Polynomial/resultant.cpp b/Polynomial/test/Polynomial/resultant.cpp
index bc62327857d..3cf23ba70b1 100644
--- a/Polynomial/test/Polynomial/resultant.cpp
+++ b/Polynomial/test/Polynomial/resultant.cpp
@@ -1,129 +1,123 @@
 
-// #define CGAL_RESULTANT_NUSE_MODULAR_ARITHMETIC 1
+// currently using modular arithmetic is far to slow,
+// that is, we use interpolation for multivariate resultants 
+// but no modular arithmetic. 
+
+// These are the defaults 
+//#define CGAL_RESULTANT_USE_MODULAR_ARITHMETIC 0
+//#define CGAL_RESULTANT_USE_DECOMPOSE 1
+//#define CGAL_RESULTANT_USE_INTERPOLATE 1 
 
 #include <CGAL/basic.h>
 
-#include <CGAL/Timer.h>
-CGAL::Timer timer1, timer2, timer3; 
-
 #include <cassert>
 #include <iostream>
 
+#include <CGAL/Polynomial/resultant.h>
+
 #include <CGAL/Arithmetic_kernel.h>
 #include <CGAL/Modular.h>
+#include <CGAL/Sqrt_extension.h>
 #include <CGAL/Polynomial.h>
 #include <CGAL/Polynomial_traits_d.h>
+
 #include <CGAL/Random.h>
+#include <CGAL/Timer.h>
 
 #include <CGAL/gen_sparse_polynomial.h>
-#include <CGAL/new_resultant.h>
-
 
 static CGAL::Random my_rnd(346); // some seed 
 
 template<class Polynomial_d> 
 void test_resultant(){
-    
-    typedef typename CGAL::Polynomial_traits_d<Polynomial_d> PT;
-    typedef typename PT::Coefficient Coeff; 
-    for (int k = 0; k < 1; k++){
-        Polynomial_d F2 = 
-            CGAL::generate_sparse_random_polynomial<Polynomial_d>(my_rnd,12); 
-        Polynomial_d F1 = 
-            CGAL::generate_sparse_random_polynomial<Polynomial_d>(my_rnd,12);
-        Polynomial_d G1 = typename PT::Move()(F1,PT::d-1,0);
-        Polynomial_d G2 = typename PT::Move()(F2,PT::d-1,0);
-        
-        CGAL::Timer timer_new;
-        timer_new.start();
-        Coeff F = CGAL::new_resultant(F1,F2,PT::d-1);
-        // Coeff G = CGAL::new_resultant(G1,G2,0);
-        timer_new.stop();
-
-        assert(F==F);
-        std::cout <<"timer1: "<< timer1.time() << std::endl;
-        std::cout <<"timer2: "<< timer2.time() << std::endl;
-        std::cout <<"timer3: "<< timer3.time() << std::endl;
-        timer1.reset();
-        timer2.reset();
-        timer3.reset();
-        std::cout <<"new time: "<< timer_new.time() << std::endl;
-    }
-    std::cout << "end resultant: " << PT::d << std::endl;
-}
+  typedef typename CGAL::Polynomial_traits_d<Polynomial_d> PT;
+  typedef typename PT::Coefficient Coeff; 
+  typename PT::Resultant resultant;
    
-
-template<class Coeff> 
-void test_multi_resultant(){
+  {
+    Polynomial_d A(0);
+    Polynomial_d B(0);
+    assert(resultant(A,B)==Coeff(0));
+  }{
+    Polynomial_d A(4);
+    Polynomial_d B(8);
+    assert(resultant(A,B)==Coeff(1));
+  }{
+    Polynomial_d f(Coeff(2),Coeff(7),Coeff(1),Coeff(8),Coeff(1),Coeff(8));
+    Polynomial_d g(Coeff(3),Coeff(1),Coeff(4),Coeff(1),Coeff(5),Coeff(9));
+    assert(resultant(f,g) == Coeff(230664271L)); // Maple
+        
+    Polynomial_d h(Coeff(3),Coeff(4),Coeff(7),Coeff(7));
+    Polynomial_d fh(f*h);
+    Polynomial_d gh(g*h);
+    assert(resultant(fh,gh) == Coeff(0) );
+  } 
     
-    typedef CGAL::Polynomial<Coeff>        Polynomial_1;
-    typedef CGAL::Polynomial<Polynomial_1> Polynomial_2;
-    typedef CGAL::Polynomial<Polynomial_2> Polynomial_3;
 
-    Polynomial_3 Q  =
-        CGAL::generate_polynomial_degree_each<Polynomial_3>(my_rnd,4,100); 
-    Polynomial_2 B1 =
-        CGAL::generate_polynomial_degree_each<Polynomial_2>(my_rnd,6,80); 
-    Polynomial_2 B2 = 
-        CGAL::generate_polynomial_degree_each<Polynomial_2>(my_rnd,6,80); 
-    std::cout << "Q:  "<< Q   << std::endl;
-    std::cout << "B1: "<< B1 << std::endl;
-    std::cout << "B2: "<< B2 <<std::endl;
-    
+  for (int k = 0; k < 1; k++){
+    Polynomial_d F2 = 
+      CGAL::generate_sparse_random_polynomial<Polynomial_d>(my_rnd,5); 
+    Polynomial_d F1 = 
+      CGAL::generate_sparse_random_polynomial<Polynomial_d>(my_rnd,5);
+    Polynomial_d G1 = typename PT::Move()(F1,PT::d-1,0);
+    Polynomial_d G2 = typename PT::Move()(F2,PT::d-1,0);
+        
     CGAL::Timer timer_new;
     timer_new.start();
-    Polynomial_2 R1 = CGAL::new_resultant(Q,Polynomial_3(B1),0);
-    Polynomial_1 R2 = CGAL::new_resultant(R1,B2,0);
-    std::cout << std::endl;
-    std::cout << R1.degree() << std::endl;
-    std::cout << R2.degree() << std::endl;
+    Coeff F = resultant(F1,F2,PT::d-1);
+    Coeff G = resultant(G1,G2,0);
     timer_new.stop();
-    
-    std::cout <<"timer1: "<< timer1.time() << std::endl;
-    std::cout <<"timer2: "<< timer2.time() << std::endl;
-    std::cout <<"timer3: "<< timer3.time() << std::endl;
-    timer1.reset();
-    timer2.reset();
-    timer3.reset();
 
-    std::cout <<"new time: "<< timer_new.time() << std::endl;
-    std::cout << "end multi resultant: " << std::endl;
+    assert(F==G);
+  }
+  std::cout << "end resultant: " << PT::d << std::endl;
 }
    
 
 int main(){
 
-    //CGAL::force_ieee_double_precision();
-    CGAL::set_pretty_mode(std::cout);
+  CGAL::Timer timer; 
+  timer.start();
+  //CGAL::force_ieee_double_precision();
+  CGAL::set_pretty_mode(std::cout);
+  
+  typedef CGAL::Arithmetic_kernel AK;
+  typedef AK::Integer Integer; 
+  typedef CGAL::Polynomial<Integer>      Polynomial_1;
+  typedef CGAL::Polynomial<Polynomial_1> Polynomial_2;
+  typedef CGAL::Polynomial<Polynomial_2> Polynomial_3;
 
-    typedef CGAL::Arithmetic_kernel AK;
-    typedef AK::Integer Integer; 
-    typedef CGAL::Polynomial<Integer>      Polynomial_1;
-    typedef CGAL::Polynomial<Polynomial_1> Polynomial_2;
-    typedef CGAL::Polynomial<Polynomial_2> Polynomial_3;
+  test_resultant<Polynomial_1>();
+  test_resultant<Polynomial_2>();
+  test_resultant<Polynomial_3>();
+   
+  typedef CGAL::Polynomial<CGAL::Modular>  MPolynomial_1;
+  typedef CGAL::Polynomial<MPolynomial_1>  MPolynomial_2;
+  typedef CGAL::Polynomial<MPolynomial_2>  MPolynomial_3;
 
-    test_multi_resultant<Integer>();
-    
-    // test_resultant<Polynomial_1>();
-    // test_resultant<Polynomial_2>();
-    //test_resultant<Polynomial_3>();
+  test_resultant<MPolynomial_1>();
+  test_resultant<MPolynomial_2>();
+  test_resultant<MPolynomial_3>();
+
+  typedef CGAL::Polynomial<AK::Rational>   RPolynomial_1;
+  typedef CGAL::Polynomial<RPolynomial_1>  RPolynomial_2;
+  typedef CGAL::Polynomial<RPolynomial_2>  RPolynomial_3;
+
+  test_resultant<RPolynomial_1>();
+  test_resultant<RPolynomial_2>();
+  test_resultant<RPolynomial_3>();  
     
+  typedef CGAL::Sqrt_extension<AK::Integer, AK::Integer> EXT;
+  typedef CGAL::Polynomial<EXT>            EPolynomial_1;
+  typedef CGAL::Polynomial<EPolynomial_1>  EPolynomial_2;
+  typedef CGAL::Polynomial<EPolynomial_2>  EPolynomial_3;
 
-    typedef CGAL::Polynomial<CGAL::Modular>  MPolynomial_1;
-    typedef CGAL::Polynomial<MPolynomial_1>  MPolynomial_2;
-    typedef CGAL::Polynomial<MPolynomial_2>  MPolynomial_3;
-
-    //test_resultant<MPolynomial_1>();
-    //test_resultant<MPolynomial_2>();
-    //test_resultant<MPolynomial_3>();
-
-    typedef CGAL::Polynomial<AK::Rational>   RPolynomial_1;
-    typedef CGAL::Polynomial<RPolynomial_1>  RPolynomial_2;
-    typedef CGAL::Polynomial<RPolynomial_2>  RPolynomial_3;
-
-    //test_resultant<RPolynomial_1>();
-    //test_resultant<RPolynomial_2>();
-    //test_resultant<RPolynomial_3>();
+  test_resultant<EPolynomial_1>();
+  test_resultant<EPolynomial_2>();
+  test_resultant<EPolynomial_3>();
+  timer.stop();
+  
+  std::cout <<" TOTAL TIME: " << timer.time();
 }