非対称ばね-自由度力学系の非線形振動に関する研究

九州大学学術情報リポジトリ

Kyushu University Institutional Repository

非対称ばね-自由度力学系の非線形振動に関する研究

松崎, 健一郎

Graduate School of Engineering, Kyushu University

https://doi.org/10.11501/3081200

出版情報：Kyushu University, 1994, 博士（工学）, 課程博士 バージョン：

権利関係：

非対称ばねー自由度力学系の非線形振動に関する研究

平成7年 1月

松 崎 健 一 郎

目 次

第 1章 序 論 1 

1 .  

1 非対称ばね非線形系に関する研究の歴史と現状 1 

1 .2 

本論文の研究目的と構成

3 

第

2

章 運 動 方 程 式 と 自 由 振 動

7 

2 .   1 

運動方程式とその無次元化

7 

2 .   2 

ばね関数とポテンシャルエネルギー

1 0   2 .   3 

無減衰自由振動

1 0   2 .  4 

二次曲線ばね系自由振動厳密解の数値化

1 3   2 .   4 .   1 

Jacobiのだ円関数cnのq展開に関する基礎的関係式

1 6   2 . 4 . 2  

厳密解の表示

1 7   2 . 4 . 3  

最大変位または最小変位が与えられる問題

1 8   2 . 4 . 4  

振動数が与えられる問題

2 0  

2 . 4 . 5

級数の打ち切り誤差 別

2 .   5 

まとめ

2 1  

第

3

章 数 値 解 析 法

2 3  

3 .   1 

定常解の数値解法

2 3  

3 .   1 .   1 

基礎式

2 3   3 . 1 . 2  

初期条件および周期の推定法

2 4   3 . 1 . 3  

解の安定判別法

2 5  

3. 

1 .  

4 運動方程式の変形 26 

3 .   2 

ポアンカレ写像

3 .   3 

分岐図

3 .   4 

スペクトル解析

3 .   5 

_{リアプノフ指数}

3 .   6 

不変多様体

2 8  

2 9  

3 0  

3 1  

3 2  

3 .  7 

まとめ

3 3  

第4章 二次曲線ばね系の解析結果

3 4  

4 . 1  

振幅特性

3 4  

4 .   2 

安定定常解領域図

4 0  

4 .   3 

分岐現象

4 4  

4 . 4  

振動波形

5 2  

4 .   5 

ポアンカレ写像

6 2  

4 .   6 

スペクトル解析

6 2  

4 .  7 

リアプノフ指数

7 5  

4 .   8 

不変多様体

7 5  

4 .   9 

まとめ

7 5  

第5章 空気ばね系の解析結果

8 1  

5 .   1 

振幅特性

8 1  

5 .   2 

安定定常解領域図

8 6  

5 .   3 

分岐現象

9 0  

5 . 4 

振動波形

9 0  

5 .   5 

ポアンカレ写像

1 0 4   5 .   6 

スペクトル解析

1 0 4   5 .   7 

リアプノフ指数

1 1 5  

5 .   8 

不変多様体

1 1 5  

5 .   9 

まとめ

1 1 5  

第6章 結 論

1 2 0  

参 考 文 献

1 2 2  

謝 辞

1 2 8  

付録A 二次曲線ばね系自由振動厳密解の数値化プログラム

1 2 9  

付録

B

定常振動問題の解と安定半Ij別のための数値解法プログラム

1 4 5  

第1章 序 論

1 .   1 

非対称ばね非線形系に関する研究の歴史と現状

非線形振動系の研究において非線形ばねにより構成される非線形力学系は最 も基本的な問題であり，また非線形振動の研究の多くは非線形ばね系を取り上 げている.非線形ばねの分類の方法としては種々あると考えられるが，対称性 の非線形ばねと非対称性の非線形ばねの二つに分類する方法は最も基本的な分 類の仕方であると考えられる.非線形ばね特性は一般に変位に関するべき級数 で表示でき，とくに線形項に加えて二次，三次の低次の非線形項で表される二 次曲線形および三次曲線形のばねは最も基本的であり，重要であると考えられ る.また，最近では空気ばね，磁気ばね，ゴムばねなどに示される，変位のべ き乗に反比例する項と定数項との差で定義されるばねが注目されてきている.

これらのばねは変位および復原力の両方に飽和を有することが特徴的な非対称 ばねであり，低次のべき級数で近似することは難しい.

さて，非線形ばね系を上述の二つの方法で 分類するとその研究の多くは対称 ばね系に集中しているのが現状である.とくに線形ばねに三次の非線形項を付 加したDuffing型ばね系が多く取り扱われ，自由振動の厳密解から強制振動の カオス様振動までその性質はほとんど解明されているといえよう.自由振動に 関しては，

1 9 1 8

年に， Duffing (1)がWeierstrassのだ円関数によりその自由振動の 厳密解を示したが，その後さらに扱いやすいJacobiのだ円関数による表示(2)，(3) 

も示され，厳密解を具体的に高精度で数値化する方法(4)，(5)， (6)も完成された 強制振動に関しても，その性質は定常振動からカオス振動まで調べられており，

上回(7)，(8)は一自由度のDu伍ng系においてカオス振動の特徴的な図であるスト レンジアトラクタを初めて発見し，さらにそのパラメータの組み合わせにより どのような種類の周期解および非周期解が存在するのかを表す詳細な安定定常 解の地図を作成した.その他にも津田らく9)によるカオス振動の研究などがある.

また， Du伍ng系は多くの非線形振動の入門的な書籍(10)，(11)に摂動法や平均法

‑2‑

などの近似解法やカオス振動の例題として多く取り上げられ，その認知度も高

一方，非対称ばね系はその研究が少なく ，まだ未解決の問題が多く取り残さ れている.本論文で取り上げる非対称ばねの一つである線形項に二次の非線形 項を付加した形で表される二次曲線ばねにより構成される非対称非線形ばね振 動系については現在までに以下に示すような研究がある.

二次曲線ばね系に関する研究は， Helmholtz(2)が耳の鼓膜が周期的に変動す る音圧を受けたときの振動をモデル化したときに取り上げたのが初めての検討 と思われる.Helmholtzは摂動法によって系が二項励振を受けたとき調和振動に 加えて結合共振が発生することを確認し後にこれはMinorsk)川〉によって再現 された.Du伍ng(4)はDu伍ng反復法によって三項近似解を求め， 自由振動に対 しでも Duffing系と同様に Weierstrassのだ円関数によりその自由振動の厳密解 が表示できることを示した.後に，田村ら(5)ぺ16)，(7)は二次曲線ばね一自由度 系の自由振動の厳密解をJacobiのだ円関数で表示し，準線形振動から強非線形 振動まで高精度に数値化する方法を完成している.Pedersen(18)やReuter(l9)は微 小な減衰を有する場合の 112次 分 数 調 波 振 動 に つ い て 解 析 を 行 っ た 植 木 ら (20)， (21)，  (22)は非対称復原力を有する振動系の例として二次曲線ばね系を取り げ，その強制振動を理論的に解析する手法として，逐次近似法や調和ノくランス 法などの理論的な解析手法の検討を行い，それらの妥当性について研究した.

ChauhanとAshwell(23)，(24)はわずかに湾曲した構造物のS形の復原力と変位の関 係を有する単純な構造物の低周波ドラミングや大振幅分数調波応答について，

周期的な外力が与えられる場合と周期的な変位が与えられる場合の両方につし てアナログ，デジタル計算機を用いて解析を行った.また，調和ノくランス法を 用いて同じ系の 112次分数調波振動について解析した(24). Nayfeh(25)はChauhan

らと同じ方程式をマルチプルスケールズ法を用いて解き， Chauhanらの論文の 間違いを指摘した.砂原ら(26)はメルニコフ関数(27)ベ28)，(29)を用いた不変多様体 の理論的解析によりカオス振動の発生する条件を導出し，シミュレーション実

‑3 ‑

験によりこの条件式の妥当性について確認した.Thompsonら(30)，(31)，  (32). (33)は 強制振動系について位相平面上の吸引領域を中心に解析を行った.

次に，本論文で取り上げるもう一つの非対称ばねは空気ばねである.空気ば ねは，最近注目されている磁気浮上物体の振動(34).(35). (36)， (37)のばねと良く似 た特徴を持ち，鉄道車両の台車(38)や自動車のサスペンションとして良く用いら れているが，これについては現在までに以下のような研究が行われている.

小田ら(39)は空気ばね懸架についての問題に対応する力学モデルを示し，これ に基づく動特性を固有振動数の立場から解析を行い，空気ばね設計諸元と動特 性の関係について検討を行った.小柳(40).(41). (42)は小田らが示した振動モデ、ル を用いて，空気ばねによる防振系の最適設計法を示し，さらに非線形性の影響 について調べた.高上ら(43)は空気圧制御による能動的防振装置に関する理論的 および実験的研究を行った.村田ら(44).(45). (46)は空気ばねを持つ強制振動系に カタストロフィ理論を適用してその特性を調べ，さらにその正しさをダイヤフ ラム型空気ばねの実験により確かめた.柏崎ら(47)は空気ばねと積層ゴムとを組 み合わせた免振・除振床システムの研究を行った.石原(48)は空気ばねの非線形 性と車体の弾性を考慮、した鉄道車両の上下・ピッチング振動を解析する手法を 用いて空気ばねの非線形性などの影響を明らかにした.Jettappa(49)は空気ばね の厳密な運動方程式を導出し，片側および両側ばねの両方について考察した.

1 . 2 

本論文の研究目的と構成

前節で示したように，非対称ばね系は対称ばね系に比べてその解析が困難な ためにその研究も少ないが，非対称の特性を有するばねは数多く存在し，対称 なばねでも大振幅の振動においては非対称の特性を示すことが少なくない.こ れらのことを考慮すると非対称ばね系の振動を解析することは大変意義深いと 思われる.

本論文では非対称ばねとして二次曲線ばねと空気ばねを取り上げる.

二次曲線ばね系，空気ばね系共に過去行われた研究は多いとはいえないが，

‑4‑

両方のばねは非対称ばねの典型として重要であると考えられる.二次曲線ばね は，そのモデ、ルの現実性は少ないが，ばね関数は線形項と二次の非線形項のみ の単純な形であり，非対称ばねの最も基本的な形状と考えることができる.ま た，空気ばねは，先に述べたように，最近注目されている磁気反発力と重力に よって構成される磁気ばねと同様に変位のべき乗に反比例する項と定数項との 差として与えられるばねとしては最も単純な形をしており，実際に存在する非 対称ばねの典型的なモデルとしてとらえることができる.これらのことを考慮、

すると二次曲線ばねおよび空気ばねを非対称ばねの典型として取り上げ，これ らのばねによって構成される非対称ばね非線形系の解析を行うことはそれ以外 の非対称ばね系の研究に対しでも大変意義深いと考えられる.

本論文を構成する各章の内容は以下のとおりである.

第2章では，まず二次曲線ばね系と空気ばね系の運動方程式を示しその無次 元化を行う.さらに，ばね関数およびポテンシャルエネルギーの特性と無減裳 自由振動の振動波形および振幅特性について解析する.二次曲線ばね系では非 線形性が強くなるとその振動波形はパルス列的な形状になり，空気ばね系にお いては上に凸な放物線を並べたような波形となり質点が弾性床をはねる様子に 似た衝突系に近い現象を示す.振幅特性の解析により，二つの系とも共通して 軟性ばね系の特徴を有すること，その最大・最小振幅および調和振幅は ω→ O

の場合ωのべき乗に， ω→1の場合 1‑ω のべき乗に反比例または比例するこ とが示される.また，二次曲線ばねを有する一自由度系の自由振動厳密解の計 算法とそのアルゴリズムについて論じている.

第3章では，強制振動の解析に用いる数値解析法について述べる.まず定常 周期解を効率よく求めるための手法として，微分方程式の数値解法と，そのた めの初期値を，解の周期条件が満足されるよう， Newton法を用いた逐次近似に より決定する操作を組み合わせて速やかに定常周期解に到達する方法をまとめ，

その安定判別法について触れる.そして，その方法に当てはまる形に第2章で 示した運動方程式の変形を行う.特に空気ばね系の場合そのばね関数に極が存

‑5 ‑

在し，その近傍では復原力が急増する.これは，数値計算法を適用するにあ たって計算精度の面から非常に問題となるので，その障害を克服するために運 動方程式に変数変換を施す手法を検討した.さらに，主にカオス振動の解析に 用いられるポアンカレ写像，分岐図，スペクトル解析， リアプノフ指数，不変 多様体などの解析手法やその計算方法について述べる.

第4章では， 二次曲線ばね系の強制振動について数値的解析を行う.調和振 動および 1/2次分数調波振動の振幅特性の解析より，特性曲線はすべての場合 でスケルトンカーブに沿って左に傾き軟性ばねの傾向を示し，励振振幅のある 値を境に曲線の形状が閉曲線から閉曲線に変化すること， 1/2次分数調波振動 は周期倍分岐により発生するものと島状に現れるものが存在すること， 定常周 期解が安定・不安定にかかわらず全く発生しない振動数範囲が存在することな どが示される.さらに，パラメータの組み合わせに対して如何なる安定状態か 存在するのかを明らかにするため安定定常解領域図を描き，それが複雑なフラ クタル構造を持つことなどを示す.そして，分岐図により振動数の緩やかな変 化に伴う振動状態、の移り変わりを調べることにより，ジャンプ現象，周期倍分 岐現象，非周期振動の存在などが確認される.また，分岐過程における振動状 態を，振動波形，ポアンカレ写像，スペクトル解析， リアプノフ指数，不変多 様体などによって解析し，非周期振動はカオス振動であると証明される.

第5章では，空気ばね系の強制振動について第4章と同様の解析を行う.振 幅特性は二次曲線ばね系と同様にスケルトンカーブに沿って左に傾き軟性ばね の特徴を持つこと， 112次分数調波振動は調和振動から周期倍分岐して発生し，

113次分数調波振動は独立して島状に現れることなどが示される.安定定常解 領域図は，特に減衰が小さい場合に励振の振幅が大きく励振の振動数が小さい 範囲で複雑になり，その構造には一定のパターンが確認されフラクタルの特徴 を有することが示される.分岐図では，ジャンプ現象，周期倍分岐現象，非周 期振動の存在などに加えて掃引の方向によって現象が異なるヒステリシスが確 認される.また，分岐過程の振動状態の解析により，非周期振動がカオス振動

‑6‑

であることの証明やカオスの危機の前後で振動状態が大きく変化することなど が示される.

第6章では，本論文で得られた結果のまとめを示す.

付録として第2章で取り扱った二次曲線ばねー自由度力学系の自由振動厳密 解の数値計算と第3章で論じた定常周期解の数値解法のサブルーチンプログラ ム(FORTRAN77)を収録している.

ー7‑

第2章 運 動 方 程 式 と 自 由 振 動

2 .   1 

運動方程式とその無次元化

図

2 . 1

に解析に用いたモデルを示す.図に示されたように，ばね関数 j(x)の 非線形ばねと減衰係数cの線形減衰器により支えられた質量mの質点が周期的 外力

F a

cosωtにより励振される一自由度系の強制振動を考える.この系の運動 方程式は次のようになる.

d²x  dx 

nl~_dt +C一+j(X) = 

F a  

coswt 

l.  dt  ^{J "}   V 

(2. 1) 

互 い

r4 t 

一波減界臨びよおω 数動振有国

︐K 数

ηR  

f t

 

ねギ ﹄ ︒

ザG

1 1 l

﹀

￨

￨ j

お ふ ノ パ 問

ザ

﹂ / ハ

J

FL

M・

f? FJ

動

3 4

一 惚 杯 人

川 凡 さ 〆

iL

一 川 訴

一 )

川 一 味

で 式

h w 瓦 ﹁ 仰

強 次 / ヘ 必

} む ま

= o

=

トlckωcc

F l  

x

ト女

( 2 .  2 )  

式

( 2 .1 )

の両辺を kで除し，式

( 2 .2 )

を用いて変形すると次のようになる.

1 d²x  2c  dx. j(x)  九 一^{̲ . 4}

一 一 一 一 一 一…^J…，

ω;d2ωめ dt k  k 

…… 

^(2. 3) 

次に，次式で示される変数変換を行う.

(2. 4) 

Fo coswt 

nonlinear spring  L./ 

f ( x )  

図

2 .1 

解析モデル

‑8 ‑

ここで， αは変位の無次元化定数である.式(2.4)を用いて式(2.3)を変形する と次式を得る.

d2~ ^2c 

d S   a f (

と/α)

a F O 一

ω^一

一 一 一 一

dr;2 ・c dτ k  k ~~~ωo ^w 

最終的に無次元化された運動方程式は次のように示される.

竺

^{+y生+中(と)= P}

^∞

^sη， 

1 

cfcl.  ・ dτ ￨  

(2.5) 

a f (

と

/ α )

2c....  aF;

ω r 

(2. 6)  ただし，中(と)= ，  J =‑， P=」 v=^ー

￨ 

k ・ Cc k ω

。 J

ここで，中(と)，y， Pおよびvはそれぞれ無次元化されたばね関数，減衰係数， 外 力の振￨幅および振動数である.

以下に二次曲線ばね系と空気ばね系の無次元化運動方程式を示す.

(1)二次曲線ばね系

二次曲線ばねのばね関数は次式で表される.

j ( x )  

^=μ

̲Qx

²  L，Q:定数

この場合， α=Q/Lとすると無次元化運動方程式は次のようになる.

d²と

a S 、

~+y ~+ とーと2 ⁼ Pcosvτ，  i 

cfcl.  . dτ ￨  

C  Ur^  1m  ( 

ただし， ^y =一戸= P=三二旦 v=ω ，/~二 ￨ 

• .J 

mL  '  L L  ' 

^V 

L 

^I 

(2)空気ばね系

(2.7) 

( 2 . 8 ) 

図2.2(a)にピストン型空気ばねのモデルを示す.図中のPaおよびpはそれぞ れ大気圧とシリンダの内圧を表す.また， ピストンの断面積を A，平衡状態に おけるシリンダ長を fとする.このときばね関数は次のように表される.

j ( x )  

⁼ (Pa ‑p)A  (2.9)  シリンダ内部の気体の状態変化がポリトロープ変化であり，そのポリトロープ 指数を nとすると，次の関係式が成立する.

P a ( A I Y  

⁼ 

P [  A ( l  

+ 

x ) r   ( 2 .  1 0 ) 

図2.2(b)にタイヤ型空気ばねのモデルを示す.この場合，外側の膨らみ部分の

‑9 ‑

影響により，シリンダ、内部の気体の状態変化が等温変化(n= 1)であるピストン 型空気ばねと等価であると考えると，式

( 2 . 9 )

，

( 2 . 1 0 )

よりタイヤ形空気ばねの ばね関数は次のようになる.

j{x) 

=PaA~

l+X 

この場合， α⁼1/1とすると無次元化運動方程式は次のようになる.

d2S  dと と ^{n ̲ ̲̲} 

‑‑，，‑‑‑ v  ‑‑‑，‑ ̲  .  し ー ‑

dr，2 ・ 'dτ 1 +と ‑‑‑̲..‑‑

c  ~  ，F I ml  ただし y=  P=‑ L， V=ω̲1二二一

， JmPaA/1' ‑ PaA' . W~ PaA 

( 2 .   1 1 ) 

( 2 .   1 2 ) 

なお，これ以降では簡便のために，式

( 2 . 8 )

，

( 2 .   1 2 )

におけると，τ，vを^x，t，ω  にそれぞれ置き換えて以下のように表すことにする.

d²x  dx ^今

一 て 一+γ一 一+x‑x.l= Pcosωf 

dt.l  . dt 

( 2 .   8 a )  

d²x  dx 

一て一+v一 一+← 一 一=Pcosωt 

dt.l  . dt  1 +x 

( 2 .   1 2 a )  

また，無次元化されたばね関数中(と)も以下j{x)に置き換えて表す.

A 

1 

+ 

x 

(a) ピストン型 _(b)タイヤ型

図

2 . 2 

空気ばね

2 .   2 

ばね関数とポテンシャルエネルギー (1)二次曲線ばね

‑10‑

文

1 2 . 3 ( a )

および(b)にばね関数とポテンシャルエネルギーを示す.ばね関数 の図からわかるように xが大きくなる方向では軟性ばねの傾向を，xが小さく なる方向では硬性ばねの傾向を持つ.また，x=oに安定な平衡点を x= 1には 不安定な平衡点を有し，安定と不安定の二つの平衡点が隣接して存在するのが 大きな特徴である.

(2)空気ばね

ばね関数およびポテンシャルエネルギーをそれぞれ図2.4(a)， (b)に示す. 一 次曲線ばねの場合と同様に，xが大きくなる方向では軟性ばねの傾向を，xが 小さくなる方向では硬性ばねの傾向を持つ.また，

x

→ ∞ で は

f ( x )

→1と復 原力の飽和が起こり，

x

→‑1 +0では

f ( x )

→一∞と変位の限界が存在する.

2 .   3 

無減衰自由振動

系の特性を調べる上で自由振動は系の本質的な特性を表し，解析の基礎およ び強制振動の解析を行うときの大きな参考となるものである.本節では上述し た二つの系の無減衰自由振動の振動波形と振幅特性について調べた.二次曲線 ばね系については第1章で述べたように，その無減衰自由振動がJacobiのだ円 関数により記述できることがわかっておりその数値化も実現されているので，

その方法により計算を行った空気ばね系については， RKG法を用いた数値シ ミュレーションによった.

二次曲線ばね系と空気ばね系の無減衰自由振動波形をそれぞれ図

2 . 5

および 図

2 .6

に示す.図

2 .5

の二次曲線ばね系については ω=O. L 0.3. 0.5. 0.7. 0.9. 0.98  の六種を，図

2 .6

の空気ばね系については ω=0.7， 0.8， 0.9， 0.95， 0.98の五種を描 いており，上下の図はそれぞれ変位と速度の波形を表している.

二次曲線ばね系の場合， ω=0.98では波形は正弦波にごく近く，準線形的な 振動となっているのに比べて， ωが小さくなるにつれて非線形性を増し，変位

1 

( H ) L﹁ 仁o h 0 H  斗4

bl} 

.口

ロ

g ‑ 1  

Co 己

、

. 

、、

.  _‑

 

•  

.  _. 

 

.  . 

_•

_. 

 

•

.  .  .  .  . 

 

ウ 山

τ ‑

x 

l M  

C m α a 

D A  

0

D  ・ 路 2 

(a)ばね関数

‑11 ‑

〆‑、、

.~.

;:S  1 

h凶

しむ

ω口

.  

 . 

•

・﹃

n u

•

 

‑同 宮口

ω臼O

円四

.  .  . 

  

_{.  .}  . 

唱 ‑

o 

1   

Displaccmcnt x  2  (b)ポテンシャルエネルギー

図

2 . 3

二次曲線ばねの特性

1 

f同h、

、‑./H 

、』、

~ 0 

H 0 

与4

bl} ロ

. ロ

号

‑1

to 己

o 

¹ 

Displacement x  ²  (a)ばね関数

.>{. 

:5  1 

k (凶 し

H ω

ロω

n u  

一 何 一 山

口 内向︼Oω

噌 ‑

o 

1 

Displaccment x  2  (b)ポテンシャルエネルギー

図

2 . 4

空気ばねの特性

8卜 ω= 0.7  1

1 ^{? ̲ ， . ‑ : t : 〈 y} V ^{i7 ン シ 〈 ¥}

同 同

0.5  ₄ 

。 。

‑0.5 

4 

き

4 S

05t

、 阪弘、

^遣さ

ⁱ o~XミhJJ5J4 JK\ヌミミー司 i 。

一0.5

。

0.5  1  ‑4  ωt/2π 

図

2 . 5 

自由振動波形(二次曲線ばね系〉

。

1  ωt/2π  図

2 .6 

自由振動波形(空気ばね系〉

︼M E

‑13 ‑

の正の方向に膨らんだ波形となり， 最終的にはパノレス列的な強非線形の振動波 斉三になる.

空気ばね系の場合， 二次曲線ばね系と同様に ω^田 1では準線形的な振動とな るが， ωが小さくなると上に凸な放物線を並べたような波形となり，弾性床を はねる質点の運動に似た衝突系に近いものとなっている.

図

2 .7

， 

2 .   8

にそれぞれ二次曲線ばね系と空気ばね系における無減衰自由振動 の振幅特性(スケルトンカーブ)を示す.中央の図は線形のグラフで， 左右両側 の図は横軸にそれぞれωと1‑ω をとり両対数で描いたものである.図中の破 線は最大振幅と最小振幅の大きさを，実線はO次から 7次までの調和振I隔を示 す.

二次曲線ばね系の場合，最大振幅および最小振l隔は ω→0でそれぞれ1.‑0.5  に収束し， ω→1では

J τ

瓦に比例して小さくなる.調和振幅は ω→0でO 次は 1に収束するが， 1次以上は ωに比例する形となり互いに漸近する.

ω→1ではO次は 1‑ω に， 1次以上は (1‑ω

r

^{/2(n:次数〉に比例する.}

空気ばね系の場合， ω→ 0では最大振幅は ω2に反比例する形で大きくなり，

最小振幅は ‑1に漸近する.調和振幅はすべての次数が最大振幅と同様に ω2に 反比例し，それぞれの比は一定となる.これより， ω→ 0では解の波形の形状 はほとんど変化せずその振幅のみ ω2に反比例して大きくなることがわかる.

ω→ lのときはxは微小な値となり

x j ( l + x )

^信 Xーがと近似できることから推察 できるように， 二次曲線ばね系とよく似た傾向を示す.

2 .   4 

二次曲線ばね系自由振動厳密解の数値化(50)

二次曲線ばね系の無減衰自由振動の厳密解はダフィング系(51)，(52)と同様に Jacobiのだ円関数によって表示できることがわかっている.本節ではその数値 化について述べる.

Al  1 

0.1 

0.01 

0.1 

Complcmentary frequcncy  1‑ω 

t

ム‑

'

l  / / 

~m1D__ ̲ ̲ ' " ノ

IH 

︑ ー

︑ ︑ ︑

︑

m¥

¥1  

X¥dA ¥ ¥ 

︑

A7 A6 As， 

0.5  Frequency  ω 

。

0.01 

自由振動の振幅特性(二次曲線ばね系〉 図

2 .7 

1 

0.5 

‑0.5  1 

0.1 

唱EA

n u  n u 

む ℃ 2 一 一 色 ESE

︒

E Z

匂5

25

E8

2門 町

一 句 日 コ

EEEdコE

一持

活

ー15

(川町何脅しな販制﹀

HC oo

‑o 

31 Hh υ5

コσ ω

﹄M

bg cω Eω

E一 仏

ou

H0 .0  

日0

.0

~

CC H 

組幹阿号附将Q掻隠恒担

3 k

n υ

ロωコ

σ ω H h

凶

問.

0

︑

.E

E

・

E

・ ‑ ‑ . . .

 

︐

.E E‑

‑

︐

.E S‑

‑a .  .

 

︐ 

.• ︐

.E

E 

Eg N 

"， 

"， 

"， 

"， 

... 

"， 

...... 

!l _E  ̲ "  

~...-

l>'(  _-~ ̲!:? 

，戸・F 司:

‑、

・岡国，

ー‑

‑ー・"

̲ . ・ ‑ ‑ ‑

ー回・・ ‑

∞ . ∞   図

3h

υ

ロωコσ ω

点

Cコ

~

H. 0 

o 

~

o  ~ ~

2 0.   

Cコ ロp町!]d^山l?:>!UO山 町pueSlU~山泊町ds!p 山nllI!U!山‘山nwp(宮内

Cコ

.‑< 

Cコ

N 

‑16 ‑

2 .   4 .   1 

Jacobiのだ円関数cnのq展開に関する基礎的関係式

Jacobiのだ円関数 cnuのフーリエ展開表示(q展開(53))は次のように表される.

cnu =  π 

京亙 m J

竺一一 x 

q 吋

^I1Jtuj2K(m)] 

K(m)

ぷ 自

1+qn

…

^{2K(m)一K(m)fni}

6 

ここで， Q， q， K(m)は次式で定義される.

Q == K(1‑m)/[ 2K(m)j

π ]  

ー，一 ‑ q三 e‑"" 

K(m)

イ

^{2 dz} 

m 

( =  

0 ~ 1)  :母数ノマラメータ

( 2 .   1 3 ) 

( 2 .   1 4 )  

式

( 2 . 1 3 )

における Q，qは誤差増殖防止のために，m ~

1 /

²では q，nl^注 1/2では Qによる表示を用いる.次に，定義式と必要となる q展開に関する基礎的関係 式について列挙しておく.

μ

三 。

^nl， μ)== 1‑n1 

Ko^三 2K(nl)jπ，K)^三 2K(1‑m)jn;

Q == rrK)/2~。

q == qo = exp( ‑2Q)， q)三 exp(̲n;2 /2Q) 

κ =  [ 1 + 2( qj + qj4 + q: +... +qjnxn +. . .) 

r 

^= 1 + 4 

  ^{; p}

^= Hj2 

I~ イ +qf+q112+… +q?l)+… l4

^1 r ̲ 

r  ( 

^1 + 

b ; )  

^2l 

μi=16q￨/ 、

I

= 160^， 

I 

V 'Vj} 

I 

11 +2(qj +qj4 +q: +

…

+qjnxn +

… )  

^{1 引 1} Kj  ¹ 

、.‑、，・，、... 

、‑...̲ ^V'‑

Hr=1+2pi， 

p

;

 

⁼ pj(1 + p;)，以 =qj+

ι 

時= 2どへ久 =2qFn+1)

( 2 .   1 5 )   ( 2 .   1 6 )   ( 2 .   1 7 )   ( 2 .   1 8 )   ( 2 . 1 9 )   ( 2 . 2 0 )  

( 2 . 2 1 )  

ただし，パま0または 1を表す(以下同様).さらに次の関数を定義しておく.

‑17 ‑

三

: ( q ; )

^三 1‑

収 i J b =‑ J6‑1

ヌ ; ー l l K ; . ) p ) '

x 

[ ( 1

一仏

)(φg+d+

円

) + { 2

p

2 [ 3 + 8 ( 1 +

p')p']p'

+μ~} p ， ]  

ただし，

φ ;

== qj ‑1‑lj 

/ 1 6  

= 8qj ( 

1 中 一 J [ ^{1 ‑2 ( 1} ー り 川 中 』

中^j== [pj ̲(qj2 +

イ ++dn+1 ) +) 戸 ] 角

p:式(2.27)参照

( 2 .   2 2 )  

2 .  4 .   2 

厳密解の表示

無減衰自由振動の運動方程式は，式(2.8a)より次のように表される.

ハU

X +  x  x  一 一

ご'

'‑

d

一︐d

( 2 .   2 3 ) 

式

( 2 .2 3 )

の解を次のように仮定する.

x(t) = D ‑

C 

cn²(Ko Qt/2) 

( 2 .   2 4 ) 

式

( 2 .2 4 )

を式

( 2 .2 3 )

に代入し， cnに関する等べき項の等置から次の関係式を得 る.

C=3

μ 。 / 2 p

^，

¹ 

D 

=    ‑ p ( 1  +  μ 2 ^。 ) / 2 p ^f 

Q = l/Ko

1 P  

ただし

p==~l-(l-μ。 )μ。 =~1-(1-μJμ1 =  ‑ ) 1 ‑ μ

^州

これらから， x(t)の最大値，最小値は，それぞれ次のようになる. max(t) x(t) ⁼ D， 

1 

min(  x)t(t

 

)

=

 

=

^{‑A = D‑C  [} 

また，

A

とDの聞には次の関係式が成立する.

( 2 .  2 5 )   ( 2 .  2 6 )  

( 2 .  2 7 )  

( 2 .   2 8 ) 

( ; ) = j ( ? ) l l + J 2 ( 1 + ( 7 ) )   _{( 2 . 2 9 )}  

ー18‑

Q o c 

dA るd℃創 れ

=

キC 2寸

Il

l1

﹂示

¥パ ソ 表心一2 う

f 吋

よ c n D G  

恥wd

℃ ' 判

ω

斗B6fd￨￨L

I / Q  

4t

fo 

課J 一ぺ/ 刊

hJ

/υ 和 x

最

( 2 .  3 0 ) 

ここで， Gnおよび Aは次式で表される. a 

=J2 g/(1+q~I) ， ^{( q o}

^→

^{0 )}

^， 

1

1 /

coshnQ， (Q→ 0) 

( 2 . 3 1 ) 

い D-3Q2~a~

31d ‑uφo 

+q O [W O 

^{‑Q)(l+呪)‑μ}

。 ( 2+q o ) / 1 6 J  

2(1 +μ

o + p ) p ‑ " "   ( 1 + q o ) 2 ( 1 +

μo + 

p )  

̲3Q2α;( +α;+α;+

一 )

ただし， ^l民主μ

。

(μ

。 + p ) / 2 ( 1 + p )

Q)三 (1‑Q)[2‑(1‑Q)]

( 2 .  3 2 ) 

︑la

E﹃﹄

' at_﹀ ' '' ' a 't l a J

︑ ︑ ︐ ︐ ︐ ︐ ︐

F

︑

J

吋3

唱E E A

〆 n  一 一

' aa︑

¥1111ノn + 

'K  

G n ・KG 

⁝ ℃

2 + 

白

0 4 2  

/1 11 11 1¥  

今4Q 

f o  

‑

‑ 一

n A 

(2. 

3 3 )  

fo

 

A

斗 ・ 今''M

n  一 一

ι E

向

︒ ︒

す

Q  釘

fh u 

‑ ‑ 一

A 

2 . 4 . 3  

最大変位または最小変位が与えられる問題 式

( 2 .2 5 )

， 

( 2 . 2 8 )

を整理すると次の関係式を得る.

ハ ド

0

・

pμ

寸仰o‑

+A 刊 一+μ・一

+

一

吋3一1﹁

︺ 一

t一

+ 帆 一+

フt

11

句︑

J一/

tt

‑

一今ム一今4

d D

 

( 2 .   3 4 ) 

1  A  3

μ ?  

2ι4(1‑μ

) / 2 + p ) p '   1-D=~(μ)

^+ p)μl 

2(1 +μ) 

+ p ) p  

式

( 2

，

3 4 )

，

( 2 . 3 5 )

を μ。または μlについて解くと次の結果を得る.

( 2 . 3 5 ) 

μ

内o‑1 

+α+吋~(1ト一α司州)(泊(1

^+ 3

^{知 刈}

^α

^吋

^)=1川+C+

^{吋 ( れ}

^1ト一2

^叫

^ω

功)布可，

^D

ただし

ι

，a

4

三

; ( ^れ

¹⁺

^ル(杓

^{0‑1)白 A =州叫=ベ=}

^杓

^(00‑

^{1 /}

C三 ;(1‑D)D=(0

判。

D=(O‑l

刈

2

Jb 

2c 

μ^{j ‑}

Jb+J4 三五

‑1+c+(2D‑l)

. J ( . 五万'

ただし

bペ三 ~~Gトド川+叫吋Aサ)G 一 Aサ十←)~ド=ベ(付0-1) 寸 白 j

^{一A =叫山=ベや(い0‑1/}

C

三 ( い

0‑1/3)^日 1‑D = (0 ‑1/2) 

(a)μ。が得られる場合

Oく μo~ 0.5に対しては次の範囲が対応する.

。 くD~ 0.5， 0く A~ 0.366 

‑19 ‑

( 2 .  3 6 )  

(2. 

3 7 )  

式

( 2 .3 6 )

により AまたはDから μ。が求まる.次に，式

( 2 .2 0 )

の逆関数として μ。から q。を，さらに式 (2.31) より q。から α" が求まる .Q については，~

( 2 . 2 6 )

に式

( 2 . 2 2 )

を適用して次式を得る.

1‑Q = 

s ( q o )   ( 2 . 3 8 )  

これから直接 l‑Qを特定することができ，式

( 2 . 3 2 )

，

( 2 . 3 3 )

より

4

が得られ る.

(b)μlが得られる場合

Oく μj_~ 0.5に対しては次の範囲が対応する.

0.5~D く 1 ， 0.366 ~ Aく 1/2

式

( 2 .3 7 )

により

1 /

2‑Aまたは l‑Dから μlが求まる.次に，式

( 2 .2 0 )

の逆関 数として μlから qjを得る.式(2.

1 8 )

を変形することにより得られる関係式

Q=

一 ( ポ / 2 ) / l n q j ( 2 . 3 9 )  

より qjから Qを，式

( 2 .3 1 )

より Qから a_nを得る. Qを式

( 2 .2 6 )

によって求め るときに必要な K。は，式(2.

1 7 )

を変形して得られる次式によって計算される.

‑20‑

Ko = 

π (

/2Q)K₁， 

( 2 . 4 0 )  

ここで，

K )

は式

( 2 .

19) によって求めることができる •

A

は式

( 2 .3 2 )

， 

( 2 .  3 3 )

に よって計算される.

2 . 4 . 4  

振動数が与えられる問題

(a) qo (μ。)が得られる場合

Oく μo^:s; 0.5に対しては次の範囲が対応する.

0.91:s;Q<1 

Qの桁落ちを防止するためにトQが与えられるとすると，式

( 2 .3 8 )

の逆関数 として l‑Qより q。を得る.そして，式

( 2 .3 1 )

より q。から α を，式

( 2 . 2 0 )

よ りμ。，式

( 2 . 3 4 )

より A，D，式

( 2 . 3 2 )

，

( 2 .   3 3 )

より A を得る.

(b) 

Q [ q )  ( μ J ]

^{が得られる場合}

O^く μ1:s; 0.5に対しては次の範囲が対応する.

。

く

^Q:s;0.91 

式

( 2 . 2 2 )

，

( 2 . 2 6 )

， 

( 2 . 4 0 )

より次式を得る.

1‑JtQ/2Q = 

S ( q J

， Q = 

πQ ( / 2 ) / ( 1 ‑ S ( q ) ) )  

式

( 2 . 1 8 )

，

( 2 . 4 1 )

より次式が得られる.

( 2 . 4 1 )  

q) = 

[ 仰 ( ー π / Q ) r ‑

^3(qt)

( 2 . 4 2 )  

上式を試行解法により解くことにより Qより

S ( q ) )

を得，三

( q J

から Q，

q l '

μl  を求めると， an^， An' A， Dなどが得られる.

2 .   4 .   5 

級数の打ち切り誤差

解は最終的にフーリエ級数の形で得られるが，フーリエ級数は有限項で打ち 切らざるを得なし、から，それに伴う誤差を検討しておく必要がある.まず， p，  Opを次のように定義する.

‑ 1 2   ‑

P 三 ~cos;nQ'

p

三 月 三

^cos;

( 2 . 4 3 )  

ー ノ る

﹀ 付 れ 必

っ か 今 ぜ っ か

( 表 (

M

什

ょ

vnノ

す 切 れ 旬 表 幻

y

ザ﹂

tt

今 ︐

h

y v E /

刈

・

1

j o

差

W

U川 町 誤

Q

庁U

7 X 4 1 6

次

A Y

相 仁

﹁ 山 門 川 す 働 ぴ ア 川 (

=

︐ ょ

2.日

21

11

ぉ 料

︑寸山内¥︑

EE tE ノJ

︑ E1 f

叶た2イ 'h v

ト

つ ぱ 一

2

川 / し

切

f

叶

¥ 勧 ー ぴ ちα

差 刊

打

Grl

誤 川 川

! で 引

℃ ん 坊 対 的

項h

い M ‑

絶 べ 限 ず の ぺ 刀

有

μ ω

解 i州

非 乞 一 円 以

LI

除d

問 料

如 刻

︑ 叶 ゾ

' i H V

エ / 川 き 川 リ

% と

&

一 の

￨

￨ フ こ

￨

ご : f j t ) ト l z ( 2 ‑ ? ) l z 2 ?

^ただし

^ろ

_p 

^E

^{く くl}

^{( 2 . 4 6 )}  

p， opはEuler‑Mac1aurinの公式(54)により次のように近似される.

p=1/2ω h  Q + tan ‑1 ( e ^‑Q ) / Q ， 

1 

/  ( Af  ，..¥rl ¥ I 

r  (  2 .   4  7  ) 

op 

=  I j

2cosh(N +2)Q + tan‑

取引は

)Q)/Q

I 

特

マ︿

V

4c  

で

ふ 川

N

ん

﹂

・

々C る れ 司 令

4D  

蓑

マ

o

t

夜

見 小 川 仇

討 課 似 誤 れ 次

h y

‑

‑

斗 ・ふ

切ηパ

l i

l ‑

‑ J

で

ち あ

2 ︐

O E

ま 川 村

↓ の で 仰

D

数 リ

. e L

‑ F︑ rt

﹄F I l l 日 恥 Z

‑ 1 i

白 川

↓

︑ ト +

j D 4

山 川

df

} d 1 人 ︐

e

︐ ・

4

!

川 川 叫

h

小山 川

J勺

U

リ 州 油 相

ω

つ い

4ト 沿

o e O

M

ト 川

φ

↓

↓

広 ん い 寸

z

ヤ Q

4

∞ 小 山

i m y p

J 1 1 '  

α

↓

﹂ { { た

式

Q

ま

( 2 . 4 8 )  

p=ta^{加 紅}n正

巾

‑1ⁱ

⑪ =  tan‑¹(

汁吋 / Q

^，

op/ p 

= 

tan‑¹( e‑{N+2)Q) /ぽ (e‑Q)

( 2 . 4 9 )  

2 .   5 

まとめ

本章では二次曲線ばね系および空気ばね系の運動方程式の導出を行い，両ば ねの特性および無減衰自由振動の振動波形と振幅特性を調べた.非線形性が強 くなると二次曲線ばね系で、はパルス列的な波形になり，空気ばね系では放物線

‑22‑

を並べたような波形になるのが特徴的である.また振幅特性により， 二つの系 とも共通して軟性ばね系の特徴を有すること，すべての特性曲線は ω→0でω の， ω→ lで 1‑ω のべき乗に反比例または比例することが判明した.さらに，

二次曲線ばね系の無減衰自由振動について，その厳密解の Jacobiのだ円関数に よる表示と数値化について述べた.

‑23 ‑

第3章 数値解析法

3 . 1  

定常解の数値解法(55)

非線形系の定常振動を解析的に求めることは非常に困難であり，ほとんど不 可能である場合が多い.そこで，摂動法，調和ノくランス法などに代表される近 似解法や計算機を利用した数値的な解法を用いらざるを得ない.数値的な方法 による場合， Runge‑Kutta‑GilI法などの微分方程式の数値解法を用いるが，これ らの方法は常微分方程式の初期値問題を解く方法であり定常振動解を求めるた めには過渡振動が収まるまで計算を続ける必要があり，減衰がない保存系では 過渡振動は永久に収まらないため定常解の計算は不可能である.そこで，より 効率よく定常解を求める方法が必要となる. AprilleとTrick(56)， (57)は定常振動 解の初期値を変数と考え， Newton法を適用して解く方法を提案した.牛田(58)

は数値微分法に基づいたNewton法と Secant法を用いて定常振動解を求める方 法を提案し，複数の定常振動解が存在するような場合Secant法とパラメータ切 換法(59)をうまく組み合わせることにより複数個の定常振動解を順次連続的に追 跡できることを見いだした小鹿(60)，(61)， (62)は二点境界値問題の数値解を求め るための計算手法として同様の初期値推定法を提案し，さらにこの方法を最も 一般的な線形の二点境界条件に拡張した本節では，これらの方法を整理し，

自律系・非自律系にかかわらず定常振動解を求めることが可能な方法を示し，

さらに解の安定判別法について述べる.

3 .   1 .   1 

基礎式

次の一階連立常微分方程式を取り扱う.

一dx 一

=x

dt 

ただし， t.独立変数(実スカラ)，

X

^三

( X

₁

X 2

^…

X J T  

¹ 

~ : n次元ベクトル X三 (X_J X₂

… x J T f  

(3.  1)