QEの計算アルゴリズムとその応用―数式処理による最適化

QEの計算アルゴリズムとその応用―数式処理による最適化

  • 作者: 穴井宏和,横山和弘
  • 出版社/メーカー: 東京大学出版会
  • 発売日: 2011/08/25
  • メディア: 単行本
  • 購入: 1人 クリック: 7回
  • この商品を含むブログを見る
 

  

にほんブログ村 科学ブログ 数学へ
にほんブログ村

 

更に続けて、上記の本のp37に紹介されている京大の入試問題。線形計画法と(線形ではない)計画法の問題となります。

「座標平面上で点P(x,y)が ( 4,x+y<=0, x+2,y>=4, 2,x-3,y>=-6 ) の範囲を動く時 ( 2,x+y, x^2+y^2 ) のそれぞれの最大値、最小値を求めよ。」

早速qepmaxを読み込みます。また条件にRという名前を付けます。

(%i1) load(qepmax)$

(%i2) R:4*x+y<=9 %and x+2*y>=4 %and 2*x-3*y>=-6;

$$ ag{%o2} left(2,x-3,ygeq -6 ight) wedge left(y+4,xleq 9 ight) wedge left(2,y+xgeq 4 ight) $$

この手の定番ですが、(x,y)がRに存在する時k1を自由変数として(k1=2,x+y)、とおいて、Q.E.するとk1の条件(範囲)を求めることができます。

(%i3) qe([[E,x],[E,y]],R %and k1=2*x+y);

$$ ag{%o3} left(mathrm{k1}-6leq 0 ight) wedge left(mathrm{k1}-2geq 0 ight) $$

同値な条件ですから、k1の最小値は2、最大値は6と求まります。ここではついでにk1=2の時のx,yを求めてみます。今度はk1の値を2と固定して、そうなるx,yの条件を求めるため、x,yは自由変数とします。

(%i4) qe([ ], R %and k1=2*x+y ),k1=2;

$$ ag{%o4} left(x=0 ight) wedge left(3,y-2,x-6=0 ight) $$

これは連立方程式ですから、Maximaで解くことができます。(目視でも簡単ですが)。

(%i5) solve([part(%,1),part(%,2)],[x,y]);

$$ ag{%o5} left[ left[ x=0 , y=2 ight] ight] $$

 

では同様に(x^2+y^2)の最大、最小も求めてみましょう。やり方は同じです。

(%i6) qe([[E,x],[E,y]],R %and k2=x^2+y^2);

$$ ag{%o6} left(4,mathrm{k2}-45leq 0 ight) wedge left(5,mathrm{k2}-16geq 0 ight) $$

これで最小値(frac{16}{5})、最大値(frac{45}{4})が求まります。ついでに最大値を与える(x,y)を求めてみます。これもやり方は同じです。

(%i7) qe([ ],R %and k2=x^2+y^2),k2=45/4;

$$ ag{%o7} left(2,x-3=0 ight) wedge left(3,y-2,x-6=0 ight) $$

(%i8) solve([part(%,1),part(%,2)],[x,y]);

$$ ag{%o8} left[ left[ x=frac{3}{2} , y=3 ight]   ight]  $$

連立方程式を解いて、最大値(frac{45}{4})を与えるx,yの値が求まりました。

(%i9) [%,45/4],numer;

$$ ag{%o9} left[ left[ left[ x=1.5 , y=3 ight]   ight]  , 11.25 ight]  $$

あとの計算と比較するために答えを浮動小数点に変換しました。

 

ところで、(2,x+y)の最大、最小を求めるのはまさに線形計画法の問題です。Maximaには線形計画法のパッケージsimplexがありますから、それで解いてみましょう。まずはsimplexパッケージを読み込みます。

(%i10) load(simplex)$

Rの3つの線形不等式を取り出して制約として指定し、minimize_lp()を使って(2,x+y)の最小値を求めます。

(%i11) minimize_lp(2*x+y, [part(R,1),part(R,2),part(R,3)]);

$$ ag{%o11} left[ 2 , left[ y=2 , x=0 ight]   ight]  $$

先ほど求めた結果と同じですね。

 

(線形ではない)計画法の問題も別の方法で解いてみます。制約条件を満たしながらどんな式でも最小値を見つける数値計算のパッケージcobylaを使います。

(%i12) load(fmin_cobyla)$

関数fmin_cobyla()は最小値しか見つけられません。(x^2+y^2)の最大値を見つけてみたいのですが、、、(-x^2-y^2)の最大値を見つけてもらえば良い訳です。

(%i13) fmin_cobyla(-(x^2+y^2),[x,y],[1,1],constraints=[part(R,1),part(R,2),part(R,3)]);

$$ ag{%o13} left[ left[ x=1.5 , y=3.0 ight]  , -11.25 , 20 , 0 ight]  $$

これもさきほどQ.E.で求めた結果と同じになります。