python解决物理问题_利用Python科学计算处理物理问题（和物理告个别）

背景：

2019 年初由于尚未学习量子力学相关知识，所以处于自学阶段。浅显的学习了曾谨言的量子力学一卷和格里菲斯编写的量子力学教材。注重将量子力学的一些基本概念了解并理解。同时老师向我们推荐了 Quantum Computation and Quantum Information 这本教材，了解了量子信息相关知识。

2019 年暑假开始量子力学课程的学习，在导师的推荐下，从 APS (美国物理学会)和 AIP (美国物理联合会)下载了与量子纠缠( Quantum Discord )相关的著名的文献和会议报告，了解了量子信息的发展历程和一些杰出的理论。其中 Unified View of Quantum and Classical Correlations 和 Quantum Discord :A Measure of Quantumness of Correlations 两篇文章影响最为深刻。对量子信息领域有了初步认识。

我也参加了相关的量子相关的报告，譬如 12 月 18 日陆朝阳教授的量子光学与量子计算背景和进程介绍， 2019 年 10 月 9 日郭光灿院士的《量子之问》，这些讲座都激发了我对量子计算、量子通信的兴趣。

我也利用空闲时间自学了 python ，掌握了实验编程所需要的基本技能，强化了自己在编程方面的知识，也学会部分 LATEX 进行论文编写。

在参加项目过程中，虽然对投身于人类探索未知及其佩服，但终究自知穷极一生也极难在物理基础领域做出杰出贡献，就此转向计算机，愿尽绵薄之力，用技术为社会做一些有价值有意义之事。

科学计算，利用python进行相关图像整理。学会了基本的3d图像，numpy，matpolib绘图工具。老实说一些底层的原理并不清楚，但还是可以总结一些经验教训的。

1. 绘制圆柱体，利用小技巧q/q，看似是1，实际上是生成了1的一个数组。

如上图，我们希望将这样的二维图像纵向拉伸变成类似于圆柱体一样的三维图像，与另外的图像进行对比，但没有找 到相关方 法， 因为这是二维图像，想要实现三维立体图就需要有两个变量作为基底，而pAB是一个单变量函数，想要加入一个变量却不改变它最终的函数值似乎是不可能实现的。

我尝试用1去直接作为第二个变量，但执行无法通过，当时没有想明白。(后面会说明)

灵机一动，我尝试用q2的具体函数值去代替q2取遍0-1之间所有常数，直接乘个q2试试？然后就有了下图。但不行啊，想要q2取遍0-1所有值，又要q2在函数中显示，不好搞。

那，要不再除个q2试试？？？然后就卧槽了，成功了。？？？。当时就很迷惑，这个和直接乘1有什么区别？ *q2/q2，不就是*1吗？然后分析发现，发现通过*q2/q2操作，实际上是生成了一个全是1的数组，这样就达到了形成二维基底的基本要求，从而画出了三维图像。

import numpy as np

from matplotlib import pyplot as plt

from mpl_toolkits.mplot3d import Axes3D

q1 = np.arange(0.01, 1, 0.01)

q2 = np.arange(0.01, 1 , 0.01) //生成一位基底

q1, q2 = np.meshgrid(q1, q2) // 混合成二维数组，形成二维基底

pCDa = (1-q1)

pCDb = (np.sqrt((1-q1)**2+q1**2)-q1)

alpha = (pCDa + pCDb) / (pCDa + 4 * pCDb)

beta = pCDb / (pCDa + 4 * pCDb)

pCDp00 = ( q1* pCDa ** 2 ) / ( pCDa*pCDa +pCDb*pCDb)

pCDp10 = ( q1* pCDb ** 2 ) / ( pCDa*pCDa +pCDb*pCDb)

pCDp01 = ( (1-q1) / 2 ) * ( pCDa + pCDb ) ** 2 / ( pCDa*pCDa +pCDb*pCDb)

pCDp11 = ( (1-q1) / 2 ) * ( pCDa - pCDb ) ** 2 / ( pCDa*pCDa +pCDb*pCDb)

s_x_pCD= -pCDp00 * np.log2(pCDp00) - pCDp01 * np.log2(pCDp01) - pCDp10 * np.log2(pCDp10) - pCDp11 * np.log2(pCDp11)

s_pCD = -q1* np.log2(q1) - (1-q1) * np.log2(1-q1)

Q_MID1 = (s_x_pCD - s_pCD) *q2 /q2 #AB或CD的关联值

fig = plt.figure()

ax = Axes3D(fig)

ax.plot_surface(q1,q2,Q_MID1)

ax.set_xlabel('value of q2')

ax.set_ylabel('value of q1')

ax.set_zlabel('the value of Q_MID1(pCD)')

plt.show()

2. 寻找图像交点，并没有使用到复杂的方法，而是简单的使用了精度为0.01周围四个点函数值相加取平均值。但数值分析里面一定说过相关的方法。

下面简单介绍下 数值分析 基本内容。

1、误差分析

模型误差(理论到模型，忽略次要因素必然带来误差)

观测误差

截断误差(无限到有限)

舍入误差(无理数到有理数)

2、求非线性方程组 f(x)=0的解

x = g(x) 的迭代法(fixed point)

牛顿-拉夫森法(Newton-Raphson)和割线法(Secant Methods)

3、求线性方程组 AX=B 的数值解法

高斯消元法(LU分解)

迭代法

雅各比行列式

高斯赛德尔

4、插值和多项式逼近

拉格朗日逼近

牛顿多项式

切比雪夫多项式

帕德逼近

5、曲线拟合

最小二乘法拟合曲线

样条函数插值

贝赛尔曲线

6、数值积分

组合梯形公式和辛普森公式

递归公式和龙贝格积分

高斯勒让德积分

7、微分方程求解

欧拉方法

休恩方法

我们希望求这个交点，算是求解非线性方程 f(x) = 0。

尝试用牛顿拉普森方法求解。(其实就是将试位法斜率用 f ' (x) 代替)

看起来还不错，我们试着放到项目上

割线法

3. 一些细节处理，比如将图像内置边框，图像标注等

plt.text(0,1,r'$max\ is\ 1.34$') /