Correctness and arithmetic resource bound for buffered phase approximation

First compute rational sine and cosine approximations at any rational $m\in [0,2]$ to absolute error $ϵ>0$. For each series, retain at least the terms through index $2$ and use the recurrence of successive terms. Subsequent absolute-term ratios are at most $1\mathrm{/}2$ and decrease: for cosine they are $m^2\mathrm{/}((2k+1)(2k+2))$, and for sine they are $m^2\mathrm{/}((2k+2)(2k+3))$. Twice the first omitted absolute term bounds the tail. Increment the index until this rational bound is at most $ϵ$. The geometric decrease gives $O(1+{\log }^{+}(1\mathrm{/}ϵ))$ terms and arithmetic operations, uniformly in $m$. The returned partial sums $c,s$ therefore satisfy $\mathrm{\mid }c-\cos m\mathrm{\mid },\mathrm{\mid }s-\sin m\mathrm{\mid }\le ϵ$.

The integral formula and $\arctan 1=\pi \mathrm{/}4$ in [F4] give $1\mathrm{/}2<\pi \mathrm{/}4<1$: on the two cells of the partition $0,1\mathrm{/}2,1$, the decreasing integrand has lower sum $\frac{1}{2}(4\mathrm{/}5+1\mathrm{/}2)=13\mathrm{/}20>1\mathrm{/}2$ and upper sum $\frac{1}{2}(1+4\mathrm{/}5)=9\mathrm{/}10<1$. Thus $2<\pi <4$. First check the coordinates for exact zero; an axis input takes instruction 2 of the algorithm. Otherwise set $a=\mathrm{\mid }u\mathrm{\mid }$, $b=\mathrm{\mid }v\mathrm{\mid }$ by exact comparisons and signs, and initialize $[L,U]=[0,2]$. It contains $\theta$ by [F1] and $\pi \mathrm{/}2<2$. While its width $w=U-L$ exceeds $h\mathrm{/}8$, set $m=(L+U)\mathrm{/}2$, $ϵ=w\mathrm{/}96$, and compute $c,s$ as in step 1.1. Form the exact real expression $q=bc-as$. If $q>2ϵ$, replace $L$ by $m$. If $q<-2ϵ$, replace $U$ by $m$. Otherwise replace the bracket by $[L,U]\cap [m-12ϵ,m+12ϵ]$. These are two exact comparisons against rational buffers; equality belongs to the third branch. No comparison with an exact transcendental value is made.

Let $D=b\cos m-a\sin m=\sin (\theta -m)$ by [F3]. Since $a,b\in [0,1]$, $\mathrm{\mid }D-q\mathrm{\mid }\le (a+b)ϵ\le 2ϵ$. Also $\mathrm{\mid }\theta -m\mathrm{\mid }\le 2$ because $\theta ,m\in [0,2]$. By [F2], the sign of $D$ is the sign of $\theta -m$, and $\mathrm{\mid }\theta -m\mathrm{\mid }\le 3\mathrm{\mid }D\mathrm{\mid }$. Thus either definite comparison in step 2.1 retains $\theta$ and halves the bracket. In the third branch, $\mathrm{\mid }D\mathrm{\mid }\le 4ϵ$ implies $\mathrm{\mid }\theta -m\mathrm{\mid }\le 12ϵ$, so the intersection also retains $\theta$ and has width at most $24ϵ=w\mathrm{/}4$. It is nonempty. Every iteration contracts by at least a factor $1\mathrm{/}2$, independently of how close the unknown phase is to a midpoint.

Obtain a rational enclosure $I_{\pi }$ for $\pi$ of width at most $h\mathrm{/}8$ by Machin's identity $\pi =16\arctan (1\mathrm{/}5)-4\arctan (1\mathrm{/}239)$. For completeness, put $\alpha =\arctan (1\mathrm{/}5)$, $\beta =\arctan (1\mathrm{/}239)$. The integral in [F4] gives $0<\beta <\alpha <1\mathrm{/}5$, so $0<4\alpha -\beta <4\mathrm{/}5<1<\pi \mathrm{/}2$ by step 2.1. Sine/cosine addition and their quotient give $\tan (2\alpha )=5\mathrm{/}12$, $\tan (4\alpha )=120\mathrm{/}119$, and $\tan (4\alpha -\beta )=(120\cdot 239-119)\mathrm{/}(119\cdot 239+120)=1$. Injectivity of the principal tangent gives $4\alpha -\beta =\pi \mathrm{/}4$. For each arctangent the positive terms $r^{2j+1}\mathrm{/}(2j+1)$ decrease to zero for $r=1\mathrm{/}5,1\mathrm{/}239$, so [F7] shows that consecutive alternating partial sums enclose its value, and their separation is at most $5^{-(2N+1)}$ after the first $N$ terms. Choose $N$ by incrementing until $20\cdot 5^{-(2N+1)}\le h\mathrm{/}8$ and take the outward linear combination of the two enclosures. This uses $O(1+\log (1\mathrm{/}h))$ rational operations. All displayed denominators are positive; all angles in these tangent identities are in its nonsingular principal range.

The loop terminates in $O(1+\log (1\mathrm{/}h))$ iterations, returning rational $L,U$ containing $\theta$ with $U-L\le h\mathrm{/}8$. Its comparison count is linear in the number of iterations. In an iteration with width $w>h\mathrm{/}8$, step 1.1 needs $O(1+\log (1\mathrm{/}h))$ rational operations. The total is therefore $O((1+\log (1\mathrm{/}h){)}^2)$. The exact real expressions compared are only rational linear combinations of the fixed coordinates $a,b$. No algebraic-coordinate approximation, root isolation, normalized-vector division, or sensitivity bound occurs in these instructions.

On an axis, the explicit OP0 branches build $J$ as follows: if $v=0,u>0$, use $[0,0]$; if $v=0,u<0$, use $I_{\pi }$; if $u=0,v>0$, use $I_{\pi }\mathrm{/}2$; if $u=0,v<0$, use $-I_{\pi }\mathrm{/}2$. The zero vector is excluded by the unit promise (the general atan2 procedure returns INVALID there). No angular bisection is needed on an axis. Otherwise use the exact signs of $u,v$ to build a rational interval $J$ for a phase $\varphi$ of $u+iv$: if $u>0,v>0$ use $[L,U]$; if $u<0,v>0$ use $I_{\pi }-[L,U]$; if $u<0,v<0$ use $-I_{\pi }+[L,U]$; if $u\ge 0,v<0$ use $-[L,U]$. Reflection, quarter-turn and parity identities give $u+iv=e^{i\varphi }$, and $\mathrm{width}(J)\le h\mathrm{/}4$. Let $t_0$ be its rational midpoint. Choose the least nonnegative integer $k$ with $2^{-k}\le h\mathrm{/}16$ by repeated halving, and round $t_0$ to a multiple of $2^{-k}$, obtaining $t$ with $\mathrm{\mid }t-t_0\mathrm{\mid }\le h\mathrm{/}16$. Thus $\mathrm{\mid }t-\varphi \mathrm{\mid }\le 3h\mathrm{/}16\le 3\mathrm{/}16$. Put $\mathrm{\Delta }=t-\varphi$. By [F8], the exponential series converges absolutely and $\mathrm{\mid }{e}^{i\mathrm{\Delta }}-1\mathrm{\mid }\le {\sum }_{n\ge 1}\mathrm{\mid }\mathrm{\Delta }{\mathrm{\mid }}^n\mathrm{/}n!\le {\sum }_{n\ge 1}\mathrm{\mid }\mathrm{\Delta }{\mathrm{\mid }}^n=\mathrm{\mid }\mathrm{\Delta }\mathrm{\mid }\mathrm{/}(1-\mathrm{\mid }\mathrm{\Delta }\mathrm{\mid })\le 3h\mathrm{/}13<h$. By [F5], $e^{it}=e^{i\varphi }{e}^{i\mathrm{\Delta }}$ and $\mathrm{\mid }{e}^{i\varphi }\mathrm{\mid }=1$, so $\mathrm{\mid }{e}^{it}-e^{i\varphi }\mathrm{\mid }<h$. The phase lies in a fixed bounded interval, so this dyadic $t$ has $O(1+\log (1\mathrm{/}h))$ bits. The constructed interval $J$ and rounding bound are the certificate.

Put $z=e^{it}-\omega$. By [F5] and [F6], $Z(t)-Z(\phi )$ has diagonal entries $z,\bar{z}$. On a vector with orthonormal coordinates $x_0,x_1$, its squared output norm is $\mathrm{\mid }z{\mathrm{\mid }}^2(\mathrm{\mid }{x}_0{\mathrm{\mid }}^2+\mathrm{\mid }{x}_1{\mathrm{\mid }}^2)$. Thus every unit vector has output norm $\mathrm{\mid }z\mathrm{\mid }$, and the operator norm is exactly $\mathrm{\mid }z\mathrm{\mid }<h$. Together with the operation counts and termination above, this proves all claims.