f2py,返回数组的Python函数(向量值函数)
发布时间:2020-12-16 21:36:36 所属栏目:Python 来源:网络整理
导读:在下面的 Python中,我有五个函数包含在func返回的数组中,我必须集成它.代码调用使用f2py生成的外部Fortran模块: import numpy as npfrom numpy import cos,sin,expfrom trapzdv import trapzdvdef func(x): return np.array([x**2,x**3,cos(x),sin(x),exp(x
在下面的
Python中,我有五个函数包含在func返回的数组中,我必须集成它.代码调用使用f2py生成的外部Fortran模块:
import numpy as np from numpy import cos,sin,exp from trapzdv import trapzdv def func(x): return np.array([x**2,x**3,cos(x),sin(x),exp(x)]) if __name__ == '__main__': xs = np.linspace(0.,20.,100) ans = trapzdv(func,xs,5) print 'from Fortran:',ans print 'exact:',np.array([20**3/3.,20**4/4.,sin(20.),-cos(20.),exp(20.)]) Fortran例程是: subroutine trapzdv(f,nf,nxs,result) integer :: I double precision :: x1,x2 integer,intent(in) :: nf,nxs double precision,dimension(nf) :: fx1,fx2 double precision,intent(in),dimension(nxs) :: xs double precision,intent(out),dimension(nf) :: result external :: f result = 0.0 do I = 2,nxs x1 = xs(I-1) x2 = xs(I) fx1 = f(x1) fx2 = f(x2) result = result + (fx1+fx2)*(x2-x1)/2 enddo return end 问题是Fortran只在func(x)中集成了第一个函数. from Fortran: [ 2666.80270721 2666.80270721 2666.80270721 2666.80270721 2666.80270721] exact: [ 2.66666667e+03 4.00000000e+04 9.12945251e-01 -4.08082062e-01 4.85165195e+08] workarond的一种方法是修改func(x)以返回给定的值 def func(x,i): return np.array([x**2,exp(x)])[i-1] 然后更改Fortran例程以使用两个参数调用该函数: subroutine trapzdv(f,x2,fx1,fx2 integer,nxs x1 = xs(I-1) x2 = xs(I) do J = 1,nf fx1 = f(x1,J) fx2 = f(x2,J) result(J) = result(J) + (fx1+fx2)*(x2-x1)/2 enddo enddo return end 哪个有效: from Fortran: [ 2.66680271e+03 4.00040812e+04 9.09838195e-01 5.89903440e-01 4.86814128e+08] exact: [ 2.66666667e+03 4.00000000e+04 9.12945251e-01 -4.08082062e-01 4.85165195e+08] 但是这里func的调用次数是必要的5倍(在实际情况下是func >有没有人知道一个更好的解决方案,使Fortran识别func(x)返回的所有数组?换句话说,使Fortran将fx1 = f(x1)构建为一个数组,其中5个元素对应于func(x)中的函数. OBS:我正在使用f2py编译-c –compiler = mingw32 -m trapzdv trapzdv.f90 解决方法
不幸的是,你无法将数组从python函数返回到Fortran.你需要一个子程序(意味着它用call语句调用),这就是f2py不允许你做的事情.
在Fortran 90中,您可以创建返回数组的函数,但这也不是f2py可以执行的操作,尤其是因为您的函数不是Fortran函数. 您唯一的选择是使用循环解决方法,或重新设计您希望python和Fortran如何交互. (编辑:李大同) 【声明】本站内容均来自网络,其相关言论仅代表作者个人观点,不代表本站立场。若无意侵犯到您的权利,请及时与联系站长删除相关内容! |