From: SourceForge.net <no...@so...> - 2008-04-30 14:11:43
|
Bugs item #1954846, was opened at 2008-04-30 10:11 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=1954846&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: Open Resolution: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: bessel_i(1/2,0) -> divide by zero error Initial Comment: bessel_i(1/2,0) produces an error. I think it should be 0 since bessel_i(1/2,x) is sqrt(2/%pi)*sinh(x)/sqrt(x). ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&group_id=4933 |
From: SourceForge.net <no...@so...> - 2008-05-02 18:05:36
|
Bugs item #1954846, was opened at 2008-04-30 09:11 Message generated for change (Comment added) made by willisbl You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&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: Open Resolution: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: bessel_i(1/2,0) -> divide by zero error Initial Comment: bessel_i(1/2,0) produces an error. I think it should be 0 since bessel_i(1/2,x) is sqrt(2/%pi)*sinh(x)/sqrt(x). ---------------------------------------------------------------------- >Comment By: Barton Willis (willisbl) Date: 2008-05-02 13:05 Message: Logged In: YES user_id=895922 Originator: NO See A&S 9.1.7: http://www.convertit.com/Go/Convertit/Reference/AMS55.ASP?Res=150&Page=360 For n > 0, bessel_j(n,0) = 0 (n needn't be integer). Maxima doesn't use this simplification for noninteger n. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&group_id=4933 |
From: SourceForge.net <no...@so...> - 2008-06-22 09:45:24
|
Bugs item #1954846, was opened at 2008-04-30 10:11 Message generated for change (Comment added) made by rtoy You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&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: Open Resolution: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: bessel_i(1/2,0) -> divide by zero error Initial Comment: bessel_i(1/2,0) produces an error. I think it should be 0 since bessel_i(1/2,x) is sqrt(2/%pi)*sinh(x)/sqrt(x). ---------------------------------------------------------------------- >Comment By: Raymond Toy (rtoy) Date: 2008-06-20 15:07 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 15:06 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Barton Willis (willisbl) Date: 2008-05-02 14:05 Message: Logged In: YES user_id=895922 Originator: NO See A&S 9.1.7: http://www.convertit.com/Go/Convertit/Reference/AMS55.ASP?Res=150&Page=360 For n > 0, bessel_j(n,0) = 0 (n needn't be integer). Maxima doesn't use this simplification for noninteger n. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&group_id=4933 |
From: SourceForge.net <no...@so...> - 2008-06-22 10:10:11
|
Bugs item #1954846, was opened at 2008-04-30 10:11 Message generated for change (Comment added) made by rtoy You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&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: Open Resolution: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: bessel_i(1/2,0) -> divide by zero error Initial Comment: bessel_i(1/2,0) produces an error. I think it should be 0 since bessel_i(1/2,x) is sqrt(2/%pi)*sinh(x)/sqrt(x). ---------------------------------------------------------------------- >Comment By: Raymond Toy (rtoy) Date: 2008-06-20 15:06 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Barton Willis (willisbl) Date: 2008-05-02 14:05 Message: Logged In: YES user_id=895922 Originator: NO See A&S 9.1.7: http://www.convertit.com/Go/Convertit/Reference/AMS55.ASP?Res=150&Page=360 For n > 0, bessel_j(n,0) = 0 (n needn't be integer). Maxima doesn't use this simplification for noninteger n. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&group_id=4933 |
From: SourceForge.net <no...@so...> - 2008-06-22 10:23:42
|
Bugs item #1954846, was opened at 2008-04-30 10:11 Message generated for change (Comment added) made by rtoy You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&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: Open Resolution: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: bessel_i(1/2,0) -> divide by zero error Initial Comment: bessel_i(1/2,0) produces an error. I think it should be 0 since bessel_i(1/2,x) is sqrt(2/%pi)*sinh(x)/sqrt(x). ---------------------------------------------------------------------- >Comment By: Raymond Toy (rtoy) Date: 2008-06-20 15:07 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 15:07 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 15:06 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Barton Willis (willisbl) Date: 2008-05-02 14:05 Message: Logged In: YES user_id=895922 Originator: NO See A&S 9.1.7: http://www.convertit.com/Go/Convertit/Reference/AMS55.ASP?Res=150&Page=360 For n > 0, bessel_j(n,0) = 0 (n needn't be integer). Maxima doesn't use this simplification for noninteger n. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&group_id=4933 |
From: SourceForge.net <no...@so...> - 2008-07-04 21:10:14
|
Bugs item #1954846, was opened at 2008-04-30 16:11 Message generated for change (Comment added) made by crategus You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&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: Open Resolution: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: bessel_i(1/2,0) -> divide by zero error Initial Comment: bessel_i(1/2,0) produces an error. I think it should be 0 since bessel_i(1/2,x) is sqrt(2/%pi)*sinh(x)/sqrt(x). ---------------------------------------------------------------------- >Comment By: Dieter Kaiser (crategus) Date: 2008-07-04 23:10 Message: Logged In: YES user_id=2039760 Originator: NO If the gloabal flag BESSELEXPAND is FALSE we get for all Bessel functions with a Rational order and a zero as argument a noun form with the CVS code 1.58. (In previous versions a lot of the zero checks were not implemented.) But, if we set the global flag BESSELEXPAND to TRUE the Bessel function will be expanded and then evaluated for the argument zero. We get for all Bessel functions the described error. I think it is a question of convention and the desired behaviour of Maxima what we should do. The code to catch the zero argument is only implemented for CL numbers. The last time the code was extended to look more carefully for zero arguments didn't extend this to Maxima Rational, Bigfloat numbers or even expression which might simplify. In my opinion the reason to look only for CL numbers is, that only for such numbers we do the numerical evaluation. We check for a zero argument to prevent Lisp Errors in the numerical slatec routines. At this point we need to do the simplification for the Bessel function with zero argument. In the case of Rational numbers it depends on the flag BESSELEXPAND what we should do. In my opinion it is the best to check for the zero argument before we actually do the expansion. If we try to simplify all Bessel functions with zero argument in advance we get a lot of problems and questions, too. The reason is that the result of the simplification with an argument of zero depends on the sign of the order and is different for complex and real order of the Bessel function. Here an example what we can do in simp-bessel-i: ((and $besselexpand (setq rat-order (max-numeric-ratio-p order 2))) (cond ((= arg 0) ;; We don't do the expansion for a zero argument. (if (> rat-order 0) 0 '$infinity)) (t ;; When order is a fraction with a denominator of 2, we ;; can express the result in terms of elementary ;; functions. ;; ;; I[1/2](z) = sqrt(2/%pi/z)*sinh(z) ;; I[-1/2](z) = sqrt(2/%pi/z)*cosh(z) (bessel-i-half-order rat-order arg)))) We no longer get the divison by zero error. Here the results: (%i7) besselexpand:false; (%o7) false (%i8) bessel_i(1/2,0); (%o8) bessel_i(1/2,0) (%i9) bessel_i(-1/2,0); (%o9) bessel_i(-1/2,0) (%i10) besselexpand:true; (%o10) true (%i11) bessel_i(1/2,0); (%o11) 0 (%i12) bessel_i(-1/2,0); (%o12) infinity This check for a zero argument before we expand the Bessel function could be implemented for all other Bessel functions too. Dieter Kaiser ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 21:07 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 21:07 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 21:06 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Barton Willis (willisbl) Date: 2008-05-02 20:05 Message: Logged In: YES user_id=895922 Originator: NO See A&S 9.1.7: http://www.convertit.com/Go/Convertit/Reference/AMS55.ASP?Res=150&Page=360 For n > 0, bessel_j(n,0) = 0 (n needn't be integer). Maxima doesn't use this simplification for noninteger n. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&group_id=4933 |
From: SourceForge.net <no...@so...> - 2008-07-04 21:46:03
|
Bugs item #1954846, was opened at 2008-04-30 16:11 Message generated for change (Comment added) made by crategus You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&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: Open Resolution: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: bessel_i(1/2,0) -> divide by zero error Initial Comment: bessel_i(1/2,0) produces an error. I think it should be 0 since bessel_i(1/2,x) is sqrt(2/%pi)*sinh(x)/sqrt(x). ---------------------------------------------------------------------- >Comment By: Dieter Kaiser (crategus) Date: 2008-07-04 23:45 Message: Logged In: YES user_id=2039760 Originator: NO Sorry, but I have omitted to check arg to be an CL number. The correct code would be: ((and $besselexpand (setq rat-order (max-numeric-ratio-p order 2))) (cond ((and (numberp arg) (= arg 0)) (if (> rat-order 0) 0 '$infinity)) ( t ;; When order is a fraction with a denominator of 2, we ;; can express the result in terms of elementary ;; functions. ;; ;; I[1/2](z) = sqrt(2/%pi/z)*sinh(z) ;; I[-1/2](z) = sqrt(2/%pi/z)*cosh(z) (bessel-i-half-order rat-order arg)))) Dieter Kaiser ---------------------------------------------------------------------- Comment By: Dieter Kaiser (crategus) Date: 2008-07-04 23:10 Message: Logged In: YES user_id=2039760 Originator: NO If the gloabal flag BESSELEXPAND is FALSE we get for all Bessel functions with a Rational order and a zero as argument a noun form with the CVS code 1.58. (In previous versions a lot of the zero checks were not implemented.) But, if we set the global flag BESSELEXPAND to TRUE the Bessel function will be expanded and then evaluated for the argument zero. We get for all Bessel functions the described error. I think it is a question of convention and the desired behaviour of Maxima what we should do. The code to catch the zero argument is only implemented for CL numbers. The last time the code was extended to look more carefully for zero arguments didn't extend this to Maxima Rational, Bigfloat numbers or even expression which might simplify. In my opinion the reason to look only for CL numbers is, that only for such numbers we do the numerical evaluation. We check for a zero argument to prevent Lisp Errors in the numerical slatec routines. At this point we need to do the simplification for the Bessel function with zero argument. In the case of Rational numbers it depends on the flag BESSELEXPAND what we should do. In my opinion it is the best to check for the zero argument before we actually do the expansion. If we try to simplify all Bessel functions with zero argument in advance we get a lot of problems and questions, too. The reason is that the result of the simplification with an argument of zero depends on the sign of the order and is different for complex and real order of the Bessel function. Here an example what we can do in simp-bessel-i: ((and $besselexpand (setq rat-order (max-numeric-ratio-p order 2))) (cond ((= arg 0) ;; We don't do the expansion for a zero argument. (if (> rat-order 0) 0 '$infinity)) (t ;; When order is a fraction with a denominator of 2, we ;; can express the result in terms of elementary ;; functions. ;; ;; I[1/2](z) = sqrt(2/%pi/z)*sinh(z) ;; I[-1/2](z) = sqrt(2/%pi/z)*cosh(z) (bessel-i-half-order rat-order arg)))) We no longer get the divison by zero error. Here the results: (%i7) besselexpand:false; (%o7) false (%i8) bessel_i(1/2,0); (%o8) bessel_i(1/2,0) (%i9) bessel_i(-1/2,0); (%o9) bessel_i(-1/2,0) (%i10) besselexpand:true; (%o10) true (%i11) bessel_i(1/2,0); (%o11) 0 (%i12) bessel_i(-1/2,0); (%o12) infinity This check for a zero argument before we expand the Bessel function could be implemented for all other Bessel functions too. Dieter Kaiser ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 21:07 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 21:07 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 21:06 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Barton Willis (willisbl) Date: 2008-05-02 20:05 Message: Logged In: YES user_id=895922 Originator: NO See A&S 9.1.7: http://www.convertit.com/Go/Convertit/Reference/AMS55.ASP?Res=150&Page=360 For n > 0, bessel_j(n,0) = 0 (n needn't be integer). Maxima doesn't use this simplification for noninteger n. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&group_id=4933 |
From: SourceForge.net <no...@so...> - 2008-07-09 19:19:59
|
Bugs item #1954846, was opened at 2008-04-30 16:11 Message generated for change (Comment added) made by crategus You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&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: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: bessel_i(1/2,0) -> divide by zero error Initial Comment: bessel_i(1/2,0) produces an error. I think it should be 0 since bessel_i(1/2,x) is sqrt(2/%pi)*sinh(x)/sqrt(x). ---------------------------------------------------------------------- >Comment By: Dieter Kaiser (crategus) Date: 2008-07-09 21:19 Message: Logged In: YES user_id=2039760 Originator: NO For all Bessel functions a test has been added, which check for a zero argument before the expansion in terms of elementary functions. Example for Bessel I: besselexpand:false; bessel_i(1/2,0) -> bessel_i(1/2,0) besselexpand:true; bessel_i(1/2,0) -> 0 Close the bug. ---------------------------------------------------------------------- Comment By: Dieter Kaiser (crategus) Date: 2008-07-04 23:45 Message: Logged In: YES user_id=2039760 Originator: NO Sorry, but I have omitted to check arg to be an CL number. The correct code would be: ((and $besselexpand (setq rat-order (max-numeric-ratio-p order 2))) (cond ((and (numberp arg) (= arg 0)) (if (> rat-order 0) 0 '$infinity)) ( t ;; When order is a fraction with a denominator of 2, we ;; can express the result in terms of elementary ;; functions. ;; ;; I[1/2](z) = sqrt(2/%pi/z)*sinh(z) ;; I[-1/2](z) = sqrt(2/%pi/z)*cosh(z) (bessel-i-half-order rat-order arg)))) Dieter Kaiser ---------------------------------------------------------------------- Comment By: Dieter Kaiser (crategus) Date: 2008-07-04 23:10 Message: Logged In: YES user_id=2039760 Originator: NO If the gloabal flag BESSELEXPAND is FALSE we get for all Bessel functions with a Rational order and a zero as argument a noun form with the CVS code 1.58. (In previous versions a lot of the zero checks were not implemented.) But, if we set the global flag BESSELEXPAND to TRUE the Bessel function will be expanded and then evaluated for the argument zero. We get for all Bessel functions the described error. I think it is a question of convention and the desired behaviour of Maxima what we should do. The code to catch the zero argument is only implemented for CL numbers. The last time the code was extended to look more carefully for zero arguments didn't extend this to Maxima Rational, Bigfloat numbers or even expression which might simplify. In my opinion the reason to look only for CL numbers is, that only for such numbers we do the numerical evaluation. We check for a zero argument to prevent Lisp Errors in the numerical slatec routines. At this point we need to do the simplification for the Bessel function with zero argument. In the case of Rational numbers it depends on the flag BESSELEXPAND what we should do. In my opinion it is the best to check for the zero argument before we actually do the expansion. If we try to simplify all Bessel functions with zero argument in advance we get a lot of problems and questions, too. The reason is that the result of the simplification with an argument of zero depends on the sign of the order and is different for complex and real order of the Bessel function. Here an example what we can do in simp-bessel-i: ((and $besselexpand (setq rat-order (max-numeric-ratio-p order 2))) (cond ((= arg 0) ;; We don't do the expansion for a zero argument. (if (> rat-order 0) 0 '$infinity)) (t ;; When order is a fraction with a denominator of 2, we ;; can express the result in terms of elementary ;; functions. ;; ;; I[1/2](z) = sqrt(2/%pi/z)*sinh(z) ;; I[-1/2](z) = sqrt(2/%pi/z)*cosh(z) (bessel-i-half-order rat-order arg)))) We no longer get the divison by zero error. Here the results: (%i7) besselexpand:false; (%o7) false (%i8) bessel_i(1/2,0); (%o8) bessel_i(1/2,0) (%i9) bessel_i(-1/2,0); (%o9) bessel_i(-1/2,0) (%i10) besselexpand:true; (%o10) true (%i11) bessel_i(1/2,0); (%o11) 0 (%i12) bessel_i(-1/2,0); (%o12) infinity This check for a zero argument before we expand the Bessel function could be implemented for all other Bessel functions too. Dieter Kaiser ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 21:07 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 21:07 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2008-06-20 21:06 Message: Logged In: YES user_id=28849 Originator: YES Doesn't seem get a divide by zero error anymore in latest CVS. It returns the noun form. But bessel_i(v,x) ~ (1/z*x)^v/gamma(v+1) for small x, it should be 0. I think the bug is in simp-bessel-i in the lines: ((or (and (numberp order) (> order 0)) (integerp order)) 0) ((and (numberp order) (< order 0)) '$infinity) Perhaps numberp should be changed to mnump and (> order 0) should be (great order 0) ---------------------------------------------------------------------- Comment By: Barton Willis (willisbl) Date: 2008-05-02 20:05 Message: Logged In: YES user_id=895922 Originator: NO See A&S 9.1.7: http://www.convertit.com/Go/Convertit/Reference/AMS55.ASP?Res=150&Page=360 For n > 0, bessel_j(n,0) = 0 (n needn't be integer). Maxima doesn't use this simplification for noninteger n. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1954846&group_id=4933 |