基于 Sage 的坐标变换(符号运算)

前日做作业时要写一些坐标变换的题目, 深感手算无望, 便使用 Sage 的符号运算功能高效地完成了所有的计算, 于是把主要的代码贴到此博客中, 以后再用相关的计算, 就不用重复输入这些代码了.

  (xx, yy, zz, yz, zx, xy) = var('xx yy zz yz zx xy')
  (r, theta, phi, z) = var('r theta phi z') 
  
  T_Polar = Matrix([              ## matrices must be defined this way.
      [cos(theta), sin(theta), 0], 
      [-sin(theta), cos(theta), 0], 
      [0, 0, 1]])
  
  T_Sophere = Matrix([
      [sin(theta)*cos(phi), sin(theta)*sin(phi), cos(theta)],
      [cos(theta)*cos(phi), cos(theta)*sin(phi), -sin(theta)],
      [-sin(phi), cos(phi), 0]])
         
  
  C = Matrix([[xx, xy, zx],
      [xy, yy, yz],
      [zx, yz, zz]])
  
  R_Polar = T_Polar * C * T_Polar.T
  
  Result_Polar = R_Polar(xx=r^2*cos(theta)^2*z, yy=r^2*sin(theta)^2*z,
                         xy=r^2*sin(theta)*cos(theta)*z, zz = z^3, 
                         yz = r*z^2*sin(theta), zx=r*z^2*cos(theta)).simplify_trig()
  
  R_Sophere = T_Sophere * C * T_Sophere.T
  
  Result_Sophere = R_Sophere(xx=sin(theta)^2*cos(phi)^2*cos(theta), 
                             yy=sin(theta)^2*sin(phi)^2*cos(theta),
                             zz=cos(theta)^3,
                             yz=sin(theta)*sin(phi)*cos(theta)^2,
                             zx=sin(theta)*cos(phi)*cos(theta)^2,
                             xy=sin(theta)^2*cos(theta)*sin(phi)*cos(phi)).simplify_trig()

后面算 Result_Sophere 和 Result_Polar 是当时计算的题目, 不是普适的程序, 贴在这里吐槽老师留作业太随意, 这尼玛是人算的吗?

结果老师下一次作业的计算是下面这个程序完成的:

x, y, z, A = var('x', 'y', 'z', 'A')
r = var('r')
var('a', latex_name='\lambda')
b = var('mu')
r(x, y, z) = sqrt(x^2 + y^2 + z^2)
u(x, y, z) = A*x*z/r^3
uxx = (u.diff(x)).diff(x)
uxx.show()
v(x, y, z) = A*y*z/r^3
vxx = (v.diff(x)).diff(x)
vxx.show()
w(x, y, z) = A*(z^2/r^3 + ((a+3*b)/(a+b))*(1/r))
wxx = (w.diff(x)).diff(x)
wxx.show()
uyy = (u.diff(y)).diff(y)
uyy.show()
uzz = (u.diff(z)).diff(z)
uzz.show()
vxy = (v.diff(x)).diff(y)
vxy.show()
wxz = (w.diff(x)).diff(z)
wxz.show()
(uxx + uyy + uzz + (a+b)/b * (uxx + vxy + wxz)).simplify_radical().show()
m = vector([u, v, w])
Gamma = 1/2 * (m.diff() + m.diff().transpose())
Gamma.show()

后面一段坐标变换的代码就不再贴出来了....

第二个程序有很多化简步骤的输出有很大的问题，下一篇博客里会仔细研究 Sage 的化简表达式的技巧。