Fun with Fourier series, Clojure lazy sequences and Incanter

Jun 6, 2010

The fourier series expansion of a square wave is given by:

sin(x) + (1/3)sin(3x) + (1/5)sin(5x) + ...

Let's try to generate (and plot) this series in Clojure.

We need to install a program called `Incanter' to generate plots; Incanter is an R-like statistical computing environment for Clojure; you can download it from here. The downloaded zip archive is self-contained - you need not install anything else to run Incanter. Just unzip the package, change into the `incanter' directory and start the Clojure REPL by running a script called `script/repl'. You can test the installation by executing this code fragment:

user=> (use '(incanter core charts stats))
user=> (def x (range 0 6.28 0.01))
user=> (def y (map sin x))
user=> (view (xy-plot x y))

We have plotted the first term of the fourier expansion - a simple sin curve!

What does the following code fragment do?

(defn seq-mul [n s]
  (map (fn [x] (* n x)) s))

It takes a number `n' and a sequence `s' and returns another sequence each element of which is multiplied by `n'.

Here is a function which will add two sequences:

(defn seq-add [a b]
  (map + a b))

We need a list of odd numbers:


(def odds (iterate #(+ % 2) 1))

Now we can generate the sequence: (x, 3x, 5x ... )


(def odd-multiples
   (map seq-mul odds (repeat x)))

`odd-multiples' is a sequence each element of which is again a sequence. We use `odd-multiples' to generate the sequence (sin(x), sin(3x), sin(5x) ...):


(def odd-harmonics
   (map (fn [x] (map sin x)) odd-multiples))

Now we need to `scale' each value in the sequence - we first generate the scale values 1, 1/3, 1/5, 1/7 ...


(def scale-values
   (map (fn [x] (/ 1.0 x)) odds))

We are now ready to generate the sequence (sin(x), (1/3)sin(3x), (1/5)sin(5x) ...):


(def fourier-series
   (map seq-mul scale-values odd-harmonics))

`fourier-series' represents an infinite series. We can take say the first 10 terms and add them up by:


(reduce seq-add (take 10 fourier))

Let's wrap this up in a convenient function:


(defn plot-fourier [n]
  (let [y (reduce seq-add (take n fourier))]
     (view (xy-plot x y))))

`plot-fourier' plots the sum of the first `n' terms of the series. Have fun by trying out different values for `n'!

[Go to Code Clojure home] [Follow me on Twitter] [Go to pramode.net home]