Method and system for testing impact sensitivity of initiating explosive device
By combining sequential optimal experimental design and Gaussian process regression model, the impact sensitivity test parameters of pyrotechnics are dynamically adjusted to achieve adaptive optimization, which solves the problems of strong test path pre-prediction and low efficiency in traditional test methods and improves test accuracy and reliability.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-12-10
- Publication Date
- 2026-03-24
AI Technical Summary
Existing methods for testing the impact sensitivity of pyrotechnics suffer from problems such as strong pre-defined test paths, low efficiency, and poor model adaptability. They cannot achieve adaptive model selection and accurate calculation of sensitivity points across the entire probability range, resulting in a waste of valuable samples and testing resources.
By employing a sequential optimal experimental design method combined with a Gaussian process regression model, and through dynamic selection of kernel functions and multi-threshold convergence criteria, the experimental results are monitored in real time and parameters are dynamically adjusted to achieve adaptive optimization and intelligent decision-making in the experimental process.
It significantly reduces the number of tests, decreases sample consumption, and improves the accuracy and reliability of test results, solving the problems of strong pre-defined test paths and low efficiency in traditional test methods.
Smart Images
Figure CN121725932A_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of pyrotechnic impact sensitivity testing technology, and in particular relates to a method and system for testing the impact sensitivity of pyrotechnics. Background Technology
[0002] As energetic materials, such as explosives, pyrotechnics are most frequently and easily triggered by mechanical impacts during their research, production, storage, and use. Accidental detonation of pyrotechnics due to impacts is common. Therefore, impact sensitivity, as a key safety indicator for measuring the strength of an energetic material's response to mechanical energy, has become a primary consideration in the safety performance evaluation of pyrotechnics.
[0003] However, traditional impact sensitivity testing methods generally suffer from the following problems: First, existing test paths are highly pre-defined, with most tests typically having a fixed number of trials or simple termination conditions, lacking a real-time feedback mechanism based on test results and unable to dynamically adjust test parameters according to the current response. This leads to unnecessary tests continuing even when the sensitivity distribution has fully converged, resulting in a significant waste of valuable samples and testing resources. Second, test efficiency is low, with some methods naturally concentrating their test points around the 50% ignition probability, requiring numerous tests to fit low-probability response points, resulting in high sample consumption. Third, model selection relies on prior experience, exhibiting poor adaptability. Whether using traditional parametric model fitting or emerging machine learning methods, the selection of models or kernel functions largely depends on the operator's experience. Therefore, there is an urgent need for an impact sensitivity testing method that can adaptively select models based on the characteristics of the current test data, accurately calculate sensitivity points across the entire probability range, and intelligently decide on the termination timing. Summary of the Invention
[0004] In view of this, the present invention aims to overcome the shortcomings of the above-mentioned problems in the prior art and proposes a method and system for testing the impact sensitivity of pyrotechnic items.
[0005] To achieve the above objectives, the technical solution of the present invention is implemented as follows:
[0006] The first aspect of this invention provides a method for testing the impact sensitivity of pyrotechnic items, comprising the following steps:
[0007] S1. Impact sensitivity tests were conducted on the pyrotechnics sequentially according to the sequential optimal experimental design method. The test data were recorded and acquired in real time. The test data included the number of tests, the stimulus amount, and the binary response variable generated by the pyrotechnics, denoted as y∈{0,1}. In the binary response variable, 1 represents the pyrotechnics igniting and 0 represents the pyrotechnics not igniting. The stimulus amount was the impact energy applied to the pyrotechnics. The stimulus amount was preprocessed to obtain the logarithm of the stimulus amount.
[0008] S2. Once the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area, the kernel function dynamic selection process is executed, and the kernel function dynamic selection process is as follows:
[0009] A candidate set containing multiple kernel functions is pre-defined, including the quadratic exponential kernel function and the Matrn kernel function. Gaussian process regression models are trained using each kernel function in the candidate set based on the currently acquired historical experimental data. For each trained Gaussian process regression model, the leave-one-out cross-validation negative log-likelihood criterion and the maximum marginal likelihood criterion are calculated in parallel, and a weighted score is calculated for each kernel function. The kernel function with the highest weighted score is selected as the optimal kernel function and locked for use in fitting the Gaussian process regression model for all subsequent monitoring periods of this experiment. The mixed result region is defined as the experimental result region where, in the stimulus quantity space explored using the sequential optimal experimental design method, within a stimulus quantity interval consisting of the maximum stimulus quantity in the all-no-fire zone as the lower limit and the minimum stimulus quantity in the all-fire zone as the upper limit, at least one fire and at least one no-fire occur simultaneously among three or more consecutive adjacent stimulus quantities.
[0010] S3. After the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area, the following steps S4-S6 are executed after each preset number of test cycles are completed.
[0011] S4. Using the logarithm of the stimulus amount in all historical test data as input and the corresponding binary response variable as output, train the Gaussian process regression model using the selected optimal kernel function to obtain the predicted ignition probability value and ignition probability prediction variance of the logarithm of the stimulus amount at the new test point, and obtain the fitting curve of the ignition probability distribution of pyrotechnics in the current test cycle.
[0012] S5. Based on the fitting curve of the ignition probability distribution obtained in step S4, calculate the estimated value of the logarithm of the stimulus quantity corresponding to a set of preset key ignition probability points by piecewise linear interpolation, and obtain the corresponding stimulus quantity estimate.
[0013] S6. Compare the estimated stimulus value of the critical ignition probability point calculated in step S5 with the estimated stimulus value of the critical ignition probability point corresponding to the previous test cycle, calculate the relative rate of change of the estimated value of each critical ignition probability point, and preset a personalized convergence threshold for each critical ignition probability point. If the relative rate of change of all critical ignition probability points is lower than their respective personalized convergence thresholds multiple times consecutively, generate a test termination suggestion; otherwise, return to step S4.
[0014] Furthermore, in step S1, the sequential optimal experimental design method first uses the quartile convergence method to explore the stimulus quantity space and locate the mixed result region where firing and non-firing outcomes coexist. After entering the mixed result region, it switches to the measurement stage based on the Fisher information matrix determinant optimal criterion, and determines the stimulus quantity for the next experimental point by solving the following formula based on the Fisher information matrix determinant optimal criterion.
[0015]
[0016] Where, x n+1 The logarithm of the stimulus intensity at the (n+1)th test point. and The sensitivity distribution parameters are obtained by maximum likelihood estimation based on the current n trial data. Let be the Fisher information matrix, and its expression be .
[0017]
[0018] in, σ is the standard deviation parameter of the sensitivity distribution obtained by maximum likelihood estimation based on historical experimental data. Let be the probability density function of the standard normal distribution. This is the cumulative distribution function of the standard normal distribution.
[0019] Furthermore, in step S2, during the dynamic selection of kernel functions, the preset set of candidate kernel functions is denoted as K, K = {k SE (x i ,x j ),k Matern32 (x i ,x j ),k Matern52 (x i ,x j )}, where k SE (x i ,x j ) is the squared exponential kernel function, used to measure the logarithm x of the stimulus magnitude in the i-th and j-th trials. i and x j The similarity between them, k Matern32 (x i ,x j ), k Matern52 (x i ,x j The following are the Matrn kernel functions for smoothness parameters ν = 3 / 2 and ν = 5 / 2, respectively.
[0020]
[0021] in, is the signal variance hyperparameter; l is the length scale hyperparameter.
[0022] When the smoothness parameter ν = 3 / 2,
[0023]
[0024] When the smoothness parameter ν = 5 / 2,
[0025]
[0026] Furthermore, the weighted score calculation process in step S2 is as follows:
[0027] For each candidate kernel function k∈K, the leave-one-out cross-validation negative log-likelihood criterion and the maximum marginal likelihood criterion are computed in parallel. The formula for calculating the leave-one-out cross-validation negative log-likelihood criterion is as follows:
[0028]
[0029] Among them, X -i y -i To exclude the input set after the i-th data point, y i The output set after excluding the i-th data point; p(y i |X -i ,y -i (k) represents the Gaussian process regression model trained with the remaining data after excluding the i-th data point, applied to y. i The predicted probability;
[0030] The formula for calculating the maximum marginal likelihood criterion is as follows:
[0031] MLL(k) = logp(y|X,k),
[0032]
[0033] Where X = [x1, x2, ..., x n ] T It is a vector consisting of the logarithms of the stimulus quantities from n historical trials, y = [y1,...,y2]. n ] T For the corresponding binary response, K XX It is an n×n covariance matrix, K XX =k(x i ,x j I is calculated from the optimal kernel function selected in step S2. n It is an n×n identity matrix. It is the noise variance hyperparameter of the Gaussian process regression model;
[0034] The leave-one-out cross-validation negative log-likelihood criterion and maximum marginal likelihood criterion for all candidate kernel functions are normalized and mapped to the [0,1] interval. The results are then standardized so that larger values indicate better performance of the Gaussian process regression model, and a weighted average is applied.
[0035] The normalized formula for the leave-one-out cross-validation negative log-likelihood criterion is:
[0036]
[0037] Where min(NLL_LOO) is the minimum value of the leave-one-out cross-validation negative log-likelihood criterion for all candidate kernel functions, and max(NLL_LOO) is the maximum value of the leave-one-out cross-validation negative log-likelihood criterion for all candidate kernel functions.
[0038] The normalized formula for the maximum marginal likelihood criterion is:
[0039]
[0040] Where MLL(k) is the maximum marginal likelihood criterion for the corresponding kernel function k, min(MLL) is the minimum value of the maximum marginal likelihood criterion for all candidate kernel functions, and max(MLL) is the maximum value of the maximum marginal likelihood criterion for all candidate kernel functions.
[0041] S(k) = w·NLL_LOO (norm) (k)+(1-w)·MLL (norm) (k),
[0042] Where w is a preset weighting coefficient, with a value range of 0.3-0.7.
[0043] Furthermore, in step S4, based on the Gaussian process regression model, the predicted ignition probability μ(x) of the logarithm of the stimulus at the new test point is calculated. * ) and the predicted variance of the probability of ignition σ 2 (x * ),
[0044]
[0045] Where, x * The logarithm of the stimulus amount at the new test site.
[0046] y = [y1, y2, ..., y n ] T It is a vector composed of the binary response values of the corresponding n historical impact sensitivity tests.
[0047] k * =[k(x * ,x1),k(x* ,x2),...,k(x * ,x n )] T It is the covariance vector between the test point x and all historical training points X.
[0048] k(x * ,x * ) is x * Its own covariance,
[0049] I n It is an n×n identity matrix.
[0050] Using the cumulative distribution function Φ of the standard normal distribution as the link function, the ignition probability of the logarithm of the stimulus at the new test point is obtained as follows: This serves as a fitted curve for the ignition probability distribution of pyrotechnics during the current test cycle.
[0051] Furthermore, in step S5, the critical ignition probability point includes the stimulus amount E with an ignition probability of 50%. 50 Stimulus level E with a 99.9% probability of triggering an ignition. 99.9 And the stimulus E with a 0.01% probability of igniting. 0.01 .
[0052] Furthermore, the specific process in step S5 is as follows: Based on the logarithmic point sequence of dense stimulus quantities predicted by the Gaussian process regression model (x... i ,P i In this process, the ignition probability corresponding to the preset critical ignition probability point is taken as the target ignition probability, and the value P of the ignition probability is found. i and P i+1 Make the value of the target firing probability p targ Given the probability of ignition P i and P i+1 Between, and then through Perform linear interpolation to obtain the target firing probability value p. target The corresponding logarithm of the stimulus quantity, and calculated using exponential operation E target =exp(x target Transform the logarithm of the stimulus amount back to the actual stimulus amount, where x target p represents the probability of the target firing. target The logarithm of the stimulus amount, E target This represents the actual target stimulus amount obtained through exponential mapping.
[0053] Furthermore, in step S6, for each critical ignition probability point p, the relative rate of change Δ p The calculation formula is in, These are the estimated stimulus values for the key ignition probability points in this test cycle. This is the estimated stimulus amount for the critical ignition probability point in the previous test cycle.
[0054] Furthermore, in step S6, after each preset number of impact sensitivity tests are completed, the test data is updated, the Gaussian process regression model is retrained, and a new fitting curve for the ignition probability distribution and an estimate of the stimulus amount at key ignition probability points are generated.
[0055] A second aspect of the present invention provides a pyrotechnic impact sensitivity testing system for implementing the above-mentioned pyrotechnic impact sensitivity testing method, comprising:
[0056] The data acquisition module is used to conduct impact sensitivity tests on pyrotechnics in sequence according to the sequential optimal experimental design method, record and acquire test data in real time, and preprocess the stimulus amount in the test data to obtain the logarithm of the stimulus amount.
[0057] The first processing module is used to execute the kernel function dynamic selection process when the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area.
[0058] The second processing module is used to sequentially call the third, fourth, and fifth processing modules to perform corresponding data processing after the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area, for each preset number of test cycles completed.
[0059] The third processing module is used to take the logarithm of the stimulus amount in all historical test data as input, the corresponding binary response variable as output, and use the selected optimal kernel function to train the Gaussian process regression model to obtain the predicted value of the ignition probability and the predicted variance of the ignition probability of the logarithm of the stimulus amount at the new test point, as well as the fitting curve of the ignition probability distribution of the pyrotechnics in the current test cycle.
[0060] The fourth processing module is used to calculate the estimated value of the logarithm of the stimulus quantity corresponding to a set of preset key ignition probability points by using piecewise linear interpolation based on the fitted curve of the ignition probability distribution, and to obtain the corresponding stimulus quantity estimate.
[0061] The fifth processing module compares the estimated stimulus value of the critical ignition probability point obtained in this calculation with the estimated stimulus value of the critical ignition probability point corresponding to the previous test cycle, calculates the relative rate of change of the estimated value of each critical ignition probability point, and presets a personalized convergence threshold for each critical ignition probability point. If the relative rate of change of all critical ignition probability points is lower than their respective personalized convergence thresholds for multiple consecutive times, a test termination suggestion is generated; otherwise, the process returns to the third processing module.
[0062] Compared with the prior art, the present invention has the following advantages:
[0063] The pyrotechnic impact sensitivity testing method described in this invention introduces a Gaussian process regression model as a non-parametric sensitivity distribution fitting tool, and employs a dynamic kernel function selection strategy and multi-threshold convergence criteria to achieve adaptive optimization and intelligent decision-making throughout the entire testing process. It solves the problems of traditional impact sensitivity testing methods, such as strong pre-defined test paths, low test efficiency, and large test quantities required for extremely high or low response probability points. Through real-time feedback of the impact sensitivity testing process and intelligent monitoring of test results, it dynamically adjusts test parameters, thereby significantly reducing the number of tests, reducing sample consumption, guiding test point settings, and improving the accuracy and reliability of test results. Attached Figure Description
[0064] The accompanying drawings, which form part of this invention, are used to provide a further understanding of the invention. The illustrative embodiments of the invention and their descriptions are used to explain the invention and do not constitute an undue limitation of the invention. In the drawings:
[0065] Figure 1 This is a flowchart of the pyrotechnic impact sensitivity test method according to Embodiment 1 of the present invention;
[0066] Figure 2 This is a graph showing the trend of key probability point estimates as a function of the number of tests and suggested test termination points in the pyrotechnic impact sensitivity test method described in Embodiment 1 of the present invention. Detailed Implementation
[0067] It should be noted that, unless otherwise specified, the embodiments and features described in the present invention can be combined with each other.
[0068] In the description of this invention, it should be understood that these descriptions are merely exemplary and not intended to limit the scope of this application. In the following detailed description, numerous specific details are set forth to provide a thorough understanding of the embodiments of this application for ease of explanation. However, it will be apparent that one or more embodiments may be practiced without these specific details. Furthermore, descriptions of well-known structures and techniques are omitted in the following description to avoid unnecessarily obscuring the concepts of this application.
[0069] The terminology used herein is for the purpose of describing particular embodiments only and is not intended to limit the scope of this application. The terms “comprising,” “including,” etc., as used herein indicate the presence of the stated features, steps, operations, and / or components, but do not exclude the presence or addition of one or more other features, steps, operations, or components.
[0070] All terms used herein (including technical and scientific terms) have the meanings commonly understood by those skilled in the art, unless otherwise defined. It should be noted that the terms used herein are to be interpreted in a manner consistent with the context of this specification, and not in an idealized or overly rigid way.
[0071] The present invention will now be described in detail with reference to the accompanying drawings and embodiments.
[0072] Example 1
[0073] like Figure 1 As shown, a method for testing the impact sensitivity of pyrotechnic items includes the following steps:
[0074] S1. Impact sensitivity tests are conducted on pyrotechnics sequentially according to the sequential optimal experimental design method. The test data are recorded and acquired in real time. The test data includes the number of tests, the stimulus amount, and the binary response variable generated by the pyrotechnics, denoted as y∈{0,1}. In the binary response variable, 1 represents the pyrotechnics igniting and 0 represents the pyrotechnics not igniting. The stimulus amount is the impact energy applied to the pyrotechnics. The stimulus amount is preprocessed to obtain the logarithm of the stimulus amount, which is used as the input variable for subsequent processes.
[0075] S2. Once the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area, the kernel function dynamic selection process is executed, and the kernel function dynamic selection process is as follows:
[0076] A candidate set of kernel functions is pre-defined, including the quadratic exponential kernel and the Matrn kernel. A Gaussian process regression model is trained using each kernel function from the candidate set based on the currently acquired historical experimental data. For each trained Gaussian process regression model, the leave-one-out cross-validation negative log-likelihood criterion (NLL_LOO) and the maximum marginal likelihood criterion (MLL) are computed in parallel, and a weighted score is calculated for each kernel function. The kernel function with the highest weighted score is selected as the optimal kernel function and locked for use in all subsequent monitoring tests of this experiment. The Gaussian process regression model for the measurement period is fitted; among them, the mixed result area is the experimental result area in which at least one ignition (response value 1) and at least one non-ignition (response value 0) occur simultaneously in three or more consecutive adjacent stimulus points within the stimulus quantity interval formed by the maximum stimulus quantity in the non-ignition zone as the lower limit and the minimum stimulus quantity in the ignition zone as the upper limit when exploring the stimulus quantity space by the sequential optimal experimental design method; in this embodiment, the preset activation threshold is determined according to the type of pyrotechnics and the experimental accuracy requirements, and it must be ensured that when the threshold is reached, it can be preliminarily determined whether the experiment has entered the mixed result area.
[0077] S3. After the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area, the following steps S4-S6 are executed after each preset number of test cycles are completed.
[0078] S4. Using the logarithm of the stimulus amount in all historical test data as input and the corresponding binary response variable as output, train a Gaussian process regression model using the selected optimal kernel function to obtain the predicted ignition probability value and ignition probability prediction variance of the logarithm of the stimulus amount at the new test point, and obtain the fitting curve of the ignition probability distribution of the pyrotechnics in the current test cycle; In this embodiment, all historical test data in step S4 includes the historical test data in step S2 and the test data within the test cycle completed in step S3;
[0079] S5. Based on the fitting curve of the ignition probability distribution obtained in step S4, calculate the estimated value of the logarithm of the stimulus quantity corresponding to a set of preset key ignition probability points by piecewise linear interpolation, and obtain the corresponding stimulus quantity estimate.
[0080] S6. Compare the estimated stimulus value of the critical ignition probability point calculated in step S5 with the estimated stimulus value of the critical ignition probability point corresponding to the previous experimental cycle, and calculate the relative rate of change of the estimated value of each critical ignition probability point; and preset a personalized convergence threshold for each critical ignition probability point. If the relative rate of change of all critical ignition probability points is lower than their respective personalized convergence thresholds for multiple consecutive cycles, a trial termination suggestion is generated; otherwise, return to step S4. In this embodiment, "multiple consecutive cycles" specifically refers to 2 to 4 consecutive cycle cycles.
[0081] In step S1, in order to conform to the statistical law that the sensitivity of pyrotechnics usually follows a log-normal distribution, the stimulus quantity E needs to be preprocessed by taking its logarithm and converting it into the model input feature x = ln(E), and the output is a binary response variable y ∈ {0,1}.
[0082] In step S1, the sequential optimal experimental design method first uses the quartile convergence method to quickly explore the stimulus space and locate the mixed result region where firing and non-firing outcomes coexist. After entering the mixed result region, it switches to the precision measurement stage based on the Fisher information matrix determinant optimal criterion. Based on the Fisher information matrix determinant optimal criterion, the stimulus quantity for the next experimental point is determined by solving the following formula.
[0083]
[0084] Where, x n+1 The logarithm of the stimulus intensity at the (n+1)th test point. and The sensitivity distribution parameters are obtained by maximum likelihood estimation based on the current n trial data. Let be the Fisher information matrix, and its expression be .
[0085]
[0086] in, σ is the standard deviation parameter of the sensitivity distribution obtained by maximum likelihood estimation based on historical experimental data. Let be the probability density function of the standard normal distribution. This is the cumulative distribution function of the standard normal distribution.
[0087] In step S2, during the dynamic selection of kernel functions, the preset set of candidate kernel functions is denoted as K, where K = {k SE (x i ,x j ),k Matern32 (x i ,x j ),k Matern52 (x i ,x j )}, where k SE (x i ,x j ) is the squared exponential kernel function, used to measure the logarithm x of the stimulus magnitude in the i-th and j-th trials. i and x j The similarity between them, k Matern32 (x i ,x j ), k Matern52 (x i ,x j The following are the Matrn kernel functions for smoothness parameters ν = 3 / 2 and ν = 5 / 2, respectively.
[0088]
[0089] in, is the signal variance hyperparameter; l is the length scale hyperparameter.
[0090] When the smoothness parameter ν = 3 / 2,
[0091]
[0092] When the smoothness parameter ν = 5 / 2,
[0093]
[0094] The weighted score calculation process in step S2 is as follows:
[0095] For each candidate kernel function k∈K, the leave-one-out cross-validation negative log-likelihood criterion and the maximum marginal likelihood criterion are computed in parallel. The formula for calculating the leave-one-out cross-validation negative log-likelihood criterion is as follows:
[0096]
[0097] Among them, X -i y -i To exclude the input set after the i-th data point, y i The output set after excluding the i-th data point; p(y i |X -i ,y -i (k) represents the Gaussian process regression model trained with the remaining data after excluding the i-th data point, applied to y. i The predicted probability;
[0098] The formula for calculating the maximum marginal likelihood criterion is as follows:
[0099] MLL(k) = logp(y|X,k),
[0100]
[0101] Where X = [x1, x2, ..., x n ] T It is a vector consisting of the logarithms of the stimulus quantities from n historical trials, y = [y1,...,y2]. n ] T For the corresponding binary response, K XX It is an n×n covariance matrix, K XX =k(x i ,x j I is calculated from the optimal kernel function selected in step S2. n It is an n×n identity matrix. It is the noise variance hyperparameter of the Gaussian process regression model;
[0102] The leave-one-out cross-validation negative log-likelihood criterion and maximum marginal likelihood criterion for all candidate kernel functions are normalized and mapped to the [0,1] interval. The results are then standardized so that larger values indicate better performance of the Gaussian process regression model, and a weighted average is applied.
[0103] The normalized formula for the leave-one-out cross-validation negative log-likelihood criterion is:
[0104]
[0105] Where min(NLL_LOO) is the minimum value of the leave-one-out cross-validation negative log-likelihood criterion for all candidate kernel functions, and max(NLL_LOO) is the maximum value of the leave-one-out cross-validation negative log-likelihood criterion for all candidate kernel functions.
[0106] The normalized formula for the maximum marginal likelihood criterion is:
[0107]
[0108] Where min(NLL_LOO) is the minimum value of the leave-one-out cross-validation negative log-likelihood criterion for all candidate kernel functions, and max(NLL_LOO) is the maximum value of the leave-one-out cross-validation negative log-likelihood criterion for all candidate kernel functions.
[0109] The normalized formula for the maximum marginal likelihood criterion is:
[0110]
[0111] Where MLL(k) is the maximum marginal likelihood criterion for the corresponding kernel function k, min(MLL) is the minimum value of the maximum marginal likelihood criterion for all candidate kernel functions, and max(MLL) is the maximum value of the maximum marginal likelihood criterion for all candidate kernel functions.
[0112] The weighted scoring formula for kernel function k is as follows:
[0113] S(k) = w·NLL_LOO (norm) (k)+(1-w)·MLL (norm) (k),
[0114] Where w is a preset weighting coefficient, with a value range of 0.3-0.7.
[0115] In step S4, based on the Gaussian process regression model, the predicted ignition probability μ(x) of the logarithm of the stimulus at the new test point is calculated. * ) and the predicted variance of the probability of ignition σ 2 (x * ),
[0116]
[0117]
[0118] Where, x * This represents the logarithm of the stimulus amount at the new test site, i.e., the logarithm of the stimulus amount to be selected for the next test.
[0119] y = [y1, y2, ..., y n ] TIt is a vector composed of the binary response values of the corresponding n historical impact sensitivity tests.
[0120] k * =[k(x * ,x1),k(x * ,x2),...,k(x * ,x n )] T It is the covariance vector between the test point x and all historical training points X.
[0121] k(x * ,x * ) is x * Its own covariance,
[0122] I n It is an n×n identity matrix.
[0123] Using the cumulative distribution function Φ of the standard normal distribution as the link function, the ignition probability of the logarithm of the stimulus at the new test point is obtained as follows: This serves as a fitted curve for the ignition probability distribution of pyrotechnics during the current test cycle.
[0124] In step S5, the critical ignition probability point includes the stimulus amount E with an ignition probability of 50%. 50 Stimulus level E with a 99.9% probability of triggering an ignition. 99.9 And the stimulus E with a 0.01% probability of igniting. 0.01 In this embodiment, the stimulus amount in each trial is considered as a point; therefore, the stimulus amount E with a 50% probability of ignition is... 50 This can be understood as the location corresponding to a 50% probability of ignition.
[0125] The specific process in step S5 is as follows: Based on the logarithmic point sequence of dense stimulus quantities predicted by the Gaussian process regression model (x... i ,P i In this process, the ignition probability corresponding to the preset critical ignition probability point is taken as the target ignition probability, and the value P of the ignition probability is found. i and P i+1 Make the value of the target firing probability p tar Given the probability of ignition P i and P i+1 Between, and then through Perform linear interpolation to obtain the target firing probability value p. target The corresponding logarithm of the stimulus quantity, and calculated using exponential operation E target =exp(x target Transform the logarithm of the stimulus amount back to the actual stimulus amount, where x target p represents the probability of the target firing. targetThe logarithm of the stimulus quantity, i.e., the logarithm of the target stimulus quantity, E target This represents the actual target stimulus quantity obtained through exponential mapping. In this embodiment, the logarithm of the stimulus quantity corresponding to the target ignition probability is the estimated value of the logarithm of the stimulus quantity corresponding to the critical ignition probability point, and the actual target stimulus quantity is the estimated stimulus quantity.
[0126] In step S6, for each critical ignition probability point p, the relative rate of change Δ p The calculation formula is in, These are the estimated stimulus values for the key ignition probability points in this test cycle. This is the estimated stimulus amount for the critical ignition probability point in the previous test cycle.
[0127] In step S6, after each preset number of impact sensitivity tests are completed, the test data is updated, the Gaussian process regression model is retrained, and a new fitting curve for the ignition probability distribution and the stimulus quantity estimate for key ignition probability points are generated.
[0128] The following example uses a hybrid propellant to conduct impact sensitivity tests, with a maximum of 90 tests. During the tests, the test energy points are dynamically selected using a sequential optimal design method, and intelligent termination is achieved based on Gaussian process regression and dynamic convergence criteria. In this example, the evaluation period starts from the 40th test, with a convergence evaluation performed every two tests. The candidate kernel function set includes the squared exponential kernel function and the Matrn kernel function (ν = 3 / 2 and ν = 5 / 2). The weighted scoring weights are set to 0.5 for leave-one-out cross-validation negative log-likelihood (NLL_LOO) and maximum marginal likelihood (MLL). The sensitivity distribution is fitted using a log-normal distribution. Four key ignition probability points are selected, including E 0.01 E1, E 50 E 99.9 Each key ignition probability point has a personalized convergence threshold and a number of consecutive successful responses, as detailed in the table below.
[0129] Table 1 Key ignition probability points and their personalized convergence thresholds
[0130]
[0131] (1) Following the sequential optimal experimental design method, impact sensitivity tests were conducted on the pyrotechnic samples in sequence, and relevant data were recorded and obtained, including the number of tests n, the stimulus amount E, and the corresponding binary response variable y∈{0,1}; 40 data were accumulated.
[0132] (2) Based on the data from the first 40 trials, the natural logarithm of the stimulus was used as the input feature for Gaussian process regression, x = ln(E). The preprocessed data from the first 40 trials are shown in the table below.
[0133] Table 2 Pretreatment data from the first 40 experiments
[0134]
[0135] After reviewing the data from the first 40 trials and confirming that the results have entered the mixed region, a weighted average of the kernel functions is calculated based on the preset set of candidate kernel functions K.
[0136] For each candidate kernel function, a Gaussian process regression model is trained based on the data from the first 40 trials, and the leave-one-out cross-validation negative log-likelihood criterion and the maximum marginal likelihood criterion are calculated.
[0137] The leave-one-out cross-validation negative log-likelihood criterion and maximum marginal likelihood criterion for all candidate kernel functions are normalized and mapped to the [0,1] interval, and unified to the direction of "the larger the value, the better". A weighted average is then performed, and the results are shown in the table below.
[0138] Table 3. Calculated values of evaluation metrics for each candidate kernel function.
[0139]
[0140] Based on the weighted scoring results, the squared exponential kernel function was selected as the optimal kernel function for this experiment. The squared exponential kernel function was then used to fit the Gaussian process regression model for all subsequent monitoring periods of this experiment, with the corresponding hyperparameters set as follows: signal variance... Length scale l = 0.42, noise variance σ n =0.31.
[0141] (3) First, construct a training set based on the data from the current 40 trials.
[0142]
[0143] Calculate the covariance matrix between training samples using the selected squared exponential kernel function:
[0144]
[0145] Based on a trained Gaussian process regression model, 1000 test points are generated at 0.01 intervals within the logarithmic interval of stimulus magnitude [7.0, 10.0], and the predicted firing probability for each test point is calculated. Based on any test point x... * The predicted probability of ignition μ(x) * ) and the predicted variance of the probability of ignition σ 2 (x * ), to obtain the ignition probability of the final new test point.
[0146] (4) The Gaussian process regression model yields discrete prediction points. To accurately calculate the key points, piecewise linear interpolation is used to obtain the key ignition probability points, which are the estimated values of the stimulus quantity corresponding to the key ignition probability.
[0147] (5) For each critical ignition probability point p, calculate the relative rate of change Δ p .
[0148] (6) After every two trials (i.e., the evaluation period set in this embodiment), the Gaussian process regression model is retrained. New experimental data is added to the training set, and the covariance matrix is recalculated using a locked quadratic exponential kernel function or optimized hyperparameters. The Gaussian process regression model is then retrained based on all updated data, generating a new fitting curve for the ignition probability distribution and estimated stimulus values for key ignition probability points. If the relative rate of change is less than a set threshold, the convergence counter for that key point is incremented by 1; otherwise, it is reset to zero. When the convergence counters for all four key ignition probability points have reached the preset four consecutive successive thresholds, the system generates a trial termination suggestion. Figure 2 The figure shown is a graph illustrating the trend of the estimated stimulus amount at the key ignition probability point obtained by the pyrotechnic impact sensitivity test method in this example as a function of the number of tests and the recommended test termination point.
[0149] Example 2
[0150] A pyrotechnic impact sensitivity testing system, used to implement the pyrotechnic impact sensitivity testing method described in Example 1, includes:
[0151] The data acquisition module is used to conduct impact sensitivity tests on pyrotechnics in sequence according to the sequential optimal experimental design method, record and acquire test data in real time, and preprocess the stimulus amount in the test data to obtain the logarithm of the stimulus amount.
[0152] The first processing module is used to execute the kernel function dynamic selection process when the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area.
[0153] The second processing module is used to sequentially call the third, fourth, and fifth processing modules to perform corresponding data processing after the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area, for each preset number of test cycles completed.
[0154] The third processing module is used to take the logarithm of the stimulus amount in all historical test data as input, the corresponding binary response variable as output, and use the selected optimal kernel function to train the Gaussian process regression model to obtain the predicted value of the ignition probability and the predicted variance of the ignition probability of the logarithm of the stimulus amount at the new test point, as well as the fitting curve of the ignition probability distribution of the pyrotechnics in the current test cycle.
[0155] The fourth processing module is used to calculate the estimated value of the logarithm of the stimulus quantity corresponding to a set of preset key ignition probability points by using piecewise linear interpolation based on the fitted curve of the ignition probability distribution, and to obtain the corresponding stimulus quantity estimate.
[0156] The fifth processing module compares the estimated stimulus value of the critical ignition probability point obtained in this calculation with the estimated stimulus value of the critical ignition probability point corresponding to the previous test cycle, calculates the relative rate of change of the estimated value of each critical ignition probability point, and presets a personalized convergence threshold for each critical ignition probability point. If the relative rate of change of all critical ignition probability points is lower than their respective personalized convergence thresholds for multiple consecutive times, a test termination suggestion is generated; otherwise, the process returns to the third processing module.
[0157] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.
Claims
1. A method for testing the impact sensitivity of pyrotechnic items, characterized in that, Includes the following steps: S1. Impact sensitivity tests were conducted on the pyrotechnics sequentially according to the sequential optimal experimental design method. The test data were recorded and acquired in real time. The test data included the number of tests, the stimulus amount, and the binary response variable generated by the pyrotechnics, denoted as y∈{0,1}. In the binary response variable, 1 represents the pyrotechnics igniting and 0 represents the pyrotechnics not igniting. The stimulus amount was the impact energy applied to the pyrotechnics. The stimulus amount was preprocessed to obtain the logarithm of the stimulus amount. S2. Once the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area, the kernel function dynamic selection process is executed, and the kernel function dynamic selection process is as follows: A candidate set containing multiple kernel functions is pre-defined, including the quadratic exponential kernel function and the Matrn kernel function. Gaussian process regression models are trained using each kernel function in the candidate set based on the currently acquired historical experimental data. For each trained Gaussian process regression model, the leave-one-out cross-validation negative log-likelihood criterion and the maximum marginal likelihood criterion are calculated in parallel, and a weighted score is calculated for each kernel function. The kernel function with the highest weighted score is selected as the optimal kernel function and locked for use in fitting the Gaussian process regression model for all subsequent monitoring periods of this experiment. The mixed result region is defined as the experimental result region where, in the stimulus quantity space explored using the sequential optimal experimental design method, within a stimulus quantity interval consisting of the maximum stimulus quantity in the all-no-fire zone as the lower limit and the minimum stimulus quantity in the all-fire zone as the upper limit, at least one fire and at least one no-fire occur simultaneously among three or more consecutive adjacent stimulus quantities. S3. After the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area, the following steps S4-S6 are executed after each preset number of test cycles are completed. S4. Using the logarithm of the stimulus amount in all historical test data as input and the corresponding binary response variable as output, train the Gaussian process regression model using the selected optimal kernel function to obtain the predicted ignition probability value and ignition probability prediction variance of the logarithm of the stimulus amount at the new test point, and obtain the fitting curve of the ignition probability distribution of pyrotechnics in the current test cycle. S5. Based on the fitting curve of the ignition probability distribution obtained in step S4, calculate the estimated value of the logarithm of the stimulus quantity corresponding to a set of preset key ignition probability points by piecewise linear interpolation, and obtain the corresponding stimulus quantity estimate. S6. Compare the estimated stimulus value of the critical ignition probability point calculated in step S5 with the estimated stimulus value of the critical ignition probability point corresponding to the previous test cycle, calculate the relative rate of change of the estimated value of each critical ignition probability point, and preset a personalized convergence threshold for each critical ignition probability point. If the relative rate of change of all critical ignition probability points is lower than their respective personalized convergence thresholds multiple times consecutively, generate a test termination suggestion; otherwise, return to step S4.
2. The method for testing the impact sensitivity of pyrotechnic items according to claim 1, characterized in that: In step S1, the sequential optimal experimental design method first uses the quartile convergence method to explore the stimulus space and locate the mixed result region where firing and non-firing outcomes coexist. After entering the mixed result region, it switches to the measurement stage based on the Fisher information matrix determinant optimality criterion. Based on the Fisher information matrix determinant optimality criterion, the stimulus quantity for the next experimental point is determined by solving the following formula. Where, x n+1 The logarithm of the stimulus intensity at the (n+1)th test point. and The sensitivity distribution parameters are obtained by maximum likelihood estimation based on the current n trial data. Let be the Fisher information matrix, and its expression be . in, σ is the standard deviation parameter of the sensitivity distribution obtained by maximum likelihood estimation based on historical experimental data. Let be the probability density function of the standard normal distribution. This is the cumulative distribution function of the standard normal distribution.
3. The method for testing the impact sensitivity of pyrotechnic items according to claim 1, characterized in that: In step S2, during the dynamic selection of kernel functions, the preset set of candidate kernel functions is denoted as K, where K = {k SE (x i ,x j ),k Matern32 (x i ,x j ),k Matern52 (x i ,x j )}, where k SE (x i ,x j ) is the squared exponential kernel function, used to measure the logarithm x of the stimulus magnitude in the i-th and j-th trials. i and x j The similarity between them, k Matern32 (x i ,x j ), k Matern52 (x i ,x j The following are the Matrn kernel functions for smoothness parameters ν = 3 / 2 and ν = 5 / 2, respectively. in, Let l be the signal variance hyperparameter and l be the length scale hyperparameter. When the smoothness parameter ν = 3 / 2, When the smoothness parameter ν = 5 / 2, 4. The method for testing the impact sensitivity of pyrotechnic items according to claim 1, characterized in that, The weighted score calculation process in step S2 is as follows: For each candidate kernel function k∈K, the leave-one-out cross-validation negative log-likelihood criterion and the maximum marginal likelihood criterion are computed in parallel. The formula for calculating the leave-one-out cross-validation negative log-likelihood criterion is as follows: Among them, X -i y -i To exclude the input set after the i-th data point, y i The output set after excluding the i-th data point; p(y i |X -i ,y -i (k) represents the Gaussian process regression model trained with the remaining data after excluding the i-th data point, applied to y. i The predicted probability; The formula for calculating the maximum marginal likelihood criterion is as follows: MLL(k) = logp(y|X,k), Where X = [x1, x2, ..., x n ] T It is a vector consisting of the logarithms of the stimulus quantities from n historical trials, y = [y1,...,y2]. n ] T For the corresponding binary response, K XX It is an n×n covariance matrix, K XX =k(x i ,x j I is calculated from the optimal kernel function selected in step S2. n It is an n×n identity matrix. It is the noise variance hyperparameter of the Gaussian process regression model; The leave-one-out cross-validation negative log-likelihood criterion and maximum marginal likelihood criterion for all candidate kernel functions are normalized and mapped to the [0,1] interval. The results are then standardized so that larger values indicate better performance of the Gaussian process regression model, and a weighted average is applied. The normalized formula for the leave-one-out cross-validation negative log-likelihood criterion is: Where min(NLL_LOO) is the minimum value of the leave-one-out cross-validation negative log-likelihood criterion for all candidate kernel functions, and max(NLL_LOO) is the maximum value of the leave-one-out cross-validation negative log-likelihood criterion for all candidate kernel functions. The normalized formula for the maximum marginal likelihood criterion is: Where MLL(k) is the maximum marginal likelihood criterion for the corresponding kernel function k, min(MLL) is the minimum value of the maximum marginal likelihood criterion for all candidate kernel functions, and max(MLL) is the maximum value of the maximum marginal likelihood criterion for all candidate kernel functions. The weighted scoring formula for kernel function k is as follows: S(k)=w·NLL_LOO (norm) (k)+(1-w)·MLL (norm) (k), Where w is a preset weighting coefficient, with a value range of 0.3-0.
7.
5. The method for testing the impact sensitivity of pyrotechnic items according to claim 4, characterized in that: In step S4, based on the Gaussian process regression model, the predicted ignition probability μ(x) of the logarithm of the stimulus at the new test point is calculated. * ) and the predicted variance of the probability of ignition σ 2 (x * ), Where, x * The logarithm of the stimulus amount at the new test site. y = [y1, y2, ..., y n ] T It is a vector composed of the binary response values of the corresponding n historical impact sensitivity tests. k * =[k(x * ,x1),k(x * ,x2),...,k(x * ,x n )] T It is the covariance vector between the test point x and all historical training points X. k(x * ,x * ) is x * Its own covariance, I n It is an n×n identity matrix. Using the cumulative distribution function Φ of the standard normal distribution as the link function, the ignition probability of the logarithm of the stimulus at the new test point is obtained as follows: This serves as a fitted curve for the ignition probability distribution of pyrotechnics during the current test cycle.
6. The method for testing the impact sensitivity of pyrotechnic items according to claim 1, characterized in that: In step S5, the critical ignition probability point includes the stimulus amount E with an ignition probability of 50%. 50 Stimulus level E with a 99.9% probability of triggering an ignition. 99.9 And the stimulus E with a 0.01% probability of igniting. 0.01 .
7. The method for testing the impact sensitivity of pyrotechnic items according to claim 1, characterized in that, The specific process in step S5 is as follows: Based on the logarithmic point sequence of dense stimulus quantities predicted by the Gaussian process regression model (x... i ,P i In this process, the ignition probability corresponding to the preset critical ignition probability point is taken as the target ignition probability, and the value P of the ignition probability is found. i and P i+1 Make the value of the target firing probability p target Given the probability of ignition P i and P i+1 Between, and then through Perform linear interpolation to obtain the target firing probability value p. target The corresponding logarithm of the stimulus quantity, and calculated using exponential operation E target =exp(x target Transform the logarithm of the stimulus amount back to the actual stimulus amount, where x target p represents the probability of the target firing. target The logarithm of the stimulus amount, E target This represents the actual target stimulus amount obtained through exponential mapping.
8. The method for testing the impact sensitivity of pyrotechnic items according to claim 1, characterized in that: In step S6, for each critical ignition probability point p, the relative rate of change Δ p The calculation formula is in, These are the estimated stimulus values for the key ignition probability points in this test cycle. This is the estimated stimulus amount for the critical ignition probability point in the previous test cycle.
9. The method for testing the impact sensitivity of pyrotechnic items according to claim 1, characterized in that: In step S6, after each preset number of impact sensitivity tests are completed, the test data is updated, the Gaussian process regression model is retrained, and a new fitting curve for the ignition probability distribution and the stimulus quantity estimate for key ignition probability points are generated.
10. A pyrotechnic impact sensitivity testing system for implementing the pyrotechnic impact sensitivity testing method as described in any one of claims 1-9, characterized in that, include: The data acquisition module is used to conduct impact sensitivity tests on pyrotechnics in sequence according to the sequential optimal experimental design method, record and acquire test data in real time, and preprocess the stimulus amount in the test data to obtain the logarithm of the stimulus amount. The first processing module is used to execute the kernel function dynamic selection process when the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area. The second processing module is used to sequentially call the third, fourth, and fifth processing modules to perform corresponding data processing after the number of tests reaches the preset activation threshold and the impact sensitivity test has entered the mixed result area, for each preset number of test cycles completed. The third processing module is used to take the logarithm of the stimulus amount in all historical test data as input, the corresponding binary response variable as output, and use the selected optimal kernel function to train the Gaussian process regression model to obtain the predicted value of the ignition probability and the predicted variance of the ignition probability of the logarithm of the stimulus amount at the new test point, as well as the fitting curve of the ignition probability distribution of the pyrotechnics in the current test cycle. The fourth processing module is used to calculate the estimated value of the logarithm of the stimulus quantity corresponding to a set of preset key ignition probability points by using piecewise linear interpolation based on the fitted curve of the ignition probability distribution, and to obtain the corresponding stimulus quantity estimate. The fifth processing module compares the estimated stimulus value of the critical ignition probability point obtained in this calculation with the estimated stimulus value of the critical ignition probability point corresponding to the previous test cycle, calculates the relative rate of change of the estimated value of each critical ignition probability point, and presets a personalized convergence threshold for each critical ignition probability point. If the relative rate of change of all critical ignition probability points is lower than their respective personalized convergence thresholds for multiple consecutive times, a test termination suggestion is generated; otherwise, the process returns to the third processing module.