ここでは、これまでつかってきた解析モデルに対して、OpenSeesPyにより固有値解析を実施してみます。
モデルの固有値解析結果を確認することでその動的特性を把握することができます。
動的解析においては、解析モデルの動的な挙動は、固有値解析により得られる動的特性に支配されます。
そのため、それらを把握することは動的解析により得られた結果の妥当性を検証する上で重要となります。
それでは、実際にOpenSeesPyを用いて上記の解析モデルの固有値解析を実行していきます。
OpenSeesPyでは定義したモデルに対して、新たにeiganコマンドを適用することでそのモデルの固有値を簡単に得ることができます。
まず、モデルの節点、質量、部材特性などをコマンドによりあらかじめ設定します。
なお、以後固有値解析結果をより確認しやすくするため、前回までに検討していた解析モデルに対して、その各節点はX方向のみに動くように制約条件を新たに加えて解析を実行するようにします。
#解析モデルの設定
op.wipe()
op.model('basic', '-ndm', 2) #2D3dof
kN_pond=0.22480894387096#kN→kips
g = 386.0885827#acceleration of gravity(in/s2)
# 節点座標の定義
m_inch=39.3700787#単位系m→inchへの変換
op.node(1,0,0*m_inch)#節点1の座標設定(節点番号，x座標，y座標)
op.node(2,0,4*m_inch)#節点2の座標設定
op.node(3,0,8*m_inch)#節点3の座標設定
op.node(4,0,12*m_inch)#節点4の座標設定
# op.node(5,1,0*m_inch)#節点1の隣接節点
#boundary 基礎固定
op.fix(1,1,1,1)
op.fix(2,0,1,1)
op.fix(3,0,1,1)
op.fix(4,0,1,1)
op.mass(2,200*kN_pond/g,200*kN_pond/g,0)
#節点2の節点質量設定（対象節点番号,鉛直方向の質量,水平方向の質量,回転方向の質量)
op.mass(3,200*kN_pond/g,200*kN_pond/g,0)#節点3の節点質量設定
op.mass(4,200*kN_pond/g,200*kN_pond/g,0)#節点4の節点質量設定
# 梁要素の定義
E_mod=205000000*kN_pond/((m_ft)**2);G_mod=78846153.8461538*kN_pond/((m_ft)**2)#鋼材を想定した材料定数の設定
op.geomTransf("Linear",1)#要素座標系から全体座標系への変換手法の定義（線形変換）
op.element('ElasticTimoshenkoBeam', 1, 1, 2, E_mod, G_mod, 0.306305284*(m_ft)**2, 0.145686451*(m_ft)**4,0.5*0.306305284*(m_ft)**2,1)#梁要素1の諸元設定
op.element('ElasticTimoshenkoBeam', 2, 2, 3, E_mod, G_mod, 0.306305284*(m_ft)**2, 0.145686451*(m_ft)**4,0.5*0.306305284*(m_ft)**2,1)#梁要素2の諸元設定
op.element('ElasticTimoshenkoBeam', 3, 3, 4, E_mod, G_mod, 0.306305284*(m_ft)**2, 0.145686451*(m_ft)**4,0.5*0.306305284*(m_ft)**2,1)#梁要素3の諸元設定次に固有値解析の実行コマンドを記述します。
コマンドにおいては、解析モデルが3つの可動節点を有し、各節点はX方向のみの自由度としているので、全部で3つの自由度があり、固有値はモデルの自由度の数だけ得られることが想定されることから、その求めたい固有値の数は3と設定します。
また、解析により得られる固有値を確認するため、print関数で得られた固有値を表示させます。
固有値解析を実施(ソルバー,求めたい固有値の数)
lamda=op.eigen('-fullGenLapack',3)
print(lamda)実行すると、以下のようなリストが表示されます。
これがOpenSeesPyにより得られた解析モデルの各次のモードの固有値λとなります。
[19048.909172331816, 149550.22525504694, 312282.70466847887]
建築構造物の特性をみるうえで重要な固有円振動数ω、固有周期Tは、それぞれ以下の式で得られます。
そこで、以下のようにmap関数を用いて全体の固有円振動数、固有周期を求めてみます。
#固有円振動数を求める。
natural_frequency= list(map(lambda x:x**0.5,lamda))
#固有周期を求める。
natural_period= list(map(lambda x:2*3.14/x,natural_frequency))
#結果を表示する。
print("Natural_frequency(Hz):",natural_frequency)
print("Natural_period(s):",natural_period)すると、以下のように各次の固有円振動数、固有周期が得られます。
Natural_frequency(Hz): [138.01778571014611, 386.7172419934841, 558.8226057242842] Natural_period(s): [0.045501382069617846, 0.01623925524403128, 0.011237913312151281]
OpenSeesPyには、動的特性を包括的に確認、出力するためのコマンドmodalProperties Commandがあります。
このコマンドをeigenコマンド実行後に表記すると、固有値解析の結果をまとめてテキストファイルなどに出力することができます。
そこで、このコマンドを実行することで結果を見てみます。
#固有値解析結果をテキストに出力する。
op.modalProperties('-print', '-file', 'ModalReport.txt', '-unorm')すると、プログラムソースと同じディレクトリ上にModalReport.txtが生成されます。
テキストファイルの中身を確認してみると、以下のように各次モードの動的特性が得られています。
特に固有円振動数（OMEGA）、固有周期（PERIOD）を確認すると、先ほど計算した値と概ね一致した値が得られていることが分かります。
# MODAL ANALYSIS REPORT
* 1. DOMAIN SIZE:
# This is the size of the problem: 2 for 2D problems, 3 for 3D problems.
2
* 2. EIGENVALUE ANALYSIS:
# MODE LAMBDA OMEGA FREQUENCY PERIOD
# ------------- ------------- ------------- ------------- -------------
1 19048.9 138.018 21.9662 0.0455245
2 149550 386.717 61.548 0.0162475
3 312283 558.823 88.9394 0.0112436また同じファイルの中には、モデルの各方向の刺激係数（MODAL PARTICIPATION FACTORS）有効質量（MODAL PARTICIPATION MASSES）、有効質量比（MODAL PARTICIPATION MASS RATIOS (%)）がそれぞれ記載されており、以下の通り確認することができます。
* 6. MODAL PARTICIPATION FACTORS:
# The participation factor for a certain mode 'a' in a certain direction 'i'
# indicates how strongly displacement along (or rotation about)
# the global axes is represented in the eigenvector of that mode.
# MODE MX MY RMZ
# ------------- ------------- ------------- -------------
1 1.22041 0 -431.849
2 -0.349292 0 -44.1118
3 -0.134143 0 11.7234
* 7. MODAL PARTICIPATION MASSES:
# The modal participation masses for each mode.
# MODE MX MY RMZ
# ------------- ------------- ------------- -------------
1 0.319346 0 39986.4
2 0.0261593 0 417.215
3 0.00385821 0 29.4685
* 8. MODAL PARTICIPATION MASSES (cumulative):
# The cumulative modal participation masses for each mode.
# MODE MX MY RMZ
# ------------- ------------- ------------- -------------
1 0.319346 0 39986.4
2 0.345506 0 40403.6
3 0.349364 0 40433.1
* 9. MODAL PARTICIPATION MASS RATIOS (%):
# The modal participation mass ratios (%) for each mode.
# MODE MX MY RMZ
# ------------- ------------- ------------- -------------
1 91.4079 0 98.8953
2 7.4877 0 1.03186
3 1.10435 0 0.0728822
* 10. MODAL PARTICIPATION MASS RATIOS (%) (cumulative):
# The cumulative modal participation mass ratios (%) for each mode.
# MODE MX MY RMZ
# ------------- ------------- ------------- -------------
1 91.4079 0 98.8953
2 98.8956 0 99.9271
3 100 0 100特に、今回検討しているX方向に着目すると、有効質量比については、3次まででぴったり100%となっています。
よって、モデルの固有モードは3つのモードのみで構成されていることが確認できます。
これは、3つの自由度を有する解析モデルの想定に沿った結果といえます。
また、X方向の3次までの有効質量は、モデル各節点に入力した総重量の値（200+200+200)×0.22480894387096（kN→kipsへの換算係数）/386.0885827（重力加速度（in/s2））=0.349363779と同等となっていることも確認できます。
したがって、対象とする多質点モデルの固有値解析は、検討方向（X方向）に関しては、入力値に対して正しく行われていそうなことが確認できます。
次回はこの固有値解析結果をもう少し詳しく見ていきたいと思います。