Title: Help about implementing the perturbation theory
Post by: CFJH on July 02, 2013, 08:35:52 PM
Hello,
I'm reading very interested the topic about sft with the significant speedups by using the perturbation theory.
http://www.fractalforums.com/announcements-and-news/superfractalthing-arbitrary-precision-mandelbrot-set-rendering-in-java/ (http://www.fractalforums.com/announcements-and-news/superfractalthing-arbitrary-precision-mandelbrot-set-rendering-in-java/)
Now I want to implement this in my own plain fractal generator written in C. But I havn't really understand how the perturbation theory is working.
Reviewing the java-code of sft is difficult to me because I'm not a java developer and some used structures are unknown to me.
Can someone explain me, how is it working or give a abstract (programming language indepent) code fragment ?
Jürgen
I'm reading very interested the topic about sft with the significant speedups by using the perturbation theory.
http://www.fractalforums.com/announcements-and-news/superfractalthing-arbitrary-precision-mandelbrot-set-rendering-in-java/ (http://www.fractalforums.com/announcements-and-news/superfractalthing-arbitrary-precision-mandelbrot-set-rendering-in-java/)
Now I want to implement this in my own plain fractal generator written in C. But I havn't really understand how the perturbation theory is working.
Reviewing the java-code of sft is difficult to me because I'm not a java developer and some used structures are unknown to me.
Can someone explain me, how is it working or give a abstract (programming language indepent) code fragment ?
Jürgen
Title: Re: Help about implementing the perturbation theory
Post by: Kalles Fraktaler on July 12, 2013, 01:16:49 PM
I recommend that you take the time and read these two threads through and the linked pdf:
http://www.fractalforums.com/announcements-and-news/superfractalthing-arbitrary-precision-mandelbrot-set-rendering-in-java/
http://www.fractalforums.com/announcements-and-news/kalles-fraktaler-2/
Even if the Java code might be difficult to digest for a C/C++ programmer all info needed is there.
http://www.fractalforums.com/announcements-and-news/superfractalthing-arbitrary-precision-mandelbrot-set-rendering-in-java/
http://www.fractalforums.com/announcements-and-news/kalles-fraktaler-2/
Even if the Java code might be difficult to digest for a C/C++ programmer all info needed is there.
Title: Re: Help about implementing the perturbation theory
Post by: Gluecker on July 13, 2013, 01:07:28 AM
the most significant speed boost will come from equation 1 (see sft_maths.pdf from the SuperFractalThing thread). Since the java code provided has additional tweaks, its not really easy to fully understand.
I will try to explain what I got so far:
- calculate one point (the reference X) near or inside your field of view at full precision and store the real and imaginary part of every calculation at lower precision. This point should be inside the M-set, or at least very close to it.
- compute the initial difference
between the point Y you want to render and the reference (using a lower precision)
- iterate over eq1 using and the stored values for X. Its important to start with the correct index for X. (Eg. I stored the starting coordinates of X as X[0], my first calculation is then X[1], so my iteration counter is initially set to 1)
- as usual the iteration ends when the iteration counter hits the iteration limit. You also want to check if the distance
gets bigger than a certain limit, since the assumption that "all the numbers are 'small'" was made. ^2+Im(\bigtriangleup_n)^2<1E-5)
works fine for me, but I am not done testing yet.
- when gets too big, we can calculate back to the current value of 
and go on with the standart mandelbrot calculation. This also seems to work with lower precision, but again, not sure about that.
hope that helps.
cheers
I will try to explain what I got so far:
- calculate one point (the reference X) near or inside your field of view at full precision and store the real and imaginary part of every calculation at lower precision. This point should be inside the M-set, or at least very close to it.
- compute the initial difference
- iterate over eq1 using
- as usual the iteration ends when the iteration counter hits the iteration limit. You also want to check if the distance
- when
hope that helps.
cheers
Title: Re: Help about implementing the perturbation theory
Post by: CFJH on August 07, 2013, 11:32:00 PM
Hallo
in the past days I had a litte time to implememt something about the pertubation theory.
Now I've an (nearly) working function which is currently using the equation (1) only.
But my problem is the bailout-value, where the interartion exits. The standard für a mandelbrot is 4 but I've to use a value 64 to get an unbroken image in all cases.
the attached Images shows
test1s : calculated with the classical method
test2s : image with bailout of 64
test3s : broken image with standard bailout of 4
Why is the image broken with standard bailout value ?
I can't find where is the problem in my code. Can somebody give me a hint ?
in the past days I had a litte time to implememt something about the pertubation theory.
Now I've an (nearly) working function which is currently using the equation (1) only.
But my problem is the bailout-value, where the interartion exits. The standard für a mandelbrot is 4 but I've to use a value 64 to get an unbroken image in all cases.
the attached Images shows
test1s : calculated with the classical method
test2s : image with bailout of 64
test3s : broken image with standard bailout of 4
Code:
int tMandel::CalcPointPertub(int ix, int iy) {
double D0r,D0i; // Delta 0
double DNr,DNi; // Delta n
double Xr,Xi; // aktuelles Xn aus dem Array
double Yr,Yi; // berechneter Yn-Wert
double Ti,Tr; // temporär
unsigned long pIter;
int Inside=1;
mpf_set_d(a,ix-0.5*SizeX);
mpf_mul(a,dx,a);
D0r=mpf_get_d(a);
mpf_set_d(a,iy-0.5*SizeY);
mpf_mul(a,dy,a);
D0i=mpf_get_d(a);
// Reference point in the center
// ∆n+1 = 2Xn ∆n + ∆n^2 + ∆0
DNr=D0r;
DNi=D0i;
for (pIter=0;pIter<MaxIter;pIter++) {
Xr=arx[pIter]; // reference values real part
Xi=ary[pIter]; // reference values imaginary part
// Calculation 2Xn ∆n
Tr = 2.0*(Xr*DNr - Xi*DNi);
Ti = 2.0*(Xr*DNi + Xi*DNr);
// + ∆n^2
Tr = Tr + (DNr*DNr - DNi*DNi);
Ti = Ti + 2.0*DNr*DNi;
// + ∆0
Tr += D0r;
Ti += D0i;
DNr = Tr;
DNi = Ti;
// ∆n = Yn − Xn -> Yn = ∆n - Xn
Yr = DNr - Xr;
Yi = DNi - Xi;
//Abort
if ((Yr*Yr+Yi*Yi)>64.0) { // <----- problem !!
Inside=0;
break;
}
if (Aborted) {
Calculation=FALSE;
return -1;
}
}
if (Inside)
return 0;
else
return pIter;
}
Why is the image broken with standard bailout value ?
I can't find where is the problem in my code. Can somebody give me a hint ?
Title: Re: Help about implementing the perturbation theory
Post by: knighty on August 08, 2013, 03:36:28 AM
// ∆n = Yn − Xn -> Yn = Xn + ∆n
Title: Re: Help about implementing the perturbation theory
Post by: Gluecker on August 10, 2013, 01:41:33 AM
yea, Y = X+delta.
other than that, higher bailouts dont hurt since the assumption "|Z|> 2 escapes" also implies |Z|> 8 does.
other than that, higher bailouts dont hurt since the assumption "|Z|> 2 escapes" also implies |Z|> 8 does.
Title: Re: Help about implementing the perturbation theory
Post by: knighty on August 10, 2013, 01:56:25 PM
Any zooms? :D
Title: Re: Help about implementing the perturbation theory
Post by: CFJH on August 10, 2013, 06:38:22 PM
Thanks for the hint,
but this solves the problem not completly. Now I can decrease the bailout value to 16 (instrad of 64) but the rendered image still looks different against an image rendered with standard method. A lower value as 16 produces a clipped image.
Is this normal or have I another bug ?
Jürgen
but this solves the problem not completly. Now I can decrease the bailout value to 16 (instrad of 64) but the rendered image still looks different against an image rendered with standard method. A lower value as 16 produces a clipped image.
Is this normal or have I another bug ?
Jürgen
Title: Re: Help about implementing the perturbation theory
Post by: knighty on August 10, 2013, 06:58:48 PM
Don't know. If it's approximately the same as the 3rd picture then I would say there still is an error in Yn calculation.
Title: Re: Help about implementing the perturbation theory
Post by: Kalles Fraktaler on August 11, 2013, 01:08:52 PM
Thanks for the hint,
but this solves the problem not completly. Now I can decrease the bailout value to 16 (instrad of 64) but the rendered image still looks different against an image rendered with standard method. A lower value as 16 produces a clipped image.
Is this normal or have I another bug ?
Jürgen
No it is NOT normal and yes you have a bug...but this solves the problem not completly. Now I can decrease the bailout value to 16 (instrad of 64) but the rendered image still looks different against an image rendered with standard method. A lower value as 16 produces a clipped image.
Is this normal or have I another bug ?
Jürgen
Sorry but I cannot find any error in you code either.
Here is my bailout validation, which looks very much like yours:
Code:
Di=D0i;
Dr=D0r;
for(antal=0;antal<m_nMaxIter;antal++){
Dnr = (m_db_dxr[antal]*Dr - m_db_dxi[antal]*Di)*2 + Dr*Dr - Di*Di + D0r;
Dni = (m_db_dxr[antal]*Di + m_db_dxi[antal]*Dr + Dr*Di)*2 + D0i;
yr=m_db_dxr[antal] + Dr;
yi=m_db_dxi[antal] + Di;
Di = Dni;
Dr = Dnr;
if(yr*yr+yi*yi > 4)
break;
}
Are you trying to do some coloring with the value of Yn? Because I have noticed that the value of Yn is unstable. It is accurate for calculating the bailout iteration, but it is not usable to do further coloring with the value since it get too far away for some locations.
Title: Re: Help about implementing the perturbation theory
Post by: CFJH on August 12, 2013, 11:39:58 PM
Hello,
now I could find my problem in the calculation :)
My mistake was, that I compute the Yn values one step too late
correct:
wrong:
Thanks to all for the help!
Jürgen
now I could find my problem in the calculation :)
My mistake was, that I compute the Yn values one step too late
correct:
Code:
Tr = 2.0*Xr*DNr - 2.0*Xi*DNi + DNr*DNr - DNi*DNi + D0r;
Ti = 2.0*Xr*DNi + 2.0*Xi*DNr + 2.0*DNr*DNi + D0i;
Yr=DNr+Xr;
Yi=DNi+Xi;
DNr = Tr;
DNi = Ti;
wrong:
Code:
Tr = 2.0*Xr*DNr - 2.0*Xi*DNi + DNr*DNr - DNi*DNi + D0r;
Ti = 2.0*Xr*DNi + 2.0*Xi*DNr + 2.0*DNr*DNi + D0i;
DNr = Tr;
DNi = Ti;
Yr=DNr+Xr;
Yi=DNi+Xi;
Thanks to all for the help!
Jürgen
Title: Re: Help about implementing the perturbation theory
Post by: laser blaster on August 15, 2013, 06:59:49 PM
Wow, this perturbation theory looks brilliant! I'm still struggling to wrap my head around the math, but if I understand it correctly, it should enable us to only have to iterate one point, store it's values at each iteration, and then discretely solve for the values of all other points in that area, without having to iterate them. That's huge! Even better, it seems like it could be applied to gpu rendering. I might be misinterpreting it, though.
Title: Re: Help about implementing the perturbation theory
Post by: cKleinhuis on August 15, 2013, 07:13:18 PM
You still need to iterate, but single precision is sufficient then
Title: Re: Help about implementing the perturbation theory
Post by: Kalles Fraktaler on August 15, 2013, 08:26:01 PM
You still need to iterate, but single precision is sufficient then
No single precision is not sufficient, at least I was not able to succeed with thatTitle: Re: Help about implementing the perturbation theory
Post by: Gluecker on August 17, 2013, 01:28:16 AM
At a certain magnification delta gets too small (if the reference is near the screen center) to be calculated with single precision (minValue for float: 1.18 × 10^−38).
One could try to find good references further away, but thats not trivial and I guess there are also other limits...
One could try to find good references further away, but thats not trivial and I guess there are also other limits...
Title: Re: Help about implementing the perturbation theory
Post by: laser blaster on August 18, 2013, 06:11:28 PM
You still need to iterate, but single precision is sufficient then
Hmm... well I've been re-reading the pdf, and it seems to indicate that you can skip iterating for at least some of the points.Now we can apply equations (3), (4) and (5) iteratively to calculate the coefficients for equation (2). Equation (2) can then be used to calculate the value for the nth iteration for all the points around X0. The approximation should be good as long as the (symbol)3 term has a magnitude signi ficantly smaller then the (symbol)2 term.
It seems like the coefficients for equation 2 only need to be calculated iteratively for one pixel. Then you can use equation 2 to estimate the position of all the pixels in the surrounding areas, but only when the symbol(3) term is significantly smaller than the (symbol)2 term.
Am I misunderstanding this?
Title: Re: Help about implementing the perturbation theory
Post by: Gluecker on August 21, 2013, 03:45:08 AM
You can "skip" the first x iterations, where x is the minimum iteration count for the image to render (minus some tolerance maybe).
alternatively, using the Mariani/Silver algorithm (as described here (http://mrob.com/pub/muency/marianisilveralgorithm.html)) could be a efficient way to find x for subsections of the image.
broadly speaking: 2) can be used to approximate delta for any iteration (n), but (n) obviously has to be less than the number of iterations where the point is actually escaping to get the correct image, so you need to guess a good (n).
alternatively, using the Mariani/Silver algorithm (as described here (http://mrob.com/pub/muency/marianisilveralgorithm.html)) could be a efficient way to find x for subsections of the image.
broadly speaking: 2) can be used to approximate delta for any iteration (n), but (n) obviously has to be less than the number of iterations where the point is actually escaping to get the correct image, so you need to guess a good (n).
Title: Re: Help about implementing the perturbation theory
Post by: therror on April 27, 2015, 12:57:33 PM
Hello! I am trying to implement perturbation theory. I have written a code similar to example by Kalles Fraktaler 08.11.2013.
Reference point is calculated as double, other as float.
Picture is almost the same as usual set. But speed decreases about 2 times... IMHO, perturbation theory is ineffective at such magnification (4.93E-8). Isn't it? If I'm wrong, why speed decreases?
Code:
dnr=d0r=spix[l].p-Pmin;
dni=d0i=spix[l].q-Qmin;
yp=yq=n=0;
while (yp*yp+yq*yq<4&&n<nmax)
{
dtr=2*(xref[n].p*dnr-xref[n].q*dni)+dnr*dnr-dni*dni+d0r;
dti=2*(xref[n].p*dni+xref[n].q*dnr+dnr*dni)+d0i;
yp=xref[n].p+dnr;
yq=xref[n++].q+dni;
dnr=dtr;
dni=dti;
}
Picture is almost the same as usual set. But speed decreases about 2 times... IMHO, perturbation theory is ineffective at such magnification (4.93E-8). Isn't it? If I'm wrong, why speed decreases?
Title: Re: Help about implementing the perturbation theory
Post by: youhn on April 27, 2015, 01:56:36 PM
I do not fully understand, but will try to help.
Do you really mean that the calculation of the same image (same location, same resolution) is actually FASTER without perturbation ?! In that case, I think you've made a mistake somewhere. But please give some more insight into what you do, so we might be able to explain.
Question to help:
What is the location? (can you provide a KFR file?)
What is the resolution?
How much time does calculation WITHOUT perturbation take?
How much time does calculation WITH perturbation take in Kalles Fraktaler?
How much time does calculation WITH perturbation take in your software?
Do you really mean that the calculation of the same image (same location, same resolution) is actually FASTER without perturbation ?! In that case, I think you've made a mistake somewhere. But please give some more insight into what you do, so we might be able to explain.
Question to help:
What is the location? (can you provide a KFR file?)
What is the resolution?
How much time does calculation WITHOUT perturbation take?
How much time does calculation WITH perturbation take in Kalles Fraktaler?
How much time does calculation WITH perturbation take in your software?
Title: Re: Help about implementing the perturbation theory
Post by: therror on April 27, 2015, 02:46:50 PM
All conditions is the same. Number of iterations, location, zoom, resolution. The only difference - few lines for perturbation and float instead double. (I never use bignum arithmetic)
Location:
range=4.9304E-8;
xcenter=-0.743643900055;
ycenter=0.131825890901.
Resolution 1280х1024. At 320x 256 usual calculation faster too.
w/o perturbation 221 second.
with perturbation 380 second in my program.
To compare with other software - good idea ;D
Location:
range=4.9304E-8;
xcenter=-0.743643900055;
ycenter=0.131825890901.
Resolution 1280х1024. At 320x 256 usual calculation faster too.
w/o perturbation 221 second.
with perturbation 380 second in my program.
To compare with other software - good idea ;D
Title: Re: Help about implementing the perturbation theory
Post by: 3dickulus on April 27, 2015, 03:26:52 PM
you can find a C++ port of SFT here (http://www.digilanti.org/cudabrot/) and some chatter about it here (http://www.fractalforums.com/programming/cuda-y-a-m-z/60/)
this was ported directly from the java source and knighty helped (a lot, tnx btw) with a few bugs, it uses ArPrec lib for the full precision part.
the deepest I went with this code was E-2346 ish ;)
this was ported directly from the java source and knighty helped (a lot, tnx btw) with a few bugs, it uses ArPrec lib for the full precision part.
the deepest I went with this code was E-2346 ish ;)
Title: Re: Help about implementing the perturbation theory
Post by: quaz0r on April 29, 2015, 12:45:55 AM
he only said hes using perturbation, not series approximation. also hes staying within the limits of hardware floating point. so of course the extra overhead and unnecessary screwing around of perturbation is going to be slower than regular iteration. and even if hes also using series approximation, lots of low-magnification locations will not benefit much from it, and the overhead of all that could still be slower than regular iteration.
Title: Re: Help about implementing the perturbation theory
Post by: cKleinhuis on April 29, 2015, 01:06:51 AM
he only said hes using perturbation, not series approximation. also hes staying within the limits of hardware floating point. so of course the extra overhead and unnecessary screwing around of perturbation is going to be slower than regular iteration. and even if hes also using series approximation, lots of low-magnification locations will not benefit much from it, and the overhead of all that could still be slower than regular iteration.
this is what is second as well, the pertubation kicks in at ridiculous depth ;) it might appear that it improves double->float as well, use max double range and use float for the rest, the thing here is that modern cpus do not have any significant double precision calculation time differencies anymore regarding float/double !!!
Title: Re: Help about implementing the perturbation theory
Post by: therror on May 14, 2015, 06:40:57 PM
Thanks, dear friends! I understand that perturbation don't help me.. But I have implemented border tracing(with rectangular borders)! I made skipped points visible, it is not error. Speed is much more!