(With some prodding from Jason H) The Newton iteration for e^(ix)=-1 is (c184) %i*%e^-(%i*x)+x+%i - %i x (d184) %i %e + x + %i which indeed converges quadratically (c185) subst(3.b0,x,%) (d185) 1.0007503399554542729b-2 %i + 3.1411200080598672221b0 (c186) subst(%,x,d184) (d186) 3.1415974073206700026b0 - 5.0129702705366286305b-5 %i (c187) subst(%,x,d184) (d187) 3.141592653828090399b0 - 1.2451741388209226812b-9 %i (c188) subst(%,x,d184) (d188) 3.1415926535897932388b0 - 7.4708305947829289906b-19 %i (c189) subst(%,x,d184) (d189) 3.1415926535897932385b0 (interesting that it got imagpart exact. Macsyma doesn't, e.g., force 3.1... + .69b-105*%i to real.) But note that the convergence is *cubic* if we just subst realparts: (c190) subst(realpart(d185),x,d184) (d190) 1.1169689640015603699b-7 %i + 3.1415926535721955587b0 (c191) subst(realpart(%),x,d184) (d191) 3.1415926535897932385b0 In fact, this has become Salamin's x+sin(x) iteration, which is simply Newton's method for cot(x/2). But Jason objected to using e, so we cheat: x 1 (d245) (%i - -) (----------- + 1) + x n %i x n (---- + 1) n (c253) (bfprecision : 2,subst([x = 3.0b0,n = 2^4],d245)) (d253) 3.0b-1 %i + 3.1b0 (c254) (bfprecision : 4,subst([x = realpart(%),n = 2^8],d245)) (d254) 2.455b-2 %i + 3.14b0 (c255) (bfprecision : 8,subst([x = realpart(%),n = 2^16],d245)) (d255) 1.6573659b-4 %i + 3.1415924b0 (c256) (bfprecision : 16,subst([x = realpart(%),n = 2^32],d245)) (d256) 3.141592653589794b0 - 1.244747560102321b-9 %i (c257) (bfprecision : 32,subst([x = realpart(%),n = 2^64],d245)) (d257) 3.1415926535897932384626433832795b0 - 1.8571745299313218459235921249323b-17 %i I carefully wrote a recursive (1+epsilon)^2^n - 1 function. Why is it needless? --rwg