数组上的逐元素操作#
在使用子程序和函数时,有三种方法可以对数组执行逐元素操作
elemental
过程显式形状数组
为向量实现操作并编写简单的包装子程序(在内部使用
reshape
)用于每个数组形状
在第一种方法中,使用elemental
关键字创建如下函数
real(dp) elemental function nroot(n, x) result(y)
integer, intent(in) :: n
real(dp), intent(in) :: x
y = x**(1._dp / n)
end function
所有参数(输入和输出)必须是标量。然后可以使用此函数与任何(兼容)形状的数组一起使用,例如
print *, nroot(2, 9._dp)
print *, nroot(2, [1._dp, 4._dp, 9._dp, 10._dp])
print *, nroot(2, reshape([1._dp, 4._dp, 9._dp, 10._dp], [2, 2]))
print *, nroot([2, 3, 4, 5], [1._dp, 4._dp, 9._dp, 10._dp])
print *, nroot([2, 3, 4, 5], 4._dp)
输出将是
3.0000000000000000
1.0000000000000000 2.0000000000000000 3.0000000000000000 3.1622776601683795
1.0000000000000000 2.0000000000000000 3.0000000000000000 3.1622776601683795
1.0000000000000000 1.5874010519681994 1.7320508075688772 1.5848931924611136
2.0000000000000000 1.5874010519681994 1.4142135623730951 1.3195079107728942
在上述情况下,通常n
是一个参数,x
是任意形状的数组,但正如您所看到的,只要最终操作有意义(如果一个参数是数组,则其他参数必须是相同形状的数组或标量),Fortran 并不关心。如果它没有意义,您将收到编译器错误。
elemental
关键字隐含pure
关键字,因此过程必须是纯的。结果是elemental procedures
只能使用pure
过程并且没有副作用。
如果可以在内部使用数组操作使基本过程算法更快,或者由于某些原因参数必须是不兼容形状的数组,则应使用其他两种方法。可以使nroot
对向量进行操作并为其他数组形状编写简单的包装器,例如
function nroot(n, x) result(y)
integer, intent(in) :: n
real(dp), intent(in) :: x(:)
real(dp) :: y(size(x))
y = x**(1._dp / n)
end function
function nroot_0d(n, x) result(y)
integer, intent(in) :: n
real(dp), intent(in) :: x
real(dp) :: y
real(dp) :: tmp(1)
tmp = nroot(n, [x])
y = tmp(1)
end function
function nroot_2d(n, x) result(y)
integer, intent(in) :: n
real(dp), intent(in) :: x(:, :)
real(dp) :: y(size(x, 1), size(x, 2))
y = reshape(nroot(n, reshape(x, [size(x)])), [size(x, 1), size(x, 2)])
end function
并按如下方式使用
print *, nroot_0d(2, 9._dp)
print *, nroot(2, [1._dp, 4._dp, 9._dp, 10._dp])
print *, nroot_2d(2, reshape([1._dp, 4._dp, 9._dp, 10._dp], [2, 2]))
这将打印
3.0000000000000000
1.0000000000000000 2.0000000000000000 3.0000000000000000 3.1622776601683795
1.0000000000000000 2.0000000000000000 3.0000000000000000 3.1622776601683795
或者可以使用显式形状数组,如下所示
function nroot(n, k, x) result(y)
integer, intent(in) :: n, k
real(dp), intent(in) :: x(k)
real(dp) :: y(k)
y = x**(1._dp / n)
end function
按如下方式使用
print *, nroot(2, 1, [9._dp])
print *, nroot(2, 4, [1._dp, 4._dp, 9._dp, 10._dp])
print *, nroot(2, 4, reshape([1._dp, 4._dp, 9._dp, 10._dp], [2, 2]))
输出与之前相同
3.0000000000000000
1.0000000000000000 2.0000000000000000 3.0000000000000000 3.1622776601683795
1.0000000000000000 2.0000000000000000 3.0000000000000000 3.1622776601683795