Monoistに掲載されていた構造物の熱伝導解析の記事をSalomeMecaで非定常熱伝導解析を行いました。
シーズヒータを内蔵したプレートををシーズヒータで１８００秒(３０分)加熱したときのプレートの表面温度を計算しています。記事から読み取れない点は推定でモデル化を行っています。
・記事名：実例で学ぶステップアップ設計者CAE（3）
URL: 設計者はどんな視点で設計者CAEを進めていくべきか【ケース2：構造物の熱伝導解析】：実例で学ぶステップアップ設計者CAE（3）（1/3 ページ） – MONOist (itmedia.co.jp)
設計者はどんな視点で設計者CAEを進めていくべきか【ケース2：構造物の熱伝導解析】
・記事の筆者：土橋美博氏
・掲載メディア：MONOist
・解析ソフトウェア：SolidWorks（バージョンは不明）
＜動作確認環境＞
Salome-Meca 2019 (CAELinux2020Lite＋VirtualBox)
１．解析対象
下図に解析モデルを示します。
１．１ モデル寸法・物性値
ヒータプレートは、プレートとシーズヒータで構成されたアセンブリ品です。
シーズヒータは、中空のパイプに対してパイプ中央部に電熱線（ニクロム線）を配置して中空部にアルミナを充填した発熱体で、両端部は電熱線がなく、端面はガイシなどの絶縁体で封止をしています。
シーズヒータは、シーズヒータのモデルは記事に詳細がないことから、本ブログでは概略計算の位置づけでシーズヒータはステンレス製の棒状発熱体とみなしました。
シーズヒータとプレートのすきまは、記事と同様にないものとしました。
・プレート：150ｍｍｘ150ｍｍｘ30ｍｍ（ヒータ穴あり）
・シーズヒータ：φ10ｘ120ｍｍ(200W) 4本
・物性値は下表を参照ください。熱伝達率は記事と同一数値を使用しました。
ステンレスの物性値は、http://www.jssa.gr.jp/contents/faq-article/q6の数値から参照しました。
| 比熱容量kg/Kmm3 | 熱伝導率 W/mmK | 熱伝達率W/Kmm2 |
| 0.00397 | 0.0163 | 9.556×10-6 |
１．２ 解析モデル
- メッシュ：四面体１次要素
- 拘束条件：プレート上面のみ自然対流による熱伝達（周囲は断熱）
発熱源はシーズヒータ全体で発熱させる - 発熱量：シーズヒータ１本あたり200Wにて30分加熱
SalomeMeca上では単位体積の発熱量0.021W/m3として入力 - 初期温度：室温20℃、シーズヒータ１本あたり200Wにて30分加熱
- 解析：非定常熱伝導解析
- 出力：プレート上面の温度
２．モデルの作成
SalomeMecaでのモデル作成は、Geometoryにて形状を作成してから、Meshにて解析用のMeshを作成します。
Geometryの場合、フィーチャーベースのモデラではないことから、一度作成したモデルを寸法変更して形状変更ができません。（SalomeMeca2019からはShaperがあるのでこちらだとできるのかもしれません）
今回は、プレートをFreeCADで作成してSalomeMecaにインポートする方法でモデルを作成して、市販の3DCADとの連携がどの程度見込めるかを確かめました。
２．１ プレートのモデル作成
plateの形状をFreeCADで作成しエクスポートからstep方式で出力します。
出力は、メニューバーからのファイル(F)をクリックし、エクスポートをクリックします。
ダイアログボックスがでますので、名前を付けて保存します。
SalomeMecaからモデルをインポートします。
インポートは、メニューバーからFileをクリックし、Import→STEPをクリックします。
クリックすると、ダイアログボックスがでるので対象のモデルを選択します。
２．２ 形状作成
(1)Geometry上でヒータの形状を作成します。
(2)Plate上に熱伝達面やヒータとの接触面をExplodeやPartationなどを用いて、Face（面）として作成します。
この操作は、SalomeMecaでモデルを作成した場合と同じです。
２．３ メッシュの作成
Mesh上にてメッシュを作成します。
CreateMeshを用いて、Plateとheaterのメッシュを作成します。Plateは、Create submeshを用いてヒータと接触する面をヒータのメッシュサイズと同一にします。
メッシュのアセンブリコマンド（Create a Compound)を用いて、プレートとヒータを一つのメッシュにします。このときに、プレートとヒータの接合面のメッシュ（接点）を接合する設定にします。
初期値１e-05に対して許容差を緩くしています。トライアンドエラーで行いました。
３．解析結果
ヒータプレート上面の任意１ケ所赤〇部の温度を確認しました。温度結果は以下です。
- 469秒：179℃
- 900秒：307℃
- 1800秒: 558℃
表面温度は、記事の解析では表面温度が200℃を基準にシーズヒータをON-OFF制御するのに対して、本解析はヒータが常にONであることから表面温度が上がり続ける結果になりました。
表面の温度分布は、記事と同様にヒータの発熱体が入っている位置の表面温度が上昇していることから意図した解析になっていると考えます。
469℃の時の表面温度は179℃であり、プレートとシーズヒータが均一に469秒間800W加熱された場合141.8℃上昇し温度161.8℃に対して近いことから、単位や入力する発熱量のオーダミスはないものと考えます。この表面温度が均一加熱した温度に対して高いのはシーズヒータに近いためと推察します。
$${温度上昇(℃)} = \frac{熱容量(J)}{熱容量(J/K）}=\frac{ヒータの熱量(W) \times 時間(s)}{比熱容量(J/mm^3・K) \times 体積(mm^3)}$$
$$=\frac{800 \times 469}{0.00397\times 665580}= 141.8℃（温度は161.8℃）$$
４．まとめ、できたこと、できなかったこと
実務設計を想定したアセンブリ部品での非定常熱伝導解析を行うことができました。
本解析のモデルは、FreeCADで作成したモデルをstep形式へエクスポートしてからSalomeMecaへインポートして作成しました。step形式へエクスポートできるCADソフトがあれば、設計したモデルからSalomeMecaへインポートを行い解析を行うことが可能です。Fusion360やCATIA、SolidWorksなどたいていの3D-CADソフトにはこのエクスポート機能がありますので使用可能です。
今後やっていきたいことは、発熱量の制御で、記事にあった任意の位置の温度で発熱量の関数を作れるかを折をみて調べたいと考えています。
ヒータとプレートの接触面の接点は、接触面のメッシュ長さを同じにして作成し、Toreranceで接点をコンピュータの自動計算でつなぎました。シーズヒータの場合、パイプと充填剤のアルミナと発熱体と合成することから、ほかの方法でやれないか調べたいと考えています。
補足１：コマンドSOURCEの動作確認
sourceというコマンドを使って解析をしたので、動作を確認しました。
発熱量は、発熱量/発熱体の体積、単位体積当たりの発熱量で入力すると解説がありました。
動作を確かめるために、単純に水1gを1秒間に1℃ずつ上昇をさせた挙動を確認しました。
水を1g温度上昇させるのに4.182Jの熱量が必要なため、1分間に水に加える発熱量は4.182/1(J/s)=4.182Wです。
(1)計算条件
- 寸法：10mm(0.01m)の立方体
- 比熱容量：4182 J/kg C = 4182000 J/m3 C
比熱：4182J/kg ℃、水1kgの体積 0.001m3より計算 - 熱伝導率：0.602 W/m K
- 発熱量： 4182000 W/m3 C
水に加える発熱量4.182W、水1gの体積10-6m3より計算 - メッシュ：六面体１次要素
- 拘束条件：なし（断熱状態）
水全体に熱を加える。 - 初期温度：20℃、シーズヒータ１本あたり200Wにて50秒加熱
- 解析：非定常熱伝導解析
(2)結果
立方体の温度は、２５秒で45℃、50秒で70℃にであり立方体内での温度分布もないことから、Sourceは立方体全体に均一に加熱させているのわかりました。
よって、Sourceは単位体積当たりの発熱量で入力します。
補足２：コマンドファイル
ヒータプレートの伝熱解析
DEBUT(identifier='0:1',
LANG='EN')
mesh0 = LIRE_MAILLAGE(identifier='1:1',
FORMAT='MED',
UNITE=2)
model = AFFE_MODELE(identifier='2:1',
AFFE=_F(MODELISATION=('3D', ),
PHENOMENE='THERMIQUE',
TOUT='OUI'),
MAILLAGE=mesh0)
almina = DEFI_MATERIAU(identifier='3:1',
THER=_F(LAMBDA=0.0036,
RHO_CP=6.364e-06))
sus304 = DEFI_MATERIAU(identifier='4:1',
THER=_F(LAMBDA=0.016,
RHO_CP=0.00397))
fieldmat = AFFE_MATERIAU(identifier='5:1',
AFFE=_F(MATER=(sus304, ),
TOUT='OUI'),
MODELE=model)
timestep = DEFI_LIST_REEL(identifier='6:1',
DEBUT=0.0,
INTERVALLE=(_F(JUSQU_A=180.0,
PAS=18),
_F(JUSQU_A=900.0,
NOMBRE=30),
_F(JUSQU_A=1800.0,
NOMBRE=20)))
load0 = AFFE_CHAR_THER(identifier='7:1',
ECHANGE=_F(COEF_H=9.557e-06,
COEF_H_INF=9.557e-06,
COEF_H_SUP=9.557e-06,
GROUP_MA=('outer', ),
TEMP_EXT=273.0,
TEMP_EXT_INF=0.0,
TEMP_EXT_SUP=600.0),
MODELE=model,
SOURCE=_F(GROUP_MA=('GrMesh_2_Volumes', 'GrMesh_3_Volumes', 'GrMesh_4_Volumes', 'GrMesh_5_Volumes'),
SOUR=0.021))
resther = THER_LINEAIRE(identifier='8:1',
CHAM_MATER=fieldmat,
ETAT_INIT=_F(VALE=20.0),
EXCIT=_F(CHARGE=load0),
INCREMENT=_F(LIST_INST=timestep),
MODELE=model)
unnamed = CALC_CHAMP(identifier='9:1',
RESULTAT=resther,
THERMIQUE=('FLUX_ELGA', 'FLUX_ELNO', 'FLUX_NOEU'))
IMPR_RESU(identifier='10:1',
RESU=(_F(RESULTAT=resther),
_F(NOM_CHAM=('FLUX_NOEU', ),
RESULTAT=unnamed)),
UNITE=80)
FIN(identifier='11:1',
)
コマンドSOURCEの動作確認
DEBUT(identifier='0:1',
LANG='EN')
mesh0 = LIRE_MAILLAGE(identifier='1:1',
FORMAT='MED',
UNITE=3)
model = AFFE_MODELE(identifier='2:1',
AFFE=_F(MODELISATION=('3D', ),
PHENOMENE='THERMIQUE',
TOUT='OUI'),
MAILLAGE=mesh0)
almi = DEFI_MATERIAU(identifier='3:1',
THER=_F(LAMBDA=200.0,
RHO_CP=2430000.0))
almina = DEFI_MATERIAU(identifier='4:1',
THER=_F(LAMBDA=30.0,
RHO_CP=3120000.0))
water = DEFI_MATERIAU(identifier='5:1',
THER=_F(LAMBDA=0.602,
RHO_CP=4182000.0))
fieldmat = AFFE_MATERIAU(identifier='6:1',
AFFE=_F(MATER=(water, ),
TOUT='OUI'),
MODELE=model)
timestep = DEFI_LIST_REEL(identifier='7:1',
DEBUT=0.0,
INTERVALLE=_F(JUSQU_A=50.0,
NOMBRE=50))
load0 = AFFE_CHAR_THER(identifier='8:1',
MODELE=model,
SOURCE=_F(GROUP_MA=('Group_Of_All_Volumes', ),
SOUR=4182000.0))
resther = THER_LINEAIRE(identifier='9:1',
CHAM_MATER=fieldmat,
ETAT_INIT=_F(VALE=20.0),
EXCIT=_F(CHARGE=load0),
INCREMENT=_F(LIST_INST=timestep),
MODELE=model)
unnamed = CALC_CHAMP(identifier='10:1',
RESULTAT=resther,
THERMIQUE=('FLUX_ELGA', 'FLUX_ELNO', 'FLUX_NOEU'))
IMPR_RESU(identifier='11:1',
RESU=(_F(RESULTAT=resther),
_F(NOM_CHAM=('FLUX_NOEU', ),
RESULTAT=unnamed)),
UNITE=80)
FIN(identifier='12:1',
)