在Python中使用三维系数数组在点(x,y)处评估二维Legendre级数

介绍

二维Legendre多项式是数学中经常用到的多项式之一，可以用于描述2D平面中的一些特性，如电场，磁场，光学等物理现象。二维Legendre多项式也可以被看做是二维版本的勒让德多项式。在本文中，我们将学习如何使用Python中的三维系数数组在给定点(x, y)处评估二维Legendre级数的值。

二维Legendre多项式

二维Legendre多项式是在单位圆内定义的多项式。给定两个整数l和m， $P_{l,m}(x, y)$ 定义如下：

$P_{l,m}(x,y) = \frac{(-1)^m}{2^ll!}\sqrt{\frac{(l+m)!(l-m)!}{(2l+1)!}}(1-x^2)^{\frac{m}{2}}\frac{\partial^l}{\partial x^l}(y^2-1)^l\frac{\partial^m}{\partial y^m}(x^2-1)^m$

其中， $l$ 和 $m$ 分别是二维Legendre多项式的次数，在给定点 $(x, y)$ 处评估二维Legendre多项式的值也称作“求解(Compute) Legendre多项式”。二维Legendre多项式可以被用于描述2D平面中的一些特性，如电场，磁场，光学等物理现象。

三维系数数组

由于我们需要在给定的点(x, y)处评估二维Legendre多项式的值，我们需要先计算一个系数数组。这个系数数组可以通过递归公式计算得到如下：

$\begin{aligned}P_{0,0}(x,y)&=1\P_{1,0}(x,y)&=x\P_{1,1}(x,y)&=y\(l-m)P_{l,m}(x,y)&=(2l-1)xP_{l-1,m}(x,y)-(l+m-1)P_{l-2,m}(x,y)\\sqrt{\frac{(2l+1)(l-m)!}{2(l+m)!}}P_{l,m}(x,y)&=(2l-1) \cdot \frac{xP_{l-1,m}(x,y)-(l+m-1)P_{l-2,m}(x,y)}{l-m}\end{aligned}$

我们可以看到，系数数组可以被递归地计算，只要我们有了初始的值 $P_{0,0}(x,y),P_{1,0}(x,y),P_{1,1}(x,y)$ 就可以利用递归公式计算所有的系数了。根据递归公式，我们可以使用一个三维数组来存储这些系数。在这个数组中，每一个元素都对应一个(l, m)的二维Legendre多项式。

下面就是如何在Python中使用递归公式生成该三维系数数组的代码，判断语言为Python：

import numpy as np

def legendre_coeffs(n):
    """
    Computes the coefficients for Legendre expansion in a 3D array. 
    Returns:
        p : ndarray [l, m, k]
           3D array containing polynomial coefficients for degree l and order m at point (x, y).
    """
    p = np.zeros((n+1, n+1, n+1))
    p[0, 0, 0] = 1.0
    p[1, 0, 0] = 1.0
    p[1, 1, 0] = 0.0
    p[1, 1, 1] = 1.0

    for l in range(2, n+1):
        for m in range(l+1):
            for k in range(n-l+1):
                t1 = (2*l-1)/(l-m)
                p[l, m, k+l] = t1*x*p[l-1, m, k+l-1] - ((l+m-1)*t1)*p[l-2, m, k+l-1]

                t2 = np.sqrt((2*l+1)*np.math.factorial(l-m)/(2*np.math.factorial(l+m)))
                p[l, m, k+l] *= t2

    return p

这个函数接受一个参数 $n$ ，表示我们需要构造一个 $n \times n \times n$ 的三维数组。例如，想要一个20次的Legendre二维多项式，我们可以这样调用这个函数：

p = legendre_coeffs(20)

这个函数会返回一个20x20x20的三维数组 $p$ ，其中的每个元素 $p_{l,m, k+l}$ 表示次数为 $l$ ，序数为 $m$ 在 $x^2 + y^2 = k^2$ 的条件下计算的多项式系数。我们将在下一节中使用这个三维数组来评估二维Legendre多项式。

在点(x,y)处评估二维Legendre级数

现在我们已经有了一个包含所有系数的三维数组，我们可以使用这个数组在椭圆上评估二维Legendre级数。首先，我们需要计算x和y的值：

x = 0.5
y = 0.5

然后，我们还需要指定一个特定的半径（ $r$ ）。在这个半径上，我们将为确保二维Legendre多项式比在其他点更可靠。例如，让我们选择半径为10：

r = 10

接下来，我们将使用上一节生成的系数数组来计算二维Legendre级数的值。我们将使用以下代码：

value = 0.0
for l in range(p.shape[0]):
    for m in range(l+1):
        t1 = p[l, m, int(np.round(np.sqrt(r**2 - x**2 - y**2)))]
        t2 = t1 * ((-1)**m * np.sqrt((2*l+1)*np.math.factorial(l-m)/(4*np.pi*np.math.factorial(l+m))))
        t3 = t2 * x**l * y**m
        value += t3

这个代码使用了上一节中的三维系数数组 $p$ 和我们指定的点 $(x, y)$ 和半径 $r$ 。我们将计算循环中的每个 $(l, m)$ 号系数，这将生成一个用于评估二维Legendre级数的项。我们将从 $l=0$ 开始计算，并在每次循环中递增 $l$ 和 $m$ 。t1计算了 $p_{l,m,k+l}$ ，其中 $k$ 是满足 $x^2 + y^2 + k^2 = r^2$ 的最大整数。t2将观测值 $t1$ 乘以系数，并t2将使用上述公式计算。最后，t3将使用观测值添加到之前的值中。