% 表題  動径座標のスペクトル法
%
% 履歴 2008/01/23 竹広真一 
%      2008/03/26 竹広真一 
%      2008/07/20 竹広真一 
%
%
\documentclass[a4j,12pt]{jarticle}
\usepackage{Dennou6}

\Dparskip

\Dtitle{動径座標のスペクトル法}

\Dauthor{竹広 真一}
\Ddate{2008 年 7 月 20 日}

\begin{document}

\maketitle

2 次元極座標, 3 次元円筒座標, および 3 次元球座標での
原点で $C^\infty$ 級の微分可能な滑らかな関数を
スペクトル法により表現するための動径座標の展開関数について述べる. 

\section{直交多項式系}

  円筒座標 $(r,\phi)$ での $C^\infty$ 級の関数 $f(r,\phi)$ を
  表現するスペクトル法をさがす. 方位角方向に Fourier 展開して
  %
  \begin{equation}
    f(r,\phi) = \sum_{m=-\infty}^\infty f_m(r)e^{im\phi}. 
  \end{equation}
  %
  このとき $f_m(r)$ は極での微分可能条件から
  $r\rightarrow 0$ で $O(r^{|m|+2p})$ と振る舞わなければならない
  ことがわかる
  \footnote{
    Matsushima and Marcus の参考文献[6] を参照すべし. 
    球座標では
    \begin{displaymath}
      g(\lambda,\varphi, r) 
      = \sum_{n=0}^\infty\sum_{m=-n}^n g_n^m(r)Y_n^m(\lambda,\varphi)
    \end{displaymath}
    なる関数について $g_n^m(r)$ は
    $r\rightarrow 0$ で $O(r^{|n|+2p})$ と振る舞わなければならないはず. 
  }. 
  ここで $p$ は非負の整数である. 

  さて, 極での条件を満たしつつ, 
  展開関数が境界で Gibbs の現象が生じないよう
  定義となる微分方程式が境界で特異であり
  \footnote{Matsushima and Marcus (1995) の参考文献11 を参照すべし}, 
  かつ重み関数が円筒および球の幾何学的形状に都合のよい
  $f_m(r)$ の表現を求めたい. 
  また, 3 項間の漸化式を満たし取りあつかいが簡単であるので
  直交多項式が望ましい. 

  この 3 つの要請を満たす定義となる微分方程式として, 
  $0\ge r \ge 1$ で定義される
  特異な Strum-Liouville 型方程式
  %
  \begin{equation}
    \Deqlab{微分方程式}
    \frac{(1-r^2)^{1-\alpha}}{r^\beta}
        \DD{}{r}\left((1-r^2)^\alpha r^\beta \DD{y}{r}\right)
   -\frac{|m|(|m|+\beta-1)}{r^2}y 
   + n(n+2\alpha+\beta-1) y = 0, 
  \end{equation}
  %
  ただし $n,m$ は整数で $0 \le |m| \le n$ であり, 
  $0<\alpha \le 1$, $\beta$ は正の整数である. 

  式\Deqref{微分方程式}は $n+m$ が偶数のとき $n$ 次の多項式の解を持ち
  $r=0$ に関するフロベニウスの級数を用いて解 $y=Q_n^m(\alpha,\beta;r)$
  を閉じた形で書くことができる.
  %
  \begin{equation}
    \Deqlab{級数解}
    Q_n^m(\alpha,\beta;r) = \sum_{p=0}^{(n-|m|)/2}
    \frac{\Ddsty (-1)^{p+(n-|m|)/2}
                 \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                 \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
         {\Ddsty p! \left(\frac{n-|m|}{2}-p\right)!
                 \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                 \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}
         r^{|m|+2p},
  \end{equation}
  %
  ただし $\gamma\equiv 2\alpha+\beta$ である
  \footnote{
    式\Deqref{微分方程式}の解を $r$ の $q$ 次から始まる級数 
    $y=\sum_{l=0}^\infty a_l r^{q+l}$ の形で求める. 
    これを式\Deqref{微分方程式}に代入すると
    %
    \begin{displaymath}
        \frac{(1-r^2)^{1-\alpha}}{r^\beta}
            \DD{}{r}\left((1-r^2)^\alpha r^\beta \DD{}{r}
              \sum_{l=0}^\infty a_l r^{q+l}\right)
       -|m|(|m|+\beta-1)\sum_{l=0}^\infty a_l r^{q+l-2}
       + n(n+2\alpha+\beta-1) \sum_{l=0}^\infty a_l r^{q+l} = 0, 
    \end{displaymath}
    %
    第 1 項目を整理すると
    \begin{eqnarray*}
      & & \frac{(1-r^2)^{1-\alpha}}{r^\beta}
            \DD{}{r}\left((1-r^2)^\alpha r^\beta \DD{}{r}
              \sum_{l=0}^\infty a_l r^{q+l}\right)
         =\frac{(1-r^2)^{1-\alpha}}{r^\beta}
            \DD{}{r}(1-r^2)^\alpha r^\beta
              \sum_{l=0}^\infty (q+l)a_l r^{q+l-1}\\
      &=& (1-r^2)\sum_{l=0}^\infty (q+l)a_l \DD{}{r} r^{q+l-1}
         + \frac{(1-r^2)^{1-\alpha}}{r^\beta}
            (-2\alpha r)(1-r^2)^{\alpha-1}r^\beta
            \sum_{l=0}^\infty (q+l)a_l r^{q+l-1} \\
      & & + \frac{(1-r^2)^{1-\alpha}}{r^\beta}
            (1-r^2)^{\alpha}\beta r^{\beta-1}
            \sum_{l=0}^\infty (q+l)a_l r^{q+l-1}\\
      &=& (1-r^2)\sum_{l=0}^\infty (q+l)(q+l-1) a_l r^{q+l-2}
         + (-2\alpha r)\sum_{l=0}^\infty (q+l)a_l r^{q+l-1}
         + \frac{(1-r^2)}{r}\beta 
           \sum_{l=0}^\infty (q+l)a_l r^{q+l-1}\\
      &=& (1-r^2)\sum_{l=0}^\infty (q+l)(q+l-1) a_l r^{q+l-2}
         - \sum_{l=0}^\infty 2\alpha (q+l)a_l r^{q+l}
         + (1-r^2)\sum_{l=0}^\infty \beta (q+l)a_l r^{q+l-2}\\
      &=& (1-r^2)\sum_{l=0}^\infty (q+l)(q+l+\beta-1) a_l r^{q+l-2}
         - \sum_{l=0}^\infty 2\alpha (q+l)a_l r^{q+l}\\
      &=& \sum_{l=0}^\infty (q+l)(q+l+\beta-1) a_l r^{q+l-2}
         -\sum_{l=0}^\infty (q+l)(q+l+\beta-1) a_l r^{q+l}
         - \sum_{l=0}^\infty 2\alpha (q+l)a_l r^{q+l}\\
      &=& \sum_{l=0}^\infty (q+l)(q+l+\beta-1) a_l r^{q+l-2}
         -\sum_{l=0}^\infty (q+l)(q+l+2\alpha+\beta-1) a_l r^{q+l}
    \end{eqnarray*}
    %
    これを元の式に代入して整理すると
    \begin{eqnarray*}
      & &
         \sum_{l=0}^\infty (q+l)(q+l+\beta-1) a_l r^{q+l-2}
           -\sum_{l=0}^\infty (q+l)(q+l+2\alpha+\beta-1) a_l r^{q+l}\\
      & & -|m|(|m|+\beta-1)\sum_{l=0}^\infty a_l r^{q+l-2}
          + n(n+2\alpha+\beta-1) \sum_{l=0}^\infty a_l r^{q+l} = 0, \\
      & &
         \sum_{l=0}^\infty [(q+l)(q+l+\beta-1) -|m|(|m|+\beta-1)]a_lr^{q+l-2}\\
      & & + \sum_{l=0}^\infty [n(n+2\alpha+\beta-1)-(q+l)(q+l+2\alpha+\beta-1)]
            a_l r^{q+l} = 0, \\
      & &
         \sum_{l=0}^\infty(q+l-|m|)(q+l+|m|+\beta-1)a_l r^{q+l-2}
       + \sum_{l=0}^\infty (n-q-l)(q+l+n+2\alpha+\beta-1)a_l r^{q+l} = 0, \\
    \end{eqnarray*}
    %
    最低次の項から次数 $q$ が定まる. 第 1 項目の $l=0$ から
    %
    \begin{displaymath}
      (q-|m|)(q+|m|+\beta-1)a_0 = 0, \quad i.e. \quad q = |m|.
    \end{displaymath}
    %
    $l=1$ からは
    %
    \begin{displaymath}
      (|m|+1-|m|)(|m|+|m|+\beta-1)a_1 = 0, \quad i.e. \quad a_1 = 0.
    \end{displaymath}
    %
    $q=|m|$ を代入して $a_l$ の漸化式を導こう. 第 1 項目の $l=0,1$ は
    0 となっているので
    %
    \begin{eqnarray*}
      & &
         \sum_{l=2}^\infty l(l+2|m|+\beta-1)a_l r^{|m|+l-2}
       + \sum_{l=0}^\infty (n-|m|-l)(l+n+|m|+ \gamma -1)a_l r^{q+l} = 0,
       \\
      & &
         \sum_{l=0}^\infty (l+2)(l+2|m|+\beta+1)a_{l+2} r^{|m|+l}
       + \sum_{l=0}^\infty (n-|m|-l)(l+n+|m|+ \gamma -1)a_l r^{|m|+l} = 0, \\
    \end{eqnarray*}
    %
    ただし $\gamma=2\alpha+\beta$ である. したがって
    %
    \begin{displaymath}
      a_{l+2} = -\frac{(n-|m|-l)(l+n+|m|+ \gamma -1)}{(l+2)(l+2|m|+\beta+1)}a_l
    \end{displaymath}
    %
    $a_1=0$ より $l$ が奇数の項は 0 である. そこで $l=2p-2$ と置き換えて
    %
    \begin{displaymath}
      a_{2p} = -\frac{(n-|m|-2p+2)(2p+n+|m|+ \gamma -3)}
                     {2p(2p+2|m|+\beta-1)}a_{2p-2}
    \end{displaymath}
    %
    これより $n-|m|$ が偶数ならば $p=(n-|m|)/2+1$ 番目の項は分子の
    一つめのカッコの中味が 0 となるので $p=(n-|m|)/2$ までで級数が閉じる. 
    さらに
    %
    \begin{eqnarray*}
      a_{2p} &=& -\frac{(n-|m|-2p+2)(2p+n+|m|+ \gamma -3)}
                     {2p(2p+2|m|+\beta-1)}a_{2p-2}\\
             &=& -\frac{\Ddsty \left(\frac{n-|m|}{2}-p+1\right)
                               \left(\frac{n+|m|+ \gamma-3}{2}+p\right)}
                       {\Ddsty p\left(\frac{2|m|+\beta-1}{2}+p\right)}a_{2p-2}\\
             &=& (-1)^2\frac{\Ddsty \left(\frac{n-|m|}{2}-p+1\right)
                               \left(\frac{n-|m|}{2}-p+2\right)
                               \left(\frac{n+|m|+ \gamma-3}{2}+p\right)
                               \left(\frac{n+|m|+ \gamma-3}{2}+p-1\right)}
                       {\Ddsty p(p-1)
                               \left(\frac{2|m|+\beta-1}{2}+p\right)
                               \left(\frac{2|m|+\beta-1}{2}+p-1\right)}a_{2p-4}\\
             &=& (-1)^{p}
                       \frac{\Ddsty \left(\frac{n-|m|}{2}-p+1\right)
                               \left(\frac{n-|m|}{2}-p+2\right)
                               \cdots
                               \left(\frac{n-|m|}{2}\right)}
                       {\Ddsty p(p-1)\cdots 1}\\
             & &       \frac{\Ddsty \left(\frac{n+|m|+ \gamma-3}{2}+p\right)
                               \left(\frac{n+|m|+ \gamma-3}{2}+p-1\right)
                               \cdots
                               \left(\frac{n+|m|+ \gamma-3}{2}+1\right)}
                            {\Ddsty
                               \left(\frac{2|m|+\beta-1}{2}+p\right)
                               \left(\frac{2|m|+\beta-1}{2}+p-1\right)
                               \cdots
                               \left(\frac{2|m|+\beta-1}{2}+1\right)}a_{0}\\
              &=& (-1)^{p}
                        \frac{\Ddsty \left(\frac{n-|m|}{2}\right)!
                                \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                                \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
                             {\Ddsty p! \left(\frac{n-|m|}{2}-p\right)!
                                \Gamma\left(\frac{n+|m|+ \gamma-1}{2}\right)
                                \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)}a_{0}
    \end{eqnarray*}
    %
    ここで $a_0$ として
    \begin{displaymath}
      a_0 = \frac{\Ddsty (-1)^{(n-|m|)/2}
                         \Gamma\left(\frac{n+|m|+ \gamma-1}{2}\right)}
                 {\Ddsty \left(\frac{n-|m|}{2}\right)! 
                         \Gamma\left(\frac{2|m|+ \gamma-1}{2}\right)}
    \end{displaymath}
    %
    を選ぶと
    \begin{displaymath}
      a_{2p} = (-1)^{p+(n-|m|)/2}
                \frac{\Ddsty  \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                              \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
                     {\Ddsty p! \left(\frac{n-|m|}{2}-p\right)!
                                \Gamma\left(\frac{2|m|+ \gamma-1}{2}\right)
                                \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)}
    \end{displaymath}
  }. 
  $m$ が偶数であれば $Q_n^m(\alpha,\beta;r)$ は $r$ の偶関数, 
  奇数であれば奇関数である. 
  $r\rightarrow 0$ で  $Q_n^m(\alpha,\beta;r)$ は $O(r^{|m|})$ で
  振る舞うので極座標 $(r,\phi)$ において $Q_n^m(\alpha,\beta;r)e^{im\phi}$
  は極の条件を正確に満たすことになる. 
  $Q_n^m(\alpha,\beta;r)$ は完全系をなし, 重み関数
  %
  \begin{equation}
    \Deqlab{重み関数}
    w(\alpha,\beta;r) = \frac{r^\beta}{(1-r^2)^{1-\alpha}}
  \end{equation}
  %
  に関して直交している. すなわち
  \footnote{
    式 \Deqref{微分方程式} に $Q_{n'}^m(\alpha,\beta;r)w(\alpha,\beta;r)$ をかけて
    領域積分すると
    %
    \begin{displaymath}
      \int_0^1 Q_{n'}^m \DD{}{r}\left((1-r^2)^\alpha r^\beta 
               \DD{Q_{n}^m}{r}\right) dr
     -\int_0^1 \frac{|m|(|m|+\beta-1)}{r^2}Q_{n'}^m Q_{n}^m w dr
     + n(n+2\alpha+\beta-1)\int_0^1 Q_{n'}^m Q_{n}^m w dr = 0, 
    \end{displaymath}
    %
    $\alpha>0$ に注意して第 1 項目を 2 度部分積分して入れ換えると, 
    \begin{eqnarray*}
      & & \int_0^1 Q_{n'}^m \DD{}{r}\left((1-r^2)^\alpha r^\beta 
               \DD{Q_{n}^m}{r}\right) dr
         = \left[Q_{n'}^m(1-r^2)^\alpha r^\beta \DD{Q_{n}^m}{r}\right]_0^1
         - \int_0^1 (1-r^2)^\alpha r^\beta \DD{Q_{n'}^m}{r}\DD{Q_{n}^m}{r}\\
      &=& -\left[Q_{n}^m (1-r^2)^\alpha r^\beta \DD{Q_{n'}^m}{r}\right]_0^1
         + \int_0^1 Q_{n}^m \DD{}{r}\left((1-r^2)^\alpha r^\beta 
               \DD{Q_{n'}^m}{r}\right) dr
          =\int_0^1 Q_{n}^m \DD{}{r}\left((1-r^2)^\alpha r^\beta 
               \DD{Q_{n'}^m}{r}\right) dr
    \end{eqnarray*}
    %
    ここで $Q_{n'}^m$ の満たす微分方程式
    %
    \begin{displaymath}
      \frac{(1-r^2)^{1-\alpha}}{r^\beta}
        \DD{}{r}\left((1-r^2)^\alpha r^\beta \DD{Q_{n'}^m}{r}\right)
        -\frac{|m|(|m|+\beta-1)}{r^2}Q_{n'}^m 
        + n'(n'+2\alpha+\beta-1)Q_{n'}^m = 0, 
    \end{displaymath}
    %
    を用いると, 
    \begin{eqnarray*}
        \int_0^1 Q_{n}^m \DD{}{r}\left((1-r^2)^\alpha r^\beta 
               \DD{Q_{n'}^m}{r}\right) dr
       &=& \int_0^1 \frac{|m|(|m|+\beta-1)}{r^2}Q_{n'}^m Q_{n}^m w dr\\
       & &  - n'(n'+2\alpha+\beta-1)\int_0^1 Q_{n'}^m Q_{n}^m w dr,
    \end{eqnarray*}
    %
    これを最初の式に代入して整理すると $|m|$ が係数に含まれる項は
    キャンセルして
    \begin{displaymath}
       [n'(n'+2\alpha+\beta-1)- n(n+2\alpha+\beta-1)]
       \int_0^1 Q_{n'}^m Q_{n}^m w dr = 0. 
    \end{displaymath}
    %
    したがって $n\neq n'$ ならば $\Ddsty\int_0^1 Q_{n'}^m Q_{n}^m w dr = 0$
    となる. 
  }
  %
  \begin{equation}
    \Deqlab{直交関係}
    \int_0^1 Q_n^m(\alpha,\beta;r) Q_{n'}^m(\alpha,\beta;r)
             w(\alpha,\beta;r) dr
   = I_n^m(\alpha,\beta)\delta_{n,n'},
  \end{equation}
  %
  積分定数 $I_n^m(\alpha,\beta)$ の漸化式が後に示される. 
  与えられた $I_n^m(\alpha,\beta)$ を用いて直交関数系を規格化して
  %
  \begin{equation}
    \Phi_n^m(\alpha,\beta;r) 
    = (I_n^m(\alpha,\beta))^{-1/2} Q_n^m(\alpha,\beta;r)
  \end{equation}
  %
  とする. $\Phi_n^m(\alpha,\beta;r)$ は極座標で与えられる単位円盤 
  $D \equiv \{(r,\phi)| 0 \le r \le 1, 0 \le \phi \le 2\pi\}$ 上での
  スカラー関数を展開するのに用いることができて
  %
  \begin{equation}
    \Deqlab{極座標スペクトル展開}
    f(r,\phi) = \sum_{n=0}^\infty \sum_{m=-n,m+n\ even}^n
                  f_n^m \Phi_n^m(\alpha,\beta;r) e^{im\phi},
  \end{equation}
  %
  ここで $f_n^m$ は複素展開係数である. 
  $\Phi_n^m(\alpha,\beta;r)$ の極での正しい振舞により
  $f(r,\phi)$ は円盤 $D$ 上で $C^\infty$ 級であることに注意されたい. 
  係数 $f_n^m$ は
  %
  \begin{equation}
    \Deqlab{極座標スペクトル係数計算}
    f_n^m = \frac{1}{2\pi}\int_0^{2\pi}\int_0^1
             f(r,\phi)\Phi_n^m(\alpha,\beta;r)e^{-im\phi}w(\alpha,\beta;r) 
             dr d\phi, 
  \end{equation}
  %
  により求められる. $C^\infty$ 級に対する展開式 (7) が
  スペクトル収束(すなわち代数的収束より速い)することは, 
  微分方程式 \Deqref{微分方程式} と (8) を用いることで示すことができる. 
  $r$ に対して部分積分をすることにより
  \footnote{
    \Deqref{微分方程式} を用いると
    \begin{displaymath}
      \Phi_n^m(\alpha,\beta;r) 
        = \frac{1}{n(n+2\alpha+\beta -1)}\left[
               - \frac{(1-r^2)^{1-\alpha}}{r^\beta}
                  \DD{}{r}\left((1-r^2)^\alpha r^\beta 
                    \DD{\Phi_n^m}{r}\right)
                +\frac{|m|(|m|+\beta-1)}{r^2}\Phi_n^m. 
                  \right]
    \end{displaymath}
    %
    これを(8)に代入して部分積分を 2 度行うと (9) を得る. 
  }, 
  %
  \begin{equation}
    f_n^m = \frac{1}{2\pi n(n+2\alpha+\beta -1)}\int_0^{2\pi}\int_0^1
             h(r,\phi)\Phi_n^m(\alpha,\beta;r)e^{-im\phi}w(\alpha,\beta;r) 
             dr d\phi, 
  \end{equation}
  %
  ただし
  %
  \begin{equation}
    h(r,\phi) = - \frac{(1-r^2)^{1-\alpha}}{r^\beta}
                  \DD{}{r}\left((1-r^2)^\alpha r^\beta \DD{f}{r}\right)
                +\frac{|m|(|m|+\beta-1)}{r^2}f. 
  \end{equation}
  %
  $f$ に $r^{|m|} e^{im\phi}$ を代入することにより, 
  $f(r,\phi)$ が極での条件を満たせば $h(r,\phi)$ もまた極での条件を満たす
  ことがわかる
  \footnote{
    実際に $r^{|m|}$ を代入すると
    \begin{eqnarray*}
      h(r,\phi) &=& - \frac{(1-r^2)^{1-\alpha}}{r^\beta}
                       \DD{}{r}\left((1-r^2)^\alpha |m|r^{|m|+\beta-1}\right)
                    +  |m|(|m|+\beta-1)r^{|m|-2}\\
                &=& - (1-r^2)|m|(|m|+\beta-1)r^{|m|-2}
                    - (1-r^2)^{1-\alpha} (1-r^2)^{\alpha-1}
                      (-2\alpha r)|m|r^{|m|-1}\\
                & &  + |m|(|m|+\beta-1)r^{|m|-2}\\
                &=& |m|(|m|+\beta-1)r^{|m|} + 2\alpha |m| r^{|m|}
                  = |m|(|m|+2\alpha+\beta-1)r^{|m|}
    \end{eqnarray*}
    %
    したがって $h(r,\phi)$ も原点で $r^{|m|}$ の振舞をする. 
  }. 
  部分積分を繰り返すことにより式(7)の $f(r,\phi)$ へのスペクトル収束を
  示すことができる
  \footnote{
    $\Phi_n^m(\alpha,\beta;r)$ を \Deqref{微分方程式} で変形して
    部分積分を $2l$ 回繰り返すと, 境界 $r=0,1$ で 
    $(1-r^2)^\alpha r^\beta=0$ となるので
    積分の寄与だけが残り, 
    \begin{displaymath}
      f_n^m = \frac{1}{2\pi n^l(n+2\alpha+\beta-1)^l}\int_0^{2\pi}\int_0^1
               {\cal D}_l f(r,\phi)
                   \Phi_n^m(\alpha,\beta;r)e^{-im\phi}w(\alpha,\beta;r) 
             dr d\phi, 
    \end{displaymath}
    %
    ただし
    %
    \begin{displaymath}
      {\cal D}_l = \left[- \frac{(1-r^2)^{1-\alpha}nh}{r^\beta}
                           \DD{}{r}\left((1-r^2)^\alpha r^\beta \DD{}{r}\right)
                        +\frac{|m|(|m|+\beta-1)}{r^2}\right]^l 
    \end{displaymath}
    %
    である. $f(r,\phi)$ が円盤上で $C^{\infty}$ 級であり, 
    さらに先の一つ目の式変形から ${\cal D}_l f$ は極での条件を満たしている.
    したがって ${\cal D}_l f$ も $C^{\infty}$ 級であるので
    $n$ によらず有界な値で積分値が制限される. 
    したがって $f_n^m$ はどのような代数的収束 $1/n^{2l}$ よりも速く収束する.
 }. 


\section{数値計算の手順}

  ここでは実際に $\Phi_n^m(\alpha,\beta;r)$ と 式(7) の $f_n^m$ を
  数値的に見積もるための手順を議論する. 
  簡単のため以下では $Q_n^m(\alpha,\beta;r)$, 
  $\Phi(\alpha,\beta;r)$, $I_n^m(\alpha,\beta)$, $w(\alpha,\beta;r)$
  の引数 $\alpha,\beta$ を省略する. 
  まずいくつかの有用な関係式を導出する.  $n>|m|$ に対して $Q_n^m$ は
  次の関係式を満たす. 
  %
  \begin{equation}
    \Deqlab{微分漸化式1}
    r\DD{}{r}(Q_n^m(r)- Q_{n-2}^m(r)) 
    = n Q_n^m(r)+ (n+\gamma-3) Q_{n-2}^m(r).
  \end{equation}
  %
  ここで $Q_{|m|-2}^m(r)\equiv 0$ である. 
  式\Deqref{微分漸化式1}は\Deqref{級数解}を
  代入することで容易に確かめられる
  \footnote{
    \begin{eqnarray*}
      & &
      r\DD{}{r}(Q_n^m(r)- Q_{n-2}^m(r)) \\
      & & = \sum_{p=0}^{(n-|m|)/2}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2}(|m|+2p)
                       \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &- \sum_{p=0}^{(n-|m|)/2-1}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2-1}(|m|+2p)
                       \Gamma\left(\frac{n+|m|+\gamma-3}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-1-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &
      = \sum_{p=0}^{(n-|m|)/2}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2}[n+(-n+|m|+2p)]
                       \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &- \sum_{p=0}^{(n-|m|)/2-1}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2-1}(|m|+2p)
                       \Gamma\left(\frac{n+|m|+\gamma-3}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-1-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &                    
      = n\sum_{p=0}^{(n-|m|)/2}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2}
                       \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &- \sum_{p=0}^{(n-|m|)/2}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2}(n-|m|-2p)
                       \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &- \sum_{p=0}^{(n-|m|)/2-1}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2-1}(|m|+2p)
                       \Gamma\left(\frac{n+|m|+\gamma-3}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-1-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &
      = nQ_n^m \\
      & &- \sum_{p=0}^{(n-|m|)/2-1}
         \frac{\Ddsty (-1)^{p+(n-|m|)/2}(n-|m|-2p)
                       \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &- \sum_{p=0}^{(n-|m|)/2-1}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2-1}(|m|+2p)
                       \Gamma\left(\frac{n+|m|+\gamma-3}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-1-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &
       = nQ_n^m\\
      & & -\sum_{p=0}^{(n-|m|)/2-1}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2}2\left(\frac{n+|m|+\gamma-3}{2}+p\right)
                       \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-1-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &- \sum_{p=0}^{(n-|m|)/2-1}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2-1}(|m|+2p)
                       \Gamma\left(\frac{n+|m|+\gamma-3}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-1-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &
       = nQ_n^m \\
      & &+ \sum_{p=0}^{(n-|m|)/2-1}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2-1}\left(n+\gamma-3+|m|+2p\right)
                       \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-1-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &- \sum_{p=0}^{(n-|m|)/2-1}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2-1}(|m|+2p)
                       \Gamma\left(\frac{n+|m|+\gamma-3}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-1-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &
       = nQ_n^m+ \sum_{p=0}^{(n-|m|)/2-1}
          \frac{\Ddsty (-1)^{p+(n-|m|)/2-1}(n+\gamma-3)
                       \Gamma\left(\frac{n+|m|+\gamma-1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
               {\Ddsty p! \left(\frac{n-|m|}{2}-1-p\right)!
                       \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                       \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}r^{|m|+2p}\\
      & &
       = nQ_n^m+ (n+\gamma-3)Q_{n-2}^m
    \end{eqnarray*}
  }. 
  式\Deqref{微分漸化式1}を用いて式\Deqref{微分方程式}の
  階数を一つ下げることができて
  \footnote{
    $Q_n^m,Q_{n-2}^m$ に対する微分方程式は\Deqref{微分方程式}より
    \begin{eqnarray*}
      & &
      \frac{(1-r^2)^{1-\alpha}}{r^\beta}
          \DD{}{r}\left((1-r^2)^\alpha r^\beta \DD{Q_n^m}{r}\right)
     -\frac{|m|(|m|+\beta-1)}{r^2}Q_n^m 
     + n(n+2\alpha+\beta-1)Q_n^m  = 0, \\
     & &
      \frac{(1-r^2)^{1-\alpha}}{r^\beta}
          \DD{}{r}\left((1-r^2)^\alpha r^\beta \DD{Q_{n-2}^m}{r}\right)
     -\frac{|m|(|m|+\beta-1)}{r^2}Q_{n-2}^m 
     + (n-2)(n+2\alpha+\beta-3) Q_{n-2}^m = 0,
    \end{eqnarray*}
    %
    2 つの式を引き算すると
    \begin{eqnarray*}
      & &
      \frac{(1-r^2)^{1-\alpha}}{r^\beta}
          \DD{}{r}\left[(1-r^2)^\alpha r^\beta \DD{}{r}(Q_n^m-Q_{n-2}^m)\right]
      -\frac{|m|(|m|+\beta-1)}{r^2}(Q_n^m - Q_{n-2}^m )\\
      & &
      + n(n+2\alpha+\beta-1)Q_n^m 
      - (n-2)(n+2\alpha+\beta-3) Q_{n-2}^m = 0.
    \end{eqnarray*}
    %
    第 1 項に\Deqref{微分漸化式1}を適用すると
    \begin{eqnarray*}
      & &
      \frac{(1-r^2)^{1-\alpha}}{r^\beta}
          \DD{}{r}\left[(1-r^2)^\alpha r^{\beta-1} 
            (n Q_n^m+ (n+\gamma-3) Q_{n-2}^m)\right]
      -\frac{|m|(|m|+\beta-1)}{r^2}(Q_n^m - Q_{n-2}^m )\\
      & &
      + n(n+2\alpha+\beta-1)Q_n^m 
      - (n-2)(n+2\alpha+\beta-3) Q_{n-2}^m = 0,
      \\
      & &
        \frac{(1-r^2)}{r}n\DD{Q_n^m}{r}
      + \frac{(1-r^2)}{r}(n+\gamma-3)\DD{Q_{n-2}^m}{r}\\
      & &
      + \frac{(1-r^2)^{1-\alpha}}{r^\beta}
          [(1-r^2)^{\alpha-1}(-2\alpha r)r^{\beta-1} 
           +(1-r^2)^\alpha (\beta-1)r^{\beta-2}]
            [n Q_n^m+ (n+\gamma-3) Q_{n-2}^m]\\
      & &
      -\frac{|m|(|m|+\beta-1)}{r^2}(Q_n^m - Q_{n-2}^m )
      + n(n+2\alpha+\beta-1)Q_n^m 
      - (n-2)(n+2\alpha+\beta-3) Q_{n-2}^m = 0,
      \\
      & &
        \frac{(1-r^2)}{r}n\DD{Q_n^m}{r}
      + \frac{(1-r^2)}{r}(n+\gamma-3)\DD{Q_{n-2}^m}{r}\\
      & &
      + \left(-2\alpha +(\beta-1)\frac{(1-r^2)}{r^2}\right)
            [n Q_n^m+ (n+\gamma-3) Q_{n-2}^m]\\
      & &
      -\frac{|m|(|m|+\beta-1)}{r^2}(Q_n^m - Q_{n-2}^m )
      + n(n+2\alpha+\beta-1)Q_n^m 
      - (n-2)(n+2\alpha+\beta-3) Q_{n-2}^m = 0,
    \end{eqnarray*}
    %
    \Deqref{微分漸化式1}を再度用いて $\Ddsty \DD{Q_{n-2}^m}{r}$ を
    消去すると
    %
    \begin{eqnarray*}
      & &
        \frac{(1-r^2)}{r}n\DD{Q_n^m}{r}
      + \frac{(1-r^2)}{r^2}(n+\gamma-3)
        \left\{r\DD{Q_n^m}{r} - [n Q_n^m+ (n+\gamma-3) Q_{n-2}^m]\right\}\\
      & &
      + \left(-2\alpha +(\beta-1)\frac{(1-r^2)}{r^2}\right)
            [n Q_n^m+ (n+\gamma-3) Q_{n-2}^m]\\
      & &
      -\frac{|m|(|m|+\beta-1)}{r^2}(Q_n^m - Q_{n-2}^m )
      + n(n+2\alpha+\beta-1)Q_n^m 
      - (n-2)(n+2\alpha+\beta-3) Q_{n-2}^m = 0,
    \end{eqnarray*}
    %
    $\Ddsty \DD{Q_n^m}{r}$ の項の係数は
    %
    \begin{displaymath}
       \frac{(1-r^2)}{r}n + \frac{(1-r^2)}{r^2}(n+\gamma-3)r
      =\frac{(1-r^2)}{r}(2n+\gamma-3).
    \end{displaymath}
    %
    $\Ddsty Q_n^m$ の項の係数は
    %
    \begin{eqnarray*}
      & & 
          -\frac{(1-r^2)}{r^2}(n+\gamma-3)n 
          +\left(-2\alpha +(\beta-1)\frac{(1-r^2)}{r^2}\right)n
          -\frac{|m|(|m|+\beta-1)}{r^2}
          + n(n+2\alpha+\beta-1)\\
      &=&  n\left[-\frac{(1-r^2)}{r^2}(n+\gamma-3)
                  +\left(-2\alpha +(\beta-1)\frac{(1-r^2)}{r^2}\right)
                  +(n+2\alpha+\beta-1)\right]
            -\frac{|m|(|m|+\beta-1)}{r^2}\\
      &=&   n\left[-\frac{(1-r^2)}{r^2}(n+\gamma-\beta-2)
                  +(n+\beta-1)\right]
            -\frac{|m|(|m|+\beta-1)}{r^2}\\
      &=&   n\left[-\frac{(n+\gamma-\beta-2)}{r^2}
                  +(n+\gamma-\beta-2+n+\beta-1)\right]
            -\frac{|m|(|m|+\beta-1)}{r^2}\\
      &=&   n\left[-\frac{(n+\gamma-\beta-2)}{r^2}
                  +(2n+\gamma-3)\right]
            -\frac{|m|(|m|+\beta-1)}{r^2}\\
      &=&   n(2n+\gamma-3)
            -\frac{|m|(|m|+\beta-1)+n(n+\gamma-\beta-2)}{r^2}.
    \end{eqnarray*}
    %
    $\Ddsty Q_{n-2}^m$ の項の係数は
    %
    \begin{eqnarray*}
      & & -\frac{(1-r^2)}{r^2}(n+\gamma-3)(n+\gamma-3)
           \left(-2\alpha +(\beta-1)\frac{(1-r^2)}{r^2}\right)(n+\gamma-3)\\
      & & +\frac{|m|(|m|+\beta-1)}{r^2}
          - (n-2)(n+2\alpha+\beta-3) \\
      &=& (n+\gamma-3)^2
          [-2\alpha -(\beta-1)](n+\gamma-3)
          - (n-2)(n+\gamma-3)\\
      & &  + \frac{-(n+\gamma-3)^2+(\beta-1)(n+\gamma-3)
                  +|m|(|m|+\beta-1)}{r^2}\\
      &=& (n+\gamma-3)[(n+\gamma-3)-2\alpha -\beta+1 - (n-2)]\\
      & &  + \frac{-(n+\gamma-3)^2+(\beta-1)(n+\gamma-3)
                  +|m|^2+|m|(\beta-1)}{r^2}\\
      &=& \frac{-(n+\gamma-|m|-3)(n+\gamma+|m|-3)
                  +(\beta-1)(n+\gamma+|m|-3)}{r^2}\\
      &=& -\frac{(n+\gamma+|m|-3)(n+\gamma-|m|-\beta-2)}{r^2}.
    \end{eqnarray*}
    %
    まとめると
    %
    \begin{eqnarray*}
      & &
      \frac{(1-r^2)}{r}(2n+\gamma-3)\DD{Q_n^m}{r}
      +\left[n(2n+\gamma-3)
            -\frac{|m|(|m|+\beta-1)+n(n+\gamma-\beta-2)}{r^2}\right]Q_n^m\\
      & &
      -\frac{(n+\gamma-|m|-3)(n+\gamma+|m|-\beta-2)}{r^2}Q_{n-2}^m =0,
      \\
      & &
        (1-r^2)r\DD{Q_n^m}{r} =
        \left[-nr^2
              +\frac{|m|(|m|+\beta-1)+n(n+\gamma-\beta-2)}{2n+\gamma-3}
        \right]Q_n^m\\
      & &+\frac{(n+\gamma+|m|-3)(n+\gamma-|m|-\beta-2)}{2n+\gamma-3}Q_{n-2}^m,
    \end{eqnarray*}
  }, 
  %
  \begin{eqnarray}
    \Deqlab{微分漸化式2}
    (1-r^2)r\DD{}{r}Q_n^m(r) 
    &=& \left(-nr^2 + \frac{|m|(|m|+\beta-1)+n(n+\gamma-\beta-2)}
                         {2n+\gamma-3}\right)Q_n^m(r)\nonumber\\
    & &+ \frac{(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)}{2n+\gamma-3}Q_{n-2}^m. 
  \end{eqnarray}
  %
  \Deqref{微分漸化式1}を再度用いて微分項を消去すると
  \footnote{
    \Deqref{微分漸化式2}を $n\rightarrow n+2$ とした式
    %
    \begin{eqnarray*}
      (1-r^2)r\DD{}{r}Q_{n+2}^m(r) 
      &=& \left(-(n+2)r^2 + \frac{|m|(|m|+\beta-1)+(n+2)(n+\gamma-\beta)}
                           {2n+\gamma+1}\right)Q_{n+2}^m\\
      & &+ \frac{(n-|m|+\gamma-\beta)(n+|m|+\gamma-1)}{2n+\gamma+1}Q_n^m. 
    \end{eqnarray*}
    %
    と\Deqref{微分漸化式2}を引き算して
    %
    \begin{eqnarray*}
      (1-r^2)r\DD{}{r}(Q_{n+2}^m-Q_n^m)
      &=&  \left(-(n+2)r^2 + \frac{|m|(|m|+\beta-1)+(n+2)(n+\gamma-\beta)}
                                  {2n+\gamma+1}\right)Q_{n+2}^m(r)\\
      & & + \frac{(n-|m|+\gamma-\beta)(n+|m|+\gamma-1)}{2n+\gamma+1}Q_n^m\\
      & & - \left(-nr^2 + \frac{|m|(|m|+\beta-1)+n(n+\gamma-\beta-2)}
                               {2n+\gamma-3}\right)Q_n^m(r)\\
      & & - \frac{(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)}
                 {2n+\gamma-3}Q_{n-2}^m. 
    \end{eqnarray*}
    %
    左辺に\Deqref{微分漸化式1}を適用すると
    %
    \begin{eqnarray*}
      & &
      (1-r^2)[(n+2)Q_{n+2}^m+(n+\gamma-1)Q_n^m]\\
      & &
       =   \left(-(n+2)r^2 + \frac{|m|(|m|+\beta-1)+(n+2)(n+\gamma-\beta)}
                                  {2n+\gamma+1}\right)Q_{n+2}^m
          + \frac{(n-|m|+\gamma-\beta)(n+|m|+\gamma-1)}{2n+\gamma+1}Q_n^m\\
      & & - \left(-nr^2 + \frac{|m|(|m|+\beta-1)+n(n+\gamma-\beta-2)}
                         {2n+\gamma-3}\right)Q_n^m(r)
         - \frac{(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)}
                {2n+\gamma-3}Q_{n-2}^m, \\
      & &
             \left(-(1-r^2)(n+2)-(n+2)r^2 
                 + \frac{|m|(|m|+\beta-1)+(n+2)(n+\gamma-\beta)}
                        {2n+\gamma+1}\right)Q_{n+2}^m\\
      & &  + \left[-(n+\gamma-1)(1-r^2) 
           + \frac{(n-|m|+\gamma-\beta)(n+|m|+\gamma-1)}
                  {2n+\gamma+1}
                     +nr^2 \right.\\
      & &
               \left.  - \frac{|m|(|m|+\beta-1)+n(n+\gamma-\beta-2)}
                              {2n+\gamma-3}\right]Q_n^m
           - \frac{(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)}
                  {2n+\gamma-3}Q_{n-2}^m. 
    \end{eqnarray*}
    %
    $Q_{n+2}^m$ の係数は
    %
    \begin{eqnarray*}
      & &
         -(1-r^2)(n+2)-(n+2)r^2 
        + \frac{|m|(|m|+\beta-1)+(n+2)(n+\gamma-\beta)}{2n+\gamma+1}\\
      &=& -(n+2)+ \frac{|m|(|m|+\beta-1)+(n+2)(n+\gamma-\beta)}{2n+\gamma+1}\\
      &=& \frac{-(n+2)(2n+\gamma+1)+ |m|(|m|+\beta-1)+(n+2)(n+\gamma-\beta)}
               {2n+\gamma+1}\\
      &=& \frac{-(n+2)(n+\beta+1)+ |m|(|m|+\beta-1)}{2n+\gamma+1}
       =  \frac{-(n+2)^2-(n+2)(\beta-1)+ |m|^2+ |m|(\beta-1)}{2n+\gamma+1}\\
      &=& \frac{-(n+|m|+2)(n-|m|+2)-(n-|m|+2)(\beta-1)}{2n+\gamma+1}
       = \frac{-(n+|m|+\beta+1)(n-|m|+2)}{2n+\gamma+1}
    \end{eqnarray*}
    %
    $Q_n^m$ の係数は
    %
    \begin{eqnarray*}
      & &  -(n+\gamma-1)(1-r^2) + \frac{(n-|m|+\gamma-\beta)(n+|m|+\gamma-1)}
                                       {2n+\gamma+1}
                     +nr^2 - \frac{|m|(|m|+\beta-1)+n(n+\gamma-\beta-2)}
                                  {2n+\gamma-3}\\
      &=&  (2n+\gamma-1)r^2
           -(n+\gamma-1) + \frac{(n-|m|+\gamma-\beta)(n+|m|+\gamma-1)}
                                {2n+\gamma+1}
                         - \frac{|m|(|m|+\beta-1)+n(n+\gamma-\beta-2)}
                                  {2n+\gamma-3}\\
    \end{eqnarray*}
    %
    $r^2$ の項以外の項を通分した分子は
    %
    \begin{eqnarray*}
      & &
           -(n+\gamma-1)(2n+\gamma+1)(2n+\gamma-3)
           +(n-|m|+\gamma-\beta)(n+|m|+\gamma-1)(2n+\gamma-3)\\
      & &  -|m|(|m|+\beta-1)(2n+\gamma+1)
           -n(n+\gamma-\beta-2)(2n+\gamma+1)\\
      &=&
           -(n+\gamma-1)(2n+\gamma+1)(2n+\gamma-3)
          +\{(n+\gamma-\beta)(n+\gamma-1)-(\beta-1)|m|-|m|^2\}(2n+\gamma-3) \\
      & &  -|m|(|m|+\beta-1)(2n+\gamma+1)
           -n(n+\gamma-\beta-2)(2n+\gamma+1)\\
      &=&
           -(n+\gamma-1)(2n+\gamma+1)(2n+\gamma-3)
          +\{(n+\gamma-\beta)(n+\gamma-1)-|m|(|m|+\beta-1)\}(2n+\gamma-3)\\
      & &  -|m|(|m|+\beta-1)(2n+\gamma+1)
           -n(n+\gamma-\beta-2)(2n+\gamma+1)\\
      &=&
           -(n+\gamma-1)(2n+\gamma+1)(2n+\gamma-3)
          + (n+\gamma-\beta)(n+\gamma-1)(2n+\gamma-3)\\
      & &  -2|m|(|m|+\beta-1)(2n+\gamma-1)
           -n(n+\gamma-\beta-2)(2n+\gamma+1)\\
      &=&
           -(n+\gamma-1)(2n+\gamma+1)(2n+\gamma-3)
          + (n+\gamma)(n+\gamma-1)(2n+\gamma-3)-\beta(n+\gamma-1)(2n+\gamma-3)\\
      & &  -2|m|(|m|+\beta-1)(2n+\gamma-1)
           -n(n+\gamma-2)(2n+\gamma+1)+\beta n(2n+\gamma+1)\\
      &=&
           -(n+\gamma-1)(2n+\gamma+1)(2n+\gamma-3)
          + (n+\gamma)(n+\gamma-1)(2n+\gamma-3)\\
      & &  -2|m|(|m|+\beta-1)(2n+\gamma-1)
           -n(n+\gamma-2)(2n+\gamma+1)
           +\beta[n(2n+\gamma+1)-(n+\gamma-1)(2n+\gamma-3)]\\
      &=&
           -(n+\gamma-1)(2n+\gamma+1)(2n+\gamma-3)
          + (n+\gamma)(n+\gamma-1)(2n+\gamma-3)\\
      & &  -2|m|(|m|+\beta-1)(2n+\gamma-1)
           -n(n+\gamma-2)(2n+\gamma+1)
           -\beta(\gamma-3)(2n+\gamma-1)\\
      &=&  
           -(n+\gamma-1)(2n+\gamma-3)(n+1)
           -n(n+\gamma-2)(2n+\gamma+1)\\
      & &  -2|m|(|m|+\beta-1)(2n+\gamma-1)
           -\beta(\gamma-3)(2n+\gamma-1)\\
      &=&  
           -(2n+1)\gamma^2 -(6n^2-2n-4)\gamma-(4n^3-6n^2-4n-3)\\
      & &  -2|m|(|m|+\beta-1)(2n+\gamma-1)
           -\beta(\gamma-3)(2n+\gamma-1)\\
      &=&  
           -(2n+\gamma-1)[(2n+1)\gamma + (2n^3-2n-3)]\\
      & &  -2|m|(|m|+\beta-1)(2n+\gamma-1)
           -\beta(\gamma-3)(2n+\gamma-1)\\
      &=&  
           -(2n+\gamma-1)[2n(n+\gamma+1)+(\gamma-3)]
           -2|m|(|m|+\beta-1)(2n+\gamma-1)
           -\beta(\gamma-3)(2n+\gamma-1)\\
      &=&  
           -(2n+\gamma-1)2n(n+\gamma+1)
           -2|m|(|m|+\beta-1)(2n+\gamma-1)
           -(\beta+1)(\gamma-3)(2n+\gamma-1)
    \end{eqnarray*}
    %
    よって, 全ての項をまとめると
    %
    \begin{eqnarray*}
      & &
           \frac{-(n+|m|+\beta+1)(n-|m|+2)}{2n+\gamma+1}Q_{n+2}^m\\
      & & +(2n+\gamma-1)\left[r^2       
           -\frac{2n(n+\gamma+1)+2|m|(|m|+\beta-1)+(\beta+1)(\gamma-3)}
                 {(2n+\gamma+1)(2n+\gamma-3)}\right]Q_n^m \\
      & &  - \frac{(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)}
                  {2n+\gamma-3}Q_{n-2}^m, 
    \end{eqnarray*}
    %
    分母をはらって
    \begin{eqnarray*}
      & &
           -(n+|m|+\beta+1)(n-|m|+2)(2n+\gamma-3)Q_{n+2}^m\\
      & & +(2n+\gamma-1)\left[(2n+\gamma+1)(2n+\gamma-3)r^2       
                              -2n(n+\gamma+1)-2|m|(|m|+\beta-1)
                              -(\beta+1)(\gamma-3)\right]Q_n^m \\
      & &  -(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)(2n+\gamma+1)Q_{n-2}^m.
    \end{eqnarray*}
  }, 
  %
  \begin{eqnarray}
    \Deqlab{漸化式}
    & &
    -(n-|m|+2)(n+|m|+\beta+1)(2n+\gamma-3)Q_{n+2}^m(r) \nonumber\\
    & &
    +(2n+\gamma-1)\{ (2n+\gamma-3)(2n+\gamma+1)r^2 -2n(n+\gamma-1)\nonumber\\
    & &
                    -2|m|(|m|+\beta-1)
                    -(\gamma-3)(\beta+1)\}Q_n^m(r)\nonumber\\
    & &                    
    -(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)(2n+\gamma+1)Q_{n-2}^m(r)=0
  \end{eqnarray}

  式 (3) において $n=|m|,|m|+2$ と置くと
  \footnote{
    \begin{eqnarray*}
       Q_{|m|}^{|m|} &=& \sum_{p=0}^0 
            \frac{\Ddsty (-1)^{p}
                         \Gamma\left(\frac{2|m|+\gamma-1}{2}+p\right)
                         \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
                 {\Ddsty p! \left(-p\right)!
                        \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                        \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}
                r^{|m|+2p}\\
          &=&
            \frac{\Ddsty \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)
                         \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
                 {\Ddsty 0! \left(0\right)!
                        \Gamma\left(\frac{2|m|+\beta+1}{2}\right)
                        \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}
                r^{|m|}=r^{|m|},\\
       Q_{|m|+2}^{|m|} &=& \sum_{p=0}^1 
            \frac{\Ddsty (-1)^{p+1}
                         \Gamma\left(\frac{2|m|+\gamma+1}{2}+p\right)
                         \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
                 {\Ddsty p! \left(1-p\right)!
                         \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                         \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}
                 r^{|m|+2p}\\
           &=&
            \frac{\Ddsty (-1)
                         \Gamma\left(\frac{2|m|+\gamma+1}{2}\right)
                         \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
                 {\Ddsty 0! \left(1\right)!
                         \Gamma\left(\frac{2|m|+\beta+1}{2}+p\right)
                         \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}
                 r^{|m|}\\
           & &+\frac{\Ddsty (-1)^{2}
                         \Gamma\left(\frac{2|m|+\gamma+1}{2}+1\right)
                         \Gamma\left(\frac{2|m|+\beta+1}{2}\right)}
                 {\Ddsty 1! \left(0\right)!
                         \Gamma\left(\frac{2|m|+\beta+1}{2}+1\right)
                         \Gamma\left(\frac{2|m|+\gamma-1}{2}\right)}
                 r^{|m|+2}\\
           &=&
              -\frac{2|m|+\gamma-1}{2}r^{|m|}
              +\frac{\Ddsty \frac{2|m|+\gamma+1}{2}
                            \frac{2|m|+\gamma-1}{2}}
                    {\Ddsty \frac{2|m|+\beta+1}{2}}
                 r^{|m|+2}\\
           &=& \frac{2|m|+\gamma-1}{2}
               \left(\frac{2|m|+\gamma+1}{2|m|+\beta+1}r^2-1 \right)r^{|m|}
    \end{eqnarray*}
  }
  %
  \begin{eqnarray}
    \Deqlab{mm次関数}
     & & Q_{|m|}^{|m|} = r^{|m|}, \\
    \Deqlab{mm-2次関数}
     & & Q_{|m|+2}^{|m|} 
        = \frac{2|m|+\gamma-1}{2}
             \left(\frac{2|m|+\gamma+1}{2|m|+\beta+1}r^2-1\right)r^{|m|}
  \end{eqnarray}
  %
  これらの値からスタートして高次の $Q_n^m(r)$ の値を
  漸化式\Deqref{漸化式}によって求めることができる. 
  正規化係数 $I_n^m$ の漸化式は\Deqref{漸化式}に直交関係\Deqref{直交関係}を
  適用することにより得られる
  \footnote{
    \Deqref{漸化式}に $Q_{n+2}w(r)$ をかけて積分すると, 
    %
    \begin{eqnarray*}
      & &
      -(n-|m|+2)(n+|m|+\beta+1)(2n+\gamma-3)I_{n+2}\\
      & &+(2n+\gamma-1)(2n+\gamma-3)(2n+\gamma+1) 
      \int_0^1 r^2 Q_{n+2}^mQ_n^m w(r)dr = 0, \\
      & &
      -(n-|m|+2)(n+|m|+\beta+1)I_{n+2}
      +(2n+\gamma-1)(2n+\gamma+1)\int_0^1 r^2 Q_{n+2}^mQ_n^m w(r)dr = 0, \\
    \end{eqnarray*}
    %
    $n\rightarrow n-2$ と次数を下げると
    %
    \begin{displaymath}
      -(n-|m|)(n+|m|+\beta-1)I_n
      +(2n+\gamma-5)(2n+\gamma-3)\int_0^1 r^2 Q_n^mQ_{n-2}^m w(r)dr
    \end{displaymath}
    %
    再度\Deqref{漸化式}に $Q_{n-2}w(r)$ をかけて積分すると, 
    \begin{eqnarray*}
      & &
      (2n+\gamma-1)(2n+\gamma-3)(2n+\gamma+1)
        \int_0^1 r^2 Q_n^mQ_{n-2}^m w(r)dr\\
      & &-(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)(2n+\gamma+1)I_{n-2}=0\\
      & &
      (2n+\gamma-1)(2n+\gamma-3)\int_0^1 r^2 Q_n^mQ_{n-2}^m w(r)dr
        -(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)I_{n-2}=0
    \end{eqnarray*}
    %
    先の式とあわせて積分の項を消去できて
    %
    \begin{eqnarray*}
      & &
      -(n-|m|)(n+|m|+\beta-1)(2n+\gamma-1) I_n
      +(2n+\gamma-5)(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)I_{n-2}=0, \\
      & &
      I_n = \frac{(2n+\gamma-5)(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)}
                   {(n-|m|)(n+|m|+\beta-1)(2n+\gamma-1)}I_{n-2}
    \end{eqnarray*}
  }. 
  \begin{equation}
    \Deqlab{正規化係数漸化式}
    I_n^m = \frac{(2n+\gamma-5)(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)}
                 {(n-|m|)(n+|m|+\beta-1)(2n+\gamma-1)}I_{n-2}^m,
  \end{equation}
  %
  スタートの値 $I_{|m|}^m$ は\Deqref{mm次関数}を\Deqref{直交関係}に
  代入して
  \footnote{
    \begin{eqnarray*}
      I_{|m|}^m &=& \int_0^1 r^{2|m|+\beta}(1-r^2)^{\alpha-1}dr
                 = \frac{1}{2}\int_0^1 t^{|m|+(\beta-1)/2}(1-t)^{\alpha-1}dt\\
                &=& \frac{1}{2}B\left(\alpha,|m|+\frac{\beta+1}{2}\right)
                = \frac{\Ddsty \Gamma(\alpha)
                              \Gamma\left(|m|+\frac{\beta+1}{2}\right)}
                       {\Ddsty 2\Gamma\left(\alpha+|m|+\frac{\beta+1}{2}\right)}
                = \frac{\Ddsty \Gamma\left(\frac{\gamma-\beta}{2}\right)
                              \Gamma\left(|m|+\frac{\beta+1}{2}\right)}
                       {\Ddsty 2\Gamma\left(|m|+\frac{\gamma+1}{2}\right)}
    \end{eqnarray*}
    %
    ただしガンマ関数, ベータ関数の公式
    \begin{displaymath}
      B(p,q) = \int_0^1(1-t)^{p-1}t^{q-1}dt, \quad
      B(p,q) = \frac{\Gamma(p)\Gamma(q)}{\Gamma(p+q)},
    \end{displaymath}
    %
    を用いた
  }
  \begin{equation}
    I_{|m|}^m= \frac{\Ddsty \Gamma\left(\frac{\gamma-\beta}{2}\right)
                            \Gamma\left(|m|+\frac{\beta+1}{2}\right)}
                    {\Ddsty 2\Gamma\left(|m|+\frac{\gamma+1}{2}\right)}.
  \end{equation}

  次に\Deqref{極座標スペクトル展開}でのスペクトル係数 $f_n^m$ を
  評価する. 偶数 $M$ に対して $f(r,\phi)$ の部分和近似 $f_M(r,\phi)$, 
  %
  \begin{equation}
    \Deqlab{スペクトル逆変換}
    f_m(r) = \sum_{n=|m|,n+m \mbox{ even}}^{\hat{M}}f_n^m\Phi_n^m(r), 
    \qquad
    f_M(r,\phi) = \sum_{m=-M+1}^{M-1} f_m(r) e^{im\phi} \sim f(r,\phi),
  \end{equation}
  %
  を考える. ここで $\hat{M}$ は $m$ が奇数のとき $\hat{M}=M-1$, 
  偶数のとき $\hat{M}=M-2$ である. 
  $M=2^p$ を選べば $f_m(r)$ を求めるのに FFT が使える. 
  \Deqref{スペクトル逆変換}の逆変換は
  %
  \begin{equation}
    \Deqlab{スペクトル正変換}
    f_n^m = \int_0^1 f_m(r)\Phi_n^m(r)w(r) dr.
  \end{equation}
  %
  ここでわれわれは $f_n^m$ を効率的に計算するための
  \Deqref{スペクトル正変換}に対するガウス積分を導出する. 
  まず関数の積 $f_m(r)\Phi_n^m(r)$ がせいぜい $2M-2$ 次の
  偶関数であることに注意しよう. 
  \Deqref{スペクトル正変換}のラグランジュ積分公式は(補遺A 参照), 
  %
  \begin{equation}
    \Deqlab{スペクトル正変換ガウス積分}
    \int_0^1 f_m(r)\Phi_n^m(r)w(r) dr 
    = \sum_{i=1}^{M/2} f_m(r_i)\Phi_n^m(r_i)w_i + E,
  \end{equation}
  %
  ここで $r_i$ は積分座標の格子点, $E$ は公式の誤差である. 
  ここでは $ 0\le r_1 < r_2 < \ldots < r_{M/2} \le 1$ を仮定する. 
  重み $w_i$ は
  %
  \begin{equation}
    \Deqlab{ガウス積分重み係数定義}
    w_i \equiv \frac{1}{\Ddsty\left.\frac{\pi'(r)}{2r}\right|_{r_i}}
                 \int_0^1 \frac{\pi(r)}{r^2-r_i^2} w(r) dr,
  \end{equation}
  %
  ただし, 
  %
  \begin{equation}
    \pi(r) = \prod_{j=1}^{M/2}(r^2 - r_j^2).
  \end{equation}
  %
  公式\Deqref{スペクトル正変換ガウス積分}は関数の積 $f_m(r)\Phi_n^m(r)$ が
  $M-2$ 次と同じかそれ以下であるときに正確になる. 
  さて, $f_m(r)\Phi_n^m(r)$ を $\pi(r)$ で割り, 
  $f_m(r)\Phi_n^m(r) = s(r) + \pi(r)q(r)$ と書き直すことができる. 
  ただし $s(r)$ と $q(r)$ はせいぜい $M-2$ 次の偶関数である. 
  これを\Deqref{スペクトル正変換ガウス積分}に代入し, 
  $1\le i \le M/2$ に対して $s(r_i) = f_m(r_i)\Phi_n^m(r_i)$ であることを
  用いると
  \footnote{
    実際に代入すると
    \begin{displaymath}
      \int_0^1 f_m(r)\Phi_n^m(r)w(r) dr 
      = \int_0^1 s(r)w(r) dr + \int_0^1 \pi(r)q(r)w(r) dr
      = \sum_{i=1}^{M/2} s(r_i)w_i + E,
    \end{displaymath}
    %
    第 1 項どうしが等しいので $E$ の表式が得られる. 
  }
  %
  \begin{equation}
    E = \int_0^1 \pi(r)q(r)w(r)dr.
  \end{equation}
  %
  したがって, もしわれわれが $\pi(r)=Q_M^0(r)$ と選べば, 
  直交関係から $E$ は 0 となり
  \footnote{
    $Q_M^0$ は $M$ 次の多項式であり $M$ 次より小さい多項式とは直交する. 
    $q(r)$ はせいぜい $M-2$ 次なので誤差の積分は 0 となる. 
  }, 
  $f_M(r,\phi)$ が $f(r,\phi)$ と等しいときには
  係数 $f_n^m$ を
  \Deqref{スペクトル正変換}と\Deqref{スペクトル正変換ガウス積分}から
  正確に計算できる. 
  よってガウス積分の格子点は $Q_M^0$ の正の零点に対応する. 
  $Q_M^0$ の零点は数値的に計算しなければならない. 

  ガウス積分の重み係数 $w_i$ は以下のようにして見積もられる. 
  まず\Deqref{漸化式}を書き直して
  \footnote{
    \Deqref{漸化式}で $m=0$ とした式は
    \begin{eqnarray*}
      & &
      -(n+2)(n+\beta+1)(2n+\gamma-3)Q_{n+2}^0(r) \\
      & &
      +(2n+\gamma-1)\{ (2n+\gamma-3)(2n+\gamma+1)r^2 -2n(n+\gamma-1)
                      -(\gamma-3)(\beta+1)\}Q_n^m(r)\\
      & &                    
      -(n+\gamma-\beta-2)(n+\gamma-3)(2n+\gamma+1)Q_{n-2}^m(r)=0, \\
   & &
      (2n+\gamma-1)(2n+\gamma-3)(2n+\gamma+1)r^2Q_n^m(r) 
      = (n+2)(n+\beta+1)(2n+\gamma-3)Q_{n+2}^0(r)\\
      & &
      +(2n+\gamma-1)[2n(n+\gamma-1)+(\gamma-3)(\beta+1)]Q_n^m(r) \\
      & &
      +(n+\gamma-\beta-2)(n+\gamma-3)(2n+\gamma+1)Q_{n-2}^m(r)=0, \\
   & &
      r^2Q_n^m(r) = \frac{(n+2)(n+\beta+1)}
                         {(2n+\gamma-1)(2n+\gamma+1)}Q_{n+2}^0(r) \\
      & &
      +\frac{2n(n+\gamma-1)+(\gamma-3)(\beta+1)}
            {(2n+\gamma-3)(2n+\gamma+1)}Q_n^m(r)
      +\frac{(n+\gamma-\beta-2)(n+\gamma-3)}
            {(2n+\gamma-1)(2n+\gamma-3)}Q_{n-2}^m(r)=0, \\
    \end{eqnarray*}
    %
    ここで
    $\Ddsty C_n \equiv \frac{(n+2)(n+\beta+1)}{(2n+\gamma-1)(2n+\gamma+1)}$ を
    用いると, 
    $\Ddsty C_{n-2} \equiv \frac{n(n+\beta-1)}{(2n+\gamma-5)(2n+\gamma-3)}$ で
    あるから,
    %
    \begin{eqnarray*}
    & &
      r^2Q_n^m(r) = C_nQ_{n+2}^0(r) \\
      & &
      +\frac{2n(n+\gamma-1)+(\gamma-3)(\beta+1)}
            {(2n+\gamma-3)(2n+\gamma+1)}Q_n^m(r)
      + C_{n-2}\frac{(2n+\gamma-5)(n+\gamma-\beta-2)(n+\gamma-3)}
                    {n(n+\beta-1)(2n+\gamma-1)}Q_{n-2}^m(r)=0, \\
    \end{eqnarray*}
    %
    $I_n^0$ で割ると
    %
    \begin{displaymath}
      \frac{r^2Q_n^m(r)}{I_n^0} = \frac{C_n}{I_n^0}Q_{n+2}^0(r)
      + D_n Q_n^m(r)
      + \frac{C_{n-2}}{I_n^0}
        \frac{(2n+\gamma-5)(n+\gamma-\beta-2)(n+\gamma-3)}
             {n(n+\beta-1)(2n+\gamma-1)}Q_{n-2}^m(r)=0.
    \end{displaymath}
    %
    ただし $\Ddsty D_n = \frac{2n(n+\gamma-1)+(\gamma-3)(\beta+1)}
                              {I_n^0(2n+\gamma-3)(2n+\gamma+1)}$ である. 
    \Deqref{正規化係数漸化式}で $m=0$ とした式から
    %
    \begin{displaymath}
      I_n^0 = \frac{(2n+\gamma-5)(n+\gamma-\beta-2)(n+\gamma-3)}
                    {n(n+\beta-1)(2n+\gamma-1)}I_{n-2}^0,
    \end{displaymath}
    %
    したがって, 
    %
    \begin{displaymath}
      \frac{r^2Q_n^m(r)}{I_n^0} = \frac{C_n}{I_n^0}Q_{n+2}^0(r)
      + D_n Q_n^m(r)
      + \frac{C_{n-2}}{I_{n-2}^0}Q_{n-2}^m(r)=0.
    \end{displaymath}
  },
  \begin{equation}
    \Deqlab{m0漸化式}
    \frac{r^2Q_n^0(r)}{I_n^0}= \frac{C_n}{I_n^0}Q_{n+2}^0(r)
                             +  D_n Q_n^0(r)
                             + \frac{C_{n-2}}{I_{n-2}^0}Q_{n-2}^0(r),
  \end{equation}
  %
  ただし
  %
  \begin{equation}
    C_n \equiv \frac{(n+2)(n+\beta+1)}{(2n+\gamma-1)(2n+\gamma+1)}. 
  \end{equation}
  %
  $D_n$ の詳細は必要ないので定義を述べない. 
  \Deqref{m0漸化式}に $Q_n^0(\rho)$ をかけて $r$ と $\rho$ を入れ換えた式を
  引き算し, $n=0$ から $M$ までの偶数を足しあわせると
  クリストッフェル--ダブルー公式が得られる
  \footnote{
    \Deqref{m0漸化式}に $Q_n^0(\rho)$ をかけた式は
    %
    \begin{displaymath}
      \frac{r^2Q_n^0(r)Q_n^0(\rho)}{I_n^0}
             = \frac{C_n}{I_n^0}Q_{n+2}^0(r)Q_n^0(\rho)
             +  D_n Q_n^0(r)Q_n^0(\rho)
             + \frac{C_{n-2}}{I_{n-2}^0}Q_{n-2}^0(r)Q_n^0(\rho),
    \end{displaymath}
    %
    $r$ と $\rho$ を入れ換えた式は
    %
    \begin{displaymath}
      \frac{\rho^2Q_n^0(r)Q_n^0(\rho)}{I_n^0}
             = \frac{C_n}{I_n^0}Q_{n+2}^0(\rho)Q_n^0(r)
             +  D_n Q_n^0(r)Q_n^0(\rho)
             + \frac{C_{n-2}}{I_{n-2}^0}Q_{n-2}^0(\rho)Q_n^0(r),
    \end{displaymath}
    %
    辺々引き算して
    %
    \begin{eqnarray*}
      & &
      \frac{(r^2-\rho^2)Q_n^0(r)Q_n^0(\rho)}{I_n^0}
             = \frac{C_n}{I_n^0}[Q_{n+2}^0(r)Q_n^0(\rho)
                                 Q_{n+2}^0(\rho)Q_n^0(r)]
             + \frac{C_{n-2}}{I_{n-2}^0}
               [Q_{n-2}^0(r)Q_n^0(\rho)-Q_{n-2}^0(\rho)Q_n^0(r)]\\
      & &
             = \frac{C_n}{I_n^0}[Q_{n+2}^0(r)Q_n^0(\rho)
                                -Q_{n+2}^0(\rho)Q_n^0(r)]
             - \frac{C_{n-2}}{I_{n-2}^0}
               [Q_n^0(r)Q_{n-2}^0(\rho)-Q_n^0(\rho)Q_{n-2}^0(r)]
    \end{eqnarray*}
    %
    $n=0$ から $M$ までの偶数に関して足しあわせると $Q_{-2}^0=0$ に注意して
    %
    \begin{displaymath}
      (r^2-\rho^2)\sum_{p=0,\mbox{ even}}^M
             \frac{Q_p^0(r)Q_p^0(\rho)}{I_p^0}
             = \frac{C_M}{I_M^0}[Q_{M+2}^0(r)Q_M^0(\rho)
                                -Q_{M+2}^0(\rho)Q_M^0(r)].
    \end{displaymath}
  }. 
  %
  \begin{equation}
    \Deqlab{クリストッフェル--ダブルー公式}
    \sum_{p=0,\mbox{even}}^{M}\frac{Q_p^0(r)Q_p^0(\rho)}{I_p^0}
      = \frac{C_M}{I_M^0}\left(
        \frac{Q_{M+2}^0(r)Q_M^0(\rho)-Q_{M+2}^0(\rho)Q_M^0(r)}
             {r^2-\rho^2} \right).
  \end{equation}
  %
  $\rho$ を $r_i$ に置き換えて $Q_M^0(r_i)=0$ であることを用い, 
  両辺に $Q_0^0(r)w(r)$ をかけて 0 から 1 まで積分し, 
  さらに\Deqref{m0漸化式}で $r=r_i$ した式を用いると
  \footnote{
    $\rho$ を $r_i$ に置き換えて $Q_M^0(r_i)=0$ であることを用いると
    \begin{displaymath}
       \sum_{p=0,\mbox{even}}^{M}\frac{Q_p^0(r_i)Q_p^0(r)}{I_p^0}
      = -\frac{C_M}{I_M^0}\frac{Q_{M+2}^0(r_i)Q_M^0(r)}{r^2-r_i^2}.
    \end{displaymath}
    %
    $Q_0^0(r)w(r)$ をかけて 0 から 1 まで積分すると, $Q_p^0Q_0^0$ の
    積分は直交関係から $p=0$ の項だけ残る. さらに $Q_0^0(r)=1$ であるから
    %
    \begin{displaymath}
       \frac{Q_0^0(r_i)}{I_0^0}\int_0^1Q_0^0(r)^2 w(r) dr
      = -\frac{C_MQ_{M+2}^0(r_i)}{I_M^0}
    \end{displaymath}
    %
    $\Ddsty \int_0^1Q_0^0(r)^2 w(r) dr=I_0^0$ なので左辺は 1 となる. 
    したがって
    %
    \begin{displaymath}
      \int_0^1 \frac{Q_M^0(r)}{r^2-r_i^2}w(r)dr
        =-\frac{I_M^0}{C_M Q_{M+2}^0(r_i)}.
    \end{displaymath}
    %
    ここで\Deqref{m0漸化式}において $n=M,r=r_i$ を適用すると 
    $Q_M^0(r_i)=0$ なので
    %
    \begin{displaymath}
       0 = \frac{C_M}{I_M^0} Q_{M+2}^0(r_i)
          + \frac{C_{M-2}}{I_{M-2}^0}Q_{M-2}^0(r_i),
    \end{displaymath}
    %
    となる. したがって
    %
    \begin{displaymath}
      \int_0^1\frac{Q_M^0(r)}{r^2-r_i^2}w(r)dr
      =\frac{I_{M-2}^0}{C_{M-2} Q_{M-2}^0(r_i)}.
    \end{displaymath}
  }
  \begin{equation}
    \Deqlab{-2乗積分}
    \int_0^1 \frac{Q_M^0}{r^2-r_i^2}w(r)dr 
     = \frac{I_{M-2}^0}{C_{M-2}^0 Q_{M-2}^0(r_i)}.
  \end{equation}
  %
  この式と\Deqref{ガウス積分重み係数定義}を $\pi(r)=Q_M^0(r)$ とともに
  用い, \Deqref{微分漸化式2}から $w_i$ の式が得られる
  \footnote{
    \Deqref{微分漸化式2}において $r=r_i,n=M, m=0$ とした式は, 
    $Q_M^0(r_i)=0$ であるから
    %
    \begin{displaymath}
      (1-r_i^2)r_i\DD{}{r}Q_M^0(r_i) 
        =\frac{(M+\gamma-\beta-2)(M+\gamma-3)}{2M+\gamma-3}Q_{M-2}^0. 
    \end{displaymath}
    %
    したがって
    %
    \begin{displaymath}
      \left.\frac{\pi'(r)}{2r}\right|_{r_i} 
      = \frac{(M+\gamma-\beta-2)(M+\gamma-3)}
             {2(2M+\gamma-3)r_i^2(1-r_i^2)}Q_{M-2}^0,
    \end{displaymath}
    %
    これと\Deqref{-2乗積分}を\Deqref{ガウス積分重み係数定義}に代入すると
    \begin{displaymath}
      w_i = \frac{2(2M+\gamma-3)r_i^2(1-r_i^2)}
                 {(M+\gamma-\beta-2)(M+\gamma-3)Q_{M-2}^0}
            \frac{I_{M-2}^0}{C_{M-2} Q_{M-2}^0(r_i)}
          =  \frac{2(2M+\gamma-3)r_i^2(1-r_i^2)}
                 {(M+\gamma-\beta-2)(M+\gamma-3)C_M[\Phi_{M-2}^0]^2}
    \end{displaymath}
  }. 
  %
  \begin{equation}
    \Deqlab{ガウス積分重み係数}
    w_i = \frac{2(2M+\gamma-3)r_i^2(1-r_i^2)}
               {(M+\gamma-\beta-2)(M+\gamma-3)C_{M-2}[\Phi_M^0(r_i)]^2}.
  \end{equation}

  選点法では積分格子点は境界条件を適用するために, 境界上に格子点が
  ある方が良い. その場合ガウス--ラダウ(Gauss-Radau)積分法を用いる
  ことができる. 
  $\pi(r)$ の適切な選び方は
  %
  \begin{equation}
    \Deqlab{ガウスラダウπ}
    \pi(r)=Q_M^0(r) - \frac{(M+\gamma-\beta-2)(M+\gamma-3)}
                           {M(M+\beta-1)}Q_{M-2}^0(r).
  \end{equation}
  %
  ただし $\pi(r)$ は $\pi(1)=0$ となるようにとってある
  \footnote{
    \Deqref{微分漸化式2}において $n=M, m=0, r=1$ と選ぶと
    \begin{eqnarray*}
      & &
       0 = \left(-M + \frac{M(M+\gamma-\beta-2)}
                             {2M+\gamma-3}\right)Q_M^m(1)
        + \frac{(M+\gamma-\beta-2)(M+\gamma-3)}{2M+\gamma-3}Q_{n-2}^m,\\
      & &
       0 = M[-(2M+\gamma-3)+(M+\gamma-\beta-2)]Q_M^m(1)
         + (M+\gamma-\beta-2)(M+\gamma-3)Q_{n-2}^m,\\
      & &
       0 = M[-(M+\beta-1)]Q_M^m(1)
         + (M+\gamma-\beta-2)(M+\gamma-3)Q_{n-2}^m,\\
      & &
       Q_M^m(1)= 
          \frac{(M+\gamma-\beta-2)(M+\gamma-3)}{M(M+\beta-1)}Q_{n-2}^m.
    \end{eqnarray*}
    %
    よって\Deqref{ガウスラダウπ}は $r=1$ で 0 となる. 
  }.
  $1\le i \le M/2-1$ の重み係数は次のように与えられる. 
  %
  \begin{equation}
    w_i = \frac{2(2M+\gamma-5)r_i^2}
               {(M+\gamma-\beta-2)(M+\gamma-3)[\Phi_{M-2}^0(r_i)]^2}
  \end{equation}
  %
  \Deqref{ガウス積分重み係数}を求めたのと同じように
  \Deqref{クリストッフェル--ダブルー公式}と $\pi(r_i)=0$ を用いて
  求めることができる
  \footnote{
    \Deqref{クリストッフェル--ダブルー公式}を $M-2$ まで適用して
    \begin{displaymath}
          \sum_{p=0,\mbox{even}}^{M-2}\frac{Q_p^0(r)Q_p^0(\rho)}{I_p^0}
      = \frac{C_{M-2}}{I_{M-2}^0}\left(
        \frac{Q_{M}^0(r)Q_{M-2}^0(\rho)-Q_M^0(\rho)Q_{M-2}^0(r)}
             {r^2-\rho^2} \right).
    \end{displaymath}
    %
    $\pi(r) = Q_M^0(r) - E_M Q_{M-2}^0(r)$ と書き表すことにすると
    $\pi(r_i)=0$ より $Q_M^0(r_i) = E_M Q_{M-2}^0(r_i)$. 
    $\rho=r_i$ ととることにより
    %
    \begin{eqnarray*}
          \sum_{p=0,\mbox{even}}^{M-2}\frac{Q_p^0(r)Q_p^0(r_i)}{I_p^0}
      &=& \frac{C_{M-2}}{I_{M-2}^0}\left(
        \frac{Q_{M}^0(r)Q_{M-2}^0(r_i)-Q_M^0(r_i)Q_{M-2}^0(r)}
             {r^2-r_i^2} \right)\\
      &=& \frac{C_{M-2}}{I_{M-2}^0}\left(
          \frac{Q_{M}^0(r)Q_{M-2}^0(r_i)-E_MQ_{M-2}^0(r_i)Q_{M-2}^0(r)}
             {r^2-r_i^2} \right)\\
      &=&\frac{C_{M-2}}{I_{M-2}^0}{I_M^0}\left(
         \frac{Q_{M}^0(r)-E_MQ_{M-2}^0(r)}
             {r^2-r_i^2} \right)
      = \frac{C_{M-2} Q_{M-2}^0(r_i)}{I_{M-2}^0}\frac{\pi(r)}{r^2-r_i^2}.
    \end{eqnarray*}
    %
    $\Ddsty \int_0^1 dr Q_0^0(r)w(r)$ を作用させると
    \Deqref{ガウス積分重み係数}を求めたのと同じように左辺が 1 となる. 
    したがって, 
    %
    \begin{eqnarray*}
      1 = \frac{C_{M-2} Q_{M-2}^0(r_i)}
               {I_{M-2}^0}\int_0^1\frac{\pi(r)}{r^2-r_i^2}dr,
      \qquad
      \int_0^1\frac{\pi(r)}{r^2-r_i^2}dr
      =\frac{I_{M-2}^0}{C_{M-2} Q_{M-2}^0(r_i)}.
    \end{eqnarray*}
    %
    一方 $\Ddsty \left.\frac{\pi'(r)}{2r}\right|_{r_i}$ の方は
    %
    \begin{eqnarray*}
      \pi'(r) = Q_M'(r_i)-E_M Q_{M-2}'(r_i)
    \end{eqnarray*}
    %
    \Deqref{微分漸化式1}を用いると
    $ r_i{Q'}_M^0(r_i)- r_i{Q'}_{M-2}(r)
    = M Q_M^0(r_i)+ (M+\gamma-3) Q_{M-2}^0(r_i)$
    であるから, 
    %
    \begin{eqnarray*}
      & &
      \pi'(r) = {Q'}_M^0(r_i)
      - E_M \left[   {Q'}_M^0(r_i)- \frac{M}{r_i} Q_M^0(r_i)
                  - \frac{M+\gamma-3)}{r_i} Q_{M-2}^0(r_i) \right] \\
      & &
      = (1-E_M){Q'}_M^0(r_i)
          + \frac{E_M}{r_i}\left[ ME_M - (M+\gamma-3) \right] Q_{M-2}^0(r_i)\\
      & &
      = (1-E_M){Q'}_M^0(r_i)
             + \frac{E_M}{r_i}
             \left[ \frac{(M+\gamma-\beta-2)(M+\gamma-3)}
                         {M+\beta-1}- (M+\gamma-3) \right] Q_{M-2}^0(r_i)\\
      & &
      = (1-E_M){Q'}_M^0(r_i)
             + \frac{E_M(M+\gamma-3)}{r_i}
             \left[ \frac{(M+\gamma-\beta-2)}
                         {M+\beta-1}- 1\right] Q_{M-2}^0(r_i)\\
      & &
      = (1-E_M){Q'}_M^0(r_i)
             + \frac{E_M(M+\gamma-3)(2M+\gamma-3)}{r_i(M+\beta-1)}
             Q_{M-2}^0(r_i).
    \end{eqnarray*}
    %
    ここで\Deqref{微分漸化式2}から $r_i\neq 1$ について
    %
    \begin{eqnarray*}
      & &
          (1-r_i^2)r_i{Q'}_M^0(r) \\
      & &
       = \left(-Mr_i^2 + \frac{M(M+\gamma-\beta-2)}
                           {2M+\gamma-3}\right)Q_M^0(r)
         + \frac{(M+\gamma-\beta-2)(M+\gamma-3)}{2M+\gamma-3}Q_{M-2}^0 \\
      & & 
       = \left(-Mr_i^2 + \frac{M(M+\gamma-\beta-2)}
                           {2M+\gamma-3}\right)E_MQ_{M-2}^0(r)
         + \frac{(M+\gamma-\beta-2)(M+\gamma-3)}{2M+\gamma-3}Q_{M-2}^0 \\
      & &
       = -r_i^2 M E_M Q_{M-2}^0(r)
        + \left[ \frac{M(M+\gamma-\beta-2)}
                      {2M+\gamma-3}E_M
                +\frac{(M+\gamma-\beta-2)(M+\gamma-3)}
                      {2M+\gamma-3}                    \right]Q_{M-2}^0 \\
      & &
       = -r_i^2 M E_M Q_{M-2}^0(r)
        + \frac{M+\gamma-\beta-2}{2M+\gamma-3}
                \left[ M E_M + M+\gamma-3 \right]Q_{M-2}^0 \\
      & &
       = -r_i^2 M E_M Q_{M-2}^0(r)
        + \frac{M+\gamma-\beta-2}{2M+\gamma-3}
                \left[ M \frac{(M+\gamma-\beta-2)(M+\gamma-3)}
                              {M(M+\beta-1)} + M+\gamma-3 \right]Q_{M-2}^0 \\
      & &
       = -r_i^2 M E_M Q_{M-2}^0(r)
        + \frac{(M+\gamma-\beta-2)(M+\gamma-3)}{2M+\gamma-3}
                \left[ \frac{(M+\gamma-\beta-2)}
                              {M+\beta-1} + 1 \right]Q_{M-2}^0 \\
      & & 
       = -r_i^2 M E_M Q_{M-2}^0(r)
        + \frac{(M+\gamma-\beta-2)(M+\gamma-3)(2M+\gamma-3)}
               {(2M+\gamma-3)(M+\beta-1)} Q_{M-2}^0 \\
      & &
       = -r_i^2 M E_M Q_{M-2}^0(r)
        + \frac{(M+\gamma-\beta-2)(M+\gamma-3)}
               {M+\beta-1} Q_{M-2}^0 \\
      &=& -r_i^2 M E_M Q_{M-2}^0(r)+ ME_M Q_{M-2}^0 
       =  (1-r_i^2) M E_M Q_{M-2}^0(r),\\
      & &
      {Q'}_M^0(r) = \frac{M E_M}{r_i} Q_{M-2}^0(r).
    \end{eqnarray*}
    %
    が成り立つので, 
    \begin{eqnarray*}
      & &
      \pi'(r) = (1-E_M)\frac{M E_M}{r_i} Q_{M-2}^0(r)
             + \frac{E_M(M+\gamma-3)(2M+\gamma-3)}{r_i(M+\beta-1)}
             Q_{M-2}^0(r_i) \\
      & &
         = \frac{E_M}{r_i}Q_{M-2}^0(r)
               \left[ (1-E_M)M
                  + \frac{(M+\gamma-3)(2M+\gamma-3)}{M+\beta-1}
                \right] \\
      & &
         = \frac{E_M}{r_i}Q_{M-2}^0(r)
               \left[ M-\frac{(M+\gamma-\beta-2)(M+\gamma-3)}{M+\beta-1}
                 + \frac{(M+\gamma-3)(2M+\gamma-3)}{M+\beta-1}
                \right]\\
      & &
         = \frac{E_M}{r_i}Q_{M-2}^0(r)
               \left[ M+\frac{M+\gamma-3}{M+\beta-1}[-(M+\gamma-\beta-2)
                 + (2M+\gamma-3)]
                \right]\\
      & &
         = \frac{E_M}{r_i}Q_{M-2}^0(r)
               \left[ M+\frac{(M+\gamma-3)(M+\beta-1)}{M+\beta-1}
                \right]\\
      & &
         = \frac{E_M(2M+\gamma-3)}{r_i}Q_{M-2}^0(r). 
    \end{eqnarray*}
    %
    よってガウスの重み係数は\Deqref{ガウス積分重み係数定義}より
    %
    \begin{eqnarray*}
      & &
      w_i \equiv \frac{1}{\Ddsty\left.\frac{\pi'(r)}{2r}\right|_{r_i}}
                 \int_0^1 \frac{\pi(r)}{r^2-r_i^2} w(r) dr 
      = \frac{2r_i^2}{E_M(2M+\gamma-3)Q_{M-2}^0(r)}
        \frac{I_{M-2}^0}{C_{M-2} Q_{M-2}^0(r_i)} \\
      & &
      = \frac{2r_i^2}{E_M C_{M-2}(2M+\gamma-3)}
        \frac{I_{M-2}^0}{[Q_{M-2}^0(r_i)]^2} \\
      & &
      = \frac{2r_i^2}{2M+\gamma-3}
        \frac{M(M+\beta-1)}{(M+\gamma-\beta-2)(M+\gamma-3)}
        \frac{(2M+\gamma-5)(2M+\gamma-3)}{M(M+\beta-1)}
        \frac{1}{[\Phi_{M-2}^0(r_i)]^2}\\
      & &
      = \frac{2r_i^2(2M+\gamma-5)}
             {(M+\gamma-\beta-2)(M+\gamma-3)}
        \frac{1}{[\Phi_{M-2}^0(r_i)]^2}.
    \end{eqnarray*}
  }
  $r_{M/2}=1$ に対応する重み係数 $w_{M/2}$ は関係式
  %
  \begin{equation}
    w_{M/2} = I_0^0 - \sum_{i=1}^{M/2-1}w_i,
  \end{equation}
  %
  から計算される
  \footnote{
    \Deqref{スペクトル正変換ガウス積分}において, $n=m=0$, 
    $f_m(r) = Q_0^0(r)=1$ と選ぶと
    \begin{displaymath}
      \int_0^1 Q_0^0\Phi_0^0(r)w(r) dr 
      = \sum_{i=1}^{M/2} Q_0^0(r_i)\Phi_0^0(r_i)w_i, 
      \qquad
      \sqrt{I_0^0}= \sum_{i=1}^{M/2} \frac{w_i}{\sqrt{I_0^0}}, \qquad
      I_0^0 = \sum_{i=1}^{M/2} w_i. 
    \end{displaymath}
  }. 

  $f_m(r)\Phi_n^m(r)$ の次数が $\pi(r)q(r)$ と同じであるので
  \Deqref{ガウスラダウπ}の形からガウス--ラダウ積分は
  $f_m(r)\Phi_n^m(r)$ の次数が $2M-4$ 以下であるときに正確となる. 
  $2M-2$ の程度の精度をのぞむためには格子点の数を 1 つ増やせばよい. 
  すると積分の精度は $2M$ となる. 
  
  2 つの関数の積はガウス積分の $M/2$ 個の格子点あるいはガウス--ラダウ積分の
  $M/2+1$ 個の格子点を選点に選ぶことで効率的に計算できる. 
  非エイリアシング化は 3/2 則を用いて実現できる. 

  数値積分を実行するためにはそれぞれの選点における $Q_n^m(r)$ の値を
  保管しておく必要がある. このことは $O(M^3)$ のメモリーを必要とし
  典型的な 2 次元計算ではこれが最大の必要メモリーとなる. 
  しかしながら同じ状況が球面調和函数のルジャンドル陪函数を用いるときに
  生じる. \Deqref{スペクトル逆変換}, \Deqref{スペクトル正変換ガウス積分}
  の変換は本質的に行列の演算であり FFT ほど効率的でない. 
  しかしながら, 行列が十分に大きければ我々の基底関数に対して
  高速変換を適用できる. 擬スペクトル法に関しては複極展開の適用が
  もう一つの選択肢である. 

\section{作用素}

  微分やその他の演算を楽な方法で計算するためには
  漸化関係式が望ましい. 
  ここでは基本的な演算である $r^2$, $r(d/dr)$ とその逆演算の
  漸化関係式を示す.(そしてヘルムホルツ演算子とその逆演算の手順を示す.)
  さて
  %
  \begin{equation}
    \Deqlab{元関数}
    g_m(r) = \sum_{n=|m|,n+m even}^\infty a_n^m Q_n^m(r),
  \end{equation}
  %
  と, ある線形作用素 $L$ があって, 
  %
  \begin{equation}
    \Deqlab{作用関数}
    Lg_m(r) = \sum_{n=|m|,n+m even}^\infty b_n^m Q_n^m(r).
  \end{equation}
  %
  ここで後に便利になるように, $n> \hat{M}$ に対して $a_n^m \equiv 0$ を
  仮定する. ただし $\hat{M}$ は\Deqref{スペクトル逆変換}で定義したものと
  同様である. 

  $L=r(d/dr)$ を考える. \Deqref{元関数}を\Deqlab{作用関数}に代入し, 
  \Deqref{微分漸化式1}を繰りかえし用いると
  \footnote{
    \Deqref{微分漸化式1}より
    %
    \begin{eqnarray*}
      & &
      \sum_{n=|m|,n+m even}^\infty b_n^m Q_n^m
      = r\DD{g}{r} = \sum_{n=|m|,n+m even}^\infty a_n^m r\DD{Q_n^m}{r}\\
      & &
      = \sum_{n=|m|,n+m even}^\infty a_n^m 
             \left(r\DD{Q_{n-2}^m}{r} + nQ_n^m + (n+\gamma-3)Q_{n-2}^m\right)\\
      & &
      = \sum_{n=|m|,n+m even}^\infty a_n^m [nQ_n^m + (n+\gamma-3)Q_{n-2}^m]\\
      & &
        +\sum_{n=|m|,n+m even}^\infty a_n^m 
             \left(r\DD{Q_{n-4}^m}{r} + (n-2)Q_{n-2}^m 
                   + (n+\gamma-5)Q_{n-4}^m\right)\\
      & &
      = \sum_{n=|m|,n+m even}^\infty a_n^m [nQ_n^m + (n+\gamma-3)Q_{n-2}^m]\\
      & &
      + \sum_{n=|m|,n+m even}^\infty a_n^m 
            [(n-2)Q_{n-2}^m + (n+\gamma-5)Q_{n-4}^m\\
      & &
        +\sum_{n=|m|,n+m even}^\infty a_n^m 
             \left(r\DD{Q_{n-6}^m}{r} + (n-4)Q_{n-4}^m 
                   + (n+\gamma-7)Q_{n-6}^m\right)\\
      & &
        + \ldots\\
      & &
      = \sum_{n=|m|,n+m even}^\infty a_n^m 
           [nQ_n^m + (2n+\gamma-5)Q_{n-2}^m + (2n+\gamma-9)Q_{n-4}^m
            + (2n+\gamma-13)Q_{n-6}^m + \ldots]\\
      & &
      = \sum_{n=|m|,n+m even}^\infty a_n^m nQ_n^m 
      + \sum_{n=|m|-2,n+m even}^\infty (2n+\gamma-1)a_{n+2}^m Q_n^m \\
      & &
      + \sum_{n=|m|-4,n+m even}^\infty (2n+\gamma-1)a_{n+4}Q_n^m
      + \sum_{n=|m|-6,n+m even}^\infty (2n+\gamma-1)a_{n+6}Q_n^m
      + \ldots
    \end{eqnarray*}
    %
    両辺に $Q_{n'}^m w(r)$ をかけて積分することにより
    %
    \begin{displaymath}
      b_{n'}^m = n a_{n'}^m + (2n'+\gamma -1)\sum_{p=n'+2,even}^\infty a_p.
    \end{displaymath}
  }
  %
  \begin{equation}
    b_n^m = n a_n^m + (2n + \gamma -1)\sum_{p=n+2,p+n=even}^\infty a_p^m.
  \end{equation}
  %
  $b_{n+2}^m$ を得るために添え字 $n$ をシフトし, 和の項を取り除くと
  $L=r(d/dr)$ についての漸化式が得られる
  \footnote{
    \begin{displaymath}
      b_n^m = n a_n^m + (2n+\gamma-1)\sum_{p=n+2,p+n=even}^\infty a_p^m.
    \end{displaymath}
    の $n$ を $n+2$ に置き換えた式は
    \begin{displaymath}
      b_{n+2}^m 
      = (n+2) a_{n+2}^m + (2n+\gamma+3)\sum_{p=n+4,p+n=even}^\infty a_p^m.
    \end{displaymath}
    %
    上の式に $2n+\gamma+3$,  下の式に $2n+\gamma-1$ をかけて引き算すれば
    %
    \begin{eqnarray*}
      & &
      (2n+\gamma+3)b_n^m - (2n+\gamma-1)b_{n+2}^m \\
      & &
      = n (2n+\gamma+3) a_n^m - (2n+\gamma-1)(n+2) a_{n+2}^m
        + (2n+\gamma-1)(2n+\gamma+3)a_{n+2}^m,\\
      & &
      (2n+\gamma+3)b_n^m = (2n+\gamma-1)b_{n+2}^m 
       + n (2n+\gamma+3) a_n^m + (2n+\gamma-1)(n+\gamma+1)a_{n+2}^m,\\
      & &
      b_n^m = \frac{2n+\gamma-1}{2n+\gamma+3}b_{n+2}^m 
       + na_n^m + \frac{(2n+\gamma-1)(n+\gamma+1)}{2n+\gamma+3}a_{n+2}^m.
    \end{eqnarray*}
  }. 
  %
  \begin{equation}
    \Deqlab{r(d/dr)漸化式}
    b_n^m = \frac{2n+\gamma-1}{2n+\gamma+3}b_{n+2}^m
          + \frac{(2n+\gamma-1)(n+\gamma+1)}{2n+\gamma+3}a_{n+2}^m
          + n a_n^m.
  \end{equation}
  %
  $b_{\hat{M}+2}^m \equiv 0$ を初期値として逆方向に数値的に
  安定に計算して $b_n^m$ について解くことができる. 
  $r(d/dr)$ の逆演算のための漸化式は,
  \Deqref{r(d/dr)漸化式}を $a_n^m$ について逆に解くことにより得られる. 
  この漸化式は定めることのできない $a_0^0$ を除いて, 
  $a_{\hat{M}+2}^m \equiv 0$ を初期値として安定に計算できる. 
  $a_0^0$ の項は $r(d/dr)^{-1}$ の積分定数として供される. 

  もしも $L=r^2$ ならば, \Deqref{作用関数}とともに
  \Deqref{漸化式}と\Deqref{元関数}を用いて, $n \ge |m|$ に対して
  \footnote{
    $g_m(r) = \sum_{n=|m|,n+m even}^\infty a_n^m Q_n^m(r)$ に対して
    $r^2$ を作用させると
    %
    \begin{eqnarray*}
      & &
      Lg_m(r) = \sum_{n=|m|,n+m even}^\infty b_n^m Q_n^m(r)
      = \sum_{n=|m|,n+m even}^\infty a_n^m r^2Q_n^m(r) \\
      & &
      = \sum_{n=|m|,n+m even}^\infty a_n^m 
        \left[  \frac{(n-|m|+2)(n+|m|+\beta+1)}
                     {(2n+\gamma-1)(2n+\gamma+1)}Q_{n+2}^m \right.\\
      & &\left.
              + \frac{2n(n+\gamma-1)+2|m|(|m|+\beta-1)+(\gamma-3)(\beta+1)}
                     {(2n+\gamma-3)(2n+\gamma+1)}Q_n^m 
              + \frac{(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)}
                     {(2n+\gamma-1)(2n+\gamma-3)}Q_{n-2}^m \right] \\
      & &
      = \sum_{n=|m|,n+m even}^\infty 
        \left[  \frac{(n-|m|+2)(n+|m|+\beta+1)}
                     {(2n+\gamma-1)(2n+\gamma+1)}a_n^m Q_{n+2}^m \right.\\
      & &
              + \frac{2n(n+\gamma-1)+2|m|(|m|+\beta-1)+(\gamma-3)(\beta+1)}
                     {(2n+\gamma-3)(2n+\gamma+1)}a_n^m Q_n^m \\
      & &\left.
              + \frac{(n-|m|+\gamma-\beta-2)(n+|m|+\gamma-3)}
                     {(2n+\gamma-1)(2n+\gamma-3)}a_n^m Q_{n-2}^m \right]
    \end{eqnarray*}
    %
    両辺に $Q_{n'}^m w(r)$ をかけて 0 から 1 まで積分すると, 
    直交関係から
    \begin{eqnarray*}
      & &
      b_n^m =  \frac{(n-|m|)(n+|m|+\beta-1)}
                    {(2n+\gamma-5)(2n+\gamma-3)}a_{n-2}^m \\
      & &
              + \frac{2n(n+\gamma-1)+2|m|(|m|+\beta-1)+(\gamma-3)(\beta+1)}
                     {(2n+\gamma-3)(2n+\gamma+1)}a_n^m \\
      & &
              + \frac{(n-|m|+\gamma-\beta)(n+|m|+\gamma-1)}
                     {(2n+\gamma+3)(2n+\gamma+1)}a_{n+2}^m
    \end{eqnarray*}
  }, 
  %
  \begin{eqnarray}
    \Deqlab{r^2漸化式}
    b_n^m &=& 
          \frac{(n-|m|)(n+|m|+\beta-1)}
               {(2n+\gamma-5)(2n+\gamma-3)}a_{n-2}^m \nonumber\\
          & &+ \frac{2n(n+\gamma-1)+2|m|(|m|+\beta-1)+(\gamma-3)(\beta+1)}
                    {(2n+\gamma+1)(2n+\gamma-3)} a_n^m \nonumber\\
          & &+ \frac{(n-|m|+\gamma-\beta)(n+|m|+\gamma-1)}
                    {(2n+\gamma+3)(2n+\gamma+1)} a_{n+2}^m,
  \end{eqnarray}
  %
  ただし $a_{|m|-2}^m\equiv 0$ である. 
  もしも $\gamma=3$ で $n=m=0$ ならば, 右辺の第 2 項は 
  $(\beta+1)/(\gamma+1)a_n^m$ とする. 
  \Deqref{r^2漸化式}によって $|m|\le n \le \hat{M}+2$ に対する $b_n^m$ を
  計算することができる. 
  $1/r^2$ 作用素に対する逆方向の漸化式は\Deqref{r^2漸化式}を $a_{n-2}^m$ に
  ついて解くことにより得られる. 
  初期値は $a_{M+2}^m \equiv 0$, $ a_{M+4}^m\equiv 0$ である. 
  しかしながらこの漸化式は数値的に不安定である. 
  $a_n^m$ は\Deqref{r^2漸化式}をピボットなしの 3 重対角行列の LU 分解して
  解くことにより安定に, かつ安価に求めることができる. 
  ここで大事なことは $b_n^m$ で表わされる関数が $r\rightarrow 0$ にて
  $O(r^{|m|+2p+2})$ で振る舞うことである. ここで $p$ は整数である. 
  でなければ結果の $g_m(r)$ が座標の特異点のまわりで不完全な振る舞いをする. 

\appendix

\section{ラグランジュ補間と積分}

  区間 [0,1] の中の格子点 $[r_1,\ldots,r_n]$ 上での関数 $f(r)$ の値
  $f(r_0),f(r_1),\ldots, f(r_n)$ が与えられているときに 
  任意の点 $r$ での $f(r)$ の値を $r^{2n-2}$ までの $r^2$ の
  多項式で表し, ガウス積分の評価を行う. 

  まず, Hilderbrand (1974) にしたがって
  格子点 $[x_0,x_1,\ldots,x_n]$ 上での関数 $f(x)$ の値
  $f(x_0),f(x_1),\ldots, f(x_n)$ が与えられているときに 
  任意の点 $x$ での $f(x)$ の値を $x^n$ までの多項式で
  補間することを考える(ラグランジュ補間法). 
  すなわち $n$ 次の多項式 $l_k(x)$ があって
  %
  \begin{equation}
    f(x) = \sum_{k=0}^n l_k(x) f(x_k)
  \end{equation}
  %
  多項式 $l_k(x)$ を見出すには, 各格子点での値が正確であることから
  求められる. すなわち
  %
  \begin{equation}
    l_k(x_k) = 1, \quad k=0,\ldots,n, \qquad
    l_k(x_l) = 0, \quad l\neq k, k=0,\ldots,n.
  \end{equation}
  %
  あるいはクロネッカーのデルタ記号を用いて
  %
  \begin{equation}
    l_k(x_l) = \delta_{kl}.
  \end{equation}
  %
  $l_k(x)$ がせいぜい $x$ の $n$ 次多項式であることと $k\neq l$ の
  条件から
  %
  \begin{equation}
    l_k(x) = C_k (x-x_0)\cdots(x-x_{k-1})(x-x_{k+1})\cdots(x-x_n)
           = C_k \frac{\Ddsty\prod_{l=0}^n (x-x_l)}{x-x_k}
  \end{equation}
  %
  と書ける. ただし $C_k$ は定数である. 
  $C_k$ を求めるには $l_k(x_k)=1$ から求められる. 
  %
  \begin{eqnarray*}
    & &
    l_k(x_k) = C_k(x_k-x_0)\cdots(x_k-x_{k-1})(x_k-x_{k+1})\cdots(x_k-x_n)=1,\\
    & &
    C_k =\frac{1}{(x_k-x_0)\cdots(x_k-x_{k-1})(x_k-x_{k+1})\cdots(x_k-x_n)}.
  \end{eqnarray*}
  %
  よって
  \begin{equation}
    l_k(x) =\frac{(x-x_0)\cdots(x-x_{k-1})(x-x_{k+1})\cdots(x-x_n)}
                 {(x_k-x_0)\cdots(x_k-x_{k-1})(x_k-x_{k+1})\cdots(x_k-x_n)}.
  \end{equation}
  ここで
  %
  \begin{displaymath}
    \pi(x) = (x-x_0)\cdots(x-x_k)\cdots(x-x_n)=\prod_{l=0}^n (x-x_l)
  \end{displaymath}
  %
  を用いると
  \begin{displaymath}
    \pi'(x) = \sum_{j=0}^n \frac{\pi(x)}{x-x_j}
  \end{displaymath}
  %
  であるから
  %
  \begin{displaymath}
    \pi'(x_k) = \left.\frac{\pi(x)}{x-x_k}\right|_{x=x_k}
              = (x_k-x_0)\cdots(x_k-x_{k-1})(x_k-x_{k+1})\cdots(x_k-x_n).
  \end{displaymath}
  %
  よって
  %
  \begin{equation}
    C_k =\frac{1}{\pi'(x_k)},
  \end{equation}
  %
  と書くことができて, 結局
  %
  \begin{equation}
    l_k(x) = \frac{\pi(x)}{(x-x_k)\pi'(x_k)}, 
    \quad \pi(x)=\prod_{l=0}^n (x-x_l).
  \end{equation}

  $f(r)$ が $r^2$ の多項式のときには, $x\rightarrow r^2$ に置き換えて
  %
  \begin{eqnarray}
    & &
    f(r) = \sum_{k=0}^n l_k(r) f(r_k), \\
    & &
    l_k(r) = \frac{(r^2-r_0^2)\cdots(r^2-r_{k-1}^2)
                   (r^2-r_{k+1}^2)\cdots(r^2-r_n^2)}
                  {(r_k^2-r_0^2)\cdots(r_k^2-r_{k-1}^2)
                   (r_k^2-r_{k+1}^2)\cdots(r_k^2-r_n^2)}
  \end{eqnarray}
  %
  ここで $\pi(r)$ を
  %
  \begin{equation}
    \pi(r) = (r^2-r_0^2)\cdots(r^2-r_k^2)\cdots(r^2-r_n^2)
           = \prod_{k=0}^n (r^2-r_k^2)
  \end{equation}
  %
  と定義すると, この場合微分は
  %
  \begin{displaymath}
    \pi'(r) = 2r\sum_{j=0}^n \frac{\pi(r)}{r^2-r_j^2}
  \end{displaymath}
  %
  であるから
  %
  \begin{displaymath}
    \left.\frac{\pi'(r)}{2r}\right|_{r_k} 
           = \left.\frac{\pi(r)}{r^2-r_k^2}\right|_{r=r_k}
           = (r_k^2-r_0^2)\cdots(r_k^2-r_{k-1}^2)
             (r_k^2-r_{k+1}^2)\cdots(r_k^2-r_n^2).
  \end{displaymath}
  %
  したがって
  %
  \begin{equation}
    l_k(r) = \frac{\pi(r)}
                  {\Ddsty(r^2-r_k^2)\left.\frac{\pi'(r)}{2r}\right|_{r_k}}.
  \end{equation}

\begin{Dreference}
  \item Matsushima, T., Marcus, P. S., 1995:
    A spectral method for polar coordinates.
    {\it J. Comp. Phys.} {\bf 120}, 365--374

  \item Hilderbrand, F. B., 1974:
    Introduction to numerical analysis. 
    {\it McGraw-Hill}, 669pp.

  \item 石岡圭一, 2004:
    スペクトル法による数値計算入門, 
    東京大学出版会, 230pp.
\end{Dreference}

\end{document}

%%% Local Variables: 
%%% mode: latex
%%% TeX-master: t
%%% End: