(* ========================================================================== *)
(* FLYSPECK - BOOK FORMALIZATION                                              *)
(*                                                                            *)
(* Chapter: Nonlinear                                                  *)
(* Author: Thomas C. Hales                                                    *)
(* Date: 2012-11-30                                                           *)
(* ========================================================================== *)

(*
Combine definitions for nonlinear inequalities here.
*)


module Nonlin_def = struct

(* From lemma.hl *)


let NONLIN = new_definition `NONLIN = T`;;


let sqrt_x1 = define `sqrt_x1 x1 x2 x3 x4 x5 x6 = sqrt x1`;;


let sqrt_x2 = define `sqrt_x2 x1 x2 x3 x4 x5 x6 = sqrt x2`;;


let sqrt_x3 = define `sqrt_x3 x1 x2 x3 x4 x5 x6 = sqrt x3`;;


let sqrt_x4 = define `sqrt_x4 x1 x2 x3 x4 x5 x6 = sqrt x4`;;


let sqrt_x5 = define `sqrt_x5 x1 x2 x3 x4 x5 x6 = sqrt x5`;;


let sqrt_x6 = define `sqrt_x6 x1 x2 x3 x4 x5 x6 = sqrt x6`;;


let halfbump_x = new_definition `halfbump_x x = bump (sqrt x / &2)`;;


let halfbump_x1 = new_definition `halfbump_x1 x1 x2 x3 x4 x5 x6 = halfbump_x x1`;;


let halfbump_x4 = new_definition `halfbump_x4 x1 x2 x3 x4 x5 x6 = halfbump_x x4`;;


let unit6 = define `unit6 x1 x2 x3 x4 x5 x6 = &1`;;


let proj_x1 = define `proj_x1 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x1`;;


let proj_x2 = define `proj_x2 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x2`;;


let proj_x3 = define `proj_x3 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x3`;;


let proj_x4 = define `proj_x4 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x4`;;


let proj_x5 = define `proj_x5 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x5`;;


let proj_x6 = define `proj_x6 (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = x6`;;


let promote = define `promote f (x1:A) (x2:B) (x3:C) (x4:D) (x5:E) (x6:F) = f x1`;;


let unit0 = define `unit0 = &1`;;
 


let pow1 = new_definition `pow1 y = y pow 1`;;


let pow2 = new_definition `pow2 y = y pow 2`;;


let pow3 = new_definition `pow3 y = y pow 3`;;


let pow4 = new_definition `pow4 y = y pow 4`;;


let promote_pow2 = new_definition `promote_pow2 x1 (x2:A) (x3:B) (x4:C) (x5:D) (x6:E) = x1 pow 2`;;


let promote_pow3 = new_definition `promote_pow3 x1 (x2:A) (x3:B) (x4:C) (x5:D) (x6:E) = x1 pow 3`;;

(*
let promote_pow3r = INST_TYPE [(`:real`,`:A`);(`:real`,`:B`);(`:real`,`:C`);(`:real`,`:D`);(`:real`,`:E`)] promote_pow3;;
*)


let compose6 = new_definition `compose6 f p1 p2 p3 p4 p5 p6 
  (x1:real) (x2:real) (x3:real) (x4:real) (x5:real) (x6:real) =
  (f:real->real->real->real->real->real->real)
  (p1 x1 x2 x3 x4 x5 x6)
    (p2 x1 x2 x3 x4 x5 x6)
    (p3 x1 x2 x3 x4 x5 x6)
    (p4 x1 x2 x3 x4 x5 x6)
    (p5 x1 x2 x3 x4 x5 x6)
    (p6 x1 x2 x3 x4 x5 x6)`;;


let scale6 = new_definition `scale6 f
   (r:real)   (x1:real) (x2:real) (x3:real) (x4:real) (x5:real) (x6:real) = 
   (f x1 x2 x3 x4 x5 x6) * r`;;




let quadratic_root_plus_curry = 
  new_definition `quadratic_root_plus_curry a b c = quadratic_root_plus (a,b,c)`;;

(* ----  *)


let sq_pow2 = 
prove_by_refinement(
  `!a x. a pow 2 <= x ==> (sqrt x * sqrt x = x)`,

  (* {{{ proof *)
  [
  REWRITE_TAC[GSYM REAL_POW_2;SQRT_POW2];
  MESON_TAC[REAL_LE_TRANS;Collect_geom.REAL_LE_SQUARE_POW];
  ]);;
  (* }}} *)


let sqrt2_sqrt2 = 
prove_by_refinement(
  `sqrt2 * sqrt2 = &2`,

  (* {{{ proof *)
  [
  REWRITE_TAC[Sphere.sqrt2];
    MATCH_MP_TAC sq_pow2;
  EXISTS_TAC`&0`;
  REAL_ARITH_TAC;
  ]);;
  (* }}} *)

(* This is valid when a > 0 *)


let sol_euler_x_div_sqrtdelta =  new_definition 
 `sol_euler_x_div_sqrtdelta x1 x2 x3 x4 x5 x6 = 
  (let a = sqrt(x1*x2*x3) + sqrt( x1)*(x2 + x3 - x4)/ &2 + 
     sqrt(x2)*(x1 + x3 - x5)/ &2 + sqrt(x3)*(x1 + x2 - x6)/ &2 in
     (matan ((delta_x x1 x2 x3 x4 x5 x6)/(&4 * a pow 2 )))/( a))`;;


let sol_euler246_x_div_sqrtdelta =  new_definition 
 `sol_euler246_x_div_sqrtdelta = rotate4 sol_euler_x_div_sqrtdelta`;;


let sol_euler345_x_div_sqrtdelta =  new_definition 
 `sol_euler345_x_div_sqrtdelta = rotate5 sol_euler_x_div_sqrtdelta`;;


let sol_euler156_x_div_sqrtdelta =  new_definition 
 `sol_euler156_x_div_sqrtdelta = rotate6 sol_euler_x_div_sqrtdelta`;;

(* This is a valid formula for dih_x / sqrt(delta) when delta >0 and d_x4 !=0 *)

let dih_x_div_sqrtdelta = `dih_x_div_sqrtdelta x1 x2 x3 x4 x5 x6 =
  ( let d_x4 = delta_x4 x1 x2 x3 x4 x5 x6 in
          let d = delta_x x1 x2 x3 x4 x5 x6 in
	  let pi_shift = if (d_x4 < &0) then pi else (&0) in
	    pi_shift + (sqrt(&4 * x1) / d_x4) * matan((&4 * x1 * d)/(d_x4 pow 2)))`;;


let dih_x_div_sqrtdelta_posbranch = new_definition
 `dih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
  ( let d_x4 = delta_x4 x1 x2 x3 x4 x5 x6 in
          let d = delta_x x1 x2 x3 x4 x5 x6 in
	     (sqrt(&4 * x1) / d_x4) * matan((&4 * x1 * d)/(d_x4 pow 2)))`;;


let ldih_x_div_sqrtdelta_posbranch = new_definition
 `ldih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
   lfun(sqrt(x1) / &2) * dih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;


let ldih2_x_div_sqrtdelta_posbranch = new_definition
 `ldih2_x_div_sqrtdelta_posbranch =  rotate2 ldih_x_div_sqrtdelta_posbranch`;;


let ldih3_x_div_sqrtdelta_posbranch = new_definition
 `ldih3_x_div_sqrtdelta_posbranch  = rotate3 ldih_x_div_sqrtdelta_posbranch`;;


let ldih5_x_div_sqrtdelta_posbranch = new_definition
 `ldih5_x_div_sqrtdelta_posbranch  = rotate5 ldih_x_div_sqrtdelta_posbranch`;;


let ldih6_x_div_sqrtdelta_posbranch = new_definition
 `ldih6_x_div_sqrtdelta_posbranch  = rotate6 ldih_x_div_sqrtdelta_posbranch`;;


let dih3_x_div_sqrtdelta_posbranch = new_definition
 `dih3_x_div_sqrtdelta_posbranch  = rotate3 dih_x_div_sqrtdelta_posbranch`;;


let dih4_x_div_sqrtdelta_posbranch = new_definition
 `dih4_x_div_sqrtdelta_posbranch  = rotate4 dih_x_div_sqrtdelta_posbranch`;;


let dih5_x_div_sqrtdelta_posbranch = new_definition
 `dih5_x_div_sqrtdelta_posbranch  = rotate5 dih_x_div_sqrtdelta_posbranch`;;


let rhazim_x_div_sqrtdelta_posbranch = new_definition
 `rhazim_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
   rho(sqrt(x1)) * dih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;


let rhazim2_x_div_sqrtdelta_posbranch = new_definition
 `rhazim2_x_div_sqrtdelta_posbranch =  
  rotate2 rhazim_x_div_sqrtdelta_posbranch`;;


let rhazim3_x_div_sqrtdelta_posbranch = new_definition
 `rhazim3_x_div_sqrtdelta_posbranch  = 
  rotate3 rhazim_x_div_sqrtdelta_posbranch`;;


let tau_residual_x  = new_definition
 `tau_residual_x x1 x2 x3 x4 x5 x6 = 
     rhazim_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6
+     rhazim2_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6
+     rhazim3_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;



let lmdih_x_div_sqrtdelta_posbranch = new_definition
 `lmdih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
   lmfun(sqrt(x1) / &2) * dih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;


let lmdih2_x_div_sqrtdelta_posbranch = new_definition 
  `lmdih2_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
         rotate2 lmdih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;


let lmdih3_x_div_sqrtdelta_posbranch = new_definition
 `lmdih3_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
  rotate3 lmdih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 `;;


let lmdih5_x_div_sqrtdelta_posbranch = new_definition
 `lmdih5_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
  rotate5 lmdih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 `;;


let lmdih6_x_div_sqrtdelta_posbranch = new_definition 
  `lmdih6_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6 =
         rotate6 lmdih_x_div_sqrtdelta_posbranch x1 x2 x3 x4 x5 x6`;;


let vol3f_sqrt2_lmplus = new_definition 
  `vol3f_sqrt2_lmplus y1 y2 (y3:real) (y4:real) (y5:real) y6 = 
    (&2 * mm1 / pi) *
             (&2 * dih_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
              &2 * dih2_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
              &2 * dih6_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
              dih3_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
              dih4_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
              dih5_y y1 y2 sqrt2 sqrt2 sqrt2 y6 - &3 * pi) -
             (&8 * mm2 / pi) *
             (
              lfun (y2 / &2) * dih2_y y1 y2 sqrt2 sqrt2 sqrt2 y6 +
              lfun (y6 / &2) * dih6_y y1 y2 sqrt2 sqrt2 sqrt2 y6)`;;


let vol3f_x_sqrt2_lmplus = new_definition
  `vol3f_x_sqrt2_lmplus x1 x2 x3 x4 x5 x6 =
   vol3f_sqrt2_lmplus (sqrt x1) (sqrt x2) (sqrt x3) (sqrt x4) (sqrt x5)      (sqrt x6)`;;


let vol3f_x_lfun = new_definition 
 `vol3f_x_lfun x1 x2 (x3:real) (x4:real) (x5:real) x6 = vol3f (sqrt x1) (sqrt x2) (sqrt x6)  sqrt2 lfun `;;


let vol3_x_sqrt = new_definition
 `vol3_x_sqrt x1 x2 (x3:real) (x4:real) (x5:real) x6  = vol_y sqrt2 sqrt2 sqrt2 (sqrt x1) (sqrt x2) (sqrt x6) `;;



let gamma23f = new_definition `gamma23f y1 y2 y3 y4 y5 y6 w1 w2 r f =
      (gamma3f y1 y2 y6 r f / &w1 + gamma3f y1 y3 y5 r f / &w2
      + (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)) `;;


let monomial = new_definition `monomial n1 n2 n3 n4 n5 n6 y1 y2 y3 y4 y5 y6 = 
   (y1 pow n1) * (y2 pow n2) * (y3 pow n3) * (y4 pow n4) * (y5 pow n5) * (y6 pow n6)`;;


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)`;;


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)`;;

(* Functional Equation *)


let uni = new_definition `uni (f,x) x1 x2 x3 x4 x5 x6 = 
  (f:A->B) ( x x1 x2 x3 x4 x5 x6)`;;


let constant6 = new_definition `constant6 c x1 x2 x3 x4 x5 x6 = c`;;


let promote3_to_6 = new_definition
    `promote3_to_6 f x1 x2 x3 x4 x5 x6 = f x1 x2 x3`;;


let promote1_to_6 = new_definition
  `promote1_to_6 f x1 x2 x3 x4 x5 x6 = f x1`;;


let rh0 = new_definition `rh0 = &1/(h0 - &1)`;;


let two6 = new_definition `two6 = constant6 ( &2)`;;


let zero6 = new_definition `zero6 = constant6 ( &0)`;;


let dummy6 = new_definition `dummy6 = constant6 ( &0)`;;


let four6 = new_definition `four6 = constant6 ( &4)`;;


let mk_126 = new_definition `mk_126 f =
  compose6 f proj_x1 proj_x2 two6 two6 two6 proj_x6`;;


let mk_456 = new_definition `mk_456 f =
  compose6 f two6 two6 two6 proj_x4 proj_x5 proj_x6`;;


let mk_135 = new_definition `mk_135 f = 
  compose6 f proj_x1 two6 proj_x3 two6 proj_x5 two6`;;


let add6 = new_definition `add6 f g x1 x2 x3 x4 x5 x6 = 
  f x1 x2 x3 x4 x5 x6 + g x1 x2 x3 x4 x5 x6`;;


let scalar6 = new_definition `scalar6 f r x1 x2 x3 x4 x5 x6 = 
  (f x1 x2 x3 x4 x5 x6) * (r:real)`;;


let mul6 = new_definition `mul6 f g x1 x2 x3 x4 x5 x6 = 
  f x1 x2 x3 x4 x5 x6 * g x1 x2 x3 x4 x5 x6`;;


let div6 = new_definition `div6 f g x1 x2 x3 x4 x5 x6 = 
  f x1 x2 x3 x4 x5 x6 / g x1 x2 x3 x4 x5 x6`;;


let sub6 = new_definition `sub6 f g x1 x2 x3 x4 x5 x6 = 
  f x1 x2 x3 x4 x5 x6  -  g x1 x2 x3 x4 x5 x6`;;


let proj_y1 = new_definition `proj_y1 x1 x2 x3 x4 x5 x6 = 
  sqrt(x1)`;;


let proj_y2 = new_definition `proj_y2 x1 x2 x3 x4 x5 x6 = 
  sqrt(x2)`;;


let proj_y3 = new_definition `proj_y3 x1 x2 x3 x4 x5 x6 = 
  sqrt(x3)`;;


let proj_y4 = new_definition `proj_y4 x1 x2 x3 x4 x5 x6 = 
  sqrt(x4)`;;


let proj_y5 = new_definition `proj_y5 x1 x2 x3 x4 x5 x6 = 
  sqrt(x5)`;;


let proj_y6 = new_definition `proj_y6 x1 x2 x3 x4 x5 x6 = 
  sqrt(x6)`;;



let domain6 = new_definition `domain6 h f g = 
  (!x1 x2 x3 x4 x5 x6. h x1 x2 x3 x4 x5 x6 ==>
      (f x1 x2 x3 x4 x5 x6 = g x1 x2 x3 x4 x5 x6))`;;

(* should be called something different, because we project out 3 coords *)


let rotate234 = new_definition `rotate234 f = 
  compose6 f proj_x2 proj_x3 proj_x4 unit6 unit6 unit6`;;


let rotate126 = new_definition `rotate126 f = 
  compose6 f proj_x1 proj_x2 proj_x6 unit6 unit6 unit6`;;


let rotate345 = new_definition `rotate345 f = 
  compose6 f proj_x3 proj_x4 proj_x5 unit6 unit6 unit6`;;

let functional_overload() = (
  overload_interface ("+",`add6`);
  overload_interface ("-",`sub6`);
  overload_interface ("/",`div6`);
  overload_interface ("*",`mul6`));;

let _ = functional_overload();;


let h0cut = new_definition `h0cut y = if (y <= &2 * h0) then &1 else &0`;;


let cos797 = new_definition `cos797 = cos(#0.797)`;;

(*  =gamma2= deprecated, Jan 23, 2013. *)

(* formula corrected Jan 23, 2013.  See GAMMAX_GAMMA2_X *)


let gamma2_x_div_azim_v2 =  
  new_definition `gamma2_x_div_azim_v2 m x = (&8 - x)* sqrt x / (&24) -
  ( &2 * (&2 * mm1/ pi) * (&1 - sqrt x / sqrt8) - 
      (&8 * mm2 / pi) * m * lfun (sqrt x / &2))`;;


let gamma2_x1_div_a_v2 = new_definition `gamma2_x1_div_a_v2 m = 
  promote1_to_6 (gamma2_x_div_azim_v2 m)`;;


let gamma3f_x_div_sqrtdelta = new_definition `gamma3f_x_div_sqrtdelta m4 m5 m6  = 
  constant6 (&1/ &12)  - 
    (scalar6 ( mk_456 (rotate5 (sol_euler_x_div_sqrtdelta)) +
		mk_456 (rotate6 (sol_euler_x_div_sqrtdelta)) + 
		mk_456 (rotate4 (sol_euler_x_div_sqrtdelta)) 
	       )  (&2 * mm1/ pi)    
    -       
      scalar6 (
	(scalar6 (uni(lfun,(scalar6 proj_y4  #0.5))) m4) * 
	  mk_456 (rotate4 (dih_x_div_sqrtdelta_posbranch)) +
	(scalar6 (uni(lfun,(scalar6 proj_y5  #0.5))) m5) * 
	  mk_456 (rotate5 (dih_x_div_sqrtdelta_posbranch)) +
	(scalar6 (uni(lfun,(scalar6 proj_y6  #0.5))) m6) * 
	  mk_456 (rotate6 (dih_x_div_sqrtdelta_posbranch))
      )  (&8 * mm2 / pi))`;;


let vol3f_456 = new_definition `vol3f_456 m4 = 
  scalar6 ( mk_456 (rotate5 sol_x) +
		mk_456 (rotate6 sol_x) + 
		mk_456 (rotate4 sol_x) 
	       )  (&2 * mm1/ pi)    
    -       
      scalar6 (
	(scalar6 (uni(lfun,(scalar6 proj_y4  #0.5))) m4) * 
	  mk_456 (rotate4 dih_x) +
	(uni(lfun,(scalar6 proj_y5  #0.5))) * 
	  mk_456 (rotate5 dih_x) +
	(uni(lfun,(scalar6 proj_y6  #0.5))) * 
	  mk_456 (rotate6 dih_x)
      )  (&8 * mm2 / pi)`;;


let gamma3_x = new_definition `gamma3_x m4  = 
  mk_456 vol_x  -     vol3f_456 m4`;;


let gamma23_full8_x = new_definition `gamma23_full8_x m1 =
     (compose6 (gamma3_x m1) 
       dummy6 dummy6 dummy6 proj_x1 proj_x2 proj_x6) + 
     (compose6 (gamma3_x m1) 
       dummy6 dummy6 dummy6 proj_x1 proj_x3 proj_x5) +
       scalar6 (dih_x -(mk_126 dih_x + mk_135 dih_x)) (#0.008)`;;


let gamma23_keep135_x = new_definition `gamma23_keep135_x m1 =
     (compose6 (gamma3_x m1) 
       dummy6 dummy6 dummy6 proj_x1 proj_x3 proj_x5) +
       scalar6 (dih_x - mk_135 dih_x) (#0.008)`;;

(* added 2013-05-20 *)

(* new 2013-May *)

    
let flat_term = new_definition `flat_term y = sol0 * (y - &2 * h0)/(&2 * h0 - &2)`;;

    
let mu_y = new_definition `mu_y y1 y2 y3 = 
     (#0.012 + #0.07 * (#2.52 - y1) + #0.01 * (#2.52 * &2 - y2 - y3 ))`;;

    
let mu6_x = new_definition `(mu6_x:real->real->real->real->real->real->real) = 
      (constant6 (#0.012) + constant6 (#0.07) * (constant6 (#2.52) - proj_y1)
	  + constant6 (#0.01) * ((constant6 (#2.52 * &2) - proj_y2 - proj_y3)))`;;

    
let taud = new_definition `taud y1 y2 y3 y4 y5 y6 = 
       sol0 * (y1 - &2 * h0)/(&2 * h0 - &2) +
	sqrt(delta_y y1 y2 y3 y4 y5 y6) * 
	 (#0.012 + #0.07 * (#2.52 - y1) + #0.01 * (#2.52 * &2 - y2 - y3 ))`;;

    
let taud_x = new_definition `taud_x  = 
      add6 (promote1_to_6 flat_term_x)
	(mul6 (uni(sqrt,delta_x))
	(compose6 (promote3_to_6 mu_y) proj_y1 proj_y2 proj_y3 dummy6 dummy6 dummy6))`;;

    
let delta_x1 = new_definition 
      `delta_x1 x1 x2 x3 x4 x5 x6 = 
      -- x1 * x4 + x2 * x5 - x3 * x5 - x2 * x6 + x3 * x6 + 
	x4 * (-- x1 + x2 + x3 - x4 + x5 + x6)`;;

    
let delta_x2 = new_definition
      `delta_x2 x1 x2 x3 x4 x5 x6 = x1 *x4 - x3* x4 - x2* x5 - x1* x6 + x3* x6 + 
      x5* (x1 - x2 + x3 + x4 - x5 + x6)`;;

    
let delta_x3 = new_definition
      `delta_x3 x1 x2 x3 x4 x5 x6 = x1 *x4 - x2 *x4 - x1* x5 + x2* x5 - 
      x3* x6 + (x1 + x2 - x3 + x4 + x5 - x6)* x6`;;

    
let delta_x4= new_definition(`delta_x4 x1 x2 x3 x4 x5 x6
				 =  -- x2* x3 -  x1* x4 + x2* x5
        + x3* x6 -  x5* x6 + x1* (-- x1 +  x2 +  x3 -  x4 +  x5 +  x6)`);;

    
let delta_x5 = new_definition
      `delta_x5 x1 x2 x3 x4 x5 x6 = -- x1* x3 + x1 *x4 - x2* x5 + x3* x6 - x4* x6 + 
      x2* (x1 - x2 + x3 + x4 - x5 + x6)`;;

    
let delta_x6 = new_definition(`delta_x6 x1 x2 x3 x4 x5 x6
        = -- x1 * x2 - x3*x6 + x1 * x4
        + x2* x5 - x4* x5 + x3*(-- x3 + x1 + x2 - x6 + x4 + x5)`);;

    
let dnum1 = new_definition `dnum1 e1 e2 e3 x4 x5 x6 = 
      (&16 - &2 * x4) * e1 + (x5 - &8) * e2 + (x6 - &8) * e3`;;

    
let ups_126 = new_definition
      `ups_126 = compose6 (promote3_to_6 ups_x) proj_x1 proj_x2 proj_x6 dummy6 dummy6 dummy6`;;

    
let taud_v4_D2_num = new_definition `taud_v4_D2_num y1 y2 y3 y4 y5 y6 = 
      (let d = delta_y y1 y2 y3 y4 y5 y6 in
     let delta' = y_of_x delta_x1 y1 y2 y3 y4 y5 y6 * (&2 * y1) in
     let delta'' = -- &8 * (y1 * y4) pow 2  + y_of_x delta_x1 y1 y2 y3 y4 y5 y6 * (&2) in
     let mu = (#0.012 + #0.07 * (#2.52 - y1) + #0.01 * (#2.52 * &2 - y2 - y3 )) in 
     let mu' = -- #0.07 in
       (mu' * d * delta' - (&1 / &4) * mu * (delta' pow 2) + (&1 / &2) * mu * d * delta'' ))`;;

    
let taud_D1_num_x = new_definition `taud_D1_num_x = 
      (let d = delta_x in
       let delta' = delta_x1 * (constant6 (&2)) * proj_y1 in
       let mu = mu6_x in
       let mu' = constant6 (-- #0.07) in
       let ft' = sol0 / (#0.52) in
	 (mu' * d + constant6(&1 / &2) * mu * delta' + constant6 ft' * uni(sqrt ,d)))`;;

    
let taud_D2_num_x = new_definition `taud_D2_num_x =
      (let d = delta_x in
       let delta' = delta_x1 * (constant6 (&2)) * proj_y1 in
       let delta'' = constant6 (-- &8) * proj_x1 * proj_x4 + delta_x1 * constant6 (&2) in
       let mu = mu6_x in
       let mu' = constant6 (-- #0.07) in
       (mu' * d * delta' - constant6(&1 / &4) * mu * (delta' *delta') + constant6(&1 / &2) * mu * d * delta'' ))`;;

    
let edge2_flatD_x1 = new_definition `edge2_flatD_x1 d x2 x3 x4 x5 x6 = 
      quadratic_root_plus (abc_of_quadratic (
	\x1.  d - delta_x x1 x2 x3 x4 x5 x6))`;;

    
let edge2_126_x = new_definition `edge2_126_x d x4 x5 = 
      compose6 (edge2_flatD_x1) (constant6 d) proj_x1 proj_x2 proj_x6 (constant6 x4) (constant6 x5)`;;

    
let edge2_135_x = new_definition `edge2_135_x d x4 x6 =
      compose6 (edge2_flatD_x1) (constant6 d) proj_x1 proj_x3 proj_x5 (constant6 x4) (constant6 x6)`;;

    
let edge2_234_x = new_definition `edge2_234_x d x5 x6 =
      compose6 (edge2_flatD_x1) (constant6 d) proj_x2 proj_x3 proj_x4 (constant6 x5) (constant6 x6)`;;

    
let flat_term2_126_x = new_definition `flat_term2_126_x d x4 x5 = 
      uni(flat_term_x,edge2_126_x d x4 x5)`;;

    
let flat_term2_135_x = new_definition `flat_term2_135_x d x4 x6 = 
      uni(flat_term_x,edge2_135_x d x4 x6)`;;

    
let flat_term2_234_x = new_definition `flat_term2_234_x d x5 x6 = 
      uni(flat_term_x,edge2_234_x d x5 x6)`;;

    
let mud_135_x = new_definition `mud_135_x_v1 y2 y4 y6 = 
      mul6 (compose6 (promote3_to_6 mu_y) (constant6 y2) proj_y1 proj_y3 dummy6 dummy6 dummy6)
	(uni(sqrt,(delta_135_x (y2*y2) (y4*y4) (y6*y6))))`;;

    
let mud_126_x = new_definition `mud_126_x_v1 y3 y4 y5 = 
      mul6 (compose6 (promote3_to_6 mu_y) (constant6 y3) proj_y1 proj_y2 dummy6 dummy6 dummy6)
	(uni(sqrt,(delta_126_x (y3*y3) (y4*y4) (y5*y5))))`;;

    
let mud_234_x = new_definition `mud_234_x_v1 y1 y5 y6 = 
      mul6 (compose6 (promote3_to_6 mu_y) (constant6 y1) proj_y2 proj_y3 dummy6 dummy6 dummy6)
	(uni(sqrt,(delta_234_x (y1*y1) (y5*y5) (y6*y6))))`;;

    
let mudLs_234_x = new_definition `mudLs_234_x d1s d2s y1 y5 y6 = 
      (mul6 (compose6 (mu6_x) (constant6 (y1*y1)) proj_x2 proj_x3 dummy6 dummy6 dummy6)
	(constant6 (&1/ (d1s+d2s)) * (delta_234_x (y1*y1) (y5*y5) (y6*y6) - constant6(d1s *d1s)) + constant6(d1s)))`;;

    
let mudLs_126_x = new_definition `mudLs_126_x d1s d2s y3 y4 y5 = 
      (mul6 (compose6 (mu6_x) (constant6 (y3*y3)) proj_x1 proj_x2 dummy6 dummy6 dummy6)
	(constant6 (&1/ (d1s+d2s)) * (delta_126_x (y3*y3) (y4*y4) (y5*y5) - constant6(d1s *d1s)) + constant6(d1s)))`;;

    
let mudLs_135_x = new_definition `mudLs_135_x d1s d2s y2 y4 y6 = 
      (mul6 (compose6 (mu6_x) (constant6 (y2*y2)) proj_x1 proj_x3 dummy6 dummy6 dummy6)
	(constant6 (&1/ (d1s+d2s)) * (delta_135_x (y2*y2) (y4*y4) (y6*y6) - constant6(d1s *d1s)) + constant6(d1s)))`;;




end;;