(57)【要約】【課題】予測精度を向上できる粉体の挙動予測方法を提供する。【解決手段】粒子２が塑性変形したときの反発係数ηを算出し、反発係数ηの算出結果に基づいて粒子２の反発力Ｆ１を算出するステップと、粒子２に発生した応力に応じて摩擦係数μを算出し、摩擦係数μの算出結果に基づいて粒子２が擦れて移動するときの摩擦力Ｆ２を算出するステップと、粒子同士の衝突によって発生する静電気力に基づいて、互いに離間する二つの粒子間の付着力Ｆ３を算出するステップと、の少なくともいずれか一つを行うことにより、個々の粒子２に対して作用する外力Ｆを算出する。【選択図】図１
2014154045.tif 000002
【特許請求の範囲】
【請求項１】
個々の粒子に対して作用する外力を算出し、前記外力を基に立てた運動方程式を解くことで、前記個々の粒子の速度、変位、および位置を算出するステップを繰り返し行うことにより、多数の前記粒子からなる粉体の挙動を予測する粉体の挙動予測方法であって、
前記粒子が塑性変形したときの反発係数を算出し、前記反発係数の算出結果に基づいて前記粒子の反発力を算出するステップと、
前記粒子に発生した応力に応じて摩擦係数を算出し、前記摩擦係数の算出結果に基づいて前記粒子が擦れて移動するときの摩擦力を算出するステップと、
前記粒子同士の衝突によって発生する静電気力に基づいて、互いに離間する二つの前記粒子間の付着力を算出するステップと、
の少なくともいずれか一つを行うことにより、前記個々の粒子に対して作用する外力を算出する、
粉体の挙動予測方法。
【請求項２】
以下の数式１を用いて前記粒子が塑性変形したときの反発係数を算出する、
請求項１に記載の粉体の挙動予測方法。
【数１】
2014154045.tif 000011
前記数式１において、ηは前記反発係数であり、ａ1〜ａ3は前記粉体の種類毎に設定される固定値であり、δは前記粒子の変形量であり、ｄは前記粒子の径である。
【請求項３】
以下の数式２を用いて前記粒子に発生した応力に応じた摩擦係数を算出する、
請求項１に記載の粉体の挙動予測方法。
【数２】
2014154045.tif 000012
前記数式２において、μは前記摩擦係数であり、ｂ1〜ｂ3は前記粉体の種類毎に設定される固定値であり、δは前記粒子の変形量であり、ｄは前記粒子の径である。
【請求項４】
以下の数式３を用いて前記互いに離間する二つの粒子の電荷をそれぞれ算出するとともに、以下の数式４を用いて前記二つの粒子間の付着力を算出する、
請求項１に記載の粉体の挙動予測方法。
【数３】
2014154045.tif 000013
前記数式３において、Ｑnは前記粒子の電荷、Ｑ0は前記粒子の元の電荷、ｃは前記粉体の種類毎に設定される固定値、Ｎは前記粒子同士の衝突回数、ｄは前記粒子の径である。
【数４】
2014154045.tif 000014
前記数式４において、Ｆ3は前記粒子間の付着力、ε0は真空の誘電率、εrは比誘電率、Ｑ1・Ｑ2は前記二つの粒子の電荷、ｒは前記二つの粒子の中心間の距離である。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、多数の粒子からなる粉体の挙動を予測する粉体の挙動予測方法に関する。
【背景技術】
【０００２】
従来から、多数の粒子からなる粉体を扱う工程、例えば、焼結体の製造工程等においては、粉体の挙動を把握するために、作業者等がカン、コツにて製造条件を調整しながら実機で焼結体を繰り返し製造している。この場合には、粉体の挙動を把握するために多大な工数が必要となる。
そこで、離散要素法を用いて粉体の挙動を予測する取組がなされている。
【０００３】
離散要素法では、粒子に作用する力（例えば、粒子の反発力および粒子が擦れて移動するときの摩擦力等）を基に運動方程式を立てて、粒子の加速度を求める。そして、離散要素法では、所定の時間増分に対して加速度を数値積分することで、粒子の速度および変位等を求める。
個別要素法では、このような計算処理を繰り返し行うことで、粉体の挙動を予測する。
【０００４】
特許文献１に開示される技術では、解析領域内に解析対象粒子を配置するとともに、解析領域の境界に運動を拘束した仮想粒子を配置する。特許文献１に開示される技術では、解析領域の境界、つまり、壁面等の表面粗さに応じて仮想粒子の径および配置間隔を設定する。
そして、特許文献１に開示される技術では、転写電界による静電力、付着力、重力、および空気抵抗力等を粒子に作用する力として運動方程式を立て、離散要素法を用いて粉体の挙動を予測する。これにより、特許文献１に開示される技術では、壁面等の表面粗さを考慮して、粉体の挙動を予測する。
【０００５】
ここで、焼結体の製造工程においては、粉末に対して荷重をかけて圧粉体を成形する成形工程が行われる。すなわち、成形工程では、前記荷重によって粒子が大きく変形し、粒子に大きな反発力が発生する。
離散要素法では、一般的に粒子が弾性変形していることを想定してこのような粒子の反発力を算出する。つまり、この場合、粒子の反発力は、粉体の種類（粉種）毎に固定の定数として設定される反発係数に基づいて算出される。
しかし、成形工程において、多くの粒子は塑性変形する。従って、成形工程において、前記反発係数は、粒子が降伏する分だけ変動する。つまり、固定の反発係数に基づいて反発力を算出するだけでは、成形工程における粉体の挙動を正確に再現できない可能性がある。
【０００６】
また、離散要素法では、一般的に、粉体の種類毎に固定の定数として設定される摩擦係数に基づいて、一の粒子が他の粒子等に擦れて移動するときの摩擦力を算出する。
しかし、成形工程のように粉体に対して荷重をかける場合には、粒子に発生する応力によって摩擦係数が変動する。つまり、固定の摩擦係数に基づいて摩擦力を算出するだけでは、粒子が擦れて移動するときの粉体の挙動を正確に再現できない可能性がある。
【０００７】
ここで、粉体を搬送する場合（例えば、粉体を貯溜する粉箱を移動させる場合等）には、粒子同士が衝突して静電気力が発生し、その結果、粉体の流動性が悪化する。
特許文献１に開示される技術では、このような流動性の悪化を考慮していないため、粉体を搬送するときの粉体の挙動を正確に再現できない可能性がある。
【０００８】
以上のように、従来技術においては、粉体の挙動を予測するときに把握している現象が少ないため、粉体の挙動を正確に再現できず、その結果、予測精度が低くなってしまう可能性があった。
この場合には、実用化できる程度まで予測精度を向上させることができない可能性があった。つまり、この場合には、実機による試験等で粉体の挙動を把握しなければならない可能性があった。
【先行技術文献】
【特許文献】
【０００９】
【特許文献１】特開２００６−２５９９１１号公報
【発明の概要】
【発明が解決しようとする課題】
【００１０】
本発明は、以上の如き状況を鑑みてなされたものであり、予測精度を向上できる粉体の挙動予測方法を提供するものである。
【課題を解決するための手段】
【００１１】
請求項１においては、個々の粒子に対して作用する外力を算出し、前記外力を基に立てた運動方程式を解くことで、前記個々の粒子の速度、変位、および位置を算出するステップを繰り返し行うことにより、多数の前記粒子からなる粉体の挙動を予測する粉体の挙動予測方法であって、前記粒子が塑性変形したときの反発係数を算出し、前記反発係数の算出結果に基づいて前記粒子の反発力を算出するステップと、前記粒子に発生した応力に応じて摩擦係数を算出し、前記摩擦係数の算出結果に基づいて前記粒子が擦れて移動するときの摩擦力を算出するステップと、前記粒子同士の衝突によって発生する静電気力に基づいて、互いに離間する二つの前記粒子間の付着力を算出するステップと、の少なくともいずれか一つを行うことにより、前記個々の粒子に対して作用する外力を算出する、ものである。
【００１２】
請求項２においては、以下の数式１を用いて前記粒子が塑性変形したときの反発係数を算出する、ものである。
【数１】
2014154045.tif 000003
前記数式１において、ηは前記反発係数であり、ａ1〜ａ3は前記粉体の種類毎に設定される固定値であり、δは前記粒子の変形量であり、ｄは前記粒子の径である。
【００１３】
請求項３においては、以下の数式２を用いて前記粒子に発生した応力に応じた摩擦係数を算出する、ものである。
【数２】
2014154045.tif 000004
前記数式２において、μは前記摩擦係数であり、ｂ1〜ｂ3は前記粉体の種類毎に設定される固定値であり、δは前記粒子の変形量であり、ｄは前記粒子の径である。
【００１４】
請求項４においては、以下の数式３を用いて前記互いに離間する二つの粒子の電荷をそれぞれ算出するとともに、以下の数式４を用いて前記二つの粒子間の付着力を算出する、ものである。
【数３】
2014154045.tif 000005
前記数式３において、Ｑnは前記粒子の電荷、Ｑ0は前記粒子の元の電荷、ｃは前記粉体の種類毎に設定される固定値、Ｎは前記粒子同士の衝突回数、ｄは前記粒子の径である。
【数４】
2014154045.tif 000006
前記数式４において、Ｆ3は前記粒子間の付着力、ε0は真空の誘電率、εrは比誘電率、Ｑ1・Ｑ2は前記二つの粒子の電荷、ｒは前記二つの粒子の中心間の距離である。
【発明の効果】
【００１５】
本発明は、予測精度を向上できる、という効果を奏する。
【図面の簡単な説明】
【００１６】
【図１】粉体の挙動予測方法の手順を示す説明図。
【図２】仮想空間にて粒子が接触している状態を示す説明図。
【図３】粒子の反発力を示す説明図。
【図４】ひずみおよび押込み荷重の関係を取得する様子を示す説明図。
【図５】ひずみおよび押込み荷重の関係を取得した結果を示す説明図。
【図６】圧粉体の部位を示す斜視図。
【図７】粒子の体積含有率から密度分布を算出する様子を示す説明図。
【図８】粒子が擦れて移動するときの摩擦力を示す説明図。
【図９】粒子の摩擦抵抗を測定する測定装置を示す説明図。
【図１０】粒子の摩擦抵抗を測定する様子を示す説明図。
【図１１】粒子の摩擦抵抗を測定した結果を示す説明図。
【図１２】粉体の流動性の悪化を検証する検証実験の様子を示す説明図。
【図１３】容器を振った回数と粒子の電圧値との関係を示す説明図。
【図１４】粒子間の付着力を示す説明図。
【図１５】圧粉体の部位毎の成形密度を予測した結果および実機にて測定した結果を示す説明図。
【図１６】焼結体の製造工程を示す説明図。（ａ）充填工程を示す図。（ｂ）成形工程を示す図。
【発明を実施するための形態】
【００１７】
以下では、第一実施形態の粉体の挙動予測方法について説明する。
【００１８】
図１に示すように、第一実施形態の粉体の挙動予測方法は、焼結体を製造する工程において、多数の粒子２からなる粉体１の挙動を予測するものである。
すなわち、第一実施形態において、粉体１は、一種類以上の金属粉末である。
【００１９】
まず、焼結体を製造する工程について説明する。
焼結体の製造工程では、充填工程および成形工程等を行う。
【００２０】
図１６（ａ）に示すように、充填工程は、下パンチ１０２を底面とするダイ１０１の内部空間に粉体１を充填する工程である。
ダイ１０１の内側面および下パンチ１０２の上面は、焼結体の形状に対応する形状に形成される。下パンチ１０２は、ダイ１０１の内側に配置され、上面がダイ１０１の上面よりも低い位置に位置する。
【００２１】
充填工程では、内部に粉体１が貯溜されるとともに下側が開口する粉箱１０３がダイ１０１の上面に配置され、ダイ１０１の上面に配置される粉箱１０３を水平方向に移動させる（図１６（ａ）に示す矢印Ｍ３参照）。このとき、粉箱１０３は、ダイ１０１の内部空間を繰り返し通過する。
これにより、充填工程では、粉箱１０３内の粉体１をダイ１０１の内部空間に落下させ、ダイ１０１の内部空間に粉体１を充填する。
粉体１は、ホッパー等の貯溜部からシュートを介して搬送され、粉箱１０３の内部に供給される。
【００２２】
図１６（ｂ）に示すように、成形工程は、ダイ１０１の内部空間に充填される粉体１に対して上パンチ１０４によって荷重をかける工程である。
上パンチ１０４の下面は、焼結体の形状に対応する形状に形成される。ダイ１０１および上下パンチ１０２・１０４は、シリンダのロッド等に連結され、それぞれ独立して上下方向に移動可能に構成されている。
【００２３】
成形工程では、ダイ１０１および上下パンチ１０２・１０４をそれぞれ下方向に移動させる。このとき、成形工程では、ダイ１０１および下パンチ１０２よりも速い速度で上パンチ１０４を下方向に移動させ、上パンチ１０４によって粉体１を押圧する（図１６（ｂ）に示す矢印Ｍ１・Ｍ２・Ｍ４参照）。
これにより、成形工程では、粉体１に対して荷重をかけて粉体１を加圧成形し、圧粉体Ｓ（図６参照）を成形する。
【００２４】
焼結体の製造工程では、成形工程を行った後で圧粉体Ｓをダイ１０１の内部空間から取り出して、圧粉体Ｓを焼き固めて焼結体を製造する。
【００２５】
このようにして製造される焼結体は、圧粉体Ｓの密度分布がばらついた場合において、密度が低い部分の強度が弱くなってしまう。つまり、この場合、焼結体は、部分的に強度が弱くなってしまう。
従って、焼結体の製造工程においては、焼結体の強度、つまり品質を向上させるために、密度分布のばらつきが小さい圧粉体Ｓを成形する必要がある。
【００２６】
圧粉体Ｓの密度分布は、充填工程および成形工程における粉体１の挙動が変わることで変動する。粉体１の挙動は、焼結体の製造条件を変更することで変わる。
すなわち、圧粉体Ｓの密度分布は、例えば、粉箱１０３の移動速度を速くすることでばらつきが小さくなったり、下パンチ１０２の移動速度を遅くすることでばらつきが大きくなったりする。
【００２７】
そこで、焼結体の製造工程においては、実際に焼結体を製造する前に粉体１の挙動を予測することで、製造条件を予め評価することが重要である。
【００２８】
図１に示すように、第一実施形態の粉体の挙動予測方法は、このような製造条件の評価のために行われ、充填工程および成形工程での粉体１の挙動を予測する。
【００２９】
次に、第一実施形態の粉体の挙動予測方法の手順について説明する。
【００３０】
粉体の挙動予測方法は、パーソナルコンピュータによって構成されるＣＡＥ（Ｃｏｍｐｕｔｅｒ  Ａｉｄｅｄ  Ｅｄｕｃａｔｉｏｎ）を用いて行われる。
【００３１】
まず、粉体の挙動予測方法では、焼結体の製造条件、つまり、粒子情報、型３Ｄモデル、充填条件、および成形条件をＣＡＥに入力する。
【００３２】
粒子情報は、例えば、粒子２の形状、粒子２の径ｄ（図２参照）、粉体１の種類（粉種）、および粒子２の質量ｍ等である。
【００３３】
なお、本実施形態において、粒子２の形状は、球状であるものとする。仮に、粒子が球状でない場合、粉体の挙動予測方法では、粒子の形状に対応するように、複数の球体を組み合わせた球体群を、粒子の形状としてＣＡＥに入力する。
【００３４】
型３Ｄモデルは、ダイ１０１および上下パンチ１０２・１０４の形状を３Ｄモデル化したものである。つまり、型３Ｄモデルは、例えば、ダイ１０１および上下パンチ１０２・１０４の形状に基づいて作成されたＣＡＤデータ等である。
【００３５】
充填条件は、例えば、粉箱１０３の形状を３Ｄモデル化したＣＡＤデータ、粉箱１０３の移動速度等である。
【００３６】
成形条件は、例えば、ダイ１０１および上下パンチ１０２・１０４の移動速度、成形工程時に粉体１にかける荷重の大きさ等である。
【００３７】
粒子情報等をＣＡＥに入力した後で、粉体の挙動予測方法では、型３Ｄモデル等に基づいてダイ１０１および粉箱１０３等に対応する領域を仮想空間上に形成するとともに、粉箱１０３に対応する領域内に粒子２を配置する。
そして、粉体の挙動予測方法では、離散要素法（個別要素法）を用いて粉体１の挙動を予測する。
【００３８】
離散要素法は、個々の粒子２に対して以下に示す運動方程式を解いて粉体１の挙動を予測するものである。
Ｆ＝ｍａ
ここで、前記運動方程式において、Ｆは個々の粒子２に対して作用する外力、ｍは挙動を予測する粒子２の質量、ａは挙動を予測する粒子２の加速度である。
【００３９】
粉体の挙動予測方法では、個々の粒子２に対して作用する外力Ｆを算出する。このような外力Ｆとしては、接触力、重力、および空気抵抗力等がある。
【００４０】
接触力は、挙動を予測する粒子２が他の粒子２と接触することで発生する力、および挙動を予測する粒子２がダイ１０１、上下パンチ１０２・１０４、および粉箱１０３と接触することで発生する力である。
すなわち、接触力は、粒子２の反発力Ｆ1（図３参照）および粒子２が擦れて移動するときの摩擦力（図８に示す矢印Ｆ2参照）等である。
【００４１】
粉体の挙動予測方法では、以下のようにして接触力を算出する。
すなわち、粉体の挙動予測方法では、仮想空間に配置した粒子２の位置（仮想空間上における座標情報）、およびＣＡＥに入力された粒子情報等に基づいて、挙動を予測する粒子２が他の粒子２と接触しているかどうかを判定する。
【００４２】
図２に示すように、粉体の挙動予測方法では、仮想空間にて二つの粒子２をオーバーラップさせて各粒子２が接触している状態を再現する。
このため、粉体の挙動予測方法では、例えば、挙動を予測する粒子２と他の粒子２との中心間の距離ｒが、前記各粒子２の半径の合計値よりも小さい場合に、挙動を予測する粒子２が他の粒子２と接触していると判定する。
【００４３】
仮に、挙動を予測する粒子２が他の粒子２と接触している場合、粉体の挙動予測方法では、接触力を算出する。
このとき、粉体の挙動予測方法では、挙動を予測する粒子２と他の粒子２との半径の合計値から各粒子２の中心間の距離ｒを減算して粒子２の変形量δを算出する。そして、粉体の挙動予測方法では、粒子２の変形量δ、反発係数η、および摩擦係数μ等に基づいて、粒子２の反発力Ｆ1および摩擦力を算出する。
【００４４】
また、粉体の挙動予測方法では、このような粒子間の接触力を算出する場合と同じ要領で、粒子２がダイ１０１、上下パンチ１０２・１０４、および粉箱１０３に対して接触しているかどうかも判定する。
すなわち、粉体の挙動予測方法では、粒子２がダイ１０１等に対応する領域に対してオーバーラップしているかどうかを確認し、粒子２がダイ１０１等に対して接触している場合に、粒子２とダイ１０１等との接触力（粒子２の反発力Ｆ1および摩擦力等）を算出する。
この場合、粉体の挙動予測方法では、粒子２の半径から、粒子２の中心からダイ１０１等に対応する領域までの距離を減算して粒子２の変形量δを算出する。
【００４５】
接触力を算出した後で、粉体の挙動予測方法では、接触の有無に関わらず粒子２に働く力、つまり、重力、および空気抵抗力等を算出する。
【００４６】
粉体の挙動予測方法では、接触力、重力、および空気抵抗力の算出結果を合算して、外力Ｆを算出する。
また、成形工程における粉体１の挙動を予測するとき、粉体の挙動予測方法では、前記接触力等と成形工程時に粉体１にかける荷重とを合算した結果を外力Ｆとして算出する。
【００４７】
外力Ｆを算出した後で、粉体の挙動予測方法では、外力Ｆの算出結果およびＣＡＥに入力された粒子２の質量ｍを運動方程式に代入して、加速度ａを算出する。
そして、粉体の挙動予測方法では、所定の時間増分に対して加速度ａを数値積分することで挙動を予測する粒子２の速度、変位、および位置を算出する。
【００４８】
粉体の挙動予測方法では、このような計算を解析対象となる全ての粒子２に対して行なうことによって、所定時間経過後の全ての粒子２の速度、変位、および位置を算出する。
粉体の挙動予測方法では、算出した全ての粒子２の速度、変位、および位置に基づいて、仮想空間上に配置した全ての粒子２を動かす。
【００４９】
粉体の挙動予測方法では、一定の時間が経過するまでこのような処理を繰り返し行うことにより、成形工程終了時の、つまり、圧粉体成形時の粒子２の位置を算出する。
粉体の挙動予測方法では、このようにして充填工程および成形工程時の粉体１の挙動を予測する。
【００５０】
このように、粉体の挙動予測方法は、外力Ｆを基に立てた運動方程式を解くことで、個々の粒子２の速度、変位、および位置を算出するステップを繰り返し行う。
【００５１】
ここで、成形工程では、上パンチ１０４によって粉体１に対して荷重をかけ、粉体１を加圧成形する（図１６（ｂ）参照）。粒子２は、このときの荷重によって大きく変形する。
【００５２】
離散要素法では、一般的に粒子２が弾性変形していることを想定してこのような場合における粒子２の反発力Ｆ1を算出する。
より詳細には、離散要素法では、一般的に粉体１の種類毎に予め設定される固定の反発係数ηおよび粒子２の変形量δを積算して、粒子２の反発力Ｆ1を算出する。
【００５３】
成形工程では、ある程度大きな荷重をかけて圧粉体Ｓを成形する。このため、図３に示すように、粒子２の多くは、成形工程時に降伏して塑性変形する。
この場合、粒子２の反発力Ｆ1は、粒子２が塑性変形した分だけ小さくなる。つまり、成形工程において、粒子２の反発力Ｆ1は、粒子２の変形量δに比例して（粒子２の変形量δに対して直線的に）発生するものではない。
【００５４】
従って、粉体１の種類毎に固定の反発係数ηを設定するだけでは、粒子２の反発力Ｆ1を正確に算出できないため、粉体１の挙動を正確に再現できない可能性がある（図５に二点鎖線で示す従来のグラフ参照）。
【００５５】
そこで、第一実施形態の粉体の挙動予測方法では、成形工程時に粒子２が塑性変形することを想定して、粒子２の反発力Ｆ1を算出する。
【００５６】
具体的には、第一実施形態の粉体の挙動予測方法では、前記時間増分量から成形工程が始まっているかどうかを判断する。
仮に成形工程が始まっている場合、第一実施形態の粉体の挙動予測方法では、以下の数式１を用いて粒子２が塑性変形したときの反発係数ηを算出する。
【００５７】
【数１】
2014154045.tif 000007
【００５８】
前記数式１において、ηは反発係数であり、ａ1〜ａ3は粉体１の種類毎に設定される固定値であり、δは粒子２の変形量であり、ｄは粒子２の径である。
【００５９】
第一実施形態の粉体の挙動予測方法では、実際に粉体１の挙動を予測する前に、以下のようにして固定値ａ1〜ａ3を具体的な数値を予め取得する。
【００６０】
すなわち、図４に示すように、第一実施形態の粉体の挙動予測方法では、加圧ヘッド１１によって粒子２を押圧し、粒子２を塑性変形させる。そして、第一実施形態の粉体の挙動予測方法では、前記押込み荷重を除荷した後の粒子２の径、つまり、塑性変形後の粒子２の径を測定する。
第一実施形態の粉体の挙動予測方法では、塑性変形前後の粒子２の径に基づいて、粒子２のひずみを算出する。
第一実施形態の粉体の挙動予測方法では、加圧ヘッド１１によって押込み荷重を徐々に大きくしながら、このようなひずみの算出を繰り返し行う。
【００６１】
これにより、図５に示すように、第一実施形態の粉体の挙動予測方法では、ひずみおよび押込み荷重の関係を取得する（図５に示す実測のグラフ参照）。
【００６２】
第一実施形態の粉体の挙動予測方法では、取得したひずみと押込み荷重との関係をｎ乗硬化則を用いて近似して、前記数式１の固定値ａ1〜ａ3の具体的な数値を取得する（図５に実線で示す近似結果のグラフ参照）。
また、第一実施形態の粉体の挙動予測方法では、焼結体の製造に用いる粉体１の種類毎にこのようなひずみおよび押込み荷重の関係を取得することで、粉体１の種類毎に前記数式１の固定値ａ1〜ａ3の具体的な数値を予め取得する。
【００６３】
粉体の挙動予測方法では、前記数式１に予め取得している固定値ａ1〜ａ3、粒子２の変形量δの算出結果、および粒子２の径ｄを代入して反発係数ηを算出する。
そして、粉体の挙動予測方法では、反発係数ηと粒子２の変形量δとの算出結果を積算して、粒子２の反発力Ｆ1を算出する。
【００６４】
前述のように、粒子２が塑性変形する場合には、仮想空間において、挙動を予測する粒子２と他の粒子２等とのオーバーラップ量が増大し、前記数式１に代入する粒子２の変形量δが増大する（図２参照）。
第一実施形態の粉体の挙動予測方法は、このような粒子２の変形量δに応じて反発係数ηを算出できるため、成形工程時に粒子２が塑性変形したときの反発係数ηの変動を考慮できる。つまり、図５に示すように、第一実施形態の粉体の挙動予測方法は、粒子２の反発力Ｆ1を正確に算出できる。
【００６５】
すなわち、第一実施形態の粉体の挙動予測方法は、実際に焼結体を製造するときに発生する、粒子２の塑性変形に起因する粒子２の反発係数ηの変動を把握した上で、粒子２の挙動を再現できる。
【００６６】
従って、第一実施形態の粉体の挙動予測方法は、粉体１の挙動の予測精度を向上できる。
【００６７】
このように、粉体の挙動予測方法は、粒子２が塑性変形したときの反発係数ηを算出し、反発係数ηの算出結果に基づいて粒子２の反発力Ｆ1を算出するステップを行う。
【００６８】
図１および図６に示すように、粒子２の挙動を予測した後で、粉体の挙動予測方法では、板面における開口部分Ｓ１や足部分Ｓ５といった圧粉体Ｓの各部位Ｓ１〜Ｓ５に対応する領域内に、何個の粒子２が位置しているかをカウントする。
そして、粉体の挙動予測方法では、圧粉体Ｓの各部位Ｓ１〜Ｓ５毎に粒子２の体積含有率、つまり、単位体積当たりに占める粒子２の個数を算出する。
【００６９】
図１および図７に示すように、粒子２の体積含有率を算出した後で、粉体の挙動予測方法では、粒子２の体積含有率の算出結果および粒子２の体積に基づいて、圧粉体Ｓの密度分布を算出する。
【００７０】
圧粉体Ｓの密度分布を算出した後で、粉体の挙動予測方法では、密度分布の算出結果に基づいてＣＡＥに入力された充填条件等の製造条件を評価する。
すなわち、粉体の挙動予測方法では、例えば、密度分布のばらつき（密度が最も高い部位と密度が最も低い部位との差）が、許容されるばらつきの範囲を超えた場合等において、入力された製造条件を微調整しながら離散要素法を用いて粉体１の挙動の予測を行う。
そして、粉体の挙動予測方法では、許容できる程度に密度分布のばらつきが小さくなったときの製造条件を、パーソナルコンピュータのディスプレイ等に出力する。
【００７１】
第一実施形態の粉体の挙動予測方法は、粒子２の塑性変形を考慮して粒子２の反発力Ｆ1を算出することで、高い精度で粉体１の挙動を予測できるため、より高い精度で（実用化できる程度に）製造条件を評価できる。
つまり、第一実施形態の粉体の挙動予測方法は、実際に焼結体を製造する前に粉体１の挙動を把握できるため、焼結体の品質を向上できる。
【００７２】
次に、第二実施形態の粉体の挙動予測方法について説明する。
【００７３】
第二実施形態の粉体の挙動予測方法は、粉体１の挙動を予測するときに、粒子２の塑性変形とは異なる現象を考慮する点が、第一実施形態の粉体の挙動予測方法と比較して異なる点である。
このため、以下では、第二実施形態で考慮する現象について詳細に説明し、他の説明を省略する。
【００７４】
図８に示すように、粉体の挙動予測方法では、外力Ｆを算出するときに、接触力として、粒子２が擦れて移動するときの摩擦力Ｆ2を算出する。
このとき、第二実施形態の粉体の挙動予測方法では、摩擦係数μ、粉体１の種類毎に設定される固定値、および粒子２の変形量δを積算して、粒子２が擦れて移動するときの摩擦力Ｆ2を算出する。
【００７５】
成形工程においては、上パンチ１０４によって粉体１に対して荷重がかかり、粒子２に大きな応力が発生している（図１６（ｂ）参照）。すなわち、成形工程において、粒子２は、大きな応力が発生した状態で、接触する他の粒子２あるいは上下パンチ１０２・１０４等に対して、擦れて移動する。
第二実施形態の粉体の挙動予測方法では、このような粒子２が擦れて移動するときの摩擦係数μに着目し、図９に示すような測定装置２０を用いて応力と摩擦係数μとの関係を確認した。
【００７６】
まず、測定装置２０の構成について説明する。
図９に示すように、測定装置２０は、上下の貯溜部材２１・２２、サーボシリンダ２３、リニアアクチュエータ２４、および二つのロードセル２５・２６を具備する。
【００７７】
上側の貯溜部材２１は、上下両側が上下方向に沿って貫通し、内部に粉体１が貯溜される。
下側の貯溜部材２２は、例えば、略直方体状の部材の上面の一部を切り欠いたような形状に形成され、前記切欠部分に粉体１が貯溜される。
上下の貯溜部材２１・２２は、水平方向における位置が合わされて、貯溜している粉体１が互いに接触する。
【００７８】
サーボシリンダ２３は、上下方向に移動可能なロッドの先端部に押圧ヘッド２３ａが取り付けられる。押圧ヘッド２３ａは、上側の貯溜部材２１の内側の形状と略同一の形状に形成され、上側の貯溜部材２１の内側に配置される。押圧ヘッド２３ａの下面は、粉体１に接触する。
サーボシリンダ２３は、ロッドを下方向に移動させることにより、上側の貯溜部材２１に貯溜される粉体１を下方向に押圧する。
【００７９】
リニアアクチュエータ２４は、水平方向に移動可能なロッドが下側の貯溜部材２２に連結される。リニアアクチュエータ２４は、ロッドを水平方向に移動させることにより、下側の貯溜部材２２を水平方向に移動させる。
【００８０】
上側のロードセル２５は、上側の貯溜部材２１の上方に配置され、サーボシリンダ２３による下方向への荷重を測定する。
下側のロードセル２６は、下側の貯溜部材２２の右方に配置され、リニアアクチュエータ２４による水平方向への荷重を測定する。
【００８１】
次に、応力と摩擦係数μとの関係を取得する手順について説明する。
【００８２】
図１０に示すように、第二実施形態の粉体の挙動予測方法では、サーボシリンダ２３を駆動させて上側の貯溜部材２１に貯溜される粉体１を下方向に押圧するとともに、リニアアクチュエータ２４を駆動させて下側の貯溜部材２２を水平方向にスライドさせる（図１０に示す矢印参照）。
各ロードセル２５・２６は、このように下側の貯溜部材２２がスライドしたときのサーボシリンダ２３およびリニアアクチュエータ２４の荷重を測定する。
【００８３】
これにより、第二実施形態の粉体の挙動予測方法では、垂直応力と剪断応力との関係を取得する（図１０の下側に示すグラフ参照）。
【００８４】
剪断応力は、垂直応力および摩擦係数μの積によって求められる。従って、図１０に示すような縦軸に剪断応力をとり、横軸に垂直応力をとったグラフの傾きは、摩擦係数μに相当する。
【００８５】
そこで、第二実施形態の粉体の挙動予測方法では、一定のピッチ毎に前記グラフを区切って、前記グラフの傾きを算出する。
これにより、図１１に示すように、第二実施形態の粉体の挙動予測方法では、垂直応力に応じた摩擦係数μを算出する。
第二実施形態の粉体の挙動予測方法では、このような摩擦係数μの算出を、焼結体の製造に用いる粉体１の種類毎に行う（図１１に示す粉種Ｐ１〜Ｐ３参照）。
【００８６】
図１１からも明らかなように、摩擦係数μは、垂直応力および粉体１の種類に応じて異なる値となっている。つまり、成形工程時の粉体１の挙動を予測する場合には、粒子２に発生した応力に応じて最適な摩擦係数μを算出することで、粒子２が擦れて移動するときの摩擦力Ｆ2をより正確に算出できる。
【００８７】
そこで、第二実施形態の粉体の挙動予測方法では、図１１に示す摩擦係数μの測定結果から導き出した以下に示す数式２を用いて、粒子２に発生した応力に応じた摩擦係数μを算出する。
【００８８】
【数２】
2014154045.tif 000008
【００８９】
前記数式２において、μは摩擦係数であり、ｂ1〜ｂ3は粉体１の種類毎に設定される固定値であり、δは粒子２の変形量であり、ｄは粒子２の径である。
【００９０】
粉体の挙動予測方法では、実際に粉体１の挙動を予測する前に、測定装置２０を用いて図１０に示すような垂直応力および摩擦係数μの関係を取得することで、前記数式２の固定値ｂ1〜ｂ3の具体的な数値を予め取得する。
【００９１】
粉体の挙動予測方法では、前記数式２に固定値ｂ1〜ｂ3、粒子２の変形量δ、および粒子２の径ｄを代入して摩擦係数μを算出する。
粉体の挙動予測方法では、摩擦係数μの算出結果と粉体１の種類毎に設定される固定値と粒子２の変形量δの算出結果とを積算して、粒子２の摩擦力Ｆ2を算出する。
なお、粒子２の摩擦力Ｆ2を算出するときに用いる固定値は、予め実験等によって算出されている。
【００９２】
前述のように、粒子２に大きな応力が発生した場合には、仮想空間において、挙動を予測する粒子２と他の粒子２等とのオーバーラップ量が増大し、前記数式２に代入する粒子２の変形量δが増大する（図２参照）。
第二実施形態の粉体の挙動予測方法は、このような粒子２の変形量δに応じて摩擦係数μを算出できるため、成形工程時に粒子２に大きな応力が発生した場合の摩擦係数μの変動を考慮できる。つまり、図１１に示すように、第二実施形態の粉体の挙動予測方法は、粉体１の種類毎に固定の摩擦係数μを用いて粒子２の摩擦力Ｆ2を算出する場合と比較して、粒子２の摩擦力Ｆ2を正確に算出できる。
【００９３】
すなわち、第二実施形態の粉体の挙動予測方法は、実際に焼結体を製造するときに発生する、粒子２が擦れて移動するときの摩擦係数μの変動を把握した上で、粒子２の挙動を再現できる。
【００９４】
従って、第二実施形態の粉体の挙動予測方法は、粉体１の挙動の予測精度を向上できる。
【００９５】
このように、第二実施形態の粉体の挙動予測方法は、粒子２に発生した応力に応じて摩擦係数μを算出し、摩擦係数μの算出結果に基づいて粒子２が擦れて移動するときの摩擦力Ｆ2を算出するステップを行う。
【００９６】
次に、第三実施形態の粉体の挙動予測方法について説明する。
【００９７】
第三実施形態の粉体の挙動予測方法は、粉体１の挙動を予測するときに、粒子２の塑性変形とは異なる現象を考慮する点が、第一実施形態の粉体の挙動予測方法と比較して異なる点である。
このため、以下では、第三実施形態で考慮する現象について詳細に説明し、他の説明を省略する。
【００９８】
第三実施形態の粉体の挙動予測方法では、粉体１の流動性に着目し、以下のような検証実験を行った。
【００９９】
図１２（ａ）に示すように、検証実験では、所定量の粉体１を入れて密閉した容器３１をモータの駆動軸等に連結し、水平方向を回転軸として容器３１を所定の角度の範囲で一定回数回転させた（図１２（ａ）に二点鎖線および点線で示す矢印参照）。つまり、検証実験では、一定回数容器３１を振った。
そして、図１２（ｂ）に示すように、検証実験では、蓋体３２を容器３１から取り外して粉体１を容器３１から漏斗３３に移し、粉体１が詰まることなく漏斗３３を通過するかどうかを確認した（図１２（ｂ）に矢印で示す粉体１参照）。
【０１００】
検証実験では、容器３１を振る回数を変更しながら、このような確認を繰り返し行った。また、検証実験では、容器３１を振らずに粉体１を漏斗３３に移し、粉体１が詰まることなく漏斗３３を通過するかどうかを確認した。
【０１０１】
このような検証実験の結果、容器３１を振って漏斗３３に移した粉体１は、漏斗３３に詰まってしまった。
一方、容器３１を振らずに漏斗３３に移した粉体１は、詰まることなく漏斗３３を通過した。
【０１０２】
以上より、容器３１を振ることによって粉体１の流動性が悪化していることがわかる。
【０１０３】
次に、検証実験では、容器３１を一定回数振った後で粉体１を容器３１から速やかに取り出して、帯電量測定装置によって粉体１の電圧値を測定した。また、検証実験では、容器３１を振らずに容器３１を粉体１から取り出して、粉体１の電圧値を測定した。
【０１０４】
図１３に示すように、容器３１を振らなかった場合の粉体１の電圧値は、ややプラスの電圧値であった（図１３に示す電圧値＋Ｖ1参照）。
一方、容器３１を一定回数振った場合の粉体１の電圧値は、マイナスの電圧値であった（図１３に示す電圧値−Ｖ2〜−Ｖ4参照）。また、容器３１を振る回数が多くした場合の粉体１の電圧値は、容器３１を振る回数が少ない場合と比較して、より大きくなった、つまり、０Ｖより大きく離れた値となった（図１３に示す振り回数Ｎ3および電圧値−Ｖ4参照）。
【０１０５】
これは、容器３１を振ることで粒子同士が衝突し、粒子２に静電気力が発生していることによるものであると考えられる。
【０１０６】
このような静電気力が発生している場合、一定時間放置することで放電現象が起こるため、粒子２の電圧値は、容器３１を振る前の電圧値に戻る。
【０１０７】
そこで、検証実験では、容器３１を一定回数振った後で容器３１を一定時間放置した。その後、検証実験では、蓋体３２を容器３１から取り外して漏斗３３に粉体１を移し、粉体１が詰まることなく漏斗３３を通過するかどうかを確認した（図１２参照）。
【０１０８】
この場合、粉体１は、詰まることなく漏斗３３を通過した。
【０１０９】
以上の検証実験より、図１４に示すように、粉同士が衝突することで粒子２に静電気力が発生し、粒子間に付着力Ｆ3、より詳細には、互いに離間する二つの粒子間に引力が発生して粉体１の流動性が悪化することがわかった。
【０１１０】
このような粉同士の衝突による粉体１の流動性の悪化は、粉体１を搬送するとき、つまり、充填工程において頻繁に発生することが想定される。
【０１１１】
そこで、第三実施形態の粉体の挙動予測方法では、充填工程時に挙動を予測する粒子２と他の粒子２との間に付着力Ｆ3が発生することを考慮して、外力Ｆを算出する。
【０１１２】
具体的には、第三実施形態の粉体の挙動予測方法では、前記時間増分量から充填工程が始まっているかどうかを判断する。
仮に充填工程が始まっている場合、第三実施形態の粉体の挙動予測方法では、挙動を予測する粒子２と他の粒子２とが接触しているかどうかを判定する。
そして、第三実施形態の粉体の挙動予測方法では、挙動を予測する粒子２と他の粒子２とが接触していない場合（各粒子２が離間している場合）に、各粒子２の間で発生する付着力Ｆ3を算出する。
【０１１３】
このとき、第三実施形態の粉体の挙動予測方法では、まず、以下の数式３を用いて、挙動を予測する粒子２と、挙動を予測する粒子２に対して離間する他の粒子２（互いに離間する二つの粒子２）の電荷Ｑ1・Ｑ2をそれぞれ算出する。
【０１１４】
【数３】
2014154045.tif 000009
【０１１５】
前記数式３において、Ｑnは粒子２の電荷、Ｑ0は元の電荷（つまり、粒子２が元々持っている電荷）、ｃは粉体１の種類毎に設定される固定値、Ｎは粒子同士の衝突回数、ｄは粒子２の径である。
【０１１６】
第三実施形態の粉体の挙動予測方法では、例えば、粉体１を一定期間静置したときの粒子２の電荷を、帯電量測定装置等を用いて測定することで、実際に粉体１の挙動を予測する前に、元の電荷Ｑ0の具体的な数値を取得する。
第三実施形態の粉体の挙動予測方法では、焼結体の製造に用いる粉体１の種類毎にこのような元の電荷Ｑ0の具体的な数値を取得する。
【０１１７】
また、第三実施形態の粉体の挙動予測方法では、粒子同士を衝突させたときの電荷を測定する実験を行うことで、実際に粉体１の挙動を予測する前に、固定値ｃを算出する。
【０１１８】
そして、第三実施形態の粉体の挙動予測方法では、粒子２が他の粒子２と接触した回数をカウントすることで、衝突回数Ｎを算出する。
【０１１９】
第三実施形態の粉体の挙動予測方法では、元の電荷Ｑ0、固定値ｃ、衝突回数Ｎの算出結果、および粒子２の径ｄを前記数式３に代入することで、挙動を予測する粒子２と、挙動を予測する粒子２に対して離間する他の粒子２との電荷Ｑ1・Ｑ2を算出する。
そして、第三実施形態の粉体の挙動予測方法では、以下の数式４を用いて粒子間の付着力Ｆ3を算出する。
【０１２０】
【数４】
2014154045.tif 000010
【０１２１】
前記数式４において、Ｆ3は粒子間の付着力、ε0は真空の誘電率、εrは比誘電率、Ｑ1・Ｑ2は各粒子２の電荷、ｒは各粒子２の中心間の距離である。
【０１２２】
第三実施形態の粉体の挙動予測方法では、固定値である真空の誘電率ε0、粉体１の種類毎に設定される比誘電率εr、各粒子２の電荷Ｑ1・Ｑ2の算出結果、各粒子２の中心間の距離ｒを前記数式４に代入することで、粒子間の付着力Ｆ3を算出する。
【０１２３】
第三実施形態の粉体の挙動予測方法では、粒子間の付着力Ｆ3の算出結果、粒子間の接触力、重力、および空気抵抗力等を合算して、外力Ｆを算出する。
そして、第三実施形態の粉体の挙動予測方法では、第一実施形態の粉体の挙動予測方法と同じ要領で粉体１の挙動を予測する。
【０１２４】
これによれば、第三実施形態の粉体の挙動予測方法は、粒子２の衝突によって発生する粉体１の流動性の悪化、つまり、粒子間の付着力Ｆ3の影響を考慮して、外力Ｆを算出できる。
すなわち、第三実施形態の粉体の挙動予測方法は、実際に焼結体を製造するときに発生する付着力Ｆ3の影響を把握した上で、粒子２の挙動を再現できる。
【０１２５】
従って、第三実施形態の粉体の挙動予測方法は、粉体１の挙動の予測精度を向上できる。
【０１２６】
このように、第三実施形態の粉体の挙動予測方法は、粒子同士の衝突によって発生する静電気力に基づいて、互いに離間する二つの粒子間の付着力Ｆ3を算出するステップを行う。
【０１２７】
なお、本発明に係る粉体の挙動予測方法は、第一実施形態にあるような反発力Ｆ1、第二実施形態にあるような摩擦力Ｆ2、および第三実施形態にあるような付着力Ｆ3の少なくともいずれか一つを算出すればよい。
つまり、本発明に係る粉体の挙動予測方法は、反発力Ｆ1および摩擦力Ｆ2を算出する、反発力Ｆ1および付着力Ｆ3を算出するといった具合に、反発力Ｆ1、摩擦力Ｆ2、および付着力Ｆ3を二つ以上算出しても構わない。
【０１２８】
図１５に示す本発明のグラフは、第一実施形態にあるような反発力Ｆ1、第二実施形態にあるような摩擦力Ｆ2、および第三実施形態にあるような付着力Ｆ3を全て算出して粉体１の挙動を予測して、図６に示す圧粉体Ｓの各部位Ｓ１〜Ｓ５の密度分布を算出したものである。
【０１２９】
図１５に示す従来のグラフは、第一実施形態にあるような反発力Ｆ1、第二実施形態にあるような摩擦力Ｆ2、および第三実施形態にあるような付着力Ｆ3を全て算出せずに粉体１の挙動を予測して、図６に示す圧粉体Ｓの各部位Ｓ１〜Ｓ５の密度分布を算出したものである。
すなわち、図１５に示す従来のグラフは、粉体１の種類毎に設定される固定の反発係数ηおよび摩擦係数μに基づいて粒子２の反発力Ｆ1および摩擦力Ｆ2を算出するとともに、粒子間の付着力Ｆ3を算出せずに粉体１の挙動を予測して、各部位Ｓ１〜Ｓ５の密度分布を算出したものである。
【０１３０】
図１５に示す実測のグラフは、ダイ１０１および上下パンチ１０２・１０４等を用いて実際に成形した圧粉体Ｓを部位Ｓ１〜Ｓ５毎に切断し、水に浮かべて密度を測定することで成形密度（密度分布）を測定したものである。
【０１３１】
このような粉体１の挙動の予測および圧粉体Ｓの成形は、それぞれ同じ製造条件で行ったものである。
【０１３２】
図１５に示すように、本発明のグラフは、部位Ｓ４を除いて従来のグラフよりも実測のグラフに近付けることができた。また、本発明のグラフは、部位Ｓ４においても、実測のグラフに対して比較的近い成形密度となった。
つまり、本発明に係る粉体の挙動予測方法は、第一実施形態にあるような反発力Ｆ1、第二実施形態にあるような摩擦力Ｆ2、および第三実施形態にあるような付着力Ｆ3を考慮することで、実測に対する誤差を全体的に低減できる。
【０１３３】
これは、予測時に実際に起こる現象を把握した上で、粉体１の挙動を再現できることによるものである。
従って、本発明に係る粉体の挙動予測方法では、予測精度を向上できる。
【０１３４】
このように、粉体の挙動予測方法は、粒子２が塑性変形したときの反発係数ηに基づいて粒子２の反発力Ｆ1を算出するステップと、応力に応じて算出した摩擦係数μに基づいて摩擦力Ｆ2を算出するステップと、粒子同士の衝突によって発生する静電気力に基づいて付着力Ｆ3を算出するステップと、の少なくともいずれか一つを行うことにより、個々の粒子２に対して作用する外力Ｆを算出する。
【符号の説明】
【０１３５】
１      粉体
２      粒子
Ｆ      外力
Ｆ1    反発力
Ｆ2    摩擦力
Ｆ3    付着力
η      反発係数
μ      摩擦係数
【図１】
2014154045.tif 000015
【図２】
2014154045.tif 000016
【図３】
2014154045.tif 000017
【図４】
2014154045.tif 000018
【図５】
2014154045.tif 000019
【図６】
2014154045.tif 000020
【図７】
2014154045.tif 000021
【図８】
2014154045.tif 000022
【図９】
2014154045.tif 000023
【図１０】
2014154045.tif 000024
【図１１】
2014154045.tif 000025
【図１２】
2014154045.tif 000026
【図１３】
2014154045.tif 000027
【図１４】
2014154045.tif 000028
【図１５】
2014154045.tif 000029
【図１６】
2014154045.tif 000030