CN115183969A - Method and system for estimating BWBN model parameters - Google Patents

Method and system for estimating BWBN model parameters Download PDF

Info

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
Application number
CN202210770940.3A
Other languages
Chinese (zh)
Inventor
赵露
丁孙玮
涂田刚
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shanghai Institute of Materials
Original Assignee
Shanghai Institute of Materials
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Shanghai Institute of Materials filed Critical Shanghai Institute of Materials
Priority to CN202210770940.3A priority Critical patent/CN115183969A/en
Publication of CN115183969A publication Critical patent/CN115183969A/en
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01MTESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
    • G01M7/00Vibration-testing of structures; Shock-testing of structures
    • G01M7/02Vibration-testing by means of a shake table
    • G01M7/022Vibration control arrangements, e.g. for generating random vibrations
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/11Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/13Differential equations
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/16Matrix 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

Method and system for estimating BWBN model parameters
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:
Figure BDA0003724195220000021
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:
Figure BDA0003724195220000022
wherein the content of the first and second substances,
Figure BDA0003724195220000023
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:
Figure BDA0003724195220000024
Figure BDA0003724195220000025
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:
Figure BDA0003724195220000031
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:
Figure BDA0003724195220000032
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,
Figure BDA0003724195220000033
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;
from
Figure BDA0003724195220000034
Generating index l at random according to probability
Figure BDA0003724195220000035
From
Figure BDA0003724195220000036
In selecting a sample
Figure BDA0003724195220000037
Will be provided with
Figure BDA0003724195220000038
Centered 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 is
Figure BDA0003724195220000039
Where r is a sample in a uniform distribution from 0 to 1, then let
Figure BDA00037241952200000310
And 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:
Figure BDA00037241952200000311
Figure BDA00037241952200000312
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 optimization
Figure BDA00037241952200000313
Number 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:
Figure BDA00037241952200000314
Figure BDA00037241952200000315
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,
Figure BDA0003724195220000041
and
Figure BDA0003724195220000042
the first and second derivatives of x are indicated,
Figure BDA0003724195220000043
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:
Figure BDA0003724195220000044
K 1 =f(x n ,y n )
Figure BDA0003724195220000045
Figure BDA0003724195220000046
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:
Figure BDA0003724195220000047
Figure BDA0003724195220000048
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,
Figure BDA0003724195220000049
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:
Figure BDA00037241952200000410
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:
Figure BDA0003724195220000061
Figure BDA0003724195220000062
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,
Figure BDA0003724195220000063
and
Figure BDA0003724195220000064
the first and second derivatives of x are indicated,
Figure BDA0003724195220000065
representing the first derivative of z.
Figure BDA0003724195220000066
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:
Figure BDA0003724195220000071
K 1 =f(x n ,y n )
Figure BDA0003724195220000072
Figure BDA0003724195220000073
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
Figure BDA0003724195220000074
Figure BDA0003724195220000081
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:
Figure BDA0003724195220000082
Figure BDA0003724195220000083
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,
Figure BDA0003724195220000084
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:
Figure BDA0003724195220000085
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:
Figure BDA0003724195220000091
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:
Figure BDA0003724195220000092
Figure BDA0003724195220000093
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:
Figure BDA0003724195220000094
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:
Figure BDA0003724195220000095
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,
Figure BDA0003724195220000096
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) From
Figure BDA0003724195220000101
Generating index l at random according to probability
Figure BDA0003724195220000102
From
Figure BDA0003724195220000103
In selecting samples
Figure BDA0003724195220000104
Will be provided with
Figure BDA0003724195220000105
Gaussian 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 is
Figure BDA0003724195220000106
Where r is a sample in a uniform distribution from 0 to 1, then let
Figure BDA0003724195220000107
And 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.
Figure BDA0003724195220000108
The chain indicating this sampling is erroneous and can be discarded and the index l regenerated if
Figure BDA0003724195220000109
Indicate 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:
Figure BDA00037241952200001010
Figure BDA00037241952200001011
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 optimization
Figure BDA00037241952200001012
Number 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
Figure BDA00037241952200001013
Figure BDA0003724195220000111
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:
Figure FDA0003724195210000011
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:
Figure FDA0003724195210000012
wherein, the first and the second end of the pipe are connected with each other,
Figure FDA0003724195210000013
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:
Figure FDA0003724195210000014
Figure FDA0003724195210000015
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:
Figure FDA0003724195210000016
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:
Figure FDA0003724195210000021
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,
Figure FDA0003724195210000022
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;
from
Figure FDA0003724195210000023
Generating index l at random according to probability
Figure FDA0003724195210000024
From
Figure FDA0003724195210000025
In selecting samples
Figure FDA0003724195210000026
Will be provided with
Figure FDA0003724195210000027
Centered 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 is
Figure FDA0003724195210000028
Where r is a sample in a uniform distribution from 0 to 1, then let
Figure FDA0003724195210000029
And 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:
Figure FDA00037241952100000210
Figure FDA00037241952100000211
wherein W represents the total number of parameters to be identified, p r Denotes sample acceptance rate, t r Indicating the target acceptance rate, na indicating the number of chains to be adjusted, i.e. in the j-th optimization
Figure FDA00037241952100000212
Number of markov chains.
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:
Figure FDA00037241952100000213
Figure FDA00037241952100000214
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,
Figure FDA0003724195210000031
and
Figure FDA0003724195210000032
the first and second derivatives of x are indicated,
Figure FDA0003724195210000033
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:
Figure FDA0003724195210000034
K 1 =f(x n ,y n )
Figure FDA0003724195210000035
Figure FDA0003724195210000036
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:
Figure FDA0003724195210000037
Figure FDA0003724195210000038
lnP 3 =lnP 1 +lnP 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,
Figure FDA0003724195210000039
represents the average value of x, and N () represents a gaussian distribution.
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:
Figure FDA00037241952100000310
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.
CN202210770940.3A 2022-06-30 2022-06-30 Method and system for estimating BWBN model parameters Pending CN115183969A (en)

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)

* Cited by examiner, † Cited by third party
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

Cited By (2)

* Cited by examiner, † Cited by third party
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
CN113259325B (en) Network security situation prediction method for optimizing Bi-LSTM based on sparrow search algorithm
CN113391211B (en) Method for predicting remaining life of lithium battery under small sample condition
CN109779791B (en) Intelligent diagnosis method for abnormal data in solid rocket engine
CN113687433B (en) Bi-LSTM-based magnetotelluric signal denoising method and system
CN111950711A (en) Second-order hybrid construction method and system of complex-valued forward neural network
CN115183969A (en) Method and system for estimating BWBN model parameters
CN110677297A (en) Combined network flow prediction method based on autoregressive moving average model and extreme learning machine
CN112086100B (en) Quantization error entropy based urban noise identification method of multilayer random neural network
CN105931130A (en) Improved ensemble Calman filter estimation method considering measurement signal loss
CN114819054A (en) Power electronic system state monitoring method based on physical information neural network
CN115982141A (en) Characteristic optimization method for time series data prediction
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
CN109217844B (en) Hyper-parameter optimization method based on pre-training random Fourier feature kernel LMS
CN116227324B (en) Fractional order memristor neural network estimation method under variance limitation
CN116304940A (en) Analog circuit fault diagnosis method based on long-short-term memory neural network
CN109088770B (en) Electromechanical system interactive network modeling method based on self-adaptive symbol transfer entropy
CN109474258B (en) Nuclear parameter optimization method of random Fourier feature kernel LMS (least mean square) based on nuclear polarization strategy
CN113642029B (en) Method and system for measuring correlation between data sample and model decision boundary
CN114372618A (en) Student score prediction method and system, computer equipment and storage medium
CN114048837A (en) Deep neural network model reinforcement method based on distributed brain-like map
CN111462479A (en) Traffic flow prediction method based on Fourier-recurrent neural network
CN112488295A (en) Method for optimizing storage life prediction of LSTM network relay by cross validation algorithm

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