Index: factory/facFqBivar.cc
===================================================================
--- factory/facFqBivar.cc	(revision 0415f923fdb78a69504a400113e65cd371cb2150)
+++ factory/facFqBivar.cc	(revision f876a663240524d1f86ed3519e45b68d1ae877b3)
@@ -336,5 +336,5 @@
             else
             {
-              if (!isInExtension (buf2, delta, k))
+              if (!isInExtension (buf2, gamma, k, delta, source, dest))
               {
                 buf /= g;
@@ -668,5 +668,5 @@
           else
           {
-            if (!isInExtension (buf2, delta, k))
+            if (!isInExtension (buf2, gamma, k, delta, source, dest))
             {
               appendTestMapDown (result, buf2, info, source, dest);
Index: factory/facFqBivar.h
===================================================================
--- factory/facFqBivar.h	(revision 0415f923fdb78a69504a400113e65cd371cb2150)
+++ factory/facFqBivar.h	(revision f876a663240524d1f86ed3519e45b68d1ae877b3)
@@ -30,4 +30,6 @@
 #include "cf_util.h"
 #include "facFqSquarefree.h"
+#include "cf_map.h"
+#include "cfNewtonPolygon.h"
 
 static const double log2exp= 1.442695041;
@@ -53,11 +55,45 @@
 ///         element is the leading coefficient.
 /// @sa FqBiSqrfFactorize(), GFBiSqrfFactorize()
+#ifdef HAVE_NTL
 inline
-CFList FpBiSqrfFactorize (const CanonicalForm & F ///< [in] a bivariate poly
+CFList FpBiSqrfFactorize (const CanonicalForm & G ///< [in] a bivariate poly
                          )
 {
   ExtensionInfo info= ExtensionInfo (false);
+  CFMap N;
+  CanonicalForm F= compress (G, N);
+  CanonicalForm contentX= content (F, 1);
+  CanonicalForm contentY= content (F, 2);
+  F /= (contentX*contentY);
+  CFFList contentXFactors, contentYFactors;
+  contentXFactors= factorize (contentX);
+  contentYFactors= factorize (contentY);
+  if (contentXFactors.getFirst().factor().inCoeffDomain())
+    contentXFactors.removeFirst();
+  if (contentYFactors.getFirst().factor().inCoeffDomain())
+    contentYFactors.removeFirst();
+  if (F.inCoeffDomain())
+  {
+    CFList result;
+    for (CFFListIterator i= contentXFactors; i.hasItem(); i++)
+      result.append (N (i.getItem().factor()));
+    for (CFFListIterator i= contentYFactors; i.hasItem(); i++)
+      result.append (N (i.getItem().factor()));
+    normalize (result);
+    result.insert (Lc (G));
+    return result;
+  }
+  mat_ZZ M;
+  vec_ZZ S;
+  F= compress (F, M, S);
   CFList result= biFactorize (F, info);
-  result.insert (Lc(F));
+  for (CFListIterator i= result; i.hasItem(); i++)
+    i.getItem()= N (decompress (i.getItem(), M, S));
+  for (CFFListIterator i= contentXFactors; i.hasItem(); i++)
+    result.append (N(i.getItem().factor()));
+  for (CFFListIterator i= contentYFactors; i.hasItem(); i++)
+    result.append (N (i.getItem().factor()));
+  normalize (result);
+  result.insert (Lc(G));
   return result;
 }
@@ -69,11 +105,44 @@
 /// @sa FpBiSqrfFactorize(), GFBiSqrfFactorize()
 inline
-CFList FqBiSqrfFactorize (const CanonicalForm & F, ///< [in] a bivariate poly
+CFList FqBiSqrfFactorize (const CanonicalForm & G, ///< [in] a bivariate poly
                           const Variable& alpha    ///< [in] algebraic variable
                          )
 {
   ExtensionInfo info= ExtensionInfo (alpha, false);
+  CFMap N;
+  CanonicalForm F= compress (G, N);
+  CanonicalForm contentX= content (F, 1);
+  CanonicalForm contentY= content (F, 2);
+  F /= (contentX*contentY);
+  CFFList contentXFactors, contentYFactors;
+  contentXFactors= factorize (contentX);
+  contentYFactors= factorize (contentY);
+  if (contentXFactors.getFirst().factor().inCoeffDomain())
+    contentXFactors.removeFirst();
+  if (contentYFactors.getFirst().factor().inCoeffDomain())
+    contentYFactors.removeFirst();
+  if (F.inCoeffDomain())
+  {
+    CFList result;
+    for (CFFListIterator i= contentXFactors; i.hasItem(); i++)
+      result.append (N (i.getItem().factor()));
+    for (CFFListIterator i= contentYFactors; i.hasItem(); i++)
+      result.append (N (i.getItem().factor()));
+    normalize (result);
+    result.insert (Lc (G));
+    return result;
+  }
+  mat_ZZ M;
+  vec_ZZ S;
+  F= compress (F, M, S);
   CFList result= biFactorize (F, info);
-  result.insert (Lc(F));
+  for (CFListIterator i= result; i.hasItem(); i++)
+    i.getItem()= N (decompress (i.getItem(), M, S));
+  for (CFFListIterator i= contentXFactors; i.hasItem(); i++)
+    result.append (N(i.getItem().factor()));
+  for (CFFListIterator i= contentYFactors; i.hasItem(); i++)
+    result.append (N (i.getItem().factor()));
+  normalize (result);
+  result.insert (Lc(G));
   return result;
 }
@@ -85,5 +154,5 @@
 /// @sa FpBiSqrfFactorize(), FqBiSqrfFactorize()
 inline
-CFList GFBiSqrfFactorize (const CanonicalForm & F ///< [in] a bivariate poly
+CFList GFBiSqrfFactorize (const CanonicalForm & G ///< [in] a bivariate poly
                          )
 {
@@ -91,6 +160,39 @@
           "GF as base field expected");
   ExtensionInfo info= ExtensionInfo (getGFDegree(), gf_name, false);
+  CFMap N;
+  CanonicalForm F= compress (G, N);
+  CanonicalForm contentX= content (F, 1);
+  CanonicalForm contentY= content (F, 2);
+  F /= (contentX*contentY);
+  CFList contentXFactors, contentYFactors;
+  contentXFactors= biFactorize (contentX, info);
+  contentYFactors= biFactorize (contentY, info);
+  if (contentXFactors.getFirst().inCoeffDomain())
+    contentXFactors.removeFirst();
+  if (contentYFactors.getFirst().inCoeffDomain())
+    contentYFactors.removeFirst();
+  if (F.inCoeffDomain())
+  {
+    CFList result;
+    for (CFListIterator i= contentXFactors; i.hasItem(); i++)
+      result.append (N (i.getItem()));
+    for (CFListIterator i= contentYFactors; i.hasItem(); i++)
+      result.append (N (i.getItem()));
+    normalize (result);
+    result.insert (Lc (G));
+    return result;
+  }
+  mat_ZZ M;
+  vec_ZZ S;
+  F= compress (F, M, S);
   CFList result= biFactorize (F, info);
-  result.insert (Lc(F));
+  for (CFListIterator i= result; i.hasItem(); i++)
+    i.getItem()= N (decompress (i.getItem(), M, S));
+  for (CFListIterator i= contentXFactors; i.hasItem(); i++)
+    result.append (N(i.getItem()));
+  for (CFListIterator i= contentYFactors; i.hasItem(); i++)
+    result.append (N (i.getItem()));
+  normalize (result);
+  result.insert (Lc(G));
   return result;
 }
@@ -102,14 +204,38 @@
 /// @sa FqBiFactorize(), GFBiFactorize()
 inline
-CFFList FpBiFactorize (const CanonicalForm & F ///< [in] a bivariate poly
+CFFList FpBiFactorize (const CanonicalForm & G ///< [in] a bivariate poly
                       )
 {
   ExtensionInfo info= ExtensionInfo (false);
   bool GF= false;
+  CFMap N;
+  CanonicalForm F= compress (G, N);
   CanonicalForm LcF= Lc (F);
+  CanonicalForm contentX= content (F, 1);
+  CanonicalForm contentY= content (F, 2);
+  F /= (contentX*contentY);
+  CFFList contentXFactors, contentYFactors;
+  contentXFactors= factorize (contentX);
+  contentYFactors= factorize (contentY);
+  if (contentXFactors.getFirst().factor().inCoeffDomain())
+    contentXFactors.removeFirst();
+  if (contentYFactors.getFirst().factor().inCoeffDomain())
+    contentYFactors.removeFirst();
+  decompress (contentXFactors, N);
+  decompress (contentYFactors, N);
+  CFFList result, resultRoot;
+  if (F.inCoeffDomain())
+  {
+    result= Union (contentXFactors, contentYFactors);
+    normalize (result);
+    result.insert (CFFactor (LcF, 1));
+    return result;
+  }
+  mat_ZZ M;
+  vec_ZZ S;
+  F= compress (F, M, S);
   CanonicalForm pthRoot, A;
   CanonicalForm sqrfP= sqrfPart (F/Lc(F), pthRoot, info.getAlpha());
   CFList buf, bufRoot;
-  CFFList result, resultRoot;
   int p= getCharacteristic();
   int l;
@@ -120,5 +246,9 @@
     result.removeFirst();
     for (CFFListIterator i= result; i.hasItem(); i++)
-      i.getItem()= CFFactor(i.getItem().factor(),i.getItem().exp()*ipower(p,l));
+      i.getItem()= CFFactor (N (decompress (i.getItem().factor(), M, S)),
+                             i.getItem().exp()*ipower (p,l));
+    result= Union (result, contentXFactors);
+    result= Union (result, contentYFactors);
+    normalize (result);
     result.insert (CFFactor (LcF, 1));
     return result;
@@ -129,4 +259,7 @@
     A= F/LcF;
     result= multiplicity (A, buf);
+    for (CFFListIterator i= result; i.hasItem(); i++)
+      i.getItem()= CFFactor (N (decompress (i.getItem().factor(), M, S)),
+                             i.getItem().exp());
   }
   if (degree (A) > 0)
@@ -134,6 +267,12 @@
     resultRoot= FpBiFactorize (A);
     resultRoot.removeFirst();
+    for (CFFListIterator i= resultRoot; i.hasItem(); i++)
+      i.getItem()= CFFactor (N (decompress (i.getItem().factor(), M, S)),
+                             i.getItem().exp());
     result= Union (result, resultRoot);
   }
+  result= Union (result, contentXFactors);
+  result= Union (result, contentYFactors);
+  normalize (result);
   result.insert (CFFactor (LcF, 1));
   return result;
@@ -146,5 +285,5 @@
 /// @sa FpBiFactorize(), FqBiFactorize()
 inline
-CFFList FqBiFactorize (const CanonicalForm & F, ///< [in] a bivariate poly
+CFFList FqBiFactorize (const CanonicalForm & G, ///< [in] a bivariate poly
                        const Variable & alpha   ///< [in] algebraic variable
                       )
@@ -152,9 +291,34 @@
   ExtensionInfo info= ExtensionInfo (alpha, false);
   bool GF= false;
+  CFMap N;
+  CanonicalForm F= compress (G, N);
   CanonicalForm LcF= Lc (F);
-  CanonicalForm pthRoot, A;
+  CanonicalForm contentX= content (F, 1);
+  CanonicalForm contentY= content (F, 2);
+  F /= (contentX*contentY);
+  CFFList contentXFactors, contentYFactors;
+  contentXFactors= factorize (contentX);
+  contentYFactors= factorize (contentY);
+  if (contentXFactors.getFirst().factor().inCoeffDomain())
+    contentXFactors.removeFirst();
+  if (contentYFactors.getFirst().factor().inCoeffDomain())
+    contentYFactors.removeFirst();
+  decompress (contentXFactors, N);
+  decompress (contentYFactors, N);
+  CFFList result, resultRoot;
+  if (F.inCoeffDomain())
+  {
+    result= Union (contentXFactors, contentYFactors);
+    normalize (result);
+    result.insert (CFFactor (LcF, 1));
+    return result;
+  }
+  mat_ZZ M;
+  vec_ZZ S;
+  CanonicalForm oldF= F;
+  F= compress (F, M, S);
+  CanonicalForm pthRoot, A, tmp;
   CanonicalForm sqrfP= sqrfPart (F/Lc(F), pthRoot, alpha);
   CFList buf, bufRoot;
-  CFFList result, resultRoot;
   int p= getCharacteristic();
   int q= ipower (p, degree (getMipo (alpha)));
@@ -166,5 +330,9 @@
     result.removeFirst();
     for (CFFListIterator i= result; i.hasItem(); i++)
-      i.getItem()= CFFactor(i.getItem().factor(),i.getItem().exp()*ipower(p,l));
+      i.getItem()= CFFactor (N (decompress (i.getItem().factor(), M, S)),
+                             i.getItem().exp()*ipower (p,l));
+    result= Union (result, contentXFactors);
+    result= Union (result, contentYFactors);
+    normalize (result);
     result.insert (CFFactor (LcF, 1));
     return result;
@@ -175,4 +343,7 @@
     A= F/LcF;
     result= multiplicity (A, buf);
+    for (CFFListIterator i= result; i.hasItem(); i++)
+      i.getItem()= CFFactor (N (decompress (i.getItem().factor(), M, S)),
+                             i.getItem().exp());
   }
   if (degree (A) > 0)
@@ -180,6 +351,12 @@
     resultRoot= FqBiFactorize (A, alpha);
     resultRoot.removeFirst();
+    for (CFFListIterator i= resultRoot; i.hasItem(); i++)
+      i.getItem()= CFFactor (N (decompress (i.getItem().factor(), M, S)),
+                             i.getItem().exp());
     result= Union (result, resultRoot);
   }
+  result= Union (result, contentXFactors);
+  result= Union (result, contentYFactors);
+  normalize (result);
   result.insert (CFFactor (LcF, 1));
   return result;
@@ -192,5 +369,5 @@
 /// @sa FpBiFactorize(), FqBiFactorize()
 inline
-CFFList GFBiFactorize (const CanonicalForm & F ///< [in] a bivariate poly
+CFFList GFBiFactorize (const CanonicalForm & G ///< [in] a bivariate poly
                       )
 {
@@ -199,9 +376,33 @@
   ExtensionInfo info= ExtensionInfo (getGFDegree(), gf_name, false);
   bool GF= true;
+  CFMap N;
+  CanonicalForm F= compress (G, N);
   CanonicalForm LcF= Lc (F);
+  CanonicalForm contentX= content (F, 1);
+  CanonicalForm contentY= content (F, 2);
+  F /= (contentX*contentY);
+  CFFList contentXFactors, contentYFactors;
+  contentXFactors= factorize (contentX);
+  contentYFactors= factorize (contentY);
+  if (contentXFactors.getFirst().factor().inCoeffDomain())
+    contentXFactors.removeFirst();
+  if (contentYFactors.getFirst().factor().inCoeffDomain())
+    contentYFactors.removeFirst();
+  decompress (contentXFactors, N);
+  decompress (contentYFactors, N);
+  CFFList result, resultRoot;
+  if (F.inCoeffDomain())
+  {
+    result= Union (contentXFactors, contentYFactors);
+    normalize (result);
+    result.insert (CFFactor (LcF, 1));
+    return result;
+  }
+  mat_ZZ M;
+  vec_ZZ S;
+  F= compress (F, M, S);
   CanonicalForm pthRoot, A;
   CanonicalForm sqrfP= sqrfPart (F/LcF, pthRoot, info.getAlpha());
   CFList buf;
-  CFFList result, resultRoot;
   int p= getCharacteristic();
   int q= ipower (p, getGFDegree());
@@ -213,5 +414,9 @@
     result.removeFirst();
     for (CFFListIterator i= result; i.hasItem(); i++)
-      i.getItem()= CFFactor(i.getItem().factor(),i.getItem().exp()*ipower(p,l));
+      i.getItem()= CFFactor (N (decompress (i.getItem().factor(), M, S)),
+                             i.getItem().exp()*ipower (p,l));
+    result= Union (result, contentXFactors);
+    result= Union (result, contentYFactors);
+    normalize (result);
     result.insert (CFFactor (LcF, 1));
     return result;
@@ -222,4 +427,7 @@
     A= F/LcF;
     result= multiplicity (A, buf);
+    for (CFFListIterator i= result; i.hasItem(); i++)
+      i.getItem()= CFFactor (N (decompress (i.getItem().factor(), M, S)),
+                             i.getItem().exp());
   }
   if (degree (A) > 0)
@@ -227,9 +435,17 @@
     resultRoot= GFBiFactorize (A);
     resultRoot.removeFirst();
+    for (CFFListIterator i= resultRoot; i.hasItem(); i++)
+      i.getItem()= CFFactor (N (decompress (i.getItem().factor(), M, S)),
+                             i.getItem().exp());
     result= Union (result, resultRoot);
   }
+  result= Union (result, contentXFactors);
+  result= Union (result, contentYFactors);
+  normalize (result);
   result.insert (CFFactor (LcF, 1));
   return result;
 }
+
+#endif
 
 /// \f$ \prod_{f\in L} {f (0, x)} \ mod\ M \f$ via divide-and-conquer
Index: factory/facFqBivarUtil.cc
===================================================================
--- factory/facFqBivarUtil.cc	(revision 0415f923fdb78a69504a400113e65cd371cb2150)
+++ factory/facFqBivarUtil.cc	(revision f876a663240524d1f86ed3519e45b68d1ae877b3)
@@ -25,4 +25,6 @@
 #include "cf_iter.h"
 #include "facFqBivarUtil.h"
+#include "cfNewtonPolygon.h"
+#include "facHensel.h"
 
 
@@ -41,5 +43,10 @@
   for (CFListIterator i= factors; i.hasItem(); i++)
     i.getItem()= N (i.getItem());
-  return;
+}
+
+void decompress (CFFList& factors, const CFMap& N)
+{
+  for (CFFListIterator i= factors; i.hasItem(); i++)
+    i.getItem()= CFFactor (N (i.getItem().factor()), i.getItem().exp());
 }
 
@@ -113,5 +120,7 @@
 
 static inline
-bool FqInExtensionHelper (const CanonicalForm& F, const CanonicalForm& gamma)
+bool FqInExtensionHelper (const CanonicalForm& F, const CanonicalForm& gamma,
+                          const CanonicalForm& delta, CFList& source,
+                          CFList& dest)
 {
   bool result= false;
@@ -123,5 +132,24 @@
       return true;
     else
-      return result;
+    {
+      int pos= findItem (source, F);
+      if (pos > 0)
+        return false;
+      Variable a;
+      hasFirstAlgVar (F, a);
+      int bound= ipower (getCharacteristic(), degree (getMipo (a)));
+      CanonicalForm buf= 1;
+      for (int i= 1; i < bound; i++)
+      {
+        buf *= gamma;
+        if (buf == F)
+        {
+          source.append (buf);
+          dest.append (power (delta, i));
+          return false;
+        }
+      }
+      return true;
+    }
   }
   else
@@ -129,5 +157,5 @@
     for (CFIterator i= F; i.hasTerms(); i++)
     {
-      result= FqInExtensionHelper (i.coeff(), gamma);
+      result= FqInExtensionHelper (i.coeff(), gamma, delta, source, dest);
       if (result == true)
         return result;
@@ -138,5 +166,6 @@
 
 bool isInExtension (const CanonicalForm& F, const CanonicalForm& gamma,
-                    const int k)
+                    const int k, const CanonicalForm& delta,
+                    CFList& source, CFList& dest)
 {
   bool result;
@@ -152,5 +181,5 @@
   else
   {
-    result= FqInExtensionHelper (F, gamma);
+    result= FqInExtensionHelper (F, gamma, delta, source, dest);
     return result;
   }
@@ -191,5 +220,5 @@
   if (k > 1)
   {
-    if (!isInExtension (g, delta, k))
+    if (!isInExtension (g, gamma, k, delta, source, dest))
     {
       g= GFMapDown (g, k);
@@ -199,5 +228,5 @@
   else if (k == 1)
   {
-    if (!isInExtension (g, delta, k));
+    if (!isInExtension (g, gamma, k, delta, source, dest));
       factors.append (g);
   }
@@ -209,5 +238,5 @@
   else if (!k && beta != Variable (1))
   {
-    if (!isInExtension (g, delta, k))
+    if (!isInExtension (g, gamma, k, delta, source, dest))
     {
       g= mapDown (g, delta, gamma, alpha, source, dest);
@@ -245,4 +274,12 @@
 }
 
+void normalize (CFFList& factors)
+{
+  for (CFFListIterator i= factors; i.hasItem(); i++) 
+    i.getItem()= CFFactor (i.getItem().factor()/Lc(i.getItem().factor()),
+                           i.getItem().exp());
+  return;
+}
+
 CFList subset (int index [], const int& s, const CFArray& elements,
                bool& noSubset)
@@ -392,4 +429,272 @@
 }
 
+CFArray
+logarithmicDerivative (const CanonicalForm& F, const CanonicalForm& G, int l,
+                       CanonicalForm& Q
+                      )
+{
+  Variable x= Variable (2);
+  Variable y= Variable (1);
+  CanonicalForm xToL= power (x, l);
+  CanonicalForm q,r;
+  CanonicalForm logDeriv;
+
+  q= newtonDiv (F, G, xToL);
+
+  logDeriv= mulMod2 (q, deriv (G, y), xToL);
+
+  logDeriv= swapvar (logDeriv, x, y);
+  int j= degree (logDeriv) + 1;
+  CFArray result= CFArray (j);
+  for (CFIterator i= logDeriv; i.hasTerms() && !logDeriv.isZero(); i++)
+  {
+    if (i.exp() == j - 1)
+    {
+      result [j - 1]= swapvar (i.coeff(), x, y);
+      j--;
+    }
+    else
+    {
+      for (; i.exp() < j - 1; j--)
+        result [j - 1]= 0;
+      result [j - 1]= swapvar (i.coeff(), x, y);
+      j--;
+    }
+  }
+  for (; j > 0; j--)
+    result [j - 1]= 0;
+  Q= q;
+  return result;
+}
+
+CFArray
+logarithmicDerivative (const CanonicalForm& F, const CanonicalForm& G, int l,
+                       int oldL, const CanonicalForm& oldQ, CanonicalForm& Q
+                      )
+{
+  Variable x= Variable (2);
+  Variable y= Variable (1);
+  CanonicalForm xToL= power (x, l);
+  CanonicalForm xToOldL= power (x, oldL);
+  CanonicalForm xToLOldL= power (x, l-oldL);
+  CanonicalForm q,r;
+  CanonicalForm logDeriv;
+
+  CanonicalForm bufF= mod (F, xToL);
+  CanonicalForm oldF= mulMod2 (G, oldQ, xToL);
+  bufF -= oldF;
+  bufF= div (bufF, xToOldL);
+
+  q= newtonDiv (bufF, G, xToLOldL);
+  q *= xToOldL;
+  q += oldQ;
+
+  logDeriv= mulMod2 (q, deriv (G, y), xToL);
+
+  logDeriv= swapvar (logDeriv, x, y);
+  int j= degree (logDeriv) + 1;
+  CFArray result= CFArray (j);
+  for (CFIterator i= logDeriv; i.hasTerms() && !logDeriv.isZero(); i++)
+  {
+    if (i.exp() == j - 1)
+    {
+      result [j - 1]= swapvar (i.coeff(), x, y);
+      j--;
+    }
+    else
+    {
+      for (; i.exp() < j - 1; j--)
+        result [j - 1]= 0;
+      result [j - 1]= swapvar (i.coeff(), x, y);
+      j--;
+    }
+  }
+  for (; j > 0; j--)
+    result [j - 1]= 0;
+  Q= q;
+  return result;
+}
+
+void
+writeInMatrix (CFMatrix& M, const CFArray& A, const int column,
+               const int startIndex
+              )
+{
+  ASSERT (A.size () - startIndex > 0, "wrong starting index");
+  ASSERT (A.size () - startIndex == M.rows(), "wrong starting index");
+  ASSERT (column > 0 && column <= M.columns(), "wrong column");
+  if (A.size() - startIndex <= 0) return;
+  int j= 1;
+  for (int i= startIndex; i < A.size(); i++, j++)
+    M (j, column)= A [i];
+}
+
+CFArray getCoeffs (const CanonicalForm& F, const int k)
+{
+  ASSERT (F.isUnivariate(), "univariate input expected");
+  if (degree (F, 2) < k)
+    return CFArray();
+
+  CFArray result= CFArray (degree (F) - k + 1);
+  CFIterator j= F;
+  for (int i= degree (F); i >= k; i--)
+  {
+    if (j.exp() == i)
+    {
+      result [i - k]= j.coeff();
+      j++; 
+      if (!j.hasTerms())
+        return result;
+    }
+    else
+      result[i - k]= 0;
+  }
+  return result;
+}
+
+CFArray getCoeffs (const CanonicalForm& F, const int k, const Variable& alpha)
+{
+  ASSERT (F.isUnivariate(), "univariate input expected");
+  if (degree (F, 2) < k)
+    return CFArray ();
+
+  int d= degree (getMipo (alpha));
+  CFArray result= CFArray ((degree (F) - k + 1)*d);
+  CFIterator j= F;
+  CanonicalForm buf;
+  CFIterator iter;
+  for (int i= degree (F); i >= k; i--)
+  {
+    if (j.exp() == i)
+    {
+      iter= j.coeff();
+      for (int l= degree (j.coeff(), alpha); l >= 0; l--)
+      {
+        if (iter.exp() == l)
+        {
+          result [(i - k)*d + l]= iter.coeff();
+          iter++;
+          if (!iter.hasTerms())
+            break;
+        }
+      }
+      j++;
+      if (!j.hasTerms())
+        return result;
+    }
+    else
+    {
+      for (int l= 0; l < d; l++)
+        result[(i - k)*d + l]= 0;
+    }
+  }
+  return result;
+}
+
+#ifdef HAVE_NTL
+CFArray
+getCoeffs (const CanonicalForm& G, const int k, const int l, const int degMipo,
+           const Variable& alpha, const CanonicalForm& evaluation,
+           const mat_zz_p& M)
+{
+  ASSERT (G.isUnivariate(), "univariate input expected");
+  CanonicalForm F= G (G.mvar() - evaluation, G.mvar());
+  if (F.isZero())
+    return CFArray ();
+
+  Variable y= Variable (2);
+  F= F (power (y, degMipo), y);
+  F= F (y, alpha);
+  zz_pX NTLF= convertFacCF2NTLzzpX (F);
+  NTLF.rep.SetLength (l*degMipo);
+  NTLF.rep= M*NTLF.rep;
+  NTLF.normalize();
+  F= convertNTLzzpX2CF (NTLF, y);
+  int d= degMipo;
+
+  if (degree (F, 2) < k)
+    return CFArray();
+
+  CFArray result= CFArray (degree (F) - k + 1);
+
+  CFIterator j= F;
+  for (int i= degree (F); i >= k; i--)
+  {
+    if (j.exp() == i)
+    {
+      result [i - k]= j.coeff();
+      j++; 
+      if (!j.hasTerms())
+        return result;
+    }
+    else
+      result[i - k]= 0;
+  }
+  return result;
+}
+#endif
+
+int * computeBounds (const CanonicalForm& F, int& n)
+{
+  n= degree (F, 1);
+  int* result= new int [n];
+  int sizeOfNewtonPolygon;
+  int** newtonPolyg= newtonPolygon (F, sizeOfNewtonPolygon);
+
+  int minXIndex= 0, minYIndex= 0, maxXIndex= 0, maxYIndex= 0;
+  int minX, minY, maxX, maxY;
+  minX= newtonPolyg [0] [0];
+  minY= newtonPolyg [0] [1];
+  maxX= minX;
+  maxY= minY;
+  for (int i= 1; i < sizeOfNewtonPolygon; i++)
+  {
+    if (minX > newtonPolyg [i] [0])
+    {
+      minX= newtonPolyg [i] [0];
+      minXIndex= i;
+    }
+    if (maxX < newtonPolyg [i] [0])
+    {
+      maxX= newtonPolyg [i] [0];
+      maxXIndex= i;
+    }
+    if (minY > newtonPolyg [i] [1])
+    {
+      minY= newtonPolyg [i] [1];
+      minYIndex= i;
+    }
+    if (maxY < newtonPolyg [i] [1])
+    {
+      maxY= newtonPolyg [i] [1];
+      maxYIndex= i;
+    }
+  }
+
+  int k= maxX;
+  bool hitMaxX= false;
+  for (int i= 0; i < n; i++)
+  {
+    if (i + 1 > maxY || i + 1 < minY)
+    {
+      result [i]= 0;
+      continue;
+    }
+    int* point= new int [2];
+    point [0]= k;
+    point [1]= i + 1;
+    while (!isInPolygon (newtonPolyg, sizeOfNewtonPolygon, point) && k > 0)
+    {
+      k--;
+      point [0]= k;
+    }
+    result [i]= k;
+    k= maxX;
+    delete [] point;
+  }
+
+  return result;
+}
+
 int
 substituteCheck (const CanonicalForm& F, const Variable& x)
Index: factory/facFqBivarUtil.h
===================================================================
--- factory/facFqBivarUtil.h	(revision 0415f923fdb78a69504a400113e65cd371cb2150)
+++ factory/facFqBivarUtil.h	(revision f876a663240524d1f86ed3519e45b68d1ae877b3)
@@ -21,4 +21,8 @@
 #include "ExtensionInfo.h"
 
+#ifdef HAVE_NTL
+#include "NTLconvert.h"
+#endif
+
 /// append @a factors2 on @a factors1
 void append (CFList& factors1,       ///< [in,out] a list of polys
@@ -29,4 +33,9 @@
 void decompress (CFList& factors, ///< [in,out] a list of polys
                  const CFMap& N   ///< [in] a map
+                );
+
+/// as above
+void decompress (CFFList& factors, ///< [in,out] a list of factors
+                 const CFMap& N    ///< [in] a map
                 );
 
@@ -61,7 +70,10 @@
                                                  ///< Fq if we are not over some
                                                  ///< GF
-                    const int k                  ///< some int k such that k
+                    const int k,                 ///< some int k such that k 
                                                  ///< divides l if we are not
                                                  ///< over some Fq
+                    const CanonicalForm& delta,  ///< [in] image of gamma
+                    CFList& source,              ///< [in,out] list of preimages
+                    CFList& dest                 ///< [in,out] list of images
                    );
 
@@ -111,4 +123,8 @@
                );
 
+/// as above
+void normalize (CFFList& factors ///< [in,out] a list of factors
+               );
+
 /// extract a subset given by @a index of size @a s from @a elements, if there
 /// is no subset we have not yet considered noSubset is set to true. @a index
@@ -158,4 +174,92 @@
                      );
 
+/// compute the coefficients of the logarithmic derivative of G mod
+/// Variable (2)^l over Fq
+///
+/// @return an array of coefficients of the logarithmic derivative of G mod
+///         Variable (2)^l
+CFArray logarithmicDerivative (const CanonicalForm& F,///<[in] a bivariate poly
+                              const CanonicalForm& G, ///<[in] a factor of F
+                              int l,                  ///<[in] lifting precision
+                              CanonicalForm& Q        ///<[in,out] F/G
+                              );
+
+/// compute the coefficients of the logarithmic derivative of G mod
+/// Variable (2)^l over Fq with oldQ= F/G mod Variable (2)^oldL
+///
+/// @return an array of coefficients of the logarithmic derivative of G mod
+///         Variable (2)^l
+CFArray
+logarithmicDerivative (const CanonicalForm& F,   ///< [in] bivariate poly
+                       const CanonicalForm& G,   ///< [in] a factor of F
+                       int l,                    ///< [in] lifting precision
+                       int oldL,                 ///< [in] old precision
+                       const CanonicalForm& oldQ,///< [in] F/G mod
+                                                 ///< Variable(2)^oldL
+                       CanonicalForm& Q          ///< [in, out] F/G
+                      );
+
+/// compute bounds for logarithmic derivative as described in K. Belabas, 
+/// M. van Hoeij, J. KlÃŒners, and A. Steel, Factoring polynomials over global 
+/// fields
+///
+/// @return @a computeBounds returns bounds as described above
+int *
+computeBounds (const CanonicalForm& F,///< [in] compressed bivariate polynomial
+               int& n                 ///< [in,out] length of output
+              );
+
+/// extract coefficients of \f$ x^i \f$ for \f$i\geq k\f$ where \f$ x \f$ is
+/// a variable of level 1
+///
+/// @return coefficients of \f$ x^i \f$ for \f$i\geq k\f$
+/// @sa {getCoeffs (const CanonicalForm&, const int, const Variable&),
+/// getCoeffs (const CanonicalForm&, const int, const int, const int, 
+/// const Variable&, const CanonicalForm&, const mat_zz_p&)}
+CFArray
+getCoeffs (const CanonicalForm& F,///< [in] compressed bivariate poly over F_p
+           const int k            ///< [in] some int
+          );
+
+/// extract coefficients of \f$ x^i \f$ for \f$i\geq k\f$ where \f$ x \f$ is
+/// a variable of level 1
+///
+/// @return coefficients of \f$ x^i \f$ for \f$i\geq k\f$
+/// @sa {getCoeffs (const CanonicalForm&, const int),
+/// getCoeffs (const CanonicalForm&, const int, const int, const int, 
+/// const Variable&, const CanonicalForm&, const mat_zz_p&)}
+CFArray
+getCoeffs (const CanonicalForm& F,///< [in] compressed bivariate poly over 
+                                  ///< F_p(alpha)
+           const int k,           ///< [in] some int
+           const Variable& alpha  ///< [in] algebraic variable
+          );
+
+#ifdef HAVE_NTL
+/// extract coefficients of \f$ x^i \f$ for \f$i\geq k\f$ where \f$ x \f$ is
+/// a variable of level 1
+///
+/// @return coefficients of \f$ x^i \f$ for \f$i\geq k\f$
+/// @sa {getCoeffs (const CanonicalForm&, const int, const Variable& alpha),
+/// getCoeffs (const CanonicalForm&, const int)}
+CFArray
+getCoeffs (const CanonicalForm& F,         ///< [in] compressed bivariate poly
+           const int k,                    ///< [in] some int
+           const int l,                    ///< [in] precision
+           const int degMipo,              ///< [in] degree of minimal poly
+           const Variable& alpha,          ///< [in] algebraic variable
+           const CanonicalForm& evaluation,///< [in] evaluation point
+           const mat_zz_p& M               ///< [in] bases change matrix
+          );
+#endif
+
+/// write A into M starting at row startIndex
+void
+writeInMatrix (CFMatrix& M,        ///< [in,out] some matrix
+               const CFArray& A,   ///< [in] array of polys
+               const int column,   ///< [in] column in which A is written
+               const int startIndex///< [in] row in which to start
+              );
+
 /// checks if a substitution x^n->x is possible
 ///
Index: factory/facFqFactorize.cc
===================================================================
--- factory/facFqFactorize.cc	(revision 0415f923fdb78a69504a400113e65cd371cb2150)
+++ factory/facFqFactorize.cc	(revision f876a663240524d1f86ed3519e45b68d1ae877b3)
@@ -696,5 +696,5 @@
         else
         {
-          if (!isInExtension (buf2, delta, k))
+          if (!isInExtension (buf2, gamma, k, delta, source, dest))
           {
             appendTestMapDown (result, buf2, info, source, dest);
@@ -924,4 +924,5 @@
     degMipoBeta= degree (getMipo (beta));
 
+  CFList source, dest;
   for (CFListIterator i= factors; i.hasItem(); i++)
   {
@@ -945,5 +946,5 @@
       else
       {
-        if (!isInExtension (gg, delta, k))
+        if (!isInExtension (gg, gamma, k, delta, source, dest))
         {
           buf /= g;
@@ -1096,5 +1097,5 @@
       else
       {
-        if (!isInExtension (gg, delta, k))
+        if (!isInExtension (gg, gamma, k, delta, source, dest))
         {
           appendTestMapDown (result, gg, info, source, dest);
