Fortran - 逻辑索引 [英] Fortran - Logical Indexing

问题描述

假设我有一个矩阵 A ，它是（mxn）和一个向量 B ，（mx 1）。此向量 B 是0和1的向量。

Suppose I have a matrix A which is (m x n) and a vector B which is (m x 1). This vector B is a vector of zeros and ones.

同样让标量 s 为 B 。

我想创建一个矩阵 C ，这是 sxn 对应于 B 等于1的行，以及向量 D ，即sx 1 ，这些元素的位置在 A 。

I want to create a matrix C which is s x n corresponding to the rows of B which equal 1, and a vector D which is s x 1, with the position of those elements in A.

Take as an example: 

      A = [1, 2, 3; 
           4, 5, 6; 
           7, 8, 9; 
           10, 11, 12; 
           13, 14, 15 ]

        B = [0; 1; 0; 1; 1]

    Then:
C = [ 4, 5, 6; 
     10, 11, 12; 
     13, 14, 15 ]
and
D = [2; 
     4
     5]

截至目前，我的fortran代码如下： / p>

As of now my fortran code looks like:

PROGRAM test1
IMPLICIT NONE

 REAL, DIMENSION(5,3) :: A
INTEGER, DIMENSION(5,1) :: B = 0
INTEGER :: i, j, k

k = 1

!Create A matrix
do i=1,5
    do j=1,3
        A(i,j) = k
        k = k+1
    end do
end do

!Create B matrix
B(2,1) = 1
B(4,1) = 1
B(5,1) = 1


end program

在matlab中我可以创建C简单地通过以类似的方式制作：C = A（逻辑（B），:)和D.

In matlab I could create C simply by making: C = A(logical(B),:), and D in a similar way.

我怎么能在fortran中做，避免循环？我查看了 where 和 forall 语句，但它们并不是我想要的。

How can I do it in fortran, avoiding loops? I have looked into the where and forall statements, but they are not exactly what I am looking for.

推荐答案

正如@francescalus建议的那样， PACK 内在函数是Fortran相当于Matlab逻辑切片，但仅限于部分。有两种类型的Matlab逻辑索引：

As @francescalus suggested, the PACK intrinsic function is the Fortran equivalent of Matlab logical slicing, but only in part. There are two types of Matlab logical indexing:

• 整个数组逻辑索引：掩码必须具有与数组相同的形状和返回值的等级为1（Matlab中的向量）。这是因为trues的位置是任意的，因此你不能保证基本上将孔戳到矩阵中的结果将是矩形的。 这是PACK在Fortran中的作用。

• Whole array logical indexing: the mask must have the same shape as the array and the returned value is of rank 1 (a vector in Matlab). This is so because the position of the trues is arbitrary, and thus you cannot guarantee that the result of, basically, poking holes into a matrix will be rectangular. This is what PACK does in Fortran.

c = a(a < 1); % Matlab: a(m,n) -> c(s)
C = PACK(A, A < 1) ! Fortran: A(m,n) -> C(s)

单维逻辑索引： mask必须是1-D，用于在数组的单个维度内选择。其他尺寸可以整体选择，也可以自己编制索引。 这就是您想要的。

Single-dimension logical indexing: the mask must be 1-D, and is used to select inside a single dimension of the array. Other dimensions can be selected whole or with indexing of their own. This is what you want.

b = a(:,1) > 0; % a(m,n) gives logical b(m)
c = a(b,:); % a(m,n) selected with b(m) -> c(s,n)

然而，PACK不直接承认此用法。 @ high-performance-mark的解决方案向您展示了如何执行此专长： SPREAD 基本上是 repmat ，所以他的解决方案将是在Matlab中看起来像这样：

However, PACK does not admit this usage directly. @high-performance-mark's solution shows you how to perform this feat: SPREAD is basically repmat, so his solution would look like this in Matlab:

b = a(:,1) > 0; % a(m,n) gives logical b(m)
bMat = repmat(b, 1, size(a,2)); % n copies of b(m) -> bMat(m,n)
c = reshape(a(bMat), [sum(b), size(a,2)]); % a(m,n) -> c(s,n)

需要进行最终重塑，因为 a（ bMat）不是矩阵，而是由于使用整个矩阵选择形式而导致大小为s * n的向量。另一个答案的Fortran代码正是如此：

Where the final reshape is required because a(bMat) is not a matrix, but a vector of size s*n due to the usage of the whole-matrix selection form. The other answer's Fortran code is doing exactly that:

c = RESHAPE( &
        PACK(a, & ! The matrix we are selecting from
            SPREAD(b==1, ! The == 1 is because you were using an INTEGER array, not LOGICAL
                   dim=2, ncopies=SIZE(a,2)) & ! This is the repmat-like code
        ), & ! The result of this PACK is cVec(s*n)
        [COUNT(b==1),SIZE(a,2)]) ! Reshape into (s,n)

分配给D的代码非常相似，但不是使用我们从一个即时生成的数组中选择包含从1到A的第一个维度长度的数字（与B的大小相同）。这次不需要掩码的任何SPREAD或RESHAPE到结果，因为源数组是1-D。

The code that assigns to D is very similar, but instead of using A we are selecting from an array that is generated on-the-fly to contains numbers from 1 to the length of the first dimension of A (which is the same as the size of B). This time there does not need to be any SPREAD of the mask or RESHAPE to the result because the source array is 1-D.

顺便说一下，Fortran 做直接支持整数向量索引，因此您可以使用代码生成D向量（带有B的真实元素的索引），然后然后执行：

Incidentally, Fortran does support integer-vector indexing directly, so you could use the code to generate the D vector (with the indices of the true elements of B) first and then do:

C = A(D,:) ! Works the same in Fortran!

节省自己的传播和重塑。

saving yourself the spreading and reshaping.

这篇关于Fortran - 逻辑索引的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！