日本機械学会［No.167-1］北陸信越支部 第 53 期総会・講演会 講演論文集 [2016.3.5 長野県長野市]
［No.167-1］日本機械学会北陸信越支部 第 53 期総会・講演会 講演論文集 [2016.3.5 長野県長野市]
307
マイクロチャネル内の液滴生成に対する有限要素解析
田井 茂俊*1，倉橋 貴彦*2 Finite element analysis for droplet generation in micro-channel
Shigetoshi TAI*1 and Takahiko KURAHASHI
*1 Nagaoka University of technology,
Kamitomioka 1603-1, Nagaoka, Niigata, 940-2137 Japan
In recent years, various technology developments are carried out using a micro chemical chip. For instance, there is a micro droplet generation technology. In the medical fields, the micro-droplet generation technology used to drug delivery systems(DDS). In the DDS, it is important to generate the monodisperse capsule shaped particles such as targeting of the drug to the lesion and sustained-release control in the drug. In this study, bubble function FEM using the VOF method to fluid separation is carried out. However, from the numerical result, it is difficult to separate two liquids by VOF method. For this reason, the FEM using the fictitious domain method is applied to simulate a droplet shape in the flow field. This method represent a moving object by the foreground mesh, and show the flow field by background mesh. In this paper, verification of foreground mesh deformation by the interfacial tension force is carried out. In addition, numerical simulation of the droplet shape deformation in the flow field in the microchannel is carried out.
Key Words : Microchannel, FEM, VOF method, Fictitious domain method, Bubble function, Doroplet generation
1. 緒 言
近年，マイクロ化学チップ（図 1）を利用した様々な液滴生成技術が注目されており，例えば，輸送にマイ クロ液滴を適切に注入する方法(1)(2)に関する研究や医療分野におけるドラッグデリバリーシステム（DDS）では病 巣への薬剤の標的指向などで薬物を包んでいるカプセル形成微粒子の単分散化が重要な点として挙げられている．
マイクロスケール流路内において生成される液滴サイズおよび液滴間隔は，流入する二流体と流量条件の組み合 わせにより変化する．そこで，事前に予測できるツールの開発を目指して有限要素解析による検討を行った．検 討では流体解析領域における各計算格子を占める流体の体積率𝜙を0から1の間の値により表すVOF法を流体の 判別に適用し，数値解析法としてはSUPG法の安定化パラメータを使用した安定化気泡関数有限要素法を適用し た．数値解析の結果は図2のように示される．図 2 は，赤色で示される酢酸エチル（𝜙 = 1.0）と青色の蒸留水
（𝜙 = 0.0）で液滴を生成させる数値実験の結果である．液滴が作成される二流体の挙動を再現できることを確認
したが，生成された液滴（酢酸エチルの球）を体積率𝜙の0から1の範囲で二流体に判別することは今回の解析 条件の場合，困難な結果となった．
そこで，本論文では図3に示すような移動物体を前景領域，流れ場を背景領域として表す仮想領域有限
要素法(3)(4)(5)による液滴生成後のマイクロスケール流路における流れ解析について数値解析による検討を行
う．検討では，液滴を表す前景領域に与える界面張力により，前景領域のメッシュの形状変形を行う．前 景領域のメッシュの形状変形では，無重力状態での液滴形状の変形を考え，表面全体の自由エネルギーが
*1 学生員，長岡技術科学大学 大学院（〒940-2137 新潟県長岡市上富岡町1603-1）
*2 正員，長岡技術科学大学 大学院 E-mail: [email protected]
（4）
（5）
最小になるときの前景領域の形状を確認した．また，以下に示す支配方程式に対してマイクロスケール流 路内の流れ解析を液滴の形状変形を踏まえて行う．
2. 仮想領域有限要素法による離散化
計算領域のイメージ図を図4に示す．本検討では流れ場を再現するための支配方程式として，背景領域 Ωに対して，式（1）に示す非圧縮粘性流体に対するナビエ・ストークスの運動方程式と式（2）に示す連 続式を用いる．式（1）,（2）は総和規約による表示をしており，本研究では平面二次元を対象とする．
𝑉𝑖̇ + 𝑉𝑗𝑉𝑖,𝑗+1
𝜌𝑃,𝑖− 𝜐(𝑉𝑖,𝑗+ 𝑉𝑗,𝑖)
,𝑗= 𝑓𝑖 （1）
𝑉𝑖,𝑗= 0 （2）
ここで𝑉𝑖は流速，𝑃は圧力である．また，𝜐は動粘性係数，𝜌は密度，𝑓𝑖 =1𝜌𝐹𝑖を表し，𝐹𝑖は外力(チップを 水平に置いているため，体積力としての重力は無視し，表面(界面)張力のみ考慮する)を表す．マイクロス ケール流路内における流体の表面(界面)張力を表すため，界面上における表面張力(6)の式（3）を主法線ベ クトル・接線ベクトルおよび曲率の幾何学的関係から書き換えたものを適用し，体積力𝐹𝑖は式（4）の様に 表される．初期条件および境界条件は式（5），（6）のように定義する．また，前景領域ωにおいては，流 速𝑣𝑖𝑎が既知であるものとし，式（7）のように与える．離散化手法としては，気泡関数要素を用いた仮想 領域有限要素法を適用し，時間方向の離散化にはクランク・ニコルソン法を適用する．
Fig.4 Image diagram of computational domain
𝐹𝑖= 𝜎𝜅𝑛̅ （3）
𝐹𝑖= 𝜎𝑡Γ́
𝑉𝑖(𝑡0) = 𝑉̂ in Ω 𝑖
Fig.1 Micro-chemical chip Fig.2 Computational result by the VOF method Fig.3 Flow analysis around the cylinder by fictitious domain method
（6）
（7）
（9）
𝑉𝑖 = 𝑉̂ on Γ𝑖 1
𝑉𝑥= 0, 𝑉𝑦 = 0 on Γ2 𝑃 = 0 on Γ3
𝑉𝑖= 𝑉𝑖𝑎 in ω
ここで，𝜎は表面（界面）張力係数であり，𝜅は表面（界面）曲率，𝑛̅は主法線ベクトル，𝑡Γは接線ベクト ルである．
流速𝑉𝑖，圧力𝑃に対する補間関数を式（8），（9）に示す．気泡関数要素を流速の補間に適用し，圧力の 補間には三角形一次要素を用いる．
𝑉𝑖= 𝑁1𝑉𝑖1+ 𝑁2𝑉𝑖2+ 𝑁3𝑉𝑖3+ 𝑁4𝑉̃𝑖4= {𝑁𝐵}𝑇{𝑉𝑖} （8）
𝑃 = 𝑁1𝑃1+ 𝑁2𝑃2+ 𝑁3𝑃3= {𝑁}𝑇{𝑃}
ここで𝑁は面積座標（𝜂1,𝜂2,𝜂3）により表される面積座標を示しており，𝑁1=𝜂1，𝑁2=𝜂2，𝑁3=𝜂3， 𝑁4=27𝜂1𝜂2𝜂3と表される．
3. 仮想領域有限要素法によるマイクロスケール流れ場の解析 3・1 液滴を表す前景領域におけるメッシュの形状変形
表面張力による前景領域のメッシュの変形について考察した．前景領域は，境界上の表面張力から得られる 移動量を考え，ラプラス方程式を解くことで各節点の移動量を求め，変形させている．今回の検討では，無重力 状態を考え，表面張力以外の力が作用しないものとした．この時，液体の表面は表面自由エネルギーという単位 面積当たり一定のエネルギーを持つため，表面積を小さくすることで，表面全体のもつ自由エネルギーが最小に なるような液面形状（球状）となる．そこで本検討では，表面全体の自由エネルギーが最小となるよう，表面積 が一定の条件で，円周が最小になったときに液体の変形を停止する処理を取り入れた．計算条件は表 1，前景領 域のメッシュ形状の初期状態とその際に発生する表面張力はそれぞれ図5，6に示す，また，各時間経過後の液滴
（前景領域メッシュ）変形の様子を図7から図10に示す．
Table.1 Computational conditions of foreground mesh ω
ω
∆𝑡 (sec.)
𝜎 (N/mm) 6.3 × 10−6
Fig.5 Deformation of foreground mesh (Initial condition T=0 sec)
Fig.7 Deformation of foreground mesh (T=0.005 sec)
Fig.6 Surface force vector (Initial condition T=0 sec)
331 600 0.0001 Time steps 120
Interfacial tension coefficient Total number of nodes in foreground domain Total number of elements in foreground domain
Time increment
3・2 マイクロスケール内の流れ場における液滴変形について
3・1項の結果に基づき，図11に示すモデルに対して，液滴の変形と合わせて，液滴の周りの流れ解析を行っ た．今回の検討では，液滴は同一箇所に留まっているものとし，流路軸方向には移動しないものとしている．計 算条件を表2に示し，流れ場と液滴（前景領域）の初期状態をそれぞれ図12，13に示す．また，各時間経過後の 流れ場および液滴形状の様子を図14から図17に示し，流れ場においては流速の大きさをコンター図にして表示 する．この結果から，背景領域と前景領域において，楕円状の液滴は表面全体の自由エネルギーが小さくなるよ うに変形されると共に，液滴周りの流れ場の計算が行えることを確認した．したがって，マイクロスケール流路 において液滴形状と流れ場を考慮した数値解析が可能であることがわかった．
Fig.8 Deformation of foreground mesh (T=0.007 sec)
Fig.9 Deformation of foreground mesh (T=0.009 sec)
Fig.10 Deformation of foreground mesh (T=0.012 sec)
331 600 1,907 3,716 0.0001
100 Total number of nodes in foreground domain
Total number of elements in foreground domain
Time increment Time steps Interfacial tension coefficient Total number of nodes in background domain Total number of elements in background domain
ω ω Ω Ω
∆𝑡 (sec.)
𝜎 (N/mm) 6.3 × 10−6 Table.2 Computational conditions of flow analysis
Fig.12 Deformation of foreground mesh (Initial condition T=0 sec)
Fig.13 Flow field of background mesh (Initial condition T=0 sec)
Fig.14 Deformation of foreground mesh (T=0.003 sec)
Fig.15 Flow field of background mesh (T=0.003 sec)
Fig.11 Calculation model and boundary conditions
1.2(mm)
0.8(mm) p = 0.0[Pa]
[mm/s]
[mm/s]
Fictitious domain FEM in domain Center of ellipse
=
[mm/s]
[mm/s]
4. 結 語
本研究では，非圧縮粘性流体に対し，移動物体を前景領域，流れ場を背景領域として表す仮想領域有限要素法 による液滴生成後のマイクロスケール流路における流れ解析の検討を行った．以下に得られた知見を示す．
・液滴（前景領域）の形状は，表面張力から境界上の移動量を算出し，ラプラス方程式により全節点の 移動量を求めることで，変形させることができることを確認した．
・仮想領域有限要素法を用いることにより，楕円状の液滴を変形させ，それに伴った流れ場の計算を行え ることがわかった．
謝 辞
本研究を実施するにあたり，マイクロ化学技研株式会社の菊谷善国氏には多大なる助言を頂いた．ここ に記して謝意を表す．
文 献
(1) 原田隆宏，鈴木幸人，塚越誠一，荒川貴博，庄司習一，“MPS法によるマイクロ流体の液滴生成シミュレーション”， 日本機械学会論文集B編，Vol. 71, No. 711 (2005), pp. 25-29.
(2) 荒川貴博，庄司習一，“混相流マイクロ流体システムを用いた微小マイクロ液滴の生成と応用”，表面科学，Vol.29, No. 6 (2008), pp. 356-363.
(3) Glowinski,R. Pan,T.W., Hesla,T.I., Joseph,D.D, and Periaux,J., “A Fictitious Domain Method with Distributed Lagrange Multipliers for the Numerical simulation of Particulate Flow”, Contemp. Math., Vol.218, (1998), pp.121-137.
(4) Glowinski,R. Pan,T.W., Hesla,T.I., Joseph,D.D, and Periaux,J., “A Fictitious Domain Approach to the Direct Numerical Simulation of Incompressible Viscous Flow”, Commun. Comput. Phys., Vol.169, (2001), pp.363-426.
(5) Glowinski,R., and He,Q., “A Least-Squares/Fictitious Domain Method for Linear Elliptic Problems with Robin Boundary Conditions”, Commun. Comput. Phys., Vol.9, No3,(2011), pp.587-606.
(6) J.U. Brackbill, D.B.Kothe and C.Zemach, “A Continuum method for mokiling surface tension” , Journal of Computational physics, 100(1992), pp.335-354.
Fig.16 Deformation of foreground mesh (T=0.01 sec)
Fig.17 Flow field of background mesh (T=0.01 sec)
