CN116244877B - Three-dimensional magnetic field numerical simulation method and system based on 3D Fourier transform - Google Patents

Three-dimensional magnetic field numerical simulation method and system based on 3D Fourier transform Download PDF

Info

Publication number
CN116244877B
CN116244877B CN202211094934.7A CN202211094934A CN116244877B CN 116244877 B CN116244877 B CN 116244877B CN 202211094934 A CN202211094934 A CN 202211094934A CN 116244877 B CN116244877 B CN 116244877B
Authority
CN
China
Prior art keywords
dimensional
magnetic field
fourier transform
node
domain
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.)
Active
Application number
CN202211094934.7A
Other languages
Chinese (zh)
Other versions
CN116244877A (en
Inventor
戴世坤
张莹
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Central South University
Original Assignee
Central South University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Central South University filed Critical Central South University
Priority to CN202211094934.7A priority Critical patent/CN116244877B/en
Publication of CN116244877A publication Critical patent/CN116244877A/en
Application granted granted Critical
Publication of CN116244877B publication Critical patent/CN116244877B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/14Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
    • G06F17/141Discrete Fourier transforms
    • G06F17/142Fast Fourier transforms, e.g. using a Cooley-Tukey type algorithm
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00Details relating to CAD techniques
    • G06F2111/10Numerical modelling
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/02Reliability analysis or reliability optimisation; Failure analysis, e.g. worst case scenario performance, failure mode and effects analysis [FMEA]
    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02ATECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/30Assessment of water resources
    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02EREDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
    • Y02E60/00Enabling technologies; Technologies with a potential or indirect contribution to GHG emissions mitigation

Abstract

The application 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 of3D 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

Three-dimensional magnetic field numerical simulation method and system based on 3D Fourier transform
Technical Field
The application relates to the field of gravity field numerical simulation, in particular to a three-dimensional magnetic field numerical simulation method and system based on 3D Fourier transform.
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 ofgravityand magnetic anomaly in a mixed space-wavenumber domain, geophysics,2018,Shikun Dai,Dongdong Zhao,etal) 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 of3D 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 application provides a three-dimensional magnetic field numerical simulation method and system based on 3D Fourier transform, which are used for solving the technical problem that the efficiency and accuracy of the existing three-dimensional magnetic field numerical simulation method are not compatible.
In order to solve the technical problems, the technical scheme provided by the application is as follows:
a three-dimensional magnetic field numerical simulation method based on 3D Fourier transform 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 domain 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 is represented, F is a spatial domain function, and F is a wave number spectrum;
for F x (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 Q 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 theIs 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 includeIt 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 samplingThe 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 bitsAnd 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:
k in x 、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:wherein U is 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 For the magnetic susceptibility of the jth node around the corresponding node of the first region, n 1 The number of nodes around the first area; 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 For the magnetic susceptibility of the jth node around the corresponding node of the second region, n 2 The number of nodes around the first area; 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 application has the following beneficial effects:
the application 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 application 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 application has other objects, features and advantages. The application 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 application and are incorporated in and constitute a part of this specification, illustrate embodiments of the application and together with the description serve to explain the application. 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 application;
FIG. 2 is a schematic diagram of a model subdivision in a preferred embodiment of the present application, 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 application;
FIG. 4 is a graph showing a comparison of a sphere numerical solution and an analytical solution in a preferred embodiment of the present application, 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 application, 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 three-dimensional magnetic field numerical simulation method based on 3D fourier transform in a preferred embodiment of the present application.
Detailed Description
Embodiments of the application are described in detail below with reference to the attached drawings, but the application 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 application 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 application 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 k is min 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 for the magnetic susceptibility of the jth node around the corresponding node of the first region, n 1 The number of nodes around the first area; 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 For the magnetic susceptibility of the jth node around the corresponding node of the second region, n 2 The number of nodes around the first area; 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
The upper part (11) is unfolded into
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
K in x 、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 k is x Representing wavenumber, F (x) is a spatial domain function, F (k) x ) Is wave number spectrum.
Discretizing the positive transform integral in equation (14)
Wherein Q 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,
the formula (15) can be written as
Order theIs intra-unit FourierTransforming node coefficients, then formula (18) is abbreviated as
When the wave number is k x If not 0, substituting formula (18) into 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 includeIt 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 carried out to obtain the Fourier transform node coefficient under zero wave number as
And accumulating the analysis expressions of different units to obtain a final one-dimensional Fourier positive transformation result. It is clear that when the spatial and frequency domains are not splitVariable time, fourier transform node coefficient W 1 、W 2 、W 3 And W is 1 0The 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 samplingThe 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 The relation of (2) is thatWherein->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 application, and the arbitrary sampling can be ensured when the gravity numerical simulation of the application is inversely transformed back to the space domain, thereby improving the precision and the efficiency.
The three-dimensional Fourier inverse transformation formula is
K in x 、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 mu 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 application 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 strength 50000nT and magnetic inclination angle 90 DEGThe magnetic declination angle is 4.5 degrees. 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 application and is not intended to limit the present application, but various modifications and variations can be made to the present application by those skilled in the art. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present application should be included in the protection scope of the present application.

Claims (9)

1. The three-dimensional magnetic field numerical simulation method based on the 3D Fourier transform 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; 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 domain magnetic field intensity and the magnetic induction intensity;
the three-dimensional Fourier positive transformation of the arbitrary sampling is three-dimensional arbitrary Fourier positive transformation adopting interpolation based on a shape function.
2. The 3D fourier transform-based three-dimensional magnetic field numerical simulation method as recited in 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, F is a wave number spectrum, and i is an imaginary number;
for F x (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 fourier transform-based three-dimensional magnetic field numerical simulation method according to claim 2, wherein the one-dimensional fourier 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 Q 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 includeIt 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.
4. A three-dimensional magnetic field numerical simulation method based on 3D fourier transform according to claim 3, wherein the wave number domain magnetic bits are obtained by three-dimensional fourier transform of arbitrary samplingThe 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 +.>
5. The 3D fourier transform-based three-dimensional magnetic field numerical simulation method as recited in claim 4, wherein the wave number domain magnetic bitsAnd wave number domain magnetic field strength->The relation of (2) is:
wherein,magnetic field strength vector +.>Components in the x, y, z directions.
6. The 3D fourier transform-based three-dimensional magnetic field numerical simulation method as recited in claim 5, wherein the three-dimensional arbitrary sampling fourier transform formula is:
k in x The wave number, k, representing the x-direction y The wavenumber, k, representing the y-direction z The wavenumber in the z-direction, F (x, y, z) being a spatial domain function, F (k) x ,k y ,k z ) Representing the wavenumber spectrum.
7. The 3D fourier transform-based three-dimensional magnetic field numerical simulation method of any of claims 1-6, wherein the three-dimensional poisson equation is:wherein U is 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.
8. The 3D fourier transform-based three-dimensional magnetic field numerical simulation method according to any one of claims 1 to 6, wherein a 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 For the magnetic susceptibility of the jth node around the corresponding node of the first region, n 1 The number of nodes around the first area; 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 For the magnetic susceptibility of the jth node around the corresponding node of the second region, n 2 The number of nodes around the second area; omega is weight, and the value range is (0, 1);
mode three: the spatial domain samples uniformly in all three directions.
9. 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 8 when the computer program is executed.
CN202211094934.7A 2022-09-05 2022-09-05 Three-dimensional magnetic field numerical simulation method and system based on 3D Fourier transform Active CN116244877B (en)

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 CN116244877A (en) 2023-06-09
CN116244877B true 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)

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

Patent Citations (11)

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

* Cited by examiner, † Cited by third party
Title
ShiKun Dai,et al..Three-Dimensional DC Anisotropic Resistivity Modeling Using a Method in the Mixed Space-Wavenumber Domin.Pure and Applied Geophysics.2022,第2183-2200页. *
各向异性介质空间-波数混合域直流电阻率法三维数值模拟研究;戴世坤等;地球物理学报;第65卷(第7期);第2729-2740页 *
重力异常场空间-波数混合域三维数值模拟;戴世坤等;地球物理学报;第63卷(第5期);第2107-2119页 *

Also Published As

Publication number Publication date
CN116244877A (en) 2023-06-09

Similar Documents

Publication Publication Date Title
Cornwell et al. The noncoplanar baselines effect in radio interferometry: The W-projection algorithm
Mantoglou et al. Simulation of random fields with the turning bands method
Bates et al. Null field approach to scalar diffraction I. General method
CN109856689B (en) Noise suppression processing method and system for superconducting aeromagnetic gradient tensor data
CN106646645B (en) A kind of gravity forward modeling accelerated method
Ergul et al. Broadband multilevel fast multipole algorithm based on an approximate diagonalization of the Green’s function
Fu et al. High‐sensitivity moment magnetometry with the quantum diamond microscope
CN112287534A (en) NUFFT-based two-dimensional magnetic anomaly fast forward modeling method and device
CN114065586B (en) Three-dimensional magnetotelluric space-wavenumber domain finite element numerical simulation method
CN109254327B (en) Exploration method and exploration system of three-dimensional ferromagnetic body
CN115659579B (en) Three-dimensional gravity field numerical simulation method and system based on 3D AS-FT
CN108984939A (en) Three-dimensional Gravity field of force forward modeling method based on 3D Gauss-FFT
Roux et al. Analysis of numerically induced oscillations in two-dimensional finite-element shallow-water models part II: Free planetary waves
Vavilova et al. The Voronoi tessellation method in astronomy
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
Papadakis et al. Spatial filtering in a 6D hybrid-Vlasov scheme to alleviate adaptive mesh refinement artifacts: a case study with Vlasiator (versions 5.0, 5.1, and 5.2. 1)
El Majouti et al. Numerical solution of two-dimensional Fredholm–Hammerstein integral equations on 2D irregular domains by using modified moving least-square method
JP7121111B2 (en) Magnetic field analysis device, analysis method, and program
CN115795231B (en) Space wave number mixed domain three-dimensional strong magnetic field iteration method and system
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
Jasper et al. Development and application of the DIRECT algorithm for leak detection in water distribution systems
CN111665550A (en) Underground medium density information inversion method

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