Update from HH
[Flyspeck/.git] / text_formalization / nonlinear / nonlin_def.hl
1 (* ========================================================================== *)
2 (* FLYSPECK - BOOK FORMALIZATION                                              *)
3 (*                                                                            *)
4 (* Chapter: Nonlinear                                                  *)
5 (* Author: Thomas C. Hales                                                    *)
6 (* Date: 2012-11-30                                                           *)
7 (* ========================================================================== *)
8
9 (*
10 Combine definitions for nonlinear inequalities here.
11 *)
12
13
14 module Nonlin_def = struct
15
16 (* From lemma.hl *)
17
18 let NONLIN = new_definition `NONLIN = T`;;
19
20 let sqrt_x1 = define `sqrt_x1 x1 x2 x3 x4 x5 x6 = sqrt x1`;;
21
22 let sqrt_x2 = define `sqrt_x2 x1 x2 x3 x4 x5 x6 = sqrt x2`;;
23
24 let sqrt_x3 = define `sqrt_x3 x1 x2 x3 x4 x5 x6 = sqrt x3`;;
25
26 let sqrt_x4 = define `sqrt_x4 x1 x2 x3 x4 x5 x6 = sqrt x4`;;
27
28 let sqrt_x5 = define `sqrt_x5 x1 x2 x3 x4 x5 x6 = sqrt x5`;;
29
30 let sqrt_x6 = define `sqrt_x6 x1 x2 x3 x4 x5 x6 = sqrt x6`;;
31
32 let halfbump_x = new_definition `halfbump_x x = bump (sqrt x / &2)`;;
33
34 let halfbump_x1 = new_definition `halfbump_x1 x1 x2 x3 x4 x5 x6 = halfbump_x x1`;;
35
36 let halfbump_x4 = new_definition `halfbump_x4 x1 x2 x3 x4 x5 x6 = halfbump_x x4`;;
37
38 let unit6 = define `unit6 x1 x2 x3 x4 x5 x6 = &1`;;
39
40 let proj_x1 = define `proj_x1 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x1`;;
41
42 let proj_x2 = define `proj_x2 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x2`;;
43
44 let proj_x3 = define `proj_x3 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x3`;;
45
46 let proj_x4 = define `proj_x4 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x4`;;
47
48 let proj_x5 = define `proj_x5 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x5`;;
49
50 let proj_x6 = define `proj_x6 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x6`;;
51
52 let promote = define `promote f (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = f x1`;;
53
54 let unit0 = define `unit0 = &1`;; 
55
56 let pow1 = new_definition `pow1 y = y pow 1`;;
57
58 let pow2 = new_definition `pow2 y = y pow 2`;;
59
60 let pow3 = new_definition `pow3 y = y pow 3`;;
61
62 let pow4 = new_definition `pow4 y = y pow 4`;;
63
64 let promote_pow2 = new_definition `promote_pow2 x1 (x2:A) (x3:B) (x4:C) (x5:D) (x6:E) = x1 pow 2`;;
65
66 let promote_pow3 = new_definition `promote_pow3 x1 (x2:A) (x3:B) (x4:C) (x5:D) (x6:E) = x1 pow 3`;;
67
68 (*
69 let promote_pow3r = INST_TYPE [(`:real`,`:A`);(`:real`,`:B`);(`:real`,`:C`);(`:real`,`:D`);(`:real`,`:E`)] promote_pow3;;
70 *)
71
72 let compose6 = new_definition `compose6 f p1 p2 p3 p4 p5 p6 
73   (x1:real) (x2:real) (x3:real) (x4:real) (x5:real) (x6:real) =
74   (f:real->real->real->real->real->real->real)
75   (p1 x1 x2 x3 x4 x5 x6)
76     (p2 x1 x2 x3 x4 x5 x6)
77     (p3 x1 x2 x3 x4 x5 x6)
78     (p4 x1 x2 x3 x4 x5 x6)
79     (p5 x1 x2 x3 x4 x5 x6)
80     (p6 x1 x2 x3 x4 x5 x6)`;;
81
82 let scale6 = new_definition `scale6 f
83    (r:real)   (x1:real) (x2:real) (x3:real) (x4:real) (x5:real) (x6:real) = 
84    (f x1 x2 x3 x4 x5 x6) * r`;;
85
86
87
88 let quadratic_root_plus_curry = 
89   new_definition `quadratic_root_plus_curry a b c = quadratic_root_plus (a,b,c)`;;
90
91 (* ----  *)
92
93 let sq_pow2 = prove_by_refinement(
94   `!a x. a pow 2 <= x ==> (sqrt x * sqrt x = x)`,
95   (* {{{ proof *)
96   [
97   REWRITE_TAC[GSYM REAL_POW_2;SQRT_POW2];
98   MESON_TAC[REAL_LE_TRANS;Collect_geom.REAL_LE_SQUARE_POW];
99   ]);;
100   (* }}} *)
101
102 let sqrt2_sqrt2 = prove_by_refinement(
103   `sqrt2 * sqrt2 = &2`,
104   (* {{{ proof *)
105   [
106   REWRITE_TAC[Sphere.sqrt2];
107     MATCH_MP_TAC sq_pow2;
108   EXISTS_TAC`&0`;
109   REAL_ARITH_TAC;
110   ]);;
111   (* }}} *)
112
113 (* This is valid when a > 0 *)
114
115 let sol_euler_x_div_sqrtdelta =  new_definition 
116  `sol_euler_x_div_sqrtdelta x1 x2 x3 x4 x5 x6 = 
117   (let a = sqrt(x1*x2*x3) + sqrt( x1)*(x2 + x3 - x4)/ &2 + 
118      sqrt(x2)*(x1 + x3 - x5)/ &2 + sqrt(x3)*(x1 + x2 - x6)/ &2 in
119      (matan ((delta_x x1 x2 x3 x4 x5 x6)/(&4 * a pow 2 )))/( a))`;;
120
121 let sol_euler246_x_div_sqrtdelta =  new_definition 
122  `sol_euler246_x_div_sqrtdelta = rotate4 sol_euler_x_div_sqrtdelta`;;
123
124 let sol_euler345_x_div_sqrtdelta =  new_definition 
125  `sol_euler345_x_div_sqrtdelta = rotate5 sol_euler_x_div_sqrtdelta`;;
126
127 let sol_euler156_x_div_sqrtdelta =  new_definition 
128  `sol_euler156_x_div_sqrtdelta = rotate6 sol_euler_x_div_sqrtdelta`;;
129
130 (* This is a valid formula for dih_x / sqrt(delta) when delta >0 and d_x4 !=0 *)
131
132 let dih_x_div_sqrtdelta = `dih_x_div_sqrtdelta x1 x2 x3 x4 x5 x6 =
133   ( let d_x4 = delta_x4 x1 x2 x3 x4 x5 x6 in
134           let d = delta_x x1 x2 x3 x4 x5 x6 in
135           let pi_shift = if (d_x4 < &0) then pi else (&0) in
136             pi_shift + (sqrt(&4 * x1) / d_x4) * matan((&4 * x1 * d)/(d_x4 pow 2)))`;;
137
138 let dih_x_div_sqrtdelta_posbranch = new_definition
139  `dih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
140   ( let d_x4 = delta_x4 x1 x2 x3 x4 x5 x6 in
141           let d = delta_x x1 x2 x3 x4 x5 x6 in
142              (sqrt(&4 * x1) / d_x4) * matan((&4 * x1 * d)/(d_x4 pow 2)))`;;
143
144 let ldih_x_div_sqrtdelta_posbranch = new_definition
145  `ldih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
146    lfun(sqrt(x1) / &2) * dih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;
147
148 let ldih2_x_div_sqrtdelta_posbranch = new_definition
149  `ldih2_x_div_sqrtdelta_posbranch =  rotate2 ldih_x_div_sqrtdelta_posbranch`;;
150
151 let ldih3_x_div_sqrtdelta_posbranch = new_definition
152  `ldih3_x_div_sqrtdelta_posbranch  = rotate3 ldih_x_div_sqrtdelta_posbranch`;;
153
154 let ldih5_x_div_sqrtdelta_posbranch = new_definition
155  `ldih5_x_div_sqrtdelta_posbranch  = rotate5 ldih_x_div_sqrtdelta_posbranch`;;
156
157 let ldih6_x_div_sqrtdelta_posbranch = new_definition
158  `ldih6_x_div_sqrtdelta_posbranch  = rotate6 ldih_x_div_sqrtdelta_posbranch`;;
159
160 let dih3_x_div_sqrtdelta_posbranch = new_definition
161  `dih3_x_div_sqrtdelta_posbranch  = rotate3 dih_x_div_sqrtdelta_posbranch`;;
162
163 let dih4_x_div_sqrtdelta_posbranch = new_definition
164  `dih4_x_div_sqrtdelta_posbranch  = rotate4 dih_x_div_sqrtdelta_posbranch`;;
165
166 let dih5_x_div_sqrtdelta_posbranch = new_definition
167  `dih5_x_div_sqrtdelta_posbranch  = rotate5 dih_x_div_sqrtdelta_posbranch`;;
168
169 let rhazim_x_div_sqrtdelta_posbranch = new_definition
170  `rhazim_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
171    rho(sqrt(x1)) * dih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;
172
173 let rhazim2_x_div_sqrtdelta_posbranch = new_definition
174  `rhazim2_x_div_sqrtdelta_posbranch =  
175   rotate2 rhazim_x_div_sqrtdelta_posbranch`;;
176
177 let rhazim3_x_div_sqrtdelta_posbranch = new_definition
178  `rhazim3_x_div_sqrtdelta_posbranch  = 
179   rotate3 rhazim_x_div_sqrtdelta_posbranch`;;
180
181 let tau_residual_x  = new_definition
182  `tau_residual_x x1 x2 x3 x4 x5 x6 = 
183      rhazim_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6
184 +     rhazim2_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6
185 +     rhazim3_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;
186
187
188 let lmdih_x_div_sqrtdelta_posbranch = new_definition
189  `lmdih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
190    lmfun(sqrt(x1) / &2) * dih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;
191
192 let lmdih2_x_div_sqrtdelta_posbranch = new_definition 
193   `lmdih2_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
194          rotate2 lmdih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;
195
196 let lmdih3_x_div_sqrtdelta_posbranch = new_definition
197  `lmdih3_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
198   rotate3 lmdih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 `;;
199
200 let lmdih5_x_div_sqrtdelta_posbranch = new_definition
201  `lmdih5_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
202   rotate5 lmdih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 `;;
203
204 let lmdih6_x_div_sqrtdelta_posbranch = new_definition 
205   `lmdih6_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
206          rotate6 lmdih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;
207
208 let vol3f_sqrt2_lmplus = new_definition 
209   `vol3f_sqrt2_lmplus y1 y2 (y3:real) (y4:real) (y5:real) y6 = 
210     (&2 * mm1 / pi) *
211              (&2 * dih_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
212               &2 * dih2_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
213               &2 * dih6_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
214               dih3_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
215               dih4_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
216               dih5_y y1 y2 sqrt2 sqrt2 sqrt2 y6 - &3 * pi) -
217              (&8 * mm2 / pi) *
218              (
219               lfun (y2 / &2) * dih2_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
220               lfun (y6 / &2) * dih6_y y1 y2 sqrt2 sqrt2 sqrt2 y6)`;;
221
222 let vol3f_x_sqrt2_lmplus = new_definition
223   `vol3f_x_sqrt2_lmplus x1 x2 x3 x4 x5 x6 =
224    vol3f_sqrt2_lmplus (sqrt x1) (sqrt x2) (sqrt x3) (sqrt x4) (sqrt x5)      (sqrt x6)`;;
225
226 let vol3f_x_lfun = new_definition 
227  `vol3f_x_lfun x1 x2 (x3:real) (x4:real) (x5:real) x6 = vol3f (sqrt x1) (sqrt x2) (sqrt x6)  sqrt2 lfun `;;
228
229 let vol3_x_sqrt = new_definition
230  `vol3_x_sqrt x1 x2 (x3:real) (x4:real) (x5:real) x6  = vol_y sqrt2 sqrt2 sqrt2 (sqrt x1) (sqrt x2) (sqrt x6) `;;
231
232
233 let gamma23f = new_definition `gamma23f y1 y2 y3 y4 y5 y6 w1 w2 r f =
234       (gamma3f y1 y2 y6 r f / &w1 + gamma3f y1 y3 y5 r f / &w2
235       + (dih_y y1 y2 y3 y4 y5 y6 - dih_y y1 y2 r r r y6 - dih_y y1 y3 r r r y5) * (vol2r y1 r - vol2f y1 r f)/(&2 * pi)) `;;
236
237 let monomial = new_definition `monomial n1 n2 n3 n4 n5 n6 y1 y2 y3 y4 y5 y6 = 
238    (y1 pow n1) * (y2 pow n2) * (y3 pow n3) * (y4 pow n4) * (y5 pow n5) * (y6 pow n6)`;;
239
240 let arclength_x_234 = new_definition `arclength_x_234  (x1:real) (x2:real) (x3:real) (x4:real) (x5:real) (x6:real) = arclength (sqrt x2) (sqrt x3) (sqrt x4)`;;
241
242 let arclength_x_126 = new_definition `arclength_x_126  (x1:real) (x2:real) (x3:real) (x4:real) (x5:real) (x6:real) = arclength (sqrt x1) (sqrt x2) (sqrt x6)`;;
243
244 (* Functional Equation *)
245
246 let uni = new_definition `uni (f,x) x1 x2 x3 x4 x5 x6 = 
247   (f:A->B) ( x x1 x2 x3 x4 x5 x6)`;;
248
249 let constant6 = new_definition `constant6 c x1 x2 x3 x4 x5 x6 = c`;;
250
251 let promote3_to_6 = new_definition
252     `promote3_to_6 f x1 x2 x3 x4 x5 x6 = f x1 x2 x3`;;
253
254 let promote1_to_6 = new_definition
255   `promote1_to_6 f x1 x2 x3 x4 x5 x6 = f x1`;;
256
257 let rh0 = new_definition `rh0 = &1/(h0 - &1)`;;
258
259 let two6 = new_definition `two6 = constant6 ( &2)`;;
260
261 let zero6 = new_definition `zero6 = constant6 ( &0)`;;
262
263 let dummy6 = new_definition `dummy6 = constant6 ( &0)`;;
264
265 let four6 = new_definition `four6 = constant6 ( &4)`;;
266
267 let mk_126 = new_definition `mk_126 f =
268   compose6 f proj_x1 proj_x2 two6 two6 two6 proj_x6`;;
269
270 let mk_456 = new_definition `mk_456 f =
271   compose6 f two6 two6 two6 proj_x4 proj_x5 proj_x6`;;
272
273 let mk_135 = new_definition `mk_135 f = 
274   compose6 f proj_x1 two6 proj_x3 two6 proj_x5 two6`;;
275
276 let add6 = new_definition `add6 f g x1 x2 x3 x4 x5 x6 = 
277   f x1 x2 x3 x4 x5 x6 + g x1 x2 x3 x4 x5 x6`;;
278
279 let scalar6 = new_definition `scalar6 f r x1 x2 x3 x4 x5 x6 = 
280   (f x1 x2 x3 x4 x5 x6) * (r:real)`;;
281
282 let mul6 = new_definition `mul6 f g x1 x2 x3 x4 x5 x6 = 
283   f x1 x2 x3 x4 x5 x6 * g x1 x2 x3 x4 x5 x6`;;
284
285 let div6 = new_definition `div6 f g x1 x2 x3 x4 x5 x6 = 
286   f x1 x2 x3 x4 x5 x6 / g x1 x2 x3 x4 x5 x6`;;
287
288 let sub6 = new_definition `sub6 f g x1 x2 x3 x4 x5 x6 = 
289   f x1 x2 x3 x4 x5 x6  -  g x1 x2 x3 x4 x5 x6`;;
290
291 let proj_y1 = new_definition `proj_y1 x1 x2 x3 x4 x5 x6 = 
292   sqrt(x1)`;;
293
294 let proj_y2 = new_definition `proj_y2 x1 x2 x3 x4 x5 x6 = 
295   sqrt(x2)`;;
296
297 let proj_y3 = new_definition `proj_y3 x1 x2 x3 x4 x5 x6 = 
298   sqrt(x3)`;;
299
300 let proj_y4 = new_definition `proj_y4 x1 x2 x3 x4 x5 x6 = 
301   sqrt(x4)`;;
302
303 let proj_y5 = new_definition `proj_y5 x1 x2 x3 x4 x5 x6 = 
304   sqrt(x5)`;;
305
306 let proj_y6 = new_definition `proj_y6 x1 x2 x3 x4 x5 x6 = 
307   sqrt(x6)`;;
308
309
310 let domain6 = new_definition `domain6 h f g = 
311   (!x1 x2 x3 x4 x5 x6. h x1 x2 x3 x4 x5 x6 ==>
312       (f x1 x2 x3 x4 x5 x6 = g x1 x2 x3 x4 x5 x6))`;;
313
314 (* should be called something different, because we project out 3 coords *)
315
316 let rotate234 = new_definition `rotate234 f = 
317   compose6 f proj_x2 proj_x3 proj_x4 unit6 unit6 unit6`;;
318
319 let rotate126 = new_definition `rotate126 f = 
320   compose6 f proj_x1 proj_x2 proj_x6 unit6 unit6 unit6`;;
321
322 let rotate345 = new_definition `rotate345 f = 
323   compose6 f proj_x3 proj_x4 proj_x5 unit6 unit6 unit6`;;
324
325 let functional_overload() = (
326   overload_interface ("+",`add6`);
327   overload_interface ("-",`sub6`);
328   overload_interface ("/",`div6`);
329   overload_interface ("*",`mul6`));;
330
331 let _ = functional_overload();;
332
333 let h0cut = new_definition `h0cut y = if (y <= &2 * h0) then &1 else &0`;;
334
335 let cos797 = new_definition `cos797 = cos(#0.797)`;;
336
337 (*  =gamma2= deprecated, Jan 23, 2013. *)
338
339 (* formula corrected Jan 23, 2013.  See GAMMAX_GAMMA2_X *)
340
341 let gamma2_x_div_azim_v2 =  
342   new_definition `gamma2_x_div_azim_v2 m x = (&8 - x)* sqrt x / (&24) -
343   ( &2 * (&2 * mm1/ pi) * (&1 - sqrt x / sqrt8) - 
344       (&8 * mm2 / pi) * m * lfun (sqrt x / &2))`;;
345
346 let gamma2_x1_div_a_v2 = new_definition `gamma2_x1_div_a_v2 m = 
347   promote1_to_6 (gamma2_x_div_azim_v2 m)`;;
348
349 let gamma3f_x_div_sqrtdelta = new_definition `gamma3f_x_div_sqrtdelta m4 m5 m6  = 
350   constant6 (&1/ &12)  - 
351     (scalar6 ( mk_456 (rotate5 (sol_euler_x_div_sqrtdelta)) +
352                 mk_456 (rotate6 (sol_euler_x_div_sqrtdelta)) + 
353                 mk_456 (rotate4 (sol_euler_x_div_sqrtdelta)) 
354                )  (&2 * mm1/ pi)    
355     -       
356       scalar6 (
357         (scalar6 (uni(lfun,(scalar6 proj_y4  #0.5))) m4) * 
358           mk_456 (rotate4 (dih_x_div_sqrtdelta_posbranch)) +
359         (scalar6 (uni(lfun,(scalar6 proj_y5  #0.5))) m5) * 
360           mk_456 (rotate5 (dih_x_div_sqrtdelta_posbranch)) +
361         (scalar6 (uni(lfun,(scalar6 proj_y6  #0.5))) m6) * 
362           mk_456 (rotate6 (dih_x_div_sqrtdelta_posbranch))
363       )  (&8 * mm2 / pi))`;;
364
365 let vol3f_456 = new_definition `vol3f_456 m4 = 
366   scalar6 ( mk_456 (rotate5 sol_x) +
367                 mk_456 (rotate6 sol_x) + 
368                 mk_456 (rotate4 sol_x) 
369                )  (&2 * mm1/ pi)    
370     -       
371       scalar6 (
372         (scalar6 (uni(lfun,(scalar6 proj_y4  #0.5))) m4) * 
373           mk_456 (rotate4 dih_x) +
374         (uni(lfun,(scalar6 proj_y5  #0.5))) * 
375           mk_456 (rotate5 dih_x) +
376         (uni(lfun,(scalar6 proj_y6  #0.5))) * 
377           mk_456 (rotate6 dih_x)
378       )  (&8 * mm2 / pi)`;;
379
380 let gamma3_x = new_definition `gamma3_x m4  = 
381   mk_456 vol_x  -     vol3f_456 m4`;;
382
383 let gamma23_full8_x = new_definition `gamma23_full8_x m1 =
384      (compose6 (gamma3_x m1) 
385        dummy6 dummy6 dummy6 proj_x1 proj_x2 proj_x6) + 
386      (compose6 (gamma3_x m1) 
387        dummy6 dummy6 dummy6 proj_x1 proj_x3 proj_x5) +
388        scalar6 (dih_x -(mk_126 dih_x + mk_135 dih_x)) (#0.008)`;;
389
390 let gamma23_keep135_x = new_definition `gamma23_keep135_x m1 =
391      (compose6 (gamma3_x m1) 
392        dummy6 dummy6 dummy6 proj_x1 proj_x3 proj_x5) +
393        scalar6 (dih_x - mk_135 dih_x) (#0.008)`;;
394
395 (* added 2013-05-20 *)
396
397 (* new 2013-May *)
398
399     let flat_term = new_definition `flat_term y = sol0 * (y - &2 * h0)/(&2 * h0 - &2)`;;
400
401     let mu_y = new_definition `mu_y y1 y2 y3 = 
402      (#0.012 + #0.07 * (#2.52 - y1) + #0.01 * (#2.52 * &2 - y2 - y3 ))`;;
403
404     let mu6_x = new_definition `(mu6_x:real->real->real->real->real->real->real) = 
405       (constant6 (#0.012) + constant6 (#0.07) * (constant6 (#2.52) - proj_y1)
406           + constant6 (#0.01) * ((constant6 (#2.52 * &2) - proj_y2 - proj_y3)))`;;
407
408     let taud = new_definition `taud y1 y2 y3 y4 y5 y6 = 
409        sol0 * (y1 - &2 * h0)/(&2 * h0 - &2) +
410         sqrt(delta_y y1 y2 y3 y4 y5 y6) * 
411          (#0.012 + #0.07 * (#2.52 - y1) + #0.01 * (#2.52 * &2 - y2 - y3 ))`;;
412
413     let taud_x = new_definition `taud_x  = 
414       add6 (promote1_to_6 flat_term_x)
415         (mul6 (uni(sqrt,delta_x))
416         (compose6 (promote3_to_6 mu_y) proj_y1 proj_y2 proj_y3 dummy6 dummy6 dummy6))`;;
417
418     let delta_x1 = new_definition 
419       `delta_x1 x1 x2 x3 x4 x5 x6 = 
420       -- x1 * x4 + x2 * x5 - x3 * x5 - x2 * x6 + x3 * x6 + 
421         x4 * (-- x1 + x2 + x3 - x4 + x5 + x6)`;;
422
423     let delta_x2 = new_definition
424       `delta_x2 x1 x2 x3 x4 x5 x6 = x1 *x4 - x3* x4 - x2* x5 - x1* x6 + x3* x6 + 
425       x5* (x1 - x2 + x3 + x4 - x5 + x6)`;;
426
427     let delta_x3 = new_definition
428       `delta_x3 x1 x2 x3 x4 x5 x6 = x1 *x4 - x2 *x4 - x1* x5 + x2* x5 - 
429       x3* x6 + (x1 + x2 - x3 + x4 + x5 - x6)* x6`;;
430
431     let delta_x4= new_definition(`delta_x4 x1 x2 x3 x4 x5 x6
432                                  =  -- x2* x3 -  x1* x4 + x2* x5
433         + x3* x6 -  x5* x6 + x1* (-- x1 +  x2 +  x3 -  x4 +  x5 +  x6)`);;
434
435     let delta_x5 = new_definition
436       `delta_x5 x1 x2 x3 x4 x5 x6 = -- x1* x3 + x1 *x4 - x2* x5 + x3* x6 - x4* x6 + 
437       x2* (x1 - x2 + x3 + x4 - x5 + x6)`;;
438
439     let delta_x6 = new_definition(`delta_x6 x1 x2 x3 x4 x5 x6
440         = -- x1 * x2 - x3*x6 + x1 * x4
441         + x2* x5 - x4* x5 + x3*(-- x3 + x1 + x2 - x6 + x4 + x5)`);;
442
443     let dnum1 = new_definition `dnum1 e1 e2 e3 x4 x5 x6 = 
444       (&16 - &2 * x4) * e1 + (x5 - &8) * e2 + (x6 - &8) * e3`;;
445
446     let ups_126 = new_definition
447       `ups_126 = compose6 (promote3_to_6 ups_x) proj_x1 proj_x2 proj_x6 dummy6 dummy6 dummy6`;;
448
449     let taud_v4_D2_num = new_definition `taud_v4_D2_num y1 y2 y3 y4 y5 y6 = 
450       (let d = delta_y y1 y2 y3 y4 y5 y6 in
451      let delta' = y_of_x delta_x1 y1 y2 y3 y4 y5 y6 * (&2 * y1) in
452      let delta'' = -- &8 * (y1 * y4) pow 2  + y_of_x delta_x1 y1 y2 y3 y4 y5 y6 * (&2) in
453      let mu = (#0.012 + #0.07 * (#2.52 - y1) + #0.01 * (#2.52 * &2 - y2 - y3 )) in 
454      let mu' = -- #0.07 in
455        (mu' * d * delta' - (&1 / &4) * mu * (delta' pow 2) + (&1 / &2) * mu * d * delta'' ))`;;
456
457     let taud_D1_num_x = new_definition `taud_D1_num_x = 
458       (let d = delta_x in
459        let delta' = delta_x1 * (constant6 (&2)) * proj_y1 in
460        let mu = mu6_x in
461        let mu' = constant6 (-- #0.07) in
462        let ft' = sol0 / (#0.52) in
463          (mu' * d + constant6(&1 / &2) * mu * delta' + constant6 ft' * uni(sqrt ,d)))`;;
464
465     let taud_D2_num_x = new_definition `taud_D2_num_x =
466       (let d = delta_x in
467        let delta' = delta_x1 * (constant6 (&2)) * proj_y1 in
468        let delta'' = constant6 (-- &8) * proj_x1 * proj_x4 + delta_x1 * constant6 (&2) in
469        let mu = mu6_x in
470        let mu' = constant6 (-- #0.07) in
471        (mu' * d * delta' - constant6(&1 / &4) * mu * (delta' *delta') + constant6(&1 / &2) * mu * d * delta'' ))`;;
472
473     let edge2_flatD_x1 = new_definition `edge2_flatD_x1 d x2 x3 x4 x5 x6 = 
474       quadratic_root_plus (abc_of_quadratic (
475         \x1.  d - delta_x x1 x2 x3 x4 x5 x6))`;;
476
477     let edge2_126_x = new_definition `edge2_126_x d x4 x5 = 
478       compose6 (edge2_flatD_x1) (constant6 d) proj_x1 proj_x2 proj_x6 (constant6 x4) (constant6 x5)`;;
479
480     let edge2_135_x = new_definition `edge2_135_x d x4 x6 =
481       compose6 (edge2_flatD_x1) (constant6 d) proj_x1 proj_x3 proj_x5 (constant6 x4) (constant6 x6)`;;
482
483     let edge2_234_x = new_definition `edge2_234_x d x5 x6 =
484       compose6 (edge2_flatD_x1) (constant6 d) proj_x2 proj_x3 proj_x4 (constant6 x5) (constant6 x6)`;;
485
486     let flat_term2_126_x = new_definition `flat_term2_126_x d x4 x5 = 
487       uni(flat_term_x,edge2_126_x d x4 x5)`;;
488
489     let flat_term2_135_x = new_definition `flat_term2_135_x d x4 x6 = 
490       uni(flat_term_x,edge2_135_x d x4 x6)`;;
491
492     let flat_term2_234_x = new_definition `flat_term2_234_x d x5 x6 = 
493       uni(flat_term_x,edge2_234_x d x5 x6)`;;
494
495     let mud_135_x = new_definition `mud_135_x_v1 y2 y4 y6 = 
496       mul6 (compose6 (promote3_to_6 mu_y) (constant6 y2) proj_y1 proj_y3 dummy6 dummy6 dummy6)
497         (uni(sqrt,(delta_135_x (y2*y2) (y4*y4) (y6*y6))))`;;
498
499     let mud_126_x = new_definition `mud_126_x_v1 y3 y4 y5 = 
500       mul6 (compose6 (promote3_to_6 mu_y) (constant6 y3) proj_y1 proj_y2 dummy6 dummy6 dummy6)
501         (uni(sqrt,(delta_126_x (y3*y3) (y4*y4) (y5*y5))))`;;
502
503     let mud_234_x = new_definition `mud_234_x_v1 y1 y5 y6 = 
504       mul6 (compose6 (promote3_to_6 mu_y) (constant6 y1) proj_y2 proj_y3 dummy6 dummy6 dummy6)
505         (uni(sqrt,(delta_234_x (y1*y1) (y5*y5) (y6*y6))))`;;
506
507     let mudLs_234_x = new_definition `mudLs_234_x d1s d2s y1 y5 y6 = 
508       (mul6 (compose6 (mu6_x) (constant6 (y1*y1)) proj_x2 proj_x3 dummy6 dummy6 dummy6)
509         (constant6 (&1/ (d1s+d2s)) * (delta_234_x (y1*y1) (y5*y5) (y6*y6) - constant6(d1s *d1s)) + constant6(d1s)))`;;
510
511     let mudLs_126_x = new_definition `mudLs_126_x d1s d2s y3 y4 y5 = 
512       (mul6 (compose6 (mu6_x) (constant6 (y3*y3)) proj_x1 proj_x2 dummy6 dummy6 dummy6)
513         (constant6 (&1/ (d1s+d2s)) * (delta_126_x (y3*y3) (y4*y4) (y5*y5) - constant6(d1s *d1s)) + constant6(d1s)))`;;
514
515     let mudLs_135_x = new_definition `mudLs_135_x d1s d2s y2 y4 y6 = 
516       (mul6 (compose6 (mu6_x) (constant6 (y2*y2)) proj_x1 proj_x3 dummy6 dummy6 dummy6)
517         (constant6 (&1/ (d1s+d2s)) * (delta_135_x (y2*y2) (y4*y4) (y6*y6) - constant6(d1s *d1s)) + constant6(d1s)))`;;
518
519
520
521
522 end;;