在扁平的fortran矩阵上循环(looping on flattened fortran matrices)

在扁平的fortran矩阵上循环(looping on flattened fortran matrices)

我有一些看起来像这样的代码：

DO I=0,500 arg1((I*54+1):(I*54+54)) = premultz*sinphi(I+1) ENDDO

简而言之，我有一个维54的数组前置。我有一个维数为501的数组sinphi。我想把sinult的第一个值乘以premultz的所有条目并将它存储在arg1的前54个条目中，然后是sinphi的第二个值乘以premultz的所有条目并将其存储在arg1的second54条目中，依此类推。

这些是扁平矩阵。 我为了速度而将它们弄平了，因为这个项目的主要目标之一是非常快速的代码。

我的问题是：在Fortran90中有更有效的编码这种计算方法吗？ 我知道Fortran有很多漂亮的数组操作可以做，我还没有完全清楚。

提前致谢。

I have a bit of code that looks like this:

DO I=0,500 arg1((I*54+1):(I*54+54)) = premultz*sinphi(I+1) ENDDO

In short, I have an array premultz of dimension 54. I have an array sinphi of dimension 501. I want to take the first value of of sinphi times all the entries of premultz and store it in the first 54 entries of arg1, then the second value of of sinphi times all the entries of premultz and store it in the second54 entries of arg1, and so on.

These are flattened matrices. I have flattened them in the interest of speed, as one of the primary goals of this project is very fast code.

My question is this: is there a more efficient way of coding this sort of calculation in Fortran90? I know that Fortran has a lot of nifty array operations that can be done that I'm not fully aware of.

Thanks in advance.

最满意答案

如果我有正确的话，这个表达式应该在一个语句中创建arg1

arg1 = reshape(spread(premultz,dim=2,ncopies=501)*& &spread(sinphi,dim=1,ncopies=54),[1,54*501])

我已经在这里修改了尺寸，这可能适合您的目的，也可能不适合您的目的。 内部表达式产生了premultz和sinphi的外积，然后将其重新整形为矢量。 您可能会发现需要重塑外部产品的转置，我没有仔细检查过。

但是，基于我对Fortran阵列内在函数的这种巧妙使用的经验，我怀疑这个或者Fortran的数组内在函数的其他大多数巧妙用法都会胜过你已经拥有的直接循环实现。 对于许多这些操作，编译器将生成数组的副本，并且复制数据相对昂贵。 当然，这是您可能想要测试的断言。

我会留给你来决定单行是否比循环更容易理解。 有时，数组语法的表达性在性能上是可接受的成本，有时则不然。

This expression, if I've got things right, ought to create arg1 in one statement

arg1 = reshape(spread(premultz,dim=2,ncopies=501)*& &spread(sinphi,dim=1,ncopies=54),[1,54*501])

I've hardwired the dimensions here, that may or may not suit your purposes. The inner expression generates the outer product of premultz and sinphi, which is then reshaped into a vector. You may find you need to reshape the transpose of the outer product, I haven't checked things very carefully.

However, based on my experience with this sort of clever use of Fortran's array intrinsics I doubt that this, or most other clever uses of Fortran's array intrinsics, will outperform the straightforward loop implementation you already have. For many of these operations the compiler is going to generate copies of arrays, and copying data is relatively expensive. Of course, this is an assertion you may want to test.

I'll leave it to you to decide if the one-liner is more comprehensible than the loops. Sometimes the expressivity of the array syntax comes at an acceptable cost in performance, sometimes it doesn't.