geometrical_optics.ml

   1 (* ========================================================================= *)
   2 (*                                                                           *)
   3 (*                        Geometrical Optics Library                         *)
   4 (*                                                                           *)
   5 (*        (c) Copyright, Muhammad Umair Siddique, Vincent Aravantinos, 2012  *)
   6 (*                       Hardware Verification Group,                        *)
   7 (*                       Concordia University                                *)
   8 (*                                                                           *)
   9 (*     Contact: <muh_sidd@ece.concordia.ca>, <vincent@ece.concordia.ca>      *)
  10 (*                                                                           *)
  11 (* Last update: Feb 16, 2012                                                 *)
  12 (*                                                                           *)
  13 (* ========================================================================= *)
  14
  15
  16 needs "Multivariate/make_complex.ml";;
  17 needs "utils.ml";;
  18 (* ------------------------------------------------------------------------- *)
  19 (* Preliminaries                                                             *)
  20 (* ------------------------------------------------------------------------- *)
  21
  22 let LET_PAIR = prove
  23   (`!P t. (let x,y = t in P x y) = P (FST t) (SND t)`,
  24   ONCE_REWRITE_TAC[GSYM PAIR] THEN PURE_REWRITE_TAC[LET_DEF;LET_END_DEF;FST;SND]
  25   THEN REWRITE_TAC[]);;
  26
  27 let ALL_REVERSE = prove
  28   (`!l P. ALL P (REVERSE l) = ALL P l`,
  29   LIST_INDUCT_TAC THEN ASM_REWRITE_TAC[REVERSE;ALL;ALL_APPEND] THEN ITAUT_TAC);;
  30
  31 let BOUNDED_VECTOR_COMPONENT = prove
  32   (`!v:real^N. (?B:real^N. !i. 1 <= i /\ i <= dimindex(:N) ==> v$i <= B$i)
  33       <=> (?b. !i. 1 <= i /\ i <= dimindex(:N) ==> v$i <= b)`,
  34   REPEAT (STRIP_TAC ORELSE EQ_TAC) THENL [
  35     MESON_TAC[REAL_LE_TRANS;REAL_ABS_LE;COMPONENT_LE_INFNORM];
  36     EXISTS_TAC `b % vec 1:real^N`
  37     THEN ASM_REWRITE_TAC[VEC_COMPONENT;VECTOR_MUL_COMPONENT;REAL_MUL_RID]
  38   ]);;
  39
  40 let BOUNDED_MATRIX_BOUNDED_VECTOR_MUL = prove
  41   (`!M:num->real^M^N.
  42       (?b. !n i j. 1 <= i /\ i <= dimindex(:N) /\ 1 <= j /\ j <= dimindex(:M)
  43           ==> abs ((M n)$i$j) <= b)
  44       ==> (!v:real^M.
  45         ?b. !n i. 1 <= i /\ i <= dimindex(:N) ==> abs ((M n ** v)$i) <= b)`,
  46   REPEAT STRIP_TAC
  47   THEN SIMP_TAC[MATRIX_VECTOR_MUL_COMPONENT;dot]
  48   THEN Pa.EXISTS_TAC "sum (1..dimindex(:M)) (\i. b * abs (v$i))"
  49   THEN REPEAT GEN_TAC THEN DISCH_TAC THEN MATCH_MP_TAC SUM_ABS_LE
  50   THEN REWRITE_TAC[FINITE_NUMSEG;IN_NUMSEG;REAL_ABS_MUL]
  51   THEN GEN_TAC THEN Pa.ASM_CASES_TAC "v$x = &0"
  52   THEN ASM_REWRITE_TAC[REAL_ABS_0;REAL_MUL_RZERO;REAL_LE_REFL]
  53   THEN MP_REWRITE_TAC REAL_LE_RMUL_EQ
  54   THEN ASM_MESON_TAC[REAL_LT_LE;REAL_ABS_ZERO;REAL_ABS_POS]);;
  55
  56
  57 (* ------------------------------------------------------------------------- *)
  58 (* System formalization                                                      *)
  59 (* ------------------------------------------------------------------------- *)
  60
  61 new_type_abbrev ("free_space",`:real#real`);;
  62
  63 let FORALL_FREE_SPACE_THM = prove
  64   (`!P. (!fs:free_space. P fs) <=> !n d. P (n,d)`,
  65   REWRITE_TAC[FORALL_PAIR_THM]);;
  66
  67 let EXISTS_FREE_SPACE_THM = prove
  68   (`!P. (?fs:free_space. P fs) <=> ?n d. P (n,d)`,
  69   REWRITE_TAC[EXISTS_PAIR_THM]);;
  70
  71 let interface_IND,interface_REC = define_type
  72   "interface = plane | spherical real";;
  73
  74 let interface_kind_IND,interface_kind_REC = define_type
  75   "interface_kind = transmitted | reflected";;
  76
  77 new_type_abbrev("optical_component",`:free_space#interface#interface_kind`);;
  78
  79 let FORALL_OPTICAL_COMPONENT_THM = prove
  80   (`!P. (!c:optical_component. P c) <=> !fs i ik. P (fs,i,ik)`,
  81   REWRITE_TAC[FORALL_PAIR_THM]);;
  82
  83 let EXISTS_OPTICAL_COMPONENT_THM = prove
  84   (`!P. (?c:optical_component. P c) <=> ?fs i ik. P (fs,i,ik)`,
  85   REWRITE_TAC[EXISTS_PAIR_THM]);;
  86
  87 new_type_abbrev("optical_system",`:optical_component list#free_space`);;
  88
  89 let FORALL_OPTICAL_SYSTEM_THM = prove
  90   (`!P. (!os:optical_system. P os) <=> !cs fs. P (cs,fs)`,
  91   REWRITE_TAC[FORALL_PAIR_THM]);;
  92
  93 let is_valid_free_space = define
  94   `is_valid_free_space ((n,d):free_space) <=> &0 < n /\ &0 <= d`;;
  95
  96 let is_valid_interface = define
  97   `(is_valid_interface (plane) <=> T) /\
  98    (is_valid_interface (spherical R) <=> ~(&0 = R))`;;
  99
 100 let is_valid_optical_component = new_definition
 101   `is_valid_optical_component (fs,i,ik) <=>
 102     is_valid_free_space fs /\ is_valid_interface i`;;
 103
 104 let is_valid_optical_system = new_definition
 105   `is_valid_optical_system (cs,fs) <=>
 106     ALL is_valid_optical_component cs /\ is_valid_free_space fs`;;
 107
 108 let head_index = define
 109   `head_index ([],(n,d)) = n /\
 110   head_index (CONS ((n,d),i) cs, (nt,dt)) = n`;;
 111
 112 let HEAD_INDEX_VALID_OPTICAL_SYSTEM = prove
 113   (`!os:optical_system. is_valid_optical_system os ==> &0 < head_index os`,
 114   REWRITE_TAC[FORALL_OPTICAL_SYSTEM_THM;is_valid_optical_system]
 115   THEN MATCH_MP_TAC list_INDUCT
 116   THEN SIMP_TAC[FORALL_OPTICAL_COMPONENT_THM;FORALL_FREE_SPACE_THM;
 117     head_index;is_valid_free_space;ALL;is_valid_optical_component]);;
 118
 119 let HEAD_INDEX_APPEND = prove
 120   (`!cs1:optical_component list cs2 fs a.
 121       head_index (APPEND cs1 cs2, fs)
 122         = head_index (cs1, (head_index (cs2,fs),a))`,
 123     MATCH_MP_TAC list_INDUCT
 124     THEN REWRITE_TAC[APPEND;head_index;FORALL_OPTICAL_COMPONENT_THM;
 125       FORALL_FREE_SPACE_THM]);;
 126
 127 let HEAD_INDEX_CONS_ALT = prove
 128   (`!fs i cs:optical_component list fs'.
 129       head_index (CONS (fs,i) cs,fs')
 130         = head_index ([]:optical_component list,fs)`,
 131   REWRITE_TAC[FORALL_FREE_SPACE_THM;head_index]);;
 132
 133 let HEAD_INDEX_CONS_EQ = prove
 134   (`!n1 n2 d1 d2 i1 i2 cs1 cs2 fs1 fs2.
 135       head_index (CONS ((n1,d1),i1) cs1,fs1)
 136         = head_index (CONS ((n2,d2),i2) cs2,fs2)
 137           <=> n1 = n2`,
 138   REWRITE_TAC[FORALL_PAIR_THM;head_index;PAIR_EQ]);;
 139
 140 let HEAD_INDEX_CONS_EQ_ALT = prove
 141   (`!c1 c2 cs1 cs2 fs1 fs2.
 142       c1 = c2 ==> head_index (CONS c1 cs1,fs1) = head_index (CONS c2 cs2,fs2)`,
 143   SIMP_TAC[FORALL_PAIR_THM;head_index;PAIR_EQ]);;
 144
 145 (* ------------------------------------------------------------------------- *)
 146 (* Ray formalization                                                         *)
 147 (* ------------------------------------------------------------------------- *)
 148
 149 new_type_abbrev("single_ray",`:real#real`);;
 150
 151 let FORALL_SINGLE_RAY_THM = prove
 152   (`!P. (!sr:single_ray. P sr) <=> !y theta. P (y,theta)`,
 153   REWRITE_TAC[FORALL_PAIR_THM]);;
 154
 155 new_type_abbrev("ray",
 156   `:single_ray # single_ray # (single_ray # single_ray) list`);;
 157
 158 let FORALL_RAY_THM = prove
 159   (`!P. (!r:ray. P r) <=> !sr1 sr2 srs. P (sr1,sr2,srs)`,
 160   REWRITE_TAC[FORALL_PAIR_THM]);;
 161
 162 let is_valid_ray_in_free_space = define
 163   `is_valid_ray_in_free_space (y0,theta0) (y1,theta1) ((n0,d):free_space) <=>
 164     y1 = y0 + d * theta0 /\ theta0 = theta1`;;
 165
 166 let is_valid_ray_at_interface = define
 167   `(is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 plane transmitted
 168     <=> y1 = y0 /\ n0 * theta0 =  n1 * theta1) /\
 169   (is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 (spherical R)
 170     transmitted <=>
 171       let phi_i = theta0 + y1 / R and phi_t = theta1 + y1 / R in
 172       y1 = y0 /\ n0 * phi_i = n1 * phi_t) /\
 173   (is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 plane reflected <=>
 174     y1 = y0 /\ n0 * theta0 = n0 * theta1) /\
 175   (is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 (spherical R)
 176     reflected <=>
 177       let phi = y1 / R - theta0  in
 178       y1 = y0 /\ theta1 = --(theta0 + &2 * phi))`;;
 179
 180 let last_single_ray = define
 181   `last_single_ray ((sr1,sr2,[]):ray) = sr2 /\
 182   last_single_ray ((sr1,sr2,CONS sr3 t):ray) = SND (LAST (CONS sr3 t))`;;
 183
 184 let fst_single_ray = define
 185   `fst_single_ray ((sr1,sr2,rs):ray )= sr1`;;
 186
 187
 188 let EXISTS_IS_VALID_RAY_IN_SYSTEM = prove(
 189   `?f. !sr1 sr2 h h' fs cs rs i ik y0 theta0 y1 theta1 y2 theta2 y3 theta3 n d n' d'.
 190   (f (sr1,sr2,[]) (CONS h cs,fs) <=> F) /\
 191   (f (sr1,sr2, CONS h' rs) ([],fs) <=> F) /\
 192   (f ((y0,theta0),(y1,theta1), []) ([],(n,d)) <=>
 193     is_valid_ray_in_free_space (y0,theta0) (y1,theta1) (n,d)) /\
 194   (f ((y0,theta0),(y1,theta1), CONS ((y2,theta2),(y3,theta3)) rs)
 195     ((CONS ((n',d'),i,ik) cs),(n,d)) <=>
 196   (is_valid_ray_in_free_space (y0,theta0) (y1,theta1) (n',d') /\
 197   (is_valid_ray_at_interface (y1,theta1) (y2,theta2) n' (head_index(cs,(n,d))) i ik))
 198   /\ f ((y2,theta2),(y3,theta3),rs) (cs,(n,d)))`,
 199   SUBGOAL_THEN `!y0 theta0 y1 theta1 n d. ?f1. (f1 ([]:(free_space#interface#interface_kind)list) = is_valid_ray_in_free_space (y0,theta0) (y1,theta1) (n,d)
 200     /\ !h t. f1 (CONS h t) = F)` (STRIP_ASSUME_TAC o REWRITE_RULE[SKOLEM_THM])
 201   THENL [
 202     REPEAT GEN_TAC THEN MATCH_ACCEPT_TAC list_RECURSION;
 203     SUBGOAL_THEN
 204       `!y2 theta2 self y0 theta0 y1 theta1 (fs:free_space). ?f2. (f2 [] = F /\ !h t.
 205       f2 (CONS h t) <=>
 206       let (n,d),i,ik = h in
 207       is_valid_ray_in_free_space (y0,theta0) (y1,theta1) (n,d) /\
 208       is_valid_ray_at_interface (y1,theta1) (y2,theta2) n (head_index(t,fs)) i ik /\
 209       self t)`
 210     (STRIP_ASSUME_TAC o REWRITE_RULE[SKOLEM_THM])
 211     THENL [
 212       REPEAT GEN_TAC THEN MATCH_ACCEPT_TAC list_RECURSION;
 213       SUBGOAL_THEN
 214         `!n:real d:real. ?g. (g [] = (\y0:real theta0:real y1:real theta1:real. f1 y0 theta0 y1 theta1 n d) /\ !h t. g (CONS h t) =
 215         \y0:real theta0:real y1:real theta1:real.
 216         let (y2,theta2),(y3,theta3) = (h:single_ray#single_ray) in
 217         (f2 y2 theta2 (g t y2 theta2 y3 theta3) y0 theta0 y1 theta1 (n,d)):(free_space#interface#interface_kind)list ->bool)`
 218         (STRIP_ASSUME_TAC o REWRITE_RULE [SKOLEM_THM])
 219       THENL [
 220         REPEAT GEN_TAC THEN MATCH_ACCEPT_TAC list_RECURSION;
 221         EXISTS_TAC `\(((y0,theta0),(y1,theta1),rs):ray) ((cs,(n,d)):optical_system). (g n d rs y0 theta0 y1 theta1 cs):bool`
 222         THEN ASM_REWRITE_TAC[LET_DEF;LET_END_DEF;CONJ_ACI;FORALL_PAIR_THM]
 223       ]
 224     ]
 225   ]);;
 226
 227 let is_valid_ray_in_system = new_specification ["is_valid_ray_in_system"] EXISTS_IS_VALID_RAY_IN_SYSTEM;;
 228
 229
 230 (* ------------------------------------------------------------------------- *)
 231 (* 2x2 Matrices                                                              *)
 232 (* ------------------------------------------------------------------------- *)
 233
 234 let mat_pow = define
 235   `mat_pow (M:real^N^N) (0:num) :real^N^N = mat 1 /\
 236    mat_pow M (SUC n) = M ** mat_pow M n`;;
 237
 238 overload_interface("pow",`mat_pow:real^N^N->num->real^N^N`);;
 239
 240 let MAT_POW_ALT = prove
 241   (`!n M. mat_pow M 0 = mat 1 /\ mat_pow M (SUC n) = mat_pow M n ** M`,
 242   REWRITE_TAC[mat_pow] THEN INDUCT_TAC
 243   THEN REWRITE_TAC[mat_pow;MATRIX_MUL_LID;MATRIX_MUL_RID]
 244   THEN ASM_REWRITE_TAC[GSYM MATRIX_MUL_ASSOC]);;
 245
 246 let mat2x2 = new_definition
 247   `mat2x2 A B C D :A^2^2 = vector[ vector[A;B]; vector[C;D] ]`;;
 248
 249 let COMMON_TAC xs =
 250   SIMP_TAC (xs @ [mat2x2;CART_EQ;VECTOR_2;DIMINDEX_2;FORALL_2;DOT_2;SUM_2]);;
 251
 252 let MAT2X2_VECTOR_MUL = prove
 253   (`!A B C D X Y.
 254     mat2x2 A B C D ** vector [X;Y] = vector[A*X+B*Y;C*X+D*Y]`,
 255   COMMON_TAC [MATRIX_VECTOR_MUL_COMPONENT]);;
 256
 257 let MAT2X2_VECTOR_MUL_GEN = prove
 258   (`!(M:real^2^2) (X:real) (Y:real).
 259     M ** vector [X;Y] = vector[(M$1$1)*X+ (M$1$2)*Y;(M$2$1*X)+ (M$2$2)*Y]`,
 260   COMMON_TAC [MATRIX_VECTOR_MUL_COMPONENT]);;
 261
 262 let MAT2X2_EQ = prove
 263   (`!A B C D P Q R S. mat2x2 A B C D = mat2x2 P Q R S <=>
 264     A = P /\ B = Q /\ C = R /\ D = S`,
 265   COMMON_TAC [] THEN CONV_TAC TAUT);;
 266
 267 let VECTOR2_EQ = prove
 268   (`!x y z t. vector [x;y] :A^2 = vector [z;t] <=> x=z /\ y=t`,
 269   COMMON_TAC[]);;
 270
 271 let MAT2X2_MUL = prove
 272   (`!A B C D P Q R S.
 273     mat2x2 A B C D ** mat2x2 P Q R S =
 274       mat2x2 (A*P + B*R) (A*Q + B*S) (C*P + D*R) (C*Q + D*S)`,
 275   COMMON_TAC[matrix_mul;MATRIX_VECTOR_MUL_COMPONENT;LAMBDA_BETA]);;
 276
 277 let MAT2X2_MAT1 = prove
 278   (`mat2x2 (&1) (&0) (&0) (&1) = mat 1`,
 279   COMMON_TAC[mat;LAMBDA_BETA] THEN CONV_TAC NUM_REDUCE_CONV);;
 280
 281 let MAT2X2_LID = prove
 282   (`!m:real^2^2. mat2x2 (&1) (&0) (&0) (&1) ** m = m`,
 283   REWRITE_TAC[MAT2X2_MAT1;MATRIX_MUL_LID]);;
 284
 285 let MAT2X2_RID = prove
 286   (`!m:real^2^2. m ** mat2x2 (&1) (&0) (&0) (&1) = m`,
 287   REWRITE_TAC[MAT2X2_MAT1;MATRIX_MUL_RID]);;
 288
 289 let MAT2X2_ID = CONJ MAT2X2_LID MAT2X2_RID;;
 290
 291 let MAT2X2_SCALAR_MUL = prove
 292   (`!A B C D k. k %% mat2x2 A B C D = mat2x2 (k*A) (k*B) (k*C) (k*D)`,
 293   COMMON_TAC [MATRIX_CMUL_COMPONENT]);;
 294
 295 let MAT2X2_SCALAR_MUL_ASSOC = prove
 296   (`!A B C D L M N P k.
 297     mat2x2 L M N P ** (k %% mat2x2 A B C D)
 298       = k %% (mat2x2 L M N P ** mat2x2 A B C D)`,
 299   SIMP_TAC[MAT2X2_SCALAR_MUL;MAT2X2_MUL;MAT2X2_EQ] THEN CONV_TAC REAL_FIELD);;
 300
 301 let MAT2X2_EQ_SCALAR_LMUL = prove
 302   (`!A B C D P Q R S. ~(k = &0) ==>
 303     (k %% mat2x2 A B C D = k %% mat2x2 P Q R S
 304       <=> A = P /\ B = Q /\ C = R /\ D = S)`,
 305   COMMON_TAC[MAT2X2_SCALAR_MUL] THEN CONV_TAC REAL_FIELD);;
 306
 307 let MAT2X2_DET = prove
 308   (`!A B C D. det (mat2x2 A B C D) = A * D - B * C`,
 309   COMMON_TAC[DET_2]);;
 310
 311 let FORALL_MATRIX_2X2 = prove
 312   (`!P. (!M:real^2^2. P M) <=> ! a b c d.P (mat2x2 a b c d )`,
 313  REWRITE_TAC[FORALL_VECTOR_2;mat2x2] );;
 314
 315 let MAT2X2_VECTOR_MUL_ALT = prove
 316  (`!(M:real^2^2) (x:real) y x' y'.  (vector[x';y'] = M ** vector[x;y]) =
 317                    (x' =  (M ** vector[x;y])$1 /\ y' = (M ** vector[x;y])$2)`,
 318 COMMON_TAC [ FORALL_MATRIX_2X2;VECTOR2_EQ;MAT2X2_VECTOR_MUL]);;
 319
 320
 321 let MAT2_MUL = prove (
 322  ` !(M:real^2^2). M**M =  mat2x2((M$1$1 * M$1$1) +  (M$1$2 * M$2$1))
 323                                 ((M$1$1 * M$1$2) +  (M$1$2 * M$2$2))
 324                                 ((M$2$1 * M$1$1) +  (M$2$2 * M$2$1))
 325                                 ((M$2$1 * M$1$2) +  (M$2$2 * M$2$2))`,
 326  COMMON_TAC[matrix_mul;MATRIX_VECTOR_MUL_COMPONENT;LAMBDA_BETA]);;
 327
 328
 329 let MAT_POW_ADD = prove (
 330 `!(M:real^a^a) m n. (M pow (m + n)) = ((M pow m)**(M pow n))`,
 331 GEN_TAC THEN INDUCT_TAC THEN
 332 REWRITE_TAC[mat_pow;ADD_CLAUSES;MATRIX_MUL_LID] THEN
 333 REWRITE_TAC[mat_pow;ADD_CLAUSES;MATRIX_MUL_LID] THEN
 334 REWRITE_TAC[GSYM MATRIX_MUL_ASSOC] THEN ASM_SIMP_TAC[]);;
 335
 336
 337 let MAT2X2_POW2_POW = prove (
 338  ` !N. ((mat2x2 A B C D) pow 2) pow N = (mat2x2 A B C D) pow (2*N)`,
 339 INDUCT_TAC THEN REWRITE_TAC[mat_pow;MULT_0] THEN
 340 REWRITE_TAC[mat_pow;ADD1;LEFT_ADD_DISTRIB;MULT_CLAUSES;MAT_POW_ADD] THEN
 341 POP_ASSUM MP_TAC THEN SIMP_TAC[EQ_SYM_EQ]);;
 342
 343 let MAT2X2_POW2_POW_GEN = prove
 344 (` !n (M:real^k^k). (M pow 2) pow n = M pow (2*n)`,
 345 INDUCT_TAC THEN REWRITE_TAC[mat_pow;MULT_0] THEN
 346 REWRITE_TAC[mat_pow;ADD1;LEFT_ADD_DISTRIB;MULT_CLAUSES;MAT_POW_ADD] THEN
 347 POP_ASSUM MP_TAC THEN SIMP_TAC[EQ_SYM_EQ]);;
 348
 349 let MAT_POW2 = prove (`
 350   !(M:real^a^a). M pow 2 = M**M`,
 351 REWRITE_TAC[mat_pow;TWO;ONE] THEN
 352 REWRITE_TAC[GSYM ONE; MATRIX_MUL_RID]);;
 353
 354 let TRACE_MAT2X2 = prove( `!A B C D. trace (mat2x2 A B C D) = A+D`,
 355  SIMP_TAC[SUM_2;trace;mat2x2;CART_EQ;DIMINDEX_2;FORALL_2;VECTOR_2]);;
 356
 357 let DET_MAT2X2 = prove( `!A B C D. det (mat2x2 A B C D) = A*D - B*C`,
 358  SIMP_TAC[SUM_2;DET_2;mat2x2;CART_EQ;DIMINDEX_2;FORALL_2;VECTOR_2]);;
 359
 360
 361 (* ------------------------------------------------------------------------- *)
 362 (* Matrix formalism                                                          *)
 363 (* ------------------------------------------------------------------------- *)
 364
 365 let free_space_matrix = new_definition
 366   `free_space_matrix (n,d) = mat2x2 (&1) d (&0) (&1)`;;
 367
 368 let interface_matrix = define
 369   `interface_matrix n0 n1 plane transmitted = mat2x2 (&1) (&0) (&0) (n0/n1) /\
 370   interface_matrix n0 n1 (spherical R) transmitted
 371     = mat2x2 (&1) (&0) ((n0-n1)/n1 * &1/R) (n0/n1)/\
 372   interface_matrix n0 n1 plane reflected = mat2x2 (&1) (&0) (&0) (&1) /\
 373   interface_matrix n0 n1 (spherical R) reflected
 374     = mat2x2 (&1) (&0) (-- &2/R) (&1)`;;
 375
 376
 377 (*-------------------------------------------------------------------------
 378 ---Verification of matrix relationships for free space and interfaces------
 379 -- ** PI= Plane Interface *** SI = Spherical Interface ** FS = Free Space--
 380 --------------------------------------------------------------------------*)
 381
 382 let common_prove t = prove (t,
 383   SIMP_TAC[free_space_matrix;FORALL_FREE_SPACE_THM;MAT2X2_VECTOR_MUL;VECTOR2_EQ;
 384     LET_END_DEF;is_valid_ray_in_free_space;REAL_MUL_LZERO;REAL_MUL_LID;LET_DEF;
 385     REAL_ADD_LID;interface_matrix;is_valid_ray_at_interface;REAL_ADD_RID]
 386   THEN CONV_TAC REAL_FIELD);;
 387
 388 let FS_MATRIX = common_prove
 389   `!fs y0 theta0 y1 theta1.
 390     is_valid_free_space fs
 391     /\ is_valid_ray_in_free_space (y0,theta0) (y1,theta1) fs
 392     ==> vector[y1;theta1] = free_space_matrix fs ** vector[y0;theta0]`;;
 393
 394 let PI_MATRIX_TRANSMITTED = common_prove
 395   `!n0 n1 y0 theta0 y1 theta1.
 396     is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 plane transmitted
 397     /\ &0 < n0 /\ &0 < n1
 398       ==> vector [y1;theta1]
 399         = interface_matrix n0 n1 plane transmitted ** vector [y0; theta0]`;;
 400
 401 let PI_MATRIX_REFLECTED = common_prove
 402   `!n0 y0 theta0 y1 theta1.
 403     is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 plane reflected
 404     /\ &0 < n0 /\ &0 < n1
 405       ==> vector [y1;theta1]
 406         = interface_matrix n0 n1 plane reflected ** vector [y0;theta0]`;;
 407
 408 let PI_MATRIX = prove
 409   (`!ik n0 n1 y0 theta0 y1 theta1.
 410      is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 plane ik
 411      /\ &0 < n0 /\ &0 < n1
 412        ==> vector [y1;theta1]
 413         = interface_matrix n0 n1 plane ik ** vector [y0;theta0]`,
 414   MATCH_MP_TAC interface_kind_IND
 415   THEN SIMP_TAC[PI_MATRIX_TRANSMITTED;PI_MATRIX_REFLECTED]);;
 416
 417 let SI_MATRIX_TRANSMITTED = common_prove
 418   `!n0 n1 y0 theta0 y1 theta1 R.
 419    is_valid_interface (spherical R) /\ &0 < n0 /\  &0 < n1
 420    /\ is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 (spherical R)
 421     transmitted
 422      ==> vector [y1;theta1]
 423       = interface_matrix n0 n1 (spherical R) transmitted ** vector [y0;theta0]`;;
 424
 425 let SI_MATRIX_REFLECTED = common_prove
 426   `!n0 y0 theta0 y1 theta1 R.
 427    is_valid_interface (spherical R) /\ &0 < n0
 428    /\ is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 (spherical R)
 429     reflected
 430      ==> vector [y1;theta1]
 431       = interface_matrix n0 n1 (spherical R) reflected ** vector [y0; theta0]`;;
 432
 433 let SI_MATRIX = prove
 434  (`!ik n0 n1 y0 theta0 y1 theta1.
 435   is_valid_interface (spherical R)
 436   /\ is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 (spherical R) ik
 437   /\ &0 < n0 /\ &0 < n1
 438    ==> vector [y1; theta1] = interface_matrix n0 n1 (spherical R) ik ** vector [y0; theta0]`,
 439   MATCH_MP_TAC interface_kind_IND
 440   THEN SIMP_TAC[SI_MATRIX_TRANSMITTED;SI_MATRIX_REFLECTED ]);;
 441
 442
 443 (*-------General Relationship for any Interface i --------------------*)
 444
 445
 446 let INTERFACE_MATRIX = prove
 447   (`!n0 n1 y0 theta0 y1 theta1 i ik.
 448   is_valid_interface i /\ &0 < n0 /\ &0 < n1 /\
 449   is_valid_ray_at_interface (y0,theta0) (y1,theta1) n0 n1 i ik
 450     ==> vector [y1;theta1] = interface_matrix n0 n1 i ik ** vector [y0;theta0]`,
 451   REPEAT (MATCH_MP_TAC interface_IND ORELSE GEN_TAC)
 452   THEN MESON_TAC[PI_MATRIX;SI_MATRIX]);;
 453
 454
 455 (*-----------------------------------------------*)
 456
 457 let system_composition = define
 458   `system_composition ([],fs) = free_space_matrix fs /\
 459    system_composition (CONS ((n',d'),i,ik) cs,n,d) =
 460      system_composition (cs,n,d)
 461      ** interface_matrix n' (head_index (cs,n,d)) i ik
 462      ** free_space_matrix (n',d')`;;
 463
 464 let IS_VALID_OPTICAL_SYSTEM_REC = prove
 465   (`!fs fs' i ik t.
 466     is_valid_optical_system ([],fs) = is_valid_free_space fs /\
 467     (is_valid_optical_system (CONS (fs,i,ik) t,fs') <=> is_valid_free_space fs
 468     /\ is_valid_interface i /\ is_valid_optical_system (t,fs'))`,
 469   REWRITE_TAC[is_valid_optical_system;ALL;is_valid_optical_component]
 470   THEN MESON_TAC[]);;
 471
 472 let LAST_SINGLE_RAY_REC = prove
 473   (`!sr1 sr2 sr3 sr4 t.
 474     last_single_ray (sr1,sr2,[]) = sr2 /\
 475     last_single_ray (sr1,sr2,CONS (sr3,sr4) t) = last_single_ray (sr3,sr4,t)`,
 476   REPEAT (LIST_INDUCT_TAC ORELSE GEN_TAC)
 477   THEN SIMP_TAC[last_single_ray;LAST;NOT_CONS_NIL]);;
 478
 479 let IS_VALID_OPTICAL_SYSTEM_HEAD_INDEX = prove
 480   (`!os:optical_system. is_valid_optical_system os ==> &0 < head_index os`,
 481   REWRITE_TAC[FORALL_OPTICAL_SYSTEM_THM] THEN MATCH_MP_TAC list_INDUCT
 482   THEN SIMP_TAC[FORALL_FREE_SPACE_THM;IS_VALID_OPTICAL_SYSTEM_REC;ALL;
 483     is_valid_free_space;head_index;FORALL_OPTICAL_COMPONENT_THM]);;
 484
 485 let SYSTEM_MATRIX = prove
 486   (`!sys ray.
 487     is_valid_optical_system sys /\ is_valid_ray_in_system ray sys ==>
 488       let ((y0,theta0),(y1,theta1),rs) = ray in
 489       let (y_n,theta_n) = last_single_ray ray in
 490         vector [y_n;theta_n] = system_composition sys ** vector [y0;theta0]`,
 491   REWRITE_TAC[FORALL_OPTICAL_SYSTEM_THM;FORALL_RAY_THM]
 492   THEN MATCH_MP_TAC list_INDUCT THEN CONJ_TAC
 493   THEN REWRITE_TAC[FORALL_SINGLE_RAY_THM;FORALL_OPTICAL_COMPONENT_THM]
 494   THENL [ALL_TAC;REPEAT GEN_TAC THEN DISCH_TAC]
 495   THEN REPEAT (MATCH_MP_TAC list_INDUCT ORELSE GEN_TAC) THEN CONJ_TAC
 496   THEN REWRITE_TAC[IS_VALID_OPTICAL_SYSTEM_REC;is_valid_ray_in_system;ALL;
 497     system_composition;head_index]
 498   THENL [
 499     REWRITE_TAC[last_single_ray;LET_DEF;LET_END_DEF;FS_MATRIX];
 500     let lemma = prove(`!P. (!a. P a) <=> !y1 th1 y2 th2. P((y1,th1),y2,th2)`,
 501       REWRITE_TAC[FORALL_PAIR_THM]) in
 502     REWRITE_TAC[lemma;is_valid_ray_in_system;GSYM MATRIX_VECTOR_MUL_ASSOC;
 503       LAST_SINGLE_RAY_REC;is_valid_free_space]
 504     THEN REPEAT STRIP_TAC THEN REPEAT LET_TAC
 505     THEN ASSUM_LIST (REWRITE_TAC o map (REWRITE_RULE[PAIR_EQ] o GSYM))
 506     THEN REPEAT_GTCL IMP_RES_THEN (wrap REWRITE_TAC)
 507       (REWRITE_RULE[IMP_CONJ;is_valid_free_space;FORALL_FREE_SPACE_THM]
 508         (GSYM FS_MATRIX))
 509     THEN IMP_RES_THEN ASSUME_TAC IS_VALID_OPTICAL_SYSTEM_HEAD_INDEX
 510     THEN REPEAT_GTCL IMP_RES_THEN (wrap REWRITE_TAC)
 511       (REWRITE_RULE[IMP_CONJ] (GSYM INTERFACE_MATRIX))
 512     THEN EVERY_ASSUM (REPEAT_GTCL IMP_RES_THEN MP_TAC o REWRITE_RULE[IMP_CONJ])
 513     THEN ASM_REWRITE_TAC[LET_DEF;LET_END_DEF] THEN SIMP_TAC[]
 514   ]);;
 515
 516