OpenRadioss

【OpenRadioss入門】熱源移動による熱解析の設定方法

こんにちは(@t_kun_kamakiri

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

今回解説するのはこちらの熱源が移動する解析です。
設定としては適当で「とりあえず計算できればいいや!」くらいの感覚ですので実現象を再現しうる精度にはなっていないでしょうが、設定方法は参考になると思います。

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

WSL2

モデルはこちらです。

github

解析設定の手順

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

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

モデル作成

モデル作成はCADソフトを使います。
無償で使えるCADソフトはぱっと思いつくのは以下があります。

FreeCADはその名の通りCADソフトですが、メッシャーもついており出力形式もSTEP、igs、stlなど他のソフトでも読み込める形式で出力が可能です。

Salome-mecaは、shaperというCAD機能もありますしメッシャーもあり、メッシュ生成後に節点を接合したり修正したり非常に多くの機能があります。FreeCADと同じように様々な出力形式に対応しています。

Fusion360は、個人的に使ったことないのでホームページを見てください_(._.)_

今回はFreeCADを採用してモデルだけを作成してSTEPファイルをとして出力することを考えます。
※FreeCADの詳しい解説はここではしません。
メッシャーは別のソフトを使います。

作成するモデルは以下の寸法の100×100の正方形の面と、5×5の正方形の面です。
解析設定で両者の厚みを1mmに設定しますが、モデルは厚みのない面(サーフェス)として作成します。
ですので、5×5の正方形の面は1mm離してモデル作成をします。

FreeCADを起動してX-Y平面にスケッチを開始します。2つの面を作成します。作成した2つ面を選択してSTEPファイルとして出力します。

今回モデル作成はFreeCADで行いましたが、STEPファイルが出力できるソフトであれば何でもよいです。

メッシュ作成

続いてメッシュ作成に移ります。
メッシュ作成もFreeCADや別のソフトでもできるのですが、今回はGmshを使います。

Gmshを使う理由はOpenRadiossの形式として出力に対応しているからです。

先ほど出力したSTEPファイルをGmshで読み込んで、グループを作ります。

メッシュ作成を行います。
メッシュのサイズ感は適当に「Min = 1mm」「Max = 10mm」とします。
以下のようにRadiossの出力形式(.rad)で出力します。
先ほど作成したグループを出力してくとOpenRadiossでの節点グループに使えるので出力するようにしておきます。ただし、節点グループはLS-DYNA?の形式で「*SET_NODE_LIST」のキーワードでカンマ区切りで出力されます。
これをOpenRadiossの10カラム区切りに編集する必要があります。
以下の流れで10カラム区切りにします。以上の作成した節点グループは使うかわからないですが、「/GRNOD/NODE/」として書いておくと良いでしょう。

Gmshから出力したファイルを確認すると以下のようになっています。

「/GRNOD/NODE」の節点グループは後ほど拘束条件や接触条件に使います。

節点と要素

モデルはメッシュという要素に分割されて作られています。

今回は扱うモデルは厚みのないシェル要素で、下の絵のように要素は3つもしくは4つの節点から構成されています。

  • 節点は節点のIDと$x,y,z$
     の座標を指定
  • 要素は3つもしくは4つの節点IDを指定

節点が3つで構成される要素は「/SH3N」のキーワードを使い、節点が4つの場合は「/SHELL」というキーワードを使います。

また、PART IDが以下のように対応していることを確認します。

001thermal_0000.rad

Gmshでメッシュ作成からRadioss形式(.rad)に出力した際にシェル要素には自動的にPART IDと関連付けされています。
100×100の正方形の面にはPART ID = 2000001が、5×5の正方形の面にはPART ID = 2000002がそれぞれ自動で付けられています。
こちらはこのままのIDでも良いですし、本記事のように後ほど開設するように変更しても良いです。

節点、シェル要素(ソリッド要素)がPARTと紐づいていることを確認したら、次に行うのはPARTに対するプロパティと材料特性の設定です。

材料特性の設定

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

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

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

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

マニュアルを見ながら必要な項目を埋めていきます。
今回は使う材料カードは温度特性も解くことができるものですので、併せて「/HEAT/MAT」のキーワードも使います。

他にも/THERM_STRESS/MATという熱膨張を扱うキーワードも使えるのですが、ここでは使用しません。

MAT ID = 1は100×100の正方形の面の材料特性に使い、MAT ID = 2は5×5の正方形の面の材料特性に使うようにします。

5×5の正方形が熱源なので800[K]としています。
100×100の正方形は初期は398[K]という設定です。

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

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

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

  • SIG_Y = 0.09025:降伏応力
  • UTS = 0.175:工学応力UTS > σy
  • EUTS = 0.24:UTSにおける工学ひずみ
  • EPS_p_max = 0.75:破壊塑性ひずみ
  • ICC(設定なし) :デフォルトは1なのでひずみ速度効果あり

つまり、蘇生ひずみεpが0.75を超えるとシェル要素が削除されるというモデル設定です。
さらに詳しく知りたい方はマニュアルを参照ください。

テキスト設定の入手方法
GUIが無い場合はテキストによる設定の記述が必要となり、マニュアルを見ながら正確に設定を書くのは正直困難です。
設定のテンプレートがあれば楽ですよね。

その場合は、マニュアルの材料設定のページにいくと、以下のように材料設定が使われている例題の紹介があるので、どれかひとつ参考になるモデルからテキスト設定をコピペすることを考えます。

例題集のモデルはモデルファイルへのアクセスから入手できるので、必要なキーワードはここから入手していけばよいでしょう。

以上で材料設定が終わりです。

※ここでは材料特性の設定を記述しただけで、PARTとは紐づいていません。

プロパティの設定

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

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

  • Ishell = 1の場合:/DEF_SHELLが定義されていない場合のデフォルト。Q4、変形モードと剛体モードに直交する粘弾性アワグラスモード(Belytschko)
  • Thick = 1.0:厚み設定

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

※ここではプロパティの設定を記述しただけで、PARTとは紐づいていません。

PARTの作成

ここで先ほど記述したプロパティと材料特性をPARTに紐づけます。

PARTは要素(ビーム要素、シェル要素、ソリッド要素)で構成されており、PARTへはプロパティのIDと材料のIDを指定します。

設定は10カラム(10行)ごとに区切られているので、記述する位置を間違えないようにしましょう。「#」で始まる行はコメント行になるので、設定のメモを残しておくとわかりやすいです。

先ほどGmshでメッシュ作成からRadioss形式(.rad)に出力した際にシェル要素には自動的にPART IDと関連付けされていましたが、以下のように変更しておきましょう。

するとPART IDはそれぞれ100×100の正方形の面にはPART ID = 1が、5×5の正方形の面にはPART ID =2となります。
ここにプロパティのIDと材料のIDを書いてPARTと関連付けます。

これでPARTにプロパティと材料が紐づけられました。

解析設定としてPARTにプロパティと材料を指定すると使用されるため、例えば使う可能性があるプロパティや材料を設定の記述だけをテキストに書いておくこともできます。
その際はプロパティや材料はそれぞれIDが被らないようにしておきましょう。

熱源の拘束条件

今回の解析は5×5の正方形が熱源となりy軸方向には自由に動けるようにしますが、その他の軸と回転は拘束するようにします。
拘束条件を与えるキーワードして以下があります。

  • /BCS
    並進および回転移動に対する節点グループ上の境界条件を定義します。

こちらは節点グループに強制速度を与えるため節点グループに「/GRNOD」を使う必要があります。しかし、節点グループはGmshで出力後に編集したように既に作成済みです。
「/GRNOD/NODE/1」の1がIDになっているので、それを使えば良いでしょう。

拘束条件のキーワード「/BC」の「grnod_ID = 2」と節点グループのIDを与えることで拘束条件が設定されます。

  • Tra:x,y,zの拘束設定(0がフリー、1が拘束)
  • rot:ωx、ωy、ωzの拘束設定(0がフリー、1が拘束)

なので節点グループID = 2はy方向のみ自由に動くことができる設定になっています。

熱源の強制速度

今回の解析は5×5の正方形が熱源となりy軸方向に0.5[m/s]で移動するものとします。
強制移動をさせる場合のキーワードして以下があります。

今回は/IMPVELを使います。

こちらも「/IMPVEL」の「grnod_ID = 2」と節点グループのIDを与えることで拘束条件が設定されます。

funct_IDTには強制速度の時間変化を与える関数テーブルをしているする必要があります。

これにより関数は常に値を1.0を取るようになっていますが、「/IMPVEL」キーワードの「Scale_y = 0.5」とスケールしているため、熱源移動のy方向の速度は常に0.5[m/s]にしています。

補足ですが、「/GRNOD/NODE/」のように節点グループの節点を全て記述するという方法もありますが、PARTに属するすべての節点を指定していることと同じなので以下のように書いても良いです。

/GRNOD/PART

「item_ID1 = 2」とPART ID を指定すると、PART ID = 2の全ての節点がグループ化されます。こちらの方が設定自体は楽ですよね。

PARTの1つの節点群をグループにしたい場合は、「/GRNOD/NODE/」を使うなどモデルによって使い分けると良いでしょう。

100×100の正方形の面の拘束

100×100の正方形の面は完全に拘束します。

こちらも「grnod_ID」には節点グループの指定が必要ですが既に以下のように節点グループは作っています。

  • Tra:x,y,zの拘束設定(0がフリー、1が拘束)
  • rot:ωx、ωy、ωzの拘束設定(0がフリー、1が拘束)

なので節点グループID = 1は完全に拘束の設定です。

接触条件の設定

100×100の正方形と5×5の正方形の熱源(800K)との接触を定義します。

接触条件の設定は解析で理解が難しいものの1つなので、自分も理解が乏しいです。
接触の定義は大別するとラグランジュ法とペナルティ法があります。

今回はペナルティ法を使うため、OpenRadiossで設定できるペナルティ法をマニュアルから引用します。

熱伝導などを扱える接触は/INTER/TYPE7なのでこちらを使います。

こちらはサーフェスと節点群との接触を定義する必要があります。

  • スレーブ:Slav_idにはサーフェス
  • マスター:Mast_idには節点グループ

それぞれPART IDを指定することで、グループは簡単に作っています。
基本的な設定は以下のようになっています。

接触剛性$K$、温度依存する摩擦係数$\mu$、熱摩擦、熱交換、放射など、非常に多くのパターンの設定が可能ですので、詳しくはマニュアルを見ることをお勧めします。

熱伝達の設定

オプションとして外気との熱伝達(自然対流、強制対流)の考慮も行うことができます。
それには「/CONVEC」キーワードを使います。

  • SURF_IDでサーフェスを指定。
    「/SURF/PART/」でPARTを指定。
  • FUNCT_IDで参照温度の時間変化の関数IDを指定
    「/FUNCT」で関数テーブルを作成
  • Hで熱伝達率を設定

熱伝達は以下のように計算されます。

\begin{align*}
\dot{Q}=H\big(T-T_{inf}(T)\big)
\end{align*}

設定は記述しましたが、熱伝達率が大きすぎるとあまり熱伝達の様子が見れなかったので今回はH=0としています。つまり、熱伝達は起こっていません。

しかし、現実には外気との熱伝達があるため設定するほうが良いです。

節点の履歴

必要に応じて出力するようにしておくと便利です。
時刻歴応答「/TH」からはじまるキーワードを使います。

  • /TH/ACCEL:加速度
  • /TH/INTER:接触力など
  • /TH/NODE:節点の変位、速度、加速度
  • /TH/MONVOL:エアバッグの体積
  • /TH/PART:重心位置、慣性モーメントなど
  • /TH/RBODY:剛体の荷重
  • /TH/RWALL:剛体壁の接触力
  • /TH/SH3N:3節点のシェル要素の応力、ひずみ
  • /TH/SHEL:4節点のシェル要素の応力、ひずみ

以下は節点座標の時刻歴を出力する設定です。

以下は剛体の荷重の出力設定です。

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

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

Engine設定

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

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

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

001thermal_0001.rad

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

  • /ANIM/VECT/CONT:接触力
  • /ANIM/NODA/TEMP:節点の温度
  • /ANIM/SHEL/TEMP:シェル要素の温度

解析実行

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

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

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

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

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

これで計算投入用のGUIが立ち上がるので、アニメーションとデータ出力の項目にチェック(②)を入れて計算を実行します。
ここにチェックを入れておくことでA-ファイル(アニメーション用ファイル)はParaViewでも閲覧できるvtk形式に変換し、T-ファイル(データ用ファイル)はcsvファイルに変換してくれます。

OpenRadiossのインストール先は下記のスクリプトで自身の環境に合わせる必要があります。

openradioss_run_script_ps.sh

問題なく計算が進んでいれば大量の「A-ファイル」が出力されているはずです。

結果処理

先ほど計算実行時に「Anim-vtk」にチェックを入れたのでA-ファイルはvtk形式に変換されているので、それらをParaViewで読み込みます。

まとめ

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

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

関連記事もどうぞ

COMMENT

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