Content from 2012-08
For a while now, I've wanted to make some of the great functionality of things like Octave (and it's parent Matlab) to Common Lisp. Mostly, I want things like the slice assignments (arrays with ranges in the indices) and matrix math. Also, it never hurts to make higher mathematical functions available in simple ways. So, I started writing Map, which I hope will do just that.
Lately, I've been working in singular optics, largely with Spatial Light Modulators, which are like very tiny LCD screens in front of a mirror. They're cool because they can transform the wavefront of a laser beam into any sort of beam you'd like! We use a fair bit of math to pump out the right images to put onto the SLM, and I'd been wanting to make a more general interface for that for a while. I'm not there yet, but I'd like to share the few things I've done so far.
In Octave, you can type
[1:10]and it will automatically be expanded to
[1 2 3 4 5 6 7 8 9 10]. This is great notation and can be used both for generating a vector or setting a range. In map, I've implemented those with
#r(1 10), which creates a 'range' class that just keeps track of the numbers, and #v(1 10), which will expand to
#(1 2 3 4 5 6 7 8 9 10), like the Octave code. This means you can do things like
(reduce #'+ #v(1 100))
To add up the numbers 1 to 100. That, of course, isn't hard, especially when Gauss came up with a much easier way, but you get the picture.
One further thing that's special, is that you can apply a function to the vector you're generating by passing another argument (you can also pass a step size by entering three numbers).
(reduce #' #v(1 .1 10 #'sin))
will give you the sum of the sine of the numbers 1, 1.1, ..., 9.9, 10! Moving on...
This is a much more meaty section, and was really fun for me, since I hadn't really implemented a ODE solver before (that I remember). Map now contains a Runge-Kutta Feldberg 45 solver, which is supposed to have an adaptive step size, and currently doesn't seem to, so I didn't do everything right, but it works to a certain accuracy.
Here's a demonstration of finding the angular position and velocity ($\theta, \omega$) of a damped, driven pendulum, taken from some of my old computational physics notes.
(defun damped-pendulum (ti yi)
(declare (type double-float ti)
(type (array double-float 2))
(vector (+ (* -1d0 (sin (elt yi 1))) (* -1d0 .5d0 (elt yi 0)))
(elt yi 0)))
(map::mute (multiple-value-setq (ts ys)
#r(0d0 1d-7 4d1)
#(0d0 1d0) 1d-10)))
Now, ts holds the time values, and ys is a 2x(length ts) array of position and velocity values! Once I implement array slices, it'll be easy to graph them, but for now, you need to manually extract one of the rows of ys to plot against ts. There's still a lot of work to be done. Anyway, the resulting graphs are:
Thanks for reading, I'll probably post more once I've got slices sorted out.
This blog covers vegan, ice cream, Work, Vegan, Uncategorized, Rust, Resistance, Research, Repair, Recipes, Physics, Nest, Music, Lisp, Linux, Links, Link, Life, Language, Information, Image, Gnome, Games, Gallery, Food, Fascism, Esperanto, Emacs, Electronics, Education, Dog, Did You Know That?, DRM, Community, Comics, Colgate, Code, Buddhism, Books, Aside, Anarchism, 3D Printing
View content from 2009-04, 2009-06, 2009-07, 2009-08, 2009-10, 2009-11, 2009-12, 2010-05, 2010-06, 2010-07, 2010-08, 2010-09, 2011-01, 2011-02, 2011-03, 2011-04, 2011-06, 2011-07, 2011-08, 2011-09, 2011-10, 2011-11, 2011-12, 2012-01, 2012-02, 2012-03, 2012-04, 2012-07, 2012-08, 2012-09, 2012-11, 2013-02, 2013-03, 2013-04, 2013-06, 2013-07, 2013-09, 2013-11, 2014-01, 2014-02, 2014-03, 2014-07, 2014-11, 2014-12, 2015-01, 2015-02, 2015-07, 2015-10, 2015-11, 2015-12, 2016-01, 2016-05, 2016-08, 2016-10, 2016-11, 2016-12, 2017-02, 2017-06, 2018-04, 2019-03, 2009-03