流体中における繊毛集団モデルのダイナミクス

1(指導教員：郡宏)
0 < θ < π2 ) を用いてシミュレーションを行った. なお
A は回転速度の最大値と最小値の差を表すパラメータで, θ は円軌道上のどこで速度最大となるかを決定するパラメータである.

はじめに大腸菌やゾウリムシなど, 水中に棲息する単細胞生物の多くは鞭毛や繊毛の集団運動を利用して遊泳している. 例えばゾウリムシは, 細胞表面が多数の繊毛で覆われているが, これらが協同的に波打つことで流れを作り泳いでいる. これらの微生物は流体の慣性がほとんど効かない低レイノルズ数 (∼ 10−3 ) の環境下で運動しており, その運動には流体を介した相互作用が働いていると考えられている.[1]
本研究ではこのような流体中での繊毛集団の運動をモデル化し, そのダイナミクスを観察した. またその運動が作り出す流れの流速を測定し, ダイナミクスのパターンによって流れの方向に現れる変化を調べた.

2

まず始めに, 回転子を一列に並べたモデルについてシミュレーションを行った.

(a)t = 700

このモデルは隣り合う粒子の位相が少しづつずれて回転し, 波を作るような挙動を示し, 実際の微生物が作る繊毛波 (メタクロナルウェーブ) の様なものが観察できた. 図 1 は, t = 700 と t = 720 のときのモデルのスナップショットである. 粒子に左から順に 1,2,...8 と番号をつけると, 8 → 1 の方向に繊毛波が伝わっている様子が分かる.
次に, このモデルに壁の効果を導入することを考えた. 回転子が一列のモデルでは回転子が作った流れは無限遠に流れていってしまうが, 実際の微生物には身体があり, 身体があることで周りの流れに変化が生じるだろうと考えたからである. 壁には大きく分けて, 壁表面での水平方向・鉛直方向の流れを許さない固定壁[4] と壁表面で水平方向の流れのみ許す自由壁の 2 種類がある. 今回は自由壁の効果を鏡像を使用して導入し,壁が 1 つの場合と 2 つの場合 (ナノポンプ) についてシミュレーションを行った. これらはどちらも A = 0.2,
θ = π4 の場合, 回転子が一列の場合と同様に前から後ろへ繊毛波を伝えるような挙動を示した.
またこれら 3 つのモデルに対して, A を −1 < A < 1,
θ を 0 < θ < π2 の範囲で変化させたときの回転子が作る流れの流速 V の長時間平均を測定した. なお, 流速
V は左向きの流れを正として測定している.

本研究では, 先行研究で用いられている円軌道上を粒子が回転する回転子モデル [2] を一列に並べることで繊毛集団を表現する. 3 次元空間のストークス流体における, ある平面上での運動を考える. 粒子の位相は φ で, 駆動力 F によって反時計回りに回転している
とき, 各粒子の位相速度 φ̇ は軌道の接線方向の単位ベクトルt, 軌道R, 他の粒子から受ける流れ場v を用いて次のように表される.

v(r i )
|R0 (φi )|(1)

ただし ωi (φi ) は固有の位相速度で, 以下のように表せる.
F (φi )ωi (φi ) =(2)ζ|R0 (φi )|なお抵抗係数 ζ は粒子の半径 a, 流体の粘性率 µ として Stokes の式 ζ = 6πµa で与えられる.
各粒子が速度U で回転するとき, それぞれの粒子は流体から抵抗力g=ζ[v(R(φi ))−U ] を受け, 周りの流体はその反作用力 −g を受ける. これにより粒子の周りには流れが生じるが, この流れが粒子間に相互作用を引き起こす. 粒子 i の座標を (xi ,yi ), 距離を l として,粒子 i が粒子 j に与える相互作用の強さを表すオセー
ンテンソルは次のように書くことができる.
)((xi −xj )(yi −yj )(xi −xj )211l +l3
l3(3)
Gij =
(y −y )2i −yj )18πµ (xi −xj )(y+ i 3j3ll(b)t = 720

図 1: 回転子を一列に並べたモデル. A = 0.2, θ = π4 .

モデル

φ̇ = ωi (φi ) + t ·シミュレーション結果

30.1

0.16
0.14

0.09

0.12
<V>

0.08

0.1
<V>

吉田理紗子

0.07

0.08
0.06

0.06

0.04
0.05

0.02

0.04

00l

0.2 0.4 0.6 0.81

1.2 1.4 1.6

θ(a)hV i の θ 依存性

式 (3) を用いて, 粒子 i が粒子 j のところに作る流れ場
vj はvj =−Gij gi と表せる. これらの式より時間刻みごとに各粒子の位置にできる流れ場を計算し, それを用いて位相速度を計算することで粒子を回転させている.
パラメータ値は粒子の半径 a = 1.0, 軌道の半径 r =5.0, 軌道間距離 d = 15.0, 時間刻み幅 dt = 0.01, 回転子の個数 n = 8 とし, 初期の位相差は 0.1π である. また駆動力 F には, 先行研究 [3] で用いられていた F = F0 [1.0 + A sin(φ + θ)] (F0 > 0, −1 < A < 1,-1 -0.8-0.6-0.4-0.2 0 0.2 0.4 0.6 0.8 1
A(b)hV i の A 依存性

図 2: パラメータによる hV i の変化. 赤:回転子を一列に並べた (壁が無い) モデル, 緑:壁が 1 つのモデル, 青:
壁が 2 つのモデル. (a)A = 0.5, (b)θ = π4 .
まず図 2(a) を見ると, どのモデルも θ の値が小さいときに hV i の値が大きく, θ の値が大きいところではhV i の値に変化は見られないということが分かった.

5987654321sin φsin φ

また図 2(b) より, hV i の値は壁の有無に関わらず Aが −0.1 程度のときに最大となることが分かった. 壁 2つのモデルの値が非常に不安定なのは, 2 つの壁の間に複雑な流れが出来ていて, その流れが回転子が繊毛波を作るのを阻害してしまっていることが原因の 1 つとして挙げられるが, 明確な理由はまだ分かっていない.
またこのグラフは, A の値で大きく 3 つの範囲に分けることが出来る. hV i の値が大きい傾向にある A < 0の範囲, A が 0 付近の範囲, hV i の値が小さい傾向にある A > 0 の範囲の 3 つである. これら 3 つの範囲に対して, モデルがどのような挙動を示しているのかをシミュレーションを行い確認した. 図 3 はそれぞれの Aの値に対する, 各粒子の sin φ の時空プロットである.
波の方向は, 粒子 8 を前と仮定している.

0250

3 つのモデルに対して, 流速は粒子 1 の左側の位置で計測している. このとき, 相互作用は距離の関数なので, 回転子番号の小さい粒子 (つまり, 流速の測定位置に近い粒子) の作る流れが回転子番号の大きい粒子に比べて流速の値に大きな影響を及ぼしていると考えられる. また反時計回りに回転する粒子が任意の場所に作る水平方向の流れには, 0 < φ < π のときに出来る左向きの流れと, π < φ < 2π のときに出来る右向きの流れがある.
これらをふまえて A > 0 のときのスナップショット(図 1(a)) と A < 0 のときのスナップショット (図 4(b))
を見比べると, 図 1(a) では回転子番号が 1,2,3... と増えるにつれて, 粒子の位相が先行していることが分かる. つまり, 流速に一番大きな影響を及ぼす回転子 1 の左向きの流れを回転子 2,3,4... が妨げる方向に運動している. 一方図 4(b) では回転子番号が増えるにつれて粒子の位相が後行している. つまり, 回転子 1 の左向きの流れを回転子 2,3,4... が後押しする方向に運動している. 実際に図 4(a) を見ると, 左向きの流速の値はA > 0(緑) よりも A < 0(青) の方が大きくなっている.
つまり, 繊毛波の流れる向きによって流速は変化し,また正の方向とは逆向きの繊毛波が流れているときに流速が大きくなると考えられる.

876543210250tt(a)A = −0.5 のとき(b)A = 0.0 のとき

sin φ8765432140

本研究では繊毛集団を表すモデルを作成し, 壁を導入したり駆動力のパラメータを変えたりすることによって現れるダイナミクスや流速の変化について考察した.
モデルは A の値によって 3 種類のダイナミクスを示し, また流速は θ = 0, A = −1 程度で最大となることが分かった. しかし壁 2 つのモデル (ナノポンプ) について, シミュレーション前は流速が一番速くなると予想していたが, 実際はかなり不安定な値を取るモデルとなってしまった. この原因がまだ分かっていないので, ポンプ内の流れ場を見るなどして考察を深めていく必要がある.
今後は回転子集団が作る流れの輸送効率を求め, 効率が最大となるような回転子の条件を見つけたい. ま
た今回はモデルに自由壁の効果を取り入れたので, 固定壁の効果を取り入れたモデルを作成し自由壁との比較を行いたい.

250t(c)A = 0.5 のとき
図 3: sinφ の時空プロット.θ = π4 で, 壁が 2 つのモデルの値を使用.

A < 0 のとき, 繊毛波を後ろから前の方向に伝えるような挙動を示した.(図 3(a)) 次に A が 0 付近のとき, 時間がいくら経っても回転子は同期しなかった.(図 3(b))
最後に A > 0 のとき, 繊毛波を前から後ろの方向に伝えるような挙動を示した.(図 3(c))

V

いま流速 V は左向きを正としているので, 前→後ろの方向に波が出来るとき (A > 0 のとき)hV i の値が大きいと予想できるが, 実際は後ろ→前 (A < 0 のとき)の方が hV i の値が大きいという結果になった. この原因を探るため, A > 0, A < 0 の場合についての流速の変化及びそのときのモデルの挙動を細かく見てみることにした (図 4).

参考文献[1] 内田就也. 流体力学的同期現象: 鞭毛や繊毛の集団運動の理解に向けて (最近の研究から). 日本物理学会誌, Vol. 67, No. 11, pp. 754–757, 2012.

21.5
10.5
0-0.5-1-1.5
100 110 120 130 140 150 160[2] Nariya Uchida and Ramin Golestanian. Genericconditions for hydrodynamic synchronization.
Physical Review Letters, Vol. 106, No. 5, p.
058104, 2011.

t(a) 流速の変化

まとめと今後の課題[3] 低レイノルズ数における流体力学的同期現象: ミ
ニマルアプローチ (生物流体力学及び関連する問題の研究).

(b)A = −0.2, t = 2500

図 4: (a) 流速の変化. 赤が V = 0, 緑が A = 0.2, 青が A = −0.2 で, 壁の無いモデルの値を使用. (b)A =[4] JR Blake and AT Chwang. Fundamental singularities of viscous flow. Journal of Engineering
Mathematics, Vol. 8, No. 1, pp. 23–29, 1974.

−0.2 のときのスナップショット.

60