(* Title: test/../rational2

   Author: Diana Meindl

   Copyright (c) Diana Meindl 2011

12345678901234567890123456789012345678901234567890123456789012345678901234567890

        10        20        30        40        50        60        70        80

*)

(*fun nth _ []      = error "nth _ []" (*Isabelle2002, still saved the hours of update*)

    | nth 1 (x::_) = x

    | nth n (_::xs) = nth (n-1) xs;*)

  fun nth xs i = List.nth (xs, i);     (*recent Isabelle: TODO update all isac code   *)

"--------------------------------------------------------";

"table of contents --------------------------------------";

"--------------------------------------------------------";

"----------- fun div2 -- --------------------------------";

"----------- fun mod_inv --------------------------------";

"----------- fun mod_div --------------------------------";

"----------- fun chinese_remainder ----------------------";

"----------- fun is_prime -------------------------------";

"----------- fun make_primes ----------------------------";

"----------- fun primes_upto ----------------------------";

"----------- fun next_prime_not_dvd ---------------------";

"----------- fun lcoeff_up ------------------------------";

"----------- fun deg_up ---------------------------------";

"----------- fun drop_lc0_up ----------------------------";

"----------- fun normalise ------------------------------";

(* order until here: (GCD_Poly_FP =) GCD_Poly = gcd_poly *)

"----------- fun %* -------------------------------------";

"----------- fun %/ -------------------------------------";

"----------- fun %-% ------------------------------------";

"----------- fun %|% ------------------------------------";

(*"----------- fun %+% ------------------------------------";*)

(*"----------- fun centr_up -----------------------------";*)

"----------- fun LANDAU_MIGNOTTE_bound ------------------";

"----------- fun chinese_remainder_up -------------------";

"----------- fun mod_up ---------------------------------";

"----------- fun mod_up_gcd -----------------------------";

"----------- fun primitive_up ---------------------------";

"----------- fun try_new_prime_up -----------------------";

"----------- fun gcd_up ---------------------------------";

"=========== Multivariate Case ==========================";

"----------- fun lcoeff ---------------------------------";

"----------- fun lexp -----------------------------------";

"----------- fun add_variable ---------------------------";

"----------- fun zero_poly ------------------------------";

"----------- drop_lc0 -----------------------------------";

"----------- fun add_monoms -----------------------------";

"----------- fun lex_ord --------------------------------";

"----------- fun order ----------------------------------";

(*"----------- fun lcoeff_up_poly ---------------------------";*)

"----------- fun lexp -----------------------------------";

"----------- fun max_deg_var ----------------------------";

"----------- infix %%/ ----------------------------------";

"----------- infix %%*%% --------------------------------";

"----------- fun gcd_coeff ------------------------------";

(*"----------- fun gcd_coeff_poly -------------------------";*)

"----------- fun eval -----------------------------------";

"----------- fun multi_to_uni ---------------------------";

"----------- fun uni_to_multi ---------------------------";

"----------- fun find_new_r  ----------------------------";

"----------- fun mult_with_new_var  ---------------------";

(*"----------- fun greater_deg  -----------------";*)

"----------- infix %%+%%  -------------------------------";

"----------- infix %%-%%  -------------------------------";

"----------- infix %%*%%  -------------------------------";

"----------- infix %%/%%  -------------------------------";

"----------- fun newton_first  --------------------------";

"----------- fun NEWTON  --------------------------------";

"----------- fun all_geq  -------------------------------";

"----------- fun greater_deg  ---------------------------";

"----------- infix %%|%% --------------------------------";

"----------- fun gcd_poly'  -----------------------------";

"----------- fun gcd_poly  ------------------------------";

"----------- fun gcd_poly downto division by zero -------";

"=========== prep Lucas-Interpretation ==================";

"----------- fun for_quot_regarding ---------------------";

"----------- fun mult_to_deg - --------------------------";

"----------- fun fact_quot ------------------------------";

"----------- fun %+% ------------------------------------";

"----------- fun %*% ------------------------------------";

"----------- fun %*/% -----------------------------------";

"----------- fun %*/% Subscript raised (line 509 of lib--";

"----------- fun EUCLID_naive_up ------------------------";

"----------- last step in EUCLID ------------------------";

" ========== END ======================================= ";

"----------- fun div2 -- --------------------------------";

"----------- fun div2 -- --------------------------------";

"----------- fun div2 -- --------------------------------";

if 5 div 2 = 2 then () else error "5 div 2 = 2 changed";

if (5 div2 2) = 2 then () else error "5 div2 2 = 2 changed";

if ~5 div 2 = ~3 then () else error "~5 div 2 = ~3 changed";

if (~5 div2 2) = ~2 then () else error "~5 div2 2 = ~2 changed";

if ~5 div ~2 = 2 then () else error "~5 div ~2 = 2 changed";

if (~5 div2 ~2) = 2 then () else error "~5 div2 ~2 = 2 changed";

if 5 div ~2 = ~3 then () else error "5 div ~2 = ~3 changed";

if (5 div2 ~2) = ~2 then () else error "5 div2 ~2 = ~2 changed";

"----------- fun mod_inv --------------------------------";

"----------- fun mod_inv --------------------------------";

"----------- fun mod_inv--- -----------------------------";

"~~~~~ fun mod_inv, args:"; val (r, m) = (5,7);

"~~~~~ fun modi, args:"; val (r, rold, m, rinv) = (r, r, m, 1);

rinv < m; (*=true*)

r mod m = 1; (*=false*)

"~~~~~ fun modi, args:"; val (r, rold, m, rinv) = (rold * (rinv + 1), rold, m, rinv + 1);

rinv < m;(*=true*)

r mod m = 1; (*=false*)

"~~~~~ fun modi, args:"; val (r, rold, m, rinv) = (rold * (rinv + 1), rold, m, rinv + 1);

rinv < m;(*=true*)

r mod m = 1; (*=true*)

"~~~~~ to mod_inv return val:"; val (rinv) = (rinv);

if mod_inv 5 7 = 3 then () else error "mod_inv 5 7 = 3 changed";

if mod_inv 3 7 = 5 then () else error "mod_inv 3 7 = 5 changed";

if mod_inv 4 339 = 85 then () else error "mod_inv 4 339 = 85 changed";

"----------- fun mod_div --------------------------------";

"----------- fun mod_div --------------------------------";

"----------- fun mod_div --------------------------------";

if mod_div 21 3 5 = 2 then () else error "mod_div changed";  21 mod 5 = (3 * 2) mod 5;

if mod_div 22 3 5 = 4 then () else error "mod_div changed";  22 mod 5 = (3 * 4) mod 5;

if mod_div 23 3 5 = 1 then () else error "mod_div changed";  23 mod 5 = (3 * 1) mod 5;

if mod_div 24 3 5 = 3 then () else error "mod_div changed";  24 mod 5 = (3 * 3) mod 5;

if mod_div 25 3 5 = 0 then () else error "mod_div changed";  25 mod 5 = (3 * 0) mod 5;

if mod_div 21 4 5 = 4 then () else error "mod_div 21 4 5 = 4 changed";

if mod_div 1 4 5 = 4 then () else error "mod_div 1 4 5 = 4 changed";

if mod_div 0 4 5 = 0 then () else error "mod_div 0 4 5 = 0 changed";

"----------- fun chinese_remainder ----------------------";

"----------- fun chinese_remainder ----------------------";

"----------- fun chinese_remainder ----------------------";

if  chinese_remainder(17,9) (3,4) = 35 then () else error "chinese_remainder(17,9)(3,4)=35 changed"; 

if  chinese_remainder(7,2) (6,11) = 17 then () else error "chinese_remainder(7,2)(6,11)=17 changed";  

"----------- fun is_prime -------------------------------";

"----------- fun is_prime -------------------------------";

"----------- fun is_prime -------------------------------";

if is_prime [2, 3] 2  = false then () else error "is_prime changed";

if is_prime [2, 3] 3  = false then () else error "is_prime changed";

if is_prime [2, 3] 4  = false then () else error "is_prime changed";

if is_prime [2, 3] 5  = true then () else error "is_prime changed";

if is_prime [2, 3, 5] 5  = false then () else error "is_prime changed";

if is_prime [2, 3] 6  = false then () else error "is_prime changed";

if is_prime [2, 3] 7  = true then () else error "is_prime changed";

if is_prime [2, 3] 25  = true then () else error "is_prime changed";

"----------- fun make_primes ----------------------------";

"----------- fun make_primes ----------------------------";

"----------- fun make_primes ----------------------------";

if make_primes [2, 3] 3 4 = [2, 3, 5] then () else error "make_primes changed";

if make_primes [2, 3] 3 5 = [2, 3, 5] then () else error "make_primes changed";

if make_primes [2, 3] 3 6 = [2, 3, 5, 7] then () else error "make_primes changed";

if make_primes [2, 3] 3 7 = [2, 3, 5, 7] then () else error "make_primes changed";

if make_primes [2, 3] 3 8 = [2, 3, 5, 7, 11] then () else error "make_primes changed";

if make_primes [2, 3] 3 9 = [2, 3, 5, 7, 11] then () else error "make_primes changed";

if make_primes [2, 3] 3 10 = [2, 3, 5, 7, 11] then () else error "make_primes changed";

if make_primes [2, 3] 3 11 = [2, 3, 5, 7, 11] then () else error "make_primes changed";

if make_primes [2, 3] 3 12 = [2, 3, 5, 7, 11, 13] then () else error "make_primes changed";

if make_primes [2, 3] 3 13 = [2, 3, 5, 7, 11, 13] then () else error "make_primes changed";

if make_primes [2, 3] 3 14 = [2, 3, 5, 7, 11, 13, 17] then () else error "make_primes changed";

if make_primes [2, 3] 3 15 = [2, 3, 5, 7, 11, 13, 17] then () else error "make_primes changed";

if make_primes [2, 3] 3 16 = [2, 3, 5, 7, 11, 13, 17] then () else error "make_primes changed";

if make_primes [2, 3] 3 17 = [2, 3, 5, 7, 11, 13, 17] then () else error "make_primes changed";

if make_primes [2, 3] 3 18 = [2, 3, 5, 7, 11, 13, 17, 19] then () else error "make_primes changed";

if make_primes [2, 3] 3 19 = [2, 3, 5, 7, 11, 13, 17, 19] then () else error "make_primes changed";

if make_primes [2, 3, 5, 7] 7 4 = [2, 3, 5, 7] then () else error "make_primes changed";

"----------- fun primes_upto ----------------------------";

"----------- fun primes_upto ----------------------------";

"----------- fun primes_upto ----------------------------";

if primes_upto 1 = [2] then () else error " primes_upto 1 = [2] changed";

if primes_upto 3 = [2,3] then () else error " primes_upto 3 = [2,3] changed";

if primes_upto 6 = [2,3,5,7] then () else error " primes_upto 6 = [2,3,5,7] changed";

if primes_upto 7 = [2,3,5,7] then () else error " primes_upto 7 = [2,3,5,7] changed";

"----------- fun next_prime_not_dvd ----------------";

"----------- fun next_prime_not_dvd ----------------";

"----------- fun next_prime_not_dvd ----------------";

if ( 1 next_prime_not_dvd 15) =  2 then () else error "next_prime_not_dvd ..changed";

if ( 2 next_prime_not_dvd 15) =  7 then () else error "next_prime_not_dvd ..changed";

if ( 3 next_prime_not_dvd 15) =  7 then () else error "next_prime_not_dvd ..changed";

if ( 4 next_prime_not_dvd 15) =  7 then () else error "next_prime_not_dvd ..changed";

if ( 5 next_prime_not_dvd 15) =  7 then () else error "next_prime_not_dvd ..changed";

if ( 6 next_prime_not_dvd 15) =  7 then () else error "next_prime_not_dvd ..changed";

if ( 7 next_prime_not_dvd 15) = 11 then () else error "next_prime_not_dvd ..changed";

"----------- fun lcoeff_up ------------------------------";

"----------- fun lcoeff_up ------------------------------";

"----------- fun lcoeff_up ------------------------------";

if lcoeff_up [3,4,5,6] = 6 then () else error "lcoeff_up (3,4,5,6) = 6 changed" ;

if lcoeff_up [3,4,5,6,0] = 6 then () else error "lcoeff_up (3,4,5,6,0) = 6 changed"  ;

"----------- fun deg_up ---------------------------------";

"----------- fun deg_up ---------------------------------";

"----------- fun deg_up ---------------------------------";

if deg_up [3,4,5,6] = 3 then () else error "deg_up [3,4,5,6] = 3 changed" ;

if deg_up [3,4,5,6,0] = 3 then () else error "deg_up [3,4,5,6,0] = 3 changed";

"----------- fun drop_lc0_up ----------------------------";

"----------- fun drop_lc0_up ----------------------------";

"----------- fun drop_lc0_up ----------------------------";

if drop_lc0_up [0, 1, 2, 3, 4, 5, 0, 0] = [0, 1, 2, 3, 4, 5] then ()

else error "drop_lc0_up [0, 1, 2, 3, 4, 5, 0, 0] = [0, 1, 2, 3, 4, 5] changed";

if drop_lc0_up [0, 1, 2, 3, 4, 5] = [0, 1, 2, 3, 4, 5] then () 

else error "drop_lc0_up [0, 1, 2, 3, 4, 5]=[0, 1, 2, 3, 4, 5] changed";

"----------- fun normalise ------------------------------";

"----------- fun normalise ------------------------------";

"----------- fun normalise ------------------------------";

if (normalise [~18, ~15, ~20, 12, 20, ~13, 2] 5) = [1, 0, 0, 1, 0, 1, 1 ] then ()

else error "(normalise [~18, ~15, ~20, 12, 20, ~13, 2] 5) = [1, 0, 0, 1, 0, 1, 1 ] changed"

"----------- fun %* ------------------------------";

"----------- fun %* ------------------------------";

"----------- fun %* ------------------------------";

if ([5, 4, 7, 8, 1] %* 5) = [25, 20, 35, 40, 5] then () 

else error "[5, 4, 7, 8, 1] %* 5 = [25, 20, 35, 40, 5] changed";

if ([5, 4, ~7, 8, ~1] %* 5) = [25, 20, ~35, 40, ~5] then () 

else error "[5, 4, ~7, 8, ~1] %* 5 = [25, 20, ~35, 40, ~5] changed";

"----------- fun %/ -------------------------------";

"----------- fun %/ -------------------------------";

"----------- fun %/ -------------------------------";

if ([4, 3, 2, 5, 6] %/ 3) = [1, 1, 0, 1, 2] then ()

else error "%/ [4, 3, 2, 5, 6] 3 = [1, 1, 0, 1, 2] changed";

if ([4, 3, 2, 0] %/ 3) = [1, 1, 0, 0] then () 

else error "%/ [4, 3, 2, 0] 3 = [1, 1, 0, 0] changed";

"----------- fun %-% ------------------------------";

"----------- fun %-% ------------------------------";

"----------- fun %-% -----------------------------";

if ([8, ~7, 0, 1] %-% [~2, 2, 3, 0]) = [10, ~9, ~3, 1] then () 

else error "([8, ~7, 0, 1] %-% [~2, 2, 3, 0]) = [10, ~9, ~3, 1] changed";

if ([8, 7, 6, 5, 4] %-% [2, 2, 3, 1, 1]) = [6, 5, 3, 4, 3] then ()

else error "([8, 7, 6, 5, 4] %-% [2, 2, 3, 1, 1]) = [6, 5, 3, 4, 3] changed";

"----------- fun %|% -------------------------------";

"----------- fun %|% -------------------------------";

"----------- fun %|% -------------------------------";

if ([2, 3] %|% [8, 22, 15]) = true then () else error "[2, 3] %|% [8, 22, 15] = true changed";

if ([4] %|% [6]) = false then () else error "[4] %|% [6] = false changed";

if ([12] %|% [6]) = false then () else error "[12] %|% [6] = false changed";

if ([8] %|% [16, 0]) = true then () else error "[8] %|% [16, 0] = true changed";

if ([3, 2] %|% [0, 0, 9, 12, 4] ) = true then () else error "[3, 2] %|% [0,0,9..]  = true changed";

if ([8, 0] %|% [16]) = true then () else error "[8,0] %|% [16] = true changed";

"----------- fun centr_up -----------------------------";

"----------- fun centr_up -----------------------------";

"----------- fun centr_up -----------------------------";

if centr_up [7, 3, 5, 8, 1, 3] 10 = [~3, 3, 5, ~2, 1, 3] then ()

else error "centr_up [7, 3, 5, 8, 1, 3] 10 = [~3, 3, 5, ~2, 1, 3] changed";

if centr_up [1, 2, 3, 4, 5] 2  = [1, 0, 1, 2, 3] then () else error "centr_up [1,.] 2 changed";

if centr_up [1, 2, 3, 4, 5] 3  = [1, ~1, 0, 1, 2] then () else error "centr_up [1,.] 3 changed";

if centr_up [1, 2, 3, 4, 5] 4  = [1, 2, ~1, 0, 1] then () else error "centr_up [1,.] 4 changed";

if centr_up [1, 2, 3, 4, 5] 5  = [1, 2, 3, ~1, 0] then () else error "centr_up [1,.] 5 changed";

if centr_up [1, 2, 3, 4, 5] 6  = [1, 2, 3, ~2, ~1] then () else error "centr_up [1,.]6 changed";

if centr_up [1, 2, 3, 4, 5] 7  = [1, 2, 3, 4, ~2] then () else error "centr_up [1,.] 7 changed";

if centr_up [1, 2, 3, 4, 5] 8  = [1, 2, 3, 4, ~3] then () else error "centr_up [1,.] 8 changed";

if centr_up [1, 2, 3, 4, 5] 9  = [1, 2, 3, 4, 5] then () else error "centr_up [1,.] 9 changed";

if centr_up [1, 2, 3, 4, 5] 10 = [1, 2, 3, 4, 5] then () else error "centr_up [1,.] 10 changed";

if centr_up [1, 2, 3, 4, 5] 11 = [1, 2, 3, 4, 5] then () else error "centr_up [1,.] 11 changed";

"----------- fun sum_lmb --------------------------------";

"----------- fun sum_lmb --------------------------------";

"----------- fun sum_lmb --------------------------------";

if sum_lmb [1, 2, 3, 4, 5] 1 = 15 then () else error "sum_lmb [1, 2, 3, 4, 5] 1 changed";

if sum_lmb [1, 2, 3, 4, 5] 2 = 55 then () else error "sum_lmb [1, 2, 3, 4, 5] 2 changed";

if sum_lmb [1, 2, 3, 4, 5] 3 = 225 then () else error "sum_lmb [1, 2, 3, 4, 5] 3 changed";

if sum_lmb [1, 2, 3, 4, 5] 4 = 979 then () else error "sum_lmb [1, 2, 3, 4, 5] 4 changed";

if sum_lmb [1, 2, 3, 4, 5] 5 = 4425 then () else error "sum_lmb [1, 2, 3, 4, 5] 5 changed";

if sum_lmb [1, 2, 3, 4, 5] 6 = 20515 then () else error "sum_lmb [1, 2, 3, 4, 5] 6 changed";

if sum_lmb [~1, 2, ~3, 4, ~5] 1 = ~3 then () else error "sum_lmb [~1, 2, ~3, 4, ~5] 1 changed";

if sum_lmb [~1, 2, ~3, 4, ~5] 2 = 55 then () else error "sum_lmb [~1, 2, ~3, 4, ~5] 2 changed";

if sum_lmb [~1, 2, ~3, 4, ~5] 3 = ~81 then () else error "sum_lmb [~1, 2, ~3, 4, ~5] 3 changed";

if sum_lmb [~1, 2, ~3, 4, ~5] 4 = 979 then () else error "sum_lmb [~1, 2, ~3, 4, ~5] 4 changed";

if sum_lmb [~1, 2, ~3, 4, ~5] 5 = ~2313 then () else error "sum_lmb [~1, 2, ~3, 4, ~5] 5 changed";

if sum_lmb [~1, 2, ~3, 4, ~5] 6 = 20515 then () else error "sum_lmb [~1, 2, ~3, 4, ~5] 6 changed";

"----------- fun LANDAU_MIGNOTTE_bound ------------------";

"----------- fun LANDAU_MIGNOTTE_bound ------------------";

"----------- fun LANDAU_MIGNOTTE_bound ------------------";

if LANDAU_MIGNOTTE_bound [1] [4, 5] = 1 then () else error "LANDAU_MIGNOTTE_bound 1";

if LANDAU_MIGNOTTE_bound [1, 2] [4, 5] = 2 then () else error "LANDAU_MIGNOTTE_bound 2";

if LANDAU_MIGNOTTE_bound [1, 2, 3] [4, 5] = 2 then () else error "LANDAU_MIGNOTTE_bound 3";

if LANDAU_MIGNOTTE_bound [1, 2, 3] [4] = 1 then () else error "LANDAU_MIGNOTTE_bound 4";

if LANDAU_MIGNOTTE_bound [1, 2, 3] [4, 5] = 2 then () else error "LANDAU_MIGNOTTE_bound 5";

if LANDAU_MIGNOTTE_bound [1, 2, 3] [4, 5, 6] = 12 then () else error "LANDAU_MIGNOTTE_bound 6";

if LANDAU_MIGNOTTE_bound [1, 2, 3] [4, 5] = LANDAU_MIGNOTTE_bound [4, 5] [1, 2, 3] then () 

  else error "LANDAU_MIGNOTTE_bound 7";

if LANDAU_MIGNOTTE_bound [~1] [4, 5] = 1 then () else error "LANDAU_MIGNOTTE_bound 11";

if LANDAU_MIGNOTTE_bound [~1, 2] [4, 5] = 2 then () else error "LANDAU_MIGNOTTE_bound 12";

if LANDAU_MIGNOTTE_bound [~1, 2, ~3] [4, ~5] = 2 then () else error "LANDAU_MIGNOTTE_bound 13";

if LANDAU_MIGNOTTE_bound [~1, 2, ~3] [4] = 1 then () else error "LANDAU_MIGNOTTE_bound 14";

if LANDAU_MIGNOTTE_bound [~1, 2, ~3] [4, ~5] = 2 then () else error "LANDAU_MIGNOTTE_bound 15";

if LANDAU_MIGNOTTE_bound [~1, 2, ~3] [4, ~5, 6] = 12 then () else error "LANDAU_MIGNOTTE_bound 16";

if LANDAU_MIGNOTTE_bound [~1, 2, ~3] [4, ~5] = LANDAU_MIGNOTTE_bound [4, 5] [1, 2, 3] then () 

  else error "LANDAU_MIGNOTTE_bound 17" 

"----------- fun chinese_remainder_up -------------------";

"----------- fun chinese_remainder_up -------------------";

"----------- fun chinese_remainder_up -------------------";

if (chinese_remainder_up (5, 7) ([2, 2, 4, 3], [3, 2, 3, 5])) = [17, 2, 24, 33] then ()

else error "chinese_remainder_up (5, 7)... changed";

"----------- fun mod_up ---------------------------------";

"----------- fun mod_up ---------------------------------";

"----------- fun mod_up ---------------------------------";

if ([5, 4, 7, 8, 1] mod_up 5) = [0, 4, 2, 3, 1] then () 

else error "[5, 4, 7, 8, 1] mod_up 5 = [0, 4, 2, 3, 1] changed" ;

if ([5, 4, ~7, 8, ~1] mod_up 5) = [0, 4, 3, 3, 4] then () 

else error "[5, 4, ~7, 8, ~1] mod_up 5  = [0, 4, 3, 3, 4] changed" ;

"----------- fun mod_up_gcd -----------------------------";

"----------- fun mod_up_gcd -----------------------------";

"----------- fun mod_up_gcd -----------------------------";

if mod_up_gcd [~18, ~15, ~20, 12, 20, ~13, 2] [8, 28, 22, ~11, ~14, 1, 2] 7 =

  [2, 6, 0, 2, 6] then () 

else error "( mod_up_gcd [~18, ~15, ~20, ...] [8, 28, 22, ...] 7 = [5, 1, 0, ..] changed";

if mod_up_gcd [8, 28, 22, ~11, ~14, 1, 2] [2, 6, 0, 2, 6] 7 = [2, 6, 0, 2, 6] then ()

else error "mod_up_gcd [8, 28, 22, ...] [2, 6, 0, ...] 7 = [1, 3, 0, 1, 3] changed";

if mod_up_gcd [20, 15, 8, 6] [8, ~2, ~2, 3] 2 = [0, 1] then ()

else error " mod_up_gcd [20, 15, 8, 6] [8, ~2, ~2, 3] 2 = [0, 1] changed";

"~~~~~ fun mod_up_gcd , args:"; 

val (a,b,m) = ([ ~13, 2,12],[ ~14, 2],5);

val moda = a mod_up m;

val modb = (replicate (length a - length(drop_lc0_up b)) 0) @ (drop_lc0_up b mod_up m) ;

val c =  mod_div (lcoeff_up moda) (lcoeff_up modb) m;

val rest = drop_lc0_up (moda %-% (modb %* c) mod_up m) mod_up m;

rest = []; (*=false*)

length rest < length b; (*=false*)

"~~~~~ fun  mod_up_gcd , args:";

val (a,b,m) = (rest, b , m );

val moda = a mod_up m

val modb = (replicate (length a - length(drop_lc0_up b)) 0) @ (drop_lc0_up b mod_up m) ;

val c =  mod_div (lcoeff_up moda) (lcoeff_up modb) m

val rest = drop_lc0_up ((moda %-% (modb %*  c)) mod_up m);

rest = [];(*=flase*)

length rest < length b; (*=true*)

"~~~~~ fun  mod_up_gcd , args:";

val (a,b,m) = (b, rest, m);

val moda = a mod_up m

val modb = (replicate (length a - length(drop_lc0_up b)) 0) @ (drop_lc0_up b mod_up m) ;

val c =  mod_div (lcoeff_up moda) (lcoeff_up modb) m

val rest = drop_lc0_up ((moda %-% (modb %*  c)) mod_up m);

rest = [];(*=flase*)

length rest < length b; (*=false*)

"~~~~~ fun  mod_up_gcd , args:";

val (a,b,m) = (rest,b, m);

val moda = a mod_up m

val modb = (replicate (length a - length(drop_lc0_up b)) 0) @ (drop_lc0_up b mod_up m) ;

val c =  mod_div (lcoeff_up moda) (lcoeff_up modb) m

val rest = drop_lc0_up ((moda %-% (modb %*  c)) mod_up m);

rest = [];(*=true*)

"~~~~~ to  return val:"; val b = b;

"----------- fun primitive_up ---------------------------";

"----------- fun primitive_up ---------------------------";

"----------- fun primitive_up ---------------------------";

if primitive_up [12, 16, 32, 44] = [3, 4, 8, 11] then () 

  else error "primitive_up [12, 16, 32, 44] = [3, 4, 8, 11] changed";

if primitive_up [4, 5, 12] =  [4, 5, 12] then () else error "primitive_up [4, 5, 12] =[4, 5, 12] changed";

if primitive_up [0] = [0] then () else error "primitive_up [0] = [0] changed";

if primitive_up [6] = [1] then () else error "primitive_up [6] = [1] changed";

"----------- fun try_new_prime_up -----------------------";

"----------- fun try_new_prime_up -----------------------";

"----------- fun try_new_prime_up -----------------------";

"~~~~~ fun gcd_up, args:"; val (a, b) = ([~18, ~15, ~20, 12, 20, ~13, 2], 

                                              [8, 28, 22, ~11, ~14, 1, 2]);

"~~~~~ fun HENSEL_lifting_up, args:"; val (a, b, d, M, p) = 

  (a, b, (abs (Integer.gcd (lcoeff_up a) (lcoeff_up b))),  

    (2 * (abs (Integer.gcd (lcoeff_up a) (lcoeff_up b))) * LANDAU_MIGNOTTE_bound a b), 1);

      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;

deg_up g_p = 0 (* = false*);

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

"~~~~~ fun try_new_prime_up, args:"; val (a, b, d, M, P, g, p) = (a, b, d, M, p, g_p, p);

P > M (* = false*);

        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 (a, b, d, M, P, g, p) = (0, 0, 0, 0, 0, 0, 0) (*isolate test below*)

val (a, b) = ([~18, ~15, ~20, 12, 20, ~13, 2], [8, 28, 22, ~11, ~14, 1, 2])

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

val M = 2 * d * LANDAU_MIGNOTTE_bound a b

val p = 1 (*-> p in head of HENSEL_lifting_up*) 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 (*see above*)

val (a, b, d, M, P, g, p) = (a, b, d, M, p, g_p, p);

if try_new_prime_up a b d M P g p = [~1, 1, ~1] then () else error "try_new_prime_up changed";

(* ---- output from changeset cf55b1438731

try_new_prime_up: a = [~18, ~15, ~20, 12, 20, ~13, 2], b = [8, 28, 22, ~11, ~14, 1, 2], 

d = 2, M = 10240, P = 12597, g = [~1, 1, ~1], p = 19 

try_new_prime_up 1 -----> [~1, 1, ~1] *)

val (a, b, d, M, P, g, p) = 

  ([~18, ~15, ~20, 12, 20, ~13, 2], [8, 28, 22, ~11, ~14, 1, 2],

    2, 10240, 3, [~1, 1, ~1], 3);

"-----";

if try_new_prime_up a b d M P g p = [~1, 1, ~1] 

  then () else error "try_new_prime_up I [~18, ~15, ... changed";

(* ---- output from changeset cf55b1438731

try_new_prime_up: a = [~18, ~15, ~20, 12, 20, ~13, 2], b = [8, 28, 22, ~11, ~14, 1, 2], 

d = 2, M = 10240, P = 96577, g = [~4, ~2, 2], p = 23 

try_new_prime_up 1 -----> [~4, ~2, 2] *)

val (a, b, d, M, P, g, p) = 

  ([~18, ~15, ~20, 12, 20, ~13, 2], [8, 28, 22, ~11, ~14, 1, 2],

    2, 10240, 5, [2, 2, 2, 2], 5);

if try_new_prime_up a b d M P g p = [~4, ~2, 2]

  then () else error "try_new_prime_up II [~18, ~15, ... changed";

(* output for "gcd_up [~1, 0 ,1] [0, 1, 1] = [1, 1]" from changeset cf55b1438731

try_new_prime_up:  a = [~1, 0, 1], b = [0, 1, 1], d = 1, M = 16, P = 2, g = [1, 1], p = 2 

try_new_prime_up:  a = [~1, 0, 1], b = [0, 1, 1], d = 1, M = 16, P = 6, g = [1, 1], p = 3 

try_new_prime_up:  a = [~1, 0, 1], b = [0, 1, 1], d = 1, M = 16, P = 30, g = [1, 1], p = 5 

try_new_prime_up 1 -----> [1, 1]*)

val (a, b, d, M, P, g, p) = ([~1, 0, 1], [0, 1, 1], 1, 16, 2, [1, 1], 2);

if try_new_prime_up a b d M P g p = [1, 1]

  then () else error "try_new_prime_up [~1, 0, ... changed";

"----------- fun gcd_up ---------------------------------";

"----------- fun gcd_up ---------------------------------";

"----------- fun gcd_up ---------------------------------";

if gcd_up [~18, ~15, ~20, 12, 20, ~13, 2] [8, 28, 22, ~11, ~14, 1, 2] = [~2, ~1, 1]

  then () else error "gcd_up [~1, 0 ... changed";

if gcd_up [~1, 0 ,1] [0, 1, 1] = [1, 1] 

  then () else error "gcd_up [~1, 0 ... changed";

"=========== Multivariate Case ==========================";

"=========== Multivariate Case ==========================";

"=========== Multivariate Case ==========================";

"----------- fun lcoeff ---------------------------------";

"----------- fun lcoeff ---------------------------------";

"----------- fun lcoeff ---------------------------------";

if lcoeff [(1,[1,1]),(2,[1,2]),(3,[1,3]),(4,[1,4]),(5,[1,5])] = 5 then () else 

error "lcoeff [(1,[1,1]),(2,[1,2]),(3,[1,3]),(4,[1,4]),(5,[1,5])] = 5 changed";

"----------- fun lexp -----------------------------------";

"----------- fun lexp -----------------------------------";

"----------- fun lexp -----------------------------------";

if lexp [(1,[1,1]),(2,[1,2]),(3,[1,3]),(4,[1,4]),(5,[1,5])] = [1,5] then () else 

error "get_coeff [(1,[1,1]),(2,[1,2]),(3,[1,3]),(4,[1,4]),(5,[1,5])] = [1,5] changed";

"----------- fun add_variable ---------------------------";

"----------- fun add_variable ---------------------------";

"----------- fun add_variable ---------------------------";

if add_variable 0

 [(1,[1,2]),(2,[2,3])] = [(1, [0, 1, 2]), (2, [0, 2, 3])] then () else 

error "add_variable [(1,[1,2]),(2,[2,3])] 0 = [(1, [0, 1, 2]), (2, [0, 2, 3])] changed";

if add_variable 1

 [(1,[1,2]),(2,[2,3])] = [(1, [1, 0,  2]), (2, [2, 0, 3])] then () else 

error "add_variable [(1,[1,2]),(2,[2,3])] 1 = [(1, [1, 0,  2]), (2, [2, 0, 3]) changed";

if add_variable 2

 [(1,[1,2]),(2,[2,3])] = [(1, [1, 2, 0]), (2, [2, 3, 0])] then () else 

error "add_variable [(1,[1,2]),(2,[2,3])] 2 = [(1, [1, 2, 0]), (2, [2, 3, 0])] changed";

"----------- fun zero_poly ------------------------------";

"----------- fun zero_poly ------------------------------";

"----------- fun zero_poly ------------------------------";

if zero_poly 1 = [(0, [0])] then () else 

error "zero_poly 1 = [(0, [0])] changed";

if zero_poly 5 = [(0, [0, 0, 0, 0, 0])] then () else 

error "zero_poly 1 = [(0, [0, 0, 0, 0, 0])] changed";

"----------- drop_lc0 -----------------------------------";

"----------- drop_lc0 -----------------------------------";

"----------- drop_lc0 -----------------------------------";

if drop_lc0 [(4,[2]),(3,[3])] = [(4,[2]),(3,[3])] then () 

else error "drop_lc0 [(4,[2]),(3,[3])] = [(4,[2]),(3,[3])] changed";

if drop_lc0 [(4,[2]),(0,[3])] = [(4,[2])] then () 

else error "drop_lc0 [(4,[2]),(3,[3])] = [(4,[2])] changed";

if drop_lc0 [(0,[2]),(3,[3])] = [(0,[2]),(3,[3])] then () 

else error "drop_lc0 [(0,[2]),(3,[3])] = [(0,[2]),(3,[3])] changed";

if drop_lc0 [(0,[2]),(0,[3])] = [(0,[0])] then () 

else error "drop_lc0  [(0,[2]),(0,[3])] = [(0,[0])] changed";

"----------- fun add_monoms -----------------------------";

"----------- fun add_monoms -----------------------------";

"----------- fun add_monoms -----------------------------";

if add_monoms [(~3,[0,0]),(4,[0,0]),(3,[1,1]),(~3,[1,1]),(2,[1,2]),(~3,[1,2])] = 

  [(1, [0, 0]), (~1, [1, 2])] 

then () else error ("add_monoms [(~3,[0,0]),(4,[0,0]),(3,[1,1]),(~3,[1,1]),(2,[1,2]),(~3,[1,2])] " ^

  "= [(1, [0, 0]), (~1, [1, 2])] changed")

"----------- fun lex_ord --------------------------------";

"----------- fun lex_ord --------------------------------";

"----------- fun lex_ord --------------------------------";

if lex_ord (000, [3]) (000, [4])

then error " lex_ord (000, [3]) (000, [4]) changed = false"  else ();

if lex_ord (000, [1,2]) (000, [0,3])

then error " lex_ord [1,2] [0,3] changed = false"  else ();

if lex_ord  (000, [1,2]) (000, [3,0])

then ()  else error " lex_ord  [1,2] [3,0] = true changed";

if lex_ord (000, [1,2]) (000, [1,2])

then error " lex_ord [1,2] [1,2] changed = false"  else ();

"----------- fun order ----------------------------------";

"----------- fun order ----------------------------------";

"----------- fun order ----------------------------------";

if order [(3,[1,1]),(2,[1,2]),(~3,[0,0]),(~3,[1,1]),(~3,[1,2]),(4,[0,0])] = 

  [(1, [0, 0]), (~1, [1, 2])] 

then () else error ("order [(3,[1,1]),(2,[1,2]),(~3,[0,0]),(~3,[1,1]),(~3,[1,2]),(4,[0,0])]" ^

  " = [(1, [0, 0]), (~1, [1, 2])] changed")

(*"----------- fun lcoeff ---------------------------";

"----------- fun lcoeff ---------------------------";

"----------- fun lcoeff ---------------------------";

if lcoeff [(~3,[0,0]),(4,[0,0]),(3,[1,1]),(~3,[1,1]),(2,[1,2]),(~3,[1,2])] = ~1 

then () else error "lcoeff [(~3,[0,0]),(4,[0,0]),(3,[1,1]),(~3,[1,1]),(2,[1,2]),(~3,[1,2])] = ~1  changed";

if lcoeff [(3,[1,1]),(2,[1,2]),(~3,[0,0])] = 2 

then () else error "lcoeff [(3,[1,1]),(2,[1,2]),(~3,[0,0])] = 2   changed";*)

"----------- fun lexp -----------------------------------";

"----------- fun lexp -----------------------------------";

"----------- fun lexp -----------------------------------";

if lexp [(~3,[0,0]),(4,[0,0]),(3,[1,1]),(~3,[1,1]),(2,[1,2])] = [1,2]

then () else error "lexp [(~3,[0,0]),(4,[0,0]),(3,[1,1]),(~3,[1,1]),(2,[1,2])] = ~1  changed";

"----------- fun max_deg_var ----------------------------";

"----------- fun max_deg_var ----------------------------";

"----------- fun max_deg_var ----------------------------";

if max_deg_var [(1,[1,2]),(2,[2,3]),(3,[5,4]),(1,[0,0])] 0 = 5 then ()

else error " max_deg_var [(1,[1,2]),(2,[2,3]),(3,[5,4]),(1,[0,0])] 0 = 5 changed";

if max_deg_var [(1,[1,2]),(2,[2,3]),(3,[5,4]),(1,[0,0])] 1 = 4 then ()

else error " max_deg_var [(1,[1,2]),(2,[2,3]),(3,[5,4]),(1,[0,0])] 1 = 4 changed";

"----------- infix %%/ ----------------------------------";

"----------- infix %%/ ----------------------------------";

"----------- infix %%/ ----------------------------------";

if ([(~3, [2, 0]), (~6, [1, 1]), (~1, [3, 1]),(1, [5, 0]), (2, [4, 1])] %%/ 2) = 

    [(~1, [2, 0]), (~3, [1, 1]), (1, [4, 1])] then ()

else error ("[(~3, [2, 0]), (~6, [1, 1]), (~1, [3, 1]),(1, [5, 0]), (2, [4, 1])] %%/ 2  =  " ^

  "[(~1, [2, 0]), (~3, [1, 1]), (1, [4, 1])] changed");

"----------- infix %%*%% --------------------------------";

"----------- infix %%*%% --------------------------------";

"----------- infix %%*%% --------------------------------";

val p1 = [(1,[1,0]), (1,[0,1])];     (* a   + b   *)

val p2 = [(1,[1,0]), (~1,[0,1])];    (* a   - b   *)

val p = [(1, [2, 0]), (~1, [0, 2])]; (* a^2 - b^2 *)

if (p1 %%*%% p2) = p then () else error "%%*%% changed";

val p1 = order [(1,[1,2,3]), (2,[2,2,2]), (3,[3,2,1])]: poly;

val p2 = order [(10,[1,2,3]), (20,[1,1,1]), (30,[3,2,1])]: poly;

if (p1 %%*%% p2) = 

  [(60, [4, 3, 2]), (90, [6, 4, 2]), (40, [3, 3, 3]), (60, [5, 4, 3]), 

   (20, [2, 3, 4]), (60, [4, 4, 4]), (20, [3, 4, 5]), (10, [2, 4, 6])] then ()

else error "%%*%% changed"

"----------- fun gcd_coeff ------------------------------";

"----------- fun gcd_coeff ------------------------------";

"----------- fun gcd_coeff ------------------------------";

if gcd_coeff [4,8,12,~2,0,~4] = 2 then () else error " gcd_coeff [4,8,12,~2,0,~4] = 2 changed";

if gcd_coeff [3,6,9,19] = 1 then () else error " gcd_coeff [3,6,9,19] = 1 changed";

(*"----------- fun gcd_coeff_poly -------------------------";

"----------- fun gcd_coeff_poly -------------------------";

"----------- fun gcd_coeff_poly -------------------------";

if gcd_coeff_poly [(3,[1,2,3,1])] = 3 then () else error " gcd_coeff_poly [(3,[1,2,3,1])] = 3 changed";

if gcd_coeff_poly [(~2, [2, 0]), (4, [5, 0]),(~4, [3, 2]), (2, [2, 3])] = 2 then () 

else error "gcd_coeff_poly [(~2, [2, 0]), (4, [5, 0]),(~4, [3, 2]), (2, [2, 3])] = 2 changed";

if gcd_coeff_poly [(~2, [2, 0]), (5, [5, 0]),(~4, [3, 2]), (2, [2, 3])] = 1 then ()

else error "gcd_coeff_poly [(~2, [2, 0]), (5, [5, 0]),(~4, [3, 2]), (2, [2, 3])] = 1 changed";*)

"----------- fun eval -----------------------------------";

"----------- fun eval -----------------------------------";

"----------- fun eval -----------------------------------";

if eval [(~3, [2, 1]), (2, [0, 1]), (1,[2, 0])] 0 1 = [(1, [0]), (~1, [1])] then () 

  else error " eval [(~3, [2, 1]), (2, [0, 1]),(1,[2,0])] 0 1 = [(1, [0]), (~1, [1])] changed";

if eval [(~3, [2, 1]), (2, [0, 1]), (1,[2, 0])] 1 1 = [(2, [0]), (~2, [2])] then () 

  else error " eval [(~3, [2, 1]), (2, [0, 1]),(1,[2,0])] 1 1 =[(2, [0]), (~2, [2])] changed";

if eval [(~3, [2, 1]), (2, [0, 1]), (1,[2, 0])] 0 3 = [(9, [0]), (~25, [1])] then () 

  else error " eval [(~3, [2, 1]), (2, [0, 1]),(1,[2,0])] 0 3 = [(9, [0]), (~25, [1])] changed";

if eval [(~3, [2, 1]), (2, [0, 1]), (1,[2, 0])] 1 3 = [(6, [0]), (~8, [2])] then () 

  else error " eval [(~3, [2, 1]), (2, [0, 1]),(1,[2,0])] 1 3 = [ (6, [0]), (~8, [2])] changed";

"----------- fun multi_to_uni ---------------------------";

"----------- fun multi_to_uni ---------------------------";

"----------- fun multi_to_uni ---------------------------";

if multi_to_uni [(~3, [1]), (2, [1]), (1, [0])] = [1, ~1]

then () else error " multi_to_uni [(~3, [1]), (2, [1]), (1, [0])] = [1, ~1] changed";

"----------- fun uni_to_multi ---------------------------";

"----------- fun uni_to_multi ---------------------------";

"----------- fun uni_to_multi ---------------------------";

if uni_to_multi [~11, 11, 22, 33, 44, 55] = 

  [(~11, [0]), (11, [1]), (22, [2]), (33, [3]), (44, [4]), (55, [5])] then ()

else error "uni_to_multi changed";

if uni_to_multi [1,2,1,1,3,4] = [(1, [0]), (2, [1]), (1, [2]), (1, [3]), (3, [4]), (4, [5])] then ()

else error "uni_to_multi [1,2,1,1,3,4] = [(1, [0]), (2, [1]), (1, [2]), (1, [3]), (3, [4]), (4, [5])] changed";

if uni_to_multi [1,2,0,0,5,6] = [(1, [0]), (2, [1]), (5, [4]), (6, [5])] then ()

else error "uni_to_multi [1,2,0,0,5,6] = [(1, [0]), (2, [1]), (5, [4]), (6, [5])] changed";

"----------- fun find_new_r  ----------------------------";

"----------- fun find_new_r  ----------------------------";

"----------- fun find_new_r  ----------------------------";

if find_new_r [(2,[0,0]),(1,[1,2])] [(4,[2,2]),(1,[1,2])] ~1 = 1 then ()

else error " find_new_r [(2,[0,0]),(1,[1,2])] [(4,[2,2]),(1,[1,2])] ~1 = 1 changed";

if find_new_r [(2,[0,0]),(1,[1,2])] [(4,[2,2]),(1,[1,2])] 1 = 2 then ()

else error "find_new_r [(2,[0,0]),(1,[1,2])] [(4,[2,2]),(1,[1,2])] 1 = 2 changed";

if find_new_r [(2,[1,2]),(~4,[0,2]),(2,[1,1])] [(1,[1,3]),(~2,[0,3])] 1 = 3 then ()

else error "find_new_r [(2,[1,2]),(~4,[0,2]),(2,[1,1])] [(1,[1,3]),(~2,[0,3])] 1 = 3 changed";

"----------- fun mult_with_new_var  ---------------------";

"----------- fun mult_with_new_var  ---------------------";

"----------- fun mult_with_new_var  ---------------------";

if mult_with_new_var [(3,[0]),(2,[1])] [~1,1] 0 = [(~3, [0, 0]), (3, [1, 0]), (~2, [0, 1]), (2, [1, 1])] then ()

else error "mult_with_new_var [(3,[0]),(2,[1])] [~1,1] 0 = [(~3, [0, 0]), (3, [1, 0]), (~2, [0, 1]), (2, [1, 1])] changed";

(*"----------- fun greater_deg  -----------------";

"----------- fun greater_deg  -----------------";475

"----------- fun greater_deg  -----------------";

greater_deg [1,2,3] [1,2,4] =2;*)

"----------- infix %%+%%  -------------------------------";

"----------- infix %%+%%  -------------------------------";

"----------- infix %%+%%  -------------------------------";

if ([(3,[0,0]),(2,[3,2]),(~3,[1,3])] %%+%% [(3,[3,2]),(~1,[0,0]),(2,[1,1])]) = [(2, [0, 0]), (2, [1, 1]), (5, [3, 2]), (~3, [1, 3])]

then () else error "([(3,[0,0]),(2,[3,2]),(~3,[1,3])] %%+%% [(3,[3,2]),(~1,[0,0]),(2,[1,1])]) = [(2, [0, 0]), (2, [1, 1]), (5, [3, 2]), (~3, [1, 3])] changed";

"----------- infix %%-%%  -------------------------------";

"----------- infix %%-%%  -------------------------------";

"----------- infix %%-%%  -------------------------------";

if ([(2, [0, 0]), (2, [1, 1]), (5, [3, 2]), (~3, [1, 3])] %%-%% [(3,[3,2]),(~1,[0,0]),(2,[1,1])]) = [(3, [0, 0]), (2, [3, 2]), (~3, [1, 3])]

then () else error "([(2, [0, 0]), (2, [1, 1]), (5, [3, 2]), (~3, [1, 3])] %%-%% [(3,[3,2]),(~1,[0,0]),(2,[1,1])]) = [(3, [0, 0]), (2, [3, 2]), (~3, [1, 3])] changed";

"----------- infix %%*%%  -------------------------------";

"----------- infix %%*%%  -------------------------------";

"----------- infix %%*%%  -------------------------------";

if ([(~1,[0]),(1,[1])] %%*%% [(~2,[0]),(1,[1])]) = [(2, [0]), (~3, [1]), (1, [2])]

then () else error "([(~1,[0]),(1,[1])] %%*%% [(~2,[0]),(1,[1])]) = [(2, [0]), (~3, [1]), (1, [2])] changed";

if ([(~2,[0]),(1,[1])] %%*%% [(1,[1]),(2,[0])]) = [(~4, [0]), (1, [2])] then ()

else error "([(~2,[0]),(1,[1])] %%*%% [(1,[1]),(2,[0])]) = [(~4, [0]), (1, [2])] changed";

if ([(3,[0,0]),(2,[0,1])] %%*%% [(1,[1,0]),(~1,[0,0])]) = mult_with_new_var [(3,[0]),(2,[1])] [~1,1] 0

then () else error 

"([(3,[0,0]),(2,[0,1])] %%*%% [(1,[1,0]),(~1,[0,0])]) = mult_with_new_var [(3,[0]),(2,[1])] [~1,1] 0 changed";

"----------- infix %%/%%  -------------------------------";

"----------- infix %%/%%  -------------------------------";

"----------- infix %%/%%  -------------------------------";

val a = [(1,[2, 0]), (1,[1, 1]), (~1,[0, 2])]; (* (x^2 + x - y) / (x - y) = *)

val b = [(~1, [1, 0]), (1, [0, 1])];

val (c, r) = (a %%/%% b);

if (c, r) = ([(~1, [0, 1])], [(1, [2, 0])]) then () else error "%%/%% changed";

if ((c %%*%% b) %%+%% r) = a then () else error "%%/%% changed";

"----------- fun newton_first  --------------------------";

"----------- fun newton_first  --------------------------";

"----------- fun newton_first  --------------------------";

if newton_first [1, 2] [[(1,[1,1])], [(4,[1,1])]] 1 = 

([(~2, [1, 0, 1]), (3, [1, 1, 1])], [~1, 1], [[(3, [1, 1])]]) 

then () else error 

"newton_first [1, 2] [[(1,[1,1])], [(4,[1,1])]] 1 = "

"([(~2, [1, 0, 1]), (3, [1, 1, 1])], [~1, 1], [[(3, [1, 1])]]) changed";

"----------- fun NEWTON  --------------------------------";

"----------- fun NEWTON  --------------------------------";

"----------- fun NEWTON  --------------------------------";

if NEWTON [1, 2 ] [[(1,[1,1])], [(4,[1,1])]] [] [1] (zero_poly 2) 1 = 

([(~2, [1, 0, 1]), (3, [1, 1, 1])], [~1, 1], [[(3, [1, 1])]])

then () else error

"NEWTON [1, 2 ] [[(1,[1,1])], [(4,[1,1])]] [] [1] (zero_poly 2) 1 = "

"([(~2, [1, 0, 1]), (3, [1, 1, 1])], [~1, 1], [[(3, [1, 1])]]) changed";

if NEWTON [1, 2, 3 ] [[(4,[1,1])], [(9,[1,1])]] [[(3, [1, 1])]] [~1, 1] [(~2, [1, 0, 1]), (3, [1, 1, 1])] 1 = 

([(1, [1, 2, 1])], [2, ~3, 1], [[(5, [1, 1])], [(1, [1, 1])]])

then () else error

"NEWTON [1, 2, 3 ] [[(4,[1,1])], [(9,[1,1])]] [[(3, [1, 1])]] [~1, 1] [(~2, [1, 0, 1]), (3, [1, 1, 1])] 1 ="

" ([(1, [1, 2, 1])], [2, ~3, 1], [[(5, [1, 1])], [(1, [1, 1])]]) changed";

"~~~~~ fun NEWTON, args:"; val (x,f,steps,t,p,ord) = ([1, 2, 3, 4],[[(9, [0]), (5, [1])], [(16, [0]), (7, [1])]], [[(5, [0]), (2, [1])], [(1, [0])]], [2, ~3, 1], [(1, [2, 0]), (~1, [0, 1]), (2, [1, 1])], 0  );

length x = 2; (* false *)

val new_value_poly = multi_to_uni((uni_to_multi t)  %%*%% (uni_to_multi [(nth x (length x -2) )* ~1, 1]));

val new_steps = [((nth f (length f -1)) %%/ ((nth x (length x - 1)) - (nth x (length x - 2)))) %%-%% ((nth f (length f -2)))];

"~~~~~ fun next_step, args:"; val (steps,new_steps, x') = (steps, new_steps, x);

steps = []; (*false*)

"~~~~~ fun next_step, args:"; val (steps,new_steps, x') = (nth_drop 0 steps, new_steps @ [(((nth new_steps (length new_steps -1)) %%-%%(nth steps 0))) %%/

                                  ((nth x' (length x' - 1)) - (nth x' (length x' - 3)))], nth_drop (length x' -2) x');steps = []; (*false*)

"~~~~~ fun next_step, args:"; val (steps,new_steps, x') = (nth_drop 0 steps, new_steps @ [(((nth new_steps (length new_steps -1)) %%-%%(nth steps 0))) %%/

                                  ((nth x' (length x' - 1)) - (nth x' (length x' - 3)))], nth_drop (length x' -2) x');

steps = []; (*true*)

val steps = new_steps;

val polynom' =  p %%+%% (mult_with_new_var (nth steps (length steps -1)) new_value_poly ord);

"----------- fun all_geq  -------------------------------";

"----------- fun all_geq  -------------------------------";

"----------- fun all_geq  -------------------------------";

if all_geq  [1,2,3,4,5] [1,2,3,4,5] = true then ()

else error " all_geq  [1,2,3,4,5] [1,2,3,4,5] = true changed";

if all_geq [1,2,3,4] [1,2,3,5] = false then () 

else error "all_geq [1,2,3,4] [1,2,3,5] = false changed"  ;

if all_geq [1,4,5,4] [0,3,4,5] = false then ()

else error "all_geq [1,2,3,4] [0,3,4,5] = false changed ";

"----------- fun greater_deg  ---------------------------";

"----------- fun greater_deg  ---------------------------";

"----------- fun greater_deg  ---------------------------";

if greater_deg  [1,2,3,4,5] [1,2,3,4,5] = 0 then ()

else error " greater_deg  [1,2,3,4,5] [1,2,3,4,5] = 0 changed";

if greater_deg [1,2,3,4] [5,2,8,1] = 1 then () 

else error "greater_deg [1,2,3,4] [1,2,3,5] = 1 changed"  ;

if greater_deg [1,4,5,4] [0,3,4,5] = 2 then ()

else error "greater_deg [1,2,3,4] [0,3,4,5] = 2 changed ";

"----------- infix %%|%% --------------------------------";

"----------- infix %%|%% --------------------------------";

"----------- infix %%|%% --------------------------------";

if [(1, [0, 0]), (~1, [0, 1])] %%|%% [(1, [0, 0]), (~1, [0, 1])] then ()

(*   1            -       y      |     1            -       y*)

else error "[(1, [0, 0]), (~1, [0, 1])] %%|%% [(1, [0, 0]), (~1, [0, 1])] = true changed";

if [(3,[1,0])] %%|%% [(9,[1,1]),(12,[2,1]),(~3,[1,2])] then ()

(*   3  x        |     9  x y  + 12 x^2 y  - 3 x y^2         *)

else error "[(3,[1,0])] %%|%% [(9,[1,1]),(12,[2,1]),(~3,[1,2])] = true changed";

if [(3,[2,1])] %%|%% [(9,[1,1]),(12,[2,1]),(~3,[1,2])] 

(*   3 x^2 y     |     9  x y  + 12 x^2 y  - 3 x y^2         *)

then error "[(3,[2,1])] %%|%% [(9,[1,1]),(12,[2,1]),(~3,[1,2])] = false changed" else ();

"----------- fun gcd_poly'  ------------------------------";

"----------- fun gcd_poly'  ------------------------------";

"----------- fun gcd_poly'  ------------------------------";

if gcd_poly' [(~3,[2,0]),(1,[5,0]),(3,[0,1]),(~6,[1,1]),(~1,[3,1]),(2,[4,1]),(1,[3,2]),(~1,[1,3]),(2,[2,3])]

[(2,[2,0]),(~2,[0,1]),(4,[1,1]),(~1,[3,1]),(1,[1,2]),(~1,[2,2]),(~1,[0,3]),(2,[1,3])] 2 0 

= [(1, [2, 0]), (~1, [0, 1]), (2, [1, 1])] then () else error 

"gcd_poly' [(~3,[2,0]),(1,[5,0]),(3,[0,1]),(~6,[1,1]),(~1,[3,1]),(2,[4,1]),(1,[3,2]),(~1,[1,3]),(2,[2,3])]"

"[(2,[2,0]),(~2,[0,1]),(4,[1,1]),(~1,[3,1]),(1,[1,2]),(~1,[2,2]),(~1,[0,3]),(2,[1,3])] 2 0 "

"= [(1, [2, 0]), (~1, [0, 1]), (2, [1, 1])] changed";

(* -xy +xy^2z+yz - 1*)(* xy +1*) (*=*) (*xy -1*)

if gcd_poly' [(~1,[0,0,0]),(1,[0,1,1]),(1,[1,2,1]),(~1,[1,1,0])] [(1,[0,0,0]),(1,[1,1,0])] 3 0 

 = [(~1, [0, 0, 0]), (~1, [1, 1, 0])] then () else error

"gcd_poly' [(~1,[0,0,0]),(1,[0,1,1]),(1,[1,2,1]),(~1,[1,1,0])] [(1,[0,0,0]),(1,[1,1,0])] 3 0 "

"= [(1, [0, 0, 0]), (1, [1, 1, 0])] changed";

"----------- fun gcd_poly  ------------------------------";

"----------- fun gcd_poly  ------------------------------";

"----------- fun gcd_poly  ------------------------------";

if gcd_poly 

  [(~3,[2,0]),(1,[5,0]),(3,[0,1]),(~6,[1,1]),(~1,[3,1]),(2,[4,1]),(1,[3,2]),(~1,[1,3]),(2,[2,3])]

 (* ~3 x^2  +     x^5 +  3 y     -  6 x y     - x^3 y +  2 x^4 y +  x^3 y^2  - x y^3 +  2 x^2 y^3*)

  [(2,[2,0]),(~2,[0,1]),(4,[1,1]),(~1,[3,1]),(1,[1,2]),(~1,[2,2]),(~1,[0,3]),(2,[1,3])]

 (* 2 x^2   -  2 y    +  4  x y  -   x^3 y  +    x y^2  - x^2 y^2  -   y^3  + 2 x y^3*) 

= (([(~3, [0, 0]), (1, [3, 0]), (1, [1, 2])], [(2, [0, 0]), (~1, [1, 1]), (1, [0, 2])]), 

 (*(  3          +    x^3     +      x y^2  ,,  2         -      x y    +       y^2   ),,*)

    [(1, [2, 0]), (~1, [0, 1]), (2, [1, 1])])

 (*     x^2     -          y     2   x  y*)

then () else error "gcd_poly ex1 changed";

if gcd_poly 

  [(~1,[0,0,0]),(1,[0,1,1]),(1,[1,2,1]),(~1,[1,1,0])]

 (* ~1         +      y z  +  x y^2 z  -     x y *)

  [(1,[0,0,0]),(1,[1,1,0])]

 (* 1         +    x y *)

  = (([(~1, [0,0,0]), (1, [0,1,1])], [(1, [0,0,0])]), [(1, [0,0,0]), (1, [1,1,0])]) 

 (*  (  ~1          +        y z   ,,  1            ),, 1          +      x y *)

then () else error "gcd_poly ex2 changed";

(* example from unipoly: *)

val a = uni_to_multi [~18, ~15, ~20, 12, 20, ~13, 2];

val b = uni_to_multi [8, 28, 22, ~11, ~14, 1, 2];

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

a' = [(9, [0]), (3, [1]), (13, [2]), (~11, [3]), (2, [4])];

b' = [(~4, [0]), (~12, [1]), (~7, [2]), (3, [3]), (2, [4])];

c = [(~2, [0]), (~1, [1]), (1, [2])];

(* checking the postcondition: *)

if a = (a' %%*%% c) andalso b = (b' %%*%% c) then () else error "gcd_poly ex-unipoly changed";

(* \<forall>c'. (c' dvd\<^isub>p  a \<and> c' dvd\<^isub>p b) \<longrightarrow> c' \<le>\<^isub>p c) can NOT be checked this way, of course *)

val cu = gcd_up [~18, ~15, ~20, 12, 20, ~13, 2] [8, 28, 22, ~11, ~14, 1, 2];

if cu = multi_to_uni c then () else error "gcd_up <> gcd_poly";

(* example from mail to Tobias Nipkos *)

val a = [(1,[2, 0]), (~1,[0, 2])]; (* x^2 - y^2 *)

val b = [(1,[2, 0]), (~1,[1, 1])]; (* x^2 - x y *)

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

a' = [(1, [1, 0]), (1, [0, 1])]; (* x + y  *)

b' = [(1, [1, 0])];              (* x      *)

c = [(1, [1, 0]), (~1, [0, 1])]; (* x - y  *)

(* checking the postcondition: *)

if a = (a' %%*%% c) andalso b = (b' %%*%% c) then () else error "gcd_poly mail changed";

(* \<forall>c'. (c' dvd\<^isub>p  a \<and> c' dvd\<^isub>p b) \<longrightarrow> c' \<le>\<^isub>p c) can NOT be checked this way, of course *)

(* regression test for --- fun gcd_poly downto division by zero ---*)

val a = [(9, [5, 2, 4])]: poly

val b = [(15, [6, 3, 1])]: poly

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

if a = (a' %%*%% c) andalso b = (b' %%*%% c) then () else error "regression test changed";

(* regression test for --- integration lev.1 -- lev.5 in rational.sml *)

val a = [(1, [0, 1, 1, 0, 1]), (1, [1, 0, 0, 1, 1])]: poly

val b = [(1, [0, 1, 0, 1, 1])]: poly

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

a' = [(1, [0, 1, 1, 0, 0]), (1, [1, 0, 0, 1, 0])];

b' = [(1, [0, 1, 0, 1, 0])];

c = [(1, [0, 0, 0, 0, 1])];

if a = (a' %%*%% c) andalso b = (b' %%*%% c) then () else error "gcd_monom changed";

"----------- fun gcd_poly downto division by zero -------";

"----------- fun gcd_poly downto division by zero -------";

"----------- fun gcd_poly downto division by zero -------";

"-------- example 187a: exception Div raised...  REGRSSION TEST FOR REMOVAL OF THIS BUG IS ABOVE";

val a = [(12, [1, 1])]: poly

val b = [(8, [0, 2])]: poly;

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

exception Div raised*)

"~~~~~ fun gcd_poly, args:"; val ((a as (_, es)::_ : poly), b) = (a, b);

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

exception Div raised*)

"~~~~~ fun gcd_poly', args:"; val ((a as (_, es1)::_ : poly), (b as (_, es2)::_ : poly), n, r) =

  (a, b, (length es), 0);

lex_ord (lmonom b) (lmonom a)       = true;

(*gcd_poly' b a n r

exception Div raised*)

"~~~~~ fun gcd_poly', args:"; val ((a as (_, es1)::_ : poly), (b as (_, es2)::_ : poly), n, r) =

  (b, a, n, r);

lex_ord (lmonom b) (lmonom a)       = false;

n > 1                               = true;

 val M = 1 + Int.min (max_deg_var a (length es1 - 2), max_deg_var b (length es2 - 2)); 

(*try_new_r a b n M r [] []

exception Div raised*)

"~~~~~ and try_new_r, args:"; val (a, b, n, M, r, r_list, steps) = (a, b, n, M, r, [], []);

        val r = find_new_r a b r;

        val r_list = r_list @ [r];

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

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

        val steps = steps @ [g_r];

(*HENSEL_lifting a b n M 1 r r_list steps g_r ( zero_poly n ) [1]

exception Div raised*)

"~~~~~ and HENSEL_lifting, args:"; val (a, b, n, M, m, r, r_list, steps, g, g_n, mult) =

  (a, b, n, M, 1, r, r_list, steps, g_r, (zero_poly n), [1]);

m > M                               = false;

        val r = find_new_r a b r; 

        val r_list = r_list @ [r];

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

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

lex_ord (lmonom g) (lmonom g_r)     = false;

lexp g = lexp g_r                   = true;

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

nth steps (length steps - 1) = zero_poly (n - 1)    = false;

(*HENSEL_lifting a b n M (m + 1) r r_list steps g_r g_n new

exception Div raised*)

"~~~~~ and HENSEL_lifting, args:"; val (a, b, n, M, m, r, r_list, steps, g, g_n, mult) =

  (a, b, n, M, m + 1, r, r_list, steps, g_r, g_n, new);

m > M                               = true;

(g_n %%|%% a andalso (g_n %%|%% b)) = false;

(*try_new_r a b n M r r_list steps

exception Div raised*)

"~~~~~ and try_new_r, args:"; val (a, b, n, M, r, r_list, steps) = (a, b, n, M, r, r_list, steps);

        val r = find_new_r a b r;

        val r_list = r_list @ [r];

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

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

        val steps = steps @ [g_r];

(*HENSEL_lifting a b n M 1 r r_list steps g_r (zero_poly n) [1]

exception Div raised*)

"~~~~~ and HENSEL_lifting, args:"; val (a, b, n, M, m, r, r_list, steps, g, g_n, mult) =

  (a, b, n, M, 1, r, r_list, steps, g_r, zero_poly n, [1]);

m > M                               = false;

        val r = find_new_r a b r; 

        val r_list = r_list @ [r];

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

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

lex_ord (lmonom g) (lmonom g_r)     = false;

lexp g = lexp g_r                   = true;

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

(nth steps (length steps - 1) = zero_poly (n - 1))  = false;

(*HENSEL_lifting a b n M (m + 1) r r_list steps g_r g_n new

exception Div raised*)

"~~~~~ and HENSEL_lifting, args:"; val (a, b, n, M, m, r, r_list, steps, g, g_n, mult) =

  (a, b, n, M, m + 1, r, r_list, steps, g_r, g_n, new);

m > M                               = true;

(g_n %%|%% a andalso g_n %%|%% b)   = false;

(*try_new_r a b n M r r_list steps

exception Div raised*)

"~~~~~ and try_new_r, args:"; val (a, b, n, M, r, r_list, steps) = (a, b, n, M, r, r_list, steps);

        val r = find_new_r a b r;

        val r_list = r_list @ [r];

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

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

        val steps = steps @ [g_r];

(*HENSEL_lifting a b n M 1 r r_list steps g_r (zero_poly n) [1]

exception Div raised*)

"~~~~~ and HENSEL_lifting, args:"; val (a, b, n, M, m, r, r_list, steps, g, g_n, mult) =

  (a, b, n, M, 1, r, r_list, steps, g_r, zero_poly n, [1]);

m > M                               = false;

        val r = find_new_r a b r; 

        val r_list = r_list @ [r];

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

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

lex_ord (lmonom g) (lmonom g_r)     = false;

lexp g = lexp g_r                   = true;

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

nth steps (length steps - 1) = zero_poly (n - 1)     = true;

(*HENSEL_lifting a b n M (M + 1) r r_list steps g_r g_n new

exception Div raised*)

"~~~~~ and HENSEL_lifting, args:"; val (a, b, n, M, m, r, r_list, steps, g, g_n, mult) =

  (a, b, n, M, m + 1, r, r_list, steps, g_r, g_n, new);

m > M                               = true;

(*(g_n %%|%% a andalso g_n %%|%% b);

exception Div raised:                !!!!! g_n = [(0, [0, 0])] !!!!!*)

"=========== prep Lucas-Interpretation ==================";

"----------- fun for_quot_regarding ---------------------";

"----------- fun for_quot_regarding ---------------------";

"----------- fun for_quot_regarding ---------------------";

(* -> quot' = [] @ [2] @ (3 * [])*)

"~~~~~ fun for_quot_regarding, args:"; val (p1, p2, p1', quot, remd) =

  ([2, 2, 2, 2, 2], [1, 2, 3], [2, 2, 2, 2, 2], [], [6, 6, 4, 2]);

       val zeros = deg_up p1' - deg_up remd - 1(*length [q]*)  (* = 0*)

       val max_qot_deg = deg_up p1 - deg_up p2                 (* = 2*)

;

zeros + 1(*length [q]*) + deg_up quot > max_qot_deg   = false

;

if for_quot_regarding p1 p2 p1' quot remd = 0 then () 

else error "for_quot_regarding 1 changed"

(* -> quot' = [0, 0, 0] @ [1] @ (1 * [])*)

"~~~~~ fun for_quot_regarding, args:"; val (p1, p2, p1', quot, remd) =

  ([3, 2, 2, 2, 1, 1, 1, 1], [1, 1, 1, 1], [3, 2, 2, 2, 1, 1, 1, 1], [], [3, 2, 2, 2]);

       val zeros = deg_up p1' - deg_up remd - 1(*length [q]*)  (* = 3*)

       val max_qot_deg = deg_up p1 - deg_up p2                 (* = 4*)

;

zeros + 1(*length [q]*) + deg_up quot > max_qot_deg   = false

;

if for_quot_regarding p1 p2 p1' quot remd = 3 then () 

else error "for_quot_regarding 2 changed"

(* -> quot' = [] @ [2] @ (1 * [0, 0, 0, 1])*)

"~~~~~ fun for_quot_regarding, args:"; val (p1, p2, p1', quot, remd) =

  ([3, 2, 2, 2, 1, 1, 1, 1], [1, 1, 1, 1], [3, 2, 2, 2], [0, 0, 0, 1], [1]);

       val zeros = deg_up p1' - deg_up remd - 1(*length [q]*)  (* = 2*)

       val max_qot_deg = deg_up p1 - deg_up p2                 (* = 4*)

;

zeros + 1(*length [q]*) + deg_up quot > max_qot_deg   = true

;

   for_quot_regarding p1 p2 p1' quot remd;

if for_quot_regarding p1 p2 p1' quot remd = 0 then () 

else error "for_quot_regarding 3 changed";

"----------- fun mult_to_deg - --------------------------";

"----------- fun mult_to_deg - --------------------------";

"----------- fun mult_to_deg - --------------------------";

  val p1 = [2,2,2, 2,2,2];

  val p2 = [7,8,6];

if mult_to_deg (deg_up p1) p2 = [0, 0, 0, 7, 8, 6] then () else error "fun mult_to_deg changed";

"----------- fun fact_quot ------------------------------";

"----------- fun fact_quot ------------------------------";

"----------- fun fact_quot ------------------------------";

  fact_quot 18 8 = (4, 9);

  fact_quot 2 6 = (3, 1);

  fact_quot ~2 6 = (3, ~1);

  fact_quot 5 6 = (6, 5);

if fact_quot 110 6 = (3, 55) then () else error "fun fact_quot changed 1";

(* this works, if (several) position(s) in the dividend are 0: *)

if fact_quot 0 3 = (1, 0) then () else error "fun fact_quot changed 2";

"----------- fun %+% ------------------------------------";

"----------- fun %+% ------------------------------------";

"----------- fun %+% ------------------------------------";

if (p1 %+% p2) = [9, 10, 8, 2, 2, 2] andalso (p2 %+% p1) = [9, 10, 8, 2, 2, 2]

  then () else error "fun %+% changed";

"----------- fun %*% ------------------------------------";

"----------- fun %*% ------------------------------------";

"----------- fun %*% ------------------------------------";