A port of the scmutils package for symbolic mathematics from mit scheme to guile, in an effort to make scmutils and the examples from the Structure and Interpretation of Classical Mechanics available on a wider variety of architectures/operating systems.

After loading the package, call (set-current-module generic-environment) to use the generic operators for +, *, etc that can take functions, vectors, etc as operands. If you do not want to replace the standard scheme operators, so you can call generic operators with g:+, g:*, etc.

Requires: guile 2.0 or higher. Plotting is implemented through calls to gnuplot, which must be installed separately. Version 1.1 now supports guile 3.0.

Functionality not available in the port:

Example session:

% guile
guile> (load "load.scm")
guile> (set-current-module generic-environment)
guile> (define D derivative)
guile> (define f (literal-function 'f))
guile> (define f^2 (expt f 2))
guile> (pe ((D f^2) 't))
(* ((derivative f) t) 2 (f t))

guile> (pe ((D sin) 's))
(cos s)

guile> (pe ((partial-derivative (lambda (x y) (* x y)) 0) 's 't ))

guile> (pe ((partial-derivative (lambda (x y) (* x y)) 1) 's 't ))

guile> (pp (expression
  (let ((k (literal-number 'k)) (m (literal-number 'm)))
      (lambda (v)
       (let ((t (s:ref v 0))
             (q (s:ref v 1))
             (p (s:ref v 2)))
         (+ (/ (square p)
                   (* 2 m))
            (* 1/2 k (square q))
            (sin t)))))
     (up (literal-number 't)
        (up (literal-number 'x)
            (literal-number 'y))
        (down (literal-number 'px)
              (literal-number 'py)))))))
(down (cos t)
      (down (* 0.5 k (+ x x)) (* 0.5 k (+ y y)))
      (up (* (+ px px) (/ 1 (* 2 m)))
          (* (+ py py) (/ 1 (* 2 m)))))

; can't apply vector as function.
; must define q as a lambda expr instead of as a vector of literal functions.
guile> (define q (lambda (t) (up ((literal-function 'x) t)  
			  ((literal-function 'y) t)  
			  ((literal-function 'z) t))))

guile> (define ((L-free-particle mass) local)
  (let ((v (ref local 2)))
    (* 1/2 mass (square v))))

guile> (define ((Gamma q) t)
  (up t
      (q t)
      ((D q) t)))

guile> (define* ((Lagrange-equations Lagrangian #:optional dissipation-function) q)
  (let ((state-path (Gamma q)))
    (if (default-object? dissipation-function)
	(- (D (compose ((partial 2) Lagrangian) state-path))
	   (compose ((partial 1) Lagrangian) state-path))
	(- (D (compose ((partial 2) Lagrangian) state-path))
	   (compose ((partial 1) Lagrangian) state-path)
	   (- (compose ((partial 2) dissipation-function) state-path))))))

guile> (define (test-path t)
  (up (+ (* 'a t) 'a0)
      (+ (* 'b t) 'b0)
      (+ (* 'c t) 'c0)))

guile> (pp (expression ((Gamma q) 't)))
(up t
    (up (x t) (y t) (z t))
    (up ((derivative x) t)
        ((derivative y) t)
        ((derivative z) t)))

guile> (pp (expression (((Lagrange-equations (L-free-particle 'm)) test-path) 't)))
(down 0 0 0)

January 8, 2024