応答解析装置、方法及びプログラム
【課題】波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現する。
【解決手段】波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルの時刻歴応答解析を行うにあたり、外部地盤を含むエネルギー伝達境界の伝達境界マトリクスを演算し(48)、インパルス応答を規定する数式として、物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を少なくとも含む数式を用い、振動が各周波数のときの前記伝達境界マトリクスの値に基づいて外部地盤を含むエネルギー伝達境界の時間領域におけるインパルス応答を演算し(50)、演算されたインパルス応答を用いて解析対象モデルに対して時間領域での応答解析を行う(54)。
【解決手段】波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルの時刻歴応答解析を行うにあたり、外部地盤を含むエネルギー伝達境界の伝達境界マトリクスを演算し(48)、インパルス応答を規定する数式として、物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を少なくとも含む数式を用い、振動が各周波数のときの前記伝達境界マトリクスの値に基づいて外部地盤を含むエネルギー伝達境界の時間領域におけるインパルス応答を演算し(50)、演算されたインパルス応答を用いて解析対象モデルに対して時間領域での応答解析を行う(54)。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は応答解析装置、方法及びプログラムに係り、特に、波動境界モデルを含む解析対象モデルに対して時間領域での振動問題の応答解析を行う応答解析装置、該応答解析装置に適用可能な応答解析方法、及び、コンピュータを前記応答解析装置として機能させるための応答解析プログラムに関する。
【背景技術】
【０００２】
地震時の建物の挙動や損傷を精度良く評価するためには、建物と地盤との相互作用効果を考慮して地震応答解析を行う必要がある。また地震応答解析は、周波数領域で応答解析を行う周波数応答解析と、時間領域で応答解析を行う時刻歴応答解析とに大別されるが、大地震時には地盤・建物とも大きく塑性化するので、建物の挙動等を精度良く評価するためには、非線形の解析が可能な時刻歴応答解析を行うことが望ましい。これに関連して本願発明者は、時刻歴応答解析の実現を目的として、周波数領域の複素関数である地盤の動的剛性を時間領域で表されるインパルス応答へ変換する技術を既に提案している（例えば非特許文献１や特許文献１を参照）。
【０００３】
ところで、地震応答解析では地盤を有限要素法(ＦＥＭ)の非線形ソリッド要素を用いてモデル化することが多い。地盤は本来半無限的な広がりを持つのに対し、ＦＥＭ解析では解析領域を有限としてモデル化する必要があるため、解析対象モデルの側面や底面には波動境界モデルが設定される。建物で励起された波動は主として、下方へは実体波として、側方へは表面波として伝播するが、表面波の方がより遠方まで伝播するため、解析対象モデルの側面及び底面のうち、側面に設定する波動境界モデルが解析結果により大きな影響を及ぼす。時間領域の解析に適用可能な側面境界モデルとしては、例えば繰返し境界、粘性境界、無反射境界が知られている。このうち、繰返し境界は、モデル化は極めて容易であるものの波動境界としての精度が低く、モデル化領域(解析対象モデル)を大きくとる必要があるので解析負荷が大きいという問題がある。また、粘性境界は斜め方向からの波動に対しての精度が低く、解析精度を上げるため同様にモデル化領域を大きくする必要がある。更に、無反射境界は現実の問題への適用に課題が多く、あまり用いられていない。
【０００４】
一方、より高精度な波動境界モデルとしてエネルギー伝達境界(Energy Transmitting Boundary、又は、Consistent Transmitting Boundary)が知られている（例えば非特許文献２を参照）。エネルギー伝達境界は、剛基盤上に平行成層をなす地盤の外端に設置され、水平方向には厳密で、上下方向は要素の変位仮定に従う（線形一次要素であれば、変位は直線的に変化する）高精度の境界であり、任意方向からの波動をほぼ完全に吸収する特性を有している。このため、解析対象モデルの側面の波動境界モデルとして上記のエネルギー伝達境界を用いれば、モデル化領域を大幅に小さくすることで、同等の解析精度での応答解析を、より小さな解析負荷で実現できるものと期待される。
【特許文献１】特許３８７８６２６号公報
【非特許文献１】中村尚弘，「地盤インピーダンスの時間領域変換による成層地盤に埋込まれた構造物の地震応答解析」，日本建築学会構造系論文集，２００３年５月，第５６７号，ｐ．６３−７０
【非特許文献２】川本眺万，林正夫，「地盤工学における有限要素解析 土質力学と岩盤力学へのアプローチ」，培風館，１９７８年，ｐ．３４１−３６０
【発明の開示】
【発明が解決しようとする課題】
【０００５】
しかしながら、エネルギー伝達境界の特性は、伝達境界マトリクスと補正力ベクトルによって表されるが、このうち伝達境界マトリクスは、周波数領域の複素数であり強い振動数依存性を有しているので、時間領域に変換することは困難であった。このため、側面境界モデルとしてエネルギー伝達境界を適用した解析対象モデルによる応答解析は、周波数領域での等価線形解析に限られているのが実情であった。
【０００６】
本発明は上記事実を考慮して成されたもので、波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現する応答解析装置、応答解析方法及び応答解析プログラムを得ることが目的である。
【課題を解決するための手段】
【０００７】
以下にエネルギー伝達境界の概要を説明する。地盤をモデル化すると共に、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界を設け、モデルの中心位置に関して逆対称条件を設けることにより１／２モデル化する。右側の１／２モデルを図１に示す。図１に示す１／２モデルでは、モデルの右側にのみエネルギー伝達境界が存在し、その左側には内部領域が位置し、その右側には外部領域(自由地盤ともいう)が位置している。この１／２モデルにおいて、２次元面内問題での自由地盤の波動方程式は次の(１)式で表される。
([Ａ]ｋ２＋i[Ｂ]ｋ＋[Ｇ]−ω２[Ｍ]){ｕ＊(ω)}＝{0} …(１)
ここで、ｋは波数、{ｕ＊(ω)}は自由地盤の変位ベクトルである。また、[Ａ],[Ｂ],[Ｇ],[Ｍ]は２ｎ×２ｎ（ｎは節点数）のマトリクスで、以下の要素毎のサブマトリクス[Ａ]j,[Ｂ]j,[Ｇ]j,[Ｍ]jの重ね合せで表される。各サブマトリクスは、各要素のラメ定数Ｇj,λjを用いて次の(２)〜(５)式で表される。ｈj,ρjは各要素の高さと密度である。
【０００８】
【数１】
【０００９】
ここで、
[Ｃ]＝[Ｇ]−ω２[Ｍ] …(６)
とし、次の(７)式の波数ｋに関する固有値問題を解く。
|[Ａ]ｋ２＋i[Ｂ]ｋ＋[Ｃ]|＝{0} …(７)
得られた４ｎ個の固有モードより、右方向に伝播する成分２ｎ個を抽出し、モードマトリクス[Ｖ]とする。これより、１／２モデルの右側に存在するエネルギー伝達境界の伝達境界マトリクスは次の(８)式の[Ｒ]で表される。但し、[Ｄ]は下記の(９)式で表されるサブマトリクス[Ｄ]jの重ね合せで表される。
[Ｒ]＝i・[Ａ][Ｖ][Ｋ][Ｖ]−１＋[Ｄ] …(８)
【００１０】
【数２】
【００１１】
なお、先の(８)式ではマトリクスが非対称となるが、以下では対称位置の項を平均することによりマトリクスを対称化する。これよりモデル全体の運動方程式は次の(10)式で表される。
(−ω２[ＭＩ]＋[ＫＩ]＋[Ｒ]){ｕ(ω)}＝−ｙ"(ω)[ＭＩ]{1}＋{Ｆ(ω)} …(10)
但し、 [ＭＩ],[ＫＩ],{ｕ(ω)}は各々、内部領域の質量マトリクス、剛性マトリクス、変位ベクトルである。また{Ｆ(ω)}は境界力ベクトルであり、次の(11)式で表せる。
{Ｆ(ω)}＝([Ｒ]−[Ｑ]){ｕ＊(ω)} …(11)
上記の(11)式のうち、[Ｑ]{ｕ＊(ω)}は補正力（切欠き力ともいう）ベクトルであり、地震動が鉛直下方より入射する場合は[Ｑ]＝[Ｄ]、それ以外の場合は[Ｑ]＝ｉｋ[Ａ]＋[Ｄ]となる。但しｋは波数である。また、[Ｒ],{ｕ＊(ω)},{Ｆ(ω)}は、(８)式までは境界部の自由度数２ｎで定義したが、内部領域の自由度に重ね合わせるため、(10)式以降では内部領域の自由度数に拡張して用いている。
【００１２】
ここで、上述したエネルギー伝達境界の伝達境界マトリクスは周波数領域の複素数であり、強い振動数依存性を有しているが、本願発明者はエネルギー伝達境界の伝達境界マトリクスについて検討・考察を重ね、エネルギー伝達境界の伝達境界マトリクスが、本質的には、同じく周波数領域の複素数で表され強い振動数依存性を有している地盤等の物体の動的剛性と類似する物理量であることに想到した。そして本願発明者は、伝達境界マトリクスが本質的には物体の動的剛性と類似する物理量であるならば、非特許文献１や特許文献１で提案している技術を適用することで、時間領域のインパルス応答へ変換できる可能性があることに想到し、当該事項を確認するために以下の解析検討を行った。
【００１３】
この解析検討では、図１に示した地盤モデル(１／２モデル)を用い、演算に用いる地盤の物性値を、剪断波速度Ｖｓ＝400(ｍ／ｓ)、ポアソン比ν＝0.4、密度ρ＝2.0(ｔ／ｍ３)、減衰定数ｈ＝0.02とした。また、伝達境界マトリクスから時間領域のインパルス応答の演算には、本願発明者が非特許文献１で提案した演算方法（便宜的にＡ法と称する）、本願発明者が特許文献１で提案した演算方法（便宜的にＢ法と称する）、当該Ｂ法に対し後述する請求項７記載の発明に従って同時成分ｍ０,ｋ０の修正を更に行う演算方法（便宜的にＣ法と称する）を各々用いると共に、時間領域への変換（インパルス応答の演算）の条件として次の表１に示す条件を適用し、表１に示す各振動数における伝達境界マトリクスの値から各演算方法によってインパルス応答を各々演算した。
【００１４】
【表１】
【００１５】
また、上記演算の精度を確認するために、上記演算によって得られたインパルス応答を周波数領域の値へ再変換し、元のデータ（伝達境界マトリクスの元の値）との偏差に基づいてインパルス応答の演算精度を各演算方法毎に評価した。なお、周波数領域の値への再変換には、インパルス応答の各成分のうち、同時成分ｋ０,ｃ０(,ｍ０)と時間遅れ成分ｋ１〜ｋ１０,ｃ１〜ｃ１０を用いた。インパルス応答の演算方法としてＣ法を適用した場合の結果を図２に示す。図２に示すように、インパルス応答の演算方法としてＣ法を適用した場合の再現値（インパルス応答を周波数領域の値へ再変換して得られる値）は、ほぼ全てのデータ点において、元のデータと一致しているか、又は元のデータとの偏差が非常に小さく、エネルギー伝達境界の伝達境界マトリクスからエネルギー伝達境界のインパルス応答を精度良く演算できることが明らかとなった。また、図示は省略するが、Ａ法及びＢ法についても、演算精度はＣ法＞Ｂ法＞Ａ法となるものの、エネルギー伝達境界の伝達境界マトリクスから、応答解析に用いるのに十分な精度でエネルギー伝達境界のインパルス応答を演算できることが確認された。
【００１６】
上記に基づき請求項１記載の発明に係る応答解析装置は、地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算するマトリクス演算手段と、物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算するインパルス応答演算手段と、前記インパルス応答演算手段によって演算された前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行う解析手段と、を含んで構成されている。
【００１７】
請求項１記載の発明では、地盤単独又は地盤と建物から成る解析対象物体に対応し、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルを解析対象としており、マトリクス演算手段は、エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算する。また、インパルス応答演算手段は、物体を振動させる外力と物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式（すなわち、前出のＡ法〜Ｃ法の何れかに対応する演算式）を用い、振動が各周波数のときの伝達境界マトリクスの値に基づいてエネルギー伝達境界のインパルス応答を演算する。これにより、先にも説明したように、エネルギー伝達境界のインパルス応答を、応答解析に用いるのに十分な精度で得ることができる。そして解析手段は、インパルス応答演算手段によって演算されたエネルギー伝達境界のインパルス応答を用いて、解析対象モデルに対して時間領域での応答解析を行う。これにより、波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現することができる。
【００１８】
なお、請求項１記載の発明におけるエネルギー伝達境界の伝達境界マトリクスは、先に説明し、かつ請求項２にも記載したように、内部領域の質量マトリクスを[ＭＩ]、剛性マトリクスを[ＫＩ]、変位ベクトルを{ｕ(ω)}、境界力ベクトルを{Ｆ(ω)}、補正力ベクトルを[Ｑ]{ｕ＊(ω)}とし、
(−ω２[ＭＩ]＋[ＫＩ]＋[Ｒ]){ｕ(ω)}＝−ｙ"(ω)[ＭＩ]{1}＋{Ｆ(ω)} …(10)
{Ｆ(ω)}＝([Ｒ]−[Ｑ]){ｕ＊(ω)} …(11)
解析対象モデル全体の運動方程式を上記の(10),(11)式で表したときに、上記の(10),(11)式の[Ｒ]で表すことができる。
【００１９】
また、詳細は後述するが、本願発明者は、本発明が２次元面内問題、２次元面外問題及び軸対称問題の何れにも有効であることを確認しており、請求項１記載の発明に係る解析対象モデルに前記波動境界モデルとして設けられているエネルギー伝達境界としては、請求項３に記載したように、２次元面内問題、２次元面外問題、及び、軸対称問題の何れかに対応するモデルを適用することができ、２次元面内問題、２次元面外問題及び軸対称問題の何れに対応するモデルであってもよい。
【００２０】
また、請求項１記載の発明において、インパルス応答演算手段は、例えば請求項４に記載したように、物体の変位に依存するインパルス応答の同時成分をｋ０、物体の速度に依存するインパルス応答の同時成分をｃ０、物体の変位に依存するインパルス応答のΔt刻みの時間遅れ成分をｋj、物体の速度に依存するインパルス応答のΔt刻みの時間遅れ成分をｃj（但しｊは自然数でtj＝Δt・j）、時間領域での物体の変位をu(t)、速度をu'(t)、加速度をu"(t)としたときに、インパルス応答ＦＡ(t)を規定する数式として、
【００２１】
【数３】
【００２２】
上記(12)式を用い、前記振動がＮ種（Ｎ＝ｎ＋１）の周波数のときの伝達境界マトリクスの値に基づいてインパルス応答を求めるように構成することができる。上記の(12)式は前出のＡ法に対応しており、エネルギー伝達境界のインパルス応答を、応答解析に用いるのに十分な精度で得ることができる。
【００２３】
また、請求項１記載の発明において、インパルス応答演算手段は、例えば請求項５に記載したように、インパルス応答を規定する数式として、剛性項及び減衰項に加えて、物体の加速度に依存し少なくとも同時成分を含んで成る質量項も含む数式を用いて、エネルギー伝達境界のインパルス応答を演算するように構成してもよい。
【００２４】
より詳しくは、インパルス応答演算手段は、例えば請求項６に記載したように、物体の変位に依存するインパルス応答の同時成分をｋ０、物体の速度に依存するインパルス応答の同時成分をｃ０、物体の加速度に依存するインパルス応答の同時成分をｍ０、物体の変位に依存するインパルス応答のΔt刻みの時間遅れ成分をｋj、物体の速度に依存するインパルス応答のΔt刻みの時間遅れ成分をｃj（但しｊは自然数でtj＝Δt・j）、時間領域での物体の変位をu(t)、速度をu'(t)、加速度をu"(t)としたときに、前記質量項も含む数式として、前記インパルス応答ＦB(t)を規定する数式である、
【００２５】
【数４】
【００２６】
上記(13)式を用い、前記振動がＮ種（Ｎ＝ｎ＋１）の周波数のときの伝達境界マトリクスの値に基づいてインパルス応答を演算するように構成してもよい。上記の(13)式は前出のＢ法に対応しており、エネルギー伝達境界のインパルス応答をＡ法によって演算する場合（インパルス応答演算手段を請求項４記載の構成とする場合）よりも、エネルギー伝達境界のインパルス応答としてより高精度な値を得ることができ、解析対象モデルに対する時間領域での応答解析の解析精度を向上させることができる。
【００２７】
また、請求項６記載の発明において、インパルス応答演算手段は、例えば請求項７に記載したように、前記(13)式における物体の加速度に依存するインパルス応答の同時成分ｍ０に対する修正値Δｍ0、物体の変位に依存するインパルス応答の同時成分ｋ０に対する修正値Δｋ0を下記の(14)式によって演算し、
【００２８】
【数５】
【００２９】
（但し、上記(14)式において、
【００３０】
【数６】
【００３１】
であり、上記(15)式におけるＲe(ＳＢ(ωｉ))は、前記(13)式を用いて演算した前記エネルギー伝達境界の前記インパルス応答から再現された、下記の(16)式で表される前記伝達境界マトリクスＳＢ(ω)の振動数ωｉでの実部の値を表し、
【００３２】
【数７】
【００３３】
上記(15)式におけるＲe(Ｄ(ωｉ))は、(13)式に基づくエネルギー伝達境界のインパルス応答の演算に用いた、振動数ωｉでの伝達境界マトリクスのデータＤ(ωｉ)のうちの実部の値を表す）、演算した修正値Δｍ0,Δｋ0を用いて同時成分ｍ０,ｋ０を修正することが好ましい。
【００３４】
本願発明者による前述の解析検討において、Ｂ法を適用して算出したエネルギー伝達境界のインパルス応答を上記の(16)式によって周波数領域の値へ再変換した場合、虚部については元のデータと精度良く一致しているものの、実部については元のデータとの偏差が比較的大きい、という結果となった(Ａ法についても同様)。これは、エネルギー伝達境界の伝達境界マトリクスが因果的ではないので時間領域への厳密な変換が不可能であり、Ｂ法によるインパルス応答の演算が近似的な演算であることが原因と推察される。このため、本願発明者は、Ｂ法を適用して算出したエネルギー伝達境界のインパルス応答のうち、実部の同時成分であるｍ０,ｋ０に対し、先の(14),(15)式に従って修正値Δｍ0,Δｋ0を演算し、演算した修正値Δｍ0,Δｋ0を用いて同時成分ｍ０,ｋ０を修正した（Ｃ法）後に、Ｃ法によって得られたエネルギー伝達境界のインパルス応答を周波数領域の値へ再変換し、元のデータと比較する解析検討を行った。その結果、前出の図２にも示すように、実部についても元のデータと精度良く一致していることが確認された。
【００３５】
請求項７記載の発明は、先の(14),(15)式に従って修正値Δｍ0,Δｋ0を演算し、演算した修正値Δｍ0,Δｋ0を用いて同時成分ｍ０,ｋ０を修正するので、上述した解析検討の結果からも明らかなように、エネルギー伝達境界のインパルス応答をＡ法やＢ法によって演算する場合（インパルス応答演算手段を請求項４又は請求項６記載の構成とする場合）よりも、エネルギー伝達境界のインパルス応答として更に高精度な値を得ることができ、解析対象モデルに対する時間領域での応答解析の解析精度を更に向上させることができる。
【００３６】
また、請求項７記載の発明において、インパルス応答演算手段は、例えば請求項８に記載したように、前記(13)式における物体の速度に依存するインパルス応答の同時成分をｃ０に対する修正値Δｃ0も下記の(17)式によって演算し、
Δｃ0＝−Ｅ／Ｂ …(17)
（但し、上記(17)式において、
【００３７】
【数８】
【００３８】
であり、上記(18)式におけるＩm(ＳＢ(ωｉ))は、前記(13)式を用いて演算したエネルギー伝達境界のインパルス応答から再現された、前記(16)式で表される伝達境界マトリクスＳＢ(ω)の振動数ωｉでの虚部の値を表し、上記(17)式におけるＩm(ＳＢ(ωｉ))は、(13)式に基づくエネルギー伝達境界のインパルス応答の演算に用いた、振動数ωｉでの伝達境界マトリクスのデータＤ(ωｉ)のうちの虚部の値を表す）、演算した修正値Δｃ0を用いて同時成分ｃ０も修正することが好ましい。
【００３９】
本願発明者による前述の解析検討では、Ｃ法によるエネルギー伝達境界のインパルス応答の演算を様々な条件で行ったが、この解析検討により、例えばインパルス応答の演算に用いるデータ点の振動数軸上における位置が等間隔でない等の場合には、演算したインパルス応答のうちの虚部についての精度が低下する（演算したインパルス応答を周波数領域の値へ再変換した場合の元のデータとの偏差が虚部について比較的大きくなる）ことが判明した。このため、本願発明者はＣ法を適用して算出したエネルギー伝達境界のインパルス応答のうち、虚部についての精度が低下していたインパルス応答に対し、虚部の同時成分であるｃ０に対し、先の(18)式に従って修正値Δｃ0を演算し、演算した修正値Δｃ0を用いて同時成分ｃ０を修正した後に、当該修正を経たエネルギー伝達境界のインパルス応答を周波数領域の値へ再変換し、元のデータと比較する解析検討を行った。その結果、虚部の精度が向上することが確認された。
【００４０】
請求項８記載の発明は、先の(18)式に従って修正値Δｃ0を演算し、演算した修正値Δｃ0を用いて同時成分ｃ０を修正するので、例えばインパルス応答の演算に用いるデータ点の振動数軸上における位置が等間隔でない等のように、演算したインパルス応答のうちの虚部についての精度が低下し易い条件であっても、エネルギー伝達境界のインパルス応答として更に高精度な値を得ることができ、解析対象モデルに対する時間領域での応答解析の解析精度を更に向上させることができる。なお、請求項８記載の発明に係る同時成分ｃ0の修正は常時行う必要はなく、例えばインパルス応答の演算に用いるデータ点の振動数軸上における位置が等間隔でない等のように、インパルス応答のうちの虚部についての精度が低下することが予想できる条件のときにのみ行うようにしてもよい。
【００４１】
請求項９記載の発明に係る応答解析方法は、コンピュータに、地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算させ、物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算させ、前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行わせるので、請求項１記載の発明と同様に、波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現することができる。
【００４２】
請求項１０記載の発明に係る応答解析プログラムは、コンピュータを、地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算するマトリクス演算手段、物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算するインパルス応答演算手段、及び、前記インパルス応答演算手段によって演算された前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行う解析手段として機能させる。
【００４３】
請求項１０記載の発明に係る応答解析プログラムは、コンピュータを、上記のマトリクス演算手段、インパルス応答演算手段及び解析手段として機能させるためのプログラムであるので、コンピュータが請求項１０記載の発明に係る応答解析プログラムを実行することで、コンピュータが請求項１０に記載の応答解析装置として機能することになり、請求項１記載の発明と同様に、波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現することができる。
【発明の効果】
【００４４】
以上説明したように本発明は、波動境界モデルとして設けられたエネルギー伝達境界の伝達境界マトリクスを演算し、インパルス応答を規定する数式として、少なくとも、物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの伝達境界マトリクスの値に基づいてエネルギー伝達境界のインパルス応答を演算し、エネルギー伝達境界のインパルス応答を用いて、解析対象モデルに対して時間領域での応答解析を行うので、波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現できる、という優れた効果を有する。
【発明を実施するための最良の形態】
【００４５】
以下、図面を参照して本発明の実施形態の一例を詳細に説明する。図３には本発明を適用可能なコンピュータ１０が示されている。コンピュータ１０は、ＣＰＵ１０Ａ、ＲＯＭやＲＡＭ等から成るメモリ１０Ｂ、ＨＤＤ(Hard Disk Drive)やフラッシュメモリ等から成る不揮発性の記憶部１０Ｃを備えており、ＣＲＴ又はＬＣＤから成るディスプレイ１２、キーボード１４、マウス１６が各々接続されている。
【００４６】
コンピュータ１０の記憶部１０Ｃには、後述する地震応答解析処理を行うための地震応答解析プログラムがインストールされている。この地震応答解析プログラムは、請求項１０に記載の応答解析プログラムに対応している。また、コンピュータ１０は請求項９及び請求項１０に記載のコンピュータに対応しており、ＣＰＵ１０Ａが地震応答解析プログラムを実行することで、請求項１等に記載の発明に係る応答解析装置として機能する。なお、コンピュータ１０としてはパーソナル・コンピュータ(ＰＣ)が好適であるが、これに限られるものではなく、例えばワークステーションであってもよいし、汎用の大型コンピュータであってもよい。
【００４７】
次に本実施形態の作用として、ＣＰＵ１０Ａが地震応答解析プログラムを実行することで実現される地震応答解析処理について、図４(及び図５)を参照して説明する。なお、以下では図６に示すように、地盤と当該地盤に埋め込まれた地下部分を有する建物が一体となった解析対象物体をモデル化すると共に、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた構成の解析対象モデルに対して時間領域での応答解析を行う場合を例に説明する。
【００４８】
地震応答解析処理では、まずステップ３０において、内部地盤（解析対象モデルの内部領域）及び建物の動特性マトリクスを演算するためのデータを取得する。なお、本実施形態において、動特性マトリクスは、地盤や建物の応答解析を行うのに必要な特性を表すマトリクスを意味し、具体的には、質量マトリクス、減衰マトリクス及び剛性マトリクスから構成され、ステップ３０で取得されるデータとしては、例えば内部地盤の各メッシュ領域毎の物性値（例えば剪断波速度Ｖｓ、ポアソン比ν、密度ρ、減衰定数ｈ等）や、建物の各階毎の定数（質量、剪断剛性、回転慣性、減衰率ｈ等）等が挙げられる。次のステップ３２では、ステップ３０で取得したデータに基づき、有限要素法(ＦＥＭ)による解析を行うことで内部地盤及び建物の動特性マトリクスを演算する。これにより、前出の(10)式に示す内部領域の質量マトリクス[ＭＩ]、内部領域の剛性マトリクス[ＫＩ]、前出の(６)式に示す減衰マトリクス[Ｃ]が各々算出される。なお、本実施形態では材料減衰の振動数非依存性を表す減衰モデルとして因果的履歴減衰モデルを用いており、上記の減衰マトリクス[Ｃ]は因果的履歴減衰モデルを用いて算出される。
【００４９】
ステップ３４では入力地震動のデータを取得し、次のステップ３６では外部地盤(無限遠方地盤)の動特性マトリクス及び地震応答を演算するためのデータとして地盤各層の物性値（例えば剪断波速度Ｖｓ、ポアソン比ν、密度ρ、減衰定数ｈ等）を取得する。ステップ３８では、先のステップ３４で取得した入力地震動のデータと、先のステップ３６で取得した外部地盤の地震応答を演算するためのデータに基づいて、外部地盤の地震応答を演算する(図５も参照)。これにより、外部地盤の地震応答として、先の(１)式における自由地盤の変位ベクトル{ｕ＊(ω)}が算出される。またステップ４０では、ステップ３８の演算によって得られた外部地盤の地震応答に対して逆フーリエ変換を行うことで、時間領域における外部地盤の地震応答を演算する(図５も参照)。
【００５０】
次のステップ４２では、ステップ３８の演算によって得られた外部地盤の地震応答と先のステップ３６で取得したデータに基づいて、外部地盤の補正力ベクトルを演算する(図５も参照)。これにより、外部地盤の補正力ベクトルとして、先の(11)式における補正力ベクトル"[Ｑ]{ｕ＊(ω)}"が算出される。ステップ４４では、ステップ４２の演算によって得られた外部地盤の補正力ベクトルに対して逆フーリエ変換を行うことで、時間領域における外部地盤の補正力ベクトルを演算する。
【００５１】
またステップ４６では、先のステップ３６で取得した外部地盤の地震応答を演算するためのデータに基づいて、外部地盤の動特性マトリクスを演算する(図５も参照)。これにより、外部地盤の動特性マトリクスとして、先の(1)式におけるマトリクス[Ａ],[Ｂ],[Ｇ],[Ｍ]やモードマトリクス[Ｖ]、先の(９)式で表されるサブマトリクス[Ｄ]jの重ね合せで表されるマトリクス[Ｄ]が算出される。次のステップ４８では、ステップ４６の演算によって得られた外部地盤の動特性マトリクスを先の(８)式に代入することで、外部地盤を含むエネルギー伝達境界の伝達境界マトリクスを演算する(図５も参照)。これにより、外部地盤を含むエネルギー伝達境界の伝達境界マトリクスとして、先の(８)式で表される伝達境界マトリクス[Ｒ]が算出される。
【００５２】
ステップ５０では、ステップ４８の演算によって得られた外部地盤を含むエネルギー伝達境界の伝達境界マトリクスを時間領域へ変換することで、外部地盤を含むエネルギー伝達境界のインパルス応答を算出する。エネルギー伝達境界のインパルス応答の演算は、具体的には以下のようにして行うことができる。
【００５３】
すなわち、まずステップ４８の演算によって得られた伝達境界マトリクスのデータから、予め設定された演算対象の周波数範囲内のＮ種の周波数（Ｎ種の角振動数ω1,…,ωN）における伝達境界マトリクスの値を表すＮ個の複素データＳB(ω1),…,ＳB(ωN)を各々抽出する。なお、演算対象の周波数範囲としては、例えば０〜20(Hz)の範囲を適用することができる。次に、インパルス応答を規定する数式として先の(13)式を、伝達境界マトリクスを規定する数式として先の(16)式を用い、(13),(16)式に基づいて導出した次の(19),(20)式の連立方程式にＮ個の複素データＳ(ω1),…,Ｓ(ωN)を代入し、この連立方程式の解を求めることで、外部地盤を含むエネルギー伝達境界のインパルス応答を、予め設定されたΔｔ刻みで演算する。
【００５４】
【数９】
【００５５】
この演算により、外部地盤を含むエネルギー伝達境界のインパルス応答として、物体の変位に依存するインパルス応答の同時成分ｋ0、物体の速度に依存するインパルス応答の同時成分ｃ0、物体の加速度に依存するインパルス応答の同時成分ｍ0のデータが得られると共に、物体の変位に依存するインパルス応答の時間遅れ成分ｋjのデータがΔt刻みでｎ個（ｎ＝Ｎ−１）得られ、物体の速度に依存するインパルス応答の時間遅れ成分ｃjのデータがΔt刻みでｎ−１個得られる。
【００５６】
続いて、先の(14),(15)式に従って修正値Δｍ0,Δｋ0を演算し、演算した修正値Δｍ0,Δｋ0を用いてインパルス応答の同時成分ｍ０,ｋ０を修正する。更に、例えばインパルス応答の演算に用いるデータ点の振動数軸上における位置が等間隔でない等のように、演算したインパルス応答のうちの虚部についての精度が低下し易い条件である場合には、先の(18)式に従って修正値Δｃ0も演算し、演算した修正値Δｃ0を用いて同時成分ｃ０を修正する。上記の演算により、外部地盤を含むエネルギー伝達境界のインパルス応答が得られ、解析対象モデルに対して時間領域での応答解析を行うことが可能となる。
【００５７】
ステップ５２では、先のステップ３２の演算によって得られた内部地盤及び建物の動特性マトリクス、先のステップ４０の演算によって得られた時間領域における外部地盤の地震応答、先のステップ４４の演算によって得られた時間領域における外部地盤の補正力ベクトル、先のステップ５０の演算によって得られた時間領域における外部地盤を含むエネルギー伝達境界の動特性マトリクス（外部地盤を含むエネルギー伝達境界のインパルス応答）を結合することで、解析対象モデルに対して時間領域の応答解析を行うための時間領域の運動方程式(次の(22)式参照)を求める(図５も参照)。
[ＭI]{ｕ"(t)}＋[ＫI]{ｕ(t)}＋{Ｒ0}＝−ｙ"(t)[ＭI]{1}＋{Ｆ(t)}＋{Ｒf(t)} …(21)
但し、
【００５８】
【数１０】
【００５９】
そしてステップ５４では、先のステップ３４で取得した入力地震動のデータが表す各時刻における地震動をステップ５２の処理によって得られた運動方程式(上記の(22),(23)式)に順に代入し、解析対象モデルに対して時間領域の地震応答解析(時刻歴応答解析)を行う。これにより、内部地盤、建物及び外部地盤が一体となった解析対象物体に地震動が入力された場合の解析対象物体の挙動を、地盤・建物の塑性化も考慮して精度良く評価することができる。
【００６０】
なお、上記ではエネルギー伝達境界のインパルス応答を、先の(13),(16)式に基づいて導出した先の (19),(20)式の連立方程式に複素データを代入することで演算する態様を説明したが、本発明に係るＢ法はこれに限定されるものではなく、(19),(20)式から算出される時間遅れ成分ｋj,ｃjのうち最初のｎ’個の時間遅れ成分が有意な値であり、それ以降の値が無視できるのであれば、先の(13),(16)式におけるΣの上限をＮ−２からｎ’に低減することも可能である（Ａ法の演算式である(12)式におけるΣの上限も同様にＮ−１からｎ’に低減可能）。これを先の(13),(16)式に適用した場合、次の(23),(24)式が得られる。
【００６１】
【数１１】
【００６２】
なお、(23),(24)式におけるｎ’＜Ｎ−１である。この(23),(24)式に基づいて導出した連立方程式に複素データを代入することで、エネルギー伝達境界のインパルス応答を演算するようにしてもよい。
【００６３】
また、上記ではＣ法を適用してエネルギー伝達境界のインパルス応答を演算する態様を説明したが、本発明はこれに限定されるものではなく、Ｂ法やＡ法を適用してエネルギー伝達境界のインパルス応答を演算することも権利範囲に含むものである。
【００６４】
また、上記では本発明を２次元面内問題に適用した態様を説明したが、これに限定されるものではなく、後述する第２実施例からも明らかなように本発明は２次元面外問題や軸対称問題にも有効であり、本発明は２次元面外問題や軸対称問題に適用してもよいことは言うまでもない。
【００６５】
更に、上記では本発明に係る応答解析プログラムに対応する地震応答解析プログラムがコンピュータ１０の記憶部１０Ｃに予め記憶（インストール）されている態様を説明したが、本発明に係る応答解析プログラムは、ＣＤ−ＲＯＭやＤＶＤ−ＲＯＭ等の記録媒体に記録されている形態で提供することも可能である。
【実施例】
【００６６】
次に、本発明の効果を確認するために本願発明者が行った解析検討の詳細を説明する。
【００６７】
〔第１実施例〕
この解析検討では、伝達境界の時間領域化とそれを用いた時刻歴応答解析の妥当性及び有効性を確認することを目的として、まず、地盤モデルに対して周波数領域の伝達境界マトリクスを作成し、時間領域変換の性状を評価する。更に、前記地盤モデルが表す地盤に埋込まれた建物を想定し、建物−地盤連成系の時刻歴地震応答解析を行った。
【００６８】
(１)伝達境界マトリクスの時間領域変換
この解析検討では、先に説明したように、図１に示した地盤モデル(１／２モデル)を用い、演算に用いる地盤の物性値を、剪断波速度Ｖｓ＝400(ｍ／ｓ)、ポアソン比ν＝0.4、密度ρ＝2.0(ｔ／ｍ３)、減衰定数ｈ＝0.02とした。なお、面外方向は単位厚さとし、面外粘性境界は未考慮とした。まず、周波数領域で(８)式で表される伝達境界マトリクス[Ｒ]及び(11)式で表される境界力ベクトル{Ｆ(ω)}を演算によって求めた。この場合、伝達境界マトリクスは水平と上下の自由度から成る６行×６列のマトリクスとなる。なお、算定する振動数の範囲は何れも０〜20Ｈzとした。
【００６９】
次に、Ａ法、Ｂ法及びＣ法を各々適用すると共に、時間領域への変換（インパルス応答の演算）の条件として先の表１に示す条件を適用し、表１に示す各振動数における伝達境界マトリクスの値から各演算方法によってインパルス応答を各々演算した。更に、上記演算によって得られたインパルス応答を周波数領域の値へ再変換し、元のデータ（伝達境界マトリクスの元の値）との偏差に基づいてインパルス応答の演算精度を各演算方法毎に評価した。なお、周波数領域の値への再変換には、インパルス応答の各成分のうち、同時成分ｋ０,ｃ０(,ｍ０)と時間遅れ成分ｋ１〜ｋ１０,ｃ１〜ｃ１０を用いた。これにより、インパルス応答の演算方法としてＣ法を適用した場合は、前出の図２に示すように、ほぼ全てのデータ点において、再現値が元のデータと一致するか、又は元のデータに非常に近い値を示す結果が得られた。
【００７０】
比較のため、Ｂ法を適用した場合及びＣ法を適用した場合のマトリクス要素Ｋ55についての再現値を図７に各々示す。図７からも明らかなように、Ｃ法を適用した場合は全データ点において再現値が元のデータに非常に近い値を示しているのに対し、Ｂ法を適用した場合、実部については元のデータに対して再現値が若干相違している。以上の結果から、Ｂ法とＣ法の結果の差異はＣ法における同時成分ｍ０,ｋ０の修正に起因しており、Ｃ法における同時成分ｍ０,ｋ０の修正が精度向上に有効であることが理解できる。
【００７１】
(２)地震応答解析
次に、本願発明者は地盤−建物連成系モデルに地震動が入力された場合の応答解析を行った。この応答解析では、解析対象モデルとして、前述した物性の地盤に埋込みを有する建物を付加した、図６に示すモデルを用いた。建物は20ｍ×20ｍの底面を有する鉄筋コンクリート造の地上３階地下１階の構造物とし、地下部分は剛体と見做した。建物の物性を次の表２に示す。
【００７２】
【表２】
【００７３】
エネルギー伝達境界マトリクスは図２と同一であるが、本解析検討では面外方向の厚さを20ｍとしたため、マトリクス各要素の値も20倍して用いた。また、入力地震動としてはEl Centro 1940 ＮＳ波（Δｔ＝0.005秒、継続時間10秒、最大加速度323.4Ｇａｌ）を用いた。
【００７４】
(2-1)材料減衰モデルの検討
本願発明者は、応答解析の実施に先立ち、まず図６の解析対象モデルのうちの内部領域モデル（外端は自由境界とする）を用いて、材料減衰の振動数非依存性を表す減衰モデル（時刻歴応答解析に用いる減衰モデル）の検討を行った。この検討では、剛性比例型減衰モデル、歪エネルギー比例型減衰モデル、因果的履歴減衰モデルを検討対象とした。複素減衰を用いた周波数領域の応答解析結果を基準として、検討対象の各減衰モデルを用いた時刻歴応答解析の結果を比較した。なお、剛性比例型減衰モデルについては、全体１次振動数（4.92Ｈz）で、地盤の減衰２％、建物の減衰３％となるように設定した。各減衰モデルを用いた応答解析の、中心軸上の地盤・建物節点（節点番号１〜７）及び外端部節点（節点番号34〜37）の最大応答水平加速度を次の表３に示す。
【００７５】
【表３】
【００７６】
上記の表３において、着色された欄の数字は周波数領域の応答解析結果に対する比率を示す。なお、この応答解析では、周波数領域の応答解析で高精度な解析結果が得られる線形の応答解析を行っているので、周波数領域の応答解析結果を基準として用いることができる。表３に示す結果より明らかなように、周波数領域の応答解析結果を基準にすると、剛性比例型減衰モデルを用いた応答解析結果は中心軸・外端部とも差異が最大８％程度でやや精度が低いのに対して、歪エネルギー比例型減衰モデルを用いた応答解析結果及び因果的履歴減衰モデルを用いた応答解析結果は、中心軸・外端部とも差異の最大が各々３％、４％であり、ともに概ね良好な値が得られている。これより、因果的履歴減衰モデルは、４角形要素等を用いた解析対象モデルによる時刻歴応答解析において、歪エネルギー比例型減衰モデルとほぼ同程度の良好な解析精度が得られることが明らかとなった。
【００７７】
(2-2)時間領域変換した伝達境界マトリクスを用いた応答解析
続いて本願発明者は、(2-1)では自由端としていたモデルの外端部に、波動境界モデルとしてのエネルギー伝達境界を接続した解析対象モデル(図６に示す解析対象モデル)を用い、本発明を適用して伝達境界マトリクスを時間領域へ変換することで得られたインパルス応答を用いて時刻歴応答解析を行った結果を周波数領域の応答解析結果と比較し、本発明の有効性を検討した。この検討において、内部減衰モデルとしては、(2-1)における検討結果に基づき歪エネルギー比例型減衰モデルと因果的履歴減衰モデルを用いた。本検討は、図６に示す解析対象モデルを用いることで、モデルの内部領域がエネルギー伝達境界を介して自由地盤の影響を受ける点で (2-1)と相違している。また、時刻歴応答解析には、伝達境界マトリクスを時間領域へ変換することで得られたインパルス応答の成分のうち、同時成分ｋ0,ｃ0,ｍ0及び時間遅れ成分ｋ1〜ｋ10,ｃ1〜ｃ10を用いた。これは図２に示した再現値の算出に用いた成分と同一である。各応答解析によって得られた、(2-1)と同じく中心軸及び外端部の各節点における最大応答水平加速度を次の表４に示す。
【００７８】
【表４】
【００７９】
表４から明らかなように、因果的履歴減衰モデルを用いた時刻歴応答解析の結果は、周波数領域の応答解析結果に対して36番節点で最大約８％の差異が生じているものの、殆どの節点では差異が２％以内で、全体としては概ね良好に対応している。
【００８０】
これに対し、(2-1)では因果的履歴減衰モデルとほぼ同様の精度が得られていた歪エネルギー比例型減衰モデルを用いた時刻歴応答解析の結果は、周波数領域の応答解析結果に対して最大で14％の差異が生じていると共に、半数以上の節点で５％以上の差異が生じており、精度が低下している。この精度低下は、(2-1)と同じく外端部自由の条件で内部領域の固有値解析を行い、全固有モードを用いて減衰マトリクスを算出したのに対し、内部領域に伝達境界マトリクスを接続したことで減衰マトリクス算出時と振動系が変化したことが原因と思われ、修正可能である。但し、歪エネルギー比例型減衰モデルは、全次数の固有値解析を必要とし、フルマトリクス減衰を用いる必要があるので、３次元大規模解析では計算負荷が非常に大きく有効性が低い。
【００８１】
一方、因果的履歴減衰モデルは、時刻歴応答解析に際して過去の変位や速度を記憶する必要はあるものの、固有値解析とフルマトリクス減衰はともに不要であり、計算負荷は小さくて済む。従って、本発明を適用して時刻歴応答解析を行う場合の減衰モデルとしては、因果的履歴減衰モデルが最適と考えられる。
【００８２】
(2-3)粘性境界を用いた応答解析結果との比較検討
また、本願発明者は、代表的な時間領域の波動境界モデルである粘性境界と比較することで、波動境界モデルとしてエネルギー伝達境界を用いたモデルの有効性を確認する比較検討を行った。本比較検討で用いる粘性境界を設けた解析対象モデルを図８に示す。なお、内部減衰モデルとしては因果的履歴減衰モデルを用いた。図８に示す解析対象モデルは、基礎幅(＝20ｍ)に対し、内部領域の幅を(Ａ)は1.5倍、(Ｂ)は３倍、(Ｃ)は５倍と変化させている。１／２逆対称モデルを用いているため、図８における内部領域の長さは上記の１／２となっている。このうち図８(Ａ)に示すモデルは、波動境界モデルとしてエネルギー伝達境界を設けた解析対象モデル（図６参照）と同一形状である。図８に示す各モデルの中心軸位置の最大応答水平加速度を、エネルギー伝達境界を用いた周波数域の解析結果と比較した結果を次の表５に示す。
【００８３】
【表５】
【００８４】
表５より明らかなように、図８(Ａ)に示す解析対象モデルを用いた時刻歴応答解析の結果は、周波数域の応答解析の結果に対し、建物頂部（節点番号１）で16％、基礎下端（節点番号4）で40％以上の差異を生じており、精度が非常に低い。これに対して、図８(Ｂ)の解析対象モデル、図８(Ｃ)の解析対象モデルとモデル化領域が大きくなるに従って差異が小さくなる傾向を示している。しかしながら、基礎幅の５倍のモデル化領域を有する図８(Ｃ)の解析対象モデルでも、基礎下端以下で周波数域の応答解析の結果に対して20％以上の差異が生じている。上記の結果を表４と比較しても明らかなように、波動境界モデルとしてエネルギー伝達境界を用いた解析対象モデルは、波動境界モデルとして粘性境界モデルを用いた解析対象モデルよりも明らかに高精度での時刻歴応答解析を実現することができる。
【００８５】
〔第２実施例〕
エネルギー伝達境界は２次元面内問題以外に２次元面外問題、軸対称問題についても定式化されている。このため、本願発明者は２次元面外問題、軸対称問題についても、エネルギー伝達境界の伝達境界マトリクスが、本発明によって精度良く時間領域のインパルス応答に変換されるか否かを確認する解析検討を行った。この解析検討に用いた地盤モデルを図９に、この地盤モデルにおける定数を表６に、時間領域への変換の条件を表７に各々示す。
【００８６】
【表６】
【００８７】
【表７】
【００８８】
また、時間領域への変換（インパルス応答の演算）にはＣ法を用い、算出されたインパルス応答の各成分のうち同時成分ｋ0,ｃ0,ｍ0及び時間遅れ成分ｋ1〜ｋ10,ｃ1〜ｃ10を用いて伝達境界マトリクスを再現した。
【００８９】
(1)２次元面外問題
２次元面外問題では１節点当りの解析自由度が面外方向の１自由度であり、伝達境界マトリクスは３行×３列となる。算出されたインパルス応答から再現された伝達境界マトリクスを、元のデータ点と共に図１０に示す。図１０からも明らかなように、再現された伝達境界マトリクスはほぼ全てのデータ点と良好に対応しており、本発明が２次元面外問題にも有効であることが確認された。
【００９０】
(2)軸対称問題
軸対称問題では１節点当りの解析自由度が３自由度であり、伝達境界マトリクスは９行×９列となる。算出されたインパルス応答から再現された９行×９列の伝達境界マトリクスのうちの１〜６行目×１〜６列目を、元のデータ点と共に図１１に示す。図１１においても、再現された伝達境界マトリクスはほぼ全てのデータ点と良好に対応しており、本発明が軸対称問題にも有効であることが確認された。
【図面の簡単な説明】
【００９１】
【図１】地盤モデル(１／２モデル)の一例を示す概念図である。
【図２】本願発明者による解析検討において、２次元面内問題で得られたインパルス応答から再現した伝達境界マトリクスを元のデータ点と共に示す線図である。
【図３】本実施形態に係るコンピュータの概略構成を示すブロック図である。
【図４】地震応答解析処理の内容を示すフローチャートである。
【図５】図４の地震応答解析処理における各ステップの関係（データの流れ）を示す概念図である。
【図６】解析対象モデルの一例を示す概念図である。
【図７】本願発明者による解析検討において、(Ａ)はＢ法を、(Ｂ)はＣ法を適用して各々演算したインパルス応答から再現した伝達境界マトリクスを元のデータ点と共に示す線図である。
【図８】本願発明者による解析検討に用いた、粘性境界モデルを設けた解析対象モデルを示す概念図である。
【図９】本願発明者による解析検討に用いた地盤モデルを示す概念図である。
【図１０】本願発明者による解析検討において、２次元面外問題で得られたインパルス応答から再現した伝達境界マトリクスを元のデータ点と共に示す線図である。
【図１１】本願発明者による解析検討において、軸対称問題で得られたインパルス応答から再現した伝達境界マトリクスを元のデータ点と共に示す線図である。
【符号の説明】
【００９２】
１０ コンピュータ
１０Ｃ 記憶部
１２ ディスプレイ
１４ キーボード
１６ マウス
【技術分野】
【０００１】
本発明は応答解析装置、方法及びプログラムに係り、特に、波動境界モデルを含む解析対象モデルに対して時間領域での振動問題の応答解析を行う応答解析装置、該応答解析装置に適用可能な応答解析方法、及び、コンピュータを前記応答解析装置として機能させるための応答解析プログラムに関する。
【背景技術】
【０００２】
地震時の建物の挙動や損傷を精度良く評価するためには、建物と地盤との相互作用効果を考慮して地震応答解析を行う必要がある。また地震応答解析は、周波数領域で応答解析を行う周波数応答解析と、時間領域で応答解析を行う時刻歴応答解析とに大別されるが、大地震時には地盤・建物とも大きく塑性化するので、建物の挙動等を精度良く評価するためには、非線形の解析が可能な時刻歴応答解析を行うことが望ましい。これに関連して本願発明者は、時刻歴応答解析の実現を目的として、周波数領域の複素関数である地盤の動的剛性を時間領域で表されるインパルス応答へ変換する技術を既に提案している（例えば非特許文献１や特許文献１を参照）。
【０００３】
ところで、地震応答解析では地盤を有限要素法(ＦＥＭ)の非線形ソリッド要素を用いてモデル化することが多い。地盤は本来半無限的な広がりを持つのに対し、ＦＥＭ解析では解析領域を有限としてモデル化する必要があるため、解析対象モデルの側面や底面には波動境界モデルが設定される。建物で励起された波動は主として、下方へは実体波として、側方へは表面波として伝播するが、表面波の方がより遠方まで伝播するため、解析対象モデルの側面及び底面のうち、側面に設定する波動境界モデルが解析結果により大きな影響を及ぼす。時間領域の解析に適用可能な側面境界モデルとしては、例えば繰返し境界、粘性境界、無反射境界が知られている。このうち、繰返し境界は、モデル化は極めて容易であるものの波動境界としての精度が低く、モデル化領域(解析対象モデル)を大きくとる必要があるので解析負荷が大きいという問題がある。また、粘性境界は斜め方向からの波動に対しての精度が低く、解析精度を上げるため同様にモデル化領域を大きくする必要がある。更に、無反射境界は現実の問題への適用に課題が多く、あまり用いられていない。
【０００４】
一方、より高精度な波動境界モデルとしてエネルギー伝達境界(Energy Transmitting Boundary、又は、Consistent Transmitting Boundary)が知られている（例えば非特許文献２を参照）。エネルギー伝達境界は、剛基盤上に平行成層をなす地盤の外端に設置され、水平方向には厳密で、上下方向は要素の変位仮定に従う（線形一次要素であれば、変位は直線的に変化する）高精度の境界であり、任意方向からの波動をほぼ完全に吸収する特性を有している。このため、解析対象モデルの側面の波動境界モデルとして上記のエネルギー伝達境界を用いれば、モデル化領域を大幅に小さくすることで、同等の解析精度での応答解析を、より小さな解析負荷で実現できるものと期待される。
【特許文献１】特許３８７８６２６号公報
【非特許文献１】中村尚弘，「地盤インピーダンスの時間領域変換による成層地盤に埋込まれた構造物の地震応答解析」，日本建築学会構造系論文集，２００３年５月，第５６７号，ｐ．６３−７０
【非特許文献２】川本眺万，林正夫，「地盤工学における有限要素解析 土質力学と岩盤力学へのアプローチ」，培風館，１９７８年，ｐ．３４１−３６０
【発明の開示】
【発明が解決しようとする課題】
【０００５】
しかしながら、エネルギー伝達境界の特性は、伝達境界マトリクスと補正力ベクトルによって表されるが、このうち伝達境界マトリクスは、周波数領域の複素数であり強い振動数依存性を有しているので、時間領域に変換することは困難であった。このため、側面境界モデルとしてエネルギー伝達境界を適用した解析対象モデルによる応答解析は、周波数領域での等価線形解析に限られているのが実情であった。
【０００６】
本発明は上記事実を考慮して成されたもので、波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現する応答解析装置、応答解析方法及び応答解析プログラムを得ることが目的である。
【課題を解決するための手段】
【０００７】
以下にエネルギー伝達境界の概要を説明する。地盤をモデル化すると共に、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界を設け、モデルの中心位置に関して逆対称条件を設けることにより１／２モデル化する。右側の１／２モデルを図１に示す。図１に示す１／２モデルでは、モデルの右側にのみエネルギー伝達境界が存在し、その左側には内部領域が位置し、その右側には外部領域(自由地盤ともいう)が位置している。この１／２モデルにおいて、２次元面内問題での自由地盤の波動方程式は次の(１)式で表される。
([Ａ]ｋ２＋i[Ｂ]ｋ＋[Ｇ]−ω２[Ｍ]){ｕ＊(ω)}＝{0} …(１)
ここで、ｋは波数、{ｕ＊(ω)}は自由地盤の変位ベクトルである。また、[Ａ],[Ｂ],[Ｇ],[Ｍ]は２ｎ×２ｎ（ｎは節点数）のマトリクスで、以下の要素毎のサブマトリクス[Ａ]j,[Ｂ]j,[Ｇ]j,[Ｍ]jの重ね合せで表される。各サブマトリクスは、各要素のラメ定数Ｇj,λjを用いて次の(２)〜(５)式で表される。ｈj,ρjは各要素の高さと密度である。
【０００８】
【数１】
【０００９】
ここで、
[Ｃ]＝[Ｇ]−ω２[Ｍ] …(６)
とし、次の(７)式の波数ｋに関する固有値問題を解く。
|[Ａ]ｋ２＋i[Ｂ]ｋ＋[Ｃ]|＝{0} …(７)
得られた４ｎ個の固有モードより、右方向に伝播する成分２ｎ個を抽出し、モードマトリクス[Ｖ]とする。これより、１／２モデルの右側に存在するエネルギー伝達境界の伝達境界マトリクスは次の(８)式の[Ｒ]で表される。但し、[Ｄ]は下記の(９)式で表されるサブマトリクス[Ｄ]jの重ね合せで表される。
[Ｒ]＝i・[Ａ][Ｖ][Ｋ][Ｖ]−１＋[Ｄ] …(８)
【００１０】
【数２】
【００１１】
なお、先の(８)式ではマトリクスが非対称となるが、以下では対称位置の項を平均することによりマトリクスを対称化する。これよりモデル全体の運動方程式は次の(10)式で表される。
(−ω２[ＭＩ]＋[ＫＩ]＋[Ｒ]){ｕ(ω)}＝−ｙ"(ω)[ＭＩ]{1}＋{Ｆ(ω)} …(10)
但し、 [ＭＩ],[ＫＩ],{ｕ(ω)}は各々、内部領域の質量マトリクス、剛性マトリクス、変位ベクトルである。また{Ｆ(ω)}は境界力ベクトルであり、次の(11)式で表せる。
{Ｆ(ω)}＝([Ｒ]−[Ｑ]){ｕ＊(ω)} …(11)
上記の(11)式のうち、[Ｑ]{ｕ＊(ω)}は補正力（切欠き力ともいう）ベクトルであり、地震動が鉛直下方より入射する場合は[Ｑ]＝[Ｄ]、それ以外の場合は[Ｑ]＝ｉｋ[Ａ]＋[Ｄ]となる。但しｋは波数である。また、[Ｒ],{ｕ＊(ω)},{Ｆ(ω)}は、(８)式までは境界部の自由度数２ｎで定義したが、内部領域の自由度に重ね合わせるため、(10)式以降では内部領域の自由度数に拡張して用いている。
【００１２】
ここで、上述したエネルギー伝達境界の伝達境界マトリクスは周波数領域の複素数であり、強い振動数依存性を有しているが、本願発明者はエネルギー伝達境界の伝達境界マトリクスについて検討・考察を重ね、エネルギー伝達境界の伝達境界マトリクスが、本質的には、同じく周波数領域の複素数で表され強い振動数依存性を有している地盤等の物体の動的剛性と類似する物理量であることに想到した。そして本願発明者は、伝達境界マトリクスが本質的には物体の動的剛性と類似する物理量であるならば、非特許文献１や特許文献１で提案している技術を適用することで、時間領域のインパルス応答へ変換できる可能性があることに想到し、当該事項を確認するために以下の解析検討を行った。
【００１３】
この解析検討では、図１に示した地盤モデル(１／２モデル)を用い、演算に用いる地盤の物性値を、剪断波速度Ｖｓ＝400(ｍ／ｓ)、ポアソン比ν＝0.4、密度ρ＝2.0(ｔ／ｍ３)、減衰定数ｈ＝0.02とした。また、伝達境界マトリクスから時間領域のインパルス応答の演算には、本願発明者が非特許文献１で提案した演算方法（便宜的にＡ法と称する）、本願発明者が特許文献１で提案した演算方法（便宜的にＢ法と称する）、当該Ｂ法に対し後述する請求項７記載の発明に従って同時成分ｍ０,ｋ０の修正を更に行う演算方法（便宜的にＣ法と称する）を各々用いると共に、時間領域への変換（インパルス応答の演算）の条件として次の表１に示す条件を適用し、表１に示す各振動数における伝達境界マトリクスの値から各演算方法によってインパルス応答を各々演算した。
【００１４】
【表１】
【００１５】
また、上記演算の精度を確認するために、上記演算によって得られたインパルス応答を周波数領域の値へ再変換し、元のデータ（伝達境界マトリクスの元の値）との偏差に基づいてインパルス応答の演算精度を各演算方法毎に評価した。なお、周波数領域の値への再変換には、インパルス応答の各成分のうち、同時成分ｋ０,ｃ０(,ｍ０)と時間遅れ成分ｋ１〜ｋ１０,ｃ１〜ｃ１０を用いた。インパルス応答の演算方法としてＣ法を適用した場合の結果を図２に示す。図２に示すように、インパルス応答の演算方法としてＣ法を適用した場合の再現値（インパルス応答を周波数領域の値へ再変換して得られる値）は、ほぼ全てのデータ点において、元のデータと一致しているか、又は元のデータとの偏差が非常に小さく、エネルギー伝達境界の伝達境界マトリクスからエネルギー伝達境界のインパルス応答を精度良く演算できることが明らかとなった。また、図示は省略するが、Ａ法及びＢ法についても、演算精度はＣ法＞Ｂ法＞Ａ法となるものの、エネルギー伝達境界の伝達境界マトリクスから、応答解析に用いるのに十分な精度でエネルギー伝達境界のインパルス応答を演算できることが確認された。
【００１６】
上記に基づき請求項１記載の発明に係る応答解析装置は、地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算するマトリクス演算手段と、物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算するインパルス応答演算手段と、前記インパルス応答演算手段によって演算された前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行う解析手段と、を含んで構成されている。
【００１７】
請求項１記載の発明では、地盤単独又は地盤と建物から成る解析対象物体に対応し、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルを解析対象としており、マトリクス演算手段は、エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算する。また、インパルス応答演算手段は、物体を振動させる外力と物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式（すなわち、前出のＡ法〜Ｃ法の何れかに対応する演算式）を用い、振動が各周波数のときの伝達境界マトリクスの値に基づいてエネルギー伝達境界のインパルス応答を演算する。これにより、先にも説明したように、エネルギー伝達境界のインパルス応答を、応答解析に用いるのに十分な精度で得ることができる。そして解析手段は、インパルス応答演算手段によって演算されたエネルギー伝達境界のインパルス応答を用いて、解析対象モデルに対して時間領域での応答解析を行う。これにより、波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現することができる。
【００１８】
なお、請求項１記載の発明におけるエネルギー伝達境界の伝達境界マトリクスは、先に説明し、かつ請求項２にも記載したように、内部領域の質量マトリクスを[ＭＩ]、剛性マトリクスを[ＫＩ]、変位ベクトルを{ｕ(ω)}、境界力ベクトルを{Ｆ(ω)}、補正力ベクトルを[Ｑ]{ｕ＊(ω)}とし、
(−ω２[ＭＩ]＋[ＫＩ]＋[Ｒ]){ｕ(ω)}＝−ｙ"(ω)[ＭＩ]{1}＋{Ｆ(ω)} …(10)
{Ｆ(ω)}＝([Ｒ]−[Ｑ]){ｕ＊(ω)} …(11)
解析対象モデル全体の運動方程式を上記の(10),(11)式で表したときに、上記の(10),(11)式の[Ｒ]で表すことができる。
【００１９】
また、詳細は後述するが、本願発明者は、本発明が２次元面内問題、２次元面外問題及び軸対称問題の何れにも有効であることを確認しており、請求項１記載の発明に係る解析対象モデルに前記波動境界モデルとして設けられているエネルギー伝達境界としては、請求項３に記載したように、２次元面内問題、２次元面外問題、及び、軸対称問題の何れかに対応するモデルを適用することができ、２次元面内問題、２次元面外問題及び軸対称問題の何れに対応するモデルであってもよい。
【００２０】
また、請求項１記載の発明において、インパルス応答演算手段は、例えば請求項４に記載したように、物体の変位に依存するインパルス応答の同時成分をｋ０、物体の速度に依存するインパルス応答の同時成分をｃ０、物体の変位に依存するインパルス応答のΔt刻みの時間遅れ成分をｋj、物体の速度に依存するインパルス応答のΔt刻みの時間遅れ成分をｃj（但しｊは自然数でtj＝Δt・j）、時間領域での物体の変位をu(t)、速度をu'(t)、加速度をu"(t)としたときに、インパルス応答ＦＡ(t)を規定する数式として、
【００２１】
【数３】
【００２２】
上記(12)式を用い、前記振動がＮ種（Ｎ＝ｎ＋１）の周波数のときの伝達境界マトリクスの値に基づいてインパルス応答を求めるように構成することができる。上記の(12)式は前出のＡ法に対応しており、エネルギー伝達境界のインパルス応答を、応答解析に用いるのに十分な精度で得ることができる。
【００２３】
また、請求項１記載の発明において、インパルス応答演算手段は、例えば請求項５に記載したように、インパルス応答を規定する数式として、剛性項及び減衰項に加えて、物体の加速度に依存し少なくとも同時成分を含んで成る質量項も含む数式を用いて、エネルギー伝達境界のインパルス応答を演算するように構成してもよい。
【００２４】
より詳しくは、インパルス応答演算手段は、例えば請求項６に記載したように、物体の変位に依存するインパルス応答の同時成分をｋ０、物体の速度に依存するインパルス応答の同時成分をｃ０、物体の加速度に依存するインパルス応答の同時成分をｍ０、物体の変位に依存するインパルス応答のΔt刻みの時間遅れ成分をｋj、物体の速度に依存するインパルス応答のΔt刻みの時間遅れ成分をｃj（但しｊは自然数でtj＝Δt・j）、時間領域での物体の変位をu(t)、速度をu'(t)、加速度をu"(t)としたときに、前記質量項も含む数式として、前記インパルス応答ＦB(t)を規定する数式である、
【００２５】
【数４】
【００２６】
上記(13)式を用い、前記振動がＮ種（Ｎ＝ｎ＋１）の周波数のときの伝達境界マトリクスの値に基づいてインパルス応答を演算するように構成してもよい。上記の(13)式は前出のＢ法に対応しており、エネルギー伝達境界のインパルス応答をＡ法によって演算する場合（インパルス応答演算手段を請求項４記載の構成とする場合）よりも、エネルギー伝達境界のインパルス応答としてより高精度な値を得ることができ、解析対象モデルに対する時間領域での応答解析の解析精度を向上させることができる。
【００２７】
また、請求項６記載の発明において、インパルス応答演算手段は、例えば請求項７に記載したように、前記(13)式における物体の加速度に依存するインパルス応答の同時成分ｍ０に対する修正値Δｍ0、物体の変位に依存するインパルス応答の同時成分ｋ０に対する修正値Δｋ0を下記の(14)式によって演算し、
【００２８】
【数５】
【００２９】
（但し、上記(14)式において、
【００３０】
【数６】
【００３１】
であり、上記(15)式におけるＲe(ＳＢ(ωｉ))は、前記(13)式を用いて演算した前記エネルギー伝達境界の前記インパルス応答から再現された、下記の(16)式で表される前記伝達境界マトリクスＳＢ(ω)の振動数ωｉでの実部の値を表し、
【００３２】
【数７】
【００３３】
上記(15)式におけるＲe(Ｄ(ωｉ))は、(13)式に基づくエネルギー伝達境界のインパルス応答の演算に用いた、振動数ωｉでの伝達境界マトリクスのデータＤ(ωｉ)のうちの実部の値を表す）、演算した修正値Δｍ0,Δｋ0を用いて同時成分ｍ０,ｋ０を修正することが好ましい。
【００３４】
本願発明者による前述の解析検討において、Ｂ法を適用して算出したエネルギー伝達境界のインパルス応答を上記の(16)式によって周波数領域の値へ再変換した場合、虚部については元のデータと精度良く一致しているものの、実部については元のデータとの偏差が比較的大きい、という結果となった(Ａ法についても同様)。これは、エネルギー伝達境界の伝達境界マトリクスが因果的ではないので時間領域への厳密な変換が不可能であり、Ｂ法によるインパルス応答の演算が近似的な演算であることが原因と推察される。このため、本願発明者は、Ｂ法を適用して算出したエネルギー伝達境界のインパルス応答のうち、実部の同時成分であるｍ０,ｋ０に対し、先の(14),(15)式に従って修正値Δｍ0,Δｋ0を演算し、演算した修正値Δｍ0,Δｋ0を用いて同時成分ｍ０,ｋ０を修正した（Ｃ法）後に、Ｃ法によって得られたエネルギー伝達境界のインパルス応答を周波数領域の値へ再変換し、元のデータと比較する解析検討を行った。その結果、前出の図２にも示すように、実部についても元のデータと精度良く一致していることが確認された。
【００３５】
請求項７記載の発明は、先の(14),(15)式に従って修正値Δｍ0,Δｋ0を演算し、演算した修正値Δｍ0,Δｋ0を用いて同時成分ｍ０,ｋ０を修正するので、上述した解析検討の結果からも明らかなように、エネルギー伝達境界のインパルス応答をＡ法やＢ法によって演算する場合（インパルス応答演算手段を請求項４又は請求項６記載の構成とする場合）よりも、エネルギー伝達境界のインパルス応答として更に高精度な値を得ることができ、解析対象モデルに対する時間領域での応答解析の解析精度を更に向上させることができる。
【００３６】
また、請求項７記載の発明において、インパルス応答演算手段は、例えば請求項８に記載したように、前記(13)式における物体の速度に依存するインパルス応答の同時成分をｃ０に対する修正値Δｃ0も下記の(17)式によって演算し、
Δｃ0＝−Ｅ／Ｂ …(17)
（但し、上記(17)式において、
【００３７】
【数８】
【００３８】
であり、上記(18)式におけるＩm(ＳＢ(ωｉ))は、前記(13)式を用いて演算したエネルギー伝達境界のインパルス応答から再現された、前記(16)式で表される伝達境界マトリクスＳＢ(ω)の振動数ωｉでの虚部の値を表し、上記(17)式におけるＩm(ＳＢ(ωｉ))は、(13)式に基づくエネルギー伝達境界のインパルス応答の演算に用いた、振動数ωｉでの伝達境界マトリクスのデータＤ(ωｉ)のうちの虚部の値を表す）、演算した修正値Δｃ0を用いて同時成分ｃ０も修正することが好ましい。
【００３９】
本願発明者による前述の解析検討では、Ｃ法によるエネルギー伝達境界のインパルス応答の演算を様々な条件で行ったが、この解析検討により、例えばインパルス応答の演算に用いるデータ点の振動数軸上における位置が等間隔でない等の場合には、演算したインパルス応答のうちの虚部についての精度が低下する（演算したインパルス応答を周波数領域の値へ再変換した場合の元のデータとの偏差が虚部について比較的大きくなる）ことが判明した。このため、本願発明者はＣ法を適用して算出したエネルギー伝達境界のインパルス応答のうち、虚部についての精度が低下していたインパルス応答に対し、虚部の同時成分であるｃ０に対し、先の(18)式に従って修正値Δｃ0を演算し、演算した修正値Δｃ0を用いて同時成分ｃ０を修正した後に、当該修正を経たエネルギー伝達境界のインパルス応答を周波数領域の値へ再変換し、元のデータと比較する解析検討を行った。その結果、虚部の精度が向上することが確認された。
【００４０】
請求項８記載の発明は、先の(18)式に従って修正値Δｃ0を演算し、演算した修正値Δｃ0を用いて同時成分ｃ０を修正するので、例えばインパルス応答の演算に用いるデータ点の振動数軸上における位置が等間隔でない等のように、演算したインパルス応答のうちの虚部についての精度が低下し易い条件であっても、エネルギー伝達境界のインパルス応答として更に高精度な値を得ることができ、解析対象モデルに対する時間領域での応答解析の解析精度を更に向上させることができる。なお、請求項８記載の発明に係る同時成分ｃ0の修正は常時行う必要はなく、例えばインパルス応答の演算に用いるデータ点の振動数軸上における位置が等間隔でない等のように、インパルス応答のうちの虚部についての精度が低下することが予想できる条件のときにのみ行うようにしてもよい。
【００４１】
請求項９記載の発明に係る応答解析方法は、コンピュータに、地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算させ、物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算させ、前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行わせるので、請求項１記載の発明と同様に、波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現することができる。
【００４２】
請求項１０記載の発明に係る応答解析プログラムは、コンピュータを、地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算するマトリクス演算手段、物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算するインパルス応答演算手段、及び、前記インパルス応答演算手段によって演算された前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行う解析手段として機能させる。
【００４３】
請求項１０記載の発明に係る応答解析プログラムは、コンピュータを、上記のマトリクス演算手段、インパルス応答演算手段及び解析手段として機能させるためのプログラムであるので、コンピュータが請求項１０記載の発明に係る応答解析プログラムを実行することで、コンピュータが請求項１０に記載の応答解析装置として機能することになり、請求項１記載の発明と同様に、波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現することができる。
【発明の効果】
【００４４】
以上説明したように本発明は、波動境界モデルとして設けられたエネルギー伝達境界の伝達境界マトリクスを演算し、インパルス応答を規定する数式として、少なくとも、物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの伝達境界マトリクスの値に基づいてエネルギー伝達境界のインパルス応答を演算し、エネルギー伝達境界のインパルス応答を用いて、解析対象モデルに対して時間領域での応答解析を行うので、波動境界モデルとしてエネルギー伝達境界が設けられた解析対象モデルに対する時間領域での応答解析を実現できる、という優れた効果を有する。
【発明を実施するための最良の形態】
【００４５】
以下、図面を参照して本発明の実施形態の一例を詳細に説明する。図３には本発明を適用可能なコンピュータ１０が示されている。コンピュータ１０は、ＣＰＵ１０Ａ、ＲＯＭやＲＡＭ等から成るメモリ１０Ｂ、ＨＤＤ(Hard Disk Drive)やフラッシュメモリ等から成る不揮発性の記憶部１０Ｃを備えており、ＣＲＴ又はＬＣＤから成るディスプレイ１２、キーボード１４、マウス１６が各々接続されている。
【００４６】
コンピュータ１０の記憶部１０Ｃには、後述する地震応答解析処理を行うための地震応答解析プログラムがインストールされている。この地震応答解析プログラムは、請求項１０に記載の応答解析プログラムに対応している。また、コンピュータ１０は請求項９及び請求項１０に記載のコンピュータに対応しており、ＣＰＵ１０Ａが地震応答解析プログラムを実行することで、請求項１等に記載の発明に係る応答解析装置として機能する。なお、コンピュータ１０としてはパーソナル・コンピュータ(ＰＣ)が好適であるが、これに限られるものではなく、例えばワークステーションであってもよいし、汎用の大型コンピュータであってもよい。
【００４７】
次に本実施形態の作用として、ＣＰＵ１０Ａが地震応答解析プログラムを実行することで実現される地震応答解析処理について、図４(及び図５)を参照して説明する。なお、以下では図６に示すように、地盤と当該地盤に埋め込まれた地下部分を有する建物が一体となった解析対象物体をモデル化すると共に、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた構成の解析対象モデルに対して時間領域での応答解析を行う場合を例に説明する。
【００４８】
地震応答解析処理では、まずステップ３０において、内部地盤（解析対象モデルの内部領域）及び建物の動特性マトリクスを演算するためのデータを取得する。なお、本実施形態において、動特性マトリクスは、地盤や建物の応答解析を行うのに必要な特性を表すマトリクスを意味し、具体的には、質量マトリクス、減衰マトリクス及び剛性マトリクスから構成され、ステップ３０で取得されるデータとしては、例えば内部地盤の各メッシュ領域毎の物性値（例えば剪断波速度Ｖｓ、ポアソン比ν、密度ρ、減衰定数ｈ等）や、建物の各階毎の定数（質量、剪断剛性、回転慣性、減衰率ｈ等）等が挙げられる。次のステップ３２では、ステップ３０で取得したデータに基づき、有限要素法(ＦＥＭ)による解析を行うことで内部地盤及び建物の動特性マトリクスを演算する。これにより、前出の(10)式に示す内部領域の質量マトリクス[ＭＩ]、内部領域の剛性マトリクス[ＫＩ]、前出の(６)式に示す減衰マトリクス[Ｃ]が各々算出される。なお、本実施形態では材料減衰の振動数非依存性を表す減衰モデルとして因果的履歴減衰モデルを用いており、上記の減衰マトリクス[Ｃ]は因果的履歴減衰モデルを用いて算出される。
【００４９】
ステップ３４では入力地震動のデータを取得し、次のステップ３６では外部地盤(無限遠方地盤)の動特性マトリクス及び地震応答を演算するためのデータとして地盤各層の物性値（例えば剪断波速度Ｖｓ、ポアソン比ν、密度ρ、減衰定数ｈ等）を取得する。ステップ３８では、先のステップ３４で取得した入力地震動のデータと、先のステップ３６で取得した外部地盤の地震応答を演算するためのデータに基づいて、外部地盤の地震応答を演算する(図５も参照)。これにより、外部地盤の地震応答として、先の(１)式における自由地盤の変位ベクトル{ｕ＊(ω)}が算出される。またステップ４０では、ステップ３８の演算によって得られた外部地盤の地震応答に対して逆フーリエ変換を行うことで、時間領域における外部地盤の地震応答を演算する(図５も参照)。
【００５０】
次のステップ４２では、ステップ３８の演算によって得られた外部地盤の地震応答と先のステップ３６で取得したデータに基づいて、外部地盤の補正力ベクトルを演算する(図５も参照)。これにより、外部地盤の補正力ベクトルとして、先の(11)式における補正力ベクトル"[Ｑ]{ｕ＊(ω)}"が算出される。ステップ４４では、ステップ４２の演算によって得られた外部地盤の補正力ベクトルに対して逆フーリエ変換を行うことで、時間領域における外部地盤の補正力ベクトルを演算する。
【００５１】
またステップ４６では、先のステップ３６で取得した外部地盤の地震応答を演算するためのデータに基づいて、外部地盤の動特性マトリクスを演算する(図５も参照)。これにより、外部地盤の動特性マトリクスとして、先の(1)式におけるマトリクス[Ａ],[Ｂ],[Ｇ],[Ｍ]やモードマトリクス[Ｖ]、先の(９)式で表されるサブマトリクス[Ｄ]jの重ね合せで表されるマトリクス[Ｄ]が算出される。次のステップ４８では、ステップ４６の演算によって得られた外部地盤の動特性マトリクスを先の(８)式に代入することで、外部地盤を含むエネルギー伝達境界の伝達境界マトリクスを演算する(図５も参照)。これにより、外部地盤を含むエネルギー伝達境界の伝達境界マトリクスとして、先の(８)式で表される伝達境界マトリクス[Ｒ]が算出される。
【００５２】
ステップ５０では、ステップ４８の演算によって得られた外部地盤を含むエネルギー伝達境界の伝達境界マトリクスを時間領域へ変換することで、外部地盤を含むエネルギー伝達境界のインパルス応答を算出する。エネルギー伝達境界のインパルス応答の演算は、具体的には以下のようにして行うことができる。
【００５３】
すなわち、まずステップ４８の演算によって得られた伝達境界マトリクスのデータから、予め設定された演算対象の周波数範囲内のＮ種の周波数（Ｎ種の角振動数ω1,…,ωN）における伝達境界マトリクスの値を表すＮ個の複素データＳB(ω1),…,ＳB(ωN)を各々抽出する。なお、演算対象の周波数範囲としては、例えば０〜20(Hz)の範囲を適用することができる。次に、インパルス応答を規定する数式として先の(13)式を、伝達境界マトリクスを規定する数式として先の(16)式を用い、(13),(16)式に基づいて導出した次の(19),(20)式の連立方程式にＮ個の複素データＳ(ω1),…,Ｓ(ωN)を代入し、この連立方程式の解を求めることで、外部地盤を含むエネルギー伝達境界のインパルス応答を、予め設定されたΔｔ刻みで演算する。
【００５４】
【数９】
【００５５】
この演算により、外部地盤を含むエネルギー伝達境界のインパルス応答として、物体の変位に依存するインパルス応答の同時成分ｋ0、物体の速度に依存するインパルス応答の同時成分ｃ0、物体の加速度に依存するインパルス応答の同時成分ｍ0のデータが得られると共に、物体の変位に依存するインパルス応答の時間遅れ成分ｋjのデータがΔt刻みでｎ個（ｎ＝Ｎ−１）得られ、物体の速度に依存するインパルス応答の時間遅れ成分ｃjのデータがΔt刻みでｎ−１個得られる。
【００５６】
続いて、先の(14),(15)式に従って修正値Δｍ0,Δｋ0を演算し、演算した修正値Δｍ0,Δｋ0を用いてインパルス応答の同時成分ｍ０,ｋ０を修正する。更に、例えばインパルス応答の演算に用いるデータ点の振動数軸上における位置が等間隔でない等のように、演算したインパルス応答のうちの虚部についての精度が低下し易い条件である場合には、先の(18)式に従って修正値Δｃ0も演算し、演算した修正値Δｃ0を用いて同時成分ｃ０を修正する。上記の演算により、外部地盤を含むエネルギー伝達境界のインパルス応答が得られ、解析対象モデルに対して時間領域での応答解析を行うことが可能となる。
【００５７】
ステップ５２では、先のステップ３２の演算によって得られた内部地盤及び建物の動特性マトリクス、先のステップ４０の演算によって得られた時間領域における外部地盤の地震応答、先のステップ４４の演算によって得られた時間領域における外部地盤の補正力ベクトル、先のステップ５０の演算によって得られた時間領域における外部地盤を含むエネルギー伝達境界の動特性マトリクス（外部地盤を含むエネルギー伝達境界のインパルス応答）を結合することで、解析対象モデルに対して時間領域の応答解析を行うための時間領域の運動方程式(次の(22)式参照)を求める(図５も参照)。
[ＭI]{ｕ"(t)}＋[ＫI]{ｕ(t)}＋{Ｒ0}＝−ｙ"(t)[ＭI]{1}＋{Ｆ(t)}＋{Ｒf(t)} …(21)
但し、
【００５８】
【数１０】
【００５９】
そしてステップ５４では、先のステップ３４で取得した入力地震動のデータが表す各時刻における地震動をステップ５２の処理によって得られた運動方程式(上記の(22),(23)式)に順に代入し、解析対象モデルに対して時間領域の地震応答解析(時刻歴応答解析)を行う。これにより、内部地盤、建物及び外部地盤が一体となった解析対象物体に地震動が入力された場合の解析対象物体の挙動を、地盤・建物の塑性化も考慮して精度良く評価することができる。
【００６０】
なお、上記ではエネルギー伝達境界のインパルス応答を、先の(13),(16)式に基づいて導出した先の (19),(20)式の連立方程式に複素データを代入することで演算する態様を説明したが、本発明に係るＢ法はこれに限定されるものではなく、(19),(20)式から算出される時間遅れ成分ｋj,ｃjのうち最初のｎ’個の時間遅れ成分が有意な値であり、それ以降の値が無視できるのであれば、先の(13),(16)式におけるΣの上限をＮ−２からｎ’に低減することも可能である（Ａ法の演算式である(12)式におけるΣの上限も同様にＮ−１からｎ’に低減可能）。これを先の(13),(16)式に適用した場合、次の(23),(24)式が得られる。
【００６１】
【数１１】
【００６２】
なお、(23),(24)式におけるｎ’＜Ｎ−１である。この(23),(24)式に基づいて導出した連立方程式に複素データを代入することで、エネルギー伝達境界のインパルス応答を演算するようにしてもよい。
【００６３】
また、上記ではＣ法を適用してエネルギー伝達境界のインパルス応答を演算する態様を説明したが、本発明はこれに限定されるものではなく、Ｂ法やＡ法を適用してエネルギー伝達境界のインパルス応答を演算することも権利範囲に含むものである。
【００６４】
また、上記では本発明を２次元面内問題に適用した態様を説明したが、これに限定されるものではなく、後述する第２実施例からも明らかなように本発明は２次元面外問題や軸対称問題にも有効であり、本発明は２次元面外問題や軸対称問題に適用してもよいことは言うまでもない。
【００６５】
更に、上記では本発明に係る応答解析プログラムに対応する地震応答解析プログラムがコンピュータ１０の記憶部１０Ｃに予め記憶（インストール）されている態様を説明したが、本発明に係る応答解析プログラムは、ＣＤ−ＲＯＭやＤＶＤ−ＲＯＭ等の記録媒体に記録されている形態で提供することも可能である。
【実施例】
【００６６】
次に、本発明の効果を確認するために本願発明者が行った解析検討の詳細を説明する。
【００６７】
〔第１実施例〕
この解析検討では、伝達境界の時間領域化とそれを用いた時刻歴応答解析の妥当性及び有効性を確認することを目的として、まず、地盤モデルに対して周波数領域の伝達境界マトリクスを作成し、時間領域変換の性状を評価する。更に、前記地盤モデルが表す地盤に埋込まれた建物を想定し、建物−地盤連成系の時刻歴地震応答解析を行った。
【００６８】
(１)伝達境界マトリクスの時間領域変換
この解析検討では、先に説明したように、図１に示した地盤モデル(１／２モデル)を用い、演算に用いる地盤の物性値を、剪断波速度Ｖｓ＝400(ｍ／ｓ)、ポアソン比ν＝0.4、密度ρ＝2.0(ｔ／ｍ３)、減衰定数ｈ＝0.02とした。なお、面外方向は単位厚さとし、面外粘性境界は未考慮とした。まず、周波数領域で(８)式で表される伝達境界マトリクス[Ｒ]及び(11)式で表される境界力ベクトル{Ｆ(ω)}を演算によって求めた。この場合、伝達境界マトリクスは水平と上下の自由度から成る６行×６列のマトリクスとなる。なお、算定する振動数の範囲は何れも０〜20Ｈzとした。
【００６９】
次に、Ａ法、Ｂ法及びＣ法を各々適用すると共に、時間領域への変換（インパルス応答の演算）の条件として先の表１に示す条件を適用し、表１に示す各振動数における伝達境界マトリクスの値から各演算方法によってインパルス応答を各々演算した。更に、上記演算によって得られたインパルス応答を周波数領域の値へ再変換し、元のデータ（伝達境界マトリクスの元の値）との偏差に基づいてインパルス応答の演算精度を各演算方法毎に評価した。なお、周波数領域の値への再変換には、インパルス応答の各成分のうち、同時成分ｋ０,ｃ０(,ｍ０)と時間遅れ成分ｋ１〜ｋ１０,ｃ１〜ｃ１０を用いた。これにより、インパルス応答の演算方法としてＣ法を適用した場合は、前出の図２に示すように、ほぼ全てのデータ点において、再現値が元のデータと一致するか、又は元のデータに非常に近い値を示す結果が得られた。
【００７０】
比較のため、Ｂ法を適用した場合及びＣ法を適用した場合のマトリクス要素Ｋ55についての再現値を図７に各々示す。図７からも明らかなように、Ｃ法を適用した場合は全データ点において再現値が元のデータに非常に近い値を示しているのに対し、Ｂ法を適用した場合、実部については元のデータに対して再現値が若干相違している。以上の結果から、Ｂ法とＣ法の結果の差異はＣ法における同時成分ｍ０,ｋ０の修正に起因しており、Ｃ法における同時成分ｍ０,ｋ０の修正が精度向上に有効であることが理解できる。
【００７１】
(２)地震応答解析
次に、本願発明者は地盤−建物連成系モデルに地震動が入力された場合の応答解析を行った。この応答解析では、解析対象モデルとして、前述した物性の地盤に埋込みを有する建物を付加した、図６に示すモデルを用いた。建物は20ｍ×20ｍの底面を有する鉄筋コンクリート造の地上３階地下１階の構造物とし、地下部分は剛体と見做した。建物の物性を次の表２に示す。
【００７２】
【表２】
【００７３】
エネルギー伝達境界マトリクスは図２と同一であるが、本解析検討では面外方向の厚さを20ｍとしたため、マトリクス各要素の値も20倍して用いた。また、入力地震動としてはEl Centro 1940 ＮＳ波（Δｔ＝0.005秒、継続時間10秒、最大加速度323.4Ｇａｌ）を用いた。
【００７４】
(2-1)材料減衰モデルの検討
本願発明者は、応答解析の実施に先立ち、まず図６の解析対象モデルのうちの内部領域モデル（外端は自由境界とする）を用いて、材料減衰の振動数非依存性を表す減衰モデル（時刻歴応答解析に用いる減衰モデル）の検討を行った。この検討では、剛性比例型減衰モデル、歪エネルギー比例型減衰モデル、因果的履歴減衰モデルを検討対象とした。複素減衰を用いた周波数領域の応答解析結果を基準として、検討対象の各減衰モデルを用いた時刻歴応答解析の結果を比較した。なお、剛性比例型減衰モデルについては、全体１次振動数（4.92Ｈz）で、地盤の減衰２％、建物の減衰３％となるように設定した。各減衰モデルを用いた応答解析の、中心軸上の地盤・建物節点（節点番号１〜７）及び外端部節点（節点番号34〜37）の最大応答水平加速度を次の表３に示す。
【００７５】
【表３】
【００７６】
上記の表３において、着色された欄の数字は周波数領域の応答解析結果に対する比率を示す。なお、この応答解析では、周波数領域の応答解析で高精度な解析結果が得られる線形の応答解析を行っているので、周波数領域の応答解析結果を基準として用いることができる。表３に示す結果より明らかなように、周波数領域の応答解析結果を基準にすると、剛性比例型減衰モデルを用いた応答解析結果は中心軸・外端部とも差異が最大８％程度でやや精度が低いのに対して、歪エネルギー比例型減衰モデルを用いた応答解析結果及び因果的履歴減衰モデルを用いた応答解析結果は、中心軸・外端部とも差異の最大が各々３％、４％であり、ともに概ね良好な値が得られている。これより、因果的履歴減衰モデルは、４角形要素等を用いた解析対象モデルによる時刻歴応答解析において、歪エネルギー比例型減衰モデルとほぼ同程度の良好な解析精度が得られることが明らかとなった。
【００７７】
(2-2)時間領域変換した伝達境界マトリクスを用いた応答解析
続いて本願発明者は、(2-1)では自由端としていたモデルの外端部に、波動境界モデルとしてのエネルギー伝達境界を接続した解析対象モデル(図６に示す解析対象モデル)を用い、本発明を適用して伝達境界マトリクスを時間領域へ変換することで得られたインパルス応答を用いて時刻歴応答解析を行った結果を周波数領域の応答解析結果と比較し、本発明の有効性を検討した。この検討において、内部減衰モデルとしては、(2-1)における検討結果に基づき歪エネルギー比例型減衰モデルと因果的履歴減衰モデルを用いた。本検討は、図６に示す解析対象モデルを用いることで、モデルの内部領域がエネルギー伝達境界を介して自由地盤の影響を受ける点で (2-1)と相違している。また、時刻歴応答解析には、伝達境界マトリクスを時間領域へ変換することで得られたインパルス応答の成分のうち、同時成分ｋ0,ｃ0,ｍ0及び時間遅れ成分ｋ1〜ｋ10,ｃ1〜ｃ10を用いた。これは図２に示した再現値の算出に用いた成分と同一である。各応答解析によって得られた、(2-1)と同じく中心軸及び外端部の各節点における最大応答水平加速度を次の表４に示す。
【００７８】
【表４】
【００７９】
表４から明らかなように、因果的履歴減衰モデルを用いた時刻歴応答解析の結果は、周波数領域の応答解析結果に対して36番節点で最大約８％の差異が生じているものの、殆どの節点では差異が２％以内で、全体としては概ね良好に対応している。
【００８０】
これに対し、(2-1)では因果的履歴減衰モデルとほぼ同様の精度が得られていた歪エネルギー比例型減衰モデルを用いた時刻歴応答解析の結果は、周波数領域の応答解析結果に対して最大で14％の差異が生じていると共に、半数以上の節点で５％以上の差異が生じており、精度が低下している。この精度低下は、(2-1)と同じく外端部自由の条件で内部領域の固有値解析を行い、全固有モードを用いて減衰マトリクスを算出したのに対し、内部領域に伝達境界マトリクスを接続したことで減衰マトリクス算出時と振動系が変化したことが原因と思われ、修正可能である。但し、歪エネルギー比例型減衰モデルは、全次数の固有値解析を必要とし、フルマトリクス減衰を用いる必要があるので、３次元大規模解析では計算負荷が非常に大きく有効性が低い。
【００８１】
一方、因果的履歴減衰モデルは、時刻歴応答解析に際して過去の変位や速度を記憶する必要はあるものの、固有値解析とフルマトリクス減衰はともに不要であり、計算負荷は小さくて済む。従って、本発明を適用して時刻歴応答解析を行う場合の減衰モデルとしては、因果的履歴減衰モデルが最適と考えられる。
【００８２】
(2-3)粘性境界を用いた応答解析結果との比較検討
また、本願発明者は、代表的な時間領域の波動境界モデルである粘性境界と比較することで、波動境界モデルとしてエネルギー伝達境界を用いたモデルの有効性を確認する比較検討を行った。本比較検討で用いる粘性境界を設けた解析対象モデルを図８に示す。なお、内部減衰モデルとしては因果的履歴減衰モデルを用いた。図８に示す解析対象モデルは、基礎幅(＝20ｍ)に対し、内部領域の幅を(Ａ)は1.5倍、(Ｂ)は３倍、(Ｃ)は５倍と変化させている。１／２逆対称モデルを用いているため、図８における内部領域の長さは上記の１／２となっている。このうち図８(Ａ)に示すモデルは、波動境界モデルとしてエネルギー伝達境界を設けた解析対象モデル（図６参照）と同一形状である。図８に示す各モデルの中心軸位置の最大応答水平加速度を、エネルギー伝達境界を用いた周波数域の解析結果と比較した結果を次の表５に示す。
【００８３】
【表５】
【００８４】
表５より明らかなように、図８(Ａ)に示す解析対象モデルを用いた時刻歴応答解析の結果は、周波数域の応答解析の結果に対し、建物頂部（節点番号１）で16％、基礎下端（節点番号4）で40％以上の差異を生じており、精度が非常に低い。これに対して、図８(Ｂ)の解析対象モデル、図８(Ｃ)の解析対象モデルとモデル化領域が大きくなるに従って差異が小さくなる傾向を示している。しかしながら、基礎幅の５倍のモデル化領域を有する図８(Ｃ)の解析対象モデルでも、基礎下端以下で周波数域の応答解析の結果に対して20％以上の差異が生じている。上記の結果を表４と比較しても明らかなように、波動境界モデルとしてエネルギー伝達境界を用いた解析対象モデルは、波動境界モデルとして粘性境界モデルを用いた解析対象モデルよりも明らかに高精度での時刻歴応答解析を実現することができる。
【００８５】
〔第２実施例〕
エネルギー伝達境界は２次元面内問題以外に２次元面外問題、軸対称問題についても定式化されている。このため、本願発明者は２次元面外問題、軸対称問題についても、エネルギー伝達境界の伝達境界マトリクスが、本発明によって精度良く時間領域のインパルス応答に変換されるか否かを確認する解析検討を行った。この解析検討に用いた地盤モデルを図９に、この地盤モデルにおける定数を表６に、時間領域への変換の条件を表７に各々示す。
【００８６】
【表６】
【００８７】
【表７】
【００８８】
また、時間領域への変換（インパルス応答の演算）にはＣ法を用い、算出されたインパルス応答の各成分のうち同時成分ｋ0,ｃ0,ｍ0及び時間遅れ成分ｋ1〜ｋ10,ｃ1〜ｃ10を用いて伝達境界マトリクスを再現した。
【００８９】
(1)２次元面外問題
２次元面外問題では１節点当りの解析自由度が面外方向の１自由度であり、伝達境界マトリクスは３行×３列となる。算出されたインパルス応答から再現された伝達境界マトリクスを、元のデータ点と共に図１０に示す。図１０からも明らかなように、再現された伝達境界マトリクスはほぼ全てのデータ点と良好に対応しており、本発明が２次元面外問題にも有効であることが確認された。
【００９０】
(2)軸対称問題
軸対称問題では１節点当りの解析自由度が３自由度であり、伝達境界マトリクスは９行×９列となる。算出されたインパルス応答から再現された９行×９列の伝達境界マトリクスのうちの１〜６行目×１〜６列目を、元のデータ点と共に図１１に示す。図１１においても、再現された伝達境界マトリクスはほぼ全てのデータ点と良好に対応しており、本発明が軸対称問題にも有効であることが確認された。
【図面の簡単な説明】
【００９１】
【図１】地盤モデル(１／２モデル)の一例を示す概念図である。
【図２】本願発明者による解析検討において、２次元面内問題で得られたインパルス応答から再現した伝達境界マトリクスを元のデータ点と共に示す線図である。
【図３】本実施形態に係るコンピュータの概略構成を示すブロック図である。
【図４】地震応答解析処理の内容を示すフローチャートである。
【図５】図４の地震応答解析処理における各ステップの関係（データの流れ）を示す概念図である。
【図６】解析対象モデルの一例を示す概念図である。
【図７】本願発明者による解析検討において、(Ａ)はＢ法を、(Ｂ)はＣ法を適用して各々演算したインパルス応答から再現した伝達境界マトリクスを元のデータ点と共に示す線図である。
【図８】本願発明者による解析検討に用いた、粘性境界モデルを設けた解析対象モデルを示す概念図である。
【図９】本願発明者による解析検討に用いた地盤モデルを示す概念図である。
【図１０】本願発明者による解析検討において、２次元面外問題で得られたインパルス応答から再現した伝達境界マトリクスを元のデータ点と共に示す線図である。
【図１１】本願発明者による解析検討において、軸対称問題で得られたインパルス応答から再現した伝達境界マトリクスを元のデータ点と共に示す線図である。
【符号の説明】
【００９２】
１０ コンピュータ
１０Ｃ 記憶部
１２ ディスプレイ
１４ キーボード
１６ マウス
【特許請求の範囲】
【請求項１】
地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算するマトリクス演算手段と、
物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算するインパルス応答演算手段と、
前記インパルス応答演算手段によって演算された前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行う解析手段と、
を含む応答解析装置。
【請求項２】
前記エネルギー伝達境界の前記伝達境界マトリクスは、内部領域の質量マトリクスを[ＭＩ]、剛性マトリクスを[ＫＩ]、変位ベクトルを{ｕ(ω)}、境界力ベクトルを{Ｆ(ω)}、補正力ベクトルを[Ｑ]{ｕ＊(ω)}とし、
(−ω２[ＭＩ]＋[ＫＩ]＋[Ｒ]){ｕ(ω)}＝−ｙ"(ω)[ＭＩ]{1}＋{Ｆ(ω)} …(１)
{Ｆ(ω)}＝([Ｒ]−[Ｑ]){ｕ＊(ω)} …(２)
前記解析対象モデル全体の運動方程式を上記の(１),(２)式で表したときに、上記の(１),(２)式の[Ｒ]で表されることを特徴とする請求項１記載の応答解析装置。
【請求項３】
前記解析対象モデルに前記波動境界モデルとして設けられているエネルギー伝達境界は、２次元面内問題、２次元面外問題、及び、軸対称問題の何れかに対応するモデルであることを特徴とする請求項１記載の応答解析装置。
【請求項４】
前記インパルス応答演算手段は、物体の変位に依存するインパルス応答の同時成分をｋ０、物体の速度に依存するインパルス応答の同時成分をｃ０、物体の変位に依存するインパルス応答のΔt刻みの時間遅れ成分をｋj、物体の速度に依存するインパルス応答のΔt刻みの時間遅れ成分をｃj（但しｊは自然数でtj＝Δt・j）、時間領域での物体の変位をu(t)、速度をu'(t)、加速度をu"(t)としたときに、前記インパルス応答ＦＡ(t)を規定する数式として、
【数１】
上記(３)式を用い、前記振動がＮ種（Ｎ＝ｎ＋１）の周波数のときの前記伝達境界マトリクスの値に基づいて前記インパルス応答を求めることを特徴とする請求項１記載の応答解析装置。
【請求項５】
前記インパルス応答演算手段は、前記インパルス応答を規定する数式として、前記剛性項及び前記減衰項に加えて、前記物体の加速度に依存し少なくとも同時成分を含んで成る質量項も含む数式を用いて、前記エネルギー伝達境界の前記インパルス応答を演算することを特徴とする請求項１記載の応答解析装置。
【請求項６】
前記インパルス応答演算手段は、物体の変位に依存するインパルス応答の同時成分をｋ０、物体の速度に依存するインパルス応答の同時成分をｃ０、物体の加速度に依存するインパルス応答の同時成分をｍ０、物体の変位に依存するインパルス応答のΔt刻みの時間遅れ成分をｋj、物体の速度に依存するインパルス応答のΔt刻みの時間遅れ成分をｃj（但しｊは自然数でtj＝Δt・j）、時間領域での物体の変位をu(t)、速度をu'(t)、加速度をu"(t)としたときに、前記質量項も含む数式として、前記インパルス応答ＦB(t)を規定する数式である、
【数２】
上記(４)式を用い、前記振動がＮ種（Ｎ＝ｎ＋１）の周波数のときの前記伝達境界マトリクスの値に基づいて前記インパルス応答を演算することを特徴とする請求項５記載の応答解析装置。
【請求項７】
前記インパルス応答演算手段は、前記(４)式における物体の加速度に依存するインパルス応答の同時成分ｍ０に対する修正値Δｍ0、物体の変位に依存するインパルス応答の同時成分ｋ０に対する修正値Δｋ0を下記の(５)式によって演算し、
【数３】
（但し、上記(５)式において、
【数４】
であり、上記(６)式におけるＲe(ＳＢ(ωｉ))は、前記(４)式を用いて演算した前記エネルギー伝達境界の前記インパルス応答から再現された、下記の(７)式で表される前記伝達境界マトリクスＳＢ(ω)の振動数ωｉでの実部の値を表し、
【数５】
上記(６)式におけるＲe(Ｄ(ωｉ))は、(４)式に基づく前記エネルギー伝達境界の前記インパルス応答の演算に用いた、振動数ωｉでの前記伝達境界マトリクスのデータＤ(ωｉ)のうちの実部の値を表す）、演算した修正値Δｍ0,Δｋ0を用いて前記同時成分ｍ０,ｋ０を修正することを特徴とする請求項６記載の応答解析装置。
【請求項８】
前記インパルス応答演算手段は、前記(４)式における物体の速度に依存するインパルス応答の同時成分をｃ０に対する修正値Δｃ0も下記の(８)式によって演算し、
Δｃ0＝−Ｅ／Ｂ …(８)
（但し、上記(８)式において、
【数６】
であり、上記(９)式におけるＩm(ＳＢ(ωｉ))は、前記(４)式を用いて演算した前記エネルギー伝達境界の前記インパルス応答から再現された、前記(７)式で表される前記伝達境界マトリクスＳＢ(ω)の振動数ωｉでの虚部の値を表し、上記(９)式におけるＩm(ＳＢ(ωｉ))は、(４)式に基づく前記エネルギー伝達境界の前記インパルス応答の演算に用いた、振動数ωｉでの前記伝達境界マトリクスのデータＤ(ωｉ)のうちの虚部の値を表す）、演算した修正値Δｃ0を用いて前記同時成分ｃ０も修正することを特徴とする請求項７記載の応答解析装置。
【請求項９】
コンピュータに、
地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算させ、
物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算させ、
前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行わせる応答解析方法。
【請求項１０】
コンピュータを、
地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算するマトリクス演算手段、
物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算するインパルス応答演算手段、
及び、前記インパルス応答演算手段によって演算された前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行う解析手段
として機能させる応答解析プログラム。
【請求項１】
地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算するマトリクス演算手段と、
物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算するインパルス応答演算手段と、
前記インパルス応答演算手段によって演算された前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行う解析手段と、
を含む応答解析装置。
【請求項２】
前記エネルギー伝達境界の前記伝達境界マトリクスは、内部領域の質量マトリクスを[ＭＩ]、剛性マトリクスを[ＫＩ]、変位ベクトルを{ｕ(ω)}、境界力ベクトルを{Ｆ(ω)}、補正力ベクトルを[Ｑ]{ｕ＊(ω)}とし、
(−ω２[ＭＩ]＋[ＫＩ]＋[Ｒ]){ｕ(ω)}＝−ｙ"(ω)[ＭＩ]{1}＋{Ｆ(ω)} …(１)
{Ｆ(ω)}＝([Ｒ]−[Ｑ]){ｕ＊(ω)} …(２)
前記解析対象モデル全体の運動方程式を上記の(１),(２)式で表したときに、上記の(１),(２)式の[Ｒ]で表されることを特徴とする請求項１記載の応答解析装置。
【請求項３】
前記解析対象モデルに前記波動境界モデルとして設けられているエネルギー伝達境界は、２次元面内問題、２次元面外問題、及び、軸対称問題の何れかに対応するモデルであることを特徴とする請求項１記載の応答解析装置。
【請求項４】
前記インパルス応答演算手段は、物体の変位に依存するインパルス応答の同時成分をｋ０、物体の速度に依存するインパルス応答の同時成分をｃ０、物体の変位に依存するインパルス応答のΔt刻みの時間遅れ成分をｋj、物体の速度に依存するインパルス応答のΔt刻みの時間遅れ成分をｃj（但しｊは自然数でtj＝Δt・j）、時間領域での物体の変位をu(t)、速度をu'(t)、加速度をu"(t)としたときに、前記インパルス応答ＦＡ(t)を規定する数式として、
【数１】
上記(３)式を用い、前記振動がＮ種（Ｎ＝ｎ＋１）の周波数のときの前記伝達境界マトリクスの値に基づいて前記インパルス応答を求めることを特徴とする請求項１記載の応答解析装置。
【請求項５】
前記インパルス応答演算手段は、前記インパルス応答を規定する数式として、前記剛性項及び前記減衰項に加えて、前記物体の加速度に依存し少なくとも同時成分を含んで成る質量項も含む数式を用いて、前記エネルギー伝達境界の前記インパルス応答を演算することを特徴とする請求項１記載の応答解析装置。
【請求項６】
前記インパルス応答演算手段は、物体の変位に依存するインパルス応答の同時成分をｋ０、物体の速度に依存するインパルス応答の同時成分をｃ０、物体の加速度に依存するインパルス応答の同時成分をｍ０、物体の変位に依存するインパルス応答のΔt刻みの時間遅れ成分をｋj、物体の速度に依存するインパルス応答のΔt刻みの時間遅れ成分をｃj（但しｊは自然数でtj＝Δt・j）、時間領域での物体の変位をu(t)、速度をu'(t)、加速度をu"(t)としたときに、前記質量項も含む数式として、前記インパルス応答ＦB(t)を規定する数式である、
【数２】
上記(４)式を用い、前記振動がＮ種（Ｎ＝ｎ＋１）の周波数のときの前記伝達境界マトリクスの値に基づいて前記インパルス応答を演算することを特徴とする請求項５記載の応答解析装置。
【請求項７】
前記インパルス応答演算手段は、前記(４)式における物体の加速度に依存するインパルス応答の同時成分ｍ０に対する修正値Δｍ0、物体の変位に依存するインパルス応答の同時成分ｋ０に対する修正値Δｋ0を下記の(５)式によって演算し、
【数３】
（但し、上記(５)式において、
【数４】
であり、上記(６)式におけるＲe(ＳＢ(ωｉ))は、前記(４)式を用いて演算した前記エネルギー伝達境界の前記インパルス応答から再現された、下記の(７)式で表される前記伝達境界マトリクスＳＢ(ω)の振動数ωｉでの実部の値を表し、
【数５】
上記(６)式におけるＲe(Ｄ(ωｉ))は、(４)式に基づく前記エネルギー伝達境界の前記インパルス応答の演算に用いた、振動数ωｉでの前記伝達境界マトリクスのデータＤ(ωｉ)のうちの実部の値を表す）、演算した修正値Δｍ0,Δｋ0を用いて前記同時成分ｍ０,ｋ０を修正することを特徴とする請求項６記載の応答解析装置。
【請求項８】
前記インパルス応答演算手段は、前記(４)式における物体の速度に依存するインパルス応答の同時成分をｃ０に対する修正値Δｃ0も下記の(８)式によって演算し、
Δｃ0＝−Ｅ／Ｂ …(８)
（但し、上記(８)式において、
【数６】
であり、上記(９)式におけるＩm(ＳＢ(ωｉ))は、前記(４)式を用いて演算した前記エネルギー伝達境界の前記インパルス応答から再現された、前記(７)式で表される前記伝達境界マトリクスＳＢ(ω)の振動数ωｉでの虚部の値を表し、上記(９)式におけるＩm(ＳＢ(ωｉ))は、(４)式に基づく前記エネルギー伝達境界の前記インパルス応答の演算に用いた、振動数ωｉでの前記伝達境界マトリクスのデータＤ(ωｉ)のうちの虚部の値を表す）、演算した修正値Δｃ0を用いて前記同時成分ｃ０も修正することを特徴とする請求項７記載の応答解析装置。
【請求項９】
コンピュータに、
地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算させ、
物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算させ、
前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行わせる応答解析方法。
【請求項１０】
コンピュータを、
地盤単独又は地盤と建物から成る解析対象物体を、モデルの内部領域と外部領域の間に波動境界モデルとしてのエネルギー伝達境界が設けられた解析対象モデルへモデル化したときの、前記エネルギー伝達境界の主要な構成要素であり、周波数領域の複素数として定義され、かつ強い振動数依存性を有する伝達境界マトリクスを演算するマトリクス演算手段、
物体を振動させる外力と前記物体の挙動との関係を時間領域で表すインパルス応答を規定する数式として、少なくとも、前記物体の変位に依存し同時成分と時間遅れ成分から成る剛性項と、前記物体の速度に依存し同時成分と時間遅れ成分から成る減衰項を含む数式を用い、前記振動が各周波数のときの前記伝達境界マトリクスの値に基づいて前記エネルギー伝達境界の前記インパルス応答を演算するインパルス応答演算手段、
及び、前記インパルス応答演算手段によって演算された前記エネルギー伝達境界のインパルス応答を用いて、前記解析対象モデルに対して時間領域での応答解析を行う解析手段
として機能させる応答解析プログラム。
【図３】
【図４】
【図５】
【図１】
【図２】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図４】
【図５】
【図１】
【図２】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【公開番号】特開２００９−２５０８０５（Ｐ２００９−２５０８０５Ａ）
【公開日】平成２１年１０月２９日（２００９．１０．２９）
【国際特許分類】
【出願番号】特願２００８−９９６７５（Ｐ２００８−９９６７５）
【出願日】平成２０年４月７日（２００８．４．７）
【出願人】（０００００３６２１）株式会社竹中工務店 (1,669)
【Ｆターム（参考）】
【公開日】平成２１年１０月２９日（２００９．１０．２９）
【国際特許分類】
【出願日】平成２０年４月７日（２００８．４．７）
【出願人】（０００００３６２１）株式会社竹中工務店 (1,669)
【Ｆターム（参考）】
[ Back to top ]