1

f2py を使用して、3 次元で単純な統合問題を実行しようとしています。

fortran コードを呼び出すpython コードは次のとおりです。

#!/Library/Frameworks/EPD64.framework/Versions/Current/bin/python
import pymods as modules
import pygauleg as gauleg
import pyint as integrator
import pylab as pl
import sys
import math
import time

############################################
# main routine #############################
############################################
zero = 0.0
one = 1.0
pi = pl.pi

Nr = 10
Nt = 10
Np = 2*Nt
r0 = zero
rf = one
NNang = Nr*Nt*Np

print 'Nr Nt Np = ', Nr, Nt, Np
print 'NNang = ', NNang
print 'r0 rf = ', r0, rf

Nx = int(math.floor( (one*NNang)**(one/3.0) ))
Ny = int(math.floor( (one*NNang)**(one/3.0) ))
Nz = int(math.floor( (one*NNang)**(one/3.0) ))

Nx = int(pl.floor(float(Nx)*1.75))
Ny = int(pl.floor(float(Ny)*1.75))
Nz = int(pl.floor(float(Nz)*1.75))

NNxyz = Nx*Ny*Nz


print 'Nx Ny Nz = ', Nx, Ny, Nz
print 'NNxyz = ', NNxyz
xyz0 = -rf
xyzf = rf

t1 = time.time()
xt = pl.zeros(Nt)
wt = pl.zeros(Nt)
gauleg.gauleg(xt, wt, 0.0, pl.pi, Nt)
print 'outside of gauleg'

fortran サブルーチンは少し長いですが、その重要な部分は最初の部分です...

  2 subroutine gauleg(x,w,x1,x2,n)
  3 !Input:   x1,x2,n
  4 !Output:  x,w
  5 !implicit none
  6 !integer, parameter :: ikind = selected_int_kind(25)
  7 !integer, parameter :: rkind = selected_real_kind(15, 307)
  8 !
  9 !real(kind = rkind), parameter :: pi = 3.14159265358979323846d00
 10 !real(kind = rkind), parameter :: one = 1.0d00
 11 !real(kind = rkind), parameter :: zero = 0.0d00
 12 use mod_gvars
 13 
 14 real(kind = rkind) :: tol = 1d-15
 15 
 17 integer :: n
 18 !!!!!f2py intent(in) n
 19 real(kind = rkind), dimension(n) :: x
 20 real(kind = rkind), dimension(n) :: w
 22 real :: x1, x2
 23 
 24 real(kind = rkind) :: z1, z, xm, xl, pp, p3, p2, p1;
 25 
 26 integer(kind = ikind) :: m
 27 integer(kind = ikind) :: i,j
 28 integer(kind = ikind) :: countmax, counter, max_counter, min_counter
 29 
 30 integer(kind = ikind) :: tenth, hundredth, thousandth
 31 
 32 print*, 'n = ', n

そして終わり...

 98 
 99 print*, 'returning'
100 
101 end subroutine

サブルーチンの先頭 (5 ~ 11 行目) のコメントは、fortran モジュールに存在する構造体ですmod_gvars。このサブルーチンが戻るまで*すべてが計画通りに進んでいるようです。出力は次のとおりです。

Nr Nt Np =  10 10 20
NNang =  2000
r0 rf =  0.0 1.0
Nx Ny Nz =  21 21 21
NNxyz =  1728
 n =           10
 m =  5
 returning
python(14167) malloc: *** error for object 0x1081f77a8: incorrect checksum for freed object - object was probably modified after being freed.
*** set a breakpoint in malloc_error_break to debug
Abort trap

サブルーチンが戻ったときにのみ問題が発生するようです。なぜこれが起こるのでしょうか?

4

1 に答える 1

3

この種の問題は通常、Python の参照カウント エラーで発生します。たとえば、f2py にバグがある場合や、Fortran で numpy 配列に割り当てられているよりも多くのメモリを上書きした場合に発生する可能性があります。これらのエラーはすべて、Fortran サブルーチンが終了した後、ランダムな時点でのみ表示されます。通常は、Python で一部のメモリの割り当てを解除したときに表示されます。

これをデバッグするには、Fortran で取得したすべての配列を出力してみてください。つまり、配列 x、w を出力して、そこにあるすべてのメモリにアクセスできることを確認してください (f2py が同じ型を使用していることなどをテストします)。 .

Fortran で境界チェックを使用していることを確認してください (少なくとも-fbounds-checkgfortran では、できれば-fcheck=allすべての問題をチェックするためだけに)。

これをデバッガーまたは valgrind で実行することもできます。問題の場所がわかる場合があります。

最後に、個人的には、Cython を使用して Fortran を直接ラップすることを好みます。その後、生成されたすべてのファイルに簡単にアクセスでき、iso_c_bindingFortran モジュールを使用して、Fortran コンパイラーが Fortran と C (Python) の間ですべての型に互換性があることを確認します。例については、こちらを参照してください。

于 2012-09-07T05:21:21.317 に答える