新しく発言する EXIT インデックスへ

多項式逆数

　　多項式逆数 しばっち 2008/03/28 20:30:14 
　　├続き しばっち 2008/03/28 20:30:59 
　　└多項式平方根 しばっち 2008/03/28 20:31:48 
　　　└続き しばっち 2008/03/28 20:32:38 

　　多項式逆数 しばっち 2008/03/28 20:30:14   ツリーへ
多項式逆数 返事を書く ノートメニュー
しばっち <dihjvcfsyu> 2008/03/28 20:30:14
F(X)=A+B*X+C*X^2+D*X^3+...
G(X)=A'+B'X+C'X^2+D'X^3+...
F(X)/G(X)=1 G(X)を求める

恒等式 (A+B*X+C*X^2+D*X^3+E*X^4+...)(A'+B'*X+C'*X^2+D'*X^3+E'*X^4+...)=1

AA'=1
AB'+BA'=0
AC'+BB'+CA'=0
AD'+BC'+CB'+DA'=0
AE'+BD'+CC'+DB'+EA'=0
:

A'=1/A
B'=-BA'/A
C'=-(BB'+CA')/A
D'=-(BC'+CB'+DA')/A
:

PUBLIC NUMERIC MAXLEVEL,MINLEVEL
LET MAXLEVEL=21
LET MINLEVEL=-2
DIM C(MINLEVEL TO MAXLEVEL)
CALL COSINE(C)
PRINT "COS(X)=";
CALL DISPLAY(C)
CALL RECIP(C)
PRINT "1/COS(X)=";
CALL DISPLAY(C)
CALL SINE(C)
PRINT "SIN(X)=";
CALL DISPLAY(C)
CALL RECIP(C)
PRINT "1/SIN(X)=";
CALL DISPLAY(C)
END

EXTERNAL SUB RECIP(C())
DIM XX(MINLEVEL TO MAXLEVEL),X(MINLEVEL TO MAXLEVEL)
LET M=DIMCHECKLOW(C)
LET X(-M)=1/C(M)
FOR K=-M+1 TO MAXLEVEL
CALL CLR(XX)
FOR I=0 TO MAXLEVEL
FOR J=-M TO MAXLEVEL-I
LET XX(I+J)=XX(I+J)+C(I)*X(J)
NEXT J
NEXT I
IF K-M<>-M THEN LET X(K-M)=-XX(K)/C(M)
NEXT K
CALL COPY(C,X)
END SUB

EXTERNAL SUB COPY(X(),Y())
FOR I=MINLEVEL TO MAXLEVEL
LET X(I)=Y(I)
NEXT I
END SUB

EXTERNAL SUB CLR(X())
FOR I=MINLEVEL TO MAXLEVEL
LET X(I)=0
NEXT I
END SUB

EXTERNAL SUB DISPLAY(A())
LET N=DIMCHECKHIGH(A)
IF N>1 THEN
IF A(N)<0 THEN PRINT "-";
IF ABS(A(N))<>1 THEN
PRINT STR$(ABS(A(N)));"*X^";STR$(N);
ELSE
PRINT "X^";STR$(N);
END IF
END IF
FOR I=N-1 TO 2 STEP -1
IF A(I)<>0 THEN
IF A(I)<0 THEN PRINT "-"; ELSE PRINT "+";
IF ABS(A(I))<>1 THEN
PRINT STR$(ABS(A(I)));"*X^";STR$(I);
ELSEIF ABS(A(I))=1 THEN
PRINT "X^";STR$(I);
END IF
END IF
NEXT I
IF A(1)<>0 THEN
IF N>1 THEN
IF A(1)<0 THEN PRINT "-"; ELSE PRINT "+";
END IF
IF ABS(A(1))<>1 THEN
PRINT STR$(ABS(A(1)));"*X";
ELSEIF ABS(A(1))=1 THEN
PRINT "X";
END IF
END IF
IF A(0)<>0 THEN
IF A(0)<0 THEN PRINT "-"; ELSE PRINT "+";
PRINT STR$(ABS(A(0)));
END IF
IF MINLEVEL<0 THEN
FOR I=-1 TO MINLEVEL STEP -1
IF A(I)<>0 THEN
IF A(I)<0 THEN PRINT "-"; ELSE PRINT "+";
IF I=-1 THEN
PRINT STR$(ABS(A(I)));"/X";
ELSE
PRINT STR$(ABS(A(I)));"/X^";STR$(ABS(I));
END IF
END IF
NEXT I
END IF
PRINT
END SUB

EXTERNAL FUNCTION DIMCHECKHIGH(X())
FOR N=MAXLEVEL TO MINLEVEL STEP -1
IF X(N)<>0 THEN EXIT FOR
NEXT N
LET DIMCHECKHIGH=N
END FUNCTION

EXTERNAL FUNCTION DIMCHECKLOW(X())
FOR N=MINLEVEL TO MAXLEVEL
IF X(N)<>0 THEN EXIT FOR
NEXT N
LET DIMCHECKLOW=N
END FUNCTION
 
　　├続き しばっち 2008/03/28 20:30:59   ツリーへ
Re: 多項式逆数 返事を書く ノートメニュー
しばっち <dihjvcfsyu> 2008/03/28 20:30:59
続き


EXTERNAL SUB SINE(X())
!'SIN(X)
CALL CLR(X)
LET X(1)=1
LET T=1
FOR I=3 TO MAXLEVEL STEP 2
LET T=-T/(I-1)/I
LET X(I)=T
NEXT I
END SUB

EXTERNAL SUB COSINE(X())
!'COS(X)
CALL CLR(X)
LET X(0)=1
LET T=1
FOR I=2 TO MAXLEVEL STEP 2
LET T=-T/(I-1)/I
LET X(I)=T
NEXT I
END SUB
 
　　└多項式平方根 しばっち 2008/03/28 20:31:48   ツリーへ
Re: 多項式逆数 返事を書く ノートメニュー
しばっち <dihjvcfsyu> 2008/03/28 20:31:48
多項式平方根


F(X)=A+B*X+C*X^2+D*X^3+...
G(X)=A'+B'*X+C'*X^2+D'*X^3...
F(X)=G(X)^2 G(X)を求める

恒等式 (A+B*X+C*X^2+D*X^3+E*X^4+...)=(A'+B'*X+C'*X^2+D'*X^3+E'*X^4+...)^2

A=A'^2
B=2A'B'
C=2A'C'+B'^2
D=2A'D'+2B'C'
:

A'=SQR(A)
B'=B/2A'
C'=(C-B'^2)/2A'
D'=(D-2B'C')/2A'
:

OPTION BASE 0
PUBLIC NUMERIC MAXLEVEL
!'1/(1-SIN(X)^2)=1/COS(X)^2=1+SIN(X)^2+SIN(X)^4+SIN(X)^6+...
LET MAXLEVEL=40
DIM C(MAXLEVEL),X(MAXLEVEL)
!'1+X^2+X^4+X^6+X^8+...=1/(1-X^2)
FOR I=0 TO MAXLEVEL/2 STEP 2
LET C(I)=1
NEXT I
CALL SINE(X)
CALL HORNER(C,X)
PRINT "1/(1-SIN(X)^2)"
CALL DISPLAY(C)
CALL SQRT(C)
PRINT "1/COS(X)"
CALL DISPLAY(C)
END

EXTERNAL SUB SQRT(C())
OPTION BASE 0
DIM XX(MAXLEVEL),X(MAXLEVEL),CC(MAXLEVEL)
FOR M=0 TO MAXLEVEL
IF C(M)<>0 THEN EXIT FOR
NEXT M
IF MOD(M,2)=1 OR C(M)<0 THEN
PRINT "計算不可"
EXIT SUB
END IF
FOR MM=MAXLEVEL TO M STEP -1
IF C(MM)<>0 THEN EXIT FOR
NEXT MM
IF MOD(MM,2)=1 THEN
PRINT "計算不可"
EXIT SUB
END IF
FOR I=0 TO MAXLEVEL-M
LET CC(I)=C(I+M)
LET C(I)=0
NEXT I
LET X(0)=SQR(CC(0))
FOR K=1 TO MM/2
CALL CLR(XX)
FOR I=0 TO MAXLEVEL
FOR J=0 TO MAXLEVEL-I
LET XX(I+J)=XX(I+J)+X(I)*X(J)
NEXT J
NEXT I
LET X(K)=(CC(K)-XX(K))/(2*X(0))
NEXT K
FOR I=0 TO MAXLEVEL-M
LET C(I+M/2)=X(I)
NEXT I
END SUB

EXTERNAL SUB COPY(X(),Y())
FOR I=0 TO MAXLEVEL
LET X(I)=Y(I)
NEXT I
END SUB

EXTERNAL SUB CLR(X())
FOR I=0 TO MAXLEVEL
LET X(I)=0
NEXT I
END SUB

EXTERNAL SUB DISPLAY(A())
LET N=DIMCHECK(A)
IF N>1 THEN
IF A(N)<0 THEN PRINT "-";
IF ABS(A(N))<>1 THEN
PRINT STR$(ABS(A(N)));"*X^";STR$(N);
ELSE
PRINT "X^";STR$(N);
END IF
END IF
FOR I=N-1 TO 2 STEP -1
IF A(I)<>0 THEN
IF A(I)<0 THEN PRINT "-"; ELSE PRINT "+";
IF ABS(A(I))<>1 THEN
PRINT STR$(ABS(A(I)));"*X^";STR$(I);
ELSEIF ABS(A(I))=1 THEN
PRINT "X^";STR$(I);
END IF
END IF
NEXT I
IF A(1)<>0 THEN
IF N>1 THEN
IF A(1)<0 THEN PRINT "-"; ELSE PRINT "+";
END IF
IF ABS(A(1))<>1 THEN
PRINT STR$(ABS(A(1)));"*X";
ELSEIF ABS(A(1))=1 THEN
PRINT "X";
END IF
END IF
IF A(0)<>0 THEN
IF A(0)<0 THEN PRINT "-"; ELSE PRINT "+";
PRINT STR$(ABS(A(0)));
END IF
PRINT
END SUB

EXTERNAL FUNCTION DIMCHECK(X())
FOR N=MAXLEVEL TO 0 STEP -1
IF X(N)<>0 THEN EXIT FOR
NEXT N
LET DIMCHECK=N
END FUNCTION

EXTERNAL SUB HORNER(F(),X())
!'F(X(XX))
OPTION BASE 0
DIM Y(MAXLEVEL)
LET N=DIMCHECK(F)
LET Y(0)=F(N)
FOR I=N-1 TO 0 STEP -1
CALL MUL(Y,X)
LET Y(0)=Y(0)+F(I)
NEXT I
CALL COPY(F,Y)
END SUB
 
　　　└続き しばっち 2008/03/28 20:32:38   ツリーへ
Re: 多項式平方根 返事を書く ノートメニュー
しばっち <dihjvcfsyu> 2008/03/28 20:32:38
続き


EXTERNAL SUB MUL(Y(),X())
OPTION BASE 0
DIM C(MAXLEVEL)
FOR J=0 TO MAXLEVEL
FOR I=0 TO MAXLEVEL-J
LET C(I+J)=C(I+J)+Y(I)*X(J)
NEXT I
NEXT J
CALL COPY(Y,C)
END SUB

EXTERNAL SUB SINE(X())
!'SIN(X)
CALL CLR(X)
LET X(1)=1
LET T=1
FOR I=3 TO MAXLEVEL STEP 2
LET T=-T/(I-1)/I
LET X(I)=T
NEXT I
END SUB

インデックスへ EXIT
新規発言を反映させるにはブラウザの更新ボタンを押してください。