XFEM解析システムと 複合材料構造解析への応用 Development of an analysis system based on the extended finite element method and application to structural analyses of composite 長嶋 利夫(上智大) Toshio NAGASHIMA (Sophia University) 2014年7月8日(火)第3回航空機CAE研究会 1 背景1 拡張有限要素法(XFEM) メッシュと独立にき裂を定義可能. 要素再分割なしにき裂進展解析が可能. 積層複合材料構造 配管 岩石 by courtesy of CRIEPI 2 背景2 アプローチ 長所 線形破壊力学 • • 比較的粗いメッシュが利用可能 線形解析 • • 初期き裂の設定が必要 エネルギー解放率、K値評価が必要 • • 応力に基づくき裂の発生 エネルギーに基づくき裂の進展 • • 詳細なメッシュが必要 (材料)非線形解析 Linear Elastic Fracture Mechanics: LEFM 結合力モデル Cohesive Zone model: CZM 短所 結合力き裂を用いた応力解析 ~ t2 ti = ti bi Γu x2 o x1 ui = ui t Imax Γσ A Γ ~ x2 ~ t1 Γ + GIC Γ ~ x1 O m2δ O m2δ t IImax u~ F 2 Γ− O GIIC (a) Mode I GIIC m1δ O m1δ F u~1 −t IImax (b) Mode II 3 背景3 NR法で収束解が得られる限界 TCT(Transverse Crack Test) 500 tc t 層間はく離:インターフェース要素 樹脂割れ:XFEMによる結合力き裂 averaged stress [MPa] 400 300 200 100 0 0 0.2 0.4 0.6 0.8 1 Displacement [mm] 非線形有限要素法解析で標準的に用いられるニュートンラプソン法による 陰的解法では収束解が得られない場合がある. 4 内容 XFEMによるき裂、損傷進展解析システムを開発し、 複合材料積層板へ適用する. • • • • • 背景、内容 解析手法の概要 開発システム 数値解析例 まとめ・今後の課題 5 本研究で用いる方法(方針) • 層間はく離など予め発生位置が分かっているき裂は通常の インターフェース要素でモデル化 • 発生位置がわかないき裂はXFEMでモデル化 • き裂面に結合力モデル(CZM)を導入可能 (解析対象によって様々なCZMを選択的に利用) • 原則としてヘビサイド関数だけを拡充し、 必要があればTIP要素を導入 • 陰解法あるいは陽解法を利用 • 二次元問題:三角形要素 準三次元問題:五面体要素 三次元問題:四面体要素 6 レベルセット関数による二次元き裂の表現 x φ >0 n + φ ψ crack tip n′ + x φ < 0 crack line ψ <0 ψ >0 Γ φ (x) = min x − x sign(n + ⋅ (x − x)) x∈Γ ψ (x) = min x − x sign(n′+ ⋅ (x − x)) x∈Γ′ 7 二次元X-FEMにおけるき裂のモデル化 ψ (x) φ (x) Enriched node with the Heaviside function Crack line Crack tip 8 レベルセット関数φによる要素の分類 normal φ + cut φ normal φ φ φ− + φ + − zero 0 φ+ φ zero − 0 φ − zero 00 − φ cut φ φ − φ− 0 0 0 cut φ+ φ φ+ − zero 0 + − φ + zero 0 φ+ φ− φ+ 0 9 レベルセット関数φ,ψによるCUT要素の分類 φ cut φ + cut + ψ ψ − φ + φ cut φ ψ − φ − φ − φ − + ψ− + φ − 0 + φ − − + ψ − φ + φ + − φ ψ+ ψ φ φ ψ+ φ− + ψ+ + − ligament ψ− − φ + ligament − φ tip ψ φ ψ + φ ligament + ψ tip ψ− − φ tip φ+ φ − ψ+ 0 ψ+ φ − 0 10 三次元XFEMにおける(内部)き裂のモデル化 x P φ P nS x ψ nC = φ ( x) min x∈S ∪ S EXT Γ S φ S EXT ψ ( x) x − x sign(n S ⋅ (x − x )) ψ= min x∈Γ x − x sign(nC ⋅ (x − x )) 11 レベルセット関数φによる テトラ要素の分類(一部) N P P P P P N Z N (2 1 1) P N N (3 0 1) P (1 0 3) N N N Z (2 1 1) N P P P P N P N (2 0 2) (2 0 2) P Z (2 1 1) P 12 半構造三次元五面体要素 Normal Cut Tip extrusion r3 x2 x1 Enriched node with the Heaviside function u h (= x, r3 ) 3 ∑ L (x)( N I I 1 B ( r3 )u I + NT ( r3 )u I +3 ) + ∑ FI (x)( N B ( r3 )a I + NT ( r3 )a I +3 ) I ∈J FI (x) ≡ LI (x)(H (φ (x)) − H (φ (x I ))) N B ( r3 ) ≡ (1 − r3 ) / 2, NT ( r3 ) ≡ (1 + r3 ) / 2 13 内挿関数 3 2D 3 u ( x) = ∑ LI u I + ∑ LI ( H (φ (x)) − H (φ (x I )))a I h = I 1= I 1 3 φ (x) = ∑ L I φ I I =1 3 u (x, r3 ) ∑ LI (x)( N B (r3 )u I + NT (r3 )u I +3 ) = Quasi-3D I =1 h 3 + ∑ LI (x)( H (φ (x)) − H (φ (x I )))( N B (r3 )a I + NT (r3 )a I +3 ) I =1 3 φ (x) = ∑ L I φ I N B ( r3 )= (1 − r3 ) / 2 NT ( r3 )= (1 + r3 ) / 2 I =1 4 3D 4 u ( x) = ∑ LI u I + ∑ LI ( H (φ (x)) − H (φ (x I )))a I h = I 1= I 1 4 φ (x) = ∑ LI φI I =1 14 解析手法 離散化式 (t ) + CU (t ) + Q(U(t )) = MU F(t ) 陰的静解析 Q(U(t )) = F(t ) K T ∆U = F − Q(U) 時間積分法:Newmark-β法(β=0.25, γ=0.5) (t + ∆= (t ) + ∆t γ U (t + ∆t ) + (1 − γ )U (t ) U t) U 陰的動解析 { } (t ) + ( ∆t )2 1 − β U (t ) + β U (t + ∆t ) = t ) U(t ) + ∆tU U(t + ∆ 2 2 − CU (K T + (1/ β ( ∆t ) )M + (γ / β∆t )C)∆U = F − Q − MU 陽的動解析 時間積分法:中央差分法 1 1 1 1 M C U F Q M U U CU(t − ∆t ) + ( t + ∆ = t ) ( t ) − ( t ) + (2 ( t ) − ( t − ∆ t )) + 2 ( ∆t )2 2∆t 2 ∆ t ( ∆t ) 対角化 条件つき安定: ∆tcrit < 2 ωmax ( 1+ ξ 2 − ξ ) 15 XFEM解析における集中質量マトリクス Uu Na ] U a 3 u ( x) = ∑ LI u I + ∑ LI (H (φ (x)) − H (φ (x I )))a I ≡ [ Nu h I ∈J =I 1 M uu ≡ ∫ ρ NTu N u dA, M ua ≡ ∫ ρ NTu N a dA, A Consistent mass matrix M cons M uu = T M ua A M ua M aa Lumped mass matrix M lump 0 diag (M uu ) = M 0 diag ( ) aa Menouillard, T., Rethore, J., Combescure A., Bung H.: Efficient explicit time stepping for the eXtended Finite Element Method (X-FEM), International Journal for Numerical Methods in Engineering, 68, (2006), pp. 911-939. 16 2D-XFEM code Code NLXT2D Latest version 3.140609 Development language ANSI-C Discretization method eXtended Finite Element Method Analysis type Static, Implicit Dynamic, Explicit Dynamic Element type Three-node triangle element (Constant Strain Triangle) Four-node interface element considering CZM Material type Isotropic, Orthotropic Enrichment type Heaviside Method use to solve the system equation Direct method: Skyline method PARDISO (Intel Math Kernel Library) Evaluation J-integral and SIF Domain Integral Method, M-integral 17 Quasi 3D-XFEM code Code NLXP3D Latest version 3.140707 Development language ANSI-C Discretization method eXtended Finite Element Method Analysis type Static, Implicit Dynamic, Explicit Dynamic Element type Six-node pentahedral element Six-node interface element considering CZM Material type Isotropic, Orthotropic Enrichment type Heaviside Method use to solve the system equation Direct method: Skyline method PARDISO (Intel Math Kernel Library) Evaluation J-integral and SIF N.A. 18 3D-XFEM code Code NLXT3D Latest version 3.140702 Development language ANSI-C Discretization method eXtended Finite Element Method Analysis type Static, Implicit Dynamic, Explicit Dynamic Element type Four-node tetrahedral element Six-node interface element considering CZM Material type Isotropic, Orthotropic Enrichment type Heaviside, Asymptotic basis Method use to solve the system equation Direct method: Skyline method PARDISO (Intel Math Kernel Library) Evaluation J-integral Domain Integral Method 19 3D-XFEMによる解析例 20 開発システム Crack definition file Finite element model file *.dat *.crk NLXT2D/NLXP3D/NLXT3D *.out Analysis results file *.log Log file *.vtf Analysis results file for visualization 21 Material properties and cohesive zone parameters Laminate Cohesive zone EL [GPa] 161 GIC [N/mm] 0.2 ET [GPa] 11.38 GIIC[N/mm] 1.0 GLT [GPa] 5.17 σImax [MPa] 60 GTT [GPa] 3.98 σIImax [MPa] 90 νLT 0.32 kI [N/mm3] 1.0 x106 νTT 0.436 kII [N/mm3] 1.0 x106 ρ [Ns2/mm4] 1.6 x 10-9 α 1.0 2 Conditions for crack initiation Conditions for crack propagation 2 < σ > τ 1 + = σ σ Im ax II m ax GI G IC α G + II G IIC α = 1 22 TCT(Transverse Crack Test) GII = σ2 tc t 4 EL t − tc σ cr = 2 GII EL Tension load t − tc ttc Tension load tc t Developments in the science and technology of composite materials,ECCM4 (1990) σ cr = 2 1.0 ×1.61×105 4−2 = 401.24 MPa 4 ×2 23 Numerical results for TCT specimens S11 [MPa] d Number of nodes 23958 Number of pentahedral elements 30720 Number of interface elements 7680 24 Numerical results for TCT specimens Displacement vs. average stress 700 2D/Static 2D/Explicit 3D/Static 3D/Explicit averaged stress [MPa] 600 500 Ref. σ cr = 401.24 MPa 400 300 200 100 0 0 0.2 0.4 0.6 0.8 1 Displacement [mm] 25 Numerical results for TCT specimens (Tetra model) S11 [MPa] d Number of nodes 23958 Number of tetrahedral elements 92160 Number of interface elements 7680 26 Numerical results for TCT specimens Displacement vs. average stress (2D/Q3D/3D) 700 2D Q3D(PENTA) 3D(TETRA) averaged stress [MPa] 600 500 400 300 200 100 0 0 0.2 0.4 0.6 0.8 1 Displacement [mm] 27 Energy Balance (MS 104) Energy [Nmm] 3 10 4 2.5 10 4 2 10 4 1.5 10 4 Uint Ukin 1 104 5 10 3 0 0 5 10-2 1 10-1 1.5 10 -1 2 10-1 time [sec] 28 Fracture analyses for composite specimens Hallett, S. R. et al, Modelling the interaction between matrix cracks and delamination damage in scaled quasiisotropic specimens, Composites Science and Technology 68 (2008) 80–89. Jiang,W-G, et al, A concise interface constitutive law for analysis of delamination and splitting in composite materials and its application to scaled notched tensile specimens, Int. J. Numer. Meth. Engng 69 (2007) 1982–1995. 29 CFRP laminate specimen under tension load Length:120 mm Width: 32 mm Thickness: 4 mm [45/90/-45/0]s Number of nodes 31944 Number of pentahedral elements 30720 Number of interface elements 23040 30 Modeling of matrix cracks by enrichment nodes 45 deg Layer 90 deg Layer -45 deg Layer 31 Stress distributions and deformation S11 [MPa] 800 Implicit/PENTA Explicit/PENTA d Average stress [MPa] 700 600 500 400 300 200 100 0 0 0.5 1 1.5 2 LPD [mm] Number of nodes 31944 Number of pentahedral elements 30720 Number of interface elements 23040 32 Comparison of numerical results 800 Implicit/PENTA Explicit/PENTA Average stress [MPa] 700 600 500 400 Hallett, S. R. et al, Modelling the interaction between matrix cracks and delamination damage in scaled quasi-isotropic specimens, Composites Science and Technology 68 (2008) 80–89. 300 200 100 0 0 0.5 1 1.5 2 LPD [mm] Iarve, E. V. et al, Mesh-independent matrix cracking and delamination modeling in laminated composite, IJNME (2011) 33 FE models with various scales m Stacking sequence dimension Thickness [mm] 1 [451/901/-451/01]s 30 x 8 1 2 [452/902/-452/02]s 60 x 16 2 4 [454/904/-454/04]s 120 x 32 4 8 [458/908/-458/08]s 240 x 64 8 m =1 m=2 m=4 m=8 Number of nodes 31944 Number of pentahedral elements 30720 Number of interface elements 23040 34 Comparison of numerical results for various scaling factors 1000 894 Averaged stress [MPa] 800 631 600 480 400 383 m=1 m=2 m=4 m=8 200 0 0 0.5 1 1.5 2 LPD [mm] MS 104 , Δt=1.0x10-7, t=0.1sec, 40mm/sec Hallett, S. R. et al, Modelling the interaction between matrix cracks and delamination damage in scaled quasi-isotropic specimens, Composites Science and Technology 68 (2008) 80–89. 35 Finite element models of composite specimen with an open hole under tension load Number of nodes 6792 Number of pentahedral elements 6280 Number of interface elements 4710 36 Stress distribution and Damage propagation S11 [MPa] d 37 FE models with various thicknesses t=4 mm t=2 mm Input data for NLXP3D t=1 mm 4 0.0 0.5 1 0.5 1 0.5 1 0.5 1 0.0 -45.0 90.0 45.0 1 1 1 1 2 2 2 2 38 Comparison of numerical results for various dimensions 500 explicitH4 explicitH2 explicitH1 averaged stress [MPa] 400 300 200 100 0 0 0.2 0.4 0.6 0.8 1 Strain [%] Jiang,W-G, et al, A concise interface constitutive law for analysis of delamination and splitting in composite materials and its application to scaled notched tensile specimens, Int. J. Numer. Meth. Engng 69 (2007) 1982–1995. 39 400 explicitH4 350 averaged stress [MPa] 300 250 200 150 100 50 0 0 0.2 0.4 0.6 0.8 1 Strain [%] 40 まとめ・今後の課題 XFEMに基づき開発した,二次元,準三次元,および 三次元き裂進展解析システムの概要と それらを用いたCFRP積層構造についての数値解析 例を示した. 今後の課題 任意の位置で発生するき裂のモデル化 繊維破断のモデル化(連続体損傷モデルの導入) リスタート機能 高速化(並列化、低減積分、陰解法と陽解法との組み合わせ、 合理的なマススケール設定) 非線形性(幾何学的非線形)の考慮 41 実問題への適用および検証 ご清聴ありがとうございました. 42
© Copyright 2024 ExpyDoc