You can subscribe to this list here.
2002 
_{Jan}

_{Feb}

_{Mar}

_{Apr}

_{May}

_{Jun}
(67) 
_{Jul}
(61) 
_{Aug}
(49) 
_{Sep}
(43) 
_{Oct}
(59) 
_{Nov}
(24) 
_{Dec}
(18) 

2003 
_{Jan}
(34) 
_{Feb}
(35) 
_{Mar}
(72) 
_{Apr}
(42) 
_{May}
(46) 
_{Jun}
(15) 
_{Jul}
(64) 
_{Aug}
(62) 
_{Sep}
(22) 
_{Oct}
(41) 
_{Nov}
(57) 
_{Dec}
(56) 
2004 
_{Jan}
(48) 
_{Feb}
(47) 
_{Mar}
(33) 
_{Apr}
(39) 
_{May}
(6) 
_{Jun}
(17) 
_{Jul}
(19) 
_{Aug}
(10) 
_{Sep}
(14) 
_{Oct}
(74) 
_{Nov}
(80) 
_{Dec}
(22) 
2005 
_{Jan}
(43) 
_{Feb}
(33) 
_{Mar}
(52) 
_{Apr}
(74) 
_{May}
(32) 
_{Jun}
(58) 
_{Jul}
(18) 
_{Aug}
(41) 
_{Sep}
(71) 
_{Oct}
(28) 
_{Nov}
(65) 
_{Dec}
(68) 
2006 
_{Jan}
(54) 
_{Feb}
(37) 
_{Mar}
(82) 
_{Apr}
(211) 
_{May}
(69) 
_{Jun}
(75) 
_{Jul}
(279) 
_{Aug}
(139) 
_{Sep}
(135) 
_{Oct}
(58) 
_{Nov}
(81) 
_{Dec}
(78) 
2007 
_{Jan}
(141) 
_{Feb}
(134) 
_{Mar}
(65) 
_{Apr}
(49) 
_{May}
(61) 
_{Jun}
(90) 
_{Jul}
(72) 
_{Aug}
(53) 
_{Sep}
(86) 
_{Oct}
(61) 
_{Nov}
(62) 
_{Dec}
(101) 
2008 
_{Jan}
(100) 
_{Feb}
(66) 
_{Mar}
(76) 
_{Apr}
(95) 
_{May}
(77) 
_{Jun}
(93) 
_{Jul}
(103) 
_{Aug}
(76) 
_{Sep}
(42) 
_{Oct}
(55) 
_{Nov}
(44) 
_{Dec}
(75) 
2009 
_{Jan}
(103) 
_{Feb}
(105) 
_{Mar}
(121) 
_{Apr}
(59) 
_{May}
(103) 
_{Jun}
(82) 
_{Jul}
(67) 
_{Aug}
(76) 
_{Sep}
(85) 
_{Oct}
(75) 
_{Nov}
(181) 
_{Dec}
(133) 
2010 
_{Jan}
(107) 
_{Feb}
(116) 
_{Mar}
(145) 
_{Apr}
(89) 
_{May}
(138) 
_{Jun}
(85) 
_{Jul}
(82) 
_{Aug}
(111) 
_{Sep}
(70) 
_{Oct}
(83) 
_{Nov}
(60) 
_{Dec}
(16) 
2011 
_{Jan}
(61) 
_{Feb}
(16) 
_{Mar}
(52) 
_{Apr}
(41) 
_{May}
(34) 
_{Jun}
(41) 
_{Jul}
(57) 
_{Aug}
(73) 
_{Sep}
(21) 
_{Oct}
(45) 
_{Nov}
(50) 
_{Dec}
(28) 
2012 
_{Jan}
(70) 
_{Feb}
(36) 
_{Mar}
(71) 
_{Apr}
(29) 
_{May}
(48) 
_{Jun}
(61) 
_{Jul}
(44) 
_{Aug}
(54) 
_{Sep}
(20) 
_{Oct}
(28) 
_{Nov}
(41) 
_{Dec}
(137) 
2013 
_{Jan}
(62) 
_{Feb}
(55) 
_{Mar}
(31) 
_{Apr}
(23) 
_{May}
(54) 
_{Jun}
(54) 
_{Jul}
(90) 
_{Aug}
(46) 
_{Sep}
(38) 
_{Oct}
(60) 
_{Nov}
(92) 
_{Dec}
(17) 
2014 
_{Jan}
(62) 
_{Feb}
(35) 
_{Mar}
(72) 
_{Apr}
(30) 
_{May}
(97) 
_{Jun}
(81) 
_{Jul}
(63) 
_{Aug}
(64) 
_{Sep}
(28) 
_{Oct}
(45) 
_{Nov}
(48) 
_{Dec}
(109) 
2015 
_{Jan}
(106) 
_{Feb}
(36) 
_{Mar}
(65) 
_{Apr}
(63) 
_{May}
(95) 
_{Jun}
(56) 
_{Jul}
(48) 
_{Aug}
(55) 
_{Sep}
(100) 
_{Oct}
(57) 
_{Nov}
(33) 
_{Dec}
(46) 
2016 
_{Jan}
(76) 
_{Feb}
(53) 
_{Mar}
(88) 
_{Apr}
(79) 
_{May}
(62) 
_{Jun}
(65) 
_{Jul}
(37) 
_{Aug}
(23) 
_{Sep}
(105) 
_{Oct}

_{Nov}

_{Dec}

S  M  T  W  T  F  S 



1
(2) 
2

3

4
(4) 
5
(9) 
6

7
(1) 
8

9
(4) 
10
(1) 
11
(4) 
12
(3) 
13
(6) 
14
(6) 
15
(3) 
16
(6) 
17
(1) 
18
(5) 
19
(4) 
20

21
(3) 
22
(1) 
23
(2) 
24
(7) 
25
(5) 
26

27
(3) 
28
(7) 
29
(1) 
30
(7) 



From: SourceForge.net <noreply@so...>  20080416 05:02:02

Bugs item #1940411, was opened at 20080411 12:03 Message generated for change (Comment added) made by nobody You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1940411&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: None Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Nobody/Anonymous (nobody) Assigned to: Nobody/Anonymous (nobody) Summary: problem with %pi in sin() Initial Comment: f(x):=sin(4*%pi*x) y:f(x) then, drawing the graph, gets a corrupted sin wave. This problem does not occur if using 2 or 3 intead of 4 in the sin argument, and also it does not occur if typing 3.14 instead of using the %pi constant. Tested using version 0.7.4 of Maxima under Win.  Comment By: Nobody/Anonymous (nobody) Date: 20080415 22:02 Message: Logged In: NO > look for "ticks" and increase it to say 500 or more. Well, I don't think it's correct to say the value of nticks is too small. The problem is that the adaptive plotting algorithm in plot2d is easily fooled into thinking that a function is constant. You can get better results by making nticks a prime number, say 7 or 11 (the default is 10). e.g. plot2d (sin(4*%pi*x), [x, 0, 10]); => spurious flat parts plot2d (sin(4*%pi*x), [x, 0, 10], [nticks, 11]); => looks right To resolve this we probably need to rethink the adaptive plotting stuff. If I'm not mistaken it has other idiosyncrasies. As a stopgap measure maybe we can change the default value of nticks to 11. Robert Dodier (not logged in)  Comment By: Barton Willis (willisbl) Date: 20080411 22:49 Message: Logged In: YES user_id=895922 Originator: NO I'm guessing that you used wxMaxima 0.74 and used the ploting > plot2d menu item. The default for the number of ticks was too small in that version. In the menu box, look for "ticks" and increase it to say 500 or more. I think there is a wxMaxima 0.74a that has a larger value for the default for ticks. Maybe you can look for that. wxMaxima 0.74sa might be part of Maxima 5.14.0 now  I think some Maxima versions 5.14.0 have wxMaxima 0.74 and others 0.74a.  Comment By: Raymond Toy (rtoy) Date: 20080411 12:49 Message: Logged In: YES user_id=28849 Originator: NO Show exactly what you did. In what way is the sine wave corrupted? What does build_info() say? Are you really using Maxima 0.7.4? Is that the wxmaxima version? I don't recall any 0.7.4 version of Maxima.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1940411&group_id=4933 
From: SourceForge.net <noreply@so...>  20080416 04:37:03

Bugs item #534874, was opened at 20020325 13:34 Message generated for change (Comment added) made by nobody You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=534874&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core Group: None Status: Closed Resolution: Fixed Priority: 5 Private: No Submitted By: Jesper Harder (harder) Assigned to: Nobody/Anonymous (nobody) Summary: Fixes for mactex.lisp Initial Comment: Here's a patch that fixes some problems in mactex.lisp: * texmexpt The patch fixes some different problems in texmexpt: tex('diff(f,x)^2); => $$%DERIVATIVE^2\left(\left(f,x,1\right)\right)$$ This is not want we want, so we need to add some more exceptions. tex(f(x)^2); => $$f^2\left(\left(x\right)\right)$$ i.e. there's one set of brackets too much. tex(f(x)^(a+b)); => $$ Error: ((MPLUS) $b $a) is not of type (OR RATIONAL LISP:FLOAT). That happens because we're doing a numeric comparison of the exponent. * texmbox \framebox doesn't work in math mode, use \boxed instead. Also,support labelled boxes. * texchoose tex(binomial(a,b)); => $$\pmatrix{a\\b}$ which doesn't typeset a proper binomial  use the builtin TeX command \choose instead.  Comment By: Nobody/Anonymous (nobody) Date: 20080415 21:37 Message: Logged In: NO Committed r1.63 src/mactex.lisp to output \choose for binomial instead of \pmatrix. Robert Dodier (not logged in)  Comment By: Robert Dodier (robert_dodier) Date: 20060325 17:24 Message: Logged In: YES user_id=501686 Closing this report, as it is already marked "fixed" and, on trying the examples mentioned in the comments, it appears that the patch or patches were applied. (Except tex(binomial(a,b)) still yields pmatrix. Oh well. If that's still a problem then we'll open a new report just for that.)  Comment By: Raymond Toy (rtoy) Date: 20020715 09:21 Message: Logged In: YES user_id=28849 Sorry for messing up your original patch. I've applied your new patch. Let me know if it's working. Also, shouldn't lsum do some instead of just printing lsum? I have no problems if you send patches to the mailing list instead of here. (I'd actually prefer that because I hate this web thing too.)  Comment By: Jesper Harder (harder) Date: 20020703 16:28 Message: Logged In: YES user_id=177224 > You didn't give an example \boxed example, so I don't know > what that does. It just replaces \framebox with \boxed, because \framebox only works in text mode not in math mode. Before the patch: tex(box(a+b)); ==> $$\framebox{b+a}$$ After the patch: tex(box(a+b)); ==> $$\boxed{b+a}$$ The texmlabox stuff adds support for labelled boxed, e.g. tex(box(a+b,foo)); ==> $$\stackrel{foo}{\boxed{b+a}}$$ > Please try this out when you can. It works fine except that you missed a part of my patch: %lsum, %limit are also execeptions, and I think there's a typo  %integral should be %integrate. Here a some examples: tex('limit(f(x),x,0)^2); ==> $$\lim ^2\left(f\left(x\right),x,0\right)$$ tex('lsum(f(x),x,0)^2); ==> $$%LSUM^2\left(f\left(x\right),x,0\right)$$ tex('integrate(f(x),x)^2); ==> $$%INTEGRATE^2\left(f\left(x\right),x\right)$$ I've attached a patch against the newest version. (BTW, I really, really dislike this awkward Sourceforge web interface  isn't there some way to submit patches via email instead?)  Comment By: Raymond Toy (rtoy) Date: 20020626 10:59 Message: Logged In: YES user_id=28849 I'm closing this as fixed. Reopen later if necessary.  Comment By: Raymond Toy (rtoy) Date: 20020626 10:56 Message: Logged In: YES user_id=28849 I'm closing this as fixed. Reopen later if necessary.  Comment By: Raymond Toy (rtoy) Date: 20020626 08:07 Message: Logged In: YES user_id=28849 I'm closing this as fixed. Reopen later if necessary.  Comment By: Raymond Toy (rtoy) Date: 20020622 10:07 Message: Logged In: YES user_id=28849 I applied your patch, with some changes, and all of your examples produce the desired result now. You didn't give an example \boxed example, so I don't know what that does. Please try this out when you can.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=534874&group_id=4933 
From: SourceForge.net <noreply@so...>  20080415 16:43:13

Bugs item #1920177, was opened at 20080319 17:07 Message generated for change (Comment added) made by rtoy You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core Group: Includes proposed fix Status: Open Resolution: None Priority: 5 Private: No Submitted By: Crategus (crategus) Assigned to: Nobody/Anonymous (nobody) Summary: Problems with the bessel functions Initial Comment: There are several problems with the bessel functions. The problems I have found can be divided in the following categories: 1. Inconsistent use of the internal arrays $besselarray, $yarray and $iarray. Example: Restart Maxima and Enter bessel_y(2,1.0). You get an Lisp error. Try first bessel_y(2,1.0) and then repeat bessel_y(2,1.0). Now you get the correct result 1.65...  %i * 3,30... Try bessel_j(2,1.0), you get 0.1149...  %i*2.8142... Next enter bessel_y(2,1.0), you get a Lisp error. The reason for the problems is that the routine bessely uses the global array $YARRAY to store values, but uses also the array $BESSELARRAY to calculate the answer. This is an error. I think the best is to eliminate the use of the global arrays completly. 2. Problematic roundoff errors: Try bessel_j(2,1.0). The result is 0,114903...  %i* 2.8142... e17. The correct result is pure real. The problem is the use of the transformation j[n]=exp(n*%pi*%i)*j[n](x) in the code which produce a small imaginary part. This is a roundoff error for an integer order. In the case of non integer values the imaginary part is correct. So you can not cut off the imaginary part in general. Here is a piece of code which will return an answer which is pure real when the order is an integer (I started to rewrite the code, introduced a function besselj and eliminated the use of the global arrays): (if (integerp order) (if (evenp order) (aref jvals n) ( (aref jvals n)) ) (let ((v (* (cis (* order pi)) (aref jvals n)))) (simplifya `((mplus) ,(realpart v) ((mtimes) $%i ,(imagpart v))) nil ) ) ) 3. Wrong mathematic: Try bessel_j(2.5,1,0). You get 0.04949... Correct is the result 2.8763... * %i Or bessel_j(3,2.0). You get 0,128943... Correct is the result 0.128943... The calculation of the bessel function j[n](x) for real argument x and negativ order n as (realpart (hankel1 order arg)) is not correct. The correct result can be obtained with the formula j[n](x) = 1/2 * (H1[n](x) + H2[n](x)). I tried the following code: (let ((result (* 0.5d0 (+ (hankel1 order arg) (hankel2 order arg))))) (complexify result) ; Problem: you get a small imaginary part ;in the case of a real result (like Problem 2) ; An alternative with a function fpround which round the result ; so that a small imaginary part will vanish ; ; (simplifya ; `((mplus) ; ,(fpround (realpart result) 14) ; ((mtimes) ; ,(fpround (imagpart result) 14) ; $%i ; ) ; ) ; nil ; ) ) This is the definition of the function fpround: (defun fpround (x &optional (digits 1)) (let ((fac (expt 10 digits))) (/ (round x (/ 1 fac)) (float fac)) ) ) I have redesigned a lot of the code for the bessel functions but the work is not finished. Is the work interesting for the project? I use GCL 2.6.6 on Windows XP for programing. Crategus  >Comment By: Raymond Toy (rtoy) Date: 20080415 12:43 Message: Logged In: YES user_id=28849 Originator: NO Thank you for the fixes! I'll fix the comment in German and I'll update the documentation too. I'll add your new tests as well. Once these are done, I'll close this bug report as fixed. If you find more issues, please open a new bug report. Thanks!  Comment By: Raymond Toy (rtoy) Date: 20080415 12:42 Message: Logged In: YES user_id=28849 Originator: NO Thank you for the fixes! I'll fix the comment in German and I'll update the documentation too. I'll add your new tests as well. Once these are done, I'll close this bug report as fixed. If you find more issues, please open a new bug report. Thanks!  Comment By: Raymond Toy (rtoy) Date: 20080415 12:42 Message: Logged In: YES user_id=28849 Originator: NO Thank you for the fixes! I'll fix the comment in German and I'll update the documentation too. I'll add your new tests as well. Once these are done, I'll close this bug report as fixed. If you find more issues, please open a new bug report. Thanks!  Comment By: Crategus (crategus) Date: 20080412 17:15 Message: Logged In: YES user_id=2039760 Originator: YES Thank you very much for your work and your interest in the changes of the code. I have seen that you have deleted the use of the global $besselarray too. I have done further numerical calculations to test the code. I think a good test to show that the different parts of the code work well is to calculate the Wronskians which are defined in A&S. I have attached a file which does the calculations in test mode. (Excuse my programing style, but I have no experience programing in Maxima. Excuse my Englisch too.) As a side effect of the extension of the Bessel functions it is now possible to get nice plots to show the symmetry of the Bessel functions. A hint: You can kill the comment in German I had added at line 465. The question asks for the sign. The sign is correct! A second hint: There a still some parts of the code which don't use the symmetries of the Bessel functions completly. I will work further on this subject. At last: The documentation has to be updated too. I would like to do this work, but my English is not good enough. Crategus File Added: test_bessel_wronskian.mac  Comment By: Raymond Toy (rtoy) Date: 20080411 10:17 Message: Logged In: YES user_id=28849 Originator: NO Changes checked in to bessel.lisp, rev 1.50.  Comment By: Raymond Toy (rtoy) Date: 20080410 12:57 Message: Logged In: YES user_id=28849 Originator: NO I have now integrated your changes into bessel.lisp. The testsuite passes and your bessel tests pass too. I changed test_bessel(a,b,d) to mean abs(ab)<10^(d), which I hope is close enough. (Perhaps the real and imaginary parts should be tested separately?) I'll check in these changes after some code clean up.  Comment By: Crategus (crategus) Date: 20080401 17:02 Message: Logged In: YES user_id=2039760 Originator: YES Again, thank you very much for your interest. I have not changed any of the derivatives. I have only changed the routines which I have listed in this thread. But I have not used the last CVS file. I have used the file which was send with the Maxima 5.14.0 release. There was one change for the derviative for bessel_k since the last release. On Sunday I have installed a tool to download the CVS files. Now I can get the really actual files. The next time I would like to optimize the routines further. Additionally, I am preparing extensive numerical tests for the functions. At the moment I have no idea how to extent the calculation to complex order. In my textbooks and in the book of Abramowitz and Stegun I can't find a relationship which help to calculate the Bessel functions for complex order. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080401 12:41 Message: Logged In: YES user_id=28849 Originator: NO Thank you very much for the updated files. (The diff should have used diff u or diff c, but no matter. The besselchanged.lisp file is perfect.) Given the imminent release of Maxima, I will wait on these changes until after the release. Then I'll try to incorporate your changes. I noticed that you changed the derivatives of some of the Bessel functions. Was that intentional? Were they wrong? (Some of the derivatives need to be in a certain form for some simplifications to be done.)  Comment By: Crategus (crategus) Date: 20080329 17:05 Message: Logged In: YES user_id=2039760 Originator: YES Next I have added a diff between the orginal and the changed file. Question: Is it possible to attach at once more than one file to a message? Crategus File Added: diff.txt  Comment By: Crategus (crategus) Date: 20080329 17:00 Message: Logged In: YES user_id=2039760 Originator: YES I have put all changes of the code for the Bessel functions in the orginal file bessel.lisp which is distributed with Maxima 5.14.0 and attached this file to this message. I have tested the code with the testsuite and the values of the mytest_bessel.mac. Crategus File Added: besselchanged.lisp  Comment By: Crategus (crategus) Date: 20080328 20:35 Message: Logged In: YES user_id=2039760 Originator: YES Thank you very much for your answer. I have redesigned the following routines: simpbesselj (old name besseljsimp) besselj (old name $bessel) simpbessely (old name besselysimp) bessely simpbesseli (old name besselisimp) besseli simpbesselk (old name besselksimp) besselk The new routines are: $hankel_1 simphankel1 $hankel_2 simphankel2 If have further introduced the nouns %hankel_1 und %hankel_2. I would prefer to call the functions $hankel_1 and $hankel_2 and not bessel_hankel. The reason is that these functions are well known as Hankel functions. Because I have collected the routines in a new file it would be a good idea to do the changes in the orginal file bessel.lisp. Than it easier to produce a diff for you. I do this the next time. But first, I would like to have a second look at the global arrays $besselarray, $yarray and $iarray. I think it is not so easy to manage these global arrays in a consistent and predictable way for the user. One reason is, that the code for the numerical calculations use the other functions too (i.e. to calculate bessely we need besselj or for the calculation of besseli we need besselj and besselk). So $besselarray can be destroyed when calculating besseli. It might be possible to prevent this effect by introducing a global flag which signals the internal call of one of the routines. A second point is, that we often dont't use the numerical slatec routines directly. In this case it is necessary to do a lot of extra calculations to obtain the values for the arrays. Perhaps it is useful to decide to fill the global arrays only for the case when we can use the values of the slatec routines directly. This is possible for positive order and argument or positiv order and complex argument of the Bessel functions. Further, an additional flag could be introduced to switch the filling of the arrays on and off. It is also to decide what we do when the user calculate a new value for a Bessel function and the calculation don't fill the global array. In this case, the user may be confused by the fact that the global arrays still hold old values. But for which order and argument? For this case it may be helpful to invalidate the array every time we calculate a new value or to store additionally the order and argument for which the values in the global array are calculated. All these problems with the global arrays arise because of the extension of the routines to negative arguments and orders. In my opinion it would be the best to use the concept of the global arrays no longer. The code becomes much more complicated and the benefit for the user might be small. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080328 12:37 Message: Logged In: YES user_id=28849 Originator: NO First, thank you very much for the fixes to the Bessel routines. I'm sure we all want the routines to be correct! Second, about $besselarray. I see the documentation for that is gone, but I think the intent for besselarray was to let the user have access to all the computed values, since some algorithms end up computing all the intermediate orders to get to the desired order. We need to think if this should go away or not. Third, it would certainly help me a lot if you actually produced a diff between your new code and the existing code. Right now, I have to go look at every single function in different places to figure out what's changed. Trying to minimize the changes would help. I know this is a lot of work for you too. Just identifying which function you actually changed would help a lot too. About the specific changes you mention in the comments. No problem with changing besselxsimp to simpbesselx. (I think most of the code uses simpbesselx, but that's too messy for my tastes.) Extending the range is very, very good. Not sure about the special tests for pure real or imaginary answers. Need to look at that some more. Definitions for Hankel functions are good. May want to name them bessel_hankel_1 to be consistent with other Bessel functions. Certainly want to add some properties like derivatives, but that can wait. All in all, I would very much like to take in your very nice changes.  Comment By: Crategus (crategus) Date: 20080327 12:50 Message: Logged In: YES user_id=2039760 Originator: YES I have tested the Bessel functions with the numerical data attached to this message. For the tests I used a function TEST_BESSEL(value, result, digits) which allow to compare the value with the result within the specified digits. Crategus File Added: mytest_bessel.mac  Comment By: Crategus (crategus) Date: 20080327 12:21 Message: Logged In: YES user_id=2039760 Originator: YES I have finished a first redesign of the code for the numerical calculation of the Bessel functions. I have extended the calculation to or changed parts of the calculation to negative orders and arguments. Perhaps the ideas and the code are interesting for the project. I have added the code to this message. A short description of the changes can be find in the header of the file. Crategus File Added: besselnew.lisp  Comment By: Crategus (crategus) Date: 20080323 13:16 Message: Logged In: YES user_id=2039760 Originator: YES I have added to this posting the code for a routine besselj which handles the above problems. The global array $BESSELARRAY is not used. Futher, I added numerical test values for the the special cases in the code. I used Maxima 5.14.0 and GCL 2.6.6 ANSI to compile the code. The testsuite runs with no unexpected errors found. 22 of the 29 numerical tests I have added give the result within a presision of 16 digits. 7 results have a precision of about 14 to 15 digits. Perhaps, the code is interesting enough for further investigation of the numerical evaluation of the Bessel functions. Crategus File Added: besselj.lisp  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 
From: SourceForge.net <noreply@so...>  20080415 16:42:56

Bugs item #1920177, was opened at 20080319 17:07 Message generated for change (Comment added) made by rtoy You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core Group: Includes proposed fix Status: Open Resolution: None Priority: 5 Private: No Submitted By: Crategus (crategus) Assigned to: Nobody/Anonymous (nobody) Summary: Problems with the bessel functions Initial Comment: There are several problems with the bessel functions. The problems I have found can be divided in the following categories: 1. Inconsistent use of the internal arrays $besselarray, $yarray and $iarray. Example: Restart Maxima and Enter bessel_y(2,1.0). You get an Lisp error. Try first bessel_y(2,1.0) and then repeat bessel_y(2,1.0). Now you get the correct result 1.65...  %i * 3,30... Try bessel_j(2,1.0), you get 0.1149...  %i*2.8142... Next enter bessel_y(2,1.0), you get a Lisp error. The reason for the problems is that the routine bessely uses the global array $YARRAY to store values, but uses also the array $BESSELARRAY to calculate the answer. This is an error. I think the best is to eliminate the use of the global arrays completly. 2. Problematic roundoff errors: Try bessel_j(2,1.0). The result is 0,114903...  %i* 2.8142... e17. The correct result is pure real. The problem is the use of the transformation j[n]=exp(n*%pi*%i)*j[n](x) in the code which produce a small imaginary part. This is a roundoff error for an integer order. In the case of non integer values the imaginary part is correct. So you can not cut off the imaginary part in general. Here is a piece of code which will return an answer which is pure real when the order is an integer (I started to rewrite the code, introduced a function besselj and eliminated the use of the global arrays): (if (integerp order) (if (evenp order) (aref jvals n) ( (aref jvals n)) ) (let ((v (* (cis (* order pi)) (aref jvals n)))) (simplifya `((mplus) ,(realpart v) ((mtimes) $%i ,(imagpart v))) nil ) ) ) 3. Wrong mathematic: Try bessel_j(2.5,1,0). You get 0.04949... Correct is the result 2.8763... * %i Or bessel_j(3,2.0). You get 0,128943... Correct is the result 0.128943... The calculation of the bessel function j[n](x) for real argument x and negativ order n as (realpart (hankel1 order arg)) is not correct. The correct result can be obtained with the formula j[n](x) = 1/2 * (H1[n](x) + H2[n](x)). I tried the following code: (let ((result (* 0.5d0 (+ (hankel1 order arg) (hankel2 order arg))))) (complexify result) ; Problem: you get a small imaginary part ;in the case of a real result (like Problem 2) ; An alternative with a function fpround which round the result ; so that a small imaginary part will vanish ; ; (simplifya ; `((mplus) ; ,(fpround (realpart result) 14) ; ((mtimes) ; ,(fpround (imagpart result) 14) ; $%i ; ) ; ) ; nil ; ) ) This is the definition of the function fpround: (defun fpround (x &optional (digits 1)) (let ((fac (expt 10 digits))) (/ (round x (/ 1 fac)) (float fac)) ) ) I have redesigned a lot of the code for the bessel functions but the work is not finished. Is the work interesting for the project? I use GCL 2.6.6 on Windows XP for programing. Crategus  >Comment By: Raymond Toy (rtoy) Date: 20080415 12:42 Message: Logged In: YES user_id=28849 Originator: NO Thank you for the fixes! I'll fix the comment in German and I'll update the documentation too. I'll add your new tests as well. Once these are done, I'll close this bug report as fixed. If you find more issues, please open a new bug report. Thanks!  Comment By: Raymond Toy (rtoy) Date: 20080415 12:42 Message: Logged In: YES user_id=28849 Originator: NO Thank you for the fixes! I'll fix the comment in German and I'll update the documentation too. I'll add your new tests as well. Once these are done, I'll close this bug report as fixed. If you find more issues, please open a new bug report. Thanks!  Comment By: Crategus (crategus) Date: 20080412 17:15 Message: Logged In: YES user_id=2039760 Originator: YES Thank you very much for your work and your interest in the changes of the code. I have seen that you have deleted the use of the global $besselarray too. I have done further numerical calculations to test the code. I think a good test to show that the different parts of the code work well is to calculate the Wronskians which are defined in A&S. I have attached a file which does the calculations in test mode. (Excuse my programing style, but I have no experience programing in Maxima. Excuse my Englisch too.) As a side effect of the extension of the Bessel functions it is now possible to get nice plots to show the symmetry of the Bessel functions. A hint: You can kill the comment in German I had added at line 465. The question asks for the sign. The sign is correct! A second hint: There a still some parts of the code which don't use the symmetries of the Bessel functions completly. I will work further on this subject. At last: The documentation has to be updated too. I would like to do this work, but my English is not good enough. Crategus File Added: test_bessel_wronskian.mac  Comment By: Raymond Toy (rtoy) Date: 20080411 10:17 Message: Logged In: YES user_id=28849 Originator: NO Changes checked in to bessel.lisp, rev 1.50.  Comment By: Raymond Toy (rtoy) Date: 20080410 12:57 Message: Logged In: YES user_id=28849 Originator: NO I have now integrated your changes into bessel.lisp. The testsuite passes and your bessel tests pass too. I changed test_bessel(a,b,d) to mean abs(ab)<10^(d), which I hope is close enough. (Perhaps the real and imaginary parts should be tested separately?) I'll check in these changes after some code clean up.  Comment By: Crategus (crategus) Date: 20080401 17:02 Message: Logged In: YES user_id=2039760 Originator: YES Again, thank you very much for your interest. I have not changed any of the derivatives. I have only changed the routines which I have listed in this thread. But I have not used the last CVS file. I have used the file which was send with the Maxima 5.14.0 release. There was one change for the derviative for bessel_k since the last release. On Sunday I have installed a tool to download the CVS files. Now I can get the really actual files. The next time I would like to optimize the routines further. Additionally, I am preparing extensive numerical tests for the functions. At the moment I have no idea how to extent the calculation to complex order. In my textbooks and in the book of Abramowitz and Stegun I can't find a relationship which help to calculate the Bessel functions for complex order. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080401 12:41 Message: Logged In: YES user_id=28849 Originator: NO Thank you very much for the updated files. (The diff should have used diff u or diff c, but no matter. The besselchanged.lisp file is perfect.) Given the imminent release of Maxima, I will wait on these changes until after the release. Then I'll try to incorporate your changes. I noticed that you changed the derivatives of some of the Bessel functions. Was that intentional? Were they wrong? (Some of the derivatives need to be in a certain form for some simplifications to be done.)  Comment By: Crategus (crategus) Date: 20080329 17:05 Message: Logged In: YES user_id=2039760 Originator: YES Next I have added a diff between the orginal and the changed file. Question: Is it possible to attach at once more than one file to a message? Crategus File Added: diff.txt  Comment By: Crategus (crategus) Date: 20080329 17:00 Message: Logged In: YES user_id=2039760 Originator: YES I have put all changes of the code for the Bessel functions in the orginal file bessel.lisp which is distributed with Maxima 5.14.0 and attached this file to this message. I have tested the code with the testsuite and the values of the mytest_bessel.mac. Crategus File Added: besselchanged.lisp  Comment By: Crategus (crategus) Date: 20080328 20:35 Message: Logged In: YES user_id=2039760 Originator: YES Thank you very much for your answer. I have redesigned the following routines: simpbesselj (old name besseljsimp) besselj (old name $bessel) simpbessely (old name besselysimp) bessely simpbesseli (old name besselisimp) besseli simpbesselk (old name besselksimp) besselk The new routines are: $hankel_1 simphankel1 $hankel_2 simphankel2 If have further introduced the nouns %hankel_1 und %hankel_2. I would prefer to call the functions $hankel_1 and $hankel_2 and not bessel_hankel. The reason is that these functions are well known as Hankel functions. Because I have collected the routines in a new file it would be a good idea to do the changes in the orginal file bessel.lisp. Than it easier to produce a diff for you. I do this the next time. But first, I would like to have a second look at the global arrays $besselarray, $yarray and $iarray. I think it is not so easy to manage these global arrays in a consistent and predictable way for the user. One reason is, that the code for the numerical calculations use the other functions too (i.e. to calculate bessely we need besselj or for the calculation of besseli we need besselj and besselk). So $besselarray can be destroyed when calculating besseli. It might be possible to prevent this effect by introducing a global flag which signals the internal call of one of the routines. A second point is, that we often dont't use the numerical slatec routines directly. In this case it is necessary to do a lot of extra calculations to obtain the values for the arrays. Perhaps it is useful to decide to fill the global arrays only for the case when we can use the values of the slatec routines directly. This is possible for positive order and argument or positiv order and complex argument of the Bessel functions. Further, an additional flag could be introduced to switch the filling of the arrays on and off. It is also to decide what we do when the user calculate a new value for a Bessel function and the calculation don't fill the global array. In this case, the user may be confused by the fact that the global arrays still hold old values. But for which order and argument? For this case it may be helpful to invalidate the array every time we calculate a new value or to store additionally the order and argument for which the values in the global array are calculated. All these problems with the global arrays arise because of the extension of the routines to negative arguments and orders. In my opinion it would be the best to use the concept of the global arrays no longer. The code becomes much more complicated and the benefit for the user might be small. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080328 12:37 Message: Logged In: YES user_id=28849 Originator: NO First, thank you very much for the fixes to the Bessel routines. I'm sure we all want the routines to be correct! Second, about $besselarray. I see the documentation for that is gone, but I think the intent for besselarray was to let the user have access to all the computed values, since some algorithms end up computing all the intermediate orders to get to the desired order. We need to think if this should go away or not. Third, it would certainly help me a lot if you actually produced a diff between your new code and the existing code. Right now, I have to go look at every single function in different places to figure out what's changed. Trying to minimize the changes would help. I know this is a lot of work for you too. Just identifying which function you actually changed would help a lot too. About the specific changes you mention in the comments. No problem with changing besselxsimp to simpbesselx. (I think most of the code uses simpbesselx, but that's too messy for my tastes.) Extending the range is very, very good. Not sure about the special tests for pure real or imaginary answers. Need to look at that some more. Definitions for Hankel functions are good. May want to name them bessel_hankel_1 to be consistent with other Bessel functions. Certainly want to add some properties like derivatives, but that can wait. All in all, I would very much like to take in your very nice changes.  Comment By: Crategus (crategus) Date: 20080327 12:50 Message: Logged In: YES user_id=2039760 Originator: YES I have tested the Bessel functions with the numerical data attached to this message. For the tests I used a function TEST_BESSEL(value, result, digits) which allow to compare the value with the result within the specified digits. Crategus File Added: mytest_bessel.mac  Comment By: Crategus (crategus) Date: 20080327 12:21 Message: Logged In: YES user_id=2039760 Originator: YES I have finished a first redesign of the code for the numerical calculation of the Bessel functions. I have extended the calculation to or changed parts of the calculation to negative orders and arguments. Perhaps the ideas and the code are interesting for the project. I have added the code to this message. A short description of the changes can be find in the header of the file. Crategus File Added: besselnew.lisp  Comment By: Crategus (crategus) Date: 20080323 13:16 Message: Logged In: YES user_id=2039760 Originator: YES I have added to this posting the code for a routine besselj which handles the above problems. The global array $BESSELARRAY is not used. Futher, I added numerical test values for the the special cases in the code. I used Maxima 5.14.0 and GCL 2.6.6 ANSI to compile the code. The testsuite runs with no unexpected errors found. 22 of the 29 numerical tests I have added give the result within a presision of 16 digits. 7 results have a precision of about 14 to 15 digits. Perhaps, the code is interesting enough for further investigation of the numerical evaluation of the Bessel functions. Crategus File Added: besselj.lisp  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 
From: SourceForge.net <noreply@so...>  20080415 16:42:49

Bugs item #1920177, was opened at 20080319 17:07 Message generated for change (Comment added) made by rtoy You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core Group: Includes proposed fix Status: Open Resolution: None Priority: 5 Private: No Submitted By: Crategus (crategus) Assigned to: Nobody/Anonymous (nobody) Summary: Problems with the bessel functions Initial Comment: There are several problems with the bessel functions. The problems I have found can be divided in the following categories: 1. Inconsistent use of the internal arrays $besselarray, $yarray and $iarray. Example: Restart Maxima and Enter bessel_y(2,1.0). You get an Lisp error. Try first bessel_y(2,1.0) and then repeat bessel_y(2,1.0). Now you get the correct result 1.65...  %i * 3,30... Try bessel_j(2,1.0), you get 0.1149...  %i*2.8142... Next enter bessel_y(2,1.0), you get a Lisp error. The reason for the problems is that the routine bessely uses the global array $YARRAY to store values, but uses also the array $BESSELARRAY to calculate the answer. This is an error. I think the best is to eliminate the use of the global arrays completly. 2. Problematic roundoff errors: Try bessel_j(2,1.0). The result is 0,114903...  %i* 2.8142... e17. The correct result is pure real. The problem is the use of the transformation j[n]=exp(n*%pi*%i)*j[n](x) in the code which produce a small imaginary part. This is a roundoff error for an integer order. In the case of non integer values the imaginary part is correct. So you can not cut off the imaginary part in general. Here is a piece of code which will return an answer which is pure real when the order is an integer (I started to rewrite the code, introduced a function besselj and eliminated the use of the global arrays): (if (integerp order) (if (evenp order) (aref jvals n) ( (aref jvals n)) ) (let ((v (* (cis (* order pi)) (aref jvals n)))) (simplifya `((mplus) ,(realpart v) ((mtimes) $%i ,(imagpart v))) nil ) ) ) 3. Wrong mathematic: Try bessel_j(2.5,1,0). You get 0.04949... Correct is the result 2.8763... * %i Or bessel_j(3,2.0). You get 0,128943... Correct is the result 0.128943... The calculation of the bessel function j[n](x) for real argument x and negativ order n as (realpart (hankel1 order arg)) is not correct. The correct result can be obtained with the formula j[n](x) = 1/2 * (H1[n](x) + H2[n](x)). I tried the following code: (let ((result (* 0.5d0 (+ (hankel1 order arg) (hankel2 order arg))))) (complexify result) ; Problem: you get a small imaginary part ;in the case of a real result (like Problem 2) ; An alternative with a function fpround which round the result ; so that a small imaginary part will vanish ; ; (simplifya ; `((mplus) ; ,(fpround (realpart result) 14) ; ((mtimes) ; ,(fpround (imagpart result) 14) ; $%i ; ) ; ) ; nil ; ) ) This is the definition of the function fpround: (defun fpround (x &optional (digits 1)) (let ((fac (expt 10 digits))) (/ (round x (/ 1 fac)) (float fac)) ) ) I have redesigned a lot of the code for the bessel functions but the work is not finished. Is the work interesting for the project? I use GCL 2.6.6 on Windows XP for programing. Crategus  >Comment By: Raymond Toy (rtoy) Date: 20080415 12:42 Message: Logged In: YES user_id=28849 Originator: NO Thank you for the fixes! I'll fix the comment in German and I'll update the documentation too. I'll add your new tests as well. Once these are done, I'll close this bug report as fixed. If you find more issues, please open a new bug report. Thanks!  Comment By: Crategus (crategus) Date: 20080412 17:15 Message: Logged In: YES user_id=2039760 Originator: YES Thank you very much for your work and your interest in the changes of the code. I have seen that you have deleted the use of the global $besselarray too. I have done further numerical calculations to test the code. I think a good test to show that the different parts of the code work well is to calculate the Wronskians which are defined in A&S. I have attached a file which does the calculations in test mode. (Excuse my programing style, but I have no experience programing in Maxima. Excuse my Englisch too.) As a side effect of the extension of the Bessel functions it is now possible to get nice plots to show the symmetry of the Bessel functions. A hint: You can kill the comment in German I had added at line 465. The question asks for the sign. The sign is correct! A second hint: There a still some parts of the code which don't use the symmetries of the Bessel functions completly. I will work further on this subject. At last: The documentation has to be updated too. I would like to do this work, but my English is not good enough. Crategus File Added: test_bessel_wronskian.mac  Comment By: Raymond Toy (rtoy) Date: 20080411 10:17 Message: Logged In: YES user_id=28849 Originator: NO Changes checked in to bessel.lisp, rev 1.50.  Comment By: Raymond Toy (rtoy) Date: 20080410 12:57 Message: Logged In: YES user_id=28849 Originator: NO I have now integrated your changes into bessel.lisp. The testsuite passes and your bessel tests pass too. I changed test_bessel(a,b,d) to mean abs(ab)<10^(d), which I hope is close enough. (Perhaps the real and imaginary parts should be tested separately?) I'll check in these changes after some code clean up.  Comment By: Crategus (crategus) Date: 20080401 17:02 Message: Logged In: YES user_id=2039760 Originator: YES Again, thank you very much for your interest. I have not changed any of the derivatives. I have only changed the routines which I have listed in this thread. But I have not used the last CVS file. I have used the file which was send with the Maxima 5.14.0 release. There was one change for the derviative for bessel_k since the last release. On Sunday I have installed a tool to download the CVS files. Now I can get the really actual files. The next time I would like to optimize the routines further. Additionally, I am preparing extensive numerical tests for the functions. At the moment I have no idea how to extent the calculation to complex order. In my textbooks and in the book of Abramowitz and Stegun I can't find a relationship which help to calculate the Bessel functions for complex order. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080401 12:41 Message: Logged In: YES user_id=28849 Originator: NO Thank you very much for the updated files. (The diff should have used diff u or diff c, but no matter. The besselchanged.lisp file is perfect.) Given the imminent release of Maxima, I will wait on these changes until after the release. Then I'll try to incorporate your changes. I noticed that you changed the derivatives of some of the Bessel functions. Was that intentional? Were they wrong? (Some of the derivatives need to be in a certain form for some simplifications to be done.)  Comment By: Crategus (crategus) Date: 20080329 17:05 Message: Logged In: YES user_id=2039760 Originator: YES Next I have added a diff between the orginal and the changed file. Question: Is it possible to attach at once more than one file to a message? Crategus File Added: diff.txt  Comment By: Crategus (crategus) Date: 20080329 17:00 Message: Logged In: YES user_id=2039760 Originator: YES I have put all changes of the code for the Bessel functions in the orginal file bessel.lisp which is distributed with Maxima 5.14.0 and attached this file to this message. I have tested the code with the testsuite and the values of the mytest_bessel.mac. Crategus File Added: besselchanged.lisp  Comment By: Crategus (crategus) Date: 20080328 20:35 Message: Logged In: YES user_id=2039760 Originator: YES Thank you very much for your answer. I have redesigned the following routines: simpbesselj (old name besseljsimp) besselj (old name $bessel) simpbessely (old name besselysimp) bessely simpbesseli (old name besselisimp) besseli simpbesselk (old name besselksimp) besselk The new routines are: $hankel_1 simphankel1 $hankel_2 simphankel2 If have further introduced the nouns %hankel_1 und %hankel_2. I would prefer to call the functions $hankel_1 and $hankel_2 and not bessel_hankel. The reason is that these functions are well known as Hankel functions. Because I have collected the routines in a new file it would be a good idea to do the changes in the orginal file bessel.lisp. Than it easier to produce a diff for you. I do this the next time. But first, I would like to have a second look at the global arrays $besselarray, $yarray and $iarray. I think it is not so easy to manage these global arrays in a consistent and predictable way for the user. One reason is, that the code for the numerical calculations use the other functions too (i.e. to calculate bessely we need besselj or for the calculation of besseli we need besselj and besselk). So $besselarray can be destroyed when calculating besseli. It might be possible to prevent this effect by introducing a global flag which signals the internal call of one of the routines. A second point is, that we often dont't use the numerical slatec routines directly. In this case it is necessary to do a lot of extra calculations to obtain the values for the arrays. Perhaps it is useful to decide to fill the global arrays only for the case when we can use the values of the slatec routines directly. This is possible for positive order and argument or positiv order and complex argument of the Bessel functions. Further, an additional flag could be introduced to switch the filling of the arrays on and off. It is also to decide what we do when the user calculate a new value for a Bessel function and the calculation don't fill the global array. In this case, the user may be confused by the fact that the global arrays still hold old values. But for which order and argument? For this case it may be helpful to invalidate the array every time we calculate a new value or to store additionally the order and argument for which the values in the global array are calculated. All these problems with the global arrays arise because of the extension of the routines to negative arguments and orders. In my opinion it would be the best to use the concept of the global arrays no longer. The code becomes much more complicated and the benefit for the user might be small. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080328 12:37 Message: Logged In: YES user_id=28849 Originator: NO First, thank you very much for the fixes to the Bessel routines. I'm sure we all want the routines to be correct! Second, about $besselarray. I see the documentation for that is gone, but I think the intent for besselarray was to let the user have access to all the computed values, since some algorithms end up computing all the intermediate orders to get to the desired order. We need to think if this should go away or not. Third, it would certainly help me a lot if you actually produced a diff between your new code and the existing code. Right now, I have to go look at every single function in different places to figure out what's changed. Trying to minimize the changes would help. I know this is a lot of work for you too. Just identifying which function you actually changed would help a lot too. About the specific changes you mention in the comments. No problem with changing besselxsimp to simpbesselx. (I think most of the code uses simpbesselx, but that's too messy for my tastes.) Extending the range is very, very good. Not sure about the special tests for pure real or imaginary answers. Need to look at that some more. Definitions for Hankel functions are good. May want to name them bessel_hankel_1 to be consistent with other Bessel functions. Certainly want to add some properties like derivatives, but that can wait. All in all, I would very much like to take in your very nice changes.  Comment By: Crategus (crategus) Date: 20080327 12:50 Message: Logged In: YES user_id=2039760 Originator: YES I have tested the Bessel functions with the numerical data attached to this message. For the tests I used a function TEST_BESSEL(value, result, digits) which allow to compare the value with the result within the specified digits. Crategus File Added: mytest_bessel.mac  Comment By: Crategus (crategus) Date: 20080327 12:21 Message: Logged In: YES user_id=2039760 Originator: YES I have finished a first redesign of the code for the numerical calculation of the Bessel functions. I have extended the calculation to or changed parts of the calculation to negative orders and arguments. Perhaps the ideas and the code are interesting for the project. I have added the code to this message. A short description of the changes can be find in the header of the file. Crategus File Added: besselnew.lisp  Comment By: Crategus (crategus) Date: 20080323 13:16 Message: Logged In: YES user_id=2039760 Originator: YES I have added to this posting the code for a routine besselj which handles the above problems. The global array $BESSELARRAY is not used. Futher, I added numerical test values for the the special cases in the code. I used Maxima 5.14.0 and GCL 2.6.6 ANSI to compile the code. The testsuite runs with no unexpected errors found. 22 of the 29 numerical tests I have added give the result within a presision of 16 digits. 7 results have a precision of about 14 to 15 digits. Perhaps, the code is interesting enough for further investigation of the numerical evaluation of the Bessel functions. Crategus File Added: besselj.lisp  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 
From: SourceForge.net <noreply@so...>  20080414 03:08:53

Bugs item #1898852, was opened at 20080221 10:46 Message generated for change (Comment added) made by macrakis You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1898852&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Simplification Group: To be reviewed Status: Open Resolution: None Priority: 5 Private: No Submitted By: David Ronis (ronis) Assigned to: Nobody/Anonymous (nobody) Summary: ev doesn't work as expected with subscripted variables Initial Comment: f(l,m):=A[l,m]+B[l,m]+A[l+1,m]; ev(f(l,m), A[l,m]=0); A[l,m] still appears in the result. However, f(l,m); ev(%, A[l,m]=0 ); works as expected.  >Comment By: Stavros Macrakis (macrakis) Date: 20080413 23:08 Message: Logged In: YES user_id=588346 Originator: NO Thanks for the bug report. EV is a very funny function, and often produces peculiar surprises like this. I hope we'll be able to fix this at some point, but in the meantime, I'd suggest you use subst to evaluate things like this: subst([a[l,m]=0],f(l,m)) instead of ev.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1898852&group_id=4933 
From: SourceForge.net <noreply@so...>  20080414 03:03:13

Bugs item #1941133, was opened at 20080412 19:31 Message generated for change (Comment added) made by macrakis You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941133&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Solving equations Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Nobody/Anonymous (nobody) Assigned to: Nobody/Anonymous (nobody) Summary: Logarithms Initial Comment: Does not always simplify equations with logarithms to the single variable.  >Comment By: Stavros Macrakis (macrakis) Date: 20080413 23:03 Message: Logged In: YES user_id=588346 Originator: NO Could you please give a more complete bug report, including an example of this issue? This report doesn't give us enough information to know what you're talking about. In the meantime, you might look at the documentation for radcan and for logcontract.  Comment By: Stavros Macrakis (macrakis) Date: 20080413 23:02 Message: Logged In: YES user_id=588346 Originator: NO Could you please give a more complete bug report, including an example of this issue? This report doesn't give us enough information to know what you're talking about. In the meantime, you might look at the documentation for radcan and for logcontract.  Comment By: Stavros Macrakis (macrakis) Date: 20080413 22:53 Message: Logged In: YES user_id=588346 Originator: NO Could you please give a more complete bug report, including an example of this issue? This report doesn't give us enough information to know what you're talking about. In the meantime, you might look at the documentation for radcan and for logcontract.  Comment By: Stavros Macrakis (macrakis) Date: 20080413 22:52 Message: Logged In: YES user_id=588346 Originator: NO Could you please give a more complete bug report, including an example of this issue? This report doesn't give us enough information to know what you're talking about. In the meantime, you might look at the documentation for radcan and for logcontract.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941133&group_id=4933 
From: SourceForge.net <noreply@so...>  20080414 03:02:05

Bugs item #1941133, was opened at 20080412 19:31 Message generated for change (Comment added) made by macrakis You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941133&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Solving equations Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Nobody/Anonymous (nobody) Assigned to: Nobody/Anonymous (nobody) Summary: Logarithms Initial Comment: Does not always simplify equations with logarithms to the single variable.  >Comment By: Stavros Macrakis (macrakis) Date: 20080413 23:02 Message: Logged In: YES user_id=588346 Originator: NO Could you please give a more complete bug report, including an example of this issue? This report doesn't give us enough information to know what you're talking about. In the meantime, you might look at the documentation for radcan and for logcontract.  Comment By: Stavros Macrakis (macrakis) Date: 20080413 22:53 Message: Logged In: YES user_id=588346 Originator: NO Could you please give a more complete bug report, including an example of this issue? This report doesn't give us enough information to know what you're talking about. In the meantime, you might look at the documentation for radcan and for logcontract.  Comment By: Stavros Macrakis (macrakis) Date: 20080413 22:52 Message: Logged In: YES user_id=588346 Originator: NO Could you please give a more complete bug report, including an example of this issue? This report doesn't give us enough information to know what you're talking about. In the meantime, you might look at the documentation for radcan and for logcontract.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941133&group_id=4933 
From: SourceForge.net <noreply@so...>  20080414 02:54:01

Bugs item #1941133, was opened at 20080412 19:31 Message generated for change (Comment added) made by macrakis You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941133&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Solving equations Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Nobody/Anonymous (nobody) Assigned to: Nobody/Anonymous (nobody) Summary: Logarithms Initial Comment: Does not always simplify equations with logarithms to the single variable.  >Comment By: Stavros Macrakis (macrakis) Date: 20080413 22:53 Message: Logged In: YES user_id=588346 Originator: NO Could you please give a more complete bug report, including an example of this issue? This report doesn't give us enough information to know what you're talking about. In the meantime, you might look at the documentation for radcan and for logcontract.  Comment By: Stavros Macrakis (macrakis) Date: 20080413 22:52 Message: Logged In: YES user_id=588346 Originator: NO Could you please give a more complete bug report, including an example of this issue? This report doesn't give us enough information to know what you're talking about. In the meantime, you might look at the documentation for radcan and for logcontract.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941133&group_id=4933 
From: SourceForge.net <noreply@so...>  20080414 02:52:55

Bugs item #1941133, was opened at 20080412 19:31 Message generated for change (Comment added) made by macrakis You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941133&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Solving equations Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Nobody/Anonymous (nobody) Assigned to: Nobody/Anonymous (nobody) Summary: Logarithms Initial Comment: Does not always simplify equations with logarithms to the single variable.  >Comment By: Stavros Macrakis (macrakis) Date: 20080413 22:52 Message: Logged In: YES user_id=588346 Originator: NO Could you please give a more complete bug report, including an example of this issue? This report doesn't give us enough information to know what you're talking about. In the meantime, you might look at the documentation for radcan and for logcontract.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941133&group_id=4933 
From: SourceForge.net <noreply@so...>  20080414 02:17:43

Bugs item #1941682, was opened at 20080413 19:17 Message generated for change (Tracker Item Submitted) made by Item Submitter You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941682&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Xmaxima or other UI Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Nobody/Anonymous (nobody) Assigned to: Nobody/Anonymous (nobody) Summary: After a while i get the message "Starting maxima timed out" Initial Comment: After a while I get the message "Starting maxima timed out. Wait longer?" and xmaxima doesn't connect to maxima. Maxima works fine, but it doesn't connect to xmaxima neither wxmaxima (this last gets the message "") I use Ubuntu 7.10, and maxima 5.10.0 and wxmaxima 0.7.1  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941682&group_id=4933 
From: SourceForge.net <noreply@so...>  20080413 18:30:03

Bugs item #1755392, was opened at 20070717 06:31 Message generated for change (Comment added) made by dgildea You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1755392&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Trigonometry Group: None >Status: Closed >Resolution: Works For Me Priority: 4 Private: No Submitted By: Barton Willis (willisbl) Assigned to: Nobody/Anonymous (nobody) Summary: trig with xxx[%pi] arguments Initial Comment: The trig functions have trouble with some arguments that involve xxx[%pi] or xxx[%pi + 1], or ... ; for example (%i1) cos(2*%pi + a[%pi]); (%o1) 1 (%i2) sin(2*%pi + a[%pi]); (%o2) 0 (%i3) tan(2*%pi + a[%pi]); (%o3) 0 (%i15) cos(2*%pi + a[%pi+1]); (%o15) 1  >Comment By: Dan Gildea (dgildea) Date: 20080413 14:30 Message: Logged In: YES user_id=1797506 Originator: NO Seems OK in current cvs. (%i42) cos(2*%pi + a[%pi]); (%o42) cos(a[%pi]) (%i43) sin(2*%pi + a[%pi]); (%o43) sin(a[%pi]) (%i44) tan(2*%pi + a[%pi]); (%o44) tan(a[%pi]) (%i45) cos(2*%pi + a[%pi+1]); (%o45) cos(a[%pi+1])  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1755392&group_id=4933 
From: SourceForge.net <noreply@so...>  20080413 18:26:54

Bugs item #1755550, was opened at 20070717 12:03 Message generated for change (Comment added) made by dgildea You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1755550&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Trigonometry Group: None >Status: Closed >Resolution: Fixed Priority: 3 Private: No Submitted By: Barton Willis (willisbl) Assigned to: Nobody/Anonymous (nobody) Summary: trig inconsistent about periodicity rules Initial Comment: OK  use perodicity to simplify to zero: (%i5) cos(2*%pi + %e^2)  cos(%e^2); (%o5) 0 Not so OK: (%i6) cos(2*%pi + %pi^2)  cos(%pi^2); (%o6) cos(%pi^2+2*%pi)cos(%pi^2) Sure, we can apply trigexpand: (%i7) trigexpand(%); (%o7) 0 This isn't a bug, I suppose. But we could do better.  >Comment By: Dan Gildea (dgildea) Date: 20080413 14:26 Message: Logged In: YES user_id=1797506 Originator: NO Fixed in trigi.lisp rev 1.30. (%i40) cos(2*%pi + %pi^2)  cos(%pi^2); (%o40) 0  Comment By: Harald Geyer (hgeyer) Date: 20070717 17:57 Message: Logged In: YES user_id=929336 Originator: NO I think the problem is, that %piargs only touches expressions that are linear in %pi. Compare this with Bug #1553866.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1755550&group_id=4933 
From: SourceForge.net <noreply@so...>  20080413 18:26:00

Bugs item #1553866, was opened at 20060907 01:13 Message generated for change (Settings changed) made by dgildea You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1553866&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Trigonometry Group: None >Status: Closed >Resolution: Fixed Priority: 3 Private: No Submitted By: Barton Willis (willisbl) Assigned to: Nobody/Anonymous (nobody) Summary: %piargs inconsistent behavior Initial Comment: (1) I can't find documenation for %piargs. (2) The way %piargs works is inconsistent and silly: (%i34) cos(a*x + %pi); (%o34) cos(a*x) (%i35) cos(%pi*x + %pi); (%o35) cos(%pi*x+%pi) < why not cos(%pi*x) ? (3) When %piargs is true, simplifying nfold compositions of trig functions takes O(3^n) time. This is due to the function linearp that does expand and forces a new simplification. Try cos(cos(....(x)))) with %piargs true then with %piargs false. Use about 14 compositions. So making %piargs less silly will also make things like cos(cos(...(x)))) much faster. Barton  >Comment By: Dan Gildea (dgildea) Date: 20080413 14:25 Message: Logged In: YES user_id=1797506 Originator: NO (1) Documentation has been added previously. (2) Fixed in trigi.lisp rev 1.30 (%i39) cos(%pi*x + %pi); (%o39) cos(%pi*x) (3) Faster now, too.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1553866&group_id=4933 
From: SourceForge.net <noreply@so...>  20080413 18:24:24

Bugs item #903190, was opened at 20040223 22:15 Message generated for change (Comment added) made by dgildea You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=903190&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Trigonometry Group: None >Status: Closed >Resolution: Fixed Priority: 5 Private: No Submitted By: Stavros Macrakis (macrakis) Assigned to: Nobody/Anonymous (nobody) Summary: sin(%pi+%pi*x) doesn't simplify Initial Comment: sin(%pi+%pi*x) doesn't simplify (cf 580721) ...similarly for other trig functions. The problem is that the %piargs code only looks for xxx+k*%pi/2 where xxx is %pifree. Even sin(f(%pi)+%pi) doesn't simplify! There are two main ways to try to fix this: 1) *Add* a purely syntactic check before or after the existing semantic check. This is appealing because a) it will clearly work and b) it can't screw up the form of the expression. 2) Replace the current semantic check with one which looks for a k*%pi/2 term among other terms. For example, sin((%pi+1)^3) has a 12*%pi term, so could simplify to sin(%pi^3+6*%pi^2+8). But is that really useful? I don't think so, and if the user wanted it with the syntactic approach, s/he could just use expand. Even less useful would be sin((%pi+1)^312*%pi). So I think solution (1) is better, even if it seems inelegant to do two independent checks.  >Comment By: Dan Gildea (dgildea) Date: 20080413 14:24 Message: Logged In: YES user_id=1797506 Originator: NO Fixed in trigi.lisp rev 1.30 (using solution 2 above, however). (%i32) sin(%pi+%pi*x); (%o32) sin(%pi*x) (%i33) sin((%pi+1)^3); (%o33) sin(%pi^3+3*%pi^2+1)  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=903190&group_id=4933 
From: SourceForge.net <noreply@so...>  20080413 18:18:08

Bugs item #580721, was opened at 20020712 14:38 Message generated for change (Comment added) made by dgildea You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=580721&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Trigonometry Group: None >Status: Closed >Resolution: Fixed Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: trigexpand bug Initial Comment: Consider this: (C1) display2d:false; (D1) FALSE (C2) tan(%pi/2+x); (D2) COT(x) (C3) tan(%pi/2+%pi*x); (D3) TAN(%PI*x+%PI/2) (C4) %,trigexpand; Division by 0  an error. Quitting. To debug this try DEBUGMODE(TRUE);) d3 should have been expanded into cot(%pi*x) instead of getting a division by zero error.  >Comment By: Dan Gildea (dgildea) Date: 20080413 14:18 Message: Logged In: YES user_id=1797506 Originator: NO Fixed in trigi.lisp rev 1.30 (%i29) tan(%pi/2+%pi*x); (%o29) cot(%pi*x) (%i30) declare(n, integer); (%o30) done (%i31) tan(n*2*%pi); (%o31) 0  Comment By: Dan Gildea (dgildea) Date: 20080325 21:29 Message: Logged In: YES user_id=1797506 Originator: NO The proposed patch disables the following simplification: (%i1) declare(n, integer); (%o1) done (%i2) tan(n*2*%pi); (%o2) 0 However, it almost fixes bugs 903190, 1553866 and 1755550.  Comment By: Rupert Swarbrick (rswarbrick) Date: 20071228 12:15 Message: Logged In: YES user_id=1673565 Originator: NO I think I've posted a patch onto the mailing list that fixes this. I can't work out how to attach files here, so the link is: http://thread.gmane.org/gmane.comp.mathematics.maxima.general/19076  Comment By: Rupert Swarbrick (rswarbrick) Date: 20071226 21:38 Message: Logged In: YES user_id=1673565 Originator: NO So the first thing that's wrong is that %piargstan/cot shouldn't return nonnil for stuff like tan(pi/2) as they aren't defined, let alone simplifiable. Here's an amended and reformatted version with variables renamed and a possible bug due to reuse of variable names eliminated too (I think) (defun %piargstan/cot (x) (displa x) (let ((coeff (linearize (coefficient x '$%pi 1))) (zlrem (coefficient x '$%pi 0)) (sinofcoeffpi) (cosofcoeffpi)) (cond ((and (zerop1 zlrem) (setq sinofcoeffpi (%piargs coeff nil)) (not (zerop1 (setq cosofcoeffpi (%piargs (cons (car coeff) (rplus 1//2 (cdr coeff))) nil))))) ;; sinofcoeffpi and cosofcoeffpi are only nonnil if they ;; are constants that %piargsoffset could compute, and we just ;; checked that cosofcoeffpi was nonzero. Thus we can just ;; return their quotient. (div sinofcoeffpi cosofcoeffpi)) ((not (mevenp (car coeff))) nil) ((integerp (setq x (mmod (cdr coeff) 2))) (consexp '%tan zlrem)) ((or (alike1 1//2 x) (alike1 '((rat) 3 2) x)) (neg (consexp '%cot zlrem))))))  Comment By: Rupert Swarbrick (rswarbrick) Date: 20071226 21:29 Message: Logged In: YES user_id=1673565 Originator: NO Ah. Of course tan(pi/2) = infinity ...  Comment By: Rupert Swarbrick (rswarbrick) Date: 20071225 20:10 Message: Logged In: YES user_id=1673565 Originator: NO Hi, a bit more information: the error's caused by a (div 1 0) which gets returned in %piargstan/cot, called by simp%tan. You can reproduce more simply by just calling tan(%pi/2), trigexpand; I'm trying to work out what the functions are _supposed_ to do. Maybe then I'll be able to fix it!  Comment By: Robert Dodier (robert_dodier) Date: 20060326 18:40 Message: Logged In: YES user_id=501686 For the record, same problem observed in maxima 5.9.3.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=580721&group_id=4933 
From: SourceForge.net <noreply@so...>  20080413 17:47:22

Bugs item #1941444, was opened at 20080413 19:47 Message generated for change (Tracker Item Submitted) made by Item Submitter You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941444&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core Group: Includes proposed fix Status: Open Resolution: None Priority: 5 Private: No Submitted By: Crategus (crategus) Assigned to: Nobody/Anonymous (nobody) Summary: Inconsistent use of global $DEBUGMODE Initial Comment: I think there is a small inconsistency with the global variable $DEBUGMODE. ********** The documentation says: Option variable: debugmode Default value: false When a Maxima error occurs, Maxima will start the debugger if debugmode is true. ... ********** Intern in the code Maxima only uses the global *MDEBUG*. When you assign the value TRUE or FALSE to $DEBUGMODE the global *MDEBUG* will be set with the assignfunction DEBUGMODE1() and all is correct. But, there is also the undocumented function $DEBUGMODE(). In the case of an error the user gets the information " an error. To debug this try debugmode(true)". So the user knows about the function, too. If you call this function the global *MDEBUG* will be set correctly as in the case above. But now the global $DEBUGMODE is not set to the new value. A program which tests the global $DEBUGMODE might not be sure if Maxima is in debugmode or not. The value of $DEBUGMODE depends on the way we have entered the debugmode. A correction in suprv1.lisp might be: (defmfun $debugmode (x) (setq $debugmode x) ; New line: Set the global $debugmode, too. (debugmode1 nil x)) Here is the assignfunction from suprv1.lisp which sets the global *MDEBUG*: (defun debugmode1 (assignvar y) (declare (ignore assignvar)) (setq *mdebug* (setq *rset y))) Crategus  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941444&group_id=4933 
From: SourceForge.net <noreply@so...>  20080412 23:31:46

Bugs item #1941133, was opened at 20080412 16:31 Message generated for change (Tracker Item Submitted) made by Item Submitter You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941133&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core  Solving equations Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Nobody/Anonymous (nobody) Assigned to: Nobody/Anonymous (nobody) Summary: Logarithms Initial Comment: Does not always simplify equations with logarithms to the single variable.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1941133&group_id=4933 
From: SourceForge.net <noreply@so...>  20080412 21:15:44

Bugs item #1920177, was opened at 20080319 22:07 Message generated for change (Comment added) made by crategus You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core Group: Includes proposed fix Status: Open Resolution: None Priority: 5 Private: No Submitted By: Crategus (crategus) Assigned to: Nobody/Anonymous (nobody) Summary: Problems with the bessel functions Initial Comment: There are several problems with the bessel functions. The problems I have found can be divided in the following categories: 1. Inconsistent use of the internal arrays $besselarray, $yarray and $iarray. Example: Restart Maxima and Enter bessel_y(2,1.0). You get an Lisp error. Try first bessel_y(2,1.0) and then repeat bessel_y(2,1.0). Now you get the correct result 1.65...  %i * 3,30... Try bessel_j(2,1.0), you get 0.1149...  %i*2.8142... Next enter bessel_y(2,1.0), you get a Lisp error. The reason for the problems is that the routine bessely uses the global array $YARRAY to store values, but uses also the array $BESSELARRAY to calculate the answer. This is an error. I think the best is to eliminate the use of the global arrays completly. 2. Problematic roundoff errors: Try bessel_j(2,1.0). The result is 0,114903...  %i* 2.8142... e17. The correct result is pure real. The problem is the use of the transformation j[n]=exp(n*%pi*%i)*j[n](x) in the code which produce a small imaginary part. This is a roundoff error for an integer order. In the case of non integer values the imaginary part is correct. So you can not cut off the imaginary part in general. Here is a piece of code which will return an answer which is pure real when the order is an integer (I started to rewrite the code, introduced a function besselj and eliminated the use of the global arrays): (if (integerp order) (if (evenp order) (aref jvals n) ( (aref jvals n)) ) (let ((v (* (cis (* order pi)) (aref jvals n)))) (simplifya `((mplus) ,(realpart v) ((mtimes) $%i ,(imagpart v))) nil ) ) ) 3. Wrong mathematic: Try bessel_j(2.5,1,0). You get 0.04949... Correct is the result 2.8763... * %i Or bessel_j(3,2.0). You get 0,128943... Correct is the result 0.128943... The calculation of the bessel function j[n](x) for real argument x and negativ order n as (realpart (hankel1 order arg)) is not correct. The correct result can be obtained with the formula j[n](x) = 1/2 * (H1[n](x) + H2[n](x)). I tried the following code: (let ((result (* 0.5d0 (+ (hankel1 order arg) (hankel2 order arg))))) (complexify result) ; Problem: you get a small imaginary part ;in the case of a real result (like Problem 2) ; An alternative with a function fpround which round the result ; so that a small imaginary part will vanish ; ; (simplifya ; `((mplus) ; ,(fpround (realpart result) 14) ; ((mtimes) ; ,(fpround (imagpart result) 14) ; $%i ; ) ; ) ; nil ; ) ) This is the definition of the function fpround: (defun fpround (x &optional (digits 1)) (let ((fac (expt 10 digits))) (/ (round x (/ 1 fac)) (float fac)) ) ) I have redesigned a lot of the code for the bessel functions but the work is not finished. Is the work interesting for the project? I use GCL 2.6.6 on Windows XP for programing. Crategus  >Comment By: Crategus (crategus) Date: 20080412 23:15 Message: Logged In: YES user_id=2039760 Originator: YES Thank you very much for your work and your interest in the changes of the code. I have seen that you have deleted the use of the global $besselarray too. I have done further numerical calculations to test the code. I think a good test to show that the different parts of the code work well is to calculate the Wronskians which are defined in A&S. I have attached a file which does the calculations in test mode. (Excuse my programing style, but I have no experience programing in Maxima. Excuse my Englisch too.) As a side effect of the extension of the Bessel functions it is now possible to get nice plots to show the symmetry of the Bessel functions. A hint: You can kill the comment in German I had added at line 465. The question asks for the sign. The sign is correct! A second hint: There a still some parts of the code which don't use the symmetries of the Bessel functions completly. I will work further on this subject. At last: The documentation has to be updated too. I would like to do this work, but my English is not good enough. Crategus File Added: test_bessel_wronskian.mac  Comment By: Raymond Toy (rtoy) Date: 20080411 16:17 Message: Logged In: YES user_id=28849 Originator: NO Changes checked in to bessel.lisp, rev 1.50.  Comment By: Raymond Toy (rtoy) Date: 20080410 18:57 Message: Logged In: YES user_id=28849 Originator: NO I have now integrated your changes into bessel.lisp. The testsuite passes and your bessel tests pass too. I changed test_bessel(a,b,d) to mean abs(ab)<10^(d), which I hope is close enough. (Perhaps the real and imaginary parts should be tested separately?) I'll check in these changes after some code clean up.  Comment By: Crategus (crategus) Date: 20080401 23:02 Message: Logged In: YES user_id=2039760 Originator: YES Again, thank you very much for your interest. I have not changed any of the derivatives. I have only changed the routines which I have listed in this thread. But I have not used the last CVS file. I have used the file which was send with the Maxima 5.14.0 release. There was one change for the derviative for bessel_k since the last release. On Sunday I have installed a tool to download the CVS files. Now I can get the really actual files. The next time I would like to optimize the routines further. Additionally, I am preparing extensive numerical tests for the functions. At the moment I have no idea how to extent the calculation to complex order. In my textbooks and in the book of Abramowitz and Stegun I can't find a relationship which help to calculate the Bessel functions for complex order. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080401 18:41 Message: Logged In: YES user_id=28849 Originator: NO Thank you very much for the updated files. (The diff should have used diff u or diff c, but no matter. The besselchanged.lisp file is perfect.) Given the imminent release of Maxima, I will wait on these changes until after the release. Then I'll try to incorporate your changes. I noticed that you changed the derivatives of some of the Bessel functions. Was that intentional? Were they wrong? (Some of the derivatives need to be in a certain form for some simplifications to be done.)  Comment By: Crategus (crategus) Date: 20080329 22:05 Message: Logged In: YES user_id=2039760 Originator: YES Next I have added a diff between the orginal and the changed file. Question: Is it possible to attach at once more than one file to a message? Crategus File Added: diff.txt  Comment By: Crategus (crategus) Date: 20080329 22:00 Message: Logged In: YES user_id=2039760 Originator: YES I have put all changes of the code for the Bessel functions in the orginal file bessel.lisp which is distributed with Maxima 5.14.0 and attached this file to this message. I have tested the code with the testsuite and the values of the mytest_bessel.mac. Crategus File Added: besselchanged.lisp  Comment By: Crategus (crategus) Date: 20080329 01:35 Message: Logged In: YES user_id=2039760 Originator: YES Thank you very much for your answer. I have redesigned the following routines: simpbesselj (old name besseljsimp) besselj (old name $bessel) simpbessely (old name besselysimp) bessely simpbesseli (old name besselisimp) besseli simpbesselk (old name besselksimp) besselk The new routines are: $hankel_1 simphankel1 $hankel_2 simphankel2 If have further introduced the nouns %hankel_1 und %hankel_2. I would prefer to call the functions $hankel_1 and $hankel_2 and not bessel_hankel. The reason is that these functions are well known as Hankel functions. Because I have collected the routines in a new file it would be a good idea to do the changes in the orginal file bessel.lisp. Than it easier to produce a diff for you. I do this the next time. But first, I would like to have a second look at the global arrays $besselarray, $yarray and $iarray. I think it is not so easy to manage these global arrays in a consistent and predictable way for the user. One reason is, that the code for the numerical calculations use the other functions too (i.e. to calculate bessely we need besselj or for the calculation of besseli we need besselj and besselk). So $besselarray can be destroyed when calculating besseli. It might be possible to prevent this effect by introducing a global flag which signals the internal call of one of the routines. A second point is, that we often dont't use the numerical slatec routines directly. In this case it is necessary to do a lot of extra calculations to obtain the values for the arrays. Perhaps it is useful to decide to fill the global arrays only for the case when we can use the values of the slatec routines directly. This is possible for positive order and argument or positiv order and complex argument of the Bessel functions. Further, an additional flag could be introduced to switch the filling of the arrays on and off. It is also to decide what we do when the user calculate a new value for a Bessel function and the calculation don't fill the global array. In this case, the user may be confused by the fact that the global arrays still hold old values. But for which order and argument? For this case it may be helpful to invalidate the array every time we calculate a new value or to store additionally the order and argument for which the values in the global array are calculated. All these problems with the global arrays arise because of the extension of the routines to negative arguments and orders. In my opinion it would be the best to use the concept of the global arrays no longer. The code becomes much more complicated and the benefit for the user might be small. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080328 17:37 Message: Logged In: YES user_id=28849 Originator: NO First, thank you very much for the fixes to the Bessel routines. I'm sure we all want the routines to be correct! Second, about $besselarray. I see the documentation for that is gone, but I think the intent for besselarray was to let the user have access to all the computed values, since some algorithms end up computing all the intermediate orders to get to the desired order. We need to think if this should go away or not. Third, it would certainly help me a lot if you actually produced a diff between your new code and the existing code. Right now, I have to go look at every single function in different places to figure out what's changed. Trying to minimize the changes would help. I know this is a lot of work for you too. Just identifying which function you actually changed would help a lot too. About the specific changes you mention in the comments. No problem with changing besselxsimp to simpbesselx. (I think most of the code uses simpbesselx, but that's too messy for my tastes.) Extending the range is very, very good. Not sure about the special tests for pure real or imaginary answers. Need to look at that some more. Definitions for Hankel functions are good. May want to name them bessel_hankel_1 to be consistent with other Bessel functions. Certainly want to add some properties like derivatives, but that can wait. All in all, I would very much like to take in your very nice changes.  Comment By: Crategus (crategus) Date: 20080327 17:50 Message: Logged In: YES user_id=2039760 Originator: YES I have tested the Bessel functions with the numerical data attached to this message. For the tests I used a function TEST_BESSEL(value, result, digits) which allow to compare the value with the result within the specified digits. Crategus File Added: mytest_bessel.mac  Comment By: Crategus (crategus) Date: 20080327 17:21 Message: Logged In: YES user_id=2039760 Originator: YES I have finished a first redesign of the code for the numerical calculation of the Bessel functions. I have extended the calculation to or changed parts of the calculation to negative orders and arguments. Perhaps the ideas and the code are interesting for the project. I have added the code to this message. A short description of the changes can be find in the header of the file. Crategus File Added: besselnew.lisp  Comment By: Crategus (crategus) Date: 20080323 18:16 Message: Logged In: YES user_id=2039760 Originator: YES I have added to this posting the code for a routine besselj which handles the above problems. The global array $BESSELARRAY is not used. Futher, I added numerical test values for the the special cases in the code. I used Maxima 5.14.0 and GCL 2.6.6 ANSI to compile the code. The testsuite runs with no unexpected errors found. 22 of the 29 numerical tests I have added give the result within a presision of 16 digits. 7 results have a precision of about 14 to 15 digits. Perhaps, the code is interesting enough for further investigation of the numerical evaluation of the Bessel functions. Crategus File Added: besselj.lisp  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 
From: SourceForge.net <noreply@so...>  20080412 05:49:49

Bugs item #1940411, was opened at 20080411 14:03 Message generated for change (Comment added) made by willisbl You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1940411&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: None Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Nobody/Anonymous (nobody) Assigned to: Nobody/Anonymous (nobody) Summary: problem with %pi in sin() Initial Comment: f(x):=sin(4*%pi*x) y:f(x) then, drawing the graph, gets a corrupted sin wave. This problem does not occur if using 2 or 3 intead of 4 in the sin argument, and also it does not occur if typing 3.14 instead of using the %pi constant. Tested using version 0.7.4 of Maxima under Win.  >Comment By: Barton Willis (willisbl) Date: 20080412 00:49 Message: Logged In: YES user_id=895922 Originator: NO I'm guessing that you used wxMaxima 0.74 and used the ploting > plot2d menu item. The default for the number of ticks was too small in that version. In the menu box, look for "ticks" and increase it to say 500 or more. I think there is a wxMaxima 0.74a that has a larger value for the default for ticks. Maybe you can look for that. wxMaxima 0.74sa might be part of Maxima 5.14.0 now  I think some Maxima versions 5.14.0 have wxMaxima 0.74 and others 0.74a.  Comment By: Raymond Toy (rtoy) Date: 20080411 14:49 Message: Logged In: YES user_id=28849 Originator: NO Show exactly what you did. In what way is the sine wave corrupted? What does build_info() say? Are you really using Maxima 0.7.4? Is that the wxmaxima version? I don't recall any 0.7.4 version of Maxima.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1940411&group_id=4933 
From: SourceForge.net <noreply@so...>  20080411 19:49:29

Bugs item #1940411, was opened at 20080411 15:03 Message generated for change (Comment added) made by rtoy You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1940411&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: None Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Nobody/Anonymous (nobody) Assigned to: Nobody/Anonymous (nobody) Summary: problem with %pi in sin() Initial Comment: f(x):=sin(4*%pi*x) y:f(x) then, drawing the graph, gets a corrupted sin wave. This problem does not occur if using 2 or 3 intead of 4 in the sin argument, and also it does not occur if typing 3.14 instead of using the %pi constant. Tested using version 0.7.4 of Maxima under Win.  >Comment By: Raymond Toy (rtoy) Date: 20080411 15:49 Message: Logged In: YES user_id=28849 Originator: NO Show exactly what you did. In what way is the sine wave corrupted? What does build_info() say? Are you really using Maxima 0.7.4? Is that the wxmaxima version? I don't recall any 0.7.4 version of Maxima.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1940411&group_id=4933 
From: SourceForge.net <noreply@so...>  20080411 19:03:34

Bugs item #1940411, was opened at 20080411 12:03 Message generated for change (Tracker Item Submitted) made by Item Submitter You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1940411&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: None Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Nobody/Anonymous (nobody) Assigned to: Nobody/Anonymous (nobody) Summary: problem with %pi in sin() Initial Comment: f(x):=sin(4*%pi*x) y:f(x) then, drawing the graph, gets a corrupted sin wave. This problem does not occur if using 2 or 3 intead of 4 in the sin argument, and also it does not occur if typing 3.14 instead of using the %pi constant. Tested using version 0.7.4 of Maxima under Win.  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1940411&group_id=4933 
From: SourceForge.net <noreply@so...>  20080411 14:17:51

Bugs item #1920177, was opened at 20080319 17:07 Message generated for change (Comment added) made by rtoy You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core Group: Includes proposed fix Status: Open Resolution: None Priority: 5 Private: No Submitted By: Crategus (crategus) Assigned to: Nobody/Anonymous (nobody) Summary: Problems with the bessel functions Initial Comment: There are several problems with the bessel functions. The problems I have found can be divided in the following categories: 1. Inconsistent use of the internal arrays $besselarray, $yarray and $iarray. Example: Restart Maxima and Enter bessel_y(2,1.0). You get an Lisp error. Try first bessel_y(2,1.0) and then repeat bessel_y(2,1.0). Now you get the correct result 1.65...  %i * 3,30... Try bessel_j(2,1.0), you get 0.1149...  %i*2.8142... Next enter bessel_y(2,1.0), you get a Lisp error. The reason for the problems is that the routine bessely uses the global array $YARRAY to store values, but uses also the array $BESSELARRAY to calculate the answer. This is an error. I think the best is to eliminate the use of the global arrays completly. 2. Problematic roundoff errors: Try bessel_j(2,1.0). The result is 0,114903...  %i* 2.8142... e17. The correct result is pure real. The problem is the use of the transformation j[n]=exp(n*%pi*%i)*j[n](x) in the code which produce a small imaginary part. This is a roundoff error for an integer order. In the case of non integer values the imaginary part is correct. So you can not cut off the imaginary part in general. Here is a piece of code which will return an answer which is pure real when the order is an integer (I started to rewrite the code, introduced a function besselj and eliminated the use of the global arrays): (if (integerp order) (if (evenp order) (aref jvals n) ( (aref jvals n)) ) (let ((v (* (cis (* order pi)) (aref jvals n)))) (simplifya `((mplus) ,(realpart v) ((mtimes) $%i ,(imagpart v))) nil ) ) ) 3. Wrong mathematic: Try bessel_j(2.5,1,0). You get 0.04949... Correct is the result 2.8763... * %i Or bessel_j(3,2.0). You get 0,128943... Correct is the result 0.128943... The calculation of the bessel function j[n](x) for real argument x and negativ order n as (realpart (hankel1 order arg)) is not correct. The correct result can be obtained with the formula j[n](x) = 1/2 * (H1[n](x) + H2[n](x)). I tried the following code: (let ((result (* 0.5d0 (+ (hankel1 order arg) (hankel2 order arg))))) (complexify result) ; Problem: you get a small imaginary part ;in the case of a real result (like Problem 2) ; An alternative with a function fpround which round the result ; so that a small imaginary part will vanish ; ; (simplifya ; `((mplus) ; ,(fpround (realpart result) 14) ; ((mtimes) ; ,(fpround (imagpart result) 14) ; $%i ; ) ; ) ; nil ; ) ) This is the definition of the function fpround: (defun fpround (x &optional (digits 1)) (let ((fac (expt 10 digits))) (/ (round x (/ 1 fac)) (float fac)) ) ) I have redesigned a lot of the code for the bessel functions but the work is not finished. Is the work interesting for the project? I use GCL 2.6.6 on Windows XP for programing. Crategus  >Comment By: Raymond Toy (rtoy) Date: 20080411 10:17 Message: Logged In: YES user_id=28849 Originator: NO Changes checked in to bessel.lisp, rev 1.50.  Comment By: Raymond Toy (rtoy) Date: 20080410 12:57 Message: Logged In: YES user_id=28849 Originator: NO I have now integrated your changes into bessel.lisp. The testsuite passes and your bessel tests pass too. I changed test_bessel(a,b,d) to mean abs(ab)<10^(d), which I hope is close enough. (Perhaps the real and imaginary parts should be tested separately?) I'll check in these changes after some code clean up.  Comment By: Crategus (crategus) Date: 20080401 17:02 Message: Logged In: YES user_id=2039760 Originator: YES Again, thank you very much for your interest. I have not changed any of the derivatives. I have only changed the routines which I have listed in this thread. But I have not used the last CVS file. I have used the file which was send with the Maxima 5.14.0 release. There was one change for the derviative for bessel_k since the last release. On Sunday I have installed a tool to download the CVS files. Now I can get the really actual files. The next time I would like to optimize the routines further. Additionally, I am preparing extensive numerical tests for the functions. At the moment I have no idea how to extent the calculation to complex order. In my textbooks and in the book of Abramowitz and Stegun I can't find a relationship which help to calculate the Bessel functions for complex order. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080401 12:41 Message: Logged In: YES user_id=28849 Originator: NO Thank you very much for the updated files. (The diff should have used diff u or diff c, but no matter. The besselchanged.lisp file is perfect.) Given the imminent release of Maxima, I will wait on these changes until after the release. Then I'll try to incorporate your changes. I noticed that you changed the derivatives of some of the Bessel functions. Was that intentional? Were they wrong? (Some of the derivatives need to be in a certain form for some simplifications to be done.)  Comment By: Crategus (crategus) Date: 20080329 17:05 Message: Logged In: YES user_id=2039760 Originator: YES Next I have added a diff between the orginal and the changed file. Question: Is it possible to attach at once more than one file to a message? Crategus File Added: diff.txt  Comment By: Crategus (crategus) Date: 20080329 17:00 Message: Logged In: YES user_id=2039760 Originator: YES I have put all changes of the code for the Bessel functions in the orginal file bessel.lisp which is distributed with Maxima 5.14.0 and attached this file to this message. I have tested the code with the testsuite and the values of the mytest_bessel.mac. Crategus File Added: besselchanged.lisp  Comment By: Crategus (crategus) Date: 20080328 20:35 Message: Logged In: YES user_id=2039760 Originator: YES Thank you very much for your answer. I have redesigned the following routines: simpbesselj (old name besseljsimp) besselj (old name $bessel) simpbessely (old name besselysimp) bessely simpbesseli (old name besselisimp) besseli simpbesselk (old name besselksimp) besselk The new routines are: $hankel_1 simphankel1 $hankel_2 simphankel2 If have further introduced the nouns %hankel_1 und %hankel_2. I would prefer to call the functions $hankel_1 and $hankel_2 and not bessel_hankel. The reason is that these functions are well known as Hankel functions. Because I have collected the routines in a new file it would be a good idea to do the changes in the orginal file bessel.lisp. Than it easier to produce a diff for you. I do this the next time. But first, I would like to have a second look at the global arrays $besselarray, $yarray and $iarray. I think it is not so easy to manage these global arrays in a consistent and predictable way for the user. One reason is, that the code for the numerical calculations use the other functions too (i.e. to calculate bessely we need besselj or for the calculation of besseli we need besselj and besselk). So $besselarray can be destroyed when calculating besseli. It might be possible to prevent this effect by introducing a global flag which signals the internal call of one of the routines. A second point is, that we often dont't use the numerical slatec routines directly. In this case it is necessary to do a lot of extra calculations to obtain the values for the arrays. Perhaps it is useful to decide to fill the global arrays only for the case when we can use the values of the slatec routines directly. This is possible for positive order and argument or positiv order and complex argument of the Bessel functions. Further, an additional flag could be introduced to switch the filling of the arrays on and off. It is also to decide what we do when the user calculate a new value for a Bessel function and the calculation don't fill the global array. In this case, the user may be confused by the fact that the global arrays still hold old values. But for which order and argument? For this case it may be helpful to invalidate the array every time we calculate a new value or to store additionally the order and argument for which the values in the global array are calculated. All these problems with the global arrays arise because of the extension of the routines to negative arguments and orders. In my opinion it would be the best to use the concept of the global arrays no longer. The code becomes much more complicated and the benefit for the user might be small. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080328 12:37 Message: Logged In: YES user_id=28849 Originator: NO First, thank you very much for the fixes to the Bessel routines. I'm sure we all want the routines to be correct! Second, about $besselarray. I see the documentation for that is gone, but I think the intent for besselarray was to let the user have access to all the computed values, since some algorithms end up computing all the intermediate orders to get to the desired order. We need to think if this should go away or not. Third, it would certainly help me a lot if you actually produced a diff between your new code and the existing code. Right now, I have to go look at every single function in different places to figure out what's changed. Trying to minimize the changes would help. I know this is a lot of work for you too. Just identifying which function you actually changed would help a lot too. About the specific changes you mention in the comments. No problem with changing besselxsimp to simpbesselx. (I think most of the code uses simpbesselx, but that's too messy for my tastes.) Extending the range is very, very good. Not sure about the special tests for pure real or imaginary answers. Need to look at that some more. Definitions for Hankel functions are good. May want to name them bessel_hankel_1 to be consistent with other Bessel functions. Certainly want to add some properties like derivatives, but that can wait. All in all, I would very much like to take in your very nice changes.  Comment By: Crategus (crategus) Date: 20080327 12:50 Message: Logged In: YES user_id=2039760 Originator: YES I have tested the Bessel functions with the numerical data attached to this message. For the tests I used a function TEST_BESSEL(value, result, digits) which allow to compare the value with the result within the specified digits. Crategus File Added: mytest_bessel.mac  Comment By: Crategus (crategus) Date: 20080327 12:21 Message: Logged In: YES user_id=2039760 Originator: YES I have finished a first redesign of the code for the numerical calculation of the Bessel functions. I have extended the calculation to or changed parts of the calculation to negative orders and arguments. Perhaps the ideas and the code are interesting for the project. I have added the code to this message. A short description of the changes can be find in the header of the file. Crategus File Added: besselnew.lisp  Comment By: Crategus (crategus) Date: 20080323 13:16 Message: Logged In: YES user_id=2039760 Originator: YES I have added to this posting the code for a routine besselj which handles the above problems. The global array $BESSELARRAY is not used. Futher, I added numerical test values for the the special cases in the code. I used Maxima 5.14.0 and GCL 2.6.6 ANSI to compile the code. The testsuite runs with no unexpected errors found. 22 of the 29 numerical tests I have added give the result within a presision of 16 digits. 7 results have a precision of about 14 to 15 digits. Perhaps, the code is interesting enough for further investigation of the numerical evaluation of the Bessel functions. Crategus File Added: besselj.lisp  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 
From: SourceForge.net <noreply@so...>  20080411 07:01:51

Bugs item #1916517, was opened at 20080317 18:40 Message generated for change (Comment added) made by satoshi_adachi You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1916517&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: None Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Satoshi Adachi (satoshi_adachi) Assigned to: Nobody/Anonymous (nobody) Summary: (F=0 and G=0) and (F=0 and FG=0) have different solutions. Initial Comment: Dear Developers of Maxima, I met a set of two albebraic equations F = 0 and G = 0 of two variables x and y such that (i) solve([F=0,G=0],[x,y]) returns no solution [] (this is wrong), whereas (ii) solve([F=0,FG=0],[x,y]) returns true solutions. This is strange, since the first set of equations F=0 and G=0 should have the same solutions of the second set of equations F=0 and FG=0. (1) The Maxima program to solve F=0 and FG=0 is  /* * solve_system_of_algebraic_equations_OK.maxima: * * F is a polynomial of x and y, which is degree 3 in x and 2 in y. * G is a polynomial of x and y, which is degree 3 in x and 2 in y. * * solve() and algsys() can solve the system equations F=0 and FG=0. */ display2d:false; F:2*x^3*y^2+2*a*x^2*y^2+3*x^2*y+2*a*x*y+x+1; G:x^3*y^2+2*a*x^2*y^2+a^2*x*y^2+x^2*y+2*a*x*y+a^2*y+1; solutions:solve([F=0,FG=0],[x,y]); /* OK: solutions are found. */  The result is  Maxima 5.14.0cvs http://maxima.sourceforge.net Using Lisp GNU Common Lisp (GCL) GCL 2.6.7 (aka GCL) Distributed under the GNU Public License. See the file COPYING. Dedicated to the memory of William Schelter. The function bug_report() provides bug reporting information. (%i1) batch(solve_system_of_algebraic_equations_OK.maxima) batching #p/Volumes/HFS+2/home/adachi/work/285/solve_system_of_algebraic_equations_OK.maxima (%i2) display2d : false (%o2) false (%i3) F:1+x+2*a*x*y+3*x^2*y+2*a*x^2*y^2+2*x^3*y^2 (%o3) 2*x^3*y^2+2*a*x^2*y^2+3*x^2*y+2*a*x*y+x+1 (%i4) G:1+a^2*y+2*a*x*y+x^2*y+a^2*x*y^2+2*a*x^2*y^2+x^3*y^2 (%o4) x^3*y^2+2*a*x^2*y^2+a^2*x*y^2+x^2*y+2*a*x*y+a^2*y+1 (%i5) solutions:solve([F = 0,FG = 0],[x,y]) (%o5) [[x = (sqrt(a)*(a^24*a)+sqrt(a4)*(2*aa^2))/(2*sqrt(a4)), y = 1/(a*sqrt(a^24*a))], [x = (sqrt(a4)*(a^22*a)+sqrt(a)*(a^24*a))/(2*sqrt(a4)), y = 1/(a*sqrt(a^24*a))]] (%o6) "solve_system_of_algebraic_equations_OK.maxima"  =============================================================================== (2) The Maxima program to solve F=0 and G=0 is  /* * solve_system_of_algebraic_equations_NG.maxima: * * F is a polynomial of x and y, which is degree 3 in x and 2 in y. * G is a polynomial of x and y, which is degree 3 in x and 2 in y. * * solve() and algsys() can NOT solve the system equations F=0 and G=0. */ display2d:false; F:2*x^3*y^2+2*a*x^2*y^2+3*x^2*y+2*a*x*y+x+1; G:x^3*y^2+2*a*x^2*y^2+a^2*x*y^2+x^2*y+2*a*x*y+a^2*y+1; solutions:solve([F=0,G=0],[x,y]); /* NG: solutions are NOT found! */  The result is  Maxima 5.14.0cvs http://maxima.sourceforge.net Using Lisp GNU Common Lisp (GCL) GCL 2.6.7 (aka GCL) Distributed under the GNU Public License. See the file COPYING. Dedicated to the memory of William Schelter. The function bug_report() provides bug reporting information. (%i1) batch(solve_system_of_algebraic_equations_NG.maxima) batching #p/Volumes/HFS+2/home/adachi/work/285/solve_system_of_algebraic_equations_NG.maxima (%i2) display2d : false (%o2) false (%i3) F:1+x+2*a*x*y+3*x^2*y+2*a*x^2*y^2+2*x^3*y^2 (%o3) 2*x^3*y^2+2*a*x^2*y^2+3*x^2*y+2*a*x*y+x+1 (%i4) G:1+a^2*y+2*a*x*y+x^2*y+a^2*x*y^2+2*a*x^2*y^2+x^3*y^2 (%o4) x^3*y^2+2*a*x^2*y^2+a^2*x*y^2+x^2*y+2*a*x*y+a^2*y+1 (%i5) solutions:solve([F = 0,G = 0],[x,y]) (%o5) [] (%o6) "solve_system_of_algebraic_equations_NG.maxima"  =============================================================================== (3) I also tried algebraic:true to solve F=0 and G=0. The program is  /* * solve_system_of_algebraic_equations_with_algebraic_NG.maxima: * * F is a polynomial of x and y, which is degree 3 in x and 2 in y. * G is a polynomial of x and y, which is degree 3 in x and 2 in y. * * solve() and algsys() can NOT solve the system equations F=0 and G=0. */ display2d:false; algebraic:true; F:2*x^3*y^2+2*a*x^2*y^2+3*x^2*y+2*a*x*y+x+1; G:x^3*y^2+2*a*x^2*y^2+a^2*x*y^2+x^2*y+2*a*x*y+a^2*y+1; solutions:solve([F=0,G=0],[x,y]); /* NG: internal error! */  The result is further strange for me as  Maxima 5.14.0cvs http://maxima.sourceforge.net Using Lisp GNU Common Lisp (GCL) GCL 2.6.7 (aka GCL) Distributed under the GNU Public License. See the file COPYING. Dedicated to the memory of William Schelter. The function bug_report() provides bug reporting information. (%i1) batch(solve_system_of_algebraic_equations_with_algebraic.maxima) batching #p/Volumes/HFS+2/home/adachi/work/285/solve_system_of_algebraic_equations_with_algebraic.maxima (%i2) display2d : false (%o2) false (%i3) algebraic:true (%o3) true (%i4) F:1+x+2*a*x*y+3*x^2*y+2*a*x^2*y^2+2*x^3*y^2 (%o4) 2*x^3*y^2+2*a*x^2*y^2+3*x^2*y+2*a*x*y+x+1 (%i5) G:1+a^2*y+2*a*x*y+x^2*y+a^2*x*y^2+2*a*x^2*y^2+x^3*y^2 (%o5) x^3*y^2+2*a*x^2*y^2+a^2*x*y^2+x^2*y+2*a*x*y+a^2*y+1 (%i6) solutions:solve([F = 0,G = 0],[x,y]) Polynomial quotient is not exact  an error. To debug this try debugmode(true);  Sincerely yours, Satoshi Adachi  >Comment By: Satoshi Adachi (satoshi_adachi) Date: 20080411 16:01 Message: Logged In: YES user_id=1953419 Originator: YES Dear Mr.Barton Willis (willisbl), Thank you very much for educating me to know how the decifiancy of algsys() that was described in my previous bug report can be overcome. Your information is very helpful! At first, I am sorry to be late in writing this letter of thanks. I had been occupied by duty works until yesterday. Today, I wrote a wrapper routine of algsys(), which is based fully on your suggestion and is named algsys_willisbl(), as follows:  /* * algsys_willisbl.maxima: * * better algsys() suggested by willisbl. * * [Remark] To use algsys_willisbl(), the following preparation is necessary: * * load("topoly"); // for elim() * algebraic:true; * * S.Adachi 2008/04/10 */ algsys_willisbl(eqns_lst, vars_lst) := block([eqns1_lst,eqns2_lst,sols_cand_lst,sols_lst], eqns1_lst:map(lambda([e], if (not atom(e) and operatorp(e,"=")) then lhs(e)rhs(e) else e), eqns_lst), eqns2_lst:first(rest(elim(eqns1_lst, vars_lst))), sols_cand_lst:algsys(eqns2_lst, vars_lst), sols_lst:listify(subset(setify(sols_cand_lst), lambda([sc], every(lambda([a], is(equal(a,0))), ratsimp(subst(sc, eqns1_lst)))))), sols_lst); /* END */  I also wrote a test program for algsys_willisbl() by using the same system of algebraic equations that was described in my previous bug report of algsys() as  /* * test_algsys_willisbl.maxima: * * test algsys_willisbl() by applying it to a system of two algebraic * equations F=0 and G=0, which cannot be solved by algsys(). * * S.Adachi 2008/04/10 */ load(topoly); load("algsys_willisbl.maxima"); algebraic:true; display2d:false; F:2*x^3*y^2+2*a*x^2*y^2+3*x^2*y+2*a*x*y+x+1; G:x^3*y^2+2*a*x^2*y^2+a^2*x*y^2+x^2*y+2*a*x*y+a^2*y+1; sols1:algsys_willisbl([F,G], [x,y]); sols2:algsys_willisbl([F=0,G=0], [x,y]); sols3:algsys_willisbl([F+1=1,G+2=2], [x,y]); /* END */  The result of executation is as follows:  [Macintosh:~/work/288:1] adachi% maxima b test_algsys_willisbl.maxima Maxima 5.14.0cvs http://maxima.sourceforge.net Using Lisp GNU Common Lisp (GCL) GCL 2.6.7 (aka GCL) Distributed under the GNU Public License. See the file COPYING. Dedicated to the memory of William Schelter. The function bug_report() provides bug reporting information. (%i1) batch(test_algsys_willisbl.maxima) batching #p/Volumes/HFS+2/home/adachi/work/288/test_algsys_willisbl.maxima (%i2) load(topoly) (%o2) /usr/local/share/maxima/5.14.0cvs/share/contrib/topoly.lisp (%i3) load(algsys_willisbl.maxima) (%o3) algsys_willisbl.maxima (%i4) algebraic : true (%o4) true (%i5) display2d : false (%o5) false (%i6) F:1+x+2*a*x*y+3*x^2*y+2*a*x^2*y^2+2*x^3*y^2 (%o6) 2*x^3*y^2+2*a*x^2*y^2+3*x^2*y+2*a*x*y+x+1 (%i7) G:1+a^2*y+2*a*x*y+x^2*y+a^2*x*y^2+2*a*x^2*y^2+x^3*y^2 (%o7) x^3*y^2+2*a*x^2*y^2+a^2*x*y^2+x^2*y+2*a*x*y+a^2*y+1 (%i8) sols1:algsys_willisbl([F,G],[x,y]) (%o8) [[x = (a*sqrt(a^24*a)a^2+2*a)/2,y = sqrt(a^24*a)/(a^34*a^2)], [x = (a*sqrt(a^24*a)+a^22*a)/2,y = sqrt(a^24*a)/(a^34*a^2)]] (%i9) sols2:algsys_willisbl([F = 0,G = 0],[x,y]) (%o9) [[x = (a*sqrt(a^24*a)a^2+2*a)/2,y = sqrt(a^24*a)/(a^34*a^2)], [x = (a*sqrt(a^24*a)+a^22*a)/2,y = sqrt(a^24*a)/(a^34*a^2)]] (%i10) sols3:algsys_willisbl([1+F = 1,2+G = 2],[x,y]) (%o10) [[x = (a*sqrt(a^24*a)a^2+2*a)/2,y = sqrt(a^24*a)/(a^34*a^2)], [x = (a*sqrt(a^24*a)+a^22*a)/2,y = sqrt(a^24*a)/(a^34*a^2)]] (%o11) "test_algsys_willisbl.maxima" [Macintosh:~/work/288:2] adachi%  I did not apply algsys_willisbl() to any other examples, yet. So, it may have still bugs. Anyway, I become now able to advance my study again without the annoyance caused by algsys(). The knowledge that you gave me is very essential to use algsys() in a better way. Thank you very much. Sincerely yours, Satoshi Adachi  Comment By: Barton Willis (willisbl) Date: 20080317 21:05 Message: Logged In: YES user_id=895922 Originator: NO Thanks for reporting this bug. The function algsys isn't nearly as good as we would like. The best I can offer is my favorite workaround (that might work); it's something like: (%i92) load(topoly)$ (%i93) algebraic : true$ Use 'elim' to triangularize the equationsthis can create spurious solutions: (%i94) elim([F,G],[x,y]); (%o94) [[],[x^4+a^2*x^3a*x^3+a^3*x^2a^2*x^2+a^3*x,x^3*y^2+2*a*x^2*y^2+a^2*x*y^2+x^2*y+2*a*x*y+a^2*y+1]] Use algsys to solve these equations: (%i95) sol : algsys(first(rest(%)),[x,y])$ Try to reject the spurious solutions: (%i96) true_sol : set()$ (%i97) for si in sol do if every(lambda([a], is(equal(a,0))), ratsimp(subst(si,[F,G]))) then true_sol : adjoin(si, true_sol); (%o97) done (%i98) true_sol; (%o98) {[x=(a*sqrt(a^24*a)a^2+2*a)/2,y=sqrt(a^24*a)/(a^34*a^2)], [x=(a*sqrt(a^24*a)+a^22*a)/2,y=sqrt(a^24*a)/(a^34*a^2)]} Two solutions *might* be spurious and two checked out OK. You'll need to inspect the two rejected solutions to see if they are really spurious. And a = 0 and a = 4 might be special cases that need to be checked. Of course, all this should be automatic; unfortunately, algsys isn't there yet. Barton  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1916517&group_id=4933 
From: SourceForge.net <noreply@so...>  20080410 16:57:52

Bugs item #1920177, was opened at 20080319 17:07 Message generated for change (Comment added) made by rtoy You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 Please note that this message will contain a full copy of the comment thread, including the initial issue submission, for this request, not just the latest update. Category: Lisp Core Group: Includes proposed fix Status: Open Resolution: None Priority: 5 Private: No Submitted By: Crategus (crategus) Assigned to: Nobody/Anonymous (nobody) Summary: Problems with the bessel functions Initial Comment: There are several problems with the bessel functions. The problems I have found can be divided in the following categories: 1. Inconsistent use of the internal arrays $besselarray, $yarray and $iarray. Example: Restart Maxima and Enter bessel_y(2,1.0). You get an Lisp error. Try first bessel_y(2,1.0) and then repeat bessel_y(2,1.0). Now you get the correct result 1.65...  %i * 3,30... Try bessel_j(2,1.0), you get 0.1149...  %i*2.8142... Next enter bessel_y(2,1.0), you get a Lisp error. The reason for the problems is that the routine bessely uses the global array $YARRAY to store values, but uses also the array $BESSELARRAY to calculate the answer. This is an error. I think the best is to eliminate the use of the global arrays completly. 2. Problematic roundoff errors: Try bessel_j(2,1.0). The result is 0,114903...  %i* 2.8142... e17. The correct result is pure real. The problem is the use of the transformation j[n]=exp(n*%pi*%i)*j[n](x) in the code which produce a small imaginary part. This is a roundoff error for an integer order. In the case of non integer values the imaginary part is correct. So you can not cut off the imaginary part in general. Here is a piece of code which will return an answer which is pure real when the order is an integer (I started to rewrite the code, introduced a function besselj and eliminated the use of the global arrays): (if (integerp order) (if (evenp order) (aref jvals n) ( (aref jvals n)) ) (let ((v (* (cis (* order pi)) (aref jvals n)))) (simplifya `((mplus) ,(realpart v) ((mtimes) $%i ,(imagpart v))) nil ) ) ) 3. Wrong mathematic: Try bessel_j(2.5,1,0). You get 0.04949... Correct is the result 2.8763... * %i Or bessel_j(3,2.0). You get 0,128943... Correct is the result 0.128943... The calculation of the bessel function j[n](x) for real argument x and negativ order n as (realpart (hankel1 order arg)) is not correct. The correct result can be obtained with the formula j[n](x) = 1/2 * (H1[n](x) + H2[n](x)). I tried the following code: (let ((result (* 0.5d0 (+ (hankel1 order arg) (hankel2 order arg))))) (complexify result) ; Problem: you get a small imaginary part ;in the case of a real result (like Problem 2) ; An alternative with a function fpround which round the result ; so that a small imaginary part will vanish ; ; (simplifya ; `((mplus) ; ,(fpround (realpart result) 14) ; ((mtimes) ; ,(fpround (imagpart result) 14) ; $%i ; ) ; ) ; nil ; ) ) This is the definition of the function fpround: (defun fpround (x &optional (digits 1)) (let ((fac (expt 10 digits))) (/ (round x (/ 1 fac)) (float fac)) ) ) I have redesigned a lot of the code for the bessel functions but the work is not finished. Is the work interesting for the project? I use GCL 2.6.6 on Windows XP for programing. Crategus  >Comment By: Raymond Toy (rtoy) Date: 20080410 12:57 Message: Logged In: YES user_id=28849 Originator: NO I have now integrated your changes into bessel.lisp. The testsuite passes and your bessel tests pass too. I changed test_bessel(a,b,d) to mean abs(ab)<10^(d), which I hope is close enough. (Perhaps the real and imaginary parts should be tested separately?) I'll check in these changes after some code clean up.  Comment By: Crategus (crategus) Date: 20080401 17:02 Message: Logged In: YES user_id=2039760 Originator: YES Again, thank you very much for your interest. I have not changed any of the derivatives. I have only changed the routines which I have listed in this thread. But I have not used the last CVS file. I have used the file which was send with the Maxima 5.14.0 release. There was one change for the derviative for bessel_k since the last release. On Sunday I have installed a tool to download the CVS files. Now I can get the really actual files. The next time I would like to optimize the routines further. Additionally, I am preparing extensive numerical tests for the functions. At the moment I have no idea how to extent the calculation to complex order. In my textbooks and in the book of Abramowitz and Stegun I can't find a relationship which help to calculate the Bessel functions for complex order. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080401 12:41 Message: Logged In: YES user_id=28849 Originator: NO Thank you very much for the updated files. (The diff should have used diff u or diff c, but no matter. The besselchanged.lisp file is perfect.) Given the imminent release of Maxima, I will wait on these changes until after the release. Then I'll try to incorporate your changes. I noticed that you changed the derivatives of some of the Bessel functions. Was that intentional? Were they wrong? (Some of the derivatives need to be in a certain form for some simplifications to be done.)  Comment By: Crategus (crategus) Date: 20080329 17:05 Message: Logged In: YES user_id=2039760 Originator: YES Next I have added a diff between the orginal and the changed file. Question: Is it possible to attach at once more than one file to a message? Crategus File Added: diff.txt  Comment By: Crategus (crategus) Date: 20080329 17:00 Message: Logged In: YES user_id=2039760 Originator: YES I have put all changes of the code for the Bessel functions in the orginal file bessel.lisp which is distributed with Maxima 5.14.0 and attached this file to this message. I have tested the code with the testsuite and the values of the mytest_bessel.mac. Crategus File Added: besselchanged.lisp  Comment By: Crategus (crategus) Date: 20080328 20:35 Message: Logged In: YES user_id=2039760 Originator: YES Thank you very much for your answer. I have redesigned the following routines: simpbesselj (old name besseljsimp) besselj (old name $bessel) simpbessely (old name besselysimp) bessely simpbesseli (old name besselisimp) besseli simpbesselk (old name besselksimp) besselk The new routines are: $hankel_1 simphankel1 $hankel_2 simphankel2 If have further introduced the nouns %hankel_1 und %hankel_2. I would prefer to call the functions $hankel_1 and $hankel_2 and not bessel_hankel. The reason is that these functions are well known as Hankel functions. Because I have collected the routines in a new file it would be a good idea to do the changes in the orginal file bessel.lisp. Than it easier to produce a diff for you. I do this the next time. But first, I would like to have a second look at the global arrays $besselarray, $yarray and $iarray. I think it is not so easy to manage these global arrays in a consistent and predictable way for the user. One reason is, that the code for the numerical calculations use the other functions too (i.e. to calculate bessely we need besselj or for the calculation of besseli we need besselj and besselk). So $besselarray can be destroyed when calculating besseli. It might be possible to prevent this effect by introducing a global flag which signals the internal call of one of the routines. A second point is, that we often dont't use the numerical slatec routines directly. In this case it is necessary to do a lot of extra calculations to obtain the values for the arrays. Perhaps it is useful to decide to fill the global arrays only for the case when we can use the values of the slatec routines directly. This is possible for positive order and argument or positiv order and complex argument of the Bessel functions. Further, an additional flag could be introduced to switch the filling of the arrays on and off. It is also to decide what we do when the user calculate a new value for a Bessel function and the calculation don't fill the global array. In this case, the user may be confused by the fact that the global arrays still hold old values. But for which order and argument? For this case it may be helpful to invalidate the array every time we calculate a new value or to store additionally the order and argument for which the values in the global array are calculated. All these problems with the global arrays arise because of the extension of the routines to negative arguments and orders. In my opinion it would be the best to use the concept of the global arrays no longer. The code becomes much more complicated and the benefit for the user might be small. Crategus  Comment By: Raymond Toy (rtoy) Date: 20080328 12:37 Message: Logged In: YES user_id=28849 Originator: NO First, thank you very much for the fixes to the Bessel routines. I'm sure we all want the routines to be correct! Second, about $besselarray. I see the documentation for that is gone, but I think the intent for besselarray was to let the user have access to all the computed values, since some algorithms end up computing all the intermediate orders to get to the desired order. We need to think if this should go away or not. Third, it would certainly help me a lot if you actually produced a diff between your new code and the existing code. Right now, I have to go look at every single function in different places to figure out what's changed. Trying to minimize the changes would help. I know this is a lot of work for you too. Just identifying which function you actually changed would help a lot too. About the specific changes you mention in the comments. No problem with changing besselxsimp to simpbesselx. (I think most of the code uses simpbesselx, but that's too messy for my tastes.) Extending the range is very, very good. Not sure about the special tests for pure real or imaginary answers. Need to look at that some more. Definitions for Hankel functions are good. May want to name them bessel_hankel_1 to be consistent with other Bessel functions. Certainly want to add some properties like derivatives, but that can wait. All in all, I would very much like to take in your very nice changes.  Comment By: Crategus (crategus) Date: 20080327 12:50 Message: Logged In: YES user_id=2039760 Originator: YES I have tested the Bessel functions with the numerical data attached to this message. For the tests I used a function TEST_BESSEL(value, result, digits) which allow to compare the value with the result within the specified digits. Crategus File Added: mytest_bessel.mac  Comment By: Crategus (crategus) Date: 20080327 12:21 Message: Logged In: YES user_id=2039760 Originator: YES I have finished a first redesign of the code for the numerical calculation of the Bessel functions. I have extended the calculation to or changed parts of the calculation to negative orders and arguments. Perhaps the ideas and the code are interesting for the project. I have added the code to this message. A short description of the changes can be find in the header of the file. Crategus File Added: besselnew.lisp  Comment By: Crategus (crategus) Date: 20080323 13:16 Message: Logged In: YES user_id=2039760 Originator: YES I have added to this posting the code for a routine besselj which handles the above problems. The global array $BESSELARRAY is not used. Futher, I added numerical test values for the the special cases in the code. I used Maxima 5.14.0 and GCL 2.6.6 ANSI to compile the code. The testsuite runs with no unexpected errors found. 22 of the 29 numerical tests I have added give the result within a presision of 16 digits. 7 results have a precision of about 14 to 15 digits. Perhaps, the code is interesting enough for further investigation of the numerical evaluation of the Bessel functions. Crategus File Added: besselj.lisp  You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1920177&group_id=4933 