Partial differential equation inverse problem solving method under sparse data
By using a combination of KNN compensation and deep neural network under sparse data conditions, the weights are dynamically adjusted and the parameter identification process is optimized, which solves the problem of insufficient recognition ability of PINN in sparse data and noise environments, and achieves high-precision and anti-noise parameter identification effect.
Patent Information
- Application Number
- CN202510159503.1
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-02-13
- Publication Date
- 2025-06-03
AI Technical Summary
Under sparse data conditions, existing PINNs show insufficient in the inverse problem of partial differential equations, especially when the observation data cannot fully cover the critical areas of parameter changes or capture the multi-scale characteristics of physical quantities, it cannot effectively identify unknown parameters and exhibit insufficient noise immunity when noise exists.
A method for solving the inverse problem of partial differential equations under sparse data is proposed. By collecting multiple observation data points, using KNN compensation technology to generate compensation data points, combining deep neural networks and dynamic weight adjustment mechanisms, the identification process of unknown parameters is optimized, and the stability and efficiency of the model are improved through priority sorting and gradual freezing parameter strategies.
High-precision parameter identification is achieved with few observation data points, and has strong noise resistance, which can effectively reduce noise interference to the identification results and improve the robustness of unknown parameter identification.
Smart Images

Figure CN120086479A_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to the technical fields of computational mathematics and artificial intelligence, and specifically to a method for solving inverse problems of partial differential equations under sparse data. Background Art
[0002] Partial differential equations (PDEs) can describe many complex phenomena in nature and engineering, including heat conduction, fluid mechanics, electromagnetic fields, elasticity, quantum mechanics, etc. Through PDEs, the laws of change of substances can be simulated and predicted, such as temperature distribution, pressure change, velocity field, wave propagation, etc. The purpose of the inverse problem of PDEs is to infer unknown parameters or source terms from limited observed data, which is one of the core tasks in engineering and scientific computing. For example, the temperature distribution can be used to infer the thermal conductivity, the relationship between pressure and flow velocity can be used to invert the viscosity of the fluid, and the measurement of the electric field and current density can be used to calculate the conductivity of the medium. However, sparse data is a common phenomenon in many inverse problems of PDEs, and there are great challenges in making full use of sparse data and combining physical constraints to improve the identification accuracy of unknown parameters.
[0003] Physics-informed neural networks (PINNs) can effectively handle high-dimensional, non-linear, and multi-scale problems by embedding physical equations into neural networks, and PINNs have become a frontier field for solving complex PDE problems. However, when PINNs are applied to the field of inverse problems of partial differential equations, especially in the case of sparse observed data, they show some significant deficiencies. Under sparse data conditions, the observation points are often limited in distribution and may not be able to fully cover the key regions of parameter changes or capture the multi-scale characteristics of physical quantities. This results in weak constraints on the physical laws in unobserved regions during the learning process of existing PINNs, and the unknown parameters cannot be effectively identified. When the observed data is sparse and accompanied by noise, the performance of PINNs further deteriorates. Since PINNs emphasize balancing data error and physical constraints through the loss function, noise will lead to the amplification of data error loss, making the model more inclined to fit the noisy data, thus destroying physical consistency. In addition, the existing PINN improvement methods have limited ability to cope with noise perturbations and are difficult to effectively distinguish noise from real physical signals. Summary of the Invention
[0004] Aiming at the deficiencies of the existing technology, the purpose of the present invention is to propose a method for solving inverse problems of partial differential equations under sparse data, including:
[0005] Step 1: Collect data of multiple observation data points in the actual industrial process;
[0006] Among them, the data of the observation data points includes the coordinate set of all observation data points and the corresponding solution set U data , the coordinate set of all observation data points A coordinate set including at least multiple initial data points and a coordinate set of N u data points inside the solution domain The coordinate set of the multiple initial data points includes the coordinates of the time and space distributions of all initial data points. The coordinate set of N u data points inside the solution domain includes the coordinates of the time and space distributions of N u data points inside the solution domain;
[0007] Step 2: Perform KNN compensation based on the coordinate set of N u data points inside the solution domain to obtain a compensated data point set
[0008] Step 3: Obtain the partial differential equation and the initial condition equation of the actual industrial process. The partial differential equation contains multiple unknown parameters. According to the partial differential equation and the initial condition equation, sort the priorities of the unknown parameters. Among them, the unknown parameters in the initial condition equation have the highest priority. For the partial differential equation, for the terms containing independent variables, the fewer the number of independent variables, the higher the priority of the corresponding unknown parameter. In the case of the same number of independent variables, the unknown parameter corresponding to the term with a lower derivative order has a higher priority. In the case of the same derivative order, the unknown parameter corresponding to the term with a lower degree of nonlinearity has a higher priority; Initialize all unknown parameters to obtain the initial values of the unknown parameters, set the initial number of iterations, use the initial number of iterations as the current number of iterations, and use the initial values of the unknown parameters as the unknown parameters of the current number of iterations;
[0009] Step 4: Calculate the physical loss loss based on the compensated data point set phsics ;
[0010] Step 5: Calculate the initial condition loss loss based on the coordinate set of multiple initial data points ic ;
[0011] Step 6: Calculate the data loss loss and the corresponding solution set U data based on the coordinate set of all observed data points data ;
[0012] Step 7: Set the loss weight w phsics of the physical loss loss phsics, set the initial condition loss ic of the loss weight w ic , set the data loss data of the loss weight w data ;
[0013] Step 8: According to the loss weight w phsics , the loss weight w ic and the loss weight w data , perform a weighted sum on the physical loss phsics , the initial condition loss ic and the data loss data to obtain the final loss value loss;
[0014] Step 9: According to the final loss value loss, update the weight matrix and bias vector in the deep neural network, and at the same time update the unknown parameters of the current iteration number to obtain the updated unknown parameters;
[0015] Step 10: According to the priority of the unknown parameters, starting from the unknown parameter with the highest priority, for this unknown parameter, judge whether this unknown parameter satisfies the convergence condition. If this unknown parameter does not satisfy the convergence condition, increment the current iteration number by one, use it as the new current iteration number, use the updated unknown parameter as the unknown parameter of the new current iteration number, and return to execute Step 4;
[0016] If this unknown parameter satisfies the convergence condition, use the unknown parameter of the current iteration number as the final value of this unknown parameter, increment the current iteration number by one, use it as the new current iteration number, keep the value of this unknown parameter unchanged, for the unknown parameters other than this unknown parameter, use the updated value as the unknown parameter of the new current iteration number, return to execute Step 4, the unknown parameters that are not updated and remain unchanged in Step 9, and in Step 10, judge whether the unknown parameter of the next priority of this unknown parameter satisfies the convergence condition, so as to obtain the values of all unknown parameters.
[0017] Optionally, Step 2 specifically includes:
[0018] Step 2.1: Obtain the coordinates of the time and space distributions of all data points on the solution domain to obtain the first set X star , the first set X star includes the coordinate set of N u data points inside the solution domain wherein, the coordinate set of N u data points inside the solution domain is expressed as:
[0019] Step 2.2: For the N u coordinate sets of data points inside the solution domain calculate the time and space coordinates x of each data point in i it, and calculate the Euclidean distance between x i and the time and space coordinates of each data point in the first set X star specifically calculated by the following formula:
[0020]
[0021] where d represents the dimension of space, y m represents the time and space coordinates of the data points in the first set X star , d(x i , y m ) represents the Euclidean distance between x i and y m , represents the coordinate value of x i in the p-th dimension, represents the coordinate value of y m in the p-th dimension;
[0022] Step 2.3: Sort the Euclidean distances between the time and space coordinates x of the data points and the time and space coordinates of each data point in the first set X i in ascending order, and obtain the time and space coordinates of the data points in the first set X corresponding to the first k + 1 Euclidean distances, getting the initial neighbor points N star (x star ) set. Remove the time and space coordinates of the data point x k from the initial neighbor point set N i (x k ) to obtain the nearest neighbor point set N' i (x i ); k (x i );
[0023] Step 2.4: Add random noise to the time and space coordinates y of each nearest neighbor point in the nearest neighbor point set N' k (x i ) to obtain the perturbed time and space coordinates of the nearest neighbor points i specifically implemented by the following formula:
[0024]
[0025] where r is a random vector sampled from a uniform distribution with a range of [-noise_scale, noise_scale], and noise_scale is the amplitude of the random noise;
[0026] Step 2.5: Determine the temporal and spatial coordinates of each perturbed nearest neighbor point to check if they satisfy the boundary conditions, which is specifically expressed by the following formula:
[0027]
[0028] where lb is the lower boundary of the solution domain, ub is the upper boundary of the solution domain, and min_distance is the minimum distance to the boundary;
[0029] For a perturbed nearest neighbor point that satisfies the boundary conditions, it is characterized as a valid point. The temporal and spatial distribution coordinates of all valid points form the set of valid neighbors V(x i );
[0030] Step 2.6: Merge the sets of valid neighbors for each data point in the coordinate set u of the N data points inside the solution domain to obtain the final set of sampled points X final , which is specifically achieved by the following formula:
[0031]
[0032] Step 2.7: Merge the final set of sampled points X final with the coordinate set u of the N data points inside the solution domain to obtain the set of compensated data points , which is specifically achieved by the following formula:
[0033]
[0034] Optionally, Step 4 is specifically achieved by the following formula:
[0035]
[0036] where N represents the number of data points in the set of compensated data points , (x r , t r ) are the temporal and spatial coordinates of the data points in the set of compensated data points , t r is the time variable, x r is the spatial variable, f(x r , tr ) is the residual of the partial differential equation, which is specifically expressed by the following formula:
[0037]
[0038] where F(·) is the partial differential operator, and u(x r , t r ) is the formula under the unknown parameters at the current iteration number. represents the derivative with respect to x r and t r . s and h respectively represent the orders of differentiation with respect to x r and t r . Both s and h are non - negative integers, and λ represents the unknown parameter at the current iteration number.
[0039] Optionally, step 5 specifically includes:
[0040] Input the coordinate sets of multiple initial data points into the deep neural network to obtain the predicted solutions of the coordinate sets of multiple initial data points , and then calculate the initial condition loss loss ic , which is specifically implemented by the following formula:
[0041]
[0042] where P represents the number of initial data points in the coordinate set . (x ic , 0) is the coordinate in the coordinate sets of multiple initial data points , u ic_pred (x ic , 0) is the predicted solution of the coordinate sets of multiple initial data points , and u ic_true (x ic , 0) is the theoretical solution calculated for (x ic , 0) based on the initial condition equation and the unknown parameters at the current iteration number.
[0043] Optionally, step 6 specifically includes:
[0044] Input the coordinate sets of all observed data points into the deep neural network to obtain the predicted solutions of the coordinate sets of all observed data points , and then calculate the data loss according to the predicted solutions of the coordinate sets of all observed data points and the corresponding solution set U data , which is specifically implemented by the following formula:
[0045]
[0046] Among them, M is the coordinate set of all observed data points The number of observed data points in, (x d , t d ) is the coordinate set of all observed data points The coordinates of the observed data points in, u true (x d , t d ) is the true solution corresponding to (x data , t d ) in U d ), u pred (x d , t d ) is the coordinate set of all observed data points The predicted solution.
[0047] Optionally, the deep neural network includes an input layer, a hidden layer, and an output layer. After the data is input into the deep neural network, it passes through the input layer to the first layer and is calculated by the following formula:
[0048] z (1) =σ(w (0) x + b (0) );
[0049] Among them, z (1) is the output of the first layer, w (0) is the weight matrix of the 0th layer, b (0) is the bias vector of the 0th layer, σ represents the activation function, and x is the coordinates of time and space input into the deep neural network;
[0050] Furthermore, the output of the first layer is calculated layer by layer through forward propagation. Specifically, the calculation from the, layer to the l+1 layer is represented by the following formula:
[0051] z (l+1) =σ (l) (w (l) z (l) + b (l) ), l = 1, 2,..., L-2;
[0052] Among them, z (l+1) is the output of the l+1 layer, σ (l) is the activation function of the, layer, w (l) is the weight matrix of the, layer, z (l) is the output of the, layer, b (l) is the bias vector of the, layer;
[0053] Furthermore, in the output layer, it is calculated by the following formula:
[0054] z(L) = w (L-1) z (L-1) + b (L-1) ;
[0055] where z (L) is the output of the L-th layer, i.e., the output of the deep neural network, w (L-1) is the weight matrix of the (L - 1)-th layer, z (L-1) is the output of the (L - 1)-th layer, b (L-1) is the bias vector of the (L - 1)-th layer.
[0056] Optionally, the parameters of the input layer, hidden layer, and output layer of the deep neural network are pre-set, and the parameters of the input layer, hidden layer, and output layer are layers = [n 0 , n 1 , n 2 ,..., n L , where n 0 is the number of neurons in the input layer, n L is the number of neurons in the output layer, n 1 , n 2 ,..., n L-1 are the number of neurons in the hidden layer.
[0057] Optionally, step 7 specifically includes:
[0058] Set the loss weight w phsics and the initial condition loss loss ic based on the magnitude of the physical loss loss phsics and the initial condition loss loss ic , which is specifically implemented through the following formula:
[0059]
[0060]
[0061] Calculate the loss weight w data according to the data loss loss data , which is specifically calculated through the following formula:
[0062]
[0063] where u pred is the set of coordinates of all observed data points of the prediction solution, is the gradient norm of u pred , which is specifically calculated through the following formula:
[0064]
[0065] Optionally, step 8 is specifically represented by the following formula:
[0066] loss = w physics ·loss physics + w ic ·loss ic + w data ·loss data .
[0067] Optionally, determining whether the unknown parameter satisfies the convergence condition in step 10 specifically includes:
[0068] Calculating the relative change rate θ between the unknown parameter of the current iteration and the unknown parameter of the previous iteration Relative_Change , which is specifically implemented by the following formula:
[0069]
[0070] where θ new is the unknown parameter of the current iteration, and θ old is the unknown parameter of the previous iteration;
[0071] Judging whether the relative change rate θ Relative_Change is less than the expected threshold θ ε . When the relative change rate θ Relative_Change is not less than the expected threshold θ ε , it indicates that the unknown parameter does not satisfy the convergence condition;
[0072] When the relative change rate θ Relative_Change is less than the expected threshold θ ε , the convergence count θ count is incremented by one, and it is judged whether the incremented convergence count reaches the expected convergence count threshold θ Design_Count . When the incremented convergence count does not reach the expected convergence count threshold θ Design_Count , it indicates that the unknown parameter does not satisfy the convergence condition;
[0073] When the incremented convergence count reaches the expected convergence count threshold θ Design_Count , it indicates that the unknown parameter satisfies the convergence condition.
[0074] The beneficial effects of adopting the above technical solutions are as follows:
[0075] 1. Low observation data requirements: High-precision parameter identification can still be achieved with fewer observation data points.
[0076] 2. Excellent anti-noise performance: Even when there is noise interference in the observation data, the parameter identification accuracy still has strong robustness.
[0077] 3. Wide applicability: Applicable to a variety of partial differential equation models, especially complex physical scenarios and multi-scale problems.
[0078] 4. Innovative optimization design: By strategies such as dynamically adjusting weights, compensating sampling points inside the solution domain, and gradually freezing parameters, the stability and efficiency of the model are improved. Description of the Drawings
[0079] Figure 1 It is a schematic flow chart of a method for solving the inverse problem of partial differential equations under sparse data in an embodiment of the present invention;
[0080] Figure 2 It is a schematic diagram of the weight distribution during the training process in an embodiment of the present invention;
[0081] Figure 3 It is for the unknown parameters β, λ 1 , λ 2 during the training process in an embodiment of the present invention; it is a schematic diagram of the convergence situation;
[0082] Figure 4 It is for the unknown parameter λ 3 during the training process in an embodiment of the present invention; it is a schematic diagram of the convergence situation. Detailed Embodiments
[0083] The following combines the drawings and embodiments to further describe in detail the specific embodiments of the present invention. The following embodiments are used to illustrate the present invention, but are not used to limit the scope of the present invention.
[0084] Aiming at the problems existing in the prior art, the present invention provides a method for solving the inverse problem of partial differential equations under sparse data, which can identify the unknown parameters contained in the PDE and the initial conditions with fewer observed data points. When the observed data contains noise, the identification accuracy of the unknown parameters also has good robustness. Specifically, the present invention provides a method for solving the inverse problem of partial differential equations under sparse data, combined with Figure 1 , it may include the following steps:
[0085] Step 1: Collect the data of multiple observed data points in the actual industrial process;
[0086] Among them, the data of the observed data points includes the coordinate set of all observed data points and the corresponding solution set U data , the coordinate set of all observed data points at least includes the coordinate set of multiple initial data points and the coordinate set of N u data points inside the solution domain The coordinate set of the multiple initial data points includes the coordinates of the time and space distributions of all the initial data points, and the coordinate set of N u data points inside the solution domain includes the coordinates of the time and space distributions of N u data points inside the solution domain;
[0087] wherein, the corresponding solution set U data at least includes the coordinate set of the initial data points corresponding solution set U initial , and the coordinate set of N u data points inside the solution domain corresponding solution set U interior ;
[0088] Furthermore, it may further include the coordinate set of multiple boundary data points including the coordinates of the time and space distributions of all the boundary data points, corresponding solution set U data may further include corresponding solution set U boundry , and the present invention is processed and calculated through U data . Furthermore, it may combine corresponding solution set U initial , corresponding solution set U interior , the coordinate set of multiple boundary data points and corresponding solution set U boundry for processing and calculation to perform more precise processing and calculation and obtain more precise results.
[0089] Step 2: Perform KNN compensation based on the coordinate set of N u data points inside the solution domain to obtain a compensated data point set
[0090] Step 2.1: Obtain the coordinates of the time and space distributions of all the data points on the solution domain to obtain a first set X star , and the first set X star includes the coordinate set of N u data points inside the solution domain wherein, the first set X star is obtained according to the set minimum sampling time and minimum sampling space, and the first set X starIt is understood as the KNN candidate sampling point set.
[0091] Among them, the coordinate set of N u data points inside the solution domain is expressed as:
[0092] Step 2.2: For the coordinates x u of each data point in the coordinate set of N data points inside the solution domain, calculate the Euclidean distance between the time and space distribution coordinates x i and the time and space distribution coordinates of each data point in the first set X i , specifically calculated by the following formula: star
[0093]
[0094]
[0094] where d represents the dimension of space, y m represents the time and space distribution coordinates of the data points in the first set X star , d(x i , y m ) represents the Euclidean distance between x i and y m , represents the coordinate value of x i in the p-th dimension, represents the coordinate value of y m in the p-th dimension;
[0095] Step 2.3: Sort the Euclidean distances between the time and space distribution coordinates x i of the data points and the time and space distribution coordinates of each data point in the first set X star from smallest to largest, and obtain the time and space distribution coordinates of the data points in the first set X star corresponding to the first k + 1 Euclidean distances, to obtain the initial neighbor set N k (x i ). In the initial neighbor set N k (x i ), remove the time and space distribution coordinates of the data point x i to obtain the nearest neighbor set N k '(x i );
[0096] where k = [10, 15, 20, 25, 30]. If the number is large, the number of k can be reduced. If the number is small, the number of k can be increased.
[0097] Step 2.4: For the set of nearest neighbor points N k '(x i ) in each of the coordinates y of the temporal and spatial distributions of the nearest neighbor points i Add random noise to obtain the coordinates of the perturbed temporal and spatial distributions of the nearest neighbor points Specifically, it is implemented through the following formula:
[0098]
[0099] where r is a random vector sampled from a uniform distribution, with a range of [-noise_scale, noise_scale], and noise_scale is the amplitude of the random noise. In a specific implementation, noise_scale = [0.01, 0.02,..., 0.1]. Among them, if the range of the solution domain is large, the value can be increased.
[0100] Step 2.5: Determine whether the coordinates of the temporal and spatial distributions of each perturbed nearest neighbor point satisfy the boundary conditions, specifically expressed by the following formula:
[0101]
[0102] where lb is the lower boundary of the solution domain, ub is the upper boundary of the solution domain, and min_distance is the minimum distance close to the boundary. In a specific implementation, min_distance = [0.01, 0.02,..., 0.1].
[0103] In the case where the perturbed nearest neighbor points satisfy the boundary conditions, it is characterized that the perturbed nearest neighbor points are valid points, and the coordinates of the temporal and spatial distributions of all valid points form the set of valid neighbors V(x i );
[0104] Step 2.6: Merge the sets of valid neighbors of each data point in the coordinate set of N u data points inside the solution domain to obtain the final sampling point set X final , specifically implemented through the following formula:
[0105]
[0106] Step 2.7: Merge the final sampling point set X final with the coordinate set of N u data points inside the solution domain to obtain the compensated data point set Specifically, it is implemented through the following formula:
[0107]
[0108] Step 3: Obtain the partial differential equation and the initial condition equation of the actual industrial process. The partial differential equation contains multiple unknown parameters. According to the partial differential equation and the initial condition equation, sort the priorities of the unknown parameters. Among them, the unknown parameters in the initial condition equation have the highest priority. For the partial differential equation, for the terms containing independent variables, the fewer the number of independent variables, the higher the priority of the corresponding unknown parameter. In the case of the same number of independent variables, the unknown parameter corresponding to the term with a lower derivative order has a higher priority. In the case of the same derivative order, the unknown parameter corresponding to the term with a lower degree of non-linearity has a higher priority. Based on this, the priorities of each term can be expressed as:
[0109] u > u 2 > u 3 > u x > uu x > u 2 u x > u xx > uu xx > u 2 u xx > u xxx > uu xxx > u 2 u xxx …;
[0110] Among them, u x represents the first derivative of u, u xx represents the second derivative of u, and so on.
[0111] If the partial differential equation is:
[0112]
[0113] Then the priority of the unknown parameters among them is:
[0114] λ 1 > λ 2 > λ 3 > λ 4 > λ 5 > λ 6 > λ 7 > λ 8 > λ 9 > λ 10 > λ 11 > λ 12 > …;
[0115] If the unknown parameter in the initial condition equation is β, then the priorities of all unknown parameters are expressed as:
[0116] β > λ1 > λ 2 > λ 3 > λ 4 > λ 5 > λ 6 > λ 7 > λ 8 > λ 9 > λ 10 > λ 11 > λ 12 > …;
[0117] Initialize all unknown parameters to obtain the initial values of the unknown parameters, set the initial iteration count, use the initial iteration count as the current iteration count, and use the initial values of the unknown parameters as the unknown parameters for the current iteration count;
[0118] Among them, steps 1 and 2 can be executed in parallel with step 3 to determine the priority of the unknown parameters.
[0119] Step 4: According to the set of compensated data points the partial differential equation, and the unknown parameters for the current iteration count, calculate the physical loss loss phsics , which is specifically implemented through the following formula:
[0120]
[0121] where N represents the number of data points in the set of compensated data points , (x r , t r ) are the coordinates of the time and space distribution of the data points in the set of compensated data points , t r is the time variable, x r is the space variable, f(x r , t r ) is the partial differential equation residual, which is specifically represented by the following formula:
[0122]
[0123] where F(·) is the partial differential operator, u(x r , t r ) is the formula under the unknown parameters for the current iteration count, represents the derivative with respect to x r and t r , s and h respectively represent the number of derivative times with respect to x r and t r , s and h are both non-negative integers, and λ represents the unknown parameters for the current iteration count.
[0124] It should be noted that λ in the formula represents the unknown parameters included in the formula. In the specific implementation process, the number of unknown parameters can be multiple.
[0125] Step 5: Based on the coordinate sets of multiple initial data points Calculate the initial condition loss loss based on the initial condition equation and the unknown parameters of the current iteration ic , specifically including:
[0126] Input the coordinate sets of multiple initial data points into the deep neural network to obtain the predicted solutions of the coordinate sets of multiple initial data points , and then calculate the initial condition loss loss ic , which is specifically implemented through the following formula:
[0127]
[0128] where P represents the number of initial data points in the coordinate set , (x ic ,0) are the coordinates of the coordinate sets of multiple initial data points in, u ic_pred (x ic ,0) are the predicted solutions of the coordinate sets of multiple initial data points , u ic_true (x ic ,0) is the theoretical solution calculated for (x ic ,0) based on the initial condition equation and the unknown parameters of the current iteration, and can be expressed by the following formula:
[0129] u ic_true (x i ,0) = g(x i ,β)
[0130] where g(x i ,β) is the initial condition equation, and β is the unknown parameter contained in g(x i ,β).
[0131] Step 6: Calculate the data loss loss based on the coordinate sets of all observed data points and the corresponding solution set U data , data specifically including:
[0132] Input the coordinate sets of all observed data points into the deep neural network to obtain the predicted solutions of the coordinate sets of all observed data points , and then, according to the coordinate sets of all observed data points The predicted solution and the corresponding solution set U data , calculate the data loss, which is specifically implemented through the following formula:
[0133]
[0134] where M is the coordinate set of all observed data points in the number of observed data points, (x d , t d ) is the coordinate set of all observed data points in the coordinates of the observed data point, u true (x d , t d ) is for U data in the true solution corresponding to (x d , t d ), u pred (x d , t d ) is the coordinate set of all observed data points of the predicted solution.
[0135] wherein, the deep neural network includes an input layer, a hidden layer and an output layer. After the data is input into the deep neural network, it passes through the input layer to the first layer and is calculated through the following formula:
[0136] z (1) =σ(w (0) x + b (0) );
[0137] where z( 1 ) is the output of the first layer, w( 0 ) is the weight matrix of the 0th layer, b( 0 ) is the bias vector of the 0th layer, σ represents the activation function, and x is the coordinates of the time and space distribution input into the deep neural network;
[0138] Furthermore, the output of the first layer is calculated layer by layer through forward propagation. Specifically, the calculation from the, layer to the l + 1 layer is represented by the following formula:
[0139] z (l+1) =σ (l) (w (l) z (l) + b (l) ), l = 1, 2,..., L - 2;
[0140] where z (l+1 ) is the output of the l + 1 layer, σ (l) is the activation function of the, layer, w (l) is the weight matrix of the, layer, z (l)is the output of the l-th layer, b(l ) is the bias vector of the l-th layer;
[0141] Furthermore, in the output layer, it is calculated by the following formula:
[0142] z (L) = w (L-1) z (L-1) + b (L-1) ;
[0143] where z (L) is the output of the L-th layer, i.e., the output of the deep neural network, w (L-1) is the weight matrix of the (L - 1)-th layer, z (L-1) is the output of the (L - 1)-th layer, b (L-1) is the bias vector of the (L - 1)-th layer.
[0144] Among them, the parameters of the input layer, hidden layer, and output layer of the deep neural network are pre-set, and the parameters of the input layer, hidden layer, and output layer are layers = [n 0 , n 1 , n 2 ,..., n L , where n 0 is the number of neurons in the input layer, n L is the number of neurons in the output layer, n 1 , n 2 ,..., n L-1 are the number of neurons in the hidden layer, and the deep neural network has a total of L layers.
[0145] Step 7: Set the loss weight W phsics of the physical loss loss phsics , set the loss weight w ic of the initial condition loss loss ic , set the loss weight W data of the data loss loss data , specifically including:
[0146] Based on the size of the initial condition loss loss phsics of the physical loss loss ic to set the loss weights W phsics and the initial condition loss loss ic , specifically implemented by the following formula:
[0147]
[0148]
[0149] According to the data loss loss data, calculate the loss weight W data , which is specifically calculated by the following formula:
[0150]
[0151] where u pred is the set of coordinates of all observed data points of the predicted solution, and pred is the gradient norm of u
[0152]
[0153] Step 8: According to the loss weight W phsics , the loss weight w ic and the loss weight W data , perform a weighted sum on the physical loss loss phsics , the initial condition loss loss ic and the data loss loss data to obtain the final loss value loss, which is specifically expressed by the following formula:
[0154] loss = w physics · loss physics + w ic · loss ic + w data · loss data .
[0155] Step 9: According to the final loss value loss, update the weight matrix and bias vector in the deep neural network, and at the same time update the unknown parameters of the current iteration number to obtain the updated unknown parameters;
[0156] Step 10: According to the priority of the unknown parameters, starting from the unknown parameter with the highest priority, for this unknown parameter, judge whether this unknown parameter meets the convergence condition. If this unknown parameter does not meet the convergence condition, increment the current iteration number by one, use it as the new current iteration number, and use the updated unknown parameter as the unknown parameter of the new current iteration number, and return to execute Step 4;
[0157] When the unknown parameter satisfies the convergence condition, the unknown parameter of the current iteration is taken as the final value of the unknown parameter, the current iteration number is incremented by one, which is taken as the new current iteration number, and the value of the unknown parameter remains unchanged. For the unknown parameters other than this unknown parameter, the updated value is taken as the unknown parameter of the new current iteration number, and step 4 is returned. The unknown parameter that remains unchanged without update in step 9 is judged in step 10 whether the unknown parameter of the next priority satisfies the convergence condition, so as to obtain the values of all unknown parameters. That is to say, if the unknown parameters are β, λ1, and λ2, then in the process of multiple iterations, first for β, it is judged whether β satisfies the convergence condition. When β satisfies the convergence condition, the value of β freezes and no longer changes. Then for λ1, it is judged whether λ1 satisfies the convergence condition. When λ1 satisfies the convergence condition, the value of λ1 freezes and no longer changes. Then for λ2, it is judged whether λ2 satisfies the convergence condition. When λ2 satisfies the convergence condition, the value of λ2 freezes and no longer changes. When all three parameters no longer change, the iteration ends.
[0158] In the specific implementation process, the Adam optimizer training realizes fast iteration. For each training epoch = 0, 1,..., nIter - 1 (nIter is the maximum number of epochs), when epoch > 2000, the parameter convergence check and gradient freezing functions are started. First, the unknown parameter β with the highest priority is judged. After β converges, the convergence of other unknown parameters is checked, and the gradients are judged and frozen step by step according to the priority; after the Adam optimizer training reaches nIter, the LBFGS optimizer training is used to realize fine optimization to improve the accuracy of parameter identification. The parameter convergence check and gradient freezing functions are always executed during the LBFGS optimizer training stage. For the unknown parameter with the lowest priority, the final identification result is determined by the value at the end of the LBFGS optimizer training.
[0159] Among them, judging whether the unknown parameter satisfies the convergence condition specifically includes:
[0160] Calculate the relative change rate θ between the unknown parameter of the current iteration and the unknown parameter of the previous iteration Relative_Change , which is specifically realized by the following formula:
[0161]
[0162] Among them, θ new is the unknown parameter of the current iteration, and θ old is the unknown parameter of the previous iteration;
[0163] It should be noted that since the process of judging unknown parameters is to judge one unknown parameter at a time, in the above formula, θ represents one of all unknown parameters, rather than the unknown parameter being θ.
[0164] Judge the relative change rate θ Relative_Change Whether it is less than the expected threshold θ ε , in specific implementation, θ ε = 1e-7. When the relative change rate θ Relative_Change is not less than the expected threshold θ ε , it indicates that this unknown parameter does not meet the convergence condition;
[0165] When the relative change rate θ Relative_Change is less than the expected threshold θ ε , the convergence count θ count increases by one. In program implementation, it can be expressed by the following formula:
[0166] θ count = θ count + 1 if θ Relative_Change < θ ε ;
[0167] Judge whether the increased convergence count reaches the expected convergence count threshold θ Design_Count , in specific implementation, θ Design_Count = [5, 10, 15, 20]. When the increased convergence count does not reach the expected convergence count threshold θ Design_Count , it indicates that this unknown parameter does not meet the convergence condition. In program implementation, it can be expressed by the following formula:
[0168] θ count >= θ Design_Count ;
[0169] When the increased convergence count reaches the expected convergence count threshold θ Design_Count , it indicates that this unknown parameter meets the convergence condition. At this time, freeze the gradient of this parameter and stop optimizing this unknown parameter. The program language can be expressed as:
[0170] θ requires_grad = False.
[0171] A method for solving the inverse problem of partial differential equations under sparse data provided by the present invention can be implemented based on a computer. Specifically, a deep neural network can be designed as a neural network module related to drugs to approximate the solution function of the PDE to be solved, providing numerical calculation support for the subsequent design of the loss function module. The input of the network is the spatial and temporal coordinates (x, t) in the solution domain, and the output is the solution function value Construct a neural network through a multi-layer fully connected network, and use the non-linear activation function Tanh to improve the non-linear fitting ability of the model; Steps 1 and 2 can be designed as an internal sampling point compensation module for the solution domain, which is used to compensate for the sparse data area in the physical solution domain, improve the coverage rate of training data by compensating sampling points, and thus enhance the prediction accuracy of the model in the whole domain; Among them, the part for calculating the loss, that is, Steps 4-8, can be designed as a loss function module. The loss function module contains loss functions of multiple components, which are used to evaluate the prediction performance of the neural network, guide the goal of model training, ensure that the model simultaneously meets the requirements of data fitting, physical constraints and initial condition constraints, and through a dynamic weight adjustment mechanism, dynamically balance the weights according to the changes of each part of the loss during the training process, avoid the dominance of a certain part of the loss in the optimization process. Steps 3, 9 and 10 can be designed as a parameter convergence check and gradient freezing module. The parameter convergence check and gradient freezing module dynamically checks the convergence of each parameter during the training process, gradually freezes the converged parameters, improves the stability of optimization, avoids the interference caused by parameter coupling, and at the same time reduces the calculation cost.
[0172] The present invention can accurately identify the unknown parameters contained in the PDE and its initial conditions when there are few observed data points. In addition, the present invention also has strong anti-noise performance. When the observed data contains a certain degree of noise, the model can effectively reduce the interference of the noise on the identification result, and thus show good robustness in the identification accuracy of the unknown parameters.
[0173] Based on the method of the present invention, the following experiments were carried out for verification:
[0174] The specific implementation case is the Allen–Cahn Equation (AC equation). The AC equation is widely used in the fields of materials science, biology, chemistry, etc., such as simulating the phase separation of multi-component alloys, the kinetic behavior in reaction-diffusion systems, etc. There are many difficulties in parameter identification of the AC equation: there is a coupling relationship between multiple parameters, different parameters have different sensitivities to the solution, and it is easy to become unbalanced during the optimization process. Solving these problems requires combining an efficient optimization strategy and a robust algorithm design. The AC equation is expressed as follows:
[0175]
[0176] Among them, the parameters in the AC formula are known. The present invention sets the parameters in the formula as unknown parameters, executes the solution of the present invention, and after obtaining the specific values of the unknown parameters after execution, compares them with the true parameters in the AC formula to verify the accuracy of the present invention. Based on this, first set the parameters in the formula as unknown parameters. Specifically, the unknown parameters to be identified for the AC equation are set as λ 1 、λ 2 、λ3 , the form with unknown parameters is: The unknown parameter to be identified in the initial condition of the AC equation is set as β, and the form with unknown parameters is: u(0,x) = βx 2 cos(πx), and set the priority of the unknown parameter.
[0177] On the public dataset, randomly collect the coordinates and solutions of the time and space distributions of 10 initial data points without noise, the coordinates and solutions of the time and space distributions of 10 boundary data points without noise, and the coordinates and solutions of the time and space distributions of 200 data points inside the solution domain without noise. Specifically,
[0178] U initial = [0.01361533, 0.03893767, 0.05425687, 0.05140531, -0.37982809, -0.37096784, -0.81483880, -0.01304558, -0.38875525, -0.21013871]
[0180]
[0181]
[0182]
[0183]
[0184] U interior=[-0.04504976,-0.99998682,-0.99756653,-0.99979940,0.41343581,0.17044241,-0.33408914,-0.40738273,-0.99998299,-0.60547817,-0.94374751,0.30930348,-0.99081770,-0.27049433,-0.03417429,-0.99995561,-0.96423673,0.09263523,-0.99304911,-0.99947628,0.93848746,-0.99973528,-0.99190750,-0.29098350,-0.97767063,-0.99870018,-0.99838377,0.94522235,0.33045712,0.97807428,0.76749687,-0.87702657,-0.97409008,0.57797891,0.83829833,0.05224407,-0.99820071,-0.74171453,-0.87643489,-0.98493964,0.06380443,0.47611510,-0.99798102,0.95169080,0.17869050,0.90317583,0.26457011,-0.26230160,0.92258924,-0.34485933,0.25941561,0.03753379,-0.93617680,-0.99998782,0.35854831,0.16227733,0.34472751,0.08707103,0.00411241,0.85975189,-0.83946491,0.77106399,-0.14497518,-0.60783067,0.78897271,-0.95080260,0.29569362,0.14526632,0.96993518,-0.91575881,-0.98160365,-0.99787833,0.07642049,-0.92445477,0.81943281,-0.99995422,-0.95794368,0.01831883,0.71382250,-0.99999851,-0.99145386,0.06462871,0.48749033,0.61365062,-0.99970100,0.00397423,-0.11334294,-0.99515189,0.07702557,-0.17916535,0.08078894,0.40095735,0.36248878,0.89525211,-0.02119549,0.11325551,-0.99856507,-0.24130411,0.00113463,-0.99966178,0.47824252,-0.99942660,0.01695351,0.44563447,0.08938043,-0.98234321,0.06759794,0.13446692,-0.52063881,0.47924384,-0.99996464,0.00256518,0.87172461,0.28013689,-0.99934899,0.21237768,-0.99357470,0.01451108,0.01336128,0.00202649,0.13984994,-0.07034347,0.01563097,0.97518663,0.13422227,0.39066506,-0.60221768,-0.87438462,0.04900193,0.02406382,0.00346158,-0.63412629,-0.72160763,-0.72160763,0.09842654,0.05058123,-0.99990086,-0.98693382,0.04514366,0.10378661,0.59153346,0.22669739,0.15528093,-0.97601685,-0.98981881,0.16887490,0.70325853,0.62891566,-0.04568577,-0.99875436,0.12850511,-0.89554122,-0.99995930,-0.95111502,0.68084814,-0.78304073,-0.39548710,-0.99571583,-0.98757378,-0.07116378,-0.98989645,-0.05490898,0.93498246,0.56725291,0.94108557,0.01068767,0.56333055,0.97000056,-0.98886501,-0.99996319,0.21655757,-0.99900500,0.11877051,-0.99842935,-0.79847068,0.12475088,0.52508538,0.11970765,0.92543041,-0.99999209,0.35013885,-0.28007962,-0.99530704,0.05973886,-0.09296339,0.08913897,-0.99999224,0.89549617,0.00652513,-0.28230378,0.42694854,-0.99764142,-0.24130411,-0.98458375,0.00065460,-0.99232849,0.17191077,-0.99360006,-0.90436637,-0.99977143].
[0186] Obtain the coordinates of the time and space distributions of all data points on the solution domain to obtain the first set X star , X star contains 102,912 data points. Specifically, X star= [(-1.0, 0.0), (-0.99609375, 0.0), (-0.9921875, 0.0), (-0.98828125, 0.0), (-0.984375, 0.0), (-0.98046875, 0.0), (-0.9765625, 0.0), (-0.97265625, 0.0), (-0.96875, 0.0), (-0.96484375, 0.0), (-0.9609375, 0.0), (-0.95703125, 0.0), (-0.953125, 0.0), (-0.94921875, 0.0), (-0.9453125, 0.0), (-0.94140625, 0.0), (-0.9375, 0.0), (-0.93359375, 0.0), (-0.9296875, 0.0), (-0.92578125, 0.0)……(0.90625, 1.0), (0.91015625, 1.0), (0.9140625, 1.0), (0.91796875, 1.0), (0.921875, 1.0), (0.92578125, 1.0), (0.9296875, 1.0), (0.93359375, 1.0), (0.9375, 1.0), (0.94140625, 1.0), (0.9453125, 1.0), (0.94921875, 1.0), (0.953125, 1.0), (0.95703125, 1.0), (0.9609375, 1.0), (0.96484375, 1.0), (0.96875, 1.0), (0.97265625, 1.0), (0.9765625, 1.0), (0.98046875, 1.0), (0.984375, 1.0), (0.98828125, 1.0), (0.9921875, 1.0), (0.99609375, 1.0)];
[0187] Calculate x i With the first set X star The Euclidean distance of the coordinates of the time and space distribution of each data point in, for the coordinates of the time and space distribution of the data point x i With the first set X star Sort the Euclidean distances of the coordinates of the time and space distribution of each data point in from smallest to largest, and obtain the coordinates of the time and space distribution of the data points in the first set X corresponding to the first 20 + 1 Euclidean distances star In, to obtain the initial neighboring points N k (x i ) set, in the initial neighboring point set Nk (x i ) to remove the coordinates of the time and space distribution of the data point x i to obtain the nearest neighbor point set N k '(x i );
[0188] For the nearest neighbor point set N k '(x i ) add random noise to the coordinates y of the time and space distribution of each nearest neighbor point i to obtain the coordinates of the time and space distribution of the perturbed nearest neighbor points Specifically, through the following formula
[0189]
[0190] Judge whether the coordinates of the time and space distribution of each perturbed nearest neighbor point satisfy the boundary conditions. Set the lower boundary lb of the solution domain as [-1 0], set the upper boundary ub of the solution domain as [0.99609375 1], and set the minimum distance min_distance close to the boundary as 0.1. Then, among the perturbed nearest neighbor points when the boundary conditions are satisfied, represent the perturbed nearest neighbor point as a valid point. The coordinates of the time and space distribution of all valid points form the set V(x of the valid neighbors of the current point i );For the N u coordinate sets of the data points inside the solution domain merge the sets of valid neighbors of each data point to obtain the final sampling point set X final , X final includes 2,851 data points. Specifically, X final=[(-0.55705201,0.18634909),(-0.52921514,0.15693076), (-0.55916010,0.18329062), (-0.47368065,0.15429404),(-0.49712708,0.16541367), (-0.50040462,0.17724587), (-0.49227157,0.14572259),(-0.51131600,0.13295304), (-0.52171903,0.15086825), (-0.51275045,0.13957872),(-0.54108221,0.13329051), (-0.57386411,0.18855477), (-0.53509578,0.13019990),(-0.56196237,0.18265541), (-0.52133277,0.13284349), (-0.57043657,0.16793872),(-0.55366137,0.10103247), (-0.51572395,0.14421102), (-0.47720053,0.18412519),(-0.55718706,0.20891115)……(-0.80171850,0.30586210),(-0.79659660,0.25445343),(-0.80090684,0.23664183), (-0.77442155,0.30803608),(-0.74928893,0.29860196),(-0.79398570,0.21455207), (-0.76282723,0.28648125),(-0.74367380,0.25174191),(-0.82893883,0.26172001), (-0.77296155,0.21523764),(-0.75289386,0.20944655),(-0.76119840,0.21679851), (-0.80267785,0.27039522),(-0.77476484,0.25545703),(-0.76875783,0.21959238), (-0.79539471,0.22014307),(-0.83844293,0.25184334),(-0.74994998,0.27158035),(-0.75629975,0.31194015),(-0.(77744735, 0.21712534);.
[0191] Combine the final set of sampling points X final with the coordinate sets of N u data points inside the solution domain to obtain a set of compensated data points which includes 3,051 data points. Specifically,[
[0192]
[0193]
[0194] Furthermore, calculate the physical loss according to step 4, where f(x r , t r ) is specifically expressed as;
[0195] f(x r , t r ) = u t - λ 1 u + λ 2 u 3 - λ 3 u xx ;
[0196] In this embodiment, for the deep neural network, the set parameters are layers = [2, 50, 50, 50, 1], that is, the number of neurons in the input layer is 2, the number of neurons in the output layer is 1, including 3 hidden layers, and the number of neurons in each hidden layer is 50. Furthermore, calculate the initial condition loss according to step 5, where the initial equation of the AC equation is expressed as:
[0197] u ic_true (x i , 0) = g(x i , β) = βx i 2 cos(πx i );
[0198] Furthermore, calculate the data loss according to steps 6 to 9 and assign weights, where Figure 2 is a schematic diagram of the weight assignment situation during the training process, where the abscissa represents the number of training times, the ordinate represents the specific value of the weight, and the weighted sum is used to obtain the final loss value loss. According to the final loss value loss, update the weight matrix and bias vector in the deep neural network, and at the same time update the unknown parameters of the current iteration number to obtain the updated unknown parameters;
[0199] According to the priority of the unknown parameters, starting from the unknown parameter with the highest priority, for this unknown parameter, determine whether this unknown parameter meets the convergence condition. Specifically, in combination with step 10, for the unknown parameter β in the initial condition equation, the relative change rate β of the current iteration value and the previous iteration value of β Relative_Change is calculated as:
[0200]
[0201] where β new is the value of the current iteration, and β old is the value of the previous iteration.
[0202] When β Relative_Change is less than the expected threshold β ε = 1e-7, it is considered that this iteration meets the convergence condition, and the convergence count β count is incremented by one:
[0203] β count = β count + 1 if β Relative_Change < β ε ;
[0204] When the convergence count reaches the expected convergence count threshold β Design_Count = 20, it is determined that this parameter has converged:
[0205] β count >= β Design_Count ;
[0206] After β converges, freeze the gradient of this parameter and stop optimizing this unknown parameter:
[0207] β requires_grad = False;
[0208] After the β gradient is frozen, for the unknown parameters λ 1 、λ 2 、λ 3 , preferentially freeze the gradient of the unknown parameter λ 1 , and the relative change rate λ 1 of the current iteration value and the previous iteration value of λ 1_Relative_Change is calculated as:
[0209]
[0210] where λ 1_new is the value of the current iteration, and λ 1_old is the value of the previous iteration.
[0211] When λ 1_Relative_Change is less than the expected threshold λ 1_ε = 1e-7, it is considered that this iteration meets the convergence condition, and the convergence count λ1_count Increment once:
[0212] λ 1_count = λ 1_count + 1 if λ 1_Relative_Change < λ 1_ε ;
[0213] The number of convergence times reaches the expected convergence times threshold λ 1_Design_Count = 20, determine that this parameter has converged:
[0214] λ 1_count >= λ 1_Design_Count ;
[0215] λ 1 After convergence, freeze the gradient of this parameter and stop optimizing this unknown parameter:
[0216] λ 1_requires_grad = False;
[0217] λ 1 After the gradient is frozen, for the unknown parameter λ 2 、λ 3 Freeze the gradient of the unknown parameter λ 2 The relative change rate λ 2 between the current iteration value and the previous iteration value of λ 2_Relative_Change is calculated as:
[0218]
[0219] where λ 2_new is the value of the current iteration, and λ 2_old is the value of the previous iteration.
[0220] When λ 2_Relative_Change is less than the expected threshold λ 2_ε = 1e - 7, it is considered that this iteration meets the convergence condition, and the number of convergence times λ 2_count increases once:
[0221] λ 2_count = λ 2_count + 1 if λ 2_Relative_Change < λ 2_ε ;
[0222] The number of convergence times reaches the expected convergence times threshold λ 2_Design_Count = 10, determine that this parameter has converged:
[0223] λ 2_count >= λ 2_Design_Count ;
[0224] λ 2 After convergence, freeze the gradient of this parameter and stop optimizing this unknown parameter:
[0225] λ 2_requires_grad = False;
[0226] The Adam optimizer training achieves fast iteration. For each training epoch where epoch = 0, 1, …, 10000 - 1 (10000 is the maximum number of epochs), when epoch > 2000, the parameter convergence check and gradient freezing functions are started. First, the unknown parameter β with the highest priority is judged. After β converges, the convergence of unknown parameters is checked, and the gradients are judged and frozen step by step according to the priority; after the Adam optimizer training reaches 10000 times, the LBFGS optimizer training is used to achieve fine optimization and improve the accuracy of parameter identification. The parameter convergence check and gradient freezing functions are always executed during the LBFGS optimizer training stage for the unknown parameter λ 3 , and the final identification result is determined by the value at the end of the LBFGS optimizer training.
[0227] At the 8535th training time, β count >= β Design_Count , the convergence value β = 0.99224716. At this time, the β gradient is frozen, and 0.99224716 is considered the final identification result of β. After the Adam optimizer training reaches 10000 times, the LBFGS optimizer continues training: at the 12230th training time, λ 1_count >= λ 1_Design_Count , the convergence value λ 1 = 4.99333715. At this time, the λ 1 gradient is frozen, and 4.99333715 is considered the 1 final identification result of λ; at the 13052nd training time, λ 2_count >= λ 2_Design_Count , the convergence value λ 2 = 4.99552059. At this time, the λ 2 gradient is frozen, and 4.99552059 is considered the 2 final identification result of λ; at the 17414th training time, the LBFGS optimizer training ends, and the convergence value λ 3 = 0.00009997, and 0.00009997 is considered the 3 final identification result of λ. In the original AC equation, β = 1, λ 1 = 5, λ 2 = 5, λ 3 = 0.0001. The relative percentage deviations of the identification results of the present invention from the original parameters are respectively: 0.77528%, 0.13326%, 0.08959%, 0.03260%, achieving accurate identification. Figure 3 shows the unknown parameter β during the training process (i.e.,Figure 3 bate) in it, λ 1 (i.e., Figure 3 lambda_1 in it), λ 2 (i.e., Figure 3 lambda_2 in it), the convergence situation of Figure 4 shows the unknown parameter λ during the training process 3 (i.e., Figure 3 lambda_3 in it), where Figure 3 and Figure 4 the abscissas of both represent the number of training times, and the ordinates of both represent the specific values of the unknown parameters.
[0228] According to the above example steps, Table 1 statistically analyzes the relative percentage deviations of the unknown parameters when 0% (the above example), 1%, 2%, 3%, 4%, and 5% times the standard deviation Gaussian noise are added to the training data of the AC equation, and Table 2 statistically analyzes the parameter values of the example.
[0229] Table 1 Relative percentage deviations of unknown parameters when 0%, 1%, 2%, 3%, 4%, and 5% times the standard deviation Gaussian noise are added
[0230] AC equation β <![CDATA[λ 1 > <![CDATA[λ 2 > <![CDATA[λ 3 > 0% 0.77528% 0.13326% 0.08959% 0.03260% 1% 0.81713% 0.03006% 0.01632% 2.06114% 2% 0.74173% 0.07901% 0.14448% 3.65281% 3% 1.17894% 0.22871% 0.38520% 4.44174% 4% 1.90691% 0.29624% 0.32673% 0.26710% 5% 2.14897% 0.55083% 0.35291% 3.68811%
[0231] Table 2 Parameter values of the example
[0232]
[0233] Thus, based on the solution of the present invention, for the implementation case of the Burgers equation:
[0234]
[0235] u(0,x) = -sin(πx),
[0236] u(t, -1) = u(t, 1) = 0,
[0237] x ∈ [-1, 1], t ∈ [0, 1]
[0238] For the unknown parameters to be identified in the Burgers equation, they are set as λ 1 and λ 2 , and the form containing the unknown parameters is: The unknown parameters to be identified in the initial conditions of the Burgers equation are set as β, and the form containing the unknown parameters is: u(0,x) = βsin(πx). On the public dataset, 10 initial point coordinates and solutions, 10 boundary point coordinates and solutions, and 200 internal point coordinates and solutions in the solution domain are randomly collected. Using the present invention for λ 1 and λ 2Identify the unknown parameters of α and β. Table 3 tabulates the relative percentage deviations of the unknown parameters when Gaussian noise with 0%, 1%, 2%, 3%, 4%, and 5% of the standard deviation is added to the training data of the Burgers equation. Table 4 tabulates the parameter values of the embodiments.
[0239] Table 3 Relative Percentage Deviations of Unknown Parameters with Gaussian Noise of 0%, 1%, 2%, 3%, 4%, and 5% of the Standard Deviation Added
[0240] Burgers equation β <![CDATA[λ 1 > <![CDATA[λ 2 > 0% 0.41296% 0.09155% 1.89574% 1% 0.46654% 0.05527% 2.87930% 2% 0.91312% 0.72419% 1.65975% 3% 1.03422% 1.47704% 5.40631% 4% 0.42570% 0.93607% 0.94435% 5% 0.45736% 1.21017% 3.23342%
[0241] Table 4 Parameter Values of the Embodiments
[0242] Burgers equation <![CDATA[β Design_Count / λ 1_Design_Count > <![CDATA[β ε / λ 1_ε > min_distance / k / noise_scale 0% 10 / 5 1.00E-07 0.01 / 20 / 0.1 1% 10 / 5 1.00E-07 0.01 / 20 / 0.1 2% 10 / 5 1.00E-07 0.01 / 25 / 0.1 3% 10 / 5 1.00E-07 0.01 / 25 / 0.1 4% 10 / 5 1.00E-07 0.01 / 30 / 0.1 5% 10 / 5 1.00E-07 0.01 / 30 / 0.1
[0243] Thus, based on the solution of the present invention, for Equation Implementation Example:
[0244]
[0245] h(0,x) = 2sech(x),
[0246] h(t, -5) = h(t, 5),
[0247]
[0248]
[0249] For the unknown parameters to be identified in the equation are set as λ 1 and λ 2 , and the form containing the unknown parameters is: For the initial conditions of the equation, the unknown parameter to be identified is set as β, and the form containing the unknown parameter is: h(0,x) = βsech(x). On the public dataset, 10 initial point coordinates and solutions, 10 boundary point coordinates and solutions, and 200 solution domain interior point coordinates and solutions are randomly collected. Using the present invention, the unknown parameters of λ 1 and λ 2 and β are identified. Table 5 tabulates the relative percentage deviations of the unknown parameters when Gaussian noise with 0%, 1%, 2%, 3%, 4%, and 5% of the standard deviation is added to the training data of the equation. Table 6 tabulates the parameter values of the embodiments.
[0250] Table 5 Relative Percentage Deviations of Unknown Parameters with Gaussian Noise of 0%, 1%, 2%, 3%, 4%, and 5% of the Standard Deviation Added
[0251]
[0252] Table 6 Parameter Values of Embodiments
[0253]
[0254] The above description is only for the preferred embodiments of the present disclosure and the explanation of the applied technical principles. Those skilled in the art should understand that the scope of the invention involved in the embodiments of the present disclosure is not limited to the technical solutions formed by the specific combination of the above technical features, but should also cover other technical solutions formed by any combination of the above technical features or their equivalent features without departing from the above inventive concept. For example, the technical solutions formed by mutually replacing the above features with the technical features (but not limited to) having similar functions disclosed in the embodiments of the present disclosure.
Claims
1. A method for solving the inverse problem of partial differential equations under sparse data, characterized in that: include: Step 1: Collect data from multiple observation points in the actual industrial process; The data of the observation data points include the coordinate set of all observation data points. and The corresponding solution set U data , the coordinate set of all observed data points At least a set of coordinates of multiple initial data points And solve for N inside the domain u The coordinates of the data points The coordinate set of the multiple initial data points Including the time and space distribution coordinates of all initial data points, the N inside the solution domain u The coordinates of the data points Including N inside the solution domain u The coordinates of the temporal and spatial distribution of the data points; Step 2: According to N inside the solution domain u The coordinates of the data points Perform KNN compensation to obtain a set of compensated data points Step 3: Obtain the partial differential equation and initial condition equation of the actual industrial process, wherein the partial differential equation contains multiple unknown parameters. According to the partial differential equation and the initial condition equation, the priorities of the unknown parameters are sorted, wherein the unknown parameters in the initial condition equation have the highest priority. For the partial differential equation, for items containing independent variables, the fewer the number of independent variables, the higher the priority of the corresponding unknown parameters. In the case of the same number of independent variables, the items with lower derivative orders have higher priorities for the unknown parameters. In the case of the same derivative orders, the items with lower nonlinearity have higher priorities for the unknown parameters. Initialize all unknown parameters to obtain initial values of the unknown parameters, set the initial number of iterations, use the initial number of iterations as the current number of iterations, and use the initial values of the unknown parameters as the unknown parameters of the current number of iterations. Step 4: Based on the compensation data point set Partial differential equation and unknown parameters of the current iteration number, calculate the physical loss loss phsics ; Step 5: Coordinate set based on multiple initial data points Initial condition equation and unknown parameters of the current iteration number, calculate the initial condition loss loss ic ; Step 6: Based on the coordinate set of all observed data points and The corresponding solution set U data , calculate data loss loss data ; Step 7: Set physical loss phsics The loss weight w phsics , set the initial condition loss loss ic The loss weight w ic , set data loss loss data The loss weight w data ; Step 8: According to the loss weight w phsics , loss weight w ic and loss weight w data , for physical loss phsics , initial condition loss loss ic and data loss data Perform weighted summation to obtain the final loss value loss; Step 9: According to the final loss value, update the weight matrix and bias vector in the deep neural network, and update the unknown parameters of the current iteration to obtain the updated unknown parameters; Step 10: According to the priority of the unknown parameters, starting from the unknown parameter with the highest priority, for the unknown parameter, determine whether the unknown parameter meets the convergence condition. If the unknown parameter does not meet the convergence condition, the current iteration number is increased by one, and it is used as the new current iteration number. The updated unknown parameter is used as the unknown parameter of the new current iteration number, and return to step 4; When the unknown parameter meets the convergence condition, the unknown parameter of the current iteration number is taken as the final value of the unknown parameter, the current iteration number is increased by 1, and it is taken as the new current iteration number, and the value of the unknown parameter is kept unchanged. For the unknown parameters other than the unknown parameter, the updated value is taken as the new unknown parameter of the current iteration number, and the process returns to execute step 4. In step 9, the unknown parameters that remain unchanged are not updated. In step 10, it is determined whether the unknown parameters of the next priority of the unknown parameter meet the convergence condition, thereby obtaining the values of all unknown parameters.
2. The method for solving the inverse problem of partial differential equations under sparse data according to claim 1, characterized in that: Step 2 specifically includes: Step 2.1: Obtain the coordinates of the time and space distribution of all data points in the solution domain and obtain the first set X star , the first set X star Including N inside the solution domain u The coordinates of the data points Among them, N inside the solution domain u The coordinates of the data points It is expressed as: Step 2.2: For N inside the solution domain u The coordinates of the data points The time and space coordinates x of each data point in i , calculate x i With the first set X star The Euclidean distance of the time and space coordinates of each data point in is calculated by the following formula: Among them, d represents the dimension of the space, y m Denotes the first set X star The time and space coordinates of the data points in , d(x i ,y m ) represents x i and m The Euclidean distance between Represents x i The coordinate value in the pth dimension, Represents y m Coordinate value in the pth dimension; Step 2.3: The time and space coordinates x of the data points i With the first set X star Sort the Euclidean distances of the time and space coordinates of each data point in from small to large, and obtain the first set X corresponding to the first k+1 Euclidean distances star The time and space coordinates of the data points in the initial neighbor point N k (x i ) set, in the initial neighboring point set N k (x i ) to remove the data point x i The time and space coordinates of the nearest neighbor point set N′ k (x i ); Step 2.4: For the nearest neighbor point set N′ k (x i ) in time and space coordinates y of each nearest neighbor point i Add random noise to obtain the time and space coordinates of the nearest neighbor after disturbance This is achieved specifically through the following formula: Where r is a random vector sampled from a uniform distribution, ranging from [-noise_scale, noise_scale], and noise_scale is the amplitude of the random noise; Step 2.5: Determine the time and space coordinates of each perturbed nearest neighbor Whether the boundary conditions are met is specifically expressed by the following formula: Among them, lb is the lower boundary of the solution domain, ub is the upper boundary of the solution domain, and min_distance is the minimum distance close to the boundary; The nearest neighbor after perturbation When the boundary conditions are met, the nearest neighbor point after the disturbance is the valid point, and the time and space coordinates of all valid points constitute the valid neighbor set V(x i ); Step 2.6: For N inside the solution domain u The coordinates of the data points The valid neighbor set of each data point in is merged to obtain the final sampling point set X final , which is specifically achieved through the following formula: Step 2.7: Set the final sampling point set X final and N inside the solution domain u The coordinates of the data points Merge to get a set of compensated data points This is achieved specifically through the following formula:
3. The method for solving the inverse problem of partial differential equations under sparse data according to claim 1, characterized in that: Step 4 is specifically implemented by the following formula: Where N represents the set of compensation data points The number of data points in (x r ,t r ) is the set of compensation data points The time and space coordinates of the data points in t r is the time variable, x r is a spatial variable, f(x r ,t r ) is the residual of the partial differential equation, which is specifically expressed by the following formula: Among them, F(·) is the partial differential operator, u(x r ,t r ) is the formula under the unknown parameters of the current iteration number, Indicates x r and t r Derivative, s and h represent the derivative of x r and t r The number of derivatives, s and h are both non-negative integers, and λ represents the unknown parameter of the current iteration number.
4. The method for solving the inverse problem of partial differential equations under sparse data according to claim 1, characterized in that: Step 5 specifically includes: The coordinates of multiple initial data points are collected Input into the deep neural network to obtain the coordinate set of multiple initial data points The predicted solution is then used to calculate the initial condition loss loss ic , which is specifically achieved through the following formula: Where P represents the coordinate set The number of initial data points in (x ic ,0) is the coordinate set of multiple initial data points The coordinates in u ic_pred (x ic ,0) is the coordinate set of multiple initial data points The predicted solution, u ic_true (x ic ,0) are unknown parameters based on the initial condition equation and the current number of iterations. ic ,0) to obtain the theoretical solution.
5. The method for solving the inverse problem of partial differential equations under sparse data according to claim 1, characterized in that: Step 6 specifically includes: The coordinates of all observation data points are collected Input into the deep neural network to obtain the coordinate set of all observed data points The predicted solution will be based on the coordinate set of all observed data points The predicted solution and The corresponding solution set U data , calculate the data loss, which is implemented by the following formula: Where M is the coordinate set of all observed data points The number of observation data points in (x d ,t d ) is the coordinate set of all observed data points The coordinates of the observation data points in true (x d ,t d ) is U data Medium (x d ,t d ) corresponds to the true solution, u pred (x d ,t d ) is the coordinate set of all observed data points The prediction solution.
6. The method for solving the inverse problem of partial differential equations under sparse data according to claim 4 or 5, characterized in that: The deep neural network includes an input layer, a hidden layer and an output layer. After the data is input into the deep neural network, it passes through the input layer to the first layer and is calculated by the following formula: z (1) =σ(w (0) x+b (0) ); Among them, z (1) is the output of the first layer, w (0) is the weight matrix of layer 0, b (0) is the bias vector of layer 0, σ represents the activation function, and x is the time and space coordinates of the input deep neural network; Then the output of the first layer is calculated layer by layer through forward propagation. Specifically, the calculation from the th layer to the l+1th layer is expressed by the following formula: With (l+1) =σ (l) (In (l) With (l) +b (l) ),l=1,2,...,L-2; Among them, z (l+1) is the output of the l+1th layer, σ (l) is the activation function of the ,th layer, w (l) is the weight matrix of the lth layer, z (l) is the output of layer l, b (l) is the bias vector of the lth layer; Then in the output layer, it is calculated by the following formula: With (L) =in (L-1) With (L-1) +b (L-1) ; Among them, z (L) is the output of the Lth layer, that is, the output of the deep neural network, w (L-1) is the weight matrix of the L-1th layer, z (L-1) is the output of the L-1 layer, b (L-1) is the bias vector of the L-1th layer.
7. The method for solving the inverse problem of partial differential equations under sparse data according to claim 6, characterized in that: The parameters of the input layer, hidden layer and output layer of the deep neural network are preset, and the parameters of the input layer, hidden layer and output layer are layers = [n0, n1, n2, ..., n L ], where n0 is the number of neurons in the input layer, n L is the number of neurons in the output layer, n1, n2, ..., n L-1 is the number of neurons in the hidden layer.
8. The method for solving the inverse problem of partial differential equations under sparse data according to claim 1, characterized in that: Step 7 specifically includes: Based on physical loss phsics Initial condition loss loss ic The size of the loss weight w is set phsics and initial condition loss ic , which is specifically achieved through the following formula: According to data loss data , calculate the loss weight w data , which is calculated by the following formula: Among them, u pred is the coordinate set of all observed data points The prediction solution of for u pred The gradient norm of is calculated by the following formula:
9. The method for solving the inverse problem of partial differential equations under sparse data according to claim 1, characterized in that: Step 8 is specifically expressed by the following formula: loss=w physics ·loss physics +w ic ·loss ic +w data ·loss data 。 10. The method for solving the inverse problem of partial differential equations under sparse data according to claim 1, characterized in that: In step 10, it is determined whether the unknown parameter meets the convergence condition, which specifically includes: Calculate the relative change rate θ between the unknown parameters of the current iteration and the unknown parameters of the previous iteration Relative_Change , which is specifically achieved through the following formula: Among them, θ new is the unknown parameter of the current iteration number, θ old is the unknown parameter of the last iteration; Determine the relative rate of change θ Relative_Change Is it less than the expected threshold θ? ε , at the relative rate of change θ Relative_Change Not less than the desired threshold θ ε When , it indicates that the unknown parameter does not meet the convergence condition; At the relative rate of change θ Relative_Change is less than the desired threshold θ ε When the convergence number θ count Increase once to determine whether the increased number of convergences reaches the expected convergence threshold θ Design_Count , after the increase, the number of convergences does not reach the expected convergence threshold θ Design_Count In the case of , it indicates that the unknown parameter does not meet the convergence condition; After the increase, the number of convergences reaches the expected convergence threshold θ Design_Count In the case of , it indicates that the unknown parameter meets the convergence condition.