CN113887102B

CN113887102B - Full-wave electromagnetic simulation method and system for integrated circuit under lossless frequency dispersion medium

Info

Publication number: CN113887102B
Application number: CN202111166771.4A
Authority: CN
Inventors: 王芬
Original assignee: Beijing Wisechip Simulation Technology Co Ltd
Current assignee: Beijing Wisechip Simulation Technology Co Ltd
Priority date: 2021-09-30
Filing date: 2021-09-30
Publication date: 2022-03-11
Anticipated expiration: 2041-09-30
Also published as: CN113887102A

Abstract

The application discloses an integrated circuit full-wave electromagnetic simulation method and system under lossless frequency dispersion medium, firstly model information of a super-large scale integrated circuit is obtained, the integrated circuit is modeled, then grid subdivision is carried out on a parallel flat plate field of the integrated circuit based on a model obtained through modeling, then a matrix equation under the lossless frequency dispersion medium is established based on a result of the grid subdivision, then a reference frequency point and a field solution of the integrated circuit are calculated through an iteration method based on the matrix equation, then a proportional coefficient between the field solution of the reference frequency point and a low-frequency field solution is obtained for a frequency point to be solved which is lower than the reference frequency point, the field solution based on the reference frequency point and the proportional coefficient obtain the field solution under the frequency point to be solved, and finally, an electromagnetic response of a full frequency band is obtained based on the field solution of all the frequency points to be solved. The method can calculate the reference frequency point of the full-wave analysis, and accurately solve the field solution of the low-frequency point based on the reference frequency point.

Description

Full-wave electromagnetic simulation method and system for integrated circuit under lossless frequency dispersion medium

Technical Field

The application relates to the technical field of electromagnetic simulation, in particular to a full-wave electromagnetic simulation method and system of an integrated circuit under a lossless frequency dispersion medium.

Background

Very large scale integrated circuits have a pronounced multi-scale structure with dimensions in the centimeter (10) range^-2m) -nm class (10)^-9m), the scale range is up to 7 orders of magnitude. On the other hand, signals transmitted by integrated circuits often have the characteristic of full-wave transmission, and the transmission frequency range of the signals covers from direct current to several GHz, which is especially problematic in integrated circuit applications of digital and mixed signal transmission. Therefore, the electromagnetic field analysis for the integrated circuit needs to be performed with a full-wave electromagnetic field analysis, and a full-wave electromagnetic field solver needs to be used to perform the full-wave electromagnetic field solution on the electromagnetic field of the integrated circuit.

However, investigation and test results show that, when a full-wave electromagnetic field is solved for the electromagnetic field problem of the integrated circuit, for the same integrated circuit model, when the test frequency is a high frequency above GHz, the solver can obtain an accurate field solution, and for the test frequency as low as MHz magnitude or lower than MHz magnitude, all solvers fail and cannot obtain an accurate field solution, and MHz and tens of MHz are just the working frequencies of many integrated circuits, so that the problem of solving the failure of the low-frequency electromagnetic field is urgently needed to be solved.

The reason for this is that the integrated circuit operates during operationThe electromagnetic field forming the integrated circuit is formed by superposing the contributions of the conduction current and the displacement current, and the contributions of the conduction current and the displacement current are equivalent at high frequency, so that the two parts are normally superposed; at low frequencies, however, the contribution of the conduction current is predominant, whereas the contribution of the displacement current is very low, the ratio of the contribution of the displacement current to the conduction current being even lower than the machine accuracy, e.g. at 10 for double-precision data storage^-16In order of magnitude, the ratio of the displacement current to the conduction current contribution at low frequencies will be comparable to or even lower than the machine accuracy, so that the ratio of the matrix elements representing the different contributions is of the order of magnitude comparable to or even lower than the machine accuracy. This fact causes that when the matrices representing different contributions are combined, the matrix elements representing the contributions of the displacement currents are completely covered by errors caused by machine precision during combination, so that the contributions of the displacement currents are negligible, but after the errors caused by machine precision are introduced, the contributions of the displacement currents are amplified by several orders of magnitude due to the errors, so that the contributions of the displacement currents are changed from negligible to resolvable, and the solution result is invalid, namely, the solved electromagnetic field is not an accurate field.

The existing method for solving the problem that the solver fails at low frequency is to combine the electromagnetic field solver based on constancy or quasi-constancy and the electromagnetic field solver based on high frequency to solve. When the test frequency is higher than a certain frequency, a high-frequency electromagnetic field solver is adopted for solving, and when the test frequency is lower than the certain frequency, a constant or quasi-constant electromagnetic field solver is adopted for solving, and then the calculation results of the two solvers are spliced.

However, this method is less accurate, firstly, because the constant or quasi-constant electromagnetic field solver involves a basic approximation, it is necessary to decouple the electric field E and the magnetic field H, forming a differential equation containing only the constant electric field or the constant magnetic field, which is correct only in the case of a strict dc field; secondly, how to set the specific frequency switched between the two solvers is unknown at present; finally, the field solved after being lower than a certain frequency is replaced by the constant field under direct current, so that the solved electromagnetic field is irrelevant to the frequency, the electromagnetic response curve of the integrated circuit of field splicing solved by the two solvers is obviously discontinuous, the electromagnetic response curve has obvious jump at the frequency point of switching of the two frequencies, and the response curve at a low frequency band is a straight line.

Therefore, in order to completely solve the problem that the existing solver is invalid in the calculation of the low-frequency field solution of the integrated circuit under the condition of no loss and dispersion media, it is necessary to accurately calculate and obtain the real solution of the full-wave maxwell equation set of the electric field E and the magnetic field H from direct current to high frequency of the integrated circuit.

Disclosure of Invention

Based on this, in order to solve the problem that the existing solver fails under the condition of low frequency of the integrated circuit under the condition of no loss and dispersion medium, and further obtain an electromagnetic field of the integrated circuit under the full frequency band including the low frequency, the application discloses the following technical scheme.

In one aspect, a full-wave electromagnetic simulation method of an integrated circuit under a lossless dispersive medium is provided, and comprises the following steps:

establishing an integrated circuit model according to layer information, layout information of each layer, via hole information and netlist information of the integrated circuit;

carrying out mesh subdivision on a parallel flat plate field of the integrated circuit by using the integrated circuit model so as to establish a matrix equation under a lossless frequency dispersion medium;

calculating a reference frequency point and a field solution thereof of the integrated circuit by an iteration method based on the matrix equation, wherein the integrated circuit field solution corresponding to the reference frequency point is an accurate solution;

obtaining a proportionality coefficient between a field solution of a reference frequency point and a low-frequency field solution, and obtaining a field solution under a frequency point to be solved based on the field solution of the reference frequency point and the proportionality coefficient, wherein a solution E (f) under the frequency point f to be solved is obtained through the following formula: e (f) ═ kE_refWhere k is a proportionality coefficient and k is a real number, E_refIs a reference frequency point f_refThe following accurate field solutions, the proportionality coefficient k being:

wherein the content of the first and second substances,

is E_refThe angular frequency omega of the electromagnetic wave is 2 pi f, b (omega) is an external excitation source of the whole finite element system, and K₂Is a dielectric constant related quality matrix of the whole finite element system;

and obtaining the electromagnetic response of the full frequency band based on the field solutions of all the frequency points to be solved.

In one possible embodiment, the establishing a matrix equation under a lossless dispersive medium includes:

establishing an electromagnetic field wave equation, and then acquiring a homogeneous equation corresponding to the electromagnetic field wave equation to obtain a functional of the homogeneous equation;

when the size of the electromagnetic field solving area reaches a set threshold value, setting the part of the functional, which is related to the electromagnetic wave, at the area boundary as 0, and dispersing the electromagnetic field solving area to obtain a discrete form of the functional:

taking a partial derivative of the functional in the discrete form and making the partial derivative be 0 to obtain a matrix equation under the following formula of the lossless dispersive medium: (K)₁-ω²K₂) E ═ b (ω), where E is the basic unit, E is the electric field, b^eAs an external excitation source for the basic cell E, E^eAn electric field vector formed by an electric field at an edge of the elementary cell e, L is the number of the elementary cells e whose entire electromagnetic field solution area is discrete,

is a stiffness matrix of the basic cell e,

is a dielectric constant dependent mass matrix of the medium of the elementary cell e,

is a basic unit eConductivity-dependent mass matrix of a medium, K₁Is the stiffness matrix of the whole finite element system, j is the imaginary unit.

In a possible implementation manner, the calculating the reference frequency point of the integrated circuit by an iterative method based on the matrix equation includes:

calculating a critical frequency point of the full-wave electromagnetic analysis of the integrated circuit based on the layout characteristic dimension of the integrated circuit and the machine precision of the simulation operation equipment, wherein the critical frequency point is a frequency point at which a solution result is credible to incredible when a matrix equation for simulating an electromagnetic field of the integrated circuit is solved;

and calculating the reference frequency point and the field solution thereof by an iterative calculation method based on the critical frequency point.

In one possible embodiment, the critical frequency points are calculated based on the following formula:

wherein f is₀Is a critical frequency point, a is the machine precision magnitude adopted in simulation operation, c is the wave velocity of electromagnetic waves in vacuum, and l is the size of a basic unit obtained by mesh division.

In a possible embodiment, the calculating the reference frequency point and the field solution thereof based on the critical frequency point by an iterative calculation method includes:

step A1, setting an iteration frequency lower limit F_minIs a critical frequency point f₀And setting an upper iteration frequency limit F_max＝Factor×f₀Wherein Factor is multiple of critical frequency point>1；

Step A2, converting the current angular frequency omega_curr＝2πF_minSubstituting the matrix equation, and solving the matrix equation to obtain omega_currField solution at angular frequency E_currCalculating the field solution E by a relative error calculation formula of the following formula_currRelative error res:

step A3, when the relative error res is less than or equal to epsilon₁Then, a reference frequency point f is obtained_ref＝ω_currPer 2 pi and field solution E thereof_ref＝E_currAnd ending; at the relative error res>ε₁Jumping to step a 4;

step A4, mixing omega_curr＝π(F_min+F_max) Substituting the matrix equation to obtain a new field solution, and calculating the relative error of the new field solution through the relative error calculation formula;

step A5, the relative error at the new field solution is less than ε₀When making F_max＝ω_currA/2 pi and jumping to step A4; the relative error at the new field solution is less than or equal to epsilon₁And is greater than or equal to epsilon₀Then, a reference frequency point f is obtained_ref＝ω_currAnd/2 pi and field solution thereof, and finishing; the relative error at the new field solution is greater than epsilon₁When making F_min＝ω_currA/2 pi and jumping to step A4;

wherein E is_currFor the current angular frequency omega_currField solution of ∈ at₀Is a preset lower error threshold value epsilon₁Is a preset upper limit of the error threshold, epsilon₀<ε₁。

On the other hand, the full-wave electromagnetic simulation method of the integrated circuit under the condition of no loss and frequency dispersion medium is also provided, and comprises the following steps:

the integrated circuit modeling module is used for establishing an integrated circuit model according to layer information, layout information of each layer, via hole information and netlist information of the integrated circuit;

the matrix equation building module is used for utilizing the integrated circuit model to carry out grid subdivision on a parallel flat plate field of the integrated circuit so as to build a matrix equation under a lossless frequency dispersion medium;

the reference frequency point calculation module is used for calculating the reference frequency point and the field solution of the integrated circuit through an iteration method based on the matrix equation, wherein the integrated circuit field solution corresponding to the reference frequency point is an accurate solution;

the to-be-solved field solution calculation module is used for obtaining a proportionality coefficient between a field solution of the reference frequency point and a low-frequency field solution, and obtaining a field solution under the to-be-solved frequency point based on the field solution of the reference frequency point and the proportionality coefficient, wherein a solution E (f) under the to-be-solved frequency point f is obtained according to the following formula: e (f) ═ kE_refWhere k is a proportionality coefficient and k is a real number, E_refIs a reference frequency point f_refThe following accurate field solutions, the proportionality coefficient k being:

wherein the content of the first and second substances,

and the electromagnetic response acquisition module is used for acquiring the electromagnetic response of the full frequency band based on the field solutions of all the frequency points to be solved.

In one possible embodiment, the matrix equation building module builds the matrix equation under the lossless dispersive medium by the following steps:

taking a partial derivative of the functional in the discrete form and making the partial derivative be 0 to obtain a matrix equation under the following formula of the lossless dispersive medium: (K)₁-ω²K₂) E ═ b (ω), where E is the basic unit, E is the electric field, b^eAs an external excitation source for the basic cell E, E^eElectric field shape of the edge of the basic cell eThe electric field vector L is the number of the elementary cells e whose whole electromagnetic field solving area is discrete,

is a stiffness matrix of the basic cell e,

a conductivity-dependent mass matrix, K, of the medium of the elementary cell e₁Is the stiffness matrix of the whole finite element system, j is the imaginary unit.

In a possible implementation manner, the reference frequency point calculating module calculates the reference frequency point of the integrated circuit by the following steps:

In a possible implementation manner, the reference frequency point calculating module calculates the reference frequency point and the field solution thereof, and specifically includes the following steps:

The application discloses a full-wave electromagnetic simulation method and a system of an integrated circuit under a lossless frequency dispersion medium, which can calculate a critical frequency point of full-wave analysis based on the characteristic dimension of the integrated circuit and the solving precision of a computer, and calculate a reference frequency point and a field solution under the reference frequency point by adopting an iteration method based on the critical frequency point, namely, by utilizing the property that the problem of generalized characteristic value is irrelevant to frequency, a reliable reference frequency point capable of obtaining an accurate solution and the field solution thereof reversely deduce the field solution under low frequency based on the reliable reference frequency point and the field solution thereof, calculate the field solution under low frequency by utilizing the expansion characteristic between the solution under low frequency to be solved and the solution under the reference frequency point, obtain the electromagnetic field of the integrated circuit under low frequency, solve the failure problem of the existing solver under the low frequency condition aiming at the integrated circuit, realize the accurate solving of the field solution of the low frequency point, and simultaneously solve the continuity problem of low frequency and high frequency response, the simulation result is more accurate, and the problem of discontinuous curves of response splicing of two solvers caused by difference of solution results at the crossed frequency points of the high frequency band and the low frequency band when different solvers are respectively adopted for the high frequency and the low frequency is avoided.

In addition, the method utilizes the generalized eigenvalue technology of the matrix to convert an equation set formed by the sparse matrix relevant to the frequency of the original problem into the generalized eigenvalue problem formed by the sparse matrix irrelevant to the frequency, and then corrects the eigenvalue of the error caused by the machine precision actually solved by a numerical calculation method according to the property of the eigenvalue of the sparse matrix, so that the error caused by the machine precision is avoided; and finally obtaining the field solution of the frequency point to be solved by utilizing the proportionality coefficient. Furthermore, the generalized eigenvalue problem formed by the finite element sparse matrix is not directly solved, the property of the generalized eigenvalue is utilized, the field solution under the low frequency is reversely deduced by taking the field solution of the reliable frequency as the reference, and the electromagnetic field of the integrated circuit under the whole frequency band from the low frequency to the high frequency is further obtained.

Drawings

The embodiments described below with reference to the drawings are exemplary and intended to be used for explaining and illustrating the present application and should not be construed as limiting the scope of the present application.

Fig. 1 is a schematic flow chart of an embodiment of a full-wave electromagnetic simulation method of an integrated circuit under a lossless dispersive medium disclosed in the present application.

Fig. 2 is a block diagram of an embodiment of an integrated circuit full-wave electromagnetic simulation system under a lossless dispersive medium disclosed in the present application.

Detailed Description

In order to make the implementation objects, technical solutions and advantages of the present application clearer, the technical solutions in the embodiments of the present application will be described in more detail below with reference to the drawings in the embodiments of the present application.

An embodiment of a full-wave electromagnetic simulation method of an integrated circuit under a lossless dispersive medium disclosed in the present application is described in detail below with reference to fig. 1.

As shown in fig. 1, the method disclosed in this embodiment includes the following steps 100 to 500.

And step 100, establishing an integrated circuit model according to layer information, layout information of each layer, via hole information and netlist information of the integrated circuit. The layer information of the integrated circuit, the layout information of each layer of the integrated circuit, the via hole information of the integrated circuit and the netlist information of the integrated circuit can be acquired from a chip design file when the electromagnetic simulation of the very large scale integrated circuit is carried out, and the layer information of the integrated circuit can comprise: the number of layers, the thickness of the layers, the medium material information of each medium layer, and the material and the thickness of each conducting layer.

After the layer information, the layout information of each layer, the via hole information and the netlist information are obtained, the integrated circuit can be modeled by using the information, and an integrated circuit model is obtained. The modeling process may include: converting layout information of the integrated circuit into discrete point information and polygon information containing layer information formed by the discrete points; converting the via hole information of the integrated circuit into a polygon prism (such as a hexagonal prism and a twelve-prism) connecting two layers; and converting the netlist information of the integrated circuit into external circuit information of the integrated circuit and topological structure information of the layout node.

In the order of the vertices of the polygon, the counterclockwise direction is positive to indicate that the inside of the polygon is filled with a conductive power supply layer, and the clockwise direction is negative to indicate that the inside of the polygon is filled with an insulating medium, that is, a hollowed polygon. The routing of the integrated circuit is converted into polygonal information in the layout of the integrated circuit.

The through holes of the integrated circuits are through holes for connecting the two layers of integrated circuits, and the two layers of integrated circuits are electrically connected at a proper position by introducing the through holes at the position. The via information of the integrated circuit may include: via position, via size, and layer number of the integrated circuit to which the upper and lower vertices of the via are connected.

The netlist information can provide circuit information and network information designed on each layer layout of the integrated circuit, wherein the circuit information provides a series of node information and connection relations thereof at different layers and different positions, the network information and the circuit information are in inclusion relation, and one network can contain a plurality of circuits. The above node information, circuit information and network information are given by unique names. The netlist information may also provide information such as devices and power supplies of the integrated circuit, and obtain topology information, i.e., the connection relationship of the external circuit.

And 200, performing mesh subdivision on a parallel flat plate field of the integrated circuit by using the integrated circuit model, and further establishing a matrix equation under the lossless frequency dispersion medium. In which, under the condition of an alternating electromagnetic field, electromagnetic waves in an integrated circuit propagate through a medium between metal layers of different layers, and the medium region between the metal layers of different layers is called a parallel flat plate field.

When electromagnetic simulation is performed on the integrated circuit by adopting an electromagnetic field numerical calculation method, a calculated field needs to be discretized based on the structure of the integrated circuit, namely mesh subdivision is performed, and then a discrete equation set is established based on the subdivided mesh for solving.

Because the layout information in the chip design file is the design file of the actual physical model, the original file may be a graphic file, and the polygons contained in the graphic file may have a large amount of redundant node information, if the simulation is directly performed on the layout polygon information in the chip design file, the input nodes of each polygon are all regarded as important nodes, so that all the input polygon nodes are inserted into the divided mesh nodes, which causes inaccurate calculation and wastes a large amount of computer resources.

Therefore, before mesh generation is carried out on the layout polygons of the multilayer integrated circuit, each layout polygon of each layer can be simplified without losing precision, so that the shape of the simplified polygons is not changed in a precision control range compared with the original input polygons, and then mesh generation is carried out on the modeling model subjected to layout simplification.

After the layout is simplified, the alignment of the polygons of the layout can be performed before the mesh generation, and then the mesh generation can be performed on the parallel flat plate field after the polygons are aligned. Specifically, if two layers of integrated circuit layouts designed have a layout polygon that is the same and completely aligned, a parallel flat field with the shape of the layout polygon will be formed between the two layers of layouts, however, since actually provided layout polygon information may generate some errors due to the processes of file format conversion and the like, the errors may cause that the finally introduced layout polygons of the two layers are not completely the same or not completely aligned, and compared with the polygon size, there may be an error of 1% or 0.1%. If the parallel flat plate field is directly identified based on the input polygon information, a plurality of parallel flat plate field fragments are generated, the fragment size is very small, the introduction of the fragment may also cause very dense mesh subdivision with poor quality, the result is wrong while the computing resource is wasted, and therefore, the layout polygons of the multilayer integrated circuit can be aligned without losing precision so as to obtain the parallel flat plate field with higher accuracy and compactness.

After the alignment of the multilayer integrated circuit layout is realized, the parallel flat plate field of the multilayer integrated circuit layout can be identified based on the layer information of the polygon of the integrated circuit layout.

In one embodiment, the matrix equation under the lossless dispersive medium is established by the following steps 210 to 250.

Step 210, establishing an electromagnetic field wave equation, and then obtaining a homogeneous equation corresponding to the electromagnetic field wave equation to obtain a functional of the homogeneous equation. Specifically, step 210 first establishes a wave equation based on the electric field E by using maxwell's equation to obtain an electromagnetic field wave equation in the following formula (1):

in the formula (1), the reaction mixture is,

is the rotation operator, mu_rIn order to obtain a relative magnetic permeability of the medium,

is the electric field vector, ω is the electromagnetic angular frequency (in rad/s), c is the wave velocity of the electromagnetic wave in vacuum, c is 3 × 10⁸m/s，ε_rIs the relative dielectric constant of the medium, j is an imaginary unit, j²＝-1，μ₀Is the permeability of a vacuum medium, mu₀＝4π×10^-7H/m, σ is the conductivity of the medium (in S/m),

current density (in A/m) for applied excitation²)。

Then, a homogeneous equation corresponding to the electromagnetic field wave equation is obtained through a variational principle, and the functional is as follows (2):

in the formula (2), V is an electromagnetic field resolving area, S is a plane surrounded by the electromagnetic field resolving area, and n is a normal vector in which any point of the plane S points outward.

And step 220, when the size of the electromagnetic field solving area reaches a set threshold value, setting the part of the functional, which is related to the electromagnetic wave, on the area boundary as 0, and dispersing the electromagnetic field solving area to obtain a discrete form of the functional. Specifically, when the electromagnetic field solution area is large enough in step 220, so that the electromagnetic wave is attenuated to approximately 0 at the boundary of the area, the functional function can be simplified as the following formula (3):

in the equation (3), the judgment method of whether the electromagnetic field solution area is large enough may be to obtain a ratio between a minimum distance from a boundary of the solution area to a source (i.e., a multilayer integrated circuit board) generating the electromagnetic wave and a wavelength of the electromagnetic wave, compare the ratio with a preset multiple, and determine that the electromagnetic field solution area is "large enough" if the ratio exceeds the preset multiple, for example, the preset multiple is 10 times, and determine that the solution area is large enough if the ratio is greater than 10.

After the functional simplification is performed to obtain the formula (3), the electromagnetic field solving area is discretized by using a sufficiently small basic unit formed by the mesh subdivision, the basic unit may be a tetrahedron, a triangular prism, a hexahedron or the like, and the electric field of any point in each discrete unit is represented by interpolating a basis function and an electric field of an edge or a surface element, as shown in the following formula (4):

in the formula (4), the reaction mixture is,

is the electric field of any point in the basic unit e, M is the number of interpolation basis functions,

the ith interpolation basis function of the basic unit E, E^eAn electric field vector formed by an electric field at an edge of the basic cell e,

the electric field value of the corresponding edge or surface element of the ith basis function of the basic unit e,

for M interpolation basis functions on an edge or bin of a basic cell

Of size M × 1, { E^eIs M interpolation basis functions on the edge or surface element of the basic unit

Corresponding electric field value

The size of the matrix form of (1) is M × 1, and T represents a transpose of the matrix.

The determination method of whether the basic unit is small enough needs to be determined according to the following two conditions, and if the following two conditions are satisfied simultaneously, the size of the basic unit is considered small enough:

1. judging that the size of the nearby grid is not larger than the characteristic size of the integrated circuit according to the size relation between the local size of the basic unit and the characteristic size of the integrated circuit;

2. the relationship between the maximum size of the basic unit and the minimum wavelength of the electromagnetic wave of the integrated circuit to be simulated, wherein the minimum wavelength of the electromagnetic wave of the integrated circuit is not less than a preset multiple, such as 10 times, of the maximum size of the basic unit.

Substituting the interpolation function into the simplified equation (3) to obtain a discrete functional expressed by the following equation (5):

in the formula (5), the reaction mixture is,

wherein the content of the first and second substances,

is the relative permeability of the medium in the region of the basic element e,

is the relative permittivity of the medium in the region of the elementary cell e,

p-th interpolation basis function of basic unit e, V^eIs an integral of the basic unit e, and L is the number of the basic units e in which the entire electromagnetic field resolving area is dispersed.

Is a stiffness matrix of the basic cell e,

is the conductivity-dependent mass matrix of the medium of the elementary cell e.

When an external stimulus is present, the formula (5) becomes the following formula (9):

in the formula (9), b^eIs an external excitation source of the basic unit e,

in the formula (10), the compound represented by the formula (10),

an external excitation current source, V, for the basic cell e^eIs an integral of the basic cell e. And step 230, taking a partial derivative of the discrete functional and setting the partial derivative to be 0 to obtain a matrix equation of the vector finite element under the lossless and frequency dispersion medium. Specifically, in step 230, according to the energy minimization principle, the solution corresponding to the electromagnetic field wave equation (1) is to solve the electric field E corresponding to the extremum value of the functional expressed by the above equation (5), so that the following equation (11) can be obtained by taking the partial derivative of the equation (5) and making the partial derivative be 0:

(K₁-ω²K₂+jωK₃)E＝0 (11)；

in the formula (11), K₁Is the stiffness matrix of the entire finite element system, K₂Is a dielectric constant-dependent mass matrix, K, of the entire finite element system₃Is the conductivity-related mass matrix of the entire finite element system.

When an external excitation source is present, a matrix equation is constructed therefrom:

(K₁-ω²K₂+jωK₃)E＝-b(ω)

wherein b (omega) is an external excitation source of the whole finite element system and is a function of the angular frequency omega of the electromagnetic wave.

Because the metal layer of the integrated circuit is an ideal conductor under the condition of no loss medium, the conductivity of the medium layer is 0; and in the presence of a dispersive medium, the dielectric constant and permeability of the dielectric layer of the integrated circuit may change due to a change in frequency, i.e., K₁、K₂And K₃The element(s) of (c) will change due to the change in frequency. Therefore, for a medium whose type is lossless and has dispersion, the stiffness matrix corresponding to the type of the medium is K in equation (12)₁The corresponding quality matrix is K in the formula (12)₂And since the conductivity is 0 at this time, K₃0. Thus, the resulting matrix equation is actually the following equation (12):

(K₁-ω²K₂)E＝-b(ω) (12)。

and step 300, calculating the reference frequency point and the field solution of the integrated circuit by an iteration method based on the matrix equation. The reference frequency point means that the field solution of the integrated circuit calculated at the reference frequency point is definitely accurate, the direct current eigenmode of the electromagnetic wave at the reference frequency point plays a main role, and the contribution of the high-order eigenmode can be ignored.

In one embodiment, step 300 includes step 310 and step 320.

And 310, calculating critical frequency points of the full-wave electromagnetic analysis of the integrated circuit based on the layout characteristic dimension of the integrated circuit and the machine precision of the simulation operation equipment. The critical frequency point is a frequency point which is determined by solving a matrix equation of electromagnetic field simulation of the integrated circuit from high frequency to low frequency in a frequency domain, wherein the error of a solving result changes from small (result credibility) to large, when the error is large to a certain degree, the solving result is not credible, and the frequency point which is determined by solving the result credible to the untrustworthy solving result is the critical frequency point.

The characteristic dimension of the layout refers to the minimum dimension and the maximum dimension of the layout, the minimum dimension refers to the dimension of the minimum functional unit of the layout, the minimum dimension may be the diameter of the minimum via hole in the layout, the width of the thinnest routing, the width of the minimum gap between the routing, or the minimum layer thickness of the multilayer integrated circuit layout, and the maximum dimension refers to the overall dimension of the whole layout plane.

In one embodiment, step 310 includes step 311 and step 312.

And 311, obtaining the magnitude ratio related to the size among different matrix elements based on the range of the layout characteristic size.

In the most advanced very large scale integrated circuits at present, the feature sizes of the layouts of different positions of the integrated circuit have dimensions in the centimeter scale (10)^-2m) -nm class (10)^-9m) if tetrahedron is used for calculation of multilayer VLSIThe field is discretized, and the smallest dimension of the discrete tetrahedral mesh is on the order of nanometers.

In addition, in the matrix expressions in expressions (6) to (8), all the basic units

Is proportional to 1/l because of the interpolation function of the basic unit

Normalization has been performed, i is the size of the elementary cell obtained by mesh division, and the volume of the elementary cell e and i³Is in direct proportion, and

and

are all constants independent of layout feature size, so K in equation (6)₁Is of the order of O (l), O (. smallcircle.) represents a comparable order of magnitude, O (l) represents a comparable order of magnitude to l, K in equation (7)₂Of the order of O (c)^-2l³) K in the formula (8)₃Of the order of O (mu)₀σl³) Thus, the matrix K in the formula (13) can be obtained₁、K₂Norm ratio of (2):

the formula (13) can embody a matrix K₁And K₂There is inevitably a difference in order of magnitude between the elements of (a) and it can also be known from equation (13) that the discrete cell size/is the only change matrix K₁And K₂And the smaller the size l, the smaller the matrix K₁And K₂The larger the element contrast difference, and the smaller the difference. It can be understood that the feature size of the layout determines the distribution of the mesh subdivision sizes, for example, if the layout is a multi-scale structure and the feature size is from centimeter level with the maximum size to nanometer level with the minimum size, the mesh subdivision sizes are determinedThe minimum size of the tetrahedral unit is nano-scale and is distributed at the position of the small-size layout; the maximum size may be in the centimeter level and is distributed at the position without the small-size layout. Therefore, step 200 is to divide the grid according to the layout information of the integrated circuit, the minimum feature size of the integrated circuit determines the minimum size of the grid unit, and the magnitude difference between the finite element system matrixes in this step is based on the minimum size of the grid unit.

Suppose that in the integrated circuit currently implementing the method, the dimension range of the characteristic dimension of the layout of the integrated circuit at different positions is centimeter level (10)^-2m) -nm class (10)^-9m), if tetrahedrons are used as basic units for meshing, i.e. discretizing, the computational field of a multi-layer VLSI, then since c is of the order of 10⁸L is of the order of 10^-9Then K is₂Norm and K of₁The ratio of the norm of (a) will also be as low as 10^-34Of the order of magnitude, i.e. a magnitude ratio of 10³⁴。

And step 312, obtaining the machine precision adopted in the simulation operation, and calculating the critical frequency point of the integrated circuit based on the machine precision and the magnitude ratio.

As can be seen from the stiffness matrix equation shown in the formula (12) and the previous discussion, all solvers fail at low frequencies, and the failure of the electromagnetic field simulation solver of the integrated circuit is caused by limited machine precision of a computer, because the matrix K in the matrix equation (12) formed by vector finite element calculation₁And K₂Is different by an order of magnitude which is inevitable when the frequency is low enough to make the frequency-dependent matrix ω in equation (12)²K₂Can lead to failure of the solver when the contribution of (a) is lost due to limited machine accuracy, since the left-hand expression of the matrix equation is now approximately equal to K₁When this happens, the solution of equation (12) solved by the solver is completely wrong, since K is now K₁Is a singular matrix.

Assuming that double-precision type data is adopted for calculation at present, the precision of a machine adopted in simulation operation is 10^-16When using double precisionData type execution ω²K₂And K₁In the subtraction of (2), if the matrix K is₁And ω²K₂The phase difference is greater than 10¹⁶The simulation computing device (e.g. computer) will directly use the matrix ω²K₂The solution of the matrix equation (12) for electromagnetic field simulation of the integrated circuit must fail to be considered zero. Even if the frequency is MHz, omega²K₂Is also compared with K₁To be smaller by 34-2 x 6-22 orders of magnitude, so when ω is executed²K₂And K₁In the subtraction, the simulation computing device (e.g. computer) will directly apply the matrix ω²K₂Considered as zero.

Therefore, the failure frequency (i.e., the critical frequency point) in the simulation of the integrated circuit needs to be calculated based on the machine precision and the magnitude ratio, and the frequency satisfying the following formula (14) is the critical frequency point:

wherein f is the frequency point to be solved which is not higher than the critical frequency point. From this formula, one can obtain:

therefore, the formula for calculating the critical frequency point is formula (15):

Assuming that the machine precision magnitude a is 16, the formula for determining whether the critical frequency point condition is satisfied is specifically: II | K₁‖/‖(2πf)²K₂‖>10¹⁶The critical frequency point is

Due to the ratio II K₁‖/‖K₂Is difficult to be directly and accurately calculated, so that O (·) operation is adopted to obtain the ratio II K₁‖/‖K₂Of the order of magnitude II, by O (c)²/l²) Instead of | K₁‖/‖K₂II this ratio is followed, but because of O (c)²/l²) Obtained is only c²/l²Rather than an exact value, c is the exact value²/l²Thus making it possible to

In one embodiment, l may be the smallest possible dimension, i.e., l ═ l_min，l_minIs the smallest dimension in the range of values for l. The value of l can be various, but for the most advanced super large scale integrated circuit, the minimum characteristic dimension of the layer structure and the layout reaches the nanometer level (10)^-9m), the size of the discrete tetrahedrons is also in the nanometer range, while if a grid size in the nanometer range can be achieved, a grid size higher than that of the nanometer process can be achieved, so l can be taken_min＝10^-9m, i.e. l 10^-9m, critical frequency point f at this moment₀160MHz, that is, the results obtained when solving the integrated circuit electromagnetic field simulation matrix equations at frequencies below 160MHz must be inaccurate.

Steps 411 and 412 can calculate critical frequency points of full-wave analysis based on the feature size of the integrated circuit and the solution accuracy of the computer, so that the electromagnetic field simulation field solution at which frequency is obtained is definitely inaccurate. After the critical frequency point is obtained, the reference frequency point can be calculated according to the critical frequency point.

And step 320, calculating the reference frequency point and the field solution thereof by an iterative calculation method based on the critical frequency point.

In one embodiment, step 320 includes steps 321 through 325.

Step 321, setting an iteration frequency lower limit F_minIs a critical frequency point f₀And setting an upper iteration frequency limit F_max＝Factor×f₀Wherein Factor is the multiple of the critical frequency point. Factor>1, and Factor may be set to 10.

Step 322, converting the current angular frequency ω_curr＝2πF_minSubstituting the matrix equation of formula (12), and solving the matrix equation to obtain omega_currField solution at angular frequency E_curr. Calculating the field solution E by a relative error calculation formula of_currTo perform checking, wherein the relative error res is used to evaluate whether the field solution of the frequency point is accurate:

in the formula, the numerator is the residual error, the denominator is the source term, the modulus of the two is the relative error of the two, E_currFor the current angular frequency omega_currAnd (4) field solution.

Step 323, when the relative error res is less than or equal to epsilon₁When the frequency is the lowest, the description shows that the precision requirement is met when the iteration of the current round is started, and the reference frequency point f is obtained_ref＝ω_currPer 2 pi and field solution E thereof_ref＝E_currAnd ending; at the relative error res>ε₁It jumps to step 324.

Step 324, let ω be_curr＝π(F_min+F_max) Substituting the matrix equation of formula (12), solving the matrix equation to obtain a new field solution E_curr'and calculating a relative error res' of the new field solution by the relative error calculation formula.

Step 325, relative error res 'at the new field solution'<ε₀When making F_max＝ω_currA/2 pi and go to step 324; the relative error res' of the new field solution is less than or equal to epsilon₁And res' is more than or equal to epsilon₀Then, a reference frequency point f is obtained_ref＝ω_currAnd/2 pi and field solution thereof, and finishing; relative error in the new field solutionPoor res'>ε₁When making F_min＝ω_currAnd/2 pi and jumps to step 324.

Wherein epsilon₁Is a preset upper limit of an error threshold value and satisfies epsilon₁>ε₀，ε₀Is a preset lower error threshold value epsilon₀Can take the value of 10^-5，ε₁Can take the value of 5 multiplied by 10^-5。

This embodiment sets up ∈₀And ε₁In step 325, the two thresholds are ended as long as the error is between the two thresholds, so that the convergence rate of the bisection iteration can be accelerated.

The error threshold is a quantitative index of whether the relative error meets the standard, if the error does not meet the standard, the frequency is still invalid (the field solution is inaccurate), if the error meets the standard, the frequency is reliable (the field solution is accurate), but the reliable frequency is not necessarily a reference frequency point, because the frequency may not be the lowest reliable frequency, the lowest reliable frequency can be obtained through iteration and used as the reference frequency point.

Steps 321 to 325 are sequentially executed steps, and are executed according to the step sequence number specified in the jump only when the jump exists in the step content, for example, there is an iteration process through the jump in step 325, as long as the relative error calculated in each iteration is less than epsilon₀Iterations occur and the field solution and relative error are recalculated, and so long as the relative error calculated for each iteration is greater than ε₁Iteration will also occur and the field solution and relative error will be recalculated until the relative error calculated after iteration is at ε₀And ε₁So as to realize the calculation of the reference frequency point through the steps 321-325.

Step 321 plus 325 can calculate the reference frequency point and the field solution at the reference frequency point by using an iterative method based on the critical frequency point, so that the electromagnetic field simulation field solution at which frequency is obtained is definitely accurate.

In the case of a lossless dielectric, the metal layer of the integrated circuit is an ideal conductor, the conductivity thereof is 0, and the matrix equation at the left end of equation (12) is equation (16):

K(ω)＝K₁-ω²K₂ (16)；

when the size of discrete tetrahedron is nano-scale, the matrix K₁And K₂By 10 of the element of³⁴Of the order of (1), so adding or subtracting two matrices directly at low frequencies results in a matrix omega²K₂Considered as zero. In order to solve the problem fundamentally, a generalized eigenvalue decomposition method is provided, and an eigenvector of the matrix is extracted, wherein the eigenvector is only related to the scale characteristic of the integrated circuit and is not related to the simulation frequency. For the matrix K₁And K₂The generalized eigenvalue problem of frequency independence is shown in formula (17):

K₁x＝λK₂x (17)；

in the formula (17), λ is an eigenvalue, x is an eigenvector, and the formula (17) has N sets of eigenvalues and solution of the eigenvector in common, that is, (λ)₁,λ₂,…,λ_N) And (x)₁,x₂,…,x_N). Due to K₁Is symmetrically semi-positive, K₂Are symmetrically positive, so that the eigenvalues λ are non-negative and real, while x and K are₁、K₂Are orthogonal.

Let Φ equal to (x)₁,x₂,…,x_N) Then obtaining formula (18) and formula (19):

Φ^TK₁Φ＝Λ (18)；

Φ^TK₂Φ＝I (19)；

where Φ is a matrix formed by all eigenvectors, where Λ is a diagonal matrix made up of eigenvalues, and I is an identity matrix.

The inverse matrix of equation (20) is obtained from equations (18) and (19):

K(ω)^-1＝Φ(Λ-ω²I)^-1Φ^T (20)；

since the eigenvalues in the formula (17) can be divided into two groups, one group is related to the physical direct current mode and the non-physical direct current mode generated by the null space of the solution, the eigenvalue is theoretically zero, and the corresponding eigenvalue and eigenvector are respectively marked as Λ₀And phi₀(ii) a The other group is related to the non-zero resonance frequency of the three-dimensional structure of the integrated circuit, is a high-order eigenmode, the eigenvalue of the eigenmode is larger than zero, and the corresponding eigenvalue and eigenvector are respectively marked as lambda_highAnd phi_highAnd due to Λ₀By changing equation (20) to equation (21) when 0:

in the formula (21), phi₀The feature vector corresponding to the zero eigenvalue, called the null space vector, Φ_highFor eigenvectors corresponding to non-zero eigenvalues, Λ_highIs a diagonal matrix of non-zero eigenvalues, phi_highAnd Λ_highI.e. higher order modes associated with non-zero resonant frequencies.

To the right of equation (21), all other expressions are frequency independent except for the expression explicitly containing ω, and the frequency dependence of the solution of the matrix equation shown in equation (12) is clearly expressed. With such a frequency continuous function, a field solution of any low frequency of the integrated circuit from a high frequency to a full wave range including a direct current can be strictly obtained without the problem that the solver fails at the low frequency. It can therefore be seen that, given an arbitrary frequency ω, the field solution is a superposition of several three-dimensional eigenmodes. For the direct current eigenmodes, i.e., eigenvectors corresponding to zero eigenvalues, their weights in the field solution are (1/ω)²) Is in direct proportion; for higher order eigenmodes, their weight in the field solution for the ith eigenmode is proportional to 1/(λ ^ k)_i-ω²) Wherein λ is_iIs the characteristic value corresponding to the ith eigenmode, and i is the serial number of the eigenmode.

At low frequencies, when the higher eigenmodes are weighted significantly less than the dc eigenmodes, i.e., [1/(λ ^ k)_i-ω²)]＜＜1/ω²Time, high order eigenmode pairThe contribution of the field solution is negligible. Therefore, expression (21) may become expression (22):

from this, the solution of the integrated circuit full-wave electromagnetic simulation matrix equation shown in equation (12) below without loss and with dispersion medium is equation (23):

and if the eigenvectors corresponding to all the zero eigenvalues of the generalized eigenvalue problem shown in the equation (17) are solved, the field under the low frequency can be accurately solved.

For a three-dimensional multi-layer VLSI, although the number of physical DC modes may be small, the null space mixes the physical DC mode and the non-physical DC mode, and the linear combination of the two modes exists in the same null space. Therefore, the physical dc mode and the non-physical dc mode cannot be distinguished from only the null space vector.

Furthermore, the size of the nullspaces cannot be reduced by discarding subsets of the nullspaces in order to increase the computation speed during the solution process, because the subsets are linearly independent of each other and each of them is essential for establishing a complete null, so that the remainder after discarding any subset of the null is incomplete, or, given an excitation vector, it can have a projection on all the null vectors, so that each null contributes to the field solution.

However, for the most advanced current very large scale integrated circuits, which have a multi-scale structure with dimensions from centimeter to nanometer, and thus the grid cells for which field division is calculated are in the order of tens of millions or even hundreds of millions, the matrix equation shown in equation (15) thus formed is in the order of tens of millions to hundreds of millions, if all the null-space vectors are solved and stored, the calculation cost is high, and a great memory storage is required, because the null space of the matrix equation is large at this time and grows linearly with the matrix size. Therefore, how to deal with the increased null space becomes the key to quickly and accurately solve the response of the very large scale integrated circuit at low frequency.

The method for solving the problem is proposed by utilizing the characteristic that the generalized eigenvalue problem shown in the formula (17) has a plurality of zero eigenvalues, and all the zero space vectors share the same zero eigenvalue due to the plurality of identical zero eigenvalues, although the eigenvectors of the zero space vectors are completely different.

Based on this, the present application uses the vector at the right end of equation (12) (the excitation vector) to reduce the dimension of the space where the field solution is located, so the right-hand source vector is always known, i.e., for a given right-hand term b (ω) solving equation (12), all null-space vectors are effectively combined together to form the field solution of the matrix equation at low frequencies, as shown in equation (23). Further, the contribution of all null space vectors representing low frequencies can be represented by a vector E₀To represent that the vector covers all low frequency solutions, as shown in equation (24):

due to the adoption of the reference frequency point f_refThe resulting field solution vector is represented by E₀The solution space is formed, so that the solution of the representation field at the reference frequency point needs to have an expression form as shown in the formula (24), that is, f_refThe solution of the field should be dominated by the DC eigenmode and satisfy [1/(λ)_i-ω²)]＜＜1/ω²I.e. by

The second term at the right end of this equation (21) is negligible, i.e., the contribution of the higher order eigenmodes is negligible. λ here_iIs any non-zero eigenvalue, which satisfies lambda_min≤λ_i≤λ_maxWherein λ is_minIs the minimum non-zero eigenvalue, λ_maxIs the largest non-zero eigenvalue, λ_minAnd λ_maxCorresponding to the lowest resonant frequency f of the three-dimensional multi-scale structure integrated circuit_minAnd the highest resonance frequency f_maxLowest resonance frequency f_minMaximum resonance frequency f corresponding to maximum size of VLSI_maxThen corresponding to the minimum mesh size of the VLSI after mesh generation, therefore f_maxAnd f_minThe ratio of the maximum dimension to the minimum dimension of the subdivision grid is the ratio of the maximum dimension to the minimum dimension of the VLSI, and the ratio of the maximum non-zero eigenvalue to the minimum non-zero eigenvalue is the square of the frequency or dimension ratio.

Let λ_maxAnd λ_minThe order of magnitude of the ratio is m, the critical frequency point f₀The corresponding characteristic value is recorded as lambda₀Then the selected reference frequency point f_refIt should satisfy: f. of₀<f_ref＜＜f_minSince the resonance frequency corresponds to a value whose generalized eigenvalue is not zero, λ_iSince the non-zero eigenvalue has a large value and a high resonance frequency, even the lowest resonance frequency is much larger than the reference frequency point, and the problem of the generalized eigenvalue shown in the formula (17) is not actually solved, the lowest resonance frequency f is set to be the lowest resonance frequency f_minIs an unknown quantity.

Thus, the reference frequency point f is determined_refThe principle of (a) is to ensure that the solver does not fail under the reference frequency point and ensure that the reference frequency point f_refAs far as possible below the lowest resonance frequency f_minThis ensures that the second term to the right of equation (21) is negligible, i.e., the contribution of the higher order eigenmodes is negligible.

Step 400, for the frequency point to be solved which is lower than the reference frequency point, according to the property of the finite element system matrix eigenvalue described above, the proportionality coefficient between the field solution of the reference frequency point and the low-frequency field solution is obtained, and the field solution under the frequency point to be solved is obtained based on the field solution of the reference frequency point and the proportionality coefficient.

By substituting formula (11) for formula (24), the following formula (25) can be obtained, and it can be understood that the reference frequency point is not only a frequency point at which the field solution can be accurately obtained by the finite element method, but also belongs to the low frequency range, that is, the field solution E (ω) at the reference frequency point satisfies formula (25):

however, when the dielectric constant and the permeability of the dielectric layers forming the multi-layered LSI have dispersion characteristics, the matrix K₁And K₂Is no longer frequency independent, i.e. the eigenvector Φ corresponding to the zero eigenvalue of the above equation₀The frequency is not related to the frequency any more, and at the moment, the solution of the frequency to be solved cannot be directly solved according to the solution under the reference frequency point. However, the dimension of the space in which the field solution is located can still be reduced by using the vector of the right-hand term b (ω) (excitation vector), i.e. for a given right-hand term, all low-frequency solutions can be covered by one vector, i.e. the contribution of all null-space vectors representing low frequencies can be covered by one vector E₀Thus, in the case of dielectric constant and permeability of dielectric layers forming a multi-layer LSI circuit having dispersion characteristics, the low frequency solution can still use a vector E representing the contribution of all the zero-space vectors of the low frequency₀Thus, the low frequency solution has a scaling relationship with the solution of the reference frequency point, and in one embodiment, the scaling relationship is shown in equation (26):

E(f)＝kE_ref (26)；

where k is a proportionality coefficient and k is a real number, E_refIs a reference frequency point f_refThe accurate field solution below.

Formula (27) can be obtained by substituting formula (26) for formula (12):

k(K₁-ω²K₂)E_ref＝b(ω) (27)；

since the matrix K is in the nanometer range of discrete tetrahedral size in the very advanced VLSI circuits₁And K₂The phase difference is of the order of 10³⁴Directly performing ω at low frequency²K₂And K₁In the subtraction of (2), the computer directly uses the matrix omega²K₂Considered to be zero, the value of k is still not solved accurately. But note that this is due toWhen E is_refIs representative of a low frequency solution, thus E_refThe method is the combination of all zero space vectors representing low frequency, and the corresponding characteristic value of the zero space vector is 0, namely K is obtained₁Φ₀0, to give K₁E_ref0. Will K₁E_refSubstitution of formula (27) with 0 yields formula (28):

k(ω²K₂)E_ref＝b(ω) (28)；

both ends of the pair of equation (28) are multiplied by non-zero vectors

To give formula (29):

the formula for calculating the scaling factor k can be obtained from the formula (29), that is, the formula (30):

as can be seen from equation (30), at low frequencies, [1/(λ) ]_i-ω²)]＜＜1/ω²When the time is short, namely the contribution of the high-order eigenmode to the field solution can be ignored, the field solution of the matrix equation under the lossless dispersive medium can pass through the reference frequency point f_refIs accurate solution E_refCalculated based on the reference frequency point f_refIs accurate solution E_refCalculating the frequency f (f) to be solved<f_ref) The following calculation formula of the field solution is formula (31):

the above formula shows that, in order to calculate the electromagnetic response of the very large scale integrated circuit under the low frequency, the reference frequency point and the accurate field solution under the reference frequency point can be calculated first, and then the electromagnetic response of the very large scale integrated circuit under the reference frequency point can be obtained by analyzing the formula (31).

The reference frequency point f_refField solution E of_refThe field solution is obtained in step 320, and in the process of calculating the reference frequency point in each iteration in step 320, the corresponding field solution is calculated in each iteration by calculating one frequency point, for example, an angular frequency ω is obtained in step 324_currField solution E of_currIf the angular frequency is determined to be the reference frequency point in step 325, the corresponding field solution E is determined_curr' is the field solution E of the reference frequency point_ref。

Therefore, the scheme of the integrated circuit full-wave electromagnetic simulation under the condition of no loss and dispersive medium can be obtained: firstly, the reference frequency point of the integrated circuit to be simulated is calculated through the step 321-325, and for the frequency band higher than the reference frequency point, the solution can be directly carried out by adopting a full three-dimensional electromagnetic field numerical calculation method; for the frequency bands lower than the reference frequency point, the solution of the frequency point to be solved is obtained through the field solution of the reference frequency point and the proportionality coefficient between the field solution and the low frequency solution, so that the problem of continuity of low-frequency and high-frequency responses is solved, the simulation result is more accurate, and the problem of discontinuous response splicing curves of two solvers caused by the difference of solution results at the frequency points where the high-frequency band and the low-frequency band intersect when different solvers are respectively adopted for the high frequency and the low frequency is avoided.

And 500, acquiring the full-band electromagnetic response required by the user based on the field solutions of all the frequency points to be solved.

After the field solutions of the frequency points in the low-frequency band and the high-frequency band are obtained, post-processing operation can be performed based on the discrete field quantity and the calculation information required by the user, and the method specifically includes at least one of the following implementation items: 1. calculating to obtain at least one distribution item of potential distribution, current distribution, power consumption distribution, loss distribution and heat distribution of different layers based on the calculated discrete field quantity, and performing graph drawing on the obtained distribution item; 2. further calculating user-provided port parameters based on the calculated discrete field quantities, including full-band frequency response characteristics of at least one of multi-port S parameters, multi-port impedance matrices, and equivalent circuit model parameters of the integrated circuit; 3. further calculating electromagnetic radiation of the integrated circuit and/or electromagnetic interference of a location of the integrated circuit based on the calculated discrete field quantities; 4. the behavioral characteristics of the components of the integrated circuit, which may be voltage-current characteristics or other characteristics/characteristics, are further calculated based on the calculated discrete field quantities, thereby further extracting the IBIS model of the integrated circuit. After the calculation of the term to be calculated, the full-wave electromagnetic simulation of the integrated circuit is completed aiming at the lossless dispersive medium.

At present, when a schematic diagram of a super-large scale circuit is designed and is realized through an integrated circuit, a lead is formed into a copper-clad integrated circuit layout through a wire, the conductivity of the copper is a finite value, the wire has a certain resistance, the thinner and longer the wire is, the larger the resistance generated on the wire is, and the larger the voltage drop caused by the wire is. Because the signal transmission of the integrated circuit actually transmits 0 and 1 signals, the 0 and 1 signals are realized by the jump of high and low levels, and the integrated circuit component identifies that the high and low levels are judged by the threshold value of the level. Since the high and low levels of transmission are applied to the reference voltage of the transmission line, and the voltage drop on the trace is an important factor causing the instability of the reference voltage, the instability of the reference voltage directly causes errors of the transmission result.

At present, when a schematic diagram of a super-large scale circuit is designed, influence of electromagnetic interference between lines on signal transmission is not considered, but when a high-frequency alternating current is transmitted through a physical transmission line, electromagnetic waves exist around the transmission line, and the electromagnetic waves influence the transmission line around the transmission line in an electromagnetic induction mode to form electromagnetic interference. If the influence of the electromagnetic interference is equivalent to or even larger than the signal transmitted by the transmission line, the interference is superimposed on the signal of the transmission line and transmitted as a transmission signal, thereby changing the originally transmitted signal and destroying the work of the integrated circuit. Therefore, calculating the electromagnetic field and the voltage and current distribution at all positions on the integrated circuit layout is an important means for analyzing the voltage drop of the integrated circuit layout through simulation calculation of the electromagnetic field.

In the aspect of performing full-wave electromagnetic simulation, as the design of large-scale integrated circuits is developed to higher frequency and the circuit complexity is increased, the simulation of high-frequency electromagnetic field usually ignores the error caused by high-order propagation mode.

In addition, in the traditional mode, when an equivalent circuit method is adopted for analysis, equivalent elements such as capacitance and inductance do not consider the change of the elements along with the frequency, so that errors are caused. For example, for a microstrip line or an interconnection line structure, due to the characteristics of complex intersection, step, bend, open circuit, slot, etc. of the microstrip line or the stripline, signal delay, distortion, reflection, etc. effects, etc. are generated on the microstrip line or the interconnection line, and crosstalk between adjacent lines occurs, when electromagnetic waves are transmitted in multiple modes. For this reason, full-wave electromagnetic simulation techniques are typically required to analyze the integrated circuit to obtain accurate S-parameters in discontinuous mode, since full-wave analysis takes into account all possible field components and boundary conditions.

An embodiment of the integrated circuit full-wave electromagnetic simulation system under the lossless dispersive medium disclosed in the present application is described in detail below with reference to fig. 2. The embodiment is a system for implementing the embodiment of the integrated circuit full-wave electromagnetic simulation method under the lossless dispersive medium.

As shown in fig. 2, the system disclosed in this embodiment mainly includes: the device comprises an integrated circuit modeling module, a mesh generation module, a matrix equation construction module, a reference frequency point calculation module, a to-be-solved field solution calculation module and an electromagnetic response acquisition module.

wherein the content of the first and second substances,

In one embodiment, the matrix equation building block builds the matrix equation under the lossless dispersive medium by:

taking a partial derivative of the functional in the discrete form and making the partial derivative be 0 to obtain a matrix equation under the following formula of the lossless dispersive medium: (K)₁-ω²K₂) E ═ b (ω), where E is the basic unit, E is the electric field, b^eAs an external excitation source for the basic cell E, E^eIs based onThe electric field vector formed by the electric field of the edge of the element e, L is the number of the discrete basic elements e in the whole electromagnetic field solving area,

is a stiffness matrix of the basic cell e,

In one embodiment, the reference frequency point calculating module calculates the reference frequency point of the integrated circuit by:

In one embodiment, the critical frequency points are calculated based on the following formula:

In one embodiment, the reference frequency point calculating module calculates the reference frequency point and the field solution thereof, and specifically includes the following steps:

where ω is the angular frequency of the electromagnetic wave, E is the electric field, E_currFor the current angular frequency omega_currField solution of ∈ at₀Is a preset lower error threshold value epsilon₁Is a preset upper limit of the error threshold, epsilon₀<ε₁。

The above description is only for the specific embodiments of the present application, but the scope of the present application is not limited thereto, and any changes or substitutions that can be easily conceived by those skilled in the art within the technical scope of the present application should be covered within the scope of the present application. Therefore, the protection scope of the present application shall be subject to the protection scope of the claims.

Claims

1. A full-wave electromagnetic simulation method of an integrated circuit under a lossless dispersive medium is characterized by comprising the following steps:

wherein the content of the first and second substances,

2. The method for full-wave electromagnetic simulation of an integrated circuit under a lossless dispersive medium according to claim 1, wherein said establishing a matrix equation under a lossless dispersive medium comprises:

is a stiffness matrix of the basic cell e,

3. The method for full-wave electromagnetic simulation of an integrated circuit under a lossless dispersive medium according to claim 2, wherein said calculating the reference frequency points of the integrated circuit by an iterative method based on said matrix equation comprises:

4. The method for full-wave electromagnetic simulation of an integrated circuit under a lossless dispersive medium according to claim 3, wherein the critical frequency points are calculated based on the following formula:

5. The method for full-wave electromagnetic simulation of an integrated circuit under a lossless dispersive medium according to claim 3 or 4, wherein the calculating the reference frequency point and the field solution thereof based on the critical frequency point by an iterative calculation method comprises:

6. An integrated circuit full-wave electromagnetic simulation system under the condition of no loss and frequency dispersion medium, which is characterized by comprising:

wherein the content of the first and second substances,

7. The system for full-wave electromagnetic simulation of an integrated circuit under a lossless dispersive medium of claim 6, wherein the matrix equation building block builds the matrix equations under a lossless dispersive medium by:

taking partial derivative of the discrete functional and making the partial derivative be 0 to obtain the following formula lossless dispersive mediumThe following matrix equation: (K)₁-ω²K₂) E ═ b (ω), where E is the basic unit, E is the electric field, b^eAs an external excitation source for the basic cell E, E^eAn electric field vector formed by an electric field at an edge of the elementary cell e, L is the number of the elementary cells e whose entire electromagnetic field solution area is discrete,

is a stiffness matrix of the basic cell e,

8. The system for full-wave electromagnetic simulation of an integrated circuit under lossless dispersive medium of claim 7, wherein the reference frequency point calculating module calculates the reference frequency point of the integrated circuit by:

9. The full-wave electromagnetic simulation system of an integrated circuit under lossless dispersive medium of claim 8, wherein the critical frequency points are calculated based on the following formula:

10. The full-wave electromagnetic simulation system of an integrated circuit under lossless dispersive medium according to claim 8 or 9, wherein the reference frequency point calculating module calculates the reference frequency point and the field solution thereof, and comprises the following steps: