bclib 2.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246
  1. scale = 50
  2. define e(x) {
  3. auto a, b, c, d, e, g, w, y, t, r
  4. r = ibase
  5. ibase = A
  6. t = scale
  7. scale = t + .434*x + 1
  8. w = 0
  9. if(x<0) {
  10. x = -x
  11. w = 1
  12. }
  13. y = 0
  14. while(x>2) {
  15. x /= 2
  16. y++
  17. }
  18. a = 1
  19. b = 1
  20. c = b
  21. d = 1
  22. e = 1
  23. for(a=1; 1; a++) {
  24. b *= x
  25. c = c*a+b
  26. d *= a
  27. g = c/d
  28. if(g == e) {
  29. g = g/1
  30. while(y--) {
  31. g *= g
  32. }
  33. scale = t
  34. if(w==1) {
  35. ibase = r
  36. return 1/g
  37. }
  38. ibase = r
  39. return g/1
  40. }
  41. e = g
  42. }
  43. }
  44. define l(x) {
  45. auto a, b, c, d, e, f, g, u, s, t, r, z
  46. r = ibase
  47. ibase = A
  48. if(x <= 0) {
  49. z = 1-10^scale
  50. ibase = r
  51. return z
  52. }
  53. t = scale
  54. f = 1
  55. scale += scale(x) - length(x) + 1
  56. s = scale
  57. while(x > 2) {
  58. s += (length(x)-scale(x))/2 + 1
  59. if(s>0) {
  60. scale = s
  61. }
  62. x = sqrt(x)
  63. f *= 2
  64. }
  65. while(x < .5) {
  66. s += (length(x)-scale(x))/2 + 1
  67. if(s>0) {
  68. scale = s
  69. }
  70. x = sqrt(x)
  71. f *= 2
  72. }
  73. scale = t + length(f) - scale(f) + 1
  74. u = (x-1)/(x+1)
  75. scale += 1.1*length(t) - 1.1*scale(t)
  76. s = u*u
  77. b = 2*f
  78. c = b
  79. d = 1
  80. e = 1
  81. for(a=3; 1; a=a+2){
  82. b *= s
  83. c = c*a + d*b
  84. d *= a
  85. g = c/d
  86. if(g==e) {
  87. scale = t
  88. ibase = r
  89. return u*c/d
  90. }
  91. e = g
  92. }
  93. }
  94. define s(x) {
  95. auto a, b, c, s, t, y, p, n, i, r
  96. r = ibase
  97. ibase = A
  98. t = scale
  99. y = x/.7853
  100. s = t + length(y) - scale(y)
  101. if(s<t) {
  102. s = t
  103. }
  104. scale = s
  105. p = a(1)
  106. scale = 0
  107. if(x>=0) {
  108. n = (x/(2*p)+1)/2
  109. }
  110. if(x<0) {
  111. n = (x/(2*p)-1)/2
  112. }
  113. x -= 4*n*p
  114. if(n%2 != 0) {
  115. x = -x
  116. }
  117. scale = t + length(1.2*t) - scale(1.2*t)
  118. y = -x*x
  119. a = x
  120. b = 1
  121. s = x
  122. for(i=3; 1; i+=2) {
  123. a *= y
  124. b *= i*(i-1)
  125. c = a/b
  126. if(c==0){
  127. scale = t
  128. ibase = r
  129. return s/1
  130. }
  131. s += c
  132. }
  133. }
  134. define c(x) {
  135. auto t, r
  136. r = ibase
  137. ibase = A
  138. t = scale
  139. scale = scale+1
  140. x = s(x + 2*a(1))
  141. scale = t
  142. ibase = r
  143. return x/1
  144. }
  145. define a(x) {
  146. auto a, b, c, d, e, f, g, s, t, r, z
  147. r = ibase
  148. ibase = A
  149. if(x==0) {
  150. return 0
  151. }
  152. if(x==1) {
  153. z = .7853981633974483096156608458198757210492923498437764/1
  154. ibase = r
  155. if(scale<52) {
  156. return z
  157. }
  158. }
  159. t = scale
  160. f = 1
  161. while(x > .5) {
  162. scale++
  163. x = -(1 - sqrt(1.+x*x))/x
  164. f *= 2
  165. }
  166. while(x < -.5) {
  167. scale++
  168. x = -(1 - sqrt(1.+x*x))/x
  169. f *= 2
  170. }
  171. s = -x*x
  172. b = f
  173. c = f
  174. d = 1
  175. e = 1
  176. for(a=3; 1; a+=2) {
  177. b *= s
  178. c = c*a + d*b
  179. d *= a
  180. g = c/d
  181. if(g==e) {
  182. scale = t
  183. ibase = r
  184. return x*c/d
  185. }
  186. e = g
  187. }
  188. }
  189. define j(n,x) {
  190. auto a,b,c,d,e,g,i,s,k,t,r
  191. r = ibase
  192. ibase = A
  193. t = scale
  194. k = 1.36*x + 1.16*t - n
  195. k = length(k) - scale(k)
  196. if(k>0) {
  197. scale += k
  198. }
  199. s = -x*x/4
  200. if(n<0) {
  201. n = -n
  202. x = -x
  203. }
  204. a = 1
  205. c = 1
  206. for(i=1; i<=n; i++) {
  207. a *= x
  208. c *= 2*i
  209. }
  210. b = a
  211. d = 1
  212. e = 1
  213. for(i=1; 1; i++) {
  214. a *= s
  215. b = b*i*(n+i) + a
  216. c *= i*(n+i)
  217. g = b/c
  218. if(g==e) {
  219. scale = t
  220. ibase = r
  221. return g/1
  222. }
  223. e = g
  224. }
  225. }