CN106682277A - Fast simulation method for non-stationary random process - Google Patents
Fast simulation method for non-stationary random process Download PDFInfo
- Publication number
- CN106682277A CN106682277A CN201611106931.5A CN201611106931A CN106682277A CN 106682277 A CN106682277 A CN 106682277A CN 201611106931 A CN201611106931 A CN 201611106931A CN 106682277 A CN106682277 A CN 106682277A
- Authority
- CN
- China
- Prior art keywords
- formula
- simulation
- random
- wave
- frequency
- 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
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
The invention discloses a fast simulation method for a non-stationary random process. The fast simulation method includes the steps: firstly, establishing the relationship between the multi-point non-stationary random process and a one-dimensional non-stationary random wave, namely, acquiring the two-dimensional evolutionary power spectrum density of a transformation random wave according to the evolutionary power spectrum and the coherence function of a given target random process; secondly, decoupling 2D EPSD by a proper orthogonal decomposition method; finally, simulating the transformation random wave by two-dimensional FFT (fast Fourier transform) to obtain the corresponding target non-stationary random process. Cholesky decomposition is omitted, decoupling time based on POD (proper orthogonal decomposition) is shortened, and the two-dimensional FFT technology is used, so that simulation efficiency is greatly improved. In addition, when parameters are properly selected, precision is still is high. In conclusion, the fast simulation method has the advantages of easiness in use and high precision and simulation efficiency, and can effectively solve the problem that the simulation efficiency of a spectral representation method is low when non-stationary random process simulation points are large.
Description
Technical field
The invention belongs to stochastic signal simulation field, and in particular to it is non-flat that a kind of mixing random wave and Proper Orthogonal are decomposed
Steady stochastic process rapid simulation method.
Background technology
Fluctuating wind speed and earthquake motion typically exhibit non-stationary characteristic, therefore available nonstationary random process description.1965
The evolution power spectrum concept (evolution PSD or EPSD) that year Priestley is proposed is widely used in describing these arbitrary excitations (goes
Except average assay).Compared to time series models and the time varying spectrum based on wavelet transformation, EPSD reflects spectrum energy with the time
Change, it is often more important that with more obvious physical significance.In addition, EPSD makes frequency domain response analyze and nonstationary random process
Simulation become to be more prone to.
Compared to time domain approach such as time-varying time series models and random trigonometric function multinomial, spectral representation method (SRM)
Because its accuracy and simplicity are widely used in the simulation of stochastic process.However, the time-varying characteristics of EPSD make to simulate
Two difficulties are showed.On the one hand, spectral power matrix will carry out a Cholesky decomposition at each moment and frequency.2013
Huang etc. constant when having separated the phase component of spectrum matrix and being derived delay coherence matrix, spectrum is greatly reduced for this
The calculating demand of matrix decomposition.On the other hand, because the simulation that can not use fast Fourier transform, nonstationary random process is imitated
Rate is more low.For this purpose, some scholars it is also proposed the non-stationary analogy method that can use fft algorithm.Li and Kareem exist
The method for proposing Stochastic Decomposition in 1991.Subsequently, they propose a kind of discrete Fourier transform and the mixing of digital filtering
Method.2014, Huang is proposed can use the small echo of FFT and the mixed method of SRM.Recently, it is yellow to be decomposed based on Proper Orthogonal
(POD) non-stationary rapid simulation method, not only efficiency and precision are all very high for the method, and are easy to use.
Although the simulation precision of nonstationary random process has obtained certain raising, to work as when simulation counts very big still
Need further to improve simulation precision.Wherein typical case is then to hit under simulation Loads of Long-span Bridges and the such as thunderstorm on power transmission line
The non-stationary wind fields such as stream, mountain area typhoon and typhoon.At present, what the main span of Loads of Long-span Bridges was built is increasingly longer.For example, build
Japanese alum strait Bridge and the main span of Messina strait Bridge built of plan respectively reach 1991 and 3300m.For many
Across power transmission line, its length is up to 8000m.In these cases, in order in the wind of time domain response analysis and wind Vehicle-Bridge Coupling System
The gentle dynamic non-line of structure is considered in vibration response analysis, then needs to simulate hundreds of simulation point.
It is well known that multipoint random process and space-time random field or random wave have closely contact.It is random for space-time
For field or random wave, time and spectrum change are all taken into account.In this case, can using steady homogenizing random wave or
The analogy method simulation multiple spot stationary random process of space-time random field.Due to the use avoided with Two-dimensional FFT that spectrum matrix decomposes,
This method can greatly improve simulation precision.Equally, said method can expand to the simulation of multiple spot nonstationary random process.
Deodatis and Shinozuka proposes the analogy method based on the non-stationary random wave of SRM within 1989.However, due to needing
Multiple spot nonstationary random process is transformed into into non-stationary random wave, the analogy method can not be directly applied for the random mistake of non-stationary
The simulation of journey.Further, since two dimension EPSD is over time and distance change, therefore above-mentioned analogy method can not use Two-dimensional FFT.
The content of the invention
The technical problem to be solved is to provide a kind of nonstationary random process rapid simulation method, which obviates
Cholesky decomposes, reduces the time based on POD decouplings and used 2D FFT techniques, therefore simulation precision to obtain greatly
Raising.
To solve above-mentioned technical problem, the technical solution used in the present invention is:
A kind of nonstationary random process rapid simulation method, comprises the following steps:
Step 1:The target power spectrum of given multiple spot nonstationary random process and coherent function
IfFor the multiple spot nonstationary random process of zero-mean, T is represented and turned
Put, its jth item elementIt is shown below:
In formula:ω is circular frequency;Aj(ω t) is gradual modulation function;Zj(ω) it is complex orthogonal increment that average is zero
Stochastic process;
The evolution spectral power matrix of the nonstationary random process is shown below:
Wherein matrix each element is expressed as:
In formula:ForCorresponding steady bilateral power spectrum function, γjk(ω) it isWithIt
Between plural coherent function, and with following property
Wherein * represents complex conjugation;Accordingly,WithCross-correlation functionIt is expressed as:
In formula:τ is time delay;
Step 2:Obtain the two-dimentional evolution spectrum density of conversion random wave
Along multiple spot nonstationary random process p of spatial axes x distribution0T () is considered along spatial axes x1,x2,…,xnPoint
The one-dimensional discrete non-stationary random wave [f of cloth0(x1,t),f0(x2,t),…,f0(xn,t)]T, according to formula (5), γ (ξ, ω) has
Following characteristic:
γ (ξ, ω)=γ*(-ξ,ω);γ (ξ, ω)=γ*(ξ,-ω) (7)
It is provided with following Fourier transform pairs to exist:
In formula:β (κ, ω) is the Fourier transformation of γ (ξ, ω), while being also non-negative and the function with regard to origin symmetry;
Formula (9) is substituted into into formula (6), then:
For conversion random wave f0(x, t), its two-dimensional autocorrelation function representation is:
Wherein:To convert the 2D EPSD of random wave;Show that the 2D EPSD for converting random wave are:
Step 3:Two-dimentional evolution power spectrum based on POD is decoupled
First by G0(x, ω, t) carry out discrete in space, frequency domain and time domain, and are reconstructed into following matrix:
In formula:xj=j Δs x corresponds to locus and j=1,2 ..., n is Spatial Cable argument;It is frequency
Vector and n2=0,1 ..., N2- 1 is frequency indices number;It is time and nt=1,2 ..., NtIt is time index number;Square
Battle array It is respectively 1 × N with the dimension of Qt, 1 × nNtAnd N2×nNt;
Orthogonal basiss are determined using Eigenvalues Decomposition method, is shown below:
RΦq=λqΦq, q=1,2 ..., N2 (16)
In formula:R is the covariance matrix of matrix Q and is expressed as form:
Wherein:λqAnd ΦqThe respectively eigen vector of q items, and satisfaction
It is by matrix Q approximate representations:
In formula:NqIt is to be superimposed the number of item and meet Nq< N2;hqIt is q item principal coordinates;The orthogonal spy of feature based vector
Property, hqIt is expressed as:
Formula (18) is expressed as again from the angle of element:
In formula:WithIt is ΦqN-th2Individual element and hqR-th element;Wherein r=(j-1) Nt
+nt;
The conitnuous forms of formula (20) are expressed as:
In formula:aq(x, t) is the real-number function of room and time;Φq(ω) be frequency even function;
Formula (21) is substituted into into formula (12), random wave decomposition two dimension EPSD is converted and is expressed as:
In formula:
Step 4:The Two-dimensional FFT efficient simulation equidistantly put
By two-dimentional evolution power spectral densityThe simulation of one-dimensional non-stationary random wave be converted into one and be
Arrange steady homogenizing random waveSimulation;Correspond toSimulation formula be shown below:
Therefore, sample f (x, t) of original non-stationary random wave is determined by following formula:
Formula (23) is expressed as again:
In formula:Re represents treating excess syndrome number;p1=0,1 ..., m × M1- 1 and p2=0,1 ..., m × M2-1;q1And q2Respectively
p1/M1And p2/M2Remainder;WithRespectively:
Wherein:WithIt is expressed as:
Δ x and Δ t should meet Δ x Δs κ=2 π/M1With Δ t Δs ω=2 π/M2, while meeting:
M1≥2N1;M2≥2N2。 (30)
Further, in step 2, during drawing the 2D EPSD of conversion random wave, cut-off wave number and frequency need to cover
The discretization points of the main energetic of two-dimentional evolution power spectral density, wave number and frequency can catch pulse spectrum.
Further, in step 2, when pulse spectrum is located at origin, simulation precision is improved using frequency shift (FS), is used
Frequency shift (FS), its simulation formula be:
In formula:It is frequency shift parameters;Take positive integer;n1=0,1 ..., N1-1;n2=0,1 ..., N2-1;WithIt is separate and is uniformly distributed in the random of [0,2 π]
Phase angle;WhenWithAfter determining respectively, obtain converting the sample Equation f of random wave(i)(x,t)。
Further, when improving simulation precision using frequency shift (FS), space increments Δ x and incremental time Δ t needs to meetWith Δ t≤2 π/(2 ωu)。
Compared with prior art, the invention has the beneficial effects as follows:Cholesky commonly used in traditional method is not needed
Decompose, and the time decoupled based on POD also obtains a certain degree of reduction.Furthermore it is also possible to great using Two-dimensional FFT technology
Improve simulation precision.The numerical example analysis shows, the inventive method has that easy to use, precision is higher and simulation precision is very high
The characteristics of.The inventive method can effectively solving nonstationary random process simulation points it is very big when, spectral representation method exist simulation
The problem of inefficiency.
Description of the drawings
Fig. 1 is the decomposition EPSD in example 1.
Fig. 2 is the two-dimentional Evolutionary Spectra reconstruct (G in example 1 based on POD0(x200, ω, t)) comparison diagram.
Fig. 3 is the two-dimentional Evolutionary Spectra reconstruct (G in example 1 based on POD0(x201, ω, t)) comparison diagram.
Fig. 4 is the sample time-histories of simulation point 200 in example 1
Fig. 5 is the sample time-histories of simulation point 201 in example 1
Fig. 6 is the comparison diagram of the auto-correlation function of simulation point 200 and theoretical value in example 1.
Fig. 7 is the comparison diagram of the auto-correlation function of simulation point 201 and theoretical value in example 1.
Fig. 8 is the comparison diagram of the cross-correlation function in example 1 between simulation point 200 and 201 and theoretical value.
Specific embodiment
With reference to the accompanying drawings and detailed description the present invention is further detailed explanation.The present invention is used for non-flat
Steady stochastic process carries out efficient, high-precision Fast simulation.Initially set up multiple spot nonstationary random process and one-dimensional non-stationary with
Contact between machine ripple, i.e., according to the evolution power spectrum (EPSD) and coherent function of given target stochastic process, draw conversion with
The two-dimentional evolution power spectral density (2D EPSD) of machine ripple;Then the method for decomposing (POD) using Proper Orthogonal is entered to 2D EPSD
Row decoupling;Finally it is simulated to converting random wave using two-dimensional fast fourier transform (2D FFT), so as to obtain corresponding mesh
The simulation time-histories of mark nonstationary random process, it is concretely comprised the following steps:
1) the target power spectrum and coherent function of multiple spot nonstationary random process are given
IfFor the multiple spot nonstationary random process of zero-mean, wherein T is represented
Transposition its jth item elementIt is shown below:
In formula:ω is circular frequency;Aj(ω t) is gradual modulation function;Zj(ω) it is complex orthogonal increment that average is zero
Stochastic process.For simplicity, it is assumed that Aj(ω t) is real number.
The evolution spectral power matrix of the nonstationary random process is shown below:
Wherein matrix each element is expressed as:
In formula:ForCorresponding steady bilateral power spectrum function;γjk(ω) it isWithIt
Between plural coherent function and with following property
Wherein * represents complex conjugation, accordingly,WithCross-correlation functionIt is expressed as:
In formula:τ is time delay.
2) the two-dimentional evolution spectrum density of conversion random wave is obtained
When stochastic process be located at straight line on when, above-mentioned multiple spot nonstationary random process and one-dimensional non-stationary random wave it
Between there is close contact.Specifically, along multiple spot nonstationary random process p of spatial axes x distribution0T () can be considered
Along spatial axes x1,x2,…,xnThe one-dimensional discrete non-stationary random wave [f of distribution0(x1,t),f0(x2,t),…,f0(xn,t)]T。
Therefore, the simulation along the multiple spot nonstationary random process of straight line can be converted into the mould of corresponding one-dimensional non-stationary random wave
Intend.
In actual applications, coherent function is often expressed as the function apart from ξ and frequencies omega.Now, it is determined that conversion random wave
2DEPSD become to be relatively easy to.According to formula (5), γ (ξ, ω) has following characteristic:
γ (ξ, ω)=γ*(-ξ,ω);γ (ξ, ω)=γ*(ξ,-ω) (7)
It is provided with following Fourier transform pairs to exist:
In formula:β (κ, ω) is the Fourier transformation of γ (ξ, ω), while being also non-negative and the function with regard to origin symmetry.
Formula (9) is substituted into formula (6), above-mentioned formula is expressed as:
For conversion random wave f0(x, t), its two-dimensional autocorrelation function representation is:
Wherein:To convert the 2D EPSD of random wave.Contrast equation (10) and (11), it can be deduced that turn
Change the 2D EPSD of random wave, be shown below:
Obviously, the two-dimentional evolution power spectral density for converting random wave is also non-negative and with regard to the function of origin symmetry.
It is similar with the analogy method of the steady homogenizing stochastic process of multiple spot based on space-time random field, above-mentioned two-dimentional evolution power
The discretization of spectrum density also can bring certain impact to the simulation precision of multipoint random process.In order to improve simulation precision, cut
Only wave number and frequency should greatly arrive the main energetic that can cover two-dimentional evolution power spectral density.In addition, wave number and frequency from
Dispersion points also should greatly be arrived and can catch pulse spectrum.When pulse spectrum is located at origin, can be further using frequency shift (FS)
Improve simulation precision.If usage frequency offsets, its simulation formula is as follows:
In formula:It is frequency shift parameters;Take positive integer;n1=0,1 ..., N1-1;n2=0,1 ..., N2-1;WithIt is separate and is uniformly distributed in the random of [0,2 π]
Phase angle.WhenWithAfter determining respectively, the sample Equation f of random wave is converted(i)(x, t) is obtained.Need note
Meaning, in order to avoid the aliasing in sampling thheorem, space increments Δ x and incremental time Δ t must is fulfilled for requiringWith Δ t≤2 π/(2 ωu)。
Once simulate discrete one-dimensional random ripple [f (x1,t),f(x2,t),…,f(xn,t)]TSample, then also obtain
Sample [the p of multiple spot nonstationary random process1(t),p2(t),…,pn(t)]T.Its conversion formula is pj(t)=f (xj, t), j=
1,2,...,n.The Cholesky that can be seen that classics by above-mentioned simulation formula decomposes and can avoid.But, due to two dimension evolution
The time of spectrum density and spatial dependence, simulation formula still can not improve simulation precision using two-dimensional FFT.
3) the two-dimentional evolution power spectrum based on POD is decoupled
In order that improving simulation precision with Two-dimensional FFT, it is necessary first to the time in the two-dimentional EPSD of conversion random wave and sky
Between variable separate.When coherent function is the function away from discrete frequency, time and space variable can be carried out using POD
Separate.By formula (12) as can be seen that the Temporal-Spatial Variables separated in 2D EPSD only need to decoupling
.For concise explanation, it is considered to the situation of ω > 0.
When being decomposed using POD, first by G0(x, ω, t) carry out discrete in space, frequency domain and time domain, and are reconstructed into
Following matrix:
In formula:xj=j Δs x corresponds to locus and j=1,2 ..., n is Spatial Cable argument;It is frequency
Vector and n2=0,1 ..., N2- 1 is frequency indices number;It is time and nt=1,2 ..., NtIt is time index number.It is aobvious
So, matrix It is respectively 1 × N with the dimension of Qt, 1 × nNtAnd N2×nNt。
In order that matrix Q is in a series of ideal quadrature basesOn projection reach maximum, adopt
The method of Eigenvalues Decomposition determines these orthogonal basiss, is shown below:
RΦq=λqΦq, q=1,2 ..., N2 (17)
In formula:R is the covariance matrix of matrix Q and is expressed as form:
Wherein:λqAnd ΦqThe respectively eigen vector of q items, and satisfaction
If order restructuring of the eigenvalue to successively decrease, less eigenvalue is by comprising most energy.Now, square
Battle array Q approximate representations be:
In formula:NqIt is to be superimposed the number of item and meet Nq< N2;hqIt is q item principal coordinates.The orthogonal spy of feature based vector
Property, hqIt is expressed as:
Meanwhile, it is nN that it is a dimensiontRow vector.
Formula (18) can again be expressed as from the angle of element:
In formula:WithIt is ΦqN-th2Individual element and hqR-th element;Wherein r=(j-1) Nt
+nt.Accordingly, the conitnuous forms of formula (20) are expressed as:
In formula:aq(x, t) is the real-number function of room and time;Φq(ω) be frequency even function.
Formula (21) is substituted into formula (12), random wave decomposition two dimension EPSD is converted and is expressed as:
In formula:
The two dimension time of EPSD is separated using POD and the effectiveness of space variable depends primarily on each EPSD with regard to the time
With the function feature of frequency variable, or characteristic distributions of the different EPSD along spatial variations.When each EPSD is in frequency
Aspect overall variation similar trend is coordinated in time, and dimensionality reduction effect will be very notable.Similar, if along spatial variations
Different EPSD have identical variation tendency, and dimensionality reduction effect equally will highly significant.
In addition, the decoupling of above-mentioned evolution spectral density function equally can separate frequency variable from time and space variable
Out, wherein characteristic vector be time and space variable function and coordinate function be frequency function.Therefore, its physics meaning
Justice is consistent with the conventional covariance analysis of POD.However, it is nN that this process is needed to dimensiont×nNtMatrix carry out spy
Levy vector and decompose (usual n > 100, Nt> 1000).And it is N that the isolation illustrated in the present invention is only needed to dimension2×N2's
Matrix is decomposed (usual N2=1024 or 2048).Obviously, the decomposition of such matrix is bad realization, even so,
This decomposition is in some other cases effective.
4) the Two-dimensional FFT efficient simulation equidistantly put
Once it is determined that formula (22), two-dimentional evolution power spectral densityOne-dimensional non-stationary random wave
Simulation can be converted into a series of steady homogenizing random wavesSimulation.Correspond toSimulation it is public
Formula is shown below:
Therefore, sample f (x, t) of original non-stationary random wave is determined by following formula:
In order that with Two-dimensional FFT, formula (23) should be expressed as again:
In formula:Re represents treating excess syndrome number;p1=0,1 ..., m × M1- 1 and p2=0,1 ..., m × M2-1;q1And q2Respectively
p1/M1And p2/M2Remainder;WithIt is as follows:
Wherein:WithIt is expressed as:
In above-mentioned equation, Δ x and Δ t should meet Δ x Δs κ=2 π/M1With Δ t Δs ω=2 π/M2.In order to avoid aliasing,
Following conditions must are fulfilled for:
M1≥2N1;M2≥2N2 (31)
It should be noted that working as N1≤n1≤M1- 1 and (or) N2≤n2≤M2When -1, in formula (31)WithNeed to set up, this is because when wave number and frequency are beyond cut-off wave number and cut-off frequency, two-dimentional evolution power spectrum
Spectrum is assumed to be 0.
The present invention and effect of the present invention are further detailed below by specific example.
Example 1:Non-stationary heterogeneous body ground motion simulation
If the distance of two earthquake synthesis points is 10m, and simulates 401 points altogether.The evolution power spectrum of earthquake is adopted
Inseparable form, its formula is shown below:
In formula:S0=1.For simulation point j=1,2 ..., 200, parameter Cj=2.5 and Dj=0.15.For simulation point j=
201,202 ..., 401, parameter Cj=4 and Dj=0.2.Fig. 1 shows and decompose under four selected frequencies spectrum G0(x200, ω, t) and G0
(x201, ω, curve of cyclical fluctuations t).As can be seen that simulation point has obvious difference between 200 and 201 Evolutionary Spectra.Base
POD methods in the present invention carry out the decoupling and reconstruct of evolution power spectrum, work as NqWhen=6, reconstruct Evolutionary Spectra develops with target
Contrast between spectrum is as shown in Figure 1.As can be seen that there is very high precision based on the reconstruct of POD.
In the implementation case, using Harichandran-Vanmarcke ' the s coherency function models for considering row wave effect
In formula:A=0.626;α=0.022;υapp=1000m/s;θ (ω) is frequency dependence distance, and is expressed from the next:
θ (ω)=d [1+ (ω/ω0)b]-1/2 (34)
Wherein:D=19700m, ω0=12.692rad/s and b=3.47.
Obviously, the coherent function is only the function away from discrete frequency, and meets the conversion condition illustrated in the present invention;Accordingly
, β (κ, ± ω) is expressed as follows:
In formula:WithObviously, β (κ, ω) be with regard to
Origin symmetry but be not quadrant symmetrical function.Based on formula (27), can obtain converting the two-dimentional EPSD of random wave.Therefore,
The simulation of multi-point Ground Motion process can be realized by converting the simulation of random wave.
In simulations, following parameter has specifically been selected:Cut-off wave number κu=π/10rad/m;Cut-off frequency ωu=20 π rad/
s;Discrete total wave number N1=512 and discrete frequency number N2=1024;M1=2N1And M2=2N2;Therefore, space interval Δ x=10m;
Time interval Δ t=0.05s;Frequency shift parametersTherefore, simulation point sum is 4096, and the persistent period of analog sample
For 204.8s.
Fig. 4 and Fig. 5 respectively illustrate a sample time-histories of simulation point 200 and 201.Fig. 6, Fig. 7 and Fig. 8 are by 2000
Auto-correlation and cross-correlation function and comparison diagram of the desired value in tri- time delays in τ=0,0.1 and 0.2s that sample is estimated.
Wherein constantWithIt is used for
Normalization correlation function, its concrete numerical value is respectively 15.42,21.36 and 14.98.As can be seen that the present invention has very high mould
Intend precision.
Example 2:Simulation precision is contrasted
POD is applied to Huang in 2015 the spectral representation method of nonstationary random process, it is proposed that one kind can be adopted
The nonstationary random process rapid simulation method of FFT.In addition, the method precision is very high and is easily achieved.Contrast for convenience, should
Method is defined as based on POD methods.In the present embodiment, the present invention will be contrasted and based on POD methods in efficiency and precision methods
Diversity.
In order to contrast, while simulating the earthquake motion stochastic process in embodiment 1 using both approaches.In simulations, two
Discrete frequency number in the method for kind is set to 1024, while it is 32 to 512 to simulate points excursion.The simulation of two methods
Matlab program calculations are all adopted, and operates in 64, Intel (R) Xeon (R) E5-2609 v2 processor (frequencies
2.50GHz) and on the computer of 32GB internal memories.The total calculating time contrast of two methods is as shown in table 1.Wherein, unit of time is
Second, time ratios are defined as the ratio of the calculating time based on POD methods and the inventive method calculating time.As can be seen that working as
When simulation points are larger, the simulation precision of the inventive method is apparently higher than the computational efficiency based on POD methods.For example, simulation is worked as
When points are more than 500, it is nearly 17000 times that the inventive method can improve simulation precision.Additionally, with the increase of simulation points, mould
The raising for intending efficiency becomes apparent from.
The contrast (s) of two methods simulation precision in the example 2 of table 1.
In terms of accuracy comparison, the simulation precision based on POD methods is mainly decoupled by POD and reconstructs what is determined.If POD
Decoupling and the precision for reconstructing are very high, then the final simulation precision of the method is equally very high.However, the precision of the present invention mainly depends on
In two aspects.First aspect is whether the two-dimentional EPSD reconstruct of conversion random wave is accurate enough.When POD is reconstructed, the item for using
When number is sufficiently large, the Evolutionary Spectra of reconstruct is almost exact value.On the other hand it is the two dimension that random wave is converted in discretization process
Whether the energy of EPSD can accurately represent.For this purpose, larger cut-off wave number, cut-off frequency, discrete wave-number and discrete should be chosen
Frequency number.Because Two-dimensional FFT efficiency is very high, these parameters affect little to computational efficiency.Additionally, frequency displacement can also be used for further
Improve simulation precision.
Claims (4)
1. a kind of nonstationary random process rapid simulation method, it is characterised in that comprise the following steps:
Step 1:The target power spectrum of given multiple spot nonstationary random process and coherent function
IfFor the multiple spot nonstationary random process of zero-mean, T represents transposition, its
Jth item elementIt is shown below:
In formula:ω is circular frequency;Aj(ω t) is gradual modulation function;Zj(ω) it is random for complex orthogonal increment that average is zero
Process;
The evolution spectral power matrix of the nonstationary random process is shown below:
Wherein matrix each element is expressed as:
In formula:ForCorresponding steady bilateral power spectrum function, γjk(ω) it isWithBetween
Plural coherent function, and with following property
Wherein * represents complex conjugation;Accordingly,WithCross-correlation functionIt is expressed as:
In formula:τ is time delay;
Step 2:Obtain the two-dimentional evolution spectrum density of conversion random wave
Along multiple spot nonstationary random process p of spatial axes x distribution0T () is considered along spatial axes x1,x2,…,xnDistribution
One-dimensional discrete non-stationary random wave [f0(x1,t),f0(x2,t),…,f0(xn,t)]T, according to formula (5), γ (ξ, ω) has following
Characteristic:
γ (ξ, ω)=γ*(-ξ,ω);γ (ξ, ω)=γ*(ξ,-ω) (7)
It is provided with following Fourier transform pairs to exist:
In formula:β (κ, ω) is the Fourier transformation of γ (ξ, ω), while being also non-negative and the function with regard to origin symmetry;
Formula (9) is substituted into into formula (6), then:
For conversion random wave f0(x, t), its two-dimensional autocorrelation function representation is:
Wherein:To convert the 2D EPSD of random wave;Show that the 2D EPSD for converting random wave are:
Step 3:Two-dimentional evolution power spectrum based on POD is decoupled
First by G0(x, ω, t) carry out discrete in space, frequency domain and time domain, and are reconstructed into following matrix:
In formula:xj=j Δs x corresponds to locus and j=1,2 ..., n is Spatial Cable argument;It is frequency vector
And n2=0,1 ..., N2- 1 is frequency indices number;It is time and nt=1,2 ..., NtIt is time index number;Matrix It is respectively 1 × N with the dimension of Qt, 1 × nNtAnd N2×nNt;
Orthogonal basiss are determined using Eigenvalues Decomposition method, is shown below:
RΦq=λqΦq, q=1,2 ..., N2 (16)
In formula:R is the covariance matrix of matrix Q and is expressed as form:
Wherein:λqAnd ΦqThe respectively eigen vector of q items, and satisfaction
It is by matrix Q approximate representations:
In formula:NqIt is to be superimposed the number of item and meet Nq< N2;hqIt is q item principal coordinates;The orthogonal property of feature based vector, hq
It is expressed as:
Formula (18) is expressed as again from the angle of element:
In formula:WithIt is ΦqN-th2Individual element and hqR-th element;Wherein r=(j-1) Nt+nt;
The conitnuous forms of formula (20) are expressed as:
In formula:aq(x, t) is the real-number function of room and time;Φq(ω) be frequency even function;
Formula (21) is substituted into into formula (12), random wave decomposition two dimension EPSD is converted and is expressed as:
In formula:
Step 4:The Two-dimensional FFT efficient simulation equidistantly put
By two-dimentional evolution power spectral densityOne-dimensional non-stationary random wave simulation be converted into it is a series of steady
Homogenizing random waveSimulation;Correspond toSimulation formula be shown below:
Therefore, sample f (x, t) of original non-stationary random wave is determined by following formula:
Formula (23) is expressed as again:
In formula:Re represents treating excess syndrome number;p1=0,1 ..., m × M1- 1 and p2=0,1 ..., m × M2-1;q1And q2Respectively p1/M1
And p2/M2Remainder;WithRespectively:
Wherein:WithIt is expressed as:
Δ x and Δ t should meet Δ x Δs κ=2 π/M1With Δ t Δs ω=2 π/M2, while meeting:
M1≥2N1;M2≥2N2。 (30) 。
2. a kind of nonstationary random process rapid simulation method as claimed in claim 1, it is characterised in that in step 2, draw
During the 2D EPSD of conversion random wave, ending wave number and frequency need to cover the main energetic of two-dimentional evolution power spectral density,
The discretization points of wave number and frequency can catch pulse spectrum.
3. a kind of nonstationary random process rapid simulation method as claimed in claim 2, it is characterised in that in step 2, when
When pulse spectrum is located at origin, simulation precision, usage frequency skew are improved using frequency shift (FS), its simulation formula is:
In formula:It is frequency shift parameters;Take positive integer;n1=0,1 ..., N1-1;n2=0,1 ..., N2-1;WithIt is separate and is uniformly distributed in the random of [0,2 π]
Phase angle;WhenWithAfter determining respectively, obtain converting the sample Equation f of random wave(i)(x,t)。
4. a kind of nonstationary random process rapid simulation method as claimed in claim 3, it is characterised in that adopt frequency shift (FS)
When improving simulation precision, space increments Δ x and incremental time Δ t need to meet Δ x≤2 π/(2 κu) and Δ t≤2 π/(2 ωu)。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201611106931.5A CN106682277A (en) | 2016-12-06 | 2016-12-06 | Fast simulation method for non-stationary random process |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201611106931.5A CN106682277A (en) | 2016-12-06 | 2016-12-06 | Fast simulation method for non-stationary random process |
Publications (1)
Publication Number | Publication Date |
---|---|
CN106682277A true CN106682277A (en) | 2017-05-17 |
Family
ID=58867305
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201611106931.5A Pending CN106682277A (en) | 2016-12-06 | 2016-12-06 | Fast simulation method for non-stationary random process |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN106682277A (en) |
Cited By (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107247687A (en) * | 2017-05-18 | 2017-10-13 | 西南交通大学 | A kind of nonstationary random process rapid simulation method of feature based Orthogonal Decomposition |
CN109490626A (en) * | 2018-12-03 | 2019-03-19 | 中车青岛四方机车车辆股份有限公司 | A kind of standard PSD acquisition methods and device based on Nonstationary Random Vibration Signals |
CN110826197A (en) * | 2019-10-21 | 2020-02-21 | 西南交通大学 | Wind speed field simulation method based on improved Cholesky decomposition closed solution |
CN111458744A (en) * | 2020-04-09 | 2020-07-28 | 西南交通大学 | Spatial rotation seismic oscillation simulation method |
CN113806845A (en) * | 2021-09-16 | 2021-12-17 | 四川农业大学 | POD interpolation-based non-stationary wind field efficient simulation method |
CN114139363A (en) * | 2021-11-25 | 2022-03-04 | 四川农业大学 | Method for simulating space variation non-stationary earthquake motion time course |
CN115526125A (en) * | 2022-09-22 | 2022-12-27 | 四川农业大学 | Numerical truncation-based high-efficiency random wind speed field simulation method |
CN116186838A (en) * | 2022-12-27 | 2023-05-30 | 武汉理工大学 | Structural random response analysis method, device and storage medium based on harmonic wavelet |
-
2016
- 2016-12-06 CN CN201611106931.5A patent/CN106682277A/en active Pending
Cited By (13)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107247687B (en) * | 2017-05-18 | 2021-06-08 | 西南交通大学 | Fast simulation method for non-stationary random process based on feature orthogonal decomposition |
CN107247687A (en) * | 2017-05-18 | 2017-10-13 | 西南交通大学 | A kind of nonstationary random process rapid simulation method of feature based Orthogonal Decomposition |
CN109490626A (en) * | 2018-12-03 | 2019-03-19 | 中车青岛四方机车车辆股份有限公司 | A kind of standard PSD acquisition methods and device based on Nonstationary Random Vibration Signals |
CN109490626B (en) * | 2018-12-03 | 2021-02-02 | 中车青岛四方机车车辆股份有限公司 | Standard PSD obtaining method and device based on non-stationary random vibration signal |
CN110826197A (en) * | 2019-10-21 | 2020-02-21 | 西南交通大学 | Wind speed field simulation method based on improved Cholesky decomposition closed solution |
CN111458744B (en) * | 2020-04-09 | 2021-07-16 | 西南交通大学 | Spatial rotation seismic oscillation simulation method |
CN111458744A (en) * | 2020-04-09 | 2020-07-28 | 西南交通大学 | Spatial rotation seismic oscillation simulation method |
CN113806845A (en) * | 2021-09-16 | 2021-12-17 | 四川农业大学 | POD interpolation-based non-stationary wind field efficient simulation method |
CN114139363A (en) * | 2021-11-25 | 2022-03-04 | 四川农业大学 | Method for simulating space variation non-stationary earthquake motion time course |
CN114139363B (en) * | 2021-11-25 | 2022-07-01 | 四川农业大学 | Method for simulating space variation non-stationary earthquake motion time course |
CN115526125A (en) * | 2022-09-22 | 2022-12-27 | 四川农业大学 | Numerical truncation-based high-efficiency random wind speed field simulation method |
CN115526125B (en) * | 2022-09-22 | 2023-09-19 | 四川农业大学 | Numerical cutoff-based efficient simulation method for random wind speed field |
CN116186838A (en) * | 2022-12-27 | 2023-05-30 | 武汉理工大学 | Structural random response analysis method, device and storage medium based on harmonic wavelet |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN106682277A (en) | Fast simulation method for non-stationary random process | |
Peng et al. | Simulation of multivariate nonstationary random processes: Hybrid stochastic wave and proper orthogonal decomposition approach | |
Liu et al. | Dimension reduction of Karhunen-Loeve expansion for simulation of stochastic processes | |
CN109883532B (en) | Sound source identification and sound field forecasting method | |
Castanié | Digital spectral analysis: parametric, non-parametric and advanced methods | |
CN105426594A (en) | Method for fast simulating wind fields on basis of stability and homogeneity of time-space field and condition interpolation | |
Singh | Breaking the limits: redefining the instantaneous frequency | |
Wang et al. | A random physical model of seismic ground motion field on local engineering site | |
CN113673076B (en) | Dynamic response solving method suitable for ocean floating structure | |
Javorskyj et al. | Analysis of mean function discrete LSM-estimator for biperiodically nonstationary random signals | |
CN102435989B (en) | Field programmable gate array (FPGA)-based general wave beam forming device | |
CN104516771A (en) | Efficient non-stationary random process simulation method | |
Kozakowski et al. | Structure factor for decorated Penrose tiling in physical space | |
Bocquet | Simulation of correlated Pareto distributed sea clutter | |
CN109870669A (en) | How soon a kind of two dimension claps mesh free compression Wave beam forming identification of sound source method | |
Douglas | Multigrid methods in science and engineering | |
Park et al. | Frequency domain spectral element model for the vibration analysis of a thin plate with arbitrary boundary conditions | |
CN108614234B (en) | Direction-of-arrival estimation method based on multi-sampling snapshot co-prime array received signal fast Fourier inverse transformation | |
Mohammadzaheri et al. | DISRAY: a distributed ray tracing by map-reduce | |
CN104077479A (en) | Parametric array sound field space distribution obtaining method based on conservation upwind scheme | |
CN112285770B (en) | Seismic wave phase spectrum disturbance method based on real number wavelet transformation | |
Gibson | Fourier expansion of disk automorphisms via scattering in layered media | |
Averbuch et al. | Direct inversion of the three-dimensional pseudo-polar Fourier transform | |
CN105572472B (en) | The frequency measurement method and system of distributed generation resource environment | |
CN114139363A (en) | Method for simulating space variation non-stationary earthquake motion time course |
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 | ||
RJ01 | Rejection of invention patent application after publication | ||
RJ01 | Rejection of invention patent application after publication |
Application publication date: 20170517 |