CN113656750B - Calculation method of magnetic induction intensity of ferromagnetic medium based on spatial wavenumber domain - Google Patents

Abstract

Based on the calculation method of magnetic induction intensity of ferromagnetic medium in the space wavenumber domain, the initial prism model containing the abnormal body is divided, a series of nodes are obtained, and the magnetic susceptibility is assigned to each node. Calculate the discrete migration wavenumber, calculate the magnetic field strength of the background field in the space domain, and obtain the magnetization model; convert the equation satisfying the magnetic potential and magnetization of the anomalous field in the space domain into a one-dimensional ordinary differential equation in the space wavenumber domain; combine the set space wavenumber domain The boundary condition of the magnetic potential of the anomalous field, the boundary value problem model satisfied by the magnetic potential of the anomalous field in the space wavenumber domain is equivalently transformed into a variational problem model and solved, and the magnetic potential of the anomalous field in the space wavenumber domain and the magnetic field strength of the anomalous field in the space wavenumber domain are obtained. It is converted into the magnetic potential of the anomalous field in the space domain and the magnetic field intensity of the anomalous field in the space domain, and the magnetic induction intensity in the space domain is obtained according to the output magnetic field intensity of the anomalous field in the space domain. The invention can carry out numerical simulation of the magnetic field for the ferromagnetic medium more accurately.

Description

Calculation method of magnetic induction intensity of ferromagnetic medium based on spatial wavenumber domain

technical field

The invention belongs to the technical field of numerical simulation of strong magnets, and particularly relates to a method for calculating the magnetic induction intensity of a strong magnetic medium based on a spatial wave number domain.

Background technique

When the magnetic susceptibility value of the medium is greater than 0.1 SI, it is generally regarded as a ferromagnetic medium, and the self-demagnetizing field in the magnetic field formed by it cannot be ignored. At present, the calculation of the magnetic induction intensity of the magnetic field of the strong magnetic medium usually ignores its self-demagnetizing field, that is, it is approximated as a weak magnetic field, which leads to a large deviation between the magnetic induction intensity obtained by numerical simulation and the actual magnetic induction intensity value, which is not conducive to the measurement of magnetic data. processing and interpretation.

For the magnetic field calculation of ferromagnetic media, the research methods of related scholars are mainly divided into spatial domain and frequency domain methods. In the space domain method, the literature (Vogel A. The application of electronic computers to the calculation of effective magnetization [J]. Geophysical Prospecting. 1963.11(1):51-58.) proposed to use a vector series to represent the magnitude of each volume element inside the magnet. The effective magnetization is solved by successive approximation, but when the magnetic susceptibility is high, the vector series converges very slowly or even diverges, and the calculation efficiency is low. In the literature (Purss, M. B. J., and J. P. Cull, A new iterative method for computing the magnetic field at high magnetic susceptibilities: Geophysics, 2005. 70: 53-62.) through a segmented model defined by spherical voxels of uniform arbitrary diameter, in The magnetic field of the high susceptibility magnetic body is calculated on an iterative basis, but the error of the high susceptibility calculation results is large. In the frequency domain method, literature (Wu, L. Y., and G.Tian, High-precision Fourier forward modeling of potential fields: Geophysics, 2014. 79, no. 5, G59–G68, doi: 10.1190/geo2014-0039.1.) Based on the shift sampling method, a new Gaussian FFT method is proposed, which is successfully applied to the numerical simulation of gravity and magnetic methods, but it is still only suitable for the case of low magnetic susceptibility. The literature (Li Kun, Dai Shikun, Chen Qingrui, etc.. Three-dimensional numerical simulation of the integral solution of magnetic anomalies in the mixed domain of spatial wavenumbers. Acta Geophysics. 2019. 62(11)) uses the characteristics of the three-dimensional spatial domain integral of magnetic potential as convolution, along the Two-dimensional Fourier transform is performed in the horizontal direction, and the three-dimensional integration problem satisfied by the magnetic potential in the space domain is transformed into a vertical one-dimensional integration problem that is independent of each other between different wavenumbers, but it is only suitable for the case of low magnetic susceptibility.

Based on the current research status above, it can be seen that the current magnetic induction intensity calculation methods rarely consider strong magnetism, and some calculation methods that consider strong magnetism have low efficiency and low precision. Therefore, it is imperative to propose a high-efficiency and high-precision calculation method of magnetic induction intensity suitable for strong magnets.

SUMMARY OF THE INVENTION

Aiming at the problem that the self-demagnetization effect is less considered in the current calculation of the magnetic induction intensity of strong magnets, and the few studies that consider the self-demagnetization effect have the problems of insufficient precision and low efficiency, the present invention aims to propose a space wavenumber domain-based method for a strong magnetic medium. Magnetic induction intensity calculation method.

For realizing the above-mentioned technical purpose, the technical scheme proposed by the present invention is:

In one aspect, the present invention provides a method for calculating the magnetic induction intensity of a ferromagnetic medium based on the spatial wavenumber domain, including:

Obtain the 3D target area including the abnormal body, and establish the initial prism model including the 3D target area;

The initial prism model is divided into equal intervals along the x , y , and z directions, and a series of nodes are obtained in the x , y , and z directions; the susceptibility is assigned to each node according to the magnetic susceptibility distribution data. Scalar, represented by χ , in which the nodes in the abnormal body part are assigned to each node according to the magnetic susceptibility value of the abnormal body; the magnetic susceptibility of the nodes without abnormal parts is 0;

Calculate the discrete offset wavenumbers in the x and y directions according to the prism model and the given Gaussian parameters in the x and y directions;

According to the earth's main magnetic field model, calculate the earth's main magnetic field intensity at each node, and use it as the magnetic field intensity of the background field in the space domain;

According to the magnetic field strength of the background field in the space domain and the magnetic field strength of the abnormal field in the space domain, the calculation model of the magnetization is obtained;

Using the two-dimensional Fourier transform to transform the three-dimensional Laplace equation satisfying the magnetic potential and magnetization of the anomalous field in the space domain into a one-dimensional ordinary differential equation in the space wavenumber domain;

Based on the one-dimensional ordinary differential equation in the space wavenumber domain, combined with the set boundary conditions that the magnetic potential of the anomalous field in the space wavenumber domain needs to satisfy, the boundary value problem model satisfied by the magnetic potential of the anomalous field in the space wavenumber domain is transformed into an equivalent variational problem Model;

By solving the variational problem model, the magnetic potential of the anomalous field in the spatial wavenumber domain is obtained;

Based on the magnetic potential of the anomalous field in the spatial wavenumber domain, the magnetic field strength of the anomalous field in the spatial wavenumber domain is obtained, and the magnetic potential of the anomalous field in the spatial wavenumber domain and the magnetic potential in the spatial wavenumber domain are converted into the magnetic potential of the anomalous field in the spatial domain by inverse Fourier transform. , the magnetic field strength of the anomalous field in the space domain;

The total magnetic field strength in the space domain is the sum of the magnetic field strength of the anomalous field in the space domain and the magnetic field strength in the background field in the space domain. The compact operator is used to tighten the total magnetic field strength in the space domain to obtain the tightened total magnetic field strength in the space domain;

Judging whether the current iterative convergence conditions are met, if so, output the spatial domain anomalous field magnetic potential and the spatial domain anomalous magnetic field intensity corresponding to the currently calculated tightened spatial domain total magnetic field strength;

According to the output magnetic field intensity of the abnormal field in the space domain, the magnetic induction intensity in the space domain is obtained.

In another aspect, the present invention provides a device for calculating the magnetic induction intensity of a ferromagnetic medium based on the spatial wavenumber domain, comprising:

The first module is used to obtain the three-dimensional target area including the abnormal body, and establish an initial prism model including the three-dimensional target area;

The second module is used to divide the initial prism model at equal intervals along the x , y , and z directions, and obtain a series of nodes along the x , y and z directions; magnetize each node according to the magnetic susceptibility distribution data rate assignment, the magnetic susceptibility is a scalar, represented by χ , in which the nodes in the abnormal body part are assigned to each node according to the magnetic susceptibility value of the abnormal body; the magnetic susceptibility of the nodes without abnormal parts is 0;

The third module is used to calculate the discrete offset wavenumber in the x and y directions according to the prism model and the given Gaussian parameters in the x and y directions;

The fourth module is used to calculate the strength of the earth's main magnetic field at each node according to the earth's main magnetic field model, and use it as the magnetic field strength of the background field in the space domain;

The fifth module is used to obtain the calculation model of the magnetization according to the magnetic field strength of the background field in the space domain and the magnetic field strength of the abnormal field in the space domain;

The sixth module is used to convert the three-dimensional Laplace equation satisfying the magnetic potential and magnetization of the anomalous field in the space domain into a one-dimensional ordinary differential equation in the space wavenumber domain by using the two-dimensional Fourier transform;

The seventh module is used to transform the boundary value problem model satisfied by the magnetic potential of the anomalous field in the spatial wavenumber domain into Equivalent variational problem models;

The eighth module is used to obtain the magnetic potential of the anomalous field in the spatial wavenumber domain by solving the variational problem model;

The ninth module is used to obtain the magnetic field intensity of the abnormal field in the spatial wavenumber domain based on the magnetic potential of the abnormal field in the spatial wavenumber domain, and convert the magnetic potential of the abnormal field in the spatial wavenumber domain and the magnetic field intensity of the abnormal field in the spatial wavenumber domain into The magnetic potential of the anomalous field in the space domain and the magnetic field strength of the anomalous field in the space domain;

The tenth module is used to tighten the total magnetic field strength of the space domain by using a compact operator to obtain the total magnetic field strength of the space domain after tightening. Sum;

The eleventh module is used to judge whether the current iterative convergence condition is met, and if so, output the magnetic potential of the abnormal field in the space domain and the magnetic field intensity of the abnormal field in the space domain corresponding to the total magnetic field intensity in the space domain after the tightening currently calculated;

The twelfth module is used to obtain the magnetic induction intensity in the space domain by solving according to the output magnetic field intensity of the abnormal field in the space domain.

On the other hand, the present invention provides a computer device, comprising a memory and a processor, the memory stores a computer program, and the processor implements the following steps when executing the computer program:

Obtain the 3D target area including the abnormal body, and establish the initial prism model including the 3D target area;

The initial prism model is divided into equal intervals along the x , y , and z directions, and a series of nodes are obtained in the x , y , and z directions; the susceptibility is assigned to each node according to the magnetic susceptibility distribution data. Scalar, represented by χ , in which the nodes in the abnormal body part are assigned to each node according to the magnetic susceptibility value of the abnormal body; the magnetic susceptibility of the nodes without abnormal parts is 0;

Calculate the discrete offset wavenumbers in the x and y directions according to the prism model and the given Gaussian parameters in the x and y directions;

According to the earth's main magnetic field model, calculate the earth's main magnetic field intensity at each node, and use it as the magnetic field intensity of the background field in the space domain;

According to the magnetic field strength of the background field in the space domain and the magnetic field strength of the abnormal field in the space domain, the calculation model of the magnetization is obtained;

Using the two-dimensional Fourier transform to transform the three-dimensional Laplace equation satisfying the magnetic potential and magnetization of the anomalous field in the space domain into a one-dimensional ordinary differential equation in the space wavenumber domain;

Based on the one-dimensional ordinary differential equation in the space wavenumber domain, combined with the set boundary conditions that the magnetic potential of the anomalous field in the space wavenumber domain needs to satisfy, the boundary value problem model satisfied by the magnetic potential of the anomalous field in the space wavenumber domain is transformed into an equivalent variational problem Model;

By solving the variational problem model, the magnetic potential of the anomalous field in the spatial wavenumber domain is obtained;

Based on the magnetic potential of the anomalous field in the spatial wavenumber domain, the magnetic field strength of the anomalous field in the spatial wavenumber domain is obtained, and the magnetic potential of the anomalous field in the spatial wavenumber domain and the magnetic potential in the spatial wavenumber domain are converted into the magnetic potential of the anomalous field in the spatial domain by inverse Fourier transform. , the magnetic field strength of the anomalous field in the space domain;

The total magnetic field strength in the space domain is the sum of the magnetic field strength of the anomalous field in the space domain and the magnetic field strength in the background field in the space domain. The compact operator is used to tighten the total magnetic field strength in the space domain to obtain the tightened total magnetic field strength in the space domain;

Judging whether the current iterative convergence conditions are met, if so, output the spatial domain anomalous field magnetic potential and the spatial domain anomalous magnetic field intensity corresponding to the currently calculated tightened spatial domain total magnetic field strength;

According to the output magnetic field intensity of the abnormal field in the space domain, the magnetic induction intensity in the space domain is obtained.

In another aspect, the present invention also provides a computer-readable storage medium on which a computer program is stored, and when the computer program is executed by a processor, the following steps are implemented:

Obtain the 3D target area including the abnormal body, and establish the initial prism model including the 3D target area;

The initial prism model is divided into equal intervals along the x , y , and z directions, and a series of nodes are obtained in the x , y , and z directions; the susceptibility is assigned to each node according to the magnetic susceptibility distribution data. Scalar, represented by χ , in which the nodes in the abnormal body part are assigned to each node according to the magnetic susceptibility value of the abnormal body; the magnetic susceptibility of the nodes without abnormal parts is 0;

Calculate the discrete offset wavenumbers in the x and y directions according to the prism model and the given Gaussian parameters in the x and y directions;

According to the earth's main magnetic field model, calculate the earth's main magnetic field intensity at each node, and use it as the magnetic field intensity of the background field in the space domain;

According to the magnetic field strength of the background field in the space domain and the magnetic field strength of the abnormal field in the space domain, the calculation model of the magnetization is obtained;

Using the two-dimensional Fourier transform to transform the three-dimensional Laplace equation satisfying the magnetic potential and magnetization of the anomalous field in the space domain into a one-dimensional ordinary differential equation in the space wavenumber domain;

Based on the one-dimensional ordinary differential equation in the space wavenumber domain, combined with the set boundary conditions that the magnetic potential of the anomalous field in the space wavenumber domain needs to satisfy, the boundary value problem model satisfied by the magnetic potential of the anomalous field in the space wavenumber domain is transformed into an equivalent variational problem Model;

By solving the variational problem model, the magnetic potential of the anomalous field in the spatial wavenumber domain is obtained;

Based on the magnetic potential of the anomalous field in the spatial wavenumber domain, the magnetic field strength of the anomalous field in the spatial wavenumber domain is obtained, and the magnetic potential of the anomalous field in the spatial wavenumber domain and the magnetic potential in the spatial wavenumber domain are converted into the magnetic potential of the anomalous field in the spatial domain by inverse Fourier transform. , the magnetic field strength of the anomalous field in the space domain;

The total magnetic field strength in the space domain is the sum of the magnetic field strength of the anomalous field in the space domain and the magnetic field strength in the background field in the space domain. The compact operator is used to tighten the total magnetic field strength in the space domain to obtain the tightened total magnetic field strength in the space domain;

Judging whether the current iterative convergence conditions are met, if so, output the spatial domain anomalous field magnetic potential and the spatial domain anomalous magnetic field intensity corresponding to the currently calculated tightened spatial domain total magnetic field strength;

According to the output magnetic field intensity of the abnormal field in the space domain, the magnetic induction intensity in the space domain is obtained.

Compared with the prior art, the advantages of the present invention are:

1. Considering the self-demagnetization effect, the present invention can more accurately perform numerical simulation of the magnetic field of the ferromagnetic medium.

2. For the case where there are unknown parameters in the right end of the model, an iterative method is proposed to solve the problem, which provides a solution idea for similar problems.

3. Further, in order to reduce the truncation effect of the Fourier transform, edge expansion processing is adopted.

4. A compact operator is added to make the algorithm converge stably.

5. The present invention reduces the three-dimensional problem to one-dimensional through Fourier transform, improves the calculation accuracy and calculation efficiency, and has good algorithm parallelism.

6. The present invention is of great significance for the estimation of iron ore deposits, the inversion and interpretation of measured data.

Description of drawings

Fig. 1 is a flow chart in an embodiment of the present invention;

2 is a schematic diagram of an initial prism model in an embodiment of the present invention;

Fig. 3 is the schematic diagram that the initial prism model shown in Fig. 2 is divided;

4 is a schematic diagram of a target area and an abnormal body in an embodiment of the present invention;

5 is a schematic diagram of the numerical solution, analytical solution and absolute error of the isotropic spherical model in an embodiment of the present invention; wherein (a), ( b ), (c ) represent the numerical solution of _Bax , _respectively , and the The absolute error of the analytical solution and the _numerical solution and analytical solution of Bax ; (d), (e), (f) represent the _numerical solution of Bay , the analytical solution of _Bay and the _absolute error of the numerical solution and analytical solution of Bay respectively Error; (g), (h), (i) represent the numerical solution of B _az , the analytical solution of B _az and the absolute error of the numerical solution and analytical solution of B _az ;

FIG. 6 is an internal structure diagram of a computer device in an embodiment of the present invention.

Detailed ways

In order to make the purposes, technical solutions and advantages of the embodiments of the present invention clearer, the following will clearly illustrate the spirit of the disclosed contents of the present invention with the accompanying drawings and detailed description. , when changes and modifications can be made by the technology taught by the content of the present invention, it does not depart from the spirit and scope of the content of the present invention. The exemplary embodiments of the present invention and their descriptions are used to explain the present invention, but are not intended to limit the present invention.

1 , a method for calculating the magnetic induction intensity of a ferromagnetic medium based on the spatial wavenumber domain provided by the present invention includes:

(S1) Acquire a 3D target area including anomalous bodies, and establish an initial prism model including the 3D target area, as shown in Figure 2, where the initial prism model is a cuboid or a cube, and the anomaly body can be of any shape. The abnormal body in Figure 2 is a cuboid, and the abnormal body is located inside the initial prism model.

(S2) The initial prism model is divided at equal intervals along the x , y , and z directions, respectively, and a series of nodes are obtained in the x , y , and z directions, as shown in Figure 3.

Next, assign the magnetic susceptibility to each node according to the magnetic susceptibility distribution data. The magnetic susceptibility is a scalar, which is represented by χ and the unit is SI.

It can be understood that a person skilled in the art can select an appropriate assignment method according to existing methods. In an embodiment of the present invention, for the nodes of the abnormal body part, each node is assigned according to the magnetic susceptibility value of the abnormal body; the magnetic susceptibility of the nodes without the abnormal body part is 0;

(S3) According to the prism model and the given Gaussian parameters in the x and y directions, calculate the discrete offset wavenumbers in the x and y directions;

(S4) According to the earth's main magnetic field model, calculate the earth's main magnetic field intensity at each node, and use it as the magnetic field intensity of the background field in the space domain;

(S5) According to the magnetic field strength of the background field in the space domain and the magnetic field strength of the abnormal field in the space domain, a calculation model of the magnetization is obtained;

(S6) Using the two-dimensional Fourier transform to convert the three-dimensional Laplace equation satisfying the magnetic potential and magnetization of the anomalous field in the spatial domain into a one-dimensional ordinary differential equation in the spatial wavenumber domain;

(S7) Based on the one-dimensional ordinary differential equation in the spatial wavenumber domain, combined with the set boundary conditions that the magnetic potential of the anomalous field in the spatial wavenumber domain needs to satisfy, transform the boundary value problem model satisfied by the magnetic potential of the anomalous field in the spatial wavenumber domain into an equivalent Variational problem models;

(S8) Obtain the magnetic potential of the anomalous field in the spatial wavenumber domain by solving the variational problem model;

(S9) Based on the magnetic potential of the anomalous field in the spatial wavenumber domain, obtain the magnetic field strength of the anomalous field in the spatial wavenumber domain, and convert the magnetic potential of the anomalous field in the spatial wavenumber domain and the magnetic field strength of the anomalous field in the spatial wavenumber domain into the anomaly in the spatial domain through inverse Fourier transform Field magnetic potential, anomalous magnetic field strength in space domain;

(S10) The total magnetic field strength of the space domain is the sum of the magnetic field strength of the anomalous field of the space domain and the magnetic field strength of the background field of the space domain. The total magnetic field strength of the space domain is tightened by using a compact operator, and the tightened total magnetic field strength of the space domain is obtained;

(S11) Determine whether the current iterative convergence conditions are met, if not, return to (S8), if so, output the spatial domain anomaly field magnetic potential and the spatial domain anomaly corresponding to the currently calculated tightened spatial domain total magnetic field strength Field magnetic field strength;

(S12) Solve according to the output magnetic field intensity of the anomalous field in the space domain to obtain the magnetic induction intensity in the space domain.

It can be understood that in step ( S3 ) of the present invention, Gaussian parameters can be set with reference to the Gaussian parameter setting methods in existing methods in the art, and the discrete offset wavenumbers in the x and y directions can be calculated.

In step (S3) of an embodiment of the present invention, the number of Gaussian points N _x in the x direction is given, the Gaussian point ta and the Gaussian coefficient A _a in the interval [-1,1], where _a = 1,2 ,..., N _x ; given the number of Gaussian points N _y in the y direction, the Gaussian point t _b and the Gaussian coefficient A _b on the interval [-1,1], where b =1,2,..., N _y .

In step (S3) of an embodiment of the present invention, the discrete offset wavenumbers in the x and y directions are calculated by the following methods:

（1）

(1)

（2）

(2)

In the formula,

,

,

表示x方向基波数，NN _x 表示初始棱柱体模型沿x方向等间隔剖分得到的节点数，

表示初始棱柱体模型沿x方向等间隔剖分时采用的单位间隔长度；k _y 表示y方向的偏移波数，

表示y方向基波数，NN _y 表示初始棱柱体模型沿y方向等间隔剖分得到的节点数，

表示初始棱柱体模型沿y方向等间隔剖分时采用的单位间隔长度。where: k _x represents the offset wavenumber in the x direction,

represents the fundamental wave number in the x -direction, NN _x represents the number of nodes obtained by dividing the initial prism model at equal intervals along the x -direction,

represents the unit interval length used when the initial prism model is divided at equal intervals along the x -direction; k _y represents the migration wavenumber in the y -direction,

represents the fundamental wave number in the y direction, NN _y represents the number of nodes obtained by dividing the initial prism model at equal intervals along the y direction,

Indicates the unit interval length used when the initial prism model is equally spaced along the y -direction.

In the step (S4) of an embodiment of the present invention, it includes:

According to the earth's main magnetic field model IGRF, calculate the earth's main magnetic field strength at each node, and take it as the space domain background field magnetic field strength H ₀ at each node, which is the background field in the numerical simulation, that is, the magnetic field when there is no abnormality , the unit is A/m. The components of the magnetic field intensity H ₀ of the background field in the space domain are expressed as H _0x , H _{0 y} , and H _{0 z} , respectively:

（3）

(3)

（4）

(4)

（5）

(5)

表示H₀的L2范数，α为目标区域磁倾角，β为目标区域磁偏角。in:

Represents the L2 norm of H ₀ , α is the magnetic inclination angle of the target area, and β is the magnetic declination angle of the target area.

In the step (S5) of an embodiment of the present invention, the calculation model of the magnetization is:

（6）

(6)

Among them, H _a represents the magnetic field strength generated by the abnormal body at each node, that is, the magnetic field strength of the abnormal field in its spatial domain, which is the abnormal field in the numerical simulation, that is, the magnetic field generated by the magnetic susceptibility abnormality, and the unit is A/m, Its three components are H _ax , _Hay , and H _az , respectively. The total magnetic field intensity H in the space domain at each node is the sum of the magnetic field intensity H ₀ of the background field in the space domain and the magnetic field intensity H _a of the anomalous field in the space domain.

In the step ( S6 ) of an embodiment of the present invention, the three-dimensional Laplace equation satisfied by the magnetic potential U _a and the magnetization M of the anomalous field in the space domain

（7）

(7)

为拉普拉斯算子，展开为

，

表示对磁化强度

求散度。上式（7）展开即为：In the formula,

is the Laplacian operator, which expands to

,

Indicates the magnetization

Find divergence. The expansion of the above formula (7) is:

Perform a two-dimensional Fourier transform on formula (7), in order to reduce the truncation effect of the Fourier transform, the edge expansion process is adopted, that is, the x and y directions are expanded by a certain distance on the basis of the target area to avoid the truncation effect on the field. The impact of value calculations.

Through the two-dimensional Fourier transform, the one-dimensional ordinary differential equation in the spatial wavenumber domain is obtained, and the z direction is reserved as the spatial domain:

（8）

(8)

表示空间波数域异常场磁位，

、

、

为波数域磁化强度，k _x 、k _y 分别为x、y方向的偏移波数。公式（8）即空间波数域异常场磁位

所满足的一维常微分方程。in

represents the magnetic potential of the anomalous field in the spatial wavenumber domain,

,

,

is the magnetization in the wavenumber domain, and k _x and ky are the offset wave numbers in the x and _y directions, respectively . Formula (8) is the magnetic potential of the anomalous field in the space wavenumber domain

The one-dimensional ordinary differential equation satisfied.

In step (S7), in order to obtain a definite solution of the governing equation, appropriate boundary conditions need to be given. In the Cartesian coordinate system, the Z axis is taken vertically downward as the positive direction, the horizontal ground of the target area is taken as the upper boundary Z _min , and the ground is far enough away from the abnormal body as the lower boundary Z _max , and the upper and lower boundary conditions are satisfied:

（9）Upper boundary:

(9)

（10）Lower boundary:

(10)

。in,

.

By combining formulas (8), (9) and (10), the boundary value problem model satisfied by the magnetic potential of the anomalous field in the spatial wavenumber domain is obtained:

（11）

(11)

Using the variational method, the boundary value problem model satisfied by the magnetic potential of the anomalous field in the spatial wavenumber domain is transformed into an equivalent variational problem model:

(12)

(12)

在单元内二次变化。In the Cartesian coordinate system, the unit is divided along the z direction, and the quadratic interpolation function is used in each unit, that is, the magnetic potential of the anomalous field in the spatial wavenumber domain

Quadratic variation within the unit.

所满足的一维常微分方程的右端项中包含有背景场和异常场，而异常场未知，因此采用迭代求解，这也给同类问题提供了求解思路（即模型右端项有未知参数的情况，提出用迭代法进行求解）。Magnetic potential of anomalous field in spatial wavenumber domain

The right-hand term of the satisfied one-dimensional ordinary differential equation includes the background field and the abnormal field, and the abnormal field is unknown, so it is solved by iterative solution, which also provides a solution idea for similar problems (that is, when the right-hand term of the model has unknown parameters, It is proposed to solve it by an iterative method).

Combining equation (6) and equation (7), it can be known that the magnetization M is determined by the product of the sum of the background field H ₀ and the anomalous field H _a and the magnetic susceptibility, and H=H ₀ +H _a , while the anomalous field H _a is unknown , so the first iteration assumes that _Ha is 0, and the sum of the anomalous field and the background field in the first iteration is replaced by the background field, so that the one-dimensional partial differential equation is transformed into a one-dimensional ordinary differential equation to solve, and the first anomalous field is obtained, and then Then, the sum of the obtained anomalous field and background field is used as the total field of the right end term for the next solution.

。The problem to be solved each time is a variational problem. The variational problem is analyzed item by item and the overall synthesis is carried out to obtain a five-diagonal equation composed of all nodes. The chase method can be used to achieve fast and efficient solutions. Spatial Wavenumber Domain Anomaly Field Magnetic Potentials at Nodes

.

，通过下式求得空间波数域异常场磁场强度：In the step (S9) of an embodiment of the present invention, the magnetic potential of the anomalous field based on the spatial wavenumber domain

, the magnetic field strength of the anomalous field in the space wavenumber domain is obtained by the following formula:

（13）

(13)

表示对z求微分，

表示空间波数域异常场磁位对z求微分。where i is an imaginary number,

means to differentiate with respect to z ,

Indicates that the magnetic potential of the anomalous field in the spatial wavenumber domain is differentiated against z .

The magnetic potential of the anomalous field in the spatial wavenumber domain and the magnetic field strength of the anomalous field in the spatial wavenumber domain can be obtained through the inverse Fourier transform to obtain the magnetic potential U _a of the anomalous field in the spatial domain and the strength H _a of the anomalous field in the spatial domain.

In the step (S10) of an embodiment of the present invention, the total magnetic field intensity in the space domain is tightened by using a compact operator to obtain the total magnetic field intensity in the space domain after tightening, wherein the total magnetic field intensity in the space domain is the magnetic field intensity of the abnormal field in the space domain and the sum of the magnetic field strength of the background field in the space domain.

（14）

(14)

分别为第j+1次迭代中的空间域背景场磁场强度和空间域异常场磁场强度，H ^j 表示第j次迭代中计算得到的空间域总磁场强度。Among them, j represents the number of iterations, H ^{j + 1} is the total magnetic field strength in the space domain calculated in the j + 1 iteration,

are the magnetic field strength of the background field in the space domain and the magnetic field strength of the abnormal field in the space domain in the j +1th iteration, respectively, and H ^j represents the total magnetic field strength in the space domain calculated in the jth iteration.

It can be understood that the preset iteration termination conditions refer to preset model calculation constraints, which are used to constrain the performance calculation process of the entire model to converge, so that the model can output results satisfying the conditions. In the present invention, the iteration termination condition in step (S11) can be set as:

（15）

(15)

表示第j次迭代中计算得到的收紧后的空间域总磁场强度，

表示第j+1次迭代中计算得到的收紧后的空间域总磁场强度。The iteration stops when the above iteration convergence conditions are met, where

represents the total magnetic field strength in the space domain after tightening calculated in the jth iteration,

Represents the total magnetic field strength in the space domain after tightening calculated in the j +1th iteration.

Of course, in practical applications, those skilled in the art can also set other iteration termination conditions based on the existing technology, conventional technical means in the field or common knowledge, which is not limited to the iteration set in the above preferred embodiments of this application. Termination condition.

Finally, in step (S12), from the relationship between the magnetic induction intensity _Ba (unit is T) of the abnormal field and the magnetic field intensity H _a of the abnormal field in the space domain, the magnetic induction intensity _Ba can be obtained, and then the three components Bax of _Ba can be _obtained , _Bay , B az _.

（16）

(16)

Among them, μ is the absolute permeability of the ferromagnetic medium, and the unit is H/m.

。The relationship between absolute permeability μ and χ is shown in equation (17), where μ ₀ is the permeability in vacuum,

.

（17）

(17)

The accuracy and efficiency of the method for calculating the magnetic induction intensity of the ferromagnetic medium based on the spatial wavenumber domain provided by the present invention are examined below.

The test computer is configured as i5-4590, the main frequency is 3.30GHz, and the memory is 12GB.

The target area is a three-dimensional prism structure, its size is 1000m×1000m×1000m, and the x and y directions are expanded by 3000m on both sides. The abnormal body in the target area is an isotropic sphere with a radius of 200m, and the center of the sphere is located at the center of the model (500m, 500m, 500m). The schematic diagram is shown in Figure 4.

The anomaly is an isotropic sphere with a spherical magnetic susceptibility of 2 SI. The background magnetic field strength of the model area is 50000nT, the magnetic inclination angle of the target area is 45°, and the magnetic declination angle of the target area is 9°. The number of nodes in the x , y , and z directions is 101, 101, and 101. In the case of a single thread, the edge-expanding FFT is adopted, and it takes 7 iterations to reach the convergence condition, which takes 56s and occupies 1.38GB of memory, which is low in memory and high in efficiency. The absolute errors of B _ax , _Bay , B _az and the analytical solution are all different from the field value by more than two orders of magnitude, which meets the requirements of numerical calculation. Figure 5 shows the numerical solutions, analytical solutions and their absolute errors of the B _ax , _Bay , B _az components calculated by the model. (a), ( b ) and ( _c ) _in The absolute error of the analytical solution and the _numerical and analytical solutions of Bax ; (d), (e) and (f) in Figure 5 represent the _numerical solution of Bay , the analytical solution of _Bay , and the numerical and _analytical solutions of Bay , respectively The absolute error of the solution; (g), (h) and (i) in Figure 5 represent the _numerical solution of Baz , the analytical solution of Baz , and the absolute _errors of the numerical and analytical solutions of _Baz , respectively.

On the basis of considering the self-demagnetization effect, the standard edge-expanding FFT is used to reduce the three-dimensional problem to one dimension, only the z direction is retained, and the one-dimensional finite element method is used. Iteratively solves, and uses compact operators to ensure the convergence of the algorithm, which greatly improves the calculation accuracy and calculation efficiency, and the algorithm has good parallelism and occupies less memory. The invention can be used in the inversion algorithm of magnetic survey data to improve the precision and efficiency of inversion.

In this embodiment, a computer device is provided, and the computer device may be a server, and its internal structure diagram may be as shown in FIG. 6 . The computer device includes a processor, memory, a network interface, and a database connected by a system bus. Among them, the processor of the computer device is used to provide computing and control capabilities. The memory of the computer device includes a non-volatile storage medium, an internal memory. The nonvolatile storage medium stores an operating system, a computer program, and a database. The internal memory provides an environment for the execution of the operating system and computer programs in the non-volatile storage medium. The computer device's database is used to store sample data. The network interface of the computer device is used to communicate with an external terminal through a network connection. When the computer program is executed by the processor, the above-mentioned method for calculating the magnetic induction intensity of a ferromagnetic medium based on the space wavenumber domain is realized.

Those skilled in the art can understand that the structure shown in FIG. 6 is only a block diagram of a partial structure related to the solution of the present application, and does not constitute a limitation on the computer equipment to which the solution of the present application is applied. Include more or fewer components than shown in the figures, or combine certain components, or have a different arrangement of components.

In one embodiment, a computer device is provided, including a memory and a processor, the memory stores a computer program, and when the processor executes the computer program, the magnetic induction intensity calculation of the ferromagnetic medium based on the spatial wavenumber domain in the above-mentioned embodiment is realized steps of the method.

In one embodiment, a computer-readable storage medium is provided, on which a computer program is stored, and when the computer program is executed by a processor, the steps of the method for calculating the magnetic induction intensity of a ferromagnetic medium based on the spatial wavenumber domain in the foregoing embodiment are implemented .

Those of ordinary skill in the art can understand that all or part of the processes in the methods of the above embodiments can be implemented by instructing relevant hardware through a computer program, and the computer program can be stored in a non-volatile computer-readable storage In the medium, when the computer program is executed, it may include the processes of the above-mentioned method embodiments. Wherein, any reference to memory, storage, database or other medium used in the various embodiments provided in this application may include non-volatile and/or volatile memory. Nonvolatile memory may include read only memory (ROM), programmable ROM (PROM), electrically programmable ROM (EPROM), electrically erasable programmable ROM (EEPROM), or flash memory. Volatile memory may include random access memory (RAM) or external cache memory. By way of illustration and not limitation, RAM is available in various forms such as static RAM (SRAM), dynamic RAM (DRAM), synchronous DRAM (SDRAM), double data rate SDRAM (DDRSDRAM), enhanced SDRAM (ESDRAM), synchronous chain Road (Synchlink) DRAM (SLDRAM), memory bus (Rambus) direct RAM (RDRAM), direct memory bus dynamic RAM (DRDRAM), and memory bus dynamic RAM (RDRAM) and so on.

The technical features of the above embodiments can be combined arbitrarily. In order to make the description simple, all possible combinations of the technical features in the above embodiments are not described. However, as long as there is no contradiction in the combination of these technical features It is considered to be the range described in this specification.

The above-mentioned embodiments only represent several embodiments of the present application, and the descriptions thereof are specific and detailed, but should not be construed as a limitation on the scope of the invention patent. It should be pointed out that for those skilled in the art, without departing from the concept of the present application, several modifications and improvements can be made, which all belong to the protection scope of the present application. Therefore, the scope of protection of the patent of the present application shall be subject to the appended claims.

Claims (10)

1. the magnetic induction intensity calculation method of the ferromagnetic medium based on the space wavenumber domain, is characterized in that, comprises: Obtain the 3D target area including the abnormal body, and establish the initial prism model including the 3D target area; The initial prism model is divided into equal intervals along the x , y , and z directions, and a series of nodes are obtained in the x , y , and z directions; the susceptibility is assigned to each node according to the magnetic susceptibility distribution data. scalar, with

represents that the nodes in the abnormal body part are assigned to each node according to the magnetic susceptibility value of the abnormal body; the magnetic susceptibility of the nodes without abnormal parts is 0; Calculate the discrete offset wavenumbers in the x and y directions according to the prism model and the given Gaussian parameters in the x and y directions; According to the earth's main magnetic field model, calculate the earth's main magnetic field intensity at each node, and use it as the magnetic field intensity of the background field in the space domain; According to the magnetic field strength of the background field in the space domain and the magnetic field strength of the abnormal field in the space domain, the calculation model of the magnetization is obtained; Using the two-dimensional Fourier transform to transform the three-dimensional Laplace equation satisfying the magnetic potential and magnetization of the anomalous field in the space domain into a one-dimensional ordinary differential equation in the space wavenumber domain; Based on the one-dimensional ordinary differential equation in the space wavenumber domain, combined with the set boundary conditions that the magnetic potential of the anomalous field in the space wavenumber domain needs to satisfy, the boundary value problem model satisfied by the magnetic potential of the anomalous field in the space wavenumber domain is transformed into an equivalent variational problem Model; By solving the variational problem model, the magnetic potential of the anomalous field in the spatial wavenumber domain is obtained; Based on the magnetic potential of the anomalous field in the spatial wavenumber domain, the magnetic field strength of the anomalous field in the spatial wavenumber domain is obtained, and the magnetic potential of the anomalous field in the spatial wavenumber domain and the magnetic potential in the spatial wavenumber domain are converted into the magnetic potential of the anomalous field in the spatial domain by inverse Fourier transform. , the magnetic field strength of the anomalous field in the space domain; The total magnetic field strength in the space domain is the sum of the magnetic field strength of the anomalous field in the space domain and the magnetic field strength in the background field in the space domain. The compact operator is used to tighten the total magnetic field strength in the space domain to obtain the tightened total magnetic field strength in the space domain; Judging whether the current iterative convergence conditions are met, if so, output the spatial domain anomalous field magnetic potential and the spatial domain anomalous magnetic field intensity corresponding to the currently calculated tightened spatial domain total magnetic field strength; According to the output magnetic field intensity of the abnormal field in the space domain, the magnetic induction intensity in the space domain is obtained. 2. The method for calculating the magnetic induction intensity of a ferromagnetic medium based on a spatial wavenumber domain according to claim 1, wherein the Gaussian parameter comprises: The number of Gaussian points in the x direction N _x , the Gaussian point _ta on the interval [-1,1], the Gaussian coefficient A _a , where _a = 1,2,..., N x ; The number of Gaussian points in the y direction N _y , the Gauss point t _b on the interval [-1,1], and the Gaussian coefficient A _b , where b =1,2,..., N _y . 3. the magnetic induction intensity calculation method of the ferromagnetic medium based on the space wavenumber domain according to claim 2, is characterized in that, x and y direction discrete shift wavenumber, is calculated by the following method:

In the formula,

,

where: k _x represents the offset wavenumber in the x direction,

represents the fundamental wave number in the x -direction, NN _x represents the number of nodes obtained by dividing the initial prism model at equal intervals along the x -direction,

represents the unit interval length used when the initial prism model is divided at equal intervals along the x -direction; k _y represents the migration wavenumber in the y -direction,

represents the fundamental wave number in the y direction, NN _y represents the number of nodes obtained by dividing the initial prism model at equal intervals along the y direction,

Indicates the unit interval length used when the initial prism model is equally spaced along the y -direction. 4. The method for calculating the magnetic induction intensity of a ferromagnetic medium based on the spatial wavenumber domain according to claim 1, 2 or 3, wherein the calculation model of the magnetization is:

Among them, H _a represents the magnetic field intensity generated by the abnormal body at each node, that is, the magnetic field intensity of the abnormal field in its spatial domain, H ₀ represents the magnetic field intensity of the space domain background field at each node, and H represents the total spatial domain magnetic field intensity at each node. Magnetic field strength. 5. The method for calculating the magnetic induction intensity of the ferromagnetic medium based on the spatial wavenumber domain according to claim 4, wherein the three-dimensional Laplace equation is satisfied for the magnetic potential U _a and the magnetization M of the abnormal field in the spatial domain

Perform a two-dimensional Fourier transform to obtain a one-dimensional ordinary differential equation in the spatial wavenumber domain:

where i is an imaginary number,

is the Laplace operator,

Indicates the magnetization

seek divergence;

represents the magnetic potential of the anomalous field in the spatial wavenumber domain,

,

,

are the magnetization in the wavenumber domain in the x , y , and z directions, respectively, k _x and ky are the offset wavenumbers in the x and _y directions, respectively ,

. 6. the magnetic induction intensity calculation method of the ferromagnetic medium based on the space wavenumber domain according to claim 5, is characterized in that, the boundary condition that the space wavenumber domain abnormal field magnetic potential needs to satisfy is: Upper boundary:

Lower boundary:

. 7. The method for calculating the magnetic induction intensity of a ferromagnetic medium based on the space wavenumber domain according to claim 6, wherein the one-dimensional ordinary differential equation in the space wavenumber domain and the boundary condition that the magnetic potential of the space wavenumber domain anomaly field needs to satisfy are composed together. The boundary value problem model of the magnetic potential satisfaction of the anomalous field in the space wavenumber domain is transformed into an equivalent variational problem model by using the variational method. 8. The method for calculating the magnetic induction intensity of a ferromagnetic medium based on the spatial wavenumber domain according to claim 7, wherein the variational problem model is solved by the shape function method to obtain the magnetic potential of the abnormal field in the spatial wavenumber domain

. 9. the method for calculating the magnetic induction intensity of the ferromagnetic medium based on the space wavenumber domain according to claim 8, is characterized in that, utilizes the compact operator to tighten, obtains the space domain total magnetic field intensity after tightening, as follows:

Among them, j represents the number of iterations, H ^{j + 1} is the total magnetic field strength in the space domain calculated in the j + 1 iteration,

are the magnetic field strength of the background field in the space domain and the magnetic field strength of the abnormal field in the space domain in the j +1th iteration, respectively, and H ^j represents the total magnetic field strength in the space domain calculated in the jth iteration. 10. The method for calculating the magnetic induction intensity of a ferromagnetic medium based on the spatial wavenumber domain according to claim 1, 2, 3, 5, 6, 7, 8 or 9, wherein the iteration termination condition is set as:

in

represents the total magnetic field strength in the space domain after tightening calculated in the jth iteration,

Represents the total magnetic field strength in the space domain after tightening calculated in the j +1th iteration.

Images