<br />
<b>Deprecated</b>:  The each() function is deprecated. This message will be suppressed on further calls in <b>/home/zhenxiangba/zhenxiangba.com/public_html/phproxy-improved-master/index.php</b> on line <b>456</b><br />
Require Coq.ZArith.ZArith.
Require Coq.ZArith.Zdiv.
Require Coq.ZArith.Zpower.
Require Coq.ZArith.Zlogarithm.

Definition base : Z := `10`.

Definition digitcountpos [p:positive] : Z := `((log_inf p) + 2)/3 + 1`.
Definition digitcount [i:Z] : Z :=
  Cases i of
    ZERO => `0`
  | (POS p) => (digitcountpos p)
  | (NEG p) => (digitcountpos p)
  end.

Definition R : Set := (Z->Z)*Z.

Definition Rdigits : R->Z->Z := Fst.
Definition Rsize : R->Z := Snd.
Coercion Rdigits : R>->FUNCLASS.

Definition ok [x:R] : Prop :=
  (n:Z)`(Zabs ((x (n + 1)) - base*(x n))) < base` /\
  `(x (-(Rsize x))) = 0`.

Definition IZR [i:Z] : R :=
  ([n:Z]
     Cases n of
       ZERO => i
     | (POS p) => `i*(Zpower_pos base p)`
     | (NEG p) => `i/(Zpower_pos base p)`
     end,
   (digitcount i)).
Coercion IZR : Z>->R.

Definition Zmax [m,n:Z] : Z := `((Zabs (m + n)) + (Zabs (m - n)))/2`.
Definition Zdiv_rounded [m,n:Z] : Z := `(m + n/2)/n`.

Definition Rplus [x,y:R] : R :=
  ([n:Z]`(Zdiv_rounded ((x (n + 1)) + (y (n + 1))) base)`,
    `(Zmax (Rsize x) (Rsize y)) + 1`).
Definition Ropp [x:R] : R := ([n:Z]`-(x n)`,(Rsize x)).
Definition Rmult [x,y:R] : R :=
  ([n:Z]`(Zdiv_rounded ((x (n + (Rsize y) + 1))*(y (n + (Rsize x) + 1)))
         (Zpower base (n + (Rsize x) + (Rsize y) + 2)))`,
   `(Rsize x) + (Rsize y)`).
Definition Rinv [x:R; m:Z] : R :=
  ([n:Z]`(Zdiv_rounded (Zpower base (n + (Rsize x) + 1)) (x (Rsize x) + 1))`,
    `-m`).
Definition Rlt [x,y:R] : Prop :=
  (EX n:Z | `(y n) - (x n) > 2`).

Definition Rminus [x,y:R] : R := (Rplus x (Ropp y)).
Definition Rdiv [x,y:R; m:Z] : R := (Rmult x (Rinv y m)).
Definition Rne [x,y:R] : Prop := (Rlt x y) \/ (Rlt y x).
Definition Req [x,y:R] : Prop := (ok x) /\ (ok y) /\ ~(Rne x y).

Fixpoint Rpowernat [x:R; n:nat] : R :=
  Cases n of
    O => `1`
  | (S m) => (Rmult (Rpowernat x m) x)
  end.
Definition Rpower [x:R; m:Z; n:Z] : R :=
  Cases n of
    ZERO => `1`
  | (POS p) => (Rpowernat x (convert p))
  | (NEG p) => (Rinv (Rpowernat x (convert p)) `n*m`)
  end.

Definition lim [a:Z->R; m:Z->Z] : R :=
  ([n:Z]`(Zdiv_rounded (a (m (n + 1)) (n + 1)) base)`,`(Rsize (a 0)) + 1`).

Fixpoint partsumnat [a:Z->R; n:nat] : R :=
  Cases n of
    O => `0`
  | (S m) => `(Rplus (partsumnat a m) (a (inject_nat m)))`
  end.
Definition partsum [a:Z->R; n:Z] : R :=
  Cases n of
    ZERO => `0`
  | (POS p) => (partsumnat a (convert p))
  | (NEG p) => (Ropp (partsumnat ([n:Z]`(a (1 - n))`) (convert p)))
  end.
Definition sum [a:Z->R; m:Z->Z] := (lim ([n:Z](partsum a n)) m).

Definition arctanterm [x:R; k:Z] : R :=
  let n = `2*k + 1` in `(Rmult (Zpower (-1) k) (Rdiv (Rpower x 0 n) n 0))`.
Definition arctan [x:R; m:Z] : R :=
  (sum (arctanterm x) [n:Z]`n/m`).

Definition pi : R :=
  `(Rminus (Rmult 16 (arctan (Rdiv 1 5 0) 1))
           (Rmult 4 (arctan (Rdiv 1 239 2) 4)))`.

Eval Compute in `(pi 2)`.

(*
Extraction "R.ml" pi.
*)