Sharedshifted_Legendre.ipynbOpen in CoCalc

shifted Legendre多項式(の正規化)を計算する。 shifted Legendre多項式とは、 (p(x),q(x))=01p(x)q(x)dx(p(x),q(x))=\int_0^1p(x)q(x)dx という内積を用いて、(1,x,x2,)(1,x,x^2,\cdots)にGram-Schmidtの直交化を施して得られる直交多項式である。

通常のLegendre多項式は、 (p(x),q(x))=11p(x)q(x)dx(p(x),q(x))=\int_{-1}^1p(x)q(x)dx という内積を用いて、(1,x,x2,)(1,x,x^2,\cdots)にGram-Schmidtの直交化を施して得られる直交多項式である。(積分範囲が違うことに注意。)

P_0(x)=1

(P0(x),P0(x)(P_0(x),P_0(x)を計算する。 (P0(x),P0(x)=01dx=1(P_0(x),P_0(x)=\int_{0}^{1}dx=1 より、

p_0(x)=P_0/sqrt((P_0^2).integral(x,0,1));p_0
x |--> 1

p0(x)=P0(P0,P0)=P0=1p_0(x)=\dfrac{P_0}{(P_0,P_0)}=P_0=1とわかる。 次に、P1(x)=x(x,p0(x))p0P_1(x)=x-(x,p_0(x))p_0を計算しよう。すなわち、

P_1=x-(x*p_0).integral(x,0,1)*p_0;P_1.simplify()
x |--> x - 1/2

次に、(P1(x),P1(x))(P_1(x),P_1(x))を計算する。すなわち、 (P1(x),P1(x)=01P12dx(P_1(x),P_1(x)=\int_{0}^{1}P_1^2dx

expand(P_1^2).simplify()
x |--> x^2 - x + 1/4

より、(P1(x),P1(x)=01P12dx=01(x2x+14)dx=112(P_1(x),P_1(x)=\int_{0}^{1}P_1^2dx=\int_{0}^{1}(x^2 - x + \dfrac{1}{4})dx=\dfrac{1}{12} よって、

p_1(x)=P_1/sqrt((P_1^2).integral(x,0,1));p_1.simplify()
x |--> sqrt(3)*(2*x - 1)

p1(x)=P1(P1,P1)=3(2x1)p_1(x)=\dfrac{P_1}{(P_1,P_1)}=\sqrt{3}(2x-1) となることがわかる。

次に、P2(x)=x2(x2,p0)p0(x2,p1)p1P_2(x)=x^2-(x^2,p_0)p_0-(x^2,p_1)p_1を計算する。すなわち、

P_2(x)=x^2-(x^2*p_0).integral(x,0,1)*p_0-(x^2*p_1).integral(x,0,1)*p_1;P_2.simplify()
x |--> x^2 - x + 1/6

次に、(P2(x),P2(x))(P_2(x),P_2(x))を計算する。 P2(x)2= x42x3+43x213x+136P_2(x)^2=\ x^{4} - 2 \, x^{3} + \frac{4}{3} \, x^{2} - \frac{1}{3} \, x + \frac{1}{36}より、

(P2(x),P2(x))=01P22dx=01(x42x3+43x213x+136)dx=1180=1512232(P_2(x),P_2(x))=\int_{0}^{1}P_2^2dx=\int_{0}^{1}(x^{4} - 2 \, x^{3} + \frac{4}{3} \, x^{2} - \frac{1}{3} \, x + \frac{1}{36})dx=\dfrac{1}{180}=\dfrac{1}{5}\dfrac{1}{2^2 3^2} とわかる。p2(x)=P2(x)/(P2(x),P2(x))p_2(x)=P_2(x)/\sqrt{(P_2(x),P_2(x))}を計算すると、

p_2(x)=P_2/sqrt((P_2^2).integral(x,0,1));p_2.simplify()
x |--> sqrt(5)*(6*x^2 - 6*x + 1)

p2(x)=5(6x26x+1)p_2(x)=\sqrt{5}(6x^2-6x+1)とわかる。次に、 P3(x)=x3(x3,p0(x))p0(x3,p1(x))p1(x3,p2(x))p2P_3(x)=x^3-(x^3,p_0(x))p_0-(x^3,p_1(x))p_1-(x^3,p_2(x))p_2を計算する。すなわち、

P_3(x)=(x^3-(x^3*p_0).integral(x,0,1)*p_0-(x^3*p_1).integral(x,0,1)*p_1-(x^3*p_2).integral(x,0,1)*p_2).simplify()

P3(x)= x332x2+35x120P_3(x)=\ x^{3} - \frac{3}{2} \, x^{2} + \frac{3}{5} \, x - \frac{1}{20}

次に、(P3(x),P3(x))=01P3(x)2dx(P_3(x),P_3(x))=\int_0^1 P_3(x)^2dxを計算する。

expand(P_3^2).simplify()
x |--> x^6 - 3*x^5 + 69/20*x^4 - 19/10*x^3 + 51/100*x^2 - 3/50*x + 1/400
(P_3^2).integral(x,0,1)
1/2800
factor(2800)
2^4 * 5^2 * 7

より、P3(x)2= x63x5+6920x41910x3+51100x2350x+1400P_3(x)^2=\ x^{6} - 3 \, x^{5} + \frac{69}{20} \, x^{4} - \frac{19}{10} \, x^{3} + \frac{51}{100} \, x^{2} - \frac{3}{50} \, x + \frac{1}{400}とわかるので、(P3(x),P3(x))=01P3(x)2dx=01( x63x5+6920x41910x3+51100x2350x+1400)dx=12800=17(1225)2(P_3(x),P_3(x))=\int_0^1 P_3(x)^2dx=\int_0^1 \left(\ x^{6} - 3 \, x^{5} + \frac{69}{20} \, x^{4} - \frac{19}{10} \, x^{3} + \frac{51}{100} \, x^{2} - \frac{3}{50} \, x + \frac{1}{400}\right)dx=\dfrac{1}{2800}=\dfrac{1}{7}\left(\dfrac{1}{2^2 5}\right)^2より、

p_3(x)=P_3/sqrt((P_3^2).integral(x,0,1));latex(p_3.simplify())
x \ {\mapsto}\ \sqrt{7} {\left(20 \, x^{3} - 30 \, x^{2} + 12 \, x - 1\right)}

p3(x)=7(20x330x2+12x1)p_3(x)=\sqrt{7} {\left(20 \, x^{3} - 30 \, x^{2} + 12 \, x - 1\right)}とわかる。