需要将以下 FORTRAN 代码转换为 C++

Question

我是一个非常糟糕的程序员，他们给了我一个程序，据说可以帮助我改进空气动力学硬件。但它是用 Fortran 编写的，我正在尝试使用 MATLAB 来运行这个程序。有什么帮助可以将其转换为 matlab 可以理解的语言吗？ （最好是c++）

      program joukow
c

c   computes joukowski airfoil and finds pressure coefficient

c   currently set up for symmetric airfoil with sharp trailing edge

c   and chord length equal to one.

c   profile is written onto prof.dat and cp onto cp.dat

c      implicit real*8(a-h,o-z)

      complex z,zeta,cw
      dimension uz(100),vz(100),xi(100),eta(100),cp(100)
      dimension xout(100),yout(100)
         open(unit=8,file='prof.dat',status='unknown')
         open(unit=9,file='cp.dat',status='unknown')
      b=1.d0
      write(6,98)
      format(2x,'input the radius of the a-circle in z plane')
      read(5,99)a
      format(f10.0)
      xl=2.*a-1.+1./(2.*a-1.)

c      xl=a+1./a

c      chord=2.*xl

      chord=2.+xl
      del=a-b

c      del =0.1d0
      do 50 i=1,100
      ri=i
      theta=6.2832d0*ri/101.d0
      x=-del+a*cos(theta)
      y=a*sin(theta)
      z=cmplx(x,y)
      zeta=z+b**2/z

c

c  xi and eta are coordinates of points on airfoil

c

      xi(i)=real(zeta)

      eta(i)=aimag(zeta)

      cw=(1.-a**2/(z+del)**2)/(1.-b**2/z**2)
c

c  uz and vz are velocity components on the airfoil assuming the free-stream

c  speed is one.
c
      uz(i)=real(cw)
      vz(i)=-aimag(cw)

c

c  xout and yout are airfoil coordinates where the leading edge is at (0,0)

c  and the chordlength is one.

c

      xout(i)=(xl+xi(i))/chord
      yout(i)=eta(i)/chord
      write(8,100)xout(i),yout(i)
      format(2x,2f10.4)
      continue

c

c  now calculate the pressure coefficient cp

c

      write(6,200)
      format(2x,'pressure coefficients')
      do 70 i=1,50
      cp(i)=1.-(uz(i)**2+vz(i)**2)
      write(9,100)xout(i),cp(i)
      continue
      stop
      end

Answer 1

Matlab 可以很好地理解 Fortran——查看文档。如果这不能满足您的要求，程序中执行任何计算的大部分代码行只需很少的修改即可输入到 Matlab 控制台中。如果你是一个糟糕的程序员，我建议你花时间将程序修改为Matlab而不是C++。如果你没有得到比我现在有时间更好的帮助，我稍后会写更多。

编辑：首先，关于使用 Matlab 中的 Fortran 源文件的一些信息。如果您确实不想（或不能或出于性能原因不这样做）将 Fortran 重写为 Matlab，然后将其转换为 MEX 文件。使用 f2c（或其他任何东西，包括您自己的时间和精力）首先将 Fortran 翻译成 C 或 C++ 对我来说似乎毫无意义。

如果您不喜欢这个想法，这里有一些将 Fortran 转换为 Matlab 的想法。

首先，所有以 C 或 c 开头的行都是注释，因此您不需要翻译它们。从您的代码开始：

  complex z,zeta,cw
  dimension uz(100),vz(100),xi(100),eta(100),cp(100)
  dimension xout(100),yout(100)

这些行声明了许多变量。在 Matlab 中使用变量之前不必声明变量，但有时这样做是有充分理由的。在 Fortran 中也不必这样做，尽管现在普遍认为这是一个坏主意。您可以使用以下语句在 Matlab 中“声明”这些变量：

uz = zeros(100,1); 
vz = zeros(100,1);

通过在 Matlab 中提前声明这些变量，您可以为它们分配一次内存，并避免一些性能降低的问题。

接下来的两行：

     open(unit=8,file='prof.dat',status='unknown')
     open(unit=9,file='cp.dat',status='unknown')

打开几个文件进行输出。它们稍后会在 write 语句中使用 - 忘记它们，改写诸如 save xout 之类的 Matlab 语句。

下一行是 Fortran，但在 Matlab 中是相同的：

  b=1.d0

下一行从控制台获取半径值：

  write(6,98)
  format(2x,'input the radius of the a-circle in z plane')
  read(5,99)a
  format(f10.0)

再次，我建议您忘记这些，只需使用 Matlab 控制台设置 a 的值。更多不需要翻译的 Fortran（尽管我建议您去掉小数点而不跟上 0，或者在它们之间放置一个空格，并且后面的 * -- .* 是 Matlab 中的特定运算符）：

  xl=2.*a-1.+1./(2.*a-1.)

  chord=2.+xl
  del=a-b

Fortran do 循环是与 Matlab for 循环相同。重写：

  do 50 i=1,100

正如

for i = 1:100

其他受访者之一所指出的那样，不清楚匹配的结束语句在哪里，您必须弄清楚这一点。请注意，我只是提供 Fortran 到 Matlab 的逐行翻译。它不是写得很好的 Fortran，而且我也没有提供写得很好的 Matlab，我会把它留给你。

这部分不需要翻译：

  ri=i
  theta=6.2832d0*ri/101.d0 
  x=-del+a*cos(theta)
  y=a*sin(theta)

cmplx 是一个 Fortran 函数，它返回一个复数，该复数具有实部 x 和虚部 y：

  z=cmplx(x,y)

在 Matlab 中，这将是 z = x + y * i。 Fortran 使用 ** 进行求幂，Matlab 使用 ^

  zeta=z+b**2/z

等等。

希望有帮助。

Answer 2

我使用了 f2matlab，然后进行了一些修改。这是清理后可编译的 fortran90 代码：

program joukow
 !
 !   computes joukowski airfoil and finds pressure coefficient
 !   currently set up for symmetric airfoil with sharp trailing edge
 !   and chord length equal to one.
 !   profile is written onto prof.dat and cp onto cp.dat
 !      implicit real*8(a-h,o-z)
 complex z,zeta,cw
 dimension uz(100),vz(100),xi(100),eta(100),cp(100)
 dimension xout(100),yout(100)
 open(unit=8,file='prof.dat',status='unknown')
 open(unit=9,file='cp.dat',status='unknown')
 b=1.d0
 write(6,98)
98 format(2x,'input the radius of the a-circle in z plane')
 read(5,99)a
99 format(f10.0)
 xl=2.*a-1.+1./(2.*a-1.)
!      xl=a+1./a
!      chord=2.*xl
 chord=2.+xl
 del=a-b
!      del =0.1d0
 do i=1,100
  ri=i
  theta=6.2832d0*ri/101.d0
  x=-del+a*cos(theta)
  y=a*sin(theta)
  z=cmplx(x,y)
  zeta=z+b**2/z
  !
  !  xi and eta are coordinates of points on airfoil
  !
  xi(i)=real(zeta)
  eta(i)=aimag(zeta)
  cw=(1.-a**2/(z+del)**2)/(1.-b**2/z**2)
  !
  !  uz and vz are velocity components on the airfoil assuming the free-stream
  !  speed is one.
  !
  uz(i)=real(cw)
  vz(i)=-aimag(cw)
  !
  !  xout and yout are airfoil coordinates where the leading edge is at (0,0)
  !  and the chordlength is one.
  !
  xout(i)=(xl+xi(i))/chord
  yout(i)=eta(i)/chord
  write(8,100)xout(i),yout(i)
100 format(2x,2f10.4)
 end do
!
!  now calculate the pressure coefficient cp
!
 write(6,200)
200 format(2x,'pressure coefficients')
 do  i=1,50
  cp(i)=1.-(uz(i)**2+vz(i)**2)
  write(9,100) xout(i),cp(i)
 end do
 stop
end program joukow

这是生成的 matlab 代码：

function hw1(varargin)
%
%   computes joukowski airfoil and finds pressure coefficient
%   currently set up for symmetric airfoil with sharp trailing edge
%   and chord length equal to one.
%   profile is written onto prof.dat and cp onto cp.dat
%      implicit real*8(a-h,o-z)

format_99=['%10.0f'];
format_100=[repmat(' ',1,2),repmat('%10.4f',1,2),'\n'];
format_200=[repmat(' ',1,2),'pressure coefficients \n'];

fid_8=fopen('prof.dat','w+');
fid_9=fopen('cp.dat','w+');
b=1.0d0;
a=input('input the radius of the a-circle in z plane');
xl=2..*a-1.+1../(2..*a-1.);
%      xl=a+1./a
%      chord=2.*xl
chord=2.+xl;
del=a-b;
%      del =0.1d0
for i=1:100;
 ri=i;
 theta=6.2832d0.*ri./101.0d0;
 x=-del+a.*cos(theta);
 y=a.*sin(theta);
 z=complex(x,y);
 zeta=z+b.^2./z;
 %
 %  xi and eta are coordinates of points on airfoil
 %
 xi(i)=real(zeta);
 eta(i)=imag(zeta);
 cw=(1.-a.^2./(z+del).^2)./(1.-b.^2./z.^2);
 %
 %  uz and vz are velocity components on the airfoil assuming the free-stream
 %  speed is one.
 %
 uz(i)=real(cw);
 vz(i)=-imag(cw);
 %
 %  xout and yout are airfoil coordinates where the leading edge is at (0,0)
 %  and the chordlength is one.
 %
 xout(i)=(xl+xi(i))./chord;
 yout(i)=eta(i)./chord;
 fprintf(fid_8,format_100,xout(i),yout(i));
end; i=100+1;
%
%  now calculate the pressure coefficient cp
%
fprintf(1,format_200);
for  i=1:50;
 cp(i)=1.-(uz(i).^2+vz(i).^2);
 fprintf(fid_9,format_100, xout(i),cp(i));
end;  i=50+1;
end %program joukow

它们都为我提供了相同的结果。不过，我没有检查算法的正确性，只是转换了代码。

Answer 3

我不知道它的支持程度如何 - 但最简单的方法是 f2c，它将 fortran 直接转换为 c 代码。