SymPy学习笔记

SymPy其实是一个冷门Python库
介绍里说SymPy是个用于符号数学的Python库
现在来看依然是个大玩具，尽管看起来包罗万象(广度上对标Mathematica)
有些地方比较鸡肋，，，，食之无味，弃之可惜
但是补充了python在符号计算的短板，这个项目还是挺牛的
没打算学完

我遇到好多地方找不到对应模块或者函数。。。

目录

• 参考
• 你可以用SymPy做什么
• 数论
• 丢番图方程
• discrete
• 物理
  • • optics
    • utilities
    • waves
    • hydrogen
    • matrices
    • qho_1d
    • sho

参考

https://docs.sympy.org/ 英语文档
https://tio.run 一个在线IDE(语言选择python 3)
https://live.sympy.org/ 一个sympy的在线IDE，可是等待结果的时间太长了
https://github.com/sympy/sympy github仓库地址

你可以用SymPy做什么

范畴论
组合数学
密码学
偏微分方程
几何学
ODE
PDE
物理
。。。

数论

动态增长的埃拉托色尼筛

from sympy import sieve
sieve._reset()
sieve.extend(100)
# sieve.extend_to_no(30)
print(25 in sieve)
print(sieve._list)

数论函数，离散对数，连分数，埃及分数，等等

from sympy import *
from sympy.ntheory import *
from sympy.ntheory.factor_ import *
from sympy.ntheory.continued_fraction import *


print(mobius(13*7*5))
print(totient(30))
print(reduced_totient(30))
print(divisor_sigma(12,2))
print(udivisor_sigma(36,3))
print(core(15**11,10))
print(digits(27,b=2))
print(primenu(30))
print(primeomega(20))
print([legendre_symbol(i,7) for i in range(7)])
print(jacobi_symbol(60,121))
print(discrete_log(41,15,7))
print(continued_fraction_periodic(3,2,7))
print(egyptian_fraction(Rational(3,7)))
print(egyptian_fraction(Rational(3, 7), "Takenouchi"))

丢番图方程

from sympy import *
from sympy.solvers import *
from sympy.solvers.diophantine import *



x,y,z=symbols("x,y,z",integer=True)
print(diophantine(2*x+3*y-5))

discrete

讲了快速傅里叶变换，快速数论变换，快速Walsh-Hadamard变换，莫比乌斯变换，等等

from sympy  import *
from sympy.discrete.convolutions import *



x,y,z=symbols("x,y,z")
u,v,w=symbols("u,v,w")



####################################
tmp1=fft([1,2,3,4],dps=15)
tmp2=ifft(tmp1)
print(f'{tmp1}\n{tmp2}')


tmp3=ntt([1,2,3,4],prime=3*2**8+1)
tmp4=intt(tmp3,prime=3*2**8+1)
print(f'{tmp3}\n{tmp4}')


tmp5=fwht([4,2,2,0,0,2,-2,0])
tmp6=ifwht(tmp5)
print(f'{tmp5}\n{tmp6}')



tmp7=mobius_transform([x,y,z],subset=False)
tmp8=inverse_mobius_transform(tmp7,subset=False)
print(f'{tmp7}\n{tmp8}')




tmp9=convolution_subset([u,v],[x,y])
tmp10=covering_product([u,v],[x,y])
tmp11=intersecting_product([u,v],[x,y])
print(f'{tmp9}\n{tmp10}\n{tmp11}')

物理

optics

• 高斯光学

from sympy.physics.optics import *


p=BeamParameter(530e-9,1,w=1e-3)
print(f'p.q    {p.q}')
print(f'p.w    {p.w}')
print(f'p.w_0    {p.w}')
print(f'p.gouy    {p.gouy}')
print(f'p.radius    {p.radius}')
print(p.q.n())
print(p.w_0.n())
print(p.z_r.n())



print("#######")
fs=FreeSpace(10)
p1=fs*p
print(p.w.n())
print(p1.w.n())

• 光学媒质

from sympy.abc import *
from sympy.physics.optics import *

m1=Medium('m1')
m2=Medium('m2',epsilon,mu)
print(f'm1.intrinsic_impedance {m1.intrinsic_impedance}')
print(f'm2.intrinsic_impedance {m2.intrinsic_impedance}')
# m1.intrinsic_impedance 149896229*pi*kilogram*meter**2/(1250000*ampere**2*second**3)
# m2.intrinsic_impedance sqrt(mu/epsilon)


print()
print(f'm1.refractive_index {m1.refractive_index}')
print(f'm2.refractive_index {m2.refractive_index}')
# m1.refractive_index 1
# m2.refractive_index 299792458*meter*sqrt(epsilon*mu)/second




print()
print(f'm1.permeability {m1.permeability}')
print(f'm2.permeability {m2.permeability}')
# m1.permeability pi*kilogram*meter/(2500000*ampere**2*second**2)
# m2.permeability mu



print()
print(f'm1.permittivity {m1.permittivity}')
print(f'm2.permittivity {m2.permittivity}')
# m1.permittivity 625000*ampere**2*second**4/(22468879468420441*pi*kilogram*meter**3)
# m2.permittivity epsilon


print()
print(f'm1.speed {m1.speed}')
print(f'm2.speed {m2.speed}')
# m1.speed 299792458*meter/second
# m2.speed 1/sqrt(epsilon*mu)

utilities

包括：折射角，菲涅耳系数，brewster_angle，临界角，磨镜者公式，mirror公式，lens公式，超焦距，横向放大率

from sympy import *
from sympy.physics.optics import *
from sympy.abc import *

# 计算布儒斯特角
print(brewster_angle(1,1.33))

# 计算超焦距
print(hyperfocal_distance(f=0.5,N=8,c=0.0033))

# 1/v-1/u=1/f 是lens公式大家应该都知道吧
print(lens_formula(focal_length=f,u=u))

# 磨镜者公式
print(lens_makers_formula(1.33,1,10,-10))

# 1/v+1/u=1/f 是mirror公式大家都知道吧
print(mirror_formula(u=u,v=v))

# 横向放大率
print(transverse_magnification(si=30,so=15))

waves

TWave类，里面有若干属性成员

from sympy import *
from sympy.physics.optics import *

A1,phi1,A2,phi2,f=symbols('A1,phi1,A2,phi2,f')
W1=TWave(A1,f,phi1)
W2=TWave(A2,f,phi2)
W3=W1+W2
print(W3)
print(W3.amplitude)
print(W3.phase)
print(W3.speed)
print(W3.wavelength)
print(W3.angular_velocity)
print(W3.frequency)
print(W3.time_period)
print(W3.wavenumber)

hydrogen

氢原子波函数

from sympy.physics.hydrogen import *
from sympy.abc import *
from sympy import *

print(E_nl(n,Z))
print(E_nl_dirac(n=2,l=1,spin_up=False))
print(R_nl(n=n,l=l,r=r,Z=Z))


# 验证归一性
print(integrate(R_nl(1,0,r,Z=2)**2*r**2,(r,0,oo)))

matrices

介绍物理中常见的矩阵：狄拉克Gamma矩阵，泡利矩阵，惯量张量矩阵，等等

from sympy.physics.matrices import *
from sympy.abc import *
from sympy import *

dx,dy,dz=symbols('x,y,z')

print(mgamma(mu=1,lower=False))
print(msigma(i=1))
print(pat_matrix(m=m,dx=dx,dy=dy,dz=dz))

qho_1d

介绍了一维量子谐振子：总能量，波函数，等等

from sympy.physics.qho_1d import *
from sympy.abc import *
from sympy import *

omega=symbols('omega')


print(E_n(x,omega))
print(psi_n(n=0,x=x,m=m,omega=omega))

sho

介绍了三维量子谐振子：总能量，径向波函数，等等

from sympy.physics.sho import *
from sympy.abc import *
from sympy import *

hw=symbols('hw')

print(E_nl(n=n,l=l,hw=hw))
print(R_nl(n=0,l=0,nu=nu,r=r))