機械の研究

Transcription

1 構造解析 :V-XFEM による疲労き裂進展シミュレーション 長嶋利夫 ( 上智大, 理研 ), 今井登 ( 理研 ). はじめに計算固体力学分野においては, 構造物内部のき裂やボイドなどの不連続な変位場を有限要素メッシュと独立に定義することができる拡張有限要素法 (X-FEM)[][] が提案され, 破壊力学問題における停留き裂解析, き裂進展解析に用いられるようになってきている []-[5] [0]-[].X-FEM を用いればき裂進展に伴うメッシュ再分割が不要になるので, き裂進展シミュレーションを従来の FEM よりも効率的に実施できる. 一方, 理化学研究所で開発されている VCAD[6]-[7] は, 現状の三次元 CAD が有する三次元形状情報に加え, 物性情報をもつボリュームデータを扱うことができる. 著者らは VCAD で用いられるボリュームデータとの入出力が可能な構造解析プログラム V-XD[8] の開発を行っている.V-XD においては, ボリュームデータを利用して仮想仕事の原理に基づく支配方程式を離散化するために必要な, 体積積分, 面積積分, さらには変位境界条件を満足させるための拘束条件方程式の導出が行われる. また, プログラム実装においては VCAD-Framework を利用する. 著者らは構造物の疲労き裂進展解析を効率的に実施するために,V-XD に X-FEM によるき裂進展解析機能を追加することを検討している. 本稿では, そのための基本技術となる X-FEM によるき裂進展解析手法および V-XD による応力解析の手法の概要を紹介する.. X-FEM によるき裂進展解析ここでは, き裂形状の表現にレベルセット法 [9] を用いた X-FEM( レベルセット X-FEM) による三次元弾性体内部の平面状のき裂のき裂進展解析手法を示す.. レベルセット関数によるき裂の表現三次元等方性弾性体内部に平面状のき裂を想定する. 図 に示すようにき裂面は - 平面上にあるものとし, き裂の前縁形状をレベルセット関数で近似表示する. すなわち, き裂近傍の領域において次式に示すような符号付き距離関数 fを定義する. T f ( ) = min sign( n( ) ( )) Γ () ここに,Γはき裂前縁, n() はΓ 上の点 におけるΓに直交する単位ベクトルである. fの絶対値は, 任意の点からき裂前縁 ( - 平面 ) へ降ろした投影点とき裂前縁との距離である. また, その符号はき裂面に対する位置を表している. すなわち, 正であればき裂前縁の外側 ( リガメント ), 負であればき裂面である. 式 () で定義されるレベルセット関数 fの値を節点で定義したうえで, これらを有限要素毎

2 に定義される内挿関数で補間することによって解析領域の任意の点でのレベルセット関数値を近似的に評価することができる. また, き裂前縁形状はレベルセット関数のゼロ等値線で表現することができる. f Crack surface Evaluation point r θ f Crack front:γ 図 レベルセット関数による三次元き裂形状の定義. 内挿関数解析領域は三次元 8 節点六面体一次要素を用いて要素分割されているものとする. X-FEM においては, き裂近傍の任意の点 における変位場 u h を, 次式を用いて近似的に表わす. u f 8 4 h h k ( ) = N ( ) u + N ( ) γ k ( f ( ), ) a + N ( ) H ( = C k = J h ( ) = 8 = N ( ) f ( ) ) b (.) (.) ここに N は通常の FEM で用いられる 8 個の節点から構成される有限要素の内挿関数,u, a J,b は節点自由度,C はき裂近傍の変位場を近似する 4 個の基底関数 γ (=,,4) の自由度をエンリッチする節点の集合,J はヘビサイド関数 H() の自由度をエンリッチする節点の集合, また,C J=φである. ここでは次式に示すような二次元き裂問題におけるき裂先端近傍の変位場の漸近解を再構成することができる 4 個の基底関数 γ (=,,4) を用いる. γ θ θ θ r cos( ), γ = r sin( ), γ = r sin( )sinθ, γ 4 = r θ cos( ) sinθ = () ここに,r,θ は評価点のき裂前縁からの距離と角度であり, レベルセット関数を用いて次式で表わされる. r f + = (4.) θ = arctan( / f ) (4.)

3 . 数値積分法 X-FEM においても通常の FEM と同様に仮想仕事の原理式に基づき離散化方程式を導出するので, その過程において, 要素単位での数値積分が必要になる. き裂前縁近傍のエンリッチされた節点を含む要素における内挿関数は, 多項式より高次の関数を含んでいるため, 通常の 8 節点六面体一次要素で用いられるようなガウスの 点公式を用いた数値積分では精度的に不十分である.X-FEM 解析においてき裂面が要素内に含まれる場合には, 領域を複数の単純な立体に分割したうえで部分領域毎にガウス積分が適用される. ここでの解析では, 剛性評価の処理を簡単化するためにエンリッチされた節点を含む要素において, ガウス積分点を増加させることによって剛性評価のための数値積分を行う..4 エネルギー解放率の評価方法 X-FEM 解析で用いられる解析モデルにおいては, 一般にき裂前縁形状と要素境界とは完全に一致しないので, き裂先端あるいは前縁でのエネルギー解放率を評価するために, 従来の FEM と組み合わせて用いられている仮想き裂閉口法を直接適用することは困難である. X-FEM の解析結果からエネルギー解放率を計算する場合, 領域積分表示された定義式を用いて計算する []-[5][0]-[]. すなわち, 図 に示すように, き裂前縁の任意の評価点においてき裂進展方向が 方向となるように局所直交座標系,, とその点を囲む任意 の領域 V を設定し, その点におけるエネルギー解放率を次式により計算する. u G = ( σ + V τ u u + ( τ + τ u + τ u u q + σ dv ) q w) + ( τ u + σ u + τ u q ) (5.) G = L c G αdc (5.) ここに α は評価点で, 領域 V の表面でゼロとなるスカラー関数,L c はき裂前縁,σ ij は応力,u j は変位,w はひずみエネルギー密度である. ここでは Sukumar ら [5] の研究に倣い, 積分領域 V として, エネルギー解放率の評価点を中心とし主軸方向が 方向となる直方体領域を設定する. その直方体内部に応力, ひずみ, 変位の評価点を設ける. 重み関数 q としては, 評価点でをとり, 直方体の側面でゼロとなるような次式で表される関数を用いる. = ( q,, ) (6) L L L ここに, L, L, L は局所座標系 方向に沿った直方体の辺の長さである. 式 (5.) に示す体積積分は, 直方体領域においては多数の積分評価点を用いて評価される.

4 Evaluation point Numerical ntegration L L L ntegration point 図 領域積分法によるエネルギー解放率の評価.5 き裂進展解析手法 X-FEM によるき裂進展解析 [0]- [] の流れを図 に示す. まず, 初期き裂形状を考慮せずに汎用の FEM 用プリプロセッサーを用いて三次元有限要素モデルを作成する. また, 別途平面状の初期き裂の前縁形状を点列とそれを結ぶ直線でモデル化しておく. つぎにプログラム処理により, 作成した FEM メッシュとき裂の前縁形状に対してエンリッチ節点を自動的に設定し,X-FEM 解析の入力データを作成する. 三次元 X-FEM 解析を実施し, き裂前縁でのエネルギー解放率, あるいは応力拡大係数を評価する. さらにこれらの値から, き裂進展則によりき裂進展量を評価し, き裂形状を更新する. このような処理を繰り返すことによりき裂進展解析が実行可能となる. FEM Modeler DFEM Data Geometry of crack front Delamination Modeler X-FEM Data X-FEM Solver Updating front tip Crack propagation law X-FEM Result Displacement, Stress/Strain Energy release rate 図 X-FEM によるき裂進展解析の流れ.6 き裂進展解析例一様分布引張荷重を受ける半楕円き裂 ( 短半径 0 mm, 長半径 00mm) を有する厚肉平板 ( 厚さ 50 mm, 長さ 000 mm, 幅 000mm, ヤング率 07GPa, ポアソン比 0.) を対象として弾性き裂進展解析を実施する. 解析においては対称性を考慮して /4 領域をモデル化する [].

5 図 4 に解析対象の形状と解析モデルを示す. 解析モデルとして構造的に 分割した メッシュ ( 節点数,49, 要素数 0,000) を用いる. ここでは繰り返し負荷される 00 MPa の一様引張り荷重によるモード の疲労き裂進展を想定して, き裂進展則を次式のように設定する. da m = CΔG (7) dn ここに a はき裂長さ,N は荷重サイクル,ΔG はエネルギー解放率範囲,C, m は実験的に決定される定数である. 与えられた半楕円形状の初期き裂形状 ( 短半径 0 mm, 長半径 00mm) に対して, 静解析を実施し, き裂前縁においてエネルギー解放率を計算し, 式 (7) を用いてき裂進展量を計算する. ここでは, 式 (7) のき裂進展則のパラメータとして m=.55,c= を与え ΔN=5,000 として静応力解析を 0 ステップだけ繰り返すことによりき裂進展解析を実施した. 初期状態と最終段階 (0 ステップ後 ) における対称面上のき裂形状を図 5に示す. 図 5 に示すように,X-FEM を用いることにより, 平板内を進展するき裂の前縁形状を有限要素メッシュと独立にモデル化でき, したがってメッシュ再分割をせずにき裂進展解析が実施できることがわかる. また, き裂最深部と表面部のき裂前縁位置と繰り返し回数との関係を図 6に示す. 最深部 (90 位置 ) では mm 程度, 表面部 (0 位置 ) では mm 程度, き裂前縁が進展した結果が得られている. 図 4 半楕円表面き裂を有する厚肉平板モデル nitial Final (After 0step) 図 5 半楕円表面き裂の進展解析結果

6 5 0 Postion [mm] 0 5 Postion [mm] Number of load cycles Number of load cycles 90 degree 0 degree 図 6 き裂前縁位置と繰り返し回数との関係. VCAD データ構造をベースとした構造解析プログラム (V-XD)[8] VCAD で扱われるボリュームデータは, オクツリー構造を有するセル情報を基本とし, セルと形状定義面との切断面 切断点などの幾何情報に加えて物性情報を保持することができる. このような特長を有するボリュームデータは, 人工物だけではなく自然物のモデル化にも有効であると考えられる. 著者らは, 損傷や欠陥を含む機械構造物や分布する特性 ( 密度, ヤング率など ) を有する生体骨などを直接モデル化して応力解析を効率的に実施することを目的として,VCAD が提供するオクツリー構造を有するボリュームデータを直接利用した応力解析手法の研究開発を行っている. 現在, そのプロトタイプ版の構造解析プログラム V-XD を開発中である.V-XD の特徴として : ボリュームデータを直接利用して応力解析を行うことができること. X-FEM( 拡張有限要素法 ) に基づいた応力解析手法であること. V-CAD-Framework を介して V-CAD の機能を利用したアプリケーションであること. などをあげることができる.. 解析手法離散化対象となる支配方程式は, 通常の FEM と同様に連続体についての仮想仕事の原理式である. 支配方程式を離散化するために, ボリュームデータのセル構造を利用する. 図 7 に示すようにボリュームデータはオクツリー構造を有するが, セル分割の基本となる最低次 ( レベル 0) のセルを 8 節点直方体形状の有限要素に対応させる. 形状定義面と交差するセル ( 境界セル ) に対しては, 拡張有限要素法 (X-FEM) の解析技術を用いて剛性評価を行う. この場合, 要素に対応するレベル 0 のセルに連結する高次のセルを用いて体積積分を行う.

7 表面荷重は境界セルの切断面を利用して面積積分することによって離散化される. ボリュームデータを利用した体積積分, 面積積分の概念を図 8に示す. また, 形状定義面において変位境界条件を与えるために, 複数個のセル切断点での拘束条件式を導出しグループ化したうえで, 節点自由度に関する多点拘束条件式を導出する. 図 9にその概念を示す. 結果として得られるシステム方程式は通常の FEM と同様な形式となるので, 既存の FEM で用いられる求解処理方法が利用可能である. Level 0 nner cell & Level 0 Node generation conventional FEM Leaf of the octree under the boundary cell with Level 0 図 7 VXD におけるモデル化の概念 (a) Volume integral (b) Area integral 図 8 境界セルにおける体積積分と面積積分 Group Group Group 図 9 境界セルにおける拘束条件評価点とそのグループ化

8 . V-XD の実装 V-XD の開発のプログラム実装においては VCAD Framework を利用する.VCAD Framework が提供する AP (Application Program nterface) 群を呼び出すことでボリュームデータ構造の詳細を意識することなく, ボリュームデータの生成, 変更, 取得などの処理が可能となる. 提供されている AP 関数の機能には 形状データを入力 VCARファイルの入出力 空間のセル分割 Kitta Cubeの生成 指定したセル位置の幾何情報の取得などがあり, 目的の処理に対する AP 関数を呼び出すことで, ボリュームデータを利用したアプリケーションを効率的に開発することができる.. 解析例機械部品を想定した形状を対象として開発中の V-XD を用いて応力解析を行った. 対象は工作機械の HSK ツールホルダであり,VOBJ 形式の表面形状データを,VCAD Framework の AP 関数を用いてボリュームデータ (VCAR データ ) に変換した. その際空間を 444 分割した.V-XD は VCAR 形式のデータに入力し, 形状を表示するとともに境界条件の付与, 解析条件の設定, さらには結果の表示などを行うことができる. 図 0に V-XD の実行画面の一例を示す. 境界条件としては, ツールシャンク部のテーパと端面を完全拘束し, 工具把持部内面の下面に下方向に圧力荷重を付与した. 解析結果から得られる変形分布およびミ -ゼス応力コンタを表示した例を図 に示す. 詳細は省くが, 通常の FEM 解析結果などと比較してほぼ妥当な結果が得られている. 図 0 V-XD による境界条件の設定

9 図 V-XD による解析結果 ( 変形およびミ - ゼス応力 ) の表示 4. おわりに本稿では X-FEM による疲労進展解析事例および V-XD による応力解析事例を示した. 今後は,V-XD に X-FEM によるき裂進展解析機能を組み込む予定である, 図 にボリュームデータでモデル化されたブラケット形状を V-XD に読み込み, さらに半楕円形状の表面き裂を設定した例を示す. このようなモデルについてのき裂進展解析を実施することが今後の課題である. 図 V-XD による半楕円表面き裂のモデル化

10 参考文献 [] Belytschko T, Black T, Elastic crack growth in finite elements with minimal remeshing, nt. j. numer. methods eng.999;45: [] Moës, N, Dolbow J, Belytschko T., A finite element method for crack growth without remeshing, nt. j. numer. methods eng. 999;46: -50. [] Moës N, Gravouil A, Belytschko T,, Non-planar D crack growth by the etended finite element and level sets-part : Mechanical model, nt. j. numer. methods eng. 00;5: [4] Gravouil A, Moës N, Belytschko T, Non-planar D crack growth by the etended finite element and level sets-part : Mechanical model, nt. j. numer. methods eng. 00;5: [5] Sukumar N., Moës N., Moran B., Belytschko T.: Etended finite element method for three-dimensional crack modeling, nternational Journal for Numerical Methods in Engineering 000; 48, [6] Kase K., Teshima Y., Usami S., Ohmori H, Teodosiu C., Makinouchi A., Volume CAD, Volume Graphics 00 Eurographics / EEE TCVG Workshop Proceedings,. Fujishiro, K. Mueller, A. Kaufman (eds.), 00; [7] Kase, K. Teshima, T., Usami, S. Kato, M., Yamazaki, S. to, M., Makinouchi, A.: Volume CAD-CW-complees based approach, Computer Aided Design, 005;7: [8] 長嶋, 今井 :VCAD-AP を用いた X-FEM に基づいた構造解析プログラム (VXD) の開発, 理研シンポジウムものつくり情報統合化研究 ( 第 5 回 ), V-CAD システム研究 - ものつくり ツールの構築と 研究 の基盤ツールへの展開 - 講演予稿集 (005-6), [9] Sethian HA, Level Set Methods and Fast Marching Methods: evolving interface in computational geometry fluid mechanics computer vision and material science, Cambridge, UK: Cambridge University Press, 999. [0] 長嶋, 寺原, 三浦 :X-FEM による三次元弾塑性応力解析システムの開発, 第 9 回計算力学講演会講演論文集,06-9 (006-), [] 作増, 長嶋, 三浦 :X-FEM による残留応力場中のき裂進展解析, 計算工学講演会論文集,(007-5), [] 作増, 長嶋 : X-FEMによる三次元き裂進展解析, 日本機械学会 007 年度年次大会講演論文集 Vol., 07-0 (007-9), -. [] 長嶋, 三浦 : 対称条件を考慮した X-FEM による三次元き裂解析, 機論 7-79, A(006),