CN115183969A - Method and system for estimating BWBN model parameters - Google Patents
Method and system for estimating BWBN model parameters Download PDFInfo
- Publication number
- CN115183969A CN115183969A CN202210770940.3A CN202210770940A CN115183969A CN 115183969 A CN115183969 A CN 115183969A CN 202210770940 A CN202210770940 A CN 202210770940A CN 115183969 A CN115183969 A CN 115183969A
- Authority
- CN
- China
- Prior art keywords
- identified
- parameter
- bwbn
- model
- value
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 33
- 238000009826 distribution Methods 0.000 claims abstract description 69
- 238000005259 measurement Methods 0.000 claims abstract description 36
- 238000002474 experimental method Methods 0.000 claims abstract description 20
- 238000013016 damping Methods 0.000 claims abstract description 7
- 238000005457 optimization Methods 0.000 claims description 54
- 239000000523 sample Substances 0.000 claims description 44
- 238000006073 displacement reaction Methods 0.000 claims description 24
- 238000004364 calculation method Methods 0.000 claims description 14
- 239000011159 matrix material Substances 0.000 claims description 12
- 230000015556 catabolic process Effects 0.000 claims description 7
- 238000006731 degradation reaction Methods 0.000 claims description 7
- 238000009827 uniform distribution Methods 0.000 claims description 5
- 230000003044 adaptive effect Effects 0.000 claims description 4
- 230000000694 effects Effects 0.000 claims description 4
- 108010074864 Factor XI Proteins 0.000 claims description 3
- 238000012546 transfer Methods 0.000 claims description 3
- 238000005070 sampling Methods 0.000 abstract description 12
- 238000004422 calculation algorithm Methods 0.000 description 6
- 230000008569 process Effects 0.000 description 4
- 238000005265 energy consumption Methods 0.000 description 3
- 230000007547 defect Effects 0.000 description 2
- 238000000605 extraction Methods 0.000 description 2
- 238000010606 normalization Methods 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 238000000342 Monte Carlo simulation Methods 0.000 description 1
- 230000009471 action Effects 0.000 description 1
- 230000032683 aging Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000009795 derivation Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000018109 developmental process Effects 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 238000003754 machining Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 238000013178 mathematical model Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000005312 nonlinear dynamic Methods 0.000 description 1
- 230000009022 nonlinear effect Effects 0.000 description 1
- 230000000737 periodic effect Effects 0.000 description 1
- 230000004044 response Effects 0.000 description 1
- 230000001629 suppression Effects 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
- 239000013598 vector Substances 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01M—TESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
- G01M7/00—Vibration-testing of structures; Shock-testing of structures
- G01M7/02—Vibration-testing by means of a shake table
- G01M7/022—Vibration control arrangements, e.g. for generating random vibrations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
- G06F17/13—Differential equations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Theoretical Computer Science (AREA)
- Software Systems (AREA)
- Databases & Information Systems (AREA)
- Algebra (AREA)
- General Engineering & Computer Science (AREA)
- Computing Systems (AREA)
- Operations Research (AREA)
- Complex Calculations (AREA)
Abstract
The invention relates to a method and a system for estimating BWBN model parameters, wherein the method comprises the following steps: acquiring measurement data and physical parameters of an experiment, substituting the measurement data and the physical parameters into a BWBN model, solving, and determining parameters to be identified; initializing the parameters to be identified, determining a value range and an initial value, supposing that each parameter to be identified meets Gaussian distribution, and respectively calculating a likelihood function to determine initial prior distribution of each parameter to be identified; respectively extracting samples from the prior distribution of each parameter to be identified, updating the estimated value of each parameter to be identified based on the samples, and repeating the step until a preset termination condition is met; and outputting the estimated value of each parameter to be identified. Compared with the prior art, the method takes the measurement data as the drive, estimates the parameters of the BWBN nonlinear hysteresis model based on the iTMCMC sampling method, has higher universality, and can be applied to the identification problem that various hysteresis damping vibrations approximately meet the BWBN model.
Description
Technical Field
The present invention relates to the field of vibration control, and in particular, to a method and a system for estimating BWBN model parameters.
Background
With the development of modern industrialization and society, the suppression or at least the attenuation of undesirable vibrations that may affect systems or structures is of great significance to the safety of building bridges, the precision of machining, the operational stability of equipment, the comfort of vehicle riding. Vibration control is typically achieved using passive, semi-active or active systems, each of which has considerable hysteretic properties. The Bouc-Wen-Baber-Noori (BWBN) model can be effectively used for analyzing and describing multi-system hysteresis performance and capturing a series of hysteresis cycle modes.
The Bouc-Wen-Baber-Noori hysteresis model is further developed on the basis of the classical Bouc-Wen hysteresis model, considering strength, stiffness degradation and pinching effects. Many engineered structures will enter a non-elastic state under dynamic loading and exhibit hysteresis characteristics. Hysteresis properties are also known as elastoplasticity. Hysteresis generally results from the non-linear properties of the material, the frictional properties of the contacting surfaces, and contact deformation between the engaging surfaces, among other things. A hysteretic relationship exists between the restoring force and the displacement of the structure under the action of load, when the load is periodic, the structure enters elastic-plastic deformation, and the displacement of the plastic deformation can not be completely restored, so that curves change along different paths to form hysteretic rings. The actual hysteresis curve of the structure is very complex and is difficult to be applied to analyzing the nonlinear characteristic of the structure, so that a hysteresis model which is convenient for mathematical description and can reflect the hysteresis characteristic of the structure needs to be established to predict the nonlinear dynamic response of the engineering structure.
The BWBN model improves the problem that the Bouc-Wen model cannot describe the pinch hysteresis phenomenon, and simultaneously has the problems that the parameter dimension is large, and a plurality of parameters in the model need to be determined when the BWBN model is required to describe the hysteresis performance of a system. The existing parameter estimation methods all need to use different filters, and the traditional MCMC method is easy to fall into a single-peak trap and has failure risk when facing high-latitude parameter estimation. TMCMC is suitable for the problem of few uncertain parameters, and if the number of uncertain parameters is large, the result may have a large deviation, and meanwhile, the efficiency of statistical estimation is low.
Disclosure of Invention
The present invention is directed to a method and system for estimating BWBN model parameters, which overcome the above-mentioned shortcomings of the prior art.
The purpose of the invention can be realized by the following technical scheme:
according to a first aspect of the present invention, there is provided a parameter estimation method of a BWBN model, including the steps of:
acquiring measurement data of an experiment, wherein the measurement data comprises force and displacement, determining physical parameters of the experiment, establishing a BWBN model, substituting the measurement data and the physical parameters into the BWBN model, solving the BWBN model, and determining W parameters to be identified in the BWBN model;
initializing each parameter to be identified respectively, determining the maximum value, the minimum value and the initial value of the parameter to be identified, supposing that each parameter to be identified meets Gaussian distribution, calculating the likelihood function of each parameter to be identified respectively to obtain initial prior distribution of each parameter to be identified, taking the initial prior distribution as the prior distribution of the parameter to be identified, and taking the initial value as the estimated value of the parameter to be identified;
respectively extracting samples from the prior distribution of each parameter to be identified, updating the estimation value of each parameter to be identified based on the samples, completing a round of optimization, repeating the step until a preset termination condition is met, and stopping iteration;
and outputting the estimated value of each parameter to be identified.
Further, in the optimization of the j-th round, samples are extracted from the prior distribution of the s-th parameter to be identified, and the estimated value of the s-th parameter to be identified is updated based on the samples, wherein j =1, 2, \8230, s =1, 2, \8230, and W specifically includes:
extracting N from the prior distribution of the s-th parameter to be identified j,s One sample:N j,s extracting a value for the sample of the s-th parameter to be identified in the preset j-th optimization;
and (4) calculating the optimized gradient coefficient of the j-th round, wherein the calculation formula is as follows:
wherein,representing the sample variation coefficient of the s-th parameter to be identified in the j-th optimization;
and respectively calculating the weighting coefficient of each sample, wherein the calculation formula is as follows:
representing the weighting coefficient of the kth sample of the s-th parameter to be identified in the j-th optimization, and D represents the measurement data of the experiment;
calculating N j,s The average value of the weighting coefficients of the samples is calculated according to the following formula:
wherein, se j,s Representing the sample weighting coefficient average value of the s-th parameter to be identified in the j-th optimization;
and (3) calculating a covariance matrix, wherein the calculation formula is as follows:
therein, sigma j,s N for representing s to-be-identified parameter in j optimization j,s The covariance matrix of the individual samples is,n for representing s to-be-identified parameter in j optimization j,s Mean value of samples, T denotes matrix transpose, ξ j Representing the self-adaptive scaling factor in the preset j-th round optimization;
fromGenerating index l at random according to probabilityFromIn selecting a sampleWill be provided withCentered gaussian distribution and covariance matrix sigma j,s As a hypothesis c To obtain a proposed distribution of c A Markov chain that starts as an initial sample distribution; if it isWhere r is a sample in a uniform distribution from 0 to 1, then letAnd updating the prior distribution and the initial value of the s-th parameter to be identified, otherwise, repeating the step.
Further, at the beginning of each round of optimization, an adaptive scaling factor ξ is determined j The calculation formula is as follows:
wherein W represents the total number of parameters to be identified, p r Denotes sample acceptance rate, t r Represents the target acceptance rate, na represents the number of strands to be adjusted,i.e. in the j-th optimizationNumber of markov chains of (c).
Further, the preset termination condition is the optimized gradient coefficient q of the j' th round j And if the value is equal to 1, considering that each parameter to be identified finds the optimal estimated value, and stopping iteration.
Further, the physical parameters of the experiment include mass and initial stiffness, and the BWBN model is as follows:
v(t)=1.0+δ v ε(t)
η(t)=1.0+δ η ε(t)
wherein m represents mass, x represents displacement, c represents damping coefficient, k e Representing stiffness, z damping displacement, F (t) force, h (z) pinch effect, v (t) strength degradation, η (t) stiffness degradation, ε (t) hysteresis dissipation energy, n, a, β, γ, δ v 、δ η As the parameter to be identified, the identification information is obtained,andthe first and second derivatives of x are indicated,representing the first derivative of z.
Further, a differential equation in the BWBN restoring force measurement model is solved using a fourth order longkutta as follows:
K 1 =f(x n ,y n )
K 4 =f(x n +h,y n +hK 3 )
wherein f () represents the differential equation of the BWBN model, x n For measuring displacements in the data, y n = F (t) force in the measurement data, h is the Longgusttal time step, K 1 ,K 2 ,K 3 ,K 4 Are transfer parameters in the dragon lattice library tower.
Further, the establishment principle of the likelihood function is as follows:
ln P 3 =ln P 1 +ln P 2
wherein σ 1 Representing the unknown variance, σ, of the prediction error of the displacement 2 An unknown variance representing the total dissipated energy prediction error,represents the average value of x, and N () represents a gaussian distribution.
Further, when calculating the likelihood function of each parameter to be identified, normalization processing needs to be performed on the predicted error dissipation energy, and the calculation formula is as follows:
wherein L = ln P 3 ,θ g For the neighbourhood of the composition of the parameters to be identified, the subscript g is the index of the forces and displacements in the measured data, N j,s And extracting a sample value of the s-th parameter to be identified in the preset j-th optimization.
Further, before the likelihood function is calculated, the total dissipated energy needs to be determined from the hysteresis curve of the measurement data.
According to a second aspect of the present invention, there is provided a parameter estimation system of a BWBN model, including:
the data acquisition module is used for acquiring measurement data and physical parameters of an experiment, wherein the measurement data comprises force and displacement;
the BWBN model module is used for establishing a BWBN model, substituting the measurement data and the physical parameters into the BWBN model, solving the BWBN model and determining W parameters to be identified in the BWBN model;
the initialization module is used for initializing each parameter to be identified, determining the maximum value, the minimum value and the initial value of the parameter to be identified, respectively calculating the likelihood function of each parameter to be identified on the assumption that each parameter to be identified meets Gaussian distribution to obtain initial prior distribution of each parameter to be identified, taking the initial prior distribution as the prior distribution of the parameter to be identified, and taking the initial value as the estimated value of the parameter to be identified;
the iteration optimization module is used for executing at least one round of optimization operation, stopping iteration when a preset termination condition is met, respectively extracting samples from the prior distribution of each parameter to be identified in each round of optimization operation, and updating the estimation value of each parameter to be identified based on the samples;
and the output module is used for outputting the estimation value of each parameter to be identified.
Compared with the prior art, the invention has the following beneficial effects:
the variable self-adaptive adjustment parameters are introduced, the sample weighting coefficients are adjusted after each Markov chain to be adjusted, so that the target distribution is adjusted, the performance of an algorithm is improved, the deviation in an estimated value is reduced, the estimated deviation of the average value and the variance value of the sample is small, the searching speed is accelerated, the sampling of the problem of multimodal or extremely flat unimodal can be effectively solved by using the advantages of the transition Markov chain Monte Carlo method, and meanwhile, the system has higher universality by applying the pinch-squeeze function of a BWBN model.
Drawings
FIG. 1 is a flow chart of the present invention;
FIG. 2 is a graph of normalized energy consumption in the example;
FIG. 3 is a probability distribution diagram of each parameter to be identified in the embodiment;
fig. 4 is a hysteresis graph in the example.
Detailed Description
The invention is described in detail below with reference to the figures and specific embodiments. The present embodiment is implemented on the premise of the technical solution of the present invention, and a detailed implementation manner and a specific operation process are given, but the scope of the present invention is not limited to the following embodiments.
Example 1:
aiming at the defects of the BWBN model parameter estimation method in the prior art, the applicant points out that the iTMMMC algorithm adaptively adjusts the proposed distribution and adjusts the sample weight after each MCMC step, so the method is superior to the original TMMMC method in the aspects of average value, standard deviation and effective sample number sampling, and can be better applied to the multi-dimensional parameter estimation problem, so the iTMMMC algorithm is improved, and the parameter estimation of the BWBN model is carried out by using the iTMMMC algorithm.
According to a first aspect of the present invention, there is provided a method for estimating parameters of a BWBN model, as shown in fig. 1, comprising the steps of:
s1, obtaining measurement data of an experiment, wherein the measurement data comprises force F (t) and displacement x, determining physical parameters of the experiment, establishing a BWBN model, substituting the measurement data and the physical parameters into the BWBN model, solving the BWBN model, and determining W to-be-identified parameters in the BWBN model;
in this embodiment, the measurement data, that is, the data acquired in the hysteresis force experiment, may also use numerical simulation data, and it should be noted that the hysteresis image represented by the measurement data or the numerical simulation data needs to conform to or approximate the BWBN model, so that the BWBN model can be used to describe the measurement data.
The BWBN restoring force measurement model is a relatively mature mathematical model, and the principle and the establishing process thereof are not described in detail, and can be understood by related practitioners. The BWBN resilience measurement model is established as follows:
v(t)=1.0+δ v ε(t)
η(t)=1.0+δ η ε(t)
wherein m represents mass, x represents displacement, c represents damping coefficient, k e Representing stiffness, z damping displacement, F (t) force, h (z) embodying the pinching effect, v (t) strength degradation, η (t) stiffness degradation, ε (t) hysteresis dissipation energy, n, a, β, γ, δ v 、δ η As the parameter to be identified, the identification information is obtained,andthe first and second derivatives of x are indicated,representing the first derivative of z.
I.e. differential equations in the BWBN model, the physical parameters in the experiment includeMass, initial stiffness, etc., will not be described in detail.
The known values and measured data are substituted into the BWBN model, which can be solved by a fourth order longstotta as follows:
K 1 =f(x n ,y n )
K 4 =f(x n +h,y n +hK 3 )
wherein f () represents a differential equation of the BWBN model, x n For measuring displacements in the data, y n = F (t) force in the measurement data, h is the Longgasta time step, K 1 ,K 2 ,K 3 ,K 4 Are transfer parameters in the dragon lattice tower.
S2, initializing each parameter to be identified respectively, determining the maximum value, the minimum value and the initial value of the parameter to be identified, supposing that each parameter to be identified meets Gaussian distribution, calculating the likelihood function of each parameter to be identified respectively to obtain initial prior distribution of each parameter to be identified, taking the initial prior distribution as the prior distribution of the parameter to be identified, and taking the initial value as the estimated value of the parameter to be identified;
after the fourth-order Runge Kutta solves the differential equation, it is equivalent to obtain the representation of each parameter to be identified, specifically, there are 17 parameters to be identified in the BWBN model, and 2 parameters are introduced to describe the input noise, so 19 parameters to be identified are totally obtained. In other application scenarios, different numbers of parameters to be identified, such as 7, 21, etc., may be set as required. According to experience, the value range and the initial value are set as shown in the following table:
TABLE 1 initialization of parameters to be identified
Calculating the inherent frequency of the system, and calculating a likelihood function to obtain the initial prior distribution of the unknown parameters to be identified on the assumption that the unknown parameters to be identified meet the Gaussian distribution; the specific establishing process of the likelihood function is as follows, a Gaussian distribution N () with the displacement mean value as a mean value and the total energy consumption prediction error as a variance is established, and the following likelihood function is established:
ln P 3 =ln P 1 +ln P 2
wherein σ 1 Representing the unknown variance, σ, of the prediction error of the displacement 2 Represents the unknown variance of the total dissipated energy prediction error,represents the average value of x, and N () represents a gaussian distribution;
and then, performing partial derivation to solve the parameter values, thereby obtaining the initial prior distribution of each parameter to be identified by calculating the likelihood function of each parameter to be identified.
When the likelihood function of each parameter to be identified is calculated, normalization processing needs to be performed on the prediction error dissipation energy, fig. 2 gives an example of a normalized energy consumption graph, and the calculation formula is as follows:
wherein L = ln P 3 ,θ g For the neighbourhood of the composition of the parameters to be identified, the subscript g is the index of the forces and displacements in the measured data, N j,s And extracting a value for the sample of the s-th parameter to be identified in the preset j-th optimization.
Before the likelihood function is calculated, the total dissipated energy needs to be determined from the hysteresis curve of the measurement data, and the calculation method is common knowledge in the art and is not described in detail.
S3, respectively extracting samples from the prior distribution of each parameter to be identified, updating the estimated value of each parameter to be identified based on the samples, completing a round of optimization, repeating the step until a preset termination condition is met, and stopping iteration;
taking the example of extracting samples from the prior distribution of the s-th parameter to be identified in the j-th optimization and updating the estimated value of the s-th parameter to be identified based on the samples, an iterative optimization process is described, wherein j =1, 2, \8230, s =1, 2, \8230, and W specifically includes:
(1) Extracting N from the prior distribution of the s-th parameter to be identified j,s One sample:N j,s for the sample extraction values of the s-th parameter to be identified in the preset j-th optimization, the sample extraction values of each parameter to be identified in each optimization can be different j,s =1000;
(2) Calculating the optimized gradient coefficient of the j 'th round, wherein the preset termination condition is the optimized gradient coefficient q of the j' th round j Equal to 1, it is considered that each parameter to be identified finds an optimal estimated value, and therefore, in this step, if the calculated gradient coefficient is 0, the iteration is stopped, the optimization is completed, and step S4 is executed. Further, a maximum number of iterations may be set as a termination condition. The formula for calculating the gradient coefficient is as follows:
q 0 =0,q j =min(|CV(τ j,s )-1|,s=1、2、...、W)
wherein, CV (τ) j,s ) The sample variation coefficient, min (| CV (τ) represents the s-th parameter to be identified in the j-th optimization j,s ) -1, s =1, 2,.. And W), i.e. to calculate the sample variation coefficient of the 19 parameters to be identified in the j-th round, and select the sample variation coefficient with the minimum difference from 1 as the optimized gradient coefficient in the j-th round, where the sample variation coefficient is a commonly used technical means in the art, and the calculation formula and principle are not described again, and the related practitioner can understand it;
(3) Calculating the weighting coefficient of each sample respectively, wherein the formula is as follows:
representing the weighting coefficient of the kth sample of the s-th parameter to be identified in the j-th optimization, and D represents the measurement data of the experiment;
(4) Calculating N j,s The average value of the weighting coefficients of the samples is calculated according to the following formula:
wherein, se j,s Representing the sample weighting coefficient average value of the s-th parameter to be identified in the j-th optimization;
(5) And (3) calculating a covariance matrix, wherein the formula is as follows:
therein, sigma j,s N for representing s to-be-identified parameter in j optimization j,s The covariance matrix of the individual samples is,n for representing s to-be-identified parameter in j optimization j,s Mean of samples, T denotes matrix transpose, ξ j Representing the self-adaptive scaling factor in the preset j-th round optimization;
(6) FromGenerating index l at random according to probabilityFromIn selecting samplesWill be provided withGaussian distribution and covariance matrix Sigma centered j,s As a hypothesis c To obtain a proposed distribution of c A Markov chain starting as an initial sample distribution; if it isWhere r is a sample in a uniform distribution from 0 to 1, then letAnd updating the prior distribution and the initial value of the s-th parameter to be identified, otherwise, repeating the step.
Generating a uniform distribution from 0 to 1, mainly for evaluating whether the probability distribution of this sampling is smaller than the distribution, if the probability of the uniform distribution is larger than the sample probability of the sampling, i.e.The chain indicating this sampling is erroneous and can be discarded and the index l regenerated ifIndicate that this sample can be accepted, followed by mu c And updating the prior distribution and the initial value of the s-th parameter to be identified.
And (4) after the steps (1) to (6) are completed, updating the prior distribution and the estimated value of the s-th parameter to be identified in the optimization of the j-th round, and then enabling s +1 to perform the steps (1) to (6) again to complete the optimization of the next parameter to be identified in the optimization of the j-th round. After the optimization of all 19 parameters to be identified is completed, because the preset termination condition is not met, j +1 is made to start the next round of iterative optimization.
In step (5), the adaptive scaling factor xi j Is determined at the beginning of each round of optimization, and the calculation formula is as follows:
wherein W represents the total number of the parameters to be identified, p r Denotes sample acceptance rate, t r Indicating the target acceptance rate and Na indicating the number of strands to be adjusted, i.e.in the j-th optimizationNumber of markov chains.
And S4, outputting the estimated value of each parameter to be identified and a probability density function pdf, wherein the probability distribution of the parameter to be identified in the embodiment is shown in FIG. 3.
In this embodiment, after multiple rounds of optimization, estimated values and variance values of 19 parameters to be identified are obtained, as shown in the following table:
TABLE 2 estimated values of parameters to be identified
In this embodiment, after obtaining the estimated values of each parameter in the BWBN model, the estimated values may be replaced into the BWBN model, and the measured data and the physical parameters of the experiment are substituted into the BWBN model again to solve the hysteresis curve of the experiment, and the simulated hysteresis curve graph is shown in fig. 4.
According to a second aspect of the present invention, there is provided a parameter estimation system of a BWBN model, including:
the data acquisition module is used for acquiring measurement data and physical parameters of an experiment, wherein the measurement data comprises force and displacement;
the BWBN model module is used for establishing a BWBN model, substituting the measurement data and the physical parameters into the BWBN model, solving the BWBN model and determining W parameters to be identified in the BWBN model;
the initialization module is used for initializing each parameter to be identified, determining the maximum value, the minimum value and the initial value of the parameter to be identified, respectively calculating the likelihood function of each parameter to be identified on the assumption that each parameter to be identified meets Gaussian distribution to obtain initial prior distribution of each parameter to be identified, taking the initial prior distribution as the prior distribution of the parameter to be identified, and taking the initial value as the estimated value of the parameter to be identified;
the iteration optimization module is used for executing at least one round of optimization operation, stopping iteration when a preset termination condition is met, respectively extracting samples from the prior distribution of each parameter to be identified in each round of optimization operation, and updating the estimation value of each parameter to be identified based on the samples;
and the output module is used for outputting the estimation value of each parameter to be identified.
The iTMCMC sampling method is based on a Bayesian model updating technology, and is used for estimating the parameters of a BWBN nonlinear hysteresis model, estimating multi-dimensional parameter vectors except for a non-sensitive hysteresis amplitude parameter A, and returning parameter posterior distribution, estimated values and variances. Compared with the MCMC sampling method, the method has the advantages that the sample weight after each MCMC sampling step is adjusted to reduce the average deviation of model evidence estimation, the aging period is applied in the MCMC step to improve the posterior approximation, and the recommended target distribution interval of the MCMC algorithm is selected in a self-adaptive mode to achieve the near-optimal acceptance rate. Therefore, the application overcomes the defects of the existing algorithm and simultaneously retains the advantages of MCMC sampling.
The application of the application can bring the following advantages:
the method has the advantages of large effective number of sampled independent samples, high sampling efficiency, small estimation deviation of the mean value and the variance value, capability of effectively solving the problem of multimodal or extremely flat unimodal sampling, wide application range, capability of well simulating hysteresis images with pinch hysteresis phenomena, simplicity and convenience in implementation, capability of realizing parameter estimation with noise influence without a filter, capability of providing probability distribution of BWBN model parameters, and great significance in evaluation of model selection.
The foregoing detailed description of the preferred embodiments of the invention has been presented. It should be understood that numerous modifications and variations could be devised by those skilled in the art in light of the present teachings without departing from the inventive concepts. Therefore, the technical solutions available to those skilled in the art through logic analysis, reasoning and limited experiments based on the prior art according to the concept of the present invention should be within the scope of protection defined by the claims.
Claims (10)
1. A method for estimating parameters of a BWBN model is characterized by comprising the following steps:
acquiring measurement data of an experiment, wherein the measurement data comprises force and displacement, determining physical parameters of the experiment, establishing a BWBN model, substituting the measurement data and the physical parameters into the BWBN model, solving the BWBN model, and determining W parameters to be identified in the BWBN model;
initializing each parameter to be identified respectively, determining the maximum value, the minimum value and the initial value of each parameter to be identified, supposing that each parameter to be identified meets Gaussian distribution, calculating the likelihood function of each parameter to be identified respectively to obtain initial prior distribution of each parameter to be identified, taking the initial prior distribution as the prior distribution of the parameter to be identified, and taking the initial value as the estimated value of the parameter to be identified;
respectively extracting samples from the prior distribution of each parameter to be identified, updating the estimated value of each parameter to be identified based on the samples, completing a round of optimization, repeating the step until a preset termination condition is met, and stopping iteration;
and outputting the estimated value of each parameter to be identified.
2. The method of claim 1, wherein in the j-th optimization, samples are extracted from the prior distribution of the s-th parameter to be identified, and the estimated value of the s-th parameter to be identified is updated based on the samples, wherein j =1, 2, \8230, s =1, 2, \8230, and W is specifically:
extracting N from the prior distribution of the s-th parameter to be identified j,s One sample:N j,s extracting a value for the sample of the s-th parameter to be identified in the preset j-th optimization;
and (4) calculating the optimized gradient coefficient of the j-th round, wherein the calculation formula is as follows:
wherein,representing the sample variation coefficient of the s-th parameter to be identified in the j-th optimization;
and respectively calculating the weighting coefficient of each sample, wherein the calculation formula is as follows:
representing the weighting coefficient of the kth sample of the s-th parameter to be identified in the j-th optimization, and D represents the measurement data of the experiment;
calculating N j,s The average value of the weighting coefficients of the samples is calculated according to the following formula:
wherein, se j,s Representing the sample weighting coefficient average value of the s-th parameter to be identified in the j-th optimization;
and (3) calculating a covariance matrix, wherein the calculation formula is as follows:
therein, sigma j,s N for representing s to-be-identified parameter in j round optimization j,s The covariance matrix of the individual samples is,n for representing s to-be-identified parameter in j optimization j,s Mean value of samples, T denotes matrix transpose, ξ j Representing the preset adaptive scaling factor in the jth round of optimization;
fromGenerating index l at random according to probabilityFromIn selecting samplesWill be provided withCentered gaussian distribution and covariance matrix sigma j,s As a hypothesis c Of the proposed distribution, to obtain μ c A Markov chain that starts as an initial sample distribution; if it isWhere r is a sample in a uniform distribution from 0 to 1, then letAnd updating the prior distribution and the initial value of the s-th parameter to be identified, otherwise, repeating the step.
3. The method of claim 2, wherein the adaptive scaling factor ξ is determined at the beginning of each round of optimization j The calculation formula is as follows:
4. The method of claim 2 for estimating the parameters of a BWBN model, whereinCharacterized in that the preset termination condition is the gradient coefficient q optimized in the j' th round j And if the value is equal to 1, considering that each parameter to be identified finds the optimal estimated value, and stopping iteration.
5. The method of claim 1, wherein the physical parameters of the experiment include mass and initial stiffness, and the BWBN model is as follows:
v(t)=1.0+δ v ε(t)
η(t)=1.0+δ η ε(t)
wherein m represents mass, x represents displacement, c represents damping coefficient, k e Representing stiffness, z damping displacement, F (t) force, h (z) pinch effect, v (t) strength degradation, η (t) stiffness degradation, ε (t) hysteresis dissipation energy, n, a, β, γ, δ v 、δ η As the parameter to be identified, the identification information is obtained,andthe first and second derivatives of x are indicated,representing the first derivative of z.
6. The parameter estimation method of the BWBN model according to claim 5, wherein the differential equation in the BWBN restoring force measurement model is solved using a fourth order longgure tower as follows:
K 1 =f(x n ,y n )
K 4 =f(x n +h,y n +hK 3 )
wherein f () represents a differential equation of the BWBN model, x n For measuring displacements in the data, y n = F (t) force in the measurement data, h is the Longgasta time step, K 1 ,K 2 ,K 3 ,K 4 Are transfer parameters in the dragon lattice tower.
7. The method of claim 5 wherein the likelihood function is built based on the following principles:
lnP 3 =lnP 1 +lnP 2
8. The method of claim 7, wherein the likelihood function of each parameter to be identified is calculated by normalizing the predicted error dissipation energy, and the calculation formula is:
wherein, L = lnP 3 ,θ g For the neighbourhood of the composition of the parameters to be identified, the subscript g is the index of the forces and displacements in the measured data, N j,s And extracting a sample value of the s-th parameter to be identified in the preset j-th optimization.
9. The method of claim 8, wherein the total dissipated energy is determined from a hysteresis curve of the measured data before the likelihood function is calculated.
10. A parameter estimation system of a BWBN model, characterized in that the parameter estimation method based on the BWBN model according to any one of claims 1-9 comprises:
the data acquisition module is used for acquiring measurement data and physical parameters of an experiment, wherein the measurement data comprises force and displacement;
the BWBN model module is used for establishing a BWBN model, substituting the measurement data and the physical parameters into the BWBN model, solving the BWBN model and determining W parameters to be identified in the BWBN model;
the initialization module is used for initializing each parameter to be identified, determining the maximum value, the minimum value and the initial value of the parameter to be identified, respectively calculating the likelihood function of each parameter to be identified on the assumption that each parameter to be identified meets Gaussian distribution to obtain initial prior distribution of each parameter to be identified, taking the initial prior distribution as the prior distribution of the parameter to be identified, and taking the initial value as the estimated value of the parameter to be identified;
the iteration optimization module is used for executing at least one round of optimization operation, stopping iteration when a preset termination condition is met, respectively extracting samples from the prior distribution of each parameter to be identified in each round of optimization operation, and updating the estimation value of each parameter to be identified based on the samples;
and the output module is used for outputting the estimation value of each parameter to be identified.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210770940.3A CN115183969A (en) | 2022-06-30 | 2022-06-30 | Method and system for estimating BWBN model parameters |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210770940.3A CN115183969A (en) | 2022-06-30 | 2022-06-30 | Method and system for estimating BWBN model parameters |
Publications (1)
Publication Number | Publication Date |
---|---|
CN115183969A true CN115183969A (en) | 2022-10-14 |
Family
ID=83514573
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202210770940.3A Pending CN115183969A (en) | 2022-06-30 | 2022-06-30 | Method and system for estimating BWBN model parameters |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN115183969A (en) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117216846A (en) * | 2023-09-12 | 2023-12-12 | 华南理工大学 | Reinforced concrete member hysteresis curve prediction method, system, equipment and medium |
-
2022
- 2022-06-30 CN CN202210770940.3A patent/CN115183969A/en active Pending
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117216846A (en) * | 2023-09-12 | 2023-12-12 | 华南理工大学 | Reinforced concrete member hysteresis curve prediction method, system, equipment and medium |
CN117216846B (en) * | 2023-09-12 | 2024-04-19 | 华南理工大学 | Reinforced concrete member hysteresis curve prediction method, system, equipment and medium |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
KR101958674B1 (en) | Actually-measured marine environment data assimilation method based on sequence recursive filtering three-dimensional variation | |
CN111783362B (en) | Method and system for determining residual service life of electric gate valve | |
CN110334580A (en) | The equipment fault classification method of changeable weight combination based on integrated increment | |
CN114819054B (en) | Power electronic system state monitoring method based on physical information neural network | |
CN110286586A (en) | A kind of MR damper hybrid modeling method | |
CN108734287A (en) | Compression method and device, terminal, the storage medium of deep neural network model | |
CN113391211B (en) | Method for predicting remaining life of lithium battery under small sample condition | |
CN109088749A (en) | The method for estimating state of complex network under a kind of random communication agreement | |
CN114897144A (en) | Complex value time sequence signal prediction method based on complex value neural network | |
CN113791351B (en) | Lithium battery life prediction method based on transfer learning and difference probability distribution | |
CN111950711A (en) | Second-order hybrid construction method and system of complex-valued forward neural network | |
CN108734268A (en) | Compression method and device, terminal, the storage medium of deep neural network model | |
CN110677297A (en) | Combined network flow prediction method based on autoregressive moving average model and extreme learning machine | |
CN115183969A (en) | Method and system for estimating BWBN model parameters | |
CN116227324B (en) | Fractional order memristor neural network estimation method under variance limitation | |
CN110962828A (en) | Method and equipment for predicting brake pressure of electric automobile | |
CN115963420A (en) | Battery SOH influence factor analysis method | |
CN116665483A (en) | Novel method for predicting residual parking space | |
CN115982141A (en) | Characteristic optimization method for time series data prediction | |
CN118036809A (en) | Fault current prediction method and medium based on snow ablation optimization cyclic neural network | |
CN118193978A (en) | Automobile roadblock avoiding method based on DQN deep reinforcement learning algorithm | |
CN113537638A (en) | Short-term wind pressure prediction method and abnormal data completion method and device for high-rise building | |
CN109217844B (en) | Hyper-parameter optimization method based on pre-training random Fourier feature kernel LMS | |
CN111462479A (en) | Traffic flow prediction method based on Fourier-recurrent neural network | |
CN109474258B (en) | Nuclear parameter optimization method of random Fourier feature kernel LMS (least mean square) based on nuclear polarization strategy |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
CB02 | Change of applicant information | ||
CB02 | Change of applicant information |
Address after: 200437 No. 99, Handan Road, Shanghai, Hongkou District Applicant after: Shanghai Material Research Institute Co.,Ltd. Address before: 200437 No. 99, Handan Road, Shanghai, Hongkou District Applicant before: SHANGHAI Research Institute OF MATERIALS |