Arndt - Algorithms for Programmers, страница 55

ARITHMETICAL ALGORITHMS323andTn+m + Tn−mTn+m − Tn−mUn+m−1 + Un−m−1Un+m−1 − Un−m−1====2 Tn Tm where n ≥ m2 (x2 + 1) Un−1 Um−1 where n ≥ m2 Un−1 Tm where n > m2 Tn Um−1 where n > m(13.343)(13.344)(13.345)(13.346)The key to the fast computation of both Tn and Un via powering algorithms is given by·¸·¸·¸kTk Tk+1T0 T1 0 −1=Uk Uk+1U0 U1 1 2 x(13.347)It is instructive to do the case k = 1 by hand: The first column of the second matrix moves the secondentries of the first matrix to the left, the second column does the recursion.

The generalization for thegeneral case is straightforward. For example, a recursion an = λ1 an−1 + λ2 an−2 + λ3 an−3 leads tokAk Ak+1 Ak+20 0 λ3A0 A1 A2B0 B1 B2  1 0 λ2  = Bk Bk+1 Bk+2 (13.348)Ck Ck+1 Ck+2C0 C1 C20 1 λ1See also [27] and [63].13.16Continued fractions *Setb1x = a0 +b3a2 +a3 +For k > 0 letpkqk(13.349)b2a1 +b4a4 + · · ·be the value of the above fraction if bk+1 is set to zero (setp−1q−1:=10andp0q0:=a01 ).Thenpkqk==ak pk−1 + bk pk−2ak qk−1 + bk qk−2(13.350a)(13.350b)(Simple continued fractions are those with bk = 1 ∀k).Pseudo code for a procedure that computes the pk , qk k = −1 .

. . n of a continued fraction :procedure ratios_from_contfrac(a[0..n], b[0..n], n, p[-1..n], q[-1..n]){p[-1] := 1q[-1] := 0p[0] := a[0]q[0] := 1for k:=1 to n{p[k] := a[k] * p[k-1] + b[k] * p[k-2]q[k] := a[k] * q[k-1] + b[k] * q[k-2]}}Pseudo code for a procedure that fills the first n terms of the simple continued fraction of (the floatingpoint number) x into the array cf[]:CHAPTER 13. ARITHMETICAL ALGORITHMS324procedure continued_fraction(x, n, cf[0..n-1]){for k:=0 to n-1{xi := floor(x)cf[k] := xix := 1 / (x-xi)}}Pseudo code for a function that computes the numerical value of a number x from (the leading n termsof) its simple continued fraction representation:function number_from_contfrac(cf[0..n-1], n){x := cf[n-1]for k:=n-2 to 0 step -1{x := 1/x + cf[k]}return x}(cf.

[24], [23], [13], [17]).It is possible to rewrite a continued fraction with positive ak , bk as an alternating sumx= a0 +∞X(−1)k+1 sk(13.351)k=1Qkb1b1 b2b1 b2 b3k+1i=1 bi= a0 +−+± . . . + (−1)± ...q0 q1q1 q2q2 q3qk qk+1where q−1 = 0, q1 = 1 and qn = an qn−1 + bn qn−2 , cf.

[49].Continued fractions with bk = 1 ∀k1x = a0 +(13.352)1a1 +1a2 +a3 +1a4 + · · ·are called simple continued fractionsSee [24] for an easy introduction and also [39].13.17Some hypergeometric identities *13.17.1DefinitionLet xk̄ := x (x + 1) (x + 2) . . . (x + k) thenµ¶a, b ¯¯Fz:=¯c∞Xak̄ bk̄ z kck̄ k!(13.353)k=0Often the so called Pochhammer symbol (x)k is used instead of the ”rising factorial power” notationxk̄ used here.

The k! in the denominators is there for historical reasons. You might want to have anadditional 1 in the lower entries in mind:µ¶µ¶2, 2 ¯¯2, 2 ¯¯Fz=”Fz”(13.354)¯¯11, 1CHAPTER 13. ARITHMETICAL ALGORITHMS325is a sum of perfect squares if z is a square.Obviously,µF¶µµ¶¶a ba+1 b+1a+2 b+2a, b ¯¯= 1+z 1+z 1+z (1 + . . .)¯zc1 c2 c+13 c+2(13.355)so by formula 13.355 hypergeometric functions with rational arguments can be computed with the binarysplitting method described in section 13.13.Hypergeometric functions can have any number of entries:µ¶∞Xa1 , . .

. , am ¯¯ak̄1 . . . ak̄m z kF=¯zb1 , . . . , bnbk̄1 . . . bk̄n k!(13.356)k=0Negative integer entries in the upper row lead to polynomials:¶µ−3, 3 ¯¯= 1 − 9 z + 18 z 2 − 10 z 3F¯z1(13.357)Negative integer entries in the lower row are verboten.13.17.2TransformationsAs obvious from the definition, entries can be swapped (capitalized symbols for readability):µ¶µ¶A, B, c ¯¯B, A, c ¯¯F= F¯z¯ze, f, ge, f, g(13.358)Usually one writes the entries in ascending order. Identical elements in the lower and upper row can becanceled:µ¶µ¶a, b, C ¯¯a, b ¯¯F= F(13.359)¯z¯ze, f, Ce, fThese are true for any number of elements, the following transformations are only valid for the givenstructure.Pfaff ’s reflection law and Euler’s identity(the latter is obtained by applying the reflection on both upper entries):µ¶µ¶1a, b ¯¯ −za, c − b ¯¯F= F¯¯z(1 − z)ac 1−zcµ¶µ¶a, c − b ¯¯ −za, b ¯¯1FF=¯¯z(1 − z)ac1−zcµ¶¯1c − a, b ¯ −z=F¯(1 − z)bc1−zµ¶µ¶¯a, b ¯c − a, c − b ¯¯(c−a−b)F= (1 − z)F¯z¯zcc(13.360)(13.361)(13.362)(13.363)A transformation by C.F.Gaussµ¶2a, 2b ¯¯F¯za + b + 12µ¶a, b ¯¯Fz¯a + b + 12µ¶a, b ¯¯1= Fwhere |z| <1 ¯ 4z(1 − z)a+b+ 22√¶µ¯2a, 2b ¯ 1 − 1 − z= F¯2a + b + 12(13.364)(13.365)CHAPTER 13.

ARITHMETICAL ALGORITHMS326Note that the right hand side of relation 13.364 does not change if z is replaced by 1 − z, so it seems thatµ¶µ¶2a, 2b ¯¯2a, 2b ¯¯Fz=F1−z(13.366)¯¯a + b + 12a + b + 12However, the relation is only true for terminating series, i.e. polynomials.

Rewriting relation 13.364 forthe argument 1−z2 we get¶¶µµ2a, 2b ¯¯ 1 − za, b ¯¯2(13.367)F= F¯¯1 − z2a + b + 12a + b + 12Ã!√¶µ2a, 2b ¯¯ 1 − 1 − z 2a, b ¯¯ 2F(13.368)= F¯¯za + b + 212a + b + 12Whipple’s identityµ1F¶µ¯ ¶+ 12 , 1 − a − b − c ¯¯ −4za, b, c¯a=(1−z)F¯¯z1 + a − b, 1 + a − c(1 − z)21 + a − b, 1 + a − c12 a, 2 a(13.369)Specializing 13.369 c = (a + 1)/2 (note the symmetry between b an c so specializing for c = (b + 1)/2produces the identical relation)µ1 1¶µ¶a, b ¯¯1a, 2 a + 12 − b ¯¯ −4z2F(13.370)F=¯z¯(1 − z)a(1 − z)21+a−b1+a−b¡¢2 !√√µ¶µ¶2a Ã1− 1−z2a, a − b + 12 ¯¯a, b ¯¯2 (1 − 1 − z)FF=(13.371)¯z¯−a + b + 12zza + b + 21With c := a − b in 13.370 one gets:µ¶a, a − c ¯¯F=¯z1+c1F(1 − z)aIf b = a in relation 13.370, then³ a, a ¯ ´¯F¯z11F(1 − z)a=µ1¶− 12 a + c ¯¯ −4z¯1+c(1 − z)2(13.372)¶− 12 a ¯¯ −4z¯1(1 − z)2(13.373)12 a, 2µ112 a, 2Observing that the hypergeometric function on the right side does not change when replacing a by 1 − aone finds Euler’s transformation for a = b, c = 1.Similarly as for the relations by Gauss,µ¶a, b ¯¯1−zF¯−1+z1+a−bÃ!√a, b ¯¯1 − 1 − z2√F¯−1+a−b1 + 1 − z2¶µa, b ¯¯2F¯1 − za + b + 21µ¶a, b ¯¯ 2Fz¯a + b + 21from 13.370 and 13.371:µ¶a µ 1 1¶1+za, 2 a + 21 − b ¯¯22=F(13.374)¯1 − z21+a−b!Ã√aµ1 1¶1 + 1 − z2a, a + 21 − b ¯¯ 2=F 2 2z(13.375)¯21+a−b!µ¶2a Ã2a, a − b + 12 ¯¯21−z=(13.376)F¯−1+z1+za + b + 12!√µ¶2a Ã2a, a − b + 12 ¯¯21 − 1 − z2√√=F(13.377)¯−a + b + 121 + 1 − z21 + 1 − z2The following is due to Ramanujan.

LetÃz0:=1−µ1−z1 + 2z¶3 !1/3(13.378)CHAPTER 13. ARITHMETICAL ALGORITHMS327thenz =µ1 2¯ ¶, ¯F 3 3 ¯ z 03=11 − (1 − z 03 )1/31 + 2 (1 − z 03 )1/3µ1 2¯ ¶, ¯(1 + 2 z) F 3 3 ¯ z 31(13.379)(13.380)Clausen’s product formulasµ¶¸2a, b ¯¯=z¯a + b + 12µ1¶ µ1¯ ¶1+ a, 41 + b ¯¯4 − a, 4 − b ¯ zF 4zF=¯¯1+a+b1−a−b· µFIf in 13.381 a = b +12F2a, a + b, 2b ¯¯¯za + b + 12 , 2a + 2bµ1¶+ a − b, 12 − a + b ¯¯z¯1 + a + b, 1 − a − b12, 2F¶(13.381)(13.382)then (two entries on the right hand side cancel)· µ¶¸2b + 12 , b ¯¯F¯z2b + 1µ= F¶2b + 12 , 2b ¯¯¯z4b + 1(13.383)and the right hand side again matches the structure on the left. The corresponding function can be´−2b³ √m.

One has Gn m (z) = (Gn (z)) .identified (see. [8]) as Gb (z) = 1+ 21−zSpecializing relation 13.382 for b = −a we get· µ1¶¸2+ a, 41 − a ¯¯F 4¯z1µ1=2F¶+ 2a, 21 , 12 − 2a ¯¯¯z1, 1(13.384)A complementary relationA relation involving arguments z and 1 − z (given in [22], p.291) isµ¶µµ¶¯ ¶a, b ¯¯a, bc − a, c − b ¯¯¯c−a−bFz=αFz+β(1−z)F1−z(13.385)¯¯¯ca+b−c+1c−a−b+1whereα = γ1 Γ(a) Γ(b) Γ(c) Γ(c − a − b)β = γ1 Γ(c) Γ(c − a) Γ(c − b) Γ(a + b − c) andγ = Γ(a) Γ(b) Γ(c − a) Γ(c − b).An integral representationFor z ∈ C − [1, ∞) one hasµ¶a, b ¯¯F=¯zcΓ(c)Γ(b) Γ(c − b)Z01 b−1tc−b−1(1 − t)a(1 − t z)dt(13.386)Differential equation³ ¯ ´¯f (z) = F a,c b ¯ z is a solution of the differential equationz (1 − z)d2 fdf+ [c − (1 + a + b) z]− abf2dzdz= 0(13.387)CHAPTER 13. ARITHMETICAL ALGORITHMS328A summation formulaSee [13].µF¶a, b ¯¯=¯1cΓ(c) Γ(c − a − b)Γ(c − a) Γ(c − b)if <(c − a − b) > 0 or b ∈ N, b < 0(13.388)See [27] chapter 15, [63] and [39].13.17.3Examples: elementary functionsAs a warmup the ‘well known’ functions like exp, log, sin, .

. . are expressed as hypergeometric functions.In some cases a transform is applied to give an alternative series.Simple series, exp and log1(1 − z)a=(1 + z)a=¶∞ µ³a¯ ´Xa+k−1 k¯z¯z =kk=0µ¶a µ ¶X−a ¯¯a kFz¯ −z =kF(13.389)(13.390)k=0exp(z) =log(1 + z) =∞³¯ ´Xzk¯¯z =k!k=0µ¶∞¯X1, 1 ¯(−1)k z k+1F¯ −z =2k+1F(13.391)(13.392)k=0erf(z)=2z√ FπÃ1232¯¯¯ − z2!(13.393)Trigonometric and hyperbolic functionssin(z) =zFsinh(z) =zF=zFcos(z) =zFcosh(z) =zFµ ¯¶∞X(−1)k z 2k+1¯ −z 2=3¯4(2k + 1)!2k=0µ ¯ 2¶∞Xz 2k+1¯z=¯34(2k + 1)!2k=0Ã!√1, 1 ¯¯ 1 − 1 − z 2by 13.3653 ¯22µ ¯¶∞X(−1)k z 2k¯ −z 2=¯14(2k)!2k=0µ ¯ 2¶∞X z 2k¯z=1¯4(2k)!2k=0(13.394)(13.395)(13.396)(13.397)(13.398)CHAPTER 13.

ARITHMETICAL ALGORITHMSÃsin(a z)= FÃcos(a z)= FÃsinh(a z)= FÃcosh(a z)= F3291+a 1−a2 , 232¯¯2¯ (sin(z))+ a2 , − a2 ¯¯2¯ (sin(z))1!(13.399)!(13.400)21+a 1−a2 , 232¯¯2¯ (sinh(z))+ a2 , − a2 ¯¯2¯ (sinh(z))1!(13.401)!(13.402)2Inverse trigonometric and hyperbolic functionsÃarcsin(z) =µ(arcsin(z))2=arcsinh(z) =zFzF1 12, 232¯¯ 2¯z!(13.403)¶1, 1, 1 ¯¯ 2z¯32, 2by 13.381Ã!1 1¯p,¯222log(z + 1 + z 2 ) = z F3 ¯ −z(13.404)(13.405)2==Ã!¯ z21,1z¯√F 23 ¯by 13.3611 + z21 + z22!Ã√1, 1 ¯¯ 1 − 1 + z 2zFby 13.3653 ¯22Ãarctan(z)= zF==¯12, 1¯3 ¯2−z=∞X(−1)k z 2k+1k=02k + 1Ã!1 1¯2,zz¯√F 232¯by 13.3611 + z21 + z22µ¶z1, 1 ¯¯ z 2Fby 13.361¯31 + z21 + z22Ãarctanh(z)!2= zFarccot(z) ===1Fz¯12, 1¯ 23 ¯z2Ã!=!1− 2zÃ1 1¯, ¯F 232¯∞Xz 2k+12k + 1(13.407)(13.408)(13.409)(13.410)(13.411)k=0¯12, 1¯3 ¯2z√1 + z2(13.406)11+z22µ¶z1, 1 ¯¯ 1F3 ¯21+z1 + z22(13.412)!(13.413)(13.414)CHAPTER 13.

ARITHMETICAL ALGORITHMS13.17.4330Elliptic K and EIn order to avoid the factorπ2we let K̃ :=2Kπ ,Ẽ :=2Eπ .µ1K̃(k)= FẼ(k)= F12, 2µ1¯ ¶¯ 2¯k− 21 ,112(13.415)¯ ¶¯ 2¯k(13.416)We further setµÑ (k) =F¶− 12 , − 21 ¯¯ 2¯k1(13.417)Most of the following relations can be written in several ways by one of the identitiespk0 =1 − k2¶2¶2µµ−4k 22k1 + k2=−=1−(1 − k 2 )21 − k2k0 2µ¶22k1 − k0k2=−=−−1 − k2k0k0 2µ¶222k1 + k01 − k0(1 + k 0 )√−= −=−=−21 − k01 − k0 21 − 1 − k2(1 − k 0 )1 − k02=1+1 + k01 + k0(13.418)(13.419)(13.420)(13.421)(13.422)Elliptic Kµ1K̃(k) ==K̃(k)12, 2F1Fk0=1ï ¶¯ 2¯k(13.423)¯ µk¯¯−1k01 12, 21F1 − k0µ1¶2 !by 13.361¶¯1 + k0¯¯−11 − k012, 2(13.424)(13.425)From relation 13.229 we get:K̃(k)With z(k) :=1−k01+k0ö !1 1¯ µ0 2,21−k¯=F 2 2¯11 + k01 + k0¶ Ã1 1¯ µ¶2 !µ1 − k0, 2 ¯ 1 − k02F=1+¯1 + k01 + k01(13.426)(13.427)(and (z(k))2 =: z 2 (k)) the last relation can be written as:K̃(k) =¡¢(1 + z(k)) K̃ z 2 (k)(13.428)CHAPTER 13. ARITHMETICAL ALGORITHMS331Elliptic EµẼ(k) =Ẽ(k) =¶− 21 , 12 ¯¯ 2F¯k1õ ¶2 !1 1¯−,k¯2 2 −k0 F¯1k0(13.429)by 13.361(13.430)The following relation resembles relation 13.426:Ã!2 ¢  1√1¯21 − k2−,−1−1−k¯2√F 2¯121 + 1 − k2õ¶2 !1 + k0− 12 , − 12 ¯¯ 1 − k 0F¯211 + k0¡¢−1(1 + z(k))Ñ z 2 (k)¡Ẽ(k) ===1+√1√F1 − k2K̃(k) =µ114, 41¯¯¯−4k 2(1 − k 2 )2(13.431)(13.432)(13.433)¶(13.434)Euler’s transform on 13.434 gives:K̃(k) =1 1 + k2Fk0 1 − k2µ1K̃(k)Ñ (k) ==== F14, 41ï µ 2k ¶2¯¯−1k0 23 34, 4¶¯¯¯ (2 k k 0 )2(13.435)by 13.364!− 14 , 34 ¯¯ −4 k 2¯21(1 − k 2 )Ã!µ¶2− 14 , 34 ¯¯2kk0 F−¯1k0 2õ¶2 !p− 14 , 14 ¯¯2k21+k F¯1k2 + 1p1 − k2 F!(13.436)Ãby 13.373(13.437)(13.438)by 13.361(13.439)332APPENDIX A.