program MatrixMultiplicationALGLIB;
{$mode objfpc}{$H+}
uses
SysUtils,
Math,
DateUtils,
Classes,
xalglib;
const
N = 1000;
type
TMatrixMultiplier = class(TThread)
public
A, B, C: xalglib.TMatrix;
OffsetA, OffsetB, OffsetC: Integer;
Rows, Cols: Integer;
procedure Execute; override;
end;
var
A, B, C, CBase: xalglib.TMatrix;
A1, A2, B1, B2: xalglib.TMatrix;
C11, C12, C21, C22: xalglib.TMatrix;
i, j, halfN, restN: integer;
StartTime, EndTime: TDateTime;
Threads: array[0..3] of TThread;
procedure TMatrixMultiplier.Execute;
begin
try
if (Length(A) > 0) and (Length(B) > 0) then
RMatrixGEMM(Rows, Cols, N, 1.0,
A, OffsetA, 0, 0,
B, 0, OffsetB, 0,
0.0, C, OffsetC, 0);
except
on E: Exception do
Writeln('Thread error: ', E.ClassName, ': ', E.Message);
end;
end;
procedure InitializeAlignedMatrix(var M: xalglib.TMatrix; Rows, Cols: Integer);
var
i: Integer;
begin
SetLength(M, Rows);
for i := 0 to Rows - 1 do
begin
SetLength(M[i], Cols);
FillChar(M[i][0], Cols * SizeOf(Double), 0);
end;
end;
procedure GenerateRandomMatrix(var M: xalglib.TMatrix; Rows, Cols: Integer);
var
i, j: Integer;
begin
InitializeAlignedMatrix(M, Rows, Cols);
for i := 0 to Rows - 1 do
for j := 0 to Cols - 1 do
M[i][j] := Random(10);
end;
function MatricesEqual(const Mat1, Mat2: xalglib.TMatrix): Boolean;
var
i, j: Integer;
begin
Result := True;
for i := 0 to High(Mat1) do
for j := 0 to High(Mat1[i]) do
if not SameValue(Mat1[i][j], Mat2[i][j], 1e-9) then
Exit(False);
end;
begin
try
// Initialize matrices with alignment
GenerateRandomMatrix(A, N, N);
GenerateRandomMatrix(B, N, N);
InitializeAlignedMatrix(C, N, N);
InitializeAlignedMatrix(CBase, N, N);
// Standard multiplication
StartTime := Now;
RMatrixGEMM(N, N, N, 1.0, A, 0, 0, 0, B, 0, 0, 0, 0.0, C, 0, 0);
EndTime := Now;
Writeln('Standard: ', MilliSecondsBetween(EndTime, StartTime), ' ms');
// Parallel block multiplication
StartTime := Now;
halfN := N div 2;
restN := N - halfN;
// Split matrix A
SetLength(A1, halfN);
SetLength(A2, restN);
for i := 0 to halfN - 1 do
begin
A1[i] := Copy(A[i], 0, N);
A2[i] := Copy(A[i + halfN], 0, N);
end;
// Split matrix B
SetLength(B1, N, halfN);
SetLength(B2, N, halfN);
for i := 0 to N - 1 do
begin
B1[i] := Copy(B[i], 0, halfN);
B2[i] := Copy(B[i], halfN, halfN);
end;
// Initialize result matrices
InitializeAlignedMatrix(C11, halfN, halfN);
InitializeAlignedMatrix(C12, halfN, restN);
InitializeAlignedMatrix(C21, restN, halfN);
InitializeAlignedMatrix(C22, restN, restN);
// Configure threads
Threads[0] := TMatrixMultiplier.Create(True);
with TMatrixMultiplier(Threads[0]) do begin
A := A1; B := B1; C := C11;
OffsetA := 0; OffsetB := 0; OffsetC := 0;
Rows := halfN; Cols := halfN;
FreeOnTerminate := False;
end;
Threads[1] := TMatrixMultiplier.Create(True);
with TMatrixMultiplier(Threads[1]) do begin
A := A2; B := B2; C := C22;
OffsetA := 0; OffsetB := 0; OffsetC := 0;
Rows := restN; Cols := restN;
FreeOnTerminate := False;
end;
Threads[2] := TMatrixMultiplier.Create(True);
with TMatrixMultiplier(Threads[2]) do begin
A := A1; B := B2; C := C12;
OffsetA := 0; OffsetB := 0; OffsetC := 0;
Rows := halfN; Cols := restN;
FreeOnTerminate := False;
end;
Threads[3] := TMatrixMultiplier.Create(True);
with TMatrixMultiplier(Threads[3]) do begin
A := A2; B := B1; C := C21;
OffsetA := 0; OffsetB := 0; OffsetC := 0;
Rows := restN; Cols := halfN;
FreeOnTerminate := False;
end;
// Start threads
for i := 0 to 3 do Threads[i].Start;
for i := 0 to 3 do Threads[i].WaitFor;
// Assemble results
for i := 0 to halfN - 1 do
begin
if (halfN > 0) and (Length(C11[i]) = halfN) then
Move(C11[i][0], C[i][0], halfN * SizeOf(Double));
if (restN > 0) and (Length(C12[i]) = restN) then
Move(C12[i][0], C[i][halfN], restN * SizeOf(Double));
end;
for i := 0 to restN - 1 do
begin
if (halfN > 0) and (Length(C21[i]) = halfN) then
Move(C21[i][0], C[i + halfN][0], halfN * SizeOf(Double));
if (restN > 0) and (Length(C22[i]) = restN) then
Move(C22[i][0], C[i + halfN][halfN], restN * SizeOf(Double));
end;
EndTime := Now;
Writeln('Parallel: ', MilliSecondsBetween(EndTime, StartTime), ' ms');
// Validation
StartTime := Now;
RMatrixGEMM(N, N, N, 1.0, A, 0, 0, 0, B, 0, 0, 0, 0.0, CBase, 0, 0);
EndTime := Now;
WriteLn('Validation: ', MatricesEqual(C, CBase));
// Cleanup
for i := 0 to 3 do Threads[i].Free;
except
on E: Exception do
Writeln('Main error: ', E.ClassName, ': ', E.Message);
end;
Readln;
end.
Recent
