UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Transform based algorithms for transient analysis of nonlinear networks Agnew, David George 1974

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Notice for Google Chrome users:
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.

Item Metadata

Download

Media
831-UBC_1974 A1 A395.pdf [ 4.02MB ]
Metadata
JSON: 831-1.0105073.json
JSON-LD: 831-1.0105073-ld.json
RDF/XML (Pretty): 831-1.0105073-rdf.xml
RDF/JSON: 831-1.0105073-rdf.json
Turtle: 831-1.0105073-turtle.txt
N-Triples: 831-1.0105073-rdf-ntriples.txt
Original Record: 831-1.0105073-source.json
Full Text
831-1.0105073-fulltext.txt
Citation
831-1.0105073.ris

Full Text

TRANSFORM BASED ALGORITHMS FOR TRANSIENT ANALYSIS OF NONLINEAR NETWORKS by DAVID GEORGE AGNEW B.Sc. U n i v e r s i t y of Calgary, Calgary 1969 M.A.Sc. U n i v e r s i t y of B r i t i s h Columbia, Vancouver 1971 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY i n the Department of E l e c t r i c a l Engineering We accept t h i s thesis' as conforming to the required standard THE UNIVERSITY OF BRITISH COLUMBIA NOVEMBER 1973 In p r e s e n t i n g t h i s t h e s i s i n p a r t i a l f u l f i l m e n t o f the requirements f o r an advanced degree at the U n i v e r s i t y of B r i t i s h Columbia, I agree t h a t the L i b r a r y s h a l l make i t f r e e l y a v a i l a b l e f o r r e f e r e n c e and study. I f u r t h e r agree t h a t p e r m i s s i o n f o r e x t e n s i v e copying of t h i s t h e s i s f o r s c h o l a r l y purposes may be granted by the Head of my Department or by h i s r e p r e s e n t a t i v e s . I t i s understood t h a t copying or p u b l i c a t i o n of t h i s t h e s i s f o r f i n a n c i a l g a i n s h a l l not be allowed w i t h o u t my w r i t t e n p e r m i s s i o n . Department of C_ A A T/A A ^  A ( ^ J ^ W / i ^ The U n i v e r s i t y of B r i t i s h Columbia Vancouver 8, Canada Date- f l . f r l k - H ^3 ABSTRACT In th i s t h e s i s , computer methods f o r the transient analysis of networks are investigated. Numerical transform techniques are developed to solve the d i f f e r e n t i a l equations a r i s i n g i n network simula-t i o n . Extensions to permit i n c l u s i o n of some nonlinear elements are considered. E f f i c i e n t methods f o r implementing the techniques are developed. For the transform techniques, e r r o r estimates are derived. Using these estimates, algorithms f o r the automatic determination of s o l u t i o n parameters are developed. Advantages over other numerical transform and numerical i n t e g r a t i o n techniques are revealed. For nonlinear networks, i t i s shown that use of a Newton-Raphson scheme for s o l v i n g nonlinear algebraic equations i s d i f f i c u l t when coupled with transform methods f o r s o l v i n g d i f f e r e n t i a l equations. Instead, an a l t e r n a t i v e technique i s developed. Steps which are e a s i l y generated, but which only approximate Newton-Raphson steps, are used. The implementation of the transform techniques and the nonlinear s o l u t i o n i s considered. A program using a sparse tableau form of network equations i s discussed. The program i s i n two sections. The f i r s t reads i n the network de s c r i p t i o n s , and writes a se r i e s of Fortran subroutines f o r performing the analysis e f f i c i e n t l y . The sub-routines must be compiled, and are used by the secon 1 part of the pro-gram to perform the actual analysis. Examples which i l l u s t r a t e the performance of the various techniques are presented. i i TABLE OF CONTENTS ABSTRACT TABLE OF CONTENTS LIST OF ILLUSTRATIONS LIST OF TABLES . ACKNOWLEDGEMENT . I. INTRODUCTION 1.1 Aim of the Thesis 1.2 Motivation for Computer-Aided Design 1.3 Background 1.4 Outline of Work I I . NUMERICAL LAPLACE TRANSFORM TECHNIQUES 2.1 Motivation. 2.2 Previous Work 2.3 Q u a l i t a t i v e Basis of Proposed Method 2.4 Quantitative Basis 2.5 High Speed Convolution 2.6 Transform E f f i c i e n c y I I I . ANALYSIS OF NONLINEAR NETWORKS 3.1 Problem Formulation 3.2 Aspects of the Problem 3.3 Pseudo-Newton Solution , 3.4 Special Techniques IV. IMPLEMENTATION OF ALGORITHM 4.1 Problem Formulation ' 4.2 Other Techniques Used V. EXAMPLES 5.1 Simple Linear Network 5.2 Saturating Inductor Example 5.3 Class B Complementary Push-Pull A m p l i f i e r . 5.4 Monostable M u l t i v i b r a t o r ' VI. CONCLUSIONS AND RECOMMENDATIONS FOR FURTHER WORK 6.1 Contributions of the Thesis 6.2 Transform Methods 6.3 Analysis of Nonlinear Networks i i i Page 6.4 Implementation 90 APPENDIX I - RESULTS PROM ANALYSIS OE SIMPLE LINEAR CIRCUIT 91 APPENDIX I I - RESULTS FROM ANALYSIS OF SATURATING INDUCTOR NETWORK 95 APPENDIX I I I - RESULTS FROM ANALYSIS OF CLASS B AMPLIFIER 98 APPENDIX IV - RESULTS FROM ANALYSIS OF MONOSTABLE MULTIVIBRATOR.. 98 REFERENCES 101 i v Figure 2.1 2.2 5.3 5.5 5.6 5.9 5.10 5.11 5.12 LIST OF ILLUSTRATIONS Page Linear Network Problem 7 7 I l l u s t r a t i v e Problem 2.3 E f f e c t s of Changes of Variable 29 2.4 E r r o r Behaviour Given by e^,(t) 32 2.5 Error Behaviour Given by e^' 33 2.6 Error Behaviour of e 33 iP 2.7 Aperiodic Convolution 41 2.8 Aperiodic Convolution from P e r i o d i c Convolution 41 2.9 Example f o r Accuracy C a l c u l a t i o n 43 3.1 Separation of Nonlinear Network 46 4.1 Example of Tableau Formation 4.2 Breakup of Tableau 5.1 Simple Linear Network '5.2 Behaviour of Simple Linear Network 71 60 62 68 Saturating Inductor Example 7 2 5.4 Behaviour of Nonlinear Inductor C i r c u i t 73 76 78 Class B A m p l i f i e r • 75 Transistor Model 75 5.7 Diode Model 5.8 Output Voltage from Class B A m p l i f i e r Current Drawn By A m p l i f i e r 79 Monostable M u l t i v i b r a t o r C i r c u i t 81 Transistor Model f o r Example 4 82 Triggering Waveform 83 5.13 Monostable Performance v 85 Triggered Switching of Monostable M u l t i v i b r a t Recovery of Monostable M u l t i v i b r a t o r v i LIST OF TABLES Table P a 8 e I. Comparison of Maximum Step Sizes f o r Given Relative Error 45 I I . Computing Times f o r Example One • 70 I I I . Computing Times for Example Two 72 IV. Transistor Parameters 76 V. Analysis of Quiescent Point of A m p l i f i e r ,A5. A6. 80 A l . Analysis of Simple Linear Network with Relative Accuracy. of 0.01 92 A2. Analysis of Simple Linear Network with R e l a t i v e Accuracy of 0.001 93 A3. Analysis of Simple Linear Network with R e l a t i v e Accuracy of 0.0001 94 A4. Analysis of Saturating Inductor Network with R e l a t i v e Accuracy of 0.01 95 .Analysis .of Saturating Inductor Network using Relative ^ Accuracy of 0.001 Analysis of Saturating Inductor Network with R e l a t i v e Accuracy of 0.0001 99 A7. Analysis of Class B A m p l i f i e r A8. Analysis of Monostable M u l t i v i b r a t o r 1 0 0 v i i ACKNOWLEDGEMENT I wish to thank my supervisor, Dr. A. D. Moore, f o r h i s many h e l p f u l suggestions and continuous encouragement while i n graduate school. Suggestions and comments from other f a c u l t y members and students are also appreciated. I would l i k e to thank Miss Norma Duggan and Miss Betty Cockburn f o r assistance i n typing the t h e s i s . F i n a l l y , f i n a n c i a l support from the National Research Council i n the form of a 1967 Science Scholarship from 1969 to 1973 i s g r a t e f u l l y acknowledged. F i n a n c i a l support f o r the p r o j e c t was also provided under NRC grant A-3357. v i i i 1 I. INTRODUCTION 1,1 Aim of the Thesis Computer-aided network analysis i s becoming an accepted part of c i r c u i t design. The success of numerical network simulation has come about because e l e c t r i c a l networks can be mathematically modelled very accurately, and because there are e f f e c t i v e numerical methods for s olving network equations. However, the computational cost of many e x i s t i n g numerical methods has proven excessive, and has held back the widespread a p p l i c a t i o n of computer techniques f o r network design. P a r t i c u l a r d i f f i c u l t y has been encountered i n s o l v i n g the d i f f e r e n t i a l equations associated with network simulations. Numerical i n t e g r a t i o n , whereby-derivatives 'are replaced with f i n i t e 'differences, i s the method normally used. However, the olde r , simpler algorithms become i n t r a c t a b l e when applied to problems with widely separated eigen-values. Recently developed approaches, which avoid eigenvalue r e s t r i c -t i o n s , are unwieldy. The d i f f i c u l t i e s mentioned are not normally encountered when analyzing networks by hand. Instead of numerical i n t e g r a t i o n , the Laplace transformation i s commonly used for manual s o l u t i o n of network d i f f e r e n t i a l equations. Unfortunately, the steps involved i n applying this powerful technique have not proven easy to automate f o r computer analysis. E x i s t i n g implementations have -not proven competitive with new developments i n numerical i n t e g r a t i o n . In t h i s thesis we present numerical techniques for using transform analysis to solve network d i f f e r e n t i a l equations. The 2 techniques appear competitive i n computing cost with recently developed i n t e g r a t i o n methods, as w e l l as bringing the benefits of transform analysis to numerical c i r c u i t simulation. However, before discussing the s p e c i f i c b e n e f i t s , we present some general reasons j u s t i f y i n g further-development of computer-oriented design methods, 1.2 Motivation for Computer-Aided Design Development of computer-aided c i r c u i t design methods has been motivated by several reasons. We consider some of those here. One d i f f i c u l t problem i s the analysis of large networks, that i s , those con s i s t i n g of many branches. Power systems and large scale integrated c i r c u i t s are examples. Even straightforward tasks such as w r i t i n g network equations become slow and tedious when done manually, and simple algebra even more so when the number of equations i s l a r g e . These are n a t u r a l tasks f o r automation. Besides replacing hand c a l c u l a t i o n s , the computer has made a variet y of new techniques possible. Algorithms f o r numerical optimiza- . 2 t i o n . have resulted i n new design techniques, such as those based on the 3 ' adjoint network . Programs which make design decisions, such as choosing component tolerances, have also been developed. A t h i r d reason why general purpose design programs have been developed i s the cost of programming the more e f f i c i e n t methods of performing numerical c a l c u l a t i o n s . A w e l l - w r i t t e n program can ensure that most of the programming need only be done once f o r many problems. F i n a l l y , general purpose simulation programs have reduced the work necessary at the breadboarding stage of design. These programs, s t a r t i n g from a simple description of the network, automatically genei~ate the network equations and solve for approximations to the responses. Using these programs, unsuitable designs may often be rejec t e d and optimum component values selected before attempting to construct a prototype c i r c u i t . The numerical techniques developed i n t h i s thesis are intended f o r s o l v i n g the network equations a r i s i n g i n c i r c u i t simulation pro-grams . 1,3 Background Recent advances i n numerical methods have brought about a f a r greater degree of success i n network simulation than was possible i n the past. A summary of two of the most s i g n i f i c a n t developments follows. Matrix methods have found almost u n i v e r s a l use i n computer-aided . network a n a l y s i s , and as indi c a t e d by network s t a t i s t i c s pub-l i s h e d by J e s s e l , network matrices are almost i n v a r i a b l y sparse. Consequently, development of sparse matrix techniques^ has had a s i g n i -f i c a n t impact, and has permitted analysis of networks f a r l a r g e r than was conceivable previously. Recent developments of a promising nature 6 7 • include use of sparse tableau formulations ' f o r network equations which are guaranteed to have a very sparse nature and p a r t i c u l a r charac-t e r i s t i c s which can be exploited, and techniques f o r e x p l i c i t or i n t e r -p r e t i v e coding of routines f o r s o l u t i o n of simultaneous l i n e a r equations, When the same problem must be solved repeatedly with only some numeri-c a l values changed, an e x p l i c i t l y coded s o l u t i o n provides much greater e f f i c i e n c y f o r s e v e r a l reasons. The indexing operations are done only at f i r s t , rather than for each s o l u t i o n . Calculations which would be repeated without change need only be done once'. Redundant operations such as m u l t i p l i c a t i o n s by 1 or -1 can be removed from the process. We make use of these techniques i n Chapter IV. The other technique of growing importance i s the use of im-p l i c i t integration- routines f o r obtaining transient response of dynamic networks (networks containing energy s t o r i n g elements). Probably the most common approach to such solutions has been to use a state v a r i a b l e 9 formulation of the network equations, and an e x p l i c i t (open-ended) i n t e g r a t i o n formula. However, that method i s severely r e s t r i c t e d because step sizes small i n comparison to the smallest network time constant are required to ensure s t a b i l i t y . Consequently an inordinate amount of computation i s necessary for solutions which change slowly r e l a t i v e to the shortest time constant. 11 12 I m p l i c i t i n t e g r a t i o n formulas ' can circumvent t h i s d i f -f i c u l t y . An obstacle which has held back - t h e i r a p p l i c a t i o n has been the requirement that at each time step a system of l i n e a r equations must be solved at l e a s t once. But with the advent of sparse matrix techniques, the a d d i t i o n a l computation i s reduced to such an extent that there i s a large o v e r a l l reduction i n computation. A further bene 13 f i t pointed out by Gear i s that the state v a r i a b l e formulation no longer holds any s i g n i f i c a n t advantage. The r e s u l t has been a trend to other formulations, p a r t i c u l a r l y nodal equations, which are more re a d i l y set up and permit greater e x p l o i t a t i o n of sparse matrix techni-ques. Gear has also developed e f f i c i e n t techniques whereby the step s i z e and order of multi'-step i m p l i c i t i n t e g r a t i o n routines may be chang 14 at each time step The multi-order i m p l i c i t i n t e g r a t i o n algorithms are not w i t h -out drawbacks. The ca l c u l a t i o n s required at each time step to choose optimum order and time increment can be considerable. Also s t a b i l i t y i s s t i l l not guaranteed i n general. While single-step A-stable"*"^ methods above order two e x i s t , they have not found general a p p l i c a t i o n because of the large amount of c a l c u l a t i o n required per time step. 16 There i s a tendency to avoid the i m p l i c i t multi-step methods above order two because they are not A-stable, although they do have bett e r s t a b i l i t y properties than e x p l i c i t methods. Furthermore, as has been pointed out by the authors of ECAP I I , the algorithms which choose the optimum order at each time step seem to favour low order formulas f o r 17 nonlinear networks . A consequence of the fact that low order methods seem most u n i v e r s a l l y e f f e c t i v e i s that short and consequently many time steps must normally be taken to obtain accuracy. The transform-based methods we s h a l l present are an a l t e r n a -t i v e - t o the integration- methods •diee-usse-d. --Our techniques r e t a i n the -goo s t a b i l i t y properties of the low order i m p l i c i t i n t e g r a t i o n methods. Other benefits w i l l be discussed as the method i s developed. 1.4 Outline of Work The techniques we propose use frequency domain information to obtain transient network s o l u t i o n s . In Chapter I I , a f t e r considering the j u s t i f i c a t i o n f o r such methods and some previous work i n the same area, we introduce a new approach to the problem of obtaining transient solutions f o r l i n e a r networks. In Chapter I I I , the method i s extended to nonlinear networks. The implementation i n a general purpose net-work simulation program i s considered i n Chapter IV. F i n a l l y , the l a s t two chapters present examples of the use of the techniques, and con-clusions . 6 I I . NUMERICAL LAPLACE TRANSFORM TECHNIQUES 2.1 Motivation As a t o o l f o r s o l v i n g d i f f e r e n t i a l equations, the Laplace transform has found wide a p p l i c a t i o n . I t has been used extensively f o r the s o l u t i o n of network problems for l i n e a r , time-invariant networks. However as a numerical technique f o r computer use, i t s a p p l i c a t i o n has been l i m i t e d , because t r a d i t i o n a l methods become d i f f i c u l t to apply ex-cept for small problems. Table look-up schemes only cover small, simple transforms. And residue calculus i s too c o s t l y and d i f f i c u l t f or use on high order networks, as w e l l as being unsuitable f o r d i s t r i b u t e d elements. However, there are some compelling reasons f o r using the La-place transform, i f e f f i c i e n t numerical techniques can be found. Systems i n v o l v i n g s t i f f d i f f e r e n t i a l equations can be handled without s p e c i a l d i f f i c u l t y ; that i s , the s t a b i l i t y problem inherent i n numerical i n t e -gration methods i s avoided. The l i n e a r i t y and time-invariance of net-works can be expl o i t e d using transform techniques. For example, t r a n s -18 f e r functions or submatrix reduction can be used to advantage with transform techniques. Using numerical i n t e g r a t i o n , there i s no b e n e f i t to be gained because at each time step c a l c u l a t i o n s must be performed to update every v a r i a b l e involved i n a d i f f e r e n t i a l r e l a t i o n (normally, capacitor charges and inductor f l u x - l i n k a g e s ) . Moreover, many d i s t r i b -uted .systems' can be handled with r e l a t i v e ease i n the frequency domain as compared to s o l u t i o n i n the time domain. This aspect has been 19 investigated by Silverberg . Further benefits suggested by Nakhla, 20 Singhal, and Vlach include ease of c a l c u l a t i n g d e r i v a t i v e s of time 7 responses, r i s e times, extrema of responses, and the a b i l i t y to handle piecewise n o n l i n e a r i t i e s . A d d i t i o n a l l y , by developing s u i t a b l e techniques enabling us to go with f a c i l i t y both from frequency- to time-domain and v i c e versa, we can hope to encompass a broader class of nonlinear elements. We pursue t h i s objective i n Chapter I I I , a f t e r considering techniques f o r l i n e a r time-invariant networks i n t h i s chapter. 2.2 Previous Work The benefits considered have encouraged a number of inves-t i g a t o r s to pursue numerical techniques based on the Laplace transform. We consider some e x i s t i n g techniques as a prelude to the present work. A l l the methods we consider only accomplish the inverse Laplace t r a n s -formation. The r a t i o n a l e f o r considering j u s t t h i s facet of a general process i s that e i t h e r simple functions ( f o r example, sinusoids or ramps) with known Laplace transforms w i l l normally serve as t e s t e x c i t a -t i o n s , or the impulse or step response i s s u f f i c i e n t because i t can be used i n a numerical time-domain convolution to determine the response 19 to a general input. This l a s t approach has been used by S i l v e r b e r g to f i n d the response of a d i s t r i b u t e d network to a general input. We consider f i r s t that category of methods based on having the transform to be inverted i n the form of a r a t i o of polynomials i n 21 22 s, the Laplace v a r i a b l e . Using t h i s form, L i o u and Thompson have developed an e f f i c i e n t technique using a s p e c i a l form of the matrix 23 exponential, and Corrington has considered expanding the quotient i n powers of ^, an expansion which i s r e a d i l y inverted. However, the requirement that the transform be a r a t i o n a l meromorph precludes general use of these ideas f o r several reasons. D i s t r i b u t e d element networks may y i e l d transcendental rather than r a t i o n a l functions of s. And even f o r lumped networks, the computational cost of obtaining the quotient of polynomials i s excessive f o r a l l but small networks. Topo^ l o g i c a l methods tend to increase i n cost extremely r a p i d l y above a ce r t a i n (quite small) s i z e of problem , and matrix methods have been observed to increase i n cost as n"*, where n i s r e l a t e d to the s i z e of the 25 problem Because of these d i f f i c u l t i e s , we r e s t r i c t further consideration to methods re q u i r i n g only numerical values of the transform f o r s p e c i -f i e d values of s. The computational cost of obtaining transforms becomes reasonable using th i s expedient, f o r i t i s normally easy to determine the imtnittances of i n d i v i d u a l network branches at s p e c i f i c values of s. Then i t i s possible to use everyday tools of network a n a l y s i s , p a r t i c u l a r l y Gaussian reduction of s t r i c t l y numerical matrices, to obtain the trans-forms of responses. Most methods capable of u t i l i z i n g such information are based on the numerical evaluation of the Complex inversion i n t e g r a l : x ( t ) = 7TT X(6) e ds (2.1) The d i r e c t approach, numerical quadrature, has been considered by Hur-26 27 witz and Sanger , but i t i s numerically expensive. On top of the considerable cost of evaluating the i n t e g r a l , there i s the necessity of recomputing i t f o r each new value of time. A common device for circumventing t h i s obstaele i e the process of f i t t i n g i n t e r p o l a t i n g functions to the transform at Selected values of s. Laguerre and Legendre polynomials have been popular for t h i s . Silverberg presents an extensive summary of r e l a t e d methods 9 Recently published methods promise answers at more reasonable st cost. One of these, using Pade approximations to the exponential e° u, 20 i s being employed i n the network analysis context by Nakhla . That approach involves evaluating the transform at the poles of the Pade approximation. Unfortunately, t h i s may require transform evaluation at a number of poles for each time point required ; approximations with up to 20 poles have been employed. Otherwise, the method i s s t r a i g h t -forward to program and provides a basis f o r e r r o r c o n t r o l . Two other methods, suggested by N e i l l and by S i l v e r b e r g , embody ideas which we s h a l l use. Both are numerical approximations to the i n v e r s i o n i n t e g r a l (2.1), using a f i n i t e i n t e r v a l along the Brom-wich contour: x(t) = 1 2irj X(s) eS t ds (2.2) where fi i s r e a l and f i n i t e . Both employ s u b s t i t u t i o n s to reduce the value of required for accuracy. 29 N e i l l suggests the f o l l o w i n g s u b s t i t u t i o n : X'(s) = X(s) - X(s+a) (2.3) which transforms to x ' ( t ) , from which x(t) = X ' ( 2 . 4 ) 1 - e The s u b s t i t u t i o n (2.3) w i l l normally reduce a function, X ( s ) , which i s -1 -2 0(s ) as s -> M , to a function X'(s) which i s 0(s ) as s -> ». X'(s) therefore becomes n e g l i g i b l e at a much smaller value of fi. N e i l l also advocates s e l e c t i n g equally spaced points i n an i n t e r v a l along a Bromwich contour on the imaginary axis In the s-plane , 10 30 so that the f a s t Fourier transform (FFT) can be used to perform the numerical i n t e g r a t i o n . I f the number of sample points along the contour i s chosen as an integer power of too, then the number of m u l t i p l i c a t i o n s necessary to obtain the i n t e g r a l f o r n values of t i s , at most, 2n l o g ? n . The idea i s not without weaknesses s however. The use of (2.4) leads to low accuracy for small tf because, for small t f x' i s divi d e d by a small number and the errors i n using the FFT to determine x' are consequently magnified f o r the early part of the s o l u t i o n . N e i l l ' s suggestion that t h i s part of the s o l u t i o n i s the most e a s i l y obtained part using other methods would seem to increase the o v e r a l l complexity of the method considerably. Furthermore, the requirement of equal spacing of points i n time and frequency over the en t i r e i n t e r v a l s f o r which x'(t) and X'(s) are s i g n i f i c a n t w i l l almost c e r t a i n l y necess-i t a t e many more sample points than necessary f o r functions which are r e l a t i v e smooth over large s u b i n t e r v a l s s a common situation,. The other process we consider i s a marching process i n time 31 32 suggested by Silv e r b e r g ' . Instead of (2.3) he suggests the s u b s t i -t u t i o n Z(s) = X(s) - £^p- ' (2.5) which gives x(t) = z(t) + x(0 +) As we s h a l l show In section 2.4, t h i s s u b s t i t u t i o n (2.5) y i e l d s a func-ti o n Z(s) which has lower order than X(s) as s -»•«>. The marching process i s ca r r i e d out by using the transform inversion only to obtain the function of time for a short i n t e r v a l . The time o r i g i n i s then reset to the end of the i n t e r v a l , and f i n a l values at the end of the l a s t i n t e r v a l are used as i n i t i a l values f o r the next. The i n t e r v a l s are chosen to be short enough that values of e x c i t a t i o n 11 do not change appreciably during one i n t e r v a l , so for that i n t e r v a l they can be treated as constants. With an i n t e r v a l so short, the i n v e r -sion i n t e g r a l need be evaluated only at the f i n a l value of t . A s i g n i f i c a n t improvement i s also effected by choosing a Brom-wich contour to the r i g h t of the imaginary axis . By t h i s expedient, the transform to be inverted becomes f a r smoother and a r e l a t i v e l y small number of sample points s u f f i c e s to provide accuracy. The process i s based on approximation of the formula . a+j°° z( t ) = T^T Z(s) e S t ds , a > 0 2« at e 2 I F Z(a + ju) e3bit da) (2.6, Choosing a l a r g e r a makes Z(s) a smoother function of to, so fewer samples are required to adequately represent the function; but a l a r g e r o also o"t increases the m u l t i p l y i n g factor e , which magnifies the errors i n approximating the i n t e g r a l , except f o r small t . Si l v e r b e r g presents h e u r i s t i c guidelines for choosing o and the sample spacing, from which he has obtained reasonable accuracy using, t y p i c a l l y , ten samples of Z(o + joj) f o r each i n t e g r a l evaluation. Although N e i l l has obtained one time point for each frequency domain sample, Silverberg's method has a d e f i n i t e o v e r a l l advantage. I t i s not r e s t r i c t e d to uniform spacing of points i n time, so does not waste large amounts of computation providing unnecessary points i n the time s o l u t i o n when the s o l u t i o n only admits of slow change. Neverthe-l e s s , reconsideration of both techniques w i l l lead to a s t i l l more e f f i c i e n t implementation of frequency domain concepts. 12 2.3 Q u a l i t a t i v e Basis of the Proposed Method We propose to adopt the f a s t F o u r i e r transform as a means of evaluating the i n v e r s i o n i n t e g r a l many times with l i t t l e a d d i t i o n a l e f f o r t compared-to a s i n g l e evaluation, along with a s u b s t i t u t i o n to reduce the bandwidth required i n the frequency domain (the i n t e r v a l along the Bromwich contour which must be considered). In addition we choose a Bromwich contour to the r i g h t of the imaginary axis i n the frequency domain. This w i l l reduce the number of samples required i n the frequency domain, but the m u l t i p l y i n g f a c t o r e a t i n (2.6) w i l l l i m i t the time i n t e r v a l over which accuracy can be maintained. There are two con-sequences: f i r s t , the s o l u t i o n must be obtained i n sections; and second, Silverberg's technique of using an i n t e r v a l so short that the e x c i t a -t i o n i s approximately a constant over the i n t e r v a l w i l l no longer s u f f i c e Consequently, we must consider e f f i c i e n t methods of obtaining approxi-mations to the Laplace transform i t s e l f , as opposed to i t s inverse. For t h i s , i t would seem l o g i c a l to look f o r an inverse procedure to an ef-f i c i e n t method of numerical Laplace transform i n v e r s i o n . However, the e x i s t i n g procedures, i f examined, do not s u f f i c e ; the technique we propose does. An error analysis i s used to e s t a b l i s h the e f f i c i e n c y of the process, and to determine the values of parameters involved. We then proceed to consider the technique as a convolution process that i s both e f f i c i e n t i n that i t uses the FFT to perform the convolution with minimum c a l c u l a t i o n , and successful i n that c e r t a i n p i t f a l l s of the pro-cess, when c a r r i e d out i n the time domain, are avoided. 2.4 Quantitative Basis The problem considered i n t h i s section i s the transient analysis of l i n e a r networks. We assume the networks can be put i n the 13 form of Figo 2„1. I — 5 * Sources y (t) i s I3> 2> Linear Network [ T(p), i ( t ) ] F i g . 2.1 Linear Network Problem The network i s described by an n x n square matrix 'operator T(p), whose elements consist of real-valued terms and terms i n v o l v i n g the d i f f e r e n t d 1 r t i a l operator p, where p = — and — = J dt. The network behaviour i s characterized by.an n-vector x ( t ) . The source vector y (t) i s augmented by adding s u f f i c i e n t zero elements to form an n-vector of e x c i t a t i o n . y (t) such that the network behaviour i s characterized by T(j>) x(t) = y e ( t ) which can be Laplace transformed to give the matrix equation T(S) X ( S ) = Y E ( S ) + Y Q ( S ) £ Y(s) (2.7) (2.8) Where Y (s) accounts f o r the i n i t i a l condition terms a r i s i n g from the o transformation r e l a t i o n s between p and s. For example, i f x ( t ) i s the vector of nodal voltages v ( t ) and y (t) i s the vector of currents i ( t ) impressed at the nodes, then the network r e l a t i o n s become 14 T( P) x ( t ) a ( i r + G + pC) v ( t ) = i ( t ) = y e ( t ) o r , i n the frequency domain, T(s) X(s) = [- r + G + sC] V (s) + 1 -° = i (s) + [c v«T) - - r (2.9) v( t ) dt] (2.10) s Y (s) However, we s h a l l not r e s t r i c t ourselves to t h i s set of v a r i a b l e s . We do assume that the problem i s s u f f i c i e n t l y well-posed that there are no impulsive sources or responses, that sources and responses are contin-uous, and that the problem has a unique s o l u t i o n , so that U(s) e x i s t s , where . •/.-" v*,''*"^5 U (s) = T _ 1 (s) (2.11) 'These assumptions "impose a minimum of r e s t r i c t i o n s f o r l i n e a r networks. P r a c t i c a l sources are always continuous, and r e s u l t i n continuous res-ponses. I f the problem does not have a unique s o l u t i o n , f o r example because i t i s disconnected, then i t i s of l i t t l e p r a c t i c a l i n t e r e s t . In terms of U(s), X(s) i s given by X (s) = U (s) Y (s) (2.12) The o v e r a l l problem of f i n d i n g a time domain response to a time domain input becomes that of s o l v i n g x ( t ) = {U(s) [Y o(s) + i f ( y e ( t ) ) ] } (2.13) using numerical techniques f o r the Laplace t r a n s f o r m , ^ , and i t s inverse, The derivation w i l l be c a r r i e d out i n terms of the inverse Laplace operation because the changes needed to'go i n the opposite d i r e c -t i o n w i l l become obvious further on. F i r s t , we observe that the s o l u t i o n i s given by the inverse Fourier transform m u l t i p l i e d by e C T t: 15 . . 1 X ( t ) - 2^J J a+j = X(s) e S t ds o~3r-at 2n X (a + ja>) e-1ut da) (2.14) In a p r a c t i c a l problem, X(s) i s r e a l i f s i s r e a l . Hence [ X(c~jw) f = X(a+jai) ( 2 ' 1 5 ) where * denotes complex conjugate. Consequently, the i n t e g r a t i o n i n t e r -v a l -can be halved: at x(t) = -— Re V I T X(c + jai) e 3 W L dw (2.16) 0 As a f i r s t step towards a numerical procedure, the i n t e g r a t i o n I n t e r v a l i s reduced to a f i n i t e i n t e r v a l . Truncating at a f i n i t e frequency a) = fl gives at /> ft x(t) = — Re X(a + jw) Q3"* da) (2.17) 11 J 0 Unfortunately, the value of Q required for good accuracy turns out to be unreasonably large. However, the problem i s made tr a c t a b l e by a sub-s t i t u t i o n process suggested by Silverberg: Z(s) = X(s) - ^p-z < ;(t) = Z"1 tZ(s)] (2-18) so that + x(t) = x ( 0 + ) [X(s) - ^ U . ] (2.19) instead of -1 x (t) = ^_ ± [X(s)] - (2.20) 16 Before proceeding with the a n a l y s i s , i t i s necessary to determine exactly what i s gained by thi s s u b s t i t u t i o n , and f o r t h i s the behaviour of X(s) for large s i s needed. This suggests the i n i t i a l value theorem, which does indeed lead to the correct r e s u l t . But using the common 33 derivation , the method i s suspect, as i t depends on s growing along a ray at less than 90° off the r e a l a x i s , and the contours of i n t e g r a t i o n we use are on, or p a r a l l e l to, the imaginary a x i s . Instead, we proceed by expanding the kernel e S t i n the Laplace i n v e r s i o n : ] r a + j r o s t , ' x(t) = ~ X(s) e o t ds J a-j°° ,a+1<» <» / sk 1 f°" J „, , v (st) . k „ a+j X(s) I ds 2 l T J J a-n- k=0 fc* L. V 1- \ s k X(s) ds (2.21) 2 ^ k=o k ! 0 - j c o which, f o r an x(t) which does not contain impulses or d i s c o n t i n u i t i e s , can be w r i t t e n as CO _^ 0-f-jOO x ( t ) = x ( 0 + ) + f (0 +) t + J L J JL. J . (2.22) k=2 " a - j r o Taking the transform of (2.22) y i e l d s -{- 0 0 0"f"j 0 0 X(s) = ^ + \ f (0 +) I f s k X ( s ) d s (2.23) s k=2 s J a-j» In general, X(s) approaches ^ 2 ~ J - a g g t e c o m e s l a r g e . But, using the s u b s t i t u t i o n (2.18), Z(s) approaches ~~2^~ (which i s also ^ L^T~) a s s s s becomes large. Now z(t) i s given by at " . . Z ( t ) = — Re Z(a+ jw) e J dco 0 at pfi . 6 Re Z(a + jw) e 3 " dta I T at + - — Re 0 Z(a + ja)) e 3 W t da) (2.24) at «ft . . z(t ) = ~ Re I Z(a + jo) du> + e (t) (2.25) ^f0 The truncation error e^ i s given by at e T ( t ) - V R £ Z(a + jo) e j W t da) (2.26) ft • For good accuracy, a large value of ft i s required i n order to keep small compared to z ( t ) . Normally, increasing ft w i l l improve the accuracy, and so the accuracy desired w i l l d i c t a t e the choice of ft. In most cases, the required value of ft i s s u f f i c i e n t l y large that Z(s) i s close to i t s asymptotic value for large s, namely —^—' j • I t i s possible (a + jai) that f o r some problems a smaller value of ft might give s u f f i c i e n t ac-curacy. However, we require from here on that f o r a l l problems the value-of ft chosen i s large enough that the asymptotic behaviour i s being ap-proached. At worst, t h i s should cause some extra computation and pro-vide b e t t e r accuracy than necessary. A l s o , we require that a « ft,. Again, we are free to choose ft large enough to meet t h i s condition. But i n f a c t , i t seems i n e v i t a b l e that t h i s c r i t e r i o n w i l l be s a t i s f i e d to meet reasonable accuracy r e -quirements. For otherwise, suppose o - ft. Then a i s large enough that Z(a) i s approaching the asymptotic behaviour 'large a o and therefore max |z(a + jai)| = |z(a + jfi) | toe [0,ft] That i s , the integrand i n the err o r function e , towards the lower l i m i t of i n t e g r a t i o n , i s comparable i n magnitude to the integrand i n the function determining z(t) . Consequently, the term e^,(t) w i l l l i k e l y not be small compared to the s o l u t i o n . Therefore we expect that ft >> a i s necessary for reasonable accuracy. With the above conditions on ft and a s a t i s f i e d , we make the foll o w i n g approximation: Z(cr + jco) z(0+) e _ z(0 +) ( 2 . 2 7 ) usft>>a .2 2 (a + jco) to Now e^ can be estimated from e 0 t f°° z(0 +) j u t , e m ( t ) - - — Re - 9 - E D W 0) _ x(0 ) e f cosut d a ) ( 2„28) J ft co2 To obtain numerical values, we must f i r s t choose s u i t a b l e values of t . From (2.24), at (t) * - — Re Z(a + jai) e J dec T T J n ' 0 at ^ Re f ' Z'(a + jco) e j U t da. (2.29) J 0 where (2.30) co £ 'ft Our approximation to z ( t ) , then, r e s u l t s from the inverse Fourier trans-form of a bandlimited function. By the sampling theorem 3 4 the approxi-mation t o z ( t ) i s completely represented by values at points i n time separated by . 1 TT A t " 2F = . ft Using samples at times t = 0, At, 2At, = ^ I = 0, 1, 2, ... the errors are given by ~~ ft / B 4 , x(0 ) e f COSOJA 0 e (£At) - - ' — ft 0) x(0 )e f cosiTTx 1 , i ^ dx nft J , x2 For % = 0, the i n t e g r a l has unit value, so e (o) „ -For other A, i n t e g r a t i o n by parts i s necessary: oo cos&Trx dx C O S AlTX , s m x , - l-n dx /oo sxnx , ~ x ~ d x Air = ( - l / - Air' { «oo _ . Air J 0 X J 0 » { - l / - Air { f - S i ( A I T ) } f x s i n t where Si(x) = J — — dt 20 35 For £ = + 1, we evaluate S i using tables , y i e l d i n g a value of the i n t e g r a l of about -0.117. For lar g e r £, use of the asymptotic s e r i e s for S i gives C O S £ T T X dx , . . I . • j - I T T T (-1) - £ T T {-j ~ 2 " £lT £ T T + 4! o 4 4 £ T T 2 ( - l ) J 2„2 T T £ + e (£) where e W l < 24 c 4 4 £ T T We can summarize the error behaviour, then, as £ = 0 e T (£At) - i ( 0 + ) Tffi . . +. aAt 0.117 x(0 ) e lift £ = + 1 - 2 ( - l ) & x ( 0 + ) e q & A t 3 2 T T £ ft £ = + 2, + 3. . . . ] } (2.36) (2.37) (2.38) The a l t e r n a t i n g sign and decreasing magnitude of the e r r o r as £ increases are reminiscent of Gibb's phenomenon, a consequence of removing the high frequency p o r t i o n of the spectrum. S u r p r i s i n g l y , there i s a ben e f i t to be gained here: i f the spacing of time points i s chosen on the basis of numerical e f f i c i e n c y , and i f an i n t e r p o l a t i o n routine i s used to provide values of responses at desired points which do not coincide with calcu-l a t e d p o i n t s , then the a l t e r n a t i n g sign of the error terms provides the p o t e n t i a l f o r p a r t i a l c a n c e l l a t i o n of the errors during the i n t e r p o l a t i o n . 21 We now consider the problem of f i n d i n g the response of a network, using only transforms and inverse transforms of functions with zero i n i t i a l conditions i n the time domain. From x(t) = u(t) * y ( t ) «t 0 ,t 0 ,t u(t - T) y(x) dx = C u(t - T) [y(x) - y ( 0 + ) ] dx U A + y ( 0 + ) u(x) dx 0 we obtain z ( t ) , using the s u b s t i t u t i o n (2.18): z ( t ) = x(t) - x(0 +) ,t + u(t - x) [y(t) - y(0 )] dx 0 + y ( o + ) u(x) dx + r° - y(0 ) U ( T > d T J o Transforming to the Laplace domain, Z(s) = U(s) j£{y(t) - y(0 )}, + - U(s) y(0 +) - J U(») y ( 0 + ) s s (2.39) (2.40) (2.41) (£.42) From (2.41), x ( t ) follows d i r e c t l y : x ( t ) = _ 1 Z(s) + x(0 +) The arguments advanced so f a r j u s t i f y i n t e g r a t i o n along a f i n i t e i n t e r -v a l on the Bromwich contour i n the Laplace transform i n v e r s i o n , and the use of sampled values i n the time domain, along with the estimated 22 errors i n doing so. But before we can, with confidence, apply numerical i n t e g r a t i o n to evaluate at ~ 0, z(t) " - — Re Z(a + jw) e J 0 J t dco (2.43) 17 J 0 i t i s necessary to consider the errors incurred i n using sampled values of Z(a + jto) . Since we. intend to make use of the FFT to perform the numerical i n t e g r a t i o n , i t might seem reasonable to use error estimates based on trapezoidal i n t e g r a t i o n , which i s s i m i l a r i n form. Such an analysis i s t y p i c a l l y based on f i n d i n g the error r e s u l t i n g from the i n c l u s i o n of each sample i n the estimate of the i n t e g r a l . But i n f a c t , the errors do not a l l add i n the same sense and estimates derived i n th i s manner are too conservative. Instead, we proceed as fo l l o w s . Beginning with the input-output r e l a t i o n s h i p X(s) = U(s) Y(s) we make a change of va r i a b l e X(s + a) = U(s + a) Y(s + a) (2.44) For a system with no r i g h t - h a l f plane poles, choice of a p o s i t i v e value of a means that the abscissa of convergence, a, i n the Laplace t r a n s -form inve r s i o n can be chosen as zero, so s becomes interchangeable with jto. Introducing the s h i f t i n g parameter a i s done s o l e l y as a convenience i n notation. In the inversion process, from examining i ro+i-» s t X ( t ) = 0—7" X ( S + a) G d S 2 i r l J a-j» (a+a)t c m ' = ® _ — — 1 X(a + a + jto) e3mz d i o ' (2.45) ^ ^ _ CO i t i s evident that the effects of a and a i n determining x(t) are i n t e r -changeable, and use of the s h i f t i n g parameter a i s therefore j u s t a 23 conceptual convenience. We w i l l examine some of the consequences of using sampled values i n the frequency domain, before attempting to estimate the act u a l errors incurred. The inversion i n t e g r a l i s being estimated from the trapezoidal i n t e g r a t i o n rule.* pa £ A t o a --co£At z(£Afe) = — » Re j Z(a + ito) e J dto J 0 a£At n / \ 17/2—1 „ A i. ~ g.e { ™~i + ) [Z(a + jkAto) e J ] -, w jr£AtAw + ~ Z(a + j~Aw) e } Aw (2.46) A- z where 2 A U T f e n e a " t r Z ( a ) K/2-1 jkMtAco . z(£At) - - Re + 1 [Z(a + jkAoj) e J ] N - i^£AtAw + ~ Z(a + j~Au>) e } (2.47) I t i s evident from (2.47) that the approximation z(£At) i s a p e r i o d i c function m u l t i p l i e d by an exponential function of t . The p e r i o d i c f a c t o r i s a Fourier series with period T ~ £At, where £ i s such that £AtAto = 2nr. 7f 2?F From (2.31), A t = - = — — , so the required value of £ i s N, and the period ' Q NAw T i s NAt ~ ~ . Evidently, the approximation cannot be v a l i d for a l l £. I f i t i s a close approximation f o r £ = 0, 1, ... N-l then i t cannot be a good approximation f o r £ = N, N+l, ... unless the s o l u t i o n i s a p e r i o d i c function m u l t i p l i e d by an exponential function of t , which i s not true i n general. Now we ask whether the approximation can be close f o r £ e [0,N-1]. This w i l l be examined i n more d e t a i l s h o r t l y , but the r e s u l t s may be anticipated by rearranging (2.46) and (2.47). z(£At) - - Re f Z(a + j oo) e j W & A t dto 1 1 J o a£At N/2-1 ,> (k-fl)Aw ioj£At = — _ R e £ Z(a + jto) e J du> 11 k=0 J kAto > a f c A t 1 Z(a) N / 2 - l jkAAtAto, 1- (i) Re { + I [Z(a+'jkAu,) e J ] e 11 2 2 k=l f Z[a + j(k+l)Au)] e J ( k + 1 ) ^ ^ t } A w (2.48) The l a s t l i n e i s equal to (2.46). We can regard the approximation as being i n the step Re Z(a + ju) e duo J kAto = ^| Re { Z(a + jkAuo) e J k £ A u ) A t + Z[a + j ( k + l ) A t o ] e j ( k + 1 ) £ A w A t } ( 2 > , 9 ) We now examine the approximation f o r d i f f e r e n t values of I i n [0, N-l] with the object of showing that the approximation becomes i n c r e a s i n g l y suspect as I increases. Setting £ = 1, «(k+l)Aw . _ ry / I • \ 3 "At , Re Z(a + jto) eJ dco v kAto « * f Re { Z(a + jkAo) e 3 k A u A t + Z[a + j (k+l)A«o] e 3 ( k + 1 ) A w A t } « *| R e { z ( a + j k A w ) eJk 2 Tr/N + z [ a + j ( k + 1 ) A w ] e J ( k + l ) 2 , / N } (2.50) I f Z(a + jto) only changes s l i g h t l y i n a frequency i n t e r v a l of width Ato, then the e n t i r e integrand Z(a + ju>) e 3 U t w i l l change only s l i g h t l y i n the i n t e r v a l , and the approximation i s reasonable. Now, s e t t i n g I = N/2, 25 Re „(k+l)Aw ' 3W2 Z(a -i- 3 to) e dw IcAco jk^AtoAt j(k+l)^AwAt = |^ Re {Z[a + jkAw] e + Z[a + j(k+l)Aw] e } = ( - l ) k ~ Re { Z[a + jkAw] - Z [a + j(k+l)Aa>] > (2.51) The k e r n e l e 3 W t changes considerably i n a frequency i n t e r v a l of width Aw. Assuming the behaviour of Z(a + jw) i s quite smooth i n the same i n t e r v a l , the r e s u l t may s t i l l be a reasonable approximation; but i t i s not as sure as when 1=1, Now suppose % - N - l . Then we have o(k+l)Aw jw(N-l)At Re Z(a + jw) e J dw ^ kAw Aw — Re > , .. , 1k(N-l)AwAt , n T -j(k+l) (N-l)AwAt, •X. Z[a + jkAw] e + Z[a + j ( k l l ) A w ] e } r j k ~ - j (k+1) ™ = l^- Re( Z[a + jkAw] e ' + Z[a + j(k+l)Aw] e • } (2.52) i w t The phase change of the kernel e i n a frequency i n t e r v a l Aw i s almost 2TT, and i f Z(a + jw) changes l i t t l e i n the i n t e r v a l , we expect that the i n t e g r a l should be r e l a t i v e l y small. However, the approximation, using only the end point values of the integrand Z(a + j a ) e ~ , W t , no longer r e f l e c t s the large change i n the k e r n e l , and the r e s u l t i s hi g h l y doubtful. This general r e s u l t , that the accuracy decreases i n the l a t t e r h a l f of the i n t e r v a l [0, (N-l)At] w i l l hold true regardless of whatever other devices (such as using the s h i f t e d function Z(s + a) ) are employed. In f a c t , we s h a l l a r b i t r a r i l y choose to ignore a l l 26 N r e s u l t s from the second h a l f of the i n t e r v a l , namely when I e [— + 1, N - l ] . Otherwise, t r y i n g to maintain accuracy throughout the i n t e r v a l [0, N-l] would require excessive computation. A r e l a t e d source of e r r o r a r i s e s when i n v e r t i n g the above process to determine frequency domain samples from time domain samples. For t h i s purpose, we must e f f e c t i v e l y consider a p e r i o d i c e x c i t a t i o n , and i t i s possible that the response of the network to t h i s p e r i o d i c input may be quite d i f f e r e n t from that due to the o r i g i n a l e x c i t a t i o n . Fortunately, the s h i f t i n g parameter can be used to control the errors a r i s i n g from t h i s source. In order to examine the mechanism involved, a number of s u b s t i t u t i o n s w i l l be introduced here. F i r s t , a modified response, x ' ( t ) , i s defined by taking the inverse transform of (2.44): ., . A , . -at x 1 ( t ) = x ( t ) e = {U(s + a) JC[ y ( t ) e ~ a t ] } (2.53) We introduce y'(t) as a modified input to the system: y'(t) = y ( t ) e" a t (2.54) and then proceed to further modify i t , forming y " ( t ) : y«'(t) = {u(t) - u ( t - T/2)} y'(t) (2.55) where y(t) i s the Heaviside step function. Introducing y"(t) i n (2.54) gives x"(t) k X - 1 {U(s + a) £ [ y " ( t ) ] } (2.56) which defines a new v a r i a b l e , x " ( t ) , r e l a t e d to x'(t) by x"(t) = x * ( t ) , t e [0, T/2] - (2.57) 27 and s i m i l a r l y to x ( t ) by x(t ) = x"(t) e a t , t e [0, T/2] (2.58) Next, a p e r i o d i c e x c i t a t i o n , y ' " ( t ) , with period T i s introduced. oo y"'(t) = I y"(t + £T) (2.59) £=_co The frequency domain desc r i p t i o n of t h i s e x c i t a t i o n only e x i s t s at f r e -quencies u> = nAta n = -2, -1, 0, 1, 2, ... (2.60) where Au = 2TT (2.61) "Defining A X , M ( J u ) . = U(a + ju>) Y , n ( J " ) (2-62) re s u l t s i n a transform which i s completely described by i t s values at dis c r e t e values of u. I t remains necessary to consider the errors com-mitted i n using X'"(jto) to determine x(t) . Before doing so, i t i s i n s t r u c t i v e to examine g r a p h i c a l l y the r e s u l t s of the su b s t i t u t i o n s considered so f a r . For t h i s , we con-s i d e r the analysis of the c i r c u i t i n F i g . 2.2. -A/W ( t ) Q v (t) = ; i ( t ) s v c ( t ) v (0 +) = 0 F i g . 2.2 I l l u s t r a t i v e Problem 28 We l e t y( t ) = v g ( t ) = u(t) (2.63) and x(t) = v (t) (2.64) Sketches of the r e s u l t s of the s u b s t i t u t i o n s can be seen i n F i g . 2.3. I t can be seen i n the f i g u r e that the e f f e c t of making the e x c i t a t i o n p e r i o d i c i s to introduce errors i n t o x"(t) f o r t e [0, T/2] due to the addition of the t a i l s of the previous cycles. To examine th i s phenomenon, observe that x(t) = x'(t) e a t = x"(t) e a t , t e [0, T/2] = '{x'"(t) - x " ( t + T) - x"(t + 21) ...} t e [0, T/2] Defining the "wrap-around" error i n x'" as e ^ ( t ) , where e'(t) = x"'(t) - x " ( t ) , t e [0, T/2] 0) i t follows that e'.(t) = x" (t + T) + x" (t + 2T) + ... at - x" (t + T) , t e t0, T/2] Now" x 1 (t) = x(t) e •at -at u(t - T) y(T> dx <"(t) -at -at u(t - T) y(x) dx t £ T/2 T/2 (2.65) (2.66) (2.67) (2.68) u(t - x) y(x) dx t > T/2 (2.69) F i g . 2.3 Eff e c t s of Changes of Variable 30 This permits an evaluation of the wrap-around e r r o r . e'(t) = x" ( t + T) ,T/2 (2.70) -at-aT f . , m . . . , = e u(t + T - x) y(x) dx J 0 (t e [0, T/2]) Now the e f f e c t of e' (t) on the c a l c u l a t i o n of x ( t ) must be estimated. w s ' x( t ) = x" e a t t e [0, T/2] = (x' M - e') e a t t e [0, T/2] (2.71) Defining , . A , at e (t) = e e w w (2.72) then e f t ) = x«"(t) e w at c(t), t e [0, T/2] and the desired r e s u l t follows immediately: e (t) * e w •aT PT/2 0 u(t + T - x) y(x) dx (2.73) (2.74) t e [0, T/2] We use (2.74) to determine a. I f we take the r e l a t i v e accuracy, E, i n a numerical approximation of a function f to be given by E = Error i n f (2.75) Value of f then the r e l a t i v e accuracy of our s o l u t i o n (just considering the wrap-around error) i s -aT e f t ) T/2 u(t + T - x) y(x) dx E(t) = (2.76) (t - x) y(x) dx t e [0, T/2] 31 Assuming U(s) has only l e f t - h a l f plane poles, then normally u(t + T - T) < u(t - T) (2.77) for t e [0, T/2]. This condition does not n e c e s s a r i l y hold; i n p a r t i -c u l a r , poles close to the imaginary axis may cause the condition (2.77) to be abrogated. But i f such poles are contributing to the response, then normally short i n t e r v a l s (and therefore small T) are required f o r accuracy, so u(t + T - x) - u(t - x ) , disregarding any impulses at the o r i g i n . Consequently, the i n t e g r a l i n the numerator of (2.76) should not be s i g n i f i c a n t l y larger than the i n t e g r a l i n the denominator, and t y p i c a l l y i t w i l l be much smaller. Therefore, we can s a f e l y —aT —aT choose e equal to the r e l a t i v e accuracy desired (e = E), to keep the wrap-around error at an acceptable l e v e l . This provides a s a t i s f a c t o r y basis for choosing a, given T. The choice of T w i l l be examined s h o r t l y . One could argue that a larger value of a would reduce the e f f e c t of wrap-around errors fur t h e r , and assure that they are not s i g n i -f i c a n t . Unfortunately, increasing a also increases the truncation errors i n using a truncated Bromwich contour i n the i n v e r s i o n process. For t h i s reason, a larger value of a i s not used. Because a value of a > 0 i s necessary to c o n t r o l wrap-around e r r o r s , the e f f e c t of t h i s choice of a on the truncation errors must be examined. From ( 2 . 3 3 ) , we have . aJlAt - . .„„ . x(0 ) e f cos&Trx . ,„ -,ON e T ( £ A t ) ^ - - ^ - — 2 — dx (2.78) Defining e^(t) as the truncation error i n x'(t) (as defined i n (2.53)), 32 » ' (£At) = e T (&At) -a&At XX(.L.L ^ COSX-TTX ~ Tfft dx (2.79) whose behaviour can be sketched as follows i n F i g . 2.4. A e TUAt> X X X X X — ! — _ — 1 -1 1 1 X 4 » i F i g . 2.4 Error Behaviour Given by e^(t) The corresponding error i n x"'(t) we denote by e^'(t) e - ( t ) = I e;( t + kT) 'T . u 'T The period T i s chosen to be an even integer (N) of time steps: T = NAt (2.80) (2.81) We also sketch e^', taking the l i b e r t y of replacing the points by a smooth curve through them, although the actual behaviour of the err o r function i s not that simple; F i g . 2.5 assumes that N i s even, r e s u l t i n g i n addition of overlapping e r r o r terms. I t follows that, i n the i n t e r v a l under consideration, t e [0, T/2] ><" (j>At) - e l ( M t ) + e' [(N-2.) At] (2.82) 33 F i g . 2.5 Error behaviour given by e^'. Now the e f f e c t of t h i s e r r o r on our approximation to x ( t ) must be examined. Since i t arises using a p e r i o d i c approximation to y ( t ) , we denote i t by e T p . a£At „ _ r n ^ (2.83) e T p (£At) = e^' (£At) e a&At t e [0,T] The term e a"" L magnifies the err o r function towards the end of the i n t e r v a l [0, T], as sketched i n F i g . 2.6. F i g . 2.6 Error Behaviour of e ^ 34 The e r r o r at t = 0 does not concern us, because i n order to use the s u b s t i t u t i o n (2,19), the i n i t i a l values of x must be computed independently of the rest of the transform. Also, since the s o l u t i o n i s only considered v a l i d for t $ T/2, the error for t > T/2 can be neglected Consequently, i t i s reasonable to require that the truncation e r r o r (from (2.38)) at t = A t be equal to the truncation e r r o r at t = T/2, giving 0.117 x(0 +) e a A t _ (2) x(0 +) e. a T / 2 rrfi Z 3 N,2 ' The r i g h t hand term i s m u l t i p l i e d by two because the use of a p e r i o d i c input causes equal contributions to the truncation e r r o r at t = T/2 from the current period of the input and from the f i r s t succeeding period of the input. From (2,84)* the number of time steps needed i s N - 3.7 e 1 / 4 a T ~ 1 / 2 a A t (2.85) or, f o r a s l i g h t l y more conservative estimate, N - 3.7 e 1 / 4 a T (2.86) In p r a c t i c e , to take f u l l advantage, of the FFT, the next l a r g e r integer power of two i s chosen for N. By r e l a t i n g the error expressions to the r e l a t i v e accuracy c r i t e r i o n defined i n (2.75), i t i s possible to s i m p l i f y the error expressions. Suppose that the desired r e l a t i v e accuracy w i t h i n the i n t e r v a l t e [0, T/2] i s given by E; then Error i n s o l u t i o n v a r i a b l e i Value of s o l u t i o n v a r i a b l e i (2.87) where the norm i s taken over a l l values of i . Since a was chosen such 35 that -aT then from (2.76) N & 3:7 E -1/4 Now from (2:38) and (2.87), ft may be determined. E 0.117 £ (0 +) TTft x. (t) (2.88) (2.89) As I measure of the truncation e r r o r , we have used the fo l l o w i n g con-venient choice for the norm: (2.90) A -= max x.(0 +) a l l i With, t h i i choice of norm, the worst behaved v a r i a b l e should meet the f i ^ t i i r l m e h t : t i l i n g t h i s norm, L ( 0 + ) i (2.91) r r_ -0-;-117-h ~ ' 'in max i i ± ( 0 + ) x T o^ - T Id tfiefefore B_uj-l7 max TTE i FtiftherBorlj if6m (2:32), i ± ( 0 + ) x.(0 +) (2.92) itta SO A t = ft T = N At ^ Ns ft 5ti B i t i t u t i n g i f i t o (2:88). :=I Nif/ft E = e a = — £n NTT (2.93) 36 E s s e n t i a l l y , (2.88), (2.92) and (2.93) provide the parameters necessary f o r the s o l u t i o n process. One further d e r i v a t i o n i s of i n t e r e s t here. From (2.61), Aw = ~ = 2 0 N ••• h - I In other words, the process gives — time domain values f o r frequency domain values, w i t h prescribed e r r o r s . This i s a considerable improve-ment over other methods discussed, which require several frequency domain values for each time domain value c a l c u l a t e d . Other comments on the e f f i c i e n c y of the process are given i n s e c t i o n 2.6. Some *further""i"emarfes on the "-error estimates are'in order. F i r s t , we have assumed that, at the cut-off frequency ft, the function Z(s) has decayed to the point that i t s behaviour i s approximately given by x(0 +) • ^ • For some functions, having e s p e c i a l l y large second or higher s order derivatives at t = 0 , t h i s may not be true. Fortunately, i t i s e a s i l y checked by examining the frequency domain samples; and when and i f t h i s d i f f i c u l t y a r i s e s , i t i s only necessary to choose a l a r g e r ft and r e s t a r t the s o l u t i o n . A second consideration concerns the addition of errors from d i f f e r e n t sources. The error parameters have been chosen so as to r e s t r i c t the i n d i v i d u a l errors to a given magnitude compared to the s o l u t i o n . In a worst case, the truncation and wrap-around errors w i l l add, to give twice the desired errors. However, the truncation errors are small except at the beginning and end of the time i n t e r v a l , and i f 37 N has been chosen somewhat conservatively as suggested e a r l i e r , then only at the beginning of the i n t e r v a l i s there a problem. Our approach has been simply to l i v e with the fact that there may be a double e r r o r committed at one point i n the i n t e r v a l . I f t h i s should prove undesirable i n some application, that error should not be hard to estimate and correct, using (2.38). Another point concerns the choice of N. Equation (2 .75) gives a minimum value, and i t i s normally desirable to go to a s l i g h t l y l a r g e r value than the minimum; but no arguments have been given f o r not choosing a much l a r g e r value. One disadvantage i s that, making use of the FFT, the computation necessary f o r the transforms goes up as N I082 N, so that the computation necessary for obtaining each time point increases. We expect that, by keeping N small, most of the work w i l l be i n s o l v i n g T(s) X(s) = Y(s) rather than i n the transforms. A second disadvantage i n choosing a l a r g e r N arises because the character of the s o l u t i o n w i l l probably change with time. I f N i s kept small, and the s o l u t i o n i s r e s t a r t e d reg u l a r l y over a new i n -t e r v a l beginning at the end of the l a s t i n t e r v a l , then the time step s i z e and other parameters are more l i k e l y to remain optimum. This follows from the observation that the time step s i z e i n a given i n t e r v a l i s based on the i n i t i a l value of the d e r i v a t i v e of the s o l u t i o n i n that i n t e r v a l . We conclude t h i s section with a" summary of the steps involved i n the process. 1. Determine i n i t i a l values of variables and accuracy desired, and use these to determine the parameters of 38 the process. From the i n i t i a l values, the parameters f o l l o w from ( 2 . 7 9 ) , (2.83) and ( 2 . 8 5 ) . 2 . Calculate y;UAt) = { y e U A t ) - y e ( 0 + ) } e " a £ A t = 2, . . . S N/2) 3. Use FFT to calculate YiaAu>). 4. Solve T(a + j£Aio) X'(£Ao>) = Y' (£Ato) + Y (a + j£Aco) fo r X' at £ «= 0, 1, . .., N/2 5. Find + Z(£Aco) = X* (£Ato) * (° ni at £ = 0, 1, . . ., N/2 a + j£Aco 6. Use inverse FFT to f i n d z ( t ) . 7. Calculate x(£At) = x(0 ) + z(£At) e <£ = 1, 2, ... N/2 As w i l l be seen when the implementation of these steps i n a program i s discussed, a s u i t a b l e organization of the network equations w i l l ensure that the programming and computing necessary f o r step 1. i s a small part of the o v e r a l l process. This a b i l i t y to determine solu-t i o n parameters automatically without a great deal of extra e f f o r t i s a desirable feature of any algorithm. 39 2.5 High Speed Convolution The idea of using a numerical convolution process to evaluate the response of a l i n e a r network or subnetwork to a time domain input 19 xs not new. For example, Sil v e r b e r g has suggested using such a process as part of a technique for f i n d i n g the response of a d i s t r i b u t e d net-work containing nonlinear elements. Unfortunately, use of numerical convolution e n t a i l s c e r t a i n d i f f i c u l t i e s . S i l v e r b e r g points out that the convolution i s awkward unless one adheres to a s i n g l e time step s i z e ; and i n general, the impulse response of a l i n e a r system w i l l contain impulses which must be dealt with separately from the rest of 3 6 the c a l c u l a t i o n . In a d d i t i o n , N e i l l points out the problem that although the input and desired response to a network may be r e l a t i v e l y smooth, the impulse response need not be. This can necessitate the use of time steps much smaller than would be required using other methods of s o l u -t i o n . The proposed method may be considered as an algorithm that i s based on convolution but avoids the d i f f i c u l t i e s described. The idea 37 of using the FFT to perform high speed convolution i s a concept used i n d i g i t a l s i g n a l processing. Although the FFT inherently gives a p e r i o d i c convolution, an aperiodic convolution i s p o s s i b l e , as pointed out by Stockham . The idea i s e a s i l y i l l u s t r a t e d with an example: we f i n d the aperiodic convolution of x ^ t ) with x 2 ( t ) , where X ; L ( t ) = x 2 ( t ) = y(t) - y ( t - l ) (2.95) This i s i l l u s t r a t e d i n F i g . 2.7. The r e s u l t i s equally possible, i f the functions are made per i o d i c with s u f f i c i e n t period (and the l i m i t s of the convolution 40 F i g . 2.7 Aperiodic Convolution Example i n t e g r a l are adjusted to cover j u s t one period): F i g . 2.8 Aperiodic Convolution from P e r i o d i c Convolution I t i s important to note that the a b i l i t y to recover the aperiodic convolution requires that both functions being convolved must be zero-valued for part of the period ( t y p i c a l l y the second h a l f of the period). The advantage of using the FFT to 'fake' an aperiodic convolution i n t h i s manner i s that the computation to perform the 2 convolution goes up as 0(N log2 N ) , instead of 0(N ), where N i s the number of samples required to represent the functions. To compare the above to the network analysis approach we have < 41 taken, assume the f o l l o w i n g problem: f i n d the response x(t) of a net-work to an input y ( t ) , where x( t ) = u(t) * y ( t ) (2.96) Using r e a l t r a n s l a t i o n i n the frequency domain, i t follows that x ( t ) e " a t = u(t) e~ a t * y(t) e ~ a t (2.97) Defining y'(t) as the p e r i o d i c extension (with period T) of r -at y'{t) e 0 $ t £ T/2 0 T/2 < t < T then y'(t) meets the requirements mentioned f o r aperiodic convolution, i . e . , that the functions convolved be zero-valued f o r the second h a l f of t h e i r period. Defining also u'(t) as the p e r i o d i c extension of u ( t ) e 0 $ t $ T then u'(t) w i l l approximately meet t h i s requirement f o r s u f f i c i e n t l y large a. By these devices, the speed advantage to be achieved by using the FFT to perform high speed convolution Is made av a i l a b l e for network a n a l y s i s . Besides that, certain concomitant benefits are gained. F i r s t l y , the e r r o r estimates from the previous section are ap p l i c a b l e , making possible an algorithm which automatically performs the ca l c u l a t i o n s i n an e f f i c i e n t manner. Secondly, the term u(t) e need never be evaluated as such, but only as U(s + a), so that the time domain impulse response need never be found. Consequently, the parameters of the transforms need to depend only on the behaviour of the inputs and outputs rather than the impulse or step responses of the network 42 considered independently of the inputs or outputs. In essence we have a method for f i n d i n g network responses rela t e d to convolution, but more e f f i c i e n t than methods using numerical quadrature i n the time domain. Also, the technique i s p o t e n t i a l l y both more e f f i c i e n t and l e s s res-t r i c t e d as to inputs than other methods which can only i n v e r t Laplace transforms. 2 .6 Transform E f f i c i e n c y I t i s not easy to compare d i f f e r e n t methods of s o l u t i o n . The only sure basis of evaluation seems to be extensive experience with a v a r i e t y of approaches. However, by applying the e r r o r estimates we have derived to a simple example, a rough comparison of our method with others may be gained. Before considering an example, however, we make some general comments. As pointed out i n section 2.4, the method w i l l produce one time point f o r each frequency point evaluated, so one network s o l u t i o n f o r each time point i s required. Assuming that the bulk of the work i s normally i n the network s o l u t i o n s , then t h i s i s roughly equivalent to the work required per time point using the low order i m p l i c i t m u l t i -step i n t e g r a t i o n formulas. Admittedly the solutions required are complex-valued whereas only a r e a l network s o l u t i o n i s involved i n i n t e g r a t i o n schemes; but counter balancing that i s the advantage that we do not have to perform c a l c u l a t i o n s to update one v a r i a b l e or input for each energy s t o r i n g element at each time step. I f these values are not otherwise needed, then c a l c u l a t i n g t h e i r values at the end of each i n t e r v a l w i l l s u f f i c e to begin the next- i n t e r v a l . The method bears other s i m i l a r i t i e s to the low order m u l t i -step i m p l i c i t i n t e g r a t i o n methods. I t has the good s t a b i l i t y properties 43 of those methods and i t i s most e f f e c t i v e when moderate accuracy w i l l s u f f i c e . High order i n t e g r a t i o n methods s t i l l seem optimum when extreme accuracy i s necessary, because with high order methods, the improvement i n accuracy can be enormous when the i n t e g r a t i o n step s i z e i s only decreased moderately. In contrast, from (2.32) and ( 2 . 9 2 ) , i t i s ap-parent that the transform process developed i n t h i s chapter has a l i n e a r r e l a t i o n s h i p between step s i z e and errors. I t only remains to show that the time step sizes necessary are not unduly short. For t h i s , the t r i v i a l example of F i g . 2.9 w i l l s u f f i c e : + v (0 +) = 1 v. c "1 F. 1 ft' v F i g . 2.9 Example for Accuracy C a l c u l a t i o n The s o l u t i o n f o r t > 0 i s v(t) = e _ t (2.98) Suppose i t i s assumed that a r e l a t i v e accuracy of 0.001 i s desired; that i s , i n the f i r s t i n t e r v a l (T/2), absolute errors no larger than 0.001 are desirable. Then « = - ~ 3 7 - 2 • ( 2 - 9 9 ) and At 0.087 . ( 2- 1 0 0> 44 The required N i s given by N = 3.7 ( O i O O l ) " 1 ^ (2.101) = 21 The next larger power of two would be used> giv i n g ah N of 32i The time i n t e r v a l over which the S o l u t i o n would be found, then, would be T N - - - /- - - • ~ = | At - 1,39 (2.102) Let us compare t h i s with a f i r s t order Integration process. We assume that a step s i z e of h i s used* and that the magnitude of the f e c a l truh-,1,2 , . . . . . . . cation e r r o r Is | v(£)| ~ In going from t to t 4- hj where % e [ t j t * h] In the I n t e r v a l T/2, there w i l l be =jj=- steps takenj so the o v e r a l l e r r o r Is For ah appfoxImatlQhj we choose the average absolute Value §f v i n [0, §]: v(|) a f = |v(t')| at ( i . i o 4 ) Then 6 s | te(0) - e(T/2)j ' (1.165) Betting e = Q-.OOlj h Is appfoximately 9.BQ26. A s i m i l a r calculation. f6 f a Second order processj assuming the l o e a l tfuheatieh e r f o f i s ap-h 3 . . . r. r proximately v(£) ==• gives a step s i z e ef about 0.088* appfexlSately the time Step S i z e we are able t© "use. In f a b l e 1; a g§mpari§bh i§ made between maximum step sizes and absolute errors for the tfansfdrm and f i r s t - and second-order -Integration methods. 45 Maximum Step Sizes f o r Transform Method F i r s t Order Integration Second Order Integration 0.01 16 0.87 0.02 0.245 0.001 0.087 0.0026 0.088 0.0001 64 0.0087 0.00082 0.024 Table I Comparison of Maximum Step Sizes for Given Relative Error Because the above numbers are a l l based on many approximations, the values should not be regarded as exact. What they do i l l u s t r a t e , however, i s that i n addition to any other advantages transform methods possess, the transform process as developed i n t h i s chapter may poten-t i a l l y compete on the basis of 'numerical efficiency with ether methods i n common use. The numerical transform methods discussed o f f e r an a l t e r n a t i v e to e x i s t i n g methods for so l v i n g the d i f f e r e n t i a l equations i n l i n e a r network anal y s i s . In the next chaptar, we consider extending the tech-niques to permit analysis of nonlinear networks. 46 IIIo ANALYSIS OF NONLINEAR NETWORKS 3.1 Problem Formulation The a b i l i t y to handle nonlinear elements i n a tra n s i e n t a n a l y s i s i s a d e f i n i t e assetc We s h a l l consider an approach f o r dealing with these elements which i s compatible with the transform techniques discussed i n the l a s t chapter*, Unfortunately i t i s not generally f e a s i b l e to handle the nonlinear elements i n the frequency domain. Some simple c h a r a c t e r i s t i c s can be handled i n the frequency domain; f o r example, instead of f [ x ( t ) ] = x 2 ( t ) (5.1) one could t r y to make use of F [ x ( s ) ] = X(s)*X(s) (5.2) but the r e s u l t i n g computations would be complicated and c o s t l y , and lac k i n g i n g e n e r a l i t y . Consequently the network must be partitionedo Linear subnetwork T(p),x(t) 2> 2> /Sources and nonlinear elements f ( x , t ) Separation of Nonlinear Network F i g . 3.1 47 The d e s c r i p t i o n of the nonlinear functions f ( x , t ) must be considered. 'We assume time-invariance of the network,, Flux-linkages and charge w i l l have to be allowed as inputs to the l i n e a r subnetwork,. Thus, a capacitor may be characterized by the r e l a t i o n s q. = f c ( v ) (3.3) i = p q (3.4) The f i r s t equation would be incorporated i n the vector f u n c t i o n f ( x , t ) and the second r e l a t i o n i n the l i n e a r matrix operator T(p). S i m i l a r l y ^ f o r an inductor, X = f A ( i ) (3.5) v = p X • (3.6) The choice of the voltage- and c u r r e n t - c o n t r o l l e d r e l a t i o n s (3 .3) and (3.5) i s not mandatory; i n general, charge c o n t r o l l e d capacitors and f l u x - c o n t r o l l e d inductors may be allowed. The only requirement i s that the problem be separable so that a l l the n o n l i n e a r i t i e s are described by algebraic functions, and a l l l i n e a r c h a r a c t e r i s t i c s , i n c l u d i n g a l l d i f f e r e n t i a l r e l a t i o n s , are described by the matrix operator T(p)„ The problem reduces to one of f i n d i n g a s o l u t i o n of T(p) 3E(t) = f [ x ( t ) , t ] f (3.7) We require that,, as before, f i s a continuous function of t . Moreover, i n the work that follows,, we w i l l be developing techniques r e l a t e d to Newton's method f o r s o l v i n g nonlinear algebraic equations. Consequently, i t i s necessary that f be a continuous, d i f f e r e n t i a b l e function of x. In the derivations which follow, we assume that T(p) x 2 ( t ) = f [ x 1 ( t ) , t ] (3.8) can be solved f o r x^r given In some cases, t h i s may present a d i f f i c u l t y i n that f i s being considered as an e x c i t a t i o n . For example, 48 there can a r i s e a s i t u a t i o n i n which a node i s being driven e f f e c t i v e l y by current sources, because the current l e a v i n g a node i s no longer a function of the voltage at the node 0 This can, be avoided by considering a r e l a t e d problem, T>(p) = T(p) - C X = f ' [ x ( t ) , t ] = f [ x ( t ) , t ] - C x r * / \ • r /& f , ( t ) (3.9) where c = 0 f o r k £ j . Because we have assumed time-invariance,. a i k given element f, of f ( x , t ) can be a f u n c t i o n of x or t but not both. • I n e f f e c t , a l i n e a r component has' been subtracted out of each : nonlinear function and i t s e f f e c t considered as part of the l i n e a r operator. The system (3.9) w i l l have the same solutions as the o r i g i n a l system ( 3 . 7 ) ; ' but now : " -T'(p) x (t), = f [ x 1 ( t ) , t ] (3.10) can be solved f o r x^ i n general. A s u i t a b l e r e l a t e d system w i l l be generated automatically i n the s o l u t i o n process. We mention the possi-* b i l i t y of analyzing the r e l a t e d problem (3-9) only to avoid any concep-t u a l d i f f i c u l t y with the d e r i v a t i o n . The transform process of the l a s t chapter i s used to solve (3»8). In the frequency domain T(s) X 2 ( S ) = J K f t z ^ t ^ t ] + P 0 ( s ) x 2 ( t ) = X - 1 U 2 ( s ) ] . (3.11) 49 P (s) i s introduced to account f o r the i n i t i a l conditions i n going from T(p) to T(s). The equations are solved numerically, using the techniques previously described. The o v e r a l l process w i l l be described by an operation U : • x 2 ( t ) = U . f [ x 1 ( t ) t t ] - (3.12) U i s therefore the inverse operation to x(p)o 3,2 Aspects of the Problem We consider f i r s t a type of problem f o r which the s o l u t i o n i s e a s i l y obtained, and then the d i f f i c u l t y of obtaining a more general s o l u t i o n . An obvious approach to f i n d i n g a f i x e d point of (3»7) i s to t r y the P i c a r d i t e r a t i o n x . + ] _ ( t ) = U . f [ i . ( t ) , t ] (3.13) Because of the time-invariance of the network, f can be broken up i n t o two components: f U , t ] = f J x ] + f j t ) (3.14) Then Defining 1+1 a i o then u ( t ) = j£ 1 [ T ~ 1 ( s ) 3 (3.16) 2 i + l ( t ) = J 0 u(t-x) f a[x.(x)JdT + U f b ( t ) (3.17) This i s a nonlinear equation of the V o l t e r r a type, and standard techniques s u f f i c e to examine the convergence. For brevity,,we consider only the scalar case.. Assume 50 | U ( T ) I - m f o r t e [0,T] and l e t the f o l l o w i n g L i p s c h i t z condition he s a t i s f i e d : f(2^) - f ( x 2 ) < k x± - x 2 f o r r e a l and x^o Defining the mapping A by x i + l ( t ) = A S i ( t ) A f1" u(t - r ) f [ x . ( x ) ] dx + U f (t) J Q 3~ = A 1 x Q ( t ) then from any two s t a r t i n g points ^ ( t ) and y Q ( t ) such that |A x Q ( t ) - A y 0(-t)| < mkn dx = mknt (3.18) (3.19) (3.20) x0 "^0 < n' (3.21) from which r r r 1 .r /, \ ,r /.\I m k nt A x n ( t ) - A y n ( t ) < (3.22) The requirements of the contraction mapping theorem are s a t i s f i e d f o r a large enough r , so the i t e r a t i o n (3.15) w i l l converge to the s o l u t i o n . Most n o n l i n e a r i t i e s can be adequately modelled so as to s a t i s f y (3.19), so that many p r a c t i c a l problems-can be solved i n t h i s manner. However, there are also many problems f o r which (3.18) i s not s a t i s f i e d , , u s u a l l y because u ( t ) contains impulses. Consequently, other techniques must be considered. An obvious approach i s to attempt the use of Newton's method. Assume we have an approximate s o l u t i o n x ^ ( t ) such that T(p) x ± ( t ) = f [ x ± ( t ) , t ] + e.(t) ' (3.23) where e.(t) i s introduced to account f o r the err o r i n x (t)„ The change 51 6x. i s desired, such that T(p) x. + ]_(t) = f [ x i + 1 ( t ) , t ] (3*24) where x. (3.25) (3»26) c ± + ] ( t ) = x ± ( t ) + 6x ± ( t ) ¥e compute 6x^ approximately by l i n e a r i z i n g f: T(p) x . ^ C t ) = f [ x , ( t ) , t ] + f [ x . ( t ) ] 6x . ( t ) l r i 1 X I 1 where f ( x) 1 S the Jacobian matrix 3x^ 3x^ ^x^ 8x 2 ¥e assume that the d e r i v a t i v e matrix f [x ( t ) ] e x i s t s and i s piecewise continuous,, To solve f o r the approximation x ± + 1 ( t ) based on the approx-imation to 6x,(t) ? we proceed as follows: 1 T(p) z . + 1 ( t ) = f [ x . ( t ) s t ] + { f z [ x . ( t ) ] 6x.(t) + f x [ x . ( t ) ] x ±(t)>- f x[x.(t)] x i(t) = f [ x . ( t ) ] + f x[x.(t)] z. + 1(t) - f z [ x . ( t ) ] x.(t) {T(p) - f x [ x . ( t ) ] } x„ + 1(t) = f [ x . ( t ) ] - f x [ x . ( t ) ] x.<t) (3.27) Transforming both sides, {T(s) - J f ? [ f x ( x . ( t ) ) ] * J X . + 1 ( s ) -= < ^ { f t x i ( t ) , t ] - f x [ x . ( t ) ] x ± ( t ) } + F Q ( s ) (3.28) 52 The appearance o f the convolution term means•that the s o l u t i o n , X ± + 1 ( s ) , i s a generalized function,, While no impediment i s seen i n extending the numerical transform techniques to encompass t h i s problem, we expect that the complexity and the amount of computation necessary would be increased g r e a t l y . Instead, i n the next s e c t i o n we develop a pseudo-Newton process, which avoids the d i f f i c u l t i e s mentioned. 3.3 Pseudo-Newton Solution Beginning with the Newton i t e r a t i o n previously considered, we w i l l derive a method of i t e r a t i v e l y generating the Newton steps without the frequency domain convolution necessary i n (3.28). This i t e r a t i v e generation of Newton steps can be considered a s u b - i t e r a t i o n w i t h i n the Newton i t e r a t i o n . A l o g i c a l s e r i e s of s u b s t i t u t i o n s , combining the s u b - i t e r a t i o n with the main Newton i t e r a t i o n , leads to a s i m p l i f i e d process we c a l l the pseudo-Newton s o l u t i o n . We begin with the Newton i t e r a t i o n given by.(3.27), namely and use E . ( t ) ^ f [ x . ( t ) , t ] - f s U . ( t ) ] z.(t) (3.29) A l s o , the f o l l o w i n g notation w i l l be used: (3.30) A rearrangement of (3°30) gives { I - f x iU} T(p) x . + 1 ( t ) = E.(t) - (3.31) Premultiplying both sides by the same s e r i e s , { i + f r iu + [f x iu] 2 + ...> a - >t(P) W t } - {i + f ±u + t ^ i ^ 2 + > \W (3<.32) 53 Provided the terms [f . ^ u ] n T(p) xi+]_("0 approach zero as n approaches i n f i n i t y ^ t h i s reduces to T(p) s i + 1 ( t ) = { I + f ^ T J + [ f 2 i U ] 2 + • } E.(t) (3.33) This can be solved i t e r a t i v e l y . x±fo = X i T ^ \,1 = \™ T ( P ) x . y 2 = E.(t) {I-H-f^TI} E.(t) T(p) x = E (t) + f x - { I + f .U + [ f .U] 2 + o o . + [ f .Tj]^ }E.(t) x i + l ~ zi,» XI XI (3.34) The underlined part of (3<>34) defines the general step of the i t e r a t i o n . Unfortunately ? l i t t l e i s known about the convergence of the s e r i e s . However? there e x i s t s a s u b s t i t u t i o n which should increase the l i k e l i h o o d of convergence, and the speed of convergence, assuming convergence i s obtained. As argued i n secti o n 3 „ 1 , r e p l a c i n g T and f by T! = T(p) - C (3.35) and f \ = f . - C (3.36) x i x i x w i l l lead to the same s o l u t i o n as the o r i g i n a l problem. C. i s a p a r t i c u l a r value of the matrix C discussed i n section 3.1 » E^(t) i s unchanged by these s u b s t i t u t i o n s , so the i t e r a t i o n can be described by T i ( p ) X i + l ( t ) = ~ { I + f x i U ' + t f x i U ' ] 2 + } E i ( : t ) ( % 3 7 ) 54 x i the s o l u t i o n w i l l be computed i n sections of length T/2, we f . - C. I f o r t-:e [t„ , t, +I/2].where t We choose 0. so that f^,. i s a minimum,, Since, using the transform proce,= a c t u a l l y want C to minimize i s the value of time at the s t a r t of a given Intervale Depending on the choice of norm and the method of computation, t h i s minimization could be quite complicated; so we have adopted the expedient of choosing G i " f x ^ X i ( t 0 + T / 4 ) ^ ( 5 o 3 8 ) which, f o r smooth functions^ i n most cases should r e s u l t i n a value of the norm close to the minimum.0 A f u r t h e r m odification that w i l l s i m p l i f y the c a l c u l a t i o n consider-'; ably suggests i t s e l f * The general step of the i t e r a t i o n to determine one Newton step i s , from (3o32) and (3 .34-) , \M - + 1 - E ± ( t ) + f s ± % j ( 3 . 3 9 ) I t i s c h a r a c t e r i s t i c of Newton's method that, when -'close 1 to a s o l u t i o n , each step reduces the error by a considerable f a c t o r . Consequently* i f ( 3 . 3 7 ) provides an i t e r a t i o n which converges to a Newton step, we can assume that each step of that i t e r a t i o n provides a b e t t e r estimate of the s o l u t i o n . Then i t should be possible to accelerate the i t e r a t i o n by using the new estimate of the s o l u t i o n to c a l c u l a t e a new estimate of E ( t ) , denoted by E ^ ( t ) , which i s c l o s e r to the f i n a l value of E ( t ) y that i s , to the value of E(t) when the i t e r a t i o n has converged to a solution,, The modified i t e r a t i o n i s T. (p) x. ., = E. (t) + f . x. . i v ^ i , j + l 1 xr 1, j = f [ x . , ( t ) , t ] » f X . + f X. . -L» 3 x 1, j x 1, J = f [ x . , ( t ) t t ] - f x. . + f x . - C. x. . ^t'J x r , j x 1, j 1 55 T.'(p) x = f [ x ( t ) , t ] - 8 x. . ( 3 . 4 0 ) But with t h i s l a s t step, the douhle index has l o s t i t s value, so the i t e r a t i o n can be s i m p l i f i e d : T ±(p) = f [ x ± ( t ) , t ] - C. x ± E f * [ x ± ( t ) , t ] ' ( 3 o 4 l ) Added i n s i g h t i n t o the process may be gained by considering i t from d i f f e r e n t viewpoints 0 For one- the i t e r a t i o n has been reduced to the Pi c a r d i t e r a t i o n discussed i n section 3 ° 2 y but applied to a modified system so that the steps taken i n the i t e r a t i o n are c l o s e r to those of Newton*s methods This should speed convergence and broaden the class of problems f o r which convergence i s achieved,, Furthermore,, r e w r i t i n g ( 3 ° 3 9 ) as [T(p> - C ±] x , + 1 = f [ x . ( t ) , t ] - C. x. ( 3 . 4 2 ) and comparing i t to ( 3 ° 2 7 ) * i t i s evident that the i t e r a t i o n i s simply Newton's method but with the funct i o n of time, f [x ( t ) ] r replaced by a constant, 'average' value„ F i n a l l y , we can also consider the i t e r a t i o n as r e p l a c i n g the step ( 3 * 3 4 ) , namely with x i + l x _ = x. . ( 3 . 4 3 ) i + 1 1*1 That i s ? we assume one step convergence of the i t e r a t i o n ( 3 « 3 4 ) j so that each step of the pseudo-Newton method i s not quite a Newton step* Of course, Newton's method sometimes f a i l s to converge, and the i t e r a t i o n developed i n t h i s section only approximates Newton's method, 56 so the convergence properties may he changed, r e s u l t i n g i n a smaller region of convergence. Consequently, f o r a general purpose algorithm, i t i s desirable to consider a device to b r i n g about convergence i f the techniques discussed so f a r should f a i l . 3°4 Special Techniques Two aspects of the nonlinear s o l u t i o n merit s p e c i a l consideration. F i r s t , we introduce some devices which may be of use i f the pseudo--Newton i t e r a t i o n s f a i l to converge. Then, a device to increase the e f f i c i e n c y of the numerical transform process^ as applied to the nonlinear problem, w i l l be explained. Various options could be adopted i f the pseudo-Newton i t e r a t i o n f a i l s to converge. For example, methods have been proposed by Beywood 41 42 •.and.J5o.ore .and,by J f e i l l to force .convergence of P i c a r d i t e r a t i o n s . Another p o s s i b i l i t y would be to use the i t e r a t i o n (j.34) several times f o r each new Ih (t) c a l c u l a t e d , so that the pseudo-Newton steps become ' c l o s e r to true Newton steps, and then employ a controlled-step Newton 43 44 method such as suggested by Broyden and Branin . However, e i t h e r of these options would complicate the algorithm considerably, so we have elected the expedient of simply shortening the i n t e r v a l s over which the s o l u t i o n i s to be found. This does not n e c e s s a r i l y e n t a i l a very large decrease i n e f f i c i e n c y using the transform process,* f o r i f the time step s i z e i s decreased by four times, then h a l f as many time points i n the transforms can be used f o r about the same accuracy, so the i n t e r v a l i s shortened by eight times. The a d d i t i o n a l c a l c u l a t i o n to get past numerically d i f f i c u l t points i n the s o l u t i o n i s ininimized i n t h i s way. I f the assumptions of continuous inputs and responses hold true, and i f the d e r i v a t i v e s of the nonlinear c h a r a c t e r i s t i c s are continuous 57 i n the i n t e r v a l , then decreasing the i n t e r v a l should decrease the v a r i a t i o n i n the functions involvedo This r e s u l t s i n the pseudo-Newton steps being closer to the true Newton steps. A l s o , one i m p l i c -a t i o n of proofs of convergence of Newton's method, such as that given 45 by Ostrowski^ 9 i s that, making reasonable assumptions (see Ostrowski), convergence I s assured i f the s t a r t i n g point i s s u f f i c i e n t l y close to the s o l u t i o n . Decreasing the i n t e r v a l w i l l r e s u l t i n the s t a r t i n g point being c l o s e r to the s o l u t i o n . This device has been s u f f i c i e n t to assure convergence i n the problems we have t r i e d . Use of the transform process permits one f u r t h e r benefit to be e x p l o i t e d i n the nonlinear s o l u t i o n . I f , i n the e a r l y I t e r a t i o n s , the higher frequency components i n the transforms are not c a l c u l a t e d , the responses Of the l i n e a r subnetwork can normally s t i l l be determined to 'Within --reasonable accuracy from the lower-frequency components. I f the nonlinear functions are only evaluated at some of the points i n the i n t e r v a l , say t = 0, 2At, 4At,... , or t = 0, 4At, 8At, ... , the r e s u l t i n g samples are u s u a l l y s u f f i c i e n t to determine the lower frequency components of the frequency domain spectrum to w i t h i n reasonable accuracy. Using t h i s procedure, convergence can be checked with a small amount of work compared to using a l l the frequency domain components and a l l the time points. I f convergence i s achieved, then the i t e r a t i o n s can be completed Using the f u l l transforms, but now s t a r t i n g with a good approximation to the s o l u t i o n . In a sense, we are i n t e r p o l a t i n g using Fourier s e r i e s to get a good approximation-to the s o l u t i o n . The o v e r a l l computation may be reduced considerably using t h i s device. 58 IV. IMPLEMENTATION OP ALGORITHM 4<>1 Problem Formulation The techniques previously discussed have been Incorporated i n a general program f o r t r a n s i e n t analysis of nonlinear networks,. We s h a l l discuss some of the techniques employed to achieve an e f f i c i e n t implementation^ A modified sparse tableau formulation, as proposed by Hachtel, 46 Brayton, and Gustavson , has been used f o r the network equations. The form of the tableau i s given i n (4ol),» Kirchhoff»s Voltage Law Kirchhoff'3 Current Law Linear Branch Constitutive Relations Nonlinear Branch Constitutive Relations Time Derivatives V(s) • (b) V (s) I(s) x ( 0 + ) (4.1) 60 where A Is the Incidence matrix I Is the u n i t matrix V I s the column vector of node voltages 1 Is the eolumn vector of Branch currents v I s the eolumn vector of branch voltages X Is the eolumn vector of djoiamie v a r i a b l e s ( c a p a c i t o r charges and Ihduetor f l u x - l i n k a g e s ) Zf 1$ Bj and E aire matrices d e f i n i n g the branch and d i f f e r e n t i a l r e l a t i o n s £- and f. are the vector functions Introduced i n the l a s t chapter £ n as the Independent v a r i a b l e s In (3.7) fhe eiample ef Elgo 4&1 w i l l serve to e l a r i f y the process* Sxainple of Tableau Formation Flgi. 4*1 Ihe tableau form fei- t h i s problem Is given By (4*2). CM • CM O I o o o o o -p r CM u_i O KJ 03 HI 03 03 CO TO HQ t=> > -H 'rj -H t> !> I> <3? 1 O 1 o 1 ° 1—1 1 — KS o o 1 1 ra o o 1 \- - - o o 1 ° o H I O 1 O 1 ° o o 1 ' o H 1 O o i •-I i _ o o 1 o -] o , o | o o o 1 H 1 I o 1 1-1 !—t ' O I 1 o o O | o o o 1 o _. — — L_ — . _ — — t - -o 1 | 1 ° 1 o 1 o H 1 H o 1 " J 1 The term -c i s the term a r i s i n g from the C matrix employed i n the nonlinear a n a l y s i s . There are several advantages to use of the tableau form. Various types of c o n t r o l l e d sources, i n c l u d i n g c u r r e n t - c o n t r o l l e d voltage and current sources, and vo l t a g e - c o n t r o l l e d voltage and current sources are a l l handled without d i f f i c u l t y . S i m i l a r l y , a broad c l a s s of ,. -nonlinear: elements may be handled without s p e c i a l provisions; we permit current and voltage c o n t r o l l e d r e sistances, current and f l u x c o n t r o l l e d inductors, and voltage and charge c o n t r o l l e d capacitances. A simpler scheme than that used i n the o r i g i n a l tableau a n a l y s i s was developed f o r s o l v i n g the system of l i n e a r equations given by ( 4 . 1 ) . As i n the o r i g i n a l work, the tableau i s analyzed at the s t a r t of the procedure, and Fortran subroutines containing the e x p l i c i t coding to do "the a c t u a l a n a l y s i s are w r i t t e n . The tableau i s considered to be i n two sections, as ind i c a t e d i n F i g . 4 .2 . I I A 0 --I 0 0 0 0 0 z I B _l_ 0 z n I n B n 0 E s i Breakup of Tableau F i g . 4 .2 Part I contains only constant, r e a l elements, whereas part I I contains some r e a l elements which change, and some complex elements which change. 63 The system of equations (4„l) i s solved once, using a conventional Gaussian e l i m i n a t i o n technique with a l l the normal indexing operations associated w i t h the sparse matrix s o l u t i o n "being performed.. The el i m i n a t i o n i s c a r r i e d out f i r s t i n the rows i n part I» Because the r e s u l t s of t h i s part of the e l i m i n a t i o n w i l l always be the same, the operations are not recorded; but the r e s u l t i n g matrix, having non-zero elements only on and above the p r i n c i p a l diagonal, i s stored* The corresponding operations f o r the input vector, however, w i l l r e s u l t i n new values f o r each new input, so a subroutine to carry Out these operations i s generated as the e l i m i n a t i o n i s performed. Because a l l the numbers m u l t i p l y i n g the input values w i l l be r e a l constants, a l l the statements i n the subroutine are w r i t t e n as m u l t i p l i c a t i o n s or d i v i s i o n s by r e a l constants, and additions and subtractions are Carried out as r e a l a r i t h m e t i c a l operations. The subroutine i s then used twice f o r the actua l a n a l y s i s , once on the r e a l parts of the input data, and once on the imaginary parts,, A f t e r the e l i m i n a t i o n i n part I i s completed, a subroutine f o r the back s u b s t i t u t i o n i n t h i s part, based on the reduced matrix, i s w r i t t e n , again c o n s i s t i n g only of statements f o r r e a l a r i t h -metical operations. In s o l v i n g the system of equations'described by the tableau, a large number of m u l t i p l i c a t i o n s and d i v i s i o n s w i l l Involve +1 or -1, and as such are a c t u a l l y unnecessary operations. In the o r i g i n a l work, those operations are kept out of the f i n a l code by using an index of v a r i a b i l i t y type. We have adopted the simpler procedure of simply performing a l l these unnecessary operations when the Gaussian e l i m i n a t i o n i s f i r s t c a r r i e d out, but, each time a Fortran statement i s about to be w r i t t e n , a l l m u l t i p l i c a t i o n s and d i v i s i o n s are checked to see i f they 64 involve +1 or -1. The unnecessary operations are removed at t h i s point. The e l i m i n a t i o n i s then c a r r i e d out i n se c t i o n I I . Because some of the elements i n part I I of the tableau are complex, the subroutines w r i t t e n f o r part I I use e n t i r e l y complex ar i t h m e t i c . Also, because some of the element'values change i n successive s o l u t i o n s , a l l the ar i t h m e t i c a l operations i n the e l i m i n a t i o n process are w r i t t e n out i n the corresponding subroutine, except f o r m u l t i p l i c a t i o n s and d i v i s i o n s by +1 or -1, as noted before. Because the elements c.. r e s u l t i n g from the nonlinear a n a l y s i s may have a r b i t r a r y values, p a r t i c u l a r l y a r b i t -r a r i l y small values, i t seems reasonable to avoid using these elements as p i v o t s , as accuracy problems could a r i s e 0 To achieve t h i s , the elements c.. are assigned a very small value f o r the f i r s t a n a l y s i s , so that the routines w i l l not normally s e l e c t those elements as pi v o t s f o r the e x p l i c i t l y ' c o d e d s o l u t i o n . S i m i l a r l y , because the'time scale f o r p r a c t i c a l c i r c u i t s i s u s u a l l y such that values of s considerably l a r g e r than u n i t y are encountered, the value of s chosen f o r the f i r s t a n alysis i s large to encourage the routines to pick elements i n v o l v i n g s as piv o t s . These techniques do not have the s o p h i s t i c a t i o n implied i n the process used i n the o r i g i n a l tableau s o l u t i o n , i n which s i x types of var i a b l e s are i d e n t i f i e d and dealt with i n d i f f e r e n t manners. But, i t i s f e l t that using the techniques described, the majority of the benefit i s gained by a simpler process. 4.2 Other Techniques Used Three fu r t h e r topics w i l l be considered: the method of s t a r t i n g the algorithm; methods of storage used; and the representation of nonlinear functions. For the f i r s t two, considerable non-numerical 65 processing i s done at the beginning to keep the actual a n a l y s i s e f f i c i e n t , , To s t a r t the algorithm or to begin a new time i n t e r v a l during program execution, values of v a r i a b l e s to be transformed and the d e r i v a t i v e s of those v a r i a b l e s are necessary. In order to determine these, we assume that the next part of the s o l u t i o n i s to be found f o r the I n t e r v a l t £ [ a , a+r/2 ] . Then the values f o r the dynamic v a r i a b l e s and the independent v a r i a b l e s c o n t r o l l i n g the nonlinear functions are found from t h e i r values at the end of the l a s t i n t e r v a l , or, f o r a = 0 , the values are supplied by the user or defaulted to zero. Two very short steps of f i r s t order i n t e g r a t i o n are then taken. The coding f o r the sparse tableau s o l u t i o n can be used d i r e c t l y f o r t h i s ; instead of (4.1) we have A 0 --I 0 V(a+h) 0 0 0 A* h 0 h 0 h • i(a+h) 0 f £(a+h) 0 Z n T n B n v(a+h) f [i(a+h),v(a+h),x(a+h)] 0 E x(a+h) E x ( a ) (4.3) A l l imaginary parts are set to zero and Ignored. Newton's method, with a c o n t r o l l e d step s i z e , i s used when nonlinear functions are present. Normally, except f o r the f i r s t i n t e r v a l , a s i n g l e i t e r a t i o n s u f f i c e s because of the small step s i z e used. Because the step s i z e i s chosen very small, the i n t e g r a t i o n errors are assumed to be n e g l i g i b l e . The f i r s t step i s intended to remove any inc o n s i s t e n c i e s i n s a t i s f y i n g the 66 network equations due to the errors incurred In using the transform process. The values obtained a f t e r the second step are used as i n i t i a l values of the network v a r i a b l e s f o r the next i n t e r v a l . The differences between the v a r i a b l e s a f t e r the f i r s t and second steps are used to estimate the time d e r i v a t i v e s of v a r i a b l e s which must undergo the transform process. The i n i t i a l values and d e r i v a t i v e s are used to cal c u l a t e the parameters of the transform process. I f the f i r s t step were not taken, the small Inconsistencies i n s a t i s f y i n g the network equations due to the err o r s i n the transform process i n the l a s t i n t e r v a l might swamp the small differences used i n estimating time d e r i v a t i v e s of network v a r i a b l e s . This method of determining the numerical transform parameters i s not the only possible route, but i t does make maximum use of e x i s t i n g coding, and adds l i t t l e to the o v e r a l l computation. Storage w i t h i n fhe program i s also handled by e x p l i c i t l y coded Fortran subroutines. To compute the transformed inputs to the l i n e a r subnetwork, several values of each input v a r i a b l e must be stored before the FFT can be used. S i m i l a r l y , as each network s o l u t i o n i s completed i n the frequency domain, the frequency domain samples of the v a r i a b l e s c o n t r o l l i n g nonlinear functions must be stored f o r subsequent processing using the FFT, i f nonlinear i t e r a t i o n s are i n progress. I f the network i s l i n e a r , or i f the nonlinear i t e r a t i o n s have converged, then samples of the transformed output v a r i a b l e s must be stored. For these purposes, storage must be al l o c a t e d , and an indexing scheme f o r s t o r i n g the va r i a b l e s i s necessary. These tasks are handled by e x p l i c i t l y coded Fortran subprograms which are generated at the beginning of the problem. A f t e r reading i n the network d e s c r i p t i o n cards, and s t o r i n g the information i n a disc 67 f i l e a f t e r pre-processing, the program generates the sparse tableau; but besides that, i t generates a main program which defines storage arrays f o r a l l purposes needed, with the minimum amount of storage required f o r each. This main program then c a l l s a "master' subroutine which controls the a n a l y s i s from then on, because the coding f o r the •master' subroutine i s f i x e d so need never be recompiled. As w e l l , subroutines are generated which extract and store the necessary v a r i a b l e a f t e r each network s o l u t i o n . The indexing of array v a r i a b l e s i s thereby done by constant Fortran i n d i c e s , rather than by v a r i a b l e i n d i c e s -computed by the program, saving computation during the actual a n a l y s i s . The method of representing nonlinear functions also deserves mention. They are described by Fortran subprograms, which must be supplied by the user. When the network d e s c r i p t i o n cards are read i n , the branch d e s c r i p t i o n cards describing nonlinear branches must each contain the name of _a user-supplied Fortran f u n c t i o n subprogramo The program, a f t e r reading these cards, generates an e x p l i c i t l y coded Fortran subroutine to c a l l the appropriate functions„ For the convenience of the user, the program has the c a p a b i l i t y of s t o r i n g a disc f i l e of model descriptions. Using t h i s f a c i l i t y , t r a n s i s t o r models, f o r example, may be entered once, and the object code f o r the appropriate nonlinear functions permanently appended to the analysis routines. Thereafter, i n describing a given c i r c u i t , any occurrence of a stored model may be described by a single card entry. Examples of r e s u l t s produced using t h i s program f o l l o w i n the next chapter. 68 V. EXAMPLES 5.1 Simple Linear Network The simple l i n e a r network of P i g . 5.1 was chosen to i l l u s t r a t e that the method i s free of eigenvalue r e s t r i c t i o n s . The network i s also designed so that the approximate behaviour may be estimated by i n s p e c t i o n . IO~6H. 2 i o o n 100 o o 1 ^ 1 0 0 u(t) Simple Linear Network P i g . 5.1 The approximate behaviour i s given by the f o l l o w i n g , where v_^(t) I s the voltage f o r t > 0 at node i . v J t ) - 100 (1 - e ~ 1 0 t) v 3 ( t ) = 1 - e v 4 ( t ) - 1 - e v 5 ( t ) * 1 - e -10 _ 6t -3 -10 J t - t <5.1) The problem was solved vising three d i f f e r e n t l e v e l s of r e l a t i v e accuracy: one part i n 100 (accuracy of 1%), one part i n 1000, and one part i n 10000. The node voltages are chosen as output v a r i a b l e s . For t h i s problem, one other feature of the program was used, that of s p e c i f y i n g 69 a t o l e r a b l e l e v e l of absolute truncation e r r o r . I f t h i s i s not done, the computing time i s very l a r g e , f o r the f o l l o w i n g reason: the value of ft chosen f o r each i n t e r v a l I s proportional to maximum a l l transformed v a r i a b l e s d e r i v a t i v e of v a r i a b l e at s t a r t of i n t e r v a l value of v a r i a b l e at s t a r t of i n t e r v a l Since a l l v a r i a b l e s not associated with the source begin w i t h zero value, the program would attempt to s t a r t the s o l u t i o n using very large values of ft , and consequently extremely short time steps. As a r e s u l t the s o l u t i o n would proceed very slowly. For t h i s problem, an a l t e r n a t i v e method of avoiding the s t a r t i n g d i f f i c u l t y would be to choose a small value f o r each i n i t i a l c ondition. For example, i f i n i t i a l conditions on the capacitors and inductors were chosen so that each node had an i n i t i a l voltage equal to 0.001, t h i s would e s s e n t i a l l y be the same as f i n d i n g the s o l u t i o n with i n i t i a l conditions of zero but with an i n i t i a l absolute e r r o r of 0.001 v o l t i n each node voltage. For the a n a l y s i s , the t o l e r a b l e absolute errors were chosen to be numerically equal to the r e l a t i v e accuracy i n each case t r i e d . For example, when using a r e l a t i v e accuracy of 0.001, absolute errors of 0.001 v o l t are permitted. Selected r e s u l t s are tabulated i n appendix I . A blank l i n e i s used to separate r e s u l t s from d i f f e r e n t i n t e r v a l s . The analysis was c a r r i e d out from time t a 0 to t > 10 sec. i n each case. The computing times used f o r the analyses are as follows (using the U n i v e r s i t y of B r i t i s h Columbia's IBM 360/67 computer) : 70 Accuracy O o O l 0.001 0.0001 C P U time(sec„) 1.5 4 . 7 3 2 . Computing Times f o r Example One Table I I These times do not include the time to set up the a n a l y s i s . In a longer example (se c t i o n 5 . 3 ) ? the various times f o r a l l steps of the s o l u t i o n (set-up, compilation, and analysis) are presented. Results are graphed i n P i g . 5 . 2 . No i n t e r p o l a t i o n was done between data points. As can be seen, the s o l u t i o n s with r e l a t i v e accuracies of 0.001 and 0.0001 are v i r t u a l l y i n d i s t i n g u i s h a b l e i n t h i s graph. I n the s o l u t i o n with r e l a t i v e accuracy of 0.01, the G-ibb's type errors are sometimes evident as small o s c i l l a t i o n s i n the s o l u t i o n . As can be seen from the r e s u l t s , the algorithm i s able to choose step s i z e s s o l e l y on the basis of the rate that the s o l u t i o n i s changing. The algorithm remains stable with large step s i z e s a f t e r the s o l u t i o n has s e t t l e d . 71 Node Voltage 3.00-1 • Time(sec*) F i g . 5.2 . Behaviour of Simple Linear Network 72 5 . 2 Saturating Inductor Example The example shown i n P i g . 5 « 3 was chosen to i l l u s t r a t e some c h a r a c t e r i s t i c s of the method when applied to a simple nonlinear network. The same s o l u t i o n parameters were chosen f o r t h i s s o l u t i o n as f o r the f i r s t example. These include three runs with r e l a t i v e accuracies of 0 . 0 1 , 0 . 0 0 1 , and 0 o 0 0 0 1 , and tole r a b l e absolute truncation errors equal to the r e l a t i v e accuracy i n each case. Also, the s o l u t i o n was ca r r i e d out from t = 0 to t > 1 0 sec. In each case. The computing times are given i n Table I I I . 1 Q vW" Saturating Inductor Example F i g . 5.3 Accuracy 0 . 0 1 0 . 0 0 1 0 . 0 0 0 1 CPU time(seco) 0 o 9 3 . 0 2 0 . Computing Times f o r Example Two Table m Some of the tabulated r e s u l t s are given i n appendix I I . The r e s u l t s are also graphed on the fo l l o w i n g page ( F i g . 5 . 4 ) . Again, 74 the r e s u l t s f o r the solutions with r e l a t i v e accuracies of 0.001 and 0.0001 are almost i n d i s t i n g u i s h a b l e i n t h i s graph. Some form of i n t e r p o l a t i o n or curve sketching would be desirable f o r the s o l u t i o n with r e l a t i v e accuracy of 0.01 . As i s evident i n F i g . 5.4, simply j o i n i n g the c a l c u l a t e d points i n the s o l u t i o n with s t r a i g h t l i n e segments r e s u l t s i n apparent slope d i s c o n t i n u i t i e s which are not c h a r a c t e r i s t i c of the s o l u t i o n . No d i f f i c u l t i e s were experienced, i n obtaining convergence. This i s seen from examining the r e s u l t s i n Appendix I I . A l l i n t e r v a l s contain the f u l l number of transform points i n each case, i n d i c a t i n g that no i n t e r v a l s were shortened because of convergence problems. 5.3 Class B Complementary Push-Pull A m p l i f i e r This example was chosen as an i l l u s t r a t i o n of the a n a l y s i s of a analogue t r a n s i s t o r network. I t i s l i n e a r i n the sense that i t . i s designed to provide l i n e a r a m p l i f i c a t i o n . But, i t also e x h i b i t s strong nonlinear behaviour by v i r t u e of the class B design. The c i r c u i t i s given i n F i g . 5.5, and the t r a n s i s t o r model i n F i g . 5.6 . 75 + 20 v. ( J ) V R r 2" 200 5 7 O R. s 5 n 4 -4v 0 6 R 2 >5000Q R l g 500 Q. *6 | s i n 2ir(10 3)t 5 *3 i o a -3 10 F. R 6 > 8 Q Class B A m p l i f i e r Figo 5»5 bb B O v w G O — ° ' k v i d = k 2 e + k 3 v d Transistor Model F i g . 5.6 76 Transi s t o r rbb k l k 2 k 3 k 4 a 1 0 0 25 i o - 6 2X10~ 7 2X10" 6 0.99 1 0 25 l O - 5 2X10~ 6 4X10" 7 0.98 1 0 -25 -IO" 5 2X10~ 6 4 X 1 0 - 7 0.98 T r a n s i s t o r Parameters Table IV o I-y r 25 V „ i . = 10"° e a + 2 x 10"' a d 6 Diode Model P i g . 5.7 The t r a n s i s t o r parameters are given i n Table IV, and the diode c h a r a c t e r i s t i c i n F i g . 5.7. The choice of parameters used was a r b i t r a r y , and does not necessarily r e f l e c t the behaviour of any e x i s t i n g t r a n s i s t o r s or diodes. The analysis was c a r r i e d out with a r e l a t i v e accuracy of 0.001 . To reduce the amount of computation f o r start-up, i n i t i a l conditions were s p e c i f i e d on the d i f f u s i o n capacitances w i t h i n the t r a n s i s t o r s such that they s t a r t w i t h i n t h e i r active regions. Also, an i n i t i a l c ondition of 10 v o l t s was put on capacitor to put the c i r c u i t close to i t s operating point. The s o l u t i o n was calcula t e d f o r two cycles of the e x c i t a t i o n , that i s , f o r two msec. The analysis took 185 sec. of computing time. For 77 a problem of t h i s s i z e , the set-up time becomes considerable. The program which processes the network d e s c r i p t i o n to generate the sparse tableau and then the e x p l i c i t Fortran coding f o r the network s o l u t i o n i s quite f a s t , t a k i n g about 10 sec. f o r t h i s problem. However, the Fortran compilation of the e x p l i c i t coding of the s o l u t i o n requires about 7 0 s e c , using the IBM Fortran Q compiler on the IBM 360/67. The r e s u l t s are graphed i n F i g . 5.8 and F i g . 5.9. The f i r s t graph shows the voltage developed across the 8^ load. The second shows the current from the 20 v o l t dc supply. The e f f e c t of the cl a s s B operation i s evident from t h i s f i g u r e . The current. I s s i n u s o i d a l while t r a n s i s t o r T^ i s on during the p o s i t i v e h a l f cycle of the output voltage. Then, f o r the negative h a l f cycle, when T^ i s cut o f f , almost a l l the current supplied i s that drawn by the d r i v e r stage T « I t w i l l be noticed that there i s a-dc o f f s e t - t o the output voltage i n F i g . 5.8. This occurs because the c i r c u i t has not yet s e t t l e d to steady-state operation a f t e r 2 msec, of operation. This was ascertained by- performing a f u r t h e r analysis to determine the quiescent operating point of the c i r c u i t . For t h i s , the analysis was repeated, but with the s i g n a l voltage, V^, set equal to zero. This a n a l y s i s was c a r r i e d out f o r 45 s e c , with a r e l a t i v e accuracy of 1$. The beginning and end of the tabulated r e s u l t s are i n Table IV ( f o l l o w i n g .Fig« 5.8 and 5.9). Several f a c t s about the c i r c u i t ' s operation are evident. The emmitters of T and T_ (node j ) are at approximately 11.5 v o l t s above 2 5 ground p o t e n t i a l , rather than the design value of 10 v o l t s . I t i s the r e s u l t i n g charging current through which causes the i n i t i a l o f f s e t i n the output voltage. Also, the supply current i s approximately 90 ma., of which about 40 ma. go through and therefore through the d r i v e r 78 F i g . 5.8 Output Voltage from Class B A m p l i f i e r 79 F i g . 5.9 Current Drawn by A m p l i f i e r 80 Analysis of Quiescent Point of -Amplifier Table V -votTAtje-AT NODE " z — -VOtTAGE-AT NODE — 3 "CURRENT" THROUGH BRANCH V2 -CURRfNT-THROUGH BRANCH ~ R 5 --V0L"TA1;£-ACROSS BRANCH -R6 --0 7 2 3 2 9 01 6~2 D - 0 6 0j9 22 0" 0 e46521524D»06 0,9210 0,697522860-06 0,9202 0,929830490-06 0 ,9J9S 0 S116213810-05 0j9189 0,139444570-05 0,9184 -OTl 62675300-OS—01917?" O e 16590610D-05 0,91?« " 9,-943 " 9,918 9,903 9,895 9,892 9,891 9,897 -».-077 657Cr»t)t-«0,7611D»01 «0.8033D«Ol PIO, 7869D»01 -0,81420-01 *0,7960D«0t -* 0". 8 1910 - 01" BO, 79880-01 0-,fl86SD<0r 0, 487BD<»Oi 0 ,4884D«>01 0 ,4889D-0 1 0,48900-01 0,4890Df01 "0 ,t|ft8-8 0 ^ 0 r -*-0 .57470*01" -.0 .8 1950-01 «0,96570*01 "0,1049 «0,1085 -0 , 1 088 •"tr. 10 65 0.4887D-01 -0,1032 0,447461800-04 0,87Sl3300D-04 -OTl 3028042O-0-3-0,1730475«D«03 0,215814660-03 0,2S656t78D-03 0,301348900-03 0(344116020^03 -0r447861800-04-,0«875533000«04 0,130320420-03 - 0,17308754O-03 0,8715 0,8561 ^ «r7 8S0T-0,8491 0 ,8478 0,8462 0,8480 0,8477 — « 7 8Tt5" 0,8561 0,0501 0,8490 " 10,80 l l . i l "l-tT-24-11,26 11,29 11 ,28 11,2"? 11.30 -rOTSv" H , U 11,24 IX,SO' s.0,1482 -0,1817 ovt 9-a-s-* 0 j 1 9 6 7 *0,1990 • <- 0,1978 -0,1977 .•0,1980 -0,1817 ^0,1948 - rO, 1967 ' 0.4402D-01 0,42370^01 0,41610-01 0,41470*01 0,41490-01 0,41460-01 0 f414lD«>0i -Ov4 4 0 20-rtt-0,42370*01 0,41740-01 0,4161D»01" 0 ,8004 1,106 — r , - 2 2 r -1 ,239 1 ,259 - 1,248-1,218 1,251 -O78OOS-1,106 1,221 - 1,239--0.-T751S279D-.0-3-0 , 177217810-03 0, 179282820-03 0,18134784D«03 0,t8341285D"03 0,185477870-03 -07187542880-03-0,18960790D-03 0,310552530-03 0,431449160-03 -0,552343790-03-0,673242420-03 0,794139050-»03 0.91503568D-03 0,103593230-02 0,115682890-02 - O r&Sirt-0 ,8501 0 ,8500 0,8499 " 0,8498 0,8498 -fl^8«97-0,8496 0,8483 0,8491 —0,6492-0,8496 0,8497 0,8500 0 ,8502 0,8505 11,24 11,24 11 ,24 11.25 11,25 "11,25 11,25 - ^ 0 T l 9 4t-PO.1945 H-0,1945 -0,1948 -0, 1948 -0,1950 — r 0 f l 9 5 0 -nO, 1952 0,1958 0,1939 -0,-1924--0, 1905 -0,1889 -0,1870 -0,1855 -0,1638 -0;4 17 2 D»01 ry2t-7~ 0,41720-01 0,4170D"0i 0,4l70Do01-0,41680i01 0,41680-01 -Oy t-l-6-70 "0-r-0,4167D«01 0.4152D-01 0,41510-01 -0-.-4 1490-0-1-0.4148D-01 0,41460^01 -0,41450-OJ-0.4143D-01 0,4 1420«01 1,218 1 ,220 1.221-1,222 1,223 -1 T225-1 ,226 1 ,232 1,217 -t-fZW-1 ,188 1 .175 - 1,160 1,147 1.132 -0T13544400O-O1-0,13442234D-01 0.13540067D-01 0,136379010-01 0,137357350*01 0 , 138335690-01 - 0 , 1 593 1 4030-0 1 0,140292370-01 -OT86-5O-0 ,8651 0 ,8652 0,6652 0,8653 0 ,8653 0,86,54" 0,8654 -11752— 11,52 11.53 11.53 11,53 11,53 - 1 1 , 5 3 — 11,53 -e0791580-OT" -0,91190-01 -0,90910-01 -0,90550-01 -0.9027D-01 -0 , 89920-0 1 *0 tR96aO-01 !t0 ,89320'-01 -0-,-« 0 6 20*01-0,10620-01 0,40620-01 0 , 4061D-01 0,40610-01 0 , 4061D-«01 0 ,40610-0 r 0,40600-01 -073340-0,3307 0,3276 0,3245 0,3214 0,3183 0,3153" 0,3122 81 stageo Prom t h i s observation, we conclude that the no-signal current through n;^  and (to avoid crossover d i s t o r t i o n ) i s approximately 50 ma. Part of the t a b u l a t i o n from which Pigs. 5.8 and 5°9 were prepared appears i n Appendix I I I . 5,4 Monostable M u l t i v i b r a t o r As an example of the a n a l y s i s of a nonlinear ( d i g i t a l ) type of c i r c u i t , a monostable m u l t i v i b r a t o r was analyzed. The c i r c u i t i s given i n P i g . 5.10 . 0 v i 10 V. R 3 > 10 4n R 4 $ 103Q 10 J fi G 2 10 'P. Monostable M u l t i v i b r a t o r C i r c u i t P i g . 5.10 The t r a n s i s t o r model used i s given i n F i g . 5.11 . 82 Q BO-v. Dl xbb 10ft i D2 v. D2 Dl 'D2 6 E 0.99 i 12 0.75 i Dl CD1 = °D2 = 1 0 ' Transistor Model f o r Example 4 F i g . 5.11 A simple c h a r a c t e r i s t i c was chosen f o r the diodes i n the t r a n s i s t o r model; i t i s i = 10~ 6 v^ + 10""4 ( v D - 0 „ 2 ) + + 10~ 2 ( v D - 0 . ? ) + + 1 0 " 1 ( v D - 0.4) + (5.1) where the subscript + has the f o l l o w i n g meaning: ( f ) + •= i o f < o f f > 0 (5.2) The t r i g g e r i n g source, v , i s described by F i g . 5.12. 83 Triggering Waveform P i g . 5.12 The t r i g g e r I s not s t a r t e d -until t = 0.3 msec, to allow the c i r c u i t time to s e t t l e a f t e r start-up. The an a l y s i s was c a r r i e d out using a r e l a t i v e accuracy c r i t e r i o n of 0.01, and required 57 sec. Tabulated r e s u l t s are given i n Appendix IV, The c o l l e c t o r voltages are graphed i n P i g s . 5.13, 5.14» and 5.15« In P i g . 5.13* the f u l l c i r c u i t operation can he seen, from start-up through one period of monostable ( u n t i l i t has f i n i s h e d s e t t l i n g ) . I n P i g . 5»14* the triggered switching of the monostable i s expanded to show the c o l l e c t o r waveforms more c l e a r l y . As can .tee seen, when <j> switches on, and I t s c o l l e c t o r voltage f a l l s , the negative going t r a n s i t i o n i s coupled through and Cv^ ( i n T^) to the c o l l e c t o r of T , causing i t to go negative f o r a short time* This i s an abnormal occurrence f o r a monostable c i r c u i t , and i s a consequence of using a s i m p l i f i e d t r a n s i s t o r model. F i n a l l y , i n P i g . 5«15* the time scale i s expanded about the i n s t a n t when the monostable r e v e r t s to i t s normal 84 s t a t e , at about time t = 1 msec. The t r a n s i s t o r model chosen i s somewhat l i m i t e d i n accuracy. However, as an examination of the tabulated r e s u l t s shows, i t i s quite s a t i s f a c t o r y f o r modelling a switching t r a n s i s t o r . Every mode of t r a n s i s t o r operation except reverse operation i s evident, from c u t - o f f through a c t i v e a m p l i f i c a t i o n to s a t u r a t i o n . Furthermore, although t h i s network does not i l l u s t r a t e reverse operation of the t r a n s i s t o r , the model i s also quite capable of e x h i b i t i n g that c h a r a c t e r i s t i c i f necessary. 85 lO.OO-i Time(msec.) F i g . 5.13 I'lonostable Performance 86 P i g . 5.14 . Triggered Switching of Konostable M u l t i v i b r a t o r 30.00 -4.00 Time (msec.) F i g . 5.15 • Recovery of Monostable M u l t i v i b r a t o r 88 VI CONCLUSIONS AND RECOMMENDATIONS FOR FURTHER WORK 6.1 Contributions of the Thesis Widespread acceptance of computer-aided network design methods requires programs and techniques -which are easy to use, s u f f i c i e n t l y general to handle many problems, and not too c o s t l y . In t h i s t h e s i s , we have discussed methods f o r l i n e a r network simulation which are quite general and reasonably e f f i c i e n t . Extensions to permit a n a l y s i s of nonlinear networks have been presented. And f i n a l l y , the implementation of these ideas i n an easily-used program has been in v e s t i g a t e d . Conclusions and some suggestions f o r continued work f o r these topics f o l l o w . 6.2 Transform Methods The Laplace transform has enjoyed wide success as an a n a l y t i c method of network s o l u t i o n . We have presented techniques which make i t competitive as a numerical method of s o l u t i o n . Compared to previously e x i s t i n g numerical transform methods, our approach seems more e f f i c i e n t and more general as to allowed inputs. Also, i t permits a simple extension to include nonlinear elements. The r e s u l t i n g algorithm has pr o v i s i o n f o r choosing step s i z e s and other parameters automatically. And, compared to i n t e g r a t i o n methods* the transform approach avoids the s t a b i l i t y problems which are often encountered. Various p o s s i b i l i t i e s e x i s t f o r f u r t h e r development of these techniques. Certain a r b i t r a r y choices have been made i n the deriv a t i o n s . For example, the sampling theorem was used as a .guide i n choosing the sample spacing i n the time domain^ but then an a r b i t r a r y choice f o r 89 the sampling i n s t a n t s was made,, The choice made i s probably not optimum f o r minimizing the frequency domain truncation e r r o r s . Furthermore, to c a l c u l a t e numerical approximations to transforms, i t was necessary to consider p e r i o d i c s i g n a l s with period T, but to disregard the r e s u l t s f o r part of the period. Again, the choice of which r e s u l t s to disregard ( i n our casej T/2 < t < T ) was a r b i t r a r y and not n e c e s s a r i l y optimum. I t would be worth i n v e s t i g a t i n g whether or not other choices might lead to improved r e s u l t s . Another p o s s i b i l i t y f o r f u r t h e r development involves a p p l i c a t i o n of the techniques discussed to other problems. As suggested e a r l i e r , i t should be easy to f i n d d e r i v a t i v e s of responses from frequency domain information. Also, i t should not be as d i f f i c u l t to incorporate p r o v i s i o n f o r d i s t r i b u t e d parameter elements, such as transmission l i n e s , as i t would be using time domain i n t e g r a t i o n . These p o s s i b i l i t i e are worth f u r t h e r consideration. 6.3 Analysis of Nonlinear Networks The transform techniques discussed In t h i s t h e s i s go beyond others i n existence i n that they are e a s i l y extended to analysis of nonlinear networks* B a s i c a l l y , the approach we have taken i s to separate the l i n e a r and nonlinear elements i n t o two separate subnetworks, and then search f o r a set of sign a l s such that the terminal behaviour of the two subnetworks i s compatible. A method of performing t h i s search i s presented. However, as suggested i n section 3.4, other i n v e s t i g a t o r s have developed methods f o r s i m i l a r problems, and perhaps one of those, or a combination of techniques, might be more e f f e c t i v e than the method we have used. This aspect deserves f u r t h e r i n v e s t i g a t i o n . 9Q 6 o 4 Implementation The techniques discussed were implemented i n a general purpose network simulation program. The sparse tableau form of network equations was adopted, and e x p l i c i t l y coded Fortran subroutines f o r the network s o l u t i o n were used. In general, these are h i g h l y e f f i c i e n t methods of performing the necessary algebra. However, the sparse tableau was o r i g i n a l l y intended f o r use with a numerical i n t e g r a t i o n scheme, and i t i s f e l t that other forms of s o l u t i o n might more f u l l y u t i l i z e the p o t e n t i a l of the transform approach,, For example, every s o l u t i o n of the equations described by the tableau r e s u l t s i n c a l c u l a t i o n of a l l branch currents and voltages, and a l l node voltages. Using transform methods., i t should be possible to work with t r a n s f e r matrices, and j u s t deal with the source si g n a l s and the desired response signalst, This i s one more item suggested f o r "further investigation,, Taken together, the methods presented i n t h i s t h e s i s are a c o n t r i -bution to the a v a i l a b l e algorithms f o r computer-aided network a n a l y s i s . Appendix I 91 Results from Analysis of Simple Linear C i r c u i t In the f o l l o w i n g pages, some of 'the tabulated r e s u l t s from the analysis of the network described i n section 5.1 are presented. In Table A l g the f i r s t page of the analysis with r e l a t i v e accuracy of 0.01 i s given. Then, In Tables A2 and A3, the same i s done f o r the analyses with r e l a t i v e accuracies of 0.001 and 0.0001 „ 92 Table A l . Analysis of Simple Linear Network with Relative Accuracy of 0.01 -vntTAGr-AT NODE 2 "VOLT AG fr AT NODE 3 -vorrAGt-AT NODE ' - — 4 -VOtTA-GT-AT NODE S -0T&5tr8878-Euwrt5 OT97S6" 0 , 170 17 1560-09 1 ,836 0,255254350-09 2,673 0,340337130-09 3 ,507" 0,425419910-09 4,331 0,510502690-09 5, 153 -'OyS95S8548u>*"09 5 ,'9 0,680668260-09 6,775 0,26050-05 0,151ID^OS 0,71490-05 O s1047DP04 0,14520^04 -OTtv24D*>a-4-0,24660-04 "0T5 3 5-2 D^T2 0 T 2 0 4 21" 0 , 700 1D->12 0 f 9978D<=I2 -0 , 14890-1 1-0,22340-11 0 , 32920-11 -0-,'f7-23D<'li"-0,65850-11 0,25610-21 0.3272D-21 -Q, 4 314 D-21" 0,58770^21 0,82030<-2i -0Tti5©Dp-2tr-0 , 16360*20 0,129580700-08 12,84 0,190974580-08 18,16 -0-.-252 3 6 8 46 D^O-S tV^iT 0,313762340-08 27,72 0.37515622D-08 32,07 "0,436550090-08 36,19" 0.49794397D-08 40,01 0,559337850-08 43,70 0,12770*03 0.2246DH-03 -0-,-35 08P~03-0,50760^03 0 ,6908D«-03 -0-,9008Dt*03 0,11340-02 0.139JD-02 0,22610-0 t3335Ds -0T50 8-3O-0,77090-0, 1 1370* O, 1625D" 0,22480' 0,3023D-09 0,59430-18 09 0,76330-18 09 0Tf01-8Uw-i7-09 08 08-08 08 0 , 1406D*-17 0, 1 9850-17 -0;-2827D«t-7-0,4Q09D»17 0,56190-17 0 , 121 060J50-07 -0,136174520^0? 0,'25128089D"07 0.31640325D-07 -0-T3 8-1-517 620^0?" 0,446631990*07 0,511746350-07 0,576060720-07 70 e86 85, 13 92,65 96,29 —98T4-9-99,41 99,93 10 0,5 0,61330-02 0,11320-0 1 0 , 17020-01 0,23120^01 -0-,2929Df 01-0,35580-01 0,41610*01 O e48UD*-0i 0 ,80670^-07 0 , 1 3 7 6 0 - 0 6 0,22940-06 0.3604D-06 -0y5-3 08O--0-6-0.74230-06 0.9941D-06 0, 12870^05-0,2247D»S4 0,294lD"i4-0,41150*14 0,6014D»14 -07 8-8930-1-4-0 , 130205-13 0 81865D"13 0,26050-13--0-,599252C80-'Oft—"lOOy? 0r45i>-0-0,114057810-05 99,89 0,6805 0 , 16819041D-05 1 00 , 1 0 ,8130 0.22232301D-05 -99,92 0 ,8863 0,27645561D-»05 100 , 1 0,9318 0.33058G21D-05 99,94 0,9557 0,-3 84720810" 0-5 10 Orl 0y9 6 93-0,438853410*05 i 00 e 0 0,9819 0,22400-03--0,5418D"03 0,94480-03 0.1410D-02 0, 18990-02 0,24130-02 -0-,-29 30Di»02-0,3459D°02 0.35291098D-04 100,0 0,9843 0,36860-01 0,66145662D*04 100, 0 0 ,9791 0 P66270»Oi -^7970002260-04 rOOyO 0v9835 0F*?4 200-'O1-0,127854790-03 100,0 -0,9802 0,1215 0,158709350«03 100,0 0,9836 0,1478 0, 189563920-03 10 0, 0 0 ,961 1 ' 0,1735 0.22041848D-03 100,0 0,9834 0,1982 0,251273050*03 100,0 0,9828 0,2226 0 , 137 1 1660D"02 0,249101100-02 0,361085600-02 0,473070100-02 100,0 100,0 100,0 100,0 0,9876 0 ,9893 0 ,9899 0 ,9900 -0T3-85-20P-U9-0,5872D"09 0,98610-09 0, 1622D-08 0,2516D^08 0 83683D«*08 -0,5129On08-0,66580^08 0, 17630P-05 0,3379D"05 -0>5843D«05-0,91860-05 0, 1333D-04 0 , 18300-04 0,2403D»0« 0,3053D"04 0,7398 0,9031 0,9662 0,9748 0,78620*03 0,17470-02 0,27840-02 0,38840-02 93 Table A2. Analysis of Simple Linear Network with. Relative Accuracy of 0.001 -yOL'TAGE' AT NODE 2 • V OUTAGE" AT NODE • 3 VOLTAGE" AT NODE a -VOL-TAGE-AT NODE 3 -07851-4278 2 0-*t1T 0,170225560*10 0,255308350-10 0,340391130-10 0,42547391D*10 0,510556690-10 -07595 639 4 8 Dw tO" 0,680722260*10 0,765805040-10 "0,850887820*10 0,935970610-10 0,102105340*09 -OTftOfrfS 6 2 09-0,119121900-09 0,i2763Q17D«09 ' 0, 136138450«09--ort-rr7-T35 2o-«-o-9-0,159402590-09 0,171031660*09 0, 18266073D"09-0, 194289790*09 0? 2059188605-09 -0-,-?i7.547930"0&-0,229i7700Ds>09 0,240&0607Df09 0,25243514!)-<09 0,264064200^09 0,275693270-09 -0-.287 322340-09-0,298951410^09 0,310580480^09 0,32220955D«09 -01-86-9113^-0, 1733 0,2578 0 e3429 0,4277 0,5124 -0T5"^ 7-2 0,6817 0,7663 0 ,8508 0,9352 i 8 020 —1TV0 4 1 , 188 1,272 "1,357 -1T4-6-5--1,582 1,696 -t~,81-0-1,924 2,039 -2yr53-2 S266 2,380 -2,494" 2,607 2,721 -2T8-33-2,947 3,059 -3,173--0l-34999545D»09-0,377775350-09 0,4055SS25D<»09 0,433335150«09 0,461115050-09 0,«88e9495D°'09 -0-.-51 6674860-09-0,544454760-09 0,572234660«09 0,600014560-09 0,627794460*09 0,655574360*09 -Oy6833S4260-09-0,711134170-09 0,738914070-09 0,766693970*09 -3,436-3,708 3 ?974 4,241-4,507 4,772 -SyO-3 7-5,300 5,564 5,826 6,087 6,348 -6-J.60 8-6,868 7,125 7,386 "07? 5120-0 8 0,18720-07 0 ,3696D"c07 0,62600-07 0,95300-07 0, 1 3540-06 -0718 250*06-0,23700-06 0,29850*06 -0 ,36740-06 0,4433D*>06 0,52650-06 -076168D«06-0.7143D-06 0,81890*06 "0,93070-06 -J8-240--0 , 4 8 8 8 0 » 0,72050" 0,11390-0,1805D-0,2782D» in-4f29D---0,59080-0,81810*-0,11010-0,14450-0,18570-0v23430-^-0,29090-0,35610* 0~, 4304D« 15 07263-7D-<r23-15 0.3000D-25 15 0,35Q4D*2S 14"-0*4 2790*25" 14 0,55120*25 14 0,74390-25 14 07fO-35D»2*-14 0,14590-24 14 0,20540^24 13—0.2867D-24-13 0,39450-24 13 0,53450-24 13 0-t-7126D>"24-13 0,93550-24 13 0,12100-23 t~3—0, 15440023--07!tft 5Dw05 0-.-5 697 D-l"* 0„ 1285D-05 0, 14740*05 -0, 1679D"05-0,18960*05 0,21270-05 -V72-3 70OT051-Q,2628D»05 0,28980*05 -0.31810-05 0,34780*05 0 t3788D"05 -0 ,-4 11-1 D«-05-0,44480-05 0,47960*05 0,51600-05 .,71060 0,87020 -0.1054D 0 S12620 0 . 10960-- O 7 I - 7 5/O 0,20460 0,23690 - 0 82722D 0,31090 0,35320 -0,-39910 0,44890 0,5026D 0,56060 • 13 13 12 12 12 -rt-2-»12 »12 -12" *12 *12 M-2-"12 -12 -12 - 0 7 2 4 T 0 D P 2 3 -0,3 I 630-23 0.4 07 6D-23 0 , 5 1 9 6 0 ^ 2 3 -0 , 6 5 3 8 0 ^ 2 3 0,6.1 4 1 0 - 2 3 -0-,-tOt) 3D---2-2-0, 1 2 2 4 0 * 2 2 0, 14000-22 0,1776D*22-0,21150*22 0 , 2 5 0 1 0 - 2 2 -0T2938D^2"2-0 , 3 4 3 1 D * 2 2 0 , 3 9 8 3 0 8 2 2 - 0 , 4 6 0 2 0 - 2 2 " -0r6-H 90-05— 0,71260-05 0,81890^05 0.9333D-05 0,i055D*"04 0,11840-04 -0-,-i-i20D<*0 4— 0,1464D-04 0, 16150»04 0,17730-04 0 , 1938D-04 0,21il0*0« -0-,-22910-Ofl— 0,24790-04 0.2672D-04 0,26750-04 -^i-f-4 3-70-' 0,93030-0,1 142D-0 8 13860-0,16620-0 , 19730--0y23200-0,2707D-0 ,31340'' 0,36050-0,41200-ft,,46830-» -0-,-52940-0,59570-0,66720-0,74 440" p-»-t2 0 12 11 11 11 11 1-1—0 11 0 n 11 11 11 l i -i i 0 11 0 11 0 7-7-3-09-0^ 22-.9661D-22 ,12520-21 ,16030*21 ,2025D«21 ,2530D"21 ,31250-21-,38230-21 ,4633D-21 ,5S69D-21 ,66410-21 , 78630*21 ,92470-21-,10810-20 ,12560-20 .1452D-20 94 Table A3.- Analysis of Simple Linear Network with R e l a t i v e Accuracy of 0.0001 -TIME- VOLTAGE-AT NODE 2 '— —"VOLTAGE" AT NODE - 3 -V0tT*6C-AT NODE a "VOtTA SC-AT NODE S if -OT856""&278^0i--r2-0 , 170765560^11 0,255848350-11 "0,340931 130*1 1 0,426013910-11 0,511096690-11 -073 9 677 948 0-11-0,681262260-1 1 0 976634504Ds-l 1 0,851427820-11 0,936510610-11 0,102159340BJO -OTtt" 06 67*20 «rj"0-0,119175900-10 0,127604170-10 0,136192450-10 0,14470073D"10 0,i5320901D«iO "0 ,"85050-02 0 0,17150-01 0 0,25600-01 0 0,34130-01 0 0.42630-01 0 0,51140-0 1 0 -075 964 0-"07 0 0,68140*01 0,76650-01 0,85150-01 0,93650-01 0,1021 - o y t - r o ? 0,1191 0,1277 -0,1361 0„1446 . 0,1531 "Oi 1"6 17 17 29 0^ 1-0 0 r l 6 t 6 0,170225560-10 0,1701 0,178733840-10 0.1786 . ••0, i-872,421£O^t0 0,187-1 0J19575O40D-10 0,1956 0,204258680-10 0,2041 -vr2t-27"6 696 O^t 0 0T2 1-2* 0,221275230^10 0,2211 0,229783510-10 0,2296 -0,23829179Ds10--0 ,2381 0,246800070^10 0,2465 0,255308350*10 0,2551 -0726 381-6630" It) 072655 0,272324900*10 0 s 2 7 2 l 0,280916280-10 0,289447650^10 -0729797902DM-0-0 ,3065.10390-10 0t3150ft-i77D»<10 0,323573140-10 0,33210451D-10 0,34063588D"10 -07349167260^10-0,35769063D-10 0,366230000-10 0,374761370-10 0.38329275D-10 0,391824120-10 -0,-400355490^10-0,408886870-10 0,417<U824D-10 0,425949610-10 0,434480900-10 0 .4430 1236DM0 0 ,2792 0 ,2879 - 0T2 9-6-3— 0 ,3049 0,3134 0.3219 0 ,3304 0,3389 -O73474-0,3559 0,3644 0,3729 0,3814 0,3899 -0 ,398 4-0,4069 0,4154 0,4239 0 e«324 0,4409 ;51"410«10-,16160-09 e5428D-09 .5976D-09 ,92350-09 ,13230-08 717940^08" ,2338D»<08 ,29530>»08 ,36420-08 ,44020-08 ,52360*06 •,614 0Df«0fi-.7119D-08 ,81680^08 ,92910-08 ,1048D-07 , U 7 5 D « 0 7 Tf309D'^7 ,1450D-07 ,1599D-?07 8-t-7!?it0-0.7 ,19170-07 ,20870-07 T2 264 D 5-0 7 ,24490-07 ,26400-07 ,28400*07 ,30460^07 ,32590-07 7-3 4 8-0 D-0-7 ,37070*07 -0729520-18-0,38070-18 0,59030-18 0,98520-18" 0,16270-17 0,25780-17 "0-,38 990-17-0,56510-17 0,76970-17 0,10700016-0 siailD«16 0 , 1 8 2 1 D • 1 6 -0 , 23 0 4R*r6-0,28680-16 0,35180-16 -0 9426 0D-16-0.5101D-16 0,60460-16 ,J71030-1 6 0,82760-16 0 ,95720-16 -..OJ. •100D-15-0,1256D«15 0,14260-15 —O7I-61-1-0^ -1-5" 0 , 1812D-15 0,2028D"15 — 0 ,22610-! 15 0,25120-15 0,27800-15 — 0,-30660-15 0,33720-15 —0739430^29— 0,4224D-29 0,46260^29 —0",52810«'29-0,63730-29 0,81300*29 —07f086"D"«28-0,14890-28 0,20620-23 - 0 ,28490-28 -0 „39000«28 0,52700-28 —0770-200^28-0 ,92 140"20 0,11920-27 — 0,-15230^27" 0 , 19200-27 0,23940-27 —07-295-20^2-7-0,36060-27 0,43640-27 --0,-52380-2-7-0,62390*27 0,73790*27 —Oy8670D-^27-0,10130-26 0,1176D«26 --0,13580«26-0,15610-26 0,17860-26 —0y2035Ow-£6-0,23060-26 0,39430-07 0,41890-07 -<ITW-37O^0-7-0,46940-07 0 ,49570-07 0.5229D-07 0 ,55070^07 0.5792D-07 -0-,-6085D-07-0,63850-07 0,66920-07 0 ,70070-0-7 0 .7329D-07 0,76580-07 - 0779940^07-0,83370-07 0,8688D*07 0,90460-07 0 , 9 411D-0 7 0,97840-07 0,37050-15 0,40570-15 -0 ,-44230-1-5-0,46140-15 0,52250-15 0,5660D«>15 0,61170-15 0,66000-15 -0,-7 1060-1 5-0 .7638D-15 0 ,81960-15 0.8780D-15 0,93920-15 0, 10030*14 -0,10700-14-0 , 1 1400-1 4 0. 12120-14 0, 12880-14 0.1367D-14 0, 14480-14 0,26270-26 0,29630*26 -0733230-26-0,3718D"26 0,41450-26 0,46100^26 0,5il2D*26 0 ,56540-26 - 0 ,6239D«<26-0.6867D-26 0,75430-26 0,8267D-26 0 ,90420*26 0 .9870D-26 -OtI0750*25-0,11700-25 0.1270D-25 0,13770-25 0,14900"25 0,16100-25 95 Appendix I I Results from Analysis of Saturating Inductor Network Some of the tabulated r e s u l t s f o r the problem given i n section 5.2 fo l l o w . In Table A4, the f u l l tabulated r e s u l t s f o r the analysis w i t h r e l a t i v e accuracy of 0.01 are given. For comparison. Tables A5 and A6 give the f i r s t pages f o r the analyses with r e l a t i v e accuracies of 0.001 and 0.0001 r e s p e c t i v e l y . Table A4 Analysis of Saturating Inductor Network with Relative Accuracy of 0.01 r T I M E VOLTAGE CURRENT AT THROUGH NODE BRANCH 2 R l > 0,85085197 i,7016799 -2,5525079 3,4033359 4,2541638 _S-t4JX49-9-L-a. 5,9558198 6,8066477 0,5619 0,2025 •»0. 1696--. -.0 ,5581 -0 ,5879 0,3414 0,5650 0,1900 0 ,7889 -0,7252 0 5 2 9 9 3 -0,3089 -0,6630 «0,6510D-01 7,7890528 £^7X1.433 a 9,7538148 10,736196 11,7185 7 7 12,700958 13.663339 14,665 720. 0,5282 0,4697 -0*217-7 u.,8 255--0.5166 0, 193"4 -0,56 5 4 -0,4 0 12 . 0, 1 157 . ... .-0,065b 0,4 39 4 • -0,305 2 0.6282 0,270b -0, 1233D?QJ 0, 876/4. 96 T T j ^ V O L T A G E CURRENT 1 4 A T THROUGH N 0 O F. B R A N C H 2 RJ 0,851067970-01 0,17018959 -0,25527239 0,34035519 0,42543796 -O^Si-O-Sitt-7-e 0,59560358 0,68068637 -0 ,765 76.9 t-7-0,85085197 0,93593477 t-^ i}.24.44-7-6-1,106100U 1,1911832 1 ,2762660-1,3613487 1,4950002 1 ,6286276 .1,7622551-.1 ,8958825 2,0295100 2,29676(19 2,4303925 — 2,5640197-2,6976472 2,8312746 3,0985295 3,2321570 3,3657844. 3,4994119 3,5902724 3,6811089 3.77.19454-3,0627819 3,9536184 4 ,J3-4i4..4 5-4-9— 4,1352914 4,2261279 4,3169644-4,4078009 4,4986374 /L, 589 413.9-4,6803104 4,7711469 4,8619834. 4,9528199 0 B8286D*»0 1 0,1565 .-0,2238 0,2831 0,3356 _0_r37-3.S» 0,4171 0,4454 _0,«662-0,4774 0 e«812 _0-,il /.46 0 ,4620 0 ,4381 ...0,4133 -.-0,3698 0,3098 0,2389 —0,.l671 0 , 9 6 1 0 0 - 0 1 0 , 2 7 3 3 0 * 0 1 -0, 1002 -0,1585 ..-0,2154 -0,2685 • »0,3226 _-JUJS-Xia-•.0,4252 -0,4692 -TO,524 0 -••0,5503 -0,5752 -0,5985 _-0.,6189. -0,6324 -0,6396 -0,6286 -0,6090 —.-0 ,5808 . -0,5428 -0 ,4981 -IU4 4 5: -0 .3902 -0 , 3288 —-0,2721.— -0.2042 0,2 1430-02 0 , 12820-01 _0 ,2873.aOi 0 ,5076D*01 0,77150^01 _o_^inaJz 0, 1439 0,1839 . 0,2268 0,2744 0 ,3240 -0 ,37 8.1 Table A5 Analysis of Saturating 0,4320 0,4907 0,5436-0,6083 Inductor Network using .Relative-Accuracy of O.OOl 0,6873 0,7594 ..0,814 0-0,8515 0 ,8693 0 ,8481 0,8113 .0,7614. 0,6980 0,6279 _0-,-x4_-6-0,4683 0 ,3788 -0,3017-0 ,2001 0,1414 0,84820-01 0.2942D-01--0 ,27930-01 »0,8610D"01 ___4-lA--9 -0,2094 -0,2751 -0,3420 « 0 , H 1 2 «0 ,4792" ___,54 7_1 -0,6092 -0,6694 .-0.7168 — ->0 ,7670 97 TIME VOLTAGE CURRENT AT THROUGH NODE BRANCH 2 Rl 0,853227970-02 0, 170405590-01 -fl,255488390*01-0>340571190-01 0,425653980-OS -0-^ J-07-3-6X8D»JXS-0,595819580-01 0 ,680902370-0 i 0,765985t7D?<01-0,85106797D»0l 0,936150770-01 1-024-23-3---0 0,11063164 0,11913992 0,12764820 0,13615647 0,14466475 _0-J-534_2341i_ 0,16168131 0,17018959 0,17869787 0, i'8720615 0,19571443 XL,-2.Q_4222__L 0,21273099 0,22123927 0,22974755 0,23825563 0,24676411 Q.,.?55 2 7-___-0,26378067 0,27228895 0,28114671 0-,285-9-8Jl4fl-0,29881424 0,30764800 0 ,31648177. 0,32531553 0,33414929 0-..3iL2_La3L_ii-0,35181682 0,36065058 0.36948434-0,3783181 1 0,38715187 0,3.9598563-0,40481940 0.41365316 ..0,42248692. 0,43132069 0,44015445 . 0,4 48__.a_-21-0,85090-02 0,16900-01 -0„2523D»01-0,3349D«01 0,41670-01 _0-,-4S-7840-"-0.l-0,5782D-0i 0,65780-01 -0.73680^01-0,81490-01 0,89240-01 _D-,-9-6-9.20jEil4-0,1045 0,1121 0,1195 0,1269 0,1342 __l,-t-a-t-j 0,23280-04 0 , 1 3540-03 0,31 1 9 D- 0 3 0,56410-03 0,8Bi7D-'03 -^-^4-2-30^02-0 , 17290-02 0,22570-02 -0,28.480^02 -0,3510002 0,42340-02 ..0., 5 0 2-9O 02 0,58830*02 0,68080-02 0,77900-02 0 ,88420-02 0 ,99490-02 JX^ JL4JL3fl___U_ 0„!486 0 ,1557 0,1627 0,1697 0,1766 _-U_L8-3_-_ 0. 1901 0,1968 .0,2034 0 ,2099 0,2163 -JU--22Z-0,1236D»Oi 0,1366D«Ql 0.1501D-01 0, 16430-*01 0,17890-01 __X)_-_1_ 0,2290 0,2351 0,21020-01 0 ,22680-0 1 0,24370-01 0,26150«01 0,27930-01 __U-L9___60_Q-L-0,31700-01 0 ,3386D«>Q1 0,2415 _jQ-.-2_4.7-L. 0,35950^01 ___,3iii--n__ai_ 0,2540 0,2602 _0.,2662. 0,2722 0,2781 _-L,2.8.4-l-0 ,2897 0 ,2954 .0.3010. 0 ,3065 0,3119 _J_,3173-0 ,3225 0,3277 . 0.3328 0,3378 0,3428 _JU3 4.7.6_ 0,4037D«01 0,42660-01 .0,44990-01 0,47390-01 0,4983D-01 __U_523..4D-_Jli_ 0 ,54690»0 1 0 ,57500-01 .0,60160-01 0,62880^0 f 0,6564O->01 -.0 , 6846D-Q 1_ 0 , 7 1 330-0 1 0 , 74240-01 0 ,77210-0 1. 0,80230-01 0,83290-01 JQ..8641D-Q1-Table A6 Analysis of riaturafdng Inductor Network with R e l a t i v e Accuracy of O.OOOl 98 Appendix I I I Results from A n a l y s i s of Class B A m p l i f i e r Table A 7 gives the f i r s t page of the tabulated r e s u l t s from the analysis of the a m p l i f i e r described i n sect i o n 5 » 3 . Appendix IV Results from Analysis of Monostable M u l t i v i b r a t o r Table A 8 gives the f i r s t page from the a n a l y s i s of the monostable m u l t i v i b r a t o r c i r c u i t described i n section 5 . 4 . 99 Table A7. Anal y s i s of Class B A m p l i f i e r -vorntGt-AT NODE " 2 -VOtTAGt-AT NODE - - 3 -CURRENT-THROUGH BRANCH V2 -CURREAT-THROUGH BRANCH Rb • -— -votTAGE-ACR033 BRANCH R6 '0,23 8 30762 D «0 7—0,9 231— 9", 9 T4~ 0,470615240-07 0,9230 9,970 0.702922B6D-07 0,9229 9,966 0,93523 0 090-07 — 0 ,9228 9,962 0,11675361D»06 0B9227 9,958 0,139984570^06 0,9226 9,955 -07l6"32r5"3*rD->*Tr6—0793-25 97-952-0,186446100^06 0,9224 9,948 -w0-,SS0^0-01 0,4856D»01"—»0y26?7D-01 -0,57760-01 BO,60010^01 »0, 61790-01 -0.6346D-01 -0,64770-01 -WOT66T2D--0T-*0,67060-01 0t4857D«»Oi 0,48580-01 0 r«859D-01 0.4860D-01 0,4862O»0t -err4 66-3D^Dt" 0 , 4 8 6 4 D » 0 1 0,30380-01 -0.3426D-01 -0,38020*01 -0,41640-01 -0. 45120-01 --«",-4 8q?D«01 0,51690-01 0.81U14456D-06 0,9204 9,899 -Qt7900D»01 0.4886D-01 -0,1006 0,144124300-05 0,9192 9,892 -0,79350-01 0 ,48900-01 -0 , 1079 "0T2 0 6"8-3 415 O--0 5 07<? t62 9,"9 01 : -0,8093 D»01 07*885 D-0l""-0-.90 960-01-0,269543990-05 0.9175 9,918 -0,80770-01 0,48760-01 -0,81650-01 0,3322538UD-05 0,9166 9,936 -0,82260-01 0.4866D-01 -0,64040*01 0 ,394963690-05 "0,9162 9,956 -0,82190-01' 0. 48550-501 -0,13740-01 0»45767353D-05 0,9155 9,975 -0,83830-01 0,4846D*0S -0,25280-01 0,520383380-05 0,9150 9,995 -0,83550-01 0,48350*01 -0,51910-02 0,830802790-05 0,9123 0,114116220*04 0,9101" 0, 145152160-04 0 ,9083 0,176188100*04 0,9070 "0,-207 22 4 Crtl 0 «• 04 ©79060-0,23825996D-04 0,9054 0.26929S920-04 0,9051 0,30 0 33166D-04—-07-9051-10,08 10,16 10,23 10,30 -1-07-3-5-10,40 10,44 -10,47--0,88390-01 -0,92600-01 -0,97130-01 -0,1013 -»07-f05e -0,1093 -0,1129 -0,1158 0,478BD-01 0.4745D-01 0,«707D»01 0,4674D»01 "01-46-450 -01 0,46190-01 0,4597D-0i . A iir -? ri r\ A * 0,83520*01 0,1639 0,2336 0,2956 -OT^SO-O— 0,3972 0,4382 0-,4579D-01—0,-4728--O7-5 6 4T4 9-360^ 0 #-0,467966860«04 0,551 184360*04 0,634401860-04 0,717619360-04 0 ,800836860-04 -07884054350-04-0,967271850-04 -OT-9063-0,9091 0,9151 0,9182 0,9239 0,9302 -0T9 368-0,9438 -STOTS^ -0,1-22-3 0,106379690-03 0,ll991220D-03 -Orl-3l«-44-70O-'03-0,142977210-03 0 , 15450972D*03 0,16604222D«03 0 , 17757473D«03 0 , 189107240*03 0,9533 0,9633 —07 972-9-0,9822 0,9910 0,9992 1 ,007 1,014 0,201364510*03 1,020 0,213501790*03 1,025 0,225639070-03 1,030 0,237776350-03 1,033 10,58 10,59 10,58-10,55 10,51 -t-0-74-6-10,40 10,31 10,21 -10i12-•10,01 9,916 9,819 9,729 9,642 9,564 9 ,488 9,424 9,370 -0,1256 *0,1266 -0,1253 -0,1228 -0,1169 ---0x1-1-4-4--0,1094 -0-7*5-41-0^ 01-0,45220-01 0,45170-01 0.4524D-0J 0,45380*01 0.4560D-01 -0745-8-7D--01-0,46190-01 -0T5W-2*-0,5769 0,5864 0,5735 0,5462 0,5040 -0T454t-0,3941 -0,1027 -0.9565D-01 --0-,69770-Ot-• 0, 847 10-01 fO,807?D*-01 -0.77620-01 s-0 , 75300-0 1 -0,73460-01 0,4664D»0t 0,47160-01 77 10-01 0.4825D-01 0,3090 0,2084 -0-.-1-097-0 ,75270-02 0f4877D"01 -0,90390-01 0,49290-01 -0 , 1873 0,«9760»01 -0,2770 0,50210-01 -0,3631 -0,72140-01 -»0,7108Di»01 »0,7032D-01 -0 ,69780«*01 0,50620-01 -0,4402 0,51020-01 -0,5158 0,51350-01 -0,5787 0,51630-01 -0,6326 Table A8 Analysis of Monostable .Tultivibrator 100 TIME VOLTAGE AT NODE 4 VOLTAGE AT NODE -5 VOLTAGE AT NODE 6 VCLTAGE AT NODE 7 0.966935060-03 3, 1631 0, 1680 0.7*150-0! 0, 2325 0.193357010-07 J,2966 0.2961 0.1284 0.4602 _0.2°00'705'0-07 0-4222 0.4216 0. 1794 . 0. 6740 "0.386684020-07 3,5461 0.5452 0.2.74 0.9930 0.54Q363060-07 3.6172 0.6152 0.2425 0.3183 0.712012U 0-07 0. 6371 C. 6340 0. 2349 0. 7148 _Q. 87*66)150-07 0. 6527 0,6483 ... . 0.2266.. 0. 5869 0,Yo3731020-06 3.6587 C.6529 0.2141 0.4465 0.129078040-06 3.6303 0.6221 0.154422050-06 3 , 7118 0. 7018 0. 17976607 0-06.__.0. 7 56 5....C. 7446 0.205110090-06 0.8040 0.7901 0. 1P22 Oo 2003 0.2010 0.2027 0. 98730- 01 0.1371 0. 37190-01 0.1027 0..363B554? 3-06- • .3.-9156 .C.336.4-0.5*? .97760-06 3, 9533 0.9101 _0,6_96340100-06 0. 9774 . 0.9193 0.860082440-06 0.9.38 0.9210 0.102382480-05 1.011 0.9233 n„1 1 3756710-05 1 , 025 0.92?> 0,135130940-05 0.151505180-05. 1 . 03 3 1. 054 0,74)549100-05 0.13 31 53300-04 0.192151700-04 0.25115 0090-04 0. 31014B4 80-04_ 0.369)4687 0-04 0,423145270-04 Oo 487143660-04 1 . 59 0 2, 058 2.49 4 2 ,9 37 3, 29S "3.66 3 4, 01 5 4. 340 0.9225 0, 9223 0.8973 0.8754 C, 85 70 Oo 33 74 0. 82 07 Oo 8033 0.7384 Oo 7672 0. 156 3 0.119 7 0„ 1026 0.94910-01 0.92360-0! 0.91420-01 3. 9221D-01 Oo 92250-01 0.1351 Oo 1712 0.2025 0. 2309 0.2558 0.27 8 2 0.2975 0. 31 37 Oo 92190-01. 0. 97050-01 0. 93 8 5 0-01 0. 95970-01 0. 92100-01 0. 95370-01 0. °290D-01 0.93850-01 0,9404D-01 Oo 9404D-01 Oo 93 8 2 0-01 Oo 93 790-01 0. 93 600-01 Oo "9356D-01" Oo 93400-01 Oo 9328D-01 _0, 63 0 32 785 0.-04_ 5 , 03 5__ 0o774506043-64 5,618 0.918184230-04 6.111 0,1 J?j_]jL6_?A0r_Li_ 3 A 2 _C. 726! 0. 6939 Co 6700 C. 6509 Oo 3366 3. 3497" Oo 3544 0.3534 Oc92990-01 o, 9,: 7:' ::-oi 0, 92580-01 Oo 92440-01 ~ 0 ~ 1 2 0 5 5 4 0 6 0 - 0 * 6 , 92 3 . C. 63 53 Oo 3493 0,92320-01 0.1 * 4 9 2 1 8 3 0-03 7,264 0.6232 0 , 3 4 3 1 0.92220-01 0 1 4 O 2 3 - 3 7 0 0-03 7, 567 3 , 613! 0,3358 0. 9213D-01 O.V63657520-33 ^ "7o 351 0. 60 3-, 0. 3263 0.9205C-01 O o ) R 5 0 0 6 O 1 3 - 0 3 8,133 0.5946 0 , . : 4 8 0.206234500-0* 3 . 5 0 9 0.5851 0,3079 0, 2 ? 7 4 6 2 9 9 0 - 33 3 . 7 3 4 _ 0 . 5 7 7 3 0. 2 9 ? 9 _ "0 . T 4 6 6 91490-33 9 , 0 4 0 0 . 5702 0. 2842 0,9)970-01 0, 9} 890-01 _0. 9! 820-01 Co 91 76 0-01 101 References 1. " E l e c t r o n i c s guide to CAD programs". E l e c t r o n i c s , A p r i l 13, 1970,. ppo 109-112. 2o Temes, G. C , "Optimization Methods i n C i r c u i t Design". In Computer-Oriented C i r c u i t Design, P 0 P. Kuo and ¥. G. Magnuson, eds., Englewood C l i f f s , N. J . : P r e n t i c e - H a l l , 1969, pp. 191-249. 3c D i r e c t o r , S.W.-and Rohrer, R 0A«, "The Generalized Adjoint Network and Network S e n s i t i v i t i e s " , IEEE Trans. C i r c u i t Theory;, Vol. v o l . CT-16, pp. 318-323, Aug. 1969. 4. J e s s e l , G.P., "Network S t a t i s t i c s f o r Computer-Aided Network An a l y s i s " , Proc 0 1973 I n t e r n a t ' l Symp. on C i r c u i t Theory, Toronto, Ontario, pp„ 283-286. 5. Tinney, W.F. and Ogbuobiri, E.C„, "Spa r s i t y Techniques: Theory and P r a c t i c e " , presented at IEEE S i x t h Region Conference, S e a t t l e , Washo, June 1970. 6o Hachtel, G.D., Brayton, R.K. and Gustavson,. F.G., "The Sparse Tableau Approach to Network A n a l y s i s and Design". IEEE Trans. C i r c u i t Theory, vole CT-18, pp„ 101-113, Jan. 1971. 7« Weeks, ¥.T., Jiminez, A.Jo, Mahoney, G.W., Qassemzadeh, H., and Scott, T.R., "Network Analysis Using a Sparse Tableau w i t h Tree Se l e c t i o n to Increase Sparseness", Proc. 1973 I n t e r n a t ' l Symp» on C i r c u i t Theory, Toronto, Ontario, pp. 165-168„ 8o Gustavson, F0G<>, L i n i g e r , Wo, and Willoughby, R.A., "Symbolic Generation of an Optimal Crout Algorithm f o r Sparse Systems of Linear Equations", J . Ass» Comput. Macho, volo 17, pp. 87-109, Jan. 1970o 9. Balabanian, N o , B i c k a r t , T . A o , and Seshu, S., E l e c t r i c a l Network Theory, New York: John Wiley and Sons, Inc., 1969, p» 229. 10o Lapidus, L., D i g i t a l Computation f o r Chemical Engineers, New Tork: M CGraw-Hill, 1962, p c 82o 11. Gear, C0 W„, "The Automatic Integration of Ordinary D i f f e r e n t i a l Equations", C. Ass 0 Comput. Macho, v o l . 14, pp. 176-180, Mar. 1970. 12. Brayton, R.K,, Gustavson, F.G., and Hachtel, G.D., "A New E f f i c i e n t Algorithm f o r Solving D i f f e r e n t i a l Algebraic Systems Using I m p l i c i t Backward D i f f e r e n t i a t i o n Formulas", Proc. IEEE, v o l , 60, pp. 98-108, Jan. 1972. 102 13. Gear, C.W., "Simultaneous Numerical Solution of D i f f e r e n t i a l -Algebraic Equations", IEEE Trans. C i r c u i t Theory, v o l . CT-18, pp. 89-95, Jan, 1971. 14. Gear, C. W,, "The Automatic Integration of Large Systems of Ordinary D i f f e r e n t i a l Equations", The Digest Record of the 1969 J o i n t Conference of Mathematical and Computer Aids to Design, Anaheim, C a l i f . , Oct. 1969, pp. 27-58. 15. Dahlquist, G,, "A Special S t a b i l i t y Problem f o r Linear Multistep Methods", BIT, v o l . 3, p. 27, 1963. 16. Weeks, W.T., Jiminez, A.J., Mahoney, G.W., Qassemzadeh, H., and Scott, T.R., op. c i t , , p. 168. 17. Branin, F.H., Hogsett, G,R., Lunde, R.L., and Kugel, L.E., "Ecap I I -An E l e c t r o n i c C i r c u i t A n alysis Program", IEEE Spectrum, v o l . 8 , pp. 18-25, June 1971. 18. Morris, A.G., "Submatrix Reduction and i t s Use i n the Simulation of Systems Having 1000 Nodes", Computer-Aided Design, pp. 19-22, Autumn, 1970. 19. S i l v e r b e r g , M., "Numerical So l u t i o n of D i s t r i b u t e d Networks Containing Nonlinear Elements", Ph. D. t h e s i s , Columbia Univ., 1967. 20. Nakhla, M., Singhal, R., and Vlach, J . , "Simple Method f o r Numerical Inversion of the Laplace Transform", Proc. Sixteenth Midwest Symposium on C i r c u i t Theory, Waterloo,. Ontario, A p r i l 1973. 21. Liou, M.L., "A Novel Method of Evaluating Transient Responses", Proc. IEEE, v o l . 54, pp. 20-23, Jan. 1966. 22. Thompson, W.E., "Evaluation of Transient Response", Proc. IEEE, v o l . 54, p. 1584, Nov. 1966. 23. Cbrrington, M., ''Simplified Calculation of Transient Response". Proc. IEEE, v o l . 53 , pp. 287-292, Mar, 1965, 24. Dawson, D. P., "The Topological Approach to Computer-Aided A n a l y s i s " . In Computer-Oriented C i r c u i t Design, P.P. Kuo and W.G. Magnuson, eds., Englewood C l i f f s , N.J.: P r e n t i c e - H a l l , 1969, p. 59. 25. N e i l l , T.B.M., "Techniques f o r C i r c u i t A n a l y s i s ( p a r t two)", Computer-Aided Design, pp, 35-43, Winter 1970. 26. Hurwitz, H., and Zweifel P.F., "Numerical Quadrature of Fourier Transform I n t e g r a l s " , MTAC, v o l . 10, pp. 140-149, 1956. 27. Saenger, A., "On Numerical Quadrature of Fourier Transforms", J . Math. Anal. Appl., v o l . 8 , pp. 1-3, 1964* 103 28. S i l v e r b e r g , Mo, op. c i t . , . p p . 22-24. 29. N e i l l , T.B.M., op. c i t . , p. 4 0 . 30. Gold, B., and Rader, C M . , D i g i t a l Processing of Signals. New York: McGraw-Hill, 1969, pp. 159-201. 31. S i l v e r b e r g , M., "A new Method of Solving State Variable Equations Permitting Large Step Siz e s " . Proc. IEEE, v o l . 56, pp. 1352-1353, Aug. 1968. 32. S i l v e r b e r g , M., " E f f i c i e n t Time-Domain Solutions Using Nodal State Variables"„ IEEE Trans. C i r c u i t Theory, v o l . CT-17* pp. 82-36, Feb. 1970. 33. Van Valkenburg, M.E., Network A n a l y s i s . Englewood C l i f f s , N.J.: P r e n t i c e - H a l l , 1965, p. 207. 34. L a t h i , B.P., Communications Systems. New York: John Wiley, 1968, p. 89. 35. Handbook of Mathematical Functions. M. Abromowitz and J . Stegun, eds., U. S. Govt. P r i n t i n g O f f i c e , Washington, D.C., 1964, p. 244. 36. N e i l l , T.B.M., "Techniques f o r C i r c u i t A n alysis (part t h r e e ) " , Computer-Aided Design. Summer 1970, p. 16. 37. Gold, B., and Rader, C M . , op» c i t . , pp. 203-232. 38. Stockham, T.G., "High Speed Convolution and C o r r e l a t i o n " , 1966 Spring J o i n t Computer Conference, AFIPS P r o c , v o l . 28, p. 230, 39. Davis, H.T., Introduction to Nonlinear D i f f e r e n t i a l and I n t e g r a l Equations. New York: Dover P u b l i c a t i o n s , 1962,, p. 415» 40. Maddox, I . J . , Elements of Functional A n a l y s i s . Cambridge, England: Cambridge U n i v e r s i t y Press, 1970, p. 56. 41. Heywood, D.R. and Moore, A.D., "Comment on an Improved Method of Analyzing Nonlinear E l e c t r i c a l Networks", Electr o n . L e t t . , v o l . 5, p. 269, 1969. 42. N e i l l , T.B.M., op. c i t . , p. 20. 43„ Broyden, C.G., "A Class of Methods f o r Solving Nonlinear Simultaneous Equations", Math. Comp., v o l . 19, Oct. 1965, p. 577* 44. Branin, F. H., "A Controlled-Step Newton Method f o r Solving Nonlinear Equations", The Digest Record of the 1969 J o i n t Conference on Mathematical and Computer Aids to Design, Anaheim, C a l i f . , p. 378, Oct. 1969. 45. Ostrowski, A.M., Solution of Equations and Systems of Equations. New York: Academic Press, 1966. Hachtel, G.D. Brayton, R.K., and Gustavson, E.G., "The Sparse Tableau Approach to Network Analysis and Design", IEEE Trans. C i r c u i t Theory, v o l . CT-18, pp. 101-113, Jan. 1971. 

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            data-media="{[{embed.selectedMedia}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0105073/manifest

Comment

Related Items