1import math 2import pytest 3from mpmath import * 4 5def test_bessel(): 6 mp.dps = 15 7 assert j0(1).ae(0.765197686557966551) 8 assert j0(pi).ae(-0.304242177644093864) 9 assert j0(1000).ae(0.0247866861524201746) 10 assert j0(-25).ae(0.0962667832759581162) 11 assert j1(1).ae(0.440050585744933516) 12 assert j1(pi).ae(0.284615343179752757) 13 assert j1(1000).ae(0.00472831190708952392) 14 assert j1(-25).ae(0.125350249580289905) 15 assert besselj(5,1).ae(0.000249757730211234431) 16 assert besselj(5+0j,1).ae(0.000249757730211234431) 17 assert besselj(5,pi).ae(0.0521411843671184747) 18 assert besselj(5,1000).ae(0.00502540694523318607) 19 assert besselj(5,-25).ae(0.0660079953984229934) 20 assert besselj(-3,2).ae(-0.128943249474402051) 21 assert besselj(-4,2).ae(0.0339957198075684341) 22 assert besselj(3,3+2j).ae(0.424718794929639595942 + 0.625665327745785804812j) 23 assert besselj(0.25,4).ae(-0.374760630804249715) 24 assert besselj(1+2j,3+4j).ae(0.319247428741872131 - 0.669557748880365678j) 25 assert (besselj(3, 10**10) * 10**5).ae(0.76765081748139204023) 26 assert bessely(-0.5, 0) == 0 27 assert bessely(0.5, 0) == -inf 28 assert bessely(1.5, 0) == -inf 29 assert bessely(0,0) == -inf 30 assert bessely(-0.4, 0) == -inf 31 assert bessely(-0.6, 0) == inf 32 assert bessely(-1, 0) == inf 33 assert bessely(-1.4, 0) == inf 34 assert bessely(-1.6, 0) == -inf 35 assert bessely(-1, 0) == inf 36 assert bessely(-2, 0) == -inf 37 assert bessely(-3, 0) == inf 38 assert bessely(0.5, 0) == -inf 39 assert bessely(1, 0) == -inf 40 assert bessely(1.5, 0) == -inf 41 assert bessely(2, 0) == -inf 42 assert bessely(2.5, 0) == -inf 43 assert bessely(3, 0) == -inf 44 assert bessely(0,0.5).ae(-0.44451873350670655715) 45 assert bessely(1,0.5).ae(-1.4714723926702430692) 46 assert bessely(-1,0.5).ae(1.4714723926702430692) 47 assert bessely(3.5,0.5).ae(-138.86400867242488443) 48 assert bessely(0,3+4j).ae(4.6047596915010138655-8.8110771408232264208j) 49 assert bessely(0,j).ae(-0.26803248203398854876+1.26606587775200833560j) 50 assert (bessely(3, 10**10) * 10**5).ae(0.21755917537013204058) 51 assert besseli(0,0) == 1 52 assert besseli(1,0) == 0 53 assert besseli(2,0) == 0 54 assert besseli(-1,0) == 0 55 assert besseli(-2,0) == 0 56 assert besseli(0,0.5).ae(1.0634833707413235193) 57 assert besseli(1,0.5).ae(0.25789430539089631636) 58 assert besseli(-1,0.5).ae(0.25789430539089631636) 59 assert besseli(3.5,0.5).ae(0.00068103597085793815863) 60 assert besseli(0,3+4j).ae(-3.3924877882755196097-1.3239458916287264815j) 61 assert besseli(0,j).ae(besselj(0,1)) 62 assert (besseli(3, 10**10) * mpf(10)**(-4342944813)).ae(4.2996028505491271875) 63 assert besselk(0,0) == inf 64 assert besselk(1,0) == inf 65 assert besselk(2,0) == inf 66 assert besselk(-1,0) == inf 67 assert besselk(-2,0) == inf 68 assert besselk(0,0.5).ae(0.92441907122766586178) 69 assert besselk(1,0.5).ae(1.6564411200033008937) 70 assert besselk(-1,0.5).ae(1.6564411200033008937) 71 assert besselk(3.5,0.5).ae(207.48418747548460607) 72 assert besselk(0,3+4j).ae(-0.007239051213570155013+0.026510418350267677215j) 73 assert besselk(0,j).ae(-0.13863371520405399968-1.20196971531720649914j) 74 assert (besselk(3, 10**10) * mpf(10)**4342944824).ae(1.1628981033356187851) 75 # test for issue 331, bug reported by Michael Hartmann 76 for n in range(10,100,10): 77 mp.dps = n 78 assert besseli(91.5,24.7708).ae("4.00830632138673963619656140653537080438462342928377020695738635559218797348548092636896796324190271316137982810144874264e-41") 79 80def test_bessel_zeros(): 81 mp.dps = 15 82 assert besseljzero(0,1).ae(2.40482555769577276869) 83 assert besseljzero(2,1).ae(5.1356223018406825563) 84 assert besseljzero(1,50).ae(157.86265540193029781) 85 assert besseljzero(10,1).ae(14.475500686554541220) 86 assert besseljzero(0.5,3).ae(9.4247779607693797153) 87 assert besseljzero(2,1,1).ae(3.0542369282271403228) 88 assert besselyzero(0,1).ae(0.89357696627916752158) 89 assert besselyzero(2,1).ae(3.3842417671495934727) 90 assert besselyzero(1,50).ae(156.29183520147840108) 91 assert besselyzero(10,1).ae(12.128927704415439387) 92 assert besselyzero(0.5,3).ae(7.8539816339744830962) 93 assert besselyzero(2,1,1).ae(5.0025829314460639452) 94 95def test_hankel(): 96 mp.dps = 15 97 assert hankel1(0,0.5).ae(0.93846980724081290423-0.44451873350670655715j) 98 assert hankel1(1,0.5).ae(0.2422684576748738864-1.4714723926702430692j) 99 assert hankel1(-1,0.5).ae(-0.2422684576748738864+1.4714723926702430692j) 100 assert hankel1(1.5,0.5).ae(0.0917016996256513026-2.5214655504213378514j) 101 assert hankel1(1.5,3+4j).ae(0.0066806866476728165382-0.0036684231610839127106j) 102 assert hankel2(0,0.5).ae(0.93846980724081290423+0.44451873350670655715j) 103 assert hankel2(1,0.5).ae(0.2422684576748738864+1.4714723926702430692j) 104 assert hankel2(-1,0.5).ae(-0.2422684576748738864-1.4714723926702430692j) 105 assert hankel2(1.5,0.5).ae(0.0917016996256513026+2.5214655504213378514j) 106 assert hankel2(1.5,3+4j).ae(14.783528526098567526-7.397390270853446512j) 107 108def test_struve(): 109 mp.dps = 15 110 assert struveh(2,3).ae(0.74238666967748318564) 111 assert struveh(-2.5,3).ae(0.41271003220971599344) 112 assert struvel(2,3).ae(1.7476573277362782744) 113 assert struvel(-2.5,3).ae(1.5153394466819651377) 114 115def test_whittaker(): 116 mp.dps = 15 117 assert whitm(2,3,4).ae(49.753745589025246591) 118 assert whitw(2,3,4).ae(14.111656223052932215) 119 120def test_kelvin(): 121 mp.dps = 15 122 assert ber(2,3).ae(0.80836846563726819091) 123 assert ber(3,4).ae(-0.28262680167242600233) 124 assert ber(-3,2).ae(-0.085611448496796363669) 125 assert bei(2,3).ae(-0.89102236377977331571) 126 assert bei(-3,2).ae(-0.14420994155731828415) 127 assert ker(2,3).ae(0.12839126695733458928) 128 assert ker(-3,2).ae(-0.29802153400559142783) 129 assert ker(0.5,3).ae(-0.085662378535217097524) 130 assert kei(2,3).ae(0.036804426134164634000) 131 assert kei(-3,2).ae(0.88682069845786731114) 132 assert kei(0.5,3).ae(0.013633041571314302948) 133 134def test_hyper_misc(): 135 mp.dps = 15 136 assert hyp0f1(1,0) == 1 137 assert hyp1f1(1,2,0) == 1 138 assert hyp1f2(1,2,3,0) == 1 139 assert hyp2f1(1,2,3,0) == 1 140 assert hyp2f2(1,2,3,4,0) == 1 141 assert hyp2f3(1,2,3,4,5,0) == 1 142 # Degenerate case: 0F0 143 assert hyper([],[],0) == 1 144 assert hyper([],[],-2).ae(exp(-2)) 145 # Degenerate case: 1F0 146 assert hyper([2],[],1.5) == 4 147 # 148 assert hyp2f1((1,3),(2,3),(5,6),mpf(27)/32).ae(1.6) 149 assert hyp2f1((1,4),(1,2),(3,4),mpf(80)/81).ae(1.8) 150 assert hyp2f1((2,3),(1,1),(3,2),(2+j)/3).ae(1.327531603558679093+0.439585080092769253j) 151 mp.dps = 25 152 v = mpc('1.2282306665029814734863026', '-0.1225033830118305184672133') 153 assert hyper([(3,4),2+j,1],[1,5,j/3],mpf(1)/5+j/8).ae(v) 154 mp.dps = 15 155 156def test_elliptic_integrals(): 157 mp.dps = 15 158 assert ellipk(0).ae(pi/2) 159 assert ellipk(0.5).ae(gamma(0.25)**2/(4*sqrt(pi))) 160 assert ellipk(1) == inf 161 assert ellipk(1+0j) == inf 162 assert ellipk(-1).ae('1.3110287771460599052') 163 assert ellipk(-2).ae('1.1714200841467698589') 164 assert isinstance(ellipk(-2), mpf) 165 assert isinstance(ellipe(-2), mpf) 166 assert ellipk(-50).ae('0.47103424540873331679') 167 mp.dps = 30 168 n1 = +fraction(99999,100000) 169 n2 = +fraction(100001,100000) 170 mp.dps = 15 171 assert ellipk(n1).ae('7.1427724505817781901') 172 assert ellipk(n2).ae(mpc('7.1427417367963090109', '-1.5707923998261688019')) 173 assert ellipe(n1).ae('1.0000332138990829170') 174 v = ellipe(n2) 175 assert v.real.ae('0.999966786328145474069137') 176 assert (v.imag*10**6).ae('7.853952181727432') 177 assert ellipk(2).ae(mpc('1.3110287771460599052', '-1.3110287771460599052')) 178 assert ellipk(50).ae(mpc('0.22326753950210985451', '-0.47434723226254522087')) 179 assert ellipk(3+4j).ae(mpc('0.91119556380496500866', '0.63133428324134524388')) 180 assert ellipk(3-4j).ae(mpc('0.91119556380496500866', '-0.63133428324134524388')) 181 assert ellipk(-3+4j).ae(mpc('0.95357894880405122483', '0.23093044503746114444')) 182 assert ellipk(-3-4j).ae(mpc('0.95357894880405122483', '-0.23093044503746114444')) 183 assert isnan(ellipk(nan)) 184 assert isnan(ellipe(nan)) 185 assert ellipk(inf) == 0 186 assert isinstance(ellipk(inf), mpc) 187 assert ellipk(-inf) == 0 188 assert ellipk(1+0j) == inf 189 assert ellipe(0).ae(pi/2) 190 assert ellipe(0.5).ae(pi**(mpf(3)/2)/gamma(0.25)**2 +gamma(0.25)**2/(8*sqrt(pi))) 191 assert ellipe(1) == 1 192 assert ellipe(1+0j) == 1 193 assert ellipe(inf) == mpc(0,inf) 194 assert ellipe(-inf) == inf 195 assert ellipe(3+4j).ae(1.4995535209333469543-1.5778790079127582745j) 196 assert ellipe(3-4j).ae(1.4995535209333469543+1.5778790079127582745j) 197 assert ellipe(-3+4j).ae(2.5804237855343377803-0.8306096791000413778j) 198 assert ellipe(-3-4j).ae(2.5804237855343377803+0.8306096791000413778j) 199 assert ellipe(2).ae(0.59907011736779610372+0.59907011736779610372j) 200 assert ellipe('1e-1000000000').ae(pi/2) 201 assert ellipk('1e-1000000000').ae(pi/2) 202 assert ellipe(-pi).ae(2.4535865983838923) 203 mp.dps = 50 204 assert ellipk(1/pi).ae('1.724756270009501831744438120951614673874904182624739673') 205 assert ellipe(1/pi).ae('1.437129808135123030101542922290970050337425479058225712') 206 assert ellipk(-10*pi).ae('0.5519067523886233967683646782286965823151896970015484512') 207 assert ellipe(-10*pi).ae('5.926192483740483797854383268707108012328213431657645509') 208 v = ellipk(pi) 209 assert v.real.ae('0.973089521698042334840454592642137667227167622330325225') 210 assert v.imag.ae('-1.156151296372835303836814390793087600271609993858798016') 211 v = ellipe(pi) 212 assert v.real.ae('0.4632848917264710404078033487934663562998345622611263332') 213 assert v.imag.ae('1.0637961621753130852473300451583414489944099504180510966') 214 mp.dps = 15 215 216def test_exp_integrals(): 217 mp.dps = 15 218 x = +e 219 z = e + sqrt(3)*j 220 assert ei(x).ae(8.21168165538361560) 221 assert li(x).ae(1.89511781635593676) 222 assert si(x).ae(1.82104026914756705) 223 assert ci(x).ae(0.213958001340379779) 224 assert shi(x).ae(4.11520706247846193) 225 assert chi(x).ae(4.09647459290515367) 226 assert fresnels(x).ae(0.437189718149787643) 227 assert fresnelc(x).ae(0.401777759590243012) 228 assert airyai(x).ae(0.0108502401568586681) 229 assert airybi(x).ae(8.98245748585468627) 230 assert ei(z).ae(3.72597969491314951 + 7.34213212314224421j) 231 assert li(z).ae(2.28662658112562502 + 1.50427225297269364j) 232 assert si(z).ae(2.48122029237669054 + 0.12684703275254834j) 233 assert ci(z).ae(0.169255590269456633 - 0.892020751420780353j) 234 assert shi(z).ae(1.85810366559344468 + 3.66435842914920263j) 235 assert chi(z).ae(1.86787602931970484 + 3.67777369399304159j) 236 assert fresnels(z/3).ae(0.034534397197008182 + 0.754859844188218737j) 237 assert fresnelc(z/3).ae(1.261581645990027372 + 0.417949198775061893j) 238 assert airyai(z).ae(-0.0162552579839056062 - 0.0018045715700210556j) 239 assert airybi(z).ae(-4.98856113282883371 + 2.08558537872180623j) 240 assert li(0) == 0.0 241 assert li(1) == -inf 242 assert li(inf) == inf 243 assert isinstance(li(0.7), mpf) 244 assert si(inf).ae(pi/2) 245 assert si(-inf).ae(-pi/2) 246 assert ci(inf) == 0 247 assert ci(0) == -inf 248 assert isinstance(ei(-0.7), mpf) 249 assert airyai(inf) == 0 250 assert airybi(inf) == inf 251 assert airyai(-inf) == 0 252 assert airybi(-inf) == 0 253 assert fresnels(inf) == 0.5 254 assert fresnelc(inf) == 0.5 255 assert fresnels(-inf) == -0.5 256 assert fresnelc(-inf) == -0.5 257 assert shi(0) == 0 258 assert shi(inf) == inf 259 assert shi(-inf) == -inf 260 assert chi(0) == -inf 261 assert chi(inf) == inf 262 263def test_ei(): 264 mp.dps = 15 265 assert ei(0) == -inf 266 assert ei(inf) == inf 267 assert ei(-inf) == -0.0 268 assert ei(20+70j).ae(6.1041351911152984397e6 - 2.7324109310519928872e6j) 269 # tests for the asymptotic expansion 270 # values checked with Mathematica ExpIntegralEi 271 mp.dps = 50 272 r = ei(20000) 273 s = '3.8781962825045010930273870085501819470698476975019e+8681' 274 assert str(r) == s 275 r = ei(-200) 276 s = '-6.8852261063076355977108174824557929738368086933303e-90' 277 assert str(r) == s 278 r =ei(20000 + 10*j) 279 sre = '-3.255138234032069402493850638874410725961401274106e+8681' 280 sim = '-2.1081929993474403520785942429469187647767369645423e+8681' 281 assert str(r.real) == sre and str(r.imag) == sim 282 mp.dps = 15 283 # More asymptotic expansions 284 assert chi(-10**6+100j).ae('1.3077239389562548386e+434288 + 7.6808956999707408158e+434287j') 285 assert shi(-10**6+100j).ae('-1.3077239389562548386e+434288 - 7.6808956999707408158e+434287j') 286 mp.dps = 15 287 assert ei(10j).ae(-0.0454564330044553726+3.2291439210137706686j) 288 assert ei(100j).ae(-0.0051488251426104921+3.1330217936839529126j) 289 u = ei(fmul(10**20, j, exact=True)) 290 assert u.real.ae(-6.4525128526578084421345e-21, abs_eps=0, rel_eps=8*eps) 291 assert u.imag.ae(pi) 292 assert ei(-10j).ae(-0.0454564330044553726-3.2291439210137706686j) 293 assert ei(-100j).ae(-0.0051488251426104921-3.1330217936839529126j) 294 u = ei(fmul(-10**20, j, exact=True)) 295 assert u.real.ae(-6.4525128526578084421345e-21, abs_eps=0, rel_eps=8*eps) 296 assert u.imag.ae(-pi) 297 assert ei(10+10j).ae(-1576.1504265768517448+436.9192317011328140j) 298 u = ei(-10+10j) 299 assert u.real.ae(7.6698978415553488362543e-7, abs_eps=0, rel_eps=8*eps) 300 assert u.imag.ae(3.141595611735621062025) 301 302def test_e1(): 303 mp.dps = 15 304 assert e1(0) == inf 305 assert e1(inf) == 0 306 assert e1(-inf) == mpc(-inf, -pi) 307 assert e1(10j).ae(0.045456433004455372635 + 0.087551267423977430100j) 308 assert e1(100j).ae(0.0051488251426104921444 - 0.0085708599058403258790j) 309 assert e1(fmul(10**20, j, exact=True)).ae(6.4525128526578084421e-21 - 7.6397040444172830039e-21j, abs_eps=0, rel_eps=8*eps) 310 assert e1(-10j).ae(0.045456433004455372635 - 0.087551267423977430100j) 311 assert e1(-100j).ae(0.0051488251426104921444 + 0.0085708599058403258790j) 312 assert e1(fmul(-10**20, j, exact=True)).ae(6.4525128526578084421e-21 + 7.6397040444172830039e-21j, abs_eps=0, rel_eps=8*eps) 313 314def test_expint(): 315 mp.dps = 15 316 assert expint(0,0) == inf 317 assert expint(0,1).ae(1/e) 318 assert expint(0,1.5).ae(2/exp(1.5)/3) 319 assert expint(1,1).ae(-ei(-1)) 320 assert expint(2,0).ae(1) 321 assert expint(3,0).ae(1/2.) 322 assert expint(4,0).ae(1/3.) 323 assert expint(-2, 0.5).ae(26/sqrt(e)) 324 assert expint(-1,-1) == 0 325 assert expint(-2,-1).ae(-e) 326 assert expint(5.5, 0).ae(2/9.) 327 assert expint(2.00000001,0).ae(100000000./100000001) 328 assert expint(2+3j,4-j).ae(0.0023461179581675065414+0.0020395540604713669262j) 329 assert expint('1.01', '1e-1000').ae(99.9999999899412802) 330 assert expint('1.000000000001', 3.5).ae(0.00697013985754701819446) 331 assert expint(2,3).ae(3*ei(-3)+exp(-3)) 332 assert (expint(10,20)*10**10).ae(0.694439055541231353) 333 assert expint(3,inf) == 0 334 assert expint(3.2,inf) == 0 335 assert expint(3.2+2j,inf) == 0 336 assert expint(1,3j).ae(-0.11962978600800032763 + 0.27785620120457163717j) 337 assert expint(1,3).ae(0.013048381094197037413) 338 assert expint(1,-3).ae(-ei(3)-pi*j) 339 #assert expint(3) == expint(1,3) 340 assert expint(1,-20).ae(-25615652.66405658882 - 3.1415926535897932385j) 341 assert expint(1000000,0).ae(1./999999) 342 assert expint(0,2+3j).ae(-0.025019798357114678171 + 0.027980439405104419040j) 343 assert expint(-1,2+3j).ae(-0.022411973626262070419 + 0.038058922011377716932j) 344 assert expint(-1.5,0) == inf 345 346def test_trig_integrals(): 347 mp.dps = 30 348 assert si(mpf(1)/1000000).ae('0.000000999999999999944444444444446111') 349 assert ci(mpf(1)/1000000).ae('-13.2382948930629912435014366276') 350 assert si(10**10).ae('1.5707963267075846569685111517747537') 351 assert ci(10**10).ae('-4.87506025174822653785729773959e-11') 352 assert si(10**100).ae(pi/2) 353 assert (ci(10**100)*10**100).ae('-0.372376123661276688262086695553') 354 assert si(-3) == -si(3) 355 assert ci(-3).ae(ci(3) + pi*j) 356 # Test complex structure 357 mp.dps = 15 358 assert mp.ci(50).ae(-0.0056283863241163054402) 359 assert mp.ci(50+2j).ae(-0.018378282946133067149+0.070352808023688336193j) 360 assert mp.ci(20j).ae(1.28078263320282943611e7+1.5707963267949j) 361 assert mp.ci(-2+20j).ae(-4.050116856873293505e6+1.207476188206989909e7j) 362 assert mp.ci(-50+2j).ae(-0.0183782829461330671+3.0712398455661049023j) 363 assert mp.ci(-50).ae(-0.0056283863241163054+3.1415926535897932385j) 364 assert mp.ci(-50-2j).ae(-0.0183782829461330671-3.0712398455661049023j) 365 assert mp.ci(-2-20j).ae(-4.050116856873293505e6-1.207476188206989909e7j) 366 assert mp.ci(-20j).ae(1.28078263320282943611e7-1.5707963267949j) 367 assert mp.ci(50-2j).ae(-0.018378282946133067149-0.070352808023688336193j) 368 assert mp.si(50).ae(1.5516170724859358947) 369 assert mp.si(50+2j).ae(1.497884414277228461-0.017515007378437448j) 370 assert mp.si(20j).ae(1.2807826332028294459e7j) 371 assert mp.si(-2+20j).ae(-1.20747603112735722103e7-4.050116856873293554e6j) 372 assert mp.si(-50+2j).ae(-1.497884414277228461-0.017515007378437448j) 373 assert mp.si(-50).ae(-1.5516170724859358947) 374 assert mp.si(-50-2j).ae(-1.497884414277228461+0.017515007378437448j) 375 assert mp.si(-2-20j).ae(-1.20747603112735722103e7+4.050116856873293554e6j) 376 assert mp.si(-20j).ae(-1.2807826332028294459e7j) 377 assert mp.si(50-2j).ae(1.497884414277228461+0.017515007378437448j) 378 assert mp.chi(50j).ae(-0.0056283863241163054+1.5707963267948966192j) 379 assert mp.chi(-2+50j).ae(-0.0183782829461330671+1.6411491348185849554j) 380 assert mp.chi(-20).ae(1.28078263320282943611e7+3.1415926535898j) 381 assert mp.chi(-20-2j).ae(-4.050116856873293505e6+1.20747571696809187053e7j) 382 assert mp.chi(-2-50j).ae(-0.0183782829461330671-1.6411491348185849554j) 383 assert mp.chi(-50j).ae(-0.0056283863241163054-1.5707963267948966192j) 384 assert mp.chi(2-50j).ae(-0.0183782829461330671-1.500443518771208283j) 385 assert mp.chi(20-2j).ae(-4.050116856873293505e6-1.20747603112735722951e7j) 386 assert mp.chi(20).ae(1.2807826332028294361e7) 387 assert mp.chi(2+50j).ae(-0.0183782829461330671+1.500443518771208283j) 388 assert mp.shi(50j).ae(1.5516170724859358947j) 389 assert mp.shi(-2+50j).ae(0.017515007378437448+1.497884414277228461j) 390 assert mp.shi(-20).ae(-1.2807826332028294459e7) 391 assert mp.shi(-20-2j).ae(4.050116856873293554e6-1.20747603112735722103e7j) 392 assert mp.shi(-2-50j).ae(0.017515007378437448-1.497884414277228461j) 393 assert mp.shi(-50j).ae(-1.5516170724859358947j) 394 assert mp.shi(2-50j).ae(-0.017515007378437448-1.497884414277228461j) 395 assert mp.shi(20-2j).ae(-4.050116856873293554e6-1.20747603112735722103e7j) 396 assert mp.shi(20).ae(1.2807826332028294459e7) 397 assert mp.shi(2+50j).ae(-0.017515007378437448+1.497884414277228461j) 398 def ae(x,y,tol=1e-12): 399 return abs(x-y) <= abs(y)*tol 400 assert fp.ci(fp.inf) == 0 401 assert ae(fp.ci(fp.ninf), fp.pi*1j) 402 assert ae(fp.si(fp.inf), fp.pi/2) 403 assert ae(fp.si(fp.ninf), -fp.pi/2) 404 assert fp.si(0) == 0 405 assert ae(fp.ci(50), -0.0056283863241163054402) 406 assert ae(fp.ci(50+2j), -0.018378282946133067149+0.070352808023688336193j) 407 assert ae(fp.ci(20j), 1.28078263320282943611e7+1.5707963267949j) 408 assert ae(fp.ci(-2+20j), -4.050116856873293505e6+1.207476188206989909e7j) 409 assert ae(fp.ci(-50+2j), -0.0183782829461330671+3.0712398455661049023j) 410 assert ae(fp.ci(-50), -0.0056283863241163054+3.1415926535897932385j) 411 assert ae(fp.ci(-50-2j), -0.0183782829461330671-3.0712398455661049023j) 412 assert ae(fp.ci(-2-20j), -4.050116856873293505e6-1.207476188206989909e7j) 413 assert ae(fp.ci(-20j), 1.28078263320282943611e7-1.5707963267949j) 414 assert ae(fp.ci(50-2j), -0.018378282946133067149-0.070352808023688336193j) 415 assert ae(fp.si(50), 1.5516170724859358947) 416 assert ae(fp.si(50+2j), 1.497884414277228461-0.017515007378437448j) 417 assert ae(fp.si(20j), 1.2807826332028294459e7j) 418 assert ae(fp.si(-2+20j), -1.20747603112735722103e7-4.050116856873293554e6j) 419 assert ae(fp.si(-50+2j), -1.497884414277228461-0.017515007378437448j) 420 assert ae(fp.si(-50), -1.5516170724859358947) 421 assert ae(fp.si(-50-2j), -1.497884414277228461+0.017515007378437448j) 422 assert ae(fp.si(-2-20j), -1.20747603112735722103e7+4.050116856873293554e6j) 423 assert ae(fp.si(-20j), -1.2807826332028294459e7j) 424 assert ae(fp.si(50-2j), 1.497884414277228461+0.017515007378437448j) 425 assert ae(fp.chi(50j), -0.0056283863241163054+1.5707963267948966192j) 426 assert ae(fp.chi(-2+50j), -0.0183782829461330671+1.6411491348185849554j) 427 assert ae(fp.chi(-20), 1.28078263320282943611e7+3.1415926535898j) 428 assert ae(fp.chi(-20-2j), -4.050116856873293505e6+1.20747571696809187053e7j) 429 assert ae(fp.chi(-2-50j), -0.0183782829461330671-1.6411491348185849554j) 430 assert ae(fp.chi(-50j), -0.0056283863241163054-1.5707963267948966192j) 431 assert ae(fp.chi(2-50j), -0.0183782829461330671-1.500443518771208283j) 432 assert ae(fp.chi(20-2j), -4.050116856873293505e6-1.20747603112735722951e7j) 433 assert ae(fp.chi(20), 1.2807826332028294361e7) 434 assert ae(fp.chi(2+50j), -0.0183782829461330671+1.500443518771208283j) 435 assert ae(fp.shi(50j), 1.5516170724859358947j) 436 assert ae(fp.shi(-2+50j), 0.017515007378437448+1.497884414277228461j) 437 assert ae(fp.shi(-20), -1.2807826332028294459e7) 438 assert ae(fp.shi(-20-2j), 4.050116856873293554e6-1.20747603112735722103e7j) 439 assert ae(fp.shi(-2-50j), 0.017515007378437448-1.497884414277228461j) 440 assert ae(fp.shi(-50j), -1.5516170724859358947j) 441 assert ae(fp.shi(2-50j), -0.017515007378437448-1.497884414277228461j) 442 assert ae(fp.shi(20-2j), -4.050116856873293554e6-1.20747603112735722103e7j) 443 assert ae(fp.shi(20), 1.2807826332028294459e7) 444 assert ae(fp.shi(2+50j), -0.017515007378437448+1.497884414277228461j) 445 446def test_airy(): 447 mp.dps = 15 448 assert (airyai(10)*10**10).ae(1.1047532552898687) 449 assert (airybi(10)/10**9).ae(0.45564115354822515) 450 assert (airyai(1000)*10**9158).ae(9.306933063179556004) 451 assert (airybi(1000)/10**9154).ae(5.4077118391949465477) 452 assert airyai(-1000).ae(0.055971895773019918842) 453 assert airybi(-1000).ae(-0.083264574117080633012) 454 assert (airyai(100+100j)*10**188).ae(2.9099582462207032076 + 2.353013591706178756j) 455 assert (airybi(100+100j)/10**185).ae(1.7086751714463652039 - 3.1416590020830804578j) 456 457def test_hyper_0f1(): 458 mp.dps = 15 459 v = 8.63911136507950465 460 assert hyper([],[(1,3)],1.5).ae(v) 461 assert hyper([],[1/3.],1.5).ae(v) 462 assert hyp0f1(1/3.,1.5).ae(v) 463 assert hyp0f1((1,3),1.5).ae(v) 464 # Asymptotic expansion 465 assert hyp0f1(3,1e9).ae('4.9679055380347771271e+27455') 466 assert hyp0f1(3,1e9j).ae('-2.1222788784457702157e+19410 + 5.0840597555401854116e+19410j') 467 468def test_hyper_1f1(): 469 mp.dps = 15 470 v = 1.2917526488617656673 471 assert hyper([(1,2)],[(3,2)],0.7).ae(v) 472 assert hyper([(1,2)],[(3,2)],0.7+0j).ae(v) 473 assert hyper([0.5],[(3,2)],0.7).ae(v) 474 assert hyper([0.5],[1.5],0.7).ae(v) 475 assert hyper([0.5],[(3,2)],0.7+0j).ae(v) 476 assert hyper([0.5],[1.5],0.7+0j).ae(v) 477 assert hyper([(1,2)],[1.5+0j],0.7).ae(v) 478 assert hyper([0.5+0j],[1.5],0.7).ae(v) 479 assert hyper([0.5+0j],[1.5+0j],0.7+0j).ae(v) 480 assert hyp1f1(0.5,1.5,0.7).ae(v) 481 assert hyp1f1((1,2),1.5,0.7).ae(v) 482 # Asymptotic expansion 483 assert hyp1f1(2,3,1e10).ae('2.1555012157015796988e+4342944809') 484 assert (hyp1f1(2,3,1e10j)*10**10).ae(-0.97501205020039745852 - 1.7462392454512132074j) 485 # Shouldn't use asymptotic expansion 486 assert hyp1f1(-2, 1, 10000).ae(49980001) 487 # Bug 488 assert hyp1f1(1j,fraction(1,3),0.415-69.739j).ae(25.857588206024346592 + 15.738060264515292063j) 489 490def test_hyper_2f1(): 491 mp.dps = 15 492 v = 1.0652207633823291032 493 assert hyper([(1,2), (3,4)], [2], 0.3).ae(v) 494 assert hyper([(1,2), 0.75], [2], 0.3).ae(v) 495 assert hyper([0.5, 0.75], [2.0], 0.3).ae(v) 496 assert hyper([0.5, 0.75], [2.0], 0.3+0j).ae(v) 497 assert hyper([0.5+0j, (3,4)], [2.0], 0.3+0j).ae(v) 498 assert hyper([0.5+0j, (3,4)], [2.0], 0.3).ae(v) 499 assert hyper([0.5, (3,4)], [2.0+0j], 0.3).ae(v) 500 assert hyper([0.5+0j, 0.75+0j], [2.0+0j], 0.3+0j).ae(v) 501 v = 1.09234681096223231717 + 0.18104859169479360380j 502 assert hyper([(1,2),0.75+j], [2], 0.5).ae(v) 503 assert hyper([0.5,0.75+j], [2.0], 0.5).ae(v) 504 assert hyper([0.5,0.75+j], [2.0], 0.5+0j).ae(v) 505 assert hyper([0.5,0.75+j], [2.0+0j], 0.5+0j).ae(v) 506 v = 0.9625 - 0.125j 507 assert hyper([(3,2),-1],[4], 0.1+j/3).ae(v) 508 assert hyper([1.5,-1.0],[4], 0.1+j/3).ae(v) 509 assert hyper([1.5,-1.0],[4+0j], 0.1+j/3).ae(v) 510 assert hyper([1.5+0j,-1.0+0j],[4+0j], 0.1+j/3).ae(v) 511 v = 1.02111069501693445001 - 0.50402252613466859521j 512 assert hyper([(2,10),(3,10)],[(4,10)],1.5).ae(v) 513 assert hyper([0.2,(3,10)],[0.4+0j],1.5).ae(v) 514 assert hyper([0.2,(3,10)],[0.4+0j],1.5+0j).ae(v) 515 v = 0.76922501362865848528 + 0.32640579593235886194j 516 assert hyper([(2,10),(3,10)],[(4,10)],4+2j).ae(v) 517 assert hyper([0.2,(3,10)],[0.4+0j],4+2j).ae(v) 518 assert hyper([0.2,(3,10)],[(4,10)],4+2j).ae(v) 519 520def test_hyper_2f1_hard(): 521 mp.dps = 15 522 # Singular cases 523 assert hyp2f1(2,-1,-1,3).ae(7) 524 assert hyp2f1(2,-1,-1,3,eliminate_all=True).ae(0.25) 525 assert hyp2f1(2,-2,-2,3).ae(34) 526 assert hyp2f1(2,-2,-2,3,eliminate_all=True).ae(0.25) 527 assert hyp2f1(2,-2,-3,3) == 14 528 assert hyp2f1(2,-3,-2,3) == inf 529 assert hyp2f1(2,-1.5,-1.5,3) == 0.25 530 assert hyp2f1(1,2,3,0) == 1 531 assert hyp2f1(0,1,0,0) == 1 532 assert hyp2f1(0,0,0,0) == 1 533 assert isnan(hyp2f1(1,1,0,0)) 534 assert hyp2f1(2,-1,-5, 0.25+0.25j).ae(1.1+0.1j) 535 assert hyp2f1(2,-5,-5, 0.25+0.25j, eliminate=False).ae(163./128 + 125./128*j) 536 assert hyp2f1(0.7235, -1, -5, 0.3).ae(1.04341) 537 assert hyp2f1(0.7235, -5, -5, 0.3, eliminate=False).ae(1.2939225017815903812) 538 assert hyp2f1(-1,-2,4,1) == 1.5 539 assert hyp2f1(1,2,-3,1) == inf 540 assert hyp2f1(-2,-2,1,1) == 6 541 assert hyp2f1(1,-2,-4,1).ae(5./3) 542 assert hyp2f1(0,-6,-4,1) == 1 543 assert hyp2f1(0,-3,-4,1) == 1 544 assert hyp2f1(0,0,0,1) == 1 545 assert hyp2f1(1,0,0,1,eliminate=False) == 1 546 assert hyp2f1(1,1,0,1) == inf 547 assert hyp2f1(1,-6,-4,1) == inf 548 assert hyp2f1(-7.2,-0.5,-4.5,1) == 0 549 assert hyp2f1(-7.2,-1,-2,1).ae(-2.6) 550 assert hyp2f1(1,-0.5,-4.5, 1) == inf 551 assert hyp2f1(1,0.5,-4.5, 1) == -inf 552 # Check evaluation on / close to unit circle 553 z = exp(j*pi/3) 554 w = (nthroot(2,3)+1)*exp(j*pi/12)/nthroot(3,4)**3 555 assert hyp2f1('1/2','1/6','1/3', z).ae(w) 556 assert hyp2f1('1/2','1/6','1/3', z.conjugate()).ae(w.conjugate()) 557 assert hyp2f1(0.25, (1,3), 2, '0.999').ae(1.06826449496030635) 558 assert hyp2f1(0.25, (1,3), 2, '1.001').ae(1.06867299254830309446-0.00001446586793975874j) 559 assert hyp2f1(0.25, (1,3), 2, -1).ae(0.96656584492524351673) 560 assert hyp2f1(0.25, (1,3), 2, j).ae(0.99041766248982072266+0.03777135604180735522j) 561 assert hyp2f1(2,3,5,'0.99').ae(27.699347904322690602) 562 assert hyp2f1((3,2),-0.5,3,'0.99').ae(0.68403036843911661388) 563 assert hyp2f1(2,3,5,1j).ae(0.37290667145974386127+0.59210004902748285917j) 564 assert fsum([hyp2f1((7,10),(2,3),(-1,2), 0.95*exp(j*k)) for k in range(1,15)]).ae(52.851400204289452922+6.244285013912953225j) 565 assert fsum([hyp2f1((7,10),(2,3),(-1,2), 1.05*exp(j*k)) for k in range(1,15)]).ae(54.506013786220655330-3.000118813413217097j) 566 assert fsum([hyp2f1((7,10),(2,3),(-1,2), exp(j*k)) for k in range(1,15)]).ae(55.792077935955314887+1.731986485778500241j) 567 assert hyp2f1(2,2.5,-3.25,0.999).ae(218373932801217082543180041.33) 568 # Branches 569 assert hyp2f1(1,1,2,1.01).ae(4.5595744415723676911-3.1104877758314784539j) 570 assert hyp2f1(1,1,2,1.01+0.1j).ae(2.4149427480552782484+1.4148224796836938829j) 571 assert hyp2f1(1,1,2,3+4j).ae(0.14576709331407297807+0.48379185417980360773j) 572 assert hyp2f1(1,1,2,4).ae(-0.27465307216702742285 - 0.78539816339744830962j) 573 assert hyp2f1(1,1,2,-4).ae(0.40235947810852509365) 574 # Other: 575 # Cancellation with a large parameter involved (bug reported on sage-devel) 576 assert hyp2f1(112, (51,10), (-9,10), -0.99999).ae(-1.6241361047970862961e-24, abs_eps=0, rel_eps=eps*16) 577 578def test_hyper_3f2_etc(): 579 assert hyper([1,2,3],[1.5,8],-1).ae(0.67108992351533333030) 580 assert hyper([1,2,3,4],[5,6,7], -1).ae(0.90232988035425506008) 581 assert hyper([1,2,3],[1.25,5], 1).ae(28.924181329701905701) 582 assert hyper([1,2,3,4],[5,6,7],5).ae(1.5192307344006649499-1.1529845225075537461j) 583 assert hyper([1,2,3,4,5],[6,7,8,9],-1).ae(0.96288759462882357253) 584 assert hyper([1,2,3,4,5],[6,7,8,9],1).ae(1.0428697385885855841) 585 assert hyper([1,2,3,4,5],[6,7,8,9],5).ae(1.33980653631074769423-0.07143405251029226699j) 586 assert hyper([1,2.79,3.08,4.37],[5.2,6.1,7.3],5).ae(1.0996321464692607231-1.7748052293979985001j) 587 assert hyper([1,1,1],[1,2],1) == inf 588 assert hyper([1,1,1],[2,(101,100)],1).ae(100.01621213528313220) 589 # slow -- covered by doctests 590 #assert hyper([1,1,1],[2,3],0.9999).ae(1.2897972005319693905) 591 592def test_hyper_u(): 593 mp.dps = 15 594 assert hyperu(2,-3,0).ae(0.05) 595 assert hyperu(2,-3.5,0).ae(4./99) 596 assert hyperu(2,0,0) == 0.5 597 assert hyperu(-5,1,0) == -120 598 assert hyperu(-5,2,0) == inf 599 assert hyperu(-5,-2,0) == 0 600 assert hyperu(7,7,3).ae(0.00014681269365593503986) #exp(3)*gammainc(-6,3) 601 assert hyperu(2,-3,4).ae(0.011836478100271995559) 602 assert hyperu(3,4,5).ae(1./125) 603 assert hyperu(2,3,0.0625) == 256 604 assert hyperu(-1,2,0.25+0.5j) == -1.75+0.5j 605 assert hyperu(0.5,1.5,7.25).ae(2/sqrt(29)) 606 assert hyperu(2,6,pi).ae(0.55804439825913399130) 607 assert (hyperu((3,2),8,100+201j)*10**4).ae(-0.3797318333856738798 - 2.9974928453561707782j) 608 assert (hyperu((5,2),(-1,2),-5000)*10**10).ae(-5.6681877926881664678j) 609 # XXX: fails because of undetected cancellation in low level series code 610 # Alternatively: could use asymptotic series here, if convergence test 611 # tweaked back to recognize this one 612 #assert (hyperu((5,2),(-1,2),-500)*10**7).ae(-1.82526906001593252847j) 613 614def test_hyper_2f0(): 615 mp.dps = 15 616 assert hyper([1,2],[],3) == hyp2f0(1,2,3) 617 assert hyp2f0(2,3,7).ae(0.0116108068639728714668 - 0.0073727413865865802130j) 618 assert hyp2f0(2,3,0) == 1 619 assert hyp2f0(0,0,0) == 1 620 assert hyp2f0(-1,-1,1).ae(2) 621 assert hyp2f0(-4,1,1.5).ae(62.5) 622 assert hyp2f0(-4,1,50).ae(147029801) 623 assert hyp2f0(-4,1,0.0001).ae(0.99960011997600240000) 624 assert hyp2f0(0.5,0.25,0.001).ae(1.0001251174078538115) 625 assert hyp2f0(0.5,0.25,3+4j).ae(0.85548875824755163518 + 0.21636041283392292973j) 626 # Important: cancellation check 627 assert hyp2f0((1,6),(5,6),-0.02371708245126284498).ae(0.996785723120804309) 628 # Should be exact; polynomial case 629 assert hyp2f0(-2,1,0.5+0.5j,zeroprec=200) == 0 630 assert hyp2f0(1,-2,0.5+0.5j,zeroprec=200) == 0 631 # There used to be a bug in thresholds that made one of the following hang 632 for d in [15, 50, 80]: 633 mp.dps = d 634 assert hyp2f0(1.5, 0.5, 0.009).ae('1.006867007239309717945323585695344927904000945829843527398772456281301440034218290443367270629519483 + 1.238277162240704919639384945859073461954721356062919829456053965502443570466701567100438048602352623e-46j') 635 636def test_hyper_1f2(): 637 mp.dps = 15 638 assert hyper([1],[2,3],4) == hyp1f2(1,2,3,4) 639 a1,b1,b2 = (1,10),(2,3),1./16 640 assert hyp1f2(a1,b1,b2,10).ae(298.7482725554557568) 641 assert hyp1f2(a1,b1,b2,100).ae(224128961.48602947604) 642 assert hyp1f2(a1,b1,b2,1000).ae(1.1669528298622675109e+27) 643 assert hyp1f2(a1,b1,b2,10000).ae(2.4780514622487212192e+86) 644 assert hyp1f2(a1,b1,b2,100000).ae(1.3885391458871523997e+274) 645 assert hyp1f2(a1,b1,b2,1000000).ae('9.8851796978960318255e+867') 646 assert hyp1f2(a1,b1,b2,10**7).ae('1.1505659189516303646e+2746') 647 assert hyp1f2(a1,b1,b2,10**8).ae('1.4672005404314334081e+8685') 648 assert hyp1f2(a1,b1,b2,10**20).ae('3.6888217332150976493e+8685889636') 649 assert hyp1f2(a1,b1,b2,10*j).ae(-16.163252524618572878 - 44.321567896480184312j) 650 assert hyp1f2(a1,b1,b2,100*j).ae(61938.155294517848171 + 637349.45215942348739j) 651 assert hyp1f2(a1,b1,b2,1000*j).ae(8455057657257695958.7 + 6261969266997571510.6j) 652 assert hyp1f2(a1,b1,b2,10000*j).ae(-8.9771211184008593089e+60 + 4.6550528111731631456e+59j) 653 assert hyp1f2(a1,b1,b2,100000*j).ae(2.6398091437239324225e+193 + 4.1658080666870618332e+193j) 654 assert hyp1f2(a1,b1,b2,1000000*j).ae('3.5999042951925965458e+613 + 1.5026014707128947992e+613j') 655 assert hyp1f2(a1,b1,b2,10**7*j).ae('-8.3208715051623234801e+1939 - 3.6752883490851869429e+1941j') 656 assert hyp1f2(a1,b1,b2,10**8*j).ae('2.0724195707891484454e+6140 - 1.3276619482724266387e+6141j') 657 assert hyp1f2(a1,b1,b2,10**20*j).ae('-1.1734497974795488504e+6141851462 + 1.1498106965385471542e+6141851462j') 658 659def test_hyper_2f3(): 660 mp.dps = 15 661 assert hyper([1,2],[3,4,5],6) == hyp2f3(1,2,3,4,5,6) 662 a1,a2,b1,b2,b3 = (1,10),(2,3),(3,10), 2, 1./16 663 # Check asymptotic expansion 664 assert hyp2f3(a1,a2,b1,b2,b3,10).ae(128.98207160698659976) 665 assert hyp2f3(a1,a2,b1,b2,b3,1000).ae(6.6309632883131273141e25) 666 assert hyp2f3(a1,a2,b1,b2,b3,10000).ae(4.6863639362713340539e84) 667 assert hyp2f3(a1,a2,b1,b2,b3,100000).ae(8.6632451236103084119e271) 668 assert hyp2f3(a1,a2,b1,b2,b3,10**6).ae('2.0291718386574980641e865') 669 assert hyp2f3(a1,a2,b1,b2,b3,10**7).ae('7.7639836665710030977e2742') 670 assert hyp2f3(a1,a2,b1,b2,b3,10**8).ae('3.2537462584071268759e8681') 671 assert hyp2f3(a1,a2,b1,b2,b3,10**20).ae('1.2966030542911614163e+8685889627') 672 assert hyp2f3(a1,a2,b1,b2,b3,10*j).ae(-18.551602185587547854 - 13.348031097874113552j) 673 assert hyp2f3(a1,a2,b1,b2,b3,100*j).ae(78634.359124504488695 + 74459.535945281973996j) 674 assert hyp2f3(a1,a2,b1,b2,b3,1000*j).ae(597682550276527901.59 - 65136194809352613.078j) 675 assert hyp2f3(a1,a2,b1,b2,b3,10000*j).ae(-1.1779696326238582496e+59 + 1.2297607505213133872e+59j) 676 assert hyp2f3(a1,a2,b1,b2,b3,100000*j).ae(2.9844228969804380301e+191 + 7.5587163231490273296e+190j) 677 assert hyp2f3(a1,a2,b1,b2,b3,1000000*j).ae('7.4859161049322370311e+610 - 2.8467477015940090189e+610j') 678 assert hyp2f3(a1,a2,b1,b2,b3,10**7*j).ae('-1.7477645579418800826e+1938 - 1.7606522995808116405e+1938j') 679 assert hyp2f3(a1,a2,b1,b2,b3,10**8*j).ae('-1.6932731942958401784e+6137 - 2.4521909113114629368e+6137j') 680 assert hyp2f3(a1,a2,b1,b2,b3,10**20*j).ae('-2.0988815677627225449e+6141851451 + 5.7708223542739208681e+6141851452j') 681 682def test_hyper_2f2(): 683 mp.dps = 15 684 assert hyper([1,2],[3,4],5) == hyp2f2(1,2,3,4,5) 685 a1,a2,b1,b2 = (3,10),4,(1,2),1./16 686 assert hyp2f2(a1,a2,b1,b2,10).ae(448225936.3377556696) 687 assert hyp2f2(a1,a2,b1,b2,10000).ae('1.2012553712966636711e+4358') 688 assert hyp2f2(a1,a2,b1,b2,-20000).ae(-0.04182343755661214626) 689 assert hyp2f2(a1,a2,b1,b2,10**20).ae('1.1148680024303263661e+43429448190325182840') 690 691def test_orthpoly(): 692 mp.dps = 15 693 assert jacobi(-4,2,3,0.7).ae(22800./4913) 694 assert jacobi(3,2,4,5.5) == 4133.125 695 assert jacobi(1.5,5/6.,4,0).ae(-1.0851951434075508417) 696 assert jacobi(-2, 1, 2, 4).ae(-0.16) 697 assert jacobi(2, -1, 2.5, 4).ae(34.59375) 698 #assert jacobi(2, -1, 2, 4) == 28.5 699 assert legendre(5, 7) == 129367 700 assert legendre(0.5,0).ae(0.53935260118837935667) 701 assert legendre(-1,-1) == 1 702 assert legendre(0,-1) == 1 703 assert legendre(0, 1) == 1 704 assert legendre(1, -1) == -1 705 assert legendre(7, 1) == 1 706 assert legendre(7, -1) == -1 707 assert legendre(8,1.5).ae(15457523./32768) 708 assert legendre(j,-j).ae(2.4448182735671431011 + 0.6928881737669934843j) 709 assert chebyu(5,1) == 6 710 assert chebyt(3,2) == 26 711 assert legendre(3.5,-1) == inf 712 assert legendre(4.5,-1) == -inf 713 assert legendre(3.5+1j,-1) == mpc(inf,inf) 714 assert legendre(4.5+1j,-1) == mpc(-inf,-inf) 715 assert laguerre(4, -2, 3).ae(-1.125) 716 assert laguerre(3, 1+j, 0.5).ae(0.2291666666666666667 + 2.5416666666666666667j) 717 718def test_hermite(): 719 mp.dps = 15 720 assert hermite(-2, 0).ae(0.5) 721 assert hermite(-1, 0).ae(0.88622692545275801365) 722 assert hermite(0, 0).ae(1) 723 assert hermite(1, 0) == 0 724 assert hermite(2, 0).ae(-2) 725 assert hermite(0, 2).ae(1) 726 assert hermite(1, 2).ae(4) 727 assert hermite(1, -2).ae(-4) 728 assert hermite(2, -2).ae(14) 729 assert hermite(0.5, 0).ae(0.69136733903629335053) 730 assert hermite(9, 0) == 0 731 assert hermite(4,4).ae(3340) 732 assert hermite(3,4).ae(464) 733 assert hermite(-4,4).ae(0.00018623860287512396181) 734 assert hermite(-3,4).ae(0.0016540169879668766270) 735 assert hermite(9, 2.5j).ae(13638725j) 736 assert hermite(9, -2.5j).ae(-13638725j) 737 assert hermite(9, 100).ae(511078883759363024000) 738 assert hermite(9, -100).ae(-511078883759363024000) 739 assert hermite(9, 100j).ae(512922083920643024000j) 740 assert hermite(9, -100j).ae(-512922083920643024000j) 741 assert hermite(-9.5, 2.5j).ae(-2.9004951258126778174e-6 + 1.7601372934039951100e-6j) 742 assert hermite(-9.5, -2.5j).ae(-2.9004951258126778174e-6 - 1.7601372934039951100e-6j) 743 assert hermite(-9.5, 100).ae(1.3776300722767084162e-22, abs_eps=0, rel_eps=eps) 744 assert hermite(-9.5, -100).ae('1.3106082028470671626e4355') 745 assert hermite(-9.5, 100j).ae(-9.7900218581864768430e-23 - 9.7900218581864768430e-23j, abs_eps=0, rel_eps=eps) 746 assert hermite(-9.5, -100j).ae(-9.7900218581864768430e-23 + 9.7900218581864768430e-23j, abs_eps=0, rel_eps=eps) 747 assert hermite(2+3j, -1-j).ae(851.3677063883687676 - 1496.4373467871007997j) 748 749def test_gegenbauer(): 750 mp.dps = 15 751 assert gegenbauer(1,2,3).ae(12) 752 assert gegenbauer(2,3,4).ae(381) 753 assert gegenbauer(0,0,0) == 0 754 assert gegenbauer(2,-1,3) == 0 755 assert gegenbauer(-7, 0.5, 3).ae(8989) 756 assert gegenbauer(1, -0.5, 3).ae(-3) 757 assert gegenbauer(1, -1.5, 3).ae(-9) 758 assert gegenbauer(1, -0.5, 3).ae(-3) 759 assert gegenbauer(-0.5, -0.5, 3).ae(-2.6383553159023906245) 760 assert gegenbauer(2+3j, 1-j, 3+4j).ae(14.880536623203696780 + 20.022029711598032898j) 761 #assert gegenbauer(-2, -0.5, 3).ae(-12) 762 763def test_legenp(): 764 mp.dps = 15 765 assert legenp(2,0,4) == legendre(2,4) 766 assert legenp(-2, -1, 0.5).ae(0.43301270189221932338) 767 assert legenp(-2, -1, 0.5, type=3).ae(0.43301270189221932338j) 768 assert legenp(-2, 1, 0.5).ae(-0.86602540378443864676) 769 assert legenp(2+j, 3+4j, -j).ae(134742.98773236786148 + 429782.72924463851745j) 770 assert legenp(2+j, 3+4j, -j, type=3).ae(802.59463394152268507 - 251.62481308942906447j) 771 assert legenp(2,4,3).ae(0) 772 assert legenp(2,4,3,type=3).ae(0) 773 assert legenp(2,1,0.5).ae(-1.2990381056766579701) 774 assert legenp(2,1,0.5,type=3).ae(1.2990381056766579701j) 775 assert legenp(3,2,3).ae(-360) 776 assert legenp(3,3,3).ae(240j*2**0.5) 777 assert legenp(3,4,3).ae(0) 778 assert legenp(0,0.5,2).ae(0.52503756790433198939 - 0.52503756790433198939j) 779 assert legenp(-1,-0.5,2).ae(0.60626116232846498110 + 0.60626116232846498110j) 780 assert legenp(-2,0.5,2).ae(1.5751127037129959682 - 1.5751127037129959682j) 781 assert legenp(-2,0.5,-0.5).ae(-0.85738275810499171286) 782 783def test_legenq(): 784 mp.dps = 15 785 f = legenq 786 # Evaluation at poles 787 assert isnan(f(3,2,1)) 788 assert isnan(f(3,2,-1)) 789 assert isnan(f(3,2,1,type=3)) 790 assert isnan(f(3,2,-1,type=3)) 791 # Evaluation at 0 792 assert f(0,1,0,type=2).ae(-1) 793 assert f(-2,2,0,type=2,zeroprec=200).ae(0) 794 assert f(1.5,3,0,type=2).ae(-2.2239343475841951023) 795 assert f(0,1,0,type=3).ae(j) 796 assert f(-2,2,0,type=3,zeroprec=200).ae(0) 797 assert f(1.5,3,0,type=3).ae(2.2239343475841951022*(1-1j)) 798 # Standard case, degree 0 799 assert f(0,0,-1.5).ae(-0.8047189562170501873 + 1.5707963267948966192j) 800 assert f(0,0,-0.5).ae(-0.54930614433405484570) 801 assert f(0,0,0,zeroprec=200).ae(0) 802 assert f(0,0,0.5).ae(0.54930614433405484570) 803 assert f(0,0,1.5).ae(0.8047189562170501873 - 1.5707963267948966192j) 804 assert f(0,0,-1.5,type=3).ae(-0.80471895621705018730) 805 assert f(0,0,-0.5,type=3).ae(-0.5493061443340548457 - 1.5707963267948966192j) 806 assert f(0,0,0,type=3).ae(-1.5707963267948966192j) 807 assert f(0,0,0.5,type=3).ae(0.5493061443340548457 - 1.5707963267948966192j) 808 assert f(0,0,1.5,type=3).ae(0.80471895621705018730) 809 # Standard case, degree 1 810 assert f(1,0,-1.5).ae(0.2070784343255752810 - 2.3561944901923449288j) 811 assert f(1,0,-0.5).ae(-0.72534692783297257715) 812 assert f(1,0,0).ae(-1) 813 assert f(1,0,0.5).ae(-0.72534692783297257715) 814 assert f(1,0,1.5).ae(0.2070784343255752810 - 2.3561944901923449288j) 815 # Standard case, degree 2 816 assert f(2,0,-1.5).ae(-0.0635669991240192885 + 4.5160394395353277803j) 817 assert f(2,0,-0.5).ae(0.81866326804175685571) 818 assert f(2,0,0,zeroprec=200).ae(0) 819 assert f(2,0,0.5).ae(-0.81866326804175685571) 820 assert f(2,0,1.5).ae(0.0635669991240192885 - 4.5160394395353277803j) 821 # Misc orders and degrees 822 assert f(2,3,1.5,type=2).ae(-5.7243340223994616228j) 823 assert f(2,3,1.5,type=3).ae(-5.7243340223994616228) 824 assert f(2,3,0.5,type=2).ae(-12.316805742712016310) 825 assert f(2,3,0.5,type=3).ae(-12.316805742712016310j) 826 assert f(2,3,-1.5,type=2).ae(-5.7243340223994616228j) 827 assert f(2,3,-1.5,type=3).ae(5.7243340223994616228) 828 assert f(2,3,-0.5,type=2).ae(-12.316805742712016310) 829 assert f(2,3,-0.5,type=3).ae(-12.316805742712016310j) 830 assert f(2+3j, 3+4j, 0.5, type=3).ae(0.0016119404873235186807 - 0.0005885900510718119836j) 831 assert f(2+3j, 3+4j, -1.5, type=3).ae(0.008451400254138808670 + 0.020645193304593235298j) 832 assert f(-2.5,1,-1.5).ae(3.9553395527435335749j) 833 assert f(-2.5,1,-0.5).ae(1.9290561746445456908) 834 assert f(-2.5,1,0).ae(1.2708196271909686299) 835 assert f(-2.5,1,0.5).ae(-0.31584812990742202869) 836 assert f(-2.5,1,1.5).ae(-3.9553395527435335742 + 0.2993235655044701706j) 837 assert f(-2.5,1,-1.5,type=3).ae(0.29932356550447017254j) 838 assert f(-2.5,1,-0.5,type=3).ae(-0.3158481299074220287 - 1.9290561746445456908j) 839 assert f(-2.5,1,0,type=3).ae(1.2708196271909686292 - 1.2708196271909686299j) 840 assert f(-2.5,1,0.5,type=3).ae(1.9290561746445456907 + 0.3158481299074220287j) 841 assert f(-2.5,1,1.5,type=3).ae(-0.29932356550447017254) 842 843def test_agm(): 844 mp.dps = 15 845 assert agm(0,0) == 0 846 assert agm(0,1) == 0 847 assert agm(1,1) == 1 848 assert agm(7,7) == 7 849 assert agm(j,j) == j 850 assert (1/agm(1,sqrt(2))).ae(0.834626841674073186) 851 assert agm(1,2).ae(1.4567910310469068692) 852 assert agm(1,3).ae(1.8636167832448965424) 853 assert agm(1,j).ae(0.599070117367796104+0.599070117367796104j) 854 assert agm(2) == agm(1,2) 855 assert agm(-3,4).ae(0.63468509766550907+1.3443087080896272j) 856 857def test_gammainc(): 858 mp.dps = 15 859 assert gammainc(2,5).ae(6*exp(-5)) 860 assert gammainc(2,0,5).ae(1-6*exp(-5)) 861 assert gammainc(2,3,5).ae(-6*exp(-5)+4*exp(-3)) 862 assert gammainc(-2.5,-0.5).ae(-0.9453087204829418812-5.3164237738936178621j) 863 assert gammainc(0,2,4).ae(0.045121158298212213088) 864 assert gammainc(0,3).ae(0.013048381094197037413) 865 assert gammainc(0,2+j,1-j).ae(0.00910653685850304839-0.22378752918074432574j) 866 assert gammainc(0,1-j).ae(0.00028162445198141833+0.17932453503935894015j) 867 assert gammainc(3,4,5,True).ae(0.11345128607046320253) 868 assert gammainc(3.5,0,inf).ae(gamma(3.5)) 869 assert gammainc(-150.5,500).ae('6.9825435345798951153e-627') 870 assert gammainc(-150.5,800).ae('4.6885137549474089431e-788') 871 assert gammainc(-3.5, -20.5).ae(0.27008820585226911 - 1310.31447140574997636j) 872 assert gammainc(-3.5, -200.5).ae(0.27008820585226911 - 5.3264597096208368435e76j) # XXX real part 873 assert gammainc(0,0,2) == inf 874 assert gammainc(1,b=1).ae(0.6321205588285576784) 875 assert gammainc(3,2,2) == 0 876 assert gammainc(2,3+j,3-j).ae(-0.28135485191849314194j) 877 assert gammainc(4+0j,1).ae(5.8860710587430771455) 878 # GH issue #301 879 assert gammainc(-1,-1).ae(-0.8231640121031084799 + 3.1415926535897932385j) 880 assert gammainc(-2,-1).ae(1.7707229202810768576 - 1.5707963267948966192j) 881 assert gammainc(-3,-1).ae(-1.4963349162467073643 + 0.5235987755982988731j) 882 assert gammainc(-4,-1).ae(1.05365418617643814992 - 0.13089969389957471827j) 883 # Regularized upper gamma 884 assert isnan(gammainc(0, 0, regularized=True)) 885 assert gammainc(-1, 0, regularized=True) == inf 886 assert gammainc(1, 0, regularized=True) == 1 887 assert gammainc(0, 5, regularized=True) == 0 888 assert gammainc(0, 2+3j, regularized=True) == 0 889 assert gammainc(0, 5000, regularized=True) == 0 890 assert gammainc(0, 10**30, regularized=True) == 0 891 assert gammainc(-1, 5, regularized=True) == 0 892 assert gammainc(-1, 5000, regularized=True) == 0 893 assert gammainc(-1, 10**30, regularized=True) == 0 894 assert gammainc(-1, -5, regularized=True) == 0 895 assert gammainc(-1, -5000, regularized=True) == 0 896 assert gammainc(-1, -10**30, regularized=True) == 0 897 assert gammainc(-1, 3+4j, regularized=True) == 0 898 assert gammainc(1, 5, regularized=True).ae(exp(-5)) 899 assert gammainc(1, 5000, regularized=True).ae(exp(-5000)) 900 assert gammainc(1, 10**30, regularized=True).ae(exp(-10**30)) 901 assert gammainc(1, 3+4j, regularized=True).ae(exp(-3-4j)) 902 assert gammainc(-1000000,2).ae('1.3669297209397347754e-301037', abs_eps=0, rel_eps=8*eps) 903 assert gammainc(-1000000,2,regularized=True) == 0 904 assert gammainc(-1000000,3+4j).ae('-1.322575609404222361e-698979 - 4.9274570591854533273e-698978j', abs_eps=0, rel_eps=8*eps) 905 assert gammainc(-1000000,3+4j,regularized=True) == 0 906 assert gammainc(2+3j, 4+5j, regularized=True).ae(0.085422013530993285774-0.052595379150390078503j) 907 assert gammainc(1000j, 1000j, regularized=True).ae(0.49702647628921131761 + 0.00297355675013575341j) 908 # Generalized 909 assert gammainc(3,4,2) == -gammainc(3,2,4) 910 assert gammainc(4, 2, 3).ae(1.2593494302978947396) 911 assert gammainc(4, 2, 3, regularized=True).ae(0.20989157171631578993) 912 assert gammainc(0, 2, 3).ae(0.035852129613864082155) 913 assert gammainc(0, 2, 3, regularized=True) == 0 914 assert gammainc(-1, 2, 3).ae(0.015219822548487616132) 915 assert gammainc(-1, 2, 3, regularized=True) == 0 916 assert gammainc(0, 2, 3).ae(0.035852129613864082155) 917 assert gammainc(0, 2, 3, regularized=True) == 0 918 # Should use upper gammas 919 assert gammainc(5, 10000, 12000).ae('1.1359381951461801687e-4327', abs_eps=0, rel_eps=8*eps) 920 # Should use lower gammas 921 assert gammainc(10000, 2, 3).ae('8.1244514125995785934e4765') 922 # GH issue 306 923 assert gammainc(3,-1-1j) == 0 924 assert gammainc(3,-1+1j) == 0 925 assert gammainc(2,-1) == 0 926 assert gammainc(2,-1+0j) == 0 927 assert gammainc(2+0j,-1) == 0 928 929def test_gammainc_expint_n(): 930 # These tests are intended to check all cases of the low-level code 931 # for upper gamma and expint with small integer index. 932 # Need to cover positive/negative arguments; small/large/huge arguments 933 # for both positive and negative indices, as well as indices 0 and 1 934 # which may be special-cased 935 mp.dps = 15 936 assert expint(-3,3.5).ae(0.021456366563296693987) 937 assert expint(-2,3.5).ae(0.014966633183073309405) 938 assert expint(-1,3.5).ae(0.011092916359219041088) 939 assert expint(0,3.5).ae(0.0086278238349481430685) 940 assert expint(1,3.5).ae(0.0069701398575483929193) 941 assert expint(2,3.5).ae(0.0058018939208991255223) 942 assert expint(3,3.5).ae(0.0049453773495857807058) 943 assert expint(-3,-3.5).ae(-4.6618170604073311319) 944 assert expint(-2,-3.5).ae(-5.5996974157555515963) 945 assert expint(-1,-3.5).ae(-6.7582555017739415818) 946 assert expint(0,-3.5).ae(-9.4615577024835182145) 947 assert expint(1,-3.5).ae(-13.925353995152335292 - 3.1415926535897932385j) 948 assert expint(2,-3.5).ae(-15.62328702434085977 - 10.995574287564276335j) 949 assert expint(3,-3.5).ae(-10.783026313250347722 - 19.242255003237483586j) 950 assert expint(-3,350).ae(2.8614825451252838069e-155, abs_eps=0, rel_eps=8*eps) 951 assert expint(-2,350).ae(2.8532837224504675901e-155, abs_eps=0, rel_eps=8*eps) 952 assert expint(-1,350).ae(2.8451316155828634555e-155, abs_eps=0, rel_eps=8*eps) 953 assert expint(0,350).ae(2.8370258275042797989e-155, abs_eps=0, rel_eps=8*eps) 954 assert expint(1,350).ae(2.8289659656701459404e-155, abs_eps=0, rel_eps=8*eps) 955 assert expint(2,350).ae(2.8209516419468505006e-155, abs_eps=0, rel_eps=8*eps) 956 assert expint(3,350).ae(2.8129824725501272171e-155, abs_eps=0, rel_eps=8*eps) 957 assert expint(-3,-350).ae(-2.8528796154044839443e+149) 958 assert expint(-2,-350).ae(-2.8610072121701264351e+149) 959 assert expint(-1,-350).ae(-2.8691813842677537647e+149) 960 assert expint(0,-350).ae(-2.8774025343659421709e+149) 961 u = expint(1,-350) 962 assert u.ae(-2.8856710698020863568e+149) 963 assert u.imag.ae(-3.1415926535897932385) 964 u = expint(2,-350) 965 assert u.ae(-2.8939874026504650534e+149) 966 assert u.imag.ae(-1099.5574287564276335) 967 u = expint(3,-350) 968 assert u.ae(-2.9023519497915044349e+149) 969 assert u.imag.ae(-192422.55003237483586) 970 assert expint(-3,350000000000000000000000).ae('2.1592908471792544286e-152003068666138139677919', abs_eps=0, rel_eps=8*eps) 971 assert expint(-2,350000000000000000000000).ae('2.1592908471792544286e-152003068666138139677919', abs_eps=0, rel_eps=8*eps) 972 assert expint(-1,350000000000000000000000).ae('2.1592908471792544286e-152003068666138139677919', abs_eps=0, rel_eps=8*eps) 973 assert expint(0,350000000000000000000000).ae('2.1592908471792544286e-152003068666138139677919', abs_eps=0, rel_eps=8*eps) 974 assert expint(1,350000000000000000000000).ae('2.1592908471792544286e-152003068666138139677919', abs_eps=0, rel_eps=8*eps) 975 assert expint(2,350000000000000000000000).ae('2.1592908471792544286e-152003068666138139677919', abs_eps=0, rel_eps=8*eps) 976 assert expint(3,350000000000000000000000).ae('2.1592908471792544286e-152003068666138139677919', abs_eps=0, rel_eps=8*eps) 977 assert expint(-3,-350000000000000000000000).ae('-3.7805306852415755699e+152003068666138139677871') 978 assert expint(-2,-350000000000000000000000).ae('-3.7805306852415755699e+152003068666138139677871') 979 assert expint(-1,-350000000000000000000000).ae('-3.7805306852415755699e+152003068666138139677871') 980 assert expint(0,-350000000000000000000000).ae('-3.7805306852415755699e+152003068666138139677871') 981 u = expint(1,-350000000000000000000000) 982 assert u.ae('-3.7805306852415755699e+152003068666138139677871') 983 assert u.imag.ae(-3.1415926535897932385) 984 u = expint(2,-350000000000000000000000) 985 assert u.imag.ae(-1.0995574287564276335e+24) 986 assert u.ae('-3.7805306852415755699e+152003068666138139677871') 987 u = expint(3,-350000000000000000000000) 988 assert u.imag.ae(-1.9242255003237483586e+47) 989 assert u.ae('-3.7805306852415755699e+152003068666138139677871') 990 # Small case; no branch cut 991 assert gammainc(-3,3.5).ae(0.00010020262545203707109) 992 assert gammainc(-2,3.5).ae(0.00040370427343557393517) 993 assert gammainc(-1,3.5).ae(0.0016576839773997501492) 994 assert gammainc(0,3.5).ae(0.0069701398575483929193) 995 assert gammainc(1,3.5).ae(0.03019738342231850074) 996 assert gammainc(2,3.5).ae(0.13588822540043325333) 997 assert gammainc(3,3.5).ae(0.64169439772426814072) 998 # Small case; with branch cut 999 assert gammainc(-3,-3.5).ae(0.03595832954467563286 + 0.52359877559829887308j) 1000 assert gammainc(-2,-3.5).ae(-0.88024704597962022221 - 1.5707963267948966192j) 1001 assert gammainc(-1,-3.5).ae(4.4637962926688170771 + 3.1415926535897932385j) 1002 assert gammainc(0,-3.5).ae(-13.925353995152335292 - 3.1415926535897932385j) 1003 assert gammainc(1,-3.5).ae(33.115451958692313751) 1004 assert gammainc(2,-3.5).ae(-82.788629896730784377) 1005 assert gammainc(3,-3.5).ae(240.08702670051927469) 1006 # Asymptotic case; no branch cut 1007 assert gammainc(-3,350).ae(6.5424095113340358813e-163, abs_eps=0, rel_eps=8*eps) 1008 assert gammainc(-2,350).ae(2.296312222489899769e-160, abs_eps=0, rel_eps=8*eps) 1009 assert gammainc(-1,350).ae(8.059861834133858573e-158, abs_eps=0, rel_eps=8*eps) 1010 assert gammainc(0,350).ae(2.8289659656701459404e-155, abs_eps=0, rel_eps=8*eps) 1011 assert gammainc(1,350).ae(9.9295903962649792963e-153, abs_eps=0, rel_eps=8*eps) 1012 assert gammainc(2,350).ae(3.485286229089007733e-150, abs_eps=0, rel_eps=8*eps) 1013 assert gammainc(3,350).ae(1.2233453960006379793e-147, abs_eps=0, rel_eps=8*eps) 1014 # Asymptotic case; branch cut 1015 u = gammainc(-3,-350) 1016 assert u.ae(6.7889565783842895085e+141) 1017 assert u.imag.ae(0.52359877559829887308) 1018 u = gammainc(-2,-350) 1019 assert u.ae(-2.3692668977889832121e+144) 1020 assert u.imag.ae(-1.5707963267948966192) 1021 u = gammainc(-1,-350) 1022 assert u.ae(8.2685354361441858669e+146) 1023 assert u.imag.ae(3.1415926535897932385) 1024 u = gammainc(0,-350) 1025 assert u.ae(-2.8856710698020863568e+149) 1026 assert u.imag.ae(-3.1415926535897932385) 1027 u = gammainc(1,-350) 1028 assert u.ae(1.0070908870280797598e+152) 1029 assert u.imag == 0 1030 u = gammainc(2,-350) 1031 assert u.ae(-3.5147471957279983618e+154) 1032 assert u.imag == 0 1033 u = gammainc(3,-350) 1034 assert u.ae(1.2266568422179417091e+157) 1035 assert u.imag == 0 1036 # Extreme asymptotic case 1037 assert gammainc(-3,350000000000000000000000).ae('5.0362468738874738859e-152003068666138139677990', abs_eps=0, rel_eps=8*eps) 1038 assert gammainc(-2,350000000000000000000000).ae('1.7626864058606158601e-152003068666138139677966', abs_eps=0, rel_eps=8*eps) 1039 assert gammainc(-1,350000000000000000000000).ae('6.1694024205121555102e-152003068666138139677943', abs_eps=0, rel_eps=8*eps) 1040 assert gammainc(0,350000000000000000000000).ae('2.1592908471792544286e-152003068666138139677919', abs_eps=0, rel_eps=8*eps) 1041 assert gammainc(1,350000000000000000000000).ae('7.5575179651273905e-152003068666138139677896', abs_eps=0, rel_eps=8*eps) 1042 assert gammainc(2,350000000000000000000000).ae('2.645131287794586675e-152003068666138139677872', abs_eps=0, rel_eps=8*eps) 1043 assert gammainc(3,350000000000000000000000).ae('9.2579595072810533625e-152003068666138139677849', abs_eps=0, rel_eps=8*eps) 1044 u = gammainc(-3,-350000000000000000000000) 1045 assert u.ae('8.8175642804468234866e+152003068666138139677800') 1046 assert u.imag.ae(0.52359877559829887308) 1047 u = gammainc(-2,-350000000000000000000000) 1048 assert u.ae('-3.0861474981563882203e+152003068666138139677824') 1049 assert u.imag.ae(-1.5707963267948966192) 1050 u = gammainc(-1,-350000000000000000000000) 1051 assert u.ae('1.0801516243547358771e+152003068666138139677848') 1052 assert u.imag.ae(3.1415926535897932385) 1053 u = gammainc(0,-350000000000000000000000) 1054 assert u.ae('-3.7805306852415755699e+152003068666138139677871') 1055 assert u.imag.ae(-3.1415926535897932385) 1056 assert gammainc(1,-350000000000000000000000).ae('1.3231857398345514495e+152003068666138139677895') 1057 assert gammainc(2,-350000000000000000000000).ae('-4.6311500894209300731e+152003068666138139677918') 1058 assert gammainc(3,-350000000000000000000000).ae('1.6209025312973255256e+152003068666138139677942') 1059 1060def test_incomplete_beta(): 1061 mp.dps = 15 1062 assert betainc(-2,-3,0.5,0.75).ae(63.4305673311255413583969) 1063 assert betainc(4.5,0.5+2j,2.5,6).ae(0.2628801146130621387903065 + 0.5162565234467020592855378j) 1064 assert betainc(4,5,0,6).ae(90747.77142857142857142857) 1065 1066def test_erf(): 1067 mp.dps = 15 1068 assert erf(0) == 0 1069 assert erf(1).ae(0.84270079294971486934) 1070 assert erf(3+4j).ae(-120.186991395079444098 - 27.750337293623902498j) 1071 assert erf(-4-3j).ae(-0.99991066178539168236 + 0.00004972026054496604j) 1072 assert erf(pi).ae(0.99999112385363235839) 1073 assert erf(1j).ae(1.6504257587975428760j) 1074 assert erf(-1j).ae(-1.6504257587975428760j) 1075 assert isinstance(erf(1), mpf) 1076 assert isinstance(erf(-1), mpf) 1077 assert isinstance(erf(0), mpf) 1078 assert isinstance(erf(0j), mpc) 1079 assert erf(inf) == 1 1080 assert erf(-inf) == -1 1081 assert erfi(0) == 0 1082 assert erfi(1/pi).ae(0.371682698493894314) 1083 assert erfi(inf) == inf 1084 assert erfi(-inf) == -inf 1085 assert erf(1+0j) == erf(1) 1086 assert erfc(1+0j) == erfc(1) 1087 assert erf(0.2+0.5j).ae(1 - erfc(0.2+0.5j)) 1088 assert erfc(0) == 1 1089 assert erfc(1).ae(1-erf(1)) 1090 assert erfc(-1).ae(1-erf(-1)) 1091 assert erfc(1/pi).ae(1-erf(1/pi)) 1092 assert erfc(-10) == 2 1093 assert erfc(-1000000) == 2 1094 assert erfc(-inf) == 2 1095 assert erfc(inf) == 0 1096 assert isnan(erfc(nan)) 1097 assert (erfc(10**4)*mpf(10)**43429453).ae('3.63998738656420') 1098 assert erf(8+9j).ae(-1072004.2525062051158 + 364149.91954310255423j) 1099 assert erfc(8+9j).ae(1072005.2525062051158 - 364149.91954310255423j) 1100 assert erfc(-8-9j).ae(-1072003.2525062051158 + 364149.91954310255423j) 1101 mp.dps = 50 1102 # This one does not use the asymptotic series 1103 assert (erfc(10)*10**45).ae('2.0884875837625447570007862949577886115608181193212') 1104 # This one does 1105 assert (erfc(50)*10**1088).ae('2.0709207788416560484484478751657887929322509209954') 1106 mp.dps = 15 1107 assert str(erfc(10**50)) == '3.66744826532555e-4342944819032518276511289189166050822943970058036665661144537831658646492088707747292249493384317534' 1108 assert erfinv(0) == 0 1109 assert erfinv(0.5).ae(0.47693627620446987338) 1110 assert erfinv(-0.5).ae(-0.47693627620446987338) 1111 assert erfinv(1) == inf 1112 assert erfinv(-1) == -inf 1113 assert erf(erfinv(0.95)).ae(0.95) 1114 assert erf(erfinv(0.999999999995)).ae(0.999999999995) 1115 assert erf(erfinv(-0.999999999995)).ae(-0.999999999995) 1116 mp.dps = 50 1117 assert erf(erfinv('0.99999999999999999999999999999995')).ae('0.99999999999999999999999999999995') 1118 assert erf(erfinv('0.999999999999999999999999999999995')).ae('0.999999999999999999999999999999995') 1119 assert erf(erfinv('-0.999999999999999999999999999999995')).ae('-0.999999999999999999999999999999995') 1120 mp.dps = 15 1121 # Complex asymptotic expansions 1122 v = erfc(50j) 1123 assert v.real == 1 1124 assert v.imag.ae('-6.1481820666053078736e+1083') 1125 assert erfc(-100+5j).ae(2) 1126 assert (erfc(100+5j)*10**4335).ae(2.3973567853824133572 - 3.9339259530609420597j) 1127 assert erfc(100+100j).ae(0.00065234366376857698698 - 0.0039357263629214118437j) 1128 1129def test_pdf(): 1130 mp.dps = 15 1131 assert npdf(-inf) == 0 1132 assert npdf(inf) == 0 1133 assert npdf(5,0,2).ae(npdf(5+4,4,2)) 1134 assert quadts(lambda x: npdf(x,-0.5,0.8), [-inf, inf]) == 1 1135 assert ncdf(0) == 0.5 1136 assert ncdf(3,3) == 0.5 1137 assert ncdf(-inf) == 0 1138 assert ncdf(inf) == 1 1139 assert ncdf(10) == 1 1140 # Verify that this is computed accurately 1141 assert (ncdf(-10)*10**24).ae(7.619853024160526) 1142 1143def test_lambertw(): 1144 mp.dps = 15 1145 assert lambertw(0) == 0 1146 assert lambertw(0+0j) == 0 1147 assert lambertw(inf) == inf 1148 assert isnan(lambertw(nan)) 1149 assert lambertw(inf,1).real == inf 1150 assert lambertw(inf,1).imag.ae(2*pi) 1151 assert lambertw(-inf,1).real == inf 1152 assert lambertw(-inf,1).imag.ae(3*pi) 1153 assert lambertw(0,-1) == -inf 1154 assert lambertw(0,1) == -inf 1155 assert lambertw(0,3) == -inf 1156 assert lambertw(e).ae(1) 1157 assert lambertw(1).ae(0.567143290409783873) 1158 assert lambertw(-pi/2).ae(j*pi/2) 1159 assert lambertw(-log(2)/2).ae(-log(2)) 1160 assert lambertw(0.25).ae(0.203888354702240164) 1161 assert lambertw(-0.25).ae(-0.357402956181388903) 1162 assert lambertw(-1./10000,0).ae(-0.000100010001500266719) 1163 assert lambertw(-0.25,-1).ae(-2.15329236411034965) 1164 assert lambertw(0.25,-1).ae(-3.00899800997004620-4.07652978899159763j) 1165 assert lambertw(-0.25,-1).ae(-2.15329236411034965) 1166 assert lambertw(0.25,1).ae(-3.00899800997004620+4.07652978899159763j) 1167 assert lambertw(-0.25,1).ae(-3.48973228422959210+7.41405453009603664j) 1168 assert lambertw(-4).ae(0.67881197132094523+1.91195078174339937j) 1169 assert lambertw(-4,1).ae(-0.66743107129800988+7.76827456802783084j) 1170 assert lambertw(-4,-1).ae(0.67881197132094523-1.91195078174339937j) 1171 assert lambertw(1000).ae(5.24960285240159623) 1172 assert lambertw(1000,1).ae(4.91492239981054535+5.44652615979447070j) 1173 assert lambertw(1000,-1).ae(4.91492239981054535-5.44652615979447070j) 1174 assert lambertw(1000,5).ae(3.5010625305312892+29.9614548941181328j) 1175 assert lambertw(3+4j).ae(1.281561806123775878+0.533095222020971071j) 1176 assert lambertw(-0.4+0.4j).ae(-0.10396515323290657+0.61899273315171632j) 1177 assert lambertw(3+4j,1).ae(-0.11691092896595324+5.61888039871282334j) 1178 assert lambertw(3+4j,-1).ae(0.25856740686699742-3.85211668616143559j) 1179 assert lambertw(-0.5,-1).ae(-0.794023632344689368-0.770111750510379110j) 1180 assert lambertw(-1./10000,1).ae(-11.82350837248724344+6.80546081842002101j) 1181 assert lambertw(-1./10000,-1).ae(-11.6671145325663544) 1182 assert lambertw(-1./10000,-2).ae(-11.82350837248724344-6.80546081842002101j) 1183 assert lambertw(-1./100000,4).ae(-14.9186890769540539+26.1856750178782046j) 1184 assert lambertw(-1./100000,5).ae(-15.0931437726379218666+32.5525721210262290086j) 1185 assert lambertw((2+j)/10).ae(0.173704503762911669+0.071781336752835511j) 1186 assert lambertw((2+j)/10,1).ae(-3.21746028349820063+4.56175438896292539j) 1187 assert lambertw((2+j)/10,-1).ae(-3.03781405002993088-3.53946629633505737j) 1188 assert lambertw((2+j)/10,4).ae(-4.6878509692773249+23.8313630697683291j) 1189 assert lambertw(-(2+j)/10).ae(-0.226933772515757933-0.164986470020154580j) 1190 assert lambertw(-(2+j)/10,1).ae(-2.43569517046110001+0.76974067544756289j) 1191 assert lambertw(-(2+j)/10,-1).ae(-3.54858738151989450-6.91627921869943589j) 1192 assert lambertw(-(2+j)/10,4).ae(-4.5500846928118151+20.6672982215434637j) 1193 mp.dps = 50 1194 assert lambertw(pi).ae('1.073658194796149172092178407024821347547745350410314531') 1195 mp.dps = 15 1196 # Former bug in generated branch 1197 assert lambertw(-0.5+0.002j).ae(-0.78917138132659918344 + 0.76743539379990327749j) 1198 assert lambertw(-0.5-0.002j).ae(-0.78917138132659918344 - 0.76743539379990327749j) 1199 assert lambertw(-0.448+0.4j).ae(-0.11855133765652382241 + 0.66570534313583423116j) 1200 assert lambertw(-0.448-0.4j).ae(-0.11855133765652382241 - 0.66570534313583423116j) 1201 assert lambertw(-0.65475+0.0001j).ae(-0.61053421111385310898+1.0396534993944097723803j) 1202 # Huge branch index 1203 w = lambertw(1,10**20) 1204 assert w.real.ae(-47.889578926290259164) 1205 assert w.imag.ae(6.2831853071795864769e+20) 1206 1207def test_lambertw_hard(): 1208 def check(x,y): 1209 y = convert(y) 1210 type_ok = True 1211 if isinstance(y, mpf): 1212 type_ok = isinstance(x, mpf) 1213 real_ok = abs(x.real-y.real) <= abs(y.real)*8*eps 1214 imag_ok = abs(x.imag-y.imag) <= abs(y.imag)*8*eps 1215 #print x, y, abs(x.real-y.real), abs(x.imag-y.imag) 1216 return real_ok and imag_ok 1217 # Evaluation near 0 1218 mp.dps = 15 1219 assert check(lambertw(1e-10), 9.999999999000000000e-11) 1220 assert check(lambertw(-1e-10), -1.000000000100000000e-10) 1221 assert check(lambertw(1e-10j), 9.999999999999999999733e-21 + 9.99999999999999999985e-11j) 1222 assert check(lambertw(-1e-10j), 9.999999999999999999733e-21 - 9.99999999999999999985e-11j) 1223 assert check(lambertw(1e-10,1), -26.303186778379041559 + 3.265093911703828397j) 1224 assert check(lambertw(-1e-10,1), -26.326236166739163892 + 6.526183280686333315j) 1225 assert check(lambertw(1e-10j,1), -26.312931726911421551 + 4.896366881798013421j) 1226 assert check(lambertw(-1e-10j,1), -26.297238779529035066 + 1.632807161345576513j) 1227 assert check(lambertw(1e-10,-1), -26.303186778379041559 - 3.265093911703828397j) 1228 assert check(lambertw(-1e-10,-1), -26.295238819246925694) 1229 assert check(lambertw(1e-10j,-1), -26.297238779529035028 - 1.6328071613455765135j) 1230 assert check(lambertw(-1e-10j,-1), -26.312931726911421551 - 4.896366881798013421j) 1231 # Test evaluation very close to the branch point -1/e 1232 # on the -1, 0, and 1 branches 1233 add = lambda x, y: fadd(x,y,exact=True) 1234 sub = lambda x, y: fsub(x,y,exact=True) 1235 addj = lambda x, y: fadd(x,fmul(y,1j,exact=True),exact=True) 1236 subj = lambda x, y: fadd(x,fmul(y,-1j,exact=True),exact=True) 1237 mp.dps = 1500 1238 a = -1/e + 10*eps 1239 d3 = mpf('1e-3') 1240 d10 = mpf('1e-10') 1241 d20 = mpf('1e-20') 1242 d40 = mpf('1e-40') 1243 d80 = mpf('1e-80') 1244 d300 = mpf('1e-300') 1245 d1000 = mpf('1e-1000') 1246 mp.dps = 15 1247 # ---- Branch 0 ---- 1248 # -1/e + eps 1249 assert check(lambertw(add(a,d3)), -0.92802015005456704876) 1250 assert check(lambertw(add(a,d10)), -0.99997668374140088071) 1251 assert check(lambertw(add(a,d20)), -0.99999999976683560186) 1252 assert lambertw(add(a,d40)) == -1 1253 assert lambertw(add(a,d80)) == -1 1254 assert lambertw(add(a,d300)) == -1 1255 assert lambertw(add(a,d1000)) == -1 1256 # -1/e - eps 1257 assert check(lambertw(sub(a,d3)), -0.99819016149860989001+0.07367191188934638577j) 1258 assert check(lambertw(sub(a,d10)), -0.9999999998187812114595992+0.0000233164398140346109194j) 1259 assert check(lambertw(sub(a,d20)), -0.99999999999999999998187+2.331643981597124203344e-10j) 1260 assert check(lambertw(sub(a,d40)), -1.0+2.33164398159712420336e-20j) 1261 assert check(lambertw(sub(a,d80)), -1.0+2.33164398159712420336e-40j) 1262 assert check(lambertw(sub(a,d300)), -1.0+2.33164398159712420336e-150j) 1263 assert check(lambertw(sub(a,d1000)), mpc(-1,'2.33164398159712420336e-500')) 1264 # -1/e + eps*j 1265 assert check(lambertw(addj(a,d3)), -0.94790387486938526634+0.05036819639190132490j) 1266 assert check(lambertw(addj(a,d10)), -0.9999835127872943680999899+0.0000164870314895821225256j) 1267 assert check(lambertw(addj(a,d20)), -0.999999999835127872929987+1.64872127051890935830e-10j) 1268 assert check(lambertw(addj(a,d40)), -0.9999999999999999999835+1.6487212707001281468305e-20j) 1269 assert check(lambertw(addj(a,d80)), -1.0 + 1.64872127070012814684865e-40j) 1270 assert check(lambertw(addj(a,d300)), -1.0 + 1.64872127070012814684865e-150j) 1271 assert check(lambertw(addj(a,d1000)), mpc(-1.0,'1.64872127070012814684865e-500')) 1272 # -1/e - eps*j 1273 assert check(lambertw(subj(a,d3)), -0.94790387486938526634-0.05036819639190132490j) 1274 assert check(lambertw(subj(a,d10)), -0.9999835127872943680999899-0.0000164870314895821225256j) 1275 assert check(lambertw(subj(a,d20)), -0.999999999835127872929987-1.64872127051890935830e-10j) 1276 assert check(lambertw(subj(a,d40)), -0.9999999999999999999835-1.6487212707001281468305e-20j) 1277 assert check(lambertw(subj(a,d80)), -1.0 - 1.64872127070012814684865e-40j) 1278 assert check(lambertw(subj(a,d300)), -1.0 - 1.64872127070012814684865e-150j) 1279 assert check(lambertw(subj(a,d1000)), mpc(-1.0,'-1.64872127070012814684865e-500')) 1280 # ---- Branch 1 ---- 1281 assert check(lambertw(addj(a,d3),1), -3.088501303219933378005990 + 7.458676867597474813950098j) 1282 assert check(lambertw(addj(a,d80),1), -3.088843015613043855957087 + 7.461489285654254556906117j) 1283 assert check(lambertw(addj(a,d300),1), -3.088843015613043855957087 + 7.461489285654254556906117j) 1284 assert check(lambertw(addj(a,d1000),1), -3.088843015613043855957087 + 7.461489285654254556906117j) 1285 assert check(lambertw(subj(a,d3),1), -1.0520914180450129534365906 + 0.0539925638125450525673175j) 1286 assert check(lambertw(subj(a,d10),1), -1.0000164872127056318529390 + 0.000016487393927159250398333077j) 1287 assert check(lambertw(subj(a,d20),1), -1.0000000001648721270700128 + 1.64872127088134693542628e-10j) 1288 assert check(lambertw(subj(a,d40),1), -1.000000000000000000016487 + 1.64872127070012814686677e-20j) 1289 assert check(lambertw(subj(a,d80),1), -1.0 + 1.64872127070012814684865e-40j) 1290 assert check(lambertw(subj(a,d300),1), -1.0 + 1.64872127070012814684865e-150j) 1291 assert check(lambertw(subj(a,d1000),1), mpc(-1.0, '1.64872127070012814684865e-500')) 1292 # ---- Branch -1 ---- 1293 # -1/e + eps 1294 assert check(lambertw(add(a,d3),-1), -1.075608941186624989414945) 1295 assert check(lambertw(add(a,d10),-1), -1.000023316621036696460620) 1296 assert check(lambertw(add(a,d20),-1), -1.000000000233164398177834) 1297 assert lambertw(add(a,d40),-1) == -1 1298 assert lambertw(add(a,d80),-1) == -1 1299 assert lambertw(add(a,d300),-1) == -1 1300 assert lambertw(add(a,d1000),-1) == -1 1301 # -1/e - eps 1302 assert check(lambertw(sub(a,d3),-1), -0.99819016149860989001-0.07367191188934638577j) 1303 assert check(lambertw(sub(a,d10),-1), -0.9999999998187812114595992-0.0000233164398140346109194j) 1304 assert check(lambertw(sub(a,d20),-1), -0.99999999999999999998187-2.331643981597124203344e-10j) 1305 assert check(lambertw(sub(a,d40),-1), -1.0-2.33164398159712420336e-20j) 1306 assert check(lambertw(sub(a,d80),-1), -1.0-2.33164398159712420336e-40j) 1307 assert check(lambertw(sub(a,d300),-1), -1.0-2.33164398159712420336e-150j) 1308 assert check(lambertw(sub(a,d1000),-1), mpc(-1,'-2.33164398159712420336e-500')) 1309 # -1/e + eps*j 1310 assert check(lambertw(addj(a,d3),-1), -1.0520914180450129534365906 - 0.0539925638125450525673175j) 1311 assert check(lambertw(addj(a,d10),-1), -1.0000164872127056318529390 - 0.0000164873939271592503983j) 1312 assert check(lambertw(addj(a,d20),-1), -1.0000000001648721270700 - 1.64872127088134693542628e-10j) 1313 assert check(lambertw(addj(a,d40),-1), -1.00000000000000000001648 - 1.6487212707001281468667726e-20j) 1314 assert check(lambertw(addj(a,d80),-1), -1.0 - 1.64872127070012814684865e-40j) 1315 assert check(lambertw(addj(a,d300),-1), -1.0 - 1.64872127070012814684865e-150j) 1316 assert check(lambertw(addj(a,d1000),-1), mpc(-1.0,'-1.64872127070012814684865e-500')) 1317 # -1/e - eps*j 1318 assert check(lambertw(subj(a,d3),-1), -3.088501303219933378005990-7.458676867597474813950098j) 1319 assert check(lambertw(subj(a,d10),-1), -3.088843015579260686911033-7.461489285372968780020716j) 1320 assert check(lambertw(subj(a,d20),-1), -3.088843015613043855953708-7.461489285654254556877988j) 1321 assert check(lambertw(subj(a,d40),-1), -3.088843015613043855957087-7.461489285654254556906117j) 1322 assert check(lambertw(subj(a,d80),-1), -3.088843015613043855957087 - 7.461489285654254556906117j) 1323 assert check(lambertw(subj(a,d300),-1), -3.088843015613043855957087 - 7.461489285654254556906117j) 1324 assert check(lambertw(subj(a,d1000),-1), -3.088843015613043855957087 - 7.461489285654254556906117j) 1325 # One more case, testing higher precision 1326 mp.dps = 500 1327 x = -1/e + mpf('1e-13') 1328 ans = "-0.99999926266961377166355784455394913638782494543377383"\ 1329 "744978844374498153493943725364881490261187530235150668593869563"\ 1330 "168276697689459394902153960200361935311512317183678882" 1331 mp.dps = 15 1332 assert lambertw(x).ae(ans) 1333 mp.dps = 50 1334 assert lambertw(x).ae(ans) 1335 mp.dps = 150 1336 assert lambertw(x).ae(ans) 1337 1338def test_meijerg(): 1339 mp.dps = 15 1340 assert meijerg([[2,3],[1]],[[0.5,2],[3,4]], 2.5).ae(4.2181028074787439386) 1341 assert meijerg([[],[1+j]],[[1],[1]], 3+4j).ae(271.46290321152464592 - 703.03330399954820169j) 1342 assert meijerg([[0.25],[1]],[[0.5],[2]],0) == 0 1343 assert meijerg([[0],[]],[[0,0,'1/3','2/3'], []], '2/27').ae(2.2019391389653314120) 1344 # Verify 1/z series being used 1345 assert meijerg([[-3],[-0.5]], [[-1],[-2.5]], -0.5).ae(-1.338096165935754898687431) 1346 assert meijerg([[1-(-1)],[1-(-2.5)]], [[1-(-3)],[1-(-0.5)]], -2.0).ae(-1.338096165935754898687431) 1347 assert meijerg([[-3],[-0.5]], [[-1],[-2.5]], -1).ae(-(pi+4)/(4*pi)) 1348 a = 2.5 1349 b = 1.25 1350 for z in [mpf(0.25), mpf(2)]: 1351 x1 = hyp1f1(a,b,z) 1352 x2 = gamma(b)/gamma(a)*meijerg([[1-a],[]],[[0],[1-b]],-z) 1353 x3 = gamma(b)/gamma(a)*meijerg([[1-0],[1-(1-b)]],[[1-(1-a)],[]],-1/z) 1354 assert x1.ae(x2) 1355 assert x1.ae(x3) 1356 1357def test_appellf1(): 1358 mp.dps = 15 1359 assert appellf1(2,-2,1,1,2,3).ae(-1.75) 1360 assert appellf1(2,1,-2,1,2,3).ae(-8) 1361 assert appellf1(2,1,-2,1,0.5,0.25).ae(1.5) 1362 assert appellf1(-2,1,3,2,3,3).ae(19) 1363 assert appellf1(1,2,3,4,0.5,0.125).ae( 1.53843285792549786518) 1364 1365def test_coulomb(): 1366 # Note: most tests are doctests 1367 # Test for a bug: 1368 mp.dps = 15 1369 assert coulombg(mpc(-5,0),2,3).ae(20.087729487721430394) 1370 1371def test_hyper_param_accuracy(): 1372 mp.dps = 15 1373 As = [n+1e-10 for n in range(-5,-1)] 1374 Bs = [n+1e-10 for n in range(-12,-5)] 1375 assert hyper(As,Bs,10).ae(-381757055858.652671927) 1376 assert legenp(0.5, 100, 0.25).ae(-2.4124576567211311755e+144) 1377 assert (hyp1f1(1000,1,-100)*10**24).ae(5.2589445437370169113) 1378 assert (hyp2f1(10, -900, 10.5, 0.99)*10**24).ae(1.9185370579660768203) 1379 assert (hyp2f1(1000,1.5,-3.5,-1.5)*10**385).ae(-2.7367529051334000764) 1380 assert hyp2f1(-5, 10, 3, 0.5, zeroprec=500) == 0 1381 assert (hyp1f1(-10000, 1000, 100)*10**424).ae(-3.1046080515824859974) 1382 assert (hyp2f1(1000,1.5,-3.5,-0.75,maxterms=100000)*10**231).ae(-4.0534790813913998643) 1383 assert legenp(2, 3, 0.25) == 0 1384 pytest.raises(ValueError, lambda: hypercomb(lambda a: [([],[],[],[],[a],[-a],0.5)], [3])) 1385 assert hypercomb(lambda a: [([],[],[],[],[a],[-a],0.5)], [3], infprec=200) == inf 1386 assert meijerg([[],[]],[[0,0,0,0],[]],0.1).ae(1.5680822343832351418) 1387 assert (besselk(400,400)*10**94).ae(1.4387057277018550583) 1388 mp.dps = 5 1389 (hyp1f1(-5000.5, 1500, 100)*10**185).ae(8.5185229673381935522) 1390 (hyp1f1(-5000, 1500, 100)*10**185).ae(9.1501213424563944311) 1391 mp.dps = 15 1392 (hyp1f1(-5000.5, 1500, 100)*10**185).ae(8.5185229673381935522) 1393 (hyp1f1(-5000, 1500, 100)*10**185).ae(9.1501213424563944311) 1394 assert hyp0f1(fadd(-20,'1e-100',exact=True), 0.25).ae(1.85014429040102783e+49) 1395 assert hyp0f1((-20*10**100+1, 10**100), 0.25).ae(1.85014429040102783e+49) 1396 1397def test_hypercomb_zero_pow(): 1398 # check that 0^0 = 1 1399 assert hypercomb(lambda a: (([0],[a],[],[],[],[],0),), [0]) == 1 1400 assert meijerg([[-1.5],[]],[[0],[-0.75]],0).ae(1.4464090846320771425) 1401 1402def test_spherharm(): 1403 mp.dps = 15 1404 t = 0.5; r = 0.25 1405 assert spherharm(0,0,t,r).ae(0.28209479177387814347) 1406 assert spherharm(1,-1,t,r).ae(0.16048941205971996369 - 0.04097967481096344271j) 1407 assert spherharm(1,0,t,r).ae(0.42878904414183579379) 1408 assert spherharm(1,1,t,r).ae(-0.16048941205971996369 - 0.04097967481096344271j) 1409 assert spherharm(2,-2,t,r).ae(0.077915886919031181734 - 0.042565643022253962264j) 1410 assert spherharm(2,-1,t,r).ae(0.31493387233497459884 - 0.08041582001959297689j) 1411 assert spherharm(2,0,t,r).ae(0.41330596756220761898) 1412 assert spherharm(2,1,t,r).ae(-0.31493387233497459884 - 0.08041582001959297689j) 1413 assert spherharm(2,2,t,r).ae(0.077915886919031181734 + 0.042565643022253962264j) 1414 assert spherharm(3,-3,t,r).ae(0.033640236589690881646 - 0.031339125318637082197j) 1415 assert spherharm(3,-2,t,r).ae(0.18091018743101461963 - 0.09883168583167010241j) 1416 assert spherharm(3,-1,t,r).ae(0.42796713930907320351 - 0.10927795157064962317j) 1417 assert spherharm(3,0,t,r).ae(0.27861659336351639787) 1418 assert spherharm(3,1,t,r).ae(-0.42796713930907320351 - 0.10927795157064962317j) 1419 assert spherharm(3,2,t,r).ae(0.18091018743101461963 + 0.09883168583167010241j) 1420 assert spherharm(3,3,t,r).ae(-0.033640236589690881646 - 0.031339125318637082197j) 1421 assert spherharm(0,-1,t,r) == 0 1422 assert spherharm(0,-2,t,r) == 0 1423 assert spherharm(0,1,t,r) == 0 1424 assert spherharm(0,2,t,r) == 0 1425 assert spherharm(1,2,t,r) == 0 1426 assert spherharm(1,3,t,r) == 0 1427 assert spherharm(1,-2,t,r) == 0 1428 assert spherharm(1,-3,t,r) == 0 1429 assert spherharm(2,3,t,r) == 0 1430 assert spherharm(2,4,t,r) == 0 1431 assert spherharm(2,-3,t,r) == 0 1432 assert spherharm(2,-4,t,r) == 0 1433 assert spherharm(3,4.5,0.5,0.25).ae(-22.831053442240790148 + 10.910526059510013757j) 1434 assert spherharm(2+3j, 1-j, 1+j, 3+4j).ae(-2.6582752037810116935 - 1.0909214905642160211j) 1435 assert spherharm(-6,2.5,t,r).ae(0.39383644983851448178 + 0.28414687085358299021j) 1436 assert spherharm(-3.5, 3, 0.5, 0.25).ae(0.014516852987544698924 - 0.015582769591477628495j) 1437 assert spherharm(-3, 3, 0.5, 0.25) == 0 1438 assert spherharm(-6, 3, 0.5, 0.25).ae(-0.16544349818782275459 - 0.15412657723253924562j) 1439 assert spherharm(-6, 1.5, 0.5, 0.25).ae(0.032208193499767402477 + 0.012678000924063664921j) 1440 assert spherharm(3,0,0,1).ae(0.74635266518023078283) 1441 assert spherharm(3,-2,0,1) == 0 1442 assert spherharm(3,-2,1,1).ae(-0.16270707338254028971 - 0.35552144137546777097j) 1443 1444def test_qfunctions(): 1445 mp.dps = 15 1446 assert qp(2,3,100).ae('2.7291482267247332183e2391') 1447 1448def test_issue_239(): 1449 mp.prec = 150 1450 x = ldexp(2476979795053773,-52) 1451 assert betainc(206, 385, 0, 0.55, 1).ae('0.99999999999999999999996570910644857895771110649954') 1452 mp.dps = 15 1453 pytest.raises(ValueError, lambda: hyp2f1(-5,5,0.5,0.5)) 1454 1455# Extra stress testing for Bessel functions 1456# Reference zeros generated with the aid of scipy.special 1457# jn_zero, jnp_zero, yn_zero, ynp_zero 1458 1459V = 15 1460M = 15 1461 1462jn_small_zeros = \ 1463[[2.4048255576957728, 1464 5.5200781102863106, 1465 8.6537279129110122, 1466 11.791534439014282, 1467 14.930917708487786, 1468 18.071063967910923, 1469 21.211636629879259, 1470 24.352471530749303, 1471 27.493479132040255, 1472 30.634606468431975, 1473 33.775820213573569, 1474 36.917098353664044, 1475 40.058425764628239, 1476 43.19979171317673, 1477 46.341188371661814], 1478 [3.8317059702075123, 1479 7.0155866698156188, 1480 10.173468135062722, 1481 13.323691936314223, 1482 16.470630050877633, 1483 19.615858510468242, 1484 22.760084380592772, 1485 25.903672087618383, 1486 29.046828534916855, 1487 32.189679910974404, 1488 35.332307550083865, 1489 38.474766234771615, 1490 41.617094212814451, 1491 44.759318997652822, 1492 47.901460887185447], 1493 [5.1356223018406826, 1494 8.4172441403998649, 1495 11.619841172149059, 1496 14.795951782351261, 1497 17.959819494987826, 1498 21.116997053021846, 1499 24.270112313573103, 1500 27.420573549984557, 1501 30.569204495516397, 1502 33.7165195092227, 1503 36.86285651128381, 1504 40.008446733478192, 1505 43.153453778371463, 1506 46.297996677236919, 1507 49.442164110416873], 1508 [6.3801618959239835, 1509 9.7610231299816697, 1510 13.015200721698434, 1511 16.223466160318768, 1512 19.409415226435012, 1513 22.582729593104442, 1514 25.748166699294978, 1515 28.908350780921758, 1516 32.064852407097709, 1517 35.218670738610115, 1518 38.370472434756944, 1519 41.520719670406776, 1520 44.669743116617253, 1521 47.817785691533302, 1522 50.965029906205183], 1523 [7.5883424345038044, 1524 11.064709488501185, 1525 14.37253667161759, 1526 17.615966049804833, 1527 20.826932956962388, 1528 24.01901952477111, 1529 27.199087765981251, 1530 30.371007667117247, 1531 33.537137711819223, 1532 36.699001128744649, 1533 39.857627302180889, 1534 43.01373772335443, 1535 46.167853512924375, 1536 49.320360686390272, 1537 52.471551398458023], 1538 [8.771483815959954, 1539 12.338604197466944, 1540 15.700174079711671, 1541 18.980133875179921, 1542 22.217799896561268, 1543 25.430341154222704, 1544 28.626618307291138, 1545 31.811716724047763, 1546 34.988781294559295, 1547 38.159868561967132, 1548 41.326383254047406, 1549 44.489319123219673, 1550 47.649399806697054, 1551 50.80716520300633, 1552 53.963026558378149], 1553 [9.9361095242176849, 1554 13.589290170541217, 1555 17.003819667816014, 1556 20.320789213566506, 1557 23.58608443558139, 1558 26.820151983411405, 1559 30.033722386570469, 1560 33.233041762847123, 1561 36.422019668258457, 1562 39.603239416075404, 1563 42.778481613199507, 1564 45.949015998042603, 1565 49.11577372476426, 1566 52.279453903601052, 1567 55.440592068853149], 1568 [11.086370019245084, 1569 14.821268727013171, 1570 18.287582832481726, 1571 21.641541019848401, 1572 24.934927887673022, 1573 28.191188459483199, 1574 31.42279419226558, 1575 34.637089352069324, 1576 37.838717382853611, 1577 41.030773691585537, 1578 44.21540850526126, 1579 47.394165755570512, 1580 50.568184679795566, 1581 53.738325371963291, 1582 56.905249991978781], 1583 [12.225092264004655, 1584 16.037774190887709, 1585 19.554536430997055, 1586 22.94517313187462, 1587 26.266814641176644, 1588 29.54565967099855, 1589 32.795800037341462, 1590 36.025615063869571, 1591 39.240447995178135, 1592 42.443887743273558, 1593 45.638444182199141, 1594 48.825930381553857, 1595 52.007691456686903, 1596 55.184747939289049, 1597 58.357889025269694], 1598 [13.354300477435331, 1599 17.241220382489128, 1600 20.807047789264107, 1601 24.233885257750552, 1602 27.583748963573006, 1603 30.885378967696675, 1604 34.154377923855096, 1605 37.400099977156589, 1606 40.628553718964528, 1607 43.843801420337347, 1608 47.048700737654032, 1609 50.245326955305383, 1610 53.435227157042058, 1611 56.619580266508436, 1612 59.799301630960228], 1613 [14.475500686554541, 1614 18.433463666966583, 1615 22.046985364697802, 1616 25.509450554182826, 1617 28.887375063530457, 1618 32.211856199712731, 1619 35.499909205373851, 1620 38.761807017881651, 1621 42.004190236671805, 1622 45.231574103535045, 1623 48.447151387269394, 1624 51.653251668165858, 1625 54.851619075963349, 1626 58.043587928232478, 1627 61.230197977292681], 1628 [15.589847884455485, 1629 19.61596690396692, 1630 23.275853726263409, 1631 26.773322545509539, 1632 30.17906117878486, 1633 33.526364075588624, 1634 36.833571341894905, 1635 40.111823270954241, 1636 43.368360947521711, 1637 46.608132676274944, 1638 49.834653510396724, 1639 53.050498959135054, 1640 56.257604715114484, 1641 59.457456908388002, 1642 62.651217388202912], 1643 [16.698249933848246, 1644 20.789906360078443, 1645 24.494885043881354, 1646 28.026709949973129, 1647 31.45996003531804, 1648 34.829986990290238, 1649 38.156377504681354, 1650 41.451092307939681, 1651 44.721943543191147, 1652 47.974293531269048, 1653 51.211967004101068, 1654 54.437776928325074, 1655 57.653844811906946, 1656 60.8618046824805, 1657 64.062937824850136], 1658 [17.801435153282442, 1659 21.95624406783631, 1660 25.705103053924724, 1661 29.270630441874802, 1662 32.731053310978403, 1663 36.123657666448762, 1664 39.469206825243883, 1665 42.780439265447158, 1666 46.06571091157561, 1667 49.330780096443524, 1668 52.579769064383396, 1669 55.815719876305778, 1670 59.040934037249271, 1671 62.257189393731728, 1672 65.465883797232125], 1673 [18.899997953174024, 1674 23.115778347252756, 1675 26.907368976182104, 1676 30.505950163896036, 1677 33.993184984781542, 1678 37.408185128639695, 1679 40.772827853501868, 1680 44.100590565798301, 1681 47.400347780543231, 1682 50.678236946479898, 1683 53.93866620912693, 1684 57.184898598119301, 1685 60.419409852130297, 1686 63.644117508962281, 1687 66.860533012260103]] 1688 1689jnp_small_zeros = \ 1690[[0.0, 1691 3.8317059702075123, 1692 7.0155866698156188, 1693 10.173468135062722, 1694 13.323691936314223, 1695 16.470630050877633, 1696 19.615858510468242, 1697 22.760084380592772, 1698 25.903672087618383, 1699 29.046828534916855, 1700 32.189679910974404, 1701 35.332307550083865, 1702 38.474766234771615, 1703 41.617094212814451, 1704 44.759318997652822], 1705 [1.8411837813406593, 1706 5.3314427735250326, 1707 8.5363163663462858, 1708 11.706004902592064, 1709 14.863588633909033, 1710 18.015527862681804, 1711 21.16436985918879, 1712 24.311326857210776, 1713 27.457050571059246, 1714 30.601922972669094, 1715 33.746182898667383, 1716 36.889987409236811, 1717 40.033444053350675, 1718 43.176628965448822, 1719 46.319597561173912], 1720 [3.0542369282271403, 1721 6.7061331941584591, 1722 9.9694678230875958, 1723 13.170370856016123, 1724 16.347522318321783, 1725 19.512912782488205, 1726 22.671581772477426, 1727 25.826037141785263, 1728 28.977672772993679, 1729 32.127327020443474, 1730 35.275535050674691, 1731 38.422654817555906, 1732 41.568934936074314, 1733 44.714553532819734, 1734 47.859641607992093], 1735 [4.2011889412105285, 1736 8.0152365983759522, 1737 11.345924310743006, 1738 14.585848286167028, 1739 17.78874786606647, 1740 20.9724769365377, 1741 24.144897432909265, 1742 27.310057930204349, 1743 30.470268806290424, 1744 33.626949182796679, 1745 36.781020675464386, 1746 39.933108623659488, 1747 43.083652662375079, 1748 46.232971081836478, 1749 49.381300092370349], 1750 [5.3175531260839944, 1751 9.2823962852416123, 1752 12.681908442638891, 1753 15.964107037731551, 1754 19.196028800048905, 1755 22.401032267689004, 1756 25.589759681386733, 1757 28.767836217666503, 1758 31.938539340972783, 1759 35.103916677346764, 1760 38.265316987088158, 1761 41.423666498500732, 1762 44.579623137359257, 1763 47.733667523865744, 1764 50.886159153182682], 1765 [6.4156163757002403, 1766 10.519860873772308, 1767 13.9871886301403, 1768 17.312842487884625, 1769 20.575514521386888, 1770 23.803581476593863, 1771 27.01030789777772, 1772 30.20284907898166, 1773 33.385443901010121, 1774 36.560777686880356, 1775 39.730640230067416, 1776 42.896273163494417, 1777 46.058566273567043, 1778 49.218174614666636, 1779 52.375591529563596], 1780 [7.501266144684147, 1781 11.734935953042708, 1782 15.268181461097873, 1783 18.637443009666202, 1784 21.931715017802236, 1785 25.183925599499626, 1786 28.409776362510085, 1787 31.617875716105035, 1788 34.81339298429743, 1789 37.999640897715301, 1790 41.178849474321413, 1791 44.352579199070217, 1792 47.521956905768113, 1793 50.687817781723741, 1794 53.85079463676896], 1795 [8.5778364897140741, 1796 12.932386237089576, 1797 16.529365884366944, 1798 19.941853366527342, 1799 23.268052926457571, 1800 26.545032061823576, 1801 29.790748583196614, 1802 33.015178641375142, 1803 36.224380548787162, 1804 39.422274578939259, 1805 42.611522172286684, 1806 45.793999658055002, 1807 48.971070951900596, 1808 52.143752969301988, 1809 55.312820330403446], 1810 [9.6474216519972168, 1811 14.115518907894618, 1812 17.774012366915256, 1813 21.229062622853124, 1814 24.587197486317681, 1815 27.889269427955092, 1816 31.155326556188325, 1817 34.39662855427218, 1818 37.620078044197086, 1819 40.830178681822041, 1820 44.030010337966153, 1821 47.221758471887113, 1822 50.407020967034367, 1823 53.586995435398319, 1824 56.762598475105272], 1825 [10.711433970699945, 1826 15.28673766733295, 1827 19.004593537946053, 1828 22.501398726777283, 1829 25.891277276839136, 1830 29.218563499936081, 1831 32.505247352375523, 1832 35.763792928808799, 1833 39.001902811514218, 1834 42.224638430753279, 1835 45.435483097475542, 1836 48.636922645305525, 1837 51.830783925834728, 1838 55.01844255063594, 1839 58.200955824859509], 1840 [11.770876674955582, 1841 16.447852748486498, 1842 20.223031412681701, 1843 23.760715860327448, 1844 27.182021527190532, 1845 30.534504754007074, 1846 33.841965775135715, 1847 37.118000423665604, 1848 40.371068905333891, 1849 43.606764901379516, 1850 46.828959446564562, 1851 50.040428970943456, 1852 53.243223214220535, 1853 56.438892058982552, 1854 59.628631306921512], 1855 [12.826491228033465, 1856 17.600266557468326, 1857 21.430854238060294, 1858 25.008518704644261, 1859 28.460857279654847, 1860 31.838424458616998, 1861 35.166714427392629, 1862 38.460388720328256, 1863 41.728625562624312, 1864 44.977526250903469, 1865 48.211333836373288, 1866 51.433105171422278, 1867 54.645106240447105, 1868 57.849056857839799, 1869 61.046288512821078], 1870 [13.878843069697276, 1871 18.745090916814406, 1872 22.629300302835503, 1873 26.246047773946584, 1874 29.72897816891134, 1875 33.131449953571661, 1876 36.480548302231658, 1877 39.791940718940855, 1878 43.075486800191012, 1879 46.337772104541405, 1880 49.583396417633095, 1881 52.815686826850452, 1882 56.037118687012179, 1883 59.249577075517968, 1884 62.454525995970462], 1885 [14.928374492964716, 1886 19.88322436109951, 1887 23.81938909003628, 1888 27.474339750968247, 1889 30.987394331665278, 1890 34.414545662167183, 1891 37.784378506209499, 1892 41.113512376883377, 1893 44.412454519229281, 1894 47.688252845993366, 1895 50.945849245830813, 1896 54.188831071035124, 1897 57.419876154678179, 1898 60.641030026538746, 1899 63.853885828967512], 1900 [15.975438807484321, 1901 21.015404934568315, 1902 25.001971500138194, 1903 28.694271223110755, 1904 32.236969407878118, 1905 35.688544091185301, 1906 39.078998185245057, 1907 42.425854432866141, 1908 45.740236776624833, 1909 49.029635055514276, 1910 52.299319390331728, 1911 55.553127779547459, 1912 58.793933759028134, 1913 62.02393848337554, 1914 65.244860767043859]] 1915 1916yn_small_zeros = \ 1917[[0.89357696627916752, 1918 3.9576784193148579, 1919 7.0860510603017727, 1920 10.222345043496417, 1921 13.361097473872763, 1922 16.500922441528091, 1923 19.64130970088794, 1924 22.782028047291559, 1925 25.922957653180923, 1926 29.064030252728398, 1927 32.205204116493281, 1928 35.346452305214321, 1929 38.487756653081537, 1930 41.629104466213808, 1931 44.770486607221993], 1932 [2.197141326031017, 1933 5.4296810407941351, 1934 8.5960058683311689, 1935 11.749154830839881, 1936 14.897442128336725, 1937 18.043402276727856, 1938 21.188068934142213, 1939 24.331942571356912, 1940 27.475294980449224, 1941 30.618286491641115, 1942 33.761017796109326, 1943 36.90355531614295, 1944 40.045944640266876, 1945 43.188218097393211, 1946 46.330399250701687], 1947 [3.3842417671495935, 1948 6.7938075132682675, 1949 10.023477979360038, 1950 13.209986710206416, 1951 16.378966558947457, 1952 19.539039990286384, 1953 22.69395593890929, 1954 25.845613720902269, 1955 28.995080395650151, 1956 32.143002257627551, 1957 35.289793869635804, 1958 38.435733485446343, 1959 41.581014867297885, 1960 44.725777117640461, 1961 47.870122696676504], 1962 [4.5270246611496439, 1963 8.0975537628604907, 1964 11.396466739595867, 1965 14.623077742393873, 1966 17.81845523294552, 1967 20.997284754187761, 1968 24.166235758581828, 1969 27.328799850405162, 1970 30.486989604098659, 1971 33.642049384702463, 1972 36.794791029185579, 1973 39.945767226378749, 1974 43.095367507846703, 1975 46.2438744334407, 1976 49.391498015725107], 1977 [5.6451478942208959, 1978 9.3616206152445429, 1979 12.730144474090465, 1980 15.999627085382479, 1981 19.22442895931681, 1982 22.424810599698521, 1983 25.610267054939328, 1984 28.785893657666548, 1985 31.954686680031668, 1986 35.118529525584828, 1987 38.278668089521758, 1988 41.435960629910073, 1989 44.591018225353424, 1990 47.744288086361052, 1991 50.896105199722123], 1992 [6.7471838248710219, 1993 10.597176726782031, 1994 14.033804104911233, 1995 17.347086393228382, 1996 20.602899017175335, 1997 23.826536030287532, 1998 27.030134937138834, 1999 30.220335654231385, 2000 33.401105611047908, 2001 36.574972486670962, 2002 39.743627733020277, 2003 42.908248189569535, 2004 46.069679073215439, 2005 49.228543693445843, 2006 52.385312123112282], 2007 [7.8377378223268716, 2008 11.811037107609447, 2009 15.313615118517857, 2010 18.670704965906724, 2011 21.958290897126571, 2012 25.206207715021249, 2013 28.429037095235496, 2014 31.634879502950644, 2015 34.828638524084437, 2016 38.013473399691765, 2017 41.19151880917741, 2018 44.364272633271975, 2019 47.53281875312084, 2020 50.697961822183806, 2021 53.860312300118388], 2022 [8.919605734873789, 2023 13.007711435388313, 2024 16.573915129085334, 2025 19.974342312352426, 2026 23.293972585596648, 2027 26.5667563757203, 2028 29.809531451608321, 2029 33.031769327150685, 2030 36.239265816598239, 2031 39.435790312675323, 2032 42.623910919472727, 2033 45.805442883111651, 2034 48.981708325514764, 2035 52.153694518185572, 2036 55.322154420959698], 2037 [9.9946283820824834, 2038 14.190361295800141, 2039 17.817887841179873, 2040 21.26093227125945, 2041 24.612576377421522, 2042 27.910524883974868, 2043 31.173701563441602, 2044 34.412862242025045, 2045 37.634648706110989, 2046 40.843415321050884, 2047 44.04214994542435, 2048 47.232978012841169, 2049 50.417456447370186, 2050 53.596753874948731, 2051 56.771765754432457], 2052 [11.064090256031013, 2053 15.361301343575925, 2054 19.047949646361388, 2055 22.532765416313869, 2056 25.91620496332662, 2057 29.2394205079349, 2058 32.523270869465881, 2059 35.779715464475261, 2060 39.016196664616095, 2061 42.237627509803703, 2062 45.4474001519274, 2063 48.647941127433196, 2064 51.841036928216499, 2065 55.028034667184916, 2066 58.209970905250097], 2067 [12.128927704415439, 2068 16.522284394784426, 2069 20.265984501212254, 2070 23.791669719454272, 2071 27.206568881574774, 2072 30.555020011020762, 2073 33.859683872746356, 2074 37.133649760307504, 2075 40.385117593813002, 2076 43.619533085646856, 2077 46.840676630553575, 2078 50.051265851897857, 2079 53.253310556711732, 2080 56.448332488918971, 2081 59.637507005589829], 2082 [13.189846995683845, 2083 17.674674253171487, 2084 21.473493977824902, 2085 25.03913093040942, 2086 28.485081336558058, 2087 31.858644293774859, 2088 35.184165245422787, 2089 38.475796636190897, 2090 41.742455848758449, 2091 44.990096293791186, 2092 48.222870660068338, 2093 51.443777308699826, 2094 54.655042589416311, 2095 57.858358441436511, 2096 61.055036135780528], 2097 [14.247395665073945, 2098 18.819555894710682, 2099 22.671697117872794, 2100 26.276375544903892, 2101 29.752925495549038, 2102 33.151412708998983, 2103 36.497763772987645, 2104 39.807134090704376, 2105 43.089121522203808, 2106 46.350163579538652, 2107 49.594769786270069, 2108 52.82620892320143, 2109 56.046916910756961, 2110 59.258751140598783, 2111 62.463155567737854], 2112 [15.30200785858925, 2113 19.957808654258601, 2114 23.861599172945054, 2115 27.504429642227545, 2116 31.011103429019229, 2117 34.434283425782942, 2118 37.801385632318459, 2119 41.128514139788358, 2120 44.425913324440663, 2121 47.700482714581842, 2122 50.957073905278458, 2123 54.199216028087261, 2124 57.429547607017405, 2125 60.65008661807661, 2126 63.862406280068586], 2127 [16.354034360047551, 2128 21.090156519983806, 2129 25.044040298785627, 2130 28.724161640881914, 2131 32.260472459522644, 2132 35.708083982611664, 2133 39.095820003878235, 2134 42.440684315990936, 2135 45.75353669045622, 2136 49.041718113283529, 2137 52.310408280968073, 2138 55.56338698149062, 2139 58.803488508906895, 2140 62.032886550960831, 2141 65.253280088312461]] 2142 2143ynp_small_zeros = \ 2144[[2.197141326031017, 2145 5.4296810407941351, 2146 8.5960058683311689, 2147 11.749154830839881, 2148 14.897442128336725, 2149 18.043402276727856, 2150 21.188068934142213, 2151 24.331942571356912, 2152 27.475294980449224, 2153 30.618286491641115, 2154 33.761017796109326, 2155 36.90355531614295, 2156 40.045944640266876, 2157 43.188218097393211, 2158 46.330399250701687], 2159 [3.6830228565851777, 2160 6.9414999536541757, 2161 10.123404655436613, 2162 13.285758156782854, 2163 16.440058007293282, 2164 19.590241756629495, 2165 22.738034717396327, 2166 25.884314618788867, 2167 29.029575819372535, 2168 32.174118233366201, 2169 35.318134458192094, 2170 38.461753870997549, 2171 41.605066618873108, 2172 44.74813744908079, 2173 47.891014070791065], 2174 [5.0025829314460639, 2175 8.3507247014130795, 2176 11.574195465217647, 2177 14.760909306207676, 2178 17.931285939466855, 2179 21.092894504412739, 2180 24.249231678519058, 2181 27.402145837145258, 2182 30.552708880564553, 2183 33.70158627151572, 2184 36.849213419846257, 2185 39.995887376143356, 2186 43.141817835750686, 2187 46.287157097544201, 2188 49.432018469138281], 2189 [6.2536332084598136, 2190 9.6987879841487711, 2191 12.972409052292216, 2192 16.19044719506921, 2193 19.38238844973613, 2194 22.559791857764261, 2195 25.728213194724094, 2196 28.890678419054777, 2197 32.048984005266337, 2198 35.204266606440635, 2199 38.357281675961019, 2200 41.508551443818436, 2201 44.658448731963676, 2202 47.807246956681162, 2203 50.95515126455207], 2204 [7.4649217367571329, 2205 11.005169149809189, 2206 14.3317235192331, 2207 17.58443601710272, 2208 20.801062338411128, 2209 23.997004122902644, 2210 27.179886689853435, 2211 30.353960608554323, 2212 33.521797098666792, 2213 36.685048382072301, 2214 39.844826969405863, 2215 43.001910515625288, 2216 46.15685955107263, 2217 49.310088614282257, 2218 52.461911043685864], 2219 [8.6495562436971983, 2220 12.280868725807848, 2221 15.660799304540377, 2222 18.949739756016503, 2223 22.192841809428241, 2224 25.409072788867674, 2225 28.608039283077593, 2226 31.795195353138159, 2227 34.973890634255288, 2228 38.14630522169358, 2229 41.313923188794905, 2230 44.477791768537617, 2231 47.638672065035628, 2232 50.797131066967842, 2233 53.953600129601663], 2234 [9.8147970120105779, 2235 13.532811875789828, 2236 16.965526446046053, 2237 20.291285512443867, 2238 23.56186260680065, 2239 26.799499736027237, 2240 30.015665481543419, 2241 33.216968050039509, 2242 36.407516858984748, 2243 39.590015243560459, 2244 42.766320595957378, 2245 45.937754257017323, 2246 49.105283450953203, 2247 52.269633324547373, 2248 55.431358715604255], 2249 [10.965152105242974, 2250 14.765687379508912, 2251 18.250123150217555, 2252 21.612750053384621, 2253 24.911310600813573, 2254 28.171051927637585, 2255 31.40518108895689, 2256 34.621401012564177, 2257 37.824552065973114, 2258 41.017847386464902, 2259 44.203512240871601, 2260 47.3831408366063, 2261 50.557907466622796, 2262 53.728697478957026, 2263 56.896191727313342], 2264 [12.103641941939539, 2265 15.982840905145284, 2266 19.517731005559611, 2267 22.916962141504605, 2268 26.243700855690533, 2269 29.525960140695407, 2270 32.778568197561124, 2271 36.010261572392516, 2272 39.226578757802172, 2273 42.43122493258747, 2274 45.626783824134354, 2275 48.815117837929515, 2276 51.997606404328863, 2277 55.175294723956816, 2278 58.348990221754937], 2279 [13.232403808592215, 2280 17.186756572616758, 2281 20.770762917490496, 2282 24.206152448722253, 2283 27.561059462697153, 2284 30.866053571250639, 2285 34.137476603379774, 2286 37.385039772270268, 2287 40.614946085165892, 2288 43.831373184731238, 2289 47.037251786726299, 2290 50.234705848765229, 2291 53.425316228549359, 2292 56.610286079882087, 2293 59.790548623216652], 2294 [14.35301374369987, 2295 18.379337301642568, 2296 22.011118775283494, 2297 25.482116178696707, 2298 28.865046588695164, 2299 32.192853922166294, 2300 35.483296655830277, 2301 38.747005493021857, 2302 41.990815194320955, 2303 45.219355876831731, 2304 48.435892856078888, 2305 51.642803925173029, 2306 54.84186659475857, 2307 58.034439083840155, 2308 61.221578745109862], 2309 [15.466672066554263, 2310 19.562077985759503, 2311 23.240325531101082, 2312 26.746322986645901, 2313 30.157042415639891, 2314 33.507642948240263, 2315 36.817212798512775, 2316 40.097251300178642, 2317 43.355193847719752, 2318 46.596103410173672, 2319 49.823567279972794, 2320 53.040208868780832, 2321 56.247996968470062, 2322 59.448441365714251, 2323 62.642721301357187], 2324 [16.574317035530872, 2325 20.73617763753932, 2326 24.459631728238804, 2327 27.999993668839644, 2328 31.438208790267783, 2329 34.811512070805535, 2330 38.140243708611251, 2331 41.436725143893739, 2332 44.708963264433333, 2333 47.962435051891027, 2334 51.201037321915983, 2335 54.427630745992975, 2336 57.644369734615238, 2337 60.852911791989989, 2338 64.054555435720397], 2339 [17.676697936439624, 2340 21.9026148697762, 2341 25.670073356263225, 2342 29.244155124266438, 2343 32.709534477396028, 2344 36.105399554497548, 2345 39.453272918267025, 2346 42.766255701958017, 2347 46.052899215578358, 2348 49.319076602061401, 2349 52.568982147952547, 2350 55.805705507386287, 2351 59.031580956740466, 2352 62.248409689597653, 2353 65.457606670836759], 2354 [18.774423978290318, 2355 23.06220035979272, 2356 26.872520985976736, 2357 30.479680663499762, 2358 33.971869047372436, 2359 37.390118854896324, 2360 40.757072537673599, 2361 44.086572292170345, 2362 47.387688809191869, 2363 50.66667461073936, 2364 53.928009929563275, 2365 57.175005343085052, 2366 60.410169281219877, 2367 63.635442539153021, 2368 66.85235358587768]] 2369 2370@pytest.mark.slow 2371def test_bessel_zeros_extra(): 2372 mp.dps = 15 2373 for v in range(V): 2374 for m in range(1,M+1): 2375 print(v, m, "of", V, M) 2376 # Twice to test cache (if used) 2377 assert besseljzero(v,m).ae(jn_small_zeros[v][m-1]) 2378 assert besseljzero(v,m).ae(jn_small_zeros[v][m-1]) 2379 assert besseljzero(v,m,1).ae(jnp_small_zeros[v][m-1]) 2380 assert besseljzero(v,m,1).ae(jnp_small_zeros[v][m-1]) 2381 assert besselyzero(v,m).ae(yn_small_zeros[v][m-1]) 2382 assert besselyzero(v,m).ae(yn_small_zeros[v][m-1]) 2383 assert besselyzero(v,m,1).ae(ynp_small_zeros[v][m-1]) 2384 assert besselyzero(v,m,1).ae(ynp_small_zeros[v][m-1]) 2385