力学的モデル

本技法では力学的モデルにより，職場を移動させ相対位置関係を得る．これまでに様々な力学的モデル の諸変形が提案されているが，どのモデルにも共通しているのは職場同士が離れていれば引力を発生さ せ，近づきすぎていれば斥力を発生させる事で最適な配置を得ようとする構造である．本研究では後に引 き続くカオス軌道を発生させるため，差分化勾配系と呼ばれるモデルに合わせた形で構成する力学的モデ ルを提案する．

以下に引力及び斥力の導出を示す．まずは一般的なエネルギー最小化問題に対する最急降下勾配系,及 び差分化勾配系について説明を行う(最急降下勾配系の説明は相吉,安田[50]を引用)．その後,本研究で 提案する「職場の円近似」「配置領域制約の無視」を組み込んだFLPの近似問題に対するエネルギー関数 として,^{目的関数と}l2型ペナルティ関数を用いたエネルギー最小化問題を定義する．そしてこのエネル ギー最小化の勾配系より式(4.28)-(4.35)で定義される引力及び斥力を示す．

最急降下勾配系[50],差分化勾配系

連続微分可能な関数E :Rⁿ →R^{を最小化する問題}

min. E(x) (4.11)

を考える．ただし,x= [x1, x2,· · ·, xn]^T とする．この問題を解くために，任意の初期点x(0) ∈R^{n から} 出発し,この問題の局所的最適解の１つx^∗に収束する軌道を与える力学系モデルとして,最も基本的な最 急降下勾配系モデルを与える．

任意の時刻t^{における位置}x(t)^{が与えられたとき},^その速度x(t)˙ ^方向の∆T ^{秒後の位置}x(t) + ∆Tx(t)˙ における関数E^の値E(x(t) + ∆Tx(t))˙ を最小にするような速度x(t)˙ ^{与える問題は},

min E(x(t) + ∆Tx(t))˙ (4.12)

で与えられる．ここで,E(x(t) + ∆Tx(t))˙ を変分∆Tx(t)˙ に関して一次近似し, E(x(t) + ∆Tx(t))˙

'E(x(t)) + ∆T∇E(x(t)) ˙x(t)) (4.13) とおくと,^式(4.13)^{は右辺の第}2項のみの最小化と等価であり,

min. ∇E(x(t)) ˙x(t)) (4.14)

s.t. ˙x(t)^TMx(t) = 1˙ (4.15)

と書き換えることができる．ただし,∇E(x(t))^は関数E^{の勾配とし},^{制約条件式}(4.15)^は有限のx(t)˙ ^を 求めるために課した単位ノルム制約であり,Mは正定値対称行列でこのノルムを定義する計量行列である．

式(4.14)をラグランジュの未定乗数法により解くと,^{一次近似の意味で関数}Eを最小にする速度として

∇x(t)) =˙ −1

φM⁻¹∇E(x(t))^T (4.16)

が得られる．ただし,φはラグランジュ乗数であり

φ= (∇E(x(t))M⁻¹+∇E(x(t))^T)^1/2 (4.17)

第4章 カオス最適化と線形計画法を用いたヒューリスティクスの提案 71 である．速度はその方向が重要であるから,1/φ=c(一定)として,式(4.16)を微分方程式として書き換え ると,

dx(t)

dt =−cM⁻¹∇E(x(t))^t, c >0 (4.18) が得られ,M =I(^単位行列)^{とした連続時間力学系}

dx(t)

dt =−c∇E(x(t))^t, c >0 (4.19) が現時刻tにおける勾配∇E(x(t))で駆動する最急降下系である．これをEuler差分化し∆T =h, c= 1 とおくと式(4.20)^{の様になる．}

x(t) =x(t+ 1)−h∇E(x(t)) (4.20)

エネルギー関数

(4.1)-(4.9)^{で示される}FLP^に対し,^{「職場の円近似」}「配置領域制約の無視」を組み込んだ近似問題と

して下記の最小化問題(4.21)(4.22)^{を定義する．}

minimize

∑N i=1

∑N j>i

fijdij (4.21)

subject to Rij−Dij ≤0 (4.22) 但し,距離dij マンハッタン距離,Dij はユークリッド距離を示す．ここで制約式(4.22)を目的関数

(4.21)に組み込んだ下記の問題を考える．

minimize E(r) (4.23)

E(r) =

∑N i=1

∑N j>i

fijdij+Gij(r) (4.24)

Gij(r) = (p/2)× {max(gij(r),0)}² (4.25)

gij(r) =Rij−Dij (4.26)

ここで,式(4.25)はl2型ペナルティ関数と呼ばれp→ ∞^{で最適化問題}(4.21)(4.22)と等価な問題と なる．ペナルティ関数(4.25)は制約を違反した場合のみ正の値をとり,そうでない場合は0となる．ま た,2乗の形になっているのは微分可能な形にするための処理である．

ここで,E(r)^{の勾配ベクトル}

∇E(r) =(

∂E(r)

∂x1 ,· · · ,^∂E(r)_∂xN ,^∂E(r)_∂y1 ,· · ·,^∂E(r)_∂yN )T

に対し,

−∇E(r) = ∆r

(4.23) E(r) xi ,

∂E(r)

∂xi

= ∂

∂xi

∑N j=1

(

fijdij+Gij(r) )

=

∑N j=1

( ∂

∂xi

fijdij+ ∂

∂xi

Gij(r) )

yiについても同様なので,それぞれ第1項,第2項を以下の様に定義できる．

[ ∂

∂xifijdij

∂

∂yifijdij

]

=−∆r^A_ij, [ ∂

∂xiGij(r)

∂

∂yiGij(r) ]

=−∆r_ij^R

ここで第1項−∆r^A_ijは職場同士が引合う方向を表すため「引力」，第2項−∆r_ij^Rは職場同士の重なり を排除する方向を表すため「斥力」と解釈する事ができる．以降より,それぞれの項について詳細に計算 を行う．

引力の導出

fijdij のxiによる偏微分を示す．Chain Ruleより,

∑N j=1

∂

∂xi

fijdij =

∑N j=1

∂( fijdij

)

∂dij

∂dij

∂xi

=

∑N j=1

fij

( ∂

∂xi

|xi−xj|)

ここで,|xi−xj|^はxi=xj 付近で微分不可能であるため下記の近似を行う．

|xi−xj| '√

(xi−xj)²+µ²−µ ここでµ^{は十分小さい数}(e.g.10⁻⁵)^{とする．これより},

∂

∂xi

(√

(xi−xj)²+µ²−µ )

= xi−xj

√(xi−xj)²+µ² xi−xj がµよりも十分大きい場合を考えると,µ→0^{となるので},

xi−xj

√(xi−xj)²+µ² → xi−xj

√(xi−xj)²

= xi−xj

|xi−xj|

= sign(xi−xj)

が導き出される．yiについても同様なので以下の様に引力を導出する事が出来る．

∆r^A_ij =fij×sign

([xj −xi

yj −yi

])

第4章 カオス最適化と線形計画法を用いたヒューリスティクスの提案 73

斥力の導出

Gij(r)^のxiによる偏微分を考える．gij(r)≤0^ならばGij(r) = 0^{より偏微分も}0^{となる．従ってこ} こではgij(r)>0の場合のみ考える．この時,

∂Gij(r)

∂xi

= ∂

∂xi

p

2(Rij−Dij)²

=−p(Dij−Rij)∂Dij

∂xi

(4.27) ここで,

∂Dij

∂xi

= ∂

∂xi

√

(xi−xj)²+ (yi−yj)²

= √ xi−xj

(xi−xj)²+ (yi−yj)²

= (xi−xj)/Dij

従って,^式(4.27)^に代入し,gij(r)≤0の場合とまとめると下記の様になる．

∂Gij(r)

∂xi

=

{ −p(^R_D^ij_ij −1)(xi−xj) (if gij(r)>0)

0 (if gij(r)≤0)

ここでgij(r)の正負を考えると,Dij >0より

gij(r)≤0 ⇐⇒ Rij−Dij ≤0

⇐⇒ (Rij

Dij −1)≤0

となる事を利用し,yiについてもまとめると以下の様に斥力を定義する事ができる．

∆r^R_ij =p×max(Rij

Dij −1,0)

[xi−xj

yi−yj

]

引力・斥力の定義

相異なる２職場i, j に対し，引力，及び斥力を下記の様に定義する．またその様子を図4.3^に示す．

∆r^A_ij =fij×sign

([xj −xi

yj −yi

])

(4.28)

∆r^R_ij =p×max(Rij

Dij −1,0)

[xi−xj

yi−yj

]

(4.29) sign(u) =

{ 1 (if u >0)

−1 (if u <0) (4.30)

Dij =

√

(xi−xj)²+ (yi−yj)² (4.31)

Rij = (Ri+Rj) (4.32)

• ∆r_ij^A,∆r_ij^R : ^職場i^が職場j^{から受ける引力},^斥力

• Dij : 職場ij間のユークリッド距離

• p : 重なりに対するペナルティ

式(4.28)の引力の導出は目的関数の1^階微分,(4.29)^{の斥力は職場}ij 間の重なりのペナルティ関数の1 階微分により得られている．

式(4.28)^{の引力に関しては},^係数fij は「バネ係数」と解釈でき，物流量の強い職場同士により強い

引力が発生し引合う効果がある．sign(u) は符号関数と呼ばれる入力の正負を判別する関数で,例えば (xi> xj, yi> yj)の時,引力は∆r^A_ij = [fij, fij]^T となる．

式 (4.29)^{の斥力に関しては},^{職場間距離} Dij が両者の半径の和 (Ri+Rj) ^{を上回った場合は係数} (_D^Rij

ij −1)が負となり引力と逆向きの斥力を発生させる．

これらを用いて，各職場i= 1· · ·N が全職場から受ける引力及び斥力を下記の様に定義する．

∆ri=

∑N j=1

(∆r^A_ij+ ∆r^R_ij) (4.33)

但し,∆ri : ^職場iが全職場から受ける引力/^{斥力の和とする．}

職場の移動（ダイナミクス）

上述の様に定義された引力・斥力を用いて下記の方程式に従い職場の移動を行う．

r(t+ 1) =r(t) +h∆r(t) (4.34)

∆r(t) = [∆r1,∆r2,· · ·,∆rN]^T (4.35)

• r(t) = [x,y]^T : t回目の移動後の職場の位置

• ∆r(t) : t回目の移動後の職場の移動方向

• h : 移動のステップ幅

第4章 カオス最適化と線形計画法を用いたヒューリスティクスの提案 75

ドキュメント内 カオス最適化を用いた職場レイアウト技法 A Chaos-based Optimization for Solving (ページ 75-80)