我需要使用 pthreads 实现以下高斯消除算法的并行版本。
procedure GAUSSIAN ELIMINATION (A, b, y)
begin
for k := 0 to n − 1 do /* Outer loop */
begin
for j := k + 1 to n − 1 do
A[k, j] := A[k, j]/A[k, k]; /* Division step */
y[k] := b[k]/A[k, k];
A[k, k] := 1;
for i := k + 1 to n − 1 do
begin
for j := k + 1 to n − 1 do
A[i, j] := A[i, j] − A[i, k] × A[k, j]; /* Elimination step */
b[i] := b[i] − A[i, k] × y[k];
A[i, k] := 0;
endfor; /* Line 9 */
endfor; /* Line 3 */
end GAUSSIAN ELIMINATION
我理解顺序实现和 pthreads,但没有得到关于如何为它的并行版本构建逻辑的单一提示(线程工作分配标准,循环并行化等)。任何线索或起点都会帮助我继续。
最佳答案
矩阵中每一行的工作需要完成所有前面的行,所以你不能那样划分工作。
但是,在一行中,每一列都可以并行处理(需要注意的是第 k
列的原始值必须保存并用于其他列的计算) .这对应于您的 j
值。
我相信你可以重新安排算法使这更容易,这样在 j
上只有一个循环:
procedure GAUSSIAN ELIMINATION (A, b, y)
begin
for k := 0 to n − 1 do /* Outer loop */
begin
for j := k + 1 to n − 1 do
begin
A[k, j] := A[k, j]/A[k, k]; /* Division step */
for i := k + 1 to n − 1 do
A[i, j] := A[i, j] − A[i, k] × A[k, j]; /* Elimination step */
endfor;
y[k] := b[k]/A[k, k];
A[k, k] := 1;
for i := k + 1 to n − 1 do
begin
b[i] := b[i] − A[i, k] × y[k];
A[i, k] := 0;
endfor;
endfor;
end GAUSSIAN ELIMINATION
请注意,j
上的循环体仅访问 j
和 k
列中的值 - 这是可以完成的循环在平行下。然后您可以进一步注意到外循环的第二部分不依赖于第一部分,因此您可以将外循环分成两部分:
procedure GAUSSIAN ELIMINATION (A, b, y)
begin
for k := 0 to n − 1 do /* Outer loop */
begin
for j := k + 1 to n − 1 do
begin
A[k, j] := A[k, j]/A[k, k]; /* Division step */
for i := k + 1 to n − 1 do
A[i, j] := A[i, j] − A[i, k] × A[k, j]; /* Elimination step */
endfor;
endfor;
for k := 0 to n − 1 do /* Outer loop, second pass */
begin
y[k] := b[k]/A[k, k];
A[k, k] := 1;
for i := k + 1 to n − 1 do
begin
b[i] := b[i] − A[i, k] × y[k];
A[i, k] := 0;
endfor;
endfor;
end GAUSSIAN ELIMINATION
您可以预先创建线程,每个线程负责对从 0
到 j
值的子集执行 j
循环n - 1
。这些线程在处理每一行后将需要一个同步屏障,因为处理下一行需要前一行的所有列的结果。为此,您可以使用 pthread_barrier_wait()
。
您会注意到并非每一列(j
的值)都有相同的工作。列 j
被该循环处理 j
次(因此列 0 执行 0 次,而列 n - 1
执行 n - 1
次)。您需要为线程分配列号,以便每个线程对每一行的工作量尽可能接近 - 这可以通过以循环方式分配列来完成。例如。如果你有三个线程 A、B 和 C 以及从 0 到 9 的 10 列,你将分配它们:
Thread A: 0, 3, 6, 9
Thread B: 1, 4, 7,
Thread C: 2, 5, 8,
线程函数看起来像这样:
for k := 0 to n − 1 do /* Outer loop */
begin
call pthread_barrier_wait(row_barrier);
for j := k + 1 + thread_number to n − 1 step n_threads do
begin
A[k, j] := A[k, j]/A[k, k]; /* Division step */
for i := k + 1 to n − 1 do
A[i, j] := A[i, j] − A[i, k] × A[k, j]; /* Elimination step */
endfor;
endfor;
call pthread_barrier_wait(phase1_barrier);
主要功能看起来像这样:
procedure GAUSSIAN ELIMINATION (A, b, y)
begin
call start_threads;
call pthread_barrier_wait(phase1_barrier);
for k := 0 to n − 1 do /* Outer loop, second pass */
begin
y[k] := b[k]/A[k, k];
A[k, k] := 1;
for i := k + 1 to n − 1 do
begin
b[i] := b[i] − A[i, k] × y[k];
A[i, k] := 0;
endfor;
endfor;
end GAUSSIAN ELIMINATION
关于c - 使用 pthreads 并行实现高斯消元,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/31697173/