CN110058307B - Full waveform inversion method based on fast quasi-Newton method - Google Patents
Full waveform inversion method based on fast quasi-Newton method Download PDFInfo
- Publication number
- CN110058307B CN110058307B CN201910367245.0A CN201910367245A CN110058307B CN 110058307 B CN110058307 B CN 110058307B CN 201910367245 A CN201910367245 A CN 201910367245A CN 110058307 B CN110058307 B CN 110058307B
- Authority
- CN
- China
- Prior art keywords
- inversion
- iteration
- data
- parameter model
- model
- 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.)
- Expired - Fee Related
Links
- 238000000034 method Methods 0.000 title claims abstract description 74
- 238000004364 calculation method Methods 0.000 claims abstract description 20
- 230000006870 function Effects 0.000 claims description 24
- 238000004088 simulation Methods 0.000 claims description 7
- 238000004422 calculation algorithm Methods 0.000 claims description 4
- 238000013501 data transformation Methods 0.000 claims description 2
- 230000001629 suppression Effects 0.000 claims description 2
- 238000005070 sampling Methods 0.000 description 3
- 230000005540 biological transmission Effects 0.000 description 2
- 239000011159 matrix material Substances 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000009933 burial Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000002708 enhancing effect Effects 0.000 description 1
- 238000001914 filtration Methods 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
Images
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
- G01V1/282—Application of seismic models, synthetic seismograms
-
- 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/30—Analysis
- G01V1/303—Analysis for determining velocity profiles or travel times
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/62—Physical property of subsurface
- G01V2210/622—Velocity, density or impedance
- G01V2210/6222—Velocity; travel time
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/66—Subsurface modeling
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
Abstract
The invention discloses a full waveform inversion method based on a fast quasi-Newton method, belongs to the technical field of seismic exploration speed modeling, and aims to improve the convergence speed of the quasi-Newton method, further improve the calculation efficiency of seismic full waveform inversion and reduce the calculation amount. The method comprises the following steps: defining a group of weighting coefficient sequences for constructing a parameter model iterative formula; introducing an intermediate iteration parameter model, and calculating a quasi-Newton iteration direction and an optimal step length based on the parameter model; and updating the original parameter model and the intermediate iteration parameter model by using the weighting coefficient, the iteration direction and the optimal step length. The invention obtains a parameter model closer to the optimal solution by calculation, and the parameter model is used for calculating the next iteration direction and the optimal step length, thereby improving the convergence speed of the target function. In each iteration calculation, the method does not increase excessive calculation amount, only has a small amount of multiplication and addition operation with low complexity, and effectively reduces the calculation amount of full waveform inversion.
Description
Technical Field
The invention belongs to the technical field of seismic exploration speed modeling, and relates to a full waveform inversion method based on a fast quasi-Newton method.
Background
The full waveform inversion technology theoretically has the capability of revealing geological structure and reservoir parameters under complex geological conditions, and the full waveform inversion of a frequency domain can realize high-precision parameter modeling by only using data of a few discrete frequencies. For solving the full waveform inversion problem, a nonlinear local optimization algorithm, namely a gradient method and a Newton method, is mainly adopted. Even if the gradient method is simple in structure and easy to implement, only gradient information is utilized in the calculation process, and the method has first-order convergence. The Newton method considers the gradient information of the target function and the Hessian matrix at the same time, and has second-order convergence, high convergence speed and high calculation precision. Although the quasi-Newton method avoids direct calculation of the Hessian matrix and effectively reduces time and memory consumption, the full waveform inversion technology is always limited to the application under the two-dimensional condition due to the limitation of calculation resources. Therefore, the research on the method for improving the convergence speed of the full waveform inversion and reducing the calculation amount and the memory requirement has important theoretical value and practical significance for enhancing the practicability of the full waveform inversion.
Disclosure of Invention
The invention provides a high-efficiency full waveform inversion method of a quasi-Newton method, aiming at the problem of huge calculated amount and storage amount of full waveform inversion of seismic reflection waves. A group of weighting coefficient sequences and intermediate iteration parameter models are introduced into a conventional quasi-Newton method, namely a parameter model which is closer to an optimal solution than the last time is calculated during each iteration, so that the number of elements in the iteration parameter model sequences is as small as possible, the target function value is reduced as much as possible, and the convergence speed of the algorithm is further accelerated. The technical scheme adopted by the invention is as follows:
step 1: data is input and processed. The input data comprises actual field seismic observation data, observation system parameters, an initial parameter model and the like. The data processing mainly comprises the steps of carrying out noise suppression, seismic channel interpolation, amplitude equalization, data transformation and the like on the seismic observation data.
Step 2: and solving the data residual error and the inversion objective function. Based on the parameters and the parameter model of the observation system, the wave equation seismic numerical simulation method is adopted to calculate and store the seismic wave field in the parameter model, namely the forward wave field, and calculation data are obtained. Subtracting the calculated data from the observed data to obtain a data residual, and then using a two-norm of the data residual as an inversion target function, namely
Wherein J denotes an inversion target function, Δ d denotes a data residual vector, and "+" denotes a conjugate transpose operation. The selection of the wave equation type needs to be determined according to actual inversion parameters, if only the inversion speed parameters adopt scalar sound wave equations, and the inversion density and the speed parameters adopt first-order stress-speed equations.
And step 3: the gradient of the inverted objective function is calculated. Based on the principle of an adjoint state method, conjugate back-propagation simulation is carried out on the data residual error in a parameter model by taking the data residual error as a seismic source to obtain a back-propagation wave field, and zero-delay cross-correlation operation is carried out on the forward-propagation wave field and the back-propagation wave field to obtain the gradient of a target function.
And 4, step 4: the iteration direction of the model update is calculated. Based on the algorithm principle of a finite-memory quasi-Newton method (L-BFGS), the gradient residual error and model residual error information in the last M (usually taken as 3-20) iteration steps are utilized to obtain the iteration direction of the model updating.
And 5: and calculating the optimal step length in the iteration direction by adopting a linear search method. The linear search method includes a non-precise line search method and a precise line search method. Although the non-precise line search method is a more common optimal step size search method, the full-waveform inversion mainly adopts a two-point parabola precise line search method, so that each iterative operation only needs to be performed once forward calculation, and the calculation amount is the minimum compared with other linear search methods.
Step 6: defining a set of weighting coefficient sequences, i.e.
Where i is the number of iterations, ηiAnd λiRespectively representing the intermediate variables and weighting coefficients of the ith iteration；
And 7: introducing an intermediate iteration parameter model yiAnd order y0＝v0(ii) a Updating the parametric model according to the following iterative formula
Where i is the number of iterations, viInverse parametric model representing the ith iteration, αiDenotes the step size of the ith iteration and d denotes the direction of the ith iteration.
In summary, due to the adoption of the technical scheme, the invention has the beneficial effects that:
the invention uses the parameter model updated by two times of iteration to carry out a special weighting treatment, obtains the parameter model closer to the optimal solution by calculation, and is used for calculating the next iteration direction and the optimal step length, thereby improving the convergence rate of the objective function. Compared with the conventional quasi-Newton method, the method of the invention does not increase excessive calculation amount in each iteration calculation, and only has a small amount of multiplication and addition operations with lower complexity. Under the condition of the same precision, the method effectively reduces the calculated amount of full waveform inversion.
Drawings
In order to more clearly illustrate the present invention, the drawings used in the embodiments will be briefly described below. It is to be expressly understood, however, that the drawings are for the purpose of providing a further understanding of the invention, and are incorporated in and constitute a part of this specification.
FIG. 1 is a flow chart of frequency domain multi-scale acoustic full waveform inversion according to the present invention;
FIG. 2 is a partial Marmousi velocity model;
FIG. 3 is an inverted initial model;
FIG. 4 is a normalized objective function descent curve obtained by single frequency (2.5Hz) inversion using the method of the present invention and a conventional quasi-Newton method (L-BFGS).
Fig. 5a is the result of a multi-scale full-waveform inversion (30 single-frequency iterations) using a conventional quasi-newton method (L-BFGS).
Fig. 5b is the result of a multi-scale full-waveform inversion (20 single frequency iterations) using a conventional quasi-newton method (L-BFGS).
FIG. 5c shows the result of performing multi-scale full-waveform inversion (single frequency iteration 20 times) using the method of the present invention.
Fig. 6 is a comparison graph of the velocity profile at lateral position X1540 m.
Detailed Description
In order to make the objects, steps and features of the present invention clearer, embodiments and applications of the present invention will be described in detail with reference to the accompanying drawings. The description information provided in the specific embodiments and applications is only an example for explaining the present invention and is not a limitation of the present invention. Various extensions and combinations of the embodiments described herein will be apparent to those skilled in the art, and the general principles defined herein may be applied to other embodiments and applications without departing from the spirit and scope of the invention.
FIG. 1 is a flow chart of frequency domain multi-scale acoustic wave full waveform inversion according to the present invention. As shown in fig. 1, a full waveform inversion method based on a fast quasi-newton method according to an embodiment of the present invention includes the following steps:
1. and inputting and processing data, wherein the input data comprises actual field earthquake observation data, observation system parameters, an initial velocity parameter model and the like. The forward record based on the real model is taken as the field seismic observation data in the embodiment.
2. Obtaining a discrete inversion frequency sequence f ═ f { f by adopting a multi-scale inversion frequency selection strategy according to input datalow,…,fhigh}. And successively inverting from low frequency to high frequency, and taking the result of low frequency inversion as an initial model of high frequency inversion.
3. And solving an objective function of the data residual error and the full waveform inversion. Introducing an intermediate iteration parameter model yiAnd order y0＝v0(ii) a The method comprises the steps of numerically simulating a sound wave equation based on a frequency domain, calculating and storing a seismic wave field u (forward wave field) in a parameter model, and obtaining the seismic wave field uCalculating data dcalAnd observation data dobsMaking difference to obtain data residual error delta d ═ dcal-dobs. The goal of full waveform inversion is to minimize the difference between the calculated and observed data, i.e., the objective function in the least squares sense can be expressed as
Where e (y) represents the inverse objective function and "+" represents the conjugate transpose operation.
4. The gradient of the objective function is calculated. Respectively carrying out conjugate back transmission simulation on the data residual errors corresponding to each single shot point based on the adjoint state method principle to obtain corresponding back transmission wave fields
Multiplying with forward wave field of corresponding shot point, and stacking according to shot set to obtain gradient at single frequency
Wherein g (y) represents the gradient of the objective function relative to the speed parameter model, omega is angular frequency, y is the parameter model, and s represents a shot set;
and R (delta d) is a source term for conjugate back propagation simulation, namely, the data residual error of the detector position is expanded to the whole model space to serve as the conjugate back propagation simulation.
5. Solving iteration direction d (y) of parameter model updating based on finite memory quasi-Newton method (L-BFGS)i) And a two-point parabola precise line search method is adopted to calculate the optimal step length alphai。
6. Calculating a weighting coefficient sequence lambdaiAnd updating the speed parameter model.
The weighting coefficient sequence is calculated by
Where i is the number of iterations, ηiAn intermediate variable representing the ith iteration;
the updating formula of the parameter model is
Where i is the number of iterations and yiAnd viRespectively representing the intermediate velocity parametric model and the inversion velocity parametric model of the ith iteration, alphaiDenotes the step size of the ith iteration and d denotes the direction of the ith iteration.
7. Judging whether a target function convergence condition is met or the maximum iteration number is reached, and if not, circularly executing the steps 3-6; if yes, terminating iteration, judging whether the current frequency is the maximum inversion frequency, if not, updating the inversion frequency, taking the result of the current low-frequency inversion as an initial model for updating the high-frequency inversion, circularly executing the steps 3-6, and if yes, terminating iteration and outputting the result of the multi-scale full-waveform inversion.
In order to more clearly demonstrate the effect of the full waveform inversion method, the present invention is described below with reference to specific embodiments.
And (3) performing a frequency domain multi-scale full-waveform inversion test on a part of the Marmousi velocity model (shown in figure 2) by adopting the full-waveform inversion method. The size of the model is 4300m × 1170m, the vertical and horizontal sampling grid is 431 × 235, and the sampling interval is 10m × 5 m. The observation system is provided with 54 seismic sources in total, the seismic sources are arranged in the range of 30-4270 m of the ground surface, the distance between the shot points is 80m, and the burial depth of the shot points is 0 m; 901 receivers are received by each shot, wherein 431 receivers are arranged on the ground surface at a track pitch of 10m (the buried depth is the same as that of the seismic source), and other 470 receivers are arranged in wells with horizontal positions of 10m and 4290m at a track pitch of 5m respectively. The seismic source function adopts a Rake wavelet with the peak frequency of 10Hz and the inversion frequency range of 1.5 and 20]Hz. Recording time3s, and a time sampling rate of 1 ms. The inversion frequency sequence obtained by adopting the multi-scale inversion frequency selection strategy is fi1.5,2,3,4,6,8,11,16, 20. FIG. 3 shows an initial velocity model v0The width of the Gaussian window function is 30 and the variance is 10, which are obtained after the Gaussian filtering is carried out on the real model. Under the condition that other parameters are the same, a conventional quasi-Newton method (L-BFGS) and the method of the invention are respectively adopted to carry out full waveform inversion, and the inversion results are shown in figures 4-6. Fig. 4 is a normalized objective function descent curve obtained by single frequency (2.5Hz) inversion using two methods, and it can be seen that the convergence rate of the objective function of the method of the present invention is faster. Fig. 5 is a velocity model obtained by performing multi-scale full-waveform inversion by using two methods, and for comparison, a velocity curve extracted from a position where the lateral position X is 1540m is shown in fig. 6. Comparing the inversion results of the method of the invention with single-frequency iteration for 20 times with the inversion results of the conventional quasi-Newton method with single-frequency iteration for 30 times, the inversion results of the method of the invention are closer to the true values; compared with the inversion results of the conventional quasi-Newton method with different single-frequency iteration times, the inversion accuracy of the conventional quasi-Newton method is not obviously improved by only increasing the single-frequency iteration times. Therefore, the invention can obviously reduce the calculation amount of the full waveform inversion and improve the calculation efficiency of the full waveform inversion on the premise of ensuring the precision.
Claims (4)
1. A full waveform inversion method based on a fast quasi-Newton method is characterized by comprising the following steps:
step 1: inputting and processing data;
the input data comprises actual field earthquake observation data, observation system parameters and an initial parameter model; the data processing mainly comprises the steps of carrying out noise suppression, seismic channel interpolation, amplitude equalization and data transformation on seismic observation data;
obtaining a discrete inversion frequency sequence f ═ f { f by adopting a multi-scale inversion frequency selection strategy according to input datalow，…，fhigh}; carrying out inversion from low frequency to high frequency one by one, and taking the result of the low frequency inversion as an initial model of the high frequency inversion;
step 2: solving a data residual error and an inversion target function;
calculating and storing a seismic wave field in the parameter model, namely a forward wave field, by adopting a wave equation seismic numerical simulation method based on the parameters and the parameter model of the observation system to obtain calculation data; subtracting the calculated data from the observed data to obtain a data residual, and then using a two-norm of the data residual as an inversion target function, namely
Wherein J represents an inversion target function, Δ d is a data residual vector, and "+" represents a conjugate transpose operation;
and step 3: calculating the gradient of an inversion target function;
based on the principle of an adjoint state method, performing conjugate back-propagation simulation on a data residual error serving as a seismic source in a parameter model to obtain a back-propagation wave field, and performing zero-delay cross-correlation operation on a forward-propagation wave field and the back-propagation wave field to obtain the gradient of a target function;
and 4, step 4: calculating the iteration direction of the model update;
based on the algorithm principle of a finite memory quasi-Newton method, obtaining the iteration direction of the model updating by using the gradient residual error and the model residual error information in the last M iteration steps;
and 5: calculating the optimal step length in the iteration direction by adopting a linear search method;
step 6: defining a set of weighting coefficient sequences, i.e.
η0＝0，
And is
Where i is the number of iterations, ηiAnd λiRespectively representing the intermediate variable and the weighting coefficient of the ith iteration;
and 7: introduction ofAn intermediate iteration parameter model yiLet y0＝v0Updating the parametric model according to the following iterative formula
Where i is the number of iterations, viInverse parametric model representing the ith iteration, αiDenotes the step size of the ith iteration and d denotes the direction of the ith iteration.
2. The full waveform inversion method based on the fast quasi-Newton method as claimed in claim 1, wherein in step 2, the wave equation type is selected according to the actual inversion parameters, if only the velocity parameters are inverted, a scalar sound wave equation is used, and the inversion density and the velocity parameters are first order stress-velocity equations.
3. The full waveform inversion method based on the fast quasi-Newton method is characterized in that M is 3-20 in step 4.
4. The full waveform inversion method based on the fast quasi-Newton method as claimed in claim 1, wherein the linear search method in step 5 comprises a non-exact line search method and an exact line search method.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910367245.0A CN110058307B (en) | 2019-05-05 | 2019-05-05 | Full waveform inversion method based on fast quasi-Newton method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910367245.0A CN110058307B (en) | 2019-05-05 | 2019-05-05 | Full waveform inversion method based on fast quasi-Newton method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110058307A CN110058307A (en) | 2019-07-26 |
| CN110058307B true CN110058307B (en) | 2020-12-18 |
Family
ID=67322156
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910367245.0A Expired - Fee Related CN110058307B (en) | 2019-05-05 | 2019-05-05 | Full waveform inversion method based on fast quasi-Newton method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110058307B (en) |
Families Citing this family (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110766628B (en) * | 2019-10-16 | 2020-12-11 | 哈尔滨工程大学 | Target edge inversion method based on multiband self-adaptive regularization iteration |
| CN111580163B (en) * | 2020-05-28 | 2021-08-10 | 中国科学院地质与地球物理研究所 | Full waveform inversion method and system based on non-monotonic search technology |
| CN111722287B (en) * | 2020-06-19 | 2021-10-08 | 南京大学 | Seismic phase characteristic identification waveform inversion method based on progressive data assimilation method |
| CN111856574A (en) * | 2020-07-14 | 2020-10-30 | 中国海洋大学 | Method for constructing velocity component based on seismic exploration pressure component of marine streamer |
| US11971513B2 (en) | 2021-05-21 | 2024-04-30 | Saudi Arabian Oil Company | System and method for forming a seismic velocity model and imaging a subterranean region |
| CN114325829B (en) * | 2021-12-21 | 2023-03-28 | 同济大学 | Full waveform inversion method based on double-difference idea |
Family Cites Families (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103207409B (en) * | 2013-04-17 | 2016-01-06 | 中国海洋石油总公司 | A kind of frequency domain full-waveform inversion seismic velocity modeling method |
| CN105005076B (en) * | 2015-06-02 | 2017-05-03 | 中国海洋石油总公司 | Seismic wave full waveform inversion method based on least square gradient update speed model |
| CN105549079A (en) * | 2016-01-12 | 2016-05-04 | 中国矿业大学(北京) | Method and device for establishing full-waveform inversion model for geophysics parameters |
| CN107505654B (en) * | 2017-06-23 | 2019-01-29 | 中国海洋大学 | Full waveform inversion method based on earthquake record integral |
| CN108873066B (en) * | 2018-06-26 | 2020-03-10 | 中国石油大学(华东) | Elastic medium wave equation reflected wave travel time inversion method |
-
2019
- 2019-05-05 CN CN201910367245.0A patent/CN110058307B/en not_active Expired - Fee Related
Also Published As
| Publication number | Publication date |
|---|---|
| CN110058307A (en) | 2019-07-26 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110058307B (en) | Full waveform inversion method based on fast quasi-Newton method | |
| CN106526674B (en) | Three-dimensional full waveform inversion energy weighting gradient preprocessing method | |
| CN106908835B (en) | Band limit Green's function filters multiple dimensioned full waveform inversion method | |
| KR20140021584A (en) | Convergence rate of full wavefield inversion using spectral shaping | |
| CN110058302A (en) | A kind of full waveform inversion method based on pre-conditional conjugate gradient accelerating algorithm | |
| CN109459789B (en) | Time-domain full waveform inversion method based on amplitude decaying and linear interpolation | |
| CN110954945B (en) | Full waveform inversion method based on dynamic random seismic source coding | |
| CN110187382B (en) | Traveling time inversion method for wave equation of reverse wave and reflected wave | |
| CN108318921A (en) | A kind of quick earthquake stochastic inversion methods based on lateral confinement | |
| CN105549080A (en) | Undulating ground surface waveform inversion method based on auxiliary coordinate system | |
| CN111580163B (en) | Full waveform inversion method and system based on non-monotonic search technology | |
| CN109946742A (en) | The pure rolling land qP shakes digital simulation method in a kind of TTI medium | |
| CN108680968B (en) | Evaluation method and device for seismic exploration data acquisition observation system in complex structural area | |
| CN117331119A (en) | Rapid earthquake wave travel time calculation method for tunnel detection | |
| CN114624766B (en) | Elastic wave least square reverse time migration gradient solving method based on traveling wave separation | |
| CN112630830A (en) | Reflected wave full-waveform inversion method and system based on Gaussian weighting | |
| CN108680957B (en) | Local cross-correlation time-frequency domain Phase-retrieval method based on weighting | |
| Wang et al. | Improved hybrid iterative optimization method for seismic full waveform inversion | |
| CN111914609B (en) | Well-seismic combined prestack geostatistical elastic parameter inversion method and device | |
| CN111208568B (en) | Time domain multi-scale full waveform inversion method and system | |
| CN115166827A (en) | Least square offset imaging method and equipment based on deconvolution imaging conditions and storage medium | |
| CN111175822B (en) | Strong scattering medium inversion method for improving direct envelope inversion and disturbance decomposition | |
| CN114154111A (en) | Bit field data downward continuation method for frequency domain continuous-fraction expansion | |
| CN113866827A (en) | Method, system, medium and device for explanatory velocity modeling seismic imaging | |
| CN110376642B (en) | Three-dimensional seismic velocity inversion method based on conical surface waves |
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 |
Address after: 610032 No. 119, Xiqing Road, Chengdu, Sichuan Applicant after: Sichuan Geological Engineering Exploration Institute Group Co.,Ltd. Address before: 610032 No. 119, Xiqing Road, Chengdu, Sichuan Applicant before: SICHUAN INSTITUTE OF GEOLOGICAL ENGINEERING INVESTIGATION |
|
| CB02 | Change of applicant information | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20201218 Termination date: 20210505 |
|
| CF01 | Termination of patent right due to non-payment of annual fee |