From da2f7b525eaef8fe78dded9e22e9f5502cd96ea7 Mon Sep 17 00:00:00 2001
From: Cezary Kaliszyk <cek@colo12-c703.uibk.ac.at>
Date: Sun, 24 Aug 2014 15:43:15 +0200
Subject: [PATCH] Update from HH
---
cvectors.ml | 2032 +++++++++++++++++++++++++++++++++++++++++++++++
em_model.ml | 391 +++++++++
frequency_equalities.ml | 311 ++++++++
make.ml | 2 +
primitive_rules.ml | 754 ++++++++++++++++++
tacticlib.ml | 81 ++
top.ml | 11 +
vectors_ext.ml | 835 +++++++++++++++++++
8 files changed, 4417 insertions(+), 0 deletions(-)
create mode 100644 cvectors.ml
create mode 100644 em_model.ml
create mode 100644 frequency_equalities.ml
create mode 100644 make.ml
create mode 100644 primitive_rules.ml
create mode 100644 tacticlib.ml
create mode 100644 top.ml
create mode 100644 vectors_ext.ml
diff --git a/cvectors.ml b/cvectors.ml
new file mode 100644
index 0000000..f6288d4
--- /dev/null
+++ b/cvectors.ml
@@ -0,0 +1,2032 @@
+(* Upadted for the latest version of HOL Light (JULY 2014) *)
+
+(* ========================================================================= *)
+(* A library for vectors of complex numbers. *)
+(* Much inspired from HOL-Light real vector library <"vectors.ml">. *)
+(* *)
+(* (c) Copyright, Sanaz Khan Afshar & Vincent Aravantinos 2011-13 *)
+(* Hardware Verification Group, *)
+(* Concordia University *)
+(* *)
+(* Contact: <s_khanaf@encs.concordia.ca> *)
+(* <vincent@encs.concordia.ca> *)
+(* *)
+(* Acknowledgements: *)
+(* - Harsh Singhal: n-dimensional dot product, utility theorems *)
+(* *)
+(* Last update: July 2013 *)
+(* *)
+(* ========================================================================= *)
+
+
+(* ========================================================================= *)
+(* ADDITIONS TO THE BASE LIBRARY *)
+(* ========================================================================= *)
+
+ (* ----------------------------------------------------------------------- *)
+ (* Additional tacticals *)
+ (* ----------------------------------------------------------------------- *)
+
+
+let SINGLE f x = f [x];;
+
+let distrib fs x = map (fun f -> f x) fs;;
+
+let DISTRIB ttacs x = EVERY (distrib ttacs x);;
+
+let REWRITE_TACS = MAP_EVERY (SINGLE REWRITE_TAC);;
+
+let GCONJUNCTS thm = map GEN_ALL (CONJUNCTS (SPEC_ALL thm));;
+
+ (* ----------------------------------------------------------------------- *)
+ (* Additions to the vectors library *)
+ (* ----------------------------------------------------------------------- *)
+
+let COMPONENT_LE_NORM_ALT = prove
+ (`!x:real^N i. 1 <= i /\ i <= dimindex (:N) ==> x$i <= norm x`,
+ MESON_TAC [REAL_ABS_LE;COMPONENT_LE_NORM;REAL_LE_TRANS]);;
+
+ (* ----------------------------------------------------------------------- *)
+ (* Additions to the library of complex numbers *)
+ (* ----------------------------------------------------------------------- *)
+
+(* Lemmas *)
+let RE_IM_NORM = prove
+ (`!x. Re x <= norm x /\ Im x <= norm x /\ abs(Re x) <= norm x
+ /\ abs(Im x) <= norm x`,
+ REWRITE_TAC[RE_DEF;IM_DEF] THEN GEN_TAC THEN REPEAT CONJ_TAC
+ THEN ((MATCH_MP_TAC COMPONENT_LE_NORM_ALT
+ THEN REWRITE_TAC[DIMINDEX_2] THEN ARITH_TAC) ORELSE SIMP_TAC [COMPONENT_LE_NORM]));;
+
+let [RE_NORM;IM_NORM;ABS_RE_NORM;ABS_IM_NORM] = GCONJUNCTS RE_IM_NORM;;
+
+let NORM_RE = prove
+ (`!x. &0 <= norm x + Re x /\ &0 <= norm x - Re x`,
+ GEN_TAC THEN MP_TAC (SPEC_ALL ABS_RE_NORM) THEN REAL_ARITH_TAC);;
+
+let [NORM_RE_ADD;NORM_RE_SUB] = GCONJUNCTS NORM_RE;;
+
+let NORM2_ADD_REAL = prove
+ (`!x y.
+ real x /\ real y ==> norm (x + ii * y) pow 2 = norm x pow 2 + norm y pow 2`,
+ SIMP_TAC[real;complex_norm;RE_ADD;IM_ADD;RE_MUL_II;IM_MUL_II;REAL_NEG_0;
+ REAL_ADD_LID;REAL_ADD_RID;REAL_POW_ZERO;ARITH_RULE `~(2=0)`;REAL_LE_POW_2;
+ SQRT_POW_2;REAL_LE_ADD]);;
+
+let COMPLEX_EQ_RCANCEL_IMP = GEN_ALL (MATCH_MP (MESON []
+ `(p <=> r \/ q) ==> (p /\ ~r ==> q) `) (SPEC_ALL COMPLEX_EQ_MUL_RCANCEL));;
+
+let COMPLEX_BALANCE_DIV_MUL = prove
+ (`!x y z t. ~(z=Cx(&0)) ==> (x = y/z * t <=> x*z = y * t)`,
+ REPEAT STRIP_TAC THEN POP_ASSUM (fun x ->
+ ASSUME_TAC (REWRITE_RULE[x] (SPEC_ALL COMPLEX_EQ_MUL_RCANCEL))
+ THEN ASSUME_TAC (REWRITE_RULE[x] (SPECL [`x:complex`;`z:complex`]
+ COMPLEX_DIV_RMUL)))
+ THEN SUBGOAL_THEN `x=y/z*t <=> x*z=(y/z*t)*z:complex` (SINGLE REWRITE_TAC)
+ THENL [ASM_REWRITE_TAC[];
+ REWRITE_TAC[SIMPLE_COMPLEX_ARITH `(y/z*t)*z=(y/z*z)*t:complex`]
+ THEN ASM_REWRITE_TAC[]]);;
+
+let CSQRT_MUL_LCX_LT = prove
+ (`!x y. &0 < x ==> csqrt(Cx x * y) = Cx(sqrt x) * csqrt y`,
+ REWRITE_TAC[csqrt;complex_mul;IM;RE;IM_CX;REAL_MUL_LZERO;REAL_ADD_RID;RE_CX;
+ REAL_SUB_RZERO]
+ THEN REPEAT STRIP_TAC THEN REPEAT COND_CASES_TAC
+ THEN FIRST_ASSUM (ASSUME_TAC o MATCH_MP REAL_LT_IMP_LE)
+ THEN ASM_SIMP_TAC[IM;RE;REAL_MUL_RZERO;SQRT_MUL]
+ THENL [
+ REPEAT (POP_ASSUM MP_TAC) THEN REWRITE_TAC[REAL_ENTIRE;REAL_MUL_POS_LE]
+ THEN REPEAT STRIP_TAC
+ THEN ASM_REWRITE_TAC[SQRT_0;REAL_MUL_LZERO;REAL_MUL_RZERO];
+ REPEAT (POP_ASSUM MP_TAC) THEN SIMP_TAC [REAL_ENTIRE]
+ THEN MESON_TAC [REAL_LT_IMP_NZ];
+ ASM_MESON_TAC [REAL_LE_MUL_EQ;REAL_ARITH `~(&0 <= y) = &0 > y`];
+ SIMP_TAC [REAL_NEG_RMUL] THEN REPEAT (POP_ASSUM MP_TAC)
+ THEN SIMP_TAC [REAL_ARITH `~(&0 <= y) = y < &0`]
+ THEN SIMP_TAC [GSYM REAL_NEG_GT0] THEN MESON_TAC[REAL_LT_IMP_LE;SQRT_MUL];
+ REPEAT (POP_ASSUM MP_TAC) THEN SIMP_TAC [REAL_ENTIRE]
+ THEN MESON_TAC [REAL_LT_IMP_NZ];
+ REPEAT (POP_ASSUM MP_TAC) THEN SIMP_TAC [REAL_ENTIRE]
+ THEN SIMP_TAC [DE_MORGAN_THM];
+ REPEAT (POP_ASSUM MP_TAC) THEN SIMP_TAC [REAL_ENTIRE]
+ THEN SIMP_TAC [DE_MORGAN_THM]; ALL_TAC] THENL [
+ SIMP_TAC [REAL_NEG_0;SQRT_0;REAL_MUL_RZERO];
+ ASM_MESON_TAC[REAL_ARITH `~(x<y /\ ~(x <=y))`];
+ ASM_MESON_TAC[REAL_ARITH `~(x<y /\ y<x)`];
+ ALL_TAC]
+ THEN REWRITE_TAC[GSYM (REWRITE_RULE[CX_DEF;complex_mul;RE;IM;
+ REAL_MUL_LZERO;REAL_ADD_RID;REAL_SUB_RZERO] COMPLEX_CMUL)]
+ THEN SIMP_TAC [NORM_MUL] THEN POP_ASSUM MP_TAC
+ THEN ASM_SIMP_TAC [GSYM REAL_ABS_REFL] THEN DISCH_TAC
+ THEN SIMP_TAC [REAL_ABS_MUL]
+ THEN ASM_SIMP_TAC [GSYM REAL_ABS_REFL]
+ THEN SIMP_TAC [GSYM REAL_ADD_LDISTRIB; GSYM REAL_SUB_LDISTRIB]
+ THEN SUBGOAL_THEN `(x*Im y) / (x*abs(Im y)) = Im y / abs(Im y)` ASSUME_TAC
+ THENL [
+ SIMP_TAC [real_div] THEN SIMP_TAC [REAL_INV_MUL]
+ THEN SIMP_TAC [GSYM REAL_MUL_ASSOC] THEN ONCE_REWRITE_TAC[REAL_MUL_AC]
+ THEN SUBGOAL_THEN `Im y * x * inv x * inv (abs(Im y)) =
+ Im y * (x * inv x) * inv (abs (Im y)) ` ASSUME_TAC
+ THENL [SIMP_TAC [REAL_MUL_AC]; ALL_TAC]
+ THEN ASM_SIMP_TAC[REAL_MUL_RINV;REAL_LT_IMP_NZ]
+ THEN SIMP_TAC [REAL_MUL_LID] THEN SIMP_TAC [REAL_MUL_AC];
+ ALL_TAC]
+ THEN ASM_SIMP_TAC[]
+ THEN SUBGOAL_THEN `sqrt x * Im y / abs(Im y) * sqrt ((norm y-Re y) / &2) =
+ Im y / abs (Im y) * sqrt x * sqrt ((norm y - Re y) / &2)` ASSUME_TAC
+ THENL [SIMP_TAC [REAL_MUL_AC]; ALL_TAC] THEN ASM_SIMP_TAC[]
+ THEN SUBGOAL_THEN `sqrt ((x * (norm y - Re y)) / &2) =
+ sqrt (x * (norm y - Re y)) / sqrt (&2)` ASSUME_TAC
+ THENL [
+ SIMP_TAC[SQRT_DIV] THEN CONJ_TAC THENL [
+ ASM_SIMP_TAC[REAL_LE_MUL_EQ;REAL_LT_IMP_LE] THEN SIMP_TAC[NORM_RE_SUB];
+ REAL_ARITH_TAC];
+ ALL_TAC]
+ THEN ASM_SIMP_TAC[] THEN SUBGOAL_THEN `sqrt ((norm y - Re y) / &2) =
+ sqrt (norm y - Re y) / sqrt (&2)` ASSUME_TAC
+ THENL [
+ SIMP_TAC[SQRT_DIV] THEN CONJ_TAC
+ THENL [SIMP_TAC [NORM_RE_SUB]; REAL_ARITH_TAC];
+ ALL_TAC ]
+ THEN ASM_SIMP_TAC[]
+ THEN SUBGOAL_THEN `sqrt ((x * (norm y + Re y)) / &2) =
+ sqrt (x * (norm y + Re y)) / sqrt (&2)` ASSUME_TAC
+ THENL [
+ SIMP_TAC[SQRT_DIV] THEN CONJ_TAC
+ THENL [
+ ASM_SIMP_TAC [REAL_LE_MUL_EQ;REAL_LT_IMP_LE]
+ THEN SIMP_TAC[NORM_RE_ADD];
+ REAL_ARITH_TAC];
+ ALL_TAC]
+ THEN SUBGOAL_THEN `sqrt ((norm y + Re y) / &2) =
+ sqrt (norm y + Re y) / sqrt (&2)` ASSUME_TAC
+ THENL [
+ SIMP_TAC[SQRT_DIV] THEN CONJ_TAC
+ THENL [SIMP_TAC[NORM_RE_ADD]; REAL_ARITH_TAC];
+ ALL_TAC]
+ THEN ASM_SIMP_TAC[] THEN SUBGOAL_THEN `&0 <= x` ASSUME_TAC
+ THENL [ ASM_SIMP_TAC [REAL_LT_IMP_LE]; ALL_TAC ]
+ THEN SIMP_TAC[COMPLEX_EQ] THEN SIMP_TAC[RE;IM] THEN CONJ_TAC
+ THENL [
+ SUBGOAL_THEN `sqrt x * sqrt (norm y + Re y) / sqrt (&2) =
+ (sqrt x * sqrt (norm y + Re y)) / sqrt (&2)` ASSUME_TAC
+ THENL [REAL_ARITH_TAC; ALL_TAC]
+ THEN ASM_MESON_TAC [SQRT_MUL;NORM_RE_ADD];
+ SUBGOAL_THEN `Im y/abs(Im y) * sqrt x * sqrt (norm y-Re y) / sqrt(&2) =
+ Im y/abs (Im y) * (sqrt x * sqrt (norm y - Re y))/sqrt(&2)` ASSUME_TAC
+ THENL [REAL_ARITH_TAC; ALL_TAC]
+ THEN ASM_MESON_TAC[SQRT_MUL;NORM_RE_SUB]]);;
+
+let CSQRT_MUL_LCX = prove
+ (`!x y. &0 <= x ==> csqrt(Cx x * y) = Cx(sqrt x) * csqrt y`,
+ REWRITE_TAC[REAL_LE_LT] THEN REPEAT STRIP_TAC
+ THEN ASM_SIMP_TAC[CSQRT_MUL_LCX_LT] THEN EXPAND_TAC "x"
+ THEN REWRITE_TAC[COMPLEX_MUL_LZERO;SQRT_0;CSQRT_0]);;
+
+let REAL_ADD_POW_2 = prove
+ (`!x y:real. (x+y) pow 2 = x pow 2 + y pow 2 + &2 * x * y`,
+ REAL_ARITH_TAC);;
+
+let COMPLEX_ADD_POW_2 = prove
+ (`!x y:complex. (x+y) pow 2 = x pow 2 + y pow 2 + Cx(&2) * x * y`,
+ REWRITE_TAC[COMPLEX_POW_2] THEN SIMPLE_COMPLEX_ARITH_TAC);;
+
+
+
+(* ----------------------------------------------------------------------- *)
+(* Additions to the topology library *)
+(* ----------------------------------------------------------------------- *)
+
+ prioritize_vector ();;
+
+(* Lemmas *)
+let FINITE_INTER_ENUM = prove
+ (`!s n. FINITE(s INTER (0..n))`,
+ MESON_TAC[FINITE_INTER;FINITE_NUMSEG]);;
+
+let NORM_PASTECART_GE1 = prove
+ (`!x y. norm x <= norm (pastecart x y)`,
+ MESON_TAC[FSTCART_PASTECART;NORM_FSTCART]);;
+
+let NORM_PASTECART_GE2 = prove
+ (`!x y. norm y <= norm (pastecart x y)`,
+ MESON_TAC[SNDCART_PASTECART;NORM_SNDCART]);;
+
+let LIM_PASTECART_EQ = prove
+ (`!net a b f:A->real^M g:A->real^N. (f --> a) net /\ (g --> b) net
+ <=> ((\x. pastecart (f x) (g x)) --> pastecart a b) net`,
+ REWRITE_TAC[MESON[] `(a <=> b) <=> (a ==> b) /\ (b ==> a)`;LIM_PASTECART;LIM;
+ MESON[]`(p\/q ==> (p \/ r) /\ (p \/ s)) <=> (~p /\ q ==> r /\ s)`;dist;
+ PASTECART_SUB]
+ THEN ASM_MESON_TAC[REAL_LET_TRANS;NORM_PASTECART_GE1;NORM_PASTECART_GE2]);;
+
+let SUMS_PASTECART = prove
+ (`!s f1:num->real^N f2:num->real^M l1 l2. (f1 sums l1) s /\ (f2 sums l2) s
+ <=> ((\x. pastecart (f1 x) (f2 x)) sums (pastecart l1 l2)) s`,
+ SIMP_TAC[sums;FINITE_INTER_ENUM;GSYM PASTECART_VSUM;LIM_PASTECART_EQ]);;
+
+let LINEAR_SUMS = prove(
+ `!s f l g. linear g ==> ((f sums l) s ==> ((g o f) sums (g l)) s)`,
+ SIMP_TAC[sums;FINITE_INTER_ENUM;GSYM LINEAR_VSUM;
+ REWRITE_RULE[o_DEF;CONTINUOUS_AT_SEQUENTIALLY] LINEAR_CONTINUOUS_AT]);;
+
+ (* ----------------------------------------------------------------------- *)
+ (* Embedding of reals in complex numbers *)
+ (* ----------------------------------------------------------------------- *)
+
+let real_of_complex = new_definition `real_of_complex c = @r. c = Cx r`;;
+
+let REAL_OF_COMPLEX = prove
+ (`!c. real c ==> Cx(real_of_complex c) = c`,
+ MESON_TAC[REAL;real_of_complex]);;
+
+let REAL_OF_COMPLEX_RE = prove
+ (`!c. real c ==> real_of_complex c = Re c`,
+ MESON_TAC[RE_CX;REAL_OF_COMPLEX]);;
+
+let REAL_OF_COMPLEX_CX = prove
+ (`!r. real_of_complex (Cx r) = r`,
+ SIMP_TAC[REAL_CX;REAL_OF_COMPLEX_RE;RE_CX]);;
+
+let REAL_OF_COMPLEX_NORM = prove
+ (`!c. real c ==> norm c = abs (real_of_complex c)`,
+ MESON_TAC[REAL_NORM;REAL_OF_COMPLEX_RE]);;
+
+let REAL_OF_COMPLEX_ADD = prove
+ (`!x y. real x /\ real y
+ ==> real_of_complex (x+y) = real_of_complex x + real_of_complex y`,
+ MESON_TAC[REAL_ADD;REAL_OF_COMPLEX_RE;RE_ADD]);;
+
+let REAL_MUL = prove
+ (`!x y. real x /\ real y ==> real (x*y)`,
+ REWRITE_TAC[real] THEN SIMPLE_COMPLEX_ARITH_TAC);;
+
+let REAL_OF_COMPLEX_MUL = prove(
+ `!x y. real x /\ real y
+ ==> real_of_complex (x*y) = real_of_complex x * real_of_complex y`,
+ MESON_TAC[REAL_MUL;REAL_OF_COMPLEX;CX_MUL;REAL_OF_COMPLEX_CX]);;
+
+let REAL_OF_COMPLEX_0 = prove(
+ `!x. real x ==> (real_of_complex x = &0 <=> x = Cx(&0))`,
+ REWRITE_TAC[REAL_EXISTS] THEN REPEAT STRIP_TAC
+ THEN ASM_SIMP_TAC[REAL_OF_COMPLEX_CX;CX_INJ]);;
+
+let REAL_COMPLEX_ADD_CNJ = prove(
+ `!x. real(cnj x + x) /\ real(x + cnj x)`,
+ REWRITE_TAC[COMPLEX_ADD_CNJ;REAL_CX]);;
+
+(* TODO
+ *let RE_EQ_NORM = prove(`!x. Re x = norm x <=> real x /\ &0 <= real_of_complex x`,
+ *)
+
+ (* ----------------------------------------------------------------------- *)
+ (* Additions to the vectors library *)
+ (* ----------------------------------------------------------------------- *)
+
+let vector_const = new_definition
+ `vector_const (k:A) :A^N = lambda i. k`;;
+let vector_map = new_definition
+ `vector_map (f:A->B) (v:A^N) :B^N = lambda i. f(v$i)`;;
+let vector_map2 = new_definition
+ `vector_map2 (f:A->B->C) (v1:A^N) (v2:B^N) :C^N =
+ lambda i. f (v1$i) (v2$i)`;;
+let vector_map3 = new_definition
+ `vector_map3 (f:A->B->C->D) (v1:A^N) (v2:B^N) (v3:C^N) :D^N =
+ lambda i. f (v1$i) (v2$i) (v3$i)`;;
+
+let FINITE_INDEX_INRANGE_2 = prove
+ (`!i. ?k. 1 <= k /\ k <= dimindex(:N) /\ (!x:A^N. x$i = x$k)
+ /\ (!x:B^N. x$i = x$k) /\ (!x:C^N. x$i = x$k) /\ (!x:D^N. x$i = x$k)`,
+ REWRITE_TAC[finite_index] THEN MESON_TAC[FINITE_INDEX_WORKS]);;
+
+let COMPONENT_TAC x =
+ REPEAT GEN_TAC THEN CHOOSE_TAC (SPEC_ALL FINITE_INDEX_INRANGE_2)
+ THEN ASM_SIMP_TAC[x;LAMBDA_BETA];;
+
+let VECTOR_CONST_COMPONENT = prove
+ (`!i k. ((vector_const k):A^N)$i = k`,
+ COMPONENT_TAC vector_const);;
+let VECTOR_MAP_COMPONENT = prove
+ (`!i f:A->B v:A^N. (vector_map f v)$i = f (v$i)`,
+ COMPONENT_TAC vector_map);;
+let VECTOR_MAP2_COMPONENT = prove
+ (`!i f:A->B->C v1:A^N v2. (vector_map2 f v1 v2)$i = f (v1$i) (v2$i)`,
+ COMPONENT_TAC vector_map2);;
+let VECTOR_MAP3_COMPONENT = prove(
+ `!i f:A->B->C->D v1:A^N v2 v3. (vector_map3 f v1 v2 v3)$i =
+ f (v1$i) (v2$i) (v3$i)`,
+ COMPONENT_TAC vector_map3);;
+
+let COMMON_TAC =
+ REWRITE_TAC[vector_const;vector_map;vector_map2;vector_map3]
+ THEN ONCE_REWRITE_TAC[CART_EQ] THEN SIMP_TAC[LAMBDA_BETA;o_DEF];;
+
+let VECTOR_MAP_VECTOR_CONST = prove
+ (`!f:A->B k. vector_map f ((vector_const k):A^N) = vector_const (f k)`,
+ COMMON_TAC);;
+
+let VECTOR_MAP_VECTOR_MAP = prove
+ (`!f:A->B g:C->A v:C^N.
+ vector_map f (vector_map g v) = vector_map (f o g) v`,
+ COMMON_TAC);;
+
+let VECTOR_MAP_VECTOR_MAP2 = prove
+ (`!f:A->B g:C->D->A u v:D^N.
+ vector_map f (vector_map2 g u v) = vector_map2 (\x y. f (g x y)) u v`,
+ COMMON_TAC);;
+
+let VECTOR_MAP2_LVECTOR_CONST = prove
+ (`!f:A->B->C k v:B^N.
+ vector_map2 f (vector_const k) v = vector_map (f k) v`,
+ COMMON_TAC);;
+
+let VECTOR_MAP2_RVECTOR_CONST = prove
+ (`!f:A->B->C k v:A^N.
+ vector_map2 f v (vector_const k) = vector_map (\x. f x k) v`,
+ COMMON_TAC);;
+
+let VECTOR_MAP2_LVECTOR_MAP = prove
+ (`!f:A->B->C g:D->A v1 v2:B^N.
+ vector_map2 f (vector_map g v1) v2 = vector_map2 (f o g) v1 v2`,
+ COMMON_TAC);;
+
+let VECTOR_MAP2_RVECTOR_MAP = prove
+ (`!f:A->B->C g:D->B v1 v2:D^N.
+ vector_map2 f v1 (vector_map g v2) = vector_map2 (\x y. f x (g y)) v1 v2`,
+ COMMON_TAC);;
+
+let VECTOR_MAP2_LVECTOR_MAP2 = prove
+ (`!f:A->B->C g:D->E->A v1 v2 v3:B^N.
+ vector_map2 f (vector_map2 g v1 v2) v3 =
+ vector_map3 (\x y. f (g x y)) v1 v2 v3`,
+ COMMON_TAC);;
+
+let VECTOR_MAP2_RVECTOR_MAP2 = prove(
+ `!f:A->B->C g:D->E->B v1 v2 v3:E^N.
+ vector_map2 f v1 (vector_map2 g v2 v3) =
+ vector_map3 (\x y z. f x (g y z)) v1 v2 v3`,
+ COMMON_TAC);;
+
+let VECTOR_MAP_SIMP_TAC = REWRITE_TAC[
+ VECTOR_MAP_VECTOR_CONST;VECTOR_MAP2_LVECTOR_CONST;
+ VECTOR_MAP2_RVECTOR_CONST;VECTOR_MAP_VECTOR_MAP;VECTOR_MAP2_RVECTOR_MAP;
+ VECTOR_MAP2_LVECTOR_MAP;VECTOR_MAP2_RVECTOR_MAP2;VECTOR_MAP2_LVECTOR_MAP2;
+ VECTOR_MAP_VECTOR_MAP2];;
+
+let VECTOR_MAP_PROPERTY_TAC fs prop =
+ REWRITE_TAC fs THEN VECTOR_MAP_SIMP_TAC THEN ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[VECTOR_MAP_COMPONENT;VECTOR_MAP2_COMPONENT;
+ VECTOR_MAP3_COMPONENT;VECTOR_CONST_COMPONENT;o_DEF;prop];;
+
+let VECTOR_MAP_PROPERTY thm f prop =
+ prove(thm,VECTOR_MAP_PROPERTY_TAC f prop);;
+
+let COMPLEX_VECTOR_MAP = prove
+ (`!f:complex->complex g. f = vector_map g
+ <=> !z. f z = complex (g (Re z), g (Im z))`,
+ REWRITE_TAC[vector_map;FUN_EQ_THM;complex] THEN REPEAT (GEN_TAC ORELSE EQ_TAC)
+ THEN ASM_SIMP_TAC[CART_EQ;DIMINDEX_2;FORALL_2;LAMBDA_BETA;VECTOR_2;RE_DEF;IM_DEF]);;
+
+let COMPLEX_NEG_IS_A_MAP = prove
+ (`(--):complex->complex = vector_map ((--):real->real)`,
+ REWRITE_TAC[COMPLEX_VECTOR_MAP;complex_neg]);;
+
+let VECTOR_NEG_IS_A_MAP = prove
+ (`(--):real^N->real^N = vector_map ((--):real->real)`,
+ REWRITE_TAC[FUN_EQ_THM;CART_EQ;VECTOR_NEG_COMPONENT;VECTOR_MAP_COMPONENT]);;
+
+let VECTOR_MAP_VECTOR_MAP_ALT = prove
+ (`!f:A^N->B^N g:C^N->A^N f' g'. f = vector_map f' /\ g = vector_map g' ==>
+ f o g = vector_map (f' o g')`,
+ SIMP_TAC[o_DEF;FUN_EQ_THM;VECTOR_MAP_VECTOR_MAP]);;
+
+let COMPLEX_VECTOR_MAP2 = prove
+ (`!f:complex->complex->complex g. f = vector_map2 g <=>
+ !z1 z2. f z1 z2 = complex (g (Re z1) (Re z2), g (Im z1) (Im z2))`,
+ REWRITE_TAC[vector_map2;FUN_EQ_THM;complex]
+ THEN REPEAT (GEN_TAC ORELSE EQ_TAC)
+ THEN ASM_SIMP_TAC[CART_EQ;DIMINDEX_2;FORALL_2;LAMBDA_BETA;VECTOR_2;RE_DEF;
+ IM_DEF]);;
+
+let VECTOR_MAP2_RVECTOR_MAP_ALT = prove(
+ `!f:complex->complex->complex g:complex->complex f' g'.
+ f = vector_map2 f' /\ g = vector_map g'
+ ==> (\x y. f x (g y)) = vector_map2 (\x y. f' x (g' y))`,
+ SIMP_TAC[FUN_EQ_THM;VECTOR_MAP2_RVECTOR_MAP]);;
+
+let COMPLEX_ADD_IS_A_MAP = prove
+ (`(+):complex->complex->complex = vector_map2 ((+):real->real->real)`,
+ REWRITE_TAC[COMPLEX_VECTOR_MAP2;complex_add]);;
+
+let VECTOR_ADD_IS_A_MAP = prove
+ (`(+):real^N->real^N->real^N = vector_map2 ((+):real->real->real)`,
+ REWRITE_TAC[FUN_EQ_THM;CART_EQ;VECTOR_ADD_COMPONENT;VECTOR_MAP2_COMPONENT]);;
+
+let COMPLEX_SUB_IS_A_MAP = prove
+ (`(-):complex->complex->complex = vector_map2 ((-):real->real->real)`,
+ ONCE_REWRITE_TAC[prove(`(-) = \x y:complex. x-y`,REWRITE_TAC[FUN_EQ_THM])]
+ THEN ONCE_REWRITE_TAC[prove(`(-) = \x y:real. x-y`,REWRITE_TAC[FUN_EQ_THM])]
+ THEN REWRITE_TAC[complex_sub;real_sub]
+ THEN MATCH_MP_TAC VECTOR_MAP2_RVECTOR_MAP_ALT
+ THEN REWRITE_TAC[COMPLEX_NEG_IS_A_MAP;COMPLEX_ADD_IS_A_MAP]);;
+
+let VECTOR_SUB_IS_A_MAP = prove
+ (`(-):real^N->real^N->real^N = vector_map2 ((-):real->real->real)`,
+ REWRITE_TAC[FUN_EQ_THM;CART_EQ;VECTOR_SUB_COMPONENT;VECTOR_MAP2_COMPONENT]);;
+
+let COMMON_TAC x =
+ SIMP_TAC[CART_EQ;pastecart;x;LAMBDA_BETA] THEN REPEAT STRIP_TAC
+ THEN REPEAT COND_CASES_TAC THEN POP_ASSUM MP_TAC THEN REWRITE_TAC[]
+ THEN SUBGOAL_THEN `1<= i-dimindex(:N) /\ i-dimindex(:N) <= dimindex(:M)`
+ ASSUME_TAC
+ THEN ASM_SIMP_TAC[LAMBDA_BETA]
+ THEN REPEAT (POP_ASSUM (MP_TAC o REWRITE_RULE[DIMINDEX_FINITE_SUM]))
+ THEN ARITH_TAC;;
+
+let PASTECART_VECTOR_MAP = prove
+ (`!f:A->B x:A^N y:A^M.
+ pastecart (vector_map f x) (vector_map f y) =
+ vector_map f (pastecart x y)`,
+ COMMON_TAC vector_map);;
+
+let PASTECART_VECTOR_MAP2 = prove
+ (`!f:A->B->C x1:A^N x2 y1:A^M y2.
+ pastecart (vector_map2 f x1 x2) (vector_map2 f y1 y2)
+ = vector_map2 f (pastecart x1 y1) (pastecart x2 y2)`,
+ COMMON_TAC vector_map2);;
+
+let vector_zip = new_definition
+ `vector_zip (v1:A^N) (v2:B^N) : (A#B)^N = lambda i. (v1$i,v2$i)`;;
+
+let VECTOR_ZIP_COMPONENT = prove
+ (`!i v1:A^N v2:B^N. (vector_zip v1 v2)$i = (v1$i,v2$i)`,
+ REPEAT GEN_TAC THEN CHOOSE_TAC (INST_TYPE [`:A#B`,`:C`] (SPEC_ALL
+ FINITE_INDEX_INRANGE_2)) THEN ASM_SIMP_TAC[vector_zip;LAMBDA_BETA]);;
+
+let vector_unzip = new_definition
+ `vector_unzip (v:(A#B)^N):(A^N)#(B^N) = vector_map FST v,vector_map SND v`;;
+
+let VECTOR_UNZIP_COMPONENT = prove
+ (`!i v:(A#B)^N. (FST (vector_unzip v))$i = FST (v$i)
+ /\ (SND (vector_unzip v))$i = SND (v$i)`,
+ REWRITE_TAC[vector_unzip;VECTOR_MAP_COMPONENT]);;
+
+let VECTOR_MAP2_AS_VECTOR_MAP = prove
+ (`!f:A->B->C v1:A^N v2:B^N.
+ vector_map2 f v1 v2 = vector_map (UNCURRY f) (vector_zip v1 v2)`,
+ REWRITE_TAC[CART_EQ;VECTOR_MAP2_COMPONENT;VECTOR_MAP_COMPONENT;
+ VECTOR_ZIP_COMPONENT;UNCURRY_DEF]);;
+
+
+
+(* ========================================================================= *)
+(* BASIC ARITHMETIC *)
+(* ========================================================================= *)
+
+make_overloadable "%" `:A-> B-> B`;;
+
+let prioritize_cvector () =
+ overload_interface("--",`(cvector_neg):complex^N->complex^N`);
+ overload_interface("+",`(cvector_add):complex^N->complex^N->complex^N`);
+ overload_interface("-",`(cvector_sub):complex^N->complex^N->complex^N`);
+ overload_interface("%",`(cvector_mul):complex->complex^N->complex^N`);;
+
+let cvector_zero = new_definition
+ `cvector_zero:complex^N = vector_const (Cx(&0))`;;
+
+let cvector_neg = new_definition
+ `cvector_neg :complex^N->complex^N = vector_map (--)`;;
+
+let cvector_add = new_definition
+ `cvector_add :complex^N->complex^N->complex^N = vector_map2 (+)`;;
+
+let cvector_sub = new_definition
+ `cvector_sub :complex^N->complex^N->complex^N = vector_map2 (-)`;;
+
+let cvector_mul = new_definition
+ `(cvector_mul:complex->complex^N->complex^N) a = vector_map (( * ) a)`;;
+
+overload_interface("%",`(%):real->real^N->real^N`);;
+prioritize_cvector ();;
+
+let CVECTOR_ZERO_COMPONENT = prove
+ (`!i. (cvector_zero:complex^N)$i = Cx(&0)`,
+ REWRITE_TAC[cvector_zero;VECTOR_CONST_COMPONENT]);;
+
+let CVECTOR_NON_ZERO = prove
+ (`!x:complex^N. ~(x=cvector_zero)
+ <=> ?i. 1 <= i /\ i <= dimindex(:N) /\ ~(x$i = Cx(&0))`,
+ GEN_TAC THEN GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV) [CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_ZERO_COMPONENT] THEN MESON_TAC[]);;
+
+let CVECTOR_ADD_COMPONENT = prove
+ (`!X Y:complex^N i. ((X + Y)$i = X$i + Y$i)`,
+ REWRITE_TAC[cvector_add;VECTOR_MAP2_COMPONENT]);;
+
+let CVECTOR_SUB_COMPONENT = prove
+ (`!X:complex^N Y i. ((X - Y)$i = X$i - Y$i)`,
+ REWRITE_TAC[cvector_sub;VECTOR_MAP2_COMPONENT]);;
+
+let CVECTOR_NEG_COMPONENT = prove
+ (`!X:complex^N i. ((--X)$i = --(X$i))`,
+ REWRITE_TAC[cvector_neg;VECTOR_MAP_COMPONENT]);;
+
+let CVECTOR_MUL_COMPONENT = prove
+ (`!c:complex X:complex^N i. ((c % X)$i = c * X$i)`,
+ REWRITE_TAC[cvector_mul;VECTOR_MAP_COMPONENT]);;
+
+(* Simple generic tactic adapted from VECTOR_ARITH_TAC *)
+let CVECTOR_ARITH_TAC =
+ let RENAMED_LAMBDA_BETA th =
+ if fst(dest_fun_ty(type_of(funpow 3 rand (concl th)))) = aty
+ then INST_TYPE [aty,bty; bty,aty] LAMBDA_BETA else LAMBDA_BETA
+ in
+ POP_ASSUM_LIST(K ALL_TAC) THEN
+ REPEAT(GEN_TAC ORELSE CONJ_TAC ORELSE DISCH_TAC ORELSE EQ_TAC) THEN
+ REPEAT(POP_ASSUM MP_TAC) THEN REWRITE_TAC[IMP_IMP; GSYM CONJ_ASSOC] THEN
+ GEN_REWRITE_TAC ONCE_DEPTH_CONV [CART_EQ] THEN
+ REWRITE_TAC[AND_FORALL_THM] THEN TRY EQ_TAC THEN
+ TRY(MATCH_MP_TAC MONO_FORALL) THEN GEN_TAC THEN
+ REWRITE_TAC[TAUT `(a ==> b) /\ (a ==> c) <=> a ==> b /\ c`;
+ TAUT `(a ==> b) \/ (a ==> c) <=> a ==> b \/ c`] THEN
+ TRY(MATCH_MP_TAC(TAUT `(a ==> b ==> c) ==> (a ==> b) ==> (a ==> c)`)) THEN
+ REWRITE_TAC[cvector_zero;cvector_add; cvector_sub; cvector_neg; cvector_mul; vector_map;vector_map2;vector_const] THEN
+ DISCH_THEN(fun th -> REWRITE_TAC[MATCH_MP(RENAMED_LAMBDA_BETA th) th]) THEN
+ SIMPLE_COMPLEX_ARITH_TAC;;
+
+let CVECTOR_ARITH tm = prove(tm,CVECTOR_ARITH_TAC);;
+
+(* ========================================================================= *)
+(* VECTOR SPACE AXIOMS AND ADDITIONAL BASIC RESULTS *)
+(* ========================================================================= *)
+
+let CVECTOR_MAP_PROPERTY thm =
+ VECTOR_MAP_PROPERTY thm [cvector_zero;cvector_add;cvector_sub;cvector_neg;
+ cvector_mul];;
+
+let CVECTOR_ADD_SYM = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. x + y = y + x`
+ COMPLEX_ADD_SYM;;
+
+let CVECTOR_ADD_ASSOC = CVECTOR_MAP_PROPERTY
+ `!x y z:complex^N. x + (y + z) = (x + y) + z`
+ COMPLEX_ADD_ASSOC;;
+
+let CVECTOR_ADD_ID = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. x + cvector_zero = x /\ cvector_zero + x = x`
+ (CONJ COMPLEX_ADD_RID COMPLEX_ADD_LID);;
+
+let [CVECTOR_ADD_RID;CVECTOR_ADD_LID] = GCONJUNCTS CVECTOR_ADD_ID;;
+
+let CVECTOR_ADD_INV = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. x + -- x = cvector_zero /\ --x + x = cvector_zero`
+ (CONJ COMPLEX_ADD_RINV COMPLEX_ADD_LINV);;
+
+let CVECTOR_MUL_ASSOC = CVECTOR_MAP_PROPERTY
+ `!a b x:complex^N. a % (b % x) = (a * b) % x`
+ COMPLEX_MUL_ASSOC;;
+
+let CVECTOR_SUB_LDISTRIB = CVECTOR_MAP_PROPERTY
+ `!c x y:complex^N. c % (x - y) = c % x - c % y`
+ COMPLEX_SUB_LDISTRIB;;
+
+let CVECTOR_SCALAR_RDIST = CVECTOR_MAP_PROPERTY
+ `!a b x:complex^N. (a + b) % x = a % x + b % x`
+ COMPLEX_ADD_RDISTRIB;;
+
+let CVECTOR_MUL_ID = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. Cx(&1) % x = x`
+ COMPLEX_MUL_LID;;
+
+let CVECTOR_SUB_REFL = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. x - x = cvector_zero`
+ COMPLEX_SUB_REFL;;
+
+let CVECTOR_SUB_RADD = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. x - (x + y) = --y`
+ COMPLEX_ADD_SUB2;;
+
+let CVECTOR_NEG_SUB = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. --(x - y) = y - x`
+ COMPLEX_NEG_SUB;;
+
+let CVECTOR_SUB_EQ = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. (x - y = cvector_zero) <=> (x = y)`
+ COMPLEX_SUB_0;;
+
+let CVECTOR_MUL_LZERO = CVECTOR_MAP_PROPERTY
+ `!x. Cx(&0) % x = cvector_zero`
+ COMPLEX_MUL_LZERO;;
+
+let CVECTOR_SUB_ADD = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. (x - y) + y = x`
+ COMPLEX_SUB_ADD;;
+
+let CVECTOR_SUB_ADD2 = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. y + (x - y) = x`
+ COMPLEX_SUB_ADD2;;
+
+let CVECTOR_ADD_LDISTRIB = CVECTOR_MAP_PROPERTY
+ `!c x y:complex^N. c % (x + y) = c % x + c % y`
+ COMPLEX_ADD_LDISTRIB;;
+
+let CVECTOR_ADD_RDISTRIB = CVECTOR_MAP_PROPERTY
+ `!a b x:complex^N. (a + b) % x = a % x + b % x`
+ COMPLEX_ADD_RDISTRIB;;
+
+let CVECTOR_SUB_RDISTRIB = CVECTOR_MAP_PROPERTY
+ `!a b x:complex^N. (a - b) % x = a % x - b % x`
+ COMPLEX_SUB_RDISTRIB;;
+
+let CVECTOR_ADD_SUB = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. (x + y:complex^N) - x = y`
+ COMPLEX_ADD_SUB;;
+
+let CVECTOR_EQ_ADDR = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. (x + y = x) <=> (y = cvector_zero)`
+ COMPLEX_EQ_ADD_LCANCEL_0;;
+
+let CVECTOR_SUB = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. x - y = x + --(y:complex^N)`
+ complex_sub;;
+
+let CVECTOR_SUB_RZERO = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. x - cvector_zero = x`
+ COMPLEX_SUB_RZERO;;
+
+let CVECTOR_MUL_RZERO = CVECTOR_MAP_PROPERTY
+ `!c:complex. c % cvector_zero = cvector_zero`
+ COMPLEX_MUL_RZERO;;
+
+let CVECTOR_MUL_LZERO = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. Cx(&0) % x = cvector_zero`
+ COMPLEX_MUL_LZERO;;
+
+let CVECTOR_MUL_EQ_0 = prove
+ (`!a:complex x:complex^N.
+ (a % x = cvector_zero <=> a = Cx(&0) \/ x = cvector_zero)`,
+ REPEAT STRIP_TAC THEN EQ_TAC THENL [
+ ASM_CASES_TAC `a=Cx(&0)` THENL [
+ ASM_REWRITE_TAC[];
+ GEN_REWRITE_TAC (RATOR_CONV o DEPTH_CONV) [CART_EQ]
+ THEN ASM_REWRITE_TAC[CVECTOR_MUL_COMPONENT;CVECTOR_ZERO_COMPONENT;
+ COMPLEX_ENTIRE]
+ THEN GEN_REWRITE_TAC (RAND_CONV o DEPTH_CONV) [CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_ZERO_COMPONENT];
+ ];
+ REPEAT STRIP_TAC THEN ASM_REWRITE_TAC[CVECTOR_MUL_RZERO;CVECTOR_MUL_LZERO];
+ ]);;
+
+let CVECTOR_NEG_MINUS1 = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. --x = (--(Cx(&1))) % x`
+ (GSYM COMPLEX_NEG_MINUS1);;
+
+let CVECTOR_SUB_LZERO = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. cvector_zero - x = --x`
+ COMPLEX_SUB_LZERO;;
+
+let CVECTOR_NEG_NEG = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. --(--(x:complex^N)) = x`
+ COMPLEX_NEG_NEG;;
+
+let CVECTOR_MUL_LNEG = CVECTOR_MAP_PROPERTY
+ `!c x:complex^N. --c % x = --(c % x)`
+ COMPLEX_MUL_LNEG;;
+
+let CVECTOR_MUL_RNEG = CVECTOR_MAP_PROPERTY
+ `!c x:complex^N. c % --x = --(c % x)`
+ COMPLEX_MUL_RNEG;;
+
+let CVECTOR_NEG_0 = CVECTOR_MAP_PROPERTY
+ `--cvector_zero = cvector_zero`
+ COMPLEX_NEG_0;;
+
+let CVECTOR_NEG_EQ_0 = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. --x = cvector_zero <=> x = cvector_zero`
+ COMPLEX_NEG_EQ_0;;
+
+let CVECTOR_ADD_AC = prove
+ (`!x y z:complex^N.
+ (x + y = y + x) /\ ((x + y) + z = x + y + z) /\ (x + y + z = y + x + z)`,
+ MESON_TAC[CVECTOR_ADD_SYM;CVECTOR_ADD_ASSOC]);;
+
+let CVECTOR_MUL_LCANCEL = prove
+ (`!a x y:complex^N. a % x = a % y <=> a = Cx(&0) \/ x = y`,
+ MESON_TAC[CVECTOR_MUL_EQ_0;CVECTOR_SUB_LDISTRIB;CVECTOR_SUB_EQ]);;
+
+let CVECTOR_MUL_RCANCEL = prove
+ (`!a b x:complex^N. a % x = b % x <=> a = b \/ x = cvector_zero`,
+ MESON_TAC[CVECTOR_MUL_EQ_0;CVECTOR_SUB_RDISTRIB;COMPLEX_SUB_0;CVECTOR_SUB_EQ]);;
+
+
+(* ========================================================================= *)
+(* LINEARITY *)
+(* ========================================================================= *)
+
+let clinear = new_definition
+ `clinear (f:complex^M->complex^N)
+ <=> (!x y. f(x + y) = f(x) + f(y)) /\ (!c x. f(c % x) = c % f(x))`;;
+
+let COMMON_TAC additional_thms =
+ SIMP_TAC[clinear] THEN REPEAT STRIP_TAC THEN ONCE_REWRITE_TAC[CART_EQ]
+ THEN SIMP_TAC(CVECTOR_ADD_COMPONENT::CVECTOR_MUL_COMPONENT::additional_thms)
+ THEN SIMPLE_COMPLEX_ARITH_TAC;;
+
+let CLINEAR_COMPOSE_CMUL = prove
+ (`!f:complex^M->complex^N c. clinear f ==> clinear (\x. c % f x)`,
+ COMMON_TAC[]);;
+
+let CLINEAR_COMPOSE_NEG = prove
+ (`!f:complex^M->complex^N. clinear f ==> clinear (\x. --(f x))`,
+ COMMON_TAC[CVECTOR_NEG_COMPONENT]);;
+
+let CLINEAR_COMPOSE_ADD = prove
+ (`!f:complex^M->complex^N g. clinear f /\ clinear g ==> clinear (\x. f x + g x)`,
+ COMMON_TAC[]);;
+
+let CLINEAR_COMPOSE_SUB = prove
+ (`!f:complex^M->complex^N g. clinear f /\ clinear g ==> clinear (\x. f x - g x)`,
+ COMMON_TAC[CVECTOR_SUB_COMPONENT]);;
+
+let CLINEAR_COMPOSE = prove
+ (`!f:complex^M->complex^N g. clinear f /\ clinear g ==> clinear (g o f)`,
+ SIMP_TAC[clinear;o_THM]);;
+
+let CLINEAR_ID = prove
+ (`clinear (\x:complex^N. x)`,
+ REWRITE_TAC[clinear]);;
+
+let CLINEAR_I = prove
+ (`clinear (I:complex^N->complex^N)`,
+ REWRITE_TAC[I_DEF;CLINEAR_ID]);;
+
+let CLINEAR_ZERO = prove
+ (`clinear ((\x. cvector_zero):complex^M->complex^N)`,
+ COMMON_TAC[CVECTOR_ZERO_COMPONENT]);;
+
+let CLINEAR_NEGATION = prove
+ (`clinear ((--):complex^N->complex^N)`,
+ COMMON_TAC[CVECTOR_NEG_COMPONENT]);;
+
+let CLINEAR_VMUL_COMPONENT = prove
+ (`!f:complex^M->complex^N v:complex^P k.
+ clinear f /\ 1 <= k /\ k <= dimindex(:N) ==> clinear (\x. (f x)$k % v)`,
+ COMMON_TAC[]);;
+
+let CLINEAR_0 = prove
+ (`!f:complex^M->complex^N. clinear f ==> (f cvector_zero = cvector_zero)`,
+ MESON_TAC[CVECTOR_MUL_LZERO;clinear]);;
+
+let CLINEAR_CMUL = prove
+ (`!f:complex^M->complex^N c x. clinear f ==> (f (c % x) = c % f x)`,
+ SIMP_TAC[clinear]);;
+
+let CLINEAR_NEG = prove
+ (`!f:complex^M->complex^N x. clinear f ==> (f (--x) = --(f x))`,
+ ONCE_REWRITE_TAC[CVECTOR_NEG_MINUS1] THEN SIMP_TAC[CLINEAR_CMUL]);;
+
+let CLINEAR_ADD = prove
+ (`!f:complex^M->complex^N x y. clinear f ==> (f (x + y) = f x + f y)`,
+ SIMP_TAC[clinear]);;
+
+let CLINEAR_SUB = prove
+ (`!f:complex^M->complex^N x y. clinear f ==> (f(x - y) = f x - f y)`,
+ SIMP_TAC[CVECTOR_SUB;CLINEAR_ADD;CLINEAR_NEG]);;
+
+let CLINEAR_INJECTIVE_0 = prove
+ (`!f:complex^M->complex^N.
+ clinear f
+ ==> ((!x y. f x = f y ==> x = y)
+ <=> (!x. f x = cvector_zero ==> x = cvector_zero))`,
+ ONCE_REWRITE_TAC[GSYM CVECTOR_SUB_EQ]
+ THEN SIMP_TAC[CVECTOR_SUB_RZERO;GSYM CLINEAR_SUB]
+ THEN MESON_TAC[CVECTOR_SUB_RZERO]);;
+
+
+
+(* ========================================================================= *)
+(* PASTING COMPLEX VECTORS *)
+(* ========================================================================= *)
+
+let CLINEAR_FSTCART_SNDCART = prove
+ (`clinear fstcart /\ clinear sndcart`,
+ SIMP_TAC[clinear;fstcart;sndcart;CART_EQ;LAMBDA_BETA;CVECTOR_ADD_COMPONENT;
+ CVECTOR_MUL_COMPONENT; DIMINDEX_FINITE_SUM;
+ ARITH_RULE `x <= a ==> x <= a + b:num`;
+ ARITH_RULE `x <= b ==> x + a <= a + b:num`]);;
+
+let FSTCART_CLINEAR = CONJUNCT1 CLINEAR_FSTCART_SNDCART;;
+let SNDCART_CLINEAR = CONJUNCT2 CLINEAR_FSTCART_SNDCART;;
+
+let FSTCART_SNDCART_CVECTOR_ZERO = prove
+ (`fstcart cvector_zero = cvector_zero /\ sndcart cvector_zero = cvector_zero`,
+ SIMP_TAC[CVECTOR_ZERO_COMPONENT;fstcart;sndcart;LAMBDA_BETA;CART_EQ;
+ DIMINDEX_FINITE_SUM;ARITH_RULE `x <= a ==> x <= a + b:num`;
+ ARITH_RULE `x <= b ==> x + a <= a + b:num`]);;
+
+let FSTCART_CVECTOR_ZERO = CONJUNCT1 FSTCART_SNDCART_CVECTOR_ZERO;;
+let SNDCART_CVECTOR_ZERO = CONJUNCT2 FSTCART_SNDCART_CVECTOR_ZERO;;
+
+let FSTCART_SNDCART_CVECTOR_ADD = prove
+ (`!x:complex^(M,N)finite_sum y.
+ fstcart(x + y) = fstcart(x) + fstcart(y)
+ /\ sndcart(x + y) = sndcart(x) + sndcart(y)`,
+ REWRITE_TAC[REWRITE_RULE[clinear] CLINEAR_FSTCART_SNDCART]);;
+
+let FSTCART_SNDCART_CVECTOR_MUL = prove
+ (`!x:complex^(M,N)finite_sum c.
+ fstcart(c % x) = c % fstcart(x) /\ sndcart(c % x) = c % sndcart(x)`,
+ REWRITE_TAC[REWRITE_RULE[clinear] CLINEAR_FSTCART_SNDCART]);;
+
+let PASTECART_TAC xs =
+ REWRITE_TAC(PASTECART_EQ::FSTCART_PASTECART::SNDCART_PASTECART::xs);;
+
+let PASTECART_CVECTOR_ZERO = prove
+ (`pastecart (cvector_zero:complex^N) (cvector_zero:complex^M) = cvector_zero`,
+ PASTECART_TAC[FSTCART_SNDCART_CVECTOR_ZERO]);;
+
+let PASTECART_EQ_CVECTOR_ZERO = prove
+ (`!x:complex^N y:complex^M.
+ pastecart x y = cvector_zero <=> x = cvector_zero /\ y = cvector_zero`,
+ PASTECART_TAC [FSTCART_SNDCART_CVECTOR_ZERO]);;
+
+let PASTECART_CVECTOR_ADD = prove
+ (`!x1 y2 x2:complex^N y2:complex^M.
+ pastecart x1 y1 + pastecart x2 y2 = pastecart (x1 + x2) (y1 + y2)`,
+ PASTECART_TAC [FSTCART_SNDCART_CVECTOR_ADD]);;
+
+let PASTECART_CVECTOR_MUL = prove
+ (`!x1 x2 c:complex.
+ pastecart (c % x1) (c % y1) = c % pastecart x1 y1`, PASTECART_TAC [FSTCART_SNDCART_CVECTOR_MUL]);;
+
+
+(* ========================================================================= *)
+(* REAL AND IMAGINARY VECTORS *)
+(* ========================================================================= *)
+
+let cvector_re = new_definition
+ `cvector_re :complex^N -> real^N = vector_map Re`;;
+let cvector_im = new_definition
+ `cvector_im :complex^N -> real^N = vector_map Im`;;
+let vector_to_cvector = new_definition
+ `vector_to_cvector :real^N -> complex^N = vector_map Cx`;;
+
+let CVECTOR_RE_COMPONENT = prove
+ (`!x:complex^N i. (cvector_re x)$i = Re (x$i)`,
+ REWRITE_TAC[cvector_re;VECTOR_MAP_COMPONENT]);;
+let CVECTOR_IM_COMPONENT = prove
+ (`!x:complex^N i. (cvector_im x)$i = Im (x$i)`,
+ REWRITE_TAC[cvector_im;VECTOR_MAP_COMPONENT]);;
+let VECTOR_TO_CVECTOR_COMPONENT = prove
+ (`!x:real^N i. (vector_to_cvector x)$i = Cx(x$i)`,
+ REWRITE_TAC[vector_to_cvector;VECTOR_MAP_COMPONENT]);;
+
+let VECTOR_TO_CVECTOR_ZERO = prove
+ (`vector_to_cvector (vec 0) = cvector_zero:complex^N`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[VECTOR_TO_CVECTOR_COMPONENT;CVECTOR_ZERO_COMPONENT;
+ VEC_COMPONENT]);;
+
+let VECTOR_TO_CVECTOR_ZERO_EQ = prove
+ (`!x:real^N. vector_to_cvector x = cvector_zero <=> x = vec 0`,
+ GEN_TAC THEN EQ_TAC THEN SIMP_TAC[VECTOR_TO_CVECTOR_ZERO]
+ THEN ONCE_REWRITE_TAC[CART_EQ]
+ THEN SIMP_TAC[VECTOR_TO_CVECTOR_COMPONENT;CVECTOR_ZERO_COMPONENT;
+ VEC_COMPONENT;CX_INJ]);;
+
+let CVECTOR_ZERO_VEC0 = prove
+ (`!x:complex^N. x = cvector_zero <=> cvector_re x = vec 0 /\ cvector_im x = vec 0`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_ZERO_COMPONENT;CVECTOR_RE_COMPONENT;
+ CVECTOR_IM_COMPONENT;VEC_COMPONENT;COMPLEX_EQ;RE_CX;IM_CX]
+ THEN MESON_TAC[]);;
+
+let VECTOR_TO_CVECTOR_MUL = prove
+ (`!a x:real^N. vector_to_cvector (a % x) = Cx a % vector_to_cvector x`,
+ ONCE_REWRITE_TAC[CART_EQ] THEN REWRITE_TAC[VECTOR_TO_CVECTOR_COMPONENT;CVECTOR_MUL_COMPONENT;VECTOR_MUL_COMPONENT;CX_MUL]);;
+
+let CVECTOR_EQ = prove
+ (`!x:complex^N y z.
+ x = vector_to_cvector y + ii % vector_to_cvector z
+ <=> cvector_re x = y /\ cvector_im x = z`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_ADD_COMPONENT;CVECTOR_MUL_COMPONENT;
+ CVECTOR_RE_COMPONENT;CVECTOR_IM_COMPONENT;VECTOR_TO_CVECTOR_COMPONENT]
+ THEN REWRITE_TAC[COMPLEX_EQ;RE_CX;IM_CX;RE_ADD;IM_ADD;RE_MUL_II;REAL_NEG_0;
+ REAL_ADD_RID;REAL_ADD_LID;IM_MUL_II] THEN MESON_TAC[]);;
+
+let CVECTOR_RE_VECTOR_TO_CVECTOR = prove
+ (`!x:real^N. cvector_re (vector_to_cvector x) = x`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_RE_COMPONENT;VECTOR_TO_CVECTOR_COMPONENT;RE_CX]);;
+
+let CVECTOR_IM_VECTOR_TO_CVECTOR = prove
+ (`!x:real^N. cvector_im (vector_to_cvector x) = vec 0`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_IM_COMPONENT;VECTOR_TO_CVECTOR_COMPONENT;IM_CX;
+ VEC_COMPONENT]);;
+
+let CVECTOR_IM_VECTOR_TO_CVECTOR_IM = prove
+ (`!x:real^N. cvector_im (ii % vector_to_cvector x) = x`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_IM_COMPONENT;VECTOR_TO_CVECTOR_COMPONENT;IM_CX;
+ VEC_COMPONENT;CVECTOR_MUL_COMPONENT;IM_MUL_II;RE_CX]);;
+
+let VECTOR_TO_CVECTOR_CVECTOR_RE_IM = prove
+ (`!x:complex^N.
+ vector_to_cvector (cvector_re x) + ii % vector_to_cvector (cvector_im x)
+ = x`,
+ GEN_TAC THEN MATCH_MP_TAC EQ_SYM THEN REWRITE_TAC[CVECTOR_EQ]);;
+
+let CVECTOR_IM_VECTOR_TO_CVECTOR_RE_IM = prove
+ (`!x y:real^N. cvector_im (vector_to_cvector x + ii % vector_to_cvector y) = y`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_IM_COMPONENT;CVECTOR_ADD_COMPONENT;
+ CVECTOR_MUL_COMPONENT;VECTOR_TO_CVECTOR_COMPONENT;IM_ADD;IM_CX;IM_MUL_II;
+ RE_CX;REAL_ADD_LID]);;
+
+let CVECTOR_RE_VECTOR_TO_CVECTOR_RE_IM = prove
+ (`!x y:real^N. cvector_re (vector_to_cvector x + ii % vector_to_cvector y)= x`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_RE_COMPONENT;CVECTOR_ADD_COMPONENT;
+ CVECTOR_MUL_COMPONENT;RE_ADD;VECTOR_TO_CVECTOR_COMPONENT;RE_CX;RE_MUL_CX;
+ RE_II;REAL_MUL_LZERO;REAL_ADD_RID]);;
+
+let CVECTOR_RE_ADD = prove
+ (`!x y:complex^N. cvector_re (x+y) = cvector_re x + cvector_re y`,
+ ONCE_REWRITE_TAC[CART_EQ] THEN REWRITE_TAC[CVECTOR_RE_COMPONENT;
+ VECTOR_ADD_COMPONENT;CVECTOR_ADD_COMPONENT;RE_ADD]);;
+
+let CVECTOR_IM_ADD = prove
+ (`!x y:complex^N. cvector_im (x+y) = cvector_im x + cvector_im y`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_IM_COMPONENT;VECTOR_ADD_COMPONENT;
+ CVECTOR_ADD_COMPONENT;IM_ADD]);;
+
+let CVECTOR_RE_MUL = prove
+ (`!a x:complex^N. cvector_re (Cx a % x) = a % cvector_re x`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_RE_COMPONENT;VECTOR_MUL_COMPONENT;
+ CVECTOR_MUL_COMPONENT;RE_MUL_CX]);;
+
+let CVECTOR_IM_MUL = prove
+ (`!a x:complex^N. cvector_im (Cx a % x) = a % cvector_im x`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_IM_COMPONENT;VECTOR_MUL_COMPONENT;
+ CVECTOR_MUL_COMPONENT;IM_MUL_CX]);;
+
+let CVECTOR_RE_VECTOR_MAP = prove
+ (`!f v:A^N. cvector_re (vector_map f v) = vector_map (Re o f) v`,
+ REWRITE_TAC[cvector_re;VECTOR_MAP_VECTOR_MAP]);;
+
+let CVECTOR_IM_VECTOR_MAP = prove
+ (`!f v:A^N. cvector_im (vector_map f v) = vector_map (Im o f) v`,
+ REWRITE_TAC[cvector_im;VECTOR_MAP_VECTOR_MAP]);;
+
+let VECTOR_MAP_CVECTOR_RE = prove
+ (`!f:real->A v:complex^N.
+ vector_map f (cvector_re v) = vector_map (f o Re) v`,
+ REWRITE_TAC[cvector_re;VECTOR_MAP_VECTOR_MAP]);;
+
+let VECTOR_MAP_CVECTOR_IM = prove
+ (`!f:real->A v:complex^N.
+ vector_map f (cvector_im v) = vector_map (f o Im) v`,
+ REWRITE_TAC[cvector_im;VECTOR_MAP_VECTOR_MAP]);;
+
+let CVECTOR_RE_VECTOR_MAP2 = prove
+ (`!f v1:A^N v2:B^N.
+ cvector_re (vector_map2 f v1 v2) = vector_map2 (\x y. Re (f x y)) v1 v2`,
+ REWRITE_TAC[cvector_re;VECTOR_MAP_VECTOR_MAP2]);;
+
+let CVECTOR_IM_VECTOR_MAP2 = prove
+ (`!f v1:A^N v2:B^N.
+ cvector_im (vector_map2 f v1 v2) = vector_map2 (\x y. Im (f x y)) v1 v2`,
+ REWRITE_TAC[cvector_im;VECTOR_MAP_VECTOR_MAP2]);;
+
+
+(* ========================================================================= *)
+(* FLATTENING COMPLEX VECTORS INTO REAL VECTORS *)
+(* *)
+(* Note: *)
+(* Theoretically, the following could be defined more generally for matrices *)
+(* instead of complex vectors, but this would require a "finite_prod" type *)
+(* constructor, which is not available right now, and which, at first sight, *)
+(* would probably require dependent types. *)
+(* ========================================================================= *)
+
+let cvector_flatten = new_definition
+ `cvector_flatten (v:complex^N) :real^(N,N) finite_sum =
+ pastecart (cvector_re v) (cvector_im v)`;;
+
+let FLATTEN_RE_IM_COMPONENT = prove
+ (`!v:complex^N i.
+ 1 <= i /\ i <= 2 * dimindex(:N)
+ ==> (cvector_flatten v)$i =
+ if i <= dimindex(:N)
+ then (cvector_re v)$i
+ else (cvector_im v)$(i-dimindex(:N))`,
+ SIMP_TAC[MULT_2;GSYM DIMINDEX_FINITE_SUM;cvector_flatten;pastecart;
+ LAMBDA_BETA]);;
+
+let complex_vector = new_definition
+ `complex_vector (v1,v2) :complex^N
+ = vector_map2 (\x y. Cx x + ii * Cx y) v1 v2`;;
+
+let COMPLEX_VECTOR_TRANSPOSE = prove(
+ `!v1 v2:real^N.
+ complex_vector (v1,v2) = vector_to_cvector v1 + ii % vector_to_cvector v2`,
+ ONCE_REWRITE_TAC[CART_EQ]
+ THEN SIMP_TAC[complex_vector;CVECTOR_ADD_COMPONENT;CVECTOR_MUL_COMPONENT;
+ VECTOR_TO_CVECTOR_COMPONENT;VECTOR_MAP2_COMPONENT]);;
+
+let cvector_unflatten = new_definition
+ `cvector_unflatten (v:real^(N,N) finite_sum) :complex^N
+ = complex_vector (fstcart v, sndcart v)`;;
+
+let UNFLATTEN_FLATTEN = prove
+ (`cvector_unflatten o cvector_flatten = I :complex^N -> complex^N`,
+ REWRITE_TAC[FUN_EQ_THM;o_DEF;I_DEF;cvector_flatten;cvector_unflatten;
+ FSTCART_PASTECART;SNDCART_PASTECART;COMPLEX_VECTOR_TRANSPOSE;
+ VECTOR_TO_CVECTOR_CVECTOR_RE_IM]);;
+
+let FLATTEN_UNFLATTEN = prove
+ (`cvector_flatten o cvector_unflatten =
+ I :real^(N,N) finite_sum -> real^(N,N) finite_sum`,
+ REWRITE_TAC[FUN_EQ_THM;o_DEF;I_DEF;cvector_flatten;cvector_unflatten;
+ PASTECART_FST_SND;COMPLEX_VECTOR_TRANSPOSE;
+ CVECTOR_RE_VECTOR_TO_CVECTOR_RE_IM;CVECTOR_IM_VECTOR_TO_CVECTOR_RE_IM]);;
+
+let FLATTEN_CLINEAR = prove
+ (`!f:complex^N->complex^M.
+ clinear f ==> linear (cvector_flatten o f o cvector_unflatten)`,
+ REWRITE_TAC[clinear;linear;cvector_flatten;cvector_unflatten;o_DEF;
+ FSTCART_ADD;SNDCART_ADD;PASTECART_ADD;complex_vector;GSYM PASTECART_CMUL]
+ THEN REPEAT STRIP_TAC THEN REPEAT (AP_TERM_TAC ORELSE MK_COMB_TAC)
+ THEN REWRITE_TAC(map GSYM [CVECTOR_RE_ADD;CVECTOR_IM_ADD;CVECTOR_RE_MUL;
+ CVECTOR_IM_MUL])
+ THEN AP_TERM_TAC THEN ASSUM_LIST (REWRITE_TAC o map GSYM)
+ THEN AP_TERM_TAC THEN ONCE_REWRITE_TAC[CART_EQ]
+ THEN SIMP_TAC[VECTOR_MAP2_COMPONENT;VECTOR_ADD_COMPONENT;
+ CVECTOR_ADD_COMPONENT;CX_ADD;VECTOR_MUL_COMPONENT;CVECTOR_MUL_COMPONENT;
+ FSTCART_CMUL;SNDCART_CMUL;CX_MUL]
+ THEN SIMPLE_COMPLEX_ARITH_TAC);;
+
+let FLATTEN_MAP = prove
+ (`!f g.
+ f = vector_map g
+ ==> !x:complex^N.
+ cvector_flatten (vector_map f x) = vector_map g (cvector_flatten x)`,
+ SIMP_TAC[cvector_flatten;CVECTOR_RE_VECTOR_MAP;CVECTOR_IM_VECTOR_MAP;
+ GSYM PASTECART_VECTOR_MAP;VECTOR_MAP_CVECTOR_RE;VECTOR_MAP_CVECTOR_IM;
+ o_DEF;IM_DEF;RE_DEF;VECTOR_MAP_COMPONENT]);;
+
+let FLATTEN_NEG = prove
+ (`!x:complex^N. cvector_flatten (--x) = --(cvector_flatten x)`,
+ REWRITE_TAC[cvector_neg;MATCH_MP FLATTEN_MAP COMPLEX_NEG_IS_A_MAP]
+ THEN REWRITE_TAC[VECTOR_NEG_IS_A_MAP]);;
+
+let CVECTOR_NEG_ALT = prove
+ (`!x:complex^N. --x = cvector_unflatten (--(cvector_flatten x))`,
+ REWRITE_TAC[GSYM FLATTEN_NEG;
+ REWRITE_RULE[o_DEF;FUN_EQ_THM;I_DEF] UNFLATTEN_FLATTEN]);;
+
+let FLATTEN_MAP2 = prove(
+ `!f g.
+ f = vector_map2 g ==>
+ !x y:complex^N.
+ cvector_flatten (vector_map2 f x y) =
+ vector_map2 g (cvector_flatten x) (cvector_flatten y)`,
+ SIMP_TAC[cvector_flatten;CVECTOR_RE_VECTOR_MAP2;CVECTOR_IM_VECTOR_MAP2;
+ CVECTOR_RE_VECTOR_MAP2;GSYM PASTECART_VECTOR_MAP2]
+ THEN REWRITE_TAC[cvector_re;cvector_im;VECTOR_MAP2_LVECTOR_MAP;
+ VECTOR_MAP2_RVECTOR_MAP]
+ THEN REPEAT MK_COMB_TAC
+ THEN REWRITE_TAC[FUN_EQ_THM;IM_DEF;RE_DEF;VECTOR_MAP2_COMPONENT;o_DEF]);;
+
+let FLATTEN_ADD = prove
+ (`!x y:complex^N.
+ cvector_flatten (x+y) = cvector_flatten x + cvector_flatten y`,
+ REWRITE_TAC[cvector_add;MATCH_MP FLATTEN_MAP2 COMPLEX_ADD_IS_A_MAP]
+ THEN REWRITE_TAC[VECTOR_ADD_IS_A_MAP]);;
+
+let CVECTOR_ADD_ALT = prove
+ (`!x y:complex^N.
+ x+y = cvector_unflatten (cvector_flatten x + cvector_flatten y)`,
+ REWRITE_TAC[GSYM FLATTEN_ADD;
+ REWRITE_RULE[o_DEF;FUN_EQ_THM;I_DEF] UNFLATTEN_FLATTEN]);;
+
+let FLATTEN_SUB = prove
+ (`!x y:complex^N. cvector_flatten (x-y) = cvector_flatten x - cvector_flatten y`,
+ REWRITE_TAC[cvector_sub;MATCH_MP FLATTEN_MAP2 COMPLEX_SUB_IS_A_MAP]
+ THEN REWRITE_TAC[VECTOR_SUB_IS_A_MAP]);;
+
+let CVECTOR_SUB_ALT = prove
+ (`!x y:complex^N. x-y = cvector_unflatten (cvector_flatten x - cvector_flatten y)`,
+ REWRITE_TAC[GSYM FLATTEN_SUB;
+ REWRITE_RULE[o_DEF;FUN_EQ_THM;I_DEF] UNFLATTEN_FLATTEN]);;
+
+
+(* ========================================================================= *)
+(* CONJUGATE VECTOR. *)
+(* ========================================================================= *)
+
+let cvector_cnj = new_definition
+ `cvector_cnj : complex^N->complex^N = vector_map cnj`;;
+
+let CVECTOR_MAP_PROPERTY thm =
+ VECTOR_MAP_PROPERTY thm [cvector_zero;cvector_add;cvector_sub;cvector_neg;
+ cvector_mul;cvector_cnj;cvector_re;cvector_im];;
+
+let CVECTOR_CNJ_ADD = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. cvector_cnj (x+y) = cvector_cnj x + cvector_cnj y`
+ CNJ_ADD;;
+
+let CVECTOR_CNJ_SUB = CVECTOR_MAP_PROPERTY
+ `!x y:complex^N. cvector_cnj (x-y) = cvector_cnj x - cvector_cnj y`
+ CNJ_SUB;;
+
+let CVECTOR_CNJ_NEG = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. cvector_cnj (--x) = --(cvector_cnj x)`
+ CNJ_NEG;;
+
+let CVECTOR_RE_CNJ = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. cvector_re (cvector_cnj x) = cvector_re x`
+ RE_CNJ;;
+
+let CVECTOR_IM_CNJ = prove
+ (`!x:complex^N. cvector_im (cvector_cnj x) = --(cvector_im x)`,
+ VECTOR_MAP_PROPERTY_TAC[cvector_im;cvector_cnj;VECTOR_NEG_IS_A_MAP] IM_CNJ);;
+
+let CVECTOR_CNJ_CNJ = CVECTOR_MAP_PROPERTY
+ `!x:complex^N. cvector_cnj (cvector_cnj x) = x`
+ CNJ_CNJ;;
+
+
+(* ========================================================================= *)
+(* CROSS PRODUCTS IN COMPLEX^3. *)
+(* ========================================================================= *)
+
+prioritize_vector();;
+
+parse_as_infix("ccross",(20,"right"));;
+
+let ccross = new_definition
+ `((ccross):complex^3 -> complex^3 -> complex^3) x y = vector [
+ x$2 * y$3 - x$3 * y$2;
+ x$3 * y$1 - x$1 * y$3;
+ x$1 * y$2 - x$2 * y$1
+ ]`;;
+
+let CCROSS_COMPONENT = prove
+ (`!x y:complex^3.
+ (x ccross y)$1 = x$2 * y$3 - x$3 * y$2
+ /\ (x ccross y)$2 = x$3 * y$1 - x$1 * y$3
+ /\ (x ccross y)$3 = x$1 * y$2 - x$2 * y$1`,
+ REWRITE_TAC[ccross;VECTOR_3]);;
+
+(* simple handy instantiation of CART_EQ for dimension 3*)
+let CART_EQ3 = prove
+ (`!x y:complex^3. x = y <=> x$1 = y$1 /\ x$2 = y$2 /\ x$3 = y$3`,
+ GEN_REWRITE_TAC (PATH_CONV "rbrblr") [CART_EQ]
+ THEN REWRITE_TAC[DIMINDEX_3;FORALL_3]);;
+
+let CCROSS_TAC lemmas =
+ REWRITE_TAC(CART_EQ3::CCROSS_COMPONENT::lemmas)
+ THEN SIMPLE_COMPLEX_ARITH_TAC;;
+
+let CCROSS_LZERO = prove
+ (`!x:complex^3. cvector_zero ccross x = cvector_zero`,
+ CCROSS_TAC[CVECTOR_ZERO_COMPONENT]);;
+
+let CCROSS_RZERO = prove
+ (`!x:complex^3. x ccross cvector_zero = cvector_zero`,
+ CCROSS_TAC[CVECTOR_ZERO_COMPONENT]);;
+
+let CCROSS_SKEW = prove
+ (`!x y:complex^3. (x ccross y) = --(y ccross x)`,
+ CCROSS_TAC[CVECTOR_NEG_COMPONENT]);;
+
+let CCROSS_REFL = prove
+ (`!x:complex^3. x ccross x = cvector_zero`,
+ CCROSS_TAC[CVECTOR_ZERO_COMPONENT]);;
+
+let CCROSS_LADD = prove
+ (`!x y z:complex^3. (x + y) ccross z = (x ccross z) + (y ccross z)`,
+ CCROSS_TAC[CVECTOR_ADD_COMPONENT]);;
+
+let CCROSS_RADD = prove
+ (`!x y z:complex^3. x ccross(y + z) = (x ccross y) + (x ccross z)`,
+ CCROSS_TAC[CVECTOR_ADD_COMPONENT]);;
+
+let CCROSS_LMUL = prove
+ (`!c x y:complex^3. (c % x) ccross y = c % (x ccross y)`,
+ CCROSS_TAC[CVECTOR_MUL_COMPONENT]);;
+
+let CCROSS_RMUL = prove
+ (`!c x y:complex^3. x ccross (c % y) = c % (x ccross y)`,
+ CCROSS_TAC[CVECTOR_MUL_COMPONENT]);;
+
+let CCROSS_LNEG = prove
+ (`!x y:complex^3. (--x) ccross y = --(x ccross y)`,
+ CCROSS_TAC[CVECTOR_NEG_COMPONENT]);;
+
+let CCROSS_RNEG = prove
+ (`!x y:complex^3. x ccross (--y) = --(x ccross y)`,
+ CCROSS_TAC[CVECTOR_NEG_COMPONENT]);;
+
+let CCROSS_JACOBI = prove
+ (`!(x:complex^3) y z.
+ x ccross (y ccross z) + y ccross (z ccross x) + z ccross (x ccross y) =
+ cvector_zero`,
+ REWRITE_TAC[CART_EQ3]
+ THEN REWRITE_TAC[CVECTOR_ADD_COMPONENT;CCROSS_COMPONENT;
+ CVECTOR_ZERO_COMPONENT] THEN SIMPLE_COMPLEX_ARITH_TAC);;
+
+
+(* ========================================================================= *)
+(* DOT PRODUCTS IN COMPLEX^N *)
+(* *)
+(* Only difference with the real case: *)
+(* we take the conjugate of the 2nd argument *)
+(* ========================================================================= *)
+
+prioritize_complex();;
+
+parse_as_infix("cdot",(20,"right"));;
+
+let cdot = new_definition
+ `(cdot) (x:complex^N) (y:complex^N) =
+ vsum (1..dimindex(:N)) (\i. x$i * cnj(y$i))`;;
+
+(* The dot product is symmetric MODULO the conjugate *)
+let CDOT_SYM = prove
+ (`!x:complex^N y. x cdot y = cnj (y cdot x)`,
+ REWRITE_TAC[cdot]
+ THEN REWRITE_TAC[MATCH_MP (SPEC_ALL CNJ_VSUM) (SPEC `dimindex (:N)` (GEN_ALL
+ (CONJUNCT1 (SPEC_ALL (REWRITE_RULE[HAS_SIZE] HAS_SIZE_NUMSEG_1)))))]
+ THEN REWRITE_TAC[CNJ_MUL;CNJ_CNJ;COMPLEX_MUL_SYM]);;
+
+let REAL_CDOT_SELF = prove
+ (`!x:complex^N. real(x cdot x)`,
+ REWRITE_TAC[REAL_CNJ;GSYM CDOT_SYM]);;
+
+(* The following theorems are usual axioms of the hermitian dot product, they are proved later on.
+ * let CDOT_SELF_POS = prove(`!x:complex^N. &0 <= real_of_complex (x cdot x)`, ...
+ * let CDOT_EQ_0 = prove(`!x:complex^N. x cdot x = Cx(&0) <=> x = cvector_zero`
+ *)
+
+let CDOT_LADD = prove
+ (`!x:complex^N y z. (x + y) cdot z = (x cdot z) + (y cdot z)`,
+ REWRITE_TAC[cdot]
+ THEN REWRITE_TAC[MATCH_MP (GSYM VSUM_ADD) (SPEC `dimindex (:N)` (GEN_ALL
+ (CONJUNCT1 (SPEC_ALL (REWRITE_RULE[HAS_SIZE] HAS_SIZE_NUMSEG_1)))))]
+ THEN REPEAT GEN_TAC THEN MATCH_MP_TAC VSUM_EQ THEN GEN_TAC THEN DISCH_TAC
+ THEN REWRITE_TAC[FUN_EQ_THM]
+ THEN REWRITE_TAC[SPECL [`(x:real^2^N)$(x':num)`;`(y:real^2^N)$(x':num)`;
+ `cnj ((z:real^2^N)$(x':num))`] (GSYM COMPLEX_ADD_RDISTRIB)]
+ THEN REWRITE_TAC[CVECTOR_ADD_COMPONENT]);;
+
+let CDOT_RADD = prove
+ (`!x:complex^N y z. x cdot (y + z) = (x cdot y) + (x cdot z)`,
+ REWRITE_TAC[cdot]
+ THEN REWRITE_TAC[MATCH_MP (GSYM VSUM_ADD) (SPEC `dimindex (:N)` (GEN_ALL
+ (CONJUNCT1 (SPEC_ALL (REWRITE_RULE[HAS_SIZE] HAS_SIZE_NUMSEG_1)))))]
+ THEN REPEAT GEN_TAC THEN MATCH_MP_TAC VSUM_EQ THEN GEN_TAC THEN DISCH_TAC
+ THEN REWRITE_TAC[FUN_EQ_THM]
+ THEN REWRITE_TAC[SPECL [`(x:real^2^N)$(x':num)`; `cnj((y:real^2^N)$(x':num))`;
+ `cnj ((z:real^2^N)$(x':num))`] (GSYM COMPLEX_ADD_LDISTRIB)]
+ THEN REWRITE_TAC[CNJ_ADD; CVECTOR_ADD_COMPONENT]);;
+
+let CDOT_LSUB = prove
+ (`!x:complex^N y z. (x - y) cdot z = (x cdot z) - (y cdot z)`,
+ REWRITE_TAC[cdot]
+ THEN REWRITE_TAC[MATCH_MP (GSYM VSUM_SUB) (SPEC `dimindex (:N)` (GEN_ALL
+ (CONJUNCT1 (SPEC_ALL (REWRITE_RULE[HAS_SIZE] HAS_SIZE_NUMSEG_1)))))]
+ THEN REPEAT GEN_TAC THEN MATCH_MP_TAC VSUM_EQ THEN GEN_TAC THEN DISCH_TAC
+ THEN REWRITE_TAC[FUN_EQ_THM]
+ THEN REWRITE_TAC[SPECL [`(x:real^2^N)$(x':num)`; `(y:real^2^N)$(x':num)`;
+ `cnj ((z:real^2^N)$(x':num))`] (GSYM COMPLEX_SUB_RDISTRIB)]
+ THEN REWRITE_TAC[CVECTOR_SUB_COMPONENT]);;
+
+let CDOT_RSUB = prove
+ (`!x:complex^N y z. x cdot (y - z) = (x cdot y) - (x cdot z)`,
+ REWRITE_TAC[cdot]
+ THEN REWRITE_TAC[MATCH_MP (GSYM VSUM_SUB) (SPEC `dimindex (:N)` (GEN_ALL
+ (CONJUNCT1 (SPEC_ALL (REWRITE_RULE[HAS_SIZE] HAS_SIZE_NUMSEG_1)))))]
+ THEN REPEAT GEN_TAC THEN MATCH_MP_TAC VSUM_EQ THEN GEN_TAC THEN DISCH_TAC
+ THEN REWRITE_TAC[FUN_EQ_THM]
+ THEN REWRITE_TAC[SPECL [`(x:real^2^N)$(x':num)`; `cnj((y:real^2^N)$(x':num))`;
+ `cnj ((z:real^2^N)$(x':num))`] (GSYM COMPLEX_SUB_LDISTRIB)]
+ THEN REWRITE_TAC[CNJ_SUB; CVECTOR_SUB_COMPONENT]);;
+
+let CDOT_LMUL = prove
+ (`!c:complex x:complex^N y. (c % x) cdot y = c * (x cdot y)`,
+ REWRITE_TAC[cdot]
+ THEN REWRITE_TAC[MATCH_MP (GSYM VSUM_COMPLEX_LMUL) (SPEC `dimindex (:N)`
+ (GEN_ALL (CONJUNCT1 (SPEC_ALL (REWRITE_RULE[HAS_SIZE]
+ HAS_SIZE_NUMSEG_1)))))]
+ THEN REWRITE_TAC[CVECTOR_MUL_COMPONENT; GSYM COMPLEX_MUL_ASSOC]);;
+
+let CDOT_RMUL = prove
+ (`!c:complex x:complex^N x y. x cdot (c % y) = cnj c * (x cdot y)`,
+ REWRITE_TAC[cdot]
+ THEN REWRITE_TAC[MATCH_MP (GSYM VSUM_COMPLEX_LMUL) (SPEC `dimindex (:N)`
+ (GEN_ALL (CONJUNCT1 (SPEC_ALL (REWRITE_RULE[HAS_SIZE]
+ HAS_SIZE_NUMSEG_1)))))]
+ THEN REWRITE_TAC[CVECTOR_MUL_COMPONENT; CNJ_MUL; COMPLEX_MUL_AC]);;
+
+let CDOT_LNEG = prove
+ (`!x:complex^N y. (--x) cdot y = --(x cdot y)`,
+ REWRITE_TAC[cdot] THEN ONCE_REWRITE_TAC[COMPLEX_NEG_MINUS1]
+ THEN REWRITE_TAC[MATCH_MP (GSYM VSUM_COMPLEX_LMUL) (SPEC `dimindex (:N)`
+ (GEN_ALL (CONJUNCT1 (SPEC_ALL (REWRITE_RULE[HAS_SIZE]
+ HAS_SIZE_NUMSEG_1)))))]
+ THEN REWRITE_TAC[CVECTOR_NEG_COMPONENT] THEN ONCE_REWRITE_TAC[GSYM
+ COMPLEX_NEG_MINUS1] THEN REWRITE_TAC[COMPLEX_NEG_LMUL]);;
+
+let CDOT_RNEG = prove
+ (`!x:complex^N y. x cdot (--y) = --(x cdot y)`,
+ REWRITE_TAC[cdot] THEN ONCE_REWRITE_TAC[COMPLEX_NEG_MINUS1]
+ THEN REWRITE_TAC[MATCH_MP (GSYM VSUM_COMPLEX_LMUL) (SPEC `dimindex (:N)`
+ (GEN_ALL (CONJUNCT1 (SPEC_ALL (REWRITE_RULE[HAS_SIZE]
+ HAS_SIZE_NUMSEG_1)))))]
+ THEN ONCE_REWRITE_TAC[GSYM COMPLEX_NEG_MINUS1]
+ THEN REWRITE_TAC[CVECTOR_NEG_COMPONENT; CNJ_NEG; COMPLEX_NEG_RMUL]);;
+
+let CDOT_LZERO = prove
+ (`!x:complex^N. cvector_zero cdot x = Cx (&0)`,
+ REWRITE_TAC[cdot] THEN REWRITE_TAC[CVECTOR_ZERO_COMPONENT]
+ THEN REWRITE_TAC[COMPLEX_MUL_LZERO; GSYM COMPLEX_VEC_0; VSUM_0]);;
+
+let CNJ_ZERO = prove(
+ `cnj (Cx(&0)) = Cx(&0)`,
+ REWRITE_TAC[cnj;RE_CX;IM_CX;CX_DEF;REAL_NEG_0]);;
+
+let CDOT_RZERO = prove(
+ `!x:complex^N. x cdot cvector_zero = Cx (&0)`,
+ REWRITE_TAC[cdot] THEN REWRITE_TAC[CVECTOR_ZERO_COMPONENT]
+ THEN REWRITE_TAC[CNJ_ZERO]
+ THEN REWRITE_TAC[COMPLEX_MUL_RZERO;GSYM COMPLEX_VEC_0;VSUM_0]);;
+
+(* Cauchy Schwarz inequality: proved later on
+ * let CDOT_CAUCHY_SCHWARZ = prove (`!x y:complex^N. norm (x cdot y) pow 2 <= cnorm2 x * cnorm2 y`,
+ * let CDOT_CAUCHY_SCHWARZ_EQUAL = prove(`!x y:complex^N. norm (x cdot y) pow 2 = cnorm2 x * cnorm2 y <=> collinear_cvectors x y`,
+*)
+
+let CDOT3 = prove
+ (`!x y:complex^3.
+ x cdot y = (x$1 * cnj (y$1) + x$2 * cnj (y$2) + x$3 * cnj (y$3))`,
+ REWRITE_TAC[cdot] THEN SIMP_TAC [DIMINDEX_3] THEN REWRITE_TAC[VSUM_3]);;
+
+let ADD_CDOT_SYM = prove(
+ `!x y:complex^N. x cdot y + y cdot x = Cx(&2 * Re(x cdot y))`,
+ MESON_TAC[CDOT_SYM;COMPLEX_ADD_CNJ]);;
+
+
+(* ========================================================================= *)
+(* RELATION WITH REAL DOT AND CROSS PRODUCTS *)
+(* ========================================================================= *)
+
+let CCROSS_LREAL = prove
+ (`!r c.
+ (vector_to_cvector r) ccross c =
+ vector_to_cvector (r cross (cvector_re c))
+ + ii % (vector_to_cvector (r cross (cvector_im c)))`,
+ REWRITE_TAC[CART_EQ3;CVECTOR_ADD_COMPONENT;CVECTOR_MUL_COMPONENT;
+ VECTOR_TO_CVECTOR_COMPONENT;CCROSS_COMPONENT;CROSS_COMPONENTS;
+ CVECTOR_RE_COMPONENT;CVECTOR_IM_COMPONENT;complex_mul;RE_CX;IM_CX;CX_DEF;
+ complex_sub;complex_neg;complex_add;RE;IM;RE_II;IM_II]
+ THEN REPEAT STRIP_TAC THEN AP_TERM_TAC THEN REWRITE_TAC[PAIR_EQ]
+ THEN ARITH_TAC);;
+
+let CCROSS_RREAL = prove
+ (`!r c.
+ c ccross (vector_to_cvector r) =
+ vector_to_cvector ((cvector_re c) cross r)
+ + ii % (vector_to_cvector ((cvector_im c) cross r))`,
+ REWRITE_TAC[CART_EQ3;CVECTOR_ADD_COMPONENT;CVECTOR_MUL_COMPONENT;
+ VECTOR_TO_CVECTOR_COMPONENT;CCROSS_COMPONENT;CROSS_COMPONENTS;
+ CVECTOR_RE_COMPONENT;CVECTOR_IM_COMPONENT;complex_mul;RE_CX;IM_CX;CX_DEF;
+ complex_sub;complex_neg;complex_add;RE;IM;RE_II;IM_II]
+ THEN REPEAT STRIP_TAC THEN AP_TERM_TAC THEN REWRITE_TAC[PAIR_EQ]
+ THEN ARITH_TAC);;
+
+let CDOT_LREAL = prove
+ (`!r c.
+ (vector_to_cvector r) cdot c =
+ Cx (r dot (cvector_re c)) - ii * Cx (r dot (cvector_im c))`,
+ REWRITE_TAC[cdot; dot; VECTOR_TO_CVECTOR_COMPONENT;CVECTOR_RE_COMPONENT;
+ CVECTOR_IM_COMPONENT] THEN REPEAT GEN_TAC
+ THEN GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV) [COMPLEX_EXPAND]
+ THEN REWRITE_TAC[MATCH_MP RE_VSUM (SPEC `dimindex (:N)` (GEN_ALL (CONJUNCT1
+ (SPEC_ALL (REWRITE_RULE[HAS_SIZE] HAS_SIZE_NUMSEG_1)))))]
+ THEN REWRITE_TAC[MATCH_MP (IM_VSUM) (SPEC `dimindex (:N)` (GEN_ALL
+ (CONJUNCT1 (SPEC_ALL (REWRITE_RULE[HAS_SIZE]
+ HAS_SIZE_NUMSEG_1)))))]
+ THEN REWRITE_TAC[RE_MUL_CX;RE_CNJ;IM_MUL_CX;IM_CNJ]
+ THEN REWRITE_TAC[COMPLEX_POLY_NEG_CLAUSES] THEN REWRITE_TAC[COMPLEX_MUL_AC]
+ THEN REWRITE_TAC[COMPLEX_MUL_ASSOC] THEN REWRITE_TAC[GSYM CX_MUL]
+ THEN REWRITE_TAC[GSYM SUM_LMUL]
+ THEN REWRITE_TAC[GSYM REAL_NEG_MINUS1;GSYM REAL_MUL_RNEG]);;
+
+let CDOT_RREAL = prove
+ (`!r c.
+ c cdot (vector_to_cvector r) =
+ Cx ((cvector_re c) dot r) + ii * Cx ((cvector_im c) dot r)`,
+ REWRITE_TAC[cdot; dot; VECTOR_TO_CVECTOR_COMPONENT;CVECTOR_RE_COMPONENT;
+ CVECTOR_IM_COMPONENT]
+ THEN REPEAT GEN_TAC
+ THEN GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV) [COMPLEX_EXPAND]
+ THEN REWRITE_TAC[MATCH_MP RE_VSUM (SPEC `dimindex (:N)` (GEN_ALL (CONJUNCT1
+ (SPEC_ALL (REWRITE_RULE[HAS_SIZE] HAS_SIZE_NUMSEG_1)))))]
+ THEN REWRITE_TAC[MATCH_MP IM_VSUM (SPEC `dimindex (:N)` (GEN_ALL (CONJUNCT1
+ (SPEC_ALL (REWRITE_RULE[HAS_SIZE] HAS_SIZE_NUMSEG_1)))))]
+ THEN REWRITE_TAC[CNJ_CX]
+ THEN REWRITE_TAC[RE_MUL_CX;RE_CNJ;IM_MUL_CX;IM_CNJ]);;
+
+
+(* ========================================================================= *)
+(* NORM, UNIT VECTORS. *)
+(* ========================================================================= *)
+
+let cnorm2 = new_definition
+ `cnorm2 (v:complex^N) = real_of_complex (v cdot v)`;;
+
+let CX_CNORM2 = prove
+ (`!v:complex^N. Cx(cnorm2 v) = v cdot v`,
+ SIMP_TAC[cnorm2;REAL_CDOT_SELF;REAL_OF_COMPLEX]);;
+
+let CNORM2_CVECTOR_ZERO = prove
+ (`cnorm2 (cvector_zero:complex^N) = &0`,
+ REWRITE_TAC[cnorm2;CDOT_RZERO;REAL_OF_COMPLEX_CX]);;
+
+let CNORM2_MODULUS = prove
+ (`!x:complex^N. cnorm2 x = (vector_map norm x) dot (vector_map norm x)`,
+ REWRITE_TAC[cnorm2;cdot;COMPLEX_MUL_CNJ;COMPLEX_POW_2;GSYM CX_MUL;
+ VSUM_CX_NUMSEG;dot;VECTOR_MAP_COMPONENT;REAL_OF_COMPLEX_CX]);;
+
+let CNORM2_EQ_0 = prove
+ (`!x:complex^N. cnorm2 x = &0 <=> x = cvector_zero`,
+ REWRITE_TAC[CNORM2_MODULUS;CX_INJ;DOT_EQ_0] THEN GEN_TAC
+ THEN GEN_REWRITE_TAC (RATOR_CONV o DEPTH_CONV) [CART_EQ]
+ THEN REWRITE_TAC[VEC_COMPONENT;VECTOR_MAP_COMPONENT;COMPLEX_NORM_ZERO]
+ THEN GEN_REWRITE_TAC (RAND_CONV o DEPTH_CONV) [CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_ZERO_COMPONENT]);;
+
+let CDOT_EQ_0 = prove
+ (`!x:complex^N. x cdot x = Cx(&0) <=> x = cvector_zero`,
+ SIMP_TAC[TAUT `(p<=>q) <=> ((p==>q) /\ (q==>p))`;CDOT_LZERO]
+ THEN GEN_TAC THEN DISCH_THEN (MP_TAC o MATCH_MP (MESON[REAL_OF_COMPLEX_CX]
+ `x = Cx y ==> real_of_complex x = y`))
+ THEN REWRITE_TAC[GSYM cnorm2;CNORM2_EQ_0]);;
+
+let CNORM2_POS = prove
+ (`!x:complex^N. &0 <= cnorm2 x`, REWRITE_TAC[CNORM2_MODULUS;DOT_POS_LE]);;
+
+let CDOT_SELF_POS = prove
+ (`!x:complex^N. &0 <= real_of_complex (x cdot x)`,
+ REWRITE_TAC[GSYM cnorm2;CNORM2_POS]);;
+
+let CNORM2_MUL = prove
+ (`!a x:complex^N. cnorm2 (a % x) = (norm a) pow 2 * cnorm2 x`,
+ SIMP_TAC[cnorm2;CDOT_LMUL;CDOT_RMUL;
+ SIMPLE_COMPLEX_ARITH `x * cnj x * y = (x * cnj x) * y`;COMPLEX_MUL_CNJ;
+ REAL_OF_COMPLEX_CX;REAL_OF_COMPLEX_MUL;REAL_CX;REAL_CDOT_SELF;
+ GSYM CX_POW]);;
+
+let CNORM2_NORM2_2 = prove
+ (`!x y:real^N.
+ cnorm2 (vector_to_cvector x + ii % vector_to_cvector y) =
+ norm x pow 2 + norm y pow 2`,
+ REWRITE_TAC[cnorm2;vector_norm;cdot;CVECTOR_ADD_COMPONENT;
+ CVECTOR_MUL_COMPONENT;VECTOR_TO_CVECTOR_COMPONENT;CNJ_ADD;CNJ_CX;CNJ_MUL;
+ CNJ_II;COMPLEX_ADD_RDISTRIB;COMPLEX_ADD_LDISTRIB;
+ SIMPLE_COMPLEX_ARITH
+ `(x*x+x*(--ii)*y)+(ii*y)*x+(ii*y)*(--ii)*y = x*x-(ii*ii)*y*y`]
+ THEN REWRITE_TAC[GSYM COMPLEX_POW_2;COMPLEX_POW_II_2;
+ SIMPLE_COMPLEX_ARITH `x-(--Cx(&1))*y = x+y`]
+ THEN SIMP_TAC[MESON[CARD_NUMSEG_1;HAS_SIZE_NUMSEG_1;FINITE_HAS_SIZE]
+ `FINITE (1..dimindex(:N))`;VSUM_ADD;GSYM CX_POW;VSUM_CX;GSYM dot;
+ REAL_POW_2;GSYM dot]
+ THEN SIMP_TAC[GSYM CX_ADD;REAL_OF_COMPLEX_CX;GSYM REAL_POW_2;DOT_POS_LE;
+ SQRT_POW_2]);;
+
+let CNORM2_NORM2 = prove
+ (`!v:complex^N.
+ cnorm2 v = norm (cvector_re v) pow 2 + norm (cvector_im v) pow 2`,
+ GEN_TAC THEN GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV) [GSYM
+ VECTOR_TO_CVECTOR_CVECTOR_RE_IM] THEN REWRITE_TAC[CNORM2_NORM2_2]);;
+
+let CNORM2_ALT = prove
+ (`!x:complex^N. cnorm2 x = norm (x cdot x)`,
+ SIMP_TAC[cnorm2;REAL_OF_COMPLEX_NORM;REAL_CDOT_SELF;EQ_SYM_EQ;REAL_ABS_REFL;
+ REWRITE_RULE[cnorm2] CNORM2_POS]);;
+
+let CNORM2_SUB = prove
+ (`!x y:complex^N. cnorm2 (x-y) = cnorm2 (y-x)`,
+ REWRITE_TAC[cnorm2;CDOT_LSUB;CDOT_RSUB] THEN REPEAT GEN_TAC THEN AP_TERM_TAC
+ THEN SIMPLE_COMPLEX_ARITH_TAC);;
+
+let CNORM2_VECTOR_TO_CVECTOR = prove
+ (`!x:real^N. cnorm2 (vector_to_cvector x) = norm x pow 2`,
+ REWRITE_TAC[CNORM2_ALT;CDOT_RREAL;CVECTOR_RE_VECTOR_TO_CVECTOR;
+ CVECTOR_IM_VECTOR_TO_CVECTOR;DOT_LZERO;COMPLEX_MUL_RZERO;COMPLEX_ADD_RID;
+ DOT_SQUARE_NORM;CX_POW;COMPLEX_NORM_POW;COMPLEX_NORM_CX;REAL_POW2_ABS]);;
+
+let cnorm = new_definition
+ `cnorm :complex^N->real = sqrt o cnorm2`;;
+
+overload_interface ("norm",`cnorm:complex^N->real`);;
+
+let CNORM_CVECTOR_ZERO = prove
+ (`norm (cvector_zero:complex^N) = &0`,
+ REWRITE_TAC[cnorm;CNORM2_CVECTOR_ZERO;o_DEF;SQRT_0]);;
+
+let CNORM_POW_2 = prove
+ (`!x:complex^N. norm x pow 2 = cnorm2 x`,
+ SIMP_TAC[cnorm;o_DEF;SQRT_POW_2;CNORM2_POS]);;
+
+let CNORM_NORM_2 = prove
+ (`!x y:real^N.
+ norm (vector_to_cvector x + ii % vector_to_cvector y) =
+ sqrt(norm x pow 2 + norm y pow 2)`,
+ REWRITE_TAC[cnorm;o_DEF;CNORM2_NORM2_2]);;
+
+let CNORM_NORM = prove(
+ `!v:complex^N.
+ norm v = sqrt(norm (cvector_re v) pow 2 + norm (cvector_im v) pow 2)`,
+ GEN_TAC THEN GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV) [GSYM
+ VECTOR_TO_CVECTOR_CVECTOR_RE_IM] THEN REWRITE_TAC[CNORM_NORM_2]);;
+
+let CNORM_MUL = prove
+ (`!a x:complex^N. norm (a % x) = norm a * norm x`,
+ SIMP_TAC[cnorm;o_DEF;CNORM2_MUL;REAL_LE_POW_2;SQRT_MUL;CNORM2_POS;
+ NORM_POS_LE;POW_2_SQRT]);;
+
+let CNORM_EQ_0 = prove
+ (`!x:complex^N. norm x = &0 <=> x = cvector_zero`,
+ SIMP_TAC[cnorm;o_DEF;SQRT_EQ_0;CNORM2_POS;CNORM2_EQ_0]);;
+
+let CNORM_POS = prove
+ (`!x:complex^N. &0 <= norm x`,
+ SIMP_TAC[cnorm;o_DEF;SQRT_POS_LE;CNORM2_POS]);;
+
+let CNORM_SUB = prove
+ (`!x y:complex^N. norm (x-y) = norm (y-x)`,
+ REWRITE_TAC[cnorm;o_DEF;CNORM2_SUB]);;
+
+let CNORM_VECTOR_TO_CVECTOR = prove
+ (`!x:real^N. norm (vector_to_cvector x) = norm x`,
+ SIMP_TAC[cnorm;o_DEF;CNORM2_VECTOR_TO_CVECTOR;POW_2_SQRT;NORM_POS_LE]);;
+
+let CNORM_BASIS = prove
+ (`!k. 1 <= k /\ k <= dimindex(:N)
+ ==> norm (vector_to_cvector (basis k :real^N)) = &1`,
+ SIMP_TAC[NORM_BASIS;CNORM_VECTOR_TO_CVECTOR]);;
+
+let CNORM_BASIS_1 = prove
+ (`norm(basis 1:real^N) = &1`,
+ SIMP_TAC[NORM_BASIS_1;CNORM_VECTOR_TO_CVECTOR]);;
+
+let CVECTOR_CHOOSE_SIZE = prove
+ (`!c. &0 <= c ==> ?x:complex^N. norm(x) = c`,
+ MESON_TAC[VECTOR_CHOOSE_SIZE;CNORM_VECTOR_TO_CVECTOR]);;
+
+(* Triangle inequality. Proved later on using Cauchy Schwarz inequality.
+ * let CNORM_TRIANGLE = prove(`!x y:complex^N. norm (x+y) <= norm x + norm y`, ...
+ *)
+
+let cunit = new_definition
+ `cunit (X:complex^N) = inv(Cx(norm X))% X`;;
+
+let CUNIT_CVECTOR_ZERO = prove
+ (`cunit cvector_zero = cvector_zero:complex^N`,
+ REWRITE_TAC[cunit;CNORM_CVECTOR_ZERO;COMPLEX_INV_0;CVECTOR_MUL_LZERO]);;
+
+let CDOT_CUNIT_MUL_CUNIT = prove
+ (`!x:complex^N. (cunit x cdot x) % cunit x = x`,
+ GEN_TAC THEN ASM_CASES_TAC `x = cvector_zero:complex^N`
+ THEN ASM_REWRITE_TAC[CUNIT_CVECTOR_ZERO;CDOT_LZERO;CVECTOR_MUL_LZERO]
+ THEN SIMP_TAC[cunit;CVECTOR_MUL_ASSOC;CDOT_LMUL;
+ SIMPLE_COMPLEX_ARITH `(x*y)*x=(x*x)*y`;GSYM COMPLEX_INV_MUL;GSYM CX_MUL;
+ GSYM REAL_POW_2;cnorm;o_DEF;CNORM2_POS;SQRT_POW_2]
+ THEN ASM_SIMP_TAC[cnorm2;REAL_OF_COMPLEX;REAL_CDOT_SELF;CDOT_EQ_0;
+ CNORM2_CVECTOR_ZERO;CVECTOR_MUL_RZERO;CNORM2_EQ_0;COMPLEX_MUL_LINV;
+ CVECTOR_MUL_ID]);;
+
+
+(* ========================================================================= *)
+(* COLLINEARITY *)
+(* ========================================================================= *)
+
+(* Definition of collinearity between complex vectors.
+ * Note: This is different from collinearity between points (which is the one defined in HOL-Light library)
+ *)
+let collinear_cvectors = new_definition
+ `collinear_cvectors x (y:complex^N) <=> ?a. y = a % x \/ x = a % y`;;
+
+let COLLINEAR_CVECTORS_SYM = prove
+ (`!x y:complex^N. collinear_cvectors x y <=> collinear_cvectors y x`,
+ REWRITE_TAC[collinear_cvectors] THEN MESON_TAC[]);;
+
+let COLLINEAR_CVECTORS_0 = prove
+ (`!x:complex^N. collinear_cvectors x cvector_zero`,
+ REWRITE_TAC[collinear_cvectors] THEN GEN_TAC THEN EXISTS_TAC `Cx(&0)`
+ THEN REWRITE_TAC[CVECTOR_MUL_LZERO]);;
+
+let NON_NULL_COLLINEARS = prove
+ (`!x y:complex^N.
+ collinear_cvectors x y /\ ~(x=cvector_zero) /\ ~(y=cvector_zero)
+ ==> ?a. ~(a=Cx(&0)) /\ y = a % x`,
+ REWRITE_TAC[collinear_cvectors] THEN REPEAT STRIP_TAC THENL [
+ ASM_MESON_TAC[CVECTOR_MUL_LZERO];
+ SUBGOAL_THEN `~(a=Cx(&0))` ASSUME_TAC THENL [
+ ASM_MESON_TAC[CVECTOR_MUL_LZERO];
+ EXISTS_TAC `inv a :complex`
+ THEN ASM_REWRITE_TAC[COMPLEX_INV_EQ_0;CVECTOR_MUL_ASSOC]
+ THEN ASM_SIMP_TAC[COMPLEX_MUL_LINV;CVECTOR_MUL_ID]]]);;
+
+let COLLINEAR_LNONNULL = prove(
+ `!x y:complex^N.
+ collinear_cvectors x y /\ ~(x=cvector_zero) ==> ?a. y = a % x`,
+ REPEAT STRIP_TAC THEN ASM_CASES_TAC `y=cvector_zero:complex^N` THENL [
+ ASM_REWRITE_TAC[] THEN EXISTS_TAC `Cx(&0)`
+ THEN ASM_MESON_TAC[CVECTOR_MUL_LZERO];
+ ASM_MESON_TAC[NON_NULL_COLLINEARS] ]);;
+
+let COLLINEAR_RNONNULL = prove(
+ `!x y:complex^N.
+ collinear_cvectors x y /\ ~(y=cvector_zero) ==> ?a. x = a % y`,
+ MESON_TAC[COLLINEAR_LNONNULL;COLLINEAR_CVECTORS_SYM]);;
+
+let COLLINEAR_RUNITREAL = prove(
+ `!x y:real^N.
+ collinear_cvectors x (vector_to_cvector y) /\ norm y = &1
+ ==> x = (x cdot (vector_to_cvector y)) % vector_to_cvector y`,
+ REPEAT STRIP_TAC
+ THEN POP_ASSUM (DISTRIB [ASSUME_TAC; ASSUME_TAC o REWRITE_RULE[NORM_EQ_0;
+ GSYM VECTOR_TO_CVECTOR_ZERO_EQ] o MATCH_MP
+ (REAL_ARITH `!x. x= &1 ==> ~(x= &0)`)])
+ THEN FIRST_X_ASSUM (fun x -> FIRST_X_ASSUM (fun y ->
+ CHOOSE_THEN (SINGLE ONCE_REWRITE_TAC) (MATCH_MP COLLINEAR_RNONNULL
+ (CONJ y x))))
+ THEN REWRITE_TAC[CDOT_LMUL;CDOT_LREAL;CVECTOR_RE_VECTOR_TO_CVECTOR;
+ CVECTOR_IM_VECTOR_TO_CVECTOR;DOT_RZERO;COMPLEX_MUL_RZERO;COMPLEX_SUB_RZERO]
+ THEN POP_ASSUM ((fun x ->
+ REWRITE_TAC[x;COMPLEX_MUL_RID]) o REWRITE_RULE[NORM_EQ_1]));;
+
+let CCROSS_COLLINEAR_CVECTORS = prove
+ (`!x y:complex^3. x ccross y = cvector_zero <=> collinear_cvectors x y`,
+ REWRITE_TAC[ccross;collinear_cvectors;CART_EQ3;VECTOR_3;
+ CVECTOR_ZERO_COMPONENT;COMPLEX_SUB_0;CVECTOR_MUL_COMPONENT]
+ THEN REPEAT GEN_TAC THEN EQ_TAC
+ THENL [
+ REPEAT (POP_ASSUM MP_TAC) THEN ASM_CASES_TAC `(x:complex^3)$1 = Cx(&0)`
+ THENL [
+ ASM_CASES_TAC `(x:complex^3)$2 = Cx(&0)` THENL [
+ ASM_CASES_TAC `(x:complex^3)$3 = Cx(&0)` THENL [
+ REPEAT DISCH_TAC THEN EXISTS_TAC `Cx(&0)`
+ THEN ASM_REWRITE_TAC[COMPLEX_POLY_CLAUSES];
+ REPEAT STRIP_TAC THEN EXISTS_TAC `(y:complex^3)$3/(x:complex^3)$3`
+ THEN ASM_SIMP_TAC[COMPLEX_BALANCE_DIV_MUL]
+ THEN ASM_MESON_TAC[COMPLEX_MUL_AC];];
+ REPEAT STRIP_TAC THEN EXISTS_TAC `(y:complex^3)$2/(x:complex^3)$2`
+ THEN ASM_SIMP_TAC[COMPLEX_BALANCE_DIV_MUL]
+ THEN ASM_MESON_TAC[COMPLEX_MUL_AC]; ];
+ REPEAT STRIP_TAC THEN EXISTS_TAC `(y:complex^3)$1/(x:complex^3)$1`
+ THEN ASM_SIMP_TAC[COMPLEX_BALANCE_DIV_MUL]
+ THEN ASM_MESON_TAC[COMPLEX_MUL_AC];];
+ SIMPLE_COMPLEX_ARITH_TAC ]);;
+
+let CVECTOR_MUL_INV = prove
+ (`!a x y:complex^N. ~(a = Cx(&0)) /\ x = a % y ==> y = inv a % x`,
+ REPEAT STRIP_TAC THEN ASM_SIMP_TAC[CVECTOR_MUL_ASSOC;
+ MESON[] `(p\/q) <=> (~p ==> q)`;COMPLEX_MUL_LINV;CVECTOR_MUL_ID]);;
+
+let CVECTOR_MUL_INV2 = prove
+ (`!a x y:complex^N. ~(x = cvector_zero) /\ x = a % y ==> y = inv a % x`,
+ REPEAT STRIP_TAC THEN ASM_CASES_TAC `a=Cx(&0)`
+ THEN ASM_MESON_TAC[CVECTOR_MUL_LZERO;CVECTOR_MUL_INV]);;
+
+
+
+let COLLINEAR_CVECTORS_VECTOR_TO_CVECTOR = prove(
+ `!x y:real^N.
+ collinear_cvectors (vector_to_cvector x) (vector_to_cvector y)
+ <=> collinear {vec 0,x,y}`,
+ REWRITE_TAC[COLLINEAR_LEMMA_ALT;collinear_cvectors]
+ THEN REPEAT (STRIP_TAC ORELSE EQ_TAC) THENL [
+ POP_ASSUM MP_TAC THEN ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_MUL_COMPONENT;VECTOR_TO_CVECTOR_COMPONENT;
+ VECTOR_MUL_COMPONENT;COMPLEX_EQ;RE_CX;RE_MUL_CX]
+ THEN REPEAT STRIP_TAC THEN DISJ2_TAC THEN EXISTS_TAC `Re a`
+ THEN ASM_SIMP_TAC[];
+ REWRITE_TAC[MESON[]`(p\/q) <=> (~p ==> q)`]
+ THEN REWRITE_TAC[GSYM VECTOR_TO_CVECTOR_ZERO_EQ]
+ THEN DISCH_TAC
+ THEN SUBGOAL_TAC "" `vector_to_cvector (y:real^N) =
+ inv a % vector_to_cvector x` [ASM_MESON_TAC[CVECTOR_MUL_INV2]]
+ THEN POP_ASSUM MP_TAC THEN ONCE_REWRITE_TAC[CART_EQ]
+ THEN REWRITE_TAC[CVECTOR_MUL_COMPONENT;VECTOR_TO_CVECTOR_COMPONENT;
+ VECTOR_MUL_COMPONENT;COMPLEX_EQ;RE_CX;RE_MUL_CX]
+ THEN REPEAT STRIP_TAC THEN EXISTS_TAC `Re(inv a)` THEN ASM_SIMP_TAC[];
+ EXISTS_TAC `Cx(&0)` THEN ASM_REWRITE_TAC[VECTOR_TO_CVECTOR_ZERO;
+ CVECTOR_MUL_LZERO];
+ ASM_REWRITE_TAC[VECTOR_TO_CVECTOR_MUL] THEN EXISTS_TAC `Cx c`
+ THEN REWRITE_TAC[];
+ ]);;
+
+
+(* ========================================================================= *)
+(* ORTHOGONALITY *)
+(* ========================================================================= *)
+
+let corthogonal = new_definition
+ `corthogonal (x:complex^N) y <=> x cdot y = Cx(&0)`;;
+
+let CORTHOGONAL_SYM = prove(
+ `!x y:complex^N. corthogonal x y <=> corthogonal y x`,
+ MESON_TAC[corthogonal;CDOT_SYM;CNJ_ZERO]);;
+
+let CORTHOGONAL_0 = prove(
+ `!x:complex^N. corthogonal cvector_zero x /\ corthogonal x cvector_zero`,
+ REWRITE_TAC[corthogonal;CDOT_LZERO;CDOT_RZERO]);;
+
+let [CORTHOGONAL_LZERO;CORTHOGONAL_RZERO] = GCONJUNCTS CORTHOGONAL_0;;
+
+let CORTHOGONAL_COLLINEAR_CVECTORS = prove
+ (`!x y:complex^N.
+ collinear_cvectors x y /\ ~(x=cvector_zero) /\ ~(y=cvector_zero)
+ ==> ~(corthogonal x y)`,
+ REWRITE_TAC[collinear_cvectors;corthogonal] THEN REPEAT STRIP_TAC
+ THEN POP_ASSUM MP_TAC
+ THEN ASM_REWRITE_TAC[CDOT_RMUL;CDOT_LMUL;COMPLEX_ENTIRE;GSYM cnorm2;
+ CDOT_EQ_0;CNJ_EQ_0]
+ THEN ASM_MESON_TAC[CVECTOR_MUL_LZERO]);;
+
+let CORTHOGONAL_MUL_CLAUSES = prove
+ (`!x y a.
+ (corthogonal x y ==> corthogonal x (a%y))
+ /\ (corthogonal x y \/ a = Cx(&0) <=> corthogonal x (a%y))
+ /\ (corthogonal x y ==> corthogonal (a%x) y)
+ /\ (corthogonal x y \/ a = Cx(&0) <=> corthogonal (a%x) y)`,
+ SIMP_TAC[corthogonal;CDOT_RMUL;CDOT_LMUL;COMPLEX_ENTIRE;CNJ_EQ_0]
+ THEN MESON_TAC[]);;
+
+let [CORTHOGONAL_RMUL;CORTHOGONAL_RMUL_EQ;CORTHOGONAL_LMUL;
+ CORTHOGONAL_LMUL_EQ] = GCONJUNCTS CORTHOGONAL_MUL_CLAUSES;;
+
+let CORTHOGONAL_LRMUL_CLAUSES = prove
+ (`!x y a b.
+ (corthogonal x y ==> corthogonal (a%x) (b%y))
+ /\ (corthogonal x y \/ a = Cx(&0) \/ b = Cx(&0)
+ <=> corthogonal (a%x) (b%y))`,
+ MESON_TAC[CORTHOGONAL_MUL_CLAUSES]);;
+
+let [CORTHOGONAL_LRMUL;CORTHOGONAL_LRMUL_EQ] =
+ GCONJUNCTS CORTHOGONAL_LRMUL_CLAUSES;;
+
+let CORTHOGONAL_REAL_CLAUSES = prove
+ (`!r c.
+ (corthogonal c (vector_to_cvector r)
+ <=> orthogonal (cvector_re c) r /\ orthogonal (cvector_im c) r)
+ /\ (corthogonal (vector_to_cvector r) c
+ <=> orthogonal r (cvector_re c) /\ orthogonal r (cvector_im c))`,
+ REWRITE_TAC[corthogonal;orthogonal;CDOT_LREAL;CDOT_RREAL;COMPLEX_SUB_0;
+ COMPLEX_EQ;RE_CX;IM_CX;RE_SUB;IM_SUB;RE_ADD;IM_ADD]
+ THEN REWRITE_TAC[RE_DEF;CX_DEF;IM_DEF;complex;complex_mul;VECTOR_2;ii]
+ THEN CONV_TAC REAL_FIELD);;
+
+let [CORTHOGONAL_RREAL;CORTHOGONAL_LREAL] =
+ GCONJUNCTS CORTHOGONAL_REAL_CLAUSES;;
+
+let CORTHOGONAL_UNIT = prove
+ (`!x y:complex^N.
+ (corthogonal x (cunit y) <=> corthogonal x y)
+ /\ (corthogonal (cunit x) y <=> corthogonal x y)`,
+ REWRITE_TAC[cunit;GSYM CORTHOGONAL_MUL_CLAUSES;COMPLEX_INV_EQ_0;CX_INJ;
+ CNORM_EQ_0]
+ THEN MESON_TAC[CORTHOGONAL_0]);;
+
+let [CORTHOGONAL_RUNIT;CORTHOGONAL_LUNIT] = GCONJUNCTS CORTHOGONAL_UNIT;;
+
+let CORTHOGONAL_PROJECTION = prove(
+ `!x y:complex^N. corthogonal (x - (x cdot cunit y) % cunit y) y`,
+ REPEAT GEN_TAC THEN ASM_CASES_TAC `y=cvector_zero:complex^N`
+ THEN ASM_REWRITE_TAC[corthogonal;CDOT_RZERO]
+ THEN REWRITE_TAC[CDOT_LSUB;cunit;CVECTOR_MUL_ASSOC;GSYM cnorm2;CDOT_LMUL;
+ CDOT_RMUL;REWRITE_RULE[REAL_CNJ] (MATCH_MP REAL_INV (SPEC_ALL REAL_CX))]
+ THEN REWRITE_TAC[COMPLEX_MUL_AC;GSYM COMPLEX_INV_MUL;GSYM COMPLEX_POW_2;
+ cnorm;o_DEF;CSQRT]
+ THEN SIMP_TAC[CNORM2_POS;CX_SQRT;cnorm2;REAL_CDOT_SELF;REAL_OF_COMPLEX;CSQRT]
+ THEN ASM_SIMP_TAC[CDOT_EQ_0;COMPLEX_MUL_RINV;COMPLEX_MUL_RID;
+ COMPLEX_SUB_REFL]);;
+
+let CDOT_PYTHAGOREAN = prove
+ (`!x y:complex^N. corthogonal x y ==> cnorm2 (x+y) = cnorm2 x + cnorm2 y`,
+ SIMP_TAC[corthogonal;cnorm2;CDOT_LADD;CDOT_RADD;COMPLEX_ADD_RID;
+ COMPLEX_ADD_LID;REAL_OF_COMPLEX_ADD;REAL_CDOT_SELF;
+ MESON[CDOT_SYM;CNJ_ZERO] `x cdot y = Cx (&0) ==> y cdot x = Cx(&0)`]);;
+
+let CDOT_CAUCHY_SCHWARZ_POW_2 = prove
+ (`!x y:complex^N. norm (x cdot y) pow 2 <= cnorm2 x * cnorm2 y`,
+ REPEAT GEN_TAC THEN ASM_CASES_TAC `y = cvector_zero:complex^N`
+ THEN ASM_REWRITE_TAC[CNORM2_CVECTOR_ZERO;CDOT_RZERO;COMPLEX_NORM_0;
+ REAL_POW_2;REAL_MUL_RZERO;REAL_OF_COMPLEX_CX;REAL_LE_REFL]
+ THEN ONCE_REWRITE_TAC[MATCH_MP (MESON[CVECTOR_SUB_ADD]
+ `(!x:complex^N y. p (x - f x y) y)
+ ==> cnorm2 x * z = cnorm2 (x - f x y + f x y) * z`) CORTHOGONAL_PROJECTION]
+ THEN MATCH_MP_TAC (GEN_ALL (MATCH_MP (MESON[] `(!x y. P x y ==> f x y = (g x y:real))
+ ==> P x y /\ a <= g x y * b ==> a <= f x y * b`) CDOT_PYTHAGOREAN))
+ THEN REWRITE_TAC[GSYM CORTHOGONAL_MUL_CLAUSES;CORTHOGONAL_RUNIT;
+ CORTHOGONAL_PROJECTION]
+ THEN SIMP_TAC[cnorm2;GSYM REAL_OF_COMPLEX_ADD;REAL_CDOT_SELF;REAL_ADD;
+ GSYM REAL_OF_COMPLEX_MUL]
+ THEN REWRITE_TACS[cunit;CDOT_RMUL;CVECTOR_MUL_ASSOC;REWRITE_RULE[REAL_CNJ]
+ (MATCH_MP REAL_INV (SPEC_ALL REAL_CX));COMPLEX_MUL_AC;GSYM COMPLEX_INV_MUL;
+ GSYM COMPLEX_POW_2;cnorm;o_DEF;CSQRT;COMPLEX_ADD_LDISTRIB;cnorm2;CDOT_RMUL;
+ CNJ_MUL;CDOT_LMUL;GSYM cnorm2;
+ REWRITE_RULE[REAL_CNJ] (MATCH_MP REAL_INV (SPEC_ALL REAL_CX))]
+ THEN SIMP_TAC[CX_SQRT;CNORM2_POS;CSQRT;CX_CNORM2]
+ THEN REWRITE_TAC[SIMPLE_COMPLEX_ARITH
+ `x * ((y * inv x) * x) * (z * inv x') * inv x'
+ = (y * z) * (x * inv x) * (x * inv x' * inv x'):complex`]
+ THEN ASM_SIMP_TAC[CDOT_EQ_0;COMPLEX_MUL_RINV;COMPLEX_MUL_LID;COMPLEX_MUL_CNJ;
+ GSYM COMPLEX_INV_MUL]
+ THEN ONCE_REWRITE_TAC[
+ GSYM (MATCH_MP REAL_OF_COMPLEX (SPEC_ALL REAL_CDOT_SELF))]
+ THEN SIMP_TAC[GSYM cnorm2;GSYM CX_SQRT;CNORM2_POS;GSYM CX_MUL;
+ GSYM COMPLEX_POW_2;GSYM CX_POW;SQRT_POW_2;GSYM CX_INV]
+ THEN ASM_SIMP_TAC[REAL_MUL_RINV;CNORM2_EQ_0;REAL_MUL_RID;GSYM CX_ADD;
+ REAL_OF_COMPLEX_CX;GSYM REAL_POW_2;REAL_LE_ADDL;REAL_LE_MUL;CNORM2_POS]);;
+
+let CDOT_CAUCHY_SCHWARZ = prove
+ (`!x y:complex^N. norm (x cdot y) <= norm x * norm y`,
+ REPEAT GEN_TAC THEN MATCH_MP_TAC (REWRITE_RULE[REAL_LE_SQUARE_ABS]
+ (REAL_ARITH `&0 <= x /\ &0 <= y /\ abs x <= abs y ==> x <= y`))
+ THEN SIMP_TAC[NORM_POS_LE;CNORM_POS;REAL_LE_MUL;REAL_POW_MUL;CNORM_POW_2;
+ CDOT_CAUCHY_SCHWARZ_POW_2]);;
+
+let CDOT_CAUCHY_SCHWARZ_POW_2_EQUAL = prove
+ (`!x y:complex^N.
+ norm (x cdot y) pow 2 = cnorm2 x * cnorm2 y <=> collinear_cvectors x y`,
+ REPEAT GEN_TAC THEN ASM_CASES_TAC `y = cvector_zero:complex^N`
+ THEN ASM_REWRITE_TAC[CNORM2_CVECTOR_ZERO;CDOT_RZERO;COMPLEX_NORM_0;
+ REAL_POW_2;REAL_MUL_RZERO;REAL_OF_COMPLEX_CX;COLLINEAR_CVECTORS_0]
+ THEN EQ_TAC THENL [
+ ONCE_REWRITE_TAC[MATCH_MP (MESON[CVECTOR_SUB_ADD]
+ `(!x:complex^N y. p (x - f x y) y) ==>
+ cnorm2 x * z = cnorm2 (x - f x y + f x y) * z`) CORTHOGONAL_PROJECTION]
+ THEN MATCH_MP_TAC (GEN_ALL (MATCH_MP (MESON[]
+ `(!x y. P x y ==> g x y = (f x y:real)) ==>
+ P x y /\ (a = f x y * z ==> R) ==> (a = g x y * z ==> R)`)
+ CDOT_PYTHAGOREAN))
+ THEN REWRITE_TAC[GSYM CORTHOGONAL_MUL_CLAUSES;CORTHOGONAL_RUNIT;
+ CORTHOGONAL_PROJECTION]
+ THEN SIMP_TAC[cnorm2;GSYM REAL_OF_COMPLEX_ADD;REAL_CDOT_SELF;REAL_ADD;
+ GSYM REAL_OF_COMPLEX_MUL]
+ THEN REWRITE_TACS[cunit;CDOT_RMUL;CVECTOR_MUL_ASSOC;REWRITE_RULE[REAL_CNJ]
+ (MATCH_MP REAL_INV (SPEC_ALL REAL_CX));COMPLEX_MUL_AC;
+ GSYM COMPLEX_INV_MUL;GSYM COMPLEX_POW_2;cnorm;o_DEF;CSQRT;
+ COMPLEX_ADD_LDISTRIB;cnorm2;CDOT_RMUL;CNJ_MUL;CDOT_LMUL;GSYM cnorm2;
+ REWRITE_RULE[REAL_CNJ] (MATCH_MP REAL_INV (SPEC_ALL REAL_CX))]
+ THEN SIMP_TAC[CX_SQRT;CNORM2_POS;CSQRT;CX_CNORM2]
+ THEN REWRITE_TAC[SIMPLE_COMPLEX_ARITH
+ `x * ((y * inv x) * x) * (z * inv x') * inv x' =
+ (y * z) * (x * inv x) * (x * inv x' * inv x'):complex`]
+ THEN ONCE_REWRITE_TAC[GSYM (MATCH_MP REAL_OF_COMPLEX
+ (SPEC_ALL REAL_CDOT_SELF))]
+ THEN SIMP_TAC[GSYM cnorm2;GSYM CX_SQRT;CNORM2_POS;GSYM CX_MUL;
+ GSYM COMPLEX_POW_2;GSYM CX_POW;SQRT_POW_2;GSYM CX_INV;REAL_POW_INV]
+ THEN ASM_SIMP_TAC[REAL_MUL_RINV;CNORM2_EQ_0;REAL_MUL_RID;GSYM CX_ADD;
+ REAL_OF_COMPLEX_CX;GSYM REAL_POW_2;REAL_LE_ADDL;REAL_LE_MUL;CNORM2_POS;
+ GSYM CX_POW;REAL_POW_ONE;COMPLEX_MUL_RID;COMPLEX_MUL_CNJ;
+ REAL_ARITH `x = y + x <=> y = &0`;REAL_ENTIRE;CNORM2_EQ_0;
+ CVECTOR_SUB_EQ;collinear_cvectors]
+ THEN MESON_TAC[];
+ REWRITE_TAC[collinear_cvectors] THEN REPEAT STRIP_TAC
+ THEN ASM_REWRITE_TAC[cnorm2;CDOT_LMUL;CDOT_RMUL;COMPLEX_NORM_MUL;
+ COMPLEX_MUL_ASSOC]
+ THEN SIMP_TAC[COMPLEX_MUL_CNJ;GSYM cnorm2;COMPLEX_NORM_CNJ;GSYM CX_POW;
+ REAL_OF_COMPLEX_MUL;REAL_CX;REAL_CDOT_SELF;REAL_OF_COMPLEX_CX;
+ GSYM CNORM2_ALT]
+ THEN SIMPLE_COMPLEX_ARITH_TAC
+ ]);;
+
+let CDOT_CAUCHY_SCHWARZ_EQUAL = prove
+ (`!x y:complex^N.
+ norm (x cdot y) = norm x * norm y <=> collinear_cvectors x y`,
+ ONCE_REWRITE_TAC[REWRITE_RULE[REAL_EQ_SQUARE_ABS] (REAL_ARITH
+ `x=y <=> abs x = abs y /\ (&0 <= x /\ &0 <= y \/ x < &0 /\ y < &0)`)]
+ THEN SIMP_TAC[NORM_POS_LE;CNORM_POS;REAL_LE_MUL;REAL_POW_MUL;CNORM_POW_2;
+ CDOT_CAUCHY_SCHWARZ_POW_2_EQUAL]);;
+
+let CNORM_TRIANGLE = prove
+ (`!x y:complex^N. norm (x+y) <= norm x + norm y`,
+ REPEAT GEN_TAC THEN MATCH_MP_TAC (REWRITE_RULE[REAL_LE_SQUARE_ABS]
+ (REAL_ARITH `abs x <= abs y /\ &0 <= x /\ &0 <= y ==> x <= y`))
+ THEN SIMP_TAC[CNORM_POS;REAL_LE_ADD;REAL_ADD_POW_2;CNORM_POW_2;cnorm2;
+ CDOT_LADD;CDOT_RADD;SIMPLE_COMPLEX_ARITH `(x+y)+z+t = x+(y+z)+t:complex`;
+ ADD_CDOT_SYM;REAL_OF_COMPLEX_ADD;REAL_CDOT_SELF;REAL_CX;REAL_ADD;
+ REAL_OF_COMPLEX_CX;REAL_ARITH `x+ &2*y+z<=x+z+ &2*t <=> y<=t:real`]
+ THEN MESON_TAC[CDOT_CAUCHY_SCHWARZ;RE_NORM;REAL_LE_TRANS]);;
+
+let REAL_ABS_SUB_CNORM = prove
+ (`!x y:complex^N. abs (norm x - norm y) <= norm (x-y)`,
+ let lemma =
+ REWRITE_RULE[CVECTOR_SUB_ADD2;REAL_ARITH `x<=y+z <=> x-y<=z:real`]
+ (SPECL [`x:complex^N`;`y-x:complex^N`] CNORM_TRIANGLE)
+ in
+ REPEAT GEN_TAC
+ THEN MATCH_MP_TAC (MATCH_MP (MESON[]
+ `(!x y. P x y <=> Q x y) ==> Q x y ==> P x y`) REAL_ABS_BOUNDS)
+ THEN ONCE_REWRITE_TAC[REAL_ARITH `--x <= y <=> --y <= x`]
+ THEN REWRITE_TAC[REAL_NEG_SUB]
+ THEN REWRITE_TAC[lemma;ONCE_REWRITE_RULE[CNORM_SUB] lemma]);;
+
+(* ========================================================================= *)
+(* VSUM *)
+(* ========================================================================= *)
+
+let cvsum = new_definition
+ `(cvsum:(A->bool)->(A->complex^N)->complex^N) s f = lambda i. vsum s (\x. (f x)$i)`;;
+
+
+(* ========================================================================= *)
+(* INFINITE SUM *)
+(* ========================================================================= *)
+
+let csummable = new_definition
+ `csummable (s:num->bool) (f:num->complex^N)
+ <=> summable s (cvector_re o f) /\ summable s (cvector_im o f)`;;
+
+let cinfsum = new_definition
+ `cinfsum (s:num->bool) (f:num->complex^N) :complex^N
+ = vector_to_cvector (infsum s (\x. cvector_re (f x)))
+ + ii % vector_to_cvector (infsum s (\x.cvector_im (f x)))`;;
+
+let CSUMMABLE_FLATTEN_CVECTOR = prove
+ (`!s (f:num->complex^N). csummable s f <=> summable s (cvector_flatten o f)`,
+ REWRITE_TAC[csummable;summable;cvector_flatten;o_DEF]
+ THEN REPEAT (STRIP_TAC ORELSE EQ_TAC)
+ THENL [
+ EXISTS_TAC `pastecart (l:real^N) (l':real^N)`
+ THEN ASM_SIMP_TAC[GSYM SUMS_PASTECART];
+ EXISTS_TAC `fstcart (l:real^(N,N) finite_sum)`
+ THEN MATCH_MP_TAC (GEN_ALL (MATCH_MP (TAUT `(p /\ q <=> r) ==> (r ==> p)`)
+ (INST_TYPE [`:N`,`:M`] (SPEC_ALL SUMS_PASTECART))))
+ THEN EXISTS_TAC `(cvector_im o f) :num->real^N`
+ THEN EXISTS_TAC `sndcart (l:real^(N,N) finite_sum)`
+ THEN ASM_REWRITE_TAC[ETA_AX;o_DEF;PASTECART_FST_SND];
+ EXISTS_TAC `sndcart (l:real^(N,N) finite_sum)`
+ THEN MATCH_MP_TAC (GEN_ALL (MATCH_MP (TAUT `(p /\ q <=> r) ==> (r ==> q)`)
+ (INST_TYPE [`:N`,`:M`] (SPEC_ALL SUMS_PASTECART))))
+ THEN EXISTS_TAC `(cvector_re o f) :num->real^N`
+ THEN EXISTS_TAC `fstcart (l:real^(N,N) finite_sum)`
+ THEN ASM_REWRITE_TAC[ETA_AX;o_DEF;PASTECART_FST_SND];
+ ]);;
+
+let FLATTEN_CINFSUM = prove
+ (`!s f.
+ csummable s f ==>
+ ((cinfsum s f):complex^N) =
+ cvector_unflatten (infsum s (cvector_flatten o f))`,
+ SIMP_TAC[cinfsum;cvector_unflatten;COMPLEX_VECTOR_TRANSPOSE;LINEAR_FSTCART;
+ LINEAR_SNDCART;CSUMMABLE_FLATTEN_CVECTOR;GSYM INFSUM_LINEAR;o_DEF;
+ cvector_flatten;FSTCART_PASTECART;SNDCART_PASTECART]);;
+
+let CSUMMABLE_LINEAR = prove
+ (`!f h:complex^N->complex^M s.
+ csummable s f /\ clinear h ==> csummable s (h o f)`,
+ REWRITE_TAC[CSUMMABLE_FLATTEN_CVECTOR] THEN REPEAT STRIP_TAC
+ THEN POP_ASSUM (ASSUME_TAC o MATCH_MP FLATTEN_CLINEAR)
+ THEN SUBGOAL_THEN
+ `cvector_flatten o (h:complex^N -> complex^M) o (f:num -> complex^N) =
+ \n. (cvector_flatten o h o cvector_unflatten) (cvector_flatten (f n))`
+ (SINGLE REWRITE_TAC)
+ THENL [
+ REWRITE_TAC[o_DEF;FUN_EQ_THM] THEN GEN_TAC THEN REPEAT AP_TERM_TAC
+ THEN REWRITE_TAC[REWRITE_RULE[o_DEF;I_DEF;FUN_EQ_THM] UNFLATTEN_FLATTEN];
+ MATCH_MP_TAC SUMMABLE_LINEAR THEN ASM_SIMP_TAC[GSYM o_DEF]
+ ]);;
+
+let CINFSUM_LINEAR = prove
+ (`!f (h:complex^M->complex^N) s.
+ csummable s f /\ clinear h ==> cinfsum s (h o f) = h (cinfsum s f)`,
+ REPEAT GEN_TAC
+ THEN DISCH_THEN (fun x -> MP_TAC (CONJ (MATCH_MP CSUMMABLE_LINEAR x) x))
+ THEN SIMP_TAC[FLATTEN_CINFSUM;CSUMMABLE_FLATTEN_CVECTOR]
+ THEN REPEAT STRIP_TAC THEN POP_ASSUM (ASSUME_TAC o MATCH_MP FLATTEN_CLINEAR)
+ THEN SUBGOAL_THEN
+ `cvector_flatten o (h:complex^M->complex^N) o (f:num->complex^M) =
+ \n. (cvector_flatten o h o cvector_unflatten) ((cvector_flatten o f) n)`
+ (SINGLE REWRITE_TAC)
+ THENL [
+ REWRITE_TAC[o_DEF;FUN_EQ_THM] THEN GEN_TAC THEN REPEAT AP_TERM_TAC
+ THEN REWRITE_TAC[REWRITE_RULE[o_DEF;I_DEF;FUN_EQ_THM] UNFLATTEN_FLATTEN];
+ FIRST_ASSUM (fun x -> FIRST_ASSUM (fun y -> REWRITE_TAC[MATCH_MP
+ (MATCH_MP (REWRITE_RULE[IMP_CONJ] INFSUM_LINEAR) x) y]))
+ THEN REWRITE_TAC[o_DEF;REWRITE_RULE[o_DEF;I_DEF;FUN_EQ_THM]
+ UNFLATTEN_FLATTEN]
+ ]);;
+
diff --git a/em_model.ml b/em_model.ml
new file mode 100644
index 0000000..4907a79
--- /dev/null
+++ b/em_model.ml
@@ -0,0 +1,391 @@
+
+(* Upadted for the latest version of HOL Light (JULY 2014) *)
+(* ========================================================================= *)
+(* Formalization of Electromagnetic Optics *)
+(* *)
+(* (c) Copyright, Sanaz Khan Afshar & Vincent Aravantinos 2011-13 *)
+(* Hardware Verification Group, *)
+(* Concordia University *)
+(* *)
+(* Contact: <s_khanaf@encs.concordia.ca> *)
+(* <vincent@encs.concordia.ca> *)
+(* *)
+(* This file deals with the definition and basic theorems about the *)
+(* electromagnetic model. *)
+(* ========================================================================= *)
+
+new_type_abbrev("time",`:real`);;
+
+(********************************)
+(* Electromagnetic fields *)
+(********************************)
+
+new_type_abbrev("single_field",`:point -> time -> complex^3`);;
+new_type_abbrev("emf", `:point -> time -> complex^3 # complex^3`);;
+
+let e_of_emf = new_definition
+ `e_of_emf (emf:emf) : single_field = \r t. let (e,h) = emf r t in e`;;
+let h_of_emf = new_definition
+ `h_of_emf (emf:emf) : single_field = \r t. let (e,h) = emf r t in h`;;
+let is_valid_emf = new_definition
+ `is_valid_emf emf
+ <=> !r t. corthogonal (h_of_emf emf r t) (e_of_emf emf r t)`;;
+
+let EMF_EQ = prove
+ (`!emf:emf r t e h. emf r t = e,h
+ <=> e_of_emf emf r t = e /\ h_of_emf emf r t = h`,
+ REPEAT STRIP_TAC THEN EQ_TAC
+ THEN SIMP_TAC[e_of_emf;h_of_emf;LET_DEFs;PAIR_EQ;LAMBDA_PAIR]
+ THEN MESON_TAC[PAIR]);;
+
+let emf_at_point_mul = new_definition
+ `emf_at_point_mul (f:complex) (emf:complex^3#complex^3) :complex^3#complex^3
+ = (f % FST emf, f % SND emf)`;;
+
+overload_interface("%",
+ `emf_at_point_mul :complex->complex^3#complex^3->complex^3#complex^3`);;
+
+let emf_add = new_definition
+ `emf_add (emf1:emf) (emf2:emf) :emf = \r t.
+ (e_of_emf emf1 r t + e_of_emf emf2 r t,
+ h_of_emf emf1 r t + h_of_emf emf2 r t)`;;
+
+overload_interface("+",`emf_add :emf->emf->emf`);;
+
+
+(********************************)
+(* Plane waves *)
+(********************************)
+
+let plane_wave = new_definition
+ `plane_wave k w e h : emf = \r t. cexp (--ii * Cx(k dot r - w*t)) % (e,h)`;;
+
+let is_plane_wave = new_definition
+ `is_plane_wave emf <=>
+ is_valid_emf emf /\ ?k w e h.
+ &0 < w /\ ~(k = vec 0) /\ emf = plane_wave k w e h
+ /\ corthogonal e (vector_to_cvector k)
+ /\ corthogonal h (vector_to_cvector k)`;;
+
+let k_of_wave = new_definition
+ `k_of_wave emf = @k. ?w e h . &0 < w /\ ~(k = vec 0)
+ /\ emf = plane_wave k w e h
+ /\ corthogonal e (vector_to_cvector k)
+ /\ corthogonal h (vector_to_cvector k)`;;
+
+let w_of_wave = new_definition
+ `w_of_wave emf = @w. ?e h . &0 < w /\ emf = plane_wave (k_of_wave emf) w e h
+ /\ corthogonal e (vector_to_cvector (k_of_wave emf))
+ /\ corthogonal h (vector_to_cvector (k_of_wave emf))`;;
+
+let e_of_wave = new_definition
+ `e_of_wave emf = @e. ?h . emf = plane_wave (k_of_wave emf) (w_of_wave emf) e h
+ /\ corthogonal e (vector_to_cvector (k_of_wave emf))
+ /\ corthogonal h (vector_to_cvector (k_of_wave emf))`;;
+
+let h_of_wave = new_definition
+ `h_of_wave emf = @h.
+ emf = plane_wave (k_of_wave emf) (w_of_wave emf) (e_of_wave emf) h
+ /\ corthogonal (e_of_wave emf) (vector_to_cvector (k_of_wave emf))
+ /\ corthogonal h (vector_to_cvector (k_of_wave emf))`;;
+
+let eh_of_wave = new_definition
+ `eh_of_wave emf = (e_of_wave emf,h_of_wave emf)`;;
+
+let scalar_of_wave = new_definition
+ `scalar_of_wave emf =
+ \r t. cexp (--ii * Cx(k_of_wave emf dot r - w_of_wave emf * t))`;;
+
+let IS_PLANE_WAVE = prove
+ (`!emf. is_plane_wave emf <=>
+ is_valid_emf emf /\ &0 < w_of_wave emf /\ ~(k_of_wave emf = vec 0)
+ /\ emf = plane_wave (k_of_wave emf) (w_of_wave emf) (e_of_wave emf)
+ (h_of_wave emf)
+ /\ corthogonal (e_of_wave emf) (vector_to_cvector (k_of_wave emf))
+ /\ corthogonal (h_of_wave emf) (vector_to_cvector (k_of_wave emf))`,
+ let COMMON_TAC x =
+ GEN_TAC THEN EQ_TAC THEN REWRITE_TAC[is_plane_wave;x] THEN MESON_TAC[]
+ in
+ let lemma1 = prove
+ (`!emf. is_plane_wave emf <=> is_valid_emf emf /\ ?w e h . &0 < w
+ /\ ~(k_of_wave emf = vec 0) /\ emf = plane_wave (k_of_wave emf) w e h
+ /\ corthogonal e (vector_to_cvector (k_of_wave emf))
+ /\ corthogonal h (vector_to_cvector (k_of_wave emf))`,
+ COMMON_TAC k_of_wave)
+ in
+ let lemma2 = prove
+ (`!emf. is_plane_wave emf <=>
+ is_valid_emf emf /\ ?e h . &0 < w_of_wave emf /\ ~(k_of_wave emf = vec 0)
+ /\ emf = plane_wave (k_of_wave emf) (w_of_wave emf) e h
+ /\ corthogonal e (vector_to_cvector (k_of_wave emf))
+ /\ corthogonal h (vector_to_cvector (k_of_wave emf))`,
+ REWRITE_TAC[lemma1] THEN COMMON_TAC w_of_wave)
+ in
+ let lemma3 = prove
+ (`!emf. is_plane_wave emf <=>
+ is_valid_emf emf /\ ?h . &0 < w_of_wave emf /\ ~(k_of_wave emf = vec 0) /\
+ emf = plane_wave (k_of_wave emf) (w_of_wave emf) (e_of_wave emf) h
+ /\ corthogonal (e_of_wave emf) (vector_to_cvector (k_of_wave emf))
+ /\ corthogonal h (vector_to_cvector (k_of_wave emf))`,
+ REWRITE_TAC[lemma2] THEN COMMON_TAC e_of_wave)
+ in
+ REWRITE_TAC[lemma3] THEN COMMON_TAC h_of_wave);;
+
+let EH_OF_EMF_PLANE_WAVE = prove
+ (`!k w e h.
+ (e_of_emf (plane_wave k w e h) = \r t. cexp (--ii * Cx(k dot r - w*t)) % e )
+ /\ h_of_emf (plane_wave k w e h) = \r t.
+ cexp (--ii * Cx(k dot r - w*t)) % h`,
+ REWRITE_TAC[e_of_emf;h_of_emf;plane_wave;emf_at_point_mul;LET_DEFs]);;
+
+let non_null_wave = new_definition
+ `non_null_wave emf <=>
+ is_plane_wave emf /\ ~(e_of_wave emf = cvector_zero)
+ /\ ~(h_of_wave emf = cvector_zero)`;;
+
+
+(***********************************)
+(* Plane interface between mediums *)
+(***********************************)
+
+ (* A medium is characterized by its refractive index *)
+new_type_abbrev("medium", `:real`);;
+
+(* The real^3 vector is a normal to the plane. It is needed in order to fix an
+ * orientation so that we can clearly say on which side is medium 1 or 2
+ * respectively. *)
+new_type_abbrev("interface",`:medium # medium # plane # real^3`);;
+
+let FORALL_INTERFACE_THM = prove
+ (`!P. (!(i:interface). P i) <=> (!n1 n2 p n. P (n1,n2,p,n))`,
+ MESON_TAC[PAIR_SURJECTIVE]);;
+
+let is_valid_interface = new_definition
+ `is_valid_interface (i:interface) <=>
+ let (n1,n2,p,n) = i in
+ &0 < n1 /\ &0 < n2 /\ plane p /\ is_normal_to_plane n p`;;
+
+let plane_of_interface = new_definition
+ `plane_of_interface (i:interface) = let (n1,n2,p,n) = i in p`;;
+let normal_of_interface = new_definition
+ `normal_of_interface (i:interface) = let (n1,n2,p,n) = i in n`;;
+let n1_of_interface = new_definition
+ `n1_of_interface (i:interface) = let (n1,n2,p,n) = i in n1`;;
+let n2_of_interface = new_definition
+ `n2_of_interface (i:interface) = let (n1,n2,p,n) = i in n2`;;
+
+(* Helpers *)
+let IS_VALID_INTERFACE_IS_NORMAL_TO_PLANE = prove
+ (`!i. is_valid_interface i ==>
+ is_normal_to_plane (normal_of_interface i) (plane_of_interface i)`,
+ SIMP_TAC[FORALL_INTERFACE_THM;is_valid_interface;LET_DEFs;
+ normal_of_interface;plane_of_interface]);;
+
+let NORMAL_OF_INTERFACE_NON_NULL = prove
+ (`!i. is_valid_interface i ==> ~(normal_of_interface i = vec 0)`,
+ SIMP_TAC[FORALL_INTERFACE_THM;is_valid_interface;is_normal_to_plane;
+ normal_of_interface;LET_DEFs]);;
+let IS_VALID_INTERFACE_PLANE = prove
+ (`!i. is_valid_interface i ==> plane (plane_of_interface i)`,
+ SIMP_TAC[FORALL_INTERFACE_THM;is_valid_interface;plane_of_interface;
+ LET_DEFs]);;
+
+(* [p] is the point where continuity holds; [n] is the normal to the tangent
+ * plane at that point. *)
+
+let boundary_conditions = new_definition
+ `boundary_conditions emf1 emf2 p r t <=>
+ !n. is_normal_to_plane n p ==>
+ let n_ccross = (ccross) (vector_to_cvector n) in
+ n_ccross (e_of_emf emf1 r t) = n_ccross (e_of_emf emf2 r t)
+ /\ n_ccross (h_of_emf emf1 r t) = n_ccross (h_of_emf emf2 r t)`;;
+
+let eta0 = new_specification ["eta0"]
+ (EXISTS (`?x. &0 < x`,`&1`) (REAL_ARITH `&0 < &1`));;
+let k0 = new_specification ["k0"]
+ (EXISTS (`?x. &0 < x`,`&1`) (REAL_ARITH `&0 < &1`));;
+
+let is_plane_wave_at_interface = new_definition
+ `is_plane_wave_at_interface i emf_i emf_r emf_t <=>
+ is_valid_interface i /\ non_null_wave emf_i /\ is_plane_wave emf_r
+ /\ is_plane_wave emf_t
+ /\ let (n1,n2,p,n) = i in
+ (!pt. pt IN p ==> !t. boundary_conditions (emf_i + emf_r) emf_t p pt t) /\
+ (let (k_i,k_r,k_t) = map_triple k_of_wave (emf_i, emf_r, emf_t) in
+ (k_i dot n > &0 /\ k_r dot n <= &0 /\ k_t dot n >= &0) /\
+ norm k_i = k0 * n1 /\ norm k_r = k0 * n1 /\ norm k_t = k0 * n2) /\
+ let emf_in_medium = \emf n.
+ h_of_wave emf = Cx(&1 / (eta0 * k0)) %
+ (vector_to_cvector (k_of_wave emf) ccross (e_of_wave emf))
+ in
+ emf_in_medium emf_i n1 /\ emf_in_medium emf_r n1
+ /\ emf_in_medium emf_t n2`;;
+
+let IS_PLANE_WAVE_AT_INTERFACE_MAGNITUDES_RELATION = prove(
+ `!i emf_i emf_r emf_t. is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ let (n1,n2,p,n) = i in
+ let (e_i,e_r,e_t) = map_triple (norm o e_of_wave) (emf_i,emf_r,emf_t) in
+ let (h_i,h_r,h_t) = map_triple (norm o h_of_wave) (emf_i,emf_r,emf_t) in
+ h_i = e_i * n1_of_interface i / eta0
+ /\ h_r = e_r * n1_of_interface i / eta0
+ /\ h_t = e_t * n2_of_interface i / eta0`,
+ SIMP_TAC[FORALL_INTERFACE_THM;LET_DEFs;map_triple;o_DEF;n1_of_interface;
+ n2_of_interface;is_plane_wave_at_interface;GSYM CCROSS_LMUL;
+ GSYM VECTOR_TO_CVECTOR_MUL;CCROSS_LREAL;CNORM_NORM_2;
+ REWRITE_RULE[REAL_ARITH `x+y=z <=> x=z-y:real`] NORM_CROSS_DOT;non_null_wave;
+ REAL_ARITH `(x-y)+(z-t)=(x+z)-(y+t):real`;REAL_POW_MUL;
+ GSYM REAL_ADD_LDISTRIB;IS_PLANE_WAVE;CORTHOGONAL_REAL_CLAUSES;DOT_LMUL]
+ THEN ONCE_REWRITE_TAC[ORTHOGONAL_SYM] THEN SIMP_TAC[orthogonal;REAL_POW_2;
+ REAL_MUL_RZERO;REAL_ADD_RID;REAL_SUB_RZERO]
+ THEN SIMP_TAC[GSYM REAL_POW_2;REAL_LE_POW_2;MESON[REAL_LE_ADD;REAL_LE_POW_2]
+ `!x y. &0 <= x pow 2 + y pow 2`;SQRT_MUL;CX_MUL]
+ THEN REWRITE_TAC[GSYM CNORM_NORM_2;VECTOR_TO_CVECTOR_CVECTOR_RE_IM]
+(* THEN SIMP_TAC[NORM_POS_LE;POW_2_SQRT;NORM_MUL;real_div;REAL_MUL_LID] *)
+ THEN SIMP_TAC[NORM_POS_LE;SQRT_POW_2;NORM_MUL;real_div;REAL_MUL_LID]
+ THEN REWRITE_TAC[MESON[REAL_LT_IMP_LE;REAL_ABS_REFL;REAL_LE_INV;REAL_LE_MUL;
+ eta0;k0] `abs (inv (eta0 * k0)) = inv (eta0 * k0)`]
+ THEN REWRITE_TAC[REAL_INV_MUL;REAL_ARITH `(x*y)*z*t=x*(y*z)*t:real`;MATCH_MP
+ REAL_MUL_LINV (GSYM (MATCH_MP REAL_LT_IMP_NE k0));REAL_MUL_LID]
+ THEN REWRITE_TAC[COMPLEX_MUL_SYM;REAL_MUL_SYM]);;
+
+(* Helpers *)
+let IS_PLANE_WAVE_AT_INTERFACE_BOUNDARY_CONDITIONS = prove
+ (`!i emf_i emf_r emf_t.
+ is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ !pt. pt IN (plane_of_interface i) ==> !t.
+ boundary_conditions (emf_i + emf_r) emf_t (plane_of_interface i) pt t`,
+ SIMP_TAC[FORALL_INTERFACE_THM;is_plane_wave_at_interface;plane_of_interface;
+ LET_DEFs]);;
+
+let IS_PLANE_WAVE_AT_INTERFACE_IS_NORMAL_TO_PLANE = prove
+ (`!i emf_i emf_r emf_t. is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ is_normal_to_plane (normal_of_interface i) (plane_of_interface i)`,
+ SIMP_TAC[is_plane_wave_at_interface;IS_VALID_INTERFACE_IS_NORMAL_TO_PLANE]);;
+
+let IS_PLANE_WAVE_AT_INTERFACE_PLANE = prove
+ (`!i emf_i emf_r emf_t. is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ plane (plane_of_interface i)`,
+ SIMP_TAC[is_plane_wave_at_interface;IS_VALID_INTERFACE_PLANE]);;
+
+let IS_PLANE_WAVE_AT_INTERFACE_NON_NULL_NORMAL = prove
+ (`!i emf_i emf_r emf_t. is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ ~(normal_of_interface i = vec 0)`,
+ SIMP_TAC[FORALL_INTERFACE_THM;is_plane_wave_at_interface;is_valid_interface;
+ normal_of_interface;is_normal_to_plane;LET_DEFs]);;
+
+
+(***********************************)
+(* TE and TM modes *)
+(***********************************)
+
+let is_mode_axis = new_definition
+ `is_mode_axis field (i:interface) emf_i emf_r emf_t v <=>
+ orthogonal v (normal_of_interface i) /\ norm v = &1 /\ !r t.
+ collinear_cvectors (field emf_i r t) (vector_to_cvector v)
+ /\ collinear_cvectors (field emf_r r t) (vector_to_cvector v)
+ /\ collinear_cvectors (field emf_t r t) (vector_to_cvector v)`;;
+
+let transverse_mode = new_definition
+ `transverse_mode field (i:interface) emf_i emf_r emf_t <=>
+ ?v. is_mode_axis field (i:interface) emf_i emf_r emf_t v`;;
+
+let TE_mode = new_definition `TE_mode = transverse_mode e_of_emf`;;
+let TM_mode = new_definition `TM_mode = transverse_mode h_of_emf`;;
+
+let mode_axis = new_definition `mode_axis field i emf_i emf_r emf_t =
+ @v. is_mode_axis field i emf_i emf_r emf_t v`;;
+
+let TE_mode_axis = new_definition `TE_mode_axis = mode_axis e_of_emf`;;
+let TM_mode_axis = new_definition `TM_mode_axis = mode_axis h_of_emf`;;
+
+let TRANSVERSE_MODE_THM = prove
+ (`!field i emf_i emf_r emf_t. transverse_mode field i emf_i emf_r emf_t <=>
+ is_mode_axis field i emf_i emf_r emf_t (mode_axis field i emf_i emf_r emf_t)`,
+ REWRITE_TAC[transverse_mode;mode_axis] THEN SELECT_ELIM_TAC
+ THEN REWRITE_TAC[]);;
+
+let TE_MODE_THM = prove
+ (`!i emf_i emf_r emf_t. TE_mode i emf_i emf_r emf_t <=>
+ is_mode_axis e_of_emf i emf_i emf_r emf_t
+ (TE_mode_axis i emf_i emf_r emf_t)`,
+ REWRITE_TAC[TE_mode;TE_mode_axis;TRANSVERSE_MODE_THM]);;
+
+let TM_MODE_THM = prove
+ (`!i emf_i emf_r emf_t. TM_mode i emf_i emf_r emf_t <=>
+ is_mode_axis h_of_emf i emf_i emf_r emf_t
+ (TM_mode_axis i emf_i emf_r emf_t)`,
+ REWRITE_TAC[TM_mode;TM_mode_axis;TRANSVERSE_MODE_THM]);;
+
+let MODE_AXIS_ORTHOGONAL_N = prove
+ (`!field i emf_i emf_r emf_t. transverse_mode field i emf_i emf_r emf_t ==>
+ orthogonal (mode_axis field i emf_i emf_r emf_t) (normal_of_interface i)`,
+ SIMP_TAC[TRANSVERSE_MODE_THM;is_mode_axis]);;
+
+let TE_MODE_AXIS_ORTHOGONAL_N = prove
+ (`!i emf_i emf_r emf_t. TE_mode i emf_i emf_r emf_t ==>
+ orthogonal (TE_mode_axis i emf_i emf_r emf_t) (normal_of_interface i)`,
+ REWRITE_TAC[TE_mode;TE_mode_axis;MODE_AXIS_ORTHOGONAL_N]);;
+
+let TM_MODE_AXIS_ORTHOGONAL_N = prove
+ (`!i emf_i emf_r emf_t. TM_mode i emf_i emf_r emf_t ==>
+ orthogonal (TM_mode_axis i emf_i emf_r emf_t) (normal_of_interface i)`,
+ REWRITE_TAC[TM_mode;TM_mode_axis;MODE_AXIS_ORTHOGONAL_N]);;
+
+let TRANSVERSE_MODE_PROJ = prove
+ (`!field i emf_i emf_r emf_t. transverse_mode field i emf_i emf_r emf_t ==>
+ !r t. let x = vector_to_cvector (mode_axis field i emf_i emf_r emf_t) in
+ field emf_i r t = (field emf_i r t cdot x) % x
+ /\ field emf_r r t = (field emf_r r t cdot x) % x
+ /\ field emf_t r t = (field emf_t r t cdot x) % x`,
+ REWRITE_TAC[TRANSVERSE_MODE_THM;is_mode_axis;LET_DEFs] THEN REPEAT STRIP_TAC
+ THEN POP_ASSUM (STRIP_ASSUME_TAC o SPEC_ALL)
+ THEN FIRST_ASSUM (ASSUME_TAC o REWRITE_RULE[NORM_EQ_0;GSYM
+ VECTOR_TO_CVECTOR_ZERO_EQ] o MATCH_MP
+ (REAL_ARITH `!x. x= &1 ==> ~(x= &0)`))
+ THEN ASM_SIMP_TAC[COLLINEAR_RUNITREAL]);;
+
+let TE_MODE_PROJ = prove
+ (`!i emf_i emf_r emf_t. TE_mode i emf_i emf_r emf_t ==> !r t.
+ let x = vector_to_cvector (TE_mode_axis i emf_i emf_r emf_t) in
+ e_of_emf emf_i r t = (e_of_emf emf_i r t cdot x) % x
+ /\ e_of_emf emf_r r t = (e_of_emf emf_r r t cdot x) % x
+ /\ e_of_emf emf_t r t = (e_of_emf emf_t r t cdot x) % x`,
+ REWRITE_TAC[TE_mode;TE_mode_axis;TRANSVERSE_MODE_PROJ]);;
+
+let TM_MODE_PROJ = prove
+ (`!i emf_i emf_r emf_t. TM_mode i emf_i emf_r emf_t ==> !r t.
+ let x = vector_to_cvector (TM_mode_axis i emf_i emf_r emf_t) in
+ h_of_emf emf_i r t = (h_of_emf emf_i r t cdot x) % x
+ /\ h_of_emf emf_r r t = (h_of_emf emf_r r t cdot x) % x
+ /\ h_of_emf emf_t r t = (h_of_emf emf_t r t cdot x) % x`,
+ REWRITE_TAC[TM_mode;TM_mode_axis;TRANSVERSE_MODE_PROJ]);;
+
+let TE_MODE_PLANEWAVE_PROJ = prove
+ (`!i emf_i emf_r emf_t. TE_mode i emf_i emf_r emf_t /\ is_plane_wave emf_i
+ /\ is_plane_wave emf_r /\ is_plane_wave emf_t ==>
+ let x = vector_to_cvector (TE_mode_axis i emf_i emf_r emf_t) in
+ e_of_wave emf_i = (e_of_wave emf_i cdot x) % x
+ /\ e_of_wave emf_r = (e_of_wave emf_r cdot x) % x
+ /\ e_of_wave emf_t = (e_of_wave emf_t cdot x) % x`,
+ REWRITE_TAC[IS_PLANE_WAVE] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN TE_MODE_PROJ (MP_TAC o REWRITE_RULE[LET_DEFs])
+ THEN ASSUM_LIST (REWRITE_TAC o mapfilter (REWRITE_RULE[EH_OF_EMF_PLANE_WAVE]
+ o MATCH_MP (MESON[] `x=y ==> e_of_emf x = e_of_emf y`)))
+ THEN REWRITE_TAC[CDOT_LMUL;GSYM CVECTOR_MUL_ASSOC;CVECTOR_MUL_LCANCEL;
+ CEXP_NZ;LET_DEFs]);;
+
+let TM_MODE_PLANEWAVE_PROJ = prove
+ (`!i emf_i emf_r emf_t. TM_mode i emf_i emf_r emf_t /\ is_plane_wave emf_i
+ /\ is_plane_wave emf_r /\ is_plane_wave emf_t ==>
+ let x = vector_to_cvector (TM_mode_axis i emf_i emf_r emf_t) in
+ h_of_wave emf_i = (h_of_wave emf_i cdot x) % x
+ /\ h_of_wave emf_r = (h_of_wave emf_r cdot x) % x
+ /\ h_of_wave emf_t = (h_of_wave emf_t cdot x) % x`,
+ REWRITE_TAC[IS_PLANE_WAVE] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN TM_MODE_PROJ (MP_TAC o REWRITE_RULE[LET_DEFs])
+ THEN ASSUM_LIST (REWRITE_TAC o mapfilter (REWRITE_RULE[EH_OF_EMF_PLANE_WAVE]
+ o MATCH_MP (MESON[] `x=y ==> h_of_emf x = h_of_emf y`)))
+ THEN REWRITE_TAC[CDOT_LMUL;GSYM CVECTOR_MUL_ASSOC;CVECTOR_MUL_LCANCEL;CEXP_NZ;
+ LET_DEFs]);;
+
+
diff --git a/frequency_equalities.ml b/frequency_equalities.ml
new file mode 100644
index 0000000..bd67006
--- /dev/null
+++ b/frequency_equalities.ml
@@ -0,0 +1,311 @@
+(* ========================================================================= *)
+(* Formalization of Electromagnetic Optics *)
+(* *)
+(* (c) Copyright, Sanaz Khan Afshar & Vincent Aravantinos 2011-13 *)
+(* Hardware Verification Group, *)
+(* Concordia University *)
+(* *)
+(* Contact: <s_khanaf@encs.concordia.ca> *)
+(* <vincent@encs.concordia.ca> *)
+(* *)
+(* Proving that (!x. A*e^(iax)+B*e^(ibx)=C*e^(icx)) ==> a=b=c /\ A+B=C *)
+(* Several versions of this result are proved that are useful in proofs. *)
+(* ========================================================================= *)
+
+
+
+(** DERIVATION OF FUNCTIONS WITH REAL ARGUMENTS AND COMPLEX VALUES
+ *
+ * The proof of the lemma requires derivation of functions whose argument is
+ * real but whose value is complex. The theoretical background to do this in
+ * HOL-Light is already present in the library, but the automation is not well
+ * developped. In addition, we need several intermediate results which are
+ * available only for complex derivations or only for general derivation, etc.
+ * Hence we need a bit of work before getting the results we want.
+ *)
+
+let REAL_CMUL_TO_VECTOR_MUL = prove
+ (`!a x. a * Cx x = x % a`,
+ REWRITE_TAC[complex_mul;RE_DEF;IM_DEF;CX_DEF;complex;CART_EQ;DIMINDEX_2;
+ FORALL_2;VECTOR_2;VECTOR_MUL_COMPONENT] THEN SIMPLE_COMPLEX_ARITH_TAC);;
+
+let LINEAR_REAL_CMUL = prove
+ (`!a. linear (\x. a * Cx (drop x))`,
+ REWRITE_TAC[REAL_CMUL_TO_VECTOR_MUL;linear;drop;VECTOR_ADD_COMPONENT;
+ VECTOR_MUL_COMPONENT]
+ THEN VECTOR_ARITH_TAC);;
+
+let HAS_VECTOR_DERIVATIVE_REAL_CMUL = prove
+ (`!a x. ((\x. a * Cx (drop x)) has_vector_derivative a) (at x)`,
+ SIMP_TAC[has_vector_derivative;GSYM REAL_CMUL_TO_VECTOR_MUL;
+ HAS_DERIVATIVE_LINEAR;LINEAR_REAL_CMUL]);;
+
+let HAS_COMPLEX_DERIVATIVE_MUL_CEXP = prove
+ (`!a z. ((\z. a * cexp z) has_complex_derivative (a * cexp z)) (at z)`,
+ REPEAT GEN_TAC THEN COMPLEX_DIFF_TAC THEN SIMPLE_COMPLEX_ARITH_TAC);;
+
+let HAS_VECTOR_DERIVATIVE_CEXP = prove
+ (`!A a x. ((\x. A * cexp(a * Cx(drop x))) has_vector_derivative (A * a *
+ cexp (a * Cx(drop x)))) (at x)`,
+ REWRITE_TAC[has_vector_derivative; GSYM (REWRITE_CONV[o_DEF]
+ `(\x. A * cexp x) o (\x. a * Cx(drop x))`); GSYM REAL_CMUL_TO_VECTOR_MUL]
+ THEN SUBGOAL_THEN
+ `!A a x. (\x'. (A * a * cexp (a * Cx (drop x))) * Cx (drop x')) = (\x'. A *
+ cexp(a * Cx(drop x)) * x') o (\x. a * Cx(drop x))` (SINGLE REWRITE_TAC)
+ THENL ON_FIRST_GOAL (REWRITE_TAC[o_DEF;FUN_EQ_THM]
+ THEN SIMPLE_COMPLEX_ARITH_TAC)
+ THEN REPEAT GEN_TAC THEN MATCH_MP_TAC DIFF_CHAIN_AT THEN CONJ_TAC
+ THENL [
+ REWRITE_TAC[REAL_CMUL_TO_VECTOR_MUL;GSYM has_vector_derivative]
+ THEN REWRITE_TAC[GSYM REAL_CMUL_TO_VECTOR_MUL;
+ HAS_VECTOR_DERIVATIVE_REAL_CMUL];
+ REWRITE_TAC[] THEN MATCH_ACCEPT_TAC (REWRITE_RULE[has_complex_derivative;
+ GSYM COMPLEX_MUL_ASSOC] HAS_COMPLEX_DERIVATIVE_MUL_CEXP);
+ ]);;
+
+let HAS_VECTOR_DERIVATIVE_SUM_CEXP = prove
+ (`!A B a b x. ((\x. A * cexp(a * Cx(drop x)) + B * cexp(b * Cx(drop x)))
+ has_vector_derivative (A * a * cexp (a * Cx(drop x)) + B * b * cexp (b *
+ Cx(drop x)))) (at x)`,
+ REPEAT GEN_TAC THEN MATCH_MP_TAC HAS_VECTOR_DERIVATIVE_ADD THEN CONJ_TAC
+ THEN MATCH_ACCEPT_TAC HAS_VECTOR_DERIVATIVE_CEXP);;
+
+
+(** MAIN RESULT *)
+
+let WAVE_SUM_EQ_CORE = prove
+ (`!a b c A B C. ~(A = Cx(&0)) /\ ~(B = Cx(&0)) /\ ~(C = Cx(&0))
+ /\ (!x. A * cexp (a * Cx x) + B * cexp (b * Cx x) = C * cexp (c * Cx x))
+ ==> a = b /\ b = c /\ A + B = C`,
+ REPEAT GEN_TAC THEN STRIP_TAC
+ THEN SUBGOAL_THEN
+ `!x. a * A * cexp (a*Cx x) + b * B * cexp (b*Cx x) = c * C * cexp (c*Cx x)`
+ ASSUME_TAC
+ THENL ON_FIRST_GOAL
+ (REWRITE_TAC[FORALL_DROP] THEN GEN_TAC
+ THEN MATCH_MP_TAC VECTOR_DERIVATIVE_UNIQUE_AT
+ THEN MAP_EVERY EXISTS_TAC [`\x. C * cexp (c * Cx(drop x))`;`x:real^1`]
+ THEN REWRITE_TAC
+ [MESON[COMPLEX_MUL_AC] `!A a. a * A * cexp x = A * a * cexp x`]
+ THEN CONJ_TAC THEN TRY (MATCH_ACCEPT_TAC HAS_VECTOR_DERIVATIVE_CEXP)
+ THEN POP_ASSUM (fun x -> REWRITE_TAC[GSYM x])
+ THEN MATCH_ACCEPT_TAC HAS_VECTOR_DERIVATIVE_SUM_CEXP)
+ THEN SUBGOAL_THEN `!x. a pow 2 * A * cexp (a*Cx x) + b
+ pow 2 * B * cexp (b*Cx x) = c pow 2 * C * cexp (c*Cx x)` ASSUME_TAC
+ THENL ON_FIRST_GOAL
+ (REWRITE_TAC[FORALL_DROP] THEN GEN_TAC
+ THEN MATCH_MP_TAC VECTOR_DERIVATIVE_UNIQUE_AT
+ THEN MAP_EVERY EXISTS_TAC [`\x.(C * c) * cexp (c * Cx(drop x))`;`x:real^1`]
+ THEN REWRITE_TAC[COMPLEX_POW_2]
+ THEN REWRITE_TAC
+ [MESON[COMPLEX_MUL_AC] `!A a. (a*a) * A * cexp x = (A*a)*a*cexp x`]
+ THEN CONJ_TAC THEN TRY (MATCH_ACCEPT_TAC HAS_VECTOR_DERIVATIVE_CEXP)
+ THEN
+ let assoc_on_cexp =
+ MESON[COMPLEX_MUL_AC] `!C c. (C*c) * cexp x = c*C*cexp x`
+ in
+ REWRITE_TAC[assoc_on_cexp]
+ THEN POP_ASSUM (fun x -> REWRITE_TAC[GSYM x; GSYM assoc_on_cexp])
+ THEN REWRITE_TAC
+ [MESON[COMPLEX_MUL_AC] `!A a. (a*A*a) * cexp x = (A*a)*a*cexp x`]
+ THEN MATCH_ACCEPT_TAC HAS_VECTOR_DERIVATIVE_SUM_CEXP)
+ THEN RULE_ASSUM_TAC (try_or_id (REWRITE_RULE[COMPLEX_MUL_RZERO;
+ COMPLEX_ADD_RID;CEXP_0;COMPLEX_MUL_RID] o SPEC `&0`))
+ THEN SUBGOAL_THEN
+ `(a*A+b*B) pow 2 = a pow 2 * A pow 2 + b pow 2 * B pow 2 + Cx(&2)*a*b*A*B`
+ MP_TAC
+ THENL ON_FIRST_GOAL
+ (REWRITE_TAC[COMPLEX_POW_2] THEN SIMPLE_COMPLEX_ARITH_TAC)
+ THEN SUBGOAL_THEN
+ `(a*A+b*B) pow 2 = a pow 2 * A pow 2 + b pow 2 * B pow 2 + a pow 2 * A *
+ B + b pow 2 * A * B :complex`
+ (SINGLE REWRITE_TAC)
+ THENL ON_FIRST_GOAL
+ (GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV) [COMPLEX_POW_2]
+ THEN ASM
+ (GEN_REWRITE_TAC (RATOR_CONV o RAND_CONV o RATOR_CONV o DEPTH_CONV)) []
+ THEN FIND_ASSUM (SINGLE REWRITE_TAC o GSYM) `A+B=C:complex`
+ THEN FIND_ASSUM (SINGLE REWRITE_TAC) `a * A + b * B = c * C:complex`
+ THEN SUBGOAL_TAC "" `(c * (A + B)) * c * C = (A+B) * c pow 2 * C:complex`
+ [REWRITE_TAC[COMPLEX_POW_2] THEN SIMPLE_COMPLEX_ARITH_TAC]
+ THEN POP_ASSUM (SINGLE REWRITE_TAC)
+ THEN FIND_ASSUM (SINGLE REWRITE_TAC o GSYM)
+ (`a pow 2 * A + b pow 2 * B = c pow 2 * C:complex`)
+ THEN REWRITE_TAC[COMPLEX_POW_2] THEN SIMPLE_COMPLEX_ARITH_TAC)
+ THEN REWRITE_TAC[COMPLEX_EQ_ADD_LCANCEL]
+ THEN SUBGOAL_THEN `a pow 2 * A * B + b pow 2 * A * B = Cx (&2) * a * b * A *
+ B <=> (a-b) pow 2 * A * B = Cx(&0)` (SINGLE REWRITE_TAC)
+ THENL ON_FIRST_GOAL
+ (REWRITE_TAC[COMPLEX_POW_2] THEN SIMPLE_COMPLEX_ARITH_TAC)
+ THEN REWRITE_TAC[COMPLEX_ENTIRE]
+ THEN RULE_ASSUM_TAC (SINGLE PURE_REWRITE_RULE (TAUT `~p <=> (p<=>F)`))
+ THEN ASM_REWRITE_TAC[COMPLEX_POW_2;COMPLEX_ENTIRE;COMPLEX_SUB_0]
+ THEN DISCH_THEN (fun x -> REWRITE_ASSUMPTIONS[x;GSYM COMPLEX_ADD_LDISTRIB]
+ THEN REWRITE_TAC[x])
+ THEN POP_ASSUM (K ALL_TAC) THEN POP_ASSUM MP_TAC
+ THEN ASM_REWRITE_TAC[COMPLEX_EQ_MUL_RCANCEL]);;
+
+
+(** WEAKER RESULTS
+ * - with only two exponentials (WAVE_SUM_EQ_SINGLE)
+ * - with B possibly equal to zero (WAVE_SUM_EQ_WEAK1)
+ * - with C possibly equal to zero (WAVE_SUM_EQ_WEAK2)
+ *)
+let WAVE_SUM_EQ_SINGLE = prove
+ (`!a b A B.
+ ~(A = Cx(&0)) /\ (!x. A * cexp (a * Cx x) = B * cexp (b * Cx x))
+ ==> a = b /\ A = B`,
+ REPEAT GEN_TAC THEN STRIP_TAC THEN SUBGOAL_THEN `~(B = Cx(&0))` ASSUME_TAC
+ THENL ON_FIRST_GOAL (POP_ASSUM (fun x -> ASM_REWRITE_TAC[REWRITE_RULE
+ [COMPLEX_MUL_RZERO;CEXP_0;COMPLEX_MUL_RID] (SPEC `&0` (GSYM x))]))
+ THEN SUBGOAL_THEN
+ `!x. A * cexp (a * Cx x) + A * cexp (a * Cx x)
+ = (Cx(&2)*B) * cexp (b * Cx x)`
+ ASSUME_TAC
+ THENL ON_FIRST_GOAL (GEN_TAC THEN POP_ASSUM (K ALL_TAC)
+ THEN POP_ASSUM (MP_TAC o SPEC_ALL) THEN SIMPLE_COMPLEX_ARITH_TAC)
+ THEN SUBGOAL_THEN `~(Cx(&2)*B=Cx(&0))` ASSUME_TAC
+ THENL ON_FIRST_GOAL
+ (POP_ASSUM (K ALL_TAC) THEN POP_ASSUM MP_TAC THEN SIMPLE_COMPLEX_ARITH_TAC)
+ THEN ASM_MESON_TAC
+ [WAVE_SUM_EQ_CORE;SIMPLE_COMPLEX_ARITH `A + A = Cx(&2) * B <=> A = B`]);;
+
+let WAVE_SUM_EQ_WEAK1 = prove
+ (`!a b c A B C.
+ ~(A = Cx(&0)) /\ ~(C = Cx(&0))
+ /\ (!x. A * cexp (a * Cx x) + B * cexp (b * Cx x) = C * cexp (c * Cx x))
+ ==> a = c /\ A + B = C`,
+ REPEAT GEN_TAC THEN ASM_CASES_TAC `B = Cx(&0)` THENL [
+ ASM_REWRITE_TAC[COMPLEX_MUL_LZERO;COMPLEX_ADD_RID]
+ THEN MESON_TAC[WAVE_SUM_EQ_SINGLE];
+ REPEAT STRIP_TAC
+ THENL ON_FIRST_GOAL (MATCH_MP_TAC EQ_TRANS THEN EXISTS_TAC `b:complex`)
+ THEN ASM_MESON_TAC[WAVE_SUM_EQ_CORE];
+ ]);;
+
+let WAVE_SUM_EQ_WEAK2 = prove
+ (`!a b c A B C.
+ ~(A = Cx(&0)) /\ ~(B = Cx(&0))
+ /\ (!x. A * cexp (a * Cx x) + B * cexp (b * Cx x) = C * cexp (c * Cx x))
+ ==> a = b /\ A + B = C`,
+ REPEAT GEN_TAC THEN ASM_CASES_TAC `C = Cx(&0)` THENL [
+ ASM_REWRITE_TAC[COMPLEX_MUL_LZERO;COMPLEX_LNEG_UNIQ;COMPLEX_NEG_LMUL]
+ THEN MESON_TAC[WAVE_SUM_EQ_SINGLE;COMPLEX_NEG_EQ_0];
+ ASM_MESON_TAC[WAVE_SUM_EQ_CORE];
+ ]);;
+
+
+(** VECTORIAL VERSIONS OF THE ABOVE THEOREMS **)
+let VEC_WAVE_SUM_EQ_CORE = prove
+ (`!a b c A B C:complex^N.
+ ~(A = cvector_zero) /\ ~(B = cvector_zero) /\ ~(C = cvector_zero)
+ /\ (!x. cexp (a * Cx x) % A + cexp (b * Cx x) % B = cexp (c * Cx x) % C)
+ ==> a = b /\ b = c /\ A + B = C`,
+ REWRITE_TAC[CVECTOR_NON_ZERO;IMP_CONJ;
+ MESON[] `(a ==> b/\c/\d) <=> (a==>b/\c) /\ (b/\c==>a==>d)`]
+ THEN REPEAT GEN_TAC
+ THEN MAP_EVERY (DISCH_THEN o X_CHOOSE_TAC) [`i:num`;`j:num`;`k:num`]
+ THEN CONJ_TAC
+ THENL ON_FIRST_GOAL (MAP_EVERY ASM_CASES_TAC [`i=(j:num)`;`j=(k:num)`])
+ THENL
+ let cexp_SYM = SIMPLE_COMPLEX_ARITH `cexp x * y = y * cexp x` in
+ let PROJECT_ON t =
+ MAP_ANTECEDENT (SPEC t o ONCE_REWRITE_RULE[SWAP_FORALL_THM] o
+ ONCE_REWRITE_RULE[CART_EQ])
+ THEN ASM_REWRITE_TAC[CVECTOR_MUL_COMPONENT;CVECTOR_ADD_COMPONENT;cexp_SYM]
+ in
+ let CONJ1_FIRST = MESON[] `(a ==> b/\c) <=> ((a==>b)/\(b==>a==>c))` in
+ [
+ PROJECT_ON `i:num` THEN ASM_MESON_TAC[WAVE_SUM_EQ_CORE];
+ REWRITE_TAC[CONJ1_FIRST] THEN CONJ_TAC
+ THENL [
+ PROJECT_ON `i:num` THEN ASM_MESON_TAC[WAVE_SUM_EQ_WEAK2];
+ DISCH_THEN (fun x -> REWRITE_TAC[x;GSYM CVECTOR_ADD_LDISTRIB])
+ THEN PROJECT_ON `k:num` THEN ASM_MESON_TAC[WAVE_SUM_EQ_SINGLE];
+ ];
+ REWRITE_TAC[MESON[] `(a ==> b/\c) <=> ((a==>c)/\(c==>a==>b))`]
+ THEN CONJ_TAC
+ THENL [
+ PROJECT_ON `j:num` THEN ONCE_REWRITE_TAC[COMPLEX_ADD_SYM]
+ THEN ASM_MESON_TAC[WAVE_SUM_EQ_WEAK1];
+ DISCH_THEN (fun x ->
+ REWRITE_TAC[x;CVECTOR_ARITH `!x:complex^N. x +y=z <=> x=z-y`;
+ GSYM CVECTOR_SUB_LDISTRIB])
+ THEN PROJECT_ON `i:num` THEN ASM_MESON_TAC[WAVE_SUM_EQ_SINGLE];
+ ];
+ ASM_CASES_TAC `(C:complex^N)$i = Cx(&0)`
+ THENL [
+ REWRITE_TAC[CONJ1_FIRST] THEN CONJ_TAC
+ THENL [
+ PROJECT_ON `i:num`
+ THEN REWRITE_TAC[COMPLEX_MUL_LZERO;COMPLEX_LNEG_UNIQ;
+ COMPLEX_NEG_LMUL]
+ THEN ASM_MESON_TAC[WAVE_SUM_EQ_SINGLE];
+ DISCH_THEN (fun x -> REWRITE_TAC[x;GSYM CVECTOR_ADD_LDISTRIB])
+ THEN PROJECT_ON `k:num` THEN ASM_MESON_TAC[WAVE_SUM_EQ_SINGLE];
+ ];
+ ONCE_REWRITE_TAC[MESON[] `a=b/\b=c <=> a=c /\ (a=b/\b=c)`]
+ THEN REWRITE_TAC
+ [MESON[] `(a==>b/\c/\d) <=> (a==>b) /\ (b==>a==>c/\d)`]
+ THEN CONJ_TAC
+ THENL [
+ PROJECT_ON `i:num` THEN ASM_MESON_TAC[WAVE_SUM_EQ_WEAK1];
+ DISCH_THEN (fun x ->
+ REWRITE_TAC[x;CVECTOR_ARITH `!x:complex^N. x +y=z <=> y=z-x`;
+ GSYM CVECTOR_SUB_LDISTRIB])
+ THEN PROJECT_ON `j:num` THEN ASM_MESON_TAC[WAVE_SUM_EQ_SINGLE]
+ ];
+ ];
+ REPLICATE_TAC 2 (DISCH_THEN (SINGLE REWRITE_TAC))
+ THEN REWRITE_TAC
+ [CVECTOR_ARITH `x % X + x % Y = x % Z <=> x % (X+Y-Z) = cvector_zero`;
+ CVECTOR_MUL_EQ_0;CEXP_NZ]
+ THEN CVECTOR_ARITH_TAC
+ ]);;
+
+let VEC_WAVE_SUM_EQ_SINGLE = prove
+ (`!a b A B:complex^N.
+ ~(A = cvector_zero) /\ (!x. cexp (a * Cx x) % A = cexp (b * Cx x) % B)
+ ==> a = b /\ A = B`,
+ REPEAT GEN_TAC THEN STRIP_TAC
+ THEN SUBGOAL_THEN `~(B:complex^N = cvector_zero)` ASSUME_TAC
+ THENL ON_FIRST_GOAL (POP_ASSUM (fun x -> ASM_REWRITE_TAC[
+ REWRITE_RULE[COMPLEX_MUL_RZERO;CEXP_0;CVECTOR_MUL_ID]
+ (SPEC `&0` (GSYM x))]))
+ THEN SUBGOAL_THEN
+ `!x. cexp (a * Cx x) % (A:complex^N) + cexp (a * Cx x) % A
+ = cexp (b * Cx x) % (Cx(&2) % B)`
+ ASSUME_TAC
+ THENL ON_FIRST_GOAL (GEN_TAC THEN POP_ASSUM (K ALL_TAC)
+ THEN POP_ASSUM (MP_TAC o SPEC_ALL) THEN CVECTOR_ARITH_TAC)
+ THEN SUBGOAL_THEN `~(Cx(&2) % (B:complex^N) = cvector_zero)` ASSUME_TAC
+ THENL ON_FIRST_GOAL (POP_ASSUM (K ALL_TAC) THEN POP_ASSUM MP_TAC
+ THEN REWRITE_TAC[CVECTOR_MUL_EQ_0;SIMPLE_COMPLEX_ARITH `~(Cx(&2) = Cx(&0))`])
+ THEN ASM_MESON_TAC
+ [VEC_WAVE_SUM_EQ_CORE;CVECTOR_ARITH `A + A = Cx(&2) % B <=> A = B`]);;
+
+let VEC_WAVE_SUM_EQ_WEAK1 = prove
+ (`!a b c A B C:complex^N.
+ ~(A = cvector_zero) /\ ~(C = cvector_zero)
+ /\ (!x. cexp (a * Cx x) % A + cexp (b * Cx x) % B = cexp (c * Cx x) % C)
+ ==> a = c /\ A + B = C`,
+ REPEAT GEN_TAC THEN ASM_CASES_TAC `(B:complex^N) = cvector_zero` THENL [
+ ASM_REWRITE_TAC[CVECTOR_MUL_RZERO;CVECTOR_ADD_ID]
+ THEN MESON_TAC[VEC_WAVE_SUM_EQ_SINGLE];
+ REPEAT STRIP_TAC
+ THENL ON_FIRST_GOAL (MATCH_MP_TAC EQ_TRANS THEN EXISTS_TAC `b:complex`)
+ THEN ASM_MESON_TAC[VEC_WAVE_SUM_EQ_CORE];
+ ]);;
+
+let VEC_WAVE_SUM_EQ_WEAK2 = prove
+ (`!a b c A B C:complex^N.
+ ~(A = cvector_zero) /\ ~(B = cvector_zero)
+ /\ (!x. cexp (a * Cx x) % A + cexp (b * Cx x) % B = cexp (c * Cx x) % C)
+ ==> a = b /\ A + B = C`,
+ REPEAT GEN_TAC THEN ASM_CASES_TAC `(C:complex^N) = cvector_zero`
+ THENL [
+ ASM_REWRITE_TAC[CVECTOR_MUL_RZERO;
+ CVECTOR_ARITH `x+y = cvector_zero <=> x= --y`;GSYM CVECTOR_MUL_RNEG]
+ THEN MESON_TAC[VEC_WAVE_SUM_EQ_SINGLE];
+ ASM_MESON_TAC[VEC_WAVE_SUM_EQ_CORE];
+ ]);;
diff --git a/make.ml b/make.ml
new file mode 100644
index 0000000..a2a9c9a
--- /dev/null
+++ b/make.ml
@@ -0,0 +1,2 @@
+List.iter needs ["Multivariate/make_complex.ml";"Multivariate/geom.ml"; "Multivariate/cross.ml";"Multivariate/wlog.ml"];;
+needs "top.ml";;
diff --git a/primitive_rules.ml b/primitive_rules.ml
new file mode 100644
index 0000000..04890d9
--- /dev/null
+++ b/primitive_rules.ml
@@ -0,0 +1,754 @@
+(* ========================================================================= *)
+(* Formalization of Electromagnetic Optics *)
+(* *)
+(* (c) Copyright, Sanaz Khan Afshar & Vincent Aravantinos 2011-14 *)
+(* Hardware Verification Group, *)
+(* Concordia University *)
+(* *)
+(* Contact: <s_khanaf@encs.concordia.ca> *)
+(* <vincent@encs.concordia.ca> *)
+(* *)
+(* This file deals with the non-trivial theorems called "primitive rules". *)
+(* ========================================================================= *)
+
+
+let WAVE_ORTHOGONALITY = prove
+ (`!emf. is_plane_wave emf ==> corthogonal (h_of_wave emf) (e_of_wave emf)`,
+ REWRITE_TAC[IS_PLANE_WAVE;is_valid_emf;plane_wave;emf_at_point_mul;e_of_emf;
+ h_of_emf;LET_DEF;LET_END_DEF]
+ THEN REPEAT STRIP_TAC THEN REPLICATE_TAC 2 (POP_ASSUM (K ALL_TAC))
+ THEN POP_ASSUM (fun x -> RULE_ASSUM_TAC (ONCE_REWRITE_RULE[x]))
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[BETA_THM;LET_DEF;LET_END_DEF])
+ THEN REPEAT (POP_ASSUM MP_TAC)
+ THEN SIMP_TAC[corthogonal;CDOT3;CVECTOR_MUL_COMPONENT;CNJ_MUL;
+ SIMPLE_COMPLEX_ARITH `(x*y)*cnj z*t=(x*cnj z)*(y*t):complex`;
+ GSYM COMPLEX_ADD_LDISTRIB;COMPLEX_ENTIRE;CEXP_NZ;CNJ_EQ_0]);;
+
+let NON_NULL_LEMMA = prove
+ (`!emf v. ~(v = vec 0) /\ non_null_wave emf ==>
+ let v_ccross = (ccross) (vector_to_cvector v) in
+ ~(v_ccross (e_of_wave emf) = cvector_zero)
+ \/ ~(v_ccross (h_of_wave emf) = cvector_zero)`,
+ REPEAT GEN_TAC THEN CONV_TAC (DEPTH_CONV let_CONV)
+ THEN ASM_CASES_TAC `vector_to_cvector v ccross e_of_wave emf = cvector_zero`
+ THEN ASM_REWRITE_TAC[non_null_wave]
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[CCROSS_COLLINEAR_CVECTORS])
+ THEN ASSUME_CONSEQUENCES (REWRITE_RULE[IMP_CONJ]
+ CORTHOGONAL_COLLINEAR_CVECTORS)
+ THEN REPEAT STRIP_TAC
+ THEN SUBGOAL_THEN `~(vector_to_cvector (v:real^3) = cvector_zero)` ASSUME_TAC
+ THENL ON_FIRST_GOAL
+ (REWRITE_TAC[CART_EQ3;CVECTOR_ZERO_COMPONENT;VECTOR_TO_CVECTOR_COMPONENT;
+ CX_INJ]
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[CART_EQ;VEC_COMPONENT;DIMINDEX_3;FORALL_3])
+ THEN ASM_REWRITE_TAC[])
+ THEN SUBGOAL_THEN
+ `?a . ~(a=Cx(&0)) /\ vector_to_cvector v = a % e_of_wave emf` STRIP_ASSUME_TAC
+ THENL ON_FIRST_GOAL (ASM_MESON_TAC[NON_NULL_COLLINEARS;COLLINEAR_CVECTORS_SYM])
+ THEN ASSUME_CONSEQUENCES WAVE_ORTHOGONALITY
+ THEN POP_ASSUM MP_TAC
+ THEN POP_ASSUM (fun x ->
+ RULE_ASSUM_TAC (REWRITE_RULE[x;CCROSS_LMUL;CVECTOR_MUL_EQ_0]))
+ THEN POP_ASSUM (fun x ->
+ RULE_ASSUM_TAC (REWRITE_RULE[x;CCROSS_COLLINEAR_CVECTORS]))
+ THEN ASM_MESON_TAC[CORTHOGONAL_COLLINEAR_CVECTORS;CORTHOGONAL_SYM]);;
+
+let NON_NULL_LEMMA_PASTECART = prove
+ (`!emf v. ~(v = vec 0) /\ non_null_wave emf ==>
+ let v_ccross = (ccross) (vector_to_cvector v) in
+ ~(pastecart (v_ccross (e_of_wave emf)) (v_ccross (h_of_wave emf)) =
+ cvector_zero)`,
+ REWRITE_TAC[PASTECART_EQ_CVECTOR_ZERO;DE_MORGAN_THM;NON_NULL_LEMMA]);;
+
+let BOUNDARY_CONDITIONS_FOR_PLANE_WAVES = prove
+ (`!i emf_i emf_r emf_t.
+ is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ let p = plane_of_interface i in
+ !n. is_normal_to_plane n p ==>
+ let n_ccross = (ccross) (vector_to_cvector n) in
+ !pt. pt IN p ==> !t.
+ let plane_component = \f_of_wave emf. cexp (--ii * Cx((k_of_wave emf) dot
+ pt - w_of_wave emf*t)) % n_ccross (f_of_wave emf) in
+ plane_component e_of_wave emf_i + plane_component e_of_wave emf_r
+ = plane_component e_of_wave emf_t
+ /\ plane_component h_of_wave emf_i + plane_component h_of_wave emf_r
+ = plane_component h_of_wave emf_t`,
+ REWRITE_TAC[FORALL_INTERFACE_THM;plane_of_interface;LET_DEF;LET_END_DEF;
+ is_plane_wave_at_interface;non_null_wave;IS_PLANE_WAVE;plane_wave;
+ emf_at_point_mul;e_of_emf;h_of_emf;LET_DEF;LET_END_DEF;map_triple;o_DEF]
+ THEN REPEAT STRIP_TAC
+ THEN POP_ASSUM (fun x -> FIRST_X_ASSUM (fun y -> MP_TAC (MATCH_MP y x)))
+ THEN ASM (GEN_REWRITE_TAC (RATOR_CONV o DEPTH_CONV)) []
+ THEN REWRITE_TAC[boundary_conditions;emf_add;e_of_emf;h_of_emf;LET_DEF;
+ LET_END_DEF;CCROSS_RADD;CCROSS_RMUL]
+ THEN DISCH_THEN (C ASM_CSQ_THEN STRIP_ASSUME_TAC) THEN ASM_REWRITE_TAC[]);;
+
+(* We combine the E and H field as one complex^6 vector.
+ * Convenient for some proofs. *)
+let PASTECART_BOUNDARY_CONDITIONS_FOR_PLANE_WAVES = prove
+ (`!i emf_i emf_r emf_t. is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ let p = plane_of_interface i in
+ !n. is_normal_to_plane n p ==>
+ let n_ccross = (ccross) (vector_to_cvector n) in
+ !pt. pt IN p ==> !t.
+ let plane_component = \emf.
+ cexp (--ii * Cx((k_of_wave emf) dot pt - w_of_wave emf*t)) % pastecart
+ (n_ccross (e_of_wave emf)) (n_ccross (h_of_wave emf)) in
+ plane_component emf_i + plane_component emf_r = plane_component emf_t`,
+ REWRITE_TAC[LET_DEF;LET_END_DEF] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN (REWRITE_RULE[LET_DEF;LET_END_DEF]
+ BOUNDARY_CONDITIONS_FOR_PLANE_WAVES) (C ASM_CSQ_THEN (C ASM_CSQ_THEN
+ (MP_TAC o SPEC_ALL)))
+ THEN PASTECART_TAC[GSYM PASTECART_CVECTOR_MUL;PASTECART_CVECTOR_ADD]);;
+
+let EXISTS_NORMAL_OF_PLANE_INTERFACE = prove
+ (`!i emf_i emf_r emf_t. is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ ?n. is_normal_to_plane n (plane_of_interface i)`,
+ REWRITE_TAC[FORALL_INTERFACE_THM;LET_DEF;LET_END_DEF;
+ is_plane_wave_at_interface;is_valid_interface;plane_of_interface]
+ THEN MESON_TAC[EXISTS_NORMAL_OF_PLANE]);;
+
+let FREQUENCY_CONSERVATION = prove
+ (`!i emf_i emf_r emf_t. is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ (non_null_wave emf_r ==> w_of_wave emf_r = w_of_wave emf_i)
+ /\ (non_null_wave emf_t ==> w_of_wave emf_t = w_of_wave emf_i)`,
+ REWRITE_TAC[FORALL_INTERFACE_THM] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN EXISTS_NORMAL_OF_PLANE_INTERFACE STRIP_ASSUME_TAC
+ THEN SUBGOAL_THEN `~(n' = vec 0:real^3) /\ plane (p:plane)`
+ STRIP_ASSUME_TAC THEN
+ TRY (RULE_ASSUM_TAC (REWRITE_RULE[is_plane_wave_at_interface;LET_DEF;
+ LET_END_DEF;is_valid_interface;is_normal_to_plane])
+ THEN ASM_REWRITE_TAC[] THEN NO_TAC)
+ THEN ASM_CSQ_THEN PLANE_NON_EMPTY (STRIP_ASSUME_TAC o REWRITE_RULE[GSYM
+ MEMBER_NOT_EMPTY])
+ THEN ASM_CSQ_THEN (REWRITE_RULE[LET_DEF;LET_END_DEF]
+ PASTECART_BOUNDARY_CONDITIONS_FOR_PLANE_WAVES) (C ASM_CSQ_THEN ASSUME_TAC)
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[plane_of_interface;LET_DEF;LET_END_DEF])
+ THEN POP_ASSUM (C ASM_CSQ_THEN MP_TAC)
+ THEN REWRITE_TAC[CX_SUB;CX_MUL;SIMPLE_COMPLEX_ARITH
+ `--ii*(x-y*z) = (ii*y)*z+ --ii*x`;CEXP_ADD;GSYM CVECTOR_MUL_ASSOC]
+ THEN ONCE_REWRITE_TAC[MESON[COMPLEX_EQ_MUL_LCANCEL;II_NZ;CX_INJ]
+ `x=y <=> ii * Cx x = ii * Cx y`]
+ THENL
+ (let APPLY_FREQ_EQ_TAC x =
+ DISCH_THEN (MP_TAC o MATCH_MP (REWRITE_RULE[MESON[]
+ `(A /\ B /\ C ==> D) <=> (C ==> A ==> B ==> D)`] x)) in
+ map APPLY_FREQ_EQ_TAC [ VEC_WAVE_SUM_EQ_WEAK2; VEC_WAVE_SUM_EQ_WEAK1])
+ THEN ANTS_TAC
+ THEN
+ let APPLY_NON_NULL_LEMMA = REWRITE_TAC[CVECTOR_MUL_EQ_0;CEXP_NZ]
+ THEN MATCH_MP_TAC (CONV_RULE (DEPTH_CONV let_CONV) NON_NULL_LEMMA_PASTECART) in
+ REPEAT (ANTS_TAC ORELSE APPLY_NON_NULL_LEMMA)
+ THEN ASM_MESON_TAC[is_plane_wave_at_interface]);;
+
+let IS_PLANE_WAVE_AT_INTERFACE_THMS =
+ [is_plane_wave_at_interface;LET_DEF;LET_END_DEF;map_triple;o_DEF;
+ is_valid_interface];;
+
+let K_PLANE_PROJECTION_CONSERVATION = prove
+ (`!i emf_i emf_r emf_t. is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ let n = unit (normal_of_interface i) in
+ (non_null_wave emf_r ==>
+ projection_on_hyperplane (k_of_wave emf_r) n =
+ projection_on_hyperplane (k_of_wave emf_i) n)
+ /\ (non_null_wave emf_t ==>
+ projection_on_hyperplane (k_of_wave emf_t) n =
+ projection_on_hyperplane (k_of_wave emf_i) n)`,
+ REWRITE_TAC[FORALL_INTERFACE_THM] THEN REPEAT GEN_TAC
+ THEN DISCH_THEN (DISTRIB [MP_TAC; LABEL_TAC "*" o MATCH_MP
+ PASTECART_BOUNDARY_CONDITIONS_FOR_PLANE_WAVES])
+ THEN REWRITE_TAC([normal_of_interface] @ IS_PLANE_WAVE_AT_INTERFACE_THMS)
+ THEN REPEAT STRIP_TAC
+ THEN REMOVE_THEN "*" (C ASM_CSQ_THEN ASSUME_TAC o
+ REWRITE_RULE[LET_DEF;LET_END_DEF;plane_of_interface])
+ THEN ASM_CSQ_THEN FORALL_PLANE_THM_2 STRIP_ASSUME_TAC
+ THEN POP_ASSUM (RULE_ASSUM_TAC o SINGLE REWRITE_RULE)
+ THEN FIRST_ASSUM (STRIP_ASSUME_TAC o CONV_RULE
+ (REWR_CONV is_normal_to_plane))
+ THEN ASM_SIMP_TAC[UNIT_THM;PROJECTION_ON_HYPERPLANE_THM]
+ THEN ASM_CSQ_THEN PLANE_DECOMP_DOT (C ASM_CSQ_THEN (C
+ ASM_CSQ_THEN ASSUME_TAC))
+ THEN MAP_EVERY (fun x -> REWRITE_TAC[VECTOR_ARITH x]
+ THEN ASSUM_LIST (GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV)))
+ [`a-b=c-d <=> a=c+b-d:real^N`;`a=c+b-d <=> c=a-b+d:real^N`]
+ THEN REWRITE_TAC[VECTOR_ARITH `a+b+c=(d+e+f)-f+c <=> a+b=d+e:real^N`]
+ THEN ASM_CSQ_THEN BASIS_NON_NULL (fun x -> FIRST_ASSUM
+ (STRIP_ASSUME_TAC o REWRITE_RULE[IN_INSERT;IN_SING;
+ MESON[] `(!x. x=a \/ x=b ==> p x) <=> p a /\ p b`;UNIT_EQ_ZERO]
+ o MATCH_MP x o CONJUNCT1 o CONV_RULE
+ (REWR_CONV is_orthogonal_plane_basis)))
+ THEN ASM_SIMP_TAC[UNIT_UNIT] THEN MK_COMB_TAC
+ THENL [AP_TERM_TAC;ALL_TAC;AP_TERM_TAC;ALL_TAC]
+ THEN AP_THM_TAC THEN AP_TERM_TAC
+ THEN ONCE_REWRITE_TAC[MESON[COMPLEX_EQ_MUL_LCANCEL;II_NZ;CX_INJ;
+ COMPLEX_NEG_EQ_0] `x=y <=> --ii * Cx x = --ii * Cx y`]
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[
+ VECTOR_ARITH `x dot (y+a%z+b%t) = x dot y+(x dot z)*a+(x dot t)*b`;
+ SIMPLE_COMPLEX_ARITH `--ii * Cx ((a+b*b'+c*c')-d) = (--ii * Cx b)
+ * Cx b' + (--ii * Cx c) * Cx c' + --ii * Cx a + ii * Cx d`;
+ CEXP_ADD])
+ THENL (map (fun f -> FIRST_X_ASSUM (ASSUME_TAC o f)) [
+ funpow 2 (SPEC `&0` o ONCE_REWRITE_RULE[SWAP_FORALL_THM]);
+ SPEC `&0` o ONCE_REWRITE_RULE[SWAP_FORALL_THM] o SPEC `&0`;
+ funpow 2 (SPEC `&0` o ONCE_REWRITE_RULE[SWAP_FORALL_THM]);
+ SPEC `&0` o ONCE_REWRITE_RULE[SWAP_FORALL_THM] o SPEC `&0`;
+ ])
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[COMPLEX_MUL_RZERO;CEXP_0;COMPLEX_MUL_RID;
+ REAL_MUL_RZERO;GSYM CVECTOR_MUL_ASSOC;CVECTOR_MUL_ID])
+ THENL
+ (let lemma = MESON[] `(!x y z t u v. p x y z t v u ==> q x y z /\ r t v u)
+ ==> (!x y z t u v. p x y z t v u ==> q x y z)` in
+ let APPLY_FREQ_EQ_TAC x =
+ POP_ASSUM (DISTRIB [ASSUME_TAC;MP_TAC o MATCH_MP (REWRITE_RULE[MESON[]
+ `(A /\ B /\ C ==> D) <=> (C ==> A ==> B ==> D)`](MATCH_MP lemma x))])
+ in
+ map APPLY_FREQ_EQ_TAC [VEC_WAVE_SUM_EQ_WEAK2;VEC_WAVE_SUM_EQ_WEAK2;
+ VEC_WAVE_SUM_EQ_WEAK1;VEC_WAVE_SUM_EQ_WEAK1])
+ THEN ANTS_TAC
+ THEN
+ let APPLY_NON_NULL_LEMMA =
+ REWRITE_TAC[CVECTOR_MUL_EQ_0;CEXP_NZ]
+ THEN MATCH_MP_TAC (CONV_RULE (DEPTH_CONV let_CONV)
+ NON_NULL_LEMMA_PASTECART)
+ in
+ REPEAT (ANTS_TAC ORELSE APPLY_NON_NULL_LEMMA)
+ THEN ASM_REWRITE_TAC[]);;
+
+let LAW_OF_REFLECTION = prove
+ (`!i emf_i emf_r emf_t. is_plane_wave_at_interface i emf_i emf_r emf_t ==>
+ let n = unit (normal_of_interface i) in
+ non_null_wave emf_r ==>
+ symetric_vectors_wrt (--(k_of_wave emf_r)) (k_of_wave emf_i) n`,
+ REWRITE_TAC[FORALL_INTERFACE_THM;LET_DEFs] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN (REWRITE_RULE[LET_DEFs] K_PLANE_PROJECTION_CONSERVATION)
+ (C ASM_CSQ_THEN MP_TAC o CONJUNCT1)
+ THEN ASM_CSQ_THEN (MESON[is_plane_wave_at_interface]
+ `is_plane_wave_at_interface i emf_i emf_r emf_t ==> is_valid_interface i`)
+ (DISTRIB [ASSUME_TAC; STRIP_ASSUME_TAC o
+ REWRITE_RULE[LET_DEFs;is_valid_interface]])
+ THEN ASM_SIMP_TAC[normal_of_interface;LET_DEFs;
+ PROJECTION_ON_HYPERPLANE_THM;symetric_vectors_wrt;
+ NORMAL_OF_INTERFACE_NON_NULL;UNIT_THM;UNIT_EQ_ZERO]
+ THEN DISCH_TAC
+ THEN SUBGOAL_THEN `!x:real^3. orthogonal (x - (x dot unit n) % unit n) ((x
+ dot unit n) % unit n)` ASSUME_TAC
+ THENL ON_FIRST_GOAL (ASM_MESON_TAC [ORTHOGONAL_RUNIT;ORTHOGONAL_CLAUSES;
+ SUB_UNIT_NORMAL_IS_ORTHOGONAL_TO_NORMAL])
+ THEN FIRST_ASSUM (REPEAT_TCL STRIP_THM_THEN (fun x ->
+ if contains_sub_term_name "k_of_wave" (concl x) then ASSUME_TAC x
+ else ALL_TAC)
+ o REWRITE_RULE[LET_DEFs;map_triple;o_DEF] o CONV_RULE (REWR_CONV
+ is_plane_wave_at_interface))
+ THEN REPLICATE_TAC 3 (POP_ASSUM (K ALL_TAC))
+ THEN POP_ASSUM (fun _ -> POP_ASSUM (fun x -> POP_ASSUM (fun y -> MP_TAC
+ (REWRITE_RULE[NORM_EQ] (TRANS x (GSYM y))))))
+ THEN FIRST_ASSUM (ASSUME_TAC o SPEC `k_of_wave emf`)
+ THEN ASM_CSQ_THEN ORTHOGONAL_SQR_NORM (SINGLE REWRITE_TAC)
+ THEN ASM_REWRITE_TAC[REAL_EQ_ADD_LCANCEL;DOT_LMUL;DOT_RMUL;DOT_SQUARE_NORM]
+ THEN SUBGOAL_TAC "" `~(n = vec 0 :real^3)`
+ [ASM_MESON_TAC[is_normal_to_plane]]
+ THEN ASM_SIMP_TAC[UNIT_THM;REAL_POW_2;REAL_MUL_RID]
+ THEN REWRITE_TAC[GSYM REAL_POW_2;GSYM REAL_EQ_SQUARE_ABS;real_abs]
+ THEN REPEAT COND_CASES_TAC THEN RULE_ASSUM_TAC (REWRITE_RULE
+ [real_ge;REWRITE_RULE[real_lt;MESON[] `(~p <=> ~q) <=> (p <=> q)`]
+ DOT_RUNIT_LT0])
+ THENL [
+ FIRST_X_ASSUM (fun x -> FIRST_X_ASSUM (ASSUME_TAC
+ o ONCE_REWRITE_RULE[GSYM DOT_RUNIT_EQ0] o GSYM
+ o CONV_RULE (REWR_CONV REAL_LE_ANTISYM) o CONJ x))
+ THEN ASM_REWRITE_TAC[DOT_LNEG;REAL_NEG_0;REAL_MUL_RZERO;VECTOR_MUL_LZERO;
+ VECTOR_ARITH `--x+y=vec 0 <=> x=y`]
+ THEN DISCH_THEN (RULE_ASSUM_TAC o SINGLE REWRITE_RULE o GSYM)
+ THEN POP_ASSUM (RULE_ASSUM_TAC o SINGLE REWRITE_RULE)
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[VECTOR_MUL_LZERO;VECTOR_SUB_RZERO])
+ THEN ASM_REWRITE_TAC[];
+ ASM_MESON_TAC[REAL_ARITH `x > &0 ==> &0 <= x`];
+ DISCH_THEN (ASSUME_TAC o GSYM)
+ THEN FIRST_X_ASSUM (fun x -> ignore (term_match [] `x-y=z-t` (concl x));
+ MP_TAC x)
+ THEN POP_ASSUM (fun x -> REWRITE_TAC(x::[VECTOR_MUL_LNEG;DOT_LNEG;
+ REAL_MUL_RNEG;VECTOR_ARITH `x-y=z- --y <=> --x+z = --(&2%y)`;
+ VECTOR_MUL_ASSOC]));
+ ASM_MESON_TAC[REAL_ARITH `x > &0 ==> &0 <= x`];
+ ]);;
+
+let PLANE_OF_INCIDENCE_LAW = prove
+ (`!i emf_i emf_r emf_t.
+ is_plane_wave_at_interface i emf_i emf_r emf_t /\ non_null_wave emf_r
+ /\ non_null_wave emf_t ==>
+ coplanar {vec 0, k_of_wave emf_i, k_of_wave emf_r, k_of_wave emf_t,
+ normal_of_interface i}`,
+ REWRITE_TAC[FORALL_INTERFACE_THM;LET_DEFs;normal_of_interface]
+ THEN REPEAT STRIP_TAC
+ THEN REWRITE_TAC[coplanar] THEN REPEAT STRIP_TAC
+ THEN MAP_EVERY EXISTS_TAC [`vec 0:real^3`;`k_of_wave emf_i`;`unit n:real^3`]
+ THEN ASM_CSQ_THEN (MESON[is_plane_wave_at_interface]
+ `is_plane_wave_at_interface i emf_i emf_r emf_t ==> is_valid_interface i`)
+ (DISTRIB [ASSUME_TAC; STRIP_ASSUME_TAC o
+ REWRITE_RULE[LET_DEFs;is_valid_interface]])
+ THEN ASM_CSQ_THEN (REWRITE_RULE[LET_DEFs] K_PLANE_PROJECTION_CONSERVATION)
+ (DISTRIB (map ((o) (C ASM_CSQ_THEN ASSUME_TAC)) [CONJUNCT1;CONJUNCT2]))
+ THEN REWRITE_TAC[INSERT_SUBSET;SING_SUBSET] THEN REPEAT CONJ_TAC
+ THENL
+ let IN_SET_TAC = MATCH_MP_TAC HULL_INC THEN SET_TAC[] in
+ let COMBINATION_TAC =
+ POP_ASSUM MP_TAC THEN POP_ASSUM MP_TAC
+ THEN ASM_SIMP_TAC[AFFINE_HULL_3_ZERO;IN_ELIM_THM;UNIV;normal_of_interface;
+ LET_DEFs;PROJECTION_ON_HYPERPLANE_THM;symetric_vectors_wrt;
+ NORMAL_OF_INTERFACE_NON_NULL;UNIT_THM;UNIT_EQ_ZERO;
+ VECTOR_ARITH `x-a%y=z-b%y <=> x= &1%z+(a-b)%y`]
+ THEN MESON_TAC[]
+ in
+ [ IN_SET_TAC; IN_SET_TAC; COMBINATION_TAC; COMBINATION_TAC;
+ REWRITE_TAC[AFFINE_HULL_3_ZERO;IN_ELIM_THM;UNIV]
+ THEN MAP_EVERY EXISTS_TAC [`&0`;`norm (n:real^3)`]
+ THEN REWRITE_TAC[VECTOR_MUL_LZERO;VECTOR_ADD_LID]
+ THEN MATCH_MP_TAC UNIT_INTRO THEN ASM_MESON_TAC[is_normal_to_plane]]);;
+
+let SNELL_LAW = prove
+ (`!i emf_i emf_r emf_t.
+ is_plane_wave_at_interface i emf_i emf_r emf_t /\ non_null_wave emf_t ==>
+ let theta = \emf. vector_angle (k_of_wave emf) (normal_of_interface i) in
+ n1_of_interface i * sin (theta emf_i) =
+ n2_of_interface i * sin (theta emf_t)`,
+ REWRITE_TAC[FORALL_INTERFACE_THM;LET_DEFs;n1_of_interface;n2_of_interface;
+ normal_of_interface]
+ THEN REPEAT STRIP_TAC
+ THEN FIRST_ASSUM (STRIP_ASSUME_TAC o REWRITE_RULE[LET_DEFs;map_triple;o_DEF;
+ is_valid_interface;IS_PLANE_WAVE;non_null_wave]
+ o CONV_RULE (REWR_CONV is_plane_wave_at_interface))
+ THEN SUBGOAL_TAC ""
+ `~(k_of_wave emf_i = vec 0) /\ ~(k_of_wave emf_t = vec 0)`
+ [ASM_REWRITE_TAC[]]
+ THEN SUBGOAL_TAC "" `~(n = vec 0:real^3)` [ASM_MESON_TAC[is_normal_to_plane]]
+ THEN ASSUM_LIST (SIMP_TAC o (@) [REWRITE_RULE[DE_MORGAN_THM;
+ MESON[] `(p\/q==>r) <=> (p==>r) /\ (q==>r)`] (CONV_RULE (DEPTH_CONV
+ COND_ELIM_CONV) vector_angle)]
+ o filter (fun x -> not (contains_sub_term_name "norm" (concl x))))
+ THEN SUBGOAL_THEN
+ `(k_of_wave emf_i dot unit n) / (k0*n1) =
+ (k_of_wave emf_i dot n) / (norm (k_of_wave emf_i) * norm n)
+ /\ (k_of_wave emf_t dot unit n) / (k0*n2) =
+ (k_of_wave emf_t dot n) / (norm (k_of_wave emf_t) * norm n)`
+ ASSUME_TAC
+ THENL ON_FIRST_GOAL (ASM_REWRITE_TAC[unit;DOT_RMUL;real_div;
+ REAL_ARITH `(x*y)*z=y*(z*x):real`;REAL_INV_MUL])
+ THEN SIMP_TAC[REWRITE_RULE[REAL_ABS_BOUNDS] NORM_CAUCHY_SCHWARZ_DIV;SIN_ACS;
+ GSYM REAL_EQ_SQUARE_ABS]
+ THEN ASM_CSQ_THEN (REWRITE_RULE[LET_DEFs] K_PLANE_PROJECTION_CONSERVATION)
+ (C ASM_CSQ_THEN MP_TAC o CONJUNCT2)
+ THEN ASM_SIMP_TAC[normal_of_interface;LET_DEFs;PROJECTION_ON_HYPERPLANE_THM;
+ NORMAL_OF_INTERFACE_NON_NULL;UNIT_THM;UNIT_EQ_ZERO]
+ THEN DISCH_THEN (MP_TAC o MATCH_MP (MESON[] `x=y ==> x dot x = y dot y`))
+ THEN SUBGOAL_THEN `!x:real^3. orthogonal (x - (x dot unit n) % unit n) ((x
+ dot unit n) % unit n)` ASSUME_TAC
+ THENL ON_FIRST_GOAL (ASM_MESON_TAC[ORTHOGONAL_RUNIT;ORTHOGONAL_CLAUSES;
+ SUB_UNIT_NORMAL_IS_ORTHOGONAL_TO_NORMAL])
+ THEN ASM_SIMP_TAC[REWRITE_RULE[REAL_ARITH `x=y+z <=> y=x-z:real`]
+ ORTHOGONAL_SQR_NORM;DOT_SQUARE_NORM;NORM_MUL;UNIT_THM;REAL_MUL_RID;
+ REAL_POW2_ABS]
+ THEN SUBGOAL_TAC "" `~(k0*n2 = &0) /\ ~(k0*n1 = &0)`
+ [ASM_MESON_TAC[NORM_EQ_0]]
+ THEN SUBGOAL_THEN `&0 < inv(k0*n1)` ASSUME_TAC
+ THENL ON_FIRST_GOAL (ASM_REWRITE_TAC[REAL_LT_INV_EQ;REAL_LT_LE]
+ THEN ASM_MESON_TAC[NORM_POS_LE])
+ THEN SUBGOAL_THEN `~((k0*n2) pow 2 = &0) /\ ~(inv ((k0 * n2) pow 2) = &0)`
+ STRIP_ASSUME_TAC
+ THENL ON_FIRST_GOAL (REWRITE_TAC[REAL_POW_EQ_0;REAL_INV_EQ_0]
+ THEN ASM_ARITH_TAC)
+ THEN POP_ASSUM (fun x -> DISCH_THEN (MP_TAC o
+ ONCE_REWRITE_RULE [REWRITE_RULE[x] (GENL [`x:real`;`y:real`] (SPECL
+ [`x:real`;`y:real`;`inv ((k0 * n2:real) pow 2)`] (GSYM
+ REAL_EQ_MUL_RCANCEL)))]))
+ THEN ASM_SIMP_TAC[REAL_SUB_RDISTRIB;REAL_MUL_RINV;GSYM real_div;
+ GSYM REAL_POW_DIV;REAL_DIV_REFL;REAL_POW_ONE]
+ THEN DISCH_THEN (K ALL_TAC)
+ THEN MATCH_MP_TAC (MESON[SQRT_MUL;POW_2_SQRT;REAL_LE_POW_2]
+ `!x y. &0 <= x /\ &0 <= y /\ &0 <= z /\ &0 <= t /\ sqrt(x pow 2 * y) =
+ sqrt(z pow 2*t) ==> x * sqrt y = z * sqrt t`)
+ THEN REPEAT CONJ_TAC THENL [
+ ASM_REAL_ARITH_TAC;
+ REWRITE_TAC[REAL_SUB_LE;ABS_SQUARE_LE_1]
+ THEN ASM_MESON_TAC[NORM_CAUCHY_SCHWARZ_DIV];
+ ASM_REAL_ARITH_TAC;
+ REWRITE_TAC[REAL_POW_DIV;REAL_ARITH `x/y-z/y=(x-z)/y:real`]
+ THEN MATCH_MP_TAC REAL_LE_DIV THEN CONJ_TAC THENL [
+ REWRITE_TAC[REAL_SUB_LE] THEN MATCH_MP_TAC REAL_POW_LE2 THEN CONJ_TAC
+ THENL [
+ REWRITE_TAC[unit;DOT_RMUL] THEN MATCH_MP_TAC REAL_LE_MUL
+ THEN REWRITE_TAC[REAL_LE_INV_EQ;NORM_POS_LE] THEN ASM_ARITH_TAC;
+ FIRST_ASSUM (SINGLE ONCE_REWRITE_TAC o REWRITE_RULE[GSYM real_div]
+ o MATCH_MP (GSYM REAL_LE_RMUL_EQ))
+ THEN ASM_SIMP_TAC[REAL_DIV_REFL]
+ THEN ASM_MESON_TAC[REWRITE_RULE[REAL_ABS_BOUNDS]
+ NORM_CAUCHY_SCHWARZ_DIV]
+ ];
+ MATCH_ACCEPT_TAC REAL_LE_POW_2;
+ ];
+ AP_TERM_TAC THEN SUBGOAL_TAC "" `~(k0 = &0) /\ ~(n1 = &0) /\ ~(n2 = &0)`
+ [ASM_MESON_TAC[REAL_MUL_LZERO;REAL_MUL_RZERO]]
+ THEN ASM_SIMP_TAC[real_div;REAL_INV_MUL;
+ REAL_ARITH `(x*y)*inv x*inv z=(x*inv x)*y*inv z:real`;REAL_MUL_RINV;
+ REAL_MUL_LID;REAL_SUB_LDISTRIB;REAL_MUL_RID]
+ THEN ASM_SIMP_TAC[GSYM REAL_POW_MUL;
+ REAL_ARITH `x*y*inv x=(x*inv x)*y:real`;
+ REAL_ARITH `x*y*(inv z*inv x)*t=(x*inv x)*y*inv z*t:real`;
+ REAL_ARITH `x*y*inv z*inv x=(x*inv x)*y*inv z:real`;
+ REAL_MUL_RINV;REAL_MUL_LID;REAL_INV_MUL]
+ THEN REWRITE_TAC[REAL_ARITH `x*inv y*inv z=(inv z*x)*inv y:real`;
+ GSYM DOT_RMUL;unit]
+ ]);;
+
+let phase_shift_at_plane = new_definition
+ `phase_shift_at_plane p n emf =
+ k_of_wave emf dot (@a. a % unit n IN p) % unit n`;;
+
+let PLANE_WAVE_WITH_PHASE_SHIFT_AT_PLANE = prove
+ (`!emf. is_plane_wave emf ==> !p. plane p ==> !n. is_normal_to_plane n p ==>
+ !r. r IN p ==> !t.
+ emf r t = cexp (--ii * Cx(projection_on_hyperplane (k_of_wave emf) (unit n)
+ dot r - w_of_wave emf * t + phase_shift_at_plane p n emf))
+ % eh_of_wave emf`,
+ let tmp = DISCH_ALL (GSYM (UNDISCH (ISPECL [`k_of_wave emf`;`unit n:real^3`]
+ PROJECTION_ON_HYPERPLANE_DECOMPOS))) in
+ let tmp' = UNDISCH_ALL (IMP_TRANS (SPEC `n:real^3` (UNDISCH (SPEC `p:plane`
+ NORMAL_OF_PLANE_NON_NULL))) (IMP_TRANS (ISPEC `n:real^3` UNIT_THM) tmp)) in
+ REWRITE_TAC[LET_DEF;LET_END_DEF;IS_PLANE_WAVE;plane_wave]
+ THEN REPEAT STRIP_TAC
+ THEN FIRST_ASSUM (let thm = MESON[] `f=g ==> g x y=z ==> f x y=z` in
+ fun x -> MATCH_MP_TAC (MATCH_MP thm x))
+ THEN REWRITE_TAC[emf_at_point_mul;PAIR_EQ;FUN_EQ_THM;eh_of_wave;e_of_emf;
+ h_of_emf;LET_DEF;LET_END_DEF]
+ THEN REPEAT (FIRST[CONJ_TAC;AP_THM_TAC;AP_TERM_TAC;GEN_TAC])
+ THEN REWRITE_TAC[REAL_ARITH `x-y=(z-y)+t <=> x=t+z:real`]
+ THEN GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV)
+ [ONCE_REWRITE_RULE[VECTOR_ADD_SYM] tmp']
+ THEN REWRITE_TAC[DOT_LADD;REAL_EQ_ADD_RCANCEL;phase_shift_at_plane]
+ THEN ABBREV_TAC `a = @a. a % unit n IN (p:plane)`
+ THEN GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV)
+ [VECTOR_ARITH `v dot r = v dot (r - a % unit n + a % unit n:real^3)`]
+ THEN REWRITE_TAC[DOT_RADD;DOT_RMUL;DOT_LMUL;REAL_ADD_LDISTRIB]
+ THEN ASM_REWRITE_TAC[MATCH_MP (UNIT_DOT_UNIT_SELF)
+ (STRIP_RULE NORMAL_OF_PLANE_NON_NULL);REAL_ARITH `x+y*z* &1=z*y <=> x= &0`]
+ THEN REWRITE_TAC[REAL_ENTIRE;GSYM orthogonal;ORTHOGONAL_LUNIT] THEN DISJ2_TAC
+ THEN ONCE_REWRITE_TAC[ORTHOGONAL_SYM]
+ THEN MATCH_MP_TAC (GENL [`pt1:point`;`pt2:point`]
+ (REWRITE_RULE[GSYM IMP_CONJ] (DISCH_ALL (STRIP_RULE
+ NORMAL_OF_PLANE_IS_ORTHOGONAL_TO_SEGMENT))))
+ THEN ASM_REWRITE_TAC[]
+ THEN EXPAND_TAC "a" THEN SELECT_ELIM_TAC
+ THEN ASM_SIMP_TAC[EXISTS_MULTIPLE_OF_NORMAL_IN_PLANE]);;
+
+let PLANE_WAVE_WITH_PHASE_SHIFT_AT_INTERFACE = prove
+ (`!emf. is_plane_wave emf ==> !i. is_valid_interface i ==> !r.
+ r IN plane_of_interface i ==> !t.
+ let n = normal_of_interface i in
+ emf r t = cexp (--ii * Cx(projection_on_hyperplane (k_of_wave emf) (unit n)
+ dot r - w_of_wave emf * t + phase_shift_at_plane (plane_of_interface i) n emf))
+ % eh_of_wave emf`,
+ REWRITE_TAC[LET_DEFs] THEN REPEAT STRIP_TAC THEN ASM_CSQ_THEN
+ IS_VALID_INTERFACE_IS_NORMAL_TO_PLANE ASSUME_TAC
+ THEN ASM_CSQ_THEN IS_VALID_INTERFACE_PLANE ASSUME_TAC
+ THEN ASM_CSQS_THEN PLANE_WAVE_WITH_PHASE_SHIFT_AT_PLANE ASSUME_TAC
+ THEN ASM_REWRITE_TAC[]);;
+
+let magnitude_shifted_at_plane = new_definition
+ `magnitude_shifted_at_plane p n emf =
+ cexp (--ii * Cx(phase_shift_at_plane p n emf)) % eh_of_wave emf`;;
+
+let E_PRESERVED_IN_TE_MODE = prove
+ (`!i emf_i emf_r emf_t.
+ is_plane_wave_at_interface i emf_i emf_r emf_t /\ non_null_wave emf_r
+ /\ non_null_wave emf_t /\ TE_mode i emf_i emf_r emf_t ==>
+ let magnitude =
+ \emf. FST (magnitude_shifted_at_plane (plane_of_interface i)
+ (normal_of_interface i) emf) cdot (vector_to_cvector (TE_mode_axis
+ i emf_i emf_r emf_t))
+ in
+ magnitude emf_r + magnitude emf_i = magnitude emf_t`,
+ REWRITE_TAC[FORALL_INTERFACE_THM;LET_DEFs;plane_of_interface;
+ normal_of_interface] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN (GEQ_IMP is_plane_wave_at_interface) (STRIP_ASSUME_TAC o
+ REWRITE_RULE[LET_DEFs])
+ THEN FIRST_ASSUM (let tmp = MESON[non_null_wave] `!emf. non_null_wave emf
+ ==> is_plane_wave emf` in ASSUME_TAC o MATCH_MP tmp)
+ THEN ASM_CSQ_THEN (GEN_ALL (MATCH_MP EQ_IMP (SPEC_ALL is_valid_interface)))
+ (STRIP_ASSUME_TAC o REWRITE_RULE[LET_DEFs])
+ THEN ASM_CSQ_THEN PLANE_NON_EMPTY (STRIP_ASSUME_TAC o
+ REWRITE_RULE[GSYM MEMBER_NOT_EMPTY])
+ THEN FIRST_ASSUM (fun x -> FIRST_X_ASSUM (C ASM_CSQ_THEN (MP_TAC o
+ REWRITE_RULE[LET_DEFs] o GEN_REWRITE_RULE (RATOR_CONV o ONCE_DEPTH_CONV)
+ [e_of_emf] o CONJUNCT1 o SPEC `&0`) o REWRITE_RULE[boundary_conditions;
+ emf_add;LET_DEFs] o C MATCH_MP x))
+ THEN ASSUM_LIST (MAP_EVERY (C ASM_CSQS_THEN (SINGLE REWRITE_TAC o CONJUNCT1
+ o SPEC `&0` o REWRITE_RULE[emf_at_point_mul;eh_of_wave;FST;SND;EMF_EQ]))
+ o mapfilter (REWRITE_RULE[FORALL_INTERFACE_THM;LET_DEFs;plane_of_interface;
+ normal_of_interface;]
+ o MATCH_MP PLANE_WAVE_WITH_PHASE_SHIFT_AT_INTERFACE))
+ THEN ASM_CSQ_THEN (REWRITE_RULE[LET_DEFs] K_PLANE_PROJECTION_CONSERVATION)
+ (CONJUNCTS_THEN (C ASM_CSQ_THEN (SINGLE REWRITE_TAC o
+ SIMP_RULE[normal_of_interface;LET_DEFs])))
+ THEN REWRITE_TAC[REAL_ARITH `x - y * &0 = x`;CX_ADD;COMPLEX_ADD_LDISTRIB;
+ CEXP_ADD;GSYM CVECTOR_MUL_ASSOC;GSYM CVECTOR_ADD_LDISTRIB;
+ CCROSS_RMUL;CVECTOR_MUL_LCANCEL;CEXP_NZ;magnitude_shifted_at_plane;
+ eh_of_wave;emf_at_point_mul;GSYM CDOT_LADD]
+ THEN ASM_CSQS_THEN (REWRITE_RULE[IMP_CONJ;LET_DEFs] TE_MODE_PLANEWAVE_PROJ)
+ (GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV) o CONJUNCTS)
+ THEN REWRITE_TAC[CCROSS_RADD;CCROSS_RMUL;CVECTOR_MUL_ASSOC;
+ GSYM CVECTOR_ADD_RDISTRIB;CVECTOR_MUL_RCANCEL;CCROSS_COLLINEAR_CVECTORS;
+ GSYM CDOT_LMUL;GSYM CDOT_LADD;COLLINEAR_CVECTORS_VECTOR_TO_CVECTOR]
+ THEN FIRST_X_ASSUM (STRIP_ASSUME_TAC o REWRITE_RULE[is_mode_axis;
+ normal_of_interface;LET_DEFs] o MATCH_MP (GEQ_IMP TE_MODE_THM))
+ THEN FIRST_ASSUM (ASSUME_TAC o REWRITE_RULE[NORM_EQ_0] o MATCH_MP
+ (REAL_ARITH `!x. x= &1 ==> ~(x= &0)`))
+ THEN FIRST_X_ASSUM (STRIP_ASSUME_TAC o MATCH_MP (GEQ_IMP is_normal_to_plane))
+ THEN ASM_CSQS_THEN (ONCE_REWRITE_RULE[ORTHOGONAL_SYM] (REWRITE_RULE[IMP_CONJ]
+ ORTHOGONAL_NON_COLLINEAR)) (SINGLE REWRITE_TAC)
+ THEN REWRITE_TAC[CVECTOR_ADD_SYM]);;
+
+let H_CROSS_Z_WRT_E_IN_TE_MODE = prove
+ (`!i emf_i emf_r emf_t.
+ is_plane_wave_at_interface i emf_i emf_r emf_t
+ /\ TE_mode i emf_i emf_r emf_t ==>
+ let p = plane_of_interface i in
+ let n = normal_of_interface i in
+ let h_cross_z_wrt_e =
+ \emf. (h_of_wave emf) ccross (vector_to_cvector n) =
+ Cx ((n dot k_of_wave emf) / (eta0 * k0)) % e_of_wave emf
+ in
+ h_cross_z_wrt_e emf_i /\ h_cross_z_wrt_e emf_r /\ h_cross_z_wrt_e emf_t`,
+ SIMP_TAC[FORALL_INTERFACE_THM;LET_DEFs;is_plane_wave_at_interface;
+ plane_of_interface;non_null_wave;normal_of_interface]
+ THEN REPEAT STRIP_TAC THEN ASM_CSQS_THEN (REWRITE_RULE[LET_DEFs;IMP_CONJ]
+ TE_MODE_PLANEWAVE_PROJ) (SINGLE ONCE_REWRITE_TAC)
+ THEN REWRITE_TAC[CCROSS_LMUL;CCROSS_LREAL;CDOT_RREAL;CVECTOR_ADD_RDISTRIB;
+ GSYM CVECTOR_MUL_ASSOC;GSYM VECTOR_TO_CVECTOR_MUL;
+ CVECTOR_IM_VECTOR_TO_CVECTOR_RE_IM;CVECTOR_RE_VECTOR_TO_CVECTOR_RE_IM;
+ CCROSS_LADD;CVECTOR_ADD_LDISTRIB]
+ THEN REWRITE_TAC[REWRITE_RULE[VECTOR_ARITH `--x = y <=> x = --y :real^N`]
+ (ONCE_REWRITE_RULE[CROSS_SKEW] CROSS_LAGRANGE);
+ CVECTOR_RE_VECTOR_TO_CVECTOR;
+ CVECTOR_IM_VECTOR_TO_CVECTOR;DOT_LZERO;VECTOR_MUL_LZERO;VECTOR_SUB_RZERO;
+ VECTOR_NEG_0;VECTOR_TO_CVECTOR_ZERO;CVECTOR_MUL_RZERO;CVECTOR_ADD_ID;
+ DOT_RMUL]
+ THEN ASM_CSQ_THEN TE_MODE_AXIS_ORTHOGONAL_N
+ (SINGLE REWRITE_TAC o REWRITE_RULE[orthogonal;normal_of_interface;LET_DEFs]
+ o ONCE_REWRITE_RULE[ORTHOGONAL_SYM])
+ THEN REWRITE_TAC[REAL_MUL_RZERO;VECTOR_MUL_LZERO;VECTOR_ARITH
+ `--(vec 0 - x) = x`;VECTOR_NEG_0;VECTOR_TO_CVECTOR_ZERO;CVECTOR_MUL_RZERO;
+ CVECTOR_ADD_ID]
+ THEN REWRITE_TAC[MESON[CVECTOR_MUL_ASSOC;COMPLEX_MUL_SYM]
+ `a % ii % v = ii % a % v:complex^N`;GSYM VECTOR_TO_CVECTOR_MUL;CVECTOR_EQ;
+ CVECTOR_RE_VECTOR_TO_CVECTOR_RE_IM;CVECTOR_IM_VECTOR_TO_CVECTOR_RE_IM;
+ VECTOR_MUL_ASSOC]
+ THEN CONJ_TAC THEN REPEAT (AP_THM_TAC ORELSE AP_TERM_TAC)
+ THEN REAL_ARITH_TAC);;
+
+let NON_PROJECTED_E_RELATION_IN_TE_MODE = prove
+ (`!i emf_i emf_r emf_t.
+ is_plane_wave_at_interface i emf_i emf_r emf_t /\ non_null_wave emf_r
+ /\ non_null_wave emf_t /\ TE_mode i emf_i emf_r emf_t ==>
+ let magnitude =
+ \emf. FST (magnitude_shifted_at_plane (plane_of_interface i)
+ (normal_of_interface i) emf) cdot (vector_to_cvector (TE_mode_axis i
+ emf_i emf_r emf_t))
+ in
+ let n = unit (normal_of_interface i) in
+ Cx (n dot k_of_wave emf_i) * (magnitude emf_i - magnitude emf_r) =
+ Cx (n dot k_of_wave emf_t) * magnitude emf_t`,
+ REWRITE_TAC[FORALL_INTERFACE_THM;LET_DEFs;normal_of_interface;
+ plane_of_interface] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN (GEQ_IMP is_plane_wave_at_interface) (STRIP_ASSUME_TAC o
+ REWRITE_RULE[LET_DEFs])
+ THEN FIRST_ASSUM (let tmp = MESON[non_null_wave] `!emf. non_null_wave emf
+ ==> is_plane_wave emf` in ASSUME_TAC o MATCH_MP tmp)
+ THEN ASM_CSQ_THEN (GEN_ALL (MATCH_MP EQ_IMP (SPEC_ALL is_valid_interface)))
+ (STRIP_ASSUME_TAC o REWRITE_RULE[LET_DEFs])
+ THEN ASM_CSQ_THEN PLANE_NON_EMPTY (STRIP_ASSUME_TAC
+ o REWRITE_RULE[GSYM MEMBER_NOT_EMPTY])
+ THEN FIRST_ASSUM (fun x -> FIRST_X_ASSUM (C ASM_CSQ_THEN (MP_TAC o
+ REWRITE_RULE[LET_DEFs]
+ o GEN_REWRITE_RULE (RATOR_CONV o ONCE_DEPTH_CONV) [h_of_emf]
+ o CONJUNCT2 o SPEC `&0`)
+ o REWRITE_RULE[boundary_conditions;emf_add;LET_DEFs] o C MATCH_MP x))
+ THEN ASSUM_LIST (MAP_EVERY (C ASM_CSQS_THEN (SINGLE REWRITE_TAC o CONJUNCT2
+ o SPEC `&0` o REWRITE_RULE[emf_at_point_mul;eh_of_wave;FST;SND;EMF_EQ]))
+ o mapfilter (REWRITE_RULE[FORALL_INTERFACE_THM;LET_DEFs;plane_of_interface;
+ normal_of_interface;]
+ o MATCH_MP PLANE_WAVE_WITH_PHASE_SHIFT_AT_INTERFACE))
+ THEN ASM_CSQ_THEN (REWRITE_RULE[LET_DEFs] K_PLANE_PROJECTION_CONSERVATION)
+ (CONJUNCTS_THEN (C ASM_CSQ_THEN (SINGLE REWRITE_TAC
+ o SIMP_RULE[normal_of_interface;LET_DEFs])))
+ THEN REWRITE_TAC[REAL_ARITH `x - y * &0 = x`;CX_ADD;COMPLEX_ADD_LDISTRIB;
+ CEXP_ADD;GSYM CVECTOR_MUL_ASSOC;GSYM CVECTOR_ADD_LDISTRIB;
+ CCROSS_RMUL;CVECTOR_MUL_LCANCEL;CEXP_NZ;magnitude_shifted_at_plane;
+ eh_of_wave;emf_at_point_mul;GSYM CDOT_LADD;CCROSS_RADD]
+ THEN ASM_CSQS_THEN (REWRITE_RULE[IMP_CONJ;CVECTOR_ARITH
+ `--x=a%y <=> x = (--a)%y:complex^N`] (ONCE_REWRITE_RULE[CCROSS_SKEW]
+ (REWRITE_RULE[LET_DEFs] H_CROSS_Z_WRT_E_IN_TE_MODE)))
+ (SINGLE REWRITE_TAC
+ o REWRITE_RULE[LET_DEFs;plane_of_interface;normal_of_interface])
+ THEN REWRITE_TAC[GSYM CX_NEG;real_div;REAL_NEG_LMUL;GSYM DOT_RNEG]
+ THEN FIRST_ASSUM (ASSUME_TAC o ONCE_REWRITE_RULE[DOT_SYM]
+ o REWRITE_RULE[DOT_RUNIT_EQ] o MATCH_MP SYMETRIC_VECTORS_PROJECTION_ON_AXIS
+ o MATCH_MP UNIT_THM o CONJUNCT1 o MATCH_MP (GEQ_IMP is_normal_to_plane))
+ THEN ASM_CSQS_THEN (SIMP_RULE[LET_DEFs] LAW_OF_REFLECTION) (fun x ->
+ POP_ASSUM (SINGLE REWRITE_TAC o C MATCH_MP (REWRITE_RULE
+ [normal_of_interface;LET_DEFs] x)))
+ THEN REWRITE_TAC[CX_MUL;DOT_RNEG;CX_NEG;COMPLEX_MUL_LNEG;CVECTOR_ARITH
+ `a%(--(u*v))%x+c%(u*v)%y = d%(--(u'*v))%z
+ <=> v%u%(a%x-c%y) = v%u'%d%z:complex^N`]
+ THEN REWRITE_TAC[CVECTOR_MUL_LCANCEL;CX_INJ;REAL_INV_EQ_0;REAL_ENTIRE;
+ MATCH_MP REAL_LT_IMP_NZ eta0;MATCH_MP REAL_LT_IMP_NZ k0]
+ THEN REWRITE_TAC[unit;DOT_LMUL;CX_MUL;GSYM COMPLEX_MUL_ASSOC;
+ COMPLEX_EQ_MUL_LCANCEL;CX_INJ;REAL_ENTIRE;REAL_INV_EQ_0;NORM_EQ_0]
+ THEN ASM_CSQS_THEN NORMAL_OF_PLANE_NON_NULL (SINGLE REWRITE_TAC)
+ THEN ASM_CSQS_THEN (REWRITE_RULE[IMP_CONJ;LET_DEFs] TE_MODE_PLANEWAVE_PROJ)
+ (GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV) o CONJUNCTS)
+ THEN REWRITE_TAC[CVECTOR_MUL_ASSOC;GSYM CVECTOR_SUB_RDISTRIB;
+ CVECTOR_MUL_RCANCEL;VECTOR_TO_CVECTOR_ZERO_EQ;GSYM NORM_EQ_0]
+ THEN FIRST_ASSUM (SINGLE REWRITE_TAC o MATCH_MP (REWRITE_RULE[is_mode_axis]
+ (GEQ_IMP TE_MODE_THM)))
+ THEN SIMP_TAC[REAL_ARITH `~(&1 = &0)`;CDOT_LMUL;COMPLEX_MUL_ASSOC]);;
+
+let COMPLEX_MUL_LINV2 = prove
+ (`!x y. ~(x=Cx(&0)) ==> inv x * (x * y) = y`,
+ SIMP_TAC[COMPLEX_MUL_ASSOC;COMPLEX_MUL_LINV;COMPLEX_MUL_LID]);;
+
+let COMPLEX_BALANCE_MUL_DIV = prove
+ (`!x y z. ~(x=Cx(&0)) ==> (x*y=z <=> y=z/x)`,
+ REPEAT STRIP_TAC THEN EQ_TAC THENL [
+ DISCH_THEN (fun x -> ASM_SIMP_TAC[GSYM x;complex_div;SIMPLE_COMPLEX_ARITH
+ `(x*y)*inv x=(x*inv x)*y:complex`;COMPLEX_MUL_RINV;COMPLEX_MUL_LID]);
+ DISCH_THEN (fun x -> ASM_SIMP_TAC[x;COMPLEX_DIV_LMUL]);
+ ]);;
+
+
+let FRESNEL_REFLECTION_TE_MODE = prove
+ (`!i emf_i emf_r emf_t.
+ is_plane_wave_at_interface i emf_i emf_r emf_t /\ non_null_wave emf_r
+ /\ non_null_wave emf_t /\ TE_mode i emf_i emf_r emf_t ==>
+ let magnitude =
+ \emf. FST (magnitude_shifted_at_plane (plane_of_interface i)
+ (normal_of_interface i) emf) cdot (vector_to_cvector (TE_mode_axis i
+ emf_i emf_r emf_t))
+ in
+ let theta = \emf. Cx(vector_angle (k_of_wave emf) (normal_of_interface i)) in
+ let n1 = Cx(n1_of_interface i) in
+ let n2 = Cx(n2_of_interface i) in
+ magnitude emf_r = (n1 * ccos (theta emf_i) - n2 * ccos (theta emf_t)) / (n1 *
+ ccos (theta emf_i) + n2 * ccos (theta emf_t)) * magnitude emf_i`,
+ REPEAT STRIP_TAC THEN ASM_CSQS_THEN (REWRITE_RULE[IMP_CONJ]
+ NON_PROJECTED_E_RELATION_IN_TE_MODE) MP_TAC
+ THEN ASM_CSQS_THEN (REWRITE_RULE[IMP_CONJ] E_PRESERVED_IN_TE_MODE) (MP_TAC o
+ GSYM)
+ THEN LET_TAC THEN REWRITE_TAC[LET_DEFs]
+ THEN DISCH_THEN (fun x -> REWRITE_TAC[x;SIMPLE_COMPLEX_ARITH
+ `x*(y-z)=t*(z+y) <=> (x+t)*z=(x-t)*y:complex`])
+ THEN DISCH_THEN (MP_TAC o MATCH_MP (MATCH_MP (MESON[]
+ `(!x y z. p x ==> (q x y z <=> r x y z)) ==> !x y z. q x y z ==> p x
+ ==> r x y z`) COMPLEX_BALANCE_MUL_DIV))
+ THEN ANTS_TAC THENL [
+ REWRITE_TAC[GSYM CX_ADD;CX_INJ] THEN MATCH_MP_TAC REAL_POS_NZ
+ THEN MATCH_MP_TAC REAL_LTE_ADD THEN REWRITE_TAC[DOT_LUNIT_POS;DOT_LUNIT_GE0]
+ THEN ONCE_REWRITE_TAC[DOT_SYM] THEN REPEAT (POP_ASSUM MP_TAC)
+ THEN SPEC_TAC (`i:interface`,`i:interface`)
+ THEN SIMP_TAC[FORALL_INTERFACE_THM;is_plane_wave_at_interface;
+ normal_of_interface;map_triple;LET_DEFs;real_ge;real_gt];
+ SIMP_TAC[SIMPLE_COMPLEX_ARITH `(x*y)/z = (x/z)*y`]
+ THEN DISCH_THEN (K ALL_TAC) THEN AP_THM_TAC THEN AP_TERM_TAC
+ THEN REWRITE_TAC[GSYM CX_ADD;GSYM CX_SUB;GSYM CX_COS;GSYM CX_MUL;
+ GSYM CX_DIV;CX_INJ]
+ THEN REPEAT (POP_ASSUM MP_TAC) THEN SPEC_TAC (`i:interface`,`i:interface`)
+ THEN REWRITE_TAC[FORALL_INTERFACE_THM;normal_of_interface;
+ plane_of_interface;LET_DEFs;n1_of_interface;n2_of_interface]
+ THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN (GEQ_IMP is_plane_wave_at_interface) (STRIP_ASSUME_TAC
+ o REWRITE_RULE[is_valid_interface;is_normal_to_plane;LET_DEFs;
+ map_triple])
+ THEN ASM_REWRITE_TAC[GSYM NORM_EQ_0]
+ THEN ASM_SIMP_TAC [GEN_ALL (DISCH_ALL (GSYM (MATCH_MP REAL_LT_IMP_NE
+ (UNDISCH_ALL (SPEC_ALL (MATCH_MP (REWRITE_RULE[IMP_CONJ] REAL_LT_MUL)
+ k0))))))]
+ THEN GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV) [VECTOR_ANGLE]
+ THEN ASM_SIMP_TAC[UNIT_THM;REAL_MUL_LID]
+ THEN REWRITE_TAC[GSYM REAL_MUL_ASSOC;GSYM REAL_ADD_LDISTRIB;
+ GSYM REAL_SUB_LDISTRIB;real_div;REAL_INV_MUL]
+ THEN REWRITE_TAC[REAL_ARITH `x*y*inv x*z=(x*inv x)*y*z:real`;
+ MATCH_MP REAL_MUL_RINV (GSYM (MATCH_MP REAL_LT_IMP_NE k0));REAL_MUL_LID]
+ THEN ASM_SIMP_TAC[unit;VECTOR_ANGLE_LMUL;REAL_INV_EQ_0;NORM_EQ_0;
+ REAL_LE_INV_EQ;NORM_POS_LE;VECTOR_ANGLE_SYM]]);;
+
+let FRESNEL_TRANSMISSION_TE_MODE = prove
+ (`!i emf_i emf_r emf_t.
+ is_plane_wave_at_interface i emf_i emf_r emf_t /\ non_null_wave emf_r
+ /\ non_null_wave emf_t /\ TE_mode i emf_i emf_r emf_t ==>
+ let magnitude =
+ \emf. FST (magnitude_shifted_at_plane (plane_of_interface i)
+ (normal_of_interface i) emf) cdot (vector_to_cvector (TE_mode_axis i
+ emf_i emf_r emf_t))
+ in
+ let theta =
+ \emf. Cx(vector_angle (k_of_wave emf) (normal_of_interface i))
+ in
+ let n1 = Cx(n1_of_interface i) in
+ let n2 = Cx(n2_of_interface i) in
+ magnitude emf_t = Cx(&2) * n1 * ccos (theta emf_i) / (n1 * ccos
+ (theta emf_i) + n2 * ccos (theta emf_t)) * magnitude emf_i`,
+ REPEAT STRIP_TAC THEN ASM_CSQS_THEN (REWRITE_RULE[IMP_CONJ]
+ E_PRESERVED_IN_TE_MODE) (MP_TAC o GSYM)
+ THEN ASM_CSQS_THEN (REWRITE_RULE[IMP_CONJ]
+ NON_PROJECTED_E_RELATION_IN_TE_MODE) MP_TAC THEN LET_TAC
+ THEN REWRITE_TAC[LET_DEFs]
+ THEN REWRITE_TAC[SIMPLE_COMPLEX_ARITH `x*(y-z)=t*u <=> x*z=x*y-t*u:complex`]
+ THEN DISCH_THEN (MP_TAC o MATCH_MP (MATCH_MP (MESON[]
+ `(!x y z. p x ==> (q x y z <=> r x y z)) ==> !x y z. q x y z ==> p x
+ ==> r x y z`) COMPLEX_BALANCE_MUL_DIV))
+ THEN MATCH_MP_TAC (MESON[] `p/\(p==>q==>r==>s) ==> ((p==>q)==>r==>s)`)
+ THEN CONJ_TAC THENL [
+ REWRITE_TAC[CX_INJ] THEN MATCH_MP_TAC REAL_POS_NZ
+ THEN REWRITE_TAC[DOT_LUNIT_POS;DOT_LUNIT_GE0]
+ THEN ONCE_REWRITE_TAC[DOT_SYM] THEN REPEAT (POP_ASSUM MP_TAC)
+ THEN SPEC_TAC (`i:interface`,`i:interface`)
+ THEN SIMP_TAC[FORALL_INTERFACE_THM;is_plane_wave_at_interface;
+ normal_of_interface;map_triple;LET_DEFs;real_gt];
+ DISCH_TAC THEN DISCH_THEN (fun x -> ASM_SIMP_TAC[x;SIMPLE_COMPLEX_ARITH
+ `t = (x*i-y*t)/x+i <=> (y/x+Cx(&1))*t = (x/x+Cx(&1))*i:complex`;
+ COMPLEX_DIV_REFL;SIMPLE_COMPLEX_ARITH `Cx(&1)+Cx(&1)=Cx(&2)`])
+ THEN ASM_CSQ_THEN COMPLEX_DIV_REFL (SINGLE REWRITE_TAC o GSYM)
+ THEN REWRITE_TAC[GSYM CX_ADD;GSYM CX_DIV;real_div;GSYM REAL_ADD_RDISTRIB]
+ THEN DISCH_THEN (MP_TAC o MATCH_MP (MATCH_MP (MESON[]
+ `(!x y z. p x ==> (q x y z <=> r x y z)) ==> (!x y z. q x y z ==> p x
+ ==> r x y z)`) COMPLEX_BALANCE_MUL_DIV))
+ THEN ANTS_TAC THENL [
+ ASM_REWRITE_TAC[CX_MUL;CX_INV;COMPLEX_ENTIRE;COMPLEX_INV_EQ_0]
+ THEN REWRITE_TAC[CX_INJ] THEN MATCH_MP_TAC REAL_POS_NZ
+ THEN MATCH_MP_TAC REAL_LET_ADD
+ THEN REWRITE_TAC[DOT_LUNIT_POS;DOT_LUNIT_GE0]
+ THEN ONCE_REWRITE_TAC[DOT_SYM] THEN REPEAT (POP_ASSUM MP_TAC)
+ THEN SPEC_TAC (`i:interface`,`i:interface`)
+ THEN SIMP_TAC[FORALL_INTERFACE_THM;is_plane_wave_at_interface;
+ normal_of_interface;map_triple;LET_DEFs;real_ge;real_gt];
+ SIMP_TAC[SIMPLE_COMPLEX_ARITH `(x*y)/z = x/z*y:complex`;GSYM CX_DIV;
+ real_div;REAL_INV_MUL;REAL_INV_INV] THEN DISCH_THEN (K ALL_TAC)
+ THEN REWRITE_TAC[GSYM CX_MUL;GSYM CX_COS;GSYM CX_ADD;GSYM CX_DIV;
+ COMPLEX_MUL_ASSOC;CX_INJ]
+ THEN AP_THM_TAC THEN AP_TERM_TAC
+ THEN REWRITE_TAC[CX_INJ;REAL_ARITH `&2*x = (&2*y)*z <=> x=y*z`]
+ THEN REPEAT (POP_ASSUM MP_TAC) THEN SPEC_TAC (`i:interface`,`i:interface`)
+ THEN REWRITE_TAC[FORALL_INTERFACE_THM;normal_of_interface;
+ plane_of_interface;LET_DEFs;n1_of_interface;n2_of_interface]
+ THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN (GEQ_IMP is_plane_wave_at_interface) (STRIP_ASSUME_TAC
+ o REWRITE_RULE[is_valid_interface;is_normal_to_plane;LET_DEFs;map_triple])
+ THEN GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV) [VECTOR_ANGLE]
+ THEN ASM_SIMP_TAC[UNIT_THM;REAL_MUL_LID]
+ THEN REWRITE_TAC[GSYM REAL_MUL_ASSOC;GSYM REAL_ADD_LDISTRIB;
+ GSYM REAL_SUB_LDISTRIB;real_div;REAL_INV_MUL]
+ THEN REWRITE_TAC[REAL_ARITH `inv x*y*x*z=(x*inv x)*y*z:real`;
+ MATCH_MP REAL_MUL_RINV (GSYM (MATCH_MP REAL_LT_IMP_NE k0));
+ REAL_MUL_LID]
+ THEN ASM_SIMP_TAC[unit;VECTOR_ANGLE_LMUL;REAL_INV_EQ_0;NORM_EQ_0;
+ REAL_LE_INV_EQ;NORM_POS_LE]
+ THEN GEN_REWRITE_TAC (RAND_CONV o ONCE_DEPTH_CONV) [VECTOR_ANGLE_SYM]
+ THEN REAL_ARITH_TAC]]);;
diff --git a/tacticlib.ml b/tacticlib.ml
new file mode 100644
index 0000000..ebfefc0
--- /dev/null
+++ b/tacticlib.ml
@@ -0,0 +1,81 @@
+(* ========================================================================= *)
+(* Formalization of Electromagnetic Optics *)
+(* *)
+(* (c) Copyright, Sanaz Khan Afshar & Vincent Aravantinos 2011-14 *)
+(* Hardware Verification Group, *)
+(* Concordia University *)
+(* *)
+(* Contact: <s_khanaf@encs.concordia.ca> *)
+(* <vincent@encs.concordia.ca> *)
+(* *)
+(* This file deals with various utilities. *)
+(* ========================================================================= *)
+
+Format.set_margin 154;;
+
+let CHANGED_RULE r thm =
+ let thm' = r thm in
+ if equals_thm thm thm' then failwith "CHANGED_RULE" else thm';;
+let MAP_ASSUMPTIONS f = REPEAT (POP_ASSUM (MP_TAC o f)) THEN REPEAT DISCH_TAC;;
+let REMOVE_LETS =
+ REPEAT LET_TAC THEN MAP_ASSUMPTIONS (repeat (CONV_RULE let_CONV));;
+let REWRITE_ASSUMPTIONS xs = MAP_ASSUMPTIONS (REWRITE_RULE xs);;
+let REWRITE_EVERYWHERE xs =
+ REWRITE_TAC xs THEN MAP_ASSUMPTIONS (REWRITE_RULE xs);;
+let ON_FIRST_GOAL ?(out_of=2) x =
+ let rec ( * ) n x = if n<=1 then [] else x :: (n-1) * x in
+ x :: out_of * ALL_TAC;;
+let STRIP_ASSUMPTIONS = REPEAT (POP_ASSUM MP_TAC) THEN REPEAT STRIP_TAC;;
+let DESTRUCT_PAIR_TAC p =
+ CHOOSE_THEN (CHOOSE_THEN (fun x -> PURE_REWRITE_TAC[x]))
+ (ISPEC p PAIR_SURJECTIVE);;
+let MAP_ANTECEDENT f = DISCH_THEN (MP_TAC o f);;
+
+let exn_to_bool f x = try ignore (f x); true with _ -> false;;
+
+let strip_parentheses s =
+ if s.[0] = '(' && s.[String.length s-1] = ')'
+ then String.sub s 1 (String.length s-2)
+ else s;;
+let contains_sub_term_name name t =
+ exn_to_bool (find_term (fun x -> ((=) name o strip_parentheses o
+ string_of_term) x)) t;;
+
+let REAL_SIMP_TAC thms =
+ ASM_REWRITE_TAC ([REAL_MUL_LZERO;REAL_MUL_RZERO;REAL_MUL_RID;REAL_MUL_LID;
+ REAL_ADD_LID;REAL_ADD_RID;REAL_SUB_RZERO;REAL_NEG_0] @ thms);;
+let COMPLEX_SIMPLIFY_TAC =
+ ASM_REWRITE_TAC[COMPLEX_ADD_LID;COMPLEX_MUL_LID;COMPLEX_MUL_LZERO;
+ COMPLEX_MUL_RZERO;COMPLEX_SUB_LZERO;COMPLEX_SUB_RZERO;COMPLEX_NEG_0 ];;
+
+let try_or_id f x = try f x with _ -> x;;
+
+let SINGLE f x = f [x];;
+let ASSUME_CONSEQUENCES x =
+ ASM (MAP_EVERY (fun y -> try ASSUME_TAC (MATCH_MP x y) with _ -> ALL_TAC))
+ [];;
+let STRIP_RULE =
+ repeat (CHANGED_RULE UNDISCH_ALL o try_or_id SPEC_ALL)
+ o REWRITE_RULE[IMP_CONJ];;
+
+let CSQS xs x = map (try_or_id (MATCH_MP x)) xs;;
+let CSQS_THEN xs x ttac = MAP_EVERY ttac (CSQS xs x);;
+let ASM_CSQS_THEN x = ASSUM_LIST o C (C CSQS_THEN x);;
+
+let ASM_CSQ_THEN ?(remove=false) ?(match_=true) x ttac =
+ (if remove then FIRST_X_ASSUM else FIRST_ASSUM)
+ (ttac o (if match_ then MATCH_MP else MP) x);;
+
+let ASM_CSQS_THEN x ttac =
+ (* looks absurd, eh? But needed in order to control the execution flow *)
+ let ttac x y = ttac x y in
+ REPEAT_TCL (fun y z -> ASM_CSQ_THEN z y ORELSE ttac z) ttac x;;
+
+let LET_SIMP_TAC = REWRITE_TAC[LET_DEF;LET_END_DEF];;
+
+let distrib fs x = map (fun f -> f x) fs;;
+
+let DISTRIB ttacs x = EVERY (distrib ttacs x);;
+let LET_DEFs = CONJ LET_DEF LET_END_DEF;;
+
+let GEQ_IMP x = GEN_ALL (MATCH_MP EQ_IMP (SPEC_ALL x));;
diff --git a/top.ml b/top.ml
new file mode 100644
index 0000000..f4092e2
--- /dev/null
+++ b/top.ml
@@ -0,0 +1,11 @@
+
+needs "tacticlib.ml";; (** some general Tactics; developed to facilitate the process of formalization **)
+needs "cvectors.ml";; (** Complex Vector Analysis **)
+needs "vectors_ext.ml";; (** some notion in real vector analysis we need in developing EM Optics, e.g., unit vector and projection on hyper plains. **)
+needs "frequency_equalities.ml";; (** Proving that (!x. A*e^(iax)+B*e^(ibx)=C*e^(icx)) ==> a=b=c /\ A+B=C **)
+needs "em_model.ml";; (** Formalizaing electromagnetic models, e.g., EM fields, plane waves, interfaces. **)
+needs "primitive_rules.ml";; (** Formalizing primitive rules of Optics, e.g., Snell's law and Fresnel Equations **)
+
+
+
+
diff --git a/vectors_ext.ml b/vectors_ext.ml
new file mode 100644
index 0000000..9fa6b5f
--- /dev/null
+++ b/vectors_ext.ml
@@ -0,0 +1,835 @@
+(* ========================================================================= *)
+(* Formalization of Electromagnetic Optics *)
+(* *)
+(* (c) Copyright, Sanaz Khan Afshar & Vincent Aravantinos 2011-14 *)
+(* Hardware Verification Group, *)
+(* Concordia University *)
+(* *)
+(* Contact: <s_khanaf@encs.concordia.ca> *)
+(* <vincent@encs.concordia.ca> *)
+(* *)
+(* Extentions of the vector library. *)
+(* ========================================================================= *)
+
+
+prioritize_vector ();;
+
+(** Additions to euclidean space *)
+
+let [ORTHOGONAL_RZERO;ORTHOGONAL_RMUL;ORTHOGONAL_RNEG;ORTHOGONAL_RADD;
+ ORTHOGONAL_RSUB;ORTHOGONAL_LZERO;ORTHOGONAL_LMUL;ORTHOGONAL_LNEG;
+ ORTHOGONAL_LADD;ORTHOGONAL_LSUB] =
+ CONJUNCTS ORTHOGONAL_CLAUSES;;
+
+let NON_VEC0 = prove
+ (`!x:real^N. ~(x=vec 0) <=> ?i. 1 <= i /\ i <= dimindex(:N) /\ ~(x$i = &0)`,
+ GEN_TAC THEN EQ_TAC THEN SIMP_TAC[vec;CART_EQ;LAMBDA_BETA]
+ THEN ASM_MESON_TAC[]);;
+
+let NON_VEC0_DIM3 = prove
+ (`!x:real^3. ~(x=vec 0) <=> ~(x$1= &0) \/ ~(x$2= &0) \/ ~(x$3= &0)`,
+ REWRITE_TAC[NON_VEC0;DIMINDEX_3]
+ THEN MESON_TAC[ARITH_RULE `1 <= i /\ i <= 3 <=> i = 1 \/ i = 2 \/ i = 3`]);;
+
+let ORTHOGONAL_SQR_NORM = prove
+ (`!x y:real^N. orthogonal (x-y) y ==> x dot x = (x-y) dot (x-y) + y dot y`,
+ REWRITE_TAC[orthogonal] THEN REPEAT STRIP_TAC
+ THEN GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV)
+ [VECTOR_ARITH `x dot x = (x-y) dot (x-y) + y dot y + &2 * (x-y) dot y`]
+ THEN ASM_REWRITE_TAC[REAL_MUL_RZERO;REAL_ADD_RID]);;
+
+let COLLINEAR_VEC0 = prove
+ (`!s:real^N->bool. collinear s /\ vec 0 IN s ==> ?u. !x. x IN s
+ ==> ?c. x = c %u`,
+ REWRITE_TAC[collinear] THEN REPEAT (STRIP_TAC ORELSE EXISTS_TAC `u:real^N`)
+ THEN ASM_MESON_TAC[VECTOR_SUB_RZERO]);;
+
+let CROSS_NORM_ORTHOGONAL = prove
+ (`!x y. orthogonal x y ==> norm (x cross y) = norm x * norm y`,
+ ONCE_REWRITE_TAC[GSYM (REWRITE_RULE[GSYM REAL_ABS_REFL] NORM_POS_LE)]
+ THEN SIMP_TAC[GSYM REAL_ABS_MUL;REAL_EQ_SQUARE_ABS;GSYM NORM_CROSS_DOT;
+ orthogonal;REAL_POW_2]
+ THEN REAL_ARITH_TAC);;
+
+
+let COLLINEAR_DEPENDENT = prove
+ (`!x y:real^N. ~(x=y) /\ collinear {vec 0,x,y} ==> dependent {x,y}`,
+ REWRITE_TAC[collinear;dependent;IN_INSERT;IN_SING] THEN REPEAT STRIP_TAC
+ THEN ASM_CASES_TAC `x=vec 0:real^N`
+ THENL ON_FIRST_GOAL (ASM_MESON_TAC[SPAN_0])
+ THEN ASM_CASES_TAC `y=vec 0:real^N`
+ THENL ON_FIRST_GOAL (ASM_MESON_TAC[SPAN_0])
+ THEN SUBGOAL_THEN `?c. x = c % y :real^N` STRIP_ASSUME_TAC THENL [
+ FIRST_ASSUM (fun x -> MAP_EVERY (STRIP_ASSUME_TAC o
+ REWRITE_RULE[IN_INSERT;IN_SING;VECTOR_SUB_RZERO] o C SPECL x) [
+ [`x:real^N`;`vec 0:real^N`]; [`y:real^N`;`vec 0:real^N`]])
+ THEN ASM_REWRITE_TAC[VECTOR_MUL_ASSOC;VECTOR_MUL_RCANCEL]
+ THEN SUBGOAL_THEN `~(c'= &0)` ASSUME_TAC
+ THENL [ ASM_MESON_TAC[VECTOR_MUL_LZERO]; EXISTS_TAC `c/c':real`
+ THEN ASM_SIMP_TAC[REAL_DIV_RMUL] ];
+ SUBGOAL_TAC "" `~(y=c%y:real^N)` [ASM_MESON_TAC[]]
+ THEN EXISTS_TAC `c%y:real^N`
+ THEN ASM_REWRITE_TAC[DELETE_INSERT;EMPTY_DELETE]
+ THEN MESON_TAC[SUBSPACE_MUL;SUBSET;SPAN_INC;IN_SING;SUBSPACE_SPAN]
+ ]);;
+
+let ORTHOGONAL_DIFFERENT = prove
+ (`!x y:real^N. orthogonal x y /\ ~(x=vec 0) /\ ~(y=vec 0) ==> ~(x=y)`,
+ REWRITE_TAC[orthogonal] THEN MESON_TAC[DOT_EQ_0]);;
+
+let ORTHOGONAL_NON_COLLINEAR = prove
+ (`!x y:real^N. orthogonal x y /\ ~(x=vec 0) /\ ~(y=vec 0)
+ ==> ~(collinear {vec 0,x,y})`,
+ REPEAT STRIP_TAC
+ THEN SUBGOAL_THEN `independent {x,y:real^N}` ASSUME_TAC THENL [
+ MATCH_MP_TAC PAIRWISE_ORTHOGONAL_INDEPENDENT
+ THEN REWRITE_TAC[pairwise;IN_INSERT;IN_SING]
+ THEN ASM_MESON_TAC[ORTHOGONAL_SYM];
+ POP_ASSUM (MP_TAC o REWRITE_RULE[independent])
+ THEN REWRITE_TAC[] THEN MATCH_MP_TAC COLLINEAR_DEPENDENT
+ THEN ASM_SIMP_TAC[ORTHOGONAL_DIFFERENT]
+ ]);;
+
+(* NORMALIZING VECTORS *)
+
+let unit = new_definition `unit v = inv(norm v) % v`;;
+
+let UNIT_THM = prove
+ (`!v. ~(v = vec 0) ==> norm(unit v) = &1`,
+ SIMP_TAC[unit;NORM_MUL;REAL_ABS_INV;REAL_ABS_NORM;GSYM NORM_EQ_0;
+ REAL_MUL_LINV]);;
+
+let UNIT_INTRO = prove
+ (`!v. ~(v = vec 0) ==> v = norm v % unit v`,
+ SIMP_TAC[unit;VECTOR_MUL_ASSOC;GSYM NORM_EQ_0;REAL_MUL_RINV;
+ VECTOR_MUL_LID]);;
+
+let UNIT_ELIM = GSYM UNIT_INTRO;;
+
+let UNIT_ZERO = prove
+ (`unit(vec 0) = vec 0`, REWRITE_TAC[unit;VECTOR_MUL_RZERO]);;
+
+let UNIT_EQ_ZERO = prove
+ (`!v. unit v = vec 0 <=> v = vec 0`,
+ REWRITE_TAC[unit;VECTOR_MUL_EQ_0;REAL_INV_EQ_0;NORM_EQ_0]);;
+
+let UNIT_UNIT = prove
+ (`!v. ~(v = vec 0) ==> unit (unit v) = unit v`,
+ SIMP_TAC[unit;UNIT_THM;REAL_INV_1;VECTOR_MUL_LID]);;
+
+let UNIT_DOT_UNIT_SELF = prove
+ (`!v. ~(v = vec 0) ==> unit v dot unit v = &1`,
+ SIMP_TAC[GSYM NORM_EQ_1;UNIT_THM]);;
+
+let DOT_UNIT_SELF = prove
+ (`!v. ~(v=vec 0) ==> v dot unit v = norm v`,
+ SIMP_TAC[MESON[UNIT_INTRO] `!v. ~(v = vec 0) ==>
+ v dot unit v = (norm v % unit v) dot unit v`;DOT_LMUL;UNIT_DOT_UNIT_SELF;
+ REAL_MUL_RID]);;
+
+let ORTHOGONAL_UNIT = prove
+ (`!x y:real^N. (orthogonal x (unit y) <=> orthogonal x y)
+ /\ (orthogonal (unit x) y <=> orthogonal x y)`,
+ let common h =
+ ASM_CASES_TAC h
+ THENL [
+ RULE_ASSUM_TAC (REWRITE_RULE[NORM_EQ_0]) THEN ASM_REWRITE_TAC[UNIT_ZERO];
+ ASM_SIMP_TAC[unit;orthogonal;DOT_RMUL;DOT_LMUL;REAL_ENTIRE;REAL_INV_EQ_0]
+ ]
+ in
+ REPEAT GEN_TAC THEN CONJ_TAC
+ THENL [ common `norm (y:real^N) = &0`; common `norm (x:real^N) = &0` ]);;
+
+let ORTHOGONAL_RUNIT = prove
+ (`!x y:real^N. orthogonal x (unit y) <=> orthogonal x y`,
+ MESON_TAC[ORTHOGONAL_UNIT]);;
+
+let ORTHOGONAL_LUNIT = prove
+ (`!x y:real^N. orthogonal (unit x) y <=> orthogonal x y`,
+ MESON_TAC[ORTHOGONAL_UNIT]);;
+
+let ORTHOGONAL_LRUNIT = prove
+ (`!x y:real^N. orthogonal (unit x) (unit y) <=> orthogonal x y`,
+ MESON_TAC[ORTHOGONAL_UNIT]);;
+
+let COLLINEAR_UNIT = prove
+ (`!x y:real^N. collinear {vec 0,x,y} <=> collinear {vec 0,unit x,unit y}`,
+ REPEAT (STRIP_TAC ORELSE EQ_TAC)
+ THEN ASM_CSQ_THEN (REWRITE_RULE[GSYM IMP_IMP] COLLINEAR_VEC0)
+ (STRIP_ASSUME_TAC o REWRITE_RULE[IN_INSERT;IN_SING])
+ THEN ASM_CSQ_THEN ~remove:true (MESON[] `(!x. x=a \/ x=b \/ x=c ==> p x) ==>
+ p b /\ p c`) STRIP_ASSUME_TAC
+ THENL [
+ REWRITE_TAC[unit];
+ MAP_EVERY ASM_CASES_TAC [`x=vec 0:real^N`;`y=vec 0:real^N`]
+ THEN REPEAT (FIRST_X_ASSUM SUBST_VAR_TAC)
+ THEN REWRITE_TAC[INSERT_AC;COLLINEAR_SING;COLLINEAR_2]
+ THEN SUBGOAL_THEN `collinear {x:real^N, y, vec 0}
+ <=> collinear {norm x % unit x, norm y % unit y, vec 0}`
+ (SINGLE REWRITE_TAC)
+ THENL ON_FIRST_GOAL (ASM_SIMP_TAC[GSYM UNIT_INTRO])
+ ]
+ THEN REWRITE_TAC[collinear;IN_INSERT;IN_SING]
+ THEN EXISTS_TAC `u:real^N` THEN REPEAT STRIP_TAC
+ THEN ASM_REWRITE_TAC[VECTOR_MUL_ASSOC;VECTOR_SUB_RZERO;VECTOR_SUB_LZERO;
+ VECTOR_NEG_0;GSYM VECTOR_MUL_LNEG;GSYM VECTOR_SUB_RDISTRIB]
+ THEN MESON_TAC[VECTOR_MUL_LZERO]);;
+
+let DOT_RUNIT_POS,DOT_LUNIT_POS = CONJ_PAIR(prove
+ (`(!x y. &0 < x dot (unit y) <=> &0 < x dot y)
+ /\ (!x y. &0 < (unit x) dot y <=> &0 < x dot y)`,
+ REWRITE_TAC[unit;DOT_RMUL;DOT_LMUL;REAL_MUL_POS_LT;REAL_LT_INV_EQ;
+ NORM_POS_LT;REWRITE_RULE[REAL_INV_NEG;REAL_NEG_GT0]
+ (SPEC `--x:real` REAL_LT_INV_EQ);
+ REWRITE_RULE[REWRITE_RULE[MESON[]`(p<=> ~q) <=> (q<=> ~p)`] real_lt]
+ NORM_POS_LE]
+ THEN MESON_TAC[DOT_RZERO;DOT_LZERO;REAL_LT_REFL]));;
+
+let DOT_RUNIT_EQ0,DOT_LUNIT_EQ0 = CONJ_PAIR(prove
+ (`(!x y. x dot (unit y) = &0 <=> x dot y = &0)
+ /\ (!x y. (unit x) dot y = &0 <=> x dot y = &0)`,
+ REWRITE_TAC[GSYM orthogonal;ORTHOGONAL_LUNIT;ORTHOGONAL_RUNIT]));;
+
+let DOT_RUNIT_LT0,DOT_LUNIT_LT0 = CONJ_PAIR(prove
+ (`(!x y. x dot (unit y) < &0 <=> x dot y < &0)
+ /\ (!x y. (unit x) dot y < &0 <=> x dot y < &0)`,
+ REWRITE_TAC[real_lt;MESON[] `(~p <=> ~q) <=> (p <=> q)`]
+ THEN REWRITE_TAC[REAL_LE_LT]
+ THEN MESON_TAC[DOT_RUNIT_EQ0;DOT_LUNIT_EQ0;DOT_RUNIT_POS;DOT_LUNIT_POS]));;
+
+let DOT_RUNIT_LE0,DOT_LUNIT_LE0 = CONJ_PAIR(prove
+ (`(!x y. x dot (unit y) <= &0 <=> x dot y <= &0)
+ /\ (!x y. (unit x) dot y <= &0 <=> x dot y <= &0)`,
+ MESON_TAC[REAL_LE_LT;DOT_RUNIT_LT0;DOT_LUNIT_LT0;DOT_RUNIT_EQ0;
+ DOT_LUNIT_EQ0]));;
+
+let DOT_RUNIT_GE0,DOT_LUNIT_GE0 = CONJ_PAIR(prove
+ (`(!x y. &0 <= x dot (unit y) <=> &0 <= x dot y)
+ /\ (!x y. &0 <= (unit x) dot y <=> &0 <= x dot y)`,
+ REWRITE_TAC[REAL_LE_LT]
+ THEN MESON_TAC[DOT_RUNIT_POS;DOT_LUNIT_POS;DOT_RUNIT_EQ0;DOT_LUNIT_EQ0]));;
+
+let DOT_RUNIT_EQ = prove
+ (`!x y z:real^N. x dot (unit z) = y dot (unit z) <=> x dot z = y dot z`,
+ REPEAT GEN_TAC THEN ASM_CASES_TAC `z=vec 0:real^N`
+ THEN ASM_REWRITE_TAC[unit;DOT_RMUL;REAL_EQ_MUL_LCANCEL;REAL_INV_EQ_0;
+ NORM_EQ_0;DOT_RZERO]);;
+
+let DOT_LUNIT_EQ = prove
+ (`!x y z:real^N. (unit z) dot x = (unit z) dot y <=> z dot x = z dot y`,
+ REPEAT GEN_TAC THEN ASM_CASES_TAC `z=vec 0:real^N`
+ THEN ASM_REWRITE_TAC[unit;DOT_LMUL;REAL_EQ_MUL_LCANCEL;REAL_INV_EQ_0;
+ NORM_EQ_0;DOT_LZERO]);;
+
+
+(* CONNECTION BETWEEN SPANS AND PLANES *)
+
+new_type_abbrev("point",`:real^3`);;
+new_type_abbrev("plane",`:point->bool`);;
+
+let PLANE_NON_EMPTY = prove
+ (`!p:plane. plane p ==> ~(p={})`,
+ REWRITE_TAC[plane;GSYM MEMBER_NOT_EMPTY] THEN REPEAT STRIP_TAC
+ THEN EXISTS_TAC `u:point`
+ THEN ASM_REWRITE_TAC[] THEN MATCH_MP_TAC (REWRITE_RULE[SUBSET] HULL_SUBSET)
+ THEN SET_TAC[]);;
+
+let is_plane_basis = new_definition
+ `is_plane_basis s p <=>
+ ?u v:real^3. s = {u,v} /\ ~(collinear {vec 0,u,v}) /\ !pt:point. pt IN p
+ ==> !pt':point. pt' IN p <=> ?a b:real. pt' = pt + a % u + b % v`;;
+
+(* That theorem is (one of) the last which makes use of the affine definition of `plane` *)
+let EXISTS_PLANE_BASIS = prove
+ (`!p. plane p ==> ?b. is_plane_basis b p`,
+ GEN_TAC
+ THEN DISCH_THEN (fun x ->
+ STRIP_ASSUME_TAC (REWRITE_RULE[plane;AFFINE_HULL_3] x)
+ THEN STRIP_ASSUME_TAC (MATCH_MP (REWRITE_RULE[GSYM MEMBER_NOT_EMPTY]
+ PLANE_NON_EMPTY) x))
+ THEN POP_ASSUM MP_TAC THEN ASM_REWRITE_TAC[is_plane_basis] THEN DISCH_TAC
+ THEN MAP_EVERY EXISTS_TAC [`{v-u,w-u:real^3}`;`v-u:real^3`; `w-u:real^3`]
+ THEN ASM_REWRITE_TAC[] THEN REPEAT STRIP_TAC
+ THENL [
+ ASM_MESON_TAC[INSERT_AC;COLLINEAR_3];
+ POP_ASSUM MP_TAC THEN REWRITE_TAC[IN_ELIM_THM] THEN REPEAT STRIP_TAC
+ THEN EQ_TAC THEN STRIP_TAC
+ THENL [
+ RULE_ASSUM_TAC (REWRITE_RULE[REAL_ARITH `x+y+z= &1 <=> x= &1-y-z`])
+ THEN MAP_EVERY EXISTS_TAC [`v''-v':real`; `w''-w':real`]
+ THEN ASM_REWRITE_TAC[] THEN VECTOR_ARITH_TAC;
+ MAP_EVERY EXISTS_TAC [`u'-a-b:real`;`v'+a:real`;`w'+b:real`]
+ THEN CONJ_TAC
+ THENL [
+ ASM_ARITH_TAC;
+ REPLICATE_TAC 2 (POP_ASSUM MP_TAC) THEN VECTOR_ARITH_TAC]
+ ]]);;
+
+
+let BASIS_NON_NULL = prove
+ (`!p. plane p ==> !b. is_plane_basis b p ==> !u. u IN b ==> ~(u=vec 0)`,
+ REWRITE_TAC[is_plane_basis] THEN REPEAT STRIP_TAC
+ THEN SUBGOAL_THEN `(u':real^3) IN {u,v}` MP_TAC
+ THENL ON_FIRST_GOAL (ASM_MESON_TAC[])
+ THEN REWRITE_TAC[IN_INSERT;IN_SING]
+ THEN ASM_MESON_TAC[INSERT_AC;COLLINEAR_2]);;
+
+let PAIR_SETS_EQ = prove
+ (`!x y z t. {x,y} = {z,t} ==> (x=z /\ y=t) \/ (x=t /\ y=z)`,
+ SET_TAC[]);;
+
+let NON_COLLINEAR3_IMPLIES_DIFFERENT = prove
+ (`!x y z. ~(collinear {x,y,z}) ==> ~(x=y) /\ ~(y=z) /\ ~(x=z)`,
+ MESON_TAC[COLLINEAR_2;INSERT_AC]);;
+
+let BASIS_DIFFERENT = prove
+ (`!p. plane p ==> !u v. is_plane_basis {u,v} p ==> ~(u=v)`,
+ REWRITE_TAC[is_plane_basis] THEN REPEAT STRIP_TAC
+ THEN ASM_MESON_TAC[PAIR_SETS_EQ;NON_COLLINEAR3_IMPLIES_DIFFERENT]);;
+
+let UNIT_OF_BASIS_IS_BASIS = prove
+ (`!p. plane p ==> !u v. is_plane_basis {u,v} p ==> is_plane_basis {unit u,unit v} p`,
+ REPEAT STRIP_TAC THEN REWRITE_TAC[is_plane_basis]
+ THEN MAP_EVERY EXISTS_TAC [`unit u:real^3`; `unit v:real^3`]
+ THEN REWRITE_TAC[]
+ THEN POP_ASSUM
+ (DISTRIB [ASSUME_TAC;STRIP_ASSUME_TAC o REWRITE_RULE[is_plane_basis]])
+ THEN ASM_REWRITE_TAC[GSYM COLLINEAR_UNIT]
+ THEN ASM_CSQ_THEN ~remove:true PAIR_SETS_EQ STRIP_ASSUME_TAC
+ THEN REPEAT (FIRST_X_ASSUM SUBST_VAR_TAC)
+ THEN ASM_CSQ_THEN BASIS_NON_NULL (C ASM_CSQ_THEN (STRIP_ASSUME_TAC
+ o REWRITE_RULE[IN_INSERT;IN_SING;
+ MESON[]`(!u. u=a \/ u=b ==> p u) <=> p a /\ p b`]))
+ THEN GEN_TAC
+ THEN DISCH_THEN (fun x -> FIRST_X_ASSUM (ASSUME_TAC o C MATCH_MP x))
+ THEN ASM_REWRITE_TAC[unit;VECTOR_MUL_ASSOC]
+ THEN REPEAT (STRIP_TAC ORELSE EQ_TAC) THEN ASM_REWRITE_TAC[]
+ THENL [
+ MAP_EVERY EXISTS_TAC [`a*norm(u':real^3)`;`b*norm(v':real^3)`]
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[GSYM NORM_EQ_0])
+ THEN ASM_SIMP_TAC[GSYM REAL_MUL_ASSOC;REAL_MUL_RINV] THEN VECTOR_ARITH_TAC;
+ MESON_TAC[];
+ MAP_EVERY EXISTS_TAC [`b*norm(v':real^3)`;`a*norm(u':real^3)`]
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[GSYM NORM_EQ_0])
+ THEN ASM_SIMP_TAC[GSYM REAL_MUL_ASSOC;REAL_MUL_RINV] THEN VECTOR_ARITH_TAC;
+ MESON_TAC[VECTOR_ADD_AC];
+ ]);;
+
+(* This theorem makes the connection between the (affine) notion of plane and the (vector) notion of span.
+ * From now on, we can generally rely only on this result to reason about the spanning set of a plane.
+ *)
+let PLANE_AS_SPAN = prove
+ (`!p. plane p ==> !b. is_plane_basis b p ==> !pt:point. pt IN p
+ ==> p = { pt + y | y IN span b }`,
+ REWRITE_TAC[is_plane_basis] THEN REPEAT STRIP_TAC
+ THEN ASM_REWRITE_TAC[EXTENSION;IN_ELIM_THM;SPAN_2;UNIV]
+ THEN ASM_MESON_TAC[]);;
+
+let PLANE_SPAN = prove
+ (`!p. plane p ==> !b. is_plane_basis b p ==> !pt:point. pt IN p
+ ==> { pt' - pt | pt' IN p } = span b`,
+ REPEAT STRIP_TAC THEN ASM_CSQ_THEN PLANE_AS_SPAN (C ASM_CSQ_THEN (C
+ ASM_CSQ_THEN (SINGLE (GEN_REWRITE_TAC (RATOR_CONV o ONCE_DEPTH_CONV)))))
+ THEN REWRITE_TAC[IN_ELIM_THM;EXTENSION] THEN GEN_TAC THEN EQ_TAC
+ THENL [
+ REPEAT STRIP_TAC THEN ASM_REWRITE_TAC[VECTOR_ARITH `(x + y)-x = y:real^N`];
+ DISCH_TAC THEN EXISTS_TAC `x+pt:real^3` THEN REWRITE_TAC (map VECTOR_ARITH
+ [`(x + y) - y = x:real^N`; `x+z=z+y <=> x=y:real^N`]) THEN ASM_MESON_TAC[]
+ ]);;
+
+let ALL_BASIS_SAME_SPAN = prove
+ (`!p. plane p ==> !b1. is_plane_basis b1 p ==> !b2. is_plane_basis b2 p
+ ==> span b1 = span b2`,
+ REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN PLANE_SPAN (fun x ->
+ ASM_CSQ_THEN ~match_:false (SPEC `b1:real^3->bool` x) ASSUME_TAC
+ THEN ASM_CSQ_THEN (SPEC `b2:real^3->bool` x) ASSUME_TAC)
+ THEN ASM_MESON_TAC[PLANE_NON_EMPTY;MEMBER_NOT_EMPTY]);;
+
+let PLANE_SUBSPACE = prove
+ (`!p. plane p ==> !pt:point. pt IN p ==> subspace { pt' - pt | pt' IN p }`,
+ REPEAT STRIP_TAC THEN ASM_CSQ_THEN EXISTS_PLANE_BASIS STRIP_ASSUME_TAC
+ THEN ASM_SIMP_TAC[PLANE_SPAN;SUBSPACE_SPAN]);;
+
+let PLANE_SING = prove
+ (`!x:real^N. ~(plane{x})`,
+ REWRITE_TAC[plane;NOT_EXISTS_THM;MESON[] `~(A/\B) <=> (A ==> ~B)`]
+ THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN NON_COLLINEAR3_IMPLIES_DIFFERENT ASSUME_TAC
+ THEN MATCH_MP_TAC (SET_RULE `~(u=v:real^N) /\ u IN {x} /\ v IN {x} ==> F`)
+ THEN ASM_REWRITE_TAC[]
+ THEN CONJ_TAC THEN MATCH_MP_TAC HULL_INC THEN SET_TAC[]);;
+
+let NON_COLLINEAR_INDEPENDENT = prove
+ (`!x y:real^N. ~(collinear {vec 0,x,y}) ==> independent{x,y}`,
+ MAP_EVERY (SINGLE ONCE_REWRITE_TAC) [independent;CONTRAPOS_THM;NOT_CLAUSES]
+ THEN REPEAT GEN_TAC
+ THEN ASM_CASES_TAC `x=y:real^N` THENL [
+ ASM_REWRITE_TAC[INSERT_AC;COLLINEAR_2];
+ REWRITE_TAC[dependent;IN_INSERT;IN_SING] THEN STRIP_TAC THENL
+ let common main_var exist_list =
+ POP_ASSUM MP_TAC
+ THEN ASM_REWRITE_TAC[DELETE_INSERT;EMPTY_DELETE;SPAN_SING;IN_ELIM_THM;
+ UNIV]
+ THEN STRIP_TAC
+ THEN REWRITE_TAC[collinear;IN_INSERT;IN_SING] THEN EXISTS_TAC main_var
+ THEN REPEAT STRIP_TAC
+ THENL map (fun t ->
+ EXISTS_TAC t THEN ASM_REWRITE_TAC[] THEN VECTOR_ARITH_TAC) exist_list
+ in
+ [ common `y:real^N`
+ [`&0`;`--u:real`;`-- &1`;`u:real`;`&0`;`u- &1`;`&1`;`&1-u`;`&0`];
+ common `x:real^N`
+ [`&0`;`-- &1`;`--u:real`;`&1`;`&0`;`&1-u`;`u:real`;`u- &1`;`&0`]
+ ]]);;
+
+let DIM_OF_PLANE_SPAN = prove
+ (`!p:plane. plane p ==> !b. is_plane_basis b p ==> dim (span b) = 2`,
+ REWRITE_TAC[is_plane_basis] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN NON_COLLINEAR_INDEPENDENT ASSUME_TAC
+ THEN ASM_CSQ_THEN NON_COLLINEAR3_IMPLIES_DIFFERENT ASSUME_TAC
+ THEN SUBGOAL_THEN `{u,v:real^3} HAS_SIZE (dim (span {u,v}))` MP_TAC
+ THENL [
+ MATCH_MP_TAC BASIS_HAS_SIZE_DIM THEN ASM_REWRITE_TAC[];
+ ASM_SIMP_TAC[HAS_SIZE;FINITE_INSERT;FINITE_EMPTY;CARD_CLAUSES;IN_INSERT;
+ IN_SING;NOT_IN_EMPTY] THEN ARITH_TAC
+ ]);;
+
+let DIM_OF_PLANE_SUBSPACE = prove
+ (`!p. plane p ==> !pt:point. pt IN p ==> dim { pt' - pt | pt' IN p } = 2`,
+ MESON_TAC[PLANE_SPAN;DIM_OF_PLANE_SPAN;EXISTS_PLANE_BASIS]);;
+
+let PLANE_SPAN_DECOMPOS = prove
+ (`!p. plane p ==> !u v. is_plane_basis {u,v} p ==>
+ !x. x IN span {u,v} <=> ?a b. x = a % u + b % v`,
+ REPEAT STRIP_TAC THEN POP_ASSUM (SINGLE REWRITE_TAC o GSYM)
+ THEN REWRITE_TAC[SPAN_2;IN_ELIM_THM;UNIV]);;
+
+let units = new_definition `units = { x | norm x = &1 }`;;
+
+let IN_UNITS = prove
+ (`!x. x IN units ==> norm x = &1`,
+ REWRITE_TAC[units;IN_ELIM_THM]);;
+
+let is_normal_to_plane = new_definition
+ `is_normal_to_plane (n:real^3) (p:plane) <=> ~(n=vec 0)
+ /\ !b. is_plane_basis b p ==> !v. v IN (span b) ==> orthogonal v n`;;
+
+let IS_NORMAL_TO_PLANE_UNIT = prove
+ (`!p n. is_normal_to_plane (unit n) p <=> is_normal_to_plane n p`,
+ REWRITE_TAC[is_normal_to_plane;ORTHOGONAL_RUNIT;UNIT_EQ_ZERO]);;
+
+let EXISTS_NORMAL_OF_PLANE = prove
+ (`!p:plane. plane p ==> ?n:real^3. is_normal_to_plane n p`,
+ REWRITE_TAC[is_normal_to_plane] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN EXISTS_PLANE_BASIS (CHOOSE_THEN (fun x ->
+ MAP_EVERY STRIP_ASSUME_TAC [x;REWRITE_RULE[is_plane_basis] x]))
+ THEN EXISTS_TAC `u cross v` THEN ASM_REWRITE_TAC[CROSS_EQ_0]
+ THEN REPEAT STRIP_TAC THEN POP_ASSUM MP_TAC
+ THEN ASM_CSQ_THEN ALL_BASIS_SAME_SPAN (C (ASM_CSQ_THEN ~remove:true) (C
+ ASM_CSQ_THEN (SINGLE REWRITE_TAC)))
+ THEN ASM_REWRITE_TAC[SPAN_2;IN_ELIM_THM;UNIV]
+ THEN REPEAT STRIP_TAC
+ THEN ASM_REWRITE_TAC[orthogonal;DOT_LADD;DOT_LMUL;DOT_CROSS_SELF]
+ THEN SIMPLE_COMPLEX_ARITH_TAC);;
+
+let NORMAL_OF_PLANE_IS_NORMAL_TO_BASIS = prove
+ (`!p. plane p ==> !b. is_plane_basis b p ==> !u. u IN b ==>
+ !n. is_normal_to_plane n p ==> orthogonal u n`,
+ MESON_TAC[is_normal_to_plane;SPAN_SUPERSET;]);;
+
+let NORMAL_OF_PLANE_NON_NULL = prove
+ (`!p. plane p ==> !n. is_normal_to_plane n p ==> ~(n=vec 0)`,
+ MESON_TAC[is_normal_to_plane]);;
+
+let NORMAL_OF_PLANE_IS_ORTHOGONAL_TO_SEGMENT = prove
+ (`!p. plane p ==> !n. is_normal_to_plane n p
+ ==> !pt1 pt2. pt1 IN p /\ pt2 IN p ==> orthogonal (pt1-pt2) n`,
+ REWRITE_TAC[is_normal_to_plane] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN EXISTS_PLANE_BASIS STRIP_ASSUME_TAC
+ THEN POP_ASSUM
+ (fun x -> FIRST_X_ASSUM (MATCH_MP_TAC o C MATCH_MP x) THEN ASSUME_TAC x)
+ THEN FIRST_X_ASSUM (STRIP_ASSUME_TAC o CONV_RULE (REWR_CONV is_plane_basis))
+ THEN POP_ASSUM (C (ASM_CSQ_THEN ~remove:true) (fun x ->
+ FIRST_X_ASSUM (STRIP_ASSUME_TAC o CONV_RULE (REWR_CONV x))))
+ THEN ASM_REWRITE_TAC[SPAN_2;UNIV;IN_ELIM_THM;VECTOR_ARITH
+ `(x+y+z)-x = y+z:real^N`]
+ THEN MESON_TAC[]);;
+
+let DIM_OF_NORMAL_SPAN = prove
+ (`!p:plane. plane p ==> !n. is_normal_to_plane n p ==> dim (span {n}) = 1`,
+ REWRITE_TAC[is_normal_to_plane] THEN REPEAT STRIP_TAC
+ THEN SUBGOAL_THEN `{n:real^3} HAS_SIZE dim (span {n})` MP_TAC
+ THENL [
+ MATCH_MP_TAC BASIS_HAS_SIZE_DIM
+ THEN REWRITE_TAC[independent;dependent;NOT_EXISTS_THM;
+ MESON[] `~(p/\q) <=> p ==> ~q`;IN_SING]
+ THEN GEN_TAC THEN DISCH_TAC
+ THEN ASM_REWRITE_TAC[DELETE_INSERT;EMPTY_DELETE;SPAN_EMPTY;IN_SING];
+ ASM_SIMP_TAC[HAS_SIZE;FINITE_INSERT;FINITE_EMPTY;CARD_CLAUSES;IN_INSERT;
+ IN_SING;NOT_IN_EMPTY]
+ THEN ARITH_TAC]);;
+
+let DIM_OF_ORTHOGONAL = prove
+ (`!p. plane p ==> !b. is_plane_basis b p
+ ==> dim {z:real^3 | !x. x IN span b ==> orthogonal z x} = 1`,
+ REPEAT STRIP_TAC
+ THEN SUBGOAL_THEN
+ `(:real^3) = { x+y
+ | x IN span b /\ y IN {z:real^3 | !x. x IN span b ==> orthogonal z x}}`
+ (ASSUME_TAC o REWRITE_RULE[DIM_UNIV;DIMINDEX_3] o
+ AP_TERM `dim:(real^3->bool)->num`)
+ THENL ON_FIRST_GOAL
+ (REWRITE_TAC[GSYM SUBSET_ANTISYM_EQ;SUBSET;UNIV;IN;IN_ELIM_THM]
+ THEN MESON_TAC[ORTHOGONAL_SUBSPACE_DECOMP_EXISTS;IN])
+ THEN SUBGOAL_THEN
+ `span b
+ INTER {z:real^3 | !x. x IN span b ==> orthogonal z x} SUBSET {vec 0}`
+ (ASSUME_TAC o GSYM o REWRITE_RULE[GSYM DIM_EQ_0])
+ THENL ON_FIRST_GOAL
+ (REWRITE_TAC[INTER;SUBSET;IN_SING;IN_ELIM_THM]
+ THEN MESON_TAC[ORTHOGONAL_REFL])
+ THEN ASM_CSQ_THEN DIM_OF_PLANE_SPAN (fun x ->
+ ASM_SIMP_TAC[ARITH_RULE `x=1 <=> 3+0=2+x`;GSYM x])
+ THEN MATCH_MP_TAC DIM_SUMS_INTER
+ THEN REWRITE_TAC[SUBSPACE_ORTHOGONAL_TO_VECTORS;ORTHOGONAL_SYM;
+ SUBSPACE_SPAN]
+ );;
+
+let NORMAL_SPAN_SUBSET_ORTHOGONAL_SUBSPACE = prove
+ (`!p. plane p ==> !b. is_plane_basis b p ==> !n. is_normal_to_plane n p
+ ==> span {n} SUBSET {z:real^3 | !x. x IN span b ==> orthogonal z x}`,
+ REWRITE_TAC[is_normal_to_plane] THEN REPEAT STRIP_TAC
+ THEN POP_ASSUM (C ASM_CSQ_THEN ASSUME_TAC)
+ THEN REWRITE_TAC[SUBSET;SPAN_SING;IN_ELIM_THM;UNIV]
+ THEN ASM_MESON_TAC[ORTHOGONAL_CLAUSES;ORTHOGONAL_SYM]);;
+
+let ORTHOGONAL_SUBSPACE_IS_NORMAL_SPAN = prove
+ (`!p. plane p ==> !b. is_plane_basis b p ==> !n. is_normal_to_plane n p
+ ==> {z:real^3 | !x. x IN span b ==> orthogonal z x} = span {n}`,
+ MAP_EVERY (SINGLE ONCE_REWRITE_TAC) [
+ ORTHOGONAL_SYM; MESON[SPAN_SPAN] `span{x}=span(span{x})`;
+ GSYM (REWRITE_RULE[GSYM SPAN_EQ_SELF] SUBSPACE_ORTHOGONAL_TO_VECTORS)]
+ THEN REPEAT STRIP_TAC THEN MATCH_MP_TAC EQ_SYM
+ THEN MATCH_MP_TAC DIM_EQ_SPAN THEN ONCE_REWRITE_TAC[ORTHOGONAL_SYM]
+ THEN CONJ_TAC
+ THENL ON_FIRST_GOAL (ASM_MESON_TAC[NORMAL_SPAN_SUBSET_ORTHOGONAL_SUBSPACE])
+ THEN MAP_EVERY (fun x ->
+ ASM_CSQ_THEN x (C ASM_CSQ_THEN (SINGLE REWRITE_TAC)))
+ [DIM_OF_ORTHOGONAL;DIM_OF_NORMAL_SPAN]
+ THEN ARITH_TAC);;
+
+let PLANE_DECOMP = prove
+ (`!p. plane p ==> !b. is_plane_basis b p ==> !n. is_normal_to_plane n p
+ ==> !v. ?w x. v = w + x % n /\ w IN span b`,
+ REPEAT STRIP_TAC
+ THEN MP_TAC (ISPEC `b:real^3->bool` ORTHOGONAL_SUBSPACE_DECOMP)
+ THEN ASM_CSQ_THEN ORTHOGONAL_SUBSPACE_IS_NORMAL_SPAN
+ (C ASM_CSQ_THEN (C ASM_CSQ_THEN (SINGLE REWRITE_TAC)))
+ THEN REWRITE_TAC[SPAN_SING;IN_ELIM_THM;EXISTS_UNIQUE_DEF;EXISTS_PAIRED_THM]
+ THEN MESON_TAC[]);;
+
+let NORMAL_ORTHOGONAL_IN_PLANE_SPAN = prove
+ (`!p. plane p ==> !b. is_plane_basis b p ==> !n. is_normal_to_plane n p
+ ==> !v. orthogonal v n <=> v IN span b`,
+ let EQ_IMPLY =
+ MATCH_MP (MESON[] `(!x y. p x y <=> q x y) ==> (!x y. p x y ==> q x y)`)
+ in
+ REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN PLANE_DECOMP (C ASM_CSQ_THEN (C ASM_CSQ_THEN
+ (STRIP_ASSUME_TAC o SPEC `v:real^3`)))
+ THEN ASM_CSQ_THEN (EQ_IMPLY is_normal_to_plane) STRIP_ASSUME_TAC
+ THEN POP_ASSUM (C ASM_CSQ_THEN ASSUME_TAC) THEN EQ_TAC THEN ASM_REWRITE_TAC[]
+ THEN ASM_CSQ_THEN (EQ_IMPLY is_plane_basis) STRIP_ASSUME_TAC
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[orthogonal])
+ THEN ASM_SIMP_TAC[orthogonal;DOT_LADD;DOT_LMUL;REAL_ADD_LID;REAL_ENTIRE;
+ DOT_EQ_0;VECTOR_MUL_LZERO;VECTOR_ADD_RID] THEN ASM_MESON_TAC[]);;
+
+let PLANE_DECOMP_DOT_NORMAL = prove
+ (`!p. plane p ==> !b. is_plane_basis b p ==> !n. is_normal_to_plane n p
+ ==> !v. ?w. w IN span b /\ v = w + (v dot (unit n)) % (unit n)`,
+ REPEAT STRIP_TAC
+ THEN FIRST_ASSUM (ASSUME_TAC o CONJUNCT1 o CONV_RULE (REWR_CONV
+ is_normal_to_plane))
+ THEN ASM_CSQ_THEN (GSYM NORMAL_ORTHOGONAL_IN_PLANE_SPAN) (C ASM_CSQ_THEN
+ (C ASM_CSQ_THEN ASSUME_TAC))
+ THEN EXISTS_TAC `v - (v dot unit n) % unit n :real^3`
+ THEN ASM_REWRITE_TAC[VECTOR_ARITH `x=x-y+z <=> y=z:real^N`]
+ THEN ONCE_REWRITE_TAC[GSYM ORTHOGONAL_RUNIT]
+ THEN REWRITE_TAC[orthogonal;DOT_LSUB;DOT_LMUL]
+ THEN ASM_MESON_TAC[UNIT_DOT_UNIT_SELF;REAL_MUL_RID;REAL_SUB_0]);;
+
+let SUB_UNIT_NORMAL_IS_ORTHOGONAL_TO_NORMAL = prove
+ (`!p. plane p ==> !n. is_normal_to_plane n p
+ ==> !x. orthogonal (x - (x dot unit n) % unit n) n`,
+ REPEAT STRIP_TAC THEN ASM_CSQ_THEN EXISTS_PLANE_BASIS STRIP_ASSUME_TAC
+ THEN ASM_CSQ_THEN PLANE_DECOMP_DOT_NORMAL
+ (C ASM_CSQ_THEN (C ASM_CSQ_THEN ASSUME_TAC))
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[VECTOR_ARITH `x=y+z <=> x-z=y`])
+ THEN ASM_MESON_TAC[NORMAL_ORTHOGONAL_IN_PLANE_SPAN]);;
+
+let is_orthogonal_plane_basis = new_definition
+ `is_orthogonal_plane_basis b p
+ <=> is_plane_basis b p /\ pairwise orthogonal b`;;
+
+let EXISTS_ORTHOGONAL_PLANE_BASIS = prove
+ (`!p. plane p ==> ?b. is_orthogonal_plane_basis b p`,
+ REWRITE_TAC[is_orthogonal_plane_basis;is_plane_basis] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN EXISTS_PLANE_BASIS (CHOOSE_THEN (fun x -> ASSUME_TAC x
+ THEN STRIP_ASSUME_TAC (REWRITE_RULE[is_plane_basis] x)))
+ THEN ASM_CSQ_THEN EXISTS_NORMAL_OF_PLANE (CHOOSE_THEN (fun x ->
+ let x' = REWRITE_RULE[is_normal_to_plane] x in
+ ASSUME_TAC x THEN ASM_CSQ_THEN (CONJUNCT2 x') ASSUME_TAC
+ THEN ASSUME_TAC (CONJUNCT1 x')))
+ THEN SUBGOAL_TAC "" `~((u:real^3)=vec 0)` [ASM SET_TAC[BASIS_NON_NULL]]
+ THEN SUBGOAL_TAC "" `orthogonal u (n:real^3)` [
+ FIRST_ASSUM MATCH_MP_TAC THEN MATCH_MP_TAC SPAN_SUPERSET
+ THEN ASM SET_TAC [SPAN_SUPERSET]
+ ]
+ THEN SUBGOAL_TAC "" `~(u cross n = vec 0)` [
+ ASM_SIMP_TAC[CROSS_EQ_0;ORTHOGONAL_NON_COLLINEAR] ]
+ THEN EXISTS_TAC `{u, u cross n}` THEN CONJ_TAC
+ THENL [
+ MAP_EVERY EXISTS_TAC [`u:real^3`;`u cross n`] THEN REWRITE_TAC[]
+ THEN CONJ_TAC THENL ON_FIRST_GOAL
+ (MATCH_MP_TAC ORTHOGONAL_NON_COLLINEAR
+ THEN ASM_REWRITE_TAC[orthogonal;DOT_CROSS_SELF])
+ THEN SUBGOAL_THEN `?a b. u cross n = a % u + b % v` STRIP_ASSUME_TAC
+ THENL ON_FIRST_GOAL
+ (ASM_MESON_TAC[GSYM PLANE_SPAN_DECOMPOS;
+ GSYM NORMAL_ORTHOGONAL_IN_PLANE_SPAN;orthogonal;DOT_CROSS_SELF])
+ THEN SUBGOAL_THEN `~(b'= &0)` (fun x -> POP_ASSUM MP_TAC THEN ASSUME_TAC x)
+ THENL ON_FIRST_GOAL
+ (ASM_CASES_TAC `a= &0` THENL [
+ ASM_MESON_TAC[VECTOR_MUL_LZERO;VECTOR_ADD_LID];
+ DISCH_THEN (fun x -> POP_ASSUM (fun y -> POP_ASSUM (MP_TAC
+ o REWRITE_RULE[x;VECTOR_MUL_LZERO;VECTOR_ADD_RID])
+ THEN ASSUME_TAC y))
+ THEN DISCH_THEN (MP_TAC o AP_TERM `(dot) (u:real^3)`)
+ THEN REWRITE_TAC[DOT_CROSS_SELF;DOT_RMUL]
+ THEN ONCE_REWRITE_TAC[EQ_SYM_EQ]
+ THEN ASM_REWRITE_TAC[REAL_ENTIRE;DOT_EQ_0]
+ ])
+ THEN DISCH_THEN (fun x ->
+ ASM_SIMP_TAC[x;VECTOR_ARITH `a%x+b%(c%x+d%y)=(a+b*c)%x+(b*d)%y`])
+ THEN REPEAT STRIP_TAC THEN EQ_TAC THENL [
+ REPEAT STRIP_TAC
+ THEN MAP_EVERY EXISTS_TAC [`a'-b''/b'*a:real`;`b''/b':real`]
+ THEN ASM_SIMP_TAC[REAL_SUB_ADD;VECTOR_ARITH `x+y=x+z <=> y=z`;
+ REAL_DIV_RMUL];
+ MESON_TAC[];
+ ];
+ REWRITE_TAC[pairwise;orthogonal] THEN SET_TAC [DOT_CROSS_SELF]
+ ]);;
+
+let ORTHOGONAL_PLANE_BASIS_AS_PAIR = prove
+ (`!b p. is_orthogonal_plane_basis b p ==> ?u v. is_orthogonal_plane_basis {u,v} p`,
+ MESON_TAC[is_orthogonal_plane_basis;is_plane_basis]);;
+
+let BASIS_OF_PLANE_ORTHONORMAL_DECOMPOS = prove
+ (`!p. plane p ==> !u v. is_orthogonal_plane_basis {u,v} p ==>
+ !x. x IN span {u,v}
+ <=> x = (x dot unit u) % unit u + (x dot unit v) % unit v`,
+ REWRITE_TAC [is_orthogonal_plane_basis] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN PLANE_SPAN_DECOMPOS (C ASM_CSQ_THEN (SINGLE REWRITE_TAC))
+ THEN EQ_TAC THENL [
+ REPEAT STRIP_TAC
+ THEN SUBGOAL_TAC "" `~(u = vec 0:real^3) /\ ~(v = vec 0:real^3)`
+ [ASM_MESON_TAC[BASIS_NON_NULL;IN_INSERT]]
+ THEN SUBGOAL_THEN
+ `orthogonal (u:real^3) (unit v) /\ orthogonal v (unit u)`
+ (ASSUME_TAC o REWRITE_RULE[orthogonal])
+ THENL ON_FIRST_GOAL
+ (RULE_ASSUM_TAC (REWRITE_RULE[pairwise;orthogonal])
+ THEN REWRITE_TAC[ORTHOGONAL_UNIT;orthogonal]
+ THEN SUBGOAL_TAC "" `~(u=v:real^3)` [ASM_MESON_TAC[BASIS_DIFFERENT]]
+ THEN CONJ_TAC THEN FIRST_ASSUM MATCH_MP_TAC THEN ASSUM_LIST SET_TAC)
+ THEN ASM_SIMP_TAC[DOT_LADD;DOT_LMUL;REAL_MUL_RZERO;REAL_ADD_RID;
+ REAL_ADD_LID;DOT_UNIT_SELF;GSYM VECTOR_MUL_ASSOC;GSYM UNIT_INTRO];
+ REWRITE_TAC[unit;VECTOR_MUL_ASSOC] THEN MESON_TAC[]
+ ]);;
+
+let PLANE_DECOMP_DOT = prove
+ (`!p:plane. plane p ==> !u v. is_orthogonal_plane_basis {u,v} p
+ ==> !n. is_normal_to_plane n p ==> !x:real^3.
+ x = (x dot unit u) % unit u + (x dot unit v) % unit v
+ + (x dot unit n) % unit n`,
+ REPEAT STRIP_TAC
+ THEN FIRST_ASSUM (STRIP_ASSUME_TAC o CONV_RULE (REWR_CONV
+ is_orthogonal_plane_basis))
+ THEN ASM_CSQ_THEN PLANE_DECOMP_DOT_NORMAL
+ (C ASM_CSQ_THEN (C ASM_CSQ_THEN (STRIP_ASSUME_TAC o SPEC `x:real^3`)))
+ THEN SUBGOAL_THEN
+ `(x:real^3) dot unit u = w dot unit u /\ x dot unit v = w dot unit v`
+ ASSUME_TAC
+ THENL ON_FIRST_GOAL
+ begin
+ SUBGOAL_THEN `(u:real^3) IN {u,v} /\ v IN {u,v}` STRIP_ASSUME_TAC
+ THENL ON_FIRST_GOAL (ASSUM_LIST SET_TAC)
+ THEN ASM_CSQ_THEN NORMAL_OF_PLANE_IS_NORMAL_TO_BASIS
+ (C ASM_CSQ_THEN ASSUME_TAC)
+ THEN POP_ASSUM (fun x -> ASM_CSQ_THEN ~remove:true x
+ (C ASM_CSQ_THEN ASSUME_TAC)
+ THEN ASM_CSQ_THEN x (C ASM_CSQ_THEN ASSUME_TAC))
+ THEN RULE_ASSUM_TAC (REWRITE_RULE[orthogonal]
+ o ONCE_REWRITE_RULE[GSYM ORTHOGONAL_LRUNIT]
+ o ONCE_REWRITE_RULE[ORTHOGONAL_SYM])
+ THEN ASM ONCE_REWRITE_TAC[]
+ THEN ASM_REWRITE_TAC[DOT_LADD;DOT_LMUL;REAL_MUL_RZERO;REAL_ADD_RID]
+ end
+ THEN ASM_CSQ_THEN (REWRITE_RULE[GSYM IMP_IMP]
+ BASIS_OF_PLANE_ORTHONORMAL_DECOMPOS) (C ASM_CSQ_THEN (ASSUME_TAC
+ o SPEC `w:real^3`))
+ THEN ASSUM_LIST (GEN_REWRITE_TAC (RATOR_CONV o DEPTH_CONV))
+ THEN ASM_REWRITE_TAC[VECTOR_ARITH `x+y=z+t+y <=> x=z+t`]
+ THEN ASM_MESON_TAC[]);;
+
+let UNIT_OF_ORTHOGONAL_BASIS_IS_ORTHOGONAL_BASIS = prove
+ (`!p. plane p ==> !u v. is_orthogonal_plane_basis {u,v} p
+ ==> is_orthogonal_plane_basis {unit u,unit v} p`,
+ SIMP_TAC[is_orthogonal_plane_basis;UNIT_OF_BASIS_IS_BASIS;pairwise;IN_INSERT;
+ IN_SING]
+ THEN REPEAT STRIP_TAC THEN ASM_REWRITE_TAC[]
+ THEN ASM_MESON_TAC[ORTHOGONAL_LRUNIT]);;
+
+let EXISTS_UNIT_ORTHOGONAL_PLANE_BASIS = prove
+ (`!p. plane p ==> ?u v. is_orthogonal_plane_basis {unit u,unit v} p`,
+ REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN EXISTS_ORTHOGONAL_PLANE_BASIS (CHOOSE_THEN
+ (DISTRIB (map (fun xs -> STRIP_ASSUME_TAC o REWRITE_RULE xs) [
+ [];
+ [is_orthogonal_plane_basis];
+ [is_orthogonal_plane_basis;is_plane_basis]])))
+ THEN FIRST_X_ASSUM SUBST_VAR_TAC
+ THEN ASM_MESON_TAC[UNIT_OF_ORTHOGONAL_BASIS_IS_ORTHOGONAL_BASIS]);;
+
+let FORALL_PLANE_THM = prove
+ (`!p. plane p ==>
+ ?u v. is_orthogonal_plane_basis {unit u,unit v} p /\ !P pt0. pt0 IN p
+ ==> ((!pt. pt IN p ==> P pt)
+ <=> (!a b. P (pt0 + a % unit u + b % unit v)))`,
+ REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN EXISTS_UNIT_ORTHOGONAL_PLANE_BASIS (CHOOSE_THEN
+ (CHOOSE_THEN (DISTRIB (map (fun xs -> STRIP_ASSUME_TAC o REWRITE_RULE xs) [
+ [];
+ [is_orthogonal_plane_basis];
+ [is_orthogonal_plane_basis;is_plane_basis]]))))
+ THEN MAP_EVERY EXISTS_TAC [`u:real^3`;`v:real^3`]
+ THEN ASM_REWRITE_TAC[] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN PAIR_SETS_EQ STRIP_ASSUME_TAC
+ THEN REPEAT (FIRST_X_ASSUM SUBST_VAR_TAC)
+ THEN ASM_CSQ_THEN PLANE_AS_SPAN (C ASM_CSQ_THEN (C ASM_CSQ_THEN (ASSUME_TAC
+ o REWRITE_RULE[SPAN_2;UNIV;IN_ELIM_THM])))
+ THEN ASSUM_LIST SET_TAC);;
+
+let FORALL_PLANE_THM_2 = prove
+ (`!p. plane p ==>
+ ?u v pt0. is_orthogonal_plane_basis {unit u,unit v} p /\ pt0 IN p
+ /\ !P. (!pt. pt IN p ==> P pt)
+ <=> (!a b. P (pt0 + a % unit u + b % unit v))`,
+ REPEAT STRIP_TAC THEN ASM_CSQ_THEN FORALL_PLANE_THM STRIP_ASSUME_TAC
+ THEN ASM_CSQ_THEN PLANE_NON_EMPTY (STRIP_ASSUME_TAC o REWRITE_RULE[GSYM MEMBER_NOT_EMPTY]) THEN ASM_MESON_TAC[]);;
+
+(** Is u the projection of v on the hyperplane represented by the non-null vector w? *)
+let is_projection_on_hyperplane = new_definition
+ `is_projection_on_hyperplane u v (w:real^N)
+ <=> norm w = &1 ==> u = v - (v dot w) % w`;;
+
+let projection_on_hyperplane = new_definition
+ `projection_on_hyperplane v w = @u. is_projection_on_hyperplane u v w`;;
+
+let PROJECTION_ON_HYPERPLANE_THM = prove
+ (`!v w. norm w = &1
+ ==> projection_on_hyperplane v w = v - (v dot w) % w`,
+ REWRITE_TAC[projection_on_hyperplane;is_projection_on_hyperplane]
+ THEN SELECT_ELIM_TAC THEN MESON_TAC[]);;
+
+let PROJECTION_ON_HYPERPLANE_THM_LNEG = prove
+ (`!v w. norm w = &1
+ ==> projection_on_hyperplane (--v) w = --(projection_on_hyperplane v w)`,
+ SIMP_TAC[PROJECTION_ON_HYPERPLANE_THM;VECTOR_NEG_SUB;DOT_LNEG;
+ VECTOR_MUL_LNEG]
+ THEN REPEAT STRIP_TAC THEN VECTOR_ARITH_TAC);;
+
+let PROJECTION_ON_HYPERPLANE_DECOMPOS = prove
+ (`!v w. norm w = &1 ==> projection_on_hyperplane v w + (v dot w) % w = v`,
+ SIMP_TAC[projection_on_hyperplane;is_projection_on_hyperplane]
+ THEN REPEAT STRIP_TAC THEN VECTOR_ARITH_TAC);;
+
+let symetric_vectors_wrt = new_definition
+ `symetric_vectors_wrt u v (w:real^3)
+ <=> norm w = &1 ==> u + v = (&2 * u dot w) % w`;;
+
+let SYMETRIC_VECTORS_WRT_SYM = prove
+ (`!u v w. symetric_vectors_wrt u v w <=> symetric_vectors_wrt v u w`,
+ REWRITE_TAC[symetric_vectors_wrt] THEN REPEAT GEN_TAC THEN EQ_TAC
+ THEN DISCH_THEN (fun x -> DISCH_THEN (DISTRIB (map ((o) ASSUME_TAC)
+ [MP x; REWRITE_RULE[NORM_EQ_1]])))
+ THEN ASSUM_LIST (GEN_REWRITE_TAC (RAND_CONV o DEPTH_CONV)
+ o map (REWRITE_RULE[VECTOR_ARITH `x+y=z <=> y=z-x:real^N`]))
+ THEN ASM_REWRITE_TAC[DOT_LSUB;DOT_LMUL;REAL_ARITH `(&2 * x) * &1 - x = x`]
+ THEN ONCE_REWRITE_TAC[VECTOR_ADD_SYM] THEN ASM_REWRITE_TAC[]);;
+
+let SYMETRIC_VECTORS_PROJECTION_ON_HYPERPLANE = prove
+ (`!u v w. norm w = &1 ==> symetric_vectors_wrt u v w
+ ==> --(projection_on_hyperplane v w) = projection_on_hyperplane u w`,
+ ONCE_REWRITE_TAC[MATCH_MP (MESON[] `(p <=> q) ==> (p <=> p /\ q)`)
+ (SPEC_ALL SYMETRIC_VECTORS_WRT_SYM)]
+ THEN SIMP_TAC[symetric_vectors_wrt;projection_on_hyperplane;
+ is_projection_on_hyperplane;
+ VECTOR_ARITH `x = (&2 * y) % z <=> y%z = inv(&2) % x`]
+ THEN REPEAT STRIP_TAC THEN VECTOR_ARITH_TAC);;
+
+let SYMETRIC_VECTORS_PROJECTION_ON_AXIS = prove
+ (`!u v w. norm w = &1 ==> symetric_vectors_wrt u v w ==> u dot w = v dot w`,
+ REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN (GEQ_IMP SYMETRIC_VECTORS_WRT_SYM) ASSUME_TAC
+ THEN REPEAT (POP_ASSUM MP_TAC) THEN SIMP_TAC[symetric_vectors_wrt]
+ THEN REWRITE_TAC[IMP_IMP;GSYM CONJ_ASSOC;REAL_EQ_MUL_LCANCEL;
+ VECTOR_ARITH `u+v=w /\ v+u=x <=> u+v=x /\ x=w`;VECTOR_MUL_RCANCEL]
+ THEN REPEAT STRIP_TAC THEN ASM_SIMP_TAC[DOT_RZERO] THEN POP_ASSUM MP_TAC
+ THEN REAL_ARITH_TAC);;
+
+
+(** Additions to affine spaces *)
+
+let AFFINE_HULL_3_ZERO = prove
+ (`affine hull {vec 0, a, b} = {u % a + v % b | u IN (:real) /\ v IN(:real)}`,
+ REWRITE_TAC[AFFINE_HULL_3;UNIV;EXTENSION;IN_ELIM_THM;VECTOR_MUL_RZERO;
+ VECTOR_ADD_LID;REAL_ARITH `x+y+z=t <=> x=t-y-z:real`]
+ THEN MESON_TAC[]);;
+
+let LSUB_PLANE_EQ_RSUB_PLANE = prove
+ (`!p:plane. plane p ==> !x. x IN p ==> {y - x | y IN p} = {x - y | y IN p}`,
+ let COMMON_TAC t =
+ EXISTS_TAC t THEN CONJ_TAC THEN TRY VECTOR_ARITH_TAC
+ THEN SUBGOAL_TAC "" `(x':real^3) IN span b` [ASM SET_TAC []]
+ in
+ REWRITE_TAC[GSPEC;SETSPEC;FUN_EQ_THM] THEN REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN EXISTS_PLANE_BASIS STRIP_ASSUME_TAC
+ THEN ASM_CSQS_THEN PLANE_AS_SPAN ASSUME_TAC THEN EQ_TAC THENL [
+ ASM_CSQS_THEN PLANE_SPAN ASSUME_TAC THEN STRIP_TAC
+ THEN COMMON_TAC `x-x':point` THEN REWRITE_TAC[VECTOR_SUB]
+ THEN POP_ASSUM (ASSUME_TAC o MATCH_MP SPAN_NEG)
+ THEN ASSUM_LIST (FIRST o List.rev_map (CHANGED_TAC o SINGLE REWRITE_TAC))
+ THEN ASSUM_LIST SET_TAC;
+ STRIP_TAC THEN ASM_CSQS_THEN PLANE_SPAN ASSUME_TAC
+ THEN COMMON_TAC `x+x':point` THEN ASSUM_LIST SET_TAC]);;
+
+let EXISTS_MULTIPLE_OF_NORMAL_IN_PLANE = prove
+ (`!p:plane. plane p ==> !n. is_normal_to_plane n p ==> ?a. a % unit n IN p`,
+ REPEAT STRIP_TAC
+ THEN ASM_CSQ_THEN PLANE_NON_EMPTY (STRIP_ASSUME_TAC
+ o REWRITE_RULE[GSYM MEMBER_NOT_EMPTY])
+ THEN ASM_CSQ_THEN EXISTS_PLANE_BASIS STRIP_ASSUME_TAC
+ THEN ASM_CSQS_THEN (GSYM PLANE_SPAN) (fun x ->
+ ASM_CSQS_THEN PLANE_DECOMP_DOT_NORMAL
+ (MP_TAC o SPEC `x:point` o REWRITE_RULE[x]))
+ THEN ASM_SIMP_TAC[LSUB_PLANE_EQ_RSUB_PLANE] THEN REWRITE_TAC[IN_ELIM_THM]
+ THEN STRIP_TAC THEN EXISTS_TAC `(x:point) dot unit n`
+ THEN MATCH_MP_TAC (MESON[] `(pt':point) IN p /\ pt' = y ==> y IN p`)
+ THEN ASM_REWRITE_TAC[] THEN REPLICATE_TAC 2 (POP_ASSUM MP_TAC)
+ THEN VECTOR_ARITH_TAC);;
+
+(* Additional functions *)
+let map_triple = new_definition `map_triple (f:A->B) (x1,x2,x3) = (f x1,f x2,f x3)`;;
--
1.7.1