Development of discrete element simulation method of vegetable oil press

During the process of plant oil extraction major part of the power input transforms into thermal energy caused by friction. The evolved heat has a very significant effect on the quality of the oil. Our aim is to approach the real pressure and thermal phenomenon with the support of discrete element method (DEM) simulations. In the present stage of our research we have created a model with more complex geometry which represents the next stage of the evolutionary model development. The article presents our newest results.


Introduction
The plant oil extraction is a process with high energy consumption.Major part of the delivered power transforms into heat.This general phenomenon caused by friction at the contacts of the seeds with each other and also with the mechanical parts.The evolved heat has a very significant effect on the oil.The heat load can modify the quality of the oil by starting chemical alterations inside the material.
In case of an oil press with worm shaft, from the temperature aspect we distinguish two basic methods; the cold and hot pressing methods.Using cold pressing procedure over the clarity of oil seeds there is no need for any further treating.According to experiences cold pressing makes possible to get oil with low phosphor content which does not require any further processing to produce biodiesel.In general, according to the class and type of the product there are limit values for the temperature to get the requested quality.
So the process needs the information of the heat distribution inside the press machine to be sure about the parameters of the extracted plant oil.

Discrete element method
Discrete element method (DEM) is a discontinuous numerical method based on molecular dynamics.It was developed and applied for analysing rock mechanics [1].Former research has shown that virtual DEM simulations were developed in correspondence with real tests [2].It is the most proper way to simulate the behaviour of materials consisting of seeds.With the support of DEM it is possible to simulate the mechanical and thermal relations in the oil press and even pre-indicate the quality of the final product.

Simulating the oil press
The examined OKB-1 machine is a screw type vegetable oil press which works with a rotating worm shaft.The head of the press machine (Fig. 1) where the extraction proceeds has a complex internal geometry.Basically the worm shaft has three sections with different shaped worms.Their function is to increase It is important to differentiate those sections where the geometry of the seeds would deform significantly, because in this stage of our research, the elements cannot be deformed [3].In this model every seed is simulated as a ball, so we cannot study them at their highly deformed state with the current material model.
The study of the oil press is made stepwise to validate the complete simulation process.Therefore the analysis is made for the section of the worm shaft.It develops from the simple built model into a more complex one.This method is reasoned because it is important to make sure of the realistic properties and effects of the settings.The examination of the realistic behaviour is only possible with evolutionary developments [4].The development of the models to set up new functions and inquiry program parts means complex and challenging programming tasks.
In the previous stages of the research we set up models with simple one sectioned screw geometry which presented the first, named inaugurator section of the press machine.It contained the surface models of the casing, the closing covers at the two ends, the shaft and the helicoid.There was a gap between the end of the worm and the closing cover.The pressure reached its highest value after the screw.The gap contributed to consolidate the stress of the particles at that high pressure avoiding false values.At the first stage of the research the worm was created using a single helicoid.This means the thickness of thread was zero.Developing the simple geometry we have created a new model.

The model
The initial model was improved by giving thickness to the thread (Fig. 2).This model contains the casing, the closing covers, the shaft and the helicoid as the previous one but with an extra helicoid which represents the thickness of the thread.In addition, there is a small surface to close the end of the worm.
The used geometry is an approach of the real section with simplifications applied to reduce the complexity of the surface model.The model has a square thread which means the thread angle is zero.The outer diameter of the worm is equal to the diameter of the casing.
There are different types of seeds which can be utilised in the press machine for oil extraction.In the previous models we had to define the size of the balls which simulate the granular material.The various sizes of the seeds are not represented in this stage of our research, so the choice is based on the dimensions of the geometrical model.Moreover we would like to minimise the number of particles which fill up the geometry to significantly decrease the duration of the simulation but still being accurate.Finally the size was defined according to ensure at least two balls positioned into radiant direction without overlapping each other [5].The applied material and kinematical properties are the same of the previous models (Tab.1).
To reduce computing time, the model was filled up with elements using circle and helix driven generating method.In the gap section between the end of the screw and the closing cover the elements were positioned by their central point one by one next to each other along the circles without overlapping each other or the surfaces.The planes of these circles were perpendicular to the axis of the shaft and the distance between them was equal to the size of the balls.There were two diameters for placing the balls because of the size choosing concept.At the section of the screw the balls were generated along helixes Global friction coefficient 0.1 to fill up the worm like the previous circle driven method.In the worm section with the helix driven method it was possible to avoid generating the balls into the spatial extent of the thread but getting a relatively dense positioned quantity of elements.

Test simulations
Because of the relatively big size and low stiffness of the particles, test runs were needed to ascertain there are no special relations between the geometrical dimensions including the diameter of the balls.In some cases the particles can be positioned in very exceptional way during the pressing process.The undeformable balls can stretch extremely instead of bouncing back.This phenomenon can cause false values.It was important to confirm that the chosen thickness of the thread is proper.
During the pressing process the slot at the end of the head generates reacting pressure against the stream of the material [5].In the test simulations the pressure was replaced by a static closing cover at the end.The exact aim of these simulations was to calculate the maximum pressure on the closing plate at different thread thicknesses including the original dimension.According to the calculated data it is possible to distinguish the incorrect setups.

Main simulation
The parameters of the main simulation are exactly the same as the values of the test simulations .But the function of the closing cover is changed.This model was made to be able to simulate the continuous material stream.The closing plate effects dynamic resistance on the elements by changing its velocity to keep the constant pressure.According to its operation method it is called servo wall.So it allows the particles to have volumetric flow.For this simulation the length of the casing is lengthened in the direction of movement of the closing cover.It prevents the model from opening while the servo wall is moving.With the simulated material stream it is possible to approach the real phenomenon in the press machine.
The aim of the main simulation is to determine the characteristic of distribution of pressure and heat development in the model.With the data we can approach the real phenomenon.In the simulation there are four calculation methods defined to get these values.First we would like to determine the longitudinal distribution of the pressure.Analysing the characteristic of the pressure graph computed on the surface of the casing we can approximately identify the values of mechanical load effecting on the material.The temperature is in the central of the whole analysis because of the quality of the extracted oil.At this stage of the research we focus on the heat development caused by friction.The purpose of the second calculation method is to determine the longitudinal distribution of heat evolution.According to the diagram of heat evolution and the volumetric flow it would be possible to approach the maximum temperature of the seeds.But because of the simple geometry we concentrate on the determination of the characteristics.The calculation of the exact temperature values is planned in the most developed models.In addition, for the analysis of pressure and heat we set up two extra computing methods to specify these values along the thread of the worm.It is reasoned because the worm transmits the kinetic energy to the material and effects pressure on it.So it is presumed that the highest values appear on the surface of the thread.

Calculation methods of pressure and heat development 2.6.1 Pressure
The value of pressure which effects on a surface is calculated from the sum of the normal forces at the contacts between the particles and the walls.During the simulation the program checks all of the contacts in the model in every time step and selects those which contact the surface.Then it is able to determine selected section by summarising the normal forces then dividing with the size of the area.The pressure was computed on the surface of sections.Because of the fluctuating behaviour of these values the final data were determined by averaging for 20000 time steps.

Heat development
The value of the generated surplus heat occurred from the contact frictional energy in the simulation was computed automatically from the contact force, the displacement and the frictional coefficient.To identify the longitudinal distribution of heat development the volume of the model was virtually divided into equal segments along the axis of the shaft.In every simulation cycle or rather in every unit of time we summarise the evolved heat at all contacts in each volume segment.It provides the thermal performance of friction.Similarly to the pressure, the values are fluctuating so the averaging method is reasoned using 20000 time steps.
In case of calculating the heat evolution on the front surface of the thread we analysed the particles which are contacting the front helicoid.To interpret the thermal performance on a surface the values have to be divided with the size of the area.It gives the dimension of heat flux density which is Wm −2 .

Test simulations
With the analysis of the results of the test simulations it is possible to distinguish the incorrect thread thicknesses.There can be dimensions which contribute to false outcome caused by special geometrical relations.The thread of the worm is built up by two helicoids.To specify its parameters with a simple variable, the thickness is defined by an angle which gives the degree of rotation of the second helicoid from the first one around the axis.This theory is illustrated in Fig. 4 at the cross-section.There were six different angles set up for the test runs.In Fig. 3. it is obvious that the curves of 0 • and 45 • are salient.At these angles the pressure does not consolidate in this time range about a pressure level like the other graphs.These false values were caused by the specially positioned particles which can extremely stretch during the pressing process.This phenomenon causes the higher pressure data at those angles.The curve of 22,5 • proves this idea because it is not located between the two salient curves but at the other ones.
Examining the correct curves, the averaged pressure values show the expected relation (Fig. 4).The averaging was done in the consolidated state which started at t 0 = 4 [s] and ended at   .The scale of the simulated angles is tight.The size of the particles delimits the angle values using in the simulation because at higher degrees the small amount of material in the model could not effect significant pressure.It would lead to inaccurate data.Moreover changing the size of the particles would highly increase the duration of the simulations.

Main simulation
In the main simulation the closing plate at the end of the model was functioned as a servo wall which effects constant resistance on the particles.The pressure on this closing cover was set to 4 kPa which is lower than its reachable maximum value according to the test runs.It means that this setup enables the volumetric flow of the material.This servo wall keeps the given pressure value with 5% inaccuracy.Comparing the length of the model to the velocity of closing cover we can note that the ve-locity has a very small value.But it has a significant effect on the computed data.
The pressure was calculated on the surface of the casing.It was divided into equal sections along the axis of the shaft.In Fig. 5. the part from 200 mm to 300 mm belongs to the last two thread and the section from 300 mm to 350 mm belong to the gap.Analysing the pressure diagram of its longitudinal distribution we can see that the significant increase starts just at the end of the worm, especially at its last thread (from 250 mm to 300 mm).The increase continues to the end of the last thread.Then it decreases a little under the local maximum.The last value seems very deflecting comparing it with the 4 kPa pressure of the servo wall.At the last thread the bordered pressure increment is caused by the property of the discontinuous material.At fluids the longitudinal distribution is more uniform because their press transmitting attribute is better.Moreover the fluids can induce flow at the slot between the worm and the casing.According to the longitudinal distribution of the pressure we can identify the most loaded part of the model.Probably the biggest quantity of heat is evolved in this section.
Examining the diagram of heat development (Fig. 6) we can see that the increase of the heat evolution starts at the last thread likewise the pressure.It reaches its maximum directly after the last thread but in the next sections the values are decreasing.It is observable that the most intense increase of thermal performance caused by friction is in the high pressure locations.The rotating helicoid affects the majority of the pressure on these faces and arises most of the heat.Calculating the longitudinal distribution of pressure gives inaccurate data because the model is built up by helicoids.It means that one section of the casing includes surfaces with different pressure conditions caused by the helical geometry.In actual fact the front helicoid affects the high pressure on the material.But each pressure calculating sections of the casing averaging both side of a thread.This method makes inaccurate to define the smooth increase of pressure.So the correct way is computing the pressure values on the front surface of the thread.Every 360 o rotation of the thread of the front helicoid was di- vided into six equal sections.
In Fig. 7. the angle values present the two last rotations of the thread of the front helicoid.The numbering starts from the end point with zero.Analysing the graph the increase of pressure looks more accurate and smooth than in the diagram of longitudinal distribution (Fig. 5).The significant change of the values is perceivable at the last five points.Fig. 8 Distribution of averaged heat development along the front helicoid The graph of heat evolution (Fig. 8) shows the same characteristic as the pressure graph (Fig. 7) with the equivalent numbering of angle values.The significant increase can be seen at the last four computing sections of the last thread rotation.These points can be approached by a linear function likewise at the pressure graph.

Conclusion
Using the discrete element method for the analysis of an oil press seems promising.The behaviour of the granular material was simulated successfully to calculate specified parameters.
As it was experienced there are special geometrical setups which could lead to false results.It means that in some cases the undeformable elements can be positioned in an exceptional way.That phenomenon leads to extreme stretch of the particles and it gives incorrect outcomes.The test runs proved that the original thread thickness defined by 90 • is appropriate for the main simulation.Because the averaged pressure on the closing plate belonging to this angle presented the expected value.
In the main simulation, pressure and heat development values were calculated by different methods.The characteristics of the longitudinal distribution of pressure graph shows that the largest part of the mechanical load effects in the segment of the last thread and in the gap.Moreover the diagram of longitudinal distribution of heat development presents similar outcome in heat load.What means that the analyses have to focus on the last segment of the oil press.The high pressure in the end of the model causes high thermal performance.It concentrates at the surface of the last thread rotation so the seeds in contact with the worm get the highest heat load.Probably the oil content of these seeds could reach the highest temperature.

Fig. 1 .
Fig. 1.The head of the oil press.

Fig. 3 .
Fig. 3.The change of pressure as function of time with different thread thicknesses defined by angle.

Fig. 4 .
Fig. 4. Averaged pressure on the closing cover in the consolidated state.

Fig. 5 .
Fig. 5. Longitudinal distribution of averaged pressure on the casing.

Fig. 7 .
Fig. 7. Distribution of averaged pressure development along the front helicoid.

Fig. 8 .
Fig. 8. Distribution of averaged pressure development along the front helicoid.