Python/Scipy 集成数组

Question

我正在尝试编写一个执行以下操作的程序:

• 从数组中获取 V 的值
• 将 V 值传递到关于 E 的积分
• 将积分结果输出到数组I中
• 绘制 I 对 V

方程式看起来很糟糕，但除了 V 之外，一切都是常数。Here is the equation .等式不是很重要。

我该如何解决这个问题？我的尝试(如下所示)并未计算从文件中读取的每个 V 值的积分。

from scipy import integrate #integrate.quad
from numpy import *
import pylab
import datetime
import time
import os
import math

# import V
fn = 'cooltemp.dat'
V = loadtxt(fn,unpack=True,usecols=[1])

# variables
del1, del2, R, E, fE, fEeV = 1,2,1,2,1,1
e = 1.602176565*10**-19

# eqn = dint(abc)
a = E/( math.sqrt( E**2 - del1**2 ) )
b = ( E+ e*V )/( math.sqrt( ( E + e*V )**2) - del2**2)
c = fE-fEeV
d = 1/(e*R) # integration constant
eqn = a*b*c

# integrate 
result = quad(lambda E: eqn,-inf,inf)

# current
I = result*d

# plot IV curve
pylab.plot(V,I,'-r')

## customise graph
pylab.legend(['degree '+str(n),'degree '+str(q),'data'])
pylab.axis([0,max(x),0,max(y)])
pylab.xlabel('voltage (V)')
pylab.ylabel('current (A)')
tc = datetime.datetime.fromtimestamp(os.path.getmtime(fn))
pylab.title('IV curve\n'+fn+'\n'+str(tc)+'\n'+str(datetime.datetime.now()))
pylab.grid(True)
pylab.show()

* 更新尝试:

from scipy import integrate
from numpy import *
import pylab
import datetime
import time
import os
import math

# import V
fn = 'cooltemp.dat'
V = loadtxt(fn,unpack=True,usecols=[1])
# print V

# variables
del1, del2, R, E, fE, fEeV = 1.0,2.0,1.0,2.0,1.0,1.0
e = 1.602176565*10**-19

I=[]
for n in range(len(V)):

    constant = 1/(e*R) # integration constant
    eqn = (E/( math.sqrt( E**2 - del1**2 ) ))*(( E + e*V[n] )/( math.sqrt( ( E + e*V[n] )**2) - del2**2))*(fE-fEeV)

    # integrate 
    result,error = integrate.quad(lambda E: eqn,-inf,inf)
    print result
    # current
    I.append(result*constant)

I = array(I)

# plot IV curve
pylab.plot(V,I,'-b')

Answer 1

你有几个问题:

您传递给 quad 的“函数”总是返回 eqn，这只是一个预先计算好的数字。您需要定义一个适当的函数，将 E 的给定值作为输入并返回被积函数。此函数还需要为 V 假定一个固定值。假设您提供的代码针对给定的 V 和 E 值计算出正确的数量(我没有检查，只是复制粘贴):

# import V
fn = 'cooltemp.dat'
V = loadtxt(fn,unpack=True,usecols=[1])
# print V

@np.vectorize
def result(x):
    def integrand(E):
        del1, del2, R, fE, fEeV = 1.0,2.0,1.0,1.0,1.0
        e = 1.602176565*10**-19
        a = E/( math.sqrt( E**2 - del1**2 ) )
        b = ( E+ e*x )/( math.sqrt( ( E + e*x )**2) - del2**2)
        c = fE-fEeV
        d = 1/(e*R) # integration constant
        return a * b * c
    return quad(integrand, -inf, inf)

I = result(V)

总结:

• result(v) 计算固定值 v 的完整积分(在 E 上)
• integrand(E) 在固定 E(积分变量)和固定 V(它从函数外部获取值，这就是被积函数的定义是嵌套的原因)计算被积函数在结果的定义中)
• @np.vectorize 技巧只是一个很好的便利函数，它允许您将 V 的数组传递给 result。 Numpy 将为您遍历这些值，并返回一个数组而不是标量