# Calculation Method for Grading Size and Grading Bounding Box of Virtual Aggregate Based on DEM

Authors: Zhigang Zhou; Jun Ling
Journal: Mathematical Problems in Engineering
 (2021)

---

## Abstract

The grading size, which can be defined as the side length of the smallest square hole through which an aggregate with irregular shape can directly pass, is an important morphology parameter and can be used to calculate the gradation of mixture material. The grading bounding box, which can be defined as the circumscribed cuboid with central axis length being the grading size, is an important visual tool for observing the size and direction of the aggregate. Virtual test to calculate the grading size of a virtual aggregate is environmentally friendly and efficient, but the result provided by current research is imprecise and the grading bounding box is also rarely mentioned. In this paper, the multilevel complete projection algorithm is proposed to precisely calculate the grading size of a virtual aggregate. The whole process of the algorithm can be expressed by formula after the operation of sphere discretization by converting the virtual aggregate shell into the discrete aggregate. Then, the discrete aggregate is projected onto a complete series of the plane to form several 2D figures, and then, each 2D figure is projected onto a complete series of the orthogonal axis to form orthogonal segments. The grading size can finally be obtained by comparing the length of the above orthogonal segments based on the key central axis length principle. The influencing factors of computational accuracy and efficiency are considered in the algorithm. Finally, the grading bounding box can be built by using the Rodrigues transformation according to the information obtained from the above algorithm.


---

## Body

## 1. Introduction

The parameter of aggregate gradation in mixture material, such as asphalt mixture [1], is an important indicator that affects its whole mechanical properties. The sieving size of a single aggregate is an important factor in calculating the gradation of its mixture, which generally refers to the minimum standard side length of the square-hole sieve that can be passed.A 3D virtual aggregate shell can be simulated by professional model software, such as 3D modeling engine ACIS and Mimics [2–4], by combining a series of equidistant 2D images obtained by photographing or X-ray CT tomography. The grading size of the virtual aggregate shell is difficult to be directly obtained because of its irregular shape. There are many achievements in finding the grading size of virtual aggregates shell in various papers, which are roughly classified as follows:(1)
Virtual screening experiment: based on the screening principle, a large number of 3D virtual aggregates are put into the virtual screen model for the virtual screening experiment. Several sieve models, such as multilayer square-hole sieve [5], banana screen [6], horizontal linear screen [7], and double-amplitude screen with variable aperture [8], have been used. The method needs to simulate the sieve model first, and then, the vibration effect is achieved by applying velocity to the sieve model, which is time-consuming(2)
Simple measuring method: based on the vernier caliper principle, the virtual aggregates are measured one by one, taking the point from the surface and calculating circularly. The skeleton model is a way to find the point [9, 10]. The farthest two points are the long axis, perpendicular to which the farthest two points are the central axis, which is the final size of virtual aggregate [11]. But the long axis of a rectangle aggregate in this method is its diagonal, clearly out of a common sense(3)
Graphic fitting method: based on the constant eigenvalue principle, a virtual aggregate is converted into a regular geometric figure, whose characteristic size can be taken as the size of aggregate: the central axis of the fitting ellipsoid or the short axis of the fitting ellipse [12–14], the diameter of the minimum circumscribed or maximum inscribed sphere (circle) [15], the diameter of the sphere (circle) based on the equal volume (area) conversion, and the least second moment of figure [16](4)
Grading bounding box method: based on the key central axis length principle [17], the central side length of the grading bounding box is chosen as the grading size. The grading bounding box [2, 3] is searched by dynamically changing the three control angles of a circumscribed cuboid of the virtual aggregate, usually rotating a single angle in turn, fixing the other two angles, and making the cuboid reach the optimal value (the minimum central axis length or the minimum volume of the real-time bounding box). In addition, the minimum Ferrette diameter method [18] for determining the size of 2D figures is also included in this methodThe former three methods have the drawbacks of result-imprecise, especially the wrong principle of the simple measuring method. The grading bounding box method is a relatively precise manner with a complicated searching path, which is hard to understand and easy to miss the best value. In fact, the grading bounding box is better concerning the result, not the cause of the grading size calculation.


## 2. Objectives and Scopes

After synthesizing and sublating the above algorithms for measuring or calculating the grading size of virtual aggregate, and based on the discrete element method and the critical central axis size principle in a practical screening experiment, a multilevel complete projection algorithm is proposed.Firstly, in order to effectively describe the irregular virtual aggregate shell, the method of sphere discretization is used to convert it into a discrete aggregate, using enough reasonable small particles ranked in a certain order to fill in its inner space.Secondly, all the particles in one discrete aggregate are projected onto a complete series of verification planes to form 2D projection figures, and then, each 2D projection figure is projected on a complete series of orthogonal axes in the verification plane to form a group of orthogonal projection segments. The grading size of discrete aggregate can be calculated accurately by judging the length relations of all segments. To enforce this rule, necessary regulatory means, such as reasonable distribution of verification planes and their inner coordinate system, should be prepared in advance.Thirdly, the influencing factors, such as the sphere discretization method, inclusion standard, diameter of particles, the number of particles, and self-defined angle increment, should be considered to balance the precision and efficiency of calculation.Lastly, the grading bounding box is constructed to judge the direction and shape of the virtual aggregate according to the finalized controlling angles of the verification plane and the obliquity of its inner orthogonal axes, in which the projection got the final grading size.


## 3. Sphere Discretization

Any virtual aggregate can be expressed as a closed outer surface and its internal space. The sphere particles used in the discrete element method (DEM) have accurate characteristics of location and size, which is more motorized than the mesh used in the finite element method (FEM) [19]. Based on this, the paper proposes filling the internal space of the virtual aggregate to be measured with small-diameter particles, forming the discrete aggregate model, and then, the grading size can be calculated according to all inner particles that do not make a difference. The calculation process is also formulated and flowered due to the characteristics of particles and can be applied to all kinds of nondiscrete virtual aggregates (geometric shells) mentioned in the introduction. The discrete aggregate model obtained by the discrete element method can be directly used for subsequent calculation.The aggregate shell is introduced into the DEM software to simulate the size and shape of it with a large number of particles arranged according to certain rules, as shown in Figure1. The bubble pack method [18] uses a geometric shell to determine the area of the Delaunay triangulation. The cubic method and hexagonal method [20, 21] use the geometric shell to cover uniform particles which are not overlapped but arranged in cubic or hexagonal order, and the geometric shell does not participate in the operation. The filling method uses a geometric shell as the boundary of the wall, inside which uniform or uneven particles are generated and balanced. The contour-filling method [22] uses a geometric shell to define the center of the uniform particles to be generated.Figure 1
Sphere discretization of virtual aggregate. (a) Bubble pack method, (b) cubic method, (c) hexagonal method, (d) filling method, and (e) contour-filling method.
(a)(b)(c)(d)(e)


## 4. Multilevel Complete Projection Algorithm

The grading size of a discrete aggregate model, which is composed of several particles, can be calculated by a multilevel complete projection algorithm. The principles and procedures are as follows.


### 4.1. Selection of Verification Plane

The projection of a single discrete aggregate in different verification planes is different; all possibilities should be tested, so it is called the completeness of verification plane selection. Besides, it has the same projection of discrete aggregate on parallel planes with the same or opposite unit normal vector, so the direction of the plane is the only factor to determine the shape and size of the projection. Based on spherical coordinates, all verification planes can be described as the planes passing through a fixed pointOx0,y0,z0 with its unit normal vector n varying with the angle ϕ and θ, as shown in Figure 2:(1)n=xn,yn,zn=sinϕcosθ,sinϕsinθ,cosϕ,ϕ∈0,90∘,θ∈0,360∘.Figure 2
Spherical coordinate system.Given the discrete features of the verification plane, in order to avoid tedious calculation and to ensure the valid result, the verification plane is selected equally spaced: the angleϕ is selected in every Δ1 (degree), and the angle θ in every Δ2ϕ (degree) on the latitudinal circumference determined by ϕ, as shown in Figure 3. Self-defined Δ1 evenly divides the 1/4 meridian circle of the imaginary sphere whose radius is R′ into several equidistant arc L⌢, so on the latitudinal circle whose radius is R′ϕ determined by every angle ϕ, the verification plane ϕ,θ is taken from θ = 0° at intervals of the equidistant arc L⌢, which corresponds to Δ2ϕ as follows:(2)ϕ=kΔ1,k=0,1,2,…,u;uΔ1=90∘,L⌢=R′Δ1,Δ2ϕ=L⌢R′ϕ=L⌢R′sinϕ=Δ1sinϕ,θ=tΔ2ϕ,t=0,1,2,…,v;vΔ2<360∘,v+1Δ2≥360∘.Figure 3
Equal arc length selection of the verification plane.Particularly, whenk = 0, the unique verification plane that is parallel to xOy plane is taken.


### 4.2. Selection of 2D Coordinate System

It is necessary for the subsequent calculation to establish an appropriate 2D coordinate system in a plane, which passes the pointOx0,y0,z0 and whose unit normal vector is n=xn,yn,zn, describing the point on it.Let the only line in the plane, which is not only passing through the pointO but also parallel to the xOy plane, be the plane’s x-axis, which can be written as x2D and whose direction vector is L1=xL1,yL1,zL1. Let the only line in the plane, which is not only passing through the point O but also perpendicular to x2D, be the plane’s y-axis, which can be written as y2D and whose direction vector is L2=xL2,yL2,zL2. Particularly, choose L1=1,0,0 and L2=0,1,0 for the verification plane parallel to the xOy plane.The following set of equations can be listed based on the geometric properties:(1)
x2D⊥n:(3)L1⋅n=xL1xn+yL1yn+zL1zn=0.(2)
x2D passing through O:(4)x−x0xL1=y−y0yL1=z−z0zL1.(3)
x2D‖xOy plane:(5)zL1=0.(4)
y2D⊥n:(6)L2⋅n=xL2xn+yL2yn+zL2zn=0.(5)
y2D passing through O:(7)x−x0xL2=y−y0yL2=z−z0zL2.(6)
y2D⊥x2D:(8)L1⋅L2=xL1xL2+yL1yL2+zL1zL2=0.By simplifying, the direction vectors of thex-axis and y-axis in the plane, whose unit normal vector is n=xn,yn,zn, can be described in the 3D coordinate system:(9)L1=yn,−xn,0,L2=xnynzn,yn2zn,−xn2yn−yn3.For a certain verification planeϕi,θj chosen in the above section, the inner coordinate system x2DOy2Dϕi,θj can be described as (10)L1ϕi,θj=sinϕi×sinθj,−cosθj,0,(11)L2ϕi,θj=sin2ϕisinθj×cosϕicosθj,cosϕisinθj,−sinϕi.


### 4.3. Calculation of Projection Coordinates

Each particle in the discrete aggregate to be measured is presented for the calculation of the projection coordinates, as shown in Figure4. Let the coordinates of a particle’s center in the 3D coordinate system be AxA,yA,zA. Then, the values of the coordinates L1,A,L2,A of point B, which is the projection point of A in the verification plane ϕi,θj, can be calculated by the geometric vector rule.Figure 4
Plane projection principle.As shown in Figure5, in the established coordinate system x2DOy2Dϕi,θj, the orthogonal axes L3 and L4 for verification are built with the obliquity angle ωk=L3,L1∈0,90∘. The linear equations of the two axes are y2D=x2D⋅tanωk and y2D=x2D⋅−1/tanωk, respectively. The values of the coordinates L3,ωk,A,L4,ωk,A of a certain particle’s center A in the coordinate system L3OL4ωk determined by a specific obliquity angle of ωk is transformed from its projection coordinate values L1,A,L2,A in the coordinate system x2DOy2D:(12)L3,ωk,AL4,ωk,A=cosωksinωk−sinωkcosωkL1,AL2,A.Figure 5
Orthogonal axis projection principle.


### 4.4. Grading Size Calculation of Projection

After transforming all the coordinate values of particles in the discrete aggregate on currentωk in current ϕi,θj, the range length Lωk,L3or4¯ in the orthogonal axes L3 or L4 can be found by firstly finding the coordinate values reached most upwards and downwards, respectively, in the axes L3 and L4, and then, taking a subtraction with the consideration of particle’s radius, as the bold line in axes shown in Figure 5,(13)Lωk,L3or4¯=maxL3or4¯,ωk,Am+rAm−minL3or4¯,ωk,An−rAn.The side lengthLωkϕi,θj of the smallest square, through which the projection of the discrete aggregate on the obliquity angle of ωk in the verification plane of ϕi,θj can pass, can be found by choosing the maximum range length:(14)Lωkϕi,θj=maxLωk,L3,Lωk,L4,where Am and An stand for arbitrary qualified particles inside the discrete aggregate, and rAm and rAn stand for the radius of the corresponding particle.


### 4.5. Multilevel Complete Projection

The specific obliquity angle ofωk is generalized to take the value in the range ω∈0,90∘ and calculate all Lωkϕi,θj values in every self-defined Δ3 in the specific verification plane ϕi,θj. Then, the side length of the smallest square Lϕi,θj, through which the projection of the discrete aggregate in the specific verification plane of ϕi,θj can pass, is calculated:(15)Lϕi,θj=minLω1ϕi,θj,Lω2ϕi,θj,…,.The specific verification planeϕi,θj is generalized to take the value in the range ϕ∈0,90∘, θ∈0,360∘, and calculate all Lϕi,θj values in every self-defined Δ1 and Δ2ϕ. Then, the side length of the smallest square L, through which the discrete aggregate can pass, is calculated:(16)L=minLϕ1θ1,Lϕ1θ2,…,Lϕ2θ1,…,.The calculation flowchart of grading size for a discrete aggregate is shown in Figure6.Figure 6
Flowchart of grading size calculation of discrete aggregate.


## 4.1. Selection of Verification Plane

The projection of a single discrete aggregate in different verification planes is different; all possibilities should be tested, so it is called the completeness of verification plane selection. Besides, it has the same projection of discrete aggregate on parallel planes with the same or opposite unit normal vector, so the direction of the plane is the only factor to determine the shape and size of the projection. Based on spherical coordinates, all verification planes can be described as the planes passing through a fixed pointOx0,y0,z0 with its unit normal vector n varying with the angle ϕ and θ, as shown in Figure 2:(1)n=xn,yn,zn=sinϕcosθ,sinϕsinθ,cosϕ,ϕ∈0,90∘,θ∈0,360∘.Figure 2
Spherical coordinate system.Given the discrete features of the verification plane, in order to avoid tedious calculation and to ensure the valid result, the verification plane is selected equally spaced: the angleϕ is selected in every Δ1 (degree), and the angle θ in every Δ2ϕ (degree) on the latitudinal circumference determined by ϕ, as shown in Figure 3. Self-defined Δ1 evenly divides the 1/4 meridian circle of the imaginary sphere whose radius is R′ into several equidistant arc L⌢, so on the latitudinal circle whose radius is R′ϕ determined by every angle ϕ, the verification plane ϕ,θ is taken from θ = 0° at intervals of the equidistant arc L⌢, which corresponds to Δ2ϕ as follows:(2)ϕ=kΔ1,k=0,1,2,…,u;uΔ1=90∘,L⌢=R′Δ1,Δ2ϕ=L⌢R′ϕ=L⌢R′sinϕ=Δ1sinϕ,θ=tΔ2ϕ,t=0,1,2,…,v;vΔ2<360∘,v+1Δ2≥360∘.Figure 3
Equal arc length selection of the verification plane.Particularly, whenk = 0, the unique verification plane that is parallel to xOy plane is taken.


## 4.2. Selection of 2D Coordinate System

It is necessary for the subsequent calculation to establish an appropriate 2D coordinate system in a plane, which passes the pointOx0,y0,z0 and whose unit normal vector is n=xn,yn,zn, describing the point on it.Let the only line in the plane, which is not only passing through the pointO but also parallel to the xOy plane, be the plane’s x-axis, which can be written as x2D and whose direction vector is L1=xL1,yL1,zL1. Let the only line in the plane, which is not only passing through the point O but also perpendicular to x2D, be the plane’s y-axis, which can be written as y2D and whose direction vector is L2=xL2,yL2,zL2. Particularly, choose L1=1,0,0 and L2=0,1,0 for the verification plane parallel to the xOy plane.The following set of equations can be listed based on the geometric properties:(1)
x2D⊥n:(3)L1⋅n=xL1xn+yL1yn+zL1zn=0.(2)
x2D passing through O:(4)x−x0xL1=y−y0yL1=z−z0zL1.(3)
x2D‖xOy plane:(5)zL1=0.(4)
y2D⊥n:(6)L2⋅n=xL2xn+yL2yn+zL2zn=0.(5)
y2D passing through O:(7)x−x0xL2=y−y0yL2=z−z0zL2.(6)
y2D⊥x2D:(8)L1⋅L2=xL1xL2+yL1yL2+zL1zL2=0.By simplifying, the direction vectors of thex-axis and y-axis in the plane, whose unit normal vector is n=xn,yn,zn, can be described in the 3D coordinate system:(9)L1=yn,−xn,0,L2=xnynzn,yn2zn,−xn2yn−yn3.For a certain verification planeϕi,θj chosen in the above section, the inner coordinate system x2DOy2Dϕi,θj can be described as (10)L1ϕi,θj=sinϕi×sinθj,−cosθj,0,(11)L2ϕi,θj=sin2ϕisinθj×cosϕicosθj,cosϕisinθj,−sinϕi.


## 4.3. Calculation of Projection Coordinates

Each particle in the discrete aggregate to be measured is presented for the calculation of the projection coordinates, as shown in Figure4. Let the coordinates of a particle’s center in the 3D coordinate system be AxA,yA,zA. Then, the values of the coordinates L1,A,L2,A of point B, which is the projection point of A in the verification plane ϕi,θj, can be calculated by the geometric vector rule.Figure 4
Plane projection principle.As shown in Figure5, in the established coordinate system x2DOy2Dϕi,θj, the orthogonal axes L3 and L4 for verification are built with the obliquity angle ωk=L3,L1∈0,90∘. The linear equations of the two axes are y2D=x2D⋅tanωk and y2D=x2D⋅−1/tanωk, respectively. The values of the coordinates L3,ωk,A,L4,ωk,A of a certain particle’s center A in the coordinate system L3OL4ωk determined by a specific obliquity angle of ωk is transformed from its projection coordinate values L1,A,L2,A in the coordinate system x2DOy2D:(12)L3,ωk,AL4,ωk,A=cosωksinωk−sinωkcosωkL1,AL2,A.Figure 5
Orthogonal axis projection principle.


## 4.4. Grading Size Calculation of Projection

After transforming all the coordinate values of particles in the discrete aggregate on currentωk in current ϕi,θj, the range length Lωk,L3or4¯ in the orthogonal axes L3 or L4 can be found by firstly finding the coordinate values reached most upwards and downwards, respectively, in the axes L3 and L4, and then, taking a subtraction with the consideration of particle’s radius, as the bold line in axes shown in Figure 5,(13)Lωk,L3or4¯=maxL3or4¯,ωk,Am+rAm−minL3or4¯,ωk,An−rAn.The side lengthLωkϕi,θj of the smallest square, through which the projection of the discrete aggregate on the obliquity angle of ωk in the verification plane of ϕi,θj can pass, can be found by choosing the maximum range length:(14)Lωkϕi,θj=maxLωk,L3,Lωk,L4,where Am and An stand for arbitrary qualified particles inside the discrete aggregate, and rAm and rAn stand for the radius of the corresponding particle.


## 4.5. Multilevel Complete Projection

The specific obliquity angle ofωk is generalized to take the value in the range ω∈0,90∘ and calculate all Lωkϕi,θj values in every self-defined Δ3 in the specific verification plane ϕi,θj. Then, the side length of the smallest square Lϕi,θj, through which the projection of the discrete aggregate in the specific verification plane of ϕi,θj can pass, is calculated:(15)Lϕi,θj=minLω1ϕi,θj,Lω2ϕi,θj,…,.The specific verification planeϕi,θj is generalized to take the value in the range ϕ∈0,90∘, θ∈0,360∘, and calculate all Lϕi,θj values in every self-defined Δ1 and Δ2ϕ. Then, the side length of the smallest square L, through which the discrete aggregate can pass, is calculated:(16)L=minLϕ1θ1,Lϕ1θ2,…,Lϕ2θ1,…,.The calculation flowchart of grading size for a discrete aggregate is shown in Figure6.Figure 6
Flowchart of grading size calculation of discrete aggregate.


## 5. Influencing Factors of Calculation Results and Efficiency

An example is shown in Figure7; a geometric shell of the octahedron is surrounded by eight random planes. The unit normal vector and the distance from the specified point inside the polyhedron to each plane direction are shown in Table 1.Figure 7
Three-view diagram of the polyhedral geometric shell. (a) Front view, (b) side view, and (c) top view.
(a)(b)(c)Table 1
The plane parameters of the polyhedral geometric shell.
PlaneUnit normal vectorDistance (mm)1(0.512, −0.422, 0.748)4.7912(−0.482, 0.396, 0.782)6.3363(0.020, −0.98, −0.183)5.6854(−0.922, −0.380, −0.074)6.9365(0.928, 0.347, −0.133)3.9866(0.090, 0.995, −0.047)5.1747(0.311, −0.332, −0.891)6.4638(−0.465, 0.571, −0.676)4.190The sphere discretization of the aggregate shell is represented by different conditions:(1)
Sphere discretization method: cubic method (usingball generate cubic … command) or filling method (using ball distribute porosity 0.4 … command and cycle 800 calm 100 command), as shown in Figure 8(2)
Inclusion standard (only for cubic method): center inclusion (as long as the sphere center of the particle is inside the polyhedral shell, it can be included) or whole inclusion (as long as the whole particle is located in the polyhedral shell, it can be included)(3)
The diameter of particles: 0.6 mm, 1.0 mm, or 1.5 mmFigure 8
Sphere discretization method. (a) Cubic method and (b) filling method.
(a)(b)By combining, the conditions of discrete aggregate after deleting the wall are formed, as shown in Figure9.Figure 9
Discrete angular aggregate model with different sphere discretization conditions. (a) Cubic method + center inclusion + 0.6 mm diameter. (b) Cubic method + center inclusion + 1.0 mm diameter. (c) Cubic method + center inclusion + 1.0 mm diameter. (d) Cubic method + center inclusion + 1.5 mm diameter. (e) Filling method + 0.6 mm diameter. (f) Filling method + 1.0 mm diameter. (g) Filling method + 1.5 mm diameter.
(a)(b)(c)(d)(e)(f)(g)The grading sizes of the above discrete aggregates are calculated by takingΔ1 = 5° and Δ3 = 5°, as shown in Table 2.Table 2
Grading size of discrete aggregate in different sphere discretization conditions.
ConditionL (mm)ϕaθbωcConditionL (mm)ϕaθbωcΔL (mm)(a)12.48335°104.61°5°(e)12.30230°100°10°0.181(b)12.68735°104.61°5°(f)12.44230°100°10°0.321(c)11.57135°104.61°5°——————(d)13.00725°106.48°0°(g)12.33140°116.68°85°0.676The smaller the diameter of the particle is, the more accurate the final calculation results are. Condition (a) and condition (e) can be used as more accurate models for the corresponding discretization method.The difference of grading size between center inclusion (condition (b)) and sphere inclusion (condition (c)) is about the diameter of the particle, the former is closer to the value of the condition (a), and all the angle variables are consistent, so sphere center inclusion is a relatively reasonable manner.When the particle diameter is fixed, the measurement value of the filling method is slightly less than that of the cubic method, which is also because the center inclusion is chosen so that part of the particle is out of the plane, while the wall of the filling method blocks the particles reaching out to some extent, which is shown in Figure8. Taking into account the fact that the particle number of the cubic method is slightly less and does not need to run a balance, it can be used in the rough calculation.When the discretization method is fixed, the measurement precision decreases with the increase of diameter. When the diameter reaches 1.5 mm, the measurement value of the filling method (condition (g)) decreases, but that of the cubic method (condition (d)) obviously increases, which is also due to the blocking effect of the wall in the filling method, while the particle diameter has a great influence on the cubic method. Besides, there is a deviation from the angle variable in both methods.In the rough calculation of grading size, the combination of the cubic method and 1 mm particle diameter can be used with accepted error. In fact, selecting 10% of the estimated size value as the diameter can obtain good results, especially when there are subsequent mechanical calculations, while the shape of the original polyhedron shell is required.In the fine calculation of grading size, the filling method should be adopted and the smaller diameter should be selected, but this will prolong the calculation time, which is proportional to the following conditions:(1)
The number of particlesnump: the bigger the diameter of particles, the less the total number of particles, but the rougher the shape it will eventually take, so the diameter of particles cannot be increased casually. Given that the inner particles contribute almost nothing to the final result, all the particles whose distance from the geometric shell is less than the diameter of particles can be grouped separately and the others can be deleted. The number of particles whether the entire or grouped is shown in Table 3.(2)
The number of verification planesnumv: numv is inversely proportional to Δ12. When Δ1 is taken as 1°, 3°, and 5°, numv is calculated as 20852, 2367, and 869, respectively.(3)
The number of orthogonal axes on each verification planenuma: numa is inversely proportional to the increment Δ3. When Δ3 is taken as 1°, 3°, and 5°, numa is calculated as 90, 30, and 18, respectively.Table 3
The total number of entire or grouped particles in a discrete aggregate model.
ConditionEntireGroupeda49871524b1088509c791407d325209e61871731f1364592g405244Taking condition (a) and condition (e) as an example, the discrete aggregate is run by taking different angle incrementsΔ1 and Δ3, as shown in Tables 4 and 5.Table 4
Grading size of discrete aggregate using the cubic method in different angle increments of projection.
Δ1Δ3L (mm)ϕaθbωcΔ3L (mm)ϕaθbωcΔ3L (mm)ϕaθbωc1°1°12.45230°104.00°5°3°12.45433°102.82°6°5°12.45230°104.00°5°3°1°12.46930°102.12°7°3°12.50733°110.16°0°5°12.46133°104.66°5°5°1°12.48335°104.60°5°3°12.49430°110.00°0°5°12.48335°104.60°5°Table 5
Grading size of discrete aggregate using the filling method in different angle increments of projection.
Δ1Δ3L (mm)ϕaθbωcΔ3L (mm)ϕaθbωcΔ3L (mm)ϕaθbωc1°1°12.22733°102.82°7°3°12.23732°103.79°6°5°12.25031°104.85°5°3°1°12.25230°102.00°8°3°12.30239°104.88°6°5°12.27433°104.66°5°5°1°12.29635°104.60°6°3°12.29635°104.60°6°5°12.30230°100°10°The results still confirm that the filling method is more accurate than the cubic method. Besides, the smaller the angle increment is, the more accurate the measured value is most of the time, but more computation is needed. Due to the discrete nature of particles, occasionally prominent or concave particles at a certain angle make measurement results not smooth, as occurs inΔ1 = 3° and Δ3 = 3° in the filling method. The maximum difference of grading size in the two methods is, respectively, 0.031 mm and 0.075 mm, which is much less than the particle diameter, so taking the two angle increment values as 5° can be accepted most of the time for less running time.For the convenience of observation, the 2D projection of the discrete aggregate in their corresponding verification plane in the increment ofΔ1 = 1° and Δ3 = 1° is shown in Figure 10.Figure 10
The 2D projection of discrete aggregates in the verification plane. (a) Cubic method and (b) filling method.
(a)(b)


## 6. Grading Bounding Box

The grading bounding box of a discrete aggregate can be defined as the cuboid, whose six surfaces are both circumscribed to the discrete aggregate without any part going beyond, and which has the same grading size as the discrete aggregate. Let the three side lengths of a grading bounding box beH≥L≥W in descending order; then, its central axis length L is necessarily the well-calculated grading size of its inner discrete aggregate.The grading bounding box can be used to simplify the judgment of the direction and size of the inner discrete aggregate. (1) The normal vector direction of the plane in which the two edgesLW are located represents the orientation of the discrete aggregate. (2) When W/H<0.4, the discrete aggregate can be judged as needle-like aggregate (L << H) or flake-like aggregate (L >> W).Let a discrete aggregate find its grading sizeL on the obliquity of ωc in the verification plane ϕa,θb, the short axis length and the long axis direction of its grading bounding box can be calculated by equations (17) and (18), and the short axis direction, the central axis direction, and the long axis length can be obtained by applying the Rodríguez transform:(17)W=minLωcϕa,θb,L3,Lωcϕa,θb,L4,(18)nLW=±1×sinϕacosθb,sinϕasinθb,cosϕa.


### 6.1. Rodríguez Transform

The Rodríguez transform [23–25] constructs a cross-product matrix K and rotation matrix P, and then, the vector n1 is rotated around the axis of k to the vector n2 by the degree of δ according to the right-hand rule, as shown in Figure 11:(19)K=0−kzkykz0−kx−kykx0,P=I+sinδK+1−cosδK2,n2=Pn1,where k=kx,ky,kzT is a unit vector and I is a third-order unit matrix.Figure 11
Rodríguez transform.As shown in Figure12, when L1 and L2 in equations (10) and (11) are rotated around nLW by degree of ωc according to the Rodríguez transform, L3 and L4 can be obtained in the 3D coordinate system by taking k=nLW, δ=ωc, and n1=L1 and L2, respectively:(20)L3=sinϕa×sinθbcosωc+cosϕacosθbsinωc,cosϕasinθbsinωc−cosθbcosωc,−sinϕasinωc,L4=sin2ϕasinθb×cosϕacosθbcosωc−sinθbsinωc,cosθbsinωc+cosϕasinθbcosωc,−sinϕacosωc.Figure 12
Parameters of the grading bounding box.After being unitized,L3 and L4 can be taken as the unit normal vector of the plane in which the two edges HW or HL are located: if the central axis length L is obtained on the axis L4, then the short axis direction nHL=±1×L3/L3, the central axis direction nHW=±1×L4/L4, and vice versa. The long axis length can be obtained by projecting algorithm, choosing the higher value range in the axes, either nHL or nHW, combined with nLW, while the lower one is the known L or W for verification.


### 6.2. Application

According to the process, the grading bounding box of two discrete aggregate models mentioned above can be constructed: three controlling anglesϕ = 30°, θ = 104°, and ω = 5° in the cubic method and ϕ = 33°, θ = 102.82°, and ω = 7° in the filling method are used to obtain two groups of grading bounding box parameters, as shown in Table 6.Table 6
Grading bounding box parameters of the discrete aggregate model.
ParametersCubic methodFilling methodH/L/W (mm)15.817/12.452/12.43116.187/12.227/12.190nLW± (−0.121, 0.485, 0.866)± (−0.121, 0.531, 0.839)nHL± (−0.293, 0.816, −0.498)± (−0.304, 0.785, −0.541)nHW± (0.948, 0.314, −0.044)± (0.945, 0.320, −0.066)With the help of the view tool, a targeted observation can be made: Figures13 and 14 show the position relationship between the grading bounding box and the discrete aggregate is observed from the original angle of view and from the long axis direction, respectively. The former objectively reflects the whole orientation of the discrete aggregate in the space, and the latter is corresponding to the 2D projection image obtained in Figure 10, which verifies the rationality of the whole process.Figure 13
Construction of grading bounding box. (a) Cubic method and (b) filling method.
(a)(b)Figure 14
The discrete aggregate in the grading bounding box observed in the long axis direction. (a) Cubic method and (b) filling method.
(a)(b)Figure15 shows a three-view diagram of the polyhedral shell enclosed in its grading bounding box obtained by the filling method. Because the discrete aggregate as the medium inherits the shape, size, and orientation of the polyhedral shell to a large extent, the grading bounding box is basically suitable for the polyhedral shell.Figure 15
Three-view diagram of polyhedral shell enclosed in the grading bounding box. (a) Front view, (b) side view, and (c) top view.
(a)(b)(c)


## 6.1. Rodríguez Transform

The Rodríguez transform [23–25] constructs a cross-product matrix K and rotation matrix P, and then, the vector n1 is rotated around the axis of k to the vector n2 by the degree of δ according to the right-hand rule, as shown in Figure 11:(19)K=0−kzkykz0−kx−kykx0,P=I+sinδK+1−cosδK2,n2=Pn1,where k=kx,ky,kzT is a unit vector and I is a third-order unit matrix.Figure 11
Rodríguez transform.As shown in Figure12, when L1 and L2 in equations (10) and (11) are rotated around nLW by degree of ωc according to the Rodríguez transform, L3 and L4 can be obtained in the 3D coordinate system by taking k=nLW, δ=ωc, and n1=L1 and L2, respectively:(20)L3=sinϕa×sinθbcosωc+cosϕacosθbsinωc,cosϕasinθbsinωc−cosθbcosωc,−sinϕasinωc,L4=sin2ϕasinθb×cosϕacosθbcosωc−sinθbsinωc,cosθbsinωc+cosϕasinθbcosωc,−sinϕacosωc.Figure 12
Parameters of the grading bounding box.After being unitized,L3 and L4 can be taken as the unit normal vector of the plane in which the two edges HW or HL are located: if the central axis length L is obtained on the axis L4, then the short axis direction nHL=±1×L3/L3, the central axis direction nHW=±1×L4/L4, and vice versa. The long axis length can be obtained by projecting algorithm, choosing the higher value range in the axes, either nHL or nHW, combined with nLW, while the lower one is the known L or W for verification.


## 6.2. Application

According to the process, the grading bounding box of two discrete aggregate models mentioned above can be constructed: three controlling anglesϕ = 30°, θ = 104°, and ω = 5° in the cubic method and ϕ = 33°, θ = 102.82°, and ω = 7° in the filling method are used to obtain two groups of grading bounding box parameters, as shown in Table 6.Table 6
Grading bounding box parameters of the discrete aggregate model.
ParametersCubic methodFilling methodH/L/W (mm)15.817/12.452/12.43116.187/12.227/12.190nLW± (−0.121, 0.485, 0.866)± (−0.121, 0.531, 0.839)nHL± (−0.293, 0.816, −0.498)± (−0.304, 0.785, −0.541)nHW± (0.948, 0.314, −0.044)± (0.945, 0.320, −0.066)With the help of the view tool, a targeted observation can be made: Figures13 and 14 show the position relationship between the grading bounding box and the discrete aggregate is observed from the original angle of view and from the long axis direction, respectively. The former objectively reflects the whole orientation of the discrete aggregate in the space, and the latter is corresponding to the 2D projection image obtained in Figure 10, which verifies the rationality of the whole process.Figure 13
Construction of grading bounding box. (a) Cubic method and (b) filling method.
(a)(b)Figure 14
The discrete aggregate in the grading bounding box observed in the long axis direction. (a) Cubic method and (b) filling method.
(a)(b)Figure15 shows a three-view diagram of the polyhedral shell enclosed in its grading bounding box obtained by the filling method. Because the discrete aggregate as the medium inherits the shape, size, and orientation of the polyhedral shell to a large extent, the grading bounding box is basically suitable for the polyhedral shell.Figure 15
Three-view diagram of polyhedral shell enclosed in the grading bounding box. (a) Front view, (b) side view, and (c) top view.
(a)(b)(c)


## 7. Summary and Conclusions

(1)
Sphere discretization is an effective means to tackle a virtual aggregate with irregular shapes. After transferring the whole geometric shell into discrete aggregate, which consists of lots of particles ranging in a certain order, the grading size and the grading bounding box can be calculated by tracing each particle’s exact position. The discrete aggregate inherits the shape, size, and direction of the original virtual aggregate, and the whole process can be formulated and flowed(2)
Multilevel complete projection algorithm can be used to calculate the grading size of discrete aggregate with a relatively accurate result. After a complete projection of the three dimensions, the grading size can be obtained by easily comparing the lengths of a series of line segments. Meanwhile, the following additional results are achieved: equal arc angle equation and the method of establishing coordinate system on any plane(3)
The filling method is more accurate to calculate the grading size than the cubic method but less efficient. A more accurate result can also be obtained by decreasing the particle diameter and the angle incrementΔ1 and Δ3, but both of them are time-consuming. When there is a follow-up mechanical calculation or low precision requirement, it is recommended that the particle diameter be 10% of the estimated size with a relatively precise shape and fewer particles(4)
Grading bounding box can be built to find the long axis direction and axis length of the virtual aggregate by using the data generated from the above algorithm and by realizing the vector rotation according to the Rodríguez transform



---

*Source: 6659285-2021-06-27.xml*