― 49 ―遺伝的アルゴリズムを用いた非線形振動系の解析
渡部雄

二＊・小林

義和Analysis of nonlinear vibration system using genetic algorithm
Yuji WATANABE ＊ and Yoshikazu KOBAYASHI（平成21年11月27日受理）
The nonlinear vibration system was analyzed by using the genetic algorithm, and accuracy
and analytical time of the analytical result were evaluated in the present study. When analyzing
it, the dimensionless values that depended on each parameter, composed the dynamic vibration
absorber was assumed, and changed to obtain the vibration frequency response curve ofnonlinear two degree of freedom system composed of nonlinear spring and the dynamic vibration
absorber. As a result, shortening analytical time and some accurate solution were obtained.

1.

緒言

2.

非線形問題は解が複雑になる場合があり，いわゆる準ニュートン法を用いた解法では局所的な解に収束してしまう。また，問題を解くのに多くの計算時間を要する。そこで遺伝的アルゴリズムを用いることで，複数の解を同時に探索でき，また実用に問題のないと考えられる精度の最適解を比較的短時間で求めることができるのではないかと考えた。本研究
では非線形ばねを使用した二自由度振動系の最適化を遺伝的アルゴリズム等の方法を用いて行い，その最適解の良否，また遺伝的アルゴリズムを用いた解析の有効性を検証することを目的とする。非線形ばねとは図 1 のようなばねを指し f ＝ kx では表すことができないばねである。ただし，f は復元力［N］，kはばね定数［N/m］，x は変位［m］である。今回の実験では f ＝ kx ±εx3という図 2 のような特性を持つ非線形ばねを仮定して解析を行った。

図1＊

非線形ばね例秋田高専専攻科学生

秋田高専研究紀要第45号

図2

非線形ばね特性

遺伝的アルゴリズム（GA）の原理（1）2.1 GA（遺伝的アルゴリズム）とはGA（遺伝的アルゴリズム）は，ある範囲内で定f義されている変数 x の関数 （x）の最大値あるいは最小値を与える x の値を，高速に求めるための最適化・検索アルゴリズムの一種である。GA は生物の進化の過程にヒントを得た比較的単純な基本原理を基にしているため，ほとんどあらゆる最適化・探索の問題に適用可能となっている。GA ではある情報を，遺伝子を持つ仮想的な生物集団として計算機内に設定し，あらかじめ定めておいた環境に適応している固体が子孫を残す確率が高くなるよう世代交代シミュレーションを実行し，遺伝子および生物集団を“進化”させる。このため，これらの仮想生物の進化によって，与えられた工学的課題の解が得られるように GA のプログラミングを行う必要がある。
GA は実際のプログラミングの詳細を規定していないので，各種の規則やパラメータの設定方法など，不確定要素が多い方法論であることが欠点として指摘されることが多い。しかしながら，むしろ厳しく規定されていないために応用範囲が広いということができる。
2.2 GA の基本的な考え方
GA では，探索空間中の探索点を 1 点ずつ順番に探索するのではなく，複数個の点を同時に探索する。
そして，各探索点を，遺伝子を持つ仮想的な生物で― 50 ―渡部雄二・小林義和あるとみなす。各個体に対して，あらかじめ設定されている環境との適応度が評価関数によって計算さfx を例とする場合，れる。一般的な関数 （x）と変数x を固体の遺伝子，f（x）を環境との適応度と考える。
低い適応度を持つ固体を淘汰して消滅させ，高い適応度を持つ固体を増殖させて親の形質を継承した遺伝子を持つ子孫の固体を生成する世代交代シミュ
レーションを実行する。この際，実際の生物の生殖でも生じる，遺伝子の交差，および突然変異と呼ばれる操作を行う。そして，最終的に最も高い適応度fの固体，言い替えれば最大値と考えられる （x）のx値を与える の値を求める。上述のような世代交代シミュレーションを分かりやすくするためにフローチャートにしたものが図 3 である。

（2）突然変異
各個体の遺伝子に相当する各ビットのどれかをランダムで選択し，そのビットの 0 を 1，あるいは 1を 0 に変更する操作を突然変異（mutation）と呼ぶ。
また，突然変異の生起確率を突然変異率（mutationrate）と呼ぶ。図 5 に突然変異の例を示す。

図5

突然変異例

この突然変異の操作によって，交差だけでは生じない遺伝子を持つ固体が生成され，探索の観点から見ると，現在の探索点から大きく離れた場所に探索点を発生させることになる。これは，生物集団が局所的な最適解に陥ったときに，そこから脱出する働きがある。このため進化シミュレーションに必要な操作ではあるが，突然変異率を大きくしすぎると，先に述べた遺伝子型の交差による親の形質の継承の特徴が失われ，探索空間中をランダムに探索することと同様になってしまうので，本研究では突然変異率を，0.1～5％程度の値とした。
図3

世代交代シミュレーション

3.
図 3 にある交差・突然変異は次のような作業を行う。
（1）交差
生成された固体の中から，2 つの固体のペアをランダムで選択し，それぞれに対して交差（crossover）と呼ぶ操作を実行する。ここで，交差が生じる確率を交差率（crossover rate）と呼ぶ。
交差は，2 つの固体の遺伝子型をランダムな位置で部分的に入れ換える操作である。今回は，最も基本的な 1 点交差（one-point crossover）と呼ばれる交差を用いた。図 4 に 1 点交差の例を示す。

解析モデルおよび解析方法

本研究で使用する解析モデルは図 6 のようになる。

図6

解析モデル

この解析モデルの運動方程式は次式で表される。

図4

一点交差例

m1 x1＋c（x＋c（x＋1 1－x0）
2 1－x2）3k（x＋k（x＋ε（x1－x0）＝01 1－x0）
2 1－x2）m2 x2＋c（x＋k（x＝02 2－x1）
2 2－x1）
●●

●

●

●●

●

●

●

●

（1）

平成22年２月

― 51 ―
遺伝的アルゴリズムを用いた非線形振動系の解析

ここで ki xi ci ε（i ＝1：主振動系，i ＝2：動吸振器）
はそれぞれ，ばね定数，質量の鉛直方向の変位，粘性減衰係数，非線形ばね定数とし，x0は主振動系の上端に作用する強制変位である。本研究で使用するばねはここで以下のようなパラメータを定義する。
x0＝ a sin ωt
x0＝ ωa cos ωtωt ＝ T
●

m1
μ＝m2
ε a2
E＝m1 ω2x1
ax2
X2＝ac1
C1＝m1ωc2
C2＝m2ω
X1＝［A1 cosT＋（B1－1）
sinT］3（7）3A32222～［A＋（B－1）］
cosT＋（B－1）［A＋（B－1）］
sinT
111
1＝ 1 144（1）’3
X1＋X1＋X1－K2 X2＋E
（C1＋ C2）（K1＋K2）（X1－sinT）＝ C1 cosT＋K1 sinT
（2）
C2
C2
K2
K2
X2－ X1＋ X2－ X1＋ X2＝0μμμμ●●●●●という式が得られる。さらに，X1, X2を（3）式のように仮定する。
X1 ＝ A1 cosT ＋ B1 sinT
X2 ＝ A2 cosT ＋ B2 sinT（3）仮定した（3）式を（2）式に代入し，変形すると，以下の（4）式が得られる。
－A1 cosT－B1 sinT＋（C1＋C（－A2）
1 sinT－B1 cosT）－C（＋（K1＋K（A2 －A2 sinT＋B2 cosT）
2） 1 cosT＋B1 sinT）
3－K（A＋E［A1 cosT＋（B1－1）
sinT］
22 cosT＋B2 sinT）＝C1 cosT＋K1 sinTC2－A2 cosT－B2 sinT－ （－A1 sinT－B1 cosT） （4）μC2K2＋ （－A2 sinT＋B2 cosT）－ （A1 cosT＋B1 sinT）
μμ
K2＋ （A2 cosT＋B2 sinT）＝0μ
ここで，（4）式に sin3T と cos3T の項が生じる。これらは数学公式（2）より，（－sin3T＋3sinT）
sin3T ＝4（cos3T ＋3cosT）
cos3T ＝4と変形し，高周波成分を無視すると，

秋田高専研究紀要第45号（6）（6）式のように近似することが出来る。従って，先の項は次のように表すことができる。

k1
m1ω2
k2
K2＝m2ω2
K1＝

これを（1）式に代入し，変形すると式のようになる。
●●3～ sinT
sin3T＝
4
3～ cosT
cos3T＝4（5）

ここで，Δ1＝ A1Δ2＝ B1－1（8）と置き，（4）式を cosT と sinT の項でまとめると，次のようになる。
A1－K2A2＋
B1－C2B2（K1＋K2－1）（C1＋C2）32
2＝C1－ EΔ（Δ11 ＋Δ1 ）4
A1＋C2A2－
B1－K2B2－（C1＋C2）（K1＋K2－1）
322＝K1－ EΔ（Δ21 ＋Δ2 ） （9）4
A2－C2B1＋C2B2＝0－K2 A1＋（K2－μ）
C2A1－C2A2－K2B1＋B2＝0（K2－μ）
そして，以下のように無次元量を定義する。
m1k2
c2, γ＝ , δ＝ ,k1
k1
c12
εac1
ε ＝ , C1＝k1
m1k1

Ω＝ ω（10）（10）式を用いて先のパラメータ C1, C2, K1, K2, E を変形し，（9）式に代入すると，以下の（11）式が得られる。
B1＋C1ΩδA2－γB2－C1Ω（1＋δ）＋（1＋γ－Ω2）
322＝1－ εΔ（Δ21 ＋Δ2 ）4
A1＋C1Ω
B1－γA2－C1δΩB2（1＋γ－Ω2）（1＋δ）
322 （11）＝C1Ω－ εΔ（Δ11 ＋Δ2 ）4C1δΩA1－γB1－C1δΩA2＋B2＝0（γ－μΩ2）
2A2＋C1δΩB2＝0－γA1－C1δΩB1＋（γ－μΩ ）振動モデルの主振動系の振幅：αと動吸振器の振幅：βは次の式で求められる。

― 52 ―
渡部雄二・小林義和

α＝ A12＋ B12
β＝ A22＋ B22（5）つ，つまり 2 の12乗の4096分割し，その中の 1 点をランダムで選んでいる。

A1, A2, B1, B2は（11）式の連立方程式を解くと求めることができる。しかし，（11）式の連立方程式は各無次元量が決まっているとしても，4 つの式に対して解が 6 つあるので容易に解くことができない。このような連立方程式をどのようにして解くかのフローチャートを以下の図 7 に示す。

図9

図7

計算方法フローチャート

図 7 に示すような流れで A1, A2, B1, B2を求める。
次に，この方法を用いて解析を行う GA のプログラ
ムの説明をする。GA の作業手順は以下の図 8 のように表される。

全探索分割例

図 9 の例では，－5～5 の間を4096分割した中の（4092，4090）という点を選んでいる。この点を 2進数に変換すると以下のようになる。
Δ1＝4092→111111111100
Δ2＝4090→111111111010
この 2 進数で表された情報を個体とし，この数を個体数とする。この 2 進数を以下のような方法で10進数に変換して計算を行う。
5－（－5）＝0.002441406
4096
5－（－5）
Δ1→4092→×4092－5＝4.99023
4096
5－（－5）
Δ2→4090→×4090－5＝4.985354096
1 目盛りの大きさ→以上のように，グラフの座標を定義された範囲の値に変換する。今回は個体数を15個，固体の遺伝子型のビット長を12×12＝24ビット，交差率を0.2（20％），突然変異率を通常時で0.05（5％）とした。
次に評価関数について説明する。
本研究では，評価関数に次のような指数関数を用いた。
図8GA フローチャート

図 8 で書かれている個体の初期値は，図 9 の例のような範囲内で定義されているΔ1, Δ2を12ビットず

Z ＝ exp（－0.5Z1）（13）
ただし，Z1＝2.0￨Δ1－Δ1￨＋2.0￨Δ2－Δ2￨指数関数 Z1＝ e－aZ1の重み：a を変化させたグラフは図10のようになる。
平成22年２月

― 53 ―
遺伝的アルゴリズムを用いた非線形振動系の解析

なぜこのような重みを評価関数に与えたかというと，もし重み：a ＝1.0として適応度を計算すると，図10から分かるように，Z1が比較的大きな値のときは，小さい変化では適応度：Z の値がほとんど変わらない。しかし，Z1の値が 5 付近になると，Z の値が急激に変化する。これでは非常に高い精度が出るが，近似解に収束しにくいので，重みを与えてやり，このような急激な変化が起こらないようにした。しかし，重みを小さくしすぎると，急激な変化は起こりにくくなるが，本来は適応度が低い場所でも適応度が高いと判断してしまうなど，正確な解が求まりにくくなってしまう。上述の理由と数回の計算により，バランスをとって，今回の研究では重み：a ＝0.5とした。

質量比：μは動吸振器の質量に関係しており，これが大きくなると動吸振器が振動し易くなり，主振動系の振動エネルギーを吸収し，振幅を減少させる。
また，共振振動数は低振動数側に移動する。しかし，現実的に考えると動吸振器の質量が主系の質量の数倍以上に大きくなるということは設計上難しいと考えられる。また質量比μを大きくしすぎると低い振動数側に共振点がシフトするという性質も強くなる
ので今回の解析では主振動系の質量の 1/2程度，つまりμ＝0.5程度を用いることにした。
②ばね定数比γ（k2/k1）の影響

図12

図104.

適応度例

解析結果・考察

4.1 各無次元量の影響
まず始めに，非線形ばねを付加しない状態で，主振動系，動吸振器から構成される 2 自由度系の振動について考え，各無次元量μ，γ，δが主振動系の振幅にどのような影響を与えるかについて調べた。
①質量比μ（m2/m1）の影響

図11
秋田高専研究紀要第45号

質量比の影響ばね定数比の影響

ばね定数比：γは，各共振点の振幅に関係しており，ばね定数比γが増加すると動吸振器の振動が起こりづらくなり，動吸振器が，主振動系の振動エネルギーを吸収しづらくなる。よって最大振幅が増加する。また，ばね定数比が大きくなるにつれて，共振よりも高い振動数域での応答振幅が小さくなるという性質を持っていることがわかる。しかし，共振点を超えるようなところでの振動の低減はあまり必要ないと思われるのでこの性質は重視しなくても良いと考えられる。以上のことから，設計の際は動吸振器を振動しやすくするように，ある程度の値まで，ばね定数比γは低くしたほうがよいと考えられる。
③減衰係数比δ（c2/c1）の影響

図13

減衰係数比の影響

― 54 ―
渡部雄二・小林義和

動粘性減衰定数比：δはγと同様に動吸振器の振動し易さに関係しており，δが小さくなると動吸振器が振動し易くなり，2 自由度系の振動系の状態に近づく。逆にδが大きくなると動吸振器が振動し難くなり，主振動系の質量と動吸振器の質量が一体となったようになる。よって，1 自由度系の状態に近くなり振幅が増加すると考えられる。最適設計パラメータを設定する際，δは小さくし過ぎず，δ＝0.25～0.5付近が良いのではないかと考えた。
これらの結果は，あくまでも今回解析したパラメータでの結果である。種々のパラメータ計算したところ，今回得られた結果と同様の傾向が見られた。
しかしながら非常に極端な組み合わせで計算した場合には異なった結果となる可能性もあると考えられる。
4.2 非線形ばねを用いたときの各無次元量の影響
4.1では，非線形ばねを用いないときの各無次元量の影響を調べたが，ここでは，非線形ばねを用いることで，各無次元量の効果がどのように変化するかを調べた。今回の解析では一例としてε＝0.3の非線形ばねを用い，無次元量γを変化させた。

図14

非線形あり（右）なし（左）によるδの影響

非線形ばねを用いない場合（図14a））と違い，図14b）では，最大振幅が右上に引っ張られるという非線形ばねの効果が出ているということが分かる。
さらに，ばね定数比：γが大きくなるにつれて振幅の増加率も大きくなっている。またε＝0.3の場合とε＝0 の場合を比較すると，γが振幅に及ぼす効果は絶対的な大きさは異なっているがその傾向は似通っていることが分かる。
μ，δなどを変化させても同様の結果を得られた
ことから，今回設定したε＝0.3程度の非線形ばねを付加しても無次元数の影響にそれほど大きな変化はないことがわかった。

4.3 解析方法の精度
次に非線形ばねを付加した場合の全探索と GA を用いた場合の精度について考察する。ε＝0.5の非線形ばねを付加した状態でμ, γ, δ, C1, εなどの無次元パラメータを変化させて多数の計算を行った。その一例としてμ＝0, γ＝0.13, δ＝0.37, C1＝0.56568, ε＝0.5のときの一自由度非線形振動系の結果を図15（a）に示し，またμ＝0.5として他のパラメータは同様の場合の動吸振器を付加した 2 自由度非線形振動系の結果を図15（b）に示す。

図15a）一自由度系応答曲線b）二自由度系応答曲線

結果として今回のモデルの遺伝的アルゴリズムでの解析は全探索の解析と比べ解析時間をかなり減少させることはできた。そして，一点ずつ解析していった場合と比較しても，遺伝的アルゴリズムを用いた解析に振幅の傾向を見る程度の精度はあると図19から見てとることができる。しかし図15a），特にΩ＝2.1～3.0の領域を見ると分かるように，遺伝的アルゴリズムを用いた場合，複数の解があるとき，確実に 2 つないし 3 つの解を同時に探索し，解を求めることは現時点では出来なかった。対応策としては評
価生物数の増加や突然変異率と交差率，交差法の変化などが考えられるが，どちらの対応策も解析時間の増加につながる可能性が高い。また図15の矢印の振動数では大きく値が離れていることが分かる。
これは評価関数において，誤差の範囲内にΔ1, Δ2が入っているからと考えられる。この問題は評価関数の精度を厳しくする，またはΔ1, Δ2の bit 数を大きくすることで改善できると思われる。しかし，図17（a），（b）どちらを見ても全探索の応答曲線を見ると非常に近い点を通っていると推測できる。よってこのような場合でも実用上の精度は問題ないと考えた。結論としては今回のプログラムでは一定の振動数に複数の振動点がある場合の解析は出来なかったが，解析時間の短縮は可能で，且つ精度に問題はな
平成22年２月

― 55 ―
遺伝的アルゴリズムを用いた非線形振動系の解析

いといえるので，まだ問題はあるが改良を重ねることで短時間での非線形振動系の解析も可能であるといえる。さらに非線形ばねを付加した場合の動吸振器の影響という点で図15を見ると非線形ばねの影響と考えられる共振点が右上に引き伸ばされる現象が，動吸振器を付加することでその影響がなくなっている。このことから動吸振器と非線形ばねは組み合わせるとお互いの影響を打ち消しあう性質をもっていると考えられる。
5.

結言

今回の研究から得られた結論は以下のようになる。
①質量比μ，ばね定数比γ，減衰係数比δが主振動系と動吸振器の振幅比に与える影響は非線形ばねを付加した場合と線形ばねを付加した場合との比較からその傾向に大きな変化は見られないことがわかった。
②非線形ばねと動吸振器は両方の特徴を打ち消しあう特性を持っており相性は悪い。

秋田高専研究紀要第45号③ GA を用いた非線形の解析は精度を考えると実用上の問題はないといえる。更に解析時間の短縮にもなった。しかし，複数の解がある場合に今回のプログラムでは複数同時の解析はできなかった。
今回のプログラムで複数の解析ができなかった原
因としてはプログラムの条件をある振動数において一つの振動点を解析した後，次の振動数に移行するというものにしたからである。複数の振動数を解析するように改良するとなると複数の振動数があるかの判断条件，次の振動数への移行条件などが問題となってくる。これを解決するようにプログラムを改良するとなると，明らかに解析時間の増加につながると予想できる。この問題を解決することが今後の課題となる。
参考文献
1）安居院猛・長尾智晴 共著「ジェネティックアルゴリズム」昭晃堂 pp.1-16
2）森口繁一・宇田川銈久・一松 信 著「岩波 数学公式」岩波書店 pp.190