在Python中写一个双和

Pal*_*l86 5 python numpy

我是StackOverflow的新手,我对Python非常陌生.

我的问题是这个......我需要写一个双和,如下:

方程

动机是这是对用于大地水准面的重力势的角度修正.

  1. 我写这笔钱很困难.并且,在你说"去这样一个这样的资源"或者对我不耐烦之前,这是我第一次完成编码/编程/不管这是什么.

  2. 这是一个使用"for"循环的好地方吗?

  3. 我有两个索引(n,m)和系数c_{nm}以及s_{nm}.txt文件的数据.每个项目都是一列.当我说usecols,我是0通过3还是1通过4


(上面的等式)

\begin{equation}
V(r, \phi, \lambda) = \sum_{n=2}^{360}\left(\frac{a}{r}\right)^{n}\sum_{m=0}^{n}\left[c_{nm}*\cos{(m\lambda)} + s_{nm}*\sin{(m\lambda)}\right]*\sqrt{\frac{(n-m)!}{(n+m)!}(2n + 1)(2 - \delta_{m0})}P_{nm}(\sin{\lambda})
\end{equation}
Run Code Online (Sandbox Code Playgroud)

jre*_*nie 4

(2) 是的,“for”循环就可以了。正如 @jpmc26 所指出的,生成器表达式是“for”循环的一个很好的替代方案。IMO,如果效率对你很重要,你会想要使用 numpy。

(3) 正如 @askewchan 所指出的,“usecols”指的是genfromtxt的参数;正如该文档中所指定的,列索引从 0 开始,因此您需要使用 0 到 3。

一个简单的实现可能没问题,因为较大的阶乘是分母,但如果您遇到数字问题,我不会感到惊讶。这里有一些可以帮助您入门的东西。请注意,您需要定义P()a。我不明白“0 到 3”与 c 和 s 有何关系,因为它们的索引范围更远。我假设每个(和增量)都有自己的值文件。

import math
import numpy
c = numpy.getfromtxt("the_c_file.txt")
s = numpy.getfromtxt("the_s_file.txt")
delta = numpy.getfromtxt("the_delta_file.txt")
def V(r, phi, lam):
  ret = 0
  for n in xrange(2, 361):
    for m in xrange(0, n + 1):
      inner = c[n,m]*math.cos(m*lam) + s[n,m]*math.sin(m*lam)
      inner *= math.sqrt(math.factorial(n-m)/math.factorial(n+m)*(2*n+1)*(2-delta[m,0]))
      inner *= P(n, m, math.sin(lam))
    ret += math.pow(a/r, n) * inner
  return ret
Run Code Online (Sandbox Code Playgroud)

确保编写单元测试来检查数学。请注意,“lambda”是保留字。