sum.rlg 14 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729
  1. Tue Apr 15 00:35:15 2008 run on win32
  2. % Tests of the SUM package.
  3. % Author: Fujio Kako (kako@kako.math.sci.hiroshima-u.ac.jp)
  4. % 1) Summations.
  5. sum(n,n);
  6. n*(n + 1)
  7. -----------
  8. 2
  9. for i:=2:10 do write sum(n**i,n);
  10. 2
  11. n*(2*n + 3*n + 1)
  12. --------------------
  13. 6
  14. 2 2
  15. n *(n + 2*n + 1)
  16. -------------------
  17. 4
  18. 4 3 2
  19. n*(6*n + 15*n + 10*n - 1)
  20. ------------------------------
  21. 30
  22. 2 4 3 2
  23. n *(2*n + 6*n + 5*n - 1)
  24. -----------------------------
  25. 12
  26. 6 5 4 2
  27. n*(6*n + 21*n + 21*n - 7*n + 1)
  28. -------------------------------------
  29. 42
  30. 2 6 5 4 2
  31. n *(3*n + 12*n + 14*n - 7*n + 2)
  32. --------------------------------------
  33. 24
  34. 8 7 6 4 2
  35. n*(10*n + 45*n + 60*n - 42*n + 20*n - 3)
  36. -----------------------------------------------
  37. 90
  38. 2 8 7 6 4 2
  39. n *(2*n + 10*n + 15*n - 14*n + 10*n - 3)
  40. -----------------------------------------------
  41. 20
  42. 10 9 8 6 4 2
  43. n*(6*n + 33*n + 55*n - 66*n + 66*n - 33*n + 5)
  44. -------------------------------------------------------
  45. 66
  46. sum((n+1)**3,n);
  47. 3 2
  48. n*(n + 6*n + 13*n + 12)
  49. ---------------------------
  50. 4
  51. sum(x**n,n);
  52. n
  53. x *x
  54. -------
  55. x - 1
  56. sum(n**2*x**n,n);
  57. n 2 2 2 2
  58. x *x*(n *x - 2*n *x + n - 2*n*x + 2*n + x + 1)
  59. --------------------------------------------------
  60. 3 2
  61. x - 3*x + 3*x - 1
  62. sum(1/n,n);
  63. 1
  64. sum(---,n)
  65. n
  66. sum(1/n/(n+2),n);
  67. n*(3*n + 5)
  68. ------------------
  69. 2
  70. 4*(n + 3*n + 2)
  71. sum(log (n/(n+1)),n);
  72. 1
  73. log(-------)
  74. n + 1
  75. % 2) Expressions including trigonometric functions.
  76. sum(sin(n*x),n);
  77. 2*n*x + x
  78. - cos(-----------)
  79. 2
  80. ---------------------
  81. x
  82. 2*sin(---)
  83. 2
  84. sum(n*sin(n*x),n,1,k);
  85. sin(k*x + x)*k - sin(k*x)*k - sin(k*x)
  86. ----------------------------------------
  87. 2*(cos(x) - 1)
  88. sum(cos((2*r-1)*pi/n),r);
  89. 2*pi*r
  90. sin(--------)
  91. n
  92. ---------------
  93. pi
  94. 2*sin(----)
  95. n
  96. sum(cos((2*r-1)*pi/n),r,1,n);
  97. 0
  98. sum(cos((2*r-1)*pi/(2*n+1)),r);
  99. 2*pi*r
  100. sin(---------)
  101. 2*n + 1
  102. ------------------
  103. pi
  104. 2*sin(---------)
  105. 2*n + 1
  106. sum(cos((2*r-1)*pi/(2*n+1)),r,1,n);
  107. 2*n*pi
  108. sin(---------)
  109. 2*n + 1
  110. ------------------
  111. pi
  112. 2*sin(---------)
  113. 2*n + 1
  114. sum(sin((2*r-1)*x),r,1,n);
  115. - cos(2*n*x) + 1
  116. -------------------
  117. 2*sin(x)
  118. sum(cos((2*r-1)*x),r,1,n);
  119. sin(2*n*x)
  120. ------------
  121. 2*sin(x)
  122. sum(sin(n*x)**2,n);
  123. - sin(2*n*x + x) + 2*sin(x)*n
  124. --------------------------------
  125. 4*sin(x)
  126. sum(cos(n*x)**2,n);
  127. sin(2*n*x + x) + 2*sin(x)*n
  128. -----------------------------
  129. 4*sin(x)
  130. sum(sin(n*x)*sin((n+1)*x),n);
  131. - sin(2*n*x + 2*x) + sin(2*x)*n
  132. ----------------------------------
  133. 4*sin(x)
  134. sum(sec(n*x)*sec((n+1)*x),n);
  135. sum(sec(n*x + x)*sec(n*x),n)
  136. sum(1/2**n*tan(x/2**n),n);
  137. x
  138. tan(----)
  139. n
  140. 2
  141. sum(-----------,n)
  142. n
  143. 2
  144. sum(sin(r*x)*sin((r+1)*x),r,1,n);
  145. - sin(2*n*x + 2*x) + sin(2*x)*n + sin(2*x)
  146. ---------------------------------------------
  147. 4*sin(x)
  148. sum(sec(r*x)*sec((r+1)*x),r,1,n);
  149. sum(sec(r*x + x)*sec(r*x),r,1,n)
  150. sum(1/2**r*tan(x/2**r),r,1,n);
  151. x
  152. tan(----)
  153. r
  154. 2
  155. sum(-----------,r,1,n)
  156. r
  157. 2
  158. sum(k*sin(k*x),k,1,n - 1);
  159. - sin(n*x - x)*n + sin(n*x)*n - sin(n*x)
  160. -------------------------------------------
  161. 2*(cos(x) - 1)
  162. sum(k*cos(k*x),k,1,n - 1);
  163. - cos(n*x - x)*n + cos(n*x)*n - cos(n*x) + 1
  164. -----------------------------------------------
  165. 2*(cos(x) - 1)
  166. sum(sin((2k - 1)*x),k,1,n);
  167. - cos(2*n*x) + 1
  168. -------------------
  169. 2*sin(x)
  170. sum(sin(x + k*y),k,0,n);
  171. 2*n*y + 2*x + y 2*x - y
  172. - cos(-----------------) + cos(---------)
  173. 2 2
  174. --------------------------------------------
  175. y
  176. 2*sin(---)
  177. 2
  178. sum(cos(x + k*y),k,0,n);
  179. 2*n*y + 2*x + y 2*x - y
  180. sin(-----------------) - sin(---------)
  181. 2 2
  182. -----------------------------------------
  183. y
  184. 2*sin(---)
  185. 2
  186. sum((-1)**(k - 1)*sin((2k - 1)*x),k,1,n + 1);
  187. n
  188. ( - 1) *sin(2*n*x + 2*x)
  189. --------------------------
  190. 2*cos(x)
  191. sum((-1)**(k - 1)*cos((2k - 1)*x),k,1,n + 1);
  192. n
  193. ( - 1) *cos(2*n*x + 2*x) + 1
  194. ------------------------------
  195. 2*cos(x)
  196. sum(r**k*sin(k*x),k,1,n - 1);
  197. n n
  198. - r *sin(n*x - x)*r + r *sin(n*x) - sin(x)*r
  199. -----------------------------------------------
  200. 2
  201. 2*cos(x)*r - r - 1
  202. sum(r**k*cos(k*x),k,0,n - 1);
  203. n n
  204. - r *cos(n*x - x)*r + r *cos(n*x) + cos(x)*r - 1
  205. ---------------------------------------------------
  206. 2
  207. 2*cos(x)*r - r - 1
  208. sum(sin(k*x)*sin((k + 1)*x),k,1,n);
  209. - sin(2*n*x + 2*x) + sin(2*x)*n + sin(2*x)
  210. ---------------------------------------------
  211. 4*sin(x)
  212. sum(sin(k*x)*sin((k + 2)*x),k,1,n);
  213. - sin(2*n*x + 3*x) + sin(3*x)*n + sin(3*x) - sin(x)*n
  214. --------------------------------------------------------
  215. 4*sin(x)
  216. sum(sin(k*x)*sin((2k - 1)*x),k,1,n);
  217. 6*n*x + x 2*n*x - 3*x 2*n*x - x 2*n*x + x
  218. ( - sin(-----------) + sin(-------------) + sin(-----------) + sin(-----------)
  219. 2 2 2 2
  220. 3*x x 3*x
  221. + sin(-----) + sin(---))/(4*sin(-----))
  222. 2 2 2
  223. % The next examples cannot be summed in closed form.
  224. sum(1/(cos(x/2**k)*2**k)**2,k,1,n);
  225. 1
  226. sum(-----------------,k,1,n)
  227. 2*k x 2
  228. 2 *cos(----)
  229. k
  230. 2
  231. sum((2**k*sin(x/2**k)**2)**2,k,1,n);
  232. 2*k x 4
  233. sum(2 *sin(----) ,k,1,n)
  234. k
  235. 2
  236. sum(tan(x/2**k)/2**k,k,0,n);
  237. x
  238. tan(----)
  239. k
  240. 2
  241. sum(-----------,k,0,n)
  242. k
  243. 2
  244. sum(cos(k**2*2*pi/n),k,0,n - 1);
  245. 2
  246. 2*k *pi
  247. sum(cos(---------),k,0,n - 1)
  248. n
  249. sum(sin(k*pi/n),k,1,n - 1);
  250. 2*n*pi - pi pi
  251. - cos(-------------) + cos(-----)
  252. 2*n 2*n
  253. ------------------------------------
  254. pi
  255. 2*sin(-----)
  256. 2*n
  257. % 3) Expressions including the factorial function.
  258. for all n,m such that fixp m let
  259. factorial(n+m)=if m > 0 then factorial(n+m-1)*(n+m)
  260. else factorial(n+m+1)/(n+m+1);
  261. sum(n*factorial(n),n);
  262. factorial(n)*(n + 1)
  263. sum(n/factorial(n+1),n);
  264. - 1
  265. ----------------------
  266. factorial(n)*(n + 1)
  267. sum((n**2+n-1)/factorial(n+2),n);
  268. - 1
  269. ----------------------
  270. factorial(n)*(n + 2)
  271. sum(n*2**n/factorial(n+2),n);
  272. n
  273. - 2*2
  274. -----------------------------
  275. 2
  276. factorial(n)*(n + 3*n + 2)
  277. sum(n*x**n/factorial(n+2),n);
  278. n
  279. x *n
  280. sum(-----------------------------------------------------,n)
  281. 2
  282. factorial(n)*n + 3*factorial(n)*n + 2*factorial(n)
  283. for all n,m such that fixp m and m > 3 let
  284. factorial((n+m)/2)= factorial((n+m)/2-1)*((n+m)/2),
  285. factorial((n-m)/2)= factorial((n-m)/2+1)/((n-m)/2+1);
  286. sum(factorial(n-1/2)/factorial(n+1),n);
  287. 2*n - 1
  288. factorial(---------)
  289. 2
  290. sum(-------------------------------,n)
  291. factorial(n)*n + factorial(n)
  292. for all n,m such that fixp m and m > 3 clear factorial((n+m)/2);
  293. for all n,m such that fixp m and m > 3 clear factorial((n-m)/2);
  294. % 4) Expressions including combination.
  295. operator comb;
  296. % Combination function.
  297. for all n ,m let comb(n,m)=factorial(n)/factorial(n-m)/factorial(m);
  298. sum((-1)**k*comb(n,k),k,1,m);
  299. m m
  300. ( - ( - 1) *factorial(n)*m + ( - 1) *factorial(n)*n
  301. - factorial( - m + n)*factorial(m)*n)/(factorial( - m + n)*factorial(m)*n)
  302. sum(comb(n + p,q)/comb(n + r,q + 2),n,1,m);
  303. ( - factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*m*p*q
  304. - 2*factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*m*p
  305. - factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*m*q
  306. - 2*factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*m
  307. 2
  308. + factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*p*q
  309. - factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*p*q*r
  310. + 2*factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*p*q
  311. - 2*factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*p*r
  312. 2
  313. + factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*q
  314. - factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*q*r
  315. + 2*factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*q
  316. - 2*factorial( - q + r)*factorial(m + p - q)*factorial(m + r)*factorial(p)*r
  317. 2
  318. - factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*m*q
  319. + factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*m*q*r
  320. - 2*factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*m*q
  321. + 2*factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*m*r
  322. 2
  323. - factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*p*q
  324. + factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*p*q*r
  325. - 2*factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*p*q
  326. + 2*factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*p*r
  327. 2
  328. - factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*q
  329. + factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*q*r
  330. - 2*factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*q
  331. + 2*factorial(m - q + r)*factorial(m + p)*factorial(p - q)*factorial(r)*r)/(
  332. factorial(m + p - q)*factorial(m + r)*factorial(p - q)*factorial(r)*(m*p*q
  333. 2 2 2 2 2
  334. - m*p*r - m*q*r + m*q + m*r - m*r - p*q + 2*p*q*r - p*r + q *r - q
  335. 2 3 2
  336. - 2*q*r + 2*q*r + r - r ))
  337. sum((-1)**(k + 1)*comb(n,k)/(k + 1),k,1,n);
  338. n
  339. -------
  340. n + 1
  341. for all n ,m clear comb(n,m);
  342. for all n,m such that fixp m clear factorial(n+m);
  343. % 3) Examples taken from
  344. % "Decision procedure for indefinite hypergeometric summation"
  345. % Proc. Natl. Acad. Sci. USA vol. 75, no. 1 pp.40-42 (1978)
  346. % R. William Gosper, Jr.
  347. %
  348. % n
  349. % ____ 2
  350. % f = || (b*k +c*k+d)
  351. % k=1
  352. %
  353. % n
  354. % ____ 2
  355. % g = || (b*k +c*k+e)
  356. % k=1
  357. %
  358. operator f,gg;
  359. % gg used to avoid possible conflict with high energy
  360. % physics operator.
  361. for all n,m such that fixp m let
  362. f(n+m)=if m > 0 then f(n+m-1)*(b*(n+m)**2+c*(n+m)+d)
  363. else f(n+m+1)/(b*(n+m+1)**2+c*(n+m+1)+d);
  364. for all n,m such that fixp m let
  365. gg(n+m)=if m > 0 then gg(n+m-1)*(b*(n+m)**2+c*(n+m)+e)
  366. else gg(n+m+1)/(b*(n+m+1)**2+c*(n+m+1)+e);
  367. sum(f(n-1)/gg(n),n);
  368. f(n)
  369. ---------------
  370. gg(n)*(d - e)
  371. sum(f(n-1)/gg(n+1),n);
  372. 2 2 2 2
  373. (f(n)*(2*b *n + 4*b *n + 2*b + 2*b*c*n + 2*b*c + 2*b*d*n + 3*b*d - 2*b*e*n
  374. 2 2 3 2 3 3
  375. - b*e + c*d - c*e + d - 2*d*e + e ))/(gg(n)*(b *d*n + 2*b *d*n + b *d
  376. 3 2 3 3 2 2 2 2
  377. - b *e*n - 2*b *e*n - b *e + b *c*d*n + b *c*d - b *c*e*n - b *c*e
  378. 2 2 2 2 2 2 2 2 2 2 2 2 2
  379. + 2*b *d *n + 4*b *d *n + 2*b *d + b *d*e - 2*b *e *n - 4*b *e *n
  380. 2 2 2 2 2 2 2 2 2
  381. - 3*b *e - b*c *d*n - 2*b*c *d*n - b*c *d + b*c *e*n + 2*b*c *e*n
  382. 2 2 2 2 2 3 2
  383. + b*c *e + 2*b*c*d *n + 2*b*c*d - 2*b*c*e *n - 2*b*c*e + b*d *n
  384. 3 3 2 2 2 2 2 2
  385. + 2*b*d *n + b*d - 3*b*d *e*n - 6*b*d *e*n - b*d *e + 3*b*d*e *n
  386. 2 2 3 2 3 3 3 3
  387. + 6*b*d*e *n + 3*b*d*e - b*e *n - 2*b*e *n - 3*b*e - c *d*n - c *d
  388. 3 3 2 2 2 3 3 2 2
  389. + c *e*n + c *e - c *d*e + c *e + c*d *n + c*d - 3*c*d *e*n - 3*c*d *e
  390. 2 2 3 3 3 2 2 3 4
  391. + 3*c*d*e *n + 3*c*d*e - c*e *n - c*e + d *e - 3*d *e + 3*d*e - e ))
  392. for all n,m such that fixp m clear f(n+m);
  393. for all n,m such that fixp m clear gg(n+m);
  394. clear f,gg;
  395. % 4) Products.
  396. prod(n/(n+2),n);
  397. 2
  398. --------------
  399. 2
  400. n + 3*n + 2
  401. prod(x**n,n);
  402. 2
  403. (n + n)/2
  404. x
  405. prod(e**(sin(n*x)),n);
  406. 1
  407. ----------------------------------
  408. cos((2*n*x + x)/2)/(2*sin(x/2))
  409. e
  410. end;
  411. Time for test: 187 ms, plus GC time: 2 ms