## maxima-bugs

 [Maxima-bugs] [ maxima-Bugs-1376860 ] specint(gammaincomplete(v, a*t)*exp(-p*t), t) seems wrong From: SourceForge.net - 2006-08-15 02:48:40 ```Bugs item #1376860, was opened at 2005-12-08 20:53 Message generated for change (Settings changed) made by robert_dodier You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1376860&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 - Integration Group: None Status: Open Resolution: None Priority: 5 Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: specint(gammaincomplete(v,a*t)*exp(-p*t),t) seems wrong Initial Comment: (%i12) assume(v>0,p>0); (%o12) [v>0,p>0] (%i13) specint(gammaincomplete(v,a*t)*exp(-p*t),t); SIMP2F1-WILL-CONTINUE-IN (%o13) a^v*(p+a)^(-v-1)*gamma(v+1)*%f[2,1]([v+1,3/2],[2],p/(p+a)) Compare this with (%i14) specint(gammagreek(v,a*t)*exp(-p*t),t); (%o14) a^v*p^(-v-1)*gamma(v+1)/((a/p+1)^v*v) This matches formula 34, p 179 in Tables of Transforms. Considering gammaincomplete = gamma(n)-gammagreek, the expression for gammaincomplete seems wrong. It might still be right if the hypergeometric function simplifies, but maxima can't, and I can't think of any way to simplify it either. ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2006-02-12 12:15 Message: Logged In: YES user_id=28849 The fix for transforming %w causes this return the result in terms of the associated Legendre function Q. Somewhat better, but I do not know if this is equivalent. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1376860&group_id=4933 ```
 [Maxima-bugs] [ maxima-Bugs-1376860 ] specint(gammaincomplete(v, a*t)*exp(-p*t), t) seems wrong From: SourceForge.net - 2008-06-15 16:11:54 ```Bugs item #1376860, was opened at 2005-12-09 04:53 Message generated for change (Comment added) made by crategus You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1376860&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 - Integration Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: specint(gammaincomplete(v,a*t)*exp(-p*t),t) seems wrong Initial Comment: (%i12) assume(v>0,p>0); (%o12) [v>0,p>0] (%i13) specint(gammaincomplete(v,a*t)*exp(-p*t),t); SIMP2F1-WILL-CONTINUE-IN (%o13) a^v*(p+a)^(-v-1)*gamma(v+1)*%f[2,1]([v+1,3/2],[2],p/(p+a)) Compare this with (%i14) specint(gammagreek(v,a*t)*exp(-p*t),t); (%o14) a^v*p^(-v-1)*gamma(v+1)/((a/p+1)^v*v) This matches formula 34, p 179 in Tables of Transforms. Considering gammaincomplete = gamma(n)-gammagreek, the expression for gammaincomplete seems wrong. It might still be right if the hypergeometric function simplifies, but maxima can't, and I can't think of any way to simplify it either. ---------------------------------------------------------------------- Comment By: Crategus (crategus) Date: 2008-06-15 18:11 Message: Logged In: YES user_id=2039760 Originator: NO I would like to suggest the following algorithm to get an evaluated integral for the Gammaincomplete function: (defun gammaincomplete-to-hypergeometric (a x) (if (eq (asksign a) '\$zero) ;; The representation as Whittaker W function for a=0 (mul* (power x (div (sub a 1) 2)) (power '\$%e (div x -2)) (\$whittaker_w (div (sub a 1) 2) (div a 2) x)) ; sign of (div a 2) correct? ;; In all other cases the representation as Hypergeometric 1F1 function (sub (list '(%gamma) a) (mul (power x a) (inv a) (\$hypergeometric (list '(mlist) a) (list '(mlist) (add a 1)) (mul -1 x)))))) This routine has to be replaced for GAMMAINCOMPLETETW and called in the routine FRACTEST2. The idea is to use the Hypergeometric Function 1F1 when the parameter A is not equal to zero. That is equivalent to use gamma(a) - gammagreek(a,x) for the evaluation. Now we use the transformation to Whittaker W only for a=0, but for this case we get a nice result too. Here some examples with the above routine: (%i1) assume(s>0,a>0); (%o1) [s > 0,a > 0] We get the expected simple answer for gammaincomplete: (%i2) radcan(specint(%e^(-s*t)*gammaincomplete(a,t),t)); (%o2) (a*gamma(a)*(s+1)^a-gamma(a+1))/(a*s*(s+1)^a) If we take a=0 we get a nice and correct result too: (%i3) radcan(specint(%e^(-s*t)*gammaincomplete(0,t),t)); (%o3) log(s+1)/s The function %ei uses the code of gammaincomplete with a=0: (%i4) radcan(specint(%e^(-s*t)*%ei(t),t)); (%o4) -log(1-s)/s For special values like a=1/2 or a=1 both algorithm get simple expressions. The results calulated with 1F1 are identical to the results of the orginal code which use the transformation to Whittaker W: (%i5) radcan(specint(%e^(-s*t)*gammaincomplete(1/2,t),t)); (%o5) (sqrt(%pi)*sqrt(s+1)-sqrt(%pi))/(s*sqrt(s+1)) (%i6) radcan(specint(%e^(-s*t)*gammaincomplete(1,t),t)); (%o6) 1/(s+1) In the above code I have used the routines \$hypergeometric and \$whittaker_w. I have introduced this functions to make the code more expressive. It is easy to use the alternavtives %f and %w. Furthermore \$specint has to be extended to integrate %f. The testsuite has no problems with this alogorithm. Dieter Kaiser ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2006-02-12 20:15 Message: Logged In: YES user_id=28849 The fix for transforming %w causes this return the result in terms of the associated Legendre function Q. Somewhat better, but I do not know if this is equivalent. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1376860&group_id=4933 ```
 [Maxima-bugs] [ maxima-Bugs-1376860 ] specint(gammaincomplete(v, a*t)*exp(-p*t), t) seems wrong From: SourceForge.net - 2008-07-18 22:17:10 ```Bugs item #1376860, was opened at 2005-12-09 04:53 Message generated for change (Comment added) made by crategus You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1376860&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 - Integration Group: None Status: Open Resolution: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: specint(gammaincomplete(v,a*t)*exp(-p*t),t) seems wrong Initial Comment: (%i12) assume(v>0,p>0); (%o12) [v>0,p>0] (%i13) specint(gammaincomplete(v,a*t)*exp(-p*t),t); SIMP2F1-WILL-CONTINUE-IN (%o13) a^v*(p+a)^(-v-1)*gamma(v+1)*%f[2,1]([v+1,3/2],[2],p/(p+a)) Compare this with (%i14) specint(gammagreek(v,a*t)*exp(-p*t),t); (%o14) a^v*p^(-v-1)*gamma(v+1)/((a/p+1)^v*v) This matches formula 34, p 179 in Tables of Transforms. Considering gammaincomplete = gamma(n)-gammagreek, the expression for gammaincomplete seems wrong. It might still be right if the hypergeometric function simplifies, but maxima can't, and I can't think of any way to simplify it either. ---------------------------------------------------------------------- >Comment By: Dieter Kaiser (crategus) Date: 2008-07-19 00:17 Message: Logged In: YES user_id=2039760 Originator: NO I have suggested to use a representation as a Hypergeometric 1F1 for a<>0 for the Gammaincomplete function to get more simple results for the Laplace transform of the Gammaincomplete function. Because \$specint has no code to do the Laplace transform for the Hypergeometric function itself (\$specint knows only an intern representation for the Hypergeometric function) it is more easy to transform to the Gammagreek function. The Gammagreek function will be transformed to the Hypergeometric 1F1 function. The following changes produce the correct results shown in this thread: =================================================================== RCS file: /cvsroot/maxima/maxima/src/hypgeo.lisp,v retrieving revision 1.40 diff -u -r1.40 hypgeo.lisp --- hypgeo.lisp 4 Jul 2008 14:12:52 -0000 1.40 +++ hypgeo.lisp 18 Jul 2008 22:03:29 -0000 @@ -2301,7 +2302,7 @@ ((eq flg 'kbateman) (kbatemantw i1 a1)) ((eq flg 'gammaincomplete) - (gammaincompletetw a1 i1)) + (gammaincomplete-to-gammagreek a1 i1)) ((eq flg 'kti) (kti i1 a1)) ((eq flg 'erfc) @@ -2582,6 +2583,28 @@ (power '\$%e (div x -2)) (wwhit x (div (sub a 1) 2)(div a 2)))) +;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; +;;; Only for a=0 we use the general representation as a Whittaker W function. +;;; In all other cases we transform to a Gammagreek function. The Gammagreek +;;; function will be further transformed to a Hypergeometric 1F1 representation. +;;; With this change we get more simple and correct results for the Laplace +;;; transform of the Gammaincomplete function. +;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; + +(defun gammaincomplete-to-gammagreek (a x) + (if (eq (asksign a) '\$zero) + ;; The representation as a Whittaker W function for a=0 + (mul + (power x (div (sub a 1) 2)) + (power '\$%e (div x -2)) + (wwhit x (div (sub a 1) 2) (div a 2))) + ;; In all other cases the representation as a Gammagreek function + (sub + (list '(%gamma) a) + (list '(\$gammagreek) a x)))) + +;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; + (defun distrexecinit (fun) (cond ((and (consp fun) (consp (car fun)) =================================================================== The testsuite has no problems with this algorithm. With this changes this bug report could be closed. Should I commit this change of the code? Dieter Kaiser ---------------------------------------------------------------------- Comment By: Dieter Kaiser (crategus) Date: 2008-06-15 18:11 Message: Logged In: YES user_id=2039760 Originator: NO I would like to suggest the following algorithm to get an evaluated integral for the Gammaincomplete function: (defun gammaincomplete-to-hypergeometric (a x) (if (eq (asksign a) '\$zero) ;; The representation as Whittaker W function for a=0 (mul* (power x (div (sub a 1) 2)) (power '\$%e (div x -2)) (\$whittaker_w (div (sub a 1) 2) (div a 2) x)) ; sign of (div a 2) correct? ;; In all other cases the representation as Hypergeometric 1F1 function (sub (list '(%gamma) a) (mul (power x a) (inv a) (\$hypergeometric (list '(mlist) a) (list '(mlist) (add a 1)) (mul -1 x)))))) This routine has to be replaced for GAMMAINCOMPLETETW and called in the routine FRACTEST2. The idea is to use the Hypergeometric Function 1F1 when the parameter A is not equal to zero. That is equivalent to use gamma(a) - gammagreek(a,x) for the evaluation. Now we use the transformation to Whittaker W only for a=0, but for this case we get a nice result too. Here some examples with the above routine: (%i1) assume(s>0,a>0); (%o1) [s > 0,a > 0] We get the expected simple answer for gammaincomplete: (%i2) radcan(specint(%e^(-s*t)*gammaincomplete(a,t),t)); (%o2) (a*gamma(a)*(s+1)^a-gamma(a+1))/(a*s*(s+1)^a) If we take a=0 we get a nice and correct result too: (%i3) radcan(specint(%e^(-s*t)*gammaincomplete(0,t),t)); (%o3) log(s+1)/s The function %ei uses the code of gammaincomplete with a=0: (%i4) radcan(specint(%e^(-s*t)*%ei(t),t)); (%o4) -log(1-s)/s For special values like a=1/2 or a=1 both algorithm get simple expressions. The results calulated with 1F1 are identical to the results of the orginal code which use the transformation to Whittaker W: (%i5) radcan(specint(%e^(-s*t)*gammaincomplete(1/2,t),t)); (%o5) (sqrt(%pi)*sqrt(s+1)-sqrt(%pi))/(s*sqrt(s+1)) (%i6) radcan(specint(%e^(-s*t)*gammaincomplete(1,t),t)); (%o6) 1/(s+1) In the above code I have used the routines \$hypergeometric and \$whittaker_w. I have introduced this functions to make the code more expressive. It is easy to use the alternavtives %f and %w. Furthermore \$specint has to be extended to integrate %f. The testsuite has no problems with this alogorithm. Dieter Kaiser ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2006-02-12 20:15 Message: Logged In: YES user_id=28849 The fix for transforming %w causes this return the result in terms of the associated Legendre function Q. Somewhat better, but I do not know if this is equivalent. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1376860&group_id=4933 ```
 [Maxima-bugs] [ maxima-Bugs-1376860 ] specint(gammaincomplete(v, a*t)*exp(-p*t), t) seems wrong From: SourceForge.net - 2008-08-06 11:11:37 ```Bugs item #1376860, was opened at 2005-12-09 04:53 Message generated for change (Comment added) made by crategus You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1376860&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 - Integration Group: None >Status: Closed Resolution: None Priority: 5 Private: No Submitted By: Raymond Toy (rtoy) Assigned to: Nobody/Anonymous (nobody) Summary: specint(gammaincomplete(v,a*t)*exp(-p*t),t) seems wrong Initial Comment: (%i12) assume(v>0,p>0); (%o12) [v>0,p>0] (%i13) specint(gammaincomplete(v,a*t)*exp(-p*t),t); SIMP2F1-WILL-CONTINUE-IN (%o13) a^v*(p+a)^(-v-1)*gamma(v+1)*%f[2,1]([v+1,3/2],[2],p/(p+a)) Compare this with (%i14) specint(gammagreek(v,a*t)*exp(-p*t),t); (%o14) a^v*p^(-v-1)*gamma(v+1)/((a/p+1)^v*v) This matches formula 34, p 179 in Tables of Transforms. Considering gammaincomplete = gamma(n)-gammagreek, the expression for gammaincomplete seems wrong. It might still be right if the hypergeometric function simplifies, but maxima can't, and I can't think of any way to simplify it either. ---------------------------------------------------------------------- >Comment By: Dieter Kaiser (crategus) Date: 2008-08-06 13:11 Message: Logged In: YES user_id=2039760 Originator: NO The code is changed as suggested. Closing the bug report. Dieter Kaiser ---------------------------------------------------------------------- Comment By: Dieter Kaiser (crategus) Date: 2008-07-19 00:17 Message: Logged In: YES user_id=2039760 Originator: NO I have suggested to use a representation as a Hypergeometric 1F1 for a<>0 for the Gammaincomplete function to get more simple results for the Laplace transform of the Gammaincomplete function. Because \$specint has no code to do the Laplace transform for the Hypergeometric function itself (\$specint knows only an intern representation for the Hypergeometric function) it is more easy to transform to the Gammagreek function. The Gammagreek function will be transformed to the Hypergeometric 1F1 function. The following changes produce the correct results shown in this thread: =================================================================== RCS file: /cvsroot/maxima/maxima/src/hypgeo.lisp,v retrieving revision 1.40 diff -u -r1.40 hypgeo.lisp --- hypgeo.lisp 4 Jul 2008 14:12:52 -0000 1.40 +++ hypgeo.lisp 18 Jul 2008 22:03:29 -0000 @@ -2301,7 +2302,7 @@ ((eq flg 'kbateman) (kbatemantw i1 a1)) ((eq flg 'gammaincomplete) - (gammaincompletetw a1 i1)) + (gammaincomplete-to-gammagreek a1 i1)) ((eq flg 'kti) (kti i1 a1)) ((eq flg 'erfc) @@ -2582,6 +2583,28 @@ (power '\$%e (div x -2)) (wwhit x (div (sub a 1) 2)(div a 2)))) +;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; +;;; Only for a=0 we use the general representation as a Whittaker W function. +;;; In all other cases we transform to a Gammagreek function. The Gammagreek +;;; function will be further transformed to a Hypergeometric 1F1 representation. +;;; With this change we get more simple and correct results for the Laplace +;;; transform of the Gammaincomplete function. +;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; + +(defun gammaincomplete-to-gammagreek (a x) + (if (eq (asksign a) '\$zero) + ;; The representation as a Whittaker W function for a=0 + (mul + (power x (div (sub a 1) 2)) + (power '\$%e (div x -2)) + (wwhit x (div (sub a 1) 2) (div a 2))) + ;; In all other cases the representation as a Gammagreek function + (sub + (list '(%gamma) a) + (list '(\$gammagreek) a x)))) + +;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; + (defun distrexecinit (fun) (cond ((and (consp fun) (consp (car fun)) =================================================================== The testsuite has no problems with this algorithm. With this changes this bug report could be closed. Should I commit this change of the code? Dieter Kaiser ---------------------------------------------------------------------- Comment By: Dieter Kaiser (crategus) Date: 2008-06-15 18:11 Message: Logged In: YES user_id=2039760 Originator: NO I would like to suggest the following algorithm to get an evaluated integral for the Gammaincomplete function: (defun gammaincomplete-to-hypergeometric (a x) (if (eq (asksign a) '\$zero) ;; The representation as Whittaker W function for a=0 (mul* (power x (div (sub a 1) 2)) (power '\$%e (div x -2)) (\$whittaker_w (div (sub a 1) 2) (div a 2) x)) ; sign of (div a 2) correct? ;; In all other cases the representation as Hypergeometric 1F1 function (sub (list '(%gamma) a) (mul (power x a) (inv a) (\$hypergeometric (list '(mlist) a) (list '(mlist) (add a 1)) (mul -1 x)))))) This routine has to be replaced for GAMMAINCOMPLETETW and called in the routine FRACTEST2. The idea is to use the Hypergeometric Function 1F1 when the parameter A is not equal to zero. That is equivalent to use gamma(a) - gammagreek(a,x) for the evaluation. Now we use the transformation to Whittaker W only for a=0, but for this case we get a nice result too. Here some examples with the above routine: (%i1) assume(s>0,a>0); (%o1) [s > 0,a > 0] We get the expected simple answer for gammaincomplete: (%i2) radcan(specint(%e^(-s*t)*gammaincomplete(a,t),t)); (%o2) (a*gamma(a)*(s+1)^a-gamma(a+1))/(a*s*(s+1)^a) If we take a=0 we get a nice and correct result too: (%i3) radcan(specint(%e^(-s*t)*gammaincomplete(0,t),t)); (%o3) log(s+1)/s The function %ei uses the code of gammaincomplete with a=0: (%i4) radcan(specint(%e^(-s*t)*%ei(t),t)); (%o4) -log(1-s)/s For special values like a=1/2 or a=1 both algorithm get simple expressions. The results calulated with 1F1 are identical to the results of the orginal code which use the transformation to Whittaker W: (%i5) radcan(specint(%e^(-s*t)*gammaincomplete(1/2,t),t)); (%o5) (sqrt(%pi)*sqrt(s+1)-sqrt(%pi))/(s*sqrt(s+1)) (%i6) radcan(specint(%e^(-s*t)*gammaincomplete(1,t),t)); (%o6) 1/(s+1) In the above code I have used the routines \$hypergeometric and \$whittaker_w. I have introduced this functions to make the code more expressive. It is easy to use the alternavtives %f and %w. Furthermore \$specint has to be extended to integrate %f. The testsuite has no problems with this alogorithm. Dieter Kaiser ---------------------------------------------------------------------- Comment By: Raymond Toy (rtoy) Date: 2006-02-12 20:15 Message: Logged In: YES user_id=28849 The fix for transforming %w causes this return the result in terms of the associated Legendre function Q. Somewhat better, but I do not know if this is equivalent. ---------------------------------------------------------------------- You can respond by visiting: https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1376860&group_id=4933 ```