src/HOL/Tools/Groebner_Basis/groebner.ML
author haftmann
Tue Oct 20 16:13:01 2009 +0200 (2009-10-20)
changeset 33037 b22e44496dc2
parent 30866 dd5117e2d73e
child 33038 8f9594c31de4
permissions -rw-r--r--
replaced old_style infixes eq_set, subset, union, inter and variants by generic versions
wenzelm@23252
     1
(*  Title:      HOL/Tools/Groebner_Basis/groebner.ML
wenzelm@23252
     2
    Author:     Amine Chaieb, TU Muenchen
wenzelm@23252
     3
*)
wenzelm@23252
     4
wenzelm@23252
     5
signature GROEBNER =
wenzelm@23252
     6
sig
wenzelm@23252
     7
  val ring_and_ideal_conv :
chaieb@30866
     8
    {idom: thm list, ring: cterm list * thm list, field: cterm list * thm list,
chaieb@30866
     9
     vars: cterm list, semiring: cterm list * thm list, ideal : thm list} ->
wenzelm@23487
    10
    (cterm -> Rat.rat) -> (Rat.rat -> cterm) ->
wenzelm@23487
    11
    conv ->  conv ->
chaieb@25251
    12
 {ring_conv : conv, 
chaieb@25251
    13
 simple_ideal: (cterm list -> cterm -> (cterm * cterm -> order) -> cterm list),
chaieb@25251
    14
 multi_ideal: cterm list -> cterm list -> cterm list -> (cterm * cterm) list,
chaieb@25251
    15
 poly_eq_ss: simpset, unwind_conv : conv}
wenzelm@23579
    16
    val ring_tac: thm list -> thm list -> Proof.context -> int -> tactic
chaieb@25251
    17
    val ideal_tac: thm list -> thm list -> Proof.context -> int -> tactic
chaieb@25251
    18
    val algebra_tac: thm list -> thm list -> Proof.context -> int -> tactic
wenzelm@23252
    19
end
wenzelm@23252
    20
chaieb@25251
    21
structure Groebner : GROEBNER =
wenzelm@23252
    22
struct
wenzelm@23252
    23
chaieb@25251
    24
open Conv Normalizer Drule Thm;
chaieb@25251
    25
chaieb@25251
    26
fun is_comb ct =
chaieb@25251
    27
  (case Thm.term_of ct of
chaieb@25251
    28
    _ $ _ => true
chaieb@25251
    29
  | _ => false);
chaieb@25251
    30
chaieb@25251
    31
val concl = Thm.cprop_of #> Thm.dest_arg;
chaieb@25251
    32
chaieb@25251
    33
fun is_binop ct ct' =
chaieb@25251
    34
  (case Thm.term_of ct' of
chaieb@25251
    35
    c $ _ $ _ => term_of ct aconv c
chaieb@25251
    36
  | _ => false);
chaieb@25251
    37
chaieb@25251
    38
fun dest_binary ct ct' =
chaieb@25251
    39
  if is_binop ct ct' then Thm.dest_binop ct'
chaieb@25251
    40
  else raise CTERM ("dest_binary: bad binop", [ct, ct'])
chaieb@25251
    41
chaieb@25251
    42
fun inst_thm inst = Thm.instantiate ([], inst);
wenzelm@23557
    43
wenzelm@23252
    44
val rat_0 = Rat.zero;
wenzelm@23252
    45
val rat_1 = Rat.one;
wenzelm@23252
    46
val minus_rat = Rat.neg;
wenzelm@23252
    47
val denominator_rat = Rat.quotient_of_rat #> snd #> Rat.rat_of_int;
wenzelm@23252
    48
fun int_of_rat a =
wenzelm@23252
    49
    case Rat.quotient_of_rat a of (i,1) => i | _ => error "int_of_rat: not an int";
haftmann@23514
    50
val lcm_rat = fn x => fn y => Rat.rat_of_int (Integer.lcm (int_of_rat x) (int_of_rat y));
wenzelm@23252
    51
wenzelm@23252
    52
val (eqF_intr, eqF_elim) =
wenzelm@23252
    53
  let val [th1,th2] = thms "PFalse"
wenzelm@23252
    54
  in (fn th => th COMP th2, fn th => th COMP th1) end;
wenzelm@23252
    55
wenzelm@23252
    56
val (PFalse, PFalse') =
wenzelm@23252
    57
 let val PFalse_eq = nth (thms "simp_thms") 13
wenzelm@23252
    58
 in (PFalse_eq RS iffD1, PFalse_eq RS iffD2) end;
wenzelm@23252
    59
wenzelm@23252
    60
chaieb@25251
    61
(* Type for recording history, i.e. how a polynomial was obtained. *)
wenzelm@23252
    62
wenzelm@23252
    63
datatype history =
wenzelm@24630
    64
   Start of int
wenzelm@24630
    65
 | Mmul of (Rat.rat * int list) * history
wenzelm@23252
    66
 | Add of history * history;
wenzelm@23252
    67
wenzelm@23252
    68
wenzelm@23252
    69
(* Monomial ordering. *)
wenzelm@23252
    70
wenzelm@23252
    71
fun morder_lt m1 m2=
wenzelm@23252
    72
    let fun lexorder l1 l2 =
wenzelm@23252
    73
            case (l1,l2) of
wenzelm@23252
    74
                ([],[]) => false
wenzelm@23252
    75
              | (x1::o1,x2::o2) => x1 > x2 orelse x1 = x2 andalso lexorder o1 o2
wenzelm@23252
    76
              | _ => error "morder: inconsistent monomial lengths"
wenzelm@24630
    77
        val n1 = Integer.sum m1
wenzelm@24630
    78
        val n2 = Integer.sum m2 in
wenzelm@23252
    79
    n1 < n2 orelse n1 = n2 andalso lexorder m1 m2
wenzelm@23252
    80
    end;
wenzelm@23252
    81
wenzelm@23252
    82
fun morder_le m1 m2 = morder_lt m1 m2 orelse (m1 = m2);
wenzelm@23252
    83
wenzelm@23252
    84
fun morder_gt m1 m2 = morder_lt m2 m1;
wenzelm@23252
    85
wenzelm@23252
    86
(* Arithmetic on canonical polynomials. *)
wenzelm@23252
    87
wenzelm@23252
    88
fun grob_neg l = map (fn (c,m) => (minus_rat c,m)) l;
wenzelm@23252
    89
wenzelm@23252
    90
fun grob_add l1 l2 =
wenzelm@23252
    91
  case (l1,l2) of
wenzelm@23252
    92
    ([],l2) => l2
wenzelm@23252
    93
  | (l1,[]) => l1
wenzelm@23252
    94
  | ((c1,m1)::o1,(c2,m2)::o2) =>
wenzelm@23252
    95
        if m1 = m2 then
wenzelm@23252
    96
          let val c = c1+/c2 val rest = grob_add o1 o2 in
wenzelm@23252
    97
          if c =/ rat_0 then rest else (c,m1)::rest end
wenzelm@23252
    98
        else if morder_lt m2 m1 then (c1,m1)::(grob_add o1 l2)
wenzelm@23252
    99
        else (c2,m2)::(grob_add l1 o2);
wenzelm@23252
   100
wenzelm@23252
   101
fun grob_sub l1 l2 = grob_add l1 (grob_neg l2);
wenzelm@23252
   102
wenzelm@24630
   103
fun grob_mmul (c1,m1) (c2,m2) = (c1*/c2, ListPair.map (op +) (m1, m2));
wenzelm@23252
   104
wenzelm@23252
   105
fun grob_cmul cm pol = map (grob_mmul cm) pol;
wenzelm@23252
   106
wenzelm@23252
   107
fun grob_mul l1 l2 =
wenzelm@23252
   108
  case l1 of
wenzelm@23252
   109
    [] => []
wenzelm@23252
   110
  | (h1::t1) => grob_add (grob_cmul h1 l2) (grob_mul t1 l2);
wenzelm@23252
   111
wenzelm@23252
   112
fun grob_inv l =
wenzelm@23252
   113
  case l of
wenzelm@23252
   114
    [(c,vs)] => if (forall (fn x => x = 0) vs) then
wenzelm@23252
   115
                  if (c =/ rat_0) then error "grob_inv: division by zero"
wenzelm@23252
   116
                  else [(rat_1 // c,vs)]
wenzelm@23252
   117
              else error "grob_inv: non-constant divisor polynomial"
wenzelm@23252
   118
  | _ => error "grob_inv: non-constant divisor polynomial";
wenzelm@23252
   119
wenzelm@23252
   120
fun grob_div l1 l2 =
wenzelm@23252
   121
  case l2 of
wenzelm@23252
   122
    [(c,l)] => if (forall (fn x => x = 0) l) then
wenzelm@23252
   123
                 if c =/ rat_0 then error "grob_div: division by zero"
wenzelm@23252
   124
                 else grob_cmul (rat_1 // c,l) l1
wenzelm@23252
   125
             else error "grob_div: non-constant divisor polynomial"
wenzelm@23252
   126
  | _ => error "grob_div: non-constant divisor polynomial";
wenzelm@23252
   127
wenzelm@23252
   128
fun grob_pow vars l n =
wenzelm@23252
   129
  if n < 0 then error "grob_pow: negative power"
wenzelm@23252
   130
  else if n = 0 then [(rat_1,map (fn v => 0) vars)]
wenzelm@24630
   131
  else grob_mul l (grob_pow vars l (n - 1));
wenzelm@23252
   132
wenzelm@23252
   133
fun degree vn p =
wenzelm@23252
   134
 case p of
wenzelm@23252
   135
  [] => error "Zero polynomial"
wenzelm@23252
   136
| [(c,ns)] => nth ns vn
wenzelm@24630
   137
| (c,ns)::p' => Int.max (nth ns vn, degree vn p');
wenzelm@23252
   138
wenzelm@23252
   139
fun head_deg vn p = let val d = degree vn p in
wenzelm@23252
   140
 (d,fold (fn (c,r) => fn q => grob_add q [(c, map_index (fn (i,n) => if i = vn then 0 else n) r)]) (filter (fn (c,ns) => c <>/ rat_0 andalso nth ns vn = d) p) []) end;
wenzelm@23252
   141
wenzelm@24630
   142
val is_zerop = forall (fn (c,ns) => c =/ rat_0 andalso forall (curry (op =) 0) ns);
wenzelm@23252
   143
val grob_pdiv =
wenzelm@23252
   144
 let fun pdiv_aux vn (n,a) p k s =
wenzelm@23252
   145
  if is_zerop s then (k,s) else
wenzelm@23252
   146
  let val (m,b) = head_deg vn s
wenzelm@23252
   147
  in if m < n then (k,s) else
wenzelm@24630
   148
     let val p' = grob_mul p [(rat_1, map_index (fn (i,v) => if i = vn then m - n else 0)
wenzelm@23252
   149
                                                (snd (hd s)))]
wenzelm@23252
   150
     in if a = b then pdiv_aux vn (n,a) p k (grob_sub s p')
wenzelm@24630
   151
        else pdiv_aux vn (n,a) p (k + 1) (grob_sub (grob_mul a s) (grob_mul b p'))
wenzelm@23252
   152
     end
wenzelm@23252
   153
  end
wenzelm@23252
   154
 in fn vn => fn s => fn p => pdiv_aux vn (head_deg vn p) p 0 s
wenzelm@23252
   155
 end;
wenzelm@23252
   156
wenzelm@23252
   157
(* Monomial division operation. *)
wenzelm@23252
   158
wenzelm@23252
   159
fun mdiv (c1,m1) (c2,m2) =
wenzelm@23252
   160
  (c1//c2,
wenzelm@24630
   161
   map2 (fn n1 => fn n2 => if n1 < n2 then error "mdiv" else n1 - n2) m1 m2);
wenzelm@23252
   162
wenzelm@23252
   163
(* Lowest common multiple of two monomials. *)
wenzelm@23252
   164
wenzelm@24630
   165
fun mlcm (c1,m1) (c2,m2) = (rat_1, ListPair.map Int.max (m1, m2));
wenzelm@23252
   166
wenzelm@23252
   167
(* Reduce monomial cm by polynomial pol, returning replacement for cm.  *)
wenzelm@23252
   168
wenzelm@23252
   169
fun reduce1 cm (pol,hpol) =
wenzelm@23252
   170
  case pol of
wenzelm@23252
   171
    [] => error "reduce1"
wenzelm@23252
   172
  | cm1::cms => ((let val (c,m) = mdiv cm cm1 in
wenzelm@23252
   173
                    (grob_cmul (minus_rat c,m) cms,
wenzelm@23252
   174
                     Mmul((minus_rat c,m),hpol)) end)
wenzelm@23252
   175
                handle  ERROR _ => error "reduce1");
wenzelm@23252
   176
wenzelm@23252
   177
(* Try this for all polynomials in a basis.  *)
wenzelm@23252
   178
fun tryfind f l =
wenzelm@23252
   179
    case l of
wenzelm@23252
   180
        [] => error "tryfind"
wenzelm@23252
   181
      | (h::t) => ((f h) handle ERROR _ => tryfind f t);
wenzelm@23252
   182
wenzelm@23252
   183
fun reduceb cm basis = tryfind (fn p => reduce1 cm p) basis;
wenzelm@23252
   184
wenzelm@23252
   185
(* Reduction of a polynomial (always picking largest monomial possible).     *)
wenzelm@23252
   186
wenzelm@23252
   187
fun reduce basis (pol,hist) =
wenzelm@23252
   188
  case pol of
wenzelm@23252
   189
    [] => (pol,hist)
wenzelm@23252
   190
  | cm::ptl => ((let val (q,hnew) = reduceb cm basis in
wenzelm@23252
   191
                   reduce basis (grob_add q ptl,Add(hnew,hist)) end)
wenzelm@23252
   192
               handle (ERROR _) =>
wenzelm@23252
   193
                   (let val (q,hist') = reduce basis (ptl,hist) in
wenzelm@23252
   194
                       (cm::q,hist') end));
wenzelm@23252
   195
wenzelm@23252
   196
(* Check for orthogonality w.r.t. LCM.                                       *)
wenzelm@23252
   197
wenzelm@23252
   198
fun orthogonal l p1 p2 =
wenzelm@23252
   199
  snd l = snd(grob_mmul (hd p1) (hd p2));
wenzelm@23252
   200
wenzelm@23252
   201
(* Compute S-polynomial of two polynomials.                                  *)
wenzelm@23252
   202
wenzelm@23252
   203
fun spoly cm ph1 ph2 =
wenzelm@23252
   204
  case (ph1,ph2) of
wenzelm@23252
   205
    (([],h),p) => ([],h)
wenzelm@23252
   206
  | (p,([],h)) => ([],h)
wenzelm@23252
   207
  | ((cm1::ptl1,his1),(cm2::ptl2,his2)) =>
wenzelm@23252
   208
        (grob_sub (grob_cmul (mdiv cm cm1) ptl1)
wenzelm@23252
   209
                  (grob_cmul (mdiv cm cm2) ptl2),
wenzelm@23252
   210
         Add(Mmul(mdiv cm cm1,his1),
wenzelm@23252
   211
             Mmul(mdiv (minus_rat(fst cm),snd cm) cm2,his2)));
wenzelm@23252
   212
wenzelm@23252
   213
(* Make a polynomial monic.                                                  *)
wenzelm@23252
   214
wenzelm@23252
   215
fun monic (pol,hist) =
wenzelm@23579
   216
  if null pol then (pol,hist) else
wenzelm@23252
   217
  let val (c',m') = hd pol in
wenzelm@23252
   218
  (map (fn (c,m) => (c//c',m)) pol,
wenzelm@23252
   219
   Mmul((rat_1 // c',map (K 0) m'),hist)) end;
wenzelm@23252
   220
wenzelm@23252
   221
(* The most popular heuristic is to order critical pairs by LCM monomial.    *)
wenzelm@23252
   222
wenzelm@23252
   223
fun forder ((c1,m1),_) ((c2,m2),_) = morder_lt m1 m2;
wenzelm@23252
   224
wenzelm@23252
   225
fun poly_lt  p q =
wenzelm@23252
   226
  case (p,q) of
wenzelm@23252
   227
    (p,[]) => false
wenzelm@23252
   228
  | ([],q) => true
wenzelm@23252
   229
  | ((c1,m1)::o1,(c2,m2)::o2) =>
wenzelm@23252
   230
        c1 </ c2 orelse
wenzelm@23252
   231
        c1 =/ c2 andalso ((morder_lt m1 m2) orelse m1 = m2 andalso poly_lt o1 o2);
wenzelm@23252
   232
wenzelm@23252
   233
fun align  ((p,hp),(q,hq)) =
wenzelm@23252
   234
  if poly_lt p q then ((p,hp),(q,hq)) else ((q,hq),(p,hp));
wenzelm@23252
   235
fun forall2 p l1 l2 =
wenzelm@23252
   236
  case (l1,l2) of
wenzelm@23252
   237
    ([],[]) => true
wenzelm@23252
   238
  | (h1::t1,h2::t2) => p h1 h2 andalso forall2 p t1 t2
wenzelm@23252
   239
  | _ => false;
wenzelm@23252
   240
wenzelm@23252
   241
fun poly_eq p1 p2 =
wenzelm@24630
   242
  forall2 (fn (c1,m1) => fn (c2,m2) => c1 =/ c2 andalso (m1: int list) = m2) p1 p2;
wenzelm@23252
   243
wenzelm@23252
   244
fun memx ((p1,h1),(p2,h2)) ppairs =
wenzelm@23252
   245
  not (exists (fn ((q1,_),(q2,_)) => poly_eq p1 q1 andalso poly_eq p2 q2) ppairs);
wenzelm@23252
   246
wenzelm@23252
   247
(* Buchberger's second criterion.                                            *)
wenzelm@23252
   248
wenzelm@23252
   249
fun criterion2 basis (lcm,((p1,h1),(p2,h2))) opairs =
wenzelm@23252
   250
  exists (fn g => not(poly_eq (fst g) p1) andalso not(poly_eq (fst g) p2) andalso
wenzelm@23252
   251
                   can (mdiv lcm) (hd(fst g)) andalso
wenzelm@23252
   252
                   not(memx (align (g,(p1,h1))) (map snd opairs)) andalso
wenzelm@23252
   253
                   not(memx (align (g,(p2,h2))) (map snd opairs))) basis;
wenzelm@23252
   254
wenzelm@23252
   255
(* Test for hitting constant polynomial.                                     *)
wenzelm@23252
   256
wenzelm@23252
   257
fun constant_poly p =
wenzelm@24630
   258
  length p = 1 andalso forall (fn x => x = 0) (snd(hd p));
wenzelm@23252
   259
wenzelm@23252
   260
(* Grobner basis algorithm.                                                  *)
chaieb@25251
   261
wenzelm@23252
   262
(* FIXME: try to get rid of mergesort? *)
wenzelm@23252
   263
fun merge ord l1 l2 =
wenzelm@23252
   264
 case l1 of
wenzelm@23252
   265
  [] => l2
wenzelm@23252
   266
 | h1::t1 =>
wenzelm@23252
   267
   case l2 of
wenzelm@23252
   268
    [] => l1
wenzelm@23252
   269
   | h2::t2 => if ord h1 h2 then h1::(merge ord t1 l2)
wenzelm@23252
   270
               else h2::(merge ord l1 t2);
wenzelm@23252
   271
fun mergesort ord l =
wenzelm@23252
   272
 let
wenzelm@23252
   273
 fun mergepairs l1 l2 =
wenzelm@23252
   274
  case (l1,l2) of
wenzelm@23252
   275
   ([s],[]) => s
wenzelm@23252
   276
 | (l,[]) => mergepairs [] l
wenzelm@23252
   277
 | (l,[s1]) => mergepairs (s1::l) []
wenzelm@23252
   278
 | (l,(s1::s2::ss)) => mergepairs ((merge ord s1 s2)::l) ss
wenzelm@23579
   279
 in if null l  then []  else mergepairs [] (map (fn x => [x]) l)
wenzelm@23252
   280
 end;
wenzelm@23252
   281
wenzelm@23252
   282
wenzelm@23252
   283
fun grobner_basis basis pairs =
chaieb@25251
   284
 case pairs of
chaieb@25251
   285
   [] => basis
chaieb@25251
   286
 | (l,(p1,p2))::opairs =>
chaieb@25251
   287
   let val (sph as (sp,hist)) = monic (reduce basis (spoly l p1 p2))
chaieb@25251
   288
   in 
chaieb@25251
   289
    if null sp orelse criterion2 basis (l,(p1,p2)) opairs
chaieb@25251
   290
    then grobner_basis basis opairs
chaieb@25251
   291
    else if constant_poly sp then grobner_basis (sph::basis) []
chaieb@25251
   292
    else 
chaieb@25251
   293
     let 
chaieb@25251
   294
      val rawcps = map (fn p => (mlcm (hd(fst p)) (hd sp),align(p,sph)))
chaieb@25251
   295
                              basis
chaieb@25251
   296
      val newcps = filter (fn (l,(p,q)) => not(orthogonal l (fst p) (fst q)))
chaieb@25251
   297
                        rawcps
chaieb@25251
   298
     in grobner_basis (sph::basis)
chaieb@25251
   299
                 (merge forder opairs (mergesort forder newcps))
chaieb@25251
   300
     end
chaieb@25251
   301
   end;
wenzelm@23252
   302
wenzelm@23252
   303
(* Interreduce initial polynomials.                                          *)
wenzelm@23252
   304
wenzelm@23252
   305
fun grobner_interreduce rpols ipols =
wenzelm@23252
   306
  case ipols of
wenzelm@23252
   307
    [] => map monic (rev rpols)
wenzelm@23252
   308
  | p::ps => let val p' = reduce (rpols @ ps) p in
wenzelm@23579
   309
             if null (fst p') then grobner_interreduce rpols ps
wenzelm@23252
   310
             else grobner_interreduce (p'::rpols) ps end;
wenzelm@23252
   311
wenzelm@23252
   312
(* Overall function.                                                         *)
wenzelm@23252
   313
wenzelm@23252
   314
fun grobner pols =
wenzelm@24630
   315
    let val npols = map2 (fn p => fn n => (p,Start n)) pols (0 upto (length pols - 1))
wenzelm@23579
   316
        val phists = filter (fn (p,_) => not (null p)) npols
wenzelm@23252
   317
        val bas = grobner_interreduce [] (map monic phists)
haftmann@25538
   318
        val prs0 = map_product pair bas bas
wenzelm@23252
   319
        val prs1 = filter (fn ((x,_),(y,_)) => poly_lt x y) prs0
wenzelm@23252
   320
        val prs2 = map (fn (p,q) => (mlcm (hd(fst p)) (hd(fst q)),(p,q))) prs1
wenzelm@23252
   321
        val prs3 =
wenzelm@23252
   322
            filter (fn (l,(p,q)) => not(orthogonal l (fst p) (fst q))) prs2 in
wenzelm@23252
   323
        grobner_basis bas (mergesort forder prs3) end;
wenzelm@23252
   324
wenzelm@23252
   325
(* Get proof of contradiction from Grobner basis.                            *)
chaieb@25251
   326
wenzelm@23252
   327
fun find p l =
wenzelm@23252
   328
  case l of
wenzelm@23252
   329
      [] => error "find"
wenzelm@23252
   330
    | (h::t) => if p(h) then h else find p t;
wenzelm@23252
   331
wenzelm@23252
   332
fun grobner_refute pols =
wenzelm@23252
   333
  let val gb = grobner pols in
wenzelm@23252
   334
  snd(find (fn (p,h) => length p = 1 andalso forall (fn x=> x=0) (snd(hd p))) gb)
wenzelm@23252
   335
  end;
wenzelm@23252
   336
wenzelm@23252
   337
(* Turn proof into a certificate as sum of multipliers.                      *)
wenzelm@23252
   338
(* In principle this is very inefficient: in a heavily shared proof it may   *)
wenzelm@23252
   339
(* make the same calculation many times. Could put in a cache or something.  *)
chaieb@25251
   340
wenzelm@23252
   341
fun resolve_proof vars prf =
wenzelm@23252
   342
  case prf of
wenzelm@23252
   343
    Start(~1) => []
wenzelm@23252
   344
  | Start m => [(m,[(rat_1,map (K 0) vars)])]
wenzelm@23252
   345
  | Mmul(pol,lin) =>
wenzelm@23252
   346
        let val lis = resolve_proof vars lin in
wenzelm@23252
   347
            map (fn (n,p) => (n,grob_cmul pol p)) lis end
wenzelm@23252
   348
  | Add(lin1,lin2) =>
wenzelm@23252
   349
        let val lis1 = resolve_proof vars lin1
wenzelm@23252
   350
            val lis2 = resolve_proof vars lin2
haftmann@33037
   351
            val dom = distinct (op =) (gen_union (op =) (map fst lis1, map fst lis2))
wenzelm@23252
   352
        in
wenzelm@23557
   353
            map (fn n => let val a = these (AList.lookup (op =) lis1 n)
wenzelm@23557
   354
                             val b = these (AList.lookup (op =) lis2 n)
wenzelm@23557
   355
                         in (n,grob_add a b) end) dom end;
wenzelm@23252
   356
wenzelm@23252
   357
(* Run the procedure and produce Weak Nullstellensatz certificate.           *)
chaieb@25251
   358
wenzelm@23252
   359
fun grobner_weak vars pols =
wenzelm@23252
   360
    let val cert = resolve_proof vars (grobner_refute pols)
wenzelm@23252
   361
        val l =
wenzelm@23252
   362
            fold_rev (fold_rev (lcm_rat o denominator_rat o fst) o snd) cert (rat_1) in
wenzelm@23252
   363
        (l,map (fn (i,p) => (i,map (fn (d,m) => (l*/d,m)) p)) cert) end;
wenzelm@23252
   364
wenzelm@23252
   365
(* Prove a polynomial is in ideal generated by others, using Grobner basis.  *)
wenzelm@23252
   366
wenzelm@23252
   367
fun grobner_ideal vars pols pol =
wenzelm@23252
   368
  let val (pol',h) = reduce (grobner pols) (grob_neg pol,Start(~1)) in
chaieb@24913
   369
  if not (null pol') then error "grobner_ideal: not in the ideal" else
wenzelm@23252
   370
  resolve_proof vars h end;
wenzelm@23252
   371
wenzelm@23252
   372
(* Produce Strong Nullstellensatz certificate for a power of pol.            *)
wenzelm@23252
   373
wenzelm@23252
   374
fun grobner_strong vars pols pol =
wenzelm@23252
   375
    let val vars' = @{cterm "True"}::vars
wenzelm@23252
   376
        val grob_z = [(rat_1,1::(map (fn x => 0) vars))]
wenzelm@23252
   377
        val grob_1 = [(rat_1,(map (fn x => 0) vars'))]
wenzelm@23252
   378
        fun augment p= map (fn (c,m) => (c,0::m)) p
wenzelm@23252
   379
        val pols' = map augment pols
wenzelm@23252
   380
        val pol' = augment pol
wenzelm@23252
   381
        val allpols = (grob_sub (grob_mul grob_z pol') grob_1)::pols'
wenzelm@23252
   382
        val (l,cert) = grobner_weak vars' allpols
wenzelm@24630
   383
        val d = fold_rev (fold_rev (curry Int.max o hd o snd) o snd) cert 0
wenzelm@23252
   384
        fun transform_monomial (c,m) =
wenzelm@23252
   385
            grob_cmul (c,tl m) (grob_pow vars pol (d - hd m))
wenzelm@23252
   386
        fun transform_polynomial q = fold_rev (grob_add o transform_monomial) q []
wenzelm@23252
   387
        val cert' = map (fn (c,q) => (c-1,transform_polynomial q))
wenzelm@23252
   388
                        (filter (fn (k,_) => k <> 0) cert) in
wenzelm@23252
   389
        (d,l,cert') end;
wenzelm@23252
   390
wenzelm@23252
   391
wenzelm@23252
   392
(* Overall parametrized universal procedure for (semi)rings.                 *)
wenzelm@23252
   393
(* We return an ideal_conv and the actual ring prover.                       *)
chaieb@25251
   394
wenzelm@23252
   395
fun refute_disj rfn tm =
wenzelm@23252
   396
 case term_of tm of
wenzelm@23252
   397
  Const("op |",_)$l$r =>
chaieb@25251
   398
   compose_single(refute_disj rfn (dest_arg tm),2,compose_single(refute_disj rfn (dest_arg1 tm),2,disjE))
wenzelm@23252
   399
  | _ => rfn tm ;
wenzelm@23252
   400
wenzelm@23252
   401
val notnotD = @{thm "notnotD"};
chaieb@25251
   402
fun mk_binop ct x y = capply (capply ct x) y
wenzelm@23252
   403
chaieb@25251
   404
val mk_comb = capply;
wenzelm@23252
   405
fun is_neg t =
wenzelm@23252
   406
    case term_of t of
wenzelm@23252
   407
      (Const("Not",_)$p) => true
wenzelm@23252
   408
    | _  => false;
wenzelm@23252
   409
fun is_eq t =
wenzelm@23252
   410
 case term_of t of
wenzelm@23252
   411
 (Const("op =",_)$_$_) => true
wenzelm@23252
   412
| _  => false;
wenzelm@23252
   413
wenzelm@23252
   414
fun end_itlist f l =
wenzelm@23252
   415
  case l of
wenzelm@23252
   416
        []     => error "end_itlist"
wenzelm@23252
   417
      | [x]    => x
wenzelm@23252
   418
      | (h::t) => f h (end_itlist f t);
wenzelm@23252
   419
wenzelm@23252
   420
val list_mk_binop = fn b => end_itlist (mk_binop b);
wenzelm@23252
   421
wenzelm@23252
   422
val list_dest_binop = fn b =>
wenzelm@23252
   423
 let fun h acc t =
chaieb@25251
   424
  ((let val (l,r) = dest_binary b t in h (h acc r) l end)
wenzelm@23252
   425
   handle CTERM _ => (t::acc)) (* Why had I handle _ => ? *)
wenzelm@23252
   426
 in h []
wenzelm@23252
   427
 end;
wenzelm@23252
   428
wenzelm@23252
   429
val strip_exists =
wenzelm@23252
   430
 let fun h (acc, t) =
wenzelm@23252
   431
      case (term_of t) of
chaieb@25251
   432
       Const("Ex",_)$Abs(x,T,p) => h (dest_abs NONE (dest_arg t) |>> (fn v => v::acc))
wenzelm@23252
   433
     | _ => (acc,t)
wenzelm@23252
   434
 in fn t => h ([],t)
wenzelm@23252
   435
 end;
wenzelm@23252
   436
wenzelm@23252
   437
fun is_forall t =
wenzelm@23252
   438
 case term_of t of
wenzelm@23252
   439
  (Const("All",_)$Abs(_,_,_)) => true
wenzelm@23252
   440
| _ => false;
wenzelm@23252
   441
wenzelm@23252
   442
val mk_object_eq = fn th => th COMP meta_eq_to_obj_eq;
wenzelm@23252
   443
val bool_simps = @{thms "bool_simps"};
wenzelm@23252
   444
val nnf_simps = @{thms "nnf_simps"};
wenzelm@23252
   445
val nnf_conv = Simplifier.rewrite (HOL_basic_ss addsimps bool_simps addsimps nnf_simps)
wenzelm@23557
   446
val weak_dnf_conv = Simplifier.rewrite (HOL_basic_ss addsimps @{thms "weak_dnf_simps"});
wenzelm@23252
   447
val initial_conv =
wenzelm@23252
   448
    Simplifier.rewrite
wenzelm@23252
   449
     (HOL_basic_ss addsimps nnf_simps
wenzelm@23252
   450
     addsimps [not_all, not_ex] addsimps map (fn th => th RS sym) (ex_simps @ all_simps));
wenzelm@23252
   451
wenzelm@23252
   452
val specl = fold_rev (fn x => fn th => instantiate' [] [SOME x] (th RS spec));
wenzelm@23252
   453
wenzelm@23252
   454
val cTrp = @{cterm "Trueprop"};
wenzelm@23252
   455
val cConj = @{cterm "op &"};
wenzelm@23252
   456
val (cNot,false_tm) = (@{cterm "Not"}, @{cterm "False"});
wenzelm@23557
   457
val assume_Trueprop = mk_comb cTrp #> assume;
wenzelm@23252
   458
val list_mk_conj = list_mk_binop cConj;
wenzelm@23252
   459
val conjs = list_dest_binop cConj;
wenzelm@23252
   460
val mk_neg = mk_comb cNot;
wenzelm@23252
   461
chaieb@25251
   462
fun striplist dest = 
chaieb@25251
   463
 let
chaieb@25251
   464
  fun h acc x = case try dest x of
chaieb@25251
   465
    SOME (a,b) => h (h acc b) a
chaieb@25251
   466
  | NONE => x::acc
chaieb@25251
   467
 in h [] end;
chaieb@25251
   468
fun list_mk_binop b = foldr1 (fn (s,t) => Thm.capply (Thm.capply b s) t);
chaieb@25251
   469
chaieb@25251
   470
val eq_commute = mk_meta_eq @{thm eq_commute};
chaieb@25251
   471
chaieb@25251
   472
fun sym_conv eq = 
chaieb@25251
   473
 let val (l,r) = Thm.dest_binop eq
chaieb@25251
   474
 in instantiate' [SOME (ctyp_of_term l)] [SOME l, SOME r] eq_commute
chaieb@25251
   475
 end;
chaieb@25251
   476
chaieb@25251
   477
  (* FIXME : copied from cqe.ML -- complex QE*)
chaieb@25251
   478
fun conjuncts ct =
chaieb@25251
   479
 case term_of ct of
chaieb@25251
   480
  @{term "op &"}$_$_ => (Thm.dest_arg1 ct)::(conjuncts (Thm.dest_arg ct))
chaieb@25251
   481
| _ => [ct];
chaieb@25251
   482
chaieb@25251
   483
fun fold1 f = foldr1 (uncurry f);
chaieb@25251
   484
chaieb@25251
   485
val list_conj = fold1 (fn c => fn c' => Thm.capply (Thm.capply @{cterm "op &"} c) c') ;
chaieb@25251
   486
chaieb@25251
   487
fun mk_conj_tab th = 
chaieb@25251
   488
 let fun h acc th = 
chaieb@25251
   489
   case prop_of th of
chaieb@25251
   490
   @{term "Trueprop"}$(@{term "op &"}$p$q) => 
chaieb@25251
   491
     h (h acc (th RS conjunct2)) (th RS conjunct1)
chaieb@25251
   492
  | @{term "Trueprop"}$p => (p,th)::acc
chaieb@25251
   493
in fold (Termtab.insert Thm.eq_thm) (h [] th) Termtab.empty end;
chaieb@25251
   494
chaieb@25251
   495
fun is_conj (@{term "op &"}$_$_) = true
chaieb@25251
   496
  | is_conj _ = false;
chaieb@25251
   497
chaieb@25251
   498
fun prove_conj tab cjs = 
chaieb@25251
   499
 case cjs of 
chaieb@25251
   500
   [c] => if is_conj (term_of c) then prove_conj tab (conjuncts c) else tab c
chaieb@25251
   501
 | c::cs => conjI OF [prove_conj tab [c], prove_conj tab cs];
chaieb@25251
   502
chaieb@25251
   503
fun conj_ac_rule eq = 
chaieb@25251
   504
 let 
chaieb@25251
   505
  val (l,r) = Thm.dest_equals eq
chaieb@25251
   506
  val ctabl = mk_conj_tab (assume (Thm.capply @{cterm Trueprop} l))
chaieb@25251
   507
  val ctabr = mk_conj_tab (assume (Thm.capply @{cterm Trueprop} r))
chaieb@25251
   508
  fun tabl c = valOf (Termtab.lookup ctabl (term_of c))
chaieb@25251
   509
  fun tabr c = valOf (Termtab.lookup ctabr (term_of c))
chaieb@25251
   510
  val thl  = prove_conj tabl (conjuncts r) |> implies_intr_hyps
chaieb@25251
   511
  val thr  = prove_conj tabr (conjuncts l) |> implies_intr_hyps
chaieb@25251
   512
  val eqI = instantiate' [] [SOME l, SOME r] @{thm iffI}
chaieb@25251
   513
 in implies_elim (implies_elim eqI thl) thr |> mk_meta_eq end;
chaieb@25251
   514
chaieb@25251
   515
 (* END FIXME.*)
chaieb@25251
   516
chaieb@25251
   517
   (* Conversion for the equivalence of existential statements where 
chaieb@25251
   518
      EX quantifiers are rearranged differently *)
chaieb@25251
   519
 fun ext T = cterm_rule (instantiate' [SOME T] []) @{cpat Ex}
chaieb@25251
   520
 fun mk_ex v t = Thm.capply (ext (ctyp_of_term v)) (Thm.cabs v t)
chaieb@25251
   521
chaieb@25251
   522
fun choose v th th' = case concl_of th of 
chaieb@25251
   523
  @{term Trueprop} $ (Const("Ex",_)$_) => 
chaieb@25251
   524
   let
chaieb@25251
   525
    val p = (funpow 2 Thm.dest_arg o cprop_of) th
chaieb@25251
   526
    val T = (hd o Thm.dest_ctyp o ctyp_of_term) p
chaieb@25251
   527
    val th0 = fconv_rule (Thm.beta_conversion true)
chaieb@25251
   528
        (instantiate' [SOME T] [SOME p, (SOME o Thm.dest_arg o cprop_of) th'] exE)
chaieb@25251
   529
    val pv = (Thm.rhs_of o Thm.beta_conversion true) 
chaieb@25251
   530
          (Thm.capply @{cterm Trueprop} (Thm.capply p v))
chaieb@25251
   531
    val th1 = forall_intr v (implies_intr pv th')
chaieb@25251
   532
   in implies_elim (implies_elim th0 th) th1  end
chaieb@25251
   533
| _ => error ""
chaieb@25251
   534
chaieb@25251
   535
fun simple_choose v th = 
chaieb@25251
   536
   choose v (assume ((Thm.capply @{cterm Trueprop} o mk_ex v) ((Thm.dest_arg o hd o #hyps o Thm.crep_thm) th))) th
chaieb@25251
   537
chaieb@25251
   538
chaieb@25251
   539
 fun mkexi v th = 
chaieb@25251
   540
  let 
chaieb@25251
   541
   val p = Thm.cabs v (Thm.dest_arg (Thm.cprop_of th))
chaieb@25251
   542
  in implies_elim 
chaieb@25251
   543
    (fconv_rule (Thm.beta_conversion true) (instantiate' [SOME (ctyp_of_term v)] [SOME p, SOME v] @{thm exI}))
chaieb@25251
   544
      th
chaieb@25251
   545
  end
chaieb@25251
   546
 fun ex_eq_conv t = 
chaieb@25251
   547
  let 
chaieb@25251
   548
  val (p0,q0) = Thm.dest_binop t
chaieb@25251
   549
  val (vs',P) = strip_exists p0 
chaieb@25251
   550
  val (vs,_) = strip_exists q0 
chaieb@25251
   551
   val th = assume (Thm.capply @{cterm Trueprop} P)
chaieb@25251
   552
   val th1 =  implies_intr_hyps (fold simple_choose vs' (fold mkexi vs th))
chaieb@25251
   553
   val th2 =  implies_intr_hyps (fold simple_choose vs (fold mkexi vs' th))
chaieb@25251
   554
   val p = (Thm.dest_arg o Thm.dest_arg1 o cprop_of) th1
chaieb@25251
   555
   val q = (Thm.dest_arg o Thm.dest_arg o cprop_of) th1
chaieb@25251
   556
  in implies_elim (implies_elim (instantiate' [] [SOME p, SOME q] iffI) th1) th2
chaieb@25251
   557
     |> mk_meta_eq
chaieb@25251
   558
  end;
chaieb@25251
   559
chaieb@25251
   560
chaieb@25251
   561
 fun getname v = case term_of v of 
chaieb@25251
   562
  Free(s,_) => s
chaieb@25251
   563
 | Var ((s,_),_) => s
chaieb@25251
   564
 | _ => "x"
chaieb@25251
   565
 fun mk_eq s t = Thm.capply (Thm.capply @{cterm "op == :: bool => _"} s) t
chaieb@25251
   566
 fun mkeq s t = Thm.capply @{cterm Trueprop} (Thm.capply (Thm.capply @{cterm "op = :: bool => _"} s) t)
chaieb@25251
   567
 fun mk_exists v th = arg_cong_rule (ext (ctyp_of_term v))
chaieb@25251
   568
   (Thm.abstract_rule (getname v) v th)
chaieb@25251
   569
 val simp_ex_conv = 
chaieb@25251
   570
     Simplifier.rewrite (HOL_basic_ss addsimps @{thms simp_thms(39)})
chaieb@25251
   571
chaieb@25251
   572
fun frees t = Thm.add_cterm_frees t [];
chaieb@25251
   573
fun free_in v t = member op aconvc (frees t) v;
chaieb@25251
   574
chaieb@25251
   575
val vsubst = let
chaieb@25251
   576
 fun vsubst (t,v) tm =  
chaieb@25251
   577
   (Thm.rhs_of o Thm.beta_conversion false) (Thm.capply (Thm.cabs v tm) t)
chaieb@25251
   578
in fold vsubst end;
wenzelm@23252
   579
wenzelm@23252
   580
wenzelm@23252
   581
(** main **)
wenzelm@23252
   582
wenzelm@23252
   583
fun ring_and_ideal_conv
chaieb@30866
   584
  {vars, semiring = (sr_ops, sr_rules), ring = (r_ops, r_rules), 
chaieb@30866
   585
   field = (f_ops, f_rules), idom, ideal}
wenzelm@23252
   586
  dest_const mk_const ring_eq_conv ring_normalize_conv =
wenzelm@23252
   587
let
wenzelm@23252
   588
  val [add_pat, mul_pat, pow_pat, zero_tm, one_tm] = sr_ops;
wenzelm@23252
   589
  val [ring_add_tm, ring_mul_tm, ring_pow_tm] =
chaieb@25251
   590
    map dest_fun2 [add_pat, mul_pat, pow_pat];
wenzelm@23252
   591
wenzelm@23252
   592
  val (ring_sub_tm, ring_neg_tm) =
wenzelm@23252
   593
    (case r_ops of
chaieb@30866
   594
     [sub_pat, neg_pat] => (dest_fun2 sub_pat, dest_fun neg_pat)
chaieb@30866
   595
    |_  => (@{cterm "True"}, @{cterm "True"}));
chaieb@30866
   596
chaieb@30866
   597
  val (field_div_tm, field_inv_tm) =
chaieb@30866
   598
    (case f_ops of
chaieb@30866
   599
       [div_pat, inv_pat] => (dest_fun2 div_pat, dest_fun inv_pat)
chaieb@30866
   600
     | _ => (@{cterm "True"}, @{cterm "True"}));
wenzelm@23252
   601
wenzelm@23252
   602
  val [idom_thm, neq_thm] = idom;
chaieb@25251
   603
  val [idl_sub, idl_add0] = 
chaieb@25251
   604
     if length ideal = 2 then ideal else [eq_commute, eq_commute]
chaieb@30866
   605
  fun ring_dest_neg t =
chaieb@30866
   606
    let val (l,r) = dest_comb t 
chaieb@30866
   607
    in if Term.could_unify(term_of l,term_of ring_neg_tm) then r 
chaieb@30866
   608
       else raise CTERM ("ring_dest_neg", [t])
chaieb@30866
   609
    end
wenzelm@23252
   610
wenzelm@23252
   611
 val ring_mk_neg = fn tm => mk_comb (ring_neg_tm) (tm);
chaieb@30866
   612
 fun field_dest_inv t =
chaieb@25251
   613
    let val (l,r) = dest_comb t in
chaieb@30866
   614
        if Term.could_unify(term_of l, term_of field_inv_tm) then r 
chaieb@30866
   615
        else raise CTERM ("field_dest_inv", [t])
wenzelm@23252
   616
    end
chaieb@25251
   617
 val ring_dest_add = dest_binary ring_add_tm;
wenzelm@23252
   618
 val ring_mk_add = mk_binop ring_add_tm;
chaieb@25251
   619
 val ring_dest_sub = dest_binary ring_sub_tm;
wenzelm@23252
   620
 val ring_mk_sub = mk_binop ring_sub_tm;
chaieb@25251
   621
 val ring_dest_mul = dest_binary ring_mul_tm;
wenzelm@23252
   622
 val ring_mk_mul = mk_binop ring_mul_tm;
chaieb@30866
   623
 val field_dest_div = dest_binary field_div_tm;
chaieb@30866
   624
 val field_mk_div = mk_binop field_div_tm;
chaieb@25251
   625
 val ring_dest_pow = dest_binary ring_pow_tm;
wenzelm@23252
   626
 val ring_mk_pow = mk_binop ring_pow_tm ;
wenzelm@23252
   627
 fun grobvars tm acc =
wenzelm@23252
   628
    if can dest_const tm then acc
chaieb@25251
   629
    else if can ring_dest_neg tm then grobvars (dest_arg tm) acc
chaieb@25251
   630
    else if can ring_dest_pow tm then grobvars (dest_arg1 tm) acc
wenzelm@23252
   631
    else if can ring_dest_add tm orelse can ring_dest_sub tm
wenzelm@23252
   632
            orelse can ring_dest_mul tm
chaieb@25251
   633
    then grobvars (dest_arg1 tm) (grobvars (dest_arg tm) acc)
chaieb@30866
   634
    else if can field_dest_inv tm
chaieb@30866
   635
         then
chaieb@30866
   636
          let val gvs = grobvars (dest_arg tm) [] 
chaieb@30866
   637
          in if null gvs then acc else tm::acc
chaieb@30866
   638
          end
chaieb@30866
   639
    else if can field_dest_div tm then
chaieb@30866
   640
         let val lvs = grobvars (dest_arg1 tm) acc
chaieb@30866
   641
             val gvs = grobvars (dest_arg tm) []
chaieb@30866
   642
          in if null gvs then lvs else tm::acc
chaieb@30866
   643
          end 
wenzelm@23252
   644
    else tm::acc ;
wenzelm@23252
   645
wenzelm@23252
   646
fun grobify_term vars tm =
wenzelm@23252
   647
((if not (member (op aconvc) vars tm) then raise CTERM ("Not a variable", [tm]) else
wenzelm@23252
   648
     [(rat_1,map (fn i => if i aconvc tm then 1 else 0) vars)])
wenzelm@23252
   649
handle  CTERM _ =>
wenzelm@23252
   650
 ((let val x = dest_const tm
wenzelm@23252
   651
 in if x =/ rat_0 then [] else [(x,map (fn v => 0) vars)]
wenzelm@23252
   652
 end)
wenzelm@23252
   653
 handle ERROR _ =>
wenzelm@23252
   654
  ((grob_neg(grobify_term vars (ring_dest_neg tm)))
wenzelm@23252
   655
  handle CTERM _ =>
wenzelm@23252
   656
   (
chaieb@30866
   657
   (grob_inv(grobify_term vars (field_dest_inv tm)))
chaieb@30866
   658
   handle CTERM _ => 
wenzelm@23252
   659
    ((let val (l,r) = ring_dest_add tm
wenzelm@23252
   660
    in grob_add (grobify_term vars l) (grobify_term vars r)
wenzelm@23252
   661
    end)
wenzelm@23252
   662
    handle CTERM _ =>
wenzelm@23252
   663
     ((let val (l,r) = ring_dest_sub tm
wenzelm@23252
   664
     in grob_sub (grobify_term vars l) (grobify_term vars r)
wenzelm@23252
   665
     end)
wenzelm@23252
   666
     handle  CTERM _ =>
wenzelm@23252
   667
      ((let val (l,r) = ring_dest_mul tm
wenzelm@23252
   668
      in grob_mul (grobify_term vars l) (grobify_term vars r)
wenzelm@23252
   669
      end)
wenzelm@23252
   670
       handle CTERM _ =>
chaieb@30866
   671
        (  (let val (l,r) = field_dest_div tm
wenzelm@23252
   672
          in grob_div (grobify_term vars l) (grobify_term vars r)
wenzelm@23252
   673
          end)
chaieb@30866
   674
         handle CTERM _ =>
wenzelm@23252
   675
          ((let val (l,r) = ring_dest_pow tm
wenzelm@23252
   676
          in grob_pow vars (grobify_term vars l) ((term_of #> HOLogic.dest_number #> snd) r)
wenzelm@23252
   677
          end)
wenzelm@23252
   678
           handle CTERM _ => error "grobify_term: unknown or invalid term")))))))));
chaieb@25251
   679
val eq_tm = idom_thm |> concl |> dest_arg |> dest_arg |> dest_fun2;
chaieb@25251
   680
val dest_eq = dest_binary eq_tm;
wenzelm@23252
   681
wenzelm@23252
   682
fun grobify_equation vars tm =
chaieb@25251
   683
    let val (l,r) = dest_binary eq_tm tm
wenzelm@23252
   684
    in grob_sub (grobify_term vars l) (grobify_term vars r)
wenzelm@23252
   685
    end;
wenzelm@23252
   686
wenzelm@23252
   687
fun grobify_equations tm =
wenzelm@23252
   688
 let
wenzelm@23252
   689
  val cjs = conjs tm
wenzelm@23252
   690
  val  rawvars = fold_rev (fn eq => fn a =>
chaieb@25251
   691
                                       grobvars (dest_arg1 eq) (grobvars (dest_arg eq) a)) cjs []
wenzelm@29269
   692
  val vars = sort (fn (x, y) => TermOrd.term_ord(term_of x,term_of y))
wenzelm@23252
   693
                  (distinct (op aconvc) rawvars)
wenzelm@23252
   694
 in (vars,map (grobify_equation vars) cjs)
wenzelm@23252
   695
 end;
wenzelm@23252
   696
wenzelm@23252
   697
val holify_polynomial =
wenzelm@23252
   698
 let fun holify_varpow (v,n) =
wenzelm@23579
   699
  if n = 1 then v else ring_mk_pow v (Numeral.mk_cnumber @{ctyp "nat"} n)  (* FIXME *)
wenzelm@23252
   700
 fun holify_monomial vars (c,m) =
wenzelm@24630
   701
  let val xps = map holify_varpow (filter (fn (_,n) => n <> 0) (vars ~~ m))
wenzelm@23252
   702
   in end_itlist ring_mk_mul (mk_const c :: xps)
wenzelm@23252
   703
  end
wenzelm@23252
   704
 fun holify_polynomial vars p =
wenzelm@23579
   705
     if null p then mk_const (rat_0)
wenzelm@23252
   706
     else end_itlist ring_mk_add (map (holify_monomial vars) p)
wenzelm@23252
   707
 in holify_polynomial
wenzelm@23252
   708
 end ;
wenzelm@23252
   709
val idom_rule = simplify (HOL_basic_ss addsimps [idom_thm]);
wenzelm@23252
   710
fun prove_nz n = eqF_elim
wenzelm@23252
   711
                 (ring_eq_conv(mk_binop eq_tm (mk_const n) (mk_const(rat_0))));
wenzelm@23252
   712
val neq_01 = prove_nz (rat_1);
wenzelm@23252
   713
fun neq_rule n th = [prove_nz n, th] MRS neq_thm;
chaieb@25251
   714
fun mk_add th1 = combination(arg_cong_rule ring_add_tm th1);
wenzelm@23252
   715
wenzelm@23252
   716
fun refute tm =
wenzelm@23557
   717
 if tm aconvc false_tm then assume_Trueprop tm else
chaieb@29800
   718
 ((let
wenzelm@23557
   719
   val (nths0,eths0) = List.partition (is_neg o concl) (HOLogic.conj_elims (assume_Trueprop tm))
chaieb@25251
   720
   val  nths = filter (is_eq o dest_arg o concl) nths0
wenzelm@23252
   721
   val eths = filter (is_eq o concl) eths0
wenzelm@23252
   722
  in
wenzelm@23252
   723
   if null eths then
wenzelm@23252
   724
    let
wenzelm@23557
   725
      val th1 = end_itlist (fn th1 => fn th2 => idom_rule(HOLogic.conj_intr th1 th2)) nths
wenzelm@23252
   726
      val th2 = Conv.fconv_rule
wenzelm@23252
   727
                ((arg_conv #> arg_conv)
wenzelm@23252
   728
                     (binop_conv ring_normalize_conv)) th1
chaieb@25251
   729
      val conc = th2 |> concl |> dest_arg
wenzelm@23252
   730
      val (l,r) = conc |> dest_eq
wenzelm@23252
   731
    in implies_intr (mk_comb cTrp tm)
chaieb@25251
   732
                    (equal_elim (arg_cong_rule cTrp (eqF_intr th2))
wenzelm@23252
   733
                           (reflexive l |> mk_object_eq))
wenzelm@23252
   734
    end
wenzelm@23252
   735
   else
wenzelm@23252
   736
   let
wenzelm@23252
   737
    val (vars,l,cert,noteqth) =(
wenzelm@23252
   738
     if null nths then
wenzelm@23252
   739
      let val (vars,pols) = grobify_equations(list_mk_conj(map concl eths))
wenzelm@23252
   740
          val (l,cert) = grobner_weak vars pols
wenzelm@23252
   741
      in (vars,l,cert,neq_01)
wenzelm@23252
   742
      end
wenzelm@23252
   743
     else
wenzelm@23252
   744
      let
wenzelm@23557
   745
       val nth = end_itlist (fn th1 => fn th2 => idom_rule(HOLogic.conj_intr th1 th2)) nths
wenzelm@23252
   746
       val (vars,pol::pols) =
chaieb@25251
   747
          grobify_equations(list_mk_conj(dest_arg(concl nth)::map concl eths))
wenzelm@23252
   748
       val (deg,l,cert) = grobner_strong vars pols pol
wenzelm@23252
   749
       val th1 = Conv.fconv_rule((arg_conv o arg_conv)(binop_conv ring_normalize_conv)) nth
wenzelm@24630
   750
       val th2 = funpow deg (idom_rule o HOLogic.conj_intr th1) neq_01
wenzelm@23252
   751
      in (vars,l,cert,th2)
wenzelm@23252
   752
      end)
wenzelm@23252
   753
    val cert_pos = map (fn (i,p) => (i,filter (fn (c,m) => c >/ rat_0) p)) cert
wenzelm@23252
   754
    val cert_neg = map (fn (i,p) => (i,map (fn (c,m) => (minus_rat c,m))
wenzelm@23252
   755
                                            (filter (fn (c,m) => c </ rat_0) p))) cert
wenzelm@23252
   756
    val  herts_pos = map (fn (i,p) => (i,holify_polynomial vars p)) cert_pos
wenzelm@23252
   757
    val  herts_neg = map (fn (i,p) => (i,holify_polynomial vars p)) cert_neg
wenzelm@23252
   758
    fun thm_fn pols =
wenzelm@23252
   759
        if null pols then reflexive(mk_const rat_0) else
wenzelm@23252
   760
        end_itlist mk_add
chaieb@25251
   761
            (map (fn (i,p) => arg_cong_rule (mk_comb ring_mul_tm p)
wenzelm@24630
   762
              (nth eths i |> mk_meta_eq)) pols)
wenzelm@23252
   763
    val th1 = thm_fn herts_pos
wenzelm@23252
   764
    val th2 = thm_fn herts_neg
wenzelm@23557
   765
    val th3 = HOLogic.conj_intr(mk_add (symmetric th1) th2 |> mk_object_eq) noteqth
wenzelm@23252
   766
    val th4 = Conv.fconv_rule ((arg_conv o arg_conv o binop_conv) ring_normalize_conv)
wenzelm@23252
   767
                               (neq_rule l th3)
chaieb@25251
   768
    val (l,r) = dest_eq(dest_arg(concl th4))
wenzelm@23252
   769
   in implies_intr (mk_comb cTrp tm)
chaieb@25251
   770
                        (equal_elim (arg_cong_rule cTrp (eqF_intr th4))
wenzelm@23252
   771
                   (reflexive l |> mk_object_eq))
wenzelm@23252
   772
   end
chaieb@29800
   773
  end) handle ERROR _ => raise CTERM ("Gorbner-refute: unable to refute",[tm]))
wenzelm@23252
   774
wenzelm@23252
   775
fun ring tm =
wenzelm@23252
   776
 let
wenzelm@23252
   777
  fun mk_forall x p =
chaieb@25251
   778
      mk_comb (cterm_rule (instantiate' [SOME (ctyp_of_term x)] []) @{cpat "All:: (?'a => bool) => _"}) (cabs x p)
chaieb@25251
   779
  val avs = add_cterm_frees tm []
wenzelm@23252
   780
  val P' = fold mk_forall avs tm
wenzelm@23252
   781
  val th1 = initial_conv(mk_neg P')
wenzelm@23252
   782
  val (evs,bod) = strip_exists(concl th1) in
chaieb@25251
   783
   if is_forall bod then raise CTERM("ring: non-universal formula",[tm])
wenzelm@23252
   784
   else
wenzelm@23252
   785
   let
wenzelm@23252
   786
    val th1a = weak_dnf_conv bod
wenzelm@23252
   787
    val boda = concl th1a
wenzelm@23252
   788
    val th2a = refute_disj refute boda
wenzelm@23252
   789
    val th2b = [mk_object_eq th1a, (th2a COMP notI) COMP PFalse'] MRS trans
wenzelm@23252
   790
    val th2 = fold (fn v => fn th => (forall_intr v th) COMP allI) evs (th2b RS PFalse)
wenzelm@23252
   791
    val th3 = equal_elim
wenzelm@23252
   792
                (Simplifier.rewrite (HOL_basic_ss addsimps [not_ex RS sym])
wenzelm@23252
   793
                          (th2 |> cprop_of)) th2
wenzelm@23252
   794
    in specl avs
wenzelm@23252
   795
             ([[[mk_object_eq th1, th3 RS PFalse'] MRS trans] MRS PFalse] MRS notnotD)
wenzelm@23252
   796
   end
wenzelm@23252
   797
 end
wenzelm@23252
   798
fun ideal tms tm ord =
wenzelm@23252
   799
 let
wenzelm@23252
   800
  val rawvars = fold_rev grobvars (tm::tms) []
wenzelm@23252
   801
  val vars = sort ord (distinct (fn (x,y) => (term_of x) aconv (term_of y)) rawvars)
wenzelm@23252
   802
  val pols = map (grobify_term vars) tms
wenzelm@23252
   803
  val pol = grobify_term vars tm
wenzelm@23252
   804
  val cert = grobner_ideal vars pols pol
wenzelm@23557
   805
 in map (fn n => these (AList.lookup (op =) cert n) |> holify_polynomial vars)
wenzelm@24630
   806
        (0 upto (length pols - 1))
wenzelm@23252
   807
 end
chaieb@25251
   808
chaieb@25251
   809
fun poly_eq_conv t = 
chaieb@25251
   810
 let val (a,b) = Thm.dest_binop t
chaieb@25251
   811
 in fconv_rule (arg_conv (arg1_conv ring_normalize_conv)) 
chaieb@25251
   812
     (instantiate' [] [SOME a, SOME b] idl_sub)
chaieb@25251
   813
 end
chaieb@25251
   814
 val poly_eq_simproc = 
chaieb@25251
   815
  let 
chaieb@25251
   816
   fun proc phi  ss t = 
chaieb@25251
   817
    let val th = poly_eq_conv t
chaieb@25251
   818
    in if Thm.is_reflexive th then NONE else SOME th
chaieb@25251
   819
    end
chaieb@25251
   820
   in make_simproc {lhss = [Thm.lhs_of idl_sub], 
chaieb@25251
   821
                name = "poly_eq_simproc", proc = proc, identifier = []}
chaieb@25251
   822
   end;
chaieb@25251
   823
  val poly_eq_ss = HOL_basic_ss addsimps simp_thms 
chaieb@25251
   824
                        addsimprocs [poly_eq_simproc]
chaieb@25251
   825
chaieb@25251
   826
 local
chaieb@25251
   827
  fun is_defined v t =
chaieb@25251
   828
  let 
chaieb@25251
   829
   val mons = striplist(dest_binary ring_add_tm) t 
chaieb@25251
   830
  in member (op aconvc) mons v andalso 
chaieb@25251
   831
    forall (fn m => v aconvc m 
chaieb@25251
   832
          orelse not(member (op aconvc) (Thm.add_cterm_frees m []) v)) mons
chaieb@25251
   833
  end
chaieb@25251
   834
chaieb@25251
   835
  fun isolate_variable vars tm =
chaieb@25251
   836
  let 
chaieb@25251
   837
   val th = poly_eq_conv tm
chaieb@25251
   838
   val th' = (sym_conv then_conv poly_eq_conv) tm
chaieb@25251
   839
   val (v,th1) = 
chaieb@25251
   840
   case find_first(fn v=> is_defined v (Thm.dest_arg1 (Thm.rhs_of th))) vars of
chaieb@25251
   841
    SOME v => (v,th')
chaieb@25251
   842
   | NONE => (valOf (find_first 
chaieb@25251
   843
          (fn v => is_defined v (Thm.dest_arg1 (Thm.rhs_of th'))) vars) ,th)
chaieb@25251
   844
   val th2 = transitive th1 
chaieb@25251
   845
        (instantiate' []  [(SOME o Thm.dest_arg1 o Thm.rhs_of) th1, SOME v] 
chaieb@25251
   846
          idl_add0)
chaieb@25251
   847
   in fconv_rule(funpow 2 arg_conv ring_normalize_conv) th2
chaieb@25251
   848
   end
chaieb@25251
   849
 in
chaieb@25251
   850
 fun unwind_polys_conv tm =
chaieb@25251
   851
 let 
chaieb@25251
   852
  val (vars,bod) = strip_exists tm
chaieb@25251
   853
  val cjs = striplist (dest_binary @{cterm "op &"}) bod
chaieb@25251
   854
  val th1 = (valOf (get_first (try (isolate_variable vars)) cjs) 
chaieb@25251
   855
             handle Option => raise CTERM ("unwind_polys_conv",[tm]))
chaieb@25251
   856
  val eq = Thm.lhs_of th1
chaieb@25251
   857
  val bod' = list_mk_binop @{cterm "op &"} (eq::(remove op aconvc eq cjs))
chaieb@25251
   858
  val th2 = conj_ac_rule (mk_eq bod bod')
chaieb@25251
   859
  val th3 = transitive th2 
chaieb@25251
   860
         (Drule.binop_cong_rule @{cterm "op &"} th1 
chaieb@25251
   861
                (reflexive (Thm.dest_arg (Thm.rhs_of th2))))
chaieb@25251
   862
  val v = Thm.dest_arg1(Thm.dest_arg1(Thm.rhs_of th3))
chaieb@25251
   863
  val vars' = (remove op aconvc v vars) @ [v]
chaieb@25251
   864
  val th4 = fconv_rule (arg_conv simp_ex_conv) (mk_exists v th3)
chaieb@25251
   865
  val th5 = ex_eq_conv (mk_eq tm (fold mk_ex (remove op aconvc v vars) (Thm.lhs_of th4)))
chaieb@25251
   866
 in transitive th5 (fold mk_exists (remove op aconvc v vars) th4)
chaieb@25251
   867
 end;
chaieb@25251
   868
end
chaieb@25251
   869
chaieb@25251
   870
local
chaieb@25251
   871
 fun scrub_var v m =
chaieb@25251
   872
  let 
chaieb@25251
   873
   val ps = striplist ring_dest_mul m 
chaieb@25251
   874
   val ps' = remove op aconvc v ps
chaieb@25251
   875
  in if null ps' then one_tm else fold1 ring_mk_mul ps'
chaieb@25251
   876
  end
chaieb@25251
   877
 fun find_multipliers v mons =
chaieb@25251
   878
  let 
chaieb@25251
   879
   val mons1 = filter (fn m => free_in v m) mons 
chaieb@25251
   880
   val mons2 = map (scrub_var v) mons1 
chaieb@25251
   881
   in  if null mons2 then zero_tm else fold1 ring_mk_add mons2
chaieb@25251
   882
  end
chaieb@25251
   883
chaieb@25251
   884
 fun isolate_monomials vars tm =
chaieb@25251
   885
 let 
chaieb@25251
   886
  val (cmons,vmons) =
chaieb@25251
   887
    List.partition (fn m => null (gen_inter op aconvc (frees m, vars)))
chaieb@25251
   888
                   (striplist ring_dest_add tm)
chaieb@25251
   889
  val cofactors = map (fn v => find_multipliers v vmons) vars
chaieb@25251
   890
  val cnc = if null cmons then zero_tm
chaieb@25251
   891
             else Thm.capply ring_neg_tm
chaieb@25251
   892
                    (list_mk_binop ring_add_tm cmons) 
chaieb@25251
   893
  in (cofactors,cnc)
chaieb@25251
   894
  end;
chaieb@25251
   895
chaieb@25251
   896
fun isolate_variables evs ps eq =
chaieb@25251
   897
 let 
chaieb@25251
   898
  val vars = filter (fn v => free_in v eq) evs
chaieb@25251
   899
  val (qs,p) = isolate_monomials vars eq
chaieb@25251
   900
  val rs = ideal (qs @ ps) p 
wenzelm@29269
   901
              (fn (s,t) => TermOrd.term_ord (term_of s, term_of t))
chaieb@25251
   902
 in (eq,Library.take (length qs, rs) ~~ vars)
chaieb@25251
   903
 end;
chaieb@25251
   904
 fun subst_in_poly i p = Thm.rhs_of (ring_normalize_conv (vsubst i p));
chaieb@25251
   905
in
chaieb@25251
   906
 fun solve_idealism evs ps eqs =
chaieb@25251
   907
  if null evs then [] else
chaieb@25251
   908
  let 
chaieb@25251
   909
   val (eq,cfs) = get_first (try (isolate_variables evs ps)) eqs |> valOf
chaieb@25251
   910
   val evs' = subtract op aconvc evs (map snd cfs)
chaieb@25251
   911
   val eqs' = map (subst_in_poly cfs) (remove op aconvc eq eqs)
chaieb@25251
   912
  in cfs @ solve_idealism evs' ps eqs'
chaieb@25251
   913
  end;
chaieb@25251
   914
end;
chaieb@25251
   915
chaieb@25251
   916
chaieb@25251
   917
in {ring_conv = ring, simple_ideal = ideal, multi_ideal = solve_idealism, 
chaieb@25251
   918
    poly_eq_ss = poly_eq_ss, unwind_conv = unwind_polys_conv}
wenzelm@23252
   919
end;
wenzelm@23252
   920
wenzelm@23252
   921
wenzelm@23252
   922
fun find_term bounds tm =
wenzelm@23252
   923
  (case term_of tm of
wenzelm@23252
   924
    Const ("op =", T) $ _ $ _ =>
wenzelm@23252
   925
      if domain_type T = HOLogic.boolT then find_args bounds tm
chaieb@25251
   926
      else dest_arg tm
chaieb@25251
   927
  | Const ("Not", _) $ _ => find_term bounds (dest_arg tm)
chaieb@25251
   928
  | Const ("All", _) $ _ => find_body bounds (dest_arg tm)
chaieb@25251
   929
  | Const ("Ex", _) $ _ => find_body bounds (dest_arg tm)
wenzelm@23252
   930
  | Const ("op &", _) $ _ $ _ => find_args bounds tm
wenzelm@23252
   931
  | Const ("op |", _) $ _ $ _ => find_args bounds tm
wenzelm@23252
   932
  | Const ("op -->", _) $ _ $ _ => find_args bounds tm
chaieb@25251
   933
  | @{term "op ==>"} $_$_ => find_args bounds tm
chaieb@25251
   934
  | Const("op ==",_)$_$_ => find_args bounds tm
chaieb@25251
   935
  | @{term Trueprop}$_ => find_term bounds (dest_arg tm)
wenzelm@23252
   936
  | _ => raise TERM ("find_term", []))
wenzelm@23252
   937
and find_args bounds tm =
wenzelm@23252
   938
  let val (t, u) = Thm.dest_binop tm
wenzelm@23252
   939
  in (find_term bounds t handle TERM _ => find_term bounds u) end
wenzelm@23252
   940
and find_body bounds b =
chaieb@25251
   941
  let val (_, b') = dest_abs (SOME (Name.bound bounds)) b
wenzelm@23252
   942
  in find_term (bounds + 1) b' end;
wenzelm@23252
   943
chaieb@25251
   944
chaieb@25251
   945
fun get_ring_ideal_convs ctxt form = 
chaieb@25251
   946
 case try (find_term 0) form of
chaieb@25251
   947
  NONE => NONE
chaieb@25251
   948
| SOME tm =>
chaieb@25251
   949
  (case NormalizerData.match ctxt tm of
chaieb@25251
   950
    NONE => NONE
chaieb@25251
   951
  | SOME (res as (theory, {is_const, dest_const, 
chaieb@25251
   952
          mk_const, conv = ring_eq_conv})) =>
chaieb@25251
   953
     SOME (ring_and_ideal_conv theory
chaieb@25251
   954
          dest_const (mk_const (ctyp_of_term tm)) (ring_eq_conv ctxt)
chaieb@25251
   955
          (semiring_normalize_wrapper ctxt res)))
chaieb@25251
   956
wenzelm@23579
   957
fun ring_solve ctxt form =
wenzelm@23252
   958
  (case try (find_term 0 (* FIXME !? *)) form of
wenzelm@23252
   959
    NONE => reflexive form
wenzelm@23252
   960
  | SOME tm =>
wenzelm@23252
   961
      (case NormalizerData.match ctxt tm of
wenzelm@23252
   962
        NONE => reflexive form
wenzelm@23252
   963
      | SOME (res as (theory, {is_const, dest_const, mk_const, conv = ring_eq_conv})) =>
chaieb@25251
   964
        #ring_conv (ring_and_ideal_conv theory
chaieb@25251
   965
          dest_const (mk_const (ctyp_of_term tm)) (ring_eq_conv ctxt)
chaieb@23330
   966
          (semiring_normalize_wrapper ctxt res)) form));
wenzelm@23252
   967
wenzelm@23579
   968
fun ring_tac add_ths del_ths ctxt =
wenzelm@23579
   969
  ObjectLogic.full_atomize_tac
chaieb@27671
   970
  THEN' asm_full_simp_tac (Simplifier.context ctxt (fst (NormalizerData.get ctxt)) delsimps del_ths addsimps add_ths)
wenzelm@23579
   971
  THEN' CSUBGOAL (fn (p, i) =>
chaieb@25251
   972
    rtac (let val form = (ObjectLogic.dest_judgment p)
chaieb@25251
   973
          in case get_ring_ideal_convs ctxt form of
chaieb@25251
   974
           NONE => reflexive form
chaieb@25251
   975
          | SOME thy => #ring_conv thy form
chaieb@25251
   976
          end) i
wenzelm@23579
   977
      handle TERM _ => no_tac
wenzelm@23579
   978
        | CTERM _ => no_tac
wenzelm@23579
   979
        | THM _ => no_tac);
chaieb@23334
   980
chaieb@25251
   981
local
chaieb@25251
   982
 fun lhs t = case term_of t of
chaieb@25251
   983
  Const("op =",_)$_$_ => Thm.dest_arg1 t
chaieb@25251
   984
 | _=> raise CTERM ("ideal_tac - lhs",[t])
chaieb@25251
   985
 fun exitac NONE = no_tac
chaieb@25251
   986
   | exitac (SOME y) = rtac (instantiate' [SOME (ctyp_of_term y)] [NONE,SOME y] exI) 1
chaieb@25251
   987
in 
chaieb@25251
   988
fun ideal_tac add_ths del_ths ctxt = 
chaieb@27671
   989
  asm_full_simp_tac 
chaieb@27671
   990
   (Simplifier.context ctxt (fst (NormalizerData.get ctxt)) delsimps del_ths addsimps add_ths) 
chaieb@27671
   991
 THEN'
chaieb@25251
   992
 CSUBGOAL (fn (p, i) =>
chaieb@25251
   993
  case get_ring_ideal_convs ctxt p of
chaieb@25251
   994
   NONE => no_tac
chaieb@25251
   995
 | SOME thy => 
chaieb@25251
   996
  let
chaieb@25251
   997
   fun poly_exists_tac {asms = asms, concl = concl, prems = prems,
chaieb@25251
   998
            params = params, context = ctxt, schematics = scs} = 
chaieb@25251
   999
    let
chaieb@25251
  1000
     val (evs,bod) = strip_exists (Thm.dest_arg concl)
chaieb@25251
  1001
     val ps = map_filter (try (lhs o Thm.dest_arg)) asms 
chaieb@25251
  1002
     val cfs = (map swap o #multi_ideal thy evs ps) 
chaieb@25251
  1003
                   (map Thm.dest_arg1 (conjuncts bod))
chaieb@25251
  1004
     val ws = map (exitac o AList.lookup op aconvc cfs) evs
chaieb@25251
  1005
    in EVERY (rev ws) THEN Method.insert_tac prems 1 
chaieb@25251
  1006
        THEN ring_tac add_ths del_ths ctxt 1
chaieb@25251
  1007
   end
chaieb@27671
  1008
  in  
chaieb@27671
  1009
     clarify_tac @{claset} i 
chaieb@27671
  1010
     THEN ObjectLogic.full_atomize_tac i 
chaieb@27671
  1011
     THEN asm_full_simp_tac (Simplifier.context ctxt (#poly_eq_ss thy)) i 
chaieb@27671
  1012
     THEN clarify_tac @{claset} i 
chaieb@25251
  1013
     THEN (REPEAT (CONVERSION (#unwind_conv thy) i))
chaieb@25251
  1014
     THEN SUBPROOF poly_exists_tac ctxt i
chaieb@25251
  1015
  end
chaieb@25251
  1016
 handle TERM _ => no_tac
chaieb@25251
  1017
     | CTERM _ => no_tac
chaieb@25251
  1018
     | THM _ => no_tac); 
wenzelm@23252
  1019
end;
chaieb@25251
  1020
chaieb@25251
  1021
fun algebra_tac add_ths del_ths ctxt i = 
chaieb@25251
  1022
 ring_tac add_ths del_ths ctxt i ORELSE ideal_tac add_ths del_ths ctxt i
chaieb@25251
  1023
 
chaieb@25251
  1024
 
chaieb@25251
  1025
chaieb@25251
  1026
end;