mandelbrot fractal
Vincent Diepeveen
diep at xs4all.nl
Tue Sep 8 12:37:05 UTC 2015
hi Folkert - you probably want to consider building a small piece of gpgpu
code for it in CUDA. The gpu's have a very small instruction set, yet they
do have quite some conditional instructions. So should be able to
calculate for a streamcore the next pixel even if other cores still are
busy with their current one, without 'if then else' type constructions.
AMD gpu's only have 24 bits available to you and the too limited
instruction set of OpenCL from conditional viewpoint seen with 32 bits is
basically 4x slower at AMD gpu's (and the top bits multiplication with 24
bits would also involve 4 units as there is no OpenCL equivalent for it
other than the one casted at 32x32 bits), whereas very cheap bit older
Nvidia gpu's have fast 32 bits integer calculations and you can pickup
even an older generation Tesla with 512 streamcores for just 192 dollar or
so on ebay.
If you are just interested in the 32 bits integer math then a gamers chip
would also be interesting f you underclock it a little.
That's teraflops of calculation power for under 200 euro and just a bit
of puzzling with conditional instructions to move to the next pixel if
current one has been calculated.
This is not too many lines of codes - one weekend of coding.
No CPU library can rival that.
Regards,
Vincent
BTW very interested in the program you got now to see how with more bits
precision mandelbrot looks like - where is it for download?
On Mon, 7 Sep 2015, folkert wrote:
> Hi,
>
> I'm trying to calculate the mandelbrot fractal using GMP.
> I've got it working and it is amazingly fast, even for data-types of
> 256 bits in size!
> Now what I'm curious about: can this piece of code be improved in some
> way? From a GMP-lib usage point of view I mean.
>
> // written by folkert at vanheusden.com
> // released under agpl v3.0
> #include <stdio.h>
> #include <gmpxx.h>
>
> void my_init_mpf(mpf_t *const what, const int n_bits)
> {
> mpf_init(*what);
> mpf_set_prec(*what, n_bits);
> }
>
> int main(int argc, char *argv[])
> {
> const int n_bits = 256;
>
> mpf_t x1, x2, y1, y2;
> mpf_init_set_d(x1, -2.0);
> mpf_init_set_d(y1, -2.0);
> mpf_init_set_d(x2, 2.0);
> mpf_init_set_d(y2, 2.0);
>
> const int max_it = 256;
> const int nx = 640;
> const int ny = 480;
>
> printf("P6 %d %d %d\n", nx, ny, max_it - 1);
>
> mpf_t two;
> mpf_init_set_d(two, 2.0);
>
> mpf_t nxg, nyg;
> mpf_init_set_ui(nxg, nx); // dummyx = xres
> mpf_init_set_ui(nyg, ny); // dummyy = yres
>
> mpf_t tx, ty;
> mpf_inits(tx, ty, NULL);
> mpf_sub(tx, x2, x1); // d2x = x2 - x1
> mpf_sub(ty, y2, y1); // d2y = y2 - y1
>
> mpf_t dx, dy;
> my_init_mpf(&dx, n_bits);
> my_init_mpf(&dy, n_bits);
> mpf_div(dx, tx, nxg); // dx = d2x / dummyx
> mpf_div(dy, ty, nyg); // dy = d2y / dummyy
>
> mpf_t cx, cy;
> my_init_mpf(&cx, n_bits);
> my_init_mpf(&cy, n_bits);
>
> mpf_t t, wx, wy, sqx, sqy;
> my_init_mpf(&t, n_bits);
> my_init_mpf(&wx, n_bits);
> my_init_mpf(&wy, n_bits);
> my_init_mpf(&sqx, n_bits);
> my_init_mpf(&sqy, n_bits);
>
> mpf_set(cy, y1);
> for(int y=0; y<ny; y++)
> {
> mpf_set(cx, x1);
>
> for(int x=0; x<nx; x++)
> {
> mpf_set_d(wx, 0.0);
> mpf_set_d(wy, 0.0);
>
> mpf_set_d(sqx, 0.0);
> mpf_set_d(sqy, 0.0);
>
> int it = 0;
>
> do {
> mpf_mul(t, wy, two); // y = 2 * y
> mpf_mul(t, t, wx); // * x
> mpf_add(wy, t, cy); // + cy
>
> //
>
> mpf_sub(wx, sqx, sqy); // x = x^2 - y^2
> mpf_add(wx, wx, cx); // x = z + cx
>
> //
>
> mpf_mul(sqx, wx, wx); // sqx = x^2
> mpf_mul(sqy, wy, wy); // sqy = y^2
>
> //
>
> mpf_add(t, sqx, sqy); // t = x^2 + y^2
> }
> while(mpf_cmp_d(t, 16.0) < 0 && ++it < max_it);
>
> double dwx = mpf_get_d(wx);
> double dwy = mpf_get_d(wy);
>
> if (it < max_it)
> printf("%c%c%c", int(it - dwx * 2.0), max_it - it, int(it + dwy * 2.0));
> else
> printf("%c%c%c", 0, 0, 0);
>
> mpf_add(cx, cx, dx);
> }
>
> mpf_add(cy, cy, dy);
> }
>
> return 0;
> }
>
>
> Please note that it emits a .pnm-file to stdout. So if you want to run
> it yourself, invoke it like this:
>
> ./a.out > test.pnm
>
> This test.pnm-file can be viewed with any Linux graphics viewer.
>
>
> regards,
>
> Folkert van Heusden
>
> --
> You've probably gotten really fed up with never winning in the Mega-
> Millions lottery. Well, weep no longer: www.smartwinning.info tells
> you everything that might help you deciding what numbers to choose.
> With nice graphs and pretty animations!
> _______________________________________________
> gmp-discuss mailing list
> gmp-discuss at gmplib.org
> https://gmplib.org/mailman/listinfo/gmp-discuss
>
More information about the gmp-discuss
mailing list