« 2012年6月 | トップページ | 2012年8月 »

2012年7月

2012年7月29日 (日)

Salome-Meca 2012.2

Code_Asterの公式サイトで、7月23日より「Salome-Meca 2012.2」のダウンロードが出来るようになったようです。

バージョンは、SALOMEが6.5.0、Code_Aster安定版が10.7、テスト版が11.2とのことです。

私も早速ダウンロードし、以前構築した仮想PC上のdebian6にインストールしてみました。

インストール方法はこちらの解説通りで、具体的な手順は下記の通りです。

----------------------

(1)Windows上でダウンロードした「SALOME-MECA-2012.2-LGPL.tgz」のファイルを、LAN内に設置したNAS(バッファローのLinkstation)にコピー。

2.仮想PC上のdebian6(rootでログイン済み)からNASのフォルダを表示させ、「SALOME-MECA-2012.2-LGPL.tgz」をdebian6上の /opt/にコピー。
Linux上からLAN上のNASを見に行く方法は
こちら

なお、仮想PCのLinux上から直接このファイルをダウンロード出来るのであれば、1,2の作業は不要。

3.debian6のファイルブラウザで/opt を表示させ、「SALOME-MECA-2012.2-LGPL.tgz」を右クリックし、「展開」で、/opt の下に解凍する。

4.「アプリケーション」-「アクセサリ」-「端末」を起動し、

 cd /opt/SALOME-MECA-2012.2-LGPL

で、ディレクトリを移動。

 ./postinstall.sh

を実行。

5.rootからログアウトし、普段使う一般ユーザー(私の場合は"aki")で再ログイン。

6.ファイルブラウザで「/opt/SALOME-MECA-2012.2-LGPL」ディレクトリを表示し、「runSalomeMeca」を、「Alt」キーを押しながらデスクトップ上へドラッグ&ドロップ。「ここにリンクを作成」を選択。

7.デスクトップ上のリンクをダブルクリック、「実行」でSalome-mecaが起動する。

以上です。

------------------------------

起動は特に問題ありませんでした。SALOMEの見た目もそれほど変わっていないようですが、様々な機能が追加されているのでしょう。

例えばX-FEMがWizardで扱えるようになっているようです。私は破壊解析関係には全く疎いのですが、一度試してみたいところです。

#########################################

さて、動作検証ですが、先日実施した「節点反力を集約する方法」の例題(reac04)を実施しました。
新しいフォルダを作成して、そこに「reac04.hdf」と「reac04.comm」の2つのみコピーし、Salomeを起動して「reac04.hdf」を読み込みました。

私はJOBの実行をastkで行っているのですが、今回はSalome上で行ってみました(Salome上でJOB実行する機能自体は以前から有ったようなのですが、慣れの問題から、私は使っていませんでした)。

Sutudy caseを作成するために、まずMeshのモジュール起動する必要があるようですので、一旦「Mesh」を立ち上げてから、すぐに「Aster」に切り替えます。

メニューの「Aster」から「Add Study case」で設定画面を開き、下図の通り設定しました。

Pic02

commファイルはディスク上からの指定、Meshはオブジェクトブラウザでの指定です(今回、commファイルは前回の事例実行時のコピーのままで、無編集です)。

Aster versionはSTA10で、メモリと制限時間も増やしています。

OKで設定完了後、オブジェクトブラウザの「Aster」の下に作成された「reac04」を右クリックして、「RUN」で実行です。

結果は下図のとおりで、各ファイルがオブジェクトブラウザで「Results」の下に
表示されます。

Pic01


med形式の結果は拡張子が「rmed」で出力されますが、これはParaVIS(ParaVIEW)で読めるようになっています(Salome6.3ではrmedを読もうとするとエラーになっていましたが、6.4から読めるようになっているようです)。

また、resuファイルなどのテキストファイルも、オブジェクトブラウザから右クリックで、テキストとして表示するか、編集することが出来ます。ただし、編集するためにはメニューの「File」- 「Preferences」の「Aster」で、テキストエディタを指定しておく必要があります。私は「/usr/bin/gedit」を指定しています。

以上のように、問題なく非線形解析が出来ました。

ただし、例題や検証問題などの実行については、まだしばらくはDEXCS-Salomeを使って行こうと思っています(現在はA3-64)。

2012年7月19日 (木)

時間ステップ幅の自動調整機能

非線形解析では、荷重や変位を細かく分割して徐々に与えていきます。通常は等分割で与えていくのですが、解析の過程の一部に非線形性が強いところがあると、その部分に合わせて全体の分割を細かくすることになってしまいます。分割が粗すぎると、そこで計算が収束せず停止してしまいます。長い計算の最後のほうで収束せずに停止してしまうと、ガッカリです。

そこで役に立つのが、時間ステップ幅の自動調整機能です。

これは、初期の時間ステップ幅で計算していて、計算が規定のイテレーション回数で収束しなかったり、剛性マトリックスが異常になるなどのエラーが発生した場合、そこで計算を中止するのではなく、時間ステップ幅を細かくした上で、その時点から計算を再開してくれるというものです。

また逆に、計算の収束が良好な場合、時間ステップ幅を増やして効率を上げます。

この機能は非線形解析を得意とするFEM解析ソフトウェアにはおおむね備わっていますが、Code_Asterでも「DEFI_LIST_INST」によって実現できます。

今回、以前実施した「追従荷重(その2)」の例題をそのまま用いて、この時間ステップ幅の自動調整機能を追加してみました。

ソフトウェアは、DEXCS2011-Salome-A3(64)(SALOME6.3 + Code_Aster STA10.5)の利用です。

この問題の詳細については、「追従荷重(その2)」、および「追従荷重」をご参照ください。

Eficasによる変更点は、まず下図の通り、基本の時間ステップ幅を定義した「DEFI_LIST_REEL(pas)」の下に、「DEFI_LIST_INST(pas2)」を追加しました。

Auto011


「DEFI_LIST_INST」のオプション・パラメータはかなり多くあり、それら全てを把握するのは大変です。詳細はマニュアル

「u4.34.03.pdf」
http://www.code-aster.org/V2/doc/default/en/man_u/u4/u4.34.03.pdf

に書かれていますが、ここでは私の設定について上のEficas画面に沿って説明します。
---------------------------------------------

●DEFI_LIST
・METHODE
 AUTOにすると時間ステップ幅が自動計算されます。

・LIST_INST
 DEFI_LIST_REELで定義した、基本の時間ステップ幅を指定します。ここでは「pas」になります。これが初期値になります。今回は1秒間を20分割しているので、0.05になります。

・PAS_MINI
 許容される最小の時間ステップ幅です。これを下回ると計算が停止します。

・PAS_MAXI
 許容される最大の時間ステップ幅です。これを上回ると計算が停止します。余裕を見た大きな値にしておくほうが良いと思います。

・NB_PAS_MAXI
 許容される最大のステップ数です。デフォルト値はかなり多いので、これに引っかかることはあまり無いと思います。

●ECHEC (時間ステップを細かくする)

・EVENEMENT
 ERREURにすると、イテレーション回数オーバーや、マトリックスSingularなどの場合に適用されます。

・SUBD_PAS
 適用された場合、現在値の時間ステップ幅を何分割するかを指定します。
 デフォルトは4で、今回はそのままです。

・SUBD_PAS_MINI
 時間ステップ幅が分割される際の、最小値を指定します。
 デフォルトは0なので、無制限に分割されてしまいます。なんらかの小さな値を入れるほうが良いと思います。なお、上記のPAS_MINIとの違いなのですが、マニュアルによるとPAS_MINIは新しい時間ステップに入るときのみ適用されて、SUBD_PAS_MINIは「under-sapnning」のときにも適用される(?)とあるのですが、いまひとつ良く分かりません。とりあえず、両者とも0.0001にしておきました。

・SUBD_METHODE
分割の方法ですが、デフォルトのUNIFORMEのままです。おそらく「等分割(?)」でしょうか。

その他にも、「SUBD_COEF_PAS_1(デフォルト=1)」などのパラメータもありますが、これも良く分かりません。

●ADAPTATION (時間ステップを粗くする)
・EVENEMENT
 適用を判断する方法を指定するようですが、よく分からないので、これもデフォルトのSEUIL(閾値)のままです。

・NB_INCR_SEUIL: 2
・NOM_PARA: NB_ITER_NEWTON
・CRIT_COMP: LE
・VALE_I: 5

以上の組み合わせで、「ニュートン法のイテレーション回数が 5以下(LE)」のステップ(インクリメント)が、「2回続いたとき」に適用される、という条件になるようです。

・PCENT_AUGM
条件が適用されたとき、現在の時間ステップ幅を何パーセント増やすかという設定です。デフォルトの100(パーセント)のままにしていますので、条件が適用された場合には、時間ステップ幅は2倍になるということです。

---------------------------------------------

以上のように「DEFI_LIST_INST」を設定し、「pas2」という名前にしました。

それからこのpas2を、「STAT_NON_LINE」の「INCREMENT」-「LIST_INST」に設定しました。

Auto012

なお、結果出力の分割数設定(DEFI_LIST_REELの“output”)は20分割にしており、これを「ARCHIVAGE」に設定しています。時間ステップ幅が自動で変わることにより、出力分割数との不整合が心配されたのですが、時間ステップ幅が自動で変化しても、初期設定の分割点(DEFI_LIST_REEL(pas)で指定した初期の分割時点)では、必ず計算をしてくれるようです。
つまり今回の時間ステップ幅の初期値が20分割なので、0秒~1秒の間の0.05秒刻みの点では、必ず計算が入ります。なので、出力分割も20分割のままでOKです。出力分割については、時間ステップ幅の初期値に整合させておけば良いようです。

さて、計算結果ですが、.mess ファイルを見たところ、計算の最初の時点で早速時間ステップ幅の分割が起こったようです。

Auto013_2


初期値が0.05なので、0.05秒で計算していたところ、エラーが起こりました。文字化けしていますが、マトリックスの問題でしょうか。イテレーション回数のオーバーではないようです。
とにかくエラーなので分割が起こり、4分割の0.0125秒にて再計算され、収束しています。

追従荷重(その2)」では、均等に80分割の時間ステップ幅まで細かくしないと計算が停止していましたので、今回の初期値(20分割・0.05秒)では収束しないことは分かります。

その後も、時間ステップを調整しながら計算が進み、結局、74インクリメント・約99秒で計算完了しました。元の80等分割のケースでは、80インクリメント・約80秒だったので、計算時間は長くなっています。

元々、この問題は非線形性の変化が少なく、自動ステップの例題としては適切でなかったかもしれません。しかし、分割数を設定するための試行錯誤をしなくても済むので、トータルでは効率的と思います。なにより、途中で計算が止まる「ガッカリ感」がありません。
また、接触や摩擦を含むケースなど、計算過程の一部で非線形性が強くなる問題では、より有用だと思います。

「DEFI_LIST_INST」については、今回の事例以外にも様々な設定項目があります。問題に応じて細かく詰めていくことが出来そうです。

今回のメッシュファイル(mmed)とcommファイルは下記です。

auto01.zip

2012年7月18日 (水)

節点反力を集約する方法

材料引っ張り試験の模擬などを行う際は、拘束した節点の反力を集約することが必要になるのですが、その方法がなかなか分からずにいました。

しかし、公式フォーラムを検索していると、下記のやり取りを参考にすれば出来そうと言うことが分かりましたので、試してみました。

http://www.code-aster.org/forum2/viewtopic.php?id=14243

試してみたのは、以前実施したゴム材料の引っ張り試験です(下記リンク)。
http://salome-meca.cocolog-nifty.com/blog/2012/04/post-51b0.html

このときは、本当は変位制御式で反力を出力したかったのですが、反力の集約方法が分からず、やむをえず荷重制御式にしていました。

このときの条件ですが、下図のとおり、ゴム板の一方の切断面であるfixXのX変位を拘束し、反対側のloadP点に2[N]の荷重をかけていきました(fixY,fixZはそれぞれ対称面なので、これらも面内に拘束です)。
loadP点は、load面とX変位を連結しています。

Reac041

このとき、fixX面に属する節点の反力を集約すれば、loadP点の荷重と一致するはずです。それを試してみました。

ソフトウェアは、DEXCS2011-Salome-A3(64)(SALOME6.3 + Code_Aster STA10.5)の利用です。

基本的には、こちらのモデルを流用しています。
まずMeshなのですが、集約する節点を名前で区別するため、節点グループに「fixX」を追加しました。それでもう一度メッシュファイルを書き出しました。

Reac042

次にcommファイルをEficasで編集しました。
反力を出力するため「POST_RELEVE_T」と「IMPR_TABLE」を下図の通り追加しました。

Reac043

POST_RELEVE_TのRESULTANTEは、今回必要な「DX」のみとしていますが、必要に応じて「DY」や「DZ」も指定できるようです。

また、集約する対象の節点は、GROUP_NOに設定します。SALOMEのMeshで名前を付けた節点グループ名です。

このPOST_RELEVE_Tには「reacF」という名前を付けましたが、これをIMPR_TABLEにて、TABLEで指定します。

なお、この2つの詳細については、下記のマニュアルが参考になると思いますが、私もまだちゃんと読んでいませんので、今回の方法以外にもっと効率の良い設定があるかもしれません。

POST_RELEVE_T:u4.81.21.pdf
http://www.code-aster.org/V2/doc/v10/en/man_u/u4/u4.81.21.pdf

IMPR_TABLE:u4.91.03.pdf
http://www.code-aster.org/V2/doc/v10/en/man_u/u4/u4.91.03.pdf

さて、反力の集約に関する部分については以上なのですが、その他にも下記の部分を変更しています。

●計算時間短縮のため、荷重を与えるステップ分割数を、50から20に減らしました。
(DEFI_LIST_REELのNOMBRE)

●収束性向上のため、「STAT_NON_LINE」の「b_meth_newton」-「NEWTON」に、「REAC_ITER = 1」を追加しました。
これの詳細に関しては、こちらの記事を参照ください。

●得られる反力の精度を上げるため、収束許容値を厳しくしました。
「STAT_NON_LINE」の「CONVERGENCE」で「REGI_GLOB_RELA」の値を0.000001にしました。
この値が0.1や0.01など緩めの値だと、荷重と反力に若干のズレが生じます。
それはそれで正常な動作なので別に構わないのですが、今回は関しては誤差があると気持ち悪いので厳しくしています。まあ、気分の問題です。

以上の設定にて計算を実行しました。

その結果、「reac04.resu」ファイル(resuファイル)に下記のように出力されました。

Reac044

荷重はX方向に2.0 [N]を20分割して与えていますが、反力はX方向マイナス側に同じ値が出力されています。これをエクセルに取り込めば、fixXの反力値として使えそうです。

これで、変位制御式の引っ張り試験も模擬できるようになりそうです。機会があれば、また変制御式での反力についても試してみたいと思います。

今回のメッシュファイル(mmed)と、commファイルは、下記です。

reac04.zip

2012年7月 3日 (火)

要素の曲げ特性(超弾性材料)

先日実施した「追従荷重(その2)」と同じゴムの片持ち梁モデルを使い、曲げ問題に対して、要素の種類とメッシュの細かさを変えることで、結果がどの程度変わるのか検証しました。

ソフトウェアは、DEXCS2011-Salome-A3(64)(SALOME6.3 + Code_Aster STA10.5)の利用です。

ジオメトリモデルは、幅2.5mm、長さ10mm、厚さ1mmのゴム平板です。材料はムーニー材料で、「C10=0.133333、C01=0.033333、Nu=0.499」です。

境界条件を与えるためのグループは、X=0の面にfixX、それと反対側の端面にloadという名前をつけました。

Hype031

このfixX面は完全に固定し、load面に下向き(Z方向マイナス向き)の荷重0.01[N]を与えて曲げることにしました。荷重は「FORCE_FACE(単位面積当たりの荷重)」で与えるので、load面の面積5mm^2で割って、値は「-0.002」となります。

メッシュですが、まず六面体で平板の厚み方向に2層、4層、8層となるモデルを、1次要素、2次要素それぞれで作成しました。これは、SalomeのMesh機能で、1DのMax_Sizeを指定してコントロールしています。平板の厚みが1mmなので、Max_Sizeを0.5にすると2層、0.25にすると4層、0.125にすると8層になります。

このメッシュ方法だと、厚み方向だけで無く全体にメッシュが細かくなってしまいますが、実際、一定の方向だけ極端に分割を細かくすることは、アスペクト比が悪化するので、あまりやらないと思います。なので、とりあえず今回はこのようにしました。

また、四面体2次要素メッシュも作成しました。今回の平板のような単純な形状なら六面体メッシュを容易に作成できますが、多くの場合は四面体要素を使わざるを得ないためです。
四面体メッシュなので「層」にはなりませんが、「2層相当」としてMax_Size=0.5としたもの、「4層相当」として、Max_Size=0.25としたものを作成しました。8層相当は要素数が多くなりすぎたために今回は使いませんでした。

Hype032

四面体の1次要素はほとんど使うことがありませんので、検証対象に入れていません。

解析の設定は、荷重以外は「追従荷重(その2)」とほぼ同じのため、省略します。
サンプルとして、「六面体1次メッシュ」のメッシュファイル(mmed)とcommファイルを下記に置いています。

hype03.zip

それから、六面体2次要素に関しては「低減積分」要素でも試しました。低減積分要素にするには、「AFFE_MODELE」-「AFFE」-「b_mecanique」で「MODELISATION: 3D_SI」を指定します。通常(完全積分要素)はここが「3D」になります。

Hype033

六面体1次要素に関しては、ゴム材料・大変形では低減積分が使えないようです。1次要素だからこそ低減積分で曲げ特性を改善したいところですが、残念です。

変位の評価は、下図のように、上面の角部のZ方向(上下方向)変位で行いました。

Hype034

さて、結果ですが、下表の通りになりました。

Hype035

表では分かりにくいので、厚み方向の層数と変位の関係をグラフにしました。

Hype036

おそらく、変位5.4mm~5.5mm程度が正解だと思うのですが、六面体1次要素は目立って精度が悪いのが分かります。2次要素はいずれも大差なく、おそらく4層もあれば実用上問題ないかと思います。

下図は計算時間です。

Hype037


8層になると節点数が増えて計算時間がかなり長くなります(これは、厚み方向だけで無く、幅・長さ方向にもメッシュが細分化されているせいでもあります)。

2次要素で4層程度にとどめるのが、時間と精度のバランスが取れているように思います。低減積分要素は計算時間とメモリ量を減らせますが、それほど大きな削減量ではないようです。今回は検証していませんが、応力値も変わってくると思いますので、利用するかどうかはケース・バイ・ケースと思います。

1次要素の精度が悪いのは、Code_Asterの精度が悪いというわけではなく、標準的な1次要素が適用されているからだと思います。

多くの商用ソフトでは、曲げ特性を改善するための様々な手法が自動的に(デフォルトで)知らないうちに適用されていたりすることが多く、普段はあまり何も考えなくてもそこそこの精度が得られます。至れり尽くせりな感じですが、あまりに自動化されてしまうと、自分が実施した計算が、いったいどんなオプションで計算されているのか分からなくなる時があります。

その点、Code_Asterは硬派というか、ちゃんと理解してオプションを設定していく必要がありますが、そのために色々と調べるので勉強になります。私も今回この検証をやってみて、いまひとつ理解が足りないことに気がつきました。Code_Asterの理論解説や、自分が使っている商用ソフトのマニュアルをもう一度ちゃんと読まないといけないと思いました。

この辺りの話について理解するための書籍は色々有りますが、最近出版された下記のものがお勧めです。せん断ロッキングや体積ロッキング、低減積分、非適合要素などの解説もあります。

  「解析塾秘伝 非線形構造解析の学び方!
   著:石川覚志、日刊工業新聞社、 ¥2,600 + 税

(実は私もちょっと編集に関わっているので、宣伝も兼ねています。)

« 2012年6月 | トップページ | 2012年8月 »