千家信息网

如何利用Python做科学计算

发表于:2024-11-15 作者:千家信息网编辑
千家信息网最后更新 2024年11月15日,这篇文章主要为大家展示了"如何利用Python做科学计算",内容简而易懂,条理清晰,希望能够帮助大家解决疑惑,下面让小编带领大家一起研究并学习一下"如何利用Python做科学计算"这篇文章吧。其核心要
千家信息网最后更新 2024年11月15日如何利用Python做科学计算

这篇文章主要为大家展示了"如何利用Python做科学计算",内容简而易懂,条理清晰,希望能够帮助大家解决疑惑,下面让小编带领大家一起研究并学习一下"如何利用Python做科学计算"这篇文章吧。

其核心要求是做一个函数拟合,但是被拟合函数是个积分表达式。最简单的方法是利用scipy库中的函数来做,下面是源代码。

# -*- coding: utf-8 -*-from scipy import integratefrom scipy import optimizefrom matplotlib import pyplotimport numpyimport pandasimport timeI = complex(0, 1)def tmp(E, params):    global I    Z = params[1]    delta = params[2]    Gamma = params[3]    complex_num = 0.5 + numpy.sqrt(pow(E+I*Gamma, 2)-delta*delta)/(2*(E+I*Gamma))    alpha = complex_num.real    eta = complex_num.imag    beta = 1 - alpha    gamma = numpy.sqrt( numpy.power(alpha+Z*Z*(alpha-beta), 2)          + numpy.power(eta*(2*Z*Z+1), 2) )    return alpha, beta, gamma, etadef factor(E, params):    P = params[0]    Z = params[1]    alpha, beta, gamma, eta = tmp(E, params)    numerator1 = numpy.sqrt((alpha*alpha+eta*eta)*(beta*beta+eta*eta))    denominator1 = gamma*gamma    AE = numerator1/denominator1    numerator2 = Z*Z*( numpy.power((alpha-beta)*Z-2*eta, 2)               + numpy.power(2*eta*Z+(alpha-beta), 2) )    denominator2 = gamma*gamma    BE = numerator2/denominator2    return 1+(1-P)*AE-BEdef dfdV_mod(E, V):    variable = E - V    if (variable >= 0):        exp = numpy.exp(-variable)    else:        exp = numpy.exp(variable)    numerator = -exp    denominator = numpy.power(exp+1, 2)    return numerator/denominator# 被积函数def integrand(E, V, params):    return dfdV_mod(E, V)*factor(E, params)# 积分def integral(V, params):    result = integrate.quad(integrand, -numpy.inf, numpy.inf, args=(V, params))    return result[0]# 画图时计算积分def integral_all(V, params):    result = numpy.zeros(V.size)    for i in range(0, V.size):        res = integrate.quad(integrand, -numpy.inf, numpy.inf, args=(V[i], params))        result[i] = res[0]    return result# 实验测量值与理论值的偏差def residual(params, g, V):    res = numpy.zeros(g.size)    for i in range(0, g.size):        res[i] = g[i] - integral(V[i], params)    return res# 将实验数据读入,需在实验数据中添加表头def ReadData(path):    dataframe = pandas.read_excel(path, sheet_name=0)    x = numpy.array(dataframe.iloc[:, 0])    y = numpy.array(dataframe.iloc[:, 1])    return y, xif __name__ == '__main__':    start = time.time()    # 赋初值    P = 0.5    Z = 1    delta = 0    Gamma = 0    # 读入数据    g, V = ReadData("data.xlsx")    # 最小二乘法拟合    params0 = numpy.array([P, Z, delta, Gamma])    result = optimize.least_squares(residual, params0, bounds=([-1, -5, -1, -1], [1, 5, 1, 1]), args=(g, V))    print("result: " + str(result))        # 画出结果    params = result.x    V_test = numpy.linspace(V.min(), V.max(), 100)    g_test = integral_all(V_test, params)    pyplot.plot(V, g, 'o', markersize=1, label='data')    pyplot.plot(V_test, g_test, label='fitted curve')    pyplot.xlabel('V')    pyplot.ylabel('g')    pyplot.show()        end = time.time()    print("time elapsed: " + str(end-start) + "s")

从代码来看还是很清晰的,主要用到两个函数,一个是integrate.quad,用来算积分,另一个是optimize.least_squares,利用最小二乘法给出函数参数值。

以上是"如何利用Python做科学计算"这篇文章的所有内容,感谢各位的阅读!相信大家都有了一定的了解,希望分享的内容对大家有所帮助,如果还想学习更多知识,欢迎关注行业资讯频道!

0