OpenRadioss

【OpenRadioss入門】 Simplified Carのポール衝突

こんにちは(@t_kun_kamakiri

Altairのサイトにこれまで商用利用としていたRadiossを2022年9月からオープンソース化するという内容の記事が投稿されていました。RadiossといえばLS-DYNAやPAM-Crashなどと並ぶAltairの衝突系(衝撃解析)の解析ソフトのイメージが強いですが、それが無償で使えるようになったというのはなかなかの衝撃でした。

今回解説するのはこちらのシンプル(おもちゃみたいな)車のポール衝突解析です。
設定としては適当で「とりあえず計算できればいいや!」くらいの感覚ですので実現象を再現しうる精度にはなっていないでしょうが、設定方法は参考になると思います。

Altairの例題集は大変参考になります。
しかし、モデル作成からRadiossは調べてキーワードというカードで設定を記述する必要があります。
そのため、GUI操作で順番にボタンをポチポチ押して行って設定するという感じではなく、適切な解析条件を自分で選択して解析設定を行う必要があります。
無償のOpenRadiossにはAltair社のHyper WorksのようなGUIが無いため、設定に苦労します。本記事では、できるだけテキスト編集をメインにOpenRadiossの基本的な設定を解説します。

WSL2

解析設定の手順

解析設定にはStarterファイルとEngineファイルの作成が必要です。
今回は以下のファイル名で設定ファイルを作成します。

単位系

物理量の値の大きさには単位系を意識する必要があります。

今回はmm-ton-secの単位系で解析設定を行います。

モデルのNastranからRadioss形式への変換

モデルは公式例題集からNastranとして入手することができます。

これをRadioss形式へ変換します。

NastranからRadioss形式への変換はGmshを使います。

変換の仕方は簡単でGmshを起動して、full_car.nasを読み込んで、Radioss形式としてエクスポートするだけです。

名前はfull_car.incに変えておきます。

full_car.inc

Node(節点)とShell(シェル要素)にIDが与えられています。

構成要素は以下のようになっています。

作成したモデルfull_car.incをメインのファイル002_full_car0000_0000.radでインクルードします。

002_full_car0000_0000.rad

PARTの確認

GmshでPARTの確認をしようとしましたが、どこにPART IDの記載があるのか全く分からなかったので、今回はNastran形式でも読み込めるLS-PrePostを使います。

しかし、GmshからRadioss形式にしたファイルを見るとPART IDが2000001から始まっているので、対応関係には注意が必要です。
(勝手にリナンバーされたんでしょうか)

LS-Dynaのキーワードを確認するとすでに、SECTION(シェル要素のプロパティ)とMAT(材料定義)がされていました・・・

こちらを参考にOpenRadioosの形に修正していきます。
Partが30個あるので面倒な作業ですが・・・

まずはWINDSHIELDから設定します。

WINDSHIELD

WINDSHIELDは以下の部分です。

パートの設定

PARTは以下のようにプロパティと材料のIDを指定します。

プロパティと材料はIDでPARTに紐づけるようにします。

プロパティの設定

シェル要素に対して要素の定義と厚みの設定を行います。

厚みは5mmにします。

このキーワードは、/PROP/TYPE1(SHELL)から引用すると、

3節点または4節点のシェル要素で使用されるシェルプロパティセットを記述します。Belytschko、QBAT、またはQEPHのシェル定式化が利用可能です。

ほとんどがデフォルト設定でよく、すべて設定しなければならないというわけではありません。
大事なところは以下です。

  • Ishell = 0の場合:4節点シェル要素定式化フラグ
  • Thick = 0.3:厚み設定

その他の設定は必要に応じてマニュアルを見て設定します。

材料特性の設定

Elasto-plastic Material /MAT/PLAS_JOHNS (WINDSHIELD)

  • [Rho_Initial] Initial Density = 2.5×10-9ton/mm3
  • [E] Young’s Modulus = 76000 MPa
  • [nu] Poisson’s Ratio = 0.3
  • [0] Yield Stress = 192 MPa
  • [K] Hardening Parameter = 220 MPa
  • [n] Hardening Exponent = 0.32

OpenRadiossには多くの材料特性の設定が用意されており、こちらの材料一覧から確認することができます。

今回使うのは以下の材料設定です。

則名 グループ タイプ モデル概要 適用分野 粘性 直交異方性 組み込み熱オプション** ひずみ速度 定式化
/MAT/LAW2 (PLAS_JOHNS)* 弾塑性材料 金属 Johnson-Cook材料モデルを使用した等方性弾塑性材料 金属、弾塑性材料、軟鋼、ひずみ速度依存、温度依存、等方および移動硬化、熱硬化性、低合金鋼、ST52、およびDP600 × × Lagrange

この材料則は、Johnson-Cook材料モデルを使用して等方性弾塑性材料を表します。

以下、こちらからの引用です。

1.これは、真の応力およびひずみ出力を伴うひずみ速度および温度効果を含む弾塑性材料モデルです。

2.このモデルでは、相当応力が塑性降伏応力よりも低い時、材料は線形弾性材料として挙動します。もっと高い応力値では、材料は塑性挙動で、応力は以下のように計算されます:

\begin{align*}
\sigma = \big(a+b\varepsilon_{p}^{n}\big)\bigg(1+c\ln\frac{\dot{\varepsilon}}{\dot{\varepsilon}_{0}}\bigg)\big(1-(T^{*})^{m}\big)
\end{align*}

$T^{*}=\frac{T-T_{r}}{T_{melt}-T_{r}}$

  • $\varepsilon_{p}$:塑性ひずみ
  • $\dot{\varepsilon}$:ひずみ速度
  • $T$:温度
  • $T_{r}$:雰囲気温度
  • $T_{melt}$:溶融温度

3.Iflag=1の場合、σy, UTSおよび に実験工学応力およびひずみデータを入力でき、パラメータabおよびnが計算されてStarter出力ファイルに出力されます。

が1つの積分点での値に到達してから、要素タイプに基づく場合:
シェル要素:対応するシェル要素が削除されます。

今回は以下のように設定しています。

  • RHO_I = 2.5E-9:密度
  • E = 76000:ヤング率
  • a= 192:応力の式の定数
  • b =220:応力の式の定数
  • n = 0.32:塑性硬化指数
  • SIG_max0 = 0:最大応力
  • ICC(設定なし) :デフォルト0なのでひずみ速度効果なし

さらに詳しく知りたい方はマニュアルを参照ください。

タイヤの設定

タイヤはご存じこちらです。

パートの設定

パートにプロパティと材料を紐づけます。

プロパティの設定

プロパティはWINDSHIELDと同じシェル要素です。

厚みは10mmにします。

材料特性の設定

Elasto-plastic Material /MAT/PLAS_JOHNS (RUBBER)

  • [Rho_Initial] Initial Density = 2×10-9ton/mm3
  • [E] Young’s Modulus = 200 MPa
  • [nu] Poisson’s Ratio = 0.49
  • [0] Yield Stress = 1e30MPa
  • [n] Hardening Exponent = 1

今回は以下のように設定しています。

  • RHO_I = 2.0E-9:密度
  • E = 200:ヤング率
  • a= 1e30:応力の式の定数
  • b =0:応力の式の定数
  • n = 1:塑性硬化指数
  • SIG_max0 = 0:最大応力
  • ICC(設定なし) :デフォルト0なのでひずみ速度効果なし

さらに詳しく知りたい方はマニュアルを参照ください。

それ以外

それ以外のパートを設定します。
それ以外というのは、WINDSHIELDとRUBBER以外のことです。

中にはSOLID要素も含まれます。

/BRICKの引用はこちら。

6面体ソリッド要素または8節点の厚肉シェル要素を定義します。要素タイプおよび定式化は/PROPカード上で定義されます。

パートの設定

パートにプロパティと材料を紐づけます。

シェル要素のプロパティはID = 1を使用します。
ソリッド要素のBRICKとなっているパートはプロパティ ID = 3にしています。

プロパティの設定

シェル要素に対するプロパティはWINDSHIELDと同じにします。
厚みはわからないので、妥協してすべて0.3mmにします。

ソリッド要素に対するプロパティは以下とします。

//PROP/TYPE20 (TSHELL)の引用はこちら。

このプロパティは、一般的な厚肉シェルプロパティセットの定義に使用されます。

Isolid = 14:HA8、8節点厚肉シェル要素、完全積分、さまざまな数のGauss点、共回転系定式化、せん断ロックフリー

他の例題でも使われているので確認してみしょう。

材料特性の設定

Elasto-plastic Material /MAT/PLAS_JOHNS (STEEL)

  • [Rho_Initial] Initial Density = 7.9×10-9ton/mm3
  • [E] Young’s Modulus = 210000 MPa
  • [nu] Poisson’s Ratio = 0.3
  • [0] Yield Stress = 200 MPa
  • [K] Hardening Parameter = 450 MPa
  • [n] Hardening Exponent = 0.5
  • [SIG_max] Maximum Stress = 425 MPa

今回は以下のように設定しています。

  • RHO_I = 2.0E-9:密度
  • E = 200:ヤング率
  • a= 200:応力の式の定数
  • b =450:応力の式の定数
  • n = 0.5:塑性硬化指数
  • SIG_max0 = 425:最大応力
  • ICC(設定なし) :デフォルト0なのでひずみ速度効果なし

さらに詳しく知りたい方はマニュアルを参照ください。

ポールの設定

ポールの設定は3Dモデルを作成するではなく、キーワードで作ることができます。

このキーワードは、/RWALLから引用すると、

無限平面、無限円筒、球形、および平行四辺形の剛壁を定義します。

  1. Set the Distance to search for main nodes to 1500.
    剛壁のセカンダリ節点は、節点のグループや初期段階で剛壁からDsearchほど下にある節点として定義できます。

セカンダリ検索の距離はポールとの接触する距離のこと1500mm = 1.5mは大きすぎるのでいったい何のことだろう・・・・

イメージを書くとこんな感じですかね。
ポールが車の真ん中までくることはないので、接触する距離を少し広めにしています。

その他の例題にも使われているので、確認してみましょう。

接触定義

接触定義はラジエータとエンジンの接触定義を行います。

このキーワードは、/INTER/TYPE7から引用すると、

インターフェースTYPE7は、メインサーフェスとセカンダリ節点グループ間の接触をモデル化する多用途衝撃インターフェースです。熱伝導および熱摩擦を考慮することもできます。

  • Istr = 2:インターフェース剛性は、メインおよびセカンダリ剛性の平均($K_{n}=\frac{K_{m}+K_{s}}{2}$)です。
    ここで、メインセグメント剛性$K_{m}$、セカンダリト剛性$K_{s}$
  • Fric = 0.2:Coulomb摩擦
  • Gapmin = 0.7 mm:最小ギャップ
  • Secon_idに節点グループでエンジンの表面の節点を指定します。
  • /GRNOD/SURFでサーフェスグループを指定しています。
    さらにサーフェスグループはパートサーフェース/SURF/PART/EXTから指定しています。
  • Main_idにはサーフェスグループでラジエータの表面を指定します。パートの外部パートサーフェスの指定は/SURF/PART/EXTで行います。
節点や表面(サーフェス)からグループを作ることができますが、GUIがない場合に節点グループやサーフェスグループを作るのは困難です(LS-PrePostやGmshを使う手がありますが・・・)。
しかし、上記のようにパートから節点やサーフェスグループを指定することができますので、そちらの方法を使います。
このようなIDの対応関係を理解することが非常に大事です。

初速

車の全体節点に初速を与えます。

単位系がmm-ton-secなので速度の単位はmm/secになります。
値だけで見ると速度m/sの単位の1000倍を設定する必要があるので注意しましょう。

35mph = 56.327km/h = 35*1000/60/60 = 15.64 m/s ≒ 156000 mm/sec

※35mphは35マイルという意味です。

このキーワードは、/INIVELから引用すると、

節点のグループに対する初速度を定義します。

節点グループに与えるというのが重要ですね。
ここでも全ての節点を選択してグループを作成するのではなく、/GRNOD/PARTでパートを指定して全節点グループを作成しています。
/INIVEL/TRAとすることで並進速度を初速に与えています。
フィールド 内容 SI 単位の例
type 初速度のタイプ。

  • TRA:材料の並進速度
  • ROT:材料の回転速度
  • T+G:材料の並進およびグリッド速度(ALE材料に対してのみ使用されます)
  • GRID:材料のグリッド速度(ALE材料に対してのみ使用されます)

節点のヒストリー

節点の時刻歴変化などの情報を出力するために節点を設定しておきます。

節点はLS-PrePostやGmshなど適当なGUIなどを使って

どの部分?って感じですが、ここです。

このキーワードは、/TH/NODEから引用すると、

節点の時刻歴を記述します。

DEFを指定していると「DX, DY, DZ, VX, VY, VZ」が出力されるようです。

 

必要に応じて設定を追加すれば良いでしょう。

以上までの設定を002_full_car0000_0000.radに保存します。
これがRadiossのstarter inputファイルになります。

Engine設定

以上まではRadiossのstarter inputファイルの設定でしたが、Radioss Engineファイルはシミュレーションの解析制御と出力を記述するファイルです。

ここでは最低限の設定だけを記述します。

  • /RUN:計算実行に関するキーワード。60ms (0.06秒)までの解析
  • /DT/NODA/CST:モデルの時間ステップを増やすためのマススケーリング
  • /TFILE:時刻歴ファイルであるT-ファイルの書き出し時間間隔
  • /ANIM/DT:アニメーションファイル(A-ファイル)を出力する開始時間と出力間隔(3ms)

002_full_car0000_0001.rad

さらに詳しく知りたい方はマニュアルを参考にしましょう。
今回は歪みや応力もアニメーションにしたかったので以下を追加しました。

  • /ANIM/ELEM/EPSP:塑性歪み
  • /ANIM/ELEM/VONM:フォン・ミーゼス応力
  • /ANIM/ELEM/ENER:比エネルギー密度
  • /ANIM/ELEM/HOURG:アワーグラスエネルギー
  • /ANIM/VECT/DISP:変位
  • /ANIM/VECT/FOPT:剛体、剛壁、および断面の力とモーメント
  • /ANIM/VECT/CONT:接触力
  • /ANIM/NODA/DMAS:ダンピング

解析実行

以上により解析の設定が終わりました。

解析に必要なファイルは以下の2つです。

  • 002_full_car0000_0001.rad:Starter ファイル(解析設定ファイル)
  • 002_full_car0000_0001.rad:Engine ファイル(解析制御と出力設定ファイル)

OpenRadiossの計算実行は、こちらで公開されているPythonで作られたGUIを使うと簡単に行えます。

WSLを使っているのでLinux環境での実行方法を参考にします。

4コアで計算させると1分ほどで計算が終わります。

結果の確認

節点の時刻歴

/TH/NODEで設定した節点の時刻歴が002_full_car0000T01.csvに出力されています。

追々グラフにまとめます。

ParaViewで確認

vtkファイルが出力されているのでParaViewで読み込んでアニメーションにするとこんな感じになります。

LS-PrePostで確認

d3plotに変換してアニメーションにするとこんな感じ。

まとめ

本記事ではOpenRadiossの設定をできるだけマニュアルを見ながら、テキストのみで解析実行できるように解説を行いました。

途中、節点グループの作成には困りましたが無理やりLS-PrePostなどを使ったりしましたので、そのあたりの良い方法は要検討です。
また、解析途中にモデルの確認ができないというのはかなり解析設定を困難にしているので、せめてメッシュ状態などを確認できる標準的な方法が今後の課題です。

とはいえOpenRadioss自体は機能も豊富なので、今後の発展に期待大です。

関連記事もどうぞ

COMMENT

メールアドレスが公開されることはありません。 が付いている欄は必須項目です