这篇文章主要为大家展示了“如何利用Python做科学计算”,内容简而易懂,条理清晰,希望能够帮助大家解决疑惑,下面让小编带领大家一起研究并学习一下“如何利用Python做科学计算”这篇文章吧。
成都创新互联公司是一家专注于网站制作、做网站与策划设计,化州网站建设哪家好?成都创新互联公司做网站,专注于网站建设10年,网设计领域的专业建站公司;建站业务涵盖:化州等地区。化州做网站价格咨询:13518219792
其核心要求是做一个函数拟合,但是被拟合函数是个积分表达式。最简单的方法是利用scipy库中的函数来做,下面是源代码。
# -*- coding: utf-8 -*- from scipy import integrate from scipy import optimize from matplotlib import pyplot import numpy import pandas import time I = 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, eta def 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-BE def 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, x if __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做科学计算”这篇文章的所有内容,感谢各位的阅读!相信大家都有了一定的了解,希望分享的内容对大家有所帮助,如果还想学习更多知识,欢迎关注创新互联行业资讯频道!