レイリー・リッツ法

レイリー・リッツ法では、切り捨てが大きくなるにつれて、演算子のスペクトルの実際の値に収束しない値を生成する可能性があります。これらの値はスペクトル汚染として知られています。^{[ 4 ] [ 6 ] [ 9 ]}場合によっては（シュレーディンガー方程式など）、方程式のすべての固有値を含み、かつ汚染を含まない近似は存在しません。^{[ 10 ]}

圧縮（およびそれに伴う汚染）のスペクトルは、演算子の数値範囲によって制限されます。多くの場合、それは本質的な数値範囲として知られる数値範囲の部分集合によって制限されます。^{[ 11 ] [ 12 ]}

マトリックス $${\displaystyle A={\begin{bmatrix}2&0&0\\0&2&1\\0&1&2\end{bmatrix}}}$$ 固有値を持つ${\displaystyle 1,2,3}$および対応する固有ベクトル $${\displaystyle \mathbf {x} _{\lambda =1}={\begin{bmatrix}0\\1\\-1\end{bmatrix}},\quad \mathbf {x} _{\lambda =2}={\begin{bmatrix}1\\0\\0\end{bmatrix}},\quad \mathbf {x} _{\lambda =3}={\begin{bmatrix}0\\1\\1\end{bmatrix}}.}$$ では、 $${\displaystyle V={\begin{bmatrix}0&0\\1&0\\0&1\end{bmatrix}},}$$ それから $${\displaystyle V^{*}AV={\begin{bmatrix}2&1\\1&2\end{bmatrix}}}$$ 固有値を持つ${\displaystyle 1,3}$および対応する固有ベクトル $${\displaystyle \mathbf {y} _{\mu =1}={\begin{bmatrix}1\\-1\end{bmatrix}},\quad \mathbf {y} _{\mu =3}={\begin{bmatrix}1\\1\end{bmatrix}},}$$ リッツ値は${\displaystyle 1,3}$そしてリッツベクトルは $${\displaystyle \mathbf {\tilde {x}} _{{\tilde {\lambda }}=1}={\begin{bmatrix}0\\1\\-1\end{bmatrix}},\quad \mathbf {\tilde {x}} _{{\tilde {\lambda }}=3}={\begin{bmatrix}0\\1\\1\end{bmatrix}}.}$$ リッツベクトルのそれぞれが、の固有ベクトルの1つと正確に一致することが観察される。${\displaystyle A}$与えられた${\displaystyle V}$リッツ値も同様に、3つの固有値のうち正確に2つを与えます。${\displaystyle A}$正確な近似の数学的な説明は、行列の列空間が${\displaystyle V}$偶然にも、それは2つの固有ベクトルによって張られる部分空間と全く同じである。${\displaystyle \mathbf {x} _{\lambda =1}}$そして${\displaystyle \mathbf {x} _{\lambda =3}}$この例では。

特異値の定義${\displaystyle \sigma }$対応する左特異ベクトルと右特異ベクトルは${\displaystyle Mv=\sigma u}$そして${\displaystyle M^{*}u=\sigma v}$エルミート正規行列にレイリー・リッツ法を単純に適用して、近似特異ベクトルと特異値の1組（左右どちらか）を見つけた。${\displaystyle M^{*}M\in \mathbb {C} ^{N\times N}}$または${\displaystyle MM^{*}\in \mathbb {C} ^{M\times M}}$どちらのサイズが小さいかによって、特異値で割るだけで、もう一方の左または右の特異ベクトルのセットを決定できます。${\displaystyle u=Mv/\sigma }$そして${\displaystyle v=M^{*}u/\sigma }$しかし、この除算は、特異値が小さい場合やゼロの場合に不安定になったり、失敗したりします。

別の方法として、例えば正規行列を次のように定義する。${\displaystyle A=M^{*}M\in \mathbb {C} ^{N\times N}}$サイズの${\displaystyle N\times N}$は、特定の${\displaystyle N\times m}$マトリックス${\displaystyle W\in \mathbb {C} ^{N\times m}}$直交列を持つレイリー・リッツ法の固有値問題${\displaystyle m\times m}$マトリックス $${\displaystyle W^{*}AW=W^{*}M^{*}MW=(MW)^{*}MW}$$ これは、${\displaystyle N\times m}$マトリックス${\displaystyle MW}$この解釈により、左および右の近似特異ベクトルを以下のように簡単に同時に計算できます。

1. 計算する${\displaystyle N\times m}$マトリックス${\displaystyle MW}$。
2. 細分化された、または経済的なサイズのSVDを計算する${\displaystyle MW=\mathbf {U} \Sigma \mathbf {V} _{h},}$と${\displaystyle N\times m}$マトリックス${\displaystyle \mathbf {U} }$、${\displaystyle m\times m}$対角行列${\displaystyle \Sigma }$、 そして${\displaystyle m\times m}$マトリックス${\displaystyle \mathbf {V} _{h}}$。
3. リッツ左の行列を計算します${\displaystyle U=\mathbf {U} }$そして、${\displaystyle V_{h}=\mathbf {V} _{h}W^{*}}$特異ベクトル。
4. 出力近似値${\displaystyle U,\Sigma ,V_{h}}$リッツ特異トリプレットと呼ばれる、選択された特異値と、元の行列の対応する左および右の特異ベクトル${\displaystyle M}$行列の列空間に限定された左特異ベクトルを持つ近似的な切断特異値分解（SVD）を表す。${\displaystyle W}$。

このアルゴリズムは、行列を後処理ステップとして使用できます。${\displaystyle W}$これは、 LOBPCGなどの固有値ソルバーの出力であり、正規行列の数値的に選択された固有ベクトルを近似したものです。${\displaystyle A=M^{*}M}$。

量子物理学では、ハミルトニアンのスペクトルは量子力学的システムによって許容される離散的なエネルギー準位の集合であり、レイリー・リッツ法は複雑な原子または核システムのエネルギー状態と波動関数を近似するために使用されます。^{[ 8 ]}実際、単一の水素原子よりも複雑なシステムについては、ハミルトニアンのスペクトルの正確な解は知られていません。^{[ 7 ]}

この場合、試行波動関数、${\displaystyle \Psi }$はシステム上でテストされます。この試行関数は、境界条件（およびその他の物理的制約）を満たすように選択されます。正確な関数は不明です。試行関数には1つ以上の調整可能なパラメータが含まれており、これらのパラメータを変化させて最低エネルギー構成を見つけます。

基底状態エネルギーは、${\displaystyle E_{0}}$は、不等式を満たす。 $${\displaystyle E_{0}\leq {\frac {\langle \Psi |{\hat {H}}|\Psi \rangle }{\langle \Psi |\Psi \rangle }}.}$$

つまり、基底状態のエネルギーはこの値よりも小さい。試行波動関数は常に基底エネルギー以上の期待値を与える。

試行波動関数が基底状態に直交することがわかっている場合、それは何らかの励起状態のエネルギーの境界を提供する。

リッツ仮説関数は、N個の既知の基底関数の線形結合である。${\displaystyle \left\lbrace \Psi _{i}\right\rbrace }$未知の係数によってパラメータ化される： $${\displaystyle \Psi =\sum _{i=1}^{N}c_{i}\Psi _{i}.}$$

既知のハミルトニアンを用いると、その期待値は次のように書ける。 $${\displaystyle \varepsilon ={\frac {\left\langle \displaystyle \sum _{i=1}^{N}c_{i}\Psi _{i}\right|{\hat {H}}\left|\displaystyle \sum _{i=1}^{N}c_{i}\Psi _{i}\right\rangle }{\left\langle \left.\displaystyle \sum _{i=1}^{N}c_{i}\Psi _{i}\right|\displaystyle \sum _{i=1}^{N}c_{i}\Psi _{i}\right\rangle }}={\frac {\displaystyle \sum _{i=1}^{N}\displaystyle \sum _{j=1}^{N}c_{i}^{*}c_{j}H_{ij}}{\displaystyle \sum _{i=1}^{N}\displaystyle \sum _{j=1}^{N}c_{i}^{*}c_{j}S_{ij}}}\equiv {\frac {A}{B}}.}$$

基底関数は通常直交しないため、重なり行列Sには非対角要素がゼロでない。${\displaystyle \left\lbrace c_{i}\right\rbrace }$または${\displaystyle \left\lbrace c_{i}^{*}\right\rbrace }$（最初の式の共役）は期待値を最小化するために使用できます。たとえば、偏微分を行うことで${\displaystyle \varepsilon }$以上${\displaystyle \left\lbrace c_{i}^{*}\right\rbrace }$ゼロの場合、 k = 1, 2, ..., Nのすべてに対して次の等式が得られます。 $${\displaystyle {\frac {\partial \varepsilon }{\partial c_{k}^{*}}}={\frac {\displaystyle \sum _{j=1}^{N}c_{j}(H_{kj}-\varepsilon S_{kj})}{B}}=0,}$$これは、 N個の永年方程式 のセットにつながる。 $${\displaystyle \sum _{j=1}^{N}c_{j}\left(H_{kj}-\varepsilon S_{kj}\right)=0\quad {\text{for}}\quad k=1,2,\dots ,N.}$$

上記の式において、エネルギー${\displaystyle \varepsilon }$そして係数${\displaystyle \left\lbrace c_{j}\right\rbrace }$は未知数である。cに関して、これは同次線形方程式系であり、これらの未知数の係数の行列式がゼロになるときに解が存在する。 $${\displaystyle \det \left(H-\varepsilon S\right)=0,}$$ これは、Nの値に対してのみ真である。${\displaystyle \varepsilon }$さらに、ハミルトニアンはエルミート演算子であるため、H行列もエルミートであり、${\displaystyle \varepsilon _{i}}$現実のものとなるでしょう。${\displaystyle \varepsilon _{i}}$(i=1,2,..,N)${\displaystyle \varepsilon _{0}}$は、使用する基底関数の基底状態に対する最良の近似値となります。残りのN-1 個のエネルギーは、励起状態エネルギーの推定値です。状態iの波動関数の近似値は、係数を求めることで得られます。${\displaystyle \left\lbrace c_{j}\right\rbrace }$対応する永年方程式から。

レイリー・リッツ法は、機械工学において、ばね質量系や断面が変化する軸上のフライホイールなど、多自由度系の近似的な実共振周波数を求めるためによく用いられる。これはレイリー法の拡張である。また、柱の座屈荷重や座屈後の挙動を求めるためにも使用できる。

システムの共振周波数を求めたい場合を考えてみましょう。まず、振動を次の形式で記述します。 $${\displaystyle y(x,t)=Y(x)\cos \omega t}$$ 未知のモード形状${\displaystyle Y(x)}$次に、運動エネルギー項と位置エネルギー項からなる系の全エネルギーを求めます。運動エネルギー項には、${\displaystyle y(x,t)}$そして、${\displaystyle \omega ^{2}}$したがって、システムの全エネルギーを計算し、次の形式で表すことができます。 $${\displaystyle E=T+V\equiv A[Y(x)]\omega ^{2}\sin ^{2}\omega t+B[Y(x)]\cos ^{2}\omega t}$$

エネルギー保存の法則により、平均運動エネルギーは平均位置エネルギーと等しくなければならない。したがって、 $${\displaystyle \omega ^{2}={\frac {B[Y(x)]}{A[Y(x)]}}=R[Y(x)]}$$これはレイリー商 とも呼ばれます。したがって、モード形状がわかっている場合${\displaystyle Y(x)}$計算できる${\displaystyle A[Y(x)]}$そして${\displaystyle B[Y(x)]}$そして、固有振動数が得られます。しかし、モード形状はまだわかりません。これを求めるために、近似することができます。${\displaystyle Y(x)}$いくつかの近似関数の組み合わせとして${\displaystyle Y_{i}(x)}$$${\displaystyle Y(x)=\sum _{i=1}^{N}c_{i}Y_{i}(x)}$$ どこ${\displaystyle c_{1},c_{2},\cdots ,c_{N}}$は決定すべき定数です。一般に、ランダムなセットを選択すると、${\displaystyle c_{1},c_{2},\cdots ,c_{N}}$システムの実際の固有モードの重ね合わせを記述します。しかし、もし私たちが求めるなら${\displaystyle c_{1},c_{2},\cdots ,c_{N}}$固有振動数が${\displaystyle \omega ^{2}}$が最小化されると、このセットで記述されるモードは${\displaystyle c_{1},c_{2},\cdots ,c_{N}}$これは、システムの実際の固有モードの中で最も低い値に近くなります。したがって、これにより最も低い固有振動数が得られます。この近似的に得られた最も低い固有モードに直交する固有モードを見つけることができれば、次のいくつかの固有振動数も近似的に求めることができます。

一般的に、${\displaystyle A[Y(x)]}$そして${\displaystyle B[Y(x)]}$係数に関して二次式で表される項の集合として${\displaystyle c_{i}}$: $${\displaystyle B[Y(x)]=\sum _{i}\sum _{j}c_{i}c_{j}K_{ij}=\mathbf {c} ^{\mathsf {T}}K\mathbf {c} }$$$${\displaystyle A[Y(x)]=\sum _{i}\sum _{j}c_{i}c_{j}M_{ij}=\mathbf {c} ^{\mathsf {T}}M\mathbf {c} }$$ どこ${\displaystyle K}$そして${\displaystyle M}$はそれぞれ、離散システムの剛性行列と質量行列である。

最小化${\displaystyle \omega ^{2}}$になる: $${\displaystyle {\frac {\partial \omega ^{2}}{\partial c_{i}}}={\frac {\partial }{\partial c_{i}}}{\frac {\mathbf {c} ^{\mathsf {T}}K\mathbf {c} }{\mathbf {c} ^{\mathsf {T}}M\mathbf {c} }}=0}$$

これを解決すると、 $${\displaystyle \mathbf {c} ^{\mathsf {T}}M\mathbf {c} {\frac {\partial \mathbf {c} ^{\mathsf {T}}K\mathbf {c} }{\partial \mathbf {c} }}-\mathbf {c} ^{\mathsf {T}}K\mathbf {c} {\frac {\partial \mathbf {c} ^{\mathsf {T}}M\mathbf {c} }{\partial \mathbf {c} }}=0}$$$${\displaystyle K\mathbf {c} -{\frac {\mathbf {c} ^{\mathsf {T}}K\mathbf {c} }{\mathbf {c} ^{\mathsf {T}}M\mathbf {c} }}M\mathbf {c} =\mathbf {0} }$$$${\displaystyle K\mathbf {c} -\omega ^{2}M\mathbf {c} =\mathbf {0} }$$

c の非自明な解を得るためには、c の行列係数の行列式がゼロである必要がある。 $${\displaystyle \det(K-\omega ^{2}M)=0}$$

これにより、システムの最初のN個の固有振動数と固有モードの解が得られます。ここで、N は近似関数の数です。

以下の議論では、システムが2つの集中ばねと2つの集中質量を持ち、2つのモード形状のみを仮定する最も単純なケースを使用します。したがって、M = [ m ₁ , m ₂ ]およびK = [ k ₁ , k ₂ ]となります。

システムにはモード形状が仮定され、2つの項があり、そのうちの1つは係数Bで重み付けされます。例えば、 Y = [ 1, 1] + B [1, −1] です。 単振動理論によれば、変位がゼロのときの速度は角周波数です。      ${\displaystyle \omega }$最大たわみ時のたわみ (y) 倍。この例では、各質量の運動エネルギー(KE) は次のようになります。${\textstyle {\frac {1}{2}}\omega ^{2}Y_{1}^{2}m_{1}}$など、各バネの位置エネルギー（PE）は${\textstyle {\frac {1}{2}}k_{1}Y_{1}^{2}}$等

また、減衰がない場合、最大運動エネルギーは最大位置エネルギーに等しいこともわかっています。したがって、 $${\displaystyle \sum _{i=1}^{2}\left({\frac {1}{2}}\omega ^{2}Y_{i}^{2}M_{i}\right)=\sum _{i=1}^{2}\left({\frac {1}{2}}K_{i}Y_{i}^{2}\right)}$$

モード形状の全体的な振幅は、常に両側から相殺されます。つまり、想定される実際の変位の大きさは重要ではなく、モード形状だけが重要です。

数学的な操作により、次の式が得られます。${\displaystyle \omega }$Bに関して微分して最小値を求める、つまり${\displaystyle d\omega /dB=0}$。これにより、B の値が得られます。${\displaystyle \omega }$が最低値です。これは、の上限解です。${\displaystyle \omega }$もし${\displaystyle \omega }$モード形状が仮定されているため、システムの予測基本周波数であると期待されますが、B は 2 つの仮定されたモード形状関数の最適な「混合」を見つけるために使用されるため、仮定の下でその上限の最小値を見つけました。

この方法には多くのコツがありますが、最も重要なのは、現実的な仮定モード形状を選択することです。たとえば、梁のたわみ問題の場合、予想される解と解析的に類似した変形形状を使用するのが賢明です。変形解の次数が低くても、単純な連結梁の簡単な問題のほとんどには4次式が適合する可能性があります。ばねと質量は離散的である必要はなく、連続的（または混合）でも構いません。分布運動エネルギーと位置エネルギーの項を簡単に記述できる場合、または連続要素を離散的な部分に分割できる場合は、この方法をスプレッドシートで簡単に使用して、かなり複雑な分布システムの固有振動数を求めることができます。

この方法は反復的に使用でき、以前の最良の解にモード形状を追加していくことができます。あるいは、多くの B と多くのモード形状を含む長い式を構築し、それらを部分的に微分することもできます。

クープマン演算子により、有限次元の非線形システムを無限次元の線形システムとして符号化することができます。一般に、これらの問題はどちらも解くのが難しいですが、後者についてはリッツ・ガラーキン法を使用して解を近似することができます。^{[ 14 ]}