%------------------------------------------------------------------------------
% sin/cos
%
%     Author: David Lester, Manchester University
%
%     Version 1.0            18/2/09   Initial Release Version
%------------------------------------------------------------------------------

sincosx: THEORY

BEGIN

  IMPORTING prelude_aux, atanx, trig_fnd@sincos, trig_fnd@trig_values,
            reals@factorial_props, trig_fnd@trig_approx, remx

  n,p: VAR nat
  x:   VAR real
  cx:  VAR cauchy_real

  real_3pi16: TYPE+ = {x | -3*pi/16 < x & x < 3*pi/16} CONTAINING 0

  cauchy_real_3pi16?(c:[nat->int]):bool 
                                    =EXISTS (x:real_3pi16): cauchy_prop(x,c)

  cauchy_real_3pi16: TYPE+=(cauchy_real_3pi16?) CONTAINING (LAMBDA p:0)

  nnsreal: TYPE+ = {x | 0 <= x & x < 1/2} CONTAINING 0

  cauchy_nnsreal?(c:[nat->int]):bool = EXISTS (x:nnsreal): cauchy_prop(x,c)

  cauchy_nnsreal: TYPE+=(cauchy_nnsreal?) CONTAINING (LAMBDA p:0)

  JUDGEMENT nnsreal           SUBTYPE_OF nnreal
  JUDGEMENT cauchy_nnsreal    SUBTYPE_OF cauchy_real
  JUDGEMENT cauchy_real_3pi16 SUBTYPE_OF cauchy_real

  sx:   VAR real_3pi16
  snx:  VAR nnsreal
  csx:  VAR cauchy_real_3pi16
  ssx:  VAR ssmallreal
  cssx: VAR cauchy_ssmallreal
  csnx: VAR cauchy_nnsreal

  cauchy_sin_series(n):cauchy_real
    = cauchy_div(cauchy_int((-1)^n),cauchy_int(factorial(2*n+1)))

  cauchy_cos_series(n):cauchy_real
    = cauchy_div(cauchy_int((-1)^n),cauchy_int(factorial(2*n)))

  sin_series_lemma: LEMMA cauchy_prop(sin_term(1)(n), cauchy_sin_series(n))

  cos_series_lemma: LEMMA cauchy_prop(cos_term(1)(n), cauchy_cos_series(n))

  cauchy_sin_drx(csnx): cauchy_real
    = (LAMBDA p: round(cauchy_powerseries(csnx,cauchy_sin_series,p+2)(p+2)/4))

  cauchy_cos_drx(csnx): cauchy_real
    = (LAMBDA p: round(cauchy_powerseries(csnx,cauchy_cos_series,p+2)(p+2)/4))

  sin_drx_lemma: LEMMA cauchy_prop(snx,csnx) AND snx /= 0 =>
                cauchy_prop(sin(sqrt(snx))/sqrt(snx),cauchy_sin_drx(csnx))

  cauchy_sin_dr(csx): cauchy_real
    = cauchy_mul(csx,cauchy_sin_drx(cauchy_mul(csx,csx)))

  cauchy_cos_dr(csx): cauchy_real = cauchy_cos_drx(cauchy_mul(csx,csx))

  sin_dr_lemma: LEMMA cauchy_prop(sx,csx) =>
                      cauchy_prop(sin(sx),cauchy_sin_dr(csx))

  cos_dr_lemma: LEMMA cauchy_prop(sx,csx) =>
                      cauchy_prop(cos(sx),cauchy_cos_dr(csx))

  cauchy_sin(cx): cauchy_real
    = LET p2  = cauchy_div2n(cauchy_pi,2),
          s2  = cauchy_sqrt(cauchy_div2n(cauchy_int(1),1)),
          s   = round(cauchy_div(cx,p2)(2)/4),
          n   = remx(8)(s),
          cy  = cauchy_sub(cx, cauchy_mul(p2, cauchy_int(s))),
          sin = cauchy_sin_dr(cy),
          cos = cauchy_cos_dr(cy)
      IN IF    n = 0 THEN sin
         ELSIF n = 1 THEN cauchy_mul(s2, cauchy_add(cos, sin))
         ELSIF n = 2 THEN cos
         ELSIF n = 3 THEN cauchy_mul(s2, cauchy_sub(cos, sin))
         ELSIF n = 4 THEN cauchy_neg(sin)
         ELSIF n = 5 THEN cauchy_neg(cauchy_mul(s2, cauchy_add(cos, sin)))
         ELSIF n = 6 THEN cauchy_neg(cos)
                     ELSE cauchy_neg(cauchy_mul(s2, cauchy_sub(cos, sin)))
         ENDIF

  sin_lemma: LEMMA cauchy_prop(x,cx) => cauchy_prop(sin(x),cauchy_sin(cx))

  cauchy_cos(cx): cauchy_real
    = LET p2  = cauchy_div2n(cauchy_pi,2),
          s2  = cauchy_sqrt(cauchy_div2n(cauchy_int(1),1)),
          s   = round(cauchy_div(cx,p2)(2)/4),
          n   = remx(8)(s),
          cy  = cauchy_sub(cx, cauchy_mul(p2, cauchy_int(s))),
          sin = cauchy_sin_dr(cy),
          cos = cauchy_cos_dr(cy)
      IN IF    n = 0 THEN cos
         ELSIF n = 1 THEN cauchy_mul(s2, cauchy_sub(cos, sin))
         ELSIF n = 2 THEN cauchy_neg(sin)
         ELSIF n = 3 THEN cauchy_neg(cauchy_mul(s2, cauchy_add(cos, sin)))
         ELSIF n = 4 THEN cauchy_neg(cos)
         ELSIF n = 5 THEN cauchy_neg(cauchy_mul(s2, cauchy_sub(cos, sin)))
         ELSIF n = 6 THEN sin
                     ELSE cauchy_mul(s2, cauchy_add(cos, sin))
         ENDIF

  cos_lemma: LEMMA cauchy_prop(x,cx) => cauchy_prop(cos(x),cauchy_cos(cx))

END sincosx