数值积分函数的图形

如何解决数值积分函数的图形

我需要绘制图表。我写的代码是：

import numpy as np
from scipy.integrate import quad
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt
import pandas as pd

mt=3.75*10**(-5) # kg massa total empregada da amostra
rho=5961 #kg/m**3 densidade da amostra
v0=mt/rho #m**3 volume da amostra
uim=0.453 #A*m**2 momento magnético ímã
g=9.82 #m/s**2 aceleração da gravidade na Terra
z0=0.0063 #m distância centro ímã parte inferior da amostra
z1=0.0093 #m distância centro ímã partesuperior da amostra
R=0.0049 #m raio do porta amostra cilíndrico
dis=(1/z0**4)-(1/z1**4)-(z0**2+R**2)/(3*(z0**2+R**2)**3)+(z1**2+R**2)/(3*(z0**2+R**2)**3)
Msv=470644.794 #A/m magnetização volumétrica
Dm=17.6*10**-9 #m diâmetro médio
sigmad=0.13 #fator de dispersão adimensional
u0=1.256*10**-6 #N/A**2 permeabilidade magnética
H=70823.95 #A/m campo aplicado
Kb=1.3807*10**-23 #m**2*kg/(s*K) constante de Boltzman
T=300 #K temperatura
Kef=105000
Hk_=2*Kef/(u0*Msv) #campo de anisotropia
t0=10**-9 #s fator de frequência
h_=H/Hk_
v1=6*Kb*T

def func(t):
    Multiplicador=3*u0*v0*uim/(32*np.pi*g)*dis*Msv*np.exp(2*sigmad**2)/(Dm*np.sqrt(2*np.pi)*sigmad)
    Integrando1=lambda D: (np.exp((-(np.log(D/Dm))**2)/(2*sigmad**2)))
    Integrando2=lambda D: (1/np.tanh(u0*np.pi*Msv*H*D**3/v1)-v1/(u0*np.pi*Msv*H*D**3))
    Integrando3=lambda D,_t: (1-np.exp(-_t*((1-(h_)**2)*((1-h_)*np.exp(-Kef*np.pi*D**3/v1*(1-h_)**2)+(1+h_)*np.exp(-Kef*np.pi*D**3/v1*(1+h_)**2)))/t0))
    Integrando4=lambda D: (1-np.exp(-27105.05*((1-(h_)**2)*((1-h_)*np.exp(-Kef*np.pi*D**3/v1*(1-h_)**2)+(1+h_)*np.exp(-Kef*np.pi*D**3/v1*(1+h_)**2)))/(t0)))
    fn1=lambda D: Integrando1(D)*Integrando2(D)*Integrando4(D)
    ma=np.asarray([quad(lambda D: Integrando1(D)*Integrando2(D)*Integrando3(D,_t),Dm/10,10*Dm)[0] for _t in t])
    ma1=np.asarray(quad(fn1,10*Dm)[0])
    return ma*Multiplicador

xFit=np.arange(0,27105,65)
plt.plot(xFit,func(xFit))
plt.show()

但数据是从以下文档中获得的，因此预计该图将采用类似于 Expected_graph 的值，但我得到了 getting_graph。 https://github.com/Henrique0501/MSM 怎么了？