The taylor series is a mechanism for approximating a function by summing polynomial terms.
First, some definitions:
f(x) | the function to approximate |
f^{n}(x) | the nth derivative of f (note that f^{0}(x) = f(x)) |
n! | n factorial, n! = n * (n - 1) * (n - 2) * ... * 2 |
The taylor series expansion of f(x) around a point (c) is then written:
n -> infinity +---- \ f^{n}(c) * (x - c) ^{n} f(x) = + ----------------- / n! +---- n = 0
As an example, consider the case where:
f(x) | => | exp(x) |
c | => | 0 |
Since the derivative of exp(x) is exp(x) we know that f^{n}(x) = exp(x) for all n.
Here are the first few terms of the series:
value of n | term | value of term |
---|---|---|
0 | exp(0) * (x - 0)^{0} / 0! | 1 |
1 | exp(0) * (x - 0)^{1} / 1! | x |
2 | exp(0) * (x - 0)^{2} / 2! | x^{2} / 2 |
3 | exp(0) * (x - 0)^{3} / 3! | x^{3} / 6 |
4 | exp(0) * (x - 0)^{4} / 4! | x^{4} / 24 |
We can write the expansion of exp(x) around c = 0 more concisely as follows:
n -> infinity +---- \ x^{n} f(x) = + -- / n! +---- n = 0
Now lets look at how accurate the approximation is as we add more terms. In this case, lets use x = 3.
exp(3) | n | sum | approximation | difference |
---|---|---|---|---|
20.08553692318766774092 | 0 | 1 | 1 | 19.08553692318766774092 |
20.08553692318766774092 | 1 | 1 + x | 4 | 16.08553692318766774092 |
20.08553692318766774092 | 2 | 1 + x + x^{2}/2 | 8.5 | 11.58553692318766774092 |
20.08553692318766774092 | 3 | 1 + x + x^{2}/2 + x^{3}/6 | 13 | 7.08553692318766774092 |
20.08553692318766774092 | 4 | 1 + x + x^{2}/2 + x^{3}/6 + x^{4}/24 | 16.375 | 3.71053692318766774092 |
As you can see, the approximation converges pretty rapidly on the actual value of exp(x).
Ok, so how do we apply this to generating a fractal? Firstly, all the statements about 'x' above apply equally well to 'z' where z is the complex number (x + iy).
In our case, assume 'current' is the point 'z' for which we wish to create an approximation of exp(z). Our algorithm (in pseudo code) looks like this:
set actualValue = exp(current) set maxIterations = 20 set tolerance = .00001 set count = 0 set zSum = [0, 0] set denominator = 1 while count is less than maxIterations and the magnitude of the difference between actualValue and zSum is greater than tolerance do zSum = zSum + current^{count} / denominator if count is greater than 0 denominator = denominator * (count + 1) endif count = count + 1 done
In this example, our expansion point (c) is zero, we can readily generalize our algorithm to handle other values of (c)
set actualValue = exp(current) set maxIterations = 20 set tolerance = .00001 set count = 0 set zSum = [0, 0] set denominator = 1 set c = [2, 3] set cValue = exp(c) while count is less than maxIterations and the magnitude of the difference between actualValue and zSum is greater than tolerance do zSum = zSum + cValue * (current - c)^{count} / denominator if count is greater than 0 denominator = denominator * (count + 1) endif count = count + 1 done
Finally, here is a full-on fractorama formula file demonstrating this technique
fractal { mapping { (-2.33672316384180778215, 3.41009999999999990905, 1.18621468926553697010, 8.00000000000000000000) => (150, 195) } formula { // The value we wish to approximate zvalue = exp(current); $xavg = 0; $yavg = 0; $mavg = 0; $xmax = 0; $ymax = 0; $mmax = 0; $xtot = 0; $ytot = 0; $mtot = 0; $maxIterations = 15; $tolerance = 1e-6; zSum = 0; $denominator = 1; // The point around which we'll create our expansion. We use // exp(current) here but we could just as well have used: // // point = [2, 3] // point = sin(current) // etc. // // There is no special significance to the fact that we're using // exp(current) and we're attempting to approximate exp. point = exp(current); // Loop until we've been through the loop '$maxIterations' times // or we deem that zSum is close enough to zvalue while($count < $maxIterations && ssq(zvalue - zSum) > $tolerance) { zSum += exp(point) * ((current - point) ^ $count) / $denominator; if($count > 0) { $denominator *= ($count + 1); } value = cot(zSum); $x = mag(value ^ .1); $y = mag(value ^ .2); $m = mag(value ^ .4); $xavg = ($x + $xavg * $count) / ($count + 1); $yavg = ($y + $yavg * $count) / ($count + 1); $mavg = ($m + $mavg * $count) / ($count + 1); $xdev = sqrt(($x - $xavg) * ($x - $xavg)); $ydev = sqrt(($y - $yavg) * ($y - $yavg)); $mdev = sqrt(($m - $mavg) * ($m - $mavg)); $xmax = max($xmax, $xdev); $ymax = max($ymax, $ydev); $mmax = max($mmax, $mdev); $xtot += $xdev; $ytot += $ydev; $mtot += $mdev; } $r = 0; $g = 0; $b = 0; if($count > 0) { $r = ($xtot/$count) / $xmax; $g = ($ytot/$count) / $ymax; $b = ($mtot/$count) / $mmax; $r = get_sin_color($r, 255, 1); $g = get_sin_color($g, 255, 1); $b = get_sin_color($b, 255, 1); } set_color($r, $g, $b); } }
And here is the resulting image: