(Joint Center)Library complex

(* (c) Copyright Microsoft Corporation and Inria. All rights reserved. *)
Require Import ssreflect ssrfun ssrbool eqtype ssrnat seq choice fintype.
Require Import bigop ssralg ssrint div ssrnum rat poly closed_field polyrcf.
Require Import matrix mxalgebra tuple mxpoly zmodp binomial realalg.

Import GRing.Theory Num.Theory.

Set Implicit Arguments.

Local Open Scope ring_scope.

Reserved Notation "x +i* y" (at level 40, left associativity, format "x +i* y").
Reserved Notation "x -i* y" (at level 40, left associativity, format "x -i* y").

Local Notation sgr := Num.sg.
Local Notation sqrtr := Num.sqrt.

CoInductive complex (R : Type) : Type := Complex { Re : R; Im : R }.

Definition real_complex_def (F : ringType) (phF : phant F) (x : F) := Complex x 0.
Notation real_complex F := (@real_complex_def _ (Phant F)).
Notation "x %:C" := (real_complex _ x)
  (at level 2, left associativity, format "x %:C") : ring_scope.
Notation "x +i* y" := (Complex x y) : ring_scope.
Notation "x -i* y" := (Complex x (- y)) : ring_scope.
Notation "x *i " := (Complex 0 x) (at level 8, format "x *i") : ring_scope.
Notation "''i'" := (Complex 0 1) : ring_scope.

Module ComplexEqChoice.
Section ComplexEqChoice.

Variable R : Type.
Notation C := (complex R).

Definition sqR_of_complex (x : C) := let: a +i× b := x in [::a; b].
Definition complex_of_sqR (x : seq R) :=
  if x is [:: a; b] then Some (a +i× b) else None.

Lemma complex_of_sqRK : pcancel sqR_of_complex complex_of_sqR.

End ComplexEqChoice.
End ComplexEqChoice.

Definition complex_eqMixin (R : eqType) :=
  PcanEqMixin (@ComplexEqChoice.complex_of_sqRK R).
Definition complex_choiceMixin (R : choiceType) :=
  PcanChoiceMixin (@ComplexEqChoice.complex_of_sqRK R).
Definition complex_countMixin (R : countType) :=
  PcanCountMixin (@ComplexEqChoice.complex_of_sqRK R).

Canonical Structure complex_eqType (R : eqType) :=
  EqType (complex R) (complex_eqMixin R).
Canonical Structure complex_choiceType (R : choiceType) :=
  ChoiceType (complex R) (complex_choiceMixin R).
Canonical Structure complex_countType (R : countType) :=
  CountType (complex R) (complex_countMixin R).

Lemma eq_complex : (R : eqType) (x y : complex R),
  (x == y) = (Re x == Re y) && (Im x == Im y).

Lemma complexr0 : (R : ringType) (x : R), x +i× 0 = x%:C.

Module ComplexField.
Section ComplexField.

Variable R : rcfType.
Local Notation C := (complex R).
Local Notation C0 := ((0 : R)%:C).
Local Notation C1 := ((1 : R)%:C).

Definition addc (x y : C) := let: a +i× b := x in let: c +i× d := y in
  (a + c) +i× (b + d).
Definition oppc (x : C) := let: a +i× b := x in (- a) +i× (- b).

Lemma addcC : commutative addc.
Lemma addcA : associative addc.

Lemma add0c : left_id C0 addc.

Lemma addNc : left_inverse C0 oppc addc.

Definition complex_ZmodMixin := ZmodMixin addcA addcC add0c addNc.
Canonical Structure complex_ZmodType := ZmodType C complex_ZmodMixin.

Definition mulc (x y : C) := let: a +i× b := x in let: c +i× d := y in
  ((a × c) - (b × d)) +i× ((a × d) + (b × c)).

Lemma mulcC : commutative mulc.

Lemma mulcA : associative mulc.

Definition invc (x : C) := let: a +i× b := x in let n2 := (a ^+ 2 + b ^+ 2) in
  (a / n2) -i× (b / n2).

Lemma mul1c : left_id C1 mulc.

Lemma mulc_addl : left_distributive mulc addc.

Lemma nonzero1c : C1 != C0.

Definition complex_comRingMixin :=
  ComRingMixin mulcA mulcC mul1c mulc_addl nonzero1c.
Canonical Structure complex_Ring := Eval hnf in RingType C complex_comRingMixin.
Canonical Structure complex_comRing := Eval hnf in ComRingType C mulcC.

Lemma mulVc : x, x != C0mulc (invc x) x = C1.

Lemma invc0 : invc C0 = C0.

Definition ComplexFieldUnitMixin := FieldUnitMixin mulVc invc0.
Canonical Structure complex_unitRing :=
  Eval hnf in UnitRingType C ComplexFieldUnitMixin.
Canonical Structure complex_comUnitRing := Eval hnf in [comUnitRingType of C].

Lemma field_axiom : GRing.Field.mixin_of complex_unitRing.

Definition ComplexFieldIdomainMixin := (FieldIdomainMixin field_axiom).
Canonical Structure complex_iDomain :=
  Eval hnf in IdomainType C (FieldIdomainMixin field_axiom).
Canonical Structure complex_fieldMixin := FieldType C field_axiom.

Ltac simpc := do ?
  [ rewrite -[(_ +i× _) - (_ +i× _)]/(_ +i× _)
  | rewrite -[(_ +i× _) + (_ +i× _)]/(_ +i× _)
  | rewrite -[(_ +i× _) × (_ +i× _)]/(_ +i× _)].

Lemma real_complex_is_rmorphism : rmorphism (real_complex R).

Canonical Structure real_complex_rmorphism :=
  RMorphism real_complex_is_rmorphism.
Canonical Structure real_complex_additive :=
  Additive real_complex_is_rmorphism.

Lemma Re_is_additive : additive (@Re R).

Canonical Structure Re_additive := Additive Re_is_additive.

Lemma Im_is_additive : additive (@Im R).

Canonical Structure Im_additive := Additive Im_is_additive.

Definition lec (x y : C) :=
  let: a +i× b := x in let: c +i× d := y in
    (d == b) && (a c).

Definition ltc (x y : C) :=
  let: a +i× b := x in let: c +i× d := y in
    (d == b) && (a < c).

Definition normc (x : C) : R :=
  let: a +i× b := x in sqrtr (a ^+ 2 + b ^+ 2).

Notation normC x := (normc x)%:C.

Lemma ltc0_add : x y, ltc 0 xltc 0 yltc 0 (x + y).

Lemma eq0_normc x : normc x = 0 → x = 0.

Lemma eq0_normC x : normC x = 0 → x = 0.

Lemma ge0_lec_total x y : lec 0 xlec 0 ylec x y || lec y x.

:TODO: put in ssralg ?
Lemma exprM (a b : R) : (a × b) ^+ 2 = a ^+ 2 × b ^+ 2.

Lemma normcM x y : normc (x × y) = normc x × normc y.

Lemma normCM x y : normC (x × y) = normC x × normC y.

Lemma subc_ge0 x y : lec 0 (y - x) = lec x y.

Lemma lec_def x y : lec x y = (normC (y - x) == y - x).

Lemma ltc_def x y : ltc x y = (y != x) && lec x y.

Lemma lec_normD x y : lec (normC (x + y)) (normC x + normC y).

Definition complex_POrderedMixin := NumMixin lec_normD ltc0_add eq0_normC
     ge0_lec_total normCM lec_def ltc_def.
Canonical Structure complex_numDomainType := NumDomainType C complex_POrderedMixin.

End ComplexField.
End ComplexField.

Canonical complex_ZmodType (R : rcfType) :=
  ZmodType (complex R) (ComplexField.complex_ZmodMixin R).
Canonical complex_Ring (R : rcfType) :=
  Eval hnf in RingType (complex R) (ComplexField.complex_comRingMixin R).
Canonical complex_comRing (R : rcfType) :=
  Eval hnf in ComRingType (complex R) (@ComplexField.mulcC R).
Canonical complex_unitRing (R : rcfType) :=
  Eval hnf in UnitRingType (complex R) (ComplexField.ComplexFieldUnitMixin R).
Canonical complex_comUnitRing (R : rcfType) :=
  Eval hnf in [comUnitRingType of (complex R)].
Canonical complex_iDomain (R : rcfType) :=
  Eval hnf in IdomainType (complex R) (FieldIdomainMixin (@ComplexField.field_axiom R)).
Canonical complex_fieldType (R : rcfType) :=
  FieldType (complex R) (@ComplexField.field_axiom R).
Canonical complex_numDomainType (R : rcfType) :=
  NumDomainType (complex R) (ComplexField.complex_POrderedMixin R).
Canonical complex_numFieldType (R : rcfType) :=
  [numFieldType of complex R].

Canonical ComplexField.real_complex_rmorphism.
Canonical ComplexField.real_complex_additive.
Canonical ComplexField.Re_additive.
Canonical ComplexField.Im_additive.

Definition conjc (R : ringType) (x : complex R) := let: a +i× b := x in a -i× b.
Notation "x ^*" := (conjc x) (at level 2, format "x ^*").

Ltac simpc := do ?
  [ rewrite -[(_ +i× _) - (_ +i× _)]/(_ +i× _)
  | rewrite -[(_ +i× _) + (_ +i× _)]/(_ +i× _)
  | rewrite -[(_ +i× _) × (_ +i× _)]/(_ +i× _)
  | rewrite -[(_ +i× _) (_ +i× _)]/((_ == _) && (_ _))
  | rewrite -[(_ +i× _) < (_ +i× _)]/((_ == _) && (_ < _))
  | rewrite -[`|_ +i× _|]/(sqrtr (_ + _))%:C
  | rewrite (mulrNN, mulrN, mulNr, opprB, opprD, mulr0, mul0r,
    subr0, sub0r, addr0, add0r, mulr1, mul1r, subrr, opprK, oppr0,
    eqxx) ].

Section ComplexTheory.

Variable R : rcfType.
Notation C := (complex R).

Lemma ReiNIm : x : C, Re (x × 'i) = - Im x.

Lemma ImiRe : x : C, Im (x × 'i) = Re x.

Lemma complexE x : x = (Re x)%:C + 'i × (Im x)%:C :> C.

Lemma real_complexE x : x%:C = x +i× 0 :> C.

Lemma sqr_i : 'i ^+ 2 = -1 :> C.

Lemma complexI : injective (real_complex R).

Lemma ler0c (x : R) : (0 x%:C) = (0 x).

Lemma lecE : x y : C, (x y) = (Im y == Im x) && (Re x Re y).

Lemma ltcE : x y : C, (x < y) = (Im y == Im x) && (Re x < Re y).

Lemma lecR : x y : R, (x%:C y%:C) = (x y).

Lemma ltcR : x y : R, (x%:C < y%:C) = (x < y).

Lemma conjc_is_rmorphism : rmorphism (@conjc R).

Canonical conjc_rmorphism := RMorphism conjc_is_rmorphism.
Canonical conjc_additive := Additive conjc_is_rmorphism.

Lemma conjcK : involutive (@conjc R).

Lemma mulcJ_ge0 (x : C) : 0 x × x ^*.

Lemma conjc_real (x : R) : x%:C^* = x%:C.

Lemma ReJ_add (x : C) : (Re x)%:C = (x + x^*) / 2%:R.

Lemma ImJ_sub (x : C) : (Im x)%:C = (x^* - x) / 2%:R × 'i.

Lemma ger0_Im (x : C) : 0 xIm x = 0.

Todo : extend theory of :
  • signed exponents

Lemma conj_ge0 : x : C, (0 x ^*) = (0 x).

Lemma conjc_nat : n, (n%:R : C)^* = n%:R.

Lemma conjc0 : (0 : C) ^* = 0.

Lemma conjc1 : (1 : C) ^* = 1.

Lemma conjc_eq0 : x : C, (x ^* == 0) = (x == 0).

Lemma conjc_inv: x : C, (x^-1)^* = (x^* )^-1.

Lemma complex_root_conj : (p : {poly C}) (x : C),
  root (map_poly (@conjc _) p) x = root p (x^*).

Lemma complex_algebraic_trans (T : comRingType) (toR : {rmorphism TR}) :
  integralRange toRintegralRange (real_complex R \o toR).

End ComplexTheory.

Section RcfDef.
Variable R : realFieldType. Notation C := (complex R).
Definition rcf_odd := forall (p : {poly R}), ~~odd (size p) -> {x | p. [x] = 0}. Definition rcf_square := forall x : R, {y | (0 <= y) && if 0 <= x then (y ^ 2 == x) else y == 0}.
Lemma rcf_odd_sqr_from_ivt : rcf_axiom R -> rcf_odd * rcf_square. Proof. move=> ivt. split. move=> p sp. move: (ivt p). admit. move=> x. case: (boolP (0 <= x)) (@ivt ('X^2 - x%:P) 0 (1 + x))=> px; last first. by move=> _; exists 0; rewrite lerr eqxx. case.

(Joint Center)by rewrite ler_paddr ?ler01.

(Joint Center)rewrite !horner_lin oppr_le0 px /=.

rewrite subr_ge0 (@ler_trans _ (1 + x)) //. by rewrite ler_paddl ?ler01 ?lerr. by rewrite ler_pemulr // addrC -subr_ge0 ?addrK // subr0 ler_paddl ?ler01.

(Joint Center)move=> y hy; rewrite /root !horner_lin; move/eqP.

move/(canRL (@addrNK _ _)); rewrite add0r=> <-. by exists y; case/andP: hy=> -> _; rewrite eqxx. Qed.
Lemma ivt_from_closed : GRing.ClosedField.axiom [ringType of C] -> rcf_axiom R. Proof. rewrite /GRing.ClosedField.axiom /= => hclosed. move=> p a b hab. Admitted.
Lemma closed_form_rcf_odd_sqr : rcf_odd -> rcf_square
  • > GRing.ClosedField.axiom [ringType of C].
Proof. Admitted.
Lemma closed_form_ivt : rcf_axiom R -> GRing.ClosedField.axiom [ringType of C]. Proof. move/rcf_odd_sqr_from_ivt; case. exact: closed_form_rcf_odd_sqr. Qed.
End RcfDef.

Section ComplexClosed.

Variable R : rcfType.
Local Notation C := (complex R).

Definition sqrtc (x : C) : C :=
  let: a +i× b := x in
  let sgr1 b := if b == 0 then 1 else sgr b in
  let r := sqrtr (a^+2 + b^+2) in
  (sqrtr ((r + a)/2%:R)) +i× (sgr1 b × sqrtr ((r - a)/2%:R)).

Lemma sqr_sqrtc : x, (sqrtc x) ^+ 2 = x.

Lemma sqrtc_sqrtr :
   (x : C), 0 xsqrtc x = (sqrtr (Re x))%:C.

Lemma sqrtc0 : sqrtc 0 = 0.

Lemma sqrtc1 : sqrtc 1 = 1.

Lemma sqrtN1 : sqrtc (-1) = 'i.

Lemma sqrtc_ge0 (x : C) : (0 sqrtc x) = (0 x).

Lemma sqrtc_eq0 (x : C) : (sqrtc x == 0) = (x == 0).

Lemma normcE x : `|x| = sqrtc (x × x^*).

Lemma sqr_normc (x : C) : (`|x| ^+ 2) = x × x^*.

Lemma normc_ge_Re (x : C) : `|Re x|%:C `|x|.

Lemma normcJ (x : C) : `|x^*| = `|x|.

Lemma invc_norm (x : C) : x^-1 = `|x|^-2 × x^*.

Lemma canonical_form (a b c : C) :
  a != 0 →
  let d := b ^+ 2 - 4%:R × a × c in
  let r1 := (- b - sqrtc d) / 2%:R / a in
  let r2 := (- b + sqrtc d) / 2%:R / a in
  a *: 'X^2 + b *: 'X + c%:P = a *: (('X - r1%:P) × ('X - r2%:P)).

Lemma monic_canonical_form (b c : C) :
  let d := b ^+ 2 - 4%:R × c in
  let r1 := (- b - sqrtc d) / 2%:R in
  let r2 := (- b + sqrtc d) / 2%:R in
  'X^2 + b *: 'X + c%:P = (('X - r1%:P) × ('X - r2%:P)).

Section extramx.
missing lemmas from matrix.v or mxalgebra.v

Lemma mul_mx_rowfree_eq0 (K : fieldType) (m n p: nat)
                         (W : 'M[K]_(m,n)) (V : 'M[K]_(n,p)) :
  row_free V(W ×m V == 0) = (W == 0).

Lemma sub_sums_genmxP (F : fieldType) (I : finType) (P : pred I) (m n : nat)
 (A : 'M[F]_(m, n)) (B_ : I'M_(m, n)) :
reflect ( u_ : I'M_m, A = \sum_(i | P i) u_ i ×m B_ i)
  (A \sum_(i | P i) <<B_ i>>)%MS.

Lemma mulmxP (K : fieldType) (m n : nat) (A B : 'M[K]_(m, n)) :
  reflect ( u : 'rV__, u ×m A = u ×m B) (A == B).

Section Skew.

Variable (K : numFieldType).

Implicit Types (phK : phant K) (n : nat).

Definition skew_vec n i j : 'rV[K]_(n × n) :=
   (mxvec ((delta_mx i j)) - (mxvec (delta_mx j i))).

Definition skew_def phK n : 'M[K]_(n × n) :=
  (\sum_(i | ((i.2 : 'I__) < (i.1 : 'I__))%N) <<skew_vec i.1 i.2>>)%MS.

Variable (n : nat).
Local Notation skew := (@skew_def (Phant K) n).

Lemma skew_direct_sum : mxdirect skew.
Hint Resolve skew_direct_sum.

Lemma rank_skew : \rank skew = (n × n.-1)./2.

Lemma skewP (M : 'rV_(n × n)) :
  reflect ((vec_mx M)^T = - vec_mx M) (M skew)%MS.

End Skew.

Notation skew K n := (@skew_def _ (Phant K) n).

Section Companion.

Variable (K : fieldType).

Lemma companion_subproof (p : {poly K}) :
  {M : 'M[K]_((size p).-1)| p \is monicchar_poly M = p}.

Definition companion (p : {poly K}) : 'M[K]_((size p).-1) :=
  projT1 (companion_subproof p).

Lemma companionK (p : {poly K}) : p \is monicchar_poly (companion p) = p.

End Companion.

Section Restriction.

Variable K : fieldType.
Variable m : nat.
Variables (V : 'M[K]_m).

Implicit Types f : 'M[K]_m.

Definition restrict f : 'M_(\rank V) := row_base V ×m f ×m (pinvmx (row_base V)).

Lemma stable_row_base f :
  (row_base V ×m f row_base V)%MS = (V ×m f V)%MS.

Lemma eigenspace_restrict f : (V ×m f V)%MS
   n a (W : 'M_(n, \rank V)),
  (W eigenspace (restrict f) a)%MS =
  (W ×m row_base V eigenspace f a)%MS.

Lemma eigenvalue_restrict f : (V ×m f V)%MS
  {subset eigenvalue (restrict f) eigenvalue f}.

Lemma restrictM : {in [pred f | (V ×m f V)%MS] &,
                      {morph restrict : f g / f ×m g}}.

End Restriction.

End extramx.
Notation skew K n := (@skew_def _ (Phant K) n).

Section Paper_HarmDerksen.

Following http://www.math.lsa.umich.edu/~hderksen/preprints/linalg.pdf quite literally except for Lemma5 where we don't use hermitian matrices. Instead we encode the morphism by hand in 'M[R]#_(n * n), which turns out to be very clumsy for formalizing commutation and the end of Lemma 4. Moreover, the Qed takes time, so it would be better good to formalize Herm C n and use it instead !

Implicit Types (K : fieldType).

Definition CommonEigenVec_def K (phK : phant K) (d r : nat) :=
   (m : nat) (V : 'M[K]_m), ~~ (d %| \rank V)
   (sf : seq 'M_m), size sf = r
  {in sf, f, (V ×m f V)%MS}
  {in sf &, f g, f ×m g = g ×m f}
  exists2 v : 'rV_m, (v != 0) & f, f \in sf
   a, (v eigenspace f a)%MS.
Notation CommonEigenVec K d r := (@CommonEigenVec_def _ (Phant K) d r).

Definition Eigen1Vec_def K (phK : phant K) (d : nat) :=
   (m : nat) (V : 'M[K]_m), ~~ (d %| \rank V)
   (f : 'M_m), (V ×m f V)%MS a, eigenvalue f a.
Notation Eigen1Vec K d := (@Eigen1Vec_def _ (Phant K) d).

Lemma Eigen1VecP (K : fieldType) (d : nat) :
  CommonEigenVec K d 1%N Eigen1Vec K d.

Lemma Lemma3 K d : Eigen1Vec K d r, CommonEigenVec K d r.+1.

Lemma Lemma4 r : CommonEigenVec R 2 r.+1.

Notation toC := (real_complex R).
Notation MtoC := (map_mx toC).

Lemma Lemma5 : Eigen1Vec C 2.

Lemma Lemma6 k r : CommonEigenVec C (2^k.+1) r.+1.

We enunciate a corollary of Theorem 7
Corollary Theorem7' (m : nat) (f : 'M[C]_m) : (0 < m)%N a, eigenvalue f a.

Lemma C_acf_axiom : GRing.ClosedField.axiom [ringType of C].

Definition C_decFieldMixin := closed_fields_QEMixin C_acf_axiom.
Canonical C_decField := DecFieldType C C_decFieldMixin.
Canonical C_closedField := ClosedFieldType C C_acf_axiom.

End Paper_HarmDerksen.

End ComplexClosed.

Definition complexalg := complex realalg.

Canonical complexalg_eqType := [eqType of complexalg].
Canonical complexalg_choiceType := [choiceType of complexalg].
Canonical complexalg_countype := [choiceType of complexalg].
Canonical complexalg_zmodType := [zmodType of complexalg].
Canonical complexalg_ringType := [ringType of complexalg].
Canonical complexalg_comRingType := [comRingType of complexalg].
Canonical complexalg_unitRingType := [unitRingType of complexalg].
Canonical complexalg_comUnitRingType := [comUnitRingType of complexalg].
Canonical complexalg_idomainType := [idomainType of complexalg].
Canonical complexalg_fieldType := [fieldType of complexalg].
Canonical complexalg_decDieldType := [decFieldType of complexalg].
Canonical complexalg_closedFieldType := [closedFieldType of complexalg].
Canonical complexalg_numDomainType := [numDomainType of complexalg].
Canonical complexalg_numFieldType := [numFieldType of complexalg].
Canonical complexalg_numClosedFieldType := [numClosedFieldType of complexalg].

Lemma complexalg_algebraic : integralRange (@ratr [unitRingType of complexalg]).