当我使用 mpmath正函数增加积分范围时,积分值会降低

问题描述

我在 python 上使用 mpmath 库来计算这个积分:

enter image description here

enter image description here

问题是,当我增加积分范围时,我的积分值慢慢减少到0..而被积函数为正!!我给你一些价值观:

积分范围:[-1.17,-1.11],积分值:1.28479400889196e-8

积分范围:[-1.2,0],积分值:2.87226592130464e-29

积分范围:[-10,10],积分值:0.0

除了只选择我的被积函数不是 0 的范围之外,发生了什么以及如何解决问题?

这是我的代码(费米狄拉克统计和电导率已针对 pmos 进行了调整,因此它不完全是我展示的积分)

import numpy as np
from math import pi,sqrt,exp,log
from scipy.integrate import quad
import scipy.constants as phys
import matplotlib.pyplot as plt
from decimal import Decimal,getcontext
getcontext().prec = 1000

import mpmath as mp


k = (phys.Boltzmann)              # Boltzmann constant
                                  #Parenthesis left due to previous decimal conversion
m = (9.1094e-31)     #electron mass

h = (phys.Planck)     #Planck constant
h_b = h/(2 * pi)        #reduced

q = (1.602176634e-19)  
e =   (1)      
dE = 3e-3 * e

pi = (pi)


m_e = m * (1.06)   #electron effective mass
m_h = m * 0.59    #hole effective mass
A_2D = m_h / (pi*h_b**2) 


###### INPUTS
T = (4.2)
V_d = (0.02)
chi = (4.05)*e
E_g = (1.16)*e


g = 4                    #degeneracy
μ = (1e3)              #mobility,cm**2/Vs
A = A_2D

W = 1.00E-05
L = 1.00E-07

#########
E_c = 0
E_v = E_c - E_g

Phi_f = (0.55)
E_f = ( E_c - E_g/2 - e*Phi_f)
Ef_d = ( E_f - e * V_d)



#######

def fd(E,E_f):             #Fermi-Dirac distribution for holes
    res = 1/(1 + mp.exp(((E - E_f)*q/(k*T))))
    return 1-res

def sigma(E):
    if E_v-E < 2 :
        res = q * A * μ *1e-4 *  dE*q *  log( 1 + mp.exp((E_v-E)/ dE)) 
    else:
        res = - q * A * μ *1e-4 *q *(E-E_v)
    return res

def integrand(E):
    return sigma(E)*( fd(E,Ef_d) -fd(E,E_f))

X = np.linspace(E_v-5.05,5,1000)

FD = np.vectorize(fd)

Sigma = np.vectorize(sigma)



Integrand = np.vectorize(integrand)
#mp.quad
Int = mp.quad(integrand,[-10,10])
I_c = W/(q*L)*Int

fig,ax1 = plt.subplots()
ax1.plot(X,FD(X,Ef_d) - FD(X,E_f),color = 'g',label='fd stat');

ax2 = ax1.twinx()

ax2.semilogy() ;
ax2.plot(X,Sigma(X),label='sigma');
ax2.plot(X,Integrand(X),label='integrand');
ax1.legend()

ax2.legend(loc='lower right')
print('Int = ',Int)

谢谢

解决方法

暂无找到可以解决该程序问题的有效方法,小编努力寻找整理中!

如果你已经找到好的解决方法,欢迎将解决方案带上本链接一起发送给小编。

小编邮箱:dio#foxmail.com (将#修改为@)