Interpolatory proper order decomposition of nonlinear transmission line circuits

Abstract

Purpose

The paper is concerned with interpolatory proper orthogonal decomposition (IPOD) methods for nonlinear transmission line circuits. This paper aims to examine several factors that must be considered when applying such model reduction techniques to this kind of circuit.

Design/methodology/approach

Two types of POD will be implemented. In each case, the choice of the order of the reduced model and the order of the interpolation space shall be considered. The stability of the models shall be explored.

Findings

The results indicate that the order for the reduced model to obtain accurate results depends on the chosen method when considering nonlinear transmission lines. The results also indicate that the structure of the nonlinear transmission line is crucial for determining the stability of the reduced models.

Originality/value

The work compares two IPOD methods and discusses the issues involved in achieving an accurate and stable reduced-order model for a nonlinear transmission line.

Citation

Condon, M. and Hayes, B. (2023), "Interpolatory proper order decomposition of nonlinear transmission line circuits", COMPEL - The international journal for computation and mathematics in electrical and electronic engineering, Vol. 42 No. 2, pp. 605-619. https://doi.org/10.1108/COMPEL-07-2022-0250

Introduction

With the advances in modern technologies, models of high-frequency circuits and systems are becoming ever more complex. The focus in this paper is on nonlinear transmission line circuits. Nonlinear transmission lines find application for high-power pulse generation (Bragg et al., 2013), edge sharpening (Ricketts et al., 2019; Bobreshov et al., 2022) and pulse shaping (Gardner et al., 2022). However, the models for these involve a large number of sections and when detailed models of the components are included, this can lead to large systems of differential equations. Nonlinear transmission lines are governed by the Korteweg-de-Vries equation when it is assumed that losses are negligible, the number of sections is very large and the nonlinear elements are modelled using a specific type of function (Giambo et al., 1984; Kuek et al., 2012; Nikoo and Hashemi, 2017). A circuit analysis approach is an alternative modelling approach and is not restricted by these constraints. Circuit design involves multiple simulations for optimization and analysis purposes. Hence, fast and reliable generation of results from models are essential. However, the larger and more nonlinear the systems become, the more computationally expensive the simulations become. To this end, model reduction strategies are considered. Their purpose is to identify the most significant behaviour of the system and thereby reduce the computational cost.

A vast array of linear model reduction techniques exist (Antoulas et al., 2001). Various nonlinear model reduction techniques exist such as empirical balanced truncation (Lall et al., 1999), piecewise linear approximation (Bond and Daniel, 2009), piecewise polynomial approximation (Dong and Roychowdhury, 2003, 2008; Qiu and Jiang, 2020) and many more. However, one important issue concerning forming a reduced-order model is preserving the stability of the original model. In Bond and Daniel (2009), stability is addressed and systems for which global and local stability of the models can be preserved are identified. In Qiu and Jiang (2020), conditions for which stability is preserved with piecewise polynomial model reduction are given.

Proper orthogonal decomposition (POD) is reliable and used for strongly nonlinear circuits (Nouri et al., 2017; Nouri and Nakhla, 2018; Prajna, 2003). However, the basic POD does not result in complexity reduction (Nouri and Nakhla, 2018) and because of this interpolatory POD (IPOD) has been proposed as in Chaturantabut and Sorensen (2010) and gappy POD (Everson and Sirovich, 1995; Willcox, 2006). In IPOD, interpolation is used to approximate the nonlinear functions. Gappy POD is used when experimental data is missing or corrupted and involves estimating the missing data. Peherstorfer et al. (2020) examine the stability of discrete empirical interpolation and gappy POD with randomized and deterministic sampling points. In their work, they show that certain sampling approaches lead to more stable reduced models.

This paper shall examine IPOD model reduction techniques for simulation for nonlinear transmission line circuit models. The paper shall review two implementations of this technique, that proposed by Nouri and Nakhla (2018) and that in Chaturantabut (2020). These shall be referred to as Method 1 and Method 2 in what follows. Several studies shall be performed and suggestions shall be made based on the findings and observations.

Interpolatory proper orthogonal decomposition (Nouri and Nakhla, 2018)

Consider the nonlinear system:

(1) dxdt=F(t,x)x(0)=x0

x = x(t) is an n-dimensional state vector and F: [0, ∞) × X →ℝⁿ is a nonlinear differentiable function, X ⊆ ℝⁿ.

An initial transient simulation of the system in (1) is performed. The training input should be wideband and excite as many as possible of the nonlinear modes of the given system (Nouri and Nakhla, 2018).

Let χ be a set of samples of the states at times t₁ to t_N:

(2) χ=[x(t1), x(t2),…,x(tN)] ∈ ℝn×N

To form the projection matrix, singular value decomposition (SVD) is performed.

(3) χ=VχΣχWχtVχtVχ=In×nVχ={v1,..,vn}∈ℝn×n

V_χ is a unitary matrix of dimension n × n, W_χ is a unitary matrix of dimension N × N and Σ_χ is a matrix with elements along its diagonal that are the singular values of χ.

The projection matrix is formed from the first m left singular vectors:

V={v1,..,vm}∈ℝn×m    m⪡n

The reduced state space is formed as:

(4) x(t)=Vx^(t)x^(t)∈ℝm

The resultant reduced model is:

(5) d d t x ^ ( t ) = V t F ( V x ^ ( t ) )  

where F is as defined in equation (1).

While the dimension of the resultant model is reduced from n to m, the complexity is not reduced as the full x(t) has to be formed to evaluate F.

With a view to reducing the complexity, IPOD is used.

Let:

(6) F=[F(x(t1)), F(x(t2)),…,F(x(tN))] ∈ℝn×N

where F is the vector field in (1) and x(t_i) are the snapshots in (2).

An SVD of F is performed:

(7) F=UFΣFWFt

The k left singular vectors corresponding to the k largest singular values are selected. The choice of k that is recommended from Nouri and Nakhla (2018) is:

(8) Rk≜1−∑r=1kσr[∑l=1mσl]−1≤α,    k≤m

σ_r/l are the diagonal elements of Σ_F, and α is a small threshold value set to achieve a balance between accuracy and efficiency.

This measure is based on the relative energy associated with the k selected basis functions. The relative energy is defined as ∑ r = 1 k σ r [ ∑ l = 1 n σ l ] − 1 . However, to determine this quantity requires a full SVD so that the n singular values are determined. Hence, Nouri and Nakhla (2018) recommend the formula in (8) where m singular values are required.

Let:

(9) U={u1,..,uk} ∈ℝn×k,  ui∈UF

An approximation for the nonlinear function is made:

(10) F(x(t))=Uc(t)

where c(t) is a coefficient vector.

Next a selector matrix P is formed that selects the k elements of F that most contribute to the space spanned by U. The algorithm for selecting such rows is given in Nouri and Nakhla (2018) and Chaturantabut (2020):

(11) F r ( x r ( t ) ) = P t F ( x ( t ) )

Not all states may contribute to F_r so:

(12) x r ( t ) = Υ t x ( t )

y is a selector matrix with elements that are either 0 or 1. It selects from x(t) those elements that contribute to F_r.

So from (10), (11) and (12):

(13) F r ( Υ t x ( t ) ) = P t U c ( t )

The approximation for the nonlinear function is thus:

(14) F ˜ ( x ( t ) ) = U c ( t ) = U ( P t U ) − 1 F r ( Υ t x ( t ) )

Now from the initial POD:

x(t)=Vx^(t)

and

F ^ ( x ^ ) ≈ V t F ˜ ( V x ^ )

So

(15) F^(x^)≈ΨrFr(Vrx^)

V r = Υ t V

x r ( t ) = V r x ^ ( t )

Ψ r = V t U ( P t U ) − 1

F r ( x r ( t ) ) = P t F ( Υ t x ( t ) )

The Jacobian is:

(16) Jr(xr)=∂Fr(xr)∂xr

Equations (15) and (16) are used in the numerical integration routine that provides the transient solution of the system and the complexity cost has been reduced. F^ and J_r avoid computation of the full n-dimensional system. However, stability of this method needs to be addressed.

Stability of reduced-order models and stabilized proper orthogonal decomposition.

Logarithmic norm (Dahlquist,1959)

Let A be a constant n × n real matrix. The logarithmic norm is defined as:

μ[A]=lim⁡h→0+||I+hA||−1h

where ‖ ‖ is the standard Euclidean norm.

It can also be determined as (Söderlind, 2006):

(17) μ[A]=λmax(A+At2)

Logarithmic Lipschitz constant (Söderlind, 2006)

Let F be a function F:Y→X and let u,v ϵ ℝⁿ

(18) M [ F ] = sup ⁡ u ≠ v ( u − v ) t ( F ( u ) − F ( v ) ) | | u − v | | 2 u , v ∈ ℝ n

When F=A∈ℝn×n, μ[A]=M[F].

Definition (infinitesimally contracting)

(Aminzare and Sontag, 2014)

A time-dependent vector field F : [ 0 , ∞ ) × X → ℝ n ,     X ⊆ ℝ n is infinitesimally contracting if:

(19) μ[  JF(t,x)]≤−c,  ∀x∈X,  ∀t≥0

where c > 0 is the contraction rate. J_F (t,x) is the Jacobian of F (x,t).

Chaturantabut (2020) gives a stronger condition:

The function F is infinitesimally contracting if sup⁡tϵ[0,∞)M[F(t,x)]<0.

Lemma 1 (Chaturantabut,2020)

Suppose the nonlinear vector field F of the full-order system is infinitesimally contracting.

Consider the reduced-order model:

(20) ddtx^=F^(x^)

where

F ^ ( x ^ ) = V t W F ( W t V x ^ )

and x ^ ∈ ℝ m ,   V ∈ ℝ n × m as in equation (4).

V has m orthonormal columns and W is a matrix such that W t V ∈ ℝ n × m has full column rank.

Then the nonlinear vector field F^(x^) is also infinitesimally contracting.

The proof of Lemma 1 is given in Chaturantabut (2020).

The following proposition, Lemma 2, is proven in Chaturantabut (2020).

Lemma 2

Suppose the nonlinear vector field F in the full model is infinitesimally contracting, then the reduced system in the form given in (20) preserves exponential stability.

Consequently, POD preserves stability for contracting systems as the reduced model is of the form in (20) with W = I.

However, the given IPOD model reduction method is not in the form of (20) and so a condition for ensuring its stability is necessary.

Lemma 3

Let F be the nonlinear vector field of the full-order system (1) and let F be infinitesimally contracting. Let the reduced system be given by:

(21) ddtx^=F^(x^)

where:

F ^ ( x ^ ) = V t U ( P t U ) − 1 P t F ( V x ^ )

and again as in equation (4), x^∈ℝm V∈ℝn×m.

Then the vector field F^ is infinitesimally contracting if:

M[℘F]<0

where:

℘ = U ( P t U ) − 1 P t

It would be advantageous to have an IPOD arrangement in the form of (20) as if the original system is exponentially stable, then the reduced system would preserve this property.

One scheme is proposed in Chaturantabut (2020).

Set W = HP^t. The goal is to choose H such that:

(22) min ⁡ H ∈ ℝ n × k | | P t F ( V x ^ ) − P t H P t F ( P H t V x ^ ) | |

One choice is W = PP^t where P is the n × k selector matrix in (11):

(23) ddtx^=F^(x^)

where:

F ^ ( x ^ ) = V t P P t F ( P P t V x ^ )

x^∈ℝm, V∈ℝn×m

r a n k ( P P t V ) = m

and

k≥m

The Jacobian is:

(24) JF=dF^dx^

Then the reduced system in (23) preserves exponential stability. Note that for this scheme to hold true k ≥ m where k is the order of the interpolation space and m is the order of the reduced model. However, as noted earlier, when selecting the order k, it was selected as ≤ m to avoid costly SVD of the F matrix in (7). Hence, to meet both stability and efficiency criteria, k = m is a suitable choice. In addition, while PP^t has a dimension of n × n, the complexity of forming PP^t V is O(km) because P is a selector matrix with only k nonzero entries.

However, the desire to have:

(25) F ( V x ^ ) = W F ( W t V x ^ )

is not fulfilled when W = PP^t. This leads to a reduction in accuracy. Results in the numerical study section shall examine these effects.

Numerical results

Example 1

Figure 1 shows the nonlinear transmission line given in Nouri and Nakhla (2018). This shall be used for the numerical studies in this paper. R = R_in = 1 Ω, L = 10 H, C = 1 F, I_d = I₀ (e^40vd − 1), I₀ = 1 A. The input is a pulse of amplitude 10 A with a rise and fall time of 1 s and a duration inclusive of rise and fall times of 10 s.

For the studies that follow, the number of transmission line sections is N = 10 unless otherwise stated. The simulation is for 40 s. The trapezoidal rule in conjunction with Newton’s method is used. Snapshots are taken every 80 time steps. The chosen time-step is 0.0031 s. Two measures of error shall be considered. The first is the root mean square error. It is defined as:

err=1Nt∑i=1Nt(y(t)−yr(t))2

where y(t) is the output from the full model and y_r (t) is the output from the reduced model. N_t is the number of time steps. This gives a measure of the error over the entire waveform.

The relative error is defined as:

relerr= |y(t)−yr(t)y(t)|

This gives an indication of where large errors occur in a time-domain waveshape.

The first important point to note is that the Jacobian matrices for this circuit are negative definite. However, the vector field for this circuit is not infinitesimally contracting. When the Jacobian J_F is formed, μ[J_F (t, x)] ≥ 0 for some t.

When implementing the IPOD, various factors affect the accuracy of the results.

The first factor is the training input. The training input should be a wideband signal and typical of that encountered in the regular application of the circuit. Pulses with sharp rise and fall times as employed by a pulse generator would be suitable for a nonlinear transmission line.

The next issue is the selection of the reduced-model order, m. For example, in Method 1, large errors can arise for certain values of m as seen in Figure 2. The explanation for this occurrence is that the selected rows of F(x_r) do not include sufficient states and hence, the algorithm produces an output that is in error. Note that the order (which row is selected first, second and so on) in which the rows are selected is dependent on U_F in (7).

Now consider the IPOD Method 2. In this case, a reduced-order model is guaranteed to be exponentially stable if the full model is exponentially stable. While this example is not infinitesimally contracting, noting that the Jacobians are Hurwitz, it is expected that the Jacobians of the reduced model would be Hurwitz. However, large errors arise again for certain m values because of the structure of the P matrix. The P matrix selects the most important equations in F to form the reduced model. The state variables involved in these selected equations are given by PPtVx^ [see equation (23)]. Until the P matrix is such that PPtVx^ includes all voltage states that are involved in the selected equations of F, there will be an error in the computed value of F and so the overall error will remain large. When all the voltage states are selected, there is a dramatic decrease in the error as evident from Figure 3. For the given set of snapshots and parameters, it is not until m = 19 that all voltage states are selected and at this point the error falls off. This situation is different with method 1 as in Method 1, the states that are selected are given by ΥtVx^ and hence, the voltage states that are involved in the selected equations are not determined or restricted by P. The consequence is that for Method 2 there is a minimum m for accurate results and this minimum may be significantly greater than that of Method 1. As in Method 1, the order in which the states are selected is dependent on U_F.

This feature of the results holds true for any number of transmission line sections, for example, when N = 5, the error remains large until m = 9 and when N = 20, the error remains large until m = 38.

The next issue is the selection of the snapshots in (2). The equations that are selected as important could differ when the interval between the snapshots change. However, tests carried out indicated that the snapshot interval did not have any significant effect on the results.

Observations

Based on the findings, the recommended strategy for model reduction of systems that have a contractive vector field is to firstly select Method 2 so that exponential stability of the reduced model is guaranteed. However, typical nonlinear transmission lines do not fall into this category. A lower order model is obtained with Method 1 if Method 2 is implemented as is. However, one recommendation is to decide what voltage states are required. The elements of P should be set to select these voltage states and this shall determine the value of m that is chosen. (It must also select the element of F that has the input to the circuit if this is not already included.) It is possible to pre-set P such that it selects the voltage states and this then determines the minimum order m for the IPOD. For the given transmission line, if P is such that PP^t selects all of the voltage states, then m = N + 1 where N is the number of transmission line sections. This reduction to O(N/2) is significant for a large N as occurs in practical realisations of nonlinear transmission line (NLTLs).

Figure 4 confirms this fact. It shows the error in the prediction of the end voltage of the transmission line with the pre-set P. As expected, the error falls off once all of the N + 1 = 11 voltage states are included.

The snapshot interval should be varied to check that there are no significant changes. Note that when P is pre-set, the snapshots affect V and not P.

Figures 5 and 6 show the transient results at the receiving end of the line. There is a larger error in Figure 6. This is expected owing to (25) not being always true for Method 2. A higher order m can be selected if greater accuracy is required. Figures 7 and 8 show the corresponding relative errors.

Example 2

The previous example had Jacobian matrices that were negative definite. However, this is not true in general of nonlinear transmission line circuits. For example, consider Figure 9 which is similar to that in Kuek (2012). The parameters chosen are also similar to Kuek (2012). The number of sections is 10. R_gen = 50 Ω, R_L = 0.16 Ω, R_load = 50 Ω, R_c = 2 Ω, C₀ = 816.14 pF, L = 1 μH, a = 2.137, b = 6.072 × 10⁻³. The nonlinear capacitance is given by:

C(v)=C0(b+(1−b)e−va)

The input pulse has a rise and fall time of 10 ns and a duration inclusive of the rise and fall time of 400 ns.

In this case, the Jacobian matrices are not always negative definite and μ[J_F] cannot be guaranteed to be negative. Indeed, forming a reduced order model for this type of circuit can be problematic as numerical issues arise when the Jacobian of the reduced vector field is unstable. For such cases, recommendations have been made in Bond and Daniel (2009). Stability issues also arose in Hasan et al. (2018). In this case, m should be increased until sufficiently large to avoid stability issues. Figure 10 shows the result from Method 1 when m = 19. It is not possible to obtain a stable reduced-order model for the given set of parameters using Method 2. This highlights the fact that forming a stable reduced-order model for a NLTL requires a thorough knowledge of the full model and its stability properties.

Conclusions

The benefit of or suitability of model reduction for circuits is dependent on a thorough knowledge of how various factors in the formation of the reduced-order model affect its operation. Choices of training inputs and snapshot intervals must be considered. Stability of the reduced-order model is paramount. When POD is used, stability is preserved (Chaturantabut, 2020). However, it does not lead to a reduction in complexity. Interpolated POD (Method 1) reduces the complexity but stability cannot be guaranteed to be preserved. Method 2 is less accurate at a cost of preserving stability for contractive systems. The paper has examined improvements to Method 2 for the specific application of nonlinear transmission lines. The paper has recommended a technique for selecting a minimal model order m and interpolation order k when the method is applied to nonlinear transmission line circuits that involve Jacobians that are Hurwitz. For nonlinear transmission lines where the Jacobians are not Hurwitz, neither method can guarantee stability and Method 1 is appropriate for such cases. The order m of such reduced models should be increased until a stable model is obtained.

Figures

Figure 1.

Nonlinear transmission line

Figure 2.

Error as m increases in Method 1

Figure 3.

Error in Method 2 with increasing m

Figure 4.

Error with Method 2 in the prediction of the endpoint voltage with pre-set P

Figure 5.

Result from Method 1 when m = 15

Figure 6.

Result from Method 2 with pre-set P when m = 15

Figure 7.

Relative error with Method 1 when m = 15

Figure 8.

Relative error with Method 2 with pre-set P when m = 15

Figure 9.

Nonlinear transmission line

Figure 10.

Result from Method 1 with m = 19