Electromagnetic transient simulation method for field programmable logic array

Abstract

Electromagnetic transient simulation method applicable to a field programmable gate array (FPGA), which integrates topological parameters of a circuit to be simulated into two matrix parameters in an initialization stage thereof; and voltage and current information at each simulation moment can be obtained only through simple matrix multiplication operation in a main part of the simulation cycle thereof. The method avoids complex initialization operation in the field programmable logic array; meanwhile, the flow of the main part of the simulation cycle in the FPGA is maximally compressed, and the efficiency of electromagnetic transient simulation based on the FPGA is greatly improved.

Claims

1. A method for real-time electromagnetic transient simulation based on a field programmable logic array (FPGA) for power system stability analysis and control research, comprising providing a host computer and a field programmable logic array (FPGA) for real-time electromagnetic transient simulation of a circuit of a power system, performing an initialization operation in the host computer, wherein a history current source vector I.sub.his.sup.1, an equivalent current source vector I.sub.src.sup.n, voltage coefficient matrix ?, and current coefficient matrix ? of a history current source expression are obtained through the following steps of (1) to (5): (1) sequentially numbering branches and nodes in the circuit to be simulated respectively, wherein the number of a grounding node is 0; (2) forming a correlation matrix M of the circuit to be simulated according to the following rules: (i) if a branch p is connected to a node q and a positive current direction defined by the branch p is an outflow node q, M (q, p)=1; (ii) if the branch p is connected to the node q and the positive current direction defined by the branch p is an inflow node q, M (q, p)=?1; (iii) if the branch p and the node q are not connected, M (q, p)=0; (3) forming a branch equivalent admittance vector Y.sub.eq, a node admittance matrix Y.sub.n, a voltage coefficient matrix ?, and a current coefficient matrix ? of a history current source expression of the circuit to be simulated according to the following substeps: (i) respectively replacing each resistance branch, inductance branch, capacitance branch and switch branch with an companion circuit model, wherein each companion circuit respectively comprises an equivalent admittance, and a history current source connected in parallel; (ii) an independent voltage source branch and an independent current source branch are represented by a Norton equivalent circuit, and each Norton equivalent circuit comprises an equivalent admittance, and an equivalent current source connected in parallel; (iii) forming branch equivalent admittance column vectors Y.sub.eq according to the branch numbers by equivalent admittances of all branches, forming branch history current source column vectors I.sub.his according to the branch numbers by history current sources of all branches, and forming branch equivalent current source column vectors I.sub.src according to the branch numbers by equivalent current sources of all branches, wherein elements at corresponding positions in I.sub.src of resistance, inductance, capacitance and switch branches are zero; and the elements at the corresponding positions in I.sub.his of the independent voltage source and the independent current source branch are zero; (iv) calculating a node admittance matrix Y.sub.n of the circuit to be simulated according to the equivalent admittance of each branch; (v) calculating the voltage coefficient matrix ? and the current coefficient matrix ? of the history current source expression, wherein the formula is as follows:
I.sub.his.sup.n+1=?Y.sub.eqV.sub.brn.sup.n+?I.sub.brn.sup.n wherein, I.sub.his.sup.n+1 is a history current source vector at a (n+1).sup.th simulation moment, a branch voltage vector at a n.sup.th simulation moment and a branch current vector at the n.sup.th simulation moment; (4) forming a node voltage/branch current coefficient matrix P and a history current source coefficient matrix Q according to the correlation matrix M, the branch equivalent admittance vector Y.sub.eq and the node admittance matrix Y.sub.n of the circuit to be simulated; $P=\left[\begin{array}{c}-Y_n^{-1}M\\ -Y_{\mathrm{eq}}{M}^T{Y}_n^{-1}M+I\end{array}\right]$
Q=Y.sub.eqM.sup.TY.sub.n.sup.?1M wherein, Y.sub.n.sup.?1 is an inverse matrix of the node admittance matrix, M.sup.T is a transpose matrix of the correlation matrix, I is an identity matrix with dimensions of N.sub.brn*N.sub.brn, and N.sub.brn is the number of branches of the circuit to be simulated; and (5) setting the initial history current source vector I.sub.his.sup.1 to zero and the current simulation moment n to 1 at end of the initiation operation; performing one or more cycles of simulation in the FPGA, wherein electric parameters at each simulation moment in each simulation cycle are obtained by a compressed matrix calculation in the FPGA through the following steps of (6) and (7): (6) calculating a node voltage vector V.sub.n.sup.n and a branch current vector I.sub.brn.sup.n at the current simulation moment according to the history current source vector I.sub.his.sup.n and the equivalent current source vector I.sub.src.sup.n at the current simulation moment, and updating the history current source vector I.sub.his.sup.n+1 at a next simulation moment, wherein the formula is as follows: $\left[\begin{array}{c}{V}_n^n\\ I_{\mathrm{brn}}^n\end{array}\right]=P\left(I_{\mathrm{his}}^n+I_{\mathrm{src}}^n\right)$ $I_{\mathrm{his}}^{n+1}=-\left(?+?\right)Q\left(I_{\mathrm{his}}^n+I_{\mathrm{src}}^n\right)+?\left(I_{\mathrm{his}}^n+I_{\mathrm{src}}^n\right)$ wherein, the equivalent current source vector I.sub.src.sup.n is automatically updated along with the n magnitudes of independent voltage sources and independent current sources, diagonal elements corresponding to switch branches in diagonal arrays ? and ? change instantly along with switch states, and other elements remain constant; and (7) returning to the step (6) when the next simulation moment n+1 becomes the current simulation moment, and running the steps (6) and (7) repeatedly to update the electric parameters real-time in the FPGA; displaying real-time simulation waveforms of the electric parameters obtained by the FPGA at each simulation moment in the host computer, and conducting stability analysis and control research of the power system based on the displayed real-time simulation waveforms; and terminating the simulation until a predetermined simulation moment is reached or an instruction of early termination is received.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) FIGS. 1A and 1B show the comparison of the waveform of transformer voltage and current under the method and off-line simulation, where FIG. 1A shows the off-line simulation result of a power electronic simulation software PLECS, and FIG. 1B shows the real-time simulation result adopting the method of the present invention.

(2) FIG. 2 shows the circuit diagram of the dual active bridge (DAB) in the present invention.

(3) FIG. 3 is a flow diagram showing the electromagnetic transient simulation method for the field programmable logic of the present invention.

DETAILED DESCRIPTIONS OF THE PRESENT INVENTION

(4) The present invention is described with reference to a dual active bridge (DAB) circuit as shown in FIG. 2, while it should not be construed as limiting the scope of the present invention.

(5) When the method of the present invention is used for specifically realizing the real-time simulation of DAB, the hardware is mainly a PXIe case of American National Instruments (NI), where PXIe-8135 is a PXIe controller mainly responsible for the simulation of a DAB control system, and is communicated with an upper computer via Ethernet to display a real-time simulation waveform on the host computer; and PXIe-7975R is an FPGA module mainly responsible for the DAB circuit simulation, and is connected with an external controller and an oscilloscope via an I/O port so as to preform hardware-in-the-loop simulation. Real-time communication is carried out between the two by a PXIe bus to complete real-time simulation.

(6) The host computer program, the PXIe controller program and the FPGA module program for carrying out real-time simulation through the method of the present invention are programmed and realized by the Labview development environment of American National Instruments (NI). Using the Labview development environment, the host computer can communicate with the PXIe controller, and display simulation waveforms, etc.; and the PXIe controller can communicate with the host computer, read and write data from the FPGA module, simulate the control system of DAB, etc. The above programs do not fall within the scope of the present invention, and there is a related program example on the official website of American National Instruments (NI), so the greater details are not provided here. In the present invention, it is programmed by Labview, and specifically implemented by the FPGA module, and FIG. 3 is a flow chart showing the electromagnetic transient simulation method for the field programmable logic array of the present invention.

(7) In one embodiment of the electromagnetic transient simulation method for the field programmable logic array of the present invention, the method comprises the steps of: (1) sequentially numbering branches and nodes in a circuit to be simulated respectively, wherein the number of a grounding node is 0, as shown in FIG. 2; (2) forming a correlation matrix M of a circuit to be simulated;

(8) $M=\left[\begin{array}{cccccccccccccc}1& 1& 1& 0& 1& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& -1& 1& 0& 0& 1& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& -1& 1& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& -1& 1& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& -1& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& -1& 1& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& -1& 1& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 1& 0& 1& 0& 1& 1\end{array}\right]$ (3) forming a branch equivalent admittance vector Y.sub.eq, a node admittance matrix Y.sub.n, a voltage coefficient matrix ? and a current coefficient matrix ? of a history current source expression of the circuit to be simulated according to the following substeps: (i) respectively replacing each resistance branch, inductance branch, capacitance branch and switch branch with an companion circuit model, wherein each companion circuit respectively comprises an equivalent admittance, and a history current source connected in parallel, and the equivalent model of the branch is shown in FIGS. 1A and 1B; (ii) an independent voltage source branch and an independent current source branch are represented by a Norton equivalent circuit, and each Norton equivalent circuit comprises an equivalent admittance, and an equivalent current source connected in parallel; (iii) forming branch equivalent admittance column vectors Y.sub.eq according to the branch numbers by equivalent admittances of all branches, forming branch history current source column vectors I.sub.his according to the branch numbers by history current sources of all branches, and forming branch equivalent current source column vectors I.sub.src according to the branch numbers by equivalent current sources of all branches, wherein elements at corresponding positions in I.sub.src of resistance, inductance, capacitance and switch branches are zero; and the elements at the corresponding positions in I.sub.his of the independent voltage source and the independent current source branch are zero (see step 6 for details of specific numerical values of I.sub.his and I.sub.src in calculation);
Y.sub.eq=[100 1600 2 2 2 2 0 3.3 2 2 2 2 400 0].sup.T (iv) calculating a node admittance matrix Y.sub.n of the circuit to be simulated according to the equivalent admittance of each branch (specific methods can refer to electromagnetic transient simulation textbooks);

(9) $Y_n=\left[\begin{array}{cccccccc}1.7040& -0.0020& -0.0020& 0& 0& 0& 0& 0\\ -0.0020& 0.0040& 0& 0& 0& 0& 0& 0\\ -0.0020& 0& 0.0040& 0& 0& 0& 0& 0\\ 0& 0& 0& 0.0033& -0.0033& 0& 0& 0\\ 0& 0& 0& -0.0033& 0.0033& 0& 0& 0\\ 0& 0& 0& 0& 0& 0.0040& 0& -0.0020\\ 0& 0& 0& 0& 0& 0& 0.0040& -0.0020\\ 0& 0& 0& 0& 0& -0.0020& -0.0020& 0.4040\end{array}\right]$ (v) according to a history current source expression of each branch referred to in 3.1), expressing a relationship between the history current source vector I.sub.his.sup.n+1 at a (n+1).sup.th simulation moment and the branch voltage vector V.sub.brn.sup.n and the branch current vector I.sub.brn.sup.n at a n.sup.th simulation moment into the following forms, and obtaining a voltage coefficient matrix ? and a current coefficient matrix ? of the history current source expression (detailed calculation values are shown in the step (6):
I.sub.his.sup.n+1=?Y.sub.eqV.sub.brn.sup.n+?I.sub.brn.sup.n (4) forming a node voltage/branch current coefficient matrix P and a history current source coefficient matrix Q which are directly used in the main part of the simulation cycle according to a correlation matrix M, a branch equivalent admittance vector Y.sub.eq and a node admittance matrix Y.sub.n of the circuit to be simulated;

(10) $P=\left[\begin{array}{c}-Y_n^{-1}M\\ -Y_{\mathrm{eq}}{M}^T{Y}_n^{-1}M+I\end{array}\right]$$Q=Y_{\mathrm{eq}}{M}^T{Y}_n^{-1}M$$P={\left[\begin{array}{ccccc}-5.8754?{10}^{-4}& -0.9401& \mathrm{.Math.}& 0& 0\\ -5.8754?{10}^{-4}& -0.4700& \mathrm{.Math.}& 0& 0\\ \mathrm{.Math.}& \mathrm{.Math.}& ?& \mathrm{.Math.}& \mathrm{.Math.}\\ 0& 0& \mathrm{.Math.}& 2.0296& -0.9949\\ 0& 0& \mathrm{.Math.}& -0.0398& 0.9999\end{array}\right]}_{22?22}$$Q={\left[\begin{array}{ccccc}0.0588& 0.0588& \mathrm{.Math.}& 0& 0\\ 0.9401& 0.9401& \mathrm{.Math.}& 0& 0\\ \mathrm{.Math.}& \mathrm{.Math.}& ?& \mathrm{.Math.}& \mathrm{.Math.}\\ 0& 0& \mathrm{.Math.}& 0.9949& 0.9949\\ 0& 0& \mathrm{.Math.}& 9.9493?{10}^{-5}& 9.9493?{10}^{-5}\end{array}\right]}_{22?22}$ wherein, Y.sub.n.sup.?1 is an inverse matrix of the node admittance matrix, M.sup.T is a transpose matrix of the correlation matrix, I is an identity matrix with dimensions of N.sub.brn*N.sub.brn, and N.sub.brn is the number of branches of the circuit to be simulated; (5) setting the initial history current source vector I.sub.his.sup.1 to zero, wherein the current simulation moment n is set to 1, and the initialization stage is completed; and the above steps are completed in an host computer; (6) completing the following functions mainly by the main part of the simulation cycle, and algorithm compiling in the field programmable logic array (FPGA) (taking n=10 as an example): (i) calculating a node voltage vector V.sub.n.sup.n and a branch current vector I.sub.brn.sup.n at the current simulation moment according to a history current source vector I.sub.his.sup.n and an equivalent current source vector I.sub.src.sup.n at the current simulation moment;

(11) $\left[\begin{array}{c}{V}_n^n\\ I_{\mathrm{brn}}^n\end{array}\right]=P\left(I_{\mathrm{his}}^n+I_{\mathrm{src}}^n\right)$
I.sub.his.sup.10=[0 ?180.9965 90.8266 ?304.9359 ?304.9359 90.8266 3.1455 0 ?6.6805 ?0.5579 ?0.5579 ?6.6805 ?0.0511 0].sup.T
I.sub.src.sup.10=[?40000 0 0 0 0 0 0 0 0 0 0 0 0 0].sup.T
V.sub.n.sup.10=[180.9965 182.1912 ?1.1947 1.0197 0.0761 0.3433 ?0.2921 0.0511].sup.T
I.sub.brn.sup.10=[?21900.3489 21717.0625 93.2160 90.0705 90.0705 93.2160 3.1455 3.1455 ?6.0962 0.1948 0.1948 ?6.0962 5.8993 0.0020].sup.T (ii) meanwhile, updating the history current source vector I.sub.his.sup.n+1 at a next simulation moment according to the history current source vector I.sub.his.sup.n and the equivalent current source vector I.sub.src.sup.n at the current simulation moment;
I.sub.his.sup.n+1=?(?+?)Q(I.sub.his.sup.n+I.sub.src.sup.n)+?(I.sub.his.sup.n+I.sub.src.sup.n)
?=[0 ?1 1 ?1 ?1 1 0 0 1 ?1 ?1 1 ?1 0], ?=diag(?).sub.14?14
?=[0 0 1 0.66 0.66 1 1 0 1 0.66 0.66 1 0 0], ?=diag(?).sub.14?14
I.sub.his.sup.11=[0 ?193.7770 86.3684 ?333.7840 ?333.7840 86.3684 3.7502 0 ?7.9206 ?0.6259 ?0.6259 ?7.9206 ?0.0689 0].sup.T wherein, the equivalent current source vector I.sub.src.sup.n is automatically updated along with the magnitudes of independent voltage sources and independent current sources, diagonal elements corresponding to switch branches in diagonal arrays ? and ? change along with switch states, and other elements remain constant; and (iii) repeatedly running the step (6) until a specified simulation moment or when an instruction of early termination is received, and finishing the simulation.
I.sub.his.sup.11=[0 ?193.7770 86.3684 ?333.7840 ?333.7840 86.3684 3.7502 0 ?7.9206 ?0.6259 ?0.6259 ?7.9206 ?0.0689 0].sup.T
I.sub.src.sup.11=[?40000 0 0 0 0 0 0 0 0 0 0 0 0 0].sup.T
V.sub.n.sup.11=[193.7770 194.8916 ?1.1145 1.3884 0.2633 0.3789 ?0.3100 0.0689].sup.T
I.sub.brn.sup.11=[?20622.2958 20448.8510 88.5975 84.8473 84.8473 88.5975 3.7502 3.7502 ?7.3006 0.1999 0.1999 ?7.3006 7.0979 0.0028].sup.T
?=[0 ?1 1 ?1 ?1 1 0 0 1 ?1 ?1 1 ?1 0], ?=diag(?).sub.14?14
?=[0 0 1 0.66 0.66 1 1 0 1 0.66 0.66 1 0 0], ?=diag(?).sub.14?14
I.sub.his.sup.12=[0 ?205.8112 82.2295 ?360.9527 ?360.9527 82.2295 4.3911 0 ?9.2285 ?0.6965 ?0.6965 ?9.2285 ?0.0898 0].sup.T

(12) Table 1 shows DAB circuit parameters for real-time simulation under the method of the present invention. Table 2 shows the utilization of hardware resources in real-time simulation of DAB under the method of the present invention. Table 3 shows the simulation step sizes and running time in one step size per cycle for real-time simulation of DAB under the method of the present invention.

(13) TABLE-US-00001 TABLE 1 DAB Circuit Parameters Parameters Symbols Values Input voltage V.sub.in 400 V Capacitance at a high C.sub.1 800 ?F voltage side Transformer leakage L 160 ?H inductance Conductive loop resistance R.sub.0 0.3 Ohm Transformer transformation N 2:1 ratio Capacitance at a low C.sub.2 200 ?F voltage side Resistive load R.sub.L 25 Ohm

(14) TABLE-US-00002 TABLE 2 Hardware resource utilization Hardware resources Utilization rate Total number of logical slices 36.9% (23419 of 63550) Logic slice register 15.3% (77752 of 508400) Logic slice LUT 24.8% (63107 of 254200) RAM block 11.4% (91 of 795) DSP48 13.3% (205 of 1540)

(15) TABLE-US-00003 TABLE 3 Real-time Performance Types of time Length of time Step size (ns) 500 Runtime (ns) in each step size 218.75