有限要素法による板材および管材の塑性局所化解析に関する研究

(1)

有限要素法による板材および管材の塑性局所化解析に関する研究

2004

年度

只野裕一

(2)

大ひずみを前提とした塑性変形解析は，材料加工，材料の破断や構造物の破壊の予測など，多岐にわたる応用領域が考えられる．材料の挙動を表現する構成式は数値解析の結果に大きく影響するため，高精度な解析結果を得ようとすれば適切な材料構成式を選択しなければならない．塑性という現象を数学的モデルによって表現するという試み自体は19世紀から行われており，20世紀半ばにHillによって体系化された^[1-1]．これらの研究による成果は古典塑性論として知られ，古典的ながらも実用的に十分な近似である場合も多く，特に微小変形の範囲であれば非常に良好な解を与える．このため，これらの理論は現在でも幅広い分野で用いられている．

しかしながら，塑性を考える上で避けることのできない大ひずみを生じるような問題に対しては，古典塑性論では説明のつかない現象の方が遙かに多い．大ひずみを考慮した塑性構成式の研究は20世紀後半にようやく本格化した^[1-2]．この結果，1960年代から 80 年代半ばにかけて基本的な定式化こそ確立されたものの，より現実的なモデルという観点からは未解決の問題も多く，今なお精力的な研究が行われている．

塑性を考慮した問題の中で，工業的に重要なものの一つが塑性不安定現象である．

金属における塑性不安定現象の概要を図1.1に示す．金属材料に引っ張り荷重を加えていくと，降伏応力を超えて塑性変形を生じた後もある程度までは荷重が増加するが，

ある点で荷重が極大点を迎え荷重が減少へと転じることが知られている．換言すればこの極大荷重こそが材料の耐えうる最大荷重であり，それを超えてしまうともはや荷重を減少させても復元力が働かなくなることを意味している．変形の進展が抑制できなくなってしまうことから，これを塑性不安定現象と呼んでいる^[1-3]．塑性不安定発現後の挙動は材料によって大きく異なるが，鋼材など金属材料の一軸引っ張りの場合には部材全体に緩やかなくびれが生じることが知られており，これを拡散くびれと呼ぶ．この段階では塑性変形は部材全体で進行しており，弾性除荷は発生していない．

さらに変形が進行すると部材内部における変形の一様性が崩れ，ある一部分の応力が低下して弾性除荷を生じる．この除荷域はその後急速に拡大していき，結果として部材の非常に狭い領域に塑性変形が集中することになる．このような現象は塑性局所化現象と呼ばれ^[1-4]，この時生じる局所的なくびれを局所くびれという．拡散くびれ発生の段階ではこのような塑性変形の局所化が生じていないことから，この二つのくびれ様式は異なる変形の段階であるとみなされる．局所化によって生じる塑性変形域は，

X型や引っ張り方向に対して斜めもしくは引っ張り方向と直交する帯状となる．せん断変形によって生じるX型や斜めの変形様式が，いわゆるせん断帯である．このような現象は，材料に予め含まれる微小な形状不整や変形過程において外部から混入する擾乱が原因であるとされ^[1-5]，数学的には一種の分岐問題であると考えられている^[1-3]．

(6)

室温における金属材料では，変形の局所化が発生すると直ちに局所くびれが非常に狭い領域に集中し，かつ局所くびれによる変形モードがそのまま破断のモードとなることが知られており，局所くびれの発生が破断にほぼ直結する．すなわち塑性局所化は破断の直接のトリガーとなる現象であり，その予測は極めて重要な課題であるといえる．

1.2 過去の研究

1.2.1 板材の局所化解析

板材の加工を考えた場合，曲げ変形に伴う面外せん断変形は無視できることも多く，

このような場合には問題は二次元（平面応力）として考えることが出来る．この時，

材料の破断は面内応力のみに支配されると考えられるため，板材の破断は面内二軸引っ張り状態における局所化現象の結果とみなすことが出来る．一様変形を保ったまま材料を加工できる限界ひずみのことを成形限界ひずみと呼ぶ．様々なひずみ経路で成形限界ひずみを評価できれば，ひずみ空間上にこれをプロットした線図が描くことが出来，これを成形限界線図（Forming limit diagram; FLD）と呼ぶ．成形限界線図は，

プレス加工，絞り加工などにおいて成形可能な変形範囲の一つの目安となる．このため，成形限界の評価に関する研究は，実験，理論の両面から広く行われている．

古くは Hill によって数学的な検討がなされており^[1-6]，その後も後述する尖り点形成を考慮したモデルによる解析^[1-7]，板材表面の不均一性を考慮した解析^[1-8]，移動硬化の影響に関する検討^[1-9]，降伏関数に関する検討^[1-10]などが試みられている．

このような成形限界の予測を目的とした板材の二軸荷重に対する解析モデルの一つとして知られているのが，Marciniak and Kuczynskiによって提案された，いわゆる M-Kモデルである[1-11][1-12]．材料の異方性を考慮した場合のM-Kモデルの概念を図1.2

Yielding point

Maximum loading point

Elongation

Load

Diffused necking

Localized necking

Failure

Fig. 1.1 Plastic instability and localization.

(7)

に示す^[1-13]． M-K モデルでは，均一な平板にバンド状の非均一な部分があると仮定する．図1.2では，初期形状において異方性主軸xˆ₁が引っ張り方向x₁に対してθÎ^だけ傾いているとしており，初期不整帯の方向 sÎに直交するベクトルmÎはx₁に対してψÎ 傾いている．初期不整帯の非均一性はどのようなものでも良く，具体例としては板厚，

材料定数などをわずかに変えることが考えられる．このようなモデルに対し，二軸ひずみ（速度）を負荷しバンド内外の変形を追跡する．ある程度のひずみまでは変形の一様性が保たれるが，やがてこの一様性が崩れる．この時のひずみが一様変形の限界であるとみなし，板材の成形限界と考える．これを様々なひずみ経路で行うことで，

ひずみ空間上に成形限界ひずみのプロットを描くことが出来る．M-Kモデルは非常に明快なモデル化であり，また汎用性も高いため各種の構成則への適用が可能である．

このためいくつかの材料モデルに対し，M-Kモデルを用いて成形限界線図を求めた研究が報告されている．Zhao et al.は予ひずみを加えた条件下でのFLD評価を試み，実験との比較も行っている^[1-14]．Xu and Weinmannは現象論的材料モデルに，異方性降伏関数を導入した際のFLDを^[1-15]，Wu et al.は結晶塑性モデルによるFLDを[1-16][1-17]

それぞれ求めている．Kuroda and Tvergaardは，現象論的モデルに複数の異方性降伏関数を適用した際のFLDを算出し，降伏関数のFLDへの影響を詳細に調査している^[1-13]． 1.2.2 現象論的塑性構成式

塑性構成式として古くから広く用いられているものの一つが，J2流れ則である^[1-18]．これは塑性変形の方向が降伏曲面の法線方向と一致すると考える関連流れ則の一種であり，降伏関数にMisesの降伏関数を用いたものである．元は微小ひずみ弾塑性のためのモデルとして考案されたものであるが，大ひずみ弾塑性や速度依存性を考慮した粘塑性モデルへの拡張も容易であり^[1-2]，現在でも様々な問題に適用されている．

Fig. 1.2 Illustration of M-K model.

ψ^I m^I

s^I x2

ˆ1

x ˆ2

x

fixed reference axis orthotropic axis

θ ^I

x1

(8)

数学的に取り扱いやすく材料定数の同定も簡単であるため実用性が高く，汎用有限要素解析コードにも広く利用されている．しかしながら，J2流れ則を用いて平面ひずみ状態における単軸引っ張り解析を行った場合，塑性不安定発生後の拡散くびれまでは再現できるものの，その後のせん断帯形成が表現できないなど，実現象を必ずしも再現できない場合があることが指摘されている[1-4][1-19][1-20]．拡散くびれに続くせん断帯形成においては，ひずみの経路が比例負荷方向から急激に変化するが，このとき J2 流れ則による塑性変形方向の予測が必ずしも現実とは整合していないことが，せん断帯が予測できない原因であるとされている．このため，ひずみ経路が比例負荷から大きく外れる際の塑性変形方向をより適切に表現できる構成式が模索された．

Budiansky は応力と塑性ひずみを一対一に対応づける変形理論について再検討した上で，これを時間微分した速度型の変形理論を示し^[1-21]，Stören and Riceはこの速度型変形理論を大変形に拡張したJ2変形理論を提案している^[1-7]．従来のJ2流れ則において塑性変形の方向は常に降伏曲面の法線方向と一致したのに対し，J2変形理論では降伏曲面の接線方向の成分もこれに加えている．J2変形理論は比例負荷から大きく外れるような負荷経路に対する解析には適用することが難しいため，これをさらに一般化する試みがなされた．Christoffersen and Hutchinson^[1-22]，後藤[1-23]-[1-33]はそれぞれJ2変形理論を基に比例負荷から大きく外れる負荷経路にも適用できる構成式を提案した．

これらの構成式は，降伏曲面の負荷点において尖り点が形成されることを再現可能なため，コーナー理論とも呼ばれる．このような構成式を用いることで，前述の平面ひずみ状態においてもせん断帯が形成されることが知られている[1-19][1-34]．J2 流れ則は塑性変形速度の方向が（偏差）応力の方向と一致することから，共軸形の構成式と呼ばれるが，J2変形理論やコーナー理論はこれらの方向が一致せず，非共軸形の構成式となっている．

局所化を伴う大変形塑性解析においては，非共軸構成式がより現実的な解を与えるものと考えられている．しかしながら，J2変形理論を用いて予測される局所化発生後の荷重の落ち込みが実験と比較して著しく大きい，Christoffersen and Hutchinson のコーナー理論（J2コーナー理論）では材料定数が多くその同定に困難が伴うなどの問題点も含んでいる．加えて，通常の塑性構成式に対してひずみ速度硬化則を適用することは熱力学的な矛盾を含んでいることや，非共軸構成式においては構成式の両辺の大きさが一致せず数学的に矛盾していることが指摘されている^[1-35]．このため，これらの矛盾を含まない構成式の構築および有限要素解析への適用が試みられており

[1-35][1-36]，ポリマの破断予測への応用も報告さている^[1-37]．この中で，熱力学的論理体

系に整合する粘塑性構成式は，必然的に非共軸性を有することが示されている．

最近になってKuroda and Tvergaardは，新たな現象論的塑性構成式を提案した^[1-38]．従来のコーナー理論の多くは塑性変形の方向を偏差応力の関数として記述しているが，Kuroda and Tvergaardはこれを変形速度の関数とした．このモデルを用いることで，

平面ひずみ状態でのせん断帯形成を伴う解析において，結晶塑性モデルとの良好な一致が得られることが報告されている^[1-39]．J2流れ則およびKuroda and Tvergaardによ

(9)

るモデルを用いて平面ひずみ単軸引っ張り解析を行った例を図1.3に示す．

Kuroda and Tvergaardはこの解析において，b.c.c.金属に対する多結晶塑性モデルによる解析結果を比較解としている．結晶塑性モデルは，金属材料の降伏が結晶のすべり変形であるという考えを基に構築されたモデルである．結晶塑性モデルは，根本的には連続体力学の上に構築されており無限に小さい物質点に対しても成立するため，

寸法効果は含まれておらず，異なるスケールの解析においても相似則が成り立つ．また前述の熱力学的な矛盾も含んでいることから，現象論的なモデリングの一種とみなすことができる．しかしながら，通常のマクロスコピックな現象論的塑性モデルと比較した場合には，一段階小さなスケールの物理モデルに立脚しているため，特に単結晶や f.c.c.，b.c.c.を有する多結晶金属に対して，構成式や降伏関数などの妥当性の検証において参照解として広く用いられている．また加工硬化曲線や異方性の進展予測などにおいて実験との整合性も示されており，本研究で取り扱うような単調引っ張り解析においては良好な解を与えると考えられる．結晶塑性モデルの特徴については付録Aに併せて示している．

Kuroda and Tvergaardモデルは古くから用いられている非法線則であるJ2コーナー理論を内包しており，材料パラメータの選択によって J2 コーナー理論と一致する．

これより，Kuroda-Tvergaard の非法線則は従来用いられてきた現象論的モデルと多結晶塑性モデルの特徴を兼ね備えており，本研究で取り扱うせん断帯を伴う局所化解析においては適切なモデルの一つであると言える．

J2 flow theory

Non-normality theory by Kuroda and Tvergaard

Equivalent plastic strain

2.0

1.0

0.0

Fig. 1.3 Shear band analysis.

(10)

1.2.3 一般化有限要素法

有限要素解析において解の精度の向上は古くからの大きな課題であり，現在に至るまで様々なアプローチが行われている^[1-40]．これらの一つとして，数学的見地から，

PU (Partition of Unity)条件を満たす形状関数および節点あたりの自由度として多項式

（の係数，定数項が変位自由度に相当）を用いる一般化有限要素法(Generalized Finite Element Method; GFEM)の開発が1980年代から行なわれており[1-41][1-42]，近年でも活発な研究報告が行なわれている[1-43][1-44][1-45][1-46]．一般化有限要素法は，コーナー節点のみを有する低次要素でも高次の補間が可能であるため，要素のゆがみに対しても極めてロバストで，信頼性の高い解析手法であると言える．一般化有限要素法の欠点として，要素当たりの自由度数ひいては計算コストが増大する傾向にある点，基底関数の一次独立性が損なわれる点等が挙げられる．

また解析精度向上のために古くから行われている試みの一つが，節点自由度として従来の並進自由度に加え，回転自由度を新たに加えた定式化を行うものである．

Allman^[1-47]，Bergan and Felippa^[1-48]によって節点回転自由度を導入した三角形要素が提案され，これを基にさまざまな改良が試みられている．このような要素を用いることで，はりの曲げ問題などにおいて従来の要素と比べて非常に精度の高い解が得られることが報告されている．ここに挙げた既往の研究より，節点回転自由度の導入が，変位型有限要素法において解の収束性や要素のゆがみに対するロバスト性の向上に有効であることが示されている．これらの要素の一部はある種の一般化有限要素とみなすことができ，従って一般化有限要素法の収束性やロバスト性を継承しているものと考えられる．

1.3 過去の研究の問題点

M-K モデルは成形限界を予測するための一つの有効な手法であると考えられているが，一方で実際の材料には M-K モデルのようなバンド状の不整が含まれていることはありえず，このようなモデル化は適切ではないという研究者も少なくない．M-K モデルは汎用性の高いモデリングであり，その物理的な妥当性を示すことが出来ればモデルの信頼性が高まるとともに，応用範囲の拡大が期待できる．しかしながら，二軸荷重を制御し，様々なひずみ（応力）経路における成形限界を実験的に求めるのは非常に困難であり，実験との比較による妥当性の確認も十分になされていないのが現状である．実験的検証が困難であることから，有限要素法などの汎用的な数値解析手法によってモデルの妥当性を検証することが期待されるが，板材の二軸引っ張り解析においては，二つの荷重を独立して制御することが必要となることから，従来の解析手法では境界条件を適切に再現することが出来ない．このため，M-Kモデルの数値解析による妥当性の検証もまた困難である．塑性加工において重要となる二軸荷重を受ける板材や管材の局所化予測精度向上のために，二軸荷重を取り扱うことの出来る新たな数値解析手法の提案が必要である．

(11)

一方，回転自由度を導入した要素や一般化有限要素について見ると，これらの要素による解析の多くは線形弾性問題への適用に留まっており，大変形塑性解析を含む非線形性問題に適用した例は少ない．Barros et al.^[1-49]は一般化有限要素法を材料非線形問題に，Terada et al.^[1-50]は一種の一般化有限要素法とみなせるFinite Cover Methodを不連続材料や接触問題に適用しているが，幾何学的な非線形性まで考慮した問題への適用，および解の収束性や要素のゆがみに対する検討は，著者の知る限り報告されていない．線形解析における解の収束性の議論は，非線形解析においても同様に成立する保証はないため，回転自由度の導入が幾何学的非線形問題の解析にも有効であるかの検証が必要であると思われる．また，大変形問題では変形に伴う要素形状の変化によって解析精度が損なわれることが懸念されるため，ゆがみに強い要素を用いた幾何学的非線形解析の定式化は有意義であると考えられる．

加えて，局所化問題を汎用的数値解析手法の一つである有限要素法で取り扱った場合，得られる解は材料モデリングのみならず使用する有限要素にも強く依存することが知られている．特に塑性変形を含む大変形解析においては，要素のゆがみの他に，

非圧縮変形の拘束から固い解が得られる locking 現象が生じることもあるため，問題によっては要素の選択が解に大きな影響を与えてしまう．このため，要素の特性を知り適切な要素を選択することが不可欠であるが，局所化問題における要素選択についての知見は未だ不十分であると思われる．すなわち，一般化有限要素も含めた各種有限要素による局所化解析の結果をまとめ，解析に適した要素を検討することが必要である．

1.4 本研究の目的

本研究では，塑性局所化現象の高精度かつ汎用的な解析手法の確立を目的として，

局所化解析に適した有限要素の選択に関する検討，および二軸荷重を受ける板材，管材の局所化解析を行うための数値解析手法の提案を目的とする．

まず有限要素解析においてしばしば問題となる locking に強い新たな一般化平面要素を提案する．提案要素を含めた各種有限要素を用い，局所化解析に適した要素について論じる．これを通じ，高精度な局所化予測のための数値解析手法を確立するとともに，有限要素法を用いた局所化解析の新たな知見が得られることが期待できる．

次に工業的にも重要な板材，管材の塑性加工における破断予測に必要となる，二軸荷重制御による新しい有限要素解析手法を提案する．これらの問題の境界条件は，従来の有限要素解析手法では厳密に取り扱うことが出来ない．そこで二軸応力比を新たな制御パラメータとすることで，二軸荷重を独立して制御した有限要素解析が可能となることを示し，数値解析例を通じてその妥当性，有用性を確認する．また提案手法による解析結果を M-K モデルによるそれと比較することで，有限要素解析における局所化発生時の変形モード，成形限界ひずみと M-K モデルの予測結果との関連を検討する．以上を通じて，数値解析手法を用いた薄肉金属材料の破断予測における新たな方向性を見いだすことが，本研究の大局的な目的である．

(12)

1.5 本論文の構成

第1章の序論に続き，第2章ではまず本研究を通じて必要となる弾粘塑性モデルについての数学的なモデリングの枠組み，および J2 流れ則と Kuroda-Tvergaard の非法線則の概要を述べる．また，弾粘塑性モデルを数値解析で取り扱う際に必要となる応力積分の手法についても触れる．その後，速度型有限要素法の基礎式である速度型仮想仕事の原理についてまとめ，そこから速度型有限要素法を定式化する手順について詳細に示す．

第3章では，速度型有限要素法における解析精度の向上を目的として，二種類の新たな一般化平面要素を提案し，その定式化を試みる．提案した要素を用い，幾何学的非線形性を考慮した有限要素解析を実施し，提案する要素が従来の要素と比較して解の収束性，要素のゆがみに対するロバスト性の二点において優れた性能を有していることを示すとともに，提案要素のlockingに対する性質も検討する．

第4章では，第3章で提案した一般化平面要素を含む各種の有限要素を用いた塑性局所化解析を実施し，要素の選択が解析結果に及ぼす影響を検討する．要素剛性を評価する際の数値積分についても触れ，要素の次数，積分手法が局所化解析にどのような影響を及ぼすかを調査するとともに，局所化解析に適した有限要素についても考察する．

第5章では，板材の成形限界を予測するために，直交する二方向の荷重を独立して制御し，かつ塑性不安定現象を再現できる新しい解析手法を提案する．提案手法を用いて様々な荷重条件の下での塑性局所化解析を行い，成形限界ひずみを算出する．提案手法による解析結果と M-Kモデルによる結果の比較を行うとともに，M-K モデルの物理的な意味についても考察する．

第6章では，第5章で提案した手法を軸力，内圧複合荷重を受ける薄肉管材の解析へと適用し，複合荷重下での管材の塑性不安定解析を実施する．これにより提案手法が薄肉管材の複合荷重解析へも適用可能であること，薄肉管材の塑性不安定条件を正しく評価できること示す．

第7章では，各章で得られた成果を要約し，本論文の結論と今後の課題を述べる．

(13)

第 2 章弾粘塑性有限要素法の定式化

2.1 緒言

本章では，弾粘塑性有限要素法の定式化について述べる．現象論的塑性モデルは，

変形のひずみ速度依存性を考慮しない弾塑性モデルと，時間依存性を考慮した弾粘塑性モデルに大別できる．本研究では室温における金属の塑性加工を対象としており，

かつ変形速度も非常に遅い状態を想定している．一般に金属材料は室温付近で低速度の変形を受ける場合，弾性変形，塑性変形ともに速度依存性はほとんど見られないため，ここではその影響について考慮する必要はないものと考えられる．また変形速度が十分に遅いことから，変形に伴う発熱やそれによる熱膨張，熱硬化についても無視できる．このため速度非依存の弾塑性構成式で十分であると言えるが，数値解析における取り扱いの容易さから，本研究では後者の弾粘塑性モデルを選択した．

まず，本研究で用いる現象論的塑性構成式の一般論^[2-1][2-2]および本研究に用いた構成式の具体形について概説した後，弾粘塑性構成式を数値解析で取り扱う際に有効な応力積分手法である rate-tangent modulus 法^[2-3]について述べる．次に速度型有限要素法の基礎方程式となる速度型仮想仕事の原理と，これを基にした速度型有限要素法の定式化を示す．最後に，二次元解析に用いられる連続体要素についても触れる^[2-4][2-5]．

金属の塑性変形を記述する構成式の構築は，結晶塑性モデルに代表されるように物質のミクロレベルにおける力学的挙動に立脚するアプローチと，マクロレベルにおける挙動のみに注目し，その応力，ひずみ関係をモデル化する現象論的なアプローチに大別できる．前者は物体内部における微視的な物理挙動の理論的な予測に基づいている．金属材料は通常多結晶構造を有し，その内部におけるすべり変形が塑性変形の要因であると考えられるため，結晶塑性モデルは非常に合理的なモデルであり，実現象の再現性も高いとされている．しかしながら，特に多結晶モデルでは，非常に多くの結晶粒についての情報を保持，計算する必要があることから，計算時間，所要メモリのいずれの観点からも多大な計算コストが必要となり，現状における計算機性能では適用可能な問題の規模は限られている．

一方現象論的塑性モデルでは，物体のマクロレベルでの挙動のみを考えるため，ミクロスケールの力学を考慮したモデルと比べて非常に少ない計算コストで問題を扱うことができ，有限要素法などを用いた実用的な計算における利便性は高い．欠点としては，複雑な塑性挙動のモデリングに際してはモデリングに必要なパラメータが増加する傾向にあり，結果として実験的な材料定数のフィッティングが難しくなる点が挙げられる．

本研究では二種類の現象論的塑性モデルを用いた．具体的には，従来から広く用いられている構成式であるJ2流れ則と，Kuroda and Tvergaardによって提案された，降伏曲面上に尖り点が形成することを考慮した非法線則^[2-6]の二つを導入する．後者は

(14)

せん断帯の発生を伴うような大変形塑性解析において，結晶塑性モデルによる予測の再現性が高いとされており^[2-7]，計算コストを低減しつつも結晶塑性モデルと近い解析結果が得られることが期待できる．

大ひずみ状態を表現する塑性構成式を扱うに際しては，物体の幾何学的な形状変化，

すなわち変形と，それに応じて物体内部に生じる応力の定義が必要となり，これらは連続体力学を基礎として成立する^[2-8]．連続体力学では，物体は物質点の連続的な集合としてとらえられる．物質点は三次元Euclid空間において座標と一対一に対応づけられ，その体積は無限小であるとする．物質点とは，このように無限小の体積を仮定した点であるにもかかわらず，分子，原子のレベルに至ることなくその点の近傍の平均的な諸量を有するとみなす，仮想的な概念である．本章における内容の多くは一般座標系においても成立するものであるが，ここでは特に空間に固定された直交デカルト座標系を参照するものとする．すなわち，ベクトル（一階のテンソル）a，二階のテンソルB，四階のテンソルCとは，それぞれその成分a_i，B_ij，C_ijklと直交基底ベクトルe_iを用いて，以下のように表現されるものとする．

i

a ei

a= (2.1)

j i

Bije e

B= ⊗ (2.2)

l k j i

Cijkle e e e

C= ⊗ ⊗ ⊗ (2.3)

なお⊗はテンソル積を表し，今後特に断りのない限り総和規約を適用する．また，二階のテンソルBの大きさ B と偏差量B′をそれぞれ以下のように定義する．

( )= BijBij

≡ B B

B tr ^T (2.4)

( )^B ^I

B

B tr

3

−1

′≡ (2.5)

ここでIは二階の単位テンソル，右上付きTはテンソルの転置を表す．

2.2 現象論的塑性構成式 2.2.1 一般的な枠組み

座標x，速度v を持つ物質点と，これと近接するもう一つの物質点を想定する．二つの物質点を結ぶ微小線素ベクトルが dx，二点の相対速度が dv であるとき， dv と dxの間に成立する線形変換Lを考える．

x x x

v v d d

d ⋅ = ⋅

∂

= ∂ L (2.6)

(15)

x v

∂

≡ ∂

L (2.7)

Lは速度こう配テンソルと呼ばれている．Lは一般に非対称テンソルであり，対称部分と反対称部分に加算分解することができる．

W D

L= + (2.8)

( ^T)

2 1 L L

D≡ + (2.9)

( ^T)

2 1 L L

W≡ − (2.10)

Lの対称成分Dを変形速度テンソルもしくはストレッチングテンソル，反対称成分W を連続体スピンテンソルと呼ぶ．

金属材料では，微小弾性，有限塑性変形を仮定できることが多い．この前提のもとで，変形速度テンソルDを弾性成分と塑性成分に加算分解する．

p e

p

e D D N

D

D= + = +Φ& (2.11)

上付き e，p は弾性，塑性成分を表している．本研究では塑性変形のみにひずみ速度依存性を仮定しており，いわゆる粘弾性は含めていない（付録 B 参照）．N^pはD^pの方向を表すテンソルであり，係数Φ&は後述する超過応力関数である．N^pの選択には任意性があり，具体的な N^pの形式については2.2.3，2.2.4に示す．また，同様に連続体スピンテンソルWを以下のように加算分解できるものとする．

p

p ω Ω

ω+ = +Φ&

= W

W (2.12)

ω，W^pはそれぞれsubstructureスピン，塑性スピンと呼ばれ，Ω^pはW^pの方向を表す．

ωを用いることで，スカラーζ^{，ベクトル}^a，^{（二階の）テンソル} B の共回転速度がそれぞれ次のように与えられる．

ζ

ζ^o = & (2.13)

a a

a^o = &−ω (2.14)

ω ωB B B

B^o = & − + (2.15)

ここで上付き_oは共回転微分を，上付き・は物質時間微分を表す記号である．式(2.12) における塑性スピンW^pの選択には任意性がある．せん断帯形成を伴うような塑性大変形解析においては，スピンの選択が解析結果に大きく影響する場合があることが知られているが[2-9][2-10][2-11][2-12]，本研究ではW^p ≡0すなわちω≡Wを仮定する．これは Jaumann速度として知られる共回転速度であり，現在広く用いられている共回転速度

(16)

の一つである．客観応力速度とスピンテンソルの歴史的背景については付録C，Dに示している．

塑性に関する状態変数として，Cauchy 応力σと i 個の構造変数s_iを考える．s_iは，

スカラー，ベクトル，二階のテンソルのいずれかであり，個数は特に規定しない．構造変数の例としては，後述する相当塑性ひずみ（スカラー），塑性異方性を考慮する場合の異方性主軸（ベクトル），移動硬化モデルにおいて応力空間における降伏曲面の移動量を表す背応力（二階のテンソル）などが挙げられる．構造変数 s_iの発展方程式は一般に，

i si

s^o =Φ& (2.16)

で与えられ， s_iは i

so の方向を表す．

変形速度テンソルの弾性部分に対して次式が成立するものとする．

:De

=C

σo (2.17)

σo 1:

e -

C

D = (2.18)

ここで，C は弾性変形を表す四階の構成則テンソルである．C の選択は任意であり，

たとえば弾性異方性を含めることも可能であるが，微小弾性変形を仮定する場合，等方性を前提としたいわゆるHookeの法則を用いる場合が多く，本研究でもこれに倣っている．このとき構成則テンソルCの成分は次式で定義される．

jl ik kl

ij

Cijkl =λδ δ +2µδ δ (2.19)

ここで λ，µ はLaméの定数であり，それぞれ

( ν)( ν)

λ ν

2 1 1+ −

= E

(2.20)

(^E+ )⁼^G

= ν

µ 21 (2.21)

である．E，ν^，GはそれぞれYoung率（縦弾性係数），Poisson比，せん断（横）弾性係数である．またδ _ij ^はKroneckerのデルタ記号である．式(2.11)を式(2.17)に代入することで，弾粘塑性モデルにおける構成式の一般表示が

: p

:D C N

C Φ&

o = −

σ (2.22)

で与えられる．

弾塑性モデルにおける降伏関数は，弾粘塑性モデルでは次式で置き換えられる．

(17)

( ^, ^, ) ⁽ ⁾ ⁰

0 p

p ⎟⎟⎠ =

⎜⎜ ⎞

⎝

− ⎛

=

m

i g

J

f Φ

ε Φ

ε &

&

σ s (2.23)

J は相当応力であり，速度非依存の弾塑性モデルと同じ関数を用いる．相当塑性ひずみε^pは状態変数の一つであると解釈できるが，ここでは特に他の状態変数を区別して表記している．すなわち，相当応力 J は Cauchy 応力，相当塑性ひずみおよびその他の状態変数の関数である．gは基準となるひずみ硬化関数，Φ&₀は基準となる超過応力関数である．m はひずみ速度敏感性を表す指数であり，式(2.23)は_m_→₀の極限において弾塑性モデルにおける降伏関数と一致する．なお，相当塑性ひずみε^p^{は相当塑性} ひずみ速度ε_&^pの時間積分として，次式で定義されるものとする．ε_&^pの具体形は後述する．

∫

≡ ^pdt

p ε

ε & (2.24)

式(2.23)を超過応力関数Φ&について解くことで，次式が得られる．

( i ) ^m

g J

1

p p

0 ( )

, , ⎥

⎦

⎢ ⎤

⎣

= ⎡

ε Φ ε

Φ& & ^σ ^s (2.25)

すなわち，Φ&は相当応力Jおよびひずみ硬化関数gの汎関数である．弾塑性解析においては，物質点が降伏しているか否か，すなわち負荷・除荷の判定が必要となるが，

式(2.22)および(2.25)によって表現される弾粘塑性構成式を用いる場合，基本的にこれを必要としない．これは式(2.25)の右辺括弧内がフィルタ関数として作用し，降伏応力（に相当する応力）以下の応力域では塑性変形の速度が無視できる程度に小さくなるためである．弾塑性モデルを用いた数値解析においては，しばしば負荷・除荷の判定が非線形解析の反復過程における収束性を悪化させることがあるが，ここに示した弾粘塑性モデルを用いた場合この問題は生じず，この点が弾粘塑性モデルの数値解析上での利点の一つと言える．また前述のように，式(2.22)，(2.25)においてmを十分に小さな値としたとき，その挙動は速度非依存の弾塑性構成式に近づくが，その場合でも負荷・除荷の判定は必要ないため，速度依存性を考えない場合でも本モデルを用いることで数値解析の安定化を図ることができる．

2.2.2 J2流れ則

弾塑性モデルとして最も広く用いられているものの一つが，J2流れ則である．塑性ひずみ空間における塑性変形の方向（式(2.11)の N^p）が，応力空間における降伏曲面の外向き法線方向と一致すると考える時，これを関連流れ則と呼ぶ．J2流れ則は関連流れ則の一種であり，降伏関数としてMisesの降伏関数を適用したものである．J2流れ則を弾粘塑性モデルへ適用する際の定式化を以下に示す．

応力空間において，降伏関数を満たす点が作る超曲面を降伏曲面と呼ぶ．この法線

(18)

方向 Nⁿは，相当応力Jを応力で微分することで得られる．

σ

∂

≡ ∂J

Nn (2.26)

Mises の相当応力は，降伏曲面の移動を伴わない等方硬化のみを考慮する場合，

Cauchy応力の偏差成分σ′を用いて，

( ) ²

1

2 :

3 ⎟

⎠

⎜ ⎞

⎝

⎛ ′ ′

≡ σ σ

σ

J (2.27)

で定義される．これを式(2.26)に代入することで，降伏曲面の法線方向Nⁿが得られる．

J2 流れ則は関連流れ則であることから， Nⁿが塑性変形の方向 N^pと一致する．以上から，

J

J σ

σ

= ′

∂

= ∂

≡ 2

n 3

p N

N (2.28)

なる関係が得られる．これを式(2.22)に代入することで，J2 流れ則に基づく弾粘塑性構成式が得られる．なお，相当塑性ひずみ速度ε&^pは次式で与えられるものとする．

( ) ^Φ ( ) ^Φ

ε_&^p ≡ ^p ^p = ^& ^p: ^p = ^&

3 : 2

3

2 D D N N (2.29)

ここで，式(2.28)において常に N^p =3 2であることを用いた．

2.2.3 Kuroda-Tvergaardの非法線則

塑性変形の進行とともに降伏曲面に尖り点が形成されることは，以前から結晶塑性モデルで予測され，最近になって実験的にも確認されている^[2-13]．尖り点形成を再現するための現象論的モデルは多く提案されているが，その中でもよく知られているものの一つであり，その後の塑性モデルに広く影響を与えたのが，Christoffersen and Hutchinsonによって提案されたJ2コーナー理論^[2-14]である．J2コーナー理論は，そのままの形式で数値解析に用いるにはやや複雑な形式であるため，Simo は数値解析に適用し易いようこれを簡略化したモデルを示した^[2-15]．Kuroda and Tvergaardによって提案された現象論的塑性構成式は，この Simo のモデルを基により一般化されたモデ

ルである^[2-6]．

降伏曲面の法線方向nと接線方向mを表すテンソルを，それぞれ次式で定義する．

σ

σ ∂

⎟⎟ ∂

⎠

⎜⎜ ⎞

⎝

⎛

∂

≡ ∂J J

n (2.30)

(19)

( ) (n D)n D

n D n m D

− ′

′

− ′

≡ ′

:

: (2.31)

ここで相当応力 J は特に限定されず，任意に選択することが出来る．n，m はいずれも正規化されており，n は Nⁿと同じ方向を有する．降伏曲面の接線は無数に存在するが，式(2.31)によってnおよび変形速度の偏差成分D′と同一平面内，かつD′とmがなす角が鋭角となるよう m が決定される．なお n とD′がなす角度θ ^{は次式で与えら} れる．