新しく発言する EXIT インデックスへ
ダイナミックプログラミング

  ダイナミックプログラミング 島村1243 2007/07/08 10:57:36  (修正1回)
  調べてみました。たぶん、こんな感じだと思... 山中和義 2007/07/10 11:49:26  (修正1回)
  │├つづき 山中和義 2007/07/10 11:50:30 
  │├山中様、有難う御座います。 島村1243 2007/07/10 23:37:16 
  │└Linuxでの実施報告1 島村1243 2007/07/11 22:29:19  (修正2回)
  !再帰型を試みたけれども、速度その他に問題... SECOND 2007/07/12 14:19:50 
   └!Variation-2全角変数なし。 SECOND 2007/07/12 20:29:41 
    └鮮やかなプログラミングですね! 島村1243 2007/07/12 22:22:23 
     ├私にも、わからない、どなたか説明できる方... SECOND 2007/07/13 03:41:49 
     └!これで、試して下さい。 SECOND 2007/07/13 17:48:02 
      └演算結果の比較 島村1243 2007/07/13 22:18:52  (修正2回)
       └!SIN(An)/2/k,SIN(An)/2は、山中氏の結果に... SECOND 2007/07/13 23:59:52  (修正3回)
        └たしかに 島村1243 2007/07/14 08:30:14 
         └十進BASICのバグの疑いがあります。というの... SECOND 2007/07/14 10:42:40 
          └FullBASICの定義関数名は変数ではなく,従っ... 白石 和夫 2007/07/14 11:27:18  (修正1回)
           ├できれば、この一件を、トレースして頂けれ... SECOND 2007/07/14 12:46:05 
           ├LOCAL文の不具合かも? 山中和義 2007/07/14 13:28:20 
           │└ありがとうございました。 SECOND 2007/07/14 14:30:48 
           └!問題が表面化するテストプログラムです。 SECOND 2007/07/14 14:24:49 
            ├FORv=始値TO限界STEP増分 白石 和夫 2007/07/14 16:35:52 
            │├簡単な例 白石 和夫 2007/07/14 16:43:42 
            │└こちらはown1を局所化したものです。 白石 和夫 2007/07/14 16:47:44 
            └以下のようにすると,同じ動作になります。 白石 和夫 2007/07/14 17:35:42 
             ├ありがとうございました。 SECOND 2007/07/14 17:46:32 
             └再帰式中ではForNextが使えない? 島村1243 2007/07/14 22:14:08  (修正2回)
              ├「内部関数定義,内部副プログラムは原則と... 白石 和夫 2007/07/15 07:09:57  (修正1回)
              │└内部関数定義や内部副プログラムは,外部手... 白石 和夫 2007/07/15 09:17:49 
              │ └大変良く分かりました。 島村1243 2007/07/15 22:45:21 
              └お疲れさまでした。代数的に、n-1番目の結果... SECOND 2007/07/15 11:40:44 
               └大変勉強になりました。 島村1243 2007/07/15 23:09:14  (修正3回)

  ダイナミックプログラミング 島村1243 2007/07/08 10:57:36  (修正1回) ツリーへ

ダイナミックプログラミング ノートメニュー
島村1243 <bjllmpcujp> 2007/07/08 10:57:36 ** この記事は1回修正されてます
パソコンが無い時代に、1次元配分過程の問題を解く手法としてのダイナミックプログラミング(DP)を本で勉強し
ました。問題の1例を挙げると次のとおりです。

中心点をoとする半径1、角度βの扇型を5個の角度θ1、θ2、θ3、θ4、θ5に分けたいがどの様に分ければ良いか?
<分ける為の条件>
(1)1個の角θ1を挟む2辺と扇型の円周との交点をそれぞれa1、b1とするとき、0-a1-b1を頂点とする三角形の面積を
S1とする。同様にして他の角度θ2、θ3、θ4、θ5についても面積S2、S3、S4、S5を定める。
(2)以上の場合に(S1+S2+S3+S4+S5)が最大になる事。

これをDPの漸化式で表せば
Fn(β)=Max{Sn(θn)+Fn-1(β-θn)} 、Sn(θn)=1/2*sin(θn)
F1(β)=Max{S1(θ1)} =1/2*sin(β)
となり、n=1からn=5迄計算して解きます。
解析的に解いた答えはθ1=θ2=θ3=θ4=θ5、つまり5等分すれば最大面積が得られる事になります。

本には、格子点法(グリッド)を用いる数値計算方法が紹介されていましたが、パソコンが無い時代での事なので、頭
の中でイメージしただけで実行はしておりませんが、一度プログラムを組めば関数「FnやSn」定義を書き替えるだ
けで他への応用が出来そうで、面白い問題だと思っています。
漸化式とMiniMaxの計算なのでBASICで数値計算する事が出来る訳ですが、この様な問題のプログラムを作成した方
おられましたら掲示頂けると有難い(SECONDさんの御言葉に甘えて恐縮)です。

  調べてみました。たぶん、こんな感じだと思... 山中和義 2007/07/10 11:49:26  (修正1回) ツリーへ

Re: ダイナミックプログラミング ノートメニュー
山中和義 <drdlxujciw> 2007/07/10 11:49:26 ** この記事は1回修正されてます
調べてみました。たぶん、こんな感じだと思います。

各角度(0〜β°)において、n分割した(n分割は、1分割とN−1分割)
・三角形面積の和の最大(Sn())
・そのときの1分割側の角度(THn())
として、各過程の計算結果を活かして、次々と求めています。


2以降はプログラムの構造が同じなので、2次元の配列を使って、FOR文などで簡略化できると思います。
 DIM S3(β),th3(β) ⇒ DIM S(3,β),th(3,β)


再帰処理や総当りが頭をかすめて、まとめが遅くなりました。
変数名に2バイト文字(β)を使っているため、Linuxでは動かないかもしれません。




!目的関数 S(x)=S1(θ1)+S2(θ2)+S3(θ3)+S4(θ4)+S5(θ5)
!     Sn(θn)=SIN(θn)/2(n=1〜5)
!制約条件 θ1+θ2+θ3+θ4+θ5=β、θn≧0(n=1〜5)

!最適性の原理に従って、再帰方程式に変換する
!φi(ai)=MAX[Si(θi) + φi-1(ai-θi)] i=2〜5、0≦θi≦β
!φ1(a1)=MAX[S1(θ1)]=S1(a1) 0≦θ1≦β

OPTION BASE 0


LET β=120 !角度 ※1度刻み

!再帰方程式
!φ1(a1)=MAX[S1(θ1)]=S1(a1) 0≦θ1≦βより
DIM S1(β),th1(β) !面積、角度
SUB P1(a1, S1(),th1())
FOR i=0 TO a1 !角度0〜a1で、その角度を1分割する三角形の面積を求める
LET S1(i)=SIN(RAD(i))/2 !半径1
LET th1(i)=i !角度を記録する
NEXT i
END SUB
CALL P1(β, S1,th1) !φ1(β)=S1(β)=SIN(β)/2
FOR i=0 TO β
PRINT USING "###°: ###° ##.########":i,th1(i),S1(i)
NEXT i

LET v=β
LET v0=th1(v) !1分割
PRINT "θ1=";v0

PRINT β, SIN(RAD(β))/2 !検算
PRINT



!φ2(a2)=MAX[S2(θ2) + φ1(a2-θ2)] 0≦θ2≦βより
DIM S2(β),th2(β)
SUB P2(a2, S2(),th2())
FOR i=0 TO a2 !角度0〜a2で、その角度を
FOR j=0 TO i !2分割する三角形の面積の和の中で、最大になるものを求める
LET v=S1(j)+ S1(i-j) !1分割+1分割
IF v>S2(i) THEN
LET S2(i)=v
LET th2(i)=j !角度を記録する
END IF
NEXT j
NEXT i
END SUB
CALL P2(β, S2,th2) !φ2(β)=MAX[S2(θ2) + φ1(β-θ2)] 0≦θ2≦β
FOR i=0 TO β
PRINT USING "###°: ###° ##.########":i,th2(i),S2(i)
NEXT i

LET v=β
LET v0=th2(v) !1分割
PRINT "θ2=";v0

LET v=v-v0
LET v0=th1(v) !残り1分割
PRINT "θ1=";v0

PRINT β/2, 2*SIN(RAD(β/2))/2 !検算
PRINT

  │├つづき 山中和義 2007/07/10 11:50:30  ツリーへ

Re: 調べてみました。たぶん、こんな感じだと思... ノートメニュー
山中和義 <drdlxujciw> 2007/07/10 11:50:30
つづき



!φ3(a3)=MAX[S3(θ3) + φ2(a3-θ3)] 0≦θ3≦βより
DIM S3(β),th3(β)
SUB P3(a3, S3(),th3())
FOR i=0 TO a3 !角度0〜a3で、その角度を
FOR j=0 TO i !3分割する三角形の面積の和の中で、最大になるものを求める
LET v=S1(j)+ S2(i-j) !1分割+2分割
IF v>S3(i) THEN
LET S3(i)=v
LET th3(i)=j !角度を記録する
END IF
NEXT j
NEXT i
END SUB
CALL P3(β, S3,th3) !φ3(β)=MAX[S3(θ3) + φ2(β-θ3)] 0≦θ3≦β
FOR i=0 TO β
PRINT USING "###°: ###° ##.########":i,th3(i),S3(i)
NEXT i

LET v=β
LET v0=th3(v) !1分割
PRINT "θ3=";v0

LET v=v-v0
LET v0=th2(v) !残り2分割
PRINT "θ2=";v0

LET v=v-v0
LET v0=th1(v) !たどっていく
PRINT "θ1=";v0

PRINT β/3, 3*SIN(RAD(β/3))/2 !検算
PRINT



!φ4(a4)=MAX[S4(θ4) + φ3(a4-θ4)] 0≦θ4≦βより
DIM S4(β),th4(β)
SUB P4(a4, S4(),th4())
FOR i=0 TO a4 !角度0〜a4で、その角度を
FOR j=0 TO i !4分割する三角形の面積の和の中で、最大になるものを求める
LET v=S1(j)+ S3(i-j) !1分割+3分割
IF v>S4(i) THEN
LET S4(i)=v
LET th4(i)=j !角度を記録する
END IF
NEXT j
NEXT i
END SUB
CALL P4(β, S4,th4) !φ4(β)=MAX[S4(θ4) + φ3(β-θ4)] 0≦θ4≦β
FOR i=0 TO β
PRINT USING "###°: ###° ##.########":i,th4(i),S4(i)
NEXT i

LET v=β
LET v0=th4(v) !1分割
PRINT "θ4=";v0

LET v=v-v0
LET v0=th3(v) !残り3分割
PRINT "θ3=";v0

LET v=v-v0
LET v0=th2(v) !たどっていく
PRINT "θ2=";v0

LET v=v-v0
LET v0=th1(v) !たどっていく
PRINT "θ1=";v0

PRINT β/4, 4*SIN(RAD(β/4))/2 !検算
PRINT



!φ5(a5)=MAX[S5(θ5) + φ4(a5-θ5)] 0≦θ5≦βより
DIM S5(β),th5(β)
SUB P5(a5, S5(),th5())
FOR i=0 TO a5 !角度0〜a5で、その角度を
FOR j=0 TO i !5分割する三角形の面積の和の中で、最大になるものを求める
LET v=S1(j)+ S4(i-j) !1分割+4分割
IF v>S5(i) THEN
LET S5(i)=v
LET th5(i)=j !角度を記録する
END IF
NEXT j
NEXT i
END SUB
CALL P5(β, S5,th5) !φ5(β)=MAX[S5(θ5) + φ4(β-θ5)] 0≦θ5≦β
FOR i=0 TO β
PRINT USING "###°: ###° ##.########":i,th5(i),S5(i)
NEXT i

LET v=β
LET v0=th5(v) !1分割
PRINT "θ5=";v0

LET v=v-v0
LET v0=th4(v) !残り4分割
PRINT "θ4=";v0

LET v=v-v0
LET v0=th3(v) !たどっていく
PRINT "θ3=";v0

LET v=v-v0
LET v0=th2(v) !たどっていく
PRINT "θ2=";v0

LET v=v-v0
LET v0=th1(v) !たどっていく
PRINT "θ1=";v0

PRINT β/5, 5*SIN(RAD(β/5))/2 !検算
PRINT


END

  │├山中様、有難う御座います。 島村1243 2007/07/10 23:37:16  ツリーへ

Re: 調べてみました。たぶん、こんな感じだと思... ノートメニュー
島村1243 <bjllmpcujp> 2007/07/10 23:37:16
山中様、有難う御座います。

短時間でプログラミング出来る山中様の力を感嘆しています。
プログラムの全容をすぐに理解する力が無いので、明日から楽しみに勉強させて頂きます。
理解できたら改めて御礼の報告を書かせて頂きます。

  │└Linuxでの実施報告1 島村1243 2007/07/11 22:29:19  (修正2回) ツリーへ

Re: 調べてみました。たぶん、こんな感じだと思... ノートメニュー
島村1243 <bjllmpcujp> 2007/07/11 22:29:19 ** この記事は2回修正されてます
Linuxでの実施報告1
プログラムの内容は理解できました。
Debian3.1(文字コードはeuc_JPを使用)とFedoraCore5(文字コード
はUTF-8を使用)上の十進BASIC-6.3.9版で、このプログラム(山中様の生コード)を走らせて見ました。中山様が心配されていた2バイト文字βがあるま
までrunしても正常に動作して解を得られました。Vine4.0でも問題ないと思います。以上ご報告致します。

今回の例題では、各段階での利得関数S1〜S5が同じ形の関数で与えられていましたが、異なる式で与えられた場合
にも半汎用的に使える様にするため、少し手を加えさせて頂こうと思っています。上手くいきましたら掲載させて
頂きます。

  !再帰型を試みたけれども、速度その他に問題... SECOND 2007/07/12 14:19:50  ツリーへ

Re: ダイナミックプログラミング ノートメニュー
SECOND <cszcthjjdj> 2007/07/12 14:19:50
!再帰型を試みたけれども、速度その他に問題を残した。
!どなたか、解決を。

!個々のθ(1~n)は、90度を越えられない。

OPTION ANGLE DEGREES
DEF Sn(θn)= 1*SIN(θn)/2

LET β=120
LET n=6

DIM θ(n)

FOR k=2 TO n
PRINT "Fn(";β;",";k;")=";
PRINT Fn(β,k);TAB(33);
FOR i=1 TO k
PRINT θ(i);
NEXT i
PRINT
NEXT k


FUNCTION Fn(β,k)
local i,fm
LET Fn=0
IF 0<k THEN
LET fm=0
FOR i=19 TO β !----------------- 本来は、1 TO β、速度に難点がある為のズル。
LET f=Sn(i)+ Fn(β-i,k-1) !--- 再帰 FUNCTION
IF f<fm THEN EXIT FOR !------- 極大値が、1つしかない前提のズル。
LET fm=f
NEXT i
LET Fn=fm
LET θ(k)=i-1
END IF
END FUNCTION

END

   └!Variation-2全角変数なし。 SECOND 2007/07/12 20:29:41  ツリーへ

Re: !再帰型を試みたけれども、速度その他に問題... ノートメニュー
SECOND <cszcthjjdj> 2007/07/12 20:29:41
!Variation-2 全角変数なし。

!個々のA(1~n)は、90度を越えられない。

OPTION ANGLE DEGREES
DEF Sn(An)= 1*SIN(An)/2

LET B=100 !120 のときは、FOR i=15 TO B を、19 TO B にすると速くなる。

DIM A(6)

PRINT "Fn(";B;",2)="; Fn(B, 2);TAB(32);A(1);A(2)
PRINT "Fn(";B;",3)="; Fn(B, 3);TAB(32);A(1);A(2);A(3)
PRINT "Fn(";B;",4)="; Fn(B, 4);TAB(32);A(1);A(2);A(3);A(4)
PRINT "Fn(";B;",5)="; Fn(B, 5);TAB(32);A(1);A(2);A(3);A(4);A(5)
PRINT "Fn(";B;",6)="; Fn(B, 6);TAB(32);A(1);A(2);A(3);A(4);A(5);A(6)


FUNCTION Fn(B,k)
local i,fm
LET Fn=0
IF 0<k THEN
LET fm=0
FOR i=15 TO B !----------------- 本来は、1 TO B、速度に難点がある為のゲタ。
LET f=Sn(i)+ Fn(B-i,k-1) !--- 再帰 FUNCTION
IF f<fm THEN EXIT FOR !------- 極大値が、1つしかない前提。
LET fm=f
NEXT i
LET Fn=fm
LET A(k)=i-1
END IF
END FUNCTION

END

    └鮮やかなプログラミングですね! 島村1243 2007/07/12 22:22:23  ツリーへ

Re: !Variation-2全角変数なし。 ノートメニュー
島村1243 <bjllmpcujp> 2007/07/12 22:22:23
鮮やかなプログラミングですね!
でも、ずーっと再帰関係式を見ているのですが、再帰式の中のFn(B-i,K-1)が

Fn(B,K)=Max{Sn(i)+Fn(B-i,K-1)} ,0<i<B
=Max{Sn(i)+Max{Sn(j)+Fn(B-i-j,K-2)}} ,0<j<B-i
= ・・・
と言う関係を満足して再帰を繰り返すしている(のだろう)と言う事が読めずに苦慮しています。
苦慮している事は下記のとおりです。

今回の例題は、全ての分割角度が等しいと言うのが正解であるがために、B=120、K=5の場合に

PRINT "Fn(";B;",5)="; Fn(B, 5);TAB(32);A(1);A(2);A(3);A(4);A(5)

で表示されるA(1);A(2);A(3);A(4);A(5)が正しい再帰関係の結果として
A(1)=A(2)=A(3)=A(4)=A(5)=24
となったのかのチェックが出来ない事です。

試しに各段での利得関数として
DEF Sn(An)= 1*SIN(An)/2
と有るのを、段に応じて重み付けを加えて
DEF Sn(An,k)= 1*SIN(An)/2*k
と変更し、
f=Sn(i)+ Fn(B-i,k-1)

f=Sn(i,k)+ Fn(B-i,k-1)
と変更したら、全ての角度が等しく無い、と言う結果を予想していたのですが、プログラムを走らせると角度が
等しいと言う結果は変わりません。

私の安直な予想結果が誤っているのであれば良いのですが。

     ├私にも、わからない、どなたか説明できる方... SECOND 2007/07/13 03:41:49  ツリーへ

Re: 鮮やかなプログラミングですね! ノートメニュー
SECOND <cszcthjjdj> 2007/07/13 03:41:49
私にも、わからない、どなたか説明できる方。

     └!これで、試して下さい。 SECOND 2007/07/13 17:48:02  ツリーへ

Re: 鮮やかなプログラミングですね! ノートメニュー
SECOND <cszcthjjdj> 2007/07/13 17:48:02
!これで、試して下さい。

!Fn 面積最大時に、固定されるはずの角度 A() が、
!その後のステージで、次々と重ね書きされ、破壊されるバグがありました。
!応急措置です。


OPTION ANGLE DEGREES
DEF Sn(An,k)= 1*SIN(An)/2
!DEF Sn(An,k)= 1*SIN(An)/2*k
!DEF Sn(An,k)= 1*SIN(An)/2/k

LET B=120 !B=100( FOR i=15 TO B )、B=120( FOR i=19 TO B )にすると速くなる。

DIM A(6)

PRINT "Fn(";B;",1)="; Fn(B, 1);TAB(32);A(1)
PRINT "Fn(";B;",2)="; Fn(B, 2);TAB(32);A(1);A(2)
PRINT "Fn(";B;",3)="; Fn(B, 3);TAB(32);A(1);A(2);A(3)
PRINT "Fn(";B;",4)="; Fn(B, 4);TAB(32);A(1);A(2);A(3);A(4)
PRINT "Fn(";B;",5)="; Fn(B, 5);TAB(32);A(1);A(2);A(3);A(4);A(5)
PRINT "Fn(";B;",6)="; Fn(B, 6);TAB(32);A(1);A(2);A(3);A(4);A(5);A(6)


FUNCTION Fn(B,k)
local i,fm,im(6)
IF k=1 THEN
LET Fn=Sn(B,k)
LET A(k)=B
ELSE
LET fm=0
FOR i=19 TO B !------------------- 本来 0 TO B 、速度に難点のため。
LET f=Sn(i,k)+ Fn(B-i,k-1) !--- 再帰 FUNCTION
IF fm<f THEN
LET fm=f
MAT im=A !------------------ 配列A(6)の全セーブ
LET im(k)=i
END IF
NEXT i
LET Fn=fm
MAT A=im
END IF
END FUNCTION

END

      └演算結果の比較 島村1243 2007/07/13 22:18:52  (修正2回) ツリーへ

Re: !これで、試して下さい。 ノートメニュー
島村1243 <bjllmpcujp> 2007/07/13 22:18:52 ** この記事は2回修正されてます
演算結果の比較

SECONDさん、ご教示有難う御座います。
修正版のプログラムコードの利得関数として2番目を有効
!DEF Sn(An,k)= 1*SIN(An)/2
DEF Sn(An,k)= 1*SIN(An)/2*k
!DEF Sn(An,k)= 1*SIN(An)/2/k
にし、RUNさせてみました。

その結果は
Fn( 120 ,1)= .433012701892219 120
Fn( 120 ,2)= 1.3228733657107 49 71
Fn( 120 ,3)= 2.1650635094611 0 60 60
Fn( 120 ,4)= 2.89254424358943 0 40 40 40
Fn( 120 ,5)= 3.5 0 30 30 30 30
となり、第1段目角度は0になって尤らしくなりましたが
2段以降の角度は同じ値になってしまい、再帰計算が正
常に行われていないことを示してます。

試しに山中さんの作成されたプログラムを、各段の利得関数に、下記の様に重みを付けてRUNしました。

2段目の重み付けは
LET v=S1(j)*2+ S1(i-j) !1分割+1分割
3段目の重み付けは
LET v=S1(j)*3+ S2(i-j) !1分割+2分割
4段目の重み付けは
LET v=S1(j)*4+ S3(i-j) !1分割+3分割
5段目の重み付けは
LET v=S1(j)*5+ S4(i-j) !1分割+4分割
としました。

その結果は
117°: 55° 3.90072234
118°: 55° 3.92561916
119°: 55° 3.95044590
120°: 56° 3.97515973
θ5= 56
θ4= 45
θ3= 19
θ2= 0
θ1= 0
となって正しい値を計算していました。

SECONDさんの再帰式では、
k=1の場合はFn(B-i,k-1)は「存在しない第0段」の計算を指示している事になり、また
B-i=0の場合も、残り角度0の場合の面積計算を支持して
いることになっています。

この所に正解が出ない原因があるかと思って
FOR i=19 TO B !---- 本来 0 TO B 速度に難点のため
 IF K=1 OR b-i THEN f=Sn(i,k)
 LET f=Sn(i,k)+ Fn(B-i,k-1) !--- 再帰 FUNCTION
 IF fm<f THEN
と修正を加えてRunして見ましたが、結果はNGでした。

SECONDさんのプログラムは極めて簡潔に纏まっていて、利得関数が式として与えられた問題の場合は非常に効率
的なので魅力があるので何とか解決できると良いのですが。。。

       └!SIN(An)/2/k,SIN(An)/2は、山中氏の結果に... SECOND 2007/07/13 23:59:52  (修正3回) ツリーへ

Re: 演算結果の比較 ノートメニュー
SECOND <cszcthjjdj> 2007/07/13 23:59:52 ** この記事は3回修正されてます
!SIN(An)/2/k, SIN(An)/2 は、山中氏の結果に一致するようですが、
!SIN(An)/2*k にすると、ズレる。これがわからない。

!----IF k=1 THEN の付いたリストでは、
!for~next 内は K=2までで終わります。理由は、
!その時 (K-1)で再帰、K=1で、再コールしますが、これは、ネストの底になり、
!Fn=Sn(B) を持って、それ以上の再帰をせずにリターン出来るからです。
!K=0 で入られる心配は、ありません。

!----IF k=1 THEN の付いていないリストでは、
!for~next 内は K=1 まで行きます。(K-1)で再帰、K=0 で、やって来ますので、
!Fn=0 を冒頭において、そのままリターンさせ、ネストの底としています。


OPTION ANGLE DEGREES
!DEF Sn(An,k)= 1*SIN(An)/2
!DEF Sn(An,k)= 1*SIN(An)/2*k
DEF Sn(An,k)= 1*SIN(An)/2/k

LET B=120 !B=100( FOR i=15 TO B )、B=120( FOR i=19 TO B ) は、一時中止。

DIM A(5)

PRINT "Fn(";B;",1)="; Fn(B, 1);TAB(32);A(1)
PRINT "Fn(";B;",2)="; Fn(B, 2);TAB(32);A(1);A(2)
PRINT "Fn(";B;",3)="; Fn(B, 3);TAB(32);A(1);A(2);A(3)
PRINT "Fn(";B;",4)="; Fn(B, 4);TAB(32);A(1);A(2);A(3);A(4)
PRINT "Fn(";B;",5)="; Fn(B, 5);TAB(32);A(1);A(2);A(3);A(4);A(5)


FUNCTION Fn(B,k)
local i,fm,im(5)
IF k=1 THEN
LET Fn=Sn(B,k)
LET A(k)=B
ELSE
LET fm=0
FOR i=0 TO B !------------------- 本来の 0 TO B に、戻す、速度に難点。
LET f=Sn(i,k)+ Fn(B-i,k-1) !--- 再帰 FUNCTION
IF fm<f THEN
LET fm=f
MAT im=A !------------------ 配列A(5)の全セーブ
LET im(k)=i
END IF
NEXT i
LET Fn=fm
MAT A=im
END IF
END FUNCTION

END

        └たしかに 島村1243 2007/07/14 08:30:14  ツリーへ

Re: !SIN(An)/2/k,SIN(An)/2は、山中氏の結果に... ノートメニュー
島村1243 <bjllmpcujp> 2007/07/14 08:30:14
たしかに

SECONDさん、ご回答有難う御座います。

「SIN(An)/2/k, SIN(An)/2 は、山中氏の結果に一致する」

上記の場合は一致する事を確認しました。
また、再帰式の場合は「k=0になると自動的にFn(B,k)はゼロになる」と言う知識が私にありませんでしたので
コードを付け加えてしまいましたが、SECONDさんの説明で良く分かりました。

「SIN(An)/2*k にすると、ズレる。これがわからない。」
上記についてプログラムの進行状況がどうなっているかを紙に書いて調べてみます。
正解が出ない状況が分かりましたら報告します。

         └十進BASICのバグの疑いがあります。というの... SECOND 2007/07/14 10:42:40  ツリーへ

Re: たしかに ノートメニュー
SECOND <cszcthjjdj> 2007/07/14 10:42:40
十進BASICのバグの疑いがあります。というのは,uBASIC で、
同じ文を、仕様の違いだけを、直して走らせて見ました。
以下の様になりました。

10 'asave"070714-1"
20 key 8,"system"+chr(13)
30 'def fnSn(An,K)=1*sin(#pi*An/180)/2
40 def fnSn(An,K)=1*sin(#pi*An/180)/2*K
50 'def fnSn(An,K)=1*sin(#pi*An/180)/2/K
60 '
70 B=120
80 '
90 dim A(5)
100 '
110 print "Fn(";B;",1)=";fnFFn(B,1);tab(32);A(1)
120 print "Fn(";B;",2)=";fnFFn(B,2);tab(32);A(1);A(2)
130 print "Fn(";B;",3)=";fnFFn(B,3);tab(32);A(1);A(2);A(3)
140 print "Fn(";B;",4)=";fnFFn(B,4);tab(32);A(1);A(2);A(3);A(4)
150 print "Fn(";B;",5)=";fnFFn(B,5);tab(32);A(1);A(2);A(3);A(4);A(5)
160 end
170 '
180 fnFFn(B,K)
190 local I,Fm,Im1,Im2,Im3,Im4,Im5
200 if K=1 then
210 :FFn=fnSn(B,K)
220 :A(K)=B
230 :else
240 :Fm=0
250 :for I=0 to B
260 :F=fnSn(I,K)+fnFFn(B-I,K-1)
270 :if Fm<F then
280 :Fm=F
290 :A(K)=I
300 :Im1=A(1):Im2=A(2):Im3=A(3):Im4=A(4):Im5=A(5)
310 :endif
320 :next I
330 :FFn=Fm
340 :A(1)=Im1:A(2)=Im2:A(3)=Im3:A(4)=Im4:A(5)=Im5
350 :endif
360 return(FFn)
370 '


run
Fn( 120 ,1)= 0.4330127018922193233 120
Fn( 120 ,2)= 0.8660254037844386467 60 60
Fn( 120 ,3)= 0.9641814145298089895 40 40 40
Fn( 120 ,4)= 0.9999999999999999997 30 30 30 30
Fn( 120 ,5)= 1.0168416076895005193 24 24 24 24 24
OK

run
Fn( 120 ,1)= 0.4330127018922193233 120
Fn( 120 ,2)= 1.3228733657107028092 49 71
Fn( 120 ,3)= 2.1793927902259533372 0 53 67
Fn( 120 ,4)= 3.0549519622798799463 0 10 49 61
Fn( 120 ,5)= 3.9751597254464342817 0 0 19 45 56
OK

run
Fn( 120 ,1)= 0.4330127018922193233 120
Fn( 120 ,2)= 0.6614366828553514046 71 49
Fn( 120 ,3)= 0.6614542285752208826 71 48 1
Fn( 120 ,4)= 0.6614542285752208826 71 48 1 0
Fn( 120 ,5)= 0.6614542285752208826 71 48 1 0 0
OK

          └FullBASICの定義関数名は変数ではなく,従っ... 白石 和夫 2007/07/14 11:27:18  (修正1回) ツリーへ

Re: 十進BASICのバグの疑いがあります。というの... ノートメニュー
白石 和夫 <fbdfvqwhki> 2007/07/14 11:27:18 ** この記事は1回修正されてます
Full BASICの定義関数名は変数ではなく,従って,関数呼び出しにおいて局所的でない,というのは,大丈夫ですか?
十進BASICの初期のバージョンでは定義関数名を局所変数として扱っていたのですが,最近のバージョンは規格に合うように修正されています。

           ├できれば、この一件を、トレースして頂けれ... SECOND 2007/07/14 12:46:05  ツリーへ

Re: FullBASICの定義関数名は変数ではなく,従っ... ノートメニュー
SECOND <cszcthjjdj> 2007/07/14 12:46:05
できれば、この一件を、トレースして頂ければ幸いです。
もう、限界に近いほど、確かめましたが、不明です。

           ├LOCAL文の不具合かも? 山中和義 2007/07/14 13:28:20  ツリーへ

Re: FullBASICの定義関数名は変数ではなく,従っ... ノートメニュー
山中和義 <drdlxujciw> 2007/07/14 13:28:20
LOCAL文の不具合かも?

内部手続き(Fn)と外部手続き(Fn3)で動作が違います。
外部手続きの方は、正しく計算されています。
プログラムは、n=4としています。


SUB Fn(n,B, th(),maxFn) !Max{Sn(θn)+Fn-1(β-θn)}
local i,v,vmax !※本サブルーチンをEXTERNAL SUBとしてもよい
local sav_th1,sav_th2,sav_th3,sav_th4,sav_th5

IF n=1 THEN !1分割の場合
LET th(1)=B
LET maxFn=SIN(RAD(B))/2
ELSE !n分割の場合
LET maxFn=0
FOR i=0 TO B !1°刻み
CALL Fn(n-1,B-i, th,vmax) !Fn-1
LET v=SIN(RAD(i))/2*n+vmax
IF v>maxFn THEN !Max{Sn+Fn-1}
LET th(n)=i
LET sav_th1=th(1)
LET sav_th2=th(2)
LET sav_th3=th(3)
LET sav_th4=th(4)
!LET sav_th5=th(5)
LET maxFn=v
END IF
NEXT i
LET th(1)=sav_th1
LET th(2)=sav_th2
LET th(3)=sav_th3
LET th(4)=sav_th4
!LET th(5)=sav_th5
END IF
END SUB


LET t0=TIME

LET β=120 !角度βの扇

LET N=4 !角度θ1、θ2、θ3、θ4、θ5に分ける
DIM th(N)

MAT th=ZER
CALL Fn(N,β, th,maxFn)
FOR i=1 TO N
PRINT th(i);
NEXT i
PRINT ,maxFn


MAT th=ZER
CALL Fn3(N,β, th,maxFn)
FOR i=1 TO N
PRINT th(i);
NEXT i
PRINT ,maxFn


END


EXTERNAL SUB Fn3(n,B, th(),maxFn) !Max{Sn(θn)+Fn-1(β-θn)}
IF n=1 THEN !1分割の場合
LET th(1)=B
LET maxFn=SIN(RAD(B))/2
ELSE !n分割の場合
LET maxFn=0
FOR i=0 TO B !1°刻み
CALL Fn3(n-1,B-i, th,vmax) !Fn-1
LET v=SIN(RAD(i))/2*n+vmax
IF v>maxFn THEN !Max{Sn+Fn-1}
LET th(n)=i
LET sav_th1=th(1)
LET sav_th2=th(2)
LET sav_th3=th(3)
LET sav_th4=th(4)
!LET sav_th5=th(5)
LET maxFn=v
END IF
NEXT i
LET th(1)=sav_th1
LET th(2)=sav_th2
LET th(3)=sav_th3
LET th(4)=sav_th4
!LET th(5)=sav_th5
END IF
END SUB


           │└ありがとうございました。 SECOND 2007/07/14 14:30:48  ツリーへ

Re: LOCAL文の不具合かも? ノートメニュー
SECOND <cszcthjjdj> 2007/07/14 14:30:48
ありがとうございました。

           └!問題が表面化するテストプログラムです。 SECOND 2007/07/14 14:24:49  ツリーへ

Re: FullBASICの定義関数名は変数ではなく,従っ... ノートメニュー
SECOND <cszcthjjdj> 2007/07/14 14:24:49
!問題が表面化するテストプログラムです。

!SIN()/2 ・・・一致します。
!SIN()/2*k ・・・一致しません。外部 Fnx の方が正常値です。
!SIN()/2/k ・・・一致します。



OPTION ANGLE DEGREES
LET B=120

PRINT "Fn(";B;",1)="; Fn(B, 1)
PRINT "Fn(";B;",2)="; Fn(B, 2)
PRINT "Fn(";B;",3)="; Fn(B, 3)
PRINT
PRINT "Fnx(";B;",1)="; Fnx(B, 1)
PRINT "Fnx(";B;",2)="; Fnx(B, 2)
PRINT "Fnx(";B;",3)="; Fnx(B, 3)


FUNCTION Fn(B,k)
local i,fm
IF k=1 THEN
LET Fn=SIN(B)/2*k
ELSE
LET fm=0
FOR i=0 TO B
LET f=SIN(i)/2*k+ Fn(B-i,k-1)
IF fm<f THEN
LET fm=f
END IF
NEXT i
LET Fn=fm
END IF
END FUNCTION

END


EXTERNAL FUNCTION Fnx(B,k)
OPTION ANGLE DEGREES
IF k=1 THEN
LET Fnx=SIN(B)/2*k
ELSE
LET fm=0
FOR i=0 TO B
LET f=SIN(i)/2*k+ Fnx(B-i,k-1)
IF fm<f THEN
LET fm=f
END IF
NEXT i
LET Fnx=fm
END IF
END FUNCTION

            ├FORv=始値TO限界STEP増分 白石 和夫 2007/07/14 16:35:52  ツリーへ

Re: !問題が表面化するテストプログラムです。 ノートメニュー
白石 和夫 <fbdfvqwhki> 2007/07/14 16:35:52
FOR v=始値 TO 限界 STEP 増分
(一連の文)
NEXT v
は,
LET own1=限界
LET own2=増分
LET v=始値
DO UNTIL (v - own1)*SGN(own2)>0
(一連の文)
LET v=v+own2
LOOP
で定義されています。
own1,own2はプログラム単位の変数なので,内部手続きに対して局所的ではありません。
これが原因の可能性があります。

            │├簡単な例 白石 和夫 2007/07/14 16:43:42  ツリーへ

Re: FORv=始値TO限界STEP増分 ノートメニュー
白石 和夫 <fbdfvqwhki> 2007/07/14 16:43:42
簡単な例
DECLARE EXTERNAL SUB t
SUB s(n)
local i
FOR i=1 TO n
CALL s(n-1)
PRINT n,i
NEXT i
END SUB
CALL s(3)
PRINT
CALL t(3)
END
EXTERNAL SUB t(n)
FOR i=1 TO n
CALL t(n-1)
PRINT n,i
NEXT i
END SUB

            │└こちらはown1を局所化したものです。 白石 和夫 2007/07/14 16:47:44  ツリーへ

Re: FORv=始値TO限界STEP増分 ノートメニュー
白石 和夫 <fbdfvqwhki> 2007/07/14 16:47:44
こちらはown1を局所化したものです。
DECLARE EXTERNAL SUB t
SUB s(n)
local i,own1
LET own1=n
LET i=1
DO UNTIL i>own1
CALL s(n-1)
PRINT n,i
LET i=i+1
LOOP
END SUB
CALL s(3)
PRINT
PRINT
CALL t(3)
END
EXTERNAL SUB t(n)
FOR i=1 TO n
CALL t(n-1)
PRINT n,i
NEXT i
END SUB

            └以下のようにすると,同じ動作になります。 白石 和夫 2007/07/14 17:35:42  ツリーへ

Re: !問題が表面化するテストプログラムです。 ノートメニュー
白石 和夫 <fbdfvqwhki> 2007/07/14 17:35:42
以下のようにすると,同じ動作になります。
OPTION ANGLE DEGREES
LET B=120

PRINT "Fn(";B;",1)="; Fn(B, 1)
PRINT "Fn(";B;",2)="; Fn(B, 2)
PRINT "Fn(";B;",3)="; Fn(B, 3)
PRINT
PRINT "Fnx(";B;",1)="; Fnx(B, 1)
PRINT "Fnx(";B;",2)="; Fnx(B, 2)
PRINT "Fnx(";B;",3)="; Fnx(B, 3)


FUNCTION Fn(B,k)
local i,fm,own1
IF k=1 THEN
LET Fn=SIN(B)/2*k
ELSE
LET fm=0
LET own1=B
LET i=0
DO UNTIL i>own1
LET f=SIN(i)/2*k+ Fn(B-i,k-1)
IF fm<f THEN
LET fm=f
END IF
LET i=i+1
LOOP
LET Fn=fm
END IF
END FUNCTION

END


EXTERNAL FUNCTION Fnx(B,k)
OPTION ANGLE DEGREES
IF k=1 THEN
LET Fnx=SIN(B)/2*k
ELSE
LET fm=0
FOR i=0 TO B
LET f=SIN(i)/2*k+ Fnx(B-i,k-1)
IF fm<f THEN
LET fm=f
END IF
NEXT i
LET Fnx=fm
END IF
END FUNCTION

             ├ありがとうございました。 SECOND 2007/07/14 17:46:32  ツリーへ

Re: 以下のようにすると,同じ動作になります。 ノートメニュー
SECOND <cszcthjjdj> 2007/07/14 17:46:32
ありがとうございました。

             └再帰式中ではForNextが使えない? 島村1243 2007/07/14 22:14:08  (修正2回) ツリーへ

Re: 以下のようにすると,同じ動作になります。 ノートメニュー
島村1243 <bjllmpcujp> 2007/07/14 22:14:08 ** この記事は2回修正されてます
再帰式中ではFor Nextが使えない?

白石先生が示されたFn(B,k)の中身を下記のように変更して
実行しましたら、正解が出ました。

OPTION ANGLE DEGREES
LET B=120
DEF s(x)=sin(x)/2 !<---SECONDさん流に戻した。

PRINT "Fn(";B;",1)="; Fn(B, 1)
PRINT "Fn(";B;",2)="; Fn(B, 2)
PRINT "Fn(";B;",3)="; Fn(B, 3)
PRINT

FUNCTION Fn(B,k)
local i,fm !,own1 <---own1を無効化
IF k=1 THEN
LET Fn=s(B)*k !<---ここ{SIN(B)/2}を変更
ELSE
LET fm=0
!LET own1=B <---ここを無効化
LET i=0
DO UNTIL i>B !<---ここown1をBに戻した
LET f=s(i)*k+ Fn(B-i,k-1) !<---ここ{SIN(i)/2}を変更
IF fm<f THEN
LET fm=f
END IF
LET i=i+1
LOOP
LET Fn=fm
END IF
END FUNCTION

次にDO LOOP 文を下記のようにFOR NEXT 文に変えたら結果は異なりました。

!LET i=0
for i=0 to B !>own1
LET f=s(i)*k+ Fn(B-i,k-1) !<---ここ{SIN(i)/2}を変更
IF fm<f THEN
LET fm=f
END IF
!LET i=i+1
!LOOP
Next i

皆様のご検討結果を元にしたこの試行結果から、
「十進BASICでは再帰式中においてはFor Next文は使えず、
DO LOOPを使わなければならない。DEF文は使える。」
が結論と理解しました。

              ├「内部関数定義,内部副プログラムは原則と... 白石 和夫 2007/07/15 07:09:57  (修正1回) ツリーへ

Re: 再帰式中ではForNextが使えない? ノートメニュー
白石 和夫 <fbdfvqwhki> 2007/07/15 07:09:57 ** この記事は1回修正されてます
「内部関数定義,内部副プログラムは原則として使わない」が正解です。
外部手続きとモジュールの機能を利用すれば,JIS規格外の命令であるLOCAL文を使わずに,確実に動作するプログラムが作れます。

              │└内部関数定義や内部副プログラムは,外部手... 白石 和夫 2007/07/15 09:17:49  ツリーへ

Re: 「内部関数定義,内部副プログラムは原則と... ノートメニュー
白石 和夫 <fbdfvqwhki> 2007/07/15 09:17:49
内部関数定義や内部副プログラムは,外部手続き定義とモジュールを使ったプログラムを書いていて,外部手続き中の中に手続きを書かなければならなくなったときにやむをえず使うものです。
C言語にはFull BASICの内部手続きに相当する概念がないのですが,それでも巨大なプログラムがCを使って書かれています。Full BASICのプログラムで内部関数定義,内部副プログラムを使う必要が生じることはめったにないと思います。

              │ └大変良く分かりました。 島村1243 2007/07/15 22:45:21  ツリーへ

Re: 内部関数定義や内部副プログラムは,外部手... ノートメニュー
島村1243 <bjllmpcujp> 2007/07/15 22:45:21
大変良く分かりました。

白石先生、お忙しいにも関わらず丁寧なご教示有難うございました。

              └お疲れさまでした。代数的に、n-1番目の結果... SECOND 2007/07/15 11:40:44  ツリーへ

Re: 再帰式中ではForNextが使えない? ノートメニュー
SECOND <cszcthjjdj> 2007/07/15 11:40:44
お疲れさまでした。代数的に、n-1 番目の結果から n 番目 を
考えるという概念は、とても数学的で、楽しい事だと、思います。
こういう問題で、悩める時間が、永く持てればいいですね。

               └大変勉強になりました。 島村1243 2007/07/15 23:09:14  (修正3回) ツリーへ

Re: お疲れさまでした。代数的に、n-1番目の結果... ノートメニュー
島村1243 <bjllmpcujp> 2007/07/15 23:09:14 ** この記事は3回修正されてます
大変勉強になりました。

SECONDさん、とても長い時間ご検討くださいまして本当
に有難うございました。
過去にDPの数学参考書でその意味合いを学んだつもり
だったのですが、実際に数値計算して確認した訳では無
かったので、真に理解していたとは言えない状態でした

今回、山中さんとSECONDさんのDP数値計算の具体的
プログラムを見せていただく事によって、そのテクニッ
クを勉強させて頂いただけでなく、DPそのものも良く
理解する事が出来ました。有難う御座いました。
今後も、数値計算を必要とする数学的に面白そうなテー
マが見つかりましたら投稿しますので、宜しくご指導お
願い致します。


このノートはこれ以上発言できません。
新しくノートを作成 して、続きを書いてください。

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