Input

• function f(x,y)
  • math equation ( for symbolic computations)
  • program subroutine ( for numerical computations)
• matrix of floating-point values ( scalar field )^[1][2] ( for numerical computations)

tracing

Description by David Rand^[3]

The plotting algorithm can be summarized thus:

a) Given a matrix of floating point values which are the values of a function z = f(x,y) given at the nodes of a grid of x and y values (the grid values are assumed equally spaced, although the horizontal and vertical spacing may differ), the program determines the minimum and maximum values of z and then computes a number of contour values (in this implementation, 10 values) by linear or logarithmic interpolation between the extrema.

b) The program "walks" about the grid of points looking for any segment (i.e. a line joining two adjacent nodes in the grid) which must be crossed by one of the contours because some contour value lies between the values of z at the nodes.

c) Having found such a segment, it finds the intersection point of the contour and the segment by linear interpolation between the nodes. It also stores the information that the current contour value has been located on the current segment, so that this operation will not be repeated.

d) The program then attempts to locate a neighbouring segment having a similar property - that is, crossed by the same contour. If it finds one, it determines the intersection point as in step c) and then draws a straight line segment joining the previous intersection point with the current one. This step is repeated until no such neighbour can be found, taking care to exclude any segment which has already been dealt with.

e) Steps b), c) and d) are repeated until no segment can be found whose intersection with any contour value has not already been processed.

Lisp

Lisp code from Maxima CAS

Function ploteq

 ploteq (exp, ...options...)

plots equipotential curves for exp, which should be an expression depending on two variables. The curves are obtained by integrating the differential equation that define the orthogonal trajectories to the solutions of the autonomous system obtained from the gradient of the expression given. The plot can also show the integral curves for that gradient system (option fieldlines).

This program also requires Xmaxima, even if its run from a Maxima session in a console, since the plot will be created by the Tk scripts in Xmaxima. By default, the plot region will be empty until the user clicks in a point (or gives its coordinate with in the set-up menu or via the trajectory_at option).

Most options accepted by plotdf can also be used for ploteq and the plot interface is the same that was described in plotdf.

Example:

  load("plotdf.lisp")$
  V: 900/((x+1)^2+y^2)^(1/2)-900/((x-1)^2+y^2)^(1/2)$
  ploteq(V,[x,-2,2],[y,-2,2],[fieldlines,"blue"])$

Clicking on a point will plot the equipotential curve that passes by that point (in red) and the orthogonal trajectory (in blue).

source code files:

• share/dynamics/dynamics.mac
• share/diffequations/drawdf.mac
• share/dynamics/plotdf.lisp

;; plotdf.mac - Adds a function plotdf() to Maxima, which draws a Direction
;;              Field for an ordinary 1st order differential equation,
;;              or for a system of two autonomous 1st order equations.
;;   
;; Copyright (C) 2004, 2008, 2011 Jaime E. Villate <villate@fe.up.pt>
;;


  
;; plot equipotential curves for a scalar field f(x,y)
(defun $ploteq (fun &rest options)
  
  (let (file cmd mfun (opts " ") (s1 '$x) (s2 '$y))
    (setf mfun `((mtimes) -1 ,fun))
    ;; if the variables are not x and y, their names must be given right
    ;; after the expression for the ode's
    (when
      (and (listp (car options)) (= (length (car options)) 3)
            (or (symbolp (cadar options)) ($subvarp (cadar options)))
            (or (symbolp (caddar options)) ($subvarp (caddar options))))
      (setq s1 (cadar options))
      (setq s2 (caddar options))
      (setq options (cdr options)))
    (defun subxy (expr)
      (if (listp expr)
          (mapcar #'subxy expr)
          (cond ((eq expr s1) '$x) ((eq expr s2) '$y) (t expr))))
    (setf mfun (mapcar #'subxy mfun))
;; the next two lines should take into account parameters given in the options
;;    (if (delete '$y (delete '$x (rest (mfuncall '$listofvars ode))))
;;        (merror "The equation(s) can depend only on 2 variable which must be specified!"))
    (setq cmd (concatenate 'string " -dxdt \""
			   (expr_to_str (mfuncall '$diff mfun '$x))
			   "\" -dydt \""
			   (expr_to_str (mfuncall '$diff mfun '$y)) 
			   "\" "))
    
    ;; parse options and copy them to string opts
    (cond (options
           (dolist (v options) 
             (setq opts (concatenate 'string opts " "
                                  (plotdf-option-to-tcl v s1 s2))))))

    (unless (search "-vectors " opts)
      (setq opts (concatenate 'string opts " -vectors {}")))
    (unless (search "-fieldlines " opts)
      (setq opts (concatenate 'string opts " -fieldlines {}")))
    (unless (search "-curves " opts)
      (setq opts (concatenate 'string opts " -curves {red}")))
    (unless (search "-windowtitle " opts)
      (setq opts (concatenate 'string opts " -windowtitle {Ploteq}")))
    (unless (search "-xaxislabel " opts)
      (setq opts (concatenate 'string opts " -xaxislabel " (ensure-string s1))))
    (unless (search "-yaxislabel " opts)
      (setq opts (concatenate 'string opts " -yaxislabel " (ensure-string s2))))
							      
    (setq file (plot-temp-file (format nil "maxout~d.xmaxima" (getpid))))
    (show-open-plot
     (with-output-to-string
         (st)
       (cond ($show_openplot (format st "plotdf ~a ~a~%" cmd opts))
             (t (format st "{plotdf ~a ~a}" cmd opts))))
     file)
    file))