CN113627027B - Method and system for simulating electromagnetic field value of non-trivial anisotropic medium - Google Patents

Abstract

The invention provides a numerical simulation method and a numerical simulation system for an electromagnetic field of a non-trivial anisotropic medium, which directly correct a Maxwell electromagnetic field control equation through a current density divergence constraint term to be more perfect, avoid the solution of extra unknown quantity or equation in the solution process and reduce the calculation burden. The method provided by the embodiment of the invention does not need to correct the current density divergence, can avoid the correction failure caused by error correction fields and algorithm restart after correction when the current density divergence correction is carried out on a non-trivial anisotropic medium in the prior art, thereby effectively suppressing the electromagnetic field pseudo solution in low-frequency application, obviously improving the convergence speed and the calculation efficiency of a solver, greatly shortening the solving time, ensuring the accuracy of numerical solution, and having important significance on the problems of electromagnetic response calculation of a complex medium, fine inversion imaging of an underground geological structure and the like.

Description

Numerical simulation method and system for electromagnetic field of non-trivial anisotropic medium

Technical Field

The present invention relates to the technical field of electromagnetic field numerical simulation, and in particular to a method and system for numerical simulation of electromagnetic field of non-trivial anisotropic medium.

Background Art

Spurious solutions refer to numerical solutions that have no physical meaning. They are common in electromagnetic numerical simulation problems and are usually closely related to the numerical methods and control equations used. It is generally believed that spurious solutions come from inappropriate approximations of the null space of the curl operator in Maxwell's equations or the eigenvectors in the double curl equations, but are usually caused by multiple factors. Electromagnetic spurious solutions can easily lead to non-uniqueness problems and produce invalid solutions, resulting in slow convergence of numerical algorithms and inaccurate numerical solutions. Therefore, in the face of different problems, predecessors have proposed a variety of technical means to try to suppress spurious solutions, such as staggered grids, finite element functions with continuous derivatives, solving magnetic fields to obtain electric fields or solving electric fields to obtain magnetic fields, current density divergence correction, and current density divergence constraints.

Current density divergence correction is widely used in low-frequency electromagnetic problems, especially in numerical simulation of geophysical magnetotelluric sounding (MT) and related methods. When the frequency is low, the conductivity term in the dual curl equation can be almost ignored; since the gradient space is the null space of the curl operator, the sum of the gradient of any scalar and the true solution is the solution of the equation, resulting in a pseudo solution. By introducing the current density divergence correction technology, in addition to the dual curl main control equation, an additional current continuity equation is solved with the current redundant divergence corresponding to the iterative electric field solution of the main control equation as the source, and the iterative electric field solution is corrected with the additional solution obtained, so that the pseudo solution is suppressed. This method has achieved good results in the low-frequency electromagnetic field simulation of conventional isotropic media and trivial (axial) anisotropic media, and significantly improved the convergence efficiency.

However, the current density divergence correction technique fails when dealing with non-trivial anisotropic media (i.e., media other than axial anisotropy) because: the current density divergence correction technique solves the additional correction field through the potential, and the potential gradient component is always distributed along the three coordinate axes of the rectangular coordinate system, which cannot be used to describe the electric field relationship corresponding to each element in the conductivity tensor, that is, the potential gradient at this time is an equivalent value to satisfy the current density divergence condition; however, when it is used to correct the iterative electric field solution and is multiplied by the conductivity tensor again in the main control equation, errors occur because the conductivity tensor at this time is the original one, not the equivalent diagonal conductivity tensor corresponding to the additional correction field. For this reason, it is urgent to provide a numerical simulation method for the electromagnetic field of non-trivial anisotropic media.

Summary of the invention

The present invention provides a numerical simulation method and system for electromagnetic fields of non-trivial anisotropic media, which are used to solve the defects in the prior art.

The present invention provides a method for numerical simulation of electromagnetic fields of non-trivial anisotropic media, comprising:

Determine the spatial topological locations of current density divergence constraints and conductivity weighting factors in nontrivial anisotropic media;

Determining the current density divergence constraint term based on the spatial topological position and the conductivity weighting factor;

Based on the current density divergence constraint term, the Maxwell electromagnetic field control equation is corrected to obtain a corrected equation;

The modified equation is solved to obtain a numerical solution of the electromagnetic field in the non-trivial anisotropic medium.

According to a numerical simulation method for electromagnetic fields of non-trivial anisotropic media provided by the present invention, the spatial topological position is determined based on the following method:

Construct a structured hexahedral staggered grid;

For the electromagnetic field in the non-trivial anisotropic medium, the electric field vector of the electromagnetic field is defined on the edges of the structured hexahedral staggered grid, the magnetic field vector of the electromagnetic field is defined on the surface of the structured hexahedral staggered grid, and the current density divergence of the electromagnetic field is defined on the nodes of the structured hexahedral staggered grid;

Based on the location of the current density divergence, a spatial topological location of the current density divergence constraint is determined.

According to a numerical simulation method for electromagnetic field of a non-trivial anisotropic medium provided by the present invention, the conductivity weighting factor is determined based on the following method:

For any node of the structured hexahedral staggered grid, averaging the diagonal elements in the grid conductivity tensor around the any node to the any node to obtain the diagonal elements of the conductivity tensor corresponding to the any node;

The root mean square of the diagonal elements of the conductivity tensor is calculated, and the inverse of the root mean square is used as the conductivity weighting factor.

According to a numerical simulation method for electromagnetic fields of non-trivial anisotropic media provided by the present invention, the current density divergence constraint term is determined based on the spatial topological position and the conductivity weighting factor, specifically including:

Determining an alternative current density divergence constraint term based on the conductivity weighting factor and the spatial topological position;

The current density divergence constraint term is determined based on the candidate current density divergence constraint term.

According to a numerical simulation method for electromagnetic fields of non-trivial anisotropic media provided by the present invention, the current density divergence constraint term is determined based on the alternative current density divergence constraint term, specifically including:

The candidate current density divergence constraint term is regularized to obtain the current density divergence constraint term.

According to a numerical simulation method for electromagnetic fields in non-trivial anisotropic media provided by the present invention, solving the correction equation to obtain a numerical solution of the electromagnetic field in the non-trivial anisotropic media specifically includes:

Based on the structured hexahedral staggered grid, the correction equation and the calculation domain corresponding to the correction equation are discretized by using the finite volume method to obtain a first discrete result corresponding to the correction equation and a second discrete result corresponding to the calculation domain;

Determining boundary conditions of the second discrete result, and determining a sparse linear equation group of the electromagnetic field based on the first discrete result and the boundary conditions;

The sparse linear equations are solved based on a Krylov subspace solver to obtain a numerical solution to the electromagnetic field.

According to a numerical simulation method for electromagnetic fields of non-trivial anisotropic media provided by the present invention, the sparse linear equations are solved based on a Krylov subspace solver to obtain the numerical solution of the electromagnetic field, which specifically includes:

Based on the Krylov subspace solver, the quasi-minimum residual method is adopted to iteratively solve the sparse linear equations to obtain the numerical solution of the electromagnetic field.

The present invention also provides a non-trivial anisotropic medium electromagnetic field numerical simulation system, comprising:

A first determination module is used to determine the spatial topological position of the current density divergence constraint and the conductivity weighting factor in the non-trivial anisotropic medium;

A second determination module, configured to determine the current density divergence constraint item based on the spatial topological position and the conductivity weighting factor;

A correction module, used for correcting the Maxwell electromagnetic field control equation based on the current density divergence constraint term to obtain a corrected equation;

A solution module is used to solve the modified equation to obtain a numerical solution of the electromagnetic field in the non-trivial anisotropic medium.

The present invention also provides an electronic device, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein when the processor executes the computer program, the steps of the numerical simulation method of the electromagnetic field of a non-trivial anisotropic medium as described above are implemented.

The present invention also provides a non-transitory computer-readable storage medium having a computer program stored thereon, wherein when the computer program is executed by a processor, the steps of any of the above-mentioned methods for numerically simulating electromagnetic fields of non-trivial anisotropic media are implemented.

The numerical simulation method and system of electromagnetic field of non-trivial anisotropic medium provided by the present invention directly corrects the Maxwell electromagnetic field control equation through the current density divergence constraint term, making it more perfect, avoiding the solution of additional unknown quantities or equations in the solution process, and reducing the calculation burden. The method provided in the embodiment of the present invention does not need to correct the current density divergence, which can avoid the correction failure caused by the wrong correction field and the restart of the algorithm after correction when correcting the current density divergence of non-trivial anisotropic media in the prior art, thereby effectively suppressing the electromagnetic field pseudo-solution in low-frequency applications, significantly improving the convergence speed and calculation efficiency of the solver, greatly shortening the solution time, and ensuring the accuracy of the numerical solution, which is of great significance to the electromagnetic response calculation of complex media and the fine inversion imaging of underground geological structures. In addition, the method provided in the embodiment of the present invention is relatively simpler and more direct to implement, which improves practicality.

BRIEF DESCRIPTION OF THE DRAWINGS

In order to more clearly illustrate the technical solutions in the present invention or the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

FIG1 is a schematic diagram of a flow chart of a numerical simulation method for electromagnetic fields of non-trivial anisotropic media provided by the present invention;

FIG2 is a finite volume three-dimensional grid partition diagram of the numerical simulation model of the electromagnetic field of a non-trivial anisotropic medium provided by the present invention;

FIG3 is a schematic diagram of the spatial topological position of the current density divergence constraint provided by the present invention;

FIG4 is a schematic diagram of electrical parameters of one-dimensional layered models of different types of anisotropic media provided by the present invention;

FIG5 is a schematic diagram of a three-dimensional COMMEMI 3D-2 model of different types of anisotropic media provided by the present invention;

FIG6 is a comparison diagram of the efficiency test results of the half-space model, the one-dimensional layered model, and the three-dimensional COMMEMI 3D-2 model of an arbitrary anisotropic medium provided by the present invention;

FIG7 is a comparison diagram of the accuracy test of the one-dimensional layered model response of any anisotropic medium provided by the present invention;

FIG8 is a schematic diagram of an asymmetric vein model of an arbitrary anisotropic medium provided by the present invention;

FIG9 is a comparison diagram of the accuracy test of the asymmetric vein model response of any anisotropic medium provided by the present invention;

FIG10 is a comparison diagram of the accuracy test of the three-dimensional COMMEMI 3D-2 model response of any anisotropic medium provided by the present invention;

FIG11 is a schematic diagram of the structure of a numerical simulation system for electromagnetic fields of non-trivial anisotropic media provided by the present invention;

FIG. 12 is a schematic diagram of the structure of an electronic device provided by the present invention.

DETAILED DESCRIPTION

In order to make the purpose, technical solution and advantages of the present invention clearer, the technical solution of the present invention will be clearly and completely described below in conjunction with the drawings of the present invention. Obviously, the described embodiments are part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

Since the numerical simulation of electromagnetic fields of non-trivial anisotropic media is crucial for the calculation of electromagnetic response of complex media and the fine inversion imaging of underground geological structures, it is necessary to study new technical means to deal with such problems. To this end, an embodiment of the present invention provides a method for numerical simulation of electromagnetic fields of non-trivial anisotropic media.

FIG1 is a flow chart of a method for numerically simulating electromagnetic fields of a non-trivial anisotropic medium provided in an embodiment of the present invention. As shown in FIG1 , the method includes:

S1, determine the spatial topological location of the current density divergence constraint and the conductivity weighting factor in nontrivial anisotropic media;

S2, determining the current density divergence constraint term based on the spatial topological position and the conductivity weighting factor;

S3, based on the current density divergence constraint term, modifying the Maxwell electromagnetic field control equation to obtain a modified equation;

S4, solving the modified equation to obtain a numerical solution of the electromagnetic field in the non-trivial anisotropic medium.

Specifically, the numerical simulation method of electromagnetic field of non-trivial anisotropic medium provided in the embodiment of the present invention is executed by a server, which can be a local server or a cloud server. The local server can be a computer, etc., which is not specifically limited in the embodiment of the present invention. In the embodiment of the present invention, the object is the numerical simulation method of electromagnetic field of non-trivial anisotropic medium. The non-trivial anisotropic medium can be a medium with non-trivial anisotropy. The non-trivial anisotropy can include horizontal anisotropy, inclination anisotropy and arbitrary anisotropy, etc., wherein arbitrary anisotropy means that the electrical main axis and the coordinate axis do not coincide. The physical or mechanical properties of non-trivial anisotropic medium, such as absorbance, refractive index, conductivity and tensile strength, will produce differences when measured along different directions.

First, step S1 is performed to determine the spatial topological position of the current density divergence constraint in the non-trivial anisotropic medium and the conductivity weighting factor. The current density divergence constraint refers to the constraint condition on the current density divergence in the non-trivial anisotropic medium introduced in the embodiment of the present invention to avoid the pseudo-solution problem in the numerical simulation process.

The current density divergence constraint can be determined by an electromagnetic field model in a non-trivial anisotropic medium, and the electromagnetic field model is a numerical simulation model of the electromagnetic field of a non-trivial anisotropic medium, which can be a structured hexahedral staggered grid as shown in FIG2. The current density divergence constraint can include node constraints, edge constraints, etc. Among them, the node is the node in the structured hexahedral staggered grid, and the edge is the edge in the structured hexahedral staggered grid.

In Figure 2, a three-dimensional rectangular coordinate system is established with the right as the positive direction of the x-axis, the downward as the positive direction of the z-axis, and the lower left as the positive direction of the y-axis. The coordinates of the central node of the structured hexahedral staggered grid are (i, j, k), the coordinates of the lower left node are (i-1, j+1, k+1), the coordinates of the lower right node are (i+1, j+1, k+1), and the coordinates of the central node of the lower right edge are (i+1, j, k+1). Each edge in the structured hexahedral staggered grid has a corresponding electric field vector and the direction of the electric field vector.

Current density divergence refers to the divergence of current density, and its value is a scalar. Current density refers to the product of the conductivity tensor and the electric field vector, which can be expressed as:

J＝σE

Where J is the current density, σ is the conductivity tensor, and E is the electric field vector.

The spatial topological position of the current density divergence constraint in the non-trivial anisotropic medium can be characterized by appropriate nodes or edges in the structured hexahedral staggered grid, which is not specifically limited in the embodiments of the present invention.

The conductivity of a non-trivial anisotropic medium is a tensor σ, and the conductivity tensor has non-diagonal elements, that is, it needs to be characterized by at least four electrical parameters. Since the constraint object in the embodiment of the present invention can be the current density divergence defined on the node, it is necessary to estimate the conductivity weighting factor λ from the multi-element conductivity tensor in the grid around the node. The value of the conductivity weighting factor λ is a scalar.

Then, step S2 is performed to determine the current density divergence constraint term according to the spatial topological position and the conductivity weighting factor. The current density divergence constraint term is used to correct the Maxwell electromagnetic field control equation, so it is introduced as an additional term.

Then, step S3 is executed to correct the Maxwell electromagnetic field control equation (which may be referred to as Maxwell equation for short) according to the current density divergence constraint term to obtain a corrected equation.

According to the geophysical magnetotelluric method (MT), the Maxwell electromagnetic field control equation is expressed in the specific form of the second-order partial differential equation of the frequency domain electric field without external field source, without considering displacement current, and with e ^-iωt as the time factor:

Wherein, ω is the angular frequency, μ is the magnetic permeability, and its value can be the vacuum magnetic permeability.

When the Maxwell electromagnetic field control equations are corrected according to the current density divergence constraint term, the current density divergence constraint term can be introduced into the above-mentioned Maxwell electromagnetic field control equations, that is, the corrected equation shown below is obtained, and the corrected equation can suppress the pseudo solution.

Where M is the current density divergence constraint.

Finally, step S4 is executed to solve the modified equation to obtain the numerical solution of the electromagnetic field in the non-trivial anisotropic medium. When solving the modified equation, an appropriate numerical method can be selected to discretize the modified equation and its computational domain, and then a linear equation group is constructed by processing boundary conditions or field sources. Finally, the numerical solution of the electromagnetic field in the non-trivial anisotropic medium can be obtained through the linear equation group.

The numerical simulation method of electromagnetic field of non-trivial anisotropic medium provided in the embodiment of the present invention directly corrects the Maxwell electromagnetic field control equation through the current density divergence constraint term, making it more perfect, avoiding the solution of additional unknown quantities or equations in the solution process, and reducing the computational burden. The method provided in the embodiment of the present invention does not need to correct the current density divergence, which can avoid the correction failure caused by the wrong correction field and the restart of the algorithm after correction when correcting the current density divergence of non-trivial anisotropic media in the prior art, thereby effectively suppressing the electromagnetic field pseudo-solution in low-frequency applications, significantly improving the convergence speed and calculation efficiency of the solver, greatly shortening the solution time, and ensuring the accuracy of the numerical solution, which is of great significance to the electromagnetic response calculation of complex media and the fine inversion imaging of underground geological structures. In addition, the method provided in the embodiment of the present invention is relatively simpler and more direct to implement, which improves practicality.

On the basis of the above-mentioned embodiment, in the numerical simulation method of electromagnetic field of non-trivial anisotropic medium provided in the embodiment of the present invention, the spatial topological position is determined based on the following method:

Construct a structured hexahedral staggered grid;

For the electromagnetic field in the non-trivial anisotropic medium, the electric field vector of the electromagnetic field is defined on the edges of the structured hexahedral staggered grid, the magnetic field vector of the electromagnetic field is defined on the surface of the structured hexahedral staggered grid, and the current density divergence of the electromagnetic field is defined on the nodes of the structured hexahedral staggered grid;

Based on the location of the current density divergence, a spatial topological location of the current density divergence constraint is determined.

Specifically, in an embodiment of the present invention, when determining the spatial topological position, a structured hexahedral staggered grid as shown in Figure 2 can be first constructed. Then, for the electromagnetic field in a non-trivial anisotropic medium, the electric field vector E of the electromagnetic field is defined on the edge of the structured hexahedral staggered grid, the magnetic field vector H of the electromagnetic field is defined on the surface of the structured hexahedral staggered grid, and the current density divergence of the electromagnetic field is defined on the node of the structured hexahedral staggered grid. Finally, according to the position of the current density divergence, that is, the node of the structured hexahedral staggered grid, the spatial topological position of the current density divergence constraint can be determined. In an embodiment of the present invention, a node or edge can be determined as the spatial topological position of the current density divergence constraint based on advantage analysis. For example, the current density divergence constraint can be directly on the node, or the current density divergence can be averaged to the edge before the constraint is applied, which is not specifically limited in the embodiment of the present invention.

In the embodiment of the present invention, the spatial topological position of the current density divergence constraint is determined by using a structured hexahedral staggered grid, so that the spatial topological position can be made more intuitive.

On the basis of the above-mentioned embodiment, in the numerical simulation method of electromagnetic field of non-trivial anisotropic medium provided in the embodiment of the present invention, the conductivity weighting factor is determined based on the following method:

For any node of the structured hexahedral staggered grid, averaging the diagonal elements in the grid conductivity tensor around the any node to the any node to obtain the diagonal elements of the conductivity tensor corresponding to the any node;

The root mean square of the diagonal elements in the conductivity tensor is calculated, and the inverse of the root mean square is used as the conductivity weighting factor.

Specifically, in the embodiment of the present invention, when determining the conductivity weighting factor, since the conductivity of the non-trivial anisotropic medium has non-diagonal elements, that is, it needs to be characterized by at least four electrical parameters, for any node of the structured hexahedral staggered grid, that is, any node, the volume of the diagonal elements in the grid conductivity tensor around the any node can be first averaged to the any node to obtain the diagonal elements of the conductivity tensor corresponding to the any node, which can be expressed as

Then, the root mean square (RMS) method can be used to calculate the root mean square of the diagonal elements of the conductivity tensor, and the inverse of the calculated root mean square is used as the conductivity weighting factor. That is,

Among them, λ is the conductivity weighting factor, and each node of the structured hexahedral staggered grid corresponds to a conductivity weighting factor.

In an embodiment of the present invention, the diagonal elements in the grid conductivity tensor around the node in the structured hexahedral staggered grid are averaged to the node to obtain the diagonal elements of the conductivity tensor corresponding to the node, and the conductivity weighting factor is obtained by calculating the root mean square of the diagonal elements of the conductivity tensor, so that the obtained conductivity weighting factor can be made more accurate.

On the basis of the above-mentioned embodiment, the numerical simulation method of electromagnetic field of non-trivial anisotropic medium provided in the embodiment of the present invention, the current density divergence constraint term is determined based on the spatial topological position and the conductivity weighting factor, specifically including:

Determining an alternative current density divergence constraint term based on the conductivity weighting factor and the spatial topological position;

The current density divergence constraint term is determined based on the candidate current density divergence constraint term.

Specifically, in the embodiment of the present invention, when determining the current density divergence constraint term, the candidate current density divergence constraint term at the spatial topological position can be determined according to the conductivity weighting factor. That is, there is:

Then, the current density divergence constraint term can be determined according to the alternative current density divergence constraint term. In the embodiment of the present invention, the alternative current density divergence constraint term can be directly used as the current density divergence constraint term, or the current density divergence constraint term can be obtained by processing the alternative current density divergence constraint term, which is not specifically limited in the embodiment of the present invention.

In the embodiment of the present invention, since the alternative current density divergence constraint term is equal to 0, it will not affect the Maxwell electromagnetic field control equations when it is used as the final current density divergence constraint term to correct the Maxwell electromagnetic field control equations, but it can avoid the occurrence of pseudo-solution problems, and can improve the efficiency of numerical simulation of electromagnetic fields of non-trivial anisotropic media without introducing new problems, thereby increasing reliability.

On the basis of the above-mentioned embodiment, the numerical simulation method of electromagnetic field of non-trivial anisotropic medium provided in the embodiment of the present invention, the determining of the current density divergence constraint term based on the alternative current density divergence constraint term specifically includes:

The candidate current density divergence constraint term is regularized to obtain the current density divergence constraint term.

Specifically, in the embodiment of the present invention, in order to balance the constraint term and the double curl term in the original Maxwell electromagnetic field control equation and make the sum of the two approximate the vector Laplacian term, thereby compressing the flat non-zero space, when determining the current density divergence constraint term according to the alternative current density divergence constraint term, the alternative current density divergence constraint term can be regularized, that is, the alternative current density divergence constraint term is transformed into:

That is, the current density divergence constraint term is obtained.

Accordingly, there are:

The correction equation can be expressed as:

In the modified equation, the first two terms on the left side of the equal sign compress the non-trivial null space of the solution, thereby suppressing pseudo-solutions.

而棱边约束项的表达式为

二者实现的功能相似。As shown in Figure 3, it is a schematic diagram of the topological position of the current density divergence constraint. Part (a) of Figure 3 represents the node constraint, and part (b) represents the edge constraint. In Figure 3, the expression of the node constraint term is

The expression of the edge constraint is

The functions achieved by the two are similar.

On the basis of the above-mentioned embodiment, the numerical simulation method of the electromagnetic field of a non-trivial anisotropic medium provided in the embodiment of the present invention, wherein the correction equation is solved to obtain the numerical solution of the electromagnetic field in the non-trivial anisotropic medium, specifically comprises:

Based on the structured hexahedral staggered grid, the correction equation and the calculation domain corresponding to the correction equation are discretized by using the finite volume method to obtain a first discrete result corresponding to the correction equation and a second discrete result corresponding to the calculation domain;

Determining boundary conditions of the second discrete result, and determining a sparse linear equation group of the electromagnetic field based on the first discrete result and the boundary conditions;

The sparse linear equations are solved based on a Krylov subspace solver to obtain a numerical solution to the electromagnetic field.

Specifically, in the embodiment of the present invention, when solving the modified equation to obtain the numerical solution of the electromagnetic field in the non-trivial anisotropic medium, the modified equation and the calculation domain corresponding to the modified equation can be discretized by using the finite volume method according to the structured hexahedral staggered grid shown in Figure 2. Among them, the finite volume method focuses on constructing discrete equations from a physical point of view. Each discrete equation is an expression of the conservation of a certain physical quantity on a finite volume. The physical concept of the derivation process is clear, the coefficients of the discrete equation have certain physical meanings, and can ensure that the discrete equation has conservation characteristics.

The first discrete result corresponding to the modified equation obtained through the discrete processing is a plurality of discrete equations, and the second discrete result corresponding to the calculation domain is a calculation domain corresponding to each discrete equation.

Then, the boundary conditions of the computational domain corresponding to each discrete equation can be determined. Considering the passive case, the boundary conditions can be the first-class boundary conditions, namely the Dirichlet boundary conditions. The field value on the boundary is determined by the electromagnetic response field value of the one-dimensional model of the boundary. The right-hand side of the equation can be formed by sorting. Therefore, according to each discrete equation and the boundary conditions of the computational domain corresponding to each discrete equation, a sparse linear equation group of the electromagnetic field can be constructed. This sparse linear equation group contains multiple discrete equations, so it is a large sparse linear equation group.

Finally, the sparse linear equations obtained above can be solved by the Krylov subspace solver to obtain the final numerical solution of the electromagnetic field, which is the result of the numerical simulation of the electromagnetic field of a non-trivial anisotropic medium. Among them, the main idea of the Krylov subspace solver is to recursively construct a residual vector for each iteration step, that is, the residual vector of the nth step is obtained by multiplying a polynomial of the coefficient matrix A with the first residual vector. However, it should be noted that the selection of the iterative polynomial should make the constructed residual vectors mutually orthogonal in a certain inner product sense, so as to ensure the minimum residual and achieve the purpose of rapid convergence.

In an embodiment of the present invention, when solving the correction equation, the finite volume method is used to discretize the correction equation and the calculation domain corresponding to the correction equation, and the boundary conditions of the calculation domain corresponding to each discrete equation are determined. A sparse linear equation group of the electromagnetic field is constructed through the discrete equations and the boundary conditions, and the Krylov subspace solver is used to solve it, which can greatly improve the solution speed, thereby improving the efficiency of numerical simulation of the electromagnetic field of non-trivial anisotropic media.

On the basis of the above embodiments, the numerical simulation method of the electromagnetic field of a non-trivial anisotropic medium provided in the embodiments of the present invention, the sparse linear equations are solved based on the Krylov subspace solver to obtain the numerical solution of the electromagnetic field, specifically comprising:

Based on the Krylov subspace solver, the quasi-minimum residual method is adopted to iteratively solve the sparse linear equations to obtain the numerical solution of the electromagnetic field.

Specifically, in an embodiment of the present invention, a sparse linear equation system can be iteratively solved by a Krylov subspace solver using a quasi-minimal residual (QMR) method to obtain a numerical solution to the electromagnetic field, which can further improve the solution speed and the efficiency of numerical simulation of electromagnetic fields in non-trivial anisotropic media.

The following examples are all based on the application of magnetotelluric (MT) numerical simulation in geophysics, and are all implemented with three-dimensional numerical algorithms and simulated in three-dimensional space. The first type of boundary condition used will impose a thicker air layer on the top of the model, which is set to low conductivity and is solved synchronously with the lower earth model; in the accuracy test, for MT, it is also necessary to convert the electromagnetic field obtained by numerical simulation into impedance and then calculate the apparent resistivity.

1) Efficiency test:

The embodiment of the present invention sets three different types of non-trivial anisotropic media and one ordinary anisotropic medium for the three models of half-space model, one-dimensional layered model and three-dimensional model, so as to comprehensively test the effect of the numerical simulation method of electromagnetic field of non-trivial anisotropic medium provided in the embodiment of the present invention. Among them, the three different types of non-trivial anisotropic media include arbitrary anisotropic media, horizontal anisotropic media and inclined anisotropic media, and the ordinary anisotropic medium includes axial anisotropic media. Table 1 shows the electrical parameters of different anisotropic half-space models. During the test, the spatial dimensions of the calculation domain along the X, Y and Z axes are set to 192km, 192km and 330km respectively, and are divided into 26×26×59 hexahedral grids.

Table 1 Electrical parameters of half-space models of different types of anisotropic media

Figure 4 shows the electrical parameters of the one-dimensional layered model of different types of anisotropic media, where LM1 represents arbitrary anisotropy, LM2 represents horizontal anisotropy, LM3 represents tilt anisotropy, and LM4 represents axial anisotropy. The spatial dimension and mesh division of the model calculation domain are the same as those of the above-mentioned half-space model. Figure 5 shows the electrical parameters of the three-dimensional COMMEMI 3D-2 model of different types of anisotropic media. In Figure 5, part (a) is a top view, with the horizontal axis in the east-west direction and the vertical axis in the north-south direction. In part (a), only block 1 (block1), block 2 (block2) and block 3 (block3) can be seen, and Profile is the position where the vertical axis takes a value of 0. Part (b) is a cross-sectional view of part (a) along Profile, with the horizontal axis in the east-west direction and the vertical axis in the depth direction. In part (b), block 1 (block1), block 2 (block2), block 3 (block3) and block 4 (block4) can be seen.

The electrical parameters of each block in Figure (5) are as follows:

block1 0.1S/m (isotropic medium)

block2 1/0.001/0.01S/m，[15°/45°/30°]

block3 0.01/1/0.01S/m，[30°/0°/0°]

block4 10S/m (isotropic medium)

For any anisotropic medium, keep all angles in the square brackets; for horizontal anisotropic medium: keep the Strike in the square brackets of block 3; for inclined anisotropic medium: keep the Dip angle in block 2; for axial anisotropic medium: discard all angles in the square brackets. COMMEMI 3D-2 is a standard model proposed by the international electromagnetic induction research community to test numerical algorithms. The original model is provided by the isotropic algorithm. In the embodiment of the present invention, the isotropic parameters in the original model are modified to anisotropic parameters; its spatial dimensions are 125.626km, 125.626km, 145km, respectively, and are divided into 46×46×40 hexahedral grids.

After actual testing of the algorithm, the computational efficiency of all models of different types of anisotropic media has been significantly improved, and the time used has been significantly shortened, indicating that the numerical simulation method of electromagnetic fields of non-trivial anisotropic media in the embodiment of the present invention can be applied not only to all types of non-trivial anisotropic media, but also to ordinary anisotropic media, and ordinary anisotropic media can be regarded as a special case of non-trivial anisotropic media. However, to avoid redundancy, only the most general case of non-trivial anisotropy to be solved in the embodiment of the present invention, that is, the efficiency of electromagnetic field numerical simulation of arbitrary anisotropic media, is shown here.

FIG6 shows the variation of the normalized relative residual of the half-space model, one-dimensional layered model, and three-dimensional COMMEMI3D-2 model of arbitrary anisotropic media with the number of iterations, and indicates the time taken. The test cycle is 500 seconds, and the polarization mode is XY mode. Among them, QMR+CCGD represents the numerical simulation method of the electromagnetic field of a non-trivial anisotropic medium provided in an embodiment of the present invention, QMR+CCDC represents the existing current density divergence correction method, and QMR represents a pure iterative solution without correction or constraint. Part (a) of FIG6 corresponds to the half-space model, part (b) corresponds to the one-dimensional layered model, and part (c) corresponds to the three-dimensional COMMEMI 3D-2 model.

FIG6 illustrates that the numerical simulation method for electromagnetic fields of non-trivial anisotropic media provided in an embodiment of the present invention can suppress spurious solutions well and significantly improve the convergence efficiency of numerical simulation of electromagnetic fields of non-trivial anisotropic media.

2) Accuracy test:

In the embodiment of the present invention, the accuracy test is performed on the one-dimensional layered model of any anisotropic medium. The model layers and model parameters of the one-dimensional layered model are shown in Table 2.

Table 2 Electrical parameters of one-dimensional layered model of arbitrary anisotropic media

The spatial dimensions of the model calculation domain are 402km, 402km, and 530km, and are divided into 30×30×59 hexahedral grids. The test period is 1-10000 seconds, and the response calculation point is located at the center of the calculation domain. FIG7 shows a comparison between the electromagnetic response calculated by the numerical simulation method of the electromagnetic field of a non-trivial anisotropic medium provided in an embodiment of the present invention and the electromagnetic response calculated by the existing analytical method. The electromagnetic response calculated by the numerical simulation method of the electromagnetic field of a non-trivial anisotropic medium provided in an embodiment of the present invention corresponds to the symbol annotation, and the electromagnetic response calculated by the existing analytical method corresponds to the dashed line annotation. In FIG7, the left part is the curve relationship between the apparent resistivity and the period (Period), and the right part is the curve relationship between the impedance phase (Impedancephase) and the period (Period).

As can be seen from FIG7 , by comparing the four modes xy, yx, xx, and yy, the electromagnetic responses all match well, which verifies that the numerical simulation method for the electromagnetic field of a non-trivial anisotropic medium provided in the embodiment of the present invention can accurately calculate the electromagnetic response of a one-dimensional layered model of any anisotropic medium.

3) Accuracy test:

In the embodiment of the present invention, the accuracy test is performed on the asymmetric dyke model of any anisotropic medium. The asymmetric dyke model is shown in FIG8 . The asymmetric dyke model includes three parts (a), (b), and (c), and the electrical parameters of the three parts are shown in Table 3:

Table 3 Electrical parameters of asymmetric vein model for arbitrary anisotropic media

The model calculation domain space dimensions are 60km, 60km, 148km, and are divided into 70×70×44 hexahedral grids. The test period is 10 seconds. The section is located in the middle of the X-axis and along the Y-axis from the center to the right end. FIG9 shows the comparison between the electromagnetic response calculated by the non-trivial anisotropic medium electromagnetic field numerical simulation method provided in the embodiment of the present invention and the electromagnetic response calculated by the two-dimensional conventional method. The electromagnetic response calculated by the non-trivial anisotropic medium electromagnetic field numerical simulation method provided in the embodiment of the present invention corresponds to the symbol mark, and the electromagnetic response calculated by the existing analytical method corresponds to the dash mark.

As can be seen from FIG9 , the electromagnetic responses of the given modes also match well, verifying that the non-trivial anisotropic medium electromagnetic field numerical simulation method provided in the embodiment of the present invention can accurately calculate the electromagnetic response of a more complex asymmetric dyke model of any anisotropic medium. In addition, it also illustrates that the non-trivial anisotropic medium electromagnetic field numerical simulation method provided in the embodiment of the present invention is also applicable to application cases where the frequency is not very low.

4) Accuracy test:

In the embodiment of the present invention, the accuracy test is performed on the three-dimensional COMMEMI 3D-2 model of any anisotropic medium. The schematic diagram of the three-dimensional COMMEMI 3D-2 model is shown in Figure 5. For any anisotropic medium, all angles in the square brackets of each block need to be retained. The spatial dimensions and grid division of the model calculation domain are consistent with the relevant description of this model in the above embodiment (i.e., all are 125.626km, 125.626km, 145km, 46×46×40), the test period is 1-10000 seconds, and the coordinates of the test point for calculating the response are (62500, 62500, 0).

FIG10 shows a comparison between the electromagnetic response calculated by the non-trivial anisotropic medium electromagnetic field numerical simulation method provided in the embodiment of the present invention and the electromagnetic response calculated by the conventional direct solution method. The electromagnetic response calculated by the non-trivial anisotropic medium electromagnetic field numerical simulation method provided in the embodiment of the present invention corresponds to the symbol annotation, and the electromagnetic response calculated by the existing analytical method corresponds to the dashed line annotation. It can be seen from FIG10 that the responses of the four modes can also be well matched with each other, which verifies that the non-trivial anisotropic medium electromagnetic field numerical simulation method provided in the embodiment of the present invention can accurately calculate the electromagnetic response of the three-dimensional COMMEMI 3D-2 model of complex arbitrary anisotropic media.

In summary, the embodiments in the efficiency test and the embodiments in the accuracy test strongly illustrate the effectiveness of the numerical simulation method of electromagnetic fields of non-trivial anisotropic media provided in the embodiments of the present invention in suppressing pseudo-solutions in the numerical simulation of electromagnetic fields of non-trivial anisotropic media, avoiding the defects of the existing current density divergence correction method, and can not only improve the convergence efficiency of the algorithm, but also ensure the accuracy of the solution. It has universal applicability, and provides a new and effective means for efficient numerical simulation of electromagnetic fields of non-trivial anisotropic media, thereby providing computational guarantees for problems such as electromagnetic response calculation of complex media and fine inversion imaging of underground geological structures.

As shown in FIG. 11 , based on the above embodiment, an embodiment of the present invention provides a numerical simulation system for electromagnetic fields of non-trivial anisotropic media, including:

A first determination module 111 is used to determine the spatial topological position of the current density divergence constraint and the conductivity weighting factor in the non-trivial anisotropic medium;

A second determination module 112, configured to determine the current density divergence constraint item based on the spatial topological position and the conductivity weighting factor;

A correction module 113, used to correct the Maxwell electromagnetic field control equation based on the current density divergence constraint term to obtain a corrected equation;

The solution module 114 is used to solve the correction equation to obtain a numerical solution of the electromagnetic field in the non-trivial anisotropic medium.

On the basis of the above-mentioned embodiment, in the non-trivial anisotropic medium electromagnetic field numerical simulation system provided in the embodiment of the present invention, the first determination module is specifically used for:

Construct a structured hexahedral staggered grid;

For the electromagnetic field in the non-trivial anisotropic medium, the electric field vector of the electromagnetic field is defined on the edges of the structured hexahedral staggered grid, the magnetic field vector of the electromagnetic field is defined on the surface of the structured hexahedral staggered grid, and the current density divergence of the electromagnetic field is defined on the nodes of the structured hexahedral staggered grid;

Based on the location of the current density divergence, a spatial topological location of the current density divergence constraint is determined.

On the basis of the above-mentioned embodiment, in the non-trivial anisotropic medium electromagnetic field numerical simulation system provided in the embodiment of the present invention, the first determining module is further specifically used for:

For any node of the structured hexahedral staggered grid, averaging the diagonal elements in the grid conductivity tensor around the any node to the any node to obtain the diagonal elements of the conductivity tensor corresponding to the any node;

The root mean square of the diagonal elements of the conductivity tensor is calculated, and the inverse of the root mean square is used as the conductivity weighting factor.

On the basis of the above-mentioned embodiment, in the non-trivial anisotropic medium electromagnetic field numerical simulation system provided in the embodiment of the present invention, the second determination module is specifically used for:

Determining an alternative current density divergence constraint term based on the conductivity weighting factor and the spatial topological position;

The current density divergence constraint term is determined based on the candidate current density divergence constraint term.

On the basis of the above-mentioned embodiment, in the non-trivial anisotropic medium electromagnetic field numerical simulation system provided in the embodiment of the present invention, the second determination module is further specifically used for:

The candidate current density divergence constraint term is regularized to obtain the current density divergence constraint term.

On the basis of the above-mentioned embodiments, in the non-trivial anisotropic medium electromagnetic field numerical simulation system provided in the embodiments of the present invention, the solution module is specifically used for:

Based on the structured hexahedral staggered grid, the correction equation and the calculation domain corresponding to the correction equation are discretized by using the finite volume method to obtain a first discrete result corresponding to the correction equation and a second discrete result corresponding to the calculation domain;

Determining boundary conditions of the second discrete result, and determining a sparse linear equation group of the electromagnetic field based on the first discrete result and the boundary conditions;

The sparse linear equations are solved based on a Krylov subspace solver to obtain a numerical solution to the electromagnetic field.

On the basis of the above-mentioned embodiment, in the non-trivial anisotropic medium electromagnetic field numerical simulation system provided in the embodiment of the present invention, the solution module is further specifically used for:

Based on the Krylov subspace solver, the quasi-minimum residual method is adopted to iteratively solve the sparse linear equations to obtain the numerical solution of the electromagnetic field.

Specifically, the functions of each module in the non-trivial anisotropic medium electromagnetic field numerical simulation system provided in the embodiment of the present invention correspond one-to-one to the operating procedures of each step in the above method embodiment, and the effects achieved are also consistent. Please refer to the above embodiment for details, which will not be repeated in the embodiment of the present invention.

FIG12 illustrates a schematic diagram of the physical structure of an electronic device. As shown in FIG12, the electronic device may include: a processor 1210, a communication interface 1220, a memory 1230 and a communication bus 1240, wherein the processor 1210, the communication interface 1220 and the memory 1230 communicate with each other through the communication bus 1240. The processor 1210 may call the logic instructions in the memory 1230 to execute the numerical simulation method of the electromagnetic field of the non-trivial anisotropic medium provided by the above-mentioned embodiments, the method comprising: determining the spatial topological position and the conductivity weighting factor of the current density divergence constraint in the non-trivial anisotropic medium; determining the current density divergence constraint term based on the spatial topological position and the conductivity weighting factor; modifying the Maxwell electromagnetic field control equation based on the current density divergence constraint term to obtain a modified equation; solving the modified equation to obtain the numerical solution of the electromagnetic field in the non-trivial anisotropic medium.

In addition, the logic instructions in the above-mentioned memory 1230 can be implemented in the form of software functional units and can be stored in a computer-readable storage medium when sold or used as an independent product. Based on such an understanding, the technical solution of the present invention, in essence, or the part that contributes to the prior art or the part of the technical solution, can be embodied in the form of a software product, which is stored in a storage medium and includes several instructions for a computer device (which can be a personal computer, a server, or a network device, etc.) to perform all or part of the steps of the method described in each embodiment of the present invention. The aforementioned storage medium includes: various media that can store program codes, such as a USB flash drive, a mobile hard disk, a read-only memory (ROM), a random access memory (RAM), a magnetic disk or an optical disk.

On the other hand, the present invention also provides a computer program product, which includes a computer program stored on a non-transitory computer-readable storage medium, and the computer program includes program instructions. When the program instructions are executed by a computer, the computer can execute the non-trivial anisotropic medium electromagnetic field numerical simulation method provided by the above embodiments, the method comprising: determining the spatial topological position and conductivity weighting factor of the current density divergence constraint in the non-trivial anisotropic medium; determining the current density divergence constraint term based on the spatial topological position and the conductivity weighting factor; based on the current density divergence constraint term, correcting the Maxwell electromagnetic field control equation to obtain a corrected equation; solving the corrected equation to obtain the electromagnetic field numerical solution in the non-trivial anisotropic medium.

On the other hand, the present invention also provides a non-transitory computer-readable storage medium having a computer program stored thereon. When the computer program is executed by a processor, it is implemented to execute the numerical simulation method of the electromagnetic field of a non-trivial anisotropic medium provided by the above embodiments, the method comprising: determining the spatial topological position and the conductivity weighting factor of the current density divergence constraint in the non-trivial anisotropic medium; determining the current density divergence constraint term based on the spatial topological position and the conductivity weighting factor; based on the current density divergence constraint term, correcting the Maxwell electromagnetic field control equation to obtain a corrected equation; solving the corrected equation to obtain a numerical solution of the electromagnetic field in the non-trivial anisotropic medium.

The device embodiments described above are merely illustrative, wherein the units described as separate components may or may not be physically separated, and the components displayed as units may or may not be physical units, that is, they may be located in one place, or may be distributed on multiple network units. Some or all of the modules may be selected according to actual needs to achieve the purpose of the scheme of this embodiment. Those of ordinary skill in the art may understand and implement it without creative effort.

Through the description of the above implementation methods, those skilled in the art can clearly understand that each implementation method can be implemented by means of software plus a necessary general hardware platform, and of course, it can also be implemented by hardware. Based on this understanding, the above technical solution is essentially or the part that contributes to the prior art can be embodied in the form of a software product, and the computer software product can be stored in a computer-readable storage medium, such as ROM/RAM, a disk, an optical disk, etc., including a number of instructions for a computer device (which can be a personal computer, a server, or a network device, etc.) to execute the methods described in each embodiment or some parts of the embodiments.

Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the aforementioned embodiments, or make equivalent replacements for some of the technical features therein. However, these modifications or replacements do not deviate the essence of the corresponding technical solutions from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims (8)

1. A numerical simulation method of an electromagnetic field of a non-trivial anisotropic medium is characterized by comprising the following steps:

determining a spatial topological position constrained by current density divergence in a non-trivial anisotropic medium and a conductivity weighting factor;

determining the current density divergence constraint term based on the spatial topological position and the conductivity weighting factor;

correcting a Maxwell electromagnetic field control equation based on the current density divergence constraint term to obtain a correction equation;

solving the correction equation to obtain an electromagnetic field numerical solution in the nontrivial anisotropic medium;

wherein the spatial topological position is determined based on:

constructing a structured hexahedral staggered grid;

for an electromagnetic field in the non-trivial anisotropic medium, defining an electric field vector of the electromagnetic field on edges of the structured hexahedral interleaved grid, defining a magnetic field vector of the electromagnetic field on faces of the structured hexahedral interleaved grid, defining a current density divergence of the electromagnetic field on nodes of the structured hexahedral interleaved grid;

determining a spatial topological position of the current density divergence constraint based on the position of the current density divergence;

solving the correction equation to obtain a numerical solution of the electromagnetic field in the nontrivial anisotropic medium, specifically comprising:

based on the structured hexahedral staggered grid, performing discrete processing on the correction equation and a calculation domain corresponding to the correction equation by adopting a finite volume method to obtain a first discrete result corresponding to the correction equation and a second discrete result corresponding to the calculation domain;

determining a boundary condition for the second discrete result and determining a sparse system of linear equations for the electromagnetic field based on the first discrete result and the boundary condition;

and solving the sparse linear equation set based on a Krylov subspace solver to obtain the electromagnetic field numerical solution.

2. A non-trivial anisotropic medium electromagnetic field numerical simulation method according to claim 1, wherein the conductivity weighting factor is determined based on:

for any node of the structured hexahedral staggered grid, averaging diagonal elements in the grid conductivity tensor around the any node to obtain a conductivity tensor diagonal element corresponding to the any node;

a root mean square of the diagonal elements of the conductivity tensor is computed, and an inverse of the root mean square is used as the conductivity weighting factor.

3. The method of numerical simulation of an electromagnetic field of a non-trivial anisotropic medium as set forth in claim 1, wherein said determining said current density divergence constraint term based on said spatial topological position and said conductivity weighting factor specifically comprises:

determining an alternative current density divergence constraint term based on the conductivity weighting factor and the spatial topological position;

determining the current density divergence constraint term based on the alternative current density divergence constraint term.

4. A non-trivial anisotropic medium electromagnetic field numerical simulation method according to claim 3, wherein said determining the current density divergence constraint term based on the candidate current density divergence constraint term, comprises in particular:

and carrying out regularization processing on the alternative current density divergence constraint term to obtain the current density divergence constraint term.

5. The non-trivial anisotropic medium electromagnetic field numerical simulation method of claim 1, wherein the solving the sparse linear equation set based on Krylov subspace solver to obtain the electromagnetic field numerical solution, specifically comprises:

and iteratively solving the sparse linear equation set by adopting a quasi-minimum residual error method based on a Krylov subspace solver to obtain the electromagnetic field numerical solution.

6. A numerical simulation system for electromagnetic fields of non-trivial anisotropic media, comprising:

a first determination module for determining a current density divergence constrained spatial topological position in a non-trivial anisotropic medium and a conductivity weighting factor;

a second determination module to determine the current density divergence constraint term based on the spatial topological position and the conductivity weighting factor;

the correction module is used for correcting a Maxwell electromagnetic field control equation based on the current density divergence constraint term to obtain a correction equation;

the solving module is used for solving the correction equation to obtain an electromagnetic field numerical solution in the non-trivial anisotropic medium;

wherein the spatial topological position is determined based on:

constructing a structured hexahedral staggered grid;

for an electromagnetic field in the non-trivial anisotropic medium, defining an electric field vector of the electromagnetic field on edges of the structured hexahedral interleaved grid, defining a magnetic field vector of the electromagnetic field on faces of the structured hexahedral interleaved grid, defining a current density divergence of the electromagnetic field on nodes of the structured hexahedral interleaved grid;

determining a spatial topological position of the current density divergence constraint based on the position of the current density divergence;

solving the correction equation to obtain a numerical solution of the electromagnetic field in the nontrivial anisotropic medium, specifically comprising:

based on the structured hexahedral staggered grid, performing discrete processing on the correction equation and a calculation domain corresponding to the correction equation by adopting a finite volume method to obtain a first discrete result corresponding to the correction equation and a second discrete result corresponding to the calculation domain;

determining a boundary condition for the second discrete result and determining a sparse system of linear equations for the electromagnetic field based on the first discrete result and the boundary condition;

and solving the sparse linear equation set based on a Krylov subspace solver to obtain the electromagnetic field numerical solution.

7. An electronic device comprising a memory, a processor and a computer program stored on said memory and executable on said processor, characterized in that said processor, when executing said program, implements the steps of the method for numerical simulation of electromagnetic fields of a non-trivial anisotropic medium as claimed in any one of claims 1 to 5.

8. A non-transitory computer readable storage medium having stored thereon a computer program, wherein the computer program when executed by a processor implements the steps of the method for numerical simulation of electromagnetic fields of a non-trivial anisotropic medium as claimed in any one of claims 1 to 5.

Images