樹脂流動解析による金型設計および成形技術の高度化に関する研究
生産技術部
中野太郎

田中徹

本研究は，金型設計技術および成形技術の高度化を目的とし，樹脂流動解析を駆使した効率的な製品開発支援により，成形品の高品質化・低コスト化を短期間に達成する最適化手法を開発するものである．本年度は，樹脂材料の充填状態に関わる樹脂流動解析精度を向上させるため，樹脂物性パラメータのチューニング法について検討した．その結果，Chebyshev の直交多項式を適用した応答曲面モデルを利用することで，射出成形品重量の予測精度を向上させることができた．

1. はじめに

2. 応答曲面法の利用

プラスチック関連産業における生産活動のグロー

樹脂流動解析結果とその結果に影響を与えるフィバル化に伴い，国内外における市場競争は激しさをッティングパラメータの関係を数式化するために応増し，製品開発期間の短縮と製品コストの大幅な削答曲面法を利用する．応答曲面とは n 個の予測変数減が求められている．このような厳しい受注環境のxi (i = 1…n) から推定される応答 y の関係を近似式で中，品質・コスト両面の競争力を向上させることが緊表現したものである 2, 3) ．

急課題となっており，製品開発の期間短縮や成形不良の問題解決に寄与するものとして，樹脂流動解析の適用が拡大している．
樹脂流動解析技術は近年めざましい発展を遂げ，金型設計など製品開発の初期段階から積極的に利用されるようになってきたが，実際に高い精度で事前予測を行うには，複数の実験における実測結果とシミュレーション結果が一致するよう，様々な設計因子を微調整して解析結果を合わせ込む必要がある．

回帰係数を β とし，2 次多項式を応答関数とした場合，応答曲面式は式 (1) で与えられる．

y = β0 +nβi x i +i=1

整作業には専門的な知識が要求される 1) ため，高度な技術やノウハウを有する解析技術者でなければ，樹脂流動解析システムを十分に使いこなすことができないのが現状である．

βii x2i +i=1n

βij xi xj(1)i<j未知の係数 β は，回帰分析により最小二乗法を用いて求めることができる．回帰係数 β の推定に用いる実験点の総数を n，説明変数の数を k ，誤差を  とす
ると，その回帰モデルは行列表示で式 (2) となる．

Y = Xβ + 

しかしながら，これらの調整に用いる因子（以下，フィッティングパラメータとする）の選定や実際の調n⎡⎤⎢ y1 ⎥⎡⎢1⎢⎢ ⎥⎢⎢ ⎥⎢⎢ ⎥⎢1⎢ y2 ⎥⎢⎢ ⎥Y =⎢ ⎥ X=⎢
⎢
⎢ . ⎥
⎢1
⎢ .. ⎥⎢⎢ ⎥⎢⎢ ⎥⎣⎣ ⎦yn1(2)⎤⎡⎤

⎡

⎤

⎥⎥

⎢
⎢

⎥
⎥

⎢⎢⎥
⎥

x11x12···x1k ⎥x21..
.

x22..
.

···..
.

⎢⎢ ⎥⎥
x2k ⎥
⎢ β1 ⎥⎢ 2 ⎥⎥⎢⎢ ⎥⎥⎥⎥ β=⎢⎥ =⎢ ⎥.. ⎥⎢ . ⎥⎢ . ⎥.
.⎥⎢⎥⎥. ⎥ ⎢.⎥ ⎢⎢ . ⎥xn1xn2···xnk⎥⎦⎢ β0 ⎥⎢⎣βk⎥⎦⎢ 1 ⎥⎢⎣n⎥⎦

そこで，経験の少ないオペレータでも樹脂流動解
析結果と実験計測結果の整合性を向上させる最適な

ここで，Y は応答値ベクトル，X は測定点による行合わせ込みを容易に達成できるよう，特殊なスキル

列，β は回帰係数ベクトル， は誤差ベクトルである．

を必要としない，標準化された樹脂物性パラメータ

誤差の 2 乗和 S を式 (3) とし,チューニング手法の開発を行う．
本年度は，充填・保圧工程における樹脂流動解析精S = T  = (Y − Xβ)T (Y − Xβ)(3)

度を向上させるため，溶融樹脂材料の圧力 (P)・比容

β の推定値を K とすると，S の K に関する停留条件

積 (V)・温度 (T) の相互関係を示す PVT 特性に着目から以下の式 (4) が導かれる．

し，その物性値に関わるフィッティングパラメータのチューニング方法について検討した．

∂S∂KK=−2X T (Y − XK) = 0=(X T X)−1 X T Y(4)

このとき，推定値 K の分散共分散行列は式 (5) になる．

V (K)

=

E[(K − E[K])(K − E[K])T ]=E[(K − β)(K − β)T ]=E[((X T X)−1 X T )((X T X)−1 X T )T ]=(X T X)−1 X T E[T ]X(X T X)−1=(X T X)−1 X T (σ 2 I)X(X T X)−1=(X T X)−1 σ 2(5)図1

作成したテストピース

σ 2 は応答 y の誤差分散である．
式 (5) よ り，回 帰 係 数 の 分 散 を 小 さ く す る に は ，(X T X)−1 の 対 角 成 分 を 小 さ く す れ ば 良 い こ と が分 か る ．σ 2 が 観 測 値 の 誤 差 に 依 存 す る の に 対 し ，(X T X)−1 は実験点に起因する行列であるため，実験の組み合わせを工夫することで回帰係数分散を最小化することができる．
本研究では，X T X を対角行列化するために，L27直交表を利用することにした．ただし，単純な多項式を回帰モデルに用いた場合，直交表を用いても直交実験とはならない．そこで，Chebyshev の直交多表1

樹脂流動解析条件

樹脂流動解析ソフト
ワークステーション

Autodesk 社 MoldﬂowInsight2011
HP 社 (Windows7 64bit)

解析順序
樹脂材料
金型温度
樹脂温度
保圧／保持制御
冷却時間
メッシュタイプグローバルエッジ長
要素数

充填＋保圧
SA 社製 ポリプロピレン
45 ℃
210 ℃
15MPa - 10sec
15sec
3D メッシュ2mm
53498 個項式を回帰モデルに用い，本実験計画が直交実験となるよう配慮した．因子が A と B の 2 つで等間隔水準区間の場合，Chebyshev の直交関数は式 (6) となる

成型品重量と樹脂流動解析で算出した予測重量を比4–6)較し，その重量差を特性値として評価する．

．

3.2 テストピースの作成
y=K00 + K10 (A − Am ) + K01 (B − Bm ))射出成型機（日精樹脂工業株式会社製 NEX2000，+K11 (A − Am )(B − Bm )(a2 − 1)h2A+K20 (A − Am )2 −122(b − 1)h2B+K02 (B − Bm )2 −122(a − 1)h2A(B − Bm )+K21 (A − Am )2 −12+K12 (A − Am )(B − Bm )2 −(b2 − 1)h2B

平成 16 年度日本自転車振興会補助）を用いてテストピースを作成した. 樹脂材料には SA 社製のポリプロピレンを用いた．
成形したテストピース形状を図 1 に示す. 合計 50 個のサンプルを作成し，その中から抜き取った 10 個の成形品の重量を測定したところ，テストピースの平均重量は 7.8585g(標準偏差 0.0210g) であった.

12+ ···(6)：特性値
yA, B：因子（パラメーター）A m , Bm：各因子の水準平均
a, b：各因子の水準数
ha , h b：各因子の水準間隔K∗∗：多項式係数(K00 は全データの平均値)

3.3 樹脂流動解析
解析ソフトウェア MoldﬂowInsight2011（Autodesk社製，平成 22 年度電源立地地域交付金補助）を用いて樹脂流動解析を実施する. 解析条件は実際に成形したテストピースの射出成形条件に合わせて設定した.
解析条件の詳細を表 1 に示す.

3.4 フィッティングパラメータの選定
射出成形機を用いた一般的な成形工程の初期段階

3. 実験方法において，溶融樹脂は金型内に注入され，金型キャビ

3.1 評価方法

ティ空間内に樹脂が充填される（充填工程）．その後，樹脂流動解析結果と実測値との差が小さくなるよ

溶融樹脂の冷却に伴い充填された樹脂が収縮するのう，フィッティングパラメータの調整を行い解析精度で，その収縮を補填するために溶融樹脂を圧縮しての向上を図る．その際，実測値との整合性を検証する押し込む工程に移行する（保圧工程）．

指標として成型品の重量を採用した. 実際に成形した射出成形に使用される溶融樹脂材料は，温度や圧力に対して比容積が変化する粘弾性体であるため，成形の過程において樹脂材料の密度は常に変動している．充填や補填に必要な溶融樹脂量は，工程中の圧力・樹脂温度等に強く影響を受けることになるため，射出成形における流動現象を高精度に解析するためには，樹脂材料の密度や圧縮性に関わる精度の高い物性値データが必要になる．
そこで樹脂の圧縮性を示す状態方程式である PVT特性に着目し，フィッティングパラメータの調整を行うことにした．樹脂の物性値データを高精度化することで，充填および保圧工程における解析精度が向上し，より高い精度で充填される樹脂重量を予測できるようになると考えられる．

MoldﬂowInsight2011 において，圧力 P，比容積 V，温度 T の関係を示す PVT 特性は，式 (7) で示す Teit式と呼ばれる近似式でフィッティングして計算に用いられる．

表2

フィッティングパラメータの水準表

水準＼因子b1[m3 /kg]b2[m3 /kg·k]b4[1/k]

b5
[k]

123

1.220E-03
1.250E-03
1.280E-03

8.313E-07
8.750E-07
9.188E-07

4.750E-03
5.000E-03
5.250E-03

440450460

解析精度を向上させるためのフィッティングパラメータとして利用することにした．また，データ解析結果から得られる成形品重量の推定式が 2 次近似式となるよう，各パラメータの水準数は 3 水準とした．
水準表を表 2 に示す．現行の設定値を第 2 水準に取り，水準幅は等間隔に設定した．

3.5 実験計画
特性値（成形品の計量値と樹脂流動解析結果の重量誤差）に対するフィッティングパラメータの影響度を効率的に解析するために，L27 直交表に割り付けて解析試験を行う．フィッティングパラメータを割り付けた直交表を表 3 に示す．

P)} + V1(7)B
ここで，V0 = b1 + b2 (T − b5 )，B = b3 exp{−b4 (T −27 パターンの PVT 特性の近似式が作成される．こb5 )}，V1 = b6 exp{b7 (T − b5 ) + b8 P }，C は定数，b1

れらの物性値データを用いて樹脂流動解析を実施し，

V = V0 {1 − C ln(1 +直交表に従ったパラメータ値の組み合わせにより，〜b8 は樹脂毎に異なるデータフィッティング係数であフィティングパラメータの影響度を分析する．そのる．これらの係数の中から，樹脂の充填重量に大き際，回帰モデルとして Chebyshev の直交多項式を適な影響を与える 4 つの係数 (b1 , b2 , b4 , b5 ) をとりあげ，用した．

表3

フィティングパラメータを割り付けた L27 直交表

4. 結果および考察
4.1 分散分析と応答曲面推定式の作成
実験計画により得られる 27 個の樹脂流動解析結果を用いて，フィッティングパラメータの影響度を解析する分散分析を行った．分散分析表を表 4 に示す．表4 の中で因子名についている下付きの 1，2 という数
字はそれぞれ 1 次成分，2 次成分を示している．また，F 値に付いている ∗∗，∗ は，それぞれ有意水準 1 ％，5 ％で有意差があることを示している．なお，寄与が小さい因子は誤差項にプールした．

図2

分散分析において有意差が認められた因子および成形品重量推定値

次数成分を用いて応答曲面推定式を作成した．成形品重量 W (g) の推定式は式 (8) となる．

W

=市販の表計算ソフト等を用いて高速にデータ処理す0.2075 − 6706(b1 − 0.00125)ることが可能であるため，従来は適用が困難であっ−142700(b2 − 0.000000875)たさまざまな最適化手法を容易に利用できるように+0.005289(b5 − 450)なる．本研究では Microsoft 社 Excel のソルバー機能を用い，フィッティングパラメータの合わせ込みを行+6662000000(b1 − 0.001250)うことにした．Excel のソルバー機能は制約条件付き{(b4 − 0.005)2 − 0.00000004167} (8)の最適化問題を解くためのツールであり，比較的規

樹脂流動解析結果と式 (8) により推定される成形品

模の小さい数理計画問題であれば，ごく簡単な操作

重量を図 2 に示す．一部相違がみられるが，式 (8) で示される応答曲面推定式により，高い精度で樹脂流手順で最適解を求めることができる．

Excel のソルバー機能の設定項目において，変化さ動解析結果を予測できていることがわかる．

せる値をパラメータ b1 ，b2 ，b4 ，b5 に指定し，目標値4.2 フィッティングパラメータのチューニングとして実測値と式 (8) から得られる推定値の差が 0 に

このように，本来であれば特性値の推定には膨大近づくよう項目を設定して演算を行った．解析結果な計算量を伴う樹脂流動解析が必要であるが，応答を表 5 に示す．パラメータ b2 ，b4 は初期値設定のま曲面法を利用することで簡単な二次多項式に置き換

ま b1 ，b5 のみを若干調整するだけで誤差を最小化でえることができる．二次多項式で表される推定式は，きるという結果が得られた．
チューニングの効果を確認するために，得られた

表4

分散分析表

最適値を用いて PVT 特性近似式を作成し，樹脂流動解析を実施した．表 6 に結果を示す．チューニング前

FactorSum ofsquareDegree offreedom

Variance

F Ratio

b1 1
b1 2b2 1b2 2b4 1
b4 2
b5 1
b5 2
b1 1 × b2 1
b1 1 × b2 2
b1 2 × b2 1
b1 2 × b2 2
b1 1 × b4 1
b1 1 × b4 2
b1 2 × b4 1
b1 2 × b4 2
b2 1 × b4 1
b2 1 × b4 2
b2 2 × b4 1b2 2 × b4 2error

7.286E-01
2.940E-06
7.019E-04
3.315E-04
1.233E-05
4.774E-04
5.036E-02
5.097E-04
2.708E-06
4.921E-04
4.669E-07
1.463E-04
5.201E-06
6.242E-04
2.618E-05
1.435E-04
1.470E-06
6.834E-05
1.627E-05
5.070E-04
1.575E-03

101
10111010101
110111127.286E-01
7.019E-04
3.315E-04
4.774E-04
5.036E-02
5.097E-04
4.921E-04
1.463E-04
6.242E-04
2.618E-05
1.435E-04
6.834E-05
1.627E-05
5.070E-04
1.313E-04

5550.64 ∗∗5.35 ∗
2.53
3.63383.67 ∗∗
3.88
3.75
1.114.76 ∗0.20
1.090.52
0.12
3.86

total

7.846E-01

26は樹脂流動解析結果と成形品重量実測値との誤差が0.2110g であったが，フィッティングパラメータを最適化することで誤差 0.0105g に改善された．
表5

Excel のソルバーにより求められた最適値

初期値
最適値

b1
[m3 /kg]b2[m3 /kg·k]

b4
[1/k]

b5[k]

1.250E-03
1.262E-03

8.750E-07
8.750E-07

5.000E-03
5.000E-03

450449表6

成形品重量
実測値との差樹脂流動解析結果の比較
初期値

最適値

8.0695g
0.2110g7.8690g0.0105g

ら 0.0105g に改善することができた．
さらに，作成した応答曲面推定式を用いて応答の可視化を行った．応答曲面を三次元グラフ化することで，視覚的に分かりやすく特性値の挙動を把握できるため，現象を理解する手助けになるとともに，パラメータの微調整作業においても有用なデータとして活用可能である．
今後は，更なる樹脂流動解析の精度向上を目指し，樹脂の粘度に関わる近似式のフィッティングを行う予図3

パラメータ b2 と b4 の応答曲面

4.3 応答の可視化定である．また，特性値に影響を与える因子の抽出をより効率的に行うために，フィティングパラメータのスクリーニング方法について検討を進める．

応答曲面近似法のもう一つの利点として，応答の可視化が挙げられる．得られた推定式を用いて応答曲面を作成し，因子が特性値に与える影響を視覚的に表現することができる．パラメータ b1 ，b5 を最適値に固定した状態でパラメータ b2 ，b4 を変化させた
とき，特性値がどのように応答するか三次元グラフ化したものを図 3 に示す．

5. おわりに
本年度は，樹脂材料の充填状態に関わる樹脂流動解析精度を向上させるために，樹脂の圧縮性を示す状態方程式である PVT 特性のチューニング法について検討した．その結果，Chebyshev の直交多項式を用いて応答曲面推定式を作成し，Excel のソルバー機能によりフィッティングパラメータの最適化を行うこと
で，樹脂流動解析の成形品予測重量誤差を 0.2110g か

参考文献

1) 吉野睦 他，シミュレーションとＳＱＣ，日本規格協会，2009，p.40-42.
2) 轟章，応答曲面法，http://www.todoroki-lab.net/
response/responsesurface.pdf
3) 大塚裕俊 他，高硬度金型材のエンドミル加工（第６報），平成 14 年度大分県産業科学技術センター
研究報告，2002，p.57-62.
4) 柏村孝義 他，実験計画法による非線形問題の最適化，朝倉書店，1998，p.119-124.
5) 田口玄一，第３版
p.48-50.

実験計画法

上，丸善，1976，

6) 宮川雅巳，実験計画法特論，日科技連出版社，
2006，p.55-66.