c++ - 如何从 C++ 调用 fortran 例程？

Question

我想从我的 C++ 代码中调用 fortran 例程 cbesj.f，我该如何实现？

以下是我完成的步骤:

1. 从 netlib amos 网页下载 cbesj.f 和依赖项，http://www.netlib.org/cgi-bin/netlibfiles.pl?filename=/amos/cbesj.f

2. 在源目录中，

   f2c -C++PR *.f

   g++ -c *.c

   ar cr libmydemo.a *.o

3. [test_cbesj.cpp][1] 和[mydemo.h][2]就是这样调用子程序的，

   g++ test_cbesj.cpp -lf2c -lm -L。 -lmydemo它返回错误:

test_cbesj.cpp:(.text+0xd6): undefined reference to `cbesj_(complex*, float*, long*, long*, complex*, long*, long*)'

在我的问题中引用 cbesj_ 子例程的正确方法是什么？谢谢!

感谢凯西:我认为你的方法是最好的。但是我还是有set fault，为什么？我们开始吧:

f77 -c *.f    

在 modemo.h 中

//File mydemo.h                                                                                                                              
#ifndef MYDEMO_H
#define MYDEMO_H
#include <stdio.h>      /* Standard Library of Input and Output */
#include "f2c.h"

extern"C" int cacai_(complex *z__, real *fnu, integer *kode, integer *mr,         integer *n, complex *y, integer *nz, real *rl, real *tol, real *el\
im, real *alim);
extern"C" int cairy_(complex *z__, integer *id, integer *kode, complex *ai,         integer *nz, integer *ierr);
extern"C" int casyi_(complex *z__, real *fnu, integer *kode, integer *n,     complex *y, integer *nz, real *rl, real *tol, real *elim, real     \
  *alim);
extern"C" int cbesj_(complex *z__, real *fnu, integer *kode, integer *n,     complex *cy, integer *nz, integer *ierr);
extern"C" int cbinu_(complex *z__, real *fnu, integer *kode, integer *n,     complex *cy, integer *nz, real *rl, real *fnul, real *tol, real *el\
im, real *alim);
extern"C" int cbknu_(complex *z__, real *fnu, integer *kode, integer *n, complex *y, integer *nz, real *tol, real *elim, real *alim);
extern"C" int cbuni_(complex *z__, real *fnu, integer *kode, integer *n,     complex *y, integer *nz, integer *nui, integer *nlast, real *fnul, \
real *tol, real *elim, real *alim);
extern"C" int ckscl_(complex *zr, real *fnu, integer *n, complex *y, integer *nz, complex *rz, real *ascle, real *tol, real *elim);
extern"C" int cmlri_(complex *z__, real *fnu, integer *kode, integer *n,     complex *y, integer *nz, real *tol);
extern"C" int crati_(complex *z__, real *fnu, integer *n, complex *cy, real *tol);
extern"C" int cs1s2_(complex *zr, complex *s1, complex *s2, integer *nz, real *ascle, real *alim, integer *iuf);
extern"C" int cseri_(complex *z__, real *fnu, integer *kode, integer *n, complex *y, integer *nz, real *tol, real *elim, real *alim);
extern"C" int cshch_(complex *z__, complex *csh, complex *cch);
extern"C" int cuchk_(complex *y, integer *nz, real *ascle, real *tol);
extern"C" int cunhj_(complex *z__, real *fnu, integer *ipmtr, real *tol, complex *phi, complex *arg, complex *zeta1, complex *zeta2, complex\
 *asum, complex *bsum);
extern"C" int cuni1_(complex *z__, real *fnu, integer *kode, integer *n, complex *y, integer *nz, integer *nlast, real *fnul, real *tol     \
  , real *elim, real *alim);
extern"C" int cuni2_(complex *z__, real *fnu, integer *kode, integer *n, complex *y, integer *nz, integer *nlast, real *fnul, real *tol     \
  , real *elim, real *alim);
extern"C" int cunik_(complex *zr, real *fnu, integer *ikflg, integer *ipmtr, real *tol, integer *init, complex *phi, complex *zeta1, complex\
 *zeta2, complex *sum, complex *cwrk);
extern"C" int cuoik_(complex *z__, real *fnu, integer *kode, integer *ikflg,     integer *n, complex *y, integer *nuf, real *tol, real *elim, re\
al *alim);
extern"C" int cwrsk_(complex *zr, real *fnu, integer *kode, integer *n,     complex *y, integer *nz, complex *cw, real *tol, real *elim, real *a\
lim);
extern"C" real gamln_(real *z__, integer *ierr);
extern"C" integer i1mach_(integer *i__);
extern"C" real r1mach_(integer *i__);
extern"C" int xerror_(char *mess, integer *nmess, integer *l1, integer *l2,     ftnlen mess_len);
    #endif

在test_cbesj.cpp中，

#include "mydemo.h"
#include "f2c.h"
#include <math.h>
#include <iostream>
#include <stdio.h>
#include <assert.h>
using namespace std;
int main(void)
{
  //  double x=86840.;                                                                                                                       
  //int nu=46431,j, err;                                                                                                                     
  complex *z,z__;
  z__.r = 3.0;z__.i = 2.0;z = &z__;
  cout << z->r << '\t' << z->i << endl;
  real *fnu;float fnu__ = 3.0;fnu = &fnu__;
  integer *kode ;long int kode__=1;kode=&kode__;
  integer *n    ;long int n__=1;n=&n__;
  complex *cy;
  integer *nz;
  integer *ierr;
  cbesj_(z, fnu, kode, n, cy, nz, ierr);
  cout << cy->r << '\t' << cy->i << endl;

  return 0;
}

然后，

g++ -c -g test_cbesj.cpp
g++ -o test *.o -lg2c
./test
3   2
Segmentation fault (core dumped)
gdb test 
GNU gdb (Ubuntu/Linaro 7.4-2012.04-0ubuntu2.1) 7.4-2012.04
Copyright (C) 2012 Free Software Foundation, Inc.
License GPLv3+: GNU GPL version 3 or later <http://gnu.org/licenses/gpl.html>
This is free software: you are free to change and redistribute it.
There is NO WARRANTY, to the extent permitted by law.  Type "show copying"
and "show warranty" for details.
This GDB was configured as "i686-linux-gnu".
For bug reporting instructions, please see:
<http://bugs.launchpad.net/gdb-linaro/>...
Reading symbols from /media/Downloads/amos-4/test...done.
(gdb) run
Starting program: /media/Downloads/amos-4/test 
3   2

Program received signal SIGSEGV, Segmentation fault.
0x0804b355 in cbesj_ ()
(gdb) frame 0
#0  0x0804b355 in cbesj_ ()
(gdb) frame 1
#1  0x0805a3ca in main () at test_cbesj.cpp:21
21    cbesj_(z, fnu, kode, n, cy, nz, ierr);

感谢roygvib的回复!实际上很好的建议。这是更改后的 test_cbesj.cpp:

complex z, cy;
  float fnu;
  long int kode, n, nz, ierr;

  z.r = 3.0; z.i = 2.0;
  fnu = 3.0;
  n = 1; kode = 1;
  cout.precision(16);
  cbesj_( &z, &fnu, &kode, &n, &cy, &nz, &ierr );
  cout << cy.r << '\t' << cy.i << endl;
  cout << "nz=" << nz << endl;
  cout << "ierr=" << ierr << Lendl;

不再有段错误。但出于某些原因，代码没有按预期工作:

./test
-1.343533039093018  -1.343533992767334
nz=0
ierr=4

答案是错误的，ierr 也从源代码中这样说:

C           NZ     - NUMBER OF COMPONENTS SET TO ZERO DUE TO UNDERFLOW,
C                    NZ= 0   , NORMAL RETURN
C                    NZ.GT.0 , LAST NZ COMPONENTS OF CY SET TO ZERO
C                              DUE TO UNDERFLOW, CY(I)=CMPLX(0.0,0.0),
C                              I = N-NZ+1,...,N
C           IERR   - ERROR FLAG
C                    IERR=0, NORMAL RETURN - COMPUTATION COMPLETED
C                    IERR=1, INPUT ERROR   - NO COMPUTATION
C                    IERR=2, OVERFLOW      - NO COMPUTATION, AIMAG(Z)
C                            TOO LARGE ON KODE=1
C                    IERR=3, CABS(Z) OR FNU+N-1 LARGE - COMPUTATION DONE
C                            BUT LOSSES OF SIGNIFCANCE BY ARGUMENT
C                            REDUCTION PRODUCE LESS THAN HALF OF MACHINE
C                            ACCURACY
C                    IERR=4, CABS(Z) OR FNU+N-1 TOO LARGE - NO COMPUTA-
C                            TION BECAUSE OF COMPLETE LOSSES OF SIGNIFI-
C                            CANCE BY ARGUMENT REDUCTION
C                    IERR=5, ERROR              - NO COMPUTATION,
C                            ALGORITHM TERMINATION CONDITION NOT MET

Answer 1

我从netlib 下载了cbesj(或zbesj)和相关文件。 ，但出于某种原因，它们不适用于 gfortran4.8(在 Linux x86_64 上)。更准确地说，cbesj 总是给出 ierr=4，而我无法编译 zbesj 因为 zabs 参数太多(但请参阅下面的更新)。

---

所以我放弃了使用原始AMOS代码并尝试了openspecfun那也是基于AMOS。我通过键入 make(使用 gfortran)编译了后者，生成了 libopenspecfun.a 等。然后我编写了以下代码来测试 zbesj:

#include <cstdio>

extern "C" {
    void zbesj_( double *zr, double *zi, double *fnu, int *kode, int *n,
                 double *Jr, double *Ji, int *nz, int *ierr );
}

int main()
{
    int    n, nz, ierr, kode;
    double fnu, zr, zi, Jr, Ji;

    fnu = 2.5;
    zr = 1.0; zi = 2.0;
    n = 1; kode = 1;

    zbesj_( &zr, &zi, &fnu, &kode, &n, &Jr, &Ji, &nz, &ierr );

    printf( "fnu = %20.15f\n", fnu );
    printf( "z   = %20.15f %20.15f\n", zr, zi );
    printf( "J   = %20.15f %20.15f\n", Jr, Ji );
    printf( "nz, ierr = %d %d\n", nz, ierr );

    return 0;
}

并编译为

g++ test_zbesj.cpp libopenspecfun.a -lgfortran

给出

fnu =    2.500000000000000
z   =    1.000000000000000    2.000000000000000
J   =   -0.394517225893988    0.297628229902939
nz, ierr = 0 0

因为 this site说答案是 -0.394517...+ 0.297628...i，zbesj 的结果似乎是正确的。

---

[更新]

通过更仔细地阅读原始代码并与 openspecfun 进行比较，发现用户需要取消注释 i1mach.f 和 r1mach 的适当部分.f(或 d1mach.f)取决于机器环境。对于 Linux x86_64，取消注释标记为

的部分似乎就足够了

C     MACHINE CONSTANTS FOR THE IBM PC

通过此修改，cbesj 可以在 ierr=0 下正常工作；否则它给出 ierr=4 因为一些常量默认为 0。对于 double 版本，我们还需要处理用户定义的ZABS，其定义与固有定义冲突。 Intel Fortran (ifort) 识别用户定义的 ZABS 并成功编译它们(尽管有很多警告)，而 gfortran 停止编译并出现语法错误(假设它是内在的)。 openspecfunc 通过用固有的重写所有 ZABS 避免了这个问题。无论如何，通过上述修改，cbesj 和 zbesj 都能正常工作(如预期)。

---

[更新 2]

事实证明，使用 d1mach.f、r1mach.f 和 i1mach.f 的 BLAS 版本，事情变得更加简单因为它们会自动确定与机器相关的常量，我们不需要手动修改代码。请参阅Netlib/FAQ详情页。同样的技巧也适用于 SLATEC(例如，这个 page )。

wget http://www.netlib.org/blas/i1mach.f 
wget http://www.netlib.org/blas/r1mach.f
wget http://www.netlib.org/blas/d1mach.f