A method for establishing node voltage resistance network based on Octree grid

Abstract

The present invention discloses a method for establishing a node voltage resistance network based on an Octree grid, and builds a node voltage resistance network system based on an expandable Octree grid. The system contains all current conservation equations and voltage conservation equations without losing any unit information, and its system matrix has symmetry and excellent solution efficiency. Numerical example results show that compared with traditional methods, this method can improve the accuracy of numerical solutions from 1% to 0.1%, and reduce the calculation time by 80%, while saving 78% of the number of units. Taking the bedrock survey in western Hubei as an example, the inversion effect is consistent with the changes in gravity high-order wavelets, aeromagnetic measurements and seismic results, verifying the practical application effect of this method.

Description

Technical Field

[0001] The present invention relates to the technical field of geoelectric structure inversion, and in particular to a method for establishing a node voltage resistance network based on an Octree grid. Background Art

[0002] Earth return current is the DC current generated in the loop between the transmission line and the ground during monopolar operation of the power grid. This DC current flows into the ground through the neutral point, and its propagation and distribution underground are affected by the underground electrical structure. Accurately calculating the earth return current through the neutral point to analyze the DC bias magnetic effect helps ensure safe power grid operation. Furthermore, since earth return current is affected by the underground electrical medium, we can construct an inverse problem for the earth return current process and use the surface earth return current signal to infer the electrical medium structure.

[0003] The propagation path of ground return current consists of two components: the ground and the power grid. A common approach to solving this problem is to treat the underground as an equivalent resistor and construct a resistor network (RN) containing both underground electrical and power grid resistance information. However, the ground's electrical structure is often complex and exhibits three-dimensional inhomogeneity. The equivalent earth resistance can only simulate simple layered geoelectrical structures, which is a major limitation on the accuracy of this method's numerical solution.

[0004] Differential discretization calculation methods are commonly used to accurately simulate the propagation and distribution of direct current in the earth. These methods include the finite element method, the finite difference method, and the finite volume method. However, these differential methods have difficulty discretizing long, thin, highly conductive wires, as otherwise they would construct differential system equations with a large number of degrees of freedom. Based on this, Kirchhoff's resistor network was proposed, which constructs node current and voltage conservation at each differential element. This resistor network-based differential discretization method avoids the difficulty of dissecting high-conductors in traditional differential methods and is suitable for computational problems involving the connection between long, straight conductors, steel casing, and geoelectric structures.

[0005] The traditional Kirchhoff's resistor network method directly superimposes and assembles the current conservation equation and the voltage conservation equation. The system is an overdetermined equation with more equations than the number of unknowns, and is usually a contradictory equation with no solution. In order to construct a square matrix that is equivalent to a full row rank matrix and a full column rank matrix, the Kirchhoff's resistor network method needs to discard some current conservation equations and voltage conservation equations, which leads to the loss of some valid unit information and the loss of numerical solution accuracy. Moreover, the matrix constructed by this method is asymmetric, with upper and lower diagonal elements corresponding to the current conservation equation and the voltage conservation equation, which greatly reduces the computational efficiency. Summary of the Invention

[0006] The present invention aims to provide a method for establishing a node voltage resistance network based on an Octree grid, so as to solve the problems of poor numerical solution accuracy and low computational efficiency of the existing Kirchhoff's resistance network.

[0007] In order to achieve the above object, the present invention provides the following technical solutions:

[0008] A method for establishing a node voltage resistance network based on an Octree grid comprises the following steps:

[0009] S1. Based on the Kirchhoff resistor network, the current is stipulated to flow from node p to the other six subsidiary nodes, and the equation is established:

[0010] U p Y pp -∑ j U j Y pj =i p (3)

[0011]

[0012] Among them, Y pp is the self-admittance, which indicates the self-association relationship of the p-th node, Y pi is the mutual admittance, which represents the association relationship between the p-th node and the j-th subsidiary node;

[0013] S2. Assemble the node voltage system equation according to the unit analysis. The system matrix and the right-hand side terms can be expressed as:

[0014] Y N×N U N×1 =I N×1 (5)

[0015]

[0016] Where N is the total number of unit nodes. In the node admittance matrix Y, when p = q, the matrix elements are self-admittances, and when p ≠ q, the matrix elements are mutual admittances. When the qth node is not an adjacent subsidiary node of the pth node, the mutual admittance is 0.

[0017] S3. Define the node at the center point of the unit. The lines connecting the center point of the unit and the centers of the six faces point to the X, Y, and Z directions respectively. Set two centrally symmetrical resistors in each direction. Let the size of the unit grid be dx, dy, and dz. The resistivity of the unit grid is ρ. i , the node resistance is expressed as:

[0018]

[0019] Among them, R x 、Ry and R z Represents the resistance in the X, Y, and Z directions respectively;

[0020] S4. Use Octree grid to refine the target area by level. Hanging nodes are generated on the coplanar surface of units at different levels. The hanging node resistance R g Expressed as:

[0021]

[0022] S5. Add the overhead line as a resistor to the node admittance matrix and rewrite formula (3) as:

[0023]

[0024] S6. According to step S5, use formula (7) to assemble into a system equation of the form of formula (5);

[0025] S7. Use GMRES as a solver and ILU as a preconditioner to solve the system equation generated in step S6 to obtain the voltage value of the node voltage system.

[0026] In step S4, the hanging node resistance R g The derivation process is:

[0027] When cells at different levels are coplanar, consider the hanging nodes and the voltage of the cell at the previous level is u i,j,k , the refined unit voltage is u i+1 / 2,j+1 / 2,k+1 / 2 , the two-level unit shares a resistor R through the coplanar g ;

[0028] Introducing ghost point u g , and its interpolation function is as follows:

[0029]

[0030] Among them, w is the weight and p is the interpolation order;

[0031] For second-order interpolation, the ghost point u g This can be achieved through planar interpolation:

[0032]

[0033]

[0034]

[0035] Approximately equate u2 and u3, and let Δu 21 =u2-u1, we can get:

[0036]

[0037]

[0038] Through the resistance relationship. Let the resistance:

[0039]

[0040]

[0041] Consider dx i =2dx i+1 / 2 , R i+1 / 2,j+1 / 2,k+1 / 2 Can be changed to:

[0042]

[0043] R g The value of R i+1 / 2,j+1 / 2,k+1 / 2 With R i,j,k sum:

[0044]

[0045] The beneficial effects of the present invention are:

[0046] Based on Kirchhoff's resistance network method, the present invention uses the node voltage principle to construct a symmetric matrix system that contains both the current conservation equation and the voltage conservation equation. This system not only does not lose any unit information, but also has a more extraordinary solution efficiency. In addition, in order to avoid redundant units caused by the lateral expansion of regular grid encryption, the present invention also introduces a more flexible and expandable octree grid to discretize the underground space. This grid is encrypted at the source point and near the surface to improve the grid accuracy. Numerical example results show that compared with traditional methods, this innovative method can improve the accuracy from 1% to 0.1%, reduce the calculation time by 80%, and save 78% of the number of units. Taking the bedrock survey in western Hubei as an example, it is verified that the actual inversion effect of this method is consistent with the changes in gravity high-order wavelets, aeromagnetic measurements and seismic results. BRIEF DESCRIPTION OF THE DRAWINGS

[0047] Figure 1 Schematic diagram of the unit center resistance network, (a) shows the unit center resistance, (b) shows the unit position in the network and the unit connection;

[0048] Figure 2 Schematic diagram of octree mesh refinement;

[0049] Figure 3 Schematic diagram of the hanging node in the xz longitudinal slice perspective;

[0050] Figure 4Schematic diagram of adding overhead line resistance to the resistance network;

[0051] Figure 5 Schematic diagram of the uniform half-space model;

[0052] Figure 6 Schematic diagram of mesh division, (a) shows the schematic diagram of mesh subdivision level in the near-surface area, (b) shows the schematic diagram of mesh subdivision level at the source or grounding electrode, (c) shows the mesh refinement of the grounding electrode and current source by K-RN, and (d) shows the refinement of the grounding electrode by free tetrahedron elements in FEM.

[0053] Figure 7 Schematic diagram of the main element distribution of the system matrix generated for K-RN and NV-RN, (a) shows the main element distribution of the system matrix of K-RN, (b) shows the main element distribution of the system matrix of NV-RN;

[0054] Figure 8 Schematic diagram of the decrease of residual error with the iteration steps in the solution process of K-RN and NV-RN;

[0055] Figure 9 Schematic diagram of the calculation results of RN and FEM numerical simulation of earth return flow, (a) shows the voltage response result on the measuring line, (b) shows the voltage response error on the measuring line, (c) shows the voltage response result at the grounding electrode, and (d) shows the voltage response error at the grounding electrode;

[0056] Figure 10 This is an overview of the work area in Hubei Province. The western Hubei region is marked with a dotted box as the target area.

[0057] Figure 11 This is a schematic diagram of the electrical structure of the western Hubei region based on NV-RN inversion;

[0058] Figure 12 Diagram of pseudocode for solving the system equations for node voltages. DETAILED DESCRIPTION

[0059] The present invention will be further described in detail below with reference to the accompanying drawings and embodiments:

[0060] A method for establishing a node voltage resistance network based on an Octree grid, comprising:

[0061] In Kirchhoff's resistor network, Kirchhoff's voltage law is expressed as:

[0062] ∑(I1R1+I2R2)=0 (1)

[0063] When the underground medium is discretized into a regular cubic unit body, each unit node is connected to six edges, and the current is stipulated to flow from node p to the other six subsidiary nodes. Substituting Ohm's law into formula (1) yields:

[0064]

[0065] Among them, R p represents the resistance between the pth node and its six subsidiary nodes. The above formula can be rewritten as:

[0066] U p Y pp -∑ j U i Y pj =i p (3)

[0067]

[0068] Formula (4) is the admittance relationship between the pth node and its subordinate nodes. pp is the self-admittance, which indicates the self-association relationship of the p-th node, Y pi is the mutual admittance, representing the relationship between the pth node and the jth subsidiary node. Formula (2) is a rewrite based on Kirchhoff's current law (KCL), automatically satisfying KCL. Furthermore, this formula incorporates the voltage relationship into the equation, improving the surface-based Kirchhoff's voltage law (KVL) to a volume-based KVL.

[0069] After each edge and node is numbered in sequence, we can assemble the node voltage system equation according to the unit analysis. The system matrix and the right-hand side terms can be expressed as:

[0070] Y N×N U N×1 =I N×1 (5)

[0071]

[0072] Where N is the total number of unit nodes. In the node admittance matrix Y, when p = q, the matrix elements are self-admittances, and when p ≠ q, the matrix elements are mutual admittances. When the qth node is not a neighboring subsidiary node of the pth node, the mutual admittance is 0. Therefore, the node admittance matrix is ​​a sparse symmetric matrix. When the number of units and the ordering method are determined, the width of the matrix is ​​also fixed. Compared with traditional methods, this method obtains a conjugate Hermite matrix for the system matrix, which does not require additional equation deletion. After solving for the node voltage, the current flowing through the resistance between the nodes can be calculated using Ohm's law.

[0073] Obtaining the resistance between nodes is a necessary condition for using the node element method. The present invention uses the unit center method to define the node at the center point of the unit, such as Figure 1 As shown in a. The lines connecting the center of the unit and the centers of the six faces can determine three axes, pointing to the X, Y, and Z directions respectively. Two resistors are set in each direction with central symmetry. Let the size of the unit grid be dx, dy, and dz, and the resistivity of the unit grid be ρ i , the resistance of each resistor can be determined by the axis length and cross-section.

[0074]

[0075] Among them, R x 、R y and R z Represents the resistance in the X, Y, and Z directions respectively.

[0076] The connection between the unit and other units is as follows Figure 1 As shown in Figure b, two resistors are connected between each node, each belonging to its own unit. This approach avoids ambiguous definitions when the nodes are coplanar or coaxial. Since there is no branch circuit between the two nodes, the two resistors belonging to different units can be considered to be connected in series. The directly added resistance values ​​can be used to calculate the mutual admittance between the two nodes and the node self-admittance.

[0077] The present invention introduces the octree grid to subdivide the target area into layers, such as Figure 2 As shown in Figure 2, hanging nodes will be generated on the coplanar surfaces of cells at different levels. Therefore, the analysis of current nodes and voltage nodes based on regular cells needs to be improved to ensure the stability of the node admittance matrix.

[0078] like Figure 3 As shown, the left unit is the previous level unit, and the right unit is the refined unit. Without considering the hanging nodes, the voltages corresponding to the units at each level are u i,j,k and u i+1 / 2,j+1 / 2,k+1 / 2 , the corresponding resistances are R i,j,k and R i+1 / 2,j+1 / 2,k+1 / 2 When the hanging nodes are taken into account, the voltage corresponding to each level unit is still u i,j,k and u i+1 / 2,j+1 / 2,k+1 / 2 , the resistors become shared by a resistor R across the coplanar g .

[0079] Introducing ghost point u g , and its interpolation function is as follows:

[0080]

[0081] Among them, w is the weight and p is the interpolation order.

[0082] For second-order interpolation, the ghost point u g This can be achieved through planar interpolation:

[0083]

[0084]

[0085]

[0086] Approximately equate u2 and u3, and let Δu 21 =u2-u1, we can get:

[0087]

[0088]

[0089] It can be seen from this that Δu g1 and Δu 2g There is a 2-fold relationship.

[0090] According to the resistance relationship, let the resistance:

[0091]

[0092]

[0093] Consider dx i =2dx i+1 / 2 , R i+1 / 2,j+1 / 2,k+1 / 2 Can be changed to:

[0094]

[0095] R i+1 / 2,j+1 / 2,k+1 / 2 With R i,j,k The 2-fold relationship is also satisfied. Under the condition of current conservation, the linear proportional relationship between resistance and voltage is confirmed. g The value should be R i+1 / 2,j+1 / 2,k+1 / 2 With R i,j,k sum:

[0096]

[0097] Assume that there is a hanging node on the i-th surface of the unit, and the resistance R in the unit pi (Formula 2) can be used R g Resistance representation.

[0098] For boundary conditions, the node voltage resistance network of the present invention is automatically cut off at the ground-air boundary, ensuring that current does not flow from underground into the air. For the underground boundary, a resistor with a resistance value smaller than the resistance value inside the boundary is added to the outside of the boundary. This resistor with a smaller resistance value simulates the attraction of the infinite range on the current in the calculation area. In order to avoid the difficulty in determining the resistance on the coplanar edges, our calculation unit uses the unit center unit instead of the vertex center unit. This will make it impossible to calculate the surface voltage, and the shallowest voltage can only be obtained at a depth of half the longitudinal dimension of the unit underground. However, the neutral point of the substation is buried underground, and this unit center method does not affect the calculation of the neutral point potential. Not only that, when the calculation area is large enough, the final results calculated by the unit center and vertex center methods are very different. When the calculation model is small and in laboratory conditions, a truncated boundary can be used to prevent current from flowing out of the calculation model area from the boundary.

[0099] When the currents in the positive and negative electrodes are too large, volume discretization can be used to distribute the current at a node to multiple surrounding grid nodes to ensure the harmony of the divergence field. This can avoid the source singularity problem to a certain extent and improve the calculation accuracy.

[0100] After establishing a complete node voltage resistance network, the overhead line can be directly added as a resistor into the node admittance matrix, such as Figure 4 As shown in Figure 3, when the two ends of the overhead line are connected to the pth unit and the p+kth unit in the calculation model respectively, the admittance correlation formula of a single node (Formula 3) needs to be rewritten as:

[0101]

[0102] When the grid order remains unchanged, equation (7) can still be used to assemble a system of equations similar to equation (5). When the source is not located at the overhead line grounding point, only individual elements in the system matrix change. When the source is located at the overhead line, new non-zero elements are added to the right-hand side of the system matrix.

[0103] The system matrix is ​​solved by ILU preconditioning and GMRES using Matlab, as shown in the following example: Figure 12 As shown:

[0104] The generated system matrix A is input to the Matlab function ilu for ILU decomposition. When using this function, set the 'type' parameter to 'ilutp' to perform an ILU decomposition with thresholding and pivot selection; set the 'droptol' parameter to 1e-5 to specify a drop tolerance of 1e-5 for the ILU decomposition. Other parameters follow the Matlab default settings. After processing the system matrix A with the ilu function, the upper and lower triangular matrices, L and U, are obtained.

[0105] The system matrix A and the right-hand side vector b, along with the resulting matrices L and U, are input into the Matlab function gmres for solution. When using this function, set the tolerance tol to 1e-12 and the maximum number of iterations to 3000. The voltage values ​​of the nodal voltage system are obtained after the solution.

[0106] Numerical model

[0107] The uniform half-space model is designed to be 98km×98km×98km with a resistivity of 100Ω·m. The ground current and the ground current are both 6250A. The coordinates of the current exit and entry points, i.e., the positive and negative poles of the current source, are (0,-20,0)km and (0,20,0)km respectively. The grounding points at both ends of the conductor are (-25,0,-25)km and (25,0,25)km respectively. The current source and measuring point locations are shown in the figure below. Figure 5 shown.

[0108] The DC resistance of the conductor is set to 0.05728Ω / km (usually 0.04-0.07Ω / km), and the resistance of the double-strand conductor is 0.02864Ω / km. The conductor distance is 70.71km, and the conductor resistance is 2.025Ω. To simplify the model, the multi-strand high-voltage line is treated as a single strand, and the conductor cross-section is set to 1e-2m. 2 If the cross section of the conductor is circular, the conductivity of the conductor is 3.4916E6 S / m.

[0109] The regular hexahedral mesh of the node voltage resistance network based on the Octree grid is M NV-RN , the overall unit size is 1km. In the near-surface area, the vertical height of the starting unit is 15.625m, and it increases to 1km according to the octree grid ( Figure 6 a), with 11 layers of cells. When the vertical cell size increases to 1 km, a grid with 88 layers of cells of equal thickness is constructed. The bottom cell extends to a depth of 92 km. The grid is denser at the positive and negative current poles and the ground ( Figure 6 b) The total number of cells in the core computation area is 0.98E6.

[0110] The regular hexahedral grid of Kirchhoff's resistor network is M K-RN The overall horizontal unit size is 1km, and the grid is encrypted with a unit size of 50m. The encrypted grid will expand to the boundary. The vertical size of the grid is as follows: NV-RN The number of cells in the underground grid is 4.3E6.

[0111] Two finite element methods with different meshes are used to compare the accuracy of the RN method. FEM-1 and M FEM-2 is a regular hexahedral grid, M FEM-3、M FEM-4 and M FEM-5 is a free tetrahedral mesh. In addition to the different mesh types, these meshes also have different densities. FEM-5 The grid used to calculate the standard response value is obtained through the a posteriori error adaptive algorithm.

[0112] M FEM-1 and M FEM-2 The underground structure is divided into regular hexahedron grids with a horizontal unit length of 1 km. FEM-1 and M FEM-2 In the figure, the number of underground subdivision units is 0.97E6. FEM-1 The number of meshes for the wire segmentation is 10598. FEM-2 The wire is meshed more densely, with a mesh size of 62222. FEM-1 and M FEM-2 The total number of grids is 0.98E6 and 1.04e6 respectively. FEM-3 and M FEM-4 M is an unstructured tetrahedral grid with different densities to divide the underground structure. FEM-3 The number of meshes in the whole is 0.06E6, and the number of meshes in the wire is 6996; FEM-4 The overall mesh number is 0.12E6, and the mesh number of the wire is 61191. FEM-Bench In the case of the adaptive standard solution, the underground structure is meshed with unstructured tetrahedrons and the adaptive encryption method is used at the measuring points. The overall number of meshes is 2.52E6, and the number of meshes for the conductor is 117481.

[0113] Table 1 RN and FEM mesh generation

[0114]

[0115] The system matrices generated by the two RN methods are as follows: Figure 7 As shown. The system matrix generated by K-RN is mainly composed of two diagonal matrices, which correspond to the behavior of simulating the propagation process of current using Kirchhoff current and voltage laws. NV-RN generates a symmetric diagonal system matrix without losing any unit information. Using GMRES as the solver and ILU as the preconditioner, the solution process of the two methods is compared as follows Figure 8As shown in the figure. Without the ILU preconditioner, both K-RN and NV-RN are difficult to solve. The residual error of K-RN does not decrease, while the residual error of NV-RN decreases slowly with each iteration. With the ILU preconditioner, the performance of the GMRES solver is significantly improved. The ILU solver for K-RN takes 1 hour, 21 minutes, and 30 seconds, while the solver for NV-RN takes only 4 minutes and 55 seconds. This indicates that the diagonal system matrix generated by NV-RN is more conducive to solving.

[0116] Comparing the results of K-RN and NV-RN, NV-RN achieves higher numerical accuracy using fewer cells. This is due to NV-RN's strategy of constructing basic resistor network cells around the cell center. It does not rely on the continuity of regular cells to establish a surface-based circulation equipotential surface, as the K-RN method does. Ultimately, NV-RN can use octree meshes or multi-level meshes.

[0117] The response results and errors of the survey lines calculated by various methods are as follows: Figure 9 As shown. By M FEM-1 and M FEM-2 The comparison of the calculation results shows that when the grids for dividing the underground are the same, the denser the grid for dividing the conductor, the higher the accuracy of the numerical solution. FEM-3 and M FEM-4 Comparison of the calculation results also demonstrates the importance of wire segmentation in improving the accuracy of the numerical solution. RN, a method that does not require additional wire segmentation, achieves numerical solution errors of less than 5% for grounding electrodes, demonstrating its superiority in handling long straight wires. Presumably, the greater the number and length of wires, the greater the advantage of the RN method.

[0118] The present invention also predicts the electrical structure through the earth return current and potential fluctuation of the power grid in western Hubei Province, and compares the results with other geophysical detection methods, gravity, magnetic and seismic detection, to test the effect of its NV-RN.

[0119] The neutral point location of the substation in western Hubei serves as the observation point, and the location of the DC grounding electrode serves as the source location. During single-unit operation, the neutral point current is recorded as observation data. Using the node voltage-resistance network of the present invention, the electrical structure of the target area is inverted to obtain an inferred electrical structure.

[0120] The inversion area covers a horizontal area of ​​300 km by 248 km, and its vertical depth exceeds the Moho, reaching 60 km. The density of source and measurement points near the surface of the inversion area is shown in Figure 10. The maximum surface cell size is 1 km, and the minimum surface cell size is 125 m. The cell size gradually increases to 4 km. The total number of cells is 6.2E6.

[0121] When the inversion fitting error reaches 10%, the predicted geoelectric model is as follows Figure 11 The horizontal slice coordinates of the geoelectric model are z = -1 km, z = -4 km, and z = -7 km. 500 Ω·m equipotential surfaces are depicted in the geoelectric model to infer the distribution of the underground bedrock.

[0122] The predicted geoelectric model shows a bedrock structure extending from the surface deep into the ground, gradually extending vertically from shallow to deep, toward the negative X-axis. Furthermore, bedrock is present on both sides of this structure, 7 km below the negative X-axis.

[0123] In the negative Y-axis direction, there are scattered bedrock bodies. These bedrock bodies are relatively scarce at the foreground (on the z = -1 km slice), but become increasingly dense with increasing depth. The bedrock bodies also increase in size.

[0124] In the geoelectric structure, the bedrock distribution in the Yichang and Zigui regions is accurately depicted through gravity and aeromagnetic ΔT anomaly results. This inversion is consistent with high-order gravity wavelet variations, aeromagnetic surveys, and seismic data.

[0125] The above is only an embodiment of the present invention, and the common knowledge such as the specific technical solutions and / or characteristics in the solution are not described in detail here. It should be pointed out that for those skilled in the art, without departing from the technical solution of the present invention, several variations and improvements can be made, which should also be regarded as the scope of protection of the present invention, and these will not affect the effect of the implementation of the present invention and the practicality of the patent. The scope of protection required by this application shall be based on the content of its claims, and the specific implementation methods and other records in the description can be used to interpret the content of the claims.

Claims

1. A method for establishing a node voltage resistance network based on an Octree grid, characterized in that: Including steps: S1. Based on the Kirchhoff resistor network, the current is stipulated to flow from node p to the other six subsidiary nodes, and the equation is established: (1) (2) in, Y pp is the self-admittance, indicating the p The self-association relationship of each node, Y pi is the mutual admittance, indicating the p Node and j The relationship between the affiliated nodes; S2. Assemble the node voltage system equation according to the unit analysis, namely formula (1) and formula (2), and the system matrix and the right-hand side term are expressed as: (3) in, N is the total number of unit nodes, in the node admittance matrix Y, p = q When the matrix element is the self-admittance, p ≠ q When the matrix elements are mutual admittance, the q The node is not p When it is a subsidiary node of a node, the mutual admittance is 0; S3. Define the node at the center point of the unit. The lines connecting the center point of the unit and the centers of the six faces point to the X, Y, and Z directions respectively. Set two centrally symmetrical resistors in each direction. Let the size of the unit grid be dx 、 dy and dz , the resistivity of the unit grid is ρ i， The node resistance is expressed as: （4） in, R x 、 R y and R z Represents the resistance in the X, Y, and Z directions respectively; S4. Use Octree grid to refine the target area by level, and generate hanging nodes on the coplanar surface of different level units. The hanging node resistance R g Expressed as: （5） S5. Add the overhead line as a resistor to the node admittance matrix and rewrite formula (3) as: (6) S6. According to step S5, use formula (6) to assemble the system equation of formula (3); S7. Use GMRES as a solver and ILU as a preconditioner to solve the system equation generated in step S6 to obtain the voltage value of the node voltage system.

2. The method for establishing a node voltage resistance network based on an Octree grid according to claim 1, wherein: In step S4, the hanging node resistance R g The derivation process is: When cells at different levels are coplanar, considering the hanging nodes, the voltage of the cell at the previous level is u i,j,k， The refined unit voltage is u i+1 / 2,j+1 / 2,k+1 / 2 , the two-level cells share a resistor across the coplanar R g ； Introducing ghost points u g , and its interpolation function is as follows: （7） in, w is the weight, p is the interpolation order; For second-order interpolation, ghost points u g This is achieved through planar interpolation: （8） （9） （10） Will u 2 and u 3 is equivalent, let Δ u 21 = u 2- u 1 can be obtained: （11） （12） According to the resistance relationship, let the resistance: （13） （14） consider dx i =2 dx i+1 / 2 , R i+1 / 2,j+1 / 2,k+1 / 2 Changes to: （15） R g The value of R i+1 / 2,j+1 / 2,k+1 / 2 and R i,j,k sum: （16）。

Images