Parameter decision method for pneumatic cleaning of side deep fertilizing device based on multi-objective optimization

Abstract

The application relates to the technical field of agricultural machinery, and discloses a pneumatic adhesion removal parameter decision method of a side deep fertilization device based on multi-target optimization, which comprises the following steps: preparing a saturated soil sample to extract mechanical parameters, and calculating total liquid bridge force of a mud particle interface; a cohesive contact model is constructed in discrete element software, surface energy parameters are corrected through a virtual calibration test; a critical adhesion removal stress balance equation containing air flow drag force and total liquid bridge force is established, a gas-solid two-way coupling control equation is combined to simulate an adhesion removal process, and air flow adhesion removal efficiency and system aerodynamic power consumption are extracted; a back propagation neural network is trained with an air flow speed at a nozzle inlet, an air flow incidence angle and a normal distance from the nozzle to the ground as control variables, a multi-target comprehensive evaluation function is constructed, a particle swarm optimization algorithm is used to solve an extreme value, and an optimal control variable combination is output to control the device to execute air flow adhesion removal; and the application improves the decision precision of parameters.

Description

Technical Field

[0001] This invention relates to the field of agricultural machinery technology, specifically to a multi-objective optimization-based method for determining the pneumatic cleaning parameters of a side-deep fertilization device. Background Technology

[0002] When rice side-deep fertilization devices are used in paddy fields, the high soil moisture content causes mud to easily adhere to the device walls. Mud accumulation can clog the fertilization channels, increasing operational resistance and affecting fertilizer discharge. Using airflow generated by a pneumatic system to remove the adhering mud is a common existing treatment method.

[0003] Current technologies often rely on manual adjustments when setting airflow removal parameters. Due to the presence of liquid bridging forces and viscous forces between paddy field mud particles, conventional discrete element method (DEM) simulations of soil typically fail to calculate the liquid bridging effect at the mud interface. This leads to deviations in the force simulation of the interaction between airflow and mud, resulting in simulation results that cannot accurately reflect the mud desorption state. Furthermore, pneumatic removal efficiency and pneumatic power consumption are mutually constrained. Existing parameter setting methods do not incorporate multi-objective calculations based on variables such as airflow velocity, incident angle, and nozzle position, lacking a decision-making mechanism that comprehensively considers removal efficiency and power consumption. This makes it difficult for fertilization devices to determine the optimal airflow control parameters. Summary of the Invention

[0004] To address the shortcomings of existing technologies, this invention provides a multi-objective optimization-based method for determining the pneumatic cleaning parameters of a side-deep fertilization device. This method solves the problems of mud adhesion easily occurring in existing side-deep fertilization devices when operating in saturated paddy fields, and the difficulty in accurately determining the control parameters of the pneumatic cleaning system.

[0005] To achieve the above objectives, the present invention provides the following technical solution: This invention provides a multi-objective optimization-based method for determining pneumatic cleaning parameters in side-deep fertilization devices, comprising the following steps: Prepare a saturated soil sample composed of mud particles with a continuous water film. Control the test probe to contact and separate from the saturated soil sample, extract the maximum normal tensile force and the maximum tangential sliding force, and calculate the normal and tangential adhesion specific pressure of the saturated soil sample. The total liquid bridging force at the mud particle contact interface was calculated, and a cohesive contact model with surface energy parameters was constructed in discrete element software. The surface energy parameters of the cohesive contact model were calibrated through virtual calibration experiments, and a discrete element model of saturated paddy field soil was constructed. A critical desorption force balance equation for mud particles, including gravity, airflow drag, interparticle adhesion contact force and total liquid bridge force, was established. The cleaning process of mud particles on the wall of the side-deep fertilization device was simulated by combining the gas-solid bidirectional coupling control equation, and the pneumatic cleaning efficiency and system pneumatic power consumption were extracted. The nozzle inlet airflow velocity, airflow injection angle, and nozzle normal distance from the ground of the side-deep fertilization device are used as control variables. Based on the control variables, the corresponding pneumatic cleaning efficiency, and the system aerodynamic power consumption, a backpropagation neural network is trained to establish a nonlinear prediction model and construct a multi-objective comprehensive evaluation function. Based on the nonlinear prediction model, the particle swarm optimization algorithm is used to solve the extreme value of the multi-objective comprehensive evaluation function, and the optimal combination of control variables is output as the airflow cleaning control parameters of the side-deep fertilization device to control the side-deep fertilization device to perform pneumatic cleaning operation.

[0006] In the above method, the calibration of soil contact mechanics parameters is achieved by calculating the liquid bridge force. After obtaining the solid-liquid contact angle, dynamic half-gap thickness, and liquid phase surface tension, the static capillary force at the mud particle contact interface is obtained based on the annular approximation. The dynamic viscous force is then obtained by combining the relative separation velocity and the Reynolds equation. The total liquid bridge force is obtained by linearly superimposing the static capillary force and the dynamic viscous force. The initial surface energy parameters of the adhesive contact model are determined by combining the equivalent contact radius of the mud particles. By comparing the bottom diffusion diameter and center height of a cylindrical measuring cylinder filled with actual saturated soil samples and a virtual measuring cylinder after collapsing and expanding under gravity, the surface energy parameters are continuously iteratively adjusted until the relative error between the simulation data and the actual data is simultaneously less than the set calibration error threshold, thus achieving the calibration of the discrete element material library parameters.

[0007] When determining the slurry desorption state, the force components of the slurry particle translational control equation in the wall normal direction are extracted. When the sum of the outward component of the airflow drag force and the outward component of gravity in the wall normal direction is greater than or equal to the sum of the inward component of the interparticle adhesive contact force and the inward component of the total liquid bridge force in the wall normal direction, the slurry particles are deemed to meet the critical desorption condition and transition to a dynamic detachment state. During the execution of the gas-solid bidirectional coupling control equation, the fluid dynamics software traverses the total volume of slurry particles within the fluid computation grid at each computation time step, calculates the local porosity, and performs volume reduction corrections on the airflow mass conservation equation. Simultaneously, the reaction forces of the slurry particles within the grid are volume-averaged to construct a two-phase momentum exchange source term, which is then substituted into the Navier-Stokes equations with local porosity to update the momentum conservation equation, achieving bidirectional data synchronization between the flow field and particle motion state.

[0008] In the decision optimization stage, Latin hypercube sampling is used to generate sample points within the range of control variables, constructing a multi-layer feedforward backpropagation neural network. Using nozzle inlet airflow velocity, airflow injection angle, and nozzle-to-ground normal distance as inputs, and pneumatic cleaning efficiency and system aerodynamic power consumption as outputs, a nonlinear prediction model is established by updating network parameters layer by layer using gradient descent. The prediction output is normalized dimensionlessly using the range transformation method, and a multi-objective comprehensive evaluation function is constructed by combining it with set target weight coefficients. Using the particle swarm optimization algorithm, parameter combinations are defined as particles, and the multi-objective comprehensive evaluation function is used as the fitness evaluation standard. Through real-time position updates of individual cognitive and social experience factors, the globally optimal position of the population is found as the theoretical optimal solution. The numerical values, after engineering rounding, are sent to the onboard programmable logic controller of the side-deep fertilization device, directly adjusting the air compressor output frequency and the stroke of the nozzle pose adjustment motor.

[0009] This invention provides a multi-objective optimization-based method for determining the pneumatic cleaning parameters of a side-deep fertilization device. It offers the following advantages: 1. This invention calculates the total liquid bridge force by extracting the maximum normal tensile force and maximum tangential sliding force from a saturated soil sample with a continuous water film, superimposing static capillary force and dynamic viscous force, and then using a virtual slump test to calibrate the surface energy parameters in the discrete element model. This method introduces the liquid bridge effect of the mud interface into the discrete element simulation for parameter calibration, solving the problem of large errors in conventional discrete element parameters when simulating saturated mud in paddy fields, and improving the calculation accuracy of the discrete element model for paddy field soil.

[0010] 2. This invention establishes the critical desorption force balance equations, which include gravity, airflow drag, adhesive contact force, and total liquid bridging force. Within each calculation time step, local porosity and two-phase momentum exchange source terms are introduced to perform gas-solid bidirectional coupling calculations. This calculation process establishes the force boundary conditions for mud desorption, simulates the momentum transfer and flow field changes between the airflow and mud particles, and provides a computational basis for extracting data on pneumatic removal efficiency and system aerodynamic power consumption.

[0011] 3. This invention utilizes a backpropagation neural network to establish a nonlinear prediction model between control variables and system performance, and combines a multi-objective comprehensive evaluation function with a particle swarm optimization algorithm for extreme value solving. This method replaces manual experience-based debugging with an optimization algorithm, balancing the improvement of cleaning efficiency with the reduction of aerodynamic power consumption; moreover, the output optimal solution, after engineering rounding, is directly used as a control command to be issued to the programmable logic controller, improving the efficiency and accuracy of airflow control parameter decision-making for the side-deep fertilization device. Attached Figure Description

[0012] Figure 1 This is a schematic diagram of the method flow of the present invention; Figure 2This is a schematic diagram of the system architecture of the present invention. Detailed Implementation

[0013] The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0014] Please see the appendix Figure 1 This invention provides a multi-objective optimization-based decision system for pneumatic cleaning parameters of a side-deep fertilization device, comprising: The macroscopic testing module is used to test the normal adhesion specific pressure and tangential adhesion specific pressure of saturated paddy field soil samples.

[0015] The mesoscopic characterization module is used to extract capillary and viscous force parameters at the interface of mud particles and construct a discrete element contact model for non-Newtonian fluids.

[0016] The coupled simulation module is used to construct a three-dimensional model of the side-deep fertilization device and the gas-solid bidirectional coupling control equations, and to perform numerical simulation calculations of the airflow desorption process of the adhering mud.

[0017] The parameter decision module is used to establish a nonlinear prediction model based on the dataset output by the coupled simulation module and solve the comprehensive evaluation objective function through an optimization algorithm to output the pneumatic cleaning parameters.

[0018] Please see the appendix Figure 2 This invention provides a method for determining pneumatic cleaning parameters of a side-deep fertilization device based on multi-objective optimization, comprising the following steps: S1. A saturated soil sample with a preset moisture content is prepared through a macroscopic testing module, and a continuous water film of a preset thickness is placed on the surface of the saturated soil sample. The test probe is controlled to contact and separate from the saturated soil sample at a preset loading and unloading speed. The maximum normal tensile force and the maximum tangential sliding force during the contact process of the test probe are recorded. The normal adhesion specific pressure and the tangential adhesion specific pressure are calculated based on the maximum normal tensile force, the maximum tangential sliding force and the effective contact area between the test probe and the saturated soil sample. S2. The liquid bridge force parameters between mud particles and between mud particles and soil contact components are obtained through the mesoscopic characterization module. The liquid bridge force parameters include solid-liquid contact angle, dynamic half-gap thickness and liquid phase surface tension. The total liquid bridge force at the contact interface is calculated based on the liquid bridge force parameters. Based on the total liquid bridge force, a cohesive contact model with surface energy properties is configured in the discrete element software. The parameters of the cohesive contact model are calibrated through virtual calibration test to complete the construction of the discrete element model of saturated paddy field soil. S3 establishes the critical desorption force balance equation of mud particles, which includes gravity, total liquid bridge force, airflow drag force and interparticle adhesion contact force, through the coupled simulation module. The gas field in the fluid dynamics software is set as the continuous phase and the mud particle field in the discrete element software is set as the discrete phase. A data interaction interface between the fluid dynamics software and the discrete element software is established. The momentum exchange calculation is performed according to the gas-solid bidirectional coupling control equation to simulate the dynamic evolution process of high-pressure airflow acting on the surface of the side-deep fertilization device. Soil adhesion quality and system airflow power consumption data under different parameter configurations are extracted. S4 extracts the airflow nozzle inlet velocity, nozzle structure angle, and device operating speed as input variables through the parameter decision module, and extracts soil adhesion quality and system airflow power consumption as output indicators. Based on the input variables and output indicators, a neural network prediction model is trained to establish a nonlinear mapping relationship. A comprehensive evaluation objective function guided by minimizing soil adhesion quality and minimizing system airflow power consumption is established, and weight coefficients are set. The minimum value of the comprehensive evaluation objective function is solved in the constraint domain through response surface analysis and optimization algorithm. The airflow nozzle inlet velocity and nozzle structure angle corresponding to the minimum value are output as the airflow cleaning control parameters of the side-deep fertilization device.

[0019] The process of mud sample preparation and state simulation in macroscopic mud adhesion characteristic testing includes the following:

[0020] Unspoiled paddy field soil from the target work area was obtained and pretreated. To eliminate the interference of the original soil moisture content and compaction history on the test results and to establish a unified physical testing benchmark, plant roots and stone impurities in the unspoiled paddy field soil needed to be removed. The cleaned soil was then placed in a ventilated environment to air dry naturally until its mass stabilized. The air-dried soil was then pulverized using a mechanical pulverizer, and the pulverized soil particles were passed through a standard test sieve with a 2 mm aperture to obtain the dry soil sample for the experiment. For the specific structure and operation control logic of the mechanical pulverizer, those skilled in the art can refer to conventional agricultural soil testing equipment, as its physical principles and mechanical transmission methods are well-known technologies in the field and will not be elaborated here.

[0021] Calculate and weigh a predetermined mass of dry soil sample and deionized water to prepare multiple sets of mixed samples with different moisture content gradients. The formula for calculating the mass moisture content is: ; In the formula, Moisture content (by weight); The mass of deionized water added is in grams. The mass of the dry soil sample used in the experiment is expressed in grams. The gradient range of moisture content was set based on historical statistical data of the side-deep fertilization operation environment. Specifically, the liquid limit and plastic limit of the soil in the target operation area were determined using a conventional Ateburg limit test, and the gradient range of moisture content was set to cover the range from the plastic limit to the liquid limit of the soil. To accurately reflect the consistency changes of the field slurry, this gradient range of moisture content is typically taken as 30% to 60%. The dry soil sample was mixed with deionized water using a planetary soil mixer and stirred until a homogeneous slurry state was achieved.

[0022] The well-mixed mud sample was dispensed into cylindrical standard sample boxes and leveled. The boxes were then sealed. The sealed boxes were placed in a constant temperature and humidity chamber and allowed to stand for a preset time gradient. This preset time gradient simulates different soaking periods of paddy field soil after actual tillage. During the settling process, the clay minerals within the mud fully absorb water and swell, the bound water film between particles gradually thickens and stabilizes, accompanied by minor particle settling and structural reorganization. To align with the typical agronomic waiting period from rotary tillage to transplanting and fertilization, the preset time gradient was set from 12 to 72 hours, with multiple samples taken at 12-hour intervals.

[0023] After settling, the cylindrical standard sample box was unsealed, and a water film interface was constructed on the upper surface of the mixed sample. Physically, after land preparation and settling, saturated paddy field mud typically undergoes phase separation due to gravity settling, resulting in a thin layer of free water. When the soil-contacting component of the side-deep fertilization device cuts into the mud, it first contacts and disrupts this water layer; therefore, the presence of this water layer directly affects the wetting and adhesion mechanism of the initial interface. Deionized water was uniformly applied to the upper surface of the mixed sample using a micro-liquid phase jetting device, forming a continuous water film with a thickness of 1 mm. To precisely control the thickness of the continuous water film, the volume of deionized water applied in a single application was calculated using the following formula: ; In the formula, The volume of deionized water applied is expressed in cubic centimeters. The area of ​​the internal base of the cylindrical standard sample box is expressed in square centimeters. The target thickness of the continuous water film is set at 0.1 cm. The resulting mixed sample with the continuous water film is the saturated soil sample, which serves as the physical benchmark for macroscopic interfacial adhesion mechanics testing. For the hydraulic circuit layout of the micro-liquid phase injection device, those skilled in the art can utilize conventional precision fluid control valve assemblies; its piping distribution is well-known in the field and will not be elaborated upon here.

[0024] The process of performing interfacial adhesion mechanics tests and extracting data in the macroscopic mud adhesion characteristics test includes the following:

[0025] A soil adhesion tester and a test probe are configured. The test probe is made of the same metal material as the soil-contacting parts of the side-deep fertilization device; the specific material is stainless steel or corrosion-resistant carbon steel, selected according to the actual manufacturing standards of the agricultural machinery. The contact surface of the test probe is mechanically ground to ensure its surface roughness matches that of the actual soil-contacting parts of the agricultural machinery, with the roughness range set between 1.6 micrometers and 6.3 micrometers. A prepared saturated soil sample with a continuous water film is fixed on the test platform of the soil adhesion tester. For the hardware layout of the drive module and force sensor of the soil adhesion tester, those skilled in the art can use the control architecture of a conventional universal testing machine; its data acquisition and motion control modules are well-known technologies in the field and will not be described further here.

[0026] Perform a normal adhesion force test. This test simulates the force process experienced by the fertilization component of a side-deep fertilization device during vertical insertion and withdrawal from the mud. The core mechanism involves determining the sum of the water film tension at the mud interface and the normal adsorption force of the soil particles that the metal component must overcome when detaching from the mud. The test probe is lowered vertically at a set loading speed until it contacts the water film interface on the surface of the saturated soil sample. After contact, the test probe applies a preset normal pressure force to the saturated soil sample and maintains the contact for a set time. The preset normal pressure force is determined based on the physical compressive stress generated by the side-deep fertilization device's own weight in the mud and its insertion depth, and its value ranges from 5 kPa to 20 kPa. The set contact time ranges from 2 seconds to 10 seconds, simulating the instantaneous contact state during agricultural machinery movement and allowing free water in the mud to form a stable solid-liquid adsorption interface on the probe surface. After the contact time, the test probe is raised vertically at a set unloading speed to separate it from the saturated soil sample. During the separation process, the force sensor of the measuring instrument records the mechanical change curve of the normal stretching process in real time, and the peak value on the curve is extracted as the maximum normal tensile force. The set unloading speed is determined by conversion based on the actual mechanical operation speed of the side-deep fertilization device, and the value ranges from 0.5 m / s to 1.5 m / s.

[0027] Perform a tangential adhesion test. This test simulates the adhesion and frictional coupling process of a furrow opener or guide component in a side-deep fertilization device as it slides horizontally through mud, determining the horizontal resistance required to shear the water film and mud layer on the component surface. Replace the sample with a saturated soil sample from the same batch and fix it in place. Control the test probe to descend vertically and contact the saturated soil sample, applying the same preset normal pressure. Maintaining the preset normal pressure under dynamic contact, control the test probe to move uniformly across the surface of the saturated soil sample at a set sliding speed. The sliding speed should be consistent with the previously set unloading speed. During horizontal movement, the tangential resistance change curve in the horizontal direction is recorded in real time using a force sensor. The maximum threshold before the relative sliding abrupt change in tangential resistance is extracted as the maximum tangential sliding force. To eliminate random errors caused by local inhomogeneities in the soil sample, both the normal and tangential adhesion tests should be repeated 3 to 5 times under the same test parameters. The arithmetic mean of the corresponding maximum peak values ​​is taken as the final recorded maximum normal tension and maximum tangential sliding force.

[0028] Calculate the interfacial adhesion specific pressure. Since the absolute value of the externally measured adhesion force is directly affected by the probe contact area, to establish a universal evaluation index unrestricted by physical dimensions, it is necessary to eliminate the differences caused by the probe's geometric dimensions. Based on the extracted maximum normal tensile force and maximum tangential sliding force, combined with the physical contact parameters of the test probe, calculate the normal adhesion specific pressure and tangential adhesion specific pressure. The specific calculation formulas are as follows: ; ; In the formula, Normal adhesion specific pressure, in Pascals; The maximum normal tensile force measured in the normal adhesion test is expressed in Newtons. The effective contact area between the probe and the saturated soil sample is measured in square meters. This is the tangential adhesion specific pressure, measured in Pascals. The maximum tangential sliding force measured in the tangential adhesion test is expressed in Newtons.

[0029] Effective contact area The effective contact area is determined based on the actual geometry of the bottom of the test probe. For a cylindrical flat-bottomed test probe, the formula for calculating the effective contact area is: ; In the formula, The effective contact area between the probe and the saturated soil sample is measured in square meters. The radius of the bottom surface of the cylindrical flat-bottomed test probe is measured in meters. The bottom radius is typically set between 0.01 meters and 0.05 meters to ensure that the contact surface contains a sufficient number of soil aggregates. The square of the base radius of the cylindrical flat-bottomed test probe; Pi is a constant with a value of approximately 3.14. The calculated normal and tangential adhesion specific pressures serve as fundamental physical quantities for quantitatively characterizing the macroscopic adhesion properties of saturated paddy field soil, guiding the setting and verification of contact parameters and surface energy parameters in subsequent mesoscopic discrete element contact models.

[0030] The process of characterizing and measuring nanoscale liquid bridge forces in the establishment of mesoscopic liquid bridge force models and the calibration of discrete element parameters includes the following:

[0031] A nano-multifunctional testing machine equipped with a microscopic vision acquisition module is described. In the microstructure of saturated paddy soil, soil particles are completely encapsulated by water. When agricultural machinery parts separate from the soil or relative displacement occurs between soil particles, the water at the particle interface is stretched, forming liquid bridges. The capillary and viscous forces generated by these liquid bridges are the microscopic mechanical causes of the high adhesion and non-Newtonian fluid properties of saturated mud. The purpose of obtaining these microscopic mechanical parameters is to reveal the essence of macroscopic adhesion phenomena and to provide physical quantity basis for subsequently reconstructing a virtual contact model in discrete element software that can realistically reflect the flow and aggregation characteristics of mud. The nano-multifunctional testing machine includes micro-nano manipulation probes for holding tiny particles and a charge-coupled device (CCD) microscopic observation system for observing the morphology of liquid bridges. For the specific physical implementation of the micro-displacement drive mechanism and high-precision mechanical sensors of the nano-multifunctional testing machine, those skilled in the art can use the underlying hardware architecture of conventional atomic force microscopes or micro / nano indenters. Their electromechanical control logic is well-known in the field and will not be elaborated here.

[0032] Single-particle test samples were prepared and microscopic liquid bridges were constructed. Representative soil particles with regular morphology were selected from dried soil samples, and their equivalent geometric dimensions were measured and recorded. The representative soil particles were fixed to the tip of a micro / nano manipulation probe of a nano-multifunctional testing machine. To ensure the particles did not detach during high-speed tensile testing, a very small amount of epoxy resin adhesive was applied to the tip of the micro / nano manipulation probe using a high-precision dispensing machine, bonding the representative soil particles to the probe tip and allowing them to cure at room temperature. A metal substrate made of the same material as the soil contact component of the side-deep fertilization device, or another fixed soil particle, was placed on the testing platform of the nano-multifunctional testing machine. A predetermined volume of deionized water was injected between the representative soil particle and the metal substrate, or between two soil particles, using a micro-injection pump to form an initial liquid bridge connecting the two. To ensure the initial liquid bridge's morphology was stable and would not collapse due to gravity, the injected deionized water volume was set from 0.5 μL to 5.0 μL, depending on the size of the representative soil particles.

[0033] Static contact parameters were extracted and capillary force was calculated. Capillary force is a static adsorption force caused by the surface tension of the liquid phase and the wetting characteristics of the interface; its essence is the additional pressure generated at the liquid-gas-liquid interface due to bending. A nano-multifunctional experimental machine was kept stationary, and a charge-coupled device (CCD) microscopic observation system was used to capture the lateral projection image of the initial liquid bridge. The tangent at the intersection of the liquid bridge and the solid surface was extracted using image processing algorithms, and the solid-liquid contact angle was measured. A calculation model for capillary force was established based on the annular approximation assumption; the specific calculation formula is as follows: ; In the formula, Capillary force at the contact interface, measured in Newtons; It is a constant; Pi is a constant value, approximately 3.14. The equivalent radius of a representative soil particle is expressed in meters. Based on the actual particle size distribution of paddy field soil, its value is usually between 0.0001 meters and 0.002 meters. The surface tension of the liquid phase is expressed in Newtons per meter. Under standard atmospheric pressure and room temperature, the surface tension of deionized water is approximately 0.072 Newtons per meter. The solid-liquid contact angle is measured by the image and is expressed in radians. For the interface between hydrophilic paddy soil and metal substrate, the value is usually between 0 radians and 1.57 radians. This is the cosine value of the solid-liquid contact angle measured from the image.

[0034] Dynamic separation parameters were extracted and viscous force was calculated. Viscosity is the dynamic resistance generated by viscous shearing of the fluid inside the liquid bridge when particles undergo relative motion, and it is directly related to the separation velocity. This dynamic viscous resistance increases during high-speed tumbling or cutting of mud in agricultural machinery components, and is a major microscopic factor contributing to increased mechanical power consumption. A micro-nano manipulation probe was controlled to vertically pull representative soil particles upwards at a preset relative separation velocity until the liquid bridge broke. A charge-coupled device (CCD) microscopic observation system synchronously recorded the pulling process at a set high frame rate, capturing the image of the frame preceding the moment the liquid bridge broke. A viscous force calculation model for the dynamic separation process was established based on the Reynolds equation, and the specific calculation formula is as follows: ; In the formula, The viscous force at the contact interface is expressed in Newtons. The dynamic viscosity of a fluid is expressed in Pascal-seconds. The dynamic viscosity of deionized water at room temperature is approximately 0.001 Pascal-seconds. It is a constant; Pi is a constant value, approximately 3.14. The normal relative separation velocity between particles is measured in meters per second. This velocity is calculated and set based on the mechanical kinematic characteristics of the side-deep fertilization device during field operation, and its range is set from 0.1 meters per second to 1.5 meters per second. The dynamic half-gap thickness, in meters, is extracted from the captured pre-fracture image. It is twice the thickness of the dynamic half-gap; The equivalent radius of a representative soil particle is given in meters. The equivalent radius of the representative soil particle is the square of the equivalent radius. The geometric definition of the dynamic half-gap thickness is: half the total vertical gap distance between the lowest point of the representative soil particle and the underlying metal substrate or the surface of the lower particle layer at the instant the liquid bridge breaks. Based on experimental observations, this thickness typically ranges from 1.0 × 10⁻⁶. -6 Meters up to 5.0 × 10 -4 Between meters.

[0035] Calculate the total liquid bridging force at the contact interface. Linearly superimpose the calculated static capillary force and dynamic viscous force to obtain the maximum total resistance force of the liquid bridge during the tensile process until fracture. The formula for calculating the total liquid bridging force is as follows: ; In the formula, The total hydraulic bridging force at the contact interface is expressed in Newtons. Capillary force at the contact interface, measured in Newtons; The viscous force at the contact interface is expressed in Newtons. The calculated total liquid bridging force quantifies the adhesion strength at the solid-liquid interface at the microscopic level. This calculation result serves as a baseline input for setting the aggregation characteristics of non-Newtonian fluids in discrete element simulation software, and is used to calibrate the adhesion parameters in the macroscopic model.

[0036] The process of contact model correction and calibration based on JKR theory in the establishment of mesoscopic liquid bridge force model and discrete element parameter calibration includes the following:

[0037] A particle contact mechanics model was configured using Discrete Element Method (DEM) software. Traditional Hertzian contact models are only applicable to dry granular materials and do not consider surface attraction at the solid-liquid interface, thus failing to characterize the high cohesion and overall aggregation of saturated paddy field soil particles under water encapsulation. To realistically reflect the non-Newtonian fluid characteristics of mud due to microscopic liquid bridging forces, a cohesive contact model incorporating surface energy properties—specifically, a Hertz-Mundell-Fleming contact model with JKR theory corrections—was configured in the DEM software as the basis for the mechanical calculations of mud particles. The principle of this cohesive contact model is based on conventional Hertzian elastic collision mechanics, introducing an additional normal adhesion force term driven by surface energy. This force requires particles to overcome to separate. The surface energy parameter is introduced to calculate the normal contact force and normal overlap between particles. The formula for calculating the normal contact force is as follows: ; The formula for calculating the normal overlap is as follows: ; In the formula, This represents the normal contact force between particles, measured in Newtons. The equivalent elastic modulus of the contacting particles is expressed in Pascals. Its value is obtained by combining the basic elastic modulus and Poisson's ratio of the two contacting particles or the particle and the geometry. For saturated paddy field mud particles, the basic elastic modulus is usually in the range of 1 MPa to 10 MPa. The equivalent contact radius of the contacting particles, in meters, is calculated from the actual radii of the two particles in contact. The radius of the contact surface after the particles are compressed is expressed in meters, and its value changes dynamically with the degree of compression between the particles. It is the square of the radius of the contact surface of the particle after it is compressed; It is the cube of the radius of the contact surface of the particles after they are compressed by force; Surface energy per unit contact area, expressed in joules per square meter; This represents the normal overlap between particles, expressed in meters. Pi is a constant value, approximately 3.14.

[0038] Calculate the initial surface energy parameters of the contact model. The underlying mechanical calculation engine of the Discrete Element Method (DEM) software cannot directly input the dynamic viscous force and static capillary force of the microscopic liquid bridges; therefore, the obtained total liquid bridge force needs to be converted into surface energy parameters recognizable by the adhesive contact model. According to the critical separation condition of JKR contact theory, the maximum normal pull-out force required to overcome when two particles physically separate is the measured total liquid bridge force. There is a direct mapping relationship between the maximum normal pull-out force and surface energy; the formula for calculating surface energy is as follows: ; In the formula, The surface energy per unit contact area, expressed in joules per square meter. The total liquid bridging force at the contact interface, expressed in Newtons. The equivalent contact radius of the contacting particles, in meters, is calculated from the actual radii of the two contacting particles. This formula establishes a physical relationship between the forces acting on the microscopic liquid bridge and the macroscopic discrete element input parameters. The calculated surface energy reflects the adhesion effect caused by moisture within the mud, and this surface energy is used as the initial benchmark setting for subsequent virtual calibration experiments.

[0039] Virtual calibration experiments were performed to calibrate the parameters of the adhesive contact model. In discrete element simulations, to improve system computational efficiency, it is usually necessary to scale up the actual micro-soil particles geometrically, typically by 2 to 5 times. Due to the enlarged particle size, the directly substituted initial surface energy parameters cannot fully reflect the macroscopic flow pattern of the mud, requiring parameter calibration through physical experiments and simulation comparisons. A cylindrical slump simulation model corresponding to the actual physical measurement environment was established. In the physical experiment, saturated soil samples with a preset moisture content were filled into a standard cylindrical graduated cylinder. The inner diameter and height of the standard cylindrical graduated cylinder ranged from 30 mm to 50 mm. After vertically lifting the graduated cylinder, the actual bottom diffusion diameter and actual center height of the mud after slumping and expanding under gravity were measured. In the discrete element software, a completely corresponding virtual graduated cylinder and particle bed were generated. The specific value ranges for the particle collision recovery coefficient, static friction coefficient, and rolling friction coefficient were set. For high-moisture-content mud, water film lubrication exists between particles. Therefore, the collision recovery coefficient is set to a range of 0.1 to 0.4, the static friction coefficient to a range of 0.1 to 0.5, and the rolling friction coefficient to a range of 0.01 to 0.2. Starting from the calculated initial surface energy parameters, the surface energy parameters and various friction coefficients are adjusted within the set ranges, and multiple slump simulations are performed. The virtual slump diffusion diameter and virtual slump center height are extracted after each simulation convergence.

[0040] A discrete element model of saturated paddy field soil was constructed. The relative error percentage between the actual collapse morphology parameters measured by physical experiments and the virtual collapse morphology parameters obtained by simulation calculations was calculated. The parameter adjustment iteration was stopped when the relative errors of the actual bottom diffusion diameter and the virtual collapse diffusion diameter, as well as the relative errors of the actual center height and the virtual collapse center height, were both less than a preset calibration error threshold. The calibration error threshold is used to measure the accuracy of the virtual model in reflecting real physical properties, and its value range is set to 3% to 5%. A set of surface energy parameters and friction coefficients that meet this calibration error threshold were used as the final parameters and configured and input into the global material library of the discrete element software to complete the reconstruction of highly adhesive saturated paddy field soil particles in the virtual environment. For the mesh generation and time step settings of the particle generation plant in the discrete element software, those skilled in the art can perform conventional configurations based on the minimum particle size and Rayleigh time step criterion. Mesh mapping and solver step size control are well-known techniques in the field and will not be elaborated upon here.

[0041] The process of setting the critical force condition for particle desorption in gas-solid two-phase coupled simulation includes the following: A three-dimensional geometric model of the side-deep fertilization device was extracted, and non-critical chamfers and bolt hole features that do not affect the macroscopic flow field were removed. This model was then imported into fluid dynamics software to construct a continuous phase computational domain. The physical principle of pneumatic removal is to utilize the high kinetic energy of a high-pressure airflow to penetrate the interface between the mud and agricultural machinery components. Through fluid drag and shear stress, the liquid bridges and surface adsorption states between particles are disrupted, forcing the mud to peel off from the component surface. Boundary conditions for the airflow inlet, pressure outlet, and soil-contacting component walls were set for the computational domain. To ensure the convergence and physical realism of the simulation calculations, the soil-contacting component walls were set as non-slip fixed boundaries. For mesh generation and turbulence model selection in the fluid dynamics software, those skilled in the art can use conventional hexahedral mesh generation methods and standard k-epsilon models. The flow field initialization configuration and solver iterative control are well-known techniques in this field and will not be elaborated upon here.

[0042] The critical desorption force equilibrium equation for mud particles was established. The desorption behavior of mud particles on the surface of the side-deep fertilization device is controlled by both the fluid domain and the discrete domain force fields. Based on Newton's second law, the translational control equation for a single mud particle in the flow field was established. The mechanical components of the force equilibrium equation include gravity, total liquid bridge force, airflow drag force, and interparticle adhesive contact force. The calculation formula for the translational control equation is as follows: ; In the formula, The mass of the mud particles is expressed in kilograms and is calculated from the volume and density of the mud particles. The density of saturated mud typically ranges from 1 kg per cubic meter to 1800 kg per cubic meter. This is the velocity vector of the mud particles, in meters per second. The derivative of the velocity vector of the mud particles; The derivative of time; Time, in seconds; The drag force exerted by the continuous phase gas on the mud particles is expressed in Newtons. This is the gravitational acceleration vector, with units of meters per square second. Its scalar value is taken as 9.81 at normal temperature and pressure. The capillary force at the contact interface is expressed in Newtons. This adhesive contact force is calculated from the JKR surface energy parameters obtained from the mesoscopic calibration. The total liquid bridging force at the contact interface is expressed in Newtons. In the macroscopic discrete element contact model, the water film that is partially free on the particle surface and has strong fluidity cannot be fully characterized by JKR theory. Therefore, the total liquid bridging force is introduced into the equation as an independent compensation resistance term to truly reflect the additional adhesion effect of the mud under extremely saturated conditions.

[0043] Airflow drag is the core driving force propelling mud particles away from the surface of the contact component. Its magnitude depends on the relative velocity distribution of the gas and solid phases and the hydrodynamic shape of the particles. The formula for calculating airflow drag is as follows: ; In the formula, The drag force exerted by the continuous phase gas on the mud particles is expressed in Newtons. It is a constant; The drag coefficient is dimensionless and its value is determined based on the relative Reynolds number of mud particles in the flow field. In irregular multiphase flow, its value is usually between 0.4 and 0.5. This is the density of the continuous phase gas, expressed in kilograms per cubic meter, and is approximately 1.225 kilograms per cubic meter under standard atmospheric pressure. This is the velocity vector of the mud particles, in meters per second. It represents the absolute value of the relative velocity between the velocity vector of the continuous phase airflow and the velocity vector of the mud particles; The windward cross-sectional area of ​​the mud particles, in square meters, is calculated based on the equivalent contact radius set in the discrete element model. This is the continuous phase airflow velocity vector, in meters per second. Since the volume of mud particles is much smaller than the size of the fluid computation grid, this airflow velocity vector is usually obtained by interpolation calculation using the velocity of the fluid grid node closest to the center of the mud particles.

[0044] Based on the translational control equations, the force components in the normal direction are extracted, and the critical desorption condition for mud particles to detach from the wall of the side-deep fertilization device is established. When the outward peeling force on the mud particles in the normal direction is greater than the inward adsorption resistance, the physical state of the particles changes from static adhesion to dynamic detachment. The inequality of the critical desorption force balance equation is expressed as follows: ; In the formula, This represents the component of the airflow drag force in the outward direction perpendicular to the wall surface, measured in Newtons. The mass of the mud particles is expressed in kilograms. This represents the component of gravitational acceleration in the outward direction perpendicular to the wall, expressed in meters per second squared. The component of the adhesive contact force in the inward direction perpendicular to the wall surface, expressed in Newtons; This represents the component of the total hydraulic bridge force in the inward direction perpendicular to the wall normal, expressed in Newtons.

[0045] A bidirectional coupling interface was established between fluid dynamics software and discrete element method (DEM) software. In the DEM software, the position coordinates, velocity vectors, and adhesive contact forces between particles were calculated, and this data was transmitted to the fluid dynamics software as momentum source terms. Based on the received particle state data, the fluid dynamics software updated the velocity and pressure field distributions within the fluid computational domain, calculated the airflow drag force, and fed it back to the DEM software. Through this data interaction mechanism, the normal force state of each mud particle in the flow field was monitored in the coupling interface, and the values ​​were substituted into the critical desorption force equilibrium equation. The system determined whether each target particle met the desorption conditions. Particles that met the desorption conditions updated their trajectory matrix, thereby simulating the dynamic peeling evolution process of mud adhering to the surface of a side-deep fertilization device under the action of high-pressure airflow.

[0046] The bidirectional coupling mechanism and simulation process of gas-solid two-phase flow include the following.

[0047] A bidirectional coupled computational domain based on the Eulerian-Lagrange framework is established. During pneumatic removal, the high-pressure airflow cuts into the mud layer, causing the mud particles to detach. The high-density mud particle group that detaches and is suspended in the airflow also occupies physical space and exerts a squeezing and obstructing effect on the local airflow field. To realistically reproduce this physical process of fluid-discrete particle interaction, the time step and spatial mapping relationship between the fluid dynamics software and the discrete element method (DEM) software are configured. The computation time step of the DEM software is set. To ensure the stability of particle collision mechanics calculations, the DEM time step is typically taken as 10% to 20% of the particle Rayleigh time step, specifically set within the range of 1×10⁻⁶. -6 5 x 10 seconds -6Seconds. The computation time step of the fluid dynamics software is set to satisfy the convergence criterion of a Courant number less than 1, typically set to 10 to 1 / 2 times the discrete element time step. Regarding spatial mapping, to ensure smooth momentum exchange, the characteristic size of the fluid computation grid is set to 3 to 5 times the equivalent contact radius of the mud particles.

[0048] The local porosity of the fluid domain is calculated, and the fluid governing equations are updated. The dynamic displacement of mud particles within the fluid grid causes changes in the effective gas flow volume within the grid. Since mud particles occupy part of the space within the fluid grid, gas can only flow in the pores between the particles. Therefore, porosity needs to be introduced to correct the standard continuous phase conservation equations by volume reduction, ensuring mass conservation in the two-phase fluid calculation. Within each fluid calculation time step, the system traverses the fluid calculation grid, calculates the total volume of mud particles within the grid, and calculates the local porosity. Based on the local porosity, the mass conservation equations for the continuous phase gas are calculated as follows: ; In the formula, The local porosity within the fluid computation grid is dimensionless. Its physical meaning is the ratio of the fluid volume not occupied by mud particles in a single grid to the total grid volume, and the value ranges from 0 to 1. This is the density of the continuous phase gas, expressed in kilograms per cubic meter, and is approximately 1.225 kilograms per cubic meter under standard atmospheric pressure. Time, in seconds; It is a dimensionless Hamiltonian differential operator used to calculate spatial gradients and divergences. This is the velocity vector of the continuous phase airflow, in meters per second. It is the partial derivative of the product of local porosity and the density of the continuous phase gas; For the partial derivative of time; It is the dot product operator; It is a constant.

[0049] A two-phase momentum exchange source term is constructed. Due to the acceleration of mud particles by the drag force of the airflow, according to Newton's third law, the mud particles exert an equal and opposite reaction force on the surrounding fluid. Introducing this reaction force into the fluid equations allows simulation of the momentum loss and flow field disturbance caused by the mud particles' reaction to the airflow. The reaction forces of all particles within the grid are volume-averaged and introduced as an external momentum source term into the fluid momentum conservation equation. The calculation formula for the momentum exchange source term is as follows: ; In the formula, For two-phase momentum exchange source terms, it represents the total reaction force of the particle group on the fluid per unit volume, with units of Newtons per cubic meter; The physical volume of the computational grid for a single fluid is expressed in cubic meters. The total number of mud particles contained within the grid at the current moment is calculated for a single fluid, and is a non-negative integer; For fluid to the first in the grid The airflow drag force exerted by each mud particle, measured in Newtons; This is the index value of the mud particles; To the first in the grid One to the first A summation operator for summing mud particles.

[0050] Substituting the momentum exchange source term into the Navier-Stokes equations with porosity, the momentum conservation equation for the continuous phase gas flow field is updated. The formula for calculating the momentum conservation equation is as follows: ; In the formula, The local porosity within the computational mesh for fluid calculation is dimensionless. This is the density of the continuous phase gas, expressed in kilograms per cubic meter. This is the velocity vector of the continuous phase airflow, in meters per second. Time, in seconds; It is a Hamiltonian differential operator, dimensionless; This represents the local pressure in the fluid field, measured in Pascals. This is the fluid viscous stress tensor, measured in Pascals, calculated from the dynamic viscosity and velocity gradient of the gas. The dynamic viscosity of air at normal temperature and pressure is approximately 1.789 × 10⁻⁶. -5 Pascal second; This is the gravitational acceleration vector, with units of meters per square second. Its scalar value is taken as 9.81 at normal temperature and pressure. This is the two-phase momentum exchange source term, with units of Newtons per cubic meter; It is the partial derivative of the product of local porosity, the density of the continuous phase gas, and the velocity vector of the continuous phase gas flow. It is the partial derivative of time.

[0051] Perform a coupled pneumatic removal simulation and output removal efficiency evaluation indicators. Set the airflow nozzle's incident velocity, incident angle, and purging duration, and start the two-way coupled solver. The airflow incident velocity is typically set to 50 m / s to 150 m / s, the incident angle to 30° to 60° with the soil contact component wall, and the purging duration to 0.1 s to 2 s. During the purging duration, the system alternately updates the displacement field of the discrete element particles with the pressure and velocity fields of the fluid domain in real time. After reaching the set purging duration, stop the calculation and extract the mud state data remaining on the surface of the side-deep fertilization device to calculate the pneumatic removal efficiency. The formula for calculating the pneumatic removal efficiency is as follows: ; In the formula, The efficiency of pneumatic cleaning is expressed as a percentage. The total mass of mud adhering to the surface of the side-deep fertilization device at the initial moment before airflow purging, in kilograms; The total mass of mud remaining on the surface of the contact component after the airflow purging time has been set, due to the failure to reach the critical desorption condition, is expressed in kilograms. The output pneumatic removal efficiency serves as a quantitative evaluation index for assessing the pneumatic removal effectiveness of different purging parameter combinations and nozzle structure design rationality. Regarding the data pointer transfer and parallel computing communication protocol between fluid dynamics software and discrete element software at the low-level memory level, those skilled in the art can use a conventional message passing interface architecture. Data synchronization between nodes and solver load balancing are well-known technologies in the field and will not be elaborated upon here.

[0052] The mapping and establishment of the regulation variable and evaluation index system in the process of constructing a parameter-integrated decision-making model for multi-factor and multi-objective optimization includes the following:

[0053] Identify the core control variables for the pneumatic desorption system. The pneumatic desorption process is influenced by multiple physical parameters, necessitating the extraction of control variables that significantly impact flow field distribution and mud desorption. A three-dimensional control variable vector is established, consisting of the nozzle inlet airflow velocity, airflow injection angle, and nozzle-to-ground normal distance. Its vector expression is as follows: .in, The inlet airflow velocity is measured in meters per second. This parameter directly determines the initial kinetic energy input to the system, and its value range is set from 50 meters per second to 150 meters per second. The airflow injection angle is the angle between the central axis of the airflow and the wall of the soil-contacting component of the side-deep fertilization device, measured in degrees. This parameter affects the distribution ratio of the shear force and normal impact force of the airflow on the wall, and its value range is set from 30 degrees to 60 degrees. The normal distance between the nozzle and the wall of the contacting component, in millimeters, determines the degree of energy attenuation of the airflow along its path before reaching the mud surface. Its value range is set from 10 mm to 50 mm.

[0054] A performance evaluation index system for pneumatic mud removal is constructed. In actual agricultural machinery operation, it is necessary not only to ensure mud removal but also to control equipment energy consumption to avoid increasing the additional power burden on the tractor. Therefore, a dual-objective evaluation system is established, including pneumatic mud removal efficiency and system pneumatic power consumption. Pneumatic mud removal efficiency... The degree of mud stripping is used to characterize its thoroughness, and its calculation method has been clearly defined in the aforementioned fluid-structure interaction simulation. The system aerodynamic power consumption is used to evaluate the energy consumption to maintain the set airflow velocity, and its calculation formula is as follows: ; In the formula, The system's aerodynamic power consumption is expressed in watts. This is the density of the continuous phase gas, expressed in kilograms per cubic meter, and is approximately 1.225 kilograms per cubic meter under standard atmospheric pressure. This is the effective cross-sectional area of ​​the airflow nozzle outlet, expressed in square meters. This value is determined based on the actual physical geometry of the selected nozzle, and typically ranges from 1.0 × 10⁻⁶. -4 From 5.0 × 10 square meters -4 square meters; The inlet airflow velocity is measured in meters per second (m / s). This parameter directly determines the initial kinetic energy input to the system, and its value is set to range from 50 m / s to 150 m / s. In the evaluation system, pneumatic cleaning efficiency is a benefit-oriented indicator with the goal of maximizing it; system aerodynamic power consumption is a cost-oriented indicator with the goal of minimizing it.

[0055] A mathematical mapping model between control variables and evaluation indicators was established. Since the bidirectional coupled simulation of gas-solid two-phase flow involves a large amount of mesh computation and particle contact retrieval, the computation time is long, making it difficult to perform traversal calculations across the entire parameter space. To improve optimization efficiency, an experimental design method was used to obtain a limited set of discrete sample point data, and a surrogate model was constructed based on this data to replace the time-consuming simulation. Combining the value range of the control variables, a three-factor, three-level Box-Behnken experimental design was performed within the design space, generating a parameter matrix containing 17 experimental combinations (i.e., 12 factorial experiments and 5 centrally repeated experiments). Each control variable in the parameter matrix was input into the aforementioned constructed CFD-DEM pneumatic removal coupled simulation model to perform the solution, extracting the corresponding pneumatic removal efficiency and system aerodynamic power consumption calculation results.

[0056] A second-order polynomial response surface model reflecting the mapping relationship between control variables and various evaluation indicators was constructed by fitting experimental data using a multivariate nonlinear regression algorithm. In fluid dynamics and particle dynamics coupled systems, the influence of input parameters on output results usually exhibits strong nonlinearity and coupled synergistic characteristics. The second-order polynomial model can fully reflect this nonlinear physical response in its quadratic and interaction terms, and its mathematical solution is relatively stable. The calculation formula of the established mathematical mapping model is as follows: ; In the formula, The first output of the system The predicted values ​​of each evaluation indicator The values ​​are positive integers, taking the values ​​1 and 2, which correspond to the pneumatic removal efficiency, respectively. With system aerodynamic power consumption ; Here, the total dimension of the control variables is... The value is 3; and The independent control variables are the nozzle inlet airflow velocity, airflow injection angle and nozzle ground normal distance after normalization. The normalization process linearly maps the original values ​​of the above three physical parameters to the dimensionless mathematical interval [-1,1] to eliminate the influence of dimensional differences on regression analysis. For constant terms; The coefficients of a linear first-order regression term reflect the independent influence of a single variable on the evaluation index. The regression coefficients are quadratic terms, reflecting the nonlinear surface effect of a single variable on the evaluation index. The coefficients of the interaction term reflect the coupled and synergistic effects on the evaluation index when two different control variables change simultaneously. For the least squares solution and related test analysis of the coefficients in the multiple regression analysis, those skilled in the art can use conventional statistical software and the principles of analysis of variance to complete the calculations. The solution process and error test criteria are well-known techniques in this field and will not be elaborated here. By establishing the above mathematical mapping model, the complex fluid-structure interaction physical process is transformed into a rapidly calculable explicit mathematical expression, providing a basic fitness function for subsequent multi-objective optimization algorithms.

[0057] The training process of the nonlinear prediction model based on BP neural network in the construction of a parameter integrated decision model for multi-factor and multi-objective optimization includes the following:

[0058] Extended data samples were collected for neural network training and testing. In the coupled hydrodynamics and discrete particle dynamics system of pneumatic desorption, the influence of input parameters on the output results exhibits high dimensionality and strong nonlinearity. The aforementioned constructed second-order polynomial response surface model has insufficient fitting accuracy under extreme parameter combinations. Backpropagation neural networks, relying on the backpropagation mechanism of errors, can approximate complex nonlinear function mappings with arbitrary precision by continuously adjusting the connection weights and biases between nodes within the network. To improve prediction accuracy, a backpropagation neural network was introduced to construct a high-precision nonlinear surrogate model. Using the Latin hypercube sampling method, spatially uniform parameter sample points were generated within a set range of nozzle inlet airflow velocity, airflow injection angle, and nozzle distance from the ground. To ensure that the neural network has good feature extraction capabilities and avoid overfitting, 1 to 500 sets of parameter sample points are typically generated. These parameter sample points are divided into training, validation, and test sets, with the data volume ratio typically set to 70%, 15%, and 15%, respectively. These parameter sample points are input into the CFD-DEM gas-solid two-phase flow coupled simulation model for solution, and the corresponding pneumatic removal efficiency and system aerodynamic power consumption calculation results are extracted to construct a dataset for network training and validation.

[0059] Configure the topology of the backpropagation neural network. Construct a multi-layer feedforward neural network containing an input layer, hidden layers, and an output layer. The number of nodes in the input layer is set to 3, corresponding to the three core control variables; the number of nodes in the output layer is set to 2, corresponding to the two evaluation indicators of pneumatic removal efficiency and system aerodynamic power consumption. The number of nodes in the hidden layer is determined by an empirical formula, the specific calculation formula is as follows: ; In the formula, The number of hidden layer nodes is a positive integer. This represents the number of nodes in the input layer, with a value of 3. This represents the number of output layer nodes, with a value of 2. A network structure adjustment constant is used to prevent overfitting or underfitting; its value is set to a range of 1 to 10. Based on this empirical formula, the number of hidden layer nodes... The value range is usually set to an integer between 4 and 12.

[0060] A forward propagation computation model for the neural network is established. The input layer data is weighted, summed, and mapped using an activation function before being passed to the hidden layers. The output calculation formula for the hidden layer nodes is as follows: ; In the formula, For the hidden layer The output value of each node is dimensionless. The independent control variables are the nozzle inlet airflow velocity, airflow injection angle and nozzle ground normal distance after normalization. The normalization process linearly maps the original values ​​of the above three physical parameters to the dimensionless mathematical interval [-1,1] to eliminate the influence of dimensional differences on regression analysis. For the input layer The mud particles reached the first layer of the hidden layer. The connection weights of each node are dimensionless. For the hidden layer The bias threshold of each node, dimensionless; The activation function is a nonlinear function. To handle the nonlinear mapping of fluid-structure interaction systems, the Sigmoid function is used, and its calculation formula is as follows: ; In the formula, The input argument of the activation function is dimensionless; The base of the natural logarithm is 2.718, and the constant value is approximately 2.718.

[0061] The data from the hidden layer is linearly weighted and then passed to the output layer to generate the final prediction result. The calculation formula for the output layer nodes is as follows: ; In the formula, For the output layer The predicted value of each node corresponds to the network calculation result of the pneumatic removal efficiency or the system aerodynamic power consumption. For the hidden layer The node to the output layer The connection weights of each node are dimensionless. For the output layer The bias threshold for each node, dimensionless. The values ​​are positive integers, taking the values ​​1 and 2, which correspond to the pneumatic removal efficiency, respectively. With system aerodynamic power consumption .

[0062] Define the loss function and perform backpropagation updates of the error. To quantify the deviation between the network's predictions and the simulated true values, mean squared error is used as the loss function. The formula for calculating the loss function is as follows: ; In the formula, The total prediction error of the network in the current forward propagation batch is dimensionless. The first obtained through CFD-DEM simulation calculation The true value of each evaluation indicator.

[0063] When the total prediction error of the network When the error exceeds a set error threshold, the system triggers the backpropagation mechanism. The partial derivatives of the error with respect to each weight and bias are calculated using gradient descent, and the network parameters are updated layer by layer along the reverse direction of the gradient. The formula for weight update is as follows: ; In the formula, The connection weights will be updated in the next iteration. For the current number The connection weights for each iteration; The iteration number is a positive integer. The network learning rate controls the step size of parameter updates, with a value range of 0.01 to 0.1. For error Regarding connection weight The partial derivatives of .

[0064] The network is trained iteratively until the convergence condition is met. The network continuously updates its parameters through forward and backward propagation on the training dataset. The set error threshold is typically set to a value between 1.0 × 10⁻⁶. -4 Up to 1.0×10 -3 When the total prediction error of the network is below the error threshold for several consecutive iterations, or when the total number of iterations reaches the set maximum number of training rounds (ranging from 0 to 5000), training stops and the current weight matrix and bias vector are saved. The trained BP neural network model will replace the time-consuming CFD-DEM coupled simulation as the fitness evaluation engine for subsequent multi-objective optimization algorithms. For the specific mathematical derivation of the chain rule for the backpropagation partial derivative matrix and the dynamic decay strategy of the learning rate in the neural network, those skilled in the art can implement it using the underlying logic of conventional deep learning frameworks. The analysis of partial derivatives and the computational graph differentiation are well-known techniques in this field and will not be elaborated upon here.

[0065] The response surface analysis and multi-objective weighted optimization decision-making process in the construction of a parametric integrated decision model for multi-factor and multi-objective optimization specifically include the following:

[0066] Response surface analysis was performed to investigate the interaction of control variables. In pneumatic cleaning systems, the effects of nozzle inlet airflow velocity, airflow injection angle, and nozzle-to-ground normal distance on cleaning performance are not entirely independent; physical coupling effects exist between the variables. For example, as the nozzle-to-ground normal distance increases, the kinetic energy attenuation along the airflow intensifies. In this case, a nonlinear increase in the nozzle inlet airflow velocity is required to maintain the shear stress on the target wall. Based on the previously trained BP neural network nonlinear prediction model, high-density gridded discrete data points were generated within the parameter design space. Keeping any one of the three control variables fixed at the central level, the pneumatic cleaning efficiency and system aerodynamic power consumption were calculated when the other two variables changed continuously throughout their range. The calculation results were extracted, and a three-dimensional spatial response surface plot and contour projection plot were plotted. By analyzing the gradient trend of the surface and the density of the contour lines, parameter-sensitive regions that are beneficial for improving pneumatic cleaning efficiency without increasing system aerodynamic power consumption were identified.

[0067] A normalized multi-objective comprehensive evaluation function is constructed. The optimization direction for pneumatic cleaning efficiency is maximization, while the optimization direction for system aerodynamic power consumption is minimization; these two are mutually constrained in terms of physical mechanisms. Higher high-pressure airflow velocity inevitably leads to better mud stripping effect, but at the same time, it increases the energy consumption of the tractor air compressor unit. To find the optimal balance between conflicting objectives, it is necessary to eliminate the differences in dimensions and numerical magnitudes of different evaluation indicators. The range transformation method is used to perform dimensionless normalization on the predicted pneumatic cleaning efficiency and system aerodynamic power consumption, and a single comprehensive evaluation function is constructed by combining it with weighting coefficients. The calculation formula for the comprehensive evaluation function is as follows: ; In the formula, The fitness value of the multi-objective comprehensive evaluation function is dimensionless; the larger the value, the better the overall system performance. The efficiency of pneumatic cleaning is expressed as a percentage. and These are the theoretical maximum and minimum values ​​of pneumatic removal efficiency within the current parameter design space, respectively. These values ​​are obtained by extracting the boundary of the sample dataset when constructing the surrogate model or by pre-calculating through single-objective optimization. The system's aerodynamic power consumption is expressed in watts. and These are the theoretical maximum and minimum values ​​of the system's aerodynamic power consumption within the current parameter design space, respectively, which are obtained through the boundary of the statistical sample dataset or single-objective pre-optimization. The target weighting coefficient for the efficiency of pneumatic cleaning is dimensionless. Let be the target weighting coefficient for the system's aerodynamic power consumption, which is dimensionless. The weighting coefficient satisfies... The constraints are as follows. In actual agricultural machinery operation, ensuring the cleanliness of the surfaces of agricultural machinery parts is the primary engineering task. Energy saving and consumption reduction are then considered on this basis. Therefore, the efficiency of pneumatic cleaning is set as the dominant factor, with a weighting coefficient of [missing information]. The value range is usually set to 0.6 to 0.8, corresponding to the system aerodynamic power consumption weighting coefficient. The value range is set to 0.2 to 0.4.

[0068] Global parameter optimization is performed using the particle swarm optimization algorithm. In a continuous three-dimensional space of control variables, exhaustively searching for the optimal solution through discrete grid points is inefficient and prone to missing theoretical extreme points. Particle swarm optimization is a heuristic search algorithm based on swarm intelligence. By simulating information sharing and collaboration mechanisms within a swarm, it can effectively search for the global optimum in complex nonlinear parameter spaces without relying on the gradient information of the objective function. Using particle swarm optimization, each parameter combination within the range of control variable values ​​is defined as a candidate solution in the search space, acting as a particle, and a multi-objective comprehensive evaluation function is used. As a fitness evaluation criterion for particles, the optimal combination of control variables that maximizes the fitness value is searched globally. The particle swarm size is initialized to 50 to 500 particles, and the maximum number of iterations is set to 500 to 500. During each iteration of optimization, particles update their search velocity and spatial position in real time based on their individual historical best position and the swarm's global best position. The formulas for calculating the particle velocity and position updates are as follows: ; ; In the formula, For the first The mud particle at the first The update speed of the next iteration in the dimensional search space, where is the mud particle index, is a positive integer and ranges from 1 to the set particle swarm size; The dimension index is a positive integer. Since there are three control variables, therefore... The value can be 1, 2, or 3; For the first The mud particle at the first The current dimensional search space The speed of each iteration; For the first The mud particle at the first The update position of the next iteration in the dimensional search space corresponds to a new set of control variable parameters; For the first The current location of each mud particle; This is the inertial weight, used to balance the algorithm's global exploration and local exploitation capabilities; its value is typically set between 0.4 and 0.9. This is an individual cognitive learning factor that guides particles toward their historical best position, with a value range of 1.5 to 2.0. This is a social experience learning factor that guides particles to approach the global optimal position of the population; its value is set between 1.5 and 2.0. and These are random numbers uniformly distributed between 0 and 1, used to increase the randomness of the search process; For the first The mud particles have so far reached [number] in [number] [units]. The individual's best historical position experienced by the dimensional; For the entire particle swarm to date in the 1st The globally optimal position that the group has experienced; Time is measured in seconds.

[0069] Output the integrated decision parameters and map them to the physical control terminal. The particle swarm optimization algorithm iterates continuously until the maximum number of iterations is reached, or the fitness increment falls below the set convergence accuracy requirement for several consecutive times. After stopping the calculation, extract the global optimum position of the population at this point. The corresponding control variable values. Considering the physical control precision limits of actual agricultural machinery pneumatic actuators and mechanical adjustment structures, the obtained theoretical optimal solution is rounded down for engineering purposes. The nozzle inlet airflow velocity is rounded to one decimal place, and the airflow injection angle and nozzle ground normal distance are rounded to integers. The processed optimal nozzle inlet airflow velocity, optimal airflow injection angle, and optimal nozzle ground normal distance combination is used as the final decision scheme and sent to the onboard programmable logic controller of the side-deep fertilization device. The controller adjusts the output frequency of the air compressor and the stroke of the nozzle pose adjustment motor according to the decision scheme, so that the actual pneumatic cleaning physical system operates at the optimal working point that can efficiently remove mud while saving power consumption. For the particle swarm algorithm boundary overflow handling mechanism and the controller analog-to-digital conversion protocol, those skilled in the art can use conventional boundary absorption strategies and standardized industrial communication protocols. The compilation of its underlying hardware driver instructions is a well-known technology in the field and will not be described in detail here.

Claims

1. A method for determining pneumatic cleaning parameters in a side-deep fertilization device based on multi-objective optimization, characterized in that, Includes the following steps: A saturated soil sample consisting of mud particles and a continuous water film was prepared. The test probe was controlled to contact and separate from the saturated soil sample. The maximum normal tensile force and the maximum tangential sliding force were extracted, and the normal adhesion specific pressure and the tangential adhesion specific pressure were calculated. The total liquid bridging force at the mud particle contact interface is calculated. Combining the normal adhesion specific pressure and the tangential adhesion specific pressure, an adhesive contact model is configured in the discrete element software, and the surface energy parameters are calibrated through a virtual calibration test to construct a discrete element model of saturated paddy field soil. Based on the discrete element model of the saturated paddy soil, a critical desorption force balance equation for mud particles, including gravity, the total liquid bridge force, airflow drag force and interparticle adhesive contact force, is established. The cleaning process is simulated by combining the gas-solid bidirectional coupling control equation, and the pneumatic cleaning efficiency and system aerodynamic power consumption are extracted. The nozzle inlet airflow velocity, airflow injection angle, and nozzle ground normal distance are used as control variables. A nonlinear prediction model is established by training a backpropagation neural network in combination with the control variables, the pneumatic cleaning efficiency, and the system aerodynamic power consumption. A multi-objective comprehensive evaluation function is constructed, and the extreme value is solved by the particle swarm optimization algorithm. The optimal combination of control variables is output as the airflow cleaning control parameters of the side-deep fertilization device to control the side-deep fertilization device to perform the cleaning operation.

2. The method for determining pneumatic cleaning parameters of a side-deep fertilization device based on multi-objective optimization according to claim 1, characterized in that, The process of preparing a saturated soil sample composed of mud particles and having a continuous water film, controlling the test probe to contact and separate from the saturated soil sample, extracting the maximum normal tensile force and the maximum tangential sliding force, and calculating the normal adhesion specific pressure and the tangential adhesion specific pressure specifically includes: According to the set moisture content, the dry soil sample was mixed with deionized water to form a slurry and then left to stand in a cylindrical standard sample box. After the mud has settled, a calculated volume of deionized water is applied to the surface of the mud using a micro-liquid phase injection device to construct a continuous water film and obtain the saturated soil sample. The test probe is controlled to descend vertically at a set loading speed to apply normal pressure to the saturated soil sample. After maintaining a set contact time, it is vertically ascended at a set unloading speed. The peak value on the mechanical change curve during the normal tensile process is recorded by a force sensor as the maximum normal tensile force. The test probe is controlled to move at a set sliding speed along the horizontal direction on the surface of the saturated soil sample while maintaining the dynamic contact state of the normal pressure force. The maximum resistance value before the relative sliding change of the tangential resistance is recorded as the maximum tangential sliding force. Obtain the contact area between the test probe and the saturated soil sample; The normal adhesion specific pressure is calculated based on the ratio of the maximum normal tensile force to the contact area, and the tangential adhesion specific pressure is calculated based on the ratio of the maximum tangential sliding force to the contact area.

3. The method for determining pneumatic cleaning parameters of a side-deep fertilization device based on multi-objective optimization according to claim 1, characterized in that, The process of calculating the total liquid bridging force at the contact interface of the mud particles specifically includes: The solid-liquid contact angle, dynamic half-gap thickness, and surface tension of the liquid phase at the contact interface of the mud particles are obtained. Based on the solid-liquid contact angle and the surface tension of the liquid phase, the static capillary force at the mud particle contact interface is calculated based on the annular approximation assumption. Based on the dynamic half-gap thickness and relative separation velocity, the dynamic viscous force at the mud particle contact interface is calculated using the Reynolds equation. The static capillary force and the dynamic viscous force are linearly superimposed to obtain the total liquid bridge force at the mud particle contact interface.

4. The method for determining pneumatic cleaning parameters of a side-deep fertilization device based on multi-objective optimization according to claim 1, characterized in that, The process of configuring the adhesive contact model in the discrete element software and calibrating the surface energy parameters through the virtual calibration experiment to construct the discrete element model of the saturated paddy field soil specifically includes: The initial surface energy parameters of the adhesive contact model are calculated based on the total liquid bridging force, the normal adhesion specific pressure, the tangential adhesion specific pressure, and the equivalent contact radius of the mud particles. The actual bottom diffusion diameter and actual center height of the saturated soil sample after it has collapsed and come to rest under the action of gravity are obtained by vertically lifting a standard cylindrical graduated cylinder containing the saturated soil sample. In the discrete element software, a corresponding virtual graduated cylinder and a particle bed containing the mud particles are established. Starting from the initial surface energy parameters, the surface energy parameters are adjusted and multiple slump simulation calculations are performed to extract the virtual slump diffusion diameter and the virtual slump center height. When the relative error between the actual bottom diffusion diameter and the virtual collapse diffusion diameter, and the relative error between the actual center height and the virtual collapse center height are both less than the calibration error threshold, the surface energy parameter adjustment iteration is stopped, and the surface energy parameter that meets the calibration error threshold is input as the final parameter into the global material library of the discrete element software, thereby completing the construction of the discrete element model of the saturated paddy field soil. The calibration error threshold is set based on the relative error statistics of actual physical collapse tests.

5. The method for determining pneumatic cleaning parameters of a side-deep fertilization device based on multi-objective optimization according to claim 1, characterized in that, The specific process for establishing the critical desorption force equilibrium equation for mud particles, which includes the gravity, the total liquid bridging force, the airflow drag force, and the interparticle adhesive contact force, is as follows: Extract the force components of the translational control equation of the mud particles in the normal direction of the contact wall; Based on the force components, the force equilibrium equation for the critical desorption of the mud particles is established according to the component of gravity in the outward direction normal to the contact wall, the component of airflow drag force in the outward direction normal to the contact wall, the component of interparticle adhesive contact force in the inward direction normal to the contact wall, and the component of total liquid bridge force in the inward direction normal to the contact wall. Based on the force balance equation for critical desorption of mud particles, when the sum of the component of the airflow drag force in the outward direction normal to the contact wall and the component of gravity in the outward direction normal to the contact wall is greater than or equal to the sum of the component of the interparticle adhesive contact force in the inward direction normal to the contact wall and the component of the total liquid bridge force in the inward direction normal to the contact wall, the mud particles are determined to meet the critical desorption condition and transition to a dynamic detachment state.

6. The method for determining pneumatic cleaning parameters of a side-deep fertilization device based on multi-objective optimization according to claim 1, characterized in that, The process of simulating the desorption process using the aforementioned gas-solid bidirectional coupling control equations and extracting the pneumatic desorption efficiency and the system's aerodynamic power consumption specifically includes: The gas-solid two-way coupling control equation is executed in the fluid dynamics software. In each calculation time step, the fluid calculation grid of the fluid dynamics software is traversed and the total volume of the mud particles contained in the fluid calculation grid is counted. The local porosity is calculated and the local porosity is introduced to correct the mass conservation equation of the airflow by volume reduction. The reaction forces of all mud particles within the fluid computation grid are volume averaged to construct a two-phase momentum exchange source term. The two-phase momentum exchange source term is substituted into the Navier-Stokes equation with the local porosity to update the momentum conservation equation of the airflow. The flow field of the airflow and the motion state of the mud particles are synchronized through bidirectional data interaction between the fluid dynamics software and the discrete element software to complete the simulation of the cleaning process. The pneumatic cleaning efficiency and the aerodynamic power consumption of the system are extracted based on the simulation results.

7. The method for determining pneumatic cleaning parameters of a side-deep fertilization device based on multi-objective optimization according to claim 1, characterized in that, The process of training the backpropagation neural network to establish the nonlinear prediction model specifically includes: A parameter matrix containing the control variables is generated using experimental design methods. This matrix is ​​then input into a CFD-DEM gas-solid bidirectional coupling simulation model constructed using the discrete element method and fluid dynamics software for solving. A second-order polynomial response surface model is then established based on the output data. The Latin hypercube sampling method is used to generate uniform parameter sample points within the range of the control variable, which are then divided into training set, validation set and test set. The backpropagation neural network with multi-layer feedforward is constructed, with the control variable as the input layer, the pneumatic removal efficiency and the system aerodynamic power consumption as the output layer, and the mean square error as the loss function. After the backpropagation neural network has passed forward, if the total prediction error is greater than the set target training error threshold, the backpropagation mechanism is triggered to update the parameters of the backpropagation neural network layer by layer using the gradient descent method until the convergence condition is met, and the trained nonlinear prediction model is output.

8. The method for determining pneumatic cleaning parameters of a side-deep fertilization device based on multi-objective optimization according to claim 1, characterized in that, The process of constructing the multi-objective comprehensive evaluation function specifically includes: The optimization direction for the pneumatic cleaning efficiency is set to maximize, and the optimization direction for the system's pneumatic power consumption is set to minimize. The range transformation method is used to perform dimensionless normalization on the pneumatic removal efficiency and the system aerodynamic power consumption under the current control variables obtained through the nonlinear prediction model. Target weight coefficients are assigned to the pneumatic removal efficiency and the system aerodynamic power consumption, respectively. The normalized pneumatic removal efficiency and the system aerodynamic power consumption are multiplied by their respective target weight coefficients and then summed to construct a single multi-objective comprehensive evaluation function.

9. The method for determining pneumatic cleaning parameters of a side-deep fertilization device based on multi-objective optimization according to claim 1, characterized in that, The process of solving for extreme values ​​using the particle swarm optimization algorithm specifically includes: Each parameter combination within the range of values ​​of the control variable is defined as a candidate solution in the search space as a particle, and the multi-objective comprehensive evaluation function is used as the fitness evaluation criterion for the particle. In each iteration of optimization, the search speed and spatial position of each particle are updated in real time by utilizing the individual particle's historical optimal position, the group's global optimal position, the individual cognitive learning factor, and the social experience learning factor. The process continues iterating until the maximum number of iterations is reached or the fitness increment meets the convergence accuracy. The values ​​corresponding to the global optimal position of the population are extracted as the theoretical optimal solution to complete the extreme value solution of the multi-objective comprehensive evaluation function.

10. The method for determining pneumatic cleaning parameters of a side-deep fertilization device based on multi-objective optimization according to claim 9, characterized in that, The process of outputting the optimal combination of control variables as the airflow clearing control parameters of the side-deep fertilization device specifically includes: The extracted theoretical optimal solution is rounded down in an engineering manner to obtain the optimal values ​​of the nozzle inlet airflow velocity, the airflow injection angle, and the nozzle normal distance from the ground, thereby generating the optimal combination of control variables; The optimal combination of control variables is sent as the airflow cleaning control parameters to the onboard programmable logic controller of the side-deep fertilization device, wherein the side-deep fertilization device includes the onboard programmable logic controller, an air compressor, and a nozzle position adjustment motor. The onboard programmable logic controller outputs control commands based on the airflow cleaning control parameters, directly adjusting the output frequency of the air compressor and the stroke of the nozzle position adjustment motor to control the side-deep fertilization device to perform cleaning operations.

Images