Merge branch 'rtoy-refactor-specvars-defint-ll-ul'
[maxima.git] / tests / rtest16.mac
blobb1586a5ca6ed55fdcb3d2864f8a315f838e2b2ed
1 (kill (all), values);
2 [];
4 /* make sure things work after a reset().  ID: 1986726 and ID: 2787047
5  *
6  * display2d is a resettable option variable. We save the value of display2d
7  * and restore it after the reset. This allows to run the testsuite in both
8  * display modes.
9  */
10 (save:display2d, done);
11 done$
12 (reset(),0);
14 (display2d:save, done);
15 done$
17 /* From A. Reiner.  Fixed in buildq.lisp rev. 1.4.
18    Exposed a bug caused by the dynamical scope of VARLIST. */
19 buildq([foo:sin(baz+bar)],1);
22 buildq([foo:[sin(baz+bar),sin(baz-bar)]],+splice(foo));
23 sin(baz+bar)+sin(baz-bar)$
25 /* Verify that extended functionality of rhs/lhs works as advertised */
27 (kill (x, y, aa, bb, cc), infix("@@"), 0);
30 map (lhs, [aa < bb, aa <= bb, aa = bb, aa # bb, equal (aa, bb), notequal (aa, bb), aa >= bb, aa > bb]);
31 [aa, aa, aa, aa, aa, aa, aa, aa];
33 map (rhs, [aa < bb, aa <= bb, aa = bb, aa # bb, equal (aa, bb), notequal (aa, bb), aa >= bb, aa > bb]);
34 [bb, bb, bb, bb, bb, bb, bb, bb];
36 map (lhs, [foo(x) := 2*x, bar(y) ::= 3*y, '(aa : bb), '(aa :: bb), ?marrow(aa, bb)]);
37 ['(foo(x)), '(bar(y)), aa, aa, aa];
39 map (rhs, [foo(x) := 2*x, bar(y) ::= 3*y, '(aa : bb), '(aa :: bb), ?marrow(aa, bb)]);
40 [2*x, 3*y, bb, bb, bb];
42 [lhs (aa @@ bb), lhs (aa @@ bb @@ cc), rhs (aa @@ bb), rhs (aa @@ bb @@ cc)];
43 [aa, aa @@ bb, bb, cc];
45 map (lhs, [aa + bb, aa - bb, aa * bb, aa / bb, sin(aa), log(aa)]);
46 [aa + bb, aa - bb, aa * bb, aa / bb, sin(aa), log(aa)];
48 map (rhs, [aa + bb, aa - bb, aa * bb, aa / bb, sin(aa), log(aa)]);
49 [0, 0, 0, 0, 0, 0];
51 kill ("@@");
52 done;
54 /* Verify that grind treats nouns correctly. string calls MSIZE in src/grind.lisp.
55  */
57 kill (all);
58 done;
60 string ('(integrate(f(x), x) + integrate(g(x), x, minf, inf) + diff(u, x) + sum(h(x), x, 1, n)));
61 "sum(h(x),x,1,n)+integrate(g(x),x,minf,inf)+integrate(f(x),x)+diff(u,x)";
63 string ('integrate(f(x), x) + 'integrate(g(x), x, minf, inf) + 'diff(u, x) + 'sum(h(x), x, 1, n));
64 "'sum(h(x),x,1,n)+'integrate(g(x),x,minf,inf)+'integrate(f(x),x)+'diff(u,x,1)";
66 /* GREAT puts nounified atoms before others, it appears ... */
67 string (%a%a + %b%b + nounify(%c%c) + nounify(%d%d) + %e%e + %f%f);
68 "%d%d+%c%c+%f%f+%e%e+%b%b+%a%a";
70 string (sin(x) * cos(x) + tan(x));
71 "tan(x)+cos(x)*sin(x)";
73 string ('foo(x, y, z) / bar(a, b, c) + 'baz(%pi - 'quux(%e ^ mumble(%i))));
74 "'foo(x,y,z)/bar(a,b,c)+'baz(%pi-'quux(%e^mumble(%i)))";
76 /* It's conceivable that someday nounified arithmetic operators would be treated differently by grind.
77  * If/when that happens, revise this example accordingly.
78  */
79 string ('"+"(a, b, '"."(c, d), '"^"(e, f)));
80 "'?mplus(a,b,'?mnctimes(c,d),'?mexpt(e,f))";
82 /* Bug 626721 */
83 logarc(atan2(y,x));
84 -%i*log((%i*y+x)/sqrt(x^2+y^2));
85 rectform(ev(%,x=-1,y=1));
86 3*%pi/4; 
89  * Bug [ 1661490 ] An integral gives a wrong result.
90  */
91 (assume(a>0, b>0, sqrt(sqrt(b^2+a^2)-a)*(sqrt(b^2+a^2)+a)^(3/2)-b^2>0),0);
93 radcan(integrate(exp(-(a+%i*b)*x^2),x,minf,inf)/(sqrt(%pi)/sqrt(a+%i*b)));
97  * [ 1663704 ] integrate(sin(r*x)^7/x^4,x,0,inf) -> r^3*false
98  *
99  * Should return the integral instead of producing false.
100  */
101 integrate(sin(a*x)^7/x^4,x,0,inf);
102 'integrate(sin(a*x)^7/x^4,x,0,inf);
104 /* we have assumed a>0 */
105 integrate(%e^(-a*r)*sin(k*r),r,0,inf);
106 k/(k^2+a^2);
109  * Bug [ 1854888 ] hgfred([5],[5], 1) doesn't simplify
110  */
111 hgfred([5],[5],1);
115  * Bug [ 1858964 ] hgfred([7],[-1], x) --/--> error
116  */
117 hgfred([7],[-1],x);
118 und;
121  * Bug [ 1858939 ] hgfred([-1],[-2],x) --> error
122  */
123 hgfred([-1],[-2],x);
124 /* Because of revision 1.110 of hyp.lisp gen_laguerre simplifies
125    -gen_laguerre(1,-3,x)/2; */
126 1+x/2;
129  * Tests for the :: operator
130  */
131 a:b;
134 a::3;
140 p:concat('p,1);
143 p::5;
149 kill(all);
150 done$
152 /* Bug [ 1860250 ] erf(-inf) --> -erf(inf) */
154 erf(-inf);
157 erf(inf);
160 erf(-x) + erf(x);
163 erf(a-b) + erf(b-a);
166 /* Bug [ 1950653 ] bessel_j not simplified 
167  * A few additional related tests added too.
168  */
170 bessel_j(1/2,%pi),besselexpand:true;
172 bessel_y(1/2,%pi/2),besselexpand:true;
175 /* Bug [ 2149714 ] fpprintprec does not work correctly
178 fpprec:16;
181 block([fpprintprec:5], string(1.23b0));
182 "1.23b0";
184 block([fpprintprec:5], string(1.2345b0));
185 "1.2345b0";
187 block([fpprintprec:5], string(1.23456789b0));
188 "1.2346b0";
190 block([fpprintprec:25], string(1.2345678901234567890123456789b0));
191 "1.234567890123457b0";
193 /* verify that fpprintprec behavior matches its description */
195 block ([L1 : [["1.2E-10","1.2E-9","1.2E-8","1.2E-7","1.2E-6","1.2E-5","1.2E-4","0.0012","0.012","0.12","1.2","1.2E+1","1.2E+2",
196   "1.2E+3","1.2E+4","1.2E+5","1.2E+6","1.2E+7","1.2E+8","1.2E+9","1.2E+10"],
197  ["1.23E-10","1.23E-9","1.23E-8","1.23E-7","1.23E-6","1.23E-5","1.23E-4","0.00123","0.0123","0.123","1.23","12.3",
198   "1.23E+2","1.23E+3","1.23E+4","1.23E+5","1.23E+6","1.23E+7","1.23E+8","1.23E+9","1.23E+10"],
199  ["1.234E-10","1.234E-9","1.234E-8","1.234E-7","1.234E-6","1.234E-5","1.234E-4","0.001234","0.01234","0.1234","1.234",
200   "12.34","123.4","1.234E+3","1.234E+4","1.234E+5","1.234E+6","1.234E+7","1.234E+8","1.234E+9","1.234E+10"],
201  ["1.2344E-10","1.2344E-9","1.2344E-8","1.2344E-7","1.2344E-6","1.2344E-5","1.2344E-4","0.0012344","0.012344",
202   "0.12344","1.2344","12.344","123.44","1234.4","1.2344E+4","1.2344E+5","1.2344E+6","1.2344E+7","1.2344E+8",
203   "1.2344E+9","1.2344E+10"],
204  ["1.23443E-10","1.23443E-9","1.23443E-8","1.23443E-7","1.23443E-6","1.23443E-5","1.23443E-4","0.00123443","0.0123443",
205   "0.123443","1.23443","12.3443","123.443","1234.43","12344.3","1.23443E+5","1.23443E+6","1.23443E+7","1.23443E+8",
206   "1.23443E+9","1.23443E+10"],
207  ["1.234432E-10","1.234432E-9","1.234432E-8","1.234432E-7","1.234432E-6","1.234432E-5","1.234432E-4","0.001234432",
208   "0.01234432","0.1234432","1.234432","12.34432","123.4432","1234.432","12344.32","123443.2","1.234432E+6",
209   "1.234432E+7","1.234432E+8","1.234432E+9","1.234432E+10"],
210  ["1.2344321E-10","1.2344321E-9","1.2344321E-8","1.2344321E-7","1.2344321E-6","1.2344321E-5","1.2344321E-4",
211   "0.0012344321","0.012344321","0.12344321","1.2344321","12.344321","123.44321","1234.4321","12344.321","123443.21",
212   "1234432.1","1.2344321E+7","1.2344321E+8","1.2344321E+9","1.2344321E+10"],
213  ["1.23443211E-10","1.23443211E-9","1.23443211E-8","1.23443211E-7","1.23443211E-6","1.23443211E-5","1.23443211E-4",
214   "0.00123443211","0.0123443211","0.123443211","1.23443211","12.3443211","123.443211","1234.43211","12344.3211",
215   "123443.211","1234432.11","1.23443211E+7","1.23443211E+8","1.23443211E+9","1.23443211E+10"],
216  ["1.234432112E-10","1.234432112E-9","1.234432112E-8","1.234432112E-7","1.234432112E-6","1.234432112E-5",
217   "1.234432112E-4","0.001234432112","0.01234432112","0.1234432112","1.234432112","12.34432112","123.4432112",
218   "1234.432112","12344.32112","123443.2112","1234432.112","1.234432112E+7","1.234432112E+8","1.234432112E+9",
219   "1.234432112E+10"],
220  ["1.2344321123E-10","1.2344321123E-9","1.2344321123E-8","1.2344321123E-7","1.2344321123E-6","1.2344321123E-5",
221   "1.2344321123E-4","0.0012344321123","0.012344321123","0.12344321123","1.2344321123","12.344321123","123.44321123",
222   "1234.4321123","12344.321123","123443.21123","1234432.1123","1.2344321123E+7","1.2344321123E+8","1.2344321123E+9",
223   "1.2344321123E+10"],
224  ["1.23443211234E-10","1.23443211234E-9","1.23443211234E-8","1.23443211234E-7","1.23443211234E-6","1.23443211234E-5",
225   "1.23443211234E-4","0.00123443211234","0.0123443211234","0.123443211234","1.23443211234","12.3443211234",
226   "123.443211234","1234.43211234","12344.3211234","123443.211234","1234432.11234","1.23443211234E+7",
227   "1.23443211234E+8","1.23443211234E+9","1.23443211234E+10"],
228  ["1.234432112344E-10","1.234432112344E-9","1.234432112344E-8","1.234432112344E-7","1.234432112344E-6",
229   "1.234432112344E-5","1.234432112344E-4","0.001234432112344","0.01234432112344","0.1234432112344","1.234432112344",
230   "12.34432112344","123.4432112344","1234.432112344","12344.32112344","123443.2112344","1234432.112344",
231   "1.234432112344E+7","1.234432112344E+8","1.234432112344E+9","1.234432112344E+10"],
232  ["1.2344321123443E-10","1.2344321123443E-9","1.2344321123443E-8","1.2344321123443E-7","1.2344321123443E-6",
233   "1.2344321123443E-5","1.2344321123443E-4","0.0012344321123443","0.012344321123443","0.12344321123443",
234   "1.2344321123443","12.344321123443","123.44321123443","1234.4321123443","12344.321123443","123443.21123443",
235   "1234432.1123443","1.2344321123443E+7","1.2344321123443E+8","1.2344321123443E+9","1.2344321123443E+10"],
236  ["1.23443211234432E-10","1.23443211234432E-9","1.23443211234432E-8","1.23443211234432E-7","1.23443211234432E-6",
237   "1.23443211234432E-5","1.23443211234432E-4","0.00123443211234432","0.0123443211234432","0.123443211234432",
238   "1.23443211234432","12.3443211234432","123.443211234432","1234.43211234432","12344.3211234432","123443.211234432",
239   "1234432.11234432","1.23443211234432E+7","1.23443211234432E+8","1.23443211234432E+9","1.23443211234432E+10"],
240  ["1.234432112344321E-10","1.234432112344321E-9","1.234432112344321E-8","1.234432112344321E-7","1.234432112344321E-6",
241   "1.234432112344321E-5","1.234432112344321E-4","0.001234432112344321","0.01234432112344321","0.1234432112344321",
242   "1.234432112344321","12.34432112344321","123.4432112344321","1234.432112344321","12344.32112344321",
243   "123443.2112344321","1234432.112344321","1.234432112344321E+7","1.234432112344321E+8","1.234432112344321E+9",
244   "1.234432112344321E+10"]],
245  L2 : block ([foo : 1.2344321123443211234],
246              makelist (block ([fpprintprec : m], makelist ([m, string (foo*10^n)], n, -10, 10)), m, 2, 16))],
247  map (lambda ([s1, s2], if sequalignore (s1, second (s2)) then true else ['fpprintprec = first (s2), 'expected = s1, 'actual = second (s2)]), apply (append, L1), apply (append, L2)),
248  delete (true, %%));
251 /* verify that fpprintprec = 0 or  fpprintprec > 16 prints readably */
253 block ([foo: float (%pi), L],
254        L: makelist (block ([fpprintprec: m], makelist ([m, foo*10^n, string (foo*10^n)], n, -10, 10)), m, [0, 17, 18, 19, 20]),
255        map (lambda ([L1], if parse_string (L1[3]) = L1[2] then true else ['fpprintprec = first (L1), 'expected = L1[2], 'actual = L1[3]]), apply (append, L)),
256        delete (true, %%));
259 /* SF bug #3213: "fpprintprec do not round bfloat correctly." */
261 (reset (fpprec, bftrunc),
262  fpprintprec:4,
263  string (1.23456789b0));
264 "1.235b0";
266 /* default fpprintprec => print all digits */
268 (reset (fpprintprec), 0);
271 string (bfloat (1000000000000*(1 + 8/9)));
272 "1.888888888888889b12";
274 string (bfloat ((1 + 8/9)/1000000000000)), fpprec=50;
275 "1.8888888888888888888888888888888888888888888888889b-12";
277 string (bfloat (1000000000000*(1 + 8/9))), fpprec=10;
278 "1.888888889b12";
280 string (bfloat ((1 + 8/9)/1000000000000)), fpprec=2;
281 "1.9b-12";
283 /* expect to see bigfloat ending in 5 to round to even */
285 map (string, map (bfloat, 1 + [1, 2, 3, 4, 5, 6, 7]/8)), fpprintprec=3;
286 ["1.12b0", "1.25b0", "1.38b0", "1.5b0", "1.62b0", "1.75b0", "1.88b0"];
288 map (string, map (bfloat, 100 + makelist (i/64, i, 1, 63))), fpprintprec=8;
289 ["1.0001562b2", "1.0003125b2", "1.0004688b2", "1.000625b2", "1.0007812b2", "1.0009375b2", "1.0010938b2", "1.00125b2",
290  "1.0014062b2", "1.0015625b2", "1.0017188b2", "1.001875b2", "1.0020312b2", "1.0021875b2", "1.0023438b2", "1.0025b2",
291  "1.0026562b2", "1.0028125b2", "1.0029688b2", "1.003125b2", "1.0032812b2", "1.0034375b2", "1.0035938b2", "1.00375b2",
292  "1.0039062b2", "1.0040625b2", "1.0042188b2", "1.004375b2", "1.0045312b2", "1.0046875b2", "1.0048438b2", "1.005b2",
293  "1.0051562b2", "1.0053125b2", "1.0054688b2", "1.005625b2", "1.0057812b2", "1.0059375b2", "1.0060938b2", "1.00625b2",
294  "1.0064062b2", "1.0065625b2", "1.0067188b2", "1.006875b2", "1.0070312b2", "1.0071875b2", "1.0073438b2", "1.0075b2",
295  "1.0076562b2", "1.0078125b2", "1.0079688b2", "1.008125b2", "1.0082812b2", "1.0084375b2", "1.0085938b2", "1.00875b2",
296  "1.0089062b2", "1.0090625b2", "1.0092188b2", "1.009375b2", "1.0095312b2", "1.0096875b2", "1.0098438b2"];
298 /* bftrunc=false => don't strip trailing 0's */
300 map (string, map (bfloat, 1 + [1, 2, 3, 4, 5, 6, 7]/8)), fpprintprec=3, bftrunc=false;
301 ["1.12b0", "1.25b0", "1.38b0", "1.50b0", "1.62b0", "1.75b0", "1.88b0"];
303 map (string, map (bfloat, 100 + makelist (i/64, i, 1, 63))), fpprintprec=8, bftrunc=false;
304 ["1.0001562b2", "1.0003125b2", "1.0004688b2", "1.0006250b2", "1.0007812b2", "1.0009375b2", "1.0010938b2", "1.0012500b2",
305  "1.0014062b2", "1.0015625b2", "1.0017188b2", "1.0018750b2", "1.0020312b2", "1.0021875b2", "1.0023438b2", "1.0025000b2",
306  "1.0026562b2", "1.0028125b2", "1.0029688b2", "1.0031250b2", "1.0032812b2", "1.0034375b2", "1.0035938b2", "1.0037500b2",
307  "1.0039062b2", "1.0040625b2", "1.0042188b2", "1.0043750b2", "1.0045312b2", "1.0046875b2", "1.0048438b2", "1.0050000b2",
308  "1.0051562b2", "1.0053125b2", "1.0054688b2", "1.0056250b2", "1.0057812b2", "1.0059375b2", "1.0060938b2", "1.0062500b2",
309  "1.0064062b2", "1.0065625b2", "1.0067188b2", "1.0068750b2", "1.0070312b2", "1.0071875b2", "1.0073438b2", "1.0075000b2",
310  "1.0076562b2", "1.0078125b2", "1.0079688b2", "1.0081250b2", "1.0082812b2", "1.0084375b2", "1.0085938b2", "1.0087500b2",
311  "1.0089062b2", "1.0090625b2", "1.0092188b2", "1.0093750b2", "1.0095312b2", "1.0096875b2", "1.0098438b2"];
313 /* fpprintprec <= fpprec */
315 every (lambda ([s], s="5.6b2"), map (string, makelist (bfloat (5000/9), fpprec, fpprintprec, 40))), fpprintprec=2;
316 true;
318 every (lambda ([s], s="5.56b-4"), map (string, makelist (bfloat (5/9000), fpprec, fpprintprec, 40))), fpprintprec=3;
319 true;
321 every (lambda ([s], s="5.556b2"), map (string, makelist (bfloat (5000/9), fpprec, fpprintprec, 40))), fpprintprec=4;
322 true;
324 every (lambda ([s], s="5.5556b-4"), map (string, makelist (bfloat (5/9000), fpprec, fpprintprec, 40))), fpprintprec=5;
325 true;
327 every (lambda ([s], s="5.55556b2"), map (string, makelist (bfloat (5000/9), fpprec, fpprintprec, 40))), fpprintprec=6;
328 true;
330 every (lambda ([s], s="5.555556b-4"), map (string, makelist (bfloat (5/9000), fpprec, fpprintprec, 40))), fpprintprec=7;
331 true;
333 every (lambda ([s], s="5.5555556b2"), map (string, makelist (bfloat (5000/9), fpprec, fpprintprec, 40))), fpprintprec=8;
334 true;
336 /* fpprintprec >= fpprec */
338 /* bfloat generally produces a number which is not exactly the same as its argument.
339  * When bfloat returns a bigger number, it should round up when formatted,
340  * and when it's smaller, it should round down.
341  * When bfloat returns an exact result, in this case it should round up (to 6).
342  * Consider this when figuring out the correct output for these examples.
343  */
345 (fpprec : 2, 
346  mybf : bfloat(5/9000),
347  is (rationalize (mybf) >= (5*sum(10^i, i, 0, fpprec)/10^(fpprec + 1))/1000));
348 true;
350 every (lambda ([s], s="5.6b-4"), map (string, makelist (mybf, fpprintprec, fpprec, 40)));
351 true;
353 (fpprec : 3,
354  mybf : bfloat(5000/9),
355  is (rationalize (mybf) >= (5*sum(10^i, i, 0, fpprec)/10^(fpprec + 1))*1000));
356 true;
358 every (lambda ([s], s="5.56b2"), map (string, makelist (mybf, fpprintprec, fpprec, 40)));
359 true;
361 (fpprec : 4,
362  mybf : bfloat(5/9000),
363  is (rationalize (mybf) >= (5*sum(10^i, i, 0, fpprec)/10^(fpprec + 1))/1000));
364 true;
366 every (lambda ([s], s="5.556b-4"), map (string, makelist (mybf, fpprintprec, fpprec, 40)));
367 true;
369 (fpprec : 5,
370  mybf : bfloat(5000/9),
371  is (rationalize (mybf) >= (5*sum(10^i, i, 0, fpprec)/10^(fpprec + 1))*1000));
372 false;
374 every (lambda ([s], s="5.5555b2"), map (string, makelist (mybf, fpprintprec, fpprec, 40)));
375 true;
377 (fpprec : 10,
378  mybf : bfloat(5/9000),
379  is (rationalize (mybf) >= (5*sum(10^i, i, 0, fpprec)/10^(fpprec + 1))/1000));
380 true;
382 every (lambda ([s], s="5.555555556b-4"), map (string, makelist (mybf, fpprintprec, fpprec, 40)));
383 true;
385 (fpprec : 20,
386  mybf : bfloat(5000/9),
387  is (rationalize (mybf) >= (5*sum(10^i, i, 0, fpprec)/10^(fpprec + 1))*1000));
388 true;
390 every (lambda ([s], s="5.5555555555555555556b2"), map (string, makelist (mybf, fpprintprec, fpprec, 40)));
391 true;
393 (fpprec : 40,
394  mybf : bfloat(5/9000),
395  is (rationalize (mybf) >= (5*sum(10^i, i, 0, fpprec)/10^(fpprec + 1))/1000));
396 false;
398 every (lambda ([s], s="5.555555555555555555555555555555555555555b-4"), map (string, makelist (mybf, fpprintprec, fpprec, 40)));
399 true;
401 (reset (fpprec), 0);
405  * Bug 2142758: integrate(sqrt(2-2*x^2)*(sqrt(2)*x^2+sqrt(2))/(4-4*x^2),x,0,1)
406  */
407 integrate(sqrt(2-2*x^2)*(sqrt(2)*x^2+sqrt(2))/(4-4*x^2),x,0,1);
408 3*%pi/8;
411 integrate(sqrt(1-x^2)*(x^2+1)/(2-2*x^2),x,0,1);
412 3*%pi/8;
414 integrate(sqrt(1-x^2)*(x^2+1)/(1-x^2),x,0,1);
415 3*%pi/4;
418  * Bug [ 2208303 ] Problem with jacobi_dn and elliptic_kc
419  */
420 jacobi_dn(elliptic_kc(m)*t,m);
421 jacobi_dn(elliptic_kc(m)*t,m);
424  * Bug [ 2180110 ] GCL do not signal an overflow converting bigfloat to float
425  */
426 errcatch(float(2b400));
429 errcatch(float(bfloat(2^1024)));
433  * Bug [ 2055235 ] Plot leaves range with jacobi functions
435  * Actually jacobi_cn(100, .7) is computed inaccurately.  Just check that abs(jacobi_cn(100,.7)) < 1
436  */
438 is(abs(jacobi_cn(100.0, 0.7)) < 1);
439 true$
442  * Bug [ 1658067 ] jacobi_sn(elliptic_kc(1-m)*%i/2,m) isn't simplified
444  * This test (and Maxima) used to be wrong.  This is related to the
445  * jacobi_sc(elliptic_kc(m)/2,m) test below.
446  */
447 jacobi_sn(elliptic_kc(1-m)*%i/2,m);
448 %i/m^(1/4)$
450 jacobi_sc(elliptic_kc(m)+u,m);
451 -jacobi_cs(u,m)/sqrt(1-m)$
454  * Maxima used to get this wrong by returning the reciprocal instead
455  */
456 jacobi_sc(elliptic_kc(m)/2,m);
457 1/(1-m)^(1/4);
460  * Bug [ 2505945 ] - hgfred([2,-1/2],[3],-x^2);
462  * Shouldn't signal from diff about non-variable second arg.
464  * The expected value here is computed from
465  * factor(ratsimp(subst([z=-x^2],hgfred([2,-1/2],[3],z))))
466  */
467 factor(ratsimp(hgfred([2,-1/2],[3],-x^2)));
468 4*(3*x^4*sqrt(x^2+1)+x^2*sqrt(x^2+1)-2*sqrt(x^2+1)+2)/(15*x^4)$
471  * Bug 2534420: asinh(%i*2b0) causes error
472  */
473 is(abs(asinh(%i*2b0)-expand(bfloat(asinh(%i*2)))) < 3b-17);
474 true;
477  * Bug 2543079: bfloat(gamma(3/4)/gamma(1/4)) is wrong.
478  */
479 bfloat(gamma(3/4)/gamma(1/4));
480 3.379891200336424b-1;
483  * Bug 2582034 - hgfred([a/2,-a/2],[1/2],z) causes error
484  */
485 (assume(zn<0), done);
486 done;
488 hgfred([a/2,-a/2],[1/2],zn);
489 ((%i*sqrt(zn)+sqrt(1-zn))^a+(sqrt(1-zn)-%i*sqrt(zn))^a)/2$
492  * Bug 2618401 - bfloat produces incorrect answer
493  */
494 is(abs(bfloat((sqrt(2)+2)*%pi^(3/2)/(8*gamma(3/4)^2))-float((sqrt(2)+2)*%pi^(3/2)/(8*gamma(3/4)^2))) < 1d-15);
495 true;
497 /* (-1.0b0)^(1/3) vs (-1.0d0)^(1/3) - ID: 619927 */
498 (-1b0)^(1/3);
499 -1.0b0;
501 (-1.0)^(1/3);
502 -1.0;
504 domain:complex;
505 complex;
507 (-1b0)^(1/3);
508 1.0b0*(-1)^(1/3);
510 (-1.0)^(1/3);
511 1.0*(-1)^(1/3);
513 domain:real;
514 real;
516 (-1b0)^(1/3),domain:complex,m1pbranch:true;
517 1.0b0*(sqrt(3)*%i/2+1/2);
519 (-1.0)^(1/3),domain:complex,m1pbranch:true;
520  1.0*(sqrt(3)*%i/2+1/2);
524  * Bug [ 2688847 ] float of rats rounds incorrectly
525  */
526 float((2^60-1)/2^60)-1;
527 0.0;
528 float((2^1000-1)/2^1000)-1;
529 0.0;
532  * Bug [ 2687962 ] hgfred([-3/2,1],[-1/2],-t) division by zero
534  * Solution from functions.wolfram.com
535  */
536 ratsimp(hgfred([-3/2,1],[-1/2], t));
537 1+3*t-3*t^(3/2)*atanh(sqrt(t));
540  * Bug 2793827: internal error in integrate
541  */
542 (assume(n>0),declare(n,integer),0);
545 integrate((g32475^n*(g32475*n-n-1)/(g32475-1)^2+1/(g32475-1)^2)/(1-g32475)
546  -(g32475^(2*n+1)*(g32475*n-n-1)/(g32475-1)^2+g32475^(n+1)/(g32475-1)^2)
547   /(1-g32475),g32475,0,1);
548 (n^2+2*n+1)/2-1/2$
550 kill(all);
551 done;
554  * Bug 609464 : 1+%e,numer and %e^%e,numer
556  * The simplifier has been extended to handle %e like other constants.
557  * In addition functions with arguments which involve %e simplify
558  * accordingly.
559  */
561 %e,numer;
562 2.7182818284590451;
564 %e+1,numer;
565 3.7182818284590451;
567 %e^%e,numer;
568 15.154262241479262;
570 %e^x,numer;
571 %e^x;
573 sin(%e),numer;
574 0.41078129050290885;
576 sin(%e+1),numer;
577 -0.54525155669233449;
579 /* Do not simplify, when %e is the base of an expression and %enumer FALSE*/
581 sin(%e^(2*x+1)),numer;
582 sin(%e^(2*x+1));
584 sin(%e^(%e^(2*x+1))),numer;
585 sin(%e^(%e^(2*x+1)));
587 /* Additionally simplifications when %enumer TRUE */
589 %enumer:true;
590 true;
592 sin(%e^x),numer;
593 sin(2.7182818284590451^x);
595 sin(%e^(%e^(2*x+1))),numer;
596 sin(2.7182818284590451^(2.7182818284590451^(2*x+1)));
598 %enumer:false;
599 false;
602  * Bug ID: 2797885 - "problem with integration"
604  * integrate(exp(%i*x)*sin(x),x) generates a Lisp error.
606  * This is a special case for the integrand: exp(a*x)*sin(b*x),
607  * with a^2+b^2 equal to zero.
608  */
610 /* This is the general case for an integral with exp and sin or cos */
611 integrate(exp(a*x)*sin(b*x),x);
612 %e^(a*x)*(a*sin(b*x)-b*cos(b*x))/(b^2+a^2);
614 integrate(exp(a*x)*cos(b*x),x);
615 %e^(a*x)*(b*sin(b*x)+a*cos(b*x))/(b^2+a^2);
617 /* Now the special case with a=%i and b=1 */
618 expand(integrate(exp(%i*x)*sin(x),x));
619 %i*x/2-%e^(2*%i*x)/4;
621 expand(integrate(exp(x)*sin(%i*x),x));
622 %i*%e^(2*x)/4-%i*x/2;
624 expand(integrate(exp(%i*x)*cos(x),x));
625 x/2-%i*%e^(2*%i*x)/4;
627 expand(integrate(exp(x)*cos(%i*x),x));
628 %e^(2*x)/4+x/2;
630 /* Bug ID: 932076 -  ode2( 'diff(y,x)=%i*y+sin(x), y, x) => div by 0
632  * This bug is related to the Bug ID: 2797885 - "problem with integration"
633  */
635 ode2('diff(y,x)-%i*y-sin(x),y,x);
636 y = (%c-%i*(x-%i*%e^-(2*%i*x)/2)/2)*%e^(%i*x);
639  * Bug ID: 826623 "simplifer returns %i*%i"
640  * 
641  * Some examples to show simplification of expressions of the form
642  * (a*b*...)^q*(a*b*...)^r, where q+r=1
643  */
645 sqrt(-%i)*sqrt(-%i)*%i;
648 sqrt(a*b)*sqrt(a*b)*a*b;
649 a^2*b^2;
651 (a*b*c)^(3/4)*(a*b*c)^(1/4)*c;
652 a*b*c^2;
655  * Bug ID: 2792493 "hgfred([1],[-5.2],x);"
656  */
657 hgfred([1],[-5.2],x);
658 %f[1,1]([-6.2],[-5.2],-x)*%e^x$
660 /* BUG ID: 721575 2/sqrt(2) doesn\'t simplify */
661 2/sqrt(2);
662 sqrt(2);
664 (1/2)*sqrt(2);
665 1/sqrt(2);
667 sqrt(2)*(1/2);
668 1/sqrt(2);
670 /* BUG ID 2029041 a*sqrt(2)/2 unsimplified */
672 a*sqrt(2)/2;
673 a/sqrt(2);
675 /* BUG ID 1923119 1/sqrt(8)-sqrt(8)/8 */
677 1/sqrt(8)-sqrt(8)/8;
680 /* BUG ID 1927178 integrate(sin(t),t,%pi/4,3*%pi/4) */
682 integrate(sin(t),t,%pi/4,3*%pi/4);
683 sqrt(2);
685 /* BUG ID: 1480562 2*a*2^k isn't simplified to a*2^(k+1) */
687 2*a*2^k;
688 a*2^(k+1);
690 a*2^k*2;
691 a*2^(k+1);
693 /* Some examples to show simplification of expressions
694  * with floating point and bigfloat numbers after improvement
695  * of plusin
696  */
698 (4.0*x-4.0*x);
699 0.0;
700 (4.0*x-3.0*x);
701 1.0*x;
702 (4.0*x-3.0*x)/2;
703 0.5*x;
705 (4.0b0*x-4.0*x);
706 0.0b0;
707 (4.0b0*x-3.0*x);
708 1.0b0*x;
709 (4.0b0*x-3.0*x)/2;
710 0.5b0*x;
712 /* BUG ID: 1996354 unsimplifed result from expand */
714 expand((%e^(-2*sqrt(2))*(%e^(2*sqrt(2))+2*%e^sqrt(2)+1)^2)/16
715        +(%e^(-2*sqrt(2))*(%e^(2*sqrt(2))-2*%e^sqrt(2)+1)^2)/16
716        -(%e^(-2*sqrt(2))*(%e^(2*sqrt(2))-1)^2)/8);
719 /* BUG ID: 631216 - "horner([...],x)/FIX" 
720    horner now maps over lists, matrices and equations.
721  */
723 horner(x^2+x=a*x^2+b*x);
724 x*(x+1) = x*(a*x+b);
725 horner([x^2+x,x^3+x,x^4+x]);
726 [x*(x+1),x*(x^2+1),x*(x^3+1)];
728 /* BUG ID: 2699862 "derivative of polylogarithm"
729  * The noun form is not put on the property list, but NIL. The routine 
730  * sdiffgrad generates a noun form, when the derivative is not known.
731  */
733 diff(li[n](x),n);
734 'diff(li[n](x),n);
736 diff(li[n*x](x),x);
737 'diff(li[n*x](x),x);
739 diff(li[n](x),x,1,n,1);
740 'diff(li[n-1](x),n)/x;
742 /* Not reported as a bug, but the same problem for the function psi */
744 diff(psi[n](x),n);
745 'diff(psi[n](x),n);
747 diff(psi[n*x](x),x);
748 'diff(psi[n*x](x),x);
750 diff(psi[n](x),x,1,n,1);
751 'diff(psi[n+1](x),n);
753 /* BUG ID: 2824909 " exp(%i*%pi/4) not simplified" 
754  * Check the simplification of exp(%i*%pi/4) and exp(-%i*pi/4)
755  */
757 exp(%i*%pi/4);
758 1/sqrt(2)+%i/sqrt(2);
759 exp(-%i*%pi/4);
760 1/sqrt(2)-%i/sqrt(2);
763  * Bug ID: 2831259 - bfloat() underflow bug
764  */
765 fpprec:500;
766 500;
767 float(0b0);
768 0.0;
771  * BUG ID: 2835098 - SIGN-PREP strangeness
772  */
774 block ([?limitp : true], sign (foo (x)));
775 pnz;
777 integrate(sqrt(2*m*(E[n]-U(x))),x,-x[0],x[0])=(n-1/2)*%pi*hbar;
778 sqrt(2)*'integrate(sqrt(m*(E[n]-U(x))),x,-x[0],x[0]) = %pi*hbar*(n-1/2);
780 integrate(f(x),x,x[0],x[1]);
781 'integrate(f(x),x,x[0],x[1]);
784  * BUG ID: 2840566 - defint fails to determine if one of its limit is real
785  */
787 (assume(b>0,c>0),done);
788 done;
790 integrate(x,x,0,sqrt(b^2+(b-c)^2));
791 (c^2-2*b*c+2*b^2)/2;
794  * BUG ID: 2842060 - unsimplified result from integrate
795  */
797 /* The result for a general symbol x */
798 integrate(1/x/sqrt(x^2-1),x);
799 -asin(1/abs(x));
801 (assume(x>0), done);
802 done;
804 /* abs(x) simplifies to x for x>0 */
805 integrate(1/x/sqrt(x^2-1),x);
806 -asin(1/x);
808 (forget(x>0), done);
809 done;
811 /* 
812  * Bug ID: 2820202 - rootscontract(%i/2) 
813  */
814 rootscontract(%i/2);
815 %i/2;
817 /* Bug ID: 2872738 - sign(-(1/n)*(-1)^n)
818  * We got the error because of the simplification
819  *  (-1)^n*(-1) -> (-1)^(n+1) and not -(-1)^n
820  * The other case 
821     (-1)*(-1)^n simplifies already to -(-1)^n
822  * Adding tests for both cases.
823  */
824 kill(all);
825 done;
826 sign(-(1/n)*(-1)^n);
828 (-1)*(-1)^n;
829 -(-1)^n;
830 (-1)^n*(-1);
831 -(-1)^n;
833 /* Bug ID: 2835634 - logcontract broken
834  * Bug ID: 1467368 - logcontract returns unsimplified expr
835  */
836 logcontract(log(x)-log(2));
837 log(x/2);
838 /* Check that we do not break the following again */
839 logcontract(log(%e*k)-log(%e^-1*k));
841 log(%e^2),logexpand:false;
844 /* Bug ID: 2880923 - realpart --> floating-point-overflow
845  */
846 sign(exp(2009));
847 pos;
848 realpart(sqrt(4*%e^2009-3)-1);
849 sqrt(4*%e^2009-3)-1;
850 sqrt(4*exp(2009));
851 2*%e^(2009/2);
853 /* Bug ID: 640332 - Need to specdisrep more systematically
854    Add the examples of the bug report.
855  */
856 ratdisrep(diff(rat(x),rat(x)));
858 diff(x,rat(x));
860 outofpois(diff(intopois(sin(x)),x));
861 cos(x);
862 taylor(intopois(sin(x)),x,0,3);
863 x-x^3/6;
864 ratsimp(intopois(sin(x)));
865 sin(x);
867 /* Bug ID: 627759 - Ratdisrep of aggregates 
868  */
869 ratdisrep(rat(x=y));
870 x = y;
871 ratdisrep(rat([x=a,y=b]));
872 [x = a,y = b];
873 ratdisrep(rat(matrix([a,b],[c,d])));
874 matrix([a,b],[c,d]);
876 /* Bug ID: 711885 - Rootscontract with imaginaries fails
877  */
878 (oldvalue:radexpand, radexpand:false, done);
879 done;
881 rootscontract(((sqrt(3)*%i+1)^(3/2)-4*%i)/sqrt(sqrt(3)*%i+1));
882 ((sqrt(-3)+1)^(3/2)-4*%i)/sqrt(sqrt(-3)+1);
884 /* It is a problem of the simplifier. Show that it works */
885 sqrt(1/(1+sqrt(-3)));
886 1/sqrt(sqrt(-3)+1);
888 (radexpand:oldvalue, done);
889 done;
891 /* BUG ID: 767556 - Distributing operations over =
892  * The operators "." and "^^" distribute over equations.
893  */
895 x . (a=b);
896 x . a = x . b;
898 (a=b)^^x;
899 a^^x = b^^x;
901 /* A more complicated example */
902 x . ((2*a+b . c) = x . (y + z))^^w;
903 x . (b . c+2*a)^^w = x . (x . (z+y))^^w;
905 /* Bug ID: 2914176 - Conversion of rational to bfloat is inaccurate
907  * The difference should be 1/262144, but we don't check for that.
908  */
909 (oldfpprec:fpprec, fpprec:5, done);
910 done;
911 is(bfloat((2^20+1)/(2^20-1)) - 1b0 > 0);
912 true;
914 /* Related to the fix for 2914176.  Didn't handle the ratio 0/1 */
915 is(equal(0b0, 0));
916 true;
918 (fpprec:oldfpprec, done);
919 done;
921 /* Bug ID:2933882 - Power function: 0^a not fully implemented
922  * Show some simplifications of 0^a
923  */
925 assume(a>0);
926 [a>0];
927 0^a;
929 errcatch(0^-a);
931 0^(a+%i);
933 0^(1/2+%i);
935 errcatch(0^(-1/2+%));
937 errcatch(0^%i);
938 []; 
940 forget(a>0);
941 [a>0];
943 /* Bug ID: 2938078 - Crash on attached input
944  */
946 declare(n,integer, j,noninteger);
947 done;
948 assume(equal(x,n), equal(y,j), equal(z,i));
949 [equal(x, n), equal(y, j), equal(z,i)];
951 featurep(x,integer);
952 true;
953 featurep(x,noninteger);
954 false;
955 featurep(y,integer);
956 false;
957 featurep(y,noninteger);
958 true;
960 diff(z+1,z);
963 remove(n,integer, j,noninteger);
964 done;
965 forget(equal(x,n), equal(y,j), equal(z,i));
966 [equal(x, n), equal(y, j), equal(z, i)];
968 /* Bug ID:  2948800 - integrate((1-cos(2*x)^2)^2/x^4,x,0,inf) wrong
969  */
971 integrate((1-cos(2*x)^2)^2/x^4,x,0,inf);
972 8*%pi/3;
974 assume(a>0);
975 [a>0];
977 /* The more general type with an argument a*x and a positive */
978 integrate((1-cos(a*x)^2)^2/x^4,x,0,inf);
979 %pi*a^3/3;
981 forget(a>0);
982 [a>0];
984 /* Bug ID: 777564 - subtraction "-"(a,b) should work/FIX */
986 "-"();
989 "-"(a);
992 "-"(2*a);
993 -2*a;
995 "-"(a+b);
996 -b-a;
998 "-"(a+b+c);
999 -c-b-a;
1001 "-"(100,20,10);
1004 map("-",[a,x,100],[b,y,20]);
1005 [a-b,x-y,80];
1007 map("-",[a,x,100],[b,y,20],[c,z,10]);
1008 [-c-b+a,-z-y+x,70];
1010 /* Bug ID: 910270 - 1/+3*x parses as 1/(+3*x)
1011  * Show that the "+" operator can be used as a prefix operator too.
1012  */
1013 1/+3*x;
1014 x*1/3;
1015 1/+x/3;
1016 1/(3*x);
1017 a^+b*c;
1018 c*a^b;
1020 /* Bug ID: 2961822 - sinh(0.0b0) causes Maxima to abort
1021  */
1022 sinh(0.0b0);
1023 0.0b0;
1025 /* Bug ID: 1219846 - properties of translated functions
1026  * The property noun is already present
1027  */
1028 kill(f);
1029 done;
1030 f(x):=x;
1031 f(x):=x;
1032 properties(f);
1033 [function,noun];
1034 translate(f);
1035 [f];
1036 properties(f);
1037 [transfun,function,noun];
1038 kill(f);
1039 done;
1041 /* Bug ID: 2968344 - gamma_incomplete(1.0, 4.368265444147715e+19) fails
1042  */
1043 gamma_incomplete(1.0, 4.368265444147715e+19);
1044 0.0;
1046 /* Bug ID: 643254 - orderlessp([rat(x)], [rat(x)])
1047  */
1048 orderlessp([rat(x)],[rat(x)]);
1049 false;
1051 /* Bug ID: 781657 - binomial(x,x) => 1, but binomial(-1,-1) => 0
1052  * binomial(x,x) simplifies to 1 only if the sign of x is known not to be 
1053  * negative.
1054  */
1055 is(equal(binomial(x,x),1));
1056 'unknown;
1057 is(equal(binomial(x^2,x^2),1));
1058 true;
1060 /* Bug ID:856209 - inconsistency between facts() and facts(v)
1061  * Show that facts(expr) now works more general.
1062  */
1063 assume(z+a>0,b>z);
1064 [a+z>0,b>z];
1065 facts(a);
1066 [a+z>0];
1067 facts(b);
1068 [b>z];
1069 facts(z);
1070 [a+z>0,b>z];
1071 facts(a+z);
1072 [a+z>0];
1073 forget(z+a>0,b>z);
1074 [a+z>0,b>z];
1076 /* Bug ID: 840848 - trigreduce doesn't enter unknown functions
1077  */
1078 trigexpand(f(sin(2*x)));
1079 f(2*cos(x)*sin(x));
1080 trigreduce(%);
1081 f(sin(2*x));
1083 /* Bug ID: 2954472 - rectform with large floats gives bad answer
1085 is(abs(rectform(1e160/(1e160+%i))-1) < 1e-160);
1086 true;
1088 is(abs(rectform(1e160/(1e160+3/2*%i))-1) < 1.5e-160);
1089 true;
1091 /* Bug ID: 2953369 - Definite Integration of 1/(a-b*cos(x)) wrong
1093  * For simplicity we test the equivalent integrate(1/(1-r*cos(x)),x,0,%pi).
1094  */
1095 /* These assumes are to answer the questions integrate (from routine unitcir) will ask */
1096 (assume(r>0,r<1,abs(sqrt(1-r^2)-1)/r-1 < 0, sqrt(1-r^2)-r+1>0), 0);
1098 integrate(1/(1-r*cos(x)),x,0,%pi);
1099 %pi/sqrt(1-r^2);
1101 /* Bug ID: 2907727 - Incorrect Integral with option integrate_use_rootsof
1102  * :true
1103  */
1104 %rnum:0;
1106 integrate((d*x^2+2*c*x+3*b)/(g*r*x^3+d*x^2+c*x+b), x), integrate_use_rootsof:true;
1107 lsum((%r1^2*d+2*%r1*c+3*b)*log(x-%r1)/(3*%r1^2*g*r+2*%r1*d+c),%r1,
1108              rootsof(g*r*%r1^3+d*%r1^2+c*%r1+b,%r1));
1110 /* Bug ID: 2880797 - bad answer in integrate(sqrt(sin(t)^2+cos(t)^2),t,0,2*%pi)
1112  */
1113 integrate(sqrt(sin(t)^2+cos(t)^2),t,0,2*%pi);
1114 2*%pi;
1116 /* Bug ID: 2980551 - Inconsistent simplification of exp(x*%i*%pi)
1118  * These examples show consistent simplification for x an expression which
1119  * can contain float or bigfloat values
1120  */
1122 exp(2*%i*%pi);
1125 exp((2+x)*%i*%pi);
1126 exp(x*%i*%pi);
1128 exp(2*%i*%pi+x*%i*%pi);
1129 exp(x*%i*%pi);
1131 log(exp((2+x)^2*%i*%pi));
1132 (2+x)^2*%i*%pi;
1134 exp(2.0*%i*%pi);
1135 1.0;
1137 exp((2.0+x)*%i*%pi);
1138 exp(1.0*x*%i*%pi);
1140 exp(2.0*%i*%pi+x*%i*%pi);
1141 exp(1.0*x*%i*%pi);
1143 log(exp((2.0+x)^2*%i*%pi));
1144 (2.0+x)^2*%i*%pi;
1146 exp(2.0b0*%i*%pi);
1147 1.0b0;
1149 exp((2.0b0+x)*%i*%pi);
1150 exp(1.0b0*x*%i*%pi);
1152 exp(2.0b0*%i*%pi+x*%i*%pi);
1153 exp(1.0b0*x*%i*%pi);
1155 log(exp((2.0b0+x)^2*%i*%pi));
1156 (2.0b0+x)^2*%i*%pi;
1158 exp(3/2*%pi*%i);
1159 -%i;
1160 exp(1.5*%pi*%i);
1161 -1.0*%i;
1162 exp(1.5b0*%pi*%i);
1163 -1.0b0*%i;
1165 exp((3/2+x)*%pi*%i);
1166 -%i*exp(%i*%pi*x);
1167 exp((1.5+x)*%pi*%i);
1168 -%i*exp(1.0*%i*%pi*x);
1169 exp((1.5b0+x)*%pi*%i);
1170 -%i*exp(1.0b0*%i*%pi*x);
1172 /* Bug ID: 2781127 - bfpsi0 of complex
1174  * (The result was not in rectangular form but it should be.)
1175  */
1176 bfpsi0(4.5 + %i,15);
1177 2.43845477527606b-1*%i + 1.41875534014717b0;
1179 /* Bug ID: 2988190 - atan2(1b20,-1b0); badly wrong
1180  * It's really a bug in atan for x < -1, so test both.
1181  */
1182 (fpprec:16, atan(-1b20));
1183 -1.570796326794897b0;
1184 atan2(1b20,-1b0);
1185 1.570796326794897b0;
1188  * Bug ID: 2991924 - Incorrect integration of rational functions
1189  */
1190 integrate(1/(x^4-2),x,0,1) - integrate(1/(x^2-sqrt(2))/(x^2+sqrt(2)),x,0,1);
1193 integrate(1/(x^6-4),x,0,1) - integrate(1/(x^3-2)/(x^3+2),x,0,1);
1196 /* BUG ID:  2113751 - Incomprehensible behavior of coeff()
1197  */
1198 coeff(2*%e^x, x, 0);
1199 2*%e^x;
1201 /* For numerical tests */
1202 closeto(e,tol):=block([numer:true,abse],abse:abs(e),if(abse<tol) then true else abse);
1203 closeto(e,tol):=block([numer:true,abse],abse:abs(e),if(abse<tol) then true else abse);
1205 /* Bug ID: 2997276 - zeta(3),numer; gives Lisp error 
1206  * 
1207  * Also add a test for complex rational argument, which wasn't handled
1208  * correctly either.
1210  * Some Lisp implementations fail these tests because things like
1211  * (cl:expt 2d0 3) only gives single-float accuracy (but with
1212  * double-float precision).
1213  */
1214 closeto(zeta(3)-1.202056903159594,1e-15), numer:true;
1215 true;
1216 closeto(zeta(3+%i)-(1.10721440843141 - .1482908671781754*%i), 1e-15);
1217 true;
1219  * Reported on mailing list 2011-05-22 by Thomas Dean:
1221  * plot2d(abs(zeta(1/2+x*%i)),[x,0,36]) causes a Lisp error with
1222  * clisp.
1224  */
1225 closeto(abs(zeta(1/2+.5*%i)) - 1.06534921249378, 1e-14);
1226 true;
1228 /* Bug ID: 2997401 - float(log(200!)) produces an error
1230  */
1231 closeto(float(log(200!))-863.2319871924054, 1e-15);
1232 true;
1234 closeto(float(log((1+200!)/7))-861.2860770433501, 1e-15);
1235 true;
1237 /* Additional tests */
1238 closeto(float(log(-1))-float(%pi)*%i, 1e-15);
1239 true;
1241 closeto(float(log((1+200!)/(-7))) - (3.141592653589793*%i + 861.2860770433501), 1e-14);
1242 true;
1244 closeto(float(log((1+200!)+(1+199!)*%i))- (.004999958333958322*%i + 863.2319996922491), 1e-15);
1245 true;
1247 closeto(float(log((1+200!)/7+(1+199!)/11*%i)) - (.003181807444342708*%i + 869.9736929490153), 1e-15);
1248 true;
1250 /* Bug ID: 2306402 - scalarp bug
1251  * Bug ID: 1985748 - array and scalar declarations yield inconsistent results
1252  * Examples from the bug report to show consistent behavior of scalarp
1253  */
1255 declare(x,scalar);
1256 done;
1257 scalarp(foo(x));
1258 true;
1259 scalarp(foo(1));
1260 true;
1261 scalarp(foo(x,1));
1262 true;
1264 scalarp(x);
1265 true;
1266 scalarp(x[1]);
1267 true;
1268 array(x,5);
1270 scalarp(x);
1271 true;
1272 scalarp(x[1]);
1273 true;
1274 nonscalarp(x);
1275 false;
1277 kill(x);
1278 done;
1280 /* Bug ID:1723548 - gradef for variables: not used in diff
1281  * Show that the total differential of f works in expressions too.
1282  */
1283 depends(f,[x,y]);
1284 [f(x,y)];
1285 diff(f);
1286 'diff(f,y,1)*del(y)+'diff(f,x,1)*del(x)$
1287 diff(3*f);
1288 3*'diff(f,y,1)*del(y)+3*'diff(f,x,1)*del(x)$
1289 diff(a*f);
1290 a*'diff(f,y,1)*del(y)+a*'diff(f,x,1)*del(x)+f*del(a)$
1291 remove(f,dependency);
1292 done;
1294 /* Bug ID: 1089719 addrow creates strange matrix
1295  */
1296 m:matrix([0,0]);
1297 matrix([0,0]);
1298 m:addrow(m,m);
1299 matrix([0,0],[0,0])$
1300 m[1,1]:11;
1303 matrix([11,0],[0,0])$
1304 kill(m);
1305 done;
1307 /* Bug ID: 1663385 - declare multiplicative - wrong simplification
1308  */
1309 declare(f,additive,f,multiplicative);
1310 done;
1311 f(a*b+c);
1312 f(a)*f(b)+f(c);
1313 kill(f);
1314 done;
1316 /* Bug ID: 816808 - subst(in)part of rat -- internal errs
1317  */
1318 substpart(x,2/3,2);
1319 2/x;
1320 substinpart(4,2/3,2);
1321 1/2;
1323 /* Bug ID: 1117533 - letsimp complains about assignment to %pi
1324  */
1325 matchdeclare(a,true);
1326 done;
1327 (let(%pi*a,foo(a)),done);
1328 done;
1329 letsimp(%pi*x);
1330 foo(x);
1331 remlet(%pi*a);
1332 done;
1334 /* Bug ID: 2805600 depends() partially prevents diff() to work
1335  */
1336 depends(t,x);
1337 [t(x)];
1338 diff(f(t),z);
1340 remove(t,dependency);
1341 done;
1343 /* Bug ID: 1184718 - AT needs soime basic simplifications
1344  */
1345 'at(2,x=0);
1348 /* Bug ID: 2998227 - spurious at(0,A=0)
1349  */
1350 taylor(integrate(gamma(x+1),x,0,A),A,0,3),nouns;
1351 A-%gamma*A^2/2+(6*%gamma^2+%pi^2)*A^3/36;
1353 /* Bug ID: 3010829 - numerical evaluation of elliptic_ec fails for argument > 1
1354  */
1355 closeto(elliptic_ec(2.0)-(.5990701173677959*%i+0.599070117367796), 1.5e-15);
1356 true;
1358 /* Bug ID: 1929287 - 0.0 + [0] ---> [0]
1359  */
1360 0.0+[0];
1361 [0.0];
1362 0.0b0+[0];
1363 [0.0b0];
1364 0.0+matrix([0,1/2,1,x]);
1365 matrix([0.0,0.5,1.0,x]);
1366 0.0b0+matrix([0,1/2,1,x]);
1367 matrix([0.0b0,5.0b-1,1.0b0,x]);
1369 /* Bug ID: 2996106 - at(diff(f(x,y),x,1,y,1),[x=a,y=b]) is wrong
1370  */
1371 at(diff(f(x,y),x,1,y,1),[x=a,y=b]);
1372 'at('diff(f(x,y),x,1,y,1),[x = a,y = b]);
1374 /* Bug report ID: 2556133 - "at" should do parallel substitutions
1375  */
1376 errcatch(at(atan2(y^2+1,x),[y=%i,x=0]));
1378 errcatch(at(atan2(y^2+1,x),[x=0,y=%i]));
1381 /* Bug report ID: 2014941 - compositions of 'at'
1382  */
1383 at(at(diff(f(x),x),[x=b]),[b=y]);
1384 'at('diff(f(x),x,1),[x = y]);
1386 at(diff(f(x,y),x,1,y,1),[x=a,y=b]) - at(diff(f(x,y),x,1,y,1),[y=b,x=a]);
1389 /* Bug report ID: 1677217 - composistions of 'at'
1390  */
1391 depends(y,[x,z]);
1392 [y(x,z)];
1393 at(at(diff(y,x),x=a),z=b);
1394 'at('at('diff(y,x,1),x = a),z = b);
1395 remove(y,dependency);
1396 done;
1398 /* Bug report ID: 3023978 - integrate(x^x+x,x) is wrong
1399  */
1400 integrate(x^x+x,x);
1401 'integrate(exp(x*log(x)),x)+x^2/2;
1403 /* Bug report ID: 2465066 - unsimplified result from integrate
1404  */
1405 matchdeclare(x, symbolp);
1406 done;
1407 (tellsimpafter('integrate(f(x),x), g(x)),done);
1408 done;
1409 integrate(5*f(x) + 7,x);
1410 5*g(x)+7*x;
1411 kill(rules);
1412 done;
1414 /* Bug report ID: 2789110 - solve, tan and atan depend on order of variables
1415  */
1416 solve(tan(x - atan(a/b)) = 0, x);
1417 [x = atan(a/b)];
1418 solve(tan(x - atan(b/a)) = 0, x);
1419 [x = atan(b/a)];
1421 /* Bug report ID: 1961494 - translated functions & values list
1422  */
1423 (kill(all), f():= x:2, translate(f));
1424 [f];
1425 f();
1427 values;
1428 [x];
1429 kill(f,x);
1430 done;
1431 /* The value of x has been removed. */
1435 /* Bug report ID: 3025038 - gruntz needs logexpand:true
1436  */
1437 gruntz( (x + 2^x) / 3^x, x, inf),logexpand:false;
1440 /* Bug report ID: 2977217 - maxima can not integrate x*exp(-1/2*(x-m)^2)
1441  */
1442 integrate(x*exp(-1/2*(x-m)^2),x);
1443 %i*(sqrt(2)*%i*gamma_incomplete(1,(m-x)^2/2)*(m-x)^2/(x-m)^2
1444    -%i*gamma_incomplete(1/2,(m-x)^2/2)*m*(m-x)/abs(x-m))/sqrt(2);
1446 /* Bug report ID: 2996542 - log(x) integration is incorrect
1447  */
1448 assume(a>0);
1449 [a>0];
1450 integrate(log(x),x,0,a);
1451 a*log(a)-a;
1452 forget(a>0);
1453 [a>0];
1455 /* Bug report ID: 3062883 - diff does not recognize indirect dependencies 
1456  *                          in expressions
1457  */
1459 depends([a,b],x,x,t);
1460 [a(x),b(x),x(t)];
1461 diff(-a,t);
1462 -'diff(a,x,1)*'diff(x,t,1);
1463 diff(a*b,t);
1464 a*'diff(b,x,1)*'diff(x,t,1)+'diff(a,x,1)*b*'diff(x,t,1);
1465 remove([a,b,x],dependency);
1466 done;
1468 /* Bug report ID: 3080397 - laplace(unit_step(-t),t,s) generates an error.
1469  */
1470 laplace(unit_step(-t),t,s);
1473 /* Bug report ID: 3081820 - lbfgs causes error
1475  * Still generates an error, but a different error that maxima
1476  * signals.
1477  */
1478 block([V:0.75, a:24, b:68, e],
1479   C(r) := 2*%pi*b*r^2 + 4*a*%pi*r + 2*b*V/r + a*V/(%pi*r^2),
1480   load(lbfgs),
1481   /* This should signal an error that we catch */
1482   e : errcatch(lbfgs(C(r), [r], [1], 1e-4, [1,0])),
1483   [e, error]);
1484 [[], ["Evaluation of gradient at ~M failed.  Bad initial point?~%", [0.0]]];
1486 /* Bug report ID: 875089 - defint(f(x)=g(x),x,0,1) -> false = false
1488  * We distribute defint more early in the code of bags to get a correct result.
1489  */
1490  defint(f(x)=g(x),x,0,1);
1491  'integrate(f(x),x,0,1)='integrate(g(x),x,0,1);
1493 /* Bug report ID: 2796194 - error doing a Fourier transform */
1494 (assume(equal(x,0)),done);
1495 done;
1497 errcatch(integrate(%pi*exp(-2*%pi*t)*exp(2*%pi*x*t*%i),t,minf,inf));
1500 error;
1501 ["defint: integral is divergent."];
1503 (forget(equal(x,0)),done);
1504 done;
1506 /* Bug reported on the mailing list
1507  * <http://www.math.utexas.edu/pipermail/maxima/2010/023024.html>
1508  * integrate(cos(2*x)*cos(x),x) is wrong.
1510  * Add a few more test that are similar to test the part of
1511  * monstertrig that deals with trig1(m*x)*trig2(n*x) where trig1 and
1512  * trig2 are either sin or cos.
1513  */
1514 integrate(cos(2*x)*cos(x),x);
1515 sin(3*x)/6+sin(x)/2;
1517 integrate(sin(2*x)*sin(x),x);
1518 sin(x)/2-sin(3*x)/6;
1520 integrate(cos(2*x)*sin(x),x);
1521 cos(x)/2-cos(3*x)/6;
1523 integrate(sin(x)*cos(2*x),x);
1524 cos(x)/2-cos(3*x)/6;
1526 /* Bug ID: 3111568 - subsequent calls to gradef hide variable dependencies
1527  */
1528 gradef(f,x,g);
1530 gradef(f,y,h);
1532 dependencies;
1533 [f(y,x)];
1534 kill(f);
1535 done;
1537 /* Bug ID: 3118770 - %edispflag:true causes a bug
1538  */
1539 %edispflag:true;
1540 true;
1541 integrate(x/(%e)^(2*x), x, 0, 1);
1542 1/4-3/(4*%e^2);
1543 reset(%edispflag);
1544 [%edispflag];
1546 /* Bug ID: 3067098 - The command timer for a Lisp function
1547  * Check that the function trisplit does not go away, when we collect
1548  * timing statistics for this function and call later kill(all).
1549  */
1550 timer(?trisplit);
1551 [?trisplit];
1552 kill(all);
1553 done;
1554 rectform(1+%i);
1555 1+%i;
1557 /* Bug ID: 3133916 - scanmap(minfactorial,a!) infinite loop
1558  */
1559 scanmap(minfactorial, a!);
1562 /* Bug ID: 3131324 - simplification of sqrt
1563  */
1564 sqrt(x^3)/sqrt(x^3);
1567 /* Bug ID: 1285104 - trigsimp and trigreduce & square roots
1568  */
1570 trigreduce(sqrt(r^2*sin(x)^2+r^2*cos(x)^2));
1571 abs(r);
1572 trigreduce(sqrt(r^2*sin(x)^2+r^2*cos(x)^2)),radexpand:all;
1574 radexpand:false;
1575 false;
1576 trigreduce(sqrt(r^2*sin(x)^2+r^2*cos(x)^2));
1577 sqrt(r^2);
1578 reset(radexpand);
1579 [radexpand];
1581 /* Bug ID: 917283 - Comment syntax confused
1582  * Show that nested comments work as expected.
1583  */
1584 a/*/**/*/+b;
1585 a+b;
1586 a/*/**/*/+/*/**/*/b;
1587 a+b;
1589 /* Bug ID: 3138054 -  bfloat problem / FIX -
1590  */
1591 exp(gamma(1/3)),bfloat;
1592 1.45696199392313b1;
1594 /* Bug ID: 3288989 - Lisp functions and linear display
1595  * Show that we do not get a Lisp error.
1596  */
1597 grind(?cdr([a,b,c]));
1598 done;
1600 /* Bug ID: 3291590 - Problems with fast arrays
1601  */
1602 (a:make_array(hashed), done);
1603 done;
1604 a[100]:100;
1605 100;
1606 a[x]:sin(x);
1607 sin(x);
1608 a[x*y]:x^2+y;
1609 y+x^2;
1611 /* Cutting out these two examples.
1612  * The ordering of the lists is different depending on the underlying Lisp.
1614  * arrayinfo(a);
1615  * [hash_table,1,100,x,x*y];
1616  * listarray(a);
1617  * [100,sin(x),y+x^2];
1618  */
1620 (f:make_array(functional, 'factorial, hashed), done);
1621 done;
1622 f[10];
1623 3628800;
1625 (kill(f), 0);
1628 (a: make_array(fixnum, 2, 2), done);
1629 done;
1630 listarray(a);
1631 [0, 0, 0, 0];
1633 use_fast_arrays:true;
1634 true;
1636 (array(a, any, 2, 2), done);
1637 done;
1638 arrayinfo(a);
1639 [declared, 2, [2, 2]];
1640 (array(a, fixnum, 2, 2), done);
1641 done;
1642 arrayinfo(a);
1643 [declared, 2, [2, 2]];
1644 (array(a, flonum, 2, 2), done);
1645 done;
1646 arrayinfo(a);
1647 [declared, 2, [2, 2]];
1648 (array(a, hashed), done);
1649 done;
1650 arrayinfo(a);
1651 [hash_table, 1];
1653 reset(use_fast_arrays);
1654 [use_fast_arrays];
1655 kill(a);
1656 done;
1658 /* Bug ID: 3247367 - expand returns unsimplified
1659  */
1660 sqrt(2)+sqrt(2);
1661 2^(3/2);
1662 sqrt(2)+sqrt(2)+sqrt(2);
1663 3*sqrt(2);
1664 sqrt(2)+sqrt(2)+sqrt(2)+sqrt(2);
1665 2^(5/2);
1666 sqrt(2)+sqrt(2)+sqrt(2)+sqrt(2)+sqrt(2);
1667 5*sqrt(2);
1669 2*sqrt(2)+3*sqrt(2);
1670 5*sqrt(2);
1671 3*sqrt(2)+2*sqrt(2);
1672 5*sqrt(2);
1673 3*sqrt(2)+2*sqrt(2)+sqrt(2);
1674 3*2^(3/2);
1675 sqrt(2)+3*sqrt(2)+2*sqrt(2);
1676 3*2^(3/2);
1678 sqrt(1/2)+sqrt(1/2);
1679 sqrt(2);
1680 sqrt(1/2)+sqrt(1/2)+sqrt(1/2);
1681 3/sqrt(2);
1682 sqrt(1/2)+sqrt(1/2)+sqrt(1/2)+sqrt(1/2);
1683 2^(3/2);
1684 sqrt(1/2)+sqrt(1/2)+sqrt(1/2)+sqrt(1/2)+sqrt(1/2);
1685 5/sqrt(2);
1687 2*sqrt(1/2)+3*sqrt(1/2);
1688 5/sqrt(2);
1689 3*sqrt(1/2)+2*sqrt(1/2);
1690 5/sqrt(2);
1691 3*sqrt(1/2)+2*sqrt(1/2)+sqrt(1/2);
1692 3*sqrt(2);
1693 sqrt(1/2)+3*sqrt(1/2)+2*sqrt(1/2);
1694 3*sqrt(2);
1696 2^(1/5)+2^(1/5);
1697 2^(6/5);
1698 2^(1/5)+2^(1/5)+2^(1/5);
1699 3*2^(1/5);
1700 2^(1/5)+2^(1/5)+2^(1/5)+2^(1/5);
1701 2^(11/5);
1702 2^(1/5)+2^(1/5)+2^(1/5)+2^(1/5)+2^(1/5);
1703 5*2^(1/5);
1705 (1/2)^(1/5)+(1/2)^(1/5);
1706 2^(4/5);
1707 (1/2)^(1/5)+(1/2)^(1/5)+(1/2)^(1/5);
1708 3/2^(1/5);
1709 (1/2)^(1/5)+(1/2)^(1/5)+(1/2)^(1/5)+(1/2)^(1/5);
1710 2^(9/5);
1711 (1/2)^(1/5)+(1/2)^(1/5)+(1/2)^(1/5)+(1/2)^(1/5)+(1/2)^(1/5);
1712 5/2^(1/5);
1714 2*(1/2)^(1/5)+3*(1/2)^(1/5);
1715 5/2^(1/5);
1716 3*(1/2)^(1/5)+2*(1/2)^(1/5);
1717 5/2^(1/5);
1719 2^sin(x)+2^sin(x);
1720 2^(sin(x)+1);
1721 2^sin(x)+2^sin(x)+2^sin(x);
1722 3*2^sin(x);
1723 2^sin(x)+2^sin(x)+2^sin(x)+2^sin(x);
1724 2^(sin(x)+2);
1725 2^sin(x)+2^sin(x)+2^sin(x)+2^sin(x)+2^sin(x);
1726 5*2^sin(x);
1728 (1/2)^sin(x)+(1/2)^sin(x);
1729 2^(1-sin(x));
1730 (1/2)^sin(x)+(1/2)^sin(x)+(1/2)^sin(x);
1731 3/2^sin(x);
1732 (1/2)^sin(x)+(1/2)^sin(x)+(1/2)^sin(x)+(1/2)^sin(x);
1733 2^(2-sin(x));
1734 (1/2)^sin(x)+(1/2)^sin(x)+(1/2)^sin(x)+(1/2)^sin(x)+(1/2)^sin(x);
1735 5/2^sin(x);
1737 (1-sqrt(5))^3-4*(1-sqrt(5))^2+8, expand;
1739 1/sqrt(2)+1/sqrt(2)+1/sqrt(2);
1740 3/sqrt(2);
1741 2^(9/5)+2^(4/5);
1742 3*2^(4/5);
1743 3*sqrt(2)+2*sqrt(2);
1744 5*sqrt(2);
1745 2*sqrt(2)+3*sqrt(2);
1746 5*sqrt(2);
1747 (1-sqrt(5))^3, expand;
1748 16-8*sqrt(5);
1749 p : z^3-2^(3/2)*%i*z^2-4*z^2+2^(5/2)*%i*z+2*z;
1750 z^3-2^(3/2)*%i*z^2-4*z^2+2^(5/2)*%i*z+2*z;
1751 divide(p, (z-2-sqrt(2)*%i),z);
1752 [z^2+(-sqrt(2)*%i-2)*z,0];
1753 2^a + 3*2^(a+1);
1754 7*2^a;
1755 2^(3/5)+2^(-2/5);
1756 3/2^(2/5);
1757 2^(3/5+x)+2^(-2/5+x);
1758 3*2^(x-2/5);
1760 /* -----------------------------------------------------------------------------
1761  * Bug ID: 1439566 - zerobern & bernpoly
1762  * Show that the option variable zerobern does not change the results.
1763  * -------------------------------------------------------------------------- */
1764 zerobern:false;
1765 false;
1766 bernpoly(x,3);
1767 x^3-3*x^2/2+x/2$
1768 bernpoly(x,5);
1769 x^5-5*x^4/2+5*x^3/3-x/6$
1770 reset(zerobern);
1771 [zerobern];
1773 /* Show that bern no longer fails with zerobern:false.
1775  * The compared values are from A&S p810, Table 23.2.
1776  */
1777 bern(28);
1778 -23749461029/870$
1779 bern(40);
1780 -261082718496449122051/13530$
1781 euler(16);
1782 19391512145$
1784 zerobern:false;
1785 false$
1786 bern(17);
1787 -7709321041217/510$
1788 bern(24);
1789 596451111593912163277961/282$
1790 euler(10);
1791 370371188237525$
1793 reset(zerobern);
1794 [zerobern]$
1796 /* -----------------------------------------------------------------------------
1797  * Bug ID: 2905929 - gcdex
1798  * -------------------------------------------------------------------------- */
1799 q0[2] : 6;
1801 ratdisrep(gcdex(x-7, x-8));
1802 [1,-1,1]$
1804 is(equal(gcdex(z^2-1, 0, z), [1,0,z^2-1]));
1805 true;
1807 is(equal(gcdex(0, z^2-1, z), [0, 1, z^2-1]));
1808 true;
1810 /* Examples from the Maxima Manual */
1811 is(equal(gcdex(x^2+1, x^3+4), [-(x^2+4*x-1)/17,(x+4)/17,1]));
1812 true$
1813 is(equal(gcdex(x*(y+1), y^2-1, x), [0,1/(y^2-1),1]));
1814 true$
1816 kill(q0);
1817 done;
1819 /* -----------------------------------------------------------------------------
1820  * Bug ID: 3389830 - Error break in rtest15 with linear display
1822  * Show that we do not get an error with grind for a prefix and a postfix
1823  * Operator when displaying the definition of a function or a macro.
1824  * This is a test of the function msz-mdef in grind.lisp.
1825  * -------------------------------------------------------------------------- */
1826 (postfix("f"), prefix("g"), done);
1827 done$
1828 grind("f"(x):=x);
1829 done$
1830 grind("f"(x)::=x);
1831 done$
1832 grind("g"(x):=x);
1833 done$
1834 grind("g"(x)::=x);
1835 done$
1837 kill ("f", "g");
1838 done;
1840 /* -----------------------------------------------------------------------------
1841  * Bug ID: 3396631 - equal terms produce different results
1842  * Correcting a bug in plusin revision 23.08.2011
1843  * -------------------------------------------------------------------------- */
1844 5*sqrt(5)+2*sqrt(3)+6*sqrt(5);
1845 11*sqrt(5)+2*sqrt(3)$
1846 5*sqrt(5)+4*sqrt(3)+6*sqrt(5);
1847 11*sqrt(5)+4*sqrt(3)$
1848 5*sqrt(5)+3*sqrt(5)+5*sqrt(3)+3*sqrt(3)+2*sqrt(75)+2*sqrt(45);
1849 14*sqrt(5)+2*3^(5/2)$
1850 5*sqrt(5)+5*sqrt(3)+3*sqrt(3)+2*sqrt(75)+2*sqrt(45)+3*sqrt(5);
1851 14*sqrt(5)+2*3^(5/2)$
1852 5*sqrt(5)+5*sqrt(3)+2*sqrt(75)+2*sqrt(45)+3*sqrt(5)+3*sqrt(3);
1853 14*sqrt(5)+2*3^(5/2)$
1855 /* -----------------------------------------------------------------------------
1856  * Bug ID 3437268: expand doesn't fully expand
1857  * Check the expected results after revision 14.11.2011 of simp.lisp.
1858  * -------------------------------------------------------------------------- */
1859 -3*a/sqrt(2)+sqrt(2)*a+a/sqrt(2);
1861 expand(-(atan((sqrt(2)*sin(8))/cos(8))+3*%pi)/sqrt(2)
1862         +atan((sqrt(2)*sin(8))/cos(8))/sqrt(2)
1863         +sqrt(2)*%pi
1864         +%pi/sqrt(2));
1868  * SF Bug: elliptic_e error observed - ID: 3438911
1870  * Test the bug by using the identity:
1871  * elliptic_e(%pi,m) = 2*elliptic_ec(m)
1873  * There's also a bug in elliptic_f.  We test this using the identity:
1874  * elliptic_f(%pi,m) = 2*elliptic_kc(m)
1875  */
1877 closeto(e,tol):=block([numer:true,abse],abse:abs(e),if(abse<tol) then true else abse);
1878 closeto(e,tol):=block([numer:true,abse],abse:abs(e),if(abse<tol) then true else abse);
1880 closeto(elliptic_e(float(%pi), .1) - 2*elliptic_ec(.1), 1e-15);
1881 true;
1882 closeto(elliptic_e(float(%pi), .9) - 2*elliptic_ec(.9), 1e-15);
1883 true;
1884 closeto(elliptic_e(bfloat(%pi), .1b0) - 2*elliptic_ec(.1b0), 1e-16);
1885 true;
1886 closeto(elliptic_e(bfloat(%pi), .9b0) - 2*elliptic_ec(.9b0), 1e-16);
1887 true;
1889 closeto(elliptic_f(float(%pi), .1) - 2*elliptic_kc(.1), 1e-15);
1890 true;
1891 closeto(elliptic_f(float(%pi), .9) - 2*elliptic_kc(.9), 1e-15);
1892 true;
1893 closeto(elliptic_f(bfloat(%pi), .1b0) - 2*elliptic_kc(.1b0), 1e-16);
1894 true;
1895 closeto(elliptic_f(bfloat(%pi), .9b0) - 2*elliptic_kc(.9b0), 1e-16);
1896 true;
1899  * Bug 3526111 - float erf (%i) not working
1900  */
1902 closeto(float(erf(%i)) - 1.650425758797543*%i, 1e-15);
1903 true;
1906  * Bug 3529992: Shi (sinh integral) wrong branch, integrate inconsistent
1907  */
1908 closeto(float(expintegral_shi(1/2) - 0.50699674981966719583), 3e-16);
1909 true;
1911 /* integrate changes k[0] --> k(0) - ID: 3530767 */
1912 integrate(x * (x^2 + k[0])/(1 + x^2),x);
1913 ((k[0]-1)*log(x^2+1))/2+x^2/2$
1915 /* polarform error on simple case - ID: 3517034 */
1916 polarform((a+1)/2 - a/2 - 1/2);
1919 polarform((a+1)/2 - a/2 - 0.5);
1920 0.0$
1922 polarform((a+1)/2 - a/2 - 0.5b0);
1923 0.0b0$
1926 /* #2531 Integration with inf */
1927 errcatch(integrate((1+1/x)^(1/2),x,1,inf));
1930 error;
1931 ["defint: integral is divergent."];
1936  * ID: 3440046: elliptic_f(0.5,1) signals error
1938  * Add a few more tests for invalid values.
1939  */
1940 closeto(elliptic_f(0.5,1)-elliptic_f(1/2,1), 1e-15);
1941 true;
1943 closeto(elliptic_f(0.5b0,1) - bfloat(elliptic_f(1/2,1)), 1b-16);
1944 true;
1946 errcatch(elliptic_f(2.0,1));
1949 errcatch(elliptic_f(2b0,1));
1953  * Bug 3428734:  integrate(bessel_y(1,z),z)  with ?z : 1
1954  */
1956 (?z:1, 0);
1959 integrate(bessel_y(1,z),z);
1960 -bessel_y(0,z);
1962 integrate(bessel_j(1,z),z);
1963 -bessel_j(0,z);
1965 integrate(bessel_k(1,z),z);
1966 -bessel_k(0,z);
1968 integrate(bessel_i(1,z),z);
1969 bessel_i(0,z);
1972  * Bug 3381301: log(-1.0b0) has small realpart
1973  */
1974 realpart(log(-1b0));
1978  * Bug 3559064: elliptic_f(2,1) is wrong.
1979  */
1980 errcatch(elliptic_f(2,1));
1984  * Bug 2528: A variable should be real if it is both real and complex
1985  */
1986 (declare(foo, real), declare(foo, complex), 0);
1989 [realpart(foo), imagpart(foo)];
1990 [foo, 0]$
1992 (kill(foo), 0);
1995 map(lambda([x], featurep(x, 'irrational)),[42,%pi,%phi,%e]);
1996 [false, true, true, true]$
1998 /* #2501 %pi/8 is definitely not an integer */
1999 integrate(log(cot(x)-1),x,0,%pi/4);
2000 (%i*li[2]((%i+1)/2)-%i*li[2](-((%i-1)/2)))/2 -(%i*(2*li[2](%i+1)-2*li[2](1-%i))+%pi*log(2))/4$
2002 integrate(log(cos(x)),x,0,%pi/2);
2003 %i*%pi^2/24-(6*%pi*log(4)+%i*%pi^2)/24$
2005 map(lambda([x], featurep(x, noninteger)),[sqrt(5),%pi,%pi/3, %pi/8,log(42),99/2013]);
2006 [true,true,true,true,true,true]$
2008 map(lambda([x], featurep(x, noninteger)),[0,1,2013]);
2009 [false,false,false]$
2011 /* #2583 sign error for integrate(x^(8*%i-1),x) */
2012 block([domain : 'real], integrate(x^(8*%i-1),x));
2013 -%i*x^(8*%i)/8$
2015 /* #2602: some-bfloatp and some-floatp recursed wrongly on rat expressions */
2016 ?some\-bfloatp(rat(1/2));
2017 false$
2019 /* #2594: Error in trigreduce for complicated expressions */
2020 subst(0, x, trigreduce(product(cos(k*x), k, 1, 8)));
2023 /* SF bug #2818: Problem with trigreduce */
2024 trigreduce(sin(1/8*%pi)*sin(3/8*%pi)*sin(5/8*%pi)*sin(7/8*%pi));
2025 1/8;
2027 /* #2591: Risch gives lisp error */
2028 (risch(asinh((z^2-1)/z)/z,z), 0);
2031 /* #2682: Function zeta fails numerically for large numbers */
2032 closeto(zeta(40.0) - 1, 1e-12);
2033 true$
2035 /* #2675 (1/3): Integration yields noun form with subscripted variable */
2036 integrate(exp(-(1+%i)*x[1]),x[1],0,inf);
2037 1/2-%i/2$
2039 /* #2688: %e^^A returns element by element exponent */
2040 is(%e^^matrix([1,2],[3,4]) = %e^matrix([1,2],[3,4]));
2041 false$
2043 /* #2676: Integral incorrect when variable is subscripted */
2044 integrate(x[1]*exp(x[1]), x[1]);
2045 exp(x[1])*(x[1]-1)$
2047 /* #2726: Integrate produces wrong answer for Gaussian Moments */
2048 (declare(m2726, even),
2049  block([tmp: integrate(exp(-x^2/2)/sqrt(2*%pi) * x^m2726, x, -1/4, 1/4)],
2050    sign (subst(m2726 = 4, tmp))));
2051 pos$
2053 /* # 2697: Inconsistent handling of Greek symbols */
2054 integrate(y(%theta)=sin(%theta),%theta,%theta[0], %theta[1]);
2055 'integrate(y(%theta),%theta,%theta[0],%theta[1]) = cos(%theta[0])-cos(%theta[1])$
2057 integrate(y(t)=sin(t),t,t[0], t[1]);
2058 'integrate(y(t),t,t[0],t[1]) = cos(t[0])-cos(t[1])$
2060 integrate(y(tau)=sin(tau),tau,tau[0], tau[1]);
2061 'integrate(y(tau),tau,tau[0],tau[1]) = cos(tau[0])-cos(tau[1])$
2063 integrate ([foo(x), bar(x)], x, x[1], x[2]);
2064 ['integrate (foo(x), x, x[1], x[2]), 'integrate (bar(x), x, x[1], x[2])];
2066 /* # 2738: Integrate encountered a Lisp error: The value 2 is not of type LIST. */
2068 (kill (x, y, I, J),
2069  x(t):=2*cos(t),
2070  y(t):=2*sin(t),
2071  I : (x(t)+y(t)^2)*sqrt(diff(x(t),t)^2+diff(y(t),t)^2),
2072  J : integrate (I, t));
2073 4*(t-tan(t)/(tan(t)^2+1))+4*sin(t)$
2075 trigsimp (diff (J, t) - I);
2077 /* bug #2980: Infinite recursion with (e: log(e), rectform(e)) */
2078 block ([e: log(e)], rectform(e));
2079 log(abs(e)) + %i*atan2(0, e)$
2081 /* bug #2159: integration_with_logabs ("integrate(tan(x),x);" etc. do not take "logabs" flag into account) */
2082 integrate([tan(x),csc(x),sec(x),cot(x),tanh(x),coth(x),csch(x)],x);
2083 [log(sec(x)),-log(csc(x)+cot(x)),log(sec(x)+tan(x)),log(sin(x)),log(cosh(x)),log(sinh(x)),log(tanh(x/2))]$
2084 integrate([tan(x),csc(x),sec(x),cot(x),tanh(x),coth(x),csch(x)],x),logabs;
2085 [log(abs(sec(x))),-log(abs(csc(x)+cot(x))),log(abs(sec(x)+tan(x))),log(abs(sin(x))),log(cosh(x)),log(abs(sinh(x))),log(abs(tanh(x/2)))]$
2087 /* Bug #3075: #3075 answer "3*false" from "integrate(3*asinh(x),x,-inf,inf)" */
2088 /* We can't check for the correct answer zero (yet), because Maxima can't solve integrate(asinh(x),x,-inf,inf). */
2089 is(integrate(3*asinh(x),x,-inf,inf)#3*false);
2090 true$
2093  * Bug #3056: exp(1b19) signals error that 1b19 doesn't have enough
2094  * precision to compute its integer part.  Add test for this and also
2095  * the original test from commit 576c7508.)
2096  */
2098 /* Don't really care what the answer is as long as we don't signal an error */
2099 is(errcatch(exp(1b19)) # []);
2100 true;
2102 /* Test from the commit 576c7508 */
2103 ceiling((207300647060*%e^-563501581931)/(403978495031*%e^-1098127402131));
2104 ceiling((207300647060*%e^534625820200)/403978495031);
2106 /* Bug #3098: numerical evaluation of li[3] */
2107 li[3](0.0);
2108 0.0;
2110 closeto(li[3](1.0) - zeta(3.0), 1e-15);
2111 true;
2113 closeto(li[3](0.5) - li[3](1/2), 1e-15);
2114 true;
2116 closeto(li[3](-2.0) - float(subst(z=-2.0, li[3](1/z) - log(-z)/6*(log(-z)^2+%pi^2))), 1e-15);
2117 true;
2119 closeto(abs(li[3](2.0) - expand(float(subst(z=2.0, li[3](1/z) - log(-z)/6*(log(-z)^2+%pi^2))))), 1e-15);
2120 true;
2122 li[2](-1);
2123 -%pi^2/12;
2125 li[2](1);
2126 %pi^2/6;
2128 li[2](1/2);
2129 %pi^2/12 - log(2)^2/2;
2131 closeto(li[2](0.5) - (%pi^2/12 - log(2)^2/2), 1.1103e-16), numer;
2132 true;
2134 closeto(li[2](2.0) - (%pi^2/4 - %i*%pi*log(2)), 2.513e-15), numer;
2135 true;
2137 closeto(li[2]((1-sqrt(5))/2) - (log((sqrt(5)-1)/2)^2/2-%pi^2/15), 1.1103e-16), numer;
2138 true;
2140 /* Catalan's constant: 0.915965594... */
2141 closeto(li[2](1.0*%i) - (-%pi^2/48 + %i*0.915965594177219015054603514932384110774149374281672134266), 1.3878e-15), numer;
2142 true;
2144 closeto(li[2](1.0-%i) - (%pi^2/16-%i*0.915965594177219015054603514932384110774149374281672134266 - %pi*%i*log(2)/4), 4.5e-16), numer;
2145 true;
2147 /* Make sure li[3](1/7),numer returns a float and not a bfloat */
2148 ?floatp(li[3](float(1/7)));
2149 true;
2151 ?floatp(realpart(li[2](1.0+1.0*%i)));
2152 true;
2154 closeto(li[3](0.5) - float((7*zeta(3))/8+log(2)^3/6-(%pi^2*log(2))/12), 1.1103e-16);
2155 true;
2157 closeto(float(li[3](exp(%pi*%i/3))) - float(zeta(3)/3 + %i*5*%pi^3/162), 5.6611e-16);
2158 true;
2161  * li[3] should not return a small imaginary part for this value.
2162  * (There are others, but the solution fixes those as well.
2163  */
2164 imagpart(li[3](-.862));
2168  * Likewise li[s](-1.5) was returning small imaginary parts.  In fact,
2169  * this was true for -2 < z < -1 because the log series was used to
2170  * compute the value.  Check that we return exactly zero now for
2171  * select values of s
2172  */
2174 makelist(imagpart(li[s](-1.5)), s, 4, 10);
2175 [0, 0, 0, 0, 0, 0, 0];  
2177 /* Bug 3582: numerical evaluation of li[2] */
2178 /* this gave an infinite loop */
2179 closeto(li[2]( 0.5 + %i* 0.6) - ( 0.7535498331267791 * %i + 0.4081870726952078 ), 1e-15);
2180 true;
2182 /* Bug 3112: zeta(n) for negative even n is inaccurate */
2183 zeta(-4.0);
2184 0.0;
2186 zeta(-6b0);
2187 0b0;
2189 /* Bug 3105: li[s](1.0) doesn't simplify */
2190 closeto(li[4](1.0)-li[4](1), 1.11022e-16);
2191 true;
2193 closeto(li[4](-1.0)-li[4](-1), 1.11022e-16);
2194 true;
2196 closeto(li[4](1b0)-li[4](1), 10^-fpprec);
2197 true;
2199 closeto(li[4](-1b0)-li[4](-1), 10^-fpprec);
2200 true;
2204  * More numeric test for polylog function li[s](z).
2206  * The reference values were obtained from
2207  * http://functions.wolfram.com/webMathematica/FunctionEvaluation.jsp?name=PolyLog
2208  */
2209 fpprec:24;
2212 closeto(li[4](0.25) - 0.25411619074634353405967371315352,
2213   5.5512e-17);
2214 true;
2216 closeto(li[4](0.25b0) - 0.25411619074634353405967371315352b0,
2217   2.06796b-25);
2218 true;
2220 closeto(li[4](0.75) - 0.79222102797282777952948578955736,
2221   2.2205e-16);
2222 true;
2224 closeto(li[4](0.75b0) - 0.79222102797282777952948578955736b0,
2225   6.20386b-25);
2226 true;
2228 closeto(li[4](1.5) - (1.7347570807760620737768805117515 - 0.0349027048283367002627421237287*%i),
2229   2.2205e-16);
2230 true;
2232 closeto(li[4](1.5b0) - (1.7347570807760620737768805117515b0 - 0.0349027048283367002627421237287b0*%i),
2233   4.14398b-25);
2234 true;
2236 closeto(li[4](10.0) - (9.6140263862742968515251940747859 - 6.3921313179656069159740055708257*%i),
2237   5.618e-15);
2238 true;
2240 closeto(li[4](10.0b0) - (9.6140263862742968515251940747859b0 - 6.3921313179656069159740055708257b0*%i),
2241   1.9364b-23);
2242 true;
2244 closeto(li[4](-5.0) - (-4.1064679790949702621073505378164),
2245   8.8818e-16);
2246 true;
2248 closeto(li[4](-5b0) - (-4.1064679790949702621073505378164b0),
2249   4.9631b-24);
2250 true;
2252 closeto(li[4](2.0*%i) - (-0.2113943747614829764326347997923 + 1.9289340331586646502356787803360*%i),
2253   7.22e-16);
2254 true;
2256 closeto(li[4](2b0*%i) - (-0.2113943747614829764326347997923b0 + 1.9289340331586646502356787803360b0*%i),
2257   1.9788b-24);
2258 true;
2260 closeto(li[5](0.25) - 0.25202158817857420100669519623555,
2261   5.55112e-17);
2262 true;
2264 closeto(li[5](0.25b0) - 0.25202158817857420100669519623555b0,
2265   5.16988b-25);
2266 true;
2268 closeto(li[5](0.75) - 0.76973541059975738097269173152535,
2269   1.111e-16);
2270 true;
2272 closeto(li[5](0.75b0) - 0.76973541059975738097269173152535b0,
2273   2.06796b-25);
2274 true;
2276 closeto(li[5](1.5) - (1.5961739456813534102689069143338 - 0.0035379572466222227823786042761*%i),
2277   4.33681e-19);
2278 true;
2280 closeto(li[5](1.5b0) - (1.5961739456813534102689069143338b0 - 0.0035379572466222227823786042761b0*%i),
2281   4.13594b-25);
2282 true;
2284 closeto(li[5](10.0) - (11.2390407376112991620107110964539 - 3.6796065713019972004384472107791*%i),
2285   9.93014e-15);
2286 true;
2288 closeto(li[5](10.0b0) - (11.2390407376112991620107110964539b0 - 3.6796065713019972004384472107791b0*%i),
2289   2.18382b-23);
2290 true;
2292 closeto(li[5](-5.0) - (-4.4800824065112010228046981292686),
2293   1e-16);
2294 true;
2296 closeto(li[5](-5b0) - (-4.4800824065112010228046981292686b0),
2297   6.61745b-24);
2298 true;
2300 closeto(li[5](2.0*%i) - (-0.1139660114783041974283299769126 + 1.9734617121122917650272273558071*%i),
2301   6.98e-16);
2302 true;
2304 closeto(li[5](2b0*%i) - (-0.1139660114783041974283299769126b0 + 1.9734617121122917650272273558071b0*%i),
2305   6.42086b-25);
2306 true;
2308 closeto(li[20](5.0) - (5.0000238783147176199891129814336 - 2.182054400942151422e-13*%i),
2309   2.41797e-15);
2310 true;
2312 closeto(li[20](5b0) - (5.0000238783147176199891129814336b0 - 2.182054400942151422b-13*%i),
2313   6.76366b-24);
2314 true;
2317  * Bug 3966: li[s](1) = zeta(s)
2318  */
2319 makelist(li[s+1/2](1), s, 1, 5);
2320 [zeta(3/2), zeta(5/2), zeta(7/2), zeta(9/2), zeta(11/2)];
2322 li[7/3](1);
2323 zeta(7/3);
2325 li[-5/2](1);
2326 li[-5/2](1);
2328 li[-5/2](1.0);
2329 li[-5/2](1.0);
2331 /* Bug #3194: No simplification of "tan(x+n*%pi)" and "cot(x+n*%pi)" with "n" being a declared integer */
2332 /* Also test that the "modulus 1/2" check in "%piargs-tan/cot" works correctly. */
2333 declare(i, integer, ei, even, oi, odd);
2334 done$
2335 [tan(x+i*%pi),cot(x+i*%pi)];
2336 [tan(x),cot(x)]$
2337 [tan(x+ei*%pi),cot(x+ei*%pi)];
2338 [tan(x),cot(x)]$
2339 [tan(x+oi*%pi),cot(x+oi*%pi)];
2340 [tan(x),cot(x)]$
2341 [tan(x-3/2*%pi),tan(x-1/2*%pi),tan(x+1/2*%pi),tan(x+3/2*%pi)];
2342 [-cot(x),-cot(x),-cot(x),-cot(x)]$
2343 [cot(x-3/2*%pi),cot(x-1/2*%pi),cot(x+1/2*%pi),cot(x+3/2*%pi)];
2344 [-tan(x),-tan(x),-tan(x),-tan(x)]$
2345 kill(i, ei, oi);
2346 done$
2348 /* Bug #3148: sign can't figure out sign(a - b) but it knows sign(b - a) where a and b are exponentials */
2349 assume(t>0);
2350 [t>0]$
2351 sign(2^(500000*t)-2^(500007*t));
2352 neg$
2353 sign(2^(500007*t)-2^(500000*t));
2354 pos$
2355 forget(t>0);
2356 [t>0]$
2358 /* Bug #3246: Integrating u'(x) * f(u(x) + c) fails for f any inverse trigonometric/hyperbolic function */
2359 /* Also make sure that all the antiderivatives of (inverse) trigonometric/hyperbolic functions are correct. */
2360 integrate(map(lambda([f], diff(u(x), x) * f(u(x) + c)), [asin,acos,atan,acsc,asec,acot,asinh,acosh,atanh,acsch,asech,acoth]), x);
2361 [(u(x)+c)*asin(u(x)+c)+sqrt(1-(u(x)+c)^2),(u(x)+c)*acos(u(x)+c)-sqrt(1-(u(x)+c)^2),(u(x)+c)*atan(u(x)+c)-log((u(x)+c)^2+1)/2,log(sqrt(1-1/(u(x)+c)^2)+1)/2-log(1-sqrt(1-1/(u(x)+c)^2))/2+(u(x)+c)*acsc(u(x)+c),-log(sqrt(1-1/(u(x)+c)^2)+1)/2+log(1-sqrt(1-1/(u(x)+c)^2))/2+(u(x)+c)*asec(u(x)+c),log((u(x)+c)^2+1)/2+(u(x)+c)*acot(u(x)+c),(u(x)+c)*asinh(u(x)+c)-sqrt((u(x)+c)^2+1),(u(x)+c)*acosh(u(x)+c)-sqrt((u(x)+c)^2-1),log(1-(u(x)+c)^2)/2+(u(x)+c)*atanh(u(x)+c),log(sqrt(1/(u(x)+c)^2+1)+1)/2-log(sqrt(1/(u(x)+c)^2+1)-1)/2+(u(x)+c)*acsch(u(x)+c),(u(x)+c)*asech(u(x)+c)-atan(sqrt(1/(u(x)+c)^2-1)),log(1-(u(x)+c)^2)/2+(u(x)+c)*acoth(u(x)+c)]$
2362 radcan(exponentialize(map(lambda([f], f(x) - diff(integrate(f(x), x), x)), [sin,cos,tan,csc,sec,cot,sinh,cosh,tanh,csch,sech,coth,asin,acos,atan,acsc,asec,acot,asinh,acosh,atanh,acsch,asech,acoth]))), domain:complex;
2363 [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]$
2365 /* Bug #3258: integrate's answer contains diff operator with negative order */
2366 freeof(diff, integrate(%e^(c*x-r^2/(4*x))/x^(3/2),x));
2367 true$
2369 /* Bug #3220: create_list doesn't bind variables properly */
2370 errcatch (create_list (%i, %i, [1, 2, 3]));
2373 create_list (bfloat (1 - 10^-50) - 1, fpprec, [16, 100]);
2374 [0.0b0, -1.0b-50]$
2376 create_list (bfloat (3.14), fpprec, 1, 3);
2377 [3.0b0, 3.1b0, 3.14b0]$
2379 /* Bug #3356: sign(nz * nz) = nz */
2380 assume(apos > 0, aneg < 0, apz >= 0, anz <= 0, notequal(apn, 0), bpos > 0, bneg < 0, bpz >= 0, bnz <= 0, notequal(bpn, 0));
2381 [apos > 0, aneg < 0, apz >= 0, anz <= 0, notequal(apn, 0), bpos > 0, bneg < 0, bpz >= 0, bnz <= 0, notequal(bpn, 0)]$
2382 map('sign,
2383         [
2384                 apos * bpos, apos * bneg, apos * bpz, apos * bnz, apos * bpn, apos * bpnz,
2385                 aneg * bpos, aneg * bneg, aneg * bpz, aneg * bnz, aneg * bpn, aneg * bpnz,
2386                 apz * bpos, apz * bneg, apz * bpz, apz * bnz, apz * bpn, apz * bpnz,
2387                 anz * bpos, anz * bneg, anz * bpz, anz * bnz, anz * bpn, anz * bpnz,
2388                 apn * bpos, apn * bneg, apn * bpz, apn * bnz, apn * bpn, apn * bpnz,
2389                 apnz * bpos, apnz * bneg, apnz * bpz, apnz * bnz, apnz * bpn, apnz * bpnz
2390         ]
2393         'pos, 'neg, 'pz, 'nz, 'pn, 'pnz,
2394         'neg, 'pos, 'nz, 'pz, 'pn, 'pnz,
2395         'pz, 'nz, 'pz, 'nz, 'pnz, 'pnz,
2396         'nz, 'pz, 'nz, 'pz, 'pnz, 'pnz,
2397         'pn, 'pn, 'pnz, 'pnz, 'pn, 'pnz,
2398         'pnz, 'pnz, 'pnz, 'pnz, 'pnz, 'pnz
2400 forget(apos > 0, aneg < 0, apz >= 0, anz <= 0, notequal(apn, 0), bpos > 0, bneg < 0, bpz >= 0, bnz <= 0, notequal(bpn, 0));
2401 [apos > 0, aneg < 0, apz >= 0, anz <= 0, notequal(apn, 0), bpos > 0, bneg < 0, bpz >= 0, bnz <= 0, notequal(bpn, 0)]$
2403 /* test cases mentioned in bug report #3356 */
2405 (kill (a, b), assume(a<=0, b<=0), sign(a*b));
2408 (kill (q, r), assume(q<0, r<0), sign(q*r));
2409 pos;
2411 sign(a*r);
2414 sign(a*-2);
2417 sign(a*q*r);
2420 sign(a*b*q);
2423 (forget (a <= 0, b <= 0, q < 0, r < 0), 0);
2426 /* Bug #3403: Function/lambda parameters declared constant cause error */
2427 declare(c, constant);
2428 done$
2429 emptyp(errcatch(f(c) := c));
2430 false$
2431 emptyp(errcatch(lambda([c], c)));
2432 false$
2433 emptyp(errcatch(f(%e) := %e));
2434 true$
2435 emptyp(errcatch(lambda([%e], %e)));
2436 true$
2437 kill(c, f);
2438 done$
2440 /*        Bug #3009: factoring exponentials causes MAKE-ARRAY error */
2441 /* a.k.a. Bug #3146: max() runs out of memory when comparing exponential functions */
2442 sign(((-1)-%e^-(4.075321792706671E-4*d_1))*%e^-(4.075321792706671E-4*d_2)+%e^-(8.150643585413342E-4*d_2)-%e^-(4.075321792706671E-4*d_1)+%e^-(8.150643585413342E-4*d_1)+1);
2443 pnz$
2444 factor(a^1000000+a+1);
2445 a^1000000+a+1$
2446 max(2.0325-6.9825*%e^-(492380.0*t),2.103-7.053*%e^-(406810.0*t));
2447 'max(2.0325-6.9825*%e^-(492380.0*t),2.103-7.053*%e^-(406810.0*t))$
2448 sign(exp(-492380*t)-2793);
2449 pnz$
2451 /* Bug #3147: sign of expressions with exponents crashes */
2452 assume(t>0);
2453 [t>0]$
2454 sign(2^(500005*t)-2^(500001*t));
2455 pos$
2456 kill(t);
2457 done$
2459 /* Test factor_max_degree. */
2460 block([factor_max_degree : 0], [factor(x^2+2*x+1), factor(x^2-1), factor(x^2+3*x+2), factor(x^2000+x^1999)]);
2461 [(x+1)^2, (x-1)*(x+1), (x+1)*(x+2), x^1999*(x+1)]$
2462 block([factor_max_degree : 1, factor_max_degree_print_warning : false], [factor(x^2+2*x+1), factor(x^2-1), factor(x^2+3*x+2), factor(x^2000+x^1999)]);
2463 [(x+1)^2, x^2-1, x^2+3*x+2, x^1999*(x+1)]$
2464 block([factor_max_degree : 2], [factor(x^2+2*x+1), factor(x^2-1), factor(x^2+3*x+2), factor(x^2000+x^1999)]);
2465 [(x+1)^2, (x-1)*(x+1), (x+1)*(x+2), x^1999*(x+1)]$
2467 /*        Bug #2928: Loops forever: ev(ratexpand(1/(x^(2/3)+1)), algebraic:true); */
2468 /* a.k.a. Bug #2994: Infinite loop: ratsimp(x/(sqrt(x^(2/3)+1)),x),algebraic; */
2469 /* a.k.a. Bug #3419: Endless loop in BPROG (rat(1/(x^(2/3)+1)), algebraic) */
2471 /* examples present in previous version of rtest16 */
2473 string (ev (rat (1/(x^(1/3)+1)), algebraic));
2474 "((x^(1/3))^2-x^(1/3)+1)/(x+1)";
2476 string (ev (rat (1/(x^(2/3)+1)), algebraic));
2477 "(x^(1/3)*x-(x^(1/3))^2+1)/(x^2+1)";
2479 string (ev (rat(1/((x-1)^(2/3)+1)), algebraic));
2480 "((x-1)^(1/3)*x-((x-1)^(1/3))^2-(x-1)^(1/3)+1)/(x^2-2*x+2)";
2482 string (ev (rat(x/(sqrt(x^(2/3)+1))), algebraic));
2483 "(sqrt(x^(2/3)+1)*x^(1/3)*x^2+(-(sqrt(x^(2/3)+1)*(x^(1/3))^2)+sqrt(x^(2/3)+1))*x)/(x^2+1)";
2485 /* additional examples courtesy of M. Talon */
2487 string (ev (rat(1/(x^(2/5)+1)), algebraic));
2488 "(((x^(1/5))^3-x^(1/5))*x+(x^(1/5))^4-(x^(1/5))^2+1)/(x^2+1)";
2490 /* Maxima cannot yet handle this next one (incomplete factorization, gets one factor but misses the other)
2492 rat(1/(x^(3/5)+x^(1/3))),algebraic;
2493  */
2495 string (ev (rat(1/(x^(4/5)+1)), algebraic));
2496 "(x^(1/5)*x^3-(x^(1/5))^2*x^2+(x^(1/5))^3*x-(x^(1/5))^4+1)/(x^4+1)";
2498 /* Bug #3431: error system variable holds unsimplified list, causing errors to be repeated when trying to access it */
2499 errcatch(0^0);
2501 errcatch(print(error), true);
2502 [true]$
2504 /* SF bug #3458: "does not interpret addcol as it should but instead interprets same as addrow" */
2506 addcol (matrix (), [11, 22, 33]);
2507 matrix ([11], [22], [33]);
2509 /* additional tests for addcol and addrow */
2511 addcol (matrix (), [11, 22, 33], [111, 222, 333]);
2512 matrix ([11, 111], [22, 222], [33, 333]);
2514 addcol (matrix (), [11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2515 matrix ([11, 111, 1111], [22, 222, 2222], [33, 333, 3333]);
2517 addcol (matrix ([1, 2, 3]), [11, 22, 33]);
2518 matrix ([1, 2, 3, 11, 22, 33]);
2520 addcol (matrix ([1, 2, 3]), [11, 22, 33], [111, 222, 333]);
2521 matrix ([1, 2, 3, 11, 22, 33, 111, 222, 333]);
2523 addcol (matrix ([1, 2, 3]), [11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2524 matrix ([1, 2, 3, 11, 22, 33, 111, 222, 333, 1111, 2222, 3333]);
2526 addcol (matrix ([1], [2], [3]), [11, 22, 33]);
2527 matrix ([1, 11], [2, 22], [3, 33]);
2529 addcol (matrix ([1], [2], [3]), [11, 22, 33], [111, 222, 333]);
2530 matrix ([1, 11, 111], [2, 22, 222], [3, 33, 333]);
2532 addcol (matrix ([1], [2], [3]), [11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2533 matrix ([1, 11, 111, 1111], [2, 22, 222, 2222], [3, 33, 333, 3333]);
2535 addrow (matrix (), [11, 22, 33]);
2536 matrix ([11, 22, 33]);
2538 addrow (matrix (), [11, 22, 33], [111, 222, 333]);
2539 matrix ([11, 22, 33], [111, 222, 333]);
2541 addrow (matrix (), [11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2542 matrix ([11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2544 addrow (matrix ([1, 2, 3]), [11, 22, 33]);
2545 matrix ([1, 2, 3], [11, 22, 33]);
2547 addrow (matrix ([1, 2, 3]), [11, 22, 33], [111, 222, 333]);
2548 matrix ([1, 2, 3], [11, 22, 33], [111, 222, 333]);
2550 addrow (matrix ([1, 2, 3]), [11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2551 matrix ([1, 2, 3], [11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2554  * Seems like these should succeed, but they fail
2555  * with error about "incompatible structure", although similar calls succeed.
2556  * Just omit these tests unless addrow/addcol are ever revised for greater consistency.
2558 addrow (matrix ([1], [2], [3]), [11, 22, 33]);
2559 matrix ([1], [2], [3], [11], [22], [33]);
2561 addrow (matrix ([1], [2], [3]), [11, 22, 33], [111, 222, 333]);
2562 matrix ([1], [2], [3], [11], [22], [33], [111], [222], [333]);
2564 addrow (matrix ([1], [2], [3]), [11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2565 matrix ([1], [2], [3], [11], [22], [33], [111], [222], [333], [1111], [2222], [3333]);
2567  */
2569 addcol (matrix (), matrix ([11], [22], [33]));
2570 matrix ([11], [22], [33]);
2572 addcol (matrix (), matrix ([11, 22, 33]));
2573 matrix ([11, 22, 33]);
2575 addcol (matrix (), matrix ([11], [22], [33]), matrix ([111], [222], [333]));
2576 matrix ([11, 111], [22, 222], [33, 333]);
2578 addcol (matrix (), matrix ([11, 22, 33]), matrix ([111, 222, 333]));
2579 matrix ([11, 22, 33, 111, 222, 333]);
2581 addcol (matrix (), matrix ([11], [22], [33]), matrix ([111], [222], [333]), matrix ([1111], [2222], [3333]));
2582 matrix ([11, 111, 1111], [22, 222, 2222], [33, 333, 3333]);
2584 addcol (matrix (), matrix ([11, 22, 33]), matrix ([111, 222, 333]), matrix ([1111, 2222, 3333]));
2585 matrix ([11, 22, 33, 111, 222, 333, 1111, 2222, 3333]);
2587 addcol (matrix ([1, 2, 3]), matrix ([11, 22, 33]));
2588 matrix ([1, 2, 3, 11, 22, 33]);
2590 addcol (matrix ([1, 2, 3]), matrix ([11, 22, 33]), matrix ([111, 222, 333]));
2591 matrix ([1, 2, 3, 11, 22, 33, 111, 222, 333]);
2593 addcol (matrix ([1, 2, 3]), matrix ([11, 22, 33]), matrix ([111, 222, 333]), matrix ([1111, 2222, 3333]));
2594 matrix ([1, 2, 3, 11, 22, 33, 111, 222, 333, 1111, 2222, 3333]);
2596 addrow (matrix (), matrix ([11], [22], [33]));
2597 matrix ([11], [22], [33]);
2599 addrow (matrix (), matrix ([11, 22, 33]));
2600 matrix ([11, 22, 33]);
2602 addrow (matrix (), matrix ([11], [22], [33]), matrix ([111], [222], [333]));
2603 matrix ([11], [22], [33], [111], [222], [333]);
2605 addrow (matrix (), matrix ([11, 22, 33]), matrix ([111, 222, 333]));
2606 matrix ([11, 22, 33], [111, 222, 333]);
2608 addrow (matrix (), matrix ([11], [22], [33]), matrix ([111], [222], [333]), matrix ([1111], [2222], [3333]));
2609 matrix ([11], [22], [33], [111], [222], [333], [1111], [2222], [3333]);
2611 addrow (matrix (), matrix ([11, 22, 33]), matrix ([111, 222, 333]), matrix ([1111, 2222, 3333]));
2612 matrix ([11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2614 addrow (matrix ([1, 2, 3]), [11, 22, 33]);
2615 matrix ([1, 2, 3], [11, 22, 33]);
2617 addrow (matrix ([1, 2, 3]), [11, 22, 33], [111, 222, 333]);
2618 matrix ([1, 2, 3], [11, 22, 33], [111, 222, 333]);
2620 addrow (matrix ([1, 2, 3]), [11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2621 matrix ([1, 2, 3], [11, 22, 33], [111, 222, 333], [1111, 2222, 3333]);
2623 /* addrow/addcol that should fail */
2625 errcatch (addrow (matrix ([1], [2], [3]), [11, 22, 33]));
2629  * Seems like these should fail, but they succeed.
2630  * Just omit these tests unless addrow/addcol are ever revised for greater consistency.
2632 errcatch (addrow (matrix ([1], [2], [3]), matrix ([11, 22, 33])));
2635 errcatch (addrow (matrix ([1, 2, 3]), matrix ([1], [2], [3])));
2638 errcatch (addcol (matrix ([1, 2, 3]), matrix ([1], [2], [3])));
2641 errcatch (addcol (matrix ([1], [2], [3]), matrix ([11, 22, 33])));
2644  */
2646 /* Bug #3532: Integrator doesn't use a new variable internally, causing facts on the original variable to be used for the substitution variable */
2647 is(integrate(cos(x) * abs(cos(x)), x, 0, %pi) = %pi / 2);
2648 false$
2650 /* Bug found on mailing list thread "Wrong integral, but antiderivative and limits are correct" */
2651 integrate(x^8/(5+x^9)^2, x, 0, inf);
2652 1/45$
2654 /* Test case for bug mentioned at https://sourceforge.net/p/maxima/mailman/message/8488105/ */
2655 is(abs(integrate(1/(x^4+1), x, 0, 1/2) - 0.49396) < 0.00001);
2656 true$
2658 /* Another bug found on mailing list thread "Wrong integral, but antiderivative and limits are correct" */
2659 freeof(?xz, residue(plog(x)/(x-ratsimp(rectform(%e^(%i*%pi/9)))), x, %e^(%i*%pi/9)));
2660 true$
2662 /* Bug #3562: "integrate(1/(1+tan(x)), x, 0, %pi/2) gives complex result, should be %pi/4" */
2663 integrate(1/(1+tan(x)), x, 0, %pi/2);
2664 %pi/4$
2666 /* Also from bug #3562 */
2667 limit(log((x^2+1)/x^2)-2*log((x-1)/x), x, 0, minus);
2670 /* #3787 fixnump checks in simpexpt */
2671 (declare(n,integer),0);
2674 (-1)^((2^62-1)*n);
2675 (-1)^n$
2677 (-1)^((2^62+1)*n);
2678 (-1)^n$
2680 (-1)^((2^62+2)*n);
2683 (remove(n,integer),0);
2687  * Check sech/csch don't overflow for large numbers.  Use relative
2688  *  error to test for accuracy.
2689  */
2691 closetorel(actual, expected, tol) := block([numer:true, relerr: abs(actual-expected)/expected], if (relerr < tol) then true else relerr);
2692 closetorel(actual, expected, tol) := block([numer:true, relerr: abs(actual-expected)/expected], if (relerr < tol) then true else relerr);
2694 closeto(sech(715.0), 6.033b-311);
2695 true;
2696 closeto(csch(715.0), 6.033b-311);
2697 true;
2701  * Check acsch doesn't overflow for small numbers.  But some lisps
2702  * like clisp don't support denormals, so we have to be careful.
2704  */
2706 if (1e-309 = 0.0) then
2707     closetorel(acsch(1e-300), acsch(bfloat(1e-300)), 1.1375b-16)
2708   else
2709     closetorel(acsch(1e-309), acsch(bfloat(1b-309)), 7.1559b-17);
2710 true;
2712 /* SF bug #3825: "apply('forget, facts()) gives Lisp error" */
2714 (kill (all), declare (F, increasing));
2715 done;
2717 facts ();
2718 [kind (F, increasing)];
2720 featurep (F, increasing);
2721 true;
2723 apply ('forget, facts ());
2724 [kind (F, increasing)];
2726 facts ();
2729 featurep (F, increasing);
2730 false;
2732 declare ([F, G], rational, [H, I], irrational);
2733 done;
2735 facts ();
2736 [kind (F, rational), kind (G, rational), kind (H, irrational), kind (I, irrational)];
2738 [featurep (F, rational),
2739  featurep (G, rational),
2740  featurep (H, irrational),
2741  featurep (I, irrational)];
2742 [true, true, true, true];
2744 forget (kind (G, rational), kind (I, irrational));
2745 [kind (G, rational), kind (I, irrational)];
2747 facts ();
2748 [kind (F, rational), kind (H, irrational)];
2750 forget (kind (F, rational), kind (H, irrational));
2751 [kind (F, rational), kind (H, irrational)];
2753 facts ();
2756 [featurep (F, rational),
2757  featurep (G, rational),
2758  featurep (H, irrational),
2759  featurep (I, irrational)];
2760 [false, false, false, false];
2762 /* Bug #3934: expand(1/(1+%i)^4) => (-4)^(-1) (unsimplified) */
2764 expand((1+%i)^-2);
2765 -%i/2;
2767 expand(1/(1+%i)^4);
2768 -1/4;
2770 expand(4/(%i+1)^4)+1;
2773 expand(1/(1+%i)^8-1/16);
2776 expand((sqrt(-1/2)+sqrt(1/2))^-8);
2779 /* Bug #3956: expand(1/((sqrt(2)-1)*(sqrt(2)+1))) => 1/1 (unsimplified) */
2781 expand(1/((sqrt(2)-1)*(sqrt(2)+1)));
2784 expand(1/(-(sqrt(3)-1)*(sqrt(3)+1)));
2785 -1/2;
2787 expand(1/(-(sqrt(n+1)-1)*(sqrt(n+1)+1)));
2788 -1/n;
2790 /* Bugs #4045 and #4062:
2791    - Different results for integration in Maxima 5.45.1 and 5.46.0
2792    - limit(li[3](x), x, inf) gives li[3](inf) */
2794 integrate((log(x)^2)/(1+x), x, 0, 1);
2795 3*zeta(3)/2;
2797 limit(li[2](x), x, -inf);
2798 minf;
2800 limit(li[2](x), x, inf);
2801 infinity;
2803 limit(li[3](x), x, -inf);
2804 minf;
2806 limit(li[3](x), x, inf);
2807 infinity;
2809 limit(li[n](x), x, -inf);
2810 'limit(li[n](x), x, -inf);
2812 limit(li[2](sin(x)), x, inf);
2813 'limit(li[2](sin(x)), x, inf);
2815 /* Bug #4048: "An incorrect limit" */
2817 makelist(limit(diff(log(tan(%pi/2*tanh(x))), x), x, inf), simpfun, ['identity, lambda([e], trigreduce(trigsimp(trigreduce(e)))), 'exponentialize, 'trigreduce, 'trigrat]);
2818 [2, 2, 2, 2, 2];
2820 /* Bug 4064: "Simple limit triggers Lisp error '1 is not of type LIST'" */
2822 [limit(tan(%pi*tanh(x)), x, inf), limit(tan(%pi*tanh(x)), x, minf)];
2823 [0, 0];
2825 limit(tan(%pi/2*tanh(x)+1), x, inf);
2826 -cot(1);
2828 [limit(cot(%pi*tanh(x)), x, inf), limit(cot(%pi*tanh(x)), x, minf)];
2829 [minf, inf];
2831 [limit(tan(%pi/2*tanh(x)), x, inf), limit(tan(%pi/2*tanh(x)), x, minf)];
2832 [inf, minf];
2834 [limit(cot(%pi/2*tanh(x)), x, inf), limit(cot(%pi/2*tanh(x)), x, minf)];
2835 [0, 0];
2837 closeto(e,tol):=block([numer:true,abse],abse:abs(e),if(abse<tol) then true else abse);
2838 closeto(e,tol):=block([numer:true,abse],abse:abs(e),if(abse<tol) then true else abse);
2840 closeto(float(li[2](%i))-(0.915965594177219*%i-0.20561675835602827), 1e-15);
2841 true;
2843 closeto(bfloat(li[2](%i))-(9.15965594177219015054603b-1*%i - 2.05616758356028304559052b-1), 10*10^(-fpprec));
2844 true;
2846 /* Bug #4317: "Maxima doesn't know that signum(real) can only be -1, 0 or +1" */
2848 is(signum(x)<=1);
2849 true;
2851 is(signum(x)>=-1);
2852 true;
2854 is(signum(x)<2);
2855 true;
2857 is(signum(x)>-2);
2858 true;
2860 is(equal(signum(x),1));
2861 unknown;
2863 is(equal(signum(x),-1));
2864 unknown;
2866 is(equal(signum(x),0));
2867 unknown;
2869 is(equal(signum(x),1/2));
2870 false;
2872 is(equal(signum(x),-1/2));
2873 false;
2875 declare(x,complex);
2876 done;
2878 is(signum(x)<=1);
2879 unknown;
2881 is(signum(x)>=-1);
2882 unknown;
2884 kill(x);
2885 done;