**1. Introduction**

The manufacturing sector has one of the largest impacts on worldwide energy use and natural resource consumption. Traditionally, research on manufacturing processes was mainly conducted to improve efficiency and accuracy and to lower costs [1]. One of the modern challenges in the field of designing manufacturing systems is to determine the optimal level of their flexibility from the point of view of the production tasks being performed. Whether or not a production process to be executed is capable of achieving the assumed performance parameters depends, among others, on the reliability of the machines and technological devices that make up the system under design [2]. As commonly known, the oldest manufacturing techniques used by humans are the grinding, chopper, and milling processes. Research on grinding tries to enhance economic and ecological properties and performance to extend grinding applications in the overall process chain—on the one hand, in the direction of increased material removal rates, avoiding turning and milling, and on the other hand, in the direction of fine finishing, thus making further abrasive finishing processes, such as lapping and polishing, obsolete [3].

Many factors influence the sustainability of a manufacturing process, and a combination of qualitative and quantitative methods are needed to discover them [4]. To resolve these issues, namely, the identification of the type of underlying impact factors, the uncertainty of the influencing factors, and the incompleteness of the evaluation information, the researchers used the evaluation method with analytical hierarchy process [5,6] and simulations [2,7–9] in the manufacturing sector. Sustainability encompasses the three pillars of economic, environmental, and social sustainability [10–12]. In the investigation of Ahmad [10], generally, there were approximately 44% (25/57) of highly applicable indicators. Among them, there were 28% (7/25) environmental, 28% (7/25) economic, and 44% (11/25) social indicators [10]. We can, therefore, see that economic and environmental factors are treated on an equal level. Social indicators in sustainable development are highly rated. Nevertheless, in the case of the two previous ones, they are more adequate for production processes, which have a more indirect nature on social impact. Beekaroo et al. [13] revealed an index with nine environmental, four economic, and two social indicators which were pertinent in sustainability measurement.

One of the most commonly used methods focusing on environmental sustainability is the product life cycle assessment (LCA). This method is used in many sectors of the economy [14–16] such as energy sector [17–19], transport sector [18,19], and food and agricultural industry [20,21]. The limitation of using this method is the need for specific and very detailed databases covering a lot of quantitative data.

Another method uses sustainability indicators (SI). In this method, the indicator can be defined "as a measure or aggregation of measures from which one can infer about the phenomenon of interest" [11]. This method can capture all three dimensions of sustainable development and help with the evaluation on many levels (e.g., enterprises, objects, processes, and products). In particular, for perpetrators with limited means and resources, SI provide a good method of analyzing sustainable development. Enterprises can assess their real situation using indicators, raise awareness, and set targets. The results of a literature review show that energy costs and GHG (Greenhouse Gas) indicators [10,11,13,19,22–25] are becoming the most commonly used indicators in sustainable planning. On this basis [11,26], four main directions of future research on sustainability indicators can also be mentioned: implementation of new optimization methods; adding further sustainable development indicators; extension of the model on a larger scale of the production system; and loosening certain assumptions. Finding optimal indicators in production processes is not an easy task; generally, we focus on minimizing energy consumption, waste production, improving the efficiency of machinery and equipment, as well as production and logistic processes in a given company. Sustainability indicators are based on measured and/or estimated data that have to be normalized, scaled, and aggregated consistently [27].

The evolution of design criteria for grinding and chopper machines is driven by functional requirements, general trends in machine tools, and cost [28]. The primary functional requirements, as named by Möhring et al. [29], are similar for all machine tools: high static and dynamic stiffness, fatigue strength, damping, thermal and long-term stability, and low weight of moving parts. The grinding and chopper process is carried out in many areas of food manufacturing [30–32] and agriculture [33–35], as well as in the industrial sector [36–38]. The mixing process is highly complicated, with a number of affecting parameters, such as the particle properties, the structure and performance of the mixer, the mixing process parameters, and the particle feeding order [38]. The crushing method and the machinery used for this purpose should be adapted to the type of material being ground and, in particular, to its mechanical properties [33,39].

Laboratories and modern industry require fast and effective grinding in small volumes [40]. In the field of agricultural mechanization, centrifugal chopper machines with a horizontal rotor have recently attracted increasing interest. Their principle of action lies in the acceleration of grains due to centrifugal inertia forces with their subsequent grinding (fragmentation) by impact or cutting [41,42]. This is due to the lower energy consumption for the grinding process in comparison with other choppers (crushers), which, in addition to the direct costs of grain destruction, have the energy costs for drifting whole and crushed grain using air. In centrifugal rotary chopper machines, the material can be supplied to crusher/shredder hammers or to a cutting pair by centrifugal inertia forces created by a horizontally rotating rotor with working elements. In this, one of the key roles of a centrifugal rotary material grinding machine is the distribution bowl or accelerator, which functions to provide a uniform and stable supply of material to the chopper bodies and impart to the material (particles) located on its surface the necessary linear velocity and trajectory of motion. Despite the rotary blenders relying upon the action of gravity to cause the powder to cascade and mix within a rotating vessel, the convective blenders employ a paddle, impeller, blade, or screw which stirs the powder inside a static vessel [43].

Most of the centrifugal rotary grinding machines in their design have a central feed of the material into the accelerator made in the form of a flat disc. However, with this material feed, the acceleration of a particle that is closest to the axis of rotation will be difficult because of the lack of initial velocity and centrifugal forces to overcome the friction forces, which provokes the appearance of a stagnant zone and an increase in the cost of overcoming them. To solve this problem, two approaches can be used: feeding the material with an offset from the axis of rotation, leading to design complexity, or using the reflective surface of the splitter, which is a straight conical surface located coaxially with the axis of rotation, ensuring the separation of the particles to be ground from the axis of rotation and giving them the initial velocity. In centrifugal impact grinders using the "stone-crushing-stone" principle with self-lining, the trajectory of the particles to be ground (milled) is equally important.

Figure 1 shows the rotary chopper machine analyzed in this study.

**Figure 1.** Tornado chopper machine: 1—feed-in tube; 2—case cover; 3—accelerator; 4—accelerating blade; 5—shredding element; 6—electric motor; 7—V-belt drive.

Figure 2 shows a scheme of self-lining formation in the accelerator.

**Figure 2.** Scheme of self-lining formation in the accelerator: 1—a splitter; 2—wear plate; 3—self-lining pocket; 4—carbide blade; 5—descent of the material from the accelerator; 6—accelerator case.

The main contribution of this study is the modeling (simulating) of the particle motion on an axisymmetric rotating curved surface with a vertical axis of rotation, whose calculation can be used for the study and design of machines, for example, centrifugal rotary chopper machines. A mathematical model of motion was developed on the basis of step-by-step numerical integration of the obtained closed system of differential equations with an equation of constraint describing the surface of revolution. An application based on the algorithm in C# under MS Visual Studio (Microsoft, Redmond, Washington, USA) was developed that enables graphical and numerical control of the calculation results. Moreover, the results of the calculations of the particle motion trajectory and kinematic indicators are presented. The conducted simulation tests show possibilities of improving the process efficiency and shortening the operation times, which results in economic benefits during sustainable manufacturing.

#### **2. Materials and Methods**

Hypothesis. One of the methods to achieve this goal takes into account the form of the accelerator, allowing us to set the initial kinematic parameters and particle trajectory, which in turn will affect the geometric dimensions, metal consumption, and technological modes of operation of a centrifugal rotary chopper or its nodes, for example, rotor (crusher) speed and material feed rate, which would then reduce the specific energy consumption per unit of output.

Methods. To solve this problem, we modeled the motion of a particle of the material on a rotating curved surface. To do this, we used the basic law of particle dynamics (dynamics of a material point) and considered the motion of a particle on the surface of an axisymmetric spinning bowl. We used step-by-step numerical integration, which allowed us to obtain data for analyzing the particle.

#### **3. Results and Discussion**

#### *3.1. Modeling Results*

The curved surface rotates around the vertical axis *z* (Figure 3). The equation of the surface is described in cylindrical coordinates by the equation

$$
\eta(\rho, \,\varphi, z) = 0 \tag{1}
$$

where ρ is the cylindrical radius, ϕ is the polar (vectorial) angle, and *z* is the application (*z*-axis).

The radius vector of a material point moving on a curved surface in the coordinate axes, rotating with it, is a function of these three coordinates, which can transform over time without breaking Equation (1):

$$
\overrightarrow{r} = \overrightarrow{r}(\rho, \varphi, z). \tag{2}
$$

The differential equation of the relative motion of a point on a rotating surface is written in vector form [40]

$$m\frac{d^2\overrightarrow{r}}{dt^2} = m\overrightarrow{\mathcal{g}} + \overrightarrow{N} + F\_{TR} + \overrightarrow{\Phi\_c} + \overrightarrow{\Phi\_c} \tag{3}$$

where *<sup>m</sup>* is the mass of a particle, *<sup>d</sup>*2<sup>→</sup> *r dt*<sup>2</sup> is its acceleration, *<sup>m</sup>*<sup>→</sup> *<sup>g</sup>* is the gravity force, <sup>→</sup> *N* is the normal reaction of a bowl surface, <sup>→</sup> *FTR* is the friction force from the surface directed opposite to the relative velocity of a particle, <sup>→</sup> <sup>Φ</sup>*<sup>e</sup>* is the centrifugal inertial force, and <sup>→</sup> Φ*<sup>c</sup>* is the Coriolis inertial force.

Forces and accelerations on the direction of the cylindrical axes of coordinates are projected. The gravity force is opposite to the *z* axis,

$$m\overrightarrow{\mathcal{g}} = m\overrightarrow{\mathcal{g}}(0, \ 0, \ -mg), \tag{4}$$

where <sup>→</sup> *g* is free-fall acceleration. The normal reaction of a rotating curved surface is

$$
\overrightarrow{N} = \lambda \cdot \text{grad}(\,\,\eta),
\tag{5}
$$

where λ = λ(*t*) is Lagrange's indeterminate multiplier [44], *grad*(η) is the vector gradient to the equation of surface (1), which has projections on the axial cylindrical coordinate system with the ϕ axis directed perpendicular to the ρ, *z*, axes and passing through the moving point, so that the axes ρ*,* ϕ*, z* form the right-hand system of vectors:

$$\begin{cases} \left( \operatorname{grad}(\eta) \right)\_{\rho} = \frac{\partial \eta}{\partial \rho} \\ \left( \operatorname{grad}(\eta) \right)\_{\rho} = \frac{\partial \eta}{\partial \psi} \cdot \frac{1}{\rho} . \end{cases} \tag{6}$$
 
$$\begin{cases} \left( \operatorname{grad}(\eta) \right)\_{z} = \frac{\partial \eta}{\partial z} \end{cases} \tag{6}$$

**Figure 3.** Design diagram of the motion of a particle on a rotating curved surface.

The projections of the normal reaction <sup>→</sup> *N* are:

$$\begin{cases} N\_{\rho} = \lambda \frac{\partial \eta}{\partial \rho} \\ N\_{\phi} = \lambda \frac{\partial \eta}{\partial \phi} \cdot \frac{1}{\rho} . \\ N\_{z} = \lambda \frac{\partial \eta}{\partial z} \end{cases} \tag{7}$$

The modulus of the normal reaction is found by

$$N = \left| \vec{\mathbf{N}} \right| = |\lambda| \sqrt{\left(\frac{\partial \eta}{\partial \rho}\right)^2 + \left(\frac{\partial \eta}{\partial \rho} \frac{1}{\rho}\right)^2 + \left(\frac{\partial \eta}{\partial z}\right)^2}. \tag{8}$$

At a constant angular velocity of rotation (spin rate) of the curved surface <sup>→</sup> ω, the centrifugal inertial force <sup>→</sup> Φ*<sup>e</sup>* is directed along the radius ρ in accordance with the result of vector products in its definition <sup>→</sup>

$$
\overrightarrow{\Phi\_{\varepsilon}} = -m\overrightarrow{\omega} \times \left(\overrightarrow{\omega} \times \overrightarrow{r}\right),
\tag{9}
$$

where <sup>→</sup> *r* is determined by Equation (2), then

$$\begin{cases} \Phi\_{\epsilon\wp} = m\nu^2 \rho \\ \quad \Phi\_{\epsilon\wp} = 0 \\ \quad \Phi\_{\epsilon z} = 0 \end{cases} . \tag{10}$$

Coriolis inertial force <sup>→</sup> Φ*<sup>c</sup>* is by definition equal to

$$
\overrightarrow{\Phi\_c} = -2m\overrightarrow{\omega} \times \overrightarrow{\upsilon}\_{\prime} \tag{11}
$$

where <sup>→</sup> *v* is the relative velocity of a material point in the moving axes of coordinates, and the modulus of relative velocity is associated with the cylindrical coordinates by

$$v = \sqrt{\dot{\rho}^2 + \left(\rho \dot{\phi}\right)^2 + \dot{z}^2}.\tag{12}$$

Then, the projections of the Coriolis force equal

$$\begin{cases} \Phi\_{c\rho} = 2ma\dot{\varphi}\rho\\ \Phi\_{c\rho} = -2m\dot{\rho}a \end{cases} \tag{13}$$

The friction force of a particle on the bowl surface <sup>→</sup> *FTR* is determined by Coulomb's law through a normal reaction and it is opposite in direction to the relative velocity:

$$
\overrightarrow{F\_{TR}} = -\left|\overrightarrow{\mathbf{N}}\right| f \frac{\overrightarrow{\boldsymbol{\upsilon}}}{\boldsymbol{\upsilon}} \tag{14}
$$

$$\begin{cases} \begin{aligned} F\_{TR\rho} &= -\overrightarrow{\text{N}} \Big| f \frac{\dot{\rho}}{v} \\ F\_{TR\rho} &= -\overrightarrow{\text{N}} f \frac{\dot{\varphi}\rho}{v} \\ F\_{TRz} &= -\overrightarrow{\text{N}} \Big| f \frac{\dot{z}}{v} \end{aligned} \end{cases} \tag{15}$$

Let us project Equation (3) on the axis of the moving system of cylindrical coordinates rotating together with the curved surface by using the found projections of all forces:

$$\begin{cases} m\left(\ddot{\rho} - \rho \dot{\varphi}^2\right) = \lambda \frac{\partial \eta}{\partial \rho} - Nf \frac{\dot{\rho}}{v} + m\omega^2 \rho + 2m\omega \dot{\varphi}\rho \\\ m \frac{1}{\rho} \frac{d}{dt} \{\rho^2 \dot{\varphi}\} = -Nf \frac{\dot{\varphi}\rho}{v} - 2m\dot{\rho}\omega \\\ m\ddot{z} = -mg + \lambda \frac{\partial \eta}{\partial z} - Nf \frac{\dot{z}}{v} \end{cases} \tag{16}$$

Let us transform the second equation of the system, taking into account the projection of acceleration on the ϕ axis:

$$\frac{1}{\rho} \frac{d}{dt} (\rho^2 \dot{\phi}) = 2\dot{\rho}\dot{\phi} + \rho \ddot{\phi}.\tag{17}$$

By bringing all forces to a unit mass of a particle *m* = 1, we obtain

$$\begin{cases} \ddot{\rho} = \rho \dot{\rho}^2 + \lambda \frac{\partial \eta}{\partial \rho} - Nf \frac{\dot{\rho}}{\overline{v}} + \omega^2 \rho + 2a \dot{\rho} \rho \\\ \ddot{\rho} = \frac{1}{\rho} (-2\dot{\rho}\dot{\rho} - Nf \frac{\dot{\rho}\rho}{\overline{v}} - 2\dot{\rho}\omega) \\\ \ddot{z} = -g + \lambda \frac{\partial \eta}{\partial \overline{z}} - Nf \frac{\dot{z}}{\overline{v}} \end{cases} \tag{18}$$

To obtain a closed system of equations, we supplement the system (18) with a coupling equation (of constraints) (1) (see Supplementary Material). Then, four unknown functions ρ(*t*), ϕ(*t*), *z*(*t*), λ(*t*) can be found by solving a system of three differential (18) and one algebraic (1) equations.

In the equation, various forms of execution of the surface of revolution are possible (1). Let us define the coupling equation in the form of an axisymmetric surface described by a power function passing through the given starting and final points of motion, written as

$$\frac{z}{z\_K} = a\_0 + \left(\frac{\rho}{\rho\_K}\right)^n,\tag{19}$$

where *a*<sup>0</sup> is a summand determined from the condition of belonging to a given surface of revolution of the initial point of the trajectory with the coordinates

$$
\rho(0) = \rho\_0, \ q(0) = 0, \ z(0) = 0 \tag{20}
$$

ρ*K*, *zK* are coordinates of the final point of the trajectory, where the surface described by Equation (19) ends for all values of the polar (vectorial) angle ϕ; *n* is an exponent in the equation of surface, which changes the degree of the bowl concavity and which can be varied to achieve the required parameters of the particle motion.

Thus, the origin of coordinates of the moving frame of reference, in which the particle motion is described, always corresponds along the *z* axis to the starting point of the motion trajectory, and the parabola vertex in the axial section of the surface of revolution lies on this axis below zero by the value *a*0:

$$a\_0 = -\left(\frac{\rho\_0}{\rho\_K}\right)^n. \tag{21}$$

Let us reduce Equation (19) to the form (1)

$$\eta(\rho,\rho,z) = \frac{z}{z\_{\mathbb{K}}} - a\_0 - \left(\frac{\rho}{\rho \mathbb{K}}\right)^n = 0,\tag{22}$$

and calculate the partial derivatives in the gradient expression to the constraint surface

$$\begin{cases} \frac{\partial \eta}{\partial \rho} = -n \left(\frac{\rho}{\rho \chi}\right)^{n-1} \frac{1}{\rho \chi} \\ \qquad \frac{\partial \eta}{\partial \varphi} = 0 \\ \qquad \frac{\partial \eta}{\partial z} = \frac{1}{z\_K} \end{cases} \tag{23}$$

These expressions in Equation (23) should be substituted in the projection of force <sup>→</sup> *N* , a normal reaction to the constraint surface in Equation (18).

From the last equation of the system (18) we find the Lagrange multiplier

$$
\lambda = \frac{1}{\left(\frac{\partial \eta}{\partial z}\right)} \left(\ddot{z} + g + Nf \frac{\dot{z}}{v}\right). \tag{24}
$$

Since, when moving along the surface, only two cylindrical coordinates are independent, then for step-by-step numerical integration, we take the first two differential Equations (18) as the basis and calculate the coordinate *z* and its time derivatives through coupling Equation (22):

$$z = z\_K a\_0 + z\_\mathcal{K} \left(\frac{\rho}{\rho\_\mathcal{K}}\right)^n \tag{25}$$

$$\dot{z} = nz\kappa \frac{1}{\rho\_{\mathbb{K}}} \left(\frac{\rho}{\rho\_{\mathbb{K}}}\right)^{n-1} \dot{\rho} \tag{26}$$

$$\bar{z} = nz\_{\bar{K}} \frac{1}{\rho\_{\bar{K}}} (n-1) \left(\frac{\rho}{\rho\_{\bar{K}}}\right)^{n-2} \frac{\bar{\rho}^2}{\rho\_{\bar{K}}} + nz\_{\bar{K}} \frac{1}{\rho\_{\bar{K}}} \left(\frac{\rho}{\rho\_{\bar{K}}}\right)^{n-1} \bar{\rho}. \tag{27}$$

In Equations (24)–(27), the coordinate values <sup>ρ</sup> and <sup>ϕ</sup>, as well as their time derivatives . <sup>ρ</sup>, . ϕ at the next integration step, are assumed to be equal to their value at the end of the previous integration step. The procedure of numerical integration of the system of Equation (18) was carried out by the method of averaged acceleration [45], using the original algorithmic program for "Windows" coded in C# under MS Visual Studio.

#### *3.2. Results and Discussion*

The results of the analysis of the kinematic parameters of the motion and trajectory of the particle, which initiate motion, with no initial linear velocity *v*<sup>0</sup> and with an angular velocity ω<sup>0</sup> equal in magnitude to the angular velocity of rotation ω*e*, for different curvilinear surfaces (Figure 4) of the accelerator, are presented in Figure 5. It shows that the use of a concave curvilinear surface (*na* = 2) allowed to double the value of the acquired velocity *vabs* in comparison with the flat surface (*na* = 0) and to increase it 1.7 times compared to a conical surface (*na* = −0.1) at a descent from a rotating surface, which is a consequence of the time spent by the particle on the surface, for example, when *na* = 2, the time was *t* = 2.05 s, when *na* = 0, the time was *t* = 0.84 s, and when *na* = −0.1, the time was *t* = 0.38 s. When the shape of the rotating surface changed from concave *na* = 2 to convex *na* = −0.1, the total length of the path traveled by the particle—that is, its trajectory—decreased.

**Figure 4.** Variants of generates of the studied curvilinear surfaces of the accelerator.

**Figure 5.** *Cont.*

**Figure 5.** Data for analyzing the velocity and trajectory of the particle motion, which is introduced with no initial linear velocity and with angular velocity equal in magnitude to the angular velocity of rotation for different curvilinear surfaces of the accelerator: (**a**) *na* = 2 (concave); (**b**) *na* = 0 (straight line); (**c**) *na* = −0.1 (conical curvilinear).

The calculations were carried out for different variable characteristics of the surface, initial indicators of kinematic parameters, and friction coefficient. Figure 6a presents an example of tabular data with the results of calculations and diagrams of changes in the motion trajectory, velocity, and acceleration of a particle along the cylindrical axes ρ and *z*.

Analysis of the diagrams presented in Figure 6b, c made it possible to note that the particle reached its maximum velocity in 0.313 s, being at a distance of ρ = 0.0092 m from the axis of rotation. At the same time, the subsequent changes in velocity, down to its descent from the concave surface, had a damped harmonic nature.

Thus, curvilinear surfaces of concave nature combining conical and concave surfaces are of great interest for further research. The trajectories of the particles and their kinematic indicators can be used to determine the relationship between the design factors and the process parameters of centrifugal choppers, such as the time when the particle is on a dispensing bowl, the descent velocity of the particle from the bowl, and the shape of accelerating blades.


**Figure 6.** *Cont.*

**Figure 6.** Example of the results of calculations of the particle motion trajectory and kinematic indicators for *na* = 2; (**a**) tabular data; (**b**) diagrams of particle displacement, velocity, and acceleration along the axis ρ; (**c**) diagrams of particle displacement, velocity, and acceleration along the axis *z*.

The ability to achieve the assumed performance parameters of the production process planned for implementation depends on the degree of reliability of the machines and technological devices included in the designed system [2]. Based on previous research [34], it can be confirmed that with the optimal values of the technological process (efficiency, time of a given operation—mixing or mining, coefficient of friction, and particle movement), it is possible to reduce energy consumption. Data on the possibility of improving the operation of production processes of various enterprises in Germany are presented by Steinhofel and others [46]. As we know, many factors influence the production process, for example, technical, economic, and social factors. Due to the above reasons, obtaining such large benefits in the whole production process, where various technological machines operate (technological operations), is difficult; nevertheless, it enables the introduction of sustainable production processes. The conducted simulation tests show possibilities of improving process efficiency and shortening the operation times, which results in economic benefits. The lower energy consumption and the increase in the efficiency of technological processes will also lead to the generation of less post-production waste (ecological benefits), which is in line with the principles of sustainable development.

#### **4. Conclusions**

In this study, a closed system of equations was obtained with a coupling equation (of constraints) that allows modeling a particle motion on a rotating surface; when numerically integrating them, it is possible to evaluate the use of a particular accelerator surface of a centrifugal rotary chopper machine to achieve the desired results in terms of downsizing, increasing productivity, or improving the quality of shredding. The use of numerical integration methods allows to consider the influence of physical and geometrical processes and kinematic parameters on the output result, for example, traveling time and velocity in cylindrical coordinates. Thus, the simulation of the particle motion

along an accelerator allows choosing rational geometrical dimensions of the chopper machine and its optimal operating practice.

**Supplementary Materials:** The Supplementary Materials are available online at http://www.mdpi.com/2071-1050/ 11/18/4873/s1.

**Author Contributions:** Conceptualization, J.C., P.A.S. and I.I.I.; Data curation, A.V.A.; Formal analysis, A.M., A.V.A., A.Y.I. and I.I.I.; Methodology, A.V.A., P.A.S. and A.Y.I.; Resources, J.C.; Software, P.A.S.; Writing—original draft, J.C., P.A.S. and A.Y.I.; Writing—review & editing, J.C. and I.I.I.

**Funding:** This research received no external funding.

**Acknowledgments:** We are thankful to the Editor and the reviewers for their valuable comments and detailed suggestions to improve the paper.

**Conflicts of Interest:** The authors declare no conflict of interest.
