src/HOL/Library/Sum_of_Squares/sum_of_squares.ML
author wenzelm
Wed, 01 Jun 2016 15:10:27 +0200
changeset 63205 97b721666890
parent 63201 f151704c08e4
child 63211 0bec0d1d9998
permissions -rw-r--r--
prefer rat numberals;
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
41474
60d091240485 renamed Sum_Of_Squares to Sum_of_Squares;
wenzelm
parents: 40721
diff changeset
     1
(*  Title:      HOL/Library/Sum_of_Squares/sum_of_squares.ML
32949
aa6c470a962a eliminated slightly odd get/set operations in favour of Unsynchronized.ref;
wenzelm
parents: 32839
diff changeset
     2
    Author:     Amine Chaieb, University of Cambridge
aa6c470a962a eliminated slightly odd get/set operations in favour of Unsynchronized.ref;
wenzelm
parents: 32839
diff changeset
     3
    Author:     Philipp Meyer, TU Muenchen
32268
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
     4
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
     5
A tactic for proving nonlinear inequalities.
32268
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
     6
*)
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
     7
32949
aa6c470a962a eliminated slightly odd get/set operations in favour of Unsynchronized.ref;
wenzelm
parents: 32839
diff changeset
     8
signature SUM_OF_SQUARES =
32268
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
     9
sig
32949
aa6c470a962a eliminated slightly odd get/set operations in favour of Unsynchronized.ref;
wenzelm
parents: 32839
diff changeset
    10
  datatype proof_method = Certificate of RealArith.pss_tree | Prover of string -> string
aa6c470a962a eliminated slightly odd get/set operations in favour of Unsynchronized.ref;
wenzelm
parents: 32839
diff changeset
    11
  val sos_tac: (RealArith.pss_tree -> unit) -> proof_method -> Proof.context -> int -> tactic
38805
b09d8603f865 Sum_Of_Squares: proper configuration options;
wenzelm
parents: 38795
diff changeset
    12
  val trace: bool Config.T
58631
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
    13
  val debug: bool Config.T
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
    14
  val trace_message: Proof.context -> (unit -> string) -> unit
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
    15
  val debug_message: Proof.context -> (unit -> string) -> unit
32332
bc5cec7b2be6 misc changes to SOS by Philipp Meyer:
wenzelm
parents: 32268
diff changeset
    16
  exception Failure of string;
32268
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
    17
end
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
    18
41474
60d091240485 renamed Sum_Of_Squares to Sum_of_Squares;
wenzelm
parents: 40721
diff changeset
    19
structure Sum_of_Squares: SUM_OF_SQUARES =
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    20
struct
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    21
33029
2fefe039edf1 uniform use of Integer.min/max;
wenzelm
parents: 33002
diff changeset
    22
val max = Integer.max;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    23
63201
f151704c08e4 tuned signature;
wenzelm
parents: 63198
diff changeset
    24
val denominator_rat = Rat.dest #> snd #> Rat.of_int;
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    25
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
    26
fun int_of_rat a =
63201
f151704c08e4 tuned signature;
wenzelm
parents: 63198
diff changeset
    27
  (case Rat.dest a of
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    28
    (i, 1) => i
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    29
  | _ => error "int_of_rat: not an int");
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    30
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    31
fun lcm_rat x y =
63201
f151704c08e4 tuned signature;
wenzelm
parents: 63198
diff changeset
    32
  Rat.of_int (Integer.lcm (int_of_rat x) (int_of_rat y));
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    33
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
    34
fun rat_pow r i =
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
    35
 let fun pow r i =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
    36
   if i = 0 then @1 else
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    37
   let val d = pow r (i div 2)
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
    38
   in d * d * (if i mod 2 = 0 then @1 else r)
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    39
   end
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    40
 in if i < 0 then pow (Rat.inv r) (~ i) else pow r i end;
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    41
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
    42
fun round_rat r =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    43
  let
63201
f151704c08e4 tuned signature;
wenzelm
parents: 63198
diff changeset
    44
    val (a,b) = Rat.dest (Rat.abs r)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    45
    val d = a div b
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
    46
    val s = if r < @0 then (Rat.neg o Rat.of_int) else Rat.of_int
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    47
    val x2 = 2 * (a - (b * d))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    48
  in s (if x2 >= b then d + 1 else d) end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    49
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    50
val abs_rat = Rat.abs;
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
    51
val pow2 = rat_pow @2;
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
    52
val pow10 = rat_pow @10;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    53
58631
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
    54
42616
92715b528e78 added Attrib.setup_config_XXX conveniences, with implicit setup of the background theory;
wenzelm
parents: 41490
diff changeset
    55
val trace = Attrib.setup_config_bool @{binding sos_trace} (K false);
58631
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
    56
val debug = Attrib.setup_config_bool @{binding sos_debug} (K false);
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
    57
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
    58
fun trace_message ctxt msg =
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
    59
  if Config.get ctxt trace orelse Config.get ctxt debug then tracing (msg ()) else ();
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
    60
fun debug_message ctxt msg = if Config.get ctxt debug then tracing (msg ()) else ();
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    61
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    62
exception Sanity;
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    63
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    64
exception Unsolvable;
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    65
32332
bc5cec7b2be6 misc changes to SOS by Philipp Meyer:
wenzelm
parents: 32268
diff changeset
    66
exception Failure of string;
bc5cec7b2be6 misc changes to SOS by Philipp Meyer:
wenzelm
parents: 32268
diff changeset
    67
32645
1cc5b24f5a01 sos method generates and uses proof certificates
Philipp Meyer
parents: 32332
diff changeset
    68
datatype proof_method =
1cc5b24f5a01 sos method generates and uses proof certificates
Philipp Meyer
parents: 32332
diff changeset
    69
    Certificate of RealArith.pss_tree
1cc5b24f5a01 sos method generates and uses proof certificates
Philipp Meyer
parents: 32332
diff changeset
    70
  | Prover of (string -> string)
1cc5b24f5a01 sos method generates and uses proof certificates
Philipp Meyer
parents: 32332
diff changeset
    71
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    72
(* Turn a rational into a decimal string with d sig digits.                  *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    73
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    74
local
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    75
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    76
fun normalize y =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
    77
  if abs_rat y < @1/10 then normalize (@10 * y) - 1
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
    78
  else if abs_rat y >= @1 then normalize (y / @10) + 1
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
    79
  else 0
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    80
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    81
in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    82
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    83
fun decimalize d x =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
    84
  if x = @0 then "0.0"
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    85
  else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    86
    let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    87
      val y = Rat.abs x
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    88
      val e = normalize y
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
    89
      val z = pow10(~ e) * y + @1
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
    90
      val k = int_of_rat (round_rat(pow10 d * z))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    91
    in
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
    92
      (if x < @0 then "-0." else "0.") ^
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    93
      implode (tl (raw_explode(string_of_int k))) ^
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    94
      (if e = 0 then "" else "e" ^ string_of_int e)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    95
    end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
    96
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    97
end;
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    98
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
    99
(* Iterations over numbers, and lists indexed by numbers.                    *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   100
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   101
fun itern k l f a =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   102
  (case l of
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   103
    [] => a
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   104
  | h::t => itern (k + 1) t f (f h k a));
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   105
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   106
fun iter (m,n) f a =
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   107
  if n < m then a
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   108
  else iter (m + 1, n) f (f m a);
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   109
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   110
(* The main types.                                                           *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   111
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   112
type vector = int * Rat.rat FuncUtil.Intfunc.table;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   113
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   114
type matrix = (int * int) * Rat.rat FuncUtil.Intpairfunc.table;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   115
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   116
fun iszero (_, r) = r = @0;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   117
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   118
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   119
(* Vectors. Conventionally indexed 1..n.                                     *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   120
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   121
fun vector_0 n = (n, FuncUtil.Intfunc.empty): vector;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   122
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   123
fun dim (v: vector) = fst v;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   124
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   125
fun vector_cmul c (v: vector) =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   126
  let val n = dim v in
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   127
    if c = @0 then vector_0 n
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   128
    else (n,FuncUtil.Intfunc.map (fn _ => fn x => c * x) (snd v))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   129
  end;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   130
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   131
fun vector_of_list l =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   132
  let val n = length l in
58635
010b610eb55d eliminated some exotic combinators;
wenzelm
parents: 58631
diff changeset
   133
    (n, fold_rev FuncUtil.Intfunc.update (1 upto n ~~ l) FuncUtil.Intfunc.empty): vector
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   134
  end;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   135
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   136
(* Matrices; again rows and columns indexed from 1.                          *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   137
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   138
fun dimensions (m: matrix) = fst m;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   139
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   140
fun row k (m: matrix) : vector =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   141
  let val (_, j) = dimensions m in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   142
    (j,
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   143
      FuncUtil.Intpairfunc.fold (fn ((i, j), c) => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   144
        if i = k then FuncUtil.Intfunc.update (j, c) a else a) (snd m) FuncUtil.Intfunc.empty)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   145
  end;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   146
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   147
(* Monomials.                                                                *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   148
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   149
fun monomial_eval assig m =
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   150
  FuncUtil.Ctermfunc.fold (fn (x, k) => fn a => a * rat_pow (FuncUtil.Ctermfunc.apply assig x) k)
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   151
    m @1;
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   152
32829
671eb46eb0a3 tuned FuncFun and FuncUtil structure in positivstellensatz.ML
Philipp Meyer
parents: 32828
diff changeset
   153
val monomial_1 = FuncUtil.Ctermfunc.empty;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   154
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   155
fun monomial_var x = FuncUtil.Ctermfunc.onefunc (x, 1);
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   156
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   157
val monomial_mul =
33002
f3f02f36a3e2 uniform use of Integer.add/mult/sum/prod;
wenzelm
parents: 32949
diff changeset
   158
  FuncUtil.Ctermfunc.combine Integer.add (K false);
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   159
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   160
fun monomial_multidegree m =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   161
  FuncUtil.Ctermfunc.fold (fn (_, k) => fn a => k + a) m 0;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   162
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   163
fun monomial_variables m = FuncUtil.Ctermfunc.dom m;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   164
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   165
(* Polynomials.                                                              *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   166
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   167
fun eval assig p =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   168
  FuncUtil.Monomialfunc.fold (fn (m, c) => fn a => a + c * monomial_eval assig m) p @0;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   169
32829
671eb46eb0a3 tuned FuncFun and FuncUtil structure in positivstellensatz.ML
Philipp Meyer
parents: 32828
diff changeset
   170
val poly_0 = FuncUtil.Monomialfunc.empty;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   171
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   172
fun poly_isconst p =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   173
  FuncUtil.Monomialfunc.fold (fn (m, _) => fn a => FuncUtil.Ctermfunc.is_empty m andalso a)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   174
    p true;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   175
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   176
fun poly_var x = FuncUtil.Monomialfunc.onefunc (monomial_var x, @1);
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   177
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   178
fun poly_const c =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   179
  if c = @0 then poly_0 else FuncUtil.Monomialfunc.onefunc (monomial_1, c);
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   180
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   181
fun poly_cmul c p =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   182
  if c = @0 then poly_0
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   183
  else FuncUtil.Monomialfunc.map (fn _ => fn x => c * x) p;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   184
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   185
fun poly_neg p = FuncUtil.Monomialfunc.map (K Rat.neg) p;
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   186
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   187
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   188
fun poly_add p1 p2 =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   189
  FuncUtil.Monomialfunc.combine (curry op +) (fn x => x = @0) p1 p2;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   190
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   191
fun poly_sub p1 p2 = poly_add p1 (poly_neg p2);
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   192
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   193
fun poly_cmmul (c,m) p =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   194
  if c = @0 then poly_0
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   195
  else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   196
    if FuncUtil.Ctermfunc.is_empty m
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   197
    then FuncUtil.Monomialfunc.map (fn _ => fn d => c * d) p
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   198
    else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   199
      FuncUtil.Monomialfunc.fold (fn (m', d) => fn a =>
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   200
          (FuncUtil.Monomialfunc.update (monomial_mul m m', c * d) a)) p poly_0;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   201
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   202
fun poly_mul p1 p2 =
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   203
  FuncUtil.Monomialfunc.fold (fn (m, c) => fn a => poly_add (poly_cmmul (c,m) p2) a) p1 poly_0;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   204
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   205
fun poly_square p = poly_mul p p;
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   206
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   207
fun poly_pow p k =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   208
  if k = 0 then poly_const @1
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   209
  else if k = 1 then p
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   210
  else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   211
    let val q = poly_square(poly_pow p (k div 2))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   212
    in if k mod 2 = 1 then poly_mul p q else q end;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   213
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   214
fun multidegree p =
44453
082edd97ffe1 comment out dead code to avoid compiler warnings
huffman
parents: 42616
diff changeset
   215
  FuncUtil.Monomialfunc.fold (fn (m, _) => fn a => max (monomial_multidegree m) a) p 0;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   216
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   217
fun poly_variables p =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   218
  sort FuncUtil.cterm_ord
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   219
    (FuncUtil.Monomialfunc.fold_rev
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   220
      (fn (m, _) => union (is_equal o FuncUtil.cterm_ord) (monomial_variables m)) p []);
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   221
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   222
(* Conversion from HOL term.                                                 *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   223
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   224
local
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   225
  val neg_tm = @{cterm "uminus :: real => _"}
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   226
  val add_tm = @{cterm "op + :: real => _"}
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   227
  val sub_tm = @{cterm "op - :: real => _"}
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   228
  val mul_tm = @{cterm "op * :: real => _"}
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   229
  val inv_tm = @{cterm "inverse :: real => _"}
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   230
  val div_tm = @{cterm "op / :: real => _"}
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   231
  val pow_tm = @{cterm "op ^ :: real => _"}
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   232
  val zero_tm = @{cterm "0:: real"}
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   233
  val is_numeral = can (HOLogic.dest_number o Thm.term_of)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   234
  fun poly_of_term tm =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   235
    if tm aconvc zero_tm then poly_0
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   236
    else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   237
      if RealArith.is_ratconst tm
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   238
      then poly_const(RealArith.dest_ratconst tm)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   239
      else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   240
       (let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   241
          val (lop, r) = Thm.dest_comb tm
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   242
        in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   243
          if lop aconvc neg_tm then poly_neg(poly_of_term r)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   244
          else if lop aconvc inv_tm then
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   245
            let val p = poly_of_term r in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   246
              if poly_isconst p
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   247
              then poly_const(Rat.inv (eval FuncUtil.Ctermfunc.empty p))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   248
              else error "poly_of_term: inverse of non-constant polyomial"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   249
            end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   250
          else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   251
           (let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   252
              val (opr,l) = Thm.dest_comb lop
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   253
            in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   254
              if opr aconvc pow_tm andalso is_numeral r
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   255
              then poly_pow (poly_of_term l) ((snd o HOLogic.dest_number o Thm.term_of) r)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   256
              else if opr aconvc add_tm
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   257
              then poly_add (poly_of_term l) (poly_of_term r)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   258
              else if opr aconvc sub_tm
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   259
              then poly_sub (poly_of_term l) (poly_of_term r)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   260
              else if opr aconvc mul_tm
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   261
              then poly_mul (poly_of_term l) (poly_of_term r)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   262
              else if opr aconvc div_tm
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   263
              then
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   264
                let
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   265
                  val p = poly_of_term l
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   266
                  val q = poly_of_term r
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   267
                in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   268
                  if poly_isconst q
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   269
                  then poly_cmul (Rat.inv (eval FuncUtil.Ctermfunc.empty q)) p
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   270
                  else error "poly_of_term: division by non-constant polynomial"
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   271
                end
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   272
              else poly_var tm
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   273
            end handle CTERM ("dest_comb",_) => poly_var tm)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   274
        end handle CTERM ("dest_comb",_) => poly_var tm)
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   275
in
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   276
  val poly_of_term = fn tm =>
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   277
    if type_of (Thm.term_of tm) = @{typ real}
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   278
    then poly_of_term tm
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   279
    else error "poly_of_term: term does not have real type"
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   280
end;
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   281
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   282
(* String of vector (just a list of space-separated numbers).                *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   283
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   284
fun sdpa_of_vector (v: vector) =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   285
  let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   286
    val n = dim v
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   287
    val strs =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   288
      map (decimalize 20 o (fn i => FuncUtil.Intfunc.tryapplyd (snd v) i @0)) (1 upto n)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   289
  in space_implode " " strs ^ "\n" end;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   290
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   291
fun triple_int_ord ((a, b, c), (a', b', c')) =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   292
  prod_ord int_ord (prod_ord int_ord int_ord) ((a, (b, c)), (a', (b', c')));
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   293
structure Inttriplefunc = FuncFun(type key = int * int * int val ord = triple_int_ord);
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   294
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   295
fun index_char str chr pos =
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   296
  if pos >= String.size str then ~1
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   297
  else if String.sub(str,pos) = chr then pos
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   298
  else index_char str chr (pos + 1);
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   299
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   300
fun rat_of_quotient (a,b) =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   301
  if b = 0 then @0 else Rat.make (a, b);
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   302
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   303
fun rat_of_string s =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   304
  let val n = index_char s #"/" 0 in
63201
f151704c08e4 tuned signature;
wenzelm
parents: 63198
diff changeset
   305
    if n = ~1 then s |> Int.fromString |> the |> Rat.of_int
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   306
    else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   307
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   308
        val SOME numer = Int.fromString(String.substring(s,0,n))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   309
        val SOME den = Int.fromString (String.substring(s,n+1,String.size s - n - 1))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   310
      in rat_of_quotient(numer, den) end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   311
  end;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   312
53975
22ee3fb9d596 backout c6297fa1031a -- strange parsers are required to make this work;
wenzelm
parents: 53972
diff changeset
   313
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   314
fun isnum x = member (op =) ["0", "1", "2", "3", "4", "5", "6", "7", "8", "9"] x;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   315
53975
22ee3fb9d596 backout c6297fa1031a -- strange parsers are required to make this work;
wenzelm
parents: 53972
diff changeset
   316
(* More parser basics. *)
22ee3fb9d596 backout c6297fa1031a -- strange parsers are required to make this work;
wenzelm
parents: 53972
diff changeset
   317
(* FIXME improper use of parser combinators ahead *)
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   318
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   319
val numeral = Scan.one isnum
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   320
val decimalint = Scan.repeat1 numeral >> (rat_of_string o implode)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   321
val decimalfrac = Scan.repeat1 numeral
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   322
  >> (fn s => rat_of_string(implode s) / pow10 (length s))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   323
val decimalsig =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   324
  decimalint -- Scan.option (Scan.$$ "." |-- decimalfrac)
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   325
  >> (fn (h,NONE) => h | (h,SOME x) => h + x)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   326
fun signed prs =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   327
     $$ "-" |-- prs >> Rat.neg
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   328
  || $$ "+" |-- prs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   329
  || prs;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   330
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   331
fun emptyin def xs = if null xs then (def, xs) else Scan.fail xs
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   332
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   333
val exponent = ($$ "e" || $$ "E") |-- signed decimalint;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   334
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   335
val decimal = signed decimalsig -- (emptyin @0|| exponent)
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   336
  >> (fn (h, x) => h * pow10 (int_of_rat x));
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   337
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   338
fun mkparser p s =
53975
22ee3fb9d596 backout c6297fa1031a -- strange parsers are required to make this work;
wenzelm
parents: 53972
diff changeset
   339
  let val (x,rst) = p (raw_explode s)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   340
  in if null rst then x else error "mkparser: unparsed input" end;
53975
22ee3fb9d596 backout c6297fa1031a -- strange parsers are required to make this work;
wenzelm
parents: 53972
diff changeset
   341
32332
bc5cec7b2be6 misc changes to SOS by Philipp Meyer:
wenzelm
parents: 32268
diff changeset
   342
(* Parse back csdp output.                                                      *)
53975
22ee3fb9d596 backout c6297fa1031a -- strange parsers are required to make this work;
wenzelm
parents: 53972
diff changeset
   343
(* FIXME improper use of parser combinators ahead *)
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   344
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   345
fun ignore _ = ((),[])
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   346
fun csdpoutput inp =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   347
  ((decimal -- Scan.repeat (Scan.$$ " " |-- Scan.option decimal) >>
33067
a70d5baa53ee use plain Scan.repeat (NB: Scan.bulk is for cascading sources -- mostly interna use);
wenzelm
parents: 33042
diff changeset
   348
    (fn (h,to) => map_filter I ((SOME h)::to))) --| ignore >> vector_of_list) inp
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   349
val parse_csdpoutput = mkparser csdpoutput
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   350
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   351
(* Try some apparently sensible scaling first. Note that this is purely to   *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   352
(* get a cleaner translation to floating-point, and doesn't affect any of    *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   353
(* the results, in principle. In practice it seems a lot better when there   *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   354
(* are extreme numbers in the original problem.                              *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   355
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   356
(* Version for (int*int*int) keys *)
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   357
local
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   358
  fun max_rat x y = if x < y then y else x
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   359
  fun common_denominator fld amat acc =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   360
    fld (fn (_,c) => fn a => lcm_rat (denominator_rat c) a) amat acc
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   361
  fun maximal_element fld amat acc =
44453
082edd97ffe1 comment out dead code to avoid compiler warnings
huffman
parents: 42616
diff changeset
   362
    fld (fn (_,c) => fn maxa => max_rat maxa (abs_rat c)) amat acc
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   363
  fun float_of_rat x =
63201
f151704c08e4 tuned signature;
wenzelm
parents: 63198
diff changeset
   364
    let val (a,b) = Rat.dest x
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   365
    in Real.fromInt a / Real.fromInt b end;
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   366
  fun int_of_float x = (trunc x handle Overflow => 0 | Domain => 0)
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   367
in
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   368
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   369
fun tri_scale_then solver (obj:vector) mats =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   370
  let
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   371
    val cd1 = fold_rev (common_denominator Inttriplefunc.fold) mats @1
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   372
    val cd2 = common_denominator FuncUtil.Intfunc.fold (snd obj) @1
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   373
    val mats' = map (Inttriplefunc.map (fn _ => fn x => cd1 * x)) mats
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   374
    val obj' = vector_cmul cd2 obj
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   375
    val max1 = fold_rev (maximal_element Inttriplefunc.fold) mats' @0
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   376
    val max2 = maximal_element FuncUtil.Intfunc.fold (snd obj') @0
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   377
    val scal1 = pow2 (20 - int_of_float(Math.ln (float_of_rat max1) / Math.ln 2.0))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   378
    val scal2 = pow2 (20 - int_of_float(Math.ln (float_of_rat max2) / Math.ln 2.0))
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   379
    val mats'' = map (Inttriplefunc.map (fn _ => fn x => x * scal1)) mats'
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   380
    val obj'' = vector_cmul scal2 obj'
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   381
  in solver obj'' mats'' end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   382
end;
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   383
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   384
(* Round a vector to "nice" rationals.                                       *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   385
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   386
fun nice_rational n x = round_rat (n * x) / n;
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   387
fun nice_vector n ((d,v) : vector) =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   388
  (d, FuncUtil.Intfunc.fold (fn (i,c) => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   389
      let val y = nice_rational n c in
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   390
        if c = @0 then a
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   391
        else FuncUtil.Intfunc.update (i,y) a
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   392
      end) v FuncUtil.Intfunc.empty): vector
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   393
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   394
fun dest_ord f x = is_equal (f x);
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   395
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   396
(* Stuff for "equations" ((int*int*int)->num functions).                         *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   397
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   398
fun tri_equation_cmul c eq =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   399
  if c = @0 then Inttriplefunc.empty
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   400
  else Inttriplefunc.map (fn _ => fn d => c * d) eq;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   401
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   402
fun tri_equation_add eq1 eq2 =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   403
  Inttriplefunc.combine (curry op +) (fn x => x = @0) eq1 eq2;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   404
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   405
fun tri_equation_eval assig eq =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   406
  let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   407
    fun value v = Inttriplefunc.apply assig v
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   408
  in Inttriplefunc.fold (fn (v, c) => fn a => a + value v * c) eq @0 end;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   409
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   410
(* Eliminate all variables, in an essentially arbitrary order.               *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   411
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   412
fun tri_eliminate_all_equations one =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   413
  let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   414
    fun choose_variable eq =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   415
      let val (v,_) = Inttriplefunc.choose eq
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   416
      in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   417
        if is_equal (triple_int_ord(v,one)) then
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   418
          let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   419
            val eq' = Inttriplefunc.delete_safe v eq
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   420
          in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   421
            if Inttriplefunc.is_empty eq' then error "choose_variable"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   422
            else fst (Inttriplefunc.choose eq')
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   423
          end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   424
        else v
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   425
      end
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   426
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   427
    fun eliminate dun eqs =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   428
      (case eqs of
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   429
        [] => dun
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   430
      | eq :: oeqs =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   431
          if Inttriplefunc.is_empty eq then eliminate dun oeqs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   432
          else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   433
            let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   434
              val v = choose_variable eq
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   435
              val a = Inttriplefunc.apply eq v
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   436
              val eq' =
63201
f151704c08e4 tuned signature;
wenzelm
parents: 63198
diff changeset
   437
                tri_equation_cmul ((Rat.of_int ~1) / a) (Inttriplefunc.delete_safe v eq)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   438
              fun elim e =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   439
                let val b = Inttriplefunc.tryapplyd e v @0 in
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   440
                  if b = @0 then e
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   441
                  else tri_equation_add e (tri_equation_cmul (Rat.neg b / a) eq)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   442
                end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   443
            in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   444
              eliminate (Inttriplefunc.update(v, eq') (Inttriplefunc.map (K elim) dun))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   445
                (map elim oeqs)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   446
            end)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   447
  in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   448
    fn eqs =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   449
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   450
        val assig = eliminate Inttriplefunc.empty eqs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   451
        val vs = Inttriplefunc.fold (fn (_, f) => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   452
          remove (dest_ord triple_int_ord) one (Inttriplefunc.dom f) @ a) assig []
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   453
      in (distinct (dest_ord triple_int_ord) vs,assig) end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   454
  end;
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   455
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   456
(* Multiply equation-parametrized poly by regular poly and add accumulator.  *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   457
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   458
fun tri_epoly_pmul p q acc =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   459
  FuncUtil.Monomialfunc.fold (fn (m1, c) => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   460
    FuncUtil.Monomialfunc.fold (fn (m2, e) => fn b =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   461
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   462
        val m =  monomial_mul m1 m2
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   463
        val es = FuncUtil.Monomialfunc.tryapplyd b m Inttriplefunc.empty
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   464
      in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   465
        FuncUtil.Monomialfunc.update (m,tri_equation_add (tri_equation_cmul c e) es) b
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   466
      end) q a) p acc;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   467
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   468
(* Hence produce the "relevant" monomials: those whose squares lie in the    *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   469
(* Newton polytope of the monomials in the input. (This is enough according  *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   470
(* to Reznik: "Extremal PSD forms with few terms", Duke Math. Journal,       *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   471
(* vol 45, pp. 363--374, 1978.                                               *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   472
(*                                                                           *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   473
(* These are ordered in sort of decreasing degree. In particular the         *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   474
(* constant monomial is last; this gives an order in diagonalization of the  *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   475
(* quadratic form that will tend to display constants.                       *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   476
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   477
(* Diagonalize (Cholesky/LDU) the matrix corresponding to a quadratic form.  *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   478
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   479
local
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   480
  fun diagonalize n i m =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   481
    if FuncUtil.Intpairfunc.is_empty (snd m) then []
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   482
    else
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   483
      let
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   484
        val a11 = FuncUtil.Intpairfunc.tryapplyd (snd m) (i,i) @0
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   485
      in
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   486
        if a11 < @0 then raise Failure "diagonalize: not PSD"
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   487
        else if a11 = @0 then
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   488
          if FuncUtil.Intfunc.is_empty (snd (row i m))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   489
          then diagonalize n (i + 1) m
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   490
          else raise Failure "diagonalize: not PSD ___ "
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   491
        else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   492
          let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   493
            val v = row i m
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   494
            val v' =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   495
              (fst v, FuncUtil.Intfunc.fold (fn (i, c) => fn a =>
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   496
                let val y = c / a11
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   497
                in if y = @0 then a else FuncUtil.Intfunc.update (i,y) a
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   498
                end) (snd v) FuncUtil.Intfunc.empty)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   499
            fun upt0 x y a =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   500
              if y = @0 then a
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   501
              else FuncUtil.Intpairfunc.update (x,y) a
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   502
            val m' =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   503
              ((n, n),
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   504
                iter (i + 1, n) (fn j =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   505
                  iter (i + 1, n) (fn k =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   506
                    (upt0 (j, k)
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   507
                      (FuncUtil.Intpairfunc.tryapplyd (snd m) (j, k) @0 -
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   508
                        FuncUtil.Intfunc.tryapplyd (snd v) j @0 *
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   509
                        FuncUtil.Intfunc.tryapplyd (snd v') k @0))))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   510
                    FuncUtil.Intpairfunc.empty)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   511
          in (a11, v') :: diagonalize n (i + 1) m' end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   512
      end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   513
in
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   514
  fun diag m =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   515
    let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   516
      val nn = dimensions m
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   517
      val n = fst nn
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   518
    in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   519
      if snd nn <> n then error "diagonalize: non-square matrix"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   520
      else diagonalize n 1 m
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   521
    end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   522
end;
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   523
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   524
(* Enumeration of monomials with given multidegree bound.                    *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   525
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   526
fun enumerate_monomials d vars =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   527
  if d < 0 then []
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   528
  else if d = 0 then [FuncUtil.Ctermfunc.empty]
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   529
  else if null vars then [monomial_1]
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   530
  else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   531
    let val alts =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   532
      map_range (fn k =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   533
        let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   534
          val oths = enumerate_monomials (d - k) (tl vars)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   535
        in map (fn ks => if k = 0 then ks else FuncUtil.Ctermfunc.update (hd vars, k) ks) oths end)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   536
        (d + 1)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   537
  in flat alts end;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   538
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   539
(* Enumerate products of distinct input polys with degree <= d.              *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   540
(* We ignore any constant input polynomials.                                 *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   541
(* Give the output polynomial and a record of how it was derived.            *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   542
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   543
fun enumerate_products d pols =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   544
  if d = 0 then [(poly_const @1, RealArith.Rational_lt @1)]
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   545
  else if d < 0 then []
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   546
  else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   547
    (case pols of
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   548
      [] => [(poly_const @1, RealArith.Rational_lt @1)]
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   549
    | (p, b) :: ps =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   550
        let val e = multidegree p in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   551
          if e = 0 then enumerate_products d ps
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   552
          else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   553
            enumerate_products d ps @
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   554
            map (fn (q, c) => (poly_mul p q, RealArith.Product (b, c)))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   555
              (enumerate_products (d - e) ps)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   556
        end)
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   557
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   558
(* Convert regular polynomial. Note that we treat (0,0,0) as -1.             *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   559
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   560
fun epoly_of_poly p =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   561
  FuncUtil.Monomialfunc.fold (fn (m, c) => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   562
      FuncUtil.Monomialfunc.update (m, Inttriplefunc.onefunc ((0, 0, 0), Rat.neg c)) a)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   563
    p FuncUtil.Monomialfunc.empty;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   564
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   565
(* String for block diagonal matrix numbered k.                              *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   566
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   567
fun sdpa_of_blockdiagonal k m =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   568
  let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   569
    val pfx = string_of_int k ^" "
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   570
    val ents =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   571
      Inttriplefunc.fold
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   572
        (fn ((b, i, j), c) => fn a => if i > j then a else ((b, i, j), c) :: a)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   573
        m []
59058
a78612c67ec0 renamed "pairself" to "apply2", in accordance to @{apply 2};
wenzelm
parents: 58635
diff changeset
   574
    val entss = sort (triple_int_ord o apply2 fst) ents
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   575
  in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   576
    fold_rev (fn ((b,i,j),c) => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   577
      pfx ^ string_of_int b ^ " " ^ string_of_int i ^ " " ^ string_of_int j ^
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   578
      " " ^ decimalize 20 c ^ "\n" ^ a) entss ""
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   579
  end;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   580
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   581
(* SDPA for problem using block diagonal (i.e. multiple SDPs)                *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   582
32268
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
   583
fun sdpa_of_blockproblem nblocks blocksizes obj mats =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   584
  let val m = length mats - 1
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   585
  in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   586
    string_of_int m ^ "\n" ^
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   587
    string_of_int nblocks ^ "\n" ^
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   588
    (space_implode " " (map string_of_int blocksizes)) ^
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   589
    "\n" ^
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   590
    sdpa_of_vector obj ^
58635
010b610eb55d eliminated some exotic combinators;
wenzelm
parents: 58631
diff changeset
   591
    fold_rev (fn (k, m) => fn a => sdpa_of_blockdiagonal (k - 1) m ^ a)
010b610eb55d eliminated some exotic combinators;
wenzelm
parents: 58631
diff changeset
   592
      (1 upto length mats ~~ mats) ""
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   593
  end;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   594
32268
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
   595
(* Run prover on a problem in block diagonal form.                       *)
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
   596
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   597
fun run_blockproblem prover nblocks blocksizes obj mats =
32268
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
   598
  parse_csdpoutput (prover (sdpa_of_blockproblem nblocks blocksizes obj mats))
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
   599
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   600
(* 3D versions of matrix operations to consider blocks separately.           *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   601
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   602
val bmatrix_add = Inttriplefunc.combine (curry op +) (fn x => x = @0);
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   603
fun bmatrix_cmul c bm =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   604
  if c = @0 then Inttriplefunc.empty
63198
c583ca33076a ad-hoc overloading for standard operations on type Rat.rat;
wenzelm
parents: 60867
diff changeset
   605
  else Inttriplefunc.map (fn _ => fn x => c * x) bm;
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   606
63201
f151704c08e4 tuned signature;
wenzelm
parents: 63198
diff changeset
   607
val bmatrix_neg = bmatrix_cmul (Rat.of_int ~1);
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   608
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   609
(* Smash a block matrix into components.                                     *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   610
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   611
fun blocks blocksizes bm =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   612
  map (fn (bs, b0) =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   613
    let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   614
      val m =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   615
        Inttriplefunc.fold
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   616
          (fn ((b, i, j), c) => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   617
            if b = b0 then FuncUtil.Intpairfunc.update ((i, j), c) a else a)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   618
        bm FuncUtil.Intpairfunc.empty
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   619
      val _ = FuncUtil.Intpairfunc.fold (fn ((i, j), _) => fn a => max a (max i j)) m 0
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   620
    in (((bs, bs), m): matrix) end)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   621
  (blocksizes ~~ (1 upto length blocksizes));
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   622
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   623
(* FIXME : Get rid of this !!!*)
32268
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
   624
local
44453
082edd97ffe1 comment out dead code to avoid compiler warnings
huffman
parents: 42616
diff changeset
   625
  fun tryfind_with msg _ [] = raise Failure msg
082edd97ffe1 comment out dead code to avoid compiler warnings
huffman
parents: 42616
diff changeset
   626
    | tryfind_with _ f (x::xs) = (f x handle Failure s => tryfind_with s f xs);
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   627
in
32268
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
   628
  fun tryfind f = tryfind_with "tryfind" f
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
   629
end
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
   630
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   631
(* Positiv- and Nullstellensatz. Flag "linf" forces a linear representation. *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   632
38805
b09d8603f865 Sum_Of_Squares: proper configuration options;
wenzelm
parents: 38795
diff changeset
   633
fun real_positivnullstellensatz_general ctxt prover linf d eqs leqs pol =
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   634
  let
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   635
    val vars =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   636
      fold_rev (union (op aconvc) o poly_variables)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   637
        (pol :: eqs @ map fst leqs) []
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   638
    val monoid =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   639
      if linf then
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   640
        (poly_const @1, RealArith.Rational_lt @1)::
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   641
        (filter (fn (p,_) => multidegree p <= d) leqs)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   642
      else enumerate_products d leqs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   643
    val nblocks = length monoid
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   644
    fun mk_idmultiplier k p =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   645
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   646
        val e = d - multidegree p
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   647
        val mons = enumerate_monomials e vars
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   648
        val nons = mons ~~ (1 upto length mons)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   649
      in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   650
        (mons,
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   651
          fold_rev (fn (m, n) =>
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   652
            FuncUtil.Monomialfunc.update (m, Inttriplefunc.onefunc ((~k, ~n, n), @1)))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   653
          nons FuncUtil.Monomialfunc.empty)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   654
      end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   655
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   656
    fun mk_sqmultiplier k (p,_) =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   657
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   658
        val e = (d - multidegree p) div 2
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   659
        val mons = enumerate_monomials e vars
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   660
        val nons = mons ~~ (1 upto length mons)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   661
      in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   662
        (mons,
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   663
          fold_rev (fn (m1, n1) =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   664
            fold_rev (fn (m2, n2) => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   665
              let val m = monomial_mul m1 m2 in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   666
                if n1 > n2 then a
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   667
                else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   668
                  let
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   669
                    val c = if n1 = n2 then @1 else @2
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   670
                    val e = FuncUtil.Monomialfunc.tryapplyd a m Inttriplefunc.empty
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   671
                  in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   672
                    FuncUtil.Monomialfunc.update
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   673
                      (m, tri_equation_add (Inttriplefunc.onefunc ((k, n1, n2), c)) e) a
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   674
                  end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   675
              end) nons) nons FuncUtil.Monomialfunc.empty)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   676
      end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   677
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   678
    val (sqmonlist,sqs) = split_list (map2 mk_sqmultiplier (1 upto length monoid) monoid)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   679
    val (_(*idmonlist*),ids) =  split_list (map2 mk_idmultiplier (1 upto length eqs) eqs)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   680
    val blocksizes = map length sqmonlist
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   681
    val bigsum =
58635
010b610eb55d eliminated some exotic combinators;
wenzelm
parents: 58631
diff changeset
   682
      fold_rev (fn (p, q) => fn a => tri_epoly_pmul p q a) (eqs ~~ ids)
010b610eb55d eliminated some exotic combinators;
wenzelm
parents: 58631
diff changeset
   683
        (fold_rev (fn ((p, _), s) => fn a => tri_epoly_pmul p s a) (monoid ~~ sqs)
010b610eb55d eliminated some exotic combinators;
wenzelm
parents: 58631
diff changeset
   684
          (epoly_of_poly (poly_neg pol)))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   685
    val eqns = FuncUtil.Monomialfunc.fold (fn (_, e) => fn a => e :: a) bigsum []
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   686
    val (pvs, assig) = tri_eliminate_all_equations (0, 0, 0) eqns
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   687
    val qvars = (0, 0, 0) :: pvs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   688
    val allassig =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   689
      fold_rev (fn v => Inttriplefunc.update (v, (Inttriplefunc.onefunc (v, @1)))) pvs assig
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   690
    fun mk_matrix v =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   691
      Inttriplefunc.fold (fn ((b, i, j), ass) => fn m =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   692
          if b < 0 then m
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   693
          else
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   694
            let val c = Inttriplefunc.tryapplyd ass v @0 in
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   695
              if c = @0 then m
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   696
              else Inttriplefunc.update ((b, j, i), c) (Inttriplefunc.update ((b, i, j), c) m)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   697
            end)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   698
        allassig Inttriplefunc.empty
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   699
    val diagents =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   700
      Inttriplefunc.fold
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   701
        (fn ((b, i, j), e) => fn a => if b > 0 andalso i = j then tri_equation_add e a else a)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   702
        allassig Inttriplefunc.empty
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   703
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   704
    val mats = map mk_matrix qvars
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   705
    val obj =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   706
      (length pvs,
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   707
        itern 1 pvs (fn v => fn i =>
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   708
          FuncUtil.Intfunc.updatep iszero (i,Inttriplefunc.tryapplyd diagents v @0))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   709
          FuncUtil.Intfunc.empty)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   710
    val raw_vec =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   711
      if null pvs then vector_0 0
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   712
      else tri_scale_then (run_blockproblem prover nblocks blocksizes) obj mats
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   713
    fun int_element (_, v) i = FuncUtil.Intfunc.tryapplyd v i @0
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   714
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   715
    fun find_rounding d =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   716
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   717
        val _ =
58631
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
   718
          debug_message ctxt (fn () => "Trying rounding with limit "^Rat.string_of_rat d ^ "\n")
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   719
        val vec = nice_vector d raw_vec
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   720
        val blockmat =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   721
          iter (1, dim vec)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   722
            (fn i => fn a => bmatrix_add (bmatrix_cmul (int_element vec i) (nth mats i)) a)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   723
            (bmatrix_neg (nth mats 0))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   724
        val allmats = blocks blocksizes blockmat
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   725
      in (vec, map diag allmats) end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   726
    val (vec, ratdias) =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   727
      if null pvs then find_rounding @1
63201
f151704c08e4 tuned signature;
wenzelm
parents: 63198
diff changeset
   728
      else tryfind find_rounding (map Rat.of_int (1 upto 31) @ map pow2 (5 upto 66))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   729
    val newassigs =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   730
      fold_rev (fn k => Inttriplefunc.update (nth pvs (k - 1), int_element vec k))
63201
f151704c08e4 tuned signature;
wenzelm
parents: 63198
diff changeset
   731
        (1 upto dim vec) (Inttriplefunc.onefunc ((0, 0, 0), Rat.of_int ~1))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   732
    val finalassigs =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   733
      Inttriplefunc.fold (fn (v, e) => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   734
        Inttriplefunc.update (v, tri_equation_eval newassigs e) a) allassig newassigs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   735
    fun poly_of_epoly p =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   736
      FuncUtil.Monomialfunc.fold (fn (v, e) => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   737
          FuncUtil.Monomialfunc.updatep iszero (v, tri_equation_eval finalassigs e) a)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   738
        p FuncUtil.Monomialfunc.empty
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   739
    fun mk_sos mons =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   740
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   741
        fun mk_sq (c, m) =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   742
          (c, fold_rev (fn k => fn a =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   743
              FuncUtil.Monomialfunc.updatep iszero (nth mons (k - 1), int_element m k) a)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   744
            (1 upto length mons) FuncUtil.Monomialfunc.empty)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   745
      in map mk_sq end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   746
    val sqs = map2 mk_sos sqmonlist ratdias
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   747
    val cfs = map poly_of_epoly ids
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   748
    val msq = filter (fn (_, b) => not (null b)) (map2 pair monoid sqs)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   749
    fun eval_sq sqs = fold_rev (fn (c, q) => poly_add (poly_cmul c (poly_mul q q))) sqs poly_0
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   750
    val sanity =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   751
      fold_rev (fn ((p, _), s) => poly_add (poly_mul p (eval_sq s))) msq
58635
010b610eb55d eliminated some exotic combinators;
wenzelm
parents: 58631
diff changeset
   752
        (fold_rev (fn (p, q) => poly_add (poly_mul p q)) (cfs ~~ eqs) (poly_neg pol))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   753
  in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   754
    if not(FuncUtil.Monomialfunc.is_empty sanity) then raise Sanity
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   755
    else (cfs, map (fn (a, b) => (snd a, b)) msq)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   756
  end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   757
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   758
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   759
(* Iterative deepening.                                                      *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   760
58631
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
   761
fun deepen ctxt f n =
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
   762
  (trace_message ctxt (fn () => "Searching with depth limit " ^ string_of_int n);
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
   763
    (f n handle Failure s =>
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
   764
      (trace_message ctxt (fn () => "failed with message: " ^ s); deepen ctxt f (n + 1))));
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   765
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   766
32645
1cc5b24f5a01 sos method generates and uses proof certificates
Philipp Meyer
parents: 32332
diff changeset
   767
(* Map back polynomials and their composites to a positivstellensatz.        *)
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   768
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   769
fun cterm_of_sqterm (c, p) = RealArith.Product (RealArith.Rational_lt c, RealArith.Square p);
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   770
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   771
fun cterm_of_sos (pr,sqs) =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   772
  if null sqs then pr
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   773
  else RealArith.Product (pr, foldr1 RealArith.Sum (map cterm_of_sqterm sqs));
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   774
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   775
(* Interface to HOL.                                                         *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   776
local
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   777
  open Conv
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   778
  val concl = Thm.dest_arg o Thm.cprop_of
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   779
  fun simple_cterm_ord t u = Term_Ord.fast_term_ord (Thm.term_of t, Thm.term_of u) = LESS
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   780
in
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   781
(* FIXME: Replace tryfind by get_first !! *)
32645
1cc5b24f5a01 sos method generates and uses proof certificates
Philipp Meyer
parents: 32332
diff changeset
   782
fun real_nonlinear_prover proof_method ctxt =
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   783
  let
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   784
    val {add = _, mul = _, neg = _, pow = _, sub = _, main = real_poly_conv} =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   785
      Semiring_Normalizer.semiring_normalizers_ord_wrapper ctxt
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   786
        (the (Semiring_Normalizer.match ctxt @{cterm "(0::real) + 1"}))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   787
        simple_cterm_ord
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   788
    fun mainf cert_choice translator (eqs, les, lts) =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   789
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   790
        val eq0 = map (poly_of_term o Thm.dest_arg1 o concl) eqs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   791
        val le0 = map (poly_of_term o Thm.dest_arg o concl) les
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   792
        val lt0 = map (poly_of_term o Thm.dest_arg o concl) lts
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   793
        val eqp0 = map_index (fn (i, t) => (t,RealArith.Axiom_eq i)) eq0
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   794
        val lep0 = map_index (fn (i, t) => (t,RealArith.Axiom_le i)) le0
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   795
        val ltp0 = map_index (fn (i, t) => (t,RealArith.Axiom_lt i)) lt0
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   796
        val (keq,eq) = List.partition (fn (p, _) => multidegree p = 0) eqp0
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   797
        val (klep,lep) = List.partition (fn (p, _) => multidegree p = 0) lep0
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   798
        val (kltp,ltp) = List.partition (fn (p, _) => multidegree p = 0) ltp0
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   799
        fun trivial_axiom (p, ax) =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   800
          (case ax of
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   801
            RealArith.Axiom_eq n =>
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   802
              if eval FuncUtil.Ctermfunc.empty p <> @0 then nth eqs n
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   803
              else raise Failure "trivial_axiom: Not a trivial axiom"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   804
          | RealArith.Axiom_le n =>
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   805
              if eval FuncUtil.Ctermfunc.empty p < @0 then nth les n
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   806
              else raise Failure "trivial_axiom: Not a trivial axiom"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   807
          | RealArith.Axiom_lt n =>
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   808
              if eval FuncUtil.Ctermfunc.empty p <= @0 then nth lts n
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   809
              else raise Failure "trivial_axiom: Not a trivial axiom"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   810
          | _ => error "trivial_axiom: Not a trivial axiom")
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   811
      in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   812
        let val th = tryfind trivial_axiom (keq @ klep @ kltp) in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   813
          (fconv_rule (arg_conv (arg1_conv (real_poly_conv ctxt))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   814
            then_conv Numeral_Simprocs.field_comp_conv ctxt) th,
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   815
            RealArith.Trivial)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   816
        end handle Failure _ =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   817
          let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   818
            val proof =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   819
              (case proof_method of
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   820
                Certificate certs =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   821
                  (* choose certificate *)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   822
                  let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   823
                    fun chose_cert [] (RealArith.Cert c) = c
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   824
                      | chose_cert (RealArith.Left::s) (RealArith.Branch (l, _)) = chose_cert s l
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   825
                      | chose_cert (RealArith.Right::s) (RealArith.Branch (_, r)) = chose_cert s r
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   826
                      | chose_cert _ _ = error "certificate tree in invalid form"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   827
                  in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   828
                    chose_cert cert_choice certs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   829
                  end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   830
              | Prover prover =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   831
                  (* call prover *)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   832
                  let
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   833
                    val pol = fold_rev poly_mul (map fst ltp) (poly_const @1)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   834
                    val leq = lep @ ltp
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   835
                    fun tryall d =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   836
                      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   837
                        val e = multidegree pol
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   838
                        val k = if e = 0 then 0 else d div e
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   839
                        val eq' = map fst eq
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   840
                      in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   841
                        tryfind (fn i =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   842
                            (d, i, real_positivnullstellensatz_general ctxt prover false d eq' leq
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   843
                              (poly_neg(poly_pow pol i))))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   844
                          (0 upto k)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   845
                      end
58631
41333b45bff9 clarified messages depending on options;
wenzelm
parents: 57889
diff changeset
   846
                    val (_,i,(cert_ideal,cert_cone)) = deepen ctxt tryall 0
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   847
                    val proofs_ideal =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   848
                      map2 (fn q => fn (_,ax) => RealArith.Eqmul(q,ax)) cert_ideal eq
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   849
                    val proofs_cone = map cterm_of_sos cert_cone
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   850
                    val proof_ne =
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   851
                      if null ltp then RealArith.Rational_lt @1
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   852
                      else
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   853
                        let val p = foldr1 RealArith.Product (map snd ltp) in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   854
                          funpow i (fn q => RealArith.Product (p, q))
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   855
                            (RealArith.Rational_lt @1)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   856
                        end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   857
                  in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   858
                    foldr1 RealArith.Sum (proof_ne :: proofs_ideal @ proofs_cone)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   859
                  end)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   860
          in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   861
            (translator (eqs,les,lts) proof, RealArith.Cert proof)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   862
          end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   863
      end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   864
  in mainf end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   865
end
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   866
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   867
(* FIXME : This is very bad!!!*)
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   868
fun subst_conv eqs t =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   869
  let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   870
    val t' = fold (Thm.lambda o Thm.lhs_of) eqs t
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   871
  in
60329
e85ed7a36b2f eliminated odd C combinator -- Isabelle/ML usually has canonical argument order;
wenzelm
parents: 59582
diff changeset
   872
    Conv.fconv_rule (Thm.beta_conversion true)
e85ed7a36b2f eliminated odd C combinator -- Isabelle/ML usually has canonical argument order;
wenzelm
parents: 59582
diff changeset
   873
      (fold (fn a => fn b => Thm.combination b a) eqs (Thm.reflexive t'))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   874
  end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   875
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   876
(* A wrapper that tries to substitute away variables first.                  *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   877
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   878
local
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   879
  open Conv
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   880
  fun simple_cterm_ord t u = Term_Ord.fast_term_ord (Thm.term_of t, Thm.term_of u) = LESS
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   881
  val concl = Thm.dest_arg o Thm.cprop_of
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   882
  val shuffle1 =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   883
    fconv_rule (rewr_conv @{lemma "(a + x == y) == (x == y - (a::real))"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   884
      by (atomize (full)) (simp add: field_simps)})
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   885
  val shuffle2 =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   886
    fconv_rule (rewr_conv @{lemma "(x + a == y) ==  (x == y - (a::real))"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   887
      by (atomize (full)) (simp add: field_simps)})
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   888
  fun substitutable_monomial fvs tm =
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   889
    (case Thm.term_of tm of
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   890
      Free (_, @{typ real}) =>
63205
97b721666890 prefer rat numberals;
wenzelm
parents: 63201
diff changeset
   891
        if not (member (op aconvc) fvs tm) then (@1, tm)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   892
        else raise Failure "substitutable_monomial"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   893
    | @{term "op * :: real => _"} $ _ $ (Free _) =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   894
        if RealArith.is_ratconst (Thm.dest_arg1 tm) andalso
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   895
          not (member (op aconvc) fvs (Thm.dest_arg tm))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   896
        then (RealArith.dest_ratconst (Thm.dest_arg1 tm), Thm.dest_arg tm)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   897
        else raise Failure "substitutable_monomial"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   898
    | @{term "op + :: real => _"}$_$_ =>
60818
5e11a6e2b044 more direct access to atomic cterms;
wenzelm
parents: 60801
diff changeset
   899
         (substitutable_monomial (Drule.cterm_add_frees (Thm.dest_arg tm) fvs) (Thm.dest_arg1 tm)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   900
           handle Failure _ =>
60818
5e11a6e2b044 more direct access to atomic cterms;
wenzelm
parents: 60801
diff changeset
   901
            substitutable_monomial (Drule.cterm_add_frees (Thm.dest_arg1 tm) fvs) (Thm.dest_arg tm))
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   902
    | _ => raise Failure "substitutable_monomial")
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   903
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   904
  fun isolate_variable v th =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   905
    let
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   906
      val w = Thm.dest_arg1 (Thm.cprop_of th)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   907
    in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   908
      if v aconvc w then th
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   909
      else
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   910
        (case Thm.term_of w of
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   911
          @{term "op + :: real => _"} $ _ $ _ =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   912
            if Thm.dest_arg1 w aconvc v then shuffle2 th
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   913
            else isolate_variable v (shuffle1 th)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   914
        | _ => error "isolate variable : This should not happen?")
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   915
   end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   916
in
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   917
32268
d50f0cb67578 Functionality for sum of squares to call a remote csdp prover
Philipp Meyer
parents: 32150
diff changeset
   918
fun real_nonlinear_subst_prover prover ctxt =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   919
  let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   920
    val {add = _, mul = real_poly_mul_conv, neg = _, pow = _, sub = _, main = real_poly_conv} =
44453
082edd97ffe1 comment out dead code to avoid compiler warnings
huffman
parents: 42616
diff changeset
   921
      Semiring_Normalizer.semiring_normalizers_ord_wrapper ctxt
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   922
        (the (Semiring_Normalizer.match ctxt @{cterm "(0::real) + 1"}))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   923
        simple_cterm_ord
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   924
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   925
    fun make_substitution th =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   926
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   927
        val (c,v) = substitutable_monomial [] (Thm.dest_arg1(concl th))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   928
        val th1 =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   929
          Drule.arg_cong_rule
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   930
            (Thm.apply @{cterm "op * :: real => _"} (RealArith.cterm_of_rat (Rat.inv c)))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   931
            (mk_meta_eq th)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   932
        val th2 = fconv_rule (binop_conv (real_poly_mul_conv ctxt)) th1
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   933
      in fconv_rule (arg_conv (real_poly_conv ctxt)) (isolate_variable v th2) end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   934
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   935
    fun oprconv cv ct =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   936
      let val g = Thm.dest_fun2 ct in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   937
        if g aconvc @{cterm "op <= :: real => _"} orelse g aconvc @{cterm "op < :: real => _"}
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   938
        then arg_conv cv ct else arg1_conv cv ct
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   939
      end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   940
    fun mainf cert_choice translator =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   941
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   942
        fun substfirst (eqs, les, lts) =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   943
          (let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   944
              val eth = tryfind make_substitution eqs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   945
              val modify =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   946
                fconv_rule (arg_conv (oprconv(subst_conv [eth] then_conv (real_poly_conv ctxt))))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   947
            in
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   948
              substfirst
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   949
                (filter_out
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   950
                  (fn t => (Thm.dest_arg1 o Thm.dest_arg o Thm.cprop_of) t aconvc @{cterm "0::real"})
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   951
                  (map modify eqs),
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   952
                  map modify les,
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   953
                  map modify lts)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   954
            end handle Failure  _ =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   955
              real_nonlinear_prover prover ctxt cert_choice translator (rev eqs, rev les, rev lts))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   956
      in substfirst end
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   957
  in mainf end
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   958
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   959
(* Overall function. *)
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   960
32645
1cc5b24f5a01 sos method generates and uses proof certificates
Philipp Meyer
parents: 32332
diff changeset
   961
fun real_sos prover ctxt =
32828
ad76967c703d removed opening of structures
Philipp Meyer
parents: 32740
diff changeset
   962
  RealArith.gen_prover_real_arith ctxt (real_nonlinear_subst_prover prover ctxt)
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   963
31119
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   964
end;
2532bb2d65c7 A decision method for universal multivariate real arithmetic with add
chaieb
parents:
diff changeset
   965
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   966
val known_sos_constants =
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   967
  [@{term "op ==>"}, @{term "Trueprop"},
56625
54505623a8ef sos accepts False, returns apply command
paulson <lp15@cam.ac.uk>
parents: 56536
diff changeset
   968
   @{term HOL.False}, @{term HOL.implies}, @{term HOL.conj}, @{term HOL.disj},
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   969
   @{term "Not"}, @{term "op = :: bool => _"},
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   970
   @{term "All :: (real => _) => _"}, @{term "Ex :: (real => _) => _"},
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   971
   @{term "op = :: real => _"}, @{term "op < :: real => _"},
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   972
   @{term "op <= :: real => _"},
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   973
   @{term "op + :: real => _"}, @{term "op - :: real => _"},
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   974
   @{term "op * :: real => _"}, @{term "uminus :: real => _"},
31512
27118561c2e0 Tuned sos tactic to reject non SOS goals
chaieb
parents: 31131
diff changeset
   975
   @{term "op / :: real => _"}, @{term "inverse :: real => _"},
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   976
   @{term "op ^ :: real => _"}, @{term "abs :: real => _"},
31512
27118561c2e0 Tuned sos tactic to reject non SOS goals
chaieb
parents: 31131
diff changeset
   977
   @{term "min :: real => _"}, @{term "max :: real => _"},
47108
2a1953f0d20d merged fork with new numeral representation (see NEWS)
huffman
parents: 46497
diff changeset
   978
   @{term "0::real"}, @{term "1::real"},
2a1953f0d20d merged fork with new numeral representation (see NEWS)
huffman
parents: 46497
diff changeset
   979
   @{term "numeral :: num => nat"},
2a1953f0d20d merged fork with new numeral representation (see NEWS)
huffman
parents: 46497
diff changeset
   980
   @{term "numeral :: num => real"},
2a1953f0d20d merged fork with new numeral representation (see NEWS)
huffman
parents: 46497
diff changeset
   981
   @{term "Num.Bit0"}, @{term "Num.Bit1"}, @{term "Num.One"}];
31512
27118561c2e0 Tuned sos tactic to reject non SOS goals
chaieb
parents: 31131
diff changeset
   982
32839
a007a7cd8c2f tuned header;
wenzelm
parents: 32838
diff changeset
   983
fun check_sos kcts ct =
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   984
  let
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
   985
    val t = Thm.term_of ct
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   986
    val _ =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   987
      if not (null (Term.add_tfrees t []) andalso null (Term.add_tvars t []))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   988
      then error "SOS: not sos. Additional type varables"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   989
      else ()
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   990
    val fs = Term.add_frees t []
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   991
    val _ =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   992
      if exists (fn ((_,T)) => not (T = @{typ "real"})) fs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   993
      then error "SOS: not sos. Variables with type not real"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   994
      else ()
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   995
    val vs = Term.add_vars t []
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   996
    val _ =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   997
      if exists (fn ((_,T)) => not (T = @{typ "real"})) vs
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   998
      then error "SOS: not sos. Variables with type not real"
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
   999
      else ()
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1000
    val ukcs = subtract (fn (t,p) => Const p aconv t) kcts (Term.add_consts t [])
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1001
    val _ =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1002
      if null ukcs then ()
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1003
      else error ("SOSO: Unknown constants in Subgoal:" ^ commas (map fst ukcs))
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1004
  in () end
31512
27118561c2e0 Tuned sos tactic to reject non SOS goals
chaieb
parents: 31131
diff changeset
  1005
60752
b48830b670a1 prefer tactics with explicit context;
wenzelm
parents: 60329
diff changeset
  1006
fun core_sos_tac print_cert prover = SUBPROOF (fn {concl, context = ctxt, ...} =>
32831
1352736e5727 changed core_sos_tac to use SUBPROOF
Philipp Meyer
parents: 32830
diff changeset
  1007
  let
32838
d9dfd30af9ae core_sos_tac: SUBPROOF body operates on subgoal 1;
wenzelm
parents: 32831
diff changeset
  1008
    val _ = check_sos known_sos_constants concl
60752
b48830b670a1 prefer tactics with explicit context;
wenzelm
parents: 60329
diff changeset
  1009
    val (th, certificates) = real_sos prover ctxt (Thm.dest_arg concl)
32838
d9dfd30af9ae core_sos_tac: SUBPROOF body operates on subgoal 1;
wenzelm
parents: 32831
diff changeset
  1010
    val _ = print_cert certificates
60752
b48830b670a1 prefer tactics with explicit context;
wenzelm
parents: 60329
diff changeset
  1011
  in resolve_tac ctxt [th] 1 end);
31131
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1012
44453
082edd97ffe1 comment out dead code to avoid compiler warnings
huffman
parents: 42616
diff changeset
  1013
fun default_SOME _ NONE v = SOME v
082edd97ffe1 comment out dead code to avoid compiler warnings
huffman
parents: 42616
diff changeset
  1014
  | default_SOME _ (SOME v) _ = SOME v;
31131
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1015
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1016
fun lift_SOME f NONE a = f a
44453
082edd97ffe1 comment out dead code to avoid compiler warnings
huffman
parents: 42616
diff changeset
  1017
  | lift_SOME _ (SOME a) _ = SOME a;
31131
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1018
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1019
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1020
local
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
  1021
  val is_numeral = can (HOLogic.dest_number o Thm.term_of)
31131
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1022
in
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1023
  fun get_denom b ct =
59582
0fbed69ff081 tuned signature -- prefer qualified names;
wenzelm
parents: 59058
diff changeset
  1024
    (case Thm.term_of ct of
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1025
      @{term "op / :: real => _"} $ _ $ _ =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1026
        if is_numeral (Thm.dest_arg ct)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1027
        then get_denom b (Thm.dest_arg1 ct)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1028
        else default_SOME (get_denom b) (get_denom b (Thm.dest_arg ct)) (Thm.dest_arg ct, b)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1029
    | @{term "op < :: real => _"} $ _ $ _ =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1030
        lift_SOME (get_denom true) (get_denom true (Thm.dest_arg ct)) (Thm.dest_arg1 ct)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1031
    | @{term "op <= :: real => _"} $ _ $ _ =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1032
        lift_SOME (get_denom true) (get_denom true (Thm.dest_arg ct)) (Thm.dest_arg1 ct)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1033
    | _ $ _ => lift_SOME (get_denom b) (get_denom b (Thm.dest_fun ct)) (Thm.dest_arg ct)
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1034
    | _ => NONE)
31131
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1035
end;
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1036
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1037
fun elim_one_denom_tac ctxt = CSUBGOAL (fn (P, i) =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1038
  (case get_denom false P of
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1039
    NONE => no_tac
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1040
  | SOME (d, ord) =>
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1041
      let
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1042
        val simp_ctxt =
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1043
          ctxt addsimps @{thms field_simps}
60867
86e7560e07d0 slight cleanup of lemmas
haftmann
parents: 60818
diff changeset
  1044
          addsimps [@{thm power_divide}]
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1045
        val th =
60801
7664e0916eec tuned signature;
wenzelm
parents: 60752
diff changeset
  1046
          Thm.instantiate' [] [SOME d, SOME (Thm.dest_arg P)]
55508
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1047
            (if ord then @{lemma "(d=0 --> P) & (d>0 --> P) & (d<(0::real) --> P) ==> P" by auto}
90c42b130652 tuned whitespace;
wenzelm
parents: 54742
diff changeset
  1048
             else @{lemma "(d=0 --> P) & (d ~= (0::real) --> P) ==> P" by blast})
60752
b48830b670a1 prefer tactics with explicit context;
wenzelm
parents: 60329
diff changeset
  1049
      in resolve_tac ctxt [th] i THEN Simplifier.asm_full_simp_tac simp_ctxt i end));
31131
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1050
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1051
fun elim_denom_tac ctxt i = REPEAT (elim_one_denom_tac ctxt i);
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1052
32949
aa6c470a962a eliminated slightly odd get/set operations in favour of Unsynchronized.ref;
wenzelm
parents: 32839
diff changeset
  1053
fun sos_tac print_cert prover ctxt =
56536
aefb4a8da31f made mult_nonneg_nonneg a simp rule
nipkow
parents: 55508
diff changeset
  1054
  (* The SOS prover breaks if mult_nonneg_nonneg is in the simpset *)
57889
049e13f616d4 reduced warnings;
wenzelm
parents: 56625
diff changeset
  1055
  let val ctxt' = Context_Position.set_visible false ctxt delsimps @{thms mult_nonneg_nonneg}
56536
aefb4a8da31f made mult_nonneg_nonneg a simp rule
nipkow
parents: 55508
diff changeset
  1056
  in Object_Logic.full_atomize_tac ctxt' THEN'
aefb4a8da31f made mult_nonneg_nonneg a simp rule
nipkow
parents: 55508
diff changeset
  1057
     elim_denom_tac ctxt' THEN'
aefb4a8da31f made mult_nonneg_nonneg a simp rule
nipkow
parents: 55508
diff changeset
  1058
     core_sos_tac print_cert prover ctxt'
aefb4a8da31f made mult_nonneg_nonneg a simp rule
nipkow
parents: 55508
diff changeset
  1059
  end;
31131
d9752181691a Now deals with division
chaieb
parents: 31119
diff changeset
  1060
31512
27118561c2e0 Tuned sos tactic to reject non SOS goals
chaieb
parents: 31131
diff changeset
  1061
end;