Microscale Discrete Element Method Simulation of the Carbon Black Aggregate Fracture Behavior in a Simple Shear Flow

One of the biggest challenges in the transition from fossil fuels to renewable energies is the development of suitable energy storage systems. Over the past decade, Li-ion batteries have proven to be a reliable solution for many applications. Since their introduction in the 1990s, their sales have increased thanks to the expansion of the electronics market. Over the last decade, the use of Li-ion batteries in electric vehicles has led to a rapid increase in the demand for high-capacity, high-performance batteries. These developments pushed the Li-ion battery into focus of the international research community. Since then numerous studies were performed all with the goal of improving the battery performance, capacity, and life span. Previous work states that in addition to the chosen materials, the microstructure of the electrodes is of a crucial importance for the performance of the cells. Further investigation revealed the complex relation between the manufacturing process of the cells that determines the particulate structure of the electrode and the functionality of the final batteries. These findings demanded a deeper insight into the interactions of components, especially in the early stages of production that govern the structure formation. One important element of the composite, particularly in the positive, are the conductive additives such as carbon black (CB). They compensate for the overall poor conductivity of the active material and have a substantial effect on the structural properties. Furthermore, the proper distribution of CB within the electrode layer is critical for the battery performance. However, the microstructure of CB changes due to experienced shear stress during the manufacturing process and therefore also its impact on overall battery capacity and performance. Hence, a large part of recent research is focused on the investigation of the relationships between the manufacturing process and the final electrode structure. For these reasons, the presented work will mainly focus on the positive. It is worthwhile to mention that the latter is often labeled as cathode in the cited literature as it references the discharge process. The authors maintain this terminology throughout the artcle. As large experimental parameter studies are associated with great time and material investment, simulative studies of the dispersion process are a promising alternative. These do not require the use of expensive materials or sophisticated analytical methods, but offer a detailed insight into the processes during the manufacturing process. For this reason, the aim of this publication is the simulative investigation of CB aggregates fracturing during the dispersion process of the cathode slurry in a planetary mixer. However, the reliability of simulation results depends largely on the applied models. Finding an appropriate approach with a sufficient accuracy while maintaining reasonable computing effort poses a major challenge for simulative studies. Furthermore, the present task includes relevant phenomena on several different scales. Figure 1 provides a schematic illustration of its multiscale character. First, there are macroscopic flow patterns which cause shear gradients at the micro level. The latter induce fracturing of nanoscaled CB aggregates. Therefore, the task is divided into a macroscale computational fluid dynamics (CFD) simulation of the flow patterns inside the planetary mixer and a microscale discrete element method E. Asylbekov, R. Trunk, Dr. M. J. Krause, Prof. H. Nirschl Institute of Mechanical Process Engineering and Mechanics Karlsruhe Institut of Technology Strasse am Forum 8, Bldg. 30.70, 76131 Karlsruhe, Germany E-mail: ermek.asylbekov@kit.edu


Introduction
One of the biggest challenges in the transition from fossil fuels to renewable energies is the development of suitable energy storage systems. Over the past decade, Li-ion batteries have proven to be a reliable solution for many applications. Since their introduction in the 1990s, their sales have increased thanks to the expansion of the electronics market. Over the last decade, the use of Li-ion batteries in electric vehicles has led to a rapid increase in the demand for high-capacity, high-performance batteries. [1] These developments pushed the Li-ion battery into focus of the international research community. Since then numerous studies were performed all with the goal of improving the battery performance, capacity, and life span. Previous work states that in addition to the chosen materials, the microstructure of the electrodes is of a crucial importance for the performance of the cells. [2][3][4] Further investigation revealed the complex relation between the manufacturing process of the cells that determines the particulate structure of the electrode and the functionality of the final batteries. [5][6][7] These findings demanded a deeper insight into the interactions of components, especially in the early stages of production that govern the structure formation.
One important element of the composite, particularly in the positive, are the conductive additives such as carbon black (CB). They compensate for the overall poor conductivity of the active material and have a substantial effect on the structural properties. [8,9] Furthermore, the proper distribution of CB within the electrode layer is critical for the battery performance. [10] However, the microstructure of CB changes due to experienced shear stress during the manufacturing process and therefore also its impact on overall battery capacity and performance. [11] Hence, a large part of recent research is focused on the investigation of the relationships between the manufacturing process and the final electrode structure. [12,13] For these reasons, the presented work will mainly focus on the positive. It is worthwhile to mention that the latter is often labeled as cathode in the cited literature as it references the discharge process. The authors maintain this terminology throughout the artcle.
As large experimental parameter studies are associated with great time and material investment, simulative studies of the dispersion process are a promising alternative. These do not require the use of expensive materials or sophisticated analytical methods, but offer a detailed insight into the processes during the manufacturing process. For this reason, the aim of this publication is the simulative investigation of CB aggregates fracturing during the dispersion process of the cathode slurry in a planetary mixer. However, the reliability of simulation results depends largely on the applied models. Finding an appropriate approach with a sufficient accuracy while maintaining reasonable computing effort poses a major challenge for simulative studies. Furthermore, the present task includes relevant phenomena on several different scales. Figure 1 provides a schematic illustration of its multiscale character. First, there are macroscopic flow patterns which cause shear gradients at the micro level. The latter induce fracturing of nanoscaled CB aggregates. Therefore, the task is divided into a macroscale computational fluid dynamics (CFD) simulation of the flow patterns inside the planetary mixer and a microscale discrete element method DOI: 10.1002/ente.202000850 The shear stress induced breaking behavior of carbon black (CB) aggregates during the manufacturing process of Li-ion batteries is investigated via microscale discrete element method (DEM) simulations. The relevant range of shear stress is chosen according to a planetary mixer and cathode slurries with high solid content. Aggregates of different sizes and shapes are modeled using a self-written algorithm based on the tunable dimension method. Then, suitable models are chosen for representing the solid bridges between the primary particles of the CB aggregates and relevant fluid forces. The results show a correlation between aggregate size and critical shear stress which is required to initiate aggregate fracturing. Furthermore, a change in aggregate shape is linked to applied stress and initial aggregate size and shape. Hence, a recommendation for an efficient disintegration of CB aggregates during the mixing process is made.
(DEM) simulation of the shear stress induced CB aggregates fracturing.
Although CFD simulations provide a distribution of the resulting shear stresses, this work only focuses on the effects of shear stress on CB aggregates of different sizes and shapes. Similar studies were performed by Higashitani et al. [14] Here the authors investigated the disintegration of loose agglomerates of spherical particles in shear and elongational flows. The authors found a correlation between the average number of particles in a fragment of the agglomerate and the ratio between hydrodynamic drag force and adhesive force between particles. Although Higashitani et al. speak of aggregates in their work, these studies only regarded agglomerates. Nichols et al. discussed the often misleading use of the terms "aggregate" and "agglomerate" in the literature and proposed a distinction of both terms by defining "agglomerates" as particular clusters connected by attractive forces and "aggregates" as constructs of primary particles which are connected with solid bridges. [15] Considering this strict definition, the investigation of aggregates fracturing processes remains a gap within the current research. With the performed simulations, the authors try to provide meaningful insight into aggregate behavior under shear stress of a viscous fluid.

Simulation Model
The DEM was introduced by Cundall in 1971 and was since then continuously improved and applied for various particle-based problems. [16] However, the simulation of large particle clusters especially at the microscale including multiphase systems is still a current issue. Although there exist many different ways to approach such complex interactions, some of which will also be addressed in the course of this publication; the basics of DEM remain valid for all applications.
The DEM treats a particle as a perfectly spherical element with a mass m, a radius r, and a moment of inertia I. Its translational and rotational acceleration dũ=dt and dω=dt are described by Newton's equations as a sum of all forcesF j or torquesT j acting on the center of the particles mass (Equation (1a,b)).
The translational and angular velocities are obtained by integration over a short time step dt.
In this work, DEM is used to calculate the behavior of nanoscaled CB aggregates, represented as clusters of spherical particles. The dominant forces include contact forcesF c , hydrodynamic drag forcesF D , and bond forcesF b . Due to the small mass of the CB primary particles and the occurring lifting force of the surrounding fluid in combination with high fluid viscosity, the gravitational forces are neglected. The contact forces result from the physical interaction of two particles and are based on the elastic contact model introduced by Hertz and later expanded by Mindlin. [17,18] This approach describes the contact force as a combination of normal and tangential forcesF n c andF t c , respectively. The basic model consists of a parallel interaction of an elastic spring with an elastic constant k and a damping element with a viscoelastic damping constant γ, both in normal and tangential directions. The calculation also includes the virtual normal and tangential overlaps δ n and δ t and the normal and tangential relative velocitiesũ n rel andũ t Even though this contact model is broadly used in DEM simulations, many linear and nonlinear models for the calculation of the elastic and viscoelastic damping constants have been proposed by several research groups. A great overview of the most common models for the normal force calculation was given by Kruggel-Emden et al. [19] Additional effects such as lubrication forces, which act as a damping agent during particle interaction are usually considered in the DEM with a very small coefficient of restitution.
In addition to the basic DEM, three major issues need to be addressed for a proper simulation setup: 1) Modeling the complex structure of the CB aggregates; 2) Finding and implementing a suitable model for a stress-induced breakup of the solid sinter bridges; 3) Finding and implementing a suitable model for the fluid force field.
Furthermore, the authors will explain as how these issues were approached.

Modeling the CB Aggregates
The complex structure of nanoscaled aggregates is very difficult to quantify. The same applies to the CB aggregates. The latter are built from small primary particles (d p % 10100 nm) which are connected to each other by sinter bridges. An aggregate can consist of several hundred primary particles. Due to their small size, the aggregates experience strong interparticle adhesive forces and form larger agglomerates. Due to their irregular and often widely branched shape, adhesion is also possible due to steric conditions. However, Forrest and Witten, [20] Sorensen, [21] and Kätzel et al. [22] have described the shape of the aggregates with two values using light scattering. They described the size of the aggregate by means of a gyration diameter d g and a fractal dimension D f . Both values were further correlated with the number N p and the size d p of the primary particles and an additional constant prefactor k f .
These correlations were used to model the aggregates of different size and fractal dimension using a self-written algorithm; the application of which was previously described by Mayer et al. [23] The algorithm uses the tunable dimension method (TDM) which allows adding primary particles to an existing aggregate without significantly changing its fractal dimension. In this work, the numerical approach for the TDM is shown in Figure 2. Starting with a very small aggregate consisting of two primary particles, more primary particles are added to it until the desired aggregate size is obtained. To maintain the preset fractal dimension, a set of possible positions for the next-to-be-placed primary particle j (Figure 2b green-gray) is generated around an already existing primary particle i ( Figure 2 red-gray). Next, all possible positions which would lead to a physical overlap of the primary particles inside the aggregate are neglected (Figure 2 red). Afterward, the resulting fractal dimension for the remaining positions ( Figure 2 blue) is calculated. Finally, the new primary particle j is placed, so that the resulting deviation from the target fractal dimension is minimal or below a predefined threshold (Figure 2a green).
The accuracy of this method highly depends on a homogeneous distribution of the possible particle positions around the particle i. This step is way more challenging than it seems on the first glance and has already been tackled and discussed not only by mathematicians but also by meteorologists who have been seeking for a way to uniformly distribute satellites in Earth's orbit. Swinbank and Purser presented the Fibonacci grids as a possible solution. [24] This method allows a very homogeneous distribution of an odd amount of pointsN on a unit sphere. The longitudinal lon k and latitudinal lat k coordinates of the kth point are given by Equation (4a-c).

Stress Break Model
Now that a tool was created to order the primary particles with a predefined structural parameters, a proper model for the solid sinter bridges is to be found. As the forces created by these solid bridges are usually larger than other adhesive forces, especially in a fluid medium, the latter will be neglected. A common approach to model bridge forces in DEM was presented by Potyondy and Cundall, which they used to simulate rock and concrete fracturing under compressive stress. [25] As the solid-bridge-breaking behavior is a major point of this work, a short review of a few relevant features will be given. Similar to the Hertz-Mindlin model, the solid bridge forces are described as a parallel spring damper system. This results in additional bond forces between two bonded particles. The differentiation between tangential F t b and normal bond forces F n b is also made here. The breaking criteria are introduced as a critical tensile strength σ c and shear stress τ c , exceeding which results in breaking the bond between the two regarded particles. Therefore, the actual tensile stress σ and shear stress τ needs to be calculated for every contact during the simulation. In addition to the bond forces, the calculation includes the normal and tangential torques T n b and T t b , and area A, radius R b , moment of inertia I, and polar moment of inertia J of the parallel bond cross section Figure 2. Procedure of the extended TDM. Already existing particle i (red-gray), next placed particle j (green-gray), neglected possible particle positions (red), remaining possible positions (blue). a) Calculated position with the smallest deviation ΔD f (green). b) Resulting aggregate after particle placement.
The most difficult and yet a crucial aspect of this model is the calibration of the critical tensile and shear stresses. This step is especially challenging in nanoscaled CB aggregates due to various factors. First, the model assumes that all bonds have the same break criteria. However, due to their origin, the solid bridges between the primary particles of a CB aggregate are very likely to have an amorphous instead of a crystalline structure and are not necessarily all equally strong. Second, due to their small size, any measurement or imaging of the solid bridges at least causes enormous effort, if possible at all. As CB is available on the market with various different structural and mechanical properties, a simple yet efficient calibration method is required.
For this reason, a simplified calibration method was designed. The calibration is carried out by applying a controlled shear stress on a slurry with a high solid content of CB. This can be achieved, e.g., in a rheometer via a cone-and-plate setup. By analyzing the particle size distribution within the slurry before and after shearing, a few calibration points are created. Now, using a parameter study, the critical tensile and shear strength can be adjusted.

Fluid Forces
The simulation of fluid-laden flows can usually be approached in different ways. A straight-forward approach is the direct numerical simulation (DNS). In this method, the entire domain, including all present particles, is resolved using a very fine mesh. Afterward, the flow field is calculated and the mesh is updated according to particle movement. However, the computational costs of DNS is very high, even at low Reynolds numbers. If combined with a densely packed domain containing thousands of particles and high shear forces, the calculating times exceed far beyond the limits of what is reasonable for a large parameter study. Therefore, modeling of the fluid forces is inevitable.
Fortunately, depending on the requirements of the results, there are different approaches to modeling fluid forces. They can basically be categorized according to the considered interactions between the elements inside the domain. The simplest way is to only model the fluid forces on the particle: the One-Way-Coupling. This approach is usually used for very low solid content, as it assumes that the presence of a few, especially small particles, does not affect the flow pattern. Furthermore, any particle collisions are ignored, as the probability of a collision between two particles in this case is very low. With increasing solid content, the influence of particles on the fluid flow gains relevance and cannot be neglected.
A common approach for the simulation of particle-laden flows with higher solid content is the coupling of DEM with CFD. [26] Hereby, the particle movement (DEM) and the fluid flow (CFD) calculations are executed successively according to their time step Δt. As Δt DEM is usually a few orders of magnitude smaller than Δt CFD , multiple DEM calculation steps are executed before the fluid field is updated. Hence, the interaction between DEM and CFD is not necessarily executed after every CFD calculation as the coupling interval depends on particle inertia and fluid forces. However, when dealing with nanoscaled particles, high solid contents and high fluid forces, the direct CFD-DEM coupling results in stability issues and tremendously large computational effort.
Therefore, a simplified Two-Way-Coupling was applied. Within this approach, the shear-flow-induced fluid forces are calculated as an additional force in DEM. This way the aggregatebreaking behavior can be calculated on the microscale with an acceptable computational effort. The latter allows large studies with a larger range of material and process parameters. A similar approach was proposed by Higashitani et al. by introducing a 2D modified DEM (mDEM). [14] However, the success of this approach depends on the fluid force model.
The fluid drag forcesF D are generally described as the product of dynamic fluid pressure and cross-sectional area A of the object. Hereby, the dynamic fluid pressure is proportional to fluid density ρ f and squared relative velocityũ rel . An additional dimensionless drag coefficient C D is usually a function of the Reynolds number and considers the objects shape (Equation (7a)).
For small Reynolds numbers (Re < 1), the drag coefficient of a spherical particle can be approximated as 24=Re. This results in Stokes' law (Equation (7c)) which includes the dynamic viscosity η. [27] However, this approach assumes Newtonian behavior of the fluid while the viscosity of the slurry follows the Herschel-Bulkley law. Yet, as the authors consider a shear rate depending viscosity and the shear stress is a product of viscosity and shear rate, this assumption does not lead to significant errors. Furthermore, Stokes' law neglects the presence of other particles and their influence on the flow field. Although maintaining small Reynolds numbers on the microscale is not an issue, the high solid content of the slurry does not allow neglecting the influence of other present particles on the flow field.
Another approach was introduced by Steinour as Equation (7c) overestimates the actual fluid drag forces on a particle in a dense package. [28] Therefore, a dimensionless correction factor is introduced. With experimental investigation on sedimentation of uniform spheres Steinour determined this factor to be a function of local porosity ϵ local The latter can be obtained by subtracting the quotient of the solid volume V s (Figure 3a, blue) and the total volume V tot (Figure 3a, red) from one. The radius r a of the shell space around the regarded particle is to be set as twice the particle radius r p .
Before the porosity approach by Steinour was implemented into the DEM simulation, a few small-scale simulations using the Lattice Boltzmann Method [29][30][31][32] were executed on the generated aggregates. In this setup, a simple shear flow with a no-slip www.advancedsciencenews.com www.entechnol.de condition was applied on a cubic cell with an edge length of 1 μm filled with a CB aggregate. Afterward, the drag forces on every primary particle within the aggregate were calculated and compared with the drag forces calculated by Equation (8). The results are shown in Figure 3b. The correction function by Steinour shows a sufficient approximation for local porosity values up to 70%. The corrected coefficient of determination R 2 is relatively low with a value of 0.31 due to a large scatter of the correction factor. The deviation of the model increases with larger local porosity values. One major critique is that this approach neglects the flow direction relative to the particle positions. Therefore, effects such as slipstream are only partially included. This is very likely to be one significant reason for the large scatter of results. However, overall this approach provides a simple model with an acceptable quality.

Simulation Setup
With all three major issues being tackled, the parametric studies can be setup. The microscale simulations contain a small representative volume element ( Figure 4). The shear stress is applied as a symmetrical shear field (simple shear). The fluid force field within the element represents the simple shear flow pattern with a predefined shear stress γ : The previously generated aggregates with a predefined fractal dimension and size are randomly placed within the microscaled shear domain. The solid content within the cell varies slightly between 45% and 65% depending on aggregate size and fractal dimension. The applied shear stress varies between 3 and 40 kPa according to occurring conditions in a planetary mixer with cathode slurries with high solid contents as used by Mayer et al. [23] The aggregate size range is between 400 nm and 1 μm and with the fractal dimension preset to different values between 1.8 and 2.5. Also the domain is filled with CB aggregates only to isolate the disintegration due to fluid shear forces. The binding agent is assumed to be fully solved within the organic solvent. Hence, a representation of binder particles is redundant. Furthermore, the effects of present binding agent are partially considered with the fluid viscosity and increased damping of particle contacts as a representation of lubrication forces.

Aggregate Size
The simulations on the microscale reveal that with increasing shear stress the CB aggregate size converges exponentially toward the limit of x 50 % 200 nm. This behavior was observed for every initial aggregate size and shape ( Figure 5). However, larger aggregates reach this limit at significantly smaller shear stress. The critical breaking stress for large aggregates (x 0 50 ¼ 1 μm) with a fractal dimension of 2.5 is around 4 kPa, while that of aggregates with a fractal dimension of 2.2, 2, and 1.8 is significantly higher with 7.8 kPa. The differentiation between different aggregate shapes is even more pronounced at smaller aggregates with x 0 50 ¼ 400 nm. Here the aggregates with a fractal dimension of 2 show a maximal critical shear stress of 21.4 kPa, while those with a fractal dimension of 2.5 break already at 14.3 kPa.
We can assume that after a successful mixing process all aggregates have passed through the zone of maximum shear stress. Hence, all following results are based on the simulation setups with an applied shear stress of 40 kPa. Figure 6 shows the resulting size distribution of the CB fragments after shearing at maximum shear stress. First, it is noticeable that neither the  www.advancedsciencenews.com www.entechnol.de initial aggregate size nor its initial fractal dimension seem to have an impact on the fracture size distribution. As already shown in Figure 5, the median of all distributions tends toward 200 nm. Shearing larger aggregates (Figure 6, black) results in slightly smaller median values of the fragment size distributions. However, the absolute difference of 30 nm can be neglected. Especially, as this effect only appears at small fractal dimensions. Furthermore, the distribution width increases slightly with the initial aggregate size. This effect is also more significant at smaller initial fractal dimensions.
In this work, the aggregate size is calculated as the diameter of gyration d g . The latter depends not only on the spatial expansion of the aggregate but also on its shape. This leads to a fluctuation of the median aggregate size during the disintegration of CB aggregates. Sometimes these fluctuation could suggest an increasing aggregate size at higher shear stress. A similar issue is known from size measurement of irregular-shaped particles with light-scattering methods. However, this small relative error can be neglected. The resulting aggregate size agrees with the aggregate size measurements of the stressed suspension by Mayer et al. [23] A closer look at the aggregate size distribution of CB in the microscale simulations also shows a good agreement with the particle size distributions of CB obtained by focused ion beam-scanning electron microscopy (FIB-SEM) tomography which were presented in the same publication.
The original size of the aggregates has no significant influence on the resulting aggregate size at high shear stresses. However, smaller aggregates withstand higher shear stress. Due to their size, small aggregates offer a short lever for acting shear forces. As a result, small aggregates experience lower momentum than large aggregates at the same shear rate. Hence, an efficient disintegration of CB aggregates requires shearing large aggregates at maximum shear stress. Low shear stress shearing results in either slightly smaller fragments which require significantly higher shear stress for further disintegration or no disintegration at all. Consequently, the dispersing process of the Li-ion battery cathodes demands high shear stress and a geometry which transports the unstressed aggregates directly to zones of high shear stress.   Figure 7 shows the distribution of the fractal dimension after shearing CB aggregates with an initial size of 1 μm and 400 nm and an initial fractal dimension varying between 1.8 and 2.5. It is remarkable that the median of all shown distributions tends toward 2. Same applies for aggregates with an initial size of 600 and 800 nm. However, after shearing aggregates with an initial fractal dimension of 2.5, the median increases with decreasing initial aggregate size. Here the median of aggregates with an initial aggregate size of 1 μm, 800, 600, and 400 nm were 1.85, 2.04, 2.15, and 2.22, respectively. Similar to the results of the aggregate size shown in Figure 6, shearing larger aggregates leads to slightly wider distributions.

Aggregate Shape
As the number of primary particles of an aggregate of a fixed size d g increases with the fractal dimension as the exponent (see Equation (3)), breaking aggregates with high fractal dimension results in numerous differently shaped fragments. This leads to an increasing distribution width with increasing initial fractal dimension and/or increasing aggregate size.
The effect of the median approaching a value of 2 seems to be independent of the original size and fractal dimension of the aggregates (Figure 7). Furthermore, aggregates with an initial dimension near 2 withstand higher shear stress ( Figure 5). Both can be explained by the form of stress. The symmetrical shear flow which is simulated in this work generates a 2D stress field. The moments acting on an aggregate increase with its spatial expansion in the direction of the shear gradient. However, if the expansion in this direction is too small, the critical stress at the sinter bridges is not exceeded and the aggregates will not break any further. Therefore, if an aggregate is rod-shaped (D f % 1) or flat (D f % 2) and is oriented perpendicular to the shear gradient, the stress on this aggregate will be minimal. To verify this claim, a closer look on fragment orientation is necessary.

Orientation of the Fragments
The orientation of an aggregate requires a quantity which can describe its expansion into different spatial directions.
In this work, the authors propose the introduction of a dominant expansion vectorõ dom . The latter is the sum of all primary particle positionsr i relative to the center of gravityr cg The angle between the expansion vectorõ dom and the direction of the shear gradient is the orientation angle α. A schematic representation ofõ dom and α is shown in Figure 8. This additional quantity provides information about the alignment of the fragments after shearing. Figure 9 shows the cumulative distributions of the orientation angle α for individual CB fragments after shearing aggregates of different size and fractal dimension at maximum shear stress. Both median and modal value tend toward 90 for aggregates with an initial size of 600 nm and larger (Figure 9a). Only the fragments of aggregates with an initial size of 400 nm and a fractal dimension of 1.8, and 2.2 result in a median of 65 and 118, respectively ( Figure 9b).
As shown in Figure 9a, for shear cell setups with large aggregates, the median value of α approximates 90 . This indicates that the aggregates are almost parallel to the shear plane which corresponds to the assumption of favorable  alignment within the shear field. As free rotation of the aggregates is strongly impeded due to high solid content, the aggregates and their fragments are either already aligned after disintegration or continue to break until they are flattened by the stress. Only the median values of α for fragments of aggregates with an initial size of 400 nm and a fractal dimension of 1.8 and 2.2 deviate significantly from 90 . As small aggregates disintegrate in only few very small fragments which can withstand high shear stress, the orientation of the fragments can deviate from 90 without resulting in further disintegration.

Conclusion
The influence of shear stress in the range of 3-40 kPa on CB aggregates with different initial size and fractal dimension was investigated. The fragment size for increasing shear stress converges to a limit of 200 nm. This behavior was observed for every initial aggregate size and fractal dimension. However, small aggregates can withstand higher shear stress before breaking. Based on these findings, the authors recommend applying high shear stress on large aggregates for an efficient dispersion process. All results agree with previous experimental studies of Mayer et al. [23] The resulting fractal dimension of the fragments converges to a median of 2 without a significant influence of initial size and fractal dimension. The 2D form of stress in a simple shear cell was found to be a plausible cause for this results. The orientation of the fragments after shearing supports this claim. Yet, the presented results refer to one specific formulation and combination of compounds. The application of the results on other compositions requires further investigation. Nevertheless, these findings help to fill the gap in the research of aggregate behavior and provide a valuable contribution into optimized manufacturing process design. Future works aim at modeling CB aggregate disintegration during dispersion on a macro scale. For this purpose, all findings will be converted in a population balance kinetics and included in a macro-scale CFD simulation.

Experimental Section
The DEM simulations were executed with LIGGGHTS-PUBLIC, an OpenSource software which allows including new and extending existing models. Furthermore, this software includes multiple different unit systems which allows microscale simulations. The shear cell was reduced to a cube with an edge length of 3 μm and a package density between 0.45-0.65 depending on the aggregate shape and size. The DEM time step was set to 10 13 s. The overall simulation time was limited to 10 À4 s to match the largest time step which was used in CFD simulations of the considered mixer. Although these settings result in a tremendous computational effort, the obtained results of the study can be summed up in a set of breaking equations which will significantly improve the calculation times of macro scale simulations.
(a) (b) Figure 9. A cumulative distribution of the resulting orientation angle α after shearing the aggregates of different initial size and fractal dimension with high shear stress. a) Aggregates with an initial size of 1 μm. b) Aggregates with an initial size of 400 nm.