【发布时间】:2020-11-03 21:09:25
【问题描述】:
我正在尝试将 python 函数传递给 f2py 生成的函数。当 python 函数的输入和输出是一个数组时,我无法弄清楚如何做到这一点。 (我是 f2py 的新手!)
这是我的 fortran 代码 (num_Dg.f):
subroutine num_Dg(func,U,jacob,n)
double precision, intent(in) :: U(n)
double precision, intent(out) :: jacob(n,n)
double precision Up(n), Um(n)
double precision g_valp(n),g_valm(n)
double precision :: rel_eps = 2.e-5
double precision eps
integer i, j
external func
eps = minval(U)*rel_eps
do j=1,n
do i=1,n
if (i .eq. j) then
Up(i) = U(i)+(eps)
Um(i) = U(i)-(eps)
else
Up(i) = 0
Um(i) = 0
endif
enddo
call func(Up,g_valp,n)
call func(Um,g_valm,n)
do i=1,n
jacob(i,j) = 0.0
jacob(i,j) = (g_valp(i)-g_valm(i))/(2*eps)
enddo
enddo
end subroutine
我使用 f2py 编译:f2py -c num_Dg.f -m num_Dg。它编译没有错误。但是当将函数导入 python 时,我得到了这个:
>>> from num_Dg import num_dg
>>> print(num_dg.__doc__)
jacob = num_dg(func,u,[n,func_extra_args])
Wrapper for ``num_dg``.
Parameters
----------
func : call-back function
u : input rank-1 array('d') with bounds (n)
Other Parameters
----------------
func_extra_args : input tuple, optional
Default: ()
n : input int, optional
Default: len(u)
Returns
-------
jacob : rank-2 array('d') with bounds (n,n)
Notes
-----
Call-back functions::
def func(up,gvalp,[n]): return
Required arguments:
up : input rank-1 array('d') with bounds (n)
gvalp : input rank-1 array('d') with bounds (n)
Optional arguments:
n : input int, optional
Default: len(up)
回调函数func有问题。输入应为数组U,输出应为gval。相反,U 和 gval 都是输入,没有输出。
我假设这就是以下 python 脚本返回垃圾的原因:
import num_Dg
def func(x):
return [x[0]**2,x[1],x[2]]
val = [0.2,0.1,0.1]
num_Dg.num_dg(func,val)
我应该如何更改我的 fortran 子例程以使回调函数正常工作?任何建议表示赞赏!谢谢!
【问题讨论】:
-
欢迎,我建议采取tour。您是否尝试为该功能编写接口块?当它是外部的时,函数签名并没有真正准确地定义。 F2py 可能会尝试以某种方式从调用中推断。
-
您在对
func的两次调用中都使用了Up,而看起来您应该依次使用Up和Um。另外,您为什么不将func称为gvalp = func(Up)?或者在参数中包含n?声明为外部的函数对其语法的检查为零,这可能会受到限制。 -
谢谢!你说的对。在
func的第二次调用中,Up应替换为Um。我将代码更改为g_valp = func(Up,n),其中n 是Up的维度。现在,当我执行print(num_dg.__doc__)时,我得到了一种看起来合理的回调函数格式。但是,当我从 python 调用程序时,我没有得到好的结果。从在 fortran 程序中使用 print 语句,看起来原因是因为对 func 的调用无法正常工作。 @VladimirF,我使用 f2py 为函数生成了一个接口块。