[Maxima] Two problems with laplace

Dieter Kaiser drdieterkaiser at web.de
Mon Mar 23 16:31:04 CDT 2009


I am working again an the laplace transforms and I have got two problems:

1. The laplace routines write and kill facts in the user database

The parameter of the Laplace transform is assumed to be positive. The algorithm
of the laplace routines writes this into the user database:

We do the assumption s>0:

(%i21) assume(s>0);
(%o21) [s > 0]

Now we do a Laplace transform which involves an integration routine:

(%i22) laplace(t^2*f(t),t,s);
(%o22) 'diff('laplace(f(t),t,s),s,2)

After finishing we have lost our fact:

(%i23) facts();
(%o23) []

This causes problems, when we would like to simplify the result further using
the fact s>0 and can cause subtle bugs in User code.

This would be a more correct session, when we do the integrations in a new
context. The new context is generated with the macro with-new-context:

(%i3) assume(s>0);
(%o3) [s > 0]
(%i4) laplace(t*f(t),t,s);
(%o4) -'diff('laplace(f(t),t,s),s,1)
(%i5) facts();
(%o5) [s > 0]

The testsuite has no problems. The diff to the original code is shown below.


2. When we do a lot of integrals, we get the error "Too many contexts"

I have observed this problem since a longer time, because I often repeat a lot
of tests with integrals. The integration routines $integrate and $defint do
their calculations in a new context which is generated with the macro
with-new-context.

If we leave the algorithm with a Maxima Error, the context will not be killed
correctly. After about 50 integrations which causes a Maxima Error, we have to
restart Maxima, because no more contexts are available.

When we introduce an unwind-protect in the macro with-new-context the problem
vanishes.

Index: maxmac.lisp
===================================================================

 (defmacro with-new-context (sub-context &rest forms)
   `(let ((context (context , at sub-context)))
-     (prog1 , at forms
+     (unwind-protect  
+       (progn , at forms)            ; I think progn is more correct!?
        (context-unwinder))))

I have tested the macro with this small sample routine:

(defun newcontext ()
  (with-new-context (context)
    (progn
      (format t "~&NEWCONTEXT conindex= ~A~%" conindex)
      (merror "WE LEAVE NEWCONTEXT WITH AN ERROR"))))

Without unwind-protect every call increases the value of conindex. After the
change, the context is correctly killed and we get for each call a context with
a condindex of 1.

Remark: The two following bugs

  integrate function raises "too many contexts" error - ID: 1960200
  monstertrig: too many contexts - ID: 1866077

had given this errors too, because of endless recursion. This has been fixed.
But the problem of "Too many contexts." is still present.


So would do you think. Should we do the two suggested extensions to overcome the
problems I have described?

Dieter Kaiser


P.S.:
These are the changes to use with-new-context in the laplace routines:

Index: laplac.lisp
===================================================================
RCS file: /cvsroot/maxima/maxima/src/laplac.lisp,v
retrieving revision 1.13
diff -u -r1.13 laplac.lisp
--- laplac.lisp   21 Mar 2009 14:32:24 -0000   1.13
+++ laplac.lisp   23 Mar 2009 20:52:38 -0000
@@ -140,9 +140,10 @@
        ;; declaration before an integration is done. Therefore we declare
        ;; the parameter of the Laplace transform to be positive before 
        ;; we call $specint too.
-       (meval `(($assume) ,@(list (list '(mgreaterp) parm 0))))
-       (setq res ($specint (mul fun (power '$%e (mul -1 var parm))) var))
-       (meval `(($forget) ,@(list (list '(mgreaterp) parm 0))))
+       (with-new-context (context)
+         (progn
+           (meval `(($assume) ,@(list (list '(mgreaterp) parm 0))))
+           (setq res ($specint (mul fun (power '$%e (mul -1 var parm))) var))))
        (if (or (isinop res '%specint)  ; Both symobls are possible, that is
                (isinop res '$specint)) ; not consistent! Check it! 02/2009
            ;; $specint has not found a result.
@@ -336,8 +337,13 @@
                      ;; $defint should not throw a Maxima error,
                      ;; therefore we set the flags errcatch and $errormsg.
                      ;; errset catches the error and returns nil
-                     (let ((errcatch t) ($errormsg nil))
-                       (errset ($defint f x a '$inf))))))
+                     (with-new-context (context)
+                       (progn
+                         (meval `(($assume) ,@(list (list '(mgreaterp) parm
0))))
+                         (meval `(($assume) ,@(list (list '(mgreaterp) x 0))))
+                         (meval `(($assume) ,@(list (list '(mgreaterp) a 0))))
+                         (let ((errcatch t) ($errormsg nil))
+                           (errset ($defint f x a '$inf))))))))
     (if tryint
    (car tryint)
    (list '(%integrate simp) f x a '$inf))))
@@ -529,14 +535,16 @@
      (and ($unknown fun)(go skip))
      (setq mult (simptimes (list '(mtimes) (exponentiate
                    (list '(mtimes simp) -1 var parm)) fun) 1 nil))
-     (meval `(($assume) ,@(list (list '(mgreaterp) parm 0))))
-     (setq tryint
-           ;; $defint should not throw a Maxima error.
-           ;; therefore we set the flags errcatch and errormsg.
-           ;; errset catches an error and returns nil.
-           (let ((errcatch t) ($errormsg nil))
-             (errset ($defint mult var 0 '$inf))))
-     (meval `(($forget) ,@(list (list '(mgreaterp) parm 0))))
+     (with-new-context (context)
+       (progn
+         (meval `(($assume) ,@(list (list '(mgreaterp) parm 0))))
+         (setq tryint
+               ;; $defint should not throw a Maxima error.
+               ;; therefore we set the flags errcatch and errormsg.
+               ;; errset catches an error and returns nil.
+               (let ((errcatch t) ($errormsg nil))
+                 (errset ($defint mult var 0 '$inf))))))
+;     (meval `(($forget) ,@(list (list '(mgreaterp) parm 0))))
      (and tryint (not (eq (caaar tryint) '%integrate))  (return (car tryint)))
      skip (return (list '(%laplace simp) fun var parm))))




More information about the Maxima mailing list