仲 井 賞 若 手 奨 励 賞 受 賞 論 文

流動層造粒における粒子成長の解析とそのモデリング
Analysis and Modeling of Granule Growth in Fluidized Bed Granulation

林 健太朗
Kentaro HAYASHI

シオノギファーマ株式会社

技術開発本部 生産技術部 製剤部門

Formulation Development, Production Technology Department, Technology Development Division,
Shionogi Pharma Co., Ltd.

要旨
本研究では、流動層造粒法の生産スケールでの実験結果、離散要素法（DEM）と数値流体力学（CFD）の連成による粒子運動解析、さらに新規に開発したポピュレーションバランスモデル（PBM）により、粒子の凝集・解砕現象とプ
ロセスパラメータとの関係を解明することを目的に基礎的研究を行った。PAT ツールによる粒子径測定結果から、造粒の急激な成長開始点およびその終了点が試料の可塑限界値に一致することを見出した。さらに、流動層造粒・乾燥工程中の粒子径推移を解明するために、DEM-CFD 連成解析から導出した粒子の衝突頻度関数を実装することにより、プロセスパラメータを考慮した新規 PBM を提案した。本 PBM により算出した粒度分布は、実際のプロ
セスと同程度に長い時間範囲において、生産スケールの実測値と極めて高い精度で一致した。さらに、本 PBM により、異なる運転条件においても粒度分布の時間変化を予測することができた。

Abstract
In this study, a novel population balance model (PBM) included process parameters of fluidized bed
granulation was proposed based on the experimental results at the production scale and the particle
collision analysis with discrete element method and computational fluid dynamics (DEM-CFD) coupling inorder to investigate the relationship between granule growth phenomenon and process parameters. It
was clarified that the start and end points of the rapid granule growth in the fluidized bed granulation
coincided the moisture content at the plastic limit. A novel PBM taking process parameters into
consideration was proposed by implementing the collision frequency function derived from DEM-CFD
coupling to reveal the changes in particle size during fluidize bed granulation including drying process.
The particle size distributions calculated using the PBM agreed extremely well with the actual values
at the actual production and process time scales. This PBM can also predict the changes in particle size
distribution under various operating conditions.

運用から、Process Analytical Technology（PAT）を

1 はじめに活用した科学的な工程理解に基づく運用へとアプ
近年の医薬品製剤開発において浸透してきたローチを転換することが重要である 2）。PAT ツール

Quality by Design（QbD）1）、あるいは、革新的製造の進歩により、製造プロセスの可視化が飛躍的に発展技術として注目されている連続生産を実現するためするとともに、製造中のリアルタイム制御と品質保証には、従来のベテラン作業者による勘と経験に基づくが可能となった。しかし、製造工程の科学的理解を深

■ Keywords ■ process analytical technology, population balance model, DEM-CFD coupling, fluidized bed granulation,
moisture content

13

製剤機械技術学会誌 (377)めるためには、高精度の測定機器を用いて工程をモニ

子径推移を理解し制御することは、品質リスクを軽減

タリングするだけでは不十分であり、PAT ツールにする上で極めて重要である。
医薬品製造で汎用される流動層造粒法では、精製水より取得した工程データを活用・解析する手法の開発

を結合液として用いた場合、最終製品の組成を変えるが求められている。
流動層造粒法は、他の造粒法と比較して良質な造粒

ことなく造粒の進行状況に合わせて加水量を任意に物が得られるだけでなく、原料の混合や造粒、乾燥な調整できるため、医薬品製造のレギュレーションの観ど複数の工程を同一の密閉容器内で行うことができ

点から好都合であることに加えて、PAT ツールによるため、工程数削減、交叉汚染防止ならびに GMP のるオンラインモニタリングおよびプロセス制御との観点からも医薬品製造に適している。流動層造粒工程

相性も良い。これまでにも湿式造粒におけるモニタリの製造方法を設定する際においても、工程理解を深め

ングおよびその終点管理に関する研究は多数報告さるためには、PAT ツールにより収集したデータを数れている 4-10）。湿式造粒においては、乾燥が必須であ

値解析・モデリングすることが必要不可欠である。しるため、造粒だけでなく乾燥までを一連の工程としてかしながら、造粒は多くの素過程からなるため、装置

粒子径推移を理解すべきであるが、流動層溶融造粒な内で起こっている現象は極めて複雑であり、特に生産ど数例を除き 10-12）、乾燥工程における粒子径変化を無スケールでの設備において、プロセスパラメータが粒視していることがほとんどであり、造粒から乾燥終点

子径推移に与える影響や粒子成長の現象そのものにまでを対象とした粒子径推移は十分に解明されていついて十分に理解されているとは言い難く、数値解析ない。

による検討もほとんどなされていない。

そこで本節では、水を結合液として用いる流動層造

そこで本研究では、流動層造粒法の生産スケールで粒工程において、生産スケールで画像解析装置によりの 実 験 結 果、 離 散 要 素 法（DEM）と 数 値 流 体 力 学

粒子径をオンラインモニタリングするとともに加水（CFD）の連成による粒子運動解析、さらに新規に開終点を制御した。さらに、造粒から乾燥までの含水率発したポピュレーションバランスモデル（PBM）によの時間変化が粒子成長に与える影響を解明するためり、粒子の凝集・解砕現象とプロセスパラメータとのに、取得した時系列データに基づき、流動層造粒機構

関係を解明することを目的に基礎的研究を行った。なを考察した 13）。

お本稿では、筆者の学位論文を基に整理・再構成して要点を報告する。

2.1

実験方法

2.1.1

試料

本研究に用いた試料は、工程を理解しやすい系であ2 オンライン粒子測定装置を用いた流動層造粒ることを理由に、水にほとんど溶けない無機化合物を機構の解析主成分とする顆粒剤 14）を模して調合した。その組成
医薬品製剤の製造における造粒の目的は、顆粒剤とおよび物理混合物としての試料物性を Table 1 およ

錠剤では異なり、前者では付着性、飛散性の防止、含び Table 2 にそれぞれ示す。

量均一性の向上、服用のしやすさなどを主な目的とす

次に、固・液・気の充填状態をトルクの推移としてるのに対し、後者では流動性の向上、圧縮性の向上を測 定 で き る ミ キ サ ー ト ル ク レ オ メ ー タ（MTR 2,

3）

主な目的としている 。特に顆粒剤の場合、造粒物のCaleva Process Solutions Ltd.）を用いて、最大のトル

粒子径は、直接、大きさや形状などの外観を決定する

ク値を示すときの乾量基準の含水率を可塑限界値とだけでなく、その後工程の充填・包装工程での製造性して求めた 15）。試料 40 g に対し、0.8 g または 1.6 g

や服用性にも大きな影響を与える。そのため、その粒の精製水を順次添加して測定したトルクを Fig. 1 にTable 1

Formulation of powder samples.

Component(378) Vol.29 No.4（2020）

FunctionContent［wt%］Dried Aluminum hydroxide gel

Active ingredient

41

Magnesium oxide

Active ingredient21

Hydroxypropylcellulose

Binder10

OthersFiller2814示す。実際の医薬品製造プロセスで生じる変動を想定

Properties of powder samples asphysical mixture.

Table 2Mass median diameter［µm］

して、流動層造粒機の初期缶体温度のみを変えること
により、異なる 3 種の加水速度条件（ⅰ；標準条件，ⅱ；

11

Sauter mean diameter［µm］

3.5

緩徐条件，ⅲ；急速条件）にて造粒を実施した。各条

Particle density［g/cm ］

2.3

件での吸気風速の実測値の推移を Fig. 2 に示す。加0.50

水終点を図中に下向き矢印で示した。流動層造粒工程

3

Bulk density［g/cm ］
3

Plastic limit on dry basis［%］中 の 粒 子 径 は 画 像 解 析 装 置（Image Eye, Dalton

39

Corp.）を用いてオンライン測定した。CCD カメラへの粉付着による測定精度低下を避けるために、造粒開2始約 25 分後に画像解析装置のプローブを取り付け、

1.5

加水終点を制御した。また工程中の顆粒の含水率は適Torque [Nꞏm]

Figure legend

質量中位径をオンラインモニタリングするとともに時、赤外線水分計（FD-720, Kett Electric Laboratory,

180 ℃）で測定した。

0.5
0

2.2
0

含水率の時間変化が質量中位径に与える影響

生産スケールにて異なる加水速度で造粒操作を

5 10 15 20 25 30 35 40 45 50
Dry-basis moisture content [%]行った際の、質量中位径と含水率の時間変化を Fig. 3に示す。加水停止時点を図中に矢印で示した。造粒工Gtorqueof powder
samples.
 ranulation
Fig.Fig.
1 1 Granulation
torqueof powder
samples.

程中の含水率を変化させた各条件において、質量中位径の変化と含水率との間に以下に示す特徴的な共通3(i)示す。本研究に用いた試料は、乾量基準含水率が約
Inlet air velocity [m/s]点が認められた。まず、粒子径測定を開始して以降、

39 % の時に最大トルク、すなわち可塑限界値を示し2た。この可塑限界値において、粒子群内の空間に占め

試料の可塑限界値相当の含水率約 39 % を閾値として、る液体の体積比により定義される飽和度は 1 となり、または緩やかであった（段階Ⅰ）。そして、含水率が1
これ以上の含水率では、一般に固体粒子が液体中に懸閾値を超えた時点で質量中位径は顕著に増大した（段濁するか、またはその充填物から液体が湧出する 16）。

階Ⅱ）
。この質量中位径の増加は加水を停止しても継これを超えるまで質量中位径の変化がほとんどない、続した。さらに、乾燥が進行し、含水率が閾値以下と0

2 . 1 . 2 3 実験装置および実験方法なったとき、質量中位径が最大値に達した。その後、Inlet air velocity [m/s](ii)
本研究では、生産スケールの流動層造粒機（NFLO-

質量中位径の減少が進んだ（段階Ⅲ）。このように、

300, Freund Corp.）に試料 195 kg を入れ、精製水を2
粉体層の上方より噴霧した。造粒条件を Table 3 に1

Table 3

Spraying conditions

Inlet air velocity [m/s]

3

2
Step

1

0

The number of spray nozzles［-］9Nozzle insert（i.d.）
［mm］
(iii)

2.2

Spray air pressure［Pa］

5 × 105

Water feed rate per nozzle［kg/min］

0.1 - 0.4

Inlet airtemperature［K］

Inlet airvelocity［m/s］

Water feed rate perpowder mass［（kg/min）/kg］

Step transfer conditions

Blend

298

0.6

NA

For 1 minSpray-1

298

0.6 to 0.9

0.02

For 25 min

20
Spray-2
DryingFig. 2を（Ⅰ）〜（Ⅲ）の三段階に分類することができた。

Operating conditions at commercial scale.

00

可塑限界値相当の閾値と加水速度に基づき、造粒工程

40 29860
80 0.9100
to 1.6120
Time,tot [min]298
368 1.8 to 1.3

0.01

When MMDS reached at 610 µm

NAWhen outlet air temperature reached at 321 K

Note：Fluidized bed granulator was used after warmed to ⅰ）319 K, ⅱ）337 K, ⅲ）314 K

Changes in inlet air velocity in (i) standard, (ii) slow, and (iii) rapid conditions.
15

製剤機械技術学会誌 (379)0Fig. 1

5 10 15 20 25 30 35 40 45 50
Dry-basis moisture content [%]Granulation torque of powder samples.

Inlet air velocity [m/s]3

各造粒条件で得られた造粒物の粒度分布を Fig. 4(i)に示す。粒子径をオンラインモニタリングして加水終点を制御することにより、加水速度を変化させたいず

2れの条件においても、同等の物性を有する顆粒が得ら
れた。

1

2.3

0

3
Inlet air velocity [m/s]

流動層造粒機構の考察

各実験条件において、総加水量はあまり重要でな
く、かつ、加水終点と造粒成長の停止は必ずしも一致(ii)

しなかった。造粒成長の停止点は、含水率が可塑限界

2値以下となった時間に一致しており、含水率と可塑限界値の大小関係が粒子径推移に大きく影響を与えていることが明らかとなった。

12 . 2 節で得られた結果は、Wauters らが提案した
0
100

(iii)

% Cumulative undersize

Inlet air velocity [m/s]

3

2

1

0

0

20

40
60
80
Time, t [min]100, i) Standard, ii) Slow
, iii) Rapid

80
60
40
20
0
100

120

1000
Particle size [m]Fig. 4

Cumulative particle size distribution of granules produced at commerCumulativeFig.
2 2 Changes
in in
inlet
airairvelocityin (i)
standard, (ii) slow, andFig.
(iii)4rapid
conditions.particle size distribution
inletvelocityin（ⅰ）
Fig.
Changesof granules produced at commercial
scale by conventional sieving test

Dry-basis moisture content [%] Mass median diameter [m]standard,（ⅱ）slow, and（ⅲ）rapid
conditions.

1000

(i)

(ii)

I

II

III

(iii)

I

II

III

I

II

III

800
600
400
200
0
50
40

Wthr

Wthr

Wthr

30
20
10
0

0

20

40
60
80
Time [min]

100

120

0

20

40
60
80
Time [min]

100

120

0

20

40
60
80
Time [min]

100

120

Fig.Fig.
3 3Changes
in mass median
diameterand dry-basis
moisturecontentin (i) standard,(ii) slow, content
and (iii) rapid
conditions
at commercial scale
in mass
mediandiameteranddry-basismoisture
in（ⅰ ）standard,（ⅱ (the） threshold of moistureChanges
content, Wthr is 38.9%)

slow, and（ⅲ）rapid conditions at commercial scale（the threshold of moisture content, W thr
is 38.9 %）（380） Vol.29 No.4（2020）16

Fig. 5

Schematic of PBM with DEM-CFD coupling.

攪拌造粒工程での粒子成長機構 17）を、流動層造粒工

凝集が急激に進行し、また、乾燥により生じた顆粒内

程においても適用できることを示唆する結果となっ部の空隙を起点として解砕が進行することが示唆さた。 2 . 2 節で分類した各段階は、（Ⅰ）核形成およびれた。⑷ オンラインモニタリングした粒子径に基づ誘導段階、（（Ⅱ）粒子成長段階、（Ⅲ）解砕段階と名付

いて、加水終点を制御することで、造粒物の粒子径をけて以下のとおり説明することができる。

一定範囲に管理可能であることを示した。

（Ⅰ）核形成および誘導段階；予め試料粉体に混合しておいた結合剤が水に濡れて結合性を発現するこ

3 DEM-CFD 連成解析を用いた新規 PBM によとにより、顆粒の核となる粒子（本研究では、これをる流動層造粒工程中の粒子成長モデリング

一次粒子とみなす）が形成される。本研究の試料については、この段階で質量中位径は約 440 µm まで到達

第 2 節において、PAT を活用して粒子径をオンラすると考えられ、これは実験条件に依存しなかった。

インモニタリングしながら、加水量を制御することにつづいて含水率が閾値を超えるまでの間、噴霧によりより、粒度を一定範囲に管理できることが示された。

粉体の含水率は増加するが、粒子径はほとんど変化しさらに、可塑限界値と流動層造粒機構の関係が解明さない。流動層造粒法では攪拌造粒法などと比較して圧れたことにより、これまで以上に工程を理解すること

密がかかりにくいため、粉体の体積変化は無視できるができた。しかし、まだ半経験則的であるため、前述と推察された。

の粒子成長機構に基づいた新規数理モデルの開発を（Ⅱ）粒子成長段階；加水により含水率が可塑限界

通じて、さらなる工程理解を推進する必要がある。

値を超えた場合、系全体の粉体表面には自由水が存在

数値解析のモデルには、粒子離散要素法（Discreteする。この自由水を結合液として利用して、核形成お

Element Method；DEM）、数値流体力学（Computa­

よび誘導段階において形成された粒子同士の凝集成

tional Fluid Dynamics；CFD）、ポピュレーションバ

長、すなわち二次凝集が進行すると考えられる。全体

ラ ン ス モ デ ル（Population Balance Model；PBM）なの含水率が閾値を超えない時間帯であっても、局所的どがある。本研究では、流動層造粒工程における粒子には、粉体表面に自由水が存在する状況は考えられるの挙動のシミュレーションには DEM-CFD 連成解析

が、粉体表面の水が十分ではないため、粒子成長よりを、凝集と解砕に関するモデル化については PBM をも繰り返し衝突による粉体全体への水の再分配が起

それぞれ選定し、両者を組み合わせるアプローチを採こっている可能性が示唆された。

用した。

（Ⅲ）解砕段階；加水を停止しても粒子表面に存在

流動層造粒機内の粒子の挙動のように、固気混相流

する水の結合力に起因する凝集現象は、すぐには停止の運動を記述する数値解析のモデルには、DEM とせず、乾燥が進行して再び閾値以下となった時点で凝CFD の連成解析（DEM-CFD カップリング）がよく知集は停止する。さらに、乾燥が進行すると造粒物内部られている。DEM は Cundall と Strack によって開発が多孔質になるため、壊れやすい顆粒になる。しかし、された方法 18）であり、複雑な相互作用を簡便に導入

流動層造粒においては、粉砕機のように大きな外力をできるため、有力視されている。これと局所体積平均

与えて積極的に粉砕しているわけではないので、一次法に基づく流体の支配方程式 19）を組み合わせ、Tsuji

粒子以下、すなわち、粒子成長段階開始時点の質量中

らは DEM-CFD 法を開発した 20）。計算機の性能向上

位径以下にはならないと考えられる。

にともない、計算可能な粒子数や規模は拡大していったが、それでも計算負荷はまだ大きく、実際の生産ス

2.4

小括

ケールを扱うには計算時間が問題になる。この問題を

流動層造粒工程中における含水率の時間変化が粒解決するためには、計算能力の高い計算機を利用する

子成長に与える影響を調べるために、画像解析装置を方法、または、粗視化技術の利用などが考えられる 21）。

用いて粒子径推移をオンラインモニタリングした。そDEM-PBM および DEM-CFD-PBM 連成の研究は、いの結果、以下の知見を得た。⑴ 造粒が急激に成長す

くつか検討がなされているが、実験結果との比較がな

る開始点は、含水率が可塑限界値を超えた時点と一致

されていないものや、実験室スケールでの実験結果とした。⑵ 造粒成長の終点は加水停止と一致せず、含の比較程度でしかない。実際の生産スケールでの検討

水率が可塑限界値以下となった時点と一致することや、実際の工程と同程度の長時間に及ぶ検討はまだなが明らかとなった。⑶ 水を結合液として用いて流動されていないのが現状である。

層造粒を行う場合、粉体表面に存在する自由水により

PBM の研究例として、Watano らは、Ouchiyama

17

製剤機械技術学会誌 (381)3.1

および Tanaka22-25）の提案した凝集モデルの付着確率

支配方程式に造粒中の粉体の含水率を考慮した PBM を提案し、転動流動層造粒の粒子成長に及ぼす含水率の影響を評価している

粗視化モデルにおける DEM-CFD 連成解析の

3.1.1

26）

固相

。Adetayo らは、回転円筒造粒を対現実的な計算時間で数値解析を実行するために、象としてストークス数と臨界ストークス数を考慮しSakai らの粗視化モデルを導入する 21）。粗視化粒子

27）

d CGM はオリジナル粒子 d ORG に対し、直径比で L 倍大

た PBM を提案している。Biggs らは、流動層溶融

造粒に対し、噴霧後に観測される顆粒解砕を顆粒成長

きく生成した。

の逆過程として考え、負の速度定数をもつ合一関数をたとおり、湿式造粒後の乾燥工程における粒子径変化

𝑑𝑑��� � 𝐿𝐿𝑑𝑑���

を無視していることがほとんどであるため、造粒工程

視化粒子およびオリジナル粒子を表す。この粗視化モから乾燥工程までを対象とした粒子径推移はまだ解デルでは、粗視化粒子とオリジナル粒子群との間に、明されておらず、十分に検討されていない。

質量と運動エネルギーが保存される。

提案している 11）。しかし、これらの研究は、前述し⑴

ここで、下添え字 CGM および ORG はそれぞれ粗そこで本節では、プロセスパラメータを考慮した新流動層造粒機内での粒子の運動を 3 次元 DEM で解規 PBM を構築し、プロセスパラメータと顆粒成長お析するにあたり、粒子は変形しない剛体ではあるが、よび解砕現象の関係を解明することを目的とした。本, i) Standard

𝑑𝑑��� � 𝐿𝐿𝑑𝑑���
接触の際に生じる変形は粒子間の微小な重なり（オー% Cumulative undersize

100, ii) Slow
, iii)究
Rapid
節における研概略を

80

バーラップ）で表現される Soft Sphere model を用い
Fig. 5 に 示 す。 は じ め に、𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � 𝐿𝐿� ��������
𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � �𝐿𝐿� 𝑘𝑘�𝐧��� 𝜹𝜹𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � 𝐿𝐿� 𝜂𝜂�𝐧��� 𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫
𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
𝑑𝑑��� � 𝐿𝐿𝑑𝑑た。接触力として、弾性反発力、粘性減衰力、および60
���
DEM-CFD
連成により流動層内の粒子の運動を模擬� ��������
𝑭𝑭
�
𝐿𝐿
𝑭𝑭
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
接触粒子間の摺動摩擦力を推定するために
Hertz し、粒子衝突頻度とプロセスパラメータの関係を解析
40
𝑑𝑑���
�
𝐿𝐿𝑑𝑑
𝑑𝑑
�
𝐿𝐿𝑑𝑑
(1)
���
���
29-32）
���
Mindlin
（no
slip）
を用いた
した。それにより導出した関係式を用いて、プロセス
。粗視化粒
𝐫𝐫𝐫𝐫
�
�
�
20
�𝐿𝐿�接触モデル
𝑘𝑘�𝐧��� 𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
� 𝐿𝐿𝜂𝜂�𝐧���
𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 , |𝐿𝐿 𝑘𝑘�𝐧��� 𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿 𝜂𝜂�𝐧��� 𝒗𝒗𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
パ ラ メ ー タ を 考 慮 し た 新 規 PBM を 開 発 し た。 本子に作用する接触力（弾性力および粘性減衰力）
は、
��
𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫
0,|𝐿𝐿� 𝑘𝑘�𝐧��� 𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿� 𝜂𝜂�𝐧��� 𝒗𝒗𝐫𝐫𝐫𝐫
�𝜇𝜇�𝐧��� �𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � 𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝐫𝐫𝐫𝐫𝐫𝐫𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
100
1000 2 節で報告
PBM を用いて算出した粒子径推移を、第Particle size [m]

した生産スケールでの造粒乾燥操作における結果

�𝒗𝒗𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 �
法線方向および接線方向についてそれぞれ次式で表

13）

される。

と比較検討し、本モデルによる予測精度および汎用性
𝑑𝑑���
���
Fig.
4
Cumulative particle size distribution of granules produced at commercial
scale�by𝐿𝐿𝑑𝑑
conventional
sieving test

𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � 𝐿𝐿� ��������
𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � �𝐿𝐿� 𝑘𝑘�𝐧��� 𝜹𝜹𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � 𝐿𝐿� 𝜂𝜂�𝐧��� 𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫
𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧を検証した 28）。

��������
𝑭𝑭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿� 𝑭𝑭
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭� �� ��������
�
𝐿𝐿
𝑘𝑘
𝑘𝑘
���
⑵
𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � 𝐿𝐿 𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � �𝐿𝐿���
𝑘𝑘�𝐧��� 𝜹𝜹𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � 𝐿𝐿� 𝜂𝜂�𝐧��� 𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫
𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
𝐫𝐫𝐫𝐫
�
�
�
�𝐿𝐿
𝑘𝑘
𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿� 𝜂𝜂�𝐧��� 𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫��𝐧���𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 , |𝐿𝐿 𝑘𝑘�𝐧��� 𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿 𝜂𝜂�𝐧��� 𝒗𝒗𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭� ��������
𝜂𝜂
�
𝐿𝐿
𝜂𝜂
���
���
𝑭𝑭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿 𝑭𝑭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭� � ���������
𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫���������
��
𝐿𝐿��𝐿𝐿𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧�𝜹𝜹��𝐿𝐿� 𝜂𝜂�𝐧���
𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫
𝑭𝑭𝑭𝑭
� 𝐿𝐿 𝑭𝑭�𝐧���
𝑘𝑘�𝐧���𝜹𝜹𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
�𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝐿𝐿�𝜂𝜂 𝑘𝑘,�𝐧���
𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫(2) 𝒗𝒗𝐫𝐫𝐫𝐫
𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
� �𝐿𝐿𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
� 𝐿𝐿� 𝜂𝜂�𝐧����𝜇𝜇
𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝜇𝜇� �𝑭𝑭
𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 |𝐿𝐿 𝑘𝑘�𝐧���
𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
𝐫𝐫𝐫𝐫𝐫𝐫 �𝐧���
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝜇𝜇
�
�𝐧���
�𝐧���
�𝒗𝒗
𝐫𝐫𝐫𝐫𝐫𝐫𝐫𝐫𝐫𝐫𝐫𝐫����𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 ��𝐿𝐿 𝑘𝑘�𝐧���𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭�𝜂𝜂�𝐧��� 𝒗𝒗𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 , |𝐿𝐿 𝑘𝑘�𝐧��� 𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿 𝜂𝜂�𝐧��� 𝒗𝒗𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 | � 𝜇𝜇�𝐧����𝐿𝐿���������
𝐿𝐿
𝑭𝑭
𝑭𝑭𝑭𝑭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
� 𝐿𝐿�𝑭𝑭��������
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 𝐫𝐫𝐫𝐫𝐫𝐫
��𝒗𝒗,|𝐿𝐿� 𝑘𝑘�𝐧���𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭� � 𝐿𝐿� 𝜂𝜂�𝐧��� 𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫| � 𝜇𝜇�𝐧����𝜇𝜇�𝐧��� �𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧�� 𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
�𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫
��𝐿𝐿 𝑘𝑘�𝐿𝐿� 𝜂𝜂
�𝐫𝐫𝐫𝐫𝐫𝐫 𝐫𝐫𝐫𝐫𝐫𝐫�𝒗𝒗, |𝐿𝐿 𝑘𝑘�𝐧���𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
� 𝐿𝐿� 𝜂𝜂�𝐧��� 𝒗𝒗𝐫𝐫𝐫𝐫
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
�𝐿𝐿 𝑘𝑘�𝐧���𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
� 𝐿𝐿��𝜂𝜂𝜹𝜹�𝐧���
𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫
,
|𝐿𝐿
𝑘𝑘
𝜹𝜹
�
𝐿𝐿
𝜂𝜂
�𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
��𝐧���
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
�𝐧���
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
�𝐧��� 𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
�𝐧��� 𝒗𝒗𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 | � 𝜇𝜇�𝐧���
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝐫𝐫𝐫𝐫𝐫𝐫𝐫𝐫𝐫𝐫𝐫𝐫��𝒗𝒗��(3) 𝐫𝐫𝐫𝐫
𝒗𝒗𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭�𝐫𝐫𝐫𝐫𝐫𝐫 |𝐿𝐿� 𝑘𝑘���𝑭𝑭� |𝐿𝐿𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭,
𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
�𝜇𝜇�𝑭𝑭�𝐧���,�𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
� 𝐫𝐫𝐫𝐫𝐫𝐫 � 𝐿𝐿 𝜂𝜂�𝐧��� 𝒗𝒗𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
�𝜇𝜇�𝐧���
�𝐧���
𝐫𝐫𝐫𝐫𝐫𝐫𝑘𝑘�𝐧���� � 𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿 𝜂𝜂�𝐧��� 𝒗𝒗𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 | � 𝜇𝜇�
𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐿𝐿���𝒗𝒗��������
�𝐧���
𝐫𝐫𝐫𝐫𝐫𝐫 𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
�𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝑘𝑘
𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
�𝑭𝑭
��𝒗𝒗�𝐿𝐿
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
�𝐧��� 𝜹𝜹𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � 𝐿𝐿 𝜂𝜂�𝐧��� 𝒗𝒗𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧
� ��������
� 𝐿𝐿� 𝑘𝑘
𝑘𝑘𝑭𝑭���
���
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿 𝑭𝑭
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭⑶

𝜂𝜂��� ��𝐿𝐿� 𝜂𝜂���
𝐫𝐫𝐫𝐫𝐫𝐫�rel
�𝐿𝐿 𝑘𝑘�𝐧���
𝜹𝜹、
�μ
𝐿𝐿�s、
𝜂𝜂�𝐧���
, |𝐿𝐿
𝑘𝑘�𝐧��� 𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿� 𝜂𝜂�𝐧��� 𝒗𝒗𝐫𝐫𝐫𝐫
ここで、
F 、k
η、
δ ，𝒗𝒗v𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭、上線、および、𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝑘𝑘��� � 𝐿𝐿� 𝑘𝑘���𝜇𝜇��𝐧���
�𝐧���� �n𝜇𝜇および
𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
�
�
下添え字
t
は、力、バネ定数、粘性減衰係数、,|𝐿𝐿 𝑘𝑘�𝐧��� 𝜹𝜹𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭 � 𝐿𝐿 𝜂𝜂�𝐧��� 𝒗𝒗𝐫𝐫𝐫𝐫
𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
𝜂𝜂��� � 𝐿𝐿� 𝜂𝜂��� �𝜇𝜇�𝐧��� �𝑭𝑭𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧𝐧 � �𝒗𝒗𝐫𝐫𝐫𝐫𝐫𝐫�𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭𝐭
�静止摩擦係数、変位、相対速度、平均値，法線方向お� 𝐿𝐿�𝑘𝑘����𝐿𝐿 𝑘𝑘���(4)
𝑘𝑘𝑘𝑘���
���𝜇𝜇�𝐧��� � 𝜇𝜇�𝐧����
よび接線方向をそれぞれ表す。したがって、粗視化粒𝜂𝜂𝜂𝜂���
� 𝐿𝐿�𝜂𝜂����𝐿𝐿 𝜂𝜂���
(5)
���
子のバネ定数、粘性減衰係数、および静止摩擦係数は、

� 𝜇𝜇�𝐧���𝜇𝜇𝜇𝜇�𝐧���� 𝜇𝜇�𝐧����𝐧���

オリジナル粒子の値と以下の関係が成り立つ。

𝑘𝑘��� � 𝐿𝐿� 𝑘𝑘����𝜂𝜂��� � 𝐿𝐿 𝜂𝜂���

𝜇𝜇�𝐧��� � 𝜇𝜇�𝐧���Fig.
5 5 Schematicof PBM
with with
DEM-CFD
coupling.
of PBM
DEM-CFD
coupling.
Fig.
Schematic（382） Vol.29 No.4（2020）

18

⑷

(6)

��������𝝎𝝎
𝐎𝐎𝐎𝐎𝐎𝐎��������
𝝉𝝉𝐂𝐂𝐂𝐂𝐂𝐂 � �� ������𝝉𝝉𝐎𝐎𝐎𝐎𝐎𝐎 � ��� 𝜇𝜇�𝐝��� 𝐹𝐹� �𝜇𝜇�𝐝��� 𝐹𝐹�𝐝��� 𝑅𝑅���𝐝��� 𝑅𝑅��� |𝝎𝝎��������|𝜇𝜇� 𝜇𝜇�𝐝���𝜇𝜇�𝐝���
�𝐝��� � 𝜇𝜇�𝐝���𝑘𝑘��� � 𝐿𝐿� 𝑘𝑘���� 𝐿𝐿
𝐿𝐿��𝜂𝜂𝑘𝑘���𝑘𝑘𝜂𝜂��� �
���

⑸

���

𝜂𝜂
� 𝐿𝐿� 𝜂𝜂���
𝜇𝜇����𝐧��� � 𝜇𝜇�𝐧���𝜇𝜇�𝐧��� � 𝜇𝜇�𝐧���⑹

また、Eqs. ⑵ - ⑶を満たすために、粗視化粒子とオリジナル粒子のせん断弾性率 G は次の関係が必要となる。
�.�

𝐺𝐺��� � � 𝐺𝐺���

⑺

�.�
𝐺𝐺粗視化粒子に作用する流体抗力
F d, CGM は、各流体��� � � 𝐺𝐺���

セルにおける気体・固体間の運動量交換係数ζORG を�.�
� ����
𝐺𝐺
��������
𝑭𝑭���
����𝒗𝒗𝐠𝐠 � 𝒗𝒗
𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝���� 𝐺𝐺���𝐎𝐎𝐎𝐎𝐎𝐎 ���用いて次式で計算する。
�����
�

�����

∑� ��� ����𝐝�����𝐝�
��������
𝑭𝑭𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝 �∑�.�
��������
�𝒗𝒗𝐠𝐠 � 𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎����1
�
𝜀𝜀𝐺𝐺���
� 1 � �� 𝐺𝐺���
����
�

⑻

�
𝜀𝜀𝐺𝐺
� 1 � �.���������
𝑭𝑭���
�𝒗𝒗𝐠𝐠 ��1
𝒗𝒗
𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝���� 𝐺𝐺
𝐎𝐎𝐎𝐎𝐎𝐎 ���
����
�����
����

⑼

����
�����
���
�

����
�����
���

∑ � ���
∑����𝐝�
����𝐝�
���������1
𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎� � ���
𝜀𝜀𝑭𝑭𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝
� 1 �� �� ��� �𝒗𝒗𝐠𝐠 �
�����
�����
∑� ��� ����𝐝�� ����𝐎𝐎𝐎𝐎𝐎𝐎����������∑� ��� ����𝐝�𝜇𝜇ここで、ρ
g、p およびμg は，それぞれ、気体密度、
�𝐝��� � 𝜇𝜇�𝐝���

圧力、気体の粘性係数を表す。ここで、上添え字 T は、��� � � �𝜀𝜀𝒗𝒗 � � 0𝐠𝐠
��
転置テンソルを示す。気体・固体間の運動量交換係数

����
ζ
ORG
� 𝒗𝒗には、高密度流動層内の流体抗力を表すのに汎𝐠𝐠 ��� � � �𝜀𝜀�� 𝒗𝒗𝐠𝐠 𝒗𝒗𝐠𝐠 � � �𝜀𝜀�� � � � �𝜀𝜀𝜀𝜀� ��𝒗𝒗𝐠𝐠 � �𝒗𝒗𝐠𝐠 �� � 𝜀𝜀�� � � 𝜁𝜁��
用されている
Wen-Yu のモデル 34） と Ergun の式 35）
��
��

36）
� � � �𝜀𝜀𝒗𝒗 �� �
0を組み合わせたモデルを用いた。このモデルでは、�� � � � �𝜀𝜀𝒗𝒗𝐠𝐠
𝐠𝐠 � 0��

気体・固体間の相間交換係数ζ
ORG は、ε >0.8 では����� 𝒗𝒗𝐠𝐠 �������� � 𝒗𝒗𝐠𝐠 � � � � �𝜀𝜀�� 𝒗𝒗𝐠𝐠 𝒗𝒗𝐠𝐠 � � �𝜀𝜀�� � � � �𝜀𝜀𝜀𝜀� ��𝒗𝒗𝐠𝐠 � �𝒗𝒗𝐠𝐠 � �� � 𝜀𝜀�� ���𝐠𝐠���𝜀𝜀�� � �𝜀𝜀��� � �𝜀𝜀𝜀𝜀� ��𝒗𝒗𝐠𝐠 � �𝒗𝒗𝐠𝐠 �� � 𝜀𝜀�� ��� � � の式、ε�
𝐠𝐠 𝒗𝒗𝐠𝐠では＿𝒗𝒗0.8��𝜀𝜀𝒗𝒗
�＜0Wen-YuErgun�の式によって得ら
��
��

れる。

����� 𝒗𝒗𝐠𝐠 ���� � � �𝜀𝜀�� 𝒗𝒗𝐠𝐠 𝒗𝒗𝐠𝐠 � � �𝜀𝜀�� � � � �𝜀𝜀𝜀𝜀� ��𝒗𝒗𝐠𝐠 � �𝒗𝒗𝐠𝐠 � �� � 𝜀𝜀�� �
�

𝐶𝐶
� �

��������𝒗𝒗
�������� |𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠 | ��.��
𝜀𝜀
𝐝
����

𝜀𝜀 � 0.�𝜁𝜁��� � ������� ��������� |𝒗𝒗��������𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎𝐠𝐠 |𝐝 𝜀𝜀 � 0.�
150
� � 1.�5
�����

����

⒁

��������𝒗𝒗
�������� |𝒗𝒗
�𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠 | ��.��
��������𝒗𝒗
� 𝐶𝐶 �������� |𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠 | 𝜀𝜀 ��.��𝐝𝜀𝜀𝜀𝜀 �
𝜀𝜀
𝐝
� 0.�
0.�
����
� 𝐶𝐶�
�
�
�
���
�.���
�
𝜁𝜁𝜁𝜁���
流体抵抗係数C
は粒状体間流れの粒子レイノルズ���D�������
��������𝒗𝒗
�� �
|𝒗𝒗
|
�
𝐎𝐎𝐎𝐎𝐎𝐎𝐠𝐠�
0.15�Re
�
𝐝
Re
�
1000
��� � � �1 ������
� 1.�5 ������
�� |𝒗𝒗��������𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠 | 𝐝 𝜀𝜀 � 0.�
��� 150 ����� ��� �𝐶𝐶�Re� p�に依存し、次式となる。
�
1.�5𝐠𝐠| ��.��
𝐝 𝜀𝜀 � 0.�
150
����
����� �� |𝒗𝒗
数��������𝒗𝒗
��������𝐎𝐎𝐎𝐎𝐎𝐎
����
�����
𝐶𝐶
𝜀𝜀
𝐝
𝜀𝜀 � 0.�
0.44,� �
Re
�
1000
�
�

���
𝜁𝜁��� � ������� ��������� |𝒗𝒗��������𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎𝐠𝐠 |��������𝒗𝒗
��� �150
�𝒗𝒗
�
𝐠𝐠 � 1.�5
��� 𝐎𝐎𝐎𝐎𝐎𝐎
𝐝 𝜀𝜀 � 0.�
�
�.���
��
�
Re
�
��
�
���
���� 1000
�.���� 𝐝 Re�� �1��
������������������
𝒗𝒗
� 0.15�Re� �
�
��� �
��� �
𝐎𝐎𝐎𝐎𝐎𝐎 ���
cell
�1
�
0.15�Re
�
�
𝐝
Re
�
1000
∑
∑
��������
���𝐝����𝐝� 𝝎𝝎
�
�
⒂
𝐶𝐶
�� ��
𝐎𝐎𝐎𝐎𝐎𝐎
𝐂𝐂𝐂𝐂𝐂𝐂
ここで、
、v�
は、それ
g、V �
i、V cell、および、��𝝎𝝎��� ��������
𝜀𝜀𝝉𝝉𝐂𝐂𝐂𝐂𝐂𝐂
� 1�����������
𝝉𝝉ε���𝜇𝜇1�𝐝���
𝐹𝐹� �𝐝���𝑅𝑅 N ��������|� �𝜇𝜇�𝐝��� 𝐹𝐹�𝐝���𝐶𝐶��𝑅𝑅��𝐎𝐎𝐎𝐎𝐎𝐎������������� ��� |𝝎𝝎|𝝎𝝎𝐂𝐂𝐂𝐂𝐂𝐂 |
0.44,
Re
�
1000
𝐎𝐎𝐎𝐎𝐎𝐎
�����
�
0.44,
Re� � 1000
i の体積、粒子 i が存在∑ ��� ����𝐝� ぞれ空間率、気体速度、粒子
𝐺𝐺��� ��.���.� 𝐺𝐺� ����.���
��
𝐺𝐺𝑭𝑭���
� ��𝐺𝐺�����
1� ���������𝝎𝝎𝝎𝝎𝐂𝐂𝐂𝐂𝐂𝐂� ���
𝐎𝐎𝐎𝐎𝐎𝐎
� 𝒗𝒗
�1
�
0.15�𝑅𝑅𝑅𝑅��𝐝
𝑅𝑅𝑅𝑅� � 1000��������
��������
�𝒗𝒗
�
�����
��������𝒗𝒗
𝝉𝝉𝐂𝐂𝐂𝐂𝐂𝐂� � ������𝝉𝝉𝐎𝐎𝐎𝐎𝐎𝐎� �𝜇𝜇�𝐝��� 𝐹𝐹�𝐝��� 𝑅𝑅�����
𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝
𝐠𝐠 𝜇𝜇�𝐝���𝐎𝐎𝐎𝐎𝐎𝐎 𝐹𝐹�𝐝������ 𝑅𝑅���� ���� �𝒗𝒗𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠 � �
する流体セルの体積、流体セル中の粒子数を表す。
������ ��
��������𝒗𝒗|𝝎𝝎��������|
|𝝎𝝎
| �𝒗𝒗
���
𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠 �
���
� ��
𝐶𝐶����Re
𝐎𝐎𝐎𝐎𝐎𝐎
𝐂𝐂𝐂𝐂𝐂𝐂
�.�
�
Re
�
⒃
�
0.44, ���𝑅𝑅𝑅𝑅� � 1000�� � � 𝐺𝐺���
�����
�����
粗視化粒子に作用する転がり摩擦トルク
τ は、次式∑�� ��� ����𝐝� �∑� ��� ����𝐝� 𝝎𝝎��������
𝝎𝝎𝐂𝐂𝐂𝐂𝐂𝐂𝐎𝐎𝐎𝐎𝐎𝐎��������𝝉𝝉𝐎𝐎𝐎𝐎𝐎𝐎� �𝜇𝜇�𝐝��� 𝐹𝐹�𝐝��� 𝑅𝑅���
33） � ���𝜇𝜇1
� 𝐹𝐹𝜀𝜀𝝉𝝉𝐂𝐂𝐂𝐂𝐂𝐂
� 1��� ������
�𝐝����𝐝��� 𝑅𝑅��� |𝝎𝝎��������||𝝎𝝎𝐂𝐂𝐂𝐂𝐂𝐂| ��������𝒗𝒗 �
で計算される。
�����𝐎𝐎𝐎𝐎𝐎𝐎
��
�����
����
� ���� �𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠
�
�
𝜇𝜇
(11)
𝜇𝜇
��������
𝑭𝑭�𝐝���
���𝐝��� �𝒗𝒗𝐠𝐠 � 𝒗𝒗
𝑅𝑅𝑅𝑅� �𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝 �
𝐎𝐎𝐎𝐎𝐎𝐎 ����������������
𝑭𝑭𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝𝐝
� �� ��� �𝒗𝒗
�
𝒗𝒗
3
.
1
.
3
モデル粒子
𝐠𝐠
𝐎𝐎𝐎𝐎𝐎𝐎��������������𝝎𝝎𝝎𝝎
𝐎𝐎𝐎𝐎𝐎𝐎������
�������������
����𝑅𝑅��� 𝐂𝐂𝐂𝐂𝐂𝐂
��𝐝���
������
𝝉𝝉𝐂𝐂𝐂𝐂𝐂𝐂 �������𝝉𝝉���� �𝜇𝜇�𝐝��� 𝐹𝐹�𝐝������ �𝑅𝑅���𝐎𝐎𝐎𝐎𝐎𝐎�� �� �𝜇𝜇�𝐝���∑𝐹𝐹
����
(11)
𝜇𝜇
DEM-CFDPBM のためのモデル粒∑𝜇𝜇� ������|𝝎𝝎��������|
|𝝎𝝎
|
�
�𝐝���
�𝐝���
���𝐝�
���𝐝�
𝐎𝐎𝐎𝐎𝐎𝐎
𝐂𝐂𝐂𝐂𝐂𝐂連成解析および
� �
�
��� �
��������
𝐝𝐝𝐝𝐝𝐝𝐝 � � ��� �𝒗𝒗𝐠𝐠 � 𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎��∑� ���
∑�
𝜀𝜀𝝎𝝎
� ���
1�
���𝐝�
���𝐝�𝝎𝝎� 1���������𝐎𝐎𝐎𝐎𝐎𝐎
𝐂𝐂𝐂𝐂𝐂𝐂
�
�
�
1
�
𝜀𝜀
�
1
�
��������
��������
子の作成方法を以下に示す。まず、Hounslow らの方法
𝑅𝑅��� |𝝎𝝎
� �𝜇𝜇�𝐝��� 𝐹𝐹�𝐝����𝐝��� 𝐹𝐹�𝐝��� 𝑅𝑅������������� |𝝎𝝎𝐂𝐂𝐂𝐂𝐂𝐂 |��������|⑽
𝐎𝐎𝐎𝐎𝐎𝐎
�����
�����
37）
� 𝜇𝜇�𝐝���(11)𝜇𝜇�𝐝���∑� ��� ����𝐝�∑� ��� ����𝐝�と同様に、次式に従い、粒子サイズの区分幅を作成
�1�
�1�
�����
�����
����
��������𝝎𝝎𝝎𝝎𝐂𝐂𝐂𝐂𝐂𝐂𝐎𝐎𝐎𝐎𝐎𝐎した。
ω は、転がり摩擦係数、接
ここで、
r 、R
����
��������
������
𝝉𝝉𝐂𝐂𝐂𝐂𝐂𝐂� ��μ𝝉𝝉𝐎𝐎𝐎𝐎𝐎𝐎
�、および
��� 𝜇𝜇�𝐝���𝐹𝐹� �𝜇𝜇�𝐝��� 𝐹𝐹�𝐝���𝑅𝑅 ��𝐝��� 𝑅𝑅��� |𝝎𝝎���|𝝎𝝎𝐂𝐂𝐂𝐂𝐂𝐂|�� ���
��������|𝐎𝐎𝐎𝐎𝐎𝐎�(11)𝜇𝜇�𝐝��� � 𝜇𝜇�𝐝����
触点から粒子中心までの距離、および角速度をそれぞ�����れ表す。したがって、粗視化粒子の転がり摩擦係数は、

����
(11)
� � � �𝜀𝜀𝒗𝒗𝐠𝐠 � � 0
� � 𝝎𝝎𝐂𝐂𝐂𝐂𝐂𝐂⒄
��������𝝎𝝎��
𝐎𝐎𝐎𝐎𝐎𝐎���𝝉𝝉𝐂𝐂𝐂𝐂𝐂𝐂 �� � ������𝝉𝝉𝐎𝐎𝐎𝐎𝐎𝐎 �𝜇𝜇�𝐝��� ��������
𝐹𝐹�𝐝���
� �𝜇𝜇�𝐝�����������𝝎𝝎𝝎𝝎𝐂𝐂𝐂𝐂𝐂𝐂 𝐹𝐹�𝐝����𝑅𝑅���
オリジナル粒子の値を用いて次式で表される。
𝐎𝐎𝐎𝐎𝐎𝐎 𝑅𝑅���� ����������
|𝝎𝝎��������||𝝎𝝎|������
𝝉𝝉��
�
�
𝝉𝝉
�
��
𝜇𝜇
𝐹𝐹
𝑅𝑅
�
�𝜇𝜇
𝐹𝐹
𝑅𝑅
𝐂𝐂𝐂𝐂𝐂𝐂
𝐂𝐂𝐂𝐂𝐂𝐂 𝒗𝒗 � 𝐎𝐎𝐎𝐎𝐎𝐎
�𝐝��� �𝐝��� ��� |𝝎𝝎�𝐝���𝐎𝐎𝐎𝐎𝐎𝐎
�𝐝��� ��� |𝝎𝝎 |������������|𝐎𝐎𝐎𝐎𝐎𝐎
𝐂𝐂𝐂𝐂𝐂𝐂 ��� �𝐠𝐠 � ��𝜀𝜀𝒗𝒗
� �0𝒗𝒗𝐠𝐠 𝒗𝒗𝐠𝐠 � � �𝜀𝜀��
𝐠𝐠���𝜀𝜀�������
�𝝎𝝎��������� � � �𝜀𝜀𝜀𝜀� ��𝒗𝒗𝐠𝐠 � �𝒗𝒗
�� � 𝜀𝜀�� � � 𝜁𝜁��� �𝒗𝒗
𝐠𝐠
𝐎𝐎𝐎𝐎𝐎𝐎 � 𝒗𝒗𝐠𝐠 ����𝐝���
�
𝜇𝜇
(11)
𝜇𝜇
𝝎𝝎
𝐎𝐎𝐎𝐎𝐎𝐎
𝐂𝐂𝐂𝐂𝐂𝐂�𝐝����⑾
𝝉𝝉𝐎𝐎𝐎𝐎𝐎𝐎 � ��� 𝜇𝜇�𝐝�����������
𝐹𝐹��
�
�𝜇𝜇
𝐹𝐹
𝑅𝑅
𝐂𝐂𝐂𝐂𝐂𝐂 � � ������
�𝐝��� 𝑅𝑅��� |𝝎𝝎
�𝐝���
�𝐝���
���
ここで、下添え字のiは粒子群の番号
（
i
=
1,
2,
...,
��������||𝝎𝝎𝐂𝐂𝐂𝐂𝐂𝐂 |𝐎𝐎𝐎𝐎𝐎𝐎����� 𝒗𝒗𝐠𝐠 ���� � ��
�𝜀𝜀𝒗𝒗
��𝜀𝜀�� �0𝒗𝒗𝐠𝐠 𝒗𝒗𝐠𝐠 � � �𝜀𝜀�� � � � �𝜀𝜀𝜀𝜀� ��𝒗𝒗𝐠𝐠 � �𝒗𝒗𝐠𝐠 � �� � 𝜀𝜀�� � � 𝜁𝜁��� �𝒗𝒗
������
�
�
�
𝒗𝒗
�
𝐠𝐠
𝐎𝐎𝐎𝐎𝐎𝐎𝐠𝐠�� ��
i max）、θi は各群の境界サイズである。以下の各変数

��� � 𝜇𝜇�𝐝���＿θi に属する粒子に対して定義した。つまはθi -1< θ＜
(11)
(11)
𝜇𝜇���𝐝��� � 𝜇𝜇�������� � �𝜀𝜀�� 𝒗𝒗𝐠𝐠 𝒗𝒗𝐠𝐠 � � �𝜀𝜀�� � � � �𝜀𝜀𝜀𝜀� ��𝒗𝒗𝐠𝐠 � �𝒗𝒗𝐠𝐠 � �� � 𝜀𝜀�� � � 𝜁𝜁��� �𝒗𝒗
𝒗𝒗𝐠𝐠 ��𝐝���
𝐎𝐎𝐎𝐎𝐎𝐎 �
��
�
�
�
�𝜀𝜀𝒗𝒗
�
�
0
𝐠𝐠
気体の支配方程式については、粗視化モデルを適用
り、各群の最大粒子サイズをその群の代表粒子サイズ3����
1 �. 𝒗𝒗2𝐠𝐠 ��気相𝜇𝜇�𝐝���
𝜇𝜇.�𝐝���
��

(11)となるようにHounslow の方法を修正した。θi はθ1������� � � �𝜀𝜀�� 𝒗𝒗𝐠𝐠 𝒗𝒗𝐠𝐠 � � �𝜀𝜀�� � � � �𝜀𝜀𝜀𝜀� ��𝒗𝒗𝐠𝐠 � �𝒗𝒗𝐠𝐠 � �� � 𝜀𝜀�� � � 𝜁𝜁��� �𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎 � 𝒗𝒗𝐠𝐠 �する場合においてもオリジナル粒子の系と同様とな����� 𝒗𝒗𝐠𝐠 ���る。連続相である気体には、空間で局所平均した連続���������𝒗𝒗
�������� |𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠 |を用いて次式で表せる。

��.��
������𝐝 𝐎𝐎𝐎𝐎𝐎𝐎 � 𝒗𝒗𝐠𝐠 �𝜀𝜀 � 0.�
� �𝜀𝜀�� � � � �𝜀𝜀𝜀𝜀
�𝒗𝒗 � �� � 𝜀𝜀� � � 𝜁𝜁𝜀𝜀����𝒗𝒗
�� � ��𝒗𝒗𝐠𝐠 � 𝐶𝐶の式、および、運動方程式を用いる。非圧縮性気体の� �� 𝐠𝐠� 0 ���� �� ��� ��𝜀𝜀𝒗𝒗𝐠𝐠
𝜁𝜁�����
�������������� |𝒗𝒗��������𝒗𝒗
�𝐎𝐎𝐎𝐎𝐎𝐎𝐠𝐠 |𝜃𝜃� � ���� 𝜃𝜃�⒅��������𝒗𝒗
�������
連続の式および運動方程式は、次式で表される。
�� |𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠 | ��.��
𝐝 𝜀𝜀𝜀𝜀 �
� 0.�
0.�
150
� � 1.�5 𝜀𝜀
𝐶𝐶
𝐝
�
�
��
����� 𝒗𝒗𝐠𝐠 � ��������
����������� � � �𝜀𝜀�� 𝒗𝒗
𝜀𝜀�� � � 𝜁𝜁��� �𝒗𝒗
モデル粒子の形状を球形とすれば、各群の境界粒子径
di は次式𝜁𝜁���
𝐠𝐠 𝒗𝒗𝐠𝐠 � � �𝜀𝜀�� � � � �𝜀𝜀𝜀𝜀� ��𝒗𝒗𝐠𝐠 � �𝒗𝒗𝐠𝐠 �� �𝐎𝐎𝐎𝐎𝐎𝐎 � 𝒗𝒗𝐠𝐠 ��� � �
������ ��������� |𝒗𝒗��������𝒗𝒗 |モデル粒子の形状を球形とすれば、各群の境界粒子����������𝒗𝒗
�������� |𝒗𝒗| ��.�� 𝐎𝐎𝐎𝐎𝐎𝐎 𝐠𝐠 𝐝 𝜀𝜀 � 0.�
�
𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠
�
1.�5
150
���
��0 �
�� � � � �𝜀𝜀𝒗𝒗
⑿
𝐝
𝜀𝜀 � 0.�𝜃𝜃� 𝑑𝑑�𝑑𝑑𝜃𝜃��d�i は次式となる。
����������� � � �𝜀𝜀𝒗𝒗𝐠𝐠�� 𝐶𝐶
�𝐠𝐠�0 ��������� �.���𝜀𝜀 ����径
� �1 �
𝜁𝜁����� � ��
�0.15�𝑅𝑅𝑅𝑅
�
�
𝐝
𝑅𝑅𝑅𝑅
�
1000
����� ��������� |𝒗𝒗��������𝒗𝒗
��𝐎𝐎𝐎𝐎𝐎𝐎 𝐠𝐠 |モデル粒子の形状を球形とすれば、各群の境界粒子径
di は次
�����
本節では、第2 節で報告した生産スケールの標準条件下におけ� 𝒗𝒗𝐠𝐠 ���𝐶𝐶�𝐝 𝜀𝜀 �
0.� �� � 𝒗𝒗𝐠𝐠 � � ��150
� � � �𝜀𝜀𝒗𝒗𝐠𝐠 � � 0����������� � �𝜀𝜀��1.�5
� �𝜀𝜀���� � �𝜀𝜀𝜀𝜀� � ��𝒗𝒗
��������𝒗𝒗
�������
� 𝒗𝒗𝐠𝐠��𝒗𝒗|𝒗𝒗
𝐠𝐠𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠 � �𝒗𝒗𝐠𝐠 �� � 𝜀𝜀�� � � 𝜁𝜁��� �𝒗𝒗
𝐎𝐎𝐎𝐎𝐎𝐎 � 𝒗𝒗𝐠𝐠 �
����
�����
𝐠𝐠 | ��.���� �0.44,
������
������𝜀𝜀�
𝒗𝒗
𝒗𝒗
�
�
�𝜀𝜀��
�
�
�
��𝒗𝒗
�
�𝒗𝒗
��
�
𝜀𝜀�
�
�
𝜁𝜁
�𝒗𝒗
�
𝒗𝒗
�
�𝜀𝜀𝜀𝜀
�.���𝜀𝜀 � 𝑅𝑅𝑅𝑅𝐝𝐠𝐠� � 1000
𝐶𝐶 � 𝐠𝐠 𝐠𝐠
𝜀𝜀 �� 0.� ��� 𝐎𝐎𝐎𝐎𝐎𝐎 𝑑𝑑 �𝐠𝐠
𝐠𝐠 �������� 𝑑𝑑⒆
������ � �� �1�� 0.15�𝑅𝑅𝑅𝑅
��� 𝐝 𝑅𝑅𝑅𝑅� � 1000
��
�
�
𝜁𝜁���
�����𝒗𝒗
�𝐎𝐎𝐎𝐎�
�� 𝒗𝒗𝐠𝐠 � 𝐠𝐠 | ��.��
�
�
𝐶𝐶
�
�����
������
��������𝒗𝒗
�
|𝒗𝒗
|
�
�
�
𝐎𝐎𝐎𝐎𝐎𝐎𝐠𝐠��������𝒗𝒗
��
�
�𝒗𝒗
�
������
� �𝜀𝜀� �𝜀𝜀��𝒗𝒗𝐝𝐠𝐠 𝒗𝒗𝐠𝐠 � � �𝜀𝜀��𝜀𝜀�0.44,� 𝜀𝜀�� � � 𝜁𝜁𝑅𝑅𝑅𝑅
�𝜀𝜀𝜀𝜀⒀
����0.�
𝐎𝐎𝐎𝐎𝐎𝐎�����𝒗𝒗
本節では、第 2 節で報告した生産スケールの標準条件下におけ𝐠𝐠 � ��𝒗𝒗𝐠𝐠�
𝐠𝐠 ��
��� �𝒗𝒗𝐎𝐎𝐎𝐎𝐎𝐎 � 𝒗𝒗𝐝 𝐠𝐠 �𝜀𝜀 � 0.�
1.�5��
𝑅𝑅𝑅𝑅� � ��150� �����
� � 1000
�.���
����
��
�1
�
0.15�𝑅𝑅𝑅𝑅
�
�
𝐝
𝑅𝑅𝑅𝑅
�
1000
������� |𝒗𝒗��������𝒗𝒗|
��
𝐎𝐎𝐎𝐎𝐎𝐎��𝐶𝐶�𝐠𝐠 �1.�5𝐝 𝜀𝜀� ���������𝒗𝒗
����0.�
���� �𝒗𝒗
�𝐎𝐎𝐎𝐎𝐎𝐎𝐠𝐠���� 𝑅𝑅𝑅𝑅��������𝒗𝒗
�������
|𝒗𝒗| ��.���
�
𝐎𝐎𝐎𝐎𝐎𝐎𝐠𝐠�0.44,𝐶𝐶� ��𝜀𝜀 𝑅𝑅𝑅𝑅𝐝� � 1000𝜀𝜀 � 0.�
����� �.���
�
製剤機械技術学会誌 (383)19� �1 �
𝜁𝜁��� � ��
0.15�𝑅𝑅𝑅𝑅�� 𝐝 𝑅𝑅𝑅𝑅
� 1000��
�𝐎𝐎𝐎𝐎𝐎𝐎
������
��������𝒗𝒗
��������𝒗𝒗
�������������
�𝒗𝒗
� |𝒗𝒗
𝐠𝐠 |
𝐎𝐎𝐎𝐎𝐎𝐎 �� 𝐠𝐠 �
�
�
𝐶𝐶
�� �
𝐝 𝜀𝜀 � 0.�
150
𝑅𝑅𝑅𝑅� � 1.�5
�.���
����� 1000
� ���
𝑅𝑅𝑅𝑅
� 𝐝 𝑅𝑅𝑅𝑅� � 1000 0.44, ��������������𝒗𝒗
�������� |𝒗𝒗
�𝐎𝐎𝐎𝐎𝐎𝐎
𝐠𝐠 | ��.��
������� |𝒗𝒗��������𝒗𝒗 |本節では、第 2 節で報告した生産スケールの標準条
件下における造粒実験13）𝜌𝜌� �𝑡𝑡� � 𝜌𝜌� �𝑡𝑡� � �に基づいて DEM-CFD 連成������

�

⒇

������ �������

したがって、各粒子群の質量
MM（）および粒子数
�𝑡𝑡� � 𝜌𝜌� �𝑡𝑡� � �𝜌𝜌�
i it(t)及び粒子数 NNii(t) は、次式で表したがって、各粒子群の質量
�
�����
��
������
������ �
�𝑡𝑡�
�𝑡𝑡
�
𝜌𝜌
�
𝜌𝜌
�
� �𝑡𝑡�（𝜌𝜌t）� 𝜌𝜌�� �𝑡𝑡�� � ������
� は、次式で表せる。
� ���� 𝑄𝑄� �𝑑𝑑��� , 𝑡𝑡��𝑀𝑀� �𝑡𝑡� � 𝑀𝑀�𝑡𝑡��𝑄𝑄�����
� �𝑑𝑑� , 𝑡𝑡��
したがって、各粒子群の質量Mi(t) 及び粒子数 Ni(t) は、次式で表�𝑑𝑑
𝑄𝑄� � , 𝑡𝑡� � �
N
したがって、各粒子群の質量M(t) 及び粒子数
及び粒子数 Nii(t)(t) は、次式で表は、次式で表したがって、各粒子群の質量Mi(t)𝑀𝑀� �𝑡𝑡� � 𝑀𝑀�𝑡𝑡��𝑄𝑄� �𝑑𝑑� , 𝑡𝑡� � 𝑄𝑄� �𝑑𝑑���i, 𝑡𝑡���������𝑡𝑡� �
�𝑑𝑑� ,, 𝑡𝑡�
𝑀𝑀� 𝑄𝑄� �𝑑𝑑��� , 𝑡𝑡��𝑀𝑀�𝑡𝑡��𝑄𝑄
�𝑡𝑡�
�𝑑𝑑
𝑀𝑀
�𝑑𝑑
𝑄𝑄
,�𝑡𝑡�𝜌𝜌𝑀𝑀�𝑡𝑡��𝑄𝑄��𝑡𝑡
� ��� ��������
𝑁𝑁�����𝑡𝑡�� 𝑡𝑡� � 𝑄𝑄� �𝑑𝑑��� , 𝑡𝑡����𝜌𝜌� �𝑡𝑡� � ���� ����� � ����� ���𝑄𝑄��𝑑𝑑�� ,, 𝑡𝑡�
𝑄𝑄�� �𝑑𝑑
𝑡𝑡� ������ ���𝑁𝑁� �𝑡𝑡� ����
したがって、各粒子群の質量Mi(t) 及び粒子数 Ni(t) は、次式で表� ����
��������� ��� �
�𝑡𝑡�
�
𝑁𝑁
�
�
�𝑡𝑡�
� ��� ����
𝑁𝑁�

�
�
� ����� �𝑑𝑑 , 𝑡𝑡� � 𝑄𝑄 �𝑑𝑑
𝑀𝑀� �𝑡𝑡� � ��𝑀𝑀�𝑡𝑡��𝑄𝑄� �� ��� , 𝑡𝑡��

解析を行った。標準条件下での流動層プロセスを通じて個々の粒子径は 2 mm を超えなかったため、計算の簡便さも考慮して、最大境界粒子径 d i max を 2 mm とした。また、標準条件下での流動層プロセスを通じて0.1 mm 未満の篩下積算分布は 1 % を超えなかったため、 最 小 境 界 粒 子 径 d 1 が 約 0.1 mm と な る よ う に＿ d i に属する粒子
i max=14 とした。粒子径 d が d i -1<d ＜を d i で代表することにより、粒子全体を i max 種の粒子群で表現することができる。そのため、第 2 節の画像

�𝑄𝑄� �𝑑𝑑��, ���
𝑡𝑡� �
�����
� t ）および Q（� ここで、
��M（はそれぞれ全粒子の3 d, t）
���

解析装置で記録した粒子径データを、目開きが 21/3 の������ ���
�
質量および体積基準の篩下積算分布を表す。
𝑁𝑁� �𝑡𝑡� �� ���� � �

等比数列である篩の組（篩の数 14、最大 2 mm、最小��� ������ �� � � ���
��� �
�
粗視化倍率
（t ）およ��������L は、粗視化粒子の総数 N CGM������� ��� ������� � ����� ������� ���
びオリジナル粒子の総数
N ORG（t ）を用いて次式で決定

約 0.1 mm）により再解析して、体積基準の篩下積算分
布を求めた。
ここで、含水率によって粒子体積は変化しない、まできる。

たは、含水率による粒子体積の変化は無視できると仮�

定すると、粒子密度ρ（は、乾量基準の含水率 W（t ）
p t）

����

およびそれらの初期値（t = t 0 ）を用いて次式のように推定できる。

���� �� � � ��� ����

粒子数が最少となる時間においても、一定の粒子数（1 × 104 個）以上を十分確保でき、かつ計算負荷を低減できるよう、粗視化倍率 L を 25 とした。
Table 4

The process parameters of fluidized bed granulation and drying at timepoints understandard conditions.

Time,t［min］

Inlet air velocity,
U 0（t ）
［m/s］

Dry-basis moisture content,
W（t ）
［wt%］

Total number of coarse grain model particles,
N CGM（t ）
［-］

31.2

1.1

37.6

3.1 × 105

42.5

1.3

41.6

1.8 × 105

54.0

1.6

43.3

8.6 × 104

67.0

1.8

37.6

6.4 × 104

74.5

1.7

30.9

7.5 × 104

84.7

1.7

20.9

1.2 × 105

94.7

1.6

10.2

1.8 × 105

101.0

1.1

6.0

1.9 × 105

105.3

1.2

5.1

2.2 × 105

Table 5

List of parameters used for model particles and walls.

Particle properties：

Poisson’s ratio［-］
Solids density［kg/m3］
Shear modulus［Pa］

0.25
Use Eq.（20）1e+6Particle-wall and particle-particle interactions：

Coefficient of restitution［-］Coefficient of static friction［-］Coefficient of rolling friction［-］

0.2
0.5
0.01

Wall properties：

Poisson’s ratio［-］
Solids density［kg/m3］
Shear modulus［Pa］

0.29
7.8e+3
7.6e+8（384） Vol.29 No.4（2020）

20

�
�
��

� ���

�
���
� �𝑡𝑡� � ���
�𝑡𝑡�𝑁𝑁�� �𝑡𝑡� ��� ���� 𝐶𝐶���
𝑟𝑟���
���
���
��
����
��

����

本節で DEM-CFD 連成シミュレーションを実施し子の衝突については restricted-in-space system41） をた（i）標準条件下の各時点における、粗視化粒子の総適用して Eq. を導出した。1 粒子の衝突頻度関数 C i

数、乾量基準の含水率および吸気風速を Table 4 に（t ）、および、2 粒子の衝突頻度関数 C i , （
は次式で
j t）

示す。試料の初期含水率は 4.4 % であった。その他、算出することができる。

モデル粒子の物性値を Table 5 に示す。粒子密度を����

除き、既往の研究 38, 39）の値を参考に用いた。
3.1.4

�� ���� ����𝐶𝐶�� �𝑡𝑡� � ����� ���
��

DEM-CFD 連成解析による流動層内の粒子
衝突評価



�
�
��
���

��� �
�� ��
�
��

2018 および ANSYS Fluent Release 19.1（ANSYS, Inc.）�𝑡𝑡� �𝐶𝐶������

を 用 い て、16 GB の RAM を 搭 載 し た Intel Core i7-粒子の挙動が安定し、実際の圧力損失値と計算値が6820HQ CPU プロセッサ（2.7 GHz）
（Intel Corp.）で実同程度となった 1 〜 2 s の 1 s 間の結果を用いて、流行した。本節では、Hertz - Mindlin（no-slip）接触モデ

動層のプロセスパラメータと 1 粒子の衝突頻度関数

ルおよび標準回転摩擦モデルを用いた。粗視化モデル

C（、および、2 粒子の衝突頻度関数 Ci , （
の関係
i t）
j t）

の流体抗力については、これに対応する EDEM Fluent

をそれぞれ解析した。

DEM-CFD 連成解析は、市販シミュレータ EDEM����� ��������� ������� ���
����

カップリングのスケーリング定数を使用した。Table 4
3.2に記載の各時点において、DEM-CFD 連成解析を実行

PBM

第 2 節と同様に、流動層造粒・乾燥工程を、核形成・した。
粒子径以外の物性値が同じである、i =1 から 14 まで

誘導段階、粒子成長段階、解砕段階に分類する。本研の計 14 種類の粗視化モデル粒子を、各時点での粒径

究では、最終的に得られる顆粒の粒度分布に大きく影分布の測定値に基づいて装置内にランダムに生成し、響し、工程管理の必要性が高いと考えられる粒子成長

自由落下させた。その 1 s 後の状態を流動化シミュレー

段階および解砕段階を対象として、プロセスパラメー

ション の た め の 初 期 配 置 とし た。こ の 後、EDEM-

タを考慮した PBM を考える。ポピュレーションバラン

Fluent カップリングにより粒子を 2 s 間流動させた。

ス式（PBE）を構築するために、以下の仮定を用いた。

粒子の挙動が安定した 1 〜 2 s のΔ t = 1 s 間に、群 i⑴ 系は完全混合バッチシステムである。

に属する粒子が他の粒子または容器内壁に衝突する回⑵ 試料の含水率が閾値を超えると、誘導段階から粒数 N i を求めた。さらに、1 〜 2 s のΔ t = 1 s 間に、群C

子成長段階へ遷移する。

i および j に属する 2 粒子が衝突する回数 N i , j も求めた。
C

粒子の衝突頻度は、Tan ら⑶ 凝集項では、粒子 2 個の凝集により粒子 1 個が生

40）

成される。

によって解析されてお⑷ 凝集確率は、試料の含水率と閾値の差（粒子表面り、それらを参考に、群 i に属する 1 つの粒子の衝突
頻度 r i （t ）を次式で定義した。
C��𝑟𝑟�� �𝑡𝑡� � � � 𝐶𝐶� �𝑡𝑡�𝑁𝑁� �𝑡𝑡���

自由水量に相当）に比例する。
⑸ 凝集する各粒子が小さいほど凝集確率は高く、凝集可能な粒子サイズには上限がある。

⑹ 試料の含水率が閾値以下になると、粒子成長段階

��Cは、群 i に属する粒子の衝突頻度関
�
� （
i t）
𝑟𝑟ここで、� �𝑡𝑡� � �� � 𝐶𝐶� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�から解砕段階へ遷移する。
⑺ 破壊項では、粒子 1 個の破壊により粒子 2 個が生数を表す。

� i および j に属する 2 つの粒子の衝突頻
同様に、群
����
�� �����𝑡𝑡�
� 𝐶𝐶��� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�
𝑟𝑟���
C �
��
����
度 r i , j （t ）を次式で定義した。

⑻ 破壊関数は、破壊によって生成する粒子の大きさ

⑼ 選択関数（破壊確率）は、閾値と試料の含水率の差

��
� ���
�
�𝑡𝑡� � ��� � 𝐶𝐶��� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� �
𝑟𝑟���
��
����

成される。
に依存せず、一定である。
（粒子内空隙率に相当）に比例する。

ここで、Ci , j （t ）は、群 i および j に属する 2 粒子

⑽ 破壊される粒子が大きいほど選択関数（破壊確率）の衝突頻度関数を表す。なお、流動層造粒工程においは高く、破壊可能な粒子サイズには下限がある。

��て、2 粒子衝突は空間率の大きな希薄領域ではなく、� ��� �ることが知られており 42）、これまでの PBM 研究 43）

流動層内部は平均して十分に混合された状態であ𝐶𝐶� �𝑡𝑡� � ���
主に濃厚領域で生じると考えられる。そのため、2
粒
��

��

� ��� ������𝐶𝐶� �𝑡𝑡� � ���� �� �
�
�
�𝑡𝑡�
𝐶𝐶���
�
�� ����� ������������� �� �

21

製剤機械技術学会誌 (385)∑����𝑡𝑡�𝑁𝑁� �� �∑���� 𝛽𝛽�,� �𝑡𝑡�𝑁𝑁𝛽𝛽���� ����
� 𝛽𝛽 �𝑡𝑡�𝑁𝑁��� �𝑡𝑡����,���� � �𝛽𝛽�,��,��𝑡𝑡�𝑁𝑁� �𝑡𝑡�� � ∑���������
����𝑡𝑡��,������ � ��������� �� ���,��� ������ �������� ∑���𝛽𝛽�,��� �� � ���� � 𝛽𝛽���,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � ���� �𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � ∑���
����

����
��
�
�� ������������ ��� ��� 𝐷𝐷� � �𝑡𝑡�
�𝑡𝑡�𝑁𝑁i� �𝑡𝑡�及び群2���
粒子の合一関数を表す。
ここで、β𝑆𝑆����
�𝑏𝑏�,� �𝑆𝑆�𝑡𝑡�
�𝑡𝑡�𝑁𝑁
�𝑡𝑡�𝑁𝑁𝛽𝛽 �𝑏𝑏�𝛽𝛽j に属する� � �𝑡𝑡� � �∑�𝐵𝐵𝛽𝛽�,� �𝑡𝑡�𝑁𝑁
� i,��j(t) �は、群� �𝑡𝑡�∑����,���� �𝑡𝑡�� �𝑡𝑡�𝑁𝑁
� �𝑡𝑡� ��∑∑��������� ��� � ���� ���,����� �,�仮定 (7) より、1 粒子が破壊され 2 粒子が生成する際に群 i の質量と同様に、本研究でも仮定⑴を用いる。 3 . 1 . 4 節で述べたとおり、凝集項については restricted-in-spacesystem41） を適用する。凝集と破壊を対象とした完全混 合 バ ッ チ シ ス テ ム の PBE を 数 値 的 に 解 く た め、
3 . 1 . 3 節で DEM-CFD 連成解析のためのモデル粒子を作成した方法と同様に、Eq. ⒄式に従い、粒子サイズの境界値を作成した。群 i について離散化した PBEは次式で表現できる。
��� ���
��� 𝐵𝐵� � �𝑡𝑡� � 𝐷𝐷� � �𝑡𝑡� � 𝐵𝐵� � �𝑡𝑡� � 𝐷𝐷� � �𝑡𝑡�

こ こ で、t 、N（、B i A（t ）、D i A（t ）、B i B（t ）、およi t）
び D i B（t ）は、それぞれ時間、群 i に属する粒子の個数、離散型の合一生成速度、合一消滅速度、破壊生成

式に書き換える。
関数を表す。

��������𝑡𝑡� � 𝐷𝐷� ��𝑡𝑡� � �𝑏𝑏���,�𝐵𝐵
�𝑡𝑡�𝑁𝑁��� �𝑡𝑡� � ∑������
b i𝑆𝑆, ������ 𝑏𝑏�,� 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑
j は次式で定義される。
𝐵𝐵��仮定⑻より、破壊関数
����� �𝑡𝑡� � 𝐷𝐷�� �𝑡𝑡� � �𝑏𝑏���,� 𝑆𝑆��� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡� � ∑� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑������� 𝑏𝑏�,� 𝑆𝑆
𝐵𝐵� �𝑡𝑡� � 𝐷𝐷� �𝑡𝑡� � �𝑏𝑏���,� 𝑆𝑆��� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡� � ∑�����
𝑏𝑏�,� 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑���𝑏𝑏�,� �
�
���
j に属す
����
�
� i, j(t) は、群 i に属する粒子が破壊されて新たに群��� �
𝐵𝐵ここで、b� �𝑡𝑡� � 𝐷𝐷� �𝑡𝑡� � �𝑏𝑏���,� 𝑆𝑆��� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡� � ∑����� 𝑏𝑏�,� 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑��� � � � 𝑏𝑏�,� � 𝑏𝑏�,��� � 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡
�関数を表す。Si,(t) は、群 i に属する粒子の選択関数を表す。
40）

衝突頻度について Eq. で定義したとおり、Tan ら
仮定 (8)
� より、破壊関数 bi, j は次式で定義される。
と同様に、群i および j に属する 2 粒子の凝集頻度 r i ,�
𝑏𝑏𝑏𝑏�,� �� ������� ���
𝑏𝑏��,�� �����,�
AC
�,��𝑡𝑡����（t）
および衝突頻度
r i , �j�𝑡𝑡�（t ）から、合一関数βi , （、𝑟𝑟�
𝛽𝛽�,� �𝑡𝑡�𝑁𝑁
�
j �,�
j t）
��
����
衝突頻度について
Eq.
(28)で定義したとおり、Tan
衝突頻度関数 C i , （および合一効率関数ψi , （の関 ら [40] と同様
j t）
j t）
�

C
� ���
𝑏𝑏jA�,�(t)� 及び衝突頻度
��,�r���
ri,係を導ける。
合一関数 βi, j(t)、衝突頻度関数 Ci, j(t)i, j (t) から、
��� � �,�𝜓𝜓�,� �𝑡𝑡�
� ���
�
�

��,�
��
�,� ���
���� �������,���,� � 𝛽𝛽 �𝑡𝑡�𝑁𝑁 �𝑡𝑡���
� ���
�𝑡𝑡�
𝑟𝑟𝑟𝑟�,�
�
�
�,�
� ���
� �𝑡𝑡� � �� � 𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� ����𝑟𝑟�,��,��𝑡𝑡�𝑁𝑁� ��𝑡𝑡� �����,� �𝑡𝑡� � ���� � 𝛽𝛽�,�
����
� ���
��,�
�
�

��,� ���

���������,����
���,����
�𝑡𝑡� �𝜓𝜓�� ��� �𝜓𝜓�,�� ��,��,� ���� ��,�
�,� �𝑡𝑡� �
�����
���
���
�� ���
���,�




�,�
�,�
�,�
�,�
�
��,�
∑����
�𝑡𝑡� ��� ���
𝑁𝑁�𝑡𝑡�
�
𝑁𝑁
(30)
�
�
���
��
�
�
�
�
速度、破壊消滅速度を表す。
𝑟𝑟�,� �𝑡𝑡� � � 𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡���� ��� � 𝐵𝐵�� �𝑡𝑡� � 𝐷𝐷�� �𝑡𝑡� � 𝐵𝐵�� �𝑡𝑡� � 𝐷𝐷�� �𝑡𝑡�
��
����
�𝑡𝑡�
�𝑡𝑡�𝐶𝐶
�𝑡𝑡�
�
𝜓𝜓
𝛽𝛽

�� � 𝐵𝐵� �𝑡𝑡� � 𝐷𝐷� �𝑡𝑡� � 𝐵𝐵� �𝑡𝑡� � 𝐷𝐷� �𝑡𝑡�
�,�
�,�
�,�
0,
𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓���
群 i に属する粒子の個数 N（
は、系全体の粒子数
i t）
��
��
A
�
�
�
����� ������
����� ������ � Δt の間に凝𝜓𝜓�i, �j ���
ここで、Nは、群 i 及び j に属する
2 �粒子が時間刻み
�𝑡𝑡� � 𝐵𝐵� �𝑡𝑡� � 𝐷𝐷� �𝑡𝑡��,� �𝑡𝑡�� ���,� ���
�,� �𝑊𝑊�𝑡𝑡� � 𝑓𝑓��� � �N（と以下の関係をもつ。
���
��t�）
�� � ��� � , 𝑓𝑓𝑓𝑓𝑟𝑟 𝑓𝑓�
𝜓𝜓�,� �𝑡𝑡� � � ��� �𝜓𝜓����
����
��
����
�
� 𝐵𝐵� � �𝑡𝑡� � 𝐷𝐷� � �𝑡𝑡� � 𝐵𝐵� � �𝑡𝑡� � 𝐷𝐷� � �𝑡𝑡�� ����� ��を次式で表す。
���� �
��
A
さらに仮定(2)、(4)、(5)
より、合一効率関数
ψ
�,�
i,
j
�,�
�
ここで、
N i , j は、群 i および j に属する 2 粒子が時𝐵𝐵���� �𝑡𝑡� � 𝐷𝐷����
� �𝑡𝑡�
𝑁𝑁�𝑡𝑡�
(30)�� � ∑��� 𝑁𝑁� �𝑡𝑡�
�𝑡𝑡� � 𝐵𝐵� �𝑡𝑡� � 𝐷𝐷� �𝑡𝑡�0,
𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓�
𝛽𝛽�,��𝑡𝑡� � 𝜓𝜓 �,�t �𝑡𝑡�𝐶𝐶
間刻みΔの間に凝集する個数を表す。
0,
𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓�𝑡𝑡�
�,��𝑡𝑡�0,
𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓��������� ��� �������� ���
���� ��
��� �
���
�
�
���
�� �
�
�
∑
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
∑�𝑡𝑡�𝑁𝑁
�𝑡𝑡�
�𝑡𝑡�
����
��
����
��
� � �����(30)� �𝛽𝛽���,���∑��� 𝛽𝛽�,�
�
�𝛽𝛽
�
�
𝛽𝛽
�
�
�
𝜓𝜓
∑
�𝑡𝑡�
𝑁𝑁�𝑡𝑡�
�
𝑁𝑁
����
����
��
��
�
�
�
�
�
�
�
�
�
�
�,�
�
�,�
�
�,�
���
仮定⑶より、2粒子が凝集して1粒子が生成する際
さらに仮定
（2）
、
（4）
、
（5）
より、合一効率関数ψ
,
���
���
�
�
�
�
�
�����
i
j
�,�
�𝑡𝑡�
����
��
����
��
𝜓𝜓�,�
� �𝜓𝜓
�
����
� ���
���� � � � �� ������� � � , �𝑓𝑓𝑓𝑓𝑟𝑟
�𝑊𝑊�𝑡𝑡�𝑁𝑁�𝑡𝑡�� ∑ � 𝑁𝑁� �𝑡𝑡�(30)
𝜓𝜓, 𝑓𝑓𝑓𝑓𝑟𝑟
𝑓𝑓�
𝑓𝑓
��� ���
𝑓𝑓�𝑡𝑡�
�� ����
���������
�𝑊𝑊�𝑡𝑡���� �� � �� �� � , 𝑓𝑓𝑓𝑓𝑟𝑟𝜓𝜓����𝑊𝑊�𝑡𝑡�
𝑓𝑓
�𝑓𝑓
𝑓𝑓���
��� ����
� 𝐵𝐵� �𝑡𝑡� � 𝐷𝐷� � �𝑡𝑡� � 𝐵𝐵の凝集項は次式で表さ
����
����
� �𝑡𝑡� � 𝐷𝐷� �𝑡𝑡�� ���� ����� �� ��
� ��
� ��
����
����
� ��
� ��
����
����
����
����
の質量保存に留意して、Eq.
（t ）を次式で表す。
��
� ��
� ��
���� �
���� �

� ���

��

(30)
��� ��������� ���
���� ��� �
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
∑�𝑡𝑡�𝑁𝑁𝐷𝐷� ����
𝐷𝐷
�𝑡𝑡�
�𝑡𝑡�
�𝑡𝑡�
∑∑𝐵𝐵���
�𝑡𝑡�
𝛽𝛽
�
�
�𝛽𝛽
�
�
𝛽𝛽
�
�
� 𝐵𝐵� � �𝑡𝑡� �れる。
�𝑁𝑁�𝑡𝑡��𝑁𝑁(30)
λa�は、合一効率定数及び合一粒子サイズ定数であ
ここで、ψ
�
�
���,�
�
�,�
�
�,�A 及び
� �𝑡𝑡�
���
���
����
����
��
����
0,
𝑊𝑊�𝑡𝑡� � 𝑊𝑊���
(30)
各時点において同一であり、含水率は時間のみに依存する。この合一��𝐵𝐵 � �𝑡𝑡� � 𝐷𝐷��� �𝑡𝑡�
(30)����������� � � ����� ������ � �𝜓𝜓�,� �𝑡𝑡� � �𝐵𝐵��� �𝑡𝑡� �
𝐷𝐷� �𝑡𝑡�𝜓𝜓� �𝑊𝑊�𝑡𝑡� � 𝑊𝑊���� �
�
�
𝑊𝑊���
��
� � , 𝑊𝑊�𝑡𝑡� �(W(t)
子成長段階においてのみ、含水率と閾値との差- Wthr) に比例し��������������� � ���� ���� �����
���� ��� ��� ������ ���
��� �
����
���
�
����
���
���� ��� ∑�����𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
∑�𝑡𝑡�𝑁𝑁
�𝑡𝑡�
𝛽𝛽
�
�∑
𝛽𝛽
�
�𝛽𝛽
�
�
𝛽𝛽
�
�
�
�
��
�
���
���
���
���,�
�
�,�
�
�,�
�
�,�
�
���
���𝛽𝛽 �𝑡𝑡�𝑁𝑁 �𝑡𝑡� � �𝛽𝛽 �𝑡𝑡�𝑁𝑁 �𝑡𝑡� � ∑ ����� 𝛽𝛽 �𝑡𝑡�𝑁𝑁 �𝑡𝑡��
� 𝑁𝑁∑∑𝐷𝐷���
𝑁𝑁�𝑡𝑡���
(30)
��� � 𝛽𝛽���,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � ���� �∑��������� � �𝑡𝑡��,��なわち、粉体の含水率が閾値を超えたとき、粉体表面に存在する水を
�,�
�
�,�
�
𝐵𝐵 � �𝑡𝑡�
�𝑡𝑡�
�����

� ������������ ���
���� ��
�����𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁
�𝑡𝑡�
���
�
�∑
𝛽𝛽
�
�𝛽𝛽
� ∑�����
𝛽𝛽�,� �𝑡𝑡�𝑁𝑁
表している。さらに、粒子径が最小境界粒子径d1 のとき最も凝集し
���,�
�
�,�
�
�,�
�
� �𝑡𝑡��
���
�𝑡𝑡�� 𝛽𝛽
∑
�𝑡𝑡�
���
�
�
𝑁𝑁�
(30)
����
��� ��� ����������� �������������ここで、ψ
���A およびλa は、合一効率定数および合一���� �� �𝑡𝑡�𝑁𝑁∑����𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�
∑
�𝑡𝑡�𝑁𝑁
�𝑡𝑡�
�
𝛽𝛽
�
�𝛽𝛽
�
�
𝛽𝛽
� ��������� � � ������𝛽𝛽𝛽𝛽
�� ��� �∑
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁
�𝑡𝑡�
∑
�𝑡𝑡�
�
�𝛽𝛽
�
�
����
���,�𝛽𝛽�,� �𝑡𝑡�𝑁𝑁
�
���,�
�
�,�
�
�
�,�
�
���
���
���
���,�
�
�,�
�
�
�����
�
����
��� � ��� � � 𝛽𝛽���,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � ����
����
∑�𝑡𝑡�
���� �𝛽𝛽��それ以上凝集しないように、各粒子径と合一粒子サイズ定数���λa (λa > 0
�,� �𝑡𝑡�𝑁𝑁�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡��
����𝑡𝑡�� 𝑆𝑆���� �∑𝛽𝛽���𝐵𝐵� ������
�𝑡𝑡�
�𝑡𝑡�
∑������� � �𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�∑�����
�∑𝐷𝐷���� ��𝑏𝑏�
𝑏𝑏�,� 𝑆𝑆� �𝑡𝑡�� �𝑡𝑡�𝑁𝑁� 𝛽𝛽� �𝑡𝑡�𝑁𝑁�粒子サイズ定数である。仮定⑴より、全ての粒子の含�,�� ���,� ��� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡� ������� � �� � 𝑏𝑏�,� � 𝑏𝑏�,��� � 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�
(30)
�
���
���
����
��
��� ���
����� ��
����𝑡𝑡�𝑁𝑁�𝑡𝑡�∑�
�𝑡𝑡�� �𝑡𝑡�𝛽𝛽�,� �𝑡𝑡�𝑁𝑁�𝛽𝛽 �𝑡𝑡�𝑁𝑁�∑� 𝛽𝛽 𝛽𝛽���,��𝑡𝑡�𝑁𝑁
�𝑡𝑡�𝑁𝑁
�𝑡𝑡�𝑁𝑁𝐵𝐵���� ����� 𝛽𝛽�,� �𝑡𝑡�𝑁𝑁
𝐷𝐷�𝛽𝛽
�𝑡𝑡�
�𝑡𝑡�
��∑
� �𝑡𝑡�
� �𝑡𝑡�� 0 から 1 の範囲に標準化した。
影響を� �𝑡𝑡��,�� �𝑡𝑡���
������ � ���,���� �,� 𝛽𝛽�,� � �
水率は各時点において同一であり、含水率は時間のみ
����
����
��� ��� �
���
����
�
���� ��
∑��� � � ��� � 𝛽𝛽���,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � � �𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � ∑���
�𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡��
����
����
���� � ������ �� �
破壊現象についても凝集現象と同様に考えると、衝突頻度についてに依存する。この合一効率関数は、含水率が閾値より��� ��� ��������� ���
∑���
�𝑡𝑡�𝑁𝑁� �𝑡𝑡����� ��
�
�𝑡𝑡� � 𝐷𝐷𝛽𝛽��,��𝑡𝑡�����𝑡𝑡�𝑁𝑁� �𝑡𝑡� ��∑��� 𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � ∑�����
𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡��� �� ���∑��� � � � 𝛽𝛽���,����� �𝛽𝛽 �𝑡𝑡�𝑁𝑁��𝑡𝑡� � ∑���� �� 𝛽𝛽�����𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁
�𝑡𝑡�
𝛽𝛽
�
�
�
高くなる粒子成長段階においてのみ粒子は凝集し、含

���,�
�
�,�
�
�,�
�
��
���
����
���� ��� ��� ���� ���� ������� ��
∑��� � � � 𝛽𝛽���,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � � �∑���� ∑�����
𝛽𝛽�,� �𝑡𝑡�𝑁𝑁（W（t�）− W thr）に比例して凝集確率� �𝑡𝑡�
��� 𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�水率と閾値との差
��
����
����
���
���
�
�
�
���� ����� � ����� �𝑡𝑡��
∑�𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
∑�𝑡𝑡�𝑁𝑁
𝛽𝛽
�
�
�𝛽𝛽
�
�
𝛽𝛽
�
�
�
�
�
�
���,�
�
�,�
�
�,�
�
��� ��
���������
が増大することを意味する。すなわち、粉体の含水率,�𝑡𝑡�（i および群
j に属する2 粒子���� ��
j �t ）
𝑏𝑏
� �
����
��𝑏𝑏
���
𝐷𝐷���i�𝑡𝑡�
𝐵𝐵ここで、β�𝑡𝑡�
�𝑡𝑡�𝑁𝑁
�𝑡𝑡�𝑁𝑁� �𝑡𝑡�
� 𝑏𝑏�,��� � 𝑆𝑆��𝑡𝑡�𝑁𝑁
�は、群� ∑������
�,����𝛽𝛽��� 𝑏𝑏�,� 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑������,�𝑆𝑆𝑆𝑆𝛽𝛽����𝑡𝑡�𝑁𝑁�����𝑡𝑡��𝑡𝑡������ ����𝑏𝑏𝑏𝑏�,� �� ������,��𝑡𝑡�𝑁𝑁
𝐷𝐷�����𝑡𝑡�
�𝑡𝑡��,������� 𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�𝐵𝐵�𝑏𝑏�,��� � 𝑆𝑆���𝑡𝑡�𝑁𝑁�∑�𝑏𝑏�����
𝑏𝑏 𝑆𝑆 �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑���
� �
���,�
���
��� �𝑡𝑡� � ∑������,�� �𝑡𝑡�
����� ��� ��� �の合一関数を表す。
�� ���
����� ���,� �が閾値を超えたとき、粉体表面に存在する水を結合剤∑��� � � ��� � 𝛽𝛽���,��𝑡𝑡�𝑁𝑁� �𝑡𝑡� ��𝑡𝑡�
�𝛽𝛽�,� �𝑡𝑡�𝑁𝑁�𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡��� � ∑���
�
�
���
���
���� 𝑆𝑆������𝑏𝑏���,���� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡� � ∑����� 𝑏𝑏�,� 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑��� � �� � 𝑏𝑏�,� � 𝑏𝑏�,��� � 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�として、粒子成長が進行することを表している。さら
仮定⑺より、1
粒子が破壊され 2 粒子が生成する際
��
����
� ��
�
�𝐵𝐵
𝐷𝐷
���
�𝑡𝑡�
�𝑡𝑡�
𝑆𝑆��� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡� � ∑�����
���𝑏𝑏𝑏𝑏�,� 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑���� �� � 𝑏𝑏�,� � 𝑏𝑏�,��� � 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�
�
�
�
���,�
���
∑
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�
�
�
𝛽𝛽
�
�
に、粒子径が最小境界粒子径 d 1 のとき最も凝集しや
の破壊項は次式で表さ�𝑡𝑡� � �∑���𝑏𝑏�,� � 𝑏𝑏�,������,��� � 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡���� � ��の質量保存に留意して、Eq.
��
����
�
�
�𝑡𝑡�
�𝑡𝑡�
�∑���
𝑏𝑏���,� 𝑆𝑆��� �𝑡𝑡�𝑁𝑁���𝑏𝑏�,� 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡���� �� � 𝑏𝑏�,� � 𝑏𝑏�,��� � 𝑆𝑆� �𝑡𝑡�𝑁𝑁�すく、かつ、最大境界粒子径
�,�
� ���
� � ∑��������れる。
d i max のときそれ以上凝
𝑟𝑟�,� �𝑡𝑡� �� 𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�
��
����
�
集しないように、各粒子径と合一粒子サイズ定数λ
a
��
����
𝑏𝑏�,�
� � �
𝐵𝐵
� 𝐷𝐷��� �𝑡𝑡� � �𝑏𝑏���,� 𝑆𝑆��� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡� � ∑�����
𝑏𝑏�,� 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑���
𝑏𝑏�,�
� ���
� �𝑡𝑡���� � �� � 𝑏𝑏�,� � 𝑏𝑏�,��� � 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡���� ��,� ���（λa>0）により、合一効率に与える粒子サイズの影響��,� ���
𝜓𝜓�,� �𝑡𝑡� � � ��� � ���
��
�
��,�
��,�
���
を
1 �の範囲に標準化した。
�𝑡𝑡� � 𝐷𝐷� � �𝑡𝑡� � �𝑏𝑏���,��𝑆𝑆����𝑡𝑡�𝑁𝑁����𝑡𝑡� � ∑������
𝑏𝑏�,���
𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�
𝑏𝑏�,� 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑����,� 0�から��� � �� � 𝑏𝑏
𝑏𝑏�,� ����
破壊現象についても凝集現象と同様に考えると、衝�����𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �∑���� 𝑏𝑏�,� � 𝑏𝑏�,��� � 𝑆𝑆� �𝑡𝑡�𝑁𝑁� �𝑡𝑡���� ���� 𝑏𝑏�,� 𝑆𝑆�
���
���,�
突頻度について Eq. で定義したとおり、群 i に属す� �����𝑡𝑡� � ��,��𝑡𝑡�𝑁𝑁� �𝑡𝑡�������
𝑟𝑟��,�
� 𝛽𝛽�,��𝑡𝑡�𝑁𝑁�𝑡𝑡� �
�𝑡𝑡�
�� � 𝛽𝛽�,�
����
𝑟𝑟�,�
�
る粒子の破壊頻度 r i B（t ）および衝突頻度 r i C（t ）から、
��
����
�
�� ���
𝑏𝑏�,���𝑁𝑁� �𝑡𝑡�
ここで、選択関数 S（b��i , ���（は、群i に属する粒子が破壊されて、衝突頻度関数 C（
および破壊効率
j t）
i t）
i t）
��� ��
������ ���
��,� ���
�,�
�,��������� �𝑡𝑡���,��𝑡𝑡��𝑡𝑡�𝑁𝑁�𝑡𝑡��� ��,�
𝑟𝑟
� 𝛽𝛽
𝜓𝜓
�
�,�
�
�
�,�
�,�
�
0,
𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓�𝑡𝑡�
�
𝑓𝑓
���
����
新たに群
に属する粒子が生成される確率、すなわち関数ψ（の関係を導ける。
���,� ��� � ��,� ���𝜓𝜓�,� �𝑡𝑡� �j ��
i t）
� �
���
��,� �����,����� ������ ���
�𝑡𝑡������ ������ � ������ ������ � �
𝜓𝜓
�
�
�,�
破壊関数を表す。
S（）は、群𝑁𝑁 �𝑡𝑡�i t�𝑓𝑓��� � �i に属する粒子の選択�𝜓𝜓� �𝑊𝑊�𝑡𝑡�� ��� � , 𝑓𝑓𝑓𝑓𝑟𝑟 𝑓𝑓�𝑡𝑡� � 𝑓𝑓���
���
����
��,�
����
��,� ���
����� ������ � � ����� ������ � �𝜓𝜓�,� �𝑡𝑡� � � ��� � ���
��,�
��,�
��
� ���
�
�𝑡𝑡� � �,� � 𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� �
��
𝑟𝑟�,�
����
(386) Vol.29�� No.4（2020）
22
�
��
��,�
�� ���
�� �𝑡𝑡� � �� � 𝛽𝛽�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� ������ ����,� ���
0, � �,�
𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓�𝑡𝑡� � 𝑓𝑓
𝜓𝜓�,� �𝑡𝑡� � �0,
𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓�𝑡𝑡� � 𝑓𝑓 ���
� ���
� ���

��

����� ⁄� ���

�𝑊𝑊��� � 𝑊𝑊�𝑡𝑡�� � �������⁄������� � , 𝑓𝑓𝑓𝑓𝑟𝑟 𝑊𝑊�𝑡𝑡� � 𝑊𝑊���
𝜓𝜓
� �𝑊𝑊
𝜓𝜓
𝑊𝑊�𝑡𝑡��
������� ⁄������ � , 𝑓𝑓𝑓𝑓𝑟𝑟
B
C
�𝑡𝑡�
𝜓𝜓
�
�
�
��� �C�
粒子の破壊頻度 ri (t) 及び衝突頻度 ri (t) から、選択関数 S𝜓𝜓ψi(t)𝑊𝑊�𝑡𝑡�
の関� 𝑊𝑊���
i(t)、衝突頻度関数
i(t) 及び破壊効率関数⁄
���������
� �𝑡𝑡� � ���� �� �
0,
𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓�𝑡𝑡� �
� 𝑊𝑊
𝑊𝑊���
0,
𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓�𝑡𝑡�
係を導ける。
���
�
��

����
�
� � 𝑆𝑆𝑟𝑟𝑟𝑟�����𝑡𝑡�� ���
�𝑡𝑡�𝑁𝑁����𝑡𝑡��𝑡𝑡��������𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �������� ������� ����𝑡𝑡� � � �𝜓𝜓�,� �𝑡𝑡�� ��� ∑���� � ���� � 𝜓𝜓���,� �𝑡𝑡�𝐶𝐶���,� �𝑡𝑡�𝑁𝑁�(38)(38)
���� � ���� ������� ��� ��� ��� ��������� ∑������ � ���
���,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � ���� �𝜓𝜓�,� �𝑡𝑡�𝐶𝐶�,� �𝑡𝑡�𝑁𝑁��� �� � ���� � 𝜓𝜓���,���𝑡𝑡�𝐶𝐶��� ���
���� ��� ��� �
�
��
����
�
��� ��∑ � � ���� 𝜓𝜓���,� �𝑡𝑡�𝐶𝐶���,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� ��𝜓𝜓�,� �𝑡𝑡�𝐶𝐶�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � ∑��� 𝜓𝜓�,� �𝑡𝑡�𝐶𝐶�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡��� ��� ��� ������� �����
����� ��� ��� ��������� ���
����� ��� � ����������� � ��������
�������
���
�𝑡𝑡�
���
𝜓𝜓
�
∑
�𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁
�𝑡𝑡�𝐶𝐶����𝑡𝑡�𝐶𝐶(39)𝜓𝜓��� ∑���
����
�𝑡𝑡�𝐶𝐶�𝜓𝜓 �𝑡𝑡�𝑁𝑁�𝑡𝑡� �
∑������
�
�
��
� 𝜓𝜓�
� 𝜓𝜓��,��𝑡𝑡�𝐶𝐶��������𝑡𝑡����,����,��,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�
�,� �𝑡𝑡�𝑁𝑁
�𝑡𝑡�
���
𝜓𝜓�
� ����������������� � �������
���� �𝑡𝑡�𝑁𝑁
� �𝑡𝑡

��
����
���� ����,�
�
�
� ���
� ���
�
�
�
�

�
���
�
�𝑡𝑡�𝑁𝑁��� �𝑡𝑡� �∑������� � � 𝜓𝜓��� �𝑡𝑡�𝐶𝐶����� 𝜓𝜓� �𝑡𝑡�𝐶𝐶� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � �������� ������ ��� �������
����
��
�
�
�
� �𝑡𝑡�
�������∑�𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁∑
�𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁�𝑡𝑡�𝜓𝜓� 𝑟𝑟𝑆𝑆�����
�
�
�𝜓𝜓
�
�
𝜓𝜓
�
�
�𝑡𝑡�
�𝑡𝑡�𝑁𝑁
�𝑡𝑡�
�
�
(38)
�𝑡𝑡�𝑁𝑁�𝑡𝑡�
∑�𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁�𝑡𝑡�
�𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁 �𝑡𝑡�(40) (42)
𝜓𝜓
𝜓𝜓
𝜓𝜓
�
�
�
�
��
�
�
�
�
�
�
�
�
���,�
���,�
�
�,�
�,�
�
�,�
�,�
�
���
�
�
���
���
���
���
���
�
�
�
�
�𝑡𝑡�
�𝑡𝑡�𝐶𝐶�𝑡𝑡����𝜓𝜓� � � � �
�����
��
����
����
 ��� 最終的に、体積基準の累積粒度分布
�
�
���
q3(
Q� 3(di, �t) �及び頻度分布
�
�
�
����
𝜓𝜓� �𝑡𝑡�𝐶𝐶� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � � � ����� � � � 𝜓𝜓� �𝑡𝑡�𝐶𝐶� �𝑡𝑡�𝑁𝑁� �𝑡𝑡
�� � � 𝜓𝜓��� �𝑡𝑡�𝐶𝐶��� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡� � ∑�����
�
�
�
���
�(di, t) 及び頻度分布
���
q
((d
d
最終的に、体積基準の累積粒度分布
Q
B
3
3
i
i+1
ここで、N
i は、群 i に属する粒子が時間刻み Δt の間に破壊される個数を表す。
∑���� �� ����1/2��
� 最終的に、体積基準の累積粒度分布 Q (d , t) 及び頻度分布 q ((dd),t)
は、それぞれ次式から得られる。
�� ���� �� ���
�
�
�
3 i
3
i
i+1
����
� ����
𝑄𝑄���𝑑𝑑��∑, �𝑡𝑡��𝑡𝑡� �
∑����� � � 𝜓𝜓� �𝑡𝑡�𝐶𝐶�ψ�𝑡𝑡�𝑁𝑁�𝑡𝑡� �
�𝑡𝑡�𝑁𝑁� B ���
(39)�𝑡𝑡�𝜓𝜓さらに仮定
�����
(42)
��
�����
������� 𝜓𝜓�� �𝑡𝑡�𝐶𝐶� � ��𝜓𝜓���� � ���
�� �𝑡𝑡�
� �𝑡𝑡�𝐶𝐶��� �𝑡𝑡�𝑁𝑁���� を次式で表す。
∑����� ����
���より、破壊効率関数
��Q
��q3((di� di+1)1/2, t) は、それぞれ次式から
����3(dt)��及び頻度分布
最終的に、体積基準の累積粒度分布
�i, �
� は、群 i �に属する粒子が時間刻みΔ
� ���
���
���
� ��� N i �(6)、(9)、(10)ここで、
t i(t)
𝑄𝑄� �𝑑𝑑
���
����
�
�
�𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁 �𝑡𝑡� � ��� ��� �𝜓𝜓�,� �𝑡𝑡�
�𝑡𝑡�
𝜓𝜓
�
�
�
� , 𝑡𝑡�� ��� ���� ��� ∑�������∑��� �� ������� ������ � 𝜓𝜓���,� �𝑡𝑡�𝐶𝐶���,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� �
∑�
����
∑
�
�
�
��
����
���� �𝜓𝜓�,�
�
�
�
���,�
���,�
�
� ��� ���� ����
𝑄𝑄� �𝑑𝑑� , 𝑡𝑡� � ������ 1/2
����
����
�
�
�
の間に破壊される個数を表す。
�
∑
����
�
�� ����� ,����� �����∑��� �� ����������
�
� ,��
⁄t)
�
�����及び頻度分布(40)
最終的に、体積基準の累積粒度分布Q�
3(di��,⁄���3((d𝑞𝑞i d��𝑑𝑑i+1) ,𝑑𝑑t) は、それぞれ次式から得られる。
⁄�������� ��� � , 𝑓𝑓𝑓𝑓𝑟𝑟 𝑊𝑊�𝑡𝑡�q�
𝑄𝑄� �𝑑𝑑�� , 𝑡𝑡������
�� ��� , 𝑡𝑡� ������𝜓𝜓𝑊𝑊�𝑡𝑡��
� �𝑊𝑊
��� �
���
��� �
�𝑓𝑓𝑓𝑓𝑟𝑟 𝑊𝑊�𝑡𝑡� � 𝑊𝑊
�
�
���
���
��� �
�� ����� ,��������
,��
�
�𝑊𝑊
𝜓𝜓
�
𝑊𝑊�𝑡𝑡���𝑊𝑊� � ,（
∑
����
�
�
�
���
���
さらに仮定⑹、⑼、⑽より、破壊効率関数ψ
最終的に、体積基準の累積粒度分布
t
）
を
Q（
およ
𝜓𝜓
�
�
�
⁄
�
�����
�
�
i
3� d i, t ）
�𝑡𝑡�
(41)
�
���⁄
����������� � ���𝑞𝑞� ��𝑑𝑑� 𝑑𝑑��� , 𝑡𝑡���� ���� 𝜓𝜓 ����������𝜓𝜓��� �𝑡𝑡����� �����
⁄�
���
��
����
∑
∑����
������
��� ���������� � � � 𝜓𝜓� �𝑡𝑡�𝐶𝐶� �𝑡𝑡�𝑁𝑁� �� ��� ,����� ,������� �𝑡𝑡�𝐶𝐶��� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡�
�
���
�
1/2
�𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁
�𝑡𝑡�
∑
𝜓𝜓
�
��
�
� 𝜓𝜓� �𝑡𝑡�𝐶𝐶� �𝑡𝑡�𝑁𝑁� �
�
�
�
���
��𝑑𝑑
𝑑𝑑
,
𝑡𝑡�
�
�i +1
���
���
���
� ���� 𝑞𝑞0,
𝑄𝑄� �𝑑𝑑� , 𝑡𝑡�次式で表す。
𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓�𝑡𝑡�
�
𝑊𝑊
び頻度分布
q
（
（
d
d
）
,
t
）
は、それぞれ次式により
�
�
���
�����
�
�
���
3
i
�
���
�� �����𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓𝑓�𝑡𝑡�
,����� ��� ,�� � 𝑊𝑊������∑��� ��0,��������� �����𝑞𝑞� ��𝑑𝑑� 𝑑𝑑����, 𝑡𝑡� ��得られる。
��� ��� �t) 及び頻度分布
及び頻度分布 qq3
最終的に、体積基準の累積粒度分布
��
ここで、ψ
B 及び λb は、破壊効率定数及び破壊粒子サイズ定数である。この破壊効率関数は、含水率が閾値以(dii,, t)
最終的に、体積基準の累積粒度分布 Q
Q33(d3�� ����� ,����� ��� ,��������� ⁄�� � �𝑞𝑞� ��𝑑𝑑� 𝑑𝑑��� , 𝑡𝑡� � �𝑊𝑊 ������𝜓𝜓
�𝑊𝑊�𝑡𝑡��
��
,𝑊𝑊�𝑡𝑡�
�
𝑊𝑊∑��� �� ������
���
下となる解砕段階においてのみ、閾値と含水率との差(Wthr - W(t)) に比例して破壊確率が増大することを意味す∑������ �����
(41)𝜓𝜓� �𝑡𝑡� � � � ������ �� � ��������� ⁄�����

�𝑑𝑑
𝑄𝑄

��� � ���� �𝑡𝑡� �� ∑�����
𝑄𝑄�� �𝑑𝑑�� ,, 𝑡𝑡���� �� ������
∑る。すなわち、粉体の含水率が閾値以下となったとき、粉体内部に存在する空隙に比例して、粒子破壊が進行す
0,
𝑊𝑊�𝑡𝑡� � 𝑊𝑊����������� ���� ��,����� ��� ,��
ることを表している。さらに、最大境界粒子径 dimax のとき最も破壊されやすく、かつ、粒子径が最小境界粒子径
�� ,��
𝑞𝑞𝑡𝑡� �� ��������
��� ,����
��� ,, 𝑡𝑡����� � �𝑞𝑞�� ��𝑑𝑑
��𝑑𝑑�� 𝑑𝑑
𝑑𝑑����������
�
�� �
���
ここで、ψ
およびλ
は、破壊効率定数および破d1 のときそれ以上破壊されないように、各粒子径と破壊粒子サイズ定数λb (λb > 0) により、破壊効率に与える粒Bb��

壊粒子サイズ定数である。この破壊効率関数は、含水子サイズの影響を0 から 1 の範囲に標準化した。したがって、Eq. (29) は次式となる。
率が閾値以下となる解砕段階においてのみ粒子は解今回構築した PBE は、Runge-Kutta 法を用いて数して破壊確率が増大することを意味する。すなわち、を Table 6 に示す。各時点における吸気風速、粒子������ ��� ���������
����
�
� ���
��
���
��
�����
��� �𝜓𝜓 �𝑡𝑡�𝐶𝐶 �𝑡𝑡�𝑁𝑁 �𝑡𝑡� � ∑��������� ��� � ������ � �� � 𝜓𝜓�� 𝜓𝜓�,� �𝑡𝑡�𝐶𝐶�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡��
�𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁� ����������� ��� ∑���
���
���
���,�
���,�
�,� �𝑡𝑡�𝐶𝐶�,��,� �𝑡𝑡�𝑁𝑁�� �𝑡𝑡� � ∑���
��� �
���（
���
砕され、閾値と含水率との差
W
W（���,�
t ））に比例の計算に用いたパラメータ�𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁��� �𝑡𝑡�
�𝑡𝑡� �thr −𝜓𝜓���,�
�𝜓𝜓�,�𝜓𝜓�,����� � ����
���� ∑������� 値解析的に解いた。PBE� � �����
��
����
�
���,�
���,�
�,�
�,�
�
�,� �𝑡𝑡�𝐶𝐶�,��,� �𝑡𝑡�𝑁𝑁�� �𝑡𝑡��
���
���
���
��
����
����
�
�

�
��
��
��
������
��������
粉体の含水率が閾値以下となったとき、粉体内部に存�𝑡𝑡�𝑁𝑁 �𝑡𝑡� � ��� � �����
�𝑡𝑡�𝐶𝐶�� �𝑡𝑡�𝑁𝑁
�𝑡𝑡�𝑁𝑁�� �𝑡𝑡�(42)����� � � � 𝜓𝜓�� �𝑡𝑡�𝐶𝐶��密度、ザウター平均径、および粒子総数を用いて衝突���� � � � � 𝜓𝜓�� �𝑡𝑡�𝐶𝐶��� �𝑡𝑡�𝐶𝐶��� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡�������𝑡𝑡�𝐶𝐶����𝑡𝑡�𝑁𝑁����𝑡𝑡� �
∑�����
𝜓𝜓���
�∑
(42)
�� ���� �� 𝜓𝜓������� 𝜓𝜓� �𝑡𝑡�𝐶𝐶� �𝑡𝑡�𝑁𝑁��� �𝑡𝑡� � ��� � ����� 𝜓𝜓�
���
���
���
���
���
�
� �𝑡𝑡�
�����
�
���
�
���
在する空隙に比例して、粒子破壊が進行することを表
頻度関数
C
（
t
）
を算出した。含水率の閾値には、ミキ��� ���� ������� ���
���� ��
�𝑡𝑡�(d� , t)∑����𝑡𝑡�𝑁𝑁� �𝑡𝑡��� ��� ∑����𝜓𝜓�,� �𝑡𝑡�𝐶𝐶�,� �𝑡𝑡�𝑁𝑁� �𝑡𝑡�q�3((d𝜓𝜓1/2�,� �𝑡𝑡�𝐶𝐶�,���� �� � ���� � 𝜓𝜓���,� �𝑡𝑡�𝐶𝐶���,� �𝑡𝑡�𝑁𝑁�Q
最終的に、体積基準の累積粒度分布i+1)1/2������している。さらに、最大境界粒子径
di maxQのとき最もサートルクレオメータで測定した可塑限界値を用い
33(dii,����
t) 及び頻度分布
及び頻度分布q3 ((dii ddi+1)1/2,, t)t) は、それぞれ次式から得られる。
は、それぞれ次式から得られる。
最終的に、体積基準の累積粒度分布33i

破壊されやすく、かつ、粒子径が最小境界粒子径d1∑�����
�������
�������
∑
�
�
�
���
�
𝑄𝑄 ��𝑑𝑑 �, 𝑡𝑡� �� ∑���� �� �����

����

i

i+1

た。流動層造粒工程における粒子成長段階開始時点

�

�

(43)

�

���
����
𝑄𝑄�� �𝑑𝑑������𝑡𝑡�𝐶𝐶�𝑡𝑡�𝑁𝑁�𝑡𝑡� � ∑����� � � 𝜓𝜓� �𝑡𝑡�𝐶𝐶� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� � � t1
� は、含水率が閾値に達した時点の時間とし、それ� � � � 𝜓𝜓� �𝑡𝑡�𝐶𝐶� �𝑡𝑡�𝑁𝑁� �𝑡𝑡� (42)� 𝜓𝜓∑∑∑���
� ��
����
������
� , 𝑡𝑡�����
のときそれ以上破壊されないように、各粒子径と破壊����������� ���� ����
���
���
�
���
���

�

�

����

以降の計算を実施した。初期粒度分布は、標準条件

粒子サイズ定数λ（により、破壊効率に与えb λ b>0）� ��,���� �� ,��� ��� ���,�������������,��
,����� ,����
𝑞𝑞𝑞𝑞����𝑑𝑑� ��������
�𝑑𝑑 ���,, 𝑡𝑡����� �1 の範囲に標準化した。
�る粒子サイズの影響を0���から� ��𝑑𝑑�� 𝑑𝑑���
��� 𝑡𝑡� �
����
���
�����

(44)
ⅰの粒子成長段階開始時点 t 1 における粒度分布を用��� ���� ���
したがって、Eq. は次式となる。

いた。時間刻み幅は、画像解析装置による粒度分布の測定間隔と同じ 10 s とした。時間 t における含水率∑���� �� ������
(43)
𝑄𝑄� �𝑑𝑑� , 𝑡𝑡� � �����
今回構築した∑��� ��PBE
����� は、Runge-Kutta 法を用いて数値解析的に解いた。PBE の計算に用いたパラメータを Table 6に示す。各時点における吸気風速、粒子密度、ザウター平均径、及び粒子総数を用いて衝突頻度関数
C(t) を算� �� ,���� �� ,��𝑞𝑞� ��𝑑𝑑� 𝑑𝑑��� , 𝑡𝑡� � � ��� ���� � �
(44)
���
�
出した。含水率の閾値には、ミキサートルクレオメータで測定した可塑限界値を用いた。流動層造粒工程におけ��parameters used for PBM calculation.
Table 6 List
る粒子成長段階開始時点
t1 of
は、含水率が閾値に達した時点の時間から決定し、それ以降の計算を実施した。初期

粒度分布は、標準条件(i) の粒子成長段階開始時点 t1 における粒度分布を用いた。時間刻み幅は、画像解析装置
Maximum of particle size diameterd i max［mm］
2による粒度分布の測定間隔と同じ
10classs とした。時間 t における含水率i max
W(t)Maximum of particle size［-］ は乾量基準の含水率の実測値を線形
14
Dry-basis moisture content at start under standard condition
補完した値で推定した。

W（t 0）［-］

4.4%

moisture content at start under slow condition
W（t 0）［-］
5.3%
合一効率定数Dry-basis
ψA、合一粒子サイズ定数λa、破壊効率定数 ψB、及び破壊粒子サイズ定数λb は、標準条件 (i) のDry-basis moisture content at start under rapid condition

W（t 0）［-］

Time step for PBM

Δt PBM［s］

10

Threshold of moisture content

W thr［-］

38.9%

ψA［-］

1.1 × 10-4

4.6%

実験値を用いて最小二乗法により決定した。generalized reduced gradient 法により、Microsoft Excel 2016 Solver(Microsoft corp.) を用いて、標準条件 (i) の各時点における累積分布の実験値と計算値の残差平方和 RSS の総和が最小値となるように決定した。
Coalescence efficiency constant

さらに、標準条件(i)sizeの実験値と計算値のフィッティングにより求めた各定数を用いて、Particleconstant for coalescenceλ［0.75 緩徐条件 (ii) 及び急a -］速条件 (iii) においてもPBM の計算を行い、粒度分布の実験値との比較を行った。
緩徐及び急速条件における初
2.9 × 10-4
Breakage efficiency
constant
ψB［-］
Particle size constant for breakage
λ［
2.3
b -］
期粒度分布は、いずれも標準条件の粒子成長段階開始時点
t1 における粒度分布を使用した。

23

製剤機械技術学会誌 (387)W （t ）は乾量基準の含水率の実測値を線形補完したメータと 1 粒子の衝突頻度関数 C（の関係を解析しi t）

値で推定した。

た。粒子径が衝突頻度関数に与える影響を調べた結Microsoft Excel 2016 Solver（Microsoft corp.）を用果、衝突頻度関数は粒子径にほとんど依存しなかった

いて、generalized reduced gradient 法により、標準ため、C（を C（t ）に書き換えて、以降の解析を実i t）

条件 ⅰの各時点における累積分布の実験値と計算値の施した。

残差平方和 RSS の総和が最小となるように、合一効

群 i の粒子数 N（を横軸に、群 i の粒子の衝突頻i t）

率定数ψA、合一粒子サイズ定数λa、破壊効率定数ψB、

度（群 i の粒子が 1 s 間に衝突した回数 N iC ）を縦軸におよび破壊粒子サイズ定数λb をそれぞれ決定した。

プロットした図を Fig. 6 に示す。図の破線は、線形

さらに、標準条件 ⅰの実験値と計算値のフィッティ

近似直線を示し、その傾きを各時点における C（t ）とングにより求めた各定数を用いて、緩徐条件 ⅱおよ

して求めた。Fig. 6 から、いずれの時点においてもび急速条件ⅲにおいても PBM の計算を行い、粒度分

衝突頻度と粒子数の間に比例関係が成立することが布の実験値との比較を行った。緩徐および急速条件に判明した。

おける初期粒度分布は、いずれも標準条件の粒子成長

ここで、流動層内の粒子の衝突頻度を巨視的な運動

段階開始時点 t 1 における粒度分布を使用した。

の大小として捉える。流動層内の全粒子に作用する重力と浮力の差と慣性力の比は、無次元量であるフルー

3.3
3.3.1

流動層内の粒子の運動解析結果

ド数 Fr を利用して次式で表すことができる。

1 粒子衝突頻度関数とプロセスパラメータ��� ��� �� ∑ �� ���の関係�� ��� ∑ �� ���

粒 子 の 挙 動 が 安 定 し た 1 〜 2 s の 1 s 間 の DEM-��� ��� ������� ������ ��� �

�

��

�

∙

���

Linear approximation107
t = 31.2 mint = 42.5 min

t = 54.0 min

105
104
103



ここで、d SV はザウター平均径を表す。巨視的な運CFD 連成解析結果を用いて、流動層のプロセスパラ

106��� ����� �

��

�

�
�
𝐶𝐶��� �𝑡𝑡� � �� ����� �𝑡𝑡�∙�
y 2= 113.4x
y 2=∙ 90.33x��������
R = 0.9906
R = 0.9802

y 2= 122.0x
R = 0.9836

ri, CGMC [s-1]

107
106

t = 67.0 min

t = 74.5 min

t = 84.7 min

105
104
10

3

10

7

106

y 2= 83.97x
R = 0.9718

y 2= 87.42xy 2= 90.92x
R = 0.9715
R =��0.9687
� ������ �∙𝐶𝐶��� �𝑡𝑡� � �� ����� �𝑡𝑡� ∙��

t = 94.7 min

t = 101.0 min

�

��

��

�����

t = 105.3 min

105
104
103

y 2= 96.11x
R = 0.9735

y 2= 109.7x
R = 0.9759

y 2= 109.5x
R = 0.9736

101 102 103 104 105 101 102 103 104 105 101 102 103 104 105

Ni, CGM [-]Fig.
betweennumberof coarse-grainedof size class i, Ni, CGM(t) and collision r��elationship between
numberof coarse-grained
particlesparticlesof sizeFig. 6
6 RRelationship��� class i , N i , CGM
� �� �
CC
� 28）
for a coarse-grained
particle
size class i, ri,�particle
(t)
[28].
���
���
（t ）and collision
rate for ofa coarse-grained
of
size
class
i
,
r
,
）
�
i
CGM
CGM （t�
��� .
���
���
(388) Vol.29 No.4（2020）1000�� ���

24

�

����

� �� � ��� ��� ��� �� ���� ��� ���

��� �

����

�

𝐶𝐶��� �𝑡𝑡� � �� ����� �𝑡𝑡� ∙3.3.2

動に着目するため、流動層内全体としての空間率は時

��

∙

�

�����

��

��

2 粒子衝突頻度関数とプロセスパラメータの関係間によらず一定であるとみなして、流体抗力は粒子の粒 子 の 挙 動 が 安 定し た 1 〜 2 s の 1 s 間 の DEM-

断面積と流体の運動エネルギーの積で表した。

CFD 連成解析結果を用いて、流動層のプロセスパラ

さらに、衝突頻度は衝突相手の数、つまり粗視化粒��� ��� �� ∑ �� ������ ������ ���� ��� �������� ��� � �
�子の総数にも依存すると考えられる。
�� ��� ∑N�CGM�� ���
����� そこで、� ��（t ）(45)
メータと粒子間衝突頻度関数 C i , （の関係を解析しj t）

粗視化粒子の衝突頻度関数 C CGM（t ）と，粗視化粒子のた。はじめに、前節で求めた 1 粒子の衝突頻度関数を

総数 N CGM（t ）と Eq. の右辺の積との関係を調べるた

用いて、Eq. を Eq. で除して次式を得た。

�

�

∙

めに、それぞれを両対数プロットした図を Fig. 7 に示す。図より、各時点のプロットが同一直線上に並ぶことがわかった。そのため、粗視化粒子の衝突頻度関�� �� �� ∑ � ���� �� ����� �� �
�

�

� �

�

�

�

��

�

�

�

�

∙

�
� �は、定数��c�b��を用いて次式で表せること
�����
数 C�CGM
（
c a、
��
� ��
� t∑）

がわかった。
107

6
10
�𝑡𝑡�

𝐶𝐶���

105

t = 31.2 min� �� ����� �𝑡𝑡� ∙

��

∙

�

������

�

R2 = 0.9836

�

����

��
���

� �� �� �� � ��� ��� ��� �� ���� �������� ���

��� �

����


�

(45)

Eq. の右辺の分子と分母を縦軸と横軸にとりプ
Linear approximationt = 54.0 min(46)/C（t ） は時
プロットした図を Fig. 8 に示す。C i , （j t）
間によらず、1 とみなせることがわかった。そのため、

なお、粗視化倍率
L の影響を確認するため、最も10
y = 122.0x
y = 113.4x4�� ���

���� ���

ロットした図の傾き C i , （/C（t ）を、時間に対して
j t）

t = 42.5 min ����� ������ ����� ���

R2 = 0.9906

総粒子数が少ない時間
67.0 min において、粗視化倍103

は C（t ）と等しいと判断した。
i, （
j t）
y 2= C
90.33xR = 0.9802
𝐶𝐶��� �𝑡𝑡� � 𝐶𝐶�𝑡𝑡�

7
率 L を10
15、20、30倍に変更して同様の評価を実施しri, CGMC [s-1]t = 67.0 min
t = 74.5 min
Fig. 7 に示す。粗視化倍率
L を変 t = 84.7 min つまり、群 i の粒子が群 j の粒子と衝突する頻度は、た結果もあわせて106���� ������ ���𝐶𝐶 �𝑡𝑡� � 𝐶𝐶�𝑡𝑡�

𝐶𝐶��� �𝑡𝑡�∙�� �� ����� �𝑡𝑡� ∙化させた場合においても、各プロットは同一直線状に5�
������

���
(46)
群i の粒子の
1 粒子の衝突頻度と群 j の粒子の個数分
�𝑡𝑡� � 𝐶𝐶�𝑡𝑡�
𝐶𝐶���

化倍率 L = 1 のときの粗視化粒子の衝突頻度関数に等���
れ次式で表すことができる。
��������

��

�����

10

�

率を掛け合わせたものに等しい。

存在することが分かった。したがって、PBM
に実装
4
10

��
y 2= 83.97x ��� ������ � y 2=� 87.42xy 2= 90.92x
最終的に、1 粒子および 2 粒子の衝突頻度はそれぞするオリジナル粒子の衝突頻度関数t�）は、粗視�𝑡𝑡� ∙= 0.9718= �0.9715
R���
R = 0.9687∙CRORG（(47)
� �� ��
𝐶𝐶��� �𝑡𝑡������ ������ �3��� ���
�������
𝐶𝐶��� �𝑡𝑡� � 𝐶𝐶�𝑡𝑡�
10∙�
𝑟𝑟 � �𝑡𝑡� � 𝐶𝐶�𝑡𝑡�𝑁𝑁 �𝑡𝑡� � � �107

し い。 す な わち、 オ リ ジ ナ ル 粒t =子101.0の 衝min突 頻 度 関 数 t = 105.3 min
t = 94.7 min
610
C ORG（t ）は、Eq.  の N CGM（t ）を N ORG（t ）に 書 き 換 え
5y 2= 96.11x
R = 0.9735 ��� ������ �

3 � � ���𝑡𝑡�𝐶𝐶���10���� �𝑡𝑡� ∙

�

∙

�

��

�
� ��������� ����� ���
���� ������ ���𝑡𝑡������������
�𝐶𝐶�𝑡𝑡�𝑁𝑁� �� �� � � ����� ∙�𝑡𝑡� �
�𝑡𝑡� �
�� �𝑡𝑡�������
𝑁𝑁� �𝑡𝑡� ���� ���� 𝐶𝐶�𝑡𝑡�𝑁𝑁
𝐶𝐶�𝑡𝑡�𝑁𝑁���𝑡𝑡�
𝑟𝑟𝑟𝑟𝑟𝑟���
���
����
���
�����∙������
�𝑡𝑡�
�
∙
�
�
����
�
�
�
�
�
�
�
� ���

10
ることで得られる。
104��� ������ ����� ���

� ���

�����

�



��

������

�

�

�



�

Ni, CGM [-]��𝑡𝑡� � 𝐶𝐶�𝑡𝑡�𝑁𝑁� �𝑡𝑡�
𝑟𝑟���

�� ���� �� �

���� ��� ��

�
����
�����

��

101 102 103 104 105 10�1 102 103 104 105 101 102 103 104 105��� ������ �

∙

�

��

𝑁𝑁� �𝑡𝑡���������

(47) ��� ���

Fig. 6
for a coarse-grained particle of size class i, ri, CGMC(t) [28].
1000����� ���
����� �� � ��� ��� ��� �� ���� ��� ���

��� �

����

Eq. から、（1−c b）が Ouchiyama と Tanaka23） に

(48)

�

2

-1

CCGM (t) [s ]�� ���

�����Ci, j, CGM (t) / CCGM (t) [-]���� ���

��
���

100

���� ���

����� ����� ��� 10 ����

10

5

��
���

� �� �

y 2= 7.1944 x0.1921
R = 0.9837
�� �� � ��� �
�

�

� ���



����class i, Ni, CGM
��(t) and�����
collision rate ����� ���
Relationship between number of coarse-grained particles of size� �� �

� ���

���
�
��
�
�
�
�
��
y 2= 109.7x
y 2= 109.5x
��� ������ � ∙𝑟𝑟���� � � 𝑁𝑁 �𝑡𝑡� �� �����𝐶𝐶�𝑡𝑡�𝑁𝑁�� �� ���������𝑡𝑡� ��� ��� � ��� ����
� �𝑡𝑡� �� ���= 0.9736
R �= 0.9759∙ ������ � 𝑁𝑁�� �𝑡𝑡� ����𝑟𝑟���
� 𝐶𝐶�𝑡𝑡�𝑁𝑁� �𝑡𝑡�
∙ ����� � �� �� ��
��𝑡𝑡�
� �
𝑟𝑟 �𝑡𝑡� � 𝐶𝐶�𝑡𝑡�𝑁𝑁R �𝑡𝑡�
���

����

�� ���� ��� �� ���� � 6 �� � ���� �

10
NCGM (t) (p(t)–g) / gFr(t)2 [-]

10

7

1

(48)
0

40
60
80
100 120
Time, t [min]Fig. 7 Relation between collision frequency function of a coarse-grained particle, CCGM(t) and NCGM (ρpFig. 7 R elation between collision frequency− ρg) / ρgFr2function(for L = 25
at several time points and L = 15, 20, 30 at tFig.
=Fig.
67.0
8 8min).
Changesfrequency function betweenCtheratio
ratio ofof collisioncollision frequency
 hangesinintheof a coarse-grained particle,; L = 15,

; L = 20,

; L = 25,

; L = 30

020

and that ofa single between
coarse-grained
particle, CCGM(t) [28].
functioncoarse-grained
parti­C C G M（t ）a n d N C G M（ ρp −ρg ）/ ρg F r 2（ f o r
L =25 at several time points and L =15,
20, 30 at t =67.0 min）.

cles, Ci, j, CGM（t ）and that of a single coarsegrained particle, C CGM（t ）28）.

25

1.5
1
0.5
0

Experimental q3(d, t)Calculated q3(d, t)

t = t1 + 10 min

製剤機械技術学会誌 (389)Time, t [min]Fig. 8 Changes in the ratio of collision frequency function between coarse
and that of a single coarse-grained particle, CCGM(t) [28].

よって定義された制限度（degree of restriction）と一致することがわかった。本研究においては、1−c b=
0.8079 であり、1 に近いため、近接粒子との接触が多いことが推察される。
3.4

衝突頻度関数を実装した PBM による計算値-1

Frequency, q3(d, t) [m ]と実測値の比較
実際の医薬品製造スケールでの粒度分布の実測値と、DEM-CFD 連成解析から導出した衝突頻度関数を実装した PBM による計算値を Fig. 9 に示す。実測値と計算値は、本研究において着目している加水後半以降の時間帯において、極めて高い精度で一致した。すなわち、本 PBM を構築するために用いた仮定
によって、実際の顆粒凝集および解砕現象のメカニズ
ムを十分に説明できることを示唆している。また、本研究の試料は、粒子サイズの揃った単分散の試料ではなく、実際の医薬品組成物であるため、本モデルは高い実用性を有しているといえる。最も計算負荷の大き

1.5
1
0.5
0
1.5
1
0.5
0
1.5
1
0.5
0
1.5
1
0.5
0
1.5
1
0.5
0
1.5
1
0.5
0
1.5
1
0.5
0
10

な DEM-CFD 連成解析を実行する回数を削減しつつExperimental q3(d, t)Calculated q3(d, t)

t = t1 + 10 min

t = t1 + 20 min

t = t1 + 30 min

t = t1 + 40 min

t = t1 + 50 min

t = t1 + 60 min

t = t1 + 70 min
100

1000

Particle diameter, d [m]

10000

C
in experimentalandFig.
 hanges
Fig.
9 9 Changes
in experimental
and calculated
frequency distributions q3(dも、実際の生産プロセスに近似した製造規模およびプ

calculated frequency distributions
q 3（d , t ）in standard condition28）.

ロセス時間に対しても実用的に適用可能な PBM を構築できたことは、非常に意義がある。
標準条件 ⅰの実験値と計算値のフィッティングに-1

Frequency, q3(d, t) [m ]より求めた各定数を用いて、緩徐条件 ⅱおよび急速

1.5
1
0.5
0
1.5
1
0.5
0
1.5
1
0.5
0
1.5
1
0.5
0
1.5
1
0.5
0
1.5
1
0.5
0
1.5
1
0.5
0
10

Experimental q3(d, t)Calculated q3(d, t)Experimental q3(d, t)Calculated q3(d, t)

t = t1 + 10 min

t = t1 + 10 min
0t = t1 + 20 min

t = t1 + 20 min
0t = t1 + 30 min

t = t1 + 30 min
0t = t1 + 40 min

t = t1 + 40 min
0t = t1 + 50 min

t = t1 + 50 min
0t = t1 + 60 min

t = t1 + 60 min
0t = t1 + 70 min

t = t1 + 70 min
100

1000

Particle diameter, d [m]0
10000 10

100

1000

Particle diameter, d [m]10000

in experimental
and calculated
frequency
distributions
q 3（d(left),
Fig.
Changes
Fig.
10 10Changes
in experimental
and calculated
frequency
distributions
q3(d, t) in slow
and rapidt ）in slow（left）and rapid（right）conditions.
(right) conditions.

(390) Vol.29 No.4（2020）26

条件 ⅲにおいても PBM の計算を行い、粒度分布の実

5 謝辞

験値との比較を行った結果を Fig. 10 に示す。緩徐条
件および急速条件のいずれにおいても、計算値と実験

この度は、仲井賞若手研究者奨励賞を賜り、誠にあ

値が良好に一致した。すなわち、標準条件の結果から

りがとうございます。身に余る光栄な賞を頂き感激す

構築したモデルを用いて、運転条件を変更した場合のるとともに、これまでご指導いただいた皆様に、この粒度分布の推移を予測的に再現できることが検証され場をお借りして御礼申し上げます。この受賞論文は、た。これにより、本 PBM が未知の運転条件に対して

著者の博士論文の中核をなす研究をまとめたものでも高い汎用性を示すことが判明した。したがって、本す。本研究を遂行するにあたり、終始多くのご支援とモデルで記述される顆粒成長および解砕の機構は、少ご指導、ご鞭撻を賜りました大阪府立大学大学院工学

なくとも実用的なレベルで正しいことが検証された。

研究科 綿野 哲教授に心より感謝申し上げます。また、企業人という立場である私に本研究を行う機会を与

3.5

小括

えてくださり、ご理解とご支援、ご協力をいただきま

流動層造粒工程において、プロセスパラメータを考した塩野義製薬株式会社およびシオノギファーマ株

慮した新規 PBM を構築するために、流動層造粒機内式会社の関係者の皆様に深くお礼申し上げます。今後の粒子運動を DEM-CFD 連成解析により数値シミュも製剤分野の技術発展に貢献できるよう精進いたし

レーションした。その結果、以下の結論を得た。⑴フます。

ルード数を用いて、1 粒子の衝突頻度関数 C（t ）を導出することができた。⑵ 2 粒子間の衝突頻度関数 C i , j

参考文献（t ）は 1 粒子の衝突頻度関数 C（t ）と等しいことが判1 ）ICH Expert Working Group, ICH Harmonized明した。

Tripartite Guideline, Pharmaceutical Develo­

この衝突頻度関数と含水率を考慮した合一関数おpment Q8, Revision 2, 2009.

よび選択関数により、新規 PBM を提案した。本モデ2 ）U.S. Department of Health and Human Services,

ルを用いて算出した粒度分布を、生産スケールでの実

Food and Drug Administration, Guidance for験結果と比較検討し、モデルの精度を検証した。さらIndustry, PAT-A framework for innovativeに、異なる運転条件において予測的に利用が可能であpharmaceutical development, manufacturing,

るかを検証した。その結果、以下の結論を得た。⑶本and quality assurance, 2004.

PBM により算出した粒度分布は、実際のプロセスと3 ）大池敦夫, 寺下敬次郎,

宮南啓,

医薬品製造に同程度に長い時間範囲において、生産スケールの実測

おける撹拌造粒 , J. Soc. Powder Technol., Jpn.,値と極めて高い精度で一致した。⑷異なる運転条件に

24, 535-541 (1987).

ついても、本 PBM により粒度分布の時間変化を予測

4 ）H. Leuenberger, Granulation, New techniques,

することができた。

Pharm. Acta Helv., 57, 72-82 (1982).
5 ）P. Holm, T. Schaefer, H. G. Kristensen, Gran­ulation in high-speed mixers Part V. Power cons­

4 おわりにumption and temperature changes during
医薬品の製造工程で汎用される流動層造粒法に着granulation, Powder Technol., 43, 213-223 (1985).

目し、工程理解の推進および粒子成長機構の解明を目

6 ）P. Holm, T. Schaefer, H. G. Kristensen,

的として、実生産スケールでの実験結果を基に数値解

Granulation in high-speed mixers Part VI.

析およびモデリングに関する基礎的研究を行った。本

Effects of process conditions on power

研究を通じて、流動層造粒・乾燥工程における顆粒凝

consumption and granule growth, Powder

集・解砕現象を詳細に解明することができただけでなTechnol., 43, 225-233 (1985).

く、PAT ツールにより取得した工程データを解析す7 ）S. Schildcrout, Rheology of pharmaceutical

る手法を構築することができた。今後、構築したモデ

granulations, J. Pharm. Pharmacol., 36, 502-505

ルや解析手法を用いて、さらなる製剤設計および製造(1984).
8 ）S. Watano, T. Morikawa, K. Miyanami, Kinetics

プロセス開発の合理化および最適化を図りたい。

of granule growth in fluidized bed granulation
with moisture control, Chem. Pharm. Bull., 43,27

製剤機械技術学会誌 (391)Powder Technol., 77, 79-87 (1993).

1764-1771 (1995).
9 ）S. Watano, Direct control of wet granulation21）M. Sakai, Y. Yamada, Y. Shigeto, Numerical sim­

processes by image processing system, Powder

ulation of cohesive particles in a fluidized bed byTechnol., 117, 163-172 (2001).

the DEM coarse grain model, J. Soc. Powder
Technol., Jpn., 47, 522-530 (2010).

10）A. Burggraeve, T. Van Den Kerhof, M. Hellings,
J. P. Remon, C. Vervaet, T. De Beer, Evaluation of22）N. Ouchiyama, T. Tanaka, Estimation of thein-line spatial filter velocimetry as PAT moni­

average number of contacts between randomly

toring tool for particle growth during fluid bed

mixed solid particles, Ind. Eng. Chem. Fundam.,

granulation, Eur. J. Pharm. Biopharm., 76, 138-146

19, 338-340 (1980).
23）N. Ouchiyama, T. Tanaka, Mathematical model(2010).

in the kinetics of granulation, Ind. Eng. Chem.

11）C. Biggs, R. Boerefijin, M. Buscan, A. Salman, M.

Process Des. Develop., 13, 383-389 (1974).

Holunslow, Fluidised bed granulation : modelling
the growth and breakage kinetics using

24）N. Ouchiyama, T. Tanaka, The probability of

population balances., Proc. World Cong. Particle

coalescence in granulation kinetics, Ind. Eng.

Technol., Sydney, pp. 629-636 (2002).

Chem. Process Des. Develop., 14, 286-289 (1975).

12）H. Tan, A. Salman, M. Holunslow, Kinetics of25）N. Ouchiyama, T. Tanaka, Kinetic analysis of

fluidised bed melt granulation V: Simultaneouscontinuous pan granulation. Possible explanations

modelling of aggregation and breakage, Chem.

for conflicting experiments and several

Eng. Sci., 60, 3847-3866 (2005).

indications for practice, Ind. Eng. Chem. Process
Des. Develop., 20, 340-348 (1981).

13）K. Hayashi, S. Watano, Investigation of fluidized
bed granulation mechanism by using an online

26）S. Watano, T. Morikawa, K. Miyanami, Math­particle size monitoring device, J. Soc. Powder

ematical model in the kinetics of agitation

Technol. Jpn., 55, 316-323 (2018).

fluidized bed granulation. Effects of moisture14）T. Tanino, Control of critical quality andcontent, damping speed and operation time on

performance attributes of powder in manufactu­

granule growth rate, Chem. Pharm. Bull., 44, 409415 (1996).

ring processes of oral dosage forms and PAT27）A.A. Adetayo, J.D. Litster, S.E. Pratsinis, B. J.

application, J. Soc. Powder Technol., Jpn., 42, 638-

Ennis, Population balance modelling of drum

647 (2005).

granulation of materials with wide size15）T. Uemura, S. Uramatsu, M. Morizane, H.
Ichikawa, Y. Fukumori, Usefulness of mixer

distribution, Powder Technol., 82, 37-49 (1995).

torque rheometer for pharmaceutical design,28）K. Hayashi, H. Nakamura, S. Watano, Numerical
study on granule aggregation and breakage inPharm Tech Jpn., 28, 605-611 (2012).

fluidized bed granulation by a novel PBM with16）日本粉体工業技術協会, 造粒ハンドブック, オー

DEM-CFD coupling approach, Powder Technol.,

ム社, p.88 (1991).

360, 1321-1336 (2020).

17）P. A. Wauters, B. Scarlett, L. X. Liu, J. D. Litster,G. M. Meesters, A population balance model for

29）H. Hertz, Über die berührung fester elastischer

high shear granulation, Chem. Eng. Commun.,

Körper (On the contact of rigid elastic solids). In

190, 1309-1334 (2003).

: Miscellaneous Papers. Jones and Schott,18）P.A. Cundall, O.D.L. Strack, Discrete numerical

Editors, J. reine und angewandte Mathematik 92,

model for granular assemblies, Geotechnique, 29,

Macmillan, London (1896), p. 156 English

47-65 (1979).

translation : H. Hertz (1882).
30）R.D. Mindlin, Compliance of elastic bodies in19）T.B. Anderson, R. Jackson, Fluid mechanical des­contact, J. Appl. Mech., 16, 259-268 (1949).

cription of fluidized beds. Equations of motion,31）R.D. Mindlin, H. Deresiewicz, Elastic spheres inInd. Eng. Chem. Fundam., 6, 527-539 (1967).

contact under varying oblique forces, ASME,20）Y. Tsuji, T. Kawaguchi, T. Tanaka, Discrete par­September, 327-344 (1953).

ticle simulation of two-dimensional fluidized bed,(392) Vol.29 No.4（2020）2832）Y. Tsuji, T. Tanaka, T. Ishida, Lagrangian nume­continuous blending processes, 1 DEM based

rical simulation of plug flow of cohesionlessmethods., Macromol. Mater. Eng., 296, 290-307

particles in a horizontal pipe, Powder Technol.,(2011).
39）M. Sen, D. Barrasso, R. Singh, R. Ramachandran,71, 239-250 (1992).
33）E. Sakaguchi, E. Ozaki, T. Igarashi, Plugging ofA multi-scale hybrid CFD-DEM-PBM description

the flow of granular materials during theof fluid-bed granulation process, Processes, 2, 89-

discharge from a silo, Int. J. Mod. Phys. B, 7,

111 (2014).
40）H. Tan, M. Goldschmidt, R. Boerefijn, M. Hou­

1949-1963 (1993).
34）C.-Y. Wen, Y.H. Yu, Mechanics of fluidization,

nslow, A. Salman, J. Kuipers, Building population

Chem. Eng. Prog. Symp., 62, 100-111 (1966).

balances for fluidized bed granulation : lessons35）S. Ergun, Fluid flow through packed columns,

from kinetic theory of granular flow, Powder
Technol., 142, 103-109 (2004).

Chem. Eng. Prog., 48, 89-94 (1952).
36）D. Gidaspow, R. Bezburuah, J. Ding, Hydro­

41）P.C. Kapur, D.W. Fuerstenau, Coalescence model

dynamics of circulating fluidized beds, kinetic

for granulation, Ind. Eng. Chem. Proc. Design

theory approach., In Fluidization VII, Proc. 7th

Dev., 8, 56-62 (1969).
42）S. Cooper, C.J. Coronella, CFD simulations ofEng. Fund. Conf. on Fluidization, 75-82 (1992).

particle mixing in a binary fluidized bed, Powder37）M.J. Hounslow, R.L. Ryall, V.R. Marshall, A

Technol., 151, 27-36 (2005).

discretized population balance for nucleation,43）K. Hayashi, S. Watano, Novel population balance

growth, and aggregation, AIChE Journal, 34,

model for granule aggregation and breakage in

1821-1832 (1988).

fluidized bed granulation and drying, Powder38）A. Dubey, A. Sarkar, M.G. Ierapetritou, C.R.

Technol., 342, 664-675 (2019).

Wassgren, M.J. Muzzio, Computational app­roaches for studying the granular dynamics of

紙面づくりに参加しませんか！（原稿等の募集）「製剤機械技術学会誌」編集委員会では、紙面のさらなる充実を図ることを目的として会員各位から以下の原稿を募集しています。一般社会動向や本会関連動向、製剤機械技術に関連して日ごろ感じていること、本誌を見た感想、また紙面を飾るイラスト・写真等々、なんでも結構です。どしどしお寄せください。寄せられた原稿は下記の指針に従い編集委員会で採否を検討し、採用原稿は適宜本誌に掲載させていただきます。
【募集内容】

1 ．オピニオン（会員の声）：従来から募集しているものです。日ごろの皆様の意見をお寄せください。刷り上がり 1/2 から 1 頁以内。タイトルを付けて提出してください。
2 ．ショートコーナー（新設）：本会に関連する話題（本誌記事に関するコメントも含む）、思いなどを簡潔にまとめたもので、刷り上がり 1/4 から 1/2 頁以内。専門技術用語の解説、Ｑ＆Ａ集なども含む。
3 ．Book review：本会会員に有益と思われる新書、成書の紹介、読後感等。刷り上がり 1/4 から 1 頁程度。
4 ．会告：本会主催以外の関連学会やセミナーなどの案内。刷り上がり 1 頁以内。
5 ．イラスト・写真：本会会員オリジナルのイラスト、絵、写真など。印刷は白黒になります。
たとえば歴史的な製剤機械や建築物、人物の写真にタイトルまたは逸話を添えて。
必ずしも製剤機械技術に関連しない趣味の写真でも可。

【投稿規定】

1 ．投稿資格者：本会会員であること。
2 ．文字数またはスペース：刷り上がり 1/4 頁とは 400 字程度、1/2 頁とは 800 字程度に相当する。
3 ．採用の可否：原稿は随時受付、採用掲載の可否及び掲載箇所・時期は編集委員会が決定する。
投稿者の氏名は希望によりイニシャルも可。編集委員会により一部修正またはタイトルを付加することがある。採用された場合には、投稿者にその旨を連絡する。
4 ．原稿料及び著作権：採用された場合、薄謝を進呈する（会告を除く）。また著作権は製剤技術学会に属する。
5 ．投稿先：製剤機械技術学会誌 編集委員会： e-mail：info@seikiken.or.jp

29

製剤機械技術学会誌 (393)