Subject: Integrator loops endlessly for lambert_w(1/x)
From: Dieter Kaiser
Date: Fri, 26 Dec 2008 05:29:41 +0100
Today I have recognized that David Billinghurst has added a test to
rtest_lambert_w.mac which fails:
integrate(lambert_w(1/x),x)
The integrator loops endlessly.
I had a long search for the problem. The integrator tries to do an
integration-by-parts and generates the following integrand:
%e^(-lambert_w(1/x))/(x*(1+lambert_w(1/x)))
This is the expression we get when we type it in at the Maxima prompt
(The numbers enumerate the different parts of the expression):
((MTIMES SIMP)
1* ((MEXPT SIMP) $X -1)
2* ((MEXPT SIMP) $%E
((MTIMES SIMP) -1 ((%LAMBERT_W SIMP) ((MEXPT SIMP) $X -1))))
3* ((MEXPT SIMP)
((MPLUS SIMP) 1 ((%LAMBERT_W SIMP) ((MEXPT SIMP) $X -1))) -1))
Then I have looked into the integrator. In the integrator the expression
is a bit different ordered. The expression is generated in the routine
partial-integration and transfered to the integrator. This is the
expression the integrator gets and tries to solve:
(%i7) :lisp (setq $form
'((MTIMES SIMP)
3* ((MEXPT SIMP)
((MPLUS SIMP) 1 ((%LAMBERT_W SIMP) ((MEXPT SIMP) $X -1))) -1)
1* ((MEXPT SIMP) $X -1)
2* ((MEXPT SIMP) $%E
((MTIMES SIMP) -1 ((%LAMBERT_W SIMP) ((MEXPT SIMP) $X -1))))))
If the routine intform finds a mexpt expression with a complicated
exponent it tries the following code:
(let* (($%e_to_numlog t)
(nexp (resimplify exp))) ; resimplify the whole integrand
(cond ((alike1 exp nexp) nil) ; no change, return nil
(t
(when *debug-integrate*
(format t "~& in INTFORM: NEW EXP is ~A~%" nexp))
(intform (setq exp nexp)))))))) ; try again with new exp
This code is involved doing the integral which fails. Now what happens
is that we get for every resimplify a new expression. The condition
(alike1 exp nexp) is never true and the integrator loops endlessly.
I have tried it directly at the Maxima prompt. Here are the results:
This is a first resimplify for the expression from above. Look at the
numbers:
(%i7) :lisp (resimplify $form)
((MTIMES SIMP)
2* ((MEXPT SIMP) $%E
((MTIMES SIMP) -1 ((%LAMBERT_W SIMP) ((MEXPT SIMP) $X -1))))
3* ((MEXPT SIMP)
((MPLUS SIMP) 1 ((%LAMBERT_W SIMP) ((MEXPT SIMP) $X -1))) -1)
1* ((MEXPT SIMP) $X -1))
We resimplify two times:
(%i7) :lisp (resimplify (resimplify $form))
((MTIMES SIMP)
1* ((MEXPT SIMP) $X -1)
2* ((MEXPT SIMP) $%E
((MTIMES SIMP) -1 ((%LAMBERT_W SIMP) ((MEXPT SIMP) $X -1))))
3* ((MEXPT SIMP)
((MPLUS SIMP) 1 ((%LAMBERT_W SIMP) ((MEXPT SIMP) $X -1))) -1))
And again:
(%i7) :lisp (resimplify (resimplify (resimplify $form)))
((MTIMES SIMP)
3* ((MEXPT SIMP)
((MPLUS SIMP) 1 ((%LAMBERT_W SIMP) ((MEXPT SIMP) $X -1))) -1)
1* ((MEXPT SIMP) $X -1)
2* ((MEXPT SIMP) $%E
((MTIMES SIMP) -1 ((%LAMBERT_W SIMP) ((MEXPT SIMP) $X -1)))))
It seems to me that we have got by accident an expression which does not
simplify as expected and change the form every time we do a further
resimplify.
At the moment I have no idea how to continue.
Dieter Kaiser