val P = P * p

               val P = P * p

               val g = centr_up ((chinese_remainder_up (P, p) (g, g_p)) mod_up P) P

               val g = centr_up ((chinese_remainder_up (P, p) (g, g_p)) mod_up P) P

             in try_new_prime_up a b d M P g p end

             in try_new_prime_up a b d M P g p end

       end

       end

        assumes d = gcd (lcoeff_up p1, lcoeff_up p2) \<and> 

        assumes d = gcd (lcoeff_up p1, lcoeff_up p2) \<and> 

                M = LANDAU_MIGNOTTE_bound \<and> p = prime  \<and>  p ~| d  \<and>

                M = LANDAU_MIGNOTTE_bound \<and> p = prime  \<and>  p ~| d  \<and>

                p1 is primitive  \<and>  p2 is primitive

                p1 is primitive  \<and>  p2 is primitive

        yields  gcd | p1  \<and>  gcd | p2  \<and>  gcd is primitive

        yields  gcd | p1  \<and>  gcd | p2  \<and>  gcd is primitive

     argumentns "a b d M p" named according to [1] p.93 *)

     argumentns "a b d M p" named according to [1] p.93 *)

     let

     let

       val p = p next_prime_not_dvd d 

       val p = p next_prime_not_dvd d 

       val g_p = centr_up (((normalise (mod_up_gcd a b p) p) %* (d mod p)) mod_up p) p 

       val g_p = centr_up (((normalise (mod_up_gcd a b p) p) %* (d mod p)) mod_up p) p 

         (*.. nesting of functions see above*)

         (*.. nesting of functions see above*)

     in

     in

       if deg_up g_p = 0 then [1] else

       if deg_up g_p = 0 then [1] else

         let 

         let 

           val g = primitive_up (try_new_prime_up a b d M p g_p p)

           val g = primitive_up (try_new_prime_up a b d M p g_p p)

         in

         in

         end

         end

     end

     end

   (* gcd_up :: unipoly \<Rightarrow> unipoly \<Rightarrow> unipoly

   (* gcd_up :: unipoly \<Rightarrow> unipoly \<Rightarrow> unipoly

      gcd_up a b = c

      gcd_up a b = c

      yields c dvd a \<and> c dvd b \<and> (\<forall>c'. (c' dvd a \<and> c' dvd b) \<longrightarrow> c' \<le> c) *)

      yields c dvd a \<and> c dvd b \<and> (\<forall>c'. (c' dvd a \<and> c' dvd b) \<longrightarrow> c' \<le> c) *)

   fun gcd_up a b =

   fun gcd_up a b =

     let val d = abs (Integer.gcd (lcoeff_up a) (lcoeff_up b))

     let val d = abs (Integer.gcd (lcoeff_up a) (lcoeff_up b))

     in 

     in 

       if a = b then a else

       if a = b then a else

     end;

     end;

 *}

 *}

 section {* gcd for multivariate polynomials *}

 section {* gcd for multivariate polynomials *}

 ML {*

 ML {*

             (* a b are made primitive and this is compensated by multiplying with gcd of divisors *)

             (* a b are made primitive and this is compensated by multiplying with gcd of divisors *)

             uni_to_multi ((gcd_up (primitive_up a') (primitive_up b')) %* 

             uni_to_multi ((gcd_up (primitive_up a') (primitive_up b')) %* 

                           (abs (Integer.gcd (gcd_coeff a') (gcd_coeff b'))))

                           (abs (Integer.gcd (gcd_coeff a') (gcd_coeff b'))))

           end

           end

      assumes length a > 1  \<and>  length b > 1

      assumes length a > 1  \<and>  length b > 1

   *)

   *)

   and try_new_r a b n M r r_list steps = 

   and try_new_r a b n M r r_list steps = 

      let 

      let 

         val r = find_new_r a b r;

         val r = find_new_r a b r;

         val r_list = r_list @ [r];

         val r_list = r_list @ [r];

         val g_r = gcd_poly' (order (eval a (n - 2) r)) 

         val g_r = gcd_poly' (order (eval a (n - 2) r)) 

                             (order (eval b (n - 2) r)) (n - 1) 0

                             (order (eval b (n - 2) r)) (n - 1) 0

         val steps = steps @ [g_r];

         val steps = steps @ [g_r];

      assumes length a > 1  \<and>  length b > 1

      assumes length a > 1  \<and>  length b > 1

   *)

   *)

     if m > M then if g_n %%|%% a andalso g_n %%|%% b then g_n else

     if m > M then if g_n %%|%% a andalso g_n %%|%% b then g_n else

     else

     else

       let 

       let 

         val r = find_new_r a b r; 

         val r = find_new_r a b r; 

         val r_list = r_list @ [r];

         val r_list = r_list @ [r];

         val g_r = gcd_poly' (order (eval a (n - 2) r)) 

         val g_r = gcd_poly' (order (eval a (n - 2) r)) 

                             (order (eval b (n - 2) r)) (n - 1) 0

                             (order (eval b (n - 2) r)) (n - 1) 0

       in  

       in  

         if lex_ord (lmonom g) (lmonom g_r)

         if lex_ord (lmonom g) (lmonom g_r)

           then

           then

             let val (g_n, new, steps) = NEWTON r_list [g, g_r] steps mult g_n (n - 2)

             let val (g_n, new, steps) = NEWTON r_list [g, g_r] steps mult g_n (n - 2)

             in 

             in 

             end

             end

           else (* \<not> lex_ord (lmonom g) (lmonom g_r) *)

           else (* \<not> lex_ord (lmonom g) (lmonom g_r) *)

       end

       end

     let val cs = (coeffs a) @ (coeffs b) @ (coeffs c)

     let val cs = (coeffs a) @ (coeffs b) @ (coeffs c)

     in length (filter (curry op < 0) cs) < length cs div 2 end

     in length (filter (curry op < 0) cs) < length cs div 2 end

   (* gcd_poly :: poly \<Rightarrow> poly \<Rightarrow> poly

   (* gcd_poly :: poly \<Rightarrow> poly \<Rightarrow> poly

      gcd_poly a b = ((a', b'), c)

      gcd_poly a b = ((a', b'), c)

   fun gcd_poly (a as (_, es)::_ : poly) b = 

   fun gcd_poly (a as (_, es)::_ : poly) b = 

     let val c = gcd_poly' a b (length es) 0

     let val c = gcd_poly' a b (length es) 0

         val (a': poly, _) = a %%/%% c

         val (a': poly, _) = a %%/%% c

         val (b': poly, _) = b %%/%% c

         val (b': poly, _) = b %%/%% c

     in 

     in 

       then ((a' %%* ~1, b' %%* ~1), c %%* ~1) (*makes results look nicer*)

       then ((a' %%* ~1, b' %%* ~1), c %%* ~1) (*makes results look nicer*)

       else ((a', b'), c)

       else ((a', b'), c)

     end

     end

 *}

 *}