Crystal size distribution soft measurement method integrated with mechanism constraint

Abstract

The invention discloses a crystal size distribution soft measurement method integrated with mechanism constraints. The method comprises the following steps: 1) data acquisition and integration; 2) modeling training and model testing: constructing a crystal size distribution soft measurement model integrated with mechanism constraint, training a training data set, injecting crystallization process dynamics knowledge into the Kolmogorov-Arnod network in a loss function form, improving the nonlinear modeling capability and credibility of the model, and performing model testing on the Kolmogorov-Arnod network; and verifying test set data by using the trained model. According to the method, the knowledge of the crystallization dynamic process is injected into the Kolmogorov-Arnod network in the form of the loss function, so that the prediction capability of the model on the crystal size distribution in the crystallization process and the credibility of the model are effectively improved.

Description

Technical Field

[0001] The present invention belongs to the technical field of soft sensor modeling in the process industry, and particularly relates to a soft measurement method for crystal size distribution incorporating mechanism constraints. Background Art

[0002] Batch crystallization, as a key separation and purification technology in the process industry, plays an important role in fields such as pharmaceuticals, chemicals, and food. By controlling conditions such as solution supersaturation, temperature, and agitation, it can efficiently achieve high-purity separation of target substances and regulation of specific crystal morphologies. Crystal size distribution (CSD), as one of the core indicators of the batch crystallization process, directly affects the quality, performance of the product, and the efficiency of subsequent processing. Accurately predicting the crystal size distribution is of great significance for optimizing the crystallization process, improving product consistency, and reducing production costs. However, CSD is difficult to directly measure online, and the measurement results often have significant delays, making it difficult to meet the requirements of real-time control. Therefore, soft sensor technology, as an effective alternative, has become a key means to solve this problem by using easily measurable process variables to model CSD. Due to the highly nonlinear and dynamic characteristics of the batch crystallization process and its complex internal mechanism, the construction of traditional soft measurement models based on process mechanisms faces many challenges. Data-driven soft sensors have gradually become the mainstream direction for studying the batch crystallization process because they can achieve modeling and prediction of complex processes without relying on in-depth understanding of internal mechanisms.

[0003] In the research of data-driven soft sensors, nonlinear modeling and interpretability are two key challenges that are interrelated and jointly affect the practicality and reliability of the model. First, the challenge of nonlinear modeling stems from the highly nonlinear and dynamic characteristics of the batch crystallization process. To accurately describe this complexity, the model needs to have stronger expressive power, usually requiring more parameters to capture the complex mapping relationship between inputs and outputs, and relying on a large amount of high-quality data. However, in the actual industrial environment, it is often very difficult to obtain sufficient and representative data. Insufficient or poor-quality data can lead to overfitting or underfitting of the model, thus affecting its prediction accuracy and generalization ability. Second, interpretability directly affects the credibility of the model prediction results. If the model prediction lacks a clear basis, engineers may doubt the results when operating conditions change, leading to decision-making errors. These two challenges are closely related: Solving the nonlinear problem can improve the model's ability to capture physical mechanisms, reduce prediction errors, and enhance engineers' trust; at the same time, improving interpretability helps the model identify key driving factors and guide the improvement of nonlinear modeling. Collaboratively solving these two problems can improve the model accuracy and reliability and promote the practical application of soft sensor technology.

[0004] In soft-sensing tasks, non-linear modeling has always been a long-standing difficult problem. Especially in complex systems, simple linear models often cannot meet the actual needs. Although traditional methods such as Convolutional Neural Networks (CNNs) and Gated Recurrent Units (GRUs) have been widely used to solve this difficult problem, these methods usually rely on a large amount of labeled data and are limited in performance when data is scarce. For complex industrial processes such as the batch crystallization process, due to high experimental costs and difficult data acquisition, there is often a lack of a large amount of labeled data. Therefore, these methods face significant challenges in practical applications. In contrast, the Kolmogorov-Arnold Network (KAN) can achieve excellent non-linear learning using relatively less data in regression tasks. However, when dealing with complex physical systems such as complex industrial processes like the batch crystallization process, the KAN model is still difficult to effectively incorporate domain knowledge, resulting in limited credibility and interpretability of its prediction results.

[0005] To effectively address the challenge of interpretability, the Physics-Informed Neural Network (PINN) provides an effective solution. The core of PINN lies in embedding physical laws (such as partial differential equations, conservation laws, etc.) as constraints into the training process of the neural network. Different from traditional data-driven models, PINN not only relies on data but also uses physical equations to guide the learning process of the model, thus ensuring that the output of the model conforms to known physical laws. This effectively improves the credibility and accuracy of the model. For example, Jia et al. proposed physical-anchored graph learning for predicting key process indicators, effectively improving the prediction of key quality variables in the wastewater treatment process. Franklin et al. proposed a method for well flow measurement oriented to physics-informed neural networks, effectively improving well flow prediction. Therefore, embedding the physical constraints of crystal distribution in the batch crystallization process into KAN can effectively enhance the interpretability and accuracy of soft sensing.

[0006] Based on the above discussion, a KAN soft-sensing model incorporating mechanism constraints (PIKAN) is proposed to address the challenges of non-linearity and interpretability in modeling the crystal size distribution of the batch crystallization process. First, KAN is used to model the crystal size distribution. Then, the population balance equation (PBE) of the crystallization process is introduced into the network structure to construct a loss function, thereby imposing physical constraints on the process state. The effectiveness of PIKAN is demonstrated through a case study of the batch crystallization process. Summary of the Invention

[0007] To address the challenges of insufficient nonlinearity and interpretability in the modeling of crystal size distribution in existing batch crystallization processes, a KAN soft sensor model incorporating mechanism constraints is proposed to solve the challenges of nonlinearity and interpretability in the modeling of crystal size distribution in batch crystallization processes. First, KAN is used to model the crystal size distribution. Then, the population balance equation (PBE) of the crystallization process is introduced into the network structure to construct a loss function, thereby imposing physical constraints on the process state. In the case of batch crystallization processes, the effectiveness of the proposed PIKAN model is demonstrated by comparison with traditional methods.

[0008] The specific technical solutions are as follows:

[0009] A soft sensing method for crystal size distribution incorporating mechanism constraints, comprising the following steps:

[0010] 1) Acquisition and integration of data

[0011] First, CSD data under N different operating conditions are obtained as a data set from the batch crystallization process simulation data platform. The data set is divided into a training set, a validation set, and a test set, and the data set is normalized;

[0012] Step 1.1: Divide the collected CSD data into a training set, a validation set, and a test set according to a set ratio;

[0013] Step 1.2: Perform normalization processing on the training data set, the validation data set, and the test data set. The data standardization formula is

[0014]

[0015] where X' is the data after normalization processing, X is the original data without normalization processing, μ is the mean of X, and X max is the maximum value of X.

[0016] 2) Modeling training and model testing

[0017] Construct a soft sensing model for crystal size distribution incorporating mechanism constraints to train the training data set. The model structure of PIKAN is as Figure 1 shown. After training, predict the test data set;

[0018] Step 2.1: The crystal size distribution in the batch crystallization process changes with the operation. As Figure 1 shown, under the m-th operating condition, at time t corresponding to the k-th product region, the data in this region can be expressed as where K is the number of variables and T is the total number of time steps. Then, the data set of the distribution under the m-th operating condition is expressed as Where G represents the number of corresponding product regions. The KAN model encodes the data of adjacent regions. Under the m-th operating condition, the data p of the corresponding product region m,t,k After encoding, the prediction y m,t,k ' is obtained. The specific process is as follows:

[0019]

[0020] Where Loss α represents the loss in the KAN prediction process, and P represents the total number of initial operating conditions.

[0021] Step 2.2: Inject the kinetic knowledge of the batch crystallization process into the model iteration process. In the batch crystallization process, the changes between variables follow the population balance model (PBM), as follows:

[0022]

[0023] Where f(L, t) represents the crystal number density, which is the key prediction index in this study, t represents time, L represents crystal size, M represents growth rate, B represents nucleation rate. δ is the Dirac function, which is used to describe the minimum crystal size L0 required for crystal nucleation.

[0024] The above partial differential equation is complex to solve. To balance the simplicity of calculation and calculation efficiency, the forward Euler method is used to obtain the numerical solution of the partial differential equation. The forward Euler method approximately calculates the crystal size distribution y m,t+1,k value at t + 1 under the m-th operating condition:

[0025]

[0026] Where Δs is the time step, ΔL is the crystal size interval step. Assuming that nucleation occurs at size L0, the corresponding interval number is k0. When k = k0, δ = 1, otherwise δ = 0.

[0027] Using the y m,t,k and y m,t,k-1 predicted by the KAN network and substituting them into the forward Euler formula, the y m,t+1,k after physical constraints can be obtained. And the crystal size distribution y m,t+1,k ' at t + 1 under the m-th operating condition predicted by the KAN network. Thus, the physical constraint loss Loss p can be calculated:

[0028]

[0029] To further balance the calculation complexity and calculation efficiency, the key kinetics of the crystallization process can be represented by the moments of the required distribution. The j-th moment of the distribution is defined as:

[0030]

[0031] where L j represents the crystal size on the j-th moment, and f(L, t) is the crystal size function, representing the concentration of crystals with size L at time t. Therefore, the characteristics of CSD can be represented by the following set of moment equations:

[0032]

[0033] Specifically, μ0 is the zero-th moment, representing the total number of crystals with sizes in the range from 0 to L, μ1 is the first moment, representing the total specific length of the given crystals, μ2 is the second moment, which is proportional to the cumulative surface area of the crystals, and μ3 is the third moment, which is proportional to the cumulative volume of the crystals in the batch.

[0034] Therefore, the moments obtained from equations (8)-(9) can be used to estimate the global characteristics of the particle size distribution. The driving force during the crystallization process is the supersaturation of the solution, and the change in supersaturation during the crystallization process determines the final CSD. The relative supersaturation S is defined as follows:

[0035]

[0036]

[0037] where z n1 、z n2 、z n3 are empirical parameters describing the nucleation kinetics, and o n1 、o n2 、o n3 are crystal growth kinetics parameters. R is the molar gas constant. The change in solute concentration can be expressed by the following formula:

[0038]

[0039] where C s represents the solvent concentration in the solvent feed, ρ c is the crystal density, and o v is the shape factor. Here, according to common practice, the shape factor is regarded as a constant.

[0040] The solubility of the crystal is a function of temperature as follows:

[0041]

[0042] where T r represents the crystallizer temperature, and C sat represents the saturation concentration of the solute at temperature T r .

[0043] The heat transfer in the crystallizer can be described by the following ordinary differential equation:

[0044]

[0045] where T j is the jacket temperature, U is the overall heat transfer coefficient, A is the heat transfer area of the crystallizer, m r is the mass of the crystallization mixture, and c p,r is the specific heat capacity of the crystallization mixture.

[0046] According to Equation (9), we can obtain

[0047]

[0048] The preliminary prediction of f(L, t) is realized by using KAN. Therefore, at time t, f(L, t) under the m-th operating condition can be discretely expressed as At this time, by using Equations (18) to (21), the predicted values μ0(m, t)', μ1(m, t)', μ2(m, t)', and μ3(m, t)' of the zero-order moment, first-order moment, second-order moment, and third-order moment of the crystal size distribution under the m-th operating condition at time t can be obtained respectively, and are written as follows:

[0049] The preliminary prediction of f(L, t) is realized by using KAN. Therefore, at time t, f(L, t) under the m-th operating condition can be discretely expressed as At this time, by using Equations (18) to (21), the predicted values μ0(m, t)', μ1(m, t)', μ2(m, t)', and μ3(m, t)' of the zero-order moment, first-order moment, second-order moment, and third-order moment of the crystal size distribution under the m-th operating condition at time t can be obtained respectively, and are written as follows:

[0050] μ0(m, t)′ = y m,t,1 + ··· + y m,t,G (19)

[0051] μ1(m, t)′ = L1 × y m,t,1 + ··· + L G × y m,t,G (20)

[0052] μ2(m, t)′ = L1 2 × y m,t,1 + ··· + L G 2 × y m,t,G (21)

[0053] μ3(m, t)′ = L1 3 × y m,t,1 + ··· + LG 3 ×y m,t,G (22)

[0054] Among them, L1 represents the crystal size of the first product area, and L G represents the Gth product area. Use automatic differentiation to calculate dμ0(m,t)' / dt, dμ1(m,t)' / dt, dμ2(m,t)' / dt, dμ3(m,t)' / dt.

[0055] According to the zero-order moment of the distribution given by formula (8), formula (9) can obtain the derivatives of the first-order moment, second-order moment, and third-order moment of the distribution with respect to time:

[0056]

[0057]

[0058] Therefore, the loss functions under different operating conditions and at different times can be expressed as follows:

[0059]

[0060] Thus, the loss of the PIKAN model can be expressed by the following formula:

[0061]

[0062] Among them, λ i , i = 1, …, 8 represents the weight coefficient. Loss PIKAN includes the loss predicted by KAN and the physical constraint loss. Loss α represents the predicted loss of KAN, and Loss p emphasizes that the prediction should follow the knowledge of the kinetics of the crystallization process, represents the loss function used to limit the zero-order moment, to represents the loss functions used to balance the first-order moment to the third-order moment, represents the loss function used to balance the solute concentration, represents the loss function used to balance the crystallizer concentration.

[0063] Step 2.3: Input the test data set into the trained PIKAN model for testing, and at the same time use the trained comparison model to test the test data set;

[0064] The beneficial effects of the present invention are mainly manifested as follows: The present invention proposes a soft measurement method for crystal size distribution incorporating mechanism constraints, uses the Kolmogorov - Arnold network to achieve non - linear modeling of crystal size distribution, and injects the physical knowledge of crystallization process dynamics into the model in the form of a loss function. While improving the model prediction accuracy, it helps to enhance the transparency and interpretability of the model. BRIEF DESCRIPTION OF THE DRAWINGS

[0065] Figure 1 is the crystal size distribution diagram at the m - th operating condition and the t - th moment of the present invention; (a) 3 - D view of the crystal size distribution changing with time; (b) 2 - D view of the crystal size distribution at the t - th moment;

[0066] Figure 2 is the framework diagram of the soft measurement of crystal size distribution based on physical information constraints of the present invention;

[0067] Figure 3 is the diagram of the crystal size distribution predicted by different models of the present invention changing with time; (a) Truevalue; (b) PIKAN model; (c) KAN model; (d) LSTM model; (e) SVR model;

[0068] Figure 4 is the result diagram of different models predicted under different noise conditions of the present invention; (a) No noise; (b) 5% noise; (c) 10% noise. DETAILED DESCRIPTION OF THE INVENTION

[0069] The present invention will be further described below in conjunction with the drawings of the specification and embodiments, but the protection scope of the present invention is not limited thereto.

[0070] A soft measurement method for crystal size distribution incorporating mechanism constraints includes the following steps:

[0071] 1) Acquisition and integration of data

[0072] The batch intermittent crystallization simulation data platform is used to obtain the required CSD data. The simulation platform is built based on the dynamics of the intermittent crystallization process, and cooling crystallization is adopted, that is, the experimental data is obtained by controlling the reaction temperature to gradually decrease. First, the batch intermittent crystallization simulation data platform simulated the CSD data under 10 different operating conditions, and a total of 2800 groups of data were obtained. Every 280 groups represent a batch. 1680 groups of data in 6 batches were used as the training set, 560 groups of data in 2 batches were used as the validation set, and 560 groups of data in 2 batches were used as the test set;

[0073] Step 1.1: Divide the collected CSD data into a training data set, a validation data set, and a test data set according to the ratio of 5:3:2;

[0074] Step 1.2: Normalize the training dataset, validation dataset, and test dataset. The data standardization formula is shown in Formula (1).

[0075] 2) Conduct PIKAN model training

[0076] Step 2.1: The crystal size distribution data in the batch crystallization process is defined as where m is the number of operating conditions, K is the number of variables, and T is the total number of time steps. Then, the dataset of the distribution under the m-th operating condition is expressed as where G represents the number of corresponding product regions. The KAN model encodes the data of adjacent regions. Under the m-th operating condition, the data p of the corresponding product region m,t,k is encoded to obtain the prediction y m,t,k ', and the specific process is shown in Formulas (2) and (3);

[0077] Step 2.2: Inject the kinetic knowledge of the batch crystallization process into the model iteration process. In the batch crystallization process, the changes between variables follow the population balance model (PBM), as shown in Formulas (4) to (31);

[0078] Step 2.3: Input the test dataset into the trained PIKAN model for testing, and at the same time use the trained comparison model to test the test dataset.

[0079] Figure 3 shows the 3D predictions of the CSD change over time in the batch crystallization process by different models. It can be found that the proposed model's prediction of the crystal size change is closest to the real distribution. Table 1 and Figure 4 show the RMSE and R of the baseline model and PIKAN in the CSD prediction of the batch crystallization process 2 . In the prediction tasks of all baseline models, the RMSE of SVR is the highest, 172.3% higher than the proposed model, indicating the limitations of machine learning in non-linear modeling. The RMSE of LSTM is 146.8% higher than the proposed model. LSTM shows some improvements compared to SVR, benefiting from the effective extraction of time dependence. KAN is 21.3% higher than the proposed model. KAN improves the accuracy by more precisely characterizing the non-linear relationship. In addition, by injecting the physical knowledge of the batch crystallization process into the model, PIKAN further improves the accuracy of the KAN model in predicting the crystal size distribution.

[0080]

[0081] In addition, to demonstrate the robustness of the model, 5% and 10% noise was added to the CSD data of the batch crystallization process. As shown in Table 2, in the 5% noise environment, the RMSE of PIKAN was 28.6% - 104.8% lower than that of the baseline model. In the 10% noise environment, the RMSE of PIKAN was 36.9% - 101.5% lower than that of the baseline model. This benefits from the advantages of KAN in nonlinear modeling and the guiding role of physical constraints in model prediction. However, neither of them is specifically tailored for noise. In short, noise will not seriously affect the performance of the model.

[0082] Table 2 RMSE of Four Soft Sensor Models under 5% and 10% Noise

[0083]

[0084] The soft sensor method for crystal size distribution incorporating mechanism constraints proposed in the present invention has been fully verified on the crystal size distribution dataset of the batch crystallization process, improving the accuracy of crystal size distribution prediction in the batch crystallization process, and having universality and generality.

[0085] The content described in the embodiments of this specification is only an enumeration of the implementation forms of the inventive concept. The protection scope of the present invention should not be regarded as limited to the specific forms stated in the embodiments. The protection scope of the present invention also covers equivalent technical means that can be conceived by those skilled in the art based on the inventive concept.

Claims

1. A soft measurement method for crystal size distribution incorporating mechanism constraints, characterized in that, It includes the following steps: 1) Data acquisition and integration: First, obtain the CSD data under N different working conditions as a data set from the batch crystallization process simulation data platform, divide the data set into a training set, a validation set, and a test set, and normalize the data set; 2) Model training and model testing: Construct a soft sensor model PIKAN for crystal size distribution based on incorporating mechanism constraints, train it with the training set data, and predict the test data set after training.

2. The soft measurement method for crystal size distribution incorporating mechanism constraints according to claim 1, characterized in that, The specific process of step 1) includes the following steps: Step 1.1: Divide the collected CSD data into a training set, a validation set, and a test set according to a ratio; Step 1.2: Perform normalization processing on the training data set, the validation data set, and the test data set. The data standardization formula is Among them, X' is the data after standardization processing, X is the original data without standardization processing, μ is the mean of X, and X max is the maximum value of X.

3. The soft measurement method for crystal size distribution incorporating mechanism constraints as described in claim 2, characterized in that The specific process of step 2) includes the following steps: Step 2.1: The crystal size distribution in the batch crystallization process changes with the operation. Under the m-th operating condition, at time t corresponding to the k-th product region, the data in this region is expressed as where K is the number of variables and T is the total number of time steps; the data set of the distribution under the m-th operating condition is expressed as where G represents the number of corresponding product regions, the KAN model encodes the data of adjacent regions, and under the m-th operating condition, the data p of the corresponding product region m,t,k after encoding, the prediction y m,t,k ' is obtained. The specific process is as follows: Among them, Loss α represents the loss during the KAN prediction process, and P represents the total number of initial operating conditions; Step 2.2: Inject the kinetic knowledge of the batch crystallization process into the model iteration process. In the batch crystallization process, the changes between variables follow the population balance model PBM, as follows: Among them, f(L,t) represents the crystal number density, t represents time, L represents crystal size, M represents growth rate, B represents nucleation rate, δ is the Dirac function, used to describe the minimum crystal size L0 required for crystal nucleation; To balance the simplicity of calculation and computational efficiency, the forward Euler method is adopted to obtain the numerical solution of the partial differential equation. The forward Euler method approximately calculates the crystal size distribution y m,t+1,k value at time t + 1 under the m-th operating condition: Among them, Δs is the time step, ΔL is the crystal size interval step. Assume that nucleation occurs at size L0, and the corresponding interval number is k0. When k = k0, δ = 1, otherwise δ = 0; y predicted using the KAN network m,t,k and y m,t,k-1 , substituting into the forward Euler formula to obtain y after physical constraints m,t+1,k , while the KAN network predicts the crystal size distribution y m,t+1,k ' at time t + 1 under the m-th operating condition, and the physical constraint loss Loss p is calculated therefrom as follows: To further balance the computational complexity and computational efficiency, the kinetics of the crystallization process is represented by the moments of the desired distribution. The jth moment of the distribution is defined as: where L j represents the crystal size on the j-th moment, and f(L, t) is the crystal size function, representing the concentration of crystal size L at time t; therefore, the characteristics of CSD are represented by the following set of moment equations: μ0 is the zero-order moment, representing the total number of crystals with sizes ranging from 0 to L. μ1 is the first-order moment, representing the specific total length of the crystals. μ2 is the second-order moment, which is proportional to the cumulative surface area of the crystals. μ3 is the third-order moment, which is proportional to the cumulative volume of the crystals in the batch; The moments obtained from equations (8) and (9) can be used to estimate the global characteristics of the particle size distribution; the driving force during the crystallization process is the supersaturation of the solution, and the change of supersaturation during the crystallization process determines the final CSD; the relative supersaturation S is defined as follows: where z n1 , z n2 , z n3 are empirical parameters describing nucleation kinetics, o n1 , o n2 , o n3 are crystal growth kinetics parameters, and R is the molar gas constant; The change of solute concentration is expressed by the following formula: where C s represents the solvent concentration in the solvent feed, ρ c is the crystal density, and o v is the shape factor; The solubility of the crystal is a function of temperature, as follows: where T r represents the temperature of the crystallization agent, in °C sat represents the saturation concentration of the solute at temperature T r in g / 100 g of solvent; The heat transfer in the crystallizer is described by the following ordinary differential equation: Among them, T j is the jacket temperature, U is the overall heat transfer coefficient, A is the heat transfer area of the crystallizer, m r is the mass of the crystallization mixture, c p,r is the specific heat capacity of the crystallization mixture; According to formula (9), we get: The preliminary prediction of f(L,t) is realized by using KAN. Therefore, at time t, f(L,t) under the m-th operating condition is discretely expressed as At this time, the predicted values μ0(m,t)', μ1(m,t)', μ2(m,t)', and μ3(m,t)' of the zero-order moment, first-order moment, second-order moment, and third-order moment of the crystal size distribution under the m-th operating condition at time t are obtained by using formulas (18) to (21) respectively, and are written as follows: μ0(m,t)′ = y m,t,1 + ··· + y m,t,G (19) μ1(m,t)′ = L1 × y m,t,1 + ··· + L G × y m,t,G (20) μ2(m,t)′ = L1 2 ×y m,t,1 +···+L G 2 ×y m,t,G (21) μ3(m,t)′ = L1 3 × y m,t,1 + ··· + L G 3 × y m,t,G (22) Among them, L1 represents the crystal size of the first product area, and L G represents the Gth product area. Use automatic differentiation to calculate dμ0(m,t)' / dt, dμ1(m,t)' / dt, dμ2(m,t)' / dt, dμ3(m,t)' / dt; According to the zero-order moment of the distribution given by formula (8), the derivatives of the first-order moment, second-order moment, and third-order moment of the distribution with respect to time are obtained from formula (9): Therefore, the loss function under different operating conditions and at different times is expressed as follows: Thus, the loss of the PIKAN model is expressed by the following formula: Among them, λ i , i = 1, …, 8 represents the weight coefficient, Loss PIKAN includes the loss predicted using KAN and the physical constraint loss, Loss α represents the predicted loss of KAN, Loss p emphasizes that the prediction should follow the knowledge of the kinetics of the crystallization process, represents the loss function for restricting the zero-order moment, to represents the loss function for balancing the first-order moment to the third-order moment, represents the loss function for balancing the solute concentration, represents the loss function for balancing the crystallizer concentration; Step 2.3: Input the test set data into the trained PIKAN model for testing, and at the same time use the trained comparison model to test the test set data.

Images