CN112989735A - Iterative solution method and device for interlayer coupling of multilayer very large scale integrated circuit - Google Patents

Iterative solution method and device for interlayer coupling of multilayer very large scale integrated circuit Download PDF

Info

Publication number
CN112989735A
CN112989735A CN202110425193.5A CN202110425193A CN112989735A CN 112989735 A CN112989735 A CN 112989735A CN 202110425193 A CN202110425193 A CN 202110425193A CN 112989735 A CN112989735 A CN 112989735A
Authority
CN
China
Prior art keywords
layer
source
integrated circuit
layers
field
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202110425193.5A
Other languages
Chinese (zh)
Other versions
CN112989735B (en
Inventor
唐章宏
邹军
王芬
汲亚飞
黄承清
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Wisechip Simulation Technology Co Ltd
Original Assignee
Beijing Wisechip Simulation Technology Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing Wisechip Simulation Technology Co Ltd filed Critical Beijing Wisechip Simulation Technology Co Ltd
Priority to CN202110425193.5A priority Critical patent/CN112989735B/en
Publication of CN112989735A publication Critical patent/CN112989735A/en
Application granted granted Critical
Publication of CN112989735B publication Critical patent/CN112989735B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/30—Circuit design
    • G06F30/32—Circuit design at the digital level
    • G06F30/33—Design verification, e.g. functional simulation or model checking
    • G06F30/3308—Design verification, e.g. functional simulation or model checking using simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
    • G06T17/20—Finite element generation, e.g. wire-frame surface description, tesselation

Abstract

The invention provides an iterative solution method and a device for interlayer coupling of a multilayer very large scale integrated circuit, wherein the iterative solution method comprises the following steps: firstly, set up themThe initial active layer of the source layer is all layers of the integrated circuit; secondly, to the firstmThe source layer is iterated, and the second step is calculated through a dyadic Green function in the iteration processmCurrent pair of layer distributionlInfluence of the layerG ml And update the firstlSource item of layer, tolApplying two-dimensional finite element to the layer to calculate its field distribution so as to update the field and current distribution of the layer, and obtaining the change of the field of the layer compared with the previous iteration resultdE ml Comparison ofG ml Determining a negligible layer with the effective influence value of the dynamically calculated dyadic Green function, and modifying the layermActive layer range of the source layer; through repeated iterations of the source layers until the influence change of all the source layersSo that the change in the field of the affected layer is less than the specified threshold and the iteration ends. According to the method and the device, the complexity of the three-dimensional problem and the occupied memory can be reduced under the condition that the calculation precision is not reduced.

Description

Iterative solution method and device for interlayer coupling of multilayer very large scale integrated circuit
Technical Field
The invention relates to the technical field of integrated circuits, in particular to an iterative solution method and device for interlayer coupling of a multilayer ultra-large scale integrated circuit.
Background
When the integrated circuit works, a high-frequency alternating electromagnetic field can be formed on a multilayer layout of the integrated circuit due to the transmission of high-speed signals, and meanwhile, in order to improve the performance of electronic equipment, reduce the volume and reduce the cost, transistors, other components and circuits are integrated on a small semiconductor substrate. In order to realize more functions, the ultra-large scale integrated circuit has a structure from tens of layers to hundreds of layers, each layer of structure is extremely complex, millions or even tens of millions of transistors are integrated, and the ultra-large scale integrated circuit has a multi-scale structure from a centimeter level to the latest nanometer level at present. In order to ensure that the integrated circuit can normally work and realize the function designed in advance, the power integrity and the signal integrity of the integrated circuit need to be ensured firstly, so that the power integrity and the signal integrity of the integrated circuit with a multi-scale structure of tens of layers and hundreds of layers need to be accurately analyzed by adopting an electromagnetic field analysis method, which is a great problem of the electromagnetic field analysis of the ultra-large scale integrated circuit.
A conventional method of analyzing the electromagnetic response of three-dimensional very large scale integrated circuits is a three-dimensional electromagnetic field numerical calculation method, such as a three-dimensional finite element method. When the electromagnetic response of the three-dimensional very large scale integrated circuit is calculated by adopting a traditional numerical calculation method, after a truncation error of a certain region is set, the whole three-dimensional integrated circuit and a limited region outside the integrated circuit are determined as a calculation region, then the whole calculation region is subjected to grid division, the electromagnetic field distribution of the whole calculation region is calculated, and the electromagnetic response such as the electromagnetic field distribution, the current voltage of a designated port and the like of each layer of the integrated circuit is further calculated. However, the characteristic dimensions of the via holes, the wires and the like of the integrated circuit are nano-scale, the dimension of the whole integrated circuit is centimeter-scale, the calculation area determined according to the truncation error is decimeter-scale and meter-scale, and hundreds of millions of grids and unknown quantities can be generated by carrying out uniform grid subdivision on the multi-scale space and then analyzing the space electromagnetic radiation of the multi-scale space, so that the hardware (memory) cost and the CPU time cost are overlarge. Therefore, the electromagnetic response of the three-dimensional large-scale integrated circuit can be calculated by adopting a method combining a finite element method and a moment method. In the three-dimensional large-scale integrated circuit area, a finite element method is adopted; in a large-scale area outside the integrated circuit, a moment method is adopted; the finite element method and the moment method are coupled at the interface of the integrated circuit and the external space. Because the moment method only integrates aiming at the interface, a large number of grid units and unknowns can be reduced, but because the scale range of the integrated circuit is from nano-scale to centimeter-scale, the finite element method directly used for solving the integrated circuit can generate a huge sparse matrix, and because the finite element method and the moment method are coupled, the formed coupling matrix is a dense matrix at the interface, the non-zero element number of the whole sparse matrix and the solving complexity of the sparse matrix are greatly increased, and the calculation time is still long.
Disclosure of Invention
Objects of the invention
Based on the problems, the invention provides an iterative solution method and device for interlayer coupling of a multilayer very large scale integrated circuit. The starting point of the invention is that the interlayer coupling of a multilayer very large scale integrated circuit is regarded as the external excitation of the integrated circuit layer, the interlayer coupling of the multilayer integrated circuit to a certain layer of integrated circuit can be regarded as the superposition of the external excitation of a plurality of layer couplings, thus the electromagnetic field distribution of the multilayer integrated circuit does not need to be considered at one time, only the single-layer electromagnetic field distribution needs to be analyzed, the coupling of other layers to the layer is regarded as a plurality of external excitations for superposition, and the coupling between the layers is continuously corrected through iteration until convergence is reached; meanwhile, the attenuation rules of the electromagnetic field and the electromagnetic wave in the space can be known, the influence of the point source on any point in the space is weakened along with the increase of the distance between the point source and the point (specifically, the influence value is inversely proportional to the distance, and the electromagnetic wave is more quickly weakened from the source point to the field point in the space due to the reflection of the layer interface). Based on this fact, when designing the iterative solution method, only the influence exerted on the layers adjacent to the affected layer is considered, the influence of the layers beyond the affected layer is not considered, and the influence of all other layers on a certain layer is not always considered, which greatly accelerates the iterative solution time.
(II) technical scheme
As a first aspect of the invention, the invention discloses an iterative solution method for interlayer coupling of a multilayer very large scale integrated circuit, which comprises the following steps:
step S100, the large scale integrated circuit is summarizedN+1 layer, each layer numbered
Figure 443859DEST_PATH_IMAGE001
When considering the second aspect of LSImWhen the current source of the layer is in the current source of the layer, the layer is called the first layermA source layer provided with a secondmActive layer of source layer
Figure 95420DEST_PATH_IMAGE002
Is divided bymOthers of the source layerNLayers of the integrated circuit, note
Figure 297994DEST_PATH_IMAGE003
I.e. firstmThe farthest distance of influence of the source layer is
Figure 427624DEST_PATH_IMAGE002
A layer; the 0 th layer is a bottom layer; the portion on which the other source layers among the source items of all the source layers affect is set to 0.
Step S200, settingm=0。
Step S300, for the secondmSource layer, using dyadic Green function to calculatemSource layer to source layerlInfluence of the layer, isG ml Then based onG ml Update the firstlFirst among source items of a layermThe portion of the source layer affecting it, tolLayer-applied two-dimensional finite element calculation of its electromagnetic field distributionUpdating the electromagnetic field and current distribution of the layer, and calculating the change amount of the electromagnetic field of the layerdE ml Wherein
Figure 163368DEST_PATH_IMAGE004
(ii) a Is provided withm=m+1, ifmIf not more than N, repeatedly executing the step 300; otherwise, step S400 is executed.
Step S400, if
Figure 720251DEST_PATH_IMAGE005
After the iteration is finished, the electromagnetic field of each layer is output, wherein
Figure 659388DEST_PATH_IMAGE006
The iteration precision is preset; otherwise, step S500 is executed.
Step S500, calculate allG ml Maximum value ofG maxAnd minimum valueG minCalculating the effective influence value of the dyadic Green function
Figure 343441DEST_PATH_IMAGE007
WhereinthredsholdThe method is a preset discarding threshold value influenced by the dyadic Green function.
Step S600, all the satisfaction values are selectedG ml <GOf the conditionG ml Is marked asG thredshold Calculate allG thredshold Middle distance layermNearest layerl near Number of layers of
Figure 684424DEST_PATH_IMAGE008
Memory for recording
Figure 412209DEST_PATH_IMAGE009
Update
Figure 353489DEST_PATH_IMAGE010
Is composed of
Figure 762604DEST_PATH_IMAGE008
Average value of (i), i.e.
Figure 595362DEST_PATH_IMAGE011
The process proceeds to step S200.
Further, the step S300 further includes a step oflIteration of layers:
s310, in the case of the last iteration, the electromagnetic field of each layer of PCB is recorded as
Figure 697310DEST_PATH_IMAGE012
；
S320, setting
Figure 673356DEST_PATH_IMAGE013
；
S330, calculatingmCurrent source pair with distributed layerslInfluence of the layerG ml Based onG ml Update the firstlFirst among source items of a layermThe portion of the source layer affecting it, tolApplying two-dimensional finite element to the layer to calculate the electromagnetic field distribution so as to update the electromagnetic field and current distribution of the layer, and calculating the change of the electromagnetic field of the layerdE ml At this time it islElectromagnetic field of the layer becomesE l =E+dE ml ；
S340、l=l+1, return to step S330 until
Figure 869851DEST_PATH_IMAGE014
。
Further, depending on the particular structure of the integrated circuit hierarchy, secondmSource layer at the secondlInfluence of layer generationG ml Can be decomposed to be located at the secondmSource layer point current source atlSuperposition of electric field expressions generated by layers, the firstmSource layer point current source atlThe electric field expression generated by the layers is a special analytical expression given by using a dyadic Green function, and the current sources of the multilayer integrated circuit are distributed in a layered mode, namely the current density distributed on each metal layer of the integrated circuit layout with a complex shape is only equal to that of the current density distributed on each metal layer of the integrated circuit layout with a complex shapexAndythe axial direction is related tozAxial direction independence, current density distribution of onlyx，yAs a function of (c).
Further, said first stepmSource layer at the secondlInfluence of layer generationG ml Can be decomposed to be located at the secondmSource layer point current source atlThe specific method for superposing the electric field expressions generated by the layers is as follows: will be located at the firstmSource layer point current source atlThe electric field expression generated by the layer is used as the integrand of two-dimensional Gaussian integration, and the first calculation is based on the linear superposition principle of the fieldmField generated at the same position by simple polygonal planar current source of source layer, second in two-dimensional plane SmSource layer at the secondlThe layer-generated field is calculated by the two-dimensional gaussian integral:
Figure 982164DEST_PATH_IMAGE015
；
wherein the content of the first and second substances,E(x,y,z) At any point in space for the current source in the two-dimensional plane Sx,y,z) The field that is generated is,
Figure 51751DEST_PATH_IMAGE016
is an arbitrary position within the two-dimensional surface S: (u,v) At any point in space (x,y,z) An expression of the dyadic green function of the generated field, (u p ,v q ) Representing a gaussian integration point corresponding to a two-dimensional gaussian integration in the two-dimensional plane S,p,qrespectively representu,vIn the first directionpA first, aqThe number of the Gaussian integration points is equal to the total number of the points,
Figure 203509DEST_PATH_IMAGE017
is the weighting factor corresponding to the gaussian integration point.
Further, influence values according to the dyadic Green function in the iteration processGDetermines a layer that can be ignored, adaptively adjustsmThe source layer adjacent theretolExtent of influence exerted by the layer
Figure 751165DEST_PATH_IMAGE010
。
On the other hand, the disclosed iteration device for interlayer coupling of the multilayer very large scale integrated circuit comprises an action layer iteration module, a source item updating module, an electromagnetic field change quantity updating module andN+1 LSI layer, each layer numbered
Figure 983563DEST_PATH_IMAGE018
。
The action layer iteration module is used for iteratively updating the action layer of the source layer
Figure 676581DEST_PATH_IMAGE010
And is provided with the firstmActive layer of source layer
Figure 361641DEST_PATH_IMAGE010
Is divided bymAll the other N layers of the source layer, i.e. the integrated circuit
Figure 712987DEST_PATH_IMAGE019
。
The source layer iteration module is used for updating an iteration source layer.
The source item updating module is used for calculating the influence G of the updated source layer on all other layers by utilizing the dyadic Green function when the source layer is updated ml 。
The change amount update module of the electromagnetic field is based on G ml Update the firstlFirst among source items of a layermThe portion of the source layer that affects it is calculated by two-dimensional finite elementslThe electromagnetic field distribution of the layer is updated to update the electromagnetic field and current distribution of the layer, thereby calculating the change amount of the electromagnetic field of the layerdE ml 。
Further, the action layer iteration module selects all the satisfied calculation unitsG ml |<GOf the conditionG ml Is marked asG thredshold WhereinGCalculating effective influence values of the dyadic Green functionAll ofG thredshold Middle distance layermNearest layerl near Number of layers of
Figure 488308DEST_PATH_IMAGE008
Is marked as
Figure 165277DEST_PATH_IMAGE020
Update
Figure 524583DEST_PATH_IMAGE021
Is composed of
Figure 679621DEST_PATH_IMAGE008
Average value of (i), i.e.
Figure 824294DEST_PATH_IMAGE022
。
Further, the effective influence value of the dyadic Green functionGThe solution of (c) is as follows: selectingG ml Maximum value ofG maxAnd minimum valueG minCalculating the effective influence value of the dyadic Green function
Figure 157318DEST_PATH_IMAGE023
WhereinthredsholdThe method is a preset discarding threshold value influenced by the dyadic Green function.
Further, depending on the particular structure of the integrated circuit hierarchy, secondmSource layer at the secondlInfluence of layer generationG ml Can be decomposed to be located at the secondmSource layer point current source atlSuperposition of electric field expressions generated by layers, the firstmSource layer point current source atlThe electric field expression generated by the layers is a special analytical expression given by using a dyadic Green function, and the current sources of the multilayer integrated circuit are distributed in a layered mode, namely the current density distributed on each metal layer of the integrated circuit layout with a complex shape is only equal to that of the current density distributed on each metal layer of the integrated circuit layout with a complex shapexAndythe axial direction is related tozAxial direction independence, current density distribution of onlyx，yAs a function of (c).
Further, said first stepmSource layer at the secondlInfluence of layer G ml Can be decomposed to be located at the secondmSource layer point current source atlThe specific method for superposing the electric field expressions generated by the layers is as follows: will be located at the firstmSource layer point current source atlThe electric field expression generated by the layer is used as the integrand of two-dimensional Gaussian integration, and the first calculation is based on the linear superposition principle of the fieldmField generated at the same position by simple polygonal planar current source of source layer, second in two-dimensional plane SmSource layer at the secondlThe layer-generated field is calculated by the two-dimensional gaussian integral:
Figure 816969DEST_PATH_IMAGE024
；
wherein the content of the first and second substances,E(x,y,z) At any point in space for the current source in the two-dimensional plane Sx,y,z) The field that is generated is,
Figure 713381DEST_PATH_IMAGE025
is at any position in the two-dimensional surface S
Figure 492987DEST_PATH_IMAGE026
At any point in space (x,y,z) The expression of the dyadic green function of the generated field,
Figure 511758DEST_PATH_IMAGE027
representing a gaussian integration point corresponding to a two-dimensional gaussian integration in the two-dimensional plane S,p,qrespectively representu,vIn the first directionpA first, aqThe number of the Gaussian integration points is equal to the total number of the points,
Figure 596389DEST_PATH_IMAGE017
is the weighting factor corresponding to the gaussian integration point.
(III) advantageous effects
According to the iterative solution method and device for the interlayer coupling of the multilayer very large scale integrated circuit, the electromagnetic field distribution and the current distribution of the affected layer are updated immediately every time the influence of the source layer on other layers is calculated, so that the source layer corresponding to the affected layer is the latest when the influence of the affected layer on other layers is calculated. The approximate solution is updated through multiple iterations to enable the final result to approach the true value, so that the complexity of the three-dimensional problem is reduced, and the time occupied by a CPU and the memory occupied by the CPU are reduced under the condition of not reducing the calculation precision.
Drawings
The embodiments described below with reference to the drawings are exemplary and intended to be used for explaining and illustrating the present invention and should not be construed as limiting the scope of the present invention.
FIG. 1 is a block diagram of the main steps of a first embodiment of the present invention;
FIG. 2 is a logic execution block diagram of a first embodiment of the present invention;
FIG. 3 is a block diagram of the modules of a second embodiment of the present invention;
fig. 4 is an exploded schematic view of the effect of a point source on a field point in the present invention.
Detailed Description
In order to make the implementation objects, technical solutions and advantages of the present invention clearer, the technical solutions in the embodiments of the present invention will be described in more detail below with reference to the accompanying drawings in the embodiments of the present invention.
It should be noted that: in the drawings, the same or similar reference numerals denote the same or similar elements or elements having the same or similar functions throughout. The embodiments described are some embodiments of the present invention, not all embodiments, and features in embodiments and embodiments in the present application may be combined with each other without conflict. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.
In the description of the present invention, it is to be understood that the terms "central," "longitudinal," "lateral," "front," "rear," "left," "right," "vertical," "horizontal," "top," "bottom," "inner," "outer," and the like are used in the orientation or positional relationship indicated in the drawings, which are used for convenience in describing the invention and for simplicity in description, and are not intended to indicate or imply that the referenced device or element must have a particular orientation, be constructed and operated in a particular orientation, and are not to be considered limiting of the scope of the invention.
The following describes in detail a first embodiment of the iterative solution method and apparatus for interlayer coupling in a multi-layer very large scale integrated circuit according to the present invention with reference to fig. 1, 2 and 4. The iterative solution method for interlayer coupling of a multilayer very large scale integrated circuit provided by the embodiment comprises the following steps:
step S100, the large scale integrated circuit is summarizedN+1 layer, each layer numbered
Figure 109541DEST_PATH_IMAGE028
When considering the second aspect of LSImWhen the current source of the layer is in the current source of the layer, the layer is called the first layermA source layer provided with a secondmActive layer of source layer
Figure 759965DEST_PATH_IMAGE010
Is divided bymOthers of the source layerNLayers of the integrated circuit, note
Figure 621742DEST_PATH_IMAGE019
I.e. firstmThe farthest distance of influence of the source layer is
Figure 239674DEST_PATH_IMAGE010
A layer; the 0 th layer is a bottom layer; setting a portion to which other source layers among source items of all source layers affect to 0;
step S200, settingm=0；
Step S300, for the secondmSource layer, using dyadic Green function to calculatemSource layer to source layerlInfluence of the layer, isG ml Then based onG ml Update the firstlFirst among source items of a layermThe portion of the source layer affecting it, tolThe electromagnetic field distribution of the layer is calculated by applying a two-dimensional finite element to the layer, so that the electromagnetic field and the current distribution of the layer are updated, and the change quantity of the electromagnetic field of the layer is calculateddE ml Wherein
Figure 477889DEST_PATH_IMAGE029
(ii) a Is provided withm=m+1, ifm≤NRepeating the step 300; otherwise, executing step S400;
further, the step S300 is the steplThe iterative process of the layers is:
s310, in the case of the last iteration, the electromagnetic field of each layer of PCB is recorded as E = E l (l=0,1,2,…,N)；
S320, setting
Figure 405656DEST_PATH_IMAGE030
；
S330, calculatingmCurrent source pair with distributed layerslInfluence of the layerG ml Based onG ml Update the firstlFirst among source items of a layermThe portion of the source layer affecting it, tolApplying two-dimensional finite element to the layer to calculate the electromagnetic field distribution so as to update the electromagnetic field and current distribution of the layer, and calculating the change dE of the electromagnetic field of the layer ml At this time it islElectromagnetic field of the layer becomes
Figure 484339DEST_PATH_IMAGE031
；
S340、l=l+1, return to step S330 until
Figure 543562DEST_PATH_IMAGE032
。
Step S400, if
Figure 647784DEST_PATH_IMAGE033
After the iteration is finished, the electromagnetic field of each layer is output, wherein
Figure 757954DEST_PATH_IMAGE034
The iteration precision is preset; otherwise, executing step S500;
step S500, calculate allG ml Maximum value ofG maxAnd minimum valueG minCalculating the effective influence value of the dyadic Green function
Figure 492691DEST_PATH_IMAGE035
WhereinthredsholdA discarding threshold value for the influence of a preset dyadic Green function;
step S600, all the satisfaction values are selectedG ml <GOf the conditionG ml Is marked asG thredshold Calculate allG thredshold Middle distance layermNearest layerl near Number of layers of
Figure 491740DEST_PATH_IMAGE008
Memory for recording
Figure 868495DEST_PATH_IMAGE036
Update
Figure 302013DEST_PATH_IMAGE010
Is composed of
Figure 473231DEST_PATH_IMAGE008
Average value of (i), i.e.
Figure 959576DEST_PATH_IMAGE037
The process proceeds to step S200.
Further, as shown in FIG. 4, according to the particular structure of the integrated circuit hierarchy, the secondmSource layer at the secondlInfluence of layer generationG ml Can be decomposed to be located at the secondmSource layer point current source atlSuperposition of electric field expressions generated by layers, the firstmSource layer point current source atlThe electric field expression generated by the layer is a special analytical expression given by using a dyadic Green function, and the specific expression of the analytical expression is as follows: aiming at the frequency domain electromagnetic field of the multilayer integrated circuit layout, the electric field intensity generated by a point source at any layer field point is calculated by adopting a parallel vector Green function, and any one of the multilayer integrated circuit layout can be solved through the following formulaThe electric field strength of any point of the layer in nine directions.
The electric field generated by the point current source at the field point is expressed as:
Figure 347877DEST_PATH_IMAGE038
Figure 150748DEST_PATH_IMAGE039
Figure 24026DEST_PATH_IMAGE040
Figure 794405DEST_PATH_IMAGE041
Figure 716224DEST_PATH_IMAGE042
Figure 186651DEST_PATH_IMAGE043
Figure 168513DEST_PATH_IMAGE044
Figure 426188DEST_PATH_IMAGE045
wherein the content of the first and second substances,
Figure 214016DEST_PATH_IMAGE046
Figure 725900DEST_PATH_IMAGE047
Figure 691713DEST_PATH_IMAGE048
Figure 921837DEST_PATH_IMAGE049
Figure 700306DEST_PATH_IMAGE050
Figure 129013DEST_PATH_IMAGE051
iis the unit of an imaginary number,i 2=-1；
Figure 452678DEST_PATH_IMAGE052
representing a Bessel function of order 0;
Figure 717569DEST_PATH_IMAGE053
representing a Bessel function of order 1;
Figure 847199DEST_PATH_IMAGE054
expressed as a function of the Bessel integral coefficient,
Figure 68096DEST_PATH_IMAGE055
；x, y, zthe coordinates of the field points are represented,
Figure 139826DEST_PATH_IMAGE056
,
Figure 78963DEST_PATH_IMAGE057
,
Figure 746705DEST_PATH_IMAGE058
representing source point coordinates; angular frequency
Figure 635157DEST_PATH_IMAGE059
，
Figure 566204DEST_PATH_IMAGE060
Represents a frequency;
Figure 789375DEST_PATH_IMAGE061
indicating that the site is at the second
Figure 510075DEST_PATH_IMAGE061
A layer of a material selected from the group consisting of,
Figure 705564DEST_PATH_IMAGE062
is as follows
Figure 823824DEST_PATH_IMAGE061
At layer boundarieszCoordinates;
Figure 534291DEST_PATH_IMAGE063
,
Figure 12677DEST_PATH_IMAGE064
respectively represent
Figure 577520DEST_PATH_IMAGE065
The number of complex waves in the horizontal and vertical directions of the layer;
Figure 850369DEST_PATH_IMAGE066
respectively represent
Figure 64444DEST_PATH_IMAGE065
A layer horizontal dielectric constant, a vertical dielectric constant;
Figure 815362DEST_PATH_IMAGE067
,
Figure 500290DEST_PATH_IMAGE068
respectively representlHorizontal magnetic conductivity and vertical magnetic conductivity of the layer;
Figure 944041DEST_PATH_IMAGE069
is shown aslThe anisotropy coefficient of the layer;
Figure 894680DEST_PATH_IMAGE070
,
Figure 465600DEST_PATH_IMAGE071
respectively representlIntegral coefficients of complex wave numbers of the horizontal and vertical layers;
Figure 755768DEST_PATH_IMAGE072
respectively representlThe undetermined coefficient of a layer,A l , B l the following linear equation is solved:
Figure 619687DEST_PATH_IMAGE073
T1is 2n×2nThe complex number matrix of (a) is,
Figure 995305DEST_PATH_IMAGE074
is of length 2nA complex vector of (a);
Figure 116056DEST_PATH_IMAGE075
Figure 323047DEST_PATH_IMAGE076
Figure 29971DEST_PATH_IMAGE077
Figure 689623DEST_PATH_IMAGE078
Figure 336767DEST_PATH_IMAGE079
the following linear equation is solved:
Figure 601526DEST_PATH_IMAGE080
T2is 2n×2nThe complex number matrix of (a) is,
Figure 869566DEST_PATH_IMAGE081
is of length 2nA complex vector of (a);
Figure 954196DEST_PATH_IMAGE082
Figure 405031DEST_PATH_IMAGE083
Figure 321035DEST_PATH_IMAGE084
Figure 713970DEST_PATH_IMAGE085
the following linear equation is solved:
Figure 535164DEST_PATH_IMAGE086
T3is 2n×2nThe complex number matrix of (a) is,
Figure 38958DEST_PATH_IMAGE087
is of length 2nA complex vector of (a);
Figure 560200DEST_PATH_IMAGE088
Figure 124037DEST_PATH_IMAGE089
Figure 432527DEST_PATH_IMAGE090
Figure 474433DEST_PATH_IMAGE091
Figure 365028DEST_PATH_IMAGE092
to representxOriented electric dipole in the second placelOf said electric field generated by said field points of the layerxA component;
Figure 584919DEST_PATH_IMAGE093
to representxOriented electric dipole in the second placelOf said electric field generated by said field points of the layeryA component;
Figure 928176DEST_PATH_IMAGE094
to representxOriented electric dipole in the second placelOf said electric field generated by said field points of the layerzA component;
Figure 23040DEST_PATH_IMAGE095
to representyOriented electric dipole in the second placelOf said electric field generated by said field points of the layerxA component;
Figure 971404DEST_PATH_IMAGE096
to representyOriented electric dipole in the second placelOf said electric field generated by said field points of the layeryA component;
Figure 831038DEST_PATH_IMAGE097
to representyOriented electric dipole in the second placelSaid electricity generated by said field points of the layerOf fieldszA component;
Figure 317383DEST_PATH_IMAGE098
to representzOriented electric dipole in the second placelOf said electric field generated by said field points of the layerxA component;
Figure 497829DEST_PATH_IMAGE099
to representzOriented electric dipole in the second placelOf said electric field generated by said field points of the layeryA component;
Figure 97437DEST_PATH_IMAGE100
to representzOriented electric dipole in the second placelOf said electric field generated by said field points of the layerzAnd (4) components.
The current sources of the multi-layer integrated circuit are distributed in a layered manner, namely the current density distributed on each metal layer of the integrated circuit layout with a complex shape is only equal to that of the current sourcexAndythe axial direction is related tozAxial direction independence, current density distribution of onlyx，yAs a function of (c).
Further, said first stepmSource layer at the secondlInfluence of layer generationG ml Can be decomposed to be located at the secondmSource layer point current source atlThe specific method for superposing the electric field expressions generated by the layers is as follows: will be located at the firstmSource layer point current source atlThe electric field expression generated by the layer is used as the integrand of two-dimensional Gaussian integration, and the first calculation is based on the linear superposition principle of the fieldmField generated at the same position by simple polygonal planar current source of source layer, second in two-dimensional plane SmSource layer at the secondlThe layer-generated field is calculated by the two-dimensional gaussian integral:
Figure 659131DEST_PATH_IMAGE101
；
wherein the content of the first and second substances,E(x,y,z) Is the two-dimensional surfaceThe current source in S is at any point in space (x,y,z) The field that is generated is,
Figure 976980DEST_PATH_IMAGE102
is an arbitrary position within the two-dimensional surface S: (u,v) At any point in space (x,y,z) The expression of the dyadic green function of the generated field,
Figure 695537DEST_PATH_IMAGE103
representing a gaussian integration point corresponding to a two-dimensional gaussian integration in the two-dimensional plane S,p,qrespectively representu,vIn the first directionpA first, aqThe number of the Gaussian integration points is equal to the total number of the points,
Figure 602182DEST_PATH_IMAGE017
is the weighting factor corresponding to the gaussian integration point.
Calculating fields generated by the current on the simple-shaped polygon at different positions of other layers of the integrated circuit, determining the fields generated by the current on the simple-shaped polygon on the layout of other layers of the integrated circuit based on the linear superposition principle of the fields, and determining the first field based on the linear superposition principle of the fieldsmSource layer at the secondlInfluence of layer generationG ml 。
Further, the specific method for calculating the two-dimensional finite element comprises the following steps:
for the direct current electric field model, the three-dimensional model of the multilayer integrated circuit refers to the conductivity in the direct current electric field model
Figure 849624DEST_PATH_IMAGE104
Potential of the electrodeuAll the distributions of (A) and (B) are three-dimensional space coordinatesx,y,z) I.e.:
Figure 405501DEST_PATH_IMAGE105
，
Figure 131012DEST_PATH_IMAGE106
the function of the three-dimensional model satisfies the following equation（1）：
Figure 705213DEST_PATH_IMAGE107
In the equation (1),
and boundary condition (2):
Figure 372823DEST_PATH_IMAGE108
in the formula
Figure 337368DEST_PATH_IMAGE109
Is a boundary of the first type and is,nis normal to the boundary of the second type,
Figure 667900DEST_PATH_IMAGE110
represents a potentialuAt the first kind boundary
Figure 299870DEST_PATH_IMAGE109
Value of above, using
Figure 872802DEST_PATH_IMAGE111
It is shown that,
Figure 652540DEST_PATH_IMAGE112
bulk current density for external circuits;
the dimension of an actual PCB or a chip packaged board in the multilayer super large scale integrated circuit is far larger than the thickness of the metal layer, so that the three-dimensional direct current field problem of the multilayer integrated circuit is simplified into a two-dimensional direct current field problem;
the field solving equation set established by the finite element method for the two-dimensional model is an equation set (3):
Figure 985432DEST_PATH_IMAGE114
in the formula (I), theI(u) In order to be a functional function,tis the thickness of the metal layer or layers,
Figure 957061DEST_PATH_IMAGE115
as a grid celleThe electrical conductivity of (a) a (b),
Figure 779524DEST_PATH_IMAGE116
as a grid celleThe potential of (a) is set to be,
Figure 233508DEST_PATH_IMAGE117
as a grid celleThe area of (a) is,
Figure 838933DEST_PATH_IMAGE118
as the density of the surface current, the current density,
Figure 242232DEST_PATH_IMAGE119
representing grid cellseThe edge of (1);
for the alternating electromagnetic field model, the three-dimensional model of the multilayer integrated circuit refers to the dielectric constant in the three-dimensional model of the electromagnetic response characteristic in the frequency domain simulation of the multilayer VLSI
Figure 924011DEST_PATH_IMAGE120
Magnetic permeability of
Figure 147182DEST_PATH_IMAGE121
Electric field intensityEMagnetic field intensityHAll the distributions of (A) and (B) are three-dimensional space coordinatesx,y,z) I.e.:
Figure 618615DEST_PATH_IMAGE122
,
Figure 63372DEST_PATH_IMAGE123
,
Figure 962057DEST_PATH_IMAGE124
，
Figure 672525DEST_PATH_IMAGE125
the function of the three-dimensional model satisfies the following equation:
Figure 636064DEST_PATH_IMAGE127
in the formulaJFor the purpose of the applied current density distribution,
Figure 748376DEST_PATH_IMAGE128
for the angular frequency simulated for the integrated circuit,
Figure 270493DEST_PATH_IMAGE129
indicating the strength of the magnetic fieldHThe degree of rotation of the screw is reduced,
Figure 202677DEST_PATH_IMAGE130
indicates the electric field intensityEThe degree of rotation of the screw is reduced,jis the unit of an imaginary number,j 2=-1；
the board size of the actual PCB or chip package in the multilayer VLSI is far larger than the metal layer spacing, the three-dimensional model of the electromagnetic response characteristics in the frequency domain simulation of the multilayer VLSI is simplified into a two-dimensional model, and the dielectric constant in the model is at the moment
Figure 281492DEST_PATH_IMAGE120
Magnetic permeability of
Figure 202305DEST_PATH_IMAGE121
Electric field intensityEMagnetic field intensityHAll the distributions are two-dimensional plane coordinates (x,y) I.e.:
Figure 442794DEST_PATH_IMAGE131
，
Figure 127853DEST_PATH_IMAGE132
，
Figure 931730DEST_PATH_IMAGE133
，
Figure 753055DEST_PATH_IMAGE134
distribution thereof andzindependent of and potential in the fielduAnd surface current densityJ sSatisfies the following conditions:
Figure 164445DEST_PATH_IMAGE135
in the formula (I), the compound is shown in the specification,
Figure 25216DEST_PATH_IMAGE136
respectively representx, y, zThe unit vector of the direction is,E zof electric field strengthzThe direction component of the light beam is,H xandH yrespectively of magnetic field strengthxAndythe direction component of the light beam is,his the metal layer spacing;
through the simplification from the three-dimensional model to the two-dimensional model, the two-dimensional finite element functional extreme value formula corresponding to the two-dimensional model is obtained as follows:
Figure 445833DEST_PATH_IMAGE137
in the formula (I), the compound is shown in the specification,
Figure 121665DEST_PATH_IMAGE138
in order to be a functional function,
Figure 156486DEST_PATH_IMAGE139
it is shown that the extreme value is taken for the functional,
Figure 816137DEST_PATH_IMAGE140
as a grid celliThe surface admittance of the first and second electrodes,
Figure 509287DEST_PATH_IMAGE141
is a boundary
Figure 259199DEST_PATH_IMAGE142
The boundary condition of the opening of (a),u kis a boundary
Figure 277971DEST_PATH_IMAGE143
The distribution of the electric potential on the upper side,
Figure 893760DEST_PATH_IMAGE144
indicating a position to the right of the boundary and infinitely close to the boundary,
Figure 843130DEST_PATH_IMAGE145
indicating a position to the left of the boundary and infinitely close to the boundary,
Figure 759134DEST_PATH_IMAGE146
representing grid cellsiThe area of (a) is,
Figure 683227DEST_PATH_IMAGE147
as a grid celliThe current density of (a) is,
Figure 5887DEST_PATH_IMAGE148
as a grid celliThe surface resistance of the glass substrate is higher than the surface resistance of the glass substrate,
Figure 306418DEST_PATH_IMAGE149
as a grid celliThe potential of (a) is set to be,kis referred to askAnd (4) a boundary.
Further, determining a negligible layer according to the magnitude of the influence value G of each layer of the dyadic Green function in the iteration process, and adaptively adjustingmSource layer to itlExtent of influence exerted by the layer
Figure 14611DEST_PATH_IMAGE010
. Since the attenuation law of the electromagnetic field and the electromagnetic wave in the space can be known, the influence of the point source on any point in the space is weakened along with the increase of the distance between the point source and the point (specifically, the influence value is inversely proportional to the distance, and the electromagnetic wave is more quickly weakened from the source point to the field point in the space due to the reflection of the layer interface), therefore, when the influence of the point source on the space point is calculated by using the parallel vector green function, the influence of the point source on the space point can be considered to be negligible when the distance between the space point and the point source is greater than a certain degree, or after the number of the medium layers separated from the space point is greater than a certain degree. Based on this fact, when designing the iterative solution method, only the influence exerted on the layers adjacent to the point source is considered, and the influence is not considered in the layers beyond the layers, which will beGreatly speeding up the iterative solution time.
Particularly, when the voltage drop and the current distribution of a power supply layer of the integrated circuit are analyzed, the working frequency is low frequency, the direct current field model is adopted for analysis, no space coupling exists between the integrated circuit layers at the moment, only physical coupling exists, namely, the layers which are connected with each other through the through hole and the external circuit are mutually coupled, at the moment, the mutual influence layers between the integrated circuit layers are determined, and iteration is not needed to influence the influence range
Figure 358873DEST_PATH_IMAGE010
And (6) correcting.
As can be seen from the above iteration steps, in the iteration process, according to the magnitude of the influence value of the dyadic green function of each layer, the range of the influence exerted by each source layer on other layers is adaptively adjusted, instead of exerting the influence of the source on other layers on all other layers every time, so that the iterative computation is accelerated. The advantage of the above iterative method is that the electromagnetic field distribution and the current distribution of the affected layer are updated immediately each time the influence of the source layer on other layers is calculated, thereby ensuring that the source layer corresponding to the affected layer is up to date when the influence of the affected layer on other layers is calculated.
A second embodiment of the iterative method and apparatus for interlayer coupling in a multi-layered very large scale integrated circuit according to the present invention is described in detail with reference to fig. 3 and 4. As shown in fig. 3 and 4, the iterative apparatus for interlayer coupling of a multi-layer very large scale integrated circuit provided in this embodiment includes an active layer iteration module, a source item update module, an electromagnetic field variation update module, a,N+1 LSI layer, each layer numbered
Figure 214834DEST_PATH_IMAGE018
。
The action layer iteration module is used for iteratively updating the action layer of the source layer
Figure 256739DEST_PATH_IMAGE010
And is provided with the firstmActive layer of source layer
Figure 644207DEST_PATH_IMAGE010
Is divided bymOthers of the source layerNLayers of integrated circuits, i.e.
Figure 175682DEST_PATH_IMAGE150
(ii) a Since the attenuation law of the electromagnetic field and the electromagnetic wave in the space can be known, the influence of the point source on any point in the space is weakened along with the increase of the distance between the point source and the point (specifically, the influence value is inversely proportional to the distance, and the electromagnetic wave is more quickly weakened from the source point to the field point in the space due to the reflection of the layer interface), therefore, when the influence of the point source on the space point is calculated by using the parallel vector green function, the influence of the point source on the space point can be considered to be negligible when the distance between the space point and the point source is greater than a certain degree, or after the number of the medium layers separated from the space point is greater than a certain degree. Based on this fact, when designing the iterative solution method, only the influence exerted on the layers adjacent to the point source is considered, and the influence is not considered in the layers beyond the layers, which greatly accelerates the iterative solution time.
The source layer iteration module is used for updatingmA source layer;
the source item updating module is used for calculating the updated second source item by utilizing a dyadic Green function when the source layer is updatedmSource layer to all the secondlInfluence of the layer, isG ml Based onG ml Update the firstlFirst among source items of a layermA portion on which the source layer affects;
the change amount updating module of the electromagnetic field is based onG ml Update the firstlFirst among source items of a layermThe portion of the source layer that affects it is calculated by two-dimensional finite elementmAfter influence of the source layerlThe electromagnetic field distribution of the layer is updated to update the electromagnetic field and current distribution of the layer, thereby calculating the change amount of the electromagnetic field of the layerdE ml 。
Further, the action layer iteration module selects all the satisfied calculation unitsG ml |<GOf the conditionG ml Is marked asG thredshold WhereinGCalculating all the effective influence values of the dyadic Green functionG thredshold Middle distance layermNearest layerl near Number of layers of
Figure 456622DEST_PATH_IMAGE008
Is marked as
Figure 82644DEST_PATH_IMAGE151
Update
Figure 31009DEST_PATH_IMAGE010
Is composed of
Figure 733386DEST_PATH_IMAGE008
Average value of (i), i.e.
Figure 252354DEST_PATH_IMAGE152
。
Further, the effective influence value of the dyadic Green functionGThe solution of (c) is as follows: selectingG ml Maximum value ofG maxAnd minimum valueG minCalculating the effective influence value of the dyadic Green function
Figure 901641DEST_PATH_IMAGE153
WhereinthredsholdThe method is a preset discarding threshold value influenced by the dyadic Green function.
Further, as shown in FIG. 4, according to the particular structure of the integrated circuit hierarchy, the secondmSource layer at the secondlInfluence of layer G ml Can be decomposed to be located at the secondmSource layer point current source atlSuperposition of electric field expressions generated by layers, the firstmSource layer point current source atlThe electric field expression generated by the layer is a special analytical expression given by using a dyadic Green function, and the specific expression of the analytical expression is as follows: aiming at the frequency domain electromagnetic field of the multilayer integrated circuit layout, the electric field intensity generated by a point source at any layer field point is calculated by adopting a parallel vector Green function, and any layer of the multilayer integrated circuit layout can be solved through the following formulaThe electric field strength of the nine orientations of the point.
The electric field generated by the point current source at the field point is expressed as:
Figure 16097DEST_PATH_IMAGE038
Figure 827058DEST_PATH_IMAGE039
Figure 98901DEST_PATH_IMAGE040
Figure 755142DEST_PATH_IMAGE154
Figure 927366DEST_PATH_IMAGE042
Figure 440387DEST_PATH_IMAGE043
Figure 933947DEST_PATH_IMAGE044
Figure 925037DEST_PATH_IMAGE045
wherein the content of the first and second substances,
Figure 155030DEST_PATH_IMAGE155
Figure 324106DEST_PATH_IMAGE156
Figure 819809DEST_PATH_IMAGE157
Figure 598278DEST_PATH_IMAGE158
Figure 495827DEST_PATH_IMAGE159
Figure 835804DEST_PATH_IMAGE160
iis the unit of an imaginary number,i 2=-1；
Figure 209016DEST_PATH_IMAGE052
representing a Bessel function of order 0;
Figure 338646DEST_PATH_IMAGE053
representing a Bessel function of order 1;
Figure 575855DEST_PATH_IMAGE054
expressed as a function of the Bessel integral coefficient,
Figure 132738DEST_PATH_IMAGE055
；x, y, zthe coordinates of the field points are represented,
Figure 71875DEST_PATH_IMAGE056
,
Figure 988884DEST_PATH_IMAGE057
,
Figure 392184DEST_PATH_IMAGE058
representing source point coordinates; angular frequency
Figure 57652DEST_PATH_IMAGE059
，
Figure 297134DEST_PATH_IMAGE060
Represents a frequency;
Figure 971829DEST_PATH_IMAGE061
indicating that the site is at the second
Figure 682165DEST_PATH_IMAGE061
A layer of a material selected from the group consisting of,
Figure 580851DEST_PATH_IMAGE062
is as follows
Figure 494580DEST_PATH_IMAGE061
At layer boundarieszCoordinates;
Figure 508718DEST_PATH_IMAGE063
,
Figure 824292DEST_PATH_IMAGE064
respectively represent
Figure 346409DEST_PATH_IMAGE061
The number of complex waves in the horizontal and vertical directions of the layer;
Figure 544173DEST_PATH_IMAGE066
respectively represent
Figure 826249DEST_PATH_IMAGE061
A layer horizontal dielectric constant, a vertical dielectric constant;
Figure 481484DEST_PATH_IMAGE067
,
Figure 987552DEST_PATH_IMAGE068
respectively representlHorizontal magnetic conductivity and vertical magnetic conductivity of the layer;
Figure 859562DEST_PATH_IMAGE069
is shown aslThe anisotropy coefficient of the layer;
Figure 476488DEST_PATH_IMAGE070
,
Figure 501075DEST_PATH_IMAGE071
respectively representlIntegral coefficients of complex wave numbers of the horizontal and vertical layers;
Figure 132039DEST_PATH_IMAGE072
respectively representlThe undetermined coefficient of a layer,A l , B l the following linear equation is solved:
Figure 304395DEST_PATH_IMAGE073
T1is 2n×2nThe complex number matrix of (a) is,
Figure 459432DEST_PATH_IMAGE074
is of length 2nA complex vector of (a);
Figure 587794DEST_PATH_IMAGE161
Figure 435665DEST_PATH_IMAGE162
Figure 783732DEST_PATH_IMAGE077
Figure 476881DEST_PATH_IMAGE163
Figure 538378DEST_PATH_IMAGE079
the following linear equation is solved:
Figure 478521DEST_PATH_IMAGE080
T2is 2n×2nThe complex number matrix of (a) is,
Figure 625469DEST_PATH_IMAGE081
is of length 2nA complex vector of (a);
Figure 122309DEST_PATH_IMAGE164
Figure 726728DEST_PATH_IMAGE083
Figure 916401DEST_PATH_IMAGE165
Figure 222748DEST_PATH_IMAGE085
the following linear equation is solved:
Figure 506968DEST_PATH_IMAGE086
T3is 2n×2nThe complex number matrix of (a) is,
Figure 543057DEST_PATH_IMAGE087
is of length 2nA complex vector of (a);
Figure 841314DEST_PATH_IMAGE167
Figure 448007DEST_PATH_IMAGE168
Figure 286650DEST_PATH_IMAGE090
Figure 114929DEST_PATH_IMAGE169
Figure 895672DEST_PATH_IMAGE092
to representxOriented electric dipole in the second placelOf said electric field generated by said field points of the layerxA component;
Figure 176612DEST_PATH_IMAGE093
to representxOriented electric dipole in the second placelOf said electric field generated by said field points of the layeryA component;
Figure 553367DEST_PATH_IMAGE094
to representxOriented electric dipole in the second placelOf said electric field generated by said field points of the layerzA component;
Figure 49201DEST_PATH_IMAGE095
to representyOriented electric dipole in the second placelOf said electric field generated by said field points of the layerxA component;
Figure 689261DEST_PATH_IMAGE096
to representyOriented electric dipole in the second placelOf said electric field generated by said field points of the layeryA component;
Figure 519814DEST_PATH_IMAGE097
to representyOriented electric dipole in the second placelSaid electricity generated by said field points of the layerOf fieldszA component;
Figure 949527DEST_PATH_IMAGE098
to representzOriented electric dipole in the second placelOf said electric field generated by said field points of the layerxA component;
Figure 486819DEST_PATH_IMAGE099
to representzOriented electric dipole in the second placelOf said electric field generated by said field points of the layeryA component;
Figure 360097DEST_PATH_IMAGE100
to representzOriented electric dipole in the second placelOf said electric field generated by said field points of the layerzAnd (4) components.
The current sources of the multi-layer integrated circuit are distributed in a layered manner, namely the current density distributed on each metal layer of the integrated circuit layout with a complex shape is only equal to that of the current sourcexAndythe axial direction is related tozAxial direction independence, current density distribution of onlyx，yAs a function of (c).
Further, said first stepmSource layer at the secondlInfluence of layer generationG ml Can be decomposed to be located at the secondmSource layer point current source atlThe specific method for superposing the electric field expressions generated by the layers is as follows: will be located at the firstmSource layer point current source atlThe electric field expression generated by the layer is used as the integrand of two-dimensional Gaussian integration, and the first calculation is based on the linear superposition principle of the fieldmField generated at the same position by simple polygonal planar current source of source layer, second in two-dimensional plane SmSource layer at the secondlThe layer-generated field is calculated by the two-dimensional gaussian integral:
Figure 366361DEST_PATH_IMAGE170
；
wherein the content of the first and second substances,E(x,y,z) Is that it isCurrent source in two-dimensional plane S at any point in space ()x,y,z) The field that is generated is,
Figure 350498DEST_PATH_IMAGE102
is an arbitrary position within the two-dimensional surface S: (u,v) At any point in space (x,y,z) The expression of the dyadic green function of the generated field,
Figure 804613DEST_PATH_IMAGE103
representing a gaussian integration point corresponding to a two-dimensional gaussian integration in the two-dimensional plane S,p,qrespectively representu,vIn the first directionpA first, aqThe number of the Gaussian integration points is equal to the total number of the points,
Figure 35743DEST_PATH_IMAGE017
is the weight factor corresponding to the gaussian integral point;
calculating fields generated by the current on the simple-shaped polygon at different positions of other layers of the integrated circuit, determining the fields generated by the current on the simple-shaped polygon on the layout of other layers of the integrated circuit based on the linear superposition principle of the fields, and determining the first field based on the linear superposition principle of the fieldsmSource layer at the secondlInfluence of layer generationG ml 。
Further, the specific method for calculating the two-dimensional finite element comprises the following steps:
for the direct current electric field model, the three-dimensional model of the multilayer integrated circuit refers to the conductivity in the direct current electric field model
Figure 575309DEST_PATH_IMAGE104
Potential of the electrodeuAll the distributions of (A) and (B) are three-dimensional space coordinatesx,y,z) I.e.:
Figure 363136DEST_PATH_IMAGE105
，
Figure 637471DEST_PATH_IMAGE106
the function of the three-dimensional model satisfies the following equationThe process (1):
Figure 586973DEST_PATH_IMAGE171
in the equation (1),
and boundary condition (2):
Figure 879414DEST_PATH_IMAGE172
in the formula
Figure 392303DEST_PATH_IMAGE109
Is a boundary of the first type and is,nis normal to the boundary of the second type,
Figure 821011DEST_PATH_IMAGE110
represents a potentialuAt the first kind boundary
Figure 206993DEST_PATH_IMAGE109
Value of above, using
Figure 409566DEST_PATH_IMAGE111
It is shown that,
Figure 539196DEST_PATH_IMAGE112
bulk current density for external circuits;
the dimension of an actual PCB or a chip packaged board in the multilayer super large scale integrated circuit is far larger than the thickness of the metal layer, so that the three-dimensional direct current field problem of the multilayer integrated circuit is simplified into a two-dimensional direct current field problem;
the field solving equation set established by the finite element method for the two-dimensional model is an equation set (3):
Figure 274940DEST_PATH_IMAGE173
in the formula (I), theI(u) In order to be a functional function,tis the thickness of the metal layer or layers,
Figure 566244DEST_PATH_IMAGE115
as a grid celleThe electrical conductivity of (a) a (b),
Figure 567698DEST_PATH_IMAGE116
as a grid celleThe potential of (a) is set to be,
Figure 189435DEST_PATH_IMAGE117
as a grid celleThe area of (a) is,
Figure 592734DEST_PATH_IMAGE118
as the density of the surface current, the current density,
Figure 258202DEST_PATH_IMAGE119
representing grid cellseThe edge of (1);
for the alternating electromagnetic field model, the three-dimensional model of the multilayer integrated circuit refers to the dielectric constant in the three-dimensional model of the electromagnetic response characteristic in the frequency domain simulation of the multilayer VLSI
Figure 996219DEST_PATH_IMAGE120
Magnetic permeability of
Figure 670914DEST_PATH_IMAGE121
Electric field intensityEMagnetic field intensityHAll the distributions of (A) and (B) are three-dimensional space coordinatesx,y,z) I.e.:
Figure 882715DEST_PATH_IMAGE174
,
Figure 984663DEST_PATH_IMAGE175
,
Figure 960710DEST_PATH_IMAGE176
，
Figure 157205DEST_PATH_IMAGE177
the function of the three-dimensional model satisfies the following equation:
Figure 269517DEST_PATH_IMAGE178
in the formulaJFor the purpose of the applied current density distribution,
Figure 339104DEST_PATH_IMAGE128
for the angular frequency simulated for the integrated circuit,
Figure 756441DEST_PATH_IMAGE129
indicating the strength of the magnetic fieldHThe degree of rotation of the screw is reduced,
Figure 225469DEST_PATH_IMAGE130
indicates the electric field intensityEThe degree of rotation of the screw is reduced,jis the unit of an imaginary number,j 2=-1；
the board size of the actual PCB or chip package in the multilayer VLSI is far larger than the metal layer spacing, the three-dimensional model of the electromagnetic response characteristics in the frequency domain simulation of the multilayer VLSI is simplified into a two-dimensional model, and the dielectric constant in the model is at the moment
Figure 457867DEST_PATH_IMAGE120
Magnetic permeability of
Figure 167197DEST_PATH_IMAGE121
Electric field intensityEMagnetic field intensityHAll the distributions are two-dimensional plane coordinates (x,y) I.e.:
Figure 806251DEST_PATH_IMAGE131
，
Figure 377172DEST_PATH_IMAGE132
，
Figure 73864DEST_PATH_IMAGE179
，
Figure 688516DEST_PATH_IMAGE180
distribution thereof andzindependent of and potential in the fielduAnd surface current densityJ sSatisfies the following conditions:
Figure 313401DEST_PATH_IMAGE181
in the formula (I), the compound is shown in the specification,
Figure 671701DEST_PATH_IMAGE136
respectively representx, y, zThe unit vector of the direction is,E zof electric field strengthzThe direction component of the light beam is,H xandH yrespectively of magnetic field strengthxAndythe direction component of the light beam is,his the metal layer spacing;
through the simplification from the three-dimensional model to the two-dimensional model, the two-dimensional finite element functional extreme value formula corresponding to the two-dimensional model is obtained as follows:
Figure 878692DEST_PATH_IMAGE137
in the formula (I), the compound is shown in the specification,
Figure 668838DEST_PATH_IMAGE182
in order to be a functional function,
Figure 266172DEST_PATH_IMAGE139
it is shown that the extreme value is taken for the functional,
Figure 474169DEST_PATH_IMAGE183
as a grid celliThe surface admittance of the first and second electrodes,
Figure 473349DEST_PATH_IMAGE141
is a boundary
Figure 180536DEST_PATH_IMAGE142
The boundary condition of the opening of (a),u kis a boundary
Figure 530746DEST_PATH_IMAGE143
The distribution of the electric potential on the upper side,
Figure 293165DEST_PATH_IMAGE144
indicating a position to the right of the boundary and infinitely close to the boundary,
Figure 396120DEST_PATH_IMAGE145
indicating a position to the left of the boundary and infinitely close to the boundary,
Figure 54634DEST_PATH_IMAGE146
representing grid cellsiThe area of (a) is,
Figure 688878DEST_PATH_IMAGE147
as a grid celliThe current density of (a) is,
Figure 677824DEST_PATH_IMAGE148
as a grid celliThe surface resistance of the glass substrate is higher than the surface resistance of the glass substrate,
Figure 713914DEST_PATH_IMAGE149
as a grid celliThe potential of (a) is set to be,kis referred to askAnd (4) a boundary.
The above shows that the device adaptively adjusts the range of the influence exerted by each source layer on other layers according to the magnitude of the influence value of the dyadic Green function of each layer in the iteration process, rather than exerting the influence of the source on other layers on all other layers every time, thereby accelerating the iterative computation. The method has the advantages that the electromagnetic field distribution and the current distribution of the affected layer are updated immediately every time the influence of the active layer on other layers is calculated, so that the active layer corresponding to the affected layer is the latest when the influence of the affected layer on other layers is calculated.
Particularly, when the voltage drop and the current distribution of a power supply layer of the integrated circuit are analyzed, the working frequency is low frequency, the direct current field model is adopted for analysis, no space coupling exists between the integrated circuit layers at the moment, only physical coupling exists, namely, the layers which are connected with each other through the through hole and the external circuit are mutually coupled, at the moment, the mutual influence layers between the integrated circuit layers are determined, and iteration is not needed to influence the influence range
Figure 995859DEST_PATH_IMAGE021
And (6) correcting.
The above description is only for the specific embodiment of the present invention, but the scope of the present invention 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 invention are included in the scope of the present invention. Therefore, the protection scope of the present invention shall be subject to the protection scope of the appended claims.

Claims (10)

1. An iterative solution method for interlayer coupling of a multilayer very large scale integrated circuit is characterized by comprising the following steps:
step S100, the large scale integrated circuit is summarizedN+1 layer, each layer numbered
Figure 406704DEST_PATH_IMAGE001
When considering the second aspect of LSImWhen the current source of the layer is in the current source of the layer, the layer is called the first layermA source layer provided with a secondmActive layer of source layer
Figure 120582DEST_PATH_IMAGE002
Is divided bymOthers of the source layerNLayers of the integrated circuit, note
Figure 634740DEST_PATH_IMAGE003
I.e. firstmThe farthest distance of influence of the source layer is
Figure 577420DEST_PATH_IMAGE002
A layer; the 0 th layer is a bottom layer; setting a portion to which other source layers among source items of all source layers affect to 0;
step S200, settingm=0；
Step S300, for the secondmSource layer, using dyadic Green function to calculatemSource layer to source layerlInfluence of the layer, denoted G ml Then based on G ml Update the firstlFirst among source items of a layermShadow of source layerPart of sound, tolThe electromagnetic field distribution of the layer is calculated by applying a two-dimensional finite element to the layer, so that the electromagnetic field and the current distribution of the layer are updated, and the change quantity of the electromagnetic field of the layer is calculateddE ml Wherein
Figure 250846DEST_PATH_IMAGE004
(ii) a Is provided withm=m+1, ifm≤NRepeating the step 300; otherwise, executing step S400;
step S400, if
Figure 119314DEST_PATH_IMAGE005
After the iteration is finished, the electromagnetic field of each layer is output, wherein
Figure 120768DEST_PATH_IMAGE006
The iteration precision is preset; otherwise, executing step S500;
step S500, calculate allG ml Maximum value ofG maxAnd minimum valueG minCalculating the effective influence value of the dyadic Green function
Figure 381985DEST_PATH_IMAGE007
WhereinthredsholdA discarding threshold value for the influence of a preset dyadic Green function;
step S600, all the satisfaction values are selectedG ml <GOf the conditionG ml Is marked asG thredshold Calculate allG thredshold Middle distance layermNearest layerl near Number of layers of
Figure 598334DEST_PATH_IMAGE008
Memory for recording
Figure 185173DEST_PATH_IMAGE009
Update
Figure 41965DEST_PATH_IMAGE002
Is composed of
Figure 919923DEST_PATH_IMAGE008
Average value of (i), i.e.
Figure 177729DEST_PATH_IMAGE010
The process proceeds to step S200.
2. The method of claim 1, wherein said step S300 further comprises the step of iteratively solving for said coupling between layers of said multilayer very large scale integrated circuitlIteration of layers:
s310, in the case of the last iteration, the electromagnetic field of each layer of PCB is recorded as
Figure 919157DEST_PATH_IMAGE011
；
S320, setting
Figure 957521DEST_PATH_IMAGE012
；
S330, calculatingmCurrent source pair with distributed layerslInfluence of the layerG ml Based onG ml Update the firstlFirst among source items of a layermThe portion of the source layer affecting it, tolApplying two-dimensional finite element to the layer to calculate the electromagnetic field distribution so as to update the electromagnetic field and current distribution of the layer, and calculating the change of the electromagnetic field of the layerdE ml At this time it islElectromagnetic field of the layer becomes
Figure 45693DEST_PATH_IMAGE013
；
S340、l=l+1, return to step S330 until
Figure 751481DEST_PATH_IMAGE014
。
3. The iterative solution method for interlayer coupling in a multilevel very large scale integrated circuit of claim 1, wherein the iterative solution method is based on a specific structure of a hierarchy of integrated circuitsmSource layer at the secondlInfluence of layer G ml Can be decomposed to be located at the secondmSource layer point current source atlSuperposition of electric field expressions generated by layers, the firstmSource layer point current source atlThe electric field expression generated by the layers is a special analytical expression given by using a dyadic Green function, and the current sources of the multilayer integrated circuit are distributed in a layered mode, namely the current density distributed on each metal layer of the integrated circuit layout with a complex shape is only equal to that of the current density distributed on each metal layer of the integrated circuit layout with a complex shapexAndythe axial direction is related tozAxial direction independence, current density distribution of onlyx，yAs a function of (c).
4. The iterative solution method for interlayer coupling in a multilevel very large scale integrated circuit of claim 3, wherein said first stepmSource layer at the secondlInfluence of layer G ml Can be decomposed to be located at the secondmSource layer point current source atlThe specific method for superposing the electric field expressions generated by the layers is as follows: will be located at the firstmSource layer point current source atlThe electric field expression generated by the layer is used as the integrand of two-dimensional Gaussian integration, and the first calculation is based on the linear superposition principle of the fieldmField generated at the same position by simple polygonal planar current source of source layer, second in two-dimensional plane SmSource layer at the secondlThe layer-generated field is calculated by the two-dimensional gaussian integral:
Figure 194970DEST_PATH_IMAGE015
；
wherein the content of the first and second substances,E(x,y,z) At any point in space for the current source in the two-dimensional plane Sx,y,z) The field that is generated is,
Figure 720629DEST_PATH_IMAGE016
is at any position in the two-dimensional surface S
Figure 533864DEST_PATH_IMAGE017
At any point in space (x,y,z) The expression of the dyadic green function of the generated field,
Figure 376049DEST_PATH_IMAGE018
representing a gaussian integration point corresponding to a two-dimensional gaussian integration in the two-dimensional plane S,p,qrespectively representu,vIn the first directionpA first, aqThe number of the Gaussian integration points is equal to the total number of the points,
Figure 944434DEST_PATH_IMAGE019
is the weighting factor corresponding to the gaussian integration point.
5. The iterative solution method for interlayer coupling in a multi-layer very large scale integrated circuit of claim 1, wherein the influence value according to the dyadic Green function is used in the iterative processGDetermines a layer that can be ignored, adaptively adjustsmThe source layer adjacent theretolExtent of influence exerted by the layer
Figure 629493DEST_PATH_IMAGE002
。
6. An iterative device for interlayer coupling of a multilayer ultra-large scale integrated circuit is characterized by comprising an action layer iteration module, a source item updating module and an electromagnetic field change quantity updating module,N+1 LSI layer, each layer numbered
Figure 354742DEST_PATH_IMAGE020
；
The action layer iteration module is used for iteratively updating the action layer of the source layer
Figure 503963DEST_PATH_IMAGE002
And is provided with the firstmActive layer of source layer
Figure 915353DEST_PATH_IMAGE002
Is divided bymOthers of the source layerNLayers of integrated circuits, i.e.
Figure 963075DEST_PATH_IMAGE003
；
The source layer iteration module is used for updating an iteration source layer;
the source item updating module is used for calculating the influence of the updated source layer on all other layers by utilizing a dyadic Green function when the source layer is updatedG ml ；
The change amount updating module of the electromagnetic field is based onG ml Update the firstlFirst among source items of a layermThe portion of the source layer that affects it is calculated by two-dimensional finite elementslThe electromagnetic field distribution of the layer is updated to update the electromagnetic field and current distribution of the layer, thereby calculating the change amount of the electromagnetic field of the layerdE ml 。
7. The iterative means for interlayer coupling in a multilevel VLSI according to claim 6, wherein said active layer iteration module selects all of the satisfied T cellsG ml |<GOf the conditionG ml Is marked asG thredshold WhereinGCalculating all the effective influence values of the dyadic Green functionG thredshold Middle distance layermNearest layerl near Number of layers of
Figure 180429DEST_PATH_IMAGE008
Is marked as
Figure 121840DEST_PATH_IMAGE021
Update
Figure 349471DEST_PATH_IMAGE022
Is composed of
Figure 743544DEST_PATH_IMAGE008
Average value of (i), i.e.
Figure 764589DEST_PATH_IMAGE023
。
8. The iterative means for interlayer coupling in a multilevel very large scale integrated circuit of claim 7, wherein the effective influence value of the dyadic Green functionGThe solution of (c) is as follows: selectingG ml Maximum value ofG maxAnd minimum valueG minCalculating the effective influence value of the dyadic Green function
Figure 170294DEST_PATH_IMAGE024
WhereinthredsholdThe method is a preset discarding threshold value influenced by the dyadic Green function.
9. The iterative means of interlayer coupling in a multilevel very large scale integrated circuit of claim 6, wherein the iterative means is based on a particular structure of a hierarchy of integrated circuitsmSource layer at the secondlInfluence of layer G ml Can be decomposed to be located at the secondmSource layer point current source atlSuperposition of electric field expressions generated by layers, the firstmSource layer point current source atlThe electric field expression generated by the layers is a special analytical expression given by using a dyadic Green function, and the current sources of the multilayer integrated circuit are distributed in a layered mode, namely the current density distributed on each metal layer of the integrated circuit layout with a complex shape is only equal to that of the current density distributed on each metal layer of the integrated circuit layout with a complex shapexAndythe axial direction is related tozAxial direction independence, current density distribution of onlyx，yAs a function of (c).
10. The iterative means of interlayer coupling in a multilevel very large scale integrated circuit of claim 9, wherein said first stepmSource layer at the secondlInfluence of layer generationG ml Can be decomposed to be located at the secondmOf the source layerThe point current source is atlThe specific method for superposing the electric field expressions generated by the layers is as follows: will be located at the firstmSource layer point current source atlThe electric field expression generated by the layer is used as the integrand of two-dimensional Gaussian integration, and the first calculation is based on the linear superposition principle of the fieldmField generated at the same position by simple polygonal planar current source of source layer, second in two-dimensional plane SmSource layer at the secondlThe layer-generated field is calculated by the two-dimensional gaussian integral:
Figure 189066DEST_PATH_IMAGE025
；
wherein the content of the first and second substances,E(x,y,z) At any point in space for the current source in the two-dimensional plane Sx,y,z) The field that is generated is,
Figure 398330DEST_PATH_IMAGE016
is at any position in the two-dimensional surface S
Figure 269072DEST_PATH_IMAGE017
At any point in space (x,y,z) The expression of the dyadic green function of the generated field,
Figure 512971DEST_PATH_IMAGE018
representing a gaussian integration point corresponding to a two-dimensional gaussian integration in the two-dimensional plane S,p,qrespectively representu,vIn the first directionpA first, aqThe number of the Gaussian integration points is equal to the total number of the points,
Figure 437065DEST_PATH_IMAGE019
is the weighting factor corresponding to the gaussian integration point.
CN202110425193.5A 2021-04-20 2021-04-20 Iterative solution method and device for interlayer coupling of multilayer very large scale integrated circuit Active CN112989735B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110425193.5A CN112989735B (en) 2021-04-20 2021-04-20 Iterative solution method and device for interlayer coupling of multilayer very large scale integrated circuit

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110425193.5A CN112989735B (en) 2021-04-20 2021-04-20 Iterative solution method and device for interlayer coupling of multilayer very large scale integrated circuit

Publications (2)

Publication Number Publication Date
CN112989735A true CN112989735A (en) 2021-06-18
CN112989735B CN112989735B (en) 2021-09-24

Family

ID=76341326

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110425193.5A Active CN112989735B (en) 2021-04-20 2021-04-20 Iterative solution method and device for interlayer coupling of multilayer very large scale integrated circuit

Country Status (1)

Country Link
CN (1) CN112989735B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114048661A (en) * 2021-11-29 2022-02-15 中南大学 Method and device for treating DC point source potential and electric field under laminar medium

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6064808A (en) * 1997-08-01 2000-05-16 Lucent Technologies Inc. Method and apparatus for designing interconnections and passive components in integrated circuits and equivalent structures by efficient parameter extraction
CN103793603A (en) * 2014-01-24 2014-05-14 同济大学 Accurate electromagnetic analysis method for chip electronic packaging structure
CN110162831A (en) * 2019-04-10 2019-08-23 华中科技大学 A kind of numeric type integrated photonic device emulation mode and system
CN111737947A (en) * 2020-08-06 2020-10-02 北京智芯仿真科技有限公司 Integrated circuit full-wave IBIS model extraction method and device based on field-circuit coupling
CN111931457A (en) * 2020-09-27 2020-11-13 北京智芯仿真科技有限公司 Multilayer integrated circuit electromagnetic field calculation method and device based on mixed order finite element

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6064808A (en) * 1997-08-01 2000-05-16 Lucent Technologies Inc. Method and apparatus for designing interconnections and passive components in integrated circuits and equivalent structures by efficient parameter extraction
CN103793603A (en) * 2014-01-24 2014-05-14 同济大学 Accurate electromagnetic analysis method for chip electronic packaging structure
CN110162831A (en) * 2019-04-10 2019-08-23 华中科技大学 A kind of numeric type integrated photonic device emulation mode and system
CN111737947A (en) * 2020-08-06 2020-10-02 北京智芯仿真科技有限公司 Integrated circuit full-wave IBIS model extraction method and device based on field-circuit coupling
CN111931457A (en) * 2020-09-27 2020-11-13 北京智芯仿真科技有限公司 Multilayer integrated circuit electromagnetic field calculation method and device based on mixed order finite element

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114048661A (en) * 2021-11-29 2022-02-15 中南大学 Method and device for treating DC point source potential and electric field under laminar medium

Also Published As

Publication number Publication date
CN112989735B (en) 2021-09-24

Similar Documents

Publication Publication Date Title
CN112989677B (en) Iteration method and device for accumulation calculation of interlayer coupling part of large-scale integrated circuit
CN112836466B (en) Two-dimensional fast iteration method and device for three-dimensional large-scale integrated circuit electromagnetic response
Reddy et al. Fast RCS computation over a frequency band using method of moments in conjunction with asymptotic waveform evaluation technique
Liu et al. Existence of ${\cal H} $-Matrix Representations of the Inverse Finite-Element Matrix of Electrodynamic Problems and ${\cal H} $-Based Fast Direct Finite-Element Solvers
Wei et al. An optimized higher order PML in domain decomposition WLP-FDTD method for time reversal analysis
CN112989735B (en) Iterative solution method and device for interlayer coupling of multilayer very large scale integrated circuit
CN112989676B (en) Iterative method and device for instantly updating current distribution of integrated circuit by interlayer coupling
CN112818633B (en) Iterative method and device for integrated circuit current distribution dynamically applied by interlayer coupling
CN112232001B (en) Self-adaptive determination method and system for ultra-wideband resonance response of integrated circuit
CN112989750B (en) Method and device for determining space electromagnetic radiation of multilayer integrated circuit
CN112257372B (en) Method and system for extracting impedance network model of integrated circuit
CN112290955A (en) Grid node coding method and system based on integrated circuit impedance network extraction
CN112989678B (en) Coarse grain parallel iteration method and device for integrated circuit interlayer coupling part accumulation
Shao et al. Signal integrity analysis of high-speed interconnects by using nonconformal domain decomposition method
CN112818585B (en) Method and device for dividing iterative computation parallel particles of integrated circuit interlayer coupling
Guo et al. Extrapolation with Range Determination of 2D Spectral Transposed Convolutional Neural Network for Advanced Packaging Problems
CN112818584B (en) Space electromagnetic radiation computing system and method for integrated circuit
Kolundzija et al. Matrix equilibration in method of moment solutions of surface integral equations
CN112989756B (en) Coarse grain parallel iteration method and device for dynamically applying coupling between layers of integrated circuit
Hollander et al. Adaptive multilevel nonuniform grid algorithm for the accelerated analysis of composite metallic–dielectric radomes
Gao et al. Efficient Full-Wave Simulation of Large-Scale Metasurfaces and Metamaterials
CN113591423B (en) Full-wave electromagnetic simulation method and system for integrated circuit under lossy frequency dispersion medium
CN112989675B (en) Coarse grain parallel iteration method and device for integrated circuit interlayer coupling instant update
Panayappan et al. A technique for handling multiscale electromagnetic problems using the finite difference time domain (FDTD) algorithm
Junkin et al. A robust 3D mesh generator for the Dey-Mittra conformal FDTD algorithm

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant