CN110441815A - Decline improved simulated annealing Rayleigh waves inversion method based on differential evolution and block coordinate - Google Patents
Decline improved simulated annealing Rayleigh waves inversion method based on differential evolution and block coordinate Download PDFInfo
- Publication number
- CN110441815A CN110441815A CN201910781681.2A CN201910781681A CN110441815A CN 110441815 A CN110441815 A CN 110441815A CN 201910781681 A CN201910781681 A CN 201910781681A CN 110441815 A CN110441815 A CN 110441815A
- Authority
- CN
- China
- Prior art keywords
- individual
- population
- inversion
- simulated annealing
- objective function
- 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.)
- Granted
Links
- 238000002922 simulated annealing Methods 0.000 title claims abstract description 50
- 238000000034 method Methods 0.000 title claims abstract description 43
- 230000007423 decrease Effects 0.000 title abstract description 6
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 32
- 238000012545 processing Methods 0.000 claims abstract description 6
- 238000005457 optimization Methods 0.000 abstract description 5
- 239000010410 layer Substances 0.000 description 21
- 239000006185 dispersion Substances 0.000 description 14
- 238000010586 diagram Methods 0.000 description 6
- 230000008901 benefit Effects 0.000 description 5
- 238000004364 calculation method Methods 0.000 description 4
- 238000004458 analytical method Methods 0.000 description 3
- 230000002068 genetic effect Effects 0.000 description 3
- 230000008569 process Effects 0.000 description 3
- 238000001514 detection method Methods 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000010438 heat treatment Methods 0.000 description 2
- 238000011835 investigation Methods 0.000 description 2
- 230000000644 propagated effect Effects 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 101000743114 Homo sapiens WASH complex subunit 4 Proteins 0.000 description 1
- 102100038143 WASH complex subunit 4 Human genes 0.000 description 1
- 238000000137 annealing Methods 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 230000000052 comparative effect Effects 0.000 description 1
- 230000004069 differentiation Effects 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 239000000284 extract Substances 0.000 description 1
- 238000009396 hybridization Methods 0.000 description 1
- 239000011229 interlayer Substances 0.000 description 1
- 230000005389 magnetism Effects 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 239000011159 matrix material Substances 0.000 description 1
- 210000003739 neck Anatomy 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000012827 research and development Methods 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/282—Application of seismic models, synthetic seismograms
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N3/00—Computing arrangements based on biological models
- G06N3/12—Computing arrangements based on biological models using genetic models
- G06N3/126—Evolutionary algorithms, e.g. genetic algorithms or genetic programming
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Remote Sensing (AREA)
- General Physics & Mathematics (AREA)
- Biophysics (AREA)
- Health & Medical Sciences (AREA)
- Evolutionary Biology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Geophysics (AREA)
- Geology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Environmental & Geological Engineering (AREA)
- Acoustics & Sound (AREA)
- Theoretical Computer Science (AREA)
- Biomedical Technology (AREA)
- Molecular Biology (AREA)
- Genetics & Genomics (AREA)
- Computational Linguistics (AREA)
- Data Mining & Analysis (AREA)
- Evolutionary Computation (AREA)
- General Health & Medical Sciences (AREA)
- Artificial Intelligence (AREA)
- Computing Systems (AREA)
- General Engineering & Computer Science (AREA)
- Mathematical Physics (AREA)
- Software Systems (AREA)
- Physiology (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
The invention discloses one kind to decline improved simulated annealing Rayleigh waves inversion method based on differential evolution and block coordinate, this method constructs inversion objective function first, differential evolution algorithm is recycled to optimize the parameter model of inversion objective function, it obtains carrying out simulated annealing processing to each individual using block coordinate descent algorithm after optimum individual, using obtained new individual as the individual of next generation population.The present invention passes through building inversion objective function, inversion objective function is optimized using differential evolution algorithm, and multiparameter problem is divided into multiple one-parameter local iterations optimization problem using block coordinate descending method, the precision of simulated annealing inverting can be effectively improved, realize simulated annealing to the controllability of precision, it to improve the stability of simulated annealing, can be applied in practical shallow layer exploration seismic data inverting, improve the precision of prediction of near surface underground structure.
Description
Technical field
The present invention relates near surface geophysics exploration techniques, are a kind of high-precision near surface Rayleigh waves inversion methods, should
Technology can be applied to near-surface velocity investigation, the necks such as geotechnical engineering, the analysis of oil-gas exploration near-surface velocity, underground engineering detection
Domain.
Background technique
Rayleigh waves be a kind of superposition mutual when being propagated by longitudinal wave and shear wave vertical component along free interface and formed, along from
The wave propagated by interface.Rayleigh waves usually along overland propagation, thickness and shear wave velocity for stratum it is very sensitive (Xia etc.,
1999), there is Dispersion (Ma Jianqing etc., 2010) in layered medium, compared to common exploration method, Rayleigh waves exploration tool
There is the advantages that high resolution, detection speed is fast (Qi Shengwen etc., 2002).Numerous scholars study Rayleigh waves, and are formed
Near surface geologic structure investigation method based on Rayleigh waves, if Rayleigh waves multi-channel analysis method is for estimating near surface shear wave speed
Degree (Park etc., 1998,1999;Xia etc., 1999).Currently, with the quickening of urbanization process, mankind's activity is near surface benefit
With gradually increasing, the Mapping The S-wave Velocity Structure inversion technique based on Rayleigh waves is as a kind of detected with high accuracy method, in Urban Underground
It is widely used (Pugin etc. and Krawczyk etc., 2013) in space exploration, has become a kind of important skill of city shallow layer exploration
Art means (summer Jiang Hai etc., 2015).Therefore research and development Rayleigh waves inversion method has important practical significance.
Current more common shear wave velocity nonlinear inversion method can be divided into two classes: local search approach
(LSMs) and global search method (Socco etc., 2010).LSMs is since calculation amount is smaller, using relatively broad in engineering.So
And such method is easily trapped into locally optimal solution.In recent years, it is contemplated that global search method can jump out locally optimal solution, disobey
Rely initial model and do not need the advantages that seeking directional derivative, more and more researchers solve frequency using global search method
Non-dramatic song line inversion problem (Lu etc., 2016) (Yu Dongkai, 2018).It is simulated annealing (SA) method and its changes using more at present
Into technology: such as Beaty carries out inverting (Beaty etc., 2002) to multi-mode Rayleigh waves dispersion curve using simulated annealing method,
Pei etc. has carried out inverting research (Pei etc., 2007) dispersion curve using simulated annealing, and Lu etc. carries out simulated annealing
It improves, proposes heating bath simulated annealing, and it is compared with L-M algorithm and Fast Simulated Annealing, show the party
Method is more suitable for inverting dispersion curve (Lu etc., 2016).But traditional analog annealing algorithm globally optimal still has and such as counts
Calculation amount is big, optimization direction is unstable, multi-parameter Simultaneous Retrieving thickness and the drawbacks such as error is larger when speed, makes it in dispersion curve
Certain restrictions is received in application in inverting.
Summary of the invention
The main purpose of the present invention is to provide a kind of precision that can be improved simulated annealing inverting in Rayleigh waves exploration
And the simulated annealing Rayleigh waves inversion method of stability, it is intended to solve the above technical problem present in existing method.
To achieve the above object, the present invention provides a kind of auspicious based on differential evolution and the improved simulated annealing of block coordinate decline
Leibo inversion method, comprising the following steps:
S1, data are calculated according to Rayleigh waves in Surface wave inversion model and observe data building inversion objective function;
S2, it is made a variation, crossover operation, is selected using population of the differential evolution algorithm to parameter in inversion objective function
The individual optimal to fitness;
S3, simulated annealing processing is carried out to each individual that step S2 is selected using block coordinate descent algorithm;
S4, the new individual for handling step S3 judge whether to meet termination condition as the individual of next-generation population;
If so, output obtains Rayleigh waves inverting optimal solution;If it is not, then return step S2.
Further, in the step S1, heat is replaced using the root mean square E for calculating the Rayleigh waves phase velocity obtained with observation
Energy in mechanics, obtains inversion objective function, is expressed as
Wherein, N indicates the stepped-frequency signal number of Rayleigh waves phase velocity, Vri obsIndicate the observation at ith sample point
Phase velocity, Vri calModel m after indicating disturbance is in the calculated phase velocity of ith sample point.
Further, the step S2 specifically include it is following step by step:
S21, the Population Size of parameter in inversion objective function is set as M, one of individual miShear wave in corresponding model
Speed is Vsi, thickness hi, individual dimension is N, and the first generation individual of population is randomly generated, is expressed as
I=1,2 ..., M, j=1,2 ..., N
Wherein,Indicate i-th of the individual speed in jth layer,Indicate thickness of i-th of individual in jth layer, (V1,
V2) and (H1,H2) value range of jth interval velocity and thickness is respectively indicated, rand () is random function;
S22, two different individual m are randomly selected in populationbAnd mc, treat the individual m of variationaMutation operation is carried out,
It is expressed as
m'i(g+1)=ma(g)+F(mb(g)-mc(g)), a ≠ b ≠ c, i=1,2 ... M
Wherein, m'i(g+1) indicate g+1 for the individual of i-th of variation in population, ma(g)、mb(g)、mc(g) g is indicated
For three Different Individuals in population, F indicates zoom factor;
S23, g+1 is become with g for i-th of individual after Population Variation for i-th of individual after Population Variation
ETTHER-OR operation is expressed as
Wherein,The corresponding shear wave velocity of individual after indicating crossover operation, Vsi(g+1) g+1 generation kind is indicated
The corresponding shear wave velocity of i-th of individual after group's variation, Vsi(g) indicate that g is corresponding for i-th of individual after Population Variation
Shear wave velocity, j indicate the jth dimension of shear wave velocity, riIt obeys (0,1) to be uniformly distributed, CR is to intersect the factor;
S24, selection operation is carried out to the individual after variation, crossover operation using greedy algorithm, be expressed as
Wherein, mi(g+1) i-th individual of the g+1 for population of selection, m are indicatedi(g) indicate the g of selection for population
I-th individual.
Further, the step S3 specifically include it is following step by step:
S31, setting initial temperature T0, final temperature Tend, temperature descent coefficient a, maximum number of iterations num, disturbance side
Formula, and parameter is divided into m block;
S32, a parameter block, and fixed other parameters block are randomly updated, judges whether to receive newly-generated parameter block, presses
This mode successively updates remaining parameter block, until reaching maximum number of iterations;
S33, temperature is reduced according to temperature descent coefficient, judges whether Current Temperatures are lower than final temperature;If so,
Obtain optimized parameter;If it is not, then return step S32.
Further, the step S32 specifically:
Current Temperatures are set as initial temperature, generate initial individuals at random, calculate inversion objective function E;
New individual is generated according to perturbation scheme, and calculates inversion objective function E ';
Calculate inversion objective function E ' and inversion objective function E difference, judge the difference whether less than zero;If so,
Receive new individual;If it is not, then receiving new individual with probability exp (- Δ E/T);
Judge whether to reach maximum number of iterations;If so, carrying out next step;If it is not, then being generated again according to disturbance
New individual.
The beneficial effects of the present invention are: the present invention passes through building inversion objective function, using differential evolution algorithm to inverting
Objective function optimizes, and multiparameter problem is divided into multiple one-parameter local iterations using block coordinate descending method and is optimized
Problem can effectively improve the precision of simulated annealing inverting, realize that simulated annealing to the controllability of precision, moves back to improve simulation
The stability of fire, can be applied in practical shallow layer exploration seismic data inverting, improve the precision of prediction of near surface underground structure.
Detailed description of the invention
Fig. 1 is of the invention to decline improved simulated annealing Rayleigh waves inversion method process based on differential evolution and block coordinate
Figure;
Fig. 2 is the contrast schematic diagram of SA and BCDESA to three model A inversion results;
Fig. 3 is the contrast schematic diagram of SA and BCDESA to three Model B inversion results;
Fig. 4 is the contrast schematic diagram of SA and BCDESA to three MODEL C inversion results;
Fig. 5 is Yellowstone seismic data and frequency dispersion energy diagram;
Fig. 6 is the inversion result schematic diagram of BCDESA and Geopsy software to Huangshi mountain park seismic data.
Specific embodiment
In order to make the objectives, technical solutions, and advantages of the present invention clearer, with reference to the accompanying drawings and embodiments, right
The present invention is further elaborated.It should be appreciated that described herein, specific examples are only used to explain the present invention, not
For limiting the present invention.
Rayleigh waves inverting is a kind of nonlinearity, multi-modal problem.Simulated annealing is main in Rayleigh waves inverting
Face two problems: first, Rayleigh waves inverting is a multiparameter problem, and as the inverting number of plies increases, inverted parameters will be at
Increase again, and simulated annealing is more suitable for solving one-parameter problem, can seriously affect final result when inverted parameters change is more.
Second, simulated annealing is using temperature as termination condition, and the direction of search is random, and when algorithm terminates, inverting may be because of parameter setting
It is unreasonable that result is caused not necessarily to reach global optimum.
The primary solutions of the embodiment of the present invention are as follows:
As shown in Figure 1, a kind of anti-based on differential evolution and improved simulated annealing (BCDESA) Rayleigh waves of block coordinate decline
Drill method, comprising the following steps:
S1, data are calculated according to Rayleigh waves in Surface wave inversion model and observe data building inversion objective function;
S2, it is made a variation, crossover operation, is selected using population of the differential evolution algorithm to parameter in inversion objective function
The individual optimal to fitness;
S3, simulated annealing processing is carried out to each individual that step S2 is selected using block coordinate descent algorithm;
S4, the new individual for handling step S3 judge whether to meet termination condition as the individual of next-generation population;
If so, output obtains Rayleigh waves inverting optimal solution;If it is not, then return step S2.
In an alternate embodiment of the present invention where, for above-mentioned steps S1 in Surface wave inversion, model m includes shear wave velocity
Vs, velocity of longitudinal wave Vp, density p and thickness h replace heating power using the root mean square E for calculating the Rayleigh waves phase velocity obtained with observation
Energy in, obtains inversion objective function, and the relationship being expressed as is
Wherein, N indicates the stepped-frequency signal number of Rayleigh waves phase velocity, Vri obsIndicate the observation at ith sample point
Phase velocity, Vri calModel m after indicating disturbance is in the calculated phase velocity of ith sample point.Shear wave velocity Vs and velocity of longitudinal wave
The relationship of Vp is
σ is Poisson's ratio.
In an alternate embodiment of the present invention where, in above-mentioned steps S2 differential evolution thought source in the proposition of early stage
Genetic algorithm simulates the hybridization in science of heredity, variation, duplication set computational algorithm.Differential evolution algorithm is relative to genetic algorithm
For, identical point be all by generating initial population at random, using the fitness of individual each in population as selection criteria, main mistake
Journey all includes variation, intersects and select three steps.The difference is that genetic algorithm is that control parent according to fitness miscellaneous
It hands over, the probability value that the filial generation generated after variation is selected, the probability that the small individual of fitness is selected in minimization problem is big,
And differential evolution algorithm variation vector be generated by parent difference vector, and with parent individuality vector intersects generation new individual to
Amount, is directly selected with its parent individuality, avoids falling into locally optimal solution to a certain extent.The size of population mainly reflects calculation
The size of species information amount in method, population is bigger, and information content is abundanter, and algorithm is more stable, but calculation amount is bigger.Differential evolution side
Method can be used to improve simulated annealing, keep simulated annealing controllable to error, improve the stability of algorithm.
The present invention makes a variation to the population of parameter in inversion objective function using differential evolution algorithm, crossover operation, choosing
Select to obtain the optimal individual of fitness, specifically include it is following step by step:
S21, the Population Size of parameter in inversion objective function is set as M, one of individual miFor a N-dimensional matrix,
Shear wave velocity is in corresponding modelThickness isIndividual dimension is N, with
Machine generates the first generation individual of population, is expressed as
I=1,2 ..., M, j=1,2 ..., N
Wherein, whereinIndicate i-th of the individual speed in jth layer,Indicate i-th of individual thickness in jth layer,
(V1,V2) and (H1,H2) value range of jth interval velocity and thickness is respectively indicated, rand () is random function;rand(V1,V2)
It indicates at random in V1To V2A random number, rand (H are generated in range1,H2) indicate at random in H1To H2One is generated in range
Random number, velocity of longitudinal wave Vp are updated by above formula, and density takes fixed value constant;The Population Size finally initialized is M, individual dimension
For N.
S22, two different individual m are randomly selected in populationbAnd mc, treat the individual m of variationaMutation operation is carried out,
It is expressed as
m'i(g+1)=ma(g)+F(mb(g)-mc(g)), a ≠ b ≠ c, i=1,2 ... M
Wherein, m'i(g+1) indicate g+1 for the individual of i-th of variation in population, ma(g)、mb(g)、mc(g) g is indicated
For three Different Individuals in population, F ∈ [0,2] indicates zoom factor, usually takes 0.5;Corresponding velocity of longitudinal wave by above formula more
Newly, density remains unchanged.
S23, g+1 is become with g for i-th of individual after Population Variation for i-th of individual after Population Variation
ETTHER-OR operation is expressed as by taking shear wave velocity Vs as an example
Wherein,The corresponding shear wave velocity of individual after indicating crossover operation, Vsi(g+1) g+1 generation kind is indicated
The corresponding shear wave velocity of i-th of individual after group's variation, Vsi(g) indicate that g is corresponding for i-th of individual after Population Variation
Shear wave velocity, j indicate the jth dimension of shear wave velocity, ri(0,1) is obeyed to be uniformly distributed, CR ∈ [0.1,0.6] is to intersect the factor, one
As take 0.3;
S24, selection operation is carried out to the individual after variation, crossover operation using greedy algorithm, be expressed as
Wherein, mi(g+1) i-th individual of the g+1 for population of selection, m are indicatedi(g) indicate the g of selection for population
I-th individual.
In an alternate embodiment of the present invention where, block coordinate descending method is a kind of non-gradient optimization in above-mentioned steps S3
Algorithm.It is continuously minimized along coordinate direction to find functional minimum value.In each iteration, algorithm passes through coordinate selection
Rule determines coordinate or coordinate block, minimizes corresponding objective function while then fixing every other coordinate or coordinate block.
Rayleigh waves inverting is a multiparameter problem, and simulated annealing can the less problem of preferable processing parameter.
Although density and velocity of longitudinal wave regards as it is known that the number of parameter is, wherein n is the number of plies in Rayleigh waves inverting.Therefore usually
Way be fixed thickness, disturb shear wave velocity, judge whether to receive the model, then fix shear wave velocity, disturb thickness, In
Judge whether to receive new thickness model, the two alternating iteration finally reaches balance or maximum number of iterations in Current Temperatures
(Pei etc., 2007).And in order to acquire optimal solution as far as possible, usually temperature is set and maximum number of iterations is larger.However work as
When the number of plies increases, it is deteriorated since the characteristics of simulated annealing will lead to efficiency of inverse process.In order to give full play to the advantage of simulated annealing, draw
Enter block coordinate descending method, multiparameter problem is divided into multiple one-parameter local iterations optimization problem, accelerating algorithm convergence.
Since differential evolution algorithm convergence rate is very fast in step S2, but with the increase of evolutionary generation, between each individual
Differentiation information be gradually reduced, convergence rate is slack-off, solve multi-solution Rayleigh waves inversion problem in, be easily trapped into part
Optimal solution;Therefore simulated annealing is added in step S3 in differential evolution, initializes one group of random population first, then to variation, friendship
Each individual after fork, selection operation carries out primary brief block coordinate and declines simulated annealing, under obtained new individual is used as
The individual of generation population, several times after iteration or after reaching termination condition, the solution for selecting fitness optimal is as output.
The present invention carries out simulated annealing processing to each individual that step S2 is selected using block coordinate descent algorithm, specifically
Including it is following step by step:
S31, setting initial temperature T0, final temperature Tend, temperature descent coefficient a, maximum number of iterations num, disturbance side
Formula, and parameter is divided into m block;
S32, a parameter block, and fixed other parameters block are randomly updated, judges whether to receive newly-generated parameter block, presses
This mode successively updates remaining parameter block, until reaching maximum number of iterations, specifically:
Current Temperatures are set as initial temperature T=T0, generate initial individuals at random, calculate inversion objective function E;
New individual is generated according to perturbation scheme, and calculates inversion objective function E ';
Calculate inversion objective function E ' and inversion objective function E difference DELTA E=E '-E, judge the difference whether less than zero
Δ E < 0;If so, receiving new individual;If it is not, then receiving new individual with probability exp (- Δ E/T);
Judge whether to reach maximum number of iterations;If so, carrying out next step;If it is not, then being generated again according to disturbance
New individual.
S33, temperature T=α T0 is reduced according to temperature descent coefficient, judges whether Current Temperatures are lower than final temperature T
≤Tend;If so, obtaining optimized parameter;If it is not, then return step S32.
The Population Variation thought that one aspect of the present invention can use differential evolution algorithm improves simulated annealing, will be in population
Termination condition of the fitness of optimum individual as algorithm realizes on the other hand control of the simulated annealing to error is also kept away
Simple differential evolution method is exempted from since individual difference information reduction causes to restrain slack-off problem, to improve global optimization method
Accuracy.
The present invention is tested using three kinds of classical theoretical models, and respectively 1 representation model A of table is incremented by for shear wave velocity
Concordant model, Model B shown in table 2 be one layer of shear wave velocity be higher than lower layer higher velocity layer model and table 3 shown in MODEL C be
One layer of shear wave velocity is lower than the low speed sandwich mould on upper layer.
Table 1, model A, three interval velocity increment type geological model parameters
Table 2, Model B, three layers of soft interlayer geological model parameter containing low speed
Table 3, MODEL C, three layers of geological model parameter of hard interbedded layer containing high speed
Using BCDESA and SA to without make an uproar theoretical model carry out inverting experiment, the inversion result of three models respectively as Fig. 2,
3, shown in 4.In figure left figure indicate respectively indicate model A, B, C theoretical dispersion curve and inverting after gained model dispersion curve
Comparative situation, right figure respectively indicate the velocity profile comparison of shear-wave velocity section figure and theoretical model that different model inversions obtain
Situation.From inversion result, it can be concluded that, in the insufficient situation of prior information, compared to SA, frequency dispersion obtained by BCDESA inverting is bent
Line with theoretical observation be more it is identical, SA inversion result does not even identify Model B.BCDESA is to model A, B, C inverting
The worst error of gained mould parameter and theoretical model parameter shear wave velocity is respectively 2.55%, 3.64% and 5.14%, SA difference
It is 3.15%, 24.38% and 13.64%,;The worst error of thickness is respectively 6.45%, 5.60% and 12.40%, SA difference
It is 8.80%, 22.47% and 33.67%.Wherein, half space velocity deviation is larger, main reason is that dispersion curve low frequency is believed
Breath is insufficient, but generally speaking, obtains underground shear wave velocity section using BCDESA and is not much different with model data, illustrates this hair
It is bright inverting to be carried out to various Rayleigh waves base rank dispersion curves well.
For the applicability for further examining BCDESA, the present invention is acquired using BCDESA in U.S. Huangshi mountain national park
Data (Pasquet etc., 2017) carry out back analysis.As shown in figure 5, left figure is one of collected earthquake record, use
10 24 channel 4.5Hz geometry earth magnetism seismic detectors, road spacing are 1m, and most trail spacing is 1m, and focus hits shake using hammer (5.4kg)
Source.Sample frequency is 0.125ms, when sampling a length of 0.75s.Right figure is the frequency dispersion energy diagram that big gun collection record extracts, from base rank frequency
Dissipate upper pickup measured data (white black circle in figure).The Rayleigh waves base rank dispersion curve picked up using BCDESA inverting, and
With the inversion result ratio of the open source software-Geopsy (www.geopsy.org) of the geophysical research and application that are called in SWIP
Compared with.Wherein Geopsy assumes that stratum is 10 layers, is 0-2.5m per thickness layer by layer, and assume that it is incremented by model for shear wave velocity.
Stratum is also assumed to be 10 layers by BCDESA, is 0-2.5m per thickness layer by layer, but the prior information of not given shear wave velocity, because above
Gross data test show that BCDESA can be with the various geological models of automatic identification.The termination error of BCDESA are set as by the present invention
2.5, other parameters with it is consistent above.Due to the search range given thickness, BCDESA restrains quickly, and in first generation population
In meet termination condition.Fig. 6 is the inversion result of BCDESA and Geopsy, and left figure is true model, BCDESA and Geopsy
Dispersion curve comparison, the error of BCDESA, Geopsy and true model is 2.5.Show the BCDESA method have compared with
Strong accuracy.Right figure is the shear-wave velocity section situation of inverting.Although the forward modeling method of the two is different, BCDESA's
Shear-wave velocity section is similar to Geopsy, and BCDESA has identified high-speed layer, therefore inversion result of the invention is more
Really.
Those of ordinary skill in the art will understand that the embodiments described herein, which is to help reader, understands this hair
Bright principle, it should be understood that protection scope of the present invention is not limited to such specific embodiments and embodiments.This field
Those of ordinary skill disclosed the technical disclosures can make according to the present invention and various not depart from the other each of essence of the invention
The specific variations and combinations of kind, these variations and combinations are still within the scope of the present invention.
Claims (5)
1. one kind declines improved simulated annealing Rayleigh waves inversion method based on differential evolution and block coordinate, which is characterized in that packet
Include following steps:
S1, data are calculated according to Rayleigh waves in Surface wave inversion model and observe data building inversion objective function;
S2, it is made a variation using population of the differential evolution algorithm to parameter in inversion objective function, crossover operation, selection is fitted
The optimal individual of response;
S3, simulated annealing processing is carried out to each individual that step S2 is selected using block coordinate descent algorithm;
S4, the new individual for handling step S3 judge whether to meet termination condition as the individual of next-generation population;If
It is that then output obtains Rayleigh waves inverting optimal solution;If it is not, then return step S2.
2. improved simulated annealing Rayleigh waves inversion method is declined based on differential evolution and block coordinate as described in claim 1,
It is characterized in that, being replaced in thermodynamics in the step S1 using the root mean square E for calculating the Rayleigh waves phase velocity obtained with observation
Energy, obtain inversion objective function, be expressed as
Wherein, N indicates the stepped-frequency signal number of Rayleigh waves phase velocity, Vri obsIndicate the observation phase velocity at ith sample point
Degree, Vri calModel m after indicating disturbance is in the calculated phase velocity of ith sample point.
3. improved simulated annealing Rayleigh waves inversion method is declined based on differential evolution and block coordinate as claimed in claim 2,
It is characterized in that, the step S2 specifically include it is following step by step:
S21, the Population Size of parameter in inversion objective function is set as M, one of individual miShear wave velocity is in corresponding model
Vsi, thickness hi, individual dimension is N, and the first generation individual of population is randomly generated, is expressed as
Wherein,Indicate i-th of the individual speed in jth layer,Indicate thickness of i-th of individual in jth layer, (V1,V2) and
(H1,H2) value range of jth interval velocity and thickness is respectively indicated, rand () is random function;
S22, two different individual m are randomly selected in populationbAnd mc, treat the individual m of variationaMutation operation is carried out, is indicated
For
m′i(g+1)=ma(g)+F(mb(g)-mc(g)), a ≠ b ≠ c, i=1,2 ... M
Wherein, m 'i(g+1) indicate g+1 for the individual of i-th of variation in population, ma(g)、mb(g)、mc(g) g generation kind is indicated
Three Different Individuals in group, F indicate zoom factor;
S23, g+1 is subjected to variation behaviour for i-th of individual after Population Variation with g for i-th of individual after Population Variation
Make, is expressed as
Wherein,The corresponding shear wave velocity of individual after indicating crossover operation, Vsi(g+1) indicate that g+1 becomes for population
The corresponding shear wave velocity of i-th of individual after different, Vsi(g) indicate g for the corresponding shear wave of i-th of individual after Population Variation
Speed, j indicate the jth dimension of shear wave velocity, riIt obeys (0,1) to be uniformly distributed, CR is to intersect the factor;
S24, selection operation is carried out to the individual after variation, crossover operation using greedy algorithm, be expressed as
Wherein, mi(g+1) i-th individual of the g+1 for population of selection, m are indicatedi(g) indicate the g of selection for the i-th of population
Individual.
4. improved simulated annealing Rayleigh waves inversion method is declined based on differential evolution and block coordinate as claimed in claim 3,
It is characterized in that, the step S3 specifically include it is following step by step:
S31, setting initial temperature T0, final temperature Tend, temperature descent coefficient a, maximum number of iterations num, perturbation scheme, and
Parameter is divided into m block;
S32, a parameter block, and fixed other parameters block are randomly updated, judges whether to receive newly-generated parameter block, by this side
Formula successively updates remaining parameter block, until reaching maximum number of iterations;
S33, temperature is reduced according to temperature descent coefficient, judges whether Current Temperatures are lower than final temperature;If so, obtaining
Optimized parameter;If it is not, then return step S32.
5. improved simulated annealing Rayleigh waves inversion method is declined based on differential evolution and block coordinate as claimed in claim 4,
It is characterized in that, the step S32 specifically:
Current Temperatures are set as initial temperature, generate initial individuals at random, calculate inversion objective function E;
New individual is generated according to perturbation scheme, and calculates inversion objective function E ';
Calculate inversion objective function E ' and inversion objective function E difference, judge the difference whether less than zero;If so, receiving
New individual;If it is not, then receiving new individual with probability exp (- Δ E/T);
Judge whether to reach maximum number of iterations;If so, carrying out next step;If it is not, then generating new according to disturbance again
Body.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910781681.2A CN110441815B (en) | 2019-08-23 | 2019-08-23 | Simulated annealing Rayleigh wave inversion method based on differential evolution and block coordinate descent |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910781681.2A CN110441815B (en) | 2019-08-23 | 2019-08-23 | Simulated annealing Rayleigh wave inversion method based on differential evolution and block coordinate descent |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110441815A true CN110441815A (en) | 2019-11-12 |
CN110441815B CN110441815B (en) | 2021-05-25 |
Family
ID=68437283
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910781681.2A Active CN110441815B (en) | 2019-08-23 | 2019-08-23 | Simulated annealing Rayleigh wave inversion method based on differential evolution and block coordinate descent |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110441815B (en) |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112182481A (en) * | 2020-10-10 | 2021-01-05 | 西安交通大学 | Seismic waveform inversion method and system based on improved differential evolution algorithm |
CN114185093A (en) * | 2021-12-07 | 2022-03-15 | 中国石油大学(北京) | Near-surface velocity model building method and device based on Rayleigh surface wave inversion |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2016140645A1 (en) * | 2015-03-02 | 2016-09-09 | Landmark Graphics Corporation | Selecting potential well locations in a reservoir grid model |
CN107909139A (en) * | 2017-11-13 | 2018-04-13 | 东南大学 | A kind of adjustable differential evolution algorithm of the crossover probability factor |
-
2019
- 2019-08-23 CN CN201910781681.2A patent/CN110441815B/en active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2016140645A1 (en) * | 2015-03-02 | 2016-09-09 | Landmark Graphics Corporation | Selecting potential well locations in a reservoir grid model |
CN107909139A (en) * | 2017-11-13 | 2018-04-13 | 东南大学 | A kind of adjustable differential evolution algorithm of the crossover probability factor |
Non-Patent Citations (3)
Title |
---|
K. S. BEATY 等: "Simulated annealing inversion of multimode Rayleigh wave dispersion curves for geological structure", 《GEOPHYSICAL JOURNAL INTERNATIONAL》 * |
YAOJUN WANG 等: "Rayleigh wave inversion based on differential evolution simulated annealing method", 《SEG TECHNICAL PROGRAM EXPANDED ABSTRACTS 2019》 * |
王峣钧 等: "基于自适应混合范数的快速多道波阻抗反演", 《CPS/SEG北京2018国际地球物理会议暨展览电子论文集》 * |
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112182481A (en) * | 2020-10-10 | 2021-01-05 | 西安交通大学 | Seismic waveform inversion method and system based on improved differential evolution algorithm |
CN112182481B (en) * | 2020-10-10 | 2022-12-09 | 西安交通大学 | Seismic waveform inversion method and system based on improved differential evolution algorithm |
CN114185093A (en) * | 2021-12-07 | 2022-03-15 | 中国石油大学(北京) | Near-surface velocity model building method and device based on Rayleigh surface wave inversion |
CN114185093B (en) * | 2021-12-07 | 2023-05-12 | 中国石油大学(北京) | Near-surface velocity model building method and device based on Rayleigh surface wave inversion |
Also Published As
Publication number | Publication date |
---|---|
CN110441815B (en) | 2021-05-25 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN109611087B (en) | Volcanic oil reservoir parameter intelligent prediction method and system | |
Xue et al. | Development of the inversion method for transient electromagnetic data | |
CN109325616A (en) | It is a kind of to be returned based on Gaussian process and the fine particle that combines of glowworm swarm algorithm is predicted and source tracing method | |
CN114723095A (en) | Missing well logging curve prediction method and device | |
CN110441815A (en) | Decline improved simulated annealing Rayleigh waves inversion method based on differential evolution and block coordinate | |
CN113568055A (en) | Aviation transient electromagnetic data retrieval method based on LSTM network | |
CN113486591B (en) | Gravity multi-parameter data density weighted inversion method for convolutional neural network result | |
CN113204054B (en) | Self-adaptive wide-area electromagnetic method induced polarization information extraction method based on reinforcement learning | |
CN113903395A (en) | BP neural network copy number variation detection method and system for improving particle swarm optimization | |
CN112035941A (en) | Prediction method for surface subsidence of deep foundation pit excavation based on BAS-BP model | |
CN114723784B (en) | Pedestrian motion trail prediction method based on domain adaptation technology | |
CN112614021A (en) | Tunnel surrounding rock geological information prediction method based on built tunnel information intelligent identification | |
CN116665067A (en) | Ore finding target area optimization system and method based on graph neural network | |
CN110988997A (en) | Hydrocarbon source rock three-dimensional space distribution quantitative prediction technology based on machine learning | |
CN109086568B (en) | Computer antibody combined mutation evolution system and method and information data processing terminal | |
CN112926251B (en) | Landslide displacement high-precision prediction method based on machine learning | |
CN113033074A (en) | Method, system and equipment for predicting porosity of policy combination mechanism fused dragonfly algorithm | |
CN116401965B (en) | Rock-soil strength parameter constraint random field simulation method based on measurement while drilling data | |
Abo-Zahhad et al. | Integrated model of DNA sequence numerical representation and artificial neural network for human donor and acceptor sites prediction | |
CN107067036A (en) | A kind of ground net corrosion rate prediction method | |
CN117371330B (en) | Magnetotelluric two-dimensional deep learning inversion method based on traditional inversion guidance | |
CN110223730A (en) | Protein and small molecule binding site prediction technique, prediction meanss | |
CN115267927B (en) | Multi-boundary curtain type geosteering method based on ant colony-gradient series algorithm | |
Wei et al. | A hybrid dung beetle optimization algorithm with simulated annealing for the numerical modeling of asymmetric wave equations | |
CN113435628B (en) | Medium-long-term runoff prediction method and system based on linear discriminant analysis and IALO-ELM |
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 | ||
GR01 | Patent grant | ||
GR01 | Patent grant |