标签: fortran

我不时会读到Fortran是或者可以比C更快进行繁重的计算.这是真的吗？我必须承认我几乎不知道Fortran,但到目前为止我见过的Fortran代码并没有表明该语言具有C所没有的功能.

如果是真的,请告诉我原因.请不要告诉我哪些语言或库对数字运算有好处,我不打算写一个app或lib来做那个,我只是很好奇.

出于好奇,我决定将我自己的矩阵乘法函数与BLAS实现进行比较......我对结果的评价最少:

自定义实现,1000x1000矩阵乘法的10次试验:

Took: 15.76542 seconds.

Run Code Online (Sandbox Code Playgroud)

BLAS实施,1000x1000矩阵乘法的10次试验:

Took: 1.32432 seconds.

Run Code Online (Sandbox Code Playgroud)

这是使用单精度浮点数.

我的实施:

template<class ValT>
void mmult(const ValT* A, int ADim1, int ADim2, const ValT* B, int BDim1, int BDim2, ValT* C)
{
    if ( ADim2!=BDim1 )
        throw std::runtime_error("Error sizes off");

    memset((void*)C,0,sizeof(ValT)*ADim1*BDim2);
    int cc2,cc1,cr1;
    for ( cc2=0 ; cc2<BDim2 ; ++cc2 )
        for ( cc1=0 ; cc1<ADim2 ; ++cc1 )
            for ( cr1=0 ; cr1<ADim1 ; ++cr1 )
                C[cc2*ADim2+cr1] += A[cc1*ADim1+cr1]*B[cc2*BDim1+cc1];
}

我有两个问题:

1. 鉴于矩阵 - 矩阵乘法说:nxm*mxn需要n*n*m次乘法,因此在1000 ^ 3或1e9次操作的情况下.我的2.6Ghz BLAS处理器如何在1.32秒内完成10*1e9操作？即使多重操作是单个操作而没有其他任何操作,也需要约4秒.
2. 为什么我的实现速度要慢得多？

我发现使用Python读取二进制文件特别困难.你能帮我个忙吗？我需要阅读这个文件,它在Fortran 90中很容易阅读

int*4 n_particles, n_groups
real*4 group_id(n_particles)
read (*) n_particles, n_groups
read (*) (group_id(j),j=1,n_particles)

具体来说,文件格式为:

Bytes 1-4 -- The integer 8.
Bytes 5-8 -- The number of particles, N.
Bytes 9-12 -- The number of groups.
Bytes 13-16 -- The integer 8.
Bytes 17-20 -- The integer 4*N.
Next many bytes -- The group ID numbers for all the particles.
Last 4 bytes -- The integer 4*N. 

我怎么用Python阅读？我尝试了一切,但从未奏效.我有没有机会在python中使用f90程序,读取这个二进制文件,然后保存我需要使用的数据？

我得到一个512 ^ 3数组,表示模拟的温度分布(用Fortran编写).该数组存储在大小约为1/2G的二进制文件中.我需要知道这个数组的最小值,最大值和平均值,因为我很快就需要了解Fortran代码,我决定试一试,并提出了以下非常简单的例程.

  integer gridsize,unit,j
  real mini,maxi
  double precision mean

  gridsize=512
  unit=40
  open(unit=unit,file='T.out',status='old',access='stream',&
       form='unformatted',action='read')
  read(unit=unit) tmp
  mini=tmp
  maxi=tmp
  mean=tmp
  do j=2,gridsize**3
      read(unit=unit) tmp
      if(tmp>maxi)then
          maxi=tmp
      elseif(tmp<mini)then
          mini=tmp
      end if
      mean=mean+tmp
  end do
  mean=mean/gridsize**3
  close(unit=unit)

在我使用的机器上,每个文件大约需要25秒.这让我觉得很长,所以我继续在Python中做了以下事情:

    import numpy

    mmap=numpy.memmap('T.out',dtype='float32',mode='r',offset=4,\
                                  shape=(512,512,512),order='F')
    mini=numpy.amin(mmap)
    maxi=numpy.amax(mmap)
    mean=numpy.mean(mmap)

现在,我预计这会更快,但我真的被吹走了.在相同条件下只需不到一秒钟.平均值偏离我的Fortran例程发现的那个(我也使用128位浮点运行,所以我不知何故更多地信任它),但仅限于第7位有效数字左右.

numpy怎么这么快？我的意思是你必须查看数组的每个条目才能找到这些值,对吧？我在Fortran程序中做了一件非常愚蠢的事情,因为它花了这么长时间吗？

编辑:

要回答评论中的问题:

• 是的,我也使用32位和64位浮点运行Fortran例程,但它对性能没有影响.
• 我用过iso_fortran_env它提供128位浮点数.
• 使用32位浮点数我的意思是相当多,所以精度确实是一个问题.
• 我以不同的顺序在不同的文件上运行这两个例程,所以缓存在比较中应该是公平的吗？
• 我实际上试过打开MP,但同时从不同位置的文件中读取.阅读完你的评论和答案后,这听起来真的很愚蠢,它使常规也需要更长的时间.我可能试一试数组操作,但也许甚至不需要.
• 文件实际上是1/2G大小,这是一个错字,谢谢.
• 我现在将尝试数组实现.

编辑2:

我实现了@Alexander Vogt和@casey在他们的答案中提出的建议,它的速度一样快,numpy但现在我有一个精确的问题,因为@Luaan指出我可能会得到.使用32位浮点数组,计算出的平均值sum为20%.干

...
real,allocatable :: tmp (:,:,:)
double precision,allocatable :: tmp2(:,:,:)
...
tmp2=tmp
mean=sum(tmp2)/size(tmp)
...

解决了这个问题但增加了计算时间(不是很多,但显着).有没有更好的方法来解决这个问题？我找不到从文件中直接读单打的方法.怎么numpy避免这个？

感谢目前为止所有的帮助.

我最近开始维护大量科学计算密集型FORTRAN代码.尽管有谷歌和两本入门级书籍,但我很难掌握四十年前语言的所有细微差别.代码充斥着"提升性能".有没有人有任何关于将FORTRAN 去优化到CS 101级别的指南或实用建议？有谁知道FORTRAN代码优化是如何运作的？是否有任何典型的FORTRAN'陷阱'可能不会发生在Java/C++/.NET引发的开发人员接管FORTRAN 77/90代码库？

在一些科学环境中,你经常离不开FORTRAN,因为大多数开发人员只知道这个习语,并且有很多遗留代码和相关经验.坦率地说,高性能编程没有很多其他跨平台选项(C++可以完成任务,但语法,零启动数组和指针与某些人不兼容).

所以,让我们假设一个新项目,必须使用Fortran 90的,但我想建立最先进的软件架构出来的,同时与最新的编译器兼容(英特尔ifort,而且还包括Sun/HP/IBM编译)

所以我想要强加一些众所周知的常识,但在我的环境中还不是标准:

• 全局变量禁止,没有getos,没有跳转标签implicit none等
• "面向对象编程"(具有数据类型和相关子例程的模块)
• 模块化/可重用的功能,文档齐全,可重用的库
• 断言/前提条件/不变量(使用预处理程序语句实现)
• 所有(大多数)子程序和"对象"的单元测试
• 一个强烈的"调试模式"(#ifdef DEBUG),带有更多检查和所有可能的英特尔编译器检查(数组边界,子程序接口等)
• 统一和强制清晰的编码风格,使用代码处理工具助手.

所有这一切的目标是拥有值得信赖,可维护和模块化的代码.然而,在许多遗留代码中,可重用性不是一个重要的目标.

我四处搜寻有关面向对象的Fortran编程契约式引用(断言/前提条件的/ etc.),发现只有丑陋的和过时的文件,通过人与没有大规模参与项目进行语法和文章,和死了的项目.

有关此主题的任何好的URL,建议,参考文件/书籍？

我读过Fortran仍然大量用于科学计算.对于已经大量投资Fortran的代码,这对我来说很有意义.

但是有理由将Fortran用于其他现代语言的新项目吗？与更流行的语言(C++,Java,Python,Ruby等)相比,Fortran中的语言设计决策是否更适合科学计算？例如,与我提到的其他语言相比,Fortran的特定语言功能是否允许编译器中的数值优化程度更高？

定义PI的动机是什么？

PI=4.D0*DATAN(1.D0)

在Fortran 77代码中？我理解它是如何工作的,但是,理由是什么？

我正在尝试学习Fortran,我看到很多不同的定义被传递,我想知道他们是否想要完成同样的事情.以下有什么区别？

• integer*4
• integer(4)
• integer(kind=4)

我正在尝试为我的领域中的一些流行的Fortran代码进行python包分发.我希望它对setup.py文件使用最标准的方法.相关的qustion有助于学习如何包装Fortran扩展.

当使用这种方法,我混合时发现一些容易混淆的行为setuptools和numpy.distutils.将两者混合是不好的做法吗？截至2015年,似乎最好setuptools尽可能多地使用.

但是,我想以Fortran一种兼容的方式构建扩展numpy. 所以我想从中导入numpy.distutils以获取Extension和setup.

我使用以下基本方法:

from setuptools.command.develop import develop
from numpy.distutils.core import Extension, setup

ext_modules=[Extension("my_package.fortran_mod", sources=['src/fortran_mod.f'])]

class MyDevelop(develop):
  def run(self):
    my_script()
    develop.run(self)

setup(
  ...
  ext_modules=ext_modules,
  cmdclass={'develop':MyDevelop})

这似乎有效,但我有疑问.

1. 混合setuptools和一般是好的做法numpy.distribute吗？
2. 我导入它们的顺序是否重要？我应该setuptools先进口吗？
3. 是否有正式的包装扩展教程numpy？也许甚至有一些讨论Fortran扩展？

一些链接

https://www.youtube.com/watch?v=R4yB-8tB0J0

http://www.fortran90.org/src/best-practices.html#interfacing-with-python