fortran中怎么定义动态数组？

Posted 2023-05-04

使用ALLOCATABLE和ALLOCATE语句
比如：INTEGER,ALLOCATABLE :: B(:)
...
ALLOCATE(B(N))
ALLOCATE(B(N))给数组B指定了N维的空间。这两个语句一定要配合使用
DEALLOCATE(B)可以解除上述定义的空间。

1.基本使用：
一维数组：
datatype name (size)
datatype: integer, real, complex, logical,或type自定义出来的类型。
integer, parameter :: students = 5
integer :: student (students)
其他定义方法：
integer a(10) !最简单的
integer, dimension (10) :: a !另外的方法
integer, dimension (10) :: a,b,c !一起定义
integer a          !f77 先声明a是整型，
dimension a(10)   ！再声明a是大小为10的数组
或type
type :: person
real ::height, weight
end type
type (person) :: a(10)
a(2)%height = 180.0


二维数组：
integer a(3, 3)
a(1,1) = 3
定义方法：
integer a(10,10)
integer, dimension (10, 10) :: a,b,c
integer a
dimension a(10,10) ! f77, 先定义类型
二维数组常用来当作矩阵使用。

多维数组：
frotran最多可声明高达七维的数组。
integer a(D1,D2,...,Dn) ！n维
一般说来，越高维数组读数据越慢。

另类声明方法：
integer a(5) 可使用的是a(1), a(2), a(3), a(4), a(5)
integer a(0:5) 下标可使用的是0～5，a(0),a(1),...,a(5)六个元素。
integer a(-3：3) 可以使用的是a(-3), a(-2), a(-1), a(0),...,a(3)七个元素。
integer a (5,0:5) 可用a(1~5, 0~5)
integer b(2:3, -1:3) 可用b(2~3, -1~3)

2.数组内容的设置
赋初值：
integer a(5)
data a /1, 2, 3, 4, 5/
data a /5*3/   !5个3    ＊表示数据重复

“隐含式”循环的功能设置初值：
integer a(5)
integer i
data(a(i),i=2,4) /2,3,4/     !a(2)=3,a(3)=3,a(4)=4
输出：write (*,*) (a(i), i=2, 4) 显示a(2),a(3),a(4)
也可以设置累加值：
(a(i),i=2,10,2) ! 循环执行5次，i分别是2,4,6,8,10。

“隐含”选还可以多层嵌套：多维数组
integer a(2,2)
integer i,j
data((a(i,j), i=1,2),j=1,2) /1,2,3,4/
里面的循环先执行，结果：
a(1,1)=1, a(2,1)=2, a(1,2)=3, a(2,2)=4

F90中可以省掉data描述，直接设置初值：
integer :: a(5) = (/1,2,3,4,5/) ! 注意：括号和除号之间不能有空格。
省掉 data直接把初值写在声明之后时，不能像data时那样，用隐式循环对数组中部分元素设初值，每个元素必须给定初值，而且写法有点不同：
integer :: i 
integer :: a(5) = (/ 1, (2, i=2,4), 5/)
!a(1) = 1
!a(2)=2,a(3)=2, a(4)=2   !(2,i=2,4) 是隐式循环
!a(5)=5
f90中的隐式循环：
integer :: a(5)= (/i, i=1,5/)
!a(1)=1, a(2)=2, a(3)=3, a(4)=4, a(5)=5

3.对整个数组的操作(数组广播)
a=5 ! a是任意维数的数组，全部元素设为5
a=(/1,2,3/) ！ 右边数字数目应该和数组的大小一样
a=b ！a,b同样维数及大小的数组，相同位置的元素赋给a
a=b+c！三个都是相同维数大小的数组，相同位置的数字相加给a
a=b-c
a=b*c！注意不等于矩阵的相乘，a(i,j)= b(i,j)*c(i,j)
a=b/c
a=sin(b)   !a(i)=sin(b(i))

4. 对部分数组的操作
a(3:5) = 5
a(3:) = 5
a(3:5) = (/3,4,5/)
a(1:3) = b(4:6)
a(1:10) = a(10:1:-1) !使用隐含循环的方法将数组a(1~10)的内容翻转。
a(:) = b(:, 2)
a(:, :) = b(:, :, 1)
!拿数组的一部分内容使用的时候，需要注意的原测：
（1）等号两边的数组元素数目要一样多。
（2）同时使用多个隐含循环的时候，较低维数循环可以看作是内层的循环。

5. WHERE: 
F95的添加功能，用来取出部分数组的内容进行设置。where命令可以经过逻辑判断来使用数组的一部分。
example:
!把数组a中小于3的元素值设置给b
where(a<3)
b=a              ！这里a,b维数大小相等
end where
也可以写成：
where (a<3) b=a    ！与if相似 
！where命令的代码简单，执行起来速度快。
！where是用来设置数组的，所以它的模块中只能出现与设置数组相关的命令，而且在它的整个程序模块中所使用的数组变量，都必须是同样维数大小的数组。
还可以配合elsewhere来处理逻辑不成立的情况：
where (a<3)
b=1
elsewhere
b=2
end where
还可作多重判断，只要elsewhere后接上逻辑判断就行了：
where (a<2)
b=1
elsewhere(a>5)
b=2
elsewhere
b=3
end where
可以嵌套的，跟循环一样可以取名字，不过在end where的时候一定要加上它的名字，用来明确要结束哪个where模块：
name: where (a<5)
b=a
end where name

where(a<5)
where(a/=2)
b=3
elsewhere
b=1
end where
elsewhere
b=0
end where

6.FORALL: F95添加
integer i
integer :: a(5)

forall (i=1:5)
a(i)=5
end forall

forall(i=1:5)
a(i)=i
end forall

forall的详细语法：
forall (triplet1 [, triplet2 [, triplet3...]], mask)
........
end forall

tripletn是用来赋值数组坐标范围的值。如forall (i=1:5)中i=1:5 就是一个triplet。

integer :: a(10, 5)
forall (i=2:10:2, j=1:5)
a(i,j) = i+j
end forall

mask是用来作条件判断的，跟where命令中使用的条件判断类似，可以用来限制forall程序模块中只作用于数组中符合条件的元素。还可以作其他限制。
forall (i=1:5, j=1:5, a(i, j)<10 )   ！只处理a中小于10的元素
a(i,j) = 1
end forall

forall (i=1:5, j=1:5, i==j)   !只处理i==j的元素
a(i,j) = 1
end forall

forall(i=1:5, j=1:5, ((i>j) .and. a(i,j)>0))
!还可赋值好几个条件，这里只处理二维矩阵的上三角部分且a(i,j)>0的元素
a(i,j)=1/a(i,j)
end forall

如果只有一行代码时候也可以省掉end forall，写在同一行：
forall (i=1:5, j=1:5, a(i,j)/=0) a(i,j)=1/a(i,j)

forall 可以多层嵌套，里面只能出现跟设置数组数值相关的程序命令，还可以在forall中使用where。不过where中不可以使用forall。
forall (i=1:5)
forall (j=1:5)
a(i,j) = 1
end forall
forall (j=6:10)
a(i,j) = 2
end forall
end forall

forall (i=1:5)
where (a(:, i) /=0)
a(:, i) = 1.0/a(:, i)
end where
end forall

7. fortran 中数组的保存规则
不管数组的形状，它所有的元素都是分布在内存中同一个连续的模块当中。
多维数组在内存中的连续模块排列情况是以一种列优先的方法排列的，数组存放在内存中的时候，会先放入第一个列中每个行的元素，然后再放入下一个列中每一行的元素。
exmaple: a(1,1)->a(2,1)->a(3,1)=>a(1,2)->a(2,2)->a(3,2)=>a(1,3)->a(2,3)->a(3,3)
也即对多维数组，会先放入维数的元素，再放入较高维数的元素。
数组元素位置的计算：
一个n维数组a(D1,D2,...,Dn)
设： sn=D1*D2*...*Dn
则a(d1,d2,d3,...,dn)在第1+(d1-1)+(d2-1)*s1+ ... +(dn-1)*sn-1个位置。
使用数组时候最好用低维的下标作为内循环，这样比较快。

8.可变大小的数组：（allocate, allocatable, deallocate）
integer :: students, error = 0
integer ,allocatable :: a(:)   !定义a是可变大小的数组
...
read(*,*) students
allocate(a(students), stat=error) ! or allocate(a(students))   ！申请数组内存空间
...
deallocate(a)    ！释放动态数组占用的内存空间
多维数组：
integer, allocatable :: a2(:, :)
integer, allocatable :: a3(:, :, :)
allocate(a2(5, 5))
allocate(a3(5, 5, 5))
或
allocate(a2(-3:3, -3:3))
相关的函数allocated可用来检查一个可变大小的数组是否已经配置内存来使用，会返回一个逻辑值。
example:
if(.not. allocated(a)) then
allocate(a(5))
end if

参考技术A 使用ALLOCATABLE和ALLOCATE语句
比如：INTEGER,ALLOCATABLE :: B(:)
...
ALLOCATE(B(N))
ALLOCATE(B(N))给数组B指定了N维的空间。这两个语句一定要配合使用
DEALLOCATE(B)可以解除上述定义的空间。

试图将一个连续的动态二维数组从 C 传递到 Fortran

【中文标题】试图将一个连续的动态二维数组从 C 传递到 Fortran

【英文标题】：Attempting to pass a contiguous dynamic 2d array from C to Fortran

【发布时间】：2016-12-09 11:58:32

【问题描述】：

我知道extra care needs to be taken 在 C 中分配内存以确保二维数组是连续的时，但是当我将它传递给 Fortran 时仍然没有得到预期的结果。下面是我尝试的玩具版本：为 2d 数组分配内存并为每个元素分配值的 main.c 文件，以及打印出 2d 数组元素的 foo.f90 文件。

#include <stdio.h>
#include <stdlib.h>

void foo_(double **,int *,int *);

int main() 
    int i, j, cols, rows;
    double *_x, **x;

    cols = 3;
    rows = 2;

//  Allocate memory
    _x = malloc(rows*cols*sizeof(double));
    x = malloc(cols*sizeof(double *));
    for(i=0;i<cols;i++)
        x[i] = &(_x[rows*i]);

//  Generate elements for the 2d array 
    for(i=0;i<cols;i++)
    for(j=0;j<rows;j++)
        x[i][j] = i*j;

//  Call Fortran subroutine foo
    foo_(x,&rows,&cols);

    return EXIT_SUCCESS;

foo.h

subroutine foo(x,rows,cols)
    use iso_c_binding
    implicit none

    integer(c_long), intent(in)                :: rows,cols
    real(c_double), intent(in), dimension(rows,cols) :: x

    integer :: i,j

    do i = 1,cols
    do j = 1,rows
        print *, j,i,x(j,i)
    end do
    end do

end subroutine

作为输出，我需要一个数组元素列表。相反，我得到以下输出

           1           1   1.1654415706619996E-316
           2           1   1.1654423611670330E-316
Segmentation fault (core dumped)

【问题讨论】：

您的 C 代码中有 no 二维数组，没有任何东西可以表示或指向一个！指针不是数组。如果您需要二维数组，请使用一个！

双指针可以表示二维数组，实际上有@Olaf，我认为问题出在 foo_(double **,int *,int *);不应该只是 foo_(double **,int,int); ?

@Olaf 你的评论要么是错误的，要么是愚蠢的，以至于在教学上毫无用处。在 C 编程语言中，存储在堆中的多维数组在运行时保留其内存，指针用作管理地址的簿记工具（cs.swarthmore.edu/~newhall/unixhelp/C_arrays.html）。

@francescalus 好点。我已编辑问题以包含输出。

@ThomasKelly：1) C 标准不强制使用特定的内存管理模型进行动态内存分配。 2）你的断言是完全错误的。锯齿状数组不是使用二维数组的正确方法。这是一个完全不同的数据结构！您可以在 C 中轻松分配 2D 数组。您可能想要更新到 17 岁的 C99，或者 - 更好 - 到 C11（也已经 5 岁）的 C 标准。告诉你使用 double ** 的人正在教古老的 C。3）这可能是你的问题的原因之一。

【参考方案1】：

你自己也是这么说的：连续！

您没有分配一个连续的数组。要分配一个，你必须写：

//  Allocate memory
double *x = malloc(rows*cols*sizeof(double));

不幸的是，您现在必须用 C 编写以下代码来索引 x：

//  Generate elements for the 2d array 
for(i=0;i<cols;i++)
    for(j=0;j<rows;j++)
        *(x+j*cols+i) = i*j;

这假设矩阵是按行优先的顺序（在内存中连续排列的一行接一行）。

注意：在 C99 中有可变长度数组，编译器在其中正确管理 x[i][j]。我没有 C99，但也许其他用户可以用 VLA 给出答案。

【讨论】：

这是类型错误。 x 是 double ** 并且您的索引是错误的。为什么不使用二维数组呢？这完全取决于 Fortran ABI 的预期。

问题标记为 C，即标准 C，即 C11。虽然 C11 使 VLA 成为可选项，但没有任何现代 C11 兼容编译器不支持它们（使现有语言功能成为可选项是委员会维护向下兼容性实践的重大突破，可能让过时编译器的维护者更好地出售它们）。跨度>

【参考方案2】：

Fortran 端的foo() 接收x 作为显式形状数组，因此我们需要传递x 的第一个元素的地址。所以，我们将原型修改为

void foo_( double *, int *, int * );

并将第一个元素的地址传递为

foo_( _x, &rows, &cols );
// or
// foo_( &( _x[0] ), &rows, &cols );
// or
// foo_( &( x[0][0] ), &rows, &cols );

另外，我们可能需要使用integer(c_int) 而不是integer(c_long) 来匹配C 端的int，这样

integer(c_int), intent(in) :: rows,cols

（在我的计算机上，使用 integer(c_long) 会导致分段错误，因为 rows 和 cols 未正确传递。）

另外，为了方便检查数组元素的对应关系，我将测试值修改为

for(i=0;i<cols;i++)
for(j=0;j<rows;j++)
    x[i][j] = (i+1) + 1000*(j+1);

并将打印语句插入到 Fortran 代码中

print *, "rows = ", rows
print *, "cols = ", cols

然后，我的计算机上的 GCC-6 (OSX 10.9) 给出了

 rows =            2
 cols =            3
           1           1   1001.0000000000000     
           2           1   2001.0000000000000     
           1           2   1002.0000000000000     
           2           2   2002.0000000000000     
           1           3   1003.0000000000000     
           2           3   2003.0000000000000

作为旁注，以下似乎也有效（而不是手动创建x）：

_x = malloc( rows * cols * sizeof(double) );

typedef double (*Array)[ rows ];
Array x = (Array) _x;

// set x[i][j] the same way

foo_( _x, &rows, &cols );
// or
// foo_( &( x[0][0] ), &rows, &cols );

但我不太确定 Array 的这种用法是否符合标准...（在 C 中）。

---

[编辑]

使用现代 C，似乎可以将 x 直接声明为具有动态大小的矩形数组，并在堆上分配内存，这样：

double (* x)[rows] = malloc( cols * rows * sizeof(double) );

// or
// double (* x)[rows] = malloc( sizeof( double[cols][rows] ) );

（请注意，变量名称“cols”和“rows”指的是 Fortran 端的变量名，因此它们在 C 端显得违反直觉。）

有了这个x，我们可以像上面一样继续调用foo_()，即，

// set x[i][j] via loops

foo_( &( x[0][0] ), &rows, &cols );

请参阅Jens Gustedts' blog（“不要使用假矩阵”）了解后一种方法的完整详细信息（感谢@Olaf 的建议）。

【讨论】：

这可能是真的 Fortran ABI 使用一维数组来连接 C 代码，但应该首先验证这一点。毕竟，这仍然与 OP 使用的锯齿状数组非常不同，也不是 2D 数组。由于它具有相同的（考虑正确的行/列顺序）布局，它应该可以工作，但它仍然具有不同的类型并且可能导致细微的错误，resp。并发症，需要其他不必要的演员..

OP 使用指向指针的指针，resp。一个指针数组。通常称为“锯齿状数组”，无论二级数组的大小是否相同。他在块中分配子数组的事实并没有改变这一点。仅 C 语言的正确方法是使用 2D 数组，但正如我所写，我不确定 Fortran 期望什么。

@Olaf 第一级数组 (_x) 是线性且连续的，对吗？然后，我的理解是x（=二级数组）提供了对_x指向的内存的矩形访问。传递给 Fortran 的是 malloc() 分配的地址，它是线性的（我假设）。另外，我的理解是“锯齿状数组”是每一行或每一列都是独立分配的，所以内存是不连续的。

但是因为 OP 使用变量名“rows”和“cols”有点令人困惑（为了匹配 Fortran 主要列并避免转置），所以乍一看代码对我来说不是很可读（对不起，只是我的意见......）

是的。 Jens 是这里的成员之一，通常可以信任 C 标准问题。也不知道他写博客。虽然 C99 现在已经 17 岁了，但这篇文章很有道理。