为什么使用 Parallel.For 的 LU 分解不起作用？

Question

我正在尝试使用 Doolittle 算法求解 LU 分解——根据这个 document。没有并行化，代码工作正常。但是，我想并行编写此代码运行 - 尝试创建一个并行的外部循环和两个内部循环。不幸的是我仍然缺少一些东西。如果我尝试首先使外部循环运行并行，我会得到一个糟糕的结果。

我试图检测可能发生碰撞的地方。后来我把那些地方都锁了，还是没有收到正确的结果。我将它们作为注释添加到复制的代码中。我做错了什么，我需要锁定其他地方吗？

外循环的正确实现是什么？

内部循环会是什么样子？

提前致谢。

算法的实现（顺序）

        //upper triangle
        var upper = new double[arr.GetLength(0), arr.GetLength(0)];
        //lower triangle
        var lower = new double[arr.GetLength(0), arr.GetLength(0)];

        //field initialization
        for (int i = 0; i < n; i++)
        {
            for (int j = i; j < n; j++)
                upper[i, j] = arr[i, j];
            for (int j = i + 1; j < n; j++)
                lower[j, i] = arr[j, i];
            lower[i, i] = 1;
        }

        for(int i=0; i<n; i++)
        {
            for (int j = i; j < n; j++)
            {
                for (int k = 0; k < i; k++)
                { 
                    upper[i, j] = upper[i, j] - (lower[i, k] * upper[k, j]);
                }
            }

            for (int j = i + 1; j < n; j++)
            {
                for (int k = 0; k < i; k++)
                {
                    lower[j, i] = lower[j, i] - (lower[j, k] * upper[k, i]);
                }
                    lower[j, i] = lower[j, i] / upper[i, i];
            }
        }

算法的实现（并行）

        //upper triangle
        var upper = new double[arr.GetLength(0), arr.GetLength(0)];
        //lower triangle
        var lower = new double[arr.GetLength(0), arr.GetLength(0)];

        //field initialization
        for (int i = 0; i < n; i++)
        {
            for (int j = i; j < n; j++)
                upper[i, j] = arr[i, j];
            for (int j = i + 1; j < n; j++)
                lower[j, i] = arr[j, i];
            lower[i, i] = 1;
        }

        //making outer loop parallel
        Parallel.For(0, n, (i, state) =>
        {
            //possibly make this loop parallel also
            for (int j = i; j < n; j++)
            {
                for (int k = 0; k < i; k++)
                {
                    //lower[i,k] is it potential problem?
                    /*
                     * I tried this solution
                     * double a;
                     * lock(lowerLock){
                     *   a = lower[i,k];
                     * }
                     * upper[i, j] = upper[i, j] - (a * upper[k, j])
                     */
                    upper[i, j] = upper[i, j] - (lower[i, k] * upper[k, j]);
                }
            }

            //possibly make this loop parallel also
            for (int j = i + 1; j < n; j++)
            {
                for (int k = 0; k < i; k++)
                {
                    //upper [k,i] is it potential problem?
                    /*
                     * I tried this solution
                     * double b;
                     * lock(upperLock){
                     *   b = upper[k, i];
                     * }
                     * lower[j, i] = lower[j, i] - (lower[j, k] * b);
                     */
                    lower[j, i] = lower[j, i] - (lower[j, k] * upper[k, i]);
                }
                    lower[j, i] = lower[j, i] / upper[i, i];
            }
        });

顺序正确的结果

Concatenation  Upper triangle  Lower triangle
 2 -1 -2       2 -1 -2          1  0  0
-2  4 -1       0  4 -1         -2  1  0
-2 -1  3       0  0  3         -2 -1  1

并行错误结果

Concatenation  Upper triangle    Lower triangle
 2 -1 -2       2 -1 -2             1  0  0
-2  4 -1       0  4 -1            -2  1  0
-2 -1  3       0  0 10 -->BAD     -2 -1  1

编辑我试图用一把锁锁定所有通往田野的方法。我知道以这种方式失去所有并行化。然而，我想至少获得正确的结果，不幸的是没有成功。

 static object mylock = new object();
            //making outer loop parallel
            Parallel.For(0, n, (i, state) =>
            {
                for (int j = i; j < n; j++)
                {
                    for (int k = 0; k < i; k++)
                    {
                        lock (mylock)
                        {
                            upper[i, j] = upper[i, j] - (lower[i, k] * upper[k, j]);
                        }
                    }
                }

                for (int j = i + 1; j < n; j++)
                {
                    for (int k = 0; k < i; k++)
                    {
                        lock (mylock)
                        {
                            lower[j, i] = lower[j, i] - (lower[j, k] * upper[k, i]);
                        }
                    }
                    lock (mylock)
                    {
                        lower[j, i] = lower[j, i] / upper[i, i];
                    }
                }
            });

Answer 1

并行循环写入同一个数组，对吧？

upper[i, j] = upper[i, j] - (lower[i, k] * upper[k, j]);

但未定义，何时哪个循环将写入您数组中的某个位置。因此，两个循环不会写入同一个索引，而是从一个索引中读取，而另一个循环可能已经写入了该索引。您不能以这种方式并行化您的算法。

为什么使用 Parallel.For 的 LU 分解不起作用？

Why does LU decomposition using Parallel.For not work?

c#

parallel-processing

linear-algebra

parallel.foreach

matrix-decomposition