[SOLVED] Calculation error when parallelized

[BarrOff]

Hello,

I recently stumbled across some Delphi code creating a visualization of the Mandelbrot set (https://en.wikipedia.org/wiki/Mandelbrot_set).
I wrote this for a math course at my school, back in the days.
After importing it into Lazarus it was still working.
Calculating every pixel of a canvas takes some time so I decided to parallelize the code.
Every point only depends on its indices in the 2d array so multiple threads should not interfere with each other.
Following the documentation about parallel procedures (https://wiki.freepascal.org/Parallel_procedures) from the wiki I was able to parallelize the code using MTProcs.
However I stumbled across a strange effect:

As soon as I use more than one thread (the number of threads to use can be set by a TEdit field) some pixels get a wrong colour (see attached image).
Rerunning the calculations the number and locations of bad pixels changes.
Using a single thread works fine though, so the basic implementation is right.
Therefore I think it must have to do something with the parallelization.

After trying to find the cause for quite some time now and not succeeding I would like to ask for help.
Can someone please help me and point out what is causing the problem after parallelization?

Thank you very much!

Here is the code:

Code: Pascal  [Select][+][-]

1. unit Mandelbrot;
2.  
3. {$IFDEF FPC}
4.   {$MODE Delphi}
5. {$ENDIF}
6.  
7. interface
8.  
9. uses
10. {$IFnDEF FPC}
11.   Windows,
12. {$ELSE}
13.   LCLIntf, LCLType,
14. {$ENDIF}
15.   SysUtils, Classes, Graphics, Controls, Forms,
16.   Dialogs, StdCtrls, ExtCtrls,
17.   MTProcs;
18.  
19. type
20.  
21.   { TForm1 }
22.  
23.   TForm1 = class(TForm)
24.     PaintBox1: TPaintBox;
25.     Edit3: TEdit;
26.     Button1: TButton;
27.     Panel1: TPanel;
28.     Edit1: TEdit;
29.     Label1: TLabel;
30.     procedure Button1Click(Sender: TObject);
31.   private
32.     { Private-Deklarationen }
33.   public
34.     { Public-Deklarationen }
35.   end;
36.  
37. var
38.   Form1: TForm1;
39.  
40. implementation
41.  
42. {$IFnDEF FPC}
43.   {$R *.dfm}
44. {$ELSE}
45.   {$R *.lfm}
46. {$ENDIF}
47.  
48. function RGB2TColor(const R, G, B: byte): integer;
49. begin
50.   Result := R + G shl 8 + B shl 16;
51. end;
52.  
53. procedure ColourSinglePixel(Index: PtrInt; Data: Pointer; Item: TMultiThreadProcItem);
54. var
55.   x, y, xx, yy, c, c1, c2: double;
56.   z, i, j, k: integer;
57.   myColor: TColor;
58.   r, g, b: byte;
59. begin
60.   c := PDouble(Data)^;
61.   z := 0;
62.   { convert index into x,y coordinates of the canvas }
63.   i := ((Index - 1) div Form1.PaintBox1.Height) + 1;
64.   j := ((Index - 1) mod Form1.PaintBox1.Height) + 1;  
65.   x := i / Form1.PaintBox1.Width;
66.   y := j / Form1.PaintBox1.Height;
67.   c1 := c * x;
68.   c2 := c * y;
69.   k := -1;
70.   while (z < 65536) and (x < 4) and (y < 4) do
71.   begin
72.     xx := x * x - y * y + c1;
73.     yy := 2 * y * x + c2;
74.     x := xx;
75.     y := yy;
76.     z := z + 1;
77.   end;
78.   if z >= 256 then
79.   begin
80.     k := z div 256;
81.     r := z div (k + 10);
82.     g := z div (k + 15);
83.     b := z div (k + 20);
84.   end
85.   else
86.   begin
87.     r := z div 2;
88.     g := z div 3 + a div 2;
89.     b := z div 5 + b div 3;
90.   end;
91.   myColor := RGB2TColor(r, g, b);
92.   Form1.PaintBox1.Canvas.Pixels[i, j] := myColor;
93. end;
94.  
95. procedure TForm1.Button1Click(Sender: TObject);
96. var
97.   c: double;
98.   pc: PDouble;
99. begin
100.   c := StrToFloat(Edit3.Text);
101.   pc := @c;
102.   ProcThreadPool.MaxThreadCount := StrToInt(Edit1.Text);
103.   ProcThreadPool.DoParallel(@ColourSinglePixel, 1, PaintBox1.Height * PaintBox1.Width, pc);
104. end;
105.  
106. end.
107.  

[Fibonacci]

If you change that to direct memory access you wont even have to use threads

Code: Pascal  [Select][+][-]

1. Form1.PaintBox1.Canvas.Pixels[i, j] := myColor;

Share full project, with form, compilable

[MarkMLl]

As soon as I use more than one thread (the number of threads to use can be set by a TEdit field) some pixels get a wrong colour (see attached image).
Rerunning the calculations the number and locations of bad pixels changes.
Using a single thread works fine though, so the basic implementation is right.
Therefore I think it must have to do something with the parallelization.

Are you trying to write to a graphical control from a background thread? If so /don't/: that's never been reliable either in Delphi or Lazarus.

I've not looked at Mandelbrot sets for around 40 years (back when an 8087 was still considered a fairly cool thing to have), but my recollection is that it's "embarrassingly parallel". If that's the case, then allocate a graphical buffer per thread, and when work is complete write it to the graphical control via Synchronize().

MarkMLl

[BarrOff]

Thank you for your answers.

@Fibonacci: Could you please elaborate a bit on

If you change that to direct memory access you wont even have to use threads

I don't understand what you mean with that.
Attached you find a mwe.

@MarkMLI: You are right, the background threads would populate the PaintBox.Canvas.Pixels.
I did not know that this is unreliable.

Could you explain/give an example or point me to some documentation about the "graphical buffers" you are talking about, please?

[Fibonacci]

PaintBox doesnt have ScanLine

Here, more info:
https://wiki.freepascal.org/Fast_direct_pixel_access

I have at least added a critical section

Code: Pascal  [Select][+][-]

1. implementation
2.  
3. {$R *.lfm}
4.  
5. var
6.   cs: TRTLCriticalSection;
7.  
8. function RGB2TColor(const R, G, B: byte): integer;
9. begin
10.   Result := R + G shl 8 + B shl 16;
11. end;
12.  
13. procedure ColourSinglePixel(Index: PtrInt; Data: Pointer; Item: TMultiThreadProcItem);
14. var
15.   x, y, xx, yy, c, c1, c2: double;
16.   z, i, j, k: integer;
17.   myColor: TColor;
18.   a, b, d: byte;
19. begin
20.   c := PDouble(Data)^;
21.   z := 0;
22.   j := ((Index - 1) mod Form1.PaintBox1.Height) + 1;
23.   i := ((Index - 1) div Form1.PaintBox1.Height) + 1;
24.   x := i / Form1.PaintBox1.Width;
25.   y := j / Form1.PaintBox1.Height;
26.   c1 := c * x;
27.   c2 := c * y;
28.   k := -1;
29.   while (z < 65536) and (x < 4) and (y < 4) do
30.   begin
31.     xx := x * x - y * y + c1;
32.     yy := 2 * y * x + c2;
33.     x := xx;
34.     y := yy;
35.     z := z + 1;
36.   end;
37.   if z >= 256 then
38.   begin
39.     k := z div 256;
40.     a := z div (k + 10);
41.     b := z div (k + 15);
42.     d := z div (k + 20);
43.   end
44.   else
45.   begin
46.     a := z div 2;
47.     b := z div 3 + a div 2;
48.     d := z div 5 + b div 3;
49.   end;
50.   myColor := RGB2TColor(a, b, d);
51.   try
52.     EnterCriticalSection(cs);
53.     Form1.PaintBox1.Canvas.Pixels[i, j] := myColor;
54.   finally
55.     LeaveCriticalSection(cs);
56.   end;
57. end;
58.  
59. procedure TForm1.Button1Click(Sender: TObject);
60. var
61.   c: double;
62.   pc: PDouble;
63. begin
64.   c := StrToFloat(Edit3.Text);
65.   pc := @c;
66.   ProcThreadPool.MaxThreadCount := StrToInt(Edit1.Text);
67.   InitCriticalSection(cs);
68.   ProcThreadPool.DoParallel(@ColourSinglePixel, 1, PaintBox1.Height * PaintBox1.Width, pc);
69.   DoneCriticalSection(cs);
70. end;

Below without and with critical section

[MarkMLl]

I have at least added a critical section

I'm not an accomplished graphical programmer, but I'd caution that there's a risk here.

Critical sections protect against a different thread attempting to access a resource, but not against the same thread.

The particular scenario that has bitten me is that:

i) A Synchronize() handler, in the context of the main thread, enters a critical section.

ii) The handler calls Application.ProcessMessages in an attempt to force an update of all visible components.

iii) The message handler loop passes control to another event handler in the context of the main thread.

iv) The event handler is allowed to enter the same critical section.

v) Havoc ensues.

That one took me several months of intermittent faults to track down, and these days I tend to use /both/ a critical section and an atomic increment/decrement to protect really critical resources.

With- again- the caveat that I don't consider myself an accomplished graphical programmer ** , in the current case I'd have thought that a simple array of pixels would have been adequate as an intermediate buffer.

** OK, I've been doing it since the 1980s but only on an occasional basis.

MarkMLl

[BarrOff]

Thanks to your comments I came to a solution which seems to work, see the code below.
Using documentation https://wiki.freepascal.org/Developing_with_Graphics#Direct_pixel_access and https://forum.lazarus.freepascal.org/index.php/topic,49855.msg362684.html?PHPSESSID=thmk0fhighki491m426042huo7#msg362684 I derived the following solution:

1. create a Bitmap with size of the canvas
2. MTProcs creates one job per line of the bitmap and distributes them over the threads
3. each thread calculates all values of a single line of the bitmap and updates it
4. once all jobs are done, draw the bitmap on the canvas

However I stumbled across another thing which threw me off:

When using multiple threads it seems that calling scanline on the bitmap sometimes returned nil.
I checked this using the debugger.
I circumvent this by retrying the scanline over and over in a while loop till it succeeds.
I could not find good documentation about the scanline function, though, so I'm not really sure why/how this happens?
Can you explain this to me or point me to the documentation for scanline?

Here is the current version:

Code: Pascal  [Select][+][-]

1. unit Mandelbrot;
2.  
3. {$IFDEF FPC}
4.   {$MODE Delphi}
5. {$ENDIF}
6.  
7. interface
8.  
9. uses
10. {$IFnDEF FPC}
11.   Windows,
12. {$ELSE}
13.   LCLIntf, LCLType,
14. {$ENDIF}
15.   SysUtils, Classes, Graphics, Controls, Forms,
16.   Dialogs, StdCtrls, ExtCtrls,
17.   MTProcs;
18.  
19. type
20.   TMandelGraphic = record
21.     c: double;
22.     bm: TBitMap;
23.   end;
24.   PMandelGraphic = ^TMandelGraphic;
25.  
26.   { TForm1 }
27.  
28.   TForm1 = class(TForm)
29.     PaintBox1: TPaintBox;
30.     Edit3: TEdit;
31.     Button1: TButton;
32.     Edit1: TEdit;
33.     Label1: TLabel;
34.     procedure Button1Click(Sender: TObject);
35.   private
36.     { Private-Deklarationen }
37.   public
38.     { Public-Deklarationen }
39.   end;
40.  
41. var
42.   Form1: TForm1;
43.  
44. implementation
45.  
46.   {$R *.lfm}
47.  
48. function RGB2TColor(const R, G, B: byte): integer;
49. begin
50.   Result := R + G shl 8 + B shl 16;
51. end;
52.  
53. procedure ColourOneLine(Index: PtrInt; Data: Pointer; Item: TMultiThreadProcItem);
54. type
55.   Trgb24 = packed record
56.     b,g,r : byte;
57.   end;
58.   Trgb24scanline = array [word] of Trgb24;
59.   Prgb24scanline = ^Trgb24scanline;
60. var
61.   x, y, xx, yy, c, c1, c2: double;
62.   z, i, k: integer;
63.   r, g, b: byte;
64.   line: Prgb24scanline;
65.   pmg: PMandelGraphic;
66. begin
67.   pmg := PMandelGraphic(Data);
68.   c := pmg^.c;
69.   line := pmg^.bm.ScanLine[Index];
70.   for i := 0 to pmg^.bm.Width - 1 do
71.   begin
72.     x := (i + 1) / Form1.PaintBox1.Width;
73.     y := (Index + 1) / Form1.PaintBox1.Height;
74.     c1 := c * x;
75.     c2 := c * y;
76.     k := -1;
77.     z := 0;
78.     while (line = nil) do
79.     begin
80.       line := pmg^.bm.ScanLine[Index];
81.       sleep(10000);
82.     end;
83.     while (z < 65536) and (x < 4) and (y < 4) do
84.     begin
85.       xx := x * x - y * y + c1;
86.       yy := 2 * y * x + c2;
87.       x := xx;
88.       y := yy;
89.       z := z + 1;
90.     end;
91.     if z >= 256 then
92.     begin
93.       k := z div 256;
94.       r := z div (k + 10);
95.       g := z div (k + 15);
96.       b := z div (k + 20);
97.     end
98.     else
99.     begin
100.       r := z div 2;
101.       g := z div 3 + a div 2;
102.       b := z div 5 + b div 3;
103.     end;
104.     line^[i].R := r;
105.     line^[i].G := g;
106.     line^[i].B := b;
107.   end;
108. end;
109.  
110. procedure TForm1.Button1Click(Sender: TObject);
111. var
112.   c: double;
113.   MandelGraphic: TMandelGraphic;
114.   pmg: PMandelGraphic;
115. begin
116.   c := StrToFloat(Edit3.Text);
117.   try
118.     MandelGraphic.c := c;
119.     MandelGraphic.bm := TBitMap.Create;
120.     MandelGraphic.bm.PixelFormat := pf24bit;
121.     MandelGraphic.bm.Width  := PaintBox1.Width;
122.     MandelGraphic.bm.Height := PaintBox1.Height;
123.     pmg := @MandelGraphic;
124.     ProcThreadPool.MaxThreadCount := StrToInt(Edit1.Text);
125.     ProcThreadPool.DoParallel(@ColourOneLine, 0, PaintBox1.Height - 1, pmg);
126.     PaintBox1.Canvas.Draw(0,0,MandelGraphic.bm);
127.   finally
128.     MandelGraphic.bm.free;
129.   end;
130. end;
131.  
132. end.
133.  

[cdbc]

Hi
The bitmap, you're trying to get a (scan)line from, in more threads  IS ALSO A SHARED OBJECT Thus you need to protect it somehow from simultaneous access, at least with a 'crit-sect'...
this line of code should be protected:

Code: Pascal  [Select][+][-]

1. begin
2.   pmg := PMandelGraphic(Data);
3.   c := pmg^.c;
4.   line := pmg^.bm.ScanLine[Index];
5.   for i := 0 to pmg^.bm.Width - 1 do
6.   begin

e.g.:

Code: Pascal  [Select][+][-]

1.   EnterCriticalSection(cs);
2.   try
3.     line := pmg^.bm.ScanLine[Index];
4.   finally
5.     LeaveCriticalSection(cs);
6.   end;

That could just maybe, mend your codes behaviour...
...and please get rid of this "while line = nil do ... sleep(10000)" atrocity.
Regards Benny

[BarrOff]

@cdbc:

The threads won't acess the same memory regions, so I thought synchronization primitves would not be necessary, even though the variable is shared.
Replacing the while construct with a critical section indeed solves this problem, thank you for pointing this out!