|
|
多項式を使って数学の問題を解く
!問題
!3次曲線 y=x^3+x^2-5x-7 が直線Lと3点A,B,Cで交わっている。
!A,Bのx座標がそれぞれ -1, 2であるとき、
!直線Lの方程式および3点A,B,Cの座標を求めよ。
OPTION ARITHMETIC RATIONAL !有理数モード
PUBLIC NUMERIC MAX_DEGREE !最高次数
LET MAX_DEGREE=20
DATA 3 !y=x^3+x^2-5x-7
DATA -7,-5,1,1
READ ff !yを読み込む
DIM F(0 TO ff)
MAT READ F
LET gg=2 !解と係数の関係より、a,bを解とする2次方程式を得る
DIM G(0 TO gg)
LET a=-1 !2点A,B
LET b=2
LET G(2)=1
LET G(1)=-(a+b)
LET G(0)=a*b
DIM Q(0 TO MAX_DEGREE),R(0 TO MAX_DEGREE)
CALL PolynomialQuotientRemainder(ff,F,gg,G, qq,Q,rr,R)
PRINT qq !点C x+2=0
MAT PRINT Q;
PRINT rr !直線L y=-x-3
MAT PRINT R;
PRINT a; PolynomialValue(a,ff,F) !3点A,B,C
PRINT b; PolynomialValue(b,ff,F)
PRINT -Q(0); PolynomialValue(-Q(0),ff,F)
END
!POLY.LIB
!f(x)=A[aa]x^aa+ … +A[2]x^2+A[1]x+A[0]、g(x)=B[bb]x^bb+ … +B[2]x^2+B[1]x+B[0]
EXTERNAL FUNCTION PolynomialValue(a,ff,F()) !多項式の値 f(a)
OPTION ARITHMETIC RATIONAL !有理数モード
LET t=F(ff)
FOR i=ff-1 TO 0 STEP -1 !ホーナー法
LET t=t*a+F(i)
NEXT i
LET PolynomialValue=t
END FUNCTION
EXTERNAL SUB PolynomialDifferential(K,ff,F(), cc,C()) !k階微分
OPTION ARITHMETIC RATIONAL !有理数モード
MAT C=ZER
FOR i=K TO ff
LET C(i-K)=PERM(i,K)*F(i)
NEXT i
LET cc=MAX(ff-K,0)
END SUB
EXTERNAL SUB PolynomialPowerN(aa,A(), N, bb,B()) !べき乗展開（多項定理による）
OPTION ARITHMETIC RATIONAL !有理数モード
DIM P(aa+1) !{p,q,…,r}の並び
MAT B=ZER
CALL PolyPowN(N,aa+1,1,P, aa,A,N, B)
LET bb=aa*N
END SUB
EXTERNAL SUB PolyPowN(m,d,s,P(), aa,A(),N, B()) !自然数mをd分割する
OPTION ARITHMETIC RATIONAL !有理数モード
IF d=1 THEN
LET P(s)=m
!{p,q,…,r}で多項式を展開する
LET c=PermFactorialM(P,aa+1) !前半部分 (p+q+ … +r)!/(p!*q!* … *r!)
LET x=0 !後半部分を加味する
FOR i=0 TO aa
LET c=c*A(i)^P(i+1) !係数
LET x=x+i*P(i+1) !べき乗
NEXT i
LET B(x)=B(x)+c !記録する
ELSEIF d>1 THEN
FOR i=0 TO m
LET P(s)=m-i
CALL PolyPowN(i,d-1,s+1,P, aa,A,N, B) !次へ
NEXT i
END IF
END SUB
EXTERNAL FUNCTION PermFactorialM(B(),M) !同じものを含む順列の「場合の数」
OPTION ARITHMETIC RATIONAL !有理数モード
LET s=B(M) !総数 r, … ,q+ … +r,p+q+ … +r
LET t=1 !組合せ comb(r,r), … ,comb(q+ … +r,q),comb(p+q+ … +r,p)
FOR i=M-1 TO 1 STEP -1
LET s=s+B(i)
LET t=t*COMB(s,B(i)) !組合せ順列
NEXT i
LET PermFactorialM=t
END FUNCTION
EXTERNAL SUB PolynomialComposition(aa,A(),bb,B(), ss,S()) !合成関数 f(g(x))
OPTION ARITHMETIC RATIONAL !有理数モード
DIM W(0 TO aa*bb)
LET S(0)=A(aa) !s=a[aa]
LET ss=0
FOR i=aa-1 TO 0 STEP -1 !ホーナー法
CALL PolynomialMultiply(ss,S,bb,B, ww,W) !s=s*X+a[i]
LET W(0)=W(0)+A(i)
MAT S=W !次へ
LET ss=ww
NEXT i
END SUB
EXTERNAL SUB PolynomialGCD(aa,A(),bb,B(), cc,C()) !最大公約数
OPTION ARITHMETIC RATIONAL !有理数モード
IF aa=0 AND bb=0 THEN !定数項のみなら
MAT C=ZER
LET C(0)=GCD(A(0),B(0))
LET cc=0
ELSE
DIM TA(0 TO MAX_DEGREE),TB(0 TO MAX_DEGREE) !作業変数
MAT TA=A
LET taa=aa
MAT TB=B
LET tbb=bb
DO WHILE NOT(tbb=0 AND TB(0)=0) !--- DO WHILE b<>0
DIM Q(0 TO MAX_DEGREE),R(0 TO MAX_DEGREE)
CALL PolynomialQuotientRemainder(taa,TA,tbb,TB, qq,Q,rr,R) !--- LET R=MOD(a,b)
MAT TA=TB !--- LET a=b
LET taa=tbb
MAT TB=R !--- LET b=R
LET tbb=rr
LOOP !--- LOOP
LET G=TA(0) !既約
FOR i=1 TO taa
LET G=GCD(G,TA(i))
NEXT i
MAT C=(1/G)*TA !--- LET GCD=a
LET cc=taa
END IF
END SUB
EXTERNAL SUB PolynomialLCM(aa,A(),bb,B(), cc,C()) !最小公倍数
OPTION ARITHMETIC RATIONAL !有理数モード
DIM G(0 TO MAX_DEGREE)
CALL PolynomialGCD(aa,A,bb,B, gg,G) !LCM=A*B/Gより
DIM Q(0 TO MAX_DEGREE),R(0 TO MAX_DEGREE)
CALL PolynomialQuotientRemainder(bb,B,gg,G, qq,Q,rr,R)
IF NOT(rr=0 AND R(0)=0) THEN !余りが0以外なら
PRINT "論理エラー"
STOP
END IF
CALL PolynomialMultiply(aa,A,qq,Q, cc,C)
END SUB
!拡張ユークリッド互除法
! f(x)=A[m]x^m+ … +A[1]x+A[0]、g(x)=B[n]x^n+ … +B[1]x+B[0]、m≧nとして、
! f(x)S(x)+g(x)T(x)=gcd(f(x),g(x))=C(x)となるS(x),T(x),C(x)を求める。
EXTERNAL SUB PolynomialExtendedGCD(aa,A(),bb,B(), ss,S(),tt,T(),cc,C()) !拡張ユークリッド互除法
OPTION ARITHMETIC RATIONAL !有理数モード
CALL Poly_ExGCD(aa,A,bb,B, ss,S,tt,T,cc,C)
LET G=C(0) !既約
FOR i=1 TO cc
LET G=GCD(G,C(i))
NEXT i
LET G=G*SGN(C(cc)) !C(x)の最高次数の係数は正とする
MAT S=(1/G)*S
MAT T=(1/G)*T
MAT C=(1/G)*C
END SUB
!f(x)S(x)+g(x)T(x)=k*gcd(f(x),g(x))=C(x)となるS(x),T(x),C(x)を求める。
EXTERNAL SUB Poly_ExGCD(aa,A(),bb,B(), ss,S(),tt,T(),cc,C()) !拡張ユークリッド互除法
OPTION ARITHMETIC RATIONAL !有理数モード
IF bb=0 AND B(0)=0 THEN !!--- IF b=0 THEN
!!--- s=1 !※f(x)*1+0*0=f(x)とする
MAT S=ZER
LET S(0)=1
LET ss=0
!!--- t=0
MAT T=ZER
LET T(0)=0
LET tt=0
!!--- c=a
MAT C=A
LET cc=aa
ELSE
!!--- q=INT(a/b), r=MOD(a,b)
DIM Q(0 TO MAX_DEGREE),R(0 TO MAX_DEGREE)
IF aa=0 AND bb=0 THEN !定数項のみ
MAT Q=ZER
LET Q(0)=INT(A(0)/B(0))
LET qq=0
MAT R=ZER
LET R(0)=MOD(A(0),B(0))
LET rr=0
ELSE
CALL PolynomialQuotientRemainder(aa,A,bb,B, qq,Q,rr,R)
END IF
!!--- CALL ExGCD(b,r, u,v,c) !k=n-1,…,3,2 まで続ける
!!--- s=v
CALL PolynomialExtendedGCD(bb,B,rr,R, tt,T,ss,S,cc,C)
!!--- t=u-v*q
DIM W(0 TO MAX_DEGREE)
CALL PolynomialMultiply(ss,S,qq,Q, ww,W)
MAT T=T-W
LET tt=ww
END IF
END SUB
!補助ルーチン
!演算関連
EXTERNAL SUB PolynomialMultiply(aa,A(),bb,B(), ss,S()) !乗算 S=A*B ※S≠A、S≠B
OPTION ARITHMETIC RATIONAL !有理数モード
MAT S=ZER
FOR i=aa TO 0 STEP -1
LET k=A(i)
IF k=0 THEN !係数が0以外なら
ELSE
FOR j=bb TO 0 STEP -1
LET S(i+j)=S(i+j)+k*B(j) !すべての係数をかける
NEXT j
END IF
NEXT i
LET ss=aa+bb!次数 ※その補正
END SUB
EXTERNAL SUB PolynomialQuotientRemainder(aa,A(),bb,B(), qq,Q(),rr,R()) !除算 ※被除数=商*除数+余り
OPTION ARITHMETIC RATIONAL !有理数モード
IF bb=0 AND B(0)=0 THEN !除数が0なら
PRINT "0で割ることはできません。"
STOP
ELSE
MAT Q=ZER !商
MAT R=A !余り
FOR t=aa TO bb STEP -1 !被除数の次数が除数のより大きいなら
IF R(t)=0 THEN !係数が0以外なら
ELSE
LET k=R(t)/B(bb) !商の係数、その次数
LET w=t-bb
LET Q(w)=k !商
FOR i=bb TO 0 STEP -1 !余り ※R=A-k*B
LET R(w+i)=R(w+i)-k*B(i)
NEXT i
END IF
NEXT t
LET qq=MAX(aa-bb,0) !次数
IF aa>=bb THEN LET t=MAX(bb-1,0) ELSE LET t=aa !次数
FOR rr=t TO 1 STEP -1 !※その補正
IF R(rr)<>0 THEN EXIT FOR
NEXT rr
END IF
END SUB
!表示関連
EXTERNAL SUB PolynomialDisplay(aa,A()) !多項式を表示する a(X)=ΣAkX^k=AnX^n+An-1X^n-1+…+A1X+A0
OPTION ARITHMETIC RATIONAL !有理数モード
CALL mono_disp(A(aa),aa) !最初の項
FOR i=aa-1 TO 0 STEP -1 !次項
LET w=A(i)
IF w>0 THEN PRINT "+";
IF w<>0 OR aa=0 THEN CALL mono_disp(w,i)
NEXT i
END SUB
EXTERNAL SUB mono_disp(ak,k) !単項式を表示する Ak*X^k
OPTION ARITHMETIC RATIONAL !有理数モード
IF k<>0 THEN !x^nで
IF ak=0 OR ak=1 THEN !係数が0,1なら
ELSEIF ak=-1 THEN !係数が－１なら
PRINT "-"; !符号
ELSE
PRINT STR$(ak);"*";
END IF
END IF
IF k=0 THEN !次数が0なら
PRINT STR$(ak);
ELSE
IF ak<>0 THEN !係数が0以外なら
PRINT "X";
IF k<>1 THEN PRINT "^";STR$(k); !次数が1以外なら
END IF
END IF
END SUB
実行結果
1
2 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1
-3 -1 0 0
-1 -2
2 -5
-2 -1
※サブルーチン部分は、POLY.LIB として共有利用してください。
●例
!問題
!2個のサイコロを投げるとき、目の和が8となる場合の数は、
!多項式 (x+x^2+x^3+x^4+x^5+x^6)^2 を展開したときの、x^8の係数を調べればよい。
OPTION ARITHMETIC RATIONAL !有理数モード
PUBLIC NUMERIC MAX_DEGREE !最高次数
LET MAX_DEGREE=20
DATA 6 !(x+x^2+x^3+x^4+x^5+x^6)^2
DATA 0,1,1,1,1,1,1
LET N=2 !べき乗
READ aa
DIM A(0 TO aa)
MAT READ A
DIM B(0 TO aa*N) !展開した多項式
CALL PolynomialPowerN(aa,A, N, bb,B)
PRINT bb
MAT PRINT B;
END
MERGE "POLY.LIB"
実行結果
12
0 0 1 2 3 4 5 6 5 4 3 2 1
●例
!問題
!3次方程式f(x)=x^3+x^2-4x+1=0の1つの解をαとする。
!他の2つの解は、-α^2-2α+2, α^2+α-3と表される。
OPTION ARITHMETIC RATIONAL !有理数モード
PUBLIC NUMERIC MAX_DEGREE !最高次数
LET MAX_DEGREE=20
DATA 3 !f(x)=x^3+x^2-4x+1=0
DATA 1,-4,1,1
DATA 2 !-α^2-2α+2
DATA 2,-2,-1
DATA 2 !α^2+α-3
DATA -3,1,1
READ ff !f(x)を読み込む
DIM F(0 TO ff)
MAT READ F
READ gg !g(α)を読み込む
DIM G(0 TO gg)
MAT READ G
DIM S(0 TO MAX_DEGREE)
CALL PolynomialComposition(ff,F,gg,G, ss,S) !f(-α^2-2α+2)を求める
PRINT ss
MAT PRINT S;
DIM Q(0 TO MAX_DEGREE),R(0 TO MAX_DEGREE)
CALL PolynomialQuotientRemainder(ss,S,ff,F, qq,Q,rr,R) !f(α)=0を考慮する
PRINT qq
MAT PRINT Q;
PRINT rr !余りが0となる
MAT PRINT R;
PRINT
!-------------------------------
READ hh !h(α)を読み込む
DIM H(0 TO hh)
MAT READ H
CALL PolynomialComposition(ff,F,hh,H, ss,S) !f(α^2+α-3)を求める
PRINT ss
MAT PRINT S;
CALL PolynomialQuotientRemainder(ss,S,ff,F, qq,Q,rr,R) !f(α)=0を考慮する
PRINT qq
MAT PRINT Q;
PRINT rr !余りが0となる
MAT PRINT R;
END
MERGE "POLY.LIB"
実行結果
6
5 -24 16 20 -5 -6 -1
3
5 -4 -5 -1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0
0 0 0 0 0 0 0
6
-5 17 9 -15 -5 3 1
3
-5 -3 2 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0
0 0 0 0 0 0 0
|
|
