• May 06, 2021, 05:34:05 PM

Login with username, password and session length

Author Topic:  How to Get More Than 3 Series Approximation Terms?  (Read 1389 times)

0 Members and 1 Guest are viewing this topic.

Offline Byte11

  • Fractal Fanatic
  • ***
  • Posts: 39
    • Fractal To Desktop
How to Get More Than 3 Series Approximation Terms?
« on: December 04, 2017, 01:16:22 AM »
I've been lurking on the forums for a while and I've managed to implement Perturbation Theory, automatic glitch fixing, and series approximation. Now I want to know how to get more than 3 terms. Is there a general formula for a given term N? What is the formula for the fourth term? The fifth?

Also, do you guy calculate series approximation terms in arbitrary precision? I found that it doesn't make a difference.

Linkback: https://fractalforums.org/index.php?topic=571.0
Fractal To Desktop - Live Render Fractals To Your Desktop!
http://fractaltodesktop.com/

Offline claude

  • 3f
  • ******
  • Posts: 1849
    • mathr.co.uk
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #1 on: December 04, 2017, 05:04:55 PM »
You can get the formulas for series approximation by using the computer algebra system maxima (or perhaps other systems, but maxima is free) and then optimizing by hand afterwards (this is more necessary for higher power multibrot sets, for power 2 it doesn't make so much difference).

Code: [Select]
$ maxima
(%i1) sumexpand: true$
(%i2) cauchysum: true$
(%i3) F(p, c, z) := z^p + c$
(%i4) D(p, c, z) := F(p, c + C, z + Z) - F(p, C, Z)$
(%i5) A(p) := niceindices(collectterms(expand(D(p,c,sum(a[i]*c^i,i,1,inf))),c))$
(%i6) tex(A(2))$
\[ \sum_{j=2}^\infty c^j \sum_{i=1}^{j−1} a_i a_{j−i} + 2 Z \sum_{i=1}^\infty c^i a_i + c \]
Then the next ai is the coefficient of ci.

The optimsation takes the form of turning \( a_i a_j + a_j a_i \) into \( 2 a_i a_j \), which saves a complex multiplication and a complex addition, and costs one scalar-complex multiplication.  For higher powers the savings are much greater.

Native double and long double do seem to have enough precision for the series coefficients, the trouble is the range, hence using a separate exponent (as in KF's floatexp, mandelbrot-perturbator's edouble, etc) to allow really big values without overflow to infinity.

Offline Byte11

  • Fractal Fanatic
  • ***
  • Posts: 39
    • Fractal To Desktop
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #2 on: December 13, 2017, 04:47:37 AM »
Claude, I've been looking at your formula for a while and I can't make sense of it. I've managed to implement a datatype that can store larger numbers.

Thankfully, I found another formula (written by you) that I found to be a lot easier to work through. Here it is: http://mathr.co.uk/blog/2016-03-06_simpler_series_approximation.html

I went through it and solved for the next two terms and they match up with the formula that Botond found for E and D here: http://www.fractalforums.com/mandel-machine/mandel-machine/msg71342/#msg71342. I'm going to keep updating this thread as I get more formulas.

I have a few questions:
In layman's terms, how does the formula you provided relate to the one I linked?

Do you hand calculate each formula or do you calculate each in a loop? If it's a loop, does it automatically optimize the formula?

Is there another version of this for higher power Mandelbrot sets? I found your stack exchange question (https://math.stackexchange.com/questions/2443377/closed-form-for-these-power-series-coefficients) so I'm pretty sure (correct me if I'm wrong) that changing the p-value will change the order of the mandelbrot set formula that maxima outputs. Problem is, I still don't know how to use that formula. Did you make a similar article for higher power mandelbrot sets?

Finally, for determining the bailout point for the series, I was using something from newman on GitHub. Basically, he takes the squared magnitude of the B term (not just B but B*DeltaSub0Squared) times a tolerancy value and once that's greater than the C term (C * DeltaSub0Cubed) it breaks. I'm using 2^-64 as the tolerancy and I'm using the pixel in the top right for my DeltaSub0. It works really well, but I've only tested it up to 2^522 (Dinkydau's Flake). Is this a sound method for determining bailout? I don't know whether it'll stop working at some point or not. How do you determine bailout? I saw knighty's post, but there were a lot of changes and disagreement on whether the formula was effective or not, so I'm not sure what to use. What do you use?

Thanks!

Offline claude

  • 3f
  • ******
  • Posts: 1849
    • mathr.co.uk
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #3 on: December 13, 2017, 06:04:14 PM »
Quote
In layman's terms, how does the formula you provided relate to the one I linked?

The "next ai is the coefficient of ci" means you have to go through the formula and separate out all the different "coefficient * ci" terms, and add all the matching power coefficients together.  At the end it works out like this:
\[ a_1 \leftarrow 2 Z a_1 + 1 \\ a_j \leftarrow 2 Z a_j + \sum_{i=1}^{j-1} a_i a_{j-i} \quad\text{ \(j > 1\) } \]
which can be optimized to
\[ a_1 \leftarrow 2 Z a_1 + 1 \\
a_j \leftarrow 2 Z a_j + 2 \sum_{i=1}^{\frac{j}{2}-1} a_i a_{j-i} + a_\frac{j}{2}^2 \quad\text{ \(j > 1\) even} \\
a_j \leftarrow 2 Z a_j + 2 \sum_{i=1}^{\frac{j-1}{2}} a_i a_{j-i} \quad\text{ \(j > 1\) odd}  \]
The optimisation is to reduce the number of complex multiplications, and is more necessary for higher powers.  A smart compiler in unsafe mode might be able to do this optimisation, but I wouldn't bet on it.
I implement the formula by a loop, means I can write one code for any order, and it's probably just as fast as unrolled code (a compiler can be told to unroll loops, too).

Quote
Is there another version of this for higher power Mandelbrot sets?
You can derive one by using the maxima code I posted and changing the argument in the last line (corresponding to p) from 2 to 3 or whatever you like.  For power 3 there will be coefficents with 3 terms multiplied together, and optimiziing this is more tricky, because there are up to 6 different orders they can occur in, rather than just 2.  For power p there will be p terms multiplied in p! orders.  One strategy might be to sort the indices of the terms you are multiplying, but then you have to count how many times that combination can occur.  Eg I think a1 a3 a3 occurs less often than a1 a2 a4, afaict, though both have indices summing to 7.  This is similar to the odd/even split for power 2.

Quote
Is this a sound method for determining bailout?
I think the answer is "not sound in general, but good enough for many cases".  You may run into overskipping (sometimes mild distortion at the edges, sometimes fantasy fractals) or underskipping (takes longer than necessary).  For sound(er) methods see knighty's truncation error formulas using ball arithmetic.  I'm not 100% sure if they are perfect yet, but they're certainly more principled than a simple magnitude check.

Quote
What do you use?
In mandelbrot-perturbator I use a similar test to yours, namely that the magnitude of all the terms is decreasing.  It's prone to underskipping so I want to replace it with knighty's formulas.

Offline Byte11

  • Fractal Fanatic
  • ***
  • Posts: 39
    • Fractal To Desktop
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #4 on: December 13, 2017, 10:49:30 PM »
In mandelbrot-perturbator I use a similar test to yours, namely that the magnitude of all the terms is decreasing.  It's prone to underskipping so I want to replace it with knighty's formulas.

Do you just store the last term that was calculated and check if it's decreased since then or do you use the derivative of the terms? Are the b terms in the article you wrote and I linked to for that purpose? What do the b terms do?

Also, when you say underskipping. How much is it underskipping by? Anything less than 500 is alright with me, but if it's skipping more than that at depths <2^1200, I'll probably take the time to learn and implement knighty's method.

Offline claude

  • 3f
  • ******
  • Posts: 1849
    • mathr.co.uk
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #5 on: December 13, 2017, 11:21:32 PM »
Do you just store the last term that was calculated and check if it's decreased since then or do you use the derivative of the terms? Are the b terms in the article you wrote and I linked to for that purpose? What do the b terms do?
I meant decreasing across the index of the terms, as in \( a_{j+1} c^{j+1} < K a_j c^j \) for all j.  Sorry for the confusion.  The b coefficients are for the derivative, but knighty pointed out you can use the derivative of the series instead of the series of the derivative, so it ends up as \( \sum j a_j c^{j-1} \) instead of needing another set of b coefficients.

Quote
How much is it underskipping by?

Impossible for me to quantify fully until I've implemented knighty's methods as comparison.  In one test knighty's SuperMB skipped around 4000 more iterations (around 300k iterations altogether) than my mandelbrot-perturbator.

Offline knighty

  • Fractal Feline
  • **
  • Posts: 198
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #6 on: December 14, 2017, 03:09:11 PM »
If you are about to implement the ball-interval-arithmetic based method, you could consider to take into account, unlike me, the rounding errors. In my implementation, the method works fine most (say 99%) of the time but it sometimes over-skips when using more than 32 coefficients. For example in the "wiggly fungus spores", with 64 coefficients, I get a glitch similar to this one.

Maybe one could use the Arb library for comparison. It uses MPFR so it may be relatively slow. I believe that using a correct implementation of ball interval arithmetic would have a side effect: The more coefficients, the more calculations, therefore the more rounding, especially when in a location with a lot of iterations. A series approximation with lot of coefficients may then under-skip more than a series with few coefficients.

Offline claude

  • 3f
  • ******
  • Posts: 1849
    • mathr.co.uk
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #7 on: December 15, 2017, 06:29:34 PM »
I had a read through some of the arb documentation, but not sure how it's polynomial series stuff works - it seems to have a polynomial with complex ball coefficients, rather than a polynomial with exact complex coefficients and a single separate radius term? maybe I'll try it over xmas...

Offline Byte11

  • Fractal Fanatic
  • ***
  • Posts: 39
    • Fractal To Desktop
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #8 on: December 18, 2017, 05:53:58 AM »
If you are about to implement the ball-interval-arithmetic based method, you could consider to take into account, unlike me, the rounding errors. In my implementation, the method works fine most (say 99%) of the time but it sometimes over-skips when using more than 32 coefficients. For example in the "wiggly fungus spores", with 64 coefficients, I get a glitch similar to this one.

At first glance, those errors are very similar to what I get when I over-skip by only a few iterations. Taking off ~10 iterations would likely solve the problem and be a lot faster than using arbitrary precision. Once I get 32 coefficients working and I update my program to actually render that deep, I'll try it out.

Offline claude

  • 3f
  • ******
  • Posts: 1849
    • mathr.co.uk
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #9 on: December 19, 2017, 08:23:40 AM »
did some more research.  acb_t is a pair of real balls aka a rectangle, not a complex disc ball, which makes using Arb likely a non-starter...

Offline knighty

  • Fractal Feline
  • **
  • Posts: 198
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #10 on: December 19, 2017, 09:40:42 AM »
Too bad. :/

BTW! for the higher exponent MB it is possible to have a series approximation by letting the computer do the calculation. In c++, one only needs to implement a class for a polynomial with the +,-,* and squaring operations. Raising to a power > 2 can be done using the multiplication and squaring operations. Now for the MB iteration, one just do:

S(t) is the polynomial used for series approximation.
Snext(t) = Pow(Sprevious(t), n);
Snext[0] = zhigh_precision; //Replace the first coefficient of Snext with the reference iterate that is computed at full precision
Snext[1] += 1; //add 1 to the second coefficient

Offline quaz0r

  • Fractal Feline
  • **
  • Posts: 153
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #11 on: December 20, 2017, 07:55:05 PM »
if using c++, how about the boost interval stuff?

http://www.boost.org/doc/libs/1_66_0/libs/numeric/interval/doc/interval.htm

its a lot of reading about a lot of stuff im not familiar with  ;D
though one good thing is being template-based it sounds like you can use your own types with it

Offline claude

  • 3f
  • ******
  • Posts: 1849
    • mathr.co.uk
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #12 on: December 20, 2017, 08:21:03 PM »
Yes reading that page I think you could use boost::multiprecision MPFR with boost::interval, though I don't know if the boilerplate (Policy) has been written already.  However this is a real two-endpoint interval, instead of a (possibly complex) midpoint-radius interval, so complex intervals may suffer from sqrt(2) too much expansion (when mutliplying by sqrt(i) for example...)

Offline Byte11

  • Fractal Fanatic
  • ***
  • Posts: 39
    • Fractal To Desktop
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #13 on: December 21, 2017, 01:04:26 AM »
I've been testing this out on DinkyDau's flake and I've been getting some weird results.

With 3 coefficients, I skip 15,765 iterations.
With 4 coefficients, the check I described above breaks at 10, 296 iterations, but I found that the approximation is actually invalid at 15,765 iterations.
With 5 coefficients, I skip 23,649 iterations.
With 6 coefficients, the check breaks at 12,086 iterations, but I found the the approximation is actually invalid at 15,765 iterations.
With 7 coefficients, I skip 15,766 iterations (yes this one skipped exactly one more iteration).

I've found that other locations follow a similar pattern (the numbers themselves are different), but 5 terms always seems to skip the most iterations.

Is series approximation really this erratic? Why does the series always skip 15,765 iterations except with 5 terms?? Am I just doing something wrong? What's weird is that it isn't just giving a completely messed up approximation with extra terms (if I change the formulas it does). It's actually giving an approximation, it's just less accurate. Is that normal?

Offline claude

  • 3f
  • ******
  • Posts: 1849
    • mathr.co.uk
Re: How to Get More Than 3 Series Approximation Terms?
« Reply #14 on: December 21, 2017, 01:26:33 AM »
By:
Quote
the check I described above
I assume you mean:
Quote
Finally, for determining the bailout point for the series, I was using something from newman on GitHub. Basically, he takes the squared magnitude of the B term (not just B but B*DeltaSub0Squared) times a tolerancy value and once that's greater than the C term (C * DeltaSub0Cubed) it breaks
(and I presume you mean less rather than greater?)

I suppose you are taking the penultimate term as B and the last one as C, instead of the second or third terms.  But maybe you should take into account more of the terms, as Botond K?sa suggested here: http://www.fractalforums.com/index.php?topic=18482.msg74200#msg74200


xx
How to implement Series Approximation in GLSL?

Started by Mr Rebooted on Programming

2 Replies
278 Views
Last post November 23, 2020, 10:45:25 PM
by sjhalayka
xx
Series approximation for the Burning Ship

Started by claude on Fractal Mathematics And New Theories

0 Replies
315 Views
Last post September 27, 2018, 06:05:07 PM
by claude
xx
Keyframe interpolation, series approximation and periods.

Started by unassigned on Programming

1 Replies
197 Views
Last post November 09, 2020, 06:27:53 PM
by claude
clip
using Misiurewicz points for perturbation and/or series approximation

Started by claude on Fractal Mathematics And New Theories

19 Replies
1538 Views
Last post November 10, 2017, 01:46:31 AM
by gerrit
clip
bad analytic DE/slopes when series approximation is disabled

Started by claude on Kalles Fraktaler

4 Replies
176 Views
Last post June 23, 2020, 07:44:35 AM
by claude