我正在尝试使用OpenmP并行化以下代码段。所示部分基本上是用于流体流动/湍流模拟中隐式求解器的Thomas矩阵算法(矩阵求逆)的一种版本。我已将其从其确切形式中删除,以使问题易于理解。
REAL(KIND=DP),INTENT(INOUT),DIMENSION(1:10) :: A,C_old,C_new,D
INTEGER :: K
!$OMP PARALLEL DO SCHEDULE(STATIC) NUM_THREADS(3)&
!$OMP SHARED(A,D)
DO K = 2,9
C_new(K) = C_old(K)/(A(K)*C_new(K-1))
END DO
!$OMP END PARALLEL
C(1) = C(10) = 0
,因为它们位于固定边界(通道的顶部和底部)。因此,无需更新它们。
很明显,C(2)
的计算需要C(1)
的信息,而C(3)
的计算需要经过修改的C(2)
的信息,依此类推。换句话说,这是前向扫描(算法在https://en.wikipedia.org/wiki/Tridiagonal_matrix_algorithm中)。因此,对于三(03)个线程,我希望线程中的迭代分布为[1,2,3,4],[4,5,6,7],[7,8,9,10]
。但是,如果我只是执行SCHEDULE STATIC
,则将无法控制它,并且分布可能是跨越3个线程的[1,3],6],10]
。由于网格点4
和7
的信息未传递到其他线程,因此会产生错误的结果。
请让我知道是否有一种方法可以强制某个线程在没有任何竞争的情况下解决索引?说,
thread #1 solve indices [1,4]
等
thread #2 solve indices [4,7]