CN116244877A - Three-dimensional magnetic field numerical simulation method and system based on 3D AS-FT - Google Patents
Three-dimensional magnetic field numerical simulation method and system based on 3D AS-FT Download PDFInfo
- Publication number
- CN116244877A CN116244877A CN202211094934.7A CN202211094934A CN116244877A CN 116244877 A CN116244877 A CN 116244877A CN 202211094934 A CN202211094934 A CN 202211094934A CN 116244877 A CN116244877 A CN 116244877A
- Authority
- CN
- China
- Prior art keywords
- dimensional
- magnetic field
- node
- domain
- wave number
- 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
- 238000000034 method Methods 0.000 title claims abstract description 46
- 238000004088 simulation Methods 0.000 title claims abstract description 30
- 238000005070 sampling Methods 0.000 claims abstract description 50
- 238000001228 spectrum Methods 0.000 claims abstract description 16
- 230000009466 transformation Effects 0.000 claims description 28
- 230000006870 function Effects 0.000 claims description 21
- 230000002159 abnormal effect Effects 0.000 claims description 20
- 230000005415 magnetization Effects 0.000 claims description 12
- 230000006698 induction Effects 0.000 claims description 8
- 238000012887 quadratic function Methods 0.000 claims description 6
- 238000004458 analytical method Methods 0.000 claims description 5
- 239000013598 vector Substances 0.000 claims description 5
- 238000004590 computer program Methods 0.000 claims description 4
- 230000010354 integration Effects 0.000 claims description 4
- 230000014509 gene expression Effects 0.000 claims description 3
- 238000012876 topography Methods 0.000 claims description 2
- 238000004364 calculation method Methods 0.000 abstract description 12
- 230000003595 spectral effect Effects 0.000 abstract description 11
- 238000004422 calculation algorithm Methods 0.000 description 5
- 238000010586 diagram Methods 0.000 description 5
- 230000008859 change Effects 0.000 description 4
- 230000000694 effects Effects 0.000 description 4
- 230000005484 gravity Effects 0.000 description 4
- 230000000737 periodic effect Effects 0.000 description 3
- 230000004048 modification Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000035699 permeability Effects 0.000 description 2
- 230000005856 abnormality Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 210000000746 body region Anatomy 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 238000000926 separation method Methods 0.000 description 1
- 238000004611 spectroscopical analysis Methods 0.000 description 1
- 238000011426 transformation method Methods 0.000 description 1
- 239000002699 waste material Substances 0.000 description 1
Images
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
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
- G06F17/141—Discrete Fourier transforms
- G06F17/142—Fast Fourier transforms, e.g. using a Cooley-Tukey type algorithm
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/02—Reliability analysis or reliability optimisation; Failure analysis, e.g. worst case scenario performance, failure mode and effects analysis [FMEA]
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/30—Assessment of water resources
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02E—REDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
- Y02E60/00—Enabling technologies; Technologies with a potential or indirect contribution to GHG emissions mitigation
Abstract
The invention discloses a three-dimensional magnetic field numerical simulation method and a system based on three-dimensional random sampling Fourier transform (Fourier transform of arbitrary sampling, abbreviated AS AS-FT), wherein a three-dimensional magnetic field differential equation is transformed into a spectral domain to be directly solved by the three-dimensional random sampling Fourier transform method, and then a solution of a spatial domain is obtained by the three-dimensional random sampling Fourier inverse transform. By means of high precision and high efficiency of 3D AS-FT, the Fourier spectrum method of the three-dimensional magnetic field can achieve higher precision, has higher calculation efficiency, can be randomly sampled, and is a great innovation for the traditional spectrum domain numerical simulation method.
Description
Technical Field
The invention relates to the field of gravity field numerical simulation, in particular to a 3D AS-FT-based three-dimensional magnetic field numerical simulation method and system.
Background
The potential field high-precision numerical simulation has important significance for geophysical interpretation and inversion, and the spectral method has wide application in magnetic field numerical simulation, and is characterized in that a solution is approximately expanded into a finite series expansion of a smooth function (generally an orthogonal polynomial), so that an approximate spectral expansion of the solution is obtained, and then an equation set of expansion coefficients is solved according to the expansion and an original equation. Spectral methods are essentially a generalization of standard separation variant techniques. The Chebyshev polynomial and the Legendre polynomial are generally chosen as the basis functions for the approximate expansion. For equations meeting periodic boundary conditions, the Fourier series and the harmonic series are simpler. The spectroscopic methods are classified as Fourier methods, chebyshev or Legendre methods. The former is applicable to periodic problems, and the latter two are applicable to non-periodic problems. The basis of the methods is to establish a space basis function, and the precision of the method directly depends on the number of the expansion terms, so that a plurality of terms are needed to achieve better precision, thereby causing large calculation amount and low efficiency.
The patent adopts a Fourier spectrum method, and the calculation efficiency of the Fourier spectrum method depends on the selected Fourier transform method. However, the accuracy of fourier forward modeling is relatively low due to inherent drawbacks such as aliasing, edge effects, added periodicity, and truncation effects when using a Discrete Fourier Transform (DFT) instead of a continuous fourier transform. If standard FFT is adopted, edge expansion processing is needed to weaken the influence of boundary effect, but the cost is increased in calculation scale; if Gauss-FFT is used for numerical simulation, the efficiency is greatly reduced. And whether standard FFT or Gauss-FFT is used, the sampling can only be uniformly performed, and the waste of calculation capacity is further caused.
Methods for performing magnetic field numerical simulation in the spectral domain by applying fourier transform have been applied in related ways. The literature (Three-dimensional numerical modeling of gravity and magnetic anomaly in a mixed space-wavenumber domain, geophysics,2018,Shikun Dai,Dongdong Zhao,et al) uses standard FFT and Gauss-FFT methods to perform numerical simulation on the magnetic field from the differential equation satisfied by the magnetic potential, and the literature (The forward modeling of 3D gravity and magnetic potential fields in space-wavenumber domains based on an integral method, geophysics,2020, dai Shikun, chen Qingrui, etc.) uses standard FFT and Gauss-FFT methods to perform numerical simulation on the magnetic field from the integral equation; however, it is pointed out that the standard FFT must spread the edges to achieve the same precision as the Gauss-FFT, whereas Gauss-FFT, although having higher precision, has smaller boundary effects, is relatively inefficient, and both standard FFT and Gauss-FFT methods can only sample unevenly, which is not flexible enough.
Disclosure of Invention
The invention provides a 3D AS-FT-based three-dimensional magnetic field numerical simulation method and a 3D AS-FT-based three-dimensional magnetic field numerical simulation system, which are used for solving the technical problem that the efficiency and the accuracy of the existing three-dimensional magnetic field numerical simulation method are incompatible.
In order to solve the technical problems, the technical scheme provided by the invention is as follows:
A3D AS-FT-based three-dimensional magnetic field numerical simulation method comprises the following steps:
constructing a three-dimensional space numerical model of the target area based on the actual topography of the target area and the actual distribution condition of the abnormal body;
performing susceptibility assignment on each node in the three-dimensional space numerical model, and calculating magnetization intensity corresponding to each node;
constructing a three-dimensional poisson equation between the magnetization intensity corresponding to each node and the space domain magnetic position thereof, carrying out three-dimensional Fourier forward transformation of random sampling on the three-dimensional poisson equation, and solving to obtain the wave number domain magnetic position of each node;
solving the wave number domain magnetic field intensity of each node according to the relation between the wave number domain magnetic bit and the wave number domain magnetic field intensity; and performing three-dimensional arbitrary sampling Fourier inverse transformation on the wave number domain magnetic field intensity obtained by solving to obtain space domain magnetic field intensity, and solving to obtain magnetic induction intensity according to the relation between the space magnetic field intensity and the magnetic induction intensity.
Preferably, the three-dimensional fourier transform of the arbitrary samples is:
carrying out x-direction one-dimensional Fourier forward transformation on the three-dimensional poisson equation;
wherein x, y, z represent three coordinate axis directions of a three-dimensional vertical coordinate system; k (k) x The wave number in the x-direction, f isA spatial domain function, F being the wavenumber spectrum;
for F (k) x Y, z) performs a one-dimensional fourier positive transform in the y-direction:
for F xy (k x ,k y Z) performing a one-dimensional fourier positive transform in the z direction:
preferably, the three-dimensional fourier positive transformation of the arbitrary sample is a three-dimensional arbitrary fourier positive transformation using shape function-based interpolation.
Preferably, the one-dimensional fourier positive transform is specifically:
let the continuous one-dimensional fourier positive transforms be expressed as:
discrete one-dimensional fourier positive transform integration is obtained:
wherein M represents the number of units, e j Represents the j-th cell, where i is an imaginary number;
interpolation of f (x) with a quadratic function:
when the quadratic interpolation shape function fitting is adopted in the units, three node coordinates in any unit are respectively set as x 1 、x 2 、x 3 ，x 2 Is the midpoint, satisfy x 1 +x 3 ＝2x 2 Each node has a value f (x 1 )、f(x 2 )、f(x 3 ) F (x) is expressed by using a quadratic function:
f(x)＝N 1 f(x 1 )+N 2 f(x 2 )+N 3 f(x 3 )
wherein ,
the above can be written as:
order the
Is the intra-unit Fourier transform node coefficient, then
The abbreviation is:
when the wave number is k x When the value is not 0, substituting the above formula into W 1 、W 2 、W 3 In (2), the intra-cell fourier transform node coefficients can be obtained:
W 1 、W 2 、W 3 the integral kernel functions all include
It is at [ x ] 1 ,x 3 ]The upper unit integral analysis solution is:
thus, k can be obtained x W when it is not 0 1 、W 2 、W 3 The semi-analytical solution is:
when the wave number is k x When the number of the organic light emitting diode is 0,
and (3) carrying out simple integration to obtain the Fourier transform node coefficient under zero wave number as follows:
and accumulating the analysis expressions of different units to obtain a final one-dimensional Fourier positive transformation result.
Preferably, the wave number domain magnetic bits are obtained by three-dimensional Fourier positive transformation of arbitrary sampling
The equation is satisfied:
wherein ,
magnetization vector in wave number domain>
The components in the x, y and z directions are solved to obtain the wave-number domain magnetic bit +.>
Preferably, the wave number domain magnetic bits
And wave number domain magnetic field strength->
The relation of (2) is:
wherein ,
magnetic field strength vector +.>
Components in the x, y, z directions.
Preferably, the three-dimensional arbitrary sampling fourier transform formula is:
in the formula kx 、k y 、k z Representing wavenumber, F (x, y, z) is a spatial domain function, F (k) x ,k y ,k z ) Representing the wavenumber spectrum.
Preferably, the three-dimensional poisson equation is: (V) 2 U a (x, y, z) =.v.m, where U a For the magnetic potential of the space domain, M is magnetization intensity, and x, y and z represent three coordinate axis directions of a three-dimensional vertical coordinate system.
Preferably, the space domain of the target model is split by adopting any one of the following modes to construct a three-dimensional space numerical model of the target area;
mode one: non-uniform subdivision is adopted for a preset first area, and encryption is carried out; wherein the first region satisfies the following formula:
χ i for the magnetic susceptibility, χ of the corresponding node of the first region j The magnetic susceptibility of the jth node around the corresponding node of the first region is given, and n is the number of the nodes around the first region; omega is weight, and the value range is (0, 1).
Mode two: sparse sampling is carried out on a preset second area, wherein the second area meets the following formula:
χ i′ for the second regionSusceptibility, χ of the corresponding node j′ The magnetic susceptibility of the jth node around the corresponding node of the second region is given, and n is the number of the nodes around the first region; omega is weight, and the value range is (0, 1);
mode three: the spatial domain samples uniformly in all three directions.
A computer system comprising a memory, a processor and a computer program stored on the memory and executable on the processor, the processor implementing the steps of the method described above when the computer program is executed.
The invention has the following beneficial effects:
the invention applies a three-dimensional random sampling Fourier transform method (3D AS-FT), transforms a three-dimensional magnetic field differential equation to a spectral domain to directly solve, and obtains a solution of a spatial domain through three-dimensional random sampling Fourier inverse transform. The three-dimensional arbitrary sampling Fourier forward transformation is applied to the transformation from the space domain to the wave number domain, so that the space domain in the x, y and z directions can be split unevenly. The abnormal body can be better fitted, and the calculation accuracy is improved. Sparse sampling is performed in the area where encryption subdivision is not needed, so that the calculation efficiency is improved. The wave number domain is inversely transformed back to the spatial domain by using a three-dimensional arbitrary sampling Fourier inverse transformation method, so that arbitrary interval sampling of the wave number is realized. The non-uniform sampling of the wave number can encrypt and split the region with strong energy and the region with weak energy and slow frequency spectrum change, thereby improving the accuracy and efficiency of the algorithm. The invention solves the abnormal magnetic bit by carrying out three-dimensional Fourier transform on poisson equation satisfied by the magnetic bit to a spectral domain and then changing the three-dimensional Fourier transform to a solution of a spatial domain. The equation is simple, and the problems of overlarge calculation amount and lower efficiency of solving the magnetic field by the traditional Fourier spectrum method are greatly improved by means of high efficiency and high precision of any Fourier transform.
In addition to the objects, features and advantages described above, the present invention has other objects, features and advantages. The invention will be described in further detail with reference to the accompanying drawings.
Drawings
The accompanying drawings, which are included to provide a further understanding of the invention and are incorporated in and constitute a part of this specification, illustrate embodiments of the invention and together with the description serve to explain the invention. In the drawings:
FIG. 1 is a schematic view of a calculation region and an abnormal body region in a preferred embodiment of the present invention;
FIG. 2 is a schematic diagram of a model subdivision in a preferred embodiment of the present invention, wherein (a) represents a model in which x, y, and z are all uniformly sampled, (b) represents a model in which x, y, and z are uniformly sampled, and (c) represents a model in which x, y, and z are all non-uniformly sampled;
FIG. 3 is a schematic diagram of a unit and nodes within the unit in a preferred embodiment of the invention;
FIG. 4 is a graph showing a comparison of a sphere numerical solution and an analytical solution in a preferred embodiment of the present invention, wherein (a) represents B ax Numerical solution, (B) graph represents B ax Analytical solution, (c) graph represents B ax Absolute error of a numerical solution compared to an analytical solution; (d) Representation B ay Numerical solution, (e) graph represents B ay Analytical solution, (f) graph represents B ay Absolute error of a numerical solution compared to an analytical solution; (g) Representation B az Numerical solution, (h) graph represents B az Analytical solution, (i) graph representation B az Absolute error of a numerical solution compared to an analytical solution;
FIG. 5 is a cross-sectional view of an absolute error in a preferred embodiment of the invention, wherein (a) represents B ax Absolute error cut-off at y=0m, (B) graph represents B ay Absolute error cut-off at x=0m, (c) represents B az Absolute error cut-off plot at y=0m;
FIG. 6 is a flow chart of a 3D AS-FT based three-dimensional magnetic field numerical simulation method in a preferred embodiment of the present invention.
Detailed Description
Embodiments of the invention are described in detail below with reference to the attached drawings, but the invention can be implemented in a number of different ways, which are defined and covered by the claims.
Embodiment one:
AS shown in fig. 6, the three-dimensional secondary magnetic field Fourier spectrum method of the 3D AS-FT based on shape function interpolation provided by the invention comprises the following steps:
step one: model building
And (5) performing geological modeling on the calculated area of the numerical simulation. The size of the whole calculation region is firstly determined, and then the distribution of abnormal bodies, which can be any complex condition and any complex shape, are determined, wherein the abnormal bodies are required to be in the calculation region. A schematic diagram of a simple model is shown in FIG. 1, in which the anomaly is a sphere.
Step two: model subdivision
Modeling in a spatial domain:
after the model is built, the model is split, and the sampling points in the x, y and z directions are Nx, ny and Nz respectively. One of the advantages of the invention is that the model subdivision is arbitrary in the x, y and z directions, the non-uniform subdivision can be adopted for encryption in a first area with a faster abnormal body change of the model, and sparse sampling can be carried out in a second area with a slower abnormal body change or without change. It is also possible that all three directions are uniformly sampled. The three-dimensional schematic diagrams of the sampling are respectively shown in fig. 2a, b and c, wherein a is uniformly sampled in three directions; b is non-uniform sampling in the x-y direction and uniform sampling in the z-direction; c is the non-uniform sampling in the x, y and z directions. Determining wavenumber k from spatial domain subdivision x ，k y ，k z The cut-off frequency (k) x ，k y ，k z Maximum positive and minimum negative of (a) and k x ，k y ，k z Is a sampling mode of the system.
The cut-off frequency is related to the minimum subdivision interval in the corresponding direction of the space domain, and the minimum subdivision interval in the x direction is set as Deltax min The minimum subdivision interval in the y direction is delta y min The z-direction minimum subdivision interval is deltaz min The corresponding cut-off frequency is
Sampling within the cut-off frequency can ensure that all spectral information is sampled. After determining the cut-off frequency, the number of samples is determined again, assuming k x ，k y ，k z The number of samples is Nkx, nky, nkz, respectively.
Can select uniform sampling, namely k x ，k y ，k z The arrangement intervals are the same; alternatively, the sampling may be uniform in the logarithmic domain.
When sampling logarithmic interval, the wave number is set to be in the range of [ -k max ,k max ]The number of the wave number domain sampling points is 2M+1, the sampling is carried out on the digital domain sampling at equal intervals, and the sampling interval is
wherein ,kmin Is a decimal, typically 10 -6 ～10 -3 。
The wave number is arranged at [ -k max ,0]Upper is
Wavenumbers are arranged in [0, k ] max ]Upper is
k x ，k y ，k z The logarithmic domain sampling of (a) can be summed by a sampling pattern, thereby giving an arrangement of the spatial domain x, y, z and the wavenumber domain kx, ky, kz.
In this embodiment, the first region satisfies the following formula:
χ i for the magnetic susceptibility, χ of the corresponding node of the first region j The magnetic susceptibility of the jth node around the corresponding node of the first region is given, and n is the number of the nodes around the first region; omega is weight, and the value range is (0, 1).
The second region satisfies the following formula:
χ i′ for the magnetic susceptibility, χ of the corresponding node of the second region j′ The magnetic susceptibility of the jth node around the corresponding node of the second region is given, and n is the number of the nodes around the first region; omega is weight, and the value range is (0, 1); in this embodiment, ω has a value of 0.2.
Step three: model susceptibility assignment
The susceptibility assignment is made to the nodes in fig. 2. The abnormal body part is assigned to each corresponding node according to the magnetic susceptibility of the abnormal body, the magnetic susceptibility on the node without the abnormal part is 0, the magnetic susceptibility is represented by χ, the scalar is represented by SI, and the unit is SI.
Step four: calculating magnetization M corresponding to the node
Calculating the earth main magnetic field intensity H at each node according to the earth main magnetic field model IGRF 0 The unit of the background field in the numerical simulation, namely the magnetic field without abnormality is A/m, and the components in three directions are respectively expressed as H 0x 、H 0y 、H 0z And this main magnetic field value is taken as the magnetic field initial value. The strength of the magnetic field generated by the abnormal body at the node is H a The unit is A/m, which is an abnormal field in numerical simulation, namely a magnetic field generated by abnormal magnetic susceptibility, and the abnormal field is ignored when the weak magnetic condition, namely the magnetic susceptibility is smaller than 0.1 SI.
The three components of the background field are calculated by the following formula, where H 0 The I represents the background field H 0 And (2) is the L2 norm of the study area, alpha is the study area dip angle, and beta is the study area declination angle.
H 0x ＝||H 0 ||·cos(α)·cos(β) (7)
H 0y ＝||H 0 ||·cos(α)·sin(β) (8)
H 0z ＝||H 0 ||·sin(α) (9)
Thus, H of each node is obtained 0 After that, the magnetization M is calculated from
M＝χH＝χ(H 0 +H a ) (10)
Wherein M is A/M and H a And neglected.
Step five: by aligning abnormal magnetic position U a Performing three-dimensional Fourier transform on the satisfied three-dimensional poisson equation to obtain wave number domain magnetic bits
Space domain magnetic potential U a And magnetization M satisfies the Poisson equation
▽ 2 U a (x,y,z)＝▽·M (11)
The upper expansion is
And performing three-dimensional Fourier forward transformation, wherein the three-dimensional Fourier forward transformation adopts three-dimensional arbitrary Fourier forward transformation based on shape function interpolation, and the three-dimensional Fourier forward transformation principle is as follows:
the three-dimensional Fourier positive transformation formula is
in the formula kx 、k y 、k z Representing wavenumber, F (x, y, z) is a spatial domain function, F (k) x ,k y ,k z ) Representing the wavenumber spectrum.
The three-dimensional transformation is completed through three one-dimensional Fourier positive transformation, and the principle of the one-dimensional Fourier positive transformation is described first.
The one-dimensional fourier positive transforms can be expressed as:
wherein kx Representing wavenumber, F (x) is a spatial domain function, F (k) x ) Is wave number spectrum.
Discrete-obtaining the positive transform integral in
Wherein M represents the number of units, e j Represents the j-th cell, where i is an imaginary number.
F (x) is interpolated by a quadratic function. When the quadratic interpolation shape function fitting is adopted in the units, three node coordinates in any unit are respectively set as x 1 、x 2 、x 3 ，x 2 Is the midpoint, satisfy x 1 +x 3 ＝2x 2 The intra-cell nodes are as shown in fig. 3:
each node has a value f (x 1 )、f(x 2 )、f(x 3 ) F (x) is obtained by using quadratic function to represent
f(x)＝N 1 f(x 1 )+N 2 f(x 2 )+N 3 f(x 3 ) (16)
wherein ,
can be written as
Order the
Is the intra-unit Fourier transform node coefficient, abbreviated as
When the wave number is k x When the value is not 0, substituting W 1 、W 2 、W 3 In the method, the intra-cell Fourier transform node coefficients can be obtained
W 1 、W 2 、W 3 The integral kernel functions all include
It is at [ x ] 1 ,x 3 ]Upper unit integral resolution
Thus, k can be obtained x W when it is not 0 1 、W 2 、W 3 Semi-analytical solution into
When the wave number is k x When the number of the organic light emitting diode is 0,
simple integration is performed to obtain Fourier transform node coefficient of +.>
And accumulating the analysis expressions of different units to obtain a final one-dimensional Fourier positive transformation result. It is easy to know that when the space domain and frequency domain subdivision is unchanged, the Fourier transform node coefficient W 1 、W 2 、W 3 And
the Fourier transform coefficients are calculated and stored in advance, so that repeated calculation can be reduced, the algorithm efficiency is improved, and the method is one of the advantages of the algorithm.
The three-dimensional Fourier transform is that the one-dimensional Fourier transform is completed on x
After that for F (k) x Y, z) performs a one-dimensional fourier transform in the y-direction:
finally to F xy (k x ,k y Z) performing a z-direction one-dimensional fourier transform:
the principle of the three-dimensional fourier transform is exactly the same as that of the process, and thus will not be described in detail.
Obtaining wave number domain magnetic bits through three-dimensional Fourier transformation of arbitrary sampling
The following equation is satisfied:
wherein ,
magnetization in the wavenumber domain>
Components in x, y, z directions.
Directly solving to obtain wave number domain magnetic position
Step six: solving the wave number domain magnetic field strength according to the relation between the magnetic potential and the magnetic field strength
Space domain magnetic potential U a And spatial domain magnetic field strength H a Is of the relationship H a ＝-▽U a, wherein />
i, j, k are unit vectors in the x, y, z directions, respectively. Thus the wave number domain magnetic bits are available>
And wave number domain magnetic field strength->
The relation of (2) is that
According to the method, the wave number domain magnetic field intensity can be obtained by solving
Step seven: the space domain magnetic field intensity H is obtained through the three-dimensional Fourier inverse transformation of arbitrary sampling a 。
The application of arbitrary sampling three-dimensional Fourier inverse transformation is also a great innovation of the invention, and the arbitrary sampling can be ensured when the gravity numerical simulation of the invention is inversely transformed back to the space domain, thereby improving the precision and the efficiency.
The three-dimensional Fourier inverse transformation formula is
in the formula kx 、k y 、k z Representing wavenumber, F (x, y, z) is a spatial domain function, F (k) x ,k y ,k z ) Representing the wavenumber spectrum. The form of the reverse transformation formula is identical to that of the forward transformation formula, the principle is also identical, and the description is omitted.
Step eight: solving to obtain magnetic induction intensity B a 。
From abnormal field magnetic induction intensity B a (in T) and the field strength H of the abnormal field a Can be used to determine the magnetic induction B a Further obtain B a Is defined by three components B ax ，B ay ，B az 。
B a ＝μH a (30)
Wherein mu is the absolute permeability of the medium and the unit is H/m. The relation between mu and χ satisfies the following formula
μ＝μ 0 (1+χ) (31)
wherein μ0 For permeability in vacuum, mu 0 ＝4π×10 -7 H/m。
The accuracy and efficiency of the three-dimensional magnetic field Fourier spectrum method based on the arbitrary sampling three-dimensional Fourier transform provided by the invention are checked.
The test computer is configured as i7-11800, the main frequency is 2.30GHz, and the memory is 32GB.
Designing a sphere model, 1000m×1000m, and ranging: x direction-500 m, y direction-500 m, z direction 0-1000 m. The abnormal body sphere model is positioned at the center (0 m,500 m), the sphere radius is 100m, the sphere magnetic susceptibility is 0.01SI, and the weak magnetic condition is adopted. The three directions of the space domain are evenly split, and 101 nodes are taken. Background magnetic field intensity 50000nT, magnetic dip angle 90 degrees and magnetic bias angle4.5 deg.. The sampling range of the wave number domain is based on the sampling theorem, and the maximum wave number is the cut-off frequency
Δx is a constant value of 10m, so the maximum wave number is pi/10, the mode of uniformly sampling the logarithmic domain is adopted, and the logarithmic minimum number is 10 -4 And the number of sampling points in three directions is 101. Calculated ground field value B ax ,B ay ,B az For example, as shown in fig. 4, the relative root mean square error is 0.41%, 0.41% and 0.23% respectively. B (B) ax The absolute error cut-off at y=0m is shown in fig. 5a, B ay The absolute error cut-off at x=0m is shown in fig. 5B, B az As shown in fig. 5c, the absolute error section diagram at y=0m shows that the absolute errors are different by more than two orders of magnitude from the analytical solutions, and the accuracy requirement is met. The memory occupies 1.02GB and takes 1.34s.
In summary, the three-dimensional magnetic field Fourier spectrum method based on the three-dimensional arbitrary sampling Fourier transform of the shape function interpolation carries out Fourier transform on the Poisson equation satisfied by the magnetic potential through the three-dimensional Fourier transform, directly obtains the magnetic potential in a spectral domain, then obtains the magnetic field intensity in the spectral domain by utilizing the relation between the magnetic potential in the spectral domain and the magnetic field intensity, and can solve the magnetic field intensity in the spatial domain by carrying out three-dimensional inverse transform on the magnetic field intensity, and finally obtains the magnetic induction intensity. The algorithm has good parallelism, occupies less memory, can randomly sample the space domain and the wave number domain, and has good guiding significance for the condition of needing non-uniform sampling. The algorithm has higher precision and efficiency by means of arbitrary sampling Fourier transform.
The above description is only of the preferred embodiments of the present invention and is not intended to limit the present invention, but various modifications and variations can be made to the present invention by those skilled in the art. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention should be included in the protection scope of the present invention.
Claims (10)
1. The 3D AS-FT-based three-dimensional magnetic field numerical simulation method is characterized by comprising the following steps of:
constructing a three-dimensional space numerical model of the target area based on the actual topography of the target area and the actual distribution condition of the abnormal body;
performing susceptibility assignment on each node in the three-dimensional space numerical model, and calculating magnetization intensity corresponding to each node;
constructing a three-dimensional poisson equation between the magnetization intensity corresponding to each node and the space domain magnetic position thereof, carrying out three-dimensional Fourier forward transformation of random sampling on the three-dimensional poisson equation, and solving to obtain the wave number domain magnetic position of each node;
solving the wave number domain magnetic field intensity of each node according to the relation between the wave number domain magnetic bit and the wave number domain magnetic field intensity; and performing three-dimensional arbitrary sampling Fourier inverse transformation on the wave number domain magnetic field intensity obtained by solving to obtain space domain magnetic field intensity, and solving to obtain magnetic induction intensity according to the relation between the space magnetic field intensity and the magnetic induction intensity.
2. The 3D AS-FT based three-dimensional magnetic field numerical simulation method of claim 1, wherein the arbitrarily sampled three-dimensional fourier transform is:
carrying out x-direction one-dimensional Fourier forward transformation on the three-dimensional poisson equation;
wherein x, y, z represent three coordinate axis directions of a three-dimensional vertical coordinate system; k (k) x The wave number in the x direction is represented, F is a spatial domain function, and F is a wave number spectrum;
for F (k) x Y, z) performs a one-dimensional fourier positive transform in the y-direction, k y Wavenumber representing y-direction:
for F xy (k x ,k y Z) performing a one-dimensional Fourier forward transform in the z direction, k z Wavenumber representing z-direction:
3. the 3D AS-FT based three-dimensional magnetic field numerical simulation method of claim 2, wherein the arbitrary sampled three-dimensional fourier transform is a three-dimensional arbitrary fourier transform employing shape function-based interpolation.
4. A 3D AS-FT based three-dimensional magnetic field numerical simulation method according to claim 3, characterized in that the one-dimensional fourier positive transform is specifically:
let the continuous one-dimensional fourier positive transforms be expressed as:
discrete the continuous one-dimensional fourier transform integral to obtain:
wherein M represents the number of units, e j Represents the j-th cell, where i is an imaginary number;
interpolation of f (x) with a quadratic function:
when the quadratic interpolation shape function fitting is adopted in the units, three node coordinates in any unit are respectively set as x 1 、x 2 、x 3 ，x 2 Is the midpoint, satisfy x 1 +x 3 ＝2x 2 Each node has a value f (x 1 )、f(x 2 )、f(x 3 ) F (x) is expressed by using a quadratic function:
f(x)＝N 1 f(x 1 )+N 2 f(x 2 )+N 3 f(x 3 )
wherein ,
the above variants are written as:
order the
For the intra-cell fourier transform node coefficients, the above abbreviation is:
when the wave number is k x When the value is not 0, substituting the above formula into W 1 、W 2 、W 3 Obtaining the intra-unit Fourier transform node coefficient:
W 1 、W 2 、W 3 the integral kernel functions all include
It is at [ x ] 1 ,x 3 ]The upper unit integral analysis solution is:
thus, get k x W when it is not 0 1 、W 2 、W 3 The semi-analytical solution is:
when the wave number is k x When the number of the organic light emitting diode is 0,
and (3) carrying out simple integration to obtain the Fourier transform node coefficient under zero wave number as follows: />
And accumulating the analysis expressions of the different units to obtain a final one-dimensional Fourier positive transformation result.
5. A 3D AS-FT based three-dimensional magnetic field numerical simulation method according to claim 3, characterized in that the wave number domain magnetic bits are obtained by arbitrary sampled three-dimensional fourier positive transformation
The equation is satisfied:
wherein ,
magnetization vector in wave number domain>
The components in the x, y and z directions are solved to obtain the wave-number domain magnetic bit +.>
6. The 3D AS-FT based three-dimensional magnetic field numerical simulation method of claim 5 wherein the wave number domain magnetic bits
And wave number domain magnetic field strength->
The relation of (2) is:
wherein ,
magnetic field strength vector +.>
Components in the x, y, z directions.
7. The 3D AS-FT based three-dimensional magnetic field numerical simulation method of claim 6, wherein the three-dimensional arbitrary sampling fourier transform formula is:
in the formula kx 、k y 、k z Representing wavenumber, F (x, y, z) is a spatial domain function, F (k) x ,k y ,k z ) Representing the wavenumber spectrum.
8. The 3D AS-FT based three-dimensional magnetic field numerical simulation method of any one of claims 1-7, wherein the three-dimensional poisson equation is: (V) 2 U a (x, y, z) =.v.m, where U a For the magnetic potential of the space domain, M is magnetization intensity, and x, y and z represent three coordinate axis directions of a three-dimensional vertical coordinate system.
9. The 3D AS-FT based three-dimensional magnetic field numerical simulation method according to any one of claims 1-7, wherein the spatial domain of the target model is split in any one of the following ways to construct a three-dimensional spatial numerical model of the target region;
mode one: non-uniform subdivision is adopted for a preset first area, and encryption is carried out; wherein the first region satisfies the following formula:
χ i for the magnetic susceptibility, χ of the corresponding node of the first region j The magnetic susceptibility of the jth node around the corresponding node of the first region is given, and n is the number of the nodes around the first region; omega is weight, and the value range is (0, 1);
mode two: sparse sampling is carried out on a preset second area, wherein the second area meets the following formula:
χ i′ for the magnetic susceptibility, χ of the corresponding node of the second region j′ The magnetic susceptibility of the jth node around the corresponding node of the second region is given, and n is the number of the nodes around the first region; omega is weight, and the value range is (0, 1);
mode three: the spatial domain samples uniformly in all three directions.
10. A computer system comprising a memory, a processor and a computer program stored on the memory and executable on the processor, characterized in that the processor implements the steps of the method of any of the preceding claims 1 to 9 when the computer program is executed.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202211094934.7A CN116244877B (en) | 2022-09-05 | 2022-09-05 | Three-dimensional magnetic field numerical simulation method and system based on 3D Fourier transform |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202211094934.7A CN116244877B (en) | 2022-09-05 | 2022-09-05 | Three-dimensional magnetic field numerical simulation method and system based on 3D Fourier transform |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN116244877A true CN116244877A (en) | 2023-06-09 |
| CN116244877B CN116244877B (en) | 2023-11-14 |
Family
ID=86631864
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202211094934.7A Active CN116244877B (en) | 2022-09-05 | 2022-09-05 | Three-dimensional magnetic field numerical simulation method and system based on 3D Fourier transform |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN116244877B (en) |
Citations (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US6414487B1 (en) * | 2000-11-22 | 2002-07-02 | Philips Medical Systems (Cleveland), Inc. | Time and memory optimized method of acquiring and reconstructing multi-shot 3D MRI data |
| JP2004093499A (en) * | 2002-09-03 | 2004-03-25 | Nippon Steel Corp | Magnetic field analyzing method, device, computer program, and computer readable storage medium |
| CN109254327A (en) * | 2018-10-30 | 2019-01-22 | 桂林理工大学 | The exploitation method and exploration system of three-dimensional ferromagnetic |
| CN111967169A (en) * | 2020-10-21 | 2020-11-20 | 中南大学 | Two-degree body weight abnormal product decomposition numerical simulation method and device |
| CN112287534A (en) * | 2020-10-21 | 2021-01-29 | 中南大学 | NUFFT-based two-dimensional magnetic anomaly fast forward modeling method and device |
| CN112949134A (en) * | 2021-03-09 | 2021-06-11 | 吉林大学 | Earth-well transient electromagnetic inversion method based on non-structural finite element method |
| CN113673163A (en) * | 2021-08-25 | 2021-11-19 | 中南大学 | Three-dimensional magnetic anisotropy constant fast forward modeling method and device and computer equipment |
| CN113962077A (en) * | 2021-10-20 | 2022-01-21 | 中南大学 | Three-dimensional anisotropic strong magnetic field numerical simulation method, device, equipment and medium |
| CN114021408A (en) * | 2021-11-05 | 2022-02-08 | 中南大学 | Two-dimensional high-intensity magnetic field numerical simulation method, device, equipment and medium |
| CN114065586A (en) * | 2021-11-22 | 2022-02-18 | 中南大学 | Three-dimensional magnetotelluric space-wavenumber domain finite element numerical simulation method |
| CN114065585A (en) * | 2021-11-22 | 2022-02-18 | 中南大学 | Three-dimensional electrical source numerical simulation method based on coulomb specification |
-
2022
- 2022-09-05 CN CN202211094934.7A patent/CN116244877B/en active Active
Patent Citations (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US6414487B1 (en) * | 2000-11-22 | 2002-07-02 | Philips Medical Systems (Cleveland), Inc. | Time and memory optimized method of acquiring and reconstructing multi-shot 3D MRI data |
| JP2004093499A (en) * | 2002-09-03 | 2004-03-25 | Nippon Steel Corp | Magnetic field analyzing method, device, computer program, and computer readable storage medium |
| CN109254327A (en) * | 2018-10-30 | 2019-01-22 | 桂林理工大学 | The exploitation method and exploration system of three-dimensional ferromagnetic |
| CN111967169A (en) * | 2020-10-21 | 2020-11-20 | 中南大学 | Two-degree body weight abnormal product decomposition numerical simulation method and device |
| CN112287534A (en) * | 2020-10-21 | 2021-01-29 | 中南大学 | NUFFT-based two-dimensional magnetic anomaly fast forward modeling method and device |
| CN112949134A (en) * | 2021-03-09 | 2021-06-11 | 吉林大学 | Earth-well transient electromagnetic inversion method based on non-structural finite element method |
| CN113673163A (en) * | 2021-08-25 | 2021-11-19 | 中南大学 | Three-dimensional magnetic anisotropy constant fast forward modeling method and device and computer equipment |
| CN113962077A (en) * | 2021-10-20 | 2022-01-21 | 中南大学 | Three-dimensional anisotropic strong magnetic field numerical simulation method, device, equipment and medium |
| CN114021408A (en) * | 2021-11-05 | 2022-02-08 | 中南大学 | Two-dimensional high-intensity magnetic field numerical simulation method, device, equipment and medium |
| CN114065586A (en) * | 2021-11-22 | 2022-02-18 | 中南大学 | Three-dimensional magnetotelluric space-wavenumber domain finite element numerical simulation method |
| CN114065585A (en) * | 2021-11-22 | 2022-02-18 | 中南大学 | Three-dimensional electrical source numerical simulation method based on coulomb specification |
Non-Patent Citations (3)
| Title |
|---|
| SHIKUN DAI, ET AL.: "Three-Dimensional DC Anisotropic Resistivity Modeling Using a Method in the Mixed Space-Wavenumber Domin", PURE AND APPLIED GEOPHYSICS, pages 2183 - 2200 * |
| 戴世坤等: "各向异性介质空间－波数混合域直流电阻率法三维数值模拟研究", 地球物理学报, vol. 65, no. 7, pages 2729 - 2740 * |
| 戴世坤等: "重力异常场空间-波数混合域三维数值模拟", 地球物理学报, vol. 63, no. 5, pages 2107 - 2119 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN116244877B (en) | 2023-11-14 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Cornwell et al. | The noncoplanar baselines effect in radio interferometry: The W-projection algorithm | |
| Sharan et al. | Application of the multiquadric method for numerical solution of elliptic partial differential equations | |
| CN112287534B (en) | NUFFT-based two-dimensional magnetic anomaly fast forward modeling method and device | |
| CN109856689B (en) | Noise suppression processing method and system for superconducting aeromagnetic gradient tensor data | |
| Fu et al. | High‐sensitivity moment magnetometry with the quantum diamond microscope | |
| CN114065586B (en) | Three-dimensional magnetotelluric space-wavenumber domain finite element numerical simulation method | |
| Lam et al. | A novel meshless approach–Local Kriging (LoKriging) method with two-dimensional structural analysis | |
| CN108984939A (en) | Three-dimensional Gravity field of force forward modeling method based on 3D Gauss-FFT | |
| CN115659579B (en) | Three-dimensional gravity field numerical simulation method and system based on 3D AS-FT | |
| Roux et al. | Analysis of numerically induced oscillations in two-dimensional finite-element shallow-water models part II: Free planetary waves | |
| CN115292973B (en) | Arbitrarily sampled space wave number domain three-dimensional magnetic field numerical simulation method and system | |
| CN116244877B (en) | Three-dimensional magnetic field numerical simulation method and system based on 3D Fourier transform | |
| Leung et al. | A level set method for three-dimensional paraxial geometrical optics with multiple point sources | |
| CN109143360B (en) | Method for determining seismic event P wave inverse azimuth angle and slowness with high resolution | |
| CN114325870A (en) | Method and system for calculating potential field gradient tensor based on cubic spline function | |
| CN115795231B (en) | Space wave number mixed domain three-dimensional strong magnetic field iteration method and system | |
| CN115267673A (en) | Sparse sound source imaging method and system considering reconstruction grid offset | |
| Guo et al. | Microwave inversion for sparse data using descent learning technique | |
| Kämmerer et al. | Computational methods for the Fourier analysis of sparse high-dimensional functions | |
| CN113267830A (en) | Two-dimensional gravity gradient and seismic data joint inversion method based on non-structural grid | |
| CN111665550A (en) | Underground medium density information inversion method | |
| CN107562690A (en) | A kind of AIM PO mixed methods based on GPU | |
| Nengem et al. | A Study of the Effects of the Shape Parameter and Type of Data Points Locations on the Accuracy of the Hermite-Based Symmetric Approach Using Positive Definite Radial Kernels | |
| Otmianowska-Mazur et al. | Magnetic field in a turbulent galactic disk | |
| CN116794735B (en) | Aviation magnetic vector gradient data equivalent source multi-component joint denoising method and device |
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 |