免震シミュレータ — 計算方法の解説

1. 解析手法の概要

免震シミュレータは、免震建物を一質点系（SDOF・上部構造を剛体とみなす）または せん断型多質点系（MDOF・上部構造を弾性とみなす）にモデル化し、 免震層の非線形復元力（バイリニア＋Masing則）を考慮した 非線形時刻歴応答解析を行うツールです。 時刻歴積分には Newmark-β法（β=1/4・平均加速度法）または Wilson-θ法（θ=1.4）による直接積分法を用います （解析モデル・積分法はツール右上の「⚙ 設定」で選択。既定は Newmark-β法）。 すべての計算はブラウザ内の JavaScript で実行されます。

解析は次の手順で進みます。

2. モデル化の前提

3. 運動方程式と単位系

地動加速度 ÿ_g を受ける一質点系の運動方程式（相対変位 x）は次のとおりです。

M·ẍ + Q(x, ẋ) = −M·ÿ_g

  M        : 質量 = ΣW / g
  ΣW       : 層重量の合計 [kN]（屋根〜1階の全層）
  g        : 重力加速度 = 980.665 cm/s²
  Q(x, ẋ) : 免震層の復元力 [kN]（バイリニア＋Masing則、変位と履歴に依存）
  ÿ_g     : 地動加速度 [cm/s²]

粘性減衰項 C·ẋ は持ちません（C = 0）。免震層のエネルギー吸収は 復元力特性の履歴ループ（履歴減衰）によってのみ表現されます。

4. 免震層の復元力特性（バイリニア骨格曲線）

免震層（積層ゴム＋ダンパー等）の水平力 Q と水平変位 D の関係を、 2本の直線で表すバイリニア型でモデル化します。 入力値は次の3つです。

骨格曲線（初載荷時の Q-D 関係）は次式で表されます。

|D| ≤ up1 のとき:  Q = K1·D
|D| > up1 のとき:  Q = sign(D) · { K1·up1 + K2·(|D| − up1) }
                     = sign(D)·Qd + K2·D　　（第2分枝）

5. Masing則による履歴ループ

地震応答中の除荷・再載荷経路は Masing則に従います。 すなわち、速度の符号が反転した点（反転点）を原点として、 骨格曲線を縦横2倍に拡大した曲線をたどります。 これにより紡錘形の履歴ループが形成され、ループ面積が吸収エネルギー（履歴減衰）に対応します。

状態変数:
  K      : 履歴分岐の深さ（1 = 骨格曲線上）
  D      : 拡大倍率（骨格曲線 = 1、分岐曲線 = 2）
  U0[·], V0[·] : 各反転点の変位・復元力（スタック）

各時間ステップの処理:
1. 速度符号が反転（ẋ_m · ẋ_m−1 < 0）したら
     反転点 (U_TN, V_TN) を前後ステップの内挿で求めてスタックに積む（K++）
     以後は D = 2 の拡大曲線に乗り移る
2. 現在の分岐の始点側の変位を追い越したら（内側ループが閉じたら）
     外側の分岐へ復帰（K −= 2。K ≤ 3 なら骨格曲線 K = 1 へ戻る）
3. 復元力の評価:
     Q = D · f( (x − U0[K−1]) / D ) + V0[K−1]
     f(·) : 4章のバイリニア骨格曲線

6. 時刻歴積分（Newmark-β法／Wilson-θ法）

運動方程式は増分形に直し、Newmark-β法（β = 1/4）または Wilson-θ法（θ = 1.4）で逐次積分します（「⚙ 設定」で選択。既定は Newmark-β法）。 Wilson-θ法は、時刻 t 〜 t+θΔt の間で加速度が線形に変化すると仮定する無条件安定の陰解法です （θ ≥ 1.37 で無条件安定）。 Newmark-β法（β = 1/4, γ = 1/2）は時間区間内の平均加速度を一定と仮定する無条件安定の陰解法で、 数値減衰を持たない点が特長です（Wilson-θ法は高振動数成分をわずかに数値減衰させます）。

6.1 積分定数

θ = 1.4,  Δt = 地震波の時間刻み

a0 = 6 / (θΔt)²      a1 = 3 / (θΔt)
a2 = 6 / (θΔt)       a3 = θΔt / 2
a4 = a0 / θ           a5 = −a2 / θ
a6 = 1 − 3/θ          a7 = Δt / 2         a8 = Δt² / 6

6.2 各時間ステップの計算

ステップ m（時刻 t_m）ごとに:

1. 有効剛性:      k̂ = k_m + a0·M
     k_m : 免震層の瞬間剛性（下記 6.3）

2. 有効増分荷重:  ΔR̂ = M·( −θ·(ÿg_m − ÿg_m−1) + a2·ẋ_m−1 + 3·ẍ_m−1 )

3. 増分変位:      Δx = ΔR̂ / k̂
     （プログラム上はコレスキー分解による連立方程式求解。
       本ツールは1自由度のためスカラー除算と等価）

4. 応答の更新:
     ẍ_m = a4·Δx + a5·ẋ_m−1 + a6·ẍ_m−1
     ẋ_m = ẋ_m−1 + a7·(ẍ_m + ẍ_m−1)
     x_m = x_m−1 + Δt·ẋ_m−1 + a8·(ẍ_m + 2·ẍ_m−1)

5. 復元力の評価:  Q_m = Masing則（5章）で x_m, ẋ_m から算定

6.3 瞬間剛性の更新

非線形性は、各ステップの変位増分と復元力増分から求めた割線剛性を 次ステップの瞬間剛性として用いることで追跡します（逐次剛性更新法）。

k_m+1 = |Q_m − Q_m−1| / |x_m − x_m−1|

  変位増分が微小（|Δx| ≤ 10⁻⁶）な場合は前ステップの剛性を据え置く

6.4 Newmark-β法（β = 1/4・平均加速度法）の場合

設定で Newmark-β法を選んだ場合は、拡大ステップ θΔt を用いず Δt 上で直接増分を解きます。 有効剛性と有効増分荷重、応答更新式は次のとおりです（γ = 1/2）。

1. 有効剛性:      k̂ = k_m + (4/Δt²)·M + (2/Δt)·C

2. 有効増分荷重:  ΔR̂ = M·( −(ÿg_m − ÿg_m−1) + (4/Δt)·ẋ_m−1 + 2·ẍ_m−1 )
                      + C·( 2·ẋ_m−1 )     ※ β=1/4 では C の加速度項係数は 0

3. 増分変位:      Δx = ΔR̂ / k̂

4. 応答の更新:
     ẍ_m = (4/Δt²)·Δx − (4/Δt)·ẋ_m−1 − ẍ_m−1
     ẋ_m = ẋ_m−1 + (Δt/2)·(ẍ_m + ẍ_m−1)
     x_m = x_m−1 + Δx

6.5 不平衡力の持ち越し（equilibrium correction）

単純な増分法では、折れ点をまたぐステップで生じた不釣合力（接線剛性による予測と 実際の復元力の差）が蓄積し、時間刻みが粗い場合に応答を過小・過大評価することがあります。 本ツールでは、前ステップの運動方程式の残差（不平衡力）を次ステップの 有効増分荷重に加算することでこれを解消しています。

res_m−1 = −M·(ẍ_m−1 + ÿg_m−1) − C·ẋ_m−1 − Q_m−1   （運動方程式の残差）

ΔR̂ ← ΔR̂ + res_m−1

6.6 多質点系への拡張

せん断型多質点モデルでは、上式の M・C・K がマトリクス（自由度 N = 層数+1）になります。 各ステップの連立方程式はコレスキー分解で解きます。 上部構造の層バネは弾性（Q = k·δ）、免震層バネのみバイリニア＋Masing則です。

M = diag(m_1, …, m_N)                     質量マトリクス（各床 m_i = W_i/g）
K = せん断型三重対角マトリクス             層バネ k_i と免震層瞬間剛性から組立
C = (2h/ω₁)·K_upper                        剛性比例型減衰（免震層の行は 0）

  h   : 上部構造の減衰定数（設定値、既定 2%）
  ω₁ : 上部構造（1F固定）の1次固有円振動数
  K_upper : 上部構造の層剛性のみで組み立てた剛性マトリクス

この C は基礎固定系の r 次モードに対し減衰定数 ξ_r = h·ω_r/ω₁（1次でちょうど h）を与えます。 免震層のエネルギー吸収は履歴減衰のみで表現する方針のため、C への免震層の寄与はゼロとしています。

固有値解析は、質量で基準化した標準固有値問題 A = M^−1/2·K·M^−1/2 を Jacobi法で解き、固有周期とモード形を求めます。 表示する固有値解析結果は上部構造のみ（1F床位置で固定した基礎固定系）のものです。 免震層は剛性が変形に応じて変化する（バイリニア）ため、固有値解析には含めません。

7. 入力地震波と振幅調整

7.1 内蔵地震波

次の5波形を内蔵しています（加速度時刻歴、単位 cm/s² = gal）。

7.2 最大速度による振幅調整

免震構造の設計では入力地震動を最大速度（Vmax）で基準化するのが慣例です。 本ツールでも、各波形の原波の最大速度に対する倍率をあらかじめ求めておき、 選択した Vmax（25・50・75・100 cm/s）に応じて加速度時刻歴全体を定数倍します。

調整倍率 = (原波の Vmax を 25 cm/s に基準化する係数) × (指定 Vmax / 25)

例: EL CENTRO-NS の 25 cm/s 基準化係数 = 0.747
    → Vmax = 50 cm/s を指定すると倍率 1.495 を全データに乗じる
「任意のVmax」を指定した場合も同じ式で任意の整数値に対応

7.3 カスタム地震波

8. 計算結果の表示

9. 検証例

サンプル「K神戸邸」を用いた計算例です。ツールで同じ条件を選ぶと再現できます。

10. 適用範囲と制限

1. 上部構造は剛体または弾性 — 剛体1質点モデルでは上部構造の層間変形・高次モードの 応答増幅は評価できません。せん断型多質点モデルでは弾性範囲の層間変形・高次モードを評価できますが、 上部構造の塑性化は考慮しません。
2. 水平1方向のみ — 上下動・2方向同時入力・ねじれは考慮しません。
3. 免震層は粘性減衰なし — オイルダンパー等の速度依存型減衰は表現できません。 免震層のエネルギー吸収はバイリニア履歴（履歴減衰）のみです （多質点時の上部構造には剛性比例型粘性減衰を考慮できます）。
4. 免震層の特性は一定 — 面圧依存・速度依存・温度依存などの特性変化、 引抜き・ハードニングは考慮しません。
5. ステップ内の平衡反復なし — Newton-Raphson 等の反復は行わず、 割線剛性の逐次更新と不平衡力の持ち越し（6.5節）で非線形性を追跡します。 折れ点近傍で微小な数値誤差が生じ得ます。
6. 学習・検討用途 — 本ツールは免震構造の挙動理解を目的としたものです。 実設計には使用できません。

11. 参考文献

• 大崎順彦『建築振動理論』彰国社

関連ページ: 免震シミュレータ  |  使用手順  |  ハイパー鋼材表  |  かんたん骨組解析