第 3 章 時間依存解析モデルの定式化と検証
3.2 若材齢挙動のモデル化
3.2.4 応力・クリープ解析フローと有限要素定式化
3.2.4.3 応力・クリープ解析の有限要素定式化
第3章 時間依存解析モデルの定式化と検証
- 44 -
図－3.11に示すように，この定式化では𝑟 = 0 時点で熱生成が生じており実現象と乖 離する。これは，断熱温度曲線から平均熱生成速度を求めたことにより，断熱温度曲線 の時刻0に対応する初期温度が反映されたためと考えられる。この手法では，r = 0に対
応するr’を初期反応度と定義し，r軸をr’軸に写像することで実現象との整合を図って
いる。Reinhardt & Blaauwendraadらが検討したコンクリート調合の範囲では，初期反応
度としてr’ = 0.1が推奨値として示されている。
第3章 時間依存解析モデルの定式化と検証
- 45 -
ここで，[𝐵]はひずみ－変位マトリクス，{𝑅̇}𝑎は外部節点力速度である。次に，材料の 構成挙動は，次式により表現される。
{𝜀̇} = {𝜀̇}
𝑒+ {𝜀̇}
𝑐+ {𝜀̇}
𝑠+ {𝜀̇}
𝑇(3.2.31) ここで，{𝜀̇}は全体ひずみ速度，{𝜀̇}𝑒は弾性ひずみ速度，{𝜀̇}𝑐はクリープひずみ速度成分，
{𝜀̇}𝑠は収縮によるひずみ速度成分，{𝜀̇}𝑇は温度膨張・収縮によるひずみ速度成分である。
すべての変数は，時間𝑡の関数である。弾性ひずみ速度と応力速度は次式により関係づ けられる。
{𝜎̇} = [𝐷]({𝜀̇} − {𝜀̇}
𝑐− {𝜀̇}
𝑠− {𝜀̇}
𝑇)
( 3 .2 . 32) 式(3.2.32)を式(3.2.30)に代入すると次式が導かれる。[𝐾]{𝛿̇} = {𝑅̇}
𝑎+ {𝑅̇}
𝜂(3.2.33) ここで，[𝐾]は要素剛性マトリクス，{𝛿̇}は未知節点変位速度，{𝑅̇}𝑎は外部節点力の速度 に対応する等価節点力ベクトル，{𝑅̇}𝜂は収縮，クリープおよび温度過程のひずみ速度に 対応する等価節点力ベクトルであり，次式により表現される。
[𝐾] = ∑ ∫ [𝐵]
𝑒 𝑣𝑒 𝑇[𝐷] [𝐵] 𝑑𝑉
(3.2.34){𝑅̇}
𝜂= ∑ ∫ [𝐵]
𝑒 𝑣𝑒 𝑇[𝐷] {𝜂̇} 𝑑𝑉
(3.2.35)ここで，{𝜂̇}は収縮ひずみ，クリープひずみおよび温度過程で生じるひずみによる速度 ベクトルである。
続いて，時間－離散化手法について説明するため，弾性挙動およびクリープ挙動は 種々の体積ひずみ変化を重ね合わせることにより求められると仮定する。ここでは，de
Borst ら［3.23］を参照して，時間𝑡における全ひずみベクトル{𝜀(𝑡)}を次式のように整
理する。
{𝜀(𝑡)} = ∫ 𝐽(𝑡, 𝜏)
0𝑡𝑪𝝈̇(𝜏)𝑑𝜏 + {𝜀}
𝑠,𝑇( 3 . 2 . 3 6)
第3章 時間依存解析モデルの定式化と検証
- 46 -
ここで，𝝈̇は応力速度テンソル，𝐽(𝑡, 𝜏)はクリープ関数，𝜏 は応力作用時のコンクリート
材齢，{𝜀}𝑠,𝑇は応力履歴の影響を受けない温度ひずみと収縮ひずみの和，𝑪はコンプライ アンスであり，各成分がポアソン比の関数として，次式により定義される。なお，太字 はマトリクスを意味している。
𝐶
𝑖𝑗𝑘𝑙=
12(1 + 𝜈)(𝛿
𝑖𝑘𝛿
𝑗𝑙+ 𝛿
𝑖𝑙𝛿
𝑗𝑘) − 𝜈𝛿
𝑖𝑗𝛿
𝑘𝑙( 3 . 2 . 3 7 )
ここで，𝜈は荷重作用時間には依存しない一定値と仮定したポアソン比，𝛿𝑖𝑗と𝛿𝑘𝑙はクロ
ネッカーデルタ(𝑖, 𝑗, 𝑘, 𝑙 = 𝑥, 𝑦, 𝑧)である。時間間隔∆tを考え，対応する応力増分とひずみ 増分を次式のように定義する。
{∆𝜎} = {𝜎(𝑡 + ∆𝑡)} − {𝜎(𝑡)}
(3.2.38){∆𝜀} = {𝜀(𝑡 + ∆𝑡)} − {𝜀(𝑡)}
(3.2.39)式(3.2.39)を考慮すると式(3.2.36)は次式で表現できる。
{∆𝜀} = [ ∫ 𝐽(𝑡 + ∆𝑡, 𝜏)
𝑡+∆𝑡
0
𝑪𝝈̇(𝜏)𝑑𝜏 + {𝜀(𝑡 + ∆𝑡)}
𝑠,𝑇] − [∫ 𝐽(𝑡, 𝜏)
𝑡
0
𝑪𝝈̇(𝜏)𝑑𝜏 + {𝜀(𝑡)}
𝑠,𝑇]
(3.2. 40)
式(3.2.40)を整理すると次式となる。
{∆𝜀} = ∫(𝐽(𝑡 + ∆𝑡, 𝜏) − 𝐽(𝑡, 𝜏))
𝑡
0
𝑪𝝈̇(𝜏)𝑑𝜏 + ∫ 𝐽(𝑡 + ∆𝑡, 𝜏)
𝑡+∆𝑡
𝑡
𝑪𝝈̇(𝜏)𝑑𝜏 + {𝜀(∆𝑡)}
𝑠,𝑇(3.2. 41) 有限要素法に組み込むために，式(3.2. 41)の応力増分がひずみ或いは応力履歴の関数 となるように再整理する。なお，応力が時間増分∆𝑡にわたり線形に変化すると仮定し
（𝝈̇ ≈ ∆𝝈/∆𝒕と仮定），一般化中点則 𝑡𝑖 ≤ 𝑡∗≤ 𝑡𝑖+1 を適用する。ここで，式(3.2.41)の右 辺第2項に着目する。この項は，現ステップでのクリープひずみ増分を表しており，こ れを∆𝜀とおいて整理すると次式となる。
第3章 時間依存解析モデルの定式化と検証
- 47 -
{∆𝜀} = ∫
𝑡𝑡+∆𝑡𝐽(𝑡 + ∆𝑡, 𝜏) 𝑪
∆𝝈∆𝒕(𝜏)𝑑𝜏
(3.2.42)式(3.2.42)は次式のように表現できる。
{∆𝜀}
𝑪
=
∆𝒕𝟏∫
𝑡𝑡+∆𝑡𝐽(𝑡 + ∆𝑡, 𝜏) (𝜏)𝑑𝜏 ∙ ∆𝝈
(3.2.43) ここで，𝐸̃
−1(𝑡
∗) =
∆𝑡1∫
𝑡𝑡+∆𝑡𝐽(𝑡 + ∆𝑡, 𝜏) 𝑑𝜏
(3.2.44)とおくと次式が得られる。
∆𝝈 = 𝑪
−𝟏∙ 𝐸̃(𝑡
∗) ∙ ∆𝜀
(3.2.45) ここで，𝑪−𝟏= 𝑫と置いて，式(3.2.41)を式(3.2.42)～式(3.2.45)に基づき∆𝜎について整理 すると最終的に次式を得る。∆𝜎 = 𝐸̃(𝑡
∗)𝑫{∆𝜀} − 𝐸̃(𝑡
∗) ∫(𝐽(𝑡 + ∆𝑡, 𝜏) − 𝐽(𝑡, 𝜏))
𝑡
0
𝝈̇(𝜏)𝑑𝜏 + 𝐸̃(𝑡
∗)𝑫{∆𝜀}
𝑠,𝑇(3.2. 46) 式(3.2.46)における𝐽(𝑡, 𝜏)は，一般的には，長期クリープシミュレーションに対して
Dirichlet級数が採用され，短期クリープシミュレーションに対してTaylor級数が採用さ
れることで実現象と良い対応を示すことが知られている［3.23］。
本論文におけるクリープひずみのモデル化は，一般的に利用されることの多い有効弾 性係数法とStep-by-step法の計算精度とコストに関する問題点を考慮して，de Borstら の二重べき乗則をTaylor級数近似する手法を採用した［3.23］。ここでは，クリープ関 数を次式のように定義する。
𝐽(𝑡, 𝜏) =
𝐸(𝜏)1[1 + 𝑞𝜏
−𝑑(𝑡 − 𝜏)
𝑝]
(3.2.47)第3章 時間依存解析モデルの定式化と検証
- 48 -
ここで，𝑞, 𝑑, 𝑝は材料定数，𝜏−𝑑は荷重の作用時間に関する係数，(𝑡 − 𝜏)𝑝はクリープ期 間に関する関数である。𝑔(𝑡 − 𝜏) = (𝑡 − 𝜏)𝑝とおいて，𝑔(𝑡 − 𝜏)を 𝑡 − 𝜏 = 𝑡𝑑（展開点）
まわりにTaylor展開して整理すると次式を得る。
𝑔(𝑡 − 𝜏) = ∑
𝑁𝑟=0ℎ
𝑟(𝑡, 𝑡
𝑑)𝜏
𝑟 (3.2.48)ここで，NはTaylor展開で使われる級数の数（本論文ではN = 5とする。），ℎ𝑟はべき数 𝑝に依存する(𝑡 − 𝑡𝑑)の関数である。式(3.2.46)右辺の第二項を{𝜎̃(𝑡)}と置き，式(3.2.48) を代入して整理すると次式が得られる。
{𝜎̃(𝑡)} = −𝐸̃(𝑡
∗)𝑞 ∑
5𝑟=0(ℎ
𝑟(𝑡 + ∆𝑡, 𝑡
𝑑) − ℎ
𝑟(𝑡, 𝑡
𝑑))𝜀̃
𝑟(𝑡)
(3.2.49)ただし，
𝜀̃
𝑟(𝑡) = ∫
0𝑡𝐸(𝜏)1𝜏
𝑟−𝑑𝜎̇(𝜏)𝑑𝜏
(3.2.50)ここで，𝜀̃𝑟項は解析時の総和により次式のように近似される。
𝜀̃
𝑟(𝑡 + ∆𝑡) = 𝜀̃
𝑟(𝑡) + ∫ 1
𝐸(𝜏) 𝜏
𝑟−𝑑𝜎̇(𝜏)𝑑𝜏
𝑡+∆𝑡
𝑡
≈ 𝜀̃
𝑟(𝑡) + (
𝐸(𝑡∆𝜎∗)∆𝑡) (
(𝑡+∆𝑡)𝑟−𝑑+1𝑟−𝑑+1−𝑡𝑟−𝑑+1)
(3.2.51)式(3.2.49)は，積分区間に実時間が現れないので，各材料点の全荷重履歴の保存を要 せず，N = 5までの近似解を保存するのみである。また，時間ステップが増しても計算 時間が長くならない特徴も有しており，計算量の節減を図ることが可能となる。