rtoy
commented
5 months ago Imported from SourceForge on 2024-07-04 21:18:21 Created by crategus on 2009-09-27 20:37:40 Original: https://sourceforge.net/p/maxima/bugs/1777/#73b8
The error is in the routine f16p217test. The Laplace transform of the Whittaker W function is only valid for an argument a*t and not a*t^v. There is no check for v=1.
This is a possible correction of the involved part of the code in the routine f16p217test:
;; Ok, we satisfy the conditions. Now extract the arg. ;; The transformation is only valid for an argument a*t. We have ;; to specialize the pattern to make sure that we satisfy the condition. (let ((l (m2-a*t a) ; more special pattern for the argument ; (m2 a ; '((mplus) ; ((coeffpt) (f hasvar) (a freevar)) ; ((coeffpp) (c zerp))) ; nil) ))
If we correct this part of code the Laplace transform for the Parabolic Cylinder D function for an argument sqrt(t) will work, but will fail for an argument t as expected.
Furthermore, we have to change the algorithm of whittest. When f16p217test fails a transformation to a Whittaker M function is done. But the algorithm for the Laplace transform of the Whittaker W function is already the general Laplace transform of the corresponding representation in terms of Whittaker M. The problem is, that we get errors when trying this transformation. So it is the best to return a noun form at this point.
This could be the change to the routine whittest:
(defun whittest (r a i1 i2) (cond ((f16p217test r a i1 i2)) (t ; The formula used in f16p217test already is already the general Laplace ; transformation. The transformation to a Whittaker M representation ; does not give any new information, but will fail with errors like undefined gamma ; functions or divisions by zero. Therefore, we return a noun form at this point. ; ;; Convert to M function and try again. ; (distrexecinit ($expand (mul (init r) ; (wtm a i1 i2)))) (setq *hyp-return-noun-flag* 'whittest-failed) )))
With these changes the Laplace transform of the Parabolic Cylinder D function will work as expected. There are no problems with the testsuite.
Dieter Kaiser
Imported from SourceForge on 2024-07-04 21:18:20 Created by crategus on 2009-09-26 23:30:08 Original: https://sourceforge.net/p/maxima/bugs/1777
A look at the implemented algorithm for the Laplace transform of a hypergeometric function shows that the integral for the Parabolic Cylinder D function should work only for an argument t^k with an exponent k < 1, e.g with an argument sqrt(t) (I have already introduced a new symbol in my sandbox for the Parabolic Cylinder D function):
(%i2) factor(ratsimp(specint(exp(-s*t)*parabolic_cylinder_d(0,sqrt(t)),t))); (%o2) 4/(4*s+1)
This is a correct result.
We can verify the result using a specific expansion:
parabolic_cylinder_d(0,z) = exp(-z^2/4)
parabolic_cylinder_d(1,z) = z*exp(-z^2/4) parabolic_cylinder_d(2,z) = (z^2-1)*exp(-z^2/4)
Now some more results with the first expansion for v=0:
(%i3) d0(z):=exp(-z^2/4)$
We can verify the result for the Parabolic Cylinder D function with an argument sqrt(t):
(%i6) factor(ratsimp(specint(exp(-s*t)*d0(sqrt(t)),t))); (%o6) 4/(4*s+1)
Now we do the integration with an argument t using the expansion. We get a result which contains the Error function:
(%i7) factor(ratsimp(specint(exp(-s*t)*d0(t),t))); (%o7) -sqrt(%pi)*%e^s^2*(erf(s)-1)
The following integral should not work, but it does and the answer is wrong:
(%i8) factor(ratsimp(specint(exp(-s*t)*parabolic_cylinder_d(0,t),t))); (%o8) 2^(3/2)*gamma(3/4)/(4*s+1)^(3/4)
The routine f16p217test does not fail as expected and produces the wrong result.
Remarks:
1. $specint uses the Parabolic Cylinder D function and its Laplace transform for the Hermite and the Erfc function. For the Hermite function the algorithm does not work as expected. For the erfc function it is more simple to integrate (1-erf(z)).
2. The Parabolic Cylinder D function is transformed in terms of the Whittaker W function. In a second step the Whittaker W function is transformed to a hypergeometric representation. This representation is integrated. But we already have a routine simpdtf which does the transformation to a hypergeometric representation in one step.
Dieter Kaiser