pianofisica

Mathematics & Physics, Maxima, a bit Python & Wolfram, and Arts

入力例で学ぶMaximaの使い方(応用:多項式の商と余り)

2つの多項式同士の割り算を考えたとき、それらの次数が低くて係数も簡単な場合には良いですが、

一般にその商と余りを手計算で求めるのはあまり現実的ではありません。

そこでMaximaなどの数式処理ソフトをコンピュータで利用する計算が大変役に立ちます。

Maximaの使い方の簡単な紹介については以下の記事を参考にしてください:

pianofisica.hatenablog.com



多項式の割り算

2つの多項式

\qquad\displaystyle{
P(x)=a_0x^p+a_1x^{p-1}+\dots+a_p \\ 
Q(x)=b_0x^q+b_1x^{q-1}+\dots+b_q
}

を考え、\ p>q としたとき、これらの割り算の商 \ H(x) と余り\ R(x)

\qquad\displaystyle{
P(x)=H(x)Q(x)+R(x)
}

ただし

\qquad\displaystyle{
H(x)=h_0x^{p-q}+h_1x^{p-q-1}+\dots+h_{p-q} \\ 
R(x)=r_0x^{q-1}+r_1x^{q-2}+\dots+r_{q-1}
}

をみたすような多項式として求まります。

これらを求める一連の関数をMaximaで構成してみたいと思います。

具体的な計算のために

\qquad\displaystyle{
P(x)=\left(\frac{x^2}{2}+3x+4\right)^5 \\ 
Q(x)=5x^3+\frac{2}{3}x+7
}

という場合で考えてみます。これらの商と余りを手計算で求めようとはなかなか思えないですよね。

以下のMaximaのソースの中では、多項式の入力を適宜変更してもらえば、

それに応じた商と余りの答えが出力されるように作ってあります。

(ただし、30次以下の多項式に限定しています)

多項式の次数

まず注意する点は、商と余りの多項式 \ H(x),\ R(x) の次数が、

入力した多項式 \ P(x),\ Q(x) の次数 \ p,\ q に応じて決定されるということです。

そこではじめに、勝手に入力した多項式について、その次数を出力する関数を定義します:

Ord(P, x):=block(n:1, while n<31 do(
c: ratcoeff(P, x, 30-n),
k[n]: if c = 0 then -1 else 30-n,
n: n+1),
lmax(makelist(k[n],n, 1,30)) )$

どのような操作をしているかというと、入力した x の多項式 P について、

その  n 次の係数がゼロでない場合に、その次数  n を要素に持つようなリスト k を作成し、

リスト k の要素の中で最大の値を出力する関数として Ord を定義しています。

多項式の入力

2つの多項式を入力し、それらの次数をみてみましょう:

P: (x^2/2+3*x+4)^5$
p:Ord(P, x);
Q: 5*x^3+2/3*x+7$
q:Ord(Q, x);

それぞれ 10、3 と出力され、正しく多項式の次数を与えていることが確認できます。

商と余りの多項式の準備

Maximaに以下のように入力して、未知の係数 h[n], r[n] を持つ形で

H:sum(h[n]*x^(p-q-n), n, 0, p-q);
R:sum(r[n]*x^(q-1-n), n, 0, q-1 );

として商と余りの多項式を準備しておきます。

商と余りの多項式の決定

条件\ P(x)=H(x)Q(x)+R(x) から、商と余りの多項式にある未知係数を決定します:

S:solve(makelist( ratcoeff(P,x,n)=ratcoeff(H*Q+R,x,n), n, 0, p),
append(makelist(h[n],n,0,p-q), makelist(r[n],n,0,q-1)));

その答えを代入することにより

H:sum(rhs(S[1][1+n])*x^(p-q-n), n, 0, p-q );
R:sum(rhs(S[1][2+p-q+n])*x^(q-1-n), n, 0, q-1);

から

\quad\displaystyle{
\begin{array}{l}
H(x)=\frac{1}{160}x^7+\frac{3}{16}x^6+\frac{2999}{1200}x^5+\frac{15573}{800}x^4+\frac{1762277}{18000}x^3+\frac{991517}{3000}x^2+\frac{403754407}{540000}x+\frac{96017993}{90000} \\ 
R(x)=\frac{2906278463}{810000}x^2-\frac{126305869}{60000}x-\frac{579965951}{90000} 
\end{array} 
}

と求まります。手計算でこれを求めるのは、ちょっと現実的でないことがわかりますね。

検算

念のために検算として確かに\ P(x)=H(x)Q(x)+R(x) が成り立っているかチェックしておきます:

is(expand(P)=expand(expand(Q*H)+R));

 {}

さいごに

さいごに、Maximaの入力を改めてひとまとめにして置いておきます:

kill(all)$

P: (x^2/2+3*x+4)^5$ 
Q: 5*x^3+2/3*x+7$

Ord(P, x):=block(n:1, while n<31 do(
c: ratcoeff(P, x, 30-n), k[n]: if c = 0 then -1 else 30-n, n: n+1),
lmax(makelist(k[n],n, 1,30)) )$
p:Ord(P, x)$ q:Ord(Q, x)$
H:sum(h[n]*x^(p-q-n), n, 0, p-q)$
R:sum(r[n]*x^(q-1-n), n, 0, q-1 )$
S:solve(makelist( ratcoeff(P,x,n)=ratcoeff(H*Q+R,x,n), n, 0, p),
append(makelist(h[n],n,0,p-q), makelist(r[n],n,0,q-1)))$
H:sum(rhs(S[1][1+n])*x^(p-q-n), n, 0, p-q );
R:sum(rhs(S[1][2+p-q+n])*x^(q-1-n), n, 0, q-1);

入力する多項式 P 、Q を変えれば、自動的にその商と余りが求まります。


以上、今回取り扱ったのは非常にシンプルな問題ですが、

コンピュータを用いる計算の有用性がわかるのではないでしょうか。


関連記事では、Maximaの使い方のまた別な応用編として、

数学や物理学のさまざまな具体的な問題の計算にMaximaを使っています。

Maximaの入力も載せているので、参考にしてみてください)


キーワードMaxima多項式、商、余り