import numpy as np from math import * def Pnm(Phi, Degree): P = np.zeros([Degree + 2, Degree + 2]) # 跨阶次正规化勒让德系数 P[1][1] = 1 P[2][1] = sin(Phi) * 3 ** 0.5 P[2][2] = sqrt(3 * (1 - sin(Phi) ** 2)) for j in range(1, 3): for i in range(3, Degree + 2): l = i - 1 m = j - 1 a = sqrt((4 * l ** 2 - 1) / (l ** 2 - m ** 2)) b = sqrt((2 * l + 1) / (2 * l - 3)) * sqrt(((l - 1) ** 2 - m ** 2) / (l ** 2 - m ** 2)) P[i][j] = a * sin(Phi) * P[i - 1][j] - b * P[i - 2][j] for j in range(3, Degree + 1): for i in range(j, j + 2): l = i - 1 m = j - 1 if (m == 2): beta = sqrt(2 * (2 * l + 1) * (l + m - 2) * (l + m - 3) / (2 * l - 3) / (l + m) / (l + m - 1)) gama = sqrt(2 * (l - m + 1) * (l - m + 2) / (l + m) / (l + m - 1)) else: beta = sqrt((2 * l + 1) * (l + m - 2) * (l + m - 3) / (2 * l - 3) / (l + m) / (l + m - 1)) gama = sqrt((l - m + 1) * (l - m + 2) / (l + m) / (l + m - 1)) P[i][j] = beta * P[i - 2][j - 2] - gama * P[i][j - 2] if ((j + 2) < Degree + 2): for i in range(j + 2, Degree + 2): l = i - 1 m = j - 1 alpha = sqrt((2 * l + 1) * (l - m) * (l - m - 1) / (2 * l - 3) / (l + m) / (l + m - 1)) if (m == 2): beta = sqrt(2 * (2 * l + 1) * (l + m - 2) * (l + m - 3) / (2 * l - 3) / (l + m) / (l + m - 1)) gama = sqrt(2 * (l - m + 1) * (l - m + 2) / (l + m) / (l + m - 1)) else: beta = sqrt((2 * l + 1) * (l + m - 2) * (l + m - 3) / (2 * l - 3) / (l + m) / (l + m - 1)) gama = sqrt((l - m + 1) * (l - m + 2) / (l + m) / (l + m - 1)) P[i][j] = alpha * P[i - 2][j] + beta * P[i - 2][j - 2] - gama * P[i][j - 2] l = Degree m = Degree beta = sqrt((2 * l + 1) * (l + m - 2) * (l + m - 3) / (2 * l - 3) / (l + m) / (l + m - 1)) gama = sqrt((l - m + 1) * (l - m + 2) / (l + m) / (l + m - 1)) P[l + 1][m + 1] = beta * P[l + 1 - 2][m + 1 - 2] - gama * P[l + 1][m + 1 - 2] return P def P_final(theta, n, m, Degree=360): Phi = pi / 2 - theta res = Pnm(Phi, Degree) return res a = P_final(radians(58), 360, 360) print(a)
时间: 2023-02-15 13:32:20 浏览: 136
for n in range(Degree): for m in range(n, Degree): if n == m: P[n, m] = sqrt((2 - 1) / 2 * factorial(n) / (4 * pi * factorial(n))) * cos(m * Phi) else: P[n, m] = sqrt((2 * n + 1) / (2 * n * (n + 1)) *
阅读全文