I have found the bug for the square of the Bessel I function. We transform to two Bessel J functions. In the transformation is missing a factor %i^-v. (In the following code I have already replaced the function 1fact and have inserted the powers of %i).

;; Laplace transform of square of Bessel I function

(cond ((setq l (onei^2 u))

(setq index1 (cdras 'v l)

arg1 (mul '$%i (cdras 'w l))

rest (mul (power '$%i (neg index1))

(power '$%i (neg index1)) ; the missing factor

(cdras 'u l)))

(return (lt1j^2 rest arg1 index1))))

Now we get:

(%i3) assume(s>0)$

We expand the Bessel function:

(%i4) specint(exp(-s*t)*bessel_i(1/2,t)^2,t),besselexpand:true;

(%o4) -log(1-4/s^2)/(2*%pi)

Now we use the hypergeometric code and get the same result:

(%i5) specint(exp(-s*t)*bessel_i(1/2,t)^2,t);

(%o5) -log(1-4/s^2)/(2*%pi)

Dieter Kaiser