式の評価、もう1つの拡張

 投稿者:山中和義  投稿日:2009年 2月20日(金)15時39分28秒
  「数の組」に対する「式の評価(計算)」ができるように次の拡張を行った。
・構文解析と演算との部分を分離
・FUNCTION文では返り値を1つのみのため、SUB文で記述

「数の組」とは
  関数、変数、定数の通常の数式の場合、「値」(1つ)
  有理数の場合、「分子と分母」(2つ)
  複素数の場合、「実数部と虚数部」(2つ)
  ベクトル、行列の場合、「要素列」(n個)
  1変数多項式の係数が整数の場合、「係数列」(n個)
 とする。


記述例 1変数多項式の係数が整数の場合

!式(中置記法)の評価 - 1変数多項式の係数が整数の場合 ※UBASIC相当

DECLARE EXTERNAL SUB expr.eval

!LET s$="(x-2)*(x-3)*(x-4)*(x-5)"
LET s$="(x^2-2*x+3)^2"
!LET s$="(x^5-5*x^3+5*x^2-1)/(x^2+3*x+1)" !商 x^3-3*x^2+3*x-1、余り 0
!LET s$="mod(2*x^3-13*x^2-26*x-15,x^2-2*x+3)" !商 2*x-9、余り -50*x+12
!LET s$="gcd(3*x^2+5*x-2,3*x^2-7*x+2)" !12*x-4 { 3*x-1 }
!LET s$="lcm(3*x^2+5*x-2,3*x^2-7*x+2)" !3/4*x^3-1/4*x^2-3*x+1 { (x+2)*(x-2)*(3*x-1) }
!LET s$="lcm(x^3-16*x,x^3-8*x^2+16*x)"


DIM a(0 TO 8) !係数 a(n)*x^n+a(n-1)*x^(n-1)+ … +a(1)*x+a(0)
CALL eval(s$, a,rc) !式

MAT PRINT a; !x^0,x^1,x^2, … の順 debug

IF rc=0 THEN CALL poly_disp(a) !結果を表示する

END



MODULE expr

PUBLIC NUMERIC p,ErrNo !共通変数

!●解析部分

!下位の共通ルーチン
EXTERNAL FUNCTION token$(s$) !1文字読み込む
   CALL EatSpace(s$)
   IF p<=LEN(s$) THEN LET token$=s$(p:p) ELSE LET token$=""
END FUNCTION

EXTERNAL SUB EatSpace(s$) !空白を読み飛ばす
   DO WHILE s$(p:p)=" " AND p<=LEN(s$)
      LET p=p+1
   LOOP
END SUB

EXTERNAL SUB CheckToken(s$,L$) !文字を確認する
   CALL EatSpace(s$)
   IF UCASE$(s$(p:p+LEN(L$)-1))<>L$ THEN CALL Error(L$&"がありません。")
   LET p=p+LEN(L$) !eat it
END SUB

EXTERNAL SUB Error(x$) !メッセージを表示する
   PRINT
   PRINT x$; p

   LET errNo=1
END SUB


!上位ルーチン
PUBLIC SUB eval
EXTERNAL SUB eval(s$, v(),rc) !式の評価
   LET errNo=0 !エラーコード
   LET p=1 !文字列へのポインタ
   CALL expression(s$, v) !計算する
   LET rc=errNo
END SUB

EXTERNAL SUB expression(s$, v()) !式
   DIM w(0 TO UBOUND(v))

   LET t$=token$(s$)
   IF t$="-" THEN !符号なら
      LET p=p+1 !eat it

      CALL term(s$, v)
      IF errNo<>0 THEN EXIT SUB

      CALL op_neg(v, v) !v=-v
   ELSE
      IF t$="+" THEN LET p=p+1 !eat it
      CALL term(s$, v)
   END IF
   IF errNo<>0 THEN EXIT SUB


   LET t$=token$(s$)
   DO WHILE t$="+" OR t$="-" !加算、減算なら
      LET p=p+1 !eat it

      CALL term(s$, w)
      IF errNo<>0 THEN EXIT SUB

      IF t$="+" THEN !計算する
         CALL op_add(v,w, v) !v=v+w
      ELSE
         CALL op_sub(v,w, v) !v=v-w
      END IF
      IF errNo<>0 THEN EXIT SUB

      LET t$=token$(s$) !次へ
   LOOP
END SUB

EXTERNAL SUB term(s$, v()) !項
   DIM w(0 TO UBOUND(v))

   CALL factor(s$,v)
   IF errNo<>0 THEN EXIT SUB

   LET t$=token$(s$)
   DO WHILE t$="*" OR t$="/" !乗算、除算なら
      LET p=p+1 !eat it

      CALL factor(s$,w)
      IF errNo<>0 THEN EXIT SUB

      IF t$="*" THEN !計算する
         CALL op_mul(v,w, v) !v=v*w
      ELSE
         CALL op_div(v,w, v) !v=v/w
      END IF
      IF errNo<>0 THEN EXIT SUB

      LET t$=token$(s$) !次へ
   LOOP
END SUB

EXTERNAL SUB factor(s$, v()) !因子
   DIM w(0 TO UBOUND(v))

   LET t$=token$(s$)
   IF t$="(" THEN !括弧なら
      LET p=p+1 !eat it

      CALL expression(s$,w) !式
      IF errNo<>0 THEN EXIT SUB

      MAT v=w

      CALL CheckToken(s$,")") !閉じ括弧か確認する
   ELSE
      CALL num(s$,v)
   END IF
   IF errNo<>0 THEN EXIT SUB


   LET t$=token$(s$)
   DO WHILE t$="^" !べき乗なら
      LET p=p+1 !eat it

      LET t$=token$(s$)
      IF t$="(" THEN !括弧なら
         LET p=p+1 !eat it

         CALL expression(s$,w) !式
         IF errNo<>0 THEN EXIT SUB

         CALL CheckToken(s$,")") !閉じ括弧か確認する
      ELSE
         CALL num(s$,w)
      END IF
      IF errNo<>0 THEN EXIT SUB


      CALL op_pow(v,w, v) !計算する v=v^w
      IF errNo<>0 THEN EXIT SUB


      LET t$=token$(s$) !次へ
   LOOP
END SUB

EXTERNAL SUB num(s$,v()) !数
   DIM w(0 TO UBOUND(v)),x(0 TO UBOUND(v))

   LET c=fnc(s$)
   IF c>0 THEN !関数なら
      LET t$=token$(s$)
      IF t$="(" THEN !括弧なら
         LET p=p+1 !eat it

         CALL expression(s$,w) !引数1
         IF errNo<>0 THEN EXIT SUB

         MAT v=w

         CALL CheckToken(s$,",") !カンマか確認する

         CALL expression(s$,w) !引数2
         IF errNo<>0 THEN EXIT SUB

         IF c=1 THEN !modpow(a,n,b)形式
            CALL CheckToken(s$,",") !カンマか確認する

            CALL expression(s$,w) !引数3
            IF errNo<>0 THEN EXIT SUB

            MAT x=w

            CALL set_fnc3(c,v,w,x, v) !v=fnc(v,w,x)
         ELSE
            CALL set_fnc2(c,v,w, v) !v=fnc(v,w)
         END IF
         IF errNo<>0 THEN EXIT SUB

         CALL CheckToken(s$,")") !閉じ括弧か確認する
      ELSE
         CALL Error("不正な文字です。")
      END IF

   ELSE
      LET c=var(s$)
      IF c>0 THEN !変数なら
         CALL set_var(c, v)
      ELSE
         LET c=number(s$)
         IF c>=0 THEN !定数(数値)なら
            CALL set_number(c, v)
         ELSE
            CALL Error("不正な文字です。")
         END IF
      END IF

   END IF
END SUB

EXTERNAL FUNCTION fnc(s$) !関数
   DATA "MODPOW","MODINV","MOD","GCD","LCM" !※文字長が大きい順
   LET k=0
   DO
      LET k=k+1
      READ IF MISSING THEN EXIT DO: d$
      IF UCASE$(s$(p:p+LEN(d$)-1))=UCASE$(d$) THEN !一致したら
         LET p=p+LEN(d$)
         LET fnc=k
         EXIT FUNCTION
      END IF
   LOOP
   LET fnc=-1
END FUNCTION

EXTERNAL FUNCTION var(s$) !変数
   LET t$=UCASE$(token$(s$)) !大文字へ
   IF "A"<=t$ AND t$<="Z" THEN
   !!!IF t$="X" THEN
      LET p=p+1 !eat it
      LET var=ORD(t$)-ORD("@") !オフセット @=0,A=1,B=2,…,Z=26
   ELSE
      LET var=-1
   END IF
END FUNCTION

EXTERNAL FUNCTION number(s$) !数値(0、正の整数)
   LET i0=p !先頭位置を記録する

   LET t$=token$(s$)
   DO WHILE t$>="0" AND t$<="9"
      LET p=p+1
      LET t$=token$(s$)
   LOOP

   LET number=-1
   IF p>i0 THEN LET number=VAL(s$(i0:p-1)) !数字列の範囲を切り取る
END FUNCTION


つづく
 

Re: 式の評価、もう1つの拡張

 投稿者:山中和義  投稿日:2009年 2月20日(金)15時42分40秒
  > No.286[元記事へ]

つづき

!●演算部分 ※「数の組」に応じて演算を定義する

EXTERNAL SUB op_neg(v1(), v()) !符号(負)
   MAT v=(-1)*v1
END SUB

EXTERNAL SUB op_add(v1(),v2(), v()) !加算
   MAT v=v1+v2
END SUB

EXTERNAL SUB op_sub(v1(),v2(), v()) !減算
   MAT v=v1-v2
END SUB

EXTERNAL SUB op_mul(v1(),v2(), v()) !乗算
   LET N=UBOUND(v)
   DIM w(0 TO N*2) !桁数は2倍になる

   MAT w=ZER
   FOR i=0 TO N !係数
      FOR j=0 TO N
         LET w(i+j)=w(i+j)+v1(i)*v2(j) !畳み込み
      NEXT j
   NEXT i
   IF poly_degree(w)>N THEN CALL Error("オーバーフロー")

   FOR i=0 TO N !下n桁をコピーする
      LET v(i)=w(i)
   NEXT i
END SUB

EXTERNAL SUB op_div(v1(),v2(), v()) !除算
   DIM Q(0 TO UBOUND(v)),R(0 TO UBOUND(v))

   IF poly_degree(v2)=0 AND v2(0)=0 THEN !定数の0なら
      CALL Error("0では割れません。")
      EXIT SUB
   END IF
   CALL poly_div(v1,v2, Q,R)
   MAT v=Q !商 ※v=INT(v1/v2)に相当
END SUB

EXTERNAL SUB op_pow(v1(),v2(), v()) !べき算
   DIM x(0 TO UBOUND(v)),T(0 TO UBOUND(v))

   MAT x=v1 !x=v1

   MAT T=ZER !定数1
   LET T(0)=1

   IF poly_degree(v1)=0 THEN !定数項のみ(定数)
      IF poly_degree(v2)>0 THEN CALL Error("べき数は整数のみ")
      IF errNo<>0 THEN EXIT SUB
      IF v1(0)=0 AND v2(0)<0 THEN CALL Error("0の負べき乗")
      IF errNo<>0 THEN EXIT SUB

      LET v(0)=v1(0)^v2(0)
   ELSE
      IF poly_degree(v2)>0 OR v2(0)<0 THEN CALL Error("べき数は非負整数のみ")
      IF errNo<>0 THEN EXIT SUB

      LET m=v2(0)
      DO UNTIL m=0
         IF MOD(m,2)=1 THEN CALL op_mul(T,x, T) !ビットが1なら、T=T*x
         IF errNo<>0 THEN EXIT SUB
         CALL op_mul(x,x, x) !x=x^2
         IF errNo<>0 THEN EXIT SUB

         LET m=INT(m/2) !2進数にする
      LOOP
      MAT v=T
   END IF
END SUB


EXTERNAL SUB set_fnc(c,v1(), v()) !関数(値)を設定する ※引数1個
!該当関数なし
END SUB

EXTERNAL SUB set_fnc2(c,v1(),v2(), v()) !関数(値)を設定する ※引数2個
   LET N=UBOUND(v)
   DIM Q(0 TO N),R(0 TO N),A(0 TO N),B(0 TO N)

   SELECT CASE c !各関数に応じて
   CASE 2 !modinv
      PRINT "modinv関数は未サポート!"
   CASE 3 !剰余
      CALL poly_div(v1,v2, Q,R)
      MAT v=R !余り
   CASE 4,5 !最大公約数、最小公倍数
      MAT A=v1
      MAT B=v2
      DO UNTIL poly_degree(B)=0 AND B(0)=0 !b=0
         CALL poly_div(A,B, Q,R) !R=MOD(a,b)
         MAT A=B !a=b
         MAT B=R !b=R
      LOOP

      IF c=5 THEN !LCM=v1*v2/GCD(v1,v2)
         CALL op_div(v1,A, Q)
         IF errNo<>0 THEN EXIT SUB
         CALL op_mul(Q,v2, v)
      ELSE !GCD
         MAT v=A
      END IF
   CASE ELSE
   END SELECT
END SUB

EXTERNAL SUB set_fnc3(c,v1(),v2(),v3(), v()) !関数(値)を設定する ※引数3個
   PRINT "modpow関数は未サポート!"
END SUB

EXTERNAL SUB set_var(c, v()) !変数(値)
   MAT v=ZER
   LET v(1)=1 !x^1の係数
END SUB

EXTERNAL SUB set_number(c, v()) !定数(数値)
   MAT v=ZER
   LET v(0)=c !x^0の係数
END SUB

END MODULE



!補助ルーチン

!演算関連

EXTERNAL SUB poly_div(A(),B(), Q(),R()) !除算 ※被除数=商*除数+余り
DIM w(0 TO UBOUND(A))

LET aa=poly_degree(A)
LET bb=poly_degree(B)

MAT Q=ZER !商、その次数
LET qq=MAX(aa-bb,0)

MAT R=A !余り、その次数
LET rr=aa

DO WHILE rr>=bb !被除数の次数が除数のより大きいなら
   IF R(rr)<>0 THEN !係数が0以外なら
      LET k=R(rr)/B(bb) !商の係数
      LET Q(rr-bb)=k !商

      MAT w=ZER !余り
      FOR i=bb TO 0 STEP -1 !R=A-k*B ※筆算参照
         LET w(rr-bb+i)=k*B(i)
      NEXT i
      MAT R=R-w
   END IF
   LET rr=rr-1 !次の次数へ
LOOP
END SUB

EXTERNAL FUNCTION poly_degree(v()) !次数を得る
FOR i=UBOUND(v) TO 1 STEP -1
   IF v(i)<>0 THEN EXIT FOR !係数が0でない最初の位置
NEXT i
LET poly_degree=i
END FUNCTION


!表示関連

EXTERNAL SUB poly_disp(A()) !多項式を表示する a(X)=ΣAkX^k=AnX^n+An-1X^n-1+…+A1X+A0
LET aa=poly_degree(A) !最初の項
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 (w=0 AND aa=0) THEN CALL mono_disp(w,i)
NEXT i
END SUB

EXTERNAL SUB mono_disp(ak,k) !単項式を表示する Ak*X^k
IF k<>0 THEN !x^nで
   IF ak=1 THEN !係数が1なら
   ELSEIF ak=-1 THEN !係数が-1なら
      PRINT "-";
   ELSE
      PRINT STR$(ak);"*";
   END IF
END IF
IF k=0 THEN !次数が0なら
   PRINT STR$(ak);
ELSEIF k=1 THEN !次数が1なら
   PRINT "X";
ELSE
   PRINT "X^";STR$(k);
END IF
END SUB


以上
 

戻る