A computational fluid dynamics prediction method for the motion state of particles in pneumatic transport

Abstract

The application discloses a kind of gas power transport in the movement state of particle computational fluid dynamics prediction method, comprising: obtaining continuous phase initial state and particle initial state in gas power transport;With continuous phase initial state and particle initial state as current time step continuous phase initial state and current time step particle initial state, current time step processing is carried out, and in current time step processing, for the continuous phase distribution function corresponding to current time step continuous phase initial state and current time step particle initial state, updated particle state and updated continuous phase distribution function are obtained by processing;According to updated particle state and updated continuous phase distribution function, the movement state of particle in gas power transport is determined.The application can accurately predict the movement state of cylindrical particle in gas power transport by combining lattice Boltzmann method and discrete element method, with high precision and low cost, which is conducive to popularization and application in food processing and pharmaceutical industry.

Description

Technical Field

[0001] This invention relates to the field of predicting the motion state of particles in fluids, and more particularly to a computational fluid dynamics method for predicting the motion state of cylindrical particles in pneumatic transport. Background Technology

[0002] The movement of solid particles along gas flow in pipelines is a common phenomenon in many industrial processes. This process is known as pneumatic conveying. Pneumatic conveying has wide applications in chemical processes, mining and metallurgy, food processing, and the pharmaceutical industry. It typically uses gas as the transport medium to transport solid particulate raw materials through pipelines to processing or storage equipment. This process enables rapid and efficient material transfer, improves production efficiency and quality control, and offers advantages such as system design flexibility and relatively low pollution risk.

[0003] Due to the interaction between discrete-phase solid particles and continuous-phase gas, pneumatic transport exhibits strong non-uniformity and multi-temporal-scale characteristics. Therefore, a deeper understanding of the interaction mechanism between particles and turbulent gas flow, and better development and utilization of the application potential of this complex flow to achieve the goals of reducing transport resistance and system energy consumption, has always been one of the important research topics in the field of fluid mechanics. Conventional experimental methods obtain the flow resistance and spatial distribution of particles in pneumatic transport systems through measurement and sampling analysis, thereby obtaining the system's energy efficiency characteristics. However, this method requires a large investment of human and material resources, and current testing techniques have high requirements for the flow environment, thus posing certain constraints on the measurement of gas-solid two-phase flow fields in practical applications.

[0004] Computational fluid dynamics (CFD), as a numerical method, physically models real-world flow systems and uses computers for numerical simulation and analysis. This method is not limited by experimental techniques and can efficiently and comprehensively obtain temporal and spatial information of gas-solid two-phase flow fields, reflecting the flow characteristics of pneumatic transport systems. However, current research in CFD for gas-solid two-phase transport mainly focuses on relatively easy-to-model spherical particles, with limited understanding of cylindrical particles commonly found in food processing and pharmaceutical industries during turbulent pneumatic transport in channels. Accurately and efficiently obtaining the influence of the anisotropic topological shape of cylindrical particles on their state distribution, including spatial orientation and aggregation preferences, and their impact on system flow resistance, presents new and significant application-oriented requirements for the development of existing CFD methods.

[0005] Therefore, existing technologies still need improvement and development. Summary of the Invention

[0006] In view of the shortcomings of the prior art, the purpose of this invention is to provide a computational fluid dynamics prediction method for the motion state of particles in pneumatic transport, so as to solve the problem that the prior art cannot predict the motion state of cylindrical particles commonly used in actual food processing and pharmaceutical industries in turbulent pneumatic transport in channels.

[0007] The technical solution of the present invention is as follows:

[0008] A computational fluid dynamics method for predicting the motion state of particles in pneumatic transport includes:

[0009] To obtain the initial state of the continuous phase and the initial state of particles in pneumatic transport;

[0010] Using the initial state of the continuous phase and the initial state of the particles as the initial state of the continuous phase and the initial state of the particles at the current time step, the current time step processing is performed, which includes:

[0011] Based on the lattice Boltzmann equation, the continuous phase distribution function is obtained using the current time step continuous phase initial state;

[0012] Calculate the interaction forces between the continuous phase and the particles, as well as the lubricating and contact forces between particles and between particles and solid walls;

[0013] Based on the interaction force, the lubrication force, and the contact force, the initial state of the particles at the current time step is updated to obtain the updated particle state.

[0014] Based on the updated particle state, the surface continuous phase distribution function of the particle is reconstructed to obtain the updated continuous phase distribution function, and the current time step processing ends.

[0015] The motion state of the particles in the pneumatic transport is determined based on the updated particle state and the updated continuous phase distribution function.

[0016] In one implementation, after ending the current time step processing, the process includes:

[0017] The updated continuous phase distribution function is used to replace the continuous phase initial state corresponding to the current time step, and the updated particle state is used to replace the particle initial state at the current time step.

[0018] Repeat the current time step processing;

[0019] After reaching the predetermined time step, the final particle state and the final continuous phase distribution function are output;

[0020] The motion state of the particles in the pneumatic transport is determined based on the final particle state and the final continuous phase distribution function.

[0021] In one embodiment, the initial state of the continuous phase includes: continuous phase flow field density, continuous phase initial velocity, and continuous phase external force term;

[0022] The initial state of the particles includes: particle size, particle volume fraction, initial particle position, particle spatial orientation angle, and initial particle velocity.

[0023] In one implementation, the continuous phase distribution function is obtained using the current time step's initial continuous phase state according to the lattice Boltzmann equation, including:

[0024] The lattice Boltzmann equations, which incorporate the initial state of the continuous phase at the current time step, are transformed into the moment space, and different moments evolve according to different relaxation times.

[0025] The evolved moments are inversely transformed to obtain the continuous phase distribution function.

[0026] In one embodiment, calculating the interaction force between the continuous phase and the particles includes:

[0027] Obtain the projection of the continuous phase distribution function onto the moment space;

[0028] By combining the projection and the continuous phase distribution function, a collision process and a migration process are performed to obtain the evolved continuous phase distribution function;

[0029] Based on the evolved continuous phase distribution function and the no-slip boundary condition, the interaction force between the continuous phase and the particles is obtained.

[0030] In one implementation, an interpolation collision bounce scheme is used to determine the no-slip boundary conditions of the continuous phase on the particles and on the solid wall.

[0031] In one implementation, updating the initial state of the particles at the current time step based on the interaction force, the lubrication force, and the contact force to obtain the updated particle state includes:

[0032] Execute multiple sub-loop steps, each of which generates a sub-loop update particle state, and uses the sub-loop update particle state as the initial state of the particle to enter the next sub-loop step.

[0033] After completing all the sub-loop steps, the updated particle state is obtained.

[0034] In one implementation, the sub-loop step includes:

[0035] Based on the initial state of the particles at the current time step, determine the first positional relationship between the particles and the second positional relationship between the particles and the solid wall.

[0036] Based on the first positional relationship and the second positional relationship, the lubrication force correction and contact force correction of the particles are determined;

[0037] Based on the lubrication force correction and the contact force correction, combined with the interaction force, the updated particle state is obtained.

[0038] In one implementation, the surface continuous phase distribution function of the particles is reconstructed based on the updated particle state to obtain the updated continuous phase distribution function, including:

[0039] Based on the updated particle state, the distance between the fluid grid point of the continuous phase and the surface of the particle, and the normal direction of the surface of the particle are determined near the surface of the particle.

[0040] Based on the distance and the normal direction, the surface continuous phase distribution function of the particle is reconstructed, and the updated continuous phase distribution function is obtained by combining the surface continuous phase distribution function and the continuous phase distribution function.

[0041] Based on the updated continuous phase distribution function, the density fluctuation value, continuous phase velocity, and continuous phase pressure of the continuous phase are calculated.

[0042] In one embodiment, the motion state of the particles includes: particle flow resistance, particle position, and particle orientation angle.

[0043] In summary, this invention discloses a computational fluid dynamics prediction method for the motion state of particles in pneumatic transport, comprising: acquiring the initial state of the continuous phase and the initial state of the particles in pneumatic transport; using the initial state of the continuous phase and the initial state of the particles as the initial state of the continuous phase and the initial state of the particles in the current time step, and performing current time step processing; determining the motion state of the particles in pneumatic transport based on the updated particle state and the updated continuous phase distribution function. The current time step processing includes: obtaining the continuous phase distribution function using the initial state of the continuous phase in the current time step according to the lattice Boltzmann equation; calculating the interaction force between the continuous phase and the particles, as well as the lubrication force and contact force between particles and between particles and solid walls; updating the initial state of the particles in the current time step according to the interaction force, the lubrication force, and the contact force, to obtain the updated particle state; reconstructing the surface continuous phase distribution function of the particles according to the updated particle state, to obtain the updated continuous phase distribution function, and ending the current time step processing. This invention combines the lattice Boltzmann method with the discrete element method to accurately predict the motion state of cylindrical particles in pneumatic transport. It is highly accurate and low in cost, which is conducive to its widespread application in the food processing and pharmaceutical industries. Attached Figure Description

[0044] To more clearly illustrate the technical solutions in the embodiments of this application, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are only some embodiments of this application. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0045] Figure 1 This is a flowchart illustrating the steps of the computational fluid dynamics prediction method for the motion state of particles in pneumatic transport as described in this invention.

[0046] Figure 2 This is a flowchart illustrating one embodiment of the computational fluid dynamics prediction method for particle motion state in pneumatic transport according to the present invention.

[0047] Figure 3 This is a schematic diagram of the pneumatic transport channel flow described in this invention.

[0048] Figure 4 This is a schematic diagram of the pneumatic transport channel flow of disc-shaped cylindrical particles with an aspect ratio of 1 / 3 in Embodiment 1 of the present invention.

[0049] Figure 5 This is a schematic diagram of the pneumatic transport channel flow of a perfectly cylindrical particle with a length-to-diameter ratio of 1 in Embodiment 2 of the present invention.

[0050] Figure 6This is a schematic diagram of the pneumatic transport channel flow of strip-shaped cylindrical particles with an aspect ratio of 3 in Embodiment 3 of the present invention.

[0051] Figure 7 This is a schematic diagram comparing the effects of different aspect ratios of cylindrical particles on the resistance of the pneumatic transport channel flow system in an embodiment of the present invention.

[0052] Figure 8 This is a schematic diagram comparing the effects of different aspect ratios of cylindrical particles on the aggregation position of flowing particles in a pneumatic transport channel, as described in an embodiment of the present invention.

[0053] Figure 9 This is a schematic diagram comparing the effects of different aspect ratios of cylindrical particles on the spatial orientation of particles flowing in a pneumatic transport channel, as described in an embodiment of the present invention.

[0054] Figure 10 This is a flowchart of the current time step processing steps in the computational fluid dynamics prediction method for particle motion state in pneumatic transport as described in this invention. Detailed Implementation

[0055] This application provides a computational fluid dynamics prediction method for the motion state of particles in pneumatic transport. To make the objectives, technical solutions, and effects of this application clearer and more explicit, the following detailed description is provided with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only for explaining this application and are not intended to limit this application.

[0056] Those skilled in the art will understand that, unless specifically stated otherwise, the singular forms “a,” “an,” “the,” and “the” used herein may also include the plural forms. It should be further understood that the term “comprising” as used in this application means the presence of the stated features, integers, steps, operations, elements, and / or components, but does not exclude the presence or addition of one or more other features, integers, steps, operations, elements, components, and / or combinations thereof. It should be understood that when we say an element is “connected” or “coupled” to another element, it can be directly connected or coupled to the other element, or there may be intermediate elements. Furthermore, “connected” or “coupled” as used herein can include wireless connections or wireless coupling. The term “and / or” as used herein includes all or any units and all combinations of one or more associated listed items.

[0057] It will be understood by those skilled in the art that, unless otherwise defined, all terms used herein (including technical and scientific terms) have the same meaning as commonly understood by one of ordinary skill in the art to which this application pertains. It should also be understood that terms such as those defined in general dictionaries should be understood to have the same meaning as in the context of the prior art, and should not be interpreted in an idealized or overly formal sense unless specifically defined as herein.

[0058] It should be understood that the sequence number and size of each step in this embodiment do not imply the order of execution. The execution order of each process is determined by its function and internal logic, and should not constitute any limitation on the implementation process of this application embodiment.

[0059] Computational fluid dynamics (CFD) can simulate the motion of particles during gas transport, thus avoiding the limitations of experimental conditions and saving human and material resources. However, existing CFD methods for gas-solid two-phase transport mainly focus on spherical particles, and the understanding of cylindrical particles in turbulent pneumatic transport in channels, which are actually encountered in food processing and pharmaceutical industries, is still very limited. This invention provides a CFD prediction method for the motion of particles in pneumatic transport based on the lattice Boltzmann-discrete element method. This method can predict the flow resistance and particle state distribution of particles in turbulent pneumatic transport in channels, especially for cylindrical particles. It can be widely used to obtain the flow resistance and particle state distribution in the pneumatic transport of cylindrical particles commonly seen in food processing and pharmaceutical industries. This method obtains the resistance and particle state distribution of turbulent pneumatic transport in channels through numerical simulation analysis of computational fluid dynamics, without the need for complex sampling analysis of actual pneumatic transport channel flow systems, avoiding the limitations of experimental conditions and thus saving a lot of human and material resources. Meanwhile, unlike conventional numerical research methods for gas-solid two-phase flow systems with spherical particles, this method accurately considers the influence of the anisotropic topological shape of cylindrical particles on flow resistance and particle state distribution. It is a new method to obtain the flow resistance and particle state distribution of turbulent pneumatic transport in channels containing cylindrical particles.

[0060] Specifically, such as Figure 1 As shown, the computational fluid dynamics prediction method for the motion state of particles in pneumatic transport according to the present invention includes the following steps:

[0061] S100: Obtain the initial state of the continuous phase and the initial state of the particles in pneumatic transport.

[0062] Furthermore, such as Figure 2As shown, the initial state of the continuous phase includes: continuous phase flow field density, continuous phase initial velocity, and continuous phase external force term; the initial state of the particles includes: particle size, particle volume fraction, particle initial position, particle spatial orientation angle, and particle initial velocity. Optionally, during the simulation prediction process, the continuous phase flow field density, the continuous phase initial velocity, and the continuous phase external force term are obtained by defining the continuous phase gas density, initial velocity, and external force term, thereby determining the initial state of the continuous phase; the particle size, particle volume fraction, particle initial position, particle spatial orientation angle, and particle initial velocity are obtained by defining the size and volume fraction of cylindrical particles, and randomly distributing the initial position and spatial orientation angle of cylindrical particles in the flow field, thereby determining the initial state of the particles. The initial velocity of the particles is set as the flow field velocity at their center of mass, and the flow inlet / outlet conditions, solid wall boundary conditions for pneumatic transport, and calculation time are simultaneously set to begin the simulation prediction process.

[0063] First, a model of a solid cylindrical particle in the discrete phase needs to be established. Unlike isotropic spherical particles, cylindrical particles have two planes with infinite radii of curvature and two rough edges. The anisotropic topological shape of the cylindrical particle introduces greater complexity to the particle modeling. To accurately represent the spatial attitude of the cylindrical particle, two coordinate systems are used, including:

[0064] Inertial reference frame<x,y,z> The inertial reference frame is a fixed spatial coordinate system used to describe the spatial translation of the continuous phase gas flow field and particles.

[0065] Satellite reference frame<x’,y’,z’> The aforementioned body reference frame is a moving coordinate system that follows the movement of the particles, with the center of mass of the cylinder as the origin and the axis of rotational symmetry as the z' axis. The rotational orientation of the cylindrical particles and the interactions between the particles are solved in this reference frame.

[0066] The inertial reference frame and the body reference frame can be transformed into each other through a rotation matrix R. The transformation equation is as follows:

[0067] R·<x′,y′,z′> =<x,y,z> R -1 ·<x,y,z> =<x′,y′,z′> .

[0068] The rotation matrix R is obtained from the quaternions [q1, q2, q3, q4]. The quaternions satisfy the normalization condition. And rotated by Euler angles φ, θ of the rigid body. express:

[0069]

[0070] The rotation matrix R is expressed as:

[0071]

[0072] Before applying interpolated collision bounce boundary conditions to cylindrical particles, it is necessary to calculate the distance q between the boundary fluid grid points and the solid wall. When constructing the distribution function for newly added fluid grid points, the local normal direction n of the solid wall also needs to be known. However, given the complex shape of cylindrical particles, calculating q and n requires first using the transformation equation between the inertial reference frame and the body reference frame to transform the fluid grid point coordinates from the inertial reference frame to the particle's body reference frame, and then calculating based on the relative positional relationship between the fluid points and the solid particles in this coordinate system.

[0073] Since cylindrical particles are rotationally symmetric, their spatial orientation can be represented by the angles between the direction vector of the axis of symmetry and the three coordinate axes of the inertial coordinate system. First, the projection of the rotational symmetry axis z' onto the inertial reference frame is obtained through coordinate transformation.<x0,y0,z0> :

[0074]

[0075] Then, the axis of symmetry and the coordinate axes are obtained through algebraic operations.<x,y,z> The included angle, i.e., the particle orientation angle <α, β, γ>:

[0076]

[0077] By establishing a model of the cylindrical particle in the above manner, the initial state of the particle can be obtained and used for subsequent prediction of the particle's motion state.

[0078] Furthermore, such as Figure 1 As shown, after step S100, the following step is also included:

[0079] S200. Using the initial state of the continuous phase and the initial state of the particles as the initial state of the continuous phase and the initial state of the particles in the current time step, perform the current time step processing.

[0080] In the current time step processing, the updated particle state and the updated continuous phase distribution function are obtained by processing the continuous phase distribution function corresponding to the initial state of the continuous phase at the current time step and the initial state of the particles at the current time step.

[0081] Specifically, such as Figure 10 As shown, the current time step processing specifically includes the following steps:

[0082] H100. Based on the lattice Boltzmann equation, the continuous phase distribution function is obtained using the initial state of the continuous phase at the current time step.

[0083] Specifically, a continuous-phase gas flow model is first established, and based on the initial state of the continuous phase in the current time step, the lattice Boltzmann equations incorporating the initial state of the continuous phase in the current time step are obtained. Specifically, based on the lattice Boltzmann method, the Boltzmann equations are discretized and projected onto the velocity directions of each particle to obtain the lattice Boltzmann equations governing the evolution of the continuous-phase gas:

[0084]

[0085] Where f i This represents the distribution function of mesoscopic particles, which is the sum of time t, spatial position x, and particle velocity e. i The function. e i The specific directions of the D3Q27 velocity set, which are discretized in two dimensions along twenty-seven velocity directions, are as follows:

[0086]

[0087] Where c = δ x / δ t δ represents the velocity per unit grid cell. x and δ t These represent the spatial step and the time step, respectively.

[0088] τ is the dimensionless relaxation time, which is related to the kinematic viscosity coefficient v of the gas:

[0089]

[0090] Where c S For the speed of sound in a lattice. f i eq The equilibrium distribution of particles is represented in incompressible flow as:

[0091]

[0092] Where ρ is the density of the continuous phase gas, u is the macroscopic velocity of the gas, I is the identity matrix, and ω i The weighting coefficients for each direction corresponding to the velocity set of D3Q27 particles are as follows:

[0093]

[0094] F i The effect of external forces on the particle distribution function is expressed as:

[0095]

[0096] Where F represents the macroscopic external force term. Since a single relaxation time τ is used in the lattice Boltzmann equations governing the evolution of continuous-phase gases, it means that the distribution function f... iThe equilibrium state approaches at the same relaxation rate in all twenty-seven velocity directions, which is inconsistent with the complete Boltzmann equation. Therefore, a multi-relaxation-time method, which better matches the predictions of the Boltzmann equation, is used to adjust the distribution function f of the lattice Boltzmann equation. i First, the process is converted to moment space. Different moments evolve according to different relaxation times. Then, the evolved moments are inversely transformed to obtain the continuous phase distribution function.

[0097] The continuous phase distribution function is written as:

[0098]

[0099] Where M is the transformation matrix, M -1 Let m be the inverse transformation matrix. i =Mf i The distribution function f of the lattice Boltzmann equation is... i The transformed moments are S = diag(s0, s1, ..., s). 26 ) is the relaxation coefficient matrix, where the relaxation coefficients related to the kinematic viscosity of the gas are as follows:

[0100]

[0101] Furthermore, such as Figure 10 As shown, in the current time step processing, after step H100, the following step is also included:

[0102] H200. Based on the continuous phase distribution function, calculate the interaction force between the continuous phase and the particles, as well as the lubricating force and contact force between particles and between particles and the solid wall.

[0103] In one implementation, such as Figure 2 As shown, after performing the "collision" and "migration" process of the continuous phase distribution function based on the lattice Boltzmann method, the continuous phase flow field information is updated, then the interpolation rebound boundary condition is performed on the particle surface, and finally the fluid force on the particle surface, i.e. the interaction force between the continuous phase and the particles, is calculated.

[0104] Specifically, the projection of the continuous phase distribution function onto the moment space is obtained, where let for Projection in moment space:

[0105]

[0106] Where u, v, and w are the three macroscopic velocity components of the gas, ρ0 is the steady part of the gas density, and δρ is the fluctuating value of the gas density. The sum of these two items constitutes the gas density ρ = ρ0 + δρ.

[0107] Then, by combining the projection with the continuous phase distribution function, a collision process and a migration process are performed to obtain the evolved continuous phase distribution function.

[0108] The "collision" process of the continuous phase distribution function at the local lattice point can be expressed as:

[0109]

[0110] Where f i * Let be the distribution function after the collision.

[0111] The "migration" process of the continuous phase distribution function to adjacent grid points can be represented as:

[0112] f i (x+e i δ t ,t+δ t )=f i * (x, t).

[0113] To accurately analyze the influence of channel walls and solid particle surfaces on fluids in pneumatic transport, the lattice Boltzmann method uses an interpolated collision bounce scheme to handle no-slip boundary conditions on solid surfaces (solid walls or particle surfaces). Specifically, the interpolated collision bounce scheme determines the no-slip boundary conditions of the continuous phase on the particle surface and the solid wall surface, thereby achieving the same second-order accuracy at solid phase (including particle and solid wall) boundaries as at internal fluid lattice points. In the interpolated collision bounce scheme, q represents the distance between the boundary fluid lattice point and the solid surface.

[0114] q=|x f -x w | / |x f -x s |,

[0115] Where x f x w and x s These represent the positions of the fluid lattice point closest to the solid surface, the solid surface itself, and the solid lattice point closest to the solid surface, respectively. When there are at least two layers of fluid lattice points outside the solid surface, different two-point interpolation methods are used to calculate the unknown distribution function f, depending on the q value. i (x f ,t+δ t When q≤0.5:

[0116]

[0117] The superscript * indicates the variable state after the collision, and the subscript... Indicates the direction of particle velocity from the fluid lattice point to the solid surface, x ff u represents the second fluid lattice point closest to the solid surface. w This indicates the current velocity of the solid surface.

[0118] When q > 0.5:

[0119]

[0120] When two solid surfaces (between particles or between a particle and a solid wall) are brought closer together to a distance of less than two grid points, x ff The two-point interpolation method no longer exists. In this case, only x is used. f Pre-collision distribution function f at lattice points i (x f (t) and the post-collision distribution function f i * (x f Construct boundary conditions with second-order accuracy for a single point (t):

[0121]

[0122] Based on the evolved continuous phase distribution function, and after processing the no-slip boundary conditions on the solid wall and particle surfaces using the interpolation collision bounce scheme, the force F exerted by the continuous phase gas on the solid particles can be calculated using the momentum exchange method. H and torque T H Thus, the interaction force between the continuous phase and the particles is obtained:

[0123]

[0124]

[0125] Among them (e) i -u w )and The relative velocities of particles on the particle surface are characterized to partially correct for Galilean invariance errors caused by the motion of solid particles. Simultaneously, the continuous transformation between fluid and solid grid points on the computational grid due to particle motion further offsets Galilean invariance errors, enabling more accurate calculation of the forces and torques exerted by the continuous phase gas on the solid particles.

[0126] Furthermore, such as Figure 10 As shown, in the current time step processing, after step H200, the following step is also included:

[0127] H300. Based on the interaction force, the lubrication force, and the contact force, update the initial state of the particles at the current time step to obtain the updated particle state.

[0128] In one implementation, step H300 consists of multiple sub-loop steps, each with identical content. Each execution of a sub-loop step generates a sub-loop updated particle state, which serves as the initial particle state for the next sub-loop step. After completing all sub-loop steps, the updated particle state is obtained. The time step of each sub-loop is much smaller than the current time step. This is because collisions between particles and between particles and solid walls are very frequent in pneumatic transport, and these collisions often occur within a very short time, causing changes in the particle's motion state. If simulation and prediction are performed using the time step of continuous phase motion, situations such as penetration or overlap, which do not conform to actual physical laws, will inevitably occur, leading to an inaccurate prediction result. By dividing the process into multiple sub-loop steps, the particle state can be simulated and calculated on a shorter time scale, ensuring a meaningful analytical solution and a calculation result that more closely matches reality.

[0129] Specifically, the sub-loop steps include:

[0130] Based on the initial state of the particles at the current time step, determine the first positional relationship between the particles and the second positional relationship between the particles and the solid wall.

[0131] Based on the first positional relationship and the second positional relationship, determine the lubrication force correction and contact force correction of the particles; and

[0132] Based on the lubrication force correction and the contact force correction, combined with the interaction force, the updated particle state is obtained.

[0133] In one implementation, such as Figure 2 As shown, firstly, the particle position and orientation are detected, and lubrication force and contact force are corrected; then, particle movement is performed, and the particle position, speed, and orientation are updated to obtain the sub-loop updated particle state. Finally, it is determined whether all sub-loop steps have been completed. If the sub-loop steps have not been completed, the obtained sub-loop updated particle state is used as the initial particle state for the next sub-loop step, and the sub-loop steps are repeated. If the sub-loop steps have been completed, the sub-loop updated particle state is output as the updated particle state.

[0134] Specifically, when solid particles in a fluid approach each other to a small distance but before contact, the fluid between particles or between particles and the solid wall thins, resulting in significant pressure exerted by the fluid on the solid surface. Because the distance between solid surfaces is small, there are often no fluid grid points between them, making it impossible to use the momentum exchange method to analyze the force exerted by the fluid on the solid surface. In this case, the influence of lubrication force on the particle state needs to be considered. If the minimum distance between two particle surfaces or between a particle surface and the wall is less than a certain threshold, and there is a relative velocity between them, contact is considered to have occurred, and a contact force is applied. In this case, the influence of the contact force on the particle state needs to be considered.

[0135] Because the characteristic timescale of the lubrication and contact force models is much smaller than that of the continuous phase gas motion, the interaction forces between the continuous phase and the particles remain unchanged in each sub-cycle step; only the contact and lubrication forces acting on the particles are updated. By dividing the simulation process into multiple sub-cycle steps, the lubrication and contact forces acting on particles or particles and walls during the "approach-contact-disengagement-movement" process can be accurately calculated, determining whether there is contact between particles or between particles and solid walls, thus ensuring stable and accurate final results.

[0136] When considering the influence of lubrication force and contact force on particle state, it is necessary to first determine the first positional relationship between the particles and the second positional relationship between the particles and the solid wall based on the initial state of the particles; then, based on the first positional relationship and the second positional relationship, determine the lubrication force correction and contact force correction of the particles; finally, based on the lubrication force correction and the contact force correction, combined with the interaction force, the updated particle state is obtained.

[0137] When the lubrication force model is introduced to correct the force on the particles, the lubrication force F L，ij writing:

[0138]

[0139] Where μ is the dynamic viscosity coefficient of the fluid, a is the equivalent radius of the solid particle, and u ij =u i- u j Let n be the relative velocity of particle i with respect to particle j. ij Let be the normal vector pointing from the centroid of particle i to the centroid of particle j, and ∈ be the ratio of the shortest distance between particles or between a particle and the wall to the particle radius a.

[0140] ∈0 represents the threshold for enabling the lubrication force model, at which point the shortest distance between particles is approximately one grid space step δ. x That is, ∈0=δ x / a. ∈1 and ∈2 represent the ratio of the distances at which the lubricating force reaches its maximum and zero, respectively. When ∈0 ≥ ∈ > ∈1, the lubricating force increases as the distance between the solid surfaces decreases; when ∈1 ≥ ∈ > ∈2, due to the surface roughness of the solid surfaces, the lubricating force will no longer increase when the distance between the solid surfaces is sufficiently small, but will remain constant; when ∈2 ≥ ∈, the effect of the lubricating force is completely replaced by the contact force, and it no longer affects the solid surface. λ is a function of ∈, written as:

[0141] Particle-particle:

[0142] Particle-wall:

[0143] When considering the influence of contact force on particle state, for cylindrical particles, since the distances from the particle's center of mass to various points on the surface are unequal, the relative orientation of the two particles needs to be considered to obtain the minimum distance. The execution process is as follows: First, calculate the distance between the centers of mass of particles i and j. When the distance between the centers of mass is less than twice the maximum value of the particle's diameter and length, transform the position vector of particle j's axis of symmetry and its center of mass relative to particle i's center of mass to the body reference frame of particle i, and determine whether the two particles are in contact; when the particle is close to the wall, determine whether they are in contact by the angle between the particle's axis of symmetry and the wall's normal direction, and the distance between the particle's center of mass and the wall.

[0144] The forces generated by the contact between two solid surfaces (particle-to-particle or particle-to-solid wall) are calculated using a soft-sphere model based on a "spring-damped" model. Assuming that minute overlapping deformations can occur between the particles, the contact force F experienced by particle i during contact with particle j is calculated. C，ij writing:

[0145] F C，ij =[-k n δ n -β n (u ij ·n ij )]n ij ,

[0146] Where δ n Let k be the overlap distance between particles i and j. n and β n These are the spring stiffness coefficient and damping coefficient, respectively:

[0147]

[0148] Where m e =m i m j / (m i +m j For effective mass, mi and m j Let e ​​represent the masses of particles i and j, respectively. d N is the coefficient of restitution. C δ represents the number of contact time steps. tp The time step for analyzing particle motion is defined. Optionally, in the sub-loop, the time step for analyzing particle motion is set to 1 / 10 of the time step for analyzing the continuous phase gas, i.e., δ. tp =0.1δ t .

[0149] When calculating the interaction force between particles and the continuous phase, the force F exerted by the continuous phase gas on the solid particles is obtained using the lattice Boltzmann method. H and torque T H The lubrication force F obtained from the discrete element method L，ij and contact force F C，ij Particle i at t+δ t The updated particle state at each moment (including translational velocity) Location angular velocity and rotation angle It can be obtained from the following equation:

[0150]

[0151]

[0152]

[0153]

[0154] Repeat the above steps until the particle completes all sub-loop steps to obtain the updated particle state.

[0155] Furthermore, such as Figure 10 As shown, in the current time step processing, after step H300, the following step is also included:

[0156] H400. Based on the updated particle state, reconstruct the surface continuous phase distribution function of the particles to obtain the updated continuous phase distribution function.

[0157] In one implementation, such as Figure 2As shown, this involves recalculating the q-values ​​of each grid point on the particle surface and the normal direction of the wall near the grid point, and reconstructing the distribution function on the new fluid grid points. Specifically, first, based on the updated particle state, the distance between the fluid grid point in the continuous phase and the particle surface, as well as the normal direction of the particle surface, are determined near the particle surface. Then, based on the distance and the normal direction, the surface continuous phase distribution function near the particle surface is reconstructed. Finally, combining the surface continuous phase distribution function and the continuous phase distribution function, the updated continuous phase distribution function is obtained.

[0158] It should be noted that the lattice Boltzmann method uses a fixed Eulerian grid to simulate moving solid particles, eliminating the need to readjust the grid for different particle positions and saving significant computational time. However, the solid grid points inside the particle will revert to fluid grid points after the particle's movement. A portion of the distribution function of the newly formed fluid grid points is obtained from the "migration" process of adjacent grid points, while other unknown distribution functions are extrapolated from the distribution functions of the surrounding fluid grid points, resulting in the evolved continuous phase distribution function, written as:

[0159] f i (x n ,t+δ t )=3f i (x+e l δ t ,t+δ t )-3f i (x+2e l δ t ,t+δ t )+f i (x+3e l δ t ,t+δ t ),

[0160] Where x n e represents the location of the newly emerging fluid grid point. l The direction of the particle velocity that makes the smallest angle with the current particle surface normal is obtained (e l ·n) / (||e l The value of ||·||n|| reaches its minimum, where n represents the current normal direction of the particle surface.

[0161] To ensure that the macroscopic velocity at the new fluid lattice points matches the particle surface velocity, the evolved continuous phase distribution function f is obtained. i (x n ,t+δ t Then, transform all distribution functions at the new fluid grid points to moment space:

[0162] m i (x n,t+δ t )=Mf i (x n ,t+δ t ),

[0163] Then, the first-order velocity moments m1, m2, and m3 are constrained to the particle surface velocities:

[0164] m1=ρu w m2=ρv w m3=ρw w ,

[0165] Where u w v w w w The current particle surface velocity u w The three components. Finally, the corrected moments are transformed back into the distribution function:

[0166] f i (x n ,t+δ t ) = M -1 m i (x n ,t+δ t ).

[0167] The evolved continuous phase distribution function f is obtained through moment space correction. i (x n ,t+δ t This not only maintains the consistency of macroscopic speed, but also helps to stabilize the numerical values ​​of the calculation process, thereby ensuring the stability of the final prediction results.

[0168] Furthermore, by solving the updated continuous phase distribution function, the macroscopic physical quantities of the continuous phase, including the density fluctuation value δρ, the continuous phase velocity u, and the continuous phase pressure p, are expressed by the following equations:

[0169]

[0170]

[0171]

[0172] Furthermore, such as Figure 1 As shown, after step S200, the following step is also included:

[0173] S300. Determine the motion state of the particles in the pneumatic transport based on the updated particle state and the updated continuous phase distribution function.

[0174] In one implementation, such as Figure 2As shown, the process first determines whether the specified time step has been completed after the current time step processing ends. If the specified time step has not been completed, the current time step processing is repeated, and the "collision" and "migration" processes of the updated continuous phase distribution function are re-executed based on the lattice Boltzmann method, and subsequent steps are executed sequentially. After confirming that the specified time step has been completed, the motion state of the particles in the pneumatic transport is output.

[0175] Specifically, firstly, the continuous phase distribution function corresponding to the initial state of the continuous phase at the current time step is replaced with the updated continuous phase distribution function, and the initial state of the particles at the current time step is replaced with the updated particle state. That is, the continuous phase distribution function and particle state output from the previous time step are used as the initial continuous phase distribution function and initial particle state for the next time step. The processing steps for the current time step are repeated until the specified time step is completed. After reaching the predetermined time step, the final particle state and the final continuous phase distribution function are output, and the motion state of the particles in the pneumatic transport is determined based on the final particle state and the final continuous phase distribution function. The motion state of the particles includes: particle flow resistance, particle position, and particle orientation angle.

[0176] Specifically, the particle position is a location. Statistical data changing over time; the particle orientation angle is a rotation angle. Statistical data changing over time; the particle flow resistance can be represented by the drag coefficient C. f express:

[0177]

[0178] Where u m u is the overall average velocity of the flow. τ The wall friction speed is written as:

[0179]

[0180] Where H is half the channel height, (-dp / dx) represents the pressure gradient, and ρ is the density of the continuous phase gas. Keeping the pressure gradient constant, the wall friction velocity u... τ Since it is a constant value, the overall average velocity u of the continuous phase gas flow is... m The larger the value, the greater the drag coefficient C. f The smaller the value, the less resistance the flow experiences, and vice versa.

[0181] The following section further illustrates the computational fluid dynamics prediction method for particle motion state in pneumatic transport, as described in this invention, using solid cylindrical particles of different geometric shapes as examples.

[0182] like Figure 3As shown, solid cylindrical particles undergo gas-solid two-phase flow within a typical pneumatic transport channel. This channel is 30 cm high and 60 cm wide. The cylindrical particles enter the gas-solid two-phase flow system from the left end of the channel along with air. The air is at ambient temperature and pressure, with a wall shear velocity of 0.018 m / s and an average velocity of 0.285 m / s. The average volume fraction of the cylindrical particles within the pneumatic transport channel is maintained at 5%, and the average flow rate is 0.577 kg / s. The shear Reynolds number is 180, and the average Reynolds number is 5500, indicating turbulent flow. An Eulerian grid with dimensions of 1 mm × 1 mm × 1 mm is used, resulting in a grid number of 600 × 300 × 300 in the length, width, and height directions of the channel. The inlet and outlet are set as periodic boundary conditions, and the channel walls and solid particle surfaces are set as no-slip boundary conditions; the time step is 0.0002s, and the calculation time is 150s. The following is a comparison of the system flow resistance and particle state distribution obtained by combining three cylindrical particles with different geometric characteristics flowing in a pneumatic transport channel, including the distribution of particle aggregation preference and spatial orientation, to reflect the influence of the anisotropic topological shape of cylindrical particles on flow resistance and particle state distribution.

[0183] Example 1

[0184] In Example 1, cylindrical particles were set with a length of 0.63 cm and a diameter of 1.89 cm along their axis of rotational symmetry, i.e., disc-shaped particles with an aspect ratio of 1 / 3. The transient position and spatial orientation distribution of the disc-shaped particles at 150 s were calculated as follows: Figure 4 As shown.

[0185] like Figure 7 The drag coefficient of the disc-shaped particle gas-solid two-phase flow system shown is 0.0180, which is the middle value among the three embodiments.

[0186] like Figure 8 As shown, the curves corresponding to Example 1 illustrate the volume fraction of the disc-shaped particles at different heights from the wall, i.e., their spatial aggregation distribution. It can be seen that the disc-shaped particles exhibit peak values ​​at 2 cm and 3.8 cm from the wall, with local volume fractions of 5.3% and 5.6%, respectively. As the wall height increases towards the center of the channel, the particle volume fraction gradually decreases to an average volume fraction of around 5%.

[0187] like Figure 9As shown in (a) and (b), the spatial orientation distribution of disc-shaped particles near the wall (0≤y≤3cm) and near the channel center (12≤y≤15cm) is represented by the probability density function of the angle between the particle's axis of symmetry and the three coordinate axes of the inertial reference frame. α represents the angle between the particle's axis of symmetry and the flow direction, β is the angle between the axis of symmetry and the wall normal, and γ is the angle between the axis of symmetry and the spanwise direction. Due to the particle's symmetry, all three angles fall within the range of [0, π / 2]. It can be seen that in the regions near the wall and near the channel center, the angle α of the disc-shaped particles has a relatively high probability of being π / 2, indicating that the axis of symmetry of the disc-shaped particles tends to be perpendicular to the flow direction, i.e., the particles are parallel to the wall, and the side with the smaller cross-sectional area faces the flow direction. This is because the resistance experienced by the disc-shaped particles is less when their axis of symmetry is perpendicular to the flow direction compared to when it is perpendicular to the wall.

[0188] Example 2

[0189] In Example 2, the particles were set to be perfectly cylindrical with a length of 1.3 cm and a diameter of 1.3 cm along their axis of rotational symmetry, i.e., an aspect ratio of 1. The transient position and spatial orientation distribution of the perfectly cylindrical particles at 150 s were calculated as follows: Figure 5 As shown.

[0190] like Figure 7 The drag coefficient of the gas-solid two-phase flow system of the cylindrical particles shown is 0.0189, which is the maximum value among the three embodiments.

[0191] like Figure 8 The curves corresponding to Example 2 illustrate the volume fraction of the cylindrical particles at different heights from the wall, i.e., their spatial aggregation distribution. It can be seen that the cylindrical particles exhibit peak values ​​at 1.5 cm and 3 cm from the wall, with local volume fractions of 5.3% and 5.4%, respectively. As the wall height increases towards the center of the channel, the particle volume fraction gradually decreases to an average volume fraction of around 5%. Simultaneously, in the region closer to the wall (0 ≤ y ≤ 1.5 cm), the cylindrical particles have the largest volume fraction compared to the disc-shaped and strip-shaped particles, which is one of the reasons why the cylindrical particles generate the highest drag coefficient.

[0192] like Figure 9 As shown in (c) and (d), the spatial orientation of the cylindrical particles is observed near the wall (0 ≤ y ≤ 3 cm) and near the channel center (12 ≤ y ≤ 15 cm). It can be seen that near the wall, the angle γ of the cylindrical particles exhibits a higher probability, meaning their axis of rotational symmetry tends to remain parallel to the flow direction near the wall; while in the mainstream region near the channel center, the axis of symmetry prefers to maintain a larger angle α with the flow direction.

[0193] Example 3

[0194] In Example 3, strip-shaped particles with a length of 2.7 cm and a diameter of 0.9 cm along their axis of rotational symmetry (i.e., an aspect ratio of 3) were used. The transient position and spatial orientation distribution of the strip-shaped particles at 150 s were calculated as follows: Figure 6 As shown.

[0195] like Figure 7 The drag coefficient of the strip-shaped particle gas-solid two-phase flow system shown is 0.0177, which is the minimum value among the three embodiments.

[0196] like Figure 8 As shown in the curve corresponding to Example 3, the volume fraction of the strip-shaped particles at different heights from the wall, i.e., their spatial aggregation distribution, is as follows. It can be seen that the strip-shaped particles only exhibit a peak at approximately 5 cm from the wall, with a local volume fraction of 5.5%. This differs from the phenomenon observed in disc-shaped and cylindrical particles, which show two peaks near the wall. As the wall height increases towards the center of the channel, the particle volume fraction gradually decreases to an average volume fraction of around 5%.

[0197] like Figure 9 As shown in (e) and (f), the spatial orientation distribution of strip-shaped particles is observed near the wall (0 ≤ y ≤ 3 cm) and near the channel center (12 ≤ y ≤ 15 cm). It can be seen that near the wall, the included angle β of the strip-shaped particles exhibits a higher probability, while the included angle α exhibits a lower probability. This indicates that the rotational symmetry axis of the strip-shaped particles tends to align with the flow direction near the wall, resulting in the smallest cross-sectional area and the least flow resistance for the particles in the flow direction. Near the channel center, the orientation angle of the strip-shaped particles does not show a significant preference. Figure 9 The orientation angle distribution shown in (e) and (f) is similar to Figure 6 The observed phenomena are consistent.

[0198] As can be seen from the above three embodiments, the computational fluid dynamics prediction method for particle motion state in pneumatic transport described in this invention can accurately model and calculate the gas-solid two-phase flow based on the anisotropic topological shape of cylindrical particles, thereby accurately predicting the flow resistance and particle state distribution of cylindrical particles in turbulent pneumatic transport in channels. This provides a theoretical basis and new technical means for better developing and utilizing the application potential of pneumatic transport systems for cylindrical particles commonly used in food processing and pharmaceutical industries, and achieving the goals of reducing transport resistance and lowering system energy consumption.

[0199] In summary, this invention discloses a computational fluid dynamics prediction method for the motion state of particles in pneumatic transport, comprising: acquiring the initial state of the continuous phase and the initial state of the particles in pneumatic transport; using the initial state of the continuous phase and the initial state of the particles as the initial state of the continuous phase and the initial state of the particles at the current time step, performing current time step processing; in the current time step processing, processing the continuous phase distribution function corresponding to the initial state of the continuous phase at the current time step and the initial state of the particles at the current time step to obtain updated particle states and updated continuous phase distribution functions; and determining the motion state of the particles in pneumatic transport based on the updated particle states and updated continuous phase distribution functions. This invention, by combining the lattice Boltzmann method and the discrete element method, can accurately predict the motion state of cylindrical particles in pneumatic transport, with high accuracy and low cost, which is beneficial for its widespread application in the food processing and pharmaceutical industries.

[0200] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of this application, and are not intended to limit them. Although this application has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of this application.

Claims

1. A computational fluid dynamics method for predicting the motion state of cylindrical particles in pneumatic transport, characterized in that, include: To obtain the initial state of the continuous phase and the initial state of particles in pneumatic transport; Using the initial state of the continuous phase and the initial state of the particles as the initial state of the continuous phase and the initial state of the particles at the current time step, the current time step processing is performed, which includes: Based on the lattice Boltzmann equation, the continuous phase distribution function is obtained using the current time step continuous phase initial state; Calculate the interaction forces between the continuous phase and the particles, as well as the lubricating and contact forces between particles and between particles and solid walls; Based on the interaction force, the lubrication force, and the contact force, the initial state of the particles at the current time step is updated to obtain the updated particle state. Based on the updated particle state, the surface continuous phase distribution function of the particle is reconstructed to obtain the updated continuous phase distribution function, and the current time step processing ends. The motion state of the particles in the pneumatic transport is determined based on the updated particle state and the updated continuous phase distribution function. The motion state of the particles includes: particle flow resistance, particle position, and particle orientation angle.

2. The computational fluid dynamics prediction method for the motion state of cylindrical particles in pneumatic transport according to claim 1, characterized in that, After the current time step processing ends, the following is included: The updated continuous phase distribution function is used to replace the continuous phase initial state corresponding to the current time step, and the updated particle state is used to replace the particle initial state at the current time step. Repeat the steps of the current time step. After reaching the predetermined time step, the final particle state and the final continuous phase distribution function are output; The motion state of the particles in the pneumatic transport is determined based on the final particle state and the final continuous phase distribution function.

3. The computational fluid dynamics prediction method for the motion state of cylindrical particles in pneumatic transport according to claim 2, characterized in that, The initial state of the continuous phase includes: the continuous phase flow field density, the initial velocity of the continuous phase, and the external force term of the continuous phase. The initial state of the particles includes: particle size, particle volume fraction, initial particle position, particle spatial orientation angle, and initial particle velocity.

4. The computational fluid dynamics prediction method for the motion state of cylindrical particles in pneumatic transport according to claim 2, characterized in that, Based on the lattice Boltzmann equation, the continuous phase distribution function is obtained using the current time step's initial state of the continuous phase, including: The lattice Boltzmann equations, which incorporate the initial state of the continuous phase at the current time step, are transformed into the moment space, and different moments evolve according to different relaxation times. The evolved moments are inversely transformed to obtain the continuous phase distribution function.

5. The computational fluid dynamics prediction method for the motion state of cylindrical particles in pneumatic transport according to claim 2, characterized in that, The calculation of the interaction force between the continuous phase and the particles includes: Obtain the projection of the continuous phase distribution function onto the moment space; By combining the projection and the continuous phase distribution function, a collision process and a migration process are performed to obtain the evolved continuous phase distribution function; Based on the evolved continuous phase distribution function and the no-slip boundary condition, the interaction force between the continuous phase and the particles is obtained.

6. The computational fluid dynamics prediction method for the motion state of cylindrical particles in pneumatic transport according to claim 5, characterized in that, The no-slip boundary conditions of the continuous phase on the particles and the solid walls are determined using an interpolated collision bounce scheme.

7. The computational fluid dynamics prediction method for the motion state of cylindrical particles in pneumatic transport according to claim 2, characterized in that, Based on the interaction force, the lubrication force, and the contact force, the initial state of the particles at the current time step is updated to obtain the updated particle state, including: Execute multiple sub-loop steps, each of which generates a sub-loop update particle state and uses the sub-loop update particle state as the initial state of the particle to enter the next sub-loop step. After completing all the sub-loop steps, the updated particle state is obtained.

8. The computational fluid dynamics prediction method for the motion state of cylindrical particles in pneumatic transport according to claim 7, characterized in that, The sub-loop steps include: Based on the initial state of the particles at the current time step, determine the first positional relationship between the particles and the second positional relationship between the particles and the solid wall. Based on the first positional relationship and the second positional relationship, the lubrication force correction and contact force correction of the particles are determined; Based on the lubrication force correction and the contact force correction, combined with the interaction force, the updated particle state is obtained.

9. The computational fluid dynamics prediction method for the motion state of cylindrical particles in pneumatic transport according to claim 2, characterized in that, Based on the updated particle state, the surface continuous phase distribution function of the particles is reconstructed to obtain the updated continuous phase distribution function, including: Based on the updated particle state, the distance between the fluid grid point of the continuous phase and the surface of the particle, and the normal direction of the surface of the particle are determined near the surface of the particle. Based on the distance and the normal direction, the surface continuous phase distribution function of the particle is reconstructed, and the updated continuous phase distribution function is obtained by combining the surface continuous phase distribution function and the continuous phase distribution function. Based on the updated continuous phase distribution function, the density fluctuation value, continuous phase velocity, and continuous phase pressure of the continuous phase are calculated.

Images