最新版を正確に読む場合には pdf ファイル をどうぞ。これは web 検索のための簡易旧版です。

• 10.5.1 振り子の運動
  • 10.5.1.1 非線形の運動方程式
  • 10.5.1.2 定性的な検討
  • 10.5.1.3 摂動法の考え方--代数方程式を例にして
    • 10.5.1.3.1 基本的な摂動法:
    • 10.5.1.3.2 べき乗選択の難しさ:
    • 10.5.1.3.3 スケールの調整--特異摂動法:
  • 10.5.1.4 摂動法による振り子の運動の検討
    • 10.5.1.4.1 素朴なアプローチ:
    • 10.5.1.4.2 特異摂動法:
  • 10.5.1.5 静的安定
  
  
• 10.5.2 リミットサイクル
• 10.5.3 係数励振
  • 10.5.3.1 変調振り子
  • 10.5.3.2 摂動法--素朴なアプローチ
  • 10.5.3.3 特異摂動法
  
  
• 10.5.4 Duffing振動子の強制振動
  • 10.5.4.1 非線形弾性体の振動
  • 10.5.4.2 スケールの変更
  • 10.5.4.3 摂動法は難しい

10.5 1自由度系の非線形振動

10.5.1 振り子の運動

10.5.1.1 非線形の運動方程式

振り子の運動方程式(10.6)を見れば明らかなように， 回転角$\theta(t)$が小さいものとしてsine関数を近似していた。 このように，現実問題は非線形であるにもかかわらず，それが線形近似の範囲で 取り扱うことができるものとして線形理論の範囲で構造力学や材料力学を 使うことの方が多い。非線形の理屈が重要になるものの代表は， 例えば章-6の座屈理論で考慮する必要がある 幾何学的非線形性と，章-9の塑性論で考慮すべき 材料学的非線形性であろう。 しかしながら，この二つの非線形性をそのまま振動問題に 組み込むことは非常に難しい。 ここでは，そういった非線形性が示すであろう振動特性について， 電気回路等で知られているような問題も含めて考察するために， 主に古典的な振動問題の代表的な非線形問題を取り上げ， その解析法としての 摂動法を紹介したい。なおこの節は，Northwestern大学Davis教授(1980年頃当時） の `Asymptotic and Perturbation Methods in Applied Mathematics' の 講義ノート^10.12を参考にした。

振り子の運動方程式を式(10.6)のようには近似しない場合には， その正しい運動方程式は

$$\ddot{\theta}+\dfrac{g}{\ell} \sin\theta=0$$ (10.254)

となる。初期条件は例えば

$$\theta(0)=a, \quad \dot{\theta}(0)=b$$ (10.255)

となる。この初期位置$a$や初速$b$が大きいと，運動方程式のsine関数は 近似できなくなる。

10.5.1.2 定性的な検討

図 10.64: 位相平面上の振り子の運動

最初は，固体摩擦振動を 検討した節-10.1.2 (3)で 用いた位相平面 上で，定性的な特性を考えてみよう。まず式(10.13)で 定義された$\omega$を用いて， 運動方程式(10.255)に$\dot{\theta}$を乗じると

$$\dot\theta \ddot\theta+\omega^2 \dot\theta \sin\theta=0 ... ...right)^2\right\}}+ \omega^2 \dot{\left(-\cos\theta\right)}=0$$

となるので，積分すれば位相平面上の運動が

$$\dfrac12\left(\dot\theta\right)^2-\omega^2 \cos\theta=C (\mbox{const.})$$ (10.256)

と表されることになる。 初期条件を代入して $C=\dfrac{b^2}{2}-\omega^2 \cos a$である。 この式(10.257)を$m\ell^2$倍すれば，それは高校の物理で 使ったことのある式になる。つまり，第1項は $\dfrac12 mv^2(t)$, ( $v(t)\equiv \ell\dot\theta(t)=\mbox{回転方向の速度}$)の運動エネルギであり， 第2項は $mg\ell\cos\theta$という 回転中心レベルを基準にした位置エネルギであるから，エネルギ保存則である。 種々の初期条件で描いたのが図-10.64である。 例えば初速が零で初期位置が$\pi$より小さい場合には，原点を中心とした 時計回りの軌跡 を描く。 初期位置が十分小さい場合には，線形理論の軌跡の式(10.28)で 示したように，ほぼ楕円状の軌跡を描く。 もちろん$\theta$が$\pi$の偶数倍の位置を 中心とする楕円軌跡も存在する。いずれも左右に振れる振り子運動である。 また後述の図-10.65にあるように， 初期位置$\pi$から始まる振動の軌跡は， 横軸上で滑らかな軌跡にならないのは興味深い。 さらに初速を大きくした場合には， 楕円のような軌跡の上下に描いたように，一方向に回転し続けて 振り子運動にはならないことを示している。

10.5.1.3 摂動法の考え方--代数方程式を例にして

10.5.1.3.1 基本的な摂動法:

次の節では，振り子運動の定量的な把握のために摂動法 というものを用いるのであるが，まず，摂動法そのものの 基本的なところを説明しよう。 そのために，摂動法を代数方程式に用いた例を示す。 以下， $0\leq \epsilon \ll 1$とする。最初の例では

$$x^2-2\epsilon x+1=0$$

を解く。正解はもちろん $x=\epsilon\pm\mbox{i}\sqrt{1-\epsilon^2}$である。 ここでは，解を$\epsilon$で漸近展開で近似し

$$x=x_0+\epsilon x_1+\epsilon^2 x_2+O(\epsilon^3)$$

と置くことにする。 記号$O$は，p.の脚注で定義したLandauの記号である。 これを上の代数方程式に代入して$\epsilon$のべき順に整理すると

$$\left(x_0^2+1\right) +\epsilon \left(2x_0 x_1-2x_0\right) +\epsilon^2 \left(x_1^2+2x_0 x_2-2x_1\right)+O(\epsilon^3)=0$$

となる。したがって，各$\epsilon$の係数を零とすることによって

$$x_0=\pm\mbox{i}, \quad x_1=1, \quad x_2=\mp \dfrac12 \mbox{i}$$

と求められることから，解は

$$x=\pm\mbox{i}+\epsilon\mp\dfrac12\epsilon^2+O(\epsilon^3)$$

である。上の正解の平方根をTaylor展開すれば明らかなように，ここまでは 正解と一致している。このように，ある鍵となるパラメータで 漸近展開して解を求める方法が摂動法である。

10.5.1.3.2 べき乗選択の難しさ:

次の例として，同じような代数方程式

$$\left(x-1\right)^2+\epsilon=0$$

を対象としてみよう。上の例と同様の漸近展開すると

$$O(1): \left(x_0-1\right)^2=0 \quad \to \quad x_0=1$$

となるが，次のべき乗は

$$O(\epsilon): 2\left(x_0-1\right) x_1+1=0 \quad \to \quad 1=0 \quad \to \quad \mbox{不能}$$

となり失敗する。しかし，これは当然だろう。というのも， 正解は $1\pm\mbox{i}\sqrt{\epsilon}$だからである。 正解を知っているなら，漸近展開は $\sqrt{\epsilon}$のべき乗で 与えなければならないことは明らかである。つまり

$$x=x_0+\sqrt{\epsilon} x_1+\epsilon x_2+O(\epsilon^{\slfrac{\mbox{\tiny$3$}}{\mbox{\tiny$2$}}})$$

と置いて，上の方程式に代入すると，$O(1)$は上と同じなので$x_0=1$であり， それ以外が

$$\begin{eqnarray*} O(\sqrt{\epsilon})&:& 2\left(x_0-1\right) x_1=0 \quad \to \q... ...x_2+1=0 \quad \to \quad x_1=\pm\mbox{i}, \quad x_2\mbox{は不定} \end{eqnarray*}$$

となった上で，高次の項をさらに演算すると，$x_k=0$, ($k\ge 2$)を得る。 これは正解である。 このように，どのような漸近展開をするか等の工夫が必要ではあるものの， べき乗の係数を求めることによって正解を得る可能性があることがわかる。

10.5.1.3.3 スケールの調整--特異摂動法:

最後に，最も高次な項に微小パラメータ$\epsilon$がある場合

$$\epsilon x^2-2x+1=0$$

を解いておこう。この問題では，$\epsilon$を無視することによって， この問題の特性を支配している最も重要な項が消えてしまうことが問題になる。 正解は

$$x=\dfrac{1}{\epsilon} \left(1\pm\sqrt{1-\epsilon}\right) =... ...dfrac12 \epsilon -\dfrac18 \epsilon^2+\cdots\right\}\right]$$

なので，二つの解は

$$x_+=\dfrac{2}{\epsilon}-\dfrac12-\dfrac18 \epsilon+\cdots, \quad x_-=\dfrac12+\dfrac18 \epsilon+\cdots$$

である。したがって，最初の例と同じように

$$x=x_0+\epsilon x_1+\epsilon^2 x_2+O(\epsilon^3)$$

と置いてしまうと$x_+$が求められないのは当然である。一応演算してみると

$$\left(-2x_0+1\right)+\epsilon \left(x_0^2-2x_1\right)+\cdots=0$$

となるので，$x_0=\dfrac12$, $x_1=\dfrac18$と求められる。 これは$x_-$に他ならないが，$x_+$は求められない。$\epsilon x^2$が 最も重要な項であるにもかかわらず，このような単純な摂動では， それ，つまり最も重要な項が 他の項に比べて無視できるほど小さいとしてしまったことになり， 失敗したと考えられる。

そこでスケールの変更が必要になる。 例えば，方程式の第1項`$\epsilon x^2$'と第3項`1'が 同じくらいの重要なスケールになって，第2項`$x$'は 無視できるようになると仮定すると，$\epsilon\to 0$のとき

$$\epsilon x^2\sim 1 \quad\to\quad \left\vert x\right\vert \sim \dfrac{1}{\sqrt{\epsilon}}\gg 1$$

となって，$x$が他の2項に比べて無視できるとした仮定と矛盾する。 そこで，方程式の第1項`$\epsilon x^2$'と第2項`$x$'が 同じくらいの重要なスケールになると仮定すると

$$\epsilon x^2\sim x \quad\to\quad \left\vert x\right\vert \sim \dfrac{1}{\epsilon}\gg 1$$

となるので，第3項の`1'は他の2項に比べて無視できることが確認でき， スケールの仮定に矛盾が無い。したがって， $\epsilon x^2\sim x$， つまり $x\sim \slfrac{1}{\epsilon}$程度と考えればいいことになる。 結局，新しい座標^10.13を

$$X\equiv \epsilon x$$ (10.257)

と定義し直（スケールの変更を）して，元の方程式を書き直すと

$$X^2-2 X+\epsilon=0$$

となる。この解を通常の漸近展開で

$$X=X_0+\epsilon X_1+\cdots$$

とした上で方程式に代入して整理すると，まず$O(1)$から

$$X_0^2-2 X_0=0 \quad\to\quad X_0=0, 2$$

を得，$O(\epsilon)$から

$$2X_0 X_1-2X_1+1=0 \quad\to\quad X_1=\left\{ \begin{array}{l... ...} \\ \slfrac{-1}{2} & X_0=2 \mbox{の場合} \end{array}\right.$$

を得る。これは

$$X=2-\dfrac{\epsilon}{2}+\cdots \quad\to\quad x=\dfrac{2}{\e... ...on}{2}+\cdots \quad\to\quad x=\dfrac12+\cdots \Rightarrow x_-$$

となり，それぞれが二つの正解$x_+$, $x_-$と一致する。 このようなスケールの変更が必要になる摂動法を特異摂動法 と呼んでいる。

10.5.1.4 摂動法による振り子の運動の検討

さて，振り子運動の非線形微分方程式に戻って， その解を摂動法によって定量評価をするための一例として， 初期位置の$a$が非線形性に直結する重要なパラメータであると考えよう。 そこで

$$\tau\equiv \omega t, \quad \Theta(\tau)\equiv\dfrac{\theta}{a}$$ (10.258)

という無次元化をしておく。こうすると，運動方程式と初期条件は

$$\ddot{\Theta}+\dfrac{1}{a} \sin\left(a\Theta\right)=0, \qquad \Theta(0)=1, \quad \dot{\Theta}(0)=B\equiv\dfrac{b}{a\omega}$$ (10.259)

と表すことができる。 ただし，この節での上付きドットは$\tau$による微分を表す。 摂動法を用い易くするために，sine関数をTaylor展開し， また簡単のために初速は$B=0$とし，解くべき式を

$$\ddot{\Theta}+\Theta=\dfrac{a^2}{6} \Theta^3+O(a^4), \qquad \Theta(0)=1, \quad \dot{\Theta}(0)=0$$ (10.260)

としておく。

10.5.1.4.1 素朴なアプローチ:

まず，代数方程式を解いたときと同様に，$a$の漸近展開で解を表しておこう。 ここには$a$の偶数乗しか存在しないから

$$\Theta(\tau)=\Theta_0+a^2 \Theta_1+a^4 \Theta_2+O(a^6)$$

と置いてもいいだろう。これを微分方程式に代入して$a$のべき乗で整理すると

$$\left(\ddot{\Theta}_0+\Theta_0\right) +a^2\left( \ddot{\Theta}_1+\Theta_1-\dfrac16 \Theta_0^3\right) +O(a^4)=0$$

となり，初期条件は

$$\Theta_0+a^2 \Theta_1+O(a^4)=1, \quad \dot{\Theta}_0+a^2 \dot{\Theta}_1+O(a^4)=0$$

となる。順にべき乗毎に支配方程式を求めると，$O(1)$は

$$\ddot{\Theta}_0+\Theta_0=0, \qquad \Theta_0(0)=1, \quad \dot{\Theta}_0(0)=0$$

となるので，一般解は $\Theta_0=A \sin\tau+B \cos\tau$で， 初期条件から$A=0$, $B=1$になるので， $\Theta_0=\cos\tau$が解である。

では，次のべき乗の$O(a^2)$の支配方程式を求めると

$$\ddot{\Theta}_1+\Theta_1=\dfrac16 \Theta_0^3, \qquad \Theta_1(0)=0, \quad \dot{\Theta}_1(0)=0$$

となる。微分方程式の右辺に，上で求められた$\Theta_0$を代入すると

$$\dfrac16 \Theta_0^3=\dfrac16 \cos^3\tau =\dfrac{1}{24}\left(3\cos\tau+\cos 3\tau\right)$$

なので，一般解は

$$\Theta_1=A_1 \sin\tau+B_1 \cos\tau -\dfrac{1}{192} \cos 3\tau+\dfrac{1}{16} \tau \sin\tau$$

と求められる。これを初期条件に代入すると$A_1=0$, $B_1=\dfrac{1}{192}$になり

$$\Theta_1=\dfrac{1}{192}\left(\cos\tau-\cos 3\tau\right) +\dfrac{1}{16} \tau \sin\tau$$

を得る。結局，ここまでのべき乗項で，解は

$$\Theta=\cos\tau+a^2\left[ \dfrac{1}{192}\left(\cos\tau-\cos 3\tau\right) +\dfrac{1}{16} \tau \sin\tau \right]+O(a^4)$$

と求められたことになる。うまくいっているように見えるが， 実は解の$a^2$の係数の中には$\tau$に比例した 項が存在することから，これは発散解になっていることがわかる。 つまり，どんなに初期位置の$a$が 小さくても，振り子のような振動解ではなく，時間とともに 回転角が大きくなるという解しか求められなかったことになり， これは前節の位相平面上の定性的な検討結果とも矛盾する。 すなわち，このアプローチは失敗なのである。

10.5.1.4.2 特異摂動法:

何がまずかったのか。実は，上のアプローチでは，振動数が$\omega$の 解を求めようとしていた。この振動数が$a$とは無関係の一定値であると していた。もしかしたら，そこがまずかったのではないか（と，賢い人が 考えたのである）。 そこで，$\omega$ではなく，$a$に依存した振動数 $\Omega(a^2) \omega$で 揺れる運動の解を探してみることにする。ここに$\Omega(a^2)$は未知の 関数である。そこで

$$s=\Omega(a^2) \omega t=\Omega(a^2) \tau$$ (10.261)

という独立変数^10.14を導入し，$\tau$ではなく$s$に関する$2\pi$周期の 解を探すことにしよう。この変数変換から

$$\D*{}{\tau}\to\Omega \D*{}{s}, \quad \D*[2]{}{\tau}\to\Omega^2 \D*[2]{}{s}$$

とすればいいので，基礎式(10.261)は

$$\Omega^2(a^2) \Theta''+\Theta=\dfrac16 a^2 \Theta^3+O(a^4), \qquad \Theta(0)=1, \quad \Theta'(0)=0$$

となる。ここのプライムは$s$による微分を表す。ここで，漸近展開を

$$\Theta=\Theta_0(s)+a^2 \Theta_1(s)+O(a^4), \quad \Omega(a^2)=\Omega_0+a^2 \Omega_1+O(a^4)$$ (10.262)

のように置くことにする。$\Omega$も展開していることに注意する。

これから

$$\Omega^2=\Omega_0^2+a^2 2 \Omega_0 \Omega_1+O(a^4)$$

であることに留意して，上の基礎方程式に代入して整理する。 まず$O(1)$の項を集めると

$$\Omega_0^2 \Theta_0''+\Theta_0=0, \qquad \Theta_0(0)=1, \quad \Theta_0'(0)=0$$

と求められる。したがって，一般解を求めて初期条件に代入すると， 結局解は

$$\Theta_0=\cos\left(\dfrac{s}{\Omega_0}\right)$$

と求められる。ここでは周期解，特に$s$については$2\pi$周期の 解を探していることから， $\slfrac{1}{\Omega_0}$は自然数でなければならない。 特に微小振動の解は$\Omega_0=1$の解に相当することから， ここではそれが答だとする。したがって

$$\Theta_0=\cos s, \quad \Omega_0=1$$ (10.263)

という解を得る。

次に，素朴なアプローチでは失敗した$O(a^2)$の項を検討しよう。 このべき乗の項を集めると

$$\Omega_0^2 \Theta_1''+\Theta_1= \dfrac16 \Theta_0^3 -2 \Omega_0 \Omega_1 \Theta_0''$$

という微分方程式を得る。初期条件を陽な形で示すのは難しいの割愛する。 この式の両辺に，$O(1)$の計算で得た結果の式(10.264)を 代入すると，結局

$$\Theta_1''+\Theta_1= \left(\dfrac18+2 \Omega_1\right) \cos s +\dfrac{1}{24} \cos 3s$$

という微分方程式になる。素朴なアプローチで発散解を発生させた 非斉次項は$\cos s$の項である。これが存在する限り，発散解が特解として 存在することになる。つまり，この項が無くなれば周期解を得ることができる。 上式の$\cos s$の項を消すためには

$$\Omega_1=-\dfrac{1}{16}$$

であればいいことになる。そのとき一般解は

$$\Theta_1=A_1 \sin s+B_1 \cos s-\dfrac{1}{192} \cos 3s$$

と求められる。

ここまでの解を集めると

$$\theta(t)=a \cos\left[\left\{1-\dfrac{1}{16} a^2+O(a^4) \right\} \omega t\right]+O(a^3)$$

と表すことができる。これは周期解であるが，振動数も 初期位置$a$に依存した値になっている。この解を，例えば

$$\theta(t)=a \cos \overline{\omega} t +O(a^3)$$

と置くと，それは，微小振動（線形系近似）のときの振動数が若干改善されて

$$\overline{\omega}\equiv \left\{1-\dfrac{1}{16} a^2+O(a^4) \right\} \omega$$ (10.264)

としたことになる。これをPoincaréの周波数シフト と呼ぶ。あるいは

$$\theta(t)=a \cos \omega \overline{t}+O(a^3)$$ (10.265)

と解釈したとすると，それは，時刻という独立変数（座標）を初期位置によって 若干伸び縮みさせたもの，つまりスケール変更したものが

$$\overline{t}\equiv \Omega(a^2) t = \left\{1-\dfrac{1}{16} a^2+O(a^4) \right\} t$$

のように定義されたと考えてもいい。 こういった考え方をmethod of strained coordinates と呼ぶことがある。 この考え方は，ミクロ・マクロのマルチスケール解析法の 一つである均質化法[113] の根拠にもなっている。 この振り子の問題のときは，`method of strained coordinates'と「Poincaréの 周波数シフト」による解が一致するが，一般には必ずしもそうではないらしい。

10.5.1.5 静的安定

図 10.65: 逆立ち振り子の運動

振り子が初速零で初期位置$\theta=\pm\pi$から 運動を始めると，図-10.65の 太線のような位相平面上の軌跡を描く。 ちょうど逆立ちした状態$\theta=\pm\pi$は分岐点のようなもので， その状態で逆転して元の方向に戻る運動の可能性と，そのまま同じ方向に 回転し続けようとする運動の可能性の二つが存在する。 位相平面上の軌跡は，$\theta=\pm\pi$で角が発生しているのは興味深い。 いずれにしても逆立ち状態は不安定である。 しかし，その状態は式(10.255)で加速度が零の場合の式から

$$\omega^2 \sin\theta=0 \quad \to \quad \theta=n \pi, \quad (n=0, \pm1, \cdots)$$

が解であることから，明らかに静的なつり合い状態である。

そこで，静的なつり合い状態が安定か不安定かの判断をするために， 動的な挙動を利用することを考えてみよう。 この手法は，剛体バネモデルの座屈を調べた節-6.2.1 (3)でも用いた。代表的な静的つり合い状態 として $\theta=\theta_0=0$と $\theta=\theta_0=\pi$の二つを対象とする。 このつり合い状態からの微小な乱れを$\Delta\theta$と記すことにすると

$$\theta=\theta_0+\Delta\theta, \quad \sin\theta=\sin\left(\theta_0+\Delta\theta\right)\simeq \Delta\theta \cos\theta_0$$

としていいから， これを運動方程式(10.255)に代入すれば，$\Delta\theta$に ついて線形化できるので

$$\Delta\ddot{\theta} +\omega^2 \cos\theta_0 \Delta\theta=0$$

を得る。したがって，ここで選んだ二つの静的つり合い状態は，それぞれ

$$\begin{eqnarray*} \mbox{正立:}\quad \theta_0=0 &\quad\to\quad& \Delta\ddot{\thet... ...Delta\theta\sim \sinh, \cosh \quad \mbox{（発散解は不安定解）} \end{eqnarray*}$$

というように微小な乱れが求められる。正立の$\theta_0=0$に 対する乱れは周期解であり， 静的つり合い位置を中心とした振り子運動をするので， もし減衰が存在すればそのうち止まることを示している。 つまり安定な状態であると判定できる。 一方，逆立ち状態の解は双曲線関数あるいは指数関数であり，発散解である。 したがって逆立ち状態は不安定であると判定できる。

10.5.2 リミットサイクル

粘性減衰振動モデルで，減衰定数が負の場合には発散する。 風荷重を受ける橋梁等では，このような動的不安定が発生すると考えられている。 似たようなもので，電気回路の振動例に一つ面白いものがある。 それはvan der Pol振動子 と呼ばれるもので，ある$\epsilon>0$に対して

$$\ddot{u}(t)+u(t)=\epsilon \left\{1-u^2(t)\right\} \dot{u}(t), \qquad u(0)=A, \quad \dot{u}(0)=0$$ (10.266)

という支配方程式でモデル化されている。容易に予想できるように，$u\gg1$で あれば正の減衰であり，その逆なら負減衰なので発散しようとする。 だから，振幅が大きくなると減衰が働き，小さい振幅になると 発散しようとする。したがってこのモデルは， 発散解も減衰解も，ある周期解に収束しようとすることが予想される。 その収束する先の周期解をリミットサイクルと呼んでいる。 この問題も摂動法で解いてみるが，リミットサイクルを求める程度なら 素朴なアプローチで十分なことがわかっている。

図 10.66: リミットサイクル

まず解を

$$u(t)=u_0(t)+\epsilon u_1(t)+\epsilon^2 u_2(t)+O(\epsilon^3)$$

と置き，式(10.267)に代入して整理すると

$$\left[\ddot{u}_0+u_0\right] +\epsilon\left[\ddot{u}_1+u_1-\left(1-u_0^2\right) \dot{u}_0\right] +O(\epsilon^2)=0$$

と，初期条件が

$$u_0(0)+\epsilon u_1(0)+O(\epsilon^2)=A, \quad \dot{u}_0(0)+\epsilon \dot{u}_1(0)+O(\epsilon^2)=0$$

になる。これより，まず$O(1)$の項は

$$\ddot{u}_0+u_0=0, \qquad u_0(0)=A, \quad \dot{u}_0(0)=0$$

となるので，解は

$$u_0(t)=A \cos t$$

と求められる。次に$O(\epsilon)$の項は

$$\ddot{u}_1+u_1=\left(1-u_0^2\right) \dot{u}_0, \qquad u_1(0)=0, \quad \dot{u}_1(0)=0$$

となる。これに上の$u_0$を代入すると，$u_1$の特解 $\left(u_1\right)_p$として

$$\left(u_1\right)_p(t)=\dfrac12 A \left(1-\dfrac14 A^2\right) t \cos t -\dfrac{1}{32} A^3 \sin 3t$$

を得る。右辺の第1項は，前節の振り子の問題と同じく発散解である。 今求めようとしているのは，リミットサイクルという周期解なので， この項は零になる必要がある。したがって$A\neq 0$だから

$$1-\dfrac14 A^2 \quad\to\quad A=\pm2$$

でなければならない。つまり，リミットサイクルは初期条件が$A=\pm2$の 振動解になる。したがって一般解は

$$u_1=\mp \dfrac14 \sin 3t+B \sin t+C \cos t$$

と求められ，初期条件から$C=0$, $B=\pm \slfrac34$となるので

$$u_1=\pm\dfrac14\left(3\sin t-\sin 3t\right)$$

を得る。以上の摂動解の範囲で

$$u=\pm2 \cos t\pm \dfrac{\epsilon}{4} \left(3\sin t-\sin 3t... ...dfrac{3\epsilon}{4} \left(\cos t-\cos 3t\right)+O(\epsilon^2)$$

と求められる。$t$を媒介変数とすれば位相平面上のリミットサイクルを 描くことができ，図-10.66にはそれを実線で示した。 これは，上で求められた$u$と$\dot{u}$を，$O(\epsilon)$までの 範囲で描いたものである。これに対し破線と一点鎖線は， 減衰定数を

$$\beta:=-\epsilon \left\{1-u^2(t)\right\}$$

と置いた上で，中央差分で数値的に求めたものである。 破線は$A=1$の初期条件から， 一点鎖線は$A=4$の初期条件から始めたものである。 実線は$\epsilon$までの摂動解なので，数値解のたどり着く先とは若干の ずれがある。しかし，発散解も減衰解もある周期解に 収束しようとしていることが明らかである。 数値解析によって，容易に非線形解析もできることも楽しいでしょ。

10.5.3 係数励振

10.5.3.1 変調振り子

図 10.67: 変調振り子

図-10.67の左側には，長さ$\ell$の振り子の 支点が上下に$v_f(t)$で強制的に変位させられ，それに伴って 左右に振動している状態を描いてある。しかし，これはおかしくないか。 本当に左右に振動するだろうか。 手で実験するのはかなり難しいが，乱れが生じないように そぉーっとゆっくり振り子の支点を上下させても， それがその運動の直角方向に振り子運動を始めるとは予想し難い。 右側の図には，同じ振り子が天井の穴を 通して設置された上で，支点を上下に$v_f(t)$で強制的に変位させた ため，振り子そのものの長さが元の長さではなく$\ell-v_f(t)$に なりながら振動している状態を描いた。 このような振り子を変調振り子 と呼んでいる。社会基盤構造の場合， 例えば吊橋や斜張橋のケーブルがその張られている面内で 上下左右に振動をしているとき，場合によってはその面外に振動を 始める現象[104]が，この変調振り子と同じものと考えられている。

まず左側の図のような運動の場合，質点の変位は

$$u=\ell \sin\theta, \quad v=\ell \left(1-\cos\theta\right)+v_f$$

であるから，加速度は

$$\ddot{u}=\ell \ddot\theta \cos\theta-\ell \left(\dot\thet... ...,\left(\dot\theta\right)^2 \cos\theta+\ddot{v}_f \eqno{(a)}$$

となる。一方，糸の張力を$T$とすると，左右上下方向の運動方程式はそれぞれ

$$-T \sin\theta=m \ddot{u}, \quad T \cos\theta-mg=m \ddot{v}$$

となるので，$T$を消去すれば

$$\ddot{u} \cos\theta+\ddot{v} \sin\theta+g \sin\theta=0 \eqno{(b)}$$

と表される。この式($b$)に式($a$)を 代入して整理すると

$$\ell \ddot\theta+\left(g+\ddot{v}_f\right) \sin\theta=0$$

のように運動方程式が求められる。ここで，支点の強制変位はcosine関数で 与えられるとして

$$v_f=-v^* \cos \omega t, \qquad \omega_0\equiv\sqrt{\dfrac{g}{\ell}}$$ (10.267)

と定義すると，上式は

$$\ddot\theta+\left\{\omega_0^2+\dfrac{v^*}{\ell} \omega^2 \cos \omega t\right\} \sin\theta=0$$

となる。もし微小な振動を考えるとすれば，sine関数を近似して

$$\ddot\theta+\left\{\omega_0^2+\dfrac{v^*}{\ell} \omega^2 \cos \omega t\right\} \theta=0$$ (10.268)

を得る。これはMathieu方程式 と呼ばれている。

図-10.67の右側の図に示した振り子の場合は

$$u=\left(\ell-v_f\right) \sin\theta, \quad v=\ell-\left(\ell-v_f\right) \cos\theta$$

となっているので，加速度が

$$\begin{eqnarray*} \ddot{u}&=&\left\{\ddot\theta \left(\ell-v_f\right) -2\dot\t... ...\dot\theta\right)^2 \left(\ell-v_f\right)\right\} \cos\theta \end{eqnarray*}$$

と表される。運動方程式は式($b$)と同じなので， それに代入して整理すると

$$\ddot\theta \left(\ell-v_f\right)-2 \dot\theta \dot{v}_f +g \sin\theta=0$$

となる。ここで $w\equiv \left(\ell-v_f\right) \theta$と定義して 上式に代入すると

$$\ddot{w}+\dfrac{\ddot{v}_f}{\ell-v_f} w +g \sin\left(\dfrac{w}{\ell-v_f}\right)=0$$

と表すことができる。さらに，変位は大きくなく，近似的に

$$\left\vert\theta\right\vert\ll 1, \quad \left\vert\dfrac{v_f}{\ell}\right\vert\ll 1$$

が成立するものとして，前の例と同様に， 式(10.268)で強制変位が与えられているとすると，上式は

$$\ddot{w}+\left(\omega_0^2 +\dfrac{v^*}{\ell} \omega^2 \cos\omega t\right) w=0$$

となる。これも式(10.269)と同じMathieu方程式である。

ここで変数変換をして

$$\tau\equiv\dfrac12 \omega t \quad\to\quad \D*{}{t}=\dfrac1... ...}{\tau}, \quad \D*[2]{}{t}=\dfrac14 \omega^2 \D*[2]{}{\tau}$$

とし， $w(t)\to u(\tau)$と置き換えると，上式は

$$\ddot{u}+\left(\dfrac{4 \omega_0^2}{\omega^2} +\dfrac{v^*}{\ell} \cos 2\tau\right) u=0$$

と変形できる。ここの上付きドットは$\tau$による微分である。そこで

$$\delta\equiv\dfrac{4 \omega_0^2}{\omega^2}, \quad \epsilon\equiv\dfrac{v^*}{\ell}$$ (10.269)

と定義すると，最終的にMathieu方程式は

$$\ddot{u}(\tau)+\left\{\delta+\epsilon \cos 2\tau\right\} u(\tau)=0$$ (10.270)

と表すことができる。

10.5.3.2 摂動法--素朴なアプローチ

最初は単純な摂動法を用いてみる。

$$u=u_0+\epsilon u_1+O(\epsilon^2)$$ (10.271)

と置いて，上式(10.271)に代入して整理すると

$$\left[\ddot{u}_0+\delta u_0\right]+\epsilon \left[ \ddot{u}_1+\delta u_1+u_0 \cos 2\tau \right]+O(\epsilon^2)=0$$

となる。これから$O(1)$のレベルの解は

$$u_0=\sin\sqrt{\delta}\tau, \quad \cos\sqrt{\delta}\tau$$

となる。ただし，$2\pi$周期の解でなければならないので

$$\sqrt{\delta}=n \quad\to\quad \delta=n^2, \quad (n=0, 1, \cdots)$$

と求められる。

次に$O(\epsilon)$のレベルの方程式は

$$\ddot{u}_1+n^2 u_1=-u_0 \cos 2\tau$$

となる。まず$n=0$の場合には

$$u_0=1, \quad \ddot{u}_1=-\cos 2\tau \quad\to\quad u_1 \sim \dfrac14 \cos 2\tau$$

のように周期解が求められるので，特に問題は生じない。しかし$n=1$の場合の 運動方程式は

$$\ddot{u}_1+u_1=\left\{ \begin{array}{l} -\sin\tau \cos 2\... ...=-\dfrac12\left(\cos\tau+\cos 3\tau\right) \end{array}\right.$$

となるので，右辺の$\sin\tau$, $\cos\tau$の項に対応して， 特解には必ず $\tau \sin\tau$, $\tau \cos\tau$の項が存在し， 発散解にしかならない。 つまり，この運動は不安定である可能性があり， 確かに実験では不安定な運動になる。

10.5.3.3 特異摂動法

そこで，特異摂動法を用いて，周期解を探してみる。 そのために，式(10.272)の摂動に加えて，$\delta$が$\epsilon$に 依存するものとして

$$\delta=\delta_0+\epsilon \delta_1+O(\epsilon^2)$$ (10.272)

という漸近展開をしておく。これを方程式(10.271)に代入して 整理すると

$$\left[\ddot{u}_0+\delta_0 u_0\right]+\epsilon \left[ \ddot... ...+\left(\delta_1+\cos 2\tau\right) u_1 \right]+O(\epsilon^3)=0$$

となる。これから$O(1)$のレベルの解は，素朴なアプローチと同じで

$$\delta_0=n^2, \quad u_0 = \sin n\tau, \quad \cos n\tau$$

と求められる。

では$n=0$, $u_0=1$の場合をまず考えてみる。$O(\epsilon)$のレベルの方程式は

$$\ddot{u}_1=-\delta_1-\cos 2\tau$$

となるので，解は

$$u_1=-\delta_1 \left(\dfrac12 \tau^2+c_1 \tau+c_2\right) +\dfrac14 \cos 2\tau+c_0$$

と求められる。周期解になるためには

$$\delta_1=0$$

でなければならないので，解は

$$u_1=\dfrac14 \cos 2\tau+c_0$$

と求められる。次に$O(\epsilon^2)$のレベルの運動方程式は

$$\ddot{u}_2=-\delta_2-c_0 \cos 2\tau -\dfrac14 \cos^2 2\tau =-\delta_2-\dfrac18 \left(1+\cos 4\tau\right)-c_0 \cos 2\tau$$

となるので，再度発散解にならない条件から

$$\delta_2=-\dfrac18$$

と求められる。ここまでの漸近展開では，$\delta$が

$$\delta=-\dfrac18 \epsilon^2+O(\epsilon^3) \eqno{(c)}$$

でないと，周期解が存在しないことがわかる。

次に$n=1$の場合を考えよう。$u_0=\sin\tau$, $\delta_0=1$の 場合の$O(\epsilon)$のレベルの運動方程式は

$$\ddot{u}_1+u_1=-\left(\delta_1+\cos 2\tau\right) \sin\tau =\left(\dfrac12-\delta_1\right) \sin\tau-\dfrac12 \sin 3\tau$$

となるので，周期解を得るには

$$\delta_1=\dfrac12,\quad u_1=\dfrac{1}{16} \sin 3\tau$$

でいいこと（斉次解はあっても無くても一緒）になる。 次に$O(\epsilon^2)$のレベルが

$$\ddot{u}_2+u_2=-\delta_2 \sin\tau -\dfrac{1}{16} \left(\d... ...t) \sin\tau -\dfrac{1}{32}\left(\sin 3\tau-\sin 5\tau\right)$$

となる。したがって，周期解を得るには $\delta_2=\slfrac{-1}{32}$で あればいいから，結局，ここまでの漸近展開では，$\delta$が

$$\delta=1+\dfrac12 \epsilon-\dfrac{1}{32} \epsilon^2+O(\epsilon^3) \eqno{(d)}$$

でないと，周期解が存在しないことがわかる。

同様の演算を$n=1$で$u_0=\cos\tau$, $\delta_0=1$の 場合で実行する。$O(\epsilon)$のレベルの運動方程式は

$$\ddot{u}_1+u_1=-\left(\delta_1+\cos 2\tau\right) \cos\tau =-\left(\dfrac12+\delta_1\right) \cos\tau-\dfrac12 \cos 3\tau$$

となるので，周期解を得るには

$$\delta_1=-\dfrac12,\quad u_1=\dfrac{1}{16} \cos 3\tau$$

でいいことになる。次に$O(\epsilon^2)$のレベルが

$$\ddot{u}_2+u_2=-\delta_2 \cos\tau -\dfrac{1}{16} \left(\c... ...) \cos\tau +\dfrac{1}{32}\left(4\cos 3\tau-\cos 5\tau\right)$$

となる。したがって，周期解を得るには $\delta_2=\slfrac{-1}{32}$で あればいいから，結局，ここまでの漸近展開では，$\delta$が

$$\delta=1-\dfrac12 \epsilon-\dfrac{1}{32} \epsilon^2+O(\epsilon^3) \eqno{(e)}$$

でないと，周期解が存在しないことがわかる。

同様の演算を$n=2$, $u_0=\sin 2\tau$, $u_0=\cos 2\tau$に対して 行うと，結局 $u_0=\sin 2\tau$に対しては

$$u_1=\dfrac{1}{24} \sin 4\tau, \qquad \delta_1=0, \quad \de... ... \delta=4-\dfrac{1}{48} \epsilon^2+O(\epsilon^3) \eqno{(f)}$$

と求められる。また $u_0=\sin 2\tau$に対しては

$$u_1=-\dfrac18+\dfrac{1}{24} \cos 4\tau, \qquad \delta_1=0, ... ... \delta=4+\dfrac{5}{48} \epsilon^2+O(\epsilon^3) \eqno{(g)}$$

でないといけないことがわかる。

図 10.68: 係数励振動--変調振り子の動的不安定

以上の式($c$)〜($g$)の関係式を 描いたのが図-10.68の左の図である。 この線上であれば解は振動解であることがわかった。 そこで，$\epsilon=0.1$として，式($d$) ($e$)で 挟まれたくさびの中に相当する$\delta=1$の場合と， くさびの外に相当する$\delta=0.8$の場合について， 中央差分法を用いてMathieu方程式を解いてみた。 その解が図-10.68の右の図である。 くさびの外の設定の解（細い実線）は振動解を示していて， 安定な解であることがわかる。くさびの中の設定（太い実線）の場合には， 解はかなり速く振幅が大きくなり，発散解であることを示唆している。 二つの解の周期も異なることに注意する。 このことから，図-10.68の左の図に示したように， 動的応答が安定解である $\epsilon-\delta$の組み合わせと， 不安定な発散解である場合が存在する。 ここでは数値的な解で安定性を確認してみたが， 数学的にはFloquetの定理を用いて，図示したような安定領域を 決定することができる。

Floquet Theorem: 周期$T$の係数を持つ線形常微分方程式の解は

$$u(t)=\exp\left(\gamma t\right) \phi(t)$$

で表される。ここに$\phi$は周期が$T$で，$\gamma$は定数で 実数とは限らない。$\Re(\gamma)$の符号で，解が発散解か 減衰解かの判断ができる。

実際には$\gamma$が零の周辺で摂動を考え

$$\phi=\phi_0+\epsilon \phi_1+O(\epsilon^2), \quad \delta=\de... ...O(\epsilon^2), \quad \gamma=0+\epsilon \gamma_1+O(\epsilon^2)$$

と置いて，Mathieu方程式を解く。例えば$\delta_0=1$の場合， くさびの中の設定では$\gamma$が実数になり発散解であることを意味し， くさびの上では$\gamma_1=0$なので周期解を意味し， くさびの外の設定では$\gamma$は純虚数になって，若干の周波数シフトを 伴った振動解になることを示すことができるそうだ。 このような動的不安定現象を係数励振 と呼んでいるが，実はこれはブランコ がこげる原理なのである。図-10.67の右の 変調振り子は，見方を変えると，振り子の重りの重心が上下運動したときの 運動と同じである。つまり，ブランコに乗った幼稚園児が足を屈伸して その重心を上下させることによって，ブランコが揺れ始めるのと同じである。 幼稚園児が，体でMathieu方程式を解くことができるというのは興味深い。 その園児が大学生になると線形の常微分方程式も解けなくなるのだから。 この現象には，1980年代にNorthwestern大学のTechnological Instituteの オープンキャンパスで，機械工学科の学生が見せてくれた模型で 初めて接した。そのときのオープンキャンパスでは，この係数励振の模型と 箱庭地盤の液状化模型には感激した。二つの類似模型は， 東北大学工学部土木工学科の元技術職員の菅原紘一さんに作ってもらって， 何回かオープンキャンパスでも見せることができた。

10.5.4 Duffing振動子の強制振動

10.5.4.1 非線形弾性体の振動

図 10.69: 非線形弾性体

社会基盤構造の場合，材料の非線形性の代表で重要なものは塑性である。 しかしながら，塑性理論は複雑で，特に除荷と載荷の変形経路が 異なることによって，動的解析に正確に組み込むことが 困難なことが多い。ここでは，もう少し簡単な非線形弾性体のような ものを対象としてみる。図-10.69に示したように， 例えばバネの抵抗力$f$がバネの伸び$u$との間に

$$f=\alpha^2 u+\epsilon \beta u^3, \quad 0< \epsilon \ll 1$$ (10.273)

のような法則が成立する場合を考えよう。$\alpha^2>0$が 線形バネの弾性係数であり，$\beta$が非線形性を代表する弾性係数で， 正の場合は硬化材料であり，負の場合は軟化材料のモデルに相当する。

このような抵抗則の線形バネに$q \cos pt$の外力が作用した場合の， 共振点付近の振動を解析したい。対応する運動方程式は

$$\ddot{u}(t)+\alpha^2 u(t)+\epsilon \beta u^3(t)=q \cos pt$$ (10.274)

と表される。初期条件は零初期条件としておく。 このような振動モデルをDuffingの振動子 と呼んでいる。この問題も摂動法を用いて解いてみよう。

$$u=u_0+\epsilon u_1+O(\epsilon^2)$$

と置き，式(10.275)に代入して整理すると

$$\left[ \ddot{u}_0+\alpha^2 u_0 - q \cos pt \right] +\epsi... ... \ddot{u}_1+\alpha^2 u_1-\beta u_0^3 \right]+O(\epsilon^2)=0$$

となる。$O(1)$のレベルからは，1自由度系の強制振動の解と同じで

$$u_0=\dfrac{q}{\alpha^2-p^2} \cos pt$$

が特解になる。上述のようにここでは，この分母が非常に小さい場合， つまり共振点付近の非線形応答を求めたい。 しかし，このまま続けて$O(\epsilon)$のレベルに代入したとしても， 分母は $\left(\alpha^2-p^2\right)$のべき乗の乗数が大きくなる だけで，さらに問題が重大になるだけである。

10.5.4.2 スケールの変更

そこで，外力レベルは小さいものとすれば，共振点付近の挙動が 求められるかもしれないと考え，当面は

$$q=\epsilon q_0 \eqno{(a)}$$

と置いた上で，時刻を外力周波数で無次元化して

$$\tau\equiv pt \quad\to\quad \D*{}{t}=p \D*{}{\tau}$$

と定義する。以下，上付きドットは$\tau$による微分を表すものとする。 上の運動方程式(10.275)を書き直すと

$$p^2 \ddot{u}+\alpha^2 u +\epsilon \beta u^3=\epsilon q_0 \cos\tau$$

となる。そこで$p$が共振点$\alpha$に近い場合を対象とすることから， 摂動展開を

$$p=\alpha+\epsilon p_1+O(\epsilon^2), \quad u=u_0+\epsilon u_1+O(\epsilon^2)$$

とすると，$p^2$は

$$p^2=\alpha^2+2\alpha p_1 \epsilon+O(\epsilon^2) \quad\to\q... ... 2\alpha p_1\simeq \dfrac{p^2-\alpha^2}{\epsilon} \eqno{(b)}$$

と表される。この展開を運動方程式に代入すると，まず$O(1)$のレベルは

$$\alpha^2 \ddot{u}_0+\alpha^2 u_0=0$$

となるので，一般解は

$$u_0=c_1 \sin\tau+c_2 \cos\tau$$

でいい。次に$O(\epsilon)$のレベルは

$$\alpha^2\left(\ddot{u}_1+u_1\right)= -2\alpha p_1 \ddot{u}_0-\beta u_0^3+q_0 \cos\tau$$

となる。

これに上の$u_0$を代入して整理すると，右辺には 問題を発生させる$\sin\tau$, $\cos\tau$の他に，$\sin 3\tau$, $\cos 3\tau$の 項が存在する。周期解を得るためには，前2者は無くなって欲しいわけである。 そこでその二つの項を集めると

$$\left[2\alpha p_1 c_1-\dfrac34 \beta c_1\left(c_1^2+c_2^... ...c34 \beta c_2\left(c_1^2+c_2^2\right)+q_0\right] \cos\tau$$

と求められる。この各係数がいずれも零になって

$$c_1\left\{2\alpha p_1-\dfrac34 \beta\left(c_1^2+c_2^2\righ... ...a p_1-\dfrac34 \beta\left(c_1^2+c_2^2\right)\right\} +q_0=0$$

が成立すれば周期解を得ることができる。この第1式の括弧が零であれば， 第2式から$q_0=0$という意味の無い解しか求められないので， 結局第1式から$c_1=0$でなければならない。したがって，第2式に代入して

$$2\alpha p_1 c_2-\dfrac34 \beta c_2^3+q_0=0$$

という関係が成立しなければならない。式($a$)の$q_0$の定義と 式($b$)の摂動展開の表現をそれぞれ，この第3, 1項に代入すると結局

$$\dfrac{p^2-\alpha^2}{\epsilon} c_2 -\dfrac34 \beta c_2^3+\dfrac{q}{\epsilon}=0$$

となることから，非線形の材料パラメータ $\left(\epsilon\beta\right)$を含んだ 振幅$c_2$を求める式が

$$\dfrac34 \left(\epsilon\beta\right) c_2^3 +\left(\alpha^2-p^2\right) c_2-q=0$$

のように求め^10.15られる。

図 10.70: Duffing振動子--硬化材料と軟化材料の共振曲線

図-10.70には， 縦軸に $\left\vert c_2\right\vert$つまり $\left\vert u\right\vert$をとり， 非線形系の場合の共振曲線を描いた。パラメータは，$\alpha=1$, $q=3$とした上で， 硬化の場合は $\epsilon\beta=0.003$とした結果であり， 軟化の場合は $\epsilon\beta=-0.0003$とした 結果である。 $p\simeq\alpha=1$付近でも有界の振幅を得る等， 線形系の共振曲線とは異なる特徴が明らかである。 特に，破線で示した振幅の振動は不安定であり，現実には生じ得ない。 したがって，共振点付近での挙動は複雑なものになることが予想される。

図 10.71: Duffing振動子--硬化の場合と軟化の場合の例

そこで，再度中央差分を用いた数値解析アプローチで，時系列挙動を 求めてみた。図-10.71は，$\alpha=1$, $q=3$として， 硬化の場合は$p=1.05$, $\epsilon\beta=0.003$とした結果であり， 軟化の場合は$p=0.95$, $\epsilon\beta=-0.0003$とした結果である。 それぞれの図の下には，位相平面上の挙動も示してある。 時系列応答はそれほど複雑ではないが，位相平面上の軌跡は面白い。 もう少し複雑な挙動を探して，いろいろなパラメータの値を試してみたが， 例えば図-10.72には$\epsilon\beta$をもう少し 大きくした場合の結果を示した。例では$\alpha=1$, $p=1$, $q=30$としている。 この挙動は複雑に見えるが， 時系列応答はさほど変わったものにはなってはいない。 各自，いろいろ試してみるといい。 このようなアプローチは，動的不安定現象にはよく用いられているようだ。 変わったところでは，心理学の分野にも力学モデルとして適用しているような 情報を，インターネットで見たことがある。

図 10.72: Duffing振動子--ちょっと面白い位相平面の軌跡（時系列応答は単純だが）

10.5.4.3 摂動法は難しい

特にDuffingの振動子の場合には，$\epsilon$によるスケールの変更を 随所で用いている。これは摂動法の使い方がとても難しいことを示唆している。 例えば，1自由度系の線形の粘性減衰自由振動の問題に戻ろう。 減衰定数に相当するパラメータを$\epsilon$とすると，運動方程式は

$$\ddot{u}+u=-\epsilon \dot{u}$$ (10.275)

となる。初期条件を

$$u(0)=1, \quad \dot{u}(0)=0$$

として，これを摂動法で解いてみる。もちろん，素朴なアプローチは 破綻するだろうと予想して，method of strained coordinatesを用い

$$t=s+\epsilon s_1(s)+\cdots, \quad u=u_0+\epsilon u_1+\cdots$$

と置いておく。以下の演算を実施すると，すぐ$O(\epsilon)$のレベルで

$$\ddot{u}_1+u_1=\sin s (1-s_1'')+\cos s (-2s_1')$$

を得る。ここにプライムは$s$による微分である。ここで既に右辺の 不要な項を消すような$s_1$が存在しなくなり，破綻してしまっている。

これはなぜか。もちろんすべての読者は既に減衰自由振動の正解を 知っているはずだが，その解は

$$u=\exp\left(-\dfrac{\epsilon}{2} t\right) \sin\left(t\sqrt{1-\slfrac{\epsilon^2}{4}}\right)$$

である。この指数関数のところとsine関数のところの時刻$t$は， それぞれ異なる速さでシフトしないといけない変数になっている。 したがって，1種類のstrained coordinateでは摂動解は求められないことになる。 これを克服するために，「速い時刻」と「遅い時刻」の二つの 時刻スケールを導入して摂動解析する方法が提案されており， それはmethod of multiple scales と呼ばれている。 これ以上は理解していないので，ここいらでやめる。

最新版を正確に読むためには pdf ファイル をどうぞ。これは web 検索のための簡易旧版です。