Multi-phase flow wear corrosion prediction method based on CFD-DEM coupling
By using the CFD-DEM coupling method, the limitations of traditional models in predicting wear and corrosion are overcome, enabling accurate prediction of wear and corrosion on industrial pipeline walls, reducing costs and increasing pipeline life.
Patent Information
- Application Number
- CN202511226831.5
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-08-29
- Publication Date
- 2026-01-23
AI Technical Summary
Traditional single CFD or DEM models have limitations when dealing with fluid-particle two-phase flow, making it difficult to accurately describe the interactions between particles and the reaction of the fluid on the particles, resulting in inaccurate prediction of wear and corrosion.
By employing the CFD-DEM coupling method, a set of fluid dynamics and particle dynamics equations is established, and combined with a user-defined API coupling interface, collaborative calculation of fluid and particles is achieved to predict wear and corrosion of industrial pipelines.
Accurately predict the degree of wear and corrosion on the walls of industrial pipelines, reduce wear and corrosion damage, extend pipeline service life, and save computing costs and memory.
Smart Images

Figure CN121389849A_ABST
Abstract
Description
TECHNICAL FIELD
[0001] The application belongs to the technical field of fluid mechanics, and particularly relates to a multiphase flow wear corrosion prediction method based on CFD-DEM coupling. BACKGROUND
[0002] Particle transport phenomenon as the core dynamics process of multiphase flow system, widely exists in the fluidized bed reactor in the field of energy and chemical industry, pneumatic conveying system in the metallurgical industry, and pollution diffusion control in environmental engineering. Through numerical simulation technology, the motion trajectory of the particle phase in the complex flow field can be analyzed, and the micro transport characteristics of the equipment interior can be deeply revealed.
[0003] The CFD-DEM (computational fluid dynamics-discrete element method) coupling method is derived from the actual needs of multiphase flow research. The traditional single CFD or DEM model has limitations in processing fluid-particle two-phase flow: CFD is good at macro fluid flow simulation, but it is difficult to accurately describe the interaction between particles and the reaction of particles to fluid; DEM can capture the motion and collision of particles in detail, but it cannot directly consider the influence of fluid on particles such as drag and pressure gradient. In order to break through this bottleneck, since the 1990s, scholars have combined CFD and DEM, and through bidirectional data transmission, the fluid and particles are calculated cooperatively: CFD solves the fluid domain, DEM calculates the particle motion, and the two are coupled through momentum and energy exchange. This method has gradually become a key tool for studying complex multiphase flow, and provides an efficient numerical means for industrial equipment design optimization and process mechanism revelation. SUMMARY
[0004] In order to solve the problems in the background art, the application provides a multiphase flow wear corrosion prediction method based on CFD-DEM coupling. The application establishes a fluid dynamics equation group considering the force of particles on fluid based on fluid field and particle field, establishes a particle dynamics equation group describing the force and motion of particles and considering the force of fluid on particles, and then realizes CFD-DEM coupling through API user-defined coupling interface.
[0005] The technical scheme adopted by the application is as follows:
[0006] The application comprises the following steps:
[0007] S1, obtaining an industrial pipeline entity structure, and then constructing a pipeline geometric model;
[0008] S2, obtaining particle parameters, and establishing a particle motion model in the pipeline geometric model according to the particle parameters;
[0009] S3, obtaining fluid parameters, and establishing a fluid motion model according to the fluid parameters;
[0010] S4, coupling calculation is performed on the particle motion model and the fluid motion model according to a preset time step, so as to obtain particle motion model states and fluid motion model states at each time step within a preset total time length;
[0011] S5, loss calculation is performed on the pipeline geometry model according to the particle motion model states and the fluid motion model states at each time step, and in combination with a wear and corrosion model, so as to obtain a wear and corrosion loss prediction result of the pipeline geometry model.
[0012] The particle parameters include a velocity of the particle, a mass of the particle, an angular velocity of the particle, a moment of inertia of the particle, a contact force of the particle, a non-contact force of the particle, an interaction force between the particle and the fluid, a volume force received by the particle, a contact torque of the particle, and a contact torque between the particle and the fluid.
[0013] The fluid parameters include a fluid density, a fluid velocity, a fluid viscosity, and a void fraction.
[0014] The particle motion model is set according to the following formula:
[0015]
[0016] F s,l +F drag +F mag +F bas +F pa +F saf
[0017] wherein v represents the velocity of the particle, m represents the mass of the particle, ω represents the angular velocity of the particle, t represents time, I represents the moment of inertia of the particle, F c,l represents the contact force of the lth particle, F nl,l represents the non-contact force of the lth particle, F s,l represents the interaction force between the lth particle and the fluid, F m,l represents the volume force received by the lth particle, M t,l represents the contact torque of the lth particle, M r,l represents the contact torque between the lth particle and the fluid.
[0018] The fluid motion model is set according to the following formula:
[0019]
[0020] F s,l +F drag +F mag +F bas +F pa +F saf
[0021] wherein ρq represents fluid density, μ eff represents fluid viscosity, g represents gravitational acceleration, α q represents porosity, F s,l represents the force between the lth particle and the fluid, p represents pressure, u i and u j respectively represent the velocity components of the fluid velocity in the i direction and the j direction, x j and x i respectively represent the coordinates of the i direction and the j direction in the calculation domain, F drag represents the drag force on the particle, F mag represents the rotational force on the particle, F bas represents the instantaneous resistance force on the particle, F pa represents the pressure gradient force on the particle, F saf represents the force generated by the pressure difference on the particle.
[0022] The step S4 specifically includes:
[0023] S4.1, calculating the contact force between the particles in the particle motion model through discrete element software according to preset particle parameters and fluid parameters;
[0024] S4.2, solving the motion velocity and position of the particle in the particle motion model according to the contact force between the particles and the force of the fluid on the particle at the last time step;
[0025] S4.3, calculating the porosity according to the motion velocity and position of the particle at a preset time step, and then solving the pressure field and velocity field of the fluid in the fluid motion model according to the porosity;
[0026] S4.4, solving the force of the fluid on the particle according to the pressure field and velocity field of the fluid;
[0027] S4.5, repeating steps S4.2 to S4.4 until the simulation time reaches the preset total time, so as to obtain the particle motion model state and the fluid motion model state at each time step within the preset total time;
[0028] The preset time step in the step S4 satisfies the following formula:
[0029]
[0030] where Δt represents the preset time step, Δt' represents the time step calculated when the particle is in a static state, G represents the shear modulus of the particle, η represents the Poisson's ratio of the particle, ρ n represents fluid density, and R represents particle radius.
[0031] The wear corrosion model is set according to the following formula:
[0032]
[0033] wherein, M s represents the mass loss of the wall surface of the abrasion corrosion, k1 represents a constant, h represents the thickness of the oxide layer of the corrosion product, p g represents the density of the oxide, d q represents the diameter of the particle, H s represents the hardness of the oxide protective film of the corrosion product, v q represents the velocity of the particle impact, a1 represents the impact angle, e represents the recovery coefficient of the oxide film, p p represents the density of the particle.
[0034] The particle motion model state of each time step includes the particle contact force, particle velocity, particle position and particle-fluid interaction force of each time step;
[0035] The fluid motion model state of each time step includes the fluid density, viscosity and velocity of each time step.
[0036] The beneficial effects of the present application are as follows:
[0037] (1) The present application innovatively proposes an abrasion corrosion model, which can accurately predict the degree of abrasion corrosion of the wall surface of an industrial pipeline, and has important theoretical value and engineering practical significance for reducing abrasion corrosion damage and improving the service life of the industrial pipeline.
[0038] (2) The present application predicts the abrasion corrosion amount and distribution of the wall surface of the industrial pipeline by means of CFD and DEM coupling, saves the calculation memory, and has excellent expansibility and adaptability;
[0039] (3) The present application is based on CFD and DEM, and through the API coupling model interface, the simulation analysis can save a large amount of time cost for experiments. BRIEF DESCRIPTION OF DRAWINGS
[0040] Figure 1 is a flow schematic diagram of the present application based on CFD-DEM coupling;
[0041] Figure 2 is a particle motion iteration schematic diagram of the present application;
[0042] Figure 3 is an abrasion corrosion diagram of the industrial pipeline obtained by the present embodiment. DETAILED DESCRIPTION
[0043] To further illustrate the embodiments, the present application provides accompanying drawings which form part of the disclosure and which are mainly used to illustrate the embodiments and can be used in conjunction with the relevant description of the specification to explain the operating principles of the embodiments. Those of ordinary skill in the art should be able to understand other possible implementations and advantages of the present application by referring to these contents.
[0044] The embodiment includes the following steps:
[0045] S1, obtaining an industrial pipeline entity structure, and then constructing a pipeline geometric model;
[0046] S2, obtaining particle parameters, and establishing a particle motion model in the pipeline geometric model according to the particle parameters;
[0047] S3, obtaining fluid parameters, and establishing a fluid motion model according to the fluid parameters;
[0048] S4, coupling calculation of the particle motion model and the fluid motion model according to a preset time step, so as to obtain particle motion model states and fluid motion model states at each time step within a preset total simulation time;
[0049] S5, loss calculation of the pipeline geometric model according to the particle motion model states and the fluid motion model states at each time step, and combining a wear and corrosion model, to obtain wear and corrosion loss of the pipeline geometric model.
[0050] Specifically, first, according to the industrial pipeline entity structure, geometric modeling is performed in FLUENT, and then the generated model is meshed. Subsequently, the meshed industrial pipeline model is imported into the discrete element software EDEM. In EDEM, initial particle parameters are set, including physical parameters and topographic parameters. Based on Newton's law of motion, the software analyzes the particle motion state by calculating the contact force between particles;
[0051] The API coupling interface between the discrete element software and the fluid simulation software is loaded in the processor, so that the particle field information in the discrete element software is imported into the fluid simulation software through the coupling interface. Then, the wear and corrosion model is loaded in the discrete element software EDEM in the processor, initial fluid information is given in the fluid simulation software, the initial information of the flow field includes the initial velocity and physical conditions of the fluid, the calculation control parameters are set, and the fluid dynamics equation set is established.
[0052] The void fraction is input into the CFD-DEM module at a predetermined time step interval to update the fluid grid information. The obtained flow field information is transmitted back to EDEM through the API coupling interface. EDEM combines the loaded corrosion and wear model and the particle motion information affected by the flow field to calculate the wear / corrosion mass loss of the industrial pipeline wall surface.
[0053] The industrial pipeline entity structure includes pipeline segments of three spatial orientations: horizontal, vertical and inclined;
[0054] The particle parameters include particle velocity, particle mass, particle angular velocity, particle moment of inertia, particle contact force, particle non-contact force, particle-fluid interaction force, particle volume force, particle contact torque and particle-fluid contact torque.
[0055] The fluid parameters include fluid density, fluid velocity, fluid viscosity and void fraction.
[0056] The particle motion model is set according to the following formula:
[0057]
[0058] F s,l +F drag +F mag +F bas +F pa +F saf
[0059] Wherein, v represents the velocity of the particle, m represents the mass of the particle, ω represents the angular velocity of the particle, t represents time, I represents the moment of inertia of the particle, F c,l represents the contact force of the lth particle, F nl,l represents the non-contact force of the lth particle, F s,l represents the interaction force of the lth particle and the fluid, F m,l represents the volume force received by the lth particle, M t,l represents the contact torque of the lth particle, M r,l represents the contact torque of the lth particle and the fluid.
[0060] As Figure 1 the red box on the left and Figure 2 shown, the particle motion is iterated in the following steps:
[0061] Step A1, the EDEM software calculates the motion of the pipeline geometry, and changes the initial distribution and initial state of the particles according to the preset conditions;
[0062] Step A2, then the particle unit retrieval is performed on the generated particles, and the contact force between the particles is calculated, the contact model API is executed at this place, and the calculation of various contact properties between the particles and the geometry is realized according to the program compiled by the user;
[0063] Step A3, then the EDEM software calculates the field force involved in the motion of the particles (such as the configured gravity field, magnetic field force, etc.);
[0064] Step A4, the volume force acting on each particle is calculated, the particle volume force API will be executed at this point, according to the user-compiled program to add different volume force model for the particle;
[0065] Step A5, after the above calculation is completed, EDEM will be based on the various forces acting on the particle to calculate the particle resultant force, and then use Newton's second law to obtain the particle velocity and acceleration information and update the particle position information.
[0066] Step A6, finally EDEM will calculate the cohesive force acting on the particle.
[0067] After the completion of a time step, if the simulation time has not reached the user-configured preset value, EDEM will start the calculation of the next time step.
[0068] In specific implementation, the discrete element software EDEM calculates and sets the angular velocity and linear velocity of the particle through its particle motion model. The particle data containing these motion states are transmitted to the fluid simulation software through the coupling interface.
[0069] The fluid motion model is set according to the following formula:
[0070]
[0071] F s,l = F drag + F mag + F bas + F pa + F saf
[0072] Where, ρ q represents the fluid density, μ eff represents the fluid viscosity, g represents the acceleration of gravity, α q represents the porosity, F s,l represents the force acting on the lth particle and the fluid, p represents the pressure, u i and u j respectively represent the velocity components of the fluid velocity in i direction and j direction, x j and x i respectively represent the coordinates of the calculation domain in i direction and j direction, F drag represents the drag force acting on the particle, F mag represents the rotational force acting on the particle, F bas represents the instantaneous resistance force acting on the particle, F pa represents the pressure gradient force acting on the particle, F saf represents the force generated by the pressure difference acting on the particle.
[0073] The calculation control parameter in the embodiment is set by using any one of the particle size, the particle shape, and the liquid inlet speed.
[0074] As shown in Figure 1 , step S4 is specifically:
[0075] S4.1, calculating the contact force between the particles in the particle motion model by the discrete element software according to the preset particle parameters and fluid parameters;
[0076] S4.2, solving the motion speed and position of the particles by the Newton's second law in the particle motion model according to the contact force between the particles and the force of the fluid on the particles at the previous time step;
[0077] S4.3, calculating the porosity according to the motion speed and position of the particles at a preset time step, and then solving the pressure field and velocity field of the fluid in the fluid motion model according to the porosity;
[0078] S4.4, solving the force of the fluid on the particles according to the pressure field and velocity field of the fluid, which is used for the calculation of the motion speed and position of the particles at the next time step;
[0079] S4.5, repeating steps S4.2 to S4.4 until the simulation time reaches the preset total time, so as to obtain the particle motion model state and the fluid motion model state at each time step within the preset total time;
[0080] The force of the fluid on the particles in the calculation of the motion speed and position of the particles at the first time step is zero.
[0081] Further, step S4 is performed as follows:
[0082] G1, giving initial flow field information and initial particle information;
[0083] G2, calculating the contact force between the particles by the discrete element software;
[0084] G3, solving the motion speed and position of the particles by the Newton's second law;
[0085] G4, transferring the speed and position information of the particles to the fluid simulation software for solving;
[0086] G5, inputting the porosity to the CFD-DEM module at a predetermined time step interval to update the fluid grid;
[0087] G6, solving the momentum equation according to the velocity field and the pressure field to update the fluid pressure field and velocity field;
[0088] G7, solving the momentum exchange source term of the particles-fluid to update the force of the fluid on the particles;
[0089] G8, repeat the above steps until the simulation reaches a pre-set time.
[0090] The time step preset in step S4 using the discrete element method satisfies the following formula:
[0091]
[0092] where Δt represents the preset time step, Δt' represents the time step calculated when the particle is in a static state, G represents the shear modulus of the particle, η represents the Poisson's ratio of the particle, ρ n represents the fluid density, and R represents the particle radius.
[0093] The erosion-corrosion model is set according to the following formula:
[0094]
[0095] where M s represents the mass loss of the erosion-corrosion wall, k1 represents a constant, h represents the oxide layer thickness of the corrosion product, ρ g represents the oxide density, d q represents the particle diameter, H s represents the hardness of the oxide protective film of the corrosion product, v q represents the particle impact velocity, α1 represents the impact angle, e represents the recovery coefficient of the oxide film, and ρ p represents the particle density.
[0096] In this embodiment, the oxide density, the constant, the impact angle, the hardness of the oxide protective film of the corrosion product, and the recovery coefficient of the oxide film are all pre-set values, the particle diameter and the particle density are specifically parameters pre-set in the particle motion model, the particle impact velocity is calculated according to the particle motion model state of each time step, the oxide layer thickness of the corrosion product is calculated according to the particle motion model state and the fluid motion model state of each time step, and the oxide layer thickness is specifically determined according to the particle contact force, the velocity, the position, and the force of the particle on the fluid in the particle motion model state of each time step, and in combination with the fluid density, the viscosity, and the velocity in the fluid motion model state, the growth and erosion removal of the oxide layer are comprehensively determined, so as to be calculated.
[0097] The particle motion model state of each time step includes the particle contact force, the particle velocity, the particle position, and the force of the particle on the fluid in each time step.
[0098] The fluid motion model state of each time step includes the fluid density, the viscosity, and the velocity in each time step.
[0099] The erosion-corrosion loss of the pipeline geometry model is specifically the pipe wall erosion-corrosion loss amount of the pipeline geometry model.
[0100] The mass loss amount graph of the wall surface of the wear corrosion is as shown in the figure. Figure 3
[0101] The API coupling interface includes a header file (.h) of interface definition and a source file (.cpp) of realizing core functions, wherein the header file structure includes: 1) copyright and version information; 2) preprocessing instructions; 3) function prototype definition; 4) class structure declaration; the source file firstly receives parameter variables generated by an EDEM simulation process or preset by a user, then performs calculation through a self-defined contact algorithm, and finally feeds back results to an EDEM main program.
[0102] Firstly, the header file includes Helper.h and CCorrosionWear.h, wherein the main function of Helper.h is to declare auxiliary functions and classes for optimizing calculation efficiency, and it needs to be added to a directory where the header file is located when compiling the API; the main function of CCorrosionWear.h is to declare functions, classes and variables of the contact model API required. In addition, the Core, ContactModels and Misc files in the EDEM installation directory must be added to the additional include directory of the compiled contact model API. Secondly, the source file includes CorrosionWear.cpp and CCorrosionWear.cpp, wherein: the main function of CorrosionWear.cpp is to realize version compatibility identification of the API and the EDEM main program; the main function of CCorrosionWear.cpp is to specifically write the functions of the contact model API declared in CCorrosionWear.h, and return the results to EDEM after calculation. Finally, after the API is compiled, a dll file will be generated, which needs to be placed together with the preconfigured txt file in the EDEM working directory, and the user can load the self-defined contact model API through the physics module.
[0103] To sum up, based on CFD and DEM and through the API coupling model interface, the simulation analysis of the application can save a lot of time cost of experiments; the application saves calculation memory and has excellent expansibility and adaptability by means of the CFD and DEM coupling means to predict the wear corrosion amount and wear corrosion distribution of the wall surface of the industrial pipeline; the application can more accurately predict the wear corrosion degree of the wall surface of the industrial pipeline, and has important theoretical value and engineering practical significance for reducing wear corrosion damage and improving the service life of the industrial pipeline.
[0104] It should be noted that the above embodiments are only used to illustrate the technical solutions of the present application, and are not intended to limit the present application; although the present application has been described in detail with reference to the above embodiments, those skilled in the art should understand that the technical solutions recorded in the above embodiments can still be modified, or some or all of the technical features can be replaced by equivalents; and these modifications or replacements do not make the essence of the corresponding technical solutions deviate from the scope of the technical solutions of the embodiments of the present application.
Claims
1. A method for predicting multiphase flow erosion corrosion based on CFD-DEM coupling, characterized in that, The method comprises the following steps: S1, obtaining an industrial pipeline entity structure, and then constructing a pipeline geometric model; S2, obtaining particle parameters, and establishing a particle motion model in the pipeline geometric model according to the particle parameters; S3, obtaining fluid parameters, and establishing a fluid motion model according to the fluid parameters; S4, coupling calculation of the particle motion model and the fluid motion model according to a preset time step, so as to obtain particle motion model states and fluid motion model states at each time step within a preset total time length; S5, loss calculation of the pipeline geometric model according to the particle motion model states and the fluid motion model states at each time step, and combining a wear and corrosion model, so as to obtain a wear and corrosion loss prediction result of the pipeline geometric model.
2. The method of claim 1, wherein the method is characterized by: The particle parameters comprise particle velocity, particle mass, particle angular velocity, particle moment of inertia, particle contact force, particle non-contact force, particle-fluid interaction force, particle volume force, particle contact torque, and particle-fluid contact torque; The fluid parameters comprise fluid density, fluid velocity, fluid viscosity, and void fraction.
3. The method of claim 1, wherein the method is characterized by: The particle motion model is set according to the following formula: F s,l = F drag + F mag + F bas + F pa + F saf where v represents the velocity of the particle, m represents the mass of the particle, ω represents the angular velocity of the particle, t represents time, I represents the moment of inertia of the particle, F c,l represents the contact force of the lth particle, F nl,l represents the non-contact force of the lth particle, F s,l represents the force of the lth particle acting on the fluid, F m,l represents the volume force experienced by the lth particle, M t,l represents the contact moment of the lth particle, M r,l represents the contact moment of the lth particle acting on the fluid.
4. The method of claim 1, wherein the method is characterized by: The fluid motion model is set according to the following formula: F s,l = F drag + F mag + F bas + F pa + F saf where p q represents fluid density, m eff represents fluid viscosity, g represents gravitational acceleration, a q represents porosity, F s,l represents the force of the lth particle and fluid, p represents pressure, u i and u j respectively represent the velocity components of the fluid velocity in i direction and j direction, x j and x i respectively represent the coordinates in i direction and j direction in the calculation domain, F drag represents the drag force on the particle, F mag represents the rotational force on the particle, F bas represents the instantaneous resistance force on the particle, F pa represents the pressure gradient force on the particle, F saf represents the force generated by the pressure difference on the particle.
5. The method of claim 1, wherein the method is characterized by: The step S4 specifically comprises: S4.1, calculating contact forces between particles in the particle motion model by discrete element software according to preset particle parameters and fluid parameters; S4.2, solving particle motion velocity and position in the particle motion model according to the contact forces between particles and fluid forces on particles at a previous time step; S4.3, calculating a void fraction according to the particle motion velocity and position at a preset time step, and then solving a fluid pressure field and velocity field in the fluid motion model according to the void fraction; S4.4, solving fluid forces on particles according to the fluid pressure field and velocity field; S4.5, repeating steps S4.2 to S4.4 until a simulation time length reaches the preset total time length, so as to obtain particle motion model states and fluid motion model states at each time step within the preset total time length.
6. The method of claim 1, wherein the method is based on CFD-DEM coupling. The preset time step in the step S4 satisfies the following formula: where Δt represents a preset time step, Δt' represents a time step calculated when the particle is in a static state, G represents a shear modulus of the particle, η represents a Poisson's ratio of the particle, p n represents a fluid density, and R represents a particle radius.
7. The method of claim 1, wherein the method is characterized by: The wear and corrosion model is set according to the following formula: where M s represents the mass loss amount of the wall surface of the wear corrosion, k1 represents a constant, h represents the thickness of the oxide layer of the corrosion product, p g represents the oxide density, d q represents the particle diameter, H s represents the hardness of the oxide protective film of the corrosion product, v q represents the particle impact velocity, a1 represents the impact angle, e represents the restitution coefficient of the oxide film, p p represents the particle density.
8. The method of claim 1, wherein the method is a CFD-DEM coupled multiphase flow erosion-corrosion prediction method. The particle motion model states at each time step comprise particle contact forces, particle velocity, particle position, and particle-fluid interaction forces at each time step; The fluid motion model states at each time step comprise fluid density, viscosity, and velocity at each time step.
Citation Information
Cited By
Method for predicting thermal wear of pump for direct coal liquefaction by considering high-temperature effect
CN121706505A
Method for predicting thermal wear of pump for direct coal liquefaction considering high temperature effect
CN121706505B
