Process optimization analysis system for vacuum concentration
Through the combination of multi-source perception, digital twin modeling and adaptive execution modules, the problems of insufficient perception accuracy and poor control system stability in the vacuum concentration process are solved, and efficient and energy-saving concentration process optimization and stable operation are achieved.
Patent Information
- Application Number
- CN202510424939.9
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-04-07
- Publication Date
- 2025-08-08
- Estimated Expiration
- Not applicable · inactive patent
AI Technical Summary
The existing vacuum concentration process has problems such as insufficient perceptual accuracy, simplified process modeling, single optimization strategy and lack of dynamic compensation and anti-interference capabilities of the control system, resulting in poor system operation stability and affecting efficiency and product quality.
The multi-source perception module is used to collect data in real time, the digital twin modeling module is used to perform high-precision modeling, and the multi-objective optimization module is used to comprehensively optimize it. The adaptive execution module is used to compensate the system lag effect and anti-interference, and a dual closed-loop adjustment mechanism is established.
Accurate perception and dynamic modeling of complex concentration processes are achieved, the stability of the concentration process and product quality are improved, the operation efficiency and service life of the equipment are improved, and the production cost is reduced.
Smart Images

Figure CN120449639A_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to the technical field of industrial automation and intelligent manufacturing, and more specifically, to a process optimization analysis system for vacuum concentration. Background Art
[0002] Vacuum concentration is a critical process in the chemical, pharmaceutical, and food industries, used to accelerate solvent evaporation by reducing pressure, thereby concentrating materials. Traditional vacuum concentration processes typically rely on manual experience or simple automated control systems, which often struggle to accurately sense and control complex process parameters in real time. For example, precise monitoring and control of parameters such as the temperature gradient within the evaporation chamber, the dynamic changes in vacuum level, and the rheological properties of the material have long been a challenge for the industry. Furthermore, existing kinetic modeling of the concentration process is often simplistic and difficult to adapt to complex industrial scenarios, resulting in limited optimization results. Regarding optimization strategies, traditional single-objective optimization methods cannot simultaneously address multiple metrics such as energy efficiency, product yield, and equipment lifespan, making them unable to meet the demands of modern industrial production for high efficiency, energy conservation, and sustainable development. Furthermore, existing systems lack effective dynamic compensation and anti-interference mechanisms for process hysteresis and external interference, further limiting their stability and reliability.
[0003] During the implementation of the present invention, the inventors discovered at least the following problems or deficiencies in the prior art: The existing vacuum concentration process suffers from insufficient sensor accuracy, making it difficult to obtain comprehensive and accurate process parameters in real time; the process modeling approach is overly simplified and unsuitable for complex industrial scenarios; the optimization strategy is single, failing to achieve simultaneous optimization of multiple objectives; and the control system lacks dynamic compensation and anti-interference capabilities, resulting in poor system operational stability. These issues have severely impacted the efficiency, product quality, and equipment lifespan of the vacuum concentration process, hindering technological advancement and development in the industry. Summary of the Invention
[0004] The present invention provides a process optimization analysis system for vacuum concentration, comprising the following multi-layer cascade modules:
[0005] The multi-source sensing module uses a distributed sensor array to collect real-time data on the temperature gradient distribution in the evaporation chamber, the dynamic waveform of the vacuum degree, the rheological properties of the material, and the latent heat parameters of the phase change;
[0006] The digital twin modeling module builds a hidden Markov chain model of process parameters and concentration dynamics based on transfer learning, and integrates historical batch data to generate an adaptive covariance matrix;
[0007] The multi-objective optimization module uses a mixed integer dynamic programming algorithm to solve the Pareto front solution set in a non-convex solution space, and simultaneously optimizes energy efficiency, product yield, and equipment life indicators;
[0008] The adaptive execution module compensates for the system hysteresis effect in real time through nonlinear model predictive control and establishes a dual closed-loop anti-interference regulation mechanism.
[0009] Furthermore, the multi-source perception module includes:
[0010] The multimodal data fusion unit uses the improved DS evidence theory to perform confidence-weighted fusion on heterogeneous data from infrared thermal imagers, ultrasonic densitometers, and micro-pressure differential sensors;
[0011] The process fingerprint extraction unit extracts the time-frequency domain feature matrix of the vacuum pulsation signal through wavelet packet decomposition, and constructs a process state fingerprint library including kurtosis factor, energy entropy and Lyapunov exponent.
[0012] Furthermore, the digital twin modeling module includes:
[0013] Variational Autoencoder (VAE), which maps high-dimensional sensor data into a latent space to generate a low-dimensional manifold representation of the process state;
[0014] The physical information neural network (PINN) embeds the Navier-Stokes equation constraints to construct a differential-data hybrid driven model of the evaporation process, and its loss function is:
[0015]
[0016] Where u is the fluid velocity field, T is the temperature field, α is the thermal diffusion coefficient, λ i To balance the weight.
[0017] Furthermore, the multi-objective optimization module includes:
[0018] Hierarchical decision-making architecture: the upper layer uses the NSGA-III algorithm to generate the global Pareto solution set, and the lower layer handles parameter uncertainty through distributed robust optimization (DRO);
[0019] Constraint processing unit, defining the dynamic feasible region:
[0020]
[0021] Where Γ1(t) is the time-varying safety threshold and Γ2(t) is the semi-positive definite matrix constraint.
[0022] Furthermore, the distributed robust optimization adopts Wasserstein fuzzy sets:
[0023]
[0024] Among them, P ∈ ={Q:W1(Q,P N )≤∈}, ξ is an uncertain parameter, and W1 is the 1-Wasserstein distance.
[0025] Furthermore, the adaptive execution module includes:
[0026] Feedforward-feedback composite controller, the feedforward channel uses Volterra series to compensate for nonlinear time lag, and the feedback channel is designed as a fractional-order PID controller:
[0027] u(t)=K p e(t)+K i D -λ e(t)+K d D μ e(t)
[0028] Where λ∈(0,1),μ∈(0,1) are fractional orders, and D is a fractional differential operator;
[0029] The actuator dynamic compensation unit establishes the Hammerstein-Wiener model of the servo motor and the pneumatic valve, and realizes dynamic linearization through inverse model decoupling.
[0030] Furthermore, the Volterra series kernel function is updated online by a recursive least squares algorithm:
[0031]
[0032] where γ k is a variable step learning rate, satisfying
[0033] Furthermore, it also includes:
[0034] The cognitive digital twin module builds a fault propagation model based on the knowledge graph and defines the impact of fault modes:
[0035]
[0036] where w i is the fault link weight, Δθ i is the parameter offset, σ i is the sensitivity coefficient;
[0037] Self-repair unit, when Ψ>Ψ th The reconstruction mechanism is triggered in real time to maintain the key functions of the system through structural adaptive reorganization.
[0038] Furthermore, the reconstruction mechanism includes:
[0039] Control strategy migration: mapping the current working condition to the nearest case library scenario and loading the pre-trained reinforcement learning strategy;
[0040] Resource reallocation, optimizing the sensor-actuator pairing relationship based on the Hungarian algorithm, and dynamically reconstructing it into multiple functionally independent virtual subsystems.
[0041] Furthermore, it also includes:
[0042] Energy and mass flow optimization module, establish Analytical model calculation of each unit Loss coefficient:
[0043]
[0044] Where T0 is the ambient temperature, m i is the mass flow rate, s i is the specific entropy, Q i is the heat flow rate;
[0045] The waste heat intelligent recovery unit decides the optimal waste heat utilization path based on the fuzzy cognitive map, giving priority to meeting the energy demand of the preheating stage.
[0046] According to the above-mentioned embodiments of the present invention, there are at least the following beneficial effects: the process optimization and analysis system for vacuum concentration of the present invention collects the temperature gradient distribution, vacuum dynamic waveform, material rheological characteristics and phase change latent heat parameters of the evaporation chamber in real time through the multi-source perception module, and combines the high-precision modeling capabilities of the digital twin modeling module to achieve accurate perception and dynamic modeling of complex concentration processes. This precise perception and modeling capability provides a solid foundation for optimization control and can effectively improve the stability of the concentration process and product quality. At the same time, the multi-objective optimization module uses a mixed integer dynamic programming algorithm to solve the Pareto front solution set in a non-convex solution space, which can simultaneously optimize the energy efficiency ratio, product yield and equipment life indicators to achieve high efficiency, energy saving and sustainable development of the process.
[0047] Furthermore, the adaptive execution module uses nonlinear model predictive control to compensate for system hysteresis in real time and establishes a dual closed-loop anti-interference regulation mechanism. This improves the system's dynamic response and anti-interference performance, ensuring stable operation of the concentration process under complex operating conditions. This intelligent control strategy not only improves equipment efficiency and service life, but also reduces production costs and enhances the overall competitiveness of the system. BRIEF DESCRIPTION OF THE DRAWINGS
[0048] The above and other objects, features and advantages of the exemplary embodiments of the present invention will become readily apparent by reading the following detailed description with reference to the accompanying drawings, in which several embodiments of the present invention are shown by way of example and not limitation, in which:
[0049] Figure 1 A schematic structural diagram of a process optimization analysis system for vacuum concentration provided by one embodiment of the present invention. DETAILED DESCRIPTION
[0050] The principles and spirit of the present invention will be described below with reference to several exemplary embodiments. It should be understood that these embodiments are provided solely to enable those skilled in the art to better understand and implement the present invention, and are not intended to limit the scope of the present invention in any way. Rather, these embodiments are provided to make the present invention more thorough and complete, and to fully convey the scope of the present invention to those skilled in the art.
[0051] Those skilled in the art will appreciate that the embodiments of the present invention may be implemented as a system, apparatus, device, method, or computer program product. Therefore, the present invention may be implemented in the following forms: entirely in hardware, entirely in software (including firmware, resident software, microcode, etc.), or in a combination of hardware and software.
[0052] It should be noted that any number of elements in the drawings is for illustration only and not for limitation, and any naming is only for distinction and does not have any limiting meaning.
[0053] Reference below Figure 1 , Figure 1 This is a schematic diagram of a process optimization analysis system for vacuum concentration provided by one embodiment of the present invention. Figure 1 As shown, a process optimization analysis system 100 for vacuum concentration includes the following multi-layer cascade modules:
[0054] The multi-source sensing module 101 uses a distributed sensor array to collect the temperature gradient distribution of the evaporation chamber, the dynamic waveform of the vacuum degree, the rheological properties of the material, and the latent heat parameters of the phase change in real time;
[0055] The digital twin modeling module 102 builds a hidden Markov chain model of process parameters and concentration dynamics based on transfer learning, and integrates historical batch data to generate an adaptive covariance matrix;
[0056] The multi-objective optimization module 103 uses a mixed integer dynamic programming algorithm to solve the Pareto front solution set in a non-convex solution space, and simultaneously optimizes energy efficiency, product yield, and equipment life indicators;
[0057] The adaptive execution module 104 compensates for the system hysteresis effect in real time through nonlinear model predictive control and establishes a dual closed-loop anti-interference regulation mechanism.
[0058] It should be noted that a process optimization and analysis system for vacuum concentration includes multiple cascaded modules. The multi-source sensing module is used to collect the temperature gradient distribution of the evaporation chamber, the dynamic waveform of the vacuum degree, the rheological properties of the material, and the latent heat parameters of the phase change in real time. The distributed sensor array is an array composed of multiple sensors of different types and is used to obtain various types of data. The digital twin modeling module constructs a hidden Markov chain model of process parameters and concentration dynamics based on transfer learning. Transfer learning is a technology that allows models to transfer knowledge between different tasks or data sets. The hidden Markov chain model is a model used to describe random processes with hidden states. At the same time, the module integrates historical batch data to generate an adaptive covariance matrix. The multi-objective optimization module uses a mixed integer dynamic programming algorithm to solve the Pareto front solution set in a non-convex solution space. The mixed integer dynamic programming algorithm is an algorithm that combines integer programming and dynamic programming. The Pareto front solution set refers to the set of solutions that cannot improve one objective without compromising other objectives in multi-objective optimization, and simultaneously optimizes energy efficiency, product yield, and equipment life indicators. The adaptive execution module compensates for the system lag effect in real time through nonlinear model predictive control. Nonlinear model predictive control is a method of predicting the future state of the system based on a nonlinear model and performing control. It also establishes a dual closed-loop anti-interference regulation mechanism.
[0059] Specifically, in a distributed sensor array, different sensors are responsible for collecting different data. For example, an infrared thermal imager can collect temperature distribution data in the evaporation chamber, an ultrasonic concentration meter can measure material concentration and thus reflect its rheological properties, and a micro-differential pressure sensor can obtain vacuum data. In the digital twin modeling module, transfer learning can transfer knowledge from existing models to new production process modeling based on the similarities between production data from different batches. The parameters of the hidden Markov chain model must be determined based on the statistical characteristics of the actual process data, such as state transition probabilities and observation probabilities. The adaptive covariance matrix is adjusted based on the fluctuations of historical batch data to better describe the model's uncertainty. In the multi-objective optimization module, the parameters of the mixed integer dynamic programming algorithm must take into account the scale and complexity of the problem, such as the range of decision variables and constraints. The energy efficiency ratio refers to the ratio of the energy effectively utilized to the total energy consumed during the concentration process. The product yield is the ratio of the actual product yield to the theoretical product yield. The equipment life indicator can be measured by factors such as the wear of key components and operating time. In the adaptive execution module, the model parameters of the nonlinear model predictive control can be determined by learning and identifying the system's historical operating data. The dual closed-loop anti-interference adjustment mechanism generally includes an inner loop and an outer loop. The inner loop is used to quickly respond to changes within the system, and the outer loop is used to resist external interference.
[0060] Preferably, during the data acquisition process of the multi-source perception module, in order to ensure the accuracy and stability of the data, the sensors can be calibrated regularly. The calibration period can be set according to the type of sensor and the use environment, generally once a week or a month. When constructing the hidden Markov chain model in the digital twin modeling module, in addition to using historical batch data, expert experience knowledge can also be introduced to optimize and correct the model. For the multi-objective optimization module, when solving the Pareto front solution set, parallel computing technology can be used to improve computing efficiency and reduce solution time. In the adaptive execution module, the inner loop and outer loop of the dual closed-loop anti-interference adjustment mechanism can adopt different control algorithms, such as the inner loop adopts the proportional-integral control algorithm and the outer loop adopts the fuzzy control algorithm to enhance the anti-interference ability of the system.
[0061] In some embodiments, the multi-source perception module includes:
[0062] The multimodal data fusion unit uses the improved DS evidence theory to perform confidence-weighted fusion on heterogeneous data from infrared thermal imagers, ultrasonic densitometers, and micro-pressure differential sensors;
[0063] The process fingerprint extraction unit extracts the time-frequency domain feature matrix of the vacuum pulsation signal through wavelet packet decomposition, and constructs a process state fingerprint library including kurtosis factor, energy entropy and Lyapunov exponent.
[0064] It should be noted that the multi-source perception module includes a multimodal data fusion unit and a process fingerprint extraction unit. The multimodal data fusion unit uses the improved DS evidence theory to weightedly fuse the different types (i.e., heterogeneous) data collected by the infrared thermal imager, ultrasonic concentration meter, and micro-pressure differential sensor according to their respective confidence levels, thereby obtaining more reliable and comprehensive data. The improved DS evidence theory here is optimized based on the traditional DS evidence theory and is used to process the fusion of uncertain information. The process fingerprint extraction unit uses the mathematical method of wavelet packet decomposition to extract the time-frequency domain feature matrix from the vacuum pulsation signal, and then constructs a process state fingerprint library containing kurtosis factor, energy entropy, and Lyapunov exponent to characterize the current process state.
[0065] Specifically, an infrared thermal imager is used to obtain image data of the temperature distribution in the evaporation chamber. An ultrasonic concentration meter measures material concentration based on the propagation characteristics of ultrasonic waves in the material. A micro-differential pressure sensor measures the pressure difference in the system, thereby reflecting the degree of vacuum. In the improved DS evidence theory, the confidence level needs to be determined based on factors such as the sensor's accuracy and stability, as well as the reliability of historical data. For example, a high-precision and stable sensor can have a higher confidence level. During wavelet packet decomposition, the number of decomposition layers can be selected based on the complexity of the signal and the accuracy of the desired features, generally between 3 and 8 layers. The kurtosis factor measures the degree to which the signal deviates from a normal distribution, the energy entropy reflects the signal's energy distribution, and the Lyapunov exponent is used to determine the chaotic characteristics of the system. When constructing a process state fingerprint library, these parameters need to be standardized to facilitate comparison and analysis.
[0066] Preferably, in the multimodal data fusion unit, a dynamic confidence update method can be used to improve the fusion effect. That is, the confidence weight of each sensor data is adjusted in real time based on the fluctuation and accuracy feedback of the real-time sensor data. In the process fingerprint extraction unit, the characteristic matrix after wavelet packet decomposition can be further reduced in dimension through principal component analysis (PCA) to remove redundant information and improve the efficiency of subsequent process state analysis. For the calculation of kurtosis factor, energy entropy, and Lyapunov exponent, parallel computing technology can be used to accelerate the calculation speed, thereby achieving real-time monitoring and rapid response to process status.
[0067] In some embodiments, the digital twin modeling module includes:
[0068] Variational Autoencoder (VAE), which maps high-dimensional sensor data into a latent space to generate a low-dimensional manifold representation of the process state;
[0069] The physical information neural network (PINN) embeds the Navier-Stokes equation constraints to construct a differential-data hybrid driven model of the evaporation process, and its loss function is:
[0070]
[0071] Where u is the fluid velocity field, T is the temperature field, α is the thermal diffusion coefficient, λ i To balance the weight.
[0072] It should be noted that the digital twin modeling module consists of a variational autoencoder (VAE) and a physical information neural network (PINN). The role of the variational autoencoder is to map high-dimensional sensor data to a latent space through a specific algorithm, thereby generating a low-dimensional manifold representation of the process state. Simply put, it converts complex high-dimensional data into a more concise and easier to analyze low-dimensional data form. The latent space is an abstract space in which data can express key features in a more compact manner. The physical information neural network (PINN) embeds the constraints of the Navier-Stokes equations to construct a differential-data hybrid driven model of the evaporation process. It combines physical laws and actual data to make the model more accurate and reliable. The loss function is used to measure the difference between the model's predicted value and the actual value, thereby optimizing the model parameters.
[0073] Specifically, high-dimensional sensor data includes data collected from various sensors in the multi-source perception module, such as the temperature gradient distribution of the evaporation chamber, the dynamic waveform of the vacuum degree, and other data. The dimensions are high and the information is complex. In the variational autoencoder, its network structure includes an encoder and a decoder. The encoder is responsible for mapping the high-dimensional data to the latent space, and the decoder restores the latent space data to high-dimensional data. During the training process, parameters such as the learning rate and the number of iterations need to be set. The learning rate is generally between 0.001 and 0.0001, and the number of iterations is set according to the amount of data and the complexity of the model, usually ranging from hundreds to thousands of times. For physical information neural networks (PINNs), the Navier-Stokes equations describe the motion laws of viscous Newtonian fluids, which play a constraint role in the model to ensure that the model conforms to physical principles. The balance weight λ in the loss function i The thermal diffusion coefficient α determines the importance of different constraints in model training. Its value needs to be determined through multiple experiments and adjustments, and is generally between 0 and 1. The thermal diffusion coefficient α depends on the physical properties of the material and the environmental conditions. It can be obtained through experimental measurement or reference to relevant literature.
[0074] Preferably, when using a variational autoencoder, a transfer learning approach can be used to utilize a model pre-trained on similar process data to speed up the training on the current process data and improve the model performance. When training a physical information neural network (PINN), in order to more accurately determine the balance weight λ i , optimization algorithms such as genetic algorithms can be used for search. When constructing a differential-data hybrid driven model for the evaporation process, in addition to embedding the Navier-Stokes equations, other relevant physical equations, such as the energy conservation equation, can also be considered to further improve the accuracy and physical significance of the model.
[0075] In some embodiments, the multi-objective optimization module includes:
[0076] Hierarchical decision-making architecture: the upper layer uses the NSGA-III algorithm to generate the global Pareto solution set, and the lower layer handles parameter uncertainty through distributed robust optimization (DRO);
[0077] Constraint processing unit, defining the dynamic feasible region:
[0078]
[0079] Where Γ1(t) is the time-varying safety threshold and Γ2(t) is the semi-positive definite matrix constraint.
[0080] It should be noted that the multi-objective optimization module mainly includes a hierarchical decision-making architecture and a constraint processing unit. In the hierarchical decision-making architecture, the upper layer uses the NSGA-III algorithm to generate a global Pareto solution set. The NSGA-III algorithm is the third generation of the non-dominated sorting genetic algorithm. It is a multi-objective evolutionary algorithm that searches for the optimal solution by simulating the natural evolution process. The global Pareto solution set refers to a set of solutions in the entire solution space that meet multiple objectives without any distinction between good and bad, and cannot improve the performance of a certain objective without reducing the performance of other objectives. The lower layer uses distributed robust optimization (DRO) to deal with parameter uncertainty. DRO is an optimization method that takes parameter uncertainty into account and can keep the optimization results robust within a certain uncertainty range. The constraint processing unit limits the scope of the solution by defining a dynamic feasible domain. The dynamic feasible domain changes over time and is expressed as a mathematical expression:
[0081] F(t)={x∈R n |g1(x,t)≤Γ1(T),g2(x,t)≤Γ2(t)}
[0082] Among them, Γ1(t) is a time-varying safety threshold used to limit certain key indicators from exceeding a specific range; ε2(t) is a semi-positive definite matrix constraint that ensures that the properties of the correlation matrix meet certain conditions.
[0083] Specifically, in the NSGA-III algorithm, the population size and genetic operation parameters must be set. The population size is typically determined based on the complexity of the problem and computing resources, generally ranging from tens to hundreds. Genetic operations include selection, crossover, and mutation. Selection can employ a tournament selection method, with the tournament size typically set to 2-5. The crossover probability is typically between 0.7 and 0.95, and the mutation probability is between 0.01 and 0.1. For distributed robust optimization, when constructing Wasserstein fuzzy sets, the value of ∈ must be determined. This controls the size of the fuzzy set, i.e., the range of uncertainty. The value of ∈ is typically set based on an estimate of parameter uncertainty and the acceptable risk level, generally between 0.01 and 0.1. Within the constraint processing unit, the time-varying safety threshold Γ1(t) is determined based on the dynamic characteristics and safety requirements of the process. For example, in a vacuum concentration process, the pressure safety threshold of the evaporation chamber varies as the concentration process progresses. The matrix elements associated with the semi-positive definite matrix constraint Γ2(t) can be set based on the physical characteristics and operating constraints of the process equipment.
[0084] Preferably, when using the NSGA-III algorithm to generate a global Pareto solution set, an elite retention strategy can be adopted to retain the excellent individuals of each generation directly to the next generation, thereby avoiding the loss of excellent solutions and accelerating the convergence speed. For distributed robust optimization, in addition to using Wasserstein fuzzy sets, other types of fuzzy sets, such as e1-fuzzy sets, e2-fuzzy sets, etc., can also be considered. The most appropriate fuzzy set form is selected according to the characteristics of parameter uncertainty in the actual problem. When defining the dynamic feasible domain, for the time-varying safety threshold Γ1(t), an adaptive adjustment method can be adopted to dynamically update the threshold according to the real-time operating status and historical data of the system to make it more in line with the actual process requirements; for the semi-positive definite matrix constraint Γ2(t), the eigenvalues of the matrix can be checked regularly to ensure that it always meets the semi-positive definite condition. If it does not meet the condition, corresponding adjustments or recalculations are made.
[0085] In some embodiments, the distributed robust optimization uses Wasserstein fuzzy sets:
[0086]
[0087] Among them, P ∈ ={Q:W1(Q,P N )≤∈}, ξ is an uncertain parameter, and W1 is the 1-Wasserstein distance.
[0088] It should be noted that Wasserstein fuzzy sets are used in the distributed robust optimization of the multi-objective optimization module. In the case of uncertain parameters ξ, it means to find the control variable u in the set U so as to minimize the expected supremum of the objective function J(u,ξ) under the fuzzy distribution Q. ∈ ={Q:W1(Q,P N )≤∈} defines the Wasserstein fuzzy set, and W1 is the 1-Wasserstein distance, which is used to measure the degree of difference between two probability distributions. N is a nominal distribution, which can be understood as a reference distribution based on existing data or experience. ∈ is a non-negative parameter that controls the size of the fuzzy set, i.e., the range of uncertainty. In this way, parameter uncertainty is taken into account during the optimization process, making the optimization results more robust.
[0089] Specifically, the uncertain parameter ξ in the vacuum concentration process may represent factors that are difficult to accurately measure or predict, such as the initial concentration of the material and ambient temperature fluctuations. N The determination of is usually based on the statistical analysis of historical process data, such as the statistical analysis of parameters such as material concentration and ambient temperature in a large number of batch production data, to obtain their probability distribution as the nominal distribution. N The calculation of ) will adopt the corresponding algorithm according to the specific distribution form. For discrete distribution, it can be obtained by calculating the sum of the product of the distance and probability between different sample points. The setting of ∈ needs to balance the conservatism and effectiveness of the optimization results. If ∈ is set to a small value, the fuzzy set P ∈ It is closer to the nominal distribution P N , the optimization result is more dependent on the nominal distribution, and may be less adaptable when facing actual uncertainty, but the calculation is relatively simple; if ∈ is set larger, the distribution range covered by the fuzzy set is wider, and the optimization result is more resistant to uncertainty, but the calculation complexity will increase. Generally, it can be adjusted within the range of 0.01-0.1 according to the actual situation. In calculating the expected E of the objective function J(u,ξ) Q [J(u,ξ)], it is necessary to use the corresponding integration or summation method to calculate according to the distribution form of Q and the specific expression of J(u,ξ).
[0090] Preferably, when determining the nominal distribution P NWhen , a dynamic update method can be used. As new process data continues to accumulate, the data is re-analyzed regularly and the nominal distribution is updated to make it more in line with the actual situation. For the selection of ∈, an adaptive adjustment strategy can be adopted. For example, in the early stage of process operation, due to the lack of understanding of uncertainty, the value of ∈ can be appropriately increased to ensure the robustness of the optimization results; as the process stabilizes and the uncertainty is gradually mastered, ∈ can be dynamically reduced according to the actual response of the system to improve the accuracy and efficiency of the optimization. When calculating the 1-Wasserstein distance W1(Q,P N ), if the amount of data is large, an approximate calculation method can be used, such as a sample-based estimation method, to reduce the amount of calculation while ensuring a certain degree of accuracy. Q [J(u,ξ)], if the distribution of Q is more complex, the Monte Carlo simulation method can be used for estimation, and the expectation can be approximated through a large number of random sampling to improve the calculation efficiency and accuracy.
[0091] In some embodiments, the adaptive execution module includes:
[0092] Feedforward-feedback composite controller, the feedforward channel uses Volterra series to compensate for nonlinear time lag, and the feedback channel is designed as a fractional-order PID controller:
[0093] u(t)=K p e(t)+K i D -λ e(t)+K d D μ e(t)
[0094] Where λ∈(0,1),μ∈(0,1) are fractional orders, and D is a fractional differential operator;
[0095] The actuator dynamic compensation unit establishes the Hammerstein-Wiener model of the servo motor and the pneumatic valve, and realizes dynamic linearization through inverse model decoupling.
[0096] It should be noted that the adaptive execution module consists of a feedforward-feedback composite controller and an actuator dynamic compensation unit. The feedforward-feedback composite controller combines the advantages of feedforward control and feedback control. The feedforward channel uses the Volterra series to compensate for the nonlinear time lag in the system. The Volterra series is a mathematical tool used to describe the input-output relationship of a nonlinear system. It can be used to model and compensate for nonlinear time lag. The feedback channel uses a fractional-order PID controller, which is expressed as
[0097] u(t)=K p e(t)+K i D -λ e(t)+Kd D μ e(t)
[0098] Here, λ∈(0,1) and μ∈(0,1) represent fractional orders, and D represents a fractional differential operator. Compared to traditional PID controllers, fractional-order PID controllers can more flexibly adjust control parameters and improve control performance. The actuator dynamic compensation unit establishes a Hammerstein-Wiener model of the servo motor and pneumatic valve. This model is a series-connected model consisting of a static nonlinear link, a linear dynamic link, and a static nonlinear link. It is used to describe the complex characteristics of the actuator. The actuator's dynamic linearization is then achieved through inverse model decoupling, ensuring that the actuator's output better meets the system control requirements.
[0099] Specifically, when using the Volterra series to compensate for nonlinear time lag, it is necessary to determine the order of the Volterra series and the kernel function. The choice of order depends on the complexity of the system nonlinearity. Generally, it can be determined by analyzing the input and output data of the system and combining trial and error. It is usually between 2 and 5 orders. The kernel function is more complicated to determine and can usually be updated online using the recursive least squares algorithm. For the fractional-order PID controller, the proportional coefficient K p , integral coefficient K i and differential coefficient K d The settings are similar to those of traditional PID controllers, but the fractional orders λ and μ need to be considered additionally. These parameters can be initially determined through empirical formulas, the Ziegler-Nichols method, etc., and then fine-tuned based on the actual response of the system. The value range of λ and μ is between (0,1). They determine the weight of the controller's integral and differential terms for the error signal. The closer the value is to 1, the stronger the effect of the corresponding term. When establishing the Hammerstein-Wiener model, it is necessary to test and analyze the static nonlinear characteristics and linear dynamic characteristics of the servo motor and pneumatic valve. For example, by inputting signals of different frequencies and amplitudes, measuring the output of the actuator, and obtaining data to fit the model parameters. The implementation of inverse model decoupling is based on the structure of the Hammerstein-Wiener model. The inverse model is obtained through mathematical transformation to eliminate the nonlinear and dynamic coupling effects of the actuator.
[0100] Preferably, when using Volterra series to compensate for nonlinear time lag, an adaptive Volterra series model can be used. According to the real-time operating status of the system, the order and kernel function of the Volterra series are dynamically adjusted to better adapt to system changes. For fractional-order PID controllers, intelligent optimization algorithms such as particle swarm optimization (PSO) or genetic algorithm (GA) can be used to optimize the parameter K. p , K i, K d , λ, and μ. These algorithms can search for optimal solutions in a wider parameter space, achieving better control results than traditional methods. When establishing the Hammerstein-Wiener model, if the actuator characteristics change over time, online modeling can be used to regularly update the model parameters to ensure model accuracy. During the inverse model decoupling process, neural network technology can be combined to leverage the powerful nonlinear mapping capabilities of neural networks to more accurately achieve dynamic linearization of the actuator, improving the system's control accuracy and response speed.
[0101] In some embodiments, the Volterra series kernel function is updated online by a recursive least squares algorithm:
[0102]
[0103] where γ k is a variable step learning rate, satisfying
[0104] It should be noted that the Volterra series kernel function is updated online through the recursive least squares algorithm. The recursive least squares algorithm is an algorithm that can gradually update the estimated parameters as new data continues to arrive. Online updating means that the Volterra series kernel function is adjusted according to the real-time data as time goes by. represents the estimated value of the n-th order Volterra series kernel function after the k+1th update, is the estimated value after the kth update. k It is a variable step learning rate, which determines the magnitude of kernel function adjustment at each update and needs to satisfy ∑γ k =∞, This ensures that the algorithm can fully explore the parameter space and ensure convergence during long-term operation. i ) is the system at t-τ i The input at time t, y(t) is the actual output of the system at time t, It is the system prediction output calculated based on the kernel function after the k-th update.
[0105] Specifically, when using the recursive least squares algorithm to update the Volterra series kernel function, we must first determine the order n of the Volterra series. Its selection depends on the complexity of the nonlinearity of the system. It can generally be determined through preliminary analysis of the system input and output data and multiple experiments. The value is usually between 2 and 5. k The initial value of can be set based on experience, for example, between 0.01 and 0.1. k =∞, Common variable step learning rate forms are
[0106]
[0107] Among them, γ0 is the initial learning rate, β is the adjustment parameter, and both γ0 and β need to be adjusted according to the convergence speed and stability of the actual system. i It represents the time delay, which is related to the response time of the system and can be determined by analyzing the dynamic characteristics of the system or by experimental measurement.
[0108]
[0109] When the system input and output data are discrete, numerical integration methods such as trapezoidal integration method or Simpson integration method can be used for approximate calculation.
[0110] Preferably, in the process of updating the Volterra series kernel function, the recursive least squares algorithm can be improved by using a forgetting factor. f (values between 0.9-0.99) can make the algorithm pay more attention to new data and gradually "forget" old data, so that it can track changes in system parameters more quickly. When the system operating environment suddenly changes, it is detected that the change in input and output data exceeds a certain threshold, and γ can be appropriately increased. k , speed up the update speed of the kernel function; when the system runs relatively stably, reduce γ k , improving the algorithm's convergence accuracy. When calculating the predicted output, for complex nonlinear systems, a numerical integration method with an adaptive integration step can be used. This dynamically adjusts the integration step size based on the changes in the integrand, improving computational efficiency and accuracy. Furthermore, to prevent error accumulation during the numerical calculation process, the kernel function can be regularly normalized.
[0111] In some embodiments, it further includes:
[0112] The cognitive digital twin module builds a fault propagation model based on the knowledge graph and defines the impact of fault modes:
[0113]
[0114] where w i is the fault link weight, Δθ i is the parameter offset, σ i is the sensitivity coefficient;
[0115] Self-repair unit, when Ψ>Ψ th The reconstruction mechanism is triggered in real time to maintain the key functions of the system through structural adaptive reorganization.
[0116] It should be noted that the system includes a cognitive digital twin module and an autonomous repair unit. The cognitive digital twin module analyzes and predicts possible system failures by building a fault propagation model based on a knowledge graph. A knowledge graph is a semantic network used to describe the relationship between entities. Here, it represents the various components of the system, the types of faults, and the relationships between them. The failure mode impact Ψ is calculated by the formula
[0117]
[0118] To calculate, it is used to measure the impact of different failure modes on the system. The autonomous repair unit will be restored when the failure mode impact Ψ is greater than the set threshold Ψ th The reconstruction mechanism is activated in real time to maintain the normal operation of the key functions of the system through structural adaptive reorganization.
[0119] Specifically, when building a fault propagation model based on a knowledge graph, it is necessary to determine the entities and relationships in the knowledge graph. Entities can include various hardware devices, software modules, fault types, etc. in the system, and relationships describe the causal relationship between them, such as a hardware fault can cause a software module to be abnormal. The fault link weight w i It reflects the importance of different fault links in the entire fault propagation process. It can be set based on statistical analysis of historical fault data and expert experience. The value range is generally between 0 and 1. The larger the weight, the greater the impact of the fault link on the system. i The difference between real-time monitoring system parameters and normal operating parameters is obtained, such as changes in temperature, pressure and other parameters. Sensitivity coefficient σ i The threshold Ψ depends on the system's sensitivity to changes in different parameters, which can also be determined based on historical data and system characteristics. th The setting of must comprehensively consider the system's fault tolerance and operational stability. An appropriate value is typically determined through extensive simulation experiments and actual operational data. During the adaptive restructuring process, it is important to identify the system's key functions and how to adjust the system structure, such as reallocating resources and adjusting control strategies.
[0120] Preferably, a dynamic update method can be used when constructing the knowledge graph. As new fault data and experience are continuously generated during the operation of the system, the entities and relationships in the knowledge graph are updated in a timely manner to make the fault propagation model more accurate. iMachine learning algorithms, such as neural networks, can be used to automatically learn and adjust weights based on a large number of failure cases, improving the accuracy and rationality of weight settings. When calculating the failure mode impact Ψ, additional factors, such as the frequency of the failure, can be added to the formula to more accurately reflect the impact of the failure on the system. When the reconstruction mechanism is triggered, in addition to adaptive structural reorganization based on pre-set rules, reinforcement learning algorithms can also be used to find the optimal solution among different reconstruction strategies, achieving more efficient system repair and maintaining critical functions.
[0121] In some embodiments, the reconstruction mechanism includes:
[0122] Control strategy migration: mapping the current working condition to the nearest case library scenario and loading the pre-trained reinforcement learning strategy;
[0123] Resource reallocation, optimizing the sensor-actuator pairing relationship based on the Hungarian algorithm, and dynamically reconstructing it into multiple functionally independent virtual subsystems.
[0124] It should be noted that the system's reconfiguration mechanism includes two parts: control strategy migration and resource reallocation. Control strategy migration involves matching the current system's operating conditions with scenarios in the case library, finding the most similar scenario, and then loading the pre-trained reinforcement learning strategy for that scenario to adjust the system's control method. The case library here is a database that stores various historical operating conditions and corresponding effective control strategies. Reinforcement learning strategies are a series of decision rules learned by the intelligent agent through continuous trial and error in the environment and based on reward feedback. Resource reallocation uses the Hungarian algorithm to optimize the pairing relationship between sensors and actuators, recombining them into multiple independent and fully functional virtual subsystems, thereby improving the overall performance and reliability of the system. The Hungarian algorithm is a combinatorial optimization algorithm for solving allocation problems and can find the optimal allocation solution.
[0125] Specifically, when migrating control strategies, an effective case library should be constructed. Each case in the case library should contain operating condition characteristic data and a corresponding reinforcement learning strategy. These operating condition characteristic data can include key process parameters such as evaporation chamber temperature, vacuum level, and material flow rate. When mapping the current operating condition to the case library scenarios, a similarity calculation method must be determined. For example, Euclidean distance or cosine similarity can be used to measure the similarity between the current operating condition and each scenario in the case library. Pre-training of the reinforcement learning strategy can use algorithms such as the Deep Q-Network (DQN) through multiple rounds of training in a simulated environment. During training, appropriate parameters such as the learning rate (e.g., 0.001) and discount factor (e.g., 0.9) should be set. For resource reallocation, when applying the Hungarian algorithm, a sensor-actuator cost matrix must first be constructed. The elements in the cost matrix represent the cost of pairing different sensors and actuators. These costs can be determined based on factors such as signal transmission loss and loss of control accuracy. The Hungarian algorithm uses this cost matrix to find the optimal pairing solution, achieving optimal resource allocation. The division of virtual subsystems should be based on the system's functional modules and control requirements to ensure that each virtual subsystem is functionally independent and can work together.
[0126] Preferably, in the control strategy migration, in order to more accurately match the current working condition with the case library scenario, the working condition characteristic data can be subjected to feature engineering processing to extract more representative features and improve the accuracy of the matching. At the same time, the case library is updated in an incremental learning manner. When new working conditions and effective control strategies appear in the system, they are added to the case library in a timely manner to continuously improve the case library. When reallocating resources, in addition to considering the cost of sensor-actuator pairing, time factors can also be introduced, such as performance changes of sensors and actuators in different time periods, and the cost matrix can be dynamically adjusted to make resource allocation more adaptable to the real-time operating status of the system. In addition, after dividing the virtual subsystems, the performance of the virtual subsystems can be evaluated regularly. If it is found that the performance of a virtual subsystem has declined, the Hungarian algorithm can be re-applied for resource reallocation to ensure that the system is always in an efficient operating state.
[0127] In some embodiments, further comprising:
[0128] Energy and mass flow optimization module, establish Analytical model calculation of each unit Loss coefficient:
[0129]
[0130] Where T0 is the ambient temperature, m i is the mass flow rate, s i is the specific entropy, Q i is the heat flow rate;
[0131] The waste heat intelligent recovery unit decides the optimal waste heat utilization path based on the fuzzy cognitive map, giving priority to meeting the energy demand of the preheating stage.
[0132] It should be noted that the system is also equipped with an energy and mass flow optimization module and a waste heat intelligent recovery unit. Analyze the model to calculate the Loss coefficient, It is a physical quantity that measures the quality of energy. Analytical models are based on thermodynamic principles and are used to assess how efficiently energy is used in a system. The loss coefficient calculation formula is:
[0133]
[0134] This coefficient reflects the energy conversion and transfer process of each unit in the system. The waste heat intelligent recovery unit uses fuzzy cognitive maps to determine the optimal solution from multiple possible waste heat utilization paths, prioritizing the energy needs of the preheating stage to improve overall energy utilization. Fuzzy cognitive maps are a soft computing method that uses concept nodes and weighted directed arcs to represent the causal relationships between various factors in the system, allowing decisions to be made through reasoning.
[0135] Specifically, in the energy and mass flow optimization module, the ambient temperature T0 can be measured and obtained by an ambient temperature sensor, and is generally based on the actual local ambient temperature. i Flow rate measurement instruments can be used to measure the material flow rate according to the material inlet and outlet conditions of different units. For example, electromagnetic flowmeters can be installed on the material conveying pipeline to measure the mass flow rate. i The value of is related to the type and state of the material (such as temperature, pressure, etc.), and can be obtained by consulting the relevant thermodynamic property table or using a specific calculation formula. i It can be calculated by measuring the temperature distribution and thermal conductivity coefficient on the surface of the equipment and combining it with the heat conduction equation, or it can be measured using a special heat flow sensor. In the waste heat intelligent recovery unit, the construction of the fuzzy cognitive map requires determining the weights of concept nodes and directed arcs. Concept nodes can include waste heat generating equipment, preheating equipment, different waste heat transmission paths, etc.; the weights of directed arcs reflect the strength of the causal relationship between different factors, and are generally determined through expert experience, historical data statistical analysis, or machine learning algorithms. When determining the best waste heat utilization path, the fuzzy cognitive map calculates the priority of each path through reasoning and selects the path with the highest priority as the best utilization path.
[0136] Preferably, in the energy flow optimization module, When calculating the loss coefficient, in order to improve the calculation accuracy, the method of dynamic parameter update can be used. Since the parameters may change during the operation of the system, the mass flow rate, heat flow rate and other parameters are monitored in real time. When the parameter change exceeds a certain threshold, the calculation is updated in time. The loss coefficient can be used to more accurately assess the energy utilization of each unit in the system. For intelligent waste heat recovery units, in addition to utilizing fuzzy cognitive maps, deep learning algorithms can also be combined to optimize waste heat utilization paths. For example, deep neural networks can be used to learn from historical waste heat utilization data and uncover potential patterns and regularities within the data, allowing for more accurate predictions of the effects of different waste heat utilization paths and further improving waste heat recovery efficiency. Furthermore, intelligent control devices can be introduced during the waste heat utilization process to dynamically adjust waste heat distribution based on the real-time energy demand of the preheating section, ensuring that the preheating section always receives an adequate energy supply.
[0137] The above-mentioned embodiments of the present invention have the following beneficial effects: by collecting a variety of process parameters in real time through the multi-source perception module, combined with the high-precision modeling capability of the digital twin modeling module, the system can achieve accurate perception and dynamic modeling of complex concentration processes, providing a solid foundation for optimization control. The multi-objective optimization module adopts a mixed integer dynamic programming algorithm, which can solve the Pareto front solution set in the non-convex solution space, and simultaneously optimize the energy efficiency ratio, product yield and equipment life indicators, thereby achieving high efficiency, energy saving and sustainable development of the process. The adaptive execution module compensates for the system lag effect in real time through nonlinear model predictive control, and establishes a dual closed-loop anti-interference adjustment mechanism, which can significantly improve the dynamic response capability and anti-interference performance of the system, and ensure the stable operation of the concentration process under complex working conditions. In addition, the system also has a cognitive digital twin module and an autonomous repair unit, which can quickly identify and handle faults through fault propagation models and reconstruction mechanisms, maintain key system functions, and further enhance the reliability and stability of the system.
[0138] In terms of energy and mass flow optimization, the system establishes Analytical model calculation of each unit The loss coefficient, combined with an intelligent waste heat recovery unit, enables efficient energy utilization and rational waste heat recovery, prioritizing the energy needs of the preheating stage and further improving the system's energy efficiency. Furthermore, the application of technologies such as distributed robust optimization and fractional-order PID controllers effectively addresses parameter uncertainty and nonlinear time lag, ensuring optimized system performance and control accuracy, providing comprehensive technical support for the intelligent upgrade of vacuum concentration processes.
[0139] Furthermore, the storage medium of the embodiment of the present application stores program instructions that can implement all the above methods, wherein the program instructions can be stored in the above storage medium in the form of a software product, including a number of instructions for causing a computer device (which can be a personal computer, server, or network device, etc.) or a processor to execute all or part of the steps of the method described in each embodiment of the present application. The aforementioned storage medium includes: various media that can store program codes, such as a USB flash drive, a mobile hard disk, a read-only memory (ROM), a random access memory (RAM), a magnetic disk or an optical disk, or a terminal device such as a computer, a server, a mobile phone, or a tablet.
[0140] The above descriptions are merely some preferred embodiments of the present invention and an illustration of the technical principles employed. Those skilled in the art should understand that the scope of the invention involved in the embodiments of the present invention is not limited to the technical solutions formed by a specific combination of the above-mentioned technical features, but should also encompass other technical solutions formed by any combination of the above-mentioned technical features or their equivalents without departing from the above-mentioned inventive concept. For example, a technical solution formed by mutually replacing the above-mentioned features with (but not limited to) technical features having similar functions disclosed in the embodiments of the present invention.
Claims
1. A process optimization analysis system for vacuum concentration, characterized in that: Includes the following multi-layer cascade modules: The multi-source sensing module uses a distributed sensor array to collect real-time data on the temperature gradient distribution in the evaporation chamber, the dynamic waveform of the vacuum degree, the rheological properties of the material, and the latent heat parameters of the phase change; The digital twin modeling module builds a hidden Markov chain model of process parameters and concentration dynamics based on transfer learning, and integrates historical batch data to generate an adaptive covariance matrix; The multi-objective optimization module uses a mixed integer dynamic programming algorithm to solve the Pareto front solution set in a non-convex solution space, and simultaneously optimizes energy efficiency, product yield, and equipment life indicators; The adaptive execution module compensates for the system hysteresis effect in real time through nonlinear model predictive control and establishes a dual closed-loop anti-interference regulation mechanism.
2. The system according to claim 1, wherein: The multi-source perception module includes: The multimodal data fusion unit uses the improved DS evidence theory to perform confidence-weighted fusion on heterogeneous data from infrared thermal imagers, ultrasonic densitometers, and micro-pressure differential sensors; The process fingerprint extraction unit extracts the time-frequency domain feature matrix of the vacuum pulsation signal through wavelet packet decomposition, and constructs a process state fingerprint library including kurtosis factor, energy entropy and Lyapunov exponent.
3. The system according to claim 1, wherein: The digital twin modeling module includes: Variational Autoencoder (VAE) maps high-dimensional sensor data into a latent space to generate a low-dimensional manifold representation of the process state; The physical information neural network PINN embeds the Navier-Stokes equation constraints to construct a differential-data hybrid driven model of the evaporation process, and its loss function is: Where u is the fluid velocity field, T is the temperature field, α is the thermal diffusion coefficient, λ i To balance the weight.
4. The system according to claim 1, wherein: The multi-objective optimization module includes: Hierarchical decision-making architecture: the upper layer uses the NSGA-III algorithm to generate the global Pareto solution set, and the lower layer handles parameter uncertainty through distributed robust optimization (DRO); Constraint processing unit, defining the dynamic feasible region: Where Γ1(t) is the time-varying safety threshold and Γ2(t) is the semi-positive definite matrix constraint.
5. The system according to claim 4, characterized in that The distributed robust optimization adopts Wasserstein fuzzy sets: Among them, P ∈ ={Q:W1(Q,P N )≤∈}, ξ is an uncertain parameter, and W1 is the 1-Wasserstein distance.
6. The system according to claim 1, wherein: The adaptive execution module includes: Feedforward-feedback composite controller, the feedforward channel uses Volterra series to compensate for nonlinear time lag, and the feedback channel is designed as a fractional-order PID controller: u(t)=K p e(t)+K i D -λ e(t)+K d D μ e(t) Where λ∈(0,1),μ∈(0,1) are fractional orders, and D is a fractional differential operator; The actuator dynamic compensation unit establishes the Hammerstein-Wiener model of the servo motor and the pneumatic valve, and realizes dynamic linearization through inverse model decoupling.
7. The system according to claim 6, characterized in that The Volterra series kernel function is updated online by the recursive least squares algorithm: Among them, γ k is a variable step learning rate, satisfying 8. The system according to claim 1, wherein: Also includes: The cognitive digital twin module builds a fault propagation model based on the knowledge graph and defines the impact of fault modes: where w i is the fault link weight, Δθ i is the parameter offset, σ i is the sensitivity coefficient; Self-repair unit, when Ψ>Ψ th The reconstruction mechanism is triggered in real time to maintain the key functions of the system through structural adaptive reorganization.
9. The system according to claim 8, characterized in that The reconstruction mechanism includes: Control strategy migration: mapping the current working condition to the nearest case library scenario and loading the pre-trained reinforcement learning strategy; Resource reallocation, optimizing the sensor-actuator pairing relationship based on the Hungarian algorithm, and dynamically reconstructing it into multiple functionally independent virtual subsystems.
10. The system according to claim 1, wherein: Also includes: Energy and mass flow optimization module, establish Analytical model calculation of each unit Loss coefficient: Where T0 is the ambient temperature, m i is the mass flow rate, s i is the specific entropy, Q i is the heat flow rate; The waste heat intelligent recovery unit decides the optimal waste heat utilization path based on the fuzzy cognitive map, giving priority to meeting the energy demand of the preheating stage.
Citation Information
Cited By
Direct modified rubber powder asphalt production tank control method
CN120871800A
Pesticide production control system and method
CN120909389A
Freeze-dried tremella fuciformis energy efficiency-quality double-excellent control system based on carrier cooking process
CN120972843A