1 Introduction

Fault diagnosis of analog circuits is a fundamental problem for design validation [8, 10, 11]. A fault can be soft if a parameter is drifted from its tolerance range, but does not lead to any topological changes, or catastrophic which is defined as cut (open-circuit) of some connecting path or short-circuit of some pair of points. They are caused by contaminated phototools, defective raw material, improper solder leveling, mechanical damage, etc. [8].

If circuit simulations take place before any testing, the diagnosis approach is termed a simulation-before-test (SBT). The results of circuit simulations are stored as patterns in a fault dictionary [1, 6, 9, 13, 18, 19]. Comparing the voltages measured at testing nodes with the information contained in the dictionary, the faulty element can be identified. During the last decades a number of methods, concepts, and techniques have been proposed to build and exploit fault dictionary, e.g., sensitivity analysis [13], neural networks [1, 2, 6, 18], and Householder formula in matrix theory [24, 26]. The fault masking by component tolerances is a difficult problem faced by SBT approach.

Diagnosing of analog circuits can be performed in DC, time and frequency domains. For example, reference [16] offers a method for parametric and catastrophic fault detection and location of linear circuits in the frequency domain, whereas the papers [14, 20, 21, 29] bring diagnostic methods using the time-domain features.

Most work devoted to SBT approach has been focused on detection and identification of a single catastrophic fault. Building a fault dictionary usually requires a very large computing power, especially in the case of multiple faults in nonlinear networks, where the amount of computations increases dramatically [24, 26].

In this paper, fault diagnosis of nonlinear circuits using SBT approach, based on DC analysis, is considered. A method is proposed for detection and identification of a single catastrophic fault in BJT and MOS circuits having multiple operating points. The circuits are used, e.g., in signal conditioning applications, function generators and switching power supplies. Even if the tested circuit does not contain any capacitors and inductors, some parasitic capacitances and inductances exist. In DC analysis they are ignored and the circuit is considered as purely resistive. The DC operating points of the resistive circuit are the equilibrium points of the original circuit. Some of them are unstable and cannot occur. Generally, to indicate the unstable operating points among all operating points, the distribution and values of parasitic capacitances and inductances are required. If a circuit has multiple DC operating points the tested output voltage may assume different values for fixed value of the input voltage. The question which of the possible values actually occurs depends on the transient state which precedes the DC steady state. Since it is unknown, all the operating points should be taken into account during fault diagnosis. According to our knowledge, no work has been focused on fault diagnosis of such class of analog circuits till now. The crucial point of the proposed approach is tracing a large number of nonlinear multivalued input-output characteristics. For this purpose, an efficient algorithm is developed, based on the theory known under the name a linear complementarity problem [5, 7, 27], considering the deviations of the circuit parameters within their tolerance ranges. To perform diagnostic test, a DC input voltage source is applied at the node accessible for excitation and output voltages are read at the nodes accessible for measurement (usually one node). This is why the problem of selecting testing points of the circuit under test is not considered in this paper. Since the method usually needs only one node accessible for measurement it can be applied both to integrated circuits and PCB-based circuitries.

2 The Theoretical Background

This section gives the theoretical background of the method proposed for tracing of input–output characteristics in nonlinear circuits comprising bipolar and MOS transistors, having multiple DC solutions. Since large number of the characteristics is required by the diagnostic method, developed in Sect. 5, a fast and efficient algorithm for tracing the characteristics is necessary. If a circuit has multiple DC solutions, the input–output characteristic is not necessarily a single-valued function of the input, usually it is multivalued or even multibranched. The characteristic can be traced using a brute-force method which is capable of finding all the solutions for different values of the input voltage. Unfortunately, this method is very time-consuming and can be applied to small-sized circuits only. There are several more efficient methods, e.g., [3, 28, 30]. However, they are rather difficult to implement and usually suffer from major shortcoming due to sharp-turning-point problem. The method proposed in [23] overcomes this problem, but it is dedicated to just BJT circuits. On the other hand, SPICE simulator usually provides incomplete characteristics and often exhibits a sharp hysteresis loop, whereas the actual characteristic is \(Z\)-type.

In this section, an efficient algorithm for tracing input–output characteristics in circuits having multiple DC solutions is proposed, based on the theory known under the name a linear complementarity problem [5, 7]. The algorithm, which extends some ideas presented in [27], is described using BJT circuits, but it can be directly adapted to MOS circuits.

Consider a circuit consisting of bipolar transistors, diodes, resistors, and voltage sources. The transistors are characterized by the Ebers–Moll model [4, 27] shown in Fig. 1. We approximate the exponential characteristics of the diodes included in the model or acting alone using piecewise-linear functions similarly as in [27]. \(N\)-segment piecewise-linear characteristic is shown in Fig. 2a. The diode specified by this characteristic can be modeled by the circuit shown in Fig. 2b including \(N-1\) ideal diodes having the characteristic shown in Fig. 3a. It is convenient to choose reversed reference direction of the voltage across the ideal diode. Then the characteristic is as shown in Fig. 3b and the diode is described by relations

$$\begin{aligned} i\ge 0,\quad v\ge 0, \quad iv=0. \end{aligned}$$
(1)
Fig. 1
figure 1

Ebers–Moll model of NPN transistor

Fig. 2
figure 2

\(N\)-segment characteristic of diode (a) and its model (b)

Fig. 3
figure 3

Characteristic of the ideal diode

To trace transfer characteristic \(v_0=f(y)\), where \(v_0\) is the output voltage and \(y\) is the input voltage, we select from the circuit all the ideal diodes, the source \(y\), and the open-circuited branch with the output voltage \(v_0\) as shown in Fig. 4. As a result an \(m\)-port is created \(({m=n+2})\) consisting of linear resistors, current-controlled current sources, and independent voltage sources.

Fig. 4
figure 4

Circuit with extracted ideal diodes, input voltage source and the output open circuited branch

To form the hybrid representation of the circuit, we replace the ideal diodes by voltage sources and connect to the output terminals a zero current source \(i_{n+2}=0\) (see Fig. 5).

Fig. 5
figure 5

Rearranged circuit of Fig. 4

Using the hybrid representation of the \(m\)-port, we write

$$\begin{aligned} \left[ {{\begin{array}{l} {\hat{\varvec{i}}}\\ {v_{n+2}}\\ \end{array}}} \right] ={\varvec{H}}\left[ {{\begin{array}{l} {\hat{\varvec{v}}}\\ {i_{n+2}}\\ \end{array}}} \right] +{\varvec{s}}, \end{aligned}$$
(2)

where

$$\begin{aligned} \hat{{\varvec{i}}}=\left[ {{\begin{array}{*{20}c} {i_1 }\\ \vdots \\ {i_n }\\ {i_{n+1}}\\ \end{array}}} \right] ,{{\hat{\varvec{v}}}}=\left[ {{\begin{array}{l} {-v_1}\\ \vdots \\ {-v_n}\\ {v_{n+1}}\\ \end{array}}} \right] ,{\varvec{s}}=\left[ {{\begin{array}{l} {s_1}\\ \vdots \\ {s_{n+2}}\\ \end{array}}} \right] \end{aligned}$$
$$\begin{aligned} {{\varvec{H}}}=\left[ {{\begin{array}{l@{\quad }l@{\quad }l@{\quad }l@{\quad }l@{\quad }l} {h_{11}}&{} \cdots &{} {h_{1n}}&{} {h_{1,n+1}}&{} {h_{1,n+2}}\\ \vdots &{} \cdots &{} \vdots &{} \vdots &{} \vdots \\ {h_{n1}}&{} \cdots &{} {h_{nn}}&{} {h_{n,n+1}}&{} {h_{n,n+2}}\\ {h_{n+1,1}}&{} \cdots &{} {h_{n+1,n}}&{} {h_{n+1,n+1}}&{} {h_{n+1,n+2} }\\ {h_{n+2,1}}&{} \cdots &{} {h_{n+2,n} }&{} {h_{n+2,n+1}}&{} {h_{n+2,n+2} }\\ \end{array}}} \right] . \end{aligned}$$

Since \(i_{n+2}=0\), we remove the column \(n+2\) of the hybrid matrix \({{\varvec{H}}}\). Moreover, we remove \(({n+1})\)-st equation of the hybrid representation because we are not interested in current \(i_{n+1}\) and extract \(({n+2})\)-nd equation. As a result, we obtain

$$\begin{aligned} {{\varvec{i}}}={{\varvec{Mv}}}+\left[ {{\begin{array}{c} {h_{1,n+1}}\\ \vdots \\ {h_{n,n+1}}\\ \end{array}}} \right] y+\left[ {{\begin{array}{l} {s_1}\\ \vdots \\ {s_n}\\ \end{array}}} \right] , \end{aligned}$$
(3)
$$\begin{aligned} v_0=v_{n+2}=-\left[ {{\begin{array}{lll} {h_{n+2,1}} &{} \cdots &{} {h_{n+2,n}}\\ \end{array}}} \right] {{\varvec{v}}}+h_{n+2,n+1} \,y+s_{n+2} , \end{aligned}$$
(4)

where

$$\begin{aligned} {{\varvec{i}}}=\left[ {{\begin{array}{l} {i_1}\\ \vdots \\ {i_n}\\ \end{array}}} \right] ,{{\varvec{v}}}=\left[ {{\begin{array}{l} {v_1}\\ \vdots \\ {v_n}\\ \end{array}}} \right] ,{{\varvec{M}}}=-\left[ {{\begin{array}{lll} {{\begin{array}{lll} {h_{11}}&{} \cdots &{} {h_{1n}}\\ \end{array}}}\\ {\cdots \cdots \cdots \cdots \cdots }\\ {{\begin{array}{lll} {h_{n1} }&{} \cdots &{} {h_{nn} }\\ \end{array}}}\\ \end{array}}} \right] . \end{aligned}$$

Letting \({{\varvec{z}}}=\left[ {z_1\cdots z_n} \right] ^{\mathrm{T}}={{\varvec{i}}},\quad {{\varvec{x}}}=\left[ {x_1\cdots x_n } \right] ^{\mathrm{T}}={{\varvec{v}}},\quad {{\varvec{b}}}=\left[ {h_{1,n+1}\cdots h_{n,n+1} } \right] ^{\mathrm{T}},\) \({{\varvec{q}}}=\left[ {s_1\cdots s_n } \right] ^{\mathrm{T}},\) where \(\text {T}\) denotes transposition, we rewrite Eq. (3) in the form

$$\begin{aligned} {{\varvec{z=q+b}}}y+{{\varvec{Mx}}}. \end{aligned}$$
(5)

Since for each ideal diode, the relations (1) are fulfilled, we can write

$$\begin{aligned} {{\varvec{x}}} \ge {{\mathbf {0}}},\quad {\varvec{z}}\ge {{\mathbf {0}}},\quad {{\varvec{z}}}^{\mathrm{T}}{{\varvec{x}}}=0. \end{aligned}$$
(6)

To trace the multibranched characteristic \(v_0=f(y)\) for \(y\in \left[ {0,Y} \right] \), we propose an algorithm sketched in Sect. 3.

3 The Proposed Algorithm

Step1 Set \(y=0\) and write on the basis of (5) and (6)

$$\begin{aligned} \begin{array}{l} {{\varvec{z}}}={{\varvec{q}}}+{{\varvec{b}}}\cdot 0+{{\varvec{Mx}}},\\ {{\varvec{x}}}\ge {{\mathbf {0}}}\,,\quad {{\varvec{z}}}\ge {{\mathbf {0}}},\,\quad {{\varvec{z}}}^{\mathrm{T}}{{\varvec{x}}}=0. \\ \end{array} \end{aligned}$$
(7)

Although \({{\varvec{b}}}\cdot 0={{\mathbf {0}}}\), we retain this term in (7). The problem specified by (7) is called a linear complementarity problem (LCP) [5, 7, 27]. To solve this problem, we choose a positive vector \({{\varvec{d}}}\in R^n\) such that \({{\varvec{d}}}+{{\varvec{q}}}>{{\mathbf {0}}}\), using the procedure described in [7], and define the homotopy [7, 27]

$$\begin{aligned} \begin{array}{l} {{\varvec{z}}}^{\mathrm{T}}{{\varvec{x}}}=0\,,\quad {{\varvec{x}}}\ge {{\mathbf {0}}}\,,\,\quad {{\varvec{z}}}\ge {{\mathbf {0}}}\,, \\ {{\varvec{z}}}={{\varvec{p}}}-\lambda {{\varvec{d}}}+{{\varvec{b}}}\cdot 0+{{\varvec{Mx}}}, \\ \end{array} \end{aligned}$$
(8)

where \({{\varvec{p}}}={{\varvec{d}}}+{{\varvec{q}}}>{{\mathbf {0}}}\) and \(\lambda \) is a variable. At \(\lambda =0\) Eq. (8) reduces to \({{\varvec{z}}}={{\varvec{p}}}+{{\varvec{Mx}}}\) and the solution \({{\varvec{x}}}={{\mathbf {0}}}\;({\varvec{z}}={\varvec{p}}>{\mathbf {0}})\) of the homotopy system is obtained. At \(\lambda =1\), we have the original LCP (7). To trace the homotopy path and find the solution, we combine the homotopy approach with Lemke’s method as described in [7, 27]. During the computation process, we execute the same operations on both terms \({{\varvec{p}}}\) and \({{\varvec{b}}}\), although the last one does not affect the solution. On any stage of the procedure, some elements of vectors \({\varvec{z}}\) and \({\varvec{x}}\) are interchanged and \({{\varvec{q}}},{\varvec{b}}\), and \({{\varvec{M}}}\) are rearranged. At the end, when \(\lambda =1\), we obtain an equation equivalent to (7)

$$\begin{aligned} {\varvec{w}}={\hat{\varvec{q}}}+{\hat{\varvec{b}}}\cdot 0+{\hat{\varvec{M}}}{\varvec{u}}, \end{aligned}$$
(9)

where vector \({\varvec{w}}\) consists of some elements of vector \({\varvec{z}}\) and some elements of vector \({{\varvec{x}}}\), whereas \({{\varvec{u}}}\) consists of their complements. Equation (9) has the solution \({{\varvec{u}}}={\mathbf {0}}\;({{\varvec{w}} ={\hat{\varvec{q}}}}>{\mathbf {0}})\). Next we use Eq. (4), where vector \({{\varvec{v}}}\) is composed of all elements \(x_i\) selected from \({\varvec{u}}\) and \({\varvec{w}}\). As a result, we find voltage \(v_0\) at \(y=0\). To find other possible solutions at \(y=0\), we continue the procedure as described in [27]. In this way, we can find several solutions at \(y=0\) and the corresponding descriptions of the form (9).

Step 2 Similarly, we find the solutions at \(y=Y\). For this purpose, we modify Eq. (5)

$$\begin{aligned} {\varvec{z}}={\tilde{\varvec{q}}}-{\varvec{b}}\tilde{y}+{{\varvec{Mx}}}, \end{aligned}$$
(10)

where \(\tilde{\varvec{q}}={\varvec{q}}+{\varvec{b}}Y\), set \(\tilde{y}=0\) and write on the basis of (10) and (6)

$$\begin{aligned} \begin{array}{l} {\varvec{z}}=\tilde{\varvec{q}}-{\varvec{b}}\cdot 0+{\varvec{Mx}}\,,\\ {\varvec{x}}\ge {\mathbf {0}}\,,\,\quad {\varvec{z}}\ge {\mathbf {0}}\,,\quad {\varvec{z}}^{\mathrm{T}}{\varvec{x}}=0\,. \\ \end{array} \end{aligned}$$
(11)

Repeating the approach described in Step 1, we find the solutions at \(\tilde{y}=0\) (or \(y=Y)\).

Step 3 Form a new homotopy corresponding to the first solution provided by Step 1 (using (9)) with new homotopy parameter \(\lambda =y\), as follows:

$$\begin{aligned} \begin{array}{l} {\varvec{w}}^{\mathrm{T}}{\varvec{u}}=0\,,\quad {\varvec{u}}\ge {\mathbf {0}}\,,\,\quad {\varvec{w}}\ge {\mathbf {0}}\,, \\ {\varvec{w}}=\hat{\varvec{q}}+\hat{\varvec{b}}\lambda +{\hat{\varvec{M}}}{\varvec{u}}. \\ \end{array} \end{aligned}$$
(12)

At \(\lambda =y=0\), the solution of (12) is known. We apply again the method being a combination of the homotopy approach and the Lemke method [7, 27] adapted to (12). Each step of this procedure leads to such value of \(\lambda =y\) that one of the complementary pairs is \(x_k=0\), \(z_k=0\). It corresponds to a breakpoint of the piecewise-linear characteristic. To find this breakpoint, we select all \(x_i \quad ( {i=1,\ldots ,n})\), create the vector \({\varvec{v}}=\left[ {x_1\cdots x_n} \right] ^{\mathrm{T}}\) and use (4) to find \(v_0 \) at this breakpoint. To find other branches of the characteristic, we repeat the described approach taking into account in succession all the solutions found in Step 1 of the algorithm. In this way, all the branches which start from \(\lambda =y=0\) are traced.

Step 4 Repeat the approach described in Step 3 for (10) with \(\lambda =\tilde{y}\) starting in succession from these solutions obtained in Step 2 which have not been found in Step 3. As \(\lambda \) is increased, \(y=Y-\tilde{y}\) decreases.

This step is applied to find the possible branches which cannot be traced starting with \(y=0\), like the branch (c) in Fig. 6.

Fig. 6
figure 6

Characteristic \(y-v_0\) consisting of three branches (a), (b), and (c)

4 Circuit Example

The proposed method has been implemented in MATLAB 2012a and tested using PC Pentium i7-2600, 4GB.

Example 1

Let us consider the BJT circuit shown in Fig. 7, being a connection of the flip–flop circuit and the line receiver. Nominal values of the resistors are indicated in this figure. The parameters of the Ebers–Moll model of the transistors are as follows: \(\alpha _\mathrm{{F}} =0.99\), \(\alpha _\mathrm{{R}}=0.5\), \(I_\mathrm{{ES}}=7.07\;\text {fA}\), \(I_\mathrm{{CS}} =14.00\;\text {fA}\), \(V_\mathrm{{T}}=25.86\;\text {mV}\), \(R_\mathrm{{E}} =\text {10}\;\varOmega \), \(R_\mathrm{{C}}=\text {10}\;\varOmega \), and \(R_\mathrm{{B}} =\text {3}\;\varOmega \). The emitter and collector diodes are modeled by the circuit shown in Fig. 2, with \(N=8\) and the following parameters

Emitter diode

$$\begin{aligned}&\hat{R}_d=25.308\;\text {M}\varOmega ;\quad \hat{R}_1=11.044\;\text {k}\varOmega ,\quad V_0^{(1)}=0.475\;\text {V};\quad \hat{R}_2=2.060\;\text {k}\varOmega ,\\&V_0^{(2)}=0.535\;\text {V}; \hat{R}_3=430\;\varOmega ,\quad V_0^{(3)}=0.576\;\text {V};\quad \hat{R}_4=77.3\;\varOmega ,\quad \\&V_0^{(4)}=0.618\;\text {V};\hat{R}_5=13.4\;\varOmega ,\quad V_0^{(5)}=0.665\;\text {V};\quad \hat{R}_6=2.34\;\varOmega ,\quad \\&V_0^{(6)}=0.711\;\text {V};\quad \hat{R}_7=0.54\;\varOmega ,\quad V_0^{(7)}=0.749\;\text {V}. \end{aligned}$$

Collector diode

$$\begin{aligned}&\!\!\hat{R}_d=19.752\;\text {M}\varOmega ;\quad \hat{R}_1=16.256\;\text {k}\varOmega ,\quad V_0^{(1)}=0.445\;\text {V};\quad \hat{R}_2=2.885\;\text {k}\varOmega ,\quad \\&\!\!V_0^{(2)}\!=\!0.510\;\text {V};\quad \hat{R}_3\!=\!428\;\varOmega ,\quad V_0^{(3)}\!=\!0.554\;\text {V};\quad \hat{R}_4\!=\!76.3\;\varOmega ,\quad V_0^{(4)}\!=\!0.600\;\text {V};\quad \\&\!\!\hat{R}_5=13.8\;\varOmega ,\quad V_0^{(5)}=0.648\;\text {V};\quad \hat{R}_6=2.26\;\varOmega ,\quad V_0^{(6)}=0.694\,\text {V};\quad \\&\!\!\hat{R}_7=0.53\;\varOmega ,\quad V_0^{(7)}=0.732\;\text {V}. \end{aligned}$$

The diode \(D_1 \) is represented by the same model as the emitter diode with the series resistance equal to \(4\;\varOmega \).

Fig. 7
figure 7

BJT circuit for Example 1

In this circuit, we trace the input–output characteristic \(v_0=f( y)\), where \(y=v_\mathrm{{in}} \). For the nominal values of the resistors, we obtain the multivalued and multibranched characteristic shown in Fig. 8. The time consumed by the method is 0.32 s.

Fig. 8
figure 8

Input–output characteristic obtained using the method described in Sect. 2

To verify this characteristic we use the brute-force method [25] enabling us to find all the solutions of the circuit for any value of the input voltage without any piecewise-linear approximation. Plot of the characteristic obtained in this way is identical as the one provided by the proposed method. On the other hand the characteristic traced by PSPICE simulator, using up and down DC sweep analyses, is depicted in Fig. 9. A comparison of the characteristics manifests that PSPICE gives a fragmentary characteristic. Figure 10 shows a family of the characteristics traced by the proposed method for 100 sets of the parameters values \(\left\{ {R_1 ,\ldots ,R_{14}}\right\} \). They are obtained by random selection from their tolerance ranges \(({\pm 5 \%})\), assuming uniform distribution.

Fig. 9
figure 9

Input–output characteristic obtained using SPICE simulator

Fig. 10
figure 10

Banded input–output characteristic

The algorithm above-discussed can be adapted to the circuits containing MOS transistors, characterized by the Shichman-Hodges model built in Level 1 of SPICE [17]. For this purpose, the equivalent model, described in [22, 27] having the structure similar to the Ebers–Moll model of bipolar transistors, should be applied.

5 The Fault Dictionary

A method of building a fault dictionary for catastrophic fault diagnosis of circuits containing bipolar and MOS transistors, having multiple DC operating points, is developed in this section. The method is based on input–output characteristics of the circuit under test, traced for all considered catastrophic faults and fault-free circuit. Every time the deviations of the circuit parameters within their tolerance ranges are considered.

For fault-free circuit and for each circuit with a single fault, we trace a family of the input-output characteristics \(v_0=f(y)\) for different values of the circuit parameters randomly selected from their tolerance ranges assuming uniform distribution. If the number of faults is \(M\), we obtain \(M+1\) families of the characteristics. Each of the families has banded branches like the characteristic depicted in Fig. 10. We choose an interval \([{y^-,y^+}]\) of \(y\), divide it into \(N\) equal subintervals and consider the points \(y_k=y^-+kh\), \(h={({y^+-y^-})}/N\), \(k=0,1,\ldots ,N\) of this interval. For any of the families of characteristics, we find and store the ranges of \(v_0 \) variation at all these points. For example, in the case shown in Fig. 10 at \(y=5\;\text {V}\), we obtain four ranges of \(v_0\) variation: \(\left[ {-0.421,-0.357} \right] \), \(\left[ {-0.245,-0.179} \right] \), \(\left[ {0.268,0.390} \right] \), and \(\left[ {0.642,0.664} \right] \), all in volts. The results are summarized in a table having the structure shown in Fig. 11, labeled F, where \(F_i \), \(i=1,\ldots ,M\), denotes \(i\)th fault and \(F_0 \) means fault-free circuit.

Fig. 11
figure 11

Structure of table F

In the place specified by \(i\)th row \(({i\in \left\{ {0,1,\ldots ,N} \right\} })\) and \(j\)th column \(({j\in \left\{ {0,1,\ldots ,M} \right\} })\), the ranges of \(v_0 \) variation at \(y_i=y^-+ih\) in the circuit with fault \(F_j \) are stored. On the basis of table F, we create \(N+1\) tables labelled \(T_0,T_1,\ldots ,T_N \) having the structure shown in Fig. 12. Each of the tables exploits the information contained in one row of table F. For example, table \(T_l \) takes into account the information contained in \(l\)th row of table F, where the ranges of \(v_0 \) variation at \(y_l=y^-+lh\) are stored for all faults and fault-free circuit. The blank spaces of \(T_l \) in Fig. 12 are filled in as follows. At the crossing of \(i\)th row, corresponding to \(F_i\quad ({i\in \left\{ {0,\ldots ,M-1} \right\} })\) and \(j\)th column corresponding to \(F_j\quad ({j\in \left\{ {1,\ldots ,M} \right\} ,j>i})\), we insert 1 if the ranges of \(v_0\) variation for faults \(F_i\) and \(F_j\) (at \(y_l=y^-+lh)\) have a common part, or 0 otherwise.

Fig. 12
figure 12

Structure of table \(T_l\)

Next we consider tables \(T_0\) and \(T_1\) and form a new table \(T_{01}\) as follows. If in table \(T_0\), the element which appears in \(i\)th row \(({i\in \left\{ {0,\ldots ,M-1} \right\} })\) and \(j\)th column \(({j\in \left\{ {1,\ldots ,M} \right\} ,j>i})\) equals 1 and in table \(T_1 \), the element located in the same place equals 0, then in table \(T_{01}\) the element equal to zero is inserted in this place. Otherwise, the element of table \(T_0\) is transferred invariable to table \(T_{01}\). Similarly, we consider other pairs of the tables: \(T_0\) and \(T_2,\ldots ,T_0\) and \(T_N\) and create tables \(T_{02},\ldots ,T_{0N}\). This procedure is continued leading to tables: \(T_{12},\ldots ,T_{1N},\ldots ,T_{N-1,N}\).

If at any stage of this procedure, a table having all entries equal to zero arises, the algorithm is terminated and a fault dictionary enabling us to identify all the considered faults is built. Let table \(T_s \quad ({s\in \left\{ {0,1,\ldots ,N} \right\} })\) have all elements equal to zero. Then at point \(y_s=y^-+sh\) (corresponding to this table) all ranges of \(v_0\) variation relating to all the faults are disconnected (no of them overlaps other one). Thus, having the value of the measured tested voltage \(v_0 =\tilde{v}_0 \) at this point, we take into account \(s\)th row of table F and find such column \(j\) that \(\tilde{v}_0\) belongs to a range located in \(s\)th row and \(j\)th column. Index \(j\) defines the fault \(F_j \). In this case, the fault dictionary consists of \(s\)th row of table F. If table \(T_\mathrm{{pr}} \) contains all elements equal to zero, the fault dictionary consists of \(p\)th and \(r\)th rows of table F. If none of the tables is satisfactory (the sum of the entries is not sufficiently small), we can repeat the described above procedure taking into account for each table \(T_{ij}\) the tables \(T_m \quad ({m\ne i,m\ne j})\) and create tables \(T_{ijm}\). This procedure can be continued as long as the sum of the entries decreases. Otherwise, the procedure terminates. Since the applied model of transistor employs the piecewise-linear approximated characteristics (Fig. 2), the obtained regions of the output voltage variation appeared in the fault dictionary are corrected. For this purpose, the original transistor model is used, with the smooth characteristics described in Fig. 1, and the Newton–Raphson algorithm applied.

Suppose that the built dictionary consists of two rows of table F. Then in order to detect and identify the actual fault, we find all regions of \(v_0\) variation at the first of the rows which include the measured voltage (at the point corresponding to this row) and form a set of the corresponding potential faults. Similarly we treat the second of the rows. As a result, we obtain two sets of potential faults and choose the common part of the sets. If the fault dictionary has been built on the basis of a table having some entries equal to one, the corresponding faults may be distinguishable or undistinguishable as explained in Example 4.

6 Numerical Examples

Example 2

Let us consider the BJT circuit shown in Fig. 7. The transistors models and the circuit parameters are as described in Example 1.

We want to diagnose fault-free circuit (\(F_{0})\) and \(M=7\) catastrophic faults: cuts of the branches AB (\(F_{1})\), CE (\(F_{3})\), FG (\(F_{4})\), KL (\(F_{5})\), and short-circuits of the pairs of points CD (\(F_{2})\), BM (\(F_{6})\), AM (\(F_{7})\). To build the fault dictionary, the procedure described in Sect. 5 has been applied. We trace eight families of input–output characteristics \(v_0=f( y)\), where \(y=v_{in} \in \left[ {0,10} \right] \;\text {V}\) with \(h=0.1\;\text {V}\), as described in Sect. 3. To compute each of the families, 100 sets of the parameters values \(\left\{ {R_1,\ldots ,R_{14}} \right\} \) are randomly selected, assuming uniform distribution, within the tolerance ranges \(\pm 5\%\). The size of table F is \(101 \times 8\). We perform the procedure described in this section and find table \(T_{ij} \), where \(i=7\), \(j=28\) having one entry equal to 1, appeared in row 0 and column 5. The index \(i=7\) corresponds to \(v_\mathrm{{in}}=0.70\;\text {V}\), whereas \(j=28\), \(v_{in}=2.80\;\text {V}\). Hence, the fault dictionary consists of the rows 7 and 28 of table F as shown in Fig. 13. It allows identifying all the discussed catastrophic faults except \(F_0 \) and \(F_5 \), which may be undistinguishable or distinguishable. To verify the proposed method, the circuit was built using 5 % tolerance resistors and several catastrophic faults were considered experimentally. Two of them are described below. In the case of fault \(F_6 \), the measured voltages were \(-0.328\;\text {V}\) for the input voltage \(y=0.7\;\text {V}\) and \(-0.230\;\text {V}\) for \(y=2.8\;\text {V}\). They allow identifying the fault correctly. In another fault \(F_5 \), the measured voltages were: \(-0.523\) and \(-0.093\;\text {V}\), respectively, leading to correct diagnosis.

Fig. 13
figure 13

Fault dictionary for Example 2

The most time-consuming part of the algorithm is tracing input–output characteristics for each of the circuit states considering a large number of sets of the parameters values. Average time of tracing one input–output characteristic is 0.32 s. Since the number of the families of the characteristics, for fault-free circuit and 7 catastrophic faults, is 8 and each of them consists of 100 characteristics corresponding to 100 sets of the parameters values, the total number of the characteristics is 800. Hence, the time consumed for their tracing is 256 s = 2.27 min. This time dominates the others, used for selecting the test voltages (62 s) and correcting the voltage ranges using the Newton–Raphson algorithm (2.8 s).

Example 3

Consider the operational amplifier-based Schmitt trigger shown in Fig. 14, where the Ebers–Moll model parameters of the transistors are as in Example 1. The emitter and collector diodes are modeled by the circuit shown in Fig. 2, with \(N=4\) and the following parameters:

Emitter diode

$$\begin{aligned}&\hat{R}_d=3250\;\text {M}\varOmega ;\quad \hat{R}_1=51.5022\;\varOmega ,\quad V_0^{(1)}=0.61\;\text {V};\quad \hat{R}_2=11.5745\;\varOmega ,\quad \\&V_0^{(2)}=0.672\;\text {V};\quad \hat{R}_3=4.7551\;\varOmega ,\quad V_0^{(3)}=0.7077\;\text {V}. \end{aligned}$$

Collector diode

$$\begin{aligned}&\hat{R}_d=2663\;\text {M}\varOmega ;\quad \hat{R}_1=54.2765\;\varOmega ,\quad V_0^{(1)}=0.6\;\text {V};\quad \hat{R}_2=11.0595\;\varOmega ,\quad \\&V_0^{(2)}=0.6651\;\text {V};\quad \hat{R}_3=4.0290\;\varOmega ,\quad V_0^{(3)}=0.7000\;\text {V}. \end{aligned}$$

We want to diagnose fault-free circuit (\(F_{0})\) and \(M=9\) catastrophic faults: cuts of the branches DJ (\(F_{1})\), BC (\(F_{2})\), and short-circuits of the pairs of points HI (\(F_{3})\), CK (\(F_{4})\), GL (\(F_{5})\), AL (\(F_{6})\), IL (\(F_{7})\), EL (\(F_{8})\), and FM (\(F_{9})\). To build the fault dictionary, the procedure described in Sect. 5 has been applied. We trace 10 families of the input–output characteristics \(v_0=f(y)\), where \(y=v_\mathrm{{in}} \in \left[ {0,10} \right] \;\text {V}\), with \(h=0.2\;\text {V}\). To create each of the families, 100 sets of the parameters values \(\left\{ {R_1 ,\ldots ,R_8} \right\} \) are randomly selected, assuming uniform distribution within the tolerance limits \(\pm 5\,\%\). The size of table F is \(51\times 10\). Performing the procedure described in Sect. 5, we find Table \(T_0\), having all elements equal to zero. The corresponding input–output voltage is \(v_\mathrm{{in}} =0\;\text {V}\). Hence, the fault dictionary consists of the row 0 of the table F, as shown in Fig. 15. It allows identifying all the aforementioned catastrophic faults.

Fig. 14
figure 14

Circuit for Example 3

Fig. 15
figure 15

Fault dictionary for Example 3

Average time of tracing one input–output characteristic of the circuit shown in Fig. 14 is 7.6 s. Since the number of the families of the characteristics, for fault-free circuit and 9 catastrophic faults, is 10 and each of them consists of 100 characteristics corresponding to 100 sets of the parameters values, the total number of the characteristics is 1,000. Hence, the time consumed for their tracing is 126.7 min. This time dominates the others, used for selecting the test voltages (42 s) and correcting the voltage ranges using the Newton–Raphson algorithm (101 s). The total time is 129 min., but it is the off-line operation.

Although the proposed diagnostic method has been explained in detail for BJT circuits, it can be directly adapted to MOS circuits characterized by the Shichman–Hodges model built in Level 1 of SPICE [17]. It can be shown [22, 27] that this model is equivalent to the circuit having the structure of the Ebers–Moll model of bipolar transistors, with the gains of the controlled sources \(\alpha _F=\alpha _R=1\) and the diodes described by the equations

$$\begin{aligned} i_1&= i_\mathrm{{EF}}=\left\{ {{\begin{array}{lll} {k({v_\mathrm{{gs}}-\left| {v_{t_0}} \right| \,})^2}&{} {\text {for}}&{}{v_\mathrm{{gs}}\ge \left| {v_{t_0}} \right| }\\ 0&{} {\text {for}}&{} {v_\mathrm{{gs}} <\left| {v_{t_0}} \right| }\\ \end{array}}} \right. {,}\\ i_2&= i_\mathrm{{CF}}=\left\{ {{\begin{array}{lll} {k({v_\mathrm{{gd}} -\left| {v_{t_0}} \right| \,})^2}&{} {\text {for}}&{} {v_\mathrm{{gd}} \ge \left| {v_{t_0}} \right| }\\ 0&{} {\text {for}}&{} {v_\mathrm{{gd}} <\left| {v_{t_0}}\right| }\\ \end{array} }} \right. . \end{aligned}$$

In the above equations \(v_{t_0}\) is the threshold voltage, \(k=\frac{K_p}{2}\frac{W}{L}\), where \(K_p\) is the transconductance parameter, \(W\) and \(L\) are the channel width and length, respectively. Using this model, we can apply the theory developed in this paper. In the case of MOS circuits, however, the circuit parameters are not resistances, but the threshold voltages \(v_{t_0}\) and the coefficients \(k\), dissipated within their tolerance ranges.

Example 4

Let us consider the waveform-reshaping circuit [12] shown in Fig. 16. The nominal values of the channel width \(W\) and length \(L\) in \(\mu \text {m}\) are indicated in Fig. 16. Nominal values of the other parameters are as follows: PMOS\(-K_p =19.485\;{\mu \text {A}}/ \mathrm{V}^{2}\), \(v_{t_0}=-0.8351\;\text {V}\), \(R_d=16.4\;\varOmega \), \(R_s=16.4\;\varOmega \), \(R_g=0\); NMOS\(-K_p =79.173\;{\mu \text {A}}/ \mathrm{V}^{2}\), \(v_{t_0}=0.5705\;\text {V}\), \(R_d=16.4\;\varOmega \), \(R_s=16.4\;\varOmega \), \(R_g=0\). Having \(W\), \(L\) and \(K_p \), the coefficient \(k=\frac{K_p}{2}\frac{W}{L}\) can be calculated. In this circuit, we wish to detect and identify \(M=8\) catastrophic faults: cuts of the branches AB (\(F_{1})\), EM (\(F_{2})\), and short-circuits of the pairs of points HK (\(F_{3})\), GA (\(F_{4})\), MC (\(F_{5})\), GH (\(F_{6})\), HA (\(F_{7})\), and LM (\(F_{8})\), as well as fault-free circuit (\(F_{0})\). MOS transistors are characterized by the model described above where the diodes are approximated by a piecewise-linear 8-segment functions (see Fig. 2). To build a fault dictionary, we trace families of the characteristics \(v_0=f( y)\), as described in Sect. 3, and consider \(y=v_{in} \in \left[ {0.8,4} \right] \;\text {V}\) with \(h=0.1\;\text {V}\). To compute each of the families, 20 sets of values of the parameters \(k\) and \(v_{t_0}\) are randomly selected assuming uniform distribution within the tolerance ranges \(\pm 5\,\% \). The family corresponding to the fault-free circuit is shown in Fig. 17.

Fig. 16
figure 16

MOS circuit for Example 4

Fig. 17
figure 17

Family of the characteristics of fault-free circuit shown in Fig. 16

The size of the table F is \(33\times 9\). We perform the procedure described in Sect. 5 and find \(T_{17}\) having two entries equal to one, appeared in row 0 and column 3 and row 1 and column 6. It enables us to identify all the discussed faults except the fault-free circuit and the fault \(F_3 \), as well as \(F_1 \) and \(F_6 \), which may be undistinguishable, on the basis of the measured output voltages at \(v_\mathrm{{in}}=0.9\;\text {V}\) and \(v_\mathrm{{in}} =1.5\;\text {V}\) corresponding to table \(T_{17}\). Hence, the fault dictionary consists of rows 1 and 7 of table F, as illustrated in Fig. 18. This dictionary shows that faults \(F_0\) and \(F_3\) are undistinguishable, whereas \(F_1\) and \(F_6\) can be undistinguishable or distinguishable depending on the measured output voltages. For example, if \(\tilde{v}_0=4.9998\;\text {V}\) at \(v_\mathrm{{in}} =0.9\;\text {V}\) then fault \(F_1\) is identified. Average time of tracing one input–output characteristic of the circuit shown in Fig. 16 is 1.4 s. Since the number of the families of the characteristics, for fault-free circuit and 8 catastrophic faults, is 9 and each of them consists of 20 characteristics corresponding to 20 sets of the parameters values, the total number of the characteristics is 180. Hence, the time consumed for their tracing is 252 s = 4.2 min. This time dominates the others, used for selecting the test voltages (12 s) and correcting the voltage ranges using the Newton–Raphson algorithm (6.8 s).

Fig. 18
figure 18

Fault dictionary for Example 4

7 Discussion and Concluding Remarks

This paper is focused on catastrophic fault diagnosis of nonlinear circuits containing bipolar and MOS transistors using simulation-before-test approach. Unlike the other works which deal with the circuits having a unique operating point, this paper, for the first time, is devoted to the circuits having multiple operating points (DC solutions).

The proposed approach to the problem of fault diagnosis of analog circuits has a general meaning, because a catastrophic fault changes the circuit topology and forms a new circuit. The faulty circuit may possess multiple DC solutions, even if the original circuit has a unique solution. For example, the universal preamplifier circuit with the test voltage source \(v_\mathrm{{in}}=0.8\;\text {V}\) connected to the node accessible for excitation, shown in Fig. 19 with the following transistors parameters: \(\alpha _F=0.9975\), \(\alpha _R=0.8\), \(I_\mathrm{{ES}}=10.22\;\text {fA}\), \(I_\mathrm{{CS}}=12.75\;\text {fA}\), \(V_\mathrm{{T}} =25.86\;\text {mV}\), \(R_\mathrm{{E}}=0.81\;\varOmega \), \(R_\mathrm{{C}}=0.33\;\varOmega \), and \(R_\mathrm{{B}} =\text {3.3}\;\varOmega \), has a unique DC solution \(v_0 =19.473\;\text {V}\) at \(v_\mathrm{{in}}=0.8\;\text {V}\). However, if the catastrophic fault, being short-circuit of points AB, occurs then the circuit possesses three DC solutions \(v_0^{(1)} =0.438\;\text {V}\), \(v_0^{(2)}=10.052\;\text {V}\), and \(v_0^{(3)} =19.473\;\text {V}\) at the same input voltage \(v_{in} =0.8\;\text {V}\). In the real circuit built for verification purpose, with 5% tolerance resistors, the measurements gave two output voltages: \(v_0^{(1)}=0.452\;\text {V}\) and \(v_0^{(3)} =19.460\;\text {V}\).

Fig. 19
figure 19

Circuit possessing a unique DC solution

The input–output characteristic \(v_0=f({v_\mathrm{{in}}})\) of the original circuit depicted in Fig. 20a is a single-valued curve, whereas the characteristic of the faulty circuit, shown in Fig. 21a is a multivalued curve. These characteristics obtained experimentally are shown in Figs. 20b and 21b. Since usually we do not know whether the circuit possesses multiple DC solutions, the fault dictionary should be built using a method which allows finding multiple solutions. Such a method works correctly also in the case of a unique solution. The approach developed in this paper fulfills this requirement and leads to reliable results. This is the main achievement of this paper.

Fig. 20
figure 20

Input–output characteristic of the circuit shown in Fig. 19 traced numerically (a) and experimentally (b)

Fig. 21
figure 21

Input–output characteristic of the faulty circuit traced numerically (a) and experimentally (b)

Fault diagnosis of the circuits possessing multiple solutions is much more complex and requires larger computing power. This is why the CPU time is much longer than in the case of circuits having a unique solution. Although the method developed in this paper for tracing input–output characteristics is very efficient, necessity of finding a great number of the characteristics makes the diagnostic approach time-consuming. Fortunately, in the applied simulation-before-test approach this is the off-line operation. Consequently, the CPU time is not as crucial as in the simulation after test approach. The method is based on the concept of linear complementarity problem and is an extension of the idea presented in reference [27] devoted to finding multiple DC operating points. Since tracing multivalued and multibranched characteristics is a basic problem of the analysis of nonlinear dynamic circuits [15] the method is a contribution to nonlinear circuit theory.

The fault dictionary, built using the proposed approach, allows detecting and identifying the faults on the basis of a diagnostic test, which usually exploits measurements at one accessible node only. The method does not make very restricted demands about accuracy of the measurements. If the obtained fault dictionary is not satisfactory another testing node should be taken and the procedure repeated.

The algorithm of the fault diagnosis, proposed in this paper, does not accept the node approach, which is commonly used in user-oriented circuit analysis programs, including SPICE. Instead, it exploits the hybrid analysis [4] to formulate the linear complementarity problem. The hybrid representation of electronic circuits is employed in almost all methods which allow finding all the DC solutions and is required by many other methods in this field.

The proposed algorithm is difficult to implementation in a computer program, what is a disadvantage. In addition the Level 1 transistors models have to be used and their piecewise-linear representations exploited. Transistors of CMOS circuits manufactured in nanometer technology should be characterized by very complex BSIM 4 or PSP 103 model. Each of these models is specified by several hundred equations, mostly nonlinear. Methods for tracing input–output characteristics applied to this class of circuits are very time-consuming and inefficient. Since the proposed approach needs huge number of the characteristics, it is not suitable for fault diagnosis of CMOS circuits manufactured in nanometer technology.

Because in the circuit having multiple DC operating points, we do not know in advance to which region of the output voltage variation belongs the measured voltage we consider all of them. In some cases, it is possible to identify the fault, even if the ranges corresponding to this fault and another one have common parts. The approach can be also useful for generating the training data for neural networks if they are used as the fault classifiers.