Numerical Analysis of Viscous Polymer Resin Mixing Processes in High-Speed Blade-Free Planetary Blender Using Smoothed Particle Hydrodynamics

Abstract

High-speed planetary mixers can rapidly and efficiently combine rheological liquids, such as polymer resins and paste materials, because of the large centrifugal forces generated by the planetary motion of the mixing vessel. Only a few attempts have been made to computationally model and analyze the intricate mixing patterns of highly viscous substances. This paper presents meshless flow simulations of the planetary mixing of polymeric fluids. This research utilized the smoothed particle hydrodynamics (SPH) approach for numerical calculations. This method has advantages over the finite-volume method, which is a grid-based computational technique, when it comes to modeling interfacial and free surface flow problems. Newtonian rheology and interfacial surface force models were used to calculate the dissipative forces in the partial differential momentum equation of fluid motion. Simulations of the flow of an uncured polyurethane resin were carried out while it was mixed in a planetary mixer, under various operating conditions. Simulations using SPH were able to accurately reproduce the intricate flow and blending pattern, providing insight into mixing mechanics and mixing index evolution characteristics according to operating conditions for the planetary mixing of polymeric fluids. The simulation results showed that the spiral band, which promotes the mixing performance, is densely and distinctively formed under high-speed operation conditions.

1. Introduction

The mixing of rheological substances with very high viscosity can be achieved with blade-free planetary mixers [1,2]. The typical configuration and operating conditions of a planetary mixer are shown schematically in Figure 1. A mixing vessel inclined at an angle to the horizontal plane simultaneously revolves and rotates in the opposite direction to induce a swirling motion in the substances. A high-speed planetary motion can achieve a rapid and effective mixing of viscous fluids such as viscous vesicular phospholipid gels [3], styrene-butadiene rubber dispersed with carbon nanotubes [4], and nanoscale fat emulsions [5]. Yamaga et al. [6] experimentally showed that the blade-free planetary mixing of powder-dispersed rheological material for denture relining had superior blending quality compared to conventional stir mixing in terms of mixing time, mixing homogeneity, antifoaming capacity, and the mechanical strength of the cured material.

As in the aforementioned studies, previous experimental research activities mainly focused on investigating the effectiveness of blade-free planetary mixers for the satisfactory blending of various substances that cannot be blended entirely by other types of mixer. However, they provided no experimental details for the mixing process through which insight into the planetary mixing mechanism could be obtained. This was mainly because the fluid flow measurement in a blade-free planetary mixer is very challenging. Conventional in vitro flow visualization techniques [7] such as molecular tagging velocimetry (MTV), particle image velocimetry (PIV), and laser Doppler velocimetry (LDV) are difficult to apply for planetary mixing cases due to the orbital motion of a fluid container. Current velocimetry systems are not compact enough to be attached to a moving container. To the author’s knowledge, empirical research on flow visualization in planetary mixers was published very recently in [8]. The authors constructed an experimental apparatus to characterize the planetary mixing of silicon oil in the low revolution speed range (less than 60 rpm) by applying three flow visualization techniques: particle concentration imaging, planar PIV, and tracer particle tracking. They concluded that an optimum rotation-to-revolution speed range exists for planetary mixing associated with flow separation due to chaotic advection.

Despite a lack of experimental data, a few numerical studies have been undertaken to investigate planetary mixing processes and the performance of powders and fluids. Their computational results have revealed that the performance of planetary mixers strongly depends on the operating conditions, i.e., the revolution and rotation speeds [9,10,11]. The computational work in [9] reported that the homogeneous blending of cohesive powder materials could only be achieved with planetary mixers at or above critical rotation and revolution speeds. It also showed that a faster speed and higher rotation-to-revolution ratio tend to demonstrate more effective blending of cohesive powders. The discrete element method used in [9], unlike the computational fluid dynamics (CFD) approach adopted in [10,11], can simulate the formation and dynamic evolution of the free surface of a partially charged particle bed. However, it cannot model the flow of a highly viscous non-particulate rheological medium, which is the principal interest of this paper. Chergui et al. [10] carried out CFD simulations to study the effect of the rotation-to-revolution speed ratio on the vortical flow structure of a Newtonian fluid under planetary motion. However, they only investigated the flow characteristics in the planetary blender without considering its mixing capabilities. Yamagata and Fujisawa [11] numerically investigated the mixing process and performance of a blade-free planetary mixer using computational fluid dynamics (CFD). They concluded that maximum mixing performance can be achieved around the rotation-to-revolution speed rate of 0.5. However, the CFD models in [10,11] ignored the presence and complex dynamic behavior of the surface of the free fluid during the mixing process by setting the entire inner space of the mixing vessel as the fluid domain. Recently, Shen et al. [12] reported both the CFD and experimental validation of the flow characteristics of a partially filled fluid in a planetary centrifugal bioreactor. Although their study was the first numerical effort to investigate free surface flow in a planetary motion with a maximum revolution rate of 160 rpm, it did not provide direct insight into the planetary mixing performance of rheological fluid conditions at high revolution speeds of around 1000 rpm, which are typical operating conditions for industrial mixing applications.

Unlike previous studies, the uniqueness of this paper lies in the numerical investigation of the flow and mixing characteristics of a highly viscous polymer resin (approximately 20,000 times greater viscosity than water) in a high-speed planetary mixer with a revolution rate of up to 1500 rpm. The numerical method called smoothed particle hydrodynamics (SPH) was used for mixing simulations due to its computational capabilities to account for highly viscous fluid flow, Lagrangian fluid particle tracking, and the free fluid surface [13]. This paper specifically investigated the influence of operating conditions such as the revolution speed and the rotation-to-revolution speed ratio on the mixing performance of two-part adhesive resins in a blade-free planetary mixer. The transport-velocity model [14] and the free surface boundary model [15] were used to improve the stability of the calculation and the accuracy of the solution in the SPH mixing simulations.

The remaining parts of this paper are organized as follows: The next section briefly introduces the SPH theories and describes the simulation setup and model parameters. The following section presents the SPH simulation results and discussions. The last section concludes this numerical study of the high-speed planetary mixing of viscous polymer resins.

2. SPH Model and Simulation Setup

This section presents a brief introduction to SPH theories as well as detailed numerical setups for viscous resin mixing simulations.

2.1. SPH Equations and Models

Smoothed particle hydrodynamics is a nonlocal Lagrangian mesh-free method that uses a smoothing kernel function for an approximation of field functions and their mathematical operations such as differentiation and integration [16,17,18,19,20]. A continuous computational domain in SPH is discretized with a finite number of spatially dispersed Lagrangian particles. Each particle carries its associated smoothing kernel, as well as physical quantities. With this particle-based kernel representation scheme, any field function $f\left(\mathbf{r}\right)$ of a position vector $\mathbf{r}$ can be approximated as $\overline{f}$ with a kernel function W by

$$\overline{f}(\mathbf{r})=\sum _{j}f\left({\mathbf{r}}_j\right)W\left(\mathbf{r}-{\mathbf{r}}_j,h\right)V_j,$$

(1)

where j is the index of a neighbor particle within the support domain associated with the smoothing length h, and ${\mathbf{r}}_j$ and $V_j$ are the position vector and volume of the particle j, respectively. Therefore, the kernel function W plays a role as a weighting function of $f\left({\mathbf{r}}_j\right)$ in the summing process for calculating $\overline{f}\left(\mathbf{r}\right)$ in Equation (1).

The continuity equation and the momentum equation for a weakly compressible SPH particle i can be expressed as

$$\frac{d{\rho }_i}{dt}={\rho }_i\sum _j\frac{m_j}{{\rho }_j}{\mathbf{v}}_{ij}·{\nabla }_i{W}_{ij}$$

(2)

and

$$\begin{array}{cc}\hfill \frac{d{\mathbf{v}}_i}{dt}=& -\frac{1}{m_i}\sum _j\left({V_i}^2+{V_j}^2\right)\left\{\left(\frac{{\rho }_j{p}_i+{\rho }_i{p}_j}{{\rho }_i+{\rho }_j}\right){\nabla }_i{W}_{ij}\right.\hfill \\ \hfill & +\left(\frac{2{\mu }_i{\mu }_j}{{\mu }_i+{\mu }_j}\right){\mathbf{v}}_{ij}\left(\frac{{\mathbf{r}}_{ij}·{\nabla }_i{W}_{ij}}{{\left|{\mathbf{r}}_{ij}\right|}^2+{\eta }^2}\right)\hfill \\ \hfill & \left.-\frac{{\alpha }_i{\kappa }_i}{V_i}\left(\frac{{\rho }_j}{{\rho }_i+{\rho }_j}{c}_i^i+\frac{{\rho }_i}{{\rho }_i+{\rho }_j}{c}_j^i\right){\nabla }_i{W}_{ij}\right\}+\mathbf{g},\hfill \end{array}$$

(3)

where $\rho$ is the density, m is the mass of the particles, ${\mathbf{v}}_{ij}={\mathbf{v}}_i-{\mathbf{v}}_j$ is the relative velocity vector, ∇ is the gradient operator, $\mu$ is the dynamic viscosity, ${\mathbf{r}}_{ij}={\mathbf{r}}_i-{\mathbf{r}}_j$ is the relative position vector, $\eta$ is a small number to avoid singularity during computation, $\alpha$ is the coefficient of liquid surface tension, $\kappa$ is the curvature of the liquid surface, $c_i^i$ and $c_j^i$ are color functions, and $\mathbf{g}$ is the acceleration of the gravity vector [16,17,18]. The color function to identify the particle at the phase interface [17] is defined by

$$c_q^p=\left\{\begin{array}{cc}1\hfill & \mathrm{if}\mathrm{particle}p\mathrm{is}\mathrm{in}\mathrm{a}\mathrm{different}\mathrm{phase}\mathrm{to}\mathrm{particle}q\hfill \\ 0\hfill & \mathrm{if}\mathrm{particle}p\mathrm{is}\mathrm{in}\mathrm{the}\mathrm{same}\mathrm{phase}\mathrm{as}\mathrm{particle}q\hfill \end{array}\right.$$

(4)

The smoothing kernel used in this study is a piecewise smooth quintic spline with a cut-off radius of $3h$ [21], as shown below.

$$W(r,h)=W(s)=\frac{7}{478\pi }\times \left\{\begin{array}{cc}{(3-s)}^5-6{(2-s)}^5+15{(1-s)}^5\hfill & 0\le s<1\hfill \\ {(3-s)}^5-6{(2-s)}^5\hfill & 1\le s<2\hfill \\ {(3-s)}^5\hfill & 2\le s<3\hfill \\ 0\hfill & s\ge 3\hfill \end{array}\right.$$

(5)

where $s=r/h$.

As proposed in [14], to mitigate a tensile instability problem that may occur while performing free surface or multiphase flow simulations with SPH [22], a transport or advection velocity ${\tilde{\mathbf{v}}}_i$ is used instead of ${\mathbf{v}}_i$ to update the particle position ${\mathbf{r}}_i$ with

$$\frac{d{\mathbf{r}}_i}{dt}={\tilde{\mathbf{v}}}_i,$$

(6)

where ${\tilde{\mathbf{v}}}_i$ can be calculated at every time step $\Delta t$ in the following discretized form [14]:

$${\tilde{\mathbf{v}}}_i(t+\Delta t)={\mathbf{v}}_i(t+\Delta t)-\frac{p_b}{m_i}\sum _j\left(V_i^2+V_j^2\right){\nabla }_i{W}_{ij}\Delta t$$

(7)

On the right-hand side of Equation (7), the second term with the background pressure $p_b$ can effectively prevent particle voids due to numerical fracture in fluid stretching situations.

An artificial equation of state was adopted as a constitutive relation between the density $\rho$ and the pressure p in the following form:

$$p=p_0\left[{\left(\frac{\rho }{{\rho }_0}\right)}^{\gamma }-1\right]+p_b$$

(8)

where $p_0=\frac{{\rho }_0{c}_0^2}{\gamma }$ is the reference pressure, ${\rho }_o$ is the reference density, $c_0$ is the speed of sound in the reference state, and $\gamma$ is the artificial adiabatic index [21].

The computation time step is determined as a minimum or lower value that satisfies four conditions, i.e., the Courant–Friedrichs–Lewy or CFL condition, the viscous condition, the body force condition, and the surface tension condition with the time criteria corresponding to each condition sequentially listed in the MIN function as below:

$$\Delta t\le \mathrm{MIN}\left\{\frac{1}{4}\frac{h}{c_{max}+\left|v_{max}\right|},\frac{1}{8}\frac{g{h}^2}{\mu },\frac{1}{4}\sqrt{\frac{h}{g}},\frac{1}{4}\sqrt{\frac{\rho h^3}{2\pi \alpha }}\right\},$$

(9)

where $c_{\max }$ is the maximum artificial speed of sound, and $v_{\max }$ is the maximum speed of fluid flow [14]. In this paper, a second-order symplectic scheme was used that considers the transport velocity in Equation (7) to integrate the governing equations in Equations (2) and (3), as in [16].

In SPH numerical simulations, particle size, material modeling, and computational time significantly impact computational accuracy in general. Particle size affects spatial resolution and solution accuracy. Smaller particles enhance numerical stability and accuracy to some extent but increase computational load [23,24]. The choice of material models influences the representation of fluid flow behaviors [25]. Advanced constitutive models enhance the precision of simulations yet increase the computational requirements. The computational time affects the scale and precision of the simulation. Longer times with more solving iterations yield more accurate results but at the cost of increased computation. Therefore, a balanced consideration of these factors is vital to achieve precise SPH results efficiently.

2.2. Simulation Setup

A planetary mixer, as its name implies, is a device that simultaneously applies orbital and rotational motion to the mixing vessel [10]. Figure 1a is a schematic diagram showing the planetary movement of the cylindrical vessel with an inner radius of r. The axis of rotation moves in a circular motion with a radius of R in the revolution plane. The revolution speed and the rotation speed are denoted by $\mathrm{\Omega }$ and $\omega$, respectively. The axis of rotation is generally inclined by an angle theta with respect to the axis of revolution. The specific mixer configuration of interest in this simulation study has antidirectional revolution and rotation axes with an inclination angle $\theta$ of 45${}^{\circ }$, which is known to exhibit optimal performance [9,26,27]. The orbital and rotational motion of the vessel starts at 0 rpm from the start time to the acceleration time $t_{\mathrm{a}}$, as schematically illustrated in Figure 1b, and then maintains the respective steady state speeds of ${\mathrm{\Omega }}_{\mathrm{ss}}$ and ${\omega }_{\mathrm{ss}}$. Here, the rotation-to-revolution speed ratio is defined as n = ${\omega }_{\mathrm{ss}}$/${\mathrm{\Omega }}_{\mathrm{ss}}$. In this paper, a numerical analysis study was conducted on the two-degree-of-freedom planetary mixing capable of an arbitrary combination of revolution and rotation speeds. Specifically, a total of 21 simulations were conducted with a combination of seven cases of revolution speeds (150 rpm, 300 rpm, 450 rpm, 600 rpm, 900 rpm, 1200 rpm, and 1500 rpm) and three ratios of rotation-to-revolution speed (1/3, 1/2, and 1) for each revolution speed.

The viscous polymer resin selected to perform blade-free planetary mixing simulations was a two-part epoxy adhesive (Scotch-Weld™ DP100 Plus Clear manufactured by 3M™, Two Harbors, MN, USA). This adhesive consists of a mercaptan-based accelerating hardener and an epoxy resin, denoted by part A and part B, respectively. The material properties of the hardener and epoxy resin are shown in

Table 1. Properties of 3M™ DP100 Plus Clear epoxy resin and mercaptan-based hardener.

Property	Part A (Hardener)	Part B (Epoxy Resin)
Density ($\mathrm{kg}/{\mathrm{m}}^3$)	1150	1170
Dynamic viscosity ($\mathrm{Pa}·\mathrm{s}$)	19.4	11.1
Admissible fluid compressibility	0.01	0.01

. Regarding fluid constitutive models, the linear Newtonian viscosity models and non-Newtonian models [9,25] such as the power law model, the Bingham plastic model, and the Herschel–Bulkley model are widely used for highly viscous and rheological fluids. The maximum shear flow rate under the operating speed conditions of this paper was 39.45 s⁻¹. A published experimental and modeling study on the rheology of uncured epoxy resin in [28] observed linearity between the shear stress and shear strain rate for a shear rate below 40 s⁻¹, and the yield stress was very low, at about 0.2 Pa. Taking into account the small yield stress value and the stress–strain rate linearity, this numerical study appropriated the modeling of the epoxy resin fluid as a Newtonian fluid in the operating velocity regime. The mixing ratio of the two parts was 1:1 in volume, and the working time of the mixture was 4 min. When preparing a small amount of epoxy adhesive, the two parts can be easily mixed using a static mixer such as a nozzle-type inline mixer [29,30]. However, when preparing a large amount, such as for the manufacturing of epoxy matrix composite materials, a motorized resin mixer, such as a blade-free high-speed blender, should be used for complex homogeneous mixing in just a few seconds, considering both the resin working time and the material handling time [3,4].

Figure 2 shows the results of the preprocessing of the mixing vessel and the two epoxy components for the simulation of planetary mixing. The opacity of the wall particles was adjusted to 0.15 to gain a view of the fluid particles inside the container. The interior radius of the cylindrical container was 4.0 cm, the height was 8.5 cm, and the capacity was 427.3 cm³. The container was filled to a capacity of 58% with a viscous fluid, with parts A and B each comprising 29%. A rectangular knob was affixed to the top of the vessel to make it easier to observe the vessel’s rotation during the simulation. The rotation arm is a structural component for implementing planetary motion. The arm was excluded from the SPH simulation because it had no direct interaction with fluid particles. Once the original meshed CAD model was converted to SPH particles with Altair SimLab™ software version 2022.2, all particle overlaps were manually corrected. The total number of modeled particles was 52,824. Parts A, B, and the vessel comprised 15,466, 15,179, and 22,179 particles, respectively. The SPH particles were placed at an average interval of 2 mm, which was dense enough to prevent fluid particles from numerically penetrating solid walls during rapid motion.

In this study, the flow and mixing of viscous polymer resins were induced by the planetary motion of the enclosed vessel without the inlet and outlet boundaries. When the steady-state revolution velocity and the rotation-to-revolution speed ratio were both set, the operating speeds during the simulation time were determined in the form shown in Figure 1b. If the revolution speed was set to accelerate at 200 rpm per second linearly, the acceleration could be calculated as $t_{\mathrm{a}}={\mathrm{\Omega }}_{\mathrm{ss}}/200$ s. Assuming that the mixing vessel was a rigid body, time-series data of the vessel’s center-of-mass position and angular velocity could be generated from the angular velocity curve. The motion data were calculated and saved as an ASCII file using an in-house Python code so that an SPH solver could look up the motion data file at each computation step and update the vessel’s spatial configuration by time interpolation and numerical integration procedures.

The SPH solver adopted in this study was Altair nanoFluidX™ (or nFX for short) software version 2022.2. The nanoFluidX™ code could start a particle simulation by reading the simulation configuration file, including the numerical, domain, phase, and motion parameters.

Table 2. SPH simulation parameters.

Parameter	Value
Simulation end time 1	Variable in range of 9.0–60.0 s
Integration time step 1	Variable in range of 2.067–20.67 $\mathsf{\mu }$s
Data saving interval	0.01 s
Reference velocity factor	1.5
Maximum number of iterations	1,000,000,000

¹ Specific parametric value for each simulation case is presented in Table 3.

shows the key simulation parameters defined in the configuration file. The mixing rate of each simulation case adjusted the simulation end time. The computation time step was determined to satisfy four criteria in Equation (9). The simulation result files for data visualization had a storage interval of 0.01 s.

Table 3. Simulation cases and running options.

Case No.	${\mathbf{\Omega }}_{\mathbf{ss}}$ (rpm)	${\mathit{\omega }}_{\mathbf{ss}}$ (rpm)	${\mathit{t}}_{\mathbf{a}}$ (s)	Simulation Time (s)	$\Delta \mathit{t}$ ($\mathsf{\mu }$s)	No. of Iterations	Computation Time 1 (h)
1	150	50	0.75	60.0	20.67	2,902,832	5.15
2	150	100	0.75	60.0	20.67	2,902,832	5.12
3	150	150	0.75	60.0	20.67	2,902,832	5.13
4	300	100	1.50	15.0	10.33	1,451,416	2.57
5	300	200	1.50	15.0	10.33	1,451,416	2.53
6	300	300	1.50	15.0	10.33	1,451,416	2.55
7	450	150	2.25	15.0	6.890	2,177,124	3.75
8	450	300	2.27	15.0	6.890	2,177,124	3.72
9	450	450	2.25	15.0	6.890	2,177,124	3.73
10	600	200	3.00	12.0	5.167	2,322,266	4.00
11	600	400	3.00	12.0	5.167	2,322,266	3.93
12	600	600	3.00	12.0	5.167	2,322,266	3.95
13	900	300	4.50	9.0	3.445	2,612,549	4.43
14	900	600	4.50	9.0	3.445	2,612,549	4.45
15	900	900	4.50	9.0	3.445	2,612,549	4.43
16	1200	400	6.00	9.0	2.584	3,483,398	5.92
17	1200	800	6.00	9.0	2.584	3,483,398	5.88
18	1200	1200	6.00	9.0	2.584	3,483,398	5.92
19	1500	500	7.50	9.0	2.067	4,354,248	7.35
20	1500	1000	7.50	9.0	2.067	4,354,248	7.38
21	1500	1500	7.50	9.0	2.067	4,354,248	7.37

¹ Simulations were performed on the GPU (model Quadro^® P5000 manufacured by Nvidia, Santa Clara, CA, USA) with 2560 CUDA cores and 16 GB memory.

shows the operating conditions (steady-state revolution speed, steady-state rotation speed, and acceleration time) and the computing information (simulation time, time step size, number of iterations, and total computation time) for the 21 simulation cases. The total parallel computing time of the GPU (graphic processing unit) was 99.26 h.

3. Simulation Results and Discussion

The particle data in the SPH solver nFX were stored as a pvtu extension, which is a parallel version of the XML-based vtu (Visualization Toolkit for Unstructured Grids) file format. These vtu and pvtu files can be post-processed by an open-source Visualization Toolikt (VTK) software package version 9.3.0.RC1. The simulation data files generated for every simulation saving interval contained information on each SPH particle’s location, speed, pressure, density, and phase identification number. The results of the SPH simulation were processed and visualized using the VTK-based open-source software ParaView version 5.11.1.

Figure 3 depicts successive post-processed images at intervals of 2 s over 60 s for Case 1 (${\mathrm{\Omega }}_{\mathrm{ss}}$ = 150 rpm, ${\omega }_{\mathrm{ss}}$ = 50 rpm) with the lowest rotation speed among the simulation cases listed in Table 3. In the same manner as in the image processing method of Figure 2b, epoxy part A and part B were represented by red and blue spherical particles, respectively, to visually confirm the particle flow and mixing patterns. During the first 2 s of mixing, the particles maintained contact with the initial boundary surface to some extent while exhibiting a rotational flow similar to that of a rigid body. Discernible flow patterns that penetrated each other appeared on the upper free surface after 2 s and lasted for about 10 s instead of occurring on the bottom surface. After 10 s, the two fluids started to swirl and create a twisted flow pattern, causing the particle groups that were initially connected to break apart. The connectivity of the clustered fluid particles deteriorated with material convection in a whirlwind manner after 20 s. A significant breakdown of the cluster of fluid particles occurred around 50 s into the observation, resulting in an increased level of disorder and blending.

The degree of mixing of materials can be qualitatively determined by visually observing the degree of disorder in the color distribution of the particles from post-processed simulation images, like in Figure 3. However, an accurate and objective mixing level can be quantified when a quantifiable mixing index is defined and used instead of a qualitative approach. This paper adopted a mixing index that was defined as follows based on the distance between the center of gravity of the two dispersed materials [9]:

$$M=1-\frac{\left|{\mathbf{r}}_A-{\mathbf{r}}_B\right|}{\left|{\mathbf{r}}_{A,0}-{\mathbf{r}}_{B,0}\right|}$$

(10)

where ${\mathbf{r}}_A$ and ${\mathbf{r}}_B$ are, respectively, the current mass center vectors of fluid part A and part B at time t, and ${\mathbf{r}}_{A,0}$ and ${\mathbf{r}}_{B,0}$ are their corresponding initial position vectors.

Before operating the mixer, the mixing index in Equation (10) initially had a value of 0. Once the two materials were perfectly mixed and their centers of gravity were ideally the same, the mixing index had a value of 1. If the distance between the centers of gravity of the two materials increased during the mixing process, the mixing index could result in a negative value due to centrifugal force. Figure 4 depicts the mixing index curve of simulation Case 1, Case 2, and Case 3 with a revolution speed of 150 rpm and a rotation-to-revolution speed ratio of 1/3, 1/2, and 1, respectively. Up to about 8 s, a higher rotation-to-revolution speed ratio led to an increased rate of change of the mixing index. For n = 1, which had the highest rotation speed among the three cases, the mixing index did not increase consistently. Instead, it fluctuated after 10 s, resulting in poorer mixing compared to the other two cases. A mixing index of 80% was achieved in 11 s in all three scenarios. However, the subsequent increase in the degree of mixing was typically much slower than the initial mixing stage. Based on this analysis, the centrifugal force generated by the 150 rpm speed setting appeared to be unable to efficiently and homogeneously stir the chosen epoxy resin.

Further simulations were conducted to verify the mixing performance with the increased centrifugal forces compared to the low-speed operating conditions of 150 rpm. Eighteen simulations (Case 4–Case 21) were carried out in a range of revolution speeds of 300–1500 rpm, taking into account the range of operation speed of conventional industrial planetary mixers. Among them, the data visualization results for six cases with n = 1/3 are shown in Figure 5 at 1.5 s intervals from 1.5 s to 9 s. The thick border lines in Figure 5 indicate that the mixing vessel was in the acceleration state before reaching the steady state. Material mixing did not make significant progress during the acceleration period when ${\mathrm{\Omega }}_{\mathrm{ss}}$ < 450 rpm. On the other hand, when ${\mathrm{\Omega }}_{\mathrm{ss}}$ > 600 rpm, a noticeable improvement in material mixing was observed during acceleration. Based on the results at t = 1.5 s and t = 3.0 s, it was evident that when the revolution speed was higher, the spiral band of the fluid particles in contact with the cylindrical wall became denser and more distinct. Figure 3 shows that particle clusters separated after 50 s at 150 rpm. On the other hand, Figure 5 shows that particle segregation occurred quickly within 10 s at 300 rpm. This led to a significant enhancement in the mixing performance, depending on the operating conditions.

Figure 6 shows the time evolution curve of the mixing index for 18 simulation cases (Case 4–Case 21) in the range of revolution speeds of 300–1500 rpm. This graph comparing the mixing indices reveals three features of the rapid planetary mixing of epoxy resin. First, as the revolution speed increases, the initial slope of the mixing index curve, which represents the initial mixing rate, also increases. Second, when ${\mathrm{\Omega }}_{\mathrm{ss}}\ge$ 900 rpm, the saturation value of the mixing index can reach nearly 100%. Third, the mixing index differences according to the rotation-to-revolution speed ratio are not significantly noticeable when ${\mathrm{\Omega }}_{\mathrm{ss}}\ge$ 600 rpm, unlike when ${\mathrm{\Omega }}_{\mathrm{ss}}\ge$ 450 rpm. Based on the qualitative analysis results presented in Figure 6, it has been determined that the optimal mixing performance is achieved through a combination of ${\mathrm{\Omega }}_{\mathrm{ss}}$ = 900 rpm and ${\omega }_{\mathrm{ss}}$ = 300 rpm. This combination of operating conditions results in a fast mixing time, a homogeneous resin mixture, and a low-energy consumption process.

Figure 7 compares the characteristics of the initial operating phase (t < 3 s) of the low-speed (${\mathrm{\Omega }}_{\mathrm{ss}}$ = 300 rpm) and high-speed (${\mathrm{\Omega }}_{\mathrm{ss}}$ = 1500 rpm) planetary mixing situations. It can be seen that the spiral band, which promoted planetary mixing performance, was densely and distinctively formed under the high-speed operation conditions. Additionally, as seen in Figure 8, the slope and area of the upper free surface were larger in the high-speed conditions than in the low-speed conditions. Viscous fluids tend to experience faster convection as a result of the lower flow resistance on the free surface compared to the cylindrical wall surface. Therefore, the mixing performance under high-speed conditions was enhanced because of the rapid and dense formation of spiral bands on the walls, while the large free surface area allowed low-resistant fluid diffusion.

The flow characteristics and mixing states of the two viscous fluids can be observed by examining the visualized particle data presented in Figure 3 and Figure 5. Furthermore, the degree of mixing can be measured and quantified from the mixing index curves shown in Figure 4 and Figure 6. An analysis of the simulation data in Figure 3, Figure 4, Figure 5 and Figure 6 can provide insight into the influence of the operating conditions on the flow patterns, mixing mechanisms, and characteristics of mixing index evolution. Based on the simulation results, the two merit figures could be defined and used to evaluate the mixing performance, namely, the rising time and the mixing time. In this paper, the 60% rising time $t_r$ measures the efficiency of mixing during the initial stage before achieving uniform stirring. It is determined by the time the mixing index first reached 60%. The 90% mixing time $t_m$ is a performance parameter that determines how quickly a homogeneous blend can be achieved. It is defined as the minimum time required for the mixing index to remain at or above 90%.

Figure 9 shows a 60% rising time bar graph extracted from the mixing index data of 21 simulation cases. At low-speed operation (${\mathrm{\Omega }}_{\mathrm{ss}}$ = 150 rpm), the rising time exceeded 6 s. However, it significantly improved to less than 3 s when operating at a high speed (${\mathrm{\Omega }}_{\mathrm{ss}}\ge$ 300 rpm). When driving at 1200 rpm or higher, it was possible to achieve a swift rising time of approximately 1 s. When operating at ${\mathrm{\Omega }}_{\mathrm{ss}}$ = 150 rpm, the rising time increased as the rotation-to-revolution speed ratio increased. However, under high-speed driving conditions of ${\mathrm{\Omega }}_{\mathrm{ss}}\ge$ 300 rpm, there was no overall consistent correlation between the rotation-to-revolution speed ratio and the rising time. It is imperative to select conditions where n = 1/3, because this would lead to significantly lower energy consumption when driving the rotary motor, as opposed to when n = 1/2 or n = 1.

The bar chart in Figure 10 displays the 90% mixing time, calculated using the mixing index data from 21 simulation cases. A minimum of 30 s of mixing was required at a low-speed setting of ${\mathrm{\Omega }}_{\mathrm{ss}}$ = 150 rpm to guarantee a uniform mixture. Under the operating conditions of ${\mathrm{\Omega }}_{\mathrm{ss}}\ge$ 300 rpm, the mixing time could be reduced by 1/10 or more compared to the case of ${\mathrm{\Omega }}_{\mathrm{ss}}$ = 150 rpm. In particular, when the revolution velocity exceeded 900 rpm, the 90% mixing process took only 4 s or less. This made it ideal for blending epoxy adhesive within a restricted time frame and ensured a consistent mixture. The correlation between the rotation-to-revolution speed ratio and the mixing time in all simulation cases was observed to be unclear. However, it should be noted that this differed from the planetary mixing of cohesive pharmaceutical powders, which demonstrates a trend towards a decreased mixing time as the rotation-to-revolution speed ratio increases [9]. Based on the analyses of Figure 9 and Figure 10, the optimal operating parameters of the two-part epoxy adhesive could be determined as ${\mathrm{\Omega }}_{\mathrm{ss}}$ = 900 rpm and n = 1/3, corresponding to ${\omega }_{\mathrm{ss}}$ = 300 rpm, taking into account the energy consumption of the process and the adequacy of the mixing time.

4. Conclusions

This paper presents the results of an SPH simulation study conducted on a two-part epoxy adhesive system. The two substances were highly viscous fluids that demanded uniform mixing and rapid processing. The scope of this numerical research extended beyond previous studies by encompassing a high-speed revolution range of up to 1500 rpm, which was necessary to account for the operational capabilities of conventional industrial mixers. It is crucial to emphasize that the mixing vessel may not be entirely filled with the working fluids in real-world applications. Therefore, it is important to consider the free surface dynamics of partially filled fluids for SPH modeling and simulation. Planetary mixing simulations were performed for 21 cases, taking into account the revolution-to-revolution ratios of 1/3, 1/2, and 1 for each revolution speed of 150 rpm, 300 rpm, 450 rpm, 600 rpm, 900 rpm, 1200 rpm, and 1500 rpm. Post-processing and data analysis resulted in visualized simulation graphics, mixing index–time curves, and bar charts on rising and mixing times. Four key conclusions were drawn by the qualitative and quantitative analysis of the post-processed images and graphs. First, it was observed that when highly viscous resins were stirred in a high-speed planetary mixer, the particle fluid flow tended to be laminar with swirling patterns rather than chaotic. Second, as the revolution speed of the mixing vessel increased, the spiral band responsible for effective mixing occurred more quickly in a denser manner. This indicated a significant improvement in both the rising and mixing times, which ultimately resulted in better mixing performance. Third, at high-speed revolution speeds above 300 rpm, there was no consistent correlation between the rotation-to-revolution speed ratio and mixing performance, except for when operating under the low-speed conditions of 150 rpm. Fourth, to ensure the optimal mixing performance of an uncured two-part epoxy adhesive in terms of rising time, mixing time, and process energy consumption, it is imperative to use a combination of a 900 rpm revolution speed and a 300 rpm rotation speed for a given mixer configuration and vessel geometry. This study is a valuable contribution to the relevant research field, as it not only amplifies our understanding of planetary mixing mechanisms, encompassing the correlation between operating conditions and mixing performance, but also presents crucial insights into the most favorable operating conditions for high-viscosity polymeric fluids.