(Joint Center)Library MathComp.closed_field

(* (c) Copyright Microsoft Corporation and Inria.                       
 You may distribute this file under the terms of the CeCILL-B license *)

Require Import ssreflect ssrfun ssrbool eqtype ssrnat seq.
Require Import bigop ssralg poly polydiv.

A proof that algebraically closed field enjoy quantifier elimination, as described in ``A formal quantifier elimination for algebraically closed fields'', proceedings of Calculemus 2010, by Cyril Cohen and Assia Mahboubi
This file constructs an instance of quantifier elimination mixin, (see the ssralg library) from the theory of polynomials with coefficients is an algebraically closed field (see the polydiv library).
This file hence deals with the transformation of formulae part, which we address by implementing one CPS style formula transformer per effective operation involved in the proof of quantifier elimination. See the paper for more details.

Set Implicit Arguments.

Import GRing.
Local Open Scope ring_scope.

Import Pdiv.Ring.
Import PreClosedField.

Section ClosedFieldQE.

Variable F : Field.type.

Variable axiom : ClosedField.axiom F.

Notation fF := (formula F).
Notation qf f := (qf_form f && rformula f).

Definition polyF := seq (term F).

Fixpoint eval_poly (e : seq F) pf :=
  if pf is c::qf then (eval_poly e qf)*'X + (eval e c)%:P else 0.

Definition rpoly (p : polyF) := all (@rterm F) p.

Fixpoint sizeT (k : natfF) (p : polyF) :=
  if p is c::q then
    sizeT (fun n
      if n is m.+1 then k m.+2 else
        GRing.If (c == 0) (k 0%N) (k 1%N)) q
   else k O%N.

Lemma sizeTP (k : natformula F) (pf : polyF) (e : seq F) :
  qf_eval e (sizeT k pf) = qf_eval e (k (size (eval_poly e pf))).

Lemma sizeT_qf (k : natformula F) (p : polyF) :
  ( n, qf (k n)) → rpoly pqf (sizeT k p).

Definition isnull (k : boolfF) (p : polyF) :=
  sizeT (fun nk (n == 0%N)) p.

Lemma isnullP (k : boolformula F) (p : polyF) (e : seq F) :
  qf_eval e (isnull k p) = qf_eval e (k (eval_poly e p == 0)).

Lemma isnull_qf (k : boolformula F) (p : polyF) :
  ( b, qf (k b)) → rpoly pqf (isnull k p).

Definition lt_sizeT (k : boolfF) (p q : polyF) : fF :=
  sizeT (fun nsizeT (fun mk (n<m)) q) p.

Definition lift (p : {poly F}) := let: q := p in map Const q.

Lemma eval_lift (e : seq F) (p : {poly F}) : eval_poly e (lift p) = p.

Fixpoint lead_coefT (k : term FfF) p :=
  if p is c::q then
    lead_coefT (fun lGRing.If (l == 0) (k c) (k l)) q
  else k (Const 0).

Lemma lead_coefTP (k : term Fformula F) :
 ( x e, qf_eval e (k x) = qf_eval e (k (Const (eval e x)))) →
   (p : polyF) (e : seq F),
  qf_eval e (lead_coefT k p) = qf_eval e (k (Const (lead_coef (eval_poly e p)))).

Lemma lead_coefT_qf (k : term Fformula F) (p : polyF) :
 ( c, rterm cqf (k c)) → rpoly pqf (lead_coefT k p).

Fixpoint amulXnT (a : term F) (n : nat) : polyF :=
  if n is n'.+1 then (Const 0) :: (amulXnT a n') else [::a].

Lemma eval_amulXnT (a : term F) (n : nat) (e : seq F) :
  eval_poly e (amulXnT a n) = (eval e a)%:P × 'X^n.

Lemma ramulXnT: a n, rterm arpoly (amulXnT a n).

Fixpoint sumpT (p q : polyF) :=
  if p is a::p' then
    if q is b::q' then (Add a b)::(sumpT p' q')
      else p
    else q.

Lemma eval_sumpT (p q : polyF) (e : seq F) :
  eval_poly e (sumpT p q) = (eval_poly e p) + (eval_poly e q).

Lemma rsumpT (p q : polyF) : rpoly prpoly qrpoly (sumpT p q).

Fixpoint mulpT (p q : polyF) :=
  if p is a :: p' then sumpT (map (Mul a) q) (Const 0::(mulpT p' q)) else [::].

Lemma eval_mulpT (p q : polyF) (e : seq F) :
  eval_poly e (mulpT p q) = (eval_poly e p) × (eval_poly e q).

Lemma rpoly_map_mul (t : term F) (p : polyF) (rt : rterm t) :
  rpoly (map (Mul t) p) = rpoly p.

Lemma rmulpT (p q : polyF) : rpoly prpoly qrpoly (mulpT p q).

Definition opppT := map (Mul (@Const F (-1))).

Lemma eval_opppT (p : polyF) (e : seq F) :
 eval_poly e (opppT p) = - eval_poly e p.

Definition natmulpT n := map (Mul (@NatConst F n)).

Lemma eval_natmulpT (p : polyF) (n : nat) (e : seq F) :
  eval_poly e (natmulpT n p) = (eval_poly e p) *+ n.

Fixpoint redivp_rec_loopT (q : polyF) sq cq (k : nat × polyF × polyFfF)
  (c : nat) (qq r : polyF) (n : nat) {struct n}:=
  sizeT (fun sr
    if sr < sq then k (c, qq, r) else
      lead_coefT (fun lr
        let m := amulXnT lr (sr - sq) in
        let qq1 := sumpT (mulpT qq [::cq]) m in
        let r1 := sumpT (mulpT r ([::cq])) (opppT (mulpT m q)) in
        if n is n1.+1 then redivp_rec_loopT q sq cq k c.+1 qq1 r1 n1
        else k (c.+1, qq1, r1)
      ) r
  ) r.

Fixpoint redivp_rec_loop (q : {poly F}) sq cq
   (k : nat) (qq r : {poly F})(n : nat) {struct n} :=
    if size r < sq then (k, qq, r) else
      let m := (lead_coef r) *: 'X^(size r - sq) in
      let qq1 := qq × cq%:P + m in
      let r1 := r × cq%:P - m × q in
      if n is n1.+1 then redivp_rec_loop q sq cq k.+1 qq1 r1 n1 else
        (k.+1, qq1, r1).

Lemma redivp_rec_loopTP (k : nat × polyF × polyFformula F) :
  ( c qq r e, qf_eval e (k (c,qq,r))
    = qf_eval e (k (c, lift (eval_poly e qq), lift (eval_poly e r))))
  → q sq cq c qq r n e
    (d := redivp_rec_loop (eval_poly e q) sq (eval e cq)
      c (eval_poly e qq) (eval_poly e r) n),
    qf_eval e (redivp_rec_loopT q sq cq k c qq r n)
    = qf_eval e (k (d.1.1, lift d.1.2, lift d.2)).

Lemma redivp_rec_loopT_qf (q : polyF) (sq : nat) (cq : term F)
 (k : nat × polyF × polyFformula F) (c : nat) (qq r : polyF) (n : nat) :
  ( r, [&& rpoly r.1.2 & rpoly r.2]qf (k r)) →
  rpoly qrterm cqrpoly qqrpoly r
    qf (redivp_rec_loopT q sq cq k c qq r n).

Definition redivpT (p : polyF) (k : nat × polyF × polyFfF)
                   (q : polyF) : fF :=
  isnull (fun b
    if b then k (0%N, [::Const 0], p) else
      sizeT (fun sq
        sizeT (fun sp
          lead_coefT (fun lq
            redivp_rec_loopT q sq lq k 0 [::Const 0] p sp
          ) q
        ) p
      ) q
  ) q.

Lemma redivp_rec_loopP (q : {poly F}) (c : nat) (qq r : {poly F}) (n : nat) :
  redivp_rec q c qq r n = redivp_rec_loop q (size q) (lead_coef q) c qq r n.

Lemma redivpTP (k : nat × polyF × polyFformula F) :
  ( c qq r e,
     qf_eval e (k (c,qq,r)) =
     qf_eval e (k (c, lift (eval_poly e qq), lift (eval_poly e r)))) →
   p q e (d := redivp (eval_poly e p) (eval_poly e q)),
    qf_eval e (redivpT p k q) = qf_eval e (k (d.1.1, lift d.1.2, lift d.2)).

Lemma redivpT_qf (p : polyF) (k : nat × polyF × polyFformula F) (q : polyF) :
  ( r, [&& rpoly r.1.2 & rpoly r.2]qf (k r)) →
   rpoly prpoly qqf (redivpT p k q).

Definition rmodpT (p : polyF) (k : polyFfF) (q : polyF) : fF :=
  redivpT p (fun dk d.2) q.
Definition rdivpT (p : polyF) (k:polyFfF) (q : polyF) : fF :=
  redivpT p (fun dk d.1.2) q.
Definition rscalpT (p : polyF) (k: natfF) (q : polyF) : fF :=
  redivpT p (fun dk d.1.1) q.
Definition rdvdpT (p : polyF) (k:boolfF) (q : polyF) : fF :=
  rmodpT p (isnull k) q.

Fixpoint rgcdp_loop n (pp qq : {poly F}) {struct n} :=
  if rmodp pp qq == 0 then qq
    else if n is n1.+1 then rgcdp_loop n1 qq (rmodp pp qq)
         else rmodp pp qq.

Fixpoint rgcdp_loopT (pp : polyF) (k : polyFformula F) n (qq : polyF) :=
  rmodpT pp (isnull
    (fun bif b then (k qq)
              else (if n is n1.+1
                    then rmodpT pp (rgcdp_loopT qq k n1) qq
                    else rmodpT pp k qq)
    )
            ) qq.

Lemma rgcdp_loopP (k : polyFformula F) :
  ( p e, qf_eval e (k p) = qf_eval e (k (lift (eval_poly e p)))) →
   n p q e,
    qf_eval e (rgcdp_loopT p k n q) =
    qf_eval e (k (lift (rgcdp_loop n (eval_poly e p) (eval_poly e q)))).

Lemma rgcdp_loopT_qf (p : polyF) (k : polyFformula F) (q : polyF) (n : nat) :
  ( r, rpoly rqf (k r)) →
  rpoly prpoly qqf (rgcdp_loopT p k n q).

Definition rgcdpT (p : polyF) k (q : polyF) : fF :=
  let aux p1 k q1 := isnull
    (fun bif b
      then (k q1)
      else (sizeT (fun n ⇒ (rgcdp_loopT p1 k n q1)) p1)) p1
    in (lt_sizeT (fun bif b then (aux q k p) else (aux p k q)) p q).

Lemma rgcdpTP (k : seq (term F) → formula F) :
  ( p e, qf_eval e (k p) = qf_eval e (k (lift (eval_poly e p)))) →
    p q e, qf_eval e (rgcdpT p k q) =
                 qf_eval e (k (lift (rgcdp (eval_poly e p) (eval_poly e q)))).

Lemma rgcdpT_qf (p : polyF) (k : polyFformula F) (q : polyF) :
  ( r, rpoly rqf (k r)) → rpoly prpoly qqf (rgcdpT p k q).

Fixpoint rgcdpTs k (ps : seq polyF) : fF :=
  if ps is p::pr then rgcdpTs (rgcdpT p k) pr else k [::Const 0].

Lemma rgcdpTsP (k : polyFformula F) :
  ( p e, qf_eval e (k p) = qf_eval e (k (lift (eval_poly e p)))) →
   ps e,
    qf_eval e (rgcdpTs k ps) =
    qf_eval e (k (lift (\big[@rgcdp _/0%:P]_(i <- ps)(eval_poly e i)))).

Definition rseq_poly ps := all rpoly ps.

Lemma rgcdpTs_qf (k : polyFformula F) (ps : seq polyF) :
  ( r, rpoly rqf (k r)) → rseq_poly psqf (rgcdpTs k ps).

Fixpoint rgdcop_recT (q : polyF) k (p : polyF) n :=
  if n is m.+1 then
    rgcdpT p (sizeT (fun sd
      if sd == 1%N then k p
      else rgcdpT p (rdivpT p (fun rrgdcop_recT q k r m)) q
    )) q
  else isnull (fun bk [::Const b%:R]) q.

Lemma rgdcop_recTP (k : polyFformula F) :
  ( p e, qf_eval e (k p) = qf_eval e (k (lift (eval_poly e p))))
  → p q n e, qf_eval e (rgdcop_recT p k q n)
    = qf_eval e (k (lift (rgdcop_rec (eval_poly e p) (eval_poly e q) n))).

Lemma rgdcop_recT_qf (p : polyF) (k : polyFformula F) (q : polyF) (n : nat) :
 ( r, rpoly rqf (k r)) →
 rpoly prpoly qqf (rgdcop_recT p k q n).

Definition rgdcopT q k p := sizeT (rgdcop_recT q k p) p.

Lemma rgdcopTP (k : polyFformula F) :
  ( p e, qf_eval e (k p) = qf_eval e (k (lift (eval_poly e p)))) →
   p q e, qf_eval e (rgdcopT p k q) =
                qf_eval e (k (lift (rgdcop (eval_poly e p) (eval_poly e q)))).

Lemma rgdcopT_qf (p : polyF) (k : polyFformula F) (q : polyF) :
  ( r, rpoly rqf (k r)) → rpoly prpoly qqf (rgdcopT p k q).

Definition ex_elim_seq (ps : seq polyF) (q : polyF) :=
  (rgcdpTs (rgdcopT q (sizeT (fun nBool (n != 1%N)))) ps).

Lemma ex_elim_seqP (ps : seq polyF) (q : polyF) (e : seq F) :
  let gp := (\big[@rgcdp _/0%:P]_(p <- ps)(eval_poly e p)) in
  qf_eval e (ex_elim_seq ps q) = (size (rgdcop (eval_poly e q) gp) != 1%N).

Lemma ex_elim_seq_qf (ps : seq polyF) (q : polyF) :
  rseq_poly psrpoly qqf (ex_elim_seq ps q).

Fixpoint abstrX (i : nat) (t : term F) :=
  match t with
    | (Var n) ⇒ if n == i then [::Const 0; Const 1] else [::t]
    | (Opp x) ⇒ opppT (abstrX i x)
    | (Add x y) ⇒ sumpT (abstrX i x) (abstrX i y)
    | (Mul x y) ⇒ mulpT (abstrX i x) (abstrX i y)
    | (NatMul x n) ⇒ natmulpT n (abstrX i x)
    | (Exp x n) ⇒ let ax := (abstrX i x) in
      iter n (mulpT ax) [::Const 1]
    | _[::t]
  end.

Lemma abstrXP (i : nat) (t : term F) (e : seq F) (x : F) :
  rterm t(eval_poly e (abstrX i t)).[x] = eval (set_nth 0 e i x) t.

Lemma rabstrX (i : nat) (t : term F) : rterm trpoly (abstrX i t).

Implicit Types tx ty : term F.

Lemma abstrX_mulM (i : nat) : {morph abstrX i : x y / Mul x y >-> mulpT x y}.

Lemma abstrX1 (i : nat) : abstrX i (Const 1) = [::Const 1].

Lemma eval_poly_mulM e : {morph eval_poly e : x y / mulpT x y >-> mul x y}.

Lemma eval_poly1 e : eval_poly e [::Const 1] = 1.

Notation abstrX_bigmul := (big_morph _ (abstrX_mulM _) (abstrX1 _)).
Notation eval_bigmul := (big_morph _ (eval_poly_mulM _) (eval_poly1 _)).
Notation bigmap_id := (big_map _ (fun _true) id).

Lemma rseq_poly_map (x : nat) (ts : seq (term F)) :
  all (@rterm _) tsrseq_poly (map (abstrX x) ts).

Definition ex_elim (x : nat) (pqs : seq (term F) × seq (term F)) :=
  ex_elim_seq (map (abstrX x) pqs.1)
  (abstrX x (\big[Mul/Const 1]_(q <- pqs.2) q)).

Lemma ex_elim_qf (x : nat) (pqs : seq (term F) × seq (term F)) :
  dnf_rterm pqsqf (ex_elim x pqs).

Lemma holds_conj : e i x ps, all (@rterm _) ps
  (holds (set_nth 0 e i x) (foldr (fun t : term FAnd (t == 0)) True ps)
   all ((@root _)^~ x) (map (eval_poly e \o abstrX i) ps)).

Lemma holds_conjn (e : seq F) (i : nat) (x : F) (ps : seq (term F)) :
  all (@rterm _) ps
  (holds (set_nth 0 e i x) (foldr (fun t : term FAnd (t != 0)) True ps)
  all (fun p~~root p x) (map (eval_poly e \o abstrX i) ps)).

Lemma holds_ex_elim : GRing.valid_QE_proj ex_elim.

Lemma wf_ex_elim : GRing.wf_QE_proj ex_elim.

Definition closed_fields_QEMixin :=
  QEdecFieldMixin wf_ex_elim holds_ex_elim.

End ClosedFieldQE.