JP2016004543A - 有限要素解析装置、該方法及び該プログラム - Google Patents
有限要素解析装置、該方法及び該プログラム Download PDFInfo
- Publication number
- JP2016004543A JP2016004543A JP2014126490A JP2014126490A JP2016004543A JP 2016004543 A JP2016004543 A JP 2016004543A JP 2014126490 A JP2014126490 A JP 2014126490A JP 2014126490 A JP2014126490 A JP 2014126490A JP 2016004543 A JP2016004543 A JP 2016004543A
- Authority
- JP
- Japan
- Prior art keywords
- finite element
- stress
- model
- equation
- strain
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Abstract
【課題】極低サイクル疲労破壊をシミュレートできる有限要素解析装置、方法及びプログラムを提供する。【解決手段】有限要素解析は、有限要素モデルの剛性行列を作成し、剛性行列と前記解析条件とで規定される剛性方程式の解に基づいて、有限要素モデルの加工硬化特性を算出する特性計算を行う応力ひずみ算出手段と、有限要素モデルの内力と外力とが釣り合うまで、応力ひずみ算出手段に特性計算を繰り返し行わせる制御手段とを備える。応力ひずみ算出手段は、剛性方程式の解に基づく試行弾性歪みを用い、繰り返し負荷を受けている部分であって、破壊が発生していない部分の応力とひずみを算出する加工硬化モデルと、破壊が発生している部分の破壊進行度と応力を算出する材料損傷モデルとにより応力とひずみとを連成させた材料損傷方程式を作成し、材料損傷方程式の解を求め、求めた解に基づいてコンシステント接線剛性を作成し、剛性行列を作成する。【選択図】図４
Description
本発明は、加工（移動）硬化モデルと損傷モデルとを用いた有限要素解析方法に関する。
近年，産業分野において機械部品の設計にＣＡＥ（Ｃｏｍｐｕｔｅｒ Ａｉｄｅｄ Ｅｎｇｉｎｅｅｒｉｎｇ）解析が行われている。この解析手法は、解析対象を要素に分割してモデル化し、物体に作用する外力により各要素に生ずるひずみや応力を逐次計算して、破壊に至る過程を計算する。
例えば、金属板から成形加工を経て製造させる金属部品のＣＡＥ解析用装置が、特許文献１に示されている。この装置は、加工部位の塑性ひずみをＣＡＥ解析の初期状態に反映させることで、成形加工時に生じる塑性変形による影響を考慮した非線形解析を行うＣＡＥ解析装置である。
ここで、機械部品の強度設計における問題の１つとして、極低サイクル疲労破壊問題がある。金属等の物体は、作用する外力がある大きさを超えると、構成材料に局部的な破壊（ボイド）が生じて、変形が徐々に大きくなる。このような物体の変形挙動（損傷）を記述する損傷モデル（力学構成式）として、例えば、Ｇｕｒｓｏｎモデルがある。また、繰り返し大負荷を受けた場合の物体の変形挙動を記述する移動（加工）硬化モデルとして、例えば、Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｃｋ則といったモデルがある。このモデルは、弾塑性変形を伴う繰返し大負荷を受けると現れるＢａｕｓｃｈｉｎｇｅｒ効果を表現することができる。
極低サイクル疲労破壊問題をシミュレートするためには、繰り返し負荷がかかる過程で生じる材料自体の損傷が、この繰り返しの過程における物体の変形挙動に及ぼす影響を考慮する必要がある。しかし、Ｇｕｒｓｏｎモデルでは、引張り負荷下で損傷が進展するが、逆に圧縮状態では損傷が回復してしまうため、繰り返し負荷を与えるような場合の損傷を表現するのに適切ではない。また、Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｃｋ則といったモデルでは、Ｂａｕｓｃｈｉｎｇｅｒ効果を表現することはできるが、繰り返しかけられる負荷による材料自体の損傷の影響を物体の変形挙動として表現することはできない。
そこで、本発明は、解析対象物の有限要素解析モデルを用いて、極低サイクル疲労破壊をシミュレートできる有限要素解析装置、有限要素解析方法及び有限要素解析プログラムを提供することを目的とする。
上述の目的を達成するために、本発明に係る有限要素解析装置は、非線形有限要素解析を行うための有限要素解析装置であって、解析対象の有限要素モデルと解析条件とを取得する取得手段と、前記有限要素モデルの剛性行列を作成し、作成した剛性行列と前記解析条件とで規定される剛性方程式の解に基づいて、前記有限要素モデルの加工硬化特性を算出する特性計算を行う応力ひずみ算出手段と、前記有限要素モデルの内力と外力とが釣り合うまで、前記応力ひずみ算出手段に前記特性計算を繰り返し行わせる制御手段とを備え、前記応力ひずみ算出手段は、前回算出した前記剛性方程式の解に基づく試行弾性歪みを用い、繰り返し負荷を受けている部分であって、破壊が発生していない部分の応力とひずみを算出する加工硬化モデルと、破壊が発生している部分の破壊進行度と応力を算出する材料損傷モデルとにより応力とひずみとを連成させた材料損傷方程式を作成し、作成した材料損傷方程式の解を求め、求めた解に基づいてコンシステント接線剛性を求め、求めたコンシステント接線剛性から前記剛性行列を作成することを特徴とする。
そして、本発明に係る有限要素解析方法は、非線形有限要素解析を行うための有限要素解析装置で用いられる有限要素解析方法であって、解析対象の有限要素モデルと解析条件とを取得する取得ステップと、前記有限要素モデルの剛性行列を作成し、作成した剛性行列と前記解析条件とで規定される剛性方程式の解に基づいて、前記有限要素モデルの加工硬化特性を算出する特性計算を行う応力ひずみ算出ステップと、前記有限要素モデルの内力と外力とが釣り合うまで、前記応力ひずみ算出ステップに前記特性計算を繰り返し行わせる制御ステップとを備え、前記応力ひずみ算出ステップは、前回算出した前記剛性方程式の解に基づく試行弾性歪みを用い、繰り返し負荷を受けている部分であって、破壊が発生していない部分の応力とひずみを算出する加工硬化モデルと、破壊が発生している部分の破壊進行度と応力を算出する材料損傷モデルとにより応力とひずみとを連成させた材料損傷方程式を作成し、作成した材料損傷方程式の解を求め、求めた解に基づいてコンシステント接線剛性を求め、求めたコンシステント接線剛性から前記剛性行列を作成することを特徴とする。
また、本発明に係る、コンピュータに実行させるための有限要素解析プログラムは、非線形有限要素解析を行うための有限要素解析装置で用いられる有限要素解析プログラムであって、解析対象の有限要素モデルと解析条件とを取得する取得手段と、前記有限要素モデルの剛性行列を作成し、作成した剛性行列と前記解析条件とで規定される剛性方程式の解に基づいて、前記有限要素モデルの加工硬化特性を算出する特性計算を行う応力ひずみ算出手段と、前記有限要素モデルの内力と外力とが釣り合うまで、前記応力ひずみ算出手段に前記特性計算を繰り返し行わせる制御手段として、コンピュータを機能させ、前記応力ひずみ算出手段は、前回算出した前記剛性方程式の解に基づく試行弾性歪みを用い、繰り返し負荷を受けている部分であって、破壊が発生していない部分の応力とひずみを算出する加工硬化モデルと、破壊が発生している部分の破壊進行度と応力を算出する材料損傷モデルとにより応力とひずみとを連成させた材料損傷方程式を作成し、作成した材料損傷方程式の解を求め、求めた解に基づいてコンシステント接線剛性を求め、求めたコンシステント接線剛性から前記剛性行列を作成することを特徴とする。
このような構成の有限要素解析装置、有限要素解析方法及び状有限要素解析プログラムによれば、加工硬化モデルと損傷モデルとを連成（強連成）したモデルを用いるので、相互に影響し合う非破壊部分と、破壊部分の剛性変化を考慮して、破壊を算出することができる。よって、サイクル負荷による破壊発生部位と破壊状況とを高精度にシミュレートすることが可能となる。
また、上述の有限要素解析装置において、前記応力ひずみ算出手段は、前記材料損傷方程式の解に基づいて、前記有限要素モデルの各要素の損傷値を求め、前記損傷値が所定の閾値を上回る要素は破壊されたと判断することが好ましい。
この構成によれば、相互に影響し合う非破壊部分と、破壊部分の剛性変化を考慮して、有限要素モデルの要素の破壊を算出することが可能となり、シミュレーションにおいて損傷値を求めることができるので、き裂が生じているかを判断することができ、破壊進展のシミュレーションを行うことが可能となる。
また、上述の有限要素解析装置において、前記加工硬化モデルは、分割背応力形Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｋ則による加工硬化モデルであることが好ましい。
この構成によれば、分割背応力モデルと損傷モデルとを連成（強連成）したモデルを用いるので、相互に影響し合う非破壊部分と、破壊部分の剛性変化を考慮して、破壊を算出することができる。
また、上述の有限要素解析装置において、前記コンシステント接線剛性は、完全陰解法により求めることが好ましい。
この構成によれば、完全陰解法による解析であるため、収束性良く、安定に計算を行うことが可能となる。
また、上述の有限要素解析装置において、前記制御手段は、前記有限要素モデルの内力と外力とが釣り合うと、前記解析条件のうちの強制変位又は外力を、所定の増分値だけ増加させて、前記応力ひずみ算出手段に前記特性計算を繰り返し行わせることが好ましい。
この構成によれば、強制変位又は外力を変えながら繰り返し前記特性計算を行わせることができるので、サイクル負荷による破壊発生部位と破壊状況を高精度に算出することが可能となる。
また、上述の有限要素解析装置において、前記有限要素モデルを表示する表示部を、更に備え、前記制御手段は、前記応力ひずみ算出手段が破壊されたと判断した要素を除いて、各要素をそれぞれの応力に応じた表示態様で、前記有限要素モデルを前記表示部に表示させることが好ましい。
この構成によれば、破壊シミュレーション結果を、表示することが可能となる。
また、上述の有限要素解析装置において、前記剛性算出手段は、作成した材料損傷方程式の解を求める前に、前記有限要素モデルの各要素の応力に対して、共回転変換を行うことが好ましい。
この構成によれば、Ｕｐｄａｔｅｄ−Ｌａｇｒａｎｇｅ形の有限変形モデルに対応することが可能となる。前記材料損傷方程式（材料損傷モデル）を有限変形モデルに対応させることで、繰り返し大負荷を受ける極低サイクル疲労破壊の中でも、変形が大きなものに対応することが可能になる。
このような構成の有限要素解析装置、有限要素解析方法及び有限要素解析プログラムでは、解析対象物の有限要素解析モデルを用いて、極低サイクル疲労破壊をシミュレートできる。
＜実施形態１＞
以下、本発明に係る実施形態を図面に基づいて説明する。
＜構成＞
図１は、有限要素解析装置の構成を示すブロック図である。図１において、有限要素解析装置１は、演算処理部１１、入力部１２、出力部１３、記憶部１４、及び、バス１９を備えて構成される。
以下、本発明に係る実施形態を図面に基づいて説明する。
＜構成＞
図１は、有限要素解析装置の構成を示すブロック図である。図１において、有限要素解析装置１は、演算処理部１１、入力部１２、出力部１３、記憶部１４、及び、バス１９を備えて構成される。
演算処理部１１は、例えば、マイクロプロセッサ等を備えて構成され、機能的に、制御部１１１、初期設定部１１２、及び、応力ひずみ解析部１１３を備える。
制御部１１１は、初期設定部１１２、及び、応力ひずみ解析部１１３を制御して有限要素解析の演算を行わせ、また、入力部１２、出力部１３、及び、記憶部１４を制御する。
初期設定部１１２は、解析対象である有限要素モデルの解析処理（シミュレーション）を行う際の初期値を設定する。初期値は、例えば、材料定数、境界条件、ステップ増分値（変位、外力）などである。
応力ひずみ解析部１１３は、ステップ増分値に応じた、応力、ひずみ、損傷値等の物理量（加工硬化特性）を算出する。具体的には、損傷モデルと加工（移動）硬化モデルとを、強連成し完全陰的に解くための方程式（アルゴリズム）を作成し、応力ひずみ挙動を求める。応力ひずみ解析部１１３の処理の詳細は、＜動作＞の項で説明する。
入力部１２は、有限要素解析装置１の演算開始指示等の各種コマンドやシステムパラメータ等の各種データを有限要素解析装置１に入力する機器であり、例えば、キーボードやマウス等である。出力部１３は、入力部１２から入力されたコマンド及びデータや本有限要素解析装置１の演算結果等を出力する機器であり、例えばＬＣＤ（Ｌｉｑｕｉｄ Ｃｒｙｓｔａｌ Ｄｉｓｐｌａｙ）、有機ＥＬ（Ｅｌｅｃｔｒｏ Ｌｕｍｉｎｅｓｅｎｃｅ）ディスプレイ等の表示装置である。
記憶部１４は、解析処理を実行する上で必要な各種プログラムや各種データを記憶する装置であり、例えば、ＲＯＭ（Ｒｅａｄ Ｏｎｌｙ Ｍｅｍｏｒｙ）やＥＥＰＲＯＭ（Ｅｌｅｃｔｒｉｃａｌｌｙ Ｅｒａｓａｂｌｅ Ｐｒｏｇｒａｍｍａｂｌｅ Ｒｅａｄ Ｏｎｌｙ Ｍｅｍｏｒｙ）等の不揮発性記憶素子、ＲＡＭ（Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）等の揮発性記憶素子およびその周辺回路等を備えて構成される。
また、記憶部１４は、有限要素モデル記憶部１４１、解析条件記憶部１４２、及び、物理量記憶部１４３を備える。
有限要素モデル記憶部１４１は、解析対象の有限要素モデルを記憶する。
解析条件記憶部１４２は、有限要素モデル記憶部１４１に記憶されている有限要素モデルのシミュレーションに必要なデータ、例えば、材料定数、ステップ増分値等を記憶する。
物理量記憶部１４３は、応力ひずみ解析部１１３が算出した物理量を記憶しておく。この物理量は、ステップ増分値が加算される毎に算出され、記憶される。この物理量に基づいて、有限要素モデルの変形挙動が画像として表示される。
これら演算処理部１１、入力部１２、出力部１３、及び、記憶部１４は、データを相互に交換することができるようにバス１９にそれぞれ接続される。
＜動作＞
次に、本実施形態の動作について説明する。
＜動作＞
次に、本実施形態の動作について説明する。
図２は、本実施形態の有限要素解析装置１の解析処理のフローチャートである。尚、フローチャート内の括弧書きの番号は、後述の式の番号を示す。例えば、（１６）は、式（１６）を示す。
まず、ユーザは、入力部１２を介して、解析処理の開始を指示するコマンド等を入力する。入力部１２を介して、解析処理の開始を指示するコマンドが入力されたことを検知した演算処理部１１の制御部１１１は、初期設定部１１２に初期値の取得を依頼する。依頼を受けた初期設定部１１２は、出力部１３に初期設定の入力を要求するメッセージを表示する。
ユーザは、解析対象の有限要素モデルを作成し、入力部１２を介して、有限要素モデル記憶部１４１に記憶させる（ステップＳ１０）。そして、ユーザは、入力部１２を介して、境界条件として、解析対象モデルに対して、どの部位をどのように固定するか、あるいはどのくらい強制的な変位を与えるか、などの拘束条件（実施形態では、境界変位ｕＢＣ）、及び、どの部位にどのような方向でどのくらい力を加えるかなどの荷重条件（実施形態では、境界外力ｆＢＣ）を、入力する（ステップＳ１１）。また、ユーザは、入力部１２を介して、解析対象モデルの材料定数を入力し（ステップＳ１２）、変位ｕ、及び、外力ｆのステップ増分値（Δｕ、Δｆ）を入力する（ステップＳ１３）。尚、実施形態では、説明の便宜上、変位ｕ及び外力ｆの両方を変化させることとして説明するが、変位ｕ及び外力ｆのうちのいずれか一方を変化させることとしてもよい。つまり、変位ｕを変化させる方法、外力ｆを変化させる方法、変位ｕ及び外力ｆを変化させる方法の３つのうちのいずれかでシミュレーションすることとすればよい。いずれを変化させるかは、実施するシミュレーションに応じて決定される。
ステップＳ１１〜Ｓ１３において入力されたデータは、解析条件記憶部１４２に記憶される。尚、ステップＳ１０〜ステップＳ１３のデータは、予め有限要素モデル記憶部１４１及び解析条件記憶部１４２に記憶させておき、ユーザが選択することとしてもよい。
初期設定部１１２による初期設定が終了すると、制御部１１１は、現在のステップｎ（最初はｎ＝１）、つまり、応力ひずみ解析部１１３に渡す変位ｕｎ、及び、外力ｆｎを算出する（ステップＳ１４）。制御部１１１は、ステップごとにステップ増分値ずつ増加させて、予め定めたステップの数だけ応力ひずみ解析部１１３に物理量を繰り返し算出させる。従って、現在のステップｎにおける変位ｕｎ、外力ｆｎは、前回の変位ｕｎ−１、外力ｆｎ−１にそれぞれΔｕ、Δｆが加算されて求められる。
制御部１１１は、変位ｕｎ、及び、外力ｆｎを応力ひずみ解析部１１３に渡して、応力ひずみ解析処理を行わせて物理量、すなわち、応力σｎ＋１、等方硬化内部変数Ｒｎ＋１、全体背応力βｎ＋１、分割背応力βｋｎ＋１、損傷値Ｄｎ＋１を求め（ステップＳ１５）、求めた物理量を物理量記憶部１４３に記憶させる（ステップＳ１６）。
制御部１１１は、変位ｕｎ、外力ｆｎそれぞれが、境界変位ｕＢＣ、境界外力ｆＢＣとなるまで（ステップＳ１７：ＮＯ）、ステップＳ１４〜Ｓ１６までの処理を繰り返し実行する（ステップＳ１７：ＹＥＳ）。
尚、極低サイクル疲労問題等のシミュレーションを行う場合は、さらに、ステップＳ１７においてＹＥＳとなった後に、ステップ増分値（Δｕ、Δｆ）を変えて、ステップＳ１４〜Ｓ１７までの処理を繰り返す。例えば、負荷を解除した場合は、Δｕ及びΔｆをマイナスとして、ステップＳ１４〜Ｓ１７までの処理を繰り返す。そして、次のサイクルでは、Δｕ及びΔｆを、プラスとしてステップＳ１４〜Ｓ１７までの処理を繰り返し、Δｕ及びΔｆを、マイナスとしてステップＳ１４〜Ｓ１７までの処理を繰り返す。
シミュレーションが終了すると、制御部１１１は、ユーザの指示に応じて、物理量記憶部１４３に記憶されている物理量に基づいて、時系列に、有限要素モデルを出力部１３に表示する。例えば、有限要素モデルの各要素を、要素の応力値に応じた色で表示するなどである。また、破壊されたと判断された要素は、描画されないようにすることで、破壊の進展の様子が表示されることになる。
次に、図３を用いて、図２のステップＳ１５の応力ひずみ解析処理について説明する。図３は、応力ひずみ解析処理のフローチャートである。
応力ひずみ解析部１１３は、制御部１１１から渡されたステップｎの変位ｕｎ、及び、外力ｆｎを用いて、以下の式（Ｇ１）で示す剛性方程式を作成する。
[Ｋ]は全体剛性マトリックス、｛Δｕ｝は節点変位ベクトル、｛Δｆ｝は節点力ベクトルを示す。応力とひずみの関係は、非線形であるため、剛性マトリックスは応力やひずみによって変化する。従って、最初のステップ１では、予め定めておいた剛性マトリックスが用いられ、ステップ２からは、後のステップＳ２８で作成された剛性マトリックスが用いられる。
次に、応力ひずみ解析部１１３は、作成した剛性方程式の解を算出し、算出した解から各要素のひずみを算出する（ステップＳ２０）。
次に、応力ひずみ解析部１１３は、各要素のＣａｕｃｈｙ応力に対して、式（１６）で示される共回転変換を行う（ステップＳ２１）。
Ｌは増分内でのスピンが一定であると仮定した時の物質点の回転を示す直交テンソルである。
次に、応力ひずみ解析部１１３は、以下の式（Ａ２０）に示す材料損傷モデルである強連成方程式を作成する（ステップＳ２２）。具体的には、損傷モデルとして連続体損傷力学に基づく損傷モデルと、加工硬化モデルとして分割背応力型Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｃｋモデルとの強連成方程式（以下、「材料損傷モデル」、「材料損傷方程式」という。）を作成する。材料損傷モデルの導出方法については、＜材料損傷モデル＞の項で説明する。
次に、応力ひずみ解析部１１３は、作成した材料損傷方程式をＮｅｗｔｏｎ−Ｒａｐｈｓｏｎ法により、求解し、塑性常数Δγを求める（ステップＳ２３、Ｓ２４）。応力ひずみ解析部１１３は、求めた塑性常数Δγを、後述する式（１４）に代入して、ひずみに応じた物理量を算出する（ステップＳ２５）。
次に、応力ひずみ解析部１１３は、算出した物理量のうちの損傷値Ｄｎ＋１を基に、破壊判定を行う（ステップＳ２６）。損傷値Ｄｎ＋１が閾値ＤＦ（材料定数）よりも大きい場合、応力ひずみ解析部１１３は、その要素は破壊されたと判断して、除去する。つまり、応力ひずみ解析部１１３によって破壊されたと判断された要素は、以後の演算の対象から外されることになり、解析対象の有限要素モデルを描画する際にも表示されないことになる。尚、実施形態では、破壊されたと判断された要素は、除去することとしているが、その要素の剛性を極端に小さくし、実質的に除去されているものと近似して計算することとしてもよい。
次に、応力ひずみ解析部１１３は、Ｃｏｎｓｉｓｔｅｎｔ接線剛性を算出する（ステップＳ２７）。尚、Ｃｏｎｓｉｓｔｅｎｔ接線剛性を算出する方法は、＜Ｃｏｎｓｉｓｔｅｎｔ接線剛性＞の項で説明する。
ここで、図４に、材料に一定の負荷をかけ続けた場合の応力ひずみ線図を示す。縦軸は、応力を示し、横軸は、ひずみを示す。曲線のグラフは、実測値を示す。横軸の実線矢印は、シミュレーションにおける変形量（ひずみ増分）を示す。この変形量は、剛性方程式の解から算出された各要素のひずみの、合計である全ひずみ増分である。
シミュレーションにおける応力ひずみ線図が実測値と近似するためには、応力ひずみ線図の傾き、つまり、接線剛性が、より実測に近く求められる必要がある。シミュレーションでは、１ステップごとに加える負荷を増加させて計算を行うので、このステップ間（２点間）においては線形と考えて、Ｃｏｎｓｉｓｔｅｎｔ接線剛性を求めるものとする。より実測に近いＣｏｎｓｉｓｔｅｎｔ接線剛性を求めるためには、図４の応力ひずみ線図で示すように、より実測に近い積分応力を求める必要がある。実施形態の材料損傷モデルは、この応力を、より実測に近く求めることができる。尚、図４の横軸の実線矢印で示される全ひずみ増分Δεは、材料損傷モデルの式（Ａ２０）におけるＭとＹｎ＋１に含まれており、材料損傷モデルにより全ひずみ増分Δεに対応する積分応力を求めることができる。Ｍは、後述の式（Ａ１８）のようにｑＥを含み、ｑＥは、式（Ａ９）のように、ηＥを含み、ηＥは、式（Ａ７）のようにσｔｒｙ ｎ＋１を含み、σｔｒｙ ｎ＋１は、式（Ａ３）のようにεｅｔｒｙ ｎ＋１を含み、εｅｔｒｙ ｎ＋１は、式（１５）のようにΔεを含む。また、Ｙｎ＋１は、式（Ａ１２）のようにεｅｔｒｙ ｎ＋１を含み、εｅｔｒｙ ｎ＋１は、式（１５）のようにΔεを含む。
Ｃｏｎｓｉｓｔｅｎｔ接線剛性Ｄｅｐ（チルド）は、以下の式（Ｇ２）で示すように、剛性マトリクス[Ｋ]を求めるために用いられる。尚、「英文字記号（チルド）」は、式中で英文字記号（変数）に「〜」が付されているものを示す。
[Ｂ]は、ひずみと変位の関係を示すＢマトリクス、ｖは要素の体積を示す。
応力ひずみ解析部１１３は、求めたＣｏｎｓｉｓｔｅｎｔ接線剛性Ｄｅｐ（チルド）を用いて、剛性マトリックス[Ｋ]を作成する（ステップＳ２８）。
応力ひずみ解析部１１３は、有限要素モデルの内力から外力を引いた値（残差力）を算出し、残差力がゼロ以下である場合（内力と外力とが釣り合った場合）は（ステップＳ２９：ＯＫ）、計算が収束したと判断して、処理を終了する。応力ひずみ解析部１１３は、この時の物理量を、現在のステップの物理量として、制御部１１１に出力し、処理を終了する。
一方、残差力がゼロより大きい場合は（ステップＳ２９：ＮＧ）、応力ひずみ解析部１１３は、ステップＳ２０からの処理を繰り返す。尚、ステップＳ２１〜Ｓ２７の処理は、有限要素モデルを構成する全要素について行われる。ステップＳ２０、Ｓ２８の処理は、ステップＳ２１〜Ｓ２７で行われた個々の要素に対する演算結果、つまり、全要素に対する処理結果を元にした、解析対象の有限要素モデル（例えば、部品）の形状全体についての振る舞いとしての処理である。
＜材料損傷モデル＞
次に、以下の式（Ａ２０）で示される方程式ｆ（Δγ）、実施形態の材料損傷モデルについて説明する。実施形態の材料損傷モデルは、分割背応力形Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｋ則による加工硬化モデルと損傷モデルとを強連成した静的、完全陰解法の計算手法である材料損傷方程式である。
Ｇはせん断弾性係数、Ｍは後述の式（Ａ１８）第２式、Ｄｎは損傷値、Ｙｎ＋１は損傷エネルギ解放率、ｒ、Ｓは損傷発展速度を決める定数である。
次に、以下の式（Ａ２０）で示される方程式ｆ（Δγ）、実施形態の材料損傷モデルについて説明する。実施形態の材料損傷モデルは、分割背応力形Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｋ則による加工硬化モデルと損傷モデルとを強連成した静的、完全陰解法の計算手法である材料損傷方程式である。
応力ひずみ解析部１１３は、材料損傷方程式を作成し、作成した方程式をＮｅｗｔｏｎ−Ｒａｐｈｓｏｎ法で解いてΔγを求める。そして、塑性常数Δγを、後述する式（１４）に代入することで、５つの物理量、応力σｎ＋１、等方硬化内部変数Ｒｎ＋１、全体背応力βｎ＋１、分割背応力βｋｎ＋１、損傷値Ｄｎ＋１を求める。
まず、分割背応力型Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｃｋモデルについて説明する。分割背応力形Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｃｋ則は、Ｃｈａｂｏｃｈｅらによって提案されたモデルであり、以下の式（１）で示される。このモデルは、Ｎ組の分割背応力βｋを重ね合わせる事で、全体背応力を表現するものである。以下、表記に関してスカラ(０階テンソル)量は細字，２階テンソル量は太字，４階テンソル量は太字にチルド記号を付けて示す。尚、以下、式中で英文字記号（変数）に「・」又は「〜」を付すものは、「英文字記号（チルド）」と記載するものとする。
応力空間における降伏曲面の移動量は、この全体背応力を用いて表現される。また各分割背応力は、以下の式（２）のように、非線形の背応力挙動を表現するＡｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｃｋ則によって表現される。
ここでａｋ、ｂｋは、各分割背応力に対する材料定数、εｐ（チルド）は塑性ひずみ、γ（チルド）は塑性乗数を示す。このように分割背応力形Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｃｋ則は全体背応力を複数の非線形特性の組み合わせにより表現する事ができるため、単一の背応力を用いた場合のＡｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｃｋ則よりも比較的複雑な背応力発展挙動を表現する事が可能となる。
また、等方硬化則として、以下の式（３）のように、Ｍ個の指数関数則を組み合わせた形式を用いる。δｙは等方硬化による相対相当応力、δｙ０は初期降伏応力、Ａｋ、Ｂｋは材料定数、Ｒは等方硬化の内部変数を示す。
分割背応力形Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｃｋ則と連続体損傷力学に基づく材料損傷を考慮した散逸ポテンシャルΨは、以下の式（４）で示される。
ここでＦｐは降伏関数、Ｆβは背応力項、ＦＤは損傷項を示す。δはＣａｕｃｈｙ応力、Ｄは材料損傷値、ｒ、Ｓは損傷発展速度を決める定数である。
式（４）中のＹは損傷発展による内部エネルギ密度の変化量を示すもので、損傷エネルギ解放率と呼ばれ、以下の式（５）のように応力三軸度関数Ｒνと関与する量である。Ｅはヤング率、νはポアソン比、ｐは静水圧、σｅｑは相当応力(Ｍｉｓｅｓ応力)を示す。
更に、弾性材料剛性テンソルをＤｅ（チルド）、偏差射影テンソルをＩｄ（チルド）、対称射影テンソルをＩｓ（チルド）、４階の恒等テンソルをＩ（チルド）、２階の等方テンソルをＩとおくと、それぞれ以下の式（６）で表記される。ただしＩの成分は、クロネッカーのデルタδと等しくなる。λ、μはＬａｍｅの定数である。
式（４）のように散逸ポテンシャルに降伏関数以外の項が含まれるため、以下の式（７）のような非関連流れ則が導かれる。
ここで有効相対応力ηを、以下の式（８）に示すようにおくと、
式（７）中の各流れベクトルは、以下の式（９）となる。
ただし、式（７）第３式には、以下の式（１０）で示されるように、分割背応力βｋに関する自由エネルギΨβに関して熱力学的応力すなわち分割背応力βｋとその内部変数Ｘｋとの関係式を用いた。
連続体損傷力学では、マイクロクラックやマイクロボイドなどの微視的な材料損傷によって内力を受け持つ有効面積の減少が、材料の巨視的な弾性材料剛性の低下を引き起こすものと考えており，構成則を式（１１）のように表現する。
式（７）に示される分割背応力の発展則をそのまま使用して応力積分を行うと、全ての分割背応力βｋが未知数となり応力積分に多大の時間を要してしまう。これを避けるためにステップｎからステップｎ＋１に掛けて塑性ひずみ増分の方向が一定であると仮定して、式（７）第３式に対してこの時間増分間での解析解を求めると、最終的に式（１２）および式（１３）が得られる。ただし添字ｎ、ｎ＋１は、そのステップにおける量を示す。式（１２）および式（３）のように時間増分間における解析解を用いる事で、各分割背応力はステップｎの既知量となり全体背応力βｎ＋１のみが背応力の未知数となる。
ただし、
である。式（７）の流れ則を後退形Ｅｕｌｅｒ表記したものと、式（１２）および弾塑性状態における降伏関数の適合条件から、式（１４）で示すＲｅｔｕｒｎ−Ｍａｐｐｉｎｇ方程式が非線形連立方程式として得られる。
ただし、εｅｔｒｙ ｎ＋１は試行弾性ひずみと呼ばれ、式（１５）のように前ステップの弾性に全ひずみ増分を加えたもので定義する。
よって、式（１４）は、応力σｎ＋１、等方硬化内部変数Ｒｎ＋１、全体背応力βｎ＋１、損傷値Ｄｎ＋１、塑性常数Δγの合計１５個の独立な未知変数成分を含む方程式となる。
さらに、式（１４）第１式より、
となる。
また有効相対応力ηｎ＋１は、式（８）に式（Ａ１）、式（１２）を用いる事で、
となる。ただし試行応力を以下の式（Ａ３）のように定義した。
さらに降伏曲面に対する法線ベクトルは式（１３）より、
となる。ただし相対相当応力を、
と定義した。これらを式（Ａ２）に代入して整理すると、
が得られる。ここで、
とする。さらに式（Ａ６）の両辺に√３／２をかけて、両辺のノルムを取ると、
となる。さらに式（１４）第５式より、
とおくと、式（Ａ８）は、
となる。式（Ａ１０）と、損傷発展式である式（１４）第４式とを連立させる事で、Ｒｅｔｕｒｎ−Ｍａｐｐｉｎｇ方程式は、
となる。この段階で、Ｒｅｔｕｒｎ−Ｍａｐｐｉｎｇ塑性方程式は、塑性乗数Δγと損傷値Ｄｎ＋１の２個のスカラ未知数にまで縮約させた。
次に、まず損傷エネルギ解放率Ｙについて考える。式（５）より、
となる。式（Ａ１２）の右辺第１項は、式（８）および式（１２）を用いて、
と表記できる。さらに式（Ａ６）より、
となる。これを式（Ａ１３）に代入すると、
となり。式（Ａ１５）右辺をｇと置く。以上より、損傷エネルギ解放率Ｙを、式（Ａ１６）のように書き直すことができる。
次に、式（Ａ１１）第1式より、
ここで、
とおくと、式（Ａ１７）から、
となる。式（Ａ１８）と（Ａ１９）とを、式（Ａ１１）第２式に代入して整理すると、式（Ａ２０）に示すような、塑性乗数Δγのみのスカラ１変数に縮約されたＲｅｔｕｒｎ−Ｍａｐｐｉｎｇ方程式が得られる。
この式（Ａ２０）で示す方程式ｆ（Δγ）が、実施形態の材料損傷方程式（モデル）である。
また有効相対応力ηｎ＋１は、式（８）に式（Ａ１）、式（１２）を用いる事で、
さらに降伏曲面に対する法線ベクトルは式（１３）より、
となる。ただし相対相当応力を、
と定義した。これらを式（Ａ２）に代入して整理すると、
が得られる。ここで、
とする。さらに式（Ａ６）の両辺に√３／２をかけて、両辺のノルムを取ると、
となる。さらに式（１４）第５式より、
とおくと、式（Ａ８）は、
となる。式（Ａ１０）と、損傷発展式である式（１４）第４式とを連立させる事で、Ｒｅｔｕｒｎ−Ｍａｐｐｉｎｇ方程式は、
となる。この段階で、Ｒｅｔｕｒｎ−Ｍａｐｐｉｎｇ塑性方程式は、塑性乗数Δγと損傷値Ｄｎ＋１の２個のスカラ未知数にまで縮約させた。
次に、まず損傷エネルギ解放率Ｙについて考える。式（５）より、
となる。式（Ａ１２）の右辺第１項は、式（８）および式（１２）を用いて、
と表記できる。さらに式（Ａ６）より、
となる。これを式（Ａ１３）に代入すると、
となり。式（Ａ１５）右辺をｇと置く。以上より、損傷エネルギ解放率Ｙを、式（Ａ１６）のように書き直すことができる。
次に、式（Ａ１１）第1式より、
ここで、
とおくと、式（Ａ１７）から、
となる。式（Ａ１８）と（Ａ１９）とを、式（Ａ１１）第２式に代入して整理すると、式（Ａ２０）に示すような、塑性乗数Δγのみのスカラ１変数に縮約されたＲｅｔｕｒｎ−Ｍａｐｐｉｎｇ方程式が得られる。
この式（Ａ２０）で示す方程式ｆ（Δγ）が、実施形態の材料損傷方程式（モデル）である。
実施形態では、さらに、応力ひずみ解析部１１３は、Ｕｐｄａｔｅｄ−Ｌａｇｒａｎｇｅ形の有限変形モデルに対応できるよう、応力積分計算（式（Ａ２０）の求解）の前にＣａｕｃｈｙ応力に対して、式（１６）で示される共回転変換を行う。
ただし、Ｌは増分内でのスピンが一定であると仮定した時の物質点の回転を示す直交テンソルである。実施形態の材料損傷方程式（材料損傷モデル）を有限変形モデルに対応させることで、繰り返し大負荷を受ける極低サイクル疲労破壊の中でも、変形が大きなものに対応することが可能になる。このような変形が大きな破壊とは、例えば、「針金を繰り返し何度も大きく曲げると数回目で折れる」、「ジュースの空き缶のフタを何度も繰り返し大きく曲げていると数回目で折れる」といった事例である。
＜Ｃｏｎｓｉｓｔｅｎｔ接線剛性＞
応力ひずみ解析部１１３が、Ｃｏｎｓｉｓｔｅｎｔ接線剛性を算出する方法について説明する。Ｃｏｎｓｉｓｔｅｎｔ接線剛性は、完全陰解法によるものである。
応力ひずみ解析部１１３が、Ｃｏｎｓｉｓｔｅｎｔ接線剛性を算出する方法について説明する。Ｃｏｎｓｉｓｔｅｎｔ接線剛性は、完全陰解法によるものである。
式（Ａ６）より、ηｎ＋１とηＥとは互いにスカラ倍されている関係に過ぎないので、両者の方向は一致する。よって、
となる。式（Ｂ１）および式（Ａ１８）を用いて、式（Ａ１）で示される構成側を以下のように変形する。
これよりσｎ＋１は内部変数を用いて、
と表記できる。以上より、Ｃｏｎｓｉｓｔｅｎｔ接線剛性Ｄｅｐ（チルド）は、σｎ＋１の各内部変数を考慮し、以下の式（Ｂ４）に示すように、
と表記できる。各微分について導出すると以下の式（Ｂ５）が得られる。
また式（Ｂ４）の塑性乗数の試行弾性ひずみによる微分は、式（Ａ２０）で示した縮約されたＲｅｔｕｒｎ−Ｍａｐｐｉｎｇ方程式を、試行弾性ひずみで微分する事により以下の式（Ｂ６）のように得られる。
右辺の分母の微分は、式（Ａ２２）から既に求まっている。そこで分子の微分について考えると、以下の式（Ｂ７）に示すように、
となる。各微分について計算すると、以下の式（Ｂ８）を得る。
以上より、陽的にＣｏｎｓｉｓｔｅｎｔ接線剛性を計算する事が可能となる。
＜シミュレーションの例＞
ここでは、図５に示す有限要素モデルを用い、上面に対してＺ方向の両振りの繰返し強制変位を与えた場合のシミュレーション結果を説明する。尚、有限要素モデルの上面以外の他の面については、単軸負荷の対称条件を満たす様に拘束するものとする。図５（ａ）は、有限要素モデルを示す図であり、図５（ｂ）は、シミュレーションで用いた材料定数の表を示す図である。シミュレーションで用いたヤング率、ポアソン比、初期降伏応力、等方硬化定数、移動硬化定数、損傷定数は、図５（ｂ）の表に示すとおりである。
となる。式（Ｂ１）および式（Ａ１８）を用いて、式（Ａ１）で示される構成側を以下のように変形する。
これよりσｎ＋１は内部変数を用いて、
と表記できる。以上より、Ｃｏｎｓｉｓｔｅｎｔ接線剛性Ｄｅｐ（チルド）は、σｎ＋１の各内部変数を考慮し、以下の式（Ｂ４）に示すように、
と表記できる。各微分について導出すると以下の式（Ｂ５）が得られる。
また式（Ｂ４）の塑性乗数の試行弾性ひずみによる微分は、式（Ａ２０）で示した縮約されたＲｅｔｕｒｎ−Ｍａｐｐｉｎｇ方程式を、試行弾性ひずみで微分する事により以下の式（Ｂ６）のように得られる。
右辺の分母の微分は、式（Ａ２２）から既に求まっている。そこで分子の微分について考えると、以下の式（Ｂ７）に示すように、
となる。各微分について計算すると、以下の式（Ｂ８）を得る。
以上より、陽的にＣｏｎｓｉｓｔｅｎｔ接線剛性を計算する事が可能となる。
＜シミュレーションの例＞
ここでは、図５に示す有限要素モデルを用い、上面に対してＺ方向の両振りの繰返し強制変位を与えた場合のシミュレーション結果を説明する。尚、有限要素モデルの上面以外の他の面については、単軸負荷の対称条件を満たす様に拘束するものとする。図５（ａ）は、有限要素モデルを示す図であり、図５（ｂ）は、シミュレーションで用いた材料定数の表を示す図である。シミュレーションで用いたヤング率、ポアソン比、初期降伏応力、等方硬化定数、移動硬化定数、損傷定数は、図５（ｂ）の表に示すとおりである。
図６は、実測値の応力ひずみ線図であり、ひずみ振幅一定で繰り返し負荷を与えた時の、応力、ひずみ挙動を示す。詳細には、ひずみは、−０．０４から＋０．０４まで変化させられ、「Ｓ５００ｓ，Δε＝８％，４５ｄａｙｓ」が示すように、入力負荷のひずみがΔ８％であり、モデルの材料はＳ５００ｓ鋼であって、４５日間腐食させた材料である。図６の応力ひずみ線図は、Ｂａｕｓｃｈｉｎｇｅｒ効果が表れたヒステリシスループとなっている。図７は、非損傷モデルによるシミュレーションの結果を示すヒステリシスループであり、図８は、実施形態の材料損傷モデルによるシミュレーションの結果を示すヒステリシスループである。
比較用の実測データの選定として、一般的な鋼材やアルミ材でも同様の傾向は見られるため差し支えないが、図６に示したように損傷に伴う変形挙動（矢印で示すような、サイクルの増加につれて生じるズレ）がより顕著に見られるため事前に腐食させたＳ５００ｓ鋼のデータを参考にした。ヤング率が一般的な鋼に対して低めとなっている理由は腐食によるものと考えられる。
また、ここで参考にした実材料の詳細な損傷定数が不明なため各材料定数として、図５（ｂ）の表に示すような現実的な金属材料の範囲で便宜的な値を入力し、定性的な変形挙動を表現できるかのみ検証した。
図７は、図５の有限要素モデルに公称ひずみ振幅４%の両振負荷を１０ｃｙｃｌｅ掛けた時のヒステリシスループであり、非損傷モデルを用いた場合の図である。図８は、実施形態の材料損傷モデルを用いた場合の図である。ヒステリシスループにおいて、横軸は強制変位量を元の長さで割った公称ひずみを示し，縦軸は引張面の荷重合計値を変形前の断面積で割った公称応力を示す。実測結果では図６のようにサイクル数の増加に伴う明瞭な加工軟化挙動が見られるが、非損傷モデルでは図７のようにヒステリシスループが塑性シェイクダウンを起こし、実測のような加工軟化挙動が見られない。一方、実施形態の材料損傷モデルでは図８のサイクル数の増加に伴う加工軟化挙動が確認できる（矢印参照）。
図９は、サイクル数に伴う各ループの引張時のピーク応力を示した図であり、縦軸は応力を示し、横軸はサイクル数を示す。白丸が非損傷モデルでの結果を示し、灰色の丸が実施形態の材料損傷モデルでの結果を示す。この図からも、実施形態の材料損傷モデルによる結果は、図６の実測と同様にサイクル数に伴う応力低下が見られ、実施形態の材料損傷モデルは、加工軟化挙動を定性的に表現できている事がわかる。
また、図１０は、実施形態の材料損傷モデルによるサイクル数に伴う損傷値の変化を示す図であり、縦軸が損傷値を示し、横軸がサイクル数を示す。図１０が示すように、実施形態の材料損傷モデルでは、サイクル数の増加に伴い損傷値が増加してゆく傾向が見られる。すなわち、実施形態の材料損傷モデルにおいて加工軟化傾向が表現されるのは、この損傷値の増加に対して、上記式（３）のように降伏関数の収縮による応力積載能力の低下が要因であると考えられる。
さらに、実測結果では図６のようにサイクル数の増加に伴い弾性領域のヤング率（剛性）が低下している様子が確認できる。図１１は、サイクル数の増加に伴うヤング率（剛性）の変化を示した図であり、縦軸がヤング率（剛性）を示し、横軸がサイクル数を示す。白丸が非損傷モデルでの結果を示し、灰色の丸が実施形態の材料損傷モデルでの結果を示す。非損傷モデルでは、弾塑性変形の履歴に関わらずヤング率（剛性）は一定値となっているが、実施形態の材料損傷モデルでは、実測結果と同様にサイクル数の増加に対してヤング率（剛性）が低下している。これは図１０のようにサイクル数の増加に伴って大きくなる損傷値に対して、式（９）で示されるように巨視的弾性材料剛性の低下が表現されるためである。
以上のように、実施形態の材料損傷モデルを用いる事で、繰返し大負荷を掛けた時の応力、ひずみ挙動を定性的に良く表現できている事がわかる。
次に、初期き裂進展挙動についてのシミュレーションを説明する。図１２に示すように、両端単純支持ばりを表す有限要素モデルに対し、矢印で示される集中負荷をかける場合をシミュレートする。有限要素モデルである梁は、Ａ５３３Ｂ鋼からなる高さＷ（＝２０ｍｍ）、幅Ｔ（＝４ｍｍ）、長さＬ（＝８０ｍｍ）の直方体であり、シミュレーションの計算は、梁の幅方向２ｍｍ分の２分の１モデルで行う。また、長さ方向の一端からＬ０（＝１２ｍｍ）の位置に、初期き裂が生じているものとする。この状態から、梁の長さ方向中央に５ｍｍの曲げ強制変位を生じさせる荷重をかけ続けた場合の結果を、図１３に示す。尚、幅方向の手前側の端面を板厚表面といい、この板厚表面から幅方向２ｍｍの位置における前記板厚表面と平行な面を板厚中央面というものとする。
次に、初期き裂進展挙動についてのシミュレーションを説明する。図１２に示すように、両端単純支持ばりを表す有限要素モデルに対し、矢印で示されるように集中荷重をかける場合をシミュレートする。
図１３（ａ）は、板厚表面のき裂進展の様子を示す図であり、画像１３１は、き裂が進展した実際の材料の写真（実験結果）であり、実線矢印は、き裂進展方向を示す。また、画像１３２は、実施形態の材料損傷モデルを用いたシミュレーション結果を示す画像である。また、図１３（ｂ）は、板厚中央面のき裂進展の様子を示す図であり、画像１３３は、き裂が進展した実際の材料の切断面の写真（実験結果）であり、破線は、き裂進展前の予き裂面を示し、実線矢印は、き裂進展方向を示す。また、画像１３４は、実施形態の材料損傷モデルを用いたシミュレーション結果を示す画像である。シミュレーション結果の画像１３２、１３４において要素（三角形）が記載されていない部分は、き裂部分を示す。
板厚表面では、図３（ａ）の実験結果の画像１３１に示すように、二手に分かれてき裂進展し、シミュレーション結果の画像１３２でも、実験結果の画像１３１と同様に、二手に分かれてき裂が進展していることが分かる。また、板厚中央面では、図３（ｂ）の実験結果の画像１３３に示すように、予き裂面から紙面右側に偏ってき裂が進展し、シミュレーション結果の画像１３４でも、実験結果の画像１３３と同様に、紙面右手に偏ってき裂が進展していることが分かる。
尚、図１３の画像１３１、１３３で示す破壊は、き裂先端近傍が大きく変形している破壊である。上述した式１６を用いて実施形態の材料損傷方程式（材料損傷モデル）を有限変形モデルに対応させたことにより、画像１３１、１３３で示すような、繰り返し大負荷を受ける極低サイクル疲労破壊の中でも変形が大きなものにも対応することが可能となり、画像１３２、１３４に示すようなシミュレーション結果を得ることが可能となっている。
実施形態の材料損傷モデルでは、損傷値が算出されるので、材料定数で示される所定の閾値を上回った場合に、き裂が生じたと判断することが可能となる。従って、損傷値が閾値を超えた要素を描画しないことで、図１３の画像１３２、画像１３４に示すようなき裂が進展した画像を表示することが可能となる。
１ 有限要素解析装置
１１ 演算処理部
１２ 入力部
１３ 出力部
１４ 記憶部
１１１ 制御部
１１２ 初期設定部
１１３ 応力ひずみ算出部
１４１ 有限要素モデル記憶部
１４２ 解析条件記憶部
１４３ 物理量記憶部
１１ 演算処理部
１２ 入力部
１３ 出力部
１４ 記憶部
１１１ 制御部
１１２ 初期設定部
１１３ 応力ひずみ算出部
１４１ 有限要素モデル記憶部
１４２ 解析条件記憶部
１４３ 物理量記憶部
Claims (9)
- 非線形有限要素解析を行うための有限要素解析装置であって、
解析対象の有限要素モデルと解析条件とを取得する取得手段と、
前記有限要素モデルの剛性行列を作成し、作成した剛性行列と前記解析条件とで規定される剛性方程式の解に基づいて、前記有限要素モデルの加工硬化特性を算出する特性計算を行う応力ひずみ算出手段と、
前記有限要素モデルの内力と外力とが釣り合うまで、前記応力ひずみ算出手段に前記特性計算を繰り返し行わせる制御手段とを備え、
前記応力ひずみ算出手段は、前回算出した前記剛性方程式の解に基づく試行弾性歪みを用い、繰り返し負荷を受けている部分であって、破壊が発生していない部分の応力とひずみを算出する加工硬化モデルと、破壊が発生している部分の破壊進行度と応力を算出する材料損傷モデルとにより応力とひずみとを連成させた材料損傷方程式を作成し、作成した材料損傷方程式の解を求め、求めた解に基づいてコンシステント接線剛性を求め、求めたコンシステント接線剛性から前記剛性行列を作成する
ことを特徴とする有限要素解析装置。 - 前記応力ひずみ算出手段は、前記材料損傷方程式の解に基づいて、前記有限要素モデルの各要素の損傷値を求め、前記損傷値が所定の閾値を上回る要素は破壊されたと判断する
ことを特徴とする請求項１に記載の有限要素解析装置。 - 前記加工硬化モデルは、分割背応力形Ａｒｍｓｔｒｏｎｇ−Ｆｒｅｄｅｒｉｋ則による加工硬化モデルである
ことを特徴とする請求項１に記載の有限要素解析装置。 - 前記コンシステント接線剛性は、完全陰解法により求める
ことを特徴とする請求項１ないし３のいずれか一項に記載の有限要素解析装置。 - 前記制御手段は、前記有限要素モデルの内力と外力とが釣り合うと、前記解析条件のうちの強制変位又は外力を、所定の増分値だけ増加させて、前記応力ひずみ算出手段に前記特性計算を繰り返し行わせる
ことを特徴とする請求項１ないし４のいずれか一項に記載の有限要素解析装置。 - 前記有限要素モデルを表示する表示部を、更に備え、
前記制御手段は、前記応力ひずみ算出手段が破壊されたと判断した要素を除いて、各要素をそれぞれの応力に応じた表示態様で、前記有限要素モデルを前記表示部に表示させる
ことを特徴とする請求項１ないし５のいずれか一項に記載の有限要素解析装置。 - 前記剛性算出手段は、作成した材料損傷方程式の解を求める前に、前記有限要素モデルの各要素の応力に対して、共回転変換を行う
ことを特徴とする請求項１ないし６のいずれか一項に記載の有限要素解析装置。 - 非線形有限要素解析を行うための有限要素解析装置で用いられる有限要素解析方法であって、
解析対象の有限要素モデルと解析条件とを取得する取得ステップと、
前記有限要素モデルの剛性行列を作成し、作成した剛性行列と前記解析条件とで規定される剛性方程式の解に基づいて、前記有限要素モデルの加工硬化特性を算出する特性計算を行う応力ひずみ算出ステップと、
前記有限要素モデルの内力と外力とが釣り合うまで、前記応力ひずみ算出ステップに前記特性計算を繰り返し行わせる制御ステップとを備え、
前記応力ひずみ算出ステップは、前回算出した前記剛性方程式の解に基づく試行弾性歪みを用い、繰り返し負荷を受けている部分であって、破壊が発生していない部分の応力とひずみを算出する加工硬化モデルと、破壊が発生している部分の破壊進行度と応力を算出する材料損傷モデルとにより応力とひずみとを連成させた材料損傷方程式を作成し、作成した材料損傷方程式の解を求め、求めた解に基づいてコンシステント接線剛性を求め、求めたコンシステント接線剛性から前記剛性行列を作成する
ことを特徴とする有限要素解析方法。 - 非線形有限要素解析を行うための有限要素解析装置で用いられる有限要素解析プログラムであって、
解析対象の有限要素モデルと解析条件とを取得する取得手段と、
前記有限要素モデルの剛性行列を作成し、作成した剛性行列と前記解析条件とで規定される剛性方程式の解に基づいて、前記有限要素モデルの加工硬化特性を算出する特性計算を行う応力ひずみ算出手段と、
前記有限要素モデルの内力と外力とが釣り合うまで、前記応力ひずみ算出手段に前記特性計算を繰り返し行わせる制御手段として、コンピュータを機能させ、
前記応力ひずみ算出手段は、前回算出した前記剛性方程式の解に基づく試行弾性歪みを用い、繰り返し負荷を受けている部分であって、破壊が発生していない部分の応力とひずみを算出する加工硬化モデルと、破壊が発生している部分の破壊進行度と応力を算出する材料損傷モデルとにより応力とひずみとを連成させた材料損傷方程式を作成し、作成した材料損傷方程式の解を求め、求めた解に基づいてコンシステント接線剛性を求め、求めたコンシステント接線剛性から前記剛性行列を作成する
ことを特徴とする有限要素解析プログラム。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2014126490A JP6135607B2 (ja) | 2014-06-19 | 2014-06-19 | 有限要素解析装置、該方法及び該プログラム |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2014126490A JP6135607B2 (ja) | 2014-06-19 | 2014-06-19 | 有限要素解析装置、該方法及び該プログラム |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JP2016004543A true JP2016004543A (ja) | 2016-01-12 |
| JP6135607B2 JP6135607B2 (ja) | 2017-05-31 |
Family
ID=55223742
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2014126490A Active JP6135607B2 (ja) | 2014-06-19 | 2014-06-19 | 有限要素解析装置、該方法及び該プログラム |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JP6135607B2 (ja) |
Cited By (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108089457A (zh) * | 2017-11-29 | 2018-05-29 | 北京航空航天大学 | 一种基于在线有限元仿真的过程质量控制方法 |
| CN108363883A (zh) * | 2018-03-07 | 2018-08-03 | 吉林大学 | 一种结构减重孔板的刚度性能计算方法和强度计算方法 |
| CN108647446A (zh) * | 2018-05-12 | 2018-10-12 | 北京交通大学 | 一种客货共线双块式无砟轨道伤损分析的方法及系统 |
| CN108920861A (zh) * | 2018-07-17 | 2018-11-30 | 暨南大学 | 一种含裂纹梁结构单元损伤因子的等效方法 |
| JP2020189310A (ja) * | 2019-05-22 | 2020-11-26 | 本田技研工業株式会社 | 加工管理システム、加工管理方法、及びプログラム |
| CN112100762A (zh) * | 2020-08-27 | 2020-12-18 | 南京航空航天大学 | 一种陶瓷基复合材料超单元结构数值仿真计算方法 |
| CN114357661A (zh) * | 2022-01-19 | 2022-04-15 | 厦门大学 | 柱塞-滑靴组件耦合动力学模型建立方法 |
| CN115482890A (zh) * | 2022-09-14 | 2022-12-16 | 东南大学 | 基于弹性应变能密度分解的局部梯度损伤模型构建方法 |
| CN117195667A (zh) * | 2023-11-08 | 2023-12-08 | 中交天津港湾工程研究院有限公司 | 一种dcm联合块石基床振夯施工仿真分析方法 |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108052717B (zh) * | 2017-12-01 | 2021-09-03 | 中国航空工业集团公司沈阳飞机设计研究所 | 一种基于局部应力-应变法的疲劳寿命校准方法 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2005141645A (ja) * | 2003-11-10 | 2005-06-02 | Ohbayashi Corp | 非線形有限要素解析装置及び方法、コンピュータプログラム、記録媒体 |
| JP2008142774A (ja) * | 2006-11-14 | 2008-06-26 | Jfe Steel Kk | 応力−ひずみ関係シミュレート方法、応力−ひずみ関係シミュレーションシステム、応力−ひずみ関係シミュレーションプログラム、及び当該プログラムを記録した記録媒体 |
| JP2012104042A (ja) * | 2010-11-12 | 2012-05-31 | Kobe Steel Ltd | 均一伸びの予測方法および均一伸びの予測プログラム |
| JP2014052909A (ja) * | 2012-09-07 | 2014-03-20 | Toyota Central R&D Labs Inc | 情報処理装置及びプログラム |
-
2014
- 2014-06-19 JP JP2014126490A patent/JP6135607B2/ja active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2005141645A (ja) * | 2003-11-10 | 2005-06-02 | Ohbayashi Corp | 非線形有限要素解析装置及び方法、コンピュータプログラム、記録媒体 |
| JP2008142774A (ja) * | 2006-11-14 | 2008-06-26 | Jfe Steel Kk | 応力−ひずみ関係シミュレート方法、応力−ひずみ関係シミュレーションシステム、応力−ひずみ関係シミュレーションプログラム、及び当該プログラムを記録した記録媒体 |
| JP2012104042A (ja) * | 2010-11-12 | 2012-05-31 | Kobe Steel Ltd | 均一伸びの予測方法および均一伸びの予測プログラム |
| JP2014052909A (ja) * | 2012-09-07 | 2014-03-20 | Toyota Central R&D Labs Inc | 情報処理装置及びプログラム |
Non-Patent Citations (2)
| Title |
|---|
| 大畑充 外４名: "繰返し荷重下での構造用鋼の延性き裂発生クライテリオン", 溶接学会論文集 [ONLINE], vol. 第２１巻，第４号, JPN6017009909, November 2003 (2003-11-01), pages 592 - 602, ISSN: 0003523351 * |
| 森本秀夫 外１名: "結晶集合組織を考慮した弾／結晶粘塑性有限要素法による深絞り加工の解析", 古河電工時報, vol. 第１０３号, JPN6017009908, 11 January 1999 (1999-01-11), pages 7 - 12, ISSN: 0003523350 * |
Cited By (14)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108089457A (zh) * | 2017-11-29 | 2018-05-29 | 北京航空航天大学 | 一种基于在线有限元仿真的过程质量控制方法 |
| CN108363883A (zh) * | 2018-03-07 | 2018-08-03 | 吉林大学 | 一种结构减重孔板的刚度性能计算方法和强度计算方法 |
| CN108647446A (zh) * | 2018-05-12 | 2018-10-12 | 北京交通大学 | 一种客货共线双块式无砟轨道伤损分析的方法及系统 |
| CN108647446B (zh) * | 2018-05-12 | 2024-04-26 | 北京交通大学 | 一种客货共线双块式无砟轨道伤损分析的方法及系统 |
| CN108920861A (zh) * | 2018-07-17 | 2018-11-30 | 暨南大学 | 一种含裂纹梁结构单元损伤因子的等效方法 |
| CN108920861B (zh) * | 2018-07-17 | 2023-01-20 | 暨南大学 | 一种含裂纹梁结构单元损伤因子的等效方法 |
| JP2020189310A (ja) * | 2019-05-22 | 2020-11-26 | 本田技研工業株式会社 | 加工管理システム、加工管理方法、及びプログラム |
| CN112100762B (zh) * | 2020-08-27 | 2024-03-08 | 南京航空航天大学 | 一种陶瓷基复合材料超单元结构数值仿真计算方法 |
| CN112100762A (zh) * | 2020-08-27 | 2020-12-18 | 南京航空航天大学 | 一种陶瓷基复合材料超单元结构数值仿真计算方法 |
| CN114357661A (zh) * | 2022-01-19 | 2022-04-15 | 厦门大学 | 柱塞-滑靴组件耦合动力学模型建立方法 |
| CN115482890B (zh) * | 2022-09-14 | 2023-10-31 | 东南大学 | 基于弹性应变能密度分解的局部梯度损伤模型构建方法 |
| CN115482890A (zh) * | 2022-09-14 | 2022-12-16 | 东南大学 | 基于弹性应变能密度分解的局部梯度损伤模型构建方法 |
| CN117195667A (zh) * | 2023-11-08 | 2023-12-08 | 中交天津港湾工程研究院有限公司 | 一种dcm联合块石基床振夯施工仿真分析方法 |
| CN117195667B (zh) * | 2023-11-08 | 2024-03-22 | 中交天津港湾工程研究院有限公司 | 一种dcm联合块石基床振夯施工仿真分析方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| JP6135607B2 (ja) | 2017-05-31 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| JP6135607B2 (ja) | 有限要素解析装置、該方法及び該プログラム | |
| Zhou et al. | Phase field modeling of quasi-static and dynamic crack propagation: COMSOL implementation and case studies | |
| De Moura et al. | Cohesive zone model for high-cycle fatigue of adhesively bonded joints under mode I loading | |
| WO2014141794A1 (ja) | 応力－ひずみ関係シミュレート方法、スプリングバック量予測方法およびスプリングバック解析装置 | |
| Nielsen et al. | Cohesive traction–separation laws for tearing of ductile metal plates | |
| Ponthot et al. | A cascade optimization methodology for automatic parameter identification and shape/process optimization in metal forming simulation | |
| JP3708067B2 (ja) | 弾塑性体の亀裂進展予測方法および変形予測方法 | |
| US20150310143A1 (en) | Deformation analysis device, deformation analysis method, and program | |
| JP2008142774A (ja) | 応力−ひずみ関係シミュレート方法、応力−ひずみ関係シミュレーションシステム、応力−ひずみ関係シミュレーションプログラム、及び当該プログラムを記録した記録媒体 | |
| JP2003270060A (ja) | 応力ひずみ解析システム及びそれに用いる応力ひずみ解析方法並びにそのプログラム | |
| Malcher et al. | Experimental program and numerical assessment for determination of stress triaxiality and J3 effects on AA6101-T4 | |
| Halama et al. | Implementation of MAKOC cyclic plasticity model with memory | |
| Sellier et al. | Orthotropic damage coupled with localized crack reclosure processing: Part II: Applications | |
| Firat | A numerical analysis of combined bending–torsion fatigue of SAE notched shaft | |
| JP6776876B2 (ja) | 高分子材料のシミュレーション方法 | |
| Jin et al. | XFEM to couple nonlocal micromechanics damage with discrete mode I cohesive fracture | |
| Pham et al. | Application of a graphical method on estimating forming limit curve of automotive sheet metals | |
| Tagiltsev et al. | Combined experimental/theoretical approach to residual stresses within multiplicative elasto-plasticity | |
| Kohansal-Vajargah et al. | Quadratic tetrahedral micropolar element for the vibration analysis of three-dimensional micro-structures | |
| Barbu et al. | Validation on large scale tests of a new hardening–softening law for the Barcelona plastic damage model | |
| Han et al. | Incorporation of gradient-enhanced microplane damage model into isogeometric analysis | |
| Bhattacharyya | A model reduction approach in space and time for fatigue damage simulation | |
| Hofman et al. | A numerical framework for simulating progressive failure in composite laminates under high-cycle fatigue loading | |
| Yuan et al. | A configurational force based anisotropic damage model for original isotropic materials | |
| Zhang et al. | A novel damage model integrated into the elastoplastic constitutive model and numerical simulations of reinforced concrete structures under cyclic loading |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20160323 |
|
| A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20170302 |
|
| TRDD | Decision of grant or rejection written | ||
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20170328 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20170410 |
|
| R150 | Certificate of patent or registration of utility model |
Ref document number: 6135607 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R150 |