CN107305592B - Numerical simulation method applied to electromagnetic transient analysis of switch circuit - Google Patents

Abstract

The invention relates to a numerical simulation method applied to electromagnetic transient analysis of a switching circuit, which comprises the following steps: writing a system state equation of the circuit in a column, and converting the system state equation into a state equation in a first order form through variable substitution; modifying a system state equation of the circuit at the moment of switching action; performing iterative solution on the first-order state equation by using a direct integration method; reinitializing the potential of the circuit node at the moment of switching action; eliminating numerical oscillation by using minimum step length backward Euler operation at the switching time of the switch state; and (4) proving the equivalence of the switching action time. The technical scheme provided by the invention can eliminate the numerical value oscillation at the switching moment, and has the characteristics of good numerical value stability and high numerical value precision (second-order precision).

Description

Numerical simulation method applied to electromagnetic transient analysis of switch circuit

Technical Field

The invention relates to a numerical simulation method, in particular to a numerical simulation method applied to electromagnetic transient analysis of a switching circuit.

Background

Electromagnetic transient analysis is an important component of power system simulation and is a necessary tool for understanding transient complex behaviors of a power system. Most of the currently used electromagnetic transient simulation software adopts a node voltage method as a bottom-layer algorithm. The node voltage method is characterized in that a numerical integration method is adopted to differentiate a characteristic equation of a dynamic element in the system to obtain an equivalent Noton equivalent circuit in the form that an equivalent calculated conductance is connected with a historical current source in parallel, a node conductance matrix of the whole circuit system is further obtained, and the instantaneous value of each node voltage in the system is obtained through solving. When a switching circuit is analyzed by using a node analysis method, there is a major problem that strong numerical oscillation occurs at the time of switching operation. Although each node analysis software adopts a certain numerical technique to eliminate oscillation, the numerical oscillation eliminating techniques are often first-order linear interpolation techniques, which may cause a decrease in numerical accuracy. It is still required to use a very small step size to ensure sufficient calculation accuracy.

Disclosure of Invention

The invention aims to provide a numerical simulation method applied to electromagnetic transient analysis of a switching circuit, which does not generate numerical oscillation at the moment of switching action and has the characteristics of good numerical stability and high numerical precision (second-order precision).

The purpose of the invention is realized by adopting the following technical scheme:

the invention provides a numerical simulation method applied to electromagnetic transient analysis of a switching circuit, which is improved in that the method comprises the following steps:

(1) writing a system state equation of the circuit in a column, and converting the system state equation into a state equation in a first order form through variable substitution;

(2) modifying a system state equation of the circuit at the moment of switching action;

(3) performing iterative solution on the first-order state equation by using a direct integration method;

(4) reinitializing the potential of the circuit node at the moment of switching action;

(5) eliminating numerical oscillation by using minimum step length backward Euler operation at the switching time of the switch state;

(6) and (4) proving the equivalence of the switching action time.

Further, in the step (1), the system state equation expression of the circuit is as follows:

wherein:

is a column vector of order n, each element value representing the potential of the corresponding node; k_CIs a capacitance coefficient matrix; k_RIs a resistivity matrix; k_LIs an inductance matrix; i is an input vector in the circuit and represents a power supply in the circuit; the three coefficient matrixes are square matrixes of n rows and n columns, and all the three coefficient matrixes meet the requirements of positive nature, symmetry and sparsity;

through variable substitution, the system state equation (r) is transformed into a first order state equation as follows:

wherein:

K₁、K₂all represent a coefficient matrix; x is a state variable; e is the unit diagonal matrix.

Further, in the step (2), if the ideal switch connected between the two nodes is operated at the switching operation time, the ideal switch is operated at the switching operation time t_dResistance coefficient matrix K in modified equation_R(ii) a And modifying the system state equation of the circuit by using an expression (R), wherein the expression (R) is as follows:

where Turn On indicates that the ideal switch is at t_dOpening at any moment; turn Off indicates that the ideal switch is at t_dDisconnection at any moment; except for K_r(i,i)＝R_s ^-1，K_r(i,j)＝-R_s ^-1，K_r(j,i)＝-R_s ^-1And K_r(j,j)＝R_s ^-1Other than, K_rAll other elements in the material are zero, 0V ideal electricityInternal resistance of voltage source R_sTends to zero; k_R,+Is a matrix of resistivity before switching action, K_R,-Is a resistance coefficient matrix after switching action; k_rIs an N-th order matrix, K_r(i, i) represents the element in the ith row and ith column in the matrix; k_r(i, j) represents the element in the ith row and jth column of the matrix; k_r(j, i) represents the element in the jth row and ith column in the matrix; k_r(j, j) represents the element in the jth row and jth column in the matrix;

further, in the step (3), for the first-order dynamic equation of the standard form, iterative solution is performed by using a direct integration method, and an iterative format of the direct integration method is represented by the following formula (c):

wherein: Δ t is the time step between the nth time instant and the time instant n + 1; beta is an optional coefficient, when the value of beta is 0, the method is called a forward difference method, the stability is stable, and the numerical precision is first-order precision; when the value of beta is 0.5, the method is called a Crank-Nicolson method, the stability is unconditionally stable, and the numerical precision is second-order precision; when the value of beta is 1, the method is called a posterior difference method, the stability is unconditionally stable, and the numerical precision is first-order precision; x is the number of_nRepresents the system state variable at the nth step, namely, the time t is n delta t; x is the number of_n+1Represents the system state variable at the n +1 th step, i.e. at the time t ═ n +1) Δ t; r_nAnd R_n+1The system input vectors of the nth step, i.e., time t ═ n Δ t and the nth +1 th step, i.e., time t ═ n +1) Δ t, are represented, respectively.

Further, the step (4) includes:

at the time t of switching action_dThe capacitor is equivalent to a voltage source, the inductor is equivalent to a current source, and the current source is based on t_d-Capacitance voltage and inductance current calculation t at time_d+Potential distribution at a time; t is t_d-Represents t_dBefore the momentAt a certain time, and t_d-t_d-Tends to zero; t is t_d+Represents t_dA certain time after the time, and t_d+-t_dTends to zero;

if the capacitance and inductance dynamic elements in the circuit are treated equivalently, only a resistance branch and a current source branch exist in the circuit to obtain a static circuit equation, and the equation in the equation (c) is solved

Completing the initialization of the circuit after the switching action;

wherein:

represents t_d+A node potential at a time; k is a coefficient matrix; i is an input vector in the circuit and represents a power supply in the circuit; the coefficient matrix K contains information K of resistance branches and voltage source resistance branches_R,d+Internal resistance information K of equivalent voltage source of sum capacitor_RC；

K＝K_R,d++K_RCWherein

For a certain branch circuit which is connected between the node I and the node j and has the capacitance C, the coefficient matrix K_rcIncluding its equivalent to an ideal voltage source, with respect to its internal resistance R_sBranch information that tends to zero; matrix K_rcMiddle removing K_rc(i,i)＝R_s ^-1，K_rc(i,j)＝-R_s ^-1，K_rc(j,i)＝-R_s ^-1，K_rc(j,j)＝R_s ^-1Except that, other elements are zero; k_rcIs an n-th order matrix, K_rc(i, i) represents a matrix K_rcRow i and column i; k_rc(i, j) represents momentMatrix K_rcRow i and column j; k_rc(j, i) represents a matrix K_rcThe element in the jth row and ith column; k_rc(j, j) represents a matrix K_rcRow j and column j;

vector I includes all current sources in the original circuit and the Nonton equivalent current source information I of the ideal voltage source_S,d+Besides, the current source also comprises information I of an inductance equivalent current source_L,d-And information I of the equivalent ideal voltage source of the capacitor_SC,d-As shown in formula xi:

I＝I_S,d++I_L,d-+I_SC,d- ⑧

wherein:

for a certain branch circuit which is connected between the node I and the node j and has inductance L, I_lIncluding information of its equivalent current source; except that I_l,d-(i)＝I_L,d-And I_l,d-(j)＝-I_L,d-Outer, I_l,d-All other elements in (A) are zero; meanwhile, for a certain branch which is connected between the node I and the node j and has the capacitance of C, I_scIncluding its equivalent norton current source information; i is_sc,d-In, except for I_sc,d-(i)＝I_SC，d-and I_sc,d-(j)＝-I_SC，d-In addition, the other elements are all zero, wherein I_SC，d-Represents t_d-The current value of the time capacitor equivalent Norton current source;

substituting the formula (c) and the formula (b) into the formula (c), and arranging the formula (c) into a format of the formula (c):

wherein:

represents t_d+The potential of each node at a time; k_R,d+Represents t_d+A resistance coefficient matrix which constantly contains information of the resistance branch and the voltage source internal resistance branch; i is_S,d+Represents t_d+Current source information vectors at time; i is_SC,d-Represents t_d-Time capacitance equivalent Norton current source information vector; l _ relative represents the inductance-related term; c _ relative represents a capacitive element-related term.

Further, the step (5) includes:

at the time t of switching action_dUsing a backward euler operation:

if at the moment t of switching action_dAt a very small time step Δ t₀Applying a backward Euler operation to the system state equation (r), where at₀＝t_d+-t_d-Using t_d-State quantity x of time_d-To find t_d+State quantity x of time_d+The calculation expression is as follows:

wherein the content of the first and second substances,

Ψ is a time integral quantity of voltage at each node, defined as

And

R_d+an input vector matrix representing the time after the switching action; k_2,d+Is a coefficient matrix; further expand formula (R) to

And formula

Taking into account Δ t₀Tends to 0, formula

Become formula

The integral quantity of the node potential along with the time is not changed before and after the switching action;

wherein: Ψ_d+And Ψ_d-Respectively representing the voltage time integral quantity of each node after the switching action and the voltage time integral quantity of each node before the switching action;

general formula

Substituted type

And arranged into

In the form of:

if formula

And formula ninthly is equivalent, stating: at the time t of switching action_dAnd one node potential reinitialization operation based on the classical theory is completed by one small-step backward Euler operation.

Further, the step (6) comprises the steps of:

<1> proof of equivalence of inductance-related terms:

for inductive branches, there are

The following holds true:

or written in a matrix format, e.g. formula

Shown in the figure:

continuing to expand to

Formula (II):

because each inductance branch has a formula

An inductance related term in the formula (ninx)

And formula

Inductance related term in

Are equivalent; i is_l,d-Represents t_d-The current value on the inductance branch at the moment; l is^-1To represent

L is the inductance value of the inductance branch; Ψ_d-(i) And Ψ_d-(j) Represents t_d-Voltage time integral quantity of the node i and the node j at the moment;

<2> proof of equivalence of capacitance-related terms:

considering the Nonton equivalent current source of the capacitance branch, the capacitance related term in the formula ninu is developed into the formula

And formula

General formula

The capacitance related term in

Writing each capacitor branch into an expansion format as shown in the formula

And

shown in the figure:

general formula

The capacitance related term in

Writing each capacitor branch into an expansion format as shown in the formula

And

shown in the figure:

due to the fact that

Represents t_d-The voltage value on the capacitor at the moment, while taking into account R_SAnd Δ t₀The fact that all tend to zero, the capacitance related term in the formula ninthly

And formula

The capacitance related term in

Equivalence; u shape_c,d-Represents t_d-At the moment, the voltage value on the capacitor branch circuit;

and

respectively representing potential values of the node i and the node j;

the switching operation uses a minimum step length backward Euler operation at the moment, and is equivalent to the operation of reinitializing the node potential according to the capacitor voltage and the inductance current in the circuit.

Compared with the closest prior art, the technical scheme provided by the invention has the following excellent effects:

because the current direction on the inductance element is always positive, when the output current of the current source passes through zero, the change rate of the inductance current jumps from negative to positive once, thereby causing the voltage at two ends of the inductance

A transition occurs every 10ms (as in fig. 4). For nodal analysis using the trapezoidal integration algorithm, numerical oscillations of the inductor voltage will inevitably occur. The PSCAD software uses an interpolation technique to eliminate numerical oscillation, and the interpolation technique is a linear interpolation technique with only first-order precision, which loses the numerical precision of calculation. In other words, although PSCAD uses a second-order precision trapezoidal method to construct equations for discrete dynamic elements, its interpolation technique at the time of switching operation reduces numerical precision to the first order. Therefore, if it is desired to maintain good numerical accuracy when the PSCAD software simulates the switching circuit, a smaller time step needs to be selected.

In the state equation method, because the post-difference operation with extremely small step length is used for carrying out reinitialization operation on the node potential at the moment of switching action, the second-order numerical accuracy can still be kept when the trapezoidal method is reused in a linear section. Therefore, compared with a node analysis method for establishing an equation based on a trapezoidal integral method, the method has the advantage of high numerical precision.

Drawings

FIG. 1 is an equivalent circuit diagram of a capacitor at the moment of switching operation provided by the present invention;

FIG. 2 is an equivalent circuit diagram of the inductor at the moment of switching operation provided by the present invention;

FIG. 3 is a circuit diagram of a bridge current source with a power frequency of 50Hz provided by the present invention;

fig. 4 is a schematic diagram showing a comparison between a simulation result of the PSCAD software provided by the present invention and a simulation result of the method of the present invention.

Detailed Description

The following describes embodiments of the present invention in further detail with reference to the accompanying drawings.

The following description and the drawings sufficiently illustrate specific embodiments of the invention to enable those skilled in the art to practice them. Other embodiments may incorporate structural, logical, electrical, process, and other changes. The examples merely typify possible variations. Individual components and functions are optional unless explicitly required, and the sequence of operations may vary. Portions and features of some embodiments may be included in or substituted for those of others. The scope of embodiments of the invention encompasses the full ambit of the claims, as well as all available equivalents of the claims. Embodiments of the invention may be referred to herein, individually or collectively, by the term "invention" merely for convenience and without intending to voluntarily limit the scope of this application to any single invention or inventive concept if more than one is in fact disclosed.

The invention relates to a numerical simulation method applied to electromagnetic transient analysis of a switching circuit, which comprises the following steps:

(1) writing a system state equation of the circuit in a column, and converting the system state equation into a state equation in a first order form through variable substitution;

the system state equation expression of a general linear circuit (a circuit having n +1 nodes, in which node 0 is defined as a reference ground node, and composed of basic elements such as a resistor, a capacitor, an inductor, and a power supply) is as follows:

wherein:

is a column vector of order n, each element value representing the potential of the corresponding node; k_CIs a capacitance coefficient matrix; k_RIs a resistivity matrix; k_LIs an inductance matrix; i is an input vector in the circuit and represents a power supply in the circuit; the three coefficient matrixes are square matrixes of n rows and n columns, and all the three coefficient matrixes meet the requirements of positive nature, symmetry and sparsity;

through variable substitution, the system state equation (r) is transformed into a first order state equation as follows:

wherein:

K₁、K₂all represent a coefficient matrix; x is a state variable; e is the unit diagonal matrix.

(2) Modifying a system state equation of the circuit at the moment of switching action;

method for rapidly modifying the matrix of state equation coefficients at the moment of a switching action when an ideal switching element in a circuit network is actuated (switched on or switched off):

if an ideal switch connected between two nodes is operated at the time of switching operation, the ideal switch is operated at the time of switching operation t₀Resistance coefficient matrix K in modified equation_R(ii) a And modifying the system state equation of the circuit by using an expression (R), wherein the expression (R) is as follows:

where Turn On indicates that the ideal switch is at t_dOpening at any moment; turn Off indicates that the ideal switch is at t_dDisconnection at any moment; except for K_r(i,i)＝R_s ^-1，K_r(i,j)＝-R_s ^-1，K_r(j,i)＝-R_s ^-1And K_r(j,j)＝R_s ^-1Other than, K_rThe other elements in the voltage source are zero, and the internal resistance R of the 0V ideal voltage source_sTends to zero; k_R,+Is a matrix of resistivity before switching action, K_R,-Is a resistance coefficient matrix after switching action; k_rIs an N-th order matrix, K_r(i, i) represents the element in the ith row and ith column in the matrix; k_r(i, j) represents the element in the ith row and jth column of the matrix; k_r(j, i) represents the element in the jth row and ith column in the matrix; k_r(j, j) represents the element in the jth row and jth column in the matrix;

(3) performing iterative solution on the first-order state equation by using a direct integration method;

for the first-order dynamic equation of the standard form II, iterative solution is carried out by using a direct integration method, and the iterative format of the direct integration method is represented by the formula V:

wherein: Δ t is the time step between the nth time instant and the time instant n + 1;x_nrepresents the system state variable at the nth step, namely, the time t is n delta t; x is the number of_n+1Represents the system state variable at the n +1 th step, i.e. at the time t ═ n +1) Δ t; r_nAnd R_n+1The system input vectors of the nth step, namely, the time t ═ n Δ t and the n +1 th step, namely, the time t ═ n +1) Δ t are respectively represented; β is an alternative coefficient, and when β takes different values, the associated names are as described in table 1 below:

TABLE 1 selection and numerical stability of beta

Among them, the trapezoidal method has a numerical accuracy of the second order and is unconditionally stable, and therefore is widely used in linear problem analysis. However, in the switching circuit analysis, the trapezoidal method generates strong numerical oscillation at the time of switching operation.

(4) Reinitializing the potential of the circuit node at the moment of switching action;

when the switching action causes the circuit topology to be at t_dWhen the time is changed, the node potential in the time is possibly at t_dJump occurs before and after the moment (as in fig. 3, the potential of the node No. 3 jumps at the moment when the diode is turned on or off), and when the trapezoidal method is used for calculation, numerical oscillation of the potential of the relevant node is caused. This requires reinitializing the potentials of the nodes of the circuit at the switching time.

According to the circuit theory, before and after the switching action, the capacitor branch has the characteristic of keeping the potential difference of two ends unchanged, and the inductor branch has the characteristic of keeping the current of the branch unchanged. Thus can be at t_dAt the moment, the capacitance is equivalent to a voltage source (as shown in FIG. 1), the inductance is equivalent to a current source (as shown in FIG. 2), and the t is the basis_d-Capacitance voltage and inductance current calculation t at time_d+The potential distribution at the moment. At the time t of switching action_dThe capacitor is equivalent to a voltage source, the inductor is equivalent to a current source, and the current source is based on t_d-Potential distribution at a time; t is t_d-Represents t_dA certain time before the time, and t_d-t_d-Tends to zero; t is t_d+Represents t_dA certain time after the time, and t_d+-t_dTends to zero;

if all the dynamic elements such as capacitance, inductance and the like in the circuit are subjected to equivalent processing, only a resistance branch and a current source branch exist in the circuit, and a static circuit equation is obtained. To solve the formula

The initialization of the circuit after the switching action can be completed.

Wherein:

represents t_d+A node potential at a time; k is a coefficient matrix; i is an input vector in the circuit and represents a power supply in the circuit; the coefficient matrix K contains information K of resistance branches and voltage source resistance branches_R,d+Internal resistance information K of equivalent voltage source of sum capacitor_RC。

(5) Eliminating numerical oscillation by using minimum step length backward Euler operation at the switching time of the switch state;

the coefficient matrix K contains information K of resistance branches and voltage source resistance branches_R,d+Internal resistance information K of equivalent voltage source of sum capacitor_RC；

For a certain branch circuit which is connected between the node I and the node j and has the capacitance C, the coefficient matrix K_rcIncluding its equivalent to an ideal voltage source, with respect to its internal resistance R_sBranch information that tends to zero; matrix K_rcMiddle removing K_rc(i,i)＝R_s ^-1，K_rc(i,j)＝-R_s ^-1，K_rc(j,i)＝-R_s ^-1，K_rc(j,j)＝R_s ^-1Except that, other elements are zero; k_rcIs an n-th order matrix, K_rc(i, i) represents a matrix K_rcRow i and column i; k_rc(i, j) represents a matrix K_rcRow i and column j; k_rc(j, i) represents a matrix K_rcThe element in the jth row and ith column; k_rc(j, j) represents a matrix K_rcRow j and column j;

vector I includes all current sources in the original circuit and the Nonton equivalent current source information I of the ideal voltage source_S,d+Besides, the current source also comprises information I of an inductance equivalent current source_L,d-And information I of the equivalent ideal voltage source of the capacitor_SC,d-As shown in formula xi:

I＝I_S,d++I_L,d-+I_SC,d- ⑧

wherein:

for a certain branch circuit which is connected between the node I and the node j and has inductance L, I_lIncluding information of its equivalent current source; except that I_l,d-(i)＝I_L,d-And I_l,d-(j)＝-I_L,d-Outer, vector I_l,d-While for a branch with a capacitance of C connected between node I and node j, I is zero_scIncluding its equivalent norton current source information; i is_SC，d-In, except for I_sc，d-(i)＝I_SC，d-And I_sc，d-(j)＝-I_SC，d-Besides, other elements are zero;

substituting the formula (c) and the formula (b) into the formula (c), and arranging the formula (c) into a format of the formula (c):

wherein:

represents t_d+The potential of each node at a time; k_R,d+Represents t_d+A resistance coefficient matrix which constantly contains information of the resistance branch and the voltage source internal resistance branch; i is_S,d+Represents t_d+Current source information vectors at time; i is_SC,d-Represents t_d-Time capacitance equivalent Norton current source information vector; l _ relative represents the inductance-related term; c _ relative represents a capacitive element-related term.

At the time t of switching action_dUsing a backward euler operation:

if at the moment t of switching action_dAt a very small time step Δ t₀Applying a backward Euler operation to the system state equation (r), where at₀＝t_d+-t_d-Using t_d-State quantity x of time_d-To find t_d+State quantity x of time_d+The calculation expression is as follows:

wherein the content of the first and second substances,

Ψ is a time integral quantity of voltage at each node, defined as

And

R_d+an input vector matrix representing the time after the switching action; k_R,d+Is a resistivity matrix after switching.

Further expand formula (R) to

And formula

Taking into account Δ t₀Tends to 0, formula

Become formula

The integral quantity of the node potential along with the time is not changed before and after the switching action;

wherein: Ψ_d+And Ψ_d+Respectively representing the voltage time integral quantity of each node after the switching action and the voltage time integral quantity of each node before the switching action;

general formula

Substituted type

And arranged into

In the form of:

if formula

And formula ninthly is equivalent, stating: at the time t of switching action_dAnd one node potential reinitialization operation based on the classical theory is completed by one small-step backward Euler operation.

(6) The equivalence of the switching action time proves that:

<1> proof of equivalence of inductance-related terms:

for inductive branches, there are

The following holds true:

or written in a matrix format, e.g. formula

Shown in the figure:

continuing to expand to

Formula (II):

because each inductance branch has a formula

An inductance related term in the formula (ninx)

And formula

Inductance related term in

Are equivalent; i is_l,d-Represents t_d-The current value on the inductance branch at the moment; l is^-1To represent

L is the inductance value of the inductance branch. Ψ_d-(i) And Ψ_d-(j) Represents t_d-And voltage time integral quantities of the node i and the node j at the moment.

<2> proof of equivalence of capacitance-related terms:

considering the Nonton equivalent current source of the capacitance branch, the capacitance related term in the formula ninu is developed into the formula

And formula

General formula

The capacitance related term in

Writing each capacitor branch into an expansion format as shown in the formula

And

shown in the figure:

general formula

The capacitance related term in

Writing each capacitor branch into an expansion format as shown in the formula

And

shown in the figure:

due to the fact that

Represents t_d-The voltage value on the capacitor at the moment, while taking into account R_SAnd Δ t₀The fact that all tend to zero, the capacitance related term in the formula ninthly

And formula

The capacitance related term in

Equivalence; u shape_c,d-Represents t_d-At the moment, the voltage value on the capacitor branch circuit;

and

respectively representing the potential values of the node i and the node j.

This proves that at the moment of switching action: the use of the minimum step-length back-difference method is equivalent to the operation of reinitializing the node potential according to the capacitor voltage and the inductive current in the circuit.

Examples

Fig. 3 shows a bridge current source rectifier circuit, and table 2 is a form in which the circuit topology is stored in a computer. The simulation analysis was performed using the patented method with the basic step size set to 2 ms.

Table 2 example leg information

Firstly, codes of a column writing method for executing a circuit network state equation on the No. 1 to No. 5 branches in the table 3 are written to obtain a coefficient matrix shown in the formula

Thereby obtaining a system equation of the circuit

And equation of state

Wherein:

e is the unit diagonal matrix.

Defining the state quantity at zero time as x⁰＝(0,0,0,0,900,0,0,0,0,0)^T。

5.1 first half cycle

At time 0+, D2 and D3 will turn on while D1 and D4 remain off due to the unidirectional turn-on characteristics of the ideal diode. Referring to the rule described in the formula, the initial state resistivity matrix K_RModified to K_R1Of formula

Wherein:

5.1.1 Single step post Difference calculation

Since D2 and D3 are turned on at time zero, it is necessary to apply a backward difference method once to obtain the state quantity at time 0+ in order to avoid numerical oscillation. The time step of this backward difference operation takes Δ t₀＝10^-7s, and will be Δ t₀And β ═ 1 substitution, giving t ═ 10^-7system state quantity at s moment

As shown in

Wherein

5.1.2 trapezoidal iterative computation

Then adjusting the time step to be 0.002-10 ^ t^-7s, according to t 10^-7system state quantity at s moment

Obtaining the state quantity x when t is 2ms by using the formula_0.002. The step length was adjusted to Δ t of 0.002s, using the formula

Sequentially meterAnd calculating the system state quantity at 4ms, 6ms, 8ms and 10 ms.

5.2 half period

When t is 10ms, D2 and D3 are turned off, D1 and D4 are turned on, and the coefficient matrix K of the first half period is determined according to the rule described in the formula₂₁Modified to K for the second half cycle₂₂Of formula

And

wherein

As in the first half cycle, at the time t of circuit topology switching is 10ms, Δ t is taken₀＝10^-7s, yielding t ═ 0.01+10^-7system state quantity at s moment

And then, after one step is adjusted, fixing the step to 2ms, and sequentially calculating to obtain the system state quantities at 12ms, 14ms, 16ms, 18ms and 20 ms.

And (4) analyzing results:

the inductor voltages obtained by the simulation are shown in table 3 and fig. 4. Meanwhile, simulation results of the PSCAD software are given for comparison, wherein the accurate solution is obtained by using the PSCAD software at the step size of 1 mu s.

TABLE 3 inductive Voltage

Although the present invention has been described in detail with reference to the above embodiments, those skilled in the art can make modifications and equivalents to the embodiments of the present invention without departing from the spirit and scope of the present invention, which is set forth in the claims of the present application.

Claims (5)

1. A numerical simulation method applied to electromagnetic transient analysis of a switching circuit is characterized by comprising the following steps:

writing a system state equation of a circuit in a row, and converting the system state equation into a state equation in a first-order form through variable substitution;

step (2) modifying a system state equation of the circuit at the moment of switching action;

step (3) iterative solution is carried out on the first-order state equation by using a direct integral method;

step (4) reinitializing the potential of the circuit node at the moment of switching action;

step (5) eliminating numerical oscillation by using minimum step length backward Euler operation at the switching time of the switch state;

step (6) carrying out equivalence certification at the moment of switching action;

in the step (1), a system state equation expression of the circuit is as follows:

wherein:

is a column vector of order N, each element value representing the potential of the corresponding node; k_CIs a capacitance coefficient matrix; k_RIs a resistivity matrix; k_LIs an inductance matrix; i is an input vector in the circuit and represents a power supply in the circuit; the three coefficient matrixes are square matrixes of N rows and N columns, and all the three coefficient matrixes meet the requirements of positive nature, symmetry and sparsity;

system state equation 1 is transformed into a first order state equation by variable substitution as follows:

wherein:

K₁、K₂all represent a coefficient matrix; x is a state variable; e is the unity diagonal matrix;

the step (6) comprises the following steps:

<1> proof of equivalence of inductance-related terms:

for the inductive branch, the existence equation (15) holds:

I_l,d-＝L^-1[Ψ_d-(i)-Ψ_d-(j)] (15)

or written in a matrix format as shown in equation (16):

continuing to expand to the following equation (17):

since equation (17) holds for each inductive branch, the inductance-related term in equation (9)

Term related to inductance in equation (14)

Are equivalent; i is_l,d-Represents t_d-The current value on the inductance branch at the moment; l is^-1To represent

L is the inductance value of the inductance branch; Ψ_d-(i) And Ψ_d-(j) Represents t_d-Voltage time integral quantity of the node i and the node j at the moment;

<2> proof of equivalence of capacitance-related terms:

considering the norton equivalent current source of the capacitive branch, the capacitance-related term in equation (9) is expanded into equation (18) and equation (19):

the capacitance in equation (14) is related to

The developed format is written according to each capacitive branch as shown in equations (20) and (21):

the capacitance in equation (14) is related to

The developed format is written according to each capacitive branch as shown in equations (22) and (23):

due to the fact that

Represents t_d-The voltage value on the capacitor at the moment, while taking into account R_SAnd Δ t₀The fact that all tend to zero, the capacitance-related term in equation (9)

Term related to capacitance in equation (14)

Equivalence; u shape_c,d-Represents t_d-At the moment, the voltage value on the capacitor branch circuit;

and

respectively representing potential values of the node i and the node j;

the switching operation uses a minimum step length backward Euler operation at the moment, and is equivalent to the operation of reinitializing the node potential according to the capacitor voltage and the inductance current in the circuit.

2. The numerical simulation method according to claim 1, wherein in the step (2), if an ideal switch connected between two nodes is operated at a switching operation time, t is at a switching operation time_dModifying the resistivity matrix K in equation (1)_R(ii) a The system state equation for the circuit is modified using equation (4), which is as follows:

where Turn On indicates that the ideal switch is at t_dOpening at any moment; turn Off indicates that the ideal switch is at t_dDisconnection at any moment; except for K_r(i,i)＝R_s ^-1，K_r(i,j)＝-R_s ^-1，K_r(j,i)＝-R_s ^-1And K_r(j,j)＝R_s ^-1Other than, K_rThe other elements in the voltage source are zero, and the internal resistance R of the 0V ideal voltage source_sTends to zero; k_R,+Is a matrix of resistivity before switching action, K_R,-Is a resistance coefficient matrix after switching action; k_rIs an N-th order matrix, K_r(i, i) represents the element in the ith row and ith column in the matrix; k_r(i, j) represents the element in the ith row and jth column of the matrix; k_r(j, i) represents the element in the jth row and ith column in the matrix; k_r(j, j) represents the element in the jth row and jth column in the matrix;

R_sis the internal resistance.

3. The numerical simulation method according to claim 1, wherein in the step (3), the iterative solution is performed by using a direct integration method for the first-order equation 2 in the standard form, and an iteration format of the direct integration method is as follows:

wherein: Δ t is the time step between the nth time instant and the time instant n + 1; beta is a coefficient, when the value of beta is 0, the method is called a forward difference method, the stability is stable, and the numerical precision is first-order precision; when the value of beta is 0.5, the method is called a Crank-Nicolson method, the stability is unconditionally stable, and the numerical precision is second-order precision; when the value of beta is 1, the method is called a posterior difference method, the stability is unconditionally stable, and the numerical precision is first-order precision; x is the number of_nRepresents the system state variable at the nth step, namely, the time t is n delta t; x is the number of_n+1Represents the system state variable at the n +1 th step, i.e. at the time t ═ n +1) Δ t; r_nAnd R_n+1The system input vectors of the nth step, i.e., time t ═ n Δ t and the nth +1 th step, i.e., time t ═ n +1) Δ t, are represented, respectively.

4. The numerical simulation method according to claim 2, wherein the step (4) includes:

at the time t of switching action_dThe capacitor is equivalent to a voltage source, the inductor is equivalent to a current source, and the current source is based on t_d-Capacitance voltage and inductance current calculation t at time_d+Potential distribution at a time; t is t_d-Represents t_dA certain time before the time, and t_d-t_d-Tends to zero; t is t_d+Represents t_dA certain time after the time, and t_d+-t_dTends to zero;

if the capacitance and inductance dynamic elements in the circuit are treated equivalently, only a resistance branch and a current source branch exist in the circuit, a static circuit equation is obtained, and the equation in the formula (6) is solved

Completing the initialization of the circuit after the switching action;

wherein:

represents t_d+A node potential at a time; k is a coefficient matrix; i is an input vector in the circuit and represents a power supply in the circuit; the coefficient matrix K contains information K of resistance branches and voltage source resistance branches_R,d+Internal resistance information K of equivalent voltage source of sum capacitor_RC；

For a certain branch circuit which is connected between the node I and the node j and has the capacitance C, the coefficient matrix K_rcIncluding its equivalent to an ideal voltage source, with respect to its internal resistance R_sBranch information that tends to zero; matrix K_rcMiddle removing K_rc(i,i)＝R_s ^-1，K_rc(i,j)＝-R_s ^-1，K_rc(j,i)＝-R_s ^-1，K_rc(j,j)＝R_s ^-1Except that, other elements are zero; k_rcIs an n-th order matrix, K_rc(i, i) represents a matrix K_rcRow i and column i; k_rc(i, j) represents a matrix K_rcRow i and column j; k_rc(j, i) represents a matrix K_rcThe element in the jth row and ith column; k_rc(j, j) represents a matrix K_rcRow j and column j;

the vector I includes all current sources and the ideal of the original circuitNorton equivalent current source information I of voltage source_S,d+Besides, the current source also comprises information I of an inductance equivalent current source_L,d-And information I of the equivalent ideal voltage source of the capacitor_SC,d-As shown in formula (8):

I＝I_S,d++I_L,d-+I_SC,d- (8)

wherein:

for a certain branch circuit which is connected between the node I and the node j and has inductance L, I_lIncluding information of its equivalent current source; except that I_l,d-(i)＝I_L,d-And I_l,d-(j)＝-I_L,d-Outer, I_l,d-All other elements in (A) are zero; meanwhile, for a certain branch which is connected between the node I and the node j and has the capacitance of C, I_scIncluding its equivalent norton current source information; i is_sc,d-In, except for I_sc,d-(i)＝I_SC，d-And I_sc,d-(j)＝-I_SC，d-In addition, the other elements are all zero, wherein I_SC，d-Represents t_d-The current value of the time capacitor equivalent Norton current source;

substituting formula (7) and formula (8) into formula (6), and arranging into a format of the following formula:

wherein:

represents t_d+The potential of each node at a time; k_R,d+Represents t_d+A resistance coefficient matrix which constantly contains information of the resistance branch and the voltage source internal resistance branch; i is_S,d+Represents t_d+Current source information vectors at time; i is_SC,d-Represents t_d-Time capacitance equivalent Norton current source information vector; l _ relative represents the inductance-related term; c _ relative represents a capacitive element-related term.

5. The numerical simulation method according to claim 2, wherein the step (5) includes:

at the time t of switching action_dUsing a backward euler operation:

if at the moment t of switching action_dAt a very small time step Δ t₀Applying a backward Euler operation to System State equation 1, where Δ t₀＝t_d+-t_d-Using t_d-State quantity x of time_d-To find t_d+State quantity x of time_d+The calculation expression is as follows:

wherein the content of the first and second substances,

Ψ is a time integral quantity of voltage at each node, defined as

And

R_d+an input vector matrix representing the time after the switching action; k_2,d+Is a coefficient matrix;

further expanding formula (10) into formula (11) and formula (12):

taking into account Δ t₀The value tends to 0, and the expression (12) is changed to the expression (13), which shows that the integral value of the node potential with time does not change before and after the switching action;

Ψ_d+＝Ψ_d- (13)

wherein: Ψ_d+And Ψ_d-Respectively representing the voltage time integral quantity of each node after the switching action and the voltage time integral quantity of each node before the switching action;

substituting formula (13) for formula (11) and arranging into the form of formula (14):

if equations (14) and (9) are equivalent, it is stated that: at the time t of switching action_dAnd one node potential reinitialization operation based on the classical theory is completed by one small-step backward Euler operation.

Images