JP3721765B2 - ３次元流れの有限要素解析法、解析装置、成形品を製造する製造方法および解析法のプログラムを記録した媒体 - Google Patents
３次元流れの有限要素解析法、解析装置、成形品を製造する製造方法および解析法のプログラムを記録した媒体 Download PDFInfo
- Publication number
- JP3721765B2 JP3721765B2 JP02789698A JP2789698A JP3721765B2 JP 3721765 B2 JP3721765 B2 JP 3721765B2 JP 02789698 A JP02789698 A JP 02789698A JP 2789698 A JP2789698 A JP 2789698A JP 3721765 B2 JP3721765 B2 JP 3721765B2
- Authority
- JP
- Japan
- Prior art keywords
- pressure
- flow
- analysis
- equation
- finite element
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
Images
Classifications
-
- B—PERFORMING OPERATIONS; TRANSPORTING
- B29—WORKING OF PLASTICS; WORKING OF SUBSTANCES IN A PLASTIC STATE IN GENERAL
- B29C—SHAPING OR JOINING OF PLASTICS; SHAPING OF MATERIAL IN A PLASTIC STATE, NOT OTHERWISE PROVIDED FOR; AFTER-TREATMENT OF THE SHAPED PRODUCTS, e.g. REPAIRING
- B29C45/00—Injection moulding, i.e. forcing the required volume of moulding material through a nozzle into a closed mould; Apparatus therefor
- B29C45/17—Component parts, details or accessories; Auxiliary operations
- B29C45/76—Measuring, controlling or regulating
- B29C45/766—Measuring, controlling or regulating the setting or resetting of moulding conditions, e.g. before starting a cycle
-
- B—PERFORMING OPERATIONS; TRANSPORTING
- B29—WORKING OF PLASTICS; WORKING OF SUBSTANCES IN A PLASTIC STATE IN GENERAL
- B29C—SHAPING OR JOINING OF PLASTICS; SHAPING OF MATERIAL IN A PLASTIC STATE, NOT OTHERWISE PROVIDED FOR; AFTER-TREATMENT OF THE SHAPED PRODUCTS, e.g. REPAIRING
- B29C45/00—Injection moulding, i.e. forcing the required volume of moulding material through a nozzle into a closed mould; Apparatus therefor
- B29C45/17—Component parts, details or accessories; Auxiliary operations
- B29C45/76—Measuring, controlling or regulating
- B29C45/7693—Measuring, controlling or regulating using rheological models of the material in the mould, e.g. finite elements method
-
- B—PERFORMING OPERATIONS; TRANSPORTING
- B29—WORKING OF PLASTICS; WORKING OF SUBSTANCES IN A PLASTIC STATE IN GENERAL
- B29C—SHAPING OR JOINING OF PLASTICS; SHAPING OF MATERIAL IN A PLASTIC STATE, NOT OTHERWISE PROVIDED FOR; AFTER-TREATMENT OF THE SHAPED PRODUCTS, e.g. REPAIRING
- B29C45/00—Injection moulding, i.e. forcing the required volume of moulding material through a nozzle into a closed mould; Apparatus therefor
- B29C45/17—Component parts, details or accessories; Auxiliary operations
- B29C45/76—Measuring, controlling or regulating
- B29C45/77—Measuring, controlling or regulating of velocity or pressure of moulding material
-
- B—PERFORMING OPERATIONS; TRANSPORTING
- B29—WORKING OF PLASTICS; WORKING OF SUBSTANCES IN A PLASTIC STATE IN GENERAL
- B29C—SHAPING OR JOINING OF PLASTICS; SHAPING OF MATERIAL IN A PLASTIC STATE, NOT OTHERWISE PROVIDED FOR; AFTER-TREATMENT OF THE SHAPED PRODUCTS, e.g. REPAIRING
- B29C45/00—Injection moulding, i.e. forcing the required volume of moulding material through a nozzle into a closed mould; Apparatus therefor
- B29C45/17—Component parts, details or accessories; Auxiliary operations
- B29C45/76—Measuring, controlling or regulating
- B29C45/78—Measuring, controlling or regulating of temperature
Landscapes
- Engineering & Computer Science (AREA)
- Manufacturing & Machinery (AREA)
- Mechanical Engineering (AREA)
- Aerodynamic Tests, Hydrodynamic Tests, Wind Tunnels, And Water Tanks (AREA)
- Injection Moulding Of Plastics Or The Like (AREA)
- Complex Calculations (AREA)
Description
【０００１】
【発明の属する技術分野】
金属の鋳造、プラスチック成形などの分野の、流動過程の有限要素解析方法および解析装置に関する。
【０００２】
【従来の技術】
従来、非圧縮性流体の流れの有限要素解析を行うばあい、安定な解をえるためには圧力と流速を求める格子点を食い違わせる必要があったこれは、Ｎａｖｉｅｒ−Ｓｔｏｋｅｓ式などの流れの式において圧力が１次微分、流速が２次微分であることに起因している。このため、有限要素法では流速の補間次数より圧力の補間次数を下げることが行われる。２次元での解析においても通常、たとえば４角形の有限要素を用いるばあい、４つの頂点で流速を求め、４角形の中心で圧力を求める必要がある。このとき、流速の補間次数は１次で圧力は０次となる。２次元解析ではＧａｌｅｒｋｉｎ法を用いて同一格子点で流速と圧力を求める方法(Computer Methods in Applied Mechanics and Engineering、 58, (1986)、 J. G. Rice and R. J. Schnipke、 An Equal-Order Velocity-Pressure Formulation that does not Exhibit Spurious Pressure Modes、 P.135-149)があるが、補間関数で重み付けを行うＧａｌｅｒｋｉｎ法は質量、運動量やエネルギーの保存を明示的に考慮した方法ではないため、安定な計算には細かな有限要素分割が必要であった。
【０００３】
３次元流れの解析法として流動コンダクタンスという量を仮定し、流速の変数を消去して未知変数を圧力のみの１つに減らして計算する方法が特開平８−９９４３１号公報、特開平９−１５０４４３号公報などに開示されているが、多孔質中の浸透流を表す近似的なダルシー流れの式にしか用いることができない。
【０００４】
これまで、３次元流れの有限要素解析で自由表面の移動を計算するばあい、流れに乗って移動する質量０のマーカー粒子を仮定して計算する方法が主流であった。また、薄板流れ（Ｈｅｌｅ−Ｓｈａｗ流れ）のように近似解析では、ＶＯＦ（Volume of Fluid）法のように格子点において流体が占める割合を計算し、この値から自由表面の位置を求めることができる。自由表面の変化が簡単なばあいは、変化した流体領域の形状に合わせて、有限要素分割をやり直したり、有限要素を追加または消去していく方法がある。
【０００５】
【発明が解決しようとする課題】
流速と圧力が求まる格子点が食い違っていると、結果を表示するばあい、このままでは、流速と圧力は異なったデータ系列となり、統一的に処理できないという問題があった。また、Ｇａｌｅｒｋｉｎ法を用いた場合は、前述のように定式化においては質量やエネルギーなどの保存則を考慮していないため、安定な解を求めるためには、有限要素分割を細かくする必要があり、結果的に計算時間の増大や解の精度の低下を招いていた。
【０００６】
また、自由表面の計算でマーカー粒子を用いるばあいは、複雑な流れのばあい、時間が経過すると粒子の分布に偏りが発生する。このため、最初から莫大な数のマーカー粒子を用いるか、適当な時間ごとにマーカー粒子の再配置を行う必要があった。これも、結果的に計算時間の増大を招く。解析領域の再有限要素分割や有限要素の追加、消去する方法は、流れが分岐したり、合流するばあいを計算するには、有限要素の形状がひずんだり、重なったりするのを避けるための多くの条件を考慮しておかねばらない。ＶＯＦ法に関しては薄板流れへの適用のみであり、３次元流れへの適用は見当たらない。
【０００７】
本発明は、このような問題点を解決するもので、３次元流れにおける流速と圧力を同一の格子点上で求め、自由表面の移動を容易に計算することができるようにすることにある。
【０００８】
【課題を解決するための手段】
本発明は、３次元流れ解析の有限要素解析法であって、立体型要素を用い、コントロールボリューム法に基づいて離散化した運動量の式と、圧力項を除いて半離散化した運動量の式を離散化した連続の式に代入したあと圧力項を離散化してえられる圧力の式と、離散化した運動量の式から導出する速度補正式とを用いて流速および圧力を求める反復計算を行い、流速と圧力を同一節点上で求めることを特徴とする解析方法である。
【０００９】
また、本発明により、３次元流れ解析の有限要素解析法に関し、立体要素を用い、コントロールボリューム法に基づいて離散化した運動量の式と、圧力項を除いて半離散化した運動量の式を離散化した連続の式に代入したあと圧力項を離散化してえられる圧力の式と、離散化した運動量の式から導出する速度補正式とを用いて流速および圧力を求める反復計算を行い、流速と圧力を同一節点上で求めることを特徴とする解析装置を提供する。
【００１０】
また、本発明は、コントロールボリュームでは流入流出する流束を見積もることにより保存則を考慮した解析となるため、本質的に計算が安定で、移流項を計算するのと同じ方法で自由表面の移動も容易に計算できる解析方法である。
【００１１】
また、本発明により、コントロールボリュームでは流入流出する流束を見積もることにより保存則を考慮した解析となるため、本質的に計算が安定で、移流項を計算するのと同じ方法で自由表面の移動も容易に計算できる解析装置を提供する。
【００１２】
また、本発明は、温度を計算するエネルギーの式と反応を計算する拡散の式を容易に組み込み、コントロールボリュームでは流入流出する流束を見積もることにより保存則を考慮した解析となるため、本質的に計算が安定で、移流項を計算するのと同じ方法で自由表面の移動も容易に計算できる解析方法である。
【００１３】
また、本発明により、温度を計算するエネルギーの式と反応を計算する拡散の式を容易に組み込み、コントロールボリュームでは流入流出する流束を見積もることにより保存則を考慮した解析となるため、本質的に計算が安定で、移流項を計算するのと同じ方法で自由表面の移動も容易に計算できる解析装置を提供する。
【００１４】
また、本発明により、えられた自由表面の移動状況、流速分布、圧力分布、温度分布および物質の濃度分布の結果をもとに成形条件を決定し、この決定した成形条件で成形品を製造する製造方法を提供する。
【００１５】
また、本発明により、３次元流れの有限要素解析法のプログラムを記録した媒体を提供する。
【００１６】
【発明の実施の形態】
以下に、添付図を参照しつつ、本発明による３次元流れの有限要素解析法とその解析装置の一例である、熱硬化性樹脂成形における樹脂流動過程の解析法および解析装置の好ましい実施の形態の例を詳細に説明する。
【００１７】
熱硬化性樹脂の自由表面が移動する３次元流れの式は以下で表される。
【００１８】
【数１】
【００１９】
【数２】
【００２０】
【数３】
【００２１】
【数４】
【００２２】
【数５】
【００２３】
【数６】
【００２４】
【数７】
【００２５】
ここでｘ、ｙ、ｚは３次元の空間座標であり、ｕ、ｖ、ｗはそれぞれ座標軸方向の樹脂の流速、ｐ、Ｔ、ｃはそれぞれ圧力、温度および熱硬化性樹脂の既反応基濃度（反応率）である。ｆは流体占有率で、ｆ＝０の位置には流体が存在しないこと、０＜ｆ＜１の位置は自由表面、ｆ＝１の位置は流体が満ちていることを示す。また、ρは密度、ηは粘度、Ｃvは比熱、ｋは熱伝導率、Ｄは拡散係数である。また、ｇx、ｇy、ｇzはそれぞれ重力加速度のｘ方向成分、ｙ方向成分、ｚ方向成分であり、Ｈは反応熱、Ｒは反応速度である。
【００２６】
解析領域を有限要素で分割し、有限要素の各節点に対して、それぞれが重ならないようにコントロールボリュームを定義する。式（１）を、コントロールボリュームで積分し離散化する。積分記号の添え字ＣＶはコントロールボリュームで積分することを表す。左辺第１項の時間微分は以下のように後退差分で表す。
【００２７】
【数８】
【００２８】
ここで、Ｎｉは節点ｉに対する補間関数、ｕiは節点ｉでのｘ方向の流速である。また、（ｎ）はｎ回目の時間ステップであることを示し、ｎ＋１ステップとｎステップの時間間隔がΔｔである。
【００２９】
式（１）の左辺第２、３、４項の移流項は非線型であり、単純に通常の補間関数を用いて積分しても解の空間的振動などの不具合が発生する。このため、本発明の一つである、流れの方向を考慮し、コントロールボリューム法にマッチした計算法を適用する。なお、この計算法は式（５）のＴ、式（６）のｃおよび式（７）の自由表面の移動に対しても同様に用いることができる。ｕ、ｖ、ｗ、Ｔ、ｃ、ｆをまとめてφで表せば、移流項のコントロールボリュームに対する積分は、ガウスの定理を適用して、面積分に変換し、近似的にφを積分の外に出してφとすると
【００３０】
【数９】
【００３１】
である。ここで、φはコントロールボリューム界面Ｓでの流束ｑsの正負に応じて、上流側の節点の値を取る。また、コントロールボリュームｉに対する流束ｑsは
【００３２】
【数１０】
【００３３】
と離散化して計算することができる。
【００３４】
右辺第１、２、３項の拡散項のコントロールボリュームに対する積分はガウスの定理により面積分に変換して補間関数で積分して離散化すれば、以下のとおりである。
【００３５】
【数１１】
【００３６】
右辺第４項の圧力項は
【００３７】
【数１２】
【００３８】
で計算する。右辺第５項の重力項の積分は以下となる。
【００３９】
【数１３】
【００４０】
以上の式（８）〜式（１３）から、ｕiに関して整理すると
【００４１】
【数１４】
【００４２】
と表すことができる。ここでｂi u,(n)はｕに関する係数で、ｎステップの時間項、圧力項、重力項が含まれる。同様にｖ、ｗに関して
【００４３】
【数１５】
【００４４】
【数１６】
【００４５】
である。
【００４６】
またＴ、ｃ、ｆに関しても同様に
【００４７】
【数１７】
【００４８】
【数１８】
【００４９】
【数１９】
【００５０】
ここで、上付き添え字のＴはＴに関する係数、上付き添え字のｃはｃに関する係数、上付き添え字のｆはｆに関する係数であることを示す。
【００５１】
式（４）の連続の式を用いて本発明の一つである圧力の式を導く。式（１）の圧力項以外の補間関数を用いて半離散化すると
【００５２】
【数２０】
【００５３】
と表される。ここでｄi u,(n)はｂi x,(n)から圧力項に関するものを除いたものとなる。圧力ｐがコントロールボリューム内で一定と仮定し、積分の外に出すと
【００５４】
【数２１】
【００５５】
ここで
【００５６】
【数２２】
【００５７】
であり
【００５８】
【数２３】
【００５９】
である。同様にｖ、ｗに関しても
【００６０】
【数２４】
【００６１】
【数２５】
【００６２】
である。節点ｉに対するコントロールボリュームで（４）式を積分して離散化すると
【００６３】
【数２６】
【００６４】
ここでＮｉは節点ｉに関する補間関数である。これに式（２１）、式（２４）、式（２５）を代入すると次式の離散化した圧力の式をえる。
【００６５】
【数２７】
【００６６】
この式をすべてのコントロールボリュームに対して重ね合せ、節点での圧力に関する全体連立方程式をえる。
【００６７】
【数２８】
【００６８】
この式は、反復法で簡単に計算でき、求めた圧力の値が空間的振動を示すことはない。
【００６９】
図１は以上の計算式を用いて行う反復計算手順の例を示したフローチャートである。流れの解析では、はじめに３次元解析領域を立体的な有限要素に分割して形状データを入力する。つぎに、流体の物性値を入力し、境界条件を入力し、初期条件を入力し、計算終了条件や反復計算の緩和係数などの解析条件を入力する。解析領域の形状のみに依存する積分を実施しておく。運動量の式（１４）〜（１６）の係数を計算し仮の圧力値で流速を求める。求めた流速で圧力の式（２８）の係数を計算し、圧力を求める。求めた圧力で次式を用いて流速を補正する。
【００７０】
【数２９】
【００７１】
【数３０】
【００７２】
【数３１】
【００７３】
必要であれば、温度の式（１７）、反応の式（１８）、自由表面の移動の式（１９）を解く。解が収束しなければ、流速、圧力を次式で更新し、必要であれば温度、反応率も次式で更新し
【００７４】
【数３２】
【００７５】
運動量の式の係数の計算に戻る。式（３２）において、φは流速または圧力または温度または反応率を表し、ｒは緩和係数で０＜ｒ＜２の範囲の値を取り、上付き添え字のＮＥＷは新しく計算された値、ＯＬＤは古い計算値であることを示す。解が収束すれば、流体占有率の式（１９）を用いて流体占有率分布を更新して自由表面を移動させ、つぎの時間のステップの計算に移る。最後の時間ステップであれば、計算を終了し、えられた解析結果をグラフィック処理して、等高線図、ベクトル分布図、グラフなどの形式で表示したプリンタなどに出力する。
【００７６】
図２は前述した３次元流れの有限要素解析に用いる解析装置の例を示す説明図である。補助記憶装置を内蔵するコンピュータ２０１に入力装置２０２、表示装置２０３が接続されている。また、プリンタ２０４、外付けの補助記憶装置２０５が接続されていても良い。入力装置２０２により３次元解析領域形状データと、流体の物性値、境界条件、初期条件を入力し補助記憶装置に保存される。このデータをもとに、コンピュータで解析を行い、えられた解析結果は通常補助記憶装置に保存され、これをたとえばグラフィック処理して表示装置３により表示する。結果の表示はプリンタ装置で行っても良い。
【００７７】
実施の形態１
図３は本発明の一実施の形態にかかわる解析対象である立方キャビティを示す斜視説明図である。図３において、３１は立方キャビティを示しており、ｗはＺ軸方向の流体の流速を示している。
【００７８】
ここでは図３に示すような一辺の長さが１で、上（ｙ＝１）のｘｙ面がｚ方向に速度１で動くばあいの立体キャビディ定常流れに適用した実施例を説明する。
【００７９】
境界条件は上（ｙ＝１）のｘｙ面でｕ＝０、ｖ＝０、ｗ＝１を与え、それ以外の壁面ではｕ＝０、ｖ＝０、ｗ＝０の条件を与えた。初期条件は、境界条件が与えられている節点以外はすべてｕ＝０、ｖ＝０、ｗ＝０、ｐ＝０とした。流体の物性値は、ρ＝１、η＝０．０１とした。また、流速および圧力の緩和係数は両者とも０．３とし、収束の判定は、相対誤差が０．００１以下になることとした。
【００８０】
解析領域を分割するために図４に示した６面体要素を用いる。図４は本発明の実施の形態１にかかわる６面体要素を示す説明図である。図４において、１、２、３、４、５、６、７および８は節点を示しており、４１は６面体要素を示している。また、９は座標軸α、βおよびγの座標軸原点である。コントロールボリュームは、たとえば節点１に対しては図４に示したように、６面体要素の中心と、節点１を共有する３つの面の中心と、節点１を共有する３つの辺の中心で定義する３つの面の内側をこの要素のサブコントロールボリュームとし、節点１を共有する要素に対するすべてのサブコントロールボリュームを加えたものとなる。解析領域は２０×２０×２０＝８０００要素に均等に分割した。なお、この要素における節点ｉに対する補間関数Ｎiとしては一般的な次式を用いた。
【００８１】
【数３３】
【００８２】
ここで、α、β、γは局所座標系で０から１までの値を取り、±は節点ｉでのＮiの値が１となるように＋か−となることを示す。全体座標系ｘ、ｙ、ｚで示されている積分は局所座標系にヤコビアンを用いて変換し、ガウスの求積法を用いて数値積分した。
【００８３】
この解析のばあい、定常流れであるため式（１）〜式（３）の時間微分は不要となり、時間ステップは１回のみとなる。また、温度の式（５）、反応の式（６）を計算する必要はない。
【００８４】
解析結果をグラフィック処理して図５の、立方キャビティでの圧力分布を示す説明図に示すような立方キャビティ５１についての圧力分布をえた。流速と圧力が同じ節点上で求まっていないときは圧力を流速が求まっている節点に補間する必要があるばあいが多いが、本発明のように流速と圧力が同じ節点上で求まっているばあいは圧力の補間などを行う必要がなく、グラフィック処理を簡単に行うことができた。この解析条件のばあいレイノルズ数は１００となるが、圧力分布はスムーズであり空間的な振動は見られない。図６は同時に計算された、立方キャビティでの流速ベクトル分布を示す説明図であり、６１は立方キャビティを示している。これにも空間的振動は見られず、良好な計算結果である。
【００８５】
なお、以上の解析は図２に示したハードウエア構成で、コンピュータにエンジニアリングワークステーションを用いて行った。収束までの反復回数は５２回で、計算時間は約３００秒とかなり短い時間で計算を完了した。
【００８６】
実施の形態２
ここでは図７に示すようなインサート物があるキャビティ（１／２モデル）の中に流体が流入する流れに適用した実施の形態を説明する。図７は、本発明にかかわる解析対象のインサート物をもつキャビティを示す説明図である。図７において、７１はキャビティであり、７２はインサート物であり、７３は対称面であり、７４はゲートである。各部の寸法は、Ｌ1が１０ｍｍ、Ｌ2が６ｍｍ、Ｌ3が１ｍｍである。
【００８７】
境界条件は入り口でｕ＝２（ｃｍ／ｓ）、ｖ＝０、ｗ＝０を与え、対称面では対称境界条件、それ以外の壁面ではｕ＝０、ｖ＝０、ｗ＝０の条件を与え、自由表面ではｐ＝０（Ｐａ）を与えた。初期条件は、境界条件が与えられている節点以外はすべてｕ＝０、ｖ＝０、ｗ＝０、ｐ＝０とした。流体の物性値は、ρ＝１（ｇ／ｃｍ3）、η＝１０（Ｐａ・ｓ）とした。また、流速および圧力の緩和係数はそれぞれ０．５と０．３とし、収束の判定は、相対誤差が０．０１以下になることとした。解析領域の分割には実施の形態１と同様な６面体要素を用い、約１００００要素に分割した。この解析のばあい、温度の式（５）、反応の式（６）は計算しない。
【００８８】
自由表面の移動の計算は以下のように移流項を計算するのと同じように式（９）を用いて行った。初期条件として、入り口節点は最初流体で満たされているとして、ｆ＝１とし、それ以外の節点ではｆ＝０とした。つぎに、ｆ＝１の節点に隣接する節点でｆ＜１の節点がｆ＝１になる時間を式（９）を用いて対応するコントロールボリュームに流入する流束を求めることによって行い、これを時間間隔Δｔとしてつぎの時間ステップの計算を行う。隣接するｆ＜０となる節点が多数あるばあいは最小のΔｔを選択する。各時間ステップの計算はｆ＞０の領域に対して行う。すべての節点でｆ＞０となるか、時間ステップが設定した最大数を超えれば計算を完了する。
【００８９】
図８は、解析結果をグラフィック処理して本発明の実施の形態２におけるインサート物をもつキャビティでの流動パターンを時間の順に（ａ）〜（ｄ）に示した説明図である。図８において、８１はキャビティであり、８３は対称面であり、８４はゲートである。キャビティ８１は図７に示したインサート物７２と同様のインサート物を含んでいる（図示せず）。これは、ｆ＝０．５の点が自由表面であるとして、時間とともに自由表面がどのように移動するかを示したものである。また、図９および図１０は、それぞれ本発明の実施の形態２にかかわるインサート物をもつキャビティでの圧力分布、流速分布を示す説明図である。９１はキャビティでの圧力分布、１０１はキャビティでの流速分布をそれぞれ示している。通常の３次元流れ解析に比較して要素分割は荒い方であるが、流動パターン、図９に示した流動途中の圧力分布はスムーズであり空間的な振動は見られない。図１０に示した流動途中の流速においても空間的振動は見られず、良好な計算結果である。これらの結果は、コントロールボリューム法を用いていることによって荒い要素分割でも安定な計算が可能となり、結果として計算時間が短くなることを示している。
【００９０】
なお、以上の解析は図２に示したハードウエア構成で、コンピュータにエンジニアリングワークステーションを用いて行った。非定常で、自由表面の移動も計算しているにもかかわらず計算時間は約１時間とかなり短い時間で計算を完了した。
【００９１】
実施の形態３
ここでは熱硬化性樹脂成形において、実施の形態２と同様に図７に示すキャビティ（１／２モデル）に樹脂が流入する過程の流れ解析の一例として行った実施の形態を説明する。
【００９２】
境界条件は入り口でｕ＝２（ｃｍ／ｓ）、ｖ＝０、ｗ＝０、Ｔ＝１２０℃、ｃ＝１を与え、対称面では対称境界条件、それ以外の壁面ではｕ＝０、ｖ＝０、ｗ＝０、Ｔ＝１８０（℃）の条件を与え、自由表面ではＰ＝０を与えた。初期条件は、境界条件が与えられている節点以外はすべてｕ＝０、ｖ＝０、ｗ＝０、ｐ＝０とした。熱硬化性樹脂にはエポキシ樹脂を用い物性値を測定した。測定した物性値は、ρ＝１．７（ｇ／ｃｍ3）、ηは非ニュートニアンであり、せん断速度γに依存し、温度と反応率にも依存するため以下の式で表され
【００９３】
【数３４】
【００９４】
ここで
【００９５】
【数３５】
【００９６】
【数３６】
【００９７】
ｎ＝０．７７、Ｂ＝４．０９×１０-6（Ｐａ・ｓ）、τ*＝２．５２×１０-9、Ｔb＝１１７９４（Ｋ）、ｃg＝０．５７、δ＝６．９２、またｃv＝０．９７（Ｊ／ｇ）、ｋ＝０．６７（Ｊ／ｍ・ｓ・Ｋ）、Ｈ＝２６０、ｇx＝０、ｇy＝０、ｇz＝９８０（ｃｍ／ｓ2）、Ｄ＝０である。反応速度Ｒには次式を用いた。
【００９８】
【数３７】
【００９９】
ここで
【０１００】
【数３８】
【０１０１】
【数３９】
【０１０２】
であり、Ａ1＝２．４４×１０9（１／ｓ）、Ｅ1＝１２１７６（Ｋ）、Ａ2＝１．４５×１０2（１／ｓ）、Ｅ2＝４１６７（Ｋ）、μ＝０．６７、ν＝０．９５とした。
【０１０３】
また、流速、圧力、温度および反応率の緩和係数はそれぞれ０．５、０．３、０．５および０．５とし、収束の判定は、相対誤差が０．０１以下になることとした。解析領域の分割には実施の形態１と同様な６面体要素を用い、約３０００要素に分割した。自由表面の移動の計算は実施の形態２のばあいと同様に行った。
【０１０４】
図１１は解析結果をグラフィック処理して本発明の実施の形態３におけるインサート物をもつキャビティでの流動パターンを時間の順に（ａ）〜（ｄ）に示した説明図である。図１１において、１１１はキャビティであり、１１３は対称面であり、１１４はゲートである。キャビティ１１１は図７に示したインサート物７２と同様のインサート物を含んでいる（図示せず）。これは、ｆ＝０．５の点が自由表面であるとして、時間とともに自由表面がどのように移動するかを示したものである。また、図１２および図１３は、それぞれ本発明の実施の形態３にかかわるインサート物をもつキャビティでの圧力分布、流速分布を示す説明図である。１２１はキャビティでの圧力分布、１３１はキャビティでの温度分布をそれぞれ示している。さらに、図１４は流動途中の反応率分布を示す説明図である。１４１はキャビティでの反応率分布を示している。図１３および図１４においても、キャビティは、図７に示したインサート物と同様のインサート物を含んでいる（図示せず）。図１２において、Ｐ1は０．９ＭＰａの等圧線を示しており、Ｐ2は０ＭＰａの等圧線を示している。また、図１３において、Ｔ1は１４０℃の等温線を示しており、Ｔ2は１８０℃の等温線を示している。さらに、図１４において、Ｒ1は０．０５の等反応率線を示している。通常の３次元流れ解析に比較して要素分割は荒い方であるが、流動パターン、図１２に示した流動途中の圧力分布を示す等圧線の形状は図１３においても、キャビティは、図７にスムーズであり空間的な振動は見られない。図１３に示した流動途中の温度においても等温線の形状に空間的振動が見られず、良好な計算結果である。図１４に示した流動途中の反応率においても等反応率線の形状に空間的振動は見られず、良好な計算結果である。これらの結果は、コントロールボリューム法を用いていることによって荒い要素分割でも安定な計算が可能となり、結果として計算時間が短くなることを示している。
【０１０５】
なお、以上の解析は図２に示したハードウエア構成で、コンピュータにエンジニアリングワークステーションを用いて行った。非定常で、自由表面の移動と温度と反応率を計算しているにもかかわらず計算時間は約３時間とかなり短い時間で計算を完了した。
【０１０６】
実施の形態４
実施の形態３のばあいと同様だが、樹脂の種類を変えて解析を行った。測定した樹脂の物性値は、粘度に関しては、ｎ＝０．８５、Ｂ＝４．９２×１０-4（Ｐａ・ｓ）、τ*＝２．８０×１０-7、Ｔｂ＝１４７９（Ｋ）、ｃg＝０．４４、δ＝３．３６であり、反応速度式に関しては、Ａ1＝６．７５×１０-1、Ｅ1＝３８２９、Ａ2＝３．８６×１０4（１／ｓ）、Ｅ2＝６４９６（Ｋ）、μ＝０．５２、ν＝０．１０１であり、その他の物性値は実施の形態３のばあいと同じである。
【０１０７】
図１５は本発明の実施の形態４におけるインサート物をもつキャビティでの圧力分布の解析結果を示す説明図である。図１４において、Ｐ3は０．５ＭＰａ等圧線を示しており、Ｐ4は０ＭＰａの等圧線を示している。実施の形態３の圧力分布の解析結果と比較して、インサート物上下の圧力差が小さく、流動圧力も小さい。
【０１０８】
えられたこれらの解析結果は、ここで用いた樹脂の方が成形性に優れ、流動過程においてインサート物を変形させたり移動させたりすることが少ないと考えられる。この樹脂を用いて成形品を製造した結果、インサート物が変形したり成形品表面から見える位置に移動することなどがない良好な成形品をえることができた。
【０１０９】
実施の形態５
本発明の３次元流れの有限要素解析法を図２の装置を用い、エディタープログラムを用いてＦＯＲＴＲＡＮ言語でコーディングして補助記憶装置にソースプログラムとして保存した。もちろんＦＯＲＴＲＡＮ以外の言語、たとえばＰＡＳＣＡＬ言語やＣ言語を用いても良い。保存したソースプログラムをＦＯＲＴＲＡＮ言語のコンパイラープログラムを用いて実行可能なロードモジュールプログラムに変換した。ソースプログラムまたはロードモジュールプログラムもしくはソースプログラムとロードモジュールプログラムをもとにＣＤ−ＲＯＭ成形用のマスター金型を作製し、このマスター金型によってソースプログラムまたはロードモジュールプログラムもしくはソースプログラムとロードモジュールプログラムを保存したＣＤ−ＲＯＭを製造した。同様に、プログラムを保存する媒体としてフロッピーディスクを選び、ソースプログラムまたはロードモジュールプログラムもしくはソースプログラムとロードモジュールプログラムを図２のコンピュータにフレキシブルディスクドライブを接続しフレキシブルディスクに保存した。プログラムを保存する媒体としてはこれらの他にもちろん、ＤＡＴテープ、カセットテープ、ＭＤ、ＤＶＤ−ＲＯＭなどを用いることができる。
【０１１０】
【発明の効果】
本発明にかかわる３次元流れの有限要素解析装置によれば、コントロールボリューム法に基づいた離散化を行い、運動量の式と連続の式から導出した解が空間的振動をおこさない圧力の式を用いることによって計算を安定に行うことができ、比較的荒い要素分割でも精度良く計算できるために、実用的な計算時間で３次元モデルを用いた解析領域形状または成形品形状に忠実な解析を行うことができる。
【０１１１】
また、コントロールボリュームに適した移流項および自由表面の移動の計算方法を用いるため、粘度が小さいばあいや流速が大きいばあいおよび複雑形状のキャビティのなかの自由表面の移動が、実用的な計算時間で精度良い解析結果をえることができる。
【０１１２】
また、圧力と流速、ならびに温度と熱硬化反応率を同じ節点で求めることができるため、圧力を求める節点が流速とずれている一般的な流れ解析法と異なり、結果をグラフィック表示するばあい、圧力を流速が求まる節点まで補間することなく容易に表示などを行うことができる。
【０１１３】
また、成形品形状の設計、成形条件の設定、材料の選定などの成形条件を決定でき、品質に優れ信頼性の高い成形品を効率的に製造することができる。
【０１１４】
また、３次元流れの有限要素解析法をプログラム化してＣＤ−ＲＯＭやフロッピー（登録商標）ディスクなどの媒体に保存すれば、容易に他の装置で利用することができる。
【図面の簡単な説明】
【図１】 本発明の３次元流れの有限要素解析法の手順の例を示すフローチャートである。
【図２】 本発明の３次元流れの有限要素解析装置のハードウエア構成の例を示す説明図である。
【図３】 本発明の一実施の形態における解析対象の立方キャビティを示す説明図である。
【図４】 本発明の一実施の形態における６面体要素を示す説明図である。
【図５】 本発明の一実施の形態における立方キャビティでの圧力分布を示す説明図である。
【図６】 本発明の一実施の形態における立方キャビティでの流速分布を示す説明図である。
【図７】 本発明の一実施の形態における解析対象のインサート物をもつキャビティの１／２モデルを示す説明図である。
【図８】 本発明の一実施の形態におけるインサート物をもつキャビティでの流動パターンを示す説明図である。
【図９】 本発明の一実施の形態におけるインサート物をもつキャビティでの圧力分布を示す説明図である。
【図１０】 本発明の一実施の形態におけるインサート物をもつキャビティでの流速分布を示す説明図である。
【図１１】 本発明の一実施の形態におけるインサート物をもつキャビティでの熱硬化性樹脂の流動パターンを示す説明図である。
【図１２】 本発明の一実施の形態におけるインサート物をもつキャビティでの熱硬化性樹脂の圧力分布を示す説明図である。
【図１３】 本発明の一実施の形態におけるインサート物をもつキャビティでの熱硬化性樹脂の温度分布を示す説明図である。
【図１４】 本発明の一実施の形態におけるインサート物をもつキャビティでの熱硬化性樹脂の反応率分布を示す説明図である。
【図１５】 本発明の一実施の形態におけるインサート物をもつキャビティでの熱硬化性樹脂の圧力分布を示す説明図である。
【符号の説明】
２０１ コンピュータ、２０２ 入力装置、２０３ 表示装置、２０４ 内蔵補助記憶装置、２０５ 外付け補助記憶装置。
【発明の属する技術分野】
金属の鋳造、プラスチック成形などの分野の、流動過程の有限要素解析方法および解析装置に関する。
【０００２】
【従来の技術】
従来、非圧縮性流体の流れの有限要素解析を行うばあい、安定な解をえるためには圧力と流速を求める格子点を食い違わせる必要があったこれは、Ｎａｖｉｅｒ−Ｓｔｏｋｅｓ式などの流れの式において圧力が１次微分、流速が２次微分であることに起因している。このため、有限要素法では流速の補間次数より圧力の補間次数を下げることが行われる。２次元での解析においても通常、たとえば４角形の有限要素を用いるばあい、４つの頂点で流速を求め、４角形の中心で圧力を求める必要がある。このとき、流速の補間次数は１次で圧力は０次となる。２次元解析ではＧａｌｅｒｋｉｎ法を用いて同一格子点で流速と圧力を求める方法(Computer Methods in Applied Mechanics and Engineering、 58, (1986)、 J. G. Rice and R. J. Schnipke、 An Equal-Order Velocity-Pressure Formulation that does not Exhibit Spurious Pressure Modes、 P.135-149)があるが、補間関数で重み付けを行うＧａｌｅｒｋｉｎ法は質量、運動量やエネルギーの保存を明示的に考慮した方法ではないため、安定な計算には細かな有限要素分割が必要であった。
【０００３】
３次元流れの解析法として流動コンダクタンスという量を仮定し、流速の変数を消去して未知変数を圧力のみの１つに減らして計算する方法が特開平８−９９４３１号公報、特開平９−１５０４４３号公報などに開示されているが、多孔質中の浸透流を表す近似的なダルシー流れの式にしか用いることができない。
【０００４】
これまで、３次元流れの有限要素解析で自由表面の移動を計算するばあい、流れに乗って移動する質量０のマーカー粒子を仮定して計算する方法が主流であった。また、薄板流れ（Ｈｅｌｅ−Ｓｈａｗ流れ）のように近似解析では、ＶＯＦ（Volume of Fluid）法のように格子点において流体が占める割合を計算し、この値から自由表面の位置を求めることができる。自由表面の変化が簡単なばあいは、変化した流体領域の形状に合わせて、有限要素分割をやり直したり、有限要素を追加または消去していく方法がある。
【０００５】
【発明が解決しようとする課題】
流速と圧力が求まる格子点が食い違っていると、結果を表示するばあい、このままでは、流速と圧力は異なったデータ系列となり、統一的に処理できないという問題があった。また、Ｇａｌｅｒｋｉｎ法を用いた場合は、前述のように定式化においては質量やエネルギーなどの保存則を考慮していないため、安定な解を求めるためには、有限要素分割を細かくする必要があり、結果的に計算時間の増大や解の精度の低下を招いていた。
【０００６】
また、自由表面の計算でマーカー粒子を用いるばあいは、複雑な流れのばあい、時間が経過すると粒子の分布に偏りが発生する。このため、最初から莫大な数のマーカー粒子を用いるか、適当な時間ごとにマーカー粒子の再配置を行う必要があった。これも、結果的に計算時間の増大を招く。解析領域の再有限要素分割や有限要素の追加、消去する方法は、流れが分岐したり、合流するばあいを計算するには、有限要素の形状がひずんだり、重なったりするのを避けるための多くの条件を考慮しておかねばらない。ＶＯＦ法に関しては薄板流れへの適用のみであり、３次元流れへの適用は見当たらない。
【０００７】
本発明は、このような問題点を解決するもので、３次元流れにおける流速と圧力を同一の格子点上で求め、自由表面の移動を容易に計算することができるようにすることにある。
【０００８】
【課題を解決するための手段】
本発明は、３次元流れ解析の有限要素解析法であって、立体型要素を用い、コントロールボリューム法に基づいて離散化した運動量の式と、圧力項を除いて半離散化した運動量の式を離散化した連続の式に代入したあと圧力項を離散化してえられる圧力の式と、離散化した運動量の式から導出する速度補正式とを用いて流速および圧力を求める反復計算を行い、流速と圧力を同一節点上で求めることを特徴とする解析方法である。
【０００９】
また、本発明により、３次元流れ解析の有限要素解析法に関し、立体要素を用い、コントロールボリューム法に基づいて離散化した運動量の式と、圧力項を除いて半離散化した運動量の式を離散化した連続の式に代入したあと圧力項を離散化してえられる圧力の式と、離散化した運動量の式から導出する速度補正式とを用いて流速および圧力を求める反復計算を行い、流速と圧力を同一節点上で求めることを特徴とする解析装置を提供する。
【００１０】
また、本発明は、コントロールボリュームでは流入流出する流束を見積もることにより保存則を考慮した解析となるため、本質的に計算が安定で、移流項を計算するのと同じ方法で自由表面の移動も容易に計算できる解析方法である。
【００１１】
また、本発明により、コントロールボリュームでは流入流出する流束を見積もることにより保存則を考慮した解析となるため、本質的に計算が安定で、移流項を計算するのと同じ方法で自由表面の移動も容易に計算できる解析装置を提供する。
【００１２】
また、本発明は、温度を計算するエネルギーの式と反応を計算する拡散の式を容易に組み込み、コントロールボリュームでは流入流出する流束を見積もることにより保存則を考慮した解析となるため、本質的に計算が安定で、移流項を計算するのと同じ方法で自由表面の移動も容易に計算できる解析方法である。
【００１３】
また、本発明により、温度を計算するエネルギーの式と反応を計算する拡散の式を容易に組み込み、コントロールボリュームでは流入流出する流束を見積もることにより保存則を考慮した解析となるため、本質的に計算が安定で、移流項を計算するのと同じ方法で自由表面の移動も容易に計算できる解析装置を提供する。
【００１４】
また、本発明により、えられた自由表面の移動状況、流速分布、圧力分布、温度分布および物質の濃度分布の結果をもとに成形条件を決定し、この決定した成形条件で成形品を製造する製造方法を提供する。
【００１５】
また、本発明により、３次元流れの有限要素解析法のプログラムを記録した媒体を提供する。
【００１６】
【発明の実施の形態】
以下に、添付図を参照しつつ、本発明による３次元流れの有限要素解析法とその解析装置の一例である、熱硬化性樹脂成形における樹脂流動過程の解析法および解析装置の好ましい実施の形態の例を詳細に説明する。
【００１７】
熱硬化性樹脂の自由表面が移動する３次元流れの式は以下で表される。
【００１８】
【数１】
【数２】
【数３】
【数４】
【数５】
【数６】
【数７】
ここでｘ、ｙ、ｚは３次元の空間座標であり、ｕ、ｖ、ｗはそれぞれ座標軸方向の樹脂の流速、ｐ、Ｔ、ｃはそれぞれ圧力、温度および熱硬化性樹脂の既反応基濃度（反応率）である。ｆは流体占有率で、ｆ＝０の位置には流体が存在しないこと、０＜ｆ＜１の位置は自由表面、ｆ＝１の位置は流体が満ちていることを示す。また、ρは密度、ηは粘度、Ｃvは比熱、ｋは熱伝導率、Ｄは拡散係数である。また、ｇx、ｇy、ｇzはそれぞれ重力加速度のｘ方向成分、ｙ方向成分、ｚ方向成分であり、Ｈは反応熱、Ｒは反応速度である。
【００２６】
解析領域を有限要素で分割し、有限要素の各節点に対して、それぞれが重ならないようにコントロールボリュームを定義する。式（１）を、コントロールボリュームで積分し離散化する。積分記号の添え字ＣＶはコントロールボリュームで積分することを表す。左辺第１項の時間微分は以下のように後退差分で表す。
【００２７】
【数８】
ここで、Ｎｉは節点ｉに対する補間関数、ｕiは節点ｉでのｘ方向の流速である。また、（ｎ）はｎ回目の時間ステップであることを示し、ｎ＋１ステップとｎステップの時間間隔がΔｔである。
【００２９】
式（１）の左辺第２、３、４項の移流項は非線型であり、単純に通常の補間関数を用いて積分しても解の空間的振動などの不具合が発生する。このため、本発明の一つである、流れの方向を考慮し、コントロールボリューム法にマッチした計算法を適用する。なお、この計算法は式（５）のＴ、式（６）のｃおよび式（７）の自由表面の移動に対しても同様に用いることができる。ｕ、ｖ、ｗ、Ｔ、ｃ、ｆをまとめてφで表せば、移流項のコントロールボリュームに対する積分は、ガウスの定理を適用して、面積分に変換し、近似的にφを積分の外に出してφとすると
【００３０】
【数９】
である。ここで、φはコントロールボリューム界面Ｓでの流束ｑsの正負に応じて、上流側の節点の値を取る。また、コントロールボリュームｉに対する流束ｑsは
【００３２】
【数１０】
と離散化して計算することができる。
【００３４】
右辺第１、２、３項の拡散項のコントロールボリュームに対する積分はガウスの定理により面積分に変換して補間関数で積分して離散化すれば、以下のとおりである。
【００３５】
【数１１】
右辺第４項の圧力項は
【００３７】
【数１２】
で計算する。右辺第５項の重力項の積分は以下となる。
【００３９】
【数１３】
以上の式（８）〜式（１３）から、ｕiに関して整理すると
【００４１】
【数１４】
と表すことができる。ここでｂi u,(n)はｕに関する係数で、ｎステップの時間項、圧力項、重力項が含まれる。同様にｖ、ｗに関して
【００４３】
【数１５】
【数１６】
である。
【００４６】
またＴ、ｃ、ｆに関しても同様に
【００４７】
【数１７】
【数１８】
【数１９】
ここで、上付き添え字のＴはＴに関する係数、上付き添え字のｃはｃに関する係数、上付き添え字のｆはｆに関する係数であることを示す。
【００５１】
式（４）の連続の式を用いて本発明の一つである圧力の式を導く。式（１）の圧力項以外の補間関数を用いて半離散化すると
【００５２】
【数２０】
と表される。ここでｄi u,(n)はｂi x,(n)から圧力項に関するものを除いたものとなる。圧力ｐがコントロールボリューム内で一定と仮定し、積分の外に出すと
【００５４】
【数２１】
ここで
【００５６】
【数２２】
であり
【００５８】
【数２３】
である。同様にｖ、ｗに関しても
【００６０】
【数２４】
【数２５】
である。節点ｉに対するコントロールボリュームで（４）式を積分して離散化すると
【００６３】
【数２６】
ここでＮｉは節点ｉに関する補間関数である。これに式（２１）、式（２４）、式（２５）を代入すると次式の離散化した圧力の式をえる。
【００６５】
【数２７】
この式をすべてのコントロールボリュームに対して重ね合せ、節点での圧力に関する全体連立方程式をえる。
【００６７】
【数２８】
この式は、反復法で簡単に計算でき、求めた圧力の値が空間的振動を示すことはない。
【００６９】
図１は以上の計算式を用いて行う反復計算手順の例を示したフローチャートである。流れの解析では、はじめに３次元解析領域を立体的な有限要素に分割して形状データを入力する。つぎに、流体の物性値を入力し、境界条件を入力し、初期条件を入力し、計算終了条件や反復計算の緩和係数などの解析条件を入力する。解析領域の形状のみに依存する積分を実施しておく。運動量の式（１４）〜（１６）の係数を計算し仮の圧力値で流速を求める。求めた流速で圧力の式（２８）の係数を計算し、圧力を求める。求めた圧力で次式を用いて流速を補正する。
【００７０】
【数２９】
【数３０】
【数３１】
必要であれば、温度の式（１７）、反応の式（１８）、自由表面の移動の式（１９）を解く。解が収束しなければ、流速、圧力を次式で更新し、必要であれば温度、反応率も次式で更新し
【００７４】
【数３２】
運動量の式の係数の計算に戻る。式（３２）において、φは流速または圧力または温度または反応率を表し、ｒは緩和係数で０＜ｒ＜２の範囲の値を取り、上付き添え字のＮＥＷは新しく計算された値、ＯＬＤは古い計算値であることを示す。解が収束すれば、流体占有率の式（１９）を用いて流体占有率分布を更新して自由表面を移動させ、つぎの時間のステップの計算に移る。最後の時間ステップであれば、計算を終了し、えられた解析結果をグラフィック処理して、等高線図、ベクトル分布図、グラフなどの形式で表示したプリンタなどに出力する。
【００７６】
図２は前述した３次元流れの有限要素解析に用いる解析装置の例を示す説明図である。補助記憶装置を内蔵するコンピュータ２０１に入力装置２０２、表示装置２０３が接続されている。また、プリンタ２０４、外付けの補助記憶装置２０５が接続されていても良い。入力装置２０２により３次元解析領域形状データと、流体の物性値、境界条件、初期条件を入力し補助記憶装置に保存される。このデータをもとに、コンピュータで解析を行い、えられた解析結果は通常補助記憶装置に保存され、これをたとえばグラフィック処理して表示装置３により表示する。結果の表示はプリンタ装置で行っても良い。
【００７７】
実施の形態１
図３は本発明の一実施の形態にかかわる解析対象である立方キャビティを示す斜視説明図である。図３において、３１は立方キャビティを示しており、ｗはＺ軸方向の流体の流速を示している。
【００７８】
ここでは図３に示すような一辺の長さが１で、上（ｙ＝１）のｘｙ面がｚ方向に速度１で動くばあいの立体キャビディ定常流れに適用した実施例を説明する。
【００７９】
境界条件は上（ｙ＝１）のｘｙ面でｕ＝０、ｖ＝０、ｗ＝１を与え、それ以外の壁面ではｕ＝０、ｖ＝０、ｗ＝０の条件を与えた。初期条件は、境界条件が与えられている節点以外はすべてｕ＝０、ｖ＝０、ｗ＝０、ｐ＝０とした。流体の物性値は、ρ＝１、η＝０．０１とした。また、流速および圧力の緩和係数は両者とも０．３とし、収束の判定は、相対誤差が０．００１以下になることとした。
【００８０】
解析領域を分割するために図４に示した６面体要素を用いる。図４は本発明の実施の形態１にかかわる６面体要素を示す説明図である。図４において、１、２、３、４、５、６、７および８は節点を示しており、４１は６面体要素を示している。また、９は座標軸α、βおよびγの座標軸原点である。コントロールボリュームは、たとえば節点１に対しては図４に示したように、６面体要素の中心と、節点１を共有する３つの面の中心と、節点１を共有する３つの辺の中心で定義する３つの面の内側をこの要素のサブコントロールボリュームとし、節点１を共有する要素に対するすべてのサブコントロールボリュームを加えたものとなる。解析領域は２０×２０×２０＝８０００要素に均等に分割した。なお、この要素における節点ｉに対する補間関数Ｎiとしては一般的な次式を用いた。
【００８１】
【数３３】
ここで、α、β、γは局所座標系で０から１までの値を取り、±は節点ｉでのＮiの値が１となるように＋か−となることを示す。全体座標系ｘ、ｙ、ｚで示されている積分は局所座標系にヤコビアンを用いて変換し、ガウスの求積法を用いて数値積分した。
【００８３】
この解析のばあい、定常流れであるため式（１）〜式（３）の時間微分は不要となり、時間ステップは１回のみとなる。また、温度の式（５）、反応の式（６）を計算する必要はない。
【００８４】
解析結果をグラフィック処理して図５の、立方キャビティでの圧力分布を示す説明図に示すような立方キャビティ５１についての圧力分布をえた。流速と圧力が同じ節点上で求まっていないときは圧力を流速が求まっている節点に補間する必要があるばあいが多いが、本発明のように流速と圧力が同じ節点上で求まっているばあいは圧力の補間などを行う必要がなく、グラフィック処理を簡単に行うことができた。この解析条件のばあいレイノルズ数は１００となるが、圧力分布はスムーズであり空間的な振動は見られない。図６は同時に計算された、立方キャビティでの流速ベクトル分布を示す説明図であり、６１は立方キャビティを示している。これにも空間的振動は見られず、良好な計算結果である。
【００８５】
なお、以上の解析は図２に示したハードウエア構成で、コンピュータにエンジニアリングワークステーションを用いて行った。収束までの反復回数は５２回で、計算時間は約３００秒とかなり短い時間で計算を完了した。
【００８６】
実施の形態２
ここでは図７に示すようなインサート物があるキャビティ（１／２モデル）の中に流体が流入する流れに適用した実施の形態を説明する。図７は、本発明にかかわる解析対象のインサート物をもつキャビティを示す説明図である。図７において、７１はキャビティであり、７２はインサート物であり、７３は対称面であり、７４はゲートである。各部の寸法は、Ｌ1が１０ｍｍ、Ｌ2が６ｍｍ、Ｌ3が１ｍｍである。
【００８７】
境界条件は入り口でｕ＝２（ｃｍ／ｓ）、ｖ＝０、ｗ＝０を与え、対称面では対称境界条件、それ以外の壁面ではｕ＝０、ｖ＝０、ｗ＝０の条件を与え、自由表面ではｐ＝０（Ｐａ）を与えた。初期条件は、境界条件が与えられている節点以外はすべてｕ＝０、ｖ＝０、ｗ＝０、ｐ＝０とした。流体の物性値は、ρ＝１（ｇ／ｃｍ3）、η＝１０（Ｐａ・ｓ）とした。また、流速および圧力の緩和係数はそれぞれ０．５と０．３とし、収束の判定は、相対誤差が０．０１以下になることとした。解析領域の分割には実施の形態１と同様な６面体要素を用い、約１００００要素に分割した。この解析のばあい、温度の式（５）、反応の式（６）は計算しない。
【００８８】
自由表面の移動の計算は以下のように移流項を計算するのと同じように式（９）を用いて行った。初期条件として、入り口節点は最初流体で満たされているとして、ｆ＝１とし、それ以外の節点ではｆ＝０とした。つぎに、ｆ＝１の節点に隣接する節点でｆ＜１の節点がｆ＝１になる時間を式（９）を用いて対応するコントロールボリュームに流入する流束を求めることによって行い、これを時間間隔Δｔとしてつぎの時間ステップの計算を行う。隣接するｆ＜０となる節点が多数あるばあいは最小のΔｔを選択する。各時間ステップの計算はｆ＞０の領域に対して行う。すべての節点でｆ＞０となるか、時間ステップが設定した最大数を超えれば計算を完了する。
【００８９】
図８は、解析結果をグラフィック処理して本発明の実施の形態２におけるインサート物をもつキャビティでの流動パターンを時間の順に（ａ）〜（ｄ）に示した説明図である。図８において、８１はキャビティであり、８３は対称面であり、８４はゲートである。キャビティ８１は図７に示したインサート物７２と同様のインサート物を含んでいる（図示せず）。これは、ｆ＝０．５の点が自由表面であるとして、時間とともに自由表面がどのように移動するかを示したものである。また、図９および図１０は、それぞれ本発明の実施の形態２にかかわるインサート物をもつキャビティでの圧力分布、流速分布を示す説明図である。９１はキャビティでの圧力分布、１０１はキャビティでの流速分布をそれぞれ示している。通常の３次元流れ解析に比較して要素分割は荒い方であるが、流動パターン、図９に示した流動途中の圧力分布はスムーズであり空間的な振動は見られない。図１０に示した流動途中の流速においても空間的振動は見られず、良好な計算結果である。これらの結果は、コントロールボリューム法を用いていることによって荒い要素分割でも安定な計算が可能となり、結果として計算時間が短くなることを示している。
【００９０】
なお、以上の解析は図２に示したハードウエア構成で、コンピュータにエンジニアリングワークステーションを用いて行った。非定常で、自由表面の移動も計算しているにもかかわらず計算時間は約１時間とかなり短い時間で計算を完了した。
【００９１】
実施の形態３
ここでは熱硬化性樹脂成形において、実施の形態２と同様に図７に示すキャビティ（１／２モデル）に樹脂が流入する過程の流れ解析の一例として行った実施の形態を説明する。
【００９２】
境界条件は入り口でｕ＝２（ｃｍ／ｓ）、ｖ＝０、ｗ＝０、Ｔ＝１２０℃、ｃ＝１を与え、対称面では対称境界条件、それ以外の壁面ではｕ＝０、ｖ＝０、ｗ＝０、Ｔ＝１８０（℃）の条件を与え、自由表面ではＰ＝０を与えた。初期条件は、境界条件が与えられている節点以外はすべてｕ＝０、ｖ＝０、ｗ＝０、ｐ＝０とした。熱硬化性樹脂にはエポキシ樹脂を用い物性値を測定した。測定した物性値は、ρ＝１．７（ｇ／ｃｍ3）、ηは非ニュートニアンであり、せん断速度γに依存し、温度と反応率にも依存するため以下の式で表され
【００９３】
【数３４】
ここで
【００９５】
【数３５】
【数３６】
ｎ＝０．７７、Ｂ＝４．０９×１０-6（Ｐａ・ｓ）、τ*＝２．５２×１０-9、Ｔb＝１１７９４（Ｋ）、ｃg＝０．５７、δ＝６．９２、またｃv＝０．９７（Ｊ／ｇ）、ｋ＝０．６７（Ｊ／ｍ・ｓ・Ｋ）、Ｈ＝２６０、ｇx＝０、ｇy＝０、ｇz＝９８０（ｃｍ／ｓ2）、Ｄ＝０である。反応速度Ｒには次式を用いた。
【００９８】
【数３７】
ここで
【０１００】
【数３８】
【数３９】
であり、Ａ1＝２．４４×１０9（１／ｓ）、Ｅ1＝１２１７６（Ｋ）、Ａ2＝１．４５×１０2（１／ｓ）、Ｅ2＝４１６７（Ｋ）、μ＝０．６７、ν＝０．９５とした。
【０１０３】
また、流速、圧力、温度および反応率の緩和係数はそれぞれ０．５、０．３、０．５および０．５とし、収束の判定は、相対誤差が０．０１以下になることとした。解析領域の分割には実施の形態１と同様な６面体要素を用い、約３０００要素に分割した。自由表面の移動の計算は実施の形態２のばあいと同様に行った。
【０１０４】
図１１は解析結果をグラフィック処理して本発明の実施の形態３におけるインサート物をもつキャビティでの流動パターンを時間の順に（ａ）〜（ｄ）に示した説明図である。図１１において、１１１はキャビティであり、１１３は対称面であり、１１４はゲートである。キャビティ１１１は図７に示したインサート物７２と同様のインサート物を含んでいる（図示せず）。これは、ｆ＝０．５の点が自由表面であるとして、時間とともに自由表面がどのように移動するかを示したものである。また、図１２および図１３は、それぞれ本発明の実施の形態３にかかわるインサート物をもつキャビティでの圧力分布、流速分布を示す説明図である。１２１はキャビティでの圧力分布、１３１はキャビティでの温度分布をそれぞれ示している。さらに、図１４は流動途中の反応率分布を示す説明図である。１４１はキャビティでの反応率分布を示している。図１３および図１４においても、キャビティは、図７に示したインサート物と同様のインサート物を含んでいる（図示せず）。図１２において、Ｐ1は０．９ＭＰａの等圧線を示しており、Ｐ2は０ＭＰａの等圧線を示している。また、図１３において、Ｔ1は１４０℃の等温線を示しており、Ｔ2は１８０℃の等温線を示している。さらに、図１４において、Ｒ1は０．０５の等反応率線を示している。通常の３次元流れ解析に比較して要素分割は荒い方であるが、流動パターン、図１２に示した流動途中の圧力分布を示す等圧線の形状は図１３においても、キャビティは、図７にスムーズであり空間的な振動は見られない。図１３に示した流動途中の温度においても等温線の形状に空間的振動が見られず、良好な計算結果である。図１４に示した流動途中の反応率においても等反応率線の形状に空間的振動は見られず、良好な計算結果である。これらの結果は、コントロールボリューム法を用いていることによって荒い要素分割でも安定な計算が可能となり、結果として計算時間が短くなることを示している。
【０１０５】
なお、以上の解析は図２に示したハードウエア構成で、コンピュータにエンジニアリングワークステーションを用いて行った。非定常で、自由表面の移動と温度と反応率を計算しているにもかかわらず計算時間は約３時間とかなり短い時間で計算を完了した。
【０１０６】
実施の形態４
実施の形態３のばあいと同様だが、樹脂の種類を変えて解析を行った。測定した樹脂の物性値は、粘度に関しては、ｎ＝０．８５、Ｂ＝４．９２×１０-4（Ｐａ・ｓ）、τ*＝２．８０×１０-7、Ｔｂ＝１４７９（Ｋ）、ｃg＝０．４４、δ＝３．３６であり、反応速度式に関しては、Ａ1＝６．７５×１０-1、Ｅ1＝３８２９、Ａ2＝３．８６×１０4（１／ｓ）、Ｅ2＝６４９６（Ｋ）、μ＝０．５２、ν＝０．１０１であり、その他の物性値は実施の形態３のばあいと同じである。
【０１０７】
図１５は本発明の実施の形態４におけるインサート物をもつキャビティでの圧力分布の解析結果を示す説明図である。図１４において、Ｐ3は０．５ＭＰａ等圧線を示しており、Ｐ4は０ＭＰａの等圧線を示している。実施の形態３の圧力分布の解析結果と比較して、インサート物上下の圧力差が小さく、流動圧力も小さい。
【０１０８】
えられたこれらの解析結果は、ここで用いた樹脂の方が成形性に優れ、流動過程においてインサート物を変形させたり移動させたりすることが少ないと考えられる。この樹脂を用いて成形品を製造した結果、インサート物が変形したり成形品表面から見える位置に移動することなどがない良好な成形品をえることができた。
【０１０９】
実施の形態５
本発明の３次元流れの有限要素解析法を図２の装置を用い、エディタープログラムを用いてＦＯＲＴＲＡＮ言語でコーディングして補助記憶装置にソースプログラムとして保存した。もちろんＦＯＲＴＲＡＮ以外の言語、たとえばＰＡＳＣＡＬ言語やＣ言語を用いても良い。保存したソースプログラムをＦＯＲＴＲＡＮ言語のコンパイラープログラムを用いて実行可能なロードモジュールプログラムに変換した。ソースプログラムまたはロードモジュールプログラムもしくはソースプログラムとロードモジュールプログラムをもとにＣＤ−ＲＯＭ成形用のマスター金型を作製し、このマスター金型によってソースプログラムまたはロードモジュールプログラムもしくはソースプログラムとロードモジュールプログラムを保存したＣＤ−ＲＯＭを製造した。同様に、プログラムを保存する媒体としてフロッピーディスクを選び、ソースプログラムまたはロードモジュールプログラムもしくはソースプログラムとロードモジュールプログラムを図２のコンピュータにフレキシブルディスクドライブを接続しフレキシブルディスクに保存した。プログラムを保存する媒体としてはこれらの他にもちろん、ＤＡＴテープ、カセットテープ、ＭＤ、ＤＶＤ−ＲＯＭなどを用いることができる。
【０１１０】
【発明の効果】
本発明にかかわる３次元流れの有限要素解析装置によれば、コントロールボリューム法に基づいた離散化を行い、運動量の式と連続の式から導出した解が空間的振動をおこさない圧力の式を用いることによって計算を安定に行うことができ、比較的荒い要素分割でも精度良く計算できるために、実用的な計算時間で３次元モデルを用いた解析領域形状または成形品形状に忠実な解析を行うことができる。
【０１１１】
また、コントロールボリュームに適した移流項および自由表面の移動の計算方法を用いるため、粘度が小さいばあいや流速が大きいばあいおよび複雑形状のキャビティのなかの自由表面の移動が、実用的な計算時間で精度良い解析結果をえることができる。
【０１１２】
また、圧力と流速、ならびに温度と熱硬化反応率を同じ節点で求めることができるため、圧力を求める節点が流速とずれている一般的な流れ解析法と異なり、結果をグラフィック表示するばあい、圧力を流速が求まる節点まで補間することなく容易に表示などを行うことができる。
【０１１３】
また、成形品形状の設計、成形条件の設定、材料の選定などの成形条件を決定でき、品質に優れ信頼性の高い成形品を効率的に製造することができる。
【０１１４】
また、３次元流れの有限要素解析法をプログラム化してＣＤ−ＲＯＭやフロッピー（登録商標）ディスクなどの媒体に保存すれば、容易に他の装置で利用することができる。
【図面の簡単な説明】
【図１】 本発明の３次元流れの有限要素解析法の手順の例を示すフローチャートである。
【図２】 本発明の３次元流れの有限要素解析装置のハードウエア構成の例を示す説明図である。
【図３】 本発明の一実施の形態における解析対象の立方キャビティを示す説明図である。
【図４】 本発明の一実施の形態における６面体要素を示す説明図である。
【図５】 本発明の一実施の形態における立方キャビティでの圧力分布を示す説明図である。
【図６】 本発明の一実施の形態における立方キャビティでの流速分布を示す説明図である。
【図７】 本発明の一実施の形態における解析対象のインサート物をもつキャビティの１／２モデルを示す説明図である。
【図８】 本発明の一実施の形態におけるインサート物をもつキャビティでの流動パターンを示す説明図である。
【図９】 本発明の一実施の形態におけるインサート物をもつキャビティでの圧力分布を示す説明図である。
【図１０】 本発明の一実施の形態におけるインサート物をもつキャビティでの流速分布を示す説明図である。
【図１１】 本発明の一実施の形態におけるインサート物をもつキャビティでの熱硬化性樹脂の流動パターンを示す説明図である。
【図１２】 本発明の一実施の形態におけるインサート物をもつキャビティでの熱硬化性樹脂の圧力分布を示す説明図である。
【図１３】 本発明の一実施の形態におけるインサート物をもつキャビティでの熱硬化性樹脂の温度分布を示す説明図である。
【図１４】 本発明の一実施の形態におけるインサート物をもつキャビティでの熱硬化性樹脂の反応率分布を示す説明図である。
【図１５】 本発明の一実施の形態におけるインサート物をもつキャビティでの熱硬化性樹脂の圧力分布を示す説明図である。
【符号の説明】
２０１ コンピュータ、２０２ 入力装置、２０３ 表示装置、２０４ 内蔵補助記憶装置、２０５ 外付け補助記憶装置。
Claims (4)
- 熱硬化性樹脂の物性値に基づく熱硬化性樹脂の３次元流れの有限要素解析装置であって、３次元解析領域を立体的な有限要素に分割した形状データと前記熱硬化性樹脂の密度、粘度、比熱、熱伝導率、拡散係数および反応速度と、境界条件および初期条件とを入力する手段と、入力した前記形状データと前記熱硬化性樹脂の密度、粘度、比熱、熱伝導率、拡散係数および反応速度と境界条件および初期条件とを記憶する補助記憶装置と、立体要素を用いるコントロールボリューム法に基づいて離散化した運動量の式
ｂ i u,(n) は運動量の式を離散化する過程で定められるｕに関する係数で、ｎステップの時間項、圧力項、重力項が含まれる。ｂ i v,(n) 、ｂ i w,(n) も同様にそれぞれｖおよびｗに関する係数である）
と、圧力項を除いてコントロールボリューム法に基づいて半離散化した運動量の式
と、コントロールボリューム法に基づいて離散化した運動量の式から導出する流速補正式
前記補助記憶装置に記憶された形状データと
前記熱硬化性樹脂の密度、粘度、比熱、熱伝導率、拡散係数および反応速度と
前記境界条件および初期条件とから、コンピュータにより流速および圧力を求める反復計算を行い、前記流速と圧力を同一節点上で解を収束させて前記流速と圧力とを求める３次元流れの有限要素解析手段と、前記同一節点上で求められた流速と圧力を記憶する補助記憶装置とを備え、前記流速と圧力により熱硬化性樹脂の３次元流れを求めることを特徴とする３次元の有限要素解析装置。 - コントロールボリュームに流入流出する流束を計算することによって自由表面の移動を計算する解析手段を含む請求項１記載の３次元流れの有限要素解析装置。
- コントロールボリューム法に基づいて離散化したエネルギーの式
- コンピュータを、請求項１記載の３次元流れの有限要素解析装置として動作させるプログラムを記録した記録媒体。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP02789698A JP3721765B2 (ja) | 1998-02-10 | 1998-02-10 | ３次元流れの有限要素解析法、解析装置、成形品を製造する製造方法および解析法のプログラムを記録した媒体 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP02789698A JP3721765B2 (ja) | 1998-02-10 | 1998-02-10 | ３次元流れの有限要素解析法、解析装置、成形品を製造する製造方法および解析法のプログラムを記録した媒体 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JPH11232250A JPH11232250A (ja) | 1999-08-27 |
| JP3721765B2 true JP3721765B2 (ja) | 2005-11-30 |
Family
ID=12233659
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP02789698A Expired - Fee Related JP3721765B2 (ja) | 1998-02-10 | 1998-02-10 | ３次元流れの有限要素解析法、解析装置、成形品を製造する製造方法および解析法のプログラムを記録した媒体 |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JP3721765B2 (ja) |
Families Citing this family (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP4754718B2 (ja) * | 2001-05-16 | 2011-08-24 | 三菱電機株式会社 | ３次元発泡流動解析方法および解析装置、ならびに解析方法のプログラムを記録した媒体 |
| TWI241949B (en) * | 2001-06-08 | 2005-10-21 | Mitsubishi Heavy Ind Ltd | Method of analyzing injection molding conditions and method for providing the analysis results thereof |
| JP2005131879A (ja) * | 2003-10-29 | 2005-05-26 | Toyo Seiki Seisakusho:Kk | 樹脂粘度特性試験システム、その方法、及びそのプログラム |
| JP4604260B2 (ja) * | 2004-09-30 | 2011-01-05 | 学校法人慶應義塾 | ＧＳＭＡＣ有限要素法による高次混合要素Ｐｏｉｓｓｏｎソルバーの数値計算手法 |
| JP4820318B2 (ja) * | 2007-03-22 | 2011-11-24 | 株式会社日立製作所 | 樹脂成形品の設計支援装置、支援方法及び支援プログラム |
| JP6164222B2 (ja) * | 2012-10-31 | 2017-07-19 | 旭硝子株式会社 | シミュレーション装置、シミュレーション方法およびプログラム |
| JP6092832B2 (ja) * | 2014-10-21 | 2017-03-08 | 株式会社日本製鋼所 | 流動挙動予測装置、その流動挙動予測方法および流動挙動予測プログラム |
| WO2016157335A1 (ja) * | 2015-03-27 | 2016-10-06 | 株式会社日立製作所 | 電子部品の樹脂流動解析方法、伝熱解析方法、熱応力解析方法、およびその装置 |
| DE102015107025B4 (de) * | 2015-05-06 | 2021-01-21 | Bt Bayern Treuhand Management & Technologie Gmbh | Verfahren und Prozessparameterwert-Ermittlungseinrichtung zur Ermittlung und Anzeige von Prozessparameterwerten in einem Spritzgussprozess, Verfahren zur Steuerung und Steuereinrichtung für eine Spritzgießanlage sowie Spritzgießanlage |
| DE102018123361A1 (de) * | 2018-09-23 | 2020-03-26 | Arburg Gmbh + Co Kg | Verfahren zur Steuerung einer Maschine zur Verarbeitung von Kunststoffen |
| CN110146252B (zh) * | 2019-04-19 | 2020-10-13 | 浙江大学 | 一种考虑环境流体作用机制的密度流实验装置和方法 |
-
1998
- 1998-02-10 JP JP02789698A patent/JP3721765B2/ja not_active Expired - Fee Related
Also Published As
| Publication number | Publication date |
|---|---|
| JPH11232250A (ja) | 1999-08-27 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Liu et al. | Homogenization of material properties in additively manufactured structures | |
| Lan et al. | An adaptive finite volume method for incompressible heat flow problems in solidification | |
| JP3721765B2 (ja) | ３次元流れの有限要素解析法、解析装置、成形品を製造する製造方法および解析法のプログラムを記録した媒体 | |
| Elliott et al. | An ALE ESFEM for solving PDEs on evolving surfaces | |
| Coward et al. | Temporal evolution of periodic disturbances in two-layer Couette flow | |
| JP2004503005A (ja) | 流体流れシミュレーションのための同次数法 | |
| Zhuo et al. | LES-based filter-matrix lattice Boltzmann model for simulating turbulent natural convection in a square cavity | |
| Chan et al. | Modeling nonisothermal impregnation of fibrous media with reactive polymer resin | |
| Shojaei et al. | Numerical simulation of three-dimensional mold filling process in resin transfer molding using quasi-steady state and partial saturation formulations | |
| CN113569450B (zh) | 一种预估与控制液滴悬浮与驻留的方法 | |
| Alphonse et al. | Directional differentiability for elliptic quasi-variational inequalities of obstacle type | |
| Karyofylli et al. | Simplex space-time meshes in thermally coupled two-phase flow simulations of mold filling | |
| Osswald et al. | A boundary element simulation of compression mold filling | |
| Soni | Fast and accurate numerical approaches for Stefan problems and crystal growth | |
| Hayes et al. | Numerical simulation of mold filling in reaction injection molding | |
| Stella et al. | The Rayleigh–Bénard problem in intermediate bounded domains | |
| Maliska et al. | An unstructured finite volume procedure for simulating flows with moving fronts | |
| Shojaei et al. | Numerical simulation of three-dimensional mold filling in resin transfer molding | |
| Béchet et al. | Adaptive mesh generation for mould filling problems in resin transfer moulding | |
| Bigot et al. | Capture of 3D moving free surfaces and material interfaces by mesh deformation | |
| Ninokata | An effective, locally exact finite-difference scheme for convection-diffusion problems | |
| Wang et al. | Adaptive grids for coupled viscous flow and transport | |
| Wang et al. | Critical evaluation of conservative and non-conservative interface treatment for Chimera grids | |
| Karyofylli et al. | Adaptive temporal refinement in injection molding | |
| Olesen | Computational fluid dynamics in microfluidic systems |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20041130 |
|
| A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20050127 |
|
| A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A821 Effective date: 20050127 |
|
| A02 | Decision of refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A02 Effective date: 20050524 |
|
| A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20050719 |
|
| A911 | Transfer to examiner for re-examination before appeal (zenchi) |
Free format text: JAPANESE INTERMEDIATE CODE: A911 Effective date: 20050728 |
|
| TRDD | Decision of grant or rejection written | ||
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20050823 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20050905 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20080922 Year of fee payment: 3 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20090922 Year of fee payment: 4 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20090922 Year of fee payment: 4 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20100922 Year of fee payment: 5 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20110922 Year of fee payment: 6 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20110922 Year of fee payment: 6 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20120922 Year of fee payment: 7 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20130922 Year of fee payment: 8 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| LAPS | Cancellation because of no payment of annual fees |